danihao123's Blog

沿用约数个数和一题的思路……可以列出此式：

\[\sigma_1(ij) = \sum_{a | i}\sum_{b | j} ab[\gcd(\frac{i}{a}, b) = 1]\]

然后我们考虑去化简原式：

\[
\begin{aligned}
\sum_{i = 1}^n\sum_{j = 1}^n\sum_{a | i}\sum_{b | j} ab\sum_{d | \frac{i}{a}, d | j}\mu(d)
\end{aligned}
\]

接下来我们考虑枚举\(d\)，但是\(a\)和\(b\)的贡献又有点麻烦了……

\(b\)的贡献还好说，直接枚举\(b\)本身是\(d\)的多少倍就是\(\sum_{b = 1}^{\lfloor\frac{n}{d}\rfloor} bd\lfloor\frac{n}{bd}\rfloor\)。那\(a\)的贡献如何考虑？

我们考虑直接去枚举\(a\)本身……可以注意到一定有\(a\le\lfloor\frac{n}{d}\rfloor\)（因为\(\frac{i}{a}\ge d\)），所以直接枚举\(a\)本身之后再考虑\(\lfloor\frac{n}{a}\rfloor\)范围内\(d\)的倍数的数量即可（相当于找\(\frac{i}{a}\)），因此贡献为\(\sum_{a = 1}^{\lfloor\frac{n}{d}\rfloor} a\lfloor\frac{n}{ad}\rfloor\)。

那么继续化简原式：

\[
\begin{aligned}
\quad&\sum_{d = 1}^n\mu(d)\sum_{b = 1}^{\lfloor\frac{n}{d}\rfloor} bd\lfloor\frac{n}{bd}\rfloor\sum_{a = 1}^{\lfloor\frac{n}{d}\rfloor} a\lfloor\frac{n}{ad}\rfloor\\
=&\sum_{d = 1}^n\mu(d)d(\sum_{i = 1}^{\lfloor\frac{n}{d}\rfloor} i\lfloor\frac{n}{id}\rfloor)^2\\
=&\sum_{d = 1}^n\mu(d)d S_1^2(\lfloor\frac{n}{d}\rfloor)
\end{aligned}
\]

此处\(S_1\)表示\(\sigma_1\)的前缀和。

接下来预处理\(\mu(d)d\)的前缀和，考虑杜教筛。我们发现该函数和\(\mathrm{id}\)卷出来就是\(\epsilon\)（证明可以考虑贝尔级数……这种方式非常有效，并且还能帮我们找到需要卷的函数，我有时间会专门撰文写一下），所以杜教筛一波即可。这部分总复杂度为\(O(n^{\frac{2}{3}})\)。

至于\(S_1\)，我们考虑不大于\(n^{\frac{2}{3}}\)可以直接预处理，剩下的用时就用\(O(\sqrt{n})\)的方法求（考虑到反演不会用到\(S_1\)的重复状态所以不需要记忆化）。用类似于杜教筛复杂度证明的方法可以证明该部分总复杂度为\(O(n^{\frac{2}{3}})\)。

代码：

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cmath>
#include <algorithm>
#include <functional>
#include <utility>
#include <unordered_map>
const int N = 1000000;
using ll = long long;
const ll ha = 1000000007LL;
const ll i2 = 500000004LL;
ll mud[N + 5], sig_1[N + 5];
int prm[N + 5]; bool vis[N + 5];
void sieve() {
  int cnt = 0;
  mud[1] = 1;
  for(int i = 2; i <= N; i ++) {
    if(!vis[i]) {
      mud[i] = (ha - i);
      prm[cnt ++] = i;
    }
    for(int j = 0; j < cnt; j ++) {
      int v = i * prm[j];
      if(v > N) break;
      vis[v] = true;
      if(i % prm[j] == 0) {
        mud[v] = 0; break;
      } else {
        mud[v] = (mud[i] * mud[prm[j]]) % ha;
      }
    }
  }
  for(int i = 1; i <= N; i ++) {
    for(int j = i; j <= N; j += i) {
      sig_1[j] = (sig_1[j] + (ll(i))) % ha;
    }
  }
  for(int i = 1; i <= N; i ++) {
    mud[i] = (mud[i - 1] + mud[i]) % ha;
    sig_1[i] = (sig_1[i - 1] + sig_1[i]) % ha;
  }
}

