Solution -「ZJOI 2014」力

Descrption

Link.

对于每一个 (i),求出:

[sum_{j=1}^{i-1}frac{a_{j}}{(i-j)^{2}}-sum_{j=i+1}^{n}frac{a_{j}}{(i-j)^{2}} ]

Solution

(f(i)=a_{i},g(i)=frac{1}{i^{2}})

然后

[sum_{j=1}^{i-1}f(j) imes g(i-j)-sum_{j=i+1}^{n}f(j) imes g(j-i) ]

可以 FFT 了。

#include<bits/stdc++.h>
using namespace std;
const double bh_pi=acos(-1);
int n;
namespace Poly
{
	typedef complex<double> comp;
	typedef vector<complex<double> > poly;
	#define len(x) (int((x).size()))
	int lim,rev[400010];
	void fft(poly &f,int op)
	{
		for(int i=0;i<lim;++i)	if(i<rev[i])	swap(f[i],f[rev[i]]);
		for(int len=2;len<=lim;len<<=1)
		{
			comp bas(cos(2*bh_pi/len),op*sin(2*bh_pi/len));
			for(int fr=0;fr<lim;fr+=len)
			{
				comp now(1,0);
				for(int ba=fr;ba<fr+(len>>1);++ba,now*=bas)
				{
					comp tmp=now*f[ba+(len>>1)];
					f[ba+(len>>1)]=f[ba]-tmp;
					f[ba]+=tmp;
				}
			}
		}
		if(op==-1)	for(int i=0;i<lim;++i)	f[i]/=lim;
	}
	poly mulPoly(poly f,poly g)
	{
		int n=len(f)+len(g)-1;
		for(lim=1;lim<=n;lim<<=1);
		for(int i=0;i<lim;++i)	rev[i]=(rev[i>>1]>>1)|((i&1)?(lim>>1):0);
		f.resize(lim),g.resize(lim);
		fft(f,1),fft(g,1);
		for(int i=0;i<lim;++i)	f[i]*=g[i];
		fft(f,-1),f.resize(n);
		return f;
	}
}using namespace Poly;
int main()
{
	poly f,g;
	scanf("%d",&n);
	f.resize(n+1),g.resize(n+1);
	for(int i=1;i<=n;++i)
	{
		double x;
		scanf("%lf",&x);
		f[i]=comp(x,0);
	}
	for(int i=1;i<=n;++i)	g[i]=comp(1.0/i/i,0);
	poly onetmp=mulPoly(f,g);
	reverse(next(f.begin()),f.end());
	poly anotmp=mulPoly(f,g);
	for(int i=1;i<=n;++i)	printf("%.3lf
",onetmp[i].real()-anotmp[n-i+1].real());
	return 0;
}
原文地址:https://www.cnblogs.com/orchid-any/p/14517246.html