Luogu5205 【模板】多项式开根（NTT+多项式求逆）

时间：2023-03-08 17:49:21

　　https://www.cnblogs.com/HocRiser/p/8207295.html 安利！

　　写NTT把i<<=1写成了i<<=2，又调了一年。发现我的日常就是数组开小调调调，变量名写错调调调，反向判if调调调，退役吧。

#include<iostream>

#include<cstdio>

#include<cmath>

#include<cstdlib>

#include<cstring>

#include<algorithm>

using namespace std;

#define ll long long

#define P 998244353

#define N 550000

char getc(){char c=getchar();while ((c<'A'||c>'Z')&&(c<'a'||c>'z')&&(c<'0'||c>'9')) c=getchar();return c;}

int gcd(int n,int m){return m==0?n:gcd(m,n%m);}

int read()

{

	int x=0,f=1;char c=getchar();

	while (c<'0'||c>'9') {if (c=='-') f=-1;c=getchar();}

	while (c>='0'&&c<='9') x=(x<<1)+(x<<3)+(c^48),c=getchar();

	return x*f;

}

int n,a[N],r[N],b[N],A[N],B[N],f[N],g[N],t;

int ksm(int a,int k)

{

	int s=1;

	for (;k;k>>=1,a=1ll*a*a%P) if (k&1) s=1ll*s*a%P;

	return s;

}

int inv(int a){return ksm(a,P-2);}

void DFT(int n,int *a,int g)

{

	for (int i=0;i<n;i++) r[i]=(r[i>>1]>>1)|(i&1)*(n>>1);

	for (int i=0;i<n;i++) if (i<r[i]) swap(a[i],a[r[i]]);

	for (int i=2;i<=n;i<<=1)

	{

		int wn=ksm(g,(P-1)/i);

		for (int j=0;j<n;j+=i)

		{

			int w=1;

			for (int k=j;k<j+(i>>1);k++,w=1ll*w*wn%P)

			{

				int x=a[k],y=1ll*w*a[k+(i>>1)]%P;

				a[k]=(x+y)%P,a[k+(i>>1)]=(x-y+P)%P;

			}

		}

	}

}

void IDFT(int *a,int n)

{

	DFT(n,a,inv(3));

	int u=inv(n);

	for (int i=0;i<n;i++) a[i]=1ll*a[i]*u%P;

}

void mul(int *a,int *b,int n)

{

	DFT(n,a,3),DFT(n,b,3);

	for (int i=0;i<n;i++) a[i]=1ll*a[i]*b[i]%P;

	IDFT(a,n);

}

void Inv(int *a,int *b,int n)

{

	if (n==1) {for (int i=0;i<t;i++) b[i]=0;b[0]=inv(a[0]);return;}

	Inv(a,b,n>>1);

	for (int i=0;i<n;i++) A[i]=a[i];

	for (int i=n;i<(n<<1);i++) A[i]=0;

	n<<=1;

	DFT(n,A,3),DFT(n,b,3);

	for (int i=0;i<n;i++) b[i]=1ll*b[i]*(P+2-1ll*A[i]*b[i]%P)%P;

	IDFT(b,n);

	n>>=1;

	for (int i=n;i<(n<<1);i++) b[i]=0;

}

void Sqrt(int *a,int *b,int n)

{

	if (n==1) {for (int i=0;i<t;i++) b[i]=0;b[0]=1;return;}

	Sqrt(a,b,n>>1);Inv(b,B,n);

	for (int i=0;i<n;i++) A[i]=a[i];

	for (int i=n;i<(n<<1);i++) A[i]=0;

	n<<=1;

	DFT(n,A,3),DFT(n,b,3),DFT(n,B,3);

	int inv2=inv(2);

	for (int i=0;i<n;i++) b[i]=(b[i]+1ll*A[i]*B[i]%P)*inv2%P;

	IDFT(b,n);

	n>>=1;

	for (int i=n;i<(n<<1);i++) b[i]=0;

}

int main()

{

#ifndef ONLINE_JUDGE

	freopen("sqrt.in","r",stdin);

	freopen("sqrt.out","w",stdout);

	const char LL[]="%I64d\n";

#else

	const char LL[]="%lld\n";

#endif

	n=read();

	for (int i=0;i<n;i++) a[i]=read();

	t=1;while (t<=n) t<<=1;

	Sqrt(a,b,t);

	for (int i=0;i<n;i++) printf("%d ",b[i]);

	return 0;

}

　　

相关文章