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

点击量：175

    理解了一些关于卷积的东西。

目录

题目描述

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

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

令\(E_i=\frac{F_i}{q_i}\)，求\(E_i\).

输入输出格式

输入格式：

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

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

输出格式：

\(n\)行，第\(i\)行输出\(E_i\)。

与标准答案误差不超过1e-2即可。

输入输出样例

输入样例#1：

5
4006373.885184
15375036.435759
1717456.469144
8514941.004912
1410681.345880

输出样例#1：

-16838672.693
3439.793
7509018.566
4595686.886
10903040.872

说明

对于\(30\%\)的数据，\(n\le 1000\)。

对于\(50\%\)的数据，\(n\le 60000\)。

对于\(100\%\)的数据，\(n\le 100000,0<q_i<1000000000\)。

题解：

    这个题对于我这种只写过FFT板子的人来说是根本看不出来的。多项式乘积的实质是卷积，即对于\(f(x)=a(x)b(x)\)，有\[f_{k}=\sum_{i=0}^{k}a_ib_{k-i}x^k\]。也就是一年级的时候老师教过：\(5\)可以分成\(1\)和\(4\)；\(5\)可以分成\(2\)和\(3\)；5可以分成3和2；5可以分成4和1。这里的分，实际上就是指指数的来源。

    看到这个题的公式要先把它化简。\[E_i=\frac{q_1}{(i-1)^2}+\frac{q_2}{(i-2)^2}+\dots+\frac{q_{i-1}}{1^2}-\frac{q_{i+1}}{1^2}-\frac{q_{i+2}}{2^2}-\dots-\frac{q_n}{(n-i)^2}\]

    对于这样一个式子，可以把分子分母联想为两个函数相乘，即令\(f(x)=q_x,g(x)=\frac{1}{x^2}\)。那么\(E_x=\sum_{i=1}^{x-1}f(i)g(x-i)-\sum_{i=x+1}^n f(i)g(i-x)\)，利用上面卷积的转化，可以发现第一个\(\sum\)是一个标准的卷积，因为下标和恒等于\(x\)。此时考虑把第二个\(\sum\)转化为卷积的形式，只要尽可能让下标和恒等于\(x\)就可以了，那么由于\(g(x)=\frac{1}{x^2}\)是个偶函数，所以\(g(x)=g(-x)\)，则上面的式子后半段变成了\(\sum_{i=x+1}^n f(i)g(x-i)\)。

    这时由于\(i>x\)，所以数组下标会出现负数，这里直接向右平移下标\(n\)个单位长度，需要重新给\(g\)数组赋值。

    再令\(A(x)=\sum_{i=1}^{x-1}f(i)g(x-i),B(x)=\sum_{i=x+1}^n f(i)g(x-i)\)。此时\(A(x)\)和\(B(x)\)均可以分别用一次卷积求出来，最后相减就得到答案了。

Code：

#include<cstdio>
#include<cstring>
#include<cmath>
const double pi=acos(-1);
struct cpl
{
    double x,y;
    cpl(double x,double y)
    {
        this->x=x;
        this->y=y;
    }
    cpl(){x=0.0,y=0.0;}
    friend cpl operator +(cpl a,cpl b)
    {return cpl(a.x+b.x,a.y+b.y);}
    friend cpl operator -(cpl a,cpl b)
    {return cpl(a.x-b.x,a.y-b.y);}
    friend cpl operator *(cpl a,cpl b)
    {return cpl(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}
}a[1<<18],b[1<<18],omg[1<<18],inv[1<<18];
int pos[1<<18],len,tot;
void init()
{
    for(int i=1;i<tot;++i)
        pos[i]=(pos[i>>1]>>1)|(i&1)<<(len-1);
    for(int i=0;i<tot;++i)
    {
        omg[i]=cpl(cos(2*pi*i/tot),sin(2*pi*i/tot));
        inv[i]=omg[i];
        inv[i].y=-inv[i].y;
    }
}
double ans[1<<18];
double q[1<<18];
void fft(cpl f[],bool flag)
{
    for(int i=0;i<tot;++i)
        if(pos[i]>i)
        {
            cpl tmp=f[i];
            f[i]=f[pos[i]];
            f[pos[i]]=tmp;
        }
    for(int bs=1;bs<=len;++bs)
        for(int i=0;i<tot;i+=(1<<bs))
        {
            int g=(1<<(bs-1));
            for(int j=0;j<g;++j)
            {
                cpl t1=f[i+j],t2;
                if(flag)
                    t2=omg[j*(1<<(len-bs))]*f[i+j+g];
                else
                    t2=inv[j*(1<<(len-bs))]*f[i+j+g];
                f[i+j]=t1+t2;
                f[i+j+g]=t1-t2;
            }
        }
}
int main()
{
    int n;
    scanf("%d",&n);
    while((1<<len)<=2*n+1)
        ++len;
    tot=(1<<len);
    init();
    for(int i=0;i<tot;++i)
        a[i]=b[i]=cpl();
    for(int i=1;i<=n;++i)
        scanf("%lf",&q[i]);
    for(int i=1;i<=n;++i)
        a[i]=cpl(q[i],0);
    for(int i=1;i<=n;++i)
        b[i]=cpl(1.0/i/i,0);
    fft(a,1);
    fft(b,1);
    for(int i=0;i<tot;++i)
        a[i]=a[i]*b[i];
    fft(a,0);
    for(int i=1;i<=n;++i)
        ans[i]=a[i].x/tot;
    for(int i=0;i<tot;++i)
        a[i]=b[i]=cpl();
    for(int i=1;i<=n;++i)
        a[i]=cpl(q[i],0);
    for(int i=1;i<=n;++i)
        b[n-i]=cpl(1.0/i/i,0);
    fft(a,1);
    fft(b,1);
    for(int i=0;i<tot;++i)
        a[i]=a[i]*b[i];
    fft(a,0);
    for(int i=1;i<=n;++i)
        ans[i]-=a[n+i].x/tot;
    for(int i=1;i<=n;++i)
        printf("%.5lf\n",ans[i]);
    return 0;
}