祖冲之的割圆术与圆周率

nyy

下面展示割圆术的计算过程

1. Clear["Global`*"];(*清除所有变量*)
   
2. a=1;(*六边形的初始边长*)
   
3. n=6;(*六边形边数*)
   
4. Do[a=Sqrt[2-Sqrt[4-a^2]];(*由n边形的边长,计算出2n边形的边长*)
   
5.    n=n*2;(*边数从n变成2n*)
   
6.    pi=a*n/2;(*周长除以直径=圆周率*)
   
7.    Print[{n,N[pi,10]}],(*打印出边数与圆周率的近似值*)
   
8. {k,1,14}]
   

复制代码

计算结果
{12,3.105828541}
{24,3.132628613}
{48,3.139350203}
{96,3.141031951}
{192,3.141452472}
{384,3.141557608}
{768,3.141583892}
{1536,3.141590463}
{3072,3.141592106}
{6144,3.141592517}
{12288,3.141592619}
{24576,3.141592645}
{49152,3.141592651}
{98304,3.141592653}

看起来没技术含量，但是在古代手算开平方根，还是很难的

圆周率差不多就等于下面的表达式
\[49152 \sqrt{2-\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{3}}}}}}}}}}}}}}}\]

其中49152=3*2^14，而上面的那个根号表达式，也有14次根号

nyy

\[a_{2n}=\sqrt{2-\sqrt{4-a_n^2}}\]
n边形的边长到2n边形的边长递推表达式

nyy

如果初始边长是根号2的正方形，那么圆周率的近似表达式

\[32768 \sqrt{2-\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2+\sqrt{2}}}}}}}}}}}}}}}\]
32768=2^15，上面也有15个根号

nyy

比如精确到小数点100位，应该多少边形？估计应该有个回归的关系

nyy

网上说祖冲之计算圆周率3.1415926，计算到了24576边形，但是根据我的计算，
似乎是只要12288就能计算到了。个人表示怀疑



祖冲之割圆术数百万次，算到了24576边形。

祖冲之认为刘徽算的结果还不算精确，他循着刘徽的方法继续算下去，排开小竹棍（算筹），从6边形算起，然后是12、24……192、384边形，边数每扩大一次，就要经过加减乘除以及开平方等11个步骤，同一运算过程要反复12次，一直算到24576边形，运算达数百万次。

这时，算筹已经从桌上摆到了地上，摆满了一屋，祖冲之还想向下计算，但已经实在无法计算了，只好就此停止。祖冲之首次将“圆周率”精算到小数第七位，即在3.1415926和3.1415927之间，他提出的“祖率”对数学的研究有重大贡献。

资料来源
https://zhidao.baidu.com/question/1763465724805890508.html

nyy

祖冲之的成果，放在今天，不到一个小时的计算就能超过他了，现代人太牛逼了

nyy

精度差不多是，3*2^n边形，精确到小数点后0.6n位

nyy

nyy 发表于 2023-12-21 09:22
精度差不多是，3*2^n边形，精确到小数点后0.6n位

鲁道夫的割圆术

1. (*鲁道夫的割圆术2^62边形*)
   
2. Clear["Global`*"];(*清除所有变量*)
   
3. a=Sqrt[2];(*正方形的初始边长*)
   
4. n=4;(*正方形边数*)
   
5. Do[a=Sqrt[2-Sqrt[4-a^2]];(*由n边形的边长,计算出2n边形的边长*)
   
6.    n=n*2;(*边数从n变成2n*)
   
7.    pi=a*n/2;(*周长除以直径=圆周率*)
   
8.    Print[{n,N[pi,10]}],(*打印出边数与圆周率的近似值*)
   
9. {k,1,60}]
   
10. N[pi-Pi,45]
    

复制代码

计算结果

1. -2.42984687184835126034142034455111983243470298*10^-37

复制代码

看这个结果，我感觉是精确到了小数点后36位，而不是35位。

大约1600年，荷兰数学家鲁道夫·范·科伊伦利用两千年前阿基米德创立的割圆术，耗尽一生心血，计算到了262 边形。你没看错就是这么多边形，得出圆周率小数点后35位，也就是3.14159 26535 89793 23846 26433 83279 50288。这个数字也叫鲁道夫数，后来他去世，墓碑上刻上了这个记录他一生成就的数字。

https://www.163.com/dy/article/EML8NCML05329IM0.html

也许他是因为计算误差，所以只到了35位，而不是36位

nyy

nyy 发表于 2023-12-22 09:34
鲁道夫的割圆术

确实是我计算错误！

1. (*鲁道夫的割圆术2^62边形*)
   
2. Clear["Global`*"];(*清除所有变量*)
   
3. a=Sqrt[2];(*正方形的初始边长*)
   
4. n=4;(*正方形边数*)
   
5. Do[a=Sqrt[2-Sqrt[4-a^2]];(*由n边形的边长,计算出2n边形的边长*)
   
6.    n=n*2;(*边数从n变成2n*)
   
7.    pi=a*n/2,(*周长除以直径=圆周率*)
   
8. {k,1,60}]
   
9. aa=N[pi,45](*鲁道夫的圆周率*)
   
10. bb=N[Pi,45](*真实的圆周率*)
    

复制代码

鲁道夫的圆周率、真实的圆周率分别如下：
3.14159265358979323846264338327950288395418471
3.14159265358979323846264338327950288419716940

确实是精确到35位，而不是36位

nyy

nyy 发表于 2023-12-22 11:28
确实是我计算错误！

鲁道夫的圆周率、真实的圆周率分别如下：

1. /^\(\S\+\)\(.*\)\n\zs\1\ze\(.*\)/e
   

复制代码

用vim命令这么搜索，就能搜索到共同部分。


被搜索内容如下：

1. 3.14159265358979323846264338327950288395418471
   
2. 3.14159265358979323846264338327950288419716940

复制代码