圆形区域上的 Poisson 方程

以 Poisson 方程为例，我们一起来看一下如何使用 FEtch 系统快速地求解偏微分方程问题。

问题描述

假设在二维区域 $\Omega$ 上函数 $u(x,y)$ 满足如下条件： $-\left(\frac{\partial^{2} u}{\partial x^{2}}+\frac{\partial^{2} u}{\partial y^{2}}\right)=xy\ \left(\text{in } \Omega\right)$

\[u(x, y) = 0 \ \left(\text{on } \partial \Omega\right)\]

其中，$\Omega$ 是位于坐标原点的一个单位圆盘， $\partial \Omega$ 是 $\Omega$ 的边界。求函数 $u(x,y)$ 的有限元解。

脚本文件

这是一个比较简单的单物理场问题。用有限元语言表述，只需要填写 4 种类型的文件，即 PDE、NFE、GCN 和 MDI 文件。

PDE 文件 (poisson.pde)，用于描述微分方程的弱形式。在有限元理论框架中，对应于单元矩阵的计算，是整个程序中最核心的部分。

defi
disp u
coor x y
shap %1 %2
gaus %3
mate ek 1.0

stif
dist = +[u/x;u/x]*ek+[u/y;u/y]*ek

load = +[u]*x*y

end

NFE 文件 (psa.nfe)，与 PDE 文件相配合，用于给出待求解的矩阵格式的线性方程组（对应于整体矩阵组装），以及对计算结果作进一步处理。

defi
stif s
load f

equation
matrix = [s]
forc = [f]

solution u
gidpost("u",1.0) u(1)

end

GCN 文件 (ps.gcn)，用于定义求解器类型。同时与 NFE 文件相配合，一起组织计算流程。

defi
a ell

#startsin a
#solvsin a

MDI 文件 (ps.mdi)，用于组织管理物理场信息。

2dxy
#a 0 1 u
pde poisson q4g2
#

实际操作中，这里的 NFE 文件可以调用系统的库文件自动生成，因此用户只需填写 PDE、GCN 和 MDI 三个文件即可，统计下来，总行数不超过 20 行！

我们通过 FEtch 客户端提交脚本文件后，会自动生成并下载有限元程序文件（ps.exe 和 run.bat）、材料参数文件（ps.mat）和 GiD 软件的前处理模板文件（ps.cnd 和 ps.bas）。编程任务至此已经顺利完成。

前处理

前处理和后处理都是通过第三方软件 GiD 来完成的。

我们首先进入 GiD 软件界面进行前处理。依次进行几何建模、施加边界条件、赋材料值和网格剖分，最终导出计算所需的数据输入文件（ps.dat)。

边界条件的施加情况如下

imag

使用的网格如下

imag

计算结果

经过计算，生成了后处理文件（ps.flavia.res）。我们再次进入 GiD 软件界面进行后处理，查看云图效果如下。

imag

至此，Poisson 方程的求解已顺利完成。

参考文献

[1] Hecht F , Pironneau O , Ohtsuka K . FreeFem++ manual. 2011.