conspire/physics/molecular/potential/morse/
mod.rs

1use crate::{math::Scalar, physics::molecular::potential::Potential};
2
3/// The Morse potential.[^1]
4/// [^1]: P.M. Morse, [Physical Review **34**, 57 (1929)](https://doi.org/10.1103/PhysRev.34.57).
5#[derive(Clone, Debug)]
6pub struct Morse {
7    /// The rest length $`x_0`$.
8    pub rest_length: Scalar,
9    /// The potential depth $`u_0`$.
10    pub depth: Scalar,
11    /// The Morse parameter $`a`$.
12    pub parameter: Scalar,
13}
14
15impl Potential for Morse {
16    /// ```math
17    /// u(x) = u_0\left[1 - e^{-a(x - x_0)}\right]^2
18    /// ```
19    fn energy(&self, length: Scalar) -> Scalar {
20        self.depth * (1.0 - (self.parameter * (self.rest_length - length)).exp()).powi(2)
21    }
22    /// ```math
23    /// f(x) = 2au_0e^{-a(x - x_0)}\left[1 - e^{-a(x - x_0)}\right]
24    /// ```
25    fn force(&self, length: Scalar) -> Scalar {
26        let exp = (self.parameter * (self.rest_length - length)).exp();
27        2.0 * self.parameter * self.depth * exp * (1.0 - exp)
28    }
29    /// ```math
30    /// k(x) = 2a^2u_0e^{-a(x - x_0)}\left[2e^{-a(x - x_0)} - 1\right]
31    /// ```
32    fn stiffness(&self, length: Scalar) -> Scalar {
33        let exp = (self.parameter * (self.rest_length - length)).exp();
34        2.0 * self.parameter.powi(2) * self.depth * exp * (2.0 * exp - 1.0)
35    }
36    /// ```math
37    /// h(x) = 2a^3u_0e^{-a(x - x_0)}\left[1 - 4e^{-a(x - x_0)}\right]
38    /// ```
39    fn anharmonicity(&self, length: Scalar) -> Scalar {
40        let exp = (self.parameter * (self.rest_length - length)).exp();
41        2.0 * self.parameter.powi(3) * self.depth * exp * (1.0 - 4.0 * exp)
42    }
43    /// ```math
44    /// \Delta x(f) = \frac{1}{a}\,\ln\left(\frac{2}{1 + \sqrt{1 - f/f_\mathrm{max}}}\right)
45    /// ```
46    fn extension(&self, force: Scalar) -> Scalar {
47        let y = force / self.peak_force();
48        if (0.0..=1.0).contains(&y) {
49            (2.0 / (1.0 + (1.0 - y).sqrt())).ln() / self.parameter
50        } else {
51            Scalar::NAN
52        }
53    }
54    /// ```math
55    /// c(f) = \frac{1}{a^2u_0}\,\frac{\left(1-f/f_\mathrm{max}\right)^{-1/2}}{1+\sqrt{1-f/f_\mathrm{max}}}
56    /// ```
57    fn compliance(&self, force: Scalar) -> Scalar {
58        let y = force / self.peak_force();
59        if (0.0..1.0).contains(&y) {
60            let s = (1.0 - y).sqrt();
61            1.0 / (self.parameter.powi(2) * self.depth) / (s * (1.0 + s))
62        } else if y == 0.0 {
63            Scalar::INFINITY
64        } else {
65            Scalar::NAN
66        }
67    }
68    /// ```math
69    /// \text{arg max }u(x) = x_0 + \frac{1}{a}\,\ln(2)
70    /// ```
71    fn peak(&self) -> Scalar {
72        self.rest_length + 2.0_f64.ln() / self.parameter
73    }
74    /// ```math
75    /// f(x_\mathrm{peak}) = \frac{au_0}{2}
76    /// ```
77    fn peak_force(&self) -> Scalar {
78        0.5 * self.parameter * self.depth
79    }
80    /// ```math
81    /// \text{arg min }u(x) = x_0
82    /// ```
83    fn rest_length(&self) -> Scalar {
84        self.rest_length
85    }
86}
conspire/physics/molecular/potential/morse/mod.rs

conspire/physics/molecular/potential/morse/
mod.rs
