题目背景
提示:原 P1829 半数集问题 已经迁移至 P1028 数的计算
题目描述
今天的数学课上,Crash小朋友学习了最小公倍数(Least Common Multiple)。对于两个正整数a和b,LCM(a, b)表示能同时整除a和b的最小正整数。例如,LCM(6, 8) = 24。
回到家后,Crash还在想着课上学的东西,为了研究最小公倍数,他画了一张NM的表格。每个格子里写了一个数字,其中第i行第j列的那个格子里写着数为LCM(i, j)。一个45的表格如下:
1 2 3 4 5
2 2 6 4 10
3 6 3 12 15
4 4 12 4 20
看着这个表格,Crash想到了很多可以思考的问题。不过他最想解决的问题却是一个十分简单的问题:这个表格中所有数的和是多少。当N和M很大时,Crash就束手无策了,因此他找到了聪明的你用程序帮他解决这个问题。由于最终结果可能会很大,Crash只想知道表格里所有数的和mod20101009的值。
输入格式
输入的第一行包含两个正整数,分别表示N和M。
输出格式
输出一个正整数,表示表格中所有数的和mod20101009的值。
输入输出样例
输入 #1
4 5
输出 #1
122
说明/提示
100%的数据满足N, M≤ 10^7。
解释:推公式呗~转大佬的来学习一下
Solution
易知原式等价于
∑ i = 1 n ∑ j = 1 m i ⋅ j gcd ( i , j ) \sum\limits_{i=1}^n\sum\limits_{j=1}^m\frac{i\cdot j}{\gcd(i,j)} i=1∑nj=1∑mgcd(i,j)i⋅j
枚举最大公因数 d,显然两个数除以 d 得到的数互质
∑ i = 1 n ∑ j = 1 m ∑ d ∣ i , d ∣ j , gcd ( i d , j d ) = 1 i ⋅ j d \sum\limits_{i=1}^n\sum\limits_{j=1}^m\sum\limits_{d|i,d|j,\gcd(\frac{i}{d},\frac{j}{d})=1}\frac{i\cdot j}{d} i=1∑nj=1∑md∣i,d∣j,gcd(di,dj)=1∑di⋅j
非常经典的 gcd \gcd gcd 式子的化法
∑ d = 1 n d ⋅ ∑ i = 1 ⌊ n d ⌋ ∑ j = 1 ⌊ m d ⌋ [ gcd ( i , j ) = 1 ] ⋅ i ⋅ j \sum\limits_{d=1}^n d\cdot\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}[\gcd(i,j)=1]\cdot i\cdot j d=1∑nd⋅i=1∑⌊dn⌋j=1∑⌊dm⌋[gcd(i,j)=1]⋅i⋅j
后半段式子中,出现了互质数对之积的和,为了让式子更简洁就把它拿出来单独计算。于是我们记
sum ( n , m ) = ∑ i = 1 n ∑ j = 1 m [ gcd ( i , j ) = 1 ] ⋅ i ⋅ j \operatorname{sum}(n,m)=\sum\limits_{i=1}^n\sum\limits_{j=1}^m [\gcd(i,j)=1]\cdot i\cdot j sum(n,m)=i=1∑nj=1∑m[gcd(i,j)=1]⋅i⋅j
接下来对 sum ( n , m ) \operatorname{sum}(n,m) sum(n,m) 进行化简。首先枚举约数,并将 [ gcd ( i , j ) = 1 ] [\gcd(i,j)=1] [gcd(i,j)=1] 表示为 ε ( gcd ( i , j ) ) ε ( g c d ( i , j ) ) ∑ d = 1 n ∑ d ∣ i n ∑ d ∣ j m μ ( d ) ⋅ i ⋅ j \varepsilon(\gcd(i,j))ε(gcd(i,j)) \sum\limits_{d=1}^n\sum\limits_{d|i}^n\sum\limits_{d|j}^m\mu(d)\cdot i\cdot j ε(gcd(i,j))ε(gcd(i,j))d=1∑nd∣i∑nd∣j∑mμ(d)⋅i⋅j
指上式中的 i,ji,j ),显然式子可以变为
∑ d = 1 n μ ( d ) ⋅ d 2 ⋅ ∑ i = 1 ⌊ n d ⌋ ∑ j = 1 ⌊ m d ⌋ i ⋅ j \sum\limits_{d=1}^n\mu(d)\cdot d^2\cdot\sum\limits_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum\limits_{j=1}^{\lfloor\frac{m}{d}\rfloor}i\cdot j d=1∑nμ(d)⋅d2⋅i=1∑⌊dn⌋j=1∑⌊dm⌋i⋅j
观察上式,前半段可以预处理前缀和;后半段又是一个范围内数对之和,记
g ( n , m ) = ∑ i = 1 n ∑ j = 1 m i ⋅ j = n ⋅ ( n + 1 ) 2 × m ⋅ ( m + 1 ) 2 g(n,m)=\sum\limits_{i=1}^n\sum\limits_{j=1}^m i\cdot j=\frac{n\cdot(n+1)}{2}\times\frac{m\cdot(m+1)}{2} g(n,m)=i=1∑nj=1∑mi⋅j=2n⋅(n+1)×2m⋅(m+1)
至此
sum ( n , m ) = ∑ d = 1 n μ ( d ) ⋅ d 2 ⋅ g ( ⌊ n d ⌋ , ⌊ m d ⌋ ) \operatorname{sum}(n,m)=\sum\limits_{d=1}^n\mu(d)\cdot d^2\cdot g(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor) sum(n,m)=d=1∑nμ(d)⋅d2⋅g(⌊dn⌋,⌊dm⌋)
我们可以 ⌊ n ⌊ n d ⌋ ⌋ \lfloor\frac{n}{\lfloor\frac{n}{d}\rfloor}\rfloor ⌊⌊dn⌋n⌋
数论分块求解 sum ( n , m ) \operatorname{sum}(n,m) sum(n,m)函数。
在求出 sum ( n , m ) \operatorname{sum}(n,m) sum(n,m) 后,回到定义 sum \operatorname{sum} sum 的地方,可得原式为
∑ d = 1 n d ⋅ sum ( ⌊ n d ⌋ , ⌊ m d ⌋ ) \sum\limits_{d=1}^n d\cdot\operatorname{sum}(\lfloor\frac{n}{d}\rfloor,\lfloor\frac{m}{d}\rfloor) d=1∑nd⋅sum(⌊dn⌋,⌊dm⌋)
可见这又是一个可以数论分块求解的式子!
最后分块就好了
#include <cstdio>
#include <algorithm>
using namespace std;
const int N=1e7;
const int mod=20101009;
int n,m,mu[N+5],p[N/10+5],sum[N+5];
bool flg[N+5];
void init() {
mu[1]=1;
int tot=0,k=min(n,m);
for(int i=2;i<=k;++i) {
if(!flg[i]) p[++tot]=i,mu[i]=-1;
for(int j=1;j<=tot&&i*p[j]<=k;++j) {
flg[i*p[j]]=1;
if(i%p[j]==0) {mu[i*p[j]]=0;break;}
mu[i*p[j]]=-mu[i];
}
}
for(int i=1;i<=k;++i) sum[i]=(sum[i-1]+1LL*i*i%mod*(mu[i]+mod))%mod;
}
int Sum(int x,int y) {
return (1LL*x*(x+1)/2%mod)*(1LL*y*(y+1)/2%mod)%mod;
}
int func(int x,int y) {
int res=0;
for(int i=1,j;i<=min(x,y);i=j+1) {
j=min(x/(x/i),y/(y/i));
res=(res+1LL*(sum[j]-sum[i-1]+mod)*Sum(x/i,y/i)%mod)%mod;
}
return res;
}
int solve(int x,int y) {
int res=0;
for(int i=1,j;i<=min(x,y);i=j+1) {
j=min(x/(x/i),y/(y/i));
res=(res+1LL*(j-i+1)*(i+j)/2%mod*func(x/i,y/i)%mod)%mod;
}
return res;
}
int main() {
scanf("%d%d",&n,&m);
init();
printf("%d\n",solve(n,m));
}