【学习笔记】简单数论-高斯消元与线性空间

发布时间 2023-08-20 21:23:33作者: The_Shadow_Dragon

友情提示

线性方程组

  • 线性方程组是由 \(M\)\(N\) 元一次方程共同构成的。
  • 线性方程组的所有系数可以写成一个 \(M\)\(N\) 列的系数矩阵,再加上每个方程等号右侧的常数,可以写成 一个 \(M\)\(N+1\) 列的增广矩阵。
  • 对增广矩阵的三类操作(初等行变化):
    • 用一个非零的数乘某一行。
    • 把其中一行的若干倍加到另一行上。
    • 交换两行的位置。
  • 用若干次初等行变化求解方程组,最后得到的矩阵被称为阶梯形矩阵,它的系数矩阵部分被称为上三角矩阵,从下往上依次带回方程组,即可得到每个未知数的解;将其进一步化简最后得到的矩阵被称为简化阶梯形矩阵,它的系数矩阵部分是一个对角矩阵,该矩阵直接给出了方程组的解(这部分可能有些难懂,建议参考蓝书)。

高斯消元(约旦消元)

  • 通过初等行列变化把增广矩阵变为简化阶梯形矩阵的线性方程组求解算法就是高斯消元算法。

  • 思想:对于每个未知数 \(x_i\) ,找到一个 \(x_i\) 的系数非零但 \(x_1 \sim x_{i-1}\) 的系数都为零的方程,然后用初等行变化把其他方程的 \(x_i\) 的系数全部消成零。

  • 特殊情形:

    • 在高斯消元过程中,可能出现 \(0=d\) 这样的方程,其中 \(d\) 是一个非零常数,表明原方程组中某些方程之间存在矛盾,有原方程组无实数解。
    • 在高斯消元过程中,有可能找不到 \(x_i\) 的系数非零但 \(x_1 \sim x_{i-1}\) 的系数都为零的方程,此时有原方程组有无穷多实数解。将存在 \(x_i\) 的系数非零但 \(x_1 \sim x_{i-1}\) 的系数都为零的方程中的 \(x_i\) 称为主元,否则称为自由元。
      • 对于每个主元,整个简化阶梯形矩阵中有且仅有一个位置 \((i,j)\) 满足该主元的系数非零,第 \(j\) 列的其他位置都为零,第 \(i\) 行的第 \(1 \sim j-1\) 列都为零。
  • 总结:

    • 高斯消元完成后,若存在系数全为零、常数不为零的行,则原方程组无实数解;若系数不全为零的行恰好有 \(n\) 个,则说明有 \(n\) 个主元,原方程组有唯一解;若系数不全为零的行有 \(k(k<n)\) 个,则说明有 \(k\) 个主元, \(n-k\) 个自由元,原方程组有无穷多实数解。
  • 例题:

    const double eps=1e-12;//考虑浮点数误差
    double a[201][201];
    int main()
    {
        int n,i,j,k,val,flag=0;
        cin>>n;
        for(i=1;i<=n;i++)
        {
            for(j=1;j<=n+1;j++)
            {
                cin>>a[i][j];
            }
        }
        for(i=1;i<=n;i++)//枚举列
        {
            val=i;
            for(j=i+1;j<=n;j++)
            {
                if(fabs(a[j][i])-fabs(a[val][i])>=eps)//选出该列最大系数
                {
                    val=j;
                }
            }
            for(j=1;j<=n+1;j++)//将第i行与第val行交换
            {
                swap(a[i][j],a[val][j]);
            }
            if(a[i][i]==0)//若最大值为0,说明该列都为0,即无实数解
            {
                flag=1;
                cout<<"No Solution"<<endl;
                break;
            }
            else
            {
                for(j=1;j<=n;j++)
                {
                    if(j!=i)//特判主元这一项
                    {
                        for(k=i+1;k<=n+1;k++)
                        {
                            a[j][k]-=a[i][k]*a[j][i]/a[i][i];//消去其他方程的x[i]的系数
                        }
                    }
                }
            }
        }
        if(flag==0)
        {
            for(i=1;i<=n;i++)//消去系数
            {
                printf("%.2lf\n",a[i][n+1]/a[i][i]);
            }
        }
        return 0;
    }
    
    const double eps=1e-12;
    double a[101][101];
    int main()
    {
        int n,i,j,k,val,flag=1;
        cin>>n;
        for(i=1;i<=n;i++)
        {
            for(j=1;j<=n+1;j++)
            {
                cin>>a[i][j];
            }
        }
        for(i=1;i<=n;i++)
        {
            val=i;
            for(j=1;j<=n;j++)//上题中是从i开始枚举的,有用式子的仅有i之后的,系数为0直接break,因为只需要判断是否有解
            //但是本题需要判断是否有解,有无数解,所以需要在系数为0后继续算下去,所以有用的式子不仅是i之后的式子,还有i之前且系数为0的式子
            {
                if(!(i>j&&fabs(a[j][j])!=0))
                {
                    if(fabs(a[j][i])-fabs(a[val][i])>=eps)
                    {
                        val=j;
                    }
                }
            }
            for(j=1;j<=n+1;j++)
            {
                swap(a[i][j],a[val][j]);
            }
            if(a[i][i]!=0)
            {
                for(j=1;j<=n;j++)
                {
                    if(j!=i)
                    {
                        for(k=i+1;k<=n+1;k++)
                        {
                            a[j][k]-=a[i][k]*a[j][i]/a[i][i];
                        }
                    }
                }
            }
        }
        for(i=1;i<=n;i++)	
        {
            if(a[i][i]==0)
            {
                if(a[i][n+1]!=0)//无解比无数解优先级更高
                {
                    flag=-1;
                    break;
                }
                else
                {
                    flag=0;
                }
            }
        }
        if(flag==1)
        {
            for(i=1;i<=n;i++)
            {
                printf("x%d=%.2lf\n",i,a[i][n+1]/a[i][i]);
            }
        }
        else
        {
            cout<<flag<<endl;
        }
        return 0;
    }