1. 程式人生 > >[bzoj2154]Crash的數字表格(mobius反演)

[bzoj2154]Crash的數字表格(mobius反演)

printf code type void include spa ray bool max

題意:$\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^m {lcm(i,j)} } $

解題關鍵:

$\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^m {lcm(i,j)} } = \sum\limits_{i = 1}^n {\sum\limits_{j = 1}^m {\frac{{i*j}}{{\gcd (i,j)}}} } $

枚舉gcd,上式化為:

$\sum\limits_{d = 1}^{\min (n,m)} {d\sum\limits_{\begin{array}{*{20}{c}}
{\gcd (i,j) = = 1}\\

{1 < = i < = n/d}\\
{1 < = j < = m/d}
\end{array}} {i*j} } $

$f(n,m,k) = \sum\limits_{\begin{array}{*{20}{c}}
{\gcd (i,j) = = k}\\
{1 < = i < = n}\\
{1 < = j < = m}
\end{array}} {i*j} $

由於

$sum(n,m) = \sum\limits_{\begin{array}{*{20}{c}}
{1 < = i < = n}\\
{1 < = j < = m}

\end{array}} {i*j} = \frac{{n(n + 1)}}{2}\frac{{m(m + 1)}}{2}$

$\begin{array}{l}
F(n,m,k) = \sum\limits_{k|d} {f(n,m,d) = \sum\limits_{\begin{array}{*{20}{c}}
{1 < = i < = n}\\
{1 < = j < = m}\\
{k|\gcd (i,j)}
\end{array}} {i*j} = {k^2}sum(\left\lfloor {\frac{n}{k}} \right\rfloor ,\left\lfloor {\frac{m}{k}} \right\rfloor )} \\

f(n,m,k) = \sum\limits_{k|d} {u(\frac{d}{k})F(n,m,d)}
\end{array}$

而此題中,$k==1$,

則,

$\begin{array}{l}
f(n,m,1) = \sum\limits_{d = 1}^{\min (n,m)} {u(d)F(n,m,d)} \\
= \sum\limits_{d = 1}^{\min (n,m)} {u(d){d^2}sum(\left\lfloor {\frac{n}{d}} \right\rfloor ,\left\lfloor {\frac{m}{d}} \right\rfloor )} \\
ans = \sum\limits_{d = 1}^{\min (n,m)} {d*f(\left\lfloor {\frac{n}{d}} \right\rfloor ,\left\lfloor {\frac{m}{d}} \right\rfloor ,1)}
\end{array}$

求解ans和$f$函數復雜度都是$O(\sqrt n )$,所以最終復雜度為$O(n)$。

 1 #include<cstdio>
 2 #include<cstring>
 3 #include<cstdlib>
 4 #include<algorithm>
 5 #include<cmath>
 6 #include<iostream>
 7 #define Sum(x,y) (x*(x+1)/2%mod*(y*(y+1)/2%mod)%mod)
 8 using namespace std;
 9 typedef long long ll;
10 const ll mod=20101009;
11 const ll maxn=10000000+200;
12 ll n,m;
13 bool vis[maxn];
14 ll prime[maxn],mu[maxn],sum1[maxn];
15 void init_mu(ll n){
16     ll cnt=0;
17     mu[1]=1;
18     for(ll i=2;i<n;i++){
19         if(!vis[i]){
20             prime[cnt++]=i;
21             mu[i]=-1;
22         }
23         for(ll j=0;j<cnt&&i*prime[j]<n;j++){
24             vis[i*prime[j]]=1;
25             if(i%prime[j]==0){mu[i*prime[j]]=0;break;}
26             else { mu[i*prime[j]]=-mu[i];}
27         }
28     }
29     sum1[0]=0;
30     for(ll i=1;i<n;i++) sum1[i]=(sum1[i-1]+1ll*mu[i]*i*i)%mod;
31 }
32 ll calf(ll n,ll m){
33     ll ans=0,pos,len=min(n,m);
34     for(ll i=1;i<=len;i=pos+1){
35         pos=min(n/(n/i),m/(m/i));
36         //ans+=(sum1[pos]-sum1[i-1])%mod*((n/i)*((n/i)+1)/2%mod*(((m/i)+1)*(m/i)/2%mod)%mod);
37         ans+=(sum1[pos]-sum1[i-1])%mod*Sum(n/i,m/i)%mod;//最好用函數寫出 
38         ans%=mod;
39         
40     }
41     return ans;
42 }
43 ll calres(ll n,ll m){
44     ll ans=0,pos,len=min(n,m);
45     for(ll i=1;i<=len;i=pos+1){
46         pos=min(n/(n/i),m/(m/i));
47         ans+=(i+pos)*(pos-i+1)/2%mod*calf(n/i,m/i)%mod;
48         ans%=mod;
49     }
50     return (ans%mod+mod)%mod;
51 }
52 
53 int main(){
54     scanf("%lld%lld",&n,&m);
55     init_mu(min(n,m)+10);
56     printf("%lld\n",calres(n,m));
57     return 0;
58 }

[bzoj2154]Crash的數字表格(mobius反演)