Skip to main content

conspire/geometry/ntree/subdivide/
mod.rs

1use crate::geometry::ntree::{
2    Orthotree,
3    node::{Kind, split::Split},
4};
5use std::{array::from_fn, ops::Add};
6
7const fn mirror_facet(facet: usize) -> usize {
8    facet ^ 1
9}
10
11pub(crate) const fn insert_bit(x: usize, axis: usize, bit: usize) -> usize {
12    let low_mask = (1usize << axis) - 1;
13    let low = x & low_mask;
14    let high = x >> axis;
15    low | (bit << axis) | (high << (axis + 1))
16}
17
18impl<const D: usize, const L: usize, const M: usize, const N: usize, T, U, V>
19    Orthotree<D, L, M, N, T, U, V>
20where
21    T: Add<Output = T> + Copy + Split + Into<usize>,
22    U: Copy + From<usize> + Into<usize>,
23    V: Copy,
24{
25    fn nodes_on_face(facet: usize) -> [usize; L] {
26        from_fn(|k| insert_bit(k, facet / 2, facet % 2))
27    }
28    fn nodes_on_other_face(face: usize) -> [usize; L] {
29        Self::nodes_on_face(mirror_facet(face))
30    }
31    pub fn subdivide(&mut self, index: U) -> Result<(), &'static str> {
32        let indices = from_fn(|n| (self.len() + n).into());
33        let mut orthants = self[index].subdivide(indices)?;
34        for (facet, node_facet) in self[index].facets.into_iter().enumerate() {
35            if let Some(facet_node) = node_facet
36                && let Some(neighbors) = self[facet_node].orthants().copied()
37            {
38                for (node, neighbor) in Self::nodes_on_face(facet)
39                    .into_iter()
40                    .zip(Self::nodes_on_other_face(facet))
41                {
42                    if orthants[node].facets[facet].is_none() {
43                        orthants[node].facets[facet] = Some(neighbors[neighbor])
44                    } else {
45                        panic!("temporary to assess need for Option<>")
46                    }
47                    if self[neighbors[neighbor]].facets[mirror_facet(facet)].is_none() {
48                        self[neighbors[neighbor]].facets[mirror_facet(facet)] = Some(indices[node])
49                    } else {
50                        panic!("temporary to assess need for Option<>")
51                    }
52                }
53            }
54        }
55        self.extend(orthants);
56        self[index].kind = Kind::Tree(indices);
57        self[index].value = None;
58        Ok(())
59    }
60}