1use crate::classification::{PointClass, SolidClassifier};
34use crate::spatial::Aabb;
35use crate::topology::BrepSolid;
36use crate::{BooleanOperation, KernelTolerances, Vec3};
37use crate::{KernelRefusal, KernelStage, OrRefuse};
38use serde::Serialize;
39
40pub const DISAGREEMENT_THRESHOLD: f64 = 0.03;
51
52const SKIP_FRACTION: f64 = 1e-4;
56
57const JITTER_FRACTION: f64 = 4e-3;
63
64const BULK_NUMERATOR: usize = 45;
67const BULK_DENOMINATOR: usize = 100;
68
69struct Rng {
73 state: u64,
74}
75
76impl Rng {
77 fn new(seed: u64) -> Self {
78 Self { state: seed }
79 }
80
81 fn next_u64(&mut self) -> u64 {
82 self.state = self.state.wrapping_add(0x9E37_79B9_7F4A_7C15);
83 let mut z = self.state;
84 z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9);
85 z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB);
86 z ^ (z >> 31)
87 }
88
89 fn unit(&mut self) -> f64 {
91 (self.next_u64() >> 11) as f64 / (1u64 << 53) as f64
92 }
93
94 fn range(&mut self, low: f64, high: f64) -> f64 {
96 low + (high - low) * self.unit()
97 }
98
99 fn unit_vector(&mut self) -> Vec3 {
101 let z = self.range(-1.0, 1.0);
102 let angle = self.range(0.0, std::f64::consts::TAU);
103 let radius = (1.0 - z * z).max(0.0).sqrt();
104 Vec3::new(radius * angle.cos(), radius * angle.sin(), z)
105 }
106}
107
108#[derive(Clone, Copy, Debug, Serialize)]
110pub struct SemanticDisagreement {
111 pub point: Vec3,
112 pub expected_in: bool,
114 pub result_in: bool,
116 pub in_a: bool,
117 pub in_b: bool,
118}
119
120#[derive(Clone, Debug, Serialize)]
122pub struct OracleReport {
123 pub sampled: usize,
125 pub on_skipped: usize,
128 pub considered: usize,
130 pub disagreements: Vec<SemanticDisagreement>,
132 pub disagreement_rate: f64,
134}
135
136impl OracleReport {
137 pub fn is_flagged(&self) -> bool {
139 self.disagreement_rate > DISAGREEMENT_THRESHOLD
140 }
141
142 pub fn sample_disagreement(&self) -> Option<SemanticDisagreement> {
144 self.disagreements.first().copied()
145 }
146}
147
148fn faces_of(solid: &BrepSolid) -> impl Iterator<Item = &crate::topology::FaceRecord> {
149 solid.shells.iter().flat_map(|shell| shell.faces.iter())
150}
151
152fn combined_bounds(first: &BrepSolid, second: &BrepSolid) -> Result<Aabb, KernelRefusal> {
153 let mut bounds = Aabb::empty();
154 for face in faces_of(first).chain(faces_of(second)) {
155 bounds.include(
156 Aabb::from_surface_controls(&face.surface)
157 .or_refuse(KernelStage::Validate, "from_surface_controls")?,
158 );
159 }
160 Ok(bounds)
161}
162
163fn expected_membership(operation: BooleanOperation, in_a: bool, in_b: bool) -> bool {
165 match operation {
166 BooleanOperation::Union => in_a || in_b,
167 BooleanOperation::Intersect => in_a && in_b,
168 BooleanOperation::Subtract => in_a && !in_b,
169 }
170}
171
172pub fn boolean_semantic_disagreement(
178 first: &BrepSolid,
179 second: &BrepSolid,
180 operation: BooleanOperation,
181 result: &BrepSolid,
182 samples: usize,
183) -> Result<OracleReport, KernelRefusal> {
184 let policy = KernelTolerances::for_pair(first, second, 1e-7);
185 let model = policy.model;
186
187 let bounds = combined_bounds(first, second)?;
188 let diagonal = bounds.diagonal();
189 if !diagonal.is_finite() || diagonal <= 0.0 || samples == 0 {
190 return Ok(OracleReport {
191 sampled: 0,
192 on_skipped: 0,
193 considered: 0,
194 disagreements: Vec::new(),
195 disagreement_rate: 0.0,
196 });
197 }
198
199 let classifier_a =
200 SolidClassifier::new(first, model).or_refuse(KernelStage::Validate, "new")?;
201 let classifier_b =
202 SolidClassifier::new(second, model).or_refuse(KernelStage::Validate, "new")?;
203 let classifier_r =
204 SolidClassifier::new(result, model).or_refuse(KernelStage::Validate, "new")?;
205
206 let skip_band = (model * 10.0).max(diagonal * SKIP_FRACTION);
207 let jitter = diagonal * JITTER_FRACTION;
208
209 let boundary_faces: Vec<&crate::topology::FaceRecord> = faces_of(first)
213 .chain(faces_of(second))
214 .chain(faces_of(result))
215 .collect();
216
217 let bulk_count = samples * BULK_NUMERATOR / BULK_DENOMINATOR;
218
219 let mut rng = Rng::new(0x0B00_00AC_1E00_5EED);
220 let mut report = OracleReport {
221 sampled: 0,
222 on_skipped: 0,
223 considered: 0,
224 disagreements: Vec::new(),
225 disagreement_rate: 0.0,
226 };
227
228 for index in 0..samples {
229 let point = if index < bulk_count || boundary_faces.is_empty() {
230 Vec3::new(
231 rng.range(bounds.minimum.x, bounds.maximum.x),
232 rng.range(bounds.minimum.y, bounds.maximum.y),
233 rng.range(bounds.minimum.z, bounds.maximum.z),
234 )
235 } else {
236 boundary_sample(&mut rng, &boundary_faces, jitter)
237 };
238 report.sampled += 1;
239
240 let near = classifier_a
244 .within_band(point, skip_band)
245 .or_refuse(KernelStage::Validate, "within_band")?
246 || classifier_b
247 .within_band(point, skip_band)
248 .or_refuse(KernelStage::Validate, "within_band")?
249 || classifier_r
250 .within_band(point, skip_band)
251 .or_refuse(KernelStage::Validate, "within_band")?;
252 if near {
253 report.on_skipped += 1;
254 continue;
255 }
256
257 let (Ok(ca), Ok(cb), Ok(cr)) = (
258 classifier_a.classify(point),
259 classifier_b.classify(point),
260 classifier_r.classify(point),
261 ) else {
262 report.on_skipped += 1;
265 continue;
266 };
267 if ca.class == PointClass::On || cb.class == PointClass::On || cr.class == PointClass::On {
270 report.on_skipped += 1;
271 continue;
272 }
273
274 let in_a = ca.class == PointClass::In;
275 let in_b = cb.class == PointClass::In;
276 let result_in = cr.class == PointClass::In;
277 let expected_in = expected_membership(operation, in_a, in_b);
278 report.considered += 1;
279 if expected_in != result_in {
280 report.disagreements.push(SemanticDisagreement {
281 point,
282 expected_in,
283 result_in,
284 in_a,
285 in_b,
286 });
287 }
288 }
289
290 report.disagreement_rate = if report.considered == 0 {
291 0.0
292 } else {
293 report.disagreements.len() as f64 / report.considered as f64
294 };
295 Ok(report)
296}
297
298fn expected_membership_nary(operation: BooleanOperation, in_operands: &[bool]) -> bool {
304 match operation {
305 BooleanOperation::Union => in_operands.iter().any(|&inside| inside),
306 BooleanOperation::Intersect => in_operands.iter().all(|&inside| inside),
307 BooleanOperation::Subtract => {
308 in_operands[0] && in_operands[1..].iter().all(|&inside| !inside)
309 }
310 }
311}
312
313pub fn boolean_semantic_disagreement_nary(
321 operands: &[BrepSolid],
322 operation: BooleanOperation,
323 result: &BrepSolid,
324 samples: usize,
325) -> Result<OracleReport, KernelRefusal> {
326 if operands.is_empty() {
327 return Err(KernelRefusal::internal(
328 KernelStage::Validate,
329 "oracle",
330 "boolean_semantic_disagreement_nary: no operands",
331 ));
332 }
333 let model = KernelTolerances::for_solid(&operands[0], 1e-7).model;
334
335 let mut bounds = Aabb::empty();
336 for operand in operands {
337 for face in faces_of(operand) {
338 bounds.include(
339 Aabb::from_surface_controls(&face.surface)
340 .or_refuse(KernelStage::Validate, "from_surface_controls")?,
341 );
342 }
343 }
344 let diagonal = bounds.diagonal();
345 if !diagonal.is_finite() || diagonal <= 0.0 || samples == 0 {
346 return Ok(OracleReport {
347 sampled: 0,
348 on_skipped: 0,
349 considered: 0,
350 disagreements: Vec::new(),
351 disagreement_rate: 0.0,
352 });
353 }
354
355 let classifiers = operands
356 .iter()
357 .map(|operand| SolidClassifier::new(operand, model))
358 .collect::<Result<Vec<_>, _>>()
359 .or_refuse(KernelStage::Validate, "csg.oracle")?;
360 let classifier_r =
361 SolidClassifier::new(result, model).or_refuse(KernelStage::Validate, "new")?;
362
363 let skip_band = (model * 10.0).max(diagonal * SKIP_FRACTION);
364 let jitter = diagonal * JITTER_FRACTION;
365
366 let boundary_faces: Vec<&crate::topology::FaceRecord> = operands
367 .iter()
368 .flat_map(faces_of)
369 .chain(faces_of(result))
370 .collect();
371
372 let bulk_count = samples * BULK_NUMERATOR / BULK_DENOMINATOR;
373
374 let mut rng = Rng::new(0x0B00_00AC_1E00_5EED);
375 let mut report = OracleReport {
376 sampled: 0,
377 on_skipped: 0,
378 considered: 0,
379 disagreements: Vec::new(),
380 disagreement_rate: 0.0,
381 };
382
383 for index in 0..samples {
384 let point = if index < bulk_count || boundary_faces.is_empty() {
385 Vec3::new(
386 rng.range(bounds.minimum.x, bounds.maximum.x),
387 rng.range(bounds.minimum.y, bounds.maximum.y),
388 rng.range(bounds.minimum.z, bounds.maximum.z),
389 )
390 } else {
391 boundary_sample(&mut rng, &boundary_faces, jitter)
392 };
393 report.sampled += 1;
394
395 let mut near = classifier_r
397 .within_band(point, skip_band)
398 .or_refuse(KernelStage::Validate, "within_band")?;
399 for classifier in &classifiers {
400 near = near
401 || classifier
402 .within_band(point, skip_band)
403 .or_refuse(KernelStage::Validate, "within_band")?;
404 }
405 if near {
406 report.on_skipped += 1;
407 continue;
408 }
409
410 let Ok(cr) = classifier_r.classify(point) else {
411 report.on_skipped += 1;
412 continue;
413 };
414 if cr.class == PointClass::On {
415 report.on_skipped += 1;
416 continue;
417 }
418 let mut in_operands = Vec::with_capacity(classifiers.len());
419 let mut ambiguous = false;
420 for classifier in &classifiers {
421 match classifier.classify(point) {
422 Ok(classification) if classification.class != PointClass::On => {
423 in_operands.push(classification.class == PointClass::In);
424 }
425 _ => {
426 ambiguous = true;
427 break;
428 }
429 }
430 }
431 if ambiguous {
432 report.on_skipped += 1;
433 continue;
434 }
435
436 let result_in = cr.class == PointClass::In;
437 let expected_in = expected_membership_nary(operation, &in_operands);
438 report.considered += 1;
439 if expected_in != result_in {
440 report.disagreements.push(SemanticDisagreement {
441 point,
442 expected_in,
443 result_in,
444 in_a: in_operands[0],
447 in_b: *in_operands.get(1).unwrap_or(&false),
448 });
449 }
450 }
451
452 report.disagreement_rate = if report.considered == 0 {
453 0.0
454 } else {
455 report.disagreements.len() as f64 / report.considered as f64
456 };
457 Ok(report)
458}
459
460fn boundary_sample(rng: &mut Rng, faces: &[&crate::topology::FaceRecord], jitter: f64) -> Vec3 {
463 let face = faces[(rng.next_u64() as usize) % faces.len()];
464 let (Ok([u0, u1]), Ok([v0, v1])) = (face.surface.domain_u(), face.surface.domain_v()) else {
465 return Vec3::default();
466 };
467 let u = rng.range(u0, u1);
468 let v = rng.range(v0, v1);
469 let base = match face.surface.evaluate(u, v) {
470 Ok(point) => point,
471 Err(_) => return Vec3::default(),
472 };
473 let direction = match face.surface.normal(u, v) {
474 Ok(normal) if normal.length() > 1e-9 => normal,
475 _ => rng.unit_vector(),
476 };
477 let sign = if rng.unit() < 0.5 { -1.0 } else { 1.0 };
478 base.add(direction.scale(sign * jitter))
479}
480
481#[derive(Clone, Copy, Debug, PartialEq, Eq, Serialize)]
483#[serde(rename_all = "snake_case")]
484pub enum ResidualKind {
485 CoincidentVertices,
487 DuplicateEdge,
490}
491
492#[derive(Clone, Debug, Serialize)]
494pub struct ResidualFusable {
495 pub kind: ResidualKind,
496 pub detail: String,
497 pub point: Vec3,
498}
499
500pub fn boolean_residual_fusables(result: &BrepSolid, tol: f64) -> Vec<ResidualFusable> {
505 let mut findings = Vec::new();
506
507 let vertices = &result.vertices;
509 for i in 0..vertices.len() {
510 for j in (i + 1)..vertices.len() {
511 let gap = vertices[i].point.sub(vertices[j].point).length();
512 if gap <= tol {
513 findings.push(ResidualFusable {
514 kind: ResidualKind::CoincidentVertices,
515 detail: format!(
516 "vertices {} and {} coincide within {gap:.3e}",
517 vertices[i].id, vertices[j].id
518 ),
519 point: vertices[i].point,
520 });
521 }
522 }
523 }
524
525 let edges = &result.edges;
530 let endpoint =
531 |edge: &crate::topology::EdgeRecord| -> Result<(Vec3, Vec3, Vec3), KernelRefusal> {
532 let start = edge
533 .curve
534 .evaluate(edge.t0)
535 .or_refuse(KernelStage::Validate, "evaluate")?;
536 let end = edge
537 .curve
538 .evaluate(edge.t1)
539 .or_refuse(KernelStage::Validate, "evaluate")?;
540 let mid = edge
541 .curve
542 .evaluate((edge.t0 + edge.t1) * 0.5)
543 .or_refuse(KernelStage::Validate, "evaluate")?;
544 Ok((start, mid, end))
545 };
546 for i in 0..edges.len() {
547 let Ok((si, mi, ei)) = endpoint(&edges[i]) else {
548 continue;
549 };
550 for j in (i + 1)..edges.len() {
551 let Ok((sj, mj, ej)) = endpoint(&edges[j]) else {
552 continue;
553 };
554 let endpoints_match = (si.sub(sj).length() <= tol && ei.sub(ej).length() <= tol)
555 || (si.sub(ej).length() <= tol && ei.sub(sj).length() <= tol);
556 if endpoints_match && mi.sub(mj).length() <= tol {
557 findings.push(ResidualFusable {
558 kind: ResidualKind::DuplicateEdge,
559 detail: format!(
560 "edges {} and {} trace the same curve within {tol:.3e}",
561 edges[i].id, edges[j].id
562 ),
563 point: mi,
564 });
565 }
566 }
567 }
568
569 findings
570}
571
572