「笔记」莫比乌斯反演

目录

• 前置知识
  • 小碎骨
  • 引理1
  • 数论分块
• 积性函数
  • 定义
  • 性质
  • 常见积性函数
• 莫比乌斯函数
  • 定义
  • 性质
  • 反演常用结论
  • 线性筛求莫比乌斯函数
• 狄利克雷(Dirichlet)卷积
  • 性质
  • 关于单位元 $e$
  • 除数函数与幂函数
  • 欧拉函数与恒等函数
  • 计算
• 莫比乌斯反演
  • 公式
  • 证明
  • 举例
• 关于一类求和式
  • 一般套路
  • 例 1
  • 例 2
  • 感悟
• 例题
  • [HAOI2011] Problem b
  • [国家集训队]Crash的数字表格
  • SP5971 LCMSUM - LCM Sum
  • [SDOI2015]约数个数和
  • P3768 简单的数学题
  • 「SDOI2017」数字表格
  • P5518 [MtOI2019]幽灵乐团 / 莫比乌斯反演基础练习题
• 写在最后

Updated on 2020.8.6
巨幅更新，对积性函数和狄利克雷卷积部分进行重构。
新增对一类特殊求和式的讲解。

Updated on 2020.9.9
添加了几个例题。

---

前置知识

小碎骨

艾佛森括号 $[P] = \begin{cases} 1 &\text{If P is true}\\ 0 &\text{Otherwise} \end{cases}$
此处 $P$ 是一个可真可假的命题。

引理1

$$\forall a,b,c\in \mathbb{Z},\left\lfloor\dfrac{a}{bc}\right\rfloor = \left\lfloor{\dfrac{\left\lfloor\dfrac{a}{b}\right\rfloor}{c}}\right\rfloor$$

证明

$$\dfrac{a}{b} = \left\lfloor{\dfrac{a}{b}}\right\rfloor + r(0\le r < 1)$$

$$\left\lfloor\dfrac{a}{bc}\right\rfloor = \left\lfloor\dfrac{a}{b}\times\dfrac{1}{c}\right\rfloor = \left\lfloor{\dfrac{1}{c}\times\left({\left\lfloor{\dfrac{a}{b} + r}\right\rfloor}\right)}\right\rfloor = \left\lfloor{\dfrac{\left\lfloor\dfrac{a}{b}\right\rfloor}{c} + \dfrac{r}{c}}\right\rfloor = \left\lfloor{\dfrac{\left\lfloor{\dfrac{a}{b}}\right\rfloor}{c}}\right\rfloor$$

数论分块

内容独立了出来，详细内容见 数论分块 - Luckyblock

对于一类含有$\left\lfloor\frac{n}{i}\right\rfloor$的求和式 ($n$ 为常数)，由于$\left\lfloor\frac{n}{i}\right\rfloor$单调不增，故存在多个区间$[l,r]$, 使得$\left\lfloor\frac{n}{i}\right\rfloor = \left\lfloor\frac{n}{j}\right\rfloor(i,j\in [l,r])$ 成立。

对于任意一个$i$，最大的满足上式的 $j=\left\lfloor{\dfrac{n}{\left\lfloor{\dfrac{n}{i}}\right\rfloor}}\right\rfloor$

---

积性函数

定义

若 $\gcd(x,y) = 1$ 且 $f(xy)=f(x)f(y)$，则 $f(n)$ 为积性函数。

性质

若 $f(x)$，$g(x)$均为积性函数，则以下函数也为积性函数：

$$\begin{aligned} & h(x) = f(x^p)\\ & h(x) = f^p(x)\\ & h(x) = f(x)g(x)\\ & h(x) = \sum_{d\mid x} f(d)g\left(\dfrac{x}{d}\right) \end{aligned}$$

常见积性函数

• 单位函数 $e(n) = [n = 1]$。
• 幂函数 $\operatorname{Id}_{k}(n) = n^k$， $\operatorname{id}_1(n)$ 通常简记为$\operatorname{id}(n)$。
• 常数函数 $1(n) = 1$。
• 因数个数 $\operatorname{d}(n) = \sum\limits_{d\mid n} 1$。
• 除数函数 $\sigma_{k}(n) = \sum\limits_{d\mid n} d^k$。
  $k=1$ 时为因数和函数，通常简记为 $\sigma(n)$。
  $k=0$ 时为因数个数函数 $\sigma_{0}(n)$。
• 欧拉函数 $\varphi(n) = \sum\limits_{i=1}^{n} [\gcd(i,n) = 1]$。
• 莫比乌斯函数 $\mu(n) = \begin{cases}1 &n=1\\0 &n\ \text{含有平方因子}\\(-1)^k &k\text{为}\ n\ \text{的本质不同质因子个数} \end{cases}$

不是很懂上面写的什么玩意？
不用深究，有个印象继续往下看就好。

---

莫比乌斯函数

定义

$\mu$ 为莫比乌斯函数，定义为

$$\mu(n) = \begin{cases} 1 &n=1\\0 &n\ \text{含有平方因子}\\(-1)^k &k\text{为}\ n\ \text{的本质不同质因子个数} \end{cases}$$

解释

令 $n = \prod\limits_{i=1}^{k} p_{i}^{c_i}$，其中$p_i$为质因子，$c_i\ge 1$。

1. $n=1$时，$\mu (n) = 1$。
2. $n\not ={1}$时 ，
   • $\exist i\in [1,k], c_i > 1$ 时，$\mu (n) = 0$。
     当某质因子出现次数大于$1$时，$\mu (n) = 0$。
   • $\forall i\in [1,k], c_i = 1$ 时，$\mu (n) = (-1)^k$。
     当每个质因子只出现一次时，即$n = \prod\limits_{i=1}^{k}p_i$，$\{p_i\}$中元素唯一。
     $\mu (n) = (-1)^k$，此处$k$为质因子的种类数。

性质

莫比乌斯函数是积性函数，且具有以下性质

$$\large \sum_{d\mid n} \mu (d) = [n=1]$$

证明，设 $n = \prod\limits_{i=1}^{k}{p_i^{c_i}}, n' = \prod\limits_{i=1}^{k}{p_i}$。

• 根据莫比乌斯函数定义，则有：$\sum\limits_{d\mid n}{\mu(d)} = \sum\limits_{d\mid n'}{\mu(d)}$。
• 若 $n'$ 的某因子 $d$，有 $\mu (d) = (-1)^i$，则它由 $i$ 个 本质不同的质因子组成。
  由于质因子总数为 $k$，满足上式的因子数为 $C_{k}^{i}$。
• 对于原求和式，转为枚举 $\mu(d)$ 的值。
  则 $\sum\limits_{d\mid n'}{\mu(d)} = \sum\limits_{i=0}^{k}{C_{k}^{i} \times (-1)^i} = \sum\limits_{i=0}^{k}{C_{k}^{i} \times (-1)^i\times 1^{k-i}}$
  根据二项式定理，上式 $= (1+(-1))^k$，
  易知该式在 $k=0$，即 $n=1$ 时为 $1$，否则为 $0$。

反演常用结论

一个反演常用结论:

$$[\gcd(i,j) = 1] = \sum\limits_{d\mid \gcd(i,j)} {\mu (d)}$$

证明 1：
设 $n = \gcd(i,j)$，则右$= \sum\limits_{d\mid n} {\mu (d)} = [n = 1] = [\gcd(i,j) = 1]=$ 左。

证明 2：
暴力展开：$[\gcd(i,j) = 1] = e(\gcd(i,j)) = \sum\limits_{d\mid \gcd(i,j)}\mu(d)$

线性筛求莫比乌斯函数

$\mu$ 为积性函数，因此可以线性筛莫比乌斯函数。

复制复制
int pnum, mu[kMaxn], p[kMaxn];
bool vis[kMaxn];
 
void Euler(int lim_) {
  vis[1] = true, mu[1] = 1ll;
  for (int i = 2; i <= lim_; ++ i) {
    if (! vis[i]) {
      p[++ pnum] = i;
      mu[i] = - 1;
    }
    for (int j = 1; j <= pnum && i * p[j] <= lim_; ++ j) {
      vis[i * p[j]] = true;
      if (i % p[j] == 0) { //勿忘平方因子
        mu[i * p[j]] = 0;
        break;
      }
      mu[i * p[j]] = - mu[i];
    }
  }
}

---

狄利克雷(Dirichlet)卷积

建议阅读 算法学习笔记(35): 狄利克雷卷积 By: Pecco。

定义两个数论函数 $f,g$ 的狄利克雷卷积为

$$\large(f\ast g) (n) = \sum_{d\mid n} f(d)g\left(\dfrac{n}{d}\right)$$

性质

1. 显然满足 交换律，结合律，分配律。
   • $f \ast g = g \ast f$
   • $(f \ast g) \ast h = f\ast (g\ast h)$
   • $f\ast (g+h) = f\ast g + f\ast h$
2. $e$ 为狄利克雷卷积的单位元，有$(f\ast e)(n) = f(n)$。
3. 若 $f,g$ 为积性函数，则 $f\ast g$ 为积性函数。

关于单位元 $e$

有：

$$e = \mu \ast 1=\sum\limits_{d\mid n} \mu (d)$$

证明

$$\begin{aligned} (f\ast e)(n) = \sum_{d\mid n} f(d)e(\dfrac{n}{d}) = \sum_{d\mid n} f(d)\left[\dfrac{n}{d} = 1\right] \end{aligned}$$

• 对于$[\dfrac{n}{d} = 1]$，当且仅当 $\dfrac{n}{d}=1$，即 $d=n$ 时为 $1$，否则为$0$。
• 则当 $d=n$ 时，$f(d)\left[\dfrac{n}{d} = 1\right] = f(n)$。
  当 $d\not ={n}$ 时，$f(d)\left[\dfrac{n}{d} = 1\right] = 0$。

综上，$(f\ast e)(n) = \sum\limits_{d\mid n} f(d)\left[\dfrac{n}{d} = 1\right] = f(n)$，满足单位元的性质。
则 $e = \mu \ast 1$ 成立。

除数函数与幂函数

幂函数 $\operatorname{Id}_{k}(n) = n^k$。
除数函数 $\sigma_{k}(n) = \sum\limits_{d\mid n} d^k$。

显然有：

$$(\operatorname{Id}_k\ast 1)(n) = \sum_{d\mid n} \operatorname{Id_k}(d) = \sum_{d\mid n} d^k = \sigma_k(n)$$

当 $k=0$ 时，$\operatorname{Id_0}=1$，$\sigma_0$ 为因数个数函数，有：

$$(1\ast 1)(n) = \sum_{d\mid n}1 = \sigma_0(n)$$

欧拉函数与恒等函数

$$\begin{aligned} \varphi \ast 1 =& \operatorname{Id}\\ \varphi =& \operatorname{Id}\ast \mu \end{aligned}$$

对于一式，当 $n=p^m$ 时（$p$ 为质数），有：

$$(\varphi \ast 1)(p^m) = \sum_{d\mid n}\varphi(d) = \varphi(1) +\sum_{i=1}^{m}\varphi(p^i) = 1 +\sum_{i=1}^{m}(p^i-p^{i-1}) = p^m$$

$p^i$ 的因子有 $p^{i-1}$ 个，为 $1\sim p^{i-1}$，故 $\varphi(p^i) = p^i-p^{i-1}$。

又 $(\varphi \ast 1)(n)$ 为积性函数，则对于任意正整数 $n$，有：

$$(\varphi \ast 1)(n) = (\varphi \ast 1)\left(\prod p^m\right) = \prod\left(\varphi \ast 1\right)(p^m) = \prod p^m = n$$

得证。

对于 2 式，在 1 式基础上两侧同时 $\ast \mu$ 即得。
左侧变为 $\varphi \ast 1\ast \mu = \varphi \ast e = \varphi$。

计算

$$(f\ast g) (n) = \sum_{d\mid n} f(d)g\left(\dfrac{n}{d}\right)$$

考虑枚举 $n$ 的因子，将贡献累加即可。
显然可以使用埃氏筛筛出所有前缀狄利克雷卷积，复杂度 $O(nk\log n)$，其中 $k$ 是计算函数值的复杂度。

---

莫比乌斯反演

反演是个啥？反演。

公式

设$f(n),g(n)$为两个数论函数。
如果有

$$\large f(n) = (g\ast 1)(n) = \sum\limits_{d\mid n}{g(d)}$$

那么有

$$\large g(n) = (\mu \ast f)(n)=\sum\limits_{d\mid n} {\mu(d)f\left(\dfrac{n}{d}\right)}$$

证明

方法一：运用卷积。

原问题为：已知 $f = g\ast 1$，证明 $g = f\ast \mu$。
易知如下转化：$f\ast \mu = g\ast 1 \ast \mu \Longrightarrow f\ast \mu = g\ast e = g$。

方法二：对原式进行数论变换。

1. 用 $\sum\limits_{d\mid n}g(d)$ 替换$f\left(\dfrac{n}{d}\right)$。

   $$\sum_{d\mid n}{\mu(d)\sum_{k\mid \frac{n}{d}}g(k)}$$

2. 变换求和顺序。

   $$\sum_{k\mid n}g(k)\sum_{d\mid \frac{n}{k}}{\mu(d)}$$

3. 依据 $\sum\limits_{d\mid n}{\mu(d)} = [n=1]$，仅当 $\dfrac{n}{k} = 1$ 时，$\sum\limits_{d\mid \frac{n}{k}}{\mu(d)} = 1$，否则为 $0$。
   此时$k=n$，故上式等价于 $\sum\limits_{k\mid n} {[n=k]\cdot g(k)} = g(n)$。

举例

从 狄利克雷(Dirichlet)卷积 部分可以知道一些积性函数的关系。
尝试对它们进行反演：

$$e = \mu \ast 1\iff\mu = \mu\ast e = \mu$$

$$\sigma_k = \operatorname{Id}_k\ast 1\iff \operatorname{Id}_k = \mu \ast \sigma_k$$

$$\operatorname{Id}=\varphi \ast 1\iff \varphi = \operatorname{Id}\ast \mu$$

---

关于一类求和式

$$\sum_{i=1}^{n}\sum_{j=1}^{m}f(\gcd(i,j))$$

一般套路

考虑构造出函数 $g$，满足下列形式：

$$f(n) = g\ast 1 = \sum_{d\mid n}g(d)$$

则求和式变为：

$$\sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{d\mid \gcd(i,j)}g(d)$$

考虑算术基本定理，发现若满足 $d\mid \gcd (i,j)$，则 $d\mid i$ 且 $d\mid j$ 成立。
考虑 $g(d)$ 在何时做出贡献，调整枚举顺序：

$$\sum_{d=1}g(d)\sum_{i=1}^n[d\mid i]\sum_{j=1}^m[d\mid j]$$

$\sum\limits_{i=1}^{n}[d\mid i]$ 等价于 $1\sim n$ 中 $d$ 的倍数的个数，则上式等价于：

$$\sum_{d=1}g(d)\left\lfloor\dfrac{n}{d}\right\rfloor\left\lfloor\dfrac{m}{d}\right\rfloor$$

数论分块求解即可。

例 1

$$\sum_{i=1}^{n}\sum_{j=1}^{m}\gcd(i,j)$$

发现此题的 $f$ 等价于 $\operatorname{Id}$，则上式等价于：

$$\begin{aligned} &\sum_{i=1}^{n}\sum_{j=1}^{m}\operatorname{Id}(\gcd(i,j))\\ =& \sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{d\mid \gcd(i,j)}\varphi(d)\\ =& \sum_{d=1}\varphi(d)\left\lfloor\dfrac{n}{d}\right\rfloor\left\lfloor\dfrac{m}{d}\right\rfloor \end{aligned}$$

例 2

$$\begin{aligned} &\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i,j) = 1]\\ =& \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}e(\gcd(i,j))\\ =& \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\sum_{d\mid \gcd(i,j)}\mu (d)\\ =& \sum_{d=1}\mu(d)\left\lfloor\dfrac{n}{d}\right\rfloor\left\lfloor\dfrac{m}{d}\right\rfloor \end{aligned}$$

