数值计算基础实验指导+部分实验源代码+复习指导+三套试题及其答案

1、数值计算基础实验指导书2010 年目录实验一直接法解线性方程组的.3实验二插值方法.12实验三数值积分.6 实验四常微分方程的数值解.8实验五迭代法解线性方程组与非线性方程.10实验一直接法解线性方程组一、实验目的掌握全选主元消去法与高斯-塞德尔法解线性方程组。二、实验内容分别写出Guass 列选主元消去法与追赶法的算法，编写程序上机调试出结果，要求所编程序适用于任何一解线性方程组问题，即能解决这一类问题，而不是某一个问题。实验中以下列数据验证程序的正确性。1、用Guass 列选主元消去法求解方程组é2.52.3-5.1ùéx1ùé3.7

2、49;ê5.39.61.5úêxú=ê3.8úêêë8.11.7-4.3úêúûêë2úx3úûê úêë5.5úû2、用追赶法求解方程组é-20000ùéx1ùé-10ùê1-2000úêxúê0 úêúê2

3、úêúê01-200úêx3ú=ê0 úêúêúêúê001-20úêx4úê0 ú0001-2x50三、实验仪器设备与材料主流微型计算机四、实验原理1、Guass 列选主元消去法对于AX=BA B是上三角矩阵。即：1）、消元过程：将（A|B）进行变换为( |)，其中AæaaLabö æ1aLaböç11121n1÷ç

4、121n1÷ça21a22La2nb2÷ç01La2nb2÷çMMMM÷®çMMMM÷ççèan1an2÷bøLannn÷ççè00÷bøLannn÷k从1 到n-1 a、列选主元选取第k列中绝对值最大元素maxaikk£i£n作为主元。b、换行akjbkc、归一化Ûaij,j=k+1,L,nÛbid、消元akj /akkbk /akk&

5、#222;akj,j=k+1,L,nÞbkaij -aikakj Þaij,i=k+1,L,n;j=k+1,L,nbi -aikbk Þbi,i=k+1,L,n2）、回代过程：由( |)解出x,x,L,x。A Bbn /ann Þxn nnn-11bk -åakjxj Þxk,k=n-1,L,2,1j=k+12、追赶法线性方程组为：æacöæx öæföç11÷ç1 ÷ç1 ÷çb2çç

6、a2c2Ob3a3c3÷çx2 ÷3÷çx ÷÷ç÷çf2 ÷3çf÷ç÷OçOOO÷çM÷=çM÷çççbn-1çOan-1cn-1÷çM÷ç÷çx÷ç÷÷÷÷n-1ç÷Mç f÷

7、1;÷÷n-1èbnan øèxn øèfn ø做LU分解为：æaöæ1böç1÷ç1÷çg2a2ççga÷ç1b2÷÷çOO÷L=ç33÷,R=ç÷ççOO÷çbOç÷çO÷çç÷çOO&

8、#247;÷1n-1÷÷ègnan øè1 ø分解公式：ìgï=aïii(i=2,3,L,n)ía1 =b1 ,ai =bi -gibi-1ï(i=2,3,L,n)iïb=ciïîai(i=1,2,L,n-1)则Ax=fÞLUx=fìLy=fÞíîUx=y回代公式：y=ìf1ï1ïíyïîïia1=fi -giyi-1ai(i=

9、2,3,L,n)ìxn =yníîxi =yi -bixi+1(i=n-1,n-2,L,1)五、实验步骤1、理解并掌握全选主元消去法与高斯-塞德尔迭代法公式；2、画出全选主元消去法与高斯-塞德尔迭代法的流程图3、使用C语言编写出相应的程序并调试验证通过六、实验报告要求1、统一使用武汉科技大学实验报告本书写，实验报告的内容要求有：实验目的、实验内容、程序流程图、源程序、运行结果及实验小结六个部分。2、源程序需打印后粘贴在实验报告册内；3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。七、实验注意事项注意如何定义数据结构以保存矩阵和解以降低算法的复杂性。八、思

