有限元法求解偏微分方程之FEniCS入门讲解（有限元法求解偏微分方程之FEniCS入门讲解）

【我已经发了同名西瓜视频，本文仅是视频中Markdown文件内容，不包括视频中讲解】

如果你希望和我一样Win10上安装FEniCS，建议采用Win10 WSL2 Ubuntu。

FEniCS是一个有限元法求解偏微分方程的开源计算平台。

• 泊松方程
• 将泊松方程转化成（弱）变分形式
• 求解微分方程边值问题（1维泊松方程）
• 求解2维泊松方程
• ParaView中查看解

泊松方程

第一个例子自然是求解著名的泊松方程。

将泊松方程转化成（弱）变分形式

首先：对泊松方程两边都乘上测试函数，然后对全域​积分。

这里需要说明一下试探函数和测试函数的概念。

测试函数，定义为：

试探函数，定义为：

然后，对上面这个变分方程进行变换

左边是双线性形式​，右边是线性形式​:

于是得到标准的变分问题：寻求​，满足

这就是和前面的泊松问题等价的变分问题。

求解微分方程边值问题（1维泊松方程）

有了前面的准备，现在可以写代码求解了。

第一步：为求解域创建网格，并定义函数空间

from dolfin import * # 单元个数 nel = 20 # 左右端点 xmin = 0 xmax = 1 # 试探/测试函数的多项式阶数 p = 2 # 创建网格 mesh = IntervalMesh(nel,xmin,xmax) # 使用连续Galerkin定义函数空间 # 每个单元上的p阶(拉格朗日)函数 V = FunctionSpace(mesh,"CG",p)

第二步：定义已提供的表达式

u0 = Constant(0) f = Constant(-1) g = Constant(1)

第三步：创建和应用Dirichlet边界条件

# 定义Dirichlet边界(x = 0) def dirichlet_boundary(x, on_boundary): return on_boundary and abs(x[0]) < DOLFIN_EPS # 在点x=0施加一个Dirichlet条件 bc = DirichletBC(V,u0,dirichlet_boundary)

第四步：定义变分问题

u = TrialFunction(V) v = TestFunction(V) a = inner(grad(u),grad(v))*dx L = f*v*dx g*v*ds

第五步：求解并绘图

# 求解 u = Function(V) solve(a == L, u, bc) # 绘图 plot(u)

求解2维泊松方程

• ​ (一个单位矩阵)
• ​(Dirichlet边界)
• ​ (Neumann边界)
• ​ (法向导数)
• ​ (源项)

类似的，也分五步求解

from dolfin import * # 第一步：为求解域创建网格，并定义函数空间 mesh = UnitSquareMesh(32, 32) V = FunctionSpace(mesh, "Lagrange", 1) # 第二步：定义已提供的表达式 u0 = Constant(0.0) f = Expression("10*exp(-(pow(x[0] - 0.5, 2) pow(x[1] - 0.5, 2)) / 0.02)", degree=2) g = Expression("sin(5*x[0])", degree=2) # 第三步：创建和应用`Dirichlet`边界条件 def boundary(x): return x[0] < DOLFIN_EPS or x[0] > 1.0 - DOLFIN_EPS bc = DirichletBC(V, u0, boundary) # 第四步：定义变分问题 u = TrialFunction(V) v = TestFunction(V) a = inner(grad(u), grad(v))*dx L = f*v*dx g*v*ds # 第五步：求解并绘图 u = Function(V) solve(a == L, u, bc) plot(u)

ParaView中查看解

# 将解保存为VTK格式 file = File("poisson.pvd") file << u

【结束】

,