數值代數–理查森外推法

實驗四

一、實驗名稱

理查森外推算法

二、實驗目的與要求：

實驗目的：掌握理查森外推算法。

實驗要求：1. 給出理查森外推算法思路，

2. 用C語言實現算法，運行環境為Microsoft Visual C++。

三、算法思路：

1. 假設函數泰勒展開式可表示為 和，將兩式相減，消去偶數項，則，整理得到下式，記L表示，表示微分形式，則有 (1)用h/2代替h，有 (2)，由(1)(2)兩式子有推廣這種方法，就是理查森外推法了。

2. 理查森外推法公式 , , 用下列公式計算，k=1,2,…,M，n=k,k+1,…，M。 則有，當n和k足夠大時D(n,k)可充分接近。

3. 上機算法 input h , M

for n=0 to M do

D(n , 0)

end do

for k=1 to M do

for n=k to M do end do

end do

output D(n , k)

四、實驗題目：

五、問題的解：

編寫程序(程序見后面附錄)，輸出結果如下：

分析得到的結果，發現在對角線附近D(n , k)的值越來越穩定，通過上面算法闡述，我們知道D(n , k)應該是越來越接近我們想求到的導數的，與實驗結果一致。

六、附錄：

實驗編程，運行環境為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);

}

double 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 D1[4][4],D2[5][5],D3[6][6];

int i,j;

for(i=0;i<=3;i++) /*第一個問題的理查森算法*/

D1[i][0]=f1(1.0/pow(2,i));

for(j=1;j<=3;j++)

for(i=j;i<=3;i++)

D1[i][j]=D1[i][j-1]+(D1[i][j-1]-D1[i-1][j-1])/(pow(4,j)-1);

printf("第一道題結果：\n");

for(i=0;i<=3;i++)

{for(j=0;j<=i;j++)

printf("%0.12f ",D1[i][j]);

printf("\n");

}

for(i=0;i<=4;i++) /*第二個問題的理查森算法*/

D2[i][0]=f2(1.0/pow(2,i));

for(j=1;j<=4;j++)

for(i=j;i<=4;i++)

D2[i][j]=D2[i][j-1]+(D2[i][j-1]-D2[i-1][j-1])/(pow(4,j)-1);

printf("第二道題結果：\n");

for(i=0;i<=4;i++)

{for(j=0;j<=i;j++)

printf("%0.12f ",D2[i][j]);

printf("\n");

}

for(i=0;i<=5;i++) /*第三個問題的理查森算法*/

D3[i][0]=f3(1.0/pow(2,i));

for(j=1;j<=5;j++)

for(i=j;i<=5;i++)

D3[i][j]=D3[i][j-1]+(D3[i][j-1]-D3[i-1][j-1])/(pow(4,j)-1);

printf("第三道題結果：\n");

for(i=0;i<=5;i++)

{for(j=0;j<=i;j++)

printf("%0.12f ",D3[i][j]);

printf("\n");

}

}

總結

以上是生活随笔為你收集整理的理查森外推法 matlab,数值代数–理查森外推法.doc的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。