泰勒-格林涡（收敛性研究）

泰勒-格林涡是具有已知精确解的经典不可压缩流动基准：在正方形 \([0, 2\pi]^2\) 上的一个周期性衰减涡，其速度与压力都可写成闭式解。脚本 examples/fluid/taylor_green/taylor_green.py 用它来对瞬态纳维-斯托克斯求解器开展定量的 h 收敛性研究——这是 流体力学 中唯一以验证而非可视化为目的的算例。

精确解

对于黏度 \(\nu\) 与时间 \(t\)，

\[\begin{split}u(x, y, t) &= -\cos(x)\, \sin(y)\, e^{-2\nu t}, \\ v(x, y, t) &= \phantom{-}\sin(x)\, \cos(y)\, e^{-2\nu t}, \\ p(x, y, t) &= -\tfrac14 \bigl(\cos(2x) + \cos(2y)\bigr)\, e^{-4\nu t}.\end{split}\]

速度以速率 \(2\nu\)（运动黏度的两倍）指数衰减；压力以速率 \(4\nu\) 衰减。代入即可验证 \((u, v, p)\) 精确满足不可压缩纳维-斯托克斯方程。

脚本在每个边界节点上施加等于当前时刻精确解的狄利克雷速度。这样既绕开了对周期性边界的需求，也将空间 / 时间离散误差单独隔离出来。

求解器

来自 圆柱绕流（涡脱落） 的瞬态 NavierStokesTransientAssembler 被原样复用——时间上采用隐式欧拉法，空间上采用带自适应 \(\tau\) 的 SUPG/PSPG。仅有的改动是：

  • 计算域为 \([0, 2\pi]^2\)，而非圆柱通道，

  • 狄利克雷速度取自精确解，而非抛物型入口，

  • 每一步都将相对于解析速度与压力的 \(L^2\) 误差写入结果文件。

列表 19 examples/fluid/taylor_green/taylor_green.py（核心部分）
def exact_velocity(points, t, nu):
    x, y = points[:, 0], points[:, 1]
    decay = math.exp(-2.0 * nu * t)
    u = -torch.cos(x) * torch.sin(y) * decay
    v =  torch.sin(x) * torch.cos(y) * decay
    return torch.stack([u, v], dim=1)

for n_grid in [10, 20, 40, 80]:
    mesh = Mesh.gen_rectangle(left=0, right=2*math.pi,
                              bottom=0, top=2*math.pi,
                              chara_length=2*math.pi/n_grid)
    # …assemble, set Dirichlet to exact velocity, time-step…
    err_u_l2 = mass_weighted_norm(u_fem - u_exact)
    err_p_l2 = mass_weighted_norm(p_fem - p_exact)

h 收敛性

脚本让网格尺寸 \(h = 2\pi / n_\text{grid}\) 遍历一组取值，并在固定的终止时刻记录 \(L^2\) 速度误差与压力误差。对于带 SUPG/PSPG 的 P1-P1 单元，预期为：

  • \(\|\mathbf{u}_h - \mathbf{u}\|_{L^2} = \mathcal{O}(h^2)\)

  • \(\|p_h - p\|_{L^2} = \mathcal{O}(h)\)

控制台输出为一张表格：

   n_grid     h       err_u_L2    rate_u    err_p_L2    rate_p
---------------------------------------------------------------
     10    0.628    1.3e-1                  4.5e-1
     20    0.314    3.4e-2     1.9          2.3e-1     1.0
     40    0.157    8.7e-3     2.0          1.2e-1     1.0
     80    0.079    2.2e-3     2.0          5.8e-2     1.0

（数值为近似值，会随所选的 \(\nu\) 与 \(\Delta t\) 略有变化）。图 taylor_green_convergence.png 在对数-对数坐标轴上展示相同数据，并附有参考的 \(h^2\) 与 \(h^1\) 斜率。

如果观测到的收敛阶与预期的 \(\mathcal{O}(h^2)\) / \(\mathcal{O}(h)\) 不符，那么求解器中必有问题——最常见的是稳定化或边界条件施加中的隐蔽 bug。

质量加权误差范数

脚本通过质量矩阵计算离散的 \(L^2\) 范数：

\[\|e_h\|_{L^2}^2 \;=\; e_h^T\, M\, e_h \;=\; \sum_K \int_K e_h^2 \,\mathrm{d}\Omega,\]

其中质量矩阵由 MassElementAssembler 装配。这是有限元误差分析的正确范数（它在各单元上积分，并以基函数加权），不同于节点误差向量的简单欧几里得范数。

输出

  • 控制台表格——每个加密步上的收敛阶。

  • ``taylor_green_convergence.png``——误差关于网格尺寸的对数-对数图，并附参考斜率。

  • ``taylor_green_results.png``——最细网格的三联面板快照：速度、压力、速度误差大小。

  • ``taylor_green.mp4``（可选）——衰减涡的动画。

t=0.5 时的泰勒-格林涡——涡量＋流线，速度＋速度矢量

图 49 taylor_green.py 在 \(t=0.5\) 时的输出。左：叠加流线的涡量场——经典的周期性 \(2\times2\) 交替符号涡阵列。右：带速度矢量的速度大小，展示了相邻卷之间特征性的“鞍”状结构。相对解析初始条件，幅值已衰减了 \(\exp(-2\nu t)\)，脚本的收敛性研究即以此为参考真值。

运行方式

cd examples/fluid/taylor_green
python taylor_green.py      # writes convergence + results pngs

下一步