【线性代数】

线代前置知识：3B1B的视频和课本
代码大量参考佳叔。

矩阵

矩阵快速幂板子。
UVA 10870 Recurrences
矩阵表示线性递推关系。

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int sz = 17;
int n,d;
int MOD;
struct mat{
    ll a[sz][sz];
    inline mat() {memset(a,0,sizeof(a));}
    inline mat operator-(const mat& T) const{
        mat res;
        for(int i=0;i<d;++i)
            for(int j=0;j<d;++j)
                res.a[i][j]=(a[i][j]-T.a[i][j]+MOD)%MOD;
        return res;
    }
    inline mat operator+(const mat& T) const{
        mat res;
        for(int i=0;i<d;++i)
            for(int j=0;j<d;++j)
                res.a[i][j]=(a[i][j]+T.a[i][j])%MOD;
        return res;
    }
    inline mat operator*(const mat& T) const{
        mat res;
        int r;
        for(int i=0;i<d;++i)
            for(int k=0;k<d;++k){
                r=a[i][k];
                for(int j=0;j<d;++j)
                    res.a[i][j]+=T.a[k][j]*r,res.a[i][j]%=MOD;
            }
        return res;
    }
    inline mat operator^(ll x) const{
        mat res,bas;
        for(int i=0;i<d;++i) res.a[i][i]=1;
        for(int i=0;i<d;++i)
            for(int j=0;j<d;++j) bas.a[i][j]=a[i][j]%MOD;
        while(x){
            if(x&1) res=res*bas;
            bas=bas*bas;
            x>>=1;
        }
        return res;
    }
};

int main(){
    while(~scanf("%d%d%d",&d,&n,&MOD),(d||n||MOD)){
        mat a;
        for(int i=0;i<d;++i){
            int tmp;scanf("%d",&tmp);
            a.a[0][i]=tmp;
        }
        for(int i=1;i<d;++i) a.a[i][i-1]=1;
        mat x;
        for(int i=d-1;i>=0;--i) {
            int tmp;scanf("%d",&tmp);
            x.a[i][0]=tmp;
        }
        if(n<=d){
            cout<<x.a[n-1][0]<<endl;
            continue;
        }
        a=a^(n-1);
        x=a*x;
        cout<<x.a[d-1][0]<<endl;
    }
}

UVA 1386 Cellular Automaton
细胞自动机。
数据范围不允许O(n^3log n)。发现转移矩阵是个circulant matrix。
两个循环矩阵乘积还是循环矩阵。
只保存一行。O(n^2log n)。最后$Ax$的时候也是用保存的那一行，推出全部，跟x以O(n^2)相乘。

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll n,m,d,k;
const int N = 504;

struct mat{
    ll v[N];
    mat() {memset(v,0,sizeof(v));}
    inline mat operator *(mat c) const{
        mat ans;
        for(int i=0;i<n;++i)
            for(int j=0;j<n;++j)
                ans.v[i]=(ans.v[i]+v[j]*c.v[(i-j+n)%n])%m;
        return ans;
    }
};

mat qpow(mat x,ll k){
    if(k==1) return x;
    mat tmp=x*x;
    mat res=qpow(tmp,k>>1);
    if(k&1) res=res*x;
    return res;
}

int b[N];

int main(){
    while(~scanf("%lld%lld%lld%lld",&n,&m,&d,&k)){
        mat a;
        for(int i=0;i<=d;++i) a.v[i]=1;
        for(int i=n-1;i>=n-d;--i) a.v[i]=1;
        a=qpow(a,k);
        for(int i=0;i<n;++i) scanf("%lld",&b[i]);
        for(int i=0;i<n;++i){
            ll ans=0;
            for(int j=0;j<n;++j)
                ans=(ans+b[j]*a.v[(j-i+n)%n])%m;
            printf("%lld%c",ans,(i==n-1?'\n':' '));
        }
    }
}

UVA 10828 Back to Kernighan-Ritchie
随机程序。
给一张概率自动机图，1号点是起点。要求每个结点的期望执行次数。
lrj思路：每个节点的答案由其pre节点得出，可列多条方程。
解多元一次方程组，用高斯约当消元法。
注意观察什么时候得出infinite，什么时候0.
结合高斯约当的特质来看。
RREF的特质就是：可逆情况下，adjacency matrix的对角线才有数字。当出现矛盾方程，会有前面为0，后面不为0的情况。当出现多余方程，会出现全为0的情况。
下图这种情况（*）也可能出现。此时，我们发现第3条方程是矛盾方程，那么$x_3$就是infinity（在此题目语意下）。那么，与他有关的第一条方程的解$x_1$也是infinity。

