Skip to main content

ravenrs_animation/
spring.rs

1/// A damped harmonic oscillator with unit mass, solved in closed form.
2///
3/// Numeric integration drifts with the frame rate and can blow up after a long
4/// frame. The analytic solution gives the same position at a given time however
5/// many frames were drawn on the way there.
6#[derive(Clone, Copy, Debug, PartialEq)]
7pub struct Spring {
8    /// Restoring force per unit of displacement.
9    pub stiffness: f32,
10    /// Braking force per unit of velocity.
11    pub damping: f32,
12}
13
14impl Spring {
15    /// A spring with this stiffness and damping.
16    pub const fn new(stiffness: f32, damping: f32) -> Self {
17        Self { stiffness, damping }
18    }
19
20    /// For things that change place or size (the `spring-snappy` token).
21    pub const fn snappy() -> Self {
22        Self::new(380.0, 30.0)
23    }
24
25    /// For panels and sheets (the `spring-gentle` token).
26    pub const fn gentle() -> Self {
27        Self::new(170.0, 26.0)
28    }
29
30    /// A spring described the way designers think about it: `response` is the
31    /// period of the undamped oscillation in seconds, and `damping_ratio` is
32    /// 1.0 for no overshoot, less for some bounce.
33    pub fn from_response(response: f32, damping_ratio: f32) -> Self {
34        let omega = std::f32::consts::TAU / response.max(1e-3);
35        Self::new(omega * omega, 2.0 * damping_ratio.max(0.0) * omega)
36    }
37
38    /// Displacement and velocity `t` seconds after starting from displacement
39    /// `x0` with velocity `v0`. The spring pulls toward zero.
40    pub fn sample(&self, x0: f32, v0: f32, t: f32) -> (f32, f32) {
41        let k = self.stiffness.max(1e-3);
42        let omega = k.sqrt();
43        let zeta = self.damping.max(0.0) / (2.0 * omega);
44        let t = t.max(0.0);
45        if (zeta - 1.0).abs() < 1e-3 {
46            let c = v0 + omega * x0;
47            let decay = (-omega * t).exp();
48            ((x0 + c * t) * decay, (c - omega * (x0 + c * t)) * decay)
49        } else if zeta < 1.0 {
50            let wd = omega * (1.0 - zeta * zeta).sqrt();
51            let a = zeta * omega;
52            let c = (v0 + a * x0) / wd;
53            let decay = (-a * t).exp();
54            let (sin, cos) = (wd * t).sin_cos();
55            let x = decay * (x0 * cos + c * sin);
56            let v = decay * ((c * wd - a * x0) * cos - (x0 * wd + a * c) * sin);
57            (x, v)
58        } else {
59            let root = (zeta * zeta - 1.0).sqrt();
60            let (r1, r2) = (-omega * (zeta - root), -omega * (zeta + root));
61            let a = (v0 - r2 * x0) / (r1 - r2);
62            let b = x0 - a;
63            let (e1, e2) = ((r1 * t).exp(), (r2 * t).exp());
64            (a * e1 + b * e2, a * r1 * e1 + b * r2 * e2)
65        }
66    }
67}