1use cellular_raza_concepts::*;
2
3use serde::{Deserialize, Serialize};
4
5#[derive(Clone, Deserialize, Serialize)]
7pub struct SimulationSetup<C, D> {
8 pub cells: Vec<C>,
11 pub domain: D,
13}
14
15#[derive(Clone, Deserialize, Serialize)]
17pub struct Settings<T, const INIT: bool> {
18 pub n_threads: core::num::NonZeroUsize,
20 pub time: T,
23 pub storage: crate::storage::StorageBuilder<INIT>,
25 pub progressbar: Option<String>,
27}
28
29impl<C, D> SimulationSetup<C, D> {
30 pub fn insert_cells<I>(&mut self, cells: I)
32 where
33 I: IntoIterator<Item = C>,
34 {
35 self.cells.extend(cells);
36 }
37
38 pub fn decompose<S>(
40 self,
41 n_subdomains: core::num::NonZeroUsize,
42 ) -> Result<DecomposedDomain<D::SubDomainIndex, S, C>, DecomposeError>
43 where
44 D: Domain<C, S>,
45 S: SubDomain,
46 {
47 self.domain.decompose(n_subdomains, self.cells)
48 }
49
50 pub fn decompose_auto_tune<S>(
56 self,
57 ) -> Result<DecomposedDomain<D::SubDomainIndex, S, C>, DecomposeError>
58 where
59 D: Domain<C, S>,
60 S: SubDomain,
61 {
62 todo!();
63 }
66}
67
68#[cfg(test)]
69mod test {
70 use super::*;
71
72 #[test]
73 fn test_config_insert_cells() {
74 let mut config = SimulationSetup {
75 cells: vec![1, 2, 3, 4],
76 domain: "hello".to_owned(),
77 };
78
79 let new_cells = vec![5, 6, 7];
80 config.insert_cells(new_cells);
81 assert_eq!(config.cells, vec![1, 2, 3, 4, 5, 6, 7]);
82 }
83
84 struct TestDomain {
85 min: f64,
86 max: f64,
87 n_voxels: usize,
88 rng_seed: u64,
89 }
90
91 #[allow(unused)]
92 struct VoxelIndex(usize);
93
94 struct TestSubDomain {
95 min: f64,
96 max: f64,
97 voxels: std::collections::BTreeMap<usize, [f64; 2]>,
98 reflect_at_boundary: (bool, bool),
99 total_voxels: usize,
100 }
101
102 impl Domain<f64, TestSubDomain> for TestDomain {
103 type SubDomainIndex = usize;
104 type VoxelIndex = VoxelIndex;
105
106 fn decompose(
107 self,
108 n_subdomains: core::num::NonZeroUsize,
109 cells: Vec<f64>,
110 ) -> Result<DecomposedDomain<Self::SubDomainIndex, TestSubDomain, f64>, DecomposeError>
111 {
112 let mut cells = cells;
113 let mut index_subdomain_cells = Vec::new();
114 let n_subdomains = n_subdomains.get();
115 let n_return_subdomains = n_subdomains.min(self.n_voxels);
116 let dx = (self.max - self.min) / (self.n_voxels as f64);
117
118 let voxel_distributions = (0..self.n_voxels)
119 .map(|i| (i, (i * n_return_subdomains).div_euclid(self.n_voxels)))
120 .fold(std::collections::BTreeMap::new(), |mut acc, x| {
121 let entry = acc.entry(x.1).or_insert(Vec::new());
122 entry.push(x.0);
123 acc
124 });
125
126 let mut n_total_voxels = 0;
127 for (subdomain_index, voxel_indices) in voxel_distributions {
128 let lower = if subdomain_index == 0 {
129 self.min
130 } else {
131 self.min + n_total_voxels as f64 * dx
132 };
133 let upper = if subdomain_index == n_subdomains - 1 {
134 self.max
135 } else {
136 self.min + (n_total_voxels + voxel_indices.len()) as f64 * dx
137 };
138 n_total_voxels += voxel_indices.len();
139 let (cells_in_subdomain, other_cells): (Vec<_>, Vec<_>) =
140 cells.into_iter().partition(|&x| lower <= x && x < upper);
141 cells = other_cells;
142
143 index_subdomain_cells.push((
144 subdomain_index,
145 TestSubDomain {
146 min: lower,
147 max: upper,
148 voxels: voxel_indices
149 .into_iter()
150 .map(|voxel_index| {
151 (
152 voxel_index,
153 [
154 self.min + voxel_index as f64 * dx,
155 self.min + (voxel_index + 1) as f64 * dx,
156 ],
157 )
158 })
159 .collect(),
160 reflect_at_boundary: (
161 subdomain_index == 0,
162 subdomain_index == n_subdomains - 1,
163 ),
164 total_voxels: self.n_voxels,
165 },
166 cells_in_subdomain,
167 ));
168 }
169
170 let n_subdomains = index_subdomain_cells.len();
171 let decomposed_domain = DecomposedDomain {
172 n_subdomains: n_subdomains.try_into().unwrap(),
173 index_subdomain_cells,
174 neighbor_map: (0..n_subdomains)
175 .map(|i| {
176 (
177 i,
178 std::collections::BTreeSet::from([
179 if i == 0 { n_subdomains } else { i - 1 },
180 i + 1,
181 ]),
182 )
183 })
184 .collect(),
185 rng_seed: self.rng_seed,
186 };
187 Ok(decomposed_domain)
188 }
189 }
190
191 impl SubDomain for TestSubDomain {
192 type VoxelIndex = usize;
193
194 fn get_all_indices(&self) -> Vec<Self::VoxelIndex> {
195 self.voxels.iter().map(|(&i, _)| i).collect()
196 }
197
198 fn get_neighbor_voxel_indices(
199 &self,
200 voxel_index: &Self::VoxelIndex,
201 ) -> Vec<Self::VoxelIndex> {
202 let mut neighbors = Vec::new();
203 if voxel_index > &0 && voxel_index <= &(self.total_voxels - 1) {
204 neighbors.push(voxel_index - 1);
205 }
206 if voxel_index >= &0 && voxel_index < &(self.total_voxels - 1) {
207 neighbors.push(voxel_index + 1);
208 }
209 neighbors
210 }
211 }
212
213 impl SortCells<f64> for TestSubDomain {
214 type VoxelIndex = usize;
215
216 fn get_voxel_index_of(
217 &self,
218 cell: &f64,
219 ) -> Result<Self::VoxelIndex, cellular_raza_concepts::BoundaryError> {
220 for (index, voxel) in self.voxels.iter() {
221 if cell >= &voxel[0] && cell <= &voxel[1] {
222 return Ok(*index);
223 }
224 }
225 Err(cellular_raza_concepts::BoundaryError(
226 "Could not find voxel which contains cell".into(),
227 ))
228 }
229 }
230
231 impl SubDomainMechanics<f64, f64> for TestSubDomain {
232 fn apply_boundary(
233 &self,
234 pos: &mut f64,
235 _vel: &mut f64,
236 ) -> Result<(), cellular_raza_concepts::BoundaryError> {
237 if self.reflect_at_boundary.0 && *pos < self.min {
238 *pos = self.min;
239 } else if self.reflect_at_boundary.1 && *pos > self.max {
240 *pos = self.max;
241 }
242 Ok(())
243 }
244 }
245
246 #[test]
247 fn test_config_to_subdomains() {
248 let min = 0.0;
250 let max = 100.0;
251 let config = SimulationSetup {
252 cells: vec![1.0, 20.0, 26.0, 41.0, 56.0, 84.0, 95.0],
253 domain: TestDomain {
254 min,
255 max,
256 n_voxels: 8,
257 rng_seed: 0,
258 },
259 };
260
261 let n_subdomains = core::num::NonZeroUsize::new(4).unwrap();
262 let decomposed_domain = config.decompose(n_subdomains).unwrap();
263 for (_, subdomain, cells) in decomposed_domain.index_subdomain_cells.into_iter() {
264 assert!(!cells.is_empty());
265 for cell in cells {
267 assert!(cell >= subdomain.min);
268 assert!(cell <= subdomain.max);
269 }
270 }
271 }
272
273 #[test]
274 fn test_apply_boundary() {
275 let min = 0.0;
277 let max = 100.0;
278 let config = SimulationSetup {
279 cells: vec![1.0, 20.0, 25.0, 50.0, 88.0],
280 domain: TestDomain {
281 min,
282 max,
283 n_voxels: 8,
284 rng_seed: 0,
285 },
286 };
287
288 let n_subdomains = core::num::NonZeroUsize::new(4).unwrap();
289 let decomposed_domain = config.decompose(n_subdomains).unwrap();
290 let mut cell_outside = -10.0;
291 for (_, subdomain, cells) in decomposed_domain.index_subdomain_cells.into_iter() {
292 for cell in cells.into_iter() {
293 let mut cell_prev = cell;
294 let mut _nothing = cell;
295 subdomain
296 .apply_boundary(&mut cell_prev, &mut _nothing)
297 .unwrap();
298 assert_eq!(cell_prev, cell);
299 }
300 let mut _nothing = cell_outside;
301 subdomain
302 .apply_boundary(&mut cell_outside, &mut _nothing)
303 .unwrap();
304 assert!(cell_outside >= min);
305 assert!(cell_outside <= max);
306 }
307 }
308
309 #[test]
310 fn test_neighbor_indices() {
311 let min = 0.0;
313 let max = 100.0;
314 let config = SimulationSetup {
315 cells: vec![1.0, 20.0, 26.0, 41.0, 56.0, 84.0, 95.0],
316 domain: TestDomain {
317 min,
318 max,
319 n_voxels: 8,
320 rng_seed: 0,
321 },
322 };
323
324 let n_subdomains = core::num::NonZeroUsize::new(4).unwrap();
325 let decomposed_domain = config.decompose(n_subdomains).unwrap();
326 for (_, subdomain, _) in decomposed_domain.index_subdomain_cells.into_iter() {
327 assert_eq!(vec![1], subdomain.get_neighbor_voxel_indices(&0));
328 assert_eq!(vec![0, 2], subdomain.get_neighbor_voxel_indices(&1));
329 assert_eq!(vec![1, 3], subdomain.get_neighbor_voxel_indices(&2));
330 assert_eq!(vec![2, 4], subdomain.get_neighbor_voxel_indices(&3));
331 assert_eq!(vec![3, 5], subdomain.get_neighbor_voxel_indices(&4));
332 assert_eq!(vec![4, 6], subdomain.get_neighbor_voxel_indices(&5));
333 assert_eq!(vec![5, 7], subdomain.get_neighbor_voxel_indices(&6));
334 assert_eq!(vec![6], subdomain.get_neighbor_voxel_indices(&7));
335 assert_eq!(vec![0_usize; 0], subdomain.get_neighbor_voxel_indices(&8));
336 assert_eq!(vec![0_usize; 0], subdomain.get_neighbor_voxel_indices(&9));
337 assert_eq!(vec![0_usize; 0], subdomain.get_neighbor_voxel_indices(&10));
338 }
339 }
340}