【BZOJ 4332】 4332: JSOI2012 分零食 (FFT+快速幂)

时间:2023-03-09 00:24:55
【BZOJ 4332】 4332: JSOI2012 分零食 (FFT+快速幂)

4332: JSOI2012 分零食

Time Limit: 10 Sec  Memory Limit: 256 MB
Submit: 119  Solved: 66

Description

这里是欢乐的进香河,这里是欢乐的幼儿园。 
今天是2月14日,星期二。在这个特殊的日子里,老师带着同学们欢乐地跳着,笑着。校长从幼儿园旁边的小吃店买了大量的零食决定分给同学们。听到这个消息,所有同学都安安静静地排好了队,大家都知道,校长不喜欢调皮的孩子。 
同学们依次排成了一列,其中有A位小朋友,有三个共同的欢乐系数O,S和U。如果有一位小朋友得到了x个糖果,那么她的欢乐程度就是f(x)=O*x2+S*x+U。 
现在校长开始分糖果了,一共有M个糖果。有些小朋友可能得不到糖果,对于那些得不到糖果的小朋友来说,欢乐程度就是1。如果一位小朋友得不到糖果,那么在她身后的小朋友们也都得不到糖果。(即这一列得不到糖果的小朋友一定是最后的连续若干位) 
所有分糖果的方案都是等概率的。现在问题是:期望情况下,所有小朋友的欢乐程度的乘积是多少?呆呆同学很快就有了一个思路,只要知道总的方案个数T和所有方案下欢乐程度乘积的总和S,就可以得到答案Ans=S/T。现在他已经求出来了T的答案,但是S怎么求呢?他就不知道了。你能告诉他么? 
因为答案很大,你只需要告诉他S对P取模后的结果。 
后记: 
虽然大家都知道,即便知道了T,知道了S对P取模后的结果,也没有办法知道期望情况下,所有小朋友欢乐程度的乘积。但是,当呆呆想到这一点的时候,已经彻底绝望了。 

Input

第一行有2个整数,分别是M和P。 
第二行有一个整数A,第三行有一个整数O。 
第四行有一个整数S,第五行有一个整数U。 

Output

一个整数S,因为答案可能很大,你只需要输出S 对P取模后的结果。 

Sample Input

4 100
4
1
0
0

Sample Output

63

样例说明
函数f(x)=x^2。一共有4份零食,4位同学。如果只有第一个同学得到,欢乐程度为16,若前两位同学得到,欢乐程度的所有可能依次为9,9,16,若有三位同学得到,欢乐程度有4,4,4,最后一种情况,每一个同学都得到了零食,欢乐程度为1。相加后得到S=63。

应上传者要求,此题不公开,如有异议,请提出.

HINT

对于100%的数据,M<=10000,P<=255,A<=108,O<=4,S<=300,U<=100。

Source

【分析】


O(n^2)做法:【实测可以过全部

 #include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
using namespace std;
int Mod;
int f[][]; int main()
{
int m,n,a,b,c;
scanf("%d%d%d%d%d%d",&m,&Mod,&n,&a,&b,&c);
if(n>m) n=m;
memset(f,,sizeof(f));
f[][]=;
int ans=,p=,ss;
for(int j=;j<=n;j++)
{
int nw=,pp=p^;
memset(f[pp],,sizeof(f[pp]));
ss=;
for(int i=;i<=m;i++)
{
if(i>=) ss=(ss+f[p][i-])%Mod;
f[pp][i]=f[pp][i-];
if(i>=) nw=(nw+ss**a)%Mod;
if(i>=) nw=(nw+f[p][i-]*(*a+a+b))%Mod;
f[pp][i]=(f[pp][i]+nw+f[p][i-]*(a+b+c))%Mod;
}p^=;
ans=(ans+f[p][m])%Mod;
}
printf("%d\n",ans);
return ;
}

具体解释见这里


O(n^2 logn)

【实测还不如暴力,不知道是不是我常数大。

【就是后面两个for改成卷积形式

 #include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cmath>
using namespace std;
const double pi=acos(-);
#define Maxn 10010*4
int Mod;
int f[Maxn],g[Maxn]; struct P
{
double x,y;
P() {x=y=;}
P(double x,double y):x(x),y(y){}
friend P operator + (P x,P y) {return P(x.x+y.x,x.y+y.y);}
friend P operator - (P x,P y) {return P(x.x-y.x,x.y-y.y);}
friend P operator * (P x,P y) {return P(x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x);}
}a[Maxn],b[Maxn]; int R[Maxn],nn,m;
void dft(P *a,int t)
{
for(int i=;i<=nn;i++) if(i<R[i]) swap(a[i],a[R[i]]);
for(int i=;i<nn;i<<=)
{
P wn(cos(pi/i),t*sin(pi/i));
for(int j=;j<nn;j+=(i<<))
{
P w(,);
for(int k=;k<i;k++,w=w*wn)
{
P x=a[j+k],y=w*a[j+k+i];
a[j+k]=x+y;a[j+k+i]=x-y;
}
}
}
if(t==-)
{
for(int i=;i<=nn;i++) a[i].x/=nn;
}
} void fft(int *X,int *Y)
{
for(int i=;i<=nn;i++)
{
a[i].x=X[i];a[i].y=;
b[i].x=Y[i];b[i].y=;
}
dft(a,);dft(b,);
for(int i=;i<=nn;i++) a[i]=a[i]*b[i];
dft(a,-);
for(int i=;i<=m;i++) X[i]=(int)(a[i].x+0.5)%Mod;
} int main()
{
int n,a,b,c;
scanf("%d%d%d%d%d%d",&m,&Mod,&n,&a,&b,&c);
if(n>m) n=m;
f[]=g[]=;
for(int i=;i<=m;i++) f[i]=g[i]=(a*i*i+b*i+c)%Mod;
int ans=; ans+=f[m]; nn=;int ll=;
while(nn<=*m) nn<<=,ll++;
for(int i=;i<=nn;i++) R[i]=(R[i>>]>>)|((i&)<<(ll-));
for(int i=;i<=n;i++)
{
fft(f,g);
ans=(ans+f[m])%Mod;
}
printf("%d\n",ans);
return ;
}