感悟

卷点什么东西，把 $g$ 卷出来。
$g$ 不一定是特殊意义的函数。

---

例题

[HAOI2011] Problem b

$n$ 组询问，每次给定参数 $a,b,c,d,k$，求：

$$\sum\limits_{x=a}^{b}\sum\limits_{y=c}^{d}[\gcd(x,y) = k]$$

$1\le n,k,a,b,c,d\le 5\times 10^4$，$a\le b,c\le d$。

令 $f(n,m) = \sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i,j) = k]$。
根据容斥原理，则原式等价于 $f(b,d) - f(a-1,d) - f(b,d-1) + f(a-1,d-1)$。
$f$ 变成了上述一类求和式的形式，考虑化简 $f$。

易知原式等价于

$$\sum\limits_{i=1}^{\left\lfloor\frac{n}{k}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor\frac{m}{k}\right\rfloor}[\gcd(i,j) = 1]$$

代入反演常用结论 $[\gcd(i,j) = 1] = \sum\limits_{d\mid \gcd(i,j)} {\mu (d)}$，上式化为

$$\sum\limits_{i=1}^{\left\lfloor\frac{n}{k}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor\frac{m}{k}\right\rfloor}\sum_{d\mid \gcd(i,j)}{\mu(d)}$$

变换求和顺序，先枚举$d\mid \gcd(i,j)$，可得

$$\sum_{d=1}\mu(d)\sum_{i=1}^{\left\lfloor\frac{n}{k}\right\rfloor}[d\mid i]\sum_{j=1}^{\left\lfloor\frac{m}{k}\right\rfloor}[d\mid j]$$

对于上式后半的解释：当$d\mid i$ 且 $d\mid j$ 时，$d\mid \gcd(i,j)$。

易知$1\sim \left\lfloor\dfrac{n}{k}\right\rfloor$ 中 $d$ 的倍数有 $\left\lfloor\dfrac{\left\lfloor\dfrac{n}{k}\right\rfloor}{d}\right\rfloor = \left\lfloor\dfrac{n}{kd}\right\rfloor$ 个（由引理 1 可知），原式变为

$$\sum_{d=1}\mu(d)\left\lfloor\dfrac{n}{kd}\right\rfloor\left\lfloor\dfrac{m}{kd}\right\rfloor$$

预处理 $\mu$ 后，显然可以数论分块求解，复杂度$O(n + T\sqrt{n})$。

---

代码

//知识点：莫比乌斯反演
/*
//By:Luckyblock
*/
#include <cstdio>
#include <cctype>
#include <algorithm>
#define ll long long
const int MARX = 6e4 + 10;
//=============================================================
int N, a, b, c, d, k;
int cnt, Mobius[MARX], Prime[MARX], Sum[MARX];
bool vis[MARX];
//=============================================================
inline int read()
{
    int f = 1, w = 0; char ch = getchar();
    for(; !isdigit(ch); ch = getchar()) if(ch == '-') f = -1;
    for(; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
    return f * w;
}
void GetMobius()
{
    Mobius[1] = 1;
    int MAX = MARX - 10;
    for(int i = 2; i <= MAX; i ++)
    {
      if(! vis[i]) Mobius[i] = - 1, Prime[++ cnt] = i;
      for(int j = 1; j <= cnt && i * Prime[j] <= MAX; j ++)
      {
        vis[i * Prime[j]] = true;
        if(i % Prime[j] == 0) break;
        Mobius[i * Prime[j]] = - Mobius[i];
      }
    }
    for(int i = 1; i <= MAX; i ++) Sum[i] = Sum[i - 1] + Mobius[i];
}
ll Calc(int x, int y)
{
    ll ans = 0ll; int max = std :: min(x, y);
    for(int l = 1, r; l <= max; l = r + 1)
      r = std :: min(x / (x / l), y / (y / l)),
      ans += (1ll * x / (1ll * l * k)) * (1ll * y / (1ll * l * k)) * (Sum[r] - Sum[l - 1]);
    return ans;
}
ll Solve()
{
    a = read(), b = read(), c = read(), d = read(), k = read();
    return Calc(b, d) - Calc(b, c - 1) - Calc(a - 1, d) + Calc(a - 1, c - 1);
}
//=============================================================
int main()
{ 
    N = read(); GetMobius();
    while(N --) printf("%lld\n", Solve());
    return 0;
}

---

[国家集训队]Crash的数字表格

给定 $n,m$ , 求：

$$\sum_{i=1}^n\sum_{j=1}^{m} \operatorname{lcm}(i,j)\bmod 20101009$$

$1\le n,m\le 10^7$。

易知原式等价于：

$$\sum_{i=1}^{n}\sum_{j=1}^{m}\dfrac{ij}{\gcd (i,j)}$$

考虑枚举 $\gcd(i,j)$，设枚举量为 $d$。
则 $d=\gcd(i,j)$ 的充要条件是满足 $d|i, d|j$ 且 $\gcd(\dfrac{i}{d},\dfrac{j}{d}) = 1$，则原式等价于：

$$\sum_{i=1}^n\sum_{j=1}^m \sum_{d=1} \dfrac{ij}{d}[d|i][d|j][\gcd(\frac{i}{d}, \frac{j}{d})=1]$$

先枚举 $d$，则原式等价于：

$$\sum_{d=1}\sum_{i=1}^{n}[d\mid i]\sum_{j=1}^m [d\mid j][\gcd(\dfrac{i}{d}, \dfrac{j}{d}=1)] \dfrac{ij}{d}$$

这个 $d$ 很烦人，把 $i,j$ 中的 $d$ 提出来，变为枚举 $\frac{i}{d}$，$\frac{j}{d}$。
消去 $d\mid i$，$d\mid j$ 的限定条件，则原式等价于：

$$\begin{aligned} &\sum_{d=1}\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}[\gcd(i,j)=1]\dfrac{id\times jd}{d}\\ =& \sum_{d=1}\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}[\gcd(i,j)=1]ijd\\ =& \sum_{d=1}d \sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}[\gcd(i,j)=1]ij \end{aligned}$$

单独考虑后半部分，设 $f(x,y) = \sum\limits_{i=1}^{x} \sum\limits_{j=1}^{y}[\gcd(i,j)=1]ij$。
发现 $f(x,y)$ 的形式与 [HAOI2011] Problem b 中的式子类似，代入 $[\gcd(i,j) = 1] = \sum\limits_{d\mid \gcd(i,j)} {\mu (d)}$ 套路一波：

