Skip to main content

conspire/geometry/ntree/dual/quadtree/
mod.rs

1#[cfg(test)]
2mod test;
3
4use crate::{
5    geometry::{
6        Coordinates,
7        mesh::{Connectivity, Mesh},
8        ntree::{
9            Quadtree,
10            dual::{Dualization, Initialize, NodeMap, Star},
11            node::split::Split,
12        },
13    },
14    math::{Scalar, TensorVec},
15};
16use std::ops::Add;
17
18const D: usize = 2;
19const N: usize = 4;
20
21impl<T, U> Dualization<D> for Quadtree<T, U>
22where
23    T: Add<Output = T> + Copy + Into<Scalar> + Into<usize> + PartialOrd + Split,
24    U: Copy + Into<usize>,
25{
26    fn dualize(&mut self) -> Mesh<D> {
27        let (center_nodes, mut coordinates, mut node_index, mut connectivity) = self.initialize();
28        let mut nodes_map = NodeMap::new();
29        edge_transition(
30            self,
31            &center_nodes,
32            &mut coordinates,
33            &mut connectivity,
34            &mut node_index,
35            &mut nodes_map,
36        );
37        self.star(&center_nodes, &mut connectivity);
38        self.rescale_coordinates(&mut coordinates);
39        (
40            vec![Connectivity::Quadrilateral(connectivity.into())],
41            coordinates,
42        )
43            .into()
44    }
45}
46
47fn edge_transition<T, U>(
48    tree: &Quadtree<T, U>,
49    center_nodes: &[usize],
50    coordinates: &mut Coordinates<D>,
51    connectivity: &mut Vec<[usize; N]>,
52    node_index: &mut usize,
53    nodes_map: &mut NodeMap<D>,
54) where
55    T: Copy + Into<Scalar> + Into<usize>,
56    U: Copy + Into<usize>,
57{
58    let mut get_or_add = |pos: [Scalar; D]| -> usize {
59        let key = pos.map(|p| (2.0 * p) as usize);
60        if let Some(&v) = nodes_map.get(&key) {
61            v
62        } else {
63            let v = *node_index;
64            coordinates.push(pos.into());
65            nodes_map.insert(key, v);
66            *node_index += 1;
67            v
68        }
69    };
70    tree.iter().for_each(|node| {
71        let node_leaves = tree.leaves(node);
72        if let Some(neighbor) = node.facets()[0]
73            && let Some(leaf_0) = node_leaves[0]
74            && let Some(leaf_2) = node_leaves[2]
75        {
76            let face_leaves = tree.orthants_leaves_on_facet(&tree[neighbor], 1);
77            if let Some([Some(g_0a), Some(g_0b)]) = face_leaves[0]
78                && let Some([Some(g_2a), Some(g_2b)]) = face_leaves[1]
79            {
80                let length: Scalar = tree[g_2a].length.into();
81                let x0: Scalar = tree[g_2a].corner[0].into();
82                let y0c: Scalar = tree[g_2a].corner[1].into();
83                let x1 = x0 + length;
84                let y0 = y0c - length * 0.5;
85                let y1 = y0 + length;
86                let new_1 = get_or_add([x1, y0]);
87                let new_2 = get_or_add([x1, y1]);
88                connectivity.push([
89                    center_nodes[g_0b.into()],
90                    new_1,
91                    new_2,
92                    center_nodes[g_2a.into()],
93                ]);
94                connectivity.push([
95                    new_1,
96                    center_nodes[leaf_0.into()],
97                    center_nodes[leaf_2.into()],
98                    new_2,
99                ]);
100                connectivity.push([
101                    center_nodes[g_2a.into()],
102                    new_2,
103                    center_nodes[leaf_2.into()],
104                    center_nodes[g_2b.into()],
105                ]);
106                connectivity.push([
107                    center_nodes[g_0a.into()],
108                    center_nodes[leaf_0.into()],
109                    new_1,
110                    center_nodes[g_0b.into()],
111                ]);
112            }
113        }
114        if let Some(neighbor) = node.facets()[1]
115            && let Some(leaf_1) = node_leaves[1]
116            && let Some(leaf_3) = node_leaves[3]
117        {
118            let face_leaves = tree.orthants_leaves_on_facet(&tree[neighbor], 0);
119            if let Some([Some(g_1a), Some(g_1b)]) = face_leaves[0]
120                && let Some([Some(g_3a), Some(g_3b)]) = face_leaves[1]
121            {
122                let length: Scalar = tree[g_3a].length.into();
123                let x0: Scalar = tree[g_3a].corner[0].into();
124                let y0c: Scalar = tree[g_3a].corner[1].into();
125                let x1 = x0;
126                let y0 = y0c - length * 0.5;
127                let y1 = y0 + length;
128                let new_1 = get_or_add([x1, y0]);
129                let new_2 = get_or_add([x1, y1]);
130                connectivity.push([
131                    new_1,
132                    center_nodes[g_1b.into()],
133                    center_nodes[g_3a.into()],
134                    new_2,
135                ]);
136                connectivity.push([
137                    new_1,
138                    new_2,
139                    center_nodes[leaf_3.into()],
140                    center_nodes[leaf_1.into()],
141                ]);
142                connectivity.push([
143                    center_nodes[g_3a.into()],
144                    center_nodes[g_3b.into()],
145                    center_nodes[leaf_3.into()],
146                    new_2,
147                ]);
148                connectivity.push([
149                    center_nodes[g_1a.into()],
150                    center_nodes[g_1b.into()],
151                    new_1,
152                    center_nodes[leaf_1.into()],
153                ]);
154            }
155        }
156        if let Some(neighbor) = node.facets()[2]
157            && let Some(leaf_0) = node_leaves[0]
158            && let Some(leaf_1) = node_leaves[1]
159        {
160            let face_leaves = tree.orthants_leaves_on_facet(&tree[neighbor], 3);
161            if let Some([Some(g_0a), Some(g_0b)]) = face_leaves[0]
162                && let Some([Some(g_1a), Some(g_1b)]) = face_leaves[1]
163            {
164                let length: Scalar = tree[g_1a].length.into();
165                let x0c: Scalar = tree[g_1a].corner[0].into();
166                let y0: Scalar = tree[g_1a].corner[1].into();
167                let y1 = y0 + length;
168                let x0 = x0c - length * 0.5;
169                let x1 = x0 + length;
170                let new_1 = get_or_add([x0, y1]);
171                let new_2 = get_or_add([x1, y1]);
172                connectivity.push([
173                    center_nodes[g_0b.into()],
174                    center_nodes[g_1a.into()],
175                    new_2,
176                    new_1,
177                ]);
178                connectivity.push([
179                    new_1,
180                    new_2,
181                    center_nodes[leaf_1.into()],
182                    center_nodes[leaf_0.into()],
183                ]);
184                connectivity.push([
185                    center_nodes[g_1a.into()],
186                    center_nodes[g_1b.into()],
187                    center_nodes[leaf_1.into()],
188                    new_2,
189                ]);
190                connectivity.push([
191                    center_nodes[g_0a.into()],
192                    center_nodes[g_0b.into()],
193                    new_1,
194                    center_nodes[leaf_0.into()],
195                ]);
196            }
197        }
198        if let Some(neighbor) = node.facets()[3]
199            && let Some(leaf_2) = node_leaves[2]
200            && let Some(leaf_3) = node_leaves[3]
201        {
202            let face_leaves = tree.orthants_leaves_on_facet(&tree[neighbor], 2);
203            if let Some([Some(g_2a), Some(g_2b)]) = face_leaves[0]
204                && let Some([Some(g_3a), Some(g_3b)]) = face_leaves[1]
205            {
206                let length: Scalar = tree[g_3a].length.into();
207                let x0c: Scalar = tree[g_3a].corner[0].into();
208                let y0c: Scalar = tree[g_3a].corner[1].into();
209                let y1 = y0c;
210                let x0 = x0c - length * 0.5;
211                let x1 = x0 + length;
212                let new_1 = get_or_add([x0, y1]);
213                let new_2 = get_or_add([x1, y1]);
214                connectivity.push([
215                    new_1,
216                    new_2,
217                    center_nodes[g_3a.into()],
218                    center_nodes[g_2b.into()],
219                ]);
220                connectivity.push([
221                    new_1,
222                    center_nodes[leaf_2.into()],
223                    center_nodes[leaf_3.into()],
224                    new_2,
225                ]);
226                connectivity.push([
227                    center_nodes[g_3a.into()],
228                    new_2,
229                    center_nodes[leaf_3.into()],
230                    center_nodes[g_3b.into()],
231                ]);
232                connectivity.push([
233                    center_nodes[g_2a.into()],
234                    center_nodes[leaf_2.into()],
235                    new_1,
236                    center_nodes[g_2b.into()],
237                ]);
238            }
239        }
240    });
241}