题解 P5277 【【模板】多项式开根（加强版）】
SSerxhs
2019-03-30 11:02:50
以下内容p均指模数998244353
前置芝士：二次剩余（模奇素数意义下）
二次剩余：求y使得$y^2=x(mod\ p)$
今年WC的zzt大佬讲课课件中提到几种求二次剩余的算法，这题应该都能用
## 法一：暴力
求出模p意义下的一个原根(其中一个是3)，bsgs求出x的指标（即求出a使得$3^a=x(mod\ p)$），然后直接指标除以二即可，复杂度$O(\sqrt{p})$
## 法二：Cipolla算法
找到一个b，使得$(b^2-x)$为模意义下的**二次非剩余**，设$(b^2-x)$为w则$y=(b+\sqrt{w})^\frac{p+1}{2}$，这里的$\sqrt{w}$不能具体求出（实际上也没有），计算时它类似于复数里的i，整个数字需要划分实部虚部并重载乘法运算。
b的求法？直接rand几个b，由于恰好有一半数字是可以当成b的，另一半不能，所以期望次数2次，总复杂度$O(log_2p)$
upd:证明如下：根据二项式定理，$y^2=b^{p+1}+\sqrt{w^{p+1}}$(其余项都含有p，可以直接消去)，根据欧拉判定可知$\sqrt{w^{p+1}}=w^{\frac{p+ 1}{2}}=-w$，根据欧拉定理可知$b^{p+1}=b^2$，所以$y^2=b^2-w=x$
为什么不讲WC的第三种做法Tonelli-Shanks 算法？~~我不会~~
以下是本题正片
## 解法一：暴力求$A^{\frac{1}{2}}$
前置芝士：[多项式幂函数（加强版）](https://www.luogu.org/problemnew/show/P5273)
把代码粘过来，注意处理当前几位为0时原来的做法会把全部值都赋值成0（因为次数太高了，原本的做法最后是把数组**右移**），这里考虑本身需要是开根号，如果后几位有k个0只要做完幂函数之后重新**左移**k/2位就可以了
另一个细节：对于数字x来说，$x^{p-1}=1(mod\ p)$，然而对于多项式f(x)来说$f^p(x)=1(mod\ p)(mod\ x^n)$，所以实际上做的是f(x)的$2^{p-2}$次方，但提取第一个非0时应该照常求其二次剩余，不能也提取${a_0}^{(2^{p-2})}$出来。
代码如下，建议吸氧使用此方法
```cpp
// luogu-judger-enable-o2
#include <stdio.h>
#include <string.h>
#include <algorithm>
#include <stdlib.h>
using namespace std;
typedef long long ll;
typedef pair<int,int> pa;
const int N=2.7e5,p=998244353,I=86583718,II=911660635;
int r[N],ig[N],yg[N],invf[N],sqrtf[N],sqrtg[N],cdqf[N],cdqg[N],lnf[N],expf[N],expg[N];
int ksmf[N],inv[N],f[N],g[N],sinf[N],sing[N],cosf[N],cosg[N],tanf[N],tang[N],cotf[N];
int cotg[N];
int n,m,i,j,l,limit,biglimit,c,w;
pa ys,yys;
inline int sj()
{
return rand()<<15|rand();
}
pa operator *(pa x,pa y)
{
ys.first=((ll)x.first*y.first+(ll)x.second*y.second%p*w)%p;
ys.second=((ll)x.first*y.second+(ll)x.second*y.first)%p;
return ys;
}
inline void dbug(int *a,int n)
{
for (int i=0;i<n;i++) printf("%d ",a[i]);puts("");
}
inline void read(int &x)
{
c=getchar();
while ((c<48)||(c>57)) c=getchar();
x=c^48;c=getchar();
while ((c>=48)&&(c<=57))
{
x=x*10+(c^48);
c=getchar();
}
}
inline int ksm(int x,int y)
{
int r=1;
while (y)
{
if (y&1) r=(ll)r*x%p;
x=(ll)x*x%p;
y>>=1;
}
return r;
}
inline pa ksm(pa x,int y)
{
yys.first=1;yys.second=0;
while (y)
{
if (y&1) yys=yys*x;
x=x*x;
y>>=1;
}
return yys;
}
inline int mosqrt(int x)
{
pa a;
int y=rand()%p;
while (ksm(w=((ll)y*y%p-x+p)%p,p-1>>1)<=1) y=rand()%p;
if (w<0) w+=p;
a.first=y;a.second=1;
a=ksm(a,p+1>>1);
return min(a.first,p-a.first);
}
inline void ycl(int l,int limit)
{
for (int i=1;i<limit;i++) r[i]=r[i>>1]>>1|(i&1)<<l;
}
inline void gg(int limit)
{
int i;
ig[limit]=ksm(yg[limit]=ksm(3,(p-1)/limit),p-2);
for (i=limit>>1;i;i>>=1)
{
yg[i]=(ll)yg[i<<1]*yg[i<<1]%p;
ig[i]=(ll)ig[i<<1]*ig[i<<1]%p;
}
inv[1]=1;
for (i=2;i<limit;i++) inv[i]=p-(ll)p/i*inv[p%i]%p;
}
inline void reverse(int limit,int *f)
{
int lim=limit>>1;
for (i=0;i<lim;i++) swap(f[i],f[limit-i-1]);
}
inline void dao(int *a,int n)
{
for (i=1;i<n;i++) a[i-1]=(ll)a[i]*i%p;a[n-1]=0;
}
inline void ji(int *a,int n)
{
for (i=n-1;~i;i--) a[i+1]=(ll)a[i]*inv[i+1]%p;
a[0]=0;
}
void dft(int *a,int xs,int limit)
{
register int i,j,k,l,w,wn,b,c;
for (i=1;i<limit;i++) if (i<r[i]) swap(a[i],a[r[i]]);
for (i=1;i<limit;i=l)
{
l=i<<1;
if (xs) wn=ig[l]; else wn=yg[l];
for (j=0;j<limit;j+=l)
{
w=1;
for (k=0;k<i;k++,w=(ll)w*wn%p)
{
b=a[j|k];c=(ll)w*a[j|k|i]%p;
a[j|k]=(b+c)%p;
a[j|k|i]=(b-c+p)%p;
}
}
}
if (xs)
{
limit=ksm(limit,p-2);
for (i=0;i<xs;i++) a[i]=(ll)a[i]*limit%p;
}
}
void polymultiply(int *f,int *g,int limit)
{
int i;
dft(f,0,limit<<1);dft(g,0,limit<<1);
for (i=0;i<limit<<1;i++) f[i]=(ll)f[i]*g[i]%p;
dft(f,limit,limit<<1);
}
void polyinv(int *f,int *g,int biglim)
{
int i,j,l=1,limit;
memset(g,0,biglim<<3);
memset(invf,0,biglim<<3);
g[0]=ksm(f[0],p-2);
for (i=2;i<=biglim;i=limit,l++)
{
limit=i<<1;
memcpy(invf,f,limit<<1);
for (j=1;j<limit;j++) r[j]=r[j>>1]>>1|(j&1)<<l;
dft(invf,0,limit);dft(g,0,limit);
for (j=0;j<limit;j++) g[j]=(ll)g[j]*(p+2-(ll)g[j]*invf[j]%p)%p;
dft(g,i,limit);
memset(g+i,0,limit<<1);
}
}
void polysqrt(int *f,int *g,int biglim)//不动f
{
memset(sqrtf,0,biglim<<3);
memset(sqrtg,0,biglim<<3);
memset(g,0,biglim<<3);
int i,j,l=1,limit;
g[0]=mosqrt(f[0]);
if ((ll)g[0]*g[0]%p!=f[0]) while (1);
for (i=2;i<=biglim;i=limit,l++)
{
limit=i<<1;
polyinv(g,sqrtg,i);
memcpy(sqrtf,f,limit<<1);
dft(sqrtf,0,limit);dft(sqrtg,0,limit);dft(g,0,limit);
for (j=0;j<limit;j++)
{
g[j]=(ll)sqrtg[j]*((ll)g[j]*g[j]%p+sqrtf[j])%p;
if (g[j]&1) g[j]=g[j]+p>>1; else g[j]>>=1;
}
dft(g,i,limit);
memset(g+i,0,limit<<1);
}
}
/*void polycdq(int l,int r,int *f,int *g)
{
if (l==r) return;
int mid=l+r>>1;
polycdq(l,mid,f,g);
int i,limit,li,siz=r-l;
limit=1;li=-1;
while (limit<=siz)
{
limit<<=1;++li;
}
for (i=l;i<=mid;i++) cdqf[i-l]=f[i];
for (i=mid-l+1;i<lim;i++) cdqf[i]=0;
for (i=1;i<=r-l;i++) cdqg[i]=g[i];G[0]=0;
for (i=r-l+1;i<lim;i++) cdqg[i]=0;
polymultiply(cdqf,cdqg,li,limit);
for (i=mid+1;i<=r;i++)
{
f[i]+=cdqf[i-l];
if (f[i]>=p) f[i]-=p;
}
polycdq(mid+1,r,f,g);
}*/
void polydivide(int *f,int *g,int *q,int *r,int limit)
{
}
void polyln(int *f,int *g,int biglim)
{
int i,limit=biglim<<1;
memcpy(lnf,f,biglim<<2);
memset(lnf+biglim,0,biglim<<2);
polyinv(f,g,biglim);
dao(lnf,biglim);
dft(lnf,0,limit);dft(g,0,limit);
for (i=0;i<limit;i++) g[i]=(ll)g[i]*lnf[i]%p;
dft(g,biglim,limit);
memset(g+biglim,0,biglim<<2);
ji(g,biglim);
}
void polyexp(int *f,int *g,int biglim)
{
memset(g,0,biglim<<3);
memset(expf,0,biglim<<3);
memset(expg,0,biglim<<3);
g[0]=1;
int i,j,l=1,limit;
for (i=2;i<=biglim;i=limit,l++)
{
limit=i<<1;
memcpy(expf,f,limit<<1);
polyln(g,expg,i);
dft(expf,0,limit);dft(g,0,limit);dft(expg,0,limit);
for (j=0;j<limit;j++) g[j]=(ll)g[j]*(1+expf[j]-expg[j]+p)%p;
dft(g,i,limit);
memset(g+i,0,limit<<1);
}
}
void polyksm1(int *f,int *g,int biglim,int cs)
{
if (cs==0)
{
g[0]=1;
return;
}
int i,j,k=1,wy=0;
if (f[0]==0)
{
for (i=1;i<biglim;i++) if (f[i]) break;
for (j=0;j<biglim-i;j++) f[j]=f[j+i];
memset(f+biglim-i,0,i<<2);
wy=cs>>1;
}
if (f[0]>1)
{
k=ksm(f[0],p-2);
for (i=1;i<biglim;i++) f[i]=(ll)f[i]*k%p;
k=mosqrt(f[0]);f[0]=1;//printf("%d\n",k);
}
polyln(f,ksmf,biglim);
for (i=0;i<biglim;i++) ksmf[i]=(ll)ksmf[i]*cs%p;
polyexp(ksmf,g,biglim);
if (k!=1)
{
for (i=0;i<biglim;i++) g[i]=(ll)g[i]*k%p;
}
if (wy)
{
for (i=wy;i<biglim;i++) g[i-wy]=g[i],g[i]=0;
}
}
void polyksm2(int *f,int *g,int biglim,int cs)
{
int limit=1,i,l=-1;
while (limit<=biglim)
{
limit<<=1;++l;
}
for (i=1;i<limit;i++) r[i]=r[i>>1]>>1|(i&1)<<l;g[0]=1;
while (cs)
{
if (cs&1)
{
memcpy(ksmf,f,biglim<<2);
memset(ksmf+biglim,0,biglim<<2);
polymultiply(g,ksmf,biglim);
memset(g+biglim,0,biglim<<2);
}
dft(f,0,limit);
for (i=0;i<limit;i++) f[i]=(ll)f[i]*f[i]%p;
dft(f,biglim,limit);
memset(f+biglim,0,biglim<<2);
cs>>=1;
}
}
void polyksm3(int *f,int *g,int biglim,int cs)
{
int i;
memset(ksmf,0,biglim<<3);
memset(g,0,biglim<<3);
polyln(f,ksmf,biglim);
for (i=0;i<biglim;i++) ksmf[i]=(ll)ksmf[i]*cs%p;
polyexp(ksmf,g,biglim);
}
void polysin(int *f,int *g,int biglim)
{
int i;
for (i=0;i<biglim;i++) g[i]=(ll)f[i]*I%p;
polyexp(g,sinf,biglim);
for (i=0;i<biglim;i++) if ((g[i]=p-(ll)f[i]*I%p)==p) g[i]-=p;
polyexp(g,sing,biglim);
for (i=0;i<biglim;i++) g[i]=(ll)(sinf[i]-sing[i]+p)*II%p;
for (i=0;i<biglim;i++) if (g[i]&1) g[i]=g[i]+p>>1; else g[i]>>=1;
}
void polycos(int *f,int *g,int biglim)
{
int i;
for (i=0;i<biglim;i++) g[i]=(ll)f[i]*I%p;
polyexp(g,cosf,biglim);
for (i=0;i<biglim;i++) if ((g[i]=p-(ll)f[i]*I%p)==p) g[i]-=p;
polyexp(g,cosg,biglim);
for (i=0;i<biglim;i++) if ((g[i]=cosf[i]+cosg[i])>=p) g[i]-=p;
for (i=0;i<biglim;i++) if (g[i]&1) g[i]=g[i]+p>>1; else g[i]>>=1;
}
void polytan(int *f,int *g,int biglim)
{
polysin(f,tanf,biglim);
polycos(f,tang,biglim);
polyinv(tang,g,biglim);
polymultiply(g,tanf,biglim);
}
void polycot(int *f,int *g,int biglim)
{
polysin(f,cotf,biglim);
polycos(f,cotg,biglim);
polyinv(cotf,g,biglim);
polymultiply(g,cotg,biglim);
}
int main()
{
read(n);
limit=1;
while (limit<n) limit<<=1;
for (i=0;i<n;i++) read(f[i]);gg(limit<<1);
polyksm1(f,g,limit,ksm(2,p-2));
for (i=0;i<n;i++) printf("%d ",g[i]);
}
```
## 解法二：推柿子
前置芝士：[多项式求逆](https://www.luogu.org/problemnew/show/P4238)、套路倍增（当然直接牛顿迭代也不是不行）
设$F^2(x)=A(x)(mod\ x^n),G^2(x)=A(x)(mod\ x^{2n})$(以下直接用A、F和G代表多项式)
则$(F-G)^2=0(mod\ x^{2n})$
拆开发现$F^2+G^2=2FG(mod\ x^{2n})$
$G^2=A(mod\ x^{2n})$所以$G=\frac{F^2+A}{2F}$
然后直接倍增即可，边界条件求一次二次剩余
代码如下
```cpp
#include <stdio.h>
#include <string.h>
#include <algorithm>
#include <stdlib.h>
using namespace std;
typedef long long ll;
typedef pair<int,int> pa;
const int N=2.7e5,p=998244353,I=86583718,II=911660635;
int r[N],ig[N],yg[N],invf[N],sqrtf[N],sqrtg[N],cdqf[N],cdqg[N],lnf[N],expf[N],expg[N];
int ksmf[N],inv[N],f[N],g[N],sinf[N],sing[N],cosf[N],cosg[N],tanf[N],tang[N],cotf[N];
int cotg[N];
int n,m,i,j,l,limit,biglimit,c,w;
pa ys,yys;
inline int sj()
{
return rand()<<15|rand();
}
pa operator *(pa x,pa y)
{
ys.first=((ll)x.first*y.first+(ll)x.second*y.second%p*w)%p;
ys.second=((ll)x.first*y.second+(ll)x.second*y.first)%p;
return ys;
}
inline void dbug(int *a,int n)
{
for (int i=0;i<n;i++) printf("%d ",a[i]);puts("");
}
inline void read(int &x)
{
c=getchar();
while ((c<48)||(c>57)) c=getchar();
x=c^48;c=getchar();
while ((c>=48)&&(c<=57))
{
x=x*10+(c^48);
c=getchar();
}
}
inline int ksm(int x,int y)
{
int r=1;
while (y)
{
if (y&1) r=(ll)r*x%p;
x=(ll)x*x%p;
y>>=1;
}
return r;
}
inline pa ksm(pa x,int y)
{
yys.first=1;yys.second=0;
while (y)
{
if (y&1) yys=yys*x;
x=x*x;
y>>=1;
}
return yys;
}
inline int mosqrt(int x)
{
pa a;
int y=rand()%p;
while (ksm(w=((ll)y*y%p-x+p)%p,p-1>>1)<=1) y=rand()%p;
if (w<0) w+=p;
a.first=y;a.second=1;
a=ksm(a,p+1>>1);
return min(a.first,p-a.first);
}
inline void ycl(int l,int limit)
{
for (int i=1;i<limit;i++) r[i]=r[i>>1]>>1|(i&1)<<l;
}
inline void gg(int limit)
{
int i;
ig[limit]=ksm(yg[limit]=ksm(3,(p-1)/limit),p-2);
for (i=limit>>1;i;i>>=1)
{
yg[i]=(ll)yg[i<<1]*yg[i<<1]%p;
ig[i]=(ll)ig[i<<1]*ig[i<<1]%p;
}
inv[1]=1;
for (i=2;i<limit;i++) inv[i]=p-(ll)p/i*inv[p%i]%p;
}
inline void reverse(int limit,int *f)
{
int lim=limit>>1;
for (i=0;i<lim;i++) swap(f[i],f[limit-i-1]);
}
inline void dao(int *a,int n)
{
for (i=1;i<n;i++) a[i-1]=(ll)a[i]*i%p;a[n-1]=0;
}
inline void ji(int *a,int n)
{
for (i=n-1;~i;i--) a[i+1]=(ll)a[i]*inv[i+1]%p;
a[0]=0;
}
void dft(int *a,int xs,int limit)
{
register int i,j,k,l,w,wn,b,c;
for (i=1;i<limit;i++) if (i<r[i]) swap(a[i],a[r[i]]);
for (i=1;i<limit;i=l)
{
l=i<<1;
if (xs) wn=ig[l]; else wn=yg[l];
for (j=0;j<limit;j+=l)
{
w=1;
for (k=0;k<i;k++,w=(ll)w*wn%p)
{
b=a[j|k];c=(ll)w*a[j|k|i]%p;
a[j|k]=(b+c)%p;
a[j|k|i]=(b-c+p)%p;
}
}
}
if (xs)
{
limit=ksm(limit,p-2);
for (i=0;i<xs;i++) a[i]=(ll)a[i]*limit%p;
}
}
void polymultiply(int *f,int *g,int limit)
{
int i;
dft(f,0,limit<<1);dft(g,0,limit<<1);
for (i=0;i<limit<<1;i++) f[i]=(ll)f[i]*g[i]%p;
dft(f,limit,limit<<1);
}
void polyinv(int *f,int *g,int biglim)
{
int i,j,l=1,limit;
memset(g,0,biglim<<3);
memset(invf,0,biglim<<3);
g[0]=ksm(f[0],p-2);
for (i=2;i<=biglim;i=limit,l++)
{
limit=i<<1;
memcpy(invf,f,limit<<1);
for (j=1;j<limit;j++) r[j]=r[j>>1]>>1|(j&1)<<l;
dft(invf,0,limit);dft(g,0,limit);
for (j=0;j<limit;j++) g[j]=(ll)g[j]*(p+2-(ll)g[j]*invf[j]%p)%p;
dft(g,i,limit);
memset(g+i,0,limit<<1);
}
}
void polysqrt(int *f,int *g,int biglim)//不动f
{
memset(sqrtf,0,biglim<<3);
memset(sqrtg,0,biglim<<3);
memset(g,0,biglim<<3);
int i,j,l=1,limit;
g[0]=mosqrt(f[0]);
if ((ll)g[0]*g[0]%p!=f[0]) while (1);
for (i=2;i<=biglim;i=limit,l++)
{
limit=i<<1;
polyinv(g,sqrtg,i);
memcpy(sqrtf,f,limit<<1);
dft(sqrtf,0,limit);dft(sqrtg,0,limit);dft(g,0,limit);
for (j=0;j<limit;j++)
{
g[j]=(ll)sqrtg[j]*((ll)g[j]*g[j]%p+sqrtf[j])%p;
if (g[j]&1) g[j]=g[j]+p>>1; else g[j]>>=1;
}
dft(g,i,limit);
memset(g+i,0,limit<<1);
}
}
/*void polycdq(int l,int r,int *f,int *g)
{
if (l==r) return;
int mid=l+r>>1;
polycdq(l,mid,f,g);
int i,limit,li,siz=r-l;
limit=1;li=-1;
while (limit<=siz)
{
limit<<=1;++li;
}
for (i=l;i<=mid;i++) cdqf[i-l]=f[i];
for (i=mid-l+1;i<lim;i++) cdqf[i]=0;
for (i=1;i<=r-l;i++) cdqg[i]=g[i];G[0]=0;
for (i=r-l+1;i<lim;i++) cdqg[i]=0;
polymultiply(cdqf,cdqg,li,limit);
for (i=mid+1;i<=r;i++)
{
f[i]+=cdqf[i-l];
if (f[i]>=p) f[i]-=p;
}
polycdq(mid+1,r,f,g);
}*/
void polydivide(int *f,int *g,int *q,int *r,int limit)
{
}
void polyln(int *f,int *g,int biglim)
{
int i,limit=biglim<<1;
memcpy(lnf,f,biglim<<2);
memset(lnf+biglim,0,biglim<<2);
polyinv(f,g,biglim);
dao(lnf,biglim);
dft(lnf,0,limit);dft(g,0,limit);
for (i=0;i<limit;i++) g[i]=(ll)g[i]*lnf[i]%p;
dft(g,biglim,limit);
memset(g+biglim,0,biglim<<2);
ji(g,biglim);
}
void polyexp(int *f,int *g,int biglim)
{
memset(g,0,biglim<<3);
memset(expf,0,biglim<<3);
memset(expg,0,biglim<<3);
g[0]=1;
int i,j,l=1,limit;
for (i=2;i<=biglim;i=limit,l++)
{
limit=i<<1;
memcpy(expf,f,limit<<1);
polyln(g,expg,i);
dft(expf,0,limit);dft(g,0,limit);dft(expg,0,limit);
for (j=0;j<limit;j++) g[j]=(ll)g[j]*(1+expf[j]-expg[j]+p)%p;
dft(g,i,limit);
memset(g+i,0,limit<<1);
}
}
void polyksm1(int *f,int *g,int biglim,int cs)
{
if (cs==0)
{
g[0]=1;
return;
}
int i,j,k=1,wy=0;
if (f[0]==0)
{
for (i=1;i<biglim;i++) if (f[i]) break;
for (j=0;j<biglim-i;j++) f[j]=f[j+i];
memset(f+biglim-i,0,i<<2);
if ((ll)i*cs>=biglim) return;
wy=i*cs;
}
if (f[0]>1)
{
k=ksm(f[0],p-2);
for (i=1;i<biglim;i++) f[i]=(ll)f[i]*k%p;
k=ksm(f[0],cs);f[0]=1;
}
polyln(f,ksmf,biglim);
for (i=0;i<biglim;i++) ksmf[i]=(ll)ksmf[i]*cs%p;
polyexp(ksmf,g,biglim);
if (k!=1)
{
for (i=0;i<biglim;i++) g[i]=(ll)g[i]*k%p;
}
if (wy)
{
for (i=biglim-1;i>=wy;i--) g[i]=g[i-wy];
memset(g,0,wy<<2);
}
}
void polyksm2(int *f,int *g,int biglim,int cs)
{
int limit=1,i,l=-1;
while (limit<=biglim)
{
limit<<=1;++l;
}
for (i=1;i<limit;i++) r[i]=r[i>>1]>>1|(i&1)<<l;g[0]=1;
while (cs)
{
if (cs&1)
{
memcpy(ksmf,f,biglim<<2);
memset(ksmf+biglim,0,biglim<<2);
polymultiply(g,ksmf,biglim);
memset(g+biglim,0,biglim<<2);
}
dft(f,0,limit);
for (i=0;i<limit;i++) f[i]=(ll)f[i]*f[i]%p;
dft(f,biglim,limit);
memset(f+biglim,0,biglim<<2);
cs>>=1;
}
}
void polyksm3(int *f,int *g,int biglim,int cs)
{
int i;
memset(ksmf,0,biglim<<3);
memset(g,0,biglim<<3);
polyln(f,ksmf,biglim);
for (i=0;i<biglim;i++) ksmf[i]=(ll)ksmf[i]*cs%p;
polyexp(ksmf,g,biglim);
}
void polysin(int *f,int *g,int biglim)
{
int i;
for (i=0;i<biglim;i++) g[i]=(ll)f[i]*I%p;
polyexp(g,sinf,biglim);
for (i=0;i<biglim;i++) if ((g[i]=p-(ll)f[i]*I%p)==p) g[i]-=p;
polyexp(g,sing,biglim);
for (i=0;i<biglim;i++) g[i]=(ll)(sinf[i]-sing[i]+p)*II%p;
for (i=0;i<biglim;i++) if (g[i]&1) g[i]=g[i]+p>>1; else g[i]>>=1;
}
void polycos(int *f,int *g,int biglim)
{
int i;
for (i=0;i<biglim;i++) g[i]=(ll)f[i]*I%p;
polyexp(g,cosf,biglim);
for (i=0;i<biglim;i++) if ((g[i]=p-(ll)f[i]*I%p)==p) g[i]-=p;
polyexp(g,cosg,biglim);
for (i=0;i<biglim;i++) if ((g[i]=cosf[i]+cosg[i])>=p) g[i]-=p;
for (i=0;i<biglim;i++) if (g[i]&1) g[i]=g[i]+p>>1; else g[i]>>=1;
}
void polytan(int *f,int *g,int biglim)
{
polysin(f,tanf,biglim);
polycos(f,tang,biglim);
polyinv(tang,g,biglim);
polymultiply(g,tanf,biglim);
}
void polycot(int *f,int *g,int biglim)
{
polysin(f,cotf,biglim);
polycos(f,cotg,biglim);
polyinv(cotf,g,biglim);
polymultiply(g,cotg,biglim);
}
int main()
{
read(n);
limit=1;
while (limit<n) limit<<=1;
for (i=0;i<n;i++) read(f[i]);gg(limit<<1);
polysqrt(f,g,limit);
for (i=0;i<n;i++) printf("%d ",g[i]);
}
```