題目背景

提示：原 P1829 半數集問題 已經遷移至 P1028 數的計算

題目描述

今天的數學課上，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只想知道表格裏所有數的和mod20101009的值。

輸入輸出格式

輸入格式：

輸入的第一行包含兩個正整數，分別表示N和M。

輸出格式：

輸出一個正整數，表示表格中所有數的和mod20101009的值。

輸入輸出樣例

4 5

122

說明

30%的數據滿足N, M≤ 10^3。

70%的數據滿足N, M≤ 10^5。

100%的數據滿足N, M≤ 10^7。

很神的一道題。。聽說反演是人類最大的本能然而我不會

公式自然還是不會寫

不過我發現以後可以用手機拍公式QAQ

先貼上代碼：

 1 #include<iostream>
 2 #include<cstdio>
 3 #include<cstring>
 4 #define mod 20101009
 5 #define M 10000010
 6 #define ll long long
 7 #define int long long
 8 using namespace
 
 std;
 9 int n,m,tot;
10 ll P[M],mul[M],pre[M];
11 bool vis[M];
12 ll ans;
13 void get_mul()
14 {
15     mul[1]=1;
16     for(int i=2;i<=n;i++)
17     {
18         if(!vis[i]) P[++tot]=i,mul[i]=-1;
19         for(int j=1;j<=tot&&P[j]*i<=n;j++)
20         {
21             vis[P[j]*i]=true;
22             if(i%P[j]==0)
23             {
24                 mul[P[j]*i]=0;
25                 break;
26             }
27             else mul[P[j]*i]=-mul[i];
28         }
29     }
30 }
31 ll cal(int x,int y)
32 {
33     ll ans=0;
34     for(int l=1,r;l<=x;l=r+1)
35     {
36         r=min(x/(x/l),y/(y/l));
37         ll sum=(ll)((1+x/l)*(x/l)/2%mod)*(ll)((1+y/l)*(y/l)/2%mod)%mod;
38         ans+=(ll)(pre[r]-pre[l-1])%mod*sum%mod;
39         ans%=mod;
40     }
41     return (ans+mod)%mod;
42 }
43 #undef int
44 int main()
45 {
46     #define int long long
47     scanf("%lld%lld",&n,&m);
48     get_mul();
49     if(n>m) swap(n,m);
50     for(int i=1;i<=n;i++) pre[i]=(pre[i-1]+(ll)i*i%mod*mul[i]%mod+mod)%mod;
51     for(int l=1,r;l<=n;l=r+1)
52     {
53         r=min(n/(n/l),m/(m/l));
54         ll sum=(ll)(l+r)*(r-l+1)/2%mod;
55         ans=(ans+(ll)cal(n/l,m/l)*sum%mod)%mod;
56     }
57     printf("%lld",ans);
58     return 0;
59 }

[BZOJ2154]Crash的數字表格