Description
设函数g(N)表示N的约数个数。现在给出一个数M,求出所有M的约数x的g(x)的K次方和。
Input
第一行输入N,K。N表示M由前N小的素数组成。接下来N行,第i+1行有一个正整数Pi,表示第Ai小的素数 有 Pi次。等式:
Output
输出一个数,表示答案。只需输出最后答案除以1000000007的余数。
Sample Input
Sample Output
解题思路:
拿到题感觉一脸不可做,反正不是反演的样子。
先来考虑一下样例:
$2^1*3^3=54$
考虑如何将答案分类。将贡献重新累加起来。
首先$g(d)=\sum\limits_{i=1}^{n}{(1+c_i)}$()其中$c_i$为第$i$项的次数)
现在考虑如何将$M$的所有因数表示出来:
以54为例
其中每个因数可以按照2的次数分类,可以是0~1次共有2种选法。
对于每种选法,对应的数中3的次数有0~3共4种选法。
而对应3的每一种选法其答案都是$[(2的次数+1)*(3的次数+1)]^k$
那么这种结果可以看做两个多项式乘积的形式,这两个多项式都是(1+2+3+...)的形式,而由于那个指数上的k是可以带入括号的。
那么答案就变成了$f(54)=(1^3+2^3)*(1^3+2^3+3^3+4^3)=9*100=900$
现在答案就可以解出来了,答案就变成了:
$f(M)=\prod\limits_{i=1}^{n}\sum\limits_{j=1}^{p_i}{j^k}$
n是可以枚举的,对于前45分的点,可以暴力快速幂求解。
剩下的怎么办。
后面的那个高次求和$\sum\limits_{i=1}^{n}{i^k}$是非常公式化的,怎么可能没有公式。
高次求和公式会形如:
$\sum\limits_{i=1}^{n}{i^k}=\sum\limits_{i=1}^{k+1}{a_i*n^i}$
所以可以这样求解(orz zwz):
$a_1x^1+a_2x^2+...+a_{k+1}x^{k+1}+(x+1)^k=a_1(x+1)+a_2(x+1)^2+...+a_{k+1}(x+1)^{k+1}$
二项式定理展开移项:
$\sum\limits_{i=1}^{k+1}{C_i^0x^0}+\sum\limits_{i=2}^{k+1}{C_i^1x^1}+...+\sum\limits_{i=k+1}^{k+1}{C_i^kx^k}=\sum\limits_{i=0}^{k}{C_k^ix^i}$
这个杨辉三角打出矩阵消元就可以解出系数了。
代码:
1 #include2 #include 3 #include 4 typedef long long lnt; 5 const lnt mod=(lnt)(1e9+7); 6 lnt n,k; 7 lnt maxs; 8 lnt C[100][100]; 9 lnt F[200001];10 lnt p[200001];11 lnt a[1001];12 lnt b[100][100];13 lnt c[100];14 lnt ksm(lnt x,lnt y)15 {16 lnt ans=1;17 while(y)18 {19 if(y&1)ans=ans*x%mod;20 x=x*x%mod;21 y=y/2;22 }23 return ans;24 }25 lnt squ(lnt x)26 {27 return x%mod*x%mod;28 }29 int main()30 {31 scanf("%lld%lld",&n,&k);32 for(int i=1;i<=n;i++)33 scanf("%lld",&p[i]),34 maxs=std::max(p[i],maxs);35 if(maxs<=100000)36 {37 lnt ans=1;38 for(int i=1;i<=maxs+1;i++)F[i]=(F[i-1]+ksm(i,k))%mod;39 for(int i=1;i<=n;i++)ans=(ans*F[p[i]+1])%mod;40 printf("%lld\n",ans);41 return 0;42 }43 if(k==3)44 {45 lnt ans=1;46 lnt Inv=ksm(4,mod-2);47 for(int i=1;i<=n;i++)48 ans=ans*squ(p[i]%mod+1)%mod*squ(p[i]%mod+2)%mod*Inv%mod;49 printf("%lld\n",ans);50 return 0;51 }52 C[0][0]=1;53 for(int i=1;i<=k+1;i++)54 {55 C[i][0]=1;56 for(int j=1;j<=i;j++)57 C[i][j]=(C[i-1][j-1]+C[i-1][j])%mod;58 }59 for(int i=1;i<=k+1;i++)60 {61 for(int j=i;j<=k+1;j++)62 {63 b[i][j]=C[j][i-1];64 }65 }66 for(int i=1;i<=k+1;i++)a[i]=C[k][i-1];67 for(int i=k+1;i>=1;i--)68 {69 a[i]=(a[i]*ksm(b[i][i],mod-2))%mod;70 for(int j=i-1;j;j--)71 {72 (a[j]-=b[j][i]*a[i]%mod)%=mod;73 }74 }75 lnt ans=1;76 for(int i=1;i<=n;i++)77 {78 lnt tmp=0;79 for(int j=1;j<=k+1;j++)80 tmp=(tmp+(a[j]*ksm(p[i]%mod+1,j))%mod)%mod;81 ans=(ans*tmp)%mod;82 }83 printf("%lld\n",(ans%mod+mod)%mod);84 return 0;85 }