bzoj3771 Triple

Posted 2020-10-25 Ren_Ivan

生成函数+FFT

之前做了几道感觉推完式子都是FFT的纯板子题，然后做这道题想到了容斥，但是不知道应该怎么上FFT，因为这个并不是纯卷积，或者说纯卷积十分累赘复杂。

FFT其实只是一个多项式点值表达式和系数表达式之间转换的一个工具，我们要做的就是理清楚各个多项式之间准确的关系，或者把卷积转化成原始的for循环理解

然而生成函数的原理是不能变的，生成函数i次方的系数就是代表对应权值为i的方案数，否则权值不对应，乘起来就乱套了。

所以我们设三个生成函数，A表示选一个权值为x的斧子，B表示选两个权值为x/2的斧子，C表示选三个权值为x/3的斧子。

那么选一个的方案数就是A，选两个的就是(A*A-B)/2，这里我们考虑原始的for循环，A*A是两层从1到n的循环，然而我们要去掉i=j的情况，就是B，然后还应该除去组合顺序，所以就是上式，选三个的也差不多，首先有三层循环，然后我们应该去掉i=j||j=k||i=k的情况，就是另外两层循环：A*B*3，但是这里i=j=k的情况多减了两次，所以我们要加上C，最后在除去组合顺序，就是(A*A*A-A*B*3+C*2)/6。

说的太多了，还是因为自己不明白。233

 1 #include <cstdio>
 2 #include <cstring>
 3 #include <iostream>
 4 #include <algorithm>
 5 #include <cmath>
 6 #define MAXN 133333
 7 using namespace std;
 8 const double pi=acos(-1.0);
 9 struct cmx{
10     double x,y;
11     cmx(){}
12     cmx(double a,double b){x=a,y=b;}
13     cmx operator + (cmx a){return cmx(x+a.x,y+a.y);}
14     cmx operator - (cmx a){return cmx(x-a.x,y-a.y);}
15     cmx operator * (cmx a){return cmx(x*a.x-y*a.y,x*a.y+y*a.x);}
16     cmx operator / (double a){return cmx(x/a,y/a);}
17     cmx operator * (double a){return cmx(x*a,y*a);}
18 }A[MAXN],B[MAXN],C[MAXN],ans[MAXN];
19 int n,N,rev[MAXN],maxn;
20 void fft(cmx *a,int o){
21     cmx dan,now,t;
22     register int i,j,k;
23     for(i=0;i<N;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
24     for(k=2;k<=N;k<<=1){
25         dan=cmx(cos(2*pi/k),o*sin(2*pi/k));
26         for(j=0;j<N;j+=k){
27             now=cmx(1,0);
28             for(i=0;i<(k>>1);i++,now=now*dan){
29                 t=now*a[i+j+(k>>1)];
30                 a[i+j+(k>>1)]=a[i+j]-t;
31                 a[i+j]=a[i+j]+t;
32             }
33         }
34     }
35     if(o==-1)for(i=0;i<N;i++)a[i].x/=N;
36 }
37 int main(){
38     scanf("%d",&n);
39     register int i;
40     register long long now;
41     for(int i=1,x;i<=n;i++){
42         scanf("%d",&x);
43         A[x].x+=1;
44         B[x*2].x+=1;
45         C[x*3].x+=1;
46         maxn=max(maxn,3*x);
47     }
48     for(N=1;N<maxn;N<<=1);
49     for(i=0;i<N;i++){
50         if(i&1)rev[i]=(rev[i>>1]>>1)|(N>>1);
51         else rev[i]=rev[i>>1]>>1;
52     }
53     fft(A,1);fft(B,1);fft(C,1);
54     for(i=0;i<N;i++)
55         ans[i]=A[i]+(A[i]*A[i]-B[i])/2+(A[i]*A[i]*A[i]-A[i]*B[i]*3+C[i]*2)/6;
56     fft(ans,-1);
57     for(i=0;i<maxn;i++){
58         now=1ll*round(ans[i].x);
59         if(now>0)printf("%d %lld\n",i,now);
60     }
61     return 0;
62 }

View Code