太原理工大学数值计算方法实验报告.docx
- 文档编号:25878043
- 上传时间:2023-06-16
- 格式:DOCX
- 页数:21
- 大小:108.04KB
太原理工大学数值计算方法实验报告.docx
《太原理工大学数值计算方法实验报告.docx》由会员分享,可在线阅读,更多相关《太原理工大学数值计算方法实验报告.docx(21页珍藏版)》请在冰豆网上搜索。
太原理工大学数值计算方法实验报告
本科实验报告
课程名称:
计算机数值方法
实验项目:
方程求根、线性方程组的直接解
法、线性方程组的迭代解法、代数插值和最
小一乘拟合多项式
实验地点:
行勉
楼
专业班级:
********^号:
*********
学生:
********
指导教师:
誌,
冬华
2016年4月8日
学生
实验成绩
实验名称
实验一方程求根
实验容和要求
熟悉使用二分法、迭代法、牛顿法、割线法等方法对给定的方程进行根的求解。
选择上
述方法中的两种方法求方程:
f(x)=x3+4x2-10=0在[1,2]的一个实根,且要求满足精度
*-5
|x-Xn|<0.5X10
(1)了解非线性方程求根的常见方法,如二分法、牛顿法、割线法。
(2)加深对方程求根方法的认识,掌握算法。
(3)会进行误差分析,并能对不同方法进行比较。
实验原理
1.二分法:
如果要求已知函数f(x)=0的根(x的解),那先要找出一个区间[a,b],
使得f(a)与f(b)异号。
根据介值定理,这个区间一定包含着方程式的根。
求该区间的中
点m=(a+b)/2,并找出f(m)的值。
若f(m)与f(a)正负号相同,则取[m,b]为新的区间,否则取[a,m]。
重复第3步和第4步,直到得到理想的精确度为止。
2.割线法是利用牛顿迭代法的思想,在根的某个领域,函数有直至二阶的连续导数,并且
不等于0,则在领域选取初值X0,X1,迭代均收敛。
(1)在区间[m,n]输入初值x0,x1.
(2)计算x2。
x2=x1-f(x1)*(x1-x0)/(f(x1)-f(x0))
(3)x0-x1,x1-x2(4)判断是否达到精度,若是输出x1,若否执行
(2)
主要仪器设备
HP计算机
实验记录
1.二分法
//方程求根(二分法).cpp:
定义控制台应用程序的入口点。
//
#include"stdafx.h"
#include"iostream"
usingnamespacestd;
classText
{
public:
floatx,y,a,b,c,n=0;
voidGetab()
{
cout<<"请输入计算区间:
(以空格隔开)"<
}
floatGetY(floatx)
{
y=x*x*x+4*x*x-10;
returny;
}
floatCalculate(floatafloatb)
{
c=(a+b)/2;
n++;
if(GetY(c)==0||((b-a)/2)<0.000005)
{
cout< } if(GetY(a)*GetY(c)<0) { returnCalculate(a,c); } if(GetY(c)*GetY(b)<0) { returnCalculate(c,b); } } }; intmain() { cout<<"方程组为: f(x)=xA3+4xA2-10=0"< Texttext; text.Getab(); a=text.a; b=text.b; text.Calculate(a,b); return0; } I QI匚^WIINDClWSX-syrtem 方程组为;代孟)円「3呜孟龙120if输入计算区同*(以空格隔开)12 1.36523为方理的禅请按任意惺雒续…. 2.割线法: //方程求根(割线法).cpp: 定义控制台应用程序的入口点。 // #include"stdafx.h" #include"iostream" usingnamespacestd; classA { public: floatx0,x1,y; floatGetY(floatx) { y=x*x*x+4*x*x-10; returny; } voidGetNumber() { cout<<"请输入两个初始近似值: (以空格隔开)"< cin>>x1; } voidCalculate(floatx0,floatx1) { floatx2; x2=x1-(GetY(x1)/(GetY(x1)-GetY(x0))*(x1-x0)); if(x2==x1) { cout< } else { cout< returnCalculate(x1,x2); } } }; intmain() { cout<<"方程组为: f(x)=xA3+4xA2-10=0"< floata,b; Atext; text.GetNumber(); a=text.xO; b=text.xl; text.Calculate(a,b); return0; SBC: \VVINDOWS\sy$tern3Cumd.exe |方程组为2f仗)=T3十4x"2-10=0请输入两个初始近似佰;(以咛格隔卄) 12 1.26316 1*33883 L36662 1.36521 L36523 1.3仙23为方程的解 请按仕意键继续…: 心得体会 使用不同的方法,可以不同程度的求得方程的解,通过二分法计算的程序实现更加了 解二分法的特点,二分法过程简单,程序容易实现,但该方法收敛比较慢一般用于求根的初始近似值,不同的方法速度不同。 面对一个复杂的问题,要学会简化处理步骤,分步骤一点一点的循序处理,只有这样,才能高效的解决一个复杂问题。 实验名称 实验二线性方程组的直接求解 实验容和要求 合理选择利用Gauss消兀法、主兀素消兀法、LU分解法、追赶法求解下列方程组: 123x,14 012x28 241X313 (1)了解线性方程组常见的直接解法,如Guass消元法、LU分解法、追赶法。 (2)加深对线性方程组求解方法的认识,掌握算法。 (3)会进行误差分析,并能对不同方法进行比较。 实验原理 1.咼斯分解法: ⑴将原方程组化为三角形方阵的方程组: lik=aik/akk aij=aij-lik*akjk=1,2,…,n-1 i=k+1,k+2,…,nj=k+1,k+2,…,n+1 ⑵由回代过程求得原方程组的解: xn=ann+1/ann xk=(akn+1-刀akjxj)/akk(k=n-1,n-2,…,2,1) 2.LU分解法: 将系数矩阵A转化为A=L*U,L为单位下三角矩阵,U为普通上三角矩阵,然后 通过解方程组l*y=b,u*x=y,来求解x. 主要仪器设备 HP计算机 实验记录 1.高斯消兀法: #include"stdio.h" #include"math.h" #include voidExchange(inti) { intj,l,k; doublemax=a0[i][i],temp;j=i; for(k=i;k<=3;k++) { if(a0[k][i]>max) { max=aO[k][i]; j=k; } } for(l=i;l<=4;l++) { temp=aO[i][l]; aO[i][l]=aO[j][l]; a0[j][l]=temp; } for(i=1;i<=3;i++) { for(j=1;j<=4;j++) { a[i][j]=a0[i][j]; } } } voiddisplayA() { inti,j;printf("\n"); for(j=1;j<=3;j++) { for(i=1;i<=4;i++) printf("%lf",a[j][i]); printf("\n"); } } voidmain() { inti,j,k; for(i=1;i<=3;i++) {for(j=1;j<=4;j++) {scanf("%lf",&a[i][j]); a0[i][j]=a[i][j]; } } displayA();printf(”列主元素消元法如下"); //消元过程 k=1; do { Exchange(k); displayA(); for(i=k+1;i<=3;i++) { l[i]=aO[i][k]/aO[k][k]; printf("l[%i][%i]=%lf",i,k,l[i]);for(j=k;j<=4;j++) {a[i][j]=aO[i][j]-l[i]*aO[k][j]; }displayA(); } k++; if(k==3)break; for(j=1;j<=3;j++) { for(i=1;i<=4;i++)a0[j][i]=a[j][i]; } }while (1); //回代过程 l[3]=a[3][4]/a[3][3]; for(k=3;k>=1;k--) { tmp=0; for(j=k+1;j<=3;j++)tmp+=a[k][j]*l[j];l[k]=(a[k][4]-tmp)/a[k][k]; } for(i=1;i<=3;i++) printf("x[%i]=%lf\n",i,l[i]); 2314 12S 4113 2.LU分解法: #include #include inti,j,k,r; doublem=0,p=0; doublea[3][3]; voidlu(doublea[3][3]) { for(i=1;i<=2;i++) { if(a[O][O]! =O)a[i][0]=a[i][0]/a[0][0]; } for(k=1;k<=2;k++) { for(j=k;j<=2;j++) { { for(r=0;r<=k-1;r++) m=m+a[k][r]*a[r][j]; }a[k][j]=a[k][j]-m; m=0; }for(i=k+1;i<=2;i++) {{for(r=0;r<=k-1;r++) P=P+a[i][r]*a[r][k]; }a[i][k]=(a[i][k]-p)/a[k][k];p=0;}}} voidmain(){ staticdoublea[3][3]={{1,2,3},{0,1,2},{2,4,1}}; staticdoubleb[3]={14,8,13}; doublec[3]; doubled[3]; doublef[3][3]; doublem=0; doublen=0; intr; inti,j; lu(a); printf("输出U的矩阵为\n"); for(i=0;i<=2;i++) { for(j=i;j<=2;j++) { f[i][j]=a[i][j]; printf("%f',f[i][j]); } printf("\n”); } printf("输出L的矩阵为\n”); for(i=0;i<=2;i++) { for(j=0;j<=i;j++) { if(i==j) { a[i][j]=1; printf("%f",a[i][j]); } else printf("%f",a[i][j]); } printf("\n"); } &[0]=b[0]; for(i=1;i<=2;i++) { for(r=0;r<=i-1;r=r+1) m=m+a[i][r]*c[r]; c[i]=b[i]-m; m=0; } d[2]=c[2]/f[2][2]; for(i=1;i>=0;i=i-1) { for(r=2;r>i;r=r-1) n=n+f[i][r]*d[r]; d[i]=(c[i]-n)/f[i][i]; n=0; } printf(”所求方程组解为x仁%f,x2=%f,x3=%f",d[0],d[1],d[2]);/*根据LU分解所得两 个矩阵及求解步骤计算所求X一组解*/ } 1.0300002.806m3.BM008 1.00(1^0 输岀L的矩阵为 1,006000 0”0000001.000090 2.0300000.0000901.000000 馬求方程组解xi=1.000000,xS=2-,x3=3.QSQSBSPressanijkeyt: ocontiniuie 心得体会 对于求解线性方程组的各种直接方法来说各有优缺点,在所有的求解方法中都应该注意其解的精度。 注意不同求解方法的不同误差求法。 编写程序的时候需要一步一步慢慢来,逐步增加自己的算法知识水平和解决问题的能力。 实验名称 实验三线性方程组的迭代求解 实验容和要求 10捲x22x37.2 x-i10x22x38.3 捲x25x34.2 使用雅可比迭代法或高斯-赛德尔迭代法对下列方程组进行求解。 实验原理 雅可比迭代法: 设线性方程组 Ax=b 的系数矩阵A可逆且主对角元素aii,a22,…,ann均不为零,令 D=diag(a11,a22,…,ann) 并将A分解成 A=(A-D)+D 从而线性方程组可写成 Dx=(D-A)x+b则有迭代公式 x(k+1)=Bx⑹+f1 其中,Bi=I-DA,f1=Db。 主要仪器设备 HP计算机 实验记录 #inelude #inelude intmain() { inti; doublex1[20],x2[20],x3[20]; doublex10,x20,x30; printf("pleaseinputx1,x2,x3: \n"); seanf("%lf%lf%lf",&x10,&x20,&x30); printf(”nx1[n]x2[n]x3[n]\n"); for(i=0;i<18;i++) { x1[0]=x10; x2[0]=x20; x3[0]=x30; x1[i+1]=0.1*x2[i]+0.2*x3[i]+0.72; x2[i+1]=0.1*x1[i]+0.2*x3[i]+0.83; x3[i+1]=0.2*x1[i]+0.2*x2[i]+0.84; printf("%5d%5lf%5lf%5lf\n",i,x1[i],x2[i],x3[i]); }return0; } Ijpleaseinputscl,x2,x2: IsB0 n xlTnl x3[nl 0 0.006060 0.000000 0.030390 1 0.720066 0.S300M 0.848908 2 fl.971686 1.3700B0 l.i&aoao 3 1・S57SSS 1・157100 1.24S298 4 1 1.185340 1.282820 5 1.095098 i・195099 1.294138 G i.S98338 1-198337 1.298039 7 1.B99442 1.199442 1,2? ? 335 8 1・099611 1.199811 1・277777 9 1.099936 1.19993C 1.299924 IM 1.899979 1.199979 1.299975 11 1.B99993 1.199993 1.299991 12 1・B9999S 1.199998 1・2? ? 9? 7 13 1.099999 1.199999 1.2999? 9 丄4 l.mm 1.200000 1.308908 15 1.1BQS0S 1.260880 1.399098 16 1・1BG80Q 1 1,3nnoaa 1? 1.160066 1.200000 1・330300 心得体会 在编写算法是不熟悉,查阅了很多资料,经过反复研究和试验后实现了题目的要求,使 用雅克比迭代法和高斯-赛德尔都可以得到方程的解,但相比之下,高斯-赛德尔的迭代次数要比雅克比的迭代次数少,能够更快的达到所求的解的精度。 实验名称 实验四代数插值和最小二乘法拟合 实验容和要求 1•学习使用拉格朗日插值法或牛顿插值法求解方法。 2•了解最小二乘法的多项式拟合的具体计算方法并且注意克服正规方程组的病态。 给定数据点(Xi,yi)如下: Xi 0 0.5 0.6 0.7 0.8 0.9 1.0 yi 1 1.75 1.96 2.19 2.44 2.71 3.00 ⑴使用扌 (2)用最丿 ⑶对比、 拉格朗日插值法或牛顿插值法,求f(0.856)的近似值.小二乘法拟合数据的(n次)多项式,求f(0.856)的近似值. 分析上两结果 实验原理 设函数在区间[a,b]上n+1互异节点xo,x1,…,xn上的函数值分别为yo,yi,…,yn, 求n次插值多项式Pn(X),满足条件 Pn(Xj)=yj,j=0,1,…,n 令 Ln(x)=yolo(x)+yili(x)+…+ynln(x)=刀yili(x) 其中10(x),l1(X),…,ln(x)为以Xo,Xl,…,xn为节点的n次插值基函数,则Ln(x)是一次数不超过n的多项式,且满足 Ln(xj)=yj,L=0,1,…,n 再由插值多项式的唯一性,得 Pn(x)=Ln(x) 主要仪器设备 HP计算机 实验记录(写出实验容中的程序代码和运行结果)(可分栏或加页) 拉格朗日插值法: #include"stdio.h" intmain() { doublem=1.0,a=0.856,l=0; inti,j; doublex[6]={0.50,0.60,0.70,0.80,0.90,1.00}; doubley[6]={1.75,1.96,2.19,2.44,2.71,3.00}; for(i=0;i<=5;i++) { for(j=0;j<=5;j++) { if(i==j)continue; m=m*((a-x[j])/(x[i]-x[j])); } l+=y[i]*m;m=1; } printf(”结果为%lf",l); return0; } 最小二乘法: #include"stdio.h" #include"math.h" intmain() {doublex[7]={0,0.5,0.6,0.7,0.8,0.9,1.0},y[7]={1,1.75,1.96,2.19,2.44,2.71,3.00},a0,a1,sum1=0,sum2=0,sum3=0,sum4=0,sum5=0,l,r; intm=6,i,k; for(i=0;i<7;i++) { sum1+=x[i]; sum2+=x[i]*x[i]; sum3+=y[i]; sum4+=x[i]*y[i]; sum5+=y[i]*y[i]; } l=sum1/(m+1); a仁(sum4-l*sum3)/(sum2-l*sum1); a0=(sum3-sum1*a1)/(m+1); doubles=sum3*a0+sum4*a1; r=sum5-s; printf("y=a0+a1*x\n"); printf("aO=%fa1=%f\t\n",a0,a1,r); doubleq=0.856,p; p=a0+a1*q; printf("y=%f\n",p); return0; } y=a0+al«xa0=0_9782Glal=1-.978261ly=2.571652 Pressanykeytocontinue 心得体会 拉格朗日插值的优点是插值多项式特别容易建立,缺点是增加节点是原有多项 式不能利用,必须重新建立,即所有基函数都要重新计算,这就造成计算量的浪费。 牛顿插值多项式的优点是增加节点时,原先的差商仍旧不变,仍可以使用。 数据拟合的具体作法是: 对给定的数据(Xi,yi)(i=0,1,…,m),在取定的函 数类中,求p(x)属于此函数类,使误差ri=p(xj-yi(i=0,1,…,m)的平方和最小,即: 刀ri2=E(Ep(xj-yi)2=min 从几何意义上讲,就是寻求与给定点(Xi,yi)(i=0,1,…,m)的距离平方和为最小 的曲线y=p(x)。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 太原 理工大学 数值 计算方法 实验 报告