1、第一章 概 述11 偏微分方程工具箱的功能偏微分方程工具箱(PDE Toolbox)提供了研究和求解空间二维偏微分方程问题的一个强大而又灵活实用的环境。PDE Toolbox 的功能包括：(1) 设置 PDE (偏微分方程)定解问题，即设置二维定解区域、边界条件以及方程的形式和系数；(2) 用有限元法 (FEM) 求解 PDE 数值解；(3) 解的可视化。无论是高级研究人员还是初学者，在使用 PDE Too1box 时都会感到非常方便。只要 PDE 定解问题的提法正确，那么，启动MATLAB 后，在 MATLAB 工作空间的命令行中键人 pdetool，系统立即产生偏微分方程工具箱(PDE T
2、oolbox)的图形用户界面(Graphical User Interface，简记为 GUI)，即 PDE 解的图形环境，这时就可以在它上面画出定解区域、设置方程和边界条件、作网格剖分、求解、作图等工作，详见 14 节中的例子。我们将在第二章详细介绍GUI 的使用，在第二章给出大量典型例子和应用实例。除了用 GUI求解 PDE 外，也可以用 M 文件的编程计算更为复杂的问题，详见第三章和第四章的内容。12 PDE Toolbox 求解的问题及其背景121 方程类型PDE Toolbox 求解的基本方程有椭圆型方程、抛物型方程、双曲型方程、特征值方程、椭圆型方程组以及非线性椭圆型方程。椭圆型方
3、程： ，(), ,cuafin椭圆型方程： i其中 是平面有界区域，c,a,f 以及未知数 u 是定义在 上的实（或复）函数。抛物型方程： (), .udcaufint双曲型方程： .2, fit特征值方程： (), ,cuadin其中 d 是定义在 上的复函数， 是待求特征值。在抛物型方程和双曲型方程中，系数 c,a,f 和 d 可以依赖于时间 t。可以求解非线性椭圆型方程：()(), ,cufuin其中 c,a,f 可以是未知函数 u 的函数。还可以求解如下 PDE 方程组；1121212()(),ccaufu利用命令行可以求解高阶方程组。对于椭圆型方程，可以用自适应网格算法，还能与非线性
4、解结合起来使用。另外，对于 Poission 方程还有一个矩形网格的快速求解器。122 边界条件（1）Dirichlet 条件 ： hur( 2 ) Neumann 条件： ()ncuqg其中 是 的边界 上的单位外法向量， 和 是定义在n,hr上的函数。对于特征值问题仅限于齐次条件： 和 。对于 0,非线性情形系数 和 可以依赖于 u；对于抛物型方程和双曲,gqhr型方程，系数可以依赖于时间 t。对于方程组情形，边界条件为( 1 ) Dirichlet 条件： 121hur212hur( 2 ) Neumann 条件： 1()()nccqg21222( 3 ) 混合边界条件为: 1hur12
5、121()()nccuqugh2 2其中 的计算要使得 Dirichlet 条件满足。在有限元法中，Dirichlet条件也称为本质边界条件，Neumann 条件称为自然边界条件。1.3 如何使用 FDE Toolbox1.3.1 定解问题的设置员简单的办法是在 PDE Tool 上直接使用图形用户界面(GUl)。设置定解问题包括三个步骤：(1)Draw 模式：使用 CSG(几何结构实体模型 )对话框画几何区域，包括矩形、圆、椭圆和多边形，也可以将它们组合使用。(2)Boundary 模式：在各个边界段上给出边界条件，(3)PDE 模式：确定方程的类型、系数 c，a，f 和 d c。也能够在不
6、同子区域上设置不同的系数(反映材料的性质) 。1.3.2 解 PDE 问题用 GUI 解 PDE 问题主要经过下面两个过程（模式）(1)Mesh 模式；生成网格自动控制网格参数。(2)Solve 模式：对于椭圆型方程还能求非线性和自适应解。对于抛物型和双曲型力程设置初始边值条件后能求出给定 t 时刻的解。对于特征值问题，能求出给定区间内的特征值；求解后可以加密网格再求解。1.3.3 使用 Toolbox 求解非标准的问题对于非标准的问题。可以用 PDE Too1box 的函数。或者用FEM(有限元法 )求解更为复杂的问题。1.3.4 计算结果的可视化从 GUI 能够使用 Plot 模式实现可视
7、化。可以使用 Color, Height和 Vector 等作图。对于抛物型和双曲型方程，还可以生成解的动画。这些操作通过命令行都很容易实现。1.3.5 应用领域在应用界面提供了丁如下应用领域结构力学平面应力问题结构力学平面应变问题静电场问题静磁场问题交流电磁场问题直流导体介质问题热传导问题9散问题这些界面都有对话框，它包括 PDE 的系数、边界条件、解的性质等。1. 4 解偏微分方程的一个例子解 Poisson 方程 ，边界条件为齐次 Dirichlet 类型。uf第一步：启动 MATLABl, 键入 pdetool，按回车键确定便可启动 GUI，然后在 Options 菜单下选择 Grid
8、 命令，打开栅格, 栅格的使用，能使用户容易确定所绘图形的大小，如图 111-1第二步：分步完成平面几何造型：R1-C1-E1+R2+C2 。用菜单或快捷工具，分别画矩形 R1、矩形 R2、椭圆 E1、圆 C1、圆 C2。画圆时，首先选中椭圆工具，按鼠标右键并拖动即可、或者在按 ctrI的同时，拖动鼠标也可绘制圆。然后在 Set formula 栏，进行编辑并用算术运算将将图形对象名称连接起来，删除默认的表达式键入R1-C1-E1+R2+C2，按等号健得到所需图形。若需要，还可进行储存形成 M 文件。选择 Boundary 菜单中 Boundary Mode 命令，进入边界模式。单击 Boun
9、dary 菜单中 Remove A11 Subdomain Borders 选项，去除子域边界。如果想将几何信息和边界信息进行存储，应选择 Boundary 菜单中的 ExPort Decomposed GeometryBoundary Conds命令，将它们分别储存于 g,b 变量中, 通过 MATLAB 形成 M 文件。第三步：选取边界单击 Boundary 菜单中 Specify Bounddy Conditions选项，打开 Boundary conditlons 对话框，输入边界条件，如图 14。本例取缺省条件。即将全部边界设为齐次 Dirichlet 条件，边界颜色显示为红色。第四
10、步：选择 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,
11、Height (3D Plot)和 Show mesh 三项。然后单击 Plot 按钮，显示三维图形解。第九步：如果要画等值线图和矢量场图，单击 Plot 菜单中Parameters选项，打开 Plot Selection 对话框选中 Contour 和Arrows 两项。然后单击 P1ot 按钮，可显示解的等值线图和矢量场图。第二章 PDE 图形用户界面2.1 PDE Toolbox 菜单File 菜单（如图 1-1）图 1-1New 新建一个几何结构实体模型（Constructive Solid Geomery,简记为 CSG） ，默认文件名为“Untitled” 。Open 从硬盘装载
12、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 菜单（如图
13、 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/squa
14、re（centered） 以中心方式画矩形/方行（Ctrl+鼠标） 。Ellipse/circle 以矩阵角点方式画椭圆/ 圆（Ctrl+ 鼠标） 。Ellipse/circle（centered ） 以中心方式画椭圆/圆（Ctrl+鼠标） 。Polygon 画多边形，单击鼠标右键可封闭多边形。Rotate 旋转已选的图形。Export Geometry Description,Set Formula，Labels将几何描述矩阵 gd、公式设置字符 sf和标识空间矩阵 ns 输出到主工作空间去。单击 Draw 菜单中 Rotate选项，可打开 Rotate 比对活框，通过输入旋转的角度，可使选
15、择的物体按输入的角度逆时针旋转。旋转中心的选择如果缺省，则为图形的质心，也可以输入旋转中心坐标。5 Boundary 菜单（如图 1-5）图 1-5Boundary Mode 进入边界模式。Specify Boundary Conditions 对于已选的边界输入条件，如果没有选择边界，则边界条件适用于所有的边界。Show Edge Labels 显示边界区域标识开关，其数据是分解几何矩阵的列数。Show Subdomain Labels 显示子区域标识开关，其数据是分解几何矩阵中的子域数值。Remove Subdomain Border 当图形进行布尔运算时，删除已选取的子域边界。Remov
16、e All Subdomain Borders 当图形进行布尔运算时，删除所有的子域边界。Export Decomposed Geometry，Boundary Conds将分解几何矩阵 g、边界条件矩阵 b 输出到主工作空间。选择 Boundary 菜单中 Specify Boundary Conditions命令可定义边界条件。在打开的 Boundary condition 对话框，可对已选的边界输入边界条件。共有如下三种不同的条件类型：NeMmann 条件 这里边界条件是由方程系数 q 和 g 确定的，在方程组的情况下(换成方程组模式)，q 是 22 矩阵，g 是 2x1 矢量。Diri
17、chlet 条件 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 Coeffici
18、ents 将当前 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
19、 Labels 显示三角形网格标识开关，三角形网格标识数据是三角形矩阵 t 的列。Parameters 修改网格生成参数。Export Mesh 输出节点矩阵 p、边界矩阵 e 和三角形矩阵 t到主工作空间。8 Solve 菜单（如图 1-8）图 1-8Solve PDE 对当前的几何结构实体 CSG、三角形网格和图形解偏微分方程。Parameters 调整 PDE 的参数。Export Solution 输出 PDE 的解矢量 u。如果可行，将计算的特征值 1 输出到主工作空间。1.1.9 Plot 菜单（如图 1-9）图 1-9Plot Solution 显示图形解。Parameters
20、打开绘图方式对话框。Export Movie 如果动画被录制了，则动画矩阵 M 将输出到主工作空间。10 Window 菜单从 Window 菜单项，可选择当前打开的所有的 MATLAB 图形窗口，被选择的窗台移至前台。11 Help 菜单Help 显示帮助信息About 显示版本信息12 PDE 工具栏主菜单下是工具栏，工具栏中喊有许多工具图标按钮，可提供快速、便捷的操作方式。从左到右 5 个按钮为绘图模式按钮，紧接着的 6 个为边界、网格、解方程和图形显示控制功能按钮，最右边的为图形缩放功能键。 （如图 1-10）图 1-10以角点方式画矩形/方行（Ctrl+鼠标） 。以中心方式画矩形/方
21、行（Ctrl+鼠标） 。以矩形角点长轴方式画椭圆/圆（Ctrl+鼠标） 。以中心方式画椭圆/圆（Ctrl+鼠标） 。画多边形，按右键可封闭多边形。进入边界模式。打开 PDE Specification（偏微分方程类型）对话框。初始化三角形网格。加密三角形网格。解偏微分方程。打开 Plot Selection 对话框，确定后给出解的三维图形。为显示缩放切换按钮。第三章 典型方程及其应用实例求解 PDE 问题主要有两种方法，一种是使用图形用户界面，另一种是采用命令行编程。前者直观简便，而后者更为灵活。21 求解椭圆方程的例子例：单位圆上的 Poisson 方程边值问题：21,|1|0uxy这一问题
22、的精确解为： 21,.4xyuy若使用图形用户界面（Graphical User Interface,简记为 GUI） ，则首先在 MATLAB 的工作窗口中键入 pdetool，按回车键确定，于是出现 PDE Toolbox 窗口。如果需要坐标网格，单击 Options 菜单下的 Grid 选项即可。下面分步进行操作。（i）画区域圆 单击工具 ，大致在（0，0）位置单击鼠标右键同时拖拉鼠标到适当位置松开，绘制圆。为了保证所绘制的圆是标准的单位圆，在所绘圆上双击，打开 Object Dialog 对话框，精确地输入圆心坐标 X-center 为 0、Y-cebter 为 0 及半径 Radiu
23、s 为 1，然后单击 OK 按钮，这样单位远已画好。（ii）设置边界条件 单击工具 ，图形边界变红，逐段双击边界，打开 Boundary Condition 对话框，输入边界条件。对于同一类型的边界，可以按 Shift 键，将多个边界同时选择，统一设置边界条件。本题选择 Dirichlet 条件，输入 h 为 1，r 为 0，然后单击OK 按钮。也可以单击 Boundary 菜单中 Specify Boundary Conditions选项，打开 Boundary Condition 对话框输入边界条件，如图 2-1。（iii ）设置方程 单击 PDE 菜单中 PDE Specificatio
24、n选项，打开 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 个，如
25、此还可继续加密。（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
26、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 选项，再单击工具 或 ，即可显示节点编号。若要输出节点坐标，只需单
27、击 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，单击 O
28、K 按钮确定。再在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(
29、p,t,u)Pdesurf(p,t,u-exact)通过命令行键入 help+命令函数，如 help pdemesh，按回车键，可以调入有关命令函数的定义、参数格式等帮助信息。22 求解抛物型方程的例子例：考虑一个带有矩形孔的金属板上的热传导问题。板的左边保持在 100 ，板的右边热量从板向环境空气定常流动，其他边及内孔c边界保持绝缘。初始 时板的温度为 0，于是概括为如下定解问题：0t00,1,|.tudtnu域 的外边界顶点坐标为（-0.5，-0.8） ， （0.5，-0.8） ， （0.5，0.8） ，（-0.5 ，0.8） 。内边界顶点坐标为（-0.005，-0.4） ， （0.05，
30、-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-
31、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）网格剖分 单击 ，或者加密网格，单击
32、 。（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选项，打开 Plo
33、t 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)23 求解双曲型方程的例子例：考虑如下二维波动方程的定解问题： ).2sin(12i3,arcto,0,0| ,1,|),(),(0xyxetutuyxut 用 GUI 求解。类
34、似前面的例子，首先作正方形区域：设置断点坐标为（-1 ，1） ， （-1，-1 ） ， （1，-1） ， （1， 1） 。在 Object Dialog 对话框中输入 Left 为-1 ， Bottom 为-1，Width 为 2，Height2，单击 OK按钮。设置边界条件：左、右边界用 Dirichlet 条件，输入 h 为 1，r为 0；上、下边界用 Neumann 条件，输入 q 为 0，g 为 0，作网格剖分。设置方程类型为 Hyperbolic（双曲型） ，键入c=1， a=0， f=0，d=1。单击 Solve 菜单中 Parameters选项，打开 Solve Paramete
35、rs 对话框，在 Time 栏中键入 linspace（0，5，31） ，设置 u 的初始值u（t0 ）为 atan(cos(pi/2*x),u的初始值 u(t0)为 3*sin(pi*x).*exp(sin(pi/2*y),Relative tolerance 为 0.01，Absolute tolerance 为0.001。最后输出图形解：单击 Plot 菜单中 Parameters选项，打开 Plot Selection 对话框，选 Color，Contour，Arrows，Height（3-D Plot）和 Show mesh 五项，单击 Plot 按钮，三维彩色图形的解如图 2-7
36、所示。图 2-7MATLAB 的动画程序：p,e,t=initmesh(squareg);x=p(1,:);y=p(2,:);u0=atan(cos(pi/2*x);ut0=3*sin(pi*x).*exp(sin(pi/2*y);n=3l;tlist=linspace(0,5,n);uu=hyperbolic(u0,ut0,tlist,squareb3,p,e,t,1,0,0,1);delta=-1:0.1:1;uxy,tn,a2,a3=tri2grid(p,t,uu(:,1),delta,delta);gp=tn;a2;a3;umax=max(max(uu);umin=min(min(uu
37、);newplot M=moviein(n);For I=1:n,pdeplot(p,e,t,xydata,uu(:,I),zdata,uu(:,I), mesh,off,xygrid, on, gridparam,gp,colorbar,off,zstyle, continuous);axis(-1 1 1 1 umin umax);caxis(umin umax);M(:,I)=getframe;EndMovie(M,10);24 求解特征值问题的例子例：混合边界条件的特征值问题：,0|43|,1|)(),(11xyxunyxu首先作区域：方形的 4 个顶点为（-1，1） ， （-1 ，-
38、1） ， （1，-1 ） ，（1，1） 。在 Object Dialog 对话框中输入 Left 为-1 ，Bottom 为-1，Width 为 2，Height2。单击 ，逐段双击边界设置边界条件：左边界选 Dirichlet 条件，输入 h 为 1，r 为 0；右边界选 Neumann 条件，输入 q 为-3/4，g 为0；上、下边界选 Neumann 条件，输入 q 为 0，g 为 0。在 PDE Specification 对话框中设置方程类型为 Eigenmodes（特征值模式） ，PDE 的系数设置为 c=1，a=0， d=1。由于右边界上 Neumann 条件，特征值可以出现负值
39、，因此在求小于 10 的特征值时应在 Solve Parameters 对话框中键入-inf 10。作网格剖分，然后单击 ，得到第一个特征值对应的特征函数图形，如图 2-8。图 2-8再在 Plot Selection 对话框的 Eigenvalue 项中选取不同的特征值，比如第二个特征值“2.06，并选 Color，Contour，Height（3-D Plot）和 Show mesh 四项，如图 2-9，可作出第二个特征值对应的特征函数的三维图形，如图 2-10。图 2-9图 2-10MATLAB 程序：p,e,t=initmesh(squareg);p,e,t=refinemesh（sq
40、uareg,p,e,t）;v,l=pdeeig(squareb2,p,e,t,1,0,1,-inf 10);pdesurf(p,t,v(:,4)25 偏微分方程在一维空间的应用已知偏微分方程初始条件时间 t 和一维空间变量 x，MATLAB自身存在对偏微分方程的解函数 pdepe。单一方程假设我们要求解热传导方程 uxt，0,）（ 1,）（ txu2,）（1.1）MATLAB 指定抛物型的偏微分方程形式如下： ),(),(),( uxuxxmtx tsutbtc ）（边界条件为： ),(),),( uxxxlll tbtqutp0rrr其中 为边界左端点， 为边界右端点，初始条件为 xl r。
41、 （我们注意到同一函数 b 在方程和边界条件中出现。 ）)(0(fu，通常情况下，每个函数都会被指定到不同的 MATLAB 文档中。也就是说，与方程有关的函数 c,b 和 s 都会被指定保存到一个MATLAB 文件中，与边界条件有关的函数 p 和 q 保存到令一个文档中（此外，我们需要注意到函数 b 是相同的只需保存一次即可） ，最后将初始函数 保存在第三个文档中。执行命令函数 pdepe 将把)(xfm 个文档结合起来进行运算，返回方程的解。在我们的例子中，有：1),(uxtcb0),(xts我们将它们保存到 MATLAB 文档中，记作 eqn1.m (m=0 以后在加详述)。Functio
42、nc,b,s=eqn1(x,t,u,DuDx)%EQN1:MATLAB function M-file that specifies%a PDE in time and one space dimension.C=1;B=DuDx;S=0;对于边界条件，我们有：P(0,t,u) = u ； q(0,t) = 0P(1,t,u) = u 1 ； q(1,t) = 0将它们保存到 MATLAB 文档中，记作 bc1.m.Function p1,q1,pr,qr=bc1(x1,u1,xr,ur,t)%BC1:MATLAB function M-file that specifies boundary conditions%for a PDE in time and one space dimension.P1=u1;Q1=0;Pr=ur-1;Qr=0;对初始条件，有：xf21）（将它保存到 MATLAB 文档中，记作 initial1.m.Function value=initial1(x)%INITIAL1:MATLAB function M-file that specifies the initial condition%for a PDE in time and one space dimension.Value=2*x/(1+x2);
