Skip to content

第25章 几何求值(Geometry Evaluation):参数化、法向量与曲率


25.1 学习目标

  • 理解曲面参数化(Surface Parametrization)的概念:CAD 几何中每个曲面自带 (u,v) 参数坐标系,掌握 getParametrization() 在物理坐标与参数坐标之间建立双向映射
  • 掌握通过 getNormal() 获取曲面上任意参数坐标处的单位法向量(Unit Normal),理解其与 (dim,tag) 实体对 API(如 getCurvature)在参数签名上的关键区别
  • 学会使用 getPrincipalCurvatures() 获取主曲率(Principal Curvatures)及其主方向(Principal Directions),理解曲率在自适应网格细化中的判据作用
  • 掌握曲线求值 API 族:getValue()(坐标)、getDerivative()(切向量)、getSecondDerivative()(二阶导),以及 reparametrizeOnSurface() 曲线-曲面联合求值
  • 了解上述 API 在 CAE 中的典型应用:压力载荷施加(需法向)、自适应网格细化(需曲率)、边界条件位置计算(需参数化)

25.2 核心概念说明

25.2.1 曲面参数化(Surface Parametrization)

在 CAD 内核中,每张曲面(Surface)都定义了一个从二维参数域 (u,v) 到三维物理空间 (x,y,z) 的映射:

S(u,v) = (x(u,v), y(u,v), z(u,v))
其中 (u,v) 属于参数定义域 [u_min, u_max] x [v_min, v_max]

参数域范围因几何类型而异:平面和裁剪曲面(Trimmed Surface)的参数域通常受裁剪边界限制,球面和圆柱面等周期曲面则可能覆盖 [0, 2*pi] 范围。

Gmsh 通过两个 API 实现物理坐标与参数坐标的双向转换:

方向 API 输入 输出
物理 -> 参数 getParametrization(dim, tag, coord, paramCoord) (x,y,z) (u,v) 或 t
参数 -> 物理 getValue(dim, tag, paramCoord, coord) (u,v) 或 t (x,y,z)

25.2.2 法向量(Normal)与曲率(Curvature)

曲面上一点的法向量和曲率是 CAE 前处理中最重要的微分几何量:

  • 法向量(Normal):垂直于曲面切平面的单位向量。在有限元分析中用于施加压力载荷(压力沿法向作用)、接触检测、壳单元局部坐标系建立。
  • 主曲率(Principal Curvatures):曲面上一点沿不同方向的弯曲程度不同,存在一对正交方向使弯曲程度达到极值,分别称为最大主曲率 k1 和最小主曲率 k2
  • 高斯曲率(Gaussian Curvature)K = k1 * k2:判断曲面局部形状(K>0 椭圆点, K<0 双曲点, K=0 抛物点)
  • 平均曲率(Mean Curvature)H = (k1 + k2) / 2:用于曲率自适应网格细化判据

API 签名关键差异:Gmsh API 中存在不一致——getNormal()getPrincipalCurvatures() 仅接收曲面标签(tag),而 getCurvature() 接收 (dim, tag) 对。这是 Gmsh API 的历史遗留问题,使用时务必注意。

25.2.3 曲线求值(Curve Evaluation)

对于一维曲线实体,Gmsh 提供三个层次的求值 API:

API 返回 数学含义 典型用途
getValue(dim, tag, t, xyz) 坐标 C(t) = (x,y,z) 边界条件位置计算
getDerivative(dim, tag, t, deriv) 切向量 C'(t) = (dx/dt, dy/dt, dz/dt) 切向载荷、流线方向
getSecondDerivative(dim, tag, t, deriv2) 二阶导 C''(t) 曲率计算、弯管应力估算

25.2.4 CAE 应用场景

  1. 压力载荷施加p * n 需要曲面上每个积分点的单位法向量 n。从网格节点获取参数坐标后,通过 getNormal() 获取精确的 CAD 法向(而非网格近似法向),避免网格粗糙时法向误差导致的载荷方向偏差。
  2. 自适应网格细化:在曲率大的区域(如倒角根部、小孔边缘)加密网格。通过 getPrincipalCurvatures() 评估局部曲率,设定曲率阈值触发局部细化。
  3. 边界条件位置计算:利用 getValue() 在曲线参数坐标上精确计算约束施加点或载荷施加点的物理坐标,适用于参数化边界条件脚本。

