欧拉计划 Project Euler65（e的有理逼近）题解

题干

e的有理逼近

2的算术平方根可以写成无限连分数的形式
$$\sqrt{2}=1+\frac{1}{2+\frac{1}{2+\frac{1}{2+\frac{1}{2+\dots }}}}$$
这个无限连分数可以简记为$\sqrt{2}=[1;(2)]$，其中$(2)$表示2无限重复。同样的，我们可以记
$$\sqrt{23}=[4;(1,3,1,8)]$$
可以证明，截取算术平方根连分数表示的一部分所组成的序列，给出了一系列最佳有理逼近值。让我们
来考虑$\sqrt{2}$的逼近值
$$1+\frac{1}{2}=\frac{3}{2}$$
$$1+\frac{1}{2+\frac{1}{2}}=\frac{3}{5}$$
$$1+\frac{1}{2+\frac{1}{2+\frac{1}{2}}}=\frac{17}{12}$$
$$1+\frac{1}{2+\frac{1}{2+\frac{1}{2+\frac{1}{2}}}}=\frac{41}{29}$$
因此$\sqrt{2}$的前十个逼近值为：
$$1,\frac{3}{2},\frac{7}{5},\frac{17}{12},\frac{41}{29},\frac{99}{70},\frac{239}{169},\frac{577}{408},\frac{1393}{985},\frac{3363}{2378},\dots$$
最令人惊讶的莫过于重要的数学常数$e$有如下连分数表示
$$e=[2;1,2,1,1,4,1,1,6,1,\dots ,1,2k,1]$$
$e$的前十个逼近值为：
$$2,3,\frac{8}{3},\frac{11}{4},\frac{19}{7},\frac{87}{32},\frac{106}{39},\frac{193}{71},\frac{1264}{465},\frac{1457}{536},\dots$$
第$10$个逼近值的分子各位数字之和为$1+4+5+7=17$
求$e$的第100个逼近值的分子各位数字之和

思路

$e$的连分数表示为：
$$e=[2;1,2,1,1,4,1,1,6,1,\dots ,1,2k,1]$$
这意味着系数序列$a_0,a_1,a_2,\dots$为：
$$a_0=2,a_1=1,a_2=2,a_3=1,a_4=1,a_5=4,a_6=1,a_7=1,a_8=6,a_9=1,\dots$$
对于$n\ge 1$，序列的模式是$1,2k,1$重复，更准确度说

• 如果$n(mod3)=1orn(mod3)=0$，则$a_n=1$
• 如果$n(mod3)=2$，则$a_n=2k$， 其中$k=(n+1)/3$

连分数的$n$个逼近值（从第$0$个开始计数）为$\frac{p_n}{q_n}$，其分子$p_n$和分母$q_n$
可以通过以下递推式进行计算

• $p_n=a_n{p}_{n-1}+p_{n-2}$
• $q_n=a_n{q}_{n-1}+q_{n-2}$

初始条件为(因为从0开始算，所以这边有负数但实际编程的时候是不需要的)

• $p_{-2}=0,p_{-1}=1$
• $q_{-2}=1,q_{-1}=0$

我们可以带几个值进去验证一下

• $p_0=a_0{p}_{-1}+p_{-2}=2$
• $p_1=a_1{p}_0+p_{-1}=3$
• $p_2=a_2{p}_1+p_0=8$
• $p_3=a_3{p}_2+p_1=11$

有了以上思路后不难写出代码，如果使用c++编程要记得模拟大数加法和乘法，不然会溢出，使用python则不需要。

code

// 272
#include <bits/stdc++.h>using namespace std;using ll = long long;// 大数加法
string add(string num1, string num2) {if (num1 == "0") return num2;if (num2 == "0") return num1;string ans = "";int i = num1.length() - 1;int j = num2.length() - 1;int carry = 0; // 进位while (i >= 0 || j >= 0 || carry) {int d1 = (i >= 0) ? (num1[i--] - '0') : 0;int d2 = (j >= 0) ? (num2[j--] - '0') : 0;int cur = d1 + d2 + carry;ans += to_string(cur % 10);carry = cur / 10;}reverse(ans.begin(), ans.end());return ans.empty() ? "0" : ans;
}// 大数乘法
// 一种显然的思路是调用加法 但是当数字很大时就会很慢
// 所以我们手动模拟一下大数乘 这里涉及到的比较小
string mul(string num, int mulp) {if (num == "0" || mulp == 0) {return "0";}if (mulp == 1) {return num;}string ans = "";int carry = 0; // 进位for (int i = num.length() - 1; i >= 0; --i) {int di = num[i] - '0';ll pr = (long long)di * mulp + carry;ans += to_string(pr % 10);carry = pr / 10;}// 处理可能剩余的最高位进位while (carry > 0) {ans += to_string(carry % 10);carry /= 10;}reverse(ans.begin(), ans.end());return ans.empty() ? "0" : ans;
}void solve() {vector<int> a(100);a[0] = 2;for (int n = 1; n < 100; ++n) {if (n % 3 == 2) {int k = (n + 1) / 3;a[n] = 2 * k;} else {a[n] = 1;}}string pmin2 = "0"; // p_{i - 2}string pmin1 = "1"; // p_{i - 1};string p_n = "0"; // 当前p_{i}for (int i = 0; i < 100; ++i) {int cur = a[i];string ter1 = mul(pmin1, cur);p_n = add(ter1, pmin2);pmin2 = pmin1;pmin1 = p_n;}long long ans = 0;for (char di : p_n) {ans += (di - '0');}cout << ans << "\n";}int main() {ios::sync_with_stdio(false);cin.tie(nullptr);int tt = 1; // cin >> tt;while (tt--) {solve();}return 0;
}