EOJ 1074多项式展开

1074多项式展开 题解

author : uniecho1

思路

其实一上来我是心存幻想,希望能不用多项式乘法切掉的…

然后简单算了一下,如果用线性求逆元预处理组合数,上二项式展开,再怎么样都是 O ( n 2 ) O(n^2) O(n2)级别…

所以肯定是要推出卷积形式上多项式乘法了

(顺便一提,这个题涉及到对 998244353 998244353 998244353取模了,必然不会是快速傅里叶变换,而是快速数论变换)

如果 A = 0 A=0 A=0的话就直接输出就行了,我们只考虑 A ! = 0 A!=0 A!=0的情况
f ( x ) = Σ i = 0 n a i ( x + A ) i = Σ i = 0 n a i Σ j = 0 i C i j x j A i − j = Σ i = 0 n x i Σ j = i n a j C j i A j − i = Σ i = 0 n x i Σ j = i n a j j ! i ! ( j − i ) ! A j − i f(x)=\Sigma_{i=0}^{n}a_i(x+A)^i=\Sigma_{i=0}^{n}a_i\Sigma_{j=0}^iC_{i}^{j}x^jA^{i-j}=\Sigma_{i=0}^nx^i\Sigma_{j=i}^na_jC_{j}^iA^{j-i}=\Sigma_{i=0}^nx^i\Sigma_{j=i}^na_j\frac{j!}{i!(j-i)!}A^{j-i} f(x)=Σi=0nai(x+A)i=Σi=0naiΣj=0iCijxjAij=Σi=0nxiΣj=inajCjiAji=Σi=0nxiΣj=inaji!(ji)!j!Aji

f ( x ) = Σ i = 0 n b i x i f(x)=\Sigma_{i=0}^nb_ix^i f(x)=Σi=0nbixi
对比一下系数,得到
b i = Σ j = i n j ! i ! ( j − i ) ! A j − i a j b_i=\Sigma_{j=i}^n\frac{j!}{i!(j-i)!}A^{j-i}a_j bi=Σj=ini!(ji)!j!Ajiaj
移项,得到
i ! ∗ b i ∗ A i = Σ j = i n j ! ( j − i ) ! A j = Σ j = 0 n − i ( n − j ) ! ( n − j − i ) ! A j a j i!*b_i*A^i=\Sigma_{j=i}^n\frac{j!}{(j-i)!}A^j=\Sigma_{j=0}^{n-i}\frac{(n-j)!}{(n-j-i)!}A^ja_j i!biAi=Σj=in(ji)!j!Aj=Σj=0ni(nji)!(nj)!Ajaj
注意到 ( n − j − i ) ! (n-j-i)! (nji)!的首项会随着 i i i的变化而不断变化,不能很好的处理,尝试用 n − i n-i ni来替换 i i i,得到
( n − i ) ! ∗ b n − i ∗ A n − i = Σ j = 0 i ( n − j ) ! A n − j ( i − j ) ! a n − j (n-i)!*b_{n-i}*A^{n-i}=\Sigma_{j=0}^i\frac{(n-j)!A^{n-j}}{(i-j)!}a_{n-j} (ni)!bniAni=Σj=0i(ij)!(nj)!Anjanj
我们依次令
p i = ( n − i ) ! b n − i A n − i p_i=(n-i)!b_{n-i}A^{n-i} pi=(ni)!bniAni
q j = ( n − j ) ! A n − j a n − j q_{j}=(n-j)!A^{n-j}a_{n-j} qj=(nj)!Anjanj
r j = 1 j ! r_j=\frac{1}{j!} rj=j!1
那么可以得到卷积形式
p i = Σ j = 0 i q j r i − j p_i=\Sigma_{j=0}^{i}q_{j}r_{i-j} pi=Σj=0iqjrij
先把这玩意儿放在一边,我们来看看多项式乘法:

给定 n n n次多项式 B ( x ) = b n x n + b n − 1 x n − 1 + . . . + b 1 x + b 0 B(x)=b_nx^n+b_{n-1}x^{n-1}+...+b_1x+b_0 B(x)=bnxn+bn1xn1+...+b1x+b0 m m m次多项式 C ( x ) = c n x n + c n − 1 x n − 1 + . . . + c 1 x + c 0 C(x)=c_nx^n+c_{n-1}x^{n-1}+...+c_1x+c0 C(x)=cnxn+cn1xn1+...+c1x+c0,求它们的乘积 n + m n+m n+m次多项式 A ( x ) = B ( x ) ∗ C ( x ) = a n + m x n + m + a n + m − 1 x n + m − 1 + . . . + a 1 x + a 0 A(x)=B(x)*C(x)=a_{n+m}x^{n+m}+a_{n+m-1}x^{n+m-1}+...+a_1x+a_0 A(x)=B(x)C(x)=an+mxn+m+an+m1xn+m1+...+a1x+a0

