偏微分方程的数值解(六)： 偏微分方程的 pdetool 解法

浅夏110

1 图形界面解法简介
对于一般的区域，任意边界条件的偏微分方程，我们可以利用 MATLAB 中 pdetool 提供的偏微分方程用户图形界面解法。 图形界面解法步骤大致上为：

  ~/ }0 C) V; ]6 C（1）定义 PDE 问题，包括二维空间范围，边界条件以及 PDE 系数等。

- J. a/ W( Q/ O（2）产生离散化之点，并将原 PDE 方程式离散化。
  \' \" d; |6 T" O* y
0 o6 C& E- G! i, w' K( K（3）利用有限元素法(finite element method；FEM)求解并显示答案。

( n3 j2 n4 Y0 A1 I在说明此解法工具之前，先介绍此 PDE 图形界面的菜单下方的功能图标(icon)按 钮。
, o0 ]# e! C: J( h9 y0 v- L, d( o
( G/ ?' \4 r+ c/ a' n$ g: t

2 图形界面解法的使用步骤
要利用 pdetool 接口求解之前，需先定义 PDE 问题，其包含三大部份：
. Y( _2 Q4 B( n0 V6 Z0 w1 Q4 M+ W
) D+ J& M4 E( i# Q+ @（1）利用绘图(draw)模式，定义 需要求解的问题的空间范围(domain)Ω 。

' [7 j: ]6 g; Y' S0 j' a（2）利用 boundary 模式，指定边界条件。

（3）利用 PDE 模式，指定 PDE 系数，即输入 c，a，f 和 d 等 PDE 模式中的系数。

9 [5 o2 K5 l6 r! z在定义 PDE 问题之后，可依以下两个步骤求解
4 y& P4 E! a6 \) P1 k* @
& Y# A# z- |* }) Q  J4 B（1）在 mesh 模式下，产生 mesh 点，以便将原问题离散化。

（2）在 solve 模式下，求解。

（3）最后，在 Plot 模式下，显示答案。


4 `* V; r' j  t& U2 E


2 [% x" k' M6 Y  n! F/ T3 K( [
& Z( v  h  a5 Y


; I/ f, ~$ x+ W0 l8 a4 r4 ?
5 }' a9 t' z1 v; O+ H! j' I* a( s( ^注意：

1. MATLAB 会以图形的方式展示结果，使用者亦可点选 plot 下之“parameters”功 能，选择适当的方式显示图形及数据。例如用 3D 方式显示求解结果。参数设置见图 10， 显示结果见图 11。


) H& T* R1 O  S% y3 M8 T( ^
9 _5 e) P3 l0 N  `2. 另外，若使用者欲将结果输出到命令窗口中，以供后续处理，可利用 solve 功能 项下之“export solution”指定变量名称来完成。
: e0 ~  o; i0 k( P
3. 如果求抛物型或双曲型方程的数值解，还需要通过“solve”菜单下的 “parameters…”选项设置初值条件。
5 Y) d4 H* b% z- M, i. f& W7 v, O7 L
+ g$ t2 i5 L2 J/ l! T4. 在上面定义边界条件和初始条件时，可以使用一些内置变量。

（1）在边界条件输入框中，可以使用如下变量： 二维坐标 x 和 y，边界线段长度参数（s s 是以箭头的方向沿边界线段从 0 增加到 1）， 外法向矢量的分量 nx 和 ny（如果需要边界的切线方向，可以通过 tx=-ny 和 ty=nx 表示）， 解 u。
" K0 Y* q0 O% m  N/ |9 _
（2）在初值条件的输入框中，也可以输入用户定义的 MATLAB 可接受变量（p， e，t，x，y）的函数。

例 11 使用 PDETOOL 重新求例 8 的数值解。

6 d4 ?9 T6 v3 ~/ n! |- {/ O' V                例8  求解正方形区域{(x, y) | −1 ≤ x, y ≤ 1}上的热传导方程



, ^$ M6 r2 W, Y8 w6 i边界条件为Dirichlet条件u = 0。

6 L8 ~/ T4 g4 h/ S8 p" y! B解 这里是抛物型方程，其中c = 1, a = 0, f = 0, d = 1。

- Y3 R6 y3 F% W3 ?1）定义 PDE 问题，包括二维空间范围，边界条件以及 PDE 系数等。我们这里就 省略了。

+ m6 {$ Q9 O8 ^3 g# e  L7 {2）区域剖分以后，通过“Mesh”菜单下的“Export Mesh…”选项可以把 p，e，t 三个参数分别输出到工作间。
9 H" g) z- d4 j
3）然后编写函数 fun1(x,y)如下：
0 M7 ^' E) R5 i/ x; u
function f=fun1(x,y);
f=zeros(length(x),1);
& o. Y& P+ ~& X6 Gix=find(x.^2+y.^2<0.16);
3 P1 S. L* h2 K' z, mf(ix)=1;
& ~, s7 \; t7 W+ t: A
其中的变量 x，y 是 MATLAB 可接受的内置变量。 设置“solve”菜单下的“parameters…”选项如下：

时间框中输入：linspace(0,0.1,20)；
8 P: D$ a$ m6 a0 V' l  L0 N0 R/ F
- J5 B% g5 w/ x6 I9 N/ c% `) _$ A初值框中输入：fun1。

4）设置“plot"菜单下的“parameters…”选项如下：选择 Height(3-D plot) 和 Animation 两项。

5）用鼠标点一下工具栏上的“＝”按钮，就可以画出数值解的 3-D 图形。
& F, y2 h- [9 R% v) Y

————————————————
版权声明：本文为CSDN博主「wamg潇潇」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
原文链接：https://blog.csdn.net/qq_29831163/article/details/89712663