Number Theoretic Function

Number Theoretic Function

「數論函數」。數論領域的經典函數。諸如因數計數函數σ₀(n)、因數總和函數σ₁(n)、質因數計數函數ω(n)、質數計數函數π(n)、互質計數函數φ(n)、質因數取捨函數μ(n)、……。維基百科整理了一份詳細列表:

https://en.wikipedia.org/wiki/Arithmetic_function

一般來說,函數輸入是正整數,函數輸出是整數。數論函數可以視作一串整數數列。

Additive Function / Multiplicative Function

數論函數通常是「加性函數」或「乘性函數」,形成遞迴公式。

additive function:       f(a + b) = f(a) × f(b)   when a ⊥ b
multiplicative function: f(a × b) = f(a) × f(b)   when a ⊥ b
加性函數:輸入相加等同輸出相乘,當輸入互質。
乘性函數:輸入相乘等同輸出相乘,當輸入互質。
completely additive function:       f(a + b) = f(a) × f(b)
completely multiplicative function: f(a × b) = f(a) × f(b)
完全加性函數:輸入相加等同輸出相乘。
完全乘性函數:輸入相乘等同輸出相乘。

演算法是篩法,或者說是Dynamic Programming。

ICPC 7837

Scalar Multiplication Function / Exponential Function

「倍率函數」與「指數函數」。

scalar multiplication function: f(a + b) = f(a) + f(b)  ⇔  f(x) = cx
exponential function:           f(a + b) = f(a) × f(b)  ⇔  f(x) = cˣ
倍率函數:輸入相加等同輸出相加。可以改寫成倍率運算。(整數系統,倍率等同乘法。)
指數函數:輸入相加等同輸出相乘。可以改寫成指數運算。

演算法是Divide-and-Conquer Method。

Linear Function

「加性函數」與「倍性函數」同時成立,稱作「線性函數」。

由於是不同領域,名稱撞名了,名稱意義不一樣。

additive function:               f(a + b) = f(a) + f(b)
homogeneous function of order 1: f(k ⋅ x) = k ⋅ f(x)
加性函數:輸入相加等同輸出相加。
倍性函數:輸入倍率等同輸出倍率。(整數系統,倍率等同乘法。)

演算法是矩陣求解。

Divisor Counting Function

List of Divisors

n的因數列表。

Divisor Counting Function σ₀(n)

n的因數數量。

Divisor Sum Function σ₁(n)

n的因數總和。

Divisor Function σₖ(n)

n的因數k次方和。k=0得到因數數量,k=1得到因數總和。

性質

σ₀(p) = 2                         when p is prime
σ₀(pⁿ) = n + 1                    when p is prime
σ₀(a × b) = σ₀(a) × σ₀(b)         when a and b are coprime
σ₀(n) = σ₀(p₁n₁) × σ₀(p₂n₂) × …   when n = p₁n₁ × p₂n₂ × … is prime factorization

希臘字母σ,唸做/ˈ sɪɡmə/,可以寫做sigma。

計算一項

試除法,時間複雜度O(sqrt(n))。

數學公式,時間複雜度仍是O(sqrt(n)),但是執行時間更短。

平方根改成立方根,時間複雜度仍是O(sqrt(n)),而且執行時間更久。弄巧成拙。負面教材。

主迴圈O(n)。最後的質數判定O(sqrt(n)) = O(n½)。

預先建立質數表,降為O(π(sqrt(n)))。π(n)是1到n的質數數量,大約是O(n/log(n))。

預先質因數分解,降為O(ω(n))。ω(n)是n的相異質因數數量,大約是O(log(n)/loglog(n))。

計算每一項

篩法。時間複雜度O(NlogN)。

GCD Distribution Function

List of GCDs

n的最大公因數列表。1到n的正整數當中,每個數與n的最大公因數。

GCD Distribution Function

n的最大公因數列表當中,某種最大公因數的數量。

