Bzoj3512 DZY Loves Math IV
Submit: 447 Solved: 223
Description
Input
Output
Sample Input
Sample Output
数据规模:
1<=n<=10^5,1<=m<=10^9,本题共4组数据。
HINT
Source
数学问题 莫比乌斯反演 杜教筛
数学题的麻烦在于,如果你拖久了不更博,就不得不再推一遍式子了。
注意到n比较小,只有1e5,这是一个可以枚举的范围,我们来看看如果枚举n会怎么样:
设$S(n,m)=\sum_{i=1}^{m} \varphi(i*n)$
有一个结论:
若$|\mu(n)|=1$,即n是无平方数,则:
$$S(n,m)=\sum_{i=1}^{m} \varphi(i) \sum_{d|gcd(i,n)} \varphi(\frac{n}{d})$$
然而博主并不知道如何能想到这个结论。手推一下发现它是对的。
http://www.cnblogs.com/ziliuziliu/p/6508562.html
↑这里有一个关于它的神奇证明
于是当n是无平方数的时候可以这么求:
$$S(n,m)=\sum_{d|n} \varphi(\frac{n}{d}) \sum_{i=1}^{\lfloor \frac{m}{d} \rfloor} \varphi(i*d)$$
$$S(n,m)=\sum_{d|n} \varphi(\frac{n}{d}) S(d,\lfloor \frac{m}{d} \rfloor)$$
记忆化一下即可。
而当$\mu(n)=0$的时候就很方便了:
我们知道
$\varphi(n)= n*\Pi_{i=1}^{k} (1-(1/p_k)) $
其中p数组是n的全部质因数。我们已经可以算出每个质因数都只出现一次时的答案,那么对于剩下的质因数,直接乘进去就行。
也就是说,找到最大的x使得$|\mu(\frac{n}{x})|=1$,
$$S(n,m)=x*S(\frac{n}{x},m)$$
后面这个S用最上面的方法求。
递归边界是$S(1,m)$,也就是 $\mu $ 的前缀和,这个可以用杜教筛算。
最后,$ans=\sum_{i=1}^{n}S(n,m)$
1 #include<iostream> 2 #include<algorithm> 3 #include<cstring> 4 #include<cstdio> 5 #include<cmath> 6 #include<map> 7 #define LL long long 8 using namespace std; 9 const int mod=1e9+7; 10 const int mxn=5000010; 11 const int inv2=500000004; 12 int read(){ 13 int x=0,f=1;char ch=getchar(); 14 while(ch<'0' || ch>'9'){if(ch=='-')f=-1;ch=getchar();} 15 while(ch>='0' && ch<='9'){x=x*10+ch-'0';ch=getchar();} 16 return x*f; 17 } 18 int pri[mxn>>1],cnt=0; 19 int phi[mxn],mx[mxn]; 20 int smm[mxn]; 21 void init(){ 22 phi[1]=1;mx[1]=1; 23 for(int i=2;i<mxn;i++){ 24 if(!phi[i]){ 25 pri[++cnt]=i; 26 phi[i]=i-1; 27 mx[i]=1; 28 } 29 for(int j=1;j<=cnt && pri[j]*i<mxn;j++){ 30 int tmp=pri[j]*i; 31 if(i%pri[j]==0){ 32 mx[tmp]=mx[i]*pri[j]; 33 phi[tmp]=phi[i]*pri[j]; 34 break; 35 } 36 phi[tmp]=phi[i]*(pri[j]-1); 37 mx[tmp]=mx[i]; 38 } 39 } 40 for(int i=1;i<mxn;i++) 41 smm[i]=((LL)smm[i-1]+phi[i])%mod; 42 return; 43 } 44 map<int,int>mphi; 45 int Calc_phi(int x){ 46 if(x<mxn)return smm[x]; 47 if(mphi.count(x))return mphi[x]; 48 int res=(LL)x*(x+1)%mod*inv2%mod; 49 for(int i=2,pos;i<=x;i=pos+1){ 50 int y=x/i; pos=x/y; 51 res=((LL)res-(pos-i+1)*Calc_phi(y))%mod; 52 } 53 res=(res+mod)%mod; 54 mphi[x]=res; 55 return res; 56 } 57 map<LL,int>mpS; 58 int solve(int n,int m){ 59 if(!m)return 0; 60 if(n==1)return Calc_phi(m); 61 LL tar=(LL)n*mod+m; 62 if(mpS.count(tar))return mpS[tar]; 63 int res=0; 64 for(int i=1;i*i<=n;i++){ 65 if(n%i==0){ 66 res=((LL)res+phi[n/i]*(LL)solve(i,m/i))%mod; 67 if(i*i!=n)res=((LL)res+phi[i]*(LL)solve(n/i,m/(n/i)))%mod; 68 } 69 } 70 mpS[tar]=res; 71 return res; 72 } 73 int n,m; 74 int main(){ 75 int i,j; 76 init(); 77 n=read();m=read(); 78 int ans=0; 79 for(i=1;i<=n;i++) 80 ans=(ans+(LL)mx[i]*solve(i/mx[i],m))%mod; 81 ans=(ans+mod)%mod; 82 printf("%d\n",ans); 83 return 0; 84 }