[ZJOI2014]力 题解
调了好长时间qwq
深刻的教训:
一个 m 次多项式和一个 n 次多项式乘起来是一个 m+n 次多项式而不是 n 次多项式……(sb错误啊啊啊啊)
题面
给定 \(n \le 10^5\), \(q_1, \cdots, q_n\)
(\(0 < q_i < 10^9\) 且为实数),
已知
\[E_i = \sum_{j=1}^{i-1}\frac{q_j}{(i-j)^2} - \sum_{j=i+1}^{n}\frac{q_j}{(j-i)^2}
\]
求
\[E_1, E_2, \cdots, E_n
\]
(当你的输出与标准答案相差不超过 10^{-2}10−2 时即被认为正确。)
首先将 \(E_i\) 拆成两项
\[A_i = \sum_{j=1}^{i-1}\frac{q_j}{(i-j)^2}
\]
\[B_i = \sum_{j=i+1}^{n}\frac{q_j}{(j-i)^2}
\]
则 \(E_i = A_i - B_i\) 。
数列 \(A\) 可以由数列 \(q\) 和数列 \(1, \frac{1}2, \frac{1}3, \cdots\) 卷积得到。
数列 \(B\) 呢?
\[B_i = \sum_{j=i+1}^{n}\frac{q_j}{(j-i)^2}
\]
\[=\sum_{j-i=1}^{n-i}\frac{q_j}{(j-i)^2}
\]
\[=\sum_{s=1}^{n-i}\frac{q_{s+i}}{s^2}
\]
山穷水尽?
看题解后,
学到许多。
一个小trick
有
\[\sum_{j=1}^{n-i} f[i+j] *g[j]
\]
设 \(f'[i] = f[n-i]\)
则
\[\sum_{j=1}^{n-i} f[i+j] *g[j]
\]
\[=\sum_{j=1}^{n-i} f'[n-(i+j)] *g[j]
\]
\[=\sum_{j=1}^{n-i} f'[n-i-j] *g[j]
\]
令 \(t = n-i\)
\[=\sum_{j=1}^{t} f'[t-j] *g[j]
\]
爽。
luogu数据AC代码
#include<bits/stdc++.h>
using namespace std;
const int maxn = 2e6+15;
const long double pi = acos(-1.0);
struct cp{
double x,y;
cp() {
x=y=0.0;
}
cp(double x, double y) {
this->x=x, this->y=y;
}
const cp operator+(const cp p)const {
return cp(x+p.x,y+p.y);
}
const cp operator-(const cp p)const {
return cp(x-p.x,y-p.y);
}
const cp operator*(const cp p)const {
return cp(x*p.x-y*p.y, x*p.y+y*p.x);
}
};
int rv[maxn];
void fft(cp *a, int n, int inv) {
for(int i=0;i<n;++i) if(i<=rv[i]) swap(a[i],a[rv[i]]);
for(int m=2;m<=n;m<<=1) {
cp w(cos(2*pi/m), inv*sin(2*pi/m));
for(int i=0;i<n;i+=m) {
cp tmp(1,0);
for(int j=0;j<(m>>1);++j) {
cp p=a[i+j], q=tmp*a[i+j+(m>>1)];
a[i+j]=p+q, a[i+j+(m>>1)]=p-q;
tmp=tmp*w;
}
}
}
if(inv==-1) {
for(int i=0;i<n;++i) a[i].x/=n;
}
}
int n,m;
long double q[maxn];
cp f[maxn], rf[maxn], g[maxn];
int main()
{
scanf("%d", &m);
for(int i=1;i<=m;++i) {
scanf("%lf", &f[i].x);
rf[m-i].x=f[i].x;
}
for(int i=1;i<=m;++i) g[i]=cp(1.0/i/i,0);
// for(int j=1;j<=m;++j) {
// double ans = 0.0;
// for(int i=1;i<=j-1;++i) ans+=q[j]*q[i]*g[j-i].x;
// for(int i=j+1;i<=m;++i) ans-=q[j]*q[i]*g[i-j].x;
// printf("%.3lf\n", ans/q[j]);
// } bf algorithm
for(n=1;n<(m<<1);n<<=1);
for(int i=0;i<n;++i) rv[i]=(rv[i>>1]>>1)|((i&1)?(n>>1):0);
fft(f,n,1), fft(rf,n,1), fft(g,n,1);
for(int i=0;i<n;++i) f[i]=f[i]*g[i];
for(int i=0;i<n;++i) rf[i]=rf[i]*g[i];
fft(f,n,-1), fft(rf,n,-1);
for(int i=1;i<=m;++i) printf("%.6lf\n", f[i].x-rf[m-i].x);
return 0;
}