线性基与高斯消元
高斯消元
作用：解方程
这里采用P2455 线性方程组作为模板题，因为这题数据比较强
模仿普通的加减消元，可以想出一种比较靠谱的方法
- 枚举主元系数最大的一个方程
- 用这个方程来消去其他所有方程的主元项
很朴素的方法，可以将最后的矩阵消为一个对角线形式，归纳一下，假设当前消到了第 \(i\) 个未知数，则其他所有方程的之前的未知数的系数都为 0，此时，我们拿着这个方程去消其他的方程，并不会让前面的未知数的系数秽土转生，那就可以顺利地消掉当前的未知数，同理，选中方程的其他未知数项会被之后的方程给消掉
要求枚举的方程系数最大一是要其本身可以去消其他的方程，二是如果这个方程的系数过小，会产生很大的精度误差
好的，如果这个方程有解的话，这种朴素的算法已经可以解决问题了，但接下来我们要解决两个十分棘手的问题：如何判断多解和无解的情况，
比较朴素的想法是这样的，如果最后的方程出现了\(0x=1\)的情况，就是无解，出现了 \(0x=0\) 的情况，就是多解，代码：
fp(i, 1, n) {
int id = i;
fp(j, i, m) if (fabsl(a[j][i]) > fabsl(a[id][i]))
id = j;
if (fabsl(a[id][i]) > eps)
fp(j, 1, n + 1) swap(a[id][j], a[i][j]);
else continue ;
fd(j, n + 1, i) a[i][j] /= a[i][i];
fp(j, 1, m)
if (i ^ j)
fd(k, n + 1, i) a[j][k] -= a[i][k] * a[j][i];
}
fp(i, 1, n)
if (fabsl(a[i][i]) < eps && fabs(a[i][n + 1]) > eps) {
cout << "0" << endl;
return 0;
}
fp(i, 1, n)
if (fabsl(a[i][i]) < eps && fabs(a[i][n + 1]) < eps) {
cout << "-1" << endl;
return 0;
}
fp(i, 1, n)
printf("x%d=%.2lf\n", i, (double)a[i][n + 1]);
于是你信誓旦旦的交了一发上去，发现：
为啥呢，我们先看一个比较明显的 Hack 数据：
3 3
0 0 1 1
0 1 0 1
0 1 0 1
这应该是一个多解的数据，但是上面的代码会出现无解
为啥呢，因为我们找不到第一个未知数系数不为 0 的方程，然后就去找第二个了，导致第一个方程和第三个是同一个，所以导致了错误判断
解决方法也很简单，如果一个方程的所有系数都为 0，而且答案不为 0，那就无解了
判断多解也比较简单了，只要一个未知数的系数为 0 ，就可以了，我们可以记录一下现在解了几个方程，在判断完无解之后再判断一下就可以了
如果方程总数大于了未知数的个数，那也不要紧，利用上面的两个方法，也可以完成判断
void guass() {
int t=1;
fp(i,1,n) {
int id=t;
fp(j,t+1,m) if(fabsl(a[j][i])>fabsl(a[id][i])) id=j;
if(fabsl(a[id][i])>eps) swap(a[id],a[t]);
else continue ;
fd(j,n+1,i) a[t][j]/=a[t][i];
fp(j,1,m)
if(t^j)
fd(k,n+1,i) a[j][k]-=a[t][k]*a[j][i];
t++;
}
fp(i,1,m) {
bool f=1;
for(int j=1;j<=n;j++) f&=fabs(a[i][j])<eps;
if(f&&fabs(a[i][n+1])>eps){
cout << "-1" << endl;
return;
}
}
if(t<=n) {
cout << "0" << endl;
return;
}
fp(i,1,n) printf("x%d=%.2Lf\n",i,a[i][n+1]);
}
注意，有些时候题目的恶心程度会使得你设的 eps 寄掉，这个时候就需要因题而论了
高斯消元的主要作用就是解方程，一般来说，如果 DP 是有后效性的，就不好从一个状态切入，这个时候就需要列出方程然后扔进高斯消元里面
P4035 球形空间产生器
打住，我知道你想 SA
这题就有点 只要把两个方程对对碰一下就可以消掉二次方项了SB了
void gauss() {
int t=1;
fp(i,1,n) {
int id=t;
fp(j,t+1,n) if(fabsl(a[j][i])>fabsl(a[id][i])) id=j;
if(fabsl(a[id][i])>eps) swap(a[id],a[t]);
else continue ;
fd(j,n+1,i) a[t][j]/=a[t][i];
fp(j,1,n)
if(t^j)
fd(k,n+1,i) a[j][k]-=a[t][k]*a[j][i];
t++;
}
fp(i,1,n) {
bool f=1;
for(int j=1;j<=n;j++) f&=fabs(a[i][j])<eps;
if(f&&fabs(a[i][n+1])>eps){
cout << "No solutions" << endl;
return;
}
}
if(t<=n){
cout << "Many solutions" << endl;return;}
fp(i,1,n) printf("%.3lf ",(double)a[i][n+1]);
}
signed main(){
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
cin>> n;
fp(i,1,n+1)
fp(j,1,n)
cin >> b[i][j];
fp(i,2,n+1)
fp(j,1,n){
a[i-1][n+1]+= b[i][j]*b[i][j]-b[i-1][j]*b[i-1][j];
a[i-1][j]=(b[i][j]-b[i-1][j])*2;
}
gauss();
return 0;
}
再来一题
UVA11542 Square
给出 \(n\) 个整数，从中选出 \(1\) 个或者多个，使得选出的整数乘积是完全平方数。一共有多少种选法？比如\(4,6,10,15\) 有 \(3\)种选法 $ {4},{6,10, 15},{4,6,10,15}$
这题是一个典型的异或方程组问题，如果我们把每一个数看成一个二进制数，则类似下图
明显，这个问题可以使用高斯消元来解决，但是由于作者实在是太弱了，无法理解高斯消元的真正原理，所以这里来讲解具有更快复杂度的线性基解法
回到题目本身，异或方程组的实际意义就是给你一堆二进制数，有多少种选择，使得这些选出来的数异或为一个二进制数
线性基的解法并不大显然，可以结合这一篇题解来理解一下
首先，对于上图这样的方程，我们改变 \(b\) ，确定 \(a\) 则不会改变解的数量（假设有解），这是因为如果我们把所有的数加入到线性基里面，则有一些数插入失败，这些数会与其中的一些数异或为 0，则我们可以选出任意组 0 ，然后再与线性基中可以的数完成组合，达到目的
答案就是二的没有成功插入的数的次方，若要求的 b 全为 0 ，就要减一
const int N=1e5+10;
int n;
int a[N];
struct linear{
bitset<100> a[100];
int cnt;
void insert(bitset<100> x){
fd(i,70,1){
if(!x[i]) continue ;
if(!a[i].any()){
a[i]=x;
cnt++;
return ;
}
x^=a[i];
}
}
}lin;
signed main(){
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
n=rd();
fp(i,1,n){
bitset<100> b;
a[i]=rd();
b.reset();
for(int j=2;j*j<=a[i];++j)
while(a[i]%j==0)
b[j]=b[j]^1,a[i]/=j;
if(a[i]!=1) b[a[i]]=b[a[i]]^1;
lin.insert(b);
}
ll ans=1;
int lim=n-lin.cnt;
fp(i,1,lim)
ans=(ans<<1)%mod;
ans=(ans-1+mod)%mod;
cout << ans << endl;
return 0;
}
最后来一道高斯消元的偏难题
P6125 有趣的游戏
小阳阳发明了一个有趣的游戏：有 \(n\) 个玩家，每个玩家都有一个长度为 \(l\) 的字母序列，任何两个玩家的字母序列不同。共有 \(m\) 种不同的字母，所有的字母序列都由这 \(m\) 种字母构成。为了方便，我们取大写字母的前 \(m\) 个字母。
例如 \(m=3,l=4,\texttt{ABAA}\) 和 \(\texttt{CBCA}\) 是两个合法的字母序列。
现在由小阳阳来操控一台神奇的机器，每个时刻机器会随机产生一个字母，其中第 \(i\) 种字母随机出来的概率为 \(\dfrac{p_i}{q_i}\) ，显然 \(\sum \limits_{k=1}^m \dfrac{p_i}{q_i}=1\)。
这样 \(T\) 个时刻后机器会产生一个长度为 \(T\) 的字母序列。
如果某个时刻某个玩家发现自己的字母序列在机器产生的字母序列中出现了，“出现”的定义是玩家的字母序列是机器产生的字母序列中连续的一段，那么我们称这个玩家获胜，游戏结束。
现在小阳阳感兴趣的一个问题是，每个玩家分别有多大的概率能获得这场游戏的胜利呢？
这题一眼先建出 ACAM ，然后列出概率的方程，然后呢？
你会发现你列出了一堆方程，但是没有除了 0 ，之外的常数，这意味着这个方程一定有无穷多组解
怎么办呢，想一想能不能凑出来一个参数，这里不知道为啥就想到了期望，对于一个终止节点的经过的期望而言，这个期望值就是其概率，而根节点则第一次一定会被经过一次，所以就有了需要的常数
复杂度\(O(n^3m^3)\)
const int N=1e3+10;
int n,m;
ld a[N][N];
ld q[N],p[N];
int ed[N];
void gauss(int n) {
int t=1;
fp(i,1,n) {
int id=t;
fp(j,t+1,n) if(fabsl(a[j][i])>fabsl(a[id][i])) id=j;
if(fabsl(a[id][i])!=0) swap(a[id],a[t]);
else continue ;
fd(j,n+1,i) a[t][j]/=a[t][i];
fp(j,1,n)
if(t^j)
fd(k,n+1,i) a[j][k]-=a[t][k]*a[j][i];
t++;
}
}
struct ac{
int son[26][N<<4],fail[N<<4],idx=1,cnt[N<<4];
void insert(string s,int id){
int n=s.length(),now=1;
s=" "+s;
fp(i,1,n){
if(!son[s[i]-'A'][now]) son[s[i]-'A'][now]=++idx;
now=son[s[i]-'A'][now];
}
ed[id]=now;
cnt[now]++;
}
void build(){
queue<int> q;
fp(i,0,m-1)
if(son[i][1]) q.push(son[i][1]),fail[son[i][1]]=1;
else son[i][1]=1;
while(!q.empty()){
int now=q.front();
q.pop();
fp(i,0,m-1)
if(!son[i][now]) son[i][now]=max(1,son[i][fail[now]]);
else fail[son[i][now]]=max(son[i][fail[now]],1),q.push(son[i][now]);
}
}
void add_edge(){
fp(i,1,idx)
fp(j,0,m-1)
if(!cnt[i])
a[son[j][i]][i]-=p[j+1]*1.0/q[j+1];
}
}acam;
signed main(){
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
cin >> n >> m >> m;
fp(i,1,m) cin >> p[i]>> q[i];
string s;
fp(i,1,n){
cin >> s;
acam.insert(s,i);
}
acam.build();
acam.add_edge();
fp(i,1,acam.idx){
a[i][i]+=1;
if(i==1) a[i][acam.idx+1]=1;
}
gauss(acam.idx);
fp(i,1,n) printf("%.2Lf\n",a[ed[i]][acam.idx+1]);
return 0;
}
当然了，对于这一类问题，其实有更好的方法来解决，实际上，我们可以只对每一个结尾设一个概率变量就可以了，这些概率变量的和一定是 1 ，而且根据这篇题解，我们可以确定一些偏序关系，从而使方程有顺序的构造
线性基
顾名思义，就是找一组数，让其能异或出来所有给出的数
模板分两种，一种是传统的，一种是压缩式的，其中传统的可以维护所有的操作，压缩式则无法实现可删除线性基
先上代码
传统：
struct linear{
ll a[101],cnt;
void insert(ll x){
fd(i,60,0)
if(x>>i&1){
if(!a[i]) return a[i]=x,++cnt,void();
x^=a[i];
}
}
}lin;
压缩：
struct liear{
ll a[101],cnt;
void insert(ll x){
fp(i,1,cnt) x=min(x,a[i]^x);
fp(i,1,cnt) a[i]=min(a[i],a[i]^x);
if(x) a[++cnt]=x;
}
}lin;
原理不太通透，先不写
基本操作：
异或最大值
传统：
ll query(ll x){
fd(i,60,0) x=max(x,x^a[i]);
return x;
}
压缩：
ll quiery(ll x){
fp(i,1,cnt) x=max(x,x^a[i]);
return x;
}
特别的，对于压缩式线性基而言，输出最大值只需要将基中所有的元素都异或起来就可以了
k小异或和
传统：
ll query(int k){
ll cut=cnt,x=0;
fd(i,60,0){
if(!a[i]) continue ;
--cur;
if( (k>>cur&1) ^ (x>>i&1)) x^=a[i];
}
return x;
}
这里解释一下判断语句，如果当前这一位要为 1 了，但是选出来的数还不为 1 ，那就要异或上，反之亦然
压缩：
sort(a+1,a+cnt+1);
ll query(int k){
ll x=0;
fp(i,1,cnt) if(k>>(i-1)&1) x^=a[i];
return x;
}
k在异或集合中的排名
经典，不难发现，对于一个二进制位，其如果在线性基中出现，则出现次数一定是一半，所以如果这一为上为1，前面比其小的位就可以随便选了，而对于没有插入成功的数而言，其能够使异或值翻一番，单独计算就可以了
压缩：
define high(x) 63-__builtin_clzll(x)
sort(a+1,a+cnt+1);
int query(int k){
ll ans=0;
fp(i,1,cnt) if(x>>high(x)&1) ans+=(1ll<<i-1);
fp(i,1,n-cnt) ans<<=1;
return ans ;
}
普通的不会，先不写了
合并
暴力插就可以了
可删除（离线）
对于每一个数做一个被删的时间戳，插入的时候，如果这一位的时间戳小于要加入的数的时间戳，很明显就可以用当前的数来代替这个数，但也不能直接把这个数给丢了，需要继续计算这个数对下面的贡献，询问的时候就把时间小于等于当前的数给加入贡献就可以了
struct linear{
ll a[50],t[50];
void insert(ll x,ll tim){
fd(i,31,0)
if((x>>i)&1ll){
if(t[i]<tim)
swap(t[i],tim),swap(a[i],x);
if(!tim) break;
x^=a[i];
}
}
void era(int tim){
fd(i,31,0) if(t[i]==tim) t[i]=a[i]=0;
}
ll query(){
ll ans=0;
fd(i,31,0) if((ans^a[i])>ans) ans^=a[i];
return ans;
}
}lin;
真正的线性基
就是高斯消元套皮，只要把普通线性基里的x^=a[i] 换成 x=x-a[i]*(x[i]/a[i][i]) 就可以了，当然这里是重载符号的，x和a[i]是向量
为了不掉精度，这里使用了有理数取模
int n,m;
vector<int> operator - (vector<int> x,vector<int> y){
vector<int> z;
fp(i,0,(int)x.size()-1)
z.push_back((mod+x[i]-y[i])%mod);
return z;
}
vector<int> operator * (vector<int> x,int y){
vector<int> z;
for(int xx:x)
z.push_back(xx*y%mod);
return z;
}
bool operator == (vector<int> x,vector<int> y){
fp(i,0,m-1)
if(x[i]^y[i]) return false;
return true;
}
ll qmi(ll x,ll p){
ll ans=1;
for(;p;p>>=1){if(p&1)(ans*=x)%=mod;(x*=x)%=mod;}
return ans;
}
struct dlinear{
vector<int> a[N];
int insert(vector <int> v){
fp(i,0,m-1){
if(!v[i]) continue ;
if(a[i].empty()){
a[i]=v;
return 1;
}
v=v-a[i]*(v[i]*qmi(a[i][i],mod-2)%mod);
}
return 0;
}
}lin;
struct node{
vector<int> v;
int c;
bool operator < (node z){
return c<z.c;
}
}a[N];
signed main(){
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
n=rd(),m=rd();
fp(i,1,n)
fp(j,1,m)
a[i].v.push_back(rd());
fp(i,1,n) a[i].c=rd();
sort(a+1,a+n+1);
ll sum=0,cnt=0;
fp(i,1,n){
int x=lin.insert(a[i].v);
sum+=x*a[i].c;
cnt+=x;
}
cout <<cnt <<' ' <<sum << endl;
return 0;
}
区间异或线性基
先给一个小问题，x和y能异或出来的东西与x和x^y能异或出来的东西一样吗
答案很明显，是一样的，这个问题还可以进行归纳，多个数一起也是一样的
问题来了，既然这两个异或出来的数是一样的，那就可以对着异或差分来建出线性基
然后就发现区间修改就没了，这题也做完了
PS:树状数组维护异或的时候请在
l处和r+1处修改
const int N=1e5+10;
int n,q;
int a[N];
struct linear{
int a[40],cnt;
void insert(int x){
fp(i,1,cnt) x=min(x,a[i]^x);
fp(i,1,cnt) a[i]=min(a[i],a[i]^x);
if(x) a[++cnt]=x;
}
linear operator + (linear x){
fp(i,1,cnt) if(a[i])x.insert(a[i]);
return x;
}
void clear(){
met(a,0),cnt=0;
}
linear(){
met(a,0),cnt=0;
}
int query(int x){
fp(i,1,cnt) x=max(x,x^a[i]);
return x;
}
};
struct segment{
linear data[N<<2];
int a[N<<2];
void modify(int now,int l,int r,int wh,int x){
if(l==r){
data[now].clear();
data[now].insert(a[now]^=x);
return ;
}
int mid((l+r)>>1);
if(wh<=mid) modify(now<<1,l,mid,wh,x);
else modify(now<<1|1,mid+1,r,wh,x);
data[now]=data[now<<1]+data[now<<1|1];
}
linear query(int now,int l,int r,int ql,int qr){
if(ql<=l&&r<=qr) return data[now];
int mid((l+r)>>1);
if(qr<=mid) return query(now<<1,l,mid,ql,qr);
else if(ql>mid) return query(now<<1|1,mid+1,r,ql,qr);
else return query(now<<1,l,mid,ql,qr)+query(now<<1|1,mid+1,r,ql,qr);
}
}seg;
struct Bit{
int c[N];
int query(int wh){
int res=0;
for(;wh;wh-=wh&-wh) res^=c[wh];
return res;
}
void modify(int l,int r,int x){
for(r++;r<=n;r+=r&-r) c[r]^=x;
for(;l<=n;l+=l&-l) c[l]^=x;
}
void add(int wh,int x){
for(;wh<=n;wh+=wh&-wh) c[wh]^=x;
}
}bit;
signed main(){
ios::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
n=rd(),q=rd();
fp(i,1,n) a[i]=rd();
fp(i,1,n) bit.modify(i,i,a[i]),seg.modify(1,1,n,i,a[i]^a[i-1]);
while(q--){
int op=rd(),l=rd(),r=rd(),x=rd();
if(op==1){
bit.modify(l,r,x);
seg.modify(1,1,n,l,x);
if(r+1<=n) seg.modify(1,1,n,r+1,x);
}
else {
linear y;
y.insert(bit.query(l));
if(r>l) y=y+seg.query(1,1,n,l+1,r);
cout << y.query(x) << endl;
}
}
return 0;
}