【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 */

 

posted @ 2017-09-26 21:52  Troywar  阅读(378)  评论(0)  编辑  收藏  举报