use crate::classification::{PointClass, SolidClassifier};
use crate::spatial::Aabb;
use crate::topology::BrepSolid;
use crate::{BooleanOperation, KernelTolerances, Vec3};
use crate::{KernelRefusal, KernelStage, OrRefuse};
use serde::Serialize;
pub const DISAGREEMENT_THRESHOLD: f64 = 0.03;
const SKIP_FRACTION: f64 = 1e-4;
const JITTER_FRACTION: f64 = 4e-3;
const BULK_NUMERATOR: usize = 45;
const BULK_DENOMINATOR: usize = 100;
struct Rng {
state: u64,
}
impl Rng {
fn new(seed: u64) -> Self {
Self { state: seed }
}
fn next_u64(&mut self) -> u64 {
self.state = self.state.wrapping_add(0x9E37_79B9_7F4A_7C15);
let mut z = self.state;
z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9);
z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB);
z ^ (z >> 31)
}
fn unit(&mut self) -> f64 {
(self.next_u64() >> 11) as f64 / (1u64 << 53) as f64
}
fn range(&mut self, low: f64, high: f64) -> f64 {
low + (high - low) * self.unit()
}
fn unit_vector(&mut self) -> Vec3 {
let z = self.range(-1.0, 1.0);
let angle = self.range(0.0, std::f64::consts::TAU);
let radius = (1.0 - z * z).max(0.0).sqrt();
Vec3::new(radius * angle.cos(), radius * angle.sin(), z)
}
}
#[derive(Clone, Copy, Debug, Serialize)]
pub struct SemanticDisagreement {
pub point: Vec3,
pub expected_in: bool,
pub result_in: bool,
pub in_a: bool,
pub in_b: bool,
}
#[derive(Clone, Debug, Serialize)]
pub struct OracleReport {
pub sampled: usize,
pub on_skipped: usize,
pub considered: usize,
pub disagreements: Vec<SemanticDisagreement>,
pub disagreement_rate: f64,
}
impl OracleReport {
pub fn is_flagged(&self) -> bool {
self.disagreement_rate > DISAGREEMENT_THRESHOLD
}
pub fn sample_disagreement(&self) -> Option<SemanticDisagreement> {
self.disagreements.first().copied()
}
}
fn faces_of(solid: &BrepSolid) -> impl Iterator<Item = &crate::topology::FaceRecord> {
solid.shells.iter().flat_map(|shell| shell.faces.iter())
}
fn combined_bounds(first: &BrepSolid, second: &BrepSolid) -> Result<Aabb, KernelRefusal> {
let mut bounds = Aabb::empty();
for face in faces_of(first).chain(faces_of(second)) {
bounds.include(
Aabb::from_surface_controls(&face.surface)
.or_refuse(KernelStage::Validate, "from_surface_controls")?,
);
}
Ok(bounds)
}
fn expected_membership(operation: BooleanOperation, in_a: bool, in_b: bool) -> bool {
match operation {
BooleanOperation::Union => in_a || in_b,
BooleanOperation::Intersect => in_a && in_b,
BooleanOperation::Subtract => in_a && !in_b,
}
}
pub fn boolean_semantic_disagreement(
first: &BrepSolid,
second: &BrepSolid,
operation: BooleanOperation,
result: &BrepSolid,
samples: usize,
) -> Result<OracleReport, KernelRefusal> {
let policy = KernelTolerances::for_pair(first, second, 1e-7);
let model = policy.model;
let bounds = combined_bounds(first, second)?;
let diagonal = bounds.diagonal();
if !diagonal.is_finite() || diagonal <= 0.0 || samples == 0 {
return Ok(OracleReport {
sampled: 0,
on_skipped: 0,
considered: 0,
disagreements: Vec::new(),
disagreement_rate: 0.0,
});
}
let classifier_a =
SolidClassifier::new(first, model).or_refuse(KernelStage::Validate, "new")?;
let classifier_b =
SolidClassifier::new(second, model).or_refuse(KernelStage::Validate, "new")?;
let classifier_r =
SolidClassifier::new(result, model).or_refuse(KernelStage::Validate, "new")?;
let skip_band = (model * 10.0).max(diagonal * SKIP_FRACTION);
let jitter = diagonal * JITTER_FRACTION;
let boundary_faces: Vec<&crate::topology::FaceRecord> = faces_of(first)
.chain(faces_of(second))
.chain(faces_of(result))
.collect();
let bulk_count = samples * BULK_NUMERATOR / BULK_DENOMINATOR;
let mut rng = Rng::new(0x0B00_00AC_1E00_5EED);
let mut report = OracleReport {
sampled: 0,
on_skipped: 0,
considered: 0,
disagreements: Vec::new(),
disagreement_rate: 0.0,
};
for index in 0..samples {
let point = if index < bulk_count || boundary_faces.is_empty() {
Vec3::new(
rng.range(bounds.minimum.x, bounds.maximum.x),
rng.range(bounds.minimum.y, bounds.maximum.y),
rng.range(bounds.minimum.z, bounds.maximum.z),
)
} else {
boundary_sample(&mut rng, &boundary_faces, jitter)
};
report.sampled += 1;
let near = classifier_a
.within_band(point, skip_band)
.or_refuse(KernelStage::Validate, "within_band")?
|| classifier_b
.within_band(point, skip_band)
.or_refuse(KernelStage::Validate, "within_band")?
|| classifier_r
.within_band(point, skip_band)
.or_refuse(KernelStage::Validate, "within_band")?;
if near {
report.on_skipped += 1;
continue;
}
let (Ok(ca), Ok(cb), Ok(cr)) = (
classifier_a.classify(point),
classifier_b.classify(point),
classifier_r.classify(point),
) else {
report.on_skipped += 1;
continue;
};
if ca.class == PointClass::On || cb.class == PointClass::On || cr.class == PointClass::On {
report.on_skipped += 1;
continue;
}
let in_a = ca.class == PointClass::In;
let in_b = cb.class == PointClass::In;
let result_in = cr.class == PointClass::In;
let expected_in = expected_membership(operation, in_a, in_b);
report.considered += 1;
if expected_in != result_in {
report.disagreements.push(SemanticDisagreement {
point,
expected_in,
result_in,
in_a,
in_b,
});
}
}
report.disagreement_rate = if report.considered == 0 {
0.0
} else {
report.disagreements.len() as f64 / report.considered as f64
};
Ok(report)
}
fn expected_membership_nary(operation: BooleanOperation, in_operands: &[bool]) -> bool {
match operation {
BooleanOperation::Union => in_operands.iter().any(|&inside| inside),
BooleanOperation::Intersect => in_operands.iter().all(|&inside| inside),
BooleanOperation::Subtract => {
in_operands[0] && in_operands[1..].iter().all(|&inside| !inside)
}
}
}
pub fn boolean_semantic_disagreement_nary(
operands: &[BrepSolid],
operation: BooleanOperation,
result: &BrepSolid,
samples: usize,
) -> Result<OracleReport, KernelRefusal> {
if operands.is_empty() {
return Err(KernelRefusal::internal(
KernelStage::Validate,
"oracle",
"boolean_semantic_disagreement_nary: no operands",
));
}
let model = KernelTolerances::for_solid(&operands[0], 1e-7).model;
let mut bounds = Aabb::empty();
for operand in operands {
for face in faces_of(operand) {
bounds.include(
Aabb::from_surface_controls(&face.surface)
.or_refuse(KernelStage::Validate, "from_surface_controls")?,
);
}
}
let diagonal = bounds.diagonal();
if !diagonal.is_finite() || diagonal <= 0.0 || samples == 0 {
return Ok(OracleReport {
sampled: 0,
on_skipped: 0,
considered: 0,
disagreements: Vec::new(),
disagreement_rate: 0.0,
});
}
let classifiers = operands
.iter()
.map(|operand| SolidClassifier::new(operand, model))
.collect::<Result<Vec<_>, _>>()
.or_refuse(KernelStage::Validate, "csg.oracle")?;
let classifier_r =
SolidClassifier::new(result, model).or_refuse(KernelStage::Validate, "new")?;
let skip_band = (model * 10.0).max(diagonal * SKIP_FRACTION);
let jitter = diagonal * JITTER_FRACTION;
let boundary_faces: Vec<&crate::topology::FaceRecord> = operands
.iter()
.flat_map(faces_of)
.chain(faces_of(result))
.collect();
let bulk_count = samples * BULK_NUMERATOR / BULK_DENOMINATOR;
let mut rng = Rng::new(0x0B00_00AC_1E00_5EED);
let mut report = OracleReport {
sampled: 0,
on_skipped: 0,
considered: 0,
disagreements: Vec::new(),
disagreement_rate: 0.0,
};
for index in 0..samples {
let point = if index < bulk_count || boundary_faces.is_empty() {
Vec3::new(
rng.range(bounds.minimum.x, bounds.maximum.x),
rng.range(bounds.minimum.y, bounds.maximum.y),
rng.range(bounds.minimum.z, bounds.maximum.z),
)
} else {
boundary_sample(&mut rng, &boundary_faces, jitter)
};
report.sampled += 1;
let mut near = classifier_r
.within_band(point, skip_band)
.or_refuse(KernelStage::Validate, "within_band")?;
for classifier in &classifiers {
near = near
|| classifier
.within_band(point, skip_band)
.or_refuse(KernelStage::Validate, "within_band")?;
}
if near {
report.on_skipped += 1;
continue;
}
let Ok(cr) = classifier_r.classify(point) else {
report.on_skipped += 1;
continue;
};
if cr.class == PointClass::On {
report.on_skipped += 1;
continue;
}
let mut in_operands = Vec::with_capacity(classifiers.len());
let mut ambiguous = false;
for classifier in &classifiers {
match classifier.classify(point) {
Ok(classification) if classification.class != PointClass::On => {
in_operands.push(classification.class == PointClass::In);
}
_ => {
ambiguous = true;
break;
}
}
}
if ambiguous {
report.on_skipped += 1;
continue;
}
let result_in = cr.class == PointClass::In;
let expected_in = expected_membership_nary(operation, &in_operands);
report.considered += 1;
if expected_in != result_in {
report.disagreements.push(SemanticDisagreement {
point,
expected_in,
result_in,
in_a: in_operands[0],
in_b: *in_operands.get(1).unwrap_or(&false),
});
}
}
report.disagreement_rate = if report.considered == 0 {
0.0
} else {
report.disagreements.len() as f64 / report.considered as f64
};
Ok(report)
}
fn boundary_sample(rng: &mut Rng, faces: &[&crate::topology::FaceRecord], jitter: f64) -> Vec3 {
let face = faces[(rng.next_u64() as usize) % faces.len()];
let (Ok([u0, u1]), Ok([v0, v1])) = (face.surface.domain_u(), face.surface.domain_v()) else {
return Vec3::default();
};
let u = rng.range(u0, u1);
let v = rng.range(v0, v1);
let base = match face.surface.evaluate(u, v) {
Ok(point) => point,
Err(_) => return Vec3::default(),
};
let direction = match face.surface.normal(u, v) {
Ok(normal) if normal.length() > 1e-9 => normal,
_ => rng.unit_vector(),
};
let sign = if rng.unit() < 0.5 { -1.0 } else { 1.0 };
base.add(direction.scale(sign * jitter))
}
#[derive(Clone, Copy, Debug, PartialEq, Eq, Serialize)]
#[serde(rename_all = "snake_case")]
pub enum ResidualKind {
CoincidentVertices,
DuplicateEdge,
}
#[derive(Clone, Debug, Serialize)]
pub struct ResidualFusable {
pub kind: ResidualKind,
pub detail: String,
pub point: Vec3,
}
pub fn boolean_residual_fusables(result: &BrepSolid, tol: f64) -> Vec<ResidualFusable> {
let mut findings = Vec::new();
let vertices = &result.vertices;
for i in 0..vertices.len() {
for j in (i + 1)..vertices.len() {
let gap = vertices[i].point.sub(vertices[j].point).length();
if gap <= tol {
findings.push(ResidualFusable {
kind: ResidualKind::CoincidentVertices,
detail: format!(
"vertices {} and {} coincide within {gap:.3e}",
vertices[i].id, vertices[j].id
),
point: vertices[i].point,
});
}
}
}
let edges = &result.edges;
let endpoint =
|edge: &crate::topology::EdgeRecord| -> Result<(Vec3, Vec3, Vec3), KernelRefusal> {
let start = edge
.curve
.evaluate(edge.t0)
.or_refuse(KernelStage::Validate, "evaluate")?;
let end = edge
.curve
.evaluate(edge.t1)
.or_refuse(KernelStage::Validate, "evaluate")?;
let mid = edge
.curve
.evaluate((edge.t0 + edge.t1) * 0.5)
.or_refuse(KernelStage::Validate, "evaluate")?;
Ok((start, mid, end))
};
for i in 0..edges.len() {
let Ok((si, mi, ei)) = endpoint(&edges[i]) else {
continue;
};
for j in (i + 1)..edges.len() {
let Ok((sj, mj, ej)) = endpoint(&edges[j]) else {
continue;
};
let endpoints_match = (si.sub(sj).length() <= tol && ei.sub(ej).length() <= tol)
|| (si.sub(ej).length() <= tol && ei.sub(sj).length() <= tol);
if endpoints_match && mi.sub(mj).length() <= tol {
findings.push(ResidualFusable {
kind: ResidualKind::DuplicateEdge,
detail: format!(
"edges {} and {} trace the same curve within {tol:.3e}",
edges[i].id, edges[j].id
),
point: mi,
});
}
}
}
findings
}