数值分析课程设计.docx
《数值分析课程设计.docx》由会员分享，可在线阅读，更多相关《数值分析课程设计.docx（15页珍藏版）》请在冰豆网上搜索。
数值分析课程设计
摘要
数值计算是具有广泛应用的数学分支，在各个学科中有最新的应用，与MATLAB软件紧密联系，揭示了一种技术或算法可以利用少量的数学知识就能在科技设计中获得巨大的回报。
本文建立在数值分析的理论基础上，能够在Matlab环境中运行，给出了理论分析、程序清单以及计算结果。
以4个主要思想为核心进行展开，首先列主元Gauss消去法求解方程组；然后对一组数据用Lagrange插值的方法求出Lagrange插值多项式；再运用Romberg求积方法，对给出定积分进行积分；最后用Runge-Kutta方法求解常微分方程数值解；在文章的最后用一道综合题结束本文。
“本书结构清晰，条理分明，实例范围较广泛。
它突出了数值分析的中心主题，给出了大量的算法，尤其是，它提供了取自现实生活的实例。
关键词：
列主元Gauss消去，Langrange插值，Romberg积分，Runge-Kutta方法
目录
实验一列主元Gauss消去法1
1实验目的与要求1
2基本原理1
3实验内容及数据来源1
4实验结论3
实验二Lagrange插值4
1实验目的与要求4
2基本原理4
3实验内容及数据来源4
4实验结论6
实验三龙贝格求积公式7
1实验目的与要求7
2基本原理7
3实验内容及数据来源7
4实验结论8
实验四Runge-Kutta方法9
1实验目的与要求9
2基本原理9
3实验内容及数据来源10
4实验结论11
实验五Gauss-Seidel迭代法实例12
1实验目的与要求12
2基本原理12
3实验内容及数据来源12
4实验结论14
心得15
参考文献16
实验一列主元Gauss消去法
1实验目的与要求
熟悉列主元Gauss消去法的基本原理；
会用matlab程序实现列主元Gauss消去法的全过程。
2基本原理
为了提高计算的数值稳定性,在消元过程中采用选择主元的方法.常采用的是列主元消去法。
给定线性方程组Ax=b,记A
（1）=A,b
（1）=b,列主元Gauss消去法的具体过程如下:
首先在增广矩阵B
（1）=（A
（1）,b
（1））的第一列元素中,取
为主元素，
。
然后进行第一步消元得增广矩阵B
（2）=（A
（2）,b
（2））。
再在矩阵B
（2）=（A
（2）,b
（2））的第二列元素中,取
为主元素，
。
然后进行第二步消元得增广矩阵B（3）=（A（3），b（3））。
按此方法继续进行下去，经过n-1步选主元和消元运算，得到增广矩阵B（n）=（A（n）,b（n））。
则方程组A（n）x=b（n）是与原方程组等价的上三角形方程组，可进行回代求解。
易证,只要|A|
0，列主元Gauss消去法就可顺利进行。
可见，列主元Gauss消去法是在每一步消元前,在主元所在的一列选取绝对值最大的元素作为主元素.而全主元Gauss消去法是在每一步消元前,在所有元素中选取绝对值最大的元素作为主元素.但由于运算量大增,实际应用中并不经常使用.
3实验内容及数据来源
下面先用Matlab编写列主元Gauss消去的函数文件如下：
%列选主元的高斯消去法
functionX=lufact_my（A,B）
%InpiutA是系数矩阵，B是右端项
%Outputx是解
[N,N]=size（A）;
X=zeros（N,1）;
Y=zeros（N,1）;
C=zeros（1,N）;
R=1:
N;
k=1;
whilek<=N-1
%求列中最大的值赋给max
[max1,j]=max（abs（A（k:
N,k）））;
%交换行
C=A（k,:
）;%C为A的第k列的值
A（k,:
）=A（j+k-1,:
）;%将A的第K列赋为最大
A（j+k-1,:
）=C;
d=R（k）;
R（k）=R（j+k-1）;
R（j+k-1）=d;
%主元为0的情况
ifA（k,k）==0
'Aissingular.nouniquesolution'
break
end
%化为上三角
form=k+1:
N
mult=A（m,k）/A（k,k）;
A（m,k）=mult;
A（m,k+1:
N）=A（m,k+1:
N）-mult*A（k,k+1:
N）;
m=m+1;
end
k=k+1;
end
%对右端项做处理，但要保证行的交换相同,要注意R（k）的作用
Y
（1）=B（R
（1））;
fork=2:
N
Y（k）=B（R（k））-A（k,1:
k-1）*Y（1:
k-1）;
end
X（N）=Y（N）/A（N,N）;
fork=N-1:
-1:
1
X（k）=（Y（k）-A（k,k+1:
N）*X（k+1:
N））/A（k,k）;
end
例题：
用列主元Gauss消去法求解方程
的解。
输入：
>>A=[223;476;789];
>>B=[4;7;9];
>>lufact（A,B）
输出：
ans=
-2.3333
-0.3333
3.1111
4实验结论
上述matlab程序可以比较好的实现列主元Gauss消去法，且经过实际例子检验出这种方法是正确的。
列主元Gauss消去法是在每一步消元前,在主元所在的一列选取绝对值最大的元素作为主元素.而全主元Gauss消去法是在每一步消元前,在所有元素中选取绝对值最大的元素作为主元素.但由于运算量大增,实际应用中并不经常使用.
实验二Lagrange插值
1实验目的与要求
熟悉Lagrange插值法的基本原理；
会用matlab程序实现Lagrange插值法；
会用Lagrange插值法求解实际问题。
2基本原理
给定（n+1）个插值节点
和对应的函数值
，利用n次拉格朗日插值多项式公式：
，其中
，可以得到插值区间任意x的函数值y为
。
从公式中可以看出，生成的多项式与用来插值的数据密切相关，数据变化则函数就要重新计算，所以当插值数据特别多的时候，计算量就会比较大。
3实验内容及数据来源
由于Matlab中并没有现成的拉格朗日插值命令，下面先编写拉格朗日的函数文件：
functions=Lagrange（x,y,x0）
nx=length（x）;
ny=length（y）;
ifnx~=ny
warning（'矢量x与y的长度应该相同'）
return
end
m=length（x0）;
fori=1:
m
t=0.0;
forj=1:
nx
u=1.0;
fork=1:
nx;
ifk~=j
u=u*（x0（i）-x（k））/（x（j）-x（k））;
end
end
t=t+u*y（j）;
end
s（m）=t;
end
Return
例题：
编好M文件后，就可以用Lagrange插值函数进行插值计算了，下面根据x和y的数据求0.9出的值。
其中，
x=[0.70.50.40.8]
y=[-0.916291-0.693147-0.356675-0.223144];
输入：
>>x=[0.70.50.40.8];
>>y=[-0.916291-0.693147-0.356675-0.223144];
>>Lagrange（x,y,0.9）
输出：
ans=
1.3930
下面重新给定数据：
求解3.6,4.9,7.3三个点上的值。
其中，
x=[578910]
y=[107641];
输入：
>>x=[578910];
>>y=[107641];
>>Lagrange（x,y,[3.44.97.3]）
输出：
ans=
006.7729
4实验结论
从上图中可以看出，拉格朗日插值的一个特点：
拟合出的多项式通过每一个测量数据点。
由于Lagrange插值的n次插值基函数
lk（x）（k=0,1,2,…n）都依赖于全部插值结点，利用公式很容易得到插值多项式，但在增加或减少结点时，插值基函数
lk（x）（k=0,1,2,…n）也随之变化，必须全部重新计算。
实验三龙贝格求积公式
1实验目的与要求
熟悉龙贝格求积公式的基本原理；
会用matlab程序实现龙贝格求积公式；
会用龙贝格求积公式求解定积分。
2基本原理
由于变步长梯形求积法算法简单，但其精度较差，收敛速度较慢，故利用变步长梯形法算法简单的优点，重新形成一个新算法，这就是龙贝格求积公式。
它是对近似值进行修正以后得到的更近似的公式，能自动改变积分步长，以使其相邻两次值的绝对误差或相对误差小于预先设定的允许误差。
龙贝格求积公式为：
3实验内容及数据来源
根据上述实验原理，可以编写龙贝格求积公式的M文件romber.m：
function[R,quad,err,h]=romber（fun,a,b,n,delta）
M=1;
h=b-a;
err=1;
J=0;
R=zeros（4,4）;
R（1,1）=h*（fun（a）+fun（b））/2;
while（（err>delta）&&（JJ=J+1;
h=h/2;
s=0;
forp=1:
M
x=a+h*（2*p-1）;
c=fun（x）;
s=s+c;
end
R（J+1,1）=R（J,1）/2+h*s;
M=2*M;
forK=1:
J
R（J+1,K+1）=R（J+1,K）+（R（J+1,K）-R（J,K））/（4.^K-1）;
end
err=abs（R（J,J）-R（J+1,K+1））;
end
quad=R（J+1,J+1）
例题：
用龙贝格算法计算定积分
在matlab的命令栏输入：
>>fun=@（x）sqrt（4*x-x.^3）;
>>romber（fun,0,1,5,10.^（-10））
输出结果：
quad=
1.2580
ans=
0.866000000
1.11771.20150000
1.20831.23851.2410000
1.24081.25161.25251.252600
1.25241.25621.25651.25661.25660
1.25651.25791.25801.25801.25801.2580
4实验结论
龙贝格求积公式算法简单，高速有效，易于编制程序，适合于计算机上操作。
采用了提高阶数与减少步长这两种提高精度的措施。
能比较准确的求出定积分的值。
但它有一个明显的缺点是：
每当把区间对分后，就要对被积函数f（x）计算它在新分点处的值，而这些函数值的个数是成倍增加的。
实验四Runge-Kutta方法
1实验目的与要求
熟悉Runge-Kutta方法的基本原理；
会用matlab程序实现Runge-Kutta方法；
会用Runge-Kutta方法求解微分方程。
2基本原理
三阶库塔方法（它的截断误差为
）:
用方程
中函数
在前一节点
上取值
的线性组合构造一个表示
的近似公式，从而避免了求
时用
的高阶导数。
该方法有多样推导方法，这里用数值积分法推导。
为此，先将微分方程
略加变形得出：
，对该式两边在相邻两节点
和
之间求积分，移项得出：
，采用不同的方法计算定积分
，便可得出数值积分的不同近似结果。
还可以构造更高阶的龙格—库塔公式，它的一般形式是：
，其中，
3实验内容及数据来源
根据上述实验原理，可以编写Runge-Kutta方法的M文件：
functionR=rk3（f,a,b,ya,N）
h=（b-a）/N;
T=zeros（1,N+1）;
Y=zeros（1,N+1）;
T=a:
h:
b;
Y
（1）=ya;
forj=1:
N
k1=feval（f,T（j）,Y（j））;
k2=feval（f,T（j）+h/3,Y（j）+h*k1/3）;
k3=feval（f,T（j）+2*h/3,Y（j）+2*h*k2/3）;
Y（j+1）=Y（j）+h*（k1+3*k3）/4;
end
R=[T'Y'];
例题：
根据给出的条件求下列微分方程的解：
建立本题数据的M文件：
functionz=f（x,y）
z=y-x/y;
以f.m文件命名。
然后，回到Matlab的命令窗口，
输入：
>>rk3（'f',0,1,1,10）
回车后得到结果：
输出：
ans=
01.0000
0.10001.1003
0.20001.2025
0.30001.3081
0.40001.4187
0.50001.5359
0.60001.6613
0.70001.7965
0.80001.9433
0.90002.1035
1.00002.2791
4实验结论
龙格—库塔方法可以解决大部分微分方程的求解问题，并且可以较好较快的解出微分方程的解。
并且，龙格—库塔方法是由Euler法改进的，故其拥有Euler法的所有优点。
实验五Gauss-Seidel迭代法实例
1实验目的与要求
熟悉Gauss-Seidel迭代法的基本原理；
会用matlab程序实现Gauss-Seidel迭代法；
会用Gauss-Seidel迭代法求解微分方程。
2基本原理
设线性方程组为AX=B，则Gauss-Seidel迭代法的迭代公式：
对应的矩阵表达式为：
其中，
称为Gauss-Seidel迭代矩阵，
上式称为Gauss-Seidel迭代法,简称为G-S迭代法
3实验内容及数据来源
下面编写Gauss-Seidel迭代法的M文件gseid.m：
functionX=gseid（A,B,P,delta,max1）
N=length（B）;
fork=1:
max1
forj=1:
N
ifj==1
X
（1）=（B
（1）-A（1,2:
N）*P（2:
N））/A（1,1）;
elseifj==N
X（N）=（B（N）-A（N,1:
N-1）*（X（1:
N-1））'）/A（N,N）;
else
X（j）=（B（j）-A（j,1:
j-1）*X（1:
j-1）-A（j,j+1:
N）*P（j+1:
N））/A（j,j）;
end
end
err=abs（norm（X'-P））;
relerr=err/（norm（X）+eps）;
P=X';
if（errbreak
end
end
X
例题：
用Gauss-Seidel迭代法求解线性方程组Ax=B,其中，
A=[4-11;4-81;-215];
B=[7-2115]'
输入：
>>A=[4-11;4-81;-215];
>>B=[7-2115]';
>>P=[000]';
>>delta=0.0000001;
>>max1=80;
>>X=gseid（A,B,P,delta,max1）
输出：
X=
2.00004.00003.0000
4实验结论
线性方程组的解法很多，如雅克比迭代法，高斯迭代法以及超松弛迭代法，较雅克比迭代法来看，高斯迭代法收敛较快.同样取精确到小数点后同一位数的近似解,高斯迭代法所需要的迭代次数次比雅克比迭代法需要的迭代次数少。
心得
短短几天的数值分析课程设计很快结束了，通过此次的训练我了解到一些简单的matlab编程方法，丰富了编程知识及常识，了解到了基础知识的重要性，以前认为的matlab都只是用来作图形的，但是课设完了我知道能利matlab解决比较复杂的数值计算问题。
通过我的努力以及遇到不会的问题主动向老师提问，我顺利的完成了这次语言课程设计，感到要想对一门课程了解深刻就必须亲自身入其中的去体验，这样才能发现服，也才能真正的学到知识，也明白了那句“纸上谈来终觉浅，深知此是要躬行自己的很多缺点，才能一一的去克”的真正意义。
纵然把结论书本背得滚瓜烂熟而不运用与实际，也不过是纸上谈兵罢了，终究在未来需要动手能力和实践人才的社会竞争中被淘汰，所以我们须要课程实习，需要工程训练，这样才能成为真正的二十一世纪合格的大学生！
参考文献
[1]石辛民.郝整清.基于Matlab的实用数值计算.清华大学出版社.北京交通大学出版社,2006
[2]苏金朋.阮沈勇.Matlab实用教程.电子工程出版社,2005
[3]张明辉.Matlab6.1..中国水利水电出版社,2008
[4]张池平.计算方法.哈尔滨工业大学出版社,2009
[5]杨万利.数值分析教程.国防工业出版社,2005
升级会员 