知其然更知其所以然中国先哲ppt课件

1、 乐经良 知其然知其然, ,更知其更知其所所以然以然. . - -中国先哲中国先哲 哪里有数哪里有数, ,哪里就哪里就有美有美. . - Proclus- Proclus 乐经良数学实验数学实验上海交大数学科学学院上海交大数学科学学院 的 计 算 乐经良你也许能写出你也许能写出 = 3.1415926535 = 3.1415926535实际问题实际问题圆周率, 我们十分熟悉的常数.用用Matlab 可以求出可以求出 到几百位到几百位 digits(100) vpa(pi)但你会计算但你会计算 的值吗？你又能用几种方法计算？的值吗？你又能用几种方法计算？ans = 3.1415926535897

2、93238462643383279502884197169399375105820974944592307816406286208998628034825342117068 乐经良刘徽割圆法刘徽割圆法从正六边形开始，逐步求边长与面积oABCD22222)(ODOCADDCADAC如图2222142)2(11)2(nnnnaaaa相应相应OAC的面积的面积4221211nnnaaADOCS设边数为设边数为62n62n的正多边形边长为的正多边形边长为 anan递推法递推法 乐经良 于是于是 的值的值 ( (刘徽计算到刘徽计算到192192边形面积，得到边形面积，得到 3.141) 3.141)nn

3、nnaS232611用用MatlabMatlab计算计算 m m文件文件function calpi(n)a(1)=1;for i=1:n-1 a(i+1)=sqrt(2-sqrt(4-a(i)2);end S=3*2(n-1)*a(n)命令窗口输入命令窗口输入 format long gcalpi(5) 乐经良如何提高精度如何提高精度提高多边形的边数提高多边形的边数 不能完全达到目的不能完全达到目的在在MatlabMatlab文件中解决文件中解决 符号运算符号运算function calpi1(n)a(1)=sym(1);for i=1:n-1 a(i+1)=sym(sqrt(2-sqrt(

4、4-a(i)2); endS=3*2(n-1)*a(n);vpa(S,60) %最后进行数值计算，60为数值计算过程中保留的有效数字 乐经良任务任务1 1 德国人鲁道夫用一生计算圆周率。德国人鲁道夫用一生计算圆周率。他同样是用圆的内接多边形逼近圆周，他同样是用圆的内接多边形逼近圆周，不过他是从正方形开始成倍增加边数。不过他是从正方形开始成倍增加边数。试推导出他计算所采用的递推公式，然试推导出他计算所采用的递推公式，然后求后求的近似值到的近似值到1010位和位和2020位位. . 乐经良121)1(5131141nn12)1(53arctan12153nxxxxxnn221422) 1(111n

5、nxxxx利用幂级数计算利用幂级数计算积分导出积分导出取取 x=1 乐经良121)1(51311 41nnnS4偶数奇数nnSnnSSnnn1211211111S ( Sn的迭代格式 ) 乐经良用用MatlabMatlab计算计算创建创建m m文件文件 calpi 2.m calpi 2.m，内容如下：，内容如下：function calpi2(n)S=0;for i=1:n if mod(i,2)=0 S=S-1/(2*i-1); else S=S+1/(2*i-1); endendS=4*S 乐经良 calpi2(1000)ans = 3.14059265383979 calpi2(100

6、00)ans = 3.14149265359003结果如何结果如何 ? ? calpi2(20000)ans = 3.14154265358982精度提高很慢精度提高很慢! ! 乐经良 能不能算得更快一点、更精确一点？能不能算得更快一点、更精确一点？431arctan21arctan12153)21(12)1()21(51)21(31214nnnMachin Machin 公式公式42391arctan51arctan4简单公式简单公式)31(12)1()31(51)31(313112153nnn 乐经良用用MatlabMatlab创建创建m m文件文件 function calpi2_1(n

7、)S=0;for i=1:n if mod(i,2)=0 S=S-1/(2*i-1)*(1/(2(2*i-1)+1/(3(2*i-1); else S=S+1/(2*i-1)*(1/(2(2*i-1)+1/(3(2*i-1); endendS=vpa(4*S,30)%观察30位有效数字 乐经良 calpi2_1(10)ans = 3.14159257960635063255949717131计算结果计算结果 calpi2_1(20)ans = 3.14159265358975625659354591335 calpi2_1(50)ans = 3.141592653589793238462643

8、38328 乐经良一个结论一个结论算法很重要算法很重要计算机速度计算机速度300300次次/ /秒秒 33.8633.861040 1040 兆兆/ /秒秒从从1950 1950 2000 2000年年104104次次 / /秒秒 1012 1012次次 / /秒，提高秒，提高1 1亿倍亿倍 算法算法 ( (解线性方程组解线性方程组 高斯消去法高斯消去法 多重网格法多重网格法) ) 计算机速度计算机速度运算次数：运算次数：10181018次次 106 106次，提高次，提高1 1万亿倍万亿倍 乐经良任务任务2 281arctan51arctan21arctan42) 2) 验证公式验证公式1)

9、 1) 用反正切函数的幂级数展开式结合有关公式用反正切函数的幂级数展开式结合有关公式 简单公式和简单公式和MachinMachin公式所用的项数公式所用的项数. .求求，若要精确到，若要精确到4040位、位、5050位数字，试比较位数字，试比较试试用此公式右端作幂级数展开完成任务试试用此公式右端作幂级数展开完成任务1)1)所需的项数所需的项数 乐经良3 3回忆在微积分中学习到的其它级数形回忆在微积分中学习到的其它级数形式是否可用来求式是否可用来求的值到的值到1010位、位、2020位、位、3030位，相应需要级数的多少项？位，相应需要级数的多少项？ 乐经良dxxA102114211)(xxy设

10、将0,1区间n等分，取xk=k/n, )(220121nnyyyyynA梯形法)(2)(61224220mmyyyyymASimpson法 利用数值积分方法利用数值积分方法yk= 1/ (1+xk2)(41231myyy 乐经良Matlab 计算创建创建 m m 文件文件梯形法梯形法function calpi3(n)x=0:1/n:1;y=1./(1+x.2);S=2*sum(y)-1-0.5;2*S/n 乐经良 calpi3(100)ans = 3.14157598692313 calpi3(500)ans = 3.14159198692313 calpi3(10000)ans = 3.1

11、4159265192314 乐经良 用数值积分计算用数值积分计算，分别用梯形法和，分别用梯形法和SimpsonSimpson法精确到法精确到1010位数字，用位数字，用SimpsonSimpson法精确到法精确到1515位数字位数字. .任务任务3 3 乐经良针与平行线相交的次数为针与平行线相交的次数为n n Monte Carlo Monte Carlo 法法从从BuffonBuffon落针实验谈起：落针实验谈起： 纸上一组平行线距离为纸上一组平行线距离为1 1，将长度为将长度为1 1的针多次地扔到的针多次地扔到纸上。若扔针次数为纸上。若扔针次数为m,m,而其中而其中Buffon指出： 的数

12、值与 m/n 有关，他由此求出求出 的近似值为的近似值为3.1423.142 乐经良设计方案设计方案 4 m/n4 m/n 计算机模拟：产生区间计算机模拟：产生区间0,10,1上数目为上数目为n n 的一的一组组 在正方形在正方形 0 x 1, 0 0 x 1, 0 y1y1上随机的投大量的点，那么上随机的投大量的点，那么落在四分之一园内的点数落在四分之一园内的点数数数m m与在正方形内的点数与在正方形内的点数n n之比之比 m/n m/n 应为这两部分图形应为这两部分图形面积之比面积之比 /4, /4, 故故随机数随机数(x, y)(x, y)，计算满足，计算满足 x2 + y2 1 x2

13、+ y2 1 的点数的点数m m 乐经良Matlab计算创建创建m m文件文件function y=calpi4(k)m=0;for n=1:k if rand(1)2+rand(1)2=1 m=m+1; end;end;y=4*m/k;观察结果，思考为什么？ 乐经良2) 2) 设计方案用计算机模拟设计方案用计算机模拟BuffonBuffon实验实验任务任务4 4看能否求得看能否求得5 5位精确数字？位精确数字？1) 1) 用用Monte Carlo Monte Carlo 法计算法计算，除了加大随机数，除了加大随机数，在随机数一定时可重复算若干次后求平均值，在随机数一定时可重复算若干次后求平

14、均值， 乐经良其他方法其他方法 1/ 的展开式 Ramanujan 公式 044396) !()263901103()!4(9801221nnnnn 算术几何平均值迭代法算术几何平均值迭代法nnnnnnnnnnbaMbabbaabalimlim,2,21, 1110020222)(211Mbannnn 乐经良利用积分利用积分dxxn20sin为奇数nnn,2! !)!1(推导公式推导公式122122543432122nnnn任务任务5 5用此公式计算用此公式计算的近似值，效果如何？的近似值，效果如何？ 乐经良有必要计算那么精确吗有必要计算那么精确吗十位小数就足以使地球周界准确到一英十位小数就足

15、以使地球周界准确到一英寸以内，三十位小数便能使整个可见宇宙周寸以内，三十位小数便能使整个可见宇宙周边准确到连最强大的显微镜都不能分辨的一个边准确到连最强大的显微镜都不能分辨的一个度量，度量， 三十五位小数的三十五位小数的值计算能把太阳系值计算能把太阳系包围起来的圆的周长，误差还不到质子直径的包围起来的圆的周长，误差还不到质子直径的百万分之一百万分之一 乐经良计算计算 的意义的意义 反映数学和计算技术发展的一个侧面 “历史上一个国家所算得的圆周率的准确程度，可以作为衡量这个国家当时数学发展水平的指标。” 3.1415926 3.1415927 113355722(领先世界900余年) 乐经良位数

16、位数 100 100万万 10.1 10.1亿亿 2061 2061亿亿 12411 12411亿亿 2.7 2.7万亿万亿 10 10万亿万亿年代年代 1973 1989 2019 2019 2019.1 2019.10 1973 1989 2019 2019 2019.1 2019.10国家国家 法国法国 美国美国 日本日本 日本日本 法国法国 日本日本人工计算：实验法人工计算：实验法 几何法几何法 分析法分析法用计算机用计算机: : 始于始于19491949 2035 2035 位位 （Von NeumannVon Neumann）最高记录：最高记录：808808位位(1948 )(19

17、48 )家用电脑家用电脑运用超级计算机运用超级计算机 已算到小数位已算到小数位6060万亿位万亿位 乐经良 引发新的概念、方法和思想引发新的概念、方法和思想 , ,产生新的问题产生新的问题 测试或检验超级计算机的各项性能测试或检验超级计算机的各项性能(Super (Super PI) PI) IntelIntel公司推出奔腾公司推出奔腾PentiumPentium时发现问题时发现问题 雅虎科技公司的研究员尼古拉斯雅虎科技公司的研究员尼古拉斯斯则斯则(Nicholas Sze)(Nicholas Sze)，采用，采用“云计算技术，利用云计算技术，利用10001000台电脑同时计算，历时台电脑同时计算，历时2323天，将圆周率精确天，将圆周率精确到小数点后到小数点后2 2千万亿位千万亿位 乐经良1) 1) 利用学习过的知