$$\begin{aligned} f(x,y) =& \sum\limits_{i=1}^{x} \sum\limits_{j=1}^{y}[\gcd(i,j)=1]ij\\ =& \sum_{i=1}^{x} \sum_{j=1}^{y}\sum_{d\mid \gcd(i,j)}\mu(d) ij\\ =& \sum_{d=1}\mu(d)\sum_{i=1}^{x}[d\mid i] \sum_{j=1}^{y}[d\mid j] ij\\ =& \sum_{d=1}\mu(d) d^2\sum_{i=1}^{\left\lfloor \frac{x}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{y}{d}\right\rfloor}ij \end{aligned}$$

前半部分 $\sum\limits_{d=1}\mu(d) d^2$，可以考虑筛出 $\mu(d)$ 后求前缀和。
后半部分是等差数列乘等差数列的形式，设 $g(p,q) = \sum\limits_{i=1}^{p} \sum\limits_{j=1}^{q}ij$，$g_{p,q}$ 可以通过下式 $O(1)$ 计算：

$$g(p,q) = \sum_{i=1}^{p} i \sum_{j=1}^{q}j= \dfrac{(1 + p)\times p}{2}\times \dfrac{(1+q)\times q}{2}$$

则对于 $f(x,y)$，有：

$$f(x,y) = \sum_{d=1}\mu(d) d^2\cdot g(\left\lfloor \frac{x}{d}\right\rfloor, \left\lfloor\frac{y}{d}\right\rfloor)$$

数论分块求解即可。

再看回原式，原式等价于：

$$\sum_{d=1}d\cdot f(\left\lfloor\frac{n}{d}\right\rfloor, \left\lfloor\frac{m}{d}\right\rfloor)$$

又是一个可以数论分块求解的形式。
线性筛预处理后 数论分块套数论分块，复杂度 $O(n + m)$，瓶颈是线性筛。

---

一些注意的点

处理时会出现求平方的运算，需要特别注意取模问题，ll 都会爆，被坑惨了。

在预处理前缀和的这个地方：

sum[i] = (sum[i - 1] + 1ll * i * i % kMod * (mu[i] + kMod)) % kMod; //仅令 mu + kMod

注意先对平方取模，在乘上 $\mu$，否则就会爆掉。
以及可以仅令 $mu + mod$。

以及这个地方：

int g(int n_, int m_) {
  return (1ll * n_ * (n_ + 1ll) / 2ll % kMod) * (1ll * m_ * (m_ + 1ll) / 2ll % kMod) % kMod;
}

平方计算，注意随时取模。

代码

//知识点：莫比乌斯反演
/*
By:Luckyblock
*/
#include <algorithm>
#include <cctype>
#include <cmath>
#include <cstdio>
#include <cstring>
#define ll long long
const ll kMod = 20101009;
const int kMaxn = 1e7 + 10;
//=============================================================
int pnum, p[kMaxn];
ll mu[kMaxn], sum[kMaxn];
bool vis[kMaxn];
//=============================================================
inline int read() {
  int f = 1, w = 0;
  char ch = getchar();
  for (; !isdigit(ch); ch = getchar())
    if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
void Getmax(int &fir_, int sec_) {
  if (sec_ > fir_) fir_ = sec_;
}
void Getmin(int &fir_, int sec_) {
  if (sec_ < fir_) fir_ = sec_;
}
void Euler(int lim_) {
  vis[1] = true, mu[1] = 1ll;
  for (int i = 2; i <= lim_; ++ i) {
    if (! vis[i]) {
      p[++ pnum] = i;
      mu[i] = - 1;
    }
    for (int j = 1; j <= pnum && i * p[j] <= lim_; ++ j) {
      vis[i * p[j]] = true;
      if (i % p[j] == 0) { //勿忘平方因子
        mu[i * p[j]] = 0;
        break;
      }
      mu[i * p[j]] = - mu[i];
    }
  }
  sum[1] = 1ll;
  for (int i = 1; i <= lim_; ++ i) {
    sum[i] = (sum[i - 1] + 1ll * i * i % kMod * (mu[i] + kMod)) % kMod; //仅令 mu + kMod
  }
}
int g(int n_, int m_) {
  return (1ll * n_ * (n_ + 1ll) / 2ll % kMod) * (1ll * m_ * (m_ + 1ll) / 2ll % kMod) % kMod;
}
int f(int n_, int m_) {
  int lim = std :: min(n_, m_), ret = 0;
  for (int l = 1, r; l <= lim; l = r + 1) {
    r = std :: min(n_ / (n_ / l), m_ / (m_ / l));
    ret = (ret + 1ll * (sum[r] - sum[l - 1] + kMod) * g(n_ / l, m_ / l) % kMod) % kMod;
  }
  return ret;
}
int Sum(ll l, ll r) {
  return (1ll * (r - l + 1ll) * (l + r) / 2ll) % kMod;
}
//=============================================================
int main() { 
  int n = read(), m = read();
  int lim = std :: min(n, m), ans = 0;
  Euler(lim);
  for (int l = 1, r; l <= lim; l = r + 1) {
    r = std :: min(n / (n / l), m / (m / l));
    ans = (ans + 1ll * Sum(l, r) * f(n / l, m / l) % kMod) % kMod;
  }
  printf("%d", ans);
  return 0;
}
/*
7718820 8445880
*/

---

SP5971 LCMSUM - LCM Sum

$T$ 次询问，每次询问给定 $n$，求：

$$\sum_{i=1}^{n}\operatorname{lcm}(i,n)$$

$1<T\le 3\times 10^5$，$1\le n\le 10^6$。

法一：无脑暴力

先拆 $\operatorname{lcm}$，原式等价于：

$$n\sum_{i=1}^{n}\dfrac{i}{\gcd(i,n)}$$

套路的枚举 $\gcd(i,n)$，调换枚举顺序，原式等价于：

$$\begin{aligned} &n\sum_{i=1}\sum_{d|i} [\gcd(i,n) = d]\dfrac{i}{d}\\ =& n\sum_{i=1}\sum_{d|i} [\gcd(\dfrac{i}{d},\dfrac{n}{d}) = 1]\dfrac{i}{d}\\ =& n\sum_{d|n}\sum_{i=1}^{n}[d|i][\gcd(\dfrac{i}{d},\dfrac{n}{d}) = 1]\dfrac{i}{d} \end{aligned}$$

把 $i,n$ 中的 $d$ 提出来，变为枚举 $\frac{i}{d}$，消去整除的条件，原式等价于：

$$n\sum_{d|n}\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[\gcd(i,\dfrac{n}{d}) = 1]i$$

代入 $[\gcd(i,j) = 1] = \sum\limits_{d\mid \gcd(i,j)} {\mu (d)}$，原式等价于：

$$n\sum_{d|n}\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor} i \sum_{t|\gcd(i,\frac{n}{d})}\mu (t)$$

值得注意的是 $t$ 的上界为 $\frac{n}{d}$，$dt\le n$。
调换枚举顺序，先枚举 $t$，原式等价于：

$$n\sum_{d|n}\sum_{t|\frac{n}{d}} \mu(t) \sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor} [t|i] i$$

套路地消去整除的条件，把 $i$ 中的 $t$ 提出来，原式等价于：

$$n\sum_{d|n}\sum_{t|\frac{n}{d}} \mu(t)t \sum_{i=1}^{\left\lfloor\frac{n}{dt}\right\rfloor} i$$

对于最后的一个求和项，设 $g(x) = \sum\limits_{i=1}^{x}i = \frac{x(x+1)}{2}$，显然可以 $O(1)$ 求解，原式等价于：

$$n\sum_{d|n}\sum_{t|\frac{n}{d}} \mu(t)t\cdot g(\left\lfloor\dfrac{n}{dt}\right\rfloor)$$

考虑枚举 $T = dt$，显然 $T\le n$。
$\mu(t)t$ 与 $d$ 无关，可以直接考虑枚举 $t|T$，原式等价于：

$$n\sum_{T=1}^{n}g(\left\lfloor\dfrac{n}{T}\right\rfloor)\sum_{t|T}\mu(t)t$$

前半块是一个数论分块的形式，可以 $O(\sqrt{n})$ 求解。
考虑后半块，设 $f(n)=\sum\limits_{d|n}\mu(d)d$，发现它是一个积性函数，可以线性筛筛出，具体地：

$$f(n)= \begin{cases} 1-n &n\in \mathrm{primes} \\ f(\frac{x}{p}) &p^2\mid n\\ f(\frac{x}{p})f(p) &p^2\nmid n \end{cases}$$

其中 $p$ 为 $n$ 的最小质因子。

此时已经可以线性筛 + 数论分块求解，复杂度 $O(n+T\sqrt{n})$，比较菜鸡，时限 500ms 过不了。
考虑筛出 $f$ 后再用埃氏筛预处理 $\sum\limits_{T=1}^{n}g(\left\lfloor\dfrac{n}{T}\right\rfloor)f(T)$，输出时乘上 $n$，复杂度变为 $O(n\log^2 n + n)$。

---

法二：

同样先拆 $\operatorname{lcm}$，枚举 $\gcd(i,n)$，调换枚举顺序，原式等价于：

$$\begin{aligned} &n\sum_{i=1}^{n}\dfrac{i}{\gcd(i,n)}\\ =& n\sum_{i=1}\sum_{d|i} [\gcd(i,n) = d]\dfrac{i}{d}\\ =& n\sum_{d|n}\sum_{i=1}^{n}[d|i][\gcd({i},{n}) = d]\dfrac{i}{d} \end{aligned}$$

把 $i,n$ 中的 $d$ 提出来，变为枚举 $\frac{i}{d}$，消去整除的条件，原式等价于：

$$\begin{aligned} &n\sum_{d|n}\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[\gcd(id,n) = d]i\\ =& n\sum_{d|n}\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[\gcd(i,\dfrac{n}{d}) = 1]i \end{aligned}$$

调整枚举对象，上式等价于：

$$n\sum_{d|n}\sum_{i=1}^{d}[\gcd(i,d) = 1]i$$

考虑 $\sum\limits_{i=1}^{d}[\gcd(i,d) = 1]i$ 的实际意义，表示 $[1,d]$ 中与 $d$ 互质的数的和。
当 $d>1$ 时，与 $d$ 互质的数总是成对存在，即若 $\gcd(i,d)=1$ 成立，则 $\gcd(d-i,d)=1$ 成立。
每对这样的数的和为 $d$，共有 $\frac{\varphi(d)}{2}$ 对这样的数。
则原式等价于：

$$n\sum_{d|n}\dfrac{\varphi(d)d}{2}$$

可以直接预处理答案。
预处理时先线性筛出 $\varphi$，再埃氏筛枚举 $i$ 的倍数，令它们的答案加上 $\frac{\varphi(i)i}{2}$，最后输出时乘上 $n$。
复杂度 $O(n\log^2 n + T)$。

---

法二代码

