conspire/math/optimize/globalization/line_search/
mod.rs

1#[cfg(test)]
2mod test;
3
4mod armijo;
5mod error;
6mod goldstein;
7mod none;
8mod residual;
9mod wolfe;
10
11pub use armijo::Armijo;
12pub use error::Error;
13pub use goldstein::Goldstein;
14pub use none::NoSearch;
15pub use residual::Residual;
16pub use wolfe::Wolfe;
17
18use crate::math::{
19    Erase, Jacobian, Quantity, Scalar, Solution, Style, StyledError, Tensor, styled_error,
20};
21use std::{
22    fmt::{Debug, Display},
23    ops::Mul,
24};
25
26/// Line searches that judge a step by a merit built from the objective.
27///
28/// The objective is the function being minimized, plus a penalty where it is
29/// constrained, so these are the searches a minimization can take and a root
30/// cannot: a root has no objective to compare the values of.
31pub trait ObjectiveMerit: Clone + Debug + Display {
32    /// Whether anything is searched at all.
33    const SEARCHES: bool = true;
34    /// Whether the value and slope of the merit are read, and so worth preparing.
35    const MERIT: bool = true;
36    /// Backtrack along a step of the objective, whose gradient is known.
37    fn search<X, J, D, W, E>(
38        &self,
39        function: impl FnMut(&X, Scalar) -> Result<Scalar, String>,
40        jacobian: impl FnMut(&X) -> Result<J, String>,
41        argument: &X,
42        jacobian0: &J,
43        decrement: &D,
44        step_size: Scalar,
45    ) -> Result<Scalar, LineSearchError>
46    where
47        J: Erase<Erased = E> + Jacobian,
48        D: Erase<Erased = E>,
49        E: Tensor,
50        X: Solution,
51        for<'a> &'a D: Mul<Quantity<W>, Output = X>;
52    /// Backtrack on a merit function of the step size alone.
53    ///
54    /// The exact penalty function is not differentiable, its norm having a kink
55    /// wherever a constraint is satisfied, so its slope along the step is
56    /// supplied rather than recovered from a gradient.
57    fn search_merit(
58        &self,
59        merit: impl FnMut(Scalar) -> Result<Scalar, String>,
60        value: Scalar,
61        slope: Scalar,
62        step_size: Scalar,
63    ) -> Result<Scalar, LineSearchError>;
64}
65
66/// Line searches that judge a step by a merit built from the residual.
67///
68/// The merit is half the squared norm of the residual being driven to zero.
69/// Along an exact Newton step its slope is the negative of twice itself, so
70/// nothing beyond the residual need be known.
71pub trait ResidualMerit: Clone + Debug + Display {
72    /// Whether anything is searched at all.
73    const SEARCHES: bool = true;
74    /// Whether the value of the merit is read, and so worth preparing.
75    const MERIT: bool = true;
76    /// Backtrack on the residual merit of the step size alone.
77    ///
78    /// The step is judged along the share `fraction` of a step, which is
79    /// one unless something shortened the step beforehand.
80    fn search_residual(
81        &self,
82        merit: impl FnMut(Scalar) -> Result<Scalar, String>,
83        value: Scalar,
84        fraction: Scalar,
85    ) -> Result<Scalar, LineSearchError>;
86}
87
88fn exhausted(search: &impl Debug, max_steps: usize) -> LineSearchError {
89    LineSearchError::MaximumStepsReached(format!("{search:?}"), max_steps)
90}
91
92/// The value and the slope at no step, once the step and direction are vetted.
93fn start<X, J, D, E>(
94    search: &impl Debug,
95    function: &mut impl FnMut(&X, Scalar) -> Result<Scalar, String>,
96    argument: &X,
97    jacobian0: &J,
98    decrement: &D,
99    step_size: Scalar,
100) -> Result<(Scalar, Scalar), LineSearchError>
101where
102    J: Erase<Erased = E> + Jacobian,
103    D: Erase<Erased = E>,
104    E: Tensor,
105{
106    if step_size <= 0.0 {
107        return Err(LineSearchError::NegativeStepSize(
108            format!("{search:?}"),
109            step_size,
110        ));
111    }
112    let Ok(value) = function(argument, 0.0) else {
113        return Err(LineSearchError::InvalidStartingPoint(format!("{search:?}")));
114    };
115    let slope = jacobian0.erase().full_contraction(decrement.erase());
116    if slope <= 0.0 {
117        return Err(LineSearchError::NotDescentDirection(format!("{search:?}")));
118    }
119    Ok((value, slope))
120}
121
122/// Cut the step back until the decrease it buys is accepted.
123///
124/// The decrease is what the objective loses, so a step that cannot be
125/// evaluated is cut back as one that buys nothing.
126#[expect(clippy::too_many_arguments)]
127fn backtrack<X, J, D, W, E>(
128    search: &impl Debug,
129    cut_back: Scalar,
130    max_steps: usize,
131    mut function: impl FnMut(&X, Scalar) -> Result<Scalar, String>,
132    argument: &X,
133    jacobian0: &J,
134    decrement: &D,
135    step_size: Scalar,
136    accepted: impl Fn(Scalar, Scalar, Scalar) -> bool,
137) -> Result<Scalar, LineSearchError>
138where
139    J: Erase<Erased = E> + Jacobian,
140    D: Erase<Erased = E>,
141    E: Tensor,
142    X: Solution,
143    for<'a> &'a D: Mul<Quantity<W>, Output = X>,
144{
145    let (value, slope) = start(
146        search,
147        &mut function,
148        argument,
149        jacobian0,
150        decrement,
151        step_size,
152    )?;
153    let mut n = step_size;
154    for _ in 0..max_steps {
155        if let Ok(trial) = function(&(decrement * Quantity::new(-n) + argument), n)
156            && accepted(value - trial, n, slope)
157        {
158            return Ok(n);
159        } else {
160            n *= cut_back
161        }
162    }
163    Err(exhausted(search, max_steps))
164}
165
166/// Cut the step back until the decrease of the merit is accepted.
167#[expect(clippy::too_many_arguments)]
168fn backtrack_merit(
169    search: &impl Debug,
170    cut_back: Scalar,
171    max_steps: usize,
172    mut merit: impl FnMut(Scalar) -> Result<Scalar, String>,
173    value: Scalar,
174    slope: Scalar,
175    step_size: Scalar,
176    accepted: impl Fn(Scalar, Scalar, Scalar) -> bool,
177) -> Result<Scalar, LineSearchError> {
178    if step_size <= 0.0 {
179        return Err(LineSearchError::NegativeStepSize(
180            format!("{search:?}"),
181            step_size,
182        ));
183    }
184    if slope <= 0.0 {
185        return Err(LineSearchError::NotDescentDirection(format!("{search:?}")));
186    }
187    let mut n = step_size;
188    for _ in 0..max_steps {
189        if let Ok(trial) = merit(n)
190            && accepted(value - trial, n, slope)
191        {
192            return Ok(n);
193        } else {
194            n *= cut_back
195        }
196    }
197    Err(exhausted(search, max_steps))
198}
199
200/// Possible errors encountered during line search.
201pub enum LineSearchError {
202    InvalidStartingPoint(String),
203    MaximumStepsReached(String, usize),
204    NegativeStepSize(String, Scalar),
205    NotDescentDirection(String),
206}
207
208impl StyledError for LineSearchError {
209    fn message(&self, style: &Style) -> String {
210        let (h, c) = (style.headline, style.frame);
211        match self {
212            Self::InvalidStartingPoint(line_search) => format!(
213                "{h}Starting point is invalid.{c}\n\
214                In line search: {line_search}."
215            ),
216            Self::MaximumStepsReached(line_search, steps) => format!(
217                "{h}Maximum number of steps ({steps}) reached.{c}\n\
218                In line search: {line_search}."
219            ),
220            Self::NegativeStepSize(line_search, step_size) => format!(
221                "{h}Negative step size ({step_size}) encountered.{c}\n\
222                In line search: {line_search}."
223            ),
224            Self::NotDescentDirection(line_search) => format!(
225                "{h}Direction is not a descent direction.{c}\n\
226                In line search: {line_search}."
227            ),
228        }
229    }
230}
231
232styled_error!(LineSearchError);
conspire/math/optimize/globalization/line_search/mod.rs

conspire/math/optimize/globalization/line_search/
mod.rs
