常微分方程数值解法在动力天文中的应用.docx
- 文档编号:2181738
- 上传时间:2022-10-27
- 格式:DOCX
- 页数:20
- 大小:226.96KB
常微分方程数值解法在动力天文中的应用.docx
《常微分方程数值解法在动力天文中的应用.docx》由会员分享,可在线阅读,更多相关《常微分方程数值解法在动力天文中的应用.docx(20页珍藏版)》请在冰豆网上搜索。
常微分方程数值解法在动力天文中的应用
青岛科技大学
学士学位论文
(2010年〜2014年)
题目:
常微分方程数值解法在动力天文中的应用
学院:
专业班级:
学生姓名:
数理学院
信计101
潘力杨
学号:
1006010129
指导教师:
王斌
起讫日期:
2014年3月
—2014年6月
常微分方程数值解法在动力天文中的应用
、尸■、亠
前言常微分方程,作为数学分支的新发展在人类认识自然,改造自然方面发挥了巨大的力量。
牛顿在研究天体力学和机械动力学时,一微分方程为工具,从而在理论上得出了行星运动的规律。
后来,法国天文学家勒维烈和英国天文学家亚当斯都各自分别利用了常微分方程,演算出了当时为止的海王星位置及其运动规律。
如今,常微分方程更是在许多学科领域内发挥着重要作用。
小行星、人造天体的运动及其规律,导弹轨道的计算,飞机空气动力学的研究,化学反应过程问题等,这些问题都可以转化为求常微分方程的解,或研究解的性质的问题。
虽然常微分方程存在众多的研究和应用领域,但事实上,只有极少数诸如线性系数常微分方程和个别典型微分方程的解能用初等函数,特殊函数或它们的级数与积分表达。
而在变系数常微分方程的解析求解上已经比较困难了,更别说一般的非线性常微分方程了。
大多数情况下,常微分方程只能用近似方法求解,近似方法可大致分为两大类:
用级数解法,逐次逼近法等求展开式的近似解析法;和在方程一些离散点上给出近似解得数值解法。
在动力天文的具体应用中,事实上由于所处的角度关系,我们很难宏观直接地去观测天体的具体运动状态,并直接以经典物理的力学方法去计算它的准确运动轨道及规律。
由于涉及到太多的不确定因素以及运动状态的复杂性,我们只能尽量地以我们的实际观测值为参数,建立模型,运用数学的方法去研究,拟合它的运动轨迹,以此为基础才能进一步作出定性的分析。
在天体运动方程的求解问题上,基本上来说也只有二体问题能给出严格解。
因此为了得出满足具体观测精度的定量解,我们不得不求助于数值解法。
利用数值解法求解常微分方程初值问题,一般只需要在满足规定精度的条件下求得解在若干相应点上的近似解,或是求出解的近似表达式就行,而无需使用复杂的方法试图去拟合解的表达式,或花费大量的计算去得出形式解。
由于它的理论简单,可操作性强且能适应各种复杂的问题条件,因此在数学建模领域有着独一无二的重要地位。
在过去技术尚未得到发展的时代,在运用数值解法求解天体运动方程的过程中,处于计算精度,方程稳定性和收敛性的需要,我们往往得做出大量繁杂重复的运算从而必须花费大量的时间,可靠性也得不到保证,这也反过来制约了数值解法的应用。
但现如今,随着计算机的普及以及编程运算技术的成熟运用,关于数值解法的应用领域也出现了新的契机,设定步长,循环次数的程式化运算大大缩短了我们的运算时间同时也提高了运算的精度,例如只要通过Matlab就可以完成简单方程求解问题的建模、迭代、作图、误差分析等一系列的工作,大大减轻了我们的负担同时也提高了运算的效率。
使得我们可以更从容地将数值解法运用到天体运动求解问题的更深邻域中去。
第一章常微分方程数值解法介绍
常微分方程可分为线性、非线性、高阶方程与方程组等类,其中线性方程包含于非线性类中,高阶方
程可化为一阶方程组,若方程组中的所有未知量看作一个向量,则方程组可写成向量形式的单个方程。
因
此研究一阶常微分方程的初值问题
dy
f(x,y),^x
dx(9-1)
y(a)=y。
的数值解法具有典型性,其中方程的解为y:
R>R。
只有保证问题(9-1)的解存在惟一的前提下,其数值解法的研究才有意义。
由常微分方程的基本理论,我们有:
定理9.1如果(9-1)中的f(x,y)满足条件
(1)f(x,y)在区域D={(x,y)a兰x兰b,-旳vy£+旳}上连续;
(2)f(x,y)在D上关于y满足Lipschitz条件,即存在常数L,使得
f(x,y)—f(x,y)兰Ly—y|
则初值问题(9-1)在区间[a,b]上存在惟一的连续解y=y(x)。
所谓数值解法,是通过常微分方程离散化而给出解在某些节点上的近似值。
对于问题(9-1),在区间〔a,bl引入若干节点
a=x0:
为:
X2:
HI:
Xn=b
m二xn1-Xn,(n=0,1川l,N-1)称为由xn到Xn1的步长,当hn=h(常数)时称为定步长,否则称为变
b—a
步长。
多数情况下,采用等步长,即h,xn=a,nh(n=1,2,|I(,N)。
记问题(9-1)的精确解为y(x),
N
记y(Xn)的近似值为yn,记f(xn,yn)为fn,『ng=1,2/,N)称为问题(9-1)的数值解。
求初值问题数值解的方法一般采用步进法,即在计算出y,^n后计算yn,,方法分为单步法和多步法。
单步法是指在计算yn1时只利用yn,而多步法在计算yn1时不仅要利用y.,还要利用前面已经计算出来的若干个yn_j,j=1,2,川,I-1°我们称要用到Yn和yn_j,j=1,2,|11,I-1的多步方法为I步方法。
不论单步法和多步法,它们都有显式方法和隐式方法之分。
显式单步法的计算公式可以写为:
yn+1=ynh(xn,yn,h),此公式右端不含ynd。
隐式单步法的计算公式可以写为:
yn+1=yn•h'(Xn,yn,yn.1,h),此公式右端含有yn-1,从而是『n计的
方程式,需要求解方程或者采用迭代法。
显式多步法的计算公式可以写为yn+1=%•h(Xn,yn,yn-1,川,丫心・1,h)。
隐式多步法的计算公式可以写为yn+1=ynh(Xn,yn!
yn,yn-1,,yn41,h),k-1-1。
显式公式与隐式公式各有特点。
显式公式的优点是使用方便,计算简单,效率高,其缺点是计算精度
低,稳定性差。
隐式公式正好与它相反,它具有计算精度高,稳定性好等优点,但求解过程很复杂,一般采用迭代法。
为了结合各自的优点,通常将显式公式与隐式公式配合使用,由显式公式提供迭代初值,再经隐式公式迭代校正。
1.1欧拉方法(显)
欧拉方法的建立可以通过以下三种方法:
(一)用差商近似导数,用向前差商y(XnQ「y(Xn)代替y(xn),则得
(n=0,1,L,N-1)
(9-2)
这样,问题(9-1)的近似解可以表示为
Yn1=Ynhf(Xn,Yn)
Y0=Y(Xo)
按式此式由初值y0经过N步迭代,可逐次算出y1,y2/yN,此方程称为差分方程。
式(9-2)称为求解
一阶常微分方程初值问题(9-1)的欧拉公式,也称显式欧拉公式。
需要说明的是,用不同的差商近似导数,将得到不同的计算公式。
(2)用Taylor多项式近似,把y(Xn』在Xn点处Taylor展开,取一次多项式近似,则得:
h2„y(Xn1)=y(Xn)hy(Xn)万y()
h2
二y(Xn)hf(Xn,y(Xn))—y()[Xn,Xn1]
设步长h的值较小,略去余项,并以yn代替y(xn),便得:
%1二ynhf(Xn,yn)
(3)将问题(9-1)中的微分方程在区间[X^Xn」上两边积分,可得:
Xn+
y(Xn卅)—y(Xn)=Jf(x,y(x))dx(n=0,1,…,N-1)
Xn
用yn1,yn分别代替y(Xn1),y(xn),若对右端积分采用取左端点的矩形公式,即:
Xn1
Jf(x,y(x))d^hf(Xn,yn)
xn
同样可得到显式欧拉公式。
1.2欧拉方法(隐)
在微分方程离散化时,用向后差商代替导数,即y(xn1)也山心°,则得到如下差分方程
hyn1fn■hf(Xn1,yn1)(n=0,1,|11,N-1)
«(9-3)
y。
二Y(x。
)
此公式称为求解问题(9-1)所得的数值解称为隐式欧拉法。
隐式欧拉法与欧拉公式(9-2)式上相似,但实际计算时却复杂得多。
欧拉公式(9-2)计算yn★的公式中
不含yn1,但是隐式欧拉法计算yn1的公式中含有yn1。
在求解yn1时,yn,Xn-1为已知,『n1是方程
yn^ynhf(Xn!
yn1)的根。
一般说来,这是一个非线性方程,当不能精确求解得到yn1的表达式时,
我们需要运用简单迭代法来求解。
迭代格式为:
■'[0]
yn1ynhf(Xn,yn)ynv'ynhf(Xni,ynk]i)(k=0,1,2,川)
由于f(x,y)满足Lipschitz条件,所以
yn羊-yn+=hf(Xn^yn』—f(X.申,y.』EhLyn申—yn申
由此可知,只要0:
:
:
hL:
:
:
1,迭代法就收敛到解yn勺。
1.3梯形公式
利用数值积分方法将微分方程离散化时,若用梯形公式计算式中右端积分,即
冷ih
ff(x,y(x))dx止二[f(Xn,y(Xn))+f(Xn+,y(Xn^))]
Xn2
并用yn,yn1代替y(Xn),y(Xni),则得计算公式
这就是求解初值问题(9-1)的梯形公式。
梯形公式也是隐式格式,一般需用迭代法求解,迭代公式为
yn°fr=yn+hf(Xn,yn)扌h-
[yn.=yn+2」(Xn,yn)+f(Xn和ynk;)](k=0,1川)
由于函数f(x,y)关于y满足Lipschitz条件,所以
[k出]
h
yn卅_yn+
=—
2
[k]
f(XnH1,yn^)-f(Xn卄y^)兰牛|yn^-Yn卅
2hL
其中L为Lipschitz常数。
因此,当0§:
:
:
1时,迭代法就收敛到解yn计。
1・4改进的欧拉法
相比于显式欧拉法和隐式欧拉法,梯形方法提高了精度,但是梯形方法,在应用迭代公式进行实际计算时,每迭代一次,都要重新计算函数f(x,y)的值,而迭代又要反复进行若干次,计算量很大。
当函数
f(x,y)比较复杂时,这个问题会变得更加突出。
为了控制计算量,我们先用欧拉公式求得一个初步的近似值yn1,称之为预测值,预测值yn1的精度可能很差,再用梯形公式将它校正一次得yn计,称为校正值。
这样的预测校正系统通常称为改进的欧拉公式。
即
预测:
ynynhf(Xn,yn)
校正:
yn1二yn■h〔f(Xn,yn)■f(Xn1,Vn1)1
2
为了便于编制程序上机,将上式改写成
yp=yn+hf(Xn,yn)
*yq=yn+hf(Xn+h,『p)(9-5)
1/、yn#=2(yp+yq)
算法实现如表
算法9.5
1输入a,b,f(x,y),h,初值y0b-a2N,n=0,x二a,y二y0h
3计算
y^yhf(x,y)y^yhf(xh,yp)
1-(ypyq)二y,xh-x输出(x,y)
4若n:
:
:
N-1,置n•1=•n,转③;否则停机。
1.5龙格—库塔(Runge-Kutta)法
Runge-Kutta法由数学家卡尔龙格和马丁威尔海姆库塔于1900年左右发明。
设初值问题(9-1)的解y=y(x)・C1
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 微分方程 数值 解法 动力 天文 中的 应用