东南大学数值分析上机报告完整版.docx

东南大学数值分析上机报告完整版.docx

《东南大学数值分析上机报告完整版.docx》由会员分享，可在线阅读，更多相关《东南大学数值分析上机报告完整版.docx（11页珍藏版）》请在冰豆网上搜索。

word格式-可编辑-感谢下载支持

数值分析上机实验报告

目录

1.chapter1舍入误差及有效数 1

2.chapter2Newton迭代法 3

3.chapter3线性代数方程组数值解法-列主元Gauss消去法 7

4.chapter3线性代数方程组数值解法-逐次超松弛迭代法 9

5.chapter4多项式插值与函数最佳逼近 10

1.chapter1舍入误差及有效数

1.1题目

设SN=j=2N1j2-1，其精确值为。

（1）编制按从大到小的顺序，计算SN的通用程序。

（2）编制按从小到大的顺序，计算SN的通用程序。

（3）按两种顺序分别计算,并指出有效位数。

（编制程序时用单精度）

（4）通过本次上机题，你明白了什么？

1.2编写相应的matlab程序

clear;

N=input（'pleaseinputN:

'）;

AValue=（（3/2-1/N-1/（N+1））/2）;

sn1=single（0）;

sn2=single（0）;

fori=2:

N

sn1=sn1+1/（i*i-1）;%从大到小相加的通用程序%

end

ep1=abs（sn1-AValue）;

forj=N:

-1:

2

sn2=sn2+1/（j*j-1）;%从小到大相加的通用程序%

end

ep2=abs（sn2-AValue）;

fprintf（'精确值为:

%f\n',AValue）;

fprintf（'从大到小的顺序累加得sn=%f\n',sn1）;

fprintf（'从大到小相加的误差ep1=%f\n',ep1）;

fprintf（'从小到大的顺序累加得sn=%f\n',sn2）;

fprintf（'从小到大相加的误差ep2=%f\n',ep2）;

disp（'================================='）;

1.3matlab运行程序结果

>>chaper1

pleaseinputN:

100

精确值为:

0.740050

从大到小的顺序累加得sn=0.740049

从大到小相加的误差ep1=0.000001

从小到大的顺序累加得sn=0.740050

从小到大相加的误差ep2=0.000000

>>chaper1

pleaseinputN:

10000

精确值为:

0.749900

从大到小的顺序累加得sn=0.749852

从大到小相加的误差ep1=0.000048

从小到大的顺序累加得sn=0.749900

从小到大相加的误差ep2=0.000000

>>chaper1

pleaseinputN:

1000000

精确值为:

0.749999

从大到小的顺序累加得sn=0.749852

从大到小相加的误差ep1=0.000147

从小到大的顺序累加得sn=0.749999

从小到大相加的误差ep2=0.000000

1.4结果分析以及感悟

按照从大到小顺序相加的有效位数为：

5,4,3。

按照从小到大顺序相加的和的有效位数为：

6,6,6。

从程序的输出误差结果可以看出，按照不同的顺序相加造成的误差限是不同的，按照从大到小相加的顺序就是一个病态问题，而按照从小到大顺序相加的误差很小，并且在从大到小顺序相加的误差随着n的增大而增大。

因此，采取从小到大的顺序累加得到的结果更加精确。

2.chapter2Newton迭代法

2.1题目

（1）给定初值及容许误差，编制牛顿法解方程f（x）=0的通用程序。

（2）给定方程,易知其有三个根

由牛顿方法的局部收敛性可知存在当时，Newton迭代序列收敛于根x2*。

试确定尽可能大的δ。

试取若干初始值，观察当时Newton序列的收敛性以及收敛于哪一个根。

（3）通过本上机题，你明白了什么？

2.2编写相应的matlab程序

定义f（x）函数

functionF=fu（x）

F=x^3/3-x;

end

定义f（x）的导函数

functionF=dfu（x）

F=x*x-1;

end

求根的通用程序

clear;

x0=input（'请输入初始值x0:

'）;

ep=input（'请输入容许误差：

'）;

flag=1;

whileflag==1

x1=x0-fu（x0）/dfu（x0）;

ifabs（x1-x0）<=ep

flag=0;

end

x0=x1;

end

fprintf（'方程的一个近似解为:

%f\n',x0）;

求sigma的通用程序

clear;

eps=input（'请输入搜索精度:

'）;

ep=input（'请输入容许误差:

'）;

flag=1;

k=0;

