代码验证斯特林公式的准确性

fliter
421 阅读6分钟

关于斯特林公式


斯特林公式（Stirling's approximation或Stirling's formula）是一个用于近似计算阶乘（n!）的公式。当要为某些极大的n求阶乘时，直接计算n!的计算量会随着n的增大而快速增长，导致计算变得不实际，尤其是在计算机程序中。斯特林公式提供了一种有效的方式来近似这种大数的阶乘，能够将求解阶乘的复杂度降低到对数级。

公式如下：

[ n!≈2πn(ne)nn! \approx \sqrt{2 \pi n} \left( \frac{n}{e} \right)^n ]

其中，(e)是自然对数的底数（约等于2.71828），(π\pi)是圆周率。

斯特林是一位苏格兰数学家,生活在大概清朝康雍乾时期. 不是现在在英超踢快乐足球的那位

但实际上这个公式最早是棣莫弗发现的,就是提出概率论中棣莫弗-拉普拉斯中心极限定理的那位,法国裔英国籍的数学家

更多可参考:

大数定律与中心极限定理

揭示函数渐近展开之奥秘——棣莫弗-拉普拉斯定理的应用与实践


作用


斯特林公式在数学、物理学、工程学和计算机科学等领域有广泛的应用，主要用于：

  • 近似计算大数的阶乘。
  • 在概率论、统计学中估算组合数和排列数。
  • 分析算法的复杂度，特别是那些涉及到阶乘计算的算法。

使用Go代码验证斯特林公式的准确性

如下编写一个简单的Go程序来计算斯特林公式的近似值，并与实际的阶乘值进行比较，以此来验证斯特林公式的准确性


package main

import (
	"fmt"
	"math"
)

// 计算阶乘 n!
func factorial(n int64) int64 {
	if n == 0 {
		return 1
	}
	return n * factorial(n-1)
}

// 斯特林公式近似计算 n!
func stirlingApproximation(n float64) float64 {
	return math.Sqrt(2*math.Pi*n) * math.Pow(n/math.E， n)
}

func main() {
	var a int64 = 5 // 示例：计算5!
	exact := factorial(a)
	approx := stirlingApproximation(float64(a))

	fmt.Printf("Exact factorial of %d is: %d\n"， a， exact)
	fmt.Printf("Stirling's approximation of %d is: %f\n"， a， approx)
	fmt.Printf("Difference: %f\n"， math.Abs(float64(exact)-approx))

	fmt.Printf("误差率: %f\n"， math.Abs(float64(exact)-approx)/float64(exact))

	fmt.Println()

	var b int64 = 10 // 示例：计算10!
	exact = factorial(b)
	approx = stirlingApproximation(float64(b))

	fmt.Printf("Exact factorial of %d is: %d\n"， b， exact)
	fmt.Printf("Stirling's approximation of %d is: %f\n"， b， approx)
	fmt.Printf("Difference: %f\n"， math.Abs(float64(exact)-approx))

	fmt.Printf("误差率: %f\n"， math.Abs(float64(exact)-approx)/float64(exact))

	fmt.Println()

	var c int64 = 20 // 示例：计算20!
	exact = factorial(c)
	approx = stirlingApproximation(float64(c))

	fmt.Printf("Exact factorial of %d is: %d\n"， c， exact)
	fmt.Printf("Stirling's approximation of %d is: %f\n"， c， approx)
	fmt.Printf("Difference: %f\n"， math.Abs(float64(exact)-approx))

	fmt.Printf("误差率: %f\n"， math.Abs(float64(exact)-approx)/float64(exact))

}

在线运行

输出:

Exact factorial of 5 is: 120
Stirling's approximation of 5 is: 118.019168
Difference: 1.980832
误差率: 0.016507

Exact factorial of 10 is: 3628800
Stirling's approximation of 10 is: 3598695.618741
Difference: 30104.381259
误差率: 0.008296

Exact factorial of 20 is: 2432902008176640000
Stirling's approximation of 20 is: 2422786846761136640.000000
Difference: 10115161415503360.000000
误差率: 0.004158

上面程序中，factorial函数递归地计算了一个数的阶乘，而stirlingApproximation函数则根据斯特林公式计算了阶乘的近似值。通过比较两者的结果，可以看到斯特林公式给出的近似值与实际阶乘值之间的差异。