25.3 C++ 代码逐段讲解

25.3.1 头文件与初始化

#include <set>
#include <vector>
#include <cmath>
#include <numeric>
#include <gmsh.h>

int main(int argc, char **argv)
{
  gmsh::initialize(argc, argv);
  gmsh::model::add("ch25");

解释<numeric> 用于 std::iota 等序列生成,<cmath> 提供 std::sqrtstd::abs 等数学函数。gmsh::model::add("ch25") 为模型命名。


25.3.2 创建融合几何体并生成曲面网格

  // ---- 1. 创建球与立方体的融合几何 ----
  int s = gmsh::model::occ::addSphere(0, 0, 0, 1);     // 单位球
  int b = gmsh::model::occ::addBox(0.5, 0, 0, 1.3, 2, 3); // 与球相交的盒子
  std::vector<std::pair<int, int>> ov;
  std::vector<std::vector<std::pair<int, int>>> ovv;
  gmsh::model::occ::fuse({{3, s}}, {{3, b}}, ov, ovv);
  gmsh::model::occ::synchronize();

  // 生成二维面网格(仅曲面)
  gmsh::model::mesh::generate(2);

解释:球和长方体做布尔融合(Boolean Union),产生一个既有球面区域又有平面区域的复合几何体。这种"曲面+平面"的组合最能体现法向量和曲率 API 的价值:球面区域的法向量处处变化且曲率非零,平面区域法向量恒定且曲率为零。mesh::generate(2) 仅生成面网格,不生成体网格,以节省计算量。


25.3.3 获取曲面网格节点及其参数坐标

  // ---- 2. 遍历所有面,获取网格节点及参数坐标 ----
  std::vector<std::pair<int, int>> entities;
  gmsh::model::getEntities(entities, 2); // dim=2: 面

  std::vector<double> allNormals;   // 用于后处理显示: 法向量
  std::vector<double> allCurvatures; // 用于后处理显示: 最大主曲率

  for(auto e : entities) {
    int surfaceTag = e.second;

    // 获取面上的网格节点
    // includeBoundary=true: 边界节点会被重新参数化到当前面上
    std::vector<std::size_t> nodeTags;
    std::vector<double> coord, paramCoord;
    gmsh::model::mesh::getNodes(nodeTags, coord, paramCoord,
                                2, surfaceTag, true);

解释getNodes() 的最后一个参数 includeBoundary=true 至关重要。内部节点(interior nodes)存储了参数坐标,但边界节点(boundary nodes)属于多个面共有——当设为 true 时,Gmsh 会将边界节点重新参数化(reparametrize)到当前面上,从而获得在该面上的 (u,v) 坐标。如果设为 false,则边界节点的 paramCoord 为空。

coord 的数据布局为 [x1,y1,z1, x2,y2,z2, ...]paramCoord 的布局为 [u1,v1, u2,v2, ...],两者按节点一一对应。


25.3.4 获取曲面法向量

    // ---- 3. 获取面上各节点参数坐标处的精确法向量 ----
    std::vector<double> normals;
    gmsh::model::getNormal(surfaceTag, paramCoord, normals);
    // 注意: getNormal() 的第一个参数是曲面标签(tag),不是(dim,tag)对
    // 这与 getCurvature(dim, tag, ...) 的签名不同!

解释getNormal(tag, parametricCoord, normals) —— 第一个参数是面标签(整数),不是 (dim, tag)。这是 Gmsh API 中的关键不一致点,也是 CAE 开发者最常见的错误来源。

normals 返回的单位法向量布局为 [nx1,ny1,nz1, nx2,ny2,nz2, ...],与输入的 paramCoord 按点一一对应。法向量的方向遵循右手定则(基于曲面的 (u,v) 参数化方向 Su x Sv)。

这里的法向量是直接从 CAD 几何求值得到的精确法向,而非从网格三角面片计算的法向近似值。在压力载荷施加场景中,使用 getNormal() 的法向可以避免网格粗糙时法向偏差导致的载荷方向错误。


25.3.5 获取曲面主曲率

    // ---- 4. 获取主曲率及主方向 ----
    std::vector<double> curvMax, curvMin, dirMax, dirMin;
    gmsh::model::getPrincipalCurvatures(surfaceTag, paramCoord,
                                         curvMax, curvMin, dirMax, dirMin);
    // getPrincipalCurvatures(tag, ...) 同样仅使用面标签,不使用(dim,tag)对

解释getPrincipalCurvatures(tag, parametricCoord, curvatureMax, curvatureMin, directionMax, directionMin) 在一次调用中返回四个向量:

输出向量 长度 (N个点) 内容
curvatureMax N 每个点的最大主曲率 k1(正值=凸,负值=凹)
curvatureMin N 每个点的最小主曲率 k2
directionMax 3N k1 对应的主方向单位向量 [dx1,dy1,dz1, dx2,dy2,dz2, ...]
directionMin 3N k2 对应的主方向单位向量

两个主方向始终正交。对于球面,k1 = k2 = 1/R(脐点,任意两个正交方向都可视为主方向)。对于平面,k1 = k2 = 0。对于圆柱面,k1 = 1/R(周向),k2 = 0(轴向)。

注意:还有一个简化 API getCurvature(dim, tag, parametricCoord, curvatures),它只返回每个点的最大曲率(单值),且使用 (dim, tag) 对作为前两个参数——与 getPrincipalCurvatures 的签名风格不同。


25.3.6 整理法向量与曲率数据用于后处理显示

    // ---- 5. 整理数据用于后处理视图 ----
    for(std::size_t i = 0; i < coord.size(); i += 3) {
      std::size_t idx = i / 3; // 点索引
      // 法向量视图: 位置 + 法向量 (VP = Vector on Points, 6个double/点)
      allNormals.insert(allNormals.end(), {
          coord[i], coord[i+1], coord[i+2],          // x, y, z
          normals[i], normals[i+1], normals[i+2]      // nx, ny, nz
      });
      // 曲率视图: 位置 + 标量值 (SP = Scalar on Points, 4个double/点)
      allCurvatures.insert(allCurvatures.end(), {
          coord[i], coord[i+1], coord[i+2],           // x, y, z
          curvMax[idx]                                 // 最大主曲率
      });
    }
  } // end for each surface

解释:将数据组织为 list-based 后处理视图格式: - "VP" (Vector on Points):每个数据条目 6 个 double — (x,y,z, vx,vy,vz) - "SP" (Scalar on Points):每个数据条目 4 个 double — (x,y,z, scalar)

这里选择了最大主曲率 curvMax 作为标量场显示;也可以选择平均曲率 (curvMax + curvMin) / 2 或高斯曲率 curvMax * curvMin 作为细化判据。


25.3.7 创建后处理视图

  // ---- 6. 创建法向量矢量视图 ----
  int vn = gmsh::view::add("Surface normals");
  gmsh::view::addListData(vn, "VP", allNormals.size() / 6, allNormals);
  gmsh::view::option::setNumber(vn, "ShowScale", 0);     // 不显示色标
  gmsh::view::option::setNumber(vn, "ArrowSizeMax", 30); // 箭头最大 30 像素
  gmsh::view::option::setNumber(vn, "ColormapNumber", 19);

  // ---- 7. 创建曲率标量视图 ----
  int vc = gmsh::view::add("Principal curvatures (max)");
  gmsh::view::addListData(vc, "SP", allCurvatures.size() / 4, allCurvatures);
  gmsh::view::option::setNumber(vc, "ShowScale", 1);
  gmsh::view::option::setNumber(vc, "IntervalsType", 2); // 连续色带

解释addListData() 创建独立于网格的 list-based 视图。"VP" 格式将法向量渲染为带箭头的矢量线;"SP" 格式将曲率值渲染为着色点云。这些视图可在 Gmsh GUI 中独立开关、叠加显示。


25.3.8 曲线参数化边界与求值

  // ---- 8. 获取某条曲线的参数化边界并等距采样求值 ----
  // 以曲线 5 为例(融合后自动生成的边)
  std::vector<double> tMin, tMax;
  gmsh::model::getParametrizationBounds(1, 5, tMin, tMax);
  // 对于曲线(dim=1): tMin/tMax 各含1个值,分别为参数 t 的下界和上界

  int N = 15; // 采样点数
  std::vector<double> tSamples, xyz1;
  for(int i = 0; i < N; i++) {
    double t = tMin[0] + i * (tMax[0] - tMin[0]) / (N - 1);
    tSamples.push_back(t);
  }
  gmsh::model::getValue(1, 5, tSamples, xyz1);
  // xyz1 布局: [x1,y1,z1, x2,y2,z2, ...], 共 N 个点

解释getParametrizationBounds(1, 5, tMin, tMax) 返回曲线参数 t 的定义域 [t_min, t_max]。对于 dim=1 的曲线,tMintMax 各为一个元素的向量;对于 dim=2 的曲面,tMin=[u_min, v_min]tMax=[u_max, v_max]

getValue(1, 5, tSamples, xyz1) 在曲线上等距采样 15 个参数点,计算对应的三维坐标。输出的 xyz1 包含 3*N 个 double。


25.3.9 曲线导数与二阶导数

  // ---- 9. 获取曲线在采样点的切向量(一阶导)和二阶导 ----
  std::vector<double> deriv1, deriv2;
  gmsh::model::getDerivative(1, 5, tSamples, deriv1);
  // deriv1 布局: [dx/dt, dy/dt, dz/dt, ...], 每点 3 个分量
  gmsh::model::getSecondDerivative(1, 5, tSamples, deriv2);
  // deriv2 布局: [d2x/dt2, d2y/dt2, d2z/dt2, ...], 每点 3 个分量

  // 计算第一个采样点处的切向量模长和曲线曲率
  double dx = deriv1[0], dy = deriv1[1], dz = deriv1[2];
  double ds = std::sqrt(dx*dx + dy*dy + dz*dz);     // 切向量模长 = |C'(t)|
  double ddx = deriv2[0], ddy = deriv2[1], ddz = deriv2[2];
  // 曲线曲率: kappa = |C' x C''| / |C'|^3
  double crossX = dy*ddz - dz*ddy;
  double crossY = dz*ddx - dx*ddz;
  double crossZ = dx*ddy - dy*ddx;
  double kappa = std::sqrt(crossX*crossX + crossY*crossY + crossZ*crossZ)
                 / (ds * ds * ds);
  gmsh::logger::write("Curve 5 curvature at t[0] = " + std::to_string(kappa));

解释getDerivative() 返回的是参数导数 dC/dt(非弧长导数)。切向量方向可直接用于局部坐标系对齐或切向载荷方向定义。结合 getDerivative()getSecondDerivative() 可计算曲线的 Frenet 曲率:

kappa = |C'(t) x C''(t)| / |C'(t)|^3

这在 CAE 中可用于:评估弯管处的应力集中系数(曲率越大,弯管应力集中越显著)、判断曲线网格加密区域。


25.3.10 曲线在曲面上的重参数化(Reparametrization)

  // ---- 10. 将曲线重参数化到曲面上,验证一致性 ----
  // 将曲线 5 上的采样点重新参数化到面 1 上
  std::vector<double> uvOnSurface, xyz2;
  gmsh::model::reparametrizeOnSurface(1, 5, tSamples, 1, uvOnSurface);
  // uvOnSurface 布局: [u1,v1, u2,v2, ...]  — 曲线采样点在面1上的(u,v)坐标

  // 用面 1 的参数坐标重新求值
  gmsh::model::getValue(2, 1, uvOnSurface, xyz2);

  // 验证两种求值方式得到相同结果
  double maxDiff = 0.0;
  for(std::size_t i = 0; i < xyz1.size(); i++)
    maxDiff = std::max(maxDiff, std::abs(xyz1[i] - xyz2[i]));
  if(maxDiff < 1e-12)
    gmsh::logger::write("Curve & surface evaluation match! (max diff < 1e-12)");
  else
    gmsh::logger::write("Mismatch: max diff = " + std::to_string(maxDiff), "warning");

解释reparametrizeOnSurface(dim, tag, paramCoord, surfaceTag, surfaceParamCoord) 实现了一条曲线(或其上的点)在指定曲面上的重新参数化。该功能仅在部分 CAD 内核中可用。

这里验证了两种求值路径的一致性: 1. 直接在曲线 5 上求值 → xyz1 2. 将曲线 5 的参数坐标重参数化到面 1 → 在面 1 上求值 → xyz2

两者应得到相同的物理坐标(因为曲线 5 是面 1 的边界)。这种验证是 CAE 工作流中的良好实践——在边界条件施加时,如果曲线和曲面求值不一致,说明 CAD 模型可能存在几何不兼容问题(如边和面之间有微小间隙)。


25.3.11 物理坐标到参数坐标的反向映射

  // ---- 11. 演示物理坐标 -> 参数坐标的反向映射 ----
  // 取第一个面(面 1)上网格的第一个节点,演示反向参数化
  std::vector<std::size_t> nTags;
  std::vector<double> nCoord, nParam;
  gmsh::model::mesh::getNodes(nTags, nCoord, nParam, 2, 1, true);
  if(nCoord.size() >= 3) {
    // 取第一个节点坐标求其参数坐标
    std::vector<double> singleCoord = {nCoord[0], nCoord[1], nCoord[2]};
    std::vector<double> fromCoord;
    gmsh::model::getParametrization(2, 1, singleCoord, fromCoord);
    // fromCoord[0] = u, fromCoord[1] = v  — 该节点在面1上的参数坐标
    gmsh::logger::write(
      "Node " + std::to_string(nTags[0]) +
      " at (" + std::to_string(nCoord[0]) + "," +
                std::to_string(nCoord[1]) + "," +
                std::to_string(nCoord[2]) + ")" +
      " -> (u=" + std::to_string(fromCoord[0]) +
      ", v=" + std::to_string(fromCoord[1]) + ")");
  }

解释getParametrization(dim, tag, coord, parametricCoord) 实现物理坐标到参数坐标的反向映射。输入物理坐标 (x,y,z),输出参数坐标 (u,v)(曲面)或 t(曲线)。这个方向通常比正向映射(参数→物理)更复杂,因为需要求解非线性反问题。对于距离曲面较远的点,可能无法正确映射或返回近似值。

CAE 应用:在接触分析中,从节点(slave node)投影到主面(master surface)上以计算穿透量和接触力时,需要反向参数化。


25.3.12 GUI 启动

  // ---- 12. 启动 GUI 查看结果 ----
  std::set<std::string> args(argv, argv + argc);
  if(!args.count("-nopopup")) gmsh::fltk::run();

  gmsh::finalize();
  return 0;
}

