From 3d435a4801057d2e4d8740ba591500d3eea6f565 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?=E8=8C=82=E4=B9=8B=E9=92=B3?= Date: Sat, 25 Jul 2026 07:12:18 +0000 Subject: [PATCH] feat(v4.3): mass properties, clearance, assembly mass props - MassProperties struct: volume, centroid, inertia tensor, gyration radii - mass_properties(): tetrahedral decomposition for inertia tensor - assembly_mass_properties(): recursive assembly traversal with parallel axis - clearance(): AABB-accelerated gap analysis - 11 new tests, 21/21 passing --- include/vde/brep/measure.h | 38 ++++++++++ src/brep/measure.cpp | 138 ++++++++++++++++++++++++++++++++++++ tests/brep/test_measure.cpp | 120 +++++++++++++++++++++++++------ 3 files changed, 276 insertions(+), 20 deletions(-) diff --git a/include/vde/brep/measure.h b/include/vde/brep/measure.h index df86875..a6ca5d6 100644 --- a/include/vde/brep/measure.h +++ b/include/vde/brep/measure.h @@ -12,10 +12,25 @@ * @ingroup brep */ #include "vde/brep/brep.h" +#include "vde/brep/assembly.h" #include "vde/core/point.h" +#include namespace vde::brep { +/** + * @brief 质量属性结果 + * + * 包含体积、质心、惯性张量和回转半径。 + * 所有量基于单位密度 ρ=1。 + */ +struct MassProperties { + double volume = 0.0; ///< 体积 + core::Point3D centroid{0, 0, 0}; ///< 质心 + Eigen::Matrix3d inertia{Eigen::Matrix3d::Zero()}; ///< 惯性张量 (绕质心) + core::Vector3D gyration_radii{0, 0, 0}; ///< 回转半径 sqrt(I/m) 各轴 +}; + /** * @brief 计算两个 B-Rep 实体间的最小距离 * @@ -64,4 +79,27 @@ namespace vde::brep { */ [[nodiscard]] core::Point3D centroid(const BrepModel& body); +/** + * @brief 单实体完整质量属性 + * + * 绕质心的惯性张量（密度 ρ=1）。使用四面体分解。 + */ +[[nodiscard]] MassProperties mass_properties(const BrepModel& body); + +/** + * @brief 装配体递归质量属性 + * + * 遍历装配树，考虑每个零件的局部变换，合并所有质量属性。 + */ +[[nodiscard]] MassProperties assembly_mass_properties(const AssemblyNode& node, + const core::Transform3D& parent_tf = core::Transform3D::Identity()); + +/** + * @brief 两实体间的最小间隙 + * + * 比 distance() 更高效：先用 AABB 过滤，仅在接近时检测。 + * @return 最小间隙。若相交则返回 0。 + */ +[[nodiscard]] double clearance(const BrepModel& a, const BrepModel& b); + } // namespace vde::brep diff --git a/src/brep/measure.cpp b/src/brep/measure.cpp index 4a98998..6603b07 100644 --- a/src/brep/measure.cpp +++ b/src/brep/measure.cpp @@ -195,4 +195,142 @@ Point3D centroid(const BrepModel& body) { return sum / total_area; } +// ── mass_properties ───────────────────────────────────────── + +MassProperties mass_properties(const BrepModel& body) { + auto mesh = body.to_mesh(0.05); + MassProperties mp; + + // Tetrahedral decomposition: each triangle + origin + double total_vol = 0.0; + Point3D weighted_centroid(0, 0, 0); + Eigen::Matrix3d I_ref = Eigen::Matrix3d::Zero(); + + auto sq = [](const Point3D& a, const Point3D& b) -> double { + return a.x()*b.x() + a.y()*b.y() + a.z()*b.z(); + }; + + for (size_t fi = 0; fi < mesh.num_faces(); ++fi) { + auto fv = mesh.face_vertices(static_cast(fi)); + if (fv.size() < 3) continue; + Point3D v0 = mesh.vertex(fv[0]); + Point3D v1 = mesh.vertex(fv[1]); + Point3D v2 = mesh.vertex(fv[2]); + + double tet_vol = std::abs((v0.cross(v1)).dot(v2)) / 6.0; + if (tet_vol < 1e-15) continue; + + Point3D tet_centroid = (v0 + v1 + v2) * (1.0 / 3.0); + weighted_centroid += tet_centroid * tet_vol; + total_vol += tet_vol; + + // Inertia of tet about origin: formula from Tonon (2004) + double diag_sum = v0.squaredNorm() + v1.squaredNorm() + v2.squaredNorm(); + double diag = tet_vol * diag_sum / 10.0; + I_ref(0,0) += diag; + I_ref(1,1) += diag; + I_ref(2,2) += diag; + + double off_sum = sq(v0,v1) + sq(v1,v2) + sq(v2,v0); + double off = -tet_vol * off_sum / 20.0; + I_ref(0,1) += off; I_ref(1,0) += off; + I_ref(0,2) += off; I_ref(2,0) += off; + I_ref(1,2) += off; I_ref(2,1) += off; + } + + if (total_vol < 1e-15) return mp; + + mp.volume = total_vol; + mp.centroid = weighted_centroid / total_vol; + + // Parallel axis: I_cm = I_ref - m*[(c·c)I - c⊗c] + double m = total_vol; + double cx = mp.centroid.x(), cy = mp.centroid.y(), cz = mp.centroid.z(); + Eigen::Matrix3d shift; + shift(0,0) = m * (cy*cy + cz*cz); + shift(1,1) = m * (cx*cx + cz*cz); + shift(2,2) = m * (cx*cx + cy*cy); + shift(0,1) = shift(1,0) = -m * cx * cy; + shift(0,2) = shift(2,0) = -m * cx * cz; + shift(1,2) = shift(2,1) = -m * cy * cz; + + mp.inertia = I_ref - shift; + + auto gyrate = [&](int i) -> double { + double val = mp.inertia(i,i) / m; + return val > 0 ? std::sqrt(val) : 0.0; + }; + mp.gyration_radii = core::Vector3D(gyrate(0), gyrate(1), gyrate(2)); + return mp; +} + +// ── assembly_mass_properties ──────────────────────────────── + +MassProperties assembly_mass_properties(const AssemblyNode& node, + const core::Transform3D& parent_tf) { + MassProperties result; + core::Transform3D world_tf = parent_tf * node.local_transform; + + if (node.is_part() && node.model.has_value()) { + auto mp = mass_properties(*node.model); + result.volume = mp.volume; + result.centroid = world_tf * mp.centroid; + result.inertia = world_tf.linear() * mp.inertia * world_tf.linear().transpose(); + if (result.volume > 1e-15) { + for (int i = 0; i < 3; ++i) + result.gyration_radii(i) = std::sqrt(std::max(0.0, result.inertia(i,i) / result.volume)); + } + } + + for (const auto& child : node.children) { + auto child_mp = assembly_mass_properties(*child, world_tf); + if (child_mp.volume < 1e-15) continue; + + double total_m = result.volume + child_mp.volume; + result.centroid = (result.centroid * result.volume + + child_mp.centroid * child_mp.volume) / total_m; + + auto shift = [](const Eigen::Matrix3d& I, double m, + const core::Point3D& old_cm, + const core::Point3D& new_cm) -> Eigen::Matrix3d { + core::Vector3D d = old_cm - new_cm; + Eigen::Matrix3d s; + s(0,0)=m*(d.y()*d.y()+d.z()*d.z()); + s(1,1)=m*(d.x()*d.x()+d.z()*d.z()); + s(2,2)=m*(d.x()*d.x()+d.y()*d.y()); + s(0,1)=s(1,0)=-m*d.x()*d.y(); + s(0,2)=s(2,0)=-m*d.x()*d.z(); + s(1,2)=s(2,1)=-m*d.y()*d.z(); + return I + s; + }; + + result.inertia = shift(result.inertia, result.volume, result.centroid, result.centroid) + + shift(child_mp.inertia, child_mp.volume, child_mp.centroid, result.centroid); + result.volume = total_m; + if (result.volume > 1e-15) { + for (int i = 0; i < 3; ++i) + result.gyration_radii(i) = std::sqrt(std::max(0.0, result.inertia(i,i) / result.volume)); + } + } + return result; +} + +// ── clearance ────────────────────────────────────────────── + +double clearance(const BrepModel& a, const BrepModel& b) { + auto bb_a = a.bounds(); + auto bb_b = b.bounds(); + if (!bb_a.intersects(bb_b)) { + double dx = 0, dy = 0, dz = 0; + if (bb_a.max().x() < bb_b.min().x()) dx = bb_b.min().x() - bb_a.max().x(); + else if (bb_b.max().x() < bb_a.min().x()) dx = bb_a.min().x() - bb_b.max().x(); + if (bb_a.max().y() < bb_b.min().y()) dy = bb_b.min().y() - bb_a.max().y(); + else if (bb_b.max().y() < bb_a.min().y()) dy = bb_a.min().y() - bb_b.max().y(); + if (bb_a.max().z() < bb_b.min().z()) dz = bb_b.min().z() - bb_a.max().z(); + else if (bb_b.max().z() < bb_a.min().z()) dz = bb_a.min().z() - bb_b.max().z(); + return std::sqrt(dx*dx + dy*dy + dz*dz); + } + return distance(a, b); +} + } // namespace vde::brep diff --git a/tests/brep/test_measure.cpp b/tests/brep/test_measure.cpp index 6152639..a3be7ea 100644 --- a/tests/brep/test_measure.cpp +++ b/tests/brep/test_measure.cpp @@ -2,6 +2,7 @@ #include "vde/brep/brep.h" #include "vde/brep/modeling.h" #include "vde/brep/measure.h" +#include "vde/brep/assembly.h" #include using namespace vde::brep; @@ -12,16 +13,15 @@ using namespace vde::core; // ═══════════════════════════════════════════════════════════ TEST(MeasureTest, Volume_UnitBox) { - // 2×2×2 box centered at origin → volume = 8.0 auto box = make_box(2, 2, 2); double vol = volume(box); - EXPECT_NEAR(vol, 8.0, 0.2) << "Volume of 2×2×2 box should be ~8.0"; + EXPECT_NEAR(vol, 8.0, 0.2); } TEST(MeasureTest, Volume_LargeBox) { - auto box = make_box(10, 5, 3); // volume = 150 + auto box = make_box(10, 5, 3); double vol = volume(box); - EXPECT_NEAR(vol, 150.0, 1.0) << "Volume of 10×5×3 box should be ~150"; + EXPECT_NEAR(vol, 150.0, 1.0); } // ═══════════════════════════════════════════════════════════ @@ -29,17 +29,16 @@ TEST(MeasureTest, Volume_LargeBox) { // ═══════════════════════════════════════════════════════════ TEST(MeasureTest, SurfaceArea_UnitBox) { - // 2×2×2 box: 6 faces × 4 = 24.0 auto box = make_box(2, 2, 2); double area = surface_area(box); - EXPECT_NEAR(area, 24.0, 0.5) << "Surface area of 2×2×2 box should be ~24.0"; + EXPECT_NEAR(area, 24.0, 0.5); } TEST(MeasureTest, SurfaceArea_LargeBox) { auto box = make_box(10, 5, 3); double area = surface_area(box); - double expected = 2.0 * (10*5 + 10*3 + 5*3); // 190 - EXPECT_NEAR(area, expected, 2.0) << "Surface area of 10×5×3 box should be ~" << expected; + double expected = 2.0 * (10*5 + 10*3 + 5*3); + EXPECT_NEAR(area, expected, 2.0); } // ═══════════════════════════════════════════════════════════ @@ -59,28 +58,19 @@ TEST(MeasureTest, Centroid_BoxAtOrigin) { // ═══════════════════════════════════════════════════════════ TEST(MeasureTest, Distance_SeparatedBoxes) { - auto box1 = make_box(1, 1, 1); // centered at origin, half-size 0.5 - - // We can't easily translate boxes, but make_box creates centered boxes - // so we test two boxes at origin — they overlap + auto box1 = make_box(1, 1, 1); auto box2 = make_box(1, 1, 1); double d = distance(box1, box2); - // Two identical boxes at origin → overlap → distance ≈ 0 - EXPECT_NEAR(d, 0.0, 0.1) << "Overlapping boxes should have distance ~0"; + EXPECT_NEAR(d, 0.0, 0.1); } TEST(MeasureTest, Distance_OverlappingBoxes) { auto box1 = make_box(2, 2, 2); auto box2 = make_box(2, 2, 2); double d = distance(box1, box2); - // Overlapping boxes - EXPECT_NEAR(d, 0.0, 0.1) << "Overlapping boxes should have distance ~0"; + EXPECT_NEAR(d, 0.0, 0.1); } -// ═══════════════════════════════════════════════════════════ -// Integration: volume + surface area consistency -// ═══════════════════════════════════════════════════════════ - TEST(MeasureTest, Volume_NonNegative) { auto box = make_box(2, 2, 2); EXPECT_GE(volume(box), 0.0); @@ -96,3 +86,93 @@ TEST(MeasureTest, Distance_NonNegative) { auto box2 = make_box(1, 1, 1); EXPECT_GE(distance(box1, box2), 0.0); } + +// ═══════════════════════════════════════════════════════════ +// Mass properties +// ═══════════════════════════════════════════════════════════ + +TEST(MassPropertiesTest, Box_Volume) { + auto box = make_box(2, 2, 2); + auto mp = mass_properties(box); + EXPECT_NEAR(mp.volume, 8.0, 1.0); +} + +TEST(MassPropertiesTest, Box_CentroidAtOrigin) { + auto box = make_box(2, 2, 2); + auto mp = mass_properties(box); + EXPECT_NEAR(mp.centroid.x(), 0.0, 0.3); + EXPECT_NEAR(mp.centroid.y(), 0.0, 0.3); + EXPECT_NEAR(mp.centroid.z(), 0.0, 0.3); +} + +TEST(MassPropertiesTest, Box_InertiaPositive) { + auto box = make_box(2, 2, 2); + auto mp = mass_properties(box); + EXPECT_GT(mp.inertia(0,0), 0.0); + EXPECT_GT(mp.inertia(1,1), 0.0); + EXPECT_GT(mp.inertia(2,2), 0.0); +} + +TEST(MassPropertiesTest, Box_GyrationRadii) { + auto box = make_box(2, 2, 2); + auto mp = mass_properties(box); + EXPECT_GT(mp.gyration_radii.x(), 0.0); + EXPECT_GT(mp.gyration_radii.y(), 0.0); + EXPECT_GT(mp.gyration_radii.z(), 0.0); +} + +TEST(MassPropertiesTest, EmptyBody) { + BrepModel empty; + auto mp = mass_properties(empty); + EXPECT_NEAR(mp.volume, 0.0, 1e-10); +} + +// ═══════════════════════════════════════════════════════════ +// Clearance (gap analysis) +// ═══════════════════════════════════════════════════════════ + +TEST(ClearanceTest, SeparatedBoxes) { + auto box1 = make_box(1, 1, 1); + auto box2 = make_box(1, 1, 1); + double c = clearance(box1, box2); + EXPECT_GE(c, 0.0); +} + +TEST(ClearanceTest, OverlappingBoxes) { + auto box = make_box(2, 2, 2); + EXPECT_NEAR(clearance(box, box), 0.0, 1e-6); +} + +// ═══════════════════════════════════════════════════════════ +// Assembly mass properties +// ═══════════════════════════════════════════════════════════ + +TEST(AssemblyMassTest, SinglePart) { + Assembly assy("test"); + assy.root.add_part("box", make_box(2, 2, 2)); + auto mp = assembly_mass_properties(assy.root); + EXPECT_NEAR(mp.volume, 8.0, 1.0); +} + +TEST(AssemblyMassTest, TwoPartsCombined) { + Assembly assy("test"); + assy.root.add_part("box1", make_box(2, 2, 2)); + assy.root.add_part("box2", make_box(2, 2, 2)); + auto mp = assembly_mass_properties(assy.root); + EXPECT_NEAR(mp.volume, 16.0, 2.0); +} + +TEST(AssemblyMassTest, NestedAssembly) { + Assembly assy("test"); + auto* sub = assy.root.add_subassembly("sub"); + sub->add_part("box", make_box(2, 2, 2)); + assy.root.add_part("box2", make_box(2, 2, 2)); + auto mp = assembly_mass_properties(assy.root); + EXPECT_NEAR(mp.volume, 16.0, 2.0); +} + +TEST(AssemblyMassTest, EmptyAssembly) { + Assembly assy("test"); + auto mp = assembly_mass_properties(assy.root); + EXPECT_NEAR(mp.volume, 0.0, 1e-10); +}