计算数值方法实验报告--太原理工大学

1、精选优质文档-倾情为你奉上计算数值方法实验报告学院：软件学院专业：软件工程班级：软件1012班学号：姓名：乔婧峰太原理工大学学生实验报告学院名称软件学院专业班级软件1012班学号学生姓名乔婧峰实验日期20124成绩课程名称数值计算方法实验题目实验一 二分法一、课题名称方程求根：熟悉使用、迭代法、牛顿法、割线法等方法对给定的方程进行根的求解。选择上述方法中的两种方法求方程：二分法f(x)=x3+4x2-10=0在1,2内的一个实根，且要求满足精度|x*-xn|<0.5×10-5迭代法：用迭代公式x=f(x)进行迭代计算，直到满足|x*-xn|<0.5×10-5 为

2、止 。割线法：x=x-f(x)/g(x),其中f(x)为给定的函数，g(x)为给定函数的导数，直到满足|x*-xn|<0.5×10-5 为止 。二、目的和意义（1）了解非线性方程求根的常见方法，如二分法、牛顿法、割线法。（2）加深对方程求根方法的认识，掌握算法。（3）会进行误差分析，并能对不同方法进行比较。三、计算公式f(x)在区间（x，y）上连续 先找到a、b属于区间（x，y），使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

3、+b)/2<0,则在区间（(a+b)/2，b)内有零点，(a+b)/2=>a，从开始继续使用 中点函数值判断。 如果f(a+b)/2>0，则在区间(a,(a+b)/2)内有零点，(a+b)/2<=b，从开始继续使用 中点函数值判断。 这样就可以不断接近零点。 通过每次把f(x)的零点所在小区间收缩一半的方法，使区间的两个端点逐步迫近函数的零点，以求得零点的近似值四、主要仪器设备Vc2008,hp五、结构程序设计迭代法: #include "stdafx.h"#include"stdio.h"#include"math.h

4、"#include"iostream"using namespace std;float main() float a; cin>>a; float t, x; x=a; do x=sqrt(10-x*x*x)/4); t=a; a=x; while(fabs(a-t)>0.5*1e-5); printf("x=%f",a); system("pause");割线法: #include "stdafx.h"#include"stdio.h"#include"

5、math.h"#include"iostream"using namespace std;float main() float c,a=1.0,b=2.0; /cin>>a>>b; 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.) break; b=c; cout<<c;六、结果讨论和分析割线法: 迭代法:分析：使用不同的方法，可以不同程度的求得方程的解，不同的方法速度不同。实验地点ZSA401指导教师李

6、志学院名称软件学院专业班级软件1012班学号学生姓名乔婧峰实验日期2012.4成绩课程名称数值计算方法实验题目实验二 线性方程组的直接解法一、课题名称线性方程组的直接解法合理利用Gauss消元法、LU分解法、追赶法求解下列方程组： （n=5,10,100）二、目的和意义（1）了解线性方程组常见的直接解法，如Guass消元法、LU分解法、追赶法。（2）加深对线性方程组求解方法的认识，掌握算法。（3）会进行误差分析，并能对不同方法进行比较。三、计算公式 高斯分解法：将原方程组化为三角形方阵的方程组：lik=aik/akk aij= aij- lik* akj k=1,2,n-1 i=k+1,k+2

7、, ,n j=k+1,k+2, ,n+1由回代过程求得原方程组的解： xn= ann+1/ ann xk=( akn+1-akj xj)/ akk (k=n-1,n-2, ,2,1) LU分解法：将系数矩阵A转化为A=L*U, L为单位下三角矩阵，U为普通上三角矩阵，然后通过解方程组l*y=b,u*x=y,来求解x.追赶法:用来求对角方程组；将系数矩阵A转化为A=L*U, L为普通下n-1对角矩阵，U为单位上n-1对角矩阵，然后通过解方程组l*y=b,u*x=y,来求解x.四、主要仪器设备Vc2008,hp五、结构程序设计 Gauss消元法：#include "stdafx.h&qu

8、ot;#include"stdio.h"#include"iostream"using namespace std; float 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

9、=%f,",i,j,aij); cout<<endl; 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=%f,",i+1,xi); LU分解法：#include "stdafx.h"#include <stdio.h>#include <math.h> #define L 30 double a L L , b L

