Problem Description

今天的数学课上，Crash小朋友学习了最小公倍数(Least Common Multiple)。对于两个正整数a和b，LCM(a, b)表示能同时被a和b整除的最小正整数。例如，LCM(6, 8) = 24。回到家后，Crash还在想着课上学的东西，为了研究最小公倍数，他画了一张N*M的表格。每个格子里写了一个数字，其中第i行第j列的那个格子里写着数为LCM(i, j)。

一个4*5的表格如下：

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只想知道表格里所有数的和mod 20101009的值。

Input

输入的第一行包含两个正整数，分别表示N和M。

100%的数据满足N, M ≤ 10^7。

Output

输出一个正整数，表示表格中所有数的和mod 20101009的值。

Examples

Input

4 5

Output

122

思路：

题目本质是给出 n、m 两个数，求 $res=\sum_{i=1}^n\sum_{j=1}^mlcm(i,j)$

由于 $lcm(i,j)=\frac{ij}{gcd(i,j)}$ ，因此 $res=\sum_{i=1}^n\sum_{j=1}^m\frac{ij}{gcd(i,j)}$

将 gcd(i,j) 提出来，有： $res=\sum_{d=1}^{min(n,m)}\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=d]\frac{ij}{d}$

把 d 提出来，有： $res=\sum_{d=1}^{min(n,m)}\sum_{i=1}^{\frac{n}{d}}\sum_{j=1}^{\frac{m}{d}}[gcd(i,j)=1]ijd$

设 $x=\frac{n}{d},y=\frac{m}{d}$ ，则有： $res=\sum_{d=1}^{min(n,m)}d\sum_{i=1}^x\sum_{j=1}^y[gcd(i,j)=1]ij$

设 $f(d)=\sum_{i=1}^x\sum_{j=1}^y[gcd(i,j)=d]ij$ ，那么， $res=\sum_{d=1}^{min(n,m)}d*f(1)$

再设 $g(d)=\sum_{i=1}^x\sum_{j=1}^y[d|gcd(i,j)]ij$

可以看出，f(d)、g(d) 符合莫比乌斯反演的倍数关系，即： $g(k)=\sum_{k|d}f(d)\Leftrightarrow f(k)=\sum_{k|d}u(\frac{d}{k})g(d)$

对 g(d) 进行化简，有： $g(d)=d^2\sum_{i=1}^{\frac{x}{d}}\sum_{j=1}^{\frac{y}{d}}[1|gcd(i,j)]ij=d^2\sum_{i=1}^{\frac{x}{d}}\sum_{j=1}^{\frac{y}{d}}ij$

那么， $f(1)=\sum_{i=1}^xu(i)g(i)=\sum_{i=1}^xu(i)i^2\sum_{i=1}^\frac{x}{d}\sum_{j=1}^{\frac{y}{d}}ij$

因此，需要预处理 u(i)*i^2 的前缀和，然后利用整除分块求出后面的 ij 的和算出 f(1) 后，再计算 res

Source Program

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<string>
#include<cstring>
#include<cmath>
#include<ctime>
#include<algorithm>
#include<utility>
#include<stack>
#include<queue>
#include<vector>
#include<set>
#include<map>
#include<bitset>
#define PI acos(-1.0)
#define INF 0x3f3f3f3f
#define LL long long
#define Pair pair<int,int>
LL quickPow(LL a,LL b){ LL res=1; while(b){if(b&1)res*=a; a*=a; b>>=1;} return res; }
LL quickModPow(LL a,LL b,LL mod){ LL res=1; a=a%mod; while(b){if(b&1)res=(a*res)%mod; a=(a*a)%mod; b>>=1;} return res; }
LL getInv(LL a,LL mod){ return quickModPow(a,mod-2,mod); }
const double EPS = 1E-10;
const int MOD = 20101009;
const int N = 12000000+5;
const int dx[] = {-1,1,0,0,-1,-1,1,1};
const int dy[] = {0,0,-1,1,-1,1,-1,1};
using namespace std;int mu[N];
int prime[N];
bool bprime[N];
int cnt;
LL sum[N];
void getMu(int n){//线性筛求莫比乌斯函数cnt=0;mu[1]=1;//根据定义，μ(1)=1memset(bprime,false,sizeof(bprime));for(int i=2;i<=n;i++){//求2~n的莫比乌斯函数if(!bprime[i]){prime[++cnt]=i;//存储质数mu[i]=-1;//i为质数时，μ(1)=-1}for(int j=1;j<=cnt&&i*prime[j]<=n;j++){//枚举i之前的素数个数bprime[i*prime[j]]=true;//不是质数if(i%prime[j])//i不是prime[j]的整数倍时，i*prime[j]就不会包含相同质因子mu[i*prime[j]]=-mu[i];//mu[k]=mu[i]*mu[prime[j]]，因为prime[j]是质数，mu值为-1else{mu[i*prime[j]]=0;break;//留到后面再筛}}}
}
LL cal(LL n,LL m){LL res=0;LL minn=min(n,m);for(LL left=1,right;left<=minn;left=right+1){right=min(n/(n/left),m/(m/left));LL G=( ((1+n/left)*(n/left)/2)%MOD ) * ( ((1+m/left)*(m/left)/2)%MOD ) %MOD;res=( res+((sum[right]-sum[left-1])%MOD*G%MOD)%MOD )%MOD;}return (res+MOD)%MOD;
}
int main(){getMu(12000000);LL n,m;scanf("%lld%lld",&n,&m);LL minn=min(n,m);for(LL i=1;i<=minn;i++){LL temp=(((i%MOD*i%MOD)%MOD)%MOD*mu[i]%MOD)%MOD;sum[i]=(sum[i-1]+temp+MOD)%MOD;}LL res=0;for(LL left=1,right;left<=minn;left=right+1){right=min(n/(n/left),m/(m/left));LL temp=(left+right)*(right-left+1)/2%MOD;res=(res+(cal(n/left,m/left)%MOD*temp%MOD)%MOD)%MOD;}printf("%lld\n",res);// system("pause");return 0;
}

