【Learning】 多项式全家桶(加减乘除逆元ln exp)

欢迎订购KFC多项式全家桶

   

约定的记号

  
   对于一个多项式\(A(x)\),若其最高次系数不为零的项是\(x^k\),则该多项式的次数\(k\).
  
​   记为\(deg(A)=k\).
  
​   对于\(x\in(k,+ \infty)\),称\(x\)都为\(A(x)\)次数界. 但一般地,我们都使用\(k+1\)作为\(A(x)\)次数界
  
  
  

取模意义下的多项式运算

  
​   这里要消除一下疑惑,两个次数界为\(n\)的多项式在取模\(x^n\)意义下相乘,虽然NTT完了之后我们会得到一个次数界为\(2n-1\)的多项式,但是\(x^n....x^{2n-1}\)的系数是没有意义的。答案多项式就是\(x^0...x^{n-1}\)的系数描述的多项式。
  
​   求分式同理,\(\frac {F(x)}{G(x)}\pmod {x^n}\)的答案多项式,就是\(F(x)*G^{-1}(x)\)的前\(n\)个系数描述的多项式。
  
  
  

多项式加法减法乘法

  
​   略去,加减为\(\mathcal O(n)\)系数运算,乘法用NTT在\(\mathcal O(n\lg n)\)内计算。
  
  
  

多项式求逆

  
​   已知次数界为\(n\)的多项式\(A(x)\),求其在模\(x^n\)意义下的逆多项式\(B(x)\),使得

\[\begin{equation} \label{eqn1} A(x)B(x)\equiv1 \pmod {x^n} \end{equation} \]

  ​ 其中,\(B(x)\)的次数小于等于\(A(x)\)的次数
  
  
  
​   采用倍增思路向上倍增求解
  
​   当\(n=1\)时,\(A(x)\)\(B(x)\)仅有一个常数项,且\(B(x)\)的常数项为\(A(x)\)常数项的逆元。多项式有无逆元也取决于这个常数是否有逆元
  