解释-nopopup 命令行参数支持无头批处理模式。在 GUI 模式下,法向量视图(箭头)和曲率视图(着色点云)可叠加在网格上查看。


25.4 完整可运行代码

// =============================================================================
//  Chapter 25: Geometry Evaluation — Parametrization, Normals & Curvatures
//
//  功能: 演示 Gmsh 几何求值 API 的完整用法,包括:
//        曲面参数化、法向量计算、主曲率提取、曲线求值与导数、
//        曲线-曲面重参数化、物理-参数坐标双向映射
//
//  编译 (需链接 Gmsh SDK):
//    g++ -o ch25_geometry_eval ch25_geometry_eval.cpp -lgmsh -std=c++11
//
//  运行:
//    ./ch25_geometry_eval              (带 GUI 交互查看法向量和曲率)
//    ./ch25_geometry_eval -nopopup     (仅运行计算,不启动 GUI)
// =============================================================================

#include <set>
#include <vector>
#include <cmath>
#include <string>
#include <algorithm>
#include <gmsh.h>

int main(int argc, char **argv)
{
  gmsh::initialize(argc, argv);
  gmsh::model::add("ch25");

  // ===== 1. 创建球与立方体的融合几何 =====
  int s = gmsh::model::occ::addSphere(0, 0, 0, 1);
  int b = gmsh::model::occ::addBox(0.5, 0, 0, 1.3, 2, 3);
  std::vector<std::pair<int, int>> ov;
  std::vector<std::vector<std::pair<int, int>>> ovv;
  gmsh::model::occ::fuse({{3, s}}, {{3, b}}, ov, ovv);
  gmsh::model::occ::synchronize();

  // 生成二维面网格
  gmsh::model::mesh::generate(2);

  // ===== 2. 遍历所有面,获取法向量和主曲率 =====
  std::vector<std::pair<int, int>> entities;
  gmsh::model::getEntities(entities, 2);

  std::vector<double> allNormals;
  std::vector<double> allCurvatures;

  for(auto e : entities) {
    int surfaceTag = e.second;

    // 获取面上的网格节点(含边界节点重参数化)
    std::vector<std::size_t> nodeTags;
    std::vector<double> coord, paramCoord;
    gmsh::model::mesh::getNodes(nodeTags, coord, paramCoord,
                                2, surfaceTag, true);

    // 获取精确法向量 (tag 仅用面标签,非 (dim,tag) 对)
    std::vector<double> normals;
    gmsh::model::getNormal(surfaceTag, paramCoord, normals);

    // 获取主曲率 (tag 仅用面标签)
    std::vector<double> curvMax, curvMin, dirMax, dirMin;
    gmsh::model::getPrincipalCurvatures(surfaceTag, paramCoord,
                                         curvMax, curvMin, dirMax, dirMin);

    // 整理为后处理视图数据格式
    for(std::size_t i = 0; i < coord.size(); i += 3) {
      std::size_t idx = i / 3;
      // VP: Vector on Points (6 doubles per point)
      allNormals.insert(allNormals.end(), {
          coord[i], coord[i+1], coord[i+2],
          normals[i], normals[i+1], normals[i+2]
      });
      // SP: Scalar on Points (4 doubles per point)
      allCurvatures.insert(allCurvatures.end(), {
          coord[i], coord[i+1], coord[i+2],
          curvMax[idx]
      });
    }
  }

  // ===== 3. 创建后处理视图 =====
  int vn = gmsh::view::add("Surface normals");
  gmsh::view::addListData(vn, "VP", allNormals.size() / 6, allNormals);
  gmsh::view::option::setNumber(vn, "ShowScale", 0);
  gmsh::view::option::setNumber(vn, "ArrowSizeMax", 30);
  gmsh::view::option::setNumber(vn, "ColormapNumber", 19);

  int vc = gmsh::view::add("Principal curvatures (max)");
  gmsh::view::addListData(vc, "SP", allCurvatures.size() / 4, allCurvatures);
  gmsh::view::option::setNumber(vc, "ShowScale", 1);
  gmsh::view::option::setNumber(vc, "IntervalsType", 2);

  // ===== 4. 曲线参数化边界与等距采样求值 =====
  std::vector<double> tMin, tMax;
  gmsh::model::getParametrizationBounds(1, 5, tMin, tMax);

  int N = 15;
  std::vector<double> tSamples, xyz1;
  for(int i = 0; i < N; i++) {
    double t = tMin[0] + i * (tMax[0] - tMin[0]) / (N - 1);
    tSamples.push_back(t);
  }
  gmsh::model::getValue(1, 5, tSamples, xyz1);

  // ===== 5. 曲线切向量与二阶导数 =====
  std::vector<double> deriv1, deriv2;
  gmsh::model::getDerivative(1, 5, tSamples, deriv1);
  gmsh::model::getSecondDerivative(1, 5, tSamples, deriv2);

  // 计算第一个采样点的 Frenet 曲率
  double dx = deriv1[0], dy = deriv1[1], dz = deriv1[2];
  double ds = std::sqrt(dx*dx + dy*dy + dz*dz);
  double ddx = deriv2[0], ddy = deriv2[1], ddz = deriv2[2];
  double cx = dy*ddz - dz*ddy;
  double cy = dz*ddx - dx*ddz;
  double cz = dx*ddy - dy*ddx;
  double kappa = std::sqrt(cx*cx + cy*cy + cz*cz) / (ds * ds * ds);
  gmsh::logger::write("Curve 5 curvature at first sample: " +
                       std::to_string(kappa));

  // ===== 6. 曲线在曲面上的重参数化与一致性验证 =====
  std::vector<double> uvOnSurface, xyz2;
  gmsh::model::reparametrizeOnSurface(1, 5, tSamples, 1, uvOnSurface);
  gmsh::model::getValue(2, 1, uvOnSurface, xyz2);

  double maxDiff = 0.0;
  for(std::size_t i = 0; i < xyz1.size(); i++)
    maxDiff = std::max(maxDiff, std::abs(xyz1[i] - xyz2[i]));
  if(maxDiff < 1e-12)
    gmsh::logger::write("Curve & surface evaluation match! (max diff < 1e-12)");
  else
    gmsh::logger::write("Mismatch: max diff = " + std::to_string(maxDiff),
                        "warning");

  // ===== 7. 物理坐标 -> 参数坐标的反向映射 =====
  std::vector<std::size_t> nTags;
  std::vector<double> nCoord, nParam;
  gmsh::model::mesh::getNodes(nTags, nCoord, nParam, 2, 1, true);
  if(nCoord.size() >= 3) {
    std::vector<double> singleCoord = {nCoord[0], nCoord[1], nCoord[2]};
    std::vector<double> fromCoord;
    gmsh::model::getParametrization(2, 1, singleCoord, fromCoord);
    gmsh::logger::write(
      "Node " + std::to_string(nTags[0]) +
      " -> (u=" + std::to_string(fromCoord[0]) +
      ", v=" + std::to_string(fromCoord[1]) + ")");
  }

  // ===== 8. 启动 GUI =====
  std::set<std::string> args(argv, argv + argc);
  if(!args.count("-nopopup")) gmsh::fltk::run();

  gmsh::finalize();
  return 0;
}

