【BZOJ】3527: [Zjoi2014]力(fft+卷积)

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

好好的一道模板题,我自己被自己坑了好久。。

首先题目看错。。。。。。。什么玩意。。。。。。。首先题目要求:

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

然后设

$$E_i=\frac{F_i}{q_i}$$

然后我将$F_j$看成$F_i$。。。。作死系列。。。。。。

然后fft写错一个地方查错。。。。

然后第二次调用fft忘记清空。。。查错。。。(完全忘记虚数也要初始化啊!!!!!!!!!!!!!!!!!!!!!!!!!!!)

虽然说调用完dft-1时,虚数部分是0,但是你知道有精度这个东西。。。。。。。。。。。。。。。。。。。

一直以为是double精度的问题?????????????????????????????????????

原来本来就是精度问题。。。。。。。。。。。。。。。。。。。。。。。。。。。。

 

首先我们来分析$E_i$

展开可得:

$$E_i=\sum_{j<i} \frac{q_j}{(j-i)^2}-\sum_{j>i} \frac{q_j}{(j-i)^2}$$

(注意i和j都变了23333333333

然后再变化一下:

$$E_i=\sum_{j=1}^{i-1} \frac{q_j}{(j-i)^2} - \sum_{j=i+1}^{n} \frac{q_j}{(j-i)^2}$$

似乎不能看出什么?但是应该根据数据范围以及i和j,我们应该可以想想:卷积。。

然后设$f[i]=q_i, g[i]=\frac{1}{i^2}$,原式可以化为

$$E_i=\sum_{j=1}^{i-1} f[j]g[i-j] - \sum_{j=i+1}^{n} f[j]g[j-i]$$

发现卷积的形式是

$$c_i=\sum_{j=0}^{i} a_j b_{i-j}$$

似乎就是j的初值不同?

两种方法:一是将读入的$q_i$按0~n-1存储,二是像我这样做。

首先学习了《具体数学》中的和式级数,发现我是多么的sb。。

首先左边的卷积显然,上下界可以去掉,我们只来说右边:

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

其实应该这样写好理解。。

$$\sum_{i<j<=n} f[j]g[j-i]$$

然后直接换元(具体数学里称交换律。。

$$\sum_{i<j+i<=n} f[j+i]g[j+i-i]$$

然后再处理一下得

$$\sum_{0<=j<=n-i} f[j+i]g[j]$$

然后再处理一下。(上式=0是因为不影响,即g[0]=0)

$$\sum_{0<=n-i-j<=n-i} f[n-i-j+i]g[n-i-j]$$

变成

$$\sum_{0<=j<=n-i} f[n-j]g[n-i-j]$$

设$t=n-i$,则

$$\sum_{0<=j<=t} f[n-j]g[t-j]$$

设$f'[j]=f[n-j]$

$$\sum_{0<=j<=t} f'[j]g[t-j]$$

囧。。。一写下来发现也很麻烦?手推的时候很多步都可以省略吧。。。

 

//以下的作废,因为太复杂了。。。。。。233

首先j=0的话,我们可以将f[0]=0就行了嘛。。。i-1我们可以g[0]=0就好了嘛。同理j=i+1也能化成j=i

那么原式变成

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

发现右边还不行?特殊的技巧。。。减法。。

于是

$$E_i=\sum_{j=0}^{i} f[j]g[i-j] - \sum_{j=0}^{n-i} f[n-j]g[(n-j)-i]$$

设$ff[j]=f[n-j]$,那么

$$E_i=\sum_{j=0}^{i} f[j]g[i-j] - \sum_{j=0}^{n-i} ff[j]g[n-i-j]$$

然后设向量A和B,那么

$$A[i]=\sum_{j=0}^{i} f[j]g[i-j]$$

$$B[i]=\sum_{j=0}^{i} ff[j]g[i-j]$$

所以

$$E_i=A[i]-B[n-i]$$

//以上的作废

然后A和B用卷积算出来即可。。

还有一个要注意的地方,数组要开大点,就是(2*n+1)*2这样。。。。。。。。否则。。。。。。。。。。。。。。。。。。

(另,这是spj)

#include <cstdio>
#include <cstring>
#include <cmath>
#include <string>
#include <iostream>
#include <algorithm>
#include <queue>
#include <set>
#include <vector>
#include <map>
using namespace std;
typedef long long ll;
#define pii pair<int, int>
#define mkpii make_pair<int, int>
#define pdi pair<double, int>
#define mkpdi make_pair<double, int>
#define pli pair<ll, int>
#define mkpli make_pair<ll, int>
#define rep(i, n) for(int i=0; i<(n); ++i)
#define for1(i,a,n) for(int i=(a);i<=(n);++i)
#define for2(i,a,n) for(int i=(a);i<(n);++i)
#define for3(i,a,n) for(int i=(a);i>=(n);--i)
#define for4(i,a,n) for(int i=(a);i>(n);--i)
#define CC(i,a) memset(i,a,sizeof(i))
#define read(a) a=getint()
#define print(a) printf("%d", a)
#define dbg(x) cout << (#x) << " = " << (x) << endl
#define error(x) (!(x)?puts("error"):0)
#define printarr2(a, b, c) for1(_, 1, b) { for1(__, 1, c) cout << a[_][__]; cout << endl; }
#define printarr1(a, b) for1(_, 1, b) cout << a[_] << '\t'; cout << endl
inline const int getint() { int r=0, k=1; char c=getchar(); for(; c<'0'||c>'9'; c=getchar()) if(c=='-') k=-1; for(; c>='0'&&c<='9'; c=getchar()) r=r*10+c-'0'; return k*r; }
inline const int max(const int &a, const int &b) { return a>b?a:b; }
inline const int min(const int &a, const int &b) { return a<b?a:b; }

struct cp {
	double r, i;
	cp(long double _r=0.0, double _i=0.0) : r(_r), i(_i) { }
	cp operator+ (const cp &x) { return cp(r + x.r, i + x.i); }
	cp operator- (const cp &x) { return cp(r - x.r, i - x.i); }
	cp operator* (const cp &x) { return cp(r*x.r - i*x.i, r*x.i + i*x.r); }
};
const int N=400005;
const double PI=acos(-1.0);
cp A[N];
int rev[N];
void dft(cp *a, int n, int flag) {
	rep(i, n) A[rev[i]]=a[i];
	rep(i, n) a[i]=A[i];
	for(int m=2; m<=n; m<<=1) {
		cp wn(cos(2.0*PI/m*flag), sin(2.0*PI/m*flag));
		for(int i=0; i<n; i+=m) {
			cp w(1.0); int mid=(m>>1);
			for(int j=0; j<mid; ++j) {
				cp u=a[i+j+mid]*w, t=a[i+j];
				a[i+j]=t+u;
				a[i+j+mid]=t-u;
				w=w*wn;
			}
		}
	}
	if(flag==-1) rep(i, n) a[i].r/=n;
}
void init(int &len) {
	int k=1, l=0;
	while(k<len) k<<=1, ++l;
	len=k;
	rep(i, len) {
		int t=i, r=0; k=l;
		while(k--) r<<=1, r|=(t&1), t>>=1;
		rev[i]=r;
	}
}
void fft(double *x, double *y, cp *a, cp *b, cp *c, int len) {
	rep(i, len) a[i].r=x[i], a[i].i=0.0; //虚数不初始化就是作死
	rep(i, len) b[i].r=y[i], b[i].i=0.0;
	dft(a, len, 1); dft(b, len, 1);
	rep(i, len) c[i]=a[i]*b[i];
	dft(c, len, -1);
}
int n, len;
double g[N], f[N], ff[N], ans[N];
cp a[N], b[N];
int main() {
	read(n); len=n*2+1;
	init(len);
	for1(i, 1, n) scanf("%lf", &f[i]);
	for1(i, 1, n) g[i]=1.0/i/i;
	for1(i, 0, n) ff[i]=f[n-i];

	fft(f, g, a, b, a, len);
	for1(i, 1, n) ans[i]=a[i].r;

	fft(ff, g, a, b, a, len);
	for1(i, 1, n) ans[i]-=a[n-i].r;

	for1(i, 1, n) printf("%.3f\n", ans[i]);
	return 0;
}

  

 


 

 

Description

给出n个数qi,给出Fj的定义如下: 
 
令Ei=Fi/qi,求Ei

Input

第一行一个整数n。 
接下来n行每行输入一个数,第i行表示qi。 

Output

n行,第i行输出Ei。 
与标准答案误差不超过1e-2即可。 

Sample Input

5
4006373.885184
15375036.435759
1717456.469144
8514941.004912
1410681.345880

 

Sample Output

 

-16838672.693
3439.793
7509018.566
4595686.886
10903040.872

 

Hint

对于30%的数据,n≤1000。 
对于50%的数据,n≤60000。 
对于100%的数据,n≤100000,0<qi<1000000000。 

Source

感谢nodgd放题

 

posted @ 2014-11-27 14:13  iwtwiioi  阅读(2879)  评论(0编辑  收藏  举报