看起来，n越大，斯特林公式计算的结果，和实际n的阶乘值之间的误差会越小。在实际应用中，通常只用斯特林公式来近似计算大数的阶乘。

如果对于非常大的n值，直接计算阶乘可能会导致整数溢出。

2的64次方: 9223372036854775807 20的阶乘: 2432902008176640000


如果n=21，就会超出int64的范围了，需要使用bigint

如下:

package main

import (
	"fmt"
	"math"
	"math/big"
)

func factorial(n int) *big.Int {
	bigN := big.NewInt(int64(n))
	if n == 0 {
		return big.NewInt(1)
	}
	result := big.NewInt(1)
	for i := big.NewInt(2); i.Cmp(bigN) <= 0; i.Add(i， big.NewInt(1)) {
		result.Mul(result， i)
	}
	return result
}

func stirlingApproximation(n float64) float64 {
	return math.Sqrt(2*math.Pi*n) * math.Pow(n/math.E， n)
}

func main() {
	var x = 50
	exact := factorial(x)
	approx := stirlingApproximation(float64(x))

	fmt.Println("Exact factorial of 50 is: "， exact.String())
	fmt.Printf("Stirling's approximation of 50 is: %f\n"， approx)

	floatExact， _ := exact.Float64()
	fmt.Printf("Difference: %f\n"， approx-floatExact)
	fmt.Printf("误差率: %f\n"， math.Abs(floatExact-approx)/floatExact)
	fmt.Println()
	fmt.Println()

	var y = 100
	exact = factorial(y)
	approx = stirlingApproximation(float64(y))

	fmt.Println("Exact factorial of 100 is: "， exact.String())
	fmt.Printf("Stirling's approximation of 100 is: %f\n"， approx)

	floatExact， _ = exact.Float64()
	fmt.Printf("Difference: %f\n"， approx-floatExact)
	fmt.Printf("误差率: %f\n"， math.Abs(floatExact-approx)/floatExact)
	fmt.Println()
	fmt.Println()

	var z = 150
	exact = factorial(z)
	approx = stirlingApproximation(float64(z))

	fmt.Println("Exact factorial of 150 is: "， exact.String())
	fmt.Printf("Stirling's approximation of 150 is: %f\n"， approx)

	floatExact， _ = exact.Float64()
	fmt.Printf("Difference: %f\n"， approx-floatExact)
	fmt.Printf("误差率: %f\n"， math.Abs(floatExact-approx)/floatExact)

}

输出:

Exact factorial of 50 is:  30414093201713378043612608166064768844377641568960512000000000000
Stirling's approximation of 50 is: 30363445939381662539822092816663010554176972669373577771666636800.000000
Difference: -50647262331712294376073382985838127571388053403738488862408704.000000
误差率: 0.001665


Exact factorial of 100 is:  93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000
Stirling's approximation of 100 is: 93248476252694086082741480126603978924814282139093327106422842133464418571195672826024908864958629150354343629008809013568909876278442701697442603499744395264.000000
Difference: -77739191250067288427407758913009195131184458378162476302178214323666797065713990911211707211836434286969363646771402787323207418418652286983705207165157376.000000
误差率: 0.000833


Exact factorial of 150 is:  57133839564458545904789328652610540031895535786011264182548375833179829124845398393126574488675311145377107878746854204162666250198684504466355949195922066574942592095735778929325357290444962472405416790722118445437122269675520000000000000000000000000000000000000
Stirling's approximation of 150 is: 57102107404794002531486784484402246707941441176612129409152808169561298568174402178138464969562539359822454467015048872053054189489998393451341243701758138113713430727594194645443385915763391969083296804782229203527249066910085527728440667512127808337057891221504.000000
Difference: -31732159664543345709716527525897500549797926077936066874981976179631513538039613253514105887862044104270936341040554986498420072620437122053086846109932370973210446316192698228766975256370780291479993145852823806527483089444575765382819081148258251866386726912.000000
误差率: 0.000555

「阶乘秘方」斯特林公式：计算大数阶乘的神奇逼近算法

avatar
文章
阅读
粉丝