10、考题若使用全主元消去法，在编程中应如何记录保存对于未知数的调换。实验二插值方法一、实验目的掌握拉格郎日插值法与牛顿插值法构造插值多项式。二、实验内容分别写出拉格郎日插值法与牛顿插值法的算法，编写程序上机调试出结果，要求所编程序适用于任何一组插值节点，即能解决这一类问题，而不是某一个问题。实验中以下列数据验证程序的正确性。已知下列函数表xi0.561600.562800.564010.56521yi0.827410.826590.825770.82495求x=0.5635时的函数值。三、实验仪器设备与材料主流微型计算机四、实验原理已知n个插值节点的函数值，则可由拉格郎日插值公式与牛顿插值公式构造

11、出插值多项式，从而由该插值多项式求出所要求点的函数值。拉格郎日插值公式与牛顿插值公式如下：1、Lagrange插值公式nLn(x)=l0(x)y0 +l1(x)y1 +.+ln(x)yn =åyklk(x)k=0(x-x)(x-x)L(x-x)(x-x)L(x-x)nx-xl(x)=01k-1 k+1 n=Õ jkk(x -x0)(xk -x1)L(xk -xk-1)(xk -xk+1)L(xk -xn)j=0j¹kxk -xj2、Newton插值公式Nn(x)=f(x0)+fx0,x1(x-x0)+fx0,x1 ,x2(x-x0)(x-x1)+L+fx0,x1

12、,Lxn(x-x0)(x-x1)L(x-xn-1)五、实验步骤1、理解并掌握拉格郎日插值法与牛顿插值法的公式；2、画出拉格郎日插值法与牛顿插值法算法的流程图；3、使用C语言编写出相应的程序并调试验证通过。六、实验报告要求1、统一使用武汉科技大学实验报告本书写，实验报告的内容要求有：实验目的、实验内容、程序流程图、源程序、运行结果及实验小结六个部分。2、源程序需打印后粘贴在实验报告册内；3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。七、实验注意事项Newton插值法在编程时应注意定义何种数据结构以保存差商。八、思考题比较Lagrange插值法与Newton插值法的异同。实验三数值积分

13、一、实验目的掌握复化梯形法与龙贝格法计算定积分。二、实验内容分别写出变步长梯形法与Romberge法计算定积分的算法，编写程序上机调试出结果，要求所编程序适用于任何类型的定积分，即能解决这一类问题，而不是某一个问题。实验中以下列数据验证程序的正确性。求1sinxdx,e£0.00001。ò0 x三、实验仪器设备与材料主流微型计算机四、实验原理通过变步长梯形法与龙贝格法，我们只要知道已知n个求积节点的函数值，则可由相应的公式求出该函数的积分值，从而不需要求该函数的原函数。变步长梯形法与龙贝格法公式如下：1、变步长梯形法å 2n-1 hTn =f(xi)+f(xi+1

14、)i=0hn-1=f(a)+2åf(xi)+f(b)2i=11hn-1T2n = Tn + åf(xi+1/2)22i=0用T2n -Tn£e来控制精度2、龙贝格法1hn-1T2n = Tn + åf(xi+1/2)22i=0Sn =41T2n - Tn33Cn =161S2n - Sn1515641Rn =C2n - Cn6363用R2n -Rn£e来控制精度五、实验步骤1、理解并掌握变步长梯形法与龙贝格法的公式；2、画出变步长梯形法与龙贝格法的流程图3、使用C语言编写出相应的程序并调试验证通过六、实验报告要求1、统一使用武汉科技大学实验报告

15、本书写，实验报告的内容要求有：实验目的、实验内容、程序流程图、源程序、运行结果及实验小结六个部分。2、源程序需打印后粘贴在实验报告册内；3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。七、实验注意事项1sinx在ò0x的定义。dx积分中，被积函数在x=0点函数值为1，对该点在程序设计中应注意对其八、思考题使用复化梯形法与复化Simpson法来计算该问题有何缺点？实验四常微分方程的数值解一、实验目的掌握改进欧拉法与四阶龙格-库塔求解一阶常微分方程的初值问题。二、实验内容分别写出改进欧拉法与四阶龙格-库塔求解的算法，编写程序上机调试出结果，要求所编程序适用于任何一阶常微分方程的

