数值实验报告Word格式.docx

数值实验报告Word格式.docx

《数值实验报告Word格式.docx》由会员分享，可在线阅读，更多相关《数值实验报告Word格式.docx（21页珍藏版）》请在冰豆网上搜索。

a[0]=1.0,b[0]=2.0;

for（i=0;

i<

11;

i++）

{

printf（"

a[i]b[i]x[i]y\n"

）;

x[i]=（b[i]+a[i]）/2;

y=x[i]*x[i]*x[i]+4*x[i]*x[i]-10;

A=a[i]*a[i]*a[i]+4*a[i]*a[i]-10;

if（（y*A）<

0）

a[i+1]=a[i],b[i+1]=x[i];

else

a[i+1]=x[i],b[i+1]=b[i];

printf（"

%f,%f,%f,%f\n"

a[i],b[i],x[i],y）;

}

}

割线法

floatmain（）

floatc,a=1.0,b=2.0;

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.000001）break;

b=c;

}

cout<

<

c;

实验结果和分析

心得体会（遇到的问题和解决方法）

使用不同的方法，可以不同程度的求得方程的解，不同的方法速度不同，求得的结果也稍有区别，当然和要求精度也有关系。

刚开始的时候用数组对二分法进行求解，发现 

循环到第二次 

就无法实现值的传递，于是换了另外一种方法代替了数组。

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

实验目的和要求

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

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

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

实验内容

合理选择利用Gauss消元法、LU分解法、追赶法求解下列方程

高斯分解法：

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

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 

（k=n-1,n-2, 

…,2,1）

LU分解法：

将系数矩阵A转化为A=L*U, 

L为单位下三角矩阵，U为普通上三角矩阵，然后通过解方程组l*y=b,u*x=y,来求解x. 

追赶法:

用来求对角方程组：

L为普通下n-1对角矩阵，U为单位上n-1对角矩阵，然后通过解方程组l*y=b,u*x=y,来求解x.

高斯消元法：

using 

namespace 

std;

float 

main（）

{ 

a[3][4]={{1,2,3,14},{0,1,2,8},{2,4,1,13}};

x[3];

sum=0;

int 

k,i,j;

for（k=0;

k<

2;

k++） 

for（i=k+1;

3;

i++） 

for（j=k+1;

j<

4;

j++） 

a[i][j]=a[i][j]-a[i][k]/a[k][k]*a[k][j];

for（j=0;

a[%d][%d]=%f,"

i,j,a[i][j]）;

x[2]=a[2][3]/a[2][2];

for（k=1;

k>

=0;

k--）

sum+=a[k][j]*x[j];

} 

x[k]=（a[k][3]-sum）/a[k][k];

printf 

（"

x[%d]=%f,"

i+1,x[i]）;

LU分解法:

#include 

stdio.h>

math.h>

#define 

L 

30 

double 

a[L][L], 

b[L], 

l[L][L], 

u[L][L], 

x[L], 

y[L];

