您的位置:首页 > 其它

Codechef Arithmetic Progressions ,分块FFT

2017-07-05 15:55 309 查看
题解网上到处都是,就不说了。

一开始块内暴力时忘了把桶清零,并且没有随下标移动更新桶内的值

然后FFT部分我一开始是这么写的

ans+=(int)(h[a[j]<<1]+0.5);


光荣地爆int了

调了下块的大小,因为FFT常数大,所以块多一点会快一些

最后1090ms(本机对于随机数据要6s多,感觉codechef的机子很强啊),似乎比卡常大爷xehoth的0.64s慢到不知哪里去。

#include<cstdio>
#include<cmath>
#include<cctype>
#include<cstring>
typedef double db;
typedef long long ll;
inline int getint(){
int x=0;
char c=getchar();
while(!isdigit(c))c=getchar();
for(;isdigit(c);c=getchar())x=x*10+c-48;
return x;
}
const db pi=acos(-1);
struct Comp{
db a,b;
Comp(db _a=0,db _b=0):a(_a),b(_b){}
Comp operator+(const Comp&x){
return Comp(a+x.a,b+x.b);
}
Comp operator-(const Comp&x){
return Comp(a-x.a,b-x.b);
}
Comp operator*(const Comp&x){
return Comp(a*x.a-b*x.b,b*x.a+a*x.b);
}
Comp operator*=(const Comp&x){
return *this=*this*x;
}
};
const int N=300010;
namespace M{
int len;
Comp x
,y
,z
,t;
inline void swap(Comp&a,Comp&b){
t=a;
a=b;
b=t;
}
inline void rader(Comp*a){
int i,j=len>>1,k;
for(i=1;i<len-1;++i){
if(i<j)swap(a[i],a[j]);
for(k=len>>1;j>=k;k>>=1)j-=k;
j+=k;
}
}
inline void FFT(Comp*a,int key){
rader(a);
int i,h,j;
Comp w,wn,u,v;
for(h=2;h<=len;h<<=1){
wn=Comp(cos(key*2*pi/h),sin(key*2*pi/h));
for(i=0;i<len;i+=h){
w=Comp(1,0);
for(j=i;j<i+(h>>1);++j){
u=a[j];
v=w*a[j+(h>>1)];
a[j]=u+v;
a[j+(h>>1)]=u-v;
w*=wn;
}
}
}
if(key==-1)for(i=0;i<len;++i)a[i].a/=len;
}
inline void cheng(db*c,db*a,db*b,int n){//a*b=c
for(len=1;len<=n<<1;len<<=1);
for(int i=0;i<len;++i)x[i]=Comp(a[i],0),y[i]=Comp(b[i],0),z[i]=Comp(0,0);
FFT(x,1);
FFT(y,1);
for(int i=0;i<len;++i)z[i]=x[i]*y[i];
FFT(z,-1);
for(int i=0;i<len;++i)c[i]=z[i].a;
}
}
int n,m,i,j,a
,cc
,be
,en
,sz,bs,k,c[30010];
db f
,g
,h
;
ll ans;
int main(){
//freopen("a.txt","r",stdin);
scanf("%d",&n);
if(n<10000)sz=sqrt(n);
else sz=7*sqrt(n);
for(i=1;i<=n;++i){
scanf("%d",a+i);
if(a[i]>m)m=a[i];
}
bs=(n-1)/sz+1;
be[1]=1;
en[1]=sz;
for(i=2;i<bs;++i)be[i]=be[i-1]+sz,en[i]=en[i-1]+sz;
be[bs]=be[bs-1]+sz;
en[bs]=n;
for(i=1;i<=bs;++i){
memset(c,0,sizeof c);
c[a[be[i]]]=1;
for(j=be[i]+1;j<en[i];++j){
for(k=j+1;k<=en[i];++k)
if(a[j]-(a[k]-a[j])>0)ans+=c[a[j]-(a[k]-a[j])];
++c[a[j]];
}
}
memset(c,0,sizeof c);
for(j=be[1];j<=en[1];++j)++c[a[j]];
for(i=2;i<=bs;++i){
for(j=be[i];j<en[i];++j)
for(k=j+1;k<=en[i];++k)
if(a[j]-(a[k]-a[j])>0)ans+=c[a[j]-(a[k]-a[j])];
for(j=be[i];j<=en[i];++j)++c[a[j]];
}
memset(c,0,sizeof c);
for(j=be[bs];j<=en[bs];++j)++c[a[j]];
for(i=bs-1;i;--i){
for(j=be[i];j<en[i];++j)
for(k=j+1;k<=en[i];++k)
if(a[k]+(a[k]-a[j])>0)ans+=c[a[k]+(a[k]-a[j])];
for(j=be[i];j<=en[i];++j)++c[a[j]];
}
for(i=be[1];i<=en[1];++i)++f[a[i]];
for(i=be[3];i<=n;++i)++g[a[i]];
for(i=2;i<bs;++i){
M::cheng(h,f,g,m);
for(j=be[i];j<=en[i];++j){
ans+=(ll)(h[a[j]<<1]+0.5);
++f[a[j]];
}
for(j=be[i+1];j<=en[i+1];++j)--g[a[j]];
}
printf("%lld\n",ans);
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: