《偏微分方程数值解PPT课件》由会员分享，可在线阅读，更多相关《偏微分方程数值解PPT课件（146页珍藏版）》请在金锄头文库上搜索。
1、1,偏微分方程数值解(Numerical Solution of Partial Differential Equations),主讲：王曰朋 ,2,参考数目,George J. Haltiner, Roger Terry Williams, Numerical Prediction and Dynamic Meteorology(2nd Edition), the United States of America, 1979.,2. Curtis F.Gerald and Patrick O., Applied Numerical Analysis, Person Education, In
2、c., 2004.,3. Eugenia Kalnay, Atmospheric Modeling, Data Assimilation and Predictability, the press Syndicate of the University of Cambridge,2003.,4. Arieh Iserles, A First Course in the Numerical Analysis of Differential Equations, Cambridge University Press,1996.,5. 李荣华，冯国忱. 微分方程数值解. 北京：人民教育出版社，198
3、0.,6. 徐长发，李红. 实用偏微分方程数值解法. 华中科技大学出版社，2003.,7. 沈桐立，田永祥等. 数值天气预报. 北京：气象出版社，2007.,3,数值天气预报PDE数值解,挪威气象学家V.Bjerknes（1904）提出数值预报的思想：通过求解一组方程的初值问题可以预报将来某个时刻的天气思想； L.F.Richardson(1922)：开创了利用数值积分进行预报天气的先例，由于一些原因（如，计算稳定性问题“Courant,1928”）并没有取得预期的效果尝试； Charney, Fjortoft, and Von Neumann(1950), 借助于Princeton大学的的计
4、算机(ENIAC)，利用一个简单的正压涡度方程（C.G.Rossby,1940）对500mb的天气形式作了24小时预报-成功；,4,The Electronic Numerical Integrator and Computer (ENIAC).,5,常微分方程的数值解,大气科学中 常微分方程和偏微分方程的关系 1. 大气行星边界层(近地面具有湍流运动特性的大气薄层,11.5km), 埃克曼（V.W.Ekman）（瑞典）螺线的导出； 2. 1963年，美国气象学家Lorenz在研究热对流的不稳定问题时，使用高截断的谱方法，由Boussinesq流体的闭合方程组得到了一个完全确定的三阶常微分方程
5、组，即著名的Lorenz系统。,6,Lorenz系统,dx / dt = a (y - x) dy / dt = x (b - z) - y dz / dt = xy - c z,其中，a=10,(Prandtl number); b=28(Rayleigh number); c=8/3; (x,y,z)_0=(0.01;0.01;1e-10),7,8,9,10,11,Franceshini 将Navier-Stokes方程截断为五维的 截谱模型如下：,12,欧拉法折线法,常微分方程能直接进行积分的是少数，而多数是借助于计算机来求常微分方程的近似解； 有限差分法是常微分方程中数值解法中通 常有
6、效的方法； 建立差分算法的两个基本的步骤： 1. 建立差分格式，包括：a. 对解的存在域剖分；b. 采用不同的算法可得到不同的逼近误差截断误差（相容性）；c.数值解对真解的精度整体截断误差（收敛性）；d.数值解收敛于真解的速度；e. 差分算法舍人误差（稳定性）.,13,2.差分格式求解 将积分方程通过差分方程转化为代数方程求解，一般常用递推算法。 在常微分方程差分法中最简单的方法是Euler方法，尽管在计算中不会使用，但从中可领悟到建立差分格式的技术路线，下面将对其作详细介绍:,14,差分方法的基本思想“就是以差商代替微商”,考虑如下两个Taylor公式：,（1）,（2）,从（1）得到：,15
7、,从（2）得到：,从（1）-（2）得到：,从（1）+（2）得到：,16,对经典的初值问题,满足Lipschitz条件,保证了方程组的初值问题有唯一解。,17,一、算法构造：,0,t,u,T,1. 在求解域上等距离分割：,2. 在 有：,微分方程的精确解,差分方程的精确解,18,3. 应用时采用如下递推方式计算：,4. 例题,对初值问题,用Euler法求解，用,即，,19,5. Euler法的几何意义,0,t,在递推的每一步，设定,过点,作 的切线，该切线的,方程为：,即：,20,二、误差分析,构造算法后，这一算法在实际中是否可行呢？也就是说是否使计算机仿真,而不失真，这还需要进一步分析。,1.
8、 局部截断误差-相容性,为了分析分析数值方法的精确度，常常在,成立的假定下，,估计误差,这种误差称为“局部截断误差”，如图。,局部截断误差是以点,的精确解,为出发值，用数值方法推进到下一个点,而产生的误差。,21,整体截断误差是以点,的初始值,为出发值，用数值方法推进i+1步到点,，所得的近似值 与精确值 的偏差：,2.整体截断误差收敛性,称为整体截断误差。,22,特例，,若不计初始误差，即,则,即,3.舍入误差稳定性,假设一个计算机仅表示4个数字（小数点后面），,那么,计算,23,我们的要求是：最初产生的小误差在以后的计算中虽然会传递下去，但不会无限制的,扩大，这就是稳定性所描述的问题。下面
9、引进稳定性的概念：,设由初值,得到精确解 ，,由初值,得到精确解 ，,若存在常数,和充分小的步长,使得,则称数值方法是稳定的。,t,u,0,24,计算例题,其解析解为：,x = 0 0.2000 0.4000 0.6000 0.8000 1.0000 y = 1.0000 1.2000 1.3733 1.5315 1.6811 1.8269,25,26,三、改进的Euler法,将微分方程,在区间,上积分，得到,用梯形法计算积分的近似值，有,于是,这是一个隐式格式，一般需要用迭代法来求 ，而用显式的Euler法提供初值。,27,为了简化计算的过程，在此基础上进一步变为如下算法：,此式称为“改进的
10、Euler法。,接下来讨论其几何意义,预估,校正,其局部截断误差为,这个问题将在下节讨论。,28,t,u,0,29,Euler法、改进的Euler法和解析解的比较,30,四、（龙格-库塔）Runge-Kutta方法,简单的Euler法是建立在Taylor级数的一项展开；,改进的Euler法是以两项Taylor级数为基础建立的，如：,如果我们截取Taylor级数的更多项会得到什么样的求解方法呢？,两个德国数学家（C.Runge & M.kutta）以这种思想为基础建立了求解微分,方程的龙格-库塔方法。它是常微分方程数值解法中使用最为广泛的方法之一。,31,一般地，一个K阶的Runge-Kutta
11、方法可用下面的公式表示：,其中，,是待定的加权系数，,是待定的系数。,Euler法就是,的R-K法。,其系数的确定如下：将,展开成,的幂级数，并与微分方程的精确解,在点,的Taylor展开式相比较，使两者的前,项相同，这样确定的R-K法，,其局部截断误差为,，根据所得关于待定系数的方程组，求出它们的值后,代入公式，就成为一个,阶R-K方法。,32,例题 以二阶R-K法为例说明上述过程,把,代入,中，有,33,经比较得到,取 为自由参数：,从而得到不同的但都是二阶的R-K方法，对应的有中点法、Heun（亨）法,以及改进的Euler法。,34,基于相同的过程，通过比较五次Taylor多项式，得到更
12、加复杂的结果，给出了包含,13个未知数的11个方程。得到多组系数，其中常用的是以下四阶R-K法：,改进的Euler法、R-K法以及解析解的比较：,35,36,五、线性多步（Linear Multistep Method）法,1. 预备知识：插值多项式,插值是离散函数逼近的重要方法，利用它可通过函数在有限个点处的取值状况， 估算出函数在其他点处的近似值。,从几何上理解：对一维而言，已知平面上n1个不同点，要寻找一条n次多项式 曲线通过这些点。插值多项式一般常见的是拉格朗日插值多项式。,2. 气象应用,不均匀站点上的气象要素数据,均匀网格点上的数据,插值,3. 拉格朗日插值多项式,拉格朗日插值多项
13、式逼近可能是求插值节点不均匀的插值多项式的最简单的方法。,实验观察结果或原始测量数据的分布通常是非均匀的。,例如，四个点可以确定一个三次多项式，其拉格朗日形式为：,37,4. Adams-Bashforth（阿达姆斯贝雪福斯）公式,首先，用以下四个点对,进行三次Langrage插值：,则,于是，有,容易算出，,例如，我们可以算得,*2,*1,38,将（*2）代入（*1）得到Adams-Bashforth公式：,基于同样的计算过程，可以得到另外一个计算公式：,这称为Adams-Moulto公式。,预估,校正,39,偏微分方程数值解,主讲：王曰朋,40,一、区域的离散,1.,2.,3.,41,则函
14、数可表示为：,二、1.（一维）一、二阶导数的有限差分近似表达式,42,2.（二维）一、二阶偏导数的有限差分近似,43,3. 抛物型方程,初条：,精确解为,（以热传导或磁扩散方程为例）,初值问题,不论初始分布如何集中，它总在瞬间影响于无穷远，虽该影响随距离按指数衰减，然而它是以无限速度传播。此乃抛物型方程解的特征。,44,三、热传导方程（抛物方程）,1. 热传导方程的介绍,2.,离散化,（1）向前差分格式：,45,计算：,这是一个显式格式（四点格式）,每一层各个节点上的值是通过一个方程组求解得到的。这可以从 下面的计算过程看出来。,46,系数矩阵为,47,计算实例：,48,49,2. 向后差分格
15、式,当知道第n层上的 时，要确定第n+1层上各点值 必须通过求解一个 线性代数方程组。,50,其矩阵表达式如下：,51,这是一个古典四点向后差分格式。计算实例,52,53,3. Crank-Nicolson格式，亦称六点对称格式,54,55,56,4. Richardson格式,这是一个五点三层差分显式格式,57,讨论：,假若由于 的作用，导致差分方程的近似解设为：,于是，我们可得到差分格式的误差方程如下：,x,t,Richardson格式是不稳定的。,58,5. 稳定性判别,Von-Neumann 稳定性,在判断有限差分近似的稳定性方法中，以Von-Neumann方法使用较为 广泛，它仅适用
16、于线性常系数的有限差分近似。其过程如下：,首先，要研究的差分方程可写为：,如，,其次，对进行变量分离：,59,最后将,代入所考察的有限差分方程。,定义为放大系数,60,下面说明，在什么条件下能使,对所有的,成立。,从上式，我们看出，,61,. 交替显隐式格式,（）显式预测隐式校正格式,在n+1/2层上用古典显式格式计算出过度值，再在n+1层上用,古典隐格式校正预测值，即：,62,(2). 跳点格式,首先将网格点（j, n）按j+n等于偶数或奇数分成两组，分别称为偶数 网点和奇数网点。从 到的计算过程中，先在偶数网格点上 用古典显式格式计算，再在奇数网点上用古典隐格式计算，即：,63,三 双曲型方程,（a）一阶常系数线性双曲型方程,（b）二阶常系数线性双曲型方程（波动方程）,其中a为常数
