【线性代数】高斯消元 - 行列式求值/生成树计数

普通代数式的消元法

考虑一个普通的二元一次方程组

\[\left \{ \begin{aligned} x_1+x_2&=6 \notag \\ x_1-x_2&=2 \notag \end{aligned} \right. \]

根据正常的数学知识，我们直到应该先消去一个未知数，解出另一个未知数，再去用另一个未知数表示消去的未知数

故对于上述\((1)\)式加上\((2)\)式，得\(2x_1=8\)，即\(x_1=4\)，代回原式得\(x_2=6-x_1=x_1-2=2\)

消元法遵循下列三个规则（懂的都懂）：

两方程互换，解不变
某一方程乘上非零数\(k\)，解不变
用某一方程加上（或减去）另外一方程，解不变

所以普通的消元法主要用于二元一次方程组的求解上，但对于多元方程组而言，普通消元法则表现出明显的力不从心~~（以前的数学题怎么都那么难）~~

线性代数的高斯消元法

德国数学家高斯分析了普通消元法并得出以下四个结论：

在消元法中，参与计算和发生改变的是方程中各变量的系数
各变量并未参与计算，且没有发生改变
可以利用系数的位置表示变量，从而省略变量
在计算中将变量简化省略，方程的解不变

根据这些结论，便能通过对增广矩阵的运用来求出多元方程组的通解，这便是高斯消元法

高斯消元法大致可分为以下五个步骤：

增广矩阵行初等行变换为行最简形
还原线性方程组
求解第一个变量
补充自由未知量
列表示方程组通解

举例说明，考虑下列多元一次方程组

\[\left \{ \begin{aligned} x_1+3x_2+x_4&=17 \notag \\ x_2+2x_4&=10 \notag \\ 2x_2+4x_4&=20 \notag \end{aligned} \right. \]

1、增广矩阵行初等行变换为行最简形

先将方程组以矩阵的方式表示

\[\left ( \begin{matrix} 1 & 3 & 0 & 1\\ 0 & 1 & 0 & 2\\ 0 & 2 & 0 & 4\\ \end{matrix} \middle | \begin{matrix} 17 \\ 10 \\ 20 \end{matrix} \right ) \]

根据增广矩阵的变换规则，将其变成行最简形

\[\xrightarrow{r_3-2r_2} \left ( \begin{matrix} 1 & 3 & 0 & 1\\ 0 & 1 & 0 & 2\\ 0 & 0 & 0 & 0\\ \end{matrix} \middle | \begin{matrix} 17 \\ 10 \\ 0 \end{matrix} \right ) \\ \xrightarrow{r_1-3r_2} \left ( \begin{matrix} 1 & 0 & 0 & -5\\ 0 & 1 & 0 & 2\\ 0 & 0 & 0 & 0\\ \end{matrix} \middle | \begin{matrix} -13 \\ 10 \\ 0 \end{matrix} \right ) \]

2、还原线性方程组

即将增广矩阵形式重新转回方程组形式

\[\left \{ \begin{aligned} x_1-5x_4&=-13 \notag \\ x_2+2x_4&=10 \notag \end{aligned} \right . \]

3、求解第一个变量

对于每个方程，将其第一个变量用其余变量表示出来

\[\left \{ \begin{aligned} x_1&=5x_4-13 \notag \\ x_2&=-2x_4+10 \notag \end{aligned} \right. \]

4、补充自由未知量

由于自由未知量不受约束，所以只能是\(x_i=x_i\)

所以将未出现的变量以上述形式表示出来

\[\left \{ \begin{aligned} x_1&=5x_4-13 \notag \\ x_2&=-2x_4+10 \notag \\ x_3&=x_3 \notag \\ x_4&=x_4 \notag \end{aligned} \right. \]

5、列表示方程组通解

\[\begin{aligned} \left ( \begin{matrix} x_1\\x_2\\x_3\\x_4 \end{matrix} \right ) =& \left ( \begin{matrix} 0\\0\\1\\0 \end{matrix} \right ) x_3 + \left ( \begin{matrix} 5\\-2\\0\\1 \end{matrix} \right ) x_4 + \left ( \begin{matrix} -13\\10\\0\\0 \end{matrix} \right ) \\ =& \left ( \begin{matrix} 0\\0\\1\\0 \end{matrix} \right ) C_1 + \left ( \begin{matrix} 5\\-2\\0\\1 \end{matrix} \right ) C_2 + \left ( \begin{matrix} -13\\10\\0\\0 \end{matrix} \right ) \end{aligned} \]

其中\(C_1,C_2\in \R\)

经过以上步骤，便得出了多元一次方程组的通解

高斯消元法的运用

行列式求值

对于行列式，存在公式

\[D=|A|=\sum_{p_1p_2\dots p_n} (-1)^r a_{(1,p_1)} a_{(2,p_2)} a_{(3,p_3)} \dots a_{(n,p_n)} \]

其中\(r\)表示\(\{p_1,p_2,p_3\dots p_n\}\)内逆序对对数

由于高斯消元法可以让我们得到行最简形矩阵，即一个对角线矩阵

此矩阵的行列式由对角线元素之积所决定

其正负符号可由交换行的数量来确定

故可以在\(O(n^3)\)时间内求出给定行列式的值

例题：NOJ 1035 - 行列式求值

#include<bits/stdc++.h>
using namespace std;

struct matrix
{
    typedef double Type;
    Type mat[1050][1050];
    int n,m;
    void init(int _n,int _m)
    {
        for(int i=1;i<=_n;i++)
            for(int j=1;j<=_m;j++)
                mat[i][j]=0;
        n=_n,m=_m;
    }
    Type getVal()
    {
        Type det=1;
        for(int i=1;i<=n;i++)
        {
            int k=i;
            for(int j=i+1;j<=n;j++)
                if(abs(mat[j][i])>abs(mat[k][i]))
                    k=j;
            if(mat[k][i]==0)
            {
                det=0;
                break;
            }
            swap(mat[i],mat[k]);
            if(i!=k)
                det=-det;
            det*=mat[i][i];
            for(int j=i+1;j<=n;j++)
                mat[i][j]/=mat[i][i];
            for(int j=1;j<=n;j++)
                if(j!=i&&mat[j][i]!=0)
                    for(int k=i+1;k<=n;k++)
                        mat[j][k]-=mat[i][k]*mat[j][i];
        }
        return det;
    }
}mat;

int main()
{
    int n;
    while(scanf("%d",&n)!=EOF)
    {
        mat.n=mat.m=n;
        for(int i=1;i<=n;i++)
            for(int j=1;j<=n;j++)
                scanf("%lf",&mat.mat[i][j]);
        printf("%.2f\n",mat.getVal());
    }
    return 0;
}

生成树计数

生成树的计数一般采用Matrix-Tree Theorem解决——

对于一张无向图，其生成树个数就等于这副图的拉普拉斯矩阵中任意一个代数余子式的值

建立拉普拉斯矩阵，通过高斯消元化成上三角型再进行求值

（行吧我也只会套模板了......）

例题：SPOJ HIGH-Highways

#include<bits/stdc++.h>
using namespace std;

struct matrix
{
    typedef long long Type;
    Type mat[15][15];
    int n,m;
    void init(int _n,int _m)
    {
        for(int i=1;i<=_n;i++)
            for(int j=1;j<=_m;j++)
                mat[i][j]=0;
        n=_n,m=_m;
    }
    void addedge(int a,int b)
    {
        mat[a][a]++;
        mat[b][b]++;
        mat[a][b]--;
        mat[b][a]--;
    }
    Type gauss()
    {
        Type res=1;
        for(int i=1;i<n;i++)
        {
            for(int j=i+1;j<n;j++)
            {
                while(mat[j][i])
                {
                    Type t=mat[i][i]/mat[j][i];
                    for(int k=i;k<n;k++)
                        mat[i][k]=mat[i][k]-t*mat[j][k];
                    swap(mat[i],mat[j]);
                    res=-res;
                }
            }
            res=res*mat[i][i];
        }
        return res;
    }
}mat;

void solve()
{
    int n,m,u,v;
    scanf("%d%d",&n,&m);
    mat.init(n,n);
    while(m--)
    {
        scanf("%d%d",&u,&v);
        mat.addedge(u,v);
    }
    printf("%lld\n",mat.gauss());
}
int main()
{
    int T;
    scanf("%d",&T);
    while(T--)
        solve();
    return 0;
}

上面用到的板子

struct matrix
{
    typedef long long Type;
    Type mat[1050][1050];
    int n,m;
    void init(int _n,int _m)
    {
        for(int i=1;i<=_n;i++)
            for(int j=1;j<=_m;j++)
                mat[i][j]=0;
        n=_n,m=_m;
    }
    Type getVal()
    {
        Type det=1;
        for(int i=1;i<=n;i++)
        {
            int k=i;
            for(int j=i+1;j<=n;j++)
                if(abs(mat[j][i])>abs(mat[k][i]))
                    k=j;
            if(mat[k][i]==0)
            {
                det=0;
                break;
            }
            swap(mat[i],mat[k]);
            if(i!=k)
                det=-det;
            det*=mat[i][i];
            for(int j=i+1;j<=n;j++)
                mat[i][j]/=mat[i][i];
            for(int j=1;j<=n;j++)
                if(j!=i&&mat[j][i]!=0)
                    for(int k=i+1;k<=n;k++)
                        mat[j][k]-=mat[i][k]*mat[j][i];
        }
        return det;
    }
    void addedge(int a,int b) //生成树计数用
    {
        mat[a][a]++;
        mat[b][b]++;
        mat[a][b]--;
        mat[b][a]--;
    }
    Type gauss() //求矩阵K的n-1阶顺序主子式
    {
        Type res=1;
        for(int i=1;i<n;i++)
        {
            for(int j=i+1;j<n;j++)
            {
                while(mat[j][i])
                {
                    Type t=mat[i][i]/mat[j][i];
                    for(int k=i;k<n;k++)
                        mat[i][k]=mat[i][k]-t*mat[j][k];
                        //mat[i][k]=(mat[i][k]-t*mat[j][k]+mod)%mod;
                    swap(mat[i],mat[j]);
                    res=-res;
                }
            }
            res=res*mat[i][i];
            //res=res*mat[i][i]%mod;
        }
        return res;
        //return (res+mod)%mod;
    }
}mat;

参考及引用：OI Wiki / Math / Gauss

posted @ 2020-08-14 17:05 StelaYuri 阅读(692) 评论(0) 编辑收藏举报

会员力量，点亮园子希望

刷新页面返回顶部

StelaYuri