$$求解\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=k]$$ $反演过程:$ $$\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=k]$$
$$\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{m}{k} \right \rfloor}\varepsilon[gcd(i,j)=1]$$
$枚举d:$
$$\sum_{d=1}^{\left \lfloor min(\left \lfloor \frac{n}{k} \right \rfloor,\left \lfloor \frac{m}{k} \right \rfloor)\right \rfloor}\mu(d)\left \lfloor \frac{n}{kd} \right \rfloor \left \lfloor \frac{m}{kd} \right \rfloor$$
$\ {先预处理\mu,然后整数分块做,复杂度为O(n + T\sqrt n)}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 5e5  + 10 ;int  mu[N]; bool  is \_prime[N];int  prime[N];int  cnt;ll sum[N];  void  Mobi () {     mu[1 ] = 1 ;     is \_prime[0 ] = is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++) {         if  (!is \_prime[i]) {             mu[i] = -1 ;             prime[++cnt] = i;         }         for  (int  j = 1 ; j <= cnt && i \* prime[j] < N; j++) {             is \_prime[i \* prime[j]] = true ;             if  (i % prime[j] == 0 ) {                 mu[i \* prime[j]] = 0 ;                 break ;             }             mu[i \* prime[j]] = -mu[i];         }     }     for (int  i = 1 ;i < N; i++) {         sum[i] = sum[i - 1 ] + mu[i];     } } ll k; ll Ans (ll n, ll m )  {    ll ans = 0 ;     n /= k, m /= k;     int  mx = min(n, m);     for (int  l = 1 , r;l <= mx; l = r + 1 ) {         r = min(n / (n / l), m / (m / l));         ans += (n / l) \* (m / l) \* (sum[r] - sum[l - 1 ]);     }     return  ans; } void  solve ()    Mobi();     int  T;     cin >> T;     while (T--) {         ll n, m;         cin >> n >> m >> k;         cout << Ans(n, m) << endl;     } } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{p\in prime}\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=p]$$ $反演过程和上面几乎一样,就多了一维p为质数的情况,直接化简为:$
$$\sum_{p\in prime}\sum_{d=1}^{\left \lfloor min(\left \lfloor \frac{n}{p} \right \rfloor,\left \lfloor \frac{m}{p} \right \rfloor)\right \rfloor}\mu(d)\left \lfloor \frac{n}{pd} \right \rfloor \left \lfloor \frac{m}{pd} \right \rfloor$$
$这里的复杂度为O(n+T\sqrt n \sqrt n),但明显不够,所以需要继续化简,这里的优化可以用T替换pd,则d=\frac{T}{p},替换后的式子为:$
$$\sum_{T=1}^{min(n,m)})\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor\sum_{pT\;p\in prime}\mu(\frac{T}{p})$$
$$令\ {F(x)=\sum_{px\;p\in prime}\mu(\frac{x}{p})}$$
$得:$
$$\sum_{p\in prime}\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=p]=\sum_{T=1}^{min(n,m)}\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor F(T)$$
$\ {O(n)处理\mu,O(nlogn)处理F,O(\sqrt n)询问分块,即复杂度为O(nlogn + T\sqrt n)}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 1e7  + 10 ;int  mu[N]; bool  is \_prime[N];int  prime[N];int  cnt;ll sum[N];  ll F[N]; void  Mobi () {     mu[1 ] = 1 ;     is \_prime[0 ] = is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++) {         if  (!is \_prime[i]) {             mu[i] = -1 ;             prime[++cnt] = i;         }         for  (int  j = 1 ; j <= cnt && i \* prime[j] < N; j++) {             is \_prime[i \* prime[j]] = true ;             if  (i % prime[j] == 0 ) {                 mu[i \* prime[j]] = 0 ;                 break ;             }             mu[i \* prime[j]] = -mu[i];         }     }     for (int  i = 1 ;i <= cnt; i++) {         for (int  j = 1 ;j \* prime[i] < N; j++) {             F[j \* prime[i]] += mu[j];         }     }     for (int  i = 1 ;i < N; i++) {         F[i] += F[i - 1 ];     } } ll k; ll Ans (ll N, ll M )  {    ll ans = 0 ;     ll n = N, m = M;     int  mx = min(n, m);     for  (int  l = 1 , r; l <= mx; l = r + 1 ) {         r = min(n / (n / l), m / (m / l));         ans += (n / l) \* (m / l) \* (F[r] - F[l - 1 ]);     }     return  ans; } void  solve ()    Mobi();     int  T;     cin >> T;     while (T--) {         ll n, m;         cin >> n >> m;         cout << Ans(n, m) << endl;     } } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=1}^n\sum_{j=i+1}^ngcd(i,j)$$
$先来看看基本形式:$ $$\sum_{i=1}^n\sum_{j=1}^ngcd(i,j)$$ $则\ {ans=\frac{Calc(n,n)-\frac{n(n+1)}{2}}{2}},首先减去gcd(i,i)的个数\frac{n(n+1)}{2},然后减去gcd(i,j)=gcd(j,i)的个数,直接除2即可,则:$ $$\sum_{i=1}^n\sum_{j=1}^n\sum_{k=1}^nk[gcd(i,j)=k]$$
$枚举k:$
$$\sum_{k=1}^nk\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=\left \lfloor \frac{i+1}{k} \right \rfloor}^{\left \lfloor \frac{n}{k} \right \rfloor}[gcd(i,j)=1]$$
$$\sum_{k=1}^nk\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{di\;dj}\mu(d)$$
$$\sum_{k=1}^nk\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}1$$
$$\sum_{k=1}^nk\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)(\left \lfloor \frac{n}{kd} \right \rfloor)^2$$
$令\ {T=kd并枚举T},得:$
$$\sum_{T=1}^n\sum_{dT}\mu(d)\frac{T}{d}(\left \lfloor \frac{n}{T} \right \rfloor)^2$$
$中间部分有\ {\sum_{dT}\mu(d)\frac{T}{d}=\varphi (T),因为\varphi =\mu * id},即得:$
$$\sum_{T=1}^n\varphi (T)(\left \lfloor \frac{n}{T} \right \rfloor)^2$$
$则最终答案为:$
$$\ {ans=\frac{\sum_{T=1}^n\varphi (T)(\left \lfloor \frac{n}{T} \right \rfloor)^2-\frac{n(n+1)}{2}}{2}}$$
$\ {O(n)预处理\varphi,\sqrt n询问分块,总时间复杂度为O(n+\sqrt n)}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 2e6  + 10 ;int  phi[N];bool  is \_prime[N];int  prime[N];int  tot;ll f[N]; void  Euler (){     phi[1 ] = 1 ;     is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++){         if (!is \_prime[i]) {             phi[i] = i - 1 ;             prime[++tot] = i;         }         for (int  j = 1 ;j <= tot && i \* prime[j] < N; j++){             is \_prime[i \* prime[j]] = true ;             if (i % prime[j]) {                 phi[i \* prime[j]] = phi[i] \* (prime[j] - 1 );             }             else {                 phi[i \* prime[j]] = phi[i] \* prime[j];                 break ;             }         }     }     for (int  i = 1 ;i < N; i++) {         f[i] = f[i - 1 ] + phi[i];     } } ll Calc (ll n )  {    ll ans = 0 ;     for (ll l = 1 , r;l <= n;l = r + 1 ) {         r = min(n, n / (n / l));         ans += (f[r] - f[l - 1 ]) \* (n / l) \* (n / l);     }     return  ans; } void  solve ()    Euler();     ll n;     cin >> n;     cout << (Calc(n) - n \* (n + 1 ) / 2 ) / 2  << endl; } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=1}^n\sum_{j=1}^mlcm (i,j)$$
$由gcd和lcm之间的关系得:$
$$\sum_{i=1}^n\sum_{j=1}^m\frac{ij}{gcd(i,j)}$$
$$\sum_{k=1}^{min(n,m)}\sum_{i=1}^n\sum_{j=1}^m\frac{ij}{k}[gcd(i,j)=k]$$
$$\sum_{k=1}^{min(n,m)}\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{m}{k} \right \rfloor}\frac{ijk^2}{k}[gcd(i,j)=1]$$
$$\sum_{k=1}^{min(n,m)}k\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{m}{k} \right \rfloor}ij\sum_{di\;dj}\mu(d)$$
$枚举d得:$
$$\sum_{k=1}^{min(n,m)}k\sum_{d=1}^{min(\left \lfloor \frac{n}{k} \right \rfloor,\left \lfloor \frac{m}{k} \right \rfloor)}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}id\sum_{j=1}^{\left \lfloor \frac{m}{kd} \right \rfloor}jd$$
$$\sum_{k=1}^{min(n,m)}k\sum_{d=1}^{min(\left \lfloor \frac{n}{k} \right \rfloor,\left \lfloor \frac{m}{k} \right \rfloor)}\mu(d)d^2\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}i\sum_{j=1}^{\left \lfloor \frac{m}{kd} \right \rfloor}j$$
$令\ {f(x)=\sum_{i=1}^{x}i,F(n,m)=\sum_{d=1}^{min(n,m)}\mu(d)d^2f(\frac{n}{d})f(\frac{m}{d})},得:$
$$\sum_{k=1}^{min(n,m)}kF(\left \lfloor\frac{n}{k}\right \rfloor,\left \lfloor\frac{m}{k}\right \rfloor)$$
$\ {O(n)处理\mu,O(1)处理f(等差数列,在计算F的过程中直接直接计算),O(\sqrt {\sqrt n})分块询问F,答案分块询问O(\sqrt n),总复杂度为O(n+n^{\frac{3}{4}})}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 1e7  + 10 ;const  ll mod = 20101009 ;int  mu[N]; bool  is \_prime[N];int  prime[N];int  cnt;ll sum[N]; void  Mobi () {     mu[1 ] = 1 ;     is \_prime[0 ] = is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++) {         if  (!is \_prime[i]) {             mu[i] = -1 ;             prime[++cnt] = i;         }         for  (int  j = 1 ; j <= cnt && i \* prime[j] < N; j++) {             is \_prime[i \* prime[j]] = true ;             if  (i % prime[j] == 0 ) {                 mu[i \* prime[j]] = 0 ;                 break ;             }             mu[i \* prime[j]] = -mu[i];         }     }     for (int  i = 1 ;i < N; i++) {         sum[i] = (sum[i - 1 ] + mu[i] \* i % mod \* i % mod) % mod;     } } ll Ans (ll n, ll m ) {     ll inv2 = (mod + 1 ) / 2 ;     ll ans = 0 ;     ll p = min(n, m);     for (int  k = 1 ;k <= p; k++) {         int  x = n / k, y = m / k, P = min(x, y);         ll Sum = 0 ;         for (int  l = 1 , r;l <= P; l = r + 1 ) {             r = min(x / (x / l), (y / (y / l)));             Sum = (Sum + (sum[r] - sum[l - 1 ] + mod) % mod \* (1  + x / l) % mod \* (x / l) % mod \* inv2 % mod \* (1  + y / l) % mod \* (y / l) % mod \* inv2 % mod) % mod;         }         ans = (ans + Sum \* k % mod) % mod;     }     return  ans; } void  solve ()    Mobi();     ll n, m;     cin >> n >> m;     cout << Ans(n, m) << endl; } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=1}^n\sum_{j=1}^nij(gcd(i,j))\;mod \;p$$ $$\sum_{k=1}^n\sum_{i=1}^n\sum_{j=1}^nijk[gcd(i,j)=k]$$
$$\sum_{k=1}^nk\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}ijk^2[gcd(i,j)=1]$$
$$\sum_{k=1}^nk^3\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}ij\sum_{di\;dj}\mu(d)$$
$$\sum_{k=1}^nk^3\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}id\sum_{j=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}jd$$
$$\sum_{k=1}^nk^3\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)d^2\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}i\sum_{j=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}j$$
$令\ {f(x)=\sum_{i=1}^xi},则式子变为:$
$$\sum_{k=1}^nk^3\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)d^2f^2(\frac{n}{kd})$$
$到这里,时间复杂度为O(n+n\sqrt n),我没尝试过,但是还可以继续化简,化简方法有替换,用T替换kd,则:$
$$\sum_{T=1}^nf^2(\frac{n}{T})T^2\sum_{kT}k*\mu(\frac{T}{k})$$
$$\sum_{T=1}^nf^2(\frac{n}{T})T^2\ {(id*\mu)T}$$
$$\sum_{T=1}^nf^2(\frac{n}{T})T^2\ {\varphi (T)}$$
$\ {可以杜教筛筛T^2\varphi (T)的前缀和,因为n最高有10^{10},明显数组存不下,所以对于5e6可以用数组直接过度求前缀和,对于后面就可以用杜教筛求前缀和,然后用map存下(STL的好处),过程如下:}$
$参考:$杜教筛 
$先看杜教筛的基本形式(与上面无关):设S(n)=\sum_{i=1}^{n}f(i),定义一个新函数为g,则有$
$$\sum_{i=1}^n(f*g)i$$
$$\sum_{i=1}^n\sum_{xy=i}f(x)g(y)$$
$$\sum_{y=1}^ng(y)\sum_{xy\leq n}f(x)$$
$$\sum_{y=1}^ng(y)\sum_{x=1}^{\left \lfloor \frac{n}{y} \right \rfloor}f(x)$$
$$\sum_{y=1}^ng(y)S(\left \lfloor \frac{n}{y} \right \rfloor)$$
$则有:$
$$\sum_{i=1}^n(f*g)i=g(1)S(n)+\sum_{y=2}^ng(y)S(\left \lfloor \frac{n}{y} \right \rfloor)$$
$$\ {S(n)=\frac{\sum_{i=1}^n(f*g)i-\sum_{y=2}^ng(y)S(\left \lfloor \frac{n}{y} \right \rfloor)}{g(1)}}$$
$对于本题,要找到合适的g函数,因为我们求的前缀和为S(n)=\sum_{i=1}^ni^2\varphi(i),f(i)=i^2\varphi(i),所以设g(n)=n^2。$
$\sum_{i=1}^ng(i)=\frac{x(x+1)(2x+1)}{6}$
$(f*g)x=\sum_{dx}f(d)g(\frac{x}{d})=\sum_{dx}d^2\varphi(d)\frac{x^2}{d^2}=x^2\sum_{dx}\varphi(d)=x^2id(x)=x^3$
$\sum_{i=1}^n(f*g)i=\sum_{i=1}^ni^3=\frac{x^2(x+1)^2}{4}$
$$\ {S(n)=\frac{n^2(n+1)^2}{4}-\sum_{y=2}^ng(y)S(\left \lfloor \frac{n}{y} \right \rfloor)}$$
$\ {O(n^{\frac{2}{3}})处理T^2\varphi (T)的前缀和,f可以O(1)算出,O(\sqrt n) 询问分块,总时间复杂度为O(n^\frac{2}{3})}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 5e6  + 10 ;ll n, mod, inv4, inv6; ll phi[N]; bool  is \_prime[N];int  prime[N];int  tot;ll sum[N]; void  Get\_phi() {    phi[1 ] = 1 ;     is \_prime[1 ] = true ;     for (int  i = 2 ;i <= N; i++){         if (!is \_prime[i]) {             phi[i] = i - 1 ;             prime[++tot] = i;         }         for (int  j = 1 ;j <= tot && i \* prime[j] < N; j++){             is \_prime[i \* prime[j]] = true ;             if (i % prime[j]) {                 phi[i \* prime[j]] = phi[i] \* (prime[j] - 1 );             }             else {                 phi[i \* prime[j]] = phi[i] \* prime[j];                 break ;             }         }     }     for (int  i = 1 ;i < N; i++) {         sum[i] = (sum[i - 1 ] + (ll)phi[i] \* i % mod \* i % mod) % mod;     } } ll quick\_pow(ll a, ll b) {     ll ans = 1 ;     while (b){         if (b & 1 ) ans = ans \* a % mod;         a = a \* a % mod;         b >>= 1 ;     }     return  ans; } ll sum1 (ll x )  {x %= mod; return  (1  + x) % mod \* x % mod \* (1  + x) % mod \* x % mod \* inv4 % mod;}ll sum2 (ll x )  {x %= mod; return  x % mod \* (1  + x) % mod \* (2  \* x + 1 ) % mod \* inv6 % mod;}map<ll, ll> mp; ll Calc (ll x )  {    if (x < N) return  sum[x];     if (mp[x]) return  mp[x];     ll num = sum1(x);     for (ll l = 2 , r;l <= x;l = r + 1 ) {         r = x / (x / l);         num = (num - Calc(x / l) % mod \* (sum2(r) - sum2(l - 1 ) + mod) % mod + mod) % mod;     }     return  mp[x] = num; } ll Ans ()  {    ll ans = 0 ;     for (ll l = 1 , r;l <= n; l = r + 1 ) {         r = n / (n / l);         ans = (ans + sum1(n / l) % mod \* (Calc(r) - Calc(l - 1 ) % mod + mod) % mod) % mod;     }     return  ans; } void  solve ()    cin >> mod >> n;     Get\_phi();     inv4 = quick\_pow(4 , mod - 2 );     inv6 = quick\_pow(6 , mod - 2 );     cout << Ans() << endl; } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=1}^n\sum_{j=1}^md (ij)$$
$$\ {约数的重要性质:(ij)=\sum_{xi}\sum_{yj}[gcd(x,y )=1]}$$
https://www.cnblogs.com/sun123zxy/p/12295533.html 
$得:$ $$\sum_{i=1}^n\sum_{j=1}^m\sum_{xi}\sum_{yj}[gcd(x,y )=1]$$
$改变枚举顺序,枚举x和y:$
$$\sum_{x=1}^n\sum_{y=1}^m\left \lfloor \frac{n}{x} \right \rfloor\left \lfloor \frac{n}{y} \right \rfloor[gcd(x,y)=1]$$
$$\sum_{i=1}^n\sum_{j=1}^m\left \lfloor \frac{n}{i} \right \rfloor\left \lfloor \frac{n}{j} \right \rfloor[gcd(i,j)=1]$$
$后面套用前面得:$
$$\sum_{i=1}^n\sum_{j=1}^m\left \lfloor \frac{n}{i} \right \rfloor\left \lfloor \frac{n}{j} \right \rfloor\sum_{di\;dj}\mu(d)$$
$枚举d得:$
$$\sum_{d=1}^{min(n,m)}\mu(d)\sum_{i=1}^n\sum_{j=1}^m\left \lfloor \frac{n}{i} \right \rfloor\left \lfloor \frac{n}{j} \right \rfloor[dgcd(i,j)]$$
$由于[dgcd(i,j)]成立的条件为d是gcd(i,j)的约数,所以可以把枚举i,j变换为枚举di,dj,从而[1gcd(i,j)]可以省略。得:$
$$\sum_{d=1}^{min(n,m)}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor}\sum_{j=1}^{\left \lfloor\frac{m}{d}\right \rfloor}\left \lfloor \frac{n}{id} \right \rfloor\left \lfloor \frac{n}{jd} \right \rfloor$$
$$\sum_{d=1}^{min(n,m)}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor}\left \lfloor \frac{n}{id} \right \rfloor\sum_{j=1}^{\left \lfloor\frac{m}{d}\right \rfloor}\left \lfloor \frac{n}{jd} \right \rfloor$$
$令\ {f(x)=\sum_{i=1}^x\left \lfloor \frac{x}{i} \right \rfloor},则:$
$$\sum_{d=1}^{min(n,m)}\mu(d)f(\frac{n}{d})f(\frac{m}{d})$$
$\ {O(n)处理\mu和其前缀和,O(n\sqrt n)处理f,O(\sqrt n)询问分块,总复杂度为O(n\sqrt n+T\sqrt n)}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 5e4  + 10 ;ll mu[N];  bool  is \_prime[N];int  prime[N];int  cnt;ll f[N]; void  Mobi () {     mu[1 ] = 1 ;     is \_prime[0 ] = is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++) {         if  (!is \_prime[i]) {             mu[i] = -1 ;             prime[++cnt] = i;         }         for  (int  j = 1 ; j <= cnt && i \* prime[j] < N; j++) {             is \_prime[i \* prime[j]] = true ;             if  (i % prime[j] == 0 ) {                 mu[i \* prime[j]] = 0 ;                 break ;             }             mu[i \* prime[j]] = -mu[i];         }     }     for (int  i = 1 ;i < N; i++) {         mu[i] += mu[i - 1 ];     }     for (int  i = 1 ;i < N; i++) {         ll ans = 0 ;         for (int  l = 1 , r;l <= i; l = r + 1 ) {             r = i / (i / l);             ans += (r - l + 1 ) \* (i / l);         }         f[i] = ans;     } } ll Ans (ll n, ll m ) {     ll ans = 0 ;     ll k = min(n, m);     for (ll l = 1 , r;l <= k; l = r + 1 ) {         r = min(n / (n / l), m / (m / l));         ans += (mu[r] - mu[l - 1 ]) \* f[m / l] \* f[n / l];     }     return  ans; } void  solve ()    int  T;     Mobi();     cin >> T;     while (T--) {         ll n, m;         cin >> n >> m;         cout << Ans(n, m) << endl;     } } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=1}^n\sum_{j=1}^m\sum_{di\;dj}d[\sum_{di\;dj}d\leq a]$$
$令\ {F(x)=\sum_{ix}i},则得:$
$$\sum_{i=1}^n\sum_{j=1}^m F(gcd(i,j))[F(gcd(i,j))\leq a]$$
$$\sum_{d=1}^{min(n,m)}\sum_{i=1}^n\sum_{j=1}^m F(d)[gcd(i,j)=d][F(d)\leq a]$$
$$\sum_{d=1}^{min(n,m)}F(d)[F(d)\leq a]\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=d]$$
$后半部分就套用前面得:$
$$\sum_{d=1}^{min(n,m)}F(d)[F(d)\leq a]\sum_{d^{‘}=1}^{min(\left \lfloor \frac{n}{d} \right \rfloor, \left \lfloor \frac{m}{d} \right \rfloor)}\mu(d^{‘})\left \lfloor \frac{n}{dd^{‘}} \right \rfloor \left \lfloor \frac{m}{dd^{‘}} \right \rfloor$$
$令T替换dd^{‘}得:$
$$\sum_{T=1}^{min(n,m)}\sum_{dT}F(d)\mu(\frac{T}{d})\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor[F(d)\leq a]$$
$$\sum_{T=1}^{min(n,m)}\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor\sum_{dT}F(d)\mu(\frac{T}{d})[F(d)\leq a]$$
$令\ {f(x)=\sum_{dx}F(d)\mu(\frac{x}{d})[F(d)\leq a]}得:$
$$\sum_{T=1}^{min(n,m)}\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor f(T)$$
$\ {O(n)处理\mu,O(nlogn)处理F,对于每次询问的a值可以离线排序,树状数组维护f,处理时间为O(nlog^{2}n),O(\sqrt n\; logn)询问分块,总时间复杂度为O(nlog^2n)+T\sqrt n\; logn}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x3f3f3f3f ##define  lc u << 1 ##define  rc u << 1  1 ##define  m (l + r) / 2 ##define  mid (t[u].l + t[u].r) / 2 ##define  lowbit(x) x & (-x) ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 1e5  + 1 ;ll mod = (1l l \* 1 ) << 31 ; struct  node {    int  id;     ll data; }F[N]; struct  Node {    int  nn, mm, id;     ll a; }q[N]; int  cnt, prime[N];bool  is \_prime[N];int  mu[N];bool  cmp\_F(node a, node b) {    return  a.data < b.data; } bool  cmp\_Q(Node a, Node b) {    return  a.a < b.a; } inline void  init () {     mu[1 ] = 1 ;     is \_prime[0 ] = is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++) {         if  (!is \_prime[i]) {             mu[i] = -1 ;             prime[++cnt] = i;         }         for  (int  j = 1 ; j <= cnt && i \* prime[j] < N; j++) {             is \_prime[i \* prime[j]] = true ;             if  (i % prime[j] == 0 ) {                 mu[i \* prime[j]] = 0 ;                 break ;             }             mu[i \* prime[j]] = -mu[i];         }     }     for (int  i = 1 ;i < N; i++) {         F[i].id = i;         for (int  j = 1 ;j \* i < N; j++) {             F[i \* j].data += i;         }     }     sort(F + 1 , F + N, cmp\_F); } ll t[N]; ll ans[N]; inline void  add (int  x, int  v  {    while (x < N) {         t[x] = (t[x] + v) % mod;         if (t[x] > mod)             t[x] %= mod;         x += lowbit(x);     } } inline ll query (int  x  {    ll ans = 0 ;     while (x) {         ans = (ans + t[x]) % mod;         if (ans > mod)             ans %= mod;         x -= lowbit(x);     }     return  ans; } inline ll Calc (int  nn, int  mm  {    if (nn > mm)         swap(nn, mm);     ll ans = 0 ;     for (int  l = 1 , r;l <= nn; l = r + 1 ) {         r = min(nn / (nn / l), mm / (mm / l));         ans = (ans + 1l l \* (nn / l) \* (mm / l) % mod \* (query(r) - query(l - 1 )) % mod) % mod;     }     return  ans; } void  solve ()    init ();     int  T;     scanf("%d" ,&T);     for (int  i = 1 ;i <= T; i++) {         scanf("%d%d%lld" ,&q[i].nn,&q[i].mm, &q[i].a);         q[i].id = i;     }     sort(q + 1 , q + T + 1 , cmp\_Q);     int  now = 1 ;     for (int  i = 1 ;i <= T; i++) {         while (now < N && F[now].data <= q[i].a) {             for (int  j = 1 ;j \* F[now].id < N; j++) {                 add (j \* F[now].id, mu[j] \* F[now].data);             }             now++;         }         ans[q[i].id] = (Calc(q[i].nn, q[i].mm) % mod + mod) % mod;     }     for (int  i = 1 ;i <= T; i++) {         printf("%lld\n" ,ans[i]);     } } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=L}^H\sum_{j=L}^H…\sum_{x=L}^H[gcd(i,j…x)=k]$$
$本题求n个数的gcd为k的个数,根据上面套路,将H/k,L/k,这样就只用求n个数的gcd为1的个数了,注意:如果L整除k,则L=L/k,否则L=L/k+1,就是向上取整,那么有一种办法非常好,就让\ {L/k+(\frac{k-1}{k})},即:$
$$\sum_{i={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}\sum_{j={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}…\sum_{x={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}[gcd(i,j…x)=1]$$
$$\sum_{i={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}\sum_{j={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}…\sum_{x={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}\sum_{di\;dj\;..\;dx}\mu(d)$$
$枚举d,令\ {l=\left \lfloor \frac{L-1}{k} \right \rfloor+1,r=\left \lfloor \frac{H}{k}\right \rfloor}:$
$$\sum_{d=1}^r\mu(d)\sum_{i=l}^r[di]\sum_{j=l}^r[dj]…\sum_{x=l}^r[dx]$$
$只有当i,j…x都是d的倍数的时候,最终答案+1,根据容斥原理,那么在[l,r]中d倍数有\ {\frac{r}{d}-\frac{l-1}{d}},对于所有的i,j…x,他们当中所有[l,r]的d的倍数的个数有\ {(\frac{r}{d}-\frac{l-1}{d})^n},即为:$
$$\sum_{d=1}^{r}\mu(d)*(\frac{r}{d}-\frac{l-1}{d})^n$$
$这里的H非常大,如果用线性筛去处理\mu,会直接Wa掉,因为没法处理到10^9次方那么大的前缀和,数组也开不了那么大。后果:$
$\ {设f(d)=\mu(d),S(r)=\sum_{d=1}^nf(d),这里的S就是要求的前缀和。}$ $我们选择\ {g=I},则\sum_{i=1}^r(f*g)i=\sum_{i=1}^r(\mu*I)i=\sum_{i=1}^r\varepsilon (i)=1。$
$枚举y:$
$$\sum_{y=1}^rg(y)\sum_{x=1}^{\left \lfloor \frac{r}{y} \right \rfloor)}f(x)$$
$$\sum_{y=1}^rg(y)S(\left \lfloor \frac{r}{y} \right \rfloor)$$
$$\sum_{i=1}^r(f*g)i =g(1)S(r)+\sum_{y=2}^rg(y)S(\left \lfloor \frac{r}{y} \right \rfloor)$$
$$S(r)=\frac{\sum_{i=1}^r(f*g)i-\sum_{y=2}^rg(y)S(\left \lfloor \frac{r}{y} \right \rfloor)}{g(1)}$$
$当g=I时,即:$
$$\ {S(r)=1-\sum_{y=2}^rS(\left \lfloor \frac{r}{y} \right \rfloor)}$$
$\ {杜教筛的复杂度为O(n\frac{2}{3}),比线性筛还要快一点,总时间复杂度为O(n^{\frac{2}{3}})。}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  ll mod = 1e9  + 7 ;const  int  N = 1e6  + 10 ;int  mu[N]; bool  is \_prime[N];int  prime[N];int  cnt;ll sum[N]; ll n, L, H, k; ll quick\_pow(ll a, ll b) {     ll ans = 1 ;     while (b) {         if (b & 1 ) ans = ans \* a % mod;         a = a \* a % mod;         b >>= 1 ;     }     return  ans; } void  Mobi () {     mu[1 ] = 1 ;     is \_prime[0 ] = is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++) {         if  (!is \_prime[i]) {             mu[i] = -1 ;             prime[++cnt] = i;         }         for  (int  j = 1 ; j <= cnt && i \* prime[j] < N; j++) {             is \_prime[i \* prime[j]] = true ;             if  (i % prime[j] == 0 ) {                 mu[i \* prime[j]] = 0 ;                 break ;             }             mu[i \* prime[j]] = -mu[i];         }     }     for (int  i = 1 ;i < N; i++) {         sum[i] = (sum[i - 1 ] + mu[i]) % mod;     } } map<ll , ll> mp; ll Calc (ll x )  {    if (x < N)         return  sum[x];     if (mp[x])         return  mp[x];     ll num = 1 ;     for (int  l = 2 , r;l <= x; l = r + 1 ) {         r = x / (x / l);         num = (num - Calc(x / l) \* (r - l + 1 ) % mod + mod) % mod;     }     return  mp[x] = num; } ll Ans ()  {    ll ans = 0 ;     for (ll l = 1 , r;l <= H;l = r + 1 ) {         r = min(H / (H / l), (L / l) ? L / (L / l) : H + 2 );         ans = (ans + (Calc(r) - Calc(l - 1 ) + mod) % mod \* quick\_pow((H / l) - (L / l), n) % mod + mod) % mod;     }     return  (ans + mod) % mod; } void  solve ()    cin >> n >> k >> L >> H;     L = (L - 1 ) / k;      H = H / k;     Mobi();     cout << Ans() << endl; } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=1}^{N}\sum_{j=1}^Nlcm(A_i,A_j)$$
$莫比乌斯反演处理的一般都是枚举变量之间的关系,所以我们要让A_i和A_j变为i和j的形式。$
$想要用i和j表示A_i和A_j并且不遗漏任何的A,所以让A等于某个i,就是让cnt_i=\ {\sum_{d=1}^N[A_d=i]},统计所有A等于某个i的个数,即cnt_i,j也同样如此,如下所示:$
$$\sum_{i=1}^{N}\sum_{j=1}^Ncnt_icnt_jlcm(i,j)$$
$这样还不够,因为i只能是[1,N]的,当Ai中有些值大于N的话,再枚举i从[1,N]就会漏掉大于N的Ai的情况,所以想要统计所有情况,即让N=Ai中最大的数即可,即\ {M=\max_{i=1}^N }。$ $所以我们需要反演的式子为:$
$$\sum_{i=1}^{M}\sum_{j=1}^Mcnt_i ⋅cnt_j ⋅lcm(i,j)$$
$$\sum_{i=1}^{M}\sum_{j=1}^Mcnt_i ⋅cnt_j ⋅\frac{ij}{gcd(i,j)}$$
$$\sum_{k=1}^M\sum_{i=1}^{M}\sum_{j=1}^Mcnt_i ⋅cnt_j ⋅\frac{ij}{k[gcd(i,j)=k]}$$
$$\sum_{k=1}^M\sum_{i=1}^{\left \lfloor \frac{M}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{M}{k} \right \rfloor}cnt_{ik} ⋅cnt_{jk} ⋅\frac{ijk^2}{k}[gcd(i,j)=1]$$
$$\sum_{k=1}^Mk\sum_{i=1}^{\left \lfloor \frac{M}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{M}{k} \right \rfloor}cnt_{ik} ⋅cnt_{jk} ⋅i ⋅j\sum_{di\;dj}\mu(d)$$
$$\sum_{k=1}^Mk\sum_{d=1}^{\left \lfloor \frac{M}{k} \right \rfloor}\mu(d)\sum_{i=1}^{\left \lfloor \frac{M}{kd} \right \rfloor}cnt_{ikd} ⋅id\sum_{j=1}^{\left \lfloor \frac{M}{kd} \right \rfloor}cnt_{jkd} ⋅jd$$
$$\sum_{k=1}^Mk\sum_{d=1}^{\left \lfloor \frac{M}{k} \right \rfloor}\mu(d)d^2(\sum_{i=1}^{\left \lfloor \frac{M}{kd} \right \rfloor}cnt_{ikd} ⋅i)^2$$
$令T=kd则:$
$$\sum_{T=1}^M(\sum_{i=1}^{\left \lfloor \frac{M}{T} \right \rfloor}cnt_{iT} ⋅i)^2\sum_{kT}\mu(k)k^2\frac{T}{k}$$
$$\sum_{T=1}^MT(\sum_{i=1}^{\left \lfloor \frac{M}{T} \right \rfloor}cnt_{iT} ⋅i)^2\sum_{kT}\mu(k)k$$
$令\ {f(x)=\sum_{ix}\mu(i)i\;,\;g(x)=\sum_{i=1}^{\left \lfloor \frac{M}{x} \right \rfloor}cnt_{ix} ⋅i=\frac{1}{x}\sum_{xt}tcnt_t},则最后的式子为:$
$$\sum_{T=1}^MTg^2(T)f(T)$$
$\ {O(n)处理A的次数,O(nlogn)处理f和g,总时间复杂度为O(nlogn)}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 5e5  + 10 ;int  mu[N]; bool  is \_prime[N];int  prime[N];int  tot;ll f[N], g[N]; int  M, n;int  cnt[N];void  Mobi () {     mu[1 ] = 1 ;     is \_prime[0 ] = is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++) {         if  (!is \_prime[i]) {             mu[i] = -1 ;             prime[++tot] = i;         }         for  (int  j = 1 ; j <= tot && i \* prime[j] < N; j++) {             is \_prime[i \* prime[j]] = true ;             if  (i % prime[j] == 0 ) {                 mu[i \* prime[j]] = 0 ;                 break ;             }             mu[i \* prime[j]] = -mu[i];         }     }     for (int  i = 1 ;i <= M; i++) {         for (int  x = i;x <= M; x += i) {             f[x] += mu[i] \* i;         }     }     for (int  x = 1 ;x <= M; x++) {         for (int  t = x;t <= M; t += x) {             g[x] += t \* cnt[t];         }         g[x] /= x;     } } void  solve ()    cin >> n;     M = n;     for (int  i = 1 ;i <= n; i++) {         int  x;         cin >> x;         cnt[x]++;         M = max(M, x);     }     Mobi();     ll ans = 0 ;     for (int  i = 1 ;i <= M; i++) {         ans += i \* g[i] \* g[i] \* f[i];     }     cout << ans << endl; } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=1}^nlcm(i,n)$$
$\left \lfloor \frac{n}{k} \right \rfloor$
$$\sum_{i=1}^n\frac{in}{gcd(i,n)}$$
$$\sum_{k=1}^n\sum_{i=1}^n\frac{in}{k[gcd(i,n)=k]}$$
$$\sum_{k=1}^n\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}in[gcd(i,\frac{n}{k})=1]$$
$$n\sum_{k=1}^n\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}i[gcd(i,\frac{n}{k})=1]$$
$$n\sum_{kn}\sum_{i=1}^ki[gcd(i,k)=1]$$
$$n\sum_{kn}\sum_{i=1}^ki\sum_{di\;dk}\mu(d)$$
$枚举d:$
$$n\sum_{kn}\sum_{dk}\mu(d)\sum_{i=1}^{\left \lfloor \frac{k}{d} \right \rfloor}id$$
$$n\sum_{kn}\sum_{dk}\mu(d)d\sum_{i=1}^{\left \lfloor \frac{k}{d} \right \rfloor}i$$
$$n\sum_{kn}\sum_{dk}\mu(d)d\frac{\left \lfloor \frac{k}{d} \right \rfloor*(\left \lfloor \frac{k}{d} \right \rfloor+1)}{2}$$
$$\frac{n}{2}\sum_{kn}{\sum_{dk}\mu(d)*\frac{k^2}{d}+\sum_{dk}\mu(d)*k}$$
$$\frac{n}{2}\sum_{kn}k{\sum_{dk}\mu(d)*\frac{k}{d}+\sum_{dk}\mu(d)}$$
$\ {这时候发现两个狄利克雷卷积,\varphi (n)=\mu*id=\sum_{dn}\mu(d)*\frac{n}{d}\;\;,\;\;\varepsilon(n) =\mu*I=\sum_{dn}\mu(d)。所以我们将其替换得:}$
$$\frac{n}{2}\sum_{kn}k[\varphi(k)+\varepsilon(k)]$$
$\ {O(n)处理\varphi,\sqrt n询问答案,总复杂度为O(n+T\sqrt n)。}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 1e6  + 10 ;int  phi[N];bool  is \_prime[N];int  prime[N];int  tot;ll f[N]; void  Euler (){     phi[1 ] = 1 ;     is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++){         if (!is \_prime[i]) {             phi[i] = i - 1 ;             prime[++tot] = i;         }         for (int  j = 1 ;j <= tot && i \* prime[j] < N; j++){             is \_prime[i \* prime[j]] = true ;             if (i % prime[j]) {                 phi[i \* prime[j]] = phi[i] \* (prime[j] - 1 );             }             else {                 phi[i \* prime[j]] = phi[i] \* prime[j];                 break ;             }         }     }     for (int  i = 1 ;i < N; i++) {         for (int  j = i;j < N; j += i) {             f[j]  += 1l l \* phi[i] \* i + (i == 1 );         }     } } void  solve ()    Euler();     int  T;     cin >> T;     while (T--) {         int  n;         cin >> n;         cout << f[n] \* n / 2  << endl;     } } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\prod_{i=1}^n\prod_{j=1}^n\frac{lcm(i,j)}{gcd(i,j)}$$
$首先还是把lcm(i,j)化简得:$
$$\prod_{i=1}^n\prod_{j=1}^n\frac{ij}{[gcd(i,j)]^2}$$
$$\prod_{i=1}^n\prod_{j=1}^nij*\prod_{i=1}^n\prod_{j=1}^n\frac{1}{[gcd(i,j)]^2}$$
$\ {我们很容易推出前半部分的数值为(n!)^{2n },后半部分的值就是\prod_{i=1}^n\prod_{j=1}^n[gcd(i,j)]^2的逆元的平方。}$
$所以现在求解:$
$$\prod_{i=1}^n\prod_{j=1}^ngcd(i,j)$$
$$\prod_{k=1}^n\prod_{i=1}^n\prod_{j=1}^nk[gcd(i,j)=k]$$
$$\prod_{k=1}^nk^{\sum_{i=1}^n\sum_{j=1}^n[gcd(i,j)=k]}$$
$这个式子的变换,本来是[1,n]中所有的gcd累乘,然后我们单独把gcd(i,j)=k的拿出来,k会有[1,n]的范围,然后我们再将这些k累乘就会达到之前的效果。$ $不过,[1,n]中的gcd不单单只有一个,可能会有多个相同的gcd,所以变换之后的式子里还需要统计gcd=k出现的次数,那只需要将他们累加即可。$
$所以没毛病。继续:$
$$\prod_{k=1}^nk^{\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}[gcd(i,j)=1]}$$
$$\prod_{k=1}^nk^{\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{di\;dj}\mu(d)}$$
$枚举d:$
$$\prod_{k=1}^nk^{\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}1}$$
$$\ {\prod_{k=1}^nk^{\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)(\left \lfloor \frac{n}{kd} \right \rfloor)^2}}$$
$最终答案为:$
$$\ {(n!)^{2n}\prod_{k=1}^nk^{\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)(\left \lfloor \frac{n}{kd} \right \rfloor)^2}}$$
$\ {O(n)处理\mu的前缀和,\sqrt n询问指数的分块,可以用欧拉降幂优化,总时间复杂度为O(n+nlog_2n}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x7f7f7f ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  int  N = 1000010 ;const  ll mod = 104857601 ;bool  is \_prime[N];int  prime[N];int  mu[N]; int  cnt;int  n;ll fac = 1 ; ll quick\_pow(ll a, ll b) {     ll ans = 1 ;     while (b) {         if (b & 1 ) ans = ans \* a % mod;         a = a \* a % mod;         b >>= 1 ;     }     return  ans; } void  Mobi () {     mu[1 ] = 1 ;     is \_prime[0 ] = is \_prime[1 ] = true ;     for (int  i = 2 ;i <= n; i++) {         if  (!is \_prime[i]) {             mu[i] = -1 ;             prime[++cnt] = i;         }         for  (int  j = 1 ; j <= cnt && i \* prime[j] <= n; j++) {             is \_prime[i \* prime[j]] = true ;             if  (i % prime[j] == 0 ) {                 mu[i \* prime[j]] = 0 ;                 break ;             }             mu[i \* prime[j]] = -mu[i];         }         mu[i] = mu[i - 1 ] + mu[i];         fac = fac \* i % mod;     } } ll Ans ()  {    ll ans = 1 ;     for (int  k = 1 ;k <= n; k++) {         ll p = 0 ;         int  t = n / k;         for (int  l = 1 , r;l <= t; l = r + 1 ) {             r = min(t, t / (t / l));             p = (p + (mu[r] - mu[l - 1 ] + mod - 1 ) \* (t / l) % (mod - 1 ) \* (t / l) % (mod - 1 )) % (mod - 1 );         }         ans = ans \* quick\_pow(k, p) % mod;     }     ll k = quick\_pow(ans, mod - 2 );     fac = quick\_pow(fac, n);     return  quick\_pow(k, 2 ) \* fac % mod \* fac % mod; } void  solve ()    cin >> n;     Mobi();     cout << Ans() << endl; } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=1}^ngcd(\left \lfloor \sqrt[3]i \right \rfloor,i)$$
$我们先将这个一重和式转化为二重和式,得:$
$$\sum_{a=1}^{\left \lfloor \sqrt[3]n \right \rfloor}\sum_{i=1}^ngcd(a,i)[\left \lfloor \sqrt[3]i \right \rfloor=a]$$
$我们来看这个式子,只有当\left \lfloor \sqrt[3]i \right \rfloor=a时,才会有对答案的贡献为gcd(a,i),所以\ {a \leq \left \lfloor \sqrt[3]i \right \rfloor < a+1,a^3 \leq i \leq (a+1)^3-1},所以上式转化为:$
$$\ {\sum_{a=1}^{\left \lfloor \sqrt[3]n \right \rfloor}\sum_{i=a^3}^{min(n,(a+1)^3-1)}gcd(a,i)}$$
$然后我们将这个二重和式拆开得:$
$${\color{Brown}\sum_{a=1}^{\left \lfloor \sqrt[3]n \right \rfloor-1}\sum_{i=a^3}^{min(n,(a+1)^3-1)}gcd(a,i)+\sum_{i=(\left \lfloor \sqrt[3]n \right \rfloor)^3}^ngcd(\left \lfloor \sqrt[3]n \right \rfloor, i)}$$
$令\sqrt[3]n=N,并使a=i,i=j(好看一点)得:$
$${\color{Brown}\sum_{i=1}^{N-1}\sum_{j=i^3}^{(i+1)^3-1}gcd(i,j)+\sum_{j=N^3}^ngcd(N, j)}$$
$先计算前半部分:$
$——————————–$
$$求解\sum_{a_1=1}^n\sum_{a_2=1}^n…\sum_{a_x=1}^n\left(\prod_{j=1}^xa_j^k\right)f(gcd(a_1,a_2…a_x))gcd(a_1,a_2…a_x)$$ $令d=gcd(a_1,a_2…a_x)得:$
$$\sum_{a_1=1}^n\sum_{a_2=1}^n…\sum_{a_x=1}^n\left(\prod_{j=1}^xa_j^k\right)f(d)d$$
$枚举d得:$
$$\sum_{d=1}^ndf(d)\sum_{a_1=1}^n\sum_{a_2=1}^n…\sum_{a_x=1}^n\left(\prod_{j=1}^xa_j^k\right)[gcd(a_1,a_2…a_x)=d]$$
$$\sum_{d=1}^nd^{kx+1}f(d)\sum_{a_1=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\sum_{a_2=1}^{\left \lfloor \frac{n}{d} \right \rfloor }…\sum_{a_x=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\left(\prod_{j=1}^xa_j^k\right)[gcd(a_1,a_2…a_x)=1]$$
$$\sum_{d=1}^nd^{kx+1}f(d)\sum_{a_1=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\sum_{a_2=1}^{\left \lfloor \frac{n}{d} \right \rfloor }…\sum_{a_x=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\left(\prod_{j=1}^xa_j^k\right)\sum_{ta_1\;ta2\;…\;ta_x}\mu(t)$$
$$\sum_{d=1}^nd^{kx+1}f(d)\left (\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor }i^k\right)^x\sum_{ta_1\;ta2\;…\;ta_x}\mu(t)$$
$枚举t得:$
$$\sum_{d=1}^nd^{kx+1}f(d)\sum_{t=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\mu(t)\left (\sum_{i=1}^{\left \lfloor \frac{n}{dt} \right \rfloor }(it)^k\right)^x$$
$$\sum_{d=1}^nd^{kx+1}f(d)\sum_{t=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\mu(t)t^{kx}\left (\sum_{i=1}^{\left \lfloor \frac{n}{dt} \right \rfloor }i^k\right)^x$$
$令T=dt并枚举T得:$
$$\sum_{T=1}^n\left (\sum_{i=1}^{\left \lfloor \frac{n}{dt} \right \rfloor }i^k\right)^xT^{kx}\sum_{dT}df(d)\mu(\frac{T}{d})$$
$令g(T)=\sum_{dT}df(d)\mu(\frac{T}{d}),则$
$$\sum_{T=1}^n\left (\sum_{i=1}^{\left \lfloor \frac{n}{T} \right \rfloor }i^k\right)^xT^{kx}g(T)$$
$\ {O(nlogn)筛f和g,并对T^{kx}g(T)和\sum_{i=1}^{\left \lfloor \frac{n}{T} \right \rfloor }i^k做前缀和,最后\sqrt n分块计算,复杂度为O(nlogn+T\sqrt n)}$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 ##include <bits/stdc++.h> using  namespace  std ;typedef long  long  ll; typedef long  double  ld; typedef pair<int , int > pdd; ##define  INF 0x3f3f3f3f ##define  lowbit(x) x & (-x) ##define  mem(a, b) memset(a , b , sizeof(a)) ##define  FOR(i, x, n) for(int i = x;i <= n; i++) const  ll mod = 1e9  + 7 ;const  int  N = 2e5  + 10 ;ll T, k, x; int  mu[N]; bool  is \_prime[N];int  prime[N];int  cnt;ll g[N], f[N]; ll sumG[N], sumi[N]; ll quick\_pow(ll a, ll b) {     ll ans = 1 ;     while (b) {         if (b & 1 ) ans = ans \* a % mod;         a = a \* a % mod;         b >>= 1 ;     }     return  ans % mod; } void  Mobi () {     f[1 ] = mu[1 ] = 1 ;     is \_prime[0 ] = is \_prime[1 ] = true ;     for (int  i = 2 ;i < N; i++) {         f[i] = 1 ;         if  (!is \_prime[i]) {             mu[i] = -1 ;             prime[++cnt] = i;         }         for  (int  j = 1 ; j <= cnt && i \* prime[j] < N; j++) {             is \_prime[i \* prime[j]] = true ;             if  (i % prime[j] == 0 ) {                 mu[i \* prime[j]] = 0 ;                 break ;             }             mu[i \* prime[j]] = -mu[i];         }     }          for (int  d = 2 ;d \* d < N; d++) {         for (int  i = d \* d;i < N; i += d \* d) {             f[i] = 0 ;         }     }          for (int  d = 1 ;d < N; d++) {         for (int  i = d;i < N; i += d) {             g[i] = (g[i] + 1l l \* d \* f[d] % mod \* mu[i / d] % mod + mod) % mod;         }     }          for (int  i = 1 ;i < N; i++) {         ll t = quick\_pow(i, k);         sumi[i] = (sumi[i - 1 ] + t) % mod;         sumG[i] = (sumG[i - 1 ] + quick\_pow(t, x) \* g[i] % mod) % mod;     } } void  solve ()    cin >> T >> k >> x;     Mobi();     while (T--) {         int  n;         cin >> n;         ll ans = 0 ;         for (int  l = 1 , r;l <= n; l = r + 1 ) {             r = min(n, n / (n / l));             ans = (ans + (sumG[r] - sumG[l - 1 ] + mod) % mod \* quick\_pow(sumi[n / l], x) % mod) % mod;         }         cout << ans << endl;     } } signed main ()  {    ios\_base::sync\_with\_stdio(false );           ##ifdef FZT\_ACM\_LOCAL     freopen("in.txt" , "r" , stdin);     freopen("out.txt" , "w" , stdout);     signed test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit(0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto start\_clock\_for\_debug = clock();         solve();         auto end\_clock\_for\_debug = clock();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug)); ##else      solve(); ##endif      return  0 ; } 
$$求解\sum_{i=1}^n\sum_{j=1}^igcd(i^a-j^a,i^b-j^b)[gcd(i,j)=1]$$
$$gcd(i^a-j^a,i^b-j^b)=i^{gcd(a,b)}-j^{gcd(a,b)}$$
$题目说gcd(a,b)=1,所以那个复杂的式子直接变成i-j.$
$$\sum_{i=1}^n\sum_{j=1}^i(i-j)[gcd(i,j)=1]$$
$$\sum_{i=1}^n\sum_{j=1}^ii[gcd(i,j)=1]-\sum_{i=1}^n\sum_{j=1}^ij[gcd(i,j)=1]$$
$左边为与i互质的个数,右边为比i小并且互质的数字和,即为:$
$$\sum_{i=1}^ni\varphi (i)-\frac{\sum_{i=1}^ni\varphi (i)+1}{2}$$
$$\frac{\sum_{i=1}^ni\varphi (i)-1}{2}$$
$因为n有1e9,所以需要杜教筛优化。$
$$f*g(n)=\sum_{dn}f(d)g(\frac{n}{d})=\sum_{dn}id(d)*\varphi (d)*g(\frac{n}{d})$$
$很显然,设g=id,则:$
$$f*g(n)=\sum_{dn}\varphi (d)id(n)=n\sum_{dn}\varphi (d)=id(n)(\varphi *I)=id^2(n)=n^2$$
$$所以得:S(n)=\sum_{i=1}^n(f*g)i-\sum_{i=2}^ng(i)S(\frac{n}{i})=\frac{n(n+1)(2n+1)}{6}-\sum_{i=2}^niS(\frac{n}{i})$$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 ##include  "bits/stdc++.h"   using  namespace  std;typedef  long  long  ll;typedef  pair<int , int > pii;typedef  pair<double , double > pdd;typedef  pair<ll, ll> pll;##define  endl "\n"   ##define  eb emplace\_back  ##define  mem(a, b) memset(a , b , sizeof(a))  const  ll INF = 0x3f3f3f3f ;const  ll mod = 1e9  + 7 ;const  double  eps = 1e-6 ;const  double  PI = acos (-1 );const  double  R = 0.57721566490153286060651209 ;template <typename  T>inline  void  read (T &a)  char  c = getchar ();T x = 0 , f = 1 ;    while  (!isdigit (c)) {if  (c == '-' )f = -1 ;c = getchar ();}     while  (isdigit (c)) {x = (x << 1 ) + (x << 3 ) + c - '0' ;c = getchar ();}a = f \* x; } const  int  N = 1e6  + 10 ;bool  is\_prime[N];int  prime[N], cnt;int  phi[N];ll sum\_i\_phi[N]; void  init ()      phi[1 ] = 1 ;     for (int  i = 2 ;i < N; i++) {         if (!is\_prime[i]) prime[++cnt] = i, phi[i] = i - 1 ;         for (int  j = 1 ;j <= cnt && i \* prime[j] < N; j++) {             is\_prime[i \* prime[j]] = 1 ;             if (i % prime[j] == 0 ) {                 phi[i \* prime[j]] = phi[i] \* prime[j];                 break ;             }             else  phi[i \* prime[j]] = phi[i] \* phi[prime[j]];         }     }     for (int  i = 1 ;i < N; i++) {         sum\_i\_phi[i] = (sum\_i\_phi[i - 1 ] + 1ll  \* i \* phi[i] % mod) % mod;     } } ll quick\_pow(ll a, ll b) {     ll ans = 1 ;     while (b) {         if (b & 1 ) ans = ans \* a % mod;          a = a \* a % mod;          b >>= 1 ;     }     return  ans % mod; } ll inv2; ll inv6; map<ll, ll> mp; ll S (ll x)   {    if (x < N) return  sum\_i\_phi[x];     else  if (mp[x]) return  mp[x];     ll ans = x \* (x + 1 ) % mod \* (2  \* x + 1 ) % mod \* inv6 % mod;     for (int  l = 2 , r;l <= x; l = r + 1 ) {         r = min (x, x / (x / l));         ans = (ans - 1ll  \* (r + l) % mod \* (r - l + 1 ) % mod \* inv2 % mod \* S (x / l)) % mod;     }     return  mp[x] = ans; } void  solve ()      init ();     inv2 = quick\_pow(2 , mod - 2 );     inv6 = quick\_pow(6 , mod - 2 );     int  \_; scanf ("%d" ,&\_);     while (\_--) {         ll n, a, b; scanf ("%lld%lld%lld" ,&n,&a,&b);         ll ans = (S (n) - 1 ) % mod;         ans = ans \* inv2 % mod;         printf ("%lld\n" ,(ans % mod + mod) % mod);     } } signed  main ()      ios\_base::sync\_with\_stdio(false );           ##ifdef  FZT\_ACM\_LOCAL      freopen ("in.txt" , "r" , stdin);     freopen ("out.txt" , "w" , stdout);     signed  test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit (0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto  start\_clock\_for\_debug = clock ();         solve ();         auto  end\_clock\_for\_debug = clock ();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback (acm\_local\_for\_debug)); ##else       solve (); ##endif       return  0 ; } 
$$求解\sum_{i=1}^A\sum_{j=1}^B\sum_{k=1}^C\phi (gcd(i,j^2,k^3))\;\;mod\;\;2^{30}$$
$首先要把\phi和gcd分开,即\phi(n)=\sum_{dn}(\phi*\mu)(d),证明如下:$ $$\phi(n)=(\mu*id)(n)$$
$$\sum_{dn}\mu(d)\frac{n}{d}$$
$$\sum_{dn}\mu(d)\sum_{d’\frac{n}{d}}\phi(d’)$$
$$\sum_{dd’n}\mu(d)\phi(d’)$$
$$\sum_{dn}\sum_{d’d}\mu(d’)\phi(\frac{d}{d’})$$
$$\sum_{dn}(\phi*\mu)(d)$$
$所以原式变为:$
$$\sum_{i=1}^A\sum_{j=1}^B\sum_{k=1}^C\sum_{di\;dj^2\;dk^3}(\phi*\mu)(d)$$
$$\sum_{d=1}^{A}(\phi*\mu)(d)\sum_{i=1\;di}^A\sum_{j=1\;dj^2}^B\sum_{k=1\;dk^3}^C1$$ \left \lceil \right \rceil $观察后面式子,对于一个x^k,若dx^k,先把d分解为\prod p_i^{a_i}.$ $则\prod p_i^{a_i}x^k,得\prod p_i^{\left \lceil \frac{a_i}{k} \right \rceil }x,所以设$
$$f_k(n)=\prod p_i^{\left \lceil \frac{a_i}{k} \right \rceil }$$
$则:$
$$ans=\sum_{d=1}^{A}(\phi * \mu)(d)\frac{A}{f_1(d)}\frac{B}{f_2(d)}\frac{C}{f_3(d)}$$
$对于f_k(n),类似分解n的形式。$ $分解质因子的过程中,记录质因子的指数,每次质因子+1时,若\%k=1时,说明该向上取整了,于是f_k(n)*=该质因子.$
$由于\phi和\mu都是积性函数,卷积之后还是积性函数,设g(d)=(\phi*\mu)(d),则$
$$\phi(n)=\sum_{dn}g(d)$$
$那么可以得:$
$$\phi(p^k)=\phi(p^{k-1})+g(p^k)$$
$$g(p^k)=\phi(p^k)-\phi(p^{k-1})$$
$$g(p^k)=(p-1)p^{k-1}-(p-1)p^{k-2}$$
$$g(p^k)=(p-1)^2p^{k-2}$$
$当我们欧拉筛的过程中:$
$$g(1)=1$$
$$g(p)=p-2$$
$$g(p^{k})=(p-1)^2p^{k-2}$$
$$g(p_1^{k_1}p_2^{k_2})=g(p_1^{k_1})g(p_2^{k_2})$$
$$……$$
Code 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 ##include  "bits/stdc++.h"   using  namespace  std;typedef  long  long  ll;const  ll mod = (1ll  << 30 );const  int  N = 1e7  + 10 ;bool  is\_prime[N];int  prime[N], cnt;int  g[N];int  f1[N], f2[N], f3[N];int  deg[N];void  init ()      f1[1 ] = f2[1 ] = f3[1 ] = g[1 ] = 1 ;     for (int  i = 2 ;i < N; i++) {         f1[i] = i;         if (!is\_prime[i]) prime[++cnt] = f2[i] = f3[i] = i, g[i] = i - 2 , deg[i] = 1 ;         for (int  j = 1 ;j <= cnt && i \* prime[j] < N; j++) {             int  now = i \* prime[j];             is\_prime[now] = 1 ;             if (i % prime[j] == 0 ) {                 deg[now] = deg[i] + 1 ;                 int  num = 1 , tmp = i;                 while (num <= 3  && tmp % prime[j] == 0 ) num++, tmp /= prime[j];                 if (num == 1 ) g[now] = g[i] \* g[prime[j]];                 else  if (num == 2 ) g[now] = g[i / prime[j]] \* (prime[j] - 1 ) \* (prime[j] - 1 );                 else  g[now] = g[i] \* prime[j];                 f2[now] = f2[i] \* (deg[now] % 2  == 1  ? prime[j] : 1 );                 f3[now] = f3[i] \* (deg[now] % 3  == 1  ? prime[j] : 1 );                 break ;             }             else  {                 deg[now] = 1 ;                 f2[now] = f2[i] \* f2[prime[j]];                 f3[now] = f3[i] \* f3[prime[j]];                 g[now] = g[i] \* g[prime[j]];             }         }     } } void  solve ()      init ();     int  \_; scanf ("%d" ,&\_);     while (\_--) {         int  A, B, C; scanf ("%d%d%d" ,&A,&B,&C);         ll ans = 0 ;         for (int  d = 1 ;d <= A; d++) {             ans = (ans + 1ll  \* g[d] \* (A / f1[d]) % mod \* (B / f2[d]) % mod \* (C / f3[d]) % mod) % mod;         }         printf ("%lld\n" ,(ans % mod + mod) % mod);     } } signed  main ()      ios\_base::sync\_with\_stdio(false );           ##ifdef  FZT\_ACM\_LOCAL      freopen ("in.txt" , "r" , stdin);     freopen ("out.txt" , "w" , stdout);     signed  test\_index\_for\_debug = 1 ;     char  acm\_local\_for\_debug = 0 ;     do  {         if  (acm\_local\_for\_debug == '$' ) exit (0 );         if  (test\_index\_for\_debug > 20 )             throw  runtime\_error("Check the stdin!!!" );         auto  start\_clock\_for\_debug = clock ();         solve ();         auto  end\_clock\_for\_debug = clock ();         cout << "Test "  << test\_index\_for\_debug << " successful"  << endl;         cerr << "Test "  << test\_index\_for\_debug++ << " Run Time: "               << double (end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s"  << endl;         cout << "--------------------------------------------------"  << endl;     } while  (cin >> acm\_local\_for\_debug && cin.putback (acm\_local\_for\_debug)); ##else       solve (); ##endif       return  0 ; } 
本文作者:jujimeizuo 本文地址https://blog.jujimeizuo.cn/2020/08/10/mobius-inversion/  本博客所有文章除特别声明外,均采用 CC BY-SA 3.0 协议。转载请注明出处!