数值计算方法实验报告(共21页).doc

精选优质文档-倾情为你奉上 本科实验报告课程名称： 计算机数值方法 实验项目： 计算机数值方法实验 实验地点： 虎峪校区致远楼B401 专业班级： 软件学院1217班 学号： xxxx 学生姓名： xxx 指导教师： xxx 2014 年 5 月 21 日太原理工大学学生实验报告学院名称软件学院专业班级1217班 学号xxxx 学生姓名xx 实验日期2014.05.21成绩课程名称数值计算方法实验题目实验一 方程求解一、实验目的和要求熟悉使用、迭代法、牛顿法、割线法等方法对给定的方程进行根的求解。选择上述方法中的两种方法求方程：二分法f(x)=x3+4x2-10=0在1,2内的一个实根，且要求满足精度|x*-xn|<0.5×10-5二、主要设备 笔记本 HP ProBook 6470b 一台 编译软件：VC+6.0三、实验内容和原理 函数f(x)在区间（x，y）上连续，先在区间（x，y）确定a与b，若f(a)，f(b)异号，说明在区间(a，b)内存在零点，然后求f(a+b)/2。假设F(a)<0,F(b)>0,a<b， 如果f(a+b)/2=0，该点即为零点； 如果f(a+b)/2<0，则区间（(a+b)/2，b)内存在零点，(a+b)/2a； 如果f(a+b)/2>0，则区间(a,(a+b)/2)内存在零点，(a+b)/2b；返回重新循环，不断接近零点。通过每次把f(x)的零点所在区间收缩一半的方法，使区间内的两个端点逐步逼近函数零点，最终求得零点近似值。四、操作方法与实验步骤 1. 二分法：#include<stdio.h>#include<stdlib.h>#include<math.h>int main() double a=1.0, b=2.0; double x,s; printf(" AnttBnttF(Xn)n"); while(1) x=(a+b)/2; s=pow(x,3)+4*x*x-10; if (-0. < s && s < 0.) break; else if(s < 0) a=x; else if(s > 0) b=x; printf("%ft%ft%fn",a,b,s); printf("X的值为：%fn",x); printf("误差：t%fn",s); return 0;2. 割线法：#include"stdio.h"#include"math.h"int main() float c,a=1.0,b=2.0; printf("每次得到的X的近似值：n"); while(1) c=b-(b*b*b+4*b*b-10)*(b-a)/(b*b*b+4*b*b-(a*a*a+4*a*a); if(fabs(b-c)<0.5*0.00001)break; b=c; printf("%fn",b); printf("X的值为：%fn",c); 五、实验结果与分析 二分法 割线法 分析： 由程序知，使用二分法和割线法均能计算出方程的根，但利用割线法要比二分法计算的次数少，并且能够较早的达到精度要求。相比之下，割线法程序代码量较少，精简明了。六、讨论、心得 本次数值计算方法程序设计实验从习题练习中跳脱出来，直接面对实用性较强的程序代码编写。效果很好，不仅加深对二分法、割线法的理解，还加强了实际用运能力。将理论知识成功地转化成实践结果。实验地点虎峪校区致远楼B401指导教师xx太原理工大学学生实验报告学院名称软件学院专业班级1217班 学号xxxx 学生姓名xx实验日期2014.05.28成绩课程名称数值计算方法实验题目实验二 线性方程组的直接解法一、实验目的和要求合理利用Gauss消元法、LU分解法、追赶法求解下列方程组： （n=5,10,100,）二、主要设备 笔记本 HP ProBook 6470b 一台 编译软件：VC+6.0三、实验内容和原理高斯消元法：将原方程组化为三角形方阵的方程组：lik=aik/akk aij= aij- lik* akj （ k=1,2,n-1 i=k+1,k+2, ,n j=k+1,k+2, ,n+1 ）由回代过程求得原方程组的解： xn= ann+1/ ann xk=( akn+1-akj xj)/ akk LU分解法：将系数矩阵A转化为A=L*U，L为单位下三角矩阵，U为普通上三角矩阵，然后通过解方程组l*y=b，u*x=y,来求解x。四、操作方法与实验步骤1. 完全主元素消元法：#include<stdio.h>#include<iostream.h>#include"math.h"float a100101;float x10;int N; void shuchu()for(int i=1;i<=N;i+)for(int j=1;j<=N+1;j+)cout<<aij<<" "<<" "cout<<endl;void shuru()cout<<"请输入矩阵阶数:"<<endl;cin>>N;cout<<"请输入矩阵各项:"<<endl;for(int i=1;i<=N;i+)for(int j=1;j<=N+1;j+)cin>>aij;cout<<endl;void main()int z10;int maxi,maxj;shuru();for(int i=1;i<=N;i+)zi=i;for(int k=1;k<N;k+)maxi=k;maxj=k;float maxv=abs(akk);for(i=k;i<=N;i+)for(int j=k;j<=N;j+)if(abs(aij)>maxv)maxv=abs(aij);maxi=i;maxj=j;if(maxi!=k) for(int j=1;j<=N+1;j+)float t=akj;akj=amaxij;amaxij=t;if(maxj!=k) for(i=1;i<=N;i+)float t=aik;aik=aimaxj;aimaxj=t;int t=zk;zk=zmaxj;zmaxj=t; for(int i=k+1;i<=N;i+) float l=aik/akk;for(int j=k;j<=N+1;j+)aij+=-l*akj;for(i=N;i>0;i-)float s=0;for(int j=i+1;j<=N;j+)s+=aij*xzj;xzi=(aiN+1-s)/aii;cout<<"完全主元素消去法之后的矩阵为："<<endl;shuchu(); for(i=1;i<=N;i+) cout<<"x"<<i<<"="<<xi<<endl;2. 列主元素消元法：#include"stdio.h"int main() float a34=1,2,3,14,0,1,2,8,2,4,1,13;float x3; float sum=0; int k,i,j; for(k=0;k<2;k+) for(i=k+1;i<3;i+) for(j=k+1;j<4;j+)aij=aij-aik/akk*akj; for(i=0;i<3;i+) for(j=0;j<4;j+)printf("a%d%d=%f ",i,j,aij); printf("n");x2=a23/a22;for(k=1;k>=0;k-)sum=0;for(j=k+1;j<3;j+)sum+=akj*xj; xk=(ak3-sum)/akk; for(i=0;i<3;i+)printf ("x%d=%fn",i+1,xi);printf("n");3. LU分解法：#include <stdio.h>#include <math.h> #define L 30double aLL,bL,lLL,uLL,xL,yL;int main() int n,i,j,k,r;printf("请输入矩阵元次：n"); scanf("%d",&n);printf("请输入矩阵各项：n"); for(i=1;i<=n;+i) for(j=1;j<=n;+j) scanf("%lf",&aij); printf("请输入方程组的常数项：n"); for(i=1;i<=n;+i) scanf("%lf",&bi); for(i=1;i<=n;+i) for(j=1;j<=n;+j) lij=0; uij=0.0; for(k=1;k<=n;+k) for(j=k;j<=n;+j) ukj=akj;for(r=1;r<k;+r) ukj-=lkr*urj; for(i=k+1;i<=n;+i) lik=aik; for(r=1;r<k;+r) lik-=lir*urk; lik/= ukk; lkk=1.0; for(i=1;i<=n;+i) yi = bi; for(j=1;j<i;+j) yi-=lij*yj; for(i=n;i>0;-i) xi = yi; for(j=i+1;j<=n;+j) xi-=uij*xj; xi/= uii; for(i=1;i<=n;+i) printf("%0.2lfn",xi); return 0;五、 实验结果与分析完全主元素消元法： 列主元素消元法： LU分解法： 分析： 对于两种高斯解方程，完全主元素跟列主元素都是先消元、再回代，由程序段可以发现，始终消去对角线下方的元素。即，为了节约内存及时效，可以不必计算出主元素下方数据。 列主元素消元法的算法设计上优于完全主元素消元法，它只需依次按列选主元素然后换行使之变到主元素位置，再进行消元即可。 列主元素消元法的耗时比完全主元素法少很多，常采用之。对于LU分解法，分解矩阵为单位下三角阵L与上三角阵U的乘积，然后解方程组Ly=b，回代，解方程组Ux=y。其中的L为n阶单位下三角阵、U为上三角阵.六、讨论、心得 本次试验中，感觉是最难的一次，完全主元素消元法程序编写过程相对来说花了好长时间。纠正各种语法、算法、思路错误。最后勉强成功，但还是有几处警告，不得解决之法。感到程序学习的不足，再加之对高斯的不甚了解。编写过程很是痛苦。 查阅各种内外部资料，这点有利有弊。突然觉得，应该再把数据结构之类的重新学习一下才行。以后多花时间在编程吧，重在理解。 必须反省一下自己的C、C+学习了，还是得多加练习，平时必须养成一种好的算法思维习惯。实验地点虎峪校区致远楼B401指导教师xx 太原理工大学学生实验报告学院名称软件学院专业班级1217班 学号xxxx 学生姓名xx实验日期2014.06.04成绩课程名称数值计算方法实验题目实验三 线性方程组的迭代解法一、实验目的和要求使用雅可比迭代法或高斯-赛德尔迭代法对下列方程组进行求解。 二、主要设备 笔记本 HP ProBook 6470b 一台 编译软件：VC+6.0三、实验内容和原理设线性方程组 Ax=b的系数矩阵A可逆，且主对角元素a11,a22,ann均不为零，令D=diag(a11,a22,ann)并将A分解成 A=(A-D)+D从而线性方程组可写成 Dx=(D-A)x+b则有迭代公式x(k+1)=B1x(k)+f1其中，B1=I-D-1A,f1=D-1b。四、操作方法与实验步骤高斯赛德尔迭代法#include "iostream"#include "iomanip"using namespace std;int main()int i,j,k=0,m,n;double t1,t2,e1,e2=0.0; cout<<"请输入精度e："cin>>e1;cout<<"请输入系数矩阵行数："cin>>m;cout<<"请输入系数矩阵列数："cin>>n;cout<<endl;double (*a)=new double *m;for(i=0;i<=m;i+)ai=new doublen;double (*b)=new double m;double (*x)=new double n;cout<<"请输入系数矩阵："<<endl;cout<<"-"<<endl;for(int num1=0;num1<m;num1+)for(int num2=0;num2<n;num2+)cin>>anum1num2;cout<<endl;cout<<"输入的系数矩阵为："<<endl;for (int num3=0;num3<m;num3+)for(int num4=0;num4<n;num4+)cout<<anum3num4<<" "cout<<endl;cout<<"请输入矩阵b："<<endl;for(int num5=0;num5<m;num5+)cin>>bnum5;cout<<"输入的矩阵b为："<<endl;for(int num6=0;num6<m;num6+)cout<<bnum6<<" " cout<<endl; for(int num7=0;num7<n;num7+)xnum7=0.0000;do cout<<"第"<<k<<"次迭代值："e2=0.0;for(i=0;i<m;i+) double sum=0.0;for(j=0;j<n;j+) if(j!=i)sum+=aij*xj;t1=xi;t2=e2;xi=(bi-sum)/aii;e2=(xi)-t1>=0?(xi)-t1:t1-(xi);e2=(e2>=t2?e2:t2);cout<<setprecision(8)<<xi<<" "cout<<endl;k+; while(e2>=e1&&k<30);cout<<"共迭代了"<<k<<"次"deletea;deleteb;deletex;return 0 ;雅克比迭代法：#include <stdio.h>#include <math.h>int main() float a33=10,-1,-2,-1,10,-2,-1,-1,5,b3=7.2,8.3,4.2;float x3=0,0,0,sum;int i,j,k,n=3;printf("tt X1tt X2tt X3n");for(k=0;k<8;k+) for(i=0;i<3;i+) sum=0; for(j=0;j<n;j+) if(i=j)continue; sum=sum+aij*xj; xi=(bi-sum)/aii; printf("第%d次迭代：t",k+1); for(i=0;i<n;i+) printf("%ft",xi);printf("n"); 五、实验结果与分析高斯赛德尔迭代法：雅克比迭代： 分析： 使用高斯-赛德尔和雅克比迭代都可以求出方程组的解，但是利用高斯-赛德尔迭代法所需的迭代次数比雅克比迭代少，能够更早的达到精度要求。 从程序中可以看出，雅克比定义的sum只有一个，而高斯赛德尔需要两个。时效性上后者要好些。六、讨论、心得 这次试验算是比较成功，要归功于授课时候的认真听讲。程序编写之前，对书本的理论知识进行了进一步的探索。预习准备工作很彻底，自然随后的一切也都很顺利。 实验地点虎峪校区致远楼B401指导教师xx太原理工大学学生实验报告学院名称软件学院专业班级1217班学号xxxx学生姓名xx实验日期2014.06.11成绩课程名称数值计算方法实验题目实验四 代数插值和最小二乘法拟合一、实验目的和要求(1)使用拉格朗日插值法或牛顿插值法求解：已知f(x)在6个点的函数值如下表所示，运用插值方法，求f(0.596)的近似值。X0.400.550.650.800.901.05f(x)0.410750.578150.696750.888111.026521.25386(2)给定数据点（xi ,yi），用最小二乘法拟合数据的多项式，并求平方误差。xi00.50.60.70.80.91.0yi11.751.962.192.442.713.00二、 主要设备 笔记本 HP ProBook 6470b 一台 编译软件：VC+6.0三、实验内容和原理 (1) 设函数在区间a,b上n+1互异节点x0,x1,xn上的函数值分别为y0,y1,yn,求n次插值多项式Pn(x),满足条件Pn(xj)=yj, j=0,1,n令Ln(x)=y0l0(x)+y1l1(x)+ynln(x)= yili(x)其中l0(x),l1(x), ln(x) 为以x0,x1,xn为节点的n次插值基函数，则Ln(x)是一次数不超过n的多项式，且满足Ln(xj)=yj, L=0,1,n再由插值多项式的唯一性，得Pn(x)Ln(x) (2 ) 建立正规方程组：（xij+k）ak=xijyi ,j=0,1,n 平方误差：I=（akxik-yi）2 对给定数据点(Xi，Yi)(i=0,1，m），在取定的函数类 中，求p(x），使误差的平方和E2最小，E2=p(Xi)-Yi2。从几何意义上讲，就是寻求与给定点 (Xi，Yi)(i=0,1，m）的距离平方和为最小的曲线y=p(x）。函数p(x）称为拟合函数或最小二乘解，求拟合函数p(x）的方法称为曲线拟合的最小二乘法。得到的两个关于a0、 a1为未知数的两个方程组，解这两个方程组得出：a0 = （Yi) / m - a1（Xi) / m a1 = mXi Yi - （Xi Yi) / mXi2 - （Xi)2 ) 即最终的拟合多项式各项系数。四、操作方法与实验步骤(1) 代数插值#include <stdio.h> #include <stdlib.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+) f0=yk; for(i=0;i<k;i+)fi+1=(fi-yi)/(xk-xi); yk=fk; return; int main() int i,n; float x20,y20,xx,yy; printf("请输入数据个数n："); scanf("%d",&n);printf("n"); for(i=0;i<=n-1;i+) printf("x%d=",i); scanf("%f",&xi);printf("y%d=",i);scanf("%f",&yi); printf("n"); printf("n"); difference(x,(float *)y,n); printf("请输入插值X："); scanf("%f",&xx); yy=y20; for(i=n-1;i>=0;i-)yy=yy*(xx-xi)+yi; printf("n近似值为：F(%f)=%fn",xx,yy); (2) 最小二乘法拟合#include<iostream.h>#include<fstream.h>#define N 15double power(double &a,int n)double b=1;for(int i=0;i<n;i+)b*=a;return b;void Gauss();double XN,YN,sumXN,sumYN,aNN,bN,lNN,xN;int main()ofstream outdata;ifstream indata;double s;int i,j,k,n,index;cout<<"请输入已知点的个数n="cin>>n;cout<<endl;cout<<"请输入X和Y:"<<endl; for(i=0;i<n;i+)cout<<"X"<<i<<"="cin>>Xi;sumX1+=Xi;cout<<"Y"<<i<<"="cin>>Yi;sumY1+=Yi;cout<<endl;cout<<"sumX1="<<sumX1<<"t"<<"sumY1="<<sumY1<<endl;cout<<"请输入拟合次数index="cin>>index;cout<<endl;i=n;sumX0=i;for(i=2;i<=2*index;i+)sumXi=0;for(j=0;j<n;j+)sumXi+=power(Xj,i);cout<<"sumX"<<i<<"="<<sumXi<<endl;for(i=2;i<=index+1;i+)sumYi=0;for(j=0;j<n;j+)sumYi+=power(Xj,i-1)*Yj;cout<<"sumY"<<i<<"="<<sumYi<<endl;for(i=1;i<=index+1;i+)for(j=1;j<=index+1;j+)aij=sumXi+j-2;bi=sumYi; k=1;dofor(j=k+1;j<=index+1;j+)ljk=ajk/akk;for(i=k+1;i<=index+1;i+)for(j=k+1;j<=index+1;j+)aij=aij-lik*akj; bi=bi-lik*bk;if(k=index+1)break;k+;while(1);xindex+1=bindex+1/aindex+1index+1;for(i=index;i>=1;i-)s=0;for(j=i+1;j<=index+1;j+)s=s+aij*xj;xi=(bi-s)/aii;cout<<"拟合系数为:"for(i=1;i<=index+1;i+)cout<<xi<<"t"double m=0;cout<<endl<<"平方误差为:"for(i=0;i<n;i+) double t=x1+x2*Xi-Yi;m=m+power(t,2);cout<<m<<endl;五、 实验结果与分析 (1) 代数插值 分析： 拉格朗日插值的优点是插值多项式特别容易建立，缺点是增加节点是原有多项式不能利用，必须重新建立，即所有基函数都要重新计算，这就造成计算量的增加。牛顿插值法则很好地避免了上述问题。 (2)最小二乘法拟合 分析： 数据拟合的具体作法是：对给定的数据（xi ,yi）（i=0,1,m），在取定的函数类中，求p(x)属于此函数类，使误差ri=p(xi)- yi (i=0,1,m)的平方和最小， 即 ri2=（p(xi)-yi）2=min从几何意义上讲，就是寻求与给定点（xi ,yi）（i=0,1,m）的距离平方和为最小的曲线y=p(x)。六、讨论、心得 本实验中，插值法有两种插值方法可以选用，最终选用了牛顿插值法。 最小二乘法拟合中， 很好地实现了最小二乘法的程序模拟。对程序编写又更进一步，分析算法、尝试编写。检查并修改其中错误。感觉收获很大。特别是对平方误差的计算模块，一次成功。最后还对程序进行优化，删去冗余以及标准化模块。实验地点虎峪校区致远楼B401指导教师xx专心-专注-专业