多项式指数函数(多项式求exp)小结

多项式指数函数

顾名思义,对于多项式 A ( x ) A(x) A(x),要求出一个多项式 G ( x ) G(x) G(x) 满足 G ( x ) ≡ e A ( x ) ( m o d x n ) G(x)\equiv e^{A(x)} \pmod {x^n} G(x)eA(x)(modxn)

两边同时求 ln ⁡ \ln ln,得到 ln ⁡ ( G ( x ) ) ≡ A ( x ) \ln(G(x))\equiv A(x) ln(G(x))A(x),即 ln ⁡ ( G ( x ) ) − A ( x ) ≡ 0 \ln(G(x))-A(x)\equiv 0 ln(G(x))A(x)0

F ( G ( x ) ) = ln ⁡ ( G ( x ) ) − A ( x ) F(G(x))=\ln(G(x))-A(x) F(G(x))=ln(G(x))A(x),那么也就是要找到 F ( x ) F(x) F(x) 的零点,即找到多项式 G ( x ) G(x) G(x) 满足 F ( G ( x ) ) ≡ 0 ( m o d x n ) F(G(x))\equiv 0 \pmod {x^n} F(G(x))0(modxn),理所当然可以把牛顿迭代拿出来用。

假设已经求出了 G 0 ( x ) G_0(x) G0(x),满足 F ( G 0 ( x ) ) ≡ 0 ( m o d x ⌈ n 2 ⌉ ) F(G_0(x))\equiv 0\pmod {x^{\lceil \frac n 2 \rceil}} F(G0(x))0(modx2n)

根据牛顿迭代,有: G ( x ) = G 0 ( x ) − F ( G 0 ( x ) ) F ′ ( G 0 ( x ) ) G(x)=G_0(x)-\dfrac {F(G_0(x))} {F'(G_0(x))} G(x)=G0(x)F(G0(x))F(G0(x))

因为我们已经知道了 A ( x ) A(x) A(x),不妨将其看成常数项,则 F ′ ( x ) = 1 x F'(x)=\dfrac 1 x F(x)=x1,带入得到 G ( x ) = G 0 ( x ) − G 0 ( x ) F ( G 0 ( x ) ) = G 0 ( x ) ( 1 − ln ⁡ ( G 0 ( x ) ) + A ( x ) ) G(x)=G_0(x)-G_0(x)F(G_0(x))=G_0(x)(1-\ln(G_0(x))+A(x)) G(x)=G0(x)G0(x)F(G0(x))=G0(x)(1ln(G0(x))+A(x))

套上求 ln ⁡ \ln ln 板子就可以了。

代码如下:

#include <cstdio>
#include <cstring>
#include <cmath>
#include <vector>
#include <algorithm>
using namespace std;
#define maxn 600010
#define mod 998244353
#define bin(x) (1<<(x))

int n;
int inv[maxn];
int ksm(int x,int y){int re=1;for(;(y&1?re=1ll*re*x%mod:0),y;y>>=1,x=1ll*x*x%mod);return re;}
#define INV(x) ksm(x,mod-2)
struct NTT{
	vector<int>w[30];NTT(){
		inv[1]=1;for(int i=2;i<=maxn-10;i++)inv[i]=1ll*(mod-mod/i)*inv[mod%i]%mod;
		for(int i=1,wn;i<=21;i++){
			w[i].resize(bin(i));w[i][0]=1;wn=ksm(3,(mod-1)/bin(i));
			for(int j=1;j<bin(i-1);j++)w[i][j]=1ll*w[i][j-1]*wn%mod;
		}
	}
	int r[maxn],limit;void dft(int *f,int lg,int type=0)
	{
		limit=bin(lg);if(type)reverse(f+1,f+limit);
		for(int i=1;i<limit;i++){r[i]=(r[i>>1]>>1)|((i&1)<<(lg-1));if(i<r[i])swap(f[i],f[r[i]]);}
		for(int mid=1,Lg=1;mid<limit;mid<<=1,Lg++)for(int j=0;j<limit;j+=(mid<<1))for(int i=0;i<mid;i++)
		{int t=1ll*w[Lg][i]*f[j+i+mid]%mod;f[j+i+mid]=(f[j+i]-t+mod)%mod;f[j+i]=(f[j+i]+t)%mod;}
	}
}ntt;
int A[maxn],B[maxn],M;
struct POLY{
	vector<int> a;int len;void rs(int N){a.resize(len=N);}POLY(){rs(M);}
	int &operator [](int x){return a[x];}
	void dft(int *A_,int lg,int ln){for(int i=0;i<bin(lg);i++)A_[i]=(i<min(ln,len)?a[i]:0);ntt.dft(A_,lg);}
	void idft(int *A_,int lg,int ln){ntt.dft(A_,lg,1);rs(ln);for(int i=0;i<ln;i++)a[i]=1ll*A_[i]*inv[bin(lg)]%mod;}
	const POLY Mul(POLY b,int ln=M){
		int lg=ceil(log2(2*ln-1)),limit=bin(lg);dft(A,lg,ln);b.dft(B,lg,ln);
		for(int i=0;i<limit;i++)B[i]=1ll*A[i]*B[i]%mod;b.idft(B,lg,ln);return b;
	}
	const POLY operator *(const POLY b){return Mul(b);}
}F,G;
void getinv(POLY &f,POLY &g,int ln=M)
{
	if(ln==1){g.rs(1);g[0]=INV(f[0]);return;}getinv(f,g,(ln+1)>>1);
	int lg=ceil(log2(2*ln-1));f.dft(A,lg,ln);g.dft(B,lg,ln);
	for(int i=0;i<bin(lg);i++)B[i]=1ll*(2ll-1ll*A[i]*B[i]%mod+mod)%mod*B[i]%mod;g.idft(B,lg,ln);
}
POLY Jifen(POLY f){f.rs(f.len+1);for(int i=f.len-1;i>0;i--)f[i]=1ll*f[i-1]*inv[i]%mod;f[0]=0;return f;}
POLY Dao(POLY f){for(int i=0;i<f.len-1;i++)f[i]=1ll*f[i+1]*(i+1)%mod;f.rs(f.len-1);return f;}
POLY getln(POLY &f,int ln=M){POLY p;getinv(f,p,ln);return Jifen(Dao(f).Mul(p,ln-1));}
void getexp(POLY &f,POLY &g,int ln=M)
{
	if(ln==1){g.rs(1);g[0]=1;return;}getexp(f,g,(ln+1)>>1);
	POLY p=getln(g,ln);for(int i=0;i<ln;i++)p[i]=(f[i]-p[i]+mod)%mod;p[0]=(p[0]+1)%mod;
	int lg=ceil(log2(2*ln-1));p.dft(A,lg,ln);g.dft(B,lg,ln);
	for(int i=0;i<bin(lg);i++)B[i]=1ll*A[i]*B[i]%mod;g.idft(B,lg,ln);
}

int main()
{
	scanf("%d",&n);F.rs(n);M=n;
	for(int i=0;i<n;i++)scanf("%d",&F[i]);
	getexp(F,G,bin((int)ceil(log2(n))));
	for(int i=0;i<n;i++)printf("%d ",G[i]);
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值