std::unordered_map<ll, ll> ma;
ll calc_id(ll x) {
  ll v1 = x, v2 = x + 1LL;
  if(x & 1LL) {
    v2 /= 2LL;
  } else {
    v1 /= 2LL;
  }
  return ((v1 * v2) % ha);
}
ll calc_mud(ll n) {
  if(n <= (ll(N))) return mud[n];
  if(ma.count(n)) return ma[n];
  ll ans = 1, las = 1;
  for(ll i = 2; i <= n;) {
    ll next = n / (n / i);
    ll nv = calc_id(next);
    ll ns = (nv - las + ha) % ha;
    ans = (ans - (ns * calc_mud(n / i)) % ha + ha) % ha;
    las = nv; i = next + 1LL;
  }
  ma[n] = ans; return ans;
}
ll calc_s1(ll n) {
  if(n <= (ll(N))) return sig_1[n];
  ll ans = 0, las = 0;
  for(ll i = 1; i <= n;) {
    ll next = n / (n / i);
    ll nv = calc_id(next);
    ll ns = (nv - las + ha) % ha;
    ans = (ans + (ns * (n / i)) % ha) % ha;
    las = nv; i = next + 1LL;
  }
  return ans;
}
ll calc(ll n) {
  ll ans = 0, las = 0;
  for(ll i = 1; i <= n;) {
    ll next = n / (n / i);
    ll nv = calc_mud(next); ll ns = (nv - las + ha) % ha;
    ll pv = calc_s1(n / i); pv = (pv * pv) % ha;
    ans = (ans + (ns * pv) % ha) % ha;
    las = nv; i = next + 1LL;
  }
  return ans;
}

int main() {
  sieve();
  ll n; scanf("%lld", &n);
  printf("%lld\n", calc(n));
  return 0;
}

首先这个题和能量采集几乎完全一致……稍微有点反演常识的人都知道答案事\(\sum_{i = 1}^n\varphi(i)\lfloor\frac{n}{i}\rfloor^2\)，然后你需要用\(\varphi\)的前缀和，这个东西直接杜教筛搞一波就行了。这样总复杂度事\(O(n^{\frac{2}{3}})\)的，下面给出证明。

先考虑杜教筛的部分。杜教筛可以看做一个转移复杂度为\(O(\sqrt{n})\)的DP，他的状态一定是通过总的\(n\)整除一个数得到的。不大于\(n^{\frac{2}{3}}\)的状态我们都事先预处理，那么我们假设所有大于\(n^{\frac{2}{3}}\)的状态都被计算了一遍，这些状态一定事通过\(n\)整除一个不大于\(n^{\frac{1}{3}}\)的数得到的。因此复杂度为：

\[\sum_{i = 1}^{n^{\frac{1}{3}}} \sqrt{\lfloor\frac{n}{i}\rfloor}\leq \int_0^{n^{\frac{1}{3}}} \sqrt{\frac{n}{x}}\mathrm{d}x = 2n^{\frac{2}{3}} = O(n^{\frac{2}{3}})\]

至于莫比乌斯反演，不考虑杜教筛的复杂度的话就只是\(O(\sqrt{n})\)而已了。

代码：

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cctype>
#include <cmath>
#include <algorithm>
#include <functional>
#include <unordered_map>
#include <utility>
using ll = long long;
constexpr ll ha = 1000000007LL;
ll pow_mod(ll a, ll b) {
  ll ans = 1, res = a;
  while(b) {
    if(1LL & b) ans = (ans * res) % ha;
    res = (res * res) % ha; b >>= 1;
  }
  return ans;
}
ll inv(ll x) {
  return pow_mod(x, ha - 2LL);
}

const int N = 10000000;
ll phi[N + 1]; int prm[N + 1]; bool vis[N + 1];
void sieve() {
  phi[1] = 1; vis[1] = true;
  int cnt = 0;
  for(int i = 2; i <= N; i ++) {
    if(!vis[i]) {
      phi[i] = i - 1; prm[cnt ++] = i;
    }
    for(int j = 0; j < cnt; j ++) {
      int v = i * prm[j]; if(v > N) break;
      vis[v] = true;
      if(i % prm[j] == 0) {
        phi[v] = phi[i] * prm[j]; break;
      } else {
        phi[v] = phi[i] * phi[prm[j]];
      }
    }
  }
  for(int i = 1; i <= N; i ++) {
    phi[i] = (phi[i] + phi[i - 1]) % ha;
  }
}

