给定(n-1)次多项式(A(x)),求一个在(mod\, x^n)意义下的多项式(B(x)),使得(B(x)^2≡A(x) (mod\, x^n)),取零次项系数最小的作为答案
系数(mod\, 998244353)
根据题意:
(B(x)^2-A(x)≡0 (mod\, x^n))
设(G(B(x))=B(x)^2-A(x)),套牛顿迭代
(B(x)≡B_0(x)-frac{G(B_0(x))}{G'(B_0(x))})
整理一下
(B(x)=frac{A(x)+B_0(x)^2}{2B_0(x)})
直接乘常数有点大,稍微变化一下
(B(x)=frac{1}{2}(A(x)B_0(x)^{-1}+B_0(x)))
递归求解
#include<bits/stdc++.h>
using namespace std;
namespace red{
#define int long long
#define y1 miao
#define eps (1e-8)
inline int read()
{
int x=0;char ch,f=1;
for(ch=getchar();(ch<'0'||ch>'9')&&ch!='-';ch=getchar());
if(ch=='-') f=0,ch=getchar();
while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+ch-'0';ch=getchar();}
return f?x:-x;
}
const int N=4e5+10,p=998244353;
int n,inv2;
int a[N],b[N],c[N],pos[N];
int invb[N];
inline int fast(int x,int k)
{
int ret=1;
while(k)
{
if(k&1) ret=ret*x%p;
x=x*x%p;
k>>=1;
}
return ret;
}
inline void ntt(int limit,int *a,int inv)
{
for(int i=0;i<limit;++i)
if(i<pos[i]) swap(a[i],a[pos[i]]);
for(int mid=1;mid<limit;mid<<=1)
{
int Wn=fast(3,(p-1)/(mid<<1));
for(int r=mid<<1,j=0;j<limit;j+=r)
{
int w=1;
for(int k=0;k<mid;++k,w=w*Wn%p)
{
int x=a[j+k],y=w*a[j+k+mid]%p;
a[j+k]=(x+y)%p;
a[j+k+mid]=(x-y+p)%p;
}
}
}
if(inv) return;
inv=fast(limit,p-2);reverse(a+1,a+limit);
for(int i=0;i<limit;++i) a[i]=a[i]*inv%p;
}
inline void poly_inv(int pw,int *a,int *b)
{
if(pw==1){b[0]=fast(a[0],p-2);return;}
poly_inv((pw+1)>>1,a,b);
int len=0,limit=1;
while(limit<(pw<<1)) limit<<=1,++len;
for(int i=0;i<limit;++i) pos[i]=(pos[i>>1]>>1)|((i&1)<<(len-1));
for(int i=0;i<pw;++i) c[i]=a[i];
for(int i=pw;i<limit;++i) c[i]=0;
ntt(limit,c,1);ntt(limit,b,1);
for(int i=0;i<limit;++i) b[i]=((2-c[i]*b[i]%p)+p)%p*b[i]%p;
ntt(limit,b,0);
for(int i=pw;i<limit;++i) b[i]=0;
}
inline void poly_sqrt(int pw,int *a,int *b)
{
if(pw==1){b[0]=1;return;}
poly_sqrt((pw+1)>>1,a,b);
for(int i=0;i<(pw<<1);++i) invb[i]=0;//记得清空数组
poly_inv(pw,b,invb);
int len=0,limit=1;
while(limit<(pw<<1)) limit<<=1,++len;
for(int i=0;i<limit;++i) pos[i]=(pos[i>>1]>>1)|((i&1)<<(len-1));
for(int i=0;i<pw;++i) c[i]=a[i];
for(int i=pw;i<limit;++i) c[i]=0;
ntt(limit,c,1);ntt(limit,invb,1);
for(int i=0;i<limit;++i) c[i]=c[i]*invb[i]%p;
ntt(limit,c,0);
for(int i=0;i<pw;++i) b[i]=(b[i]+c[i])%p*inv2%p;
for(int i=pw;i<limit;++i) b[i]=0;
}
inline void main()
{
n=read();inv2=fast(2,p-2);
for(int i=0;i<n;++i) a[i]=read();
poly_sqrt(n,a,b);
for(int i=0;i<n;++i) printf("%lld ",b[i]);
}
}
signed main()
{
red::main();
return 0;
}