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

浅夏110

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

（1）定义 PDE 问题，包括二维空间范围，边界条件以及 PDE 系数等。

（2）产生离散化之点，并将原 PDE 方程式离散化。

, M8 m) v& X9 X% s（3）利用有限元素法(finite element method；FEM)求解并显示答案。
0 I+ p9 g( _4 p3 I7 x! b+ ]
! v6 A* d$ E$ Z* B( H; U4 X- p; C1 V, c在说明此解法工具之前，先介绍此 PDE 图形界面的菜单下方的功能图标(icon)按 钮。
5 _. {9 i; y9 w# f
, |* P7 Q) k, ?# r2 P

5 `, K& H# b# p1 H2 图形界面解法的使用步骤
, }' |) U/ _2 D# p  V要利用 pdetool 接口求解之前，需先定义 PDE 问题，其包含三大部份：

0 Q+ [% h& A9 a; V, s  D) y（1）利用绘图(draw)模式，定义 需要求解的问题的空间范围(domain)Ω 。

0 S/ q* |8 |2 g: d9 G9 Z- T0 d（2）利用 boundary 模式，指定边界条件。

（3）利用 PDE 模式，指定 PDE 系数，即输入 c，a，f 和 d 等 PDE 模式中的系数。
5 o1 C7 W# ?; l$ R7 x. f
在定义 PDE 问题之后，可依以下两个步骤求解
/ q, o0 D8 [5 _
7 J4 A8 D/ J0 K5 z（1）在 mesh 模式下，产生 mesh 点，以便将原问题离散化。
- r( M( V% [0 o6 u& _) |; A
0 q1 x, F( ^! g（2）在 solve 模式下，求解。
0 H; Y0 Q5 |8 q& _( N) j( Z' U) m/ n
) }2 m9 K& ^; i8 [（3）最后，在 Plot 模式下，显示答案。
( \4 P7 b' r* b
6 T8 \' }; c' X" ]3 z1 K* z1 ?; z
3 b4 r& m' ?/ A5 A6 q; M
2 q! H$ E* Y1 e  N  b
' N1 h! {; n( B, O: o

) v+ f; l; z! N
7 H7 c1 Z9 f; s, K8 T7 B) l
  ~2 e' K' P! @& f

注意：

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

3 h* t# h& w  f5 @( Y9 t; w

2. 另外，若使用者欲将结果输出到命令窗口中，以供后续处理，可利用 solve 功能 项下之“export solution”指定变量名称来完成。
% ^. R& _4 f. L( p4 M; P: L* e: H
3. 如果求抛物型或双曲型方程的数值解，还需要通过“solve”菜单下的 “parameters…”选项设置初值条件。
# n+ Q- u3 `+ |$ s. j! R" m# Y
4. 在上面定义边界条件和初始条件时，可以使用一些内置变量。
9 R' g; o9 [- {* \& j
  f& \5 j8 T. [（1）在边界条件输入框中，可以使用如下变量： 二维坐标 x 和 y，边界线段长度参数（s s 是以箭头的方向沿边界线段从 0 增加到 1）， 外法向矢量的分量 nx 和 ny（如果需要边界的切线方向，可以通过 tx=-ny 和 ty=nx 表示）， 解 u。

（2）在初值条件的输入框中，也可以输入用户定义的 MATLAB 可接受变量（p， e，t，x，y）的函数。
$ M. |" p6 M; ?6 ?
例 11 使用 PDETOOL 重新求例 8 的数值解。

                例8  求解正方形区域{(x, y) | −1 ≤ x, y ≤ 1}上的热传导方程
4 Q2 Q9 v1 H/ Y6 `' d, S; y( a
9 }; w- \0 _5 B% A

/ r4 @# {' g; Q6 c$ L, |: s边界条件为Dirichlet条件u = 0。

解 这里是抛物型方程，其中c = 1, a = 0, f = 0, d = 1。
% ?  G6 Q( f+ w2 i3 o0 V' ~
1）定义 PDE 问题，包括二维空间范围，边界条件以及 PDE 系数等。我们这里就 省略了。

. U0 z- u$ T) x, b2）区域剖分以后，通过“Mesh”菜单下的“Export Mesh…”选项可以把 p，e，t 三个参数分别输出到工作间。
* T/ J2 |- K" t6 G  ]" Y) x
3）然后编写函数 fun1(x,y)如下：
9 j- h" f' J3 C
- k; X, W. C6 a2 N9 ^function f=fun1(x,y);
% L9 ]% n$ M0 ?0 bf=zeros(length(x),1);
; D& n$ G: u6 X1 Z, xix=find(x.^2+y.^2<0.16);
. W* z9 {6 y3 U3 Q9 k, Uf(ix)=1;

% ^+ C% w1 c: ^7 d其中的变量 x，y 是 MATLAB 可接受的内置变量。 设置“solve”菜单下的“parameters…”选项如下：
) b$ z8 q0 q. M; _$ V, @
时间框中输入：linspace(0,0.1,20)；
5 Y. |! B3 x  s9 q8 y  @& J  Q7 g; X
! o3 j4 H" K* q) L; l% s/ Q初值框中输入：fun1。

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

/ `3 k3 `5 h0 V# U: [5）用鼠标点一下工具栏上的“＝”按钮，就可以画出数值解的 3-D 图形。


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

* N8 v2 `1 H- P/ F5 d. ?% }6 y
/ N# U& F: f$ p: p% N4 K