题目大意：给定长为\(n\)的数列\(p\)，\(F_j=\sum_{i=1}^{j-1}\frac{q_iq_j}{(i-j)^2}-\sum_{i=j+1}^{n}\frac{q_iq_j}{(i-j)^2}\)，对于\(i\in[1,n]\)，求\(E_i=\frac{F_i}{q_i}\)
FFT、卷积
分析：
\[F_j=\sum_{i=1}^{j-1}\frac{q_iq_j}{(i-j)^2}-\sum_{i=j+1}^{n}\frac{q_iq_j}{(i-j)^2}
\]
首先最后要除\(q_i\)，那么我们直接约去
\[E_j=\sum_{i=1}^{j-1}\frac{q_i}{(i-j)^2}-\sum_{i=j+1}^{n}\frac{q_i}{(i-j)^2}
\]
正着刚这个式子十分困难，我们可以从侧面，考虑每一个\(q_i\)会对哪些位置的答案产生怎样的贡献
不难发现，对于式子的前半部分，\(q_i\)会对\(q_{i+k}\)产生\(\frac{q_i}{k^2}\)的贡献
把负号移进来，\(q_i\)会对\(q_{i-k}\)产生\(\frac{-q_i}{k^2}\)的贡献
观察发现，这个东西与卷积十分相似
假设我们有
\(A[i]=q_i\)
\(B[i]=\begin{cases}0 \quad i = 0 \\ \frac{1}{i^2} \quad i > 0,i\in[1,n] \\ \frac{-1}{i^2} \quad i < 0,i\in[-n,-1]\end{cases}\)
那么\(A,B\)做个卷积就可以得到答案
下标没法为负数我们直接平移一下\(B\)就可以了
卷积直接\(FFT\)
#include <iostream>
#include <iomanip>
#include <cmath>
#include <algorithm>
using namespace std;
const int maxn = 4e5 + 100;
typedef double type;
const type pi = acos(-1);
struct complex{
double x,y;
complex operator + (const complex &rhs)const{return complex{x + rhs.x,y + rhs.y};}
complex operator - (const complex &rhs)const{return complex{x - rhs.x,y - rhs.y};}
complex operator * (const complex &rhs)const{return complex{x * rhs.x - y * rhs.y,x * rhs.y + y * rhs.x};}
}A[maxn << 1],B[maxn << 1];
type q[maxn];
int n,N,tr[maxn << 1];
inline void fft(complex *f,int 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;
const complex unit{cos(2 * pi / p),sin(2 * pi / p) * flag};
for(int k = 0;k < N;k += p){
complex g{1,0};
for(int l = k;l < k + len;l++){
complex t = f[l + len] * g;
f[l + len] = f[l] - t;
f[l] = f[l] + t;
g = g * unit;
}
}
}
}
int main(){
ios::sync_with_stdio(false);
cin >> n;
for(int i = 1;i <= n;i++)cin >> A[i].x,q[i] = A[i].x;
for(int i = 1;i <= n;i++){
B[n + 1 - i].x = (type(-1) / i) / i;
B[n + 1 + i].x = (type(1) / i) / i;
}
for(N = 1;N <= 3 * n + 1;N <<= 1);
for(int i = 0;i < N;i++)
tr[i] = (tr[i >> 1] >> 1) | ((i & 1) ? (N >> 1) : 0);
fft(A,1),fft(B,1);
for(int i = 0;i < N;i++)A[i] = A[i] * B[i];
fft(A,-1);
for(int i = 0;i < N;i++)A[i].x /= type(N);
for(int i = n + 2;i <= 2 * n + 1;i++)cout << fixed << A[i].x << '\n';
return 0;
}