logの和の近似

\log{n!}の近似によくスターリングの公式が用いられます。
\log{n!} = n\log{n}-n
これで十分なことも多いですが、もう少しよい近似を考えてみましょう。

積分で近似して、
 \displaystyle S_1 = \int_2^{n+1}{\log{(x-\frac{1}{2}})dx}

この差を取ると、
 \displaystyle S - S_1 = \sum_{k=2}^n{(\log{k} - \int_k^{k+1}{\log{(x-\frac{1}{2}})dx})}

summationの中は、
 \log{k} - (k + \frac{1}{2})\log{(k + \frac{1}{2})} + (k - \frac{1}{2})\log{(k - \frac{1}{2})} + 1

 \log{(x+a)}テイラー展開すると、
 \log{(x+a)} = \log{x} + ax^{-1} - \frac{1}{2}a^2x^{-2} + \frac{1}{3}a^3x^{-3} + \cdots

 (x+a)\log{(x+a)}
 = (x+a)\log{x} + (x+a)(ax^{-1} - \frac{1}{2}a^2x^{-2} + \frac{1}{3}a^3x^{-3} + \cdots)
 = (x+a)\log{x} + a + a^2x^{-1} - \frac{1}{2}(a^2x^{-1}+a^3x^{-2}) + \frac{1}{3}(a^3x^{-2}+a^4x^{-3}) + \cdots
 = (x+a)\log{x} + a + \frac{1}{1\cdot 2}a^2x^{-1} - \frac{1}{2\cdot 3}a^3x^{-2} + \cdots

 (x+a)\log{(x+a)} - (x-a)\log{(x-a)}
 = 2a\log{x} + 2a - \frac{2}{2\cdot 3}a^3x^{-2} + \frac{2}{4\cdot 5}a^5x^{-4} + \cdots

だから、summationの中は、
 \frac{1}{2\cdot 3\cdot 2^2}k^{-2} + \frac{1}{4\cdot 5\cdot 2^4}k^{-4} + \cdots

 \displaystyle S - S_1 = \frac{1}{2\cdot 3\cdot 2^2}\sum_{k=2}^n{\frac{1}{k^2}} + \frac{1}{4\cdot 5\cdot 2^4}\sum_{k=2}^n{\frac{1}{k^4}} + \cdots

右辺第一項はリーマン・ゼータ関数

 \zeta(s) = \displaystyle \sum_{n=1}^{\infty}{\frac{1}{n^s}}

の特殊値

 \zeta(2) = \frac{\pi^2}{6}

が使えて、

 \displaystyle \sum_{k=2}^n{\frac{1}{k^2}} = \zeta(2) - 1 - \sum_{k=n+1}^{\infty}{\frac{1}{k^2}}

右辺第3項は上と同様な積分による近似が使えて、

 \displaystyle \sum_{k=n+1}^{\infty}{\frac{1}{k^2}} = \frac{1}{n+\frac{1}{2}} + O(n^{-3})

 S-S_1の右辺の第2項以降は、数値計算して、

from itertools import count

eps = 1e-22
v = []
for m in count(4, 2):
	a = 1.0 / (m*(m+1)*2**m)
	if a / 2**m < eps:
		break
	for k in count(2):
		b = a / k**m
		if b < eps:
			break
		v.append(b)

print(sum(sorted(v)))
0.000263942581577

 0.000263942581577 + O(n^{-3})

ということが分かりました。
以上をまとめると、

 S = (n+\frac{1}{2})\log{(n+\frac{1}{2})} - n + 1 - \frac{3}{2}\log{\frac{3}{2}} + \pi^2/144 + \frac{1}{24n+12} + 0.000263942581577 + O(n^{-3})

さらに定数項をまとめると、

 S = (n+\frac{1}{2})\log{(n+\frac{1}{2})} - n + \frac{1}{24n+12} + 0.46060519987133985 + O(n^{-3})

検証してみましょう。

from math import log

for n in (10, 20, 30):
	S = sum(log(m) for m in range(1, n+1))
	S2 = (n+0.5)*log(n+0.5)-n-1.0/(24*n+12)+0.41893853320467317
	print(n, S-S2)
>||

>||
10 2.0936225819667698e-06
20 2.8191458056880947e-07
30 8.56364152923561e-08

およそ-3乗で誤差が減っていっていることが分かると思います。

 n+\frac{1}{2}が気持ち悪かったら、そこも展開すると、

 S = n\log{n} - n + \frac{1}{2}\log{n} + \frac{1}{12n} + 0.91893853320467317 + O(n^{-3})

となりました。