geometry.rs

1use super::models::Cell3d;2use crate::dim::geometry;3use crate::dim::models::{dtype_for, Cell2d};4use mrlycore::cell::remap;5use mrlycore::errors::{value_error, Result};6use mrlycore::tensor::Tensor;7use std::sync::OnceLock;89pub use crate::dim::geometry::{magic, mosaic, perforate};1011/// Returns the 24 rotation triples that reach each distinct cube orientation.12pub fn orientations() -> &'static Vec<(usize, usize, usize)> {13    static TABLE: OnceLock<Vec<(usize, usize, usize)>> = OnceLock::new();14    TABLE.get_or_init(|| {15        let mut probe = Tensor::new(vec![3, 3, 3]);16        for (flat, item) in probe.bytes_mut().iter_mut().enumerate() {17            *item = flat as u8;18        }19        let mut seen: Vec<Vec<u8>> = Vec::new();20        let mut table = Vec::new();21        for a in 0..4 {22            for b in 0..4 {23                for c in 0..4 {24                    let image = probe.rot90(a, (1, 2)).rot90(b, (0, 2)).rot90(c, (0, 1));25                    if !seen.contains(&image.bytes().to_vec()) {26                        seen.push(image.bytes().to_vec());27                        table.push((a, b, c));28                    }29                }30            }31        }32        table33    })34}3536/// Merges the cells into one cube arranged width by height by depth.37pub fn merge(cells: &[Cell3d], width: usize, height: usize, depth: usize) -> Result<Cell3d> {38    geometry::merge_reps(cells, &[height, width, depth])39}4041/// Orients a copy of the cell by each mask value and merges them in the mask's shape.42pub fn special(mask: &Tensor, cell: &Cell3d) -> Result<Cell3d> {43    if mask.shape.len() != 3 {44        return value_error("special mask must be 3d.");45    }46    if mask.bytes().iter().any(|&v| v > 23) {47        return value_error("Invalid orientation value. Must be 0..23.");48    }49    let oriented: Result<Vec<Cell3d>> = mask50        .bytes()51        .iter()52        .map(|&k| cell.clone().orient(k as usize))53        .collect();54    let oriented = oriented?;55    merge(&oriented, mask.shape[1], mask.shape[0], mask.shape[2])56}5758fn slice_map(cube: &[usize], axis: usize, index: usize) -> (Vec<usize>, Vec<usize>) {59    let strides = [cube[1] * cube[2], cube[2], 1];60    let kept: Vec<usize> = (0..3).filter(|&a| a != axis).collect();61    let shape = vec![cube[kept[0]], cube[kept[1]]];62    let mut map = Vec::with_capacity(shape[0] * shape[1]);63    for row in 0..shape[0] {64        for col in 0..shape[1] {65            let mut at = [0usize; 3];66            at[axis] = index;67            at[kept[0]] = row;68            at[kept[1]] = col;69            map.push(at[0] * strides[0] + at[1] * strides[1] + at[2]);70        }71    }72    (shape, map)73}7475/// Takes the flat cell left when one axis of the cube is fixed at an index, colors and tags with it.76///77/// ```78/// let sponge = mrlymath::three::carpet(3, 2).unwrap();79/// let front = mrlymath::three::slice(&sponge, 2, 0).unwrap();80/// assert_eq!(front, mrlymath::two::carpet(3, 2).unwrap());81/// ```82pub fn slice(cell: &Cell3d, axis: usize, index: usize) -> Result<Cell2d> {83    if axis > 2 {84        return value_error("slice axis is past the cube's rank.");85    }86    if index >= cell.types().shape[axis] {87        return value_error("slice index is past the axis.");88    }89    let (shape, map) = slice_map(&cell.types().shape, axis, index);90    Ok(Cell2d {91        cell: remap(&cell.cell, &map, &shape),92    })93}9495// EXTRUDE9697fn lift_map(shape: &[usize], axis: usize, depth: usize) -> (Vec<usize>, Vec<usize>) {98    let mut lifted = shape.to_vec();99    lifted.insert(axis, depth);100    let strides = [lifted[1] * lifted[2], lifted[2], 1];101    let (height, width) = (shape[0], shape[1]);102    let mut map = vec![0usize; lifted.iter().product()];103    for plane in 0..depth {104        for y in 0..height {105            for x in 0..width {106                let mut at = vec![y, x];107                at.insert(axis, plane);108                map[at[0] * strides[0] + at[1] * strides[1] + at[2]] = y * width + x;109            }110        }111    }112    (lifted, map)113}114115/// Lifts a flat cell into a cube by repeating it depth times along a new axis, colors and tags with it.116///117/// This is the inverse of slice: every slice of the lift on that axis is the flat cell again.118///119/// ```120/// let flat = mrlymath::two::carpet(3, 2).unwrap();121/// let cube = mrlymath::three::extrude(&flat, 2, 4).unwrap();122/// assert_eq!(cube.depth(), 4);123/// assert_eq!(mrlymath::three::slice(&cube, 2, 3).unwrap(), flat);124/// ```125pub fn extrude(cell: &Cell2d, axis: usize, depth: usize) -> Result<Cell3d> {126    if axis > 2 {127        return value_error("extrude axis must be 0, 1 or 2.");128    }129    if depth == 0 {130        return value_error("extrude depth must be at least 1.");131    }132    let (shape, map) = lift_map(&cell.types().shape, axis, depth);133    Ok(Cell3d {134        cell: remap(&cell.cell, &map, &shape),135    })136}137138// LAYERS139140/// Tags every site with its Manhattan distance from the cube's center, the diamond shells.141///142/// The cell's own layers count Chebyshev shells, which are boxes; these are octahedra.143pub fn manhattan_layers(mut cell: Cell3d) -> Cell3d {144    let shape = cell.types().shape.clone();145    let strides = [shape[1] * shape[2], shape[2], 1];146    let mut rings = vec![0i64; shape.iter().product()];147    for i in 0..shape[0] {148        for j in 0..shape[1] {149            for k in 0..shape[2] {150                let at = [i, j, k];151                let reach: f64 = (0..3)152                    .map(|axis| (at[axis] as f64 - (shape[axis] as f64 - 1.0) / 2.0).abs())153                    .sum();154                rings[i * strides[0] + j * strides[1] + k] = reach.floor() as i64;155            }156        }157    }158    let peak = rings.iter().copied().max().unwrap_or(0);159    let mut tags = Tensor::typed(shape, dtype_for(peak));160    for (flat, &ring) in rings.iter().enumerate() {161        tags.put(flat, ring);162    }163    cell.cell.tags = Some(tags);164    cell165}166167#[cfg(test)]168mod tests {169    use super::*;170    use crate::three::designs;171    use crate::two;172    #[test]173    fn exactly_24_orientations() {174        let table = orientations();175        assert_eq!(table.len(), 24);176        assert_eq!(table[0], (0, 0, 0));177    }178    #[test]179    fn orientations_preserve_sum_and_shape() {180        let c = designs::carpet(3, 1).unwrap();181        for i in 0..24 {182            let o = c.clone().orient(i).unwrap();183            assert_eq!(o.types().sum(), c.types().sum());184            assert_eq!(o.types().shape, c.types().shape);185        }186        assert!(c.clone().orient(24).is_err());187    }188    #[test]189    fn orientations_are_distinct_on_chiral_design() {190        let tree = designs::xtree(3, 1).unwrap();191        let images: Vec<Vec<u8>> = (0..24)192            .map(|i| tree.clone().orient(i).unwrap().types().bytes().to_vec())193            .collect();194        let mut unique = images.clone();195        unique.sort();196        unique.dedup();197        assert!(unique.len() >= 3);198    }199    #[test]200    fn special_identity_is_tile() {201        let c = designs::carpet(3, 1).unwrap();202        let mask = Tensor::new(vec![2, 2, 2]);203        let s = special(&mask, &c).unwrap();204        assert_eq!(s, c.clone().tile(2, 2, 2));205    }206    #[test]207    fn only_the_carpet_and_two_trees_face_a_flat_name() {208        for n in [3, 5] {209            for level in [1, 2] {210                let carpet = slice(&designs::carpet(n, level).unwrap(), 2, 0).unwrap();211                let net = slice(&designs::net(n, level).unwrap(), 2, 0).unwrap();212                let xtree = slice(&designs::xtree(n, level).unwrap(), 2, 0).unwrap();213                let ytree = slice(&designs::ytree(n, level).unwrap(), 2, 0).unwrap();214                let ztree = slice(&designs::ztree(n, level).unwrap(), 2, 0).unwrap();215                let void = slice(&designs::void(n, level).unwrap(), 2, 0).unwrap();216                assert_eq!(carpet, two::carpet(n, level).unwrap());217                assert_eq!(xtree, two::vtree(n, level).unwrap());218                assert_eq!(ytree, two::htree(n, level).unwrap());219                assert_eq!(void, ztree);220                assert_ne!(net, two::net(n, level).unwrap());221                assert_ne!(void, two::void(n, level).unwrap());222                assert_ne!(ztree, two::htree(n, level).unwrap());223                assert_ne!(ztree, two::vtree(n, level).unwrap());224            }225        }226    }227    #[test]228    fn extrude_undoes_slice_on_every_axis() {229        use mrlycore::cell::mapping;230        use mrlycore::enums::Mode;231        let flat = two::carpet(3, 2)232            .unwrap()233            .layers()234            .paint(&mapping(), Mode::Index);235        for axis in 0..3 {236            let cube = extrude(&flat, axis, 4).unwrap();237            assert_eq!(cube.types().shape[axis], 4);238            for index in 0..4 {239                assert_eq!(slice(&cube, axis, index).unwrap(), flat);240            }241            assert_eq!(cube.types().sum(), 4 * flat.types().sum());242        }243        assert!(extrude(&flat, 3, 2).is_err());244        assert!(extrude(&flat, 0, 0).is_err());245        assert!(slice(&extrude(&flat, 2, 1).unwrap(), 3, 0).is_err());246        assert!(slice(&extrude(&flat, 2, 1).unwrap(), 2, 1).is_err());247    }248    #[test]249    fn extrude_lifts_a_flat_face_of_a_cube_back() {250        let sponge = designs::carpet(3, 2).unwrap();251        let face = slice(&sponge, 2, 0).unwrap();252        let column = extrude(&face, 2, sponge.depth()).unwrap();253        assert_eq!(column.types().shape, sponge.types().shape);254        assert_eq!(slice(&column, 2, 5).unwrap(), face);255    }256    #[test]257    fn extrude_carries_colors_and_tags() {258        use mrlycore::cell::mapping;259        use mrlycore::enums::Mode;260        let flat = two::carpet(3, 1)261            .unwrap()262            .layers()263            .paint(&mapping(), Mode::Type);264        let cube = extrude(&flat, 0, 3).unwrap();265        let colors = cube.cell.colors.as_ref().unwrap();266        let tags = cube.cell.tags.as_ref().unwrap();267        let flat_colors = flat.cell.colors.as_ref().unwrap();268        let flat_tags = flat.cell.tags.as_ref().unwrap();269        assert_eq!(colors.len(), 27);270        assert_eq!(tags.shape, vec![3, 3, 3]);271        for plane in 0..3 {272            for site in 0..9 {273                assert_eq!(colors[plane * 9 + site], flat_colors[site]);274                assert_eq!(tags.at(plane * 9 + site), flat_tags.at(site));275            }276        }277    }278    #[test]279    fn manhattan_layers_are_diamond_shells() {280        let cube = manhattan_layers(designs::ones(3, 1).unwrap());281        let tags = cube.cell.tags.as_ref().unwrap();282        assert_eq!(tags.get(&[1, 1, 1]), 0);283        assert_eq!(tags.get(&[0, 1, 1]), 1);284        assert_eq!(tags.get(&[0, 0, 1]), 2);285        assert_eq!(tags.get(&[0, 0, 0]), 3);286        let boxes = designs::ones(3, 1).unwrap().layers();287        assert_eq!(boxes.cell.tags.as_ref().unwrap().get(&[0, 0, 0]), 1);288    }289    #[test]290    fn manhattan_layers_widen_past_a_byte() {291        use mrlycore::tensor::Dtype;292        let small = manhattan_layers(designs::ones(3, 1).unwrap());293        assert_eq!(small.cell.tags.as_ref().unwrap().dtype(), Dtype::U8);294        let long = manhattan_layers(Cell3d::new(Tensor::full(vec![1, 1, 600], 1)));295        let tags = long.cell.tags.as_ref().unwrap();296        assert_eq!(tags.dtype(), Dtype::U16);297        assert_eq!(tags.at(0), 299);298    }299}