hoomd_md/modify/
thermalize_angular_momentum.rs

1// Copyright (c) 2024-2026 The Regents of the University of Michigan.
2// Part of hoomd-rs, released under the BSD 3-Clause License.
3
4//! Implement `ThermalizeAngularMomentum`
5
6use super::ThermalizeAngularMomentum;
7use hoomd_microstate::{
8    Body, Microstate, SiteKey, Tagged, Transform,
9    boundary::{GenerateGhosts, Wrap},
10    property::{AngularMomentum, DynamicOrientedPoint, MomentOfInertia, Position},
11};
12use hoomd_spatial::PointUpdate;
13use hoomd_vector::{Angle, Cartesian, Outer, Versor, Wedge};
14use rand_distr::{Distribution, Normal};
15
16impl<P, S, X, C> ThermalizeAngularMomentum<DynamicOrientedPoint<P, Angle>, S>
17    for Microstate<DynamicOrientedPoint<P, Angle>, S, X, C>
18where
19    P: Copy + Wedge + Outer,
20    DynamicOrientedPoint<P, Angle>: Clone + Transform<S>,
21    S: Position<Position = P> + Default,
22    X: PointUpdate<P, SiteKey>,
23    C: Wrap<DynamicOrientedPoint<P, Angle>> + Wrap<S> + GenerateGhosts<S>,
24{
25    #[inline]
26    fn thermalize_angular_momentum_with_filter<
27        F: Fn(&Tagged<Body<DynamicOrientedPoint<P, Angle>, S>>) -> bool,
28    >(
29        &mut self,
30        temperature: f64,
31        should_thermalize_body: F,
32    ) {
33        let mut rng = self.counter().make_rng();
34
35        for body_index in 0..self.bodies().len() {
36            let body = &self.bodies()[body_index];
37            if !should_thermalize_body(body) {
38                continue;
39            }
40
41            let mut body_properties = body.item.properties.clone();
42
43            let moment_of_inertia = body_properties.moment_of_inertia();
44            let sigma = (temperature * moment_of_inertia).sqrt();
45            let normal = Normal::new(0.0, sigma).expect("Normal distribution should be valid");
46
47            *body_properties.angular_momentum_mut() = normal.sample(&mut rng);
48
49            self.update_body_properties(body_index, body_properties)
50                .expect("Bodies and sites should remain in simulation boundary.");
51        }
52
53        self.increment_substep();
54    }
55}
56
57impl<P, S, X, C> ThermalizeAngularMomentum<DynamicOrientedPoint<P, Versor>, S>
58    for Microstate<DynamicOrientedPoint<P, Versor>, S, X, C>
59where
60    P: Copy + Wedge + Outer,
61    DynamicOrientedPoint<P, Versor>: Clone + Transform<S>,
62    S: Position<Position = P> + Default,
63    X: PointUpdate<P, SiteKey>,
64    C: Wrap<DynamicOrientedPoint<P, Versor>> + Wrap<S> + GenerateGhosts<S>,
65{
66    #[inline]
67    fn thermalize_angular_momentum_with_filter<
68        F: Fn(&Tagged<Body<DynamicOrientedPoint<P, Versor>, S>>) -> bool,
69    >(
70        &mut self,
71        temperature: f64,
72        should_thermalize_body: F,
73    ) {
74        let mut rng = self.counter().make_rng();
75
76        for body_index in 0..self.bodies().len() {
77            let body = &self.bodies()[body_index];
78            if !should_thermalize_body(body) {
79                continue;
80            }
81
82            let mut body_properties = body.item.properties.clone();
83
84            let moment_of_inertia = body_properties.moment_of_inertia();
85
86            let x_nonzero = moment_of_inertia[0] > 0.0;
87            let y_nonzero = moment_of_inertia[1] > 0.0;
88            let z_nonzero = moment_of_inertia[2] > 0.0;
89            let sigma_x = (temperature * moment_of_inertia[0]).sqrt();
90            let sigma_y = (temperature * moment_of_inertia[1]).sqrt();
91            let sigma_z = (temperature * moment_of_inertia[2]).sqrt();
92            let normal_x = Normal::new(0.0, sigma_x).expect("Normal distribution should be valid.");
93            let normal_y = Normal::new(0.0, sigma_y).expect("Normal distribution should be valid.");
94            let normal_z = Normal::new(0.0, sigma_z).expect("Normal distribution should be valid.");
95
96            let mut angular_momentum_new = Cartesian::<3>::default();
97
98            if x_nonzero {
99                angular_momentum_new[0] = normal_x.sample(&mut rng);
100            }
101            if y_nonzero {
102                angular_momentum_new[1] = normal_y.sample(&mut rng);
103            }
104            if z_nonzero {
105                angular_momentum_new[2] = normal_z.sample(&mut rng);
106            }
107
108            *body_properties.angular_momentum_mut() = angular_momentum_new;
109            self.update_body_properties(body_index, body_properties)
110                .expect("Bodies and sites should remain in simulation boundary.");
111        }
112
113        self.increment_substep();
114    }
115}
hoomd_md/modify/thermalize_angular_momentum.rs

hoomd_md/modify/
thermalize_angular_momentum.rs
