51nod 1237 最大公约数之和 V3(杜教筛)_const int inv6 = 166666668;是什么-优快云博客

本文链接：https://blog.youkuaiyun.com/qq_41791981/article/details/98478166

本文深入解析了杜教筛算法，一种高效处理大规模数据的数学筛选技术。通过对比传统线性筛的局限性，文章详细介绍了杜教筛的原理、构造方法及求解过程，并提供了AC代码实例，展示了其在解决特定数学问题上的优越性。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

题意：计算 $\sum_{i=1}^{n}\sum_{j=1}^{n}gcd(i,j)$

由于数据范围达到1e10，显然直接线性筛是完成不了的，所以我们选择杜教筛。

$\sum_{i=1}^{n}\sum_{j=1}^{n}gcd(i,j)=\sum_{k=1}^{n}k\sum_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum_{j=1}^{\lfloor \frac{n}{k} \rfloor}\sum_{d|i.d|j} \mu(d) =\sum_{k=1}^{n}k\sum_{d=1}^{n} \mu(d) \lfloor\frac{n}{kd} \rfloor^2=\sum_{d=1}^{n} \lfloor \frac{n}{d} \rfloor^2 \phi(d)$

对该式分块，令 $f(n)=\sum_{i=1}^{n} \phi(i)$

构造 $g(n)=\sum_{i=1}^{n} i$ ,带入 $f(n)=g(n)-\sum_{i=2}^{n}f(\lfloor \frac{n}{i} \rfloor)$ ，求解。

AC代码：

//#pragma comment(linker, “/STACK:1024000000,1024000000”
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 4e6 + 7;
const int inf = 0x3f3f3f3f;
const int mod = 1e9 + 7;
const int inv2 = 5e8 + 4;
//const double pi = acos(-1.0);
const double eps = 1e-6;
const int inv = 250000002;
const int inv6 = 166666668;
int gcd(int a,int b){return b ? gcd(b,a % b) : a;}
int prime[N],cnt;
ll g[N],phi[N];
bool isprime[N];
ll n;
ll fpow(ll a,ll b)
{
    ll res = 1;
    while(b){
        if(b & 1) res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}
inline ll gao(ll x)
{
    x %= mod;
    return (x + 1) % mod * x % mod * inv2 % mod;
}
void get_phi()
{
   phi[1] = 1;
   for(int i = 2;i < N;++i){
       if(!isprime[i]){
             prime[++cnt] = i;
             phi[i] = i - 1;
        }
       for(int j = 1;j <= cnt && i * prime[j] < N;++j)
       {
          isprime[i * prime[j]] = 1;
          if(i % prime[j] == 0){
             phi[i * prime[j]] = phi[i] * prime[j];
             break;
          }
          phi[i * prime[j]] = phi[i] * (prime[j] - 1);
       }
   }
   for(ll i = 1;i < N;++i) g[i] = (g[i - 1] + phi[i]  % mod) % mod;
}

map<ll ,ll>mp;
ll solve(ll pos)
{
    if(pos < N) return g[pos];
    if(mp[pos]) return mp[pos];
    ll res = gao(pos),last;
    for(ll i = 2;i <= pos;i = last + 1){
        last = pos / (pos / i);
        res = ((res - (last - i + 1) * solve(pos / i) % mod) % mod + mod) % mod;
    }
    mp[pos] = res;
    return res;
}
ll work(ll pos)
{
    ll last,ans = 0;
    for(ll i = 1;i <= pos;i = last + 1){
        last = pos / (pos / i);
        ans = (ans + (pos / i) % mod * ((pos / i) % mod) % mod * (solve(last) - solve(i - 1)) % mod) % mod + mod;
        ans %= mod;
    }
    return ans;
}
int main()
{
    get_phi();
    scanf("%lld",&n);
    printf("%lld\n",work(n));
    return 0;
}