数值分析实验指导.docx
《数值分析实验指导.docx》由会员分享，可在线阅读，更多相关《数值分析实验指导.docx（16页珍藏版）》请在冰豆网上搜索。
数值分析实验指导
《数值分析》实验指导
2011.09.01
前言
本实验指导是依据信息与计算机科学专业本科教学计划编写，与教材讲义配合使用，旨在引导学生使用计算机开展数值试验，掌握数值算法和程序设计的基本原理和技能。
学生通过选择算法，编写程序，分析数值结果，写数值实验报告，课堂讨论等环节的综合训练，从而逐步掌握数值试验的方法和技巧，获得多方面的计算经验。
教学计划中共设16学时实验，各实验的名称及学时如下
序号
实验项目
学时
备注
1
实验一、插值法
2
设计性
2
实验二、数值积分
（一）
2
设计性
3
实验三、三点Gauss公式
2
设计性
4
实验四、方程求根的迭代法
2
设计性
5
实验五、Choleskey方法
2
设计性
6
实验六、Gauss消去法求解线性方程组
2
设计性
7
实验七、常微分方程数值解法
2
设计性
8
实验八、线性方程组的迭代法
2
设计性
实验一插值法
实验目的
一、拉格朗日插值法
1、 掌握拉格朗日插值法的基本理论。
2、会编写拉格朗日插值法的通用子程序，解决一些基本问题。
二、牛顿插值法
1、掌握牛顿插值法的基本思路和步骤；
2、培养编程与上机调试能力。
实验内容
1、已知函数数据如下：
0.20.240.280.320.360.40
0.198670.23770.27640.31460.35230.3894
（1）计算各阶差商；
（2）利用拉格朗日插值多项式
和牛顿插值多项式
计算近似函数
，
参考实验步骤
1、计算Lagrange插值多项式在x=x0处的值
function[y0,N]=Lagrange_eval（X,Y,x0）
%X,Y是已知插值节点坐标
%x0是插值点
%y0是Lagrange插值多项式在x0处的值
%N（i）是Lagrange插值基函数
m=length（X）;
N=zeros（m,1）;
y0=0;
fori=1:
m
N（i）=1;
forj=1:
m
ifj~=i
N（i）=N（i）*（x0-X（j））/（X（i）-X（j））;
end
end
y0=y0+Y（i）*N（i）
end
2、牛顿插值法计算步骤
1．输入
值及（
；要计算的函数点
。
2．对给定的
由
计算
的值。
3．输出
。
实验二数值积分
（一）
实验目的
掌握复化梯形公式和复化Simpson公式，会编写程序，上机进行实例计算。
实验内容
编写复化梯形公式和复化Simpson公式通用子程序，计算
，步长
实验步骤
一、程序范例
1、复化梯形公式计算积分
functionS=FT（f,a,b,N）
%f表示被积函数的函数名
%a,b表示被积区间[a,b]的端点
%N表示区间个数
%S是用复化梯形公式求得的积分值
h=（b-a）/N;
fa=feval（f,a）;
fb=feval（f,b）;
S=fa+fb;
x=a;
fori=1:
N-1
x=x+h;
fx=feval（f,x）;%小区间端点函数值
S=S+2*fx;
end
S=h*S/2;
2、定步长复化Simpson公式计算积分
functionS=FSimpson（f,a,b,N）
%f表示被积函数的函数名
%a,b表示被积区间[a,b]的端点
%N表示区间个数
%S是用复化Simpson公式求得的积分值
h=（b-a）/N;
fa=feval（f,a）;
fb=feval（f,b）;
S=fa+fb;
x=a;
fori=1:
N
x=x+h/2;
fx=feval（f,x）;%小区间中点函数值
S=S+4*fx;
x=x+h/2;
fx=feval（f,x）;%小区间右端点函数值
S=S+2*fx;
end
S=h*S/6;
实验三三点Gauss公式
实验目的
利用三点Gauss公式计算被积函数
在给定区间上的积分值。
实验内容
利用三点Gauss公式计算积分
。
实验步骤
一、程序范例
三点Gauss公式
FunctionG=TGauss（f,a,b）
%f表示被积函数句柄
%a，b表示被积函数区间
的端点
%G是用三点Gauss公式求得的积分值
x1=（a+b）/2-sqrt（3/5）*（b-a）/2;
x2=（a+b）/2+sqrt（3/5）*（b-a）/2;
G=（b-a）*（5*feval（f,x1）/9+8*feval（f,（a+b）/2）/9+5*feval（f,x2）/9）/2;
思考问题：
考虑两点Gauss公式计算积分值。
实验四方程求根的迭代法
实验目的
用二分法求解非线性方程
在区间
内的根。
实验内容
求方程
在区间
内的实根，使精度达到
。
实验步骤
一、程序范例
二分法
Function[x,k]=deminimetbod（a,b,emg）
%a,b表示求解区间
的端点
%f表示所求方程的函数名
%emg是精度指标
%x表示所求近似解
%k表示循环次数
fa=feval（f,a）;
fab=feval（f,（a+b）/2）;
k=0;
x=a;
whileabs（b-a）>emg
iffab==0
x=（a+b）/2;
return;
elseiffa*fb<0
b=（a+b）/2;
else
a=（a+b）/2;
end
fa=feval（f,a）;
fab=feval（f,（a+b）/2）;
k=k+1;
end
x=（a+b）/2;
实验五Choleskey方法
实验目的
用Choleskey分解法求解对称线性方程组
。
实验内容
用Choleskey分解法求解对称线性方程组
实验步骤
一、程序范例
Choleskey方法
Functionx=Chol_decompose
%用Choleskey分解法求解对称线性方程组
%A是对称矩阵
%L是单位下三角矩阵
%D是对角阵
%对矩阵A进行三角分解：
A=LDLT
N=length（A）;
L=zeros（N,N）;D=aeros（1,N）;
fori=1:
N
L（i,i）=1;
end
D
（1）=A（1,1）;
fori=2:
N
forj=1:
i-1
ifj==1
L（i,j）=A（i,j）/D（j）;
else
sum1=0;
fork=1:
j-1
sum1=sum1+L（i,k）*D（k）*L（j,k）;
end
L（i,j）=（A（i,j）-sum1）/D（j）;
end
end
sum2=0;
fork=1:
i-1
sum2=sum2+L（I,k）^2*D（k）;
end
D（i）=A（I,i）-sum2;
end
%分别求解方程组Ly=b;LTx=y/D
y=zeros（1,N）;
y
（1）=b
（1）;
fori=2:
N
sumi=0;
fork=1:
i-1
sumi=sumi+L（i,k）*x（k）;
end
y（i）=b（i）-sumi;
end
x=zeros（i,N）;
x（N）=y（N）/D（N）;
fori=N-1:
-1:
1
sumi=0;
fork=i+1:
N
sumi=sumi+L（k,i）*x（k）;
end
x（i）=y（i）/D（i）-sumi;
end
思考问题：
考虑用追赶法求解三对角方程组。
实验六Gauss消去法求解线性方程组
实验目的
1、 使学生掌握求解线性方程组的直接解法的基本理论。
2、使学生借助于计算机这一高效计算工具和本课程中的基本数值计算方法来解决某些数学问题，获取其近似解，培养学生的科学计算能力。
实验内容
1、用Gauss消去法求解线性方程组，上机解
实验步骤
一、程序范例
1、Gauss消去法
functionx=Gauss_Elimination（A,b）
%用Gauss消去法求解线性方程组Ax=b
%x是未知向量
n=length（b）;
x=zeros（n,1）;
d1=0;
%消去过程
fork=1:
n-1;
fori=k+1:
n
c=A（i,k）/A（k,k）;
forj=k+1:
n
A（i,j）=A（i,j）-c*A（k,j）;
end
b（i）=b（i）-b（k）*c;
A（i,k）=0;
end
end
%回代求解
x（n）=b（n）/A（n,n）
fori=n-1:
-1:
1
sum=0;
forj=i+1:
n
sum=sum+A（i,j）*x（j）;
end
x（i）=（b（i）-sum）/A（i,i）
end
结果：
ans=
思考问题：
1、考虑全主元Gauss消去法求解方程组
实验七常微分方程数值解法
实验目的
掌握Euler法及改进的Euler算法，会编写程序，
实验要求：
编写Euler法及改进的Euler算法通用子程序，并用于解下列问题的近似解
实验步骤
一、程序范例
1、Euler算法求解初值问题
functionE=Euler（f,a,b,N,ya）
%f1是微分方程右端函数名
%a,b是自变量的取值区间[a,b]的端点
%N是区间等分的个数
%ya表示初值y（a）
%E=[x',y']是自变量X和解Y所组成的矩阵
h=（b-a）/N;
y=zeros（1,N+1）;
x=zeros（1,N+1）;
y
（1）=ya;
x=a:
h:
b;
fori=1:
N
y（i+1）=y（i）+h*feval（f,x（i）,y（i））;
end
E=[x',y'];
2、改进的Euler算法求解初值问题
functionE=MendEuler（f,a,b,N,ya）
%f是微分方程右端函数名
%a,b是自变量的取值区间[a,b]的端点
%N是区间等分的个数
%ya表示初值y（a）
%E=[x',y']是自变量X和解Y所组成的矩阵
h=（b-a）/N;
y=zeros（1,N+1）;
x=zeros（1,N+1）;
y
（1）=ya;
x=a:
h:
b;
fori=1:
N
y1=y（i）+h*feval（f,x（i）,y（i））;
y2=y（i）+h*feval（f,x（i+1）,y1）;
y（i+1）=（y1+y2）/2;
end
E=[x',y'];
思考问题：
1、考虑梯形公式和后退Euler法求解常微分方程
2、用4阶R-K法求解常微分方程的初值问题.
实验八线性方程组的迭代法
实验目的
1、 使学生掌握求解线性方程组的迭代法的基本理论。
2、使学生借助于计算机这一高效计算工具和本课程中的基本数值计算方法来解决某些数学问题，获取其近似解，培养学生的科学计算能力。
实验内容
分别用Jacobi迭代法和Gauss-Seidel迭代法解
实验步骤
一、程序范例
1、Jacobi迭代法
Function[x,k]=Jacobimethod（A,b,x0,N,emg）
%A是线性方程组的系数矩阵
%b是右端向量
%x0是初始估计
%N表示迭代次数上限，若迭代次数大于N，则迭代失败
%emg表示精度
%用Jacobi迭代法求解线性方程组A*x=b的解
%k表示迭代次数
%x表示用迭代法求得的线性方程组的近似解
n=length（A）;
x1=zero（n,1）;x2=zero（n,1）;
x1=x0;k=0;
r=max（abs（b-A*x1））;
whiler>emg
fori=1:
n
sum=0;
forj=1:
n
ifi~=j
sum=sum+A（i,j）*x1（j）;
end
end
x2（i）=（b（i）-sum）/A（i,i）;
end
r=max（abs（x2-x1））;
x1=x2;
k=k+1;
ifk>N
disp（‘迭代失败，返回’）;
return;
end
end
x=x1;
2、Gauss-Seidel迭代法
Function[x,k]=Gaussmethod（A,b,x0,N,emg）
%A是线性方程组的系数矩阵
%b是右端向量
%x0是初始估计
%N表示迭代次数上限，若迭代次数大于N，则迭代失败
%emg表示精度
%用Gauss-Seideli迭代法求解线性方程组A*x=b的解
%k表示迭代次数
%x表示用迭代法求得的线性方程组的近似解
n=length（A）;
x1=zero（n,1）;x2=zero（n,1）;
x1=x0;k=0;
r=max（abs（b-A*x1））;
whiler>emg
fori=1:
n
sum=0;
forj=1:
n
ifj>i
sum=sum+A（i,j）*x1（j）;
elseifjsum=sum+A（i,j）*x2（j）;
end
end
x2（i）=（b（i）-sum）/A（i,i）;
end
r=max（abs（x2-x1））;
x1=x2;
k=k+1;
ifk>N
disp（‘迭代失败，返回’）;
return;
end
end
x=x1;
思考问题：
1、考虑超松弛法求解线性方程组。
数值分析的实验日期分别写１２月２日
朱凤双（1576961540） 14:
25:
56
５日
朱凤双（1576961540） 14:
26:
02
９日
朱凤双（1576961540） 14:
26:
05
１２日
朱凤双（1576961540） 14:
26:
10
１５日
朱凤双（1576961540） 14:
26:
15
１９日
朱凤双（1576961540） 14:
26:
30
２５日
朱凤双（1576961540） 14:
26:
34
２６日
2、
升级会员 