数学研发网设为首页收藏本站
切换到窄版

 找回密码
 欢迎注册
数学研发论坛»论坛 › 【算法研发】 › 算法交流 › 一个计算圆周率任意精度的spigot算法研究 ...
12345678下一页
返回列表 发新帖
楼主: plp626

[讨论] 一个计算圆周率任意精度的spigot算法研究

[复制链接]
gracias
发表于 2011-4-24 11:07:50 | 显示全部楼层
http://user.qzone.qq.com/414353346/blog/1288275811 http://numbers.computation.free. ... grams/tinycodes.pdf
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

plp626
 楼主| 发表于 2011-4-24 11:38:11 | 显示全部楼层
11# gracias 恩，后面那个资料我看过，所以才如此感兴趣，只是两年了，对于细节我又忘记了；但大概的思路我还有； 我是要问任意整数开m次方的那个收敛级数，这个问题我已经有了想法，就是如9楼所讨论的 若不限定b=1，更一般的情况就是，解关于正整数a,p,q,非0整数b 的不定方程：
$a*q^m-(a-b)*N*p^m=0$
其中N,m为正整数常数，p,q互质，a,b互质，|b|
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

plp626
 楼主| 发表于 2011-4-25 06:53:05 | 显示全部楼层
http://user.qzone.qq.com/414353346/blog/1288275811 http://numbers.computation.free. ... grams/tinycodes.pdf gracias 发表于 2011-4-24 11:07
谢谢你提供的关于这个算法的详细资料，认真看过后有种相见恨晚的感觉； http://www.comlab.ox.ac.uk/resea ... lgprog-20031107.pdf http://www.cs.umb.edu/~offner/files/pi.pdf http://www.mathpropress.com/stan/bibliography/spigot.pdf http://numbers.computation.free. ... grams/tinycodes.pdf
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

plp626
 楼主| 发表于 2011-4-25 06:53:48 | 显示全部楼层
言归正传，继续思考这个不定方程的解法。。
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

liangbch
发表于 2011-4-25 17:47:22 | 显示全部楼层
5# plp626 谢谢楼主，有时间好好学下批处理。
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

▄︻┻═┳一‥
发表于 2011-4-26 08:58:43 | 显示全部楼层
本帖最后由 ▄︻┻═┳一‥ 于 2011-4-26 11:27 编辑 当m=2时，也就是开方的时候，用程序暴力搜索1000以内的所有整数解得到： 只是俺的这个程序算法太笨蛋了，运行了半天才得到这些
  2. +------+------+------+------+
  3. | N | a | p | q |
  4. +------+------+------+------+
  5. | 2 | 9 | 3 | 4 |
  6. | 2 | 50 | 5 | 7 |
  7. | 2 | 289 | 17 | 24 |
  8. | 3 | 4 | 2 | 3 |
  9. | 3 | 49 | 7 | 12 |
  10. | 3 | 676 | 26 | 45 |
  11. | 5 | 81 | 9 | 20 |
  12. | 6 | 25 | 5 | 12 |
  13. | 6 | 243 | 9 | 22 |
  14. | 7 | 64 | 8 | 21 |
  15. | 8 | 9 | 3 | 8 |
  16. | 8 | 50 | 5 | 14 |
  17. | 8 | 289 | 17 | 48 |
  18. | 10 | 361 | 19 | 60 |
  19. | 11 | 100 | 10 | 33 |
  20. | 12 | 49 | 7 | 24 |
  21. | 12 | 676 | 13 | 45 |
  22. | 13 | 325 | 5 | 18 |
  23. | 14 | 8 | 2 | 7 |
  24. | 14 | 225 | 15 | 56 |
  25. | 15 | 16 | 4 | 15 |
  26. | 15 | 961 | 31 | 120 |
  27. | 18 | 50 | 5 | 21 |
  28. | 18 | 289 | 17 | 72 |
  29. | 20 | 81 | 9 | 40 |
  30. | 21 | 28 | 2 | 9 |
  31. | 22 | 99 | 3 | 14 |
  32. | 23 | 576 | 24 | 115 |
  33. | 24 | 25 | 5 | 24 |
  34. | 24 | 243 | 9 | 44 |
  35. | 27 | 4 | 2 | 9 |
  36. | 27 | 49 | 7 | 36 |
  37. | 27 | 676 | 26 | 135 |
  38. | 28 | 64 | 4 | 21 |
  39. | 30 | 121 | 11 | 60 |
  40. | 32 | 9 | 3 | 16 |
  41. | 32 | 50 | 5 | 28 |
  42. | 32 | 289 | 17 | 96 |
  43. | 33 | 12 | 2 | 11 |
  44. | 33 | 529 | 23 | 132 |
  45. | 34 | 18 | 3 | 17 |
  46. | 35 | 36 | 6 | 35 |
  47. | 39 | 625 | 25 | 156 |
  48. | 40 | 361 | 19 | 120 |
  49. | 42 | 169 | 13 | 84 |
  50. | 44 | 100 | 5 | 33 |
  51. | 45 | 81 | 3 | 20 |
  52. | 48 | 49 | 7 | 48 |
  53. | 48 | 676 | 13 | 90 |
  54. | 50 | 9 | 3 | 20 |
  55. | 50 | 289 | 17 | 120 |
  56. | 52 | 325 | 5 | 36 |
  57. | 54 | 25 | 5 | 36 |
  58. | 54 | 243 | 3 | 22 |
  59. | 54 | 755 | 349 | 938 |
  60. | 55 | 45 | 3 | 22 |
  61. | 56 | 225 | 15 | 112 |
  62. | 57 | 76 | 2 | 15 |
  63. | 58 | 887 | 959 | 134 |
  64. | 60 | 16 | 2 | 15 |
