绕多个障碍物的流动

在 Colab 中打开

一个稳态纳维-斯托克斯算例，既展示了 TensorMesh 的有限元机制，也同样突出了它的网格生成能力。脚本 examples/fluid/flow_obstacles/flow_obstacles.py 求解 \(\mathrm{Re} = 150\) 下穿过一个 \(3 \times 1\) 通道的不可压缩流动，通道内含有六个大小与位置各异的圆形障碍物。网格通过 MeshGen 借助 CSG（构造实体几何）以编程方式生成：从一个矩形出发，再减去六个圆盘。

离散化采用与 顶盖驱动方腔 相同的 Taylor-Hood P2-P1 混合单元对——与 圆柱绕流（涡脱落） 类似，MeshGen 网格是线性的，因此二次速度空间由 MixedElementAssembler 在其上按拓扑方式生成。

问题

\[\rho\, (\mathbf{u} \cdot \nabla)\mathbf{u} \;=\; -\nabla p + \mu\, \Delta \mathbf{u} \quad \text{in } \Omega, \qquad \nabla \cdot \mathbf{u} = 0,\]

其中通道 \((0, 3) \times (0, 1)\) 减去位于以下位置的六个圆盘

\[ \begin{align}\begin{aligned}(0.5, 0.5)\!:\,r\!=\!0.10, \quad (1.0, 0.3)\!:\,r\!=\!0.08, \quad (1.0, 0.7)\!:\,r\!=\!0.08,\\(1.5, 0.5)\!:\,r\!=\!0.12, \quad (2.0, 0.4)\!:\,r\!=\!0.07, \quad (2.0, 0.6)\!:\,r\!=\!0.07.\end{aligned}\end{align} \]

边界条件：

  • 入口（\(x = 0\)）：抛物型分布，\(u_x(y) = 4\,y\,(1 - y)\)，

  • 壁面（\(y = 0\)、\(y = 1\)）与障碍物表面：无滑移，

  • 出口（\(x = 3\)）："do-nothing"（自然边界条件），同时固定压力基准——完全不需要施加压力边界条件。

在 \(\mathrm{Re} = 150\) 时流动是稳态的（无脱落），每个障碍物的尾迹与下一个相互作用，解中呈现出靠得很近的圆柱之间特征性的狭窄射流。

通过 CSG 生成网格

这个脚本最有意思的部分是网格生成代码块——它简短、自成体系，并展示了 MeshGen 所暴露的 add_… / remove_… 风格的 CSG：

列表 22 examples/fluid/flow_obstacles/flow_obstacles.py（核心部分）
from tensormesh import MeshGen

gen = MeshGen(chara_length=1.0/n_grid)
gen.add_rectangle(0, 0, 3.0, 1.0)
for ox, oy, orad in [
    (0.5, 0.5, 0.10), (1.0, 0.3, 0.08), (1.0, 0.7, 0.08),
    (1.5, 0.5, 0.12), (2.0, 0.4, 0.07), (2.0, 0.6, 0.07),
]:
    gen.remove_circle(ox, oy, orad)
mesh = gen.gen().double()

在底层，MeshGen 将布尔运算与三角剖分交给 Gmsh 的 OpenCASCADE 后端处理，随后再按其内部单元排序转换回 Mesh。完整的 MeshGen API 参见 网格。

求解器

与 顶盖驱动方腔 相同——使用相同的双场 NavierStokesAssembler（Picard 线性化对流，四行标量被积函数），但具有更丰富的边界掩码。由于速度自由度包含按拓扑生成的边节点，掩码直接从布局自身的坐标和边界检测构建——无需对障碍圆进行几何距离测试：

x_u = layout.points("u")                                # P2 node coords
is_boundary = layout.split(layout.boundary_mask("u"))["u"][:, 0]
is_inlet = is_boundary & (x_u[:, 0] < eps)
is_outlet = x_u[:, 0] > LENGTH - eps
no_slip = is_boundary & ~is_inlet & ~is_outlet          # walls + obstacles

bc_mask = layout.dof_mask("u", node_mask=is_inlet | no_slip)
y_in = x_u[is_inlet, 1]
bc_val[layout.dof_mask("u", node_mask=is_inlet, component=0)] = \
    4.0 * y_in * (1.0 - y_in)

# Picard loop: lagged velocity on the P2 field's own DOFs
K = assembler(field_data={"w": ("u", w)})

障碍物表面自动由拓扑分类得出：它们是既非入口也非出口的边界节点。

输出

flow_obstacles.png 是最终图像：速度大小与压力等值线叠加在障碍物轮廓之上。图中可做的合理性检查包括：

  • \(x \approx 1.0\) 处上排两个障碍物之间的高速射流，

  • 位于 \((1.5, 0.5)\) 的最大障碍物后方的宽阔回流区，

  • 沿通道近乎均匀的压降。

绕 6 个随机圆形障碍物稳态流动的速度大小与压力

图 63 flow_obstacles.py 在 Re = 150 时的输出（Taylor-Hood P2-P1）。左：速度大小——高速射流在障碍对之间挤压，并在下游扩展为尾流。右：压力——每个障碍物上游形成滞止，尾流中出现低压区，通道内沿流向压力近乎均匀下降。

运行方式

cd examples/fluid/flow_obstacles
python flow_obstacles.py     # writes flow_obstacles.png

在默认分辨率（51k 自由度）下，Picard 迭代在 8 次迭代内收敛。

下一步