zoukankan      html  css  js  c++  java
  • 【Luogu4512】多项式除法(FFT)

    题面

    洛谷

    题解

    模板题。。。
    我直接蒯我写的东西。。。

    这个除法是带余除法,所以并不能直接求逆解决。
    要求的就是给定两个多项式(A(x),B(x)),其项数为(n,m)
    求解一个(n-m)项的多项式(C(x)),以及一个小于(n-m)项的多项式(R(x))
    满足:(A(x)=B(x)*C(x)+R(x))
    定义一个操作(R),其中(R)就是(Reverse)(A^R(x)=x^nA(frac{1}{x}))。这个操作说白点就是(A(x)[x^i])对应(A^R(x)[x^{n-i}]),也就是把所有的系数给翻转过来。
    然后推式子:

    [egin{aligned} A(x)&=B(x)*C(x)+R(x)\ A(frac{1}{x})&=B(frac{1}{x})*C(frac{1}{x})+R(frac{1}{x})\ x^nA(frac{1}{x})&=(x^m*B(frac{1}{x}))*(x^{n-m}*C(frac{1}{x}))+x^nR(frac{1}{x})\ A^R(x)&=B^R(x)*C^R(x)+R^R(x)*x^{n-m+1} end{aligned}]

    到了这里我们把等号换成同余,把整个式子在模(x^{n-m+1})意义下进行。

    [egin{aligned} A^R(x)&equiv B^R(x)*C^R(x)+R^R(x)*x^{n-m+1}\ &equiv B^R(x)*C^R(x) end{aligned}]

    所以在模意义下,我们可以利用多项式求逆求解(C^R(x)=frac{A^R(x)}{B^R(x)}mod x^{n-m+1})
    那么就有(R(x)=A(x)-B(x)*C(x))

    #include<iostream>
    #include<cstdio>
    #include<algorithm>
    using namespace std;
    #define MOD 998244353
    #define MAX 300000
    inline int read()
    {
    	int x=0;bool t=false;char ch=getchar();
    	while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
    	if(ch=='-')t=true,ch=getchar();
    	while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
    	return t?-x:x;
    }
    int fpow(int a,int b)
    {
    	int s=1;
    	while(b){if(b&1)s=1ll*s*a%MOD;a=1ll*a*a%MOD;b>>=1;}
    	return s;
    }
    int r[MAX],W[MAX];
    int N,n,m,NN;
    void NTT(int *P,int opt,int N)
    {
    	int l=0;for(int i=1;i<N;i<<=1)++l;
    	for(int i=0;i<N;++i)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
    	for(int i=0;i<N;++i)if(i<r[i])swap(P[i],P[r[i]]);
    	for(int i=1;i<N;i<<=1)
    	{
    		int w=fpow(3,(MOD-1)/(i<<1));W[0]=1;
    		for(int k=1;k<i;++k)W[k]=1ll*W[k-1]*w%MOD;
    		for(int p=i<<1,j=0;j<N;j+=p)
    			for(int k=0;k<i;++k)
    			{
    				int X=P[j+k],Y=1ll*W[k]*P[i+j+k]%MOD;
    				P[j+k]=(X+Y)%MOD;P[i+j+k]=(X-Y+MOD)%MOD;
    			}
    	}
    	if(opt==-1)
    	{
    		reverse(&P[1],&P[N]);
    		for(int i=0,inv=fpow(N,MOD-2);i<N;++i)P[i]=1ll*P[i]*inv%MOD;
    	}
    }
    int A[MAX],B[MAX];
    void Inv(int *a,int *b,int len) 
    {
    	if(len==1){b[0]=fpow(a[0],MOD-2);return;}
    	Inv(a,b,len>>1);
    	for(int i=0;i<len;++i)A[i]=a[i],B[i]=b[i];
    	NTT(A,1,len<<1);NTT(B,1,len<<1);
    	for(int i=0;i<len<<1;++i)A[i]=1ll*A[i]*B[i]%MOD*B[i]%MOD;
    	NTT(A,-1,len<<1);
    	for(int i=0;i<len;++i)b[i]=(b[i]+b[i])%MOD;
    	for(int i=0;i<len;++i)b[i]=(b[i]+MOD-A[i])%MOD;
    	for(int i=0;i<len<<1;++i)A[i]=B[i]=0;
    }
    int F[MAX],G[MAX],InvG[MAX];
    int Q[MAX],R[MAX];
    int main()
    {
    	n=read();m=read();
    	for(int i=0;i<=n;++i)F[i]=read();
    	for(int i=0;i<=m;++i)G[i]=read();
    	reverse(&F[0],&F[n+1]);reverse(&G[0],&G[m+1]);
    	for(NN=1;NN<=n-m+1;NN<<=1);Inv(G,InvG,NN);
    	for(N=1;N<=n+NN;N<<=1);
    	NTT(InvG,1,N);NTT(F,1,N);
    	for(int i=0;i<N;++i)Q[i]=1ll*F[i]*InvG[i]%MOD;
    	NTT(InvG,-1,N);NTT(F,-1,N);NTT(Q,-1,N);
    	for(int i=n-m+1;i<N;++i)Q[i]=0;
    	reverse(&F[0],&F[n+1]);reverse(&G[0],&G[m+1]);reverse(&Q[0],&Q[n-m+1]);
    	for(int i=0;i<=n-m;++i)printf("%d ",Q[i]);puts("");
    	for(N=1;N<=n;N<<=1);
    	NTT(G,1,N);NTT(Q,1,N);
    	for(int i=0;i<N;++i)G[i]=1ll*G[i]*Q[i]%MOD;
    	NTT(G,-1,N);
    	for(int i=0;i<N;++i)R[i]=(F[i]+MOD-G[i])%MOD;
    	for(int i=0;i<m;++i)printf("%d ",R[i]);puts("");
    	return 0;
    
    }
    
  • 相关阅读:
    资深项目经理推荐的几款免费/开源项目管理工具
    内网穿透工具frp简单使用教程
    10部全尺度欧美宫斗剧!献给不甘平淡的你
    Spring Boot后端+Vue前端+微信小程序,完整的开源解决方案!
    搭建Keepalived + Nginx + Tomcat的高可用负载均衡架构
    集成Activiti工作流的J2EE快速开发框架
    国内5大前端团队网站,你了解多少
    5 天 4000 star 的一个爆款开源项目
    「干货」常用的10个网络DOS命令,菜鸟学了变高手
    js自定义正则表达式
  • 原文地址:https://www.cnblogs.com/cjyyb/p/10129563.html
Copyright © 2011-2022 走看看