differential_equations/
utils.rs

1//! Utility functions for the differential equation solvers
2
3use crate::{
4    error::Error,
5    traits::{Real, State},
6};
7
8/// Constrain the step size to be within the bounds of `h_min` and `h_max`.
9///
10/// If `h` is less than `h_min`, returns `h_min` with the same sign as `h`.
11/// If `h` is greater than `h_max`, returns `h_max` with the same sign as `h`.
12/// Otherwise, returns `h` unchanged.
13///
14/// # Arguments
15/// * `h` - Step size to constrain.
16/// * `h_min` - Minimum step size.
17/// * `h_max` - Maximum step size.
18///
19/// # Returns
20/// * The step size constrained to be within the bounds of `h_min` and `h_max`.
21///
22pub fn constrain_step_size<T: Real>(h: T, h_min: T, h_max: T) -> T {
23    // Determine the direction of the step size
24    let sign = h.signum();
25    // Bound the step size
26    if h.abs() < h_min {
27        sign * h_min
28    } else if h.abs() > h_max {
29        sign * h_max
30    } else {
31        h
32    }
33}
34
35/// Validate the step size parameters.
36///
37/// Checks the following:
38/// * `tf` cannot be equal to `t0`.
39/// * `h0` has the same sign as `tf - t0`.
40/// * `h_min` and `h_max` are non-negative.
41/// * `h_min` is less than or equal to `h_max`.
42/// * `|h0|` is greater than or equal to `h_min`.
43/// * `|h0|` is less than or equal to `h_max`.
44/// * `|h0|` is less than or equal to `|tf - t0|`.
45/// * `h0` is not zero.
46///
47/// If any of the checks fail, returns `Err(Error::BadInput)` with a descriptive message.
48/// Else returns `Ok(h0)` indicating the step size is valid.
49///
50/// # Arguments
51/// * `h0` - Initial step size.
52/// * `h_min` - Minimum step size.
53/// * `h_max` - Maximum step size.
54/// * `t0` - Initial time.
55/// * `tf` - Final time.
56///
57/// # Returns
58/// * `Result<Real, Error>` - Ok if all checks pass, Err if any check fails.
59///
60pub fn validate_step_size_parameters<T: Real, Y: State<T>>(
61    h0: T,
62    h_min: T,
63    h_max: T,
64    t0: T,
65    tf: T,
66) -> Result<T, Error<T, Y>> {
67    // Check if tf == t0
68    if tf == t0 {
69        return Err(Error::BadInput {
70            msg: format!(
71                "Invalid input: tf ({:?}) cannot be equal to t0 ({:?})",
72                tf, t0
73            ),
74        });
75    }
76
77    // Determine direction of the step size
78    let sign = (tf - t0).signum();
79
80    // Check h0 has same sign as tf - t0
81    if h0.signum() != sign {
82        return Err(Error::BadInput {
83            msg: format!(
84                "Invalid input: Initial step size ({:?}) must have the same sign as the integration direction (sign of tf - t0 = {:?})",
85                h0,
86                tf - t0
87            ),
88        });
89    }
90
91    // Check h_min and h_max bounds
92    if h_min < T::zero() {
93        return Err(Error::BadInput {
94            msg: format!(
95                "Invalid input: Minimum step size ({:?}) must be non-negative",
96                h_min
97            ),
98        });
99    }
100    if h_max < T::zero() {
101        return Err(Error::BadInput {
102            msg: format!(
103                "Invalid input: Maximum step size ({:?}) must be non-negative",
104                h_max
105            ),
106        });
107    }
108    if h_min > h_max {
109        return Err(Error::BadInput {
110            msg: format!(
111                "Invalid input: Minimum step size ({:?}) must be less than or equal to maximum step size ({:?})",
112                h_min, h_max
113            ),
114        });
115    }
116
117    // Check h0 bounds
118    if h0.abs() < h_min {
119        return Err(Error::BadInput {
120            msg: format!(
121                "Invalid input: Absolute value of initial step size ({:?}) must be greater than or equal to minimum step size ({:?})",
122                h0.abs(),
123                h_min
124            ),
125        });
126    }
127    if h0.abs() > h_max {
128        return Err(Error::BadInput {
129            msg: format!(
130                "Invalid input: Absolute value of initial step size ({:?}) must be less than or equal to maximum step size ({:?})",
131                h0.abs(),
132                h_max
133            ),
134        });
135    }
136
137    // Check h0 is not larger then integration interval
138    if h0.abs() > (tf - t0).abs() {
139        return Err(Error::BadInput {
140            msg: format!(
141                "Invalid input: Absolute value of initial step size ({:?}) must be less than or equal to the absolute value of the integration interval (tf - t0 = {:?})",
142                h0.abs(),
143                (tf - t0).abs()
144            ),
145        });
146    }
147
148    // Check h0 is not zero
149    if h0 == T::zero() {
150        return Err(Error::BadInput {
151            msg: format!("Invalid input: Initial step size ({:?}) cannot be zero", h0),
152        });
153    }
154
155    // Return Ok if all bounds are valid return the step size
156    Ok(h0)
157}