「學(xué)習(xí)筆記」莫比烏斯反演
點(diǎn)擊查看目錄
目錄
- 「學(xué)習(xí)筆記」莫比烏斯反演
- 前置知識(shí)
- 整除分塊
- 積性函數(shù)
- 線性篩任意積性函數(shù)
- 莫比烏斯反演
- 莫比烏斯函數(shù)
- 莫比烏斯反演公式
- 例題
- [HAOI2011] Problem b
- YY 的 GCD
- [SDOI2014] 數(shù)表
- DZY Loves Math
- [SDOI2015] 約數(shù)個(gè)數(shù)和
- [SDOI2017] 數(shù)字表格
- 于神之怒加強(qiáng)版
- [國(guó)家集訓(xùn)隊(duì)] Crash的數(shù)字表格 / JZPTAB
- [湖北省隊(duì)互測(cè)2014] 一個(gè)人的數(shù)論
前置知識(shí)
整除分塊
考慮快速求：
\[\sum_{i = 1}^{n}\left\lfloor\dfrac{n}{i}\right\rfloor
\]
發(fā)現(xiàn)連續(xù)的一段 \(\left\lfloor\dfrac{n}{i}\right\rfloor\) 取值是一樣的，而且取值最多只有 \(2\sqrt{n}\) 種，考慮從這里入手，把連續(xù)的一段統(tǒng)一處理 .
當(dāng)前段左端點(diǎn) \(l\) 可以通過(guò)上一段右端點(diǎn)加一得到，如何快速求右端點(diǎn) \(r\) ？給出結(jié)論是 \(r = \left\lfloor\dfrac{n}{\left\lfloor\frac{n}{l}\right\rfloor}\right\rfloor\)，證明如下：
對(duì)于 \(\left\lfloor\frac{n}{x}\right\rfloor = k\)，有 \(n = xk + r(0\le r\le x)\)，可推導(dǎo)出不等式 \(n\ge xk\)，移項(xiàng)得 \(x\le\left\lfloor\frac{n}{k}\right\rfloor\)，此時(shí) \(x\) 最大為 \(\left\lfloor\frac{n}{k}\right\rfloor\) .
當(dāng)前塊 \(k = \left\lfloor\frac{n}{l}\right\rfloor\)，所以右端點(diǎn)（最大值）取 \(r = \left\lfloor\dfrac{n}{\left\lfloor\frac{n}{l}\right\rfloor}\right\rfloor\) .
\(O(1)\) 跑到下一個(gè)塊，最多 \(2\sqrt{n}\) 個(gè)塊，所以時(shí)間復(fù)雜度 \(O(\sqrt{n})\)，比較厲害 .
積性函數(shù)
給定數(shù)論函數(shù) \(f(x)\)，如果對(duì)于任意一組互質(zhì)的整數(shù) \(a, b\) 存在 \(f(ab) = f(a)f(b)\)，則稱 \(f(x)\) 為「積性函數(shù)」 .
特別的，如果對(duì)于任何一組（不要求互質(zhì)）整數(shù) \(a, b\) 都存在 \(f(ab) = f(a)f(b)\)，則稱 \(f(x)\) 為「完全積性函數(shù)」 .
比較常見(jiàn)的積性函數(shù)：
- \(\varphi(x)\)：歐拉函數(shù) .
- \(\sigma_{k}(x)\)：約數(shù)函數(shù)，公式為 \(\sigma_{k}(x) = \sum_{d\mid x}d^k\) . 為方便一般把 \(\sigma_{0}(x)\) 簡(jiǎn)記為 \(\tau\)，把 \(\sigma_{1}(x)\) 簡(jiǎn)記為 \(\sigma(x)\) .
- \(\mu(x)\)：莫比烏斯函數(shù)，本文核心內(nèi)容，放在下文寫(xiě) .
比較常見(jiàn)的完全積性函數(shù)：
- \(\epsilon(x) = [x = 1]\) .
- \(\operatorname{id}_k(x)=x^k\)，\(\operatorname{id}_{1}(x)\) 通常簡(jiǎn)記作 \(\operatorname{id}(x)\) .
- \(1(x) = 1\) .
好像狄利克雷卷積會(huì)用，這里掛個(gè)名 .
線性篩任意積性函數(shù)
看名字就感覺(jué)很厲害！
所有積性函數(shù)都可以被線性篩，如果只求積性函數(shù)前綴和還有「杜教篩」這種復(fù)雜度低于線性的高級(jí)篩法 . 但沒(méi)學(xué)，目前題也不用 .
注意到線性篩中所有數(shù)都只會(huì)被它的最小質(zhì)因子篩到，那么所有只有一個(gè)質(zhì)因子的數(shù)的函數(shù)值都可以基于積性函數(shù)「\(f(ab) = f(a)f(b)\)」這一性質(zhì)進(jìn)行處理 .
考慮不止一個(gè)質(zhì)因子的數(shù)怎么辦 . 設(shè)當(dāng)前篩的函數(shù)為 \(f(x)\)，最小質(zhì)因子為 \(j\)，處理到的數(shù)為 \(i\times j\)，其中 \(i = j^k\times p_1^{c_1}\times p_2^{c_2}\times\cdots\times p_n^{c_n}\) . 把 \(i\times j\) 分解為 \(\dfrac{i}{j^k}\) 和 \(j^{k + 1}\) 這兩個(gè)互質(zhì)的數(shù)即可算出 \(f(i\times j) = f\left(\dfrac{i}{j^k}\right)f\left(j^{k + 1}\right)\) . 那你存一下對(duì)于每個(gè)數(shù)的 \(k\)（最小質(zhì)因數(shù)的指數(shù)）就行了 .
如果完全積性就隨便了 .
一份篩 \(\sigma(x)\) 的實(shí)現(xiàn)，更具一般性：
d[1] = 1;
_for (i, 2, MX) {
if (!vis[i]) {
prime.push_back (i);
low[i] = i, d[i] = i + 1;
}
far (j, prime) {
if (i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) {
low[i * j] = low[i] * j;
if (low[i] == i) d[i * j] = d[i] + i * j;
else d[i * j] = d[i / low[i]] * d[low[i] * j];
break;
}
low[i * j] = j, d[i * j] = d[i] * d[j];
}
}
一份篩 \(\mu(x)\) 的實(shí)現(xiàn)，由于定義了有完全平方因子的數(shù)的函數(shù)值所以很簡(jiǎn)單（不清楚定義先往下看）：
mu[1] = 1;
_for (i, 2, MX) {
if (!vis[i]) prime.push_back (i), mu[i] = -1;
far (j, prime) {
if (i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) { mu[i * j] = 0; break; }
mu[i * j] = -mu[i];
}
}
莫比烏斯反演
莫比烏斯函數(shù)
\[\mu(x) =
\begin{cases}
1 &x = 1\\
0 &x\ 含有完全平方因子\\
(-1)^k &x\ 不含有完全平方因子，不同質(zhì)因數(shù)個(gè)數(shù)為\ k
\end{cases}
\]
性質(zhì)：
- \(\sum_{d\mid n}\mu(d) = [n = 1]\)：最重要的性質(zhì) . 常用于轉(zhuǎn)化 \([\gcd(x, y) = 1]\)，在狄利克雷卷積中也會(huì)出現(xiàn) .
- 積性函數(shù)：意味著可以被快速篩出來(lái) .
莫比烏斯反演公式
形式一：
\[f(n) = \sum_{d\mid n}g(d) \Longleftrightarrow g(n) = \sum_{d\mid n}f(d)\mu\left(\frac{n}n5n3t3z\right)
\]
證明直接大力代入，同時(shí)使用莫比烏斯函數(shù)性質(zhì)：
\[\begin{aligned}
\sum_{d\mid n}\mu(d)f\left(\frac{n}n5n3t3z\right)
&= \sum_{d\mid n}\mu(d)\sum_{k\mid \frac{n}n5n3t3z}g(k)\\
&= \sum_{k\mid n}g(k)\sum_{d\mid \frac{n}{k}}\mu(d)\\
&= \sum_{k\mid n}\left[\frac{n}{k} = 1\right]g(k)\\
&= g(n)\\
\end{aligned}
\]
\[\tag*{$\square$}
\]
形式二：
\[f(n) = \sum_{n\mid d}g(d) \Longleftrightarrow g(n) = \sum_{n\mid d}f(d)\mu\left(\fracn5n3t3z{n}\right)
\]
這個(gè)形式好像很不常用啊 .
證明：
令 \(k = \fracn5n3t3z{n}\) .
\[\begin{aligned}
\sum_{n\mid d}f(d)\mu\left(\fracn5n3t3z{n}\right)
&= \sum_{k = 1}^{+\infty}\mu(k)g(nk)\\
&= \sum_{k = 1}^{+\infty}\mu(k)\sum_{nk\mid t}f(t)\\
&= \sum_{n\mid t}f(t)\sum_{k\mid \frac{t}{n}}\mu(k)\\
&= \sum_{n\mid t}\left[\frac{t}{n} = 1\right]f(k)\\
&= g(n)\\
\end{aligned}
\]
\[\tag*{$\square$}
\]
（事實(shí)上你會(huì)發(fā)現(xiàn)，這兩個(gè)式子完全不會(huì)也能做題，因?yàn)橐话銇?lái)說(shuō)你都可以用直接推的方式代替）
例題
會(huì)把本來(lái)想把常見(jiàn)套路單開(kāi)一個(gè)部分的，想了想還是放在例題中比較好 .
都標(biāo)的十分明顯 .
[HAOI2011] Problem b
這么典的么！
可以容斥，所以只考慮上限為 \(n, m\) 怎么做 . 下文為了方便欽定 \(n < m\) .
有式子：
\[\sum_{i = 1}^{n}\sum_{j = 1}^{m}[\gcd(i, j) = k]
\]
「套路一」：變換上界：
\[\sum_{i = 1}^{\left\lfloor n/k\right\rfloor}\sum_{j = 1}^{\left\lfloor m/k\right\rfloor}[\gcd(i, j) = 1]
\]
「套路二」：用莫比烏斯函數(shù)的性質(zhì)把 \([\gcd(i, j) = 1]\) 換掉：
\[\sum_{i = 1}^{\left\lfloor n/k\right\rfloor}\sum_{j = 1}^{\left\lfloor m/k\right\rfloor}\sum_{d\mid\gcd(i, j)}\mu(d)
\]
「套路三」：變換求和順序，枚舉 \(d\)：
\[\sum_{d = 1}^{n}\mu(d)\sum_{i = 1}^{\left\lfloor n/k\right\rfloor}\sum_{j = 1}^{\left\lfloor m/k\right\rfloor}[d\mid\gcd(i, j)]
\]
化簡(jiǎn)成：
\[\sum_{d = 1}^{n}\mu(d)\left\lfloor\frac{n}{dk}\right\rfloor\left\lfloor\frac{m}{dk}\right\rfloor
\]
復(fù)雜度是 \(O(n)\) 的，考慮使用整除分塊，預(yù)處理出 \(\mu\) 函數(shù)前綴和，即可 \(O(\sqrt{n})\) 通過(guò) .
點(diǎn)擊查看代碼
const int N = 5e4 + 10, MX = 5e4;
namespace SOLVE {
int a, b, c, d, k, mu[N], smu[N]; ll ans;
bool vis[N]; std::vector <int> prime;
inline int rnt () {
int x = 0, w = 1; char c = getchar ();
while (!isdigit (c)) { if (c == '-') w = -1; c = getchar (); }
while (isdigit (c)) x = (x << 3) + (x << 1) + (c ^ 48), c = getchar ();
return x * w;
}
inline void Pre () {
smu[1] = mu[1] = 1;
_for (i, 2, MX) {
if (!vis[i]) prime.push_back (i), mu[i] = -1;
far (j, prime) {
if (i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) { mu[i * j] = 0; break; }
mu[i * j] = -mu[i];
}
smu[i] = smu[i - 1] + mu[i];
}
return;
}
inline int F (int a, int b) {
ll sum = 0;
for (int l = 1, r = 0; l * k <= std::min (a, b); l = r + 1) {
r = std::min (a / (a / l), b / (b / l));
sum += 1ll * (smu[r] - smu[l - 1]) * (a / l / k) * (b / l / k);
}
return sum;
}
inline void In () {
a = rnt (), b = rnt (), c = rnt (), d = rnt (), k = rnt ();
return;
}
inline void Solve () {
ans = F (b, d) - F (a - 1, d) - F (b, c - 1) + F (a - 1, c - 1);
return;
}
inline void Out () {
printf ("%lld\n", ans);
return;
}
}
YY 的 GCD
考慮枚舉質(zhì)數(shù)：
\[\sum_{p\in\mathbb{P}}\sum_{i = 1}^{n}\sum_{j = 1}^{m}[\gcd(i, j) = p]
\]
使用套路一：
\[\sum_{p\in\mathbb{P}}\sum_{i = 1}^{\left\lfloor n / p\right\rfloor}\sum_{j = 1}^{\left\lfloor m / p\right\rfloor}[\gcd(i, j) = 1]
\]
使用套路二：
\[\sum_{p\in\mathbb{P}}\sum_{i = 1}^{\left\lfloor n / p\right\rfloor}\sum_{j = 1}^{\left\lfloor m / p\right\rfloor}\sum_{d\mid\gcd(i, j)}\mu(d)
\]
使用套路三：
\[\sum_{p\in\mathbb{P}}\sum_{d = 1}^{n}\mu(d)\left\lfloor\frac{n}{pd}\right\rfloor\left\lfloor\frac{m}{pd}\right\rfloor
\]
「套路四」：為了方便整除分塊，枚舉分母：
令 \(T = pd\)，則：
\[\begin{aligned}
&\sum_{p\in\mathbb{P}}\sum_{d = 1}^{n}\mu(d)\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\\
=&\sum_{t = 1}^{n}\sum_{d\mid T}\mu(d)\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\\
=&\sum_{t = 1}^{n}\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\sum_{d\mid T}\mu(d)\\
\end{aligned}
\]
預(yù)處理出每個(gè)數(shù)的所有因數(shù)的 \(\mu\) 函數(shù)之和即可 . 這個(gè)題可以直接用調(diào)和級(jí)數(shù)預(yù)處理搞，比較厲害 .
點(diǎn)擊查看代碼
const int N = 1e7 + 10, MX = 1e7;
namespace SOLVE {
int n, m, mu[N], mus[N]; ll ans;
bool vis[N]; std::vector <int> prime;
inline ll rnt () {
ll x = 0, w = 1; char c = getchar ();
while (!isdigit (c)) { if (c == '-') w = -1; c = getchar (); }
while (isdigit (c)) x = (x << 3) + (x << 1) + (c ^ 48), c = getchar ();
return x * w;
}
inline void Pre () {
mu[1] = 1;
_for (i, 2, MX) {
if (!vis[i]) prime.push_back (i), mu[i] = -1;
far (j, prime) {
if (1ll * i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) { mu[i * j] = 0; break; }
mu[i * j] = -mu[i];
}
}
far (p, prime) _for (i, 1, MX / p) mus[i * p] += mu[i];
_for (i, 1, MX) mus[i] += mus[i - 1];
return;
}
inline void In () {
n = rnt (), m = rnt ();
return;
}
inline void Solve () {
ans = 0;
for (int l = 1, r = 0; l <= std::min (n, m); l = r + 1) {
r = std::min (n / (n / l), m / (m / l));
ans += 1ll * (mus[r] - mus[l - 1]) * (n / l) * (m / l);
}
return;
}
inline void Out () {
printf ("%lld\n", ans);
return;
}
}
[SDOI2014] 數(shù)表
入門題是不是過(guò)了，那不寫(xiě)什么套路幾了 .
首先我們不考慮 \(\sigma(k)\le a\) 的限制寫(xiě)出式子然后推：
\[\begin{aligned}
\sum_{k = 1}^{n}\sigma(k)\sum_{i = 1}^{n}\sum_{j = 1}^{m}[\gcd(i, j) = k]
&=\sum_{k = 1}^{n}\sigma(k)\sum_{i = 1}^{\left\lfloor n / p\right\rfloor}\sum_{j = 1}^{\left\lfloor m / p\right\rfloor}[\gcd(i, j) = 1]\\
&=\sum_{k = 1}^{n}\sigma(k)\sum_{i = 1}^{\left\lfloor n / p\right\rfloor}\sum_{j = 1}^{\left\lfloor m / p\right\rfloor}\sum_{d\mid\gcd(i, j)}\mu(d)\\
&=\sum_{k = 1}^{n}\sigma(k)\sum_{d = 1}^{n}\mu(d)\left\lfloor\frac{n}{dk}\right\rfloor\left\lfloor\frac{m}{dk}\right\rfloor\\
\end{aligned}
\]
令 \(T = dk\)：
\[\begin{aligned}
\sum_{k = 1}^{n}\sigma(k)\sum_{d = 1}^{n}\mu(d)\left\lfloor\frac{n}{dk}\right\rfloor\left\lfloor\frac{m}{dk}\right\rfloor
&=\sum_{T = 1}^{n}\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\sum_{k\mid T}\sigma(k)\mu\left(\frac{T}{k}\right)
\end{aligned}
\]
發(fā)現(xiàn)后面是個(gè)狄利克雷卷積的形式，然后發(fā)現(xiàn)卷不了，因?yàn)檫€有個(gè) \(\sigma(k)\le a\) 的限制 .
那么考慮對(duì)每個(gè)詢問(wèn)的 \(a\) 升序排序，每次詢問(wèn)時(shí)加入滿足 \(\sigma(k)\le a\) 的 \(k\)，由于查詢的是前綴和可以直接用樹(shù)狀數(shù)組維護(hù) .
模數(shù)比較特殊是 \(2^{31}\)，直接自然溢出可以跑得飛快，實(shí)測(cè)正常取模 2.08s，自然溢出 1.15s .
點(diǎn)擊查看代碼
const ll N = 1e5 + 10, MX = 1e5, Q = 2e4 + 10, P = 1ll << 31;
namespace BIT {
class BIT {
public:
int n;
private:
int b[N];
inline int lowbit (int x) { return x & -x; }
public:
inline void Update (int x, int y) {
while (x <= n) b[x] += y, x += lowbit (x);
return;
}
inline int Query (int x) {
int sum = 0;
while (x) sum += b[x], x -= lowbit (x);
return sum;
}
};
}
namespace SOLVE {
int q, mu[N], si[N], od[N], low[N], ans[Q];
std::vector <int> prime; bool vis[N]; BIT::BIT tr;
class QU {
public:
int n, m, a, id;
friend inline bool operator < (QU a, QU b) {
return a.a < b.a;
}
} qu[Q];
inline int rnt () {
int x = 0, w = 1; char c = getchar ();
while (!isdigit (c)) { if (c == '-') w = -1; c = getchar (); }
while (isdigit (c)) x = (x << 3) + (x << 1) + (c ^ 48), c = getchar ();
return x * w;
}
inline void Pre () {
tr.n = MX;
mu[1] = 1, si[1] = 1, od[1] = 1;
_for (i, 2, MX) {
if (!vis[i]) {
prime.push_back (i);
mu[i] = -1;
low[i] = i, si[i] = i + 1;
}
far (j, prime) {
if (i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) {
mu[i * j] = 0;
low[i * j] = low[i] * j;
if (low[i] == i) si[i * j] = si[i] + i * j;
else si[i * j] = si[i / low[i]] * si[low[i] * j];
break;
}
mu[i * j] = -mu[i];
low[i * j] = j, si[i * j] = si[i] * si[j];
}
od[i] = i;
}
std::sort (od + 1, od + MX + 1, [](int i, int j) { return si[i] < si[j]; });
return;
}
inline void In () {
q = rnt ();
_for (i, 1, q) qu[i].n = rnt (), qu[i].m = rnt (), qu[i].a = rnt (), qu[i].id = i;
return;
}
inline void Solve () {
std::sort (qu + 1, qu + q + 1);
int tmp = 0;
_for (i, 1, q) {
while (tmp < MX && si[od[tmp + 1]] <= qu[i].a) {
int p = od[++tmp];
_for (i, 1, MX / p) tr.Update (i * p, si[p] * mu[i]);
}
int sum = 0, a = qu[i].n, b = qu[i].m;
for (int l = 1, r = 0; l <= std::min (a, b); l = r + 1) {
r = std::min (a / (a / l), b / (b / l));
sum = (sum + (tr.Query (r) - tr.Query (l - 1)) * (a / l) * (b / l));
}
ans[qu[i].id] = sum;
}
return;
}
inline void Out () {
_for (i, 1, q) printf ("%lld\n", ans[i] < 0 ? ans[i] + P : ans[i]);
return;
}
}
DZY Loves Math
\[\begin{aligned}
\sum_{i = 1}^{n}\sum_{j = 1}^{m}f(\gcd(i, j))
&= \sum_{d = 1}^{n}f(d)\sum_{i = 1}^{n}\sum_{j = 1}^{m}[\gcd(i, j) = d]\\
&= \sum_{d = 1}^{n}f(d)\sum_{i = 1}^{\left\lfloor n/d\right\rfloor}\sum_{j = 1}^{\left\lfloor m/d\right\rfloor}\sum_{x\mid \gcd(i, j)}\mu(x)\\
&= \sum_{d = 1}^{n}f(d)\sum_{x = 1}^{n}\mu(x)\left\lfloor\frac{n}{dx}\right\rfloor\left\lfloor\frac{m}{dx}\right\rfloor\\
\end{aligned}
\]
令 \(T = dx\)：
\[\begin{aligned}
\sum_{d = 1}^{n}f(d)\sum_{x = 1}^{n}\mu(x)\left\lfloor\frac{n}{dx}\right\rfloor\left\lfloor\frac{m}{dx}\right\rfloor
&= \sum_{d = 1}^{n}f(d)\sum_{x = 1}^{n}\mu\left(\frac{T}n5n3t3z\right)\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\\
&= \sum_{T = 1}^{n}\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\sum_{d\mid T}f(d)\mu\left(\frac{T}n5n3t3z\right)\\
\end{aligned}
\]
調(diào)和級(jí)數(shù)預(yù)處理 .
點(diǎn)擊查看代碼
const ll N = 1e7 + 10, MX = 1e7;
namespace SOLVE {
ll n, m, mu[N], f[N], low[N], cnt[N], sum[N], ans;
bool vis[N]; std::vector <ll> prime;
inline ll rnt () {
ll x = 0, w = 1; char c = getchar ();
while (!isdigit (c)) { if (c == '-') w = -1; c = getchar (); }
while (isdigit (c)) x = (x << 3) + (x << 1) + (c ^ 48), c = getchar ();
return x * w;
}
inline void Pre () {
mu[1] = 1, f[1] = 0;
_for (i, 2, MX) {
if (!vis[i]) prime.push_back (i), mu[i] = -1, f[i] = 1, low[i] = i;
far (j, prime) {
if (i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) {
mu[i * j] = 0, low[i * j] = low[i] * j;
f[i * j] = std::max (f[i], f[low[i]] + 1);
break;
}
mu[i * j] = -mu[i], f[i * j] = f[i], low[i * j] = j;
}
}
_for (i, 1, MX) _for (j, 1, MX / i) sum[i * j] += mu[i] * f[j];
_for (i, 1, MX) sum[i] += sum[i - 1];
return;
}
inline void In () {
n = rnt (), m = rnt ();
return;
}
inline void Solve () {
ans = 0;
for (ll l = 1, r = 0; l <= std::min (n, m); l = r + 1) {
r = std::min (n / (n / l), m / (m / l));
ans += (sum[r] - sum[l - 1]) * (n / l) * (m / l);
}
return;
}
inline void Out () {
printf ("%lld\n", ans);
return;
}
}
[SDOI2015] 約數(shù)個(gè)數(shù)和
給出一個(gè)人類想不出來(lái)的結(jié)論：
\[d(ij) = \sum_{x\mid i}\sum_{y\mid j}[\gcd(x, y) = 1]
\]
簡(jiǎn)單解釋：假設(shè)對(duì)于 \(ij\) 的質(zhì)因數(shù) \(p\)，其在 \(i\) 中的指數(shù)為 \(a\)，在 \(j\) 中的指數(shù)為 \(b\)，則其在 \(ij\) 中的指數(shù)為 \(a + b\)，考慮如何選出來(lái)所有可能的 \(a + b + 1\) 種指數(shù) . 我們欽定要選的指數(shù) \(k\le a\) 時(shí)從 \(i\) 取，否則把 \(i\) 取光后從 \(j\) 里取剩下的 . 為了方便，「\(i\) 全選后在 \(j\) 里選 \(k - a\) 個(gè)」變?yōu)椤?span id="n5n3t3z" class="math inline">\(i\) 里不選直接在 \(j\) 里選剩下的 \(k - a\) 個(gè)」，那就能保證選出來(lái)的兩個(gè)數(shù)是互質(zhì)的，就成了上面那個(gè)式子 .
呃呃，寫(xiě)的好抽象，感性理解罷 .
然后就直接推式子：
\[\begin{aligned}
\sum_{i = 1}^{n}\sum_{j = 1}^{m}d(ij)
&= \sum_{i = 1}^{n}\sum_{j = 1}^{m}\sum_{x\mid i}\sum_{y\mid j}[\gcd(x, y) = 1]\\
&= \sum_{i = 1}^{n}\sum_{j = 1}^{m}\sum_{x\mid i}\sum_{y\mid j}\sum_{k\mid\gcd(x, y)}\mu(k)\\
&= \sum_{x = 1}^{n}\sum_{y = 1}^{m}\bigg\lfloor\frac{n}{x}\bigg\rfloor\bigg\lfloor\frac{m}{y}\bigg\rfloor\sum_{k\mid\gcd(x, y)}\mu(k)\\
&= \sum_{k = 1}^{n}\mu(k)\sum_{x = 1}^{\lfloor n/k\rfloor}\sum_{y = 1}^{\lfloor m/k\rfloor}\bigg\lfloor\frac{n}{kx}\bigg\rfloor\bigg\lfloor\frac{m}{ky}\bigg\rfloor\\
&= \sum_{k = 1}^{n}\mu(k)(\sum_{x = 1}^{\lfloor n/k\rfloor}\bigg\lfloor\frac{n}{kx}\bigg\rfloor)(\sum_{y = 1}^{\lfloor m/k\rfloor}\bigg\lfloor\frac{m}{ky}\bigg\rfloor)\\
\end{aligned}
\]
考慮預(yù)處理出所有 \(\sum_{i = 1}^{n}\left\lfloor\frac{n}{i}\right\rfloor\)，可以用整除分塊做到單次詢問(wèn) \(O(\sqrt{n})\) .
點(diǎn)擊查看代碼
const ll N = 5e4 + 10, MX = 5e4;
namespace SOLVE {
ll n, m, mu[N], mus[N], sum[N], ans;
bool vis[N]; std::vector <ll> prime;
inline ll rnt () {
ll x = 0, w = 1; char c = getchar ();
while (!isdigit (c)) { if (c == '-') w = -1; c = getchar (); }
while (isdigit (c)) x = (x << 3) + (x << 1) + (c ^ 48), c = getchar ();
return x * w;
}
inline void Pre () {
mus[1] = mu[1] = 1;
_for (i, 2, MX) {
if (!vis[i]) prime.push_back (i), mu[i] = -1;
far (j, prime) {
if (i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) { mu[i * j] = 0; break; }
mu[i * j] = -mu[i];
}
mus[i] = mus[i - 1] + mu[i];
}
_for (i, 1, MX) {
for (ll l = 1, r = 0; l <= i; l = r + 1) {
r = i / (i / l);
sum[i] += (r - l + 1) * (i / l);
}
}
return;
}
inline void In () {
n = rnt (), m = rnt ();
return;
}
inline void Solve () {
ans = 0;
for (ll l = 1, r = 0; l <= std::min (n, m); l = r + 1) {
r = std::min (n / (n / l), m / (m / l));
ans += (mus[r] - mus[l - 1]) * sum[n / l] * sum[m / l];
}
return;
}
inline void Out () {
printf ("%lld\n", ans);
return;
}
}
[SDOI2017] 數(shù)字表格
這個(gè)是乘積，交換運(yùn)算順序要化成指數(shù).
\[\begin{aligned}
\prod_{i = 1}^{n}\prod_{j = 1}^{m}f(\gcd(i, j))
&= \prod_{k = 1}^{n}f(k)^{\sum_{i = 1}^{n}\sum_{j = 1}^{m}[\gcd(i, j) = k]}\\
&= \prod_{k = 1}^{n}f(k)^{\sum_{d = 1}^{n}\mu(d)\left\lfloor\frac{n}{dk}\right\rfloor\left\lfloor\frac{m}{dk}\right\rfloor}\\
\end{aligned}
\]
令 \(T = dk\)：
\[\begin{aligned}
\prod_{k = 1}^{n}f(k)^{\sum_{d = 1}^{n}\mu(d)\left\lfloor\frac{n}{dk}\right\rfloor\left\lfloor\frac{m}{dk}\right\rfloor}
&= \prod_{T = 1}^{n}\prod_{d\mid T}f\left(\frac{T}n5n3t3z\right)^{\mu(d)\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor}\\
&= \prod_{T = 1}^{n}\left(\prod_{d\mid T}f\left(\frac{T}n5n3t3z\right)^{\mu(d)}\right)^{\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor}\\
\end{aligned}
\]
\(\prod_{d\mid T}f\left(\frac{T}n5n3t3z\right)^{\mu(d)}\) 是可以調(diào)和級(jí)數(shù)預(yù)處理出來(lái)的，剩下部分使用整除分塊解決 .
點(diǎn)擊查看代碼
const ll N = 1e6 + 10, MX = 1e6, P = 1e9 + 7;
namespace SOLVE {
ll n, m, f[N], mu[N], pro[N], inv_f[N], inv_pro[N], ans;
bool vis[N]; std::vector <ll> prime;
inline ll rnt () {
ll x = 0, w = 1; char c = getchar ();
while (!isdigit (c)) { if (c == '-') w = -1; c = getchar (); }
while (isdigit (c)) x = (x << 3) + (x << 1) + (c ^ 48), c = getchar ();
return x * w;
}
inline ll FastPow (ll a, ll b) {
ll ans = 1;
while (b) {
if (b & 1) ans = ans * a % P;
a = a * a % P, b >>= 1;
}
return ans;
}
inline void Pre () {
f[1] = 1, inv_f[1] = 1, mu[1] = 1;
_for (i, 2, MX) {
f[i] = (f[i - 1] + f[i - 2]) % P;
inv_f[i] = FastPow (f[i], P - 2);
if (!vis[i]) prime.push_back (i), mu[i] = -1;
far (j, prime) {
if (i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) { mu[i * j] = 0; break; }
mu[i * j] = -mu[i];
}
}
pro[0] = 1, inv_pro[0] = 1;
_for (i, 1, MX) pro[i] = 1;
_for (i, 1, MX) {
_for (j, 1, MX / i) {
if (mu[j] > 0) pro[i * j] = pro[i * j] * f[i] % P;
else if (mu[j] < 0) pro[i * j] = pro[i * j] * inv_f[i] % P;
}
pro[i] = pro[i] * pro[i - 1] % P;
inv_pro[i] = FastPow (pro[i], P - 2);
}
return;
}
inline void In () {
n = rnt (), m = rnt ();
return;
}
inline void Solve () {
ans = 1;
for (ll l = 1, r = 0; l <= std::min (n, m); l = r + 1) {
r = std::min (n / (n / l), m / (m / l));
ans = ans * FastPow (pro[r] * inv_pro[l - 1] % P, (n / l) * (m / l) % (P - 1)) % P;
}
return;
}
inline void Out () {
printf ("%lld\n", (ans + P) % P);
return;
}
}
于神之怒加強(qiáng)版
推導(dǎo)部分比較平凡：
\[\begin{aligned}
\sum_{i = 1}^{n}\sum_{j = 1}^{m}\gcd(i, j)^k
&= \sum_{d = 1}^{n}d^k\sum_{i = 1}^{n}\sum_{j = 1}^{m}[\gcd(i, j) = d]\\
&= \sum_{d = 1}^{n}d^k\sum_{i = 1}^{\left\lfloor n/d\right\rfloor}\sum_{j = 1}^{\left\lfloor m/d\right\rfloor}\sum_{x\mid \gcd(i, j)}\mu(x)\\
&= \sum_{d = 1}^{n}d^k\sum_{x = 1}^{n}\mu(x)\left\lfloor\frac{n}{dx}\right\rfloor\left\lfloor\frac{m}{dx}\right\rfloor\\
\end{aligned}
\]
令 \(T = dx\)：
\[\begin{aligned}
\sum_{d = 1}^{n}d^k\sum_{x = 1}^{n}\mu(x)\left\lfloor\frac{n}{dx}\right\rfloor\left\lfloor\frac{m}{dx}\right\rfloor
&= \sum_{d = 1}^{n}d^k\sum_{x = 1}^{n}\mu\left(\frac{T}n5n3t3z\right)\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\\
&= \sum_{T = 1}^{n}\left\lfloor\frac{n}{T}\right\rfloor\left\lfloor\frac{m}{T}\right\rfloor\sum_{d\mid T}d^k\mu\left(\frac{T}n5n3t3z\right)\\
\end{aligned}
\]
令 \(g(T) = \sum_{d\mid T}^{n}d^k\mu\left(\frac{T}n5n3t3z\right)\)，這部分可以埃篩處理，雖然能過(guò)但是交一發(fā)就會(huì)喜提最劣解 .
觀察發(fā)現(xiàn)是個(gè)卷積形式，而且卷的兩個(gè)函數(shù)都是積性函數(shù)，說(shuō)明 \(g\) 也是個(gè)積性函數(shù) . 考慮如何線性篩篩出來(lái) .
如果 \(T\) 是一個(gè)質(zhì)數(shù)則有 \(g(T) = T^k - 1\)，否則根據(jù)積性函數(shù)定義把它分解質(zhì)因數(shù)然后推一下：
\[\begin{aligned}
g(T)
&=\prod_{i = 1}^{t}g\left(p_i^{c_i}\right)\\
&=\prod_{i = 1}^{t}((p_i^{c_i - 1})^k\mu(p_i) + (p_i^{c_i})^k\mu(1))\\
&=\prod_{i = 1}^{t}p_i^{k(c_i - 1)}(p_i^k - 1)\\
\end{aligned}
\]
就可以線性篩預(yù)處理了 .
點(diǎn)擊查看代碼
const ll N = 5e6 + 10, MX = 5e6, P = 1e9 + 7;
namespace SOLVE {
ll k, n, m, mu[N], dk[N], sum[N], ans;
bool vis[N]; std::vector <ll> prime;
inline ll rnt () {
ll x = 0, w = 1; char c = getchar ();
while (!isdigit (c)) { if (c == '-') w = -1; c = getchar (); }
while (isdigit (c)) x = (x << 3) + (x << 1) + (c ^ 48), c = getchar ();
return x * w;
}
inline ll FastPow (ll a, ll b) {
ll ans = 1;
while (b) {
if (b & 1) ans = ans * a % P;
a = a * a % P, b >>= 1;
}
return ans;
}
inline void Pre () {
mu[1] = 1, dk[1] = 1;
_for (i, 2, MX) {
if (!vis[i]) prime.push_back (i), mu[i] = -1;
far (j, prime) {
if (i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) { mu[i * j] = 0; break; }
mu[i * j] = -mu[i];
}
dk[i] = FastPow (i, k);
}
_for (i, 1, MX) _for (j, 1, MX / i) sum[i * j] += dk[i] * mu[j] % P;
_for (i, 1, MX) sum[i] = ((sum[i - 1] + sum[i]) % P + P) % P;
return;
}
inline void In () {
n = rnt (), m = rnt ();
return;
}
inline void Solve () {
ans = 0;
for (ll l = 1, r = 0; l <= std::min (n, m); l = r + 1) {
r = std::min (n / (n / l), m / (m / l));
ans = (ans + (sum[r] - sum[l - 1] + P) * (n / l) % P * (m / l) % P) % P;
}
return;
}
inline void Out () {
printf ("%lld\n", ans);
return;
}
}
[國(guó)家集訓(xùn)隊(duì)] Crash的數(shù)字表格 / JZPTAB
硬推：
\[\begin{aligned}
\sum_{i = 1}^{n}\sum_{j = 1}^{m}\operatorname{lcm}(i, j)
&= \sum_{i = 1}^{n}\sum_{j = 1}^{m}\frac{ij}{\gcd(i, j)}\\
&= \sum_{k = 1}^{n}\frac{1}{k}\sum_{i = 1}^{n}\sum_{j = 1}^{m}ij[\gcd(i, j) = k]\\
&= \sum_{k = 1}^{n}\frac{1}{k}\sum_{i = 1}^{\left\lfloor n/k\right\rfloor}\sum_{j = 1}^{\left\lfloor m/k\right\rfloor}ij[\gcd(i, j) = 1]\\
&= \sum_{k = 1}^{n}\frac{1}{k}\sum_{i = 1}^{\left\lfloor n/k\right\rfloor}\sum_{j = 1}^{\left\lfloor m/k\right\rfloor}k^2ij[\gcd(i, j) = 1]\\
&= \sum_{k = 1}^{n}k\sum_{i = 1}^{\left\lfloor n/k\right\rfloor}\sum_{j = 1}^{\left\lfloor m/k\right\rfloor}ij\sum_{d\mid\gcd(i, j)}\mu(d)\\
&= \sum_{k = 1}^{n}k\sum_{d = 1}^{n}\mu(d)\sum_{i = 1}^{\left\lfloor n/dk\right\rfloor}\sum_{j = 1}^{\left\lfloor m/dk\right\rfloor}d^2ij\\
&= \sum_{k = 1}^{n}k\sum_{d = 1}^{n}d^2\mu(d)\left(\sum_{i = 1}^{\left\lfloor n/dk\right\rfloor}i\right)\left(\sum_{j = 1}^{\left\lfloor m/dk\right\rfloor}j\right)\\
\end{aligned}
\]
令 \(T = dk, S(n) = \dbinom{n + 1}{2}\)：
\[\begin{aligned}
\sum_{k = 1}^{n}k\sum_{d = 1}^{n}d^2\mu(d)\left(\sum_{i = 1}^{\left\lfloor n/dk\right\rfloor}i\right)\left(\sum_{j = 1}^{\left\lfloor m/dk\right\rfloor}j\right)
&= \sum_{k = 1}^{n}k\sum_{d = 1}^{n}d^2\mu(d)S\left(\left\lfloor\frac{n}{T}\right\rfloor\right)S\left(\left\lfloor\frac{m}{T}\right\rfloor\right)\\
&= \sum_{T = 1}^{n}S\left(\left\lfloor\frac{n}{T}\right\rfloor\right)S\left(\left\lfloor\frac{m}{T}\right\rfloor\right)\sum_{d\mid T}d^2\mu(d)\frac{T}n5n3t3z\\
&= \sum_{T = 1}^{n}S\left(\left\lfloor\frac{n}{T}\right\rfloor\right)S\left(\left\lfloor\frac{m}{T}\right\rfloor\right)T\sum_{d\mid T}d\mu(d)\\
\end{aligned}
\]
考慮線性篩出 \(g(T) = \sum_{d\mid T}d\mu(d)\) 后整除分塊 .
卷的這兩個(gè)都是積性函數(shù)所以 \(g(T)\) 也是積性函數(shù)，對(duì)于 \(T = ip(p\in\mathbb{P}) \) 進(jìn)行分類討論：
- \(i = 1(T\in\mathbb{P})\)：\(g(T) = 1\mu(1) + T\mu(T) = 1 - T\) .
- \(p\) 為 \(i\) 的質(zhì)因子：\(\mu(p^2) = 0\) 所以產(chǎn)生不了任何貢獻(xiàn) .
- \(p\perp i\)：由于是積性函數(shù)所以直接算 \(g(T) = g(i)g(p)\) .
點(diǎn)擊查看代碼
const int N = 1e7 + 10, MX = 1e7, P = 1e8 + 9;
namespace SOLVE {
int n, m, g[N], sum[N], s[N], ans;
bool vis[N]; std::vector <int> prime;
inline int rnt () {
int x = 0, w = 1; char c = getchar ();
while (!isdigit (c)) { if (c == '-') w = -1; c = getchar (); }
while (isdigit (c)) x = (x << 3) + (x << 1) + (c ^ 48), c = getchar ();
return x * w;
}
inline void Pre () {
g[1] = 1;
_for (i, 2, MX) {
if (!vis[i]) prime.push_back (i), g[i] = 1 - i + P;
far (j, prime) {
if (1ll * i * j > MX) break;
vis[i * j] = true;
if (!(i % j)) { g[i * j] = g[i]; break; }
g[i * j] = 1ll * g[i] * g[j] % P;
}
}
_for (i, 1, MX) {
sum[i] = (sum[i - 1] + 1ll * g[i] * i % P) % P;
s[i] = (s[i - 1] + i) % P;
}
return;
}
inline void In () {
n = rnt (), m = rnt ();
if (n > m) std::swap (n, m);
return;
}
inline void Solve () {
ans = 0;
for (int l = 1, r = 0; l <= n; l = r + 1) {
r = std::min (n / (n / l), m / (m / l));
ans = (ans + 1ll * (sum[r] - sum[l - 1] + P) * s[n / l] % P * s[m / l] % P) % P;
}
return;
}
inline void Out () {
printf ("%d\n", ans);
return;
}
}
[湖北省隊(duì)互測(cè)2014] 一個(gè)人的數(shù)論
推式子：
\[\begin{aligned}
\sum_{i = 1}^{n}i^k[\gcd(i, j) = 1]
&= \sum_{i = 1}^{n}i^k\sum_{d\mid i, d\mid n}\mu(d)\\
&= \sum_{d\mid n}\mu(d)\sum_{d\mid i}i^k\\
&= \sum_{d\mid n}\mu(d)d^k\sum_{i = 1}^{\lfloor n/d\rfloor}i^k\\
\end{aligned}
\]
發(fā)現(xiàn)后面自然數(shù)冪和 .
喜報(bào)：我不會(huì)伯努利數(shù)！
但是之前擾動(dòng)法的閑話中寫(xiě)過(guò)（原來(lái)閑話能有這種用途）：
\[S_{k}(n) = \sum_{i = 0}^{n}i^k = \frac{(n + 1) ^ {k + 1} - \sum_{j = 0}^{k - 1}\dbinom{k + 1}{j}S_j(n)}{(k + 1)} \]
還是可以看出來(lái)它是一個(gè) \(k+1\) 次多項(xiàng)式的 .
設(shè) \(f(x) = \sum_{i = 0}^{x}i^k = \sum_{i = 0}^{k + 1}a_ix^i\) .
那么繼續(xù)推式子：
\[\begin{aligned}
\sum_{d\mid n}\mu(d)d^k\sum_{i = 1}^{\lfloor n/d\rfloor}i^k
&= \sum_{d\mid n}\mu(d)d^kf\left(\frac{n}n5n3t3z\right)\\
&= \sum_{d\mid n}\mu(d)d^k\sum_{i = 0}^{k + 1}a_{i}\left(\frac{n}n5n3t3z\right)^{i}\\
&= \sum_{i = 0}^{k + 1}a_{i}n^{i}\sum_{d\mid n}\mu(d)d^{k - i}\\
\end{aligned}
\]
令 \(g(n) = \sum_{d\mid n}\mu(d)d^{k - i}\)，是個(gè)積性函數(shù)，所以 \(g(n) = \prod_{i = 1}^{w}g(p_i^{\alpha_i})\) . 考慮計(jì)算 \(g(p_i^{\alpha_i})\)，\(mu\) 只有 \(d = 1\) 和 \(d = p_i\) 時(shí)非零，易得 \(g(p_i^{\alpha_i}) = \mu(1)1^{k - i} + \mu(p_i)p_i^{k - i} = 1 - p_i^{k - i}\)，即：文章來(lái)源：http://www.zghlxwxcb.cn/news/detail-658198.html
\[\sum_{i = 0}^{k + 1}a_{i}n^{i}\prod_{i = 1}^{w}(1 - p_i^{k - i})
\]
唯一的問(wèn)題是 \(a_i\) 怎么求，數(shù)據(jù)范圍允許 \(O(d^3)\)，所以直接考慮直接高斯消元 .文章來(lái)源地址http://www.zghlxwxcb.cn/news/detail-658198.html
點(diǎn)擊查看代碼
const ll W = 1e3 + 10, N = 110, P = 1e9 + 7;
namespace SOLVE {
ll k, w, n, p[W], alpha[W], g[N][N], a[N], ans;
inline ll rnt () {
ll x = 0, w = 1; char c = getchar ();
while (!isdigit (c)) { if (c == '-') w = -1; c = getchar (); }
while (isdigit (c)) x = (x << 3) + (x << 1) + (c ^ 48), c = getchar ();
return x * w;
}
inline ll FastPow (ll a, ll b) {
ll ans = 1;
while (b) {
if (b & 1) ans = ans * a % P;
a = a * a % P, b >>= 1;
}
return ans;
}
inline ll Inv (ll a) { return FastPow (a, P - 2); }
inline void Build () {
_for (i, 0, k + 1) {
_for (j, 0, k + 1) g[i][j] = FastPow (i + 1, j);
_for (j, 1, i + 1) g[i][k + 2] = (g[i][k + 2] + FastPow (j, k)) % P;
}
return;
}
inline void Gauss () {
_for (i, 0, k + 1) {
ll l = i;
_for (j, i + 1, k + 1) if (g[j][i] > g[l][i]) l = j;
std::swap (g[i], g[l]);
ll fm = Inv (g[i][i]);
_for (j, i, k + 2) g[i][j] = g[i][j] * fm % P;
_for (j, 0, k + 1) {
if (i == j) continue;
_for (l, i + 1, k + 2) g[j][l] = (g[j][l] + P - g[j][i] * g[i][l] % P) % P;
}
}
_for (i, 0, k + 1) a[i] = g[i][k + 2];
return;
}
inline void In () {
k = rnt (), w = rnt (), n = 1;
_for (i, 1, w) {
p[i] = rnt (), alpha[i] = rnt ();
n = n * FastPow (p[i], alpha[i]) % P;
}
return;
}
inline void Solve () {
Build (), Gauss ();
ans = 0;
_for (i, 0, k + 1) {
ll prod = 1;
_for (j, 1, w) prod = prod * (1 - FastPow (p[j], k - i + P - 1) + P) % P;
ans = (ans + a[i] * FastPow (n, i) % P * prod % P) % P;
}
return;
}
inline void Out () {
printf ("%lld\n", ans);
return;
}
}
到了這里，關(guān)于「學(xué)習(xí)筆記」莫比烏斯反演的文章就介紹完了。如果您還想了解更多內(nèi)容，請(qǐng)?jiān)谟疑辖撬阉鱐OY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章，希望大家以后多多支持TOY模板網(wǎng)！