16、数值解问题，即能解决这一类问题，而不是某一个问题。实验中以下列数据验证程序的正确性。ìy¢=-xy2求íîy(0)=2(0£x£5)步长h=0.25。三、实验仪器设备与材料主流微型计算机四、实验原理常微分方程的数值解主要采用“步进式”，即求解过程顺着节点排列次序一步一步向前推进，在单步法中改进欧拉法和四阶龙格-库塔法公式如下：1、改进欧拉法yn+1 =yn +hf(xn,yn)yn+1=yn+hf(x2n,yn)+f(xn+1,yn+1)2、四阶龙格-库塔法ìhïyn+1 =yn +6(k1 +2k2 +2k3 +

17、k4)ïïk1 =ïík2 =ïïf(xn,yn)f(x +hn2h,yn +hk1)2hïk3 =ïïîk4 =f(xn +,yn +k2)22f(xn +h,yn +hk3)五、实验步骤1、理解并掌握改进欧拉法与四阶龙格-库塔法的公式；2、画出改进欧拉法与四阶龙格-库塔法的流程图3、使用C语言编写出相应的程序并调试验证通过六、实验报告要求1、统一使用武汉科技大学实验报告本书写，实验报告的内容要求有：实验目的、实验内容、程序流程图、源程序、运行结果及实验小结六个部分。2、源程序需打印后粘贴在实

18、验报告册内；3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。七、实验注意事项ìy¢=-xy2íîy(0)=2精度的变化八、思考题(0£x£5)的精确解为y =2/(1+x2)，通过调整步长，观察结果的如何对四阶龙格-库塔法进行改进，以保证结果的精度。实验五迭代法解线性方程组与非线性方程一、实验目的掌握高斯-塞德尔迭代法求解线性方程组与牛顿迭代法求方程根。二、实验内容分别写出高斯-塞德尔迭代法与牛顿迭代法的算法，编写程序上机调试出结果，要求所编程序适用于任何一个方程的求根，即能解决这一类问题，而不是某一个问题。实验中以下列数据

19、验证程序的正确性。1、高斯-塞德尔迭代法求解线性方程组é721-2ùéx1ùé4ùê9153-2úêxúê7úêúê2ú=êúê-2-2 115úêx3úê-1úêúêúêúë13213ûëx4ûë0û2、用牛顿迭代法求方程x3 -x-1=

