diff --git a/src/algebra/csc/core.rs b/src/algebra/csc/core.rs index f6e48971..dd5be1e3 100644 --- a/src/algebra/csc/core.rs +++ b/src/algebra/csc/core.rs @@ -51,7 +51,7 @@ pub struct CscMatrix { pub n: usize, /// CSC format column pointer. /// - /// Ths field should have length `n+1`. The last entry corresponds + /// This field should have length `n+1`. The last entry corresponds /// to the the number of nonzeros and should agree with the lengths /// of the `rowval` and `nzval` fields. pub colptr: Vec, @@ -419,7 +419,7 @@ where self.size() == other.size() && self.colptr == other.colptr && self.rowval == other.rowval } - /// Same as is_equal_sparsity, but returns an error indicating the reason + /// Same as `is_equal_sparsity`, but returns an error indicating the reason /// for failure if the matrices do not have equivalent sparsity patterns. pub fn check_equal_sparsity(&self, other: &Self) -> Result<(), SparseFormatError> { if self.size() != other.size() { diff --git a/src/algebra/dense/types.rs b/src/algebra/dense/types.rs index 732b7691..4a4358ac 100644 --- a/src/algebra/dense/types.rs +++ b/src/algebra/dense/types.rs @@ -10,7 +10,7 @@ where { /// dimensions pub size: (usize, usize), - /// vector of data in column major formmat + /// vector of data in column major format pub data: S, pub(crate) phantom: std::marker::PhantomData, } diff --git a/src/algebra/error_types.rs b/src/algebra/error_types.rs index 5baa755c..7000449d 100644 --- a/src/algebra/error_types.rs +++ b/src/algebra/error_types.rs @@ -4,21 +4,27 @@ use thiserror::Error; #[derive(Error, Debug)] pub enum MatrixConcatenationError { #[error("Incompatible dimensions")] + /// Indicates inputs have incompatible dimension IncompatibleDimension, } #[derive(Error, Debug)] /// Error type returned by sparse matrix assembly operations. pub enum SparseFormatError { + /// Matrix dimension fields and/or array lengths are incompatible #[error("Matrix dimension fields and/or array lengths are incompatible")] IncompatibleDimension, + /// Data is not sorted by row index within each column #[error("Data is not sorted by row index within each column")] BadRowOrdering, #[error("Row value exceeds the matrix row dimension")] + /// Row value exceeds the matrix row dimension BadRowval, #[error("Bad column pointer values")] + /// Matrix column pointer values are defective BadColptr, #[error("sparsity pattern mismatch")] + /// Operation on matrices that have mismatching sparsity patterns SparsityMismatch, } diff --git a/src/algebra/floats.rs b/src/algebra/floats.rs index e1d53c0e..c8501a01 100644 --- a/src/algebra/floats.rs +++ b/src/algebra/floats.rs @@ -8,8 +8,10 @@ use crate::algebra::dense::BlasFloatT; /// Core traits for internal floating point values. /// /// This trait defines a subset of bounds for `FloatT`, which is preferred -/// throughout for use in the solver. When the "sdp" feature is enabled, +/// throughout for use in the solver. When the `sdp` feature is enabled, /// `FloatT` is additionally restricted to f32/f64 types supported by BLAS. +/// When the `faer-sparse` feature is enabled, `FloatT` is additionally +/// restricted to types implementing `RealField` from the `faer` crate. pub trait CoreFloatT: 'static + Send @@ -51,24 +53,32 @@ impl CoreFloatT for T where cfg_if::cfg_if! { if #[cfg(feature="sdp")] { + /// A marker trait implemented on types supported by BLAS (i.e. f32 and f64) + /// when the package is compiled with the "sdp" feature using a BLAS/LAPACK library + #[doc(hidden)] pub trait MaybeBlasFloatT : BlasFloatT {} impl MaybeBlasFloatT for T where T: BlasFloatT {} } else { + #[doc(hidden)] pub trait MaybeBlasFloatT {} impl MaybeBlasFloatT for T {} } } // if "faer" is enabled, we must add an additional bound -// to restrict compilation to Reals implementing SimpleEntity +// to restrict compilation to types implementing RealField cfg_if::cfg_if! { if #[cfg(feature="faer-sparse")] { + #[doc(hidden)] + /// A marker trait implemented on types supported by faer-rs + /// when the package is compiled with the "faer-sparse" feature pub trait MaybeFaerFloatT : faer_traits::RealField {} impl MaybeFaerFloatT for T where T: faer_traits::RealField {} } else { + #[doc(hidden)] pub trait MaybeFaerFloatT {} impl MaybeFaerFloatT for T {} } @@ -86,7 +96,7 @@ cfg_if::cfg_if! { pub trait FloatT: CoreFloatT + MaybeBlasFloatT + MaybeFaerFloatT {} impl FloatT for T where T: CoreFloatT + MaybeBlasFloatT + MaybeFaerFloatT {} -/// Trait for convering Rust primitives to [`FloatT`](crate::algebra::FloatT) +/// Trait for converting Rust primitives to [`FloatT`](crate::algebra::FloatT) /// /// This convenience trait is implemented on f32/64 and u32/64. This trait /// is required internally by the solver for converting constant primitives @@ -97,7 +107,7 @@ impl FloatT for T where T: CoreFloatT + MaybeBlasFloatT + MaybeFaerFloatT {} // NB: `AsFloatT` is a convenience trait for f32/64 and u32/64 // so that we can do things like (2.0).as_T() everywhere on // constants, rather than the awful T::from_f32(2.0).unwrap() -pub trait AsFloatT: 'static { +pub(crate) trait AsFloatT: 'static { fn as_T(&self) -> T; } diff --git a/src/algebra/math_traits.rs b/src/algebra/math_traits.rs index 002e952c..e775f961 100644 --- a/src/algebra/math_traits.rs +++ b/src/algebra/math_traits.rs @@ -63,7 +63,7 @@ pub trait VectorMath { /// Dot product fn dot(&self, y: &Self) -> T; - // computes dot(z + αdz,s + αds) without intermediate allocation + /// computes dot(z + αdz,s + αds) without intermediate allocation fn dot_shifted(z: &[T], s: &[T], dz: &[T], ds: &[T], α: T) -> T; /// Standard Euclidian or 2-norm distance from `self` to `y` @@ -93,7 +93,7 @@ pub trait VectorMath { /// Inf-norm of an elementwise scaling of `self` by `v` fn norm_one_scaled(&self, v: &Self) -> T; - // Inf-norm of vector difference + /// Inf-norm of vector difference fn norm_inf_diff(&self, b: &Self) -> T; /// Minimum value in vector @@ -152,7 +152,7 @@ pub trait MatrixMath { fn col_norms_no_reset(&self, norms: &mut [T]); /// Compute columnwise infinity norm operations on - /// a symmstric matrix + /// a symmetric matrix fn col_norms_sym(&self, norms: &mut [T]); /// Compute columnwise infinity norm operations on diff --git a/src/algebra/matrix_traits.rs b/src/algebra/matrix_traits.rs index 3dc33b44..eec5e2e3 100644 --- a/src/algebra/matrix_traits.rs +++ b/src/algebra/matrix_traits.rs @@ -39,7 +39,7 @@ pub trait BlockConcatenate: Sized { /// Errors if column dimensions are incompatible fn vcat(A: &Self, B: &Self) -> Result; - /// general block concatentation + /// general block concatenation fn hvcat(mats: &[&[&Self]]) -> Result; /// block diagonal concatenation diff --git a/src/algebra/sparsevector/mod.rs b/src/algebra/sparsevector/mod.rs index d59c5cfc..9c5aa91c 100644 --- a/src/algebra/sparsevector/mod.rs +++ b/src/algebra/sparsevector/mod.rs @@ -15,7 +15,7 @@ pub(crate) struct SparseVector { pub nzval: Vec, } -/// Creates a SparseVector from a dense slice. +/// Creates a `SparseVector` from a dense slice. impl SparseVector where T: Num + Copy, diff --git a/src/lib.rs b/src/lib.rs index 40cc3edb..75025468 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -55,6 +55,7 @@ //Rust hates greek characters #![allow(confusable_idents)] +#![warn(missing_docs)] const VERSION: &str = env!("CARGO_PKG_VERSION"); @@ -64,6 +65,9 @@ pub mod solver; pub(crate) mod stdio; pub mod timers; +pub(crate) mod utils; +pub use crate::utils::infbounds::*; + #[cfg(feature = "python")] pub mod python; @@ -79,7 +83,7 @@ macro_rules! printbuildenv { }; } -// print detailed build info to stdout +/// print detailed build configuration info to stdout #[allow(clippy::explicit_write)] pub fn buildinfo() { use std::io::Write; @@ -109,6 +113,8 @@ pub fn buildinfo() { writeln!(crate::stdio::stdout(), "no build info available").unwrap(); } +pub(crate) const _INFINITY_DEFAULT: f64 = 1e20; + #[test] fn test_buildinfo() { buildinfo(); diff --git a/src/python/mod.rs b/src/python/mod.rs index d80c74c8..755a4d46 100644 --- a/src/python/mod.rs +++ b/src/python/mod.rs @@ -1,7 +1,7 @@ //! Clarabel Python interface. //! //! This module implements a wrapper for the Rust version of Python using -//! [PYO3](https://pyo3.rs/). To build these wrappers from `cargo`, compile the crate with +//! [PyO3](https://pyo3.rs/). To build these wrappers from `cargo`, compile the crate with //! `--features python`. This module has no public API. //! //! It should not normally be necessary to compile the Python wrapper from diff --git a/src/python/module_py.rs b/src/python/module_py.rs index d4ecb00e..b75dd33e 100644 --- a/src/python/module_py.rs +++ b/src/python/module_py.rs @@ -12,15 +12,15 @@ fn force_load_blas_lapack_py() { // get/set for the solver's internal infinity limit #[pyfunction(name = "get_infinity")] fn get_infinity_py() -> f64 { - crate::solver::get_infinity() + crate::get_infinity() } #[pyfunction(name = "set_infinity")] fn set_infinity_py(v: f64) { - crate::solver::set_infinity(v); + crate::set_infinity(v); } #[pyfunction(name = "default_infinity")] fn default_infinity_py() { - crate::solver::default_infinity(); + crate::default_infinity(); } #[pyfunction(name = "buildinfo")] fn buildinfo_py() { diff --git a/src/qdldl/qdldl.rs b/src/qdldl/qdldl.rs index 7171efbc..2c55adf2 100644 --- a/src/qdldl/qdldl.rs +++ b/src/qdldl/qdldl.rs @@ -6,39 +6,56 @@ use std::iter::zip; use thiserror::Error; /// Error codes returnable from [`QDLDLFactorisation`](QDLDLFactorisation) factor operations - #[derive(Error, Debug)] pub enum QDLDLError { #[error("Matrix dimension fields are incompatible")] + /// Matrix dimension fields are incompatible IncompatibleDimension, #[error("Matrix has a zero column")] + /// Matrix has a zero column EmptyColumn, #[error("Matrix is not upper triangular")] + /// Matrix is not upper triangular NotUpperTriangular, + /// Matrix factorization produced a zero pivot #[error("Matrix factorization produced a zero pivot")] ZeroPivot, #[error("Invalid permutation vector")] + /// Invalid permutation vector supplied InvalidPermutation, } -/// Required settings for [`QDLDLFactorisation`](QDLDLFactorisation) - #[derive(Builder, Debug, Clone)] +#[allow(missing_docs)] +/// Required settings for [`QDLDLFactorisation`](QDLDLFactorisation) pub struct QDLDLSettings { + /// "dense scale" parameter for AMD ordering #[builder(default = "1.0")] - amd_dense_scale: f64, + pub amd_dense_scale: f64, + + /// optional ueser-supplied custom permutation vector for the matrix #[builder(default = "None", setter(strip_option))] - perm: Option>, + pub perm: Option>, + + /// Logical factorisation only, no numerical factorisation #[builder(default = "false")] - logical: bool, + pub logical: bool, + + /// optional user-supplied signs of the diagonal elements of D in LDL^T #[builder(default = "None", setter(strip_option))] - Dsigns: Option>, + pub Dsigns: Option>, + + /// Enable regularization during factorisation #[builder(default = "true")] - regularize_enable: bool, + pub regularize_enable: bool, + + /// Regularization epsilon parameter #[builder(default = "(1e-12).as_T()")] - regularize_eps: T, + pub regularize_eps: T, + + /// Regularization delta parameter #[builder(default = "(1e-7).as_T()")] - regularize_delta: T, + pub regularize_delta: T, } impl Default for QDLDLSettings @@ -51,29 +68,30 @@ where } /// Performs $LDL^T$ factorization of a symmetric quasidefinite matrix - #[derive(Debug)] pub struct QDLDLFactorisation { - // permutation vector + /// permutation vector pub perm: Vec, - // inverse permutation + /// inverse permutation #[allow(dead_code)] //Unused because we call ipermute in solve instead. Keep anyway. iperm: Vec, - // lower triangular factor + /// lower triangular factor L in LDL^T pub L: CscMatrix, - // D and is inverse for A = LDL^T + /// vector of diagonal elements of D in LDL^T pub D: Vec, + /// vector of reciprocal diagonal elements of D in LDL^T pub Dinv: Vec, - // workspace data + /// internal workspace data workspace: QDLDLWorkspace, - // is it logical factorisation only? - is_logical: bool, + /// true of factorisation is symbolic only + is_symbolic: bool, } impl QDLDLFactorisation where T: FloatT, { + /// create a new LDL^T factorisation pub fn new( Ain: &CscMatrix, opts: Option>, @@ -83,18 +101,21 @@ where _qdldl_new(Ain, opts) } + /// returns the number of positive eigenvalues pub fn positive_inertia(&self) -> usize { self.workspace.positive_inertia } + /// returns the number of regularisation shifts + /// that were applied during factorisation pub fn regularize_count(&self) -> usize { self.workspace.regularize_count } - // Solves Ax = b using LDL factors for A. - // Solves in place (x replaces b) + /// Solves Ax = b using LDL factors for A. + /// Solves in place (x replaces b) pub fn solve(&mut self, b: &mut [T]) { // bomb if logical factorisation only - assert!(!self.is_logical); + assert!(!self.is_symbolic); // bomb if b is the wrong size assert_eq!(b.len(), self.D.len()); @@ -116,6 +137,8 @@ where ipermute(b, tmp, &self.perm); } + /// Update a subset of the values of the matrix to be (re)factored. See [`refactor`](crate::qdldl::QDLDLFactorisation::refactor) + /// pub fn update_values(&mut self, indices: &[usize], values: &[T]) { let nzval = &mut self.workspace.triuA.nzval; // post perm internal data let AtoPAPt = &self.workspace.AtoPAPt; //mapping from input matrix entries to triuA @@ -125,6 +148,8 @@ where } } + /// Update a subset of the values of the matrix to be (re)factored. See [`refactor`](crate::qdldl::QDLDLFactorisation::refactor) + /// pub fn scale_values(&mut self, indices: &[usize], scale: T) { let nzval = &mut self.workspace.triuA.nzval; // post perm internal data let AtoPAPt = &self.workspace.AtoPAPt; //mapping from input matrix entries to triuA @@ -134,6 +159,10 @@ where } } + /// Shifts a subset of the values of the matrix to be (re)factored. The + /// values are offset by `offset`, with `signs` a vector of +/- 1 values + /// indicating the direction of shifts. See [`refactor`](crate::qdldl::QDLDLFactorisation::refactor) + /// pub fn offset_values(&mut self, indices: &[usize], offset: T, signs: &[i8]) { assert_eq!(indices.len(), signs.len()); @@ -141,22 +170,32 @@ where let AtoPAPt = &self.workspace.AtoPAPt; //mapping from input matrix entries to triuA for (&idx, &sign) in zip(indices, signs) { - let sign: T = T::from_i8(sign).unwrap(); - nzval[AtoPAPt[idx]] += offset * sign; + match sign.signum() { + 1 => { + nzval[AtoPAPt[idx]] += offset; + } + -1 => { + nzval[AtoPAPt[idx]] -= offset; + } + _ => {} + } } } + /// Refactor a matrix after its data has been modified. See [`update_values`](crate::qdldl::QDLDLFactorisation::update_values), + /// [`scale_values`](crate::qdldl::QDLDLFactorisation::scale_values) and [`offset_values`](crate::qdldl::QDLDLFactorisation::offset_values) + /// pub fn refactor(&mut self) -> Result<(), QDLDLError> { // It never makes sense to call refactor for a logical // factorization since it will always be the same. Calling // this function implies that we want a numerical factorization - self.is_logical = false; + self.is_symbolic = false; _factor( &mut self.L, &mut self.D, &mut self.Dinv, &mut self.workspace, - self.is_logical, + self.is_symbolic, ) } } @@ -241,7 +280,7 @@ fn _qdldl_new( D, Dinv, workspace, - is_logical: opts.logical, + is_symbolic: opts.logical, }) } diff --git a/src/solver/core/solver.rs b/src/solver/core/solver.rs index fe7ee8e1..88a9ef92 100644 --- a/src/solver/core/solver.rs +++ b/src/solver/core/solver.rs @@ -97,7 +97,9 @@ pub trait SolverJSONReadWrite: Sized where T: FloatT, { + /// write internal problem data to a JSON file fn write_to_file(&self, file: &mut std::fs::File) -> Result<(), std::io::Error>; + /// load problem data from a JSON file previously saved using [`write_to_file`](self::SolverJSONReadWrite::write_to_file) fn read_from_file( file: &mut std::fs::File, settings: Option>, diff --git a/src/solver/core/traits.rs b/src/solver/core/traits.rs index 5d5a309f..d3b3511b 100644 --- a/src/solver/core/traits.rs +++ b/src/solver/core/traits.rs @@ -17,8 +17,11 @@ use crate::timers::*; /// Data for a conic optimization problem. pub trait ProblemData { + /// associated variable type type V: Variables; + /// associated cone type type C: Cone; + /// associated settings type type SE: Settings; /// Equilibrate internal data before solver starts. @@ -27,9 +30,13 @@ pub trait ProblemData { /// Variables for a conic optimization problem. pub trait Variables { + /// associated problem data type type D: ProblemData; + /// associated problem residuals type type R: Residuals; + /// associated cone type type C: Cone; + /// associated settings type type SE: Settings; /// Compute the scaled duality gap. @@ -86,7 +93,9 @@ pub trait Variables { /// Residuals for a conic optimization problem. pub trait Residuals { + /// associated problem data type type D: ProblemData; + /// associated variable type type V: Variables; /// Compute residuals for the given variables. @@ -96,9 +105,13 @@ pub trait Residuals { /// KKT linear solver object. pub trait KKTSystem { + /// associated problem data type type D: ProblemData; + /// associated variable type type V: Variables; + /// associated cone type type C: Cone; + /// associated settings type type SE: Settings; /// Update the KKT system. In particular, update KKT @@ -132,8 +145,11 @@ pub trait InfoPrint where T: FloatT, { + /// associated problem data type type D: ProblemData; + /// associated cone type type C: Cone; + /// associated settings type type SE: Settings; /// Print the solver configuration, e.g. settings etc. @@ -161,7 +177,9 @@ pub trait Info: InfoPrint where T: FloatT, { + /// associated variables type type V: Variables; + /// associated problem residuals type type R: Residuals; /// Reset internal data, particularly solve timers. @@ -185,8 +203,9 @@ where /// Return `true` if termination conditions have been reached. fn check_termination(&mut self, residuals: &Self::R, settings: &Self::SE, iter: u32) -> bool; - // save and recover prior iterates + /// save a prior iterate fn save_prev_iterate(&mut self, variables: &Self::V, prev_variables: &mut Self::V); + /// restore a prior iterate fn reset_to_prev_iterate(&mut self, variables: &mut Self::V, prev_variables: &Self::V); /// Record some of the top level solver's choice of various @@ -194,16 +213,21 @@ where /// `σ = ` multiplier for the updated centering parameter. fn save_scalars(&mut self, μ: T, α: T, σ: T, iter: u32); - /// Report or update termination status + /// Report the termination status fn get_status(&self) -> SolverStatus; + /// Set the termination status fn set_status(&mut self, status: SolverStatus); } /// Solution for a conic optimization problem. pub trait Solution { + /// Associated problem data type type D: ProblemData; + /// Associated problem variable type type V: Variables; + /// Associated progress information type type I: Info; + /// Associated solver settings settings type SE: Settings; /// Compute solution from the Variables at solver termination @@ -221,7 +245,7 @@ pub trait Solution { /// Settings for a conic optimization problem. /// -/// Implementors of this trait can define any internal or problem +/// Implementers of this trait can define any internal or problem /// specific settings they wish. They must, however, also maintain /// a settings object of type [`CoreSettings`](crate::solver::core::CoreSettings) /// and return this to the solver internally. diff --git a/src/solver/implementations/default/data_updating.rs b/src/solver/implementations/default/data_updating.rs index bb3a1a58..e56ec14c 100644 --- a/src/solver/implementations/default/data_updating.rs +++ b/src/solver/implementations/default/data_updating.rs @@ -9,16 +9,20 @@ use thiserror::Error; #[derive(Error, Debug)] pub enum DataUpdateError { #[error("Data updates are not allowed when presolver is active")] + /// Data updates are not allowed when presolver is active PresolveIsActive, #[cfg(feature = "sdp")] #[error("Data updates are not allowed when chordal decomposition is active")] + /// Data updates are not allowed when chordal decomposition is active ChordalDecompositionIsActive, #[error("Data formatting error")] + /// Data formatting error. See [`SparseFormatError`] BadFormat(#[from] SparseFormatError), } /// Trait for updating problem data matrices (`P` and `A`) from various data types pub trait MatrixProblemDataUpdate { + /// Update matrix entries using associated left/right conditioners and scaling terms fn update_matrix( &self, M: &mut CscMatrix, @@ -30,6 +34,7 @@ pub trait MatrixProblemDataUpdate { /// Trait for updating problem data vectors (`q`` and `b`) from various data types pub trait VectorProblemDataUpdate { + /// Update vector entries using associated elementwise and overall scaling terms fn update_vector( &self, v: &mut [T], @@ -42,6 +47,8 @@ impl DefaultSolver where T: FloatT, { + // PJG: rustdoc fails to resolve links to `update_P`, `update_q`, `update_A`, `update_b` below + /// Overwrites internal problem data structures in a solver object with new data, avoiding new memory allocations. /// See `update_P`, `update_q`, `update_A`, `update_b` for allowable inputs. /// diff --git a/src/solver/implementations/default/equilibration.rs b/src/solver/implementations/default/equilibration.rs index 54d7ef19..813dde0c 100644 --- a/src/solver/implementations/default/equilibration.rs +++ b/src/solver/implementations/default/equilibration.rs @@ -10,12 +10,15 @@ pub struct DefaultEquilibrationData { // scaling matrices for problem data equilibration // fields d,e,dinv,einv are vectors of scaling values // to be treated as diagonal scaling data + /// Vector of variable scaling terms pub d: Vec, + /// Vector of inverse variable scaling terms pub dinv: Vec, + /// Vector of constraint scaling terms pub e: Vec, + /// Vector of inverse constraint scaling terms pub einv: Vec, - - // overall scaling for objective function + /// overall scaling for objective function pub c: T, } @@ -23,6 +26,7 @@ impl DefaultEquilibrationData where T: FloatT, { + /// creates a new equilibration object pub fn new(n: usize, m: usize) -> Self { // Left/Right diagonal scaling for problem data let d = vec![T::one(); n]; diff --git a/src/solver/implementations/default/info.rs b/src/solver/implementations/default/info.rs index 1c2deaa5..9c2abc9c 100644 --- a/src/solver/implementations/default/info.rs +++ b/src/solver/implementations/default/info.rs @@ -8,29 +8,49 @@ use crate::timers::*; #[repr(C)] #[derive(Default, Debug, Clone)] pub struct DefaultInfo { + /// interior point path parameter μ pub μ: T, + /// interior point path parameter reduction ratio σ pub sigma: T, + /// step length for the current iteration pub step_length: T, + /// number of iterations pub iterations: u32, + /// primal objective value pub cost_primal: T, + /// dual objective value pub cost_dual: T, + /// primal residual pub res_primal: T, + /// dual residual pub res_dual: T, + /// primal infeasibility residual pub res_primal_inf: T, + /// dual infeasibility residual pub res_dual_inf: T, + /// absolute duality gap pub gap_abs: T, + /// relative duality gap pub gap_rel: T, + /// κ/τ ratio pub ktratio: T, // previous iterate + /// primal object value from previous iteration prev_cost_primal: T, + /// dual objective value from previous iteration prev_cost_dual: T, + /// primal residual from previous iteration prev_res_primal: T, + /// dual residual from previous iteration prev_res_dual: T, + /// absolute duality gap from previous iteration prev_gap_abs: T, + /// relative duality gap from previous iteration prev_gap_rel: T, - + /// solve time pub solve_time: f64, + /// solver status pub status: SolverStatus, } @@ -38,6 +58,7 @@ impl DefaultInfo where T: FloatT, { + /// creates a new `DefaultInfo` object pub fn new() -> Self { Self::default() } diff --git a/src/solver/implementations/default/kktsystem.rs b/src/solver/implementations/default/kktsystem.rs index e8f0176e..b01eac4e 100644 --- a/src/solver/implementations/default/kktsystem.rs +++ b/src/solver/implementations/default/kktsystem.rs @@ -35,6 +35,7 @@ impl DefaultKKTSystem where T: FloatT, { + /// Create a new KKT system solver pub fn new( data: &DefaultProblemData, cones: &CompositeCone, diff --git a/src/solver/implementations/default/mod.rs b/src/solver/implementations/default/mod.rs index f0cea8c8..75edc814 100644 --- a/src/solver/implementations/default/mod.rs +++ b/src/solver/implementations/default/mod.rs @@ -21,7 +21,7 @@ pub use data_updating::*; pub use equilibration::*; pub use info::*; pub use kktsystem::*; -pub use presolver::*; +pub(crate) use presolver::*; pub use problemdata::*; pub use residuals::*; pub use settings::*; diff --git a/src/solver/implementations/default/presolver.rs b/src/solver/implementations/default/presolver.rs index 46ab0cf1..9f59c7ae 100644 --- a/src/solver/implementations/default/presolver.rs +++ b/src/solver/implementations/default/presolver.rs @@ -18,7 +18,7 @@ pub(crate) struct PresolverRowReductionIndex { /// Presolver data for the standard solver implementation #[derive(Debug)] -pub struct Presolver { +pub(crate) struct Presolver { // original cones of the problem // PJG: not currently used. Here for future presolver pub(crate) _init_cones: Vec>, @@ -41,13 +41,14 @@ impl Presolver where T: FloatT, { - pub fn new( + /// create a new presolver object + pub(crate) fn new( _A: &CscMatrix, b: &[T], cones: &[SupportedConeT], _settings: &DefaultSettings, ) -> Self { - let infbound = crate::solver::get_infinity(); + let infbound = crate::get_infinity(); // make copy of cones to protect from user interference let init_cones = cones.to_vec(); @@ -64,10 +65,12 @@ where } } - pub fn is_reduced(&self) -> bool { + /// true if the presolver has reduced the problem + pub(crate) fn is_reduced(&self) -> bool { self.reduce_map.is_some() } - pub fn count_reduced(&self) -> usize { + /// returns number of constraints eliminated + pub(crate) fn count_reduced(&self) -> usize { self.mfull - self.mreduced } diff --git a/src/solver/implementations/default/problemdata.rs b/src/solver/implementations/default/problemdata.rs index 27edcad7..5795d295 100644 --- a/src/solver/implementations/default/problemdata.rs +++ b/src/solver/implementations/default/problemdata.rs @@ -18,14 +18,21 @@ use crate::solver::chordal::ChordalInfo; /// Standard-form solver type implementing the [`ProblemData`](crate::solver::core::traits::ProblemData) trait pub struct DefaultProblemData { - // the main KKT residuals + /// The matrix P in the quadratic objective term pub P: CscMatrix, + /// The vector q in the quadratic objective term pub q: Vec, + /// The matrix A in the constraints pub A: CscMatrix, + /// The vector b in the constraints pub b: Vec, + /// Vector of cones in the problem pub cones: Vec>, + /// Number of variables pub n: usize, + /// Number of constraints pub m: usize, + /// Equilibration data for the problem pub equilibration: DefaultEquilibrationData, // unscaled inf norms of linear terms. Set to "None" @@ -44,6 +51,7 @@ impl DefaultProblemData where T: FloatT, { + /// Create a new `DefaultProblemData` object pub fn new( P: &CscMatrix, q: &[T], @@ -111,7 +119,7 @@ where //for inf values that were not in a reduced cone //this is not considered part of the "presolve", so //can always happen regardless of user settings - let infbound = crate::solver::get_infinity().as_T(); + let infbound = crate::get_infinity().as_T(); b_new.scalarop(|x| T::min(x, infbound)); // this ensures m is the *reduced* size m diff --git a/src/solver/implementations/default/residuals.rs b/src/solver/implementations/default/residuals.rs index 568b4a23..df9e2878 100644 --- a/src/solver/implementations/default/residuals.rs +++ b/src/solver/implementations/default/residuals.rs @@ -10,29 +10,30 @@ use crate::solver::core::traits::Residuals; /// Standard-form solver type implementing the [`Residuals`](crate::solver::core::traits::Residuals) trait pub struct DefaultResiduals { // the main KKT residuals - pub rx: Vec, - pub rz: Vec, - pub rτ: T, + pub(crate) rx: Vec, + pub(crate) rz: Vec, + pub(crate) rτ: T, // partial residuals for infeasibility checks - pub rx_inf: Vec, - pub rz_inf: Vec, + pub(crate) rx_inf: Vec, + pub(crate) rz_inf: Vec, // various inner products. // NB: these are invariant w.r.t equilibration - pub dot_qx: T, - pub dot_bz: T, - pub dot_sz: T, - pub dot_xPx: T, + pub(crate) dot_qx: T, + pub(crate) dot_bz: T, + pub(crate) dot_sz: T, + pub(crate) dot_xPx: T, // the product Px by itself. Required for infeasibilty checks - pub Px: Vec, + pub(crate) Px: Vec, } impl DefaultResiduals where T: FloatT, { + /// Create a new `DefaultResiduals` object pub fn new(n: usize, m: usize) -> Self { let rx = vec![T::zero(); n]; let rz = vec![T::zero(); m]; diff --git a/src/solver/implementations/default/settings.rs b/src/solver/implementations/default/settings.rs index 0e3dbf75..d97ae503 100644 --- a/src/solver/implementations/default/settings.rs +++ b/src/solver/implementations/default/settings.rs @@ -97,7 +97,7 @@ pub struct DefaultSettings { #[builder(default = "(1e+4).as_T()")] pub equilibrate_max_scaling: T, - ///linesearch backtracking + ///line search backtracking #[builder(default = "(0.8).as_T()")] pub linesearch_backtrack_step: T, @@ -224,8 +224,8 @@ impl DefaultSettingsBuilder where T: FloatT, { + /// check that the specified direct_solve_method is valid pub fn validate(&self) -> Result<(), String> { - // check that the direct solve method is valid if let Some(ref direct_solve_method) = self.direct_solve_method { validate_direct_solve_method(direct_solve_method.as_str())?; } @@ -252,6 +252,7 @@ impl DefaultSettings where T: FloatT, { + /// Checks that the settings are valid pub fn validate(&self) -> Result<(), String> { validate_direct_solve_method(&self.direct_solve_method)?; diff --git a/src/solver/implementations/default/solution.rs b/src/solver/implementations/default/solution.rs index 7b39b44d..be032b10 100644 --- a/src/solver/implementations/default/solution.rs +++ b/src/solver/implementations/default/solution.rs @@ -9,15 +9,25 @@ use crate::{ /// Standard-form solver type implementing the [`Solution`](crate::solver::core::traits::Solution) trait #[derive(Debug)] pub struct DefaultSolution { + /// primal solution pub x: Vec, + /// dual solution (in dual cone) pub z: Vec, + /// vector of slacks (in primal cone) pub s: Vec, + /// final solver status pub status: SolverStatus, + /// primal objective value pub obj_val: T, + /// dual objective value pub obj_val_dual: T, + /// solve time in seconds pub solve_time: f64, + /// number of iterations pub iterations: u32, + /// primal residual pub r_prim: T, + /// dual residual pub r_dual: T, } @@ -25,6 +35,7 @@ impl DefaultSolution where T: FloatT, { + /// Create a new `DefaultSolution` object pub fn new(n: usize, m: usize) -> Self { let x = vec![T::zero(); n]; let z = vec![T::zero(); m]; diff --git a/src/solver/implementations/default/variables.rs b/src/solver/implementations/default/variables.rs index f8610d6d..49dda749 100644 --- a/src/solver/implementations/default/variables.rs +++ b/src/solver/implementations/default/variables.rs @@ -38,6 +38,7 @@ impl DefaultVariables where T: FloatT, { + /// Create a new `DefaultVariables` object pub fn new(n: usize, m: usize) -> Self { let x = vec![T::zero(); n]; let s = vec![T::zero(); m]; diff --git a/src/solver/mod.rs b/src/solver/mod.rs index 1259001d..221c4ed6 100644 --- a/src/solver/mod.rs +++ b/src/solver/mod.rs @@ -11,12 +11,9 @@ //! It is also possible to implement a custom solver by defining a collection //! of custom types that together implement all of the required core //! [traits] for objects in Clarabel's core solver. - -pub(crate) const _INFINITY_DEFAULT: f64 = 1e20; // internal module structure pub(crate) mod core; pub mod implementations; -pub(crate) mod utils; // chordal decomposition only if SDPs are enabled #[cfg(feature = "sdp")] @@ -26,8 +23,6 @@ pub(crate) mod chordal; //and rearrange public modules a bit to give a more //user friendly API -pub use crate::solver::utils::infbounds::*; - //allows declaration of cone constraints pub use crate::solver::core::cones::{SupportedConeT, SupportedConeT::*}; diff --git a/src/timers/timers.rs b/src/timers/timers.rs index 1177ec26..4f5d3c38 100644 --- a/src/timers/timers.rs +++ b/src/timers/timers.rs @@ -1,3 +1,4 @@ +#![allow(missing_docs)] use std::collections::HashMap; use std::ops::{Deref, DerefMut}; diff --git a/src/solver/utils/atomic.rs b/src/utils/atomic.rs similarity index 100% rename from src/solver/utils/atomic.rs rename to src/utils/atomic.rs diff --git a/src/solver/utils/infbounds.rs b/src/utils/infbounds.rs similarity index 86% rename from src/solver/utils/infbounds.rs rename to src/utils/infbounds.rs index 0e34e6ee..10c4188d 100644 --- a/src/solver/utils/infbounds.rs +++ b/src/utils/infbounds.rs @@ -1,16 +1,16 @@ -use crate::solver::utils::atomic::{AtomicF64, Ordering}; +use crate::utils::atomic::{AtomicF64, Ordering}; use lazy_static::lazy_static; /// Constant indicating that an inequality bound is to be treated as infinite. /// /// If the setting [`presolve_enable`](crate::solver::DefaultSettings::presolve_enable) /// is `true`, any such constraints are removed. Bounds for all other cones with -/// values greather than this are capped at this value. +/// values greater than this are capped at this value. /// A custom constant for this bound can be specified using [`set_infinity`]. /// /// Setting the infinity bound to a custom constant applies at module level. /// -pub const INFINITY_DEFAULT: f64 = crate::solver::_INFINITY_DEFAULT; +pub const INFINITY_DEFAULT: f64 = crate::_INFINITY_DEFAULT; lazy_static! { static ref INFINITY: AtomicF64 = AtomicF64::new(INFINITY_DEFAULT); diff --git a/src/solver/utils/mod.rs b/src/utils/mod.rs similarity index 100% rename from src/solver/utils/mod.rs rename to src/utils/mod.rs diff --git a/tests/presolve.rs b/tests/presolve.rs index f7f54236..e2d301a2 100644 --- a/tests/presolve.rs +++ b/tests/presolve.rs @@ -1,6 +1,7 @@ #![allow(non_snake_case)] use clarabel::{algebra::*, solver::*}; +use clarabel::{default_infinity, get_infinity, set_infinity}; #[allow(clippy::type_complexity)] fn presolve_test_data() -> (