#include<bits/stdc++.h>
using namespace std;
int n;
const int N = 110;
const double eps = 1e-8;
int out[N];  // 用来算概率
vector<int> pre[N];  // 前驱

typedef double Matrix[N][N];

Matrix matrix;

void gauss_jordan(Matrix A,int n){
    int i,j,k,r;
    for(i=0;i<n;++i){
        r=i;
        for(j=i+1;j<n;++j)
            if(fabs(A[j][i])>fabs(A[r][i])) r=j;
        if(fabs(A[r][i])<eps) continue;//i是自由变量
        if(r!=i) for(j=0;j<=n;++j) swap(A[i][j],A[r][j]);
        // 变成对角矩阵
        for(k=0;k<n;++k) if(k!=i)
            for(j=n;j>=i;--j)
                A[k][j]-=A[k][i]/A[i][i]*A[i][j];
    }
}

int inf[N];  // check if is inf
int main(){
    int cas=0;
    while(~scanf("%d",&n),n){
        int u,v;
        for(int i=0;i<n;++i) pre[i].clear();
        memset(out,0,sizeof(out));
        while(~scanf("%d%d",&u,&v),(u||v)){
            u--,v--;
            pre[v].push_back(u);
            out[u]++;
        }
        memset(matrix,0,sizeof(matrix));
        for(int i=0;i<n;++i){
            matrix[i][i]=1;
            for(int j=0;j<pre[i].size();++j)
                matrix[i][pre[i][j]]-=1.0/out[pre[i][j]];
            if(i==0) matrix[i][n]=1;
        }
        gauss_jordan(matrix,n);
        memset(inf,0,sizeof(inf));
        for(int i=n-1;i>=0;--i){
            if(fabs(matrix[i][i])<eps&&fabs(matrix[i][n])>eps) inf[i]=1;
            for(int j=i+1;j<n;++j)
                if(fabs(matrix[i][j])>eps&&inf[j]) inf[i]=1;
                // 处理(*)的情况
        }
        int q;scanf("%d",&q);
        cout<<"Case #"<<++cas<<":"<<endl;
        while(q--){
            scanf("%d",&u);u--;
            if(inf[u]) cout<<"infinity"<<endl;
            else printf("%.3lf\n",fabs(matrix[u][u])<eps?0.0:matrix[u][n]/matrix[u][u]);
        }
    }
}

洛谷[SDOI2010]外星千足虫
异或方程组裸题，问是否有自由变量，有的话输出cannot determine，无就输出用了前多少个方程求出了所有解，并输出解。
玄学的时间复杂度，居然过了（。

#include<bits/stdc++.h>
using namespace std;
char ss[1002];
int A[2002][1002];

void result(int m,int n){
    int i=0,j=0,k,r,u;
    int las=0;
    while(i<m&&j<n){
        r = i;
        for(k=i;k<m;++k)
            if(A[k][j]) {r = k;break;}
        if(A[r][j]){
            las = max(las,r);
            if(r!=i) for(k=0;k<=n;++k) swap(A[r][k],A[i][k]);
            for(u=i+1;u<m;++u) if(A[u][j])
                for(k=i;k<=n;++k) A[u][k] ^= A[i][k];
            i++;
        }
        j++;
    }
    if(i<n){
        cout<<"Cannot Determine"<<endl;return;
    }// free variable
    cout<<las+1<<endl;
    for(i=n-1;i>=0;i--)
        for(j=i-1;j>=0;j--)
        if(A[j][i]) A[j][n] ^= A[i][n];
    for(i=0;i<n;++i) {
        if(A[i][n]) cout<<"?y7M#"<<endl;
        else cout<<"Earth"<<endl;
    }
}

int main() {
    int n,m;cin>>n>>m;
    int res;
    for(int i=0;i<m;++i){
        // m 条方程, n 个变量
        scanf("%s",ss);
        scanf("%d",&res);
        for(int j=0;j<n;++j) {
            A[i][j] = ss[j] - '0';
        }
        A[i][n] = res;
    }
    result(m,n);
}

线性基

一篇讲得很好的博文