Crash的数字表格（HYSBZ-2154）相关推荐

【BZOJ】【2154】Crash的数字表格
莫比乌斯反演 PoPoQQQ讲义第4题题解:http://www.cnblogs.com/jianglangcaijin/archive/2013/11/27/3446169.html 感觉两次sq ...
BZOJ 2154 [国家集训队]Crash的数字表格 / JZPTAB（莫比乌斯反演，经典好题）（Luogu P1829）
整理的算法模板合集: ACM模板点我看算法全家桶系列!!! 实际上是一个全新的精炼模板整合计划 P1829 [国家集训队]Crash的数字表格 / JZPTAB(反演,经典好题) Problem S ...
BZOJ 2154 Crash的数字表格 (莫比乌斯反演)
Crash的数字表格今天的数学课上,Crash小朋友学习了最小公倍数(Least Common Multiple).对于两个正整数a和b,LCM(a, b)表示能同时被a和b整除的最小正整数.例如, ...
bzoj2154 Crash的数字表格
2154: Crash的数字表格 Time Limit: 20 Sec Memory Limit: 259 MB Submit: 4549 Solved: 1643 [Submit][Status ...
bzoj2154: Crash的数字表格
好神的莫比乌斯函数然后O(sqrt(n)*sqrt(n))好神的优化啊. #include<cstdio> #include<cstring> #include<ccty ...
P1829 [国家集训队]Crash的数字表格（推了好久的mobius反演）
P1829 [国家集训队]Crash的数字表格 / JZPTAB 推导过程 ∑i=1n∑j=1mlcm(i,j)\sum_{i = 1} ^{n} \sum_{j = 1} ^{m} lcm(i, j ...
P1829 [国家集训队]Crash的数字表格 / JZPTAB
P1829 [国家集训队]Crash的数字表格 / JZPTAB 题意: 求∑i=1n∑j=1mlcm(i,j)\sum_{i=1}^{n}\sum_{j=1}^{m}lcm(i,j)∑i=1n∑j ...
【BZOJ2154】Crash的数字表格 [莫比乌斯反演]
Crash的数字表格 Time Limit: 20 Sec Memory Limit: 259 MB [Submit][Status][Discuss] Description 今天的数学课上,Cr ...
P1829 [国家集训队]Crash的数字表格 / JZPTAB（莫比乌斯反演）
[国家集训队]Crash的数字表格 / JZPTAB 题目描述今天的数学课上,Crash 小朋友学习了最小公倍数(Least Common Multiple).对于两个正整数 a a a 和 b b ...
Crash的数字表格
Crash的数字表格求$\sum_{i=1}^N\sum_{j=1}^Mlcm(i,j)$ 解设$N<M$,显然有 \[\sum_{i=1}^N\sum_{j=1}^M\frac{i ...

Crash的数字表格（HYSBZ-2154）

Problem Description

Input

Output

Examples

Source Program

Crash的数字表格（HYSBZ-2154）相关推荐

最新文章

热门文章