简单思考一下,发现有
a i = Σ j = 0 i b j c i − j a_i=\Sigma_{j=0}^{i}b_jc_{i-j} ai=Σj=0ibjcij
和上面的卷积形式一模一样,所以实际上我们现在要做的就是求一个多项式乘法

考虑一下做法:如果直接展开,显然还是 O ( n 2 ) O(n^2) O(n2),鉴定为纯纯的超时

这个时候就要引入新的概念:多项式的点值表达形式。具体来说,对于一个 n n n次多项式 B ( x ) B(x) B(x),如果给定了其在 n n n个不同 x x x处的值,那么我们可以唯一确定 B ( x ) B(x) B(x).回想一下线性代数的知识,这个结论还是比较容易得到的

那么用点值表达形式的多项式怎么求多项式乘法 A ( x ) = B ( x ) ∗ C ( x ) A(x)=B(x)*C(x) A(x)=B(x)C(x)呢?如果两个多项式的取值的 x x x集合相同,那么我们只需要将对应于同一 x x x B ( x ) 与 C ( x ) B(x)与C(x) B(x)C(x)相乘,然后就能得到 A ( x ) A(x) A(x) x x x处的值,进一步得到 A ( x ) A(x) A(x)的点值表达式

现在我们可以 O ( n ) O(n) O(n)进行点值表达式的求积了,然而我们还需要考虑如何将 B ( x ) B(x) B(x) C ( x ) C(x) C(x)从一般表达转换为点值表达,以及如何把 A ( x ) A(x) A(x)从点值表达转换为多项式表达

如果随便取 n n n个不同的 x x x进行多项式求值(注意这里其实已经 n 2 n^2 n2复杂度了),然后再拉格朗日插值反解系数,时间复杂度依旧是 O ( n 2 ) O(n^2) O(n2),鉴定为寄

所以我们要通过精心选取 x x x来进行求值,来减少多项式求值与插值的复杂度

所谓精心选取,放在 F F T FFT FFT里,是选取复数单位根;放在 N T T NTT NTT里,是选取原根。这俩都有一些比较神奇的性质

限于篇幅,就题论题,只说 N T T NTT NTT

  • 原根的定义:设 m m m为正整数, a a a是整数,若 a m o d    m a \mod m amodm等于 ϕ ( m ) \phi(m) ϕ(m)(欧拉函数),那么称 a a a为模 m m m的一个原根

  • 阶的定义:若对于存在一个最小正整数 d d d使得 a d m o d    m = 1 a^d \mod m=1 admodm=1,那么称 d d d a m o d    m a \mod m amodm的阶

998244353 998244353 998244353是一个素数(可以记一下,这个数经常出现),它的一个原根为 3 3 3

基于上述定义,以及现在模数是一个素数的情况(注意这一点很重要,任意模数 N T T NTT NTT会非常棘手),原根会有以下有趣的性质:

  • 对于任意 0 < i < p 0 < i < p 0<i<p g i m o d    p g^i\mod p gimodp互不相同(一句话证明:如果中间出现相同的了,说明会有更小的 d d d使得 a d m o d    m = 1 a^d \mod m=1 admodm=1,与原根定义矛盾)

有什么用呢?假如我们现在要求一个 n = 2 l n=2^l n=2l次多项式 A ( x ) A(x) A(x) n n n个不同点的值:

  1. g n = g p − 1 n g_n=g^{\frac{p-1}{n}} gn=gnp1(ps. 998244353 = 2 23 ∗ 119 + 1 998244353=2^{23}*119+1 998244353=223119+1,只要 n n n没有大于 2 23 2^{23} 223就不会有问题,对于这个题而言,显然够了)
  2. A ( x ) A(x) A(x)按项数的奇偶性拆分成两部分 A ( x ) = Σ i = 0 2 l a i x i = Σ i = 0 2 l − 1 a 2 i x 2 i + x Σ i = 0 2 l − 1 a 2 i + 1 x 2 i A(x)=\Sigma_{i=0}^{2^l}a_ix^i=\Sigma_{i=0}^{2^{l-1}}a_{2i}x^{2i}+x\Sigma_{i=0}^{2^{l-1}}a_{2i+1}x^{2i} A(x)=Σi=02laixi=Σi=02l1a2ix2i+xΣi=02l1a2i+1x2i
  3. 我们将 g n k g_n^{k} gnk g n k + n 2 g_n^{k+\frac{n}{2}} gnk+2n分别带入 A ( x ) A(x) A(x)求值(其中 0 ≤ k < n / 2 0\le k<n/2 0k<n/2),我们会得到

