寫程式算階乘，關於階乘的三兩事

正整數的階乘(Factorial)是小於或等於該數的所有正整數的乘積，若正整數為N，用N!來表示N的階乘。至於0!是階乘的特例，它被定義為1。階乘的運算具有遞迴(recursion)概念，常作為學習程式邏輯的材料，而這篇文章將會說明階乘相關的程式計算。

階乘(Factorial)函數

計算某數(正整數)的階乘，是指計算所有小於及等於該數的正整數的積。方程式表示如下：

#{{{
\begin{eqnarray}
Factorial(n) &=& n! \nonumber \\
&=& 1 \times 2 \times 3 \times \dotsb \times (n - 1) \times n \nonumber \\
&=& [1 \times 2 \times 3 \times \dotsb \times (n - 1)] \times n \nonumber \\
&=& n \times (n - 1)!
\end{eqnarray}
}}}#

且定義

#{{{
0! = 1
}}}#

我們可以很直覺地寫出如下的遞迴程式：

pub fn factorial(n: u8) -> u128 {
    if n <= 1 {
        1
    } else {
        n as u128 * factorial(n - 1)
    }
}

public static long factorial(final byte n) {
    if (n <= 1) {
        return 1;
    } else {
        return n * factorial((byte)(n - 1));
    }
}

function factorial(n: number) {
    if (n <= 1) {
        return 1;
    } else {
        return n * factorial(n - 1);
    }
}

func Factorial(n uint8) uint64 {
	if n <= 1 {
		return 1
	} else {
		return uint64(n) * Factorial(n-1)
	}
}

階乘的數值大小成長得很快，以下分別列出不同位元的無號整數所能儲存的最大階乘值：

32位元：12!
64位元：20!
128位元：34!
256位元：57!

非遞迴的階乘函數

以遞迴的方式呼叫函數，每次都需要建立新的堆疊框(Stack Frame)，堆疊框的數量太多就會導致堆疊溢出(Stack Overflow)，所以在實際應用上我們會比較偏好不去使用遞迴函數。上面的階乘遞迴函數，可以參考以下文章來轉成非遞迴的版本：

https://magiclen.org/recursive-to-iterative/

轉換後的非遞迴階乘函數如下：

pub fn factorial(n: u8) -> u128 {
    let mut accumulator = 1u128;

    for i in 2..=n as u128 {
        accumulator *= i;
    }

    accumulator
}

public static long factorial(final byte n) {
    long accumulator = 1;
 
    for (byte i = 2; i <= n; ++i) {
        accumulator *= i;
    }
 
    return accumulator;
}

function factorial(n: number) {
    let accumulator = 1;
 
    for (let i = 2;i <= n;i++) {
        accumulator *= i;
    }
 
    return accumulator;
}

func Factorial(n uint8) uint64 {
	accumulator := uint64(1)

	for i := uint8(2); i <= n; i++ {
		accumulator *= uint64(i)
	}

	return accumulator
}

階乘的大數運算

要實作階乘的大數運算，只要會乘法的大數運算即可。可以參考以下文章來了解加減乘除的大數運算：

https://magiclen.org/numeral-system-big-integer/

`N!`有幾位數？

前面提到階乘的數值大小成長得很快，那它到底會有幾位數呢？

要求出數字#{{N}}#有幾位數，我們可以用下面這個公式來直接算出：

#{{{
\lfloor {log_{10}N} \rfloor + 1
}}}#

所以，#{{N!}}#有幾位數，可以這樣計算：

#{{{
\begin{eqnarray}
&&\lfloor {logN!} \rfloor + 1 \nonumber \\
&=& \lfloor {log(1 \times 2 \times 3 \times \dotsb \times N)} \rfloor + 1 \nonumber \\
&=& \lfloor {log1 + log2 + log3 + \dotsb + logN} \rfloor + 1
\end{eqnarray}
}}}#

程式實作如下：

pub fn count_factorial_digits(n: u32) -> u32 {
    let mut sum = 0.0;

    for i in 2..=n {
        sum += (i as f64).log10();
    }

    sum.floor() as u32 + 1
}

public static int countFactorialDigits(final int n) {
    double sum = 0.0;
    
    for (int i = 2; i <= n; i++) {
        sum += Math.log10(i);
    }
    
    return (int)Math.floor(sum) + 1;
}

function countFactorialDigits(n: number) {
    let sum = 0.0;

    for (let i = 2;i <= n;i++) {
        sum += Math.log10(i);
    }
    
    return Math.floor(sum) + 1;
}

func CountFactorialDigits(n uint32) uint32 {
	sum := 0.0

	for i := uint32(2); i <= n; i++ {
		sum += math.Log10(float64(i))
	}

	return uint32(math.Floor(sum)) + 1
}

由於#{{log1 = 0}}#，所以計次迴圈從2開始。

