高斯消元

高斯消元用来求解线性方程组

构造增广矩阵,然后对增广矩阵消元

每次选取这一列绝对值最大的值作为主元,可以避免精度误差,如果发现这一列都为(0),则方程无解

然后将主元系数化为(1),矩阵化为上三角矩阵后,便可以回代求解

支持判是否多解

(code:)

bool gauss()
{
    for(int i=1;i<=n;++i)
    {
        int ma=i;
        for(int j=i+1;j<=n;++j)
            if(fabs(a[j][i])>fabs(a[ma][i]))
                ma=j;
        swap(a[ma],a[i]);
        double d=a[i][i];
        if(fabs(d)<eps) return false;
        for(int j=i;j<=n+1;++j) a[i][j]/=d;
        for(int j=i+1;j<=n;++j)
        {
            d=a[j][i];
            for(int k=i;k<=n+1;++k)
                a[j][k]-=a[i][k]*d;
        }
    }
    for(int i=n;i;--i)
        for(int j=i-1;j;--j)
            a[j][n+1]-=a[j][i]*a[i][n+1];
    return true;
}

支持判是否多解和是否无解

(code:)

void gauss()
{
    for(int i=1;i<=n;++i)
    {
        int ma=i;
        for(int j=i+1;j<=n;++j)
            if(fabs(a[j][i])>fabs(a[ma][i]))
                ma=j;
        swap(a[ma],a[i]);
        double d=a[i][i];
        if(fabs(d)<eps) continue;
        for(int j=i;j<=n+1;++j) a[i][j]/=d;
        for(int j=1;j<=n;++j)
        {
            if(i==j) continue;
            d=a[j][i];
            for(int k=i;k<=n+1;++k)
                a[j][k]-=a[i][k]*d;
        }
    }
    for(int i=1;i<=n;++i)
    {
        int num=0;
        for(int j=1;j<=n+1;++j)
            if(fabs(a[i][j])<eps) num++;
        if(num==n+1) inf_flag=true;
        if(num==n&&fabs(a[i][n+1])>=eps) no_flag=true;
    }
    if(inf_flag||no_flag) return;
    for(int i=n;i;--i)
        for(int j=i-1;j;--j)
            a[j][n+1]-=a[j][i]*a[i][n+1];
}

高斯消元还可以用来求逆矩阵

一个矩阵(A)若可以通过初等行变换转化为初等行变换,那么这个矩阵存在逆矩阵

三种初等行变换:

交换两行

将一行的所有元素乘上(k(k eq 0))

将一行的所有元素加上另一行的(k)

将每个初等变换对应初等矩阵,得(P_x...P_2P_1A=PA=I),那么(P=A^{-1})

在要求逆的(n)阶矩阵右侧加入一个(n)阶单位矩阵,用高斯消元将左侧矩阵变换成单位矩阵后,右侧矩阵即为所求的逆矩阵

(code:)

bool gauss()
{
    for(int i=1;i<=n;++i)
    {
        for(int j=i;j<=n;++j)
        {   
            if(a[j][i])
            {
                swap(a[j],a[i]);
                break;
            }
        }
        if(!a[i][i]) return false;
        ll d=inv(a[i][i]);
        for(int j=i;j<=m;++j) a[i][j]=a[i][j]*d%mod;
        for(int j=1;j<=n;++j)
        {
            if(i==j) continue;
            d=a[j][i]%mod;
            for(int k=i;k<=m;++k)
                a[j][k]=((a[j][k]-a[i][k]*d%mod)%mod+mod)%mod;
        }
    }
    return true;
}
原文地址:https://www.cnblogs.com/lhm-/p/12229653.html