138 lines
4.4 KiB
Rust
138 lines
4.4 KiB
Rust
// From Niri (GPL-3.0) https://github.com/YaLTeR/niri/tree/db49deb7fd2fbe805ceec060aa4dec65009ad7a7
|
|
use std::time::Duration;
|
|
|
|
#[derive(Debug, Clone, Copy)]
|
|
pub struct SpringParams {
|
|
pub damping: f64,
|
|
pub mass: f64,
|
|
pub stiffness: f64,
|
|
pub epsilon: f64,
|
|
}
|
|
|
|
#[derive(Debug, Clone, Copy)]
|
|
pub struct Spring {
|
|
pub from: f64,
|
|
pub to: f64,
|
|
pub initial_velocity: f64,
|
|
pub params: SpringParams,
|
|
}
|
|
|
|
impl SpringParams {
|
|
pub fn new(damping_ratio: f64, stiffness: f64, epsilon: f64) -> Self {
|
|
let damping_ratio = damping_ratio.max(0.);
|
|
let stiffness = stiffness.max(0.);
|
|
let epsilon = epsilon.max(0.);
|
|
|
|
let mass = 1.;
|
|
let critical_damping = 2. * (mass * stiffness).sqrt();
|
|
let damping = damping_ratio * critical_damping;
|
|
|
|
Self {
|
|
damping,
|
|
mass,
|
|
stiffness,
|
|
epsilon,
|
|
}
|
|
}
|
|
}
|
|
|
|
impl Spring {
|
|
pub fn value_at(&self, t: Duration) -> f64 {
|
|
self.oscillate(t.as_secs_f64())
|
|
}
|
|
|
|
// Based on libadwaita (LGPL-2.1-or-later):
|
|
// https://gitlab.gnome.org/GNOME/libadwaita/-/blob/1.4.4/src/adw-spring-animation.c,
|
|
// which itself is based on (MIT):
|
|
// https://github.com/robb/RBBAnimation/blob/master/RBBAnimation/RBBSpringAnimation.m
|
|
/// Computes and returns the duration until the spring is at rest.
|
|
pub fn duration(&self) -> Duration {
|
|
const DELTA: f64 = 0.001;
|
|
|
|
let beta = self.params.damping / (2. * self.params.mass);
|
|
|
|
if beta.abs() <= f64::EPSILON || beta < 0. {
|
|
return Duration::MAX;
|
|
}
|
|
|
|
let omega0 = (self.params.stiffness / self.params.mass).sqrt();
|
|
|
|
// As first ansatz for the overdamped solution,
|
|
// and general estimation for the oscillating ones
|
|
// we take the value of the envelope when it's < epsilon.
|
|
let mut x0 = -self.params.epsilon.ln() / beta;
|
|
|
|
// f64::EPSILON is too small for this specific comparison, so we use
|
|
// f32::EPSILON even though it's doubles.
|
|
if (beta - omega0).abs() <= f64::from(f32::EPSILON) || beta < omega0 {
|
|
return Duration::from_secs_f64(x0);
|
|
}
|
|
|
|
// Since the overdamped solution decays way slower than the envelope
|
|
// we need to use the value of the oscillation itself.
|
|
// Newton's root finding method is a good candidate in this particular case:
|
|
// https://en.wikipedia.org/wiki/Newton%27s_method
|
|
let mut y0 = self.oscillate(x0);
|
|
let m = (self.oscillate(x0 + DELTA) - y0) / DELTA;
|
|
|
|
let mut x1 = (self.to - y0 + m * x0) / m;
|
|
let mut y1 = self.oscillate(x1);
|
|
|
|
let mut i = 0;
|
|
while (self.to - y1).abs() > self.params.epsilon {
|
|
if i > 1000 {
|
|
return Duration::ZERO;
|
|
}
|
|
|
|
x0 = x1;
|
|
y0 = y1;
|
|
|
|
let m = (self.oscillate(x0 + DELTA) - y0) / DELTA;
|
|
|
|
x1 = (self.to - y0 + m * x0) / m;
|
|
y1 = self.oscillate(x1);
|
|
i += 1;
|
|
}
|
|
|
|
Duration::from_secs_f64(x1)
|
|
}
|
|
|
|
/// Returns the spring position at a given time in seconds.
|
|
fn oscillate(&self, t: f64) -> f64 {
|
|
let b = self.params.damping;
|
|
let m = self.params.mass;
|
|
let k = self.params.stiffness;
|
|
let v0 = self.initial_velocity;
|
|
|
|
let beta = b / (2. * m);
|
|
let omega0 = (k / m).sqrt();
|
|
|
|
let x0 = self.from - self.to;
|
|
|
|
let envelope = (-beta * t).exp();
|
|
|
|
// Solutions of the form C1*e^(lambda1*x) + C2*e^(lambda2*x)
|
|
// for the differential equation m*ẍ+b*ẋ+kx = 0
|
|
|
|
// f64::EPSILON is too small for this specific comparison, so we use
|
|
// f32::EPSILON even though it's doubles.
|
|
if (beta - omega0).abs() <= f64::from(f32::EPSILON) {
|
|
// Critically damped.
|
|
self.to + envelope * (x0 + (beta * x0 + v0) * t)
|
|
} else if beta < omega0 {
|
|
// Underdamped.
|
|
let omega1 = ((omega0 * omega0) - (beta * beta)).sqrt();
|
|
|
|
self.to
|
|
+ envelope
|
|
* (x0 * (omega1 * t).cos() + ((beta * x0 + v0) / omega1) * (omega1 * t).sin())
|
|
} else {
|
|
// Overdamped.
|
|
let omega2 = ((beta * beta) - (omega0 * omega0)).sqrt();
|
|
|
|
self.to
|
|
+ envelope
|
|
* (x0 * (omega2 * t).cosh() + ((beta * x0 + v0) / omega2) * (omega2 * t).sinh())
|
|
}
|
|
}
|
|
}
|