https://zybuluo.com/ysner/note/1106886
本质
模拟加减消元，先消成上三角再代入求解。
据说是小学内容？？？
列方程
有这个量的影响就赋系数。
实现（加减方程组）
- step 1：交换（将当前列系数绝对值最大的放到当前行）
- step 2：将每行第一系数化为1
- step 3：消成倒三角（在下面的方程中减去上面的方程）
- step 4：回代求解
il void Gauss()
{
fp(i,1,n)//枚举列
{
re int now=i;
fp(j,i+1,n)//枚举行
if(fabs(a[j][i])+eps>fabs(a[now][i])) now=j;
if(now^i) swap(a[i],a[now]);
fp(j,i+1,n)
fq(k,n+1,i)
a[j][k]-=a[i][k]*a[j][i]/a[i][i];
}
fq(i,n,1)
{
ans[i]=a[i][n+1];
fq(j,n,i+1) ans[i]-=ans[j]*a[i][j];
ans[i]/=a[i][i];
}
}
最关键的就是(a[j][k]-=a[i][k]*a[j][i]/a[i][i])了。
此时上面所有方程第一系数在(/a[i][i])后已经被视为(1)。
想想平时怎么解方程的。
先把上面方程乘上当前方程第一系数（即(*a[j][i]))，再同位项相减即可。
实现（异或方程组)
解这种方程不难，和解加减方程差不多，原理是不断消去倒三角左边的系数。一开始消下面的系数，往回代时消上面的。
除了continue，代码思路一致。
il void Gauss()
{
fp(i,1,n)
{
re int now=i;
fp(j,i+1,n) if(a[j][i]>a[now][i]) now=j;
if(now^i) swap(a[now],a[i]);
if(!a[i][i]) continue;
fp(j,i+1,n) if(a[j][i]) a[j]^=a[i];
}
fq(i,n,1)
fq(j,i-1,1)
if(a[j][i]) a[j]^=a[i];
}
特殊情况
- 无解：第一系数为0,常数非0
- 无数解：第一系数为0,常数为0
在开头交换后判一下系数最大值是否为0,有一个就记(flag=1)，回代时方程形式为(x=b),判b是否为0即可。
（注意先判无解，判完后若(flag=1),就是无穷解）
例题
T1 Luogu3389【模板】高斯消元法
见上
T2 hihoCoder1196 高斯消元法二
据说异或方程组高斯消元和普通高斯消元一模一样，就是把减改成异或。
若只有01建议用bitset压常数
bitset<50>a[50];
il void Gauss()
{
fp(i,1,30)
{
re int now=i;
fp(j,i+1,30)
if(a[j][i]>a[now][i]) now=j;
if(now^i) swap(a[now],a[i]);
fp(j,i+1,30) if(a[j][i]) a[j]^=a[i];
}
fq(i,30,1)
fq(j,i-1,1)
if(a[j][i]) a[j]^=a[i];
}
int main()
{
fp(i,1,5) scanf("%s",s[i]+1);
fp(i,1,5)
fp(j,1,6)
{
re int p=pos(i,j),l=pos(i,j-1),r=pos(i,j+1),u=pos(i-1,j),d=pos(i+1,j);
a[p][p]=1;a[p][31]=(s[i][j]-'0')^1;
if(i>1) a[p][u]=1;if(i<5) a[p][d]=1;if(j>1) a[p][l]=1;if(j<6) a[p][r]=1;
}
Gauss();
return 0;
}
T3Luogu2962 [USACO09NOV]灯Lights
给定一个无向联通图，开始点值全为0。你可以选择将一个点和与其相邻的点同时取反，问最少进行多少次这样的操作可以使点值全为1。
显然高斯消元解异或方程组。
而且由于方程内值只有(0/1)，可以开(bitset)使复杂度降到(O(frac{n^3}{64}))。
注意到这种方程解起来很容易得到形如(0x=0)的表达式（由于保证有解，不考虑(0x=b))。此时我们要枚举这个(x),才能继续往回代。
bitset<50>a[50];
int n,m,ans=1e9,tag[50];
il void Gauss()
{
fp(i,1,n)
{
re int now=i;
fp(j,i+1,n) if(a[j][i]>a[now][i]) now=j;
if(now^i) swap(a[now],a[i]);
if(!a[i][i]) continue;
fp(j,i+1,n) if(a[j][i]) a[j]^=a[i];
}
}
il void dfs(re int x,re int tot)
{
if(tot>=ans) return;
if(!x) {ans=tot;return;}
if(a[x][x])
{
re int t=a[x][n+1];
fp(i,x+1,n) if(a[x][i]) t^=tag[i];
tag[x]=t;
dfs(x-1,tot+t);
}
else
{
tag[x]=0;dfs(x-1,tot);
tag[x]=1;dfs(x-1,tot+1);
}
}
int main()
{
n=gi();m=gi();
fp(i,1,m)
{
re int u=gi(),v=gi();
a[u][v]=a[v][u]=1;
}
fp(i,1,n) a[i][i]=a[i][n+1]=1;
Gauss();
dfs(n,0);
printf("%d
",ans);
return 0;
}
T4Luogu2447 [SDOI2010]外星千足虫
给m个方程,最多n个未知数，解出方程并询问这最少需要要前几个方程组。
太板子了。
注意至少要(n)个方程。
il int Gauss(re int m)
{
fp(i,1,m) fp(j,1,n) a[i]=aa[i];
fp(i,1,n)
{
re int now=i;
fp(j,i+1,m) if(a[j][i]>a[now][i]) now=j;
if(now^i) swap(a[now],a[i]);
if(a[i][i]==0) return 0;
fp(j,i+1,m) if(a[j][i]) a[j]^=a[i];
}
fq(i,m,1)
fq(j,i-1,1)
if(a[j][i]) a[j]^=a[i];
if(mm==m) fp(i,1,n) ans[i]=a[i][n+1];
return 1;
}
int main()
{
n=gi();mm=m=gi();
fp(i,1,m)
{
scanf("%s",s+1);re int len=strlen(s+1);
fp(j,1,len)
if(s[j]=='1') a[i][j]=aa[i][j]=1;
a[i][n+1]=aa[i][n+1]=gi();
}
if(!Gauss(m)) {puts("Cannot Determine");return 0;}
re int l=n,r=m;
while(l<r)
{
re int mid=l+r>>1;
if(Gauss(mid)) r=mid;
else l=mid+1;
}
printf("%d
",l);
fp(i,1,n)
if(ans[i]&1) puts("?y7M#");
else puts("Earth");
return 0;
}
T5Luogu2973 [USACO10HOL]赶小猪
一个无向图，节点1有一个炸弹，在每个单位时间内，有p/q的概率在这个节点炸掉，有1-p/q的概率随机选择一条出去的路到其他的节点上。问最终炸弹在每个节点上爆炸的概率。
看到题就想起了[HNOI2013]游走。
好像思路方程都是一样的咦。
不过高斯消元解的应该是炸弹到某点的概率，炸弹爆炸概率最后再乘。
il void Gauss()
{
fp(i,1,n)//枚举列数
{
re int now=i;
fp(j,i+1,n)//枚举行
if(fabs(a[j][i])+eps>fabs(a[now][i])) now=j;
if(now^i) swap(a[i],a[now]);
fp(j,i+1,n)
fq(k,n+1,i)
a[j][k]-=a[i][k]*a[j][i]/a[i][i];
}
fq(i,n,1)
{
ans[i]=a[i][n+1];
fq(j,n,i+1) ans[i]-=ans[j]*a[i][j];
ans[i]/=a[i][i];
}
}
int main()
{
memset(h,-1,sizeof(h));
n=gi();m=gi();p=gi();q=gi();pi=(double)1.0*p/q;
fp(i,1,m)
{
re int u=gi(),v=gi();
add(u,v);
}
a[1][n+1]=1;
fp(u,1,n)
{
a[u][u]=1;
for(re int i=h[u];i+1;i=e[i].next)
{
re int v=e[i].to;
a[u][v]=-1.0*(1-pi)/in[v];
}
}
Gauss();
fp(i,1,n) printf("%.9lf
",ans[i]*pi);
return 0;
}
T6 bzoj3270博物馆
给定一个无向联通图，两个人分别在(a),(b)两点。在每一轮决策中，他们有(pi)的概率选择不动，有(1-pi)的概率等概率到任意一相邻房间。问在每个房间相遇的概率。
这题挺新奇的是吧。。。
要求每个房间的概率，我们就要求出每个状态的概率（即(a)在哪个点，(b)在哪个点）。
所以我们把一个状态的概率设为未知数。
方程？
(P_{a在i，b在j}-sum P_{ain i的邻点，bin i的邻点状态转移至此的概率}=0)
（懒得详写分类讨论)
注意不要从已相遇状态转移过来就成。
void Gauss()
{
fp(i,1,tot)
{
re int now=i;
fp(j,i+1,tot) if(fabs(a[j][i])+eps>fabs(a[now][i])) now=j;
if(now^i) swap(a[now],a[i]);
fp(j,i+1,tot)
fq(k,tot+1,i)
a[j][k]-=a[i][k]*a[j][i]/a[i][i];
}
fq(i,tot,1)
{
fq(j,tot,i+1) a[i][tot+1]-=a[j][tot+1]*a[i][j];
a[i][tot+1]/=a[i][i];
}
}
int main()
{
n=gi();m=gi();A=gi();B=gi();
fp(i,1,m)
{
re int u=gi(),v=gi();
d[u][v]=d[v][u]=1;++in[u];++in[v];
}
fp(i,1,n) scanf("%lf",&p[i]),d[i][i]=1;
fp(i,1,n) fp(j,1,n) id[i][j]=++tot;
a[id[A][B]][tot+1]=1;
fp(i,1,n)
fp(j,1,n)
{
a[id[i][j]][id[i][j]]=1;
fp(ii,1,n)
fp(jj,1,n)
if(d[ii][i]&&d[jj][j]&&ii!=jj)//ii!=jj
{
re double p1,p2;
if(i==ii) p1=p[ii];else p1=(1-p[ii])/in[ii];
if(j==jj) p2=p[jj];else p2=(1-p[jj])/in[jj];
a[id[i][j]][id[ii][jj]]+=-p1*p2;//+=???
}
}
Gauss();
fp(i,1,n) printf("%.6lf ",a[id[i][i]][tot+1]);
return 0;
}
T7 Luogu3164 [CQOI2014]和谐矩阵
我们称一个由0和1组成的矩阵是和谐的，当且仅当每个元素都有偶数个相邻的1。一个元素相邻的元素包括它本身，及他上下左右的4个元素（如果存在）。给定矩阵的行数和列数，请计算并输出一个和谐的矩阵。注意：所有元素为0的矩阵是不允许的。
直接上方程啊。
啥？输出全为(0)？把方程里的自由元赋为(1)即可.
int main()
{
n=gi();m=gi();tot=n*m;
fp(i,1,n)
fp(j,1,m)
{
re int p=pos(i,j),l=pos(i,j-1),r=pos(i,j+1),u=pos(i-1,j),d=pos(i+1,j);
a[p][p]=1;
if (i>1) a[p][u]=1;
if (j>1) a[p][l]=1;
if (i<n) a[p][d]=1;
if (j<m) a[p][r]=1;
}
fp(i,1,tot)
{
if (!a[i][i])
{
fp(j,i+1,tot)
if (a[j][i]) {swap(a[i],a[j]);break;}
}
fp(j,i+1,tot)
if (a[j][i]) a[j]^=a[i];
}
fq(i,tot,1)
{
ans[i]=a[i][tot+1];
fq(j,tot,i+1)
if (a[i][j]) ans[i]^=ans[j];
if (!a[i][i]) ans[i]=1;
}
for (int i=1;i<=n;++i,puts(""))
fp(j,1,m)
printf("%d ",ans[pos(i,j)]);
return 0;
}
T8 Luogu3265 [JLOI2015]装备购买
(n)个装备，每个装备(m)个属性，每个装备还有个价格。如果手里有的装备的每一项属性为它们分配系数（实数）后可以相加得到某件装备，则不必要买这件装备。求最多装备下的最小花费。
看到这道题可以想起线性基的性质：其中没有一个元素能被其它元素异或和求得。
那这题也类似？可以拿线性基的(insert)的方法套？
线性基里存(m)维，拿一个表达式进去时，若该表达式对应维系数为0就换维，如果该维没访问过就把这个式子整个塞进去，否则就把这个式子当前维消掉（高斯消元）。
于是，式子只有两个结果，一个是被塞进去，另一个是被消成蛋。
(竟然卡精度)
il int Insert(re int x)
{
fp(i,1,m)
{
if(fabs(a[x].s[i])<=eps) continue;
if(vis[i]) fq(j,m,i) a[x].s[j]-=p[i][j]*a[x].s[i]/p[i][i];
else {vis[i]=1;fp(j,i,m) p[i][j]=a[x].s[j];return 1;}
}
return 0;
}
int main()
{
n=gi();m=gi();
fp(i,1,n) fp(j,1,m) a[i].s[j]=gi();
fp(i,1,n) a[i].c=gi();
sort(a+1,a+1+n);
//fp(i,1,n) printf("%d %d %d
",a[i].s[1],a[i].s[2],a[i].s[3]);
fp(i,1,n)
if(Insert(i)) ++tot,ans+=a[i].c;
printf("%d %d
",tot,ans);
return 0;
}