hdu5730(FFT+多项式求逆）

（最近在班上好无聊，总感觉同学都好拽，也没有能说很多话的小姐姐，而且都很爱学习。老师也很死板，诸多要求。但我还是我，毕竟看过那么多治愈番。）

终有一天会降临，只属于我的苍蓝天使！

稳定的题面传送门

题意是这样的，长度分别为1..n的棒子各有a[i]种，每种都有无限条，问组成长度为n的大棒子有多少种方案，mod313。

设f[i]表示长度为i的方案数，根据题意f[0]=1，考虑dp。

枚举最后一条棒子的长度，就有

f[x]=∑i=0xa[i]f[x−i]

朴素就是n3的，相当于做n次卷积，若用FFT优化卷积就可以n2logn。当然很多dalao一眼就看出这是个分治+FFT的模板题。

所以下面我们考虑生成函数。

是F为f的生成函数，A为a的生成函数，就有

F=F∗A+1

+1是因为f[0]=1，a[0]=0。

然后有F=11−A

就是一个多项式逆元了，复杂度省了一个logn。

没学过？那还不去膜L-leader的博客？

这里mod的是313，FFT误差卡了我一晚，有几个能卡过这题的方法（我不会NTT怎么搞）。

1、单位根要预处理，不要累乘。

2、对于系数都是整数的卷积，把系数>>15，最后结果<<30。

3、对于多项式A*B*C，用2的方法先做A*B，再(A*B)*C。

4、对于多项式加法，在系数表示的时候做，不要在点值表示的时候做。

（NTT就不用想那么多了）

#include <iostream>
#include <fstream>
#include <algorithm>
#include <cmath>
#include <ctime>
#include <cstdio>
#include <cstdlib>
#include <cstring>

using namespace std;
#define mmst(a, b) memset(a, b, sizeof(a))
#define mmcp(a, b) memcpy(a, b, sizeof(b))

typedef long long LL;

const int N=800100,p=313;
const double pi=acos(-1);

struct yy
{
double x,y;
yy(double a=0,double b=0):x(a),y(b){}
};

yy operator +(yy a,yy b){return yy(a.x+b.x,a.y+b.y);}
yy operator -(yy a,yy b){return yy(a.x-b.x,a.y-b.y);}
yy operator *(yy a,yy b){return yy(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}

int n,rev
;
yy omega
[2];

void init(int lim)
{
int k=-1;
n=1;
while(n<lim)
n<<=1,k++;
for(int i=0;i<n;i++)
rev[i]=(rev[i>>1]>>1) | ((i&1)<<k);
for(int i=0;i<n;i++)
{
omega[i][1]=yy(cos(2*pi/n*i),sin(2*pi/n*i));
omega[i][0]=yy(cos(2*pi/n*i),-sin(2*pi/n*i));
}
}

void fft(yy *a,int ops)
{
for(int i=0;i<n;i++)
if(i<rev[i])
swap(a[i],a[rev[i]]);

for(int l=2;l<=n;l<<=1)
{
int m=l>>1;
for(int i=0;i<n;i+=l)
for(int k=0;k<m;k++)
{
yy t=a[i+k+m]*omega[n/l*k][ops];
a[i+k+m]=a[i+k]-t;
a[i+k]=a[i+k]+t;
}
}
if(!ops)
for(int i=0;i<n;i++)
a[i].x/=1.0*n;
}

int nn,a
;
LL pp=p<<15,ans
;
yy X
,Y
;

void Inv(int len)
{
if(len==1)
return;
Inv(len/2);

for(int i=0;i<len/2;i++)
X[i].x=(1.0*ans[i]/32768);
af4e

init(len+1);
fft(X,1);
for(int i=0;i<n;i++)
X[i]=X[i]*X[i];
fft(X,0);

for(int i=0;i<n;i++)
{
X[i].x=(((LL)floor(X[i].x*(1<<30ll)+0.3)%p)+p)%p;
X[i].x/=32768;
X[i].y=0.0;
}

for(int i=0;i<len;i++)
Y[i].x=(1.0*a[i]/32768);

init(len*2);
fft(X,1);
fft(Y,1);
for(int i=0;i<n;i++)
X[i]=X[i]*Y[i];
fft(X,0);

for(int i=0;i<len;i++)
ans[i]=(2*ans[i]-((LL)floor(X[i].x*(1<<30ll)+0.3)%p)+p)%p;

for(int i=0;i<n;i++)
X[i].x=X[i].y=Y[i].x=Y[i].y=0.0;
}
int main()
{
while(1)
{
cin>>nn;
if(!nn)
break;
mmst(ans,0);
ans[0]=1;
mmst(a,0);

for(int i=1;i<=nn;i++)
scanf("%d",&a[i]);

for(int i=1;i<=nn;i++)
a[i]=p-a[i]%p;
a[0]=1;

init(nn+1);
Inv(n);
cout<<(ans[nn]+pp)%p<<endl;
}

return 0;
}