数值分析课程实验报告
姓 名：
学 号：
学 院： 机 电 学 院
日 期： 2015 年 X 月X 日
目 录
实验 函数插值方法 1
实验二 函数逼曲线拟合 5
实验三 数值积分数值微分 7
实验四 线方程组直接解法 9
实验五 解线性方程组迭代法 15
实验六 非线性方程求根 19
实验七 矩阵特征值问题计算 21
实验八 常微分方程初值问题数值解法 24
实验 函数插值方法
问题提出
定元函数n+1节点值试Lagrange公式求插值项式分段二次Lagrange插值项式
数：
（1）
04
055
065
080
095
105
041075
057815
069675
090
100
125382
求五次Lagrange项式分段三次插值项式计算
值（提示：结果 ）
（2）
1
2
3
4
5
6
7
0368
0135
0050
0018
0007
0002
0001
试构造Lagrange项式计算值（提示：结果 ）
二求
1 利Lagrange插值公式
编写出插值项式程序
2 出插值项式分段三次插值项式表达式
3 根节点选取原问题（2）三点插值二点插值结果
4 插值问题Newton插值项式结果Newton插值项式：
中：
三目意义
1 学会常插值方法求函数似表达式解决实际问题
2 明确插值项式分段插值项式优缺点
3 熟悉插值方法程序编制
4 果绘出插值函数曲线观察光滑性
四实验步骤
（1）
04
055
065
080
095
105
041075
057815
069675
090
100
125382
求五次Lagrange项式分段三次插值项式计算
值（提示：结果 ）
第步：先matlab中定义lagranM文件拉格朗日函数代码：
function[cl]lagran(xy)
wlength(x)
nw1
lzeros(ww)
for k1n+1
v1
for j1n+1
if(k~j)
vconv(vpoly(x(j)))(x(k)x(j))
end
end
l(k)v
end
cy*l
end
第二步：然matlab命令窗口输入：
>>>> x[04 055 065 080095 105]y[041075 057815 069675 090 100 125382]
>> lagran(xy)
回车：
ans 1216264 4227503 5725667 3772549 1219718 150845
出求拉格朗日项式
p（x）1216264x54227503x4+5725667x33772549x2+1219718x150845
第三步：编辑窗口输入命令：
>> x[04 055 065 080095 105]
>> y1216264*x^54227503*x^4+5725667*x^33772549*x^2+1219718*x150845
>> plot(xy)
命令执行图示图形然
>> x0596
>> y1216264*x^54227503*x^4+5725667*x^33772549*x^2+1219718*x15084
y 06262
f（0596）06262
理f（099）10547
（2）
1
2
3
4
5
6
7
0368
0135
0050
0018
0007
0002
0001
试构造Lagrange项式分段三次插值项式计算值（提示：结果 ）
实验步骤：
第步定义
function[cl]lagran(xy)
wlength(x)
nw1
lzeros(ww)
for k1n+1
v1
for j1n+1
if(k~j)
vconv(vpoly(x(j)))(x(k)x(j))
end
end
l(k)v
end
cy*l
end
定义完拉格朗日M文件
第二步：然matlab命令窗口输入：
>>>> x[1 2 3 4 5 6 7] y[0368 0135 0050 0018 0007 0002 0001]
>> lagran(xy)
回车：
ans 00001 00016 00186 01175 04419 09683 09950
出求拉格朗日项式
p（x）00001x600016x5+00186x401175x3+04419x209683x+09950
第三步：编辑窗口输入命令：
>> x[1 2 3 4 5 6 7]
>> y00001*x^600016*x^5+00186*x^401175*x^3+04419*x^209683*x+09950
>> plot(xy)
命令执行图示图形然
>> x18
>> y1216264*x^54227503*x^4+5725667*x^33772549*x^2+1219718*x15084
y 01650
f（0596）06262
理f（615）23644
五实验结
插值离散数基础补插连续函数条连续曲线通全部定离散数点离散函数逼重方法利通函数限点处取值状况估算出函数点处似值
实验二 函数逼曲线拟合
问题提出
机数中找出规律性出似表达式问题生产实践科学实验中量存通常利数二法求拟合曲线
某冶炼程中根统计数含碳量时间关系试求含碳量时间t拟合曲线
t(分)
0 5 10 15 20 25 30 35 40 45 50 55
0 127 216 286 344 387 415 437 451 458 402 464
二求
1二法进行曲线拟合
2似解析表达式
3印出拟合函数印出误差
4外选取似表达式尝试拟合效果较
5* 绘制出曲线拟合图
三目意义
1掌握曲线拟合二法
2二法解超定线代数方程组
3探索拟合函数选择拟合精度间关系
四实验步骤：
第步先写出线性二法M文件
function clspoly(xym)
nlength(x)
bzeros(1m+1)
fzeros(nm+1)
for k1m+1
f(k)x^(k1)
end
af＇*f
bf＇*y＇
ca＼b
cflipud(c)
第二步命令窗口输入：
>>lspoly([0510152025303540455055][0127216286344387415437451458402464]2)
回车：
ans
00024
02037
02305
求拟合曲线y00024x2+02037x+02305
编辑窗口输入命令：
>> x[0510152025303540455055]
>> y00024*x^2+02037*x+02305
>> plot(xy)
命令执行图
五实验结
分析复杂实验数时常采分段曲线拟合方法利方法段实现佳逼段边界满足连续性导性分段函数光滑算法出相应误差分析出该方法分段曲线拟合中应方法凸轮实验数动分段拟合
实验三 数值积分数值微分
问题提出
选复合梯形公式复合Simpson公式Romberg算法计算
（1）
（2）
（3）
（4）
二求
1 编制数值积分算法程序
2 分两种算法计算积分较结果
3 分取步长试较计算结果（n 10 20等）
4 定精度求ε试变步长算法确定佳步长
三目意义
1 深刻认识数值积分法意义
2 明确数值积分精度步长关系
3 根定积分计算方法考虑二重积分计算问题
四 实验步骤
第步：编写种积分程序
复合梯形程序：
function ITX(xy)
nlength(x)mlength(y)
if n~m
error(＇The lengths of X and Y must be equal＇)
return
end
h(x(n)x(1))(n1)
a[1 2*ones(1n2) 1]
Ih2*sum(a*y)
复合Simpson程序：
function s simpr1(fabn)
h(ba)(2*n)
s10
s20
for k110
xa+h*(2*k1)
s1s1+feval(fx)
end
for k1(101)
xa+h*2*k
s2s2+feval(fx)
end
sh*(feval(fa)+feval(fb)+4*s1+2*s2)3
end
Romberg程序：
function I Romber_yang(funabep)
if nargin<4
ep1e5end
m1 hba
Ih2*(feval(funa)+feval(funb))T(11)I
while 1
N2^(m1)hh2II2
for i1N
II+h*feval(funa+(2*i1)*h)
end
T(m+11)IM2*Nk1
while M>1
T(m+1k+1)(4^k*T(m+1k)T(mk))(4^k1)
MM2kk+1
end
if abs(T(kk)T(k1k1))break
end
mm+1
end
IT(kk)
第二步：复合梯形公式复合Simpson公式Romberg公式积分进行计算
1 积分Ι0144sin2xdx梯形积分T049871011辛普森积分S049871111Romberg积分R049871111
2 积分Ι01sinXXdxf(0)1梯形积分T094607307辛普森积分S094607308Romberg积分R094607307
3 积分Ι01eX4+X2dx梯形积分T039081248辛普森积分S039081185Romberg积分R039081885
4 积分Ι01ln1+X1+X2dx梯形积分T027218912辛普森积分S027219844Romberg积分R027219827
五实验结
通实验学会复合梯形公式复合Simpson公式Romberg公式编程应掌握MATLAB提供计算积分种函数方法
实验四 线方程组直接解法
问题提出
出列类型线性方程组请适算法计算解
1 设线性方程组
2 设称正定阵系数阵线方程组
3 三角形线性方程组
二求
1 述三方程组分利Gauss序消法Gauss列元消法方根法改进方根法追赶法求解（选择）
2 应结构程序设计编出通程序
3 较计算结果分析数值解误差原
4 利相应模块输出系数矩阵三角分解式
三目意义
1通该课题实验体会模块化结构程序设计方法优点
2运学计算方法解决类线性方程组直接算法
3提高分析解决问题力做学致
4 通三角形线性方程组解法体会稀疏线性方程组解法特点
四实验步骤：
列元高斯消法matlabM文件程序
function [xdetindex]Gauss(Ab)
求线形方程组列元Gauss消法中
A方程组系数矩阵
b方程组右端项
x方程组解
det系数矩阵A行列式值
index指标变量index0表示计算失败index1表示计算成功
[nm]size(A) nblength(b)
方程组行列维数相等时停止计算输出出错信息
if n~m
error(＇The rows and columns of matrix A must be equal＇)
return
end
方程组右端项维数匹配时停止计算输出出错信息
if m~nb
error(＇The columns of A must be equal the length of b＇)
return
end
开始计算先赋初值
index1det1xzeros(n1)
for k1n1
选元
a_max0
for ikn
if abs(A(ik))>a_max
a_maxabs(A(ik))ri
end
end
if a_max<1e10
index0return
end
交换两行
if r>k
for jkn
zA(kj)A(kj)A(rj)A(rj)z
end
zb(k)b(k)b(r)b(r)zdetdet
end
消元程
for ik+1n
mA(ik)A(kk)
for jk+1n
A(ij)A(ij)m*A(kj)
end
b(i)b(i)m*b(k)
end
detdet*A(kk)
end
detdet*A(nn)
回代程
if abs(A(nn))<1e10
index0return
end
for kn11
for jk+1n
b(k)b(k)A(kj)*x(j)
end
x(k)b(k)A(kk)
end
然命令窗口输入
>> A[4 2 3 1 2 1 0 0 0 08 6 5 3 6 5 0 1 0 04 2 2 1 3 2 1 0 3 10 2 1 5 1 3 1 1 9 44 2 6 1 6 7 3 3 2 38 6 8 5 7 17 2 6 3 50 2 1 3 4 2 5 3 0 116 10 11 9 17 34 2 1 2 24 6 2 7 13 9 2 0 12 40 0 1 8 3 24 8 6 3 1]
>> b[5 12 3 2 3 46 13 38 19 21]
>> gauss(Ab)
ans
10000
10000
00000
10000
20000
00000
30000
10000
10000
20000
高斯约消法maltabM文件程序
function [xflag]Gau_Jor(Ab)
求线形方程组列元Gauss约法消法中
A方程组系数矩阵
b方程组右端项
x方程组解
[nm]size(A) nblength(b)
方程组行列维数相等时停止计算输出出错信息
if n~m
error(＇The rows and columns of matrix A must be equal＇)
return
end
方程组右端项维数匹配时停止计算输出出错信息
if m~nb
error(＇The columns of A must be equal the length of b＇)
return
end
开始计算先赋初值
flag＇ok＇xzeros(n1)
for k1n
选元
max10
for ikn
if abs(A(ik))>max1
max1abs(A(ik))ri
end
end
if max1<1e10
falg＇failure＇return
end
交换两行
if r>k
for jkn
zA(kj)A(kj)A(rj)A(rj)z
end
zb(k)b(k)b(r)b(r)z
end
消元程
b(k)b(k)A(kk)
for jk+1n
A(kj)A(kj)A(kk)
end
for i1n
if i~k
for jk+1n
A(ij)A(ij)A(ik)*A(kj)
end
b(i)b(i)A(ik)*b(k)
end
end
end
输出x
for i1n
x(i)b(i)
end
然保存命令窗口输入：
>> A[4 2 4 0 2 4 0 02 2 1 2 1 3 2 04 1 14 1 8 3 5 60 2 1 6 1 4 3 32 1 8 1 22 4 10 34 3 3 4 4 11 1 40 2 5 3 10 1 14 20 0 6 3 3 4 2 19]
>> b[0 6 20 23 9 22 15 45]
>> Gau_Jor(Ab)
ans
1211481
1401127
297515
601528
109120
267963
54259
20185
五实验结
LU法调matlab中函数lu中L三角直接计算结果计算进行行变换果进行行变b变样会麻烦
实验五 解线性方程组迭代法
问题提出
实验四列目意义线性方程组试分选Jacobi 迭代法GaussSeidel迭代法SOR方法计算解
二求
1体会迭代法求解线性方程组消法做较
2分精度求迭代次数体会该迭代法收敛快慢
3方程组23SOR方法时选取松弛子ω080911112等试算法收敛性影响找出选松弛子佳者
4出种算法设计程序计算结果
三目意义
1通机计算体会迭代法求解线性方程组特点消法较
2运学迭代法算法解决类线性方程组编出算法程序
3体会机计算时终止步骤k >（予迭代次数）迭代法敛散性意义
4 体会初始解松弛子选取计算结果影响
四实验步骤
第步编写实验需Jacobi迭代法GaussSeidel迭代法SOR迭代法程序
Jacobi迭代法：
function [xkindex]J(Abepitmax)
if nargin<4
itmax100
end
if nargin<3 ep1e5
end
nlength(A)
k0
xzeros(n1)yzeros(n1)index1
while 1
for i1n
y(i)b(i)
for j1n
if j~i
y(i)y(i)A(ij)*x(j)
end
end
if abs(A(ii))<1e10|kitmax
index0return
end
y(i)y(i)A(ii)
end
if norm(yxinf) break
end
xykk+1
end
GaussSeidel迭代法：
function [xkindex]G(Abepitmax)
if nargin<4
itmax100
end
if nargin<3 ep1e5
end
nlength(A)
k0
xzeros(n1)yzeros(n1)index1
while 1
yx
for i1n
zb(i)
for j1n
if j~i
zzA(ij)*x(j)
end
end
if abs(A(ii))<1e10|kitmax
index0return
end
zzA(ii)x(i)z
end
if norm(yxinf) break
end
kk+1
end
SOR迭代法：
function [xkindex]SOR(Abepwitmax)
if nargin<5
itmax100
end
if nargin<4 w1end
if nargin<3 ep1e5
end
nlength(A)
k0
xzeros(n1)yzeros(n1)index1
while 1
yx
for i1n
zb(i)
for j1n
if j~i
zzA(ij)*x(j)
end
end
if abs(A(ii))<1e10|kitmax
index0return
end
zzA(ii)x(i)(1w)*x(i)+w*z
end
if norm(yxinf) break
end
kk+1
end
第二步：实验出方程代入程序计算结果
1 设线性方程组
2 设称正定阵系数阵线方程组
3三角形线性方程组
五实验结
迭代法解线性方程组重实方法特适求解实际中量出现系数矩阵稀疏阵型线性方程组通次实验学会Jacobi迭代法GaussSeidel迭代法SOR迭代法程序编写掌握优缺点适条件
实验六 非线性方程求根
问题提出
设方程三实根
现采面六种计算格式求 f(x)0根
1
2
3
4
5
6
二求
1编制程序进行运算印出种迭代格式敛散情况
2事误差估计控制迭代次数印出迭代次数
3初始值选取迭代收敛影响
4分析迭代收敛发散原
三目意义
1通实验进步解方程求根算法
2认识选择计算格式重性
3掌握迭代算法精度控制
4明确迭代收敛性初值选取关系
四实验步骤
第步：编写实验需程序
function [x_starindexit]DD(funxepitmax)
if nargin<4
itmax100
end
if nargin<3
ep1e5
end
index0k1
while k x1xxfeval(funx)
if abs(xx1) index1break
end
kk+1
end
x_starxitk
第二步：分述六种形式表达式计算方程根结果
1 x10x20
2 x1穷x203473
3 x118794x218794
4 x103473x203473
5 x118794x218794
6 x118794x203473
五 实验结
非线性方程求解析解时候困难采数值方法容易求似解次实验采迭代法求非线性方程根非线性方程选迭代形式收敛程度样造成效率精确度差
实验七 矩阵特征值问题计算
问题提出
利冪法反冪法求方阵模模特征值应特征量
设矩阵A特征分布：
试求列矩阵
（1） 求
取
结果
（2） 求
取
结果：
（3） 求
取
结果
(4)
取
收敛慢例子迭代次达
结果
（5）
似特征值试幂法求应特征量改进特征值（原点移法）
取
结果
二求
1掌握冪法反冪法求矩阵部分特征值算法程序设计
2会原点移法改进算法加速收敛矩阵BAPI取P值试求效果
3试取初始量观察结果影响
4矩阵特征值分布计算
三目意义
1求矩阵部分特征值问题具重实际意义求矩阵谱半径稳定性问题求矩阵模特征值
2进步掌握冪法反冪法原点移加速法程序设计技巧
3问题中题（5）反应利原点移反冪法求矩阵特征值特征量
四实验步骤
第步：写出实验需幂法求特征值反幂法求特征值程序
幂法程序：
function [muindex]TZ(Aepitmax)
if nargin<3
itmax100
end
if nargin<2 ep1e5
end
nlength(A)
uones(n1)
index0k0m10
while k vA*u[vmaxi]max(abs(v))
mv(i)uvm
if abs(mm1) index1break
end
m1mkk+1
end
反幂法程序：
function [muindex]FTZ(Aepitmax)
if nargin<3
itmax100
end
if nargin<2 ep1e5
end
nlength(A)
uones(n1)
index0k0m10
invAinv(A)
while k vinvA*u[vmaxi]max(abs(v))
mv(i)uvm
if abs(mm1) index1break
end
m1mkk+1
end
第二步：矩阵代入程序出结果
λ164210X1（00462 03749 10000）Tλ334723x3（10000 05229 02422）T
λ1213053X1（08724 05401 09973 05644 04972 10000）Tλ616214
五实验结
求n阶方阵A特征值特征量实际中常常碰问题通次实验掌握幂法反幂法求方阵特征值特征量绝值特征值特征量
实验八 常微分方程初值问题数值解法
问题提出
科学计算中常遇微分方程（组）初值问题需利Euler法改进Euler法RungKutta方法求数值解诸问题：
（1）
分取h010204时数值解 初值问题精确解
（2）
r3Adams显式预 校式求解
取步长h01四阶标准RK方法求值
（3）
改进Euler法四阶标准RK方法求解
取步长001计算数值解参考结果
（4）利四阶标准R K方法求二阶方程初值问题数值解
（I）
(II)
(III)
(IV)
二求
1 根初值问题数值算法分选择二初值问题编程计算
2 试分取步长考察某节点处数值解误差变化情况
3 试算法求解某初值问题结果异常
4 分析算法优缺点
三目意义
1 熟悉种初值问题算法编出算法程序
2 明确种算法精度选步长密切关系
3 通计算更加解种算法优越性
四实验步骤
function [xy]euler(funx0xfinaly0n)
if nargin<5n50
end
h(xfinalx0)n
x(1)x0y(1)y0
for i1n
x(i+1)x(i)+h
y(i+1)y(i)+h*feval(funx(i)y(i))
end
实验程序分析
（Ⅰ） （1）算法程序
function E Euler_1(funx0y0xNN)
Euler前公式中
fun阶微分方程函数
x0y0初始条件
xN取值范围端点
h区间步长
N区间数
xXn构成量
yyn构成量
xzeros(1N+1)yzeros(1N+1)
x(1)x0y(1)y0
h(xNx0)N
for n1N
x(n+1)x(n)+h
y(n+1)y(n)+h*feval(funx(n)y(n))
end
T[x＇y＇]
function zf(xy)
z4*xyx*y
(2)运行程序
>> Euler_1(＇f＇03220)
结果 ：
>> Euler_1(＇f＇03220)
T
0 30000
01000 29836
02000 29517
03000 29058
04000 28481
05000 27810
06000 27073
07000 26297
08000 25511
09000 24739
10000 24004
11000 23325
12000 22714
13000 22177
14000 21717
15000 21332
16000 21017
17000 20765
18000 20567
19000 20414
20000 20299
五实验结
科学技术工程问题常微分方程形式建立数学模型微分方程求解意义绝数微分方程问题难者根解析解研究微分方程数值方法非常意义
文档香网(httpswwwxiangdangnet)户传
《香当网》用户分享的内容，不代表《香当网》观点或立场，请自行判断内容的真实性和可靠性！
该内容是文档的文本内容，更好的格式请下载文档
数值分析课程实验报告
姓 名：
学 号：
学 院： 机 电 学 院
日 期： 2015 年 X 月X 日
目 录
实验 函数插值方法 1
实验二 函数逼曲线拟合 5
实验三 数值积分数值微分 7
实验四 线方程组直接解法 9
实验五 解线性方程组迭代法 15
实验六 非线性方程求根 19
实验七 矩阵特征值问题计算 21
实验八 常微分方程初值问题数值解法 24
实验 函数插值方法
问题提出
定元函数n+1节点值试Lagrange公式求插值项式分段二次Lagrange插值项式
数：
（1）
04
055
065
080
095
105
041075
057815
069675
090
100
125382
求五次Lagrange项式分段三次插值项式计算
值（提示：结果 ）
（2）
1
2
3
4
5
6
7
0368
0135
0050
0018
0007
0002
0001
试构造Lagrange项式计算值（提示：结果 ）
二求
1 利Lagrange插值公式
编写出插值项式程序
2 出插值项式分段三次插值项式表达式
3 根节点选取原问题（2）三点插值二点插值结果
4 插值问题Newton插值项式结果Newton插值项式：
中：
三目意义
1 学会常插值方法求函数似表达式解决实际问题
2 明确插值项式分段插值项式优缺点
3 熟悉插值方法程序编制
4 果绘出插值函数曲线观察光滑性
四实验步骤
（1）
04
055
065
080
095
105
041075
057815
069675
090
100
125382
求五次Lagrange项式分段三次插值项式计算
值（提示：结果 ）
第步：先matlab中定义lagranM文件拉格朗日函数代码：
function[cl]lagran(xy)
wlength(x)
nw1
lzeros(ww)
for k1n+1
v1
for j1n+1
if(k~j)
vconv(vpoly(x(j)))(x(k)x(j))
end
end
l(k)v
end
cy*l
end
第二步：然matlab命令窗口输入：
>>>> x[04 055 065 080095 105]y[041075 057815 069675 090 100 125382]
>> lagran(xy)
回车：
ans 1216264 4227503 5725667 3772549 1219718 150845
出求拉格朗日项式
p（x）1216264x54227503x4+5725667x33772549x2+1219718x150845
第三步：编辑窗口输入命令：
>> x[04 055 065 080095 105]
>> y1216264*x^54227503*x^4+5725667*x^33772549*x^2+1219718*x150845
>> plot(xy)
命令执行图示图形然
>> x0596
>> y1216264*x^54227503*x^4+5725667*x^33772549*x^2+1219718*x15084
y 06262
f（0596）06262
理f（099）10547
（2）
1
2
3
4
5
6
7
0368
0135
0050
0018
0007
0002
0001
试构造Lagrange项式分段三次插值项式计算值（提示：结果 ）
实验步骤：
第步定义
function[cl]lagran(xy)
wlength(x)
nw1
lzeros(ww)
for k1n+1
v1
for j1n+1
if(k~j)
vconv(vpoly(x(j)))(x(k)x(j))
end
end
l(k)v
end
cy*l
end
定义完拉格朗日M文件
第二步：然matlab命令窗口输入：
>>>> x[1 2 3 4 5 6 7] y[0368 0135 0050 0018 0007 0002 0001]
>> lagran(xy)
回车：
ans 00001 00016 00186 01175 04419 09683 09950
出求拉格朗日项式
p（x）00001x600016x5+00186x401175x3+04419x209683x+09950
第三步：编辑窗口输入命令：
>> x[1 2 3 4 5 6 7]
>> y00001*x^600016*x^5+00186*x^401175*x^3+04419*x^209683*x+09950
>> plot(xy)
命令执行图示图形然
>> x18
>> y1216264*x^54227503*x^4+5725667*x^33772549*x^2+1219718*x15084
y 01650
f（0596）06262
理f（615）23644
五实验结
插值离散数基础补插连续函数条连续曲线通全部定离散数点离散函数逼重方法利通函数限点处取值状况估算出函数点处似值
实验二 函数逼曲线拟合
问题提出
机数中找出规律性出似表达式问题生产实践科学实验中量存通常利数二法求拟合曲线
某冶炼程中根统计数含碳量时间关系试求含碳量时间t拟合曲线
t(分)
0 5 10 15 20 25 30 35 40 45 50 55
0 127 216 286 344 387 415 437 451 458 402 464
二求
1二法进行曲线拟合
2似解析表达式
3印出拟合函数印出误差
4外选取似表达式尝试拟合效果较
5* 绘制出曲线拟合图
三目意义
1掌握曲线拟合二法
2二法解超定线代数方程组
3探索拟合函数选择拟合精度间关系
四实验步骤：
第步先写出线性二法M文件
function clspoly(xym)
nlength(x)
bzeros(1m+1)
fzeros(nm+1)
for k1m+1
f(k)x^(k1)
end
af＇*f
bf＇*y＇
ca＼b
cflipud(c)
第二步命令窗口输入：
>>lspoly([0510152025303540455055][0127216286344387415437451458402464]2)
回车：
ans
00024
02037
02305
求拟合曲线y00024x2+02037x+02305
编辑窗口输入命令：
>> x[0510152025303540455055]
>> y00024*x^2+02037*x+02305
>> plot(xy)
命令执行图
五实验结
分析复杂实验数时常采分段曲线拟合方法利方法段实现佳逼段边界满足连续性导性分段函数光滑算法出相应误差分析出该方法分段曲线拟合中应方法凸轮实验数动分段拟合
实验三 数值积分数值微分
问题提出
选复合梯形公式复合Simpson公式Romberg算法计算
（1）
（2）
（3）
（4）
二求
1 编制数值积分算法程序
2 分两种算法计算积分较结果
3 分取步长试较计算结果（n 10 20等）
4 定精度求ε试变步长算法确定佳步长
三目意义
1 深刻认识数值积分法意义
2 明确数值积分精度步长关系
3 根定积分计算方法考虑二重积分计算问题
四 实验步骤
第步：编写种积分程序
复合梯形程序：
function ITX(xy)
nlength(x)mlength(y)
if n~m
error(＇The lengths of X and Y must be equal＇)
return
end
h(x(n)x(1))(n1)
a[1 2*ones(1n2) 1]
Ih2*sum(a*y)
复合Simpson程序：
function s simpr1(fabn)
h(ba)(2*n)
s10
s20
for k110
xa+h*(2*k1)
s1s1+feval(fx)
end
for k1(101)
xa+h*2*k
s2s2+feval(fx)
end
sh*(feval(fa)+feval(fb)+4*s1+2*s2)3
end
Romberg程序：
function I Romber_yang(funabep)
if nargin<4
ep1e5end
m1 hba
Ih2*(feval(funa)+feval(funb))T(11)I
while 1
N2^(m1)hh2II2
for i1N
II+h*feval(funa+(2*i1)*h)
end
T(m+11)IM2*Nk1
while M>1
T(m+1k+1)(4^k*T(m+1k)T(mk))(4^k1)
MM2kk+1
end
if abs(T(kk)T(k1k1))
end
mm+1
end
IT(kk)
第二步：复合梯形公式复合Simpson公式Romberg公式积分进行计算
1 积分Ι0144sin2xdx梯形积分T049871011辛普森积分S049871111Romberg积分R049871111
2 积分Ι01sinXXdxf(0)1梯形积分T094607307辛普森积分S094607308Romberg积分R094607307
3 积分Ι01eX4+X2dx梯形积分T039081248辛普森积分S039081185Romberg积分R039081885
4 积分Ι01ln1+X1+X2dx梯形积分T027218912辛普森积分S027219844Romberg积分R027219827
五实验结
通实验学会复合梯形公式复合Simpson公式Romberg公式编程应掌握MATLAB提供计算积分种函数方法
实验四 线方程组直接解法
问题提出
出列类型线性方程组请适算法计算解
1 设线性方程组
2 设称正定阵系数阵线方程组
3 三角形线性方程组
二求
1 述三方程组分利Gauss序消法Gauss列元消法方根法改进方根法追赶法求解（选择）
2 应结构程序设计编出通程序
3 较计算结果分析数值解误差原
4 利相应模块输出系数矩阵三角分解式
三目意义
1通该课题实验体会模块化结构程序设计方法优点
2运学计算方法解决类线性方程组直接算法
3提高分析解决问题力做学致
4 通三角形线性方程组解法体会稀疏线性方程组解法特点
四实验步骤：
列元高斯消法matlabM文件程序
function [xdetindex]Gauss(Ab)
求线形方程组列元Gauss消法中
A方程组系数矩阵
b方程组右端项
x方程组解
det系数矩阵A行列式值
index指标变量index0表示计算失败index1表示计算成功
[nm]size(A) nblength(b)
方程组行列维数相等时停止计算输出出错信息
if n~m
error(＇The rows and columns of matrix A must be equal＇)
return
end
方程组右端项维数匹配时停止计算输出出错信息
if m~nb
error(＇The columns of A must be equal the length of b＇)
return
end
开始计算先赋初值
index1det1xzeros(n1)
for k1n1
选元
a_max0
for ikn
if abs(A(ik))>a_max
a_maxabs(A(ik))ri
end
end
if a_max<1e10
index0return
end
交换两行
if r>k
for jkn
zA(kj)A(kj)A(rj)A(rj)z
end
zb(k)b(k)b(r)b(r)zdetdet
end
消元程
for ik+1n
mA(ik)A(kk)
for jk+1n
A(ij)A(ij)m*A(kj)
end
b(i)b(i)m*b(k)
end
detdet*A(kk)
end
detdet*A(nn)
回代程
if abs(A(nn))<1e10
index0return
end
for kn11
for jk+1n
b(k)b(k)A(kj)*x(j)
end
x(k)b(k)A(kk)
end
然命令窗口输入
>> A[4 2 3 1 2 1 0 0 0 08 6 5 3 6 5 0 1 0 04 2 2 1 3 2 1 0 3 10 2 1 5 1 3 1 1 9 44 2 6 1 6 7 3 3 2 38 6 8 5 7 17 2 6 3 50 2 1 3 4 2 5 3 0 116 10 11 9 17 34 2 1 2 24 6 2 7 13 9 2 0 12 40 0 1 8 3 24 8 6 3 1]
>> b[5 12 3 2 3 46 13 38 19 21]
>> gauss(Ab)
ans
10000
10000
00000
10000
20000
00000
30000
10000
10000
20000
高斯约消法maltabM文件程序
function [xflag]Gau_Jor(Ab)
求线形方程组列元Gauss约法消法中
A方程组系数矩阵
b方程组右端项
x方程组解
[nm]size(A) nblength(b)
方程组行列维数相等时停止计算输出出错信息
if n~m
error(＇The rows and columns of matrix A must be equal＇)
return
end
方程组右端项维数匹配时停止计算输出出错信息
if m~nb
error(＇The columns of A must be equal the length of b＇)
return
end
开始计算先赋初值
flag＇ok＇xzeros(n1)
for k1n
选元
max10
for ikn
if abs(A(ik))>max1
max1abs(A(ik))ri
end
end
if max1<1e10
falg＇failure＇return
end
交换两行
if r>k
for jkn
zA(kj)A(kj)A(rj)A(rj)z
end
zb(k)b(k)b(r)b(r)z
end
消元程
b(k)b(k)A(kk)
for jk+1n
A(kj)A(kj)A(kk)
end
for i1n
if i~k
for jk+1n
A(ij)A(ij)A(ik)*A(kj)
end
b(i)b(i)A(ik)*b(k)
end
end
end
输出x
for i1n
x(i)b(i)
end
然保存命令窗口输入：
>> A[4 2 4 0 2 4 0 02 2 1 2 1 3 2 04 1 14 1 8 3 5 60 2 1 6 1 4 3 32 1 8 1 22 4 10 34 3 3 4 4 11 1 40 2 5 3 10 1 14 20 0 6 3 3 4 2 19]
>> b[0 6 20 23 9 22 15 45]
>> Gau_Jor(Ab)
ans
1211481
1401127
297515
601528
109120
267963
54259
20185
五实验结
LU法调matlab中函数lu中L三角直接计算结果计算进行行变换果进行行变b变样会麻烦
实验五 解线性方程组迭代法
问题提出
实验四列目意义线性方程组试分选Jacobi 迭代法GaussSeidel迭代法SOR方法计算解
二求
1体会迭代法求解线性方程组消法做较
2分精度求迭代次数体会该迭代法收敛快慢
3方程组23SOR方法时选取松弛子ω080911112等试算法收敛性影响找出选松弛子佳者
4出种算法设计程序计算结果
三目意义
1通机计算体会迭代法求解线性方程组特点消法较
2运学迭代法算法解决类线性方程组编出算法程序
3体会机计算时终止步骤k >（予迭代次数）迭代法敛散性意义
4 体会初始解松弛子选取计算结果影响
四实验步骤
第步编写实验需Jacobi迭代法GaussSeidel迭代法SOR迭代法程序
Jacobi迭代法：
function [xkindex]J(Abepitmax)
if nargin<4
itmax100
end
if nargin<3 ep1e5
end
nlength(A)
k0
xzeros(n1)yzeros(n1)index1
while 1
for i1n
y(i)b(i)
for j1n
if j~i
y(i)y(i)A(ij)*x(j)
end
end
if abs(A(ii))<1e10|kitmax
index0return
end
y(i)y(i)A(ii)
end
if norm(yxinf)
end
xykk+1
end
GaussSeidel迭代法：
function [xkindex]G(Abepitmax)
if nargin<4
itmax100
end
if nargin<3 ep1e5
end
nlength(A)
k0
xzeros(n1)yzeros(n1)index1
while 1
yx
for i1n
zb(i)
for j1n
if j~i
zzA(ij)*x(j)
end
end
if abs(A(ii))<1e10|kitmax
index0return
end
zzA(ii)x(i)z
end
if norm(yxinf)
end
kk+1
end
SOR迭代法：
function [xkindex]SOR(Abepwitmax)
if nargin<5
itmax100
end
if nargin<4 w1end
if nargin<3 ep1e5
end
nlength(A)
k0
xzeros(n1)yzeros(n1)index1
while 1
yx
for i1n
zb(i)
for j1n
if j~i
zzA(ij)*x(j)
end
end
if abs(A(ii))<1e10|kitmax
index0return
end
zzA(ii)x(i)(1w)*x(i)+w*z
end
if norm(yxinf)
end
kk+1
end
第二步：实验出方程代入程序计算结果
1 设线性方程组
2 设称正定阵系数阵线方程组
3三角形线性方程组
五实验结
迭代法解线性方程组重实方法特适求解实际中量出现系数矩阵稀疏阵型线性方程组通次实验学会Jacobi迭代法GaussSeidel迭代法SOR迭代法程序编写掌握优缺点适条件
实验六 非线性方程求根
问题提出
设方程三实根
现采面六种计算格式求 f(x)0根
1
2
3
4
5
6
二求
1编制程序进行运算印出种迭代格式敛散情况
2事误差估计控制迭代次数印出迭代次数
3初始值选取迭代收敛影响
4分析迭代收敛发散原
三目意义
1通实验进步解方程求根算法
2认识选择计算格式重性
3掌握迭代算法精度控制
4明确迭代收敛性初值选取关系
四实验步骤
第步：编写实验需程序
function [x_starindexit]DD(funxepitmax)
if nargin<4
itmax100
end
if nargin<3
ep1e5
end
index0k1
while k
if abs(xx1)
end
kk+1
end
x_starxitk
第二步：分述六种形式表达式计算方程根结果
1 x10x20
2 x1穷x203473
3 x118794x218794
4 x103473x203473
5 x118794x218794
6 x118794x203473
五 实验结
非线性方程求解析解时候困难采数值方法容易求似解次实验采迭代法求非线性方程根非线性方程选迭代形式收敛程度样造成效率精确度差
实验七 矩阵特征值问题计算
问题提出
利冪法反冪法求方阵模模特征值应特征量
设矩阵A特征分布：
试求列矩阵
（1） 求
取
结果
（2） 求
取
结果：
（3） 求
取
结果
(4)
取
收敛慢例子迭代次达
结果
（5）
似特征值试幂法求应特征量改进特征值（原点移法）
取
结果
二求
1掌握冪法反冪法求矩阵部分特征值算法程序设计
2会原点移法改进算法加速收敛矩阵BAPI取P值试求效果
3试取初始量观察结果影响
4矩阵特征值分布计算
三目意义
1求矩阵部分特征值问题具重实际意义求矩阵谱半径稳定性问题求矩阵模特征值
2进步掌握冪法反冪法原点移加速法程序设计技巧
3问题中题（5）反应利原点移反冪法求矩阵特征值特征量
四实验步骤
第步：写出实验需幂法求特征值反幂法求特征值程序
幂法程序：
function [muindex]TZ(Aepitmax)
if nargin<3
itmax100
end
if nargin<2 ep1e5
end
nlength(A)
uones(n1)
index0k0m10
while k
mv(i)uvm
if abs(mm1)
end
m1mkk+1
end
反幂法程序：
function [muindex]FTZ(Aepitmax)
if nargin<3
itmax100
end
if nargin<2 ep1e5
end
nlength(A)
uones(n1)
index0k0m10
invAinv(A)
while k
mv(i)uvm
if abs(mm1)
end
m1mkk+1
end
第二步：矩阵代入程序出结果
λ164210X1（00462 03749 10000）Tλ334723x3（10000 05229 02422）T
λ1213053X1（08724 05401 09973 05644 04972 10000）Tλ616214
五实验结
求n阶方阵A特征值特征量实际中常常碰问题通次实验掌握幂法反幂法求方阵特征值特征量绝值特征值特征量
实验八 常微分方程初值问题数值解法
问题提出
科学计算中常遇微分方程（组）初值问题需利Euler法改进Euler法RungKutta方法求数值解诸问题：
（1）
分取h010204时数值解 初值问题精确解
（2）
r3Adams显式预 校式求解
取步长h01四阶标准RK方法求值
（3）
改进Euler法四阶标准RK方法求解
取步长001计算数值解参考结果
（4）利四阶标准R K方法求二阶方程初值问题数值解
（I）
(II)
(III)
(IV)
二求
1 根初值问题数值算法分选择二初值问题编程计算
2 试分取步长考察某节点处数值解误差变化情况
3 试算法求解某初值问题结果异常
4 分析算法优缺点
三目意义
1 熟悉种初值问题算法编出算法程序
2 明确种算法精度选步长密切关系
3 通计算更加解种算法优越性
四实验步骤
function [xy]euler(funx0xfinaly0n)
if nargin<5n50
end
h(xfinalx0)n
x(1)x0y(1)y0
for i1n
x(i+1)x(i)+h
y(i+1)y(i)+h*feval(funx(i)y(i))
end
实验程序分析
（Ⅰ） （1）算法程序
function E Euler_1(funx0y0xNN)
Euler前公式中
fun阶微分方程函数
x0y0初始条件
xN取值范围端点
h区间步长
N区间数
xXn构成量
yyn构成量
xzeros(1N+1)yzeros(1N+1)
x(1)x0y(1)y0
h(xNx0)N
for n1N
x(n+1)x(n)+h
y(n+1)y(n)+h*feval(funx(n)y(n))
end
T[x＇y＇]
function zf(xy)
z4*xyx*y
(2)运行程序
>> Euler_1(＇f＇03220)
结果 ：
>> Euler_1(＇f＇03220)
T
0 30000
01000 29836
02000 29517
03000 29058
04000 28481
05000 27810
06000 27073
07000 26297
08000 25511
09000 24739
10000 24004
11000 23325
12000 22714
13000 22177
14000 21717
15000 21332
16000 21017
17000 20765
18000 20567
19000 20414
20000 20299
五实验结
科学技术工程问题常微分方程形式建立数学模型微分方程求解意义绝数微分方程问题难者根解析解研究微分方程数值方法非常意义
文档香网(httpswwwxiangdangnet)户传
《香当网》用户分享的内容，不代表《香当网》观点或立场，请自行判断内容的真实性和可靠性！
该内容是文档的文本内容，更好的格式请下载文档