【题解】[ZJOI2014] 力 FFT 洛谷 P3338

[ZJOI2014] 力

题意

给出 $n$ 个数 $q_i$,
定义
$$ F_j = \sum_{i=1}^{j-1} \frac{q_i \times q_j}{(i-j)^2} – \sum_{i=j+1}^{n} \frac{q_i \times q_j}{(i-j)^2} $$
$$ E_i = \frac{F_i}{q_i} $$

求 $E_i$ ($i \in [1,n]$).

思路

首先, 提公因式, 化简可以得到,
$$ E_j = \sum_{i=1}^{j-1} \frac{q_i}{(j-i)^2} – \sum_{i=j+1}^{n} \frac{q_i}{(j-i)^2} $$ (把 $q_j$ 提出来后约掉)

注意到式子可以分为两个部分, 这两个部分中的因子可以分成 $i$ 和 $j-i$ 两个部分, 即这两个部分的编号和为 $j$, 那么我们就找到加法卷积了.

设 $F[i] = \frac{1}{i^2}$, $G[i] = q_i$.

则
$$ E_j = \sum_{i=1}^{j-1} F[j-i] \times G[i] – \sum_{i=j+1}^{n} F[j-i] \times G[i] $$

这样, 我们就可以分两步愉快地 FFT 了.

需要注意的是, 后半部分的 $j-i$ 为负数, 所以需要将整个数组往前移动 $n$ 位, 最后统计答案的时候也要往前移动 $n$ 位 (没听懂的话可以看一下代码).

代码

#include<bits/stdc++.h>
#define _USE_MATH_DEFINES
#define ll long long
#define db double
using namespace std;
const int _=1e5+7;
const int __=3e5+7;
const db Pi=M_PI;
struct cn{
  db a,b;
  cn operator + (const cn &x) const {
    return (cn){a+x.a,b+x.b};
  }
  cn operator - (const cn &x) const {
    return (cn){a-x.a,b-x.b};
  }
  cn operator * (const cn &x) const {
    return (cn){a*x.a-b*x.b,a*x.b+b*x.a};
  }
};
int N,n,num[__];
db q[_],E[_];
cn f[__],g[__];
void FFT(cn *f,int id){
  for(int i=0;i<n;i++)
    if(i<num[i]) swap(f[i],f[num[i]]);
  for(int len=2;len<=n;len<<=1){
    int gap=len>>1;
    cn w1=(cn){cos(2*Pi/len),sin(2*Pi/len)*id};
    for(int k=0;k<n;k+=len){
      cn w=(cn){1,0};
      for(int i=k;i<k+gap;i++,w=w*w1){
    cn tmp=w*f[i+gap];
    f[i+gap]=f[i]-tmp;
    f[i]=f[i]+tmp;
      }
    }
  }
}
int main(){
  //freopen("force.in","r",stdin);
  cin>>N;
  for(int i=1;i<=N;i++) scanf("%lf",&q[i]);
  n=1; while(n<=2*N) n<<=1;
  for(ll i=1;i<=N;i++){    // 这里需要开 long long
    f[i].a=(db)1ll/(i*i);
    g[i].a=q[i];
  }
  for(int i=0;i<n;i++)
    num[i]=(num[i>>1]>>1)|((i&1) ?n>>1 :0);
  FFT(f,1);
  FFT(g,1);
  for(int i=0;i<n;i++) f[i]=f[i]*g[i];
  FFT(f,-1);
  for(int i=1;i<=N;i++) E[i]=f[i].a/n;
  memset(f,0,sizeof(f));
  memset(g,0,sizeof(g));
  for(ll i=0;i<=N;i++){
    f[i].a=(db)1ll/((i-N)*(i-N)); // 往前移动 N 位
    g[i].a=q[i];
  }
  f[N].a=0;
  FFT(f,1);
  FFT(g,1);
  for(int i=0;i<n;i++) f[i]=f[i]*g[i];
  FFT(f,-1);
  for(int i=1;i<=N;i++) printf("%.3lf\n",E[i]-f[i+N].a/n);  // 统计答案时也往前移动 N 位
  return 0;
}