forge_core/

problem.rs

//! The unified optimization problem.
//!
//! [`Problem`] is the single interface every forge optimizer consumes: a
//! box-constrained, real-valued search space plus a scalar objective. By
//! convention forge **minimizes**; a candidate whose objective is non-finite
//! (`NaN`/`±∞`) is treated as infeasible and rejected, so models that blow up
//! on degenerate parameters (a common case in hydrological calibration) need no
//! special handling.
//!
//! Two adapters cover the everyday cases without writing a `struct`:
//!
//! - [`func`] wraps a closure with explicit bounds.
//! - [`Maximize`] flips the sense for objectives that should be *maximized*
//!   (e.g. NSE/KGE in rainflow), so the minimizing core stays the single
//!   convention.
//!
//! Combinatorial problems (Anvil's bit-flip simulated annealing) are served by
//! a separate abstraction migrated in a later milestone; v0.1 is the
//! continuous, population/real-vector substrate that unblocks rainflow.

/// Inclusive search bounds for one decision variable: `lower <= x <= upper`.
pub type Bound = (f64, f64);

/// A box-constrained, real-valued minimization problem.
pub trait Problem {
    /// Number of decision variables.
    fn dim(&self) -> usize;

    /// Per-variable inclusive bounds, length [`Problem::dim`].
    fn bounds(&self) -> &[Bound];

    /// Objective value to **minimize**. A non-finite return marks the point as
    /// infeasible; optimizers will reject it rather than crash.
    fn objective(&self, x: &[f64]) -> f64;
}

/// A box-constrained, real-valued **multi-objective** minimization problem.
///
/// All objectives are minimized (negate any to be maximized). The
/// multi-objective optimizers ([`NsgaII`], [`NsgaIII`], [`SmsEmoa`]) return a
/// Pareto front rather than a single solution. As with [`Problem`], a
/// non-finite objective value marks a point as infeasible — it is dominated by
/// every feasible point.
///
/// [`NsgaII`]: crate::algo::NsgaII
/// [`NsgaIII`]: crate::algo::NsgaIII
/// [`SmsEmoa`]: crate::algo::SmsEmoa
pub trait MultiProblem {
    /// Number of decision variables.
    fn dim(&self) -> usize;

    /// Per-variable inclusive bounds, length [`MultiProblem::dim`].
    fn bounds(&self) -> &[Bound];

    /// Number of objectives (≥ 2 for a genuine multi-objective problem).
    fn n_objectives(&self) -> usize;

    /// The objective vector to minimize, length [`MultiProblem::n_objectives`].
    fn objectives(&self, x: &[f64]) -> Vec<f64>;
}

/// A [`MultiProblem`] defined inline by a closure, a bounds vector, and an
/// objective count.
///
/// ```
/// use forge_core::problem::{multi_func, MultiProblem};
/// // Schaffer N.1: minimize x² and (x−2)² on [-5, 5].
/// let sch = multi_func(vec![(-5.0, 5.0)], 2, |x| vec![x[0] * x[0], (x[0] - 2.0).powi(2)]);
/// assert_eq!(sch.n_objectives(), 2);
/// assert_eq!(sch.objectives(&[1.0]), vec![1.0, 1.0]);
/// ```
pub fn multi_func<F>(bounds: Vec<Bound>, n_objectives: usize, f: F) -> MultiFunc<F>
where
    F: Fn(&[f64]) -> Vec<f64>,
{
    MultiFunc {
        bounds,
        n_objectives,
        f,
    }
}

/// Closure-backed multi-objective problem produced by [`multi_func`].
pub struct MultiFunc<F> {
    bounds: Vec<Bound>,
    n_objectives: usize,
    f: F,
}

impl<F> MultiProblem for MultiFunc<F>
where
    F: Fn(&[f64]) -> Vec<f64>,
{
    fn dim(&self) -> usize {
        self.bounds.len()
    }
    fn bounds(&self) -> &[Bound] {
        &self.bounds
    }
    fn n_objectives(&self) -> usize {
        self.n_objectives
    }
    fn objectives(&self, x: &[f64]) -> Vec<f64> {
        (self.f)(x)
    }
}

/// Validates that a problem is well-formed: at least one variable and every
/// bound strictly ordered (`lower < upper`), rejecting `NaN` bounds.
pub fn validate(problem: &dyn Problem) -> Result<(), BoundsError> {
    let b = problem.bounds();
    if b.is_empty() || problem.dim() == 0 {
        return Err(BoundsError::Empty);
    }
    if b.len() != problem.dim() {
        return Err(BoundsError::DimMismatch {
            dim: problem.dim(),
            bounds: b.len(),
        });
    }
    for (i, &(lo, hi)) in b.iter().enumerate() {
        // `partial_cmp` rejects NaN too, unlike a plain `lo >= hi`.
        if lo.partial_cmp(&hi) != Some(std::cmp::Ordering::Less) {
            return Err(BoundsError::NotOrdered { dim: i, lo, hi });
        }
    }
    Ok(())
}

