【BZOJ】3527: [Zjoi2014]力（fft+卷积）
http://www.lydsy.com/JudgeOnline/problem.php?id=3527
好好的一道模板题，我自己被自己坑了好久。。
首先题目看错。。。。。。。什么玩意。。。。。。。首先题目要求：
$$F_j=\sum_{i<j} \frac{q_i q_j}{(i-j)^2}-\sum_{i>j} \frac{q_i q_j}{(i-j)^2}$$
然后设
$$E_i=\frac{F_i}{q_i}$$
然后我将$F_j$看成$F_i$。。。。作死系列。。。。。。
然后fft写错一个地方查错。。。。
然后第二次调用fft忘记清空。。。查错。。。（完全忘记虚数也要初始化啊！！！！！！！！！！！！！！！！！！！！！！！！！！！）
虽然说调用完dft-1时，虚数部分是0，但是你知道有精度这个东西。。。。。。。。。。。。。。。。。。。
一直以为是double精度的问题？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？？
原来本来就是精度问题。。。。。。。。。。。。。。。。。。。。。。。。。。。。
首先我们来分析$E_i$
展开可得：
$$E_i=\sum_{j<i} \frac{q_j}{(j-i)^2}-\sum_{j>i} \frac{q_j}{(j-i)^2}$$
（注意i和j都变了23333333333
然后再变化一下：
$$E_i=\sum_{j=1}^{i-1} \frac{q_j}{(j-i)^2} - \sum_{j=i+1}^{n} \frac{q_j}{(j-i)^2}$$
似乎不能看出什么？但是应该根据数据范围以及i和j，我们应该可以想想：卷积。。
然后设$f[i]=q_i, g[i]=\frac{1}{i^2}$，原式可以化为
$$E_i=\sum_{j=1}^{i-1} f[j]g[i-j] - \sum_{j=i+1}^{n} f[j]g[j-i]$$
发现卷积的形式是
$$c_i=\sum_{j=0}^{i} a_j b_{i-j}$$
似乎就是j的初值不同？
两种方法：一是将读入的$q_i$按0～n-1存储，二是像我这样做。
首先学习了《具体数学》中的和式级数，发现我是多么的sb。。
首先左边的卷积显然，上下界可以去掉，我们只来说右边：
$$\sum_{j=i+1}^{n} f[j]g[j-i]$$
其实应该这样写好理解。。
$$\sum_{i<j<=n} f[j]g[j-i]$$
然后直接换元（具体数学里称交换律。。
$$\sum_{i<j+i<=n} f[j+i]g[j+i-i]$$
然后再处理一下得
$$\sum_{0<=j<=n-i} f[j+i]g[j]$$
然后再处理一下。（上式=0是因为不影响，即g[0]=0）
$$\sum_{0<=n-i-j<=n-i} f[n-i-j+i]g[n-i-j]$$
变成
$$\sum_{0<=j<=n-i} f[n-j]g[n-i-j]$$
设$t=n-i$，则
$$\sum_{0<=j<=t} f[n-j]g[t-j]$$
设$f'[j]=f[n-j]$
则
$$\sum_{0<=j<=t} f'[j]g[t-j]$$
囧。。。一写下来发现也很麻烦？手推的时候很多步都可以省略吧。。。
//以下的作废，因为太复杂了。。。。。。233
首先j=0的话，我们可以将f[0]=0就行了嘛。。。i-1我们可以g[0]=0就好了嘛。同理j=i+1也能化成j=i
那么原式变成
$$E_i=\sum_{j=0}^{i} f[j]g[i-j] - \sum_{j=i}^{n} f[j]g[j-i]$$
发现右边还不行？特殊的技巧。。。减法。。
于是
$$E_i=\sum_{j=0}^{i} f[j]g[i-j] - \sum_{j=0}^{n-i} f[n-j]g[(n-j)-i]$$
设$ff[j]=f[n-j]$，那么
$$E_i=\sum_{j=0}^{i} f[j]g[i-j] - \sum_{j=0}^{n-i} ff[j]g[n-i-j]$$
然后设向量A和B，那么
$$A[i]=\sum_{j=0}^{i} f[j]g[i-j]$$
$$B[i]=\sum_{j=0}^{i} ff[j]g[i-j]$$
所以
$$E_i=A[i]-B[n-i]$$
//以上的作废
然后A和B用卷积算出来即可。。
还有一个要注意的地方，数组要开大点，就是(2*n+1)*2这样。。。。。。。。否则。。。。。。。。。。。。。。。。。。
（另，这是spj）
#include <cstdio>
#include <cstring>
#include <cmath>
#include <string>
#include <iostream>
#include <algorithm>
#include <queue>
#include <set>
#include <vector>
#include <map>
using namespace std;
typedef long long ll;
#define pii pair<int, int>
#define mkpii make_pair<int, int>
#define pdi pair<double, int>
#define mkpdi make_pair<double, int>
#define pli pair<ll, int>
#define mkpli make_pair<ll, int>
#define rep(i, n) for(int i=0; i<(n); ++i)
#define for1(i,a,n) for(int i=(a);i<=(n);++i)
#define for2(i,a,n) for(int i=(a);i<(n);++i)
#define for3(i,a,n) for(int i=(a);i>=(n);--i)
#define for4(i,a,n) for(int i=(a);i>(n);--i)
#define CC(i,a) memset(i,a,sizeof(i))
#define read(a) a=getint()
#define print(a) printf("%d", a)
#define dbg(x) cout << (#x) << " = " << (x) << endl
#define error(x) (!(x)?puts("error"):0)
#define printarr2(a, b, c) for1(_, 1, b) { for1(__, 1, c) cout << a[_][__]; cout << endl; }
#define printarr1(a, b) for1(_, 1, b) cout << a[_] << '\t'; cout << endl
inline const int getint() { int r=0, k=1; char c=getchar(); for(; c<'0'||c>'9'; c=getchar()) if(c=='-') k=-1; for(; c>='0'&&c<='9'; c=getchar()) r=r*10+c-'0'; return k*r; }
inline const int max(const int &a, const int &b) { return a>b?a:b; }
inline const int min(const int &a, const int &b) { return a<b?a:b; }
struct cp {
double r, i;
cp(long double _r=0.0, double _i=0.0) : r(_r), i(_i) { }
cp operator+ (const cp &x) { return cp(r + x.r, i + x.i); }
cp operator- (const cp &x) { return cp(r - x.r, i - x.i); }
cp operator* (const cp &x) { return cp(r*x.r - i*x.i, r*x.i + i*x.r); }
};
const int N=400005;
const double PI=acos(-1.0);
cp A[N];
int rev[N];
void dft(cp *a, int n, int flag) {
rep(i, n) A[rev[i]]=a[i];
rep(i, n) a[i]=A[i];
for(int m=2; m<=n; m<<=1) {
cp wn(cos(2.0*PI/m*flag), sin(2.0*PI/m*flag));
for(int i=0; i<n; i+=m) {
cp w(1.0); int mid=(m>>1);
for(int j=0; j<mid; ++j) {
cp u=a[i+j+mid]*w, t=a[i+j];
a[i+j]=t+u;
a[i+j+mid]=t-u;
w=w*wn;
}
}
}
if(flag==-1) rep(i, n) a[i].r/=n;
}
void init(int &len) {
int k=1, l=0;
while(k<len) k<<=1, ++l;
len=k;
rep(i, len) {
int t=i, r=0; k=l;
while(k--) r<<=1, r|=(t&1), t>>=1;
rev[i]=r;
}
}
void fft(double *x, double *y, cp *a, cp *b, cp *c, int len) {
rep(i, len) a[i].r=x[i], a[i].i=0.0; //虚数不初始化就是作死
rep(i, len) b[i].r=y[i], b[i].i=0.0;
dft(a, len, 1); dft(b, len, 1);
rep(i, len) c[i]=a[i]*b[i];
dft(c, len, -1);
}
int n, len;
double g[N], f[N], ff[N], ans[N];
cp a[N], b[N];
int main() {
read(n); len=n*2+1;
init(len);
for1(i, 1, n) scanf("%lf", &f[i]);
for1(i, 1, n) g[i]=1.0/i/i;
for1(i, 0, n) ff[i]=f[n-i];
fft(f, g, a, b, a, len);
for1(i, 1, n) ans[i]=a[i].r;
fft(ff, g, a, b, a, len);
for1(i, 1, n) ans[i]-=a[n-i].r;
for1(i, 1, n) printf("%.3f\n", ans[i]);
return 0;
}
Description
给出n个数qi，给出Fj的定义如下： 
令Ei=Fi/qi，求Ei.
Input
第一行一个整数n。
接下来n行每行输入一个数，第i行表示qi。
Output
n行，第i行输出Ei。
与标准答案误差不超过1e-2即可。
Sample Input
5
4006373.885184
15375036.435759
1717456.469144
8514941.004912
1410681.345880
Sample Output
-16838672.693
3439.793
7509018.566
4595686.886
10903040.872
Hint
对于30%的数据，n≤1000。
对于50%的数据，n≤60000。
对于100%的数据，n≤100000，0<qi<1000000000。
Source
感谢nodgd放题
