计算方法A上机实验报告

. / 《计算方法》上机实验报告班级：XXXXXX小组成员：XXXXXXXXXXXXXXXXXXXXXXXXXXXX任课教师：XXX二〇一八年五月二十五日前言通过进行多次的上机实验，我们结合课本上的内容以及老师对我们的指导，能够较为熟练地掌握Newton 迭代法、Jacobi 迭代法、Gauss-Seidel 迭代法、Newton 插值法、Lagrange 插值法和Gauss 求积公式等六种算法的原理和使用方法，并参考课本例题进行了MATLAB 程序的编写。

以下为本次上机实验报告，按照实验内容共分为六部分。

实验一：一、实验名称及题目： Newton 迭代法例2.7(P38):应用Newton 迭代法求在附近的数值解，并使其满足.二、解题思路：设'x 是0)(=x f 的根,选取0x 作为'x 初始近似值,过点())(,00x f x 做曲线)(x f y =的切线L ,L 的方程为))((')(000x x x f x f y -+=,求出L 与x 轴交点的横坐标)(')(0001x f x f x x -=,称1x 为'x 的一次近似值,过点))(,(11x f x 做曲线)(x f y =的切线,求该切线与x 轴的横坐标)(')(1112x f x f x x -=称2x 为'x 的二次近似值,重复以上过程,得'x 的近似值序列{}n x ,把)(')(1n n n n x f x f x x -=+称为'x 的1+n 次近似值，这种求解方法就是牛顿迭代法。

三、Matlab 程序代码：function newton_iteration(x0,tol) syms z %定义自变量 format long %定义精度 f=z*z*z-z-1;f1=diff(f);%求导 y=subs(f,z,x0);y1=subs(f1,z,x0);%向函数中代值 x1=x0-y/y1; k=1;while abs(x1-x0)>=tol x0=x1;y=subs(f,z,x0); y1=subs(f1,z,x0); x1=x0-y/y1;k=k+1; endx=double(x1) K四、运行结果：实验二：一、实验名称及题目：Jacobi 迭代法例3.7（P74）：试利用Jacobi 迭代公式求解方程组要求数值解为方程组的精确解. 二、解题思路：首先将方程组中的系数矩阵A 分解成三部分，即：U D L A ++=，D 为对角阵，L 为下三角矩阵，U 为上三角矩阵。

计算方法A上机报告学院（系）：电气工程学院学生姓名：陶然学号：授课老师：完成日期：2019年12月03日西安交通大学Xi'an Jiaotong University目录1 QR分解法求解线性方程组 (2)1.1 算法原理 (2)1.1.1 基于吉文斯变换的QR分解 (2)1.1.2 基于豪斯霍尔德变换的QR分解 (3)1.2 程序流程图 (4)1.2.1 基于吉文斯变换的QR分解流程图 (4)1.2.2 基于豪斯霍尔德变换的QR分解流程图 (5)1.3 程序使用说明 (5)1.3.1 基于吉文斯变换的QR分解程序说明 (5)1.3.2 基于豪斯霍尔德变换的QR分解程序说明 (7)1.4 算例计算结果 (8)2 共轭梯度法求解线性方程组 (10)2.1 算法原理 (10)2.2 程序流程图 (10)2.3 程序使用说明 (11)2.4 算例计算结果 (12)3 三次样条插值 (14)3.1 算法原理 (14)3.2 程序流程图 (16)3.3 程序使用说明 (17)3.4 算例计算结果 (19)4 龙贝格积分 (21)4.1 算法原理 (21)4.2 程序流程图 (22)4.3 程序使用说明 (23)4.4 算例计算结果 (24)结论 (26)1 QR 分解法求解线性方程组1.1 算法原理矩阵的QR 分解是指，可以将矩阵A 分解为一个正交矩阵Q 和一个上三角矩阵R 的乘积，实际中，QR 分解经常被用来解决线性最小二乘问题，分解情况如图1.1所示。

=⨯图1.1 QR 分解示意图本次上机学习主要进行了两个最基本的正交变换—吉文斯（Givens ）变换和豪斯霍尔德（Householder ）变换，并由此导出矩阵的QR 分解以及求解线性方程组的的方法。

1.1.1 基于吉文斯变换的QR 分解吉文斯矩阵也称初等旋转阵，如式（1.1）所示，它把n 阶单位矩阵I 的第,i j 行的对角元改为c ，将第i 行第j 列的元素改为s ，第j 行第i 列的元素改为s −后形成的矩阵。

《计算方法》实验报告指导教师：学院：班级：团队成员：一、题目例2.7应用Newton 迭代法求方程210x x --=在1x =附近的数值解k x ，并使其满足8110k k x x ---<原理：在方程()0f x =解的隔离区间[],a b 上选取合适的迭代初值0x ，过曲线()y f x =的点()()00x f x ，引切线()()()1000:'l y f x f x x x =+-其与x 轴相交于点：()()0100 'f x x x f x =-，进一步,过曲线()y f x =的点()()11x f x , 引切线()()()2111: 'l y f x f x x x =+-其与x 轴相交于点：()()1211 'f x x x f x =-如此循环往复，可得一列逼近方程()0f x =精确解*x 的点01k x x x ,,,,,其一般表达式为：()()111 'k k k k f x x x f x ---=-该公式所表述的求解方法称为Newton 迭代法或切线法。

程序：function y=f(x)%定义原函数y=x^3-x-1;endfunction y1=f1(x0)%求导函数在x0点的值syms x;t=diff(f(x),x);y1=subs(t,x,x0);endfunction newton_iteration(x0,tol)%输入初始迭代点x0及精度tol x1=x0-f(x0)/f1(x0);k=1;%调用f函数和f1函数while abs(x1-x0)>=tolx0=x1;x1=x0-f(x0)/f1(x0);k=k+1;endfprintf('满足精度要求的数值为x(%d)=%1.16g\n',k,x1);fprintf('迭代次数为k=%d\n',k);end结果：二、题目例3.7试利用Jacobi 迭代公式求解方程组1234451111101112115181111034x x x x ----⎛⎫⎛⎫⎛⎫ ⎪ ⎪ ⎪--- ⎪ ⎪ ⎪= ⎪ ⎪ ⎪--- ⎪ ⎪ ⎪---⎝⎭⎝⎭⎝⎭ 要求数值解()k X 满足()4210k --≤X X ，其中T (1,2,3,4)=X 为方程组的精确解。

计算方法上机作业报告姓名：李小盼班级：计算方法B3班学号：6230489477专业：热能工程2016年《计算方法B》上机题目一.计算机语言要求使用语言不做限制，可以使用C、C++、FORTRAN、VC、VB、C#、Matlab、PHP、JavaScript等语言。

二.实习报告内容上机题目完成后，必须交一份上机报告。

上机报告中应对每道题目包括：（1）上机题目内容;（2）详细说明实现的思想、算法依据、算法实现的结构;（3）详细完整的源程序，并附相关的注释说明;（4）给出必要的计算结果（若数据量大，则只需列出主要的数据内容），并对结果进行分析;（5）对上机中出现的问题进行分析总结；三.实习报告要求1.提供一份完整的上机报告的电子文档；然后再提供一份与电子文档对应的纸质上机报告。

2.电子文档中提供上机过程中的所有源程序、输出数据，以及可以运行的文件。

如果程序的运行环境特殊，请附上运行程序所需要的软件环境。

3.上机报告严禁抄袭，如发现有抄袭现象，所有涉及抄袭的上机报告将被以作弊处理，并按零分处理,不再另行通知。

4.上机报告电子版一、二、三班发送到邮箱:xjtujsff@,上机作业纸面作业请送到:理科楼338。

上机实习题目1.某通信公司在一次施工中，需要在水面宽度为20米的河沟底部沿直线走向铺设一条沟底光缆。

在铺设光缆之前需要对沟底的地形进行初步探测，从而估计所需光缆的长度，为工程预算提供依据。

已探测到一组等分点位置的深度数据(单位:米)如下表所示：(1)请用合适的曲线拟合所测数据点；(2)估算所需光缆长度的近似值，并作出铺设河底光缆的曲线图；解：1.1 实现思想选用曲线拟合数据点时，一方面要满足插值条件，即保证所得曲线经过以上所有数据点，另一方面也要保证曲线具有足够的“光滑性”，故这里采用三次样条插值法构造插值函数。

为构成封闭方程组，边界条件选取自然三次样条，以获得三次样条插值的关键参数M0和M n。

1.2 算法的依据以及结构三次样条插值具有较好的光滑性，可以用于对数据点的光滑拟合。

计算方法A上机实验报告姓名：苏福班级：硕4020 学号：3114161019一、上机练习目的1）复习和巩固数值计算方法的基本数学模型，全面掌握运用计算机进行数值计算的具体过程及相关问题。

2）利用计算机语言独立编写、调试数值计算方法程序，培养学生利用计算机和所学理论知识分析解决实际问题的能力。

二、上机练习任务1）利用计算机语言编写并调试一系列数值方法计算通用程序，并能正确计算给定题目，掌握调试技能。

2）掌握文件使用编程技能，如文件的各类操作，数据格式设计、通用程序运行过程中文件输入输出运行方式设计等。

3）写出上机练习报告。

三、上机题目1. 共轭梯度法求解线性方程组。

（第三章）2. 三次样条插值（第四章）3. 龙贝格积分（第六章）4. 四阶龙格－库塔法求解常微分方程的初值问题四、上机报告题目1：共轭梯度法求解线性方程组1．算法原理共轭梯度法是把求解线性方程组的问题转化为求解一个与之等价的二次函数极小值的问题。

从任意给定的初始点出发，沿一组关于矩阵A共轭的方向进行线性搜索，在无舍入误差的假定下，最多迭代n 次（其中n 为矩阵A 的阶数），就可求得二次函数的极小值，也就求得了线性方程组Ax b =的解。

定理：设A 是n 阶对称正定矩阵，则x *是方程组Ax b =的解得充分必要条件是x *是二次函数1()2TT f x x Ax b x =-的极小点，即 ()()min nx R Ax b f x f x **∈=⇔=共轭梯度法的计算公式：(0)(0)(0)()()()()(1)()()(1)(1)(1)()()()(1)(1)()k T k k k T k k k k k k k k T k k k T k k k k k d r b Ax r d d Ad xx d r b Ax r Ad d Ad d r d ααββ++++++⎧==-⎪⎪=⎪⎪=+⎪⎨=-⎪⎪⎪=-⎪⎪=+⎩2. 程序框图（1）编写共轭梯度法求解对称正定矩阵的线性方程组见附录（myge.m）：function x=myge(A,b)输入对称正定矩阵及对应的列向量，初始向量设为0，精度取为810 。

1.拉格朗日插值多项式，用于离散数据的拟合#include <stdio.h>#include <conio.h>#include <alloc.h>float lagrange(float *x,float *y,float xx,int n) /*拉格朗日插值算法*/{ int i,j;float *a,yy=0.0; /*a作为临时变量，记录拉格朗日插值多项式*/a=(float *)malloc(n*sizeof(float));for(i=0;i<=n-1;i++){ a[i]=y[i];for(j=0;j<=n-1;j++)if(j!=i) a[i]*=(xx-x[j])/(x[i]-x[j]);yy+=a[i];}free(a);return yy;}main(){ int i,n;float x[20],y[20],xx,yy;printf("Input n:");scanf("%d",&n);if(n>=20) {printf("Error!The value of n must in (0,20)."); getch();return 1;} if(n<=0) {printf("Error! The value of n must in (0,20)."); getch(); return 1;} for(i=0;i<=n-1;i++){ printf("x[%d]:",i);scanf("%f",&x[i]);}printf("\n");for(i=0;i<=n-1;i++){ printf("y[%d]:",i);scanf("%f",&y[i]);}printf("\n");printf("Input xx:");scanf("%f",&xx);yy=lagrange(x,y,xx,n);printf("x=%f,y=%f\n",xx,yy);getch();}2.牛顿插值多项式，用于离散数据的拟合#include <stdio.h>#include <conio.h>#include <alloc.h>void difference(float *x,float *y,int n){ float *f;int k,i;f=(float *)malloc(n*sizeof(float));for(k=1;k<=n;k++){ f[0]=y[k];for(i=0;i<k;i++)f[i+1]=(f[i]-y[i])/(x[k]-x[i]);y[k]=f[k];}return;}main(){ int i,n;float x[20],y[20],xx,yy;printf("Input n:");scanf("%d",&n);if(n>=20) {printf("Error! The value of n must in (0,20)."); getch(); return 1;} if(n<=0) {printf("Error! The value of n must in (0,20).");getch(); return 1;} for(i=0;i<=n-1;i++){ printf("x[%d]:",i);scanf("%f",&x[i]);}printf("\n");for(i=0;i<=n-1;i++){ printf("y[%d]:",i);scanf("%f",&y[i]);}printf("\n");difference(x,(float *)y,n);printf("Input xx:");scanf("%f",&xx);yy=y[20];for(i=n-1;i>=0;i--) yy=yy*(xx-x[i])+y[i];printf("NewtonInter(%f)=%f",xx,yy);getch();}3.高斯列主元消去法,求解其次线性方程组第一种#include<stdio.h>#include <math.h>#define N 20int main(){ int n,i,j,k;int mi,tmp,mx;float a[N][N],b[N],x[N];printf("\nInput n:");scanf("%d",&n);if(n>N){ printf("The input n should in(0,N)!\n");getch();return 1;}if(n<=0){ printf("The input n should in(0,N)!\n");getch();return 1;}printf("Now input a(i,j),i,j=0...%d:\n",n-1); for(i=0;i<n;i++){ for(j=0;j<n;j++)scanf("%f",&a[i][j]);}printf("Now input b(i),i,j=0...%d:\n",n-1);for(i=0;i<n;i++)scanf("%f",&b[i]);for(i=0;i<n-2;i++){ for(j=i+1,mi=i,mx=fabs(a[i][j]);j<n-1;j++) if(fabs(a[j][i])>mx){ mi=j;mx=fabs(a[j][i]);}if(i<mi){ tmp=b[i];b[i]=b[mi];b[mi]=tmp;for(j=i;j<n;j++){ tmp=a[i][j];a[i][j]=a[mi][j];a[mi][j]=tmp;}}for(j=i+1;j<n;j++){ tmp=-a[j][i]/a[i][i];b[j]+=b[i]*tmp;for(k=i;k<n;k++)a[j][k]+=a[i][k]*tmp;}}x[n-1]=b[n-1]/a[n-1][n-1];for(i=n-2;i>=0;i--){ x[i]=b[i];for(j=i+1;j<n;j++)x[i]-=a[i][j]*x[j];x[i]/=a[i][i];}for(i=0;i<n;i++)printf("Answer:\n x[%d]=%f\n",i,x[i]); getch();return 0;}第二种#include<math.h>#include<stdio.h>#define NUMBER 20#define Esc 0x1b#define Enter 0x0dfloat A[NUMBER][NUMBER+1] ,ark;int flag,n;exchange(int r,int k);float max(int k);message();main(){float x[NUMBER];int r,k,i,j;char celect;clrscr();printf("\n\nUse Gauss.");printf("\n\n1.Jie please press Enter."); printf("\n\n2.Exit press Esc.");celect=getch();if(celect==Esc)exit(0);printf("\n\n input n=");scanf("%d",&n);printf(" \n\nInput matrix A and B:"); for(i=1;i<=n;i++){printf("\n\nInput a%d1--a%d%d and b%d:",i,i,n,i);for(j=1;j<=n+1;j++) scanf("%f",&A[i][j]);}for(k=1;k<=n-1;k++){ark=max(k);if(ark==0){printf("\n\nIt's wrong!");message();}else if(flag!=k)exchange(flag,k);for(i=k+1;i<=n;i++)for(j=k+1;j<=n+1;j++)A[i][j]=A[i][j]-A[k][j]*A[i][k]/A[k][k];}x[n]=A[n][n+1]/A[n][n];for( k=n-1;k>=1;k--){float me=0;for(j=k+1;j<=n;j++){me=me+A[k][j]*x[j];}x[k]=(A[k][n+1]-me)/A[k][k];}for(i=1;i<=n;i++){printf(" \n\nx%d=%f",i,x[i]);}message();}exchange(int r,int k){int i;for(i=1;i<=n+1;i++)A[0][i]=A[r][i];for(i=1;i<=n+1;i++)A[r][i]=A[k][i];for(i=1;i<=n+1;i++)A[k][i]=A[0][i];}float max(int k){int i;float temp=0;for(i=k;i<=n;i++)if(fabs(A[i][k])>temp){temp=fabs(A[i][k]);flag=i;}return temp;}message(){printf("\n\n Go on Enter ,Exit press Esc!");switch(getch()){case Enter: main();case Esc: exit(0);default:{printf("\n\nInput error!");message();} }}4.龙贝格求积公式，求解定积分#include<stdio.h>#include<math.h>#define f(x) (sin(x)/x)#define N 20#define MAX 20#define a 2#define b 4#define e 0.00001float LBG(float p,float q,int n){ int i;float sum=0,h=(q-p)/n;for (i=1;i<n;i++)sum+=f(p+i*h);sum+=(f(p)+f(q))/2;return(h*sum);}void main(){ int i;int n=N,m=0;float T[MAX+1][2];T[0][1]=LBG(a,b,n);n*=2;for(m=1;m<MAX;m++){ for(i=0;i<m;i++)T[i][0]=T[i][1];T[0][1]=LBG(a,b,n);n*=2;for(i=1;i<=m;i++)T[i][1]=T[i-1][1]+(T[i-1][1]-T[i-1][0])/(pow(2,2*m)-1);if((T[m-1][1]<T[m][1]+e)&&(T[m-1][1]>T[m][1]-e)){ printf("Answer=%f\n",T[m][1]); getch();return ;}}}5.牛顿迭代公式，求方程的近似解#include<stdio.h>#include<math.h>#include<conio.h>#define N 100#define PS 1e-5#define TA 1e-5float Newton(float (*f)(float),float(*f1)(float),float x0 ) { float x1,d=0;int k=0;do{ x1= x0-f(x0)/f1(x0);if((k++>N)||(fabs(f1(x1))<PS)){ printf("\nFailed!");getch();exit();}d=(fabs(x1)<1?x1-x0:(x1-x0)/x1);x0=x1;printf("x(%d)=%f\n",k,x0);}while((fabs(d))>PS&&fabs(f(x1))>TA) ;return x1;}float f(float x){ return x*x*x+x*x-3*x-3; }float f1(float x){ return 3.0*x*x+2*x-3; }void main(){ float f(float);float f1(float);float x0,y0;printf("Input x0: ");scanf("%f",&x0);printf("x(0)=%f\n",x0);y0=Newton(f,f1,x0);printf("\nThe root is x=%f\n",y0); getch();}。

计算方法上上机实习报告在本次计算方法的上机实习中，我深入体验了数值计算的魅力和挑战，通过实际操作和实践，对计算方法有了更深刻的理解和认识。

实习的目的在于将课堂上学到的理论知识运用到实际的计算中，熟悉各种数值算法的实现过程，提高编程能力和解决实际问题的能力。

我们使用了具体编程语言和软件名称进行编程和计算。

在实习过程中，我首先接触到的是数值逼近的相关内容。

通过多项式插值和曲线拟合的练习，我明白了如何用简单的函数去近似复杂的曲线。

例如，拉格朗日插值法和牛顿插值法让我能够根据给定的离散数据点构建出一个连续的函数，从而对未知点进行预测。

在实际操作中，我需要仔细处理数据的输入和输出，以及算法中的细节，如边界条件和误差控制。

数值积分是另一个重要的部分。

通过梯形公式和辛普森公式，我学会了如何对给定的函数进行数值积分。

在编程实现时，要合理地选择积分区间和步长，以达到所需的精度。

同时，我也了解到了数值积分方法的误差来源和误差估计方法，这对于评估计算结果的可靠性非常重要。

线性方程组的求解是计算方法中的核心内容之一。

我分别使用了高斯消元法和迭代法（如雅克比迭代法和高斯赛德尔迭代法）来求解线性方程组。

在实际编程中，我深刻体会到了算法的效率和稳定性的重要性。

对于大规模的线性方程组，选择合适的算法可以大大提高计算速度和精度。

在非线性方程求根方面，我运用了二分法、牛顿法和割线法等方法。

这些方法各有特点，二分法简单但收敛速度较慢，牛顿法收敛速度快但需要计算导数。

在实际应用中，需要根据方程的特点和求解的要求选择合适的方法。

在实习中，我也遇到了不少问题和挑战。

首先是编程中的错误，如语法错误、逻辑错误等，这需要我耐心地调试和修改代码。

其次，对于一些复杂的算法，理解其原理和实现细节并不容易，需要反复查阅资料和思考。

还有就是数值计算中的误差问题，有时候由于误差的积累，导致计算结果与预期相差较大，需要通过调整算法参数或者采用更精确的算法来解决。

实习题11. 用两种不同的顺序计算644834.11000012≈∑=-n n，试分析其误差的变化解：从n=1开始累加，n 逐步增大，直到n=10000；从n=10000开始累加，n 逐步减小，直至1。

算法1的C 语言程序如下： #include<stdio.h> #include<math.h> void main() { float n=0.0; int i; for(i=1;i<=10000;i++) { n=n+1.0/(i*i); } printf("%-100f",n); printf("\n"); float m=0.0; int j; for(j=10000;j>=1;j--) { m=m+1.0/(j*j); } printf("%-7f",m); printf("\n"); }运行后结果如下：结论： 4.设∑=-=Nj N j S 2211，已知其精确值为)11123(21+--N N 。

1）编制按从大到小的顺序计算N S 的程序； 2）编制按从小到大的顺序计算N S 的程序；3）按2种顺序分别计算30000100001000,,S S S ，并指出有效位数。

解：1）从大到小的C语言算法如下：#include<stdio.h>#include<math.h>#include<iostream>using namespace std;void main(){float n=0.0;int i;int N;cout<<"Please input N"<<endl;cin>>N;for(i=N;i>1;i--){n=n+1.0/(i*i-1);N=N-1;}printf("%-100f",n);printf("\n");}执行后结果为：N=2时，运行结果为：N=3时，运行结果为：N=100时，运行结果为：N=4000时，运行结果为：2）从小到大的C语言算法如下：#include<stdio.h>#include<math.h>#include<iostream>using namespace std;void main(){float n=0.0;int i;int N;cout<<"Please input N"<<endl;cin>>N;for(i=2;i<=N;i++){n=n+1.0/(i*i-1);}printf("%-100f",n);printf("\n");}执行后结果为：N=2时，运行结果为：N=3时，运行结果为：N=100时，运行结果为：N=4000时，运行结果为：结论：通过比较可知：N 的值较小时两种算法的运算结果相差不大，但随着N 的逐渐增大，两种算法的运行结果相差越来越大。

计算方法与实习上机实验报告一、引言本文旨在介绍和展示我们在“计算方法”课程中的实习上机实验环节所完成的一些关键任务和所取得的成果。

该实验课程的目标是让我们更深入地理解和应用各种计算方法，并在实际操作中提高我们的编程和问题解决能力。

二、实验内容与目标实验的主要内容是利用各种计算方法解决实际数学问题。

我们被要求使用编程语言（如Python或Java）来实现和解决这些问题。

这些问题包括使用牛顿法求解平方根，使用蒙特卡洛方法计算圆周率，以及使用最优化方法求解函数的最小值等。

实验的目标不仅是让我们掌握计算方法的基本理论，更是要让我们能够在实际操作中运用这些方法。

我们需要在实习过程中，通过与同伴们合作，共同解决问题，提高我们的团队合作能力和问题解决能力。

三、实验过程与问题解决策略在实验过程中，我们遇到了许多问题，如编程错误、理解困难和时间压力等。

我们通过相互讨论、查阅资料和寻求教师帮助等方式，成功地解决了这些问题。

例如，在实现牛顿法求解平方根时，我们一开始对导数的计算和理解出现了一些错误。

但我们通过查阅相关资料和讨论，最终理解了导数的正确计算方法，并成功地实现了牛顿法。

四、实验结果与结论通过这次实习上机实验，我们不仅深入理解了计算方法的基本理论，还在实际操作中提高了我们的编程和问题解决能力。

我们的成果包括编写出了能有效求解平方根、计算圆周率和求解函数最小值的程序。

这次实习上机实验非常成功。

我们的团队不仅在理论学习和实践操作上取得了显著的进步，还在团队合作和问题解决方面积累了宝贵的经验。

这次实验使我们对计算方法有了更深的理解和认识，也提高了我们的编程技能和解决问题的能力。

五、反思与展望回顾这次实验，我们意识到在实验过程中，我们需要更好地管理我们的时间和压力。

在解决问题时，我们需要更有效地利用我们的知识和资源。

在未来，我们希望能够更加熟练地运用计算方法，并能够更有效地解决问题。

我们也希望能够将所学的计算方法应用到更广泛的领域中，如数据分析、科学研究和工业生产等。

计算方法第三次上机实习报告实验报告课程名称： 计算方法 指导老师： 郑太英 成绩：__________________实验名称： 第三次上机作业 实验类型： matlab 同组学生姓名： 一、实验目的和要求（必填） 二、实验内容和原理（必填）三、主要仪器设备（必填） 四、操作方法和实验步骤五、实验数据记录和处理 六、实验结果与分析（必填） 七、讨论、心得一、实验目的用龙贝格算法计算积分，要求误差不超过ε=×105二、实验原理龙贝格算法是由递推算法得来的。

由梯形公式得出辛普森公式得出柯特斯公式最后得到龙贝格公式。

设将求积区间[a ，b]分为n 个等分，则一共有n+1个等分点,k x a kh =+,0,1,b ah k n-==，n 。

这里用n T 表示复化梯形法求得的积分值，其下标n 表示等分数。

先考察下一个字段[1,k k x x +]，其中点()11212k k k x x x ++=+，在该子段上二分前后两个积分值 ()()112k k hT f x f x +=+⎡⎤⎣⎦ ()()21124k k k h T f x f x f x ++⎡⎤⎛⎫=++⎢⎥ ⎪⎢⎥⎝⎭⎣⎦显然有下列关系 2112122k h T T f x +⎛⎫=+ ⎪⎝⎭将关系式关于k 从0到n-1累加求和，即可得递推公式12102122n n n k k h T T f x -+=⎛⎫=+ ⎪⎝⎭∑需要强调指出的是，上式中的b a h n -=代表二分前的步长，而1212k x a k h +⎛⎫=++ ⎪⎝⎭根据梯形法的误差公式，积分值n T 的截断误差大致与2h 成正比，因此步长减半后误差将减至四分之一，即有21114n n T T -≈- 将上式移项整理，知 2211()3n n n T T T -≈-按上式，积分值2n T 的误差大致等于21()3n n T T -，如果用这个误差值作为2n T 的一种补偿，可以期望，所得的()222141333n n n n n T T T T T T =+-=-应当是更好的结果。