O(nlogn^2)做法:

把上面的第一个循环改成快速幂

从7点搞到了现在。。

  理解奥爷爷的代码好久啊。。。但是打的真心短。。。

  下面那个qpow是一个矩阵乘法的快速幂!!【我傻啊看了很久才看出来。。

  所以求$G^1+G^2+...G^n$

【BZOJ 4332】 4332: JSOI2012 分零食 (FFT+快速幂)【BZOJ 4332】 4332: JSOI2012 分零食 (FFT+快速幂)【BZOJ 4332】 4332: JSOI2012 分零食 (FFT+快速幂)

呵呵【我看这个看了好久

【呵呵

其实我觉得这个卷积有点迷

但是不管了,明天再说

关于$a^1+a^2+a^3...+a^n$,当然如果你只是要求数的话,可以直接套等比数列求和公式。但如果不是数呢?或者mod下没有逆元怎么求分母呢?

可以用 矩阵乘法快速幂 ,就类似上面的那个方法的。

打成数值的形式是这样:

 #include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
using namespace std; int main()
{
int x,y;
scanf("%d%d",&x,&y);
int ans=,xx=x;
while(y)
{
if(y&)
{
ans=ans*x+xx*;
}
xx=xx*x+xx*;
x=x*x;
y>>=;
}
printf("%d\n",ans);
}

上面的求卷积的幂的和,就是这样子做的。

然后是O(n logn^2)的代码

 #include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cmath>
using namespace std;
#define Maxn 10010*4
const double pi=acos(-);
int Mod; struct P
{
double x,y;
P() {x=y=;}
P(double x,double y):x(x),y(y){}
friend P operator + (P x,P y) {return P(x.x+y.x,x.y+y.y);}
friend P operator - (P x,P y) {return P(x.x-y.x,x.y-y.y);}
friend P operator * (P x,P y) {return P(x.x*y.x-x.y*y.y,x.x*y.y+x.y*y.x);}
}a[Maxn],b[Maxn]; int R[Maxn],nn,m;
void dft(P *a,int f)
{
for(int i=;i<nn;i++) if(i<R[i]) swap(a[i],a[R[i]]);
for(int i=;i<nn;i<<=)
{
P wn(cos(pi/i),f*sin(pi/i));
for(int j=;j<nn;j+=i<<)
{
P w(,);
for(int k=;k<i;k++,w=w*wn)
{
P x=a[j+k],y=w*a[j+k+i];
a[j+k]=x+y;a[j+k+i]=x-y;
}
}
}
if(f==-)
{
for(int i=;i<=nn;i++) a[i].x/=nn,a[i].y/=nn;
}
} int A[Maxn],B[Maxn],C[Maxn],nw[Maxn];
int aa,bb,cc;
void fft(int *A,int *B)
{
for(int i=;i<nn;i++)
{
a[i].x=A[i];a[i].y=;
b[i].x=B[i];b[i].y=;
}
dft(a,);dft(b,);
for(int i=;i<=nn;i++) a[i]=a[i]*b[i];
dft(a,-);
for(int i=;i<=m;i++) A[i]=((int)(a[i].x+0.5)%Mod);
} void add(int *A,int *B)
{
for(int i=;i<=m;i++) A[i]=(A[i]+B[i])%Mod;
} void qpow(int k)
{
for(int i=;i<=m;i++) A[i]=;
for(int i=;i<=m;i++) C[i]=B[i]=(aa*i*i+bb*i+cc)%Mod;
while(k)
{
if(k&)
{
fft(A,B);
add(A,C);
}
for(int i=;i<=m;i++) nw[i]=C[i];
fft(C,B);
add(C,nw);
fft(B,B);
k>>=;
}
} int main()
{
int n;
scanf("%d%d%d%d%d%d",&m,&Mod,&n,&aa,&bb,&cc);
if(n>m) n=m;
nn=;int ll=;
while(nn<=*m) nn<<=,ll++;
for(int i=;i<=nn;i++) R[i]=(R[i>>]>>)|((i&)<<(ll-));
qpow(n);
printf("%d\n",A[m]);
return ;
}

2017-04-14 21:46:24