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

时间：2020-11-12 17:08:52　来源：勤学考试网本文已影响人

TAniLLN UNBTRSIIY OF TECHNOLOGY

本科实验报告

课程名称：

计算机数值方法

实验项目：

方程求根、线性方程组的直接解

法、线性方程组的迭代解法、代数插值和最

小二乘拟合多项式

实验地点：行勉楼

专业班级：

********

学号：

*********

学生姓名：

********

指导教师：

李

誌，崔

冬华

2016年4月8日

学生姓名

实验成绩

实验名称

实验一方程求根

实验内容和要求

熟悉使用二分法、迭代法、牛顿法、割线法等方法对给定的方程进行根的求解。选择上

述方法中的两种方法求方程：f(x)=x 3+4x2-10=0在［1,2］内的一个实根，且要求满足精度

* -5

|x -Xn|<0.5 X 10

了解非线性方程求根的常见方法，如二分法、牛顿法、割线法。

加深对方程求根方法的认识，掌握算法。

(3 )会进行误差分析，并能对不同方法进行比较。

实验原理

二分法：如果要求已知函数 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步，直到得到理想的精确度为止。

割线法是利用牛顿迭代法的思想，在根的某个领域内，函数有直至二阶的连续导数，并

且不等于0,则在领域内选取初值 xo,x 1,迭代均收敛。

在区间［m ,n］内输入初值x0, x1.

计算 x2。x2=x1-f(x1)*(x1-x0)/(f(x1)-f(x0))

x0-x1,x1-x2 (4)判断是否达到精度，若是输出x1，若否执行(2)

主要仪器设备

HP计算机

实验记录

1.二分法

//方程求根(二分法).cpp :定义控制台应用程序的入口点。

#i nclude "stdafx.h"

#i nclude"iostream"

using n amespace std;

class Text

