1use core::fmt;
4
5use axiolid_contracts::ExecutionOptions;
6use axiolid_core::{Frame3, Point3, Tolerance, Vec3};
7use axiolid_mesh::TriMesh;
8
9use crate::{MeshPlaneSection, SectionLimits};
10
11#[derive(Debug, Clone, PartialEq)]
12pub enum ConformanceFailure {
13 ReturnedError(String),
14 EmptyCentralSection,
15 IncorrectEvidence,
16 NonDeterministic,
17 WrongContourCount {
19 case: &'static str,
20 expected: usize,
21 actual: usize,
22 },
23 WrongSectionArea {
25 case: &'static str,
26 expected: f64,
27 actual: f64,
28 },
29 TangentPlaneHasArea {
32 case: &'static str,
33 area: f64,
34 },
35}
36
37#[derive(Debug, Clone, Default, PartialEq)]
38pub struct ConformanceReport {
39 pub failures: Vec<ConformanceFailure>,
40}
41
42impl ConformanceReport {
43 pub fn is_success(&self) -> bool {
44 self.failures.is_empty()
45 }
46}
47
48impl fmt::Display for ConformanceReport {
49 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
50 if self.is_success() {
51 write!(f, "conformant")
52 } else {
53 write!(f, "{:?}", self.failures)
54 }
55 }
56}
57
58pub struct ConformanceSuite;
59
60impl ConformanceSuite {
61 pub fn run(provider: &dyn MeshPlaneSection) -> ConformanceReport {
62 let mut report = ConformanceReport::default();
63 let mesh = unit_cube();
64 let frame = Frame3 {
65 origin: Point3::new(0.0, 0.0, 0.5),
66 x: Vec3::X,
67 y: Vec3::Y,
68 z: Vec3::Z,
69 };
70 let limits = SectionLimits::new(8, 12, 32, 4);
71 let options = ExecutionOptions::new(Tolerance::new(1e-9, 1e-9).expect("valid tolerance"));
72 let first = provider.section(&mesh, frame, limits, &options);
73 let second = provider.section(&mesh, frame, limits, &options);
74 match (first, second) {
75 (Ok(a), Ok(b)) => {
76 if a.contours.is_empty() {
77 report
78 .failures
79 .push(ConformanceFailure::EmptyCentralSection);
80 }
81 if a.evidence.source_triangles != mesh.triangles().len()
82 || !a.evidence.is_derived_from_input_mesh()
83 {
84 report.failures.push(ConformanceFailure::IncorrectEvidence);
85 }
86 if a != b {
87 report.failures.push(ConformanceFailure::NonDeterministic);
88 }
89 }
90 (Err(error), _) | (_, Err(error)) => report
91 .failures
92 .push(ConformanceFailure::ReturnedError(error.to_string())),
93 }
94
95 let sphere = icosphere(1.0, 3);
105 for spec in &[
106 SphereCase {
108 case: "sphere central",
109 height: 0.0,
110 contours: 1,
111 tolerance: 0.01,
112 },
113 SphereCase {
116 case: "sphere h=0.5",
117 height: 0.5,
118 contours: 1,
119 tolerance: 0.02,
120 },
121 SphereCase {
123 case: "sphere above pole",
124 height: 1.5,
125 contours: 0,
126 tolerance: 0.0,
127 },
128 SphereCase {
132 case: "sphere tangent pole",
133 height: 1.0,
134 contours: 0,
135 tolerance: 0.0,
136 },
137 ] {
138 check_sphere_section(provider, &mut report, &sphere, 1.0, spec);
139 }
140
141 report
142 }
143}
144
145fn unit_cube() -> TriMesh {
146 TriMesh::new(
147 vec![
148 Point3::new(0.0, 0.0, 0.0),
149 Point3::new(1.0, 0.0, 0.0),
150 Point3::new(1.0, 1.0, 0.0),
151 Point3::new(0.0, 1.0, 0.0),
152 Point3::new(0.0, 0.0, 1.0),
153 Point3::new(1.0, 0.0, 1.0),
154 Point3::new(1.0, 1.0, 1.0),
155 Point3::new(0.0, 1.0, 1.0),
156 ],
157 vec![
158 0, 2, 1, 0, 3, 2, 4, 5, 6, 4, 6, 7, 0, 1, 5, 0, 5, 4, 1, 2, 6, 1, 6, 5, 2, 3, 7, 2, 7,
159 6, 3, 0, 4, 3, 4, 7,
160 ],
161 )
162}
163
164fn contour_area(contour: &crate::SectionContour) -> f64 {
166 let p = &contour.points;
167 if p.len() < 3 {
168 return 0.0;
169 }
170 let mut twice = 0.0;
171 for i in 0..p.len() {
172 let a = p[i];
173 let b = p[(i + 1) % p.len()];
174 twice += a.x * b.y - b.x * a.y;
175 }
176 (twice / 2.0).abs()
177}
178
179fn total_area(contours: &[crate::SectionContour]) -> f64 {
181 contours.iter().map(contour_area).sum()
182}
183
184fn icosphere(radius: f64, subdivisions: u32) -> TriMesh {
192 let t = (1.0 + 5.0_f64.sqrt()) / 2.0;
193 let mut v: Vec<[f64; 3]> = vec![
194 [-1.0, t, 0.0],
195 [1.0, t, 0.0],
196 [-1.0, -t, 0.0],
197 [1.0, -t, 0.0],
198 [0.0, -1.0, t],
199 [0.0, 1.0, t],
200 [0.0, -1.0, -t],
201 [0.0, 1.0, -t],
202 [t, 0.0, -1.0],
203 [t, 0.0, 1.0],
204 [-t, 0.0, -1.0],
205 [-t, 0.0, 1.0],
206 ];
207 let mut f: Vec<[u32; 3]> = vec![
208 [0, 11, 5],
209 [0, 5, 1],
210 [0, 1, 7],
211 [0, 7, 10],
212 [0, 10, 11],
213 [1, 5, 9],
214 [5, 11, 4],
215 [11, 10, 2],
216 [10, 7, 6],
217 [7, 1, 8],
218 [3, 9, 4],
219 [3, 4, 2],
220 [3, 2, 6],
221 [3, 6, 8],
222 [3, 8, 9],
223 [4, 9, 5],
224 [2, 4, 11],
225 [6, 2, 10],
226 [8, 6, 7],
227 [9, 8, 1],
228 ];
229 for _ in 0..subdivisions {
230 let mut mid: std::collections::HashMap<(u32, u32), u32> = std::collections::HashMap::new();
231 let mut next: Vec<[u32; 3]> = Vec::with_capacity(f.len() * 4);
232 for tri in &f {
233 let mut m = [0u32; 3];
234 for e in 0..3 {
235 let (a, b) = (tri[e], tri[(e + 1) % 3]);
236 let key = (a.min(b), a.max(b));
237 m[e] = *mid.entry(key).or_insert_with(|| {
238 let (pa, pb) = (v[a as usize], v[b as usize]);
239 v.push([
240 (pa[0] + pb[0]) * 0.5,
241 (pa[1] + pb[1]) * 0.5,
242 (pa[2] + pb[2]) * 0.5,
243 ]);
244 (v.len() - 1) as u32
245 });
246 }
247 next.push([tri[0], m[0], m[2]]);
248 next.push([tri[1], m[1], m[0]]);
249 next.push([tri[2], m[2], m[1]]);
250 next.push([m[0], m[1], m[2]]);
251 }
252 f = next;
253 }
254 let positions: Vec<Point3> = v
255 .iter()
256 .map(|p| {
257 let l = (p[0] * p[0] + p[1] * p[1] + p[2] * p[2]).sqrt();
258 let k = radius / l;
259 Point3::new(p[0] * k, p[1] * k, p[2] * k)
260 })
261 .collect();
262 let indices: Vec<u32> = f.iter().flat_map(|t| [t[0], t[1], t[2]]).collect();
263 TriMesh::new(positions, indices)
264}
265
266struct SphereCase {
268 case: &'static str,
270 height: f64,
272 contours: usize,
274 tolerance: f64,
276}
277
278fn check_sphere_section(
286 provider: &dyn MeshPlaneSection,
287 report: &mut ConformanceReport,
288 sphere: &TriMesh,
289 radius: f64,
290 spec: &SphereCase,
291) {
292 let SphereCase {
293 case,
294 height,
295 contours: expect_contours,
296 tolerance: rel_tolerance,
297 } = *spec;
298 let frame = Frame3 {
299 origin: Point3::new(0.0, 0.0, height),
300 x: Vec3::X,
301 y: Vec3::Y,
302 z: Vec3::Z,
303 };
304 let limits = SectionLimits::new(4096, 8192, 16384, 64);
305 let options = ExecutionOptions::new(Tolerance::new(1e-9, 1e-9).expect("valid tolerance"));
306 match provider.section(sphere, frame, limits, &options) {
307 Err(error) => report
308 .failures
309 .push(ConformanceFailure::ReturnedError(error.to_string())),
310 Ok(outcome) => {
311 let closed = outcome.contours.len();
312 if closed != expect_contours {
313 report.failures.push(ConformanceFailure::WrongContourCount {
314 case,
315 expected: expect_contours,
316 actual: closed,
317 });
318 }
319 let area = total_area(&outcome.contours);
320 if expect_contours == 0 {
321 if area > 1e-9 {
324 report
325 .failures
326 .push(ConformanceFailure::TangentPlaneHasArea { case, area });
327 }
328 return;
329 }
330 let expected = core::f64::consts::PI * (radius * radius - height * height);
331 if expected > 0.0 && (area - expected).abs() / expected > rel_tolerance {
332 report.failures.push(ConformanceFailure::WrongSectionArea {
333 case,
334 expected,
335 actual: area,
336 });
337 }
338 }
339 }
340}