數論函數及相關基本定義html
素數的線性篩ide
線性篩能夠在嚴格$O(n)$的時間內篩出積性函數的值,函數
它有常見的套路post
假設$n = p_1^{a_1} p_2^{a_2} \dots p_k^{a_k}$url
若是咱們能快速得出$f(p_i)$和$f(p_i^{k+1})$的取值,那麼直接套板子就好了spa
在下文中如無特殊說明,默認$p_i$表示$n$質因數分解以後第$i$個質數,$a_i$表示$p_i$的指數code
常見的有如下幾種htm
比較簡單,這也是篩其餘積性函數的基礎blog
#include<cstdio> const int MAXN = 1e4 + 10; int N, prime[MAXN], vis[MAXN], tot; void GetPrime(int N) { vis[1] = 1; for(int i = 2; i <= N; i++) { if(!vis[i]) prime[++tot] = i; for(int j = 1; j <= tot && i * prime[j] <= N; j++) { vis[i * prime[j]] = 1; if(!(i % prime[j])) break; } } } int main() { GetPrime(1e3); for(int i = 1; i <= tot; i++) printf("%d ", prime[i]); return 0; }
這個也是比較常見的get
根據莫比烏斯函數的定義
$$\mu =\begin{cases}\left( -1\right) ^{k}\left( n=p_{1}p_{2}\ldots p_{k}\right) \\ 0\left( \exists P^{2}|n\right) \\ 1\left( n=1\right) \end{cases}$$
直接篩就能夠了
#include<cstdio> const int MAXN = 1e4 + 10; int N, prime[MAXN], vis[MAXN], mu[MAXN], tot; void GetMu(int N) { vis[1] = mu[1] = 1; for(int i = 2; i <= N; i++) { if(!vis[i]) prime[++tot] = i, mu[i] = -1; for(int j = 1; j <= tot && i * prime[j] <= N; j++) { vis[i * prime[j]] = 1; if(!(i % prime[j])) { mu[i * prime[j]] = 0; break; } mu[i * prime[j]] = mu[i] * mu[prime[j]]; //根據莫比烏斯函數的定義,這裏也能夠寫爲 //mu[i * prime[j]] = -mu[i]; } } } int main() { GetMu(1e3); for(int i = 1; i <= tot; i++) printf("%d ", mu[i]); return 0; }
這個貌似更經常使用一點qwq。
我在之前的文章中也詳細的講過,這裏只放一下代碼
#include<cstdio> const int MAXN = 1e4 + 10; int N, prime[MAXN], vis[MAXN], phi[MAXN], tot; void GetPhi(int N) { vis[1] = phi[1] = 1; for(int i = 2; i <= N; i++) { if(!vis[i]) prime[++tot] = i, phi[i] = i - 1; for(int j = 1; j <= tot && i * prime[j] <= N; j++) { vis[i * prime[j]] = 1; if(!(i % prime[j])) { phi[i * prime[j]] = phi[i] * prime[j]; break; } phi[i * prime[j]] = phi[i] * phi[prime[j]]; } } } int main() { GetPhi(1e3); for(int i = 1; i <= tot; i++) printf("%d ", phi[i]); return 0; }
這個就稍微高端一點了,按照上面的套路,咱們只須要考慮最小的質因子對每一個數的貢獻
$n = p_1^{a_1} p_2^{a_2} \dots p_k^{a_k}$
設$d(i)$表示$i$的約數個數
那麼根據約數定理
$d(i) = \prod_{i = 1}^k (a_i + 1) $
記$a(i)$表示$n$的最小的質因子($a_1$)的指數。
$d(p_i) = 2$
當$i % p_j = 0$時,考慮$i * p_j$,實際上也就是$a_i$的指數多了$1$
咱們先除去原來的,再加上新的就OK了
#include<cstdio> const int MAXN = 1e4 + 10; int N, prime[MAXN], vis[MAXN], D[MAXN], a[MAXN], tot; void GetD(int N) { vis[1] = D[1] = a[1] = 1; for(int i = 2; i <= N; i++) { if(!vis[i]) prime[++tot] = i, D[i] = 2, a[i] = 1; for(int j = 1; j <= tot && i * prime[j] <= N; j++) { vis[i * prime[j]] = 1; if(!(i % prime[j])) { D[i * prime[j]] = D[i] / (a[i] + 1) * (a[i] + 2); a[i * prime[j]] = a[i] + 1; break; } D[i * prime[j]] = D[i] * D[prime[j]]; a[i * prime[j]] = 1; } } } int main() { GetD(1e3); for(int i = 1; i <= tot; i++) printf("%d ", D[i]); return 0; }
一樣,按照上面的套路考慮
$n = p_1^{a_1} p_2^{a_2} \dots p_k^{a_k}$
設$SD(i)$表示$i$的約數和
$SD(n) = \prod_{i = 1}^k (\sum_{j = 1}^{a_i} p_i^j)$
$sum(i)$表示$i$最小的質因子的貢獻,即$sum(i) = \sum_{i = 1}^{a_1}p_1^j$
$low(i)$表示$i$最小質因子的指數,$low(i) = a_1$
有了這三個咱們就能夠轉移了
一樣是考慮$i$的最小的因子對答案的貢獻,應該比較好推,看代碼吧
#include<cstdio> const int MAXN = 1e4 + 10; int N, prime[MAXN], vis[MAXN], SD[MAXN], sum[MAXN], low[MAXN], tot; void GetSumD(int N) { vis[1] = SD[1] = low[1] = sum[1] = 1; for(int i = 2; i <= N; i++) { if(!vis[i]) prime[++tot] = i, sum[i] = SD[i] = i + 1, low[i] = i; for(int j = 1; j <= tot && i * prime[j] <= N; j++) { vis[i * prime[j]] = 1; if(!(i % prime[j])) { low[i * prime[j]] = low[i] * prime[j]; sum[i * prime[j]] = sum[i] + low[i * prime[j]]; SD[i * prime[j]] = SD[i] / sum[i] * sum[i * prime[j]]; break; } low[i * prime[j]] = prime[j]; sum[i * prime[j]] = prime[j] + 1; //這裏low和sum不是積性函數 SD[i * prime[j]] = SD[i] * SD[prime[j]]; } } } int main() { GetSumD(1e3); for(int i = 1; i <= tot; i++) printf("%d ", SD[i]); return 0; }
不少狀況下咱們會遇到求兩個積性函數狄利克雷卷積的狀況
很顯然,這個函數也是積性函數,咱們考慮如何求得
爲了方便篩,咱們須要把問題無限簡化,
設$low(i)$表示$p_1^{a_1}$
考慮篩法中最關鍵的地方
$i \% p_j = 0$,
有了$low(i)$,此時咱們須要分兩種狀況討論
1. $low(i) = i$,此時$i$必定是某個素數的冪的形式(不然就會break掉)
這裏就用到了我最開始說的那個套路
若是咱們能快速的利用$f(p_i^{k})$更新出$f(p_i^{k + 1})$,那這個素數就很容易篩了
2. $low(i) \not = i$,那麼$i / low(i)$必定與$low(i) * p_j$是互質的,咱們能夠直接利用積性函數的性質去更新
C++版的僞代碼
vis[1] = low[1] = 1; H[1] = 初始化 for(int i = 2; i <= N; i++) { if(!vis[i]) prime[++tot] = i, mu[i] = -1, H[i] = 質數的狀況, low[i] = i; for(int j = 1; j <= tot && i * prime[j] <= N; j++) { vis[i * prime[j]] = 1; if(!(i % prime[j])) { low[i * prime[j]] = (low[i] * prime[j]); if(low[i] == i) H[i * prime[j]] = 特殊判斷; else H[i * prime[j]] = H[i / low[i]] * H[prime[j] * low[i]]; break; } H[i * prime[j]] = H[i] * H[prime[j]]; low[i * prime[j]] = prime[j]; } }