【bzoj 3601】一个人的数论 (莫比乌斯反演+伯努利数)
题解:
(吐槽:网上题解那个不严谨猜测真是没谁了……关键是还猜得辣么准……)
直接化简到求和那一段:
$f_{d}(n)=\sum_{t|n}\mu(t)t^{d}\sum_{i=1}^{\frac{n}{t}}i^{d}$
$设S_{d}(T)=\sum_{i=1}^{T}i^{d}$
那这个是什么呢?伯努利数(我会说我百度找到的吗……)
$s_{p}(T)=\sum_{i=1}^{p+1}\frac{(-1)^{p+1-i}C_{p+1}^{i}B_{p+1-i}}{p+1}n^{i}$(这个是百度那个公式化过来的)
然后$n^{i}$前面那一堆玩意就是网上题解的$a_{i}$。
接下来的化简我就不解释了……http://www.cnblogs.com/jianglangcaijin/p/4033399.html
然后我们只要求出C和B这题就没了。(貌似可以把d的范围再扩10倍233)
代码:
1 #include<cstdio> 2 using namespace std; 3 typedef long long ll; 4 const ll mod=1e9+7; 5 const int N=10000; 6 inline ll read(){ 7 ll s=0,k=1;char ch=getchar(); 8 while(ch<'0'|ch>'9') ch=='-'?k=-1:0,ch=getchar(); 9 while(ch>47&ch<='9') s=s*10+(ch^48),ch=getchar(); 10 return s*k; 11 } 12 inline ll powmod(ll a,ll b){ 13 ll ans=1; 14 if(b<0) 15 return powmod(powmod(a,mod-2),-b); 16 a%=mod; 17 while(b){ 18 if(b&1) ans=ans*a%mod; 19 b>>=1;a=a*a%mod; 20 }return ans; 21 } 22 ll w,d; 23 ll p[N],pk[N]; 24 ll tot=1; 25 inline ll calc(int n){ 26 ll t=powmod(tot,n); 27 for(int i=1;i<=w;i++){ 28 t=t*(1ll-powmod(p[i],d-n))%mod; 29 } 30 if(t<0) t+=mod; 31 return t; 32 } 33 ll c[105][105],b[105]; 34 int main(){ 35 d=read(),w=read(); 36 ll n=w; 37 for(int i=1;i<=n;i++){ 38 p[i]=read(),pk[i]=read(); 39 tot=tot*powmod(p[i],pk[i])%mod; 40 } 41 c[0][0]=1; 42 for(int i=1;i<=101;i++){ 43 c[i][0]=1; 44 for(int j=1;j<=i;j++) 45 c[i][j]=(c[i-1][j]+c[i-1][j-1])%mod; 46 } 47 b[0]=1; 48 for(int i=1;i<=101;i++){ 49 for(int j=0;j<i;j++) 50 b[i]=(b[i]+c[i+1][j]*b[j])%mod; 51 b[i]=b[i]*(-powmod(i+1,mod-2))%mod+mod; 52 b[i]%=mod; 53 } 54 ll ans=0; 55 ll inv=powmod(d+1,mod-2); 56 for(int i=1;i<=d+1;i++){ 57 ll temp=((d+1-i&1)?-1:1)*c[d+1][i]*b[d+1-i]%mod*inv%mod; 58 if(temp==0) 59 continue; 60 temp=temp*calc(i)%mod; 61 ans+=temp; 62 ans%=mod; 63 } 64 printf("%lld\n",(ans%mod+mod)%mod); 65 } 66 /* 67 3 2 68 2 1 69 5 1 70 */
没有什么不可能。