Spoj 8372 Triple Sums

传送门：http://www.spoj.com/problems/TSUM/

思路：先不管i<j<k这个条件

构造一个多项式A(x)=sigma x^a[i]

那么S=a[i]+a[j]+a[k]的个数就是x^(a[i]+a[j]+a[k]=S)的个数

为了让指数相加，我们把A(x)进行立方

那么个数就是A(x)^3中S次项的系数

然后进行容斥,为了方便，以下省去指数a[i]

A(x)^3=∑x*x*x(就是三个数都相等)+3∑x*x*y(就是两个数相等，系数是3是因为有3种排列方式)−6∑xyz(就是3个数都不等，系数是6同理)

现在我们要得到的是∑xyz

其中3∑x*x*y=∑x*x*∑x(有两个以上相同的)-∑x*x*x（三个都相同的）

代入原式变形

∑xyz=(A(x)^3-3*sigma(x*x)*sigma(x)+2*sigma(x*x*x))/6 （符号突然粘不进来了...）

然后就上FFT就可以了

递归太慢。。。。

#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define ll long long
const int maxn=1<<17;
const double pi=3.1415926535897932384626433832795;
using namespace std;
struct plex{double r,i;}tmp[maxn];

plex operator +(plex a,plex b){return (plex){a.r+b.r,a.i+b.i};}
plex operator -(plex a,plex b){return (plex){a.r-b.r,a.i-b.i};}
plex operator *(plex a,plex b){return (plex){a.r*b.r-a.i*b.i,a.r*b.i+a.i*b.r};}

struct DFT{
	plex a[maxn];
	void fft(int bg,int step,int size,int op){
		if (size==1) return;
		fft(bg,step<<1,size>>1,op),fft(bg+step,step<<1,size>>1,op);
		plex w=(plex){1,0},t=(plex){cos(2.0*pi/size),sin(2.0*pi*op/size)};
		int p=bg,p0=bg,p1=bg+step;
		for (int i=0;i<size/2;i++){
			tmp[p]=a[p0]+w*a[p1];
			tmp[p+size/2*step]=a[p0]-w*a[p1];
			p+=step,p0+=(step<<1),p1+=(step<<1),w=w*t;
		}
		for (int i=bg;size;size--,i+=step) a[i]=tmp[i];
	}
}A,B,C;
int n,a[maxn],b[maxn],c[maxn];

int main(){
	scanf("%d",&n);
	for (int i=0,x;i<n;i++)
		scanf("%d",&x),a[x+20000]++,b[x+x+40000]++,c[3*x+60000]++;
	for (int i=0;i<maxn;i++) A.a[i].r=a[i],B.a[i].r=b[i],C.a[i].r=c[i];
	A.fft(0,1,maxn,1),B.fft(0,1,maxn,1);
	for (int i=0;i<maxn;i++) C.a[i]=A.a[i]*(A.a[i]*A.a[i]-(plex){3.0,0.0}*B.a[i]);
	C.fft(0,1,maxn,-1);
	for (int i=0;i<maxn;i++){
		ll ans=((ll)(C.a[i].r/maxn+0.5)+2*c[i])/6;
		if (ans) printf("%d : %I64d
",i-60000,ans);
	}
	return 0;
}

(∑x)3−3(∑x2)(∑x)+2∑x36