5.2 KiB
5.2 KiB
1. Faulhaber's Formula 概述
Faulhaber's formula 用于计算 前 n 个正整数的 p 次幂之和：
S_p(n) = \sum_{k=1}^n k^p = 1^p + 2^p + 3^p + \cdots + n^p
公式的一般形式为：
S_p(n) = \frac{1}{p+1} \sum_{j=0}^p \binom{p+1}{j} B_j \cdot n^{p+1-j}
其中 B_j 是伯努利数（Bernoulli numbers）。
特例
p = 0: S_0(n) = n
p = 1: $S_1(n) = \frac{n(n+1)}{2}$（三角形数）
p = 2: $S_2(n) = \frac{n(n+1)(2n+1)}{6}$（平方金字塔数）
p = 3: $S_3(n) = \left[\frac{n(n+1)}{2}\right]^2$（有趣的事实：等于 $S_1(n)^2$）
2. 历史背景
- Johann Faulhaber（1580-1635）：德国数学家，计算了
p=1到p=17的公式，但未发现通式 - Jacob Bernoulli（1654-1705）：独立研究并发现了通用公式，引入了伯努利数
- 公式因此以 Faulhaber 命名，但现代形式主要归功于 Bernoulli
3. 伯努利数（Bernoulli Numbers）详解
3.1 定义
伯努利数由生成函数定义：
\frac{x}{e^x - 1} = \sum_{k=0}^\infty B_k \frac{x^k}{k!}
3.2 递推关系
更实用的计算方法是递推公式：
\sum_{k=0}^{m} \binom{m+1}{k} B_k = 0 \quad \text{对于} \ m \ge 1
由此可得：
B_m = -\frac{1}{m+1} \sum_{k=0}^{m-1} \binom{m+1}{k} B_k
3.3 前几个伯努利数
| n | B_n |
说明 |
|---|---|---|
| 0 | 1 | B_0 = 1 |
| 1 | -1/2 | |
| 2 | 1/6 | |
| 3 | 0 | （所有奇数大于1的都是0） |
| 4 | -1/30 | |
| 5 | 0 | |
| 6 | 1/42 | |
| 7 | 0 | |
| 8 | -1/30 | |
| 9 | 0 | |
| 10 | 5/66 | |
| 12 | -691/2730 | 这是第一个分子大于1的 |
| 14 | 7/6 |
重要性质：
- 对于所有奇数 $n > 1$，有
B_n = 0 B_1 = -1/2是唯一的非零奇数项- 符号交替变化：
B_2 > 0, B_4 < 0, B_6 > 0, \ldots
3.4 计算示例
计算 $B_2$：
使用递推式（m=2）：
\binom{3}{0}B_0 + \binom{3}{1}B_1 + \binom{3}{2}B_2 = 0
1 \cdot 1 + 3 \cdot \left(-\frac{1}{2}\right) + 3 \cdot B_2 = 0
1 - \frac{3}{2} + 3B_2 = 0
3B_2 = \frac{1}{2}
B_2 = \frac{1}{6}
计算 $B_4$：
使用递推式（m=4）：
\sum_{k=0}^4 \binom{5}{k} B_k = 0
1 \cdot 1 + 5 \cdot \left(-\frac{1}{2}\right) + 10 \cdot \frac{1}{6} + 10 \cdot 0 + 5 \cdot B_4 = 0
1 - \frac{5}{2} + \frac{10}{6} + 5B_4 = 0
-\frac{3}{2} + \frac{5}{3} + 5B_4 = 0
-\frac{9}{6} + \frac{10}{6} + 5B_4 = 0
\frac{1}{6} + 5B_4 = 0
B_4 = -\frac{1}{30}
4. 使用 Faulhaber's Formula 计算具体例子
示例：计算 S_4(n) = \sum_{k=1}^n k^4
使用公式：
S_4(n) = \frac{1}{5} \sum_{j=0}^4 \binom{5}{j} B_j \cdot n^{5-j}
展开：
S_4(n) = \frac{1}{5}\left[\binom{5}{0}B_0 n^5 + \binom{5}{1}B_1 n^4 + \binom{5}{2}B_2 n^3 + \binom{5}{3}B_3 n^2 + \binom{5}{4}B_4 n\right]
代入伯努利数：
S_4(n) = \frac{1}{5}\left[1 \cdot 1 \cdot n^5 + 5 \cdot \left(-\frac{1}{2}\right) n^4 + 10 \cdot \frac{1}{6} n^3 + 10 \cdot 0 \cdot n^2 + 5 \cdot \left(-\frac{1}{30}\right) n\right]
化简：
S_4(n) = \frac{1}{5}\left[n^5 - \frac{5}{2}n^4 + \frac{10}{6}n^3 - \frac{1}{6}n\right]
= \frac{1}{5}\left[n^5 - \frac{5}{2}n^4 + \frac{5}{3}n^3 - \frac{1}{6}n\right]
= \frac{n^5}{5} - \frac{n^4}{2} + \frac{n^3}{3} - \frac{n}{30}
通分后：
S_4(n) = \frac{6n^5 - 15n^4 + 10n^3 - n}{30}
= \frac{n(n+1)(2n+1)(3n^2+3n-1)}{30}
5. 算法实现思路
计算伯努利数的 Python 代码框架：
def bernoulli_numbers(n):
"""
计算前n个伯努利数 B_0 到 B_n
"""
B = [0] * (n + 1)
B[0] = 1 # B_0 = 1
for m in range(1, n + 1):
sum_val = 0
for k in range(m):
# 二项式系数 binomial(m+1, k)
binom = binomial_coefficient(m + 1, k)
sum_val += binom * B[k]
B[m] = -sum_val / (m + 1)
return B
计算幂和：
def faulhaber_sum(p, n):
"""
使用Faulhaber公式计算 sum_{k=1}^n k^p
"""
B = bernoulli_numbers(p)
result = 0
for j in range(p + 1):
binom = binomial_coefficient(p + 1, j)
term = binom * B[j] * (n ** (p + 1 - j))
result += term
return result / (p + 1)
6. 重要关系与推广
6.1 与黎曼ζ函数的关系
对于偶数伯努利数：
B_{2k} = (-1)^{k+1} \frac{2(2k)!}{(2\pi)^{2k}} \zeta(2k)
例如：
\zeta(2) = \sum_{n=1}^\infty \frac{1}{n^2} = \frac{\pi^2}{6}B_2 = \frac{1}{6}
6.2 与多项式的关系
S_p(n) 是关于 n 的 p+1 次多项式，且：
- 常数项为 0
- 系数为有理数
- 首项系数为
\frac{1}{p+1}
6.3 生成函数视角
幂和也有生成函数表示：
\sum_{p=0}^\infty S_p(n) \frac{t^p}{p!} = \frac{t}{e^t - 1} \cdot \frac{e^{(n+1)t} - 1}{t}
这正是伯努利数生成函数与几何级数的乘积。
7. 计算复杂度与优化
- 直接计算：$O(n \cdot p)$（暴力求和）
- Faulhaber公式：$O(p^2)$（计算伯努利数）+ $O(p)$（求值）
- 对于固定p，多次查询不同n：预处理伯努利数后每次查询
O(p) - 对于大p（如 $p > 10^6$）：需要更高级算法和模运算