ZJOI 2014力

Posted 2022-09-13 luoshuitianyi

Problem

Description

给出 \(n\) 个数 \(q_i\)，给出 \(F_j\) 的定义如下：

\[F_j=\sum_i<j \fracq_iq_j(i-j)^2 - \sum_i>j \fracq_iq_j(i-j)^2\]

令 \(E_i=F_i/q_i\)，求 \(E_i\)。

Input Format

第一行一个整数\(n\)。

接下来 \(n\) 行每行输入一个数，第 \(i\) 行表示 \(q_i\)。

Output Format

\(n\) 行，第 \(i\) 行输出 \(E_i\)。与标准答案误差不超过 \(10^-2\) 即可。

Sample

Input

5
4006373.885184
15375036.435759
1717456.469144
8514941.004912
1410681.345880

Output

-16838672.693
3439.793
7509018.566
4595686.886
10903040.872

Range

对于所有的数据，\(n\leq 100000,\ 0<q_i<10^9\)。

Algorithm

\(FFT\)

Mentality

\[ E_i=\sum_j=1^i-1\fracq_j(i-j)^2-\sum_j=i+1^n\fracq_j(j-i)^2 \]

设 \(g_i=\frac1i^2\) ，则有：
\[ E_i=\sum_j=1^i-1q_jg_i-j-\sum_j=i+1^nq_jg_j-i \]
单独算每个 \(E_i\) 的前一部分，再算所有 \(E_i\) 的后一部分。

由于 \(q_0=g_0=0\) ，则对于：
\[ \sum_j=1^i-1q_jg_i-j=\sum_j=0^iq_jg_i-j \]
直接 \(FFT\) 求得。

对于后一部分，将 \(q\) 数组翻转得到数组 \(p\) 。

则有：
\[ \sum_j=i+1^nq_jg_j-i=\sum_j=0^n-i+1p_jg_n-i+1-j \]
同样直接 \(FFT\) 计算。

Code

#include <algorithm>
#include <cmath>
#include <complex>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <map>
#include <queue>
#include <set>
#include <vector>
using namespace std;
#define cp complex<double>
const int Max_n = 3e5 + 5;
const double pi = acos(-1);
int n;
int rev[Max_n];
double dv[Max_n], ans[Max_n], q[Max_n];
cp f[Max_n], g[Max_n];
namespace FFT 
int lim, bit;
void dft(cp *f, int t) 
  for (int i = 0; i < lim; i++)
    if (rev[i] > i) swap(f[rev[i]], f[i]);
  for (int len = 1; len < lim; len <<= 1) 
    cp Wn = exp(cp(0, t * pi / len));
    for (int i = 0; i < lim; i += len << 1) 
      cp Wnk(1, 0);
      for (int k = i; k < i + len; k++, Wnk *= Wn) 
        cp x = f[k], y = Wnk * f[k + len];
        f[k] = x + y, f[k + len] = x - y;
      
    
  

void fft(double *a, double *b, int tot) 
  lim = 1, bit = 0;
  while (lim <= tot) lim <<= 1, bit++;
  for (int i = 0; i < lim; i++)
    rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (bit - 1));
  for (int i = 0; i < lim; i++) f[i] = a[i], g[i] = b[i];
  dft(f, 1), dft(g, 1);
  for (int i = 0; i < lim; i++) f[i] *= g[i];
  dft(f, -1);
  for (int i = 0; i < lim; i++) f[i] /= lim;

  // namespace FFT
using namespace FFT;
int main() 
#ifndef ONLINE_JUDGE
  freopen("3338.in", "r", stdin);
  freopen("3338.out", "w", stdout);
#endif
  cin >> n;
  for (int i = 1; i <= n; i++) scanf("%lf", &q[i]);
  for (int i = 1; i <= n; i++) dv[i] = (double)(1.0 / i / i);
  fft(q, dv, n << 1);
  for (int i = 1; i <= n; i++) ans[i] += f[i].real();
  reverse(q + 1, q + n + 1);
  fft(q, dv, n << 1);
  for (int i = 1; i <= n; i++)
    printf("%.2lf\n", (ans[i] -= f[n - i + 1].real()));