x0=0;

whileflag==1;

sigma=k*eps;

x0=sigma;

k=k+1;

m=0;

flag1=1;

whileflag1==1&&m<=10^3

x1=x0-fu（x0）/dfu（x0）;

ifabs（x1-x0）

flag1=0;

end

m=m+1;

x0=x1;

end

ifflag1==1||abs（x0）>=ep

flag=0;

end

end

fprintf（'最大的sigma值为：

%f\n',sigma）;

2.3运行结果

寻找最大的sigma值

主要是在0的基础上，不断的增加步长，带入Newton公式，验证该值是否收敛于0，不断的循环，最后得到最小的不收敛于0的sigma值，此时也为最大满足收敛于0的最大的sigma值。

改变不同的步长，分别得到不同的sigma值，取其中的最小值，即为满足条件的最大的sigma值。

程序相应的运行结果如下：

>>chapter2_2

请输入搜索精度:

10^-6

请输入容许误差:

10^-6

最大的sigma值为：

0.774597

>>chapter2_2

请输入搜索精度:

10^-4

请输入容许误差:

10^-6

最大的sigma值为：

0.774600

>>chapter2_2

请输入搜索精度:

10^-2

请输入容许误差:

10^-6

最大的sigma值为：

0.780000

运行chapter2_1程序

（1）当初值x0属于（-∞，-1）内时，程序运行结果如下，

>>chaper2_1

请输入初始值x0:

-10000

请输入容许误差：

10^-6

方程的一个近似解为:

-1.732051

>>chaper2_1

请输入初始值x0:

-100

请输入容许误差：

10^-6

方程的一个近似解为:

-1.732051

>>chaper2_1

请输入初始值x0:

-10

请输入容许误差：

10^-6

方程的一个近似解为:

-1.732051

>>chaper2_1

请输入初始值x0:

-1.1

请输入容许误差：

10^-6

方程的一个近似解为:

-1.732051

可以得出不论取何值，Newton迭代式收敛，方程的近似解都收敛于-3。

（2）当初值x0属于-1，-δ内时，程序运行结果如下，

>>chaper2_1

请输入初始值x0:

-0.9

请输入容许误差：

10^-6

方程的一个近似解为:

1.732051

>>chaper2_1

请输入初始值x0:

-0.85

请输入容许误差：

10^-6

方程的一个近似解为:

1.732051

>>chaper2_1

请输入初始值x0:

-0.774598

请输入容许误差：

10^-6

方程的一个近似解为:

1.732051

可以得出不论取何值，在此区间上Newton迭代式不收敛。

（3）当初值x0属于（-δ，δ）内时，程序运行结果如下，

>>chaper2_1

请输入初始值x0:

-0.76

请输入容许误差：

10^-6

方程的一个近似解为:

0.000000

>>chaper2_1

请输入初始值x0:

-0.5

请输入容许误差：

10^-6

方程的一个近似解为:

0.000000

>>chaper2_1

请输入初始值x0:

-0.1

请输入容许误差：

10^-6

方程的一个近似解为:

0.000000

>>chaper2_1

请输入初始值x0:

-0.01

请输入容许误差：

10^-6

方程的一个近似解为:

-0.000000

>>chaper2_1

请输入初始值x0:

0.01

请输入容许误差：

10^-6

方程的一个近似解为:

0.000000

>>chaper2_1

请输入初始值x0:

0.1

请输入容许误差：

10^-6

方程的一个近似解为:

0.000000

>>chaper2_1

请输入初始值x0:

0.5

请输入容许误差：

10^-6

方程的一个近似解为:

0.000000

>>chaper2_1

请输入初始值x0:

0.76

请输入容许误差：

10^-6

方程的一个近似解为:

0.000000

可以得出在此区间内，不论取何值，Newton迭代式收敛，方程的近似解都收敛于0。

（4）当初值x0属于（δ，1）内时，运行程序结果如下，

>>chaper2_1

请输入初始值x0:

0.774598

请输入容许误差：

10^-6

方程的一个近似解为:

-1.732051

>>chaper2_1

请输入初始值x0:

0.8

请输入容许误差：

10^-6

方程的一个近似解为:

-1.732051

>>chaper2_1

请输入初始值x0:

0.9

请输入容许误差：

10^-6

方程的一个近似解为:

-1.732051

>>chaper2_1

请输入初始值x0:

0.95

请输入容许误差：

10^-6

方程的一个近似解为:

-1.732051

