port: surface tessellation splits, iso boundaries and CSI leaf driver

Co-Authored-By: OpenAI Codex <noreply@openai.com>
This commit is contained in:
lillith
2026-10-13 12:42:52 +02:00
parent c16729604a
commit 2b6bcd6ac3
6 changed files with 3590 additions and 0 deletions
@@ -0,0 +1,318 @@
//! Differential runner for surface.surface_p9. Reference: SmSurface.cpp:11644,
//! 13887, 14877. A quadratic Bezier patch is used on both sides. CSI fixtures use
//! one planar patch and one line span, so both subdivision trees have one leaf.
use crate::{Case, CaseContext, HarnessError, Output, Rng, UnitDiff};
use serde::{Deserialize, Serialize};
use serde_json::Value;
use solidcore::curve::curve::Curve;
use solidcore::curve::curve_cache::{CurveCache, CurveCacheOptions};
use solidcore::geom::extent1d::Extent1d;
use solidcore::geom::extent2d::Extent2d;
use solidcore::geom::extent3d::Extent3d;
use solidcore::geom::global_solver::{GlobalSolver, SolverOperation};
use solidcore::geom::pseudo_box::PseudoBox;
use solidcore::geom::tree::{
AuxDataType, NodeClass, NodeData, NodeId, SplitDir, SurfaceAux, Tree, TreeNode, TreeOwner,
};
use solidcore::geom::vector3d::Vec3;
use solidcore::kernel::Kernel;
use solidcore::nurbs::b_spline_curve::BSplineCurve;
use solidcore::nurbs::b_spline_surface::{BSplineSurface, SurfParamType};
use solidcore::nurbs::crv_basic::Curve as NlCurve;
use solidcore::nurbs::math::{CPoint, NOW};
use solidcore::nurbs::srf_basic::Surface as NlSurface;
use solidcore::surface::surface::p9::CsiGlobalSolver;
use solidcore::surface::surface::Surface;
use solidcore::surface::surface_cache::{
BezierAux2d, BezierPatch, SurfaceCache, SurfaceCacheOptions,
};
use solidcore_ffi::units::surface_surface_p9 as ffi;
pub struct Runner;
#[derive(Debug, Clone, Serialize, Deserialize)]
struct Inputs {
op: i32,
net: [[f64; 3]; 9],
#[serde(default = "domain")]
domain: [f64; 4],
direction: i32,
controls: [f64; 3],
#[serde(default)]
diagonals: bool,
#[serde(default = "line")]
line: [[f64; 3]; 2],
}
fn domain() -> [f64; 4] {
[0.0, 0.0, 1.0, 1.0]
}
fn line() -> [[f64; 3]; 2] {
[[0.5, 0.5, -5.0], [0.5, 0.5, 5.0]]
}
fn harness(e: impl std::fmt::Display) -> HarnessError {
HarnessError::new(e.to_string())
}
fn surface(i: &Inputs) -> Result<Surface, HarnessError> {
let net = i
.net
.chunks(3)
.map(|row| {
row.iter()
.map(|p| CPoint::new(p[0], p[1], p[2], NOW))
.collect()
})
.collect();
let knots = vec![0.0, 0.0, 0.0, 1.0, 1.0, 1.0];
let nurb = NlSurface::new(net, 2, 2, knots.clone(), knots).map_err(harness)?;
Ok(Surface::from(BSplineSurface::from_nurb(nurb, false)))
}
fn output(flat: Vec<f64>) -> Output {
let (status, rest) = flat.split_first().unwrap_or((&-1.0, &[]));
Output::new()
.scalar("status", *status)
.values("out", rest.to_vec())
}
fn run_rust(i: &Inputs) -> Result<Output, HarnessError> {
let k = Kernel::default();
let surface = surface(i)?;
let direction =
SurfParamType::from_code(i.direction as u64).ok_or_else(|| harness("direction"))?;
let domain = Extent2d::from_bounds(i.domain[0], i.domain[1], i.domain[2], i.domain[3]);
let result = match i.op {
0 => surface
.find_tessellation_splits(
&k,
direction,
&domain,
i.controls[0],
i.controls[1],
i.controls[2],
i.diagonals,
)
.map(|v| {
let mut out = vec![0.0];
out.extend(v);
out
}),
1 => surface
.create_iso_boundaries(&k, direction, i.controls[0])
.and_then(|curves| {
let mut out = vec![0.0, if curves.is_some() { 2.0 } else { 0.0 }];
if let Some((a, b)) = curves {
for curve in [a, b] {
for n in 0..=4 {
let p = curve.evaluate_point(
curve.natural_interval().evaluate(f64::from(n) / 4.0),
)?;
out.extend([p.x, p.y, p.z]);
}
}
}
Ok(out)
}),
2 => {
let points = i
.line
.iter()
.map(|p| CPoint::new(p[0], p[1], p[2], NOW))
.collect();
let curve = Curve::from(BSplineCurve::from_nurb(
NlCurve::new(points, 1, vec![0.0, 0.0, 1.0, 1.0]).map_err(harness)?,
3,
false,
));
let interval = Extent1d::from_min_max(0.0, 1.0);
let mut cc = CurveCache::<i64>::new(
curve.clone(),
interval,
CurveCacheOptions {
make_full_tree: true,
..CurveCacheOptions::default()
},
);
cc.set_curve_handle(Some(0));
match cc.tessellate(&k) {
Err(e) => Err(e),
Ok(()) => {
let mut cache = SurfaceCache::<i64>::new(
surface,
Extent2d::from_bounds(0.0, 0.0, 1.0, 1.0),
SurfaceCacheOptions::default(),
);
cache.set_surface_handle(Some(1));
let uv = cache.uv_domain();
let mut bbox = Extent3d::new();
for p in i.net {
bbox = bbox.union(&Extent3d::from_points(
Vec3::new(p[0], p[1], p[2]),
Vec3::new(p[0], p[1], p[2]),
));
}
cache.patches_mut().push(BezierPatch {
aux: BezierAux2d {
uv_domain: uv,
owning_tree_node: Some(NodeId(0)),
..BezierAux2d::default()
},
pseudo_box: PseudoBox::from_extent3d(&bbox),
..BezierPatch::default()
});
let mut tree = Tree::new(Some(TreeOwner::SurfaceCache));
tree.add_node(TreeNode {
bbox,
aux_type: AuxDataType::BezierSurface,
data: NodeData::Surface(SurfaceAux {
split_dir: SplitDir::Unknown,
is_marked: false,
node_class: NodeClass::Unknown,
face: None,
surface: Some(1),
slot: Some(0),
}),
..TreeNode::new()
});
tree.set_top_node(Some(NodeId(0)));
cache.set_tree(Some(tree));
let curve_tree = cc.tree().ok_or_else(|| harness("curve tree"))?;
let tree = cache.tree().ok_or_else(|| harness("surface tree"))?;
let mut base = GlobalSolver::with_trees(curve_tree, tree, 3);
base.operation = SolverOperation::Intersect;
base.tolerance_3d = 1e-6;
let mut solver = CsiGlobalSolver {
base,
curve: &curve,
surface_cache: &cache,
kernel: &k,
found_coincidence: false,
interval,
uv_domain: domain,
};
solver
.local_solve_branch(&[Some(NodeId(0)), Some(NodeId(0)), None, None])
.map(|more| {
let mut out = vec![
0.0,
f64::from(u8::from(more)),
solver.base.solutions.len() as f64,
];
for solution in solver.base.solutions.as_slice() {
out.extend([
f64::from(solution.solution_type.code()),
solution.num_variables as f64,
solution.num_objects as f64,
solution.start.solution_value,
solution.end.solution_value,
]);
out.extend(solution.start.parameters.iter().take(3));
out.extend(solution.end.parameters.iter().take(3));
out.extend(
solution
.nodes
.iter()
.take(2)
.map(|v| v.map_or(-1.0, |n| n.0 as f64)),
);
out.extend(
solution
.objects
.iter()
.take(2)
.map(|v| v.map_or(-1.0, |n| n as f64)),
);
}
out
})
}
}
}
_ => return Err(harness("invalid op")),
};
Ok(match result {
Ok(v) => output(v),
Err(e) => Output::unported(&e).unwrap_or_else(|| {
let status = match e {
solidcore::SolidError::Upstream { code, .. } => f64::from(code),
_ => 1007.0,
};
output(vec![status])
}),
})
}
impl UnitDiff for Runner {
fn id(&self) -> &'static str {
"surface.surface_p9"
}
fn cases(&self, ctx: &CaseContext) -> Result<Vec<Case>, HarnessError> {
let ctx = CaseContext {
seeds: ctx.seeds.max(600),
..ctx.clone()
};
crate::runner::default_cases(self, &ctx)
}
fn run_cpp(&self, case: &Case) -> Result<Output, HarnessError> {
let i: Inputs = case.parse()?;
ffi::eval(
i.op,
&i.net,
i.domain,
i.direction,
i.controls,
i.diagonals,
i.line,
)
.map(output)
.map_err(harness)
}
fn run_rust(&self, case: &Case) -> Result<Output, HarnessError> {
run_rust(&case.parse()?)
}
fn generate(&self, rng: &mut Rng) -> Option<Value> {
let mut net = [[0.0; 3]; 9];
let scale = rng.uniform(0.5, 3.0);
let height = rng.uniform(0.0, 1.5);
for (index, p) in net.iter_mut().enumerate() {
*p = [
(index / 3) as f64 * scale,
(index % 3) as f64 * scale,
rng.uniform(-height, height),
];
}
let direction = rng.below(2) as i32;
let (op, controls, diagonals) = if rng.chance(0.25) {
(1, [0.01, 0.0, 0.0], false)
} else {
let mode = rng.below(4);
(
0,
[
if mode == 0 {
0.0
} else {
rng.uniform(0.02, 0.2)
},
if mode < 2 {
0.0
} else {
rng.uniform(5.0, 30.0)
},
if mode == 3 {
rng.uniform(0.5, 2.0)
} else {
0.0
},
],
rng.chance(0.5),
)
};
serde_json::to_value(Inputs {
op,
net,
direction,
domain: domain(),
controls,
diagonals,
line: line(),
})
.ok()
}
}
@@ -0,0 +1,110 @@
//! Safe wrappers for `surface.surface_p9`; shim: `surface_surface_p9.cpp`.
//! Upstream: `source/SMLib/src/SmSurface.cpp` (tessellation, iso boundaries and CSI).
//! Inputs are one non-rational quadratic 3x3 Bezier patch, a domain within [0,1]^2,
//! bounded controls and, for CSI, a non-degenerate line. No counts cross the FFI.
use crate::{call_with_buffer, FfiError};
mod sys {
extern "C" {
pub fn sc_surface_surface_p9_eval(
op: i32,
xyz: *const f64,
bounds: *const f64,
direction: i32,
chord: f64,
angle: f64,
max_side: f64,
diagonals: i32,
line: *const f64,
output: *mut f64,
cap: i64,
out_len: *mut i64,
) -> i32;
}
}
/// Evaluates a bounded surface case: op 0 splits, 1 iso boundaries, 2 CSI leaf solve.
///
/// # Errors
/// Invalid data is rejected before C++; shim failures return `FfiError::Status`.
pub fn eval(
op: i32,
net: &[[f64; 3]; 9],
domain: [f64; 4],
direction: i32,
controls: [f64; 3],
diagonals: bool,
line: [[f64; 3]; 2],
) -> Result<Vec<f64>, FfiError> {
let invalid = || FfiError::InvalidArgument("unsafe surface p9 input".into());
if !(0..=2).contains(&op)
|| !(0..=9).contains(&direction)
|| net
.iter()
.flatten()
.any(|v| !v.is_finite() || v.abs() > 10.0)
|| domain
.iter()
.any(|v| !v.is_finite() || !(0.0..=1.0).contains(v))
|| domain[2] - domain[0] < 1e-4
|| domain[3] - domain[1] < 1e-4
|| controls.iter().any(|v| !v.is_finite() || *v < 0.0)
|| (controls[0] != 0.0 && !(1e-3..=10.0).contains(&controls[0]))
|| (controls[1] != 0.0 && !(1.0..=90.0).contains(&controls[1]))
|| (controls[2] != 0.0 && !(0.05..=100.0).contains(&controls[2]))
|| line
.iter()
.flatten()
.any(|v| !v.is_finite() || v.abs() > 20.0)
{
return Err(invalid());
}
if op == 0 && direction > 1 {
return Err(invalid());
}
if controls[2] > 0.0 {
let distances = if direction == 0 {
[(0, 2), (3, 5), (6, 8)]
} else {
[(0, 6), (1, 7), (2, 8)]
};
for (a, b) in distances {
let pa = net.get(a).ok_or_else(invalid)?;
let pb = net.get(b).ok_or_else(invalid)?;
if pa
.iter()
.zip(pb)
.map(|(a, b)| (a - b) * (a - b))
.sum::<f64>()
< 1e-4
{
return Err(invalid());
}
}
}
if op == 2 && line[0] == line[1] {
return Err(invalid());
}
let xyz: Vec<f64> = net.iter().flatten().copied().collect();
let ray: Vec<f64> = line.iter().flatten().copied().collect();
call_with_buffer(128, |buffer, len| {
// SAFETY: fixed input lengths, finite bounded values and valid enumerators
// were checked above; the shim copies output only within buffer capacity.
unsafe {
sys::sc_surface_surface_p9_eval(
op,
xyz.as_ptr(),
domain.as_ptr(),
direction,
controls[0],
controls[1],
controls[2],
i32::from(diagonals),
ray.as_ptr(),
buffer.as_mut_ptr(),
buffer.len() as i64,
&raw mut *len,
)
}
})
}
+1
View File
@@ -8,5 +8,6 @@
//! functions).
pub mod p2;
pub mod p9;
pub use p2::{Surface, SurfaceKind};
File diff suppressed because it is too large Load Diff
File diff suppressed because it is too large Load Diff
+99
View File
@@ -0,0 +1,99 @@
// SPDX-License-Identifier: Apache-2.0
// Unit surface.surface_p9. Calls the pinned SmSurface methods and the file-local
// SmCSIGlobalSolver class (declaration reproduced from SmSurface.cpp:13781).
#include "sc_shim.h"
#include "StdAfx.h"
#include <SmBSplineSurface.h>
#include <SmBSplineCurve.h>
#include <SmContext.h>
#include <SmCurveCache.h>
#include <SmSurfaceCache.h>
#include <SmGlobalSolver.h>
#include <SmSolutionArray.h>
#include <SmTree.h>
#include <nurbs.h>
#include <vector>
class SmCSIGlobalSolver : public SmGlobalSolver {
public:
const SmSurface* m_pSurf;
SmBoolean m_bFoundCoincidence;
const SmExtent1d& m_crInterval;
const SmExtent2d& m_crUVDomain;
SmCSIGlobalSolver(SmTree*, SmTree*, SmBoolean, const SmExtent1d&, const SmExtent2d&);
virtual SmStatus LocalSolve(SmTreeNode* [SM_GS_MAX_TREES], SmBoolean&);
virtual SmSolutionPairType IdenticalSolutions(const SmSolution&, const SmSolution&, const SmSolution*&) const;
};
SC_EXPORT int32_t sc_surface_surface_p9_eval(int32_t op, const double* xyz, const double* bounds,
int32_t direction, double chord, double angle, double max_side, int32_t diagonals,
const double* line, double* output, int64_t cap, int64_t* out_len) {
SC_TRY
SmContext context;
gw_CPOINT points[3][3]; gw_CPOINT* rows[3] = {points[0],points[1],points[2]};
for (int i=0;i<3;i++) for(int j=0;j<3;j++) {
auto& p=points[i][j]; const int at=3*(3*i+j);
p.x=xyz[at];p.y=xyz[at+1];p.z=xyz[at+2];p.w=NL_NOW;
}
double knots[6]={0,0,0,1,1,1};
SmBSplineSurface* surface=new(context) SmBSplineSurface(2,2,2,2,5,5,rows,knots,knots);
SmExtent2d domain(bounds[0],bounds[1],bounds[2],bounds[3]);
std::vector<double> out;
auto status=[](SmStatus r){return r==SM_SUCCESS ? 0.0 : (double)r;};
if(op==0) {
SmTArray<double> splits;
SmStatus r=surface->FindTessellationSplits((SmSurfParamType)direction,domain,chord,angle,max_side,splits,diagonals?TRUE:FALSE);
out.push_back(status(r));
if(r==SM_SUCCESS) for(ULONG i=0;i<splits.GetSize();i++) out.push_back(splits[i]);
} else if(op==1) {
SmBSplineCurve* a=NULL;SmBSplineCurve* b=NULL;
SmStatus r=surface->CreateIsoBoundaries(context,(SmSurfParamType)direction,SmApproxTol3d(chord),a,b);
out.push_back(status(r));
if(r==SM_SUCCESS) {
out.push_back((a?1.0:0.0)+(b?1.0:0.0));
for(auto curve:{a,b}) if(curve) for(int i=0;i<=4;i++) {
SmPoint3d p;
SmStatus er=curve->EvaluatePoint(curve->GetNaturalInterval().Evaluate((double)i/4.0),p);
if(er!=SM_SUCCESS) { delete a;delete b;delete surface;return SC_ERR_ARG; }
out.push_back(p.x);out.push_back(p.y);out.push_back(p.z);
}
}
delete a;delete b;
} else if(op==2) {
gw_CPOINT cp[2]={{line[0],line[1],line[2],NL_NOW},{line[3],line[4],line[5],NL_NOW}};
double ck[4]={0,0,1,1};
SmBSplineCurve* curve=new(context) SmBSplineCurve(3,1,3,1,cp,ck);
{
SmExtent1d interval(0.0,1.0);
SmCurveCache cc(*curve,interval,TRUE);
SmSurfaceCache sc(*surface,surface->GetNaturalUVDomain());
SmStatus r=cc.Tessellate();
if(r==SM_SUCCESS) r=sc.BuildTree(NULL);
if(r==SM_SUCCESS) {
SmCSIGlobalSolver solver(cc.GetTree(),sc.GetTree(),FALSE,interval,domain);
SmSolutionArray solutions;solver.SetSolutions(&solutions);
solver.SetSolverOperation(SM_SO_INTERSECT);solver.Set3DTolerance(1e-6);solver.SetSolutionRequested(SM_SR_ALL);
SmTreeNode* branch[SM_GS_MAX_TREES]={cc.GetTree()->GetFirstLeafNode(),sc.GetTree()->GetFirstLeafNode(),NULL,NULL};
SmBoolean more=FALSE;r=solver.LocalSolve(branch,more);
out.push_back(status(r));
if(r==SM_SUCCESS) {
out.push_back(more==TRUE?1.0:0.0);out.push_back((double)solutions.GetSize());
for(ULONG i=0;i<solutions.GetSize();i++) {
const SmSolution& s=solutions[i];
out.push_back((double)s.m_eSolutionType);out.push_back((double)s.m_lNumVariables);out.push_back((double)s.m_lNumObjects);
out.push_back(s.m_vStart.m_dSolutionValue);out.push_back(s.m_vEnd.m_dSolutionValue);
for(int j=0;j<3;j++) out.push_back(s.m_vStart[j]);
for(int j=0;j<3;j++) out.push_back(s.m_vEnd[j]);
for(int j=0;j<2;j++) out.push_back(s.m_apNodes[j]==branch[j]?0.0:-1.0);
for(int j=0;j<2;j++) out.push_back(s.m_apObjects[j]==curve?0.0:(s.m_apObjects[j]==surface?1.0:-1.0));
}
}
} else out.push_back(status(r));
}
delete curve;
} else {delete surface;return SC_ERR_ARG;}
delete surface;
return sc_copy_out(out.data(),(int64_t)out.size(),output,cap,out_len);
SC_CATCH
}