10、, l L L , u L L , x L , y L ; int main() int n, i, j, k, r; scanf( "%d", &n ); for ( i = 1; i <= n; +i ) for ( j = 1; j <= n; +j ) scanf( "%lf", &a i j ); for ( i = 1; i <= n; +i ) scanf( "%lf", &b i ); for ( i = 1; i <= n; +i ) for ( j = 1; j &l

11、t;= n; +j ) l i j =0; u i j = 0.0; for ( k = 1; k <= n; +k ) for ( j = k; j <= n; +j ) u k j = a k j ; for ( r = 1; r < k; +r ) u k j -= l k r * u r j ; for ( i = k + 1; i <= n; +i ) l i k = a i k ; for ( r = 1; r < k; +r ) l i k -= l i r * u r k ; l i k /= u k k ; l k k = 1.0; for (

12、i = 1; i <= n; +i ) y i = b i ; for ( j = 1; j < i; +j ) y i -= l i j * y j ; for ( i = n; i > 0; -i ) x i = y i ; for ( j = i + 1; j <= n; +j ) x i -= u i j * x j ; x i /= u i i ; for ( i = 1; i <= n; +i ) printf( "%0.2lfn", x i ); return 0;追赶法：#include "stdafx.h"

13、#include "stdio.h"void main() FILE *f; double a15,b15,c15,d15; double t; int i,n; f=fopen("zgf.txt","r"); fscanf(f,"%d",&n); fscanf(f,"%lf%lf%lf",&b1,&c1,&d1); for(i=2;i<=n-1;i+) fscanf(f,"%lf%lf%lf%lf",&ai,&bi,&a

14、mp;ci,&di); fscanf(f,"%lf%lf%lf",&an,&bn,&dn); fclose(f); c1=c1/b1; d1=d1/b1; For(i=2;i<=n-1;i+) t=bi-ci-1*ai; ci=ci/t; di=(di-di-1*ai)/t; dn=(dn-dn-1*an)/(bn-cn-1*an); for(i=n-1;i>=1;i-) di=di-ci*di+1; printf("n*n"); for(i=1;i<=n;i+) printf("d%2d=%l

15、fn",i,di);Zgf.txt文件中的内容是：52 1 -71 2 1 -5 1 2 1 -5 1 2 1 -5 1 2 -5六、结果讨论和分析Gauss消元法：LU分解法：追赶法：分析从消元过程可以看出，对于n阶线性方程组，只要各步主元素不为零，经过n-1步消元，就可以得到一个等价的系数矩阵为上三角形阵的方程组，然后再利用回代过程可求得原方程组的解.消元过程相当于分解 A为单位下三角阵L与上三角阵U的乘积，解方程组Ly=b回代过程就是解方程组Ux=y。其中的L为n阶单位下三角阵、U为上三角阵. 在 A 的LU 分解中, L取下三角阵, U 取单位上三角阵,这样求解方程组Ax=d

16、 的方法称为追赶法.实验地点ZSA401指导教师李志学院名称软件学院专业班级软件1012学号学生姓名乔婧峰实验日期20114成绩课程名称数值计算方法实验题目实验三 线性方程组的迭代解法一、课题名称线性方程组的迭代解法使用雅可比迭代法或高斯-赛德尔迭代法对下列方程组进行求解。 二、目的和意义学习使用雅可比迭代法或高斯-赛德尔迭代法三、计算公式雅克比迭代法：设线性方程组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其中，B

17、1=I-D-1A,f1=D-1b。四、主要仪器设备Vc2008,hp五、结构程序设计高斯-赛德尔迭代法：#include "stdafx.h"#include <stdio.h>#include <math.h>void main() float a33=10,-1,-2,-1,10,-2,-1,-1,5,b3=7.2,8.3,4.2;float x3=0,0,0,sum1,sum2;int i,j,k,n=3;for (k=0;k<10;k+) for(i=0;i<n;i+) sum1=0; sum2=0;for(j=0;j<i;

18、j+) sum1=sum1+aij*xj; for(j=i+1;j<3;j+) sum2=sum2+aij*xj; xi=(bi-sum1-sum2)/aii; for(i=0;i<n;i+) printf("x%d=%f,",i+1,xi); printf("n"); 雅克比迭代：#include "stdafx.h"#include <stdio.h>#include <math.h>void main() float a33=10,-1,-2,-1,10,-2,-1,-1,5,b3=7.2,8

19、.3,4.2;float x3=0,0,0,sum1;int i,j,k,n=3;for (k=0;k<10;k+) for(i=0;i<3;i+) sum1=0; for(j=0;j<n;j+) if(i=j) continue; sum1=sum1+aij*xj; xi=(bi-sum1)/aii; for(i=0;i<n;i+) printf("x%d=%f,",i+1,xi);printf("n"); 六、实验结果与分析：高斯-赛德尔迭代法：雅克比迭代：分析：使用高斯-赛德尔和雅克比迭代都可以求出方程组的解，但是利用高斯-

20、赛德尔迭代法所需的迭代次数比雅克比迭代少，能够更早的达到精度要求。实验地点ZSA401指导教师李志学院名称软件学院专业班级软件1012班学号学生姓名乔婧峰实验日期20114成绩课程名称数值计算方法实验题目实验四 最小二乘法拟合多项式一、课题名称（1）了解矩阵特征值与特征向量问题解法，掌握幂法。（2）加深对矩阵特征值与特征向量问题求解方法的认识，掌握算法。（3）会进行误差分析。二、目的和意义 学习使用最小二乘法拟合多项式三、计算公式幂法：由已知的非零向量x0和矩阵A的乘幂构造向量序列xn以计算矩阵A的按模最大特征值及其特征向量的方法，称为幂法。迭代公式：结果可取 四、主要仪器设备Vc2008,h

