magba/sources/
source.rs

1/*
2 * Magba is licensed under The 3-Clause BSD, see LICENSE.
3 * Copyright 2025 Sira Pornsiriprasert <code@psira.me>
4 */
5
6//! Core traits and types for magnetic sources and collections of sources.
7
8use std::fmt::{Debug, Display};
9
10use nalgebra::{Point3, RealField, Translation3, UnitQuaternion, Vector3};
11
12use crate::{crate_util, geometry::Transform, Float};
13
14/// Trait shared by objects that generate magnetic field.
15#[allow(non_snake_case)]
16pub trait Field<T: RealField + Copy> {
17    /// Compute the magnetic field (B) at the given points.
18    ///
19    /// # Arguments
20    /// - `points`: Slice of observer positions (m)
21    ///
22    /// # Returns
23    /// - B-field vectors at each observer.
24    fn get_B(&self, points: &[Point3<T>]) -> Vec<Vector3<T>>;
25}
26
27/// Trait shared by magnetic sources.
28///
29/// Requires [`Transform`] and [`Field`].
30pub trait Source<T: RealField + Copy>:
31    Transform<T> + Field<T> + Debug + Send + Sync + Display
32{
33}
34
35macro_rules! impl_default {
36    () => {
37        fn default() -> Self {
38            Self {
39                position: Point3::origin(),
40                orientation: UnitQuaternion::identity(),
41                children: Vec::new(),
42            }
43        }
44    };
45}
46
47/// Implement deep Transform for objects with children.
48macro_rules! impl_transform_collection {
49    () => {
50        #[inline]
51        fn position(&self) -> Point3<T> {
52            self.position
53        }
54
55        #[inline]
56        fn orientation(&self) -> UnitQuaternion<T> {
57            self.orientation
58        }
59
60        #[inline]
61        fn set_position(&mut self, position: Point3<T>) {
62            let translation = Translation3::from(position - &self.position);
63            self.children
64                .iter_mut()
65                .for_each(|source| source.translate(&translation));
66
67            self.position = position
68        }
69
70        #[inline]
71        fn set_orientation(&mut self, orientation: UnitQuaternion<T>) {
72            let rotation = orientation * &self.orientation.inverse();
73            self.children
74                .iter_mut()
75                .for_each(|source| source.rotate_anchor(&rotation, &self.position));
76
77            self.orientation = orientation;
78        }
79
80        #[inline]
81        fn translate(&mut self, translation: &Translation3<T>) {
82            self.children
83                .iter_mut()
84                .for_each(|source| source.translate(translation));
85
86            self.position = translation.transform_point(&self.position);
87        }
88
89        #[inline]
90        fn rotate(&mut self, rotation: &UnitQuaternion<T>) {
91            #[cfg(feature = "parallel")]
92            if self.children.len() > 5000 {
93                use rayon::iter::{IntoParallelRefMutIterator, ParallelIterator};
94
95                self.children
96                    .par_iter_mut()
97                    .for_each(|source| source.rotate_anchor(rotation, &self.position));
98
99                self.orientation = rotation * self.orientation;
100                return;
101            }
102
103            self.children
104                .iter_mut()
105                .for_each(|source| source.rotate_anchor(rotation, &self.position));
106            self.orientation = rotation * self.orientation;
107        }
108
109        #[inline]
110        fn rotate_anchor(&mut self, rotation: &UnitQuaternion<T>, anchor: &Point3<T>) {
111            #[cfg(feature = "parallel")]
112            if self.children.len() > 5000 {
113                use rayon::iter::{IntoParallelRefMutIterator, ParallelIterator};
114                self.children
115                    .par_iter_mut()
116                    .for_each(|source| source.rotate_anchor(rotation, anchor));
117            } else {
118                self.children
119                    .iter_mut()
120                    .for_each(|source| source.rotate_anchor(rotation, anchor));
121            }
122
123            #[cfg(not(feature = "parallel"))]
124            {
125                self.children
126                    .iter_mut()
127                    .for_each(|source| source.rotate_anchor(rotation, anchor));
128            }
129
130            let local_position = self.position - anchor;
131            self.position = Point3::from(rotation * local_position + Vector3::from(anchor.coords));
132            self.orientation = rotation * &self.orientation;
133        }
134    };
135}
136
137/// Implement Field for source collection-like structs.
138macro_rules! impl_field_collection {
139    () => {
140        #[inline]
141        fn get_B(&self, points: &[Point3<T>]) -> Vec<Vector3<T>> {
142            let mut net_field = vec![Vector3::zeros(); points.len()];
143            #[cfg(feature = "parallel")]
144            {
145                use rayon::prelude::*;
146                let b_fields: Vec<_> = self
147                    .children
148                    .par_iter()
149                    .map(|source| source.get_B(points))
150                    .collect();
151                b_fields.iter().for_each(|child_b_field| {
152                    net_field
153                        .iter_mut()
154                        .zip(child_b_field)
155                        .for_each(|(sum, b)| *sum += b)
156                });
157            }
158            #[cfg(not(feature = "parallel"))]
159            {
160                for source in &self.children {
161                    let b_fields = source.get_B(points);
162                    net_field
163                        .iter_mut()
164                        .zip(b_fields)
165                        .for_each(|(sum, b)| *sum += b);
166                }
167            }
168            net_field
169        }
170    };
171}
172
173/// Stack-allocated collection of a single source type.
174///
175/// # Fields
176/// - `position`: Center of the collection (m), where the children reference.
177/// - `orientation`: Orientation of the collection, where the children reference.
178/// - `children`: An ordered-vec of homogeneous magnetic sources.
179///
180/// # Example
181/// ```
182/// use magba::sources::*;
183/// use nalgebra::*;
184///
185/// let magnet = CylinderMagnet::<f64>::default();
186/// let mut collection = SourceCollection::default();
187/// collection.add(magnet);
188/// ```
189#[derive(Debug)]
190pub struct SourceCollection<S: Source<T>, T: Float> {
191    /// Center of the collection (m), where the children reference
192    position: Point3<T>,
193    /// Orientation of the collection, where the children reference
194    orientation: UnitQuaternion<T>,
195
196    /// An ordered-vec of homogeneous magnetic sources
197    children: Vec<S>,
198}
199
200impl<S: Source<T> + PartialEq, T: Float> PartialEq for SourceCollection<S, T> {
201    fn eq(&self, other: &Self) -> bool {
202        self.position == other.position
203            && self.orientation == other.orientation
204            && self.children.len() == other.children.len()
205            && self.children.iter().all(|source| {
206                other
207                    .children
208                    .iter()
209                    .any(|other_source| source.eq(other_source))
210            })
211    }
212}
213
214impl<S: Source<T>, T: Float> Source<T> for SourceCollection<S, T> {}
215
216impl<S: Source<T>, T: Float> SourceCollection<S, T> {
217    /// Initialize [SourceCollection].
218    pub fn new(position: Point3<T>, orientation: UnitQuaternion<T>, sources: Vec<S>) -> Self {
219        Self {
220            position,
221            orientation,
222            children: sources,
223        }
224    }
225
226    /// Initialize [SourceCollection] from a vec of homogeneous [Source].
227    pub fn from_sources(sources: Vec<S>) -> Self {
228        Self {
229            position: Point3::origin(),
230            orientation: UnitQuaternion::identity(),
231            children: sources,
232        }
233    }
234
235    /// Add [Source] to the collection.
236    pub fn add(&mut self, source: S) {
237        self.children.push(source);
238    }
239
240    /// Add multiple [Source] to the collection.
241    pub fn add_sources(&mut self, source: &mut Vec<S>) {
242        self.children.append(source);
243    }
244}
245
246impl<S: Source<T>, T: Float> Display for SourceCollection<S, T> {
247    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
248        let len = self.children.len();
249        writeln!(
250            f,
251            "SourceCollection ({} children) at pos={}, q={}",
252            len,
253            crate_util::format_point3!(self.position),
254            crate_util::format_quat!(self.orientation)
255        )?;
256        for (i, source) in self.children.iter().enumerate() {
257            if i + 1 != len {
258                writeln!(f, "├── {}: {}", i, source)?;
259            } else {
260                write!(f, "└── {}: {}", i, source)?;
261            }
262        }
263        Ok(())
264    }
265}
266
267impl<S: Source<T>, T: Float> Default for SourceCollection<S, T> {
268    impl_default!();
269}
270
271impl<S: Source<T>, T: Float> Transform<T> for SourceCollection<S, T> {
272    impl_transform_collection!();
273}
274
275impl<S: Source<T>, T: Float> Field<T> for SourceCollection<S, T> {
276    impl_field_collection!();
277}
278
279/// Heap-allocated collection of multiple source types.
280///
281/// # Fields
282/// - `position`: Center of the collection (m), where the children reference.
283/// - `orientation`: Orientation of the collection, where the children reference.
284/// - `children`: An ordered-vec of heterogeneous magnetic sources.
285///
286/// # Example
287/// ```
288/// use magba::sources::*;
289/// use nalgebra::*;
290///
291/// let mut collection = MultiSourceCollection::<f64>::default();
292/// collection.add(Box::new(CylinderMagnet::default()));
293/// collection.add(Box::new(CuboidMagnet::default()));
294/// ```
295#[derive(Debug)]
296pub struct MultiSourceCollection<T: Float> {
297    /// Center of the collection (m), where the children reference
298    position: Point3<T>,
299    /// Orientation of the collection, where the children reference
300    orientation: UnitQuaternion<T>,
301
302    /// An ordered-vec of heterogeneous magnetic sources
303    children: Vec<Box<dyn Source<T>>>,
304}
305
306impl<T: Float> Source<T> for MultiSourceCollection<T> {}
307
308impl<T: Float> MultiSourceCollection<T> {
309    /// Initialize [`MultiSourceCollection`].
310    pub fn new(
311        position: Point3<T>,
312        orientation: UnitQuaternion<T>,
313        sources: Vec<Box<dyn Source<T>>>,
314    ) -> Self {
315        Self {
316            position,
317            orientation,
318            children: sources,
319        }
320    }
321
322    /// Initialize [MultiSourceCollection] from a vec of [Source].
323    pub fn from_sources(sources: Vec<Box<dyn Source<T>>>) -> Self {
324        Self {
325            position: Point3::origin(),
326            orientation: UnitQuaternion::identity(),
327            children: sources,
328        }
329    }
330
331    /// Add [Source] to the collection.
332    pub fn add(&mut self, source: Box<dyn Source<T>>) {
333        self.children.push(source);
334    }
335
336    /// Add multiple [Source] to the collection.
337    pub fn add_sources(&mut self, source: &mut Vec<Box<dyn Source<T>>>) {
338        self.children.append(source);
339    }
340}
341
342impl<T: Float> Default for MultiSourceCollection<T> {
343    impl_default!();
344}
345
346impl<T: Float> Display for MultiSourceCollection<T> {
347    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
348        let len = self.children.len();
349        writeln!(
350            f,
351            "MultiSourceCollection ({} children) at pos={}, q={}",
352            len,
353            crate_util::format_point3!(self.position),
354            crate_util::format_quat!(self.orientation)
355        )?;
356        for (i, source) in self.children.iter().enumerate() {
357            if i + 1 != len {
358                writeln!(f, "├── {}: {}", i, source)?;
359            } else {
360                write!(f, "└── {}: {}", i, source)?;
361            }
362        }
363        Ok(())
364    }
365}
366
367impl<T: Float> Transform<T> for MultiSourceCollection<T> {
368    impl_transform_collection!();
369}
370
371impl<T: Float> Field<T> for MultiSourceCollection<T> {
372    impl_field_collection!();
373}
374
375#[cfg(test)]
376mod base_source_collection_tests {
377    use std::f64::consts::FRAC_PI_2;
378
379    use super::*;
380    use crate::{sources::*, testing_util::quat_from_rotvec};
381
382    #[test]
383    fn test_collection_display() {
384        let magnet1 = CylinderMagnet::new(
385            Point3::new(4.0, 5.0, 6.0),
386            UnitQuaternion::identity(),
387            Vector3::new(1.0, 2.0, 3.0),
388            0.1,
389            0.3,
390        );
391        let magnet2 = CylinderMagnet::new(
392            Point3::new(10.0, 11.0, 12.0),
393            quat_from_rotvec(FRAC_PI_2, 0.0, 0.0),
394            Vector3::new(7.0, 8.0, 9.0),
395            0.1,
396            0.3,
397        );
398
399        let collection = SourceCollection::new(
400            Point3::origin(),
401            UnitQuaternion::identity(),
402            vec![magnet1, magnet2],
403        );
404
405        assert_eq!("SourceCollection (2 children) at pos=[0, 0, 0], q=[0, 0, 0, 1]
406├── 0: CylinderMagnet (pol=[1, 2, 3], r=0.1, h=0.3) at pos=[4, 5, 6], q=[0, 0, 0, 1]
407└── 1: CylinderMagnet (pol=[7, 8, 9], r=0.1, h=0.3) at pos=[10, 11, 12], q=[0.7071067811865475, 0, 0, 0.7071067811865476]", format!("{}", collection))
408    }
409}
410
411#[cfg(test)]
412mod cylinder_collection_tests {
413    use std::f64::consts::PI;
414
415    use super::*;
416    use crate::{sources::*, testing_util::*};
417
418    fn get_collection() -> SourceCollection<CylinderMagnet<f64>, f64> {
419        let mut collection = SourceCollection::default();
420        collection.add(CylinderMagnet::new(
421            Point3::new(0.009389999999999999, 0.0, -0.006),
422            quat_from_rotvec(1.2091995761561452, 1.209199576156145, 1.2091995761561452),
423            Vector3::new(1.0, 2.0, 3.0),
424            1.5e-3,
425            4e-3,
426        ));
427        collection.add(CylinderMagnet::new(
428            Point3::new(-0.004694999999999998, 0.008131978541535878, -0.006),
429            quat_from_rotvec(1.5315599088338596, 0.41038024073191587, 0.4103802407319159),
430            Vector3::new(0.4, 0.5, 0.6),
431            2e-3,
432            5e-3,
433        ));
434        collection.add(CylinderMagnet::new(
435            Point3::new(-0.004695000000000004, -0.008131978541535875, -0.006),
436            quat_from_rotvec(1.5315599088338594, -0.410380240731917, -0.41038024073191703),
437            Vector3::new(0.9, 0.8, 0.6),
438            2.5e-3,
439            6e-3,
440        ));
441        collection
442    }
443
444    #[test]
445    fn test_collection() {
446        let collection = get_collection();
447        test_B_magnet!(@small, &collection, "cylinder-collection.csv", 5e-9);
448    }
449
450    #[test]
451    fn test_collection_translate() {
452        let mut collection = get_collection();
453        let translation = Translation3::new(0.01, 0.015, 0.02);
454        collection.translate(&translation);
455        test_B_magnet!(@small, &collection, "cylinder-collection-translate.csv", 1e-8);
456
457        collection.translate(&translation.inverse());
458        collection.set_position(Point3::new(0.01, 0.015, 0.02));
459        test_B_magnet!(@small, &collection, "cylinder-collection-translate.csv", 1e-8);
460    }
461
462    #[test]
463    fn test_collection_rotate() {
464        let mut collection = get_collection();
465        let rotation = quat_from_rotvec(PI / 3.0, PI / 4.0, PI / 5.0);
466        collection.rotate(&rotation);
467        test_B_magnet!(@small, &collection, "cylinder-collection-rotate.csv", 1e-9);
468
469        collection.rotate(&rotation.inverse());
470        collection.set_orientation(rotation);
471        test_B_magnet!(@small, &collection, "cylinder-collection-rotate.csv", 5e-8);
472
473        collection.set_position(Point3::new(0.01, 0.015, 0.02));
474        test_B_magnet!(@small, &collection, "cylinder-collection-translate-rotate.csv", 5e-8);
475    }
476}
477
478#[cfg(test)]
479mod cuboid_collection_tests {
480    use std::f64::consts::{FRAC_PI_3, PI};
481
482    use super::*;
483    use crate::{sources::*, testing_util::*};
484
485    fn get_collection() -> SourceCollection<CuboidMagnet<f64>, f64> {
486        let mut collection = SourceCollection::default();
487        collection.add(CuboidMagnet::new(
488            Point3::new(0.005, 0.01, 0.015),
489            UnitQuaternion::identity(),
490            Vector3::new(0.1, 0.2, 0.3),
491            Vector3::new(0.02, 0.02, 0.03),
492        ));
493        collection.add(CuboidMagnet::new(
494            Point3::new(0.015, 0.005, 0.01),
495            quat_from_rotvec(0.0, FRAC_PI_3, 0.0),
496            Vector3::new(0.1, 0.2, 0.3),
497            Vector3::new(0.02, 0.02, 0.03),
498        ));
499        collection.add(CuboidMagnet::new(
500            Point3::new(0.01, 0.015, 0.005),
501            quat_from_rotvec(0.0, 0.0, FRAC_PI_3),
502            Vector3::new(0.1, 0.2, 0.3),
503            Vector3::new(0.02, 0.02, 0.03),
504        ));
505        collection
506    }
507
508    #[test]
509    fn test_collection() {
510        let collection = get_collection();
511        test_B_magnet!(@small, &collection, "cuboid-collection.csv", 2e-13);
512    }
513
514    #[test]
515    fn test_collection_translate() {
516        let mut collection = get_collection();
517        let translation = Translation3::new(0.01, 0.015, 0.02);
518        collection.translate(&translation);
519        test_B_magnet!(@small, &collection, "cuboid-collection-translate.csv", 2e-13);
520
521        collection.translate(&translation.inverse());
522        collection.set_position(Point3::new(0.01, 0.015, 0.02));
523        test_B_magnet!(@small, &collection, "cuboid-collection-translate.csv", 2e-13);
524    }
525
526    #[test]
527    fn test_collection_rotate() {
528        let mut collection = get_collection();
529        let rotation = quat_from_rotvec(PI / 3.0, PI / 4.0, PI / 5.0);
530        collection.rotate(&rotation);
531        test_B_magnet!(@small, &collection, "cuboid-collection-rotate.csv", 2e-13);
532
533        collection.rotate(&rotation.inverse());
534        collection.set_orientation(rotation);
535        test_B_magnet!(@small, &collection, "cuboid-collection-rotate.csv", 2e-13);
536
537        collection.set_position(Point3::new(0.01, 0.015, 0.02));
538        test_B_magnet!(@small, &collection, "cuboid-collection-translate-rotate.csv", 2e-13);
539    }
540}
541
542#[cfg(test)]
543mod multi_source_collection_tests {
544    use std::f64::consts::{FRAC_PI_3, PI};
545
546    use super::*;
547    use crate::{sources::*, testing_util::*};
548
549    fn get_collection() -> MultiSourceCollection<f64> {
550        let mut collection = MultiSourceCollection::default();
551        collection.add(Box::new(CylinderMagnet::new(
552            Point3::new(0.005, 0.01, 0.015),
553            UnitQuaternion::identity(),
554            Vector3::new(0.1, 0.2, 0.3),
555            0.02,
556            0.05,
557        )));
558        collection.add(Box::new(CuboidMagnet::new(
559            Point3::new(0.015, 0.005, 0.01),
560            quat_from_rotvec(0.0, FRAC_PI_3, 0.0),
561            Vector3::new(0.1, 0.2, 0.3),
562            Vector3::new(0.02, 0.02, 0.03),
563        )));
564        collection
565    }
566
567    #[test]
568    fn test_collection() {
569        let collection = get_collection();
570        test_B_magnet!(@small, &collection, "multi-collection.csv", 5e-11);
571    }
572
573    #[test]
574    fn test_collection_translate() {
575        let mut collection = get_collection();
576        let translation = Translation3::new(0.01, 0.015, 0.02);
577        collection.translate(&translation);
578        test_B_magnet!(@small, &collection, "multi-collection-translate.csv", 5e-11);
579
580        collection.translate(&translation.inverse());
581        collection.set_position(Point3::new(0.01, 0.015, 0.02));
582        test_B_magnet!(@small, &collection, "multi-collection-translate.csv", 5e-11);
583    }
584
585    #[test]
586    fn test_collection_rotate() {
587        let mut collection = get_collection();
588        let rotation = quat_from_rotvec(PI / 3.0, PI / 4.0, PI / 5.0);
589        collection.rotate(&rotation);
590        test_B_magnet!(@small, &collection, "multi-collection-rotate.csv", 5e-11);
591
592        collection.rotate(&rotation.inverse());
593        collection.set_orientation(rotation);
594        test_B_magnet!(@small, &collection, "multi-collection-rotate.csv", 5e-11);
595
596        collection.set_position(Point3::new(0.01, 0.015, 0.02));
597        test_B_magnet!(@small, &collection, "multi-collection-translate-rotate.csv", 5e-11);
598    }
599
600    #[test]
601    fn test_collection_display() {
602        let collection = get_collection();
603
604        println!("{}", collection);
605        assert_eq!("MultiSourceCollection (2 children) at pos=[0, 0, 0], q=[0, 0, 0, 1]
606├── 0: CylinderMagnet (pol=[0.1, 0.2, 0.3], r=0.02, h=0.05) at pos=[0.005, 0.01, 0.015], q=[0, 0, 0, 1]
607└── 1: CuboidMagnet (pol=[0.1, 0.2, 0.3], dim=[0.02, 0.02, 0.03]) at pos=[0.015, 0.005, 0.01], q=[0, 0.5, 0, 0.8660254037844386]",
608         format!("{}", collection))
609    }
610}
magba/sources/source.rs

magba/sources/
source.rs
