[数论][莫比乌斯反演][杜教筛] 51Nod 1238 最小公倍数之和 V3

最新推荐文章于 2020-12-13 21:51:42 发布

原创最新推荐文章于 2020-12-13 21:51:42 发布 · 560 阅读

1 ·

CC 4.0 BY-SA版权

数论同时被 2 个专栏收录

20 篇文章

订阅专栏

莫比乌斯反演&杜教筛

8 篇文章

订阅专栏

该博客详细介绍了如何利用数论中的莫比乌斯反演和杜教筛方法解决51Nod 1238题——计算最小公倍数之和的问题，探讨了数论在算法中的应用。

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

$Description$

求

\sum i = 1 n \sum j = 1 n l c m (i, j)

$\sum_{i=1}^n\sum_{j=1}^nlcm(i,j)$
推柿子：

A (n) A n s = = = = = = = = \sum i = 1 n l c m (i, n) n \sum d ∣ n \sum i = 1 n i d [(n d, i d) = 1] n \sum d ∣ n \sum i = 1 d i [(d, i) = 1] n \sum d ∣ n [ d = 1 ] + d φ ( d ) 2 n 2 + 1 2 \sum d ∣ n d φ (d) 2 \sum i = 1 n A (i) - n ( n + 1 ) 2 \sum i = 1 n i \sum d ∣ i d φ (d) \sum d = 1 n d 2 φ (d) \sum i = 1 ⌊ n d ⌋ i

$\begin{eqnarray} A(n)&=&\sum_{i=1}^nlcm(i,n)\\ &=&n\sum_{d\mid n}\sum_{i=1}^n{i\over d}[({n\over d}, {i\over d})=1]\\ &=&n\sum_{d\mid n}\sum_{i=1}^di[(d, i)=1]\\ &=&n\sum_{d\mid n}{{[d=1]+d\varphi(d)}\over2}\\ &=&{n\over2}+{1\over2}\sum_{d\mid n}d\varphi(d)\\ Ans&=&2\sum_{i=1}^nA(i)-{{n(n+1)}\over2}\\ &=&\sum_{i=1}^ni\sum_{d\mid i}d\varphi(d)\\ &=&\sum_{d=1}^nd^2\varphi(d)\sum_{i=1}^{\lfloor{n\over d}\rfloor}i \end{eqnarray}$
后面那个东西取值只有

O(n√) $O(\sqrt{n})$ 种。
所以杜教筛

S (n) = = = \sum i = 1 n (i d 2 \cdot φ) (i) \sum i = 1 n [(φ * 1) \cdot i d 2] (i) - \sum i = 2 n S (⌊ n i ⌋) i 2 \sum i = 1 n i 3 - \sum i = 2 n S (⌊ n i ⌋) i 2

$\begin{eqnarray} S(n)&=&\sum_{i=1}^n(id^2\cdot\varphi)(i)\\ &=&\sum_{i=1}^n\bigg[(\varphi*1)\cdot id^2\bigg](i)-\sum_{i=2}^nS(\lfloor{n\over i}\rfloor)i^2\\ &=&\sum_{i=1}^ni^3-\sum_{i=2}^nS(\lfloor{n\over i}\rfloor)i^2 \end{eqnarray}$
还有一种处理方法，从 jiry_2博客上看到的。
这里写图片描述

#include <bits/stdc++.h>
using namespace std;

typedef long long ll;
const ll N = 10101010;
const ll MOD = 1000000007;
const ll INV2 = (MOD + 1) / 2;
const ll INV6 = (MOD + 1) / 6;

inline char get(void) {
  static char buf[100000], *S = buf, *T = buf;
  if (S == T) {
    T = (S = buf) + fread(buf, 1, 100000, stdin);
    if (S == T) return EOF;
  }
  return *S++;
}
template<typename T>
inline void read(T &x) {
  static char c; x = 0; ll sgn = 0;
  for (c = get(); c < '0' || c > '9'; c = get()) if (c == '-') sgn = 1;
  for (; c >= '0' && c <= '9'; c = get()) x = x * 10 + c - '0';
  if (sgn) x = -x;
}

int phi[N], pre[N];
int vis[N], prime[N];
ll Pcnt, lim, ans;
ll n;
map<ll, int> mp;

inline void Add(ll &x, ll a) {
  x = (x + a) % MOD;
}
inline ll Mod(ll x) {
  return (x % MOD + MOD) % MOD;
}
inline void Pre(int n) {
  phi[1] = 1; int x;
  for (int i = 2; i <= n; i++) {
    if (!vis[i]) {
      phi[i] = i - 1; prime[++Pcnt] = i;
    }
    for (int j = 1; j <= Pcnt && (x = prime[j] * i) <= n; j++) {
      vis[x] = 1;
      if (i % prime[j]) {
    phi[x] = phi[prime[j]] * phi[i];
      } else {
    phi[x] = prime[j] * phi[i]; break;
      }
    }
  }
  for (int i = 1; i <= n; i++)
    pre[i] = ((ll)phi[i] * i % MOD * i % MOD + pre[i - 1]) % MOD;
}
inline ll Sum1(ll n) {
  n %= MOD;
  return n * (n + 1) / 2 % MOD;
}
inline ll Sum2(ll n) {
  n %= MOD; ll n2 = (2 * n + 1) % MOD;
  return n * (n + 1) % MOD * n2 % MOD * INV6 % MOD;
}
inline ll Sum2(ll l, ll r) {
  return Mod(Sum2(r) - Sum2(l - 1));
}
inline ll Sum3(ll n) {
  n = Sum1(n);
  return n * n % MOD;
}
inline ll S(ll n) {
  if (n <= lim) return pre[n];
  if (mp.count(n)) return mp[n];
  ll res = Sum3(n); ll pos;
  for (ll i = 2; i <= n; i = pos + 1) {
    pos = n / (n / i);
    Add(res, MOD - (ll)Sum2(i, pos) * S(n / i) % MOD);
  }
  return mp[n] = res;
}

int main(void) {
  freopen("1.in", "r", stdin);
  freopen("1.out", "w", stdout);
  read(n);
  Pre(lim = (int)pow(n, 0.7));
  ll pos;
  for (ll d = 1; d <= n; d = pos + 1) {
    pos = n / (n / d);
    Add(ans, (ll)Mod(S(pos) - S(d - 1)) * Sum1(n / d) % MOD);
  }
  cout << ans << endl;
  return 0;
}