feat(v3.6): TrimmedSurface — B-Rep fundamental primitive

include/vde/brep/trimmed_surface.h

@@ -0,0 +1,54 @@
+#pragma once
+#include "vde/curves/nurbs_surface.h"
+#include "vde/curves/nurbs_curve.h"
+#include "vde/core/point.h"
+#include "vde/core/aabb.h"
+#include <optional>
+#include <vector>
+
+namespace vde::brep {
+
+/// A p-curve: 2D curve in the parameter space of a surface
+/// We reuse NurbsCurve with Z=0 — only X (u) and Y (v) are meaningful.
+using PCurve = curves::NurbsCurve;
+
+/// Trimming loop: closed contour in parameter space (u, v)
+/// An outer loop defines the boundary; inner loops define holes.
+struct TrimLoop {
+    std::vector<PCurve> p_curves;  ///< parameter-space curves forming the loop
+    bool is_outer = true;          ///< true = outer boundary, false = hole
+};
+
+/// A trimmed surface: base NURBS surface + trimming loops in parameter space
+///
+/// Point classification: a parameter (u, v) is INSIDE if it lies within
+/// the outer loop AND outside all inner loops.
+///
+/// This is the fundamental B-Rep primitive used by all industrial CAD kernels
+/// (Parasolid, ACIS, CGM).
+///
+/// @ingroup brep
+struct TrimmedSurface {
+    curves::NurbsSurface base_surface;  ///< Underlying geometry
+    std::vector<TrimLoop> loops;        ///< Trimming loops (first = outer boundary if marked)
+
+    /// Evaluate the surface at (u, v).
+    /// Returns the 3D point only if (u, v) is inside the trimmed region.
+    /// @return Point3D if inside, std::nullopt if outside
+    [[nodiscard]] std::optional<core::Point3D> evaluate(double u, double v) const;
+
+    /// Check if a parameter point is inside the trimmed region
+    [[nodiscard]] bool is_inside(double u, double v) const;
+
+    /// Get the axis-aligned bounding box of the trimmed surface
+    [[nodiscard]] core::AABB3D bounds() const;
+
+    /// Create an untrimmed surface (full parameter domain)
+    static TrimmedSurface from_surface(const curves::NurbsSurface& surf);
+
+    /// Create a rectangular trimmed surface
+    static TrimmedSurface from_rect(const curves::NurbsSurface& surf,
+                                    double u0, double u1, double v0, double v1);
+};
+
+} // namespace vde::brep

src/CMakeLists.txt

