（编辑：针对脾气暴躁的评论，不，这不是家庭作业。我正在研究音高检测，获取一系列潜在的谐波峰值，并尝试构建基频的候选对象。所以，这实际上是一个非常实际的问题。 )
考虑（例如）pi 的最佳分数近似值，按分母递增排序：3/1、22/7、355/113、...
挑战：创建一个整洁的C 算法，该算法将为给定的浮点数生成第 n 个商近似 a/b，同时返回差异。
calcBestFrac(float frac, int n, int * a, int * b, float * err) {...}
我认为最好的技术是连分数
去掉 pi 的小数部分，得到 3
现在，余数是 0.14159... = 1/7.06251..
所以下一个最好的有理数是 3 + 1/7 = 22/7
从 7.06251 中减去 7 得到 0.06251.. 大约是 1/15.99659..
称它为 16，那么下一个最佳近似值是
3 + 1/(7 + 1/16) = 355/113
但是，将其转换为干净的 C 代码绝非易事。如果我得到一些整洁的东西，我会发布。与此同时，有人可能会喜欢它作为脑筋急转弯。
[因为您要求将此作为答案而不是评论。]
对于任何实数，其连分数的收敛点 p[k]/q[k] 始终是最佳有理逼近，但它们并非都是最佳有理逼近。要获得所有这些，您还必须采用半收敛/中间值——(p[k]+n*p[k+1])/(q[k]+n*q[k+1])某个整数 n≥1 的形式的分数。取 n=a[k+2] 得到 p[k+2]/q[k+2]，取的整数 n 是来自 floor(a[k+2]/2) 或 ceiling(a[ k+2]/2)，到 a[k+2]。维基百科上也提到了这一点。
π 的连分数是 [3; 7, 15, 1, 292, 1, 1, 1, 2, 1, 3, 1, 14, 2…]（OEIS中的序列A001203），收敛的序列是3/1, 22/7, 333/106 , 355/113, 103993/33102… ( A002485 / A002486 ), 最佳逼近序列为 3/1, 13/4, 16/5, 19/6, 22/7, 179/57… ( A063674 / A063673）。
所以算法说 π = [3; 7, 15, 1, 292, 1, 1,…] 是
3/1 = [3]
13/4 = [3; 4]
16/5 = [3; 5]
19/6 = [3; 6]
22/7 = [3; 7]
179/57 = [3; 7, 8]
201/64 = [3; 7, 9]
223/71 = [3; 7, 10]
245/78 = [3; 7, 11]
267/85 = [3; 7, 12]
289/92 = [3; 7, 13]
311/99 = [3; 7, 14]
333/106 = [3; 7, 15]
355/113 = [3; 7, 15, 1]
52163/16604 = [3; 7, 15, 1, 146]
52518/16717 = [3; 7, 15, 1, 147]
… (all the fractions from [3; 7, 15, 1, 148] to [3; 7, 15, 1, 291])…
103993/33102 = [3; 7, 15, 1, 292]
104348/33215 = [3; 7, 15, 1, 292, 1]
...
这是一个 C 程序，它给定一个正实数，生成它的连分数、它的收敛和最佳有理逼近序列。该函数find_cf找到连分数（将项放在 a[] 中，将收敛项放在 p[] 和 q[] 中——排除全局变量），并且该函数all_best打印所有最佳有理近似值。
#include <math.h>
#include <stdio.h>
#include <assert.h>
// number of terms in continued fraction.
// 15 is the max without precision errors for M_PI
#define MAX 15
#define eps 1e-9
long p[MAX], q[MAX], a[MAX], len;
void find_cf(double x) {
int i;
//The first two convergents are 0/1 and 1/0
p[0] = 0; q[0] = 1;
p[1] = 1; q[1] = 0;
//The rest of the convergents (and continued fraction)
for(i=2; i<MAX; ++i) {
a[i] = lrint(floor(x));
p[i] = a[i]*p[i-1] + p[i-2];
q[i] = a[i]*q[i-1] + q[i-2];
printf("%ld: %ld/%ld\n", a[i], p[i], q[i]);
len = i;
if(fabs(x-a[i])<eps) return;
x = 1.0/(x - a[i]);
}
}
void all_best(double x) {
find_cf(x); printf("\n");
int i, n; long cp, cq;
for(i=2; i<len; ++i) {
//Test n = a[i+1]/2. Enough to test only when a[i+1] is even, actually...
n = a[i+1]/2; cp = n*p[i]+p[i-1]; cq = n*q[i]+q[i-1];
if(fabs(x-(double)cp/cq) < fabs(x-(double)p[i]/q[i]))
printf("%ld/%ld, ", cp, cq);
//And print all the rest, no need to test
for(n = (a[i+1]+2)/2; n<=a[i+1]; ++n) {
printf("%ld/%ld, ", n*p[i]+p[i-1], n*q[i]+q[i-1]);
}
}
}
int main(int argc, char **argv) {
double x;
if(argc==1) { x = M_PI; } else { sscanf(argv[1], "%lf", &x); }
assert(x>0); printf("%.15lf\n\n", x);
all_best(x); printf("\n");
return 0;
}
对于 π，这是该程序的输出，大约在 0.003 秒内（即，它确实比遍历所有可能的分母更好！），换行以提高可读性：
% ./a.out
3.141592653589793
3: 3/1
7: 22/7
15: 333/106
1: 355/113
292: 103993/33102
1: 104348/33215
1: 208341/66317
1: 312689/99532
2: 833719/265381
1: 1146408/364913
3: 4272943/1360120
1: 5419351/1725033
14: 80143857/25510582
13/4, 16/5, 19/6, 22/7, 179/57, 201/64, 223/71, 245/78, 267/85, 289/92, 311/99,
333/106, 355/113, 52163/16604, 52518/16717, 52873/16830, 53228/16943, 53583/17056,
53938/17169, 54293/17282, 54648/17395, 55003/17508, 55358/17621, 55713/17734,
56068/17847, 56423/17960, 56778/18073, 57133/18186, 57488/18299, 57843/18412,
58198/18525, 58553/18638, 58908/18751, 59263/18864, 59618/18977, 59973/19090,
60328/19203, 60683/19316, 61038/19429, 61393/19542, 61748/19655, 62103/19768,
62458/19881, 62813/19994, 63168/20107, 63523/20220, 63878/20333, 64233/20446,
64588/20559, 64943/20672, 65298/20785, 65653/20898, 66008/21011, 66363/21124,
66718/21237, 67073/21350, 67428/21463, 67783/21576, 68138/21689, 68493/21802,
68848/21915, 69203/22028, 69558/22141, 69913/22254, 70268/22367, 70623/22480,
70978/22593, 71333/22706, 71688/22819, 72043/22932, 72398/23045, 72753/23158,
73108/23271, 73463/23384, 73818/23497, 74173/23610, 74528/23723, 74883/23836,
75238/23949, 75593/24062, 75948/24175, 76303/24288, 76658/24401, 77013/24514,
77368/24627, 77723/24740, 78078/24853, 78433/24966, 78788/25079, 79143/25192,
79498/25305, 79853/25418, 80208/25531, 80563/25644, 80918/25757, 81273/25870,
81628/25983, 81983/26096, 82338/26209, 82693/26322, 83048/26435, 83403/26548,
83758/26661, 84113/26774, 84468/26887, 84823/27000, 85178/27113, 85533/27226,
85888/27339, 86243/27452, 86598/27565, 86953/27678, 87308/27791, 87663/27904,
88018/28017, 88373/28130, 88728/28243, 89083/28356, 89438/28469, 89793/28582,
90148/28695, 90503/28808, 90858/28921, 91213/29034, 91568/29147, 91923/29260,
92278/29373, 92633/29486, 92988/29599, 93343/29712, 93698/29825, 94053/29938,
94408/30051, 94763/30164, 95118/30277, 95473/30390, 95828/30503, 96183/30616,
96538/30729, 96893/30842, 97248/30955, 97603/31068, 97958/31181, 98313/31294,
98668/31407, 99023/31520, 99378/31633, 99733/31746, 100088/31859, 100443/31972,
100798/32085, 101153/32198, 101508/32311, 101863/32424, 102218/32537, 102573/32650,
102928/32763, 103283/32876, 103638/32989, 103993/33102, 104348/33215, 208341/66317,
312689/99532, 833719/265381, 1146408/364913, 3126535/995207,
4272943/1360120, 5419351/1725033, 42208400/13435351, 47627751/15160384,
53047102/16885417, 58466453/18610450, 63885804/20335483, 69305155/22060516,
74724506/23785549, 80143857/25510582,
所有这些术语都是正确的，但如果你增加 MAX，你会因为精度而开始出错。我对仅使用 13 个收敛得到的项数印象深刻。（如您所见，有一个小错误，有时它不会打印第一个“n/1”近似值，或者打印不正确——我留给您修复！）
你可以试试 √2，它的连分数是 [1; 2, 2, 2, 2…]:
% ./a.out 1.41421356237309504880
1.414213562373095
1: 1/1
2: 3/2
2: 7/5
2: 17/12
2: 41/29
2: 99/70
2: 239/169
2: 577/408
2: 1393/985
2: 3363/2378
2: 8119/5741
2: 19601/13860
2: 47321/33461
3/2, 4/3, 7/5, 17/12, 24/17, 41/29, 99/70, 140/99, 239/169, 577/408, 816/577, 1393/985, 3363/2378, 4756/3363, 8119/5741, 19601/13860, 47321/33461,
或黄金比例 φ = (1+√5)/2 其连分数为 [1; 1, 1, 1, …]:
% ./a.out 1.61803398874989484820
1.618033988749895
1: 1/1
1: 2/1
1: 3/2
1: 5/3
1: 8/5
1: 13/8
1: 21/13
1: 34/21
1: 55/34
1: 89/55
1: 144/89
1: 233/144
1: 377/233
2/1, 3/2, 5/3, 8/5, 13/8, 21/13, 34/21, 55/34, 89/55, 144/89, 233/144, 377/233,
（见斐波那契数列？这里的收敛值都是近似值。）
或者像 4/3 = [1; 这样的有理数 3]：
% ./a.out 1.33333333333333333333
1.333333333333333
1: 1/1
3: 4/3
3/2, 4/3,
或 14/11 = [1; 3、1、2]：
% ./a.out 1.27272727272727272727
1.272727272727273
1: 1/1
3: 4/3
1: 5/4
2: 14/11
3/2, 4/3, 5/4, 9/7, 14/11,
享受！
C 程序很好，除了你不能相信对余数的检查，从计算 x*pq 也可以看出：
Iteration #1: 3: 3/1 - delta: 0.141592653589793116, rem: 0.141592653589793116
Iteration #2: 7: 22/7 - delta: -0.008851424871448188, rem: 0.062513305931051878
Iteration #3: 15: 333/106 - delta: 0.008821280518070296, rem: 0.996594406684156776
Iteration #4: 1: 355/113 - delta: -0.000030144353377892, rem: 0.003417231014946418
Iteration #5: 292: 103993/33102 - delta: 0.000019129331725765, rem: 0.634590879621879211
Iteration #6: 1: 104348/33215 - delta: -0.000011015021655680, rem: 0.575818424298580694
Iteration #7: 1: 208341/66317 - delta: 0.000008114310077190, rem: 0.736658567704091524
Iteration #8: 1: 312689/99532 - delta: -0.000002900711592702, rem: 0.357480987585133375
Iteration #9: 2: 833719/265381 - delta: 0.000002312886920208, rem: 0.797351564778957706
Iteration #10: 1: 1146408/364913 - delta: -0.000000587824615650, rem: 0.254151925163927682
Iteration #11: 3: 4272943/1360120 - delta: 0.000000549413016415, rem: 0.934654436927838420
Iteration #12: 1: 5419351/1725033 - delta: -0.000000038411599235, rem: 0.069914142051204637
Iteration #13: 14: 80143857/25510582 - delta: 0.000000011648808140, rem: 0.303257833981669641
Iteration #14: 3: 245850922/78256779 - delta: -0.000000003463355824, rem: 0.297524047014214316
Iteration #15: 3: 817696623/260280919 - delta: 0.000000001280568540, rem: 0.361072861287829440
Iteration #16: 2: 1881244168/598818617 - delta: -0.000000000931322575, rem: 0.769524124392304913
Iteration #17: 1: 2698940791/859099536 - delta: 0.000000000232830644, rem: 0.299504418772708979
Iteration #18: 3: 9978066541/3176117225 - delta: 0.000000000000000000, rem: 0.338848902789946401 ******* 'true' deviation below epsilon threshold