📚 课程目录第3章 FEniCS核心对象初识

3. FEniCS核心对象初识

在FEniCS中求解偏微分方程，本质上是在五个核心对象之间建立一条“从几何到数值解”的链条：Mesh描述计算区域，FunctionSpace定义离散函数空间，TrialFunction与TestFunction构成变分形式的未知量与检验函数，Function则承载最终求得的数值解。理解这五个对象的创建顺序与依赖关系，是读懂后续所有FEniCS代码的前提。本节我们用一个简单的一维泊松问题，在交互式Python中逐步构建这些对象，并观察它们之间的实际联系。

从几何开始：Mesh的创建

Mesh是有限元计算的几何基础，它把连续区域剖分为有限个单元。FEniCS提供了多种创建方式：UnitIntervalMesh、UnitSquareMesh、UnitCubeMesh用于规则区域，RectangleMesh可指定范围，而Mesh()配合外部网格文件则用于复杂几何。我们选择一维区间[0,1]，划分为10个单元：

from fenics import *

mesh = UnitIntervalMesh(10)
print(mesh.num_vertices())   # 11
print(mesh.num_cells())      # 10
print(mesh.coordinates()[:5])  # 前5个顶点坐标

这里num_vertices()返回顶点数11，num_cells()返回单元数10，符合一维区间10等分的拓扑结构。Mesh对象本身只包含几何与拓扑信息，不携带任何函数值，它是后续所有对象的“地基”。

在网格上定义函数空间：FunctionSpace

有了网格，下一步是选择离散函数空间。FEniCS通过FunctionSpace(mesh, family, degree)创建空间，其中family指定单元类型（如"Lagrange"、"DG"），degree指定多项式阶数。对于泊松问题，我们使用一阶拉格朗日连续元：

V = FunctionSpace(mesh, "Lagrange", 1)
print(V.dim())   # 11，与顶点数一致

一维一阶拉格朗日空间的自由度恰好对应每个顶点，因此V.dim()为11。若改为degree=2，维度会变为21，因为每个单元内部增加了自由度。FunctionSpace决定了离散解的精度与自由度分布，是连接几何与代数方程的关键桥梁。

变分形式的两个“占位符”：TrialFunction与TestFunction

FEniCS用符号计算的方式表达变分形式。对于泊松方程−Δu=f，其弱形式为：求u∈V，使得对任意v∈V，有∫∇u·∇v dx = ∫f v dx。这里u是未知的试探函数，v是任意的检验函数。在代码中，它们分别由TrialFunction(V)和TestFunction(V)创建：

u = TrialFunction(V)
v = TestFunction(V)

这两个对象并不是具体的数值函数，而是“符号占位符”，用于构造双线性形式a(u,v)和线性形式L(v)。它们必须依附于同一个FunctionSpace V，否则维度不匹配会导致组装失败。可以这样理解：TrialFunction代表“我们要求的那个解”，TestFunction代表“用来测试方程的任意函数”，二者共同把偏微分方程转化为线性代数系统。

组装与求解：从符号到数值

接下来定义源项f和边界条件。取f=1，并在两端施加零Dirichlet边界条件：

f = Constant(1.0)
a = dot(grad(u), grad(v)) * dx
L = f * v * dx

bc = DirichletBC(V, Constant(0.0), "on_boundary")
uh = Function(V)
solve(a == L, uh, bc)

这里Function(V)创建了一个属于空间V的具体函数对象，初始值为零，求解后被填充为数值解。solve内部完成了组装矩阵、施加边界条件、求解线性系统的全过程。求解完成后，uh就是一个可以在任意点求值、可以绘图、可以计算误差的真实函数。

五个对象的关系全景

把上述过程串起来，五个对象形成一条清晰的依赖链：Mesh提供几何，FunctionSpace在Mesh上定义离散空间，TrialFunction和TestFunction从FunctionSpace派生用于构造变分形式，Function同样属于FunctionSpace并承载求解结果。下面的SVG展示了这一结构与数据流向：

图中蓝色箭头表示空间依赖：Mesh被FunctionSpace使用，FunctionSpace派生出TrialFunction、TestFunction和Function。橙色箭头表示变分形式的构造与求解数据流：TrialFunction与TestFunction共同进入双线性形式a(u,v)和线性形式L(v)，组装求解后结果写入Function。这条链条中任何一环出错——比如FunctionSpace与Mesh不匹配、TrialFunction与TestFunction来自不同空间——都会在组装阶段报错。

验证与观察

求解完成后，可以用几个简单操作确认对象状态。例如在区间中点x=0.5处求值，理论解为u(0.5)=0.125（对−u''=1、u(0)=u(1)=0），一阶有限元在10单元下会给出接近但略有偏差的数值：

import numpy as np
print(uh(0.5))            # 约0.125，具体值随网格略有差异
print(uh.compute_vertex_values()[:5])  # 前5个顶点上的解值
print(assemble(uh*dx))    # 解的积分，约0.0417

这些输出说明uh已经是一个可操作的数值函数，而不再是符号占位符。若把V改为二次元重新求解，uh(0.5)会更接近0.125，这正是FunctionSpace中degree参数影响精度的直接体现。

记忆口诀：Mesh是地，FunctionSpace是尺，TrialFunction是“求什么”，TestFunction是“拿什么测”，Function是“算出来的答案”。五者顺序不可颠倒，空间必须一致。

掌握这五个对象的创建与关系后，后续无论是处理二维热传导、弹性力学还是Navier-Stokes方程，代码骨架都遵循同一模式：建网格、选空间、定试探与检验函数、写变分形式、求解到Function。本节的一维泊松例子虽小，却完整覆盖了这条主线的每一个环节，值得反复动手修改参数观察变化。