​   假设已经求得\(A(x)\)在模\(x^{\lceil\frac n 2\rceil}\)意义下的逆元\(B'(x)\)

\[\begin{equation} A(x)B'(x)=1 \pmod{x^{\lceil\frac n 2\rceil}} \end{equation} \]

  ​ 而\((1)\)放在模\(x^{\lceil\frac n 2\rceil}\)意义下同样成立,有

\[\begin{equation} \label{eqn2} A(x)B(x)\equiv1 \pmod {x^{\lceil\frac n 2\rceil}} \end{equation} \]

  ​ 将\((3)-(2)\)得到

\[\begin{equation} B(x)-B'(x)\equiv 0 \pmod {x^{\lceil\frac n 2\rceil}} \end{equation} \]

​   平方得到

\[\begin{equation} B^2(x)-2B(x)B'(x)+B'^2(x)\equiv 0 \pmod {x^n} \end{equation} \]

  ​ 模数同时平方的原因是:原本多项式模\(x^{\lceil\frac n 2\rceil}\)为0,说明\(0\)~\(\lceil\frac n 2\rceil-1\)项的系数为0,平方后由于系数相乘,这些0系数会导致0~n-1项系数为0,也即模\(x^n\)为0
  
  ​ 两边同乘\(A(x)\),消去\(B(x)\)并移项,得到

\[B(x) \equiv 2B'(x) - A(x)B'^2(x) \pmod {x^n} \]

  ​ 至此可以递归倍增求解,伪代码如下,忽略清零、取模操作:
  

void polyInv(int *a,int *b,int n){ // a是要求逆元的多项式,b是在模x^n意义下的a的逆元
	if(n==1) 令b为一个次数为1,常数项为a[0]的逆元的多项式,返回;
  	polyInv(a,b,(n+1)>>1);
  	ntt_init(dega+degb+degb); //ntt长度
  	static int A[]=a,B[]=b; //临时数组,防止影响到传入的指针
  	ntt(A);
  	ntt(B);
  	for(int i=0;i<nttlen;i++)
    	A[i]=2*B[i]-A[i]*B[i]*B[i]; //点值计算
  	ntt(A,-1);
  	b=A;
}

  
  完整代码如下:
  

void polyInv(int *a,int *b,int n){
	if(n==1){
		b[0]=fmi(a[0],MOD-2);
		return;
	}
	int m=(n+1)>>1;
	polyInv(a,b,m);
	NTT::init(2*m+n);
	static int ta[S],tb[S];
	for(int i=0;i<n;i++) ta[i]=a[i],tb[i]=b[i];
	for(int i=n;i<NTT::n;i++) ta[i]=tb[i]=0;//NTT前对后半也要清空
	NTT::ntt(ta,0);
	NTT::ntt(tb,0);
	for(int i=0;i<NTT::n;i++) 
		ta[i]=(2LL*tb[i]-1LL*ta[i]*tb[i]%MOD*tb[i]%MOD)%MOD;
	NTT::ntt(ta,1);
	for(int i=0;i<n;i++) b[i]=ta[i];
}

  
​   Tips:调用整个递归之前需要清空\(b\)数组,否则\(b\)在赋值给\(B\)准备NTT时可能有处于\([(n+1)>>1,n)\)中的杂乱数字,将影响NTT。或者每次递归完成后都清空一次也可以。总之记得清空无用数据
  
​   时间复杂度\(O(n log n)\),不过我不知道怎么证。
  
  
  
  
  

多项式除法及取模

  
  ​ 给定多项式\(A(x)\)和多项式\(B(x)\),求两个多项式\(D(x)\)\(R(x)\),使得

\[\begin{equation} \label{1} A(x)=D(x)B(x)+R(x) \end{equation} \]

  ​ 其中,\(deg(A)>=deg(B)\) , \(deg(D)\leq deg(A)-deg(B)\)\(deg(R)<deg(B)\)
  
  
  
​   除法用奇妙变化解决.
  
​   引入一个操作:翻转操作. 对于一个次数为\(n\)的多项式\(A(x)\),定义

\[A^R(x)=x^nA(\frac 1 x) \]

​   会发现\(A^R(x)\)的系数相对于\(A(x)\)完全\(reverse\)了一下
  
​   将\((1)\)\(x\)换成\(\frac 1 x\),并两边同乘\(x^n\),记\(n=deg(A)\), \(m=deg(B)\), \(deg(D)=n-m\), \(deg(R)=m-1\).

\[\begin{aligned} x^nA(\frac1x)&=x^{n-m}D(\frac1x)x^mB(\frac1x)+x^{n-m+1}x^{m-1}R(\frac1x)\\ A^R(x)&=D^R(x)B^R(x)+x^{n-m+1}R^R(x) \end{aligned} \]

​   我们发现\(D^R(x)\)的次数仍然等于\(n-m\),如果把上式放在模\(x^{n-m+1}\)意义下,我们会发现\(D^R(x)\)不会受到任何影响,而\(x^{n-m+1}R^R(x)\)被完全模掉了!
  
​   于是

\[\begin{aligned} A^R(x)&\equiv D^R(x)B^R(x)\pmod{x^{n-m+1}}\\ D^R(x)&\equiv A^R(x){B^R}^{-1}(x)\pmod{x^{n-m+1}} \end{aligned} \]

​   直接对\(A\)系数翻转,\(B\)系数翻转,求\(B\)的逆,\(A\)\(B\)一乘,把结果的系数再次翻转,就是\(D\)了!
  

void polyDiv(int *a,int *b,int *d){// a和b对应上文的A和B,d是上文的D,是答案
  	if(deg(a)<deg(b)) 令d为一个0多项式,返回;
	reverse(a);
  	reverse(b);
  	static int invb[];
  	polyInv(b,invb,deg(b)+1); //b次数界是deg(b)+1
  	d=a*invb;
  	reverse(d);
}

  
  ​ 如果你还要求\(R(x)\),带回最初的式子直接算就好,多项式乘法用*直接代替了
  

void polyMod(int *a,int *b,int *r){
  	if(deg(a)<deg(b)) 令r为a,返回;
  	static int d[];
  	polyDiv(a,b,d);
  	r=a-d*b;
}

  
  
  
  
  

多项式求\(\ln\)

  
​   已知次数界为\(n\)的多项式\(f(x)\),并知关系\(f(x)=e^{g(x)}\)。求未知多项式\(g(x)=\ln f(x)\)
  
​   两边求导,得到

\[g'(x)=\frac {f'(x)}{f(x)} \]

  ​ 所以

\[g(x)=\int\frac{f'(x)}{f(x)} \]

​   先对\(f\)求导(\(\mathcal O(n)\)),再与\(f\)的逆多项式相乘(都是\(\mathcal O(n \lg n)\)),对结果再积分(\(\mathcal O(n)\)),得到\(g(x)\)
  

void polyDeri(int *a,int *b,int n){//计算一个次数界为n的多项式a的导数,放入b
	for(int i=0;i<n-1;i++) b[i]=1LL*(i+1)*a[i+1]%MOD;
	b[n-1]=0;
}
void polyInte(int *a,int *b,int n){//计算一个次数界为n的多项式a的积分,放入b
	for(int i=n-1;i>0;i--) 
		b[i]=1LL*a[i-1]*inv[i]%MOD;
	b[0]=0;	
}
void polyLn(int *f,int *g,int n){//整体计算于模x^n意义下进行
	static int t1[S],t2[S];
	polyDeri(f,t1,n);//求导,放入t1
	memset(t2,0,sizeof t2);
	polyInv(f,t2,n);//求逆,放入t2
	NTT::init(n*2-1);//接下来是t1与t2相乘
	for(int i=n;i<NTT::n;i++) t1[i]=t2[i]=0;
	NTT::ntt(t1,0);
	NTT::ntt(t2,0);
	for(int i=0;i<NTT::n;i++) t1[i]=1LL*t1[i]*t2[i]%MOD;
	NTT::ntt(t1,1);//此时t1只有0~n-1项的系数有意义,即为上述分式
	polyInte(t1,g,n);//积分,放入g
}

  
  
  
  
  

多项式牛顿迭代法

  
​   已知多项式\(g(x)\),求\(f(x)\)使得\(g(f(x))\equiv0\pmod{x^n}\)
  
​ 做法:倍增模数\(x^n\)
  
​   当\(n=1\)时,\(f(x)=1\)(令常数项为1)。
  
​   假设已知\(f_0(x)\)使得

\[g(f_0(x))\equiv0\pmod{x^{\lceil \frac n2\rceil}} \]

​   则构造\(f(x)\)

\[f(x)\equiv f_0(x)-\frac{g(f_0(x))}{g'(f_0(x))}\pmod{x^n} \]

​   可满足\(g(f(x))\equiv0\pmod{x^n}\)
  
  
  
  
  

多项式求exp

  
​   已知次数界为\(n\)的多项式\(f(x)\),并知关系\(g(x)\equiv e^{f(x)}\pmod{x^n}\)。求未知多项式\(g(x)\)
  
​   前置技能:多项式牛顿迭代法。
  
​   先来转化一下式子:

\[\begin{aligned} g(x)-e^{f(x)}&\equiv 0\pmod{x^n}\\ \ln g(x)-f(x)&\equiv 0 \pmod{x^n} \end{aligned} \]

​   其实就是要求解多项式函数

\[H(g(x))=\ln g(x)-f(x)\pmod{x^n} \]

​   的零点多项式。其中\(g(x)\)是变量,\(f(x)\)是常量多项式。
  
​   套用多项式牛顿迭代,假设已知多项式\(g_0(x)\)满足

\[\ln g_0(x)-f(x)\equiv 0\pmod{x^{\lceil \frac n2\rceil}} \]

  ​ 则构造\(g(x)\)

\[\begin{aligned} g(x)&\equiv g_0(x)-\frac{\ln g_0(x)-f(x)}{\frac 1 {g_0(x)}}\\ &\equiv g_0(x)*(1-\ln g_0(x)+f(x))\pmod{x^n} \end{aligned} \]

  ​ 可满足\(\ln g(x)-f(x)\equiv 0\pmod{x^{n}}\)
  
  ​ 向上倍增即可。
  

void polyExp(int *f,int *g,int n){
	if(n==1){
		g[0]=1;
		return;
	}
	static int t[S];
	polyExp(f,g,(n+1)>>1);//求得g0,放在g中
	polyLn(g,t,n);//求g0的ln,放在t中
	for(int i=0;i<n;i++) t[i]=(f[i]-t[i])%MOD;//构造后面一个括号的多项式,还是放在t中
	t[0]++;
	NTT::init(n*2-1);//g0和t相乘
	for(int i=n;i<NTT::n;i++) t[i]=g[i]=0;
	NTT::ntt(t,0);
	NTT::ntt(g,0);
	for(int i=0;i<NTT::n;i++) g[i]=1LL*t[i]*g[i]%MOD;
	NTT::ntt(g,1);
	for(int i=n;i<NTT::n;i++) g[i]=0;//清空不必要的东西
}

  
​   同求逆,最好在整体调用之前清空\(g\)数组,以免不必要的错误。
    

posted @ 2018-04-05 22:21  RogerDTZ  阅读(658)  评论(6编辑  收藏  举报