码迷,mamicode.com
首页 > 其他好文 > 详细

【莫比乌斯反演】BZOJ2154 Crash的数字表格

时间:2015-06-24 12:30:46      阅读:117      评论:0      收藏:0      [点我收藏+]

标签:

Description

  求sigma lcm(x,y),x<=n,y<=m。n,m<=1e7。

 

Solution

  lcm没有什么直接做的好方法,用lcm=x*y/gcd转成gcd来做

  就是要求sigma d*f(x/d,y/d)

  f(x,y)为x和y以内gcd正好为1的对数

  F为所有对数,于是有F(x,y)=x*(x+1)/2*y*(y+1)/2

  f(x,y)=sigma (1<=i<=x) i*i*mu(i)*F(x/i,y/i) 

  f用莫比乌斯反演解决,这两个式子都套上分块优化到sqrt,于是总复杂度sqrt*sqrt=n

  分块优化具体可以见上一篇blog

 

Code

  一开始全开LL MLE了一发

  然后又WA了两发,第一次是有一地方算的时候溢出

  一开始为了避免MLE把prime数组/50,但其实只能/20的样子

 1 #include<cstdio>
 2 #include<algorithm>
 3 #include<cstring>
 4 #define ll long long
 5 using namespace std;
 6 const int maxn=1e7+5,mod=20101009;
 7 
 8 bool flag[maxn];
 9 int prime[maxn],mu[maxn],cnt;
10 int sum[maxn],s[maxn];
11 int n,m;
12 
13 void getmu(){
14     mu[1]=1;
15     for(int i=2;i<=n;i++){
16         if(!flag[i]){
17             mu[i]=-1;
18             prime[++cnt]=i;
19         }
20         for(int j=1;i*prime[j]<=n&&j<=cnt;j++){
21             flag[i*prime[j]]=1;
22             if(i%prime[j]==0){
23                 mu[i*prime[j]]=0;
24                 break;
25             }
26             mu[i*prime[j]]=-mu[i];
27         }
28     }
29     for(int i=1;i<=n;i++)
30         sum[i]=(sum[i-1]+(ll)i*i*mu[i]%mod)%mod;
31 }
32 
33 ll F(int x,int y){
34     return (ll)((ll)x*(x+1)/2%mod)*((ll)y*(y+1)/2%mod)%mod;
35 }
36 
37 ll f(int x,int y){
38     ll ret=0;
39     int p;
40     for(int i=1;i<=x;i=p+1){
41         p=min(x/(x/i),y/(y/i));
42         ret=(ret+(ll)(sum[p]-sum[i-1])*F(x/i,y/i))%mod;
43     }
44     return ret;
45 }
46 
47 int main(){
48     scanf("%d%d",&n,&m);
49     if(n>m) swap(n,m);
50     for(int i=1;i<=n;i++) s[i]=(s[i-1]+i)%mod;
51     getmu();
52     
53     int pos;
54     ll ans=0;
55     for(int i=1;i<=n;i=pos+1){
56         pos=min(n/(n/i),m/(m/i));
57         ans=(ans+(ll)(s[pos]-s[i-1])*f(n/i,m/i))%mod;
58     }
59     printf("%lld\n",(ans+mod)%mod);
60     return 0;
61 }

 

【莫比乌斯反演】BZOJ2154 Crash的数字表格

标签:

原文地址:http://www.cnblogs.com/xkui/p/4597275.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!