hoomd_interaction/univariate/harmonic.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 [`Harmonic`]
5
6use serde::{Deserialize, Serialize};
7
8use super::{UnivariateEnergy, UnivariateForce};
9
10/// Quadratic potential well centered at a given separation distance.
11///
12/// ```math
13/// U = \frac{1}{2} k (r - r_0)^2
14/// ```
15///
16/// Compute the harmonic potential and force as a function of `r` with
17/// equilibrium spring length `r_0`.
18///
19/// # Examples
20///
21/// ```
22/// use approxim::{assert_abs_diff_eq, assert_relative_eq};
23/// use hoomd_interaction::univariate::{
24/// Harmonic, UnivariateEnergy, UnivariateForce,
25/// };
26///
27/// let k = 2.0;
28/// let r_0 = 0.0;
29///
30/// let harmonic = Harmonic { k, r_0 };
31/// assert_abs_diff_eq!(harmonic.energy(0.0), 0.0);
32/// assert_relative_eq!(harmonic.energy(1.0), 1.0);
33/// assert_abs_diff_eq!(harmonic.force(1.0), -2.0, epsilon = 1e-12);
34/// ```
35///
36/// The parameters are public fields and may be accessed directly:
37///
38/// ```
39/// use hoomd_interaction::univariate::Harmonic;
40///
41/// let mut harmonic = Harmonic { k: 1.0, r_0: 0.0 };
42/// harmonic.k = 5.0;
43/// harmonic.r_0 = 1.0;
44/// ```
45#[derive(Clone, Debug, PartialEq, Serialize, Deserialize)]
46pub struct Harmonic {
47 /// Spring constant $`[\mathrm{energy}] [\mathrm{length}]^{-2}`$.
48 pub k: f64,
49 /// Equilibrium spring length $`[\mathrm{length}]`$.
50 pub r_0: f64,
51}
52
53impl UnivariateEnergy for Harmonic {
54 #[inline]
55 fn energy(&self, r: f64) -> f64 {
56 0.5 * self.k * (r - self.r_0) * (r - self.r_0)
57 }
58}
59
60impl UnivariateForce for Harmonic {
61 #[inline]
62 fn force(&self, r: f64) -> f64 {
63 -self.k * (r - self.r_0)
64 }
65}
66
67#[cfg(test)]
68mod tests {
69 use super::*;
70 use approxim::assert_relative_eq;
71 use rstest::*;
72
73 #[rstest]
74 fn zero_energy_point(#[values(1.0, 2.0, 5.0, 10.0)] k: f64, #[values(0.0, 1.0, 2.0)] r_0: f64) {
75 let harmonic = Harmonic { k, r_0 };
76
77 assert_eq!(harmonic.k, k);
78 assert_eq!(harmonic.r_0, r_0);
79
80 assert_eq!(harmonic.energy(r_0), 0.0);
81 assert_eq!(harmonic.force(r_0), 0.0);
82 }
83
84 #[rstest]
85 fn general_case(#[values(1.0, 2.0, 5.0, 10.0)] k: f64, #[values(0.0, 1.0, 2.0)] r_0: f64) {
86 let r = 5.0;
87 let harmonic = Harmonic { k, r_0 };
88
89 assert_eq!(harmonic.k, k);
90 assert_eq!(harmonic.r_0, r_0);
91
92 let expected_energy = 0.5 * k * (r - r_0) * (r - r_0);
93 let expected_force = -k * (r - r_0);
94
95 assert_relative_eq!(harmonic.energy(r), expected_energy);
96 assert_relative_eq!(harmonic.force(r), expected_force);
97 }
98}