ZJOI2014 力 推式子 FFT板题 VJ1000题祭

ZJOI2014 力 推式子 FFT板题 VJ1000题祭

题意

定义

[F_j = sum_{i = 1}^{j - 1} frac{q_i imes q_j} {(i - j)^2} - sum_{i = j+ 1}^n frac{q_i imes q_j}{(i - j)^2}\ E_i = F_i / q_i ]

对于(1 leq i leq n),求出(E_i)

[1 leq n leq 10^5\ 0 < q_i < 10^9 ]

分析

直接推式子，我们希望推成卷积的形式

[E_j = frac{F_j}{q_j}\ = sum_{i = 1}^{j} frac{q_i}{(i - j)^2} - sum_{i=j}^n frac{q_i}{(i - j)^2}\ f_i = q_i,g_i = frac{1}{i^2}\ = sum_{i=1}^j f_i g_{j-i} - sum_{i=j}^n f_i g_{i-j}\ 约定f_0 = g_0 = 0\ = sum_{i = 0}^i f_i g_{j - i} - sum_{i=0}^t f'_{t-i}g_i ]

于是对左右卷积即可

代码

#include<bits/stdc++.h>
#define re register
#define pii pair<int,int>
#define fi first
#define se second 
using namespace std;
typedef long long ll;

const double PI = acos(-1.0);
const int maxn = 4e5 + 5;

ll rd(){
	ll x = 0;
	int f = 1;
	char ch = getchar();
	while(ch < '0' || ch > '9') {
		if(ch == '-') f = -1;
		ch = getchar();
	} 
	while(ch >= '0' && ch <= '9') {
		x = x * 10 + ch - '0';
		ch = getchar();
	}
	return x * f;
}

int n,m;

struct CP{
	CP(double xx = 0,double yy = 0) {x = xx,y = yy;}
	double x,y;
	CP operator + (CP const& B) {return CP(x + B.x,y + B.y);}
	CP operator - (CP const& B) {return CP(x - B.x,y - B.y);}
	CP operator * (CP const& B) {return CP(x * B.x - y * B.y,x * B.y + y * B.x);}
}f[maxn << 1],p[maxn << 1],c[maxn << 1];

int tr[maxn << 1];

void fft(CP *f,bool flag){
	for(int i = 0;i < n;i++)
		if(i < tr[i]) swap(f[i],f[tr[i]]);
	for(int p = 2;p <= n;p <<= 1) {
		int len = p >> 1;
		CP tG(cos(2 * PI / p),sin(2 * PI / p));
		if(!flag) tG.y *= -1;
		for(int k = 0;k < n;k += p) {
			CP buf(1,0);
			for(int l = k;l < k + len;l++){
				CP tt = buf * f[len + l];
				f[len + l] = f[l] - tt;
				f[l] = f[l] + tt;
				buf = buf * tG;
			}
		}
	}
}

int main(){
	n = rd();
	int nn = n;
	m = n;
	for(int i = 1;i <= n;i++)
		scanf("%lf",&f[i].x),c[n - i].x = f[i].x,p[i].x = (double)(1.0 / i / i); 
	for(m += n,n = 1;n <= m;n <<= 1);
	for(int i = 0;i < n;i++)
		tr[i] = (tr[i >> 1] >> 1) | ((i & 1) ? n >> 1 : 0);
	fft(f,1);
	fft(p,1);
	fft(c,1);
	for(int i = 0;i < n;i++)
		f[i] = f[i] * p[i],c[i] = c[i] * p[i];
	fft(f,0);
	fft(c,0);
	for(int i = 1;i <= nn;i++){
		double res = (f[i].x - c[nn - i].x ) *1.0 / n;
		printf("%.3f
",res);
	}
}