上面程式的時間複雜度為#{{O(N)}}#，其實還可以更快，快到#{{O(1)}}#。

史特靈公式(Stirling's formula)，是用來計算階乘近似值的公式，如下：

#{{{
n! \approx {\sqrt{2 \pi n}}{({n \over e})}^n
}}}#

雖然史特靈公式只是近似公式，但要拿它來計算階乘的位數已經足夠準確了。位數的計算方式如下：

#{{{
\begin{eqnarray}
&&\lfloor {logN!} \rfloor + 1 \nonumber \\
&\approx& \lfloor {log({\sqrt{2 \pi N}}{({N \over e})}^N)} \rfloor + 1 \nonumber \\
&=& \lfloor { {1 \over 2}log(2 \pi N) + Nlog({N \over e}) } \rfloor + 1 \nonumber \\
&=& \lfloor { {1 \over 2} \times [log(2 \pi) + logN] + N \times [logN - log(e)] } \rfloor + 1
\end{eqnarray}
}}}#

程式實作如下：

pub fn count_factorial_digits(n: u32) -> u32 {
    let log_2_pi = (2.0 * std::f64::consts::PI).log10();

    let n_f64 = n as f64;
    let log_n = n_f64.log10();

    (0.5 * (log_2_pi + log_n) + n_f64 * (log_n - std::f64::consts::LOG10_E)).floor() as u32 + 1
}

private static final double LOG_E = Math.log10(Math.E);
private static final double LOG_2_PI = Math.log10(2.0 * Math.PI);

public static int countFactorialDigits(final int n) {
    final double logN = Math.log10(n);
    
    return (int)Math.floor(0.5 * (LOG_2_PI + logN) + n * (logN - LOG_E)) + 1;
}

const LOG_2_PI = Math.log10(2.0 * Math.PI);

function countFactorialDigits(n: number) {
    const logN = Math.log10(n);
    
    return Math.floor((0.5 * (LOG_2_PI + logN)) + (n * (logN - Math.LOG10E))) + 1;
}

func CountFactorialDigits(n uint32) uint32 {
	log2Pi := math.Log10(2.0 * math.Pi)

	nF64 := float64(n)
	logN := math.Log10(nF64)

	return uint32(math.Floor(0.5*(log2Pi+logN)+nF64*(logN-math.Log10E))) + 1
}

`N!`的尾端有幾個`0`？

階乘是連續的正整數相乘所得的結果，什麼情況下尾端會多出一個0？那就是#{{2 \times 5}}#的時候啦！N!的質因數中有幾組2和5，N!的尾端就會幾個0。由於2的數量肯定不會比5的數量少，所以我們只需要去計算N!在質因數分解後有幾個5，就是尾端的0的數量。

以下公式可以算出N!在質因數分解後有幾個5：

#{{{
\lfloor {N \over 5} \rfloor + \lfloor {N \over 25} \rfloor + \lfloor {N \over 125} \rfloor + \dotsb + \lfloor {N \over 5^i} \rfloor, 5^i \leq N
}}}#

例如，#{{4! = 24}}#，尾端有0個0；#{{5! = 120}}#，尾端有1個0；#{{9! = 362880}}#，尾端有1個0；#{{10! = 3628800}}#，尾端有2個0；#{{24! = 620448401733239439360000}}#，尾端有4個0；#{{25! = 15511210043330985984000000}}#，尾端有6個0。

在用程式實作之前，可以將上面的式子寫成以下這樣，會更容易實作：

#{{{
\lfloor {N \over 5} \rfloor + \lfloor {{N \over 5} \over 5} \rfloor + \lfloor {{{N \over 5} \over 5} \over 5} \rfloor + \dotsb
}}}#

當最外層的分子小於5時，就可以不用再加下去了，因為只會取整數部份來相加。

程式實作如下：

pub fn count_factorial_tail_zeros(mut n: u32) -> u32 {
    let mut sum = 0;

    while n >= 5 {
        n /= 5;
        sum += n;
    }

    sum
}

public static int countFactorialTailZeros(int n) {
    int sum = 0;
    
    while (n >= 5) {
        n /= 5;
        sum += n;
    }
    
    return sum;
}

function countFactorialTailZeros(n: number) {
    let sum = 0;

    while (n >= 5) {
        n = Math.floor(n / 5);
        sum += n;
    }
    
    return sum;
}

func CountFactorialTailZeros(n uint32) uint32 {
	sum := uint32(0)

	for n >= 5 {
		n /= 5
		sum += n
	}

	return sum
}

寫程式算階乘，關於階乘的三兩事

階乘(Factorial)函數

非遞迴的階乘函數

階乘的大數運算

N!有幾位數？

N!的尾端有幾個0？

`N!`有幾位數？

`N!`的尾端有幾個`0`？