Skip to main content

axiolid_reference/
primitive.rs

1//! Tessellation of CSG primitives into closed solids.
2//!
3//! Every variant is analytic, so each is emitted directly. A block has no
4//! curvature to sample; a sphere's is exactly known. The chord budget still
5//! decides the radial segment count so a primitive and a swept face of the
6//! same radius agree on what a tolerance means.
7
8use axiolid_contracts::{GeomError, GeomResult};
9use axiolid_core::{Point3, Scalar, Tolerance};
10use axiolid_mesh::TriMesh;
11use axiolid_primitive::Primitive;
12
13/// Radial segments needed for `radius` under a chord budget.
14///
15/// Same sagitta rule the curve flattener uses: r(1 - cos(pi/n)) <= tol.
16/// Clamped low so a degenerate tolerance cannot ask for an unbounded mesh.
17fn segments(radius: Scalar, tolerance: Scalar) -> usize {
18    if !(radius.is_finite() && radius > 0.0 && tolerance.is_finite() && tolerance > 0.0) {
19        return 3;
20    }
21    let ratio = 1.0 - (tolerance / radius).min(1.0);
22    let n = (core::f64::consts::PI / ratio.acos().max(1e-9)).ceil();
23    (n as usize).clamp(3, 4096)
24}
25
26/// Tessellate one CSG primitive into a closed, outward-wound solid.
27///
28/// Outward winding is not decoration: `axiolid-mesh-boolean-boolmesh` and the clash
29/// containment test both read signed volume, and an inverted primitive
30/// silently produces negative volume and wrong verdicts.
31pub fn tessellate_primitive(primitive: &Primitive, tolerance: Tolerance) -> GeomResult<TriMesh> {
32    let tol = tolerance.linear();
33    match primitive {
34        Primitive::Block { x, y, z } => block(*x, *y, *z),
35        Primitive::Sphere { radius } => sphere(*radius, tol),
36        Primitive::Cylinder { radius, height } => cylinder(*radius, *height, tol),
37        Primitive::Cone { radius, height } => cone(*radius, *height, tol),
38        Primitive::Pyramid { x, y, height } => pyramid(*x, *y, *height),
39        // The enum is non_exhaustive: a new family is unsupported, never
40        // silently approximated by the nearest one.
41        _ => Err(GeomError::Unsupported {
42            backend: crate::ScalarBoolean::ID,
43            operation: axiolid_contracts::Operation::Tessellation,
44        }),
45    }
46}
47
48/// Validate a positive finite extent, naming the offender.
49fn positive(value: Scalar, what: &str) -> GeomResult<Scalar> {
50    if !value.is_finite() || value <= 0.0 {
51        return Err(GeomError::InvalidInput(format!(
52            "{what} must be positive and finite, got {value}"
53        )));
54    }
55    Ok(value)
56}
57
58/// Axis-aligned block centred on the local origin.
59fn block(x: Scalar, y: Scalar, z: Scalar) -> GeomResult<TriMesh> {
60    let (hx, hy, hz) = (
61        positive(x, "block x")? / 2.0,
62        positive(y, "block y")? / 2.0,
63        positive(z, "block z")? / 2.0,
64    );
65    let p = vec![
66        Point3::new(-hx, -hy, -hz),
67        Point3::new(hx, -hy, -hz),
68        Point3::new(hx, hy, -hz),
69        Point3::new(-hx, hy, -hz),
70        Point3::new(-hx, -hy, hz),
71        Point3::new(hx, -hy, hz),
72        Point3::new(hx, hy, hz),
73        Point3::new(-hx, hy, hz),
74    ];
75    // Outward winding, verified by the positive-volume test rather than by
76    // reading the index list.
77    let i = vec![
78        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, 6,
79        3, 0, 4, 3, 4, 7,
80    ];
81    Ok(TriMesh::new(p, i))
82}
83
84/// Rectangular pyramid: base on z = 0, apex on +z.
85fn pyramid(x: Scalar, y: Scalar, height: Scalar) -> GeomResult<TriMesh> {
86    let (hx, hy) = (
87        positive(x, "pyramid x")? / 2.0,
88        positive(y, "pyramid y")? / 2.0,
89    );
90    let h = positive(height, "pyramid height")?;
91    let p = vec![
92        Point3::new(-hx, -hy, 0.0),
93        Point3::new(hx, -hy, 0.0),
94        Point3::new(hx, hy, 0.0),
95        Point3::new(-hx, hy, 0.0),
96        Point3::new(0.0, 0.0, h),
97    ];
98    let i = vec![0, 2, 1, 0, 3, 2, 0, 1, 4, 1, 2, 4, 2, 3, 4, 3, 0, 4];
99    Ok(TriMesh::new(p, i))
100}
101
102/// A ring of `n` points at `radius`, height `z`.
103fn ring(radius: Scalar, z: Scalar, n: usize) -> Vec<Point3> {
104    (0..n)
105        .map(|k| {
106            let a = core::f64::consts::TAU * (k as Scalar) / (n as Scalar);
107            Point3::new(radius * a.cos(), radius * a.sin(), z)
108        })
109        .collect()
110}
111
112/// Cylinder along +z, base on z = 0.
113fn cylinder(radius: Scalar, height: Scalar, tol: Scalar) -> GeomResult<TriMesh> {
114    let r = positive(radius, "cylinder radius")?;
115    let h = positive(height, "cylinder height")?;
116    let n = segments(r, tol);
117    let mut p = ring(r, 0.0, n);
118    p.extend(ring(r, h, n));
119    p.push(Point3::new(0.0, 0.0, 0.0));
120    p.push(Point3::new(0.0, 0.0, h));
121    let (bc, tc) = (2 * n, 2 * n + 1);
122    let mut i = Vec::with_capacity(n * 12);
123    for k in 0..n {
124        let (a, b) = (k, (k + 1) % n);
125        // Side quad, then the two caps. The base fan is wound opposite to
126        // the top so both face away from the enclosed volume.
127        i.extend([a as u32, b as u32, (b + n) as u32]);
128        i.extend([a as u32, (b + n) as u32, (a + n) as u32]);
129        i.extend([bc as u32, b as u32, a as u32]);
130        i.extend([tc as u32, (a + n) as u32, (b + n) as u32]);
131    }
132    Ok(TriMesh::new(p, i))
133}
134
135/// Cone along +z: base ring on z = 0, apex at height.
136fn cone(radius: Scalar, height: Scalar, tol: Scalar) -> GeomResult<TriMesh> {
137    let r = positive(radius, "cone radius")?;
138    let h = positive(height, "cone height")?;
139    let n = segments(r, tol);
140    let mut p = ring(r, 0.0, n);
141    p.push(Point3::new(0.0, 0.0, 0.0));
142    p.push(Point3::new(0.0, 0.0, h));
143    let (base, apex) = (n, n + 1);
144    let mut i = Vec::with_capacity(n * 6);
145    for k in 0..n {
146        let (a, b) = (k as u32, ((k + 1) % n) as u32);
147        i.extend([base as u32, b, a]);
148        i.extend([a, b, apex as u32]);
149    }
150    Ok(TriMesh::new(p, i))
151}
152
153/// Sphere centred on the local origin, as a UV mesh.
154fn sphere(radius: Scalar, tol: Scalar) -> GeomResult<TriMesh> {
155    let r = positive(radius, "sphere radius")?;
156    let n = segments(r, tol);
157    // Half as many stacks as segments: the polar direction spans PI, not TAU,
158    // so equal counts would oversample it by 2x for the same chord error.
159    let stacks = (n / 2).max(2);
160    let mut p = Vec::with_capacity((stacks + 1) * n);
161    for i in 0..=stacks {
162        let v = core::f64::consts::PI * (i as Scalar) / (stacks as Scalar);
163        for j in 0..n {
164            let u = core::f64::consts::TAU * (j as Scalar) / (n as Scalar);
165            p.push(Point3::new(
166                r * v.sin() * u.cos(),
167                r * v.sin() * u.sin(),
168                r * v.cos(),
169            ));
170        }
171    }
172    let mut idx = Vec::with_capacity(stacks * n * 6);
173    for i in 0..stacks {
174        for j in 0..n {
175            let jn = (j + 1) % n;
176            // A pole row's n entries are the same point, so they must map
177            // to ONE index or no triangle around the pole shares an edge.
178            let row = |r: usize, k: usize| -> u32 {
179                if r == 0 || r == stacks {
180                    (r * n) as u32
181                } else {
182                    (r * n + k) as u32
183                }
184            };
185            let a = row(i, j);
186            let b = row(i, jn);
187            let c = row(i + 1, j);
188            let d = row(i + 1, jn);
189            // Pole rows collapse to a point; skip the degenerate half.
190            if i > 0 {
191                idx.extend([a, c, b]);
192            }
193            if i + 1 < stacks {
194                idx.extend([b, c, d]);
195            }
196        }
197    }
198    Ok(TriMesh::new(p, idx))
199}