atomecs/laser_cooling/
doppler.rs

1//! Calculations of the Doppler shift.
2extern crate rayon;
3extern crate serde;
4use specs::prelude::*;
5
6use super::CoolingLight;
7use crate::atom::Velocity;
8use crate::laser::gaussian::GaussianBeam;
9use crate::laser::index::LaserIndex;
10use serde::Serialize;
11use specs::{Component, Join, ReadStorage, System, VecStorage, WriteStorage};
12
13const LASER_CACHE_SIZE: usize = 16;
14
15/// Represents the Dopplershift of the atom with respect to each beam due to the atom velocity
16#[derive(Clone, Copy, Serialize)]
17pub struct DopplerShiftSampler {
18    /// detuning value in rad/s
19    pub doppler_shift: f64,
20}
21
22impl Default for DopplerShiftSampler {
23    fn default() -> Self {
24        DopplerShiftSampler {
25            /// Doppler shift with respect to laser beam, in SI units of rad/s.
26            doppler_shift: f64::NAN,
27        }
28    }
29}
30
31/// This system calculates the Doppler shift for each atom in each cooling beam.
32///
33/// The result is stored in `DopplerShiftSamplers`
34pub struct CalculateDopplerShiftSystem<const N: usize>;
35
36impl<'a, const N: usize> System<'a> for CalculateDopplerShiftSystem<N> {
37    type SystemData = (
38        ReadStorage<'a, CoolingLight>,
39        ReadStorage<'a, LaserIndex>,
40        ReadStorage<'a, GaussianBeam>,
41        WriteStorage<'a, DopplerShiftSamplers<N>>,
42        ReadStorage<'a, Velocity>,
43    );
44
45    fn run(&mut self, (cooling, indices, gaussian, mut samplers, velocities): Self::SystemData) {
46        use rayon::prelude::*;
47
48        // There are typically only a small number of lasers in a simulation.
49        // For a speedup, cache the required components into thread memory,
50        // so they can be distributed to parallel workers during the atom loop.
51        type CachedLaser = (CoolingLight, LaserIndex, GaussianBeam);
52        let laser_cache: Vec<CachedLaser> = (&cooling, &indices, &gaussian)
53            .join()
54            .map(|(cooling, index, gaussian)| (*cooling, *index, *gaussian))
55            .collect();
56
57        // Perform the iteration over atoms, `LASER_CACHE_SIZE` at a time.
58        for base_index in (0..laser_cache.len()).step_by(LASER_CACHE_SIZE) {
59            let max_index = laser_cache.len().min(base_index + LASER_CACHE_SIZE);
60            let slice = &laser_cache[base_index..max_index];
61            let mut laser_array = vec![laser_cache[0]; LASER_CACHE_SIZE];
62            laser_array[..max_index].copy_from_slice(slice);
63            let number_in_iteration = slice.len();
64
65            (&mut samplers, &velocities)
66                .par_join()
67                .for_each(|(sampler, vel)| {
68                    for (cooling, index, gaussian) in laser_array.iter().take(number_in_iteration) {
69                        sampler.contents[index.index].doppler_shift = vel
70                            .vel
71                            .dot(&(gaussian.direction.normalize() * cooling.wavenumber()));
72                    }
73                })
74        }
75    }
76}
77
78/// Component that holds a list of `DopplerShiftSampler`s
79///
80/// Each list entry corresponds to the detuning with respect to a CoolingLight entity
81/// and is indext via `CoolingLightIndex`
82#[derive(Clone, Copy, Serialize)]
83pub struct DopplerShiftSamplers<const N: usize> {
84    /// List of all `DopplerShiftSampler`s
85    #[serde(with = "serde_arrays")]
86    pub contents: [DopplerShiftSampler; N],
87}
88impl<const N: usize> Component for DopplerShiftSamplers<N> {
89    type Storage = VecStorage<Self>;
90}
91
92/// This system initialises all `DopplerShiftSamplers` to a NAN value.
93///
94/// It also ensures that the size of the `DopplerShiftSamplers` components match the number of CoolingLight entities in the world.
95pub struct InitialiseDopplerShiftSamplersSystem<const N: usize>;
96
97impl<'a, const N: usize> System<'a> for InitialiseDopplerShiftSamplersSystem<N> {
98    type SystemData = (WriteStorage<'a, DopplerShiftSamplers<N>>,);
99    fn run(&mut self, (mut samplers,): Self::SystemData) {
100        use rayon::prelude::*;
101
102        (&mut samplers).par_join().for_each(|mut sampler| {
103            sampler.contents = [DopplerShiftSampler::default(); N];
104        });
105    }
106}
107
108#[cfg(test)]
109pub mod tests {
110
111    use super::*;
112    use crate::constant::PI;
113    use crate::laser_cooling::CoolingLight;
114    use assert_approx_eq::assert_approx_eq;
115    extern crate nalgebra;
116    use crate::laser::{gaussian, DEFAULT_BEAM_LIMIT};
117    use nalgebra::Vector3;
118
119    #[test]
120    fn test_calculate_doppler_shift_system() {
121        let mut test_world = World::new();
122        test_world.register::<LaserIndex>();
123        test_world.register::<CoolingLight>();
124        test_world.register::<GaussianBeam>();
125        test_world.register::<Velocity>();
126        test_world.register::<DopplerShiftSamplers<{ DEFAULT_BEAM_LIMIT }>>();
127
128        let wavelength = 780e-9;
129        test_world
130            .create_entity()
131            .with(CoolingLight {
132                polarization: 1,
133                wavelength,
134            })
135            .with(LaserIndex {
136                index: 0,
137                initiated: true,
138            })
139            .with(GaussianBeam {
140                direction: Vector3::new(1.0, 0.0, 0.0),
141                intersection: Vector3::new(0.0, 0.0, 0.0),
142                e_radius: 2.0,
143                power: 1.0,
144                rayleigh_range: gaussian::calculate_rayleigh_range(&wavelength, &2.0),
145                ellipticity: 0.0,
146            })
147            .build();
148
149        let atom_velocity = 100.0;
150        let sampler1 = test_world
151            .create_entity()
152            .with(Velocity {
153                vel: Vector3::new(atom_velocity, 0.0, 0.0),
154            })
155            .with(DopplerShiftSamplers {
156                contents: [DopplerShiftSampler::default(); crate::laser::DEFAULT_BEAM_LIMIT],
157            })
158            .build();
159
160        let mut system = CalculateDopplerShiftSystem::<{ DEFAULT_BEAM_LIMIT }>;
161        system.run_now(&test_world);
162        test_world.maintain();
163        let sampler_storage =
164            test_world.read_storage::<DopplerShiftSamplers<{ DEFAULT_BEAM_LIMIT }>>();
165
166        assert_approx_eq!(
167            sampler_storage
168                .get(sampler1)
169                .expect("entity not found")
170                .contents[0]
171                .doppler_shift,
172            2.0 * PI / wavelength * atom_velocity,
173            1e-5_f64
174        );
175    }
176}