use std::ops::Index;
use crate::errors::QlResult;
use crate::math::comparison::close_enough;
use crate::types::{Real, Size, Time};
use crate::{fail, require};
#[derive(Clone, Debug, PartialEq, Default)]
pub struct TimeGrid {
times: Vec<Time>,
dt: Vec<Time>,
mandatory_times: Vec<Time>,
}
impl TimeGrid {
#[allow(clippy::neg_cmp_op_on_partial_ord)]
pub fn new(end: Time, steps: Size) -> QlResult<Self> {
require!(end > 0.0, "negative times not allowed");
require!(steps > 0, "at least one step required");
let dt = end / steps as Real;
let times = (0..=steps).map(|i| dt * i as Real).collect();
Ok(TimeGrid {
times,
dt: vec![dt; steps],
mandatory_times: vec![end],
})
}
#[allow(clippy::float_cmp, clippy::neg_cmp_op_on_partial_ord)]
pub fn with_mandatory_times(times: &[Time], steps: Size) -> QlResult<Self> {
require!(!times.is_empty(), "empty time sequence");
let mut mandatory = times.to_vec();
mandatory.sort_by(|a, b| {
a.partial_cmp(b)
.expect("time-grid mandatory times must be totally ordered")
});
require!(mandatory[0] >= 0.0, "negative times not allowed");
mandatory.dedup_by(|a, b| close_enough(*a, *b));
let last = *mandatory
.last()
.expect("mandatory is non-empty after the guard");
let dt_max = if steps == 0 {
let mut diff = Vec::with_capacity(mandatory.len());
diff.push(mandatory[0]);
for pair in mandatory.windows(2) {
diff.push(pair[1] - pair[0]);
}
require!(
!diff.is_empty(),
"at least two distinct points required in time grid"
);
if diff.first() == Some(&0.0) {
diff.remove(0);
}
require!(!diff.is_empty(), "not enough distinct points in time grid");
diff.into_iter().fold(Real::INFINITY, Real::min)
} else {
last / steps as Real
};
let mut grid_times = vec![0.0];
let mut period_begin = 0.0;
for &period_end in &mandatory {
if period_end != 0.0 {
let n_steps = (((period_end - period_begin) / dt_max).round()).max(1.0) as Size;
let dt = (period_end - period_begin) / n_steps as Real;
for n in 1..=n_steps {
grid_times.push(period_begin + n as Real * dt);
}
}
period_begin = period_end;
}
let dt = grid_times.windows(2).map(|w| w[1] - w[0]).collect();
Ok(TimeGrid {
times: grid_times,
dt,
mandatory_times: mandatory,
})
}
pub fn dt(&self, i: Size) -> Time {
self.dt[i]
}
pub fn mandatory_times(&self) -> &[Time] {
&self.mandatory_times
}
pub fn times(&self) -> &[Time] {
&self.times
}
pub fn size(&self) -> Size {
self.times.len()
}
pub fn empty(&self) -> bool {
self.times.is_empty()
}
pub fn at(&self, i: Size) -> Option<Time> {
self.times.get(i).copied()
}
pub fn closest_index(&self, t: Time) -> Size {
let result = self.times.partition_point(|&x| x < t);
if result == 0 {
0
} else if result == self.times.len() {
self.times.len() - 1
} else {
let dt1 = self.times[result] - t;
let dt2 = t - self.times[result - 1];
if dt1 < dt2 { result } else { result - 1 }
}
}
pub fn index(&self, t: Time) -> QlResult<Size> {
let i = self.closest_index(t);
if close_enough(t, self.times[i]) {
Ok(i)
} else {
fail!(
"using inadequate time grid: no node is close enough to the \
required time t = {t} (closest node is t1 = {})",
self.times[i]
);
}
}
pub fn front(&self) -> Option<Time> {
self.times.first().copied()
}
pub fn back(&self) -> Option<Time> {
self.times.last().copied()
}
}
impl Index<Size> for TimeGrid {
type Output = Time;
fn index(&self, i: Size) -> &Time {
&self.times[i]
}
}
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn regular_grid_matches_reference() {
let grid = TimeGrid::new(1.0, 4).unwrap();
assert_eq!(grid.size(), 5);
assert_eq!(grid.times(), &[0.0, 0.25, 0.5, 0.75, 1.0]);
for i in 0..4 {
assert_eq!(grid.dt(i), 0.25);
}
assert_eq!(grid.front(), Some(0.0));
assert_eq!(grid.back(), Some(1.0));
assert_eq!(grid.mandatory_times(), &[1.0]);
}
#[test]
fn indexing_and_bounds_checks() {
let grid = TimeGrid::new(1.0, 4).unwrap();
assert_eq!(grid[2], 0.5);
assert_eq!(grid.at(4), Some(1.0));
assert_eq!(grid.at(5), None);
assert!(!grid.empty());
}
#[test]
fn default_grid_is_empty() {
let grid = TimeGrid::default();
assert!(grid.empty());
assert_eq!(grid.size(), 0);
assert_eq!(grid.front(), None);
assert_eq!(grid.back(), None);
}
#[test]
fn non_positive_end_is_rejected() {
let err = TimeGrid::new(0.0, 4).unwrap_err();
assert_eq!(err.message(), "negative times not allowed");
assert!(TimeGrid::new(-1.0, 4).is_err());
}
#[test]
fn zero_steps_is_rejected() {
let err = TimeGrid::new(1.0, 0).unwrap_err();
assert_eq!(err.message(), "at least one step required");
}
#[test]
fn index_resolves_grid_aligned_times() {
let grid = TimeGrid::new(1.0, 4).unwrap();
assert_eq!(grid.index(0.0).unwrap(), 0);
assert_eq!(grid.index(0.25).unwrap(), 1);
assert_eq!(grid.index(0.5).unwrap(), 2);
assert_eq!(grid.index(1.0).unwrap(), 4);
}
#[test]
fn index_rejects_off_grid_times() {
let grid = TimeGrid::new(1.0, 4).unwrap();
assert!(grid.index(0.3).is_err());
assert!(grid.index(1.5).is_err());
assert!(grid.index(-0.1).is_err());
}
#[test]
fn with_mandatory_times_places_every_time_on_an_exact_node() {
let grid = TimeGrid::with_mandatory_times(&[1.0, 0.5, 1.0], 2).unwrap();
assert_eq!(grid.times(), &[0.0, 0.5, 1.0]);
assert_eq!(grid.mandatory_times(), &[0.5, 1.0]);
assert_eq!(grid[grid.index(0.5).unwrap()], 0.5);
assert_eq!(grid[grid.index(1.0).unwrap()], 1.0);
}
#[test]
fn with_mandatory_times_inserts_regular_interior_points() {
let grid = TimeGrid::with_mandatory_times(&[1.0, 3.0], 4).unwrap();
assert_eq!(grid.size(), 5);
assert_eq!(grid[1], 1.0, "first mandatory time is a verbatim node");
assert_eq!(
grid.back(),
Some(3.0),
"last mandatory time is the final node"
);
assert_eq!(grid.mandatory_times(), &[1.0, 3.0]);
assert!((grid.dt(1) - 2.0 / 3.0).abs() < 1e-15);
}
#[test]
fn with_mandatory_times_auto_spaces_when_steps_is_zero() {
let grid = TimeGrid::with_mandatory_times(&[1.0, 2.0, 4.0], 0).unwrap();
assert_eq!(grid.times(), &[0.0, 1.0, 2.0, 3.0, 4.0]);
assert_eq!(grid.mandatory_times(), &[1.0, 2.0, 4.0]);
}
#[test]
fn with_mandatory_times_rejects_empty_and_negative() {
assert_eq!(
TimeGrid::with_mandatory_times(&[], 4)
.unwrap_err()
.message(),
"empty time sequence"
);
assert_eq!(
TimeGrid::with_mandatory_times(&[-1.0, 2.0], 4)
.unwrap_err()
.message(),
"negative times not allowed"
);
}
#[test]
fn closest_index_snaps_to_nearest_node() {
let grid = TimeGrid::new(1.0, 4).unwrap();
assert_eq!(grid.closest_index(0.3), 1);
assert_eq!(grid.closest_index(0.4), 2);
assert_eq!(grid.closest_index(2.0), 4);
assert_eq!(grid.closest_index(-1.0), 0);
assert_eq!(grid.closest_index(0.125), 0);
}
}