f o r m u l a 1 : A ( g n k ) = Σ i = 0 2 l − 1 a 2 i g n 2 k i + g n k Σ i = 0 2 l − 1 a 2 i + 1 g n 2 k i formula1: A(g_n^k)=\Sigma_{i=0}^{2^{l-1}}a_{2i}g_n^{2ki}+g_n^k\Sigma_{i=0}^{2^{l-1}}a_{2i+1}g_n^{2ki} formula1:A(gnk)=Σi=02l1a2ign2ki+gnkΣi=02l1a2i+1gn2ki

f o r m u l a 2 : A ( g n k + n 2 ) = Σ i = 0 2 l − 1 a 2 i g n 2 k i + n i + g n k + n 2 Σ i = 0 2 l − 1 a 2 i + 1 g n 2 k i + n i formula2:A(g_n^{k+\frac{n}{2}})=\Sigma_{i=0}^{2^{l-1}}a_{2i}g_n^{2ki+ni}+g_n^{k+\frac{n}{2}}\Sigma_{i=0}^{2^{l-1}}a_{2i+1}g_n^{2ki+ni} formula2:A(gnk+2n)=Σi=02l1a2ign2ki+ni+gnk+2nΣi=02l1a2i+1gn2ki+ni

引理:

  • g n n = 1 g_n^{n}=1 gnn=1(一句话证明: g n n = g p − 1 = g ϕ ( p ) = 1 g_n^n=g^{p-1}=g^{\phi(p)}=1 gnn=gp1=gϕ(p)=1)(欧拉定理)
  • g n k + n = g n k g_n^{k+n}=g_n^{k} gnk+n=gnk(一句话证明: g n k + n = g n k ∗ g n n = g n k g_n^{k+n}=g_n^k*g_n^n=g_n^k gnk+n=gnkgnn=gnk)
  • g n n 2 = − 1 g_n^{\frac{n}{2}}=-1 gn2n=1(一句话证明: g n n = ( g n n 2 ) 2 g_n^n=(g_n^{\frac{n}{2}})^2 gnn=(gn2n)2,所以 g n n 2 = ± 1 g_n^{\frac{n}{2}}=\pm1 gn2n=±1,根据之前互不相同的结论 g n n 2 = − 1 g_n^{\frac{n}{2}}=-1 gn2n=1)
  • g n k + n 2 = − g n k g_n^{k+\frac{n}{2}}=-g_n^{k} gnk+2n=gnk(一句话证明: g n k + n 2 = g n k ∗ g n n 2 = − g n k g_n^{k+\frac{n}{2}}=g_n^k*g_n^\frac{n}{2}=-g_n^k gnk+2n=gnkgn2n=gnk)

所以我们可以把第二个式子换成
f o r m u l a 3 : A ( g n k + n 2 ) = Σ i = 0 2 l − 1 a 2 i g n 2 k i − g n k Σ i = 0 2 l − 1 a 2 i + 1 g n 2 k i formula3:A(g_n^{k+\frac{n}{2}})=\Sigma_{i=0}^{2^{l-1}}a_{2i}g_n^{2ki}-g_n^{k}\Sigma_{i=0}^{2^{l-1}}a_{2i+1}g_n^{2ki} formula3:A(gnk+2n)=Σi=02l1a2ign2kignkΣi=02l1a2i+1gn2ki
f o r m u l a 1 formula1 formula1 f o r m u l a 3 formula3 formula3对比一下,我们会发现,唯一的区别,就是中间那个正负号

这说明,我们只需要计算一半的值,剩下的一半可以直接改正负号得到

再来观察一下 f o r m u l a 3 formula3 formula3:

