cellular_raza_concepts/
domain.rs

1use std::collections::{BTreeMap, BTreeSet};
2
3use crate::errors::{BoundaryError, DecomposeError};
4
5/// Provides an abstraction of the physical total simulation domain.
6///
7/// [cellular_raza](https://github.com/jonaspleyer/cellular_raza) uses domain-decomposition
8/// algorithms to split up the computational workload over multiple physical regions.
9/// That's why the domain itself is mostly responsible for being deconstructed
10/// into smaller [SubDomains](SubDomain) which can then be used to numerically solve our system.
11///
12/// This trait can be automatically implemented when the [SortCells], [DomainRngSeed],
13/// and [DomainCreateSubDomains] are satisfied together with a small number of trait bounds to hash
14/// and compare indices.
15pub trait Domain<C, S, Ci = Vec<C>> {
16    /// Subdomains can be identified by their unique [SubDomainIndex](Domain::SubDomainIndex).
17    /// The backend uses this property to construct a mapping (graph) between subdomains.
18    type SubDomainIndex;
19
20    /// Similarly to the [SubDomainIndex](Domain::SubDomainIndex), voxels can be accessed by
21    /// their unique index. The backend will use this information to construct a mapping
22    /// (graph) between voxels inside their respective subdomains.
23    type VoxelIndex;
24
25    /// Deconstructs the [Domain] into its respective subdomains.
26    ///
27    /// When using the blanket implementation of this function, the following steps are carried
28    /// out:
29    /// Its functionality consists of the following steps:
30    /// 1. Decompose the Domain into [Subdomains](SubDomain)
31    /// 2. Build a neighbor map between [SubDomains](SubDomain)
32    /// 3. Sort cells to their respective [SubDomain]
33    /// However, to increase performance or avoid trait bounds, one can also opt to implement this
34    /// trait directly.
35    fn decompose(
36        self,
37        n_subdomains: core::num::NonZeroUsize,
38        cells: Ci,
39    ) -> Result<DecomposedDomain<Self::SubDomainIndex, S, C>, DecomposeError>;
40}
41
42/// Manage the current rng seed of a [Domain]
43pub trait DomainRngSeed {
44    // fn set_rng_seed(&mut self, seed: u64);
45
46    /// Obtains the current rng seed
47    fn get_rng_seed(&self) -> u64;
48}
49
50/// Generate [SubDomains](SubDomain) from an existing [Domain]
51pub trait DomainCreateSubDomains<S> {
52    /// This should always be identical to [Domain::SubDomainIndex].
53    type SubDomainIndex;
54    /// This should always be identical to [Domain::VoxelIndex].
55    type VoxelIndex;
56
57    /// Generates at most `n_subdomains`. This function can also return a lower amount of
58    /// subdomains but never less than 1.
59    fn create_subdomains(
60        &self,
61        n_subdomains: core::num::NonZeroUsize,
62    ) -> Result<
63        impl IntoIterator<Item = (Self::SubDomainIndex, S, Vec<Self::VoxelIndex>)>,
64        DecomposeError,
65    >;
66}
67
68/// Generated by the [decompose](Domain::decompose) method. The backend will know how to
69/// deal with this type and crate a working simulation from it.
70pub struct DecomposedDomain<I, S, C> {
71    /// Number of spawned [SubDomains](SubDomain). This number is guaranteed to be
72    /// smaller or equal to the number may be different to the one given to the
73    /// [Domain::decompose] method.
74    /// Such behaviour can result from not being able to construct as many subdomains as desired.
75    /// Note that this function will attempt to construct more [SubDomains](SubDomain)
76    /// than available CPUs if given a larger number.
77    pub n_subdomains: core::num::NonZeroUsize,
78    /// Vector containing properties of individual [SubDomains](SubDomain).
79    /// Entries are [Domain::SubDomainIndex], [SubDomain], and a vector of cells.
80    // TODO can be use another iterator than Vec<(I, S, Vec<C>)>?
81    pub index_subdomain_cells: Vec<(I, S, Vec<C>)>,
82    /// Encapsulates how the subdomains are linked to each other.
83    /// Eg. two subdomains without any boundary will never appear in each others collection
84    /// of neighbors.
85    /// For the future, we might opt to change to an undirected graph rather than a [BTreeMap].
86    pub neighbor_map: BTreeMap<I, BTreeSet<I>>,
87    /// Initial seed of the simulation for random number generation.
88    pub rng_seed: u64,
89}
90
91/// Subdomains are produced by decomposing a [Domain] into multiple physical regions.
92///
93/// # Derivation
94/// ```
95/// # use cellular_raza_concepts::*;
96/// struct MySubDomain {
97///     x_min: f32,
98///     x_max: f32,
99///     n: usize,
100/// }
101///
102/// impl SubDomain for MySubDomain {
103///     type VoxelIndex = usize;
104///
105///     fn get_neighbor_voxel_indices(
106///         &self,
107///         voxel_index: &Self::VoxelIndex
108///     ) -> Vec<Self::VoxelIndex> {
109///         (voxel_index.saturating_sub(1)..voxel_index.saturating_add(1).min(self.n)+1)
110///             .filter(|k| k!=voxel_index)
111///             .collect()
112///     }
113///
114///     fn get_all_indices(&self) -> Vec<Self::VoxelIndex> {
115///         (0..self.n).collect()
116///     }
117/// }
118///
119/// #[derive(SubDomain)]
120/// struct MyNewSubDomain {
121///     #[Base]
122///     base: MySubDomain,
123/// }
124/// # let _my_sdm = MyNewSubDomain {
125/// #     base: MySubDomain {
126/// #         x_min: -20.0,
127/// #         x_max: -11.0,
128/// #         n: 20,
129/// #     }
130/// # };
131/// # assert_eq!(_my_sdm.get_all_indices(), (0..20).collect::<Vec<_>>());
132/// # assert_eq!(_my_sdm.get_neighbor_voxel_indices(&0), vec![1]);
133/// # assert_eq!(_my_sdm.get_neighbor_voxel_indices(&3), vec![2,4]);
134/// # assert_eq!(_my_sdm.get_neighbor_voxel_indices(&7), vec![6,8]);
135/// ```
136pub trait SubDomain {
137    /// Individual Voxels inside each subdomain can be accessed by this index.
138    type VoxelIndex;
139
140    /// Obtains the neighbor voxels of the specified voxel index. This function behaves similarly
141    /// to [SortCells::get_voxel_index_of] in that it also has to return
142    /// indices which are in other [SubDomains](SubDomain).
143    fn get_neighbor_voxel_indices(&self, voxel_index: &Self::VoxelIndex) -> Vec<Self::VoxelIndex>;
144
145    // fn apply_boundary(&self, cell: &mut C) -> Result<(), BoundaryError>;
146
147    /// Get all voxel indices of this [SubDomain].
148    fn get_all_indices(&self) -> Vec<Self::VoxelIndex>;
149}
150
151/// Assign an [VoxelIndex](SortCells::VoxelIndex) to a given cell.
152///
153/// This trait is used by the [Domain] and [SubDomain] trait to assign a [Domain::SubDomainIndex]
154/// and [SubDomain::VoxelIndex] respectively.
155///
156/// # [SubDomain]
157/// This trait is supposed to return the correct voxel index of the cell
158/// even if this index is inside another [SubDomain].
159/// This restriction might be lifted in the future but is still
160/// required now.
161pub trait SortCells<C> {
162    /// An index which determines to which next smaller unit the cell should be assigned.
163    type VoxelIndex;
164
165    /// If given a cell, we can sort this cell into the corresponding sub unit.
166    fn get_voxel_index_of(&self, cell: &C) -> Result<Self::VoxelIndex, BoundaryError>;
167}
168
169/// Apply boundary conditions to a cells position and velocity.
170///
171/// # Derivation
172/// ```
173/// # use cellular_raza_concepts::*;
174/// # use cellular_raza_concepts::BoundaryError;
175/// struct MyMechanics {
176///     x_min: f64,
177///     x_max: f64,
178/// }
179///
180/// impl SubDomainMechanics<f64, f64> for MyMechanics {
181///     fn apply_boundary(&self, pos: &mut f64, vel: &mut f64) -> Result<(), BoundaryError> {
182///         if *pos < self.x_min {
183///             *vel = vel.abs();
184///         }
185///         if *pos > self.x_max {
186///             *vel = -vel.abs();
187///         }
188///         *pos = pos.clamp(self.x_min, self.x_max);
189///         Ok(())
190///     }
191/// }
192///
193/// #[derive(SubDomain)]
194/// struct MySubDomain {
195///     #[Mechanics]
196///     mechanics: MyMechanics,
197/// }
198/// # let _my_sdm = MySubDomain {
199/// #     mechanics: MyMechanics {
200/// #         x_min: 1.0,
201/// #         x_max: 33.0,
202/// #     }
203/// # };
204/// # let mut pos = 0.0;
205/// # let mut vel = - 0.1;
206/// # _my_sdm.apply_boundary(&mut pos, &mut vel).unwrap();
207/// # assert_eq!(pos, 1.0);
208/// # assert_eq!(vel, 0.1);
209/// ```
210pub trait SubDomainMechanics<Pos, Vel> {
211    /// If the subdomain has boundary conditions, this function will enforce them onto the cells.
212    /// For the future, we plan to replace this function to additionally obtain information
213    /// about the previous and current location of the cell.
214    fn apply_boundary(&self, pos: &mut Pos, vel: &mut Vel) -> Result<(), BoundaryError>;
215}
216
217/// Apply a force on a cell depending on its position and velocity.
218///
219/// # Derivation
220/// ```
221/// # use cellular_raza_concepts::*;
222/// struct MyForce {
223///     damping: f64,
224/// }
225///
226/// impl SubDomainForce<f64, f64, f64> for MyForce {
227///     fn calculate_custom_force(&self, pos: &f64, vel: &f64) -> Result<f64, CalcError> {
228///         Ok(- self.damping * vel)
229///     }
230/// }
231///
232/// #[derive(SubDomain)]
233/// struct MySubDomain {
234///     #[Force]
235///     force: MyForce,
236/// }
237/// # let _my_sdm = MySubDomain {
238/// #     force: MyForce {
239/// #         damping: 0.1,
240/// #     }
241/// # };
242/// # let calculated_force = _my_sdm.calculate_custom_force(&0.0, &1.0).unwrap();
243/// # assert_eq!(calculated_force, -0.1);
244/// ```
245pub trait SubDomainForce<Pos, Vel, For> {
246    ///
247    fn calculate_custom_force(&self, pos: &Pos, vel: &Vel) -> Result<For, crate::CalcError>;
248}
249
250/// Describes extracellular reactions and fluid dynamics
251///
252/// # Derivation
253/// ```
254/// # use cellular_raza_concepts::*;
255///
256/// #[derive(Clone, Debug)]
257/// struct MyReactions<const N: usize> {
258///     values: Vec<f32>,
259///     pos: [f32; N],
260/// }
261///
262/// impl<const N: usize> SubDomainReactions<[f32; N], Vec<f32>, f32> for MyReactions<N> {
263///     type NeighborValue = Vec<f32>;
264///     type BorderInfo = Self;
265///
266///     fn treat_increments<I, J>(
267///         &mut self,
268///         neighbors: I,
269///         sources: J,
270///     ) -> Result<(), CalcError>
271///     where
272///         I: IntoIterator<Item = Self::NeighborValue>,
273///         J: IntoIterator<Item = ([f32; N], Vec<f32>)>,
274///     {
275///         Ok(())
276///     }
277///
278///     fn update_fluid_dynamics(&mut self, dt: f32) -> Result<(), CalcError> {
279///         Ok(())
280///     }
281///
282///     fn get_extracellular_at_pos(&self, pos: &[f32; N]) -> Result<Vec<f32>, CalcError> {
283///         Ok(self.values.clone())
284///     }
285///
286///     fn get_neighbor_value(&self, border_info: Self::BorderInfo) -> Self::NeighborValue {
287///         self.values.clone()
288///     }
289///
290///     fn get_border_info(&self) -> Self::BorderInfo {
291///         self.clone()
292///     }
293/// }
294///
295/// #[derive(SubDomain)]
296/// struct DerivedSubDomain<const N: usize> {
297///     #[Reactions]
298///     reactions: MyReactions<N>,
299/// }
300/// ```
301pub trait SubDomainReactions<Pos, Re, Float> {
302    /// Extracellular value of neighbor
303    type NeighborValue;
304    /// Exchanged information to locate neighboring subdomains.
305    type BorderInfo;
306
307    /// Combines increments which have been obtained by neighbors and cell-sources
308    fn treat_increments<I, J>(&mut self, neighbors: I, sources: J) -> Result<(), crate::CalcError>
309    where
310        I: IntoIterator<Item = Self::NeighborValue>,
311        J: IntoIterator<Item = (Pos, Re)>;
312
313    /// Main update function to calculate new values of extracellular concentrations.
314    fn update_fluid_dynamics(&mut self, dt: Float) -> Result<(), crate::CalcError>;
315
316    /// Obtain extracellular concentrations at given point.
317    fn get_extracellular_at_pos(&self, pos: &Pos) -> Result<Re, crate::CalcError>;
318    /// Obtains the [NeighborValue] which should be sent to the neighbor which has exposed the given
319    /// [BorderInfo].
320    fn get_neighbor_value(&self, border_info: Self::BorderInfo) -> Self::NeighborValue;
321    /// Obtains the [BorderInfo] used to retrieve the [NeighborValue].
322    fn get_border_info(&self) -> Self::BorderInfo;
323}
324
325/// This trait derives the different aspects of a [SubDomain].
326///
327/// It serves similarly as the [cellular_raza_concepts_derive::CellAgent] trait to quickly
328/// build new structures from already existing functionality.
329///
330/// | Attribute | Trait | Implemented |
331/// | ---  | --- |:---:|
332/// | `Base` | [SubDomain] | ✅ |
333/// | `SortCells` | [SortCells] | ✅ |
334/// | `Mechanics` | [SubDomainMechanics] | ✅ |
335/// | `Force` | [SubDomainForce] | ✅  |
336/// | `Reactions` | [SubDomainReactions] | ❌ |
337///
338/// # Example Usage
339/// ```
340/// # use cellular_raza_concepts::*;
341/// # struct MySubDomain;
342/// # impl SubDomain for MySubDomain {
343/// #     type VoxelIndex = usize;
344/// #     fn get_neighbor_voxel_indices(&self, voxel_index: &Self::VoxelIndex) -> Vec<usize> {
345/// #         Vec::new()
346/// #     }
347/// #     fn get_all_indices(&self) -> Vec<Self::VoxelIndex> {
348/// #         Vec::new()
349/// #     }
350/// # }
351/// #[derive(SubDomain)]
352/// struct MyDerivedSubDomain {
353///     #[Base]
354///     s: MySubDomain,
355/// }
356/// # let derived_subdomain = MyDerivedSubDomain {
357/// #     s: MySubDomain,
358/// # };
359/// # let all_indices = derived_subdomain.get_all_indices();
360/// # assert_eq!(all_indices.len(), 0);
361/// ```
362#[doc(inline)]
363pub use cellular_raza_concepts_derive::SubDomain;
364
365// TODO
366#[doc(inline)]
367pub use cellular_raza_concepts_derive::Domain;
cellular_raza_concepts/domain.rs

cellular_raza_concepts/
domain.rs
