bzoj3601 一个人的数论 (拉格朗日插值求系数)
Description
有一天hjy96想到了一个数论问题:
对于一个非负整数d和一个正整数n,定义fa(n)为所有小于n且与n互质的正整数的d次方之和。如\(f_3(10) = 1^3+3^3+7^3+ 9^3\)。
现给定d,n,求fa(n)的值。输出答案对\(10^9+7\)取模后的结果。
hjy96当然知道怎么做啦!但是他想考考你.....
Input
由于n可能很大,我们给出n的质因数分解式。
\[n=\prod\limits_{i=1}^wp_i^{\alpha_i}
\]
第一行有用空格隔开的一个非负整数d,一个正整数w。
接下来w行,每行有两个用空格隔开的正整数\(p_i,\alpha_i\)。(保证\(p_i\)为素数且互不相同)
Output
一行,一个非负整数表示\(f_d(n)\)对\(10^9+7\)取模后的结果。
Sample Input
3 2
2 1
5 1
Sample Output
1100
数据范围与约定
对于所有数据 \(1\le d\le 100,1 \le w \le 1000,2 \le p_i \le 10^9,1 \le \alpha_i \le 10^9.\)
Solution
由题可知
\[ans=\sum\limits_{i=1}^n[gcd(i,n)=1]i^d
\]
推一下
\[ans=\sum\limits_{i=1}^n[gcd(i,n)=1]i^d\\
=\sum\limits_{i=1}^n\sum\limits_{t|x且t|n}\mu(t)i^d\\
=\sum\limits_{t|n}\mu(t)\sum\limits_{i=1}^{\frac{n}{t}}(i\times t)^d\\
=\sum\limits_{t|n}\mu(t)t^d\sum\limits_{i=1}^{\frac{n}{t}}i^d\\
\]
\(\sum\limits_{i=1}^{\frac{n}{t}}i^d可以表示为一个关于\frac{n}{t}的d+1次多项式\sum\limits_{i=0}^{d+1}f_i(\frac{n}{t})^i\)
\(f_i\)可以用拉格朗日插值或高斯消元求出(我用的是拉格朗日)
\[\therefore ans=\sum\limits_{t|n}\mu(t)t^d\sum\limits_{i=0}^{d+1}f_i(\frac{n}{t})^i\\
=\sum\limits_{t|n}\mu(t)\sum\limits_{i=0}^{d+1}f_in^it^{d-i}\\
=\sum\limits_{i=0}^{d+1}\sum\limits_{t|n}\mu(t)f_in^it^{d-i}\\
=\sum\limits_{i=0}^{d+1}f_in^i\sum\limits_{t|n}\mu(t)t^{d-i}\\
\]
容易发现\(\sum\limits_{t|n}\mu(t)t^{d-i}\)是一个积性函数
设\(F(n)=\sum\limits_{t|n}\mu(t)t^{d-i}\)
则有\(F(n)=\prod\limits_{i=1}^{w}F(p_i^{\alpha_i})\)
而\(F(p_i^{\alpha_i})\)显然等于\(1-p_i^{d-i}\)
\[\therefore ans=
\sum\limits_{i=0}^{d+1}f_in^i\prod\limits_{i=1}^{w}(1-p_i^{d-i})\\
\]
Code
#include<bits/stdc++.h>
const int N=1005;
const int mod=1e9+7;
int Pow(int x,int f=mod-2){
int re=1;
while(f){
if(f&1)re=1ll*re*x%mod;
f>>=1;
x=1ll*x*x%mod;
}
return re;
}
namespace Lagrange{
int inv[N],f_[N],f[N],y[N],tmp[N];
void Solve(int n){
for(int i=1;i<=n;i++){
int div=1,lst=0;
for(int j=1;j<=n;j++)
if(i!=j)div=1ll*div*(mod+i-j)%mod;
div=1ll*y[i]*Pow(div)%mod;
for(int j=0;j<n;j++){
tmp[j]=1ll*(lst+mod-f_[j])*inv[i]%mod;
f[j]=(f[j]+(1ll*div*tmp[j]%mod))%mod;
lst=tmp[j];
}
}
}
void Pre(int d){
f_[0]=1;
for(int i=1;i<=d+2;std::swap(f_,tmp),i++){
tmp[0]=0,inv[i]=Pow(i);
for(int j=1;j<=i;j++)tmp[j]=f_[j-1];
for(int j=0;j<=i;j++)tmp[j]=(tmp[j]+(1ll*f_[j]*(mod-i)%mod))%mod;
}
for(int i=1;i<=d+2;i++)y[i]=(y[i-1]+Pow(i,d))%mod;
}
}
int H(int t,int p){
if(t<0)return (1+mod-Pow(Pow(p,-t)))%mod;
return (1+mod-Pow(p,t))%mod;
}
int d,w,a[N],p[N],ans=0;
signed main(){
using namespace Lagrange;
scanf("%d%d",&d,&w);
for(int i=1;i<=w;i++)scanf("%d%d",p+i,a+i);
Pre(d);
Solve(d+2);
int pown=1;
for(int i=1;i<=d+1;i++){
int tmp=1;
for(int j=1;j<=w;j++){
pown=1ll*pown*Pow(p[j],a[j])%mod;
tmp=1ll*tmp*H(d-i,p[j])%mod;
}
ans=(ans+(1ll*f[i]*pown%mod*tmp%mod))%mod;
}
printf("%d",ans);
return 0;
}