可以得出不论取何值，在此区间上Newton迭代式不收敛。

（5）当初值x0属于区间1，+∞内时，运行程序结果如下，

>>chaper2_1

请输入初始值x0:

1.1

请输入容许误差：

10^-6

方程的一个近似解为:

1.732051

>>chaper2_1

请输入初始值x0:

10

请输入容许误差：

10^-6

方程的一个近似解为:

1.732051

>>chaper2_1

请输入初始值x0:

100

请输入容许误差：

10^-6

方程的一个近似解为:

1.732051

>>chaper2_1

请输入初始值x0:

10000

请输入容许误差：

10^-6

方程的一个近似解为:

1.732051

可以得出，不管x0取何值，Newton迭代式都收敛，且收敛于根3。

3.chapter3线性代数方程组数值解法-列主元Gauss消去法

3.1题目

对于某电路的分析，归结为求解线性方程组RI=V，其中

（1）编制解n阶线性方程组Ax=b的列主元高斯消去法的通用程序；

（2）用所编程序解线性方程组RI=V，并打印出解向量，保留5位有效数字；

（3）本题编程之中，你提高了哪些编程能力？

3.2程序编写：

%通用Gauss列主元消去法

n=input（'输入线性方程组阶数:

n='）;

b=zeros（1,n）;

A=input（'输入系数矩阵：

A=\n'）;

b（1,:

）=input（'输入线性方程组右端向量：

b=\n'）;

b=b';

C=[A,b];

fori=1:

n-1

[maximum,index]=max（abs（C（i:

n,i）））;%求取C矩阵列主元以及其所在列中位置

index=index+i-1;

T=C（index,:

）;

C（index,:

）=C（i,:

）;

C（i,:

）=T;

fork=i+1:

n

ifC（k,i）~=0

C（k,:

）=C（k,:

）-C（k,i）/C（i,i）*C（i,:

）;

end

end

end

%回代求解

x=zeros（n,1）;

x（n）=C（n,n+1）/C（n,n）;

fori=n-1:

-1:

1

x（i）=（C（i,n+1）-C（i,i+1:

n）*x（i+1:

n,1））/C（i,i）;

end

disp（'该方程组的解为:

'）;

fprintf（'%.5g\n',x）;

3.3运行结果

输入线性方程组阶数:

n=9

输入系数矩阵：

A=

[31-13000-10000;-1335-90-110000;0-931-1000000;00-1079-30000-9;000-3057-70-50;0000-747-3000;00000-304100;0000-50027-2;000-9000-229]

输入线性方程组右端向量：

b=

[-1527-230-2012-7710]

该方程组的解为:

-0.28923

0.34544

-0.71281

-0.22061

-0.4304

0.15431

-0.057823

0.20105

0.29023

3.4结果分析

从程序运行结果可得，该线性方程组的解向量为：

[-0.289230.34544-0.71281-0.22061-0.43040.15431-0.0578230.201050.29023]。

通过该道题程序的编写，我加深了对Gauss列主元消去法的理解，也增强了MATLAB对矩阵、数组处理的能力。

4.chapter3线性代数方程组数值解法-逐次超松弛迭代法

4.1题目

（1）编制解n阶线性方程组Ax=b的SOR方法的通用程序（要求）；

（2）对于35题中所给的线性方程组，取松弛因子，容许误差，打印松弛因子、迭代次数、最佳松弛因子及解向量。

4.2程序编写

R=zeros（9,9）;

R（1,1）=31/2;R（1,2）=-13;R（1,6）=-10;R（2,2）=35/2;R（2,3）=-9;R（2,5）=-11;

R（3,3）=31/2;R（3,4）=-10;R（4,4）=79/2;R（4,5）=-30;R（4,9）=-9;R（5,5）=57/2;

R（5,6）=-7;R（5,8）=-5;R（6,6）=47/2;R（6,7）=-30;R（7,7）=41/2;

R（8,8）=27/2;R（8,9）=-2;R（9,9）=29/2;

R=R+R';R（6,1）=0;R（5,2）=0;

V=[-15;27;-23;0;-20;12;-7;7;10];

eps1=0.5e-5;

detaeps=ones（99,1）;

k=zeros（99,1）;

deta=zeros（9,1）;

x=zeros（9,1）;

x1=ones（9,99）;

x0=ones（9,99）;

fori=1:

99

w=i/50;

while（detaeps（i）>eps1）

forj=1:

9

sum=0;

