题面:

http://wenku.baidu.com/link?url=D2ORnA9xjgSxa2GlYLB7gGiYgBcXsy-Aw0kVYTjTE-iYhH1s7h8xXGmnaMwl32SYznVvzodyKZgHODl_ekzIFwEsO64ZOMIvQbsah-9kZiW

提交:

http://www.lydsy.com/JudgeOnline/problem.php?id=3527

给出n个数字(q1~qn),定义$$F_i=\sum_{j<i}{q_iq_j\over (i-j)^2}-\sum_{j>i}{q_iq_j\over (i-j)^2}$$

然后设$$E_i={F_i\over q_i}$$

求  \(E_i\)

 

分析

我们先把式子化简一下$$E_i=\sum_{j<i}{q_j\over (i-j)^2}-\sum_{j>i}{q_j\over (i-j)^2}$$

然后我们令$$f[i]=q_i,g[i]={1\over i^2}$$

然后发现左边好像卷积的形式$$c_i=\sum_{j=0}^ia_jb_{i-j}$$

但是没有 \(j=0\) 和 \(j=i\) 的情况.没关系,我们令 \(f[i]=0 , g[i]=0\) .

这样的话原式\[\sum_{j=1}^{i-1}{q_j\over (i-j)^2}\]就和\[\sum_{j=0}^i{q_j\over (i-j)^2}\]相等了

这样左边就是\[A_i=\sum_{j=0}^if[j]g[i-j]\]的卷积了

我们来看右边,右边也很有成为卷积的潜质啊,可是不满足\(0\le{j}\le{i}\)啊.没关系,我们发现左边的式子和右边的式子正好相反,所以我们考虑"倒过来",

于是就有原式\[\sum_{j=i+1}^nf[j]g[i-j]\]等价于\[\sum_{j=i}^nf[j]g[j-i]\]又等价于\[\sum_{j=0}^{n-i}f[n-j]g[n-j-i]\]

这个变化可以通过换元实现.我们可以设\(f'[x]=f'[n-x]\),这样的话右式就等于$$\sum_{j=0}^{n-i}f'[j]g[n-i-j]$$

这样右边就是$$B_i=\sum_{j=0}^if'[j]g[i-j]$$的卷积了

最后$$E_i=A_i+B_{n-i}$$

开心地FFT吧~

 

注意:

1.由于数组是从1开始的,预留了\(f[0]与g[0]\)的位置,所以长度应该+1,所以len=2*n-1+1+1=2*n+1.

 

 1 #include <bits/stdc++.h>
 2 using namespace std;
 3 
 4 const int maxn=263000;
 5 const double pi=acos(-1.0);
 6 int n,len;
 7 int rev[maxn];
 8 double f[maxn],f_[maxn],g[maxn],ans[maxn];
 9 struct cp{
10     double r,i;
11     cp(double r=0,double i=0):r(r),i(i){}
12     cp operator + (const cp &x) const { return cp(r+x.r,i+x.i); }
13     cp operator - (const cp &x) const { return cp(r-x.r,i-x.i); }
14     cp operator * (const cp &x) const { return cp(r*x.r-i*x.i,r*x.i+i*x.r); }
15 }a[maxn],b[maxn],c[maxn],A[maxn];
16 void brc(int &len){
17     memset(rev,-1,sizeof rev);
18     int k=1,l=0;
19     while(k<len) k<<=1, l++;
20     len=k;
21     for(int i=1;i<len-1;i++){
22         if(rev[i]!=-1) continue;
23         int x=i,y=0,m=l;
24         while(m--) y<<=1, y|=(x&1), x>>=1;
25         rev[i]=y; rev[y]=i;
26     }
27 }
28 void dft(cp *a,int n,int op){
29     for(int i=1;i<n-1;i++) A[rev[i]]=a[i];
30     for(int i=1;i<n-1;i++) a[i]=A[i];
31     for(int m=2;m<=n;m<<=1){
32         cp wn(cos(2.0*pi/m*op),sin(2.0*pi/m*op));
33         for(int i=0;i<n;i+=m){
34             cp w(1.0); int k=m>>1;
35             for(int j=0;j<k;j++){
36                 cp u=a[i+j],t=w*a[i+j+k];
37                 a[i+j]=u+t;
38                 a[i+j+k]=u-t;
39                 w=w*wn;
40             }
41         }
42     }
43     if(op==-1)for(int i=0;i<n;i++) a[i].r/=n;
44 }
45 void fft(double *x,double *y,cp *a,cp* b,cp *c,int n){
46     for(int i=0;i<len;i++) a[i]=cp(x[i]), b[i]=cp(y[i]);
47     dft(a,n,1); dft(b,n,1);
48     for(int i=0;i<n;i++) c[i]=a[i]*b[i];
49     dft(c,n,-1);
50 }
51 int main(){
52     scanf("%d",&n); len=2*n+1;
53     brc(len);
54     for(int i=1;i<=n;i++) scanf("%lf",&f[i]);
55     for(int i=1;i<=n;i++) g[i]=1.0/i/i;
56     for(int i=0;i<=n;i++) f_[i]=f[n-i];
57     fft(f,g,a,b,a,len);
58     for(int i=1;i<=n;i++) ans[i]=a[i].r;
59     fft(f_,g,a,b,a,len);
60     for(int i=1;i<=n;i++) ans[i]-=a[n-i].r;
61     for(int i=1;i<=n;i++) printf("%.3lf\n",ans[i]);
62     return 0;
63 }

View Code