【bzoj3771】Triple 【FFT 生成函数 容斥原理】

Triple

Description

我们讲一个悲伤的故事。
从前有一个贫穷的樵夫在河边砍柴。
这时候河里出现了一个水神,夺过了他的斧头,说:
“这把斧头,是不是你的?”
樵夫一看:“是啊是啊!”
水神把斧头扔在一边,又拿起一个东西问:
“这把斧头,是不是你的?”
樵夫看不清楚,但又怕真的是自己的斧头,只好又答:“是啊是啊!”
水神又把手上的东西扔在一边,拿起第三个东西问:
“这把斧头,是不是你的?”
樵夫还是看不清楚,但是他觉得再这样下去他就没法砍柴了。
于是他又一次答:“是啊是啊!真的是!”
水神看着他,哈哈大笑道:
“你看看你现在的样子,真是丑陋!”
之后就消失了。
樵夫觉得很坑爹,他今天不仅没有砍到柴,还丢了一把斧头给那个水神。
于是他准备回家换一把斧头。
回家之后他才发现真正坑爹的事情才刚开始。
水神拿着的的确是他的斧头。
但是不一定是他拿出去的那把,还有可能是水神不知道怎么偷偷从他家里拿走的。
换句话说,水神可能拿走了他的一把,两把或者三把斧头。
樵夫觉得今天真是倒霉透了,但不管怎么样日子还得过。
他想统计他的损失。
樵夫的每一把斧头都有一个价值,不同斧头的价值不同。总损失就是丢掉的斧头价值和。
他想对于每个可能的总损失,计算有几种可能的方案。
注意:如果水神拿走了两把斧头a和b,(a,b)和(b,a)视为一种方案。拿走三把斧头时,(a,b,c),(b,c,a),(c,a,b),(c,b,a),(b,a,c),(a,c,b)视为一种方案。

Input

第一行是整数N,表示有N把斧头。
接下来n行升序输入N个数字Ai,表示每把斧头的价值。

Output

若干行,按升序对于所有可能的总损失输出一行x y,x为损失值,y为方案数。

HINT

所有数据满足:Ai<=40000
题解:我们搞3个多项式A,B,C,分别代表每个物品选一个、选两个、选三个。
(这玩意好像叫做生成函数)
设当前要求价值为i的答案数。
选一把斧头的情况:答案为A[i]。
选两把斧头的情况:答案为(A[i]A[i]B[i])/2。A[i]*A[i]中包含了选相同的,所以要减去B[i],还要除掉排列数。
选三把斧头的情况:答案为(A[i]A[i]A[i]3A[i]B[i]+2C[i])/6。具体就是A[i]A[i]A[i]包含了一种选三个一样的和三种选两个 一样的,A*B包含了一种选三个一样的和一种选两个一样的。(A[i]A[i]A[i]3A[i]B[i])多减了两种选三个一样的,所以再加上两个C[i]。最后除掉排列数。手动容斥理解一下就好了。
具体见代码。

#include<cstdio>
#include<algorithm>
#include<cmath>
using namespace std;
const int N=270005;
const double pi=3.141592653589793;
int n,m,x,rev[N];
long long ans;
struct complex{
    double x,y;
    complex(){
        x=y=0;
    }
    complex(double x,double y):x(x),y(y){}
    friend complex operator + (const complex &a,const complex &b){
        return complex(a.x+b.x,a.y+b.y);
    }
    friend complex operator - (const complex &a,const complex &b){
        return complex(a.x-b.x,a.y-b.y);
    }
    friend complex operator * (const complex &a,const complex &b){
        return complex(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);
    }
    friend complex operator * (const complex &a,const double &b){
        return complex(a.x*b,a.y*b);
    }
    friend complex operator / (const complex &a,const double &b){
        return complex(a.x/b,a.y/b);
    }
}a[N],b[N],c[N];
void fft(complex *a,int dft){
    for(int i=0;i<n;i++){
        if(i<rev[i]){
            swap(a[i],a[rev[i]]);
        }
    }
    for(int i=1;i<n;i<<=1){
        complex wn=complex(cos(pi/i),dft*sin(pi/i));
        for(int j=0;j<n;j+=i<<1){
            complex w=complex(1,0),x,y;
            for(int k=j;k<j+i;k++,w=w*wn){
                x=a[k];
                y=w*a[k+i];
                a[k]=x+y;
                a[k+i]=x-y;
            }
        }
    }
    if(dft==-1){
        for(int i=0;i<n;i++){
            a[i]=a[i]/n;
        }
    }
}
int main(){
    scanf("%d",&n);
    for(int i=1;i<=n;i++){
        scanf("%d",&x);
        a[x].x++;
        b[2*x].x++;
        c[3*x].x++;
        m=max(m,3*x);
    }
    for(n=1;n<=m;n<<=1);
    for(int i=0;i<n;i++){
        rev[i]=(rev[i>>1]>>1)|((i&1)*(n>>1));
    }
    fft(a,1);
    fft(b,1);
    fft(c,1);
    for(int i=0;i<n;i++){
        a[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;
    }
    fft(a,-1);
    for(int i=0;i<n;i++){
        ans=(long long)(a[i].x+0.5);
        if(ans){
            printf("%d %lld\n",i,ans);
        }
    }
    return 0;
}
posted @ 2018-04-19 22:26  ez_2016gdgzoi471  阅读(85)  评论(0编辑  收藏  举报