算法挑战：为浮点生成连续分数
（ 编辑 ：回应脾气暴躁的评论，不，这不是功课。我正在进行音调检测，获取一系列潜在的谐波峰值，并尝试构建基频的候选者。因此，它实际上是一个非常实际的问题。 ）
考虑（例如）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]总是最佳有理逼近，但它们并非都是最合理的近似。 为了获得所有这些，你还必须采用半收敛/中介 – forms的分数(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的整数是来自任一层（a [k + 2] / 2）或上限（a [ k + 2] / 2），到[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找到连续的分数（将项放在[]中，收敛项放在p []和q []中 – 请原谅全局变量），函数all_best打印所有最合理的近似值。
#include #include #include // 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; i0); 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,
（参见Fibonacci数？这里的收敛数都是近似值。）
或者像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,
请享用！