BZOJ 3771: Triple [快速傅里叶变换 生成函数 容斥原理]

时间：2023-03-08 21:46:35

题意：n个物品，可以用1/2/3个不同的物品组成不同的价值，求每种价值有多少种方案(顺序不同算一种)

【生成函数】：

构造这么一个多项式函数g(x)，使得n次项系数为a[n]。

普通型生成函数用于解决多重集的组合问题

生成函数的x无实际意义 通常可以化为一个简单的式子

组合数的生成函数 A(x)=(1+x)^n  C(n,k)=a[k] 可以这么想，一次要么选择1要么选择x，选择x系数就会+1

生成函数系数为方案数，次数为价值

A(x) 选一个

B(x) A每项平方 选两个

C(x) A每项三次方 选三个

然后容斥原理算答案 听好想的看代码吧

注意计算的时候可以一直用点值，最后在再IDFT变系数表示

#include <iostream>

#include <cstdio>

#include <string>

#include <algorithm>

#include <cmath>

using namespace std;

const int N=+;

inline int read(){

    char c=getchar();int x=,f=;

    while(c<''||c>''){if(c=='-')f=-;c=getchar();}

    while(c>=''&&c<=''){x=x*+c-'';c=getchar();}

    return x*f;

}

const double PI=acos(-);

struct Vector{

    double x,y;

    Vector(double a=,double b=):x(a),y(b){}

};

typedef Vector CD;

Vector operator +(Vector a,Vector b){return Vector(a.x+b.x,a.y+b.y);}

Vector operator -(Vector a,Vector b){return Vector(a.x-b.x,a.y-b.y);}

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

Vector operator *(Vector a,double b){return Vector(a.x*b,a.y*b);}

Vector operator /(Vector a,double b){return Vector(a.x/b,a.y/b);}

Vector conj(Vector a){return Vector(a.x,-a.y);}

struct FastFourierTransform{

    int n,rev[N];

    CD omega[N],omegaInv[N];

    void ini(int m){

        n=;

        while(n<m) n<<=;

        for(int k=;k<n;k++)

            omega[k]=CD(cos(*PI/n*k),sin(*PI/n*k)),

            omegaInv[k]=conj(omega[k]);

        int k=;

        while((<<k)<n) k++;

        for(int i=;i<n;i++){

            int t=;

            for(int j=;j<k;j++) if(i&(<<j)) t|=(<<(k-j-));

            rev[i]=t;

        }

    }

    void transform(CD *a,CD *omega){

        for(int i=;i<n;i++) if(i<rev[i]) swap(a[i],a[rev[i]]);

        for(int l=;l<=n;l<<=){

            int m=l>>;

            for(CD *p=a;p!=a+n;p+=l)

                for(int k=;k<m;k++){

                    CD t=omega[n/l*k]*p[k+m];

                    p[k+m]=p[k]-t;

                    p[k]=p[k]+t;

                }

        }

    }

    void DFT(CD *a,int flag){

        if(flag==) transform(a,omega);

        else{

            transform(a,omegaInv);

            for(int i=;i<n;i++) a[i].x/=(double)n;

        }

    }

}fft;

int n,m,w;

CD A[N],B[N],C[N],ans[N];

int main(){

    freopen("in","r",stdin);

    n=read();

    for(int i=;i<=n;i++){

        w=read(); m=max(m,w);

        A[w].x=;

        B[w*].x=;

        C[w*].x=;

    }

    m=m*;

    fft.ini(m);

    fft.DFT(A,);fft.DFT(B,);fft.DFT(C,);

    for(int i=;i<fft.n;i++)

        ans[i]=ans[i]+A[i]+(A[i]*A[i]-B[i])/2.0+(A[i]*A[i]*A[i]-3.0*A[i]*B[i]+2.0*C[i])/6.0;

    fft.DFT(ans,-);

    for(int i=;i<m;i++) if(int(ans[i].x+0.5)) printf("%d %d\n",i,int(ans[i].x+0.5));

}

相关文章