《数值代数-理查森外推法》由会员分享,可在线阅读,更多相关《数值代数-理查森外推法(6页珍藏版)》请在金锄头文库上搜索。
1、实验四一、实验名称理查森外推算法二、实验目的与要求:实验目的:掌握理查森外推算法。实验要求:1. 给出理查森外推算法思路,2. 用 C 语言实现算法,运行环境为 Microsoft Visual C+。三、算法思路:1. 假设函数 泰勒展开式可表示为 )(xf 0)(!1)(kkxfhxf和 ,将两式相减,消去偶数项,0)(1(!)(kkxfhhf则 ,整理得到下0)12()!12()()(k kxfhxfxf式 ,记1)12( )!2()(21)( kkxfhxfhffL 表示 , 表示微分形式 ,)(xf)( )(fxf则有 (1)用 h/2 代替 h,有42hah(2) ,由(1) (2
2、)两6/)/(42式子有 推广这6/4/)(3)/(34Lhahh种方法,就是理查森外推法了。2. 理查森外推法公式 , , 用下)2/()0,(nhnDM0列公式计算 ,)1,(14),(14),( knDknDknkkk=1,2,M,n=k,k+1, ,M。 则有 ,(),(2khOL当 n 和 k 足够大时 D(n,k)可充分接近 。xf3. 上机算法 input h , Mfor n=0 to M doD(n , 0) )2/(nhend do for k=1 to M do for n=k to M do )14/(),1(),()1,(),( knDkknDend do end d
3、o output D(n , k) )0,(nkMn四、实验题目:五、问题的解:编写程序(程序见后面附录) ,输出结果如下:分析得到的结果,发现在对角线附近 D(n , k)的值越来越稳定,通过上面算法阐述,我们知道 D(n , k)应该是越来越接近我们想求到的导数 的,与实验结果一致。)(xf六、附录:实验编程,运行环境为 Microsoft Visual C+#include #include #include double f1(double x) /定义函数 f1(x)/ double y;y=(log(3.0+x)-log(3.0-x)/(2.0*x);return(y);doubl
4、e f2(double x) /定义函数 f2(x)/ double y;y=(tan(asin(0.8)+x)-tan(asin(0.8)-x)/(2.0*x);return(y);double f3(double x) /定义函数 f3(x)/ double y;y=(sin(x*x+x/3.0)-sin(x*x-x/3.0)/(2.0*x);return(y);void main()double D144,D255,D366;int i,j;for(i=0;i=3;i+) /*第一个问题的理查森算法*/D1i0=f1(1.0/pow(2,i);for(j=1;j=3;j+)for(i=j
5、;i=3;i+)D1ij=D1ij-1+(D1ij-1-D1i-1j-1)/(pow(4,j)-1);printf(第一道题结果:n);for(i=0;i=3;i+)for(j=0;j=i;j+)printf(%0.12f ,D1ij);printf(n);for(i=0;i=4;i+) /*第二个问题的理查森算法*/D2i0=f2(1.0/pow(2,i);for(j=1;j=4;j+)for(i=j;i=4;i+)D2ij=D2ij-1+(D2ij-1-D2i-1j-1)/(pow(4,j)-1);printf(第二道题结果:n);for(i=0;i=4;i+)for(j=0;j=i;j+)printf(%0.12f ,D2ij);printf(n);for(i=0;i=5;i+) /*第三个问题的理查森算法*/ D3i0=f3(1.0/pow(2,i);for(j=1;j=5;j+)for(i=j;i=5;i+)D3ij=D3ij-1+(D3ij-1-D3i-1j-1)/(pow(4,j)-1);printf(第三道题结果:n);for(i=0;i=5;i+)for(j=0;j=i;j+)printf(%0.12f ,D3ij);printf(n);