稍有常识的人都能看出dp方程
fn=2(n2)−∑i=1n−1(n−1i−1)fi2(n−i2)
稍稍整理一下就可以分治fft了。
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define LL long long
const int p=1004535809,g=3,maxn=600010;
int dp[maxn],pw[maxn],c2[maxn],inv[maxn],fac[maxn],f[maxn],h[maxn],rev[maxn],n;
int inc(int x,int y)
{
x+=y;
return x>=p?x-p:x;
}
int dec(int x,int y)
{
x-=y;
return x<0?x+p:x;
}
int pow(int b,int k)
{
int ret=1;
for (;k;k>>=1,b=(LL)b*b%p)
if (k&1) ret=(LL)ret*b%p;
return ret;
}
int powLL(int b,LL k)
{
int ret=1;
for (;k;k>>=1,b=(LL)b*b%p)
if (k&1) ret=(LL)ret*b%p;
return ret;
}
void fft(int *a,int m,int t,int flag)
{
int x,t1,t2;
for (int i=0;i<m;i++)
if (rev[i]>i) swap(a[i],a[rev[i]]);
for (int i=0;i<t;i++)
for (int j=0;j<m;j+=1<<(i+1))
{
x=0;
for (int k=j;k<j+(1<<i);k++)
{
t1=a[k];
t2=a[k+(1<<i)];
a[k]=inc(t1,(LL)t2*pw[x]%p);
a[k+(1<<i)]=dec(t1,(LL)t2*pw[x]%p);
x+=flag*(1<<(t-i-1));
if (x<0) x+=m;
}
}
if (flag==-1)
{
x=pow(m,p-2);
for (int i=0;i<m;i++) f[i]=(LL)f[i]*x%p;
}
}
void solve(int l,int r)
{
if (l==r)
{
dp[l]=dec(c2[l],(LL)fac[l-1]*dp[l]%p);
return;
}
int mid=(l+r)>>1,m=1,t=0;
solve(l,mid);
for (int i=l;i<=mid;i++) f[i-l]=(LL)dp[i]*inv[i-1]%p;
for (int i=1;i<=r-l;i++) h[i-1]=(LL)c2[i]*inv[i]%p;
while (m<=2*(r-l-1)) m<<=1,t++;
for (int i=mid-l+1;i<m;i++) f[i]=0;
for (int i=r-l;i<m;i++) h[i]=0;
pw[0]=1;
pw[1]=pow(g,(p-1)/m);
for (int i=2;i<m;i++) pw[i]=(LL)pw[i-1]*pw[1]%p;
for (int i=0;i<m;i++)
{
rev[i]=0;
for (int j=0;j<t;j++) rev[i]|=((i>>j)&1)<<(t-j-1);
}
fft(f,m,t,1);
fft(h,m,t,1);
for (int i=0;i<m;i++) f[i]=(LL)f[i]*h[i]%p;
fft(f,m,t,-1);
for (int i=mid+1;i<=r;i++) dp[i]=inc(dp[i],f[i-l-1]);
solve(mid+1,r);
}
int main()
{
scanf("%d",&n);
inv[0]=fac[0]=1;
for (int i=1;i<=n;i++) fac[i]=(LL)fac[i-1]*i%p;
inv[n]=pow(fac[n],p-2);
for (int i=n-1;i>=1;i--) inv[i]=(LL)inv[i+1]*(i+1)%p;
for (int i=1;i<=n;i++) c2[i]=powLL(2,(LL)i*(i-1)/2);
solve(1,n);
printf("%d\n",dp[n]);
}