✨ feat(project):添加欧拉项目第6题解决方案及相关依赖
📝 docs(project):添加Faulhaber公式详细文档说明 ⬆️ chore(project):添加numpy依赖以支持数学计算
This commit is contained in:
43
solutions/0006.SumSquareDifference/eular_6.py
Normal file
43
solutions/0006.SumSquareDifference/eular_6.py
Normal file
@@ -0,0 +1,43 @@
|
||||
"""
|
||||
The sum of the squares of the first ten natural numbers is 1^2 + 2^2 + ... + 10^2 = 385,
|
||||
The square of the sum of the first ten natural numbers is (1+2+...+10)^2 = 55^2 = 3025,
|
||||
Hence the difference between the sum of the squares of the first ten natural numbers and the square of the sum is 3025 - 385 = 2640.
|
||||
Find the difference between the sum of the squares of the first one hundred natural numbers and the square of the sum.
|
||||
"""
|
||||
|
||||
import time
|
||||
|
||||
import numpy as np
|
||||
|
||||
|
||||
def timer(func):
|
||||
def wrapper(*args, **kwargs):
|
||||
start = time.time()
|
||||
result = func(*args, **kwargs)
|
||||
end = time.time()
|
||||
print(f"{func.__name__} took {end - start:.6f} seconds")
|
||||
return result
|
||||
|
||||
return wrapper
|
||||
|
||||
|
||||
@timer
|
||||
def sum_square_difference(n: int) -> int:
|
||||
res = 0
|
||||
for i in range(1, n + 1):
|
||||
res += sum(np.array(list(range(1, i))) * i)
|
||||
|
||||
return 2 * res
|
||||
|
||||
|
||||
@timer
|
||||
def better_sum_square(n: int) -> int:
|
||||
"""https://en.wikipedia.org/wiki/Faulhaber%27s_formula"""
|
||||
su = ((1 + n) * n // 2) ** 2
|
||||
sq = (2 * n + 1) * (n + 1) * n // 6
|
||||
return su - sq
|
||||
|
||||
|
||||
if __name__ == "__main__":
|
||||
print(sum_square_difference(100))
|
||||
print(better_sum_square(100))
|
||||
183
solutions/0006.SumSquareDifference/readme.md
Normal file
183
solutions/0006.SumSquareDifference/readme.md
Normal file
@@ -0,0 +1,183 @@
|
||||
## 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 代码框架:
|
||||
|
||||
```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
|
||||
```
|
||||
|
||||
### 计算幂和:
|
||||
|
||||
```python
|
||||
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$):需要更高级算法和模运算
|
||||
Reference in New Issue
Block a user