●BZOJ 2154 Crash的数字表格

题链:

http://www.lydsy.com/JudgeOnline/problem.php?id=2154

题解:

莫比乌斯反演。

题意还是很清楚的,就不赘述了。

显然有

$ANS=\sum_{i=1\;j=1}^{n\;m} lcm(i,j)$

化为较为熟悉的gcd形式:

$\quad\quad=\sum_{i=1\;j=1}^{n\;m} \frac{i \times j}{gcd(i,j)}$

令$g$为gcd的值,$F(n,m)=\sum i\times j,满足1\leq i \leq \lfloor \frac{n}{g}\rfloor,1\leq j \leq \lfloor \frac{m}{g} \rfloor,且gcd(i,j)=1$

那么,$ANS=\sum_{g=1}^{min(n,m)} g\times F(n,m)$

如果已知$F(n,m)$的值,那么ANS就可以在$O\sqrt N$的复杂度内求出。


接下来看看对于一个确定的g,$F(n,m)$怎么求

即我们要求满足$1\leq i \leq \lfloor \frac{n}{g}\rfloor,1\leq j \leq \lfloor \frac{m}{g} \rfloor$,且$gcd(i,j)=1$的所有$i \times j$的和

还记得这个题么?●HDU 1695 GCD,让求的是满足上述条件的$(i,j)$的对数。

我们仍可以类似地去求:(令$x=\lfloor \frac{n}{g}\rfloor,y=\lfloor \frac{m}{g} \rfloor$)

令$t(k)为gcd(i,j)=k的i\times j的和$

$T(k)为gcd(x,y)=\lambda k的i\times j的和$

显然$T(k)=\sum_{k|d}{f(d)}$,即T为t的倍数关系和函数

考虑一下如何计算T(k)。

首先i的取值有:$k,2k,3k,\cdots,\lfloor \frac{x}{k} \rfloor k$

首先j的取值有:$k,2k,3k,\cdots,\lfloor \frac{y}{k} \rfloor k$

任意的两两组合都是T(k)的一部分。

令$sum(a,b)=\frac{(1+a)a}{2}\times\frac{(1+b)b}{2}$

所以$T(k)=k^2\times sum(\lfloor \frac{x}{k} \rfloor,\lfloor \frac{y}{k} \rfloor)$

那么由莫比乌斯反演公式的形式二(倍数关系那个式子):

$t(k)=\sum_{k|d}\mu(\frac{d}{k})T(d)$

$\quad\quad=\sum_{k|d}\mu(\frac{d}{k})d^2\times sum(\lfloor \frac{x}{d} \rfloor,\lfloor \frac{y}{d} \rfloor)$

而我们要求的是t(1),所以

$t(1)=\sum_{d=1}^{min(\lfloor \frac{n}{g}\rfloor,\lfloor \frac{m}{g}\rfloor)}\mu(d)d^2\times sum(\lfloor \frac{x}{d} \rfloor,\lfloor \frac{y}{d} \rfloor)$

显然,这个式子也可以在$O\sqrt N$的复杂度内求出。

所以综上,时间复杂度为$O(N)$

代码:

#include<cstdio>
#include<cstring>
#include<iostream>
#define MAXN 10000050
using namespace std;
const int mod=20101009;
int mu[MAXN],pmui2[MAXN];
void Sieve(int n){
static bool np[MAXN];
static int prime[MAXN],pnt;
mu[1]=pmui2[1]=1;
for(int i=2;i<=n;i++){
if(!np[i]) prime[++pnt]=i,mu[i]=-1;
for(int j=1;j<=pnt&&i<=n/prime[j];j++){
np[i*prime[j]]=1;
if(i%prime[j]) mu[i*prime[j]]=-mu[i];
else{mu[i*prime[j]]=0; break;}
}
pmui2[i]=(pmui2[i-1]+1ll*mu[i]*i%mod*i%mod)%mod;
}
}
int sum(int n,int m){
return (1ll*(1+n)*n/2%mod)*(1ll*(1+m)*m/2%mod)%mod;
}
int F(int n,int m){
int mini=min(n,m),ret=0;
for(int d=1,last;d<=mini;d=last+1){
last=min(n/(n/d),m/(m/d));
ret=(1ll*ret+1ll*(pmui2[last]-pmui2[d-1]+mod)%mod*sum(n/d,m/d)%mod)%mod;
}
return ret;
}
int main(){
int n,m,mini,ans=0;
scanf("%d%d",&n,&m); mini=min(n,m);
Sieve(mini);
for(int g=1,last;g<=mini;g=last+1){
last=min(n/(n/g),m/(m/g));
ans=(1ll*ans+1ll*(g+last)*(last-g+1)/2%mod*F(n/g,m/g)%mod)%mod;
}
printf("%d\n",(ans+mod)%mod);
return 0;
}

  

上一篇:CentOS6.5yum配置本地源


下一篇:C语言有哪些鲜为人知的特性?