一类二维抛物型方程的ADI格式

【摘要】本文针对一类二维抛物型方程，建立了一个在空间和时间方向上均具有二阶精度的ADI格式，并分析其稳定性. 比较以往算法，此格式具有精度较高，无条件稳定等优点，同时，该方法通过求解两个线性代数方程得到原问题的解，避免了非线性迭代运算，提高了计算效率.

【关键词】二维抛物型方程；ADI格式；稳定性；截断误差

1.引言

抛物型偏微分方程在研究热传导过程、一些扩散现象及电磁场传播等许多问题中都有广泛的应用，对这一类方程数值解法的研究一直是科研工作者关注的热点问题之一，其中高精度的有限差分方法更是受到了越来越多的重视. 考虑如下的初边值问题[1]：

其中，为常数.

在文献[1]中对问题（1）建立了差分格式，格式的截断误差阶为.本文将在文献[1]的基础上进一步研究问题（1）的高效差分格式，建立了一个高精度的交替方向隐式差分格式（即ADI格式），提高了时间方向上的精度，并给出相应的稳定性分析。

2.差分格式的建立

为了得到问题（1）的有限差分格式，首先将求解区域进行网格剖分，结点为. 选取正整数L和N，并令为方向上的网格步长，为方向上的网格步长，记

假定第层的已知，则由第（Ⅰ）步，通过解一个三对角线性代数方程组求出，再由第（Ⅱ）步，再解一个三对角线性代数方程组即可求出. 所以，只要利用追赶法求解两个三对角线性代数方程组即可，此时计算量与储存量都较小. 另外，在处理边界条件时，为了提高精度，采取中心差分，这样会出现虚拟值，此时，只要再与格式中的方程联立，即可消去虚拟值[2].

3. 稳定性分析

下面采用von Newmann方法[3]对上述D格式进行稳定性分析. 一般地，低阶项不影响差分格式的稳定性，所以不妨略去项，并对（3）、（5）式消去中间变量得：

利用Taylor展开式求误差，可知此处建立的D格式的截断误差阶为.

参考文献：

(完整版)大连理工大学高等数值分析抛物型方程有限差分法

