数值666 (2)

时间:2022-04-18 08:00:07 阅读: 最新文章 文档下载
说明:文章内容仅供预览,部分内容可能不全。下载后的文档,内容与下面显示的完全一致。下载之前请确认下面内容是否您想要的,是否完整无缺。
龙贝格求积分算法

4.龙贝格求积分算法

4.1算法说明



生成JK的逼近表R(J,K),并以R(J1,J1)为最终解来逼近积分













b

a

f(x)dxR(J,J)

逼近R(J,K)存在于一个特别的下三角矩阵中,0列元素R(J,0)用基于2J[a,b]子区间的连续梯形方法计算,然后利用龙贝格公式计算R(J,K)。当1KJ时,第J行的元素为













R(J,K)R(J,K1)

R(J,K1)R(J1,K1)



4K1

|R(J,J)R(J1,J1)|tol时,程序在第(J1)行结束。

4.2龙贝格求积分算法流程图

4-1 流程图



14 32


数值计算课程设计

4.3龙贝格求积分算法程序调试 求解积分方程果如图4-2

9

1

xdx为例,对龙贝格求积分算法程序进行编译和链接,执行后得到结



4.4龙贝格求积分算法程序代码

4-2 执行结果

#include #include

static double T[200][200];

double fun_x( double x ) //被积函数 {

return sqrt(x);}

double Romberg( double a, double b, double Eps ) {

int k=0; //用来记录把区间[a,b]2等分的次数 T[0][0] = (b-a)/2*( fun_x(a)+fun_x(b) ); do{ k++;

double temp=0;

for( int i=1; i<=pow(2,k-1); i++ ) {

temp += fun_x(a+(2*i-1)*(b-a)/pow(2,k));}

T[0][k] = 0.5*(T[0][k-1]+((b-a)/pow(2,k-1))*temp);

for( int m=1; m<=k; m++ ) {

T[m][k-m] = ( pow(4,m)*T[m-1][k-m+1]-T[m-1][k-m] )/( pow(4,m) - 1 );}

}while(fabs( T[k][0]-T[k-1][0] )>=Eps);

return T[k][0];} int main() {

double a,b,Eps;

cout<<"请输入积分上下限ab及精度Eps:"; cin>>a>>b>>Eps;

cout<<"所求积分的近似值为:"<return 0; }

15 32


本文来源:https://www.wddqw.com/doc/4b42bef3ae51f01dc281e53a580216fc700a531f.html