西南交通大学数值分析报告文档格式.docx
《西南交通大学数值分析报告文档格式.docx》由会员分享,可在线阅读,更多相关《西南交通大学数值分析报告文档格式.docx(20页珍藏版)》请在冰豆网上搜索。
比如一定用插值吗?
(5)逼近某个函数不用插值方式,有何变通之举?
(6)函数之间的误差如何度量,逼近的标准又是什么?
(7)如何比较好的使用插值多项式呢?
1.2计算思路
本题我选用拉格朗日插值函数,拉格朗日插值函数的构建算法是从101组数据中等距选取n组数据来构造n阶的拉格朗日插值函数,然后在[-5,5]区间选取m个插值点带入插值多项式计算出m组(x,y),本题中,我分别取了n=6,8,10,20四种不同的阶数,然后利用matlab绘图命令查看插值函数图像与原函数图像的逼近效果。
1.3计算结果
当阶数为5阶时,插值图像与原函数图象如下所示(红色曲线为插值函数曲线,蓝色为原函数曲线)
表1(迭代数据)
i
Xi
Yi
1
-5.0000
25.0000
26
-2.5000
6.2693
51
0
10.0000
71
1.0000
101
5.0000
由图象可以看出,除了在插值点附近的拟合效果还可以,其他的点都差距比较大,总的来说效果不是很好。
当阶数为7阶时图象如下:
表2(迭代数据)
17
-3.4000
11.5600
33
-1.8000
2.8832
49
-0.2000
5.2312
65
1.4000
3.0219
81
3.0000
8.9991
97
4.6000
21.1600
由此图象可以看出,利用7阶插值来进行计算,在0值附近的拟合不是很好,但总得来说比5阶函数较好。
当阶数为10阶时,图象如下:
表3(迭代数据)
12
-3.9000
15.2100
23
-2.8000
7.8405
34
-1.7000
2.5554
...
…
89
3.8000
14.4400
100
4.9000
24.0100
利用10阶函数进行插值计算时,在函数的两端会出现相比于原函数偏离比较大的点,这就是runge(龙格)现象,但不是很明显。
当阶数为20阶时,图象如下:
表4(迭代数据)
6
-4.5000
20.2500
11
-4.0000
16.0000
16
-3.5000
12.2500
91
4.0000
96
4.5000
由此图像可以看出,runge(龙格)现象已经非常明显了,由于在两端的龙格现象太明显,所以此图是在去除了两端偏离比较大的点。
1.4总结
利用拉格朗日插值函数构造插值多项式是比较常用的插值计算的方法之一,其解决了不需要求解方程组的问题,它也有一定的缺点:
(1)从实验中可以看出,并不是拉格朗日的插值阶数越高,其插值效果越好,虽然在一些局部拟合的效果较好,但是当阶数超过7阶时就容易出现runge现象,由上图中的阶数n=20时就可以看出;
(2)但是当插值阶数太少时,插值的效果又不是很好,由上图中阶数n=5时可以看出,插值效果非常不好;
(3)由上面的图象可以分析,利用拉格朗日进行插值计算时,所选取的插值点并不是在区间里平均选取最好,根据不同的函数曲线,应当用不同的方法选点。
实验二
2.1题目
松弛因子对SOR法收敛速度的影响。
用SOR法求解方程组Ax=b,其中
要求程序中不存系数矩阵A,分别对不同的阶数取w=1.1,1.2,...,1.9进行迭代,记录近似解x(k)达到||x(k)-x(k-1)||<
10-6时所用的迭代次数k,观察松弛因子对收敛速度的影响,并观察当w≤0或w≥2会有什么影响?
2.2松弛思想分析
SOR法(松弛法),是由高斯-塞德尔迭代方法演变来的。
高斯-塞德尔迭代方法的迭代格式为
,而松弛法的迭代格式为
,通过大量的实验发现,ω的取值很关键,如果取了好的ω,迭代方程的收敛速度会加快,根据理论基础可以知道,当0<
ω<
2时松弛法才能够收敛。
2.3问题的求解
分别取不同阶数系数矩阵和不同松弛因子进行计算,结果如下图所示:
当取阶数为5时,ω=1.1时:
迭代次数:
迭代过程:
x=0.82500.77690.76360.76001.0340
x=0.95610.94530.94261.01761.0014
x=0.98930.98681.00691.00051.0000
x=0.99741.00251.00011.00001.0000
x=1.00101.00011.00001.00001.0000
x=0.99991.00001.00001.00001.0000
x=1.00001.00001.00001.00001.0000
当阶数n=6,ω=1.3时:
x=0.97500.96690.96420.96340.96311.2880
x=0.99670.99730.99790.99831.10410.9474
x=1.00011.00021.00011.03440.96291.0037
x=1.00001.00001.01110.98121.00621.0009
x=1.00001.00360.99171.00501.00000.9997
x=1.00120.99661.00300.99950.99971.0000
x=0.99851.00150.99940.99991.00001.0000
x=1.00090.99971.00001.00011.00001.0000
x=0.99961.00001.00001.00001.00001.0000
x=1.00011.00001.00001.00001.00001.0000
x=1.00001.00001.00001.00001.00001.0000
当阶数n=5时,ω=-1时:
请输入阶数n:
5
请输入w的值:
-1
num=1,error=2.53e+000
num=2,error=3.93e+000
num=3,error=7.35e+000
num=4,error=1.80e+001
num=5,error=4.84e+001
num=6,error=1.32e+002
num=7,error=3.63e+002
num=10;
error=1.00e+003
num=11,error=2.77e+003
num=12,error=7.68e+003
num=13,error=2.12e+004
num=14,error=5.87e+004
num=15,error=1.62e+005
num=16,error=4.47e+005
num=17,error=1.23e+006
num=20,error=3.38e+006
num=21,error=9.27e+006
num=22,error=2.54e+007
num=23,error=6.96e+007
num=24,error=1.90e+008
num=25,error=5.20e+008
num=26,error=1.42e+009
num=27,error=3.87e+009
num=30,error=1.06e+010
num=31,error=2.88e+010
num=32,error=7.84e+010
num=33,error=2.13e+011
num=34,error=5.81e+011
num=35,error=1.58e+012
num=36,error=4.29e+012
x=
1.0e+012*
-0.30190.8370-1.53842.1285-1.9856
由数据中可以分析得出,当ω=-1时,迭代是不收敛的。
当阶数n=5时,ω=2.5时
2.5
num=1,error=1.38e+001
num=2,error=2.11e+001
num=3,error=2.58e+001
num=4,error=2.80e+001
num=5,error=4.63e+001
num=6,error=6.61e+001
num=7,error=8.74e+001
num=10,error=1.38e+002
num=11,error=1.72e+002
num=12,error=2.63e+002
num=13,error=4.26e+002
num=14,error=8.10e+002
num=15,error=1.51e+003
num=16,error=2.52e+003
num=17,error=3.58e+003
num=20,error=4.46e+003
num=21,error=6.34e+003
num=22,error=9.29e+003
num=23,error=1.25e+004
num=24,error=1.83e+004
num=25,error=2.82e+004
num=26,error=3.74e+004
num=27,error=6.20e+004
num=30,error=1.18e+005
num=31,error=1.96e+005
num=32,error=2.94e+005
num=33,error=4.09e+005
num=34,error=5.63e+005
num=35,error=7.87e+005
num=36,error=1.24e+006
1.0e+005*
0.46131.5993-1.47881.3274-2.6027
由数据中可以分析得出,当ω=2.5时,迭代也是不收敛的。
2.4总结
(1)在本题中,迭代次数和矩阵的未知数有以及ω的取值有一定的联系,不一定是未知数多迭代次数就多,也有可能在取相同ω时,未知数多迭代次数反而较少,这和SOR法的迭代原理有关。
(2)当ω≤0或ω≥2.0,SOR法不收敛,这是必要条件。
实验三
3.1题目
用Runge-Kutta4阶算法对初值问题y/=-20*y,y(0)=1按不同步长求解,用于观察稳定区间的作用,推荐两种步长h=0.1,0.2。
注:
此方程的精确解为:
y=e-20x
3.2Runge-Kutta法的基本思想
Runge-Kutta法的基本思想是通过f(x,y)某些点函数值的适当线性组合替换Euler法中的f(xk,yk),可能使得方法的精确度更高。
由于有了这种思想,因此标准的四阶Runge-Kutta的算法如下(h为求解步长):
K1=hf(xk,yk)
如果yk某种方法第k步的近似值,y(xk)是其准确值,其绝对误差为δk,即有δk=y(xk)-yk。
假定第k步之后的计算不再有舍入误差,只是由δk引起的扰动δm,都有︳δm︳<
︳δk︳,则称此方法是绝对稳定的。
数值解法的稳定性一般与步长h以及f(x,y)有关。
3.3问题的求解
给定初值x为0,我选定迭代终点为1,,分别选取h=0.1、h=0.2进行计算,并对迭代过程和迭代图像进行分析。
计算结果如下:
当h=0.1时
请输入a的值:
请输入b的值:
请输入步长h:
0.1
x[1]=0.000000y[1]=1.000000yi[1]=1.000000error[1]=0.000000
x[2]=0.100000y[2]=0.333333yi[2]=0.135335error[2]=0.197998
x[3]=0.200000y[3]=0.111111yi[3]=0.018316error[3]=0.092795
x[4]=0.300000y[4]=0.037037yi[4]=0.002479error[4]=0.034558
x[5]=0.400000y[5]=0.012346yi[5]=0.000335error[5]=0.012010
x[6]=0.500000y[6]=0.004115yi[6]=0.000045error[6]=0.004070
x[7]=0.600000y[7]=0.001372yi[7]=0.000006error[7]=0.001366
x[8]=0.700000y[8]=0.000457yi[8]=0.000001error[8]=0.000456
x[9]=0.800000y[9]=0.000152yi[9]=0.000000error[9]=0.000152
当h=0.2时
由图象和数据中可以看出,当步长为0.2时,迭代是不收敛的。
3.4问题的总结
(1)用标准的四阶Runge-Kutta算法计算常微分方程时,不同的步长算法的稳定性有影响。
不够稳定的步长下面的计算,误差会越来越大,结果失真严重。
(2)一般情况下,步长越小,标准的四阶Runge-Kutta算法的稳定性越高,精度也越高。
总结
通过此次数值分析上机实习,认识到了计算机编程在数值分析课程中的重要性。
也让我意识到了学习编程语言的重要性。
要想要将数学应用于实际工程中,特别是对于工科的学生来讲,这是一门非常重要的实践课程。
在此次报告中,首次接触了Matlab这门软件,之所以选这门,是因为很多工程中对数据的处理需要使用,这对我本身的专业是非常重要的。
本学期学习的数值分析方法,让我在对数学的理解上又有了新的认识。
在一连串看似杂乱无章的数据中,用数值分析进行处理,有了很多的收获。
在学习的过程中,意识到了这门课对我专业的重要性,让我对自己的专业有了更好的发挥。
对我写论文以及解决实际工程问题有很大的帮助。
最后,感谢老师给我这样的机会去接触这门语言,虽然只了解了皮毛,可是仍然收获颇多。
由于初次接触这门软件,在报告中仍难免会有不完善甚至错误的地方,望谅解!
附录
实验一程序
functionyy=L(xi)
f=fopen('
xk.txt'
'
r'
);
x=fscanf(f,'
%f'
g=fopen('
yk.txt'
y=fscanf(g,'
n=input('
请输入插值函数阶数:
'
s=0;
xi=[-3.5:
0.1:
3.5];
fori=1:
fix(101/(n-1)):
t=ones(1,length(xi));
forj=1:
ifj~=i
t=t.*(xi-x(j))/(x(i)-x(j));
end
s=s+t*y(i);
yy=s;
plot(x,y,'
b'
xi,yy,'
)
实验二程序
n=input('
w=input('
fori=1:
n
ifi==j
a(i,j)=-4;
elseif(i==j+1|i==j-1)
a(i,j)=1;
else
a(i,j)=0;
end
if(i==1|i==n)
b(i)=-3;
b(i)=-2;
x=zeros(1,n);
fornum=1:
30
error=0;
xb=x(i);
ifi~=j,s=s+a(i,j)*x(j);
x(i)=w*(b(i)-s)/a(i,i)+(1-w)*x(i);
error=error+abs(x(i)-xb);
fprintf('
num=%3.o,error=%7.2e\n'
num,error)
iferror/3<
0.000001,break;
x
实验三程序
a=input('
b=input('
h=input('
x=[a:
h:
b];
y
(1)=1;
yi
(1)=exp(-20*x
(1));
error
(1)=0;
n=(b-a)/h+1;
fori=2:
k1=-20*y(i-1)*h;
k2=-20*(y(i-1)+k1/2)*h;
k3=-20*(y(i-1)+k2/2)*h;
k4=-20*(y(i-1)+k3)*h;
y(i)=y(i-1)+(k1+2*k2+2*k3+k4)/6;
yi(i)=exp(-20*x(i));
error(i)=y(i)-yi(i);
x[%2.0f]=%10.6fy[%2.0f]=%10.6fyi[%2.0f]=%10.6ferror[%2.0f]=%10.6f\n'
i-1,x(i-1),i-1,y(i-1),i-1,yi(i-1),i-1,error(i-1))
plot(x,y,'
g'
x,yi,'