/// Validates that a multi-objective problem is well-formed: valid bounds (as
/// in [`validate`]) plus at least two objectives, and an `objectives()` vector
/// whose length matches `n_objectives()` at the box midpoint.
pub fn validate_multi(problem: &dyn MultiProblem) -> Result<(), BoundsError> {
    struct AsProblem<'a>(&'a dyn MultiProblem);
    impl Problem for AsProblem<'_> {
        fn dim(&self) -> usize {
            self.0.dim()
        }
        fn bounds(&self) -> &[Bound] {
            self.0.bounds()
        }
        fn objective(&self, _x: &[f64]) -> f64 {
            0.0
        }
    }
    validate(&AsProblem(problem))?;
    let m = problem.n_objectives();
    if m < 2 {
        return Err(BoundsError::BadObjectiveCount {
            declared: m,
            got: m,
        });
    }
    // Probe once at the box midpoint: a closure returning the wrong number of
    // objectives would otherwise panic deep inside the non-dominated sort.
    let mid: Vec<f64> = problem
        .bounds()
        .iter()
        .map(|&(lo, hi)| 0.5 * (lo + hi))
        .collect();
    let got = problem.objectives(&mid).len();
    if got != m {
        return Err(BoundsError::BadObjectiveCount { declared: m, got });
    }
    Ok(())
}

/// Why a problem's bounds are invalid.
#[derive(Debug, Clone, PartialEq)]
#[non_exhaustive]
pub enum BoundsError {
    /// No decision variables were given.
    Empty,
    /// `bounds().len()` disagrees with `dim()`.
    DimMismatch {
        /// What `dim()` reported.
        dim: usize,
        /// How many bound pairs `bounds()` returned.
        bounds: usize,
    },
    /// Dimension `dim` has `lower >= upper` (or a NaN bound).
    NotOrdered {
        /// Index of the offending variable.
        dim: usize,
        /// Its lower bound.
        lo: f64,
        /// Its upper bound.
        hi: f64,
    },
    /// `objectives()` length disagrees with `n_objectives()` (or fewer than 2).
    BadObjectiveCount {
        /// What `n_objectives()` declared.
        declared: usize,
        /// What `objectives()` actually returned.
        got: usize,
    },
}

impl std::fmt::Display for BoundsError {
    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
        match self {
            BoundsError::Empty => write!(f, "problem must have at least one variable"),
            BoundsError::DimMismatch { dim, bounds } => {
                write!(f, "dim() = {dim} but bounds() has {bounds} entries")
            }
            BoundsError::NotOrdered { dim, lo, hi } => {
                write!(
                    f,
                    "dimension {dim}: lower bound ({lo}) must be < upper bound ({hi})"
                )
            }
            BoundsError::BadObjectiveCount { declared, got } => {
                write!(
                    f,
                    "n_objectives() declares {declared} but objectives() returned {got} (need >= 2)"
                )
            }
        }
    }
}

impl std::error::Error for BoundsError {}

/// A [`Problem`] defined inline by a closure and a bounds vector.
///
/// ```
/// use forge_core::problem::{func, Problem};
/// let sphere = func(vec![(-5.0, 5.0); 3], |x| x.iter().map(|v| v * v).sum());
/// assert_eq!(sphere.dim(), 3);
/// assert_eq!(sphere.objective(&[0.0, 0.0, 0.0]), 0.0);
/// ```
pub fn func<F>(bounds: Vec<Bound>, f: F) -> Func<F>
where
    F: Fn(&[f64]) -> f64,
{
    Func { bounds, f }
}

/// Closure-backed problem produced by [`func`].
pub struct Func<F> {
    bounds: Vec<Bound>,
    f: F,
}

impl<F> Problem for Func<F>
where
    F: Fn(&[f64]) -> f64,
{
    fn dim(&self) -> usize {
        self.bounds.len()
    }
    fn bounds(&self) -> &[Bound] {
        &self.bounds
    }
    fn objective(&self, x: &[f64]) -> f64 {
        (self.f)(x)
    }
}

/// Wraps a problem so the engine **maximizes** its objective instead of
/// minimizing it, by negating finite values (non-finite stay rejected).
///
/// Lets a maximization client (e.g. rainflow maximizing NSE) reuse the
/// minimizing core unchanged. The reported objective in a [`Solution`] is the
/// negated, internal value; recover the original sense by negating it back, or
/// read [`crate::Report::best_value_maximized`].
///
/// [`Solution`]: crate::Solution
pub struct Maximize<P>(pub P);

impl<P: Problem> Problem for Maximize<P> {
    fn dim(&self) -> usize {
        self.0.dim()
    }
    fn bounds(&self) -> &[Bound] {
        self.0.bounds()
    }
    fn objective(&self, x: &[f64]) -> f64 {
        let v = self.0.objective(x);
        if v.is_finite() {
            -v
        } else {
            v
        }
    }
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn func_adapter_works() {
        let p = func(vec![(-1.0, 1.0), (-1.0, 1.0)], |x| x[0] + x[1]);
        assert_eq!(p.dim(), 2);
        assert_eq!(p.objective(&[0.3, 0.4]), 0.7);
    }

    #[test]
    fn maximize_negates_finite_only() {
        let p = Maximize(func(vec![(0.0, 1.0)], |x| x[0]));
        assert_eq!(p.objective(&[0.6]), -0.6);
        let bad = Maximize(func(vec![(0.0, 1.0)], |_| f64::NAN));
        assert!(bad.objective(&[0.5]).is_nan());
    }

    #[test]
    fn validate_catches_bad_bounds() {
        assert_eq!(
            validate(&func(vec![], |_| 0.0)).unwrap_err(),
            BoundsError::Empty
        );
        assert!(matches!(
            validate(&func(vec![(1.0, 0.0)], |_| 0.0)).unwrap_err(),
            BoundsError::NotOrdered { dim: 0, .. }
        ));
        assert!(matches!(
            validate(&func(vec![(0.0, f64::NAN)], |_| 0.0)).unwrap_err(),
            BoundsError::NotOrdered { .. }
        ));
        assert!(validate(&func(vec![(0.0, 1.0), (-2.0, 2.0)], |_| 0.0)).is_ok());
    }
}