20、0的近似根，e£0.00001，牛顿法的初始值为1。三、实验仪器设备与材料主流微型计算机四、实验原理二分法通过将含根区间逐步二分，从而将根的区间缩小到容许误差范围。牛顿通过迭代的方法逐步趋进于精确解，该两种方法的公式如下：1、高斯-塞德尔迭代法1）判断线性方程组是否主对角占优nåaijj=1j¹i£aii,i=1,2,L,n2）直接分离xi，即nxi =(di -åbijxj)/aii,i=1,2,L,nj=1建立高斯-塞德尔迭代格式为：(k+1)xii-1=(di -åaijxj(k+1)n-åaijxj(k)/aii,i

21、=1,2,L,nj=1j=i+13）取初值迭代求解至所要求的精度为止。2、牛顿法kx=x-f(xk)k+1kf¢(x)五、实验步骤1、理解并掌握二分法与牛顿法的公式；2、画出二分法与牛顿法的流程图3、使用C语言编写出相应的程序并调试验证通过六、实验报告要求1、统一使用武汉科技大学实验报告本书写，实验报告的内容要求有：实验目的、实验内容、程序流程图、源程序、运行结果及实验小结六个部分。2、源程序需打印后粘贴在实验报告册内；3、运行结果以屏幕截图形式保存并打印后粘贴在实验报告册内。七、实验注意事项对于二分法应注意二分后如何判断根的区间，对于牛顿法注意如何确定迭代过程的结束八、思考题若使用

22、牛顿法是发散的，如何对牛顿法进行改进以保证其收敛性。前三个实验的程序代码（C/C+）和运行结果截图#include<iostream>#include<stdlib.h>#include<math.h>usingnamespacestd;Gauss全选主元解方程组的源程序及运行结果class Matrixpublic: Matrix();Matrix();voidSetMatrix(constintn,constdoubleesp1);/构造线性方程组相应的矩阵，n为方程的未知数数目，esp1为要求的精度voidMax(constintr);/全选主元voi

23、dChangeRC(constintr);/根据主元变换矩阵的行或列voidEliminate(constintr);/处理消元工作voidResult()const;/计算方程的解voidCalculate();intGetRank()const;/返回矩阵的行数doubleGetX(consti)const;/确定方程组的第i个解(1<=i<=N)private:/指针a和b分别用于存储方程组的未知数系数和方程"="右边的常数，esp存/储精度double*a,*b,esp;/指针flag用于记录方程组解的顺序int*flag;/以下的结构体用于在全选主元中

24、记录最大主元的位置structcoordinateintrow,column;location;intN;/方程组未知数的数目;intmain()intn;doubleesp1; Matrixmatrix; docout<<"请输入方阵的阶数："cin>>n;if(n<0)n=0;/如何控制非法字符的输入？while(n=0);docout<<"请输入计算精度："cin>>esp1;if(esp1<0)esp1=0;/输入不合法的精度就把精度置0while(esp1=0);cout<<

25、;"输入线性方程组的增广矩阵:n" matrix.SetMatrix(n,esp1);/设置矩阵内的数据matrix.Calculate();/计算方程组的解/输出方程组的解cout<<"nn方程组的解如下:n"for(inti=1;i<=matrix.GetRank();+i)cout<<"X"<<i<<":"<<matrix.GetX(i)<<endl;return0;Matrix:Matrix()/将Matrix类的数据成员初始化a

26、=NULL; b=NULL; flag=NULL; location.row=0; location.column=0; esp=0;N=0;Matrix:Matrix()/释放指针a、b和flag指向的内存空间deletea;deleteb;deleteflag;voidMatrix:SetMatrix(constintn,constdoubleesp1)N=n;esp=esp1;a=newdoubleN*N; b=newdoubleN; flag=newintN;/判断是否成功分配存储区if(a=NULL|b=NULL|flag=NULL)cout<<"分配存储区失败

27、！n"exit(EXIT_FAILURE);/读取线性方程组的增广矩阵for(inti=0;i<N;+i)for(intj=0;j<N;+j) cin>>*(a+i*N+j);cin>>*(b+i);/flag中存储的值对应相应的x值，当方程的解由于列变换交换后，flag中/的值也相应交换，最后用于恢复解的顺序for(i=0;i<N;+i)*(flag+i)=i;voidMatrix:Max(constintr)doublemax=0;for(inti=r;i<N;+i)for(intj=r;j<N;+j)if(max<fa

28、bs(*(a+i*N+j)max=fabs(*(a+i*N+j);/设定最大主元的行、列location.row=i;location.column=j;/最大主元小于输入的精度时，认为方程组无解，退出程序if(max<=esp)cout<<"方程组无解！n"exit(EXIT_FAILURE);voidMatrix:ChangeRC(constintr)doubletemp;/如果最大主元所在的行不在当前行，则进行行变换if(location.row!=r)for(inti=r;i<N;+i)temp=*(a+r*N+i);*(a+r*N+i)=*

29、(a+location.row*N+i);*(a+location.row*N+i)=temp;temp=br; br=blocation.row; blocation.row=temp;/若最大主元所在的列不在当前的r列，则进行列变换if(location.column!=r)for(inti=r;i<N;+i)temp=*(a+i*N+r);*(a+i*N+r)=*(a+i*N+location.column);*(a+i*N+location.column)=temp;/交换flag中的元素来标记方程解的位置变化inttemp1;temp1=*(flag+r);*(flag+r)=

30、*(flag+location.column);*(flag+location.column)=temp1;voidMatrix:Eliminate(constintr)if(fabs(*(a+N*r+r)<=esp)cout<<"方程组无解！n"exit(EXIT_FAILURE);for(inti=r+1;i<N;+i)for(intj=r+1;j<N;+j)(*(a+i*N+j)-=(*(a+i*N+r)*(*(a+r*N+j)/(*(a+r*N+r); (*(b+i)-=(*(b+r)*(*(a+i*N+r)/(*(a+r*N+r);v

31、oidMatrix:Result()constif(fabs(*(a+N*(N-1)+N-1)<=esp)cout<<"方程组无解！n"exit(EXIT_FAILURE);doubletemp;*(b+N-1)/=(*(a+N*(N-1)+N-1);/求出XN-1/依次求出Xi(i=N-2,N-3····，1)for(inti=N-2;i>=0;-i)temp=0;for(intj=i+1;j<N;+j)temp+=(*(a+i*N+j)*(*(b+j);*(b+i)=(*(b+i)-temp)/(*(

32、a+i*N+i);/根据flag中的数据用冒泡排序法恢复方程组解的次序for(i=0;i<N-1;+i)for(intj=0;j<N-i-1;+j)if(*(flag+j)>*(flag+j+1)inttemp1;/交换解的顺序temp=*(b+j);*(b+j)=*(b+j+1);*(b+j+1)=temp;/交换用于标记的元素的顺序temp1=*(flag+j);*(flag+j)=*(flag+j+1);*(flag+j+1)=temp1;voidMatrix:Calculate()/根据矩阵行数重复进行寻找最大主元、变换行或列、消元for(inti=0;i<Ge

33、tRank()-1;+i) Max(i); ChangeRC(i); Eliminate(i);Result();intMatrix:GetRank()constreturnN;/返回矩阵的行数doubleMatrix:GetX(constinti)constreturn*(b+i-1);运行结果追赶法求解方程组的算法：1输入方程组的维数n，将主对角元素b(i)(i=0:n-1)，主对角元素左边的元素a(i)(i=0:n-2)，主对角元素右边的元素c(i)(i=0:n-2)，右端项的元素f(i)(i=0:n-1)2对方程组的系数矩阵作Crout分解, (0)=b(0),对于i=0:n-2,c(

34、i):=(i):=c(i)/b(i),b(i+1):=(i+1):= b(i+1)-a(i)* (i)3.解方程组Ly=fb(0):=y(0):=f(0)/ (0):=f(0)/b(0)对于i=1：n-1，b(i):=y(i):=f(i)-a(i-1)*y(i-1)/b(i)4.解方程组Ux=ya(n-1):=x(n-1):=b(n-1)对于i=n-2:0,a(i)=x(i):=b(i)-c(i)*a(i+1);5输出方程组的解a(0:n-1)#include<iostream>#include<stdlib.h>用追赶法求解方程组的源程序及运行结果usingnames

35、pacestd;class MatrixThrpublic: MatrixThr();MatrixThr();voidSetMatrixThr(constintn);/设置三对角矩阵的数据voidResult();/计算三对角矩阵的解doubleGetX(constinti)const;/取得第i个解,i从1开始intGetN() const;/返回未知数的数目private:intN;/N为未知数的数目/b为矩阵主对角线的元素首地址，a为主对角线左边一斜条元素的首地址，/c为主对角线右边一斜条元素首地址，f为方程组的常数首地址double*a,*b,*c,*f;intmain()Matrix

36、Thrmatrix;intn;docout<<"输入三对角方程组的变量的数目N:"cin>>n;while(n<3);cout<<"请依次输入三对角方程组每行的数据(0元素除外):n" matrix.SetMatrixThr(n);/计算方程组的解matrix.Result();/输出方程组的解cout<<"方程的解如下:n"for(inti=1;i<=matrix.GetN();+i)cout<<"X"<<i<<&qu

37、ot;:"<<matrix.GetX(i)<<endl;return0;MatrixThr:MatrixThr()/初始化相关数据N=0; a=NULL; b=NULL; c=NULL;f=NULL;MatrixThr:MatrixThr()/释放分配的内存空间delete a; delete b; delete c;delete f;voidMatrixThr:SetMatrixThr(constintn)/根据输入的未知数个数设置矩阵的数据N=n;a=new doubleN; b=new doubleN; c=new doubleN; f=newdoubl

38、eN;/若内存分配失败，退出程序if(a=NULL|b=NULL|c=NULL|f=NULL)cout<<"初始化分配存储空间失败！n"exit(EXIT_FAILURE);/依次输入三对角矩阵每行的元素cin>>*b>>*c>>*f;for(inti=1;i<N-1;+i) cin>>*(a+i-1)>>*(b+i)>>*(c+i)>>*(f+i);cin>>*(a+i-1)>>*(b+i)>>*(f+i);voidMatrixThr:

39、Result()/对系数矩阵A作Crout分解for(inti=0;i<N-1;+i)/将U中的存于指针b中，L中的存于指针c中*(c+i)/=(*(b+i);*(b+i+1)-=(*(a+i)*(*(c+i);/解方程组Ly=f，求得的y值存于指针b中*b=(*f)/(*b);for(i=1;i<N;+i)*(b+i)=(*(f+i)-(*(a+i-1)*(*(b+i-1)/(*(b+i);/解方程组Ux=y，求得的x值存于指针a中*(a+N-1)=*(b+N-1);for(i=N-2;i>=0;-i)*(a+i)=*(b+i)-(*(c+i)*(*(a+i+1);intM

40、atrixThr:GetN()constreturnN;doubleMatrixThr:GetX(constinti)constreturn*(a+i-1);运行结果#include<iostream>Lagrange插值法的源程序：#include<stdlib.h>usingnamespacestd;class Lagrangepublic: Lagrange();Lagrange();voidSetLagrange(constintn);/根据用户的输入设置Lagrange类中的插值点数据boolExist(constdoublex,constinti);/检测是

41、否输入了与前i个插值结点横坐标相同的点intGetN()const;/获取插值结点的数目voidCalculate(constdoublea);/计算横坐标a对应的函数值doubleGetResult()const;/返回计算的函数值private:intN;/插值结点的数目double*x,*y,zx,zy;/x、y分别用于存储插值点的数据，zx、zy表示所求的坐标点;intmain()intn=2; LagrangeL; doif(n<2)cout<<"用于模拟曲线的插值点数目N>2"<<endl;cout<<"

42、请输入用于模拟曲线的插值点数目N:"cin>>n;while(n<2);cout<<"请输入各插值点横、纵坐标的数据n" L.SetLagrange(n);doublea;cout<<"请输入所求点的横坐标X:"cin>>a; L.Calculate(a);cout<<"横坐标"<<a<<"-对应的函数值为:"<<L.GetResult()<<endl;return0;Lagrange:Lag

43、range()/初始化相关数据N=0;x=y=NULL;zx=zy=0;Lagrange:Lagrange()/释放分配给指针的内存空间delete x;delete y;boolLagrange:Exist(constdoublea,constinti)/遍历以输入的插值点，看是否重复插入横坐标相同的插值点for(intj=0;j<i;+j)if(a=*(x+j)returntrue;returnfalse;voidLagrange:SetLagrange(constintn)N=n;x=newdoubleN;y=newdoubleN;/判断是否成功为指针分配了内存空间if(x=NUL

44、L|y=NULL)cout<<"分配存储空间失败！"<<endl;exit(EXIT_FAILURE);/输入插值点for(inti=0;i<GetN();+i)cin>>*(x+i)>>*(y+i);/如果不是输入第一个坐标值，则会对输入的横坐标进行合法性检测while(i!=0&&Exist(*(x+i),i)=true)cout<<"输入了重复的横坐标-请重新输入第"<<i+1<<"个插入结点的数据n"cin>>

45、*(x+i)>>*(y+i);voidLagrange:Calculate(constdoublea)zx=a;for(inti=0;i<GetN();+i)/计算插值基函数doubletemp=1;for(intj=0;j<i;+j)temp*=(zx-*(x+j)/(*(x+i)-*(x+j);for(+j;j<GetN();+j)temp*=(zx-*(x+j)/(*(x+i)-*(x+j);zy+=(*(y+i)*temp);intLagrange:GetN()constreturnN;/返回插值点的数目doubleLagrange:GetResult()

46、constreturnzy;/返回求得的函数值Lagrange插值法的运行结果#include<iostream>usingnamespacestd;Newton插值法的源程序：class Matrixpublic: Matrix();Matrix();voidSetMatrix(constintn);/根据用户输入的插值点的数据设置计算结果的矩阵intGetN()const;/返回插值点的数目voidCalculate();/计算差商doubleGetResult(doublex)const;/根据输入的横坐标求函数值，返回运算结果private:double*a,*f;/a和f

47、分别用于存储插值点的横坐标和相应的函数值intN;/记录插值点的数目;intmain()Matrixmatrix;intn;docout<<"输入插值点的数目N："cin>>n;while(n<2);cout<<"输入插值点的数据："<<endl;matrix.SetMatrix(n);matrix.Calculate();doublex;cout<<"输入所求的横坐标："cin>>x;cout<<"所求的函数值为："<

48、<matrix.GetResult(x)<<endl;return1;Matrix:Matrix()/初始化数据a=f=NULL;N=0;Matrix:Matrix()/释放指针a、f指向的内存空间delete a;delete f;voidMatrix:SetMatrix(intn)N=n;/为插值点创建动态数组a=newdoubleN;f=newdoubleN;/输入插值点的数据for(inti=0;i<GetN();+i) cin>>*(a+i)>>*(f+i);intMatrix:GetN()constreturnN;voidMatrix

49、:Calculate()/将差商存储在一个一维数组内for(inti=0;i<GetN();+i)for(intj=GetN()-1;j>i;-j)*(f+j)=(*(f+j)-*(f+j-1)/(*(a+j)-*(a+j-i-1);doubleMatrix:GetResult(doublex)const/利用差商和插值点的横坐标及第一个插值点的纵坐标计算函数值double result=*f;/指针f指向第一个插值点的纵坐标for(inti=1;i<GetN();+i)doubletemp=1;for(intj=0;j<i;+j)temp*=(x-*(a+j);/计算

50、(x-x0)*(x-x1)*···*(x-x(i-1)result+=(*(f+i)*temp;/指针(f+i)指向第i号差商returnresult;Newton插值法的运行结果变步长梯形法求定积分的源程序：#include<iostream>#include<math.h>using namespacestd;doublefunction(constdoublex);/求被积函数的值并返回/accumulate()为求定积分的函数，a、b 分别为积分的上下限,默认精度为0.00001 doubleaccumulate(constdoub

51、lea,constdoubleb,constdoubleeps=0.00001);intmain()doublea,b,eps;/a,b分别为定积分的上限和下限,h为步长,eps 为要求的精度a=0; b=1; eps=0.00001;cout<<"(sinX/X)在0 到1上的积分为:"<<accumulate(a,b)<<endl;return0;doublefunction(constdoublex)if(x=0)return1;/x为0 时函数值为1 return(sin(x)/x);doubleaccumulate(constd

52、oublea,constdoubleb,constdouble eps)intn=1;doubleh=b-a;/h为步长doubleT1=h*(function(a)+function(b)/2;doubleT2=T1/2+h*function(a+h/2)/2;/计算结果不满足精度，则继续while(fabs(T2-T1)>eps)n*=2; h/=2;/步长折半T1=T2;/利用T1 计算T2 doubletemp=0;for(int i=1;i<=n;+i)temp+=function(a+(i-1.0/2)*h); T2=T1/2+h*temp/2;returnT2;/返回

53、步长h最合适时定积分的近似值1sinx变步长梯形法求定积分ò0的运行结果xRomberg算法求定积分的源程序：#include<iostream>#include<math.h>usingnamespacestd;/求被积函数的值并返回doublefun(constdoublex);/Romberge()为求定积分的函数，a、b分别为积分的上下限,默认精度为0.00001 doubleRomberg(constdoublea,constdoubleb,constdoubleeps=0.0001);/函数Tm()为T-数表的计算公式doubleTm(constd

54、oubleT1,constdoubleT2,constintm);intmain()doublea,b,eps;/a,b分别为定积分的上限和下限,h为步长,eps 为要求的精度a=0; b=1; eps=0.00001;cout<<"所求积分为:"<<Romberg(a,b,eps)<<endl;return0;doublefun(constdoublex)if(x=0)return1;/x为0 时函数值为1 return(sin(x)/x);/返回被积函数的值doubleTm(constdouble T1,constdoubleT2,c

55、onstintm)/根据T1、T2 和m计算出下一个数并返回return(pow(4,m)*T1-T2)/(pow(4,m)-1);doubleRomberg(constdoublea,constdoubleb,constdoubleeps)doubleT10;/用于存储T-数表某一行的数据int l=1;T0=(b-a)*(fun(a)+fun(b)/2;doubleT1,T2;intflag;/flag用于标记十分求得符合精度的值doflag=0;intn=pow(2,l-1);/n表示区间等分数doubleh=(b-a)/n;/h为步长doubletemp=0;for(int i=0;i<n;+i)temp+=fun(a+(i+1.0/2)*h); T1=T0/2+h*temp/2