问一个关于NTT的玄学问题
  • 板块学术版
  • 楼主2018heyuyang
  • 当前回复14
  • 已保存回复14
  • 发布时间2020/8/14 16:04
  • 上次更新2023/11/6 20:18:53
查看原帖
问一个关于NTT的玄学问题
120634
2018heyuyang楼主2020/8/14 16:04

多项式乘法逆(4238)和多项式开根(5205)这两题

#include<cstdio>
using namespace std;
inline int read()
{
	register char c=getchar();register int s=0;bool f=1;
	while(c<48){if(c=='-')f=0;c=getchar();}
	while(c>47)s=(s<<1)+(s<<3)+(c^48),c=getchar();
	return f?s:-s;
}
template<typename hyy>void swap(hyy &a,hyy &b){hyy c=a;a=b;b=c;}
const int N=1<<20;
const int mod=998244353;
const int inv2=499122177;
int suan(int a,int b)
{
	int ans=1;
	while(b)
	{
		if(b&1)ans=1ll*ans*a%mod;
		a=1ll*a*a%mod;b>>=1;
	}
	return ans;
}
const int G=3;//原根 G^(0~mod-1)为互不相同的整数 使用原根可以代替单位圆的"循环"性质 
const int invG=suan(G,mod-2);//这样就可以反着转圈圈啦 

int rev[N],qwq=-1;
void type(int n)
{
	if(qwq==n)return ;qwq=n;
	for(int i=0;i<n;i++)rev[i]=(rev[i>>1]>>1)|((i&1)?(n>>1):0);
}
void ntt(int *f,int op,int n)
{
	for(int i=0;i<n;i++)if(i<rev[i])swap(f[i],f[rev[i]]);
	int buf,tG,x,y;
	for(int p=2;p<=n;p<<=1)//p:合并的长度 
	{
		int len=p>>1;
		tG=suan(op==1?G:invG,(mod-1)/p);//乘上TG就是转1/n圈哦 
		//要注意是正向操作还是反向毒奶 
		for(int k=0;k<n;k+=p)//0~n-1
		{
			buf=1;//w(0,n)
			for(int l=k;l<k+len;l++)//扫描长度为(p>>1)的左区间 
			{
				x=f[l];y=1ll*buf*f[l+len]%mod;
				f[l+len]=(x-y+mod)%mod;
				f[l]=(x+y)%mod;
				buf=1ll*buf*tG%mod;
			}
		}
	}
	if(op!=1)
	{
		int invn=suan(n,mod-2);
		for(int i=0;i<n;i++)f[i]=1ll*f[i]*invn%mod;
	}
}
int A[N],B[N];
void cheng(int *a,int n,int *b,int m,int *c,int lim)//c为结果 
{//严格取a的0~n项 
	for(int i=0;i<n;i++)A[i]=a[i];
	for(int i=0;i<m;i++)B[i]=b[i];
	for(m+=n-2,n=1;n<=m;n<<=1);//n为2的幂次 
	
	type(n);
	ntt(A,1,n);ntt(B,1,n);//DFT
	for(int i=0;i<n;i++)c[i]=1ll*A[i]*B[i]%mod;//点乘 
	ntt(c,-1,n);//IDFT
	
	for(int i=0;i<n;i++)A[i]=B[i]=0;
	for(int i=lim;i<n;i++)c[i]=0;
}
int tmp1[N],tmp2[N];
void invF(int *a,int *c,int n)
{
	c[0]=suan(a[0],mod-2);
	int len,p;
	for(len=2;len<(n<<1);len<<=1)//len是目标长度哦 
	{
		p=len<<1;
		for(int i=0;i<len;i++)tmp1[i]=a[i],tmp2[i]=c[i];
		type(p);
		ntt(tmp1,1,p);ntt(tmp2,1,p);
		for(int i=0;i<p;i++)c[i]=(2ll-1ll*tmp1[i]*tmp2[i]%mod+mod)*tmp2[i]%mod;
		ntt(c,-1,p);
		for(int i=len;i<p;i++)c[i]=0;
	}
	for(int i=0;i<len;i++)tmp1[i]=tmp2[i]=0;
	for(int i=n;i<len;i++)c[i]=0;
}
/*
F(x)*G(x)=1(mod x^n)
F(x)*G'(x)=1(mod x^(n/2))
G(x)-G'(x)=0(mod x^(n/2))
G(x)^2-2G(x)*G'(x)+G'(x)^2=0(mod x^n)
G(x)-2G'(x)-F(x)G'(x)^2=0(mod x^n)
G(x)=2G'(x)-F(x)G'(x)^2(mod x^n)
G(0)=f[0]^(-1)

这里G(x)的目标长度为len 
但是F(x)G'(x)^2会让最高项到达x^(2*len-1) 
所以NTT的长度为p=len<<1 
*/
int tmp3[N],tmp4[N];
void sqrtF(int *a,int *c,int n)
{
	c[0]=1;
	int len,p;
	for(len=2;len<(n<<1);len<<=1)
	{
		p=len<<1;
		for(int i=0;i<len;i++)tmp3[i]=a[i];
		invF(c,tmp4,len);
		type(p);
		ntt(tmp3,1,p);ntt(tmp4,1,p);
		for(int i=0;i<p;i++)tmp3[i]=1ll*tmp3[i]*tmp4[i]%mod;
		ntt(tmp3,-1,p);
		for(int i=0;i<len;i++)c[i]=1ll*(c[i]+tmp3[i])*inv2%mod;
		for(int i=len;i<p;i++)c[i]=0;
	}
	for(int i=0;i<len;i++)tmp3[i]=tmp4[i]=0;
	for(int i=n;i<len;i++)c[i]=0;
}
/*
G(x)^2-F(x)=0(mod x^n)
G'(x)^2-F(x)=0(mod x^(n/2))
G(x)=G'(x)-(G'(x)^2-F(x))/(2G'(x))
G(x)=(F(x)+G'(x)^2)/(2G'(x))=inv2*(F(x)/G'(x)+G'(x))
*/
int a[N],b[N],c[N];
#include<cstring>
int main()
{
	int n;scanf("%d",&n);
	for(int i=0;i<n;i++)scanf("%d",&a[i]);c[1]=233333;
	sqrtF(a,c,n);
	for(int i=0;i<n;i++)printf("%d ",c[i]);
	return 0;
}

可以看到我读入完数据后面有一个c[1]=233333

但是ta不影响答案,why?

2020/8/14 16:04
加载中...