differential_equations/
utils.rs

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