[CSP-S模拟测试]:随（快速幂+数学）
题目描述
给出$n$个正整数$a_1,a_2...a_n$和一个质数mod。一个变量$x$初始为$1$。进行$m$次操作.每次在$n$个数中随机选一个$a_i$,然后$x=x\times a_i$.问$m$次操作之后$x$的取值的期望。答案一定可以表示成$\frac{a}{b}$的精确分数形式。$a$和$b$可能很大,所以只需要输出$a\times b^{{10}^9+5}$模${10}^9+7$的结果。
输入格式
第一行三个整数$n,m,mod$。
接下来一行$n$个空格隔开的正整数$a_1,a_2...a_n$。
输出格式
一行一个整数表示答案。
样例
样例输入
2 1 3
1 2
样例输出
500000005
数据范围与提示
第$1$个测试点：$mod=2$。
第$2$个测试点：$n=1$。
第$3,4,5$个测试点：$m\leqslant 1,000,1\leqslant mod \leqslant 300$。
第$6,7,8$个测试点：$1\leqslant mod \leqslant 300$。
对于全部测试点：$1\leqslant a_i <mod,mod$为质数$1\leqslant mod \leqslant 1,000,1\leqslant n\leqslant {10}^5,1\leqslant m\leqslant {10}^9$。
吕老师教你学数学
质数$P$的原根$g$满足$1\leqslant rt<P$，且$rt$的$1$次方，$2$次方$...(P-1)$次方在模$P$意义下可以取遍$1$到$(P-1)$的所有整数。
欧拉定理：对于质数$P,1\leqslant x<P$的任意$x$的$(P-1)$次方在模$P$意义下都为$1$。
显然，原根的$1$次方，$2$次方$...（P-2）$次方在模$P$意义下都不为$1$，只有$(P-1)$次方在模$P$意义下为$1$。
这也是一个数成为原根的充分必要条件。
题解
第$1$个测试点：
$mod=2$，不要慌，不用考虑$mod$完会有$0$的情况，因为$1\leqslant a_i <mod$，于是$a_i$一定是$1$，所以$x$就是$1$，直接$puts("1");return 0;$就好了。
期望得分：$10$分，共$10$分。
第$2$个测试点：
$n=1$，因为每次只能选这一个数，所以快速幂直接搞就好了。
期望得分：$10$分，共$20$分。
第$3,4,5$个测试点：
考虑$DP$，定义$dp[i][j]$表示第$i$次操作后$x$为$j$的概率。
那么我们能写出$DP$式子：$dp[i][j\times a[k]\% mod]=\sum \limits_{k=1}^ndp[i-1][[j]$。
时间复杂度：$\Theta (n\times m\times mod)$。
但是注意到$3,4,5$这三个测试点并没有对$n$的范围加以约束，所以显然会$T$飞起，那么我们考虑优化这个$DP$。
接着注意$1\leqslant a_i <mod$，所以一共最多会有$mod-1$种不同的$a_i$，试图从这里入手。
开一个数组$sum[i]$记录$i$这个数出现的次数。
那么$DP$式子即可转化为：$dp[i][j\times k\%mod]=\sum \limits_{k=0}^{mod-1}dp[i-1][j]\times sum[k]$。
时间复杂度：$\Theta (m\times {mod}^2)$。
期望得分：$30$分$(RE)$，共$50$分。
如果你说你不希望$RE$，那么你可以把上面$dp$的第一维滚起来。
期望得分：$30$分$(TLE)$，共$50$分。
第$6,7,8$个测试点：
考虑对上面的算法进行优化，$m$很大，考虑矩阵快速幂。
时间复杂度：$\Theta ({mod}^3\times \log m)$。
期望得分：$60$分，共$80$分。
第$9,10$个测试点：
考虑对第$3,4,5$个测试点的$DP$进行优化，还是从$m$下手。
根据上面的$DP$式子可得：$dp[2\times i][j\times k\% mod]=\sum \limits_{j=0}^{mod-1}\sum \limits_{k=0}^{mod-1}dp[i][j]\times dp[i][k]$。
接着使用滚动数组，快速幂求解即可。
代码时刻
$30$分解法（非滚动）：
#include<bits/stdc++.h>
using namespace std;
int n,m,mod;
long long sum[1001];
long long fz,fm,ans;
long long dp[1001][1001];
long long qpow(long long x,long long y)
{
long long res=1;
while(y)
{
if(y&1)res=(res*x)%1000000007;
x=(x*x)%1000000007;
y>>=1;
}
return res;
}
int main()
{
scanf("%d%d%d",&n,&m,&mod);
for(int i=1;i<=n;i++)
{
int x;
scanf("%d",&x);
sum[x]++;
}
dp[0][1]=1;
for(int i=1;i<=m;i++)
for(int j=0;j<mod;j++)
for(int k=0;k<mod;k++)
dp[i][j*k%mod]=(dp[i][j*k%mod]+dp[i-1][j]*sum[k]%1000000007)%1000000007;
for(int i=0;i<mod;i++)
fz=(fz+dp[m][i]*i%1000000007)%1000000007;
fm=qpow(n,m);
fm=qpow(fm,1000000005);
ans=fz*fm%1000000007;
printf("%lld",ans);
return 0;
}
$30$分解法（滚动）：
#include<bits/stdc++.h>
using namespace std;
int n,m,mod;
long long sum[1001];
long long ans;
long long dp[2][1001];
long long qpow(long long x,long long y)
{
long long res=1;
while(y)
{
if(y&1)res=(res*x)%1000000007;
x=(x*x)%1000000007;
y>>=1;
}
return res;
}
int main()
{
scanf("%d%d%d",&n,&m,&mod);
for(int i=1;i<=n;i++)
{
int x;
scanf("%d",&x);
sum[x]++;
}
dp[0][1]=1;
for(int i=1;i<=m;i++)
{
memset(dp[i&1],0,sizeof(dp[i&1]));
for(int j=0;j<mod;j++)
for(int k=0;k<mod;k++)
dp[i&1][j*k%mod]=(dp[i&1][j*k%mod]+dp[i&1^1][j]*sum[k]%1000000007)%1000000007;
}
for(int i=0;i<mod;i++)
ans=(ans+dp[m&1][i]*i%1000000007)%1000000007;
printf("%lld",ans*qpow(qpow(n,m),1000000005)%1000000007);
return 0;
}
$100$分解法：
#include<bits/stdc++.h>
using namespace std;
int n,m,mod;
long long dp[2][100001],dp03[2][100001];
bool p,q;
long long ans;
long long qpow(long long x,long long y)
{
long long res=1;
while(y)
{
if(y&1)res=(res*x)%1000000007;
x=(x*x)%1000000007;
y>>=1;
}
return res;
}
void judge()
{
q^=1;
memset(dp03[q],0,sizeof(dp03[q]));
for(int i=0;i<mod;i++)
for(int j=0;j<mod;j++)
dp03[q][1LL*i*j%mod]=(dp03[q][1LL*i*j%mod]+dp03[q^1][i]*dp03[q^1][j]%1000000007)%1000000007;
}
int main()
{
scanf("%d%d%d",&n,&m,&mod);
for(int i=1;i<=n;i++)
{
int x;
scanf("%d",&x);
dp03[0][x]++;
}
dp[0][1]=1;
long long inv=qpow(qpow(n,m),1000000005)%1000000007;
while(m)
{
if(m&1)
{
p^=1;
memset(dp[p],0,sizeof(dp[p]));
for(int i=0;i<mod;i++)
for(int j=0;j<mod;j++)
dp[p][1LL*i*j%mod]=(dp[p][1LL*i*j%mod]+dp[p^1][i]*dp03[q][j]%1000000007)%1000000007;
}
m>>=1;
judge();
}
for(int i=0;i<mod;i++)
ans=(ans+i*dp[p][i])%1000000007;
ans=(ans*inv)%1000000007;
printf("%lld",ans);
return 0;
}
rp++