25.5 关键 API 速查表

25.5.1 曲面求值(Surface Evaluation)

API 签名 说明
getParametrization (dim, tag, coord, paramCoord) 物理坐标 -> 参数坐标。dim=1 返回 [t],dim=2 返回 [u,v]
getParametrizationBounds (dim, tag, min, max) 获取参数域范围。dim=1: min/max 各 1 个 t 值;dim=2: 各 2 个值 [u,v]
getValue (dim, tag, paramCoord, coord) 参数坐标 -> 物理坐标。dim=1: paramCoord=[t1,t2,...];dim=2: paramCoord=[u1,v1,u2,v2,...]
getNormal (tag, paramCoord, normals) 获取曲面法向量。tag 仅为面标签(整数),不使用 (dim,tag)
getCurvature (dim, tag, paramCoord, curvatures) 获取最大曲率(单值/点)。使用 (dim,tag)
getPrincipalCurvatures (tag, paramCoord, curvMax, curvMin, dirMax, dirMin) 获取两个主曲率及其 3D 方向向量。tag 仅为面标签
reparametrizeOnSurface (dim, tag, paramCoord, surfaceTag, surfaceParamCoord) 将边界实体(曲线/点)的参数坐标转换到目标曲面的参数坐标

25.5.2 曲线求值(Curve Evaluation)

