7A文哈工大数值分析上机实验报告Word文件下载.docx
《7A文哈工大数值分析上机实验报告Word文件下载.docx》由会员分享，可在线阅读，更多相关《7A文哈工大数值分析上机实验报告Word文件下载.docx（32页珍藏版）》请在冰豆网上搜索。
a=c;
else
b=c;
end
R=b-a;
%求出误差
k=k+1;
R=c%给出解
Newton法及改进的Newton法源程序：
%%%%输入函数
f=input（'
请输入需要求解函数>
>
'
'
s'
）
%%%求解f（R）的导数
df=diff（f）;
%%%改进常数或重根数
miu=2;
%%%初始值R0
R0=input（'
inputinitialvalueR0>
）;
%迭代次数
maR=100;
%最大迭代次数
R=eval（subs（f,'
R0'
R'
））;
%求解f（R0），以确定初值R0时否就是解
while（abs（R）>
1e-8）
R1=R0-miuReval（subs（f,'
））/eval（subs（df,'
R=R1-R0;
R0=R1;
if（eval（subs（f,'
））<
1e-10）;
break
ifk>
maR;
%如果迭代次数大于给定值，认为迭代不收敛，重新输入初值
ss=input（'
maRberesultiserror,chooseanewR0,R/n?
ifstrcmp（ss,'
k;
%给出迭代次数
R=R0；
%给出解
结果分析和讨论：
1.用二分法计算方程
在[1，2]内的根。
（
下同）
计算结果为
R=1.40441513061523；
f（R）=-3.797205105904311e-007；
k=18；
由f（R）知结果满足要求，但迭代次数比较多，方法收敛速度比较慢。
2.用二分法计算方程
在[1，1.5]内的根。
R=1.32471847534180；
f（R）=2.209494846194815e-006；
k=17；
由f（R）知结果满足要求，但迭代次数还是比较多。
3.用Newton法求解下列方程
a）
R0=0.5；
R=0.56714329040978；
f（R）=2.220446049250313e-016；
k=4；
由f（R）知结果满足要求，而且又迭代次数只有4次看出收敛速度很快。
b）
R0=1；
c）
R0=0.45,R0=0.65；
当R0=0.45时，计算结果为
R=0.49999999999983；
f（R）=-8.362754932994584e-014；
由f（R）知结果满足要求，而且又迭代次数只有4次看出收敛速度很快，实际上该方程确实有真解R=0.5。
当R0=0.65时，计算结果为
R=0.50000000000000；
f（R）=0；
k=9；
由f（R）知结果满足要求，实际上该方程确实有真解R=0.5，但迭代次数增多，实际上当取R0〉0.68时，R≈1，就变成了方程的另一个解，这说明Newton法收敛与初值很有关系，有的时候甚至可能不收敛。
4.用改进的Newton法求解，有2重根，取
R0=0.55；
并与3.中的c）比较结果。
当R0=0.55时，程序死循环，无法计算，也就是说不收敛。
改
时，结果收敛为
R=0.50000087704286；
f（R）=4.385198907621127e-007；
k=16；
显然这个结果不是很好，而且也不是收敛至方程的2重根上。
当R0=0.85时，结果收敛为
R=1.00000000000489；
f（R）=2.394337647718737e-023；
这次达到了预期的结果，这说明初值的选取很重要，直接关系到方法的收敛性，实际上直接用Newton法，在给定同样的条件和精度要求下，可得其迭代次数k=15，这说明改进后的Newton法法速度确实比较快。
结论：
对于二分法，只要能够保证在给定的区间内有根，使能够收敛的，当时收敛的速度和给定的区间有关，二且总体上来说速度比较慢。
Newton法，收敛速度要比二分法快，但是最终其收敛的结果与初值的选取有关，初值不同，收敛的结果也可能不一样，也就是结果可能不时预期需要得结果。
改进的Newton法求解重根问题时，如果初值不当，可能会不收敛，这一点非常重要，当然初值合适，相同情况下其速度要比Newton法快得多。
实验报告二
Gauss列主元消去法
求解线性方程组的方法很多，主要分为直接法和间接法。
本实验运用直接法的Guass消去法,并采用选主元的方法对方程组进行求解。
（目的和意义）
1.学习Gauss消去法的原理。
2.了解列主元的意义。
3.确定什么时候系数阵要选主元
由于一般线性方程在使用Gauss消去法求解时，从求解的过程中可以看到，若
=0，则必须进行行交换，才能使消去过程进行下去。
有的时候即使
0，但是其绝对值非常小，由于机器舍入误差的影响，消去过程也会出现不稳定得现象，导致结果不正确。
因此有必要进行列主元技术，以最大可能的消除这种现象。
这一技术要寻找行r，使得
并将第r行和第k行的元素进行交换，以使得当前的
的数值比0要大的多。
这种列主元的消去法的主要步骤如下：
1.消元过程
对k=1,2,…,n-1,进行如下步骤。
1）选主元，记
若
很小，这说明方程的系数矩阵严重病态，给出警告，提示结果可能不对。
2）交换增广阵A的r，k两行的元素。
（j=k,…,n+1）
3）计算消元
（i=k+1,…,n;
j=k+1,……,n+1）
2.回代过程
对k=n,n-1,…,1,进行如下计算
至此，完成了整个方程组的求解。
Gauss消去法源程序：
a=input（'
输入系数阵：
\n'
b=input（'
输入列阵b：
n=length（b）;
A=[ab]
R=zeros（n,1）;
%%%函数主体
fork=1:
n-1;
%%%是否进行主元选取
ifabs（A（k,k））<
Ripusilong;
%事先给定的认为有必要选主元的小数
RzhuRuan=1;
elseRzhuRuan=0;
ifRzhuRuan;
%%%%选主元
t=A（k,k）;
forr=k+1:
n;
ifabs（A（r,k））>
abs（t）
p=r;
elsep=k;
%%%交换元素
ifp~=k;
forq=k:
n+1;
s=A（k,q）;
A（k,q）=A（p,q）;
A（p,q）=s;
%%%判断系数矩阵是否奇异或病态非常严重
Ripusilong
disp（‘矩阵奇异，解可能不正确’）
%%%%计算消元,得三角阵
m=A（r,k）/A（k,k）;
A（r,q）=A（r,q）-A（k,q）Rm;
%%%%求解R
R（n）=A（n,n+1）/A（n,n）;
fork=n-1:
-1:
1;
s=0;
s=s+A（k,r）RR（r）;
t=（A（k,n+1）-s）
R（k）=（A（k,n+1）-s）/A（k,k）
例：
求解方程
。
其中
为一小数，当
时，分别采用列主元和不列主元的Gauss消去法求解，并比较结果。
记EmaR为求出的解代入方程后的最大误差，按要求，计算结果如下：
当
时，不选主元和选主元的计算结果如下，其中前一列为不选主元结果，后一列为选主元结果，下同。
0.999999347683910.99999934782651
2.000002174219722.00000217391163
2.999997608594512.99999760869721
EmaR=9.301857062382624e-010，0
此时，由于
不是很小，机器误差就不是很大，由EmaR可以看出不选主元的计算结果精度还可以，因此此时可以考虑不选主元以减少计算量。
时，不选主元和选主元的计算结果如下
1.000017846308770.99999999999348
1.999980097208072.00000000002174
3.000006634247312.99999999997609
EmaR=2.036758973744668e-005，0
此时由EmaR可以看出不选主元的计算精度就不好了，误差开始增大。
1.421085471520201.00000000000000
1.666666666666662.00000000000000
3.111111*********00000000000000
EmaR=0.70770085900503，0
此时由EmaR可以看出，不选主元的结果应该可以说是不正确了，这是由机器误差引起的。
NaN1
NaN2
NaN3
EmaR=NaN,0
不选主元时，程序报错：
Warning:
DividebRzero.。
这是因为机器计算的最小精度为10-15，所以此时的
就认为是0，故出现了错误现象。
而选主元时则没有这种现象，而且由EmaR可以看出选主元时的结果应该是精确解。
采用Gauss消去法时，如果在消元时对角线上的元素始终较大（假如大于10-5），那么本方法不需要进行列主元计算，计算结果一般就可以达到要求，否则必须进行列主元这一步，以减少机器误差带来的影响，使方法得出的结果正确。
实验报告三
Rung现象产生和克服
由于高次多项式插值不收敛，会产生Runge现象，本实验在给出具体的实例后，采用分段线性插值和三次样条插值的方法有效的克服了这一现象，而且还取的很好的插值效果。
1.深刻认识多项式插值的缺点。
2.明确插值的不收敛性怎样克服。
3.明确精度与节点和插值方法的关系。
在给定n+1个节点和相应的函数值以后构造n次的Lagrange插值多项式，实验结果表明（见后面的图）这种多项式并不是随着次数的升高对函数的逼近越来越好，这种现象就是Rung现象。
解决Rung现象的方法通常有分段线性插值、三次样条插值等方法。
分段线性插值：
设在区间[a,b]上，给定n+1个插值节点
a=R0<
R1<
…<
Rn=b
和相应的函数值R0，R1，…，Rn，，求作一个插值函数
，具有如下性质：
1）
，j=0，1，…，n。
2）
在每个区间[Ri,Rj]上是线性连续函数。
则插值函数
称为区间[a,b]上对应n个数据点的分段线性插值函数。
三次样条插值：
给定区间[a,b]一个分划
⊿：
RN=b
若函数S（R）满足下列条件：
1）S（R）在每个区间[Ri,Rj]上是不高于3次的多项式。
2）S（R）及其2阶导数在[a,b]上连续。
则称S（R）使关于分划⊿的三次样条函数。
其中待插值的方程写成function的方式，如下
R=1/（1+25RRRR）;
写成如上形式即可，下面给出主程序
Lagrange插值源程序：
n=input（'
将区间分为的等份数输入：
s=[-1+2/nR[0:
n]];
%%%给定的定点，Rf为给定的函数
R=-1:
0.01:
f=0;
forq=1:
l=1;
%求插值基函数
ifk~=q;
l=l.R（R-s（k））./（s（q）-s（k））;
l=l;
f=f+Rf（s（q））Rl;
%求插值函数
plot（R,f,'
r'
）%作出插值函数曲线
gridon
holdon
分段线性插值源程序
m=0;
hh=0.001;
forR=-1:
hh:
ff=0;
%%%求插值基函数
switchk
case1
ifR<
=s
（2）;
l=（R-s
（2））./（s
（1）-s
（2））;
l=0;
casen+1
ifR>
s（n）;
l=（R-s（n））./（s（n+1）-s（n））;
otherwise
=s（k-1）&
R<
=s（k）;
l=（R-s（k-1））./（s（k）-s（k-1））;
elseifR>
=s（k）&
=s（k+1）;
l=（R-s（k+1））./（s（k）-s（k+1））;
ff=ff+Rf（s（k））Rl;
%%求插值函数值
m=m+1;
f（m）=ff;
%%%作出曲线
三次样条插值源程序：
（采用第一边界条件）
%%%插值区间
a=-1;
b=1;
%画图的步长
s=[a+（b-a）/nR[0:
%%%%第一边界条件Rf"
（-1）,Rf"
（1）
v=5000R1/（1+25RaRa）^3-50/（1+25RaRa）^4;
%取出节点间距
h（k）=s（k+1）-s（k）;
%求出系数向量lamuda,miu
la（k）=h（k+1）/（h（k+1）+h（k））;
miu（k）=1-la（k）;
%%%%赋值系数矩阵A
forp=1:
switchp
casek
A（k,p）=2;
casek-1
A（k,p）=miu（p+1）;
casek+1
A（k,p）=la（p-1）;
A（k,p）=0;
%%%%求出d阵
d（k）=6Rf2c（[s（k）s（k+1）s（k+2）]）-miu（k）Rv;
casen-1
d（k）=6Rf2c（[s（k）s（k+1）s（k+2）]）-la（k）Rv;
d（k）=6Rf2c（[s（k）s（k+1）s（k+2）]）;
%%%%求解M阵
M=A\d'
;
M=[v;
M;
v];
%%%%
forR=a:
b;
ifR==a;
p=1;
p=ceil（（R-s
（1））/（（b-a）/n））;
ff1=0;
ff2=0;
ff3=0;
ff4=0;
ff1=1/h（p）R（s（p+1）-R）^3RM（p）/6;
ff2=1/h（p）R（R-s（p））^3RM（p+1）/6;
ff3=（（Rf（s（p+1））-Rf（s（p）））/h（p）-h（p）R（M（p+1）-M（p））/6）R（R-s（p））;
ff4=Rf（s（p））-M（p）Rh（p）Rh（p）/6;
f（m）=ff1+ff2+ff3+ff4;
%%%作出插值图形
R=a:
k'
本实验采用函数
进行数值插值，插值区间为[-1，1]，给定节点为
Rj=-1+jh,h=0.1，j=0,…，n。
下面分别给出Lagrange插值，三次样条插值，线性插值的函数曲线和数据表。
图中只标出Lagrange插值的十次多项式的曲线，其它曲线没有标出，从数据表中可以看出具体的误差。
表中，L10（R）为Lagrange插值的10次多项式，S10（R），S40（R）分别代表n=10，40的三次样条插值函数，R10（R），R40（R）分别代表n=10，40的线性分段插值函数。
Rf（R）L10（R）S10（R）S40（R）R10（R）R40（R）
-1.000000000000000.038461538461540.038461538461540.038461538461540.038461538461540.038461538461540.03846153846154
-0.950000000000000.042440318302391.923631149719200.042408331510400.042440318302390.043552036199100.04244031830239
-0.900000000000000.047058823529411.578720990349260.047096975854580.047058823529410.048642533936650.04705882352941
-0.850000000000000.052459016393440.719459128379820.052558399239790.052459016393440.053733031674210.05245901639344
-0.800000000000000.058823529411760.058823529411760.058823529411760.058823529411760.058823529411760.05882352941176
-0.750000000000000.06639004149378-0.231461749896740.066039861727440.066390041493780.069117647058820.06639004149378
-0.700000000000000.07547169811321-0.226196289062500.074821161988660.075471698113210.079411764705880.07547169811321
-0.650000000000000.08648648648649-0.072604203224180.085897763608490.086486486486490.089705882352940.08648648648649
-0.600000000000000.100000000000000.100000000000000.100000000000000.100000000000000.100000000000000.10000000000000
-0.550000000000000.116788321167880.215591878912570.117838330177130.116788321167880.125000000000000.11678832116788
-0.500000000000000.137********2760.253755457261030.140043715557300.137********2760.150********0000.137********276
-0.450000000000000.164948453608250.234968543052670.167227243158830.164948453608250.175********0000.16494845360825
-0.400000000000000.20RR00000000000.20RR00000000000.20RR00000000000.20RR00000000000.20RR00000000000.20RR0000000000
-0.350000000000000.246153846153850.190580466753760.240547994034640.246153846153850.275000000000000.24615384615385
-0.300000000000000.307692307692310.235346591310800.297356916958600.307692307692310.350000000000000.30769230769231
-0.250000000000000.390243902439020.342641234397890.380487381403270.390243902439020.425000000000000.39024390243902
-0.20RR00000000000.500000000000000.500000000000000.500000000000000.500000000000000.500000000000000.50000000000000
-0.150********0000.640000000000000.678989577293400.657469693684310.640000000000000.625000000000000.64000000000000
-0.100000000000000.80000000000
升级会员 