Files
ViewDesignEngine/include/vde/mesh/fea_mesh.h
T
茂之钳 bb0029234f
CI / Build & Test (push) Failing after 31s
CI / Release Build (push) Failing after 38s
Build & Test / build-and-test (push) Has been cancelled
Build & Test / python-bindings (push) Has been cancelled
feat(v12): generative surfaces (Sweep/Loft/Net) + curve tools + surface analysis
M1 — Generative Surfaces:
- sweep_surface: Explicit/Spine/TwoGuides modes
- loft_surface: multi-section interpolation + guide + tangent constraints
- net_surface: bidirectional curve grid → NURBS

M3 — Curve Tools + Analysis:
- curve_tools: project_curve, parallel_curve, connect_curve(G1-G3), helix, isoparametric
- surface_analysis: inflection_lines, checker_mapping, surface_checker_report (A/B/C/D)
- 14/14 tests passed

Fixed: constraint_solver.h duplicate declaration, fea_mesh.h comment syntax

Pending: M2 surface editing (retrying)
2026-07-27 09:08:01 +08:00

545 lines
20 KiB
C++
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
#pragma once
/**
* @file fea_mesh.h
* @brief 高质量FEA有限元网格生成模块
*
* 提供四面体、边界层棱柱、六面体网格生成、自适应细化、
* 网格质量评估及CAE格式(Abaqus/Ansys/Nastran)导出。
*
* ## 功能概览
* - tetrahedral_mesh: Delaunay 3D + 质量优化四面体网格
* - boundary_layer_mesh: 沿壁面法向拉伸棱柱边界层
* - hexahedral_mesh: sweep/extrude 六面体网格
* - adaptive_refinement: 基于误差估计的局部细化
* - fea_quality_report: skewness/aspect_ratio/jacobian/orthogonality 报告
* - export_abaqus / export_ansys / export_nastran: CAE 格式导出
*
* @ingroup mesh
*/
#include "vde/mesh/halfedge_mesh.h"
#include "vde/mesh/delaunay_3d.h"
#include "vde/brep/brep.h"
#include <vector>
#include <array>
#include <string>
#include <functional>
#include <cstdio>
namespace vde::mesh {
using core::Point3D;
using core::Vector3D;
// ════════════════════════════════════════════════════════
// FEA 单元类型
// ════════════════════════════════════════════════════════
/**
* @brief FEA 单元拓扑类型
* @ingroup mesh
*/
enum class FEAElementType {
Tet4, ///< 4节点线性四面体 (C3D4)
Tet10, ///< 10节点二次四面体 (C3D10)
Hex8, ///< 8节点线性六面体 (C3D8)
Hex20, ///< 20节点二次六面体 (C3D20)
Wedge6, ///< 6节点线性三棱柱/楔形 (C3D6)
Wedge15, ///< 15节点二次三棱柱 (C3D15)
};
/// 获取单元类型对应的每单元节点数
inline int nodes_per_element(FEAElementType t) {
switch (t) {
case FEAElementType::Tet4: return 4;
case FEAElementType::Tet10: return 10;
case FEAElementType::Hex8: return 8;
case FEAElementType::Hex20: return 20;
case FEAElementType::Wedge6: return 6;
case FEAElementType::Wedge15:return 15;
}
return 0;
}
/// 获取单元类型名称
inline const char* element_type_name(FEAElementType t) {
switch (t) {
case FEAElementType::Tet4: return "Tet4";
case FEAElementType::Tet10: return "Tet10";
case FEAElementType::Hex8: return "Hex8";
case FEAElementType::Hex20: return "Hex20";
case FEAElementType::Wedge6: return "Wedge6";
case FEAElementType::Wedge15:return "Wedge15";
}
return "Unknown";
}
// ════════════════════════════════════════════════════════
// FEAMesh 数据结构
// ════════════════════════════════════════════════════════
/**
* @brief FEA有限元网格
*
* 存储顶点、单元连通性、边界面及单元类型信息。
* 支持Tet4/Tet10/Hex8/Wedge6等常见FEA单元。
*
* @ingroup mesh
*/
struct FEAMesh {
std::vector<Point3D> vertices; ///< 顶点坐标数组
std::vector<std::vector<int>> elements; ///< 各单元顶点索引(长度 = nodes_per_element
FEAElementType element_type = FEAElementType::Tet4; ///< 单元类型
std::vector<int> material_ids; ///< 各单元材料ID(可选,默认0)
std::vector<std::vector<int>> boundary_faces;///< 边界面(用于施加BC),顶点索引数组
std::vector<int> boundary_face_tags; ///< 边界面标签(与boundary_faces一一对应)
/** @brief 单元总数 */
[[nodiscard]] size_t num_elements() const { return elements.size(); }
/** @brief 顶点总数 */
[[nodiscard]] size_t num_vertices() const { return vertices.size(); }
/** @brief 边界面总数 */
[[nodiscard]] size_t num_boundary_faces() const { return boundary_faces.size(); }
/** @brief 每单元节点数 */
[[nodiscard]] int npe() const { return nodes_per_element(element_type); }
};
// ════════════════════════════════════════════════════════
// 网格生成参数
// ════════════════════════════════════════════════════════
/**
* @brief 四面体网格生成参数
* @ingroup mesh
*/
struct TetMeshParams {
double min_angle = 15.0; ///< 最小二面角(度),太小导致刚度矩阵病态
double max_aspect_ratio = 5.0; ///< 最大长宽比(外接球半径/最短边)
double max_size = 0.0; ///< 最大单元尺寸,0 = 自动(包围盒对角线的 5%)
double min_size = 0.0; ///< 最小单元尺寸,0 = max_size/10
int quality_iterations = 3; ///< 质量优化迭代次数(Laplacian光顺+边翻转)
bool keep_surface = true; ///< 是否保持原表面顶点位置
};
/**
* @brief 边界层网格参数
* @ingroup mesh
*/
struct BLPParams {
double first_cell_height = 0.001; ///< 第一层棱柱高度(壁面法向距离)
double growth_rate = 1.2; ///< 相邻层高度增长率(>1.0
int num_layers = 5; ///< 棱柱层数
double total_thickness = 0.0; ///< 边界层总厚度,0 = 自动(first * (growth^n - 1)/(growth - 1)
bool smooth_transition = true; ///< 是否平滑过渡到内部网格
};
/**
* @brief 六面体网格参数
* @ingroup mesh
*/
struct HexMeshParams {
int sweep_layers = 1; ///< sweep 方向层数
double sweep_distance = 0.0; ///< sweep 总距离,0 = 自动(包围盒最小维度)
bool use_quad_dominant = true; ///< 源面是否使用四边形主导网格
double cell_size = 0.0; ///< 目标单元尺寸,0 = 自动
};
// ════════════════════════════════════════════════════════
// 网格质量评估
// ════════════════════════════════════════════════════════
/**
* @brief FEA网格质量报告
*
* 包含 skewness(偏斜度)、aspect_ratio(长宽比)、
* jacobian(雅可比)、orthogonality(正交性)四项
* 核心FEA网格质量指标。
*
* @ingroup mesh
*/
struct FEAQualityReport {
/// 偏斜度 (skewness): 实际单元与理想单元的偏差,[0,1],0=最优
double min_skewness = 0.0;
double max_skewness = 0.0;
double avg_skewness = 0.0;
/// 长宽比 (aspect ratio): 最长边/最短边,(Fluent: 1=最佳)
double min_aspect_ratio = 0.0;
double max_aspect_ratio = 0.0;
double avg_aspect_ratio = 0.0;
/// 归一化雅可比 (scaled jacobian): [0,1]1=正单元,0=退化,<0=反转
double min_jacobian = 0.0;
double max_jacobian = 0.0;
double avg_jacobian = 0.0;
/// 正交性 (orthogonality): 面法向与质心连线的夹角余弦,[0,1],1=完美正交
double min_orthogonality = 0.0;
double max_orthogonality = 0.0;
double avg_orthogonality = 0.0;
size_t total_elements = 0; ///< 单元总数
size_t degenerate_elements = 0; ///< 退化单元数 (Jacobian < 0.01)
size_t inverted_elements = 0; ///< 反转单元数 (Jacobian < 0)
/// 分级统计 (A=优秀 B=良好 C=合格 D=差 F=不合格,按Jacobian)
size_t grade_A = 0, grade_B = 0, grade_C = 0, grade_D = 0, grade_F = 0;
};
// ════════════════════════════════════════════════════════
// 自适应细化
// ════════════════════════════════════════════════════════
/**
* @brief 误差估计函数类型
*
* 对给定单元返回估计误差值(≥0),值越大表示该单元越需要细化。
*
* @param element_idx 单元索引
* @param mesh 所属网格
* @return 误差估计值(非负)
*/
using ErrorEstimator = std::function<double(int element_idx, const FEAMesh& mesh)>;
/**
* @brief 自适应细化参数
* @ingroup mesh
*/
struct RefinementParams {
double error_threshold = 0.5; ///< 相对误差阈值(高于此值的单元被细化)
double refinement_fraction = 0.2; ///< 每次细化单元比例(0~1
int max_iterations = 3; ///< 最大细化迭代次数
bool smooth_after_refine = true; ///< 细化后是否光顺
};
// ════════════════════════════════════════════════════════
// 网格生成 API
// ════════════════════════════════════════════════════════
/**
* @brief 生成四面体FEA网格(Delaunay 3D + 质量优化)
*
* 对 BrepModel 进行 tessellation 获取表面三角网格,
* 采样内部点后进行 Delaunay 3D 四面体化,
* 再进行 Laplacian 光顺 + 边翻转优化质量。
*
* @param body B-Rep实体模型
* @param params 四面体化参数(最小角、长宽比约束等)
* @return FEAMesh 四面体网格
*
* @code{.cpp}
* auto box = vde::brep::make_box(1, 1, 1);
* TetMeshParams p;
* p.min_angle = 18.0;
* p.max_aspect_ratio = 4.0;
* auto mesh = tetrahedral_mesh(box, p);
* @endcode
* @ingroup mesh
*/
FEAMesh tetrahedral_mesh(const brep::BrepModel& body,
const TetMeshParams& params = {});
/**
* @brief 生成边界层棱柱网格
*
* 沿B-Rep模型壁面法向拉伸多层棱柱单元,
* 用于CFD边界层捕捉(y+控制)。
* 先对表面进行三角化,然后每层沿法向偏移。
*
* @param body B-Rep实体模型
* @param params 边界层参数(首层高度、增长率、层数)
* @return FEAMesh 棱柱网格 (Wedge6)
*
* @code{.cpp}
* BLPParams p;
* p.first_cell_height = 0.0001;
* p.growth_rate = 1.15;
* p.num_layers = 10;
* auto bl = boundary_layer_mesh(wing, p);
* @endcode
* @ingroup mesh
*/
FEAMesh boundary_layer_mesh(const brep::BrepModel& body,
const BLPParams& params = {});
/**
* @brief 生成六面体网格(sweep/extrude策略)
*
* 将B-Rep模型的一个面(源面)四边形化,
* 然后沿法向 sweep 生成六面体层。
* 适用于拉伸体/棱柱体等 sweepable 几何。
*
* @param body B-Rep实体模型
* @param params 六面体化参数
* @return FEAMesh 六面体网格 (Hex8)
*
* @code{.cpp}
* auto cylinder = vde::brep::make_cylinder(1, 3);
* HexMeshParams p;
* p.sweep_layers = 5;
* auto mesh = hexahedral_mesh(cylinder, p);
* @endcode
* @ingroup mesh
*/
FEAMesh hexahedral_mesh(const brep::BrepModel& body,
const HexMeshParams& params = {});
// ════════════════════════════════════════════════════════
// 高级六面体网格生成 API
// ════════════════════════════════════════════════════════
/// 六面体质量优化选项
struct HexOptimizationParams {
double min_scaled_jacobian = 0.3; ///< 最小可接受雅可比
double max_aspect_ratio = 10.0; ///< 最大长宽比
int max_iterations = 50; ///< 最大优化迭代次数
double convergence_tol = 1e-4; ///< 收敛容差
bool untangle = true; ///< 是否翻转修复(untangle inverted elements
bool laplacian_smooth = true; ///< 是否 Laplacian 光顺
bool optimization_based_smooth = true; ///< 是否基于优化的光顺
};
/// 六面体质量优化结果
struct HexOptimizationResult {
FEAMesh optimized_mesh; ///< 优化后的网格
FEAQualityReport initial_quality; ///< 初始质量报告
FEAQualityReport final_quality; ///< 最终质量报告
int iterations = 0; ///< 实际迭代次数
int inverted_fixed = 0; ///< 修复的反转单元数
int degenerate_fixed = 0; ///< 修复的退化单元数
bool converged = false; ///< 是否收敛
};
/// 面标识结构
struct FaceRegion {
std::vector<int> face_indices; ///< 面的顶点索引列表
std::string name; ///< 面名称
};
/// 多块分区定义
struct BlockRegion {
std::vector<Point3D> corner_vertices; ///< 8个角顶点(i,j,k 顺序)
int n_i = 2; ///< i方向分段数
int n_j = 2; ///< j方向分段数
int n_k = 2; ///< k方向分段数
};
/**
* @brief 映射法六面体网格生成
*
* 给定 B-Rep 实体的源面和目标面,使用映射法生成结构化六面体网格。
* 核心思路:将源面四边形网格通过坐标映射变换到目标面,
* 并在中间层线性插值,形成规则的六面体网格。
*
* 算法:
* 1. 在源面生成四边形网格(重心映射或栅格映射)
* 2. 在目标面生成对应的四边形网格
* 3. 沿 sweep 方向线性插值生成中间层
* 4. 连接各层为六面体单元
*
* @param body B-Rep实体模型
* @param source_face 源面标识
* @param target_face 目标面标识
* @param layers 层数(默认 4
* @return FEAMesh 六面体网格 (Hex8)
*
* @note 对标 Ansys ICEM CFD Blocking / Pointwise
* @ingroup mesh
*/
[[nodiscard]] FEAMesh mapped_hex_mesh(
const brep::BrepModel& body,
const FaceRegion& source_face,
const FaceRegion& target_face,
int layers = 4);
/**
* @brief 子映射法六面体网格生成
*
* 自动将复杂几何体分区为多个可映射子区域,
* 然后在每个子区域内用映射法生成六面体网格,
* 最后缝合各子区域网格。
*
* 算法:
* 1. 分析几何体的拓扑特征(角点、棱边)
* 2. 自动分区为可扫掠(sweepable)子区域
* 3. 对每个子区域执行映射法六面体化
* 4. 合并子区域网格,保证界面节点一致
*
* @param body B-Rep实体模型
* @return FEAMesh 六面体网格 (Hex8)
*
* @note 对标 Ansys Workbench MultiZone / CUBIT Submap
* @ingroup mesh
*/
[[nodiscard]] FEAMesh submapped_hex_mesh(
const brep::BrepModel& body);
/**
* @brief 多块结构化六面体网格生成
*
* 对几何体进行显式多块分解,每个块独立生成结构化六面体网格,
* 确保块间界面节点一一对应。
*
* 算法:
* 1. 根据用户定义或自动识别的块分区
* 2. 每个块内使用 transfinite interpolation (TFI) 生成节点
* 3. 块间共享面节点,保证 conformal 连接
* 4. 组装块为整体网格
*
* @param body B-Rep实体模型
* @param blocks 块分区定义列表
* @return FEAMesh 六面体网格 (Hex8)
*
* @note 对标 Ansys ICEM CFD Hexa Blocking / TrueGrid
* @ingroup mesh
*/
[[nodiscard]] FEAMesh multi_block_hex_mesh(
const brep::BrepModel& body,
const std::vector<BlockRegion>& blocks);
/**
* @brief 六面体网格质量优化与翻转修复
*
* 对六面体网格进行质量优化,包括反转变形修复(untangling)、
* Laplacian 光顺和基于优化的光顺。
*
* 算法:
* - untangling: 检测负雅可比单元,通过解局部优化问题翻转节点
* - Laplacian 光顺: 将内部节点移到邻域质心
* - 优化光顺: 最小化目标函数(扭曲度+体积变化)
*
* @param mesh 输入 FEA 六面体网格
* @param params 优化参数
* @return HexOptimizationResult 优化结果(含优化前后质量报告)
*
* @note 对标 CUBIT MeshGems / Pointwise T-Rex Quality
* @ingroup mesh
*/
[[nodiscard]] HexOptimizationResult hex_quality_optimization(
const FEAMesh& mesh,
const HexOptimizationParams& params = {});
// ════════════════════════════════════════════════════════
// 自适应细化
/**
* 基于误差估计函数对高误差区域进行局部细分。
* 支持边中点分裂细化(适用于Tet4→4×Tet4子单元)。
*
* @param mesh 输入FEA网格
* @param estimator 误差估计函数
* @param params 细化参数
* @return FEAMesh 细化后的网格
*
* @code{.cpp}
* auto estimator = [](int ei, const FEAMesh& m) -> double {
* // 基于应力梯度的误差估计
* return stress_gradient_norm(ei, m);
* };
* auto refined = adaptive_refinement(mesh, estimator);
* @endcode
* @ingroup mesh
*/
FEAMesh adaptive_refinement(const FEAMesh& mesh,
const ErrorEstimator& estimator,
const RefinementParams& params = {});
// ════════════════════════════════════════════════════════
// 质量评估 API
// ════════════════════════════════════════════════════════
/**
* @brief 评估FEA网格质量
*
* 计算四项核心质量指标:
* - skewness: 偏斜度(实际vs理想单元的偏差)
* - aspect_ratio: 长宽比
* - jacobian: 归一化雅可比行列式
* - orthogonality: 正交性(面法向与质心连线夹角)
*
* 自动根据 element_type 选择对应算法(Tet/Hex/Wedge)。
*
* @param mesh FEA网格
* @return FEAQualityReport 质量报告
*
* @ingroup mesh
*/
FEAQualityReport fea_quality_report(const FEAMesh& mesh);
// ── 单项质量指标 ──
/**
* @brief 计算单元的偏斜度
* @param verts 单元顶点数组
* @param type 单元类型
* @return skewness [0,1]0=理想单元
*/
double element_skewness(const std::vector<Point3D>& verts, FEAElementType type);
/**
* @brief 计算单元的长宽比
* @param verts 单元顶点数组
* @param type 单元类型
* @return aspect_ratio ≥ 11=理想
*/
double element_aspect_ratio(const std::vector<Point3D>& verts, FEAElementType type);
/**
* @brief 计算单元的归一化雅可比
* @param verts 单元顶点数组
* @param type 单元类型
* @return scaled_jacobian [0,1]1=正四面体/正方体,0=退化
*/
double element_scaled_jacobian(const std::vector<Point3D>& verts, FEAElementType type);
/**
* @brief 计算单元的正交性
* @param verts 单元顶点数组
* @param type 单元类型
* @return orthogonality [0,1]1=完全正交(面法向与对面质心连线平行)
*/
double element_orthogonality(const std::vector<Point3D>& verts, FEAElementType type);
// ════════════════════════════════════════════════════════
// CAE 格式导出
// ════════════════════════════════════════════════════════
/**
* @brief 导出 Abaqus 输入文件 (.inp)
*
* 生成包含 *NODE 和 *ELEMENT 卡的 Abaqus 输入文件。
* 支持 Tet4 (C3D4), Tet10 (C3D10), Hex8 (C3D8), Wedge6 (C3D6)。
*
* @param mesh FEA网格
* @param filepath 输出文件路径
* @return 成功返回 true
*
* @ingroup mesh
*/
bool export_abaqus(const FEAMesh& mesh, const std::string& filepath);
/**
* @brief 导出 Ansys CDB 文件 (.cdb)
*
* 生成包含 NBLOCK 和 EBLOCK 的 Ansys CDB 格式文件。
*
* @param mesh FEA网格
* @param filepath 输出文件路径
* @return 成功返回 true
*
* @ingroup mesh
*/
bool export_ansys(const FEAMesh& mesh, const std::string& filepath);
/**
* @brief 导出 Nastran 批量数据文件 (.bdf/.dat)
*
* 生成包含 GRID 和 CTETRA/CHEXA/CPENTA 卡的 Nastran 输入文件。
*
* @param mesh FEA网格
* @param filepath 输出文件路径
* @return 成功返回 true
*
* @ingroup mesh
*/
bool export_nastran(const FEAMesh& mesh, const std::string& filepath);
} // namespace vde::mesh