−
∇
2
u
(
x
)
=
f
(
x
)
,
x
in
Ω
,
(
1
)
-\nabla^{2} u(\boldsymbol{x})=f(\boldsymbol{x}), \quad \boldsymbol{x} \text { in } \Omega,\qquad(1)
−∇2u(x)=f(x),x in Ω,(1)
u
(
x
)
=
u
D
(
x
)
,
x
on
∂
Ω
,
(
2
)
u(\boldsymbol{x})=u_{\mathrm{D}}(\boldsymbol{x}), \quad \boldsymbol{x} \text { on } \partial \Omega,\qquad(2)
u(x)=uD(x),x on ∂Ω,(2)
其中,
u
=
u
(
x
)
u=u(\boldsymbol{x})
u=u(x) 为未知函数,
f
=
f
(
x
)
f=f(\boldsymbol{x})
f=f(x) 为规定函数,
∇
2
\nabla^{2}
∇2 为拉普拉斯算子(常写为
∇
\nabla
∇),
Ω
\Omega
Ω 为空间域,
∂
Ω
\partial \Omega
∂Ω 为
Ω
\Omega
Ω的边界。 泊松问题,包括偏微分方程
−
∇
2
u
=
f
-\nabla^{2} u=f
−∇2u=f 和
∂
Ω
\partial \Omega
∂Ω 上的边界条件
u
=
u
D
u=u_{\mathrm{D}}
u=uD,是边界值问题的一个例子,在开始使用 FEniCS 解决它之前必须精确说明。
在坐标为 x 和 y 的二维空间中,我们可以写出泊松方程为
−
∂
2
u
∂
x
2
−
∂
2
u
∂
y
2
=
f
(
x
,
y
)
,
(
3
)
-\frac{\partial^{2} u}{\partial x^{2}}-\frac{\partial^{2} u}{\partial y^{2}}=f(x, y),\qquad(3)
−∂x2∂2u−∂y2∂2u=f(x,y),(3)
未知数
u
\boldsymbol{u}
u 现在是两个变量的函数,
u
=
u
(
x
,
y
)
u=u(x, y)
u=u(x,y),在二维域
Ω
\Omega
Ω 上定义。
泊松方程出现在许多物理环境中,包括热传导、静电、物质扩散、弹性杆的扭曲、无粘性流体流动和水波。 此外,该方程出现在更复杂的偏微分方程系统的数值分裂策略中,尤其是 Navier-Stokes 方程。
求解边界值问题包括以下步骤:
有限元变分公式
FEniCS 基于有限元方法,它是 PDE 数值解的通用且高效的数学机器。 有限元方法的起点是以变分形式表示的偏微分方程。
将 PDE 转化为变分问题的基本方法是将 PDE 乘以函数
v
v
v,在域
Ω
\Omega
Ω 上对所得方程进行积分,并通过具有二阶导数的部分项进行积分。 乘以 PDE 的函数
v
v
v 称为测试函数。 要逼近的未知函数
u
u
u 称为试验函数。 术语试验和测试功能也用于 FEniCS 程序。 试验和测试函数属于某些所谓的函数空间,这些函数空间指定了函数的属性
在本例中,我们首先将泊松方程乘以测试函数
v
v
v 并在
Ω
\Omega
Ω 上积分:
−
∫
Ω
(
∇
2
u
)
v
d
x
=
∫
Ω
f
v
d
x
,
(
4
)
-\int_{\Omega}\left(\nabla^{2} u\right) v \mathrm{~d} x=\int_{\Omega} f v \mathrm{~d} x,\qquad(4)
−∫Ω(∇2u)v dx=∫Ωfv dx,(4)
我们在这里让 dx 表示域 Ω 上积分的微分元素。 我们稍后将让 ds 表示在 Ω 边界上积分的微分元素。
当我们推导出变分公式时,一个常见的规则是我们尽量保持
u
u
u 和
v
v
v 的导数的阶数尽可能小。 在这里,我们有
u
u
u 的二阶空间导数,可以通过应用部分积分技术将其转换为
u
u
u 和
v
v
v 的一阶导数。 公式是这样写的
−
∫
Ω
(
∇
2
u
)
v
d
x
=
∫
Ω
∇
u
⋅
∇
v
d
x
−
∫
∂
Ω
∂
u
∂
n
v
d
s
,
(
5
)
-\int_{\Omega}\left(\nabla^{2} u\right) v \mathrm{~d} x=\int_{\Omega} \nabla u \cdot \nabla v \mathrm{~d} x-\int_{\partial \Omega} \frac{\partial u}{\partial n} v \mathrm{~d} s,\qquad(5)
−∫Ω(∇2u)v dx=∫Ω∇u⋅∇v dx−∫∂Ω∂n∂uv ds,(5)
其中
∂
u
∂
n
=
∇
u
⋅
n
\frac{\partial u}{\partial n}=\nabla u \cdot n
∂n∂u=∇u⋅n 是
u
u
u 在边界外法线方向
n
n
n 上的导数。
变分公式的另一个特点是测试函数
v
v
v 需要在解
u
u
u 已知的边界部分消失。 在当前问题中,这意味着整个边界
∂
Ω
\partial \Omega
∂Ω 上的
v
=
0
v=0
v=0。 因此,等式(5)右边的第二项消失了。 从等式(4)和(5)可以得出
∫
Ω
∇
u
⋅
∇
v
d
x
=
∫
Ω
f
v
d
x
,
(
6
)
\int_{\Omega} \nabla u \cdot \nabla v \mathrm{~d} x=\int_{\Omega} f v \mathrm{~d} x,\qquad(6)
∫Ω∇u⋅∇v dx=∫Ωfv dx,(6)
抽象有限元变分公式
事实证明,为变分问题引入以下规范符号是很方便的:找到
u
∈
V
u \in V
u∈V 使得
a
(
u
,
v
)
=
L
(
v
)
∀
v
∈
V
^
,
(
9
)
a(u, v)=L(v) \quad \forall v \in \hat{V},\qquad(9)
a(u,v)=L(v)∀v∈V^,(9)
对于泊松方程,我们有:
a
(
u
,
v
)
=
∫
Ω
∇
u
⋅
∇
v
d
x
,
(
10
)
a(u, v)=\int_{\Omega} \nabla u \cdot \nabla v \mathrm{~d} x,\qquad(10)
a(u,v)=∫Ω∇u⋅∇v dx,(10)
L
(
v
)
=
∫
Ω
f
v
d
x
,
(
11
)
L(v)=\int_{\Omega} f v \mathrm{~d} x,\qquad(11)
L(v)=∫Ωfv dx,(11)
实现代码(Python)
from fenics import *
# Create mesh and define function space
mesh = UnitSquareMesh(8, 8)
V = FunctionSpace(mesh, ’P’, 1)
# Define boundary condition
u_D = Expression(’1 + x[0]*x[0] + 2*x[1]*x[1]’, degree=2)
def boundary(x, on_boundary):
return on_boundary
bc = DirichletBC(V, u_D, boundary)
# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
f = Constant(-6.0)
a = dot(grad(u), grad(v))*dx
L = f*v*dx
# Compute solution
u = Function(V)
solve(a == L, u, bc)
...
详情参阅 - 亚图跨际