【学习笔记】数学知识-高斯消元与线性空间
友情提示
- 本博客内部分内容因缺乏样例，可能较为晦涩难懂，建议一并参考蓝书，数论小白都能看懂的线性方程组及其解法，OI Wiki。
线性方程组
- 线性方程组是形如 \(\begin{cases}a_{1,1} x_1 + a_{1,2} x_2 + a_{1,3} x_3 + \dots + a_{1,n} x_n &=a_{1,n+1}\\a_{2,1} x_1 + a_{2,2} x_2 + a_{2,3} x_3 + \dots + a_{2,n} x_n &=a_{2,n+1}\ \\ \dots &= \dots \\a_{m,1} x_1 + a_{m,2} x_2 + a_{m,3} x_3 + \dots + a_{m,n} x_m &=a_{m,n+1}\end{cases}\) ，由 \(m\) 个 \(n\) 元一次方程共同构成的。
- 线性方程组的所有系数可以写成一个 \(m\) 行 \(n\) 列的系数矩阵，再加上每个方程等号右侧的常数，可以写成 一个 \(m\) 行 \(n+1\) 列的增广矩阵，形如： \(\begin{array}{|c c c c c|c|}a_{1,1}&a_{1,2}&a_{1,3}& \dots &a_{1,n}&a_{1,n+1}\\a_{2,1}&a_{2,2}&a_{2,3}& \dots &a_{2,n}&a_{2,n+1}\\ \dots & \dots & \dots & \dots & \dots & \dots \\a_{m,1}&a_{m,2}&a_{m,3}& \dots &a_{m,n}&a_{m,n+1}\end{array}\) 。
- 对增广矩阵的三类操作（初等行变化）：
- 用一个非零的数乘某一行。
- 把其中一行的若干倍加到另一行上。
- 交换两行的位置。
- 用若干次初等行变化求解方程组，最后得到的矩阵被称为阶梯形矩阵，它的系数矩阵部分被称为上三角矩阵，从下往上依次带回方程组，即可得到每个未知数的解；将其进一步化简最后得到的矩阵被称为简化阶梯形矩阵，它的系数矩阵部分是一个对角矩阵，该矩阵直接给出了方程组的解（这部分可能有些难懂，建议参考蓝书）。
高斯消元
定义
- 通过初等行列变化把增广矩阵变为简化阶梯形矩阵的线性方程组求解算法就是高斯消元算法。
- 思想：对于每个未知数 \(x_i\) ，找到一个 \(x_i\) 的系数非零但 \(x_1 \sim x_{i-1}\) 的系数都为零的方程，然后用初等行变化把其他方程的 \(x_i\) 的系数全部消成零。
- 特殊情形：
- 在高斯消元过程中，可能出现 \(0x_{i}=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(1 \le k<n)\) 个，则说明有 \(k\) 个主元， \(n-k\) 个自由元，原方程组存在多解。
消元
回代消元
- 将增广矩阵化作阶梯型矩阵后，自下而上减去已求出的未知数乘以其系数，求出此时要求的未知数。
约旦消元
- 在消元的过程中，对主元上、下方的式子均进行初等行变化，使得最终第 \(i\) 行只有主元的系数不为 \(0\) ，接着相除即可。
例题
luogu P3389 【模板】高斯消元法
-
板子。
点击查看回代消元代码
const double eps=1e-10; double g[110][110],ans[110]; void Gauss(int n,double ans[]) { for(int i=1;i<=n;i++) { int val=i; for(int j=i+1;j<=n;j++) { if(fabs(g[j][i])-fabs(g[val][i])>=eps) { val=j; } } for(int j=1;j<=n+1;j++) { swap(g[i][j],g[val][j]); } if(fabs(g[i][i])<eps) { cout<<"No Solution"<<endl; exit(0); } for(int j=i+1;j<=n;j++) { for(int k=i+1;k<=n+1;k++) { g[j][k]-=g[i][k]*g[j][i]/g[i][i]; } g[j][i]=0;//这一行可以删去，因为后续不会再访问 g[j][i]，且 g[j][i] 后续一直为 0 } } for(int i=n;i>=1;i--) { ans[i]=g[i][n+1]; for(int j=i+1;j<=n;j++) { ans[i]-=g[i][j]*ans[j]; } ans[i]/=g[i][i]; } } int main() { int n,i,j; cin>>n; for(i=1;i<=n;i++) { for(j=1;j<=n+1;j++) { cin>>g[i][j]; } } Gauss(n,ans); for(i=1;i<=n;i++) { printf("%.2lf\n",ans[i]); } return 0; }点击查看约旦消元代码
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; }
luogu P2455 [SDOI2006] 线性方程组
-
板子。
点击查看 hack 数据
in: 6 1 2 3 4 5 6 10 0 0 1 2 3 4 15 0 0 0 2 6 7 30 0 0 0 0 6 8 10 0 0 0 0 0 1 30 0 0 0 0 0 0 0 ans:0 -
回代消元
- 方程组中可能会存在无用方程（在消元过程中主元系数为 \(0\) 的方程），我们可以选择记录 \(r\) 表示当前正在用第 \(r\) 个方程求解。最终，若 \(r=n+1\) 说明有唯一解，否则无解或无数组解。
点击查看回代消元代码
const double eps=1e-10; double g[110][110],ans[110]; void Gauss(int n,double ans[]) { int r=1; for(int i=1;i<=n;i++) { int val=r; for(int j=r+1;j<=n;j++) { if(fabs(g[j][i])-fabs(g[val][i])>=eps) { val=j; } } for(int j=1;j<=n+1;j++) { swap(g[r][j],g[val][j]); } if(fabs(g[r][i])>=eps) { for(int j=r+1;j<=n;j++) { for(int k=i+1;k<=n+1;k++) { g[j][k]-=g[r][k]*g[j][i]/g[r][i]; } g[j][i]=0; } r++; } } if(r<=n) { for(int i=r;i<=n;i++) { if(fabs(g[i][n+1])>=eps) { cout<<"-1"<<endl; exit(0); } } cout<<"0"<<endl; exit(0); } else { for(int i=n;i>=1;i--) { ans[i]=g[i][n+1]; for(int j=i+1;j<=n;j++) { ans[i]-=g[i][j]*ans[j]; } ans[i]/=g[i][i]; } } } int main() { int n,i,j; cin>>n; for(i=1;i<=n;i++) { for(j=1;j<=n+1;j++) { cin>>g[i][j]; } } Gauss(n,ans); for(i=1;i<=n;i++) { printf("x%d=%.2lf\n",i,fabs(ans[i])<eps?0:ans[i]); } return 0; } -
约旦消元
- 上题的做法无法判断无解和无数解。消元顺序会影响答案。主要原因是因为把有用的式子放在了无用的位置，且后续不再来此判断是否有用。
- 上题中是从 \(i\) 开始枚举的，有用式子仅为 \(i\) 之后的，当遇到主元系数为 \(0\) 直接退出。但是本题需要判断是否有解和存在多解，所以需要在系数为 \(0\) 后继续算下去，因为有用的式子还包括 \(i\) 之前且系数为 \(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++) { if((i<=j||fabs(a[j][j])<eps)&&fabs(a[j][i])-fabs(a[val][i])>=eps)//系数为 0 的也应该继续算下去 { val=j; } } for(j=1;j<=n+1;j++) { swap(a[i][j],a[val][j]); } if(fabs(a[i][i])>=eps) { 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(fabs(a[i][i])<eps)//精度 { if(fabs(a[i][n+1])>=eps)//无解比无数解优先级更高 { 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; }
luogu P4035 [JSOI2008] 球形空间产生器
-
设球心的 \(n\) 维坐标分别为 \(x_1,x_2, \dots ,x_n\) ，此时有对于每一个 \(i(1 \le i \le n+1)\) 均满足 \(\sum\limits_{j=1}^{n}(a_{i,j}-x_j)^2=dis\) ，其中 \(dis\) 为常数。
-
将相邻两个方程作差，得到对于每一个 \(i(1 \le i \le n)\) 均满足 \(\sum\limits_{j=1}^{n}(a_{i,j}-x_j)^2-(a_{i+1,j}-x_j)^2=0\) ，去括号并移项得 \(\sum\limits_{j=1}^{n}2(a_{i,j}-a_{i+1,j})x_j=\sum\limits_{j=1}^{n} a_{i,j}^2-a_{i+1,j}^2\) 。接着高斯消元即可。
点击查看代码
const double eps=1e-12; double a[50][50],b[50][50]; int main() { int n,i,j,k,val; double ans=0; cin>>n; for(i=1;i<=n+1;i++) { for(j=1;j<=n;j++) { cin>>b[i][j]; } } for(i=1;i<=n;i++) { for(j=1;j<=n;j++) { a[i][j]=2*(b[i][j]-b[i+1][j]); a[i][n+1]+=b[i][j]*b[i][j]-b[i+1][j]*b[i+1][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++) { 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++) { printf("%.3lf ",a[i][n+1]/a[i][i]); } return 0; }
luogu P2973 [USACO10HOL] Driving Out the Piggies G
-
有后效性的 \(DP\) ，可以使用高斯消元求解。
-
设 \(x_i\) 表示炸弹在第 \(i\) 个点爆炸的概率， \(din_i\) 表示 \(i\) 的入度。
- 当 \(i=1\) 时，炸弹可能一开始就爆炸或在与之相连的点不爆炸的情况下，从与之相连的点转移过来，有 \(x_i=\dfrac{p}{q}+\sum\limits_{(i,j)\in E}{} x_j \times(1-\dfrac{p}{q}) \times \dfrac{1}{din_j}\) 。
- 当 \(i \ne 1\) 时，炸弹只会在与之相连的点不爆炸的情况下，从与之相连的点转移过来，有 \(x_i=\sum\limits_{(i,j)\in E}{} x_j \times(1-\dfrac{p}{q}) \times \dfrac{1}{din_j}\) 。
-
建立线性方程组：\(\begin{cases}x_1-\sum\limits_{(1,j)\in E}{} x_j \times(1-\dfrac{p}{q}) \times \dfrac{1}{din_j}&=\dfrac{p}{q}\\x_2-\sum\limits_{(2,j)\in E}{} x_j \times(1-\dfrac{p}{q}) \times \dfrac{1}{din_j}&=0\ \\ \dots &= \dots \\x_n-\sum\limits_{(n,j)\in E}{} x_j \times(1-\dfrac{p}{q}) \times \dfrac{1}{din_j}&=0\end{cases}\)
-
接着进行高斯消元即可。
点击查看代码
const double eps=1e-12; double a[305][305]; int b[305][305],din[305]; int main() { int n,m,i,j,k,u,v,val; double p,q; cin>>n>>m>>p>>q; for(i=1;i<=m;i++) { cin>>u>>v; b[u][v]=b[v][u]=1; din[u]++; din[v]++; } a[1][n+1]=p/q; for(i=1;i<=n;i++) { a[i][i]=1; for(j=1;j<=n;j++) { if(b[i][j]==1) { a[i][j]=-(1.0-p/q)*1.0/din[j]; } } } for(i=1;i<=n;i++) { val=i; for(j=i;j<=n;j++) { 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++) { printf("%.9lf\n",a[i][n+1]/a[i][i]); } return 0; }
应用
求解异或方程
定义
- 异或方程组：形如 \(\begin{cases}a_{1,1} x_1 \bigoplus a_{1,2} x_2 \bigoplus a_{1,3} x_3 \bigoplus \dots \bigoplus a_{1,n} x_n &=a_{1,n+1}\\a_{2,1} x_1 \bigoplus a_{2,2} x_2 \bigoplus a_{2,3} x_3 \bigoplus \dots \bigoplus a_{2,n} x_n &=a_{2,n+1}\ \\ \dots &= \dots \\a_{m,1} x_1 \bigoplus a_{m,2} x_2 \bigoplus a_{m,3} x_3 \bigoplus \dots \bigoplus a_{m,n} x_m &=a_{m,n+1}\end{cases}\) 的方程组。
- 异或性质
- 归零率：\(a \bigoplus a=0\)
- 恒等律：\(a \bigoplus 0=a\)
- 交换律：\(a \bigoplus b=b \bigoplus a\)
- 结合律: \(a \bigoplus b \bigoplus c=a \bigoplus (b \bigoplus c)=(a \bigoplus b) \bigoplus c\)
- 自反性(异或的逆运算为它本身): \(a \bigoplus b \bigoplus b=a\)
- 用高斯消元法解异或方程组时应使用异或消元而非加减消元，且不需要进行乘除改变系数（因为系数只有 \(1\) 和 \(0\) ）。
- 可用
bitset优化。
- 异或性质
例题
luogu P2962 [USACO09NOV] Lights G
-
板子。
点击查看代码
int a[50][50],vis[50],ans=0; void dfs(int x,int num,int n) { if(num>=ans) { return; } else { if(x==0) { ans=min(ans,num); return; } else { if(a[x][x]==1)//如果是主元 { int ls=a[x][n+1],i;//回代 for(i=x+1;i<=n;i++) { if(a[x][i]==1) { ls^=vis[i]; } } dfs(x-1,num+ls,n); } else//否则是自由元 { vis[x]=0;//不进行操作 dfs(x-1,num,n); vis[x]=1;//进行操作 dfs(x-1,num+1,n); } } } } int main() { int n,m,i,j,k,val,l,r,flag=0; cin>>n>>m; for(i=1;i<=n;i++) { a[i][n+1]=a[i][i]=1;//第i个方程的第i个系数初始化为1 //一个点最后的状态取决于最初位置 a[i][n+1]=0^1 } for(i=1;i<=m;i++) { cin>>l>>r; a[l][r]=a[r][l]=1;//若l与r有连边，则将第l个方程的第r个系数与第r个方程的第l个系数初始化为1 } for(i=1;i<=n;i++) { val=i; while(val<=n&&a[val][i]==0) { val++; } if(val!=n+1) { for(j=1;j<=n+1;j++) { swap(a[i][j],a[val][j]); } for(j=1;j<=n;j++) { if(j!=i&&a[j][i]==1)//保证必须是主元 { for(k=i+1;k<=n+1;k++) { a[j][k]^=a[i][k]; } } } } else//说明存在自由元 { flag=1; } } if(flag==0)//说明不存在自由元 { for(i=1;i<=n;i++) { ans+=a[i][n+1]; } cout<<ans<<endl; } else { ans=0x7f7f7f7f; dfs(n,0,n);//进行回带 cout<<ans<<endl; } return 0; }
luogu P2447 [SDOI2010] 外星千足虫
-
令每只地球千足虫的足的总数为 \(0\) ，每只外星千足虫的足的总数为 \(1\) ， \(a_{i,j}(1 \le i \le m,1 \le j \le n)\) 表示第 \(i\) 次称量时第 \(j\) 只外星千足虫是否被放入机器（ \(1\) 表示放入， \(0\) 表示不放入），此时有同余方程组：\(\begin{cases}a_{1,1} x_1 + a_{1,2} x_2 + a_{1,3} x_3 \dots + a_{1,n} x_n &\equiv a_{1,n+1} \pmod 2\\a_{2,1} x_1 + a_{2,2} x_2 + a_{2,3} x_3 \dots + a_{2,n} x_n &\equiv a_{2,n+1} \pmod 2\ \\ \dots &\equiv \dots \pmod 2\\a_{m,1} x_1 + a_{m,2} x_2 + a_{m,3} x_3 \dots + a_{m,n} x_n &\equiv a_{m,n+1} \pmod 2\end{cases}\) ，考虑到 \(\bmod 2\) 的特殊性质，又可有异或同余方程组：\(\begin{cases}a_{1,1} x_1 \bigoplus a_{1,2} x_2 \bigoplus a_{1,3} x_3 \dots \bigoplus a_{1,n} x_n &\equiv a_{1,n+1} \pmod 2\\a_{2,1} x_1 \bigoplus a_{2,2} x_2 \bigoplus a_{2,3} x_3 \dots \bigoplus a_{2,n} x_n &\equiv a_{2,n+1} \pmod 2\ \\ \dots &\equiv \dots \pmod 2\\a_{m,1} x_1 \bigoplus a_{m,2} x_2 \bigoplus a_{m,3} x_3 \dots \bigoplus a_{m,n} x_n &\equiv a_{m,n+1} \pmod 2\end{cases}\) ，又因为本题的 \(a_{i,j},x_i \in [0,1](1 \le i \le m,1 \le j \le n)\) ，可得到本题方程：\(\begin{cases}a_{1,1} x_1 \bigoplus a_{1,2} x_2 \bigoplus a_{1,3} x_3 \dots \bigoplus a_{1,n} x_n &= a_{1,n+1} \\a_{2,1} x_1 \bigoplus a_{2,2} x_2 \bigoplus a_{2,3} x_3 \dots \bigoplus a_{2,n} x_n &= a_{2,n+1} \ \\ \dots &= \dots \\a_{m,1} x_1 \bigoplus a_{m,2} x_2 \bigoplus a_{m,3} x_3 \dots \bigoplus a_{m,n} x_n &= a_{m,n+1} \end{cases}\) ，高斯消元即可。
- 数据可能没造满， \(n \le 1000\) 的 \(n^3\) 跑过去了 ，如果不放心就用 \(bitset\) 优化一下。
- 可能是 \(luogu\) 的测评机太快了，要是 \(AcCoders\) 的估计会被卡。
点击查看代码
bitset<2500> a[2500]; int main() { int n,m,i,j,ls,val,flag=0,ans=0; char pd; cin>>n>>m; for(i=1;i<=m;i++) { for(j=1;j<=n;j++) { cin>>pd; a[i][j]=pd-'0'; } cin>>ls; a[i][n+1]=ls; } for(i=1;i<=n;i++) { val=i; while(val<=m&&a[val][i]==0) { val++; } if(val!=m+1) { ans=max(ans,val);//二分查找ans但会T，但是有每个val均满足a[val][i]=1且val最小，取max即可 swap(a[i],a[val]); for(j=1;j<=m;j++) { if(j!=i&&a[j][i]==1) { a[j]^=a[i]; } } } else//说明有自由元，存在多解 { flag=1; cout<<"Cannot Determine"<<endl; break; } } if(flag==0) { cout<<ans<<endl; for(i=1;i<=n;i++) { if(a[i][n+1]==1) { cout<<"?y7M#"<<endl; } else { cout<<"Earth"<<endl; } } } return 0; } - 数据可能没造满， \(n \le 1000\) 的 \(n^3\) 跑过去了 ，如果不放心就用 \(bitset\) 优化一下。
luogu P3164 [CQOI2014] 和谐矩阵
-
题目要求对于任意一个点 \((i,j)\) ，其中 \(1 \le i\le n,1 \le j \le m\)，均有 \(x_{i-1,j} \times [1 \le i-1] \bigoplus x_{i,j-1} \times [1 \le j-1] \bigoplus x_{i,j} \bigoplus x_{i+1,j} \times [i+1 \le n] \bigoplus x_{i,j+1} \times [j+1 \le m]=0\) 。
-
考虑建立 \(nm\) 元线性异或方程组，\(\begin{cases}a_{1,1} x_1 \bigoplus a_{1,2} x_2 \bigoplus a_{1,3} x_3 \dots \bigoplus a_{1,nm} x_{nm} &= 0 \\a_{2,1} x_1 \bigoplus a_{2,2} x_2 \bigoplus a_{2,3} x_3 \dots \bigoplus a_{2,nm} x_{nm} &= 0 \ \\ \dots &= \dots \\a_{nm,1} x_1 \bigoplus a_{nm,2} x_2 \bigoplus a_{nm,3} x_3 \dots \bigoplus a_{nm,nm} x_{nm} &= 0 \end{cases}\) ，其中如果一个点 \((i',j')\) 与 \((i,j)\) 相邻或者 \(i=i',j=j'\) ，则有 \(a_{i',j'}=1\) ，否则 \(a_{i',j'}=0\) ，接着进行高斯消元，然后求解。
- 为了防止全零矩阵，自由元取 \(1\) 即可。
- \(n,m \le 40\) ，需要 \(bitset\) 优化。
点击查看代码
int dir[5][2]={{0,0},{0,-1},{0,1},{-1,0},{1,0}},ans[2000]; bitset<2000>a[2000]; int main() { int n,m,i,j,k,x,y,val; cin>>n>>m; for(i=1;i<=n;i++) { for(j=1;j<=m;j++) { for(k=0;k<=4;k++) { x=i+dir[k][0]; y=j+dir[k][1]; if(1<=x&&x<=n&&1<=y&&y<=m) { a[(i-1)*m+j][(x-1)*m+y]=1; } } } a[i][n*m+1]=0; } for(i=1;i<=n*m;i++) { val=i; while(val<=n*m&&a[val][i]==0) { val++; } if(val!=n*m+1) { swap(a[i],a[val]); for(j=1;j<=n*m;j++) { if(j!=i&&a[j][i]==1) { a[j]^=a[i]; } } } else//自由元直接取1 { ans[i]=1; } } for(i=n*m;i>=1;i--)//求解 { for(j=i+1;j<=n*m;j++) { ans[i]^=ans[j]*a[i][j]; } } for(i=1;i<=n*m;i++) { cout<<ans[i]<<" "; if(i%m==0) { cout<<endl; } } return 0; }
luogu P10499 开关问题
-
判断有无解，用类似 luogu P2455 [SDOI2006] 线性方程组 的方法进行处理。
-
若存在 \(num\) 个自由元，易知自由元仅可取 \(0\) 或 \(1\) ，故答案总数为 \(2^{num}\) 。
点击查看代码
int start[50],en[50],a[50][50]; int main() { int t,n,i,j,k,l,r,ii,h,val,ans; cin>>t; for(h=1;h<=t;h++) { cin>>n; ans=0; memset(a,0,sizeof(a)); for(i=1;i<=n;i++) { cin>>start[i]; } for(i=1;i<=n;i++) { cin>>en[i]; a[i][i]=1; a[i][n+1]=start[i]^en[i]; } while(cin>>l>>r) { if(l==0&&r==0) { break; } else { a[r][l]=1; } } for(i=1;i<=n;i++) { val=i; for(j=1;j<=n;j++) { if(!(i>j&&a[j][j]!=0)) { if(a[j][i]>a[val][i]) { val=j; } } } for(j=1;j<=n+1;j++) { swap(a[i][j],a[val][j]); } for(j=1;j<=n;j++) { if(j!=i&&a[j][i]==1) { for(k=i+1;k<=n+1;k++) { a[j][k]^=a[i][k]; } } } } for(i=1;i<=n;i++) { if(a[i][i]==0) { if(a[i][n+1]!=0) { ans=-1; break; } else { ans++; } } } if(ans==-1) { cout<<"Oh,it's impossible~!!"<<endl; } else { cout<<(1<<ans)<<endl; } } return 0; }
求解带状矩阵
定义
- 带状矩阵又称 Band-Matrix 。
- 对于 \(n \times n\) 的方阵 \(A\) ，若它的全部非零元素落在一个以 主对角线 为中心的带状区域内，这个带状区域包含主对角线，以及主对角线上下各 \(d\) 条对角线上的元素，则称该方阵为半带宽为 \(d\) 的带状矩阵，形如 \(A=\begin{bmatrix} a_{1,1} & a_{1,2} & 0 & 0 & 0 & 0 & 0 \\ a_{2,1} & a_{2,2} & a_{2,3} & 0 & 0 & 0 & 0 \\ 0 & a_{3,2} & a_{3,3} & a_{3,4} & 0 & 0 & 0 \\ 0 & 0 & a_{4,3} & a_{4,4} & a_{4,5} & 0 & 0 \\ \vdots & \ddots & \ddots & \ddots & \ddots & \ddots & \vdots \\ \vdots & \dots & \dots & \dots & \dots & \dots & \vdots \end{bmatrix}\) 。
- 对于半带宽为 \(d\) 的带状矩阵 \(A\) 中的任意非零元素 \(a_{i,j}\) 均有 \(|i-j| \le d\) ，即 \(a_{i,i}\) 向右和向下拓展有不超过 \(d-1\) 个非零数字。
- 对于带状矩阵 \(A\) 若存在最小的正整数 \(m\) 使得 \(|i-j| \ge m\) 时均有 \(a_{i,j}=0\) ，则称 \(w=2m-1\) 为带状矩阵 \(A\) 的带宽， \(d=2w-1\) 为带状矩阵 \(A\) 的半带宽。
消元
- 在带状矩阵中有很多位置根本不需要进行消元。
- 设（增广）矩阵为 \(G\) ，消元方法基本同普通的高斯消元（回代消元版），但当主元 \(G_{i,i}\) 为 \(0\) 时的处理方式不一样，常见的处理方式有如下两种，代码均为第一种写法。
- 仍选择交换行，即将第 \(i\) 行与它下面的主元非 \(0\) 的某行进行交换。这样的话就破坏了带状的性质，每次右边最多多出来 \(d\) 个数，保险起见每次往右消元 \(2d\) 个数即可。
- 进行换元，选择交换列。需要记录 \(from_{i}\) 表示第 \(i\) 列表示的是最开始的第 \(i\) 个元。
- 判断是否有唯一解和无解同普通高斯消元。
例题
CF24D Broken robot
线性基
定义
- 线性空间 是一个关于以下两个运算封闭的向量集合：
- 向量加法 \(a+b\) ，其中 \(a,b\) 均为向量。
- 标量乘法（数乘运算） \(k \times a\) ，其中 \(a\) 为向量， \(k\) 为常数（标量）。
- 给定若干个向量 \(a_{1},a_{2}, \dots,a_{k}\) ，若向量 \(b\) 能被 \(a_{1},a_{2}, \dots,a_{k}\) 通过向量加法和标量乘法得到，则称向量 \(b\) 能被向量 \(a_{1},a_{2}, \dots,a_{k}\) 线性表出 。另外 \(a_{1},a_{2}, \dots,a_{k}\) 能表出的所有向量构成一个线性空间， \(a_{1},a_{2}, \dots,a_{k}\) 被称为这个线性空间的 生成子集 。
- 任意选出线性空间的若干个向量，若其中存在一个向量能被其他向量线性表出，则称这些向量 线性相关 ，否则称这些向量 线性无关 。
- 线性无关的生成子集被称为线性空间的 线性基 ，简称 基 。关于 线性基 的另一个定义是线性空间的极大线性无关子集。
- 一个线性空间的所有线性基包含的向量个数都相等，这个数偶被称为线性空间的 维数 。
- 在 OI 中，一般称 \(n\) 维实线性空间 \(\mathbf{R}^{n}\) 下的线性基称为 实数线性基 ， \(n\) 维布尔域线性空间 \(\mathbf{Z}_{2}^{n}\) 下的线性基称为 异或线性基 。在无特殊说明下，线性基一般代指异或线性基。
- 以异或线性基为例，我们可以根据给定的一组序列 \(X=\{ x_{1},x_{2}, \dots ,x_{m} \}\) 构造出一组线性基 \(B=\{ b_{1},b_{2}, \dots ,b_{k} \}\) ，使得 \(B\) 满足以下性质：
- \(B\) 中任意非空子集的异或和不为 \(0\) 。
- 对于 \(X\) 中的任意元素 \(x\) ，都可以在 \(B\) 中取出若干元素使其异或和为 \(x\) 。
- 对于任意满足以上两条性质的集合 \(B'\) ，其元素个数不会小于 \(B\) 的元素个数。
- 符合线性基中 极大 的定义。
- 因此可以用来实现
- 判断一个数能否表示成某数集子集的异或和。
- 求一个数表示成某数集子集异或和的方案数。
- 求某数集子集的最大/最小/第 \(k\) 大/第 \(k\) 小异或和。
- 求一个数在某数集子集异或和中的排名。
- 以异或线性基为例，我们可以根据给定的一组序列 \(X=\{ x_{1},x_{2}, \dots ,x_{m} \}\) 构造出一组线性基 \(B=\{ b_{1},b_{2}, \dots ,b_{k} \}\) ，使得 \(B\) 满足以下性质：
- 对于一个 \(n\) 行 \(m\) 列的矩阵，可以把它的每一行看作一个长度为 \(m\) 的向量，称为 行向量。矩阵的 \(n\) 个行向量能够表出的所有向量构成一个线性空间，这个线性空间的维数被称为矩阵的 行秩 。类似地，我们可以定义 列向量 和 列秩 。时间上矩阵的行秩一定等于列秩，它们都被称为矩阵的 秩 。
- 进一步地，把这个 \(n \times m\) 的矩阵看做系数矩阵（增广矩阵的最后一列全看作零）进行高斯消元，可以得到一个简化阶梯型矩阵。因初等行变换就是行向量之间进的向量加法与标量乘法运算，所以高斯消元不改变矩阵的行变量表出的线性空间，所有非零行向量线性无关，而这些所有非零行向量就是该线性空间的一个线性基，非零行向量的个数就是矩阵的秩。
构造
- 以构造异或线性基为例；实数线性基仅需做一些简单修改，详见 luogu P3265 [JLOI2015] 装备购买 题解 。
- 贪心法
- 对原数集的每个数 \(p\) 转为二进制后从高往地位扫，对于第 \(x\) 位是 \(1\) 的，若 \(b_{x}\) 不存在，则令 \(b_{x}=p\) 并结束扫描；若存在，则令 \(p=p \bigoplus b_{x}\) 。
- 实质是若集合中存在两个元素 \(x,y\) 不妨把它们替换成 \(x, x \bigoplus y\) 。
- 贪心法构成的线性基满足以下性质：
- 线性基没有异或和为 \(0\) 的子集。
- 线性基中各数二进制最高位不同。
- 查询原数集合中最大异或和，只需要将线性基从高往低位扫，若异或上当前扫到的 \(b_{x}\) 能使答案变大，则把答案异或上 \(b_{x}\) 。
- 查询原数集合中最小异或和就是线性基集合中所有元素中最小的那个。
- 判断一个数 \(p\) 能否表示成某数集子集的异或和，只需要类似插入的过程，若最后插入的 \(p\) 被异或成了 \(0\) ，则能被异或出来；否则在中途的过程中一定存在一个位置的 \(b_{x}\) 等于 \(0\) ，说明不能异或出来。
- 代码实现
-
点击查看代码
ll d[70]; void insert(ll x) { for(ll i=60;i>=0;i--) { if((x>>i)&1) { if(d[i]==0) { d[i]=x; break; } x^=d[i]; } } } int main() { ll n,x,ans=0,i; cin>>n; for(i=1;i<=n;i++) { cin>>x; insert(x); } for(i=60;i>=0;i--) { ans=max(ans,ans^d[i]); } cout<<ans<<endl; return 0; }
-
- 对原数集的每个数 \(p\) 转为二进制后从高往地位扫，对于第 \(x\) 位是 \(1\) 的，若 \(b_{x}\) 不存在，则令 \(b_{x}=p\) 并结束扫描；若存在，则令 \(p=p \bigoplus b_{x}\) 。
- 高斯消元法
- 类似异或高斯消元的写法。
- 高斯消元后的矩阵是一个行简化阶梯形矩阵，包含了贪心法构造的线性基满足的两条基础。
- 所以查询原数集合中最大异或和时只需要将线性基内的所有元素异或起来输出即可。
- 代码实现
-
点击查看代码
ll x[70],d[70]; ll gauss(ll n) { ll row=1; for(ll i=1;i<=n;i++) { d[i]=x[i]; } for(ll col=60;col>=0&&row<=n;col--) { for(ll i=row;i<=n;i++) { if((d[i]>>col)&1) { swap(d[row],d[i]); break; } } if((d[row]>>col)&1) { for(ll i=1;i<=n;i++) { if(i!=row&&((d[i]>>col)&1)) { d[i]^=d[row]; } } row++; } } row--; return row; } int main() { ll n,row,ans=0,i; cin>>n; for(i=1;i<=n;i++) { cin>>x[i]; } row=gauss(n); for(i=1;i<=row;i++) { ans^=d[i]; } cout<<ans<<endl; return 0; }
-
- 设向量长度为 \(n\) ，总数为 \(m\) ，则构造异或线性基的时间复杂度为 \(O(nm)\) ，其中高斯消元法的常数略大；构造实数线性基的时间复杂度为 \(O(n^{2}m)\) 。
- 在总数为 \(n\) ，值域为 \(V\) 时，构造异或线性基的时间复杂度为 \(O(n \log V)\) 。
例题
luogu P3857 [TJOI2008] 彩灯
luogu P4570 [BJWC2011] 元素
luogu P4301 [CQOI2013] 新Nim游戏
luogu P4151 [WC2011] 最大XOR和路径
HDU3949 XOR
luogu P3265 [JLOI2015] 装备购买
[ABC283G] Partial Xor Enumeration
CF959F Mahmoud and Ehab and yet another xor task
BZOJ4184 shallot
luogu P5556 圣剑护符
luogu P4839 P 哥的桶
CF1100F Ivan and Burgers
luogu P3292 [SCOI2016] 幸运数字
luogu P5607 [Ynoi2013] 无力回天 NOI2017
luogu P4869 albus就是要第一个出场
BZOJ3811 玛里苟斯
luogu P10778 BZOJ3569 DZY Loves Chinese II
BZOJ2322 梦想封印
本文来自博客园，作者：hzoi_Shadow，原文链接：https://www.cnblogs.com/The-Shadow-Dragon/p/17644605.html，未经允许严禁转载。
版权声明：本作品采用 「署名-非商业性使用-相同方式共享 4.0 国际」许可协议（CC BY-NC-SA 4.0） 进行许可。
