Skip to main content

Crate ripopt

Crate ripopt 

Source
Expand description

§ripopt — Rust Interior Point Optimizer

ripopt is a primal-dual interior point solver for nonlinear programming (NLP) problems:

min  f(x)
 x
s.t. g_l ≤ g(x) ≤ g_u
     x_l ≤  x   ≤ x_u

It closely follows the algorithm of Ipopt and achieves comparable or better solve rates on standard benchmark suites (HS120, CUTEst).

§Quick Start

Implement NlpProblem for your problem, then call solve:

use ripopt::{NlpProblem, SolveResult, SolverOptions, SolveStatus, solve};

struct MyProblem;

impl NlpProblem for MyProblem {
    fn num_variables(&self) -> usize { 2 }
    fn num_constraints(&self) -> usize { 0 }
    fn bounds(&self, x_l: &mut [f64], x_u: &mut [f64]) {
        x_l.fill(f64::NEG_INFINITY);
        x_u.fill(f64::INFINITY);
    }
    fn constraint_bounds(&self, _g_l: &mut [f64], _g_u: &mut [f64]) {}
    fn initial_point(&self, x0: &mut [f64]) { x0[0] = 0.5; x0[1] = 0.5; }
    fn objective(&self, x: &[f64], _new_x: bool, obj: &mut f64) -> bool {
        *obj = (1.0 - x[0]).powi(2) + 100.0 * (x[1] - x[0].powi(2)).powi(2);
        true
    }
    fn gradient(&self, x: &[f64], _new_x: bool, grad: &mut [f64]) -> bool {
        grad[0] = -2.0 * (1.0 - x[0]) - 400.0 * x[0] * (x[1] - x[0].powi(2));
        grad[1] = 200.0 * (x[1] - x[0].powi(2));
        true
    }
    fn constraints(&self, _x: &[f64], _new_x: bool, _g: &mut [f64]) -> bool { true }
    fn jacobian_structure(&self) -> (Vec<usize>, Vec<usize>) { (vec![], vec![]) }
    fn jacobian_values(&self, _x: &[f64], _new_x: bool, _vals: &mut [f64]) -> bool { true }
    fn hessian_structure(&self) -> (Vec<usize>, Vec<usize>) {
        (vec![0, 1, 1], vec![0, 0, 1])
    }
    fn hessian_values(&self, x: &[f64], _new_x: bool, obj_factor: f64, _lambda: &[f64], vals: &mut [f64]) -> bool {
        vals[0] = obj_factor * (2.0 - 400.0 * (x[1] - x[0].powi(2)) + 800.0 * x[0].powi(2));
        vals[1] = obj_factor * (-400.0 * x[0]);
        vals[2] = obj_factor * 200.0;
        true
    }
}

let result = solve(&MyProblem, &SolverOptions::default());
assert_eq!(result.status, SolveStatus::Optimal);

§Key Types

§Algorithm

The solver implements:

  • Mehrotra predictor-corrector IPM with Gondzio corrections
  • Filter line search with second-order corrections
  • Gauss-Newton and NLP restoration phases
  • Fallback cascade: IPM → L-BFGS → Augmented Lagrangian → SQP → slack reformulation
  • Sparse (multifrontal LDL^T) and dense (Bunch-Kaufman LDL^T) linear solvers
  • Parametric sensitivity analysis (solve_with_sensitivity)

Re-exports§

pub use options::SolverOptions;
pub use options::LinearSolverChoice;
pub use problem::NlpProblem;
pub use result::SolveResult;
pub use result::SolverDiagnostics;
pub use result::SolveStatus;
pub use sensitivity::ParametricNlpProblem;
pub use sensitivity::SensitivityContext;
pub use sensitivity::SensitivityResult;

Modules§

augmented_lagrangian
Augmented Lagrangian method for constrained problems.
c_api
C API for ripopt — mirrors the Ipopt C interface.
convergence
filter
intermediate
Intermediate callback infrastructure for ripopt.
ipm
kkt
lbfgs
L-BFGS unconstrained minimizer.
linear_solver
linearity
Near-linear constraint detection.
nl
options
preprocessing
Preprocessing: eliminate fixed variables and redundant constraints before solving.
problem
restoration
restoration_nlp
result
sensitivity
Parametric sensitivity analysis for nonlinear programs.
slack_formulation
sqp
Simple SQP (Sequential Quadratic Programming) fallback solver.
trace
Per-iteration TSV tracing for direction-diff harness against Ipopt.
warmstart

Functions§

solve
Solve a nonlinear programming problem using the interior point method.
solve_with_sensitivity
Solve and retain factored KKT for parametric sensitivity analysis.