题目传送门 http://www.51nod.com/onlineJudge/questionCode.html#!problemId=1237
数学题真是做的又爽又痛苦，爽在于只要推出来公式基本上就是AC，痛苦就在于推公式。。。
题意很简单，求
$Largesumlimits_{i=1}^{n}sumlimits_{j=1}^{n}gcd(i,j)$
其中$nle 10^{10}$
这个题有很多做法，除了普及组的$O(n^2log n)$做法，还有用莫比乌斯反演+分块优化的$O(n)$做法
然而因为这个题两个维度的限制是相等的，都是$n$，所以可以用杜教筛做到$O(n^{frac{2}{3}})$
然后推一波公式
$Largesumlimits_{i=1}^{n}sumlimits_{j=1}^{n}gcd(i,j)$
按照套路，可以枚举最大公因数$g$，于是有
$Largesumlimits_{g=1}^{n}sumlimits_{g=gcd(i,j)}g$
继续变形
$Largesumlimits_{g=1}^{n}gcdot(sumlimits_{g=gcd(i,j)}1)$
$Largesumlimits_{g=1}^{n}gcdot(sumlimits_{1=gcd(i,j)}^{ilelfloorfrac{n}{g} floor,jlelfloorfrac{n}{g} floor}1)$
到这里，我们就可以用莫比乌斯反演的套路做到$O(n)$的复杂度了，但是我们换种形式继续推式子
我们先关注这一部分的变形，暂时不管前面的那一部分
$Largesumlimits_{1=gcd(i,j)}^{ilelfloorfrac{n}{g} floor,jlelfloorfrac{n}{g} floor}1$
展开，得到
$Largesumlimits_{i=1}^{lfloorfrac{n}{g} floor}sumlimits_{j=1}^{lfloorfrac{n}{g} floor}[gcd(i,j)=1]cdot 1$
其中$[gcd(i,j)=1]$的意思是，当$gcd(i,j)=1$的时候这个东西的值为$1$，否则为$0$
然后用$varphi()$函数进行化简，于是式子就变成了
$Large 2cdot(sumlimits_{i=1}^{lfloorfrac{n}{g} floor}varphi(i))-1$
注意到中间的一部分是$varphi()$函数的前缀和，于是我们可以用杜教筛算
为了方便，设函数
$Large S(n)=sumlimits_{i=1}^{n}varphi(i)$
然后再看一眼总的式子
$Largesumlimits_{g=1}^{n}gcdot(2cdot S(lfloorfrac{n}{g} floor)-1)$
于是$lfloorfrac{n}{g} floor$是整除，可以用分块优化，这样的话，除去杜教筛求$S()$的部分，复杂度是$O(n^{frac{1}{2}})$
因为杜教筛有记忆化，复杂度和分块优化的部分是分离的，所以总复杂度是
$Large O(n^{frac{1}{2}}+n^{frac{2}{3}})=O(n^{frac{2}{3}})$
单点时限$5s$，我杜教筛记忆化用的$map$，没有自己手写$hash$，于是复杂度就多了$O(log n)$，但是仍然每个点跑到了$1s$以内，还是挺不错的
上代码：
1 #include <cstdio> 2 #include <algorithm> 3 #include <cstring> 4 #include <map> 5 6 using namespace std; 7 typedef long long ll; 8 const int MOD = 1e9+7; 9 const int LIMIT = 4641589; 10 11 bool vis[LIMIT] = {0}; 12 ll prime[LIMIT], pidx = 0, phi[LIMIT] = {0}, pfphi[LIMIT] = {0}; 13 void prelude_phi() { // 线性筛预处理一部分前缀和 14 phi[1] = 1; 15 for( ll i = 2; i < LIMIT; ++i ) { 16 if( !vis[i] ) { 17 prime[pidx++] = i; 18 phi[i] = i-1; 19 } 20 for( int j = 0; j < pidx; ++j ) { 21 ll k = i * prime[j]; 22 if( k >= LIMIT ) break; 23 vis[k] = true; 24 if( i % prime[j] ) phi[k] = phi[i] * phi[prime[j]] % MOD; 25 else { 26 phi[k] = phi[i] * prime[j] % MOD; 27 break; 28 } 29 } 30 } 31 for( int i = 1; i < LIMIT; ++i ) 32 pfphi[i] = (pfphi[i-1] + phi[i]) % MOD; 33 } 34 35 ll inv2; 36 ll pow_mod( ll a, ll b ) { 37 if( !b ) return 1; 38 ll rtn = pow_mod(a,b>>1); 39 rtn = rtn * rtn % MOD; 40 if( b&1 ) rtn = rtn * a % MOD; 41 return rtn; 42 } 43 ll inv( ll x ) { 44 return pow_mod(x,MOD-2); 45 } 46 47 map<ll,ll> mp; // 记忆化用 48 ll S( ll n ) { // 杜教筛求前缀和 49 if( n < LIMIT ) return pfphi[n]; 50 if( mp.count(n) ) return mp[n]; 51 ll ans = (n % MOD) * ((n+1) % MOD) % MOD * inv2 % MOD; 52 for( ll i = 2, j; i <= n; i = j+1 ) { 53 j = n/(n/i); 54 ans = (ans - S(n/i) * ((j-i+1) % MOD) % MOD + MOD) % MOD; 55 } 56 mp[n] = ans; 57 return ans; 58 } 59 60 ll n; 61 62 void solve() { 63 ll ans = 0; 64 for( ll i = 1, j; i <= n; i = j+1 ) { 65 j = n/(n/i); // 分块优化 66 ans = (ans + (((i+j) % MOD) * ((j-i+1) % MOD) % MOD * inv2 % MOD) * ((2 * S(n/i) - 1) % MOD) % MOD ) % MOD; 67 } 68 printf( "%lld ", ans ); 69 } 70 71 int main() { 72 prelude_phi(); 73 inv2 = inv(2); 74 scanf( "%lld", &n ); 75 solve(); 76 return 0; 77 }