偏微分方程的MATLAB解法.docx

引 言偏微分方程定解问题有着广泛的应用背景。人们用偏微分方程来描述、解释或者预见各种自然现象，并用于科学和工程技术的各个领域fll。然而，对于广大应用工作者来说，从偏微分方程模型出发，使用有限元法或有限差分法求解都要耗费很大的工作量，才能得到数值解。现在，MATLAB PDEToolbox已实现对于空间二维问题高速、准确的求解过程。偏微分方程如果一个微分方程中出现的未知函数只含一个自变量，这个方程叫做常微分方程，也简称微分方程；如果一个微分方程中出现多元函数的偏导数，或者说如果未知函数和几个变量有关，而且方程中出现未知函数对几个变量的导数，那么这种微分方程就是偏微分方程。常用的方法有变分法和有限差分法。变分法是把定解问题转化成变分问题，再求变分问题的近似解；有限差分法是把定解问题转化成代数方程，然后用计算机进行计算；还有一种更有意义的模拟法，它用另一个物理的问题实验研究来代替所研究某个物理问题的定解。虽然物理现象本质不同，但是抽象地表示在数学上是同一个定解问题，如研究某个不规则形状的物体里的稳定温度分布问题，由于求解比较困难，可作相应的静电场或稳恒电流场实验研究，测定场中各处的电势，从而也解决了所研究的稳定温度场中的温度分布问题。随着物理科学所研究的现象在广度和深度两方面的扩展，偏微分方程的应用范围更广泛。从数学自身的角度看，偏微分方程的求解促使数学在函数论、变分法、级数展开、常微分方程、代数、微分几何等各方面进行发展。从这个角度说，偏微分方程变成了数学的中心。一、MATLAB方法简介及应用1.1 MATLAB简介MATLAB是美国MathWorks公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和Simulink两大部分。1.2 Matlab主要功能数值分析数值和符号计算工程与科学绘图控制系统的设计与仿真数字图像处理数字信号处理通讯系统设计与仿真财务与金融工程1.3 优势特点1) 高效的数值计算及符号计算功能，能使用户从繁杂的数学运算分析中解脱出来；2) 具有完备的图形处理功能，实现计算结果和编程的可视化；3) 友好的用户界面及接近数学表达式的自然化语言，使学者易于学习和掌握；4) 功能丰富的应用工具箱(如信号处理工具箱、通信工具箱等) ，为用户提供了大量方便实用的处理工具。1.4 MATLAB 产品族可以用来进行以下各种工作数值分析数值和符号计算工程与科学绘图控制系统的设计与仿真数字图像处理技术数字信号处理技术通讯系统设计与仿真财务与金融工程管理与调度优化计算（运筹学）MATLAB 的应用范围非常广，包括信号和图像处理、通讯MATLAB在通讯系统设计与仿真的应用、控制系统设计、测试和测量、财务建模和分析以及计算生物学等众多应用领域。附加的工具箱（单独提供的专用MATLAB函数集）扩展了MATLAB 环境，以解决这些应用领域内特定类型的问题。二、Laplacian算子简介Laplacian算子：Poisson方程(ellipptic)：Laplacian算子的特征值问题： Heat equation(parabolic):Wave equation(hyperbolic): 五点离散:Poisson方程离散： 特征值问题：热方程：波动方程：波动方程：椭圆方程： 特征值方程： 热方程： 波动方程： 热方程： 波动方程： 三、Matlab解偏微分方程解偏微分方程不是一件轻松的事情，但是偏微分方程在自然科学和工程领域应用很广，因此，研究解偏微分方程的方法、开发解偏微分方程的工具是数学和计算机领域中的一项重要工作。MATLAB提供了专门用于解二维偏微分方程的工具箱，使用这个工具箱，一方面解偏微分方程，另一方面，可以让我们学习如何把求解数学问题的过程与方法工程化。应当承认，我们国家在数学软件的开发方面还比较落后，MATLAB是当今世界上最好的数学软件之一，通过对这个软件的认识，有助于研发我们自己的数学软件。MATLAB的偏微分方程工具箱名字叫pdetool，它采用有限元法解偏微分方程。用这个工具箱可以解如下方程。椭圆方程抛物线方程双曲线方程特征值方程所有的方程都在二维平面 域上。方程中， 是Laplace算子，u是待解的未知函数，c, a, f是已知的实值标量函数，d是已知的复值函数，是未知的特征值。MATLAB提供了两种方法解决PDE问题，一是pdepe函数，它可以求解一般的PDEs，具有较大的通用性，但只支持命令行形式调用。二是PDE工具箱，可以求解特殊PDE问题，但有较大的局限性，比如只能求解二阶PDE问题，并且不能解决偏微分方程组。它提供了GUI界面，可以从繁杂的编程中解脱出来，同时还可以通过FileSave As直接生成M代码。3.1函数解法pdepe函数介绍它的调用格式为sol=pdepe(m，pdefun，pdeic，pdebc,x,t)输入参数:pdefun：是PDE的问题描述函数 pdebc：是PDE的边界条件描述函数 pdeic：是PDE的初值条件输出参数：sol：是一个三维数组，sol(:,:,i)表示ui的解，换句话说uk对应x(i)和t(j)时的解为sol(i,j，k)，通过sol，我们可以使用pdeval()直接计算某个点的函数值。实例讲解例：试求解下面的偏微分 其中， ，且满足初始条件及边界条件解：（1）对照给出的偏微分方程，根据标注形式，则原方程可以改写为 可见m=0,且目标PDE函数functionc，f，s-pdefun(x，t，u，du)c=1；1；f_0024*du(1)；017*du(2)； temp=u(1)一u(2)；s=一1；1*(exp(573*temp)一exp(-1146*temp)；（2）边界条件改写为下边界上边界边界条件函数functionpaqa，pb，qb-pdebe(xa，ua，xb，ub，t)a表示下边界，b表示上边界pa=0；ua(2)；qa=1；0 ；pb=ub(1)一1；O；qb=O；1；(3)初值条件改写为初值条件函数 function u0=pdeic(x)u0=1；0；(4)最后编写主调函数 ticm=0；x=0 0.005 0.01 0.05 0.1 0.2 0.5 0.7 0.9 0.95 0.99 0.995 1；T=0 0.005 0.01 0.05 0.1 0.5 1 1.5 2；sol=pdepe(m，pdefun，pdeic，pdebc，x，t)；U1=sol(：，：，1)；u2=sol(：，：，2)； figuresurf(x，t，u1)title(ul(x，t)xlabel(Distance x)ylabel(Time t)figuresurf(x，t，u2)tiffe(u2(x，t)xlabel(Distance x)ylabel(Time t)结果图32偏微分方程的pdetool解法3.2.1偏微分方程工具箱的功能偏微分方程工具箱(PDE Toolbox)提供了研究和求解空间二维偏微分方程问题的一个强大而又灵活实用的环境。PDE Toolbox的功能包括：(1)设置PDE (偏微分方程)定解问题，即设置二维定解区域、边界条件以及方程的形式和系数；(2)用有限元法 (FEM) 求解PDE数值解；(3)解的可视化。无论是高级研究人员还是初学者，在使用PDE Too1box时都会感到非常方便。只要PDE定解问题的提法正确，那么，启动MATLAB后，在MATLAB工作空间的命令行中键人pdetool，系统立即产生偏微分方程工具箱(PDE Toolbox)的图形用户界面(Graphical User Interface，简记为GUI)，即PDE解的图形环境，这时就可以在它上面画出定解区域、设置方程和边界条件、作网格剖分、求解、作图等工作，详见14节中的例子。我们将在第二章详细介绍GUI的使用，在第二章给出大量典型例子和应用实例。除了用GUI求解PDE外，也可以用M文件的编程计算更为复杂的问题。 3.2.2PDE Toolbox求解的问题及其背景3.2.2.1方程类型PDE Toolbox求解的基本方程有椭圆型方程、抛物型方程、双曲型方程、特征值方程、椭圆型方程组以及非线性椭圆型方程。椭圆型方程：，椭圆型方程：其中是平面有界区域，c,a,f以及未知数u是定义在上的实（或复）函数。抛物型方程：双曲型方程：.特征值方程：其中d是定义在上的复函数，是待求特征值。在抛物型方程和双曲型方程中，系数c,a,f和d可以依赖于时间t。可以求解非线性椭圆型方程：其中c,a,f可以是未知函数u的函数。还可以求解如下PDE方程组； 利用命令行可以求解高阶方程组。对于椭圆型方程，可以用自适应网格算法，还能与非线性解结合起来使用。另外，对于Poission方程还有一个矩形网格的快速求解器。3.2.2.2边界条件(1)Dirichlet条件 ：(2)Neumann 条件： 其中是的边界上的单位外法向量，和是定义在上的函数。对于特征值问题仅限于齐次条件：和。对于非线性情形系数和可以依赖于u；对于抛物型方程和双曲型方程，系数可以依赖于时间t。对于方程组情形，边界条件为(1)Dirichlet条件：(2)Neumann条件： (3)混合边界条件为: 其中的计算要使得Dirichlet条件满足。在有限元法中，Dirichlet条件也称为本质边界条件，Neumann条件称为自然边界条件。3.3 如何使用FDE Toolbox3.3.1定解问题的设置原简单的办法是在PDE Tool上直接使用图形用户界面(GUl)。设置定解问题包括三个步骤：(1)Draw模式：使用CSG(几何结构实体模型)对话框画几何区域，包括矩形、圆、椭圆和多边形，也可以将它们组合使用。(2)Boundary模式：在各个边界段上给出边界条件，(3)PDE模式：确定方程的类型、系数c，a，f和d c。也能够在不同子区域上设置不同的系数(反映材料的性质)。3.3.2解PDE问题用GUI解PDE问题主要经过下面两个过程（模式）(1)Mesh模式；生成网格自动控制网格参数。(2)Solve模式：对于椭圆型方程还能求非线性和自适应解。对于抛物型和双曲型力程设置初始边值条件后能求出给定t时刻的解。对于特征值问题，能求出给定区间内的特征值；求解后可以加密网格再求解。3.3.3使用Toolbox求解非标准的问题对于非标准的问题。可以用PDE Too1box的函数。或者用FEM(有限元法)求解更为复杂的问题。3.3.4 计算结果的可视化从GUI能够使用Plot模式实现可视化。可以使用Color, Height和Vector等作图。对于抛物型和双曲型方程，还可以生成解的动画。这些操作通过命令行都很容易实现。3.3.5 应用领域在应用界面提供了丁如下应用领域结构力学平面应力问题结构力学平面应变问题静电场问题静磁场问题交流电磁场问题直流导体介质问题热传导问题3.4解偏微分方程的一个例子解Poisson方程，边界条件为齐次Dirichlet类型。第一步：启动MATLABl, 键入pdetool，按回车键确定便可启动GUI，然后在Options菜单下选择Grid命令，打开栅格, 栅格的使用，能使用户容易确定所绘图形的大小，如图11 1-1第二步：分步完成平面几何造型：R1-C1-E1+R2+C2。用菜单或快捷工具，分别画矩形R1、矩形R2、椭圆E1、圆C1、圆C2。画圆时，首先选中椭圆工具，按鼠标右键并拖动即可、或者在按ctrI的同时，拖动鼠标也可绘制圆。然后在Set formula栏，进行编辑并用算术运算将将图形对象名称连接起来，删除默认的表达式键入R1-C1-E1+R2+C2，按等号健得到所需图形。若需要，还可进行储存形成M文件。选择Boundary菜单中Boundary Mode命令，进入边界模式。单击Boundary菜单中Remove A11 Subdomain Borders选项，去除子域边界。如果想将几何信息和边界信息进行存储，应选择Boundary菜单中的ExPort Decomposed GeometryBoundary Conds命令，将它们分别储存于g,b变量中, 通过MATLAB形成M文件。第三步：选取边界单击Boundary菜单中Specify Bounddy Conditions选项，打开Boundary conditlons对话框，输入边界条件，如图14。本例取缺省条件。即将全部边界设为齐次Dirichlet条件，边界颜色显示为红色。第四步：选择PDE菜单中PDE Mode命令，进入PDE模式。单击PDE菜单中PDE Specification选项，打开PDE对话框，设置方程类型。本例取缺省设置，类型为椭圆型，参数c,a,f分别为1，0，10。第五步：选择Mesh菜单中Initialize Mesh命令，进行网格剖分。第六步：选择Mesh菜单中Refine Mesh命令，对网格加密。第七步：选择Solve菜单中So1ve PDE命令，解偏微分方程并显示图形解。第八步：单击Plot菜单中Parameters选项，打开Plot selection对话框，选中Color, Height (3D Plot)和Show mesh三项。然后单击Plot按钮，显示三维图形解。第九步：如果要画等值线图和矢量场图，单击Plot菜单中Parameters选项，打开Plot Selection对话框选中Contour和Arrows两项。然后单击P1ot按钮，可显示解的等值线图和矢量场图。3.4.1PDE Toolbox菜单1.File菜单（如图1-1）图1-1New 新建一个几何结构实体模型（Constructive SolidGeomery,简记为CSG），默认文件名为“Untitled”。Open 从硬盘装载M文件Save 将在GUI内完成的成果储存到一个M文件中。Save As 将在GUI内完成的成果储存到另外一个M文件中。Print 将PDE工具箱完成的图形送到打印机内进行硬拷贝。Exit 退出PDE工具图形用户界面。2、Edit菜单（如图1-2）图1-2Undo 在绘制多边形时退回到上一步操作。Cut 将已选实体剪切到剪贴板上。Copy 将已选实体拷贝到剪贴板上。Paste 将剪贴板上的实体粘贴到当前几何结构实体模型中。Clear 删除已选的实体。Select All 选择当前几何结构实体造型CSG中的所有实体及其边界和字域。3 Options菜单（如图1-3）图1-3Grid 绘图时打开或关闭栅格。Grid Spacing 调整栅格的大小。Snap 打开或关闭捕捉栅格功能。Axes Limits 设置绘图轴的坐标范围。Axes Equal 打开或关闭绘图方轴。Turn off Toolbar Help 关闭工具栏按钮的帮助信息。Zoom 打开或关闭图形缩放功能。Application 选择应用的模式。Refresh 重新显示PDE工具箱中的图形实体。4、Draw菜单（如图1-4）图1-4Draw Mode 进入绘图模式。Rectangle/square 以角点方式画矩形/方行（Ctrl+鼠标）。Rectangle/square（centered） 以中心方式画矩形/方行（Ctrl+鼠标）。Ellipse/circle 以矩阵角点方式画椭圆/圆（Ctrl+鼠标）。Ellipse/circle（centered） 以中心方式画椭圆/圆（Ctrl+鼠标）。Polygon 画多边形，单击鼠标右键可封闭多边形。Rotate 旋转已选的图形。Export Geometry Description,Set Formula，Labels 将几何描述矩阵gd、公式设置字符sf和标识空间矩阵ns输出到主工作空间去。单击Draw菜单中Rotate选项，可打开Rotate比对活框，通过输入旋转的角度，可使选择的物体按输入的角度逆时针旋转。旋转中心的选择如果缺省，则为图形的质心，也可以输入旋转中心坐标。5、Boundary菜单（如图1-5）图1-5Boundary Mode 进入边界模式。Specify Boundary Conditions 对于已选的边界输入条件，如果没有选择边界，则边界条件适用于所有的边界。Show Edge Labels 显示边界区域标识开关，其数据是分解几何矩阵的列数。Show Subdomain Labels 显示子区域标识开关，其数据是分解几何矩阵中的子域数值。Remove Subdomain Border 当图形进行布尔运算时，删除已选取的子域边界。Remove All Subdomain Borders 当图形进行布尔运算时，删除所有的子域边界。Export Decomposed Geometry，Boundary Conds 将分解几何矩阵g、边界条件矩阵b输出到主工作空间。选择Boundary菜单中Specify Boundary Conditions命令可定义边界条件。在打开的Boundary condition对话框，可对已选的边界输入边界条件。共有如下三种不同的条件类型：NeMmann条件 这里边界条件是由方程系数q和g确定的，在方程组的情况下(换成方程组模式)，q是22矩阵，g是2x1矢量。Dirichlet条件 u定义在边界上，边界条件方程是价h*u=r，这里h是可以选样的权因子(通常为1)。在方程组情况下，h是2x2矩阵，r是2x l矢量，混合边界条件(仅适合于方程组情形) 它是Dirichlet和Neumann的混合边界条件，q是2x 2矩阵，g是2x1矢量，h是1x 2矢量，r是一个标量。6、PDE菜单（如图1-6）图1-6PDE Mode 进入偏微分方程模式。Show Subdomain Labels 显示子区域标识开关。PDE Specification 调整PDE参数和类型。Export PDE Coefficients 将当前PDE参数c,a,f,d输出到主工作空间，其参数变量为字符类型。7、Mesh菜单（如图1-7）图1-7Mesh Mode 输入网格模式。Initialize Mesh 建立和显示初始化三角形网格。Refine Mesh 加密当前三角型网格。Jiggle Mesh 优化网格。Undo Mesh Change 退回上一次网格操作。Display Triangle Quality 用01之间数字化的颜色显示三角形网格的质量，大于0.6的网格可接受的。Show Node Labels 显示网格节点标识开关，节点标识数据是点矩阵p的列。Show Triangle Labels 显示三角形网格标识开关，三角形网格标识数据是三角形矩阵t的列。Parameters 修改网格生成参数。Export Mesh 输出节点矩阵p、边界矩阵e和三角形矩阵t到主工作空间。8、Solve菜单（如图1-8）图1-8Solve PDE 对当前的几何结构实体CSG、三角形网格和图形解偏微分方程。Parameters 调整PDE的参数。Export Solution 输出PDE的解矢量u。如果可行，将计算的特征值1输出到主工作空间。9、Plot菜单（如图1-9）图1-9Plot Solution 显示图形解。Parameters 打开绘图方式对话框。Export Movie 如果动画被录制了，则动画矩阵M将输出到主工作空间。10、Window菜单从Window菜单项，可选择当前打开的所有的MATLAB图形窗口，被选择的窗台移至前台。11、Help菜单Help 显示帮助信息About 显示版本信息主菜单下是工具栏，工具栏中喊有许多工具图标按钮，可提供快速、便捷的操作方式。从左到右5个按钮为绘图模式按钮，紧接着的6个为边界、网格、解方程和图形显示控制功能按钮，最右边的为图形缩放功能键。（如图1-10） 图1-10 以角点方式画矩形/方行（Ctrl+鼠标）。 以中心方式画矩形/方行（Ctrl+鼠标）。 以矩形角点长轴方式画椭圆/圆（Ctrl+鼠标）。 以中心方式画椭圆/圆（Ctrl+鼠标）。 画多边形，按右键可封闭多边形。 进入边界模式。 打开PDE Specification（偏微分方程类型）对话框。 初始化三角形网格。 加密三角形网格。 解偏微分方程。 打开Plot Selection对话框，确定后给出解的三维图形。 为显示缩放切换按钮。求解PDE问题主要有两种方法，一种是使用图形用户界面，另一种是采用命令行编程。前者直观简便，而后者更为灵活。3.4.2求解椭圆方程及抛物型方程的例子例1：单位圆上的Poisson方程边值问题：这一问题的精确解为：若使用图形用户界面（Graphical User Interface,简记为GUI），则首先在MATLAB的工作窗口中键入pdetool，按回车键确定，于是出现PDE Toolbox窗口。如果需要坐标网格，单击Options菜单下的Grid选项即可。下面分步进行操作。（i）画区域圆 单击工具，大致在（0，0）位置单击鼠标右键同时拖拉鼠标到适当位置松开，绘制圆。为了保证所绘制的圆是标准的单位圆，在所绘圆上双击，打开Object Dialog对话框，精确地输入圆心坐标X-center为0、Y-cebter为0及半径Radius为1，然后单击OK按钮，这样单位远已画好。（ii）设置边界条件 单击工具，图形边界变红，逐段双击边界，打开Boundary Condition对话框，输入边界条件。对于同一类型的边界，可以按Shift键，将多个边界同时选择，统一设置边界条件。本题选择Dirichlet条件，输入h为1，r为0，然后单击OK按钮。也可以单击Boundary菜单中Specify Boundary Conditions选项，打开Boundary Condition对话框输入边界条件，如图2-1。（iii）设置方程 单击PDE菜单中PDE Specification选项，打开PDE Specification对话框，选项方程类型。本题单击Elliptic，输入c为1，a为0，f为1，然后单击OK按钮，如图2-2。图2-1图2-2（iv）网格剖分 单击工具，或者单击Mesh菜单中Initialize Mesh选项，可进行初始网格剖分，这时在PDE Toolbox窗口下方的状态栏内显示初始问网格的节点数和三角形单元数。本题节点数为144个，三角形单元数为254个。如果需要网格加密，再单击，或者单击Mesh菜单中Refine Mesh选项，这时节点数变为541个，三角形单元数为1016个，如此还可继续加密。（v）解方程 单击工具，或者单击Solve菜单中Solve菜单中Solve PDE选项，可显示方程色彩解。如果单击Plot菜单中Parameters选项，出现Plot Selection对话框，如图2-3，从中可以选择Color，Contour，Arrows,Deformed mesh，Height（3-D polt），还可以设置等值线的数目等。本例中选择Color，Contour，Height（3-D polt）和Show mesh四项，然后单击Plot按钮，方程的图形解如图2-4所示。除了作定解问题解u的图形外，也可以作|grad u|，|cgrad u|等图形。图2-3图2-4（vi）与精确解作比较 单击Plot菜单中Parameters选项，打开Plot Selection对话框，在Height（3-D plot）行Property下拉框中选user entry，且在该行的User entry输入框中键入u-（1-x.2-y.2）/4,单击lot按钮就可以看到解的绝对误差图形，如图-可见在边界处误差为。图2-5（vii）输出网格节点的编号、单元编号以及节点坐标 单击Mesh菜单中Show Node Labels选项，再单击工具或，即可显示节点编号。若要输出节点坐标，只需单击Mesh菜单中Export Mesh选项，这时打开的Export对话框中默认值为p，e，t，这里p，e，t分别表示points（点）、edges（边）、triangles（三角形）数据的变量，单击OK按钮。然后在MATLAB命令窗口键入p，按回车键确定，即可显示出节点按编号排列的坐标（二维数组）；键入e，按回车键，则显示边界线段数据矩阵（7维数组）；输入t，按回车键，则显示三角形单元数据矩阵（4维数组）。（viii）输出解的数值 单击Solve菜单中Export Solution选项，在打开的Export对话框中输入u，单击OK按钮确定。再在MATLAB命令窗口中输入u，按回车确定，即显示按节点编号排列的解的数值。我们也可以用MATLAB程序求解PDE问题，同时显示解的图形：p,e,t=initmesh(circleg,hmax,1);Error=;err=1;While err0.001,p,e,t=refinemesh(circleg,p,e,t);U=assempde(circleb1,p,e,t,1,0,1);Exact=(1-p(1,:).2-p(2,:).2)/4;Err=norm(u-exact,inf);Error=error err;EndPdemesh(p,e,t)Pdesurf(p,t,u)Pdesurf(p,t,u-exact)通过命令行键入help+命令函数，如help pdemesh，按回车键，可以调入有关命令函数的定义、参数格式等帮助信息。例2：求解抛物型方程的例子考虑一个带有矩形孔的金属板上的热传导问题。板的左边保持在100，板的右边热量从板向环境空气定常流动，其他边及内孔边界保持绝缘。初始时板的温度为0，于是概括为如下定解问题：域的外边界顶点坐标为（-0.5，-0.8），（0.5，-0.8），（0.5，0.8），（-0.5，0.8）。内边界顶点坐标为（-0.005，-0.4），（0.05，-0.4），（0.05，0.4），（-0.05，0.4）。使用GUI求解这一问题。在PDE Toolbox窗口的工具栏中选择Generic Scalar模式。（i）区域设置 单击工具，在窗口拖拉出一个矩形，双击矩形区域，在Object Dialog对话框中输入Left为-0.5，Bottom为-0.8，Width为1，Height为1.6，单击OK按钮，显示矩形区域R1。用同样方法作内孔R2，只要设置Left=-0.05，Bottom=-0.4，Width=0.1，Height=0.8即可。然后在Set formula栏中键入R1-R2。（ii）设置边界条件 单击，使边界变红色，然后分别双击每段边界，打开Boundary Cvondition对话框，设置边界条件。在左边界条件。在左边界上，选择Dirichlet条件，输入h为1，r为100；右边界上，选择Neumann条件，输入g为-1，q为0；其他边界上，选择Neumann条件，输入g为0，q为0。（iii）设置方程类型 单击，打开PDE Specification对话框，设置方程类型为Parabolic（抛物型），d=1，c=1，a=0，f=0，单击OK按钮。（iv）网格剖分 单击，或者加密网格，单击。（v）初值和误差的设置 单击Solve菜单中Parameters选项，打开Solve Parameters对话框，输入Time为0:5，u（t0）为0，Relative tolerance为0.01，Absolute tolerance为0.001，然后单击OK按钮。（vi）数值解的输出 单击Solve菜单中Export Solution选项，在Export对话框中输入u，单击OK按钮。再在MATLAB命令窗口中键入u，按回车键，这时显示出按节点编号的数值解。（vii）解的图形 单击Plot菜单中Parameters选项，打开Plot Selection对话框，选Color，Contour，Arrows，单击Plot按钮，窗口显示出Time=5时解的彩色图形，如图2-6。图2-6MATLAB程序：p,e,t=initmesh(crackg);u=parabolic(0,0:0.5:5,crackb,p,e,t,1,0,0,1);pdeplot(p,e,t,xydata,u(:,11),mesh,off,colormap,hot)结 论偏微分方程是一门实用性很强的学科。对于理论研究和实际应用中提出的偏微分方程问题，由于其边界和边界条件复杂等原因，寻求解析解是非常困难甚至不可能的，利用计算机研究相应问题的数值解法是十分必要的。编程实现从偏微分方程数值求解全过程需要很好的理论基础和编程技巧，而偏微分方程工具箱提供了研究和求解空间二维偏微分方程问题的一个强大而又灵活实用的环境借助于这个工具，我们可以从繁琐。共性的求解步骤中解脱出来而专注于问题的核心即问题的描述，定义及简化，边界条件的确定，求解方法和精度控制的选择等，大大提高了求解效率。参考文献1 孙良数值分析方法 .东北大学学报M. 北京：自然科学版，1998：16-21.2 德淑敏.求解任意多边形区域二维偏微分方程.中国农业大学学报2007:89-102.3 顾正平.周边式二次沉淀池流态的计算模拟研究环境科学研究19954 德淑敏,黄成,梅树立.中国农业大学学报M.万方数据.2007致 谢在论文完成之际，我们要特别感谢我的指导老师林爽老师的热情关怀和悉心指导。在我撰写论文的过程中，林爽老师倾注了大量的心血和汗水，无论是在论文的选题、构思和资料的收集方面，还是在论文的研究方法以及成文定稿方面，我都得到了林爽老师悉心细致的教诲和无私的帮助，特别是他广博的学识、深厚的学术素养、严谨的治学精神和一丝不苟的工作作风使我终生受益，在此表示真诚地感谢和深深的谢意。 在论文的写作过程中，也得到了许多同学的宝贵建议，同时还到许多在工作过程中许多同学的支持和帮助，在此一并致以诚挚的谢意。感谢所有关心、支持、帮助过我的良师益友。34