API 签名 说明
getValue (1, tag, t, xyz) 曲线在参数 t 处的三维坐标。t 向量: [t1,t2,...],xyz: [x1,y1,z1,...]
getDerivative (1, tag, t, deriv) 一阶参数导数(切向量方向)。deriv: [dx/dt, dy/dt, dz/dt, ...]
getSecondDerivative (1, tag, t, deriv2) 二阶参数导数。deriv2: [d2x/dt2, d2y/dt2, d2z/dt2, ...]

25.5.3 曲面求值(Surface Evaluation)— dim=2 时的导数

API 输出布局(每点) 说明
getDerivative [dx/du, dy/du, dz/du, dx/dv, dy/dv, dz/dv] 6 个分量/点:对 u 和对 v 的偏导数
getSecondDerivative [d2x/du2, d2y/du2, d2z/du2, d2x/dv2, d2y/dv2, d2z/dv2, d2x/dudv, d2y/dudv, d2z/dudv] 9 个分量/点:二阶纯导数和混合偏导数

25.5.4 API 签名风格对比(重要)

API 实体标识方式 示例
getNormal tag(面标签) getNormal(5, param, norm)
getPrincipalCurvatures tag(面标签) getPrincipalCurvatures(5, param, ...)
getCurvature (dim, tag) getCurvature(2, 5, param, curv)
getValue (dim, tag) getValue(2, 5, paramCoord, xyz)
getDerivative (dim, tag) getDerivative(1, 5, t, deriv)
getSecondDerivative (dim, tag) getSecondDerivative(1, 5, t, deriv2)
getParametrization (dim, tag) getParametrization(2, 5, coord, param)
getParametrizationBounds (dim, tag) getParametrizationBounds(1, 5, min, max)