21、p五、结构程序设计 五、结果讨论和分析分析： 幂法是一种求任意矩阵A的按模最大特征值及其对应特征向量的迭代算法。该方法的最大优点是计算简单，容易在计算机上实现，对稀疏矩阵较为适合，但有时收敛速度很慢。实验地点综合楼六层606室指导教师王峥学院名称计算机科学与技术专业班级计算机学号学生姓名某某实验日期2011-6-20成绩课程名称数值计算方法实验题目实验五 代数插值一、 课题名称使用拉格朗日插值法或牛顿插值法求解：已知f(x)在6个点的函数值如下表所示，运用插值方法，求f(0.596)的近似值。x0.400.550.650.800.901.05f(x)0.410750.578150.696750

22、.888111.026521.25386二、 目的和意义学习使用拉格朗日插值法或牛顿插值法求解三、 计算公式设函数在区间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)四、 结构程序设计#include<st

23、dio.h>#include<stdlib.h>#include<iostream.h>typedef struct data float x; float y;Data; /变量x和函数值y的结构Data d20; /最多二十组数据float f(int s,int t) /牛顿插值法，用以返回插商 if(t=s+1) return (dt.y-ds.y)/(dt.x-ds.x); else return (f(s+1,t)-f(s,t-1)/(dt.x-ds.x); float Newton(float x,int count) int n; while(1)

24、 cout<<"请输入n值(即n次插值):"/获得插值次数 cin>>n; if(n<=count-1)/ 插值次数不得大于count1次 break; else system("cls"); float t=1.0; float y=d0.y; float yt=0.0; for(int j=1;j<=n;j+) t=(x-dj-1.x)*t; yt=f(0,j)*t; y=y+yt; return y;float lagrange(float x,int count) float y=0.0; for(int k=

25、0;k<count;k+)/这儿默认为count1次插值 float p=1.0;/初始化p for(int j=0;j<count;j+) /计算p的值 if(k=j)continue;/判定是否为同一个数 p=p*(x-dj.x)/(dk.x-dj.x); y=y+p*dk.y;/求和 return y;/返回y的值void main() float x,y; int count; while(1) cout<<"请输入xi,yi的组数，不得超过20组:"/要求用户输入数据组数 cin>>count; if(count<=20)

26、 break;/检查输入的是否合法 system("cls"); /获得各组数据 for(int i=0;i<count;i+) cout<<"请输入第"<<i+1<<"组x的值:" cin>>di.x; cout<<"请输入第"<<i+1<<"组y的值:" cin>>di.y; system("cls"); cout<<"请输入x的值:"/获

27、得变量x的值 cin>>x; while(1) int choice=3; cout<<"请您选择使用哪种插值法计算："<<endl; cout<<" (0):退出"<<endl; cout<<" (1):Lagrange"<<endl; cout<<" (2):Newton"<<endl; cout<<"输入你的选择：" cin>>choice;/取得用户的选择项

28、 if(choice=2) cout<<"你选择了牛顿插值计算方法，其结果为：" y=Newton(x,count);break;/调用相应的处理函数 if(choice=1) cout<<"你选择了拉格朗日插值计算方法，其结果为：" y=lagrange(x,count);break;/调用相应的处理函数 if(choice=0) break; system("cls"); cout<<"输入错误!"<<endl; cout<<x<<&quo

29、t; , "<<y<<endl;/输出最终结果 五、结果讨论和分析分析：拉格朗日插值的优点是插值多项式特别容易建立，缺点是增加节点是原有多项式不能利用，必须重新建立，即所有基函数都要重新计算，这就造成计算量的浪费。实验地点综合楼六层606室指导教师王峥学院名称计算机科学与技术专业班级计算机学号学生姓名某某实验日期2011-6-20成绩课程名称数值计算方法实验题目实验六 最小二乘法拟合多项式一、课题名称给定数据点（xi ,yi），用最小二乘法拟合数据的多项式，并求平方误差。xi00.50.60.70.80.91.0yi11.751.962.192.442.713

30、.00二、目的和意义1熟练运用已学计算方法求解方程组2加深对计算方法技巧，选择正确的计算方法来求解各种方程组3培养使用电子计算机进行科学计算和解决问题的能力三、计算公式建立正规方程组：（xij+k）ak=xijyi ,j=0,1,n 平方误差：I=（akxik-yi）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();d

31、ouble XN,YN,sumXN,sumYN,aNN,bN,lNN,xN;void 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<<"="c

32、in>>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<<