main（）{ 

n, 

i, 

j, 

k, 

r;

scanf（ 

"

%d"

&

n 

for 

（ 

i=1;

=n;

++i 

） 

for 

j=1;

++j 

%lf"

a[ 

i 

][ 

j 

] 

（i=1;

++i）{ 

scanf（"

b[i]）;

++i）

（j=1;

++j） 

l[i][j] 

u[i][j] 

=0.0;

（k=1;

++k） 

（j=k;

u[k][j] 

= 

a[k][j];

（r=1;

r<

k;

++r）

-= 

l[k][r] 

* 

u[r][j];

（i=k+1;

++i） 

l[ 

k 

];

（r 

1;

r 

++r）{ 

l[i][k] 

l[i][r] 

u[r][k];

l[i][k] 

/= 

u[k][k];

l[k][k] 

1.0;

y[i] 

b[i];

i;

y[j];

（i=n;

i>

0;

--i） 

x[i]=y[i];

（j=i+1;

x[i] 

x[j];

x[i] 

u[i][i];

printf（ 

%0.2lf\n"

x[ 

return0

追赶法：

FILE*f;

doublea[15],b[15],c[15],d[15];

doublet;

inti,n;

/**********************************************/

f=fopen（"

zgf.txt"

"

r"

fscanf（f,"

&

n）;

%lf%lf%lf"

b[1],&

c[1],&

d[1]）;

for（i=2;

=n-1;

fscanf（f,"

%lf%lf%lf%lf"

a[i],&

b[i],&

c[i],&

d[i]）;

a[n],&

b[n],&

d[n]）;

fclose（f）;

/*********************************************/

c[1]=c[1]/b[1];

d[1]=d[1]/b[1];

t=b[i]-c[i-1]*a[i];

c[i]=c[i]/t;

d[i]=（d[i]-d[i-1]*a[i]）/t;

d[n]=（d[n]-d[n-1]*a[n]）/（b[n]-c[n-1]*a[n]）;

for（i=n-1;

=1;

i--）d[i]=d[i]-c[i]*d[i+1];

printf（"

\n********************************\n"

for（i=1;

d[%2d]=%lf\n"

i,d[i]）;

高斯消元法结果：

LU分解法结果：

追赶法结果：

从消元过程可以看出，对于n阶线性方程组，只要各步主元素不为零，经过n-1步消元，就可以得到一个等价的系数矩阵为上三角形阵的方程组，然后再利用回代过程可求得原方程组的解. 

由于列主元素法相似且优于完全主元素法 

所以省略了后者。

消元过程相当于分解 

A为单位下三角阵L与上三角阵U的乘积，解方程组Ly=b回代过程就是解方程组Ux=y。

其中的L为n阶单位下三角阵、U为上三角阵. 

在 

A 

的LU 

分解中, 

L取下三角阵, 

U 

取单位上三角阵,这样求解方程组Ax=d 

的方法称为追赶法。

另外是追赶法和其他方法求同一方程结果不一样，我多次修改源程序，也不知道原因。

再就是追赶法有很大的局限性 

还待改良

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

（1）了解线性方程组的迭代求解方法，雅可比迭代法或高斯-赛德尔迭代法等典型方法。

（2）帮助学生全面消化已学的相关课程内容，深刻理解计算数值方法课程的内涵，培养使用电子计算机进行科学计算和解决问题的能力。

（3）进行基本技能训练和巩固。

使学生得到选择算法、编写程序、分析数值结果、写数值试验报告、课堂讨论等环节的综合训练

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

实验原理：

把矩阵A分解矩阵N和P， 

其中N为非奇异矩阵，于是，Nx=Px+b 

即x=Bx+f 

据此，可以建立迭代公式x（k+1）=Bx（k）+1;

若序列｛X

（k）

｝收敛, 

lim 

x（k 

）=x* 

显然有 

x*=Bx*+f;

即,极限x*便是所求方程组的解。

#include<

main（） 

x1[20] 

x2[20],x3[20];

x10, 

x20, 

x30;

请输入x1,x2,x3的初值：

\n"

x10,&

x30）;

x1[n] 

x2[n] 

x3[n] 

18;

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;

%5d 

%5lf 

%5lf\n"

i,x1[i],x2[i],x3[i]）;

}

使用高斯-赛德尔和雅克比迭代都可以求出方程组的解，但是利用高斯-赛德尔迭代法所需的迭代次数比雅克比迭代少，能够更早的达到精度要求。

产生误差的原因，舍入误差及数的存储空间是有限的，因此出现误差。

雅可比迭代法显然其速度不如高斯-赛德尔迭代法好。

这是由其算法的设计本身来决定的。

同时高斯-赛德尔算法由于都涉及了两次循环且同时未引入新的变量，因此相比而言高斯赛德尔迭代法较优。

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

（1）了解矩阵特征值与特征向量问题解法，掌握幂法。

（2）加深对矩阵特征值与特征向量问题求解方法的认识，掌握算法。

（3）熟练运用已学计算方法求解方程组 

（4）加深对计算方法技巧，选择正确的计算方法来求解各种方程组 

（5）培养使用电子计算机进行科学计算和解决问题的能力

代数插值：

用拉格朗日插值法或牛顿插值法求解：

已知f（x）在6个点的函数值如下表所示，运用插值方法，求f（0.596）的近似值。

x

0.40

0.55

0.65

0.80

0.90

1.05

f（x）

0.41075

0.57815

0.69675

0.88811

1.02652

1.25386

最小二乘法拟合多项式：

给定数据点（xi,yi），用最小二乘法拟合数据的多项式，并求平方误差。

xi

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

原理：

设函数在区间[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） 

最小二乘法原理：

iostream.h>

I;

char 

L;

M[100][100];

x[100],y[100];

X=1,xx=0,w=1,N=0,P,R=1;

n=5;

//cout<

请输入所求均差阶数：

;

//求所有阶差

//cin>

>

n;

//for（int 

i=0;

//{ 

/*cout<

请输入x"

的值：

endl;

cin>

x[i];

cout<

请输入y"

y[i];

M[i][0]=x[i];

M[i][1]=y[i];

*///

M[0][0]=0.40;

M[0][1]=0.41075;

M[1][0]=0.55;

M[1][1]=0.57815;

M[2][0]=0.65;

M[2][1]=0.69675;

M[3][0]=0.80;

M[3][1]=0.88811;

M[4][0]=0.90;

M[4][1]=1.02652;

M[5][0]=1.05;

M[5][1]=1.25386;

for（ 

j=2;

=n+1;

j++）

for（int 

M[i][j]=（M[i][j-1]-M[i-1][j-1]）/（M[i][0]-M[i-j+1][0]）;

for（int 

其"

阶均差为：

M[i][i+1]<

请输入x的值：

x="

xx;

X*=xx-M[i][0];

N+=M[i+1][i+2]*X;

P=M[0][1]+N;

其函数值：

y="

P<

// 

fstream.h>

N 

15 

power（double 

a,int 

n） 

b=1;

b*=a;

return 

b;

void 

Gauss（）;

X[N],Y[N],sumX[N],sumY[N],a[N][N],b[N],l[N][N],x[N];

s;

i,j,k,n,index;

请输入已知点的个数n="

n=7;

请输入X和Y:

//输入给定数据 

X[0]=0.0;

Y[0]=1.00;

X[1]=0.5;

Y[1]=1.75;

X[2]=0.6;

Y[2]=1.96;

X[3]=0.7;

Y[3]=2.19;

X[4]=0.8;

Y[4]=2.44;

X[5]=0.9;

Y[5]=2.71;

X[6]=1.0;

Y[6]=3.00;

//绑定数据 

//可以解绑由下for循环 

输入任何数据 

sumX[1]+=X[i];

sumY[1]+=Y[i]

sumX[1]="

sumX[1]<

\t"

sumY[1]="

sumY[1]<

请输入拟合次数index="

index;

i=n;

sumX[0]=i;

for（i=2;

=2*index;

sumX[i]=0;

sumX[i]+=power（X[j],i）;

sumX["

]="

sumX[i]<

=index+1;

sumY[i]=0;

sumY[i]+=power（X[j],i-1）*Y[j];

sumY["

sumY[i]<

for（i=1;

//建立正规方程组 

for（j=1;

a[i][j]=sumX[i+j-2];

b[i]=sumY[i];

k=1;

//用高斯消元法解方程组 

do{ 

l[j][k]=a[j][k]/a[k][k];

for（i=k+1;

i++）{ 

a[i][j]=a[i][j]-l[i][k]*a[k][j];

b[i]=b[i]-l[i][k]*b[k];

if（k==index+1） 

break;

k++;

}while

（1）;

x[index+1]=b[index+1]/a[index+1][ind