(模板)多项式乘法对任意数取模

最新推荐文章于 2021-07-23 14:06:04 发布

yamiedie_

最新推荐文章于 2021-07-23 14:06:04 发布

阅读量2.8k

点赞数

CC 4.0 BY-SA版权

分类专栏：模板

本文链接：https://blog.youkuaiyun.com/u013654696/article/details/49593915

模板专栏收录该内容

25 篇文章

订阅专栏

该博客介绍了如何利用NTT（快速数位变换）算法进行多项式乘法，并对结果进行模操作，特别是针对模数MOD=1000000007的情况。博客中详细解释了NTT的初始化、位翻转、蝶形运算等步骤，并展示了如何结合中国剩余定理（CRT）解决多个模数下的乘法问题。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

// 多项式乘法 系数对MOD=1000000007取模， 常数巨大，慎用
// 只要选的K个素数乘积大于MOD*MOD*N,理论上MOD可以任取。
#define MOD 1000000007
#define K 3

const int m[K] = {1004535809, 998244353, 104857601};
#define G 3


int qpow(int x, int k, int P) {
	int ret = 1;
	while(k) {
		if(k & 1) ret = 1LL * ret * x % P;
		k >>= 1;
		x = 1LL * x * x % P;
	}
	return ret;
}

struct _NTT {
	int wn[25], P;

	void init(int _P) {
		P = _P;
		for(int i = 1; i <= 21; ++i) {      
			int t = 1 << i;      
			wn[i] = qpow(G, (P - 1) / t, P);      
		}
    }
	void change(int *y, int len) {
		for(int i = 1, j = len / 2; i < len - 1; ++i) {      
			if(i < j) swap(y[i], y[j]);      
			int k = len / 2;      
			while(j >= k) {      
				j -= k;      
				k /= 2;      
			}      
			j += k;      
		} 
	}
	void NTT(int *y, int len, int on) {
		change(y, len);      
		int id = 0;      
      
		for(int h = 2; h <= len; h <<= 1) {      
			++id;      
			for(int j = 0; j < len; j += h) {      
				int w = 1;      
				for(int k = j; k < j + h / 2; ++k) {      
					int u = y[k];      
					int t = 1LL * y[k+h/2] * w % P;     
					y[k] = u + t;      
					if(y[k] >= P) y[k] -= P;      
					y[k+h/2] = u - t + P;      
					if(y[k+h/2] >= P) y[k+h/2] -= P;  
					w = 1LL * w * wn[id] % P;
				}      
			}      
		}      
		if(on == -1) {      
			for(int i = 1; i < len / 2; ++i) swap(y[i], y[len-i]);      
			int inv = qpow(len, P - 2, P);      
			for(int i = 0; i < len; ++i)   
				y[i] = 1LL * y[i] * inv % P;
		}      
	}
	void mul(int A[], int B[], int len) {
		NTT(A, len, 1);
		NTT(B, len, 1);
		for(int i = 0; i < len; ++i) A[i] = 1LL * A[i] * B[i] % P;
		NTT(A, len, -1);
	}
}ntt[K];

int tmp[N][K], t1[N], t2[N];
int r[K][K];

int CRT(int a[]) {
	int x[K];
	for(int i = 0; i < K; ++i) {
		x[i] = a[i];
		for(int j = 0; j < i; ++j) {
			int t = (x[i] - x[j]) % m[i];
			if(t < 0) t += m[i];
			x[i] = 1LL * t * r[j][i] % m[i];
		}
	}
	int mul = 1, ret = x[0] % MOD;
	for(int i = 1; i < K; ++i) {
		mul = 1LL * mul * m[i-1] % MOD;
		ret += 1LL * x[i] * mul % MOD;
		if(ret >= MOD) ret -= MOD;
	}
	return ret;
}

void mul(int A[], int B[], int len) {
	for(int id = 0; id < K; ++id) {

		for(int i = 0; i < len; ++i) {
			t1[i] = A[i];
			t2[i] = B[i];
		}
		ntt[id].mul(t1, t2, len);
		for(int i = 0; i < len; ++i) 
			tmp[i][id] = t1[i];
	}
	for(int i = 0; i < len; ++i) {
		A[i] = CRT(tmp[i]);

	}
}

void init() {
	for(int i = 0; i < K; ++i) {
		for(int j = 0; j < i; ++j) {
			r[j][i] = qpow(m[j], m[i] - 2, m[i]);
		}
	}
	for(int i = 0; i < K; ++i) {
		ntt[i].init(m[i]);
	}
}