1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
|
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;
}
if (self.to - self.from).abs() <= f64::EPSILON {
return Duration::ZERO;
}
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);
// Overdamped springs have some numerical stability issues...
if !y1.is_finite() {
return Duration::from_secs_f64(x0);
}
i += 1;
}
Duration::from_secs_f64(x1)
}
/// Computes and returns the duration until the spring reaches its target position.
pub fn clamped_duration(&self) -> Option<Duration> {
let beta = self.params.damping / (2. * self.params.mass);
if beta.abs() <= f64::EPSILON || beta < 0. {
return Some(Duration::MAX);
}
if (self.to - self.from).abs() <= f64::EPSILON {
return Some(Duration::ZERO);
}
// The first frame is not that important and we avoid finding the trivial 0 for in-place
// animations.
let mut i = 1u16;
let mut y = self.oscillate(f64::from(i) / 1000.);
while (self.to - self.from > f64::EPSILON && self.to - y > self.params.epsilon)
|| (self.from - self.to > f64::EPSILON && y - self.to > self.params.epsilon)
{
if i > 3000 {
return None;
}
i += 1;
y = self.oscillate(f64::from(i) / 1000.);
}
Some(Duration::from_millis(u64::from(i)))
}
/// 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())
}
}
}
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn overdamped_spring_equal_from_to_nan() {
let spring = Spring {
from: 0.,
to: 0.,
initial_velocity: 0.,
params: SpringParams::new(1.15, 850., 0.0001),
};
let _ = spring.duration();
let _ = spring.clamped_duration();
let _ = spring.value_at(Duration::ZERO);
}
#[test]
fn overdamped_spring_duration_panic() {
let spring = Spring {
from: 0.,
to: 1.,
initial_velocity: 0.,
params: SpringParams::new(6., 1200., 0.0001),
};
let _ = spring.duration();
let _ = spring.clamped_duration();
let _ = spring.value_at(Duration::ZERO);
}
}
|