复制代码
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

▄︻┻═┳一‥
发表于 2011-4-26 09:28:49 | 显示全部楼层
上面输出表面有些N在1000内还无解 猜想某些N也许在计算机的32位整数内无解。 看来这个级数对于N开m次方不具有通用性。
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

liangbch
发表于 2011-4-26 09:43:19 | 显示全部楼层
// 对于一个给定的浮点数，求一个一定范围内的p/q,使得p/q近似等于这个浮点数，有1种2分迭代法，可迅速逼近，下面给出c语言代码
  1. #include <stdlib.h>
  2. #include <stdio.h>
  3. #include <math.h>
  5. typedef unsigned long DWORD;
  7. void print_pq(double f,DWORD limit)
  8. {
  9. DWORD low_p,low_q,high_p,high_q,mid_p,mid_q;
  11. //初始化low_q,low_q,high_p,high_q
  12. low_q=DWORD(floor(f));
  13. high_q=DWORD(ceil(f));
  14. low_p=high_p=1;
  16. mid_p=low_p+high_p;
  17. mid_q=low_q+high_q;
  19. while ( mid_p<limit )
  20. {
  21. printf("%u/%u=%.12f\n",mid_q,mid_p,double(mid_q)/double(mid_p) );
  23. if ( double(mid_q)/double(mid_p) > f)
  24. {
  25. high_p=mid_p; high_q=mid_q;
  26. }
  27. else
  28. {
  29. low_p=mid_p; low_q=mid_q;
  30. }
  31. mid_p=low_p+high_p;
  32. mid_q=low_q+high_q;
  33. }
  34. }
  36. int main(int argc, char* argv[])
  37. {
  38. print_pq(sqrt(2),1000000); //求分子在1000000以内的sqrt(2)有理数逼近
  39. return 0;
  40. }
复制代码
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

liangbch
发表于 2011-4-26 09:44:38 | 显示全部楼层
这里给出楼上代码的运行结果
  1. 3/2=1.500000000000
  2. 4/3=1.333333333333
  3. 7/5=1.400000000000
  4. 10/7=1.428571428571
  5. 17/12=1.416666666667
  6. 24/17=1.411764705882
  7. 41/29=1.413793103448
  8. 58/41=1.414634146341
  9. 99/70=1.414285714286
  10. 140/99=1.414141414141
  11. 239/169=1.414201183432
  12. 338/239=1.414225941423
  13. 577/408=1.414215686275
  14. 816/577=1.414211438475
  15. 1393/985=1.414213197970
  16. 1970/1393=1.414213926777
  17. 3363/2378=1.414213624895
  18. 4756/3363=1.414213499851
  19. 8119/5741=1.414213551646
  20. 11482/8119=1.414213573100
  21. 19601/13860=1.414213564214
  22. 27720/19601=1.414213560533
  23. 47321/33461=1.414213562057
  24. 66922/47321=1.414213562689
  25. 114243/80782=1.414213562427
  26. 161564/114243=1.414213562319
  27. 275807/195025=1.414213562364
  28. 390050/275807=1.414213562382
  29. 665857/470832=1.414213562375
  30. 941664/665857=1.414213562372
复制代码
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

▄︻┻═┳一‥
发表于 2011-4-26 09:48:30 | 显示全部楼层
18# liangbch 你说的可是最佳有理逼近？ 那个是连分数，只是当分母超过 21亿的时候，就不好算了，要自己构造大数据类型做除法了。。。。 还是不如spigot算法 省事。。。
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

12345678下一页
返回列表 发新帖
高级模式
B Color Image Link Quote Code Smilies
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2024-11-22 05:09 , Processed in 0.026838 second(s), 15 queries .

Powered by Discuz! X3.5

© 2001-2024 Discuz! Team.

快速回复 返回顶部 返回列表