GCD Sum Function

n的最大公因數列表當中,每個最大公因數的總和。

計算一項

同除以d,問題變成Coprime Counting Function。

UVa 12425

計算每一項

Dynamic Programming。時間複雜度O(NlogN)。

Coprime Counting Function

Coprime Counting Function φ(n)
(Euler's Totient Function)

1到n的正整數當中,跟n互質(最大公因數是一)的數,總共有幾個。

性質

φ(p) = p - 1                   when p is prime
φ(pⁿ) = pⁿ - pⁿ⁻¹              when p is prime
φ(a × b) = φ(a) × φ(b)         when a and b are coprime
φ(n) = φ(p₁n₁) × φ(p₂n₂) × …   when n = p₁n₁ × p₂n₂ × … is prime factorization
φ(a × b) = φ(a) × φ(b) × gcd(a,b) ÷ φ(gcd(a,b))

希臘字母φ,唸做/faɪ/,可以寫做phi。

計算一項

數學公式。不用一個一個去計算最大公因數,非常有效率!

首先是質因數分解:

n = p₁n₁ × p₂n₂ × ... × pₖnₖ   where p₁ ... pₖ are primes

以質因數計算Coprime Counting Function:

φ(n) = n × (1 - 1/p₁) × (1 - 1/p₂) × ... × (1 - 1/pₖ)

採用這個順序計算,避免溢位:

n ÷ p₁ × (p₁-1) ÷ p₂ × (p₂-1) ÷ ... ÷ pₖ × (pₖ-1)

試除法,時間複雜度O(sqrt(n))。

預先建立質數表,降為O(π(sqrt(n)))。π(n)是1到n的質數數量,大約是O(n/log(n))。

預先質因數分解,降為O(ω(n))。ω(n)是n的相異質因數數量,大約是O(log(n)/loglog(n))。

UVa 10299 10179 11064

計算每一項

篩法。時間複雜度O(NloglogN)。

線性篩法。時間複雜度O(N)。

首先運用線性篩法,求得每個數的最小質因數。然後運用倍性函數,實施Dynamic Programming。時間複雜度O(N)。

UVa 10820 10990 11327 11424 11426 12493

Möbius Function

Möbius Function μ(n)

n進行質因數分解。質因數重複為0、不重複為±1。不重複的情況下,偶數個質因數為+1、奇數個質因數為-1。

性質

μ(p¹) = -1                     when p is prime
μ(p²) = μ(p³) = ... = 0        when p is prime
μ(a × b) = μ(a) × μ(b)         when a and b are coprime
μ(n) = μ(p₁n₁) × μ(p₂n₂) × …   when n = p₁n₁ × p₂n₂ × … is prime factorization

希臘字母μ,唸做/mju:/,可以寫做mu。

計算一項

試除法,時間複雜度O(sqrt(n))。

預先建立質數表,降為O(π(sqrt(n)))。π(n)是1到n的質數數量,大約是O(n/log(n))。

預先質因數分解,降為O(ω(n))。ω(n)是n的相異質因數數量,大約是O(log(n)/loglog(n))。

計算每一項

線性篩法。時間複雜度O(N)。

首先運用線性篩法,求得每個數的最小質因數。然後運用倍性函數,實施Dynamic Programming。時間複雜度O(N)。

ICPC 2116 luogu P4318

Prime Counting Function

Prime Counting Function π(n)

1到n的正整數當中,總共有幾個質數。

計算一項

Meissel–Lehmer Algorithm。可能是O(n¾ / log(n))。

https://www.luogu.com.cn/blog/Leasier/Meissel-Lehmer

Lagarias–Miller–Odlyzko Algorithm。可能是O(n / log²(n))。

https://caojiangxia.github.io/Meissel-Lehmer/

luogu P3912 P7884

計算每一項

線性篩法。時間複雜度O(N)。

Summatory Function

Summatory Function

「累和函數」。累積和、前綴和。

數論函數的累和,通常擁有快速演算法。

Power Sum Function Sₖ(n)

「冪和函數」。正整數1到n,各自k次方的總和。

計算一項:時間複雜度O(n)。

計算一項:時間複雜度O(1)。

List of Quotients【尚無正式名稱】

「商數列表」。quo(n,1) quo(n,2) ... quo(n,n)。

正整數n,分別除以1到n,所得到的商數們。

n | list of quotients | formula
--| ------------------| --------------------------------------------
1 | 1                 | quo(1,1)
2 | 2 1               | quo(2,1) quo(2,2)
3 | 3 1 1             | quo(3,1) quo(3,2) quo(3,3)
4 | 4 2 1 1           | quo(4,1) quo(4,2) quo(4,3) quo(4,4)
5 | 5 2 1 1 1         | quo(5,1) quo(5,2) quo(5,3) quo(5,4) quo(5,5)

n的商數列表最多2sqrt(n)種數值。

分成兩種情況:
一、除數小於等於sqrt(n):除數總共sqrt(n)種,於是商數最多sqrt(n)種。
二、除數大於sqrt(n):商數小於sqrt(n),於是商數最多sqrt(n)種。

求得每種數值所在位置,時間複雜度O(sqrt(n))。

中文網路俗稱此演算法為「整除分块」。

Quotient Sum Function Q(n)【OEIS A006218

「商數總和函數」。商數列表總和。

計算一項:時間複雜度O(n)。

計算一項:援引方才的演算法,時間複雜度O(sqrt(n))。

計算一項:時間複雜度O(sqrt(n))。

Q(n) = (floor(n/1) + ... + floor(n/r)) × 2 - r²
where r = floor(sqrt(n))

計算一項:凸函數曲線下方格點,時間複雜度可能是O(n)。

List of Remainders【尚無正式名稱】

「餘數列表」。rem(n,1) rem(n,2) ... rem(n,n)。

正整數n,分別除以1到n,所得到的餘數們。

n | list of remainders | formula
--| -------------------| --------------------------------------------
1 | 0                  | rem(1,1)
2 | 0 0                | rem(2,1) rem(2,2)
3 | 0 1 0              | rem(3,1) rem(3,2) rem(3,3)
4 | 0 0 1 0            | rem(4,1) rem(4,2) rem(4,3) rem(4,4)
5 | 0 1 2 1 0          | rem(5,1) rem(5,2) rem(5,3) rem(5,4) rem(5,5)

Remainder Sum Function R(n)【OEIS A004125

「餘數總和函數」。餘數列表總和。

計算一項:時間複雜度O(n)。

餘數不易計算rem(n,1) rem(n,2) ... rem(n,n),整除部分容易計算n-rem(n,1) n-rem(n,2) ... n-rem(n,n)。此即商數列表的變種1⋅quo(n,1) 2⋅quo(n,2) ... n⋅quo(n,n)。其總和扣掉n²再變號,就是正確答案。

計算一項:援引方才的演算法,時間複雜度O(sqrt(n))。

計算一項:凸函數曲線下方格點,時間複雜度可能是O(n)。

Summatory Divisor Counting Function D₀(n)

計算一項:時間複雜度O(n sqrt(n))。

商數總和函數、因數計數函數累積和,恰好相等。

證明方式:建立因數列表、套用篩法。

計算一項:時間複雜度O(sqrt(n))。

因數計數函數,得以快速計算一項。

luogu P2424 P6028

Summatory Divisor Sum Function D₁(n)

計算一項:時間複雜度O(n sqrt(n))。

餘數總和函數、因數總和函數累積和,恰好互補,相加是n²。

計算一項:時間複雜度O(sqrt(n))。

因數總和函數,得以快速計算一項。

Coprime Counting Function

f(n)   =   ∑   [gcd(i,n) = 1]   Coprime Counting Function
         1≤i≤n                  Euler's Totient Function

f(n)   =   ∑   [gcd(i,j) = 1]   Coprime Triangular Counting Function
        1≤i<j≤n

f(n,m) =   ∑   [gcd(i,j) = 1]   Coprime Double Counting Function
         1≤i≤n
         1≤j≤m

預處理O(n)、計算一項O(1)。只能處理n = m的情況。

預處理O(n)、計算一項O(sqrt(n))。

沒有預處理、計算一項O(n log(n))。

luogu P2158

GCD Distribution Function

f(d;n)   =   ∑   [gcd(i,n) = d]   GCD Distribution Function
           1≤i≤n

f(d;n)   =   ∑   [gcd(i,j) = d]   GCD Triangular Distribution Function
          1≤i<j≤n

f(d;n,m) =   ∑   [gcd(i,j) = d]   GCD Double Distribution Function
           1≤i≤n
           1≤j≤m

luogu P4450 P2522 P2568 P2257 P3172

GCD Sum Function

F(n)   =   ∑   gcd(i,n)   GCD Sum Function
         1≤i≤n            Pillai's Arithmetical Function

F(n)   =   ∑   gcd(i,j)   GCD Triangular Sum Function
        1≤i<j≤n

F(n,m) =   ∑   gcd(i,j)   GCD Double Sum Function
         1≤i≤n
         1≤j≤m

拆開成乘性卷積。然後整除分塊。

F(n) =  ∑  gcd(i,n) = ∑ d φ(n/d) = ∑ (n/d) φ(d)
      1≤i≤n          d|n          d|n

預處理O(n)、計算一項O(sqrt(n))。

沒有預處理、計算一項O(n log(n))。

UVa 11417 11424 11426 luogu P1390 P2398 P3768

Summatory Multiplicative Convolution

乘性卷積的累和,擁有特殊演算法,快速計算一項。

sum (a∗b)(i) = sum sum a(d)b(i/d)
i≤n            i≤n d|i

sum (a∗b)(i) = sum  sum a(d)b(i) = sum a(d) sum b(i)
i≤n            d≤n i≤n/d           d≤n     i≤n/d

sum (a∗b)(i) = sum a(d)B(n/d)
i≤n            d≤n
求得乘性卷積的累和sum (a∗b),已知A與B。

按照公式計算a∗b累和,時間複雜度O(n sqrt(n))。如果a累和與b累和,可以快速計算一項,需時O(f(n))與O(g(n)),那麼運用整除分塊降至O((f(n) + g(n)) sqrt(n))。

一、因數成雙成對,窮舉一半足矣。O(n sqrt(n))。

二、運用整除分塊,分成sqrt(n)塊,每段計算一項。O((f(n) + g(n)) sqrt(n))。

Summatory Multiplicative Function

https://codeforces.com/blog/entry/54150

乘性函數的累和,擁有特殊演算法,快速計算一項。

sum (a∗b)(i) = sum a(d)B(n/d)
i≤n            d≤n

sum (a∗b)(i) = a(1)B(n/1) + sum a(d)B(n/d)
i≤n                        2≤d≤n

B(n) = (sum (a∗b)(i) - sum a(d)B(n/d)) / a(1)
        i≤n           2≤d≤n
已知sum (a∗b)與a,已知b恰是乘性函數,求得該乘性函數的累和B。

運用線性篩法計算b累和,時間複雜度O(n)。如果a∗b累和與a,可以快速計算一項,需時O(1),那麼運用分部積分降至O(n¾)。而b是乘性函數,運用線性篩法降至O(n)。

一、當a∗b累和與a,可以快速計算一項,那麼運用分部積分,拆出第一項,形成遞迴公式。快速計算一項,假定是O(1),方便計算遞迴公式的時間複雜度。O(n¾)。

二、運用線性篩法,預先計算B的前m項。O(m + n/sqrt(m))。令m = n/sqrt(m)讓時間複雜度達到最低。O(n)。

中文網路俗稱「杜教筛」。

Summatory Coprime Counting Function

回到正題。φ的累和,有兩種演算法:

一、乘性卷積的累和(正向):採用公式φ = μ ∗ ι。

sum φ(i) = sum (μ∗ι)(i) = sum μ(d)I(n/d)
i≤n        i≤n            d≤n

二、乘性函數的累和(反向):採用公式ι = 𝟏 ∗ φ。

Φ(n) = (sum ι(i) - sum 𝟏(d)Φ(n/d)) / 𝟏(1)
        i≤n       2≤d≤n

反向的優點是不必另外計算μ。

計算一項:時間複雜度O(n)。

luogu P4213

Summatory Möbius Function

反向。公式ε = 𝟏 ∗ μ。

計算一項:時間複雜度O(n)。

Summatory Prime Counting Function

https://min-25.hatenablog.com/entry/2018/11/11/172216

https://www.cnblogs.com/zjp-shadow/p/11093242.html

luogu P5493

Ramanujan's Sum(Under Construction!)

Ramanujan's Sum

[totient multiplicative]
φ(mn) = φ(m)φ(n)
與m互質的數a  與n互質的數b 與mn互質的數x
{ x = a (mod m)     m個m個排整齊剛好是一堆column
{ x = b (mod n)     n個n個排整齊剛好是一堆column
窮舉a與b實施中國餘數定理得到x

[totient divisor sum]
sum φ(d) = n
d|n
1/20 2/20 .... 20/20
約分之後,分母只可能是因數
分母是20的有φ(20)個
分母是10的有φ(10)個
分母是5的有φ(5)個
分母是4 2 1 有 φ(4) φ(2) φ(1)個

[totient gcd]
 sum gcd(i,n) = sum φ(n) * n/d  原理同[euler divisor sum]  n/d是gcd的值
i=1~n           d|n
gcd數列累積和   φ(n)數列只取因數項

[totient 取捨]
φ(n) = sum μ(d) * n/d = n sum μ(d)/d
       d|n                d|n
μ(d) = n質因數分解  質因數各一個相乘  零個偶數個相乘+1 奇數個相乘-1  其他0
順帶一提,隨便一個poset都可以定義取捨,μ(n)是定在因數關係的poset。

[Dirichlet convolution]
sum f(d) * g(n/d)
d|n

[數列變成多項式]
複數系統(z-transform)  http://en.wikipedia.org/wiki/Dirichlet_series
餘數系統(數論轉換row)  http://en.wikipedia.org/wiki/Bell_series

[滴哩卷積的計算]
令c(n)=1,補到原本的等式裡面,弄成卷積的樣子
sum μ(d) * c(n/d) = { 1 when n = 1    [mobius divisor]
d|n                 { 0 otherwise
sum φ(d) * c(n/d) = n               [totient divisor]
d|n
sum μ(d) * n/d = φ(n)              [totient 取捨]
d|n
https://math.stackexchange.com/questions/135351/

μ(n)和c(n)=1的卷積是脈衝函數。μ(n)的倒數是水平線c(n)=1。非常奇葩。
φ(n)和c(n)=1的卷積是45度斜線id(n)=n。也很奇葩。
μ(n)和id(n)=n的卷積是φ(n)。

[滴哩卷積的計算]
  id   id
μ -> φ -> gcd累積和

設函數id(n)=n
 sum gcd(i,n) = φ(n) * id(n)  [euler gcd]
i=1~n
 sum f(gcd(i,n)) = φ(n) * f(n)
i=1~n

要如何簡單表示id的倒數?這樣就可以很快用gcd累積和來算出φ(n)
http://math.stackexchange.com/questions/423317/
http://oeis.org/A023900

要如何簡單表示φ(n)的倒數?
https://math.stackexchange.com/questions/3073401/

[拉馬努金和]
http://en.wikipedia.org/wiki/Ramanujan%27s_sum
                 nk
eta_q(n) = sum  w      n個複數波,n是倍率,q是取樣數,k是索引值,通通加起來
          k=1~q  q      首先定參一下取樣數是多少

         = { 0  if q!|n 宛如餘數加法循環,各數(各等分角)剛好出現一遍,加起來是0
           { q  if q|n  取樣通通抽中0度角=1。

                nk
c_q(n) = sum   w        針對每一種波,取樣時取不可約分(互質)的那些點,φ(q)個
       (k,q)=1  q

         if q is prime  [totient φ(q)=q-1 if q is prime 的推廣]
       = { -1   if q!|n 全部的取樣和 - 可以約分的取樣和 = 不可約分的取樣和
         { q-1  if q|n      0或q         1(只有第n點)         -1或q-1

eta_q(n) = sum c_d(n)   [totient divisor sum]的推廣,從計數推廣成複數單位圓
           d|q          (需要仔細研究一下怎麼分配位置,因為循環所以可以位移吧)
sum eta_d(n) * μ(q/d) = c_q(n)  卷積倒數
d|q
sum d        * μ(q/d) = c_q(n)  照eta的公式,等於0的就刪掉,剩下有因數的
d|gcd(q,n)
sum μ(q/d)   * d      = c_q(n)  調整一下寫法
sum μ(d)     * q/d    = c_q(n)  調整一下寫法
                 [1;31m是因數的,滴哩卷積/取捨後,得到互質的 [m

φ(q) = c_q(0)           平坦波,取樣必定都是0度角=1。照定義來,取小於n的互質數
μ(q) = c_q(1)           一倍波。這很奇葩,需要想一下。
c_q(0) = sum c_d(1) * q/d   [totient 取捨]
         d|q
c_q(n) = sum c_d(1) * q/d   前面推導的式子 (n=0是特例,剛好是totient取捨)
        d|gcd(q,n)          這可以DP:比較小的取樣結果,遞推比較大的取樣結果

[傅立葉轉換]
滴哩卷積下,μ(n)和φ(d)和g(n)=1和g(n)=n經過某種傅立葉轉換之後不知道是啥鬼。

循環卷積下,gcd數列dft的第一項(不是第零項)是φ(n)
      gcd數列與一倍頻率複數波的內積是φ(n)
不知道怎麼證

[傅立葉轉換,針對特例gcd]
http://www.kurims.kyoto-u.ac.jp/EMIS/journals/INTEGERS/papers/i50/i50.pdf
注意到沒有第零項
                     ij
gcd(i,n) = 1/n sum (w  * sum c_d(j) * n/d )  n是數量,i是索引值
              j=1~n      d|n

                 ij
         = sum (w  * sum c_d(j)/d )  外面1/n跟裡面n相消
          j=1~n      d|n

注意到j=1的那一項,sum c_d(1) * n/d = c_d(0) = φ(n)
                   d|n

還可以算複數系統的最大公因數???

List of GCDs的傅立葉轉換是Coprime Counting Function

https://zhuanlan.zhihu.com/p/166530236

《The Fourier transform of functions of the greatest common divisor》

       N-1
φ(n) =  ∑ { gcd(n, k) ÷ e𝑖⋅2π⋅(k/N)⋅n }
       k=0

Fourier Ramanujan Transform

http://perso.ens-lyon.fr/patrick.flandrin/NewOrl2.pdf

基底:複數波改成Ramanujan's Sum。範數:N改成φ(n)。正交:點積為零改成互質。逆轉換需要用到無限長數列。

Fourier Transform
       N-1
X(n) =  ∑ { x(k) ÷ e𝑖⋅2π⋅(k/N)⋅n }
       k=0
Fourier Ramanujan Transform
       N          q
X(n) =  ∑ { x(q) ⋅ ∑ e𝑖⋅2π⋅(k/q)⋅n }
       q=1        k=1,k⊥q