题目地址

洛谷P3338

Solution

第一道FFT的应用AC祭!

我们要求:

\[E_j=\frac{F_j}{q_j}=\sum_{i<j}\frac{q_i}{(i-j)^2}-\sum_{i>j}\frac{q_i}{(i-j)^2}
\]

(\(q_j\) 直接在除法的时候消掉了qwq)

  • Step 0 卷积是什么?

首先我们要有明确的目标,我们要把上面的式子推成卷积的形式,我们就要来回顾一下卷积是什么。卷积的形式如下:

\[C_k=\sum_{i=0}^{k}A_i*B_{k-i}
\]

  • Step 1 直接推式子

有了目标,我们就好来推式子了(推式子真好玩),下面给出推理的重要步骤,尽量没有繁琐的步骤,读者可以自己思考一下。

\[E_j=\sum_{i<j}\frac{q_i}{(i-j)^2}-\sum_{i>j}\frac{q_i}{(i-j)^2}
\]

改变一下 \(\sum\) 的上下标表示形式,原式变成

\[\sum_{i=1}^{j-1}\frac{q_i}{(i-j)^2}-\sum_{i=j+1}^{n}\frac{q_i}{(i-j)^2}
\]

如果 \(i=j\) 也累加进去,对答案不影响,所以式子变成

\[\sum_{i=1}^{j}\frac{q_i}{(i-j)^2}-\sum_{i=j}^{n}\frac{q_i}{(i-j)^2}
\]

  • Step 2 转化

设 \(f[i]=q_i,g[i]=\frac{1}{i^2}\) ,所以原式变成

\[\sum_{i=1}^{j}f[i]*g[j-i]-\sum_{i=j}^{n}f[i]*g[i-j]
\]

令 \(f[0]=0,g[0]=0\),则原式变成

\[\sum_{i=0}^{j}f[i]*g[j-i]-\sum_{i=j}^{n}f[i]*g[i-j]
\]

这时我们发现,左边已经是一个卷积的形式,所以我们直接来推右边

将 \(\sum_{i=j}^{n}f[i]*g[i-j]\) 展开,发现:

\[f[j]*g[0]+f[j+1]*g[1]+...+f[j+(n-j)]*g[n-j]
\]

所以我们可以将原式写成

\[\sum_{i=0}^{n-j}f[j+i]*g[i]
\]

  • Step 3 继续转化

