题意简述
求
∑ i = 1 n ∑ j = 1 m l c m ( i , j ) \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}lcm(i,j) i=1∑nj=1∑mlcm(i,j)
(其中 n , m < = 1 e 7 n,m<=1e7 n,m<=1e7)
数据
输入
4 5
输出
122
思路
一看就知道是反演。。。
来推式子吧!
令
g
=
g
c
d
(
i
,
j
)
g=gcd(i,j)
g=gcd(i,j),
n
<
=
m
n<=m
n<=m
原式
=
∑
i
=
1
n
∑
j
=
1
m
i
j
g
这
个
就
是
化
了
一
下
l
c
m
(
i
,
j
)
=
∑
d
=
1
n
∑
i
=
1
n
∑
j
=
1
m
[
g
=
=
d
]
i
j
d
这
个
就
是
提
前
枚
举
g
c
d
(
i
,
j
)
,
然
后
后
面
找
哪
些
会
算
到
=
∑
d
=
1
n
∑
i
d
=
d
n
∑
j
d
=
d
m
[
g
=
=
1
]
i
d
∗
j
d
d
=
∑
d
=
1
n
∑
i
=
1
n
/
d
∑
j
=
1
m
/
d
[
g
=
=
1
]
i
j
d
i
,
j
换
成
枚
举
d
的
几
倍
,
注
意
换
成
枚
举
倍
数
后
i
变
成
i
d
,
j
变
成
j
d
,
所
以
多
乘
一
个
d
2
,
原
来
是
i
j
d
,
变
成
i
j
d
=
∑
d
=
1
n
∑
i
=
1
n
/
d
∑
j
=
1
m
/
d
∑
q
∣
g
μ
(
q
)
i
j
d
把
[
g
=
=
1
]
换
成
了
ϵ
(
g
)
,
然
后
用
了
一
下
反
演
式
=
∑
d
=
1
n
∑
q
=
1
n
/
d
∑
i
=
1
n
/
d
∑
j
=
1
m
/
d
[
q
∣
g
]
μ
(
q
)
i
j
d
提
前
枚
举
q
,
看
后
面
有
哪
些
算
到
的
=
∑
d
=
1
n
∑
q
=
1
n
/
d
∑
i
q
=
1
n
/
d
∑
j
q
=
1
m
/
d
[
1
∣
g
]
μ
(
q
)
i
q
∗
j
q
∗
d
i
,
j
有
变
成
了
枚
举
q
的
倍
数
,
注
意
i
变
成
i
q
,
j
变
成
j
q
,
所
以
多
乘
一
个
q
2
=
∑
d
=
1
n
∑
q
=
1
n
/
d
∑
i
=
1
n
/
d
q
∑
j
=
1
m
/
d
q
μ
(
q
)
i
j
q
2
d
稍
作
整
理
,
去
掉
了
没
什
么
用
的
[
1
∣
g
]
(
因
为
不
管
g
多
少
这
个
都
成
立
)
=
∑
d
=
1
n
d
∑
q
=
1
n
/
d
q
2
μ
(
q
)
∑
i
=
1
n
/
d
q
i
∑
j
=
1
m
/
d
q
j
进
一
步
的
整
理
,
把
无
关
的
提
到
最
前
=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\frac{ij}{g}\\ 这个就是化了一下lcm(i,j)\\ =\sum\limits_{d=1}^{n}\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[g==d]\frac{ij}{d}\\ 这个就是提前枚举gcd(i,j),然后后面找哪些会算到\\ =\sum\limits_{d=1}^{n}\sum\limits_{id=d}^{n}\sum\limits_{jd=d}^{m}[g==1]\frac{id*jd}{d}\\ =\sum\limits_{d=1}^{n}\sum\limits_{i=1}^{n/d}\sum\limits_{j=1}^{m/d}[g==1]ijd\\ i,j换成枚举d的几倍,注意换成枚举倍数后i变成id,j变成jd,所以多乘一个d^2,原来是\frac{ij}{d},变成ijd\\ =\sum\limits_{d=1}^{n}\sum\limits_{i=1}^{n/d}\sum\limits_{j=1}^{m/d}\sum\limits_{q|g}\mu(q)ijd\\ 把[g==1]换成了\epsilon(g),然后用了一下反演式\\ =\sum\limits_{d=1}^{n}\sum\limits_{q=1}^{n/d}\sum\limits_{i=1}^{n/d}\sum\limits_{j=1}^{m/d}[q|g]\mu(q)ijd\\ 提前枚举q,看后面有哪些算到的\\ =\sum\limits_{d=1}^{n}\sum\limits_{q=1}^{n/d}\sum\limits_{iq=1}^{n/d}\sum\limits_{jq=1}^{m/d}[1|g]\mu(q)iq*jq*d\\ i,j有变成了枚举q的倍数,注意i变成iq,j变成jq,所以多乘一个q^2\\ =\sum\limits_{d=1}^{n}\sum\limits_{q=1}^{n/d}\sum\limits_{i=1}^{n/dq}\sum\limits_{j=1}^{m/dq}\mu(q)ijq^2d\\ 稍作整理,去掉了没什么用的[1|g](因为不管g多少这个都成立)\\ =\sum\limits_{d=1}^{n}d\sum\limits_{q=1}^{n/d}q^2\mu(q)\sum\limits_{i=1}^{n/dq}i\sum\limits_{j=1}^{m/dq}j\\ 进一步的整理,把无关的提到最前\\
=i=1∑nj=1∑mgij这个就是化了一下lcm(i,j)=d=1∑ni=1∑nj=1∑m[g==d]dij这个就是提前枚举gcd(i,j),然后后面找哪些会算到=d=1∑nid=d∑njd=d∑m[g==1]did∗jd=d=1∑ni=1∑n/dj=1∑m/d[g==1]ijdi,j换成枚举d的几倍,注意换成枚举倍数后i变成id,j变成jd,所以多乘一个d2,原来是dij,变成ijd=d=1∑ni=1∑n/dj=1∑m/dq∣g∑μ(q)ijd把[g==1]换成了ϵ(g),然后用了一下反演式=d=1∑nq=1∑n/di=1∑n/dj=1∑m/d[q∣g]μ(q)ijd提前枚举q,看后面有哪些算到的=d=1∑nq=1∑n/diq=1∑n/djq=1∑m/d[1∣g]μ(q)iq∗jq∗di,j有变成了枚举q的倍数,注意i变成iq,j变成jq,所以多乘一个q2=d=1∑nq=1∑n/di=1∑n/dqj=1∑m/dqμ(q)ijq2d稍作整理,去掉了没什么用的[1∣g](因为不管g多少这个都成立)=d=1∑ndq=1∑n/dq2μ(q)i=1∑n/dqij=1∑m/dqj进一步的整理,把无关的提到最前
然后现在我们会发现,什么都好算了:
- 最后面那两个 s i g m a sigma sigma:就是等差数列求和,设 s u m ( n , m ) = ( ∑ i = 1 n i ) ( ∑ j = 1 m j ) sum(n,m)=(\sum\limits_{i=1}^{n}i)(\sum\limits_{j=1}^{m}j) sum(n,m)=(i=1∑ni)(j=1∑mj),显然根据搞死求和公式,它 = ( n ( n + 1 ) 2 ) ( m ( m + 1 ) 2 ) =(\frac{n(n+1)}{2})(\frac{m(m+1)}{2}) =(2n(n+1))(2m(m+1))。后面两个 s i g m a sigma sigma就变成了 s u m ( n / d q , m / d q ) sum(n/dq,m/dq) sum(n/dq,m/dq),就珂以 O ( 1 ) O(1) O(1)求了。
- 如果能求出 μ ( q ) \mu(q) μ(q)的所有值,那么 q 2 μ ( q ) q^2\mu(q) q2μ(q)的所有值也能求,同时也就能求 q 2 μ ( q ) q^2\mu(q) q2μ(q)的前缀和。那么 ∑ q = 1 n / d q 2 μ ( q ) ∑ i = 1 n / d q i ∑ j = 1 m / d q j = ∑ q = 1 n / d q 2 μ ( q ) s ( n / d q , m / d q ) \sum\limits_{q=1}^{n/d}q^2\mu(q)\sum\limits_{i=1}^{n/dq}i\sum\limits_{j=1}^{m/dq}j=\sum\limits_{q=1}^{n/d}q^2\mu(q)s(n/dq,m/dq) q=1∑n/dq2μ(q)i=1∑n/dqij=1∑m/dqj=q=1∑n/dq2μ(q)s(n/dq,m/dq)就珂以整除分块了,设其为 f ( n / d , m / d ) f(n/d,m/d) f(n/d,m/d)
- 设了不少函数!现在原式变为 ∑ d = 1 n d f ( n / d , m / d ) \sum\limits_{d=1}^{n}df(n/d,m/d) d=1∑ndf(n/d,m/d),这个整除分块太 t m ^{tm} tm明显了吧!
计算一下时间复杂度。一次
f
(
n
/
d
,
m
/
d
)
f(n/d,m/d)
f(n/d,m/d)最坏是
f
(
n
,
m
)
f(n,m)
f(n,m),也就是
O
(
n
)
O(\sqrt{n})
O(n)。
外层在套一层整除分块,也就是
O
(
n
)
O(\sqrt{n})
O(n)。
两层一乘,
O
(
n
)
O(n)
O(n)过了。
(而实际上并不是准确的
O
(
n
)
O(n)
O(n),稍微要小一点,应该是
O
(
∑
i
=
1
n
n
/
i
)
O(\sum\limits_{i=1}^{\sqrt{n}}\sqrt{n/i})
O(i=1∑nn/i),反正小于
O
(
n
)
O(n)
O(n)就是了。
O
(
n
)
O(n)
O(n)是能过的,我们的代码当然珂以过)
代码:
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define mod 20101009
#define N 10000010
namespace Flandle_Scarlet
{
int mu[N],primes[N];
bool notp[N];
int d2mu[N];
//d2mu[i]=i*i*μ(i)的前缀和
void Init()
{
mu[1]=notp[1]=1;
int n=10000000;
int &cnt=primes[0];
for(int i=2;i<=n;++i)//线性筛mu
{
if (!notp[i])
{
primes[++cnt]=i;
mu[i]=-1;
}
for(int j=1;j<=cnt and i*primes[j]<=n;++j)
{
int u=primes[j];
notp[i*u]=1;
if (i%u==0)
{
mu[i*u]=0;
break;
}
else
{
mu[i*u]=-mu[i];
}
}
}
for(int i=1;i<=n;++i)
{
d2mu[i]=(d2mu[i-1]+i*i%mod*(mu[i]+mod))%mod;
//这个式子十分显然,但是注意取膜时的优先级
}
}
int sum(int n,int m)
//也就是上面说的sum
//但是一定要注意优先级!!!我错了好多次!!!
{
return ((n*(n+1)/2)%mod*(m*(m+1)/2)%mod)%mod;
}
int f(int n,int m)
//也就是上面说的f
{
int ans=0;
for(int l=1,r;l<=min(n,m);l=r+1)
{
r=min(n/(n/l),m/(m/l));
//d2mu[r]-d2mu[l-1]计算中间和
//sum(n/l,m/l)是相同部分
ans=(ans+(d2mu[r]-d2mu[l-1]+mod)*sum(n/l,m/l)%mod)%mod;
}
return ans;
}
int calc(int n,int m)
//主求解函数
{
int ans=0;
for(int l=1,r;l<=min(n,m);l=r+1)
{
r=min(n/(n/l),m/(m/l));
//这个太显然了。。。不想解释
ans=(ans+(r-l+1)*(l+r)/2%mod*f(n/l,m/l)%mod)%mod;
}
return ans;
}
int n,m;
void Main()
{
Init();
scanf("%lld%lld",&n,&m);
printf("%lld\n",calc(n,m));
}
};
main()
{
Flandle_Scarlet::Main();
return 0;
}