复习笔记——一些卷积

快速傅里叶变换

把多项式在系数表示与点值表示之间 (O(nlogn)) 转换

正变换是由系数 (Rightarrow) 点值
一个 (n) 次函数可用 (n+1) 个点值唯一确定,所以重点是选 (n+1) 个点快速求值

单位复数根

(x^n=1) 在复数意义下的解为 (n) 次复数根,显然有 (n) 个解。
联想高中学过的复数三角形式 (e^{ui}=cos(u)+isin(u))
((e^{ui})^n=1) ,则 (e^{uni}=cos(nu)+isin(nu)=1),则 (cos(nu)=1,sin(nu)=0 Rightarrow nu=2kpi,kin mathbb{Z})
于是 (u=frac{2kpi}{u},kin mathbb{Z}),所以 (e^{frac{2kpi i}{n}},kin mathbb{Z}) 都是解,当这无穷个数中其实只有 (n) 个不同的数(因为 (e^{ui}=e^{(u+2pi)i})
(omega_{n}=e^{frac{2pi i}{n}}) ,则 (omega_n^k=e^{frac{2kpi i}{n}}),显然原方程的解是 ({omega_n^k|k=0,1,dots,n-1 })

(omega_n) 为单位复数根,下面是它的一些性质:

  1. (omega_n^n=1)
  2. 消去引理:(omega_{dn}^{dk}=omega_n^k)
    证:显然 (omega_{dn}^{dk}=e^{frac{2dkpi i}{dn}}=e^{frac{2kpi i}{n}}=omega_n^k)
  3. 折半引理:(n)(n) 次单位复数根的平方的集合为 (n/2)(n/2) 次单位复数根集合
    证:由消去引理 (omega_n^{2k}=omega_{n/2}^k)
    对于 (k=0,1,dots,n/2-1) ,恰好构成的是 ({omega_{n/2}^k| k=0,1,dots,n/2-1})
    而对于 (k=k=n/2,n/2+1,dots,n-1),构成的是 ({omega_{n/2}^{k+n/2}| k=0,1,dots,n/2-1})
    注意到 (omega_{n}^{k+n}=omega_n^nomega_n^k=omega_n^k) ,所以这两个集合是一样的
  4. (omega_{2n}^{n+k}=omega_{2n}^{k}omega_{2n}^{n}=-omega_{2n}^{k})

选择 ({omega_n^k|k=0,1,dots,n-1 }) ((n) 为2的幂) 为要计算点值的点

DFT

分治计算。
原多项式 (A(x)=sumlimits_{i=0}^{n-1} a_ix^i)
(A^{[0]}(x)=sumlimits_{i=0}^{n/2-1} a_{2i}x^{i},A^{[1]}(x)=sumlimits_{i=0}^{n/2-1} a_{2i+1}x^{i})
相当于把奇偶项分开了,则 (A(x)=A^{[0]}(x^2)+xA^{[1]}(x^2))

由单位复数根的性质(公式中 (k<n/2) ) :

[egin{equation*} egin{aligned} A(omega_n^k)&=A^{[0]}(omega_{n/2}^k)+omega_n^kA^{[1]}(omega_{n/2}^k) \ A(omega_n^{k+n/2})&=A^{[0]}(omega_{n/2}^{k+n/2})+omega_n^{k+n/2}A^{[1]}(omega_{n/2}^{k+n/2}) \ &=A^{[0]}(omega_{n/2}^k)-omega_n^k A^{[1]}(omega_{n/2}^k) end{aligned} end{equation*} ]

于是可以递归再合并了
当然代码里写的是非递归

IDFT

理论不太会……但代码与 (DFT) 几乎一样……

注意 (l) 的大小很玄妙……一定一定别开小了……
(r[])for(int i=1;i<l;i++) r[i]=(r[i>>1]>>1)|((i&1)*(l>>1));
(IDFT) 后的最终结果是 (int(A[i].r/l+0.5))

const double pi = 3.1415926535897932384626433832795;
struct c{ //复数
    double r,i;
    c() { r=i=0.0; }
    c(double x,double y) { r=x; i=y; }
    c operator + (const c &b) { return c(r+b.r,i+b.i); }
    c operator += (const c &b) { return *this=*this+b; }
    c operator - (const c &b) { return c(r-b.r,i-b.i); }
    c operator -= (const c &b) { return *this=*this-b; }
    c operator * (const c &b) { return c(r*b.r-i*b.i,r*b.i+i*b.r); }
    c operator *= (const c &b) { return *this=*this*b; }
}x[N];
int l,r[N];
void fft(c *A,int ty){
    for(int i=0;i<l;i++) x[r[i]]=A[i];
    for(int i=0;i<l;i++) A[i]=x[i];
    for(int i=2;i<=l;i<<=1){
        c wn(cos(pi*2/i),ty*sin(pi*2/i));
        for(int j=0;j<l;j+=i){
            c w(1,0);
            for(int k=j;k<j+i/2;k++){
                c t=A[k+i/2]*w;
                A[k+i/2]=A[k]-t;
                A[k]+=t;
                w*=wn;
            }
        }
    }
}

快速数论变换

理论跟 (FFT) 差不多,区别是取模,将单位复数根换成了原根的幂
常见原根 (g) :998244353,1004535809的原根都是3
然而具体的我也不太清楚 (qwq)

int l,r[N],X[N];
void ntt(int *A,int ty){
	for(int i=0;i<l;i++) X[r[i]]=A[i];
	for(int i=0;i<l;i++) A[i]=X[i];
	for(int i=2;i<=l;i<<=1){
		int wn=Pow_mod(3,(P-1)/i);
		if(ty==-1) wn=Pow_mod(wn,P-2);
		for(int j=0;j<l;j+=i){
			int w=1;
			for(int k=j;k<j+i/2;k++){
				int t=1ll*A[k+i/2]*w%P;
				A[k+i/2]=Minus(A[k],t);
				A[k]=Plus(A[k],t);
				w=1ll*w*wn%P;
			}
		}
	}
	if(ty==1) return;
	int Inv=Pow_mod(l,P-2);
	for(int i=0;i<l;i++) A[i]=1ll*A[i]*Inv%P;
}


快速沃尔什变换

也是类似系数表示与点值表示之间的转换

与(&)卷积

(C[i]=sumlimits_{j&k=i} A[j] imes B[k])
需要设计出数组 (A'=FWT(A)) 满足 (C'[i]=A'[i] imes B'[i])

考虑这样设计 (A'[i]=sumlimits_{j&i=i} A[j])
那么 (C'[i]=sumlimits_{j&i=i} sumlimits_{k&l=j} A[k] imes B[l]=sumlimits_{j&i=i} A[j] imes sumlimits_{j&i=i} B[j]=A'[i] imes B'[i])
验证成功。

考虑如何求 (A'[]) (即 (FWT) 过程)
(A0[]) 表示下标最高位为0的,(A1[]) 表示下标最高位为1的
(FWT(A)=FWT(A0,A1)=(FWT(A0)+FWT(A1),FWT(A1)))
(IFWT(A)=IFWT(A0,A1)=(IFWT(A0)-IFWT(A1),IFWT(A1)))

void fwt_and(int *A,int ty){
	for(int i=2;i<=n;i<<=1)
		for(int j=0;j<n;j+=i)
			for(int k=j;k<j+i/2;k++){
				if(ty==1) A[k]=Plus(A[k],A[k+i/2]);
                                else A[k]=Minus(A[k],A[k+i/2]);
                        } 
}

或(|)卷积

(C[i]=sumlimits_{j|k=i} A[j] imes B[k])
需要设计出数组 (A'=FWT(A)) 满足 (C'[i]=A'[i] imes B'[i])

考虑这样设计 (A'[i]=sumlimits_{j|i=i} A[j])
那么 (C'[i]=sumlimits_{j|i=i} sumlimits_{k|l=j} A[k] imes B[l]=sumlimits_{j|i=i} A[j] imes sumlimits_{j|i=i} B[j]=A'[i] imes B'[i])
验证成功。

考虑如何求 (A'[]) (即 (FWT) 过程)
(A0[]) 表示下标最高位为0的,(A1[]) 表示下标最高位为1的
(FWT(A)=FWT(A0,A1)=(FWT(A0),FWT(A0)+FWT(A1)))
(IFWT(A)=IFWT(A0,A1)=(IFWT(A0),IFWT(A1)-IFWT(A0)))

void fwt_or(int *A,int ty){
	for(int i=2;i<=n;i<<=1)
		for(int j=0;j<n;j+=i)
			for(int k=j;k<j+i/2;k++){
				if(ty==1) A[k+i/2]=Plus(A[k+i/2],A[k]);
                                else A[k+i/2]=Minus(A[k+i/2],A[k]);
                        } 
}

异或(^)卷积

这个有点难……
原理:令 (i&j) 中 1 的个数的奇偶性为 “(i)(j) 的奇偶性”,那么 “(i)(k) 的奇偶性” 异或 “(k)(j) 的奇偶性” 等于 “(i hat{} j)(k) 的奇偶性”

(A'[i]=sumlimits_{i&j奇偶性为0} A[j]-sumlimits_{i&j奇偶性为1} A[j])

[egin{equation*} egin{aligned} C'[i]&=sumlimits_{i&j奇偶性为0} sumlimits_{khat{}l=j} A[k]B[l]-sumlimits_{i&j奇偶性为1} sumlimits_{khat{}l=j} A[k]B[l] \ &=sumlimits_{i&j奇偶性为0,i&k奇偶性为0} A[j]B[k]+sumlimits_{i&j奇偶性为1,i&k奇偶性为1} A[j]B[k]\ & -sumlimits_{i&j奇偶性为1,i&k奇偶性为0} A[j]B[k]-sumlimits_{i&j奇偶性为0,i&k奇偶性为1} A[j]B[k] \ &=(sumlimits_{i&j奇偶性为0} A[j]-sumlimits_{i&j奇偶性为1} A[j]) imes (sumlimits_{i&j奇偶性为0} B[j]-sumlimits_{i&j奇偶性为1} B[j] )\ &= A'[i] imes B'[i] end{aligned} end{equation*} ]

验证成功。

考虑如何求 (A'[]) (即 (FWT) 过程)
(A0[]) 表示下标最高位为0的,(A1[]) 表示下标最高位为1的
(FWT(A)=FWT(A0,A1)=(FWT(A0)+FWT(A1),FWT(A0)-FWT(A1)))
(IFWT(A)=IFWT(A0,A1)=(frac{IFWT(A0)+IFWT(A1)}{2},frac{IFWT(A0)-IFWT(A1)}{2}))

void fwt_xor(int *A,int ty){
	for(int i=2;i<=n;i<<=1)
		for(int j=0;j<n;j+=i)
			for(int k=j;k<j+i/2;k++){
				int t=A[k+i/2];
				A[k+i/2]=Minus(A[k],t);
				A[k]=Plus(A[k],t);
				if(ty==-1) 
					A[k+i/2]=1ll*A[k+i/2]*Q%P,A[k]=1ll*A[k]*Q%P;//Q=1/2(mod P)
			}
}

子集卷积

(C[i]=sumlimits_{j|k=i,j&k=0} A[j] imes B[k])

如果没有 (j&k=0) 的限制,就是或卷积。
有了这个限制就要求在 (j|k=i) 的基础上满足 (|i|=|j|+|k|) ,绝对值表示1的个数
(|i|=|j|+|k|) 让我们想到了普通的卷积。

于是我们卷积套卷积,外层普通卷积,内层或卷积。
具体地,令 (f[i][j]=[|j|==i] imes A[j],g[i][j]=[|j|==i] imes B[j])
那么每个 (f[i],g[i]) 都可看成是多项式
(h[i]=sumlimits_{j+k=i} f[j] oplus g[k]) ,两个多项式 (aoplus b) 表示或卷积
最终 (C[i]=h[|i|][i])

下面的是洛谷模板题代码

#include<cstdio>
#include<iostream>
#include<algorithm>

#define P 1000000009

using namespace std;

int read(){
	int x=0;
	char ch=getchar();
	while(!isdigit(ch)) ch=getchar();
	while(isdigit(ch)) x=x*10+ch-'0',ch=getchar();
	return x;
}

const int N = 1<<20;

int n,m,a[N],b[N];

int Plus(int x,int y) { return x+y<P?x+y:x+y-P; }
int Minus(int x,int y) { return x>=y?x-y:x-y+P; }
void fwt(int *A,int ty){
	for(int i=2;i<=m;i<<=1)
		for(int j=0;j<m;j+=i)
			for(int k=j;k<j+i/2;k++){
				if(ty==1) A[k+i/2]=Plus(A[k+i/2],A[k]);
				else A[k+i/2]=Minus(A[k+i/2],A[k]); 
			}
}
int f[25][N],g[25][N],num[N];
int h[25][N];

int main()
{
	n=read();
	m=(1<<n);
	for(int i=0;i<m;i++) a[i]=read();
	for(int i=0;i<m;i++) b[i]=read();
	
	for(int i=0;i<m;i++){
		int t=0,x=i;
		for(;x;x>>=1) t+=(x&1);
		num[i]=t;
		f[t][i]=a[i];
		g[t][i]=b[i];
	}
	for(int i=0;i<=n;i++){ //pre fwt
		fwt(f[i],1); fwt(g[i],1);
	}
	
	for(int i=0;i<=n;i++)
		for(int j=0;j<=i;j++)
			for(int k=0;k<m;k++) 
				h[i][k]=Plus(h[i][k],1ll*f[j][k]*g[i-j][k]%P); //don't fwt back
	
	for(int i=0;i<=n;i++) fwt(h[i],-1);
	
	for(int i=0;i<m;i++) printf("%d ",h[num[i]][i]);
	
	return 0;
}

狄利克雷卷积

(h=f*g),即 (h[i]=sumlimits_{d|i} f[d] imes g[frac{i}{d}])

在杜教筛中会用到
主要是有两个重要的式子,这里先声明几个函数:
(mu(n)) 莫比乌斯函数,(phi(n)) 欧拉函数
(epsilon(n)=[n=1]) 元函数,(I(n)=1) 恒等函数,(id(n)=n) 单位函数

重要式子来了:

  1. (mu * I=epsilon) 这个式子确保了莫比乌斯反演的正确性
  2. (phi * I=id) 看着很美妙但不会证…

原文地址:https://www.cnblogs.com/lindalee/p/13456437.html