fort=1:

j-1

sum=sum+R（j,t）*x1（t,i）;

end

fort=j+1:

9

sum=sum+R（j,t）*x1（t,i）;

end

x1（j,i）=（1-w）*x0（j,i）+w*（V（j）-sum）/R（j,j）;

end

k（i）=k（i）+1;

forj=1:

9

deta（j）=abs（x1（j,i）-x0（j,i））;

end

detaeps（i）=max（deta）;

x0=x1;

end

end

p=min（k）;

fori=1:

9

x（i）=x0（i,p）;

end

disp（p）;

disp（x）;

4.3运行结果

最少迭代次数p=10;最佳迭代因子w=0.2；解析向量x=[-0.28920.3455-0.7128-0.2206-0.43040.1544-0.0570.20110.2902]。

4.4结果分析

通过此次的编程，对比不同的迭代因子，迭代因此不同，收敛速度也不相同，再次加强对于SOR的理解，提高MATLAB的编程能力。

5.chapter4多项式插值与函数最佳逼近

5.1题目

（1）编制求第一型3次样条插值函数的通用程序;

（2）已知汽车曲线型值点的数据如下：

0

1

2

3

4

5

6

7

8

9

10

2.51

3.30

4.04

4.70

5.22

5.54

5.78

5.40

5.57

5.70

5.80

端点条件为=0.8，=0.2。

用所编制程序求车门的3次样条插值函数S（x）,并打印出S（i+0.5）（i=0,1,…9）。

5.2程序编写

clc;

clear;

x=[0;1;2;3;4;5;6;7;8;9;10];

y=[2.51;3.3;4.04;4.7;5.22;5.54;5.78;5.4;5.57;5.7;5.8];

dy=[0.8;0.2];

h=zeros（8,1）;

u=zeros（9,1）;

nameda=zeros（9,1）;

d=zeros（11,1）;

mm=zeros（11,1）;

m=zeros（11,11）;

k=0;

fori=1:

10

h（i）=x（i+1）-x（i）;

end

fori=1:

9

u（i）=h（i）/（h（i）+h（i+1））;

nameda（i）=1-u（i）;

end

d

（1）=6*（（y

（2）-y

（1））/h

（1）-dy

（1））/h

（1）;

d（11）=6*（-（y（11）-y（10））/h（10）+dy

（2））/h（10）;

fori=2:

10

d（i）=6*（（y（i+1）-y（i））/h（i）-（y（i）-y（i-1））/h（i-1））/（x（i+1）-x（i-1））;

end

fori=2:

10

m（i,i-1）=u（i-1）;

m（i,i）=2;

m（i,i+1）=nameda（i-1）;

end

m（1,1）=2;m（1,2）=1;m（11,10）=1;m（11,11）=2;

mm=inv（m）*d;

fx=zeros（1,10）;

forj=1:

10

t=input（'请输入0到10之间的一个整数:

'）;

t=t+0.5;

i=fix（t）;

fx（j）=y（i+1）+（（y（i+2）-y（i+1））/h（i+1）-h（i+1）*（1/3*mm（i+1）+1/6*mm（i+2）））*（t-x（i+1））+0.5*mm（i+1）*（t-x（i+1））^2+1/（6*h（i+1））*（mm（i+2）-mm（i+1））*（t-x（i+1））^3;

disp（fx）;

end

sx=zeros（901,1）;

forj=0:

0.01:

9

i=fix（j）;

k=k+1;

sx（k）=y（i+1）+（（y（i+2）-y（i+1））/h（i+1）-h（i+1）*（1/3*mm（i+1）+1/6*mm（i+2）））*（j-x（i+1））+0.5*mm（i+1）*（j-x（i+1））^2+1/（6*h（i+1））*（mm（i+2）-mm（i+1））*（j-x（i+1））^3;

end

5.3运行结果

（1）i分别取0到9之间的整数时，函数输出结果

i

0

1

2

3

4

S（i+0.5）

2.9086

3.6784

4.3815

4.9882

5.3833

i

5

6

7

8

9

S（i+0.5）

5.7237

5.5944

5.4299

5.6598

5.7323

（2）y-x曲线以及拟合曲线

y-x曲线

拟合曲线

5.4结果分析

从以上的两个曲线可以看出，拟合曲线与y-x的曲线基本一致，三次样条函数较好的拟合原函数，通过本次的编程，提高了我对三次样条函数的认识，加深我对三次样条函数的理解。