hoomd_interaction/univariate/
weeks_chandler_anderson.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 [`WeeksChandlerAnderson`]
5
6use serde::{Deserialize, Serialize};
7
8use super::{LennardJones, UnivariateEnergy, UnivariateForce};
9
10/// Potential with a steep repulsive core.
11///
12/// ```math
13/// U(r) = \begin{cases}
14/// 4 \varepsilon \left[ \left( \frac{\sigma}{r} \right)^{12} - \left( \frac{\sigma}{r} \right)^{6} \right] + \varepsilon & r \lt 2^{1/6} \sigma \\
15///
16/// 0 & r \ge 2^{1/6} \sigma
17/// \end{cases}
18/// ```
19///
20/// Compute the Weeks-Chandler-Anderson (WCA) potential and force as a function of `r`.
21///
22/// # Examples
23///
24/// Basic usage:
25///
26/// ```
27/// use approxim::{assert_abs_diff_eq, assert_relative_eq};
28/// use hoomd_interaction::univariate::{
29///     UnivariateEnergy, UnivariateForce, WeeksChandlerAnderson,
30/// };
31///
32/// let epsilon = 1.5;
33/// let sigma = 2.5;
34///
35/// let wca = WeeksChandlerAnderson { epsilon, sigma };
36/// assert_relative_eq!(wca.energy(sigma), epsilon);
37/// assert_abs_diff_eq!(wca.energy(2.0 * sigma), 0.0);
38/// assert_relative_eq!(wca.energy(2.0_f64.powf(1.0 / 6.0) * sigma), 0.0);
39/// assert_abs_diff_eq!(
40///     wca.force(2.0_f64.powf(1.0 / 6.0) * sigma),
41///     0.0,
42///     epsilon = 1e-12
43/// );
44/// ```
45///
46/// The parameters are public fields and may be accessed directly:
47///
48/// ```
49/// use hoomd_interaction::univariate::WeeksChandlerAnderson;
50///
51/// let mut wca = WeeksChandlerAnderson::default();
52/// wca.epsilon = 1.5;
53/// wca.sigma = 3.0;
54/// ```
55#[derive(Clone, Debug, PartialEq, Serialize, Deserialize)]
56pub struct WeeksChandlerAnderson {
57    /// Energy scale *(\[energy\])*.
58    pub epsilon: f64,
59    /// Interaction width *(\[length\])*.
60    pub sigma: f64,
61}
62
63impl Default for WeeksChandlerAnderson {
64    /// Construct a [`WeeksChandlerAnderson`] with default parameters (epsilon=1.0, sigma=1.0)
65    ///
66    /// # Example
67    ///
68    /// ```
69    /// use hoomd_interaction::univariate::WeeksChandlerAnderson;
70    ///
71    /// let wca = WeeksChandlerAnderson::default();
72    /// ```
73    #[inline]
74    fn default() -> Self {
75        WeeksChandlerAnderson {
76            epsilon: 1.0,
77            sigma: 1.0,
78        }
79    }
80}
81
82impl UnivariateEnergy for WeeksChandlerAnderson {
83    #[inline]
84    fn energy(&self, r: f64) -> f64 {
85        if r < 2.0_f64.powf(1.0 / 6.0) * self.sigma {
86            let lj = LennardJones::<12, 6> {
87                epsilon: self.epsilon,
88                sigma: self.sigma,
89            };
90            lj.energy(r) + self.epsilon
91        } else {
92            0.0
93        }
94    }
95}
96
97impl UnivariateForce for WeeksChandlerAnderson {
98    #[inline]
99    fn force(&self, r: f64) -> f64 {
100        if r < 2.0_f64.powf(1.0 / 6.0) * self.sigma {
101            let lj = LennardJones::<12, 6> {
102                epsilon: self.epsilon,
103                sigma: self.sigma,
104            };
105            lj.force(r)
106        } else {
107            0.0
108        }
109    }
110}
111
112#[cfg(test)]
113mod tests {
114    use super::*;
115    use approxim::{assert_abs_diff_eq, assert_relative_eq};
116    use rstest::*;
117
118    #[rstest]
119    fn special_points_12_6(
120        #[values(1.0, 2.0, 12.125, 0.25)] epsilon: f64,
121        #[values(1.0, 2.0, 0.5)] sigma: f64,
122    ) {
123        let wca = WeeksChandlerAnderson { epsilon, sigma };
124
125        assert_eq!(wca.epsilon, epsilon);
126        assert_eq!(wca.sigma, sigma);
127
128        // Zero crossing (shifted by epsilon)
129        assert_relative_eq!(wca.energy(sigma), epsilon);
130        assert_relative_eq!(wca.force(sigma), 24.0 * epsilon / sigma);
131
132        // Bottom of the well
133        assert_abs_diff_eq!(wca.energy(2.0_f64.powf(1.0 / 6.0) * sigma), 0.0);
134        assert_abs_diff_eq!(
135            wca.force(2.0_f64.powf(1.0 / 6.0) * sigma),
136            0.0,
137            epsilon = 1e-12
138        );
139
140        // r = 2 sigma
141        assert_eq!(wca.energy(2.0 * sigma), 0.0);
142        assert_eq!(wca.force(2.0 * sigma), 0.0);
143    }
144}
hoomd_interaction/univariate/weeks_chandler_anderson.rs

hoomd_interaction/univariate/
weeks_chandler_anderson.rs