抛物型方程有限差分法 1. 简单差分法考虑一维模型热传导方程 (1.1) )(22x f x u a t u +??=??，T t ≤<0 其中a 为常数。)(x f 是给定的连续函数。（1.1）的定解问题分两类：第一，初值问题（Cauchy 问题）：求足够光滑的函数()t x u ,，满足方程（1.1）和初始条件：（1.2） ()()x x u ?=0,, ∞<<∞-x 第二，初边值问题（也称混合问题）：求足够光滑的函数()t x u ,，满足方程（1.1）和初始条件： ()13.1 ()()x x u ?=0,, l x l <<- 及边值条件 ()23.1 ()()0,,0==t l u t u ， T t ≤≤0 假定()x f 和()x ?在相应的区域光滑，并且于()0,0，()0,l 两点满足相容条件,则上述问题有唯一的充分光滑的解。

现在考虑边值问题（1.1），（1.3）的差分逼近取 N l h = 为空间步长，M T = τ为时间步长，其中N ，M 是自然数， jh x x j ==, ()N j ,,1,0Λ=； τ k y y k ==, ()M k ,,1,0Λ= 将矩形域G {}T t l x ≤≤≤≤=0;0分割成矩形网格。其中 ()j i y x ,表示网格节点； h G 表示网格内点（位于开矩形G 中的网格节点）的集合； h G 表示位于闭矩形G 中的网格节点的集合； h Γ表示h G -h G 网格边界点的集合。 k j u 表示定义在网点()k i t x ,处的待求近似解，N j ≤≤0，M k ≤≤0。注意到在节点()k i t x ,处的微商和差商之间的下列关系（(,)k j k j u u x t t t ????≡ ? ????）： ()() ()ττ O t u t x u t x u k j k j k j +??? ????=-+,,1 ()() ()2112,,ττ O t u t x u t x u k j k j k j +??? ????=--+ ()()()h O x u h t x u t x u k j k j k j +??? ????=-+,,1 ()() ()h O x u h t x u t x u k j k j k j +??? ????=--,,1 ()() ()2112,,h O x u h t x u t x u k j k j k j +??? ????=--+ ()()() ()2 222 11,,2,h O x u h t x u t x u t x u k j k j k j k j +???? ????=+--+ 可得到以下几种最简差分格式

抛物型方程的计算方法

分类号：O241.82 本科生毕业论文（设计）题目：一类抛物型方程的计算方法作者单位数学与信息科学学院作者姓名专业班级2011级数学与应用数学创新2班指导教师论文完成时间二〇一五年四月

一类抛物型方程的数值计算方法 (数学与信息科学学院数学与应用数学专业2011级创新2班) 指导教师摘要: 抛物型方程数值求解常用方法有差分方法、有限元方法等。差分方法是一种对方程直接进行离散化后得到的差分计算格式，有限元方法是基于抛物型方程的变分形式给出的数值计算格式.本文首先给出抛物型方程的差分计算方法，并分析了相应差分格式的收敛性、稳定性等基本理论问题.然后，给出抛物型方程的有限元计算方法及理论分析. 关键词：差分方法，有限元方法，收敛性，稳定性 Numerical computation methods for a parabolic equation Yan qian (Class 2, Grade 2011, College of Mathematics and Information Science) Advisor: Nie hua Abstract: The common methods to solve parabolic equations include differential method, finite element method etc. The main idea of differential method is to construct differential schemes by discretizing differential equations directly. Finite element scheme is based on the variational method of parabolic equations. In this article, we give some differential schemes for a parabolic equation and analyze their convergence and stability. Moreover, the finite element method and the corresponding theoretical analysis for parabolic equation are established. Key words: differential method, finite element method, convergence, stability

一维抛物线偏微分方程数值解法(附图及matlab程序)

一维抛物线偏微分方程数值解法（4）上一篇参看一维抛物线偏微分方程数值解法（3）（附图及matlab程序）解一维抛物线型方程(理论书籍可以参看孙志忠：偏微分方程数值解法） Ut-Uxx=0, 00) U(x,0)=e^x, 0<=x<=1, U(0,t)=e^t,U(1,t)=e^(1+t), 0

抛物形扩散方程的有限差分法及数值实例

偏微分方程数值解所在学院：数学与统计学院课题名称：抛物形扩散方程的有限差分法及数值实例学生姓名：向聘

抛物形扩散方程的有限差分法及数值实例 1.1抛物型扩散方程抛物型偏微分方程是一类重要的偏微分方程。考虑一维热传导方程： 22(),0u u a f x t T t x ??=+<≤?? (1.1.1) 其中a 是常数，()f x 是给定的连续函数。按照初边值条件的不同给法，可将(1.1.1)的定解分为两类：第一，初值问题（Cauchy 问题）：求足够光滑的函数()t x u ,，满足方程(1.1.1)和初始条件： ()()x x u ?=0,, ∞<<∞-x (1.1.2) 第二，初边值问题（也称混合问题）：求足够光滑的函数()t x u ,，满足方程(1.1.1)和初始条件： ()()x x u ?=0,, 0x l << (1.1.3) 及边值条件 ()()0,,0==t l u t u ， T t ≤≤0 (1.1.4) 假定()x f 和()x ?在相应的区域光滑，并且于()0,0，()0,l 两点满足相容条件,则上述问题有唯一的充分光滑的解。 1.2抛物线扩散方程的求解下面考虑如下热传导方程 22()(0.)(,)0(,0)()u u a f x t x u t u L t u x x ????=+????? ==??=??? (1.2.1) 其中，0x l <<，T t ≤≤0,a (常数)是扩散系数。取N l h = 为空间步长，M T =τ为时间步长，其中N ，M 是自然数，用两族

平行直线jh x x j ==, ()N j ,,1,0 =和k t t k τ ==, ()M k ,,1,0 =将矩形域 G {}T t l x ≤≤≤≤=0;0分割成矩形网格。其中 (),j k x t 表示网格节点；h G 表示网格内点（位于开矩形G 中的网格节点）的集合；h G 表示位于闭矩形G 中的网格节点的集合；h Γ表示h G -h G 网格边界点的集合。 k j u 表示定义在网点(),j k x t 处的待求近似解，N j ≤≤0，M k ≤≤0。现在对方程进行差分近似： (一) 向前差分格式 =-+τ k j k j u u 111 2 2(())k k k j j j j j j u u u a f f f x h +--++= (1.2.2) ()j j j x u ??==0， k u 0=k N u =0 (1.2.3) 计算后得： 111(12)k k k k j j j j j u ru r u ru f τ++-=+-++ (1.2.4) 其中，2 a r h τ = ，1,,1,0-=N j ，1,,1,0-=M k 。显然，这是一个四点显示格式，每一层各个节点上的值是通过一个方程组求解到的。方程组如下： 1000 121011000 232121000 3432310001121(12)(12)(12)(12)N N N N N u ru r u ru f u ru r u ru f u ru r u ru f u ru r u ru f ττττ----?=+-++?=+-++??=+-++? ???=+-++? (1.2.5) 若记 () T k N k k k u u u 1 21,,,-= u ，()()()()T N x x x 121,,,-=???? ，()()()()T N x f x f x f 121,,,-=τττ f 则显格式(1.2.4)可写成向量形式 10 ,0,1,,1 k k k M φ +?=+=-?=? u Au f u (1.2.6) 其中

维抛物线偏微分方程数值解法

一维抛物线偏微分方程数值解法（2）上一篇文章请参看一维抛物线偏微分方程数值解法（1）解一维抛物线型方程(理论书籍可以参看孙志忠：偏微分方程数值解法） Ut-Uxx=0, 00) U(x,0)=e^x, 0<=x<=1, U(0,t)=e^t,U(1,t)=e^(1+t), 0

ADI(交替方向隐格式)求解二维抛物方程(含matlab程序)

ADI 法求解二维抛物方程学校：中国石油大学（华东）学院：理学院：道德时间：2013.4.27 1、ADI 法介绍作为模型，考虑二维热传导方程的边值问题：（3.6.1）,0,,0(,,0)(,)(0,,)(,,)(,0,)(,,)0t xx yy u u u x y l t u x y x y u y t u l y t u x t u x l t ?=+<<>?? =??====? 取空间步长1h M ，时间步长0，作两族平行于坐标轴的网线： ,,,0,1, ,,j k x x jh y y kh j k M =====将区域0,x y l ≤≤分割成2M 个小矩形。第一个ADI 算法（交替方向隐格式）是Peaceman 和Rachford （1955）提出的。方法：由第n 层到第n+1层计算分为两步：（1）第一步： 12,1 2 n j k xx yy u +从n->n+，求u 对向后差分,u 向前差分，构造出差分格式为： 1 (3.6.1)11112222,,1,,1,,1,,1 2 2 1 222,,2-22=2 1()n n n n n n n n j k j k j k j k j k j k j k j k n n x j k y j k h h h τ δδ+ + + + +-+-+-+-+=+u u u u u u u u （ + ） u u （2）第二步：12,1 2 n j k xx yy u +从n+->n+1，求u 对向前差分,u 向后差分，构造出差分格式为： 2 (3.6.1)1111111222,,1,,1,,1,,1 2 2 1 2212,,2-22=2 1()n n n n n n n n j k j k j k j k j k j k j k j k n n x j k y j k h h h τ δδ+ + + +++++-+-++-+-+=+u u u u u u u u （ + ） u u 其中12 11 ,1,,1,0,1,2, ,()22n j k M n n n τ+=-=+=+上表表示在t=t 取值。假定第n 层的,n j k u 已求得，则由1(3.6.1)求出12 ,n j k u +，这只需按行 (1, ,1)j M =-解一些具有三对角系数矩阵的方程组；再由2(3.6.1)求出

抛物形扩散方程的有限差分法及数值实例(完整资料).doc

此文档下载后即可编辑偏微分方程数值解所在学院：数学与统计学院课题名称：抛物形扩散方程的有限差分法及数值实例学生姓名：向聘抛物形扩散方程的有限差分法及数值实例

1.1抛物型扩散方程抛物型偏微分方程是一类重要的偏微分方程。考虑一维热传导方程： 22(),0u u a f x t T t x ??=+<≤?? (1.1.1) 其中a 是常数，()f x 是给定的连续函数。按照初边值条件的不同给法，可将(1.1.1)的定解分为两类：第一，初值问题（Cauchy 问题）：求足够光滑的函数()t x u ,，满足方程(1.1.1)和初始条件： ()()x x u ?=0,, ∞<<∞-x (1.1.2) 第二，初边值问题（也称混合问题）：求足够光滑的函数()t x u ,，满足方程(1.1.1)和初始条件： ()()x x u ?=0,, 0x l << (1.1.3) 及边值条件 ()()0,,0==t l u t u ， T t ≤≤0 (1.1.4) 假定()x f 和()x ?在相应的区域光滑，并且于()0,0，()0,l 两点满足相容条件,则上述问题有唯一的充分光滑的解。 1.2抛物线扩散方程的求解下面考虑如下热传导方程

22()(0.)(,)0(,0)()u u a f x t x u t u L t u x x ????=+????? ==??=??? (1.2.1) 其中，0x l <<，T t ≤≤0,a (常数)是扩散系数。取N l h = 为空间步长，M T = τ为时间步长，其中N ，M 是自然数，用两族平行直线 jh x x j ==, ()N j ,,1,0Λ=和 k t t k τ ==, ()M k ,,1,0Λ=将矩形域G {}T t l x ≤≤≤≤=0;0分割成矩形网格。其中 (),j k x t 表示网格节点；h G 表示网格内点（位于开矩形G 中的网格节点）的集合；h G 表示位于闭矩形G 中的网格节点的集合；h Γ表示h G -h G 网格边界点的集合。 k j u 表示定义在网点(),j k x t 处的待求近似解， N j ≤≤0，M k ≤≤0。现在对方程进行差分近似： (一) 向前差分格式 =-+τ k j k j u u 111 2 2(()) k k k j j j j j j u u u a f f f x h +--++= (1.2.2) () j j j x u ??==0， k u 0 = k N u =0 (1.2.3) 计算后得： 111(12)k k k k j j j j j u ru r u ru f τ++-=+-++ (1.2.4) 其中，2 a r h τ = ，1,,1,0-=N j Λ，1,,1,0-=M k Λ。显然，这是一个四点显示格式，每一层各个节点上的值是通过一个方程组求解到的。方程组如下：

二维抛物线方程数值解法(ADI隐式交替法)方法

ADI隐式交替法三种解法及误差分析（一般的教材上只说第一种）理论部分参看孙志忠：偏微分方程数值解法注意： 1.最好不要直接看程序，中间很多公式很烦人的（一定要小心），我写了两天，终于写对了。 2.中间：例如r*(u(i-1,m1,k)+u(i+1,m1,k))形式写成分形式：r*u(i-1,m1,k)+r*u(i+1,m1,k)后面会出错，我也不是很清楚为什么，可能由于舍入误差，或者大数吃掉小数的影响。 3.下面有三个程序 4.具体理论看书，先仔细看书（孙志忠：偏微分方程数值解法）或者网上搜一些理论。 Matlab程序： 1.function [u u0 p e x y t]=ADI1(h1,h2,m1,m2,n) %ADI解二维抛物线型偏微分方程(P-R交替隐式，截断) %此程序用的是追赶法解线性方程组 %h1为空间步长，h2为时间步长 %m1,m2分别为x方向，y方向网格数，n为时间网格数 %p为精确解，u为数值解,e为误差 %定义u0(i,j,k)=u(i,j,k+1/2),因为矩阵中，i,j,k必须全为整数 x=(0:m1)*h1+0;%定义x0,y0,t0是为了f(x,t)~=0的情况% y=(0:m2)*h1+0; t=(0:n)*h2+0; t0=(0:n)*h2+1/2*h2; for k=1:n+1 for i=1:m2+1 for j=1:m1+1 f(i,j,k)=-1.5*exp(0.5*(x(j)+y(i))-t0(k)); end end end for i=1:m2+1 for j=1:m1+1 u(i,j,1)=exp(0.5*(x(j)+y(i))); end end for k=1:n+1 for i=1:m2+1 u(i,[1 m1+1],k)=[exp(0.5*y(i)-t(k)) exp(0.5*(1+y(i))-t(k))]; u0(i,[1 m1+1],k)=[exp(0.5*y(i)-t0(k)) exp(0.5*(1+y(i))-t0(k))] ; end end for k=1:n+1 for j=1:m1+1

抛物型方程求解

22 10,01,01(,0),01(0,),(1,),01 (,)x t t x t u u x t t x u x e x u t e u t e t u x t e ++??-=<<<≤??=≤≤==<≤= 运行：前向euler 法 [xx,tt,uh]=equationepaowu2('myfun','myfun1','myfun1','myfun2',1,[0,1],[0,1],[1/10,1/200]); function [xx,tt,uh]=equationepaowu2(myfun,myfun1,myfun2,myfun3,a,xxx,ttt,step) %利用差分方法求抛物型方程数值解； %myfun--方程右端f(x,t); %myfun1--u(x,0); %myfun2--u(t1,t); %myfun3--u(t2,t); %[x1,x2]--x 的取值范围； %[t1,t2]--t 的取值范围； %a-正常数 %h,tao-分别是x,t 方向的步长。 %—————————————————————— %激活函数 f=fcnchk(myfun); f1=fcnchk(myfun1); f2=fcnchk(myfun2); f3=fcnchk(myfun3); x1=xxx(1);x2=xxx(2); t1=ttt(1);t2=ttt(2); h=step(1);tao=step(2); %__________________________________ %划分网格,x1-nt+1行，nx+1列。 x=linspace(x1,x2,round((x2-x1)/h)+1); t=linspace(t1,t2,round((t2-t1)/tao)+1); nx=size(x,2); nt=size(t,2); [xx,tt]=meshgrid(x,t); %________________________________________ %赋初值及边值 size(x1) size(x) U0=zeros(size(xx));

热传导方程抛物型偏微分方程和基本知识

1. 热传导的基本概念 1.1温度场一物体或系统内部，只要各点存在温度差，热就可以从高温点向低温点传导，即产生热流。因此物体或系统内的温度分布情况决定着由热传导方式引起的传热速率（导热速率）。温度场：在任一瞬间，物体或系统内各点的温度分布总和。因此，温度场内任一点的温度为该点位置和时间的函数。〖说明〗若温度场内各点的温度随时间变化，此温度场为非稳态温度场，对应于非稳态的导热状态。若温度场内各点的温度不随时间变化，此温度场为稳态温度场，对应于稳态的导热状态。若物体内的温度仅沿一个坐标方向发生变化，且不随时间变化，此温度场为一维稳态温度场。 1.2 等温面在同一时刻，具有相同温度的各点组成的面称为等温面。因为在空间同一点不可能同时有两个不同的温度，所以温度不同的等温面不会相交。 1.3 温度梯度从任一点起沿等温面移动，温度无变化，故无热量传递；而沿和等温面相交的任一方向移动，温度发生变化，即有热量传递。温度随距离的变化程度沿法向最大。温度梯度：相邻两等温面间温差△t与其距离△n之比的极限。〖说明〗温度梯度为向量，其正方向为温度增加的方向，与传热方向相反。稳定的一维温度场，温度梯度可表示为：grad t = dt/dx

2. 热传导的基本定律——傅立叶定律物体或系统内导热速率的产生，是由于存在温度梯度的结果，且热流方向和温度降低的方向一致，即与负的温度梯度方向一致，后者称为温度降度。傅立叶定律是用以确定在物体各点存在温度差时，因热传导而产生的导热速率大小的定律。定义：通过等温面导热速率，与其等温面的面积及温度梯度成正比： q = dQ/ds = -λ·dT/dX 式中：q 是热通量（热流密度），W/m2 dQ是导热速率，W dS是等温表面的面积，m2 λ是比例系数，称为导热系数，W/m·℃ dT / dX 为垂直与等温面方向的温度梯度 “－”表示热流方向与温度梯度方向相反 3. 导热系数将傅立叶定律整理，得导热系数定义式： λ= q/(dT/dX) 物理意义：导热系数在数值上等于单位温度梯度下的热通量。因此，导热系数表征物体导热能力的大小，是物质的物性常数之一。其大小取决于物质的组成结构、状态、温度和压强等。导热系数大小由实验测定，其数值随状态变化很大。 3.1 固体的导热系数金属：35～420W/(m·℃)，非金属：0.2～3.0W/ (m·℃) 〖说明〗

抛物型方程

前言抛物型方程解的估计及其应用 1前言数学物理方程主要指从物理学及其它各门自然科学、技术科学中所产生的偏微分方程（有时也包括积分方程、微分方程等），它们反映了有关的未知变量关于时间的导数和关于空间变量的导数之间的制约关系．连续介质力学、电磁学、量子力学等等方面的基本方程属于数学物理方程的范围．它以具有物理背景的偏微分方程（组）作为研究的主要对象．它与其他数学分支及物理、化学等自然科学和工程技术的很多领域都有着广泛的联系，因此，无论在历史上还是在今天的现实生活中，它对于推动数学理论的发展，加强理论与实际的联系，帮助人们认识世界和改造世界都起着重要大的作用．微积分产生以后，人们就开始把力学中的一些问题，归结为偏微分方程进行研究．早在18世纪初，人们已经将弦线振动问题归结为弦振动方程，并探讨了它的解法．随后，人们又陆续了解了流体的运动、弹性体的平衡和振动、热传导、电磁相互作用、原子核和电子的相互作用、化学反应过程等等自然现象的基本规律，把它们写成偏微分方程的形式，并且求出了典型问题的解答，从而能通过实践，验证这些基本规律的正确性，显示了数学物理方程对于认识自然界基本规律的重要性．有了基本规律，人们还要利用这些基本规律来研究复杂的自然现象和解决复杂的工程技术问题，这就需要求出数学物理方程中许多特定问题的解答，随着计算机的出现及计算技术的发展，即使是相当复杂的问题，也可以计算出足够精确的数值来，这对于预测自然现象的变化(如气象预报)和进行各种工程设计（如机械强度的计算）都有着很重要的作用．在研究数学物理方程的同时，人们对偏微分方程的性质也了解得越来越多、越来越深入，形成数学中的一门重要分支——偏微分方程理论．它既有悠久的历史，又不断地更新着它的对象、内容和方法．它直接联系着众多自然现象和实际问题，不断地提出或产生需要解决的新课题和新方法．它所面临的数学问题多样而复杂，不断地促进着许多相关数学分支（如泛函分析、复变函数、微分几何、计算数学等）的发展，并从它们之间引进许多有力的解决问题的工具．因此，数学物理方程又是纯粹数学的

用显式格式求解二维抛物型偏微分方程

用显式格式求解二维抛物型偏微分方程 2010-05-14 10:41 function varargout=liu(varargin) T=1;a=1;h=1/32;dt=1/200; [X,T,Z]=chfenmethed(h,dt,a,T); mesh(X,T,Z(:,:,3)); shading flat; % xlabel('X','FontSize',14); % ylabel('t','FontSize',14); % zlabel('error','FontSize',14); % title('误差图'); function [X,Y,Z]=chfenmethed(h,dt,a,T); %求解下问题 %u_t-a*(u_xx+u_yy)=f(x,y,t) 0

r=a*dt/h^2; [X,Y]=meshgrid(x,y); Z=zeros(m,m,n); U=zeros(m,m,n); for i=1:m for j=1:m U(i,j,1)=d(x(i),y(j)); end end for j=2:n for k=1:m U(1,k,j)=g0(y(k),t(j)); U(m,k,j)=g1(y(k),t(j)); U(k,1,j)=h0(x(k),t(j)); U(k,m,j)=h1(x(k),t(j)); end end for k=2:n for i=2:m-1 for j=2:m-1 U(i,j,k)=U(i,j,k-1)+r*a*(U(i+1,j,k-1)+U(i-1,j,k-1)+U(i,j+1,k-1)... +U(i,j-1,k-1)-4*U(i,j,k-1))+f(x(i),y(j),t(k-1));

有限元解二维变系数抛物方程

%%%%% 真解u=sin(pi*x)*sin(pi*y)*sin(t) %%%%% 方程diff(u,t)-div(a(x,y)*grad(u))=f %%%%% a(x,y)=x^2+y^2+1 %%%%% f=sin(pi*x)*sin(pi*y)*cos(t)-2*pi^2*sin(pi*x)*sin(pi*y)*sin(t)*(x^2+y^2+1)+2*pi*sin(t)*(x*cos( pi*x)*sin(pi*y)+y*sin(pi*x)*cos(pi*y)) %clear all % clc %%%%finite element code for parabolic equation with constant coefficient %%%mesh%% node=[0,0;1,0;1,1;0,1]; elem=[2,3,1;4,1,3]; T=1; bdEdge=setboundary(node,elem,’Dirichlet’); n=input(‘Please input initial mesh:’); M=input(‘M=’); for i=1:n [node,elem,bdEdge]=uniformrefine(node,elem,bdEdge); end N=size(node,1); NT=size(elem,1); S=1/NT; r=1/M; A=zeros(N,N);; u=zeros(N,M+1); u1=zeros(N,1); f=inline(‘sin(pi*xx(1,1))*sin(pi*xx(1,2))*cos(t)-2*pi^2*sin(pi*xx(1,1))*sin(pi*xx(1,2))*sin(t)*(x x(1,1)^2+xx(1,2)^2+1)+2*pi*sin(t)*(xx(1,1)*cos(pi*xx(1,1))*sin(pi*xx(1,2))+y*sin(pi*xx(1,1))* cos(pi*xx(1,2)); a=inline(‘xx(1,1)^2+xx(1,2)^2+1’,’xx’); [lambda,weight]=quadpts(5); p=node’; T=elem’; totalEdge=[elem(:,[2,3]);elem(:,[3,1]);elem(:,[1,2])]; isBdEdge=reshape(bdEdge,3*NT,1); Dirichlet=totalEdge(isBdEdge==1),:); isBdNode=false(N,1); isBdNode(Dirichlet)=true; bdNode=find(isBdNode); freeNode=find(~isBdNode); for j=2:M+1 for i=1:NT F=zeros(N,1); F_ele=zeros(1,3);

有限差分法求解抛物型方程说明

有限差分法求解抛物型方程偏微分方程只是在一些特殊情况下，才能求得定解问题解的解析式，对比较复杂的问题要找到解的解析表达式是困难的，因此需采用数值方法来求解.有限差分法是一种发展较早且比较成熟的数值求解方法，只适用于几何形状规则的结构化网格.它在微分方程中用差商代替偏导数，得到相应的差分方程，通过解差分方程得到微分方程解的近似值.本章主要介绍有限差分法的基本思想，并给出一些具体的数值实例. §1 差分方法的基本思想有限差分法把偏微分方程的求解区域划分为有限个网格节点组成的网格，主要采用Taylor 级数展开等方法，在每个网格节点上用有限差分近似公式代替方程中的导数，从而建立以网格节点上的函数值为未知数的代数方程组. 有限差分格式，从格式的精度来划分，有一阶格式、二阶格式和高阶格式.从差分的空间形式来考虑，可分为中心格式和逆风格式.考虑时间因子的影响，差分格式还可以分为显格式、隐格式和显隐交替格式等.目前常见的差分格式，主要是上述几种格式的组合，不同的组合构成不同的差分格式. 泰勒级数展开法对有限差分格式的分类和公式的建立起着十分重要的作用.下面采用泰勒展开式导出一个自变量系统的若干有限差分表达式. 首先考虑单变量函数()u x ，如图1把区域x 离散为一批结点，记 0()(), =0,1,2,i i u x u x ih u i =+= 图1 单变量函数离散化函数()u x 在点i x 处的泰勒展开式为 23 ()()()()()2!3! i i i i i u x u x u x h u x u x h h h ''''''+=++ ++ （1）或 23 ()()()()()2!3! i i i i i u x u x u x h u x u x h h h ''''''-=-+ -+ （2）式(1)和(2)重新整理可得 2()()()()()2!3! i i i i i u x h u x u x u x u x h h h '''''+-'= --- （3）

一类二维抛物型方程的ADI格式

一类二维抛物型方程的ADI格式【摘要】本文针对一类二维抛物型方程，建立了一个在空间和时间方向上均具有二阶精度的ADI格式，并分析其稳定性. 比较以往算法，此格式具有精度较高，无条件稳定等优点，同时，该方法通过求解两个线性代数方程得到原问题的解，避免了非线性迭代运算，提高了计算效率. 【关键词】二维抛物型方程；ADI格式；稳定性；截断误差 1.引言抛物型偏微分方程在研究热传导过程、一些扩散现象及电磁场传播等许多问题中都有广泛的应用，对这一类方程数值解法的研究一直是科研工作者关注的热点问题之一，其中高精度的有限差分方法更是受到了越来越多的重视. 考虑如下的初边值问题[1]：其中，为常数. 在文献[1]中对问题（1）建立了差分格式，格式的截断误差阶为.本文将在文献[1]的基础上进一步研究问题（1）的高效差分格式，建立了一个高精度的交替方向隐式差分格式（即ADI格式），提高了时间方向上的精度，并给出相应的稳定性分析。 2.差分格式的建立为了得到问题（1）的有限差分格式，首先将求解区域进行网格剖分，结点为. 选取正整数L和N，并令为方向上的网格步长，为方向上的网格步长，记假定第层的已知，则由第（Ⅰ）步，通过解一个三对角线性代数方程组求出，再由第（Ⅱ）步，再解一个三对角线性代数方程组即可求出. 所以，只要利用追赶法求解两个三对角线性代数方程组即可，此时计算量与储存量都较小. 另外，在处理边界条件时，为了提高精度，采取中心差分，这样会出现虚拟值，此时，只要再与格式中的方程联立，即可消去虚拟值[2]. 3. 稳定性分析下面采用von Newmann方法[3]对上述D格式进行稳定性分析. 一般地，低阶项不影响差分格式的稳定性，所以不妨略去项，并对（3）、（5）式消去中间变量得：利用Taylor展开式求误差，可知此处建立的D格式的截断误差阶为. 参考文献：

扩散问题的偏微分方程模型,数学建模

第七节扩散问题的偏微分方程模型物质的扩散问题，在石油开采、环境污染、疾病流行、化学反应、新闻传播、煤矿瓦斯爆炸、农田墒情、水利工程、生态问题、房屋基建、神经传导、药物在人体内分布以及超导、液晶、燃烧等诸多自然科学与工程技术领域，十分普遍地存在着. 显然，对这些问题的研究是十分必要的，其中的数学含量极大. 事实上，凡与反应扩散有关的现象，大都能由线性或非线性抛物型偏微分方程作为数学模型来定量或定性地加以解决. ＭＣＭ的试题来自实际，是“真问题⊕数学建模⊕计算机处理”的“三合一”准科研性质的一种竞赛，对上述这种有普遍意义和数学含量高，必须用计算机处理才能得到数值解的扩散问题，当然成为试题的重要来源，例如，ＡＭＣＭ-90Ａ，就是这类试题；ＡＭＣＭ-90Ａ要研究治疗帕金森症的多巴胺（dopamine ）在人脑中的分布，此药液注射后在脑子里经历的是扩散衰减过程，可以由线性抛物型方程这一数学模型来刻划. ＡＭＣＭ-90Ａ要研究单层住宅混凝土地板中的温度变化，也属扩散（热传导）问题，其数学模型与ＡＭＣＭ-90Ａ一样，也是线性抛物型方程. 本文交代扩散问题建模的思路以及如何推导出相应的抛物型方程，如何利用积分变换求解、如何确定方程与解的表达式中的参数等关键数学过程，且以ＡＭＣＭ-90Ａ题为例，显示一个较细致的分析、建模、求解过程. §1 抛物型方程的导出设(,,,)u x y z t 是t 时刻点(,,)x y z 处一种物质的浓度. 任取一个闭曲面S ，它所围的区域是Ω，由于扩散，从t 到t t +?时刻这段时间内，通过S 流入Ω的质量为 2 221(cos cos cos )dSd t t t S u u u M a b c t x y z αβγ+????=++???? ??. 由高斯公式得 2222 221222()d d d d t t t u u u M a b c x y z t x y z +?Ω ???=++???? ???. （1）其中，222,,a b c 分别是沿,,x y z 方向的扩散系数. 由于衰减（例如吸收、代谢等），Ω内的质量减少为 2 2d d d d t t t M k u x y z t +?Ω =? ???，（2）其中2 k 是衰减系数. 由物质不灭定律，在Ω内由于扩散与衰减的合作用，积存于Ω内的质量为12M M -. 换一种角度看，Ω内由于深度之变化引起的质量增加为 3[(,,,)(,,,)]d d d d d d d . (3)t t t M u x y z t t u x y z t x y z u x y z t t Ω +?Ω =+?-?=????? ??? 显然312M M M =-，即

抛物型方程差分方法

偏微分方程数值解复习提纲一.基本内容:(1)椭圆型方程差分方法;(2)抛物型方程差分方法;(3)双曲型方程差分方法;(4)椭圆型方程的有限元方法. 二.基本概念: (1)显式和隐式差分格式,网格比和加密路径; (2)差分格式的截断误差、相容性、稳定性、收敛性、逼近精度阶和收敛阶; (3)双曲型方程(组)的特征与Riemann不变量,差分格式的依赖区域和CFL条件; (4)差分格式的增长因子和增长矩阵、振幅误差与相位误差、耗散与色散、群速度； (5)双曲守恒方程的弱解与激波传播速度； (6)守恒性与守恒型差分格式、有限体积法； (7)差分格式的Fourier分析与L2稳定性、最大值原理与L∞稳定性、实用稳定性和强稳定性、网格的P`e clet数; (8)椭圆边值问题的变分形式与弱解、强制边界条件与自然边界条件； (9)Galerkin方法与Ritz方法,协调与非协调有限元方法； (10)有限元与有限元空间,有限元插值算子与插值函数,有限元方程与有限元解； (11)有限元的仿射等价与等参等价,有限元剖分的正则性和拟一致性. 三.基本方法与技巧: (1)比较函数与利用最大值原理的误差分析; (2)Taylor展开、Fourier分析、最大值原理； (3)修正方程分析、能量法分析； (4)充分利用解的守恒性和特征,以及适当处理初始条件与边界条件； (5)Sobolev空间及其基本性质，如嵌入定理、迹定理,Poincar′e-Friedrichs不等式; (6)仿射等价、多项式不变算子、商空间与商范数、Sobolev空间半范数的关系; (7)Aubin-Nische技巧,bramble-Hilbert引理,双线性引理. 四.基本格式: (1)二维Poisson方程的五点差分格式; (2)抛物型方程的显式差分格式、隐式差分格式、Crank-Nicolson格式和θ-方法； (3)具有热守恒性质的格式； (4)ADI格式与LOD格式； (5)双曲型方程的迎风格式、Lax-Wendro?格式、盒式格式和蛙跳格式；

抛物形扩散方程的有限差分法与数值实例

偏微分方程数值解所在学院:数学与统计学院课题名称:抛物形扩散方程的有限差分法及数值实例学生：向聘

抛物形扩散方程的有限差分法及数值实例 1.1抛物型扩散方程抛物型偏微分方程是一类重要的偏微分方程。考虑一维热传导方程: 2 u a 2 x 其中a 是常数，f(x)是给定的连续函数。按照初边值条件的不同给法，可将(1.1.1) 的定解分为两类：第一，初值问题(Cauchy 问题)：求足够光滑的函数u x, t ，满足方程(1.1.1) 和初始条件： u x, 0 x , x (1.1.2) 第二，初边值问题(也称混合问题)：求足够光滑的函数u x,t ，满足方程(1.1.1) 和初始条件： u x, 0 x , 0 x l (1.1.3) 及边值条件 u0,t ul,t 0， 0 t T (1.1.4) 假定f X 和 x 在相应的区域光滑，并且于0,0，l,0两点满足相容条件，则上述问题有唯一的充分光滑的解。 1.2抛物线扩散方程的求解 F 面考虑如下热传导方程 u(0.t) u(L,t) 0 u(x,0) (x) 其中，0 x l ，0 t T, a (常数)是扩散系数。 f(x),0 t T (1.1.1) (1.2.1)

取h秸为空间步长, M为时间步长，其中N，M是自然数，用两族

平行直线x X j jh , j 0,1, , N 和t t k k , k 0,1, , M 将矩形域 G 0 x l; 0 t T 分割成矩形网格。其中 X j ,t k 表示网格节点；G h 表示网格点（位于开矩形G 中的网格节点）的集合；G h 表示位于闭矩形G 中的网格节点的集合；h 表示G h - G h 网格边界点的集合。 u k 表示定义在网点x j ,t k 处的待求近似解，0 现在对方程进行差分近似: （一）向前差分格式若记则显格式（1.2.4）可写成向量形式 U 其中 U j j X j , k k U 0 =U N =0 (1.2.3) 计算后得： 1 U j rU j 1 (1 2r)U k k r V 1 f j (1.2.4) 其中，r 各，j 0,1, , N 1，k h 0,1, ,M 1。显然, 这是一个四点显示格式，每一层各个节点上的值是通过一个方程组求 k 1 k U j U j h 2 k k k U j 1 2U j U j 1 a (f j f (X j )) 解到的。方程组如下: (1.2.2) 0 rU 2 (1 2r )U 0 0 叫 f 1 0 rU a (1 2r)U 0 rU 1 f 2 (1 0 rU 4 2r)U a rU 2 f 3 M 1 U N 1 0 「U N (1 2「)U N 1 「U N 2 f N (1.2.5) k k k k T u U 1 ,U 2, , U N 1 X 1 , X 2 , T X N 1 ， f f X ! , f X 2 , T X N 1 1 u Au 0,1,L ,M 1 (1.2.6)