[普通型母函数+容斥+FFT]BZOJ 3771:Triple

看到ZZK大佬的blog写的太好啦，所以就“转”来了。

题目概述

CHNJZ有$n$把价值不一样的斧子，ZigZagK偷走了$1$把或$2$把或$3$把，对于每个可能的总损失，计算有几种可能的方案。

解题报告

emm……显然是母函数啊？但是有数量限制。由于最多偷走三把，所以我们可以直接三种情况都讨论过去。

$1$把：母函数：$A(x)=x^{a_1}+x^{a_2}+x^{a_3}+...+x^{a_n}$ 。

$2$把：母函数： ，但是这样显然算进了同把斧子用了两次的情况，所以要去掉$B(x)=x^{2a_1}+x^{2a_2}+x^{2a_3}+...+x^{2a_n}$。又因为这样算进了排列，所以最后要除以$2$。

$3$把：母函数：$A^3(x)$ ，去掉$3A(x)B(x)$，补回$2[C(x)=x^{3a_1}+x^{3a_2}+x^{3a_3}+...+x^{3a_n}]$，除以 6 <script type="math/tex" id="MathJax-Element-14">6</script>。

因为项数比较多，所以用FFT进行多项式乘法。

解题报告

#include<cmath>
#include<cstdio>
#include<algorithm>
using namespace std;
const int maxn=300005;
const double pi=acos(-1);
struct jz{
  double x,y;
  jz(double x=0,double y=0):x(x),y(y){}
  jz operator-(const jz &b)const{return jz(x-b.x,y-b.y);}
  jz operator+(const jz &b)const{return jz(x+b.x,y+b.y);}
  jz operator*(const jz &b)const{return jz(x*b.x-y*b.y,x*b.y+y*b.x);}
}ans[maxn],A[maxn],B[maxn],AB[maxn],A2[maxn],C[maxn],A3[maxn];
int n,m,re[maxn];
inline int _read(){
  int num=0;char ch=getchar();
  while(ch<'0'||ch>'9') ch=getchar();
  while(ch>='0'&&ch<='9') num=num*10+ch-48,ch=getchar();
  return num;
}
void FFT(jz *a,int f){
  for (int i=0;i<n;i++) if (i<re[i]) swap(a[i],a[re[i]]);
  for (int i=1;i<n;i<<=1){
    jz w(1,0),wn(cos(pi/i),sin(f*pi/i)),x,y;
    for (int j=0;j<n;j+=i<<1,w=jz(1,0))
      for (int k=0;k<i;k++,w=w*wn){
    x=a[j+k];y=w*a[j+k+i];
    a[j+k]=x+y;a[j+k+i]=x-y;
      }
  }
}
int main(){
  n=_read();
  for (int i=1;i<=n;i++){
    int x=_read();m=max(x,m);
    A[x]=B[2*x]=C[3*x]=jz(1,0);
  }m*=3;
  int l=0;for (n=1;n<=(m<<1);l++,n<<=1);
  for (int i=0;i<n;i++) re[i]=((re[i>>1]>>1)|((i&1)<<(l-1)));
  for (int i=0;i<=m;i++) ans[i].x=A[i].x-B[i].x/2+C[i].x/3+0.1;
  FFT(A,1);FFT(B,1);
  for (int i=0;i<n;i++) AB[i]=A[i]*B[i],A2[i]=A[i]*A[i],A3[i]=A[i]*A[i]*A[i];
  FFT(AB,-1);FFT(A2,-1);FFT(A3,-1);
  for (int i=0;i<=m;i++){
    int w=ans[i].x+(A2[i].x/2+A3[i].x/6-AB[i].x/2)/n+0.1;
    if (w) printf("%d %d\n",i,w);
  }
}