因为 g n 2 k i = g p − 1 n ∗ 2 k i = g p − 1 2 l − 1 ∗ k i g_n^{2ki}=g^{\frac{p-1}{n}*2ki}=g^{\frac{p-1}{2^{l-1}}*ki} gn2ki=gnp12ki=g2l1p1ki

我们令 n ′ = 2 l − 1 n'=2^{l-1} n=2l1,那么我们要计算的前一个部分会变成
Σ i = 0 2 l − 1 a 2 i g n 2 k i = Σ i = 0 n ′ a 2 i g n ′ k i = A e v e n ( g n ′ k i ) \Sigma_{i=0}^{2^{l-1}}a_{2i}g_n^{2ki}=\Sigma_{i=0}^{n'}a_{2i}g_{n'}^{ki}=A^{even}(g_{n'}^{ki}) Σi=02l1a2ign2ki=Σi=0na2ignki=Aeven(gnki)
本质是对 A A A中偶数项进行 g n ′ k i g_{n'}^{ki} gnki(其中 0 ≤ k < n ′ 0 \le k<n' 0k<n)处的表达式求值,我们依旧可以使用 g n ′ k i g_{n'}^{ki} gnki g n ′ k i + n ′ 2 g_{n'}^{ki+\frac{n'}{2}} gnki+2n来成对计算,从而降低复杂度;后一个部分同理

也就是说,这个计算量减半的过程是可递归的

简单看一下复杂度: f ( n ) = 2 ∗ f ( n 2 ) + O ( n ) f(n)=2*f(\frac{n}{2})+O(n) f(n)=2f(2n)+O(n)根据主定理复杂度为 O ( n l o g n ) O(nlogn) O(nlogn)

所以,现在我们做到了以 O ( n l o g n ) O(nlogn) O(nlogn)的时间代价从多项式的一般表达式得到点值表达式

再来看一下怎么从点值表达式求得一般表达式

假如我们现在有 A ( x ) A(x) A(x) n {n} n g n k i g_{n}^{ki} gnki出的值

f o r m u l a 1 + f o r m u l a 3 2 \frac{formula1+formula3}{2} 2formula1+formula3可以直接得到 A e v e n ( x ) A^{even}(x) Aeven(x) n ′ n' n g n ′ k i g_{n'}^{ki} gnki处的值

f o r m u l a 1 − f o r m u l a 3 2 ∗ ( g n k i ) − 1 \frac{formula1-formula3}{2}*(g_n^{ki})^{-1} 2formula1formula3(gnki)1可以直接得到 A o d d ( x ) A^{odd}(x) Aodd(x) n ′ n' n g n ′ k i g_{n'}^{ki} gnki处的值

依旧是一个递归问题, e v e n even even o d d odd odd不断划分最后就得到某个系数的值了

显然这个过程的时间复杂度也是 O ( n l o g n ) O(nlog n) O(nlogn)

综上所述,我们的流程变成了
B 和 C 求 点 值 表 达 O ( n l o g n ) − > B 和 C 在 每 一 个 点 求 积 得 到 A 的 点 值 表 达 式 O ( n ) − > 把 A 转 换 为 一 般 表 达 式 O ( n l o g n ) B和C求点值表达O(nlogn)->B和C在每一个点求积得到A的点值表达式O(n)->把A转换为一般表达式O(nlog n) BCO(nlogn)>BCAO(n)>AO(nlogn)
时间复杂度成功降下来了

不过现在是递归进行快速数论变换,常数巨大,有可能卡不过去,所以我们还需要把它改成迭代版 N T T NT T NTT

[外链图片转存失败,源站可能有防盗链机制,建议将图片保存下来直接上传(img-PYAkOxDz-1653984449189)(C:\Users\uniecho1\AppData\Roaming\Typora\typora-user-images\image-20220529104427727.png)]

观察一下每次按奇偶性分类,最后的结果就是二进制表达翻转排序(其实考虑每一次都是按最后一位是 0 / 1 0/1 0/1分类的话挺好证明的,我就不证了)

所以我们可以直接先预处理排序结果,然后自底向上迭代

好了假如现在我们通过 N T T NTT NTT i ! b i A i i!b_iA^i i!biAi求出来了,那么接下来就只需要对 i ! i! i! A i A^i Ai求逆元再乘回去就行了

A i A^i Ai很简单, i ! i! i!也就是个线性求逆元的事儿, i − 1 = − ⌊ ( P / i ) ⌋ ( P m o d    i ) − 1 i^{-1}=-\lfloor (P/i)\rfloor(P \mod i)^{-1} i1=(P/i)(Pmodi)1,然后累乘就行了

完结撒花,附上完整代码:

#include<bits/stdc++.h>
#define int long long
using namespace std;

const int MAXN=1e6+5;
const int mod=998244353;
const int G=3;
const int invG=332748118;

int quickPow(int x,int k) { //经典快速幂
	if(!k)
		return 1;
	else if(k&1)
		return quickPow(x*x%mod,k>>1)*x%mod;
	else
		return quickPow(x*x%mod,k>>1);
}

struct convolve {
	int lim,len;
	int R[MAXN];
	void NTT(int* A,int flag) {
		//A:实现NTT变换的数组指针
		//flag=1:一般表达式->点值表达式
		//flag=-1:点值表达式->一般表达式
		for(int i=0; i<lim; i++)
			if(i<R[i])
				swap(A[i],A[R[i]]);//显然R[R[i]]=i,所以交换一下就完事儿了
		for(int i=1; i<lim; i<<=1) {
			int n=i<<1;//一次性处理的区间的长度
			int gn=quickPow(flag>0?G:invG,(mod-1)/n);
            //一般表达式转点值表达式的时候获得n次原根
            //点值表达式转一般表达式的获得其逆元
			for(int j=0; j<lim; j+=n) { //枚举处理哪个区间
				int g=1;
				for(int k=0; k<i; k++) {
                    //对于上一层的区间而言,A[j+k]和A[j+i+k]使用的都是g_{n/2}^{k}
                    //所以放在这里直接取用就行
					int x=A[j+k],y=g*A[j+i+k]%mod;
					A[j+k]=(x+y)%mod;
					A[j+i+k]=(x-y+mod)%mod;
                    //注意,对于点值表达式转一般表达式的情况,我们没有进行除2的操作,这个会在最后补上
					g=g*gn%mod;//获得下一个g
				}
			}
		}
	}
	void conv(int N,int* A,int* B,int* C) {
		lim=1;
		while(lim<=2*N)//添加0项使得多项式项数为2^lim,方便后面均分操作
			lim<<=1,++len;
		for(int i=0; i<lim; i++)
			R[i]=(R[i>>1]>>1)|((i&1)<<(len-1));//计算i的二进制翻转得到的值,本质上是一个动态规划
		NTT(B,1);
		NTT(C,1);
		for(int i=0; i<lim; i++)
			A[i]=B[i]*C[i]%mod;//直接求积
		NTT(A,-1);
		int invlim=quickPow(lim,mod-2);
		for(int i=0; i<=N; i++)
			A[i]=A[i]*invlim%mod;//把NTT里面没有除去的2给除了
	}
} solution;

int n,A;
int a[MAXN],p[MAXN],q[MAXN],r[MAXN],fac[MAXN],inv[MAXN],invfac[MAXN];

signed main() {
	//freopen("in.txt","r",stdin);
	ios::sync_with_stdio(false);
	cin>>n;
	for(int i=0; i<=n; i++)
		cin>>a[i];
	cin>>A;
	if(!A)
		for(int i=0; i<=n; i++)//特判一下A=0的情况,直接输出就行
			cout<<a[i]<<" ";
	else {
		fac[0]=fac[1]=1;
		for(int i=2; i<=n; i++)
			fac[i]=fac[i-1]*i%mod;//阶乘
		inv[0]=inv[1]=1;
		for(int i=2; i<=n; i++)
			inv[i]=(mod-mod/i)*inv[mod%i]%mod;//线性求逆元
		invfac[0]=invfac[1]=1;
		for(int i=2; i<=n; i++)
			invfac[i]=invfac[i-1]*inv[i]%mod;//阶乘逆元
		for(int i=0; i<=n; i++)
			q[i]=fac[i]*quickPow(A,i)%mod*a[i]%mod;
		reverse(q,q+n+1);//获得q
		for(int i=0; i<=n; i++)
			r[i]=invfac[i];//获得r
		solution.conv(n,p,q,r);
		reverse(p,p+n+1);//转成原来的顺序
		int invA=quickPow(A,mod-2);
		for(int i=0; i<=n; i++)
			p[i]=p[i]*invfac[i]%mod*quickPow(invA,i)%mod;//把多余的项给抵消掉
		for(int i=0; i<=n; i++)
			cout<<p[i]<<" ";//输出答案
        return 0;
	}
	return 0;
}
  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值