//知识点：莫比乌斯反演
/*
By:Luckyblock
*/
#include <algorithm>
#include <cctype>
#include <cstdio>
#include <cstring>
#define ll long long
const int kMaxn = 1e6;
//=============================================================
ll phi[kMaxn + 10], ans[kMaxn + 10];
int pnum, p[kMaxn + 10];
bool flag[kMaxn + 10];
//=============================================================
inline int read() {
  int f = 1, w = 0;
  char ch = getchar();
  for (; !isdigit(ch); ch = getchar())
    if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
void GetMax(int &fir, int sec) {
  if (sec > fir) fir = sec;
}
void GetMin(int &fir, int sec) {
  if (sec < fir) fir = sec;
}
void GetPrime() {
  phi[1] = 1, flag[1] = true; //注意初始化
  for (int i = 2; i <= kMaxn; ++ i) {
    if (! flag[i]) {
      p[++ pnum] = i;
      phi[i] = i - 1ll;
    }
    for (int j = 1; j <= pnum && i * p[j] <= kMaxn; ++ j) {
      flag[i * p[j]] = true;
      if (i % p[j]) {
        phi[i * p[j]] = phi[i] * phi[p[j]];
      } else {
        phi[i * p[j]] = phi[i] * p[j];
        break;
      }
    }
  }
  for (int i = 1; i <= kMaxn; ++ i) {
    for (int j = 1; i * j <= kMaxn; ++ j) {
      ans[i * j] += (i == 1 ? 1 : 1ll * phi[i] * i / 2);
    }
  }
}
//=============================================================
int main() { 
  GetPrime();
  int T = read();
  while (T --) {
    int n = read();
    printf("%lld\n", 1ll * ans[n] * n);
  }
  return 0; 
}

---

[SDOI2015]约数个数和

$T$ 次询问，每次询问给定 $n,m$。
定义 >$\operatorname{d}(i)$ 为 $i$ 的约数个数，求：

$$\sum_{i=1}^{n}\sum_{j=1}^m\operatorname{d}(ij)$$

$1<T,n\le 5\times 10^4$。

一个结论：

$$\operatorname{d}(ij) = \sum_{x|i}\sum_{y|j}[\gcd(x,y)=1]$$

证明

先考虑 $i = p^a$，$j=p^b(p\in \mathrm{primes})$ 的情况，有：

$$\operatorname{d}(p^{a+b})=\sum_{x|p^a}\sum_{y|p^b}[\gcd(x,y)=1]$$

对于等式左侧，$p^{a+b}$ 的约数个数为 $a+b+1$。
对于等式右侧，在保证 $\gcd(x,y)=1$ 成立的情况下，有贡献的数对 $(x,y)$ 只能是下列三种形式：

• $x>0,y-0$，$x$ 有 $a$ 种取值方法。
• $x=0,y>0$，$y$ 有 $b$ 种取值方法。
• $x=0,y=0$。

则等式右侧贡献次数为 $a+b+1$ 次，等于 $p^{a+b}$ 的约数个数。
则当 $i = p^a$，$j=p^b(p\in \mathrm{primes})$ 时等式成立。

又不同质因子间相互独立，上述情况可拓展到一般的情况。

---

对 $\operatorname{d}(i,j)$ 进行化简，代入 $[\gcd(i,j) = 1] = \sum\limits_{d\mid \gcd(i,j)} {\mu (d)}$，原式等价于：

$$\begin{aligned} \operatorname{d}(ij) =& \sum_{x|i}\sum_{y|j}[\gcd(x,y)=1]\\ =& \sum_{x|i}\sum_{y|j}\sum\limits_{d\mid \gcd(x,y)} {\mu (d)} \end{aligned}$$

调换枚举顺序，先枚举 $d$，原式等价于：

$$\sum_{d=1}[d|i][d|j]{\mu (d)}\sum_{x|i}[d|x]\sum_{y|j}[d|y]$$

把各项中的 $d$ 提出来，消去整除的条件，原式等价于：

$$\begin{aligned} &\sum_{d=1}[d|i][d|j]{\mu (d)}\sum_{x|\frac{i}{d}}\sum_{y|\frac{j}{d}}1\\ =& \sum_{d=1}[d|i][d|j]{\mu (d)}\cdot \operatorname{d}(\dfrac{i}{d})\operatorname{d}(\dfrac{j}{d}) \end{aligned}$$

---

将 $\operatorname{d}(ij)$ 代回原式，原式等价于：

$$\sum_{i=1}^{n}\sum_{j=1}^m \sum_{d=1}[d|i][d|j]{\mu (d)}\cdot \operatorname{d}(\dfrac{i}{d})\operatorname{d}(\dfrac{j}{d})$$

调换枚举顺序，先枚举 $d$，原式等价于：

$$\sum_{d=1}{\mu (d)}\sum_{i=1}^{n}[d|i]\sum_{j=1}^m [d|j]\cdot \operatorname{d}(\dfrac{i}{d})\operatorname{d}(\dfrac{j}{d})$$

把 $i,j$ 中的 $d$ 提出来，变为枚举 $\frac{i}{d}, \frac{j}{d}$，消去整除的条件，原式等价于：

$$\sum_{d=1}{\mu (d)}\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\operatorname{d}(i)\sum_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}\operatorname{d}(j)$$

考虑预处理 $S(x) = \sum\limits_{i=1}^{x}\operatorname{d}(i)$，则原式等价于：

$$\sum_{d=1}{\mu (d)}S\left({\left\lfloor\frac{n}{d}\right\rfloor}\right)S\left({\left\lfloor\frac{m}{d}\right\rfloor}\right)$$

线性筛预处理 $\mu,\operatorname{d}$，数论分块求解即可，复杂度 $O(n+T\sqrt{n})$。

---

代码

//知识点：莫比乌斯反演
/*
By:Luckyblock
*/
#include <algorithm>
#include <cctype>
#include <cmath>
#include <cstdio>
#include <cstring>
#define ll long long
const int kMaxn = 5e4 + 10;
//=============================================================
int pnum, p[kMaxn];
ll mu[kMaxn], num[kMaxn], d[kMaxn]; //num 为最小质因子的次数
ll summu[kMaxn], sumd[kMaxn];
bool vis[kMaxn];
//=============================================================
inline int read() {
  int f = 1, w = 0;
  char ch = getchar();
  for (; !isdigit(ch); ch = getchar())
    if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
void Getmax(int &fir_, int sec_) {
  if (sec_ > fir_) fir_ = sec_;
}
void Getmin(int &fir_, int sec_) {
  if (sec_ < fir_) fir_ = sec_;
}
void Euler(int lim_) {
  vis[1] = true;
  mu[1] = d[1] = 1ll;
  for (int i = 2; i <= lim_; ++ i) {
    if (! vis[i]) {
      p[++ pnum] = i;
      mu[i] = - 1;
      num[i] = 1;
      d[i] = 2;
    }
    for (int j = 1; j <= pnum && i * p[j] <= lim_; ++ j) {
      vis[i * p[j]] = true;
      if (i % p[j] == 0) { //勿忘平方因子
        mu[i * p[j]] = 0;
        num[i * p[j]] = num[i] + 1;
        d[i * p[j]] = 1ll * d[i] / num[i * p[j]] * (num[i * p[j]] + 1ll);
        break;
      }
      mu[i * p[j]] = - mu[i];
      num[i * p[j]] = 1;
      d[i * p[j]] = 2ll * d[i]; //
 
    }
  }
  for (int i = 1; i <= lim_; ++ i) {
    summu[i] = summu[i - 1] + mu[i];
    sumd[i] = sumd[i - 1] + d[i];
  }
}
 
//=============================================================
int main() { 
  Euler(kMaxn - 10);
  int T = read();
  while (T --) {
    int n = read(), m = read(), lim = std :: min(n, m);
    ll ans = 0ll;
    for (int l = 1, r; l <= lim; l = r + 1) {
      r = std :: min(n / (n / l), m / (m / l));
      ans += 1ll * (summu[r] - summu[l - 1]) * sumd[n / l] * sumd[m / l]; //
    }
    printf("%lld\n", ans);
  }
  return 0;
}
/*
in
1
32 43
*/
/*
out
15420
*/

---

P3768 简单的数学题

给定参数 $n$、$p$，求：

$$\left(\sum_{i=1}^n\sum_{j=1}^ni\cdot j\cdot \gcd(i,j)\right)\bmod p$$

$n\leq10^{10}$，$5\times10^8\leq p\leq1.1\times10^9$，$p\in \mathrm{primes}$。
时限 4s。

无脑套路暴力。

考虑先枚举 $\gcd(i,j)$，原式等价于：

$$\begin{aligned} &\sum_{d=1}d\sum_{i=1}^{n}\sum_{j=1}^{n}[\gcd(i,j)=d]ij\\ =& \sum_{d=1}d\sum_{i=1}^{n}[d|i]\sum_{j=1}^{n}[d|j][\gcd(\dfrac{i}{d},\dfrac{j}{d})=1]ij \end{aligned}$$

提出 $d$，变为枚举 $\frac{i}{d}$ 与 $\frac{j}{d}$，消去整除的条件，原式等价于：

$$\begin{aligned} &\sum_{d=1}d\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[\gcd(i,j)=1]id\cdot jd\\ =& \sum_{d=1} d^3\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[\gcd(i,j)=1]ij \end{aligned}$$

代入 $[\gcd(i,j) = 1] = \sum\limits_{d\mid \gcd(i,j)} {\mu (d)}$，原式等价于：

$$\sum_{d=1} d^3\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{n}{d}\right\rfloor}ij\sum_{t|\gcd(i,j)}\mu (t)$$

值得注意的是 $t$ 的上界为 $\frac{n}{d}$，$dt\le n$。
调换枚举顺序，先枚举 $t$，原式等价于：