25.6 注意事项与最佳实践

  1. getNormal vs getCurvature 的签名不一致getNormal(tag, ...) 的第一个参数是面标签(int),而 getCurvature(dim, tag, ...) 的第一个参数是维度(dim)。这是 Gmsh API 的历史遗留问题,混用时极易出错。建议在代码中以注释标注(如本节示例),或封装一层统一接口。

  2. 边界节点的参数坐标mesh::getNodes(..., includeBoundary=true) 时,边界节点会被重新参数化(reparametrized)到当前面上。相邻两个面对同一边界节点的参数化结果通常不同(因为参数域不同),这是正常现象。

  3. getValue 的 dim=0 特殊情况:对于点实体(dim=0),getValue(0, tag, {}, coord)parametricCoord 参数为空向量 {},因为点没有参数坐标。

  4. getNormal 返回的是 CAD 精确法向:它不是从网格三角面片插值得到的近似法向,而是直接从 CAD B-rep 几何求值。因此即使在非常粗糙的网格上,法向量也是精确的。在压力载荷施加时应优先使用 getNormal() 而非从网格计算法向。

  5. 主曲率与主方向的正交性getPrincipalCurvatures 返回的两个主方向在连续几何上严格正交。但在退化曲面(如球面极点)附近,方向可能不稳定——这是微分几何的固有特性(脐点处所有方向都是主方向)。

  6. getDerivative 返回的是参数导数,非弧长导数:切向量 C'(t) 的模长通常不等于 1。需要单位切向量时,应除以 |C'(t)|。同理,Frenet 曲率需要额外计算 |C' x C''| / |C'|^3

  7. reparametrizeOnSurface 的可用性:该功能仅对部分 CAD 内核可用(特别是 OpenCASCADE)。使用前应检查返回值是否为空。对于不可用的实体,备选方案是通过投影或最近点搜索近似替代。

  8. 物理坐标反向映射的收敛性getParametrization() 通过 Newton 迭代求解非线性反问题。如果提供的物理坐标离曲面太远(超出 CAD 容差),迭代可能失败或返回无意义值。建议在调用前确保坐标确实在面上(可通过 isInside() 预检)。

  9. 网格节点的参数坐标自动存储:当通过 mesh::generate() 生成网格后,每个网格节点自动存储了其在所属几何实体上的参数坐标。这些坐标可通过 mesh::getNodes() 获取,无需再次调用 getParametrization()(除非节点被移动或修改)。

  10. 后处理视图的性能:当网格节点数很大(如数百万)时,将所有法向量打包为单个 list-based view 可能耗尽内存。可考虑分批创建视图、使用 model-based view(依赖网格插值),或仅采样显示(每隔 N 个点取一个)。


25.7 本章小结

  • 曲面参数化 (u,v) <-> (x,y,z) 是 CAD 内核的核心功能,Gmsh 通过 getParametrization()(物理->参数)和 getValue()(参数->物理)实现双向映射,getParametrizationBounds() 提供定义域边界。
  • getNormal(tag, ...) 以面标签(int)为第一参数,返回 CAD 精确单位法向量,与 getCurvature(dim, tag, ...)(dim,tag) 签名不一致,是常见错误来源。
  • getPrincipalCurvatures() 在一次调用中返回两个主曲率值及其对应的 3D 主方向向量,是自适应网格细化最重要的几何判据来源。
  • 曲线求值 API 族(getValue/getDerivative/getSecondDerivative)分别提供位置、切向量和二阶导,可用于 Frenet 曲率计算和局部坐标系建立。
  • reparametrizeOnSurface() 实现曲线在曲面上的重参数化,多用于验证边界一致性或进行边界约束的坐标转换。
  • 所有几何求值 API 都是直接从 CAD B-rep 求值,精度独立于网格密度——这是 Gmsh 作为几何感知(geometry-aware)前处理工具的核心优势。