引入 \(f'[i]=f[n-i]\) ,实际上这是一种翻转的套路。则原式可写为

\[\sum_{i=0}^{n-j}f'[n-(j+i)]*g[i]
\]

\[\sum_{i=0}^{n-j}f'[n-j-i]*g[i]
\]

令 \(t = n-j\),则原式等于

\[\sum_{i=0}^{t}f'[t-i]*g[i]
\]

至此,我们成功的把两个式子化成了卷积!!!

总结一下:

\[E_j=\sum_{i=0}^{j}f[i]*g[j-i]-\sum_{i=0}^{t}f'[t-i)]*g[i]
\]

- Step 4 \(FFT\)加速卷积

设有多项式 \(A(x)=\sum_{i=0}^{n}f[i]\) , \(B(x)=\sum_{i=0}^{n}g[i]\) , \(C(x)=\sum_{i=0}^{n}f'[i]\),

我们令 \(L(x)=A(x)*B(x)\) , \(R(x)=C(x)*B(x)\)

所以 \(Ans[i]=l_i - r_{n-i}\) (\(l_i,r_i\) 分别是多项式 \(L(x)\) 和 \(R(x)\) \(x^i\) 的系数)

- Step 5 读到这里,你和暴力选手还没有差别 (逃

世界上卡你精度的办法有千千万万种。 ----By me

注意 \(g[i]=\frac{1}{i^2}\)

楼主在处理这里的时候写的是 1.0/(i*i*1.0),交上去只有 \(30pts\)。

在题解中看到大佬们处理这里时用的 (double)(1.0 / i / i),改一下后,\(100pts\)。

如果大家对处理精度问题有什么独特的见解,记得和我分享分享QwQ

Code

Talk is cheap.Show me the code.

#include<bits/stdc++.h>
using namespace std;
inline int read() {
int x=0,f=1; char ch=getchar();
while(ch<'0' || ch>'9') { if(ch=='-') f=-1; ch=getchar(); }
while(ch>='0'&&ch<='9') { x=(x<<3)+(x<<1)+(ch^48); ch=getchar(); }
return x * f;
}
const int N = 400007;
const double pi = acos(-1.0);
int n;
int rev[N];
struct CP {
double x,y;
CP operator + (CP el) { return (CP)<%x+el.x , y+el.y%>; }
CP operator - (CP el) { return (CP)<%x-el.x , y-el.y%>; }
CP operator * (CP el) { return (CP)<%x*el.x-y*el.y , x*el.y+y*el.x%>; }
}a[N],b[N],c[N];
void FFT(CP *A,int n,int flag) {
for(int i=0;i<n;++i) if(i < rev[i]) swap(A[i],A[rev[i]]);
for(int mid=1;mid<n;mid<<=1) {
CP Wn = (CP)<%cos(2*pi/(mid<<1)) , flag*sin(2*pi/(mid<<1))%>;
for(int i=0;i<n;i+=(mid<<1)) {
CP W = (CP)<%1 , 0%>;
for(int j=0;j<mid;++j,W=(W*Wn)) {
CP tmp0 = A[i+j], tmp1 = W*A[i+mid+j];
A[i+j] = tmp0 + tmp1;
A[i+mid+j] = tmp0 - tmp1;
}
}
}
if(flag == -1) {
for(int i=0;i<n;++i) A[i].x /= n;
}
}
int main()
{
//freopen("Li.txt","r",stdin);
//freopen("My.out","w",stdout);
n = read();
for(int i=1;i<=n;++i) {
scanf("%lf",&a[i].x);
c[n-i].x = a[i].x;
b[i].x = (double)(1.0 / i / i);
}
int lim = 1, L = 0; while(lim <= (n<<1)) lim <<= 1, ++L;
for(int i=0;i<lim;++i) rev[i] = (rev[i>>1]>>1) | ((i&1)<<(L-1));
FFT(a,lim,1), FFT(b,lim,1), FFT(c,lim,1);
for(int i=0;i<lim;++i) a[i] = a[i]*b[i], c[i] = c[i]*b[i];
FFT(a,lim,-1), FFT(c,lim,-1);
for(int i=1;i<=n;++i)
printf("%.3lf\n",a[i].x-c[n-i].x);
return 0;
}
/*
5
4006373.885184
15375036.435759
1717456.469144
8514941.004912
1410681.345880 -16838672.693
3439.793
7509018.566
4595686.886
10903040.872
*/

Summary

推理过程大部分是我自己的手推,和其他大佬不同请见谅,毕竟条条大路通罗马呗。

过程中的转折点就是我推不动的时候,所这个题目让我学会了:

  • 卷积 (雾 ,要把式子推成卷积形式

  • 巧设数组代替抽象的数学式子

  • 翻转序列的技巧

我多项式还是太菜了呢,赶紧去做题吧!QAQ

最新文章

  1. javascript关于立即函数
  2. Redis的三种启动方式
  3. iOS之initialize与load
  4. 将windows上的文件同步到linux上
  5. 使用WinINet和WinHTTP实现Http访问
  6. 【转】java环境配置
  7. TypeScript学习笔记(二):基本数据类型及数据转换
  8. CSS之鼠标经过字体光标形状的改变
  9. Fiddler On Linux
  10. python Hbase Thrift pycharm 及引入包
  11. NodeJs的包漏洞扫描与漏洞测试攻击
  12. Servlet做简单的ajax增删改查(分页)
  13. [js高手之路]node js系列课程-图解express+supervisor+ejs用法
  14. 改变oracle数据库归档模式_译文
  15. C++ stl 怎么打印内存内容?
  16. sql sever insert into混合嵌套插入
  17. LeetCode算法扫题系列83
  18. Golang查缺补漏(一)
  19. 运维grep语法
  20. Spring boot(四)thymeleaf使用介绍

热门文章

  1. Activity 的状态都有哪些?
  2. ceph集群添加osd
  3. 远程连接mysql要点 虚拟主机定义与分类
  4. nodejs 之简单web服务器
  5. 阶段3 2.Spring_08.面向切面编程 AOP_9 spring基于注解的AOP配置
  6. Linux监控命令之==&gt;strace
  7. 查询sq字段逗号分隔的方式
  8. 【Linux 网络编程】常用套接字类型
  9. Python笔记(二)
  10. 滑雪(dp或记忆化搜索)