偏微分方程边值问题的数值解法.docx

1、偏微分方程边值问题的数值解法求解偏微分方程的边值问题本实验学习使用MATLAB的图形用户命令pdetool来求解偏微分方程的边值问 题。这个工具是用有限元方法来求解的，而且采用三角元。我们用内个例题来说明 它的用法。一、MATLAB支持的偏微分方程类型考虑平面有界区域D上的二阶椭圆型PDE边值问题：-VHcVw) + aw = f其中(1) V = (2)叮退D上的己知函数(3) c是标赞或2x2的函数方阵未知函数为“3) (x)eO它的边界条件分为三类：(1)Direchlet 条件：加=/(2)Neumann 条件：wC(cVw) + w = g(3)混合边界条件：在边界dD上部分为Dir

2、echlet条件，另外部分为Neumann条件。其中 是定义在边界dD 的已知函数，另外C 也可以是一个2*2的函数矩阵，n是沿边界的外法线的单位向量。在使用pdetool时要向它提供这些已知参数。二、例题例题1 用pdetool求解-Am = 1 D: x? +y2 1wLd = 0解：首先在MATLAB的工作命令行中键入pdetool ,按回牟键确定，于是出 现 PDE Toolbox 窗口，选 Genenic Scalar 模式.(1 )画区域圆单击椭圆工具按钮，大致在(0, 0)位置单击鼠标右键，拖拉鼠标到适当位 置松开。为了保证所绘制的圆是标准的单位园，在所绘园上双击，打开Objec

3、t Dialog对话框，精确地输入圆心坐标X-center为0、Y-center为0及半径 Radius为1 ,然后单击0K按钮，这样单位画已画好.(2)设置边界条件单击工具边界模式按钮，图形边界变红，逐段双击边界，打开Boundary condition对话框.输入边界条件.对于同一类型的边界，可以按Sh辻t键，将多 个边界同时选择，统一设边界条件.本题选择Dirichlet条件，输入h为1 , r 为0。，然后单击0K按钮.也可以单击Boundary菜单中Spoc辻y Boundary Condition选项，打开Boundary Condition对话框输入边界条件.(3)设冒方程单击偏

4、微分方程按钮，打开PDE Spec辻ication对话框，选择方程类型本题 选Ellintic (椭圆型)，输入c为1 , &为0 , f为1，然后单击0K按钮.(4)网格剖分单击网格工具，或者单击Mesh菜单中Initialize Mesh项，可进行初始网 格剖分.这时在PDE Toolbox窗口下方的状态栏内显示出初始网格的节点数和三 角形单元数.本题节点数为144个，三角形单元数为254个(图?？)。如果要细化网格，单击细化工具，或者单击Mesh菜单中Refine Mesh选项，节点数成为541 个，三角形单元数为1016个。(5)解方程单击解方程工具，或者单击S olve菜单中Solv

5、e PDE选项，可求得方程 数值解并用彩色图形显示。单击作图工具，或者单击Plot菜单中Parameter选 项，出现Plot selection对话框.从中选择于Height ( 3D plot),然后单击 Plot按钮，方程的图形解如图?？所示。除了作定解问题解U的图形外，也可以作等图形（6）输出网格节点的编号、单元编号以及节点坐标单击Mesh菜单中Show Node Labels选项，再单击网格匸具，即可显示节 点编号（图？ ？）。若要输出节点坐标，只需单击Mesh菜单中Export Mesh 选项，这时打开的Export对话框中的默认值为p e t,这里p、e、t分别表示 point

6、（点）、edges（边）、triangles（三角形）数据变量，单击0K按钮，然后 在MATLAB命令行键入p,即可以显示按节点编号排列的坐标；键入e再回车则显示 边界数据矩阵（7维数组）；键入t按回车则显示三角形单元数据矩阵（4维数 组）。点、边、单元的部分输出为:P 二Columns 1 through 11-1. 0000 0. 0000 1. 0000 0. 0000 一0. 7071 0. 70710. 7071-0. 7071-0. 9808-0. 9239 -0 8315-0. 0000 -1. 000001.0000 -o. 7071 -0. 7071 0. 70710. 70

7、71-0 1951-0. 3827-0. 5556e 二Columns 1 through111. 00009. 0000 10. 000011.00005. 000012.000013.0000 14.00002. 0000 15. 000016. 00009. 000010. 0000 11. 00005. 000012.000013. 000014.00002. 0000 15. 0000 16. 0000 17. 00000 0. 1250 0.2500 0.3750 0.50000. 62500. 7500 0. 875000. 12500. 25000. 12500. 25000.

8、 37500. 50000. 62500.75000. 87501. 00000. 12500. 25000. 37501. 00001. 00001. 00001. 00001. 00001. 00001. 00001. 00002. 00002. 00002. 00001. 00001. 00001. 00001. 00001. 00001. 00001. 00001. 00001. 00001. 00001.000000000000000t -Columns 1through 1832 1420262917231001166 891 99451213971 23 486 7285321

9、9105 121314289 9781988492991279489 119 119 95118 11890701261 11 111 11111 11 11111(7)输出近似解数值单击Solve菜单中Export Solution。选项，在Export对话框中输入u 再单击OK按钮，再在MATLAB命令行中输入u并回车，就会显示按节点编号排列的 解u的数值。(8)近似解和准确解的比较方程的准确解为：1x2 V2飒x,y) =厂为了与准确解比较，单击Plot菜单中Parameters选项，打开Plot Selection 对话框，在 Height (3D plot)行的 Property 下拉框中选 User Entry, 并且输入u- (1 -x.2-y./4 ,单击Plot按钮，就可以看到误差曲面，其数量级为104练习1 用pdetool 11具求解课本P128的第2、3、4这儿题的解，并作岀 图形。练习2 用pdetool工具求解课本P418的第2、4这两题的解（三角单单元 形状不限），并作出图形。