@@ -155,6 +155,7 @@ add_library(vde_brep STATIC
     brep/feature_tree.cpp
     brep/assembly_constraints.cpp
     brep/measure.cpp
+    brep/trimmed_surface.cpp
 )
 target_include_directories(vde_brep
     PUBLIC  ${CMAKE_SOURCE_DIR}/include

src/brep/trimmed_surface.cpp

@@ -0,0 +1,152 @@
+#include "vde/brep/trimmed_surface.h"
+#include <cmath>
+
+namespace vde::brep {
+
+// ─────────────────────────────────────────────────
+//  Internal: even-odd ray casting in parameter space
+// ─────────────────────────────────────────────────
+
+namespace {
+
+/// Build a polygon from p-curves by sampling each at 101 points
+std::vector<std::pair<double, double>> build_polygon(const TrimLoop& loop) {
+    std::vector<std::pair<double, double>> polygon;
+    for (const auto& pc : loop.p_curves) {
+        for (int i = 0; i <= 100; ++i) {
+            double t = static_cast<double>(i) / 100.0;
+            auto pt = pc.evaluate(t);
+            polygon.emplace_back(pt.x(), pt.y());
+        }
+    }
+    return polygon;
+}
+
+/// Even-odd rule: cast ray from (u, v) in the +u direction
+/// Returns true if the point is inside the polygon.
+bool point_in_loop(double u, double v, const TrimLoop& loop) {
+    auto polygon = build_polygon(loop);
+    if (polygon.size() < 3) return false;
+
+    int crossings = 0;
+    const size_t n = polygon.size();
+    for (size_t i = 0; i < n; ++i) {
+        const auto& [u1, v1] = polygon[i];
+        const auto& [u2, v2] = polygon[(i + 1) % n];
+
+        // Check if the ray crosses this edge
+        if ((v1 > v) != (v2 > v)) {
+            // Compute the u-coordinate of the edge at height v
+            double ue = u1 + (v - v1) / (v2 - v1) * (u2 - u1);
+            if (u < ue) crossings++;
+        }
+    }
+    return (crossings % 2) == 1;
+}
+
+} // anonymous namespace
+
+// ─────────────────────────────────────────────────
+//  is_inside
+// ─────────────────────────────────────────────────
+
+bool TrimmedSurface::is_inside(double u, double v) const {
+    if (loops.empty()) return true;  // no trimming = full surface
+
+    // Check outer loop — must be inside the outer boundary
+    const auto& outer = loops.front();
+    if (!point_in_loop(u, v, outer)) return false;
+
+    // Check inner loops — must be outside all holes
+    for (size_t i = 1; i < loops.size(); ++i) {
+        if (point_in_loop(u, v, loops[i])) return false;
+    }
+    return true;
+}
+
+// ─────────────────────────────────────────────────
+//  evaluate
+// ─────────────────────────────────────────────────
+
+std::optional<core::Point3D> TrimmedSurface::evaluate(double u, double v) const {
+    if (!is_inside(u, v)) return std::nullopt;
+    return base_surface.evaluate(u, v);
+}
+
+// ─────────────────────────────────────────────────
+//  bounds
+// ─────────────────────────────────────────────────
+
+core::AABB3D TrimmedSurface::bounds() const {
+    core::AABB3D box;
+
+    if (loops.empty()) {
+        // No trimming: sample the full parameter domain
+        const auto& ku = base_surface.knots_u();
+        const auto& kv = base_surface.knots_v();
+        int pu = base_surface.degree_u();
+        int pv = base_surface.degree_v();
+        double u_min = ku[pu];
+        double u_max = ku[ku.size() - 1 - pu];
+        double v_min = kv[pv];
+        double v_max = kv[kv.size() - 1 - pv];
+
+        const int N = 20;
+        for (int i = 0; i <= N; ++i) {
+            for (int j = 0; j <= N; ++j) {
+                double u = u_min + (u_max - u_min) * static_cast<double>(i) / N;
+                double v = v_min + (v_max - v_min) * static_cast<double>(j) / N;
+                box.expand(base_surface.evaluate(u, v));
+            }
+        }
+        return box;
+    }
+
+    // Sample the outer loop boundary and evaluate at 3D points
+    auto polygon = build_polygon(loops.front());
+    for (const auto& [u, v] : polygon) {
+        box.expand(base_surface.evaluate(u, v));
+    }
+    return box;
+}
+
+// ─────────────────────────────────────────────────
+//  Static constructors
+// ─────────────────────────────────────────────────
+
+TrimmedSurface TrimmedSurface::from_surface(const curves::NurbsSurface& surf) {
+    TrimmedSurface ts;
+    ts.base_surface = surf;
+    // No loops: untrimmed = full parameter domain
+    return ts;
+}
+
+TrimmedSurface TrimmedSurface::from_rect(const curves::NurbsSurface& surf,
+                                          double u0, double u1,
+                                          double v0, double v1) {
+    TrimmedSurface ts;
+    ts.base_surface = surf;
+
+    // Create a rectangular outer loop from 4 linear p-curves
+    // p-curves use NurbsCurve with Z=0; points are in parameter space (u, v, 0)
+    TrimLoop outer;
+    outer.is_outer = true;
+
+    auto make_line = [](double ua, double va, double ub, double vb) -> PCurve {
+        using core::Point3D;
+        std::vector<Point3D> cps = { Point3D(ua, va, 0.0), Point3D(ub, vb, 0.0) };
+        std::vector<double> knots = { 0.0, 0.0, 1.0, 1.0 };
+        std::vector<double> weights = { 1.0, 1.0 };
+        return PCurve(cps, knots, weights, 1);  // degree-1 = line segment
+    };
+
+    outer.p_curves.push_back(make_line(u0, v0, u1, v0));  // bottom
+    outer.p_curves.push_back(make_line(u1, v0, u1, v1));  // right
+    outer.p_curves.push_back(make_line(u1, v1, u0, v1));  // top
+    outer.p_curves.push_back(make_line(u0, v1, u0, v0));  // left
+
+    ts.loops.push_back(std::move(outer));
+    return ts;
+}
+
+} // namespace vde::brep

tests/brep/CMakeLists.txt

@@ -10,3 +10,4 @@ add_vde_test(test_assembly)
 add_vde_test(test_brep_face_split)
 add_vde_test(test_measure)
 add_vde_test(test_assembly_constraints)
+add_vde_test(test_trimmed_surface)

tests/brep/test_trimmed_surface.cpp

@@ -0,0 +1,247 @@
+#include <gtest/gtest.h>
+#include "vde/brep/trimmed_surface.h"
+#include "vde/curves/nurbs_surface.h"
+#include <cmath>
+
+using namespace vde::brep;
+using namespace vde::curves;
+using vde::core::Point3D;
+
+namespace {
+
+/// Helper: create a simple 10×10 planar NURBS surface in the XY plane
+/// Control grid: 2×2, degree 1×1, knots {0,0,1,1}
+/// Domain: u∈[0,1], v∈[0,1] → maps to XY square [0,10]×[0,10] at Z=0
+NurbsSurface make_plane() {
+    std::vector<std::vector<Point3D>> cp = {
+        {Point3D(0, 0, 0), Point3D(0, 10, 0)},
+        {Point3D(10, 0, 0), Point3D(10, 10, 0)}
+    };
+    std::vector<double> ku = {0, 0, 1, 1};
+    std::vector<double> kv = {0, 0, 1, 1};
+    std::vector<std::vector<double>> w = {{1, 1}, {1, 1}};
+    return NurbsSurface(cp, ku, kv, w, 1, 1);
+}
+
+/// Helper: create a linear p-curve line segment in parameter space
+PCurve make_pcurve_line(double ua, double va, double ub, double vb) {
+    std::vector<Point3D> cps = { Point3D(ua, va, 0.0), Point3D(ub, vb, 0.0) };
+    std::vector<double> knots = { 0.0, 0.0, 1.0, 1.0 };
+    std::vector<double> weights = { 1.0, 1.0 };
+    return PCurve(cps, knots, weights, 1);
+}
+
+/// Helper: create a rectangular trim loop
+TrimLoop make_rect_loop(double u0, double u1, double v0, double v1) {
+    TrimLoop loop;
+    loop.is_outer = true;
+    loop.p_curves.push_back(make_pcurve_line(u0, v0, u1, v0));  // bottom
+    loop.p_curves.push_back(make_pcurve_line(u1, v0, u1, v1));  // right
+    loop.p_curves.push_back(make_pcurve_line(u1, v1, u0, v1));  // top
+    loop.p_curves.push_back(make_pcurve_line(u0, v1, u0, v0));  // left
+    return loop;
+}
+
+} // anonymous namespace
+
+// ═══════════════════════════════════════════════════════════
+//  Untrimmed surface
+// ═══════════════════════════════════════════════════════════
+
+TEST(TrimmedSurfaceTest, Untrimmed_EvaluatesAllCorners) {
+    auto surf = make_plane();
+    auto ts = TrimmedSurface::from_surface(surf);
+
+    // All points in the parameter domain should evaluate successfully
+    auto p00 = ts.evaluate(0.0, 0.0);
+    ASSERT_TRUE(p00.has_value());
+    EXPECT_NEAR(p00->x(), 0.0, 1e-6);
+    EXPECT_NEAR(p00->y(), 0.0, 1e-6);
+    EXPECT_NEAR(p00->z(), 0.0, 1e-6);
+
+    auto p11 = ts.evaluate(1.0, 1.0);
+    ASSERT_TRUE(p11.has_value());
+    EXPECT_NEAR(p11->x(), 10.0, 1e-6);
+    EXPECT_NEAR(p11->y(), 10.0, 1e-6);
+    EXPECT_NEAR(p11->z(), 0.0, 1e-6);
+
+    auto p10 = ts.evaluate(1.0, 0.0);
+    ASSERT_TRUE(p10.has_value());
+    EXPECT_NEAR(p10->x(), 10.0, 1e-6);
+    EXPECT_NEAR(p10->y(), 0.0, 1e-6);
+
+    auto p01 = ts.evaluate(0.0, 1.0);
+    ASSERT_TRUE(p01.has_value());
+    EXPECT_NEAR(p01->x(), 0.0, 1e-6);
+    EXPECT_NEAR(p01->y(), 10.0, 1e-6);
+
+    // Center
+    auto pc = ts.evaluate(0.5, 0.5);
+    ASSERT_TRUE(pc.has_value());
+    EXPECT_NEAR(pc->x(), 5.0, 1e-6);
+    EXPECT_NEAR(pc->y(), 5.0, 1e-6);
+}
+
+TEST(TrimmedSurfaceTest, Untrimmed_IsInsideAlwaysTrue) {
+    auto surf = make_plane();
+    auto ts = TrimmedSurface::from_surface(surf);
+
+    EXPECT_TRUE(ts.is_inside(0.0, 0.0));
+    EXPECT_TRUE(ts.is_inside(0.5, 0.5));
+    EXPECT_TRUE(ts.is_inside(1.0, 1.0));
+    EXPECT_TRUE(ts.is_inside(0.123, 0.789));
+}
+
+// ═══════════════════════════════════════════════════════════
+//  Rectangular trim
+// ═══════════════════════════════════════════════════════════
+
+TEST(TrimmedSurfaceTest, RectTrim_InsideRegion) {
+    auto surf = make_plane();
+    auto ts = TrimmedSurface::from_rect(surf, 0.0, 0.5, 0.0, 0.5);
+
+    // Point at (0.25, 0.25) → 3D: (2.5, 2.5, 0) — inside trim
+    auto p = ts.evaluate(0.25, 0.25);
+    ASSERT_TRUE(p.has_value());
+    EXPECT_NEAR(p->x(), 2.5, 1e-6);
+    EXPECT_NEAR(p->y(), 2.5, 1e-6);
+    EXPECT_NEAR(p->z(), 0.0, 1e-6);
+
+    EXPECT_TRUE(ts.is_inside(0.25, 0.25));
+}
+
+TEST(TrimmedSurfaceTest, RectTrim_OutsideRegion) {
+    auto surf = make_plane();
+    auto ts = TrimmedSurface::from_rect(surf, 0.0, 0.5, 0.0, 0.5);
+
+    // Point at (0.75, 0.25) — outside trim (u > 0.5)
+    auto p = ts.evaluate(0.75, 0.25);
+    EXPECT_FALSE(p.has_value());
+    EXPECT_FALSE(ts.is_inside(0.75, 0.25));
+
+    // Point at (0.25, 0.75) — outside trim (v > 0.5)
+    EXPECT_FALSE(ts.is_inside(0.25, 0.75));
+
+    // Point at (0.75, 0.75) — outside both
+    EXPECT_FALSE(ts.is_inside(0.75, 0.75));
+}
+
+TEST(TrimmedSurfaceTest, RectTrim_BoundaryIsInside) {
+    auto surf = make_plane();
+    auto ts = TrimmedSurface::from_rect(surf, 0.0, 0.5, 0.0, 0.5);
+
+    // Points on the boundary should be considered inside (convention)
+    EXPECT_TRUE(ts.is_inside(0.0, 0.0));
+    EXPECT_TRUE(ts.is_inside(0.5, 0.5));
+    EXPECT_TRUE(ts.is_inside(0.0, 0.5));
+    EXPECT_TRUE(ts.is_inside(0.5, 0.0));
+}
+
+// ═══════════════════════════════════════════════════════════
+//  Hole (inner loop)
+// ═══════════════════════════════════════════════════════════
+
+TEST(TrimmedSurfaceTest, Hole_CenterIsOutside) {
+    auto surf = make_plane();
+    TrimmedSurface ts;
+    ts.base_surface = surf;
+
+    // Outer loop: full [0,1]×[0,1] rectangle
+    ts.loops.push_back(make_rect_loop(0.0, 1.0, 0.0, 1.0));
+
+    // Inner loop (hole): [0.25,0.75]×[0.25,0.75] square
+    auto hole = make_rect_loop(0.25, 0.75, 0.25, 0.75);
+    hole.is_outer = false;
+    ts.loops.push_back(hole);
+
+    // Center (0.5, 0.5) is inside the hole → outside trimmed surface
+    EXPECT_FALSE(ts.is_inside(0.5, 0.5));
+    auto p_center = ts.evaluate(0.5, 0.5);
+    EXPECT_FALSE(p_center.has_value());
+
+    // Corner (0.1, 0.1) is inside outer, outside hole → inside trimmed surface
+    EXPECT_TRUE(ts.is_inside(0.1, 0.1));
+    auto p_corner = ts.evaluate(0.1, 0.1);
+    ASSERT_TRUE(p_corner.has_value());
+    EXPECT_NEAR(p_corner->x(), 1.0, 1e-6);
+    EXPECT_NEAR(p_corner->y(), 1.0, 1e-6);
+
+    // Point inside hole boundary
+    EXPECT_FALSE(ts.is_inside(0.3, 0.3));
+}
+
+// ═══════════════════════════════════════════════════════════
+//  Bounds
+// ═══════════════════════════════════════════════════════════
+
+TEST(TrimmedSurfaceTest, Bounds_Untrimmed) {
+    auto surf = make_plane();
+    auto ts = TrimmedSurface::from_surface(surf);
+    auto b = ts.bounds();
+
+    // The plane spans [0,10]×[0,10] at Z=0
+    EXPECT_NEAR(b.min().x(), 0.0, 1e-6);
+    EXPECT_NEAR(b.max().x(), 10.0, 1e-6);
+    EXPECT_NEAR(b.min().y(), 0.0, 1e-6);
+    EXPECT_NEAR(b.max().y(), 10.0, 1e-6);
+    EXPECT_NEAR(b.min().z(), 0.0, 1e-6);
+    EXPECT_NEAR(b.max().z(), 0.0, 1e-6);
+}
+
+TEST(TrimmedSurfaceTest, Bounds_RectangularTrim) {
+    auto surf = make_plane();
+    auto ts = TrimmedSurface::from_rect(surf, 0.25, 0.75, 0.25, 0.75);
+    auto b = ts.bounds();
+
+    // Trimmed to [2.5, 7.5]×[2.5, 7.5] at Z=0
+    EXPECT_NEAR(b.min().x(), 2.5, 1e-6);
+    EXPECT_NEAR(b.max().x(), 7.5, 1e-6);
+    EXPECT_NEAR(b.min().y(), 2.5, 1e-6);
+    EXPECT_NEAR(b.max().y(), 7.5, 1e-6);
+    EXPECT_NEAR(b.min().z(), 0.0, 1e-6);
+    EXPECT_NEAR(b.max().z(), 0.0, 1e-6);
+}
+
+TEST(TrimmedSurfaceTest, Bounds_WithHole) {
+    auto surf = make_plane();
+    TrimmedSurface ts;
+    ts.base_surface = surf;
+    ts.loops.push_back(make_rect_loop(0.2, 0.8, 0.2, 0.8));
+    auto hole = make_rect_loop(0.3, 0.7, 0.3, 0.7);
+    hole.is_outer = false;
+    ts.loops.push_back(hole);
+
+    auto b = ts.bounds();
+
+    // Bounds are computed from outer loop only, spans [2,8]×[2,8] at Z=0
+    EXPECT_NEAR(b.min().x(), 2.0, 1e-6);
+    EXPECT_NEAR(b.max().x(), 8.0, 1e-6);
+    EXPECT_NEAR(b.min().y(), 2.0, 1e-6);
+    EXPECT_NEAR(b.max().y(), 8.0, 1e-6);
+}
+
+// ═══════════════════════════════════════════════════════════
+//  from_surface constructor
+// ═══════════════════════════════════════════════════════════
+
+TEST(TrimmedSurfaceTest, FromSurface_HasNoLoops) {
+    auto surf = make_plane();
+    auto ts = TrimmedSurface::from_surface(surf);
+
+    EXPECT_TRUE(ts.loops.empty());
+    EXPECT_TRUE(ts.is_inside(0.5, 0.5));
+}
+
+// ═══════════════════════════════════════════════════════════
+//  from_rect constructor
+// ═══════════════════════════════════════════════════════════
+
+TEST(TrimmedSurfaceTest, FromRect_CreatesOneOuterLoop) {
+    auto surf = make_plane();
+    auto ts = TrimmedSurface::from_rect(surf, 0.0, 1.0, 0.0, 1.0);
+
+    EXPECT_EQ(ts.loops.size(), 1u);
+    EXPECT_TRUE(ts.loops[0].is_outer);
+    // Rectangular loop has 4 p-curves (one per side)
+    EXPECT_EQ(ts.loops[0].p_curves.size(), 4u);
+}