{

public:

float x, y, a, b, c, n - 0;

void Getab()

{

cout << "请输入计算区间：(以空格隔开)"<< endl; cin >> a >> b;

}

float GetY(float x)

{

y = x*x*x + 4 * x*x - 10;

return y;

}

float Calculate(float afloat b)

{

c = (a + b) / 2;

n++;

if (GetY(c) == 0 || ((b - a) / 2) < 0.000005)

{

cout << c <<"为方程的解"<< endl; return 0;

}

if (GetY(a)*GetY(c) < 0)

{

return Calculate(a,c);

}

if (GetY(c)*GetY(b)< 0)

{

retur n Calculate(c,b);

}

};

int mai n()

{

cout << "方程组为：f(x)=xA3+4xA2-10=0" << endl; float a, b;

Text text;

text.Getab();

a = text.a;

b = text.b;

text.Calculate(a, b);

return 0;

}

QI 匚 ^WIINDClWSX-syrtem

方程组为；代孟）円「3呜孟龙120 if输入计算区同* （以空格隔开） 1 2

1.36523为方理的禅请按任意惺雒续….

2.割线法：

//方程求根（割线法）.cpp :定义控制台应用程序的入口点。

#i nclude "stdafx.h"

#in clude"iostream"

using n amespace std;

class A

{

public:

float x0,x1,y;

float GetY(float x)

{

y= x*x*x+4*x*x-10;

return y;

}

void GetNumber()

{

cout<<"请输入两个初始近似值：(以空格隔开)"<< endl; cin >> x0;

cin >> x1;

}

void Calculate(float x0,float x1)

{

float x2;

x2 = x1 - (GetY(x1) / (GetY(x1) - GetY(x0))*(x1 - x0));

if (x2==x1)

{

cout <<x2<<"为方程的解"<< endl;

}

else

{

cout << x2 << en dl;

return Calculate(x1, x2);

}

};

int mai n()

{

cout << "方程组为：f(x)=xA3+4xA2-10=0" << endl;

float a, b;

A text;

text.GetNumber();

a = text.xO;

b = text.xl;

text.Calculate(a,b);

return 0;

SB C：\VVI N D OWS\sy $ tern 3 Cum d.exe

|方程组为2 f仗）=T3十4x"2-10=0 请输入两个初始近似佰；（以咛格隔卄）

1 2

1. 26316

1 * 33883

L 36662

1. 36521

L 36523

1. 3仙23为方程的解

请按仕意键继续…:

心得体会

使用不同的方法，可以不同程度的求得方程的解，通过二分法计算的程序实现更加了

解二分法的特点，二分法过程简单，程序容易实现，但该方法收敛比较慢一般用于求根的初始近似值，不同的方法速度不同。面对一个复杂的问题，要学会简化处理步骤，分步骤一点一点的循序处理，只有这样，才能高效的解决一个复杂问题。

实验名称

实验二线性方程组的直接求解

实验内容和要求

合理选择利用Gauss消兀法、主兀素消兀法、 LU分解法、追赶法求解下列方程组：

2 3 x, 14

0 1 2 x2 8

4 1 X3 13

了解线性方程组常见的直接解法，如 Guass消元法、LU分解法、追赶法。

加深对线性方程组求解方法的认识，掌握算法。

会进行误差分析，并能对不同方法进行比较。

实验原理

咼斯分解法：

⑴将原方程组化为三角形方阵的方程组：

l ik =aik /a kk

a ij = a ij - l ik * a kj k=1,2, …,n-1

i=k+1,k+2, …,n j=k+1,k+2, …,n+1

⑵由回代过程求得原方程组的解：

x n= a nn+1/ a nn

x k=( a kn+1-刀 akj x j)/ a kk (k=n-1,n-2, …,2,1)

LU分解法：

将系数矩阵A转化为A=L*U, L为单位下三角矩阵，U为普通上三角矩阵，然后

通过解方程组l*y=b,u*x=y, 来求解x.

主要仪器设备

HP计算机

实验记录

1.高斯消兀法：

#i nclude "stdio.h"

#i nclude "math.h"

#i nclude <stdio.h> double a[5][ 6],a0 [5][6]; double l[5],tmp;

void Excha nge(i nt i)

{

int j,l,k;

double max=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];

}

void displayA()

{

int i,j;pri ntf("\n");

for(j=1;j<=3;j++)

{

for(i=1;i<=4;i++)

prin tf("%lf ",a[j][i]);

prin tf("\n");

}

void mai n()

{

int i,j,k;

for(i=1;i<=3;i++)

{ for(j=1;j<=4;j++)

{ scan f("%lf",&a[i][j]);

a0[i][j]=a[i][j];

}

displayA(); printf(”列主元素消元法如下");

//消元过程

k=1;

{

Excha nge(k);

displayA();

for(i=k+1;i<=3;i++)

{

l[i]=aO[i][k]/aO[k][k];

prin tf("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++)

prin tf("x[%i]=%lf\n",i,l[i]);

2 3 14

1 2 S

4 1 13

2-0B000U1 060000 4 *0盹盹 0列主元索消元法如下2.G6OOO0 Q?6Q0C00 1.000000 1[21E1J-0 2?06000Q e^sesea 1.000000 1L31E1J-0 2.QH060Q

2-0B000U

1 .9839801

2.060000 4 *0盹盹 0

列主元索消元法如下

2.G6OOO0 Q?6Q0C00 1.000000 1[21E1J-0 2?06000Q e^sesea 1.000000 1L31E1J-0 2.QH060Q

4.000900 1?000盹0 s.aeena 000000

4.B0B00B 1?000900

2.0BMH0 530000

4.D000O0

0.066060

2.0U060Q

1.0000001 s. aeawa

4.00W00

1 - 0009001

0.000000 0.000000

l[3H2]=B.OB(90n

2.QU0O0Q 4.00dUUU

6.000000

1 - 9039001

3.090099

2.03丽阳

1? 900008

1.099009

2 ? 030300

3? 009300

1? aaoaaa

2 ? 090090

a ” 000000

1.00090(9

2.000900

2.^00088

1.0389km

2.090900

2.590009

1 .

2.098990

14.000030

9 _B008SB

U? 003000

L3?000000 8?000000 14.000HAR

直3?000000 B?B0000B

14.000000

13.000000

8.S00000

7?580688

13.000000

8.B00000

7.500000

6.000000 0.000000 2.530300

x[l1^1.000000

工[2]=2-盹0盹0

x[3 1=3 .盹盹00

Pre^s 殆ny key to cont inuie—

13.000000

8.S3000B

7-500000

2.LU分解法：

#i nclude<stdio.h>

#in clude<math.h>

int i,j,k,r;

double m=0,p=0;

double a[3][3];

void lu(double a[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; } } }

void mai n() {

static double a[3][3]={{1,2,3},{0,1,2},{2,4,1}};

static double b[3]={14,8,13};

double c[3];

double d[3];

double f[3][3];

double m=0;

double n=0;

int r;

int i,j;

lu(a);

printf("输出U的矩阵为\n");

for(i=0;i<=2;i++)

{

for(j=i;j<=2;j++)

{

f[i][j]=a[i][j];

prin tf(" %f',f[i][j]);

}

prin tf("\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]);

}

prin tf("\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.030000 2.806m 3.BM008

1 .00(1^0

输岀L的矩阵为

1,006000

0”000000 1.000090

2.030000 0.000090 1.000000

馬求方程组解xi =1.000000, xS =2 -, x3 =3 .QSQSBSPress an ij key t：o con tiniuie

心得体会

对于求解线性方程组的各种直接方法来说各有优缺点，在所有的求解方法中都应该注意其解的精度。注意不同求解方法的不同误差求法。编写程序的时候需要一步一步慢慢来，逐步增加自己的算法知识水平和解决问题的能力。

实验名称

实验三线性方程组的迭代求解

实验内容和要求

10捲 x2 2x3 7.2

x-i 10x2 2x3 8.3

捲 x2 5x3 4.2

使用雅可比迭代法或高斯-赛德尔迭代法对下列方程组进行求解。

实验原理

雅可比迭代法：

设线性方程组

Ax=b

的系数矩阵A可逆且主对角元素 aii,a22,…,ann均不为零，令

D=diag(a 11 ,a 22,…,a nn)

并将A分解成

A=(A-D)+D

从而线性方程组可写成

Dx=(D-A) x+b 则有迭代公式

x(k+1)=Bx⑹ +f1

其中，Bi=I-D A,f 1=D b。

主要仪器设备

HP计算机

实验记录

#inelude <stdio.h>

#in elude <math.h>

int main()

{

int i;

double x1[20],x2[20],x3[20];

double x10,x20,x30;

printf( "please in put x1,x2,x3:\n");

seanf("%lf%lf%lf"，&x10,&x20,&x30);

printf(” n x1[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]);

} return 0;

}

Ijplease input scl ,x2 ,x2 :

Is B 0

xlTnl

x3[nl

0.006060

0.000000

0.030390

0.720066

0.S300M

0.848908

fl.971686

1.3700B0

l.i&aoao

1?S57SSS

1?157100

1.24S298

1.185340

1.282820

1.095098

i?195099

1.294138

i.S98338

1-198337

1.298039

1.B99442

1.199442

1,2335

1?099611

1.199811

1?277777

1.099936

1.19993C

1.299924

1.899979

1.199979

1.299975

1.B99993

1.199993

1.299991

1?B9999S

1.199998

1?29?7

1.099999

1.199999

1.2999?9

丄4

l.mm

1.200000

1.308908

1.1BQS0S

1.260880

1.399098

1?1BG80Q

1,3nnoaa

1.160066

1.200000

1?330300

心得体会

在编写算法是不熟悉，查阅了很多资料，经过反复研究和试验后实现了题目的要求，使

用雅克比迭代法和高斯-赛德尔都可以得到方程的解，但相比之下，高斯-赛德尔的迭代次数要比雅克比的迭代次数少，能够更快的达到所求的解的精度。

实验名称

实验四代数插值和最小二乘法拟合

实验内容和要求

1?学习使用拉格朗日插值法或牛顿插值法求解方法。

2?了解最小二乘法的多项式拟合的具体计算方法并且注意克服正规方程组的病态。

　给定数据点(Xi , yi)如下：

0.5

0.6

0.7

0.8

0.9

1.0

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,x 1,…,xn上的函数值分别为 yo,yi,…,y n,

求n次插值多项式Pn(X),满足条件

Pn(Xj)=yj, j=0,1, …,n

令

Ln(x)=y ol o(x)+y il i(x)+ …+ynl n(x)= 刀 yi l i (x)

其中10(x),l 1(X)，…,l n(x)为以Xo,Xl，…,x n为节点的n次插值基函数，则Ln(x)是一次数不超过 n的多项式，且满足

Ln(xj)=yj, L=0,1, …,n

再由插值多项式的唯一性，得

Pn(x) = Ln(x)

主要仪器设备

HP计算机

实验记录(写出实验内容中的程序代码和运行结果)(可分栏或加页)

拉格朗日插值法：

#i nclude "stdio.h"

int mai n()

{

double m=1.0,a=0.856,l=0;

int i,j;

double x[6]={0.50,0.60,0.70,0.80,0.90,1.00};

double y[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) con ti nue;

m=m*((a-x[j])/(x[i]-x[j]));

}

l+=y[i]*m; m=1;

}

printf(” 结果为％lf",l); return 0;

}

芦果 ^2-588736

Pi*ess 曰n刃 l<e^ to cent inue

最小二乘法：

#i nclude "stdio.h"

#i nclude "math.h"

int mai n()

{double x[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;

int m=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);

double s=sum3*a0+sum4*a1;

r=sum5-s;

prin tf("y=a0+a1*x\n");

prin tf("aO=%f a1= %f\t\n",a0,a1,r);

double q=0.856,p;

p=a0+a1*q;

prin tf("y=%f\n",p);

return 0;

}

a0=0_878261 al=1-978261 ly=2.571652

Press any key to continue

心得体会

拉格朗日插值的优点是插值多项式特别容易建立，缺点是增加节点是原有多项

式不能利用，必须重新建立，即所有基函数都要重新计算，这就造成计算量的浪费。

牛顿插值多项式的优点是增加节点时，原先的差商仍旧不变，仍可以使用。

数据拟合的具体作法是：对给定的数据( Xi ,yi) (i=0,1,…,m),在取定的函

数类中，求p(x)属于此函数类，使误差 m=p(xi)-y i (i=0,1,…,m)的平方和最小，

即：

刀 ri2=E(E p(xj-yi) 2=min

从几何意义上讲，就是寻求与给定点( Xi , yi) (i=0,1，…,m)的距离平方和为最小

的曲线y=p(x)。

相关热词搜索： 实验报告 数值 计算方法 实验 太原理工大学数值计算方法实验报告x

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

最新文章

热门文章