bzoj 3771: Triple 快速傅里叶变换 FFT

题目大意:

给出(n)个互不相同的物品,每个物品有价值(x_i(x_i leq 40000))如果可以从中取一个或两个或三个物品。问能够组合出来的所有价值和对应的方案数,全部输出.取值时,((a,b))((b,a))算作一种

题解:

不要说什没有给(n)的范围什么的,题目表明了(n leq 40000)
我们首先想到的是直接取
我们构造出只取一个物品时的母函数

[A(x) = x + x^2 + x^3 + ... ]

然后我们取两个物品的时候的母函数就可以用乘积来表示
(A(x)*A(x))
但是我们发现这样取到了不合法的情况,每个物品只能取一遍,但是直接计算的话会出现同一个物品取了两次的情况
所以我们应该减去这个不合法的情况(即容斥原理)
什么是不合法情况呢?就是同一个物品被取到了两次.
所以我们构造出取同一个物品取了两次的母函数,即

[B(x) = x^2 + x^4 + x^6 + ... ]

所以我们取两个的方案数就是(A(x)^2 - B(x))
相应的我们构造出同一个物品取了三次的母函数,即

[C(x) = x^3 + x^6 + x^9 + ... ]

所以我们考虑在取三个的时候我们取到的有(AAB, ABA, BAA, AAA)
所以由容斥原理,取三个的时候的方案即为(A(x)^3 - 3A(x)B(x) + 2C(x))
但是我们都没有考虑顺序的问题,我们考虑顺去,除去所有元素的全排列
那么我们最终的答案就是

[frac{A(x)}{1!} + frac{A(x)^2-B(x)}{2!} + frac{A(x)^3 - 3A(x)B(x) + 2C(x)}{3!} ]

使用FFT计算即可

#include <cmath>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;
typedef long long ll;
inline void read(int &x){
	x=0;char ch;bool flag = false;
	while(ch=getchar(),ch<'!');if(ch == '-') ch=getchar(),flag = true;
	while(x=10*x+ch-'0',ch=getchar(),ch>'!');if(flag) x=-x;
}
inline int cat_max(const int &a,const int &b){return a>b ? a:b;}
inline int cat_min(const int &a,const int &b){return a<b ? a:b;}
const int maxn = 208000;
const double pi = acos(-1);
struct complex{
	double x,y;
	complex(){}
	complex(double a,double b){x=a;y=b;}
	complex operator + (const complex &r){return complex(x+r.x,y+r.y);}
	complex operator - (const complex &r){return complex(x-r.x,y-r.y);}
	complex operator * (const complex &r){return complex(x*r.x-y*r.y,x*r.y+y*r.x);}
	complex operator / (const double &r){return complex(x/r,y/r);}
};
void FFT(complex *x,int n,int p){
	for(int i=0,t=0;i<n;++i){
		if(i > t) swap(x[i],x[t]);
		for(int j=n>>1;(t^=j) < j;j>>=1);
	}
	for(int m=2;m<=n;m<<=1){
		complex wn(cos(p*2*pi/m),sin(p*2*pi/m));
		for(int i=0;i<n;i+=m){
			complex w(1,0),u;
			int k = m>>1;
			for(int j=0;j<k;++j,w=w*wn){
				u = x[i+j+k]*w;
				x[i+j+k] = x[i+j] - u;
				x[i+j] = x[i+j] + u;
			}
		}
	}
	if(p == -1) for(int i=0;i<n;++i) x[i] = x[i]/n;
}
complex a[maxn],b[maxn],c[maxn],d[maxn];
int main(){
	int n;read(n);
	int nn = 0;
	for(int i=0,x;i<n;++i){
		read(x);
		a[x].x += 1.0;
		b[x<<1].x += 1.0;
		c[x*3].x += 1.0;
		nn = cat_max(nn,x*3);
	}
	int len = 1;
	for(int i=1;(i>>1)<nn;i<<=1) len = i;
	FFT(a,len,1);FFT(b,len,1);
	for(int i=0;i<len;++i) d[i] = a[i]*a[i]*a[i]/6.0 + (a[i]*a[i] - a[i]*b[i] - b[i])/2.0 + a[i];
	FFT(d,len,-1);
	for(int i=0;i<len;++i) d[i] = d[i] + c[i]/3.0;
	for(int i=0;i<len;++i){
		int x = (int)(d[i].x + 0.5);
		if(x == 0) continue;
		printf("%d %d
",i,x);
	}
	getchar();getchar();
	return 0;
}
  
原文地址:https://www.cnblogs.com/Skyminer/p/6357342.html