zoukankan      html  css  js  c++  java
  • 线性代数学习笔记

    1 高斯消元法

    这里采用更加简洁的高斯-约旦消元法,直接消成对角矩阵。

    考虑到还要处理无解和无穷解,我们记录一个 nline,表示前 nline 个方程时有用的(交换过之后)。然后我们只需要对剩下的方程判断一下即可。

    SDOI2006 线性方程组

    #include<bits/stdc++.h>
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=59;
    const double eps=1e-7;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n;
    double a[N][N],x[N];
    
    int gauss() {
    	int nl=1;
    	rep(i,1,n) {
    		int t=nl;
    		rep(j,nl+1,n) if(fabs(a[j][i])>fabs(a[t][i])) t=j;
    		if(t!=nl) swap(a[nl],a[t]);
    		if(fabs(a[nl][i])<eps) continue;
    		rep(j,1,n) if(nl!=j) {
    			double b=a[j][i]/a[nl][i];
    			rep(k,i+1,n+1) a[j][k]-=a[nl][k]*b;
    		}
    		nl++;
    	}
    	if(nl<n+1) {
    		while(nl<n) if(a[nl++][n+1]) return -1;
    		return 0; 
    	}
    	rep(i,1,n) x[i]=fabs(a[i][n+1]/a[i][i])<eps?0:a[i][n+1]/a[i][i];
    	return 1;
    }
    
    int main() {
    	n=read();
    	rep(i,1,n) rep(j,1,n+1) a[i][j]=read();
    	int ans=gauss();
    	if(ans!=1) printf("%d
    ",ans);
    	else {
    		rep(i,1,n) printf("x%d=%.2lf
    ",i,x[i]);
    	}
    	return 0;
    }
    

    1.2 高斯消元解异或方程

    本质还是一样的,不过有一点,就是因为值只有 0 和 1,所以可以通过 bitset 去优化高斯消元。

    题目:P2447 [SDOI2010]外星千足虫

    这道题还要额外求一个最早到第几个方程可以判断解。我们用贪心的方法,每一个元都找到离它的最近且系数为 1 的方程即可。

    #include<bits/stdc++.h>
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=2009;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,m,ans;
    bitset<N>a[N];
    char c[N];
    
    bool gauss() {
    	rep(i,1,n) {
    		int t=i;
    		while(t<=m&&!a[t][i]) t++;
    		if(t==m+1) return 0;
    		ans=max(ans,t);
    		if(i!=t) swap(a[t],a[i]);
    		rep(j,1,m) if(i!=j&&a[j][i]) a[j]^=a[i];
    	}
    	return 1;
    }
    
    int main() {
    	n=read(), m=read();
    	rep(i,1,m) {
    		scanf("%s",c+1); int v=read();
    		rep(j,1,n) a[i][j]=c[j]-'0'; a[i][n+1]=v;
    	}
    	bool res=gauss();
    	if(!res) puts("Cannot Determine");
    	else {
    		printf("%d
    ",ans);
    		rep(i,1,n) puts(a[i][n+1]?"?y7M#":"Earth");
    	}
    	return 0;
    }
    

    [CQOI2014]和谐矩阵

    题目中的要求可以转换为 (a_{x}{y-1} ext{ xor } a_{x-1}{y} ext{ xor } a_{x+1}{y} ext{ xor } a_{x}{y+1}=0)

    然后我们把所有自由元赋值为 1,于是我们就可以得到一个合法的解了。

    #include<bits/stdc++.h>
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=49,K=1609;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,m,x[K],ans[K];
    bitset<K>a[K];
    vector<int>fre;
    bool ok(int x,int y) {return x>0&&x<=n&&y>0&&y<=m;}
    int hsh(int x,int y) {return m*(x-1)+y;}
    
    void gauss() {
    	rep(i,1,n*m) {
    		int t=i;
    		while(t<=n*m&&!a[t][i]) t++;
    		if(t==n*m+1) ans[i]=1;
    		if(t!=i) swap(a[t],a[i]);
    		rep(j,i+1,n*m) if(a[j][i]) a[j]^=a[i];
    	}
    	per(i,n*m,1) rep(j,i+1,n*m) if(a[i][j]) ans[i]^=ans[j];
    }
    
    int main() {
    	n=read(), m=read();
    	rep(i,1,n) rep(j,1,m) {
    		int h=hsh(i,j);
    		if(ok(i+1,j)) a[h][hsh(i+1,j)]=1;
    		if(ok(i,j+1)) a[h][hsh(i,j+1)]=1;
    		if(ok(i-1,j)) a[h][hsh(i-1,j)]=1;
    		if(ok(i,j-1)) a[h][hsh(i,j-1)]=1;
    		a[h][h]=1;
    	}
    	gauss();
    	rep(i,1,n) {
    		rep(j,1,m) printf("%d ",ans[hsh(i,j)]);
    		puts("");
    	}
    	return 0;
    }
    

    1.3 高斯消元解决DP问题

    P2973 [USACO10HOL]Driving Out the Piggies G

    期望 DP+高斯消元

    直接求概率行不通。所以我们设 (f_u) 代表未爆炸经过节点 (u) 的期望次数。

    (d_u)(u) 的度数,可以得到方程

    [f(u)=(1-frac{p}{q})sum_{v o u}frac{1}{d_v}f(v) ]

    (u) 的答案就是 (f_u imes frac{p}{q})

    然后这个式子有后效性。不过不要紧,我们有高斯消元。

    式子转化成

    [f(u)+sum_{v o u}(frac{p}{q}-1) imes d_v imes f(v)=0 ]

    由于初始状态是 (f(1)=1),所以第一个方程的等于号后面的常数项为 (1)

    #include<bits/stdc++.h>
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=309;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,m,p,q,d[N];
    bool e[N][N];
    double k,a[N][N],f[N];
    
    void gauss() {
    	rep(i,1,n) {
    		int t=i;
    		rep(j,i+1,n) if(fabs(a[j][i])>fabs(a[t][i])) t=j;
    		if(t!=i) swap(a[i],a[t]);
    		rep(j,1,n) if(i!=j) {
    			double b=a[j][i]/a[i][i];
    			rep(k,i+1,n+1) a[j][k]-=a[i][k]*b;
    		}
    	}
    	rep(i,1,n) f[i]=a[i][n+1]/a[i][i];
    }
    
    int main() {
    	n=read(), m=read(), p=read(), q=read();
    	k=1.*p/q;
    	rep(i,1,m) {
    		int u=read(), v=read();
    		e[u][v]=e[v][u]=1, d[u]++, d[v]++;
    	}
    	rep(i,1,n) {
    		a[i][i]=1, a[i][n+1]=(i==1);
    		rep(j,1,n) if(i!=j&&e[i][j]) a[i][j]=(k-1)/d[j];
    	}
    	gauss();
    	rep(i,1,n) printf("%.9lf
    ",f[i]*k);
    	return 0;
    }
    

    P3211 [HNOI2011]XOR和路径

    由于异或的期望不是期望的异或,所以我们不能直接处理。考虑将位拆开,分位求期望,这样就可以了。对于所有 (X),设 (f(u)) 表示从 (u) 走到 (n) 然后得到的答案的第 (X) 位为 1 的概率。

    [f_u=frac{1}{d_u}(sum_{w_{u,v}=0}f_v+sum_{w_{u,v}=1}1-f_v) ]

    (k_i) 为对于第 (X)(sum w_{u,v})

    可得

    [d_uf_u-sum_{w_{u,v}=0}f_v+sum_{w_{u,v}=1}f_v=k_u ]

    #include<bits/stdc++.h>
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=109,M=10009;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,m,d[N],k[N];
    double ans,a[N][N],f[N];
    vector<int>e[N],p[N];
    
    void add(int u,int v,int w){e[u].push_back(v),p[u].push_back(w);}
    bool bit(int s,int x) {return s&(1<<x);} 
    
    void gauss() {
    	rep(i,1,n) {
    		int t=i;
    		rep(j,i+1,n) if(fabs(a[j][i])>fabs(a[t][i])) t=j;
    		if(t!=i) swap(a[i],a[t]);
    		rep(j,1,n) if(i!=j) {
    			double b=a[j][i]/a[i][i];
    			rep(k,i+1,n+1) a[j][k]-=a[i][k]*b;
    		}
    	}
    	rep(i,1,n) f[i]=a[i][n+1]/a[i][i];
    }
    
    signed main() {
    	n=read(), m=read();
    	rep(i,1,m) {
    		int u=read(), v=read(), w=read();
    		add(u,v,w),d[u]++;
    		if(u!=v) add(v,u,w),d[v]++;
    	}
    	rep(h,0,30) {
    		memset(a,0,sizeof(a));
    		a[n][n]=1;
    		rep(u,1,n-1) {
    			a[u][u]=d[u];
    			for(int i=0;i<e[u].size();i++) {
    				int v=e[u][i],w=bit(p[u][i],h);
    				if(w) a[u][v]++,a[u][n+1]++;
    				else a[u][v]--;
    			}
    		}
    		gauss();
    		ans+=f[1]*(1<<h);
    	} 
    	printf("%.3lf
    ",ans);
    	return 0;
    } 
    

    2 行列式

    矩阵行列式的定义为

    [D=sum (-1)^v prod a_{i,l_i} ]

    其中 (l) 为排列,(v)(l) 的逆序对个数。

    行列式有一些性质

    • 行交换,行列式取反
    • 行减去另一行的倍数,行列式不变
    • 矩阵同时乘 (k),行列式乘 (k)

    消元之后可以得到一个三角矩阵,行列式的值即对角线之积。注意消元的时候,要用第 (j) 行减去第 (i) 行的倍数。特殊情况,当矩阵非满秩时,行列式为 0。然而很多东西是要取模的,而且 (a) 甚至还不一定在模意义下有逆元。

    由于可以相减,考虑类似于辗转相除法。对于当前行和要消的行,我们不断相减,直到被消掉。至于复杂度,还是均摊 (O(n^3)) 的。

    洛谷 行列式求值 模板

    #include<bits/stdc++.h>
    #define int long long
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=609;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,mod,a[N][N];
    
    int gauss(int det=1) {
    	rep(i,1,n) {
    		rep(j,i+1,n) {
    			if(a[j][i]>a[i][i]) det=-det, swap(a[j],a[i]);
    			while(a[j][i]) {
    				int b=a[i][i]/a[j][i];
    				rep(k,i,n) a[i][k]=(a[i][k]+(mod-b)*a[j][k])%mod;
    				det=-det, swap(a[j],a[i]);
    			}
    		}
    		det=(det*a[i][i]%mod+mod)%mod;
    	}
    	return det;
    }
    
    signed main() {
    	n=read(), mod=read();
    	rep(i,1,n) rep(j,1,n) a[i][j]=read(), a[i][j]%=mod;
    	int det=gauss();
    	printf("%lld
    ",det);
    	return 0;
    }
    

    3 矩阵树定理

    设无向无权图的基尔霍夫矩阵 (K) 为度数矩阵 (D) 减去邻接矩阵 (A)。则其生成树个数为 ( ext{det}(K')),其中 (K')(K) 去掉任意一行一列的 (n-1) 阶主子式。

    对于带权无向图,度数矩阵变成了相邻边的权值和。基尔霍夫矩阵的 (n-1) 阶主子式的行列式为所有生成树边权乘积的总和。

    对于有向图,有两种情况:根向树和叶向树。根向树定义 (D) 为出边边权总和,而叶向树定义 (D) 为入边边权总和。注意,设根为 (r),有向图情况去掉的行列必须是第 (r)(r) 列。

    洛谷 矩阵树定理 模板

    我们发现对于有向图,直接去掉一行一列很烦。由于他和标号关系不大,所以我们可以有一个小技巧,就是让他和第n个节点交换一下,这样去掉一行一列就方便多了。

    #include<bits/stdc++.h>
    #define int long long
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=309,mod=1e9+7;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int a[N][N];
    
    int gauss(int n,int det=1) {
    	rep(i,1,n) {
    		rep(j,i+1,n) {
    			if(a[j][i]>a[i][i]) det=-det, swap(a[j],a[i]);
    			while(a[j][i]) {
    				int b=a[i][i]/a[j][i];
    				rep(k,i,n) a[i][k]=(a[i][k]+(mod-b)*a[j][k])%mod;
    				det=-det, swap(a[j],a[i]);
    			}
    		}
    		det=(det*a[i][i]%mod+mod)%mod;
    	}
    	return det;
    }
    
    signed main() {
    	int n=read(), m=read(), t=read();
    	rep(i,1,m) {
    		int u=read(), v=read(), w=read();
    		if(t==0) {
    			(a[u][u]+=w)%=mod, (a[v][v]+=w)%=mod, (a[u][v]-=w)%=mod, (a[v][u]-=w)%=mod;
    		} else {
    			if(u==1) u=n; else if(u==n) u=1;
    			if(v==1) v=n; else if(v==n) v=1;
    			(a[v][v]+=w)%=mod, (a[u][v]-=w)%=mod;
    		}
    	}
    	printf("%lld
    ",gauss(n-1));
    	return 0;
    }
    

    P4455 [CQOI2018]社交网络

    有向图的Matrix-Tree定理模板。

    P3317 [SDOI2014]重建

    矩阵树求出来的是 (sum_T prod_{ein T} p_e)。但是题目中还有一个限制,就是其他边不能选择,即

    [ans=sum_T prod_{ein T}p_e prod_{e otin T} (1-p_e) ]

    考虑把后面的 (e otin T) 化成 (ein T) 方便矩阵树去求。考虑用乘法原理的补集转换,可以得到

    [ans=sum_T prod_{ein T}p_e frac{prod (1-p_e)}{、prod_{ein T}(1-p_e)} ]

    整理一下式子,即

    [ans=(prod(1-p_e)) imes sum_T{prod_{ein T}frac{p_e}{1-p_e}} ]

    用矩阵树暴力算一遍即可。

    #include<bits/stdc++.h>
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=59;
    
    double ans=1,a[N][N];
    
    double gauss(int n,double ret=1) {
    	rep(i,1,n) {
    		int t=i;
    		rep(j,i+1,n) if(fabs(a[j][i])>fabs(a[t][i])) t=j;
    		if(t!=i) swap(a[i],a[t]), ret=-ret;
    		if(!a[i][i]) return 0;
    		rep(j,i+1,n) {
    			double b=a[j][i]/a[i][i];
    			rep(k,i,n) a[j][k]-=a[i][k]*b;
    		}
    		ret*=a[i][i];
    	}
    	return ret;
    }
    
    int main() {
    	int n; scanf("%d",&n);
    	rep(i,1,n) {
    		rep(j,1,n) {
    			double w; scanf("%lf",&w);
    			if(w==1) w-=0.000001;
    			if(i<j) ans*=(1-w);
    			w=w/(1-w);
    			a[i][i]+=w, a[i][j]-=w;
    		}
    	}
    	ans*=gauss(n-1);
    	printf("%.4lf
    ",ans);
    	return 0;
    }
    

    P4336 [SHOI2016]黑暗前的幻想乡

    经典题。题目中强调每个公司只能用一遍,但实际上我们的矩阵树没法直接做到这一点。那怎么办?容斥!容斥一遍然后就可以安心矩阵树了。

    复杂度 (O(2^{n-1}n^3))

    #include<bits/stdc++.h>
    #define int long long
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=19,mod=1e9+7;
    typedef pair<int,int> pii;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,c[N],a[N][N],ans;
    vector<pii>r[N];
    
    bool bit(int s,int i) {return s&(1<<i);}
    
    int gauss(int n,int det=1) {
    	rep(i,1,n) {
    		rep(j,i+1,n) {
    			if(a[j][i]>a[i][i]) det=-det, swap(a[j],a[i]);
    			while(a[j][i]) {
    				int b=a[i][i]/a[j][i];
    				rep(k,i,n) a[i][k]=(a[i][k]+(mod-b)*a[j][k])%mod;
    				det=-det, swap(a[j],a[i]);
    			}
    		}
    		det=(det*a[i][i]%mod+mod)%mod;
    	}
    	return det;
    }
    
    signed main() {
    	n=read();
    	rep(i,1,n-1) {
    		c[i]=read();
    		if(c[i]==0) {puts("0"); return 0;}
    		rep(j,1,c[i]) r[i].push_back(make_pair(read(),read()));
    	}
    	int S=(1<<n-1)-1;
    	rep(s,1,S) {
    		memset(a,0,sizeof(a));
    		int pop=0;
    		rep(i,1,n-1) if(bit(s,i-1)) {
    			pop++;
    			rep(j,0,c[i]-1) {
    				int u=r[i][j].first, v=r[i][j].second;
    				a[u][u]++, a[v][v]++, a[u][v]--, a[v][u]--;
    			}
    		}
    		int det=gauss(n-1);
    		ans=(ans+det*((n-pop)%2?1:-1))%mod;
    	}
    	printf("%lld
    ",(ans+mod)%mod);
    	return 0;
    }
    

    P4208 [JSOI2008]最小生成树计数

    首先对于一个图的所有最小生成树有一个性质:每种权值数量的边在不同的最小生成树中是一样的,而且每种边权的边所形成的联通状态也是一样的。于是,我们可以对于每种在kruskal中出现的不同的边权做生成树计数,然后再把他们乘起来。对于计算一种边权的答案,我们考虑缩点,把边权不是目前计算的边权的边的两点缩在一起,然后再跑一遍 Matrix-Tree 即可。

    #include<bits/stdc++.h>
    #define int long long
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=1009,mod=31011;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,m,mst[N],tot,ans=1,use[N],cnt;
    struct edge {int u,v,w;} e[N];
    bool cmp(const edge &a,const edge &b) {return a.w<b.w;}
    
    int id[N],p[N],pc;
    int find(int i) {return i==id[i]?i:id[i]=find(id[i]);}
    void unite(int u,int v) {id[find(u)]=find(v);}
    
    int a[109][109];
    int gauss(int n,int det=1) {
    	rep(i,1,n) {
    		rep(j,i+1,n) {
    			if(a[j][i]>a[i][i]) det=-det, swap(a[j],a[i]);
    			while(a[j][i]) {
    				int b=a[i][i]/a[j][i];
    				rep(k,i,n) a[i][k]=(a[i][k]+(mod-b)*a[j][k])%mod;
    				det=-det, swap(a[j],a[i]);
    			}
    		}
    		det=(det*a[i][i]%mod+mod)%mod;
    	}
    	return det;
    }
    
    signed main() {
    	n=read(), m=read();
    	rep(i,1,m) {
    		int a=read(), b=read(), c=read();
    		e[i]=(edge){a,b,c};
    	}
    	sort(e+1,e+m+1,cmp);
    	rep(i,1,n) id[i]=i;
    	rep(i,1,m) if(find(e[i].u)!=find(e[i].v)) {
    		unite(e[i].u,e[i].v);
    		mst[++tot]=i;
    		if(use[cnt]!=e[i].w) use[++cnt]=e[i].w;
    	}
    	rep(i,1,cnt) {
    		int w=use[i];
    		rep(i,1,n) id[i]=i;
    		rep(i,1,tot) {
    			if(e[mst[i]].w!=w) unite(e[mst[i]].u,e[mst[i]].v);
    		}
    		memset(a,0,sizeof(a));
    		pc=0;
    		rep(i,1,n) if(find(i)==i) p[i]=++pc;
    		rep(i,1,n) p[i]=p[id[i]];
    		rep(i,1,m) if(e[i].w==w) {
    			int u=p[e[i].u], v=p[e[i].v];
    			a[u][u]++, a[v][v]++, a[u][v]--, a[v][u]--;
    		}
    		ans=ans*gauss(pc-1)%mod;
    	}
    	printf("%lld
    ",ans);
    	return 0;
    }
    

    3.2 BEST 定理

    BEST 定理:对于有向图 (G),设其以 (k) 为根的根向生成树个数为 (t(G,k)),则其欧拉回路数为

    [ec(G)=t(G,k)prod_{uin V} (deg_u-1)! ]

    BZOJ3659. Which Dreamed It

    模板题。注意最终答案要乘上 (d_1)

    #include<bits/stdc++.h>
    #define int long long
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=109,mod=1000003;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int T,n,a[N][N],d[N],fac[200009];
    int gauss(int det=1) {
    	rep(i,2,n) {
    		rep(j,i+1,n) {
    			if(a[j][i]>a[i][i]) det=-det, swap(a[j],a[i]);
    			while(a[j][i]) {
    				int b=a[i][i]/a[j][i];
    				rep(k,i,n) a[i][k]=(a[i][k]+(mod-b)*a[j][k])%mod;
    				det=-det, swap(a[j],a[i]);
    			}
    		}
    		det=(det*a[i][i]%mod+mod)%mod;
    	}
    	return det;
    }
    
    signed main() {
    	T=read();
    	fac[0]=1;
    	rep(i,1,200000) fac[i]=fac[i-1]*i%mod;
    	while(T--) {
    		n=read();
    		memset(a,0,sizeof(a));
    		rep(i,1,n) {
    			d[i]=read();
    			rep(j,1,d[i]) {
    				int v=read();
    				a[i][i]++, a[i][v]--;
    			}
    		}
    		if(n==1&&d[1]==0) {puts("1");continue;}
    		int ans=gauss();
    		rep(i,1,n) ans=ans*fac[d[i]-1]%mod;
    		ans=ans*d[1]%mod;
    		printf("%lld
    ",ans);
    	}
    	return 0;
    }
    

    4 线性基

    线性基本质是维护一些向量的基底。但是由于常见的线性基都是做异或问题的。所以……

    对于二进制来说,原集合 T 的线性基是一个集合,满足线性基中的元素可以通过异或得到 T 中的每一个元素,且这个集合最小。

    性质

    • 线性基中的元素异或算得 T 中的每一个数的方案唯一
    • 线性基中的元素的二进制最高位不同(第i个数的最高位即i)
    • 线性基中的元素互相异或,其异或集合不变

    实现

    • 插入

    令插入的数位x,设其最高位是i,有两种情况。第一种,线性基中没有最高位是i的元素,那么直接在位置i插入即可。第二种,线性基中存在最高位是i的元素,那么我们考虑让x去异或这个线性基中的元素,并继续重复以上步骤。如果未能插入,代表我们已经可以用线性基的元素表示它了,就不需要插入了。

    • 求最小值

    显然,答案为线性基中的最小值。注意要特判 0。

    • 求最大值

    按位贪心,从高到低遍历线性基。假如到第i位时,答案中的i为0,那么我们异或上线性基中的i。

    洛谷 线性基 模板

    #include<bits/stdc++.h>
    #define int long long
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=59;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,s[N],lb[N];
    
    bool bit(int s,int i) {return s&(1ll<<i);} 
    void insert(int x) {
    	per(i,52,0)
    		if(bit(x,i)) {
    			if(!lb[i]) {lb[i]=x; return;}
    			x^=lb[i];
    		}
    }
    int query(int ret=0) {
    	per(i,50,0)
    		if(!bit(ret,i)) ret^=lb[i];
    	return ret;
    }
    
    signed main() {
    	n=read();
    	rep(i,1,n) s[i]=read(), insert(s[i]);
    	printf("%lld
    ",query());
    	return 0;
    }
    

    P4151 [WC2011]最大XOR和路径

    一个路径可以拆成一条简单链加上很多环。容易证明这条链的异或和异或上所有所选的环的异或和 等于答案路径的异或和。所以现在考虑如何找环。但是显然我们不可能把所有环都拎出来。但是考虑对于环套环(即两个环公用一条边,会形成一个新的大环等情况)的情况,我们实则只需要像仙人掌一样找到能cover所有边的一些环,其他环必定可以动过这些找到的环异或出来。

    所以我们在找到所有这样的环,异或和塞到线性基里面,再找一条链即可。

    #include<bits/stdc++.h>
    #define int long long
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=100009;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    struct edge {int to,nxt,w;} e[N<<1]; int hd[N],tot;
    void add(int u,int v,int w) {e[++tot]=(edge){v,hd[u],w};hd[u]=tot;}
    
    int n,m,s[N],lb[N];
    bool vst[N];
    
    bool bit(int s,int i) {return s&(1ll<<i);}
    void insert(int x) {
    	per(i,60,0)
    		if(bit(x,i)) {
    			if(!lb[i]) {lb[i]=x; return;}
    			x^=lb[i];
    		}
    }
    int query(int ret) {
    	per(i,60,0)
    		if(!bit(ret,i)) ret^=lb[i];
    	return ret;
    }
    
    void dfs(int u,int res) {
    	vst[u]=1, s[u]=res;
    	for(int i=hd[u],v;i;i=e[i].nxt) {
    		if(!vst[v=e[i].to]) dfs(v,res^e[i].w);
    		else insert(s[v]^res^e[i].w);
    	}
    }
    
    signed main() {
    	n=read(), m=read();
    	rep(i,1,m) {
    		int s=read(), t=read(), d=read();
    		add(s,t,d), add(t,s,d);
    	}
    	dfs(1,0);
    	printf("%lld
    ",query(s[n]));
    	return 0;
    }
    

    P4570 [BJWC2011]元素

    题意:找到一个最大的元素集合,使得其不存在任何子集的ai的异或和为0,且bi之和最大。

    线性基+贪心。很明显,如果一个数在插入线性基的时候发现已经能被表示,则这个数就不能算到集合中去。我们用贪心的想法,按bi从大往小插。容易证明贪心是对的,因为假如一个bi大的数字无法插入,我们可以剔除掉线性基中已有的bi小的数,使得答案更优。

    #include<bits/stdc++.h>
    #define int long long
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=1009;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,ans;
    struct stone {int a,b;} s[N];
    bool cmp(const stone &a,const stone &b) {return a.b>b.b;}
    
    int lb[N];
    bool bit(int s,int i) {return s&(1ll<<i);}
    bool insert(int x) {
    	per(h,62,0) if(bit(x,h)) {
    		if(!lb[h]) return lb[h]=x, 1;
    		x^=lb[h];
    	}
    	return 0;
    }
    
    signed main() {
    	n=read();
    	rep(i,1,n) s[i].a=read(), s[i].b=read();
    	sort(s+1,s+n+1,cmp);
    	rep(i,1,n) if(insert(s[i].a)) ans+=s[i].b;
    	printf("%lld
    ",ans);
    	return 0;
    }
    

    P3265 [JLOI2015]装备购买

    题意:找到一组线性无关的向量,使得其价值总和最小

    贪心思路和上题类似。但是这道题线性基回到了其本来的意思,存储的是一堆向量。和高斯消元解异或方程一样,异或上的线性基无非就是消元时只需要异或一下,而普通的线性基消元时就必须学着高斯消元的办法。

    这题tm卡精度。

    #include<bits/stdc++.h>
    #define rep(i,a,b) for(register int i=(a);i<=(b);i++)
    #define per(i,a,b) for(register int i=(a);i>=(b);i--)
    using namespace std;
    const int N=509;
    const double eps=1e-4;
    
    inline int read() {
        register int x=0, f=1; register char c=getchar();
        while(c<'0'||c>'9') {if(c=='-') f=-1; c=getchar();}
        while(c>='0'&&c<='9') {x=(x<<3)+(x<<1)+c-48,c=getchar();}
        return x*f;
    }
    
    int n,m,ans1,ans2;
    struct eq {vector<long double>z; int val;} a[N];
    bool cmp(const eq &a,const eq &b) {return a.val<b.val;}
    
    long double lb[509][509];
    bool insert(vector<long double>x) {
    	rep(i,1,m) if(fabs(x[i])>eps) {
    		if(fabs(lb[i][i])<eps) {
    			rep(j,i,m) lb[i][j]=x[j];
    			return 1;
    		}
    		double k=x[i]/lb[i][i];
    		rep(j,i,m) x[j]-=k*lb[i][j];
    	}
    	return 0;
    }
    
    int main() {
    	n=read(), m=read();
    	rep(i,1,n) {
    		a[i].z.push_back(0);
    		rep(j,1,m) a[i].z.push_back(read());
    	}
    	rep(i,1,n) a[i].val=read();
    	sort(a+1,a+n+1,cmp);
    	rep(i,1,n) if(insert(a[i].z)) ans1++,ans2+=a[i].val;
    	printf("%d %d
    ",ans1,ans2);
    	return 0;
    }
    
  • 相关阅读:
    转 Linux查看版本信息及CPU内核、型号等
    freeswitch ODBC error: ODBC NOT AVAILABLE!
    asterisk 命令
    Freeswitch mod 安装
    数据库压缩备份
    IEnumreable的使用小结
    新的Layout布局系统
    前台网站开发手记
    容器服务是如何在高速增长下保持高可用性的
    Kubernetes问题排查流程
  • 原文地址:https://www.cnblogs.com/TetrisCandy/p/14162624.html
Copyright © 2011-2022 走看看