Skip to main content

conspire/geometry/grid/from/tessellation/
mod.rs

1#[cfg(test)]
2mod test;
3
4use crate::{
5    geometry::{Coordinate, grid::Voxels, mesh::Tessellation},
6    math::Quantity,
7    units::Length,
8};
9use std::{
10    array::from_fn,
11    thread::{available_parallelism, scope},
12};
13
14impl Voxels<usize> {
15    pub fn from_tessellation(tessellation: &Tessellation, size: Quantity<Length>) -> Self {
16        let mesh = tessellation.mesh();
17        let bvh = tessellation.bvh();
18        let elements: Vec<&[usize]> = mesh.connectivities().iter().flatten().collect();
19        let coordinates = mesh.coordinates();
20        let mut min = [Quantity::<Length>::new(f64::INFINITY); 3];
21        let mut max = [Quantity::<Length>::new(f64::NEG_INFINITY); 3];
22        for point in coordinates {
23            (0..3).for_each(|ax| {
24                min[ax] = min[ax].min(point[ax]);
25                max[ax] = max[ax].max(point[ax]);
26            });
27        }
28        let nel: [usize; 3] =
29            from_fn(|ax| ((((max[ax] - min[ax]) / size).ceil().value()) as usize).max(1));
30        let [nx, ny, _] = nel;
31        let layer = nx * ny;
32        let mut data = vec![0usize; layer * nel[2]];
33        let direction = Coordinate::from([1.0, 0.01, 0.001]);
34        let threads = available_parallelism().map_or(1, |threads| threads.get());
35        let chunk_size = data.len().div_ceil(threads).max(1);
36        scope(|scope| {
37            let (elements, direction) = (&elements, &direction);
38            data.chunks_mut(chunk_size)
39                .enumerate()
40                .for_each(|(chunk, voxels)| {
41                    let offset = chunk * chunk_size;
42                    scope.spawn(move || {
43                        voxels.iter_mut().enumerate().for_each(|(local, voxel)| {
44                            let flat = offset + local;
45                            let index = [flat % nx, flat / nx % ny, flat / layer];
46                            let center = Coordinate::from(from_fn::<_, 3, _>(|ax| {
47                                min[ax] + (index[ax] as f64 + 0.5) * size
48                            }));
49                            let ray = (center, direction.clone()).into();
50                            if bvh.intersections(&ray, coordinates, elements) % 2 == 1 {
51                                *voxel = 1;
52                            }
53                        });
54                    });
55                });
56        });
57        Voxels::new(data, nel)
58    }
59}