n<=100000个实数$q_i$求数组f:
$f(i)=\sum_{j<i}\frac{q_j}{(j-i)^2}-\sum _{j>i}\frac{q_j}{(j-i)^2},i\in [0,n)$
题目的i和j是反过来的。。难受,不管啦。
两部分,先来第一个:
$\sum_{j<i}\frac{q_j}{(j-i)^2}= \sum_{j=0}^{i-1}q_j\frac{1}{(i-j)^2}$,一个卷积,简单。
再来第二个:
等等你会<号不会>??把数列反过来不就完了??等等我们要优雅点,假装是自己“推”出来的:
$\sum _{j>i}\frac{q_j}{(j-i)^2}=\sum_{j=i+1}^{n-1}q_j\frac{1}{(j-i)^2}$
闪一句:令k=j-i。
$=\sum_{k=1}^{n-i-1}q_{k+i}\frac{1}{k^2}=\sum_{k=1}^{n-i-1}q‘_{n-i-1-k}\frac{1}{k^2}$
很好又一卷积。
写的时候注意下,一个数组乘完就废了。。别再变回去,精度伤不起
1 #include<cstdio> 2 #include<cstdlib> 3 #include<algorithm> 4 #include<cstdlib> 5 #include<cstring> 6 #include<complex> 7 #include<math.h> 8 //#include<iostream> 9 using namespace std; 10 11 int n,m,wei; 12 #define maxn 262222 13 const double pi=acos(-1); 14 typedef complex<double> cp; 15 cp a[maxn],b[maxn],c[maxn],d[maxn],e[maxn],f[maxn]; int rev[maxn]; 16 17 void dft(cp *a,int n,int type) 18 { 19 for (int i=0;i<n;i++) if (i<rev[i]) {cp t=a[i]; a[i]=a[rev[i]]; a[rev[i]]=t;} 20 for (int i=1;i<n;i<<=1) 21 { 22 cp base=cp(cos(pi/i),type*sin(pi/i)); 23 for (int j=0,p=(i<<1);j<n;j+=p) 24 { 25 cp t=cp(1,0); 26 for (int k=0;k<i;k++,t*=base) 27 { 28 cp tmp=t*a[j+k+i]; 29 a[j+k+i]=a[j+k]-tmp; 30 a[j+k]+=tmp; 31 } 32 } 33 } 34 } 35 36 void mul(cp *a,cp *b,cp *c) 37 { 38 if (!rev[1]) for (int i=1;i<n;i++) 39 { 40 rev[i]=0; 41 for (int j=0;j<wei;j++) rev[i]|=((i>>j)&1)<<(wei-j-1); 42 } 43 dft(a,n,1); dft(b,n,1); 44 for (int i=0;i<n;i++) c[i]=a[i]*b[i]; 45 dft(c,n,-1); 46 for (int i=0;i<n;i++) c[i]/=n; 47 } 48 49 int main() 50 { 51 scanf("%d",&n); double x; 52 for (int i=0;i<n;i++) scanf("%lf",&x),a[i]=d[n-1-i]=x; a[n]=d[n]=0; 53 b[0]=e[0]=0; for (int i=1;i<=n;i++) b[i]=e[i]=1.0/i/i; 54 m=n+n; for (n=1,wei=0;n<=m;n<<=1,wei++); mul(a,b,c); m>>=1; 55 m<<=1; mul(d,e,f); m>>=1; 56 for (int i=0;i<m;i++) printf("%.3lf\n",c[i].real()-f[m-i-1].real()); 57 return 0; 58 }