time: 2021-11-05 21:00:22
tags:
- 筆記 聽課筆記 數論/積性函數 數論/莫比烏斯反演
數論函數
定義域為正整數，陪域為實數的函數。
積性函數
定義當 \((a,b)=1\) 時滿足 \(f(ab)=f(a)f(b)\) 的函數為積性函數。而對於任意 \(a,b\)，\(f(ab)=f(a)f(b)\) 都成立的函數叫做完全積性函數。
常見的積性函數有
- 恆等函數 \(I(n)=1\)
- 冪函數 \(I_k(n)=n^k\)
- 單位函數 \(id(n)=n\)
- 元函數 \(\varepsilon (n)=\begin{cases}1,&n=1\\0,&n>1\end{cases}\)
- 因子和函數 \(\sigma(n)=\sum_{d\mid n}d\)
- 約數個數函數 \(d(n)=\sum_{d\mid n} 1\)
【例 關於積性】
如設 \(\sigma_k(n)=\sum_{d\mid n}d^k\)，則 \(\sigma_k(n)\) 也是積性函數。考慮設互質的兩個數 \(a=\prod_{i=1}^x p_i^{q_i},b=\prod_{i=1}^yh_i^{s_i}\).
那么 \(\sigma_k(a)=\sum_{i=1}^x(1+p_i+\dots+p_i^{q_i})^k\)，\(\sigma_k(b)=\sum_{i=1}^y(1+h_i+\dots+h_i^{s_i})^k\)，所以顯然有 \(\sigma_k(ab)=\sigma_k(a)\sigma_k(b)\).
線性篩積性函數
線性篩任何積性函數都要分兩種情況考慮， 一是當前數 \(i\bmod p_j\ne0\)，二是 \(i\bmod p_j=0\). \(i\bmod p_j\ne 0\) 的情況意味着 \((i,p_j)=1\)，則 \(f(ip_j)=f(i)f(p_j)\) 即可。\(i\bmod p_j=0\) 的情況意味着 \(p_j\) 是 \(i\) 最小的質因數。設 \(i=p_j^q\times k\)，則我們可以記下 \(q\) 的值，如果能 \(O(1)\) 直接計算 \(f(p_j^{q+1})\)，或者 \(O(1)\) 從 \(f(p_j^q)\) 轉移得到 \(f(p_j^{q+1})\) 的值，就可以算出 \(f(ip_j)=f(i)/f(p_j^q)\times f(p_j^{q+1})\) 或者 \(f(ip_j)=f(i/p_j^q)\times f(p_j^{q+1})\)，不管怎么算，都需要記下對於 \(i\) 來說 \(f(p_j^q)\) 的值。
不過對於一些特殊的積性函數，有特殊的計算技巧使得我們不用考慮那么多東西。
篩 \(\varphi\)
對於 \(p\mid i\) 的情況，有 \(\varphi(pi)=pi\prod_{j=1}^q\frac{s_i-1}{s_i}=p\times \varphi(i)\).
什么額外的東西都不用記。
篩 \(\mu\)
對於 \(p\mid i\) 的情況，發現 \(p\) 的冪次大於了 1，因此 \(\mu(pi)=0\).
什么額外的東西都不用記。
篩 \(d\)（約數個數）
對於 \(p\mid i\) 的情況，則 \(d(pi)=d(i)/d(p^q)\times (p^{q+1})=d(i)/(q+1)\times(q+2)\).
需要記下最小質因子的次數 \(q\)。
篩 \(\sigma\)（約數和）
對於 \(p\mid i\) 的情況，\(\sigma(pi)=\sigma(i)/\sigma(p^q)\times \sigma (p^{q+1})=\sigma(i)/\sigma(p^q)\times (\sigma(p^q)+p^{q+1})\).
需要記下 \(q\)，\(p^q\).
代碼
什么？你說 \(\sigma\) 會爆 int？這與我有什么關系
#include <bits/stdc++.h>
using namespace std;
const int N = 1e7 + 5;
int flag[N], phi[N], mu[N], d[N], sig[N];
int q[N], pq[N], p[N], v, c;
void sieve(int n) {
flag[1] = true, phi[1] = 1, mu[1] = 1, d[1] = 1, sig[1] = 1;
for(int i = 2; i <= n; i++) {
if(!flag[i]) {
p[++c] = i, phi[i] = i - 1, mu[i] = -1, d[i] = 2, sig[i] = 1 + i;
q[i] = 1, pq[i] = i;
}
for(int j = 1; j <= c && (v = i * p[j]) <= n; j++) {
flag[v] = true;
if(i % p[j]) {
phi[v] = phi[p[j]] * phi[i];
mu[v] = -mu[i]; // mu[v] = mu[p[j]] * mu[i]
d[v] = d[p[j]] * d[i];
sig[v] = sig[p[j]] * sig[i];
q[v] = 1, pq[v] = p[j];
} else {
phi[v] = p[j] * phi[i];
mu[v] = 0;
d[v] = d[i] / (q[i] + 1) * (q[i] + 2);
sig[v] = sig[i] / sig[pq[i]] * (sig[pq[i]] + pq[i] * p[j]);
q[v] = q[i] + 1, pq[v] = pq[i] * p[j];
break;
}
}
}
}
int main() {
sieve(1e7);
for(int i = 1; i <= 20; i++) {
debug(i, phi[i], mu[i], d[i], sig[i], q[i], pq[i]);
}
return 0;
}
狄利克雷卷積
定義兩個數論函數 \(f,g\) 的狄利克雷卷積 \(f\ast g=\sum_{d\mid n}f(d)g(\frac nd)\). 狄利克雷卷積具有交換律和結合律。兩個積性函數在狄利克雷卷積后得到的函數仍然是積性函數。
\(\varphi\ast I=id\)
即 \(\sum_{d\mid n}\varphi(d)=n\). 考慮兩個集合，\(S_1=\{(i,n)|1\le i\le n\}\)，\(S_2=\{(i,p)|1\le i\le p\le n,\gcd(i,p)=1,p\mid n\}\). 將 \(S_1\) 中的元素視作分數 \(\frac in\) 后構成分數集合 \(P_1\). 這些分數約分后一定在 \(S_2\) 構成的分數集合 \(P_2\) 中，所以 \(|S_1|\le|S_2|\).
而將 \(P_2\) 中的分數 \(\frac ip\) 化為 \(\frac {ig}n\)，其中 \(g=\frac np\) 后，由於 \(\frac ip\) 的值兩兩不同，所以 \(ig\) 也兩兩不同。而 \(1\le ig\le n\)，所以 \(P_2\) 中的元素也都在 \(P_1\) 中。
綜上，\(|S_1|=|S_2|\)，即 \(n=\sum_{d\mid n}\varphi(d)\).
\(\mu\ast I=\varepsilon\)
當 \(n>1\) 時
當 \(n=1\) 時顯然有 \(\mu\ast I=1\).
\(\mu\ast id=\varphi\)
第一種證法：
第二種證法：
考慮容斥計算 \(\varphi\)，比如 \(\varphi(60)=60-30-20-12+10+6+4-2=16\). 即有 0 個公共質因數的數的個數等於至少有 0 個的，減去至少有一個的，加上至少有兩個的 ...
而右邊的部分就是 \(\mu\ast id\).
整除分塊（數論分塊）
求 \(\sum_{i=1}^n\lfloor\frac ni\rfloor\)，如果直接枚舉是 \(O(n)\) 的，使用整除分塊的方法就是 \(O(\sqrt n)\).
考慮 \(\lfloor\frac ni\rfloor\) 在 \(i\le \sqrt n\) 的時候至多有 \(\sqrt n\) 種不同的取值，而在 \(i>\sqrt n\) 時，\(\lfloor\frac ni\rfloor< \sqrt n\)，至多有 \(\sqrt n\) 種不同的取值。因此 \(\lfloor\frac ni\rfloor\) 最多有 \(2\sqrt n\) 種不同的取值。
在算法實現時，先枚舉 i=1 的情況，此時 \(\lfloor \frac ni\rfloor =n\)，我們試圖 \(O(1)\) 算出 \(\lfloor \frac ni\rfloor\) 同樣等於 \(n\) 的數的右端點。
\(\lfloor\frac ni\rfloor=k\) 實際上意味着 \(ki+p=n,0\le p<i\)，而如果 \(\lfloor\frac n{i+d}\rfloor=k\)，同樣有 \(k(i+d)+p'=n,0\le p<i+d\). 於是有 \(p'=p-kd\)，所以 \(d\le\frac pk\)，\(d_{max}=\lfloor\frac pk\rfloor\).
所以
在算出這個右端點 \(i_{max}\) 后，下一個左端點就等於 \(i_{max}+1\).
莫比烏斯反演定理
【定理 莫比烏斯反演定理 形式 1】
若
則有
【證明】
【定理 莫比烏斯反演定理 形式 2】
這個其實在做題的時候更常用。
若
則
【證明】
![[files/Pasted image 20211110172030.png]]
莫比烏斯反演入門題
P2158 [SDOI2008] 儀仗隊
作為體育委員，C 君負責這次運動會儀仗隊的訓練。儀仗隊是由學生組成的 \(N \times N\) 的方陣，為了保證隊伍在行進中整齊划一，C 君會跟在儀仗隊的左后方，根據其視線所及的學生人數來判斷隊伍是否整齊（如下圖）。
現在，C 君希望你告訴他隊伍整齊時能看到的學生人數。
顯然答案為 \(2\sum_{i=1}^{n-1}\varphi(i)+2\)，但是可以用一點高妙的方法做。
答案為
將 \([\gcd(i,j)=1]\) 轉化為 \(\sum_{d\mid \gcd(i,j)}\mu(d)\) 是個常用技巧，后面也會用到。有些題解說這個就是莫比烏斯反演？？
注意 \(n=1\) 的時候答案為 0。。。
#include <bits/stdc++.h>
using namespace std;
const int N = 40005;
int n;
int flag[N], phi[N], mu[N], d[N], sig[N];
int q[N], pq[N], p[N], v, c;
void sieve(int n) {
flag[1] = true, phi[1] = 1, mu[1] = 1, d[1] = 1, sig[1] = 1;
for(int i = 2; i <= n; i++) {
if(!flag[i]) {
p[++c] = i, phi[i] = i - 1, mu[i] = -1, d[i] = 2, sig[i] = 1 + i;
q[i] = 1, pq[i] = i;
}
for(int j = 1; j <= c && (v = i * p[j]) <= n; j++) {
flag[v] = true;
if(i % p[j]) {
phi[v] = phi[p[j]] * phi[i];
mu[v] = -mu[i]; // mu[v] = mu[p[j]] * mu[i]
d[v] = d[p[j]] * d[i];
sig[v] = sig[p[j]] * sig[i];
q[v] = 1, pq[v] = p[j];
} else {
phi[v] = p[j] * phi[i];
mu[v] = 0;
d[v] = d[i] / (q[i] + 1) * (q[i] + 2);
sig[v] = sig[i] / sig[pq[i]] * (sig[pq[i]] + pq[i] * p[j]);
q[v] = q[i] + 1, pq[v] = pq[i] * p[j];
break;
}
}
}
}
int main() {
cin >> n;
sieve(n);
long long ans = 0;
for(int i = 1; i <= n - 1; i++)
ans += mu[i] * (long long)((n - 1) / i) * ((n - 1) / i);
if(n > 1) ans += 2;
cout << ans << '\n';
return 0;
}
P1829 [國家集訓隊]Crash的數字表格 / JZPTAB
今天的數學課上，Crash 小朋友學習了最小公倍數（Least Common Multiple）。對於兩個正整數 \(a\) 和 \(b\)，\(\text{lcm}(a,b)\) 表示能同時整除 \(a\) 和 \(b\) 的最小正整數。例如，\(\text{lcm}(6, 8) = 24\)。
回到家后，Crash 還在想着課上學的東西，為了研究最小公倍數，他畫了一張 \(n \times m\) 的表格。每個格子里寫了一個數字，其中第 \(i\) 行第 \(j\) 列的那個格子里寫着數為 \(\text{lcm}(i, j)\)。
看着這個表格，Crash 想到了很多可以思考的問題。不過他最想解決的問題卻是一個十分簡單的問題：這個表格中所有數的和是多少。當 \(n\) 和 \(m\) 很大時，Crash 就束手無策了，因此他找到了聰明的你用程序幫他解決這個問題。由於最終結果可能會很大，Crash 只想知道表格里所有數的和 \(\bmod 20101009\) 的值。
設 \(n\le m\)，答案為
現在將枚舉 \(i,j\) 換為枚舉 \(d\) 的倍數，得到
\(n/d\) 可以整除分塊一下，乘上一個 \(\sqrt n\)，然后 \(n/dl\) 又可以整除分塊一下，所以總復雜度是 \(O(n)\).
代碼：
記得把減法的地方都加上 MD.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 1e7 + 5, MD = 20101009;
int n, m;
int flag[N], mu[N], f[N], S[N], p[N], c;
void init(int n) {
flag[1] = true, mu[1] = 1, f[1] = 1, S[1] = 1;
for(int i = 2, v; i <= n; i++) {
if(!flag[i]) p[++c] = i, mu[i] = -1 + MD;
for(int j = 1; j <= c && (v = i * p[j]) <= n; j++) {
flag[v] = true;
if(i % p[j]) mu[v] = (-mu[i] + MD) % MD;
else {
mu[v] = 0;
break;
}
}
f[i] = (f[i - 1] + ((ll)i * i) % MD * mu[i] % MD) % MD;
S[i] = (S[i - 1] + i) % MD;
}
}
int sum(int x, int y) {
int res = 0;
for(int i = 1, j; i <= x; i = j + 1) {
j = min(x / (x / i), y / (y / i));
res += (ll)(f[j] - f[i - 1] + MD) * S[x / i] % MD * S[y / i] % MD;
res %= MD;
}
return res;
}
int calc(int n, int m) {
int res = 0;
for(int i = 1, j; i <= n; i = j + 1) {
j = min(n / (n / i), m / (m / i));
res += (ll)(S[j] - S[i - 1] + MD) * sum(n / i, m / i) % MD;
res %= MD;
}
return res;
}
int main() {
cin >> n >> m;
if(n > m) swap(n, m);
init(m);
cout << calc(n, m) << '\n';
return 0;
}
BZOJ2694 lcm
注意這一步相當於把 \(di,dj\) 換成了 \(i,j\)，相應的在后面就要乘上 \(d^2\). 然后用上面的技巧 \([\gcd(i,j)=1]=\sum_{l\mid\gcd(i,j)}\mu(l)\) 得到
令 \(T=dl,F(i)=\frac{x(x+1)}2\)，得到
如果此時先枚舉 \(T\)，再枚舉 \(d\)，那么 \(l\) 就可以直接確定，即 \(\frac Td\). 寫出來就是
針對 \(F(A/T)F(B/T)\) 的不同取值，可以整除分塊。后面的東西是個積性函數，但是線性篩不好篩，可以枚舉 \(d\)，然后對 \(d\) 的倍數計算貢獻。復雜度是 \(n(1+\frac 12+\frac 13+\dots +\frac 1n)=O(n\log n)\).
代碼：
模數是 \(2^{30}\)，可以直接利用 unsigned int 類型的自動溢出。
#include <bits/stdc++.h>
#define int unsigned int
using namespace std;
const int N = 4e6 + 5;
int T, n, m;
int flag[N], f[N], mu[N], p[N], S[N], v, c;
void sieve1(int n) { // 線性篩 mu
mu[1] = 1;
for(int i = 2; i <= n; i++) {
if(!flag[i]) mu[i] = -1, p[++c] = i;
for(int j = 1; j <= n && (v = p[j] * i) <= n; j++) {
flag[v] = true;
if(i % p[j]) mu[v] = -mu[i];
else {
mu[v] = 0;
break;
}
}
}
}
void sieve2(int n) { // n log n 篩 f
S[1] = 1;
for(int i = 1; i <= n; i++) {
for(int j = i; j <= n; j += i) {
f[j] += i * mu[i] * mu[i] * mu[j / i] * (j / i) * (j / i);
}
S[i] = S[i - 1] + i;
}
for(int i = 1; i <= n; i++) f[i] += f[i - 1]; // 對 f 做前綴和，方便整除分塊
}
int calc(int n, int m) {
int res = 0;
for(int i = 1, j; i <= n; i = j + 1) {
j = min(n / (n / i), m / (m / i));
res += (f[j] - f[i - 1]) * S[n / i] * S[m / i];
}
return res;
}
signed main() {
sieve1(4e6);
sieve2(4e6);
for(cin >> T; T; T--) {
cin >> n >> m;
if(n > m) swap(n, m);
cout << calc(n, m) % (1 << 30) << '\n';
}
return 0;
}
P3327 [SDOI2015]約數個數和
設 \(d(x)\) 為 \(x\) 的約數個數，給定 \(n,m\)，求
\[\sum_{i=1}^n\sum_{j=1}^md(ij) \]\(T,n,m\le 50000\)
【結論】
【證明】 設 \(ij\) 的一個質因子 \(p\) 的冪次為 \(c\)，那么 \(d(ij)\) 就是所有 \(c+1\) 相乘，意義是從 \(p^0,p^1,\dots,p^c\) 里面選出來一個。
那么設 \(i\) 的質因子 \(p\) 的冪次為 \(a\)，\(j\) 的質因子的冪次為 \(b\)。把「在 \(i\) 中選擇 \(p^x\)」視作「在 \(ij\) 中選擇 \(p^x\)」，把「在 \(j\) 中選擇 \(p^y\)」視作「在 \(ij\) 中選擇 \(p^{a+x}\)」. 這樣就可以組合出 \(ij\) 所有的約數，而條件是要求 \(i\) 的約數和 \(j\) 的約數互質。
【結論】
【證明】
設 \(x=Ai+p,\ (0\le p<A)\)，則
由於 \(\frac pA<1\)，所以
原式得證。
現在開始做題了
然后這一步我直接用莫比烏斯函數 \(\sum_{d\mid \gcd(x,y)}\) 替換了 \([\gcd(x,y)=1]\)，弄出來
孩子懵逼了，雖然可以算，但是復雜度不對。
回到式 (1)，其實可以用一個技巧化簡，把其中兩個求和號搞掉。
又出現了式 (1) 中求和號底下有個整除條件的形式！於是可以確定枚舉的東西的上下界，然后把條件丟到要求得東西里去。最后改變一下枚舉的東西，就可以變得更簡單。
這個技巧一直都在用，只是這里歸納一下。
把枚舉 \(x,y\) 改成枚舉 \(dx,dy\)，得到
令 \(S(x)=\sum_{i=1}^x\lfloor\frac xi\rfloor\)，考慮「枚舉每個數的倍數，做出 1 的貢獻」就等於「所有數的約數個數和」，有
所以 \(S(x)\) 可以 \(O(n)\) 線性篩。最后答案為
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N = 5e4 + 5;
int T, n, m;
int flag[N], p[N], mu[N], d[N], q[N], S[N], Smu[N], c;
void sieve(int n) {
flag[1] = true, mu[1] = 1, d[1] = 1, S[1] = 1, Smu[1] = 1;
for(int i = 2, v; i <= n; i++) {
if(!flag[i]) p[++c] = i, mu[i] = -1, d[i] = 2, q[i] = 1;
for(int j = 1; j <= c && (v = i * p[j]) <= n; j++) {
flag[v] = true;
if(i % p[j]) {
mu[v] = -mu[i];
d[v] = d[i] * d[p[j]];
q[v] = 1;
} else {
mu[v] = 0;
d[v] = d[i] / (q[i] + 1) * (q[i] + 2);
q[v] = q[i] + 1;
break;
}
}
S[i] = S[i - 1] + d[i];
Smu[i] = Smu[i - 1] + mu[i];
}
}
int calc(int n, int m) {
int res = 0;
for(int i = 1, j; i <= n; i = j + 1) {
j = min(n / (n / i), m / (m / i));
res += (Smu[j] - Smu[i - 1]) * S[n / i] * S[m / i];
}
return res;
}
signed main() {
ios::sync_with_stdio(false); cin.tie(nullptr);
sieve(5e4);
for(cin >> T; T; T--) {
cin >> n >> m;
if(n > m) swap(n, m);
cout << calc(n, m) << '\n';
}
return 0;
}