std::unordered_map<ll, ll> ma;
const ll i2 = inv(2LL);
ll phi_sum(ll n) {
  if(n <= (ll(N))) return phi[n];
  if(ma.count(n)) return ma[n];
  ll fg = (n % ha) * ((n + 1LL) % ha) % ha * i2 % ha;
  for(ll i = 2; i <= n;) {
    ll next = n / (n / i);
    ll len = (next - i + 1LL) % ha;
    fg = (fg - (len * phi_sum(n / i)) % ha + ha) % ha;
    i = next + 1LL;
  }
  ma[n] = fg; return fg;
}
ll calc(ll n) {
  ll ans = 0, las = 0;
  for(ll i = 1; i <= n;) {
    ll next = n / (n / i);
    ll b = ((n / i) % ha) * ((n / i) % ha) % ha;
    ll ns = phi_sum(next);
    ans = (ans + (((ns - las + ha) % ha) * b % ha)) % ha; las = ns;
    i = next + 1LL;
  }
  return ans;
}

int main() {
  sieve();
  ll n; scanf("%lld", &n);
  printf("%lld\n", calc(n));
  return 0;
}

MIN-MAX容斥第一题……

众所周知斐波那契数有个非常好的性质：

\[\gcd(\{f_i | i\in T\}) = f_{\gcd\{T\}}\]

但很不幸的事情是，lcm没有类似的性质。那我们就考虑怎么把lcm改成gcd。

根据MIN-MAX容斥，可以得到这个式子（可以理解为对指数的容斥。这里设\(S\)为全集）：

\[\prod_{T\subseteq S, T\neq\emptyset} f_{\gcd\{T\}}^{(-1)^{|T| + 1}}\]

gcd套在上面非常烦人，考虑用狄利克雷卷积去化解它。考虑有个函数\(g\)满足\(f_n = \prod_{d|n} g(d)\)，然后xjb化简一下之后柿子变成：

\[\prod_{d} g(d)^{\sum_{T\subseteq S, T\neq\emptyset, \forall x\in T, d | x}(-1)^{|T| + 1}}\]

再去考虑那个指数。当且仅当\(d\in S\)时，\(g(d)\)才会对答案做贡献。并且无论\(S\)中\(d\)的约数有多少个（只要大于），最后指数肯定都是1。因为其实大于零的情况那个指数也是一个所有最小值都是1的MIN-MAX容斥……这样求出来最大值肯定是1了，那个指数自然就是1。

于是乎柿子变成了：

\[\prod_{\exists x\in S, d | x} g(d)\]

然后搞出来\(S\)中所有数的约数很容易。当务之急就是求\(g\)了（这里其实就是莫比乌斯反演了）。不过球莫比乌斯函数反而让问题复杂化了……移下项可以发现：

\[g(n) = f_n\cdot (\prod_{d | n, d\neq n} g(d))^{-1}\]

然后就很愉快的搞出来了（逃

代码：

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cmath>
#include <cctype>
#include <algorithm>
#include <utility>
#include <queue>
const int maxn = 1000005;
using ll = long long;
const ll ha = 1000000007LL;
ll pow_mod(ll a, ll b) {
  ll ans = 1LL, res = a;
  while(b) {
    if(1LL & b) {
      ans = (ans * res) % ha;
    }
    res = (res * res) % ha;
    b /= 2LL;
  }
  return ans;
}
ll inv(ll x) {
  return pow_mod(x, ha - 2LL);
}

ll f[maxn], g[maxn];
void process_f() {
  int N = 1000000;
  f[0] = 0LL; f[1] = 1LL;
  for(int i = 2; i <= N; i ++) {
    f[i] = (f[i - 1] + f[i - 2]) % ha;
  }
  std::fill(g + 1, g + 1 + N, 1LL);
  for(int i = 1; i <= N; i ++) {
    g[i] = (f[i] * inv(g[i])) % ha;
    for(int j = i * 2; j <= N; j += i) {
      g[j] = (g[j] * g[i]) % ha;
    }
  }
#ifdef LOCAL
  puts("g has been processed!");
#endif
}

bool vis[maxn];
void process_vis() {
  int N = 1000000;
  for(int i = 1; i <= N; i ++) {
    for(int j = i * 2; j <= N; j += i) {
      vis[i] = (vis[i] || vis[j]);
    }
  }
}

int main() {
  process_f();
  int n; scanf("%d", &n);
  for(int i = 1; i <= n; i ++) {
    int a; scanf("%d", &a);
    vis[a] = true;
  }
  process_vis();
  ll ans = 1LL;
  for(int i = 1; i <= 1000000; i ++) {
    if(vis[i]) {
      ans = (ans * g[i]) % ha;
    }
  }
  printf("%lld\n", ans);
  return 0;
}