$$\sum_{d=1} d^3\sum_{t=1}\mu(t)\sum_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[t|i] \sum_{j=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[t|j]ij$$

和上面一样，提出 $t$，套路地消去整除的条件，原式等价于：

$$\sum_{d=1} d^3\sum_{t=1}\mu(t)t^2\sum_{i=1}^{\left\lfloor\frac{n}{dt}\right\rfloor}\sum_{j=1}^{\left\lfloor\frac{n}{dt}\right\rfloor}ij$$

发现后面两个求和是等差数列乘等差数列的形式。
设 $g(p,q) = \sum\limits_{i=1}^{p} \sum\limits_{j=1}^{q}ij$，$g_{p,q}$ 可以通过下式 $O(1)$ 计算：

$$g(p,q) = \sum_{i=1}^{p} i \sum_{j=1}^{q}j= \dfrac{(1 + p)\times p}{2}\times \dfrac{(1+q)\times q}{2}$$

代入原式，原式等价于：

$$\sum_{d=1} d^3\sum_{t=1}\mu(t)t^2\cdot g\left({\left\lfloor\frac{n}{dt}\right\rfloor},{\left\lfloor\frac{n}{dt}\right\rfloor}\right)$$

考虑枚举 $T = dt$，显然 $T\le n$。
再考虑枚举 $d|T$，即可得到 $t = \frac{T}{d}$，原式等价于：

$$\begin{aligned} &\sum_{T=1}^{n}g\left({\left\lfloor\frac{n}{T}\right\rfloor},{\left\lfloor\frac{n}{T}\right\rfloor}\right)\sum_{d|T}d^3\mu{\left(\dfrac{T}{d}\right)}\left(\dfrac{T}{d}\right)^2\\ =& \sum_{T=1}^{n}T^2 g\left({\left\lfloor\frac{n}{T}\right\rfloor},{\left\lfloor\frac{n}{T}\right\rfloor}\right)\sum_{d|T}d\cdot \mu{\left(\dfrac{T}{d}\right)} \end{aligned}$$

对于后面这一坨，用 $\sum\limits_{d|T}d\cdot \mu{\left(\frac{T}{d}\right)} = \operatorname{Id} \ast \mu(T)= \varphi(T)$ 反演，则原式等价于：

$$\sum_{T=1}^{n}T^2 \varphi(T) \cdot g\left({\left\lfloor\frac{n}{T}\right\rfloor},{\left\lfloor\frac{n}{T}\right\rfloor}\right)$$

后半块可以数论分块，考虑前半块。
发现前半段即为 $\operatorname{Id}^2(T)\varphi(T)$，又是前缀和形式，考虑杜教筛。

有：

$$f(n) = \operatorname{Id}^2\varphi(n)$$

考虑找到一个函数 $g$，构造函数 $h = f\ast g$ 使其便于求值，有：

$$h(n) = \sum_{d|n} d^2\varphi(d)\cdot g\left(\dfrac{n}{d}\right)$$

看到同时存在 $d$ 和 $\frac{n}{d}$，考虑把 $d^2$ 消去。
令 $g = \operatorname{Id}^2$，有：

$$\begin{aligned} h(n) =& \sum_{d|n} d^2\varphi(d)\cdot \left(\dfrac{n}{d}\right)^2\\ =& n^2\sum_{d|n} \varphi(d)\\ =& n^2 \cdot \varphi \ast 1(n)\\ \end{aligned}$$

又 $\varphi \ast 1 = \operatorname{Id}$，则有：

$$h(n) = n^3$$

找到了合适的 $g$，套杜教筛的公式。

$$\begin{aligned} g(1)S(n) &= \sum_{i=1}^{n}h(i) - \sum_{i=2}^{n}g(i)S\left(\left\lfloor\dfrac{n}{i}\right\rfloor\right)\\ S(n) &= \sum_{i=1}^{n} i^3 - \sum_{i=2}^{n} i^2\cdot S\left(\left\lfloor\dfrac{n}{i}\right\rfloor\right) \end{aligned}$$

前一项是自然数的立方和，有 $\sum\limits_{i=1}^{n} i^3 = (\frac{n(n+1)}{2})^2$。证明详见：自然数前n项平方和、立方和公式及证明 - 百度文库。
后一项直接等差数列求和 + 数论分块求解即可。

---

「SDOI2017」数字表格

设 $f_{i}$ 表示斐波那契数列的第 $i$ 项。
$T$ 组数据，每次给定参数 $n,m$，求：

$$\prod_{i=1}^{n}\prod_{j=1}^{m}f_{\gcd(i,j)} \pmod {10^9 + 7}$$

$1\le T\le 10^3$，$1\le n,m\le 10^6$。
5S，256MB。

以下钦定 $n\ge m$。
大力化式子，先套路地枚举 $\gcd(i,j)$，用初中知识把两个 $\prod$ 化到指数位置，原式等于：

$$\large\begin{aligned} &\prod_{d = 1}^{n}\prod_{i=1}^{n}\prod_{j=1}^{m}f_{d}[\gcd(i,j)=d]\\ =&\prod_{d = 1}^{n}f_{d}^{\left(\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i,j)=d]\right)} \end{aligned}$$

分母套路一波，有：

$$\begin{aligned} &\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}[\gcd(i,j) = d]\\ =& \sum\limits_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}[\gcd(i,j) = 1]\\ =& \sum\limits_{i=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor\frac{m}{d}\right\rfloor}\sum_{k\mid \gcd(i,j)}\mu (k)\\ =& \sum_{k=1}\mu(k)\left\lfloor\dfrac{n}{kd}\right\rfloor\left\lfloor\dfrac{m}{kd}\right\rfloor \end{aligned}$$

代回原式，原式等于：

$$\large \begin{aligned} &\prod_{d = 1}^{n}f_{d}^{\left(\sum\limits_{k=1}\mu(k)\left\lfloor\frac{n}{kd}\right\rfloor\left\lfloor\frac{m}{kd}\right\rfloor\right)}\\ =& \prod_{d = 1}^{n}\left(f_{d}^{{\sum\limits_{k=1}\mu(k)}}\right)^{\left\lfloor\frac{n}{kd}\right\rfloor\left\lfloor\frac{m}{kd}\right\rfloor} \end{aligned}$$

考虑再暴力拆一波，原式等于：

$$\large \begin{aligned} &\prod_{d = 1}^{n}\left(f_{d}^{{\sum\limits_{k=1}\mu(k)}}\right)^{\left\lfloor\frac{n}{kd}\right\rfloor\left\lfloor\frac{m}{kd}\right\rfloor}\\ =& \prod_{d = 1}^{n}\left(\prod_{k=1}^{\left\lfloor\frac{n}{d}\right\rfloor}f_{d}^{\mu(k)\left\lfloor\frac{n}{kd}\right\rfloor\left\lfloor\frac{m}{kd}\right\rfloor}\right) \end{aligned}$$

做不动了，但发现变量仅有 $k,d,kd$，考虑更换枚举对象改为枚举 $t = kd$ 与 $d$，则原式等于：

$$\large\prod_{t=1}^{n}\left(\prod_{d | t}^{n}f_{d}^{{\mu(\frac{t}{d})}}\right)^{\left\lfloor\frac{n}{t}\right\rfloor\left\lfloor\frac{m}{t}\right\rfloor}$$

枚举对象变成了约数形式。从后面的式子推前面的式子是比较显然的，可以认为这种枚举 $t=kd$ 的形式是一种逆向操作。

设：

$$\large g(t)=\prod_{d | t}^{n}f_{d}^{{\mu(\frac{t}{d})}}$$

$g(t)$ 可以用类似埃氏筛的方法 $O(n\log ^2 n)$ 地预处理出来。再把 $g$ 代回原式，原式等于：

$$\large\prod_{t=1}^{n}g(t)^{\left\lfloor\frac{n}{t}\right\rfloor\left\lfloor\frac{m}{t}\right\rfloor}$$

可以考虑预处理 $g(t)$ 的前缀积，数论分块枚举指数求解即可。

总时间复杂度 $O(n\log ^2 n + T\sqrt n)$，轻微卡常可以过。

//知识点：莫比乌斯反演 
/*
By:Luckyblock
*/
#include <algorithm>
#include <cctype>
#include <cstdio>
#include <cstring>
#define LL long long
const LL mod = 1e9 + 7;
const int kN = 1e6;
//=============================================================
LL n, m, ans;
int p_num, p[kN + 10];
bool vis[kN + 10];
LL mu[kN + 10], f[kN + 10], g[kN + 10], prod[kN + 10];
LL invf[kN + 10], invp[kN];
//=============================================================
inline int read() {
  int f = 1, w = 0;
  char ch = getchar();
  for (; !isdigit(ch); ch = getchar())
    if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) {
    w = (w << 3) + (w << 1) + (ch ^ '0');
  }
  return f * w;
}
void Chkmax(int &fir_, int sec_) {
  if (sec_ > fir_) fir_ = sec_;
}
void Chkmin(int &fir_, int sec_) {
  if (sec_ < fir_) fir_ = sec_;
}
LL QPow(LL x_, LL y_) {
  x_ %= mod;
  LL ret = 1;
  for (; y_; y_ >>= 1ll) {
    if (y_ & 1) ret = ret * x_ % mod;
    x_ = x_ * x_ % mod;
  }
  return ret;
}
void Euler() {
  vis[1] = true, mu[1] = 1;
  for (int i = 2; i <= kN; ++ i) {
    if (! vis[i]) {
      p[++ p_num] = i;
      mu[i] = -1;
    }
    for (int j = 1; j <= p_num && i * p[j] <= kN; ++ j) {
      vis[i * p[j]] = true;
      if (i % p[j] == 0) {
        mu[i * p[j]] = 0;
        break;
      }
      mu[i * p[j]] = -mu[i];
    }
  }
}
void Prepare() {
  g[1] = g[2] = 1;
  f[1] = f[2] = 1;
  invf[1] = invf[2] = 1;
  for (int i = 3; i <= kN; ++ i) {
    g[i] = 1;
    f[i] = (f[i - 1] + f[i - 2]) % mod;
    invf[i] = QPow(f[i], mod - 2);
  }
  
  Euler();
  for (int d = 1; d <= kN; ++ d) {
    for (int j = 1; d * j <= kN; ++ j) {
      if (mu[j] == 1) {
        g[d * j] = g[d * j] * f[d] % mod;
      } else if (mu[j] == -1) {
        g[d * j] = g[d * j] * invf[d] % mod;
      }
    }
  }
  invp[0] = prod[0] = 1;
  for (int i = 1; i <= kN; ++ i) {
    prod[i] = prod[i - 1] * g[i] % mod;
    invp[i] = QPow(prod[i], mod - 2);
  }
}
//=============================================================
int main() {
  Prepare();
  int T = read();
  while (T -- ){
    n = read(), m = read(), ans = 1;
    if (n < m) std::swap(n, m);
    for (LL l = 1, r = 1; l <= m; l = r + 1) {
      r = std::min(n / (n / l), m / (m / l));
      ans = (ans * QPow(prod[r] * invp[l - 1] % mod, 1ll * (n / l) * (m / l))) % mod;
    }
    printf("%lld\n", ans);
  }
  return 0;
}

---

P5518 [MtOI2019]幽灵乐团 / 莫比乌斯反演基础练习题

给定参数 $p$，有 $T$ 组数据，每次给定参数 $A,B,C$，求：

$$\prod_{i=1}^{A}\prod_{j=1}^{B}\prod_{k=1}^{C}\left(\dfrac{\operatorname{lcm}(i,j)}{\gcd(i,k)}\right)^{f(type)}$$

其中 $f(type)$ 的取值如下：

$$f(type) = \begin{cases} 1 &type = 0\\ i\times j\times k &type = 1\\ \gcd(i,j,k) &type = 2 \end{cases}$$

$1\le A,B,C\le 10^5$，$10^7\le p\le 1.05\times 10^9$，$p\in \mathbb{P}$，$T=70$。
2.5S，128MB。

先化下原式，原式等于：

$$\prod_{i=1}^{A}\prod_{j=1}^{B}\prod_{k=1}^{C}\left(\dfrac{i\times j }{\gcd(i,j)\times \gcd(i,k)}\right)^{f(type)}$$

发现每一项仅与两个变量有关，设：

$$\begin{aligned} f_1(a,b,c) &= \prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} i^{f(type)}\\ f_2(a,b,c) &= \prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} \gcd(i,j)^{f(type)} \end{aligned}$$

发现 $\prod$ 可以随意交换，则原式等价于：

$$\dfrac{f_1(a,b,c)\times f_1(b,a,c)}{f_2(a,b,c)\times f_2(a,c,b)}$$

考虑在 $type$ 取值不同时，如何快速求得 $f_1$ 与 $f_2$。
一共有 $6$ 个需要推导的式子，不妨就叫它们 $1\sim 6$ 面了（

---

type = 0

$$\begin{aligned} f_1(a,b,c) &= \prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} i\\ f_2(a,b,c) &= \prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} \gcd(i,j) \end{aligned}$$

对于 1 面，显然有：

$$\prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} i = \left(\prod_{i=1}^{a}i\right)^{b\times c}$$

预处理阶乘 + 快速幂即可，单次计算时间复杂度 $O(\log n)$。

---

再考虑 2 面，套路地枚举 $\gcd$，显然有：

$$\large \begin{aligned} &\prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} \gcd(i,j)\\ =&\left(\prod_{i=1}^{a}\prod_{j=1}^{b}\gcd(i,j)\right)^c\\ =& \left(\prod_{d=1} d^{\left(\sum\limits_{i=1}^{a}\sum\limits_{j=1}^{b}[\gcd(i,j) = d]\right)}\right)^c \end{aligned}$$

指数是个套路，可以看这里：P3455 [POI2007]ZAP-Queries。于是有：

$$\begin{aligned} &\sum\limits_{i=1}^{a}\sum\limits_{j=1}^{b}[\gcd(i,j) = d]\\ =& \sum\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor}[\gcd(i,j) = 1]\\ =& \sum\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor}\sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor}\sum_{k\mid \gcd(i,j)}\mu (k)\\ =& \sum_{k=1}\mu(k)\left\lfloor\dfrac{a}{kd}\right\rfloor\left\lfloor\dfrac{b}{kd}\right\rfloor \end{aligned}$$

代回原式，略做处理，则原式等于：

$$\large \begin{aligned} &\left(\prod_{d=1} d^{\left(\sum\limits_{k=1}\mu(k)\left\lfloor\frac{a}{kd}\right\rfloor\left\lfloor\frac{b}{kd}\right\rfloor\right)}\right)^c\\ =& \left(\prod_{d=1} \left(d^{\sum\limits_{k=1}\mu(k)}\right)^{\left\lfloor\frac{a}{kd}\right\rfloor\left\lfloor\frac{b}{kd}\right\rfloor}\right)^c\\ =& \prod_{d=1} \left(\prod_{k=1}^{\left\lfloor\frac{n}{d}\right\rfloor}d^{\left(\mu(k)\left\lfloor\frac{a}{kd}\right\rfloor\left\lfloor\frac{b}{kd}\right\rfloor\right)}\right)^c \end{aligned}$$

像「SDOI2017」数字表格 一样，考虑枚举 $t=kd$ 和 $d$，则原式等于：

$$\large \prod_{t=1}^{n}\left(\left(\prod_{d|t} d^{\mu{\left(\frac{t}{d}\right)}}\right)^{\left\lfloor\frac{a}{t}\right\rfloor\left\lfloor\frac{b}{t}\right\rfloor}\right)^c$$

设：

$$\large g_0(t) = \prod_{d|t}d^{\mu\left(\frac{t}{d}\right)}$$

线性筛预处理 $\mu$ 后，$g_0(t)$ 可以用埃氏筛预处理，时间复杂度 $O(n\log n)$。再代回原式，原式等于：

$$\large \prod_{t=1}^{a}\left(g_0(t)^{\left\lfloor\frac{a}{t}\right\rfloor\left\lfloor\frac{b}{t}\right\rfloor}\right)^c$$

预处理 $g_0(t)$ 的前缀积和前缀积的逆元，复杂度 $O(n\log n)$。
数论分块 + 快速幂计算即可，单次时间复杂度 $O(\sqrt n\log n)$。

---

type = 1

$$\begin{aligned} f_1(a,b,c) &= \prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} i^{i\times j\times k}\\ f_2(a,b,c) &= \prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} \gcd(i,j)^{i\times j\times k} \end{aligned}$$

考虑 3 面，把 $\prod k$ 扔到指数位置，有：

$$\large \prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} i^{i\times j\times k} = \prod_{i=1}^{a}\prod_{j=1}^{b}i^{\left(i\times j\times \sum\limits_{k = 1}^{c} k\right)}$$

再把 $\prod j$ 也扔到指数位置，引入 $\operatorname{sum}(n) = \sum_{i=1}^{n} i = \frac{n(n+1)}{2}$，原式等于：

$$\left(\prod_{i=1}^{a}i^i\right)^{\operatorname{sum}(b)\times \operatorname{sum}(c)}$$

预处理 $i^i$ 的前缀积，复杂度 $O(n\log n)$。
指数可以 $O(1)$ 算出，然后快速幂，单次时间复杂度 $O(\log n)$。

根据费马小定理，指数需要对 $p - 1$ 取模。注意 $p-1$ 不是质数，计算 $\operatorname{sum}$ 时不能用逆元，但乘不爆 LL，直接算就行。

---

再考虑 4 面，发现 $k$ 与 $\gcd$ 无关，则同样把 $\prod k$ 扔到指数位置，则有：

$$\large \begin{aligned} &\prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} \gcd(i,j)^{i\times j\times k}\\ =& \left(\prod_{i=1}^a\prod_{j=1}^b\gcd(i,j)^{i\times j}\right)^{\operatorname{sum}(c)} \end{aligned}$$

套路地枚举 $\gcd$，原式等于：

$$\large \left(\prod_{d=1}d^{\left(\sum\limits_{i=1}^a \sum\limits_{j=1}^b i\times j[\gcd(i,j)=d]\right)}\right)^{\operatorname{sum}(c)}$$

大力化简指数，有：

$$\large \begin{aligned} &\sum\limits_{i=1}^a \sum\limits_{j=1}^b i\times j[\gcd(i,j)=d]\\ =& d^2 \sum\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor} \sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor} i\times j[\gcd(i,j)=1\\ =& d^2 \sum\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor} i \sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor} j\sum\limits_{t|\gcd(i,j)}\mu(t)\\ =& d^2 \sum\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor} i \sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor} j\sum\limits_{k|\gcd(i,j)}\mu(k)\\ =& d^2 \sum\limits_{k=1}\mu(k)\sum\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor} i[k|i] \sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor} j[k|j]\\ =& d^2 \sum\limits_{k=1}k^2\mu(k)\sum\limits_{i=1}^{\left\lfloor\frac{a}{kd}\right\rfloor} i\sum\limits_{j=1}^{\left\lfloor\frac{b}{kd}\right\rfloor} j\\ =& d^2\sum\limits_{k=1}k^2\mu(k)\operatorname{sum}{\left(\left\lfloor\frac{a}{kd}\right\rfloor\right)} \operatorname{sum}{\left(\left\lfloor\frac{b}{kd}\right\rfloor\right)}\\ \end{aligned}$$

指数化不动了，代回原式，原式等于：

$$\large \left(\prod_{d=1}d^{\left(d^2\sum\limits_{k=1}k^2\mu(k)\operatorname{sum}{\left(\left\lfloor\frac{a}{kd}\right\rfloor\right)} \operatorname{sum}{\left(\left\lfloor\frac{b}{kd}\right\rfloor\right)}\right)}\right)^{\operatorname{sum}(c)}$$

同 2 面的情况，先展开一下，再枚举 $t=kd$ 和 $d$，原式等于：

$$\large \begin{aligned} &\left(\prod_{d=1}\left(\prod_{k=1}^{\left\lfloor\frac{n}{d}\right\rfloor}d^{\left(d^2 k^2\mu(k)\right)}\right)^{\left(\operatorname{sum}{\left(\left\lfloor\frac{a}{kd}\right\rfloor\right)} \operatorname{sum}{\left(\left\lfloor\frac{b}{kd}\right\rfloor\right)}\right)}\right)^{\operatorname{sum}(c)}\\ =& \prod_{t=1}\left(\left(\prod_{d|t}d^{\left(d^2\left(\frac{t}{d}\right)^2\mu\left(\frac{t}{d}\right)\right)}\right)^{\operatorname{sum}{\left(\left\lfloor\frac{a}{t}\right\rfloor\right)} \operatorname{sum}{\left(\left\lfloor\frac{b}{t}\right\rfloor\right)}}\right)^{\operatorname{sum}(c)}\\ =& \prod_{t=1}\left(\left(\prod_{d|t}d^{\left(t^2\mu\left(\frac{t}{d}\right)\right)}\right)^{\operatorname{sum}{\left(\left\lfloor\frac{a}{t}\right\rfloor\right)} \operatorname{sum}{\left(\left\lfloor\frac{b}{t}\right\rfloor\right)}}\right)^{\operatorname{sum}(c)} \end{aligned}$$

与二面相同，设：

$$\large g_1(t) = \prod_{d|t}d^{\left(t^2\mu\left(\frac{t}{d}\right)\right)}$$

$g_1(t)$ 可以用埃氏筛套快速幂筛出，时间复杂度 $O(n\log^2 n)$。再代回原式，原式等于：

$$\prod_{t=1}\left(g_1(t)^{\operatorname{sum}{\left(\left\lfloor\frac{a}{t}\right\rfloor\right)} \operatorname{sum}{\left(\left\lfloor\frac{b}{t}\right\rfloor\right)}}\right)^{\operatorname{sum}(c)}$$

同样预处理 $g_1(t)$ 的前缀积及其逆元，时间复杂度 $O(n\log n)$。
整除分块 + 快速幂即可，单次时间复杂度 $O(\sqrt n\log n)$。

注意指数的取模。

---

type = 2

$$\begin{aligned} f_1(a,b,c) &= \prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} i^{\gcd(i,j,k)}\\ f_2(a,b,c) &= \prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} \gcd(i,j)^{\gcd(i,j,k)} \end{aligned}$$

考虑 5 面，手段同上，大力反演化简一波，再调换枚举对象，则有：

$$\large \begin{aligned} &\prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} i^{\gcd(i,j,k)}\\ =&\prod_{d=1}\prod\limits_{i=1}^{a}i^{\left(\sum\limits_{j=1}^{b}\sum\limits_{k=1}^{c}[\gcd(i,j,k)=d]\right)}\\ =& \prod_{d=1}\prod\limits_{i=1}^{a}i^{\left(\sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor}\sum\limits_{k=1}^{\left\lfloor\frac{c}{d}\right\rfloor}[\gcd(\frac{i}{d},j,k)=1]\right)}\\ =& \prod_{d=1}\prod\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor}(id)^{\left(d\sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor}\sum\limits_{k=1}^{\left\lfloor\frac{c}{d}\right\rfloor}[\gcd(i,j,k)=1]\right)}\\ =& \prod_{d=1}\prod\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor}(id)^{\left(d\sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor}\sum\limits_{k=1}^{\left\lfloor\frac{c}{d}\right\rfloor}\sum\limits_{x|\gcd(i,j,k)}{\mu(x)}\right)}\\ =& \prod_{d=1}\prod\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor}(id)^{\left(d\sum\limits_{x=1}\mu(x)[x|i]\sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor}[x|j]\sum\limits_{k=1}^{\left\lfloor\frac{c}{d}\right\rfloor}[x|k]\right)}\\ =& \prod_{d=1}\prod_{x=1}\prod\limits_{i=1}^{\left\lfloor\frac{a}{d}\right\rfloor}(id)^{\left(d\times \mu(x)[x|i]\sum\limits_{j=1}^{\left\lfloor\frac{b}{d}\right\rfloor}[x|j]\sum\limits_{k=1}^{\left\lfloor\frac{c}{d}\right\rfloor}[x|k]\right)}\\ =& \prod_{d=1}\prod_{x=1}\prod\limits_{i=1}^{\left\lfloor\frac{a}{xd}\right\rfloor}(ixd)^{\left(d\times \mu(x){\left\lfloor\frac{b}{xd}\right\rfloor}{\left\lfloor\frac{c}{xd}\right\rfloor}\right)}\\ =& \prod_{t = 1}\prod_{d|T}\prod_{i=1}^{\left\lfloor\frac{a}{t}\right\rfloor}(it)^{\left(d\times \mu\left(\frac{t}{d}\right){\left\lfloor\frac{b}{t}\right\rfloor}{\left\lfloor\frac{c}{t}\right\rfloor}\right)}\\ =& \prod_{t = 1}\prod_{d|T}\left(t^{\left\lfloor\frac{a}{t}\right\rfloor}\prod_{i=1}^{\left\lfloor\frac{a}{t}\right\rfloor}i\right)^{d\times \mu\left(\frac{t}{d}\right){\left\lfloor\frac{b}{t}\right\rfloor}{\left\lfloor\frac{c}{t}\right\rfloor}}\\ \end{aligned}$$

引入 $\operatorname{fac}(n) = \prod_{i=1}^{n} i$，再根据枚举对象调整一下指数，原式等于：

$$\large \begin{aligned} &\prod_{t = 1}\prod_{d|t}\left(t^{\left\lfloor\frac{a}{t}\right\rfloor}\times \operatorname{fac}\left(\left\lfloor\frac{a}{t}\right\rfloor\right)\right)^{\left(d\times \mu\left(\frac{t}{d}\right){\left\lfloor\frac{b}{t}\right\rfloor}{\left\lfloor\frac{c}{t}\right\rfloor}\right)}\\ =& \prod_{t = 1}\left(\prod_{d|t}\left(t^{\left\lfloor\frac{a}{t}\right\rfloor}\times \operatorname{fac}\left(\left\lfloor\frac{a}{t}\right\rfloor\right)\right)^{d\times \mu\left(\frac{t}{d}\right)}\right)^{{\left\lfloor\frac{b}{t}\right\rfloor}{\left\lfloor\frac{c}{t}\right\rfloor}}\\ =& \prod_{t = 1}\left(\left(t^{\left\lfloor\frac{a}{t}\right\rfloor}\times \operatorname{fac}\left(\left\lfloor\frac{a}{t}\right\rfloor\right)\right)^{\sum\limits_{d|t}d\times \mu\left(\frac{t}{d}\right)}\right)^{{\left\lfloor\frac{b}{t}\right\rfloor}{\left\lfloor\frac{c}{t}\right\rfloor}} \end{aligned}$$

指数中出现了一个经典的狄利克雷卷积的形式，对其进行反演。
将 $(\operatorname{Id}\ast \mu) (n)= \varphi (n)$ 代入原式，原式等于：

$$\large \begin{aligned} &\prod_{t = 1}\left(t^{\left\lfloor\frac{a}{t}\right\rfloor}\times \operatorname{fac}\left(\left\lfloor\frac{a}{t}\right\rfloor\right)\right)^{\varphi(t){\left\lfloor\frac{b}{t}\right\rfloor}{\left\lfloor\frac{c}{t}\right\rfloor}}\\ =& \prod_{t = 1}\left(t^{\varphi(t)\left\lfloor\frac{a}{t}\right\rfloor{\left\lfloor\frac{b}{t}\right\rfloor}{\left\lfloor\frac{c}{t}\right\rfloor}}\times \operatorname{fac}\left(\left\lfloor\frac{a}{t}\right\rfloor\right)^{\varphi(t){\left\lfloor\frac{b}{t}\right\rfloor}{\left\lfloor\frac{c}{t}\right\rfloor}}\right) \end{aligned}$$

预处理 $t^{\varphi(t)}$ 的前缀积及逆元，阶乘的前缀积及阶乘逆元，$\pmod {p-1}$ 下的 $\varphi(t)$ 的前缀和（指数
），时间复杂度 $O(n\log n)$。
同样整除分块 + 快速幂即可，单次时间复杂度 $O(\sqrt n\log n)$。

---

然后是最掉 sans 的 6 面。有 $\gcd(i,j,k) = \gcd(\gcd(i,j), k)$，考虑先枚举 $\gcd(i,j)$，然后套路化式子，则有：

$$\large \begin{aligned} &\prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} \gcd(i,j)^{\gcd(i,j,k)}\\ =& \prod_{d=1}\prod_{i=1}^{a}\prod_{j=1}^{b}\prod_{k=1}^{c} [\gcd(i,j)=d] d^{\gcd(d,k)}\\ =& \prod_{d=1} \left(d^{\left(\sum\limits_{i=1}^{a}\sum\limits_{j=1}^{b} [\gcd(i,j)=d]\right)}\right)^{\sum\limits_{k=1}^{c}\gcd(d,k)} \end{aligned}$$

先考虑最外面的指数，这也是个套路，可以参考 一个例子。用 $\operatorname{Id} = \varphi \ast 1$ 反演，显然有：

$$\large \begin{aligned} &\sum\limits_{k=1}^{c}\gcd(d,k)\\ =& \sum\limits_{k=1}^{c}\sum_{x|\gcd(d,k)}\varphi(x)\\ =& \sum_{x=1}\varphi(x)[x|d]\sum_{k=1}^{c}[x|k]\\ =& \sum_{x|d}\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor \end{aligned}$$

再考虑里面的指数，发现这式子在 2 面已经推了一遍了，于是直接拿过来用，有：

$$\sum\limits_{i=1}^{a}\sum\limits_{j=1}^{b}[\gcd(i,j) = d]=\sum_{y=1}\mu(y)\left\lfloor\dfrac{a}{yd}\right\rfloor\left\lfloor\dfrac{b}{yd}\right\rfloor$$

将化简后的两个指数代入原式，原式等于：

$$\large \begin{aligned} &\prod_{d=1} \left(d^{\left(\sum\limits_{y=1}\mu(y)\left\lfloor\frac{a}{yd}\right\rfloor\left\lfloor\frac{b}{yd}\right\rfloor\right)}\right)^{\sum\limits_{x|d}\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor}\\ =& \prod_{d=1} \left(\prod\limits_{y=1}d^{\left(\mu(y)\left\lfloor\frac{a}{yd}\right\rfloor\left\lfloor\frac{b}{yd}\right\rfloor\right)}\right)^{\sum\limits_{x|d}\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor} \end{aligned}$$

与 2、4 面同样套路地，考虑枚举 $t=yd$ 和 $d$，再略作调整，原式等于：

$$\large \begin{aligned} &\prod_{d=1} \left(\prod\limits_{y=1}d^{\left(\mu(y)\left\lfloor\frac{a}{yd}\right\rfloor\left\lfloor\frac{b}{yd}\right\rfloor\right)}\right)^{\sum\limits_{x|d}\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor}\\ =& \prod_{t=1}\prod_{d|t} d^{\left(\mu\left(\frac{t}{d}\right)\left\lfloor\frac{a}{t}\right\rfloor\left\lfloor\frac{b}{t}\right\rfloor\sum\limits_{x|d}\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\right)}\\ =& \prod_{t=1}\left(\prod_{d|t} d^{\left(\mu\left(\frac{t}{d}\right)\sum\limits_{x|d}\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\right)}\right)^{\left\lfloor\frac{a}{t}\right\rfloor\left\lfloor\frac{b}{t}\right\rfloor}\\ =& \prod_{t=1}\left(\prod_{d|t} \prod_{x|d}d^{\left(\mu\left(\frac{t}{d}\right)\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\right)}\right)^{\left\lfloor\frac{a}{t}\right\rfloor\left\lfloor\frac{b}{t}\right\rfloor} \end{aligned}$$

发现要同时枚举 $d$ 和 $x$，化不动了。
从题解里学到一个比较神的技巧，考虑把 $d$ 拆成 $x$ 和 $\frac{d}{x}$ 分别计算贡献再相乘，即分别计算下两式：

$$\large \begin{aligned} &\prod_{t=1}\left(\prod_{d|t} \prod_{x|d}x^{\left(\mu\left(\frac{t}{d}\right)\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\right)}\right)^{\left\lfloor\frac{a}{t}\right\rfloor\left\lfloor\frac{b}{t}\right\rfloor}\\ &\prod_{t=1}\left(\prod_{d|t} \prod_{x|d}{\left(\frac{d}{x}\right)}^{\left(\mu\left(\frac{t}{d}\right)\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\right)}\right)^{\left\lfloor\frac{a}{t}\right\rfloor\left\lfloor\frac{b}{t}\right\rfloor} \end{aligned}$$

---

先考虑 $x$ 的情况，首先把枚举 $x$ 调整到最外层。设 $\operatorname{lim}=\max(a,b,c)$，则原式等于：

$$\large \begin{aligned} &\prod_{x=1} \prod_{t=1}^{\operatorname{lim}}[x|t]\left(\prod_{d|t} [x|d]{x}^{\left(\mu\left(\frac{t}{d}\right)\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\right)}\right)^{\left\lfloor\frac{a}{t}\right\rfloor\left\lfloor\frac{b}{t}\right\rfloor}\\ =& \prod_{x=1} \prod_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}\left(\prod_{d|t} {x}^{\left(\mu\left(\frac{tx}{dx}\right)\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\right)}\right)^{\left\lfloor\frac{a}{tx}\right\rfloor\left\lfloor\frac{b}{tx}\right\rfloor}\\ =& \prod_{x=1} \prod_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}\prod_{d|t} {x}^{\left(\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\left\lfloor\frac{a}{tx}\right\rfloor\left\lfloor\frac{b}{tx}\right\rfloor\mu\left(\frac{t}{d}\right)\right)} \end{aligned}$$

把 $\prod {t}$ 挪到指数位置，原式等于：

$$\large \begin{aligned} &\prod_{x=1} {x}^{\left(\sum\limits_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}\sum\limits_{d|t}\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\left\lfloor\frac{a}{tx}\right\rfloor\left\lfloor\frac{b}{tx}\right\rfloor\mu\left(\frac{t}{d}\right)\right)}\\ =& \prod_{x=1} {x}^{\left(\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\sum\limits_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}\left\lfloor\frac{a}{tx}\right\rfloor\left\lfloor\frac{b}{tx}\right\rfloor\sum\limits_{d|t}\mu\left(\frac{t}{d}\right)\right)} \end{aligned}$$

指数中又出现了一个经典的狄利克雷卷积的形式，对其进行反演。
将 $(\mu \ast 1) (n)= \epsilon (n)=[n=1]$ 代入原式，原式等于：

$$\large \begin{aligned} &\prod_{x=1} {x}^{\left(\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\sum\limits_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}\left\lfloor\frac{a}{tx}\right\rfloor\left\lfloor\frac{b}{tx}\right\rfloor[t=1]\right)}\\ =& \prod_{x=1} {x}^{\left(\varphi(x)\left\lfloor\frac{a}{x}\right\rfloor\left\lfloor\frac{b}{x}\right\rfloor\left\lfloor\frac{c}{x}\right\rfloor\right)} \end{aligned}$$

得到了一个非常优美的式子，而且发现这个式子是 5 面最终答案的一部分。同 5 面的做法，直接整除分块即可。

---

再考虑 $\frac{d}{x}$ 的情况，同上先把枚举 $x$ 放到最外层，并调整一下指数，则原式等于：

$$\large \begin{aligned} &\prod_{x=1} \prod_{t=1}^{\operatorname{lim}}[x|t]\left(\prod_{d|t} [x|d]{\left(\frac{d}{x}\right)}^{\left(\mu\left(\frac{t}{d}\right)\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\right)}\right)^{\left\lfloor\frac{a}{t}\right\rfloor\left\lfloor\frac{b}{t}\right\rfloor}\\ =& \prod_{x=1} \prod_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}\left(\prod_{d|tx} [x|d]{\left(\frac{d}{x}\right)}^{\left(\mu\left(\frac{tx}{d}\right)\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor\right)}\right)^{\left\lfloor\frac{a}{tx}\right\rfloor\left\lfloor\frac{b}{tx}\right\rfloor}\\ =& \prod_{x=1} \left(\prod_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}\left(\prod_{d|tx} [x|d]{\left(\frac{d}{x}\right)}^{\mu\left(\frac{tx}{d}\right)}\right)^{\left\lfloor\frac{a}{tx}\right\rfloor\left\lfloor\frac{b}{tx}\right\rfloor}\right)^{\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor} \end{aligned}$$

考虑枚举 $dx$，替换原来的 $d$，注意一下这里的倍数关系。原式等于：

$$\large \prod_{x=1} \left(\prod_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}\left(\prod_{d|t}d^{\mu\left(\frac{t}{d}\right)}\right)^{\left\lfloor\frac{a}{tx}\right\rfloor\left\lfloor\frac{b}{tx}\right\rfloor}\right)^{\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor}$$

发现最内层的式子 $\prod_{d|t}d^{\mu\left(\frac{t}{d}\right)}$，即为二面处理过的 $g_0(t)$。直接代入，原式等于：

$$\large \prod_{x=1} \left(\prod_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}g_0(t)^{\left\lfloor\frac{a}{tx}\right\rfloor\left\lfloor\frac{b}{tx}\right\rfloor}\right)^{\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor}$$

一个小结论，证明可以看 这里：

$$\forall a,b,c\in \mathbb{Z},\left\lfloor\dfrac{a}{bc}\right\rfloor = \left\lfloor{\dfrac{\left\lfloor\dfrac{a}{b}\right\rfloor}{c}}\right\rfloor$$

则原式等于：

$$\large \prod_{x=1} \left(\prod_{t=1}^{\left\lfloor\frac{\operatorname{lim}}{x}\right\rfloor}g_0(t)^{\left\lfloor\frac{\left\lfloor\frac{a}{x}\right\rfloor}{t}\right\rfloor\left\lfloor\frac{\left\lfloor\frac{b}{x}\right\rfloor}{t}\right\rfloor}\right)^{\varphi(x)\left\lfloor\frac{c}{x}\right\rfloor}$$

于是可以先对外层整除分块，再对内层整除分块。

然后就做完了，哈哈哈。

---

一些实现上的小技巧：

• 逆元能预处理就预处理。
• 注意对指数取模，模数为 $p-1$。

//知识点：莫比乌斯反演 
/*
By:Luckyblock
用了比较清晰易懂的变量名，阅读体验应该会比较好。  
vsc 的自动补全真是太好啦！
*/
#include <algorithm>
#include <cctype>
#include <cstdio>
#include <cstring>
using std::min;
using std::max;
#define LL long long
const int Lim = 1e5;
const int kN = 1e5 + 10;
//=============================================================
LL A, B, C, mod, ans;
int T, p_num, p[kN];
bool vis[kN];
LL mu[kN], phi[kN], fac[kN], g[2][kN];
LL sumphi[kN], prodt_phi[kN], prodi_i[kN], prodg[2][kN];
LL inv[kN], inv_fac[kN], inv_prodt_phi[kN], inv_prodg[2][kN];
//=============================================================
inline int read() {
  int f = 1, w = 0;
  char ch = getchar();
  for (; !isdigit(ch); ch = getchar())
    if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) {
    w = (w << 3) + (w << 1) + (ch ^ '0');
  }
  return f * w;
}
void Chkmax(int &fir_, int sec_) {
  if (sec_ > fir_) fir_ = sec_;
}
void Chkmin(int &fir_, int sec_) {
  if (sec_ < fir_) fir_ = sec_;
}
LL QPow(LL x_, LL y_) {
  x_ %= mod;
  y_ %= mod - 1;
  LL ret = 1;
  for (; y_; y_ >>= 1ll) {
    if (y_ & 1) ret = ret * x_ % mod;
    x_ = x_ * x_ % mod;
  }
  return ret;
}
LL Inv(LL x_) {
  return QPow(x_, mod - 2);
}
LL Sum(LL n_) {
  return ((n_ * (n_ + 1ll)) / 2ll) % (mod - 1);
}
void Euler() {
  vis[1] = true, mu[1] = phi[1] = 1; //初值
  for (int i = 2; i <= Lim; ++ i) {
    if (! vis[i]) {
      p[++ p_num] = i;
      mu[i] = -1;
      phi[i] = i - 1;
    }
    for (int j = 1; j <= p_num && i * p[j] <= Lim; ++ j) {
      vis[i * p[j]] = true;
      if (i % p[j] == 0) {
        mu[i * p[j]] = 0;
        phi[i * p[j]] = phi[i] * p[j];
        break;
      }
      mu[i * p[j]] = -mu[i];
      phi[i * p[j]] = phi[i] * (p[j] - 1);
    }
  }
}
void Prepare() {
  Euler();
  inv[1] = fac[0] = prodt_phi[0] = prodi_i[0] = 1;
  for (int i = 1; i <= Lim; ++ i) {
    g[0][i] = g[1][i] = 1;
    fac[i] = 1ll * fac[i - 1] * i % mod;
    sumphi[i] = (sumphi[i - 1] + phi[i]) % (mod - 1);
    prodi_i[i] = prodi_i[i - 1] * QPow(i, i) % mod;
    if (i > 1) inv[i] = (mod - mod / i) * inv[mod % i] % mod;
 
    prodt_phi[i] = prodt_phi[i - 1] * QPow(i, phi[i]) % mod;
    inv_prodt_phi[i] = Inv(prodt_phi[i]);
  }
 
  for (int d = 1; d <= Lim; ++ d) {
    for (int j = 1; d * j <= Lim; ++ j) {
      int t = d * j;
      if (mu[j] == 1) {
        g[0][t] = g[0][t] * d % mod;
        g[1][t] = g[1][t] * QPow(1ll * d, 1ll * t * t) % mod;
      } else if (mu[j] == -1) {
        g[0][t] = g[0][t] * inv[d] % mod;
        g[1][t] = g[1][t] * Inv(QPow(1ll * d, 1ll * t * t)) % mod;
      }
    }
  }
  inv_prodg[0][0] = prodg[0][0] = 1;
  inv_prodg[1][0] = prodg[1][0] = 1;
  inv_prodt_phi[0] = 1;
  for (int i = 1; i <= Lim; ++ i) {
    for (int j = 0; j <= 1; ++ j) {
      prodg[j][i] = prodg[j][i - 1] * g[j][i] % mod;
      inv_prodg[j][i] = Inv(prodg[j][i]);
    }
  }
}
LL f1(LL a_, LL b_, LL c_, int type) {
  if (! type) return QPow(fac[a_], b_ * c_);
  if (type == 1) return QPow(prodi_i[a_], Sum(b_) * Sum(c_));
  LL ret = 1, lim = min(min(a_, b_), c_);
  for (LL l = 1, r = 1; l <= lim; l = r + 1) {
    r = min(min(a_ / (a_ / l), b_ / (b_ / l)), c_ / (c_ / l));
    ret = ret * QPow(prodt_phi[r] * inv_prodt_phi[l - 1], (a_ / l) * (b_ / l) % (mod - 1) * (c_ / l)) % mod;
    ret = ret * QPow(fac[a_ / l], (sumphi[r] - sumphi[l - 1] + mod - 1) % (mod - 1) * (b_ / l) % (mod - 1) * (c_ / l)) % mod;
  }
  return ret;
}
LL f2_2(LL a_, LL b_) { 
  LL ret = 1;
  for (LL l = 1, r = 1; l <= min(a_, b_); l = r + 1) {
    r = min(a_ / (a_ / l), b_ / (b_ / l));
    ret = ret * QPow(prodg[0][r] * inv_prodg[0][l - 1], 1ll * (a_ / l) * (b_ / l)) % mod;
  }
  return ret;
}
LL f2(LL a_, LL b_, LL c_, int type) {
  LL ret = 1;
  if (! type) {
    for (LL l = 1, r = 1; l <= min(a_, b_); l = r + 1) {
      r = min(a_ / (a_ / l), b_ / (b_ / l));
      LL val = QPow(prodg[0][r] * inv_prodg[0][l - 1], 1ll * (a_ / l) * (b_ / l));
      ret = (ret * QPow(val, c_)) % mod;
    }
  } else if (type == 1) {
    for (LL l = 1, r = 1; l <= min(a_, b_); l = r + 1) {
      r = min(a_ / (a_ / l), b_ / (b_ / l));
      LL val = QPow(prodg[1][r] * inv_prodg[1][l - 1], Sum(a_ / l) * Sum(b_ / l));
      ret = (ret * QPow(val, Sum(c_))) % mod;
    }
  } else {
    LL lim = min(min(a_, b_), c_);
    for (LL l = 1, r = 1; l <= lim; l = r + 1) {
      r = min(min(a_ / (a_ / l), b_ / (b_ / l)), c_ / (c_ / l));
      ret = ret * QPow(f2_2(a_ / l, b_ / l), (sumphi[r] - sumphi[l - 1] + mod - 1) % (mod - 1) * (c_ / l)) % mod;
      ret = ret * QPow(prodt_phi[r] * inv_prodt_phi[l - 1], (a_ / l) * (b_ / l) % (mod - 1) * (c_ / l)) % mod;
    }
  }
  return ret;
}
//=============================================================
int main() {
  T = read(), mod = read();
  Prepare();
  while (T -- ) {
    A = read(), B = read(), C = read();
    for (int i = 0; i <= 2; ++ i) {
      ans = f1(A, B, C, i) * f1(B, A, C, i) % mod;
      ans = ans * Inv(f2(A, B, C, i)) % mod * Inv(f2(A, C, B, i)) % mod;
      printf("%lld ", ans);  
    }
    printf("\n");
  }
  return 0;
}

---

写在最后

参考资料：
Oi-Wiki-莫比乌斯反演
算法学习笔记(35): 狄利克雷卷积 By: Pecco
题解 SP5971 【LCMSUM - LCM Sum】 - BJpers2 的博客
题解 SP5971 【LCMSUM - LCM Sum】 - Venus 的博客
题解 P3327 【[SDOI2015]约数个数和】 - suncongbo 的博客

把 Oi-Wiki 上的内容进行了复制 整理扩充。
我是个没有感情的复读机(大雾)