[ZJOI2014]力[卷积,FFT]

最新推荐文章于 2023-08-28 23:43:59 发布

Hardict

最新推荐文章于 2023-08-28 23:43:59 发布

阅读量225

点赞数

分类专栏：多项式-FFT 多项式-卷积文章标签： FFT 卷积

本文链接：https://blog.youkuaiyun.com/Hardict/article/details/82631485

版权

多项式-FFT 同时被 2 个专栏收录

2 篇文章

订阅专栏

多项式-卷积

1 篇文章

订阅专栏

题目

给出个数 $q_i$ ，给出 F_j 的定义如下：

$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}$ ，求 E_i .

数据范围: $2\leq n \leq 100000, 0\leq q_i \leq 1000000000$

分析

先考虑 n=4 的情况,即有 q_1,q_2,q_3,q_4 四个数

$\\ F_1=0*q_1-\frac{1}{1^2}*q_2-\frac{1}{2^2}*q_3-\frac{1}{3^2}*q_4\\ F_2=\frac{1}{1^2}*q_1+0*q_2-\frac{1}{1^2}*q_3-\frac{1}{2^2}*q_4\\ F_3=\frac{1}{2^2}*q_1+\frac{1}{1^2}*q_2+0*q_3-\frac{1}{q^2}*q_4\\ F_4=\frac{1}{3^2}*q_1+\frac{1}{2^2}*q_2+\frac{1}{1^2}*q_3+0*q_4$

考虑序列 $f = \{ q_1,q_2,q_3,...,q_n \}$ , $g= \{ -\frac{1}{(n-1)^2}, -\frac{1}{(n-2)}^2,...,-\frac{1}{1^2},0,\frac{1}{1^2},...,\frac{1}{(n-1)^2} \}$ .考虑卷积 $(f\bigotimes g)(n) = \sum_{i=0}^{n}f(i)g(n-i)$ m即为上述序列相乘的形式,以序列 $f,g$ 分别为生成函数 F(x),G(x) 的前k项系数(其余项系数为0),答案即为 C(x)=F(x)*G(x) 的系数,以FFT优化即可.注意精度.

#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <complex>
using namespace std;

const int MAXN = 1e6 + 5;
const double PI = acos(-1.0);
int N;
double p[MAXN];
complex<double> F[MAXN], G[MAXN], C[MAXN];
void FFT(int len, complex<double> *a, int type);

int main(){
  ios::sync_with_stdio(false);
  cin >> N;
  int i;
  for(i = 1; i <= N; i++) cin >> p[i];
  for(i = 0; i <= N - 1; i++) F[i] = complex<double>(p[i + 1], 0);
  for(i = 0; i <= N - 2; i++) G[i] = complex<double>(-1.0 / ((N - 1.0 - i) * (N - 1.0 - i)), 0);
  for(i = N; i <= 2 * N - 2; i++) G[i] = complex<double>(1.0 / ((N - 1.0 - i) * (N - 1.0 - i)), 0);
  G[N - 1] = complex<double>(0, 0);
  int limit = 1;
  while(limit <= N + 2 * N - 1) limit <<= 1;
  FFT(limit, F, 1), FFT(limit, G, 1);
  for(i = 0; i <= limit; i++)
    C[i] = F[i] * G[i];
  FFT(limit, C, -1);
  for(i = N - 1; i <= 2 * N - 2; i++)
    printf("%.6lf\n", C[i].real() / limit);
  return 0;
}

void bit_reverse(int limit, complex<double> *a);
void FFT(int len, complex<double> *a, int type){
  int i, j, l;
  bit_reverse(len, a);
  for(l = 2; l <= len; l <<= 1){
    int mid = l >> 1;
    complex<double> wn = complex<double>(cos(2.0 * PI / l), type * sin(2.0 * PI / l));
    for(i = 0; i < len; i += l){
      complex<double> w = complex<double>(1, 0);
      for(j = 0; j < mid; j++, w = w * wn){
        complex<double> x = a[i + j], y = w * a[i + j + mid];
        a[i + j] = x + y;
        a[i + j + mid] = x - y;
      }
    }
  }
}

void bit_reverse(int limit, complex<double> *a){
  int i, j, l;
  for(i = 0, j = 0; i < limit; i++){
    if(i > j) swap(a[i], a[j]);
    for(l = limit >> 1; (j & l); j ^= l, l >>= 1);
    j ^= l;
  }
}