| Title: | R Interface to Proximal Interior Point Quadratic Programming Solver |
|---|---|
| Description: | An embedded proximal interior point quadratic programming solver, which can solve dense and sparse quadratic programs, described in Schwan, Jiang, Kuhn, and Jones (2023) <doi:10.48550/arXiv.2304.00290>. Combining an infeasible interior point method with the proximal method of multipliers, the algorithm can handle ill-conditioned convex quadratic programming problems without the need for linear independence of the constraints. The solver is written in header only 'C++ 14' leveraging the 'Eigen' library for vectorized linear algebra. For small dense problems, vectorized instructions and cache locality can be exploited more efficiently. Allocation free problem updates and re-solves are also provided. |
| Authors: | Balasubramanian Narasimhan [aut, cre], Roland Schwan [aut, cph], Yuning Jiang [aut], Daniel Kuhn [aut], Colin N. Jones [aut] |
| Maintainer: | Balasubramanian Narasimhan <[email protected]> |
| License: | BSD_2_clause + file LICENSE |
| Version: | 0.6.2 |
| Built: | 2026-05-19 09:21:05 UTC |
| Source: | https://github.com/predict-epfl/piqp-r |
Get dimensions of a PIQP model
get_dims(model, ...)get_dims(model, ...)
model |
A |
... |
not used |
a list with named elements n, p, and m
Get settings of a PIQP model
get_settings(model, ...)get_settings(model, ...)
model |
A |
... |
not used |
a list of current settings
Convert a plain matrix or simple triplet form matrix to a Matrix::dgCMatrix (implicit) form
make_csc_matrix(x)make_csc_matrix(x)
x |
a matrix or a simple triplet form matrix |
a list of row pointer, column pointer, and values corresponding to a Matrix::dgCMatrix object
Convert a plain matrix or simple triplet form matrix to a Matrix::dsCMatrix (implicit, upper) form
make_csc_symm_matrix(m)make_csc_symm_matrix(m)
m |
a matrix or a simple triplet form matrix |
a list of row pointer, column pointer, and values corresponding to a Matrix::dsCMatrix object
PIQP Solver object
piqp( P = NULL, c = NULL, A = NULL, b = NULL, G = NULL, h_l = NULL, h_u = NULL, x_l = NULL, x_u = NULL, settings = list(), backend = c("auto", "sparse", "dense") )piqp( P = NULL, c = NULL, A = NULL, b = NULL, G = NULL, h_l = NULL, h_u = NULL, x_l = NULL, x_u = NULL, settings = list(), backend = c("auto", "sparse", "dense") )
P |
dense or sparse matrix of class dgCMatrix or coercible into such, must be positive semidefinite |
c |
numeric vector |
A |
dense or sparse matrix of class dgCMatrix or coercible into such |
b |
numeric vector |
G |
dense or sparse matrix of class dgCMatrix or coercible into such |
h_l |
numeric vector of lower inequality bounds, default |
h_u |
numeric vector of upper inequality bounds, default |
x_l |
a numeric vector of lower variable bounds, default |
x_u |
a numeric vector of upper variable bounds, default |
settings |
list with optimization parameters, empty by default; see |
backend |
which backend to use, if auto and P, A or G are sparse then sparse backend is used ( |
Allows one to solve a parametric
problem with for example warm starts between updates of the parameter, c.f. the examples.
The object returned by piqp contains several methods which can be used to either update/get details of the
problem, modify the optimization settings or attempt to solve the problem.
An S7 object of class "piqp_model" with methods
solve(), update(), get_settings(), get_dims(), update_settings()
which can be used to solve the problem with updated settings / parameters.
model = piqp(P = NULL, c = NULL, A = NULL, b = NULL, G = NULL,
h_l = NULL, h_u = NULL, x_l = NULL, x_u = NULL,
settings = piqp_settings(),
backend = c("auto", "sparse", "dense"))
solve(model)
update(model, P = NULL, c = NULL, A = NULL, b = NULL, G = NULL,
h_l = NULL, h_u = NULL, x_l = NULL, x_u = NULL)
get_settings(model)
get_dims(model)
update_settings(model, new_settings = piqp_settings())
## example, adapted from PIQP documentation library(piqp) library(Matrix) P <- Matrix(c(6., 0., 0., 4.), 2, 2, sparse = TRUE) c <- c(-1., -4.) A <- Matrix(c(1., -2.), 1, 2, sparse = TRUE) b <- c(1.) G <- Matrix(c(1., 2., -1., 0.), 2, 2, sparse = TRUE) h_u <- c(0.2, -1.) x_l <- c(-1., -Inf) x_u <- c(1., Inf) settings <- list(verbose = TRUE) model <- piqp(P, c, A, b, G, h_u = h_u, x_l = x_l, x_u = x_u, settings = settings) # Solve res <- solve(model) res$x # Define new data A_new <- Matrix(c(1., -3.), 1, 2, sparse = TRUE) h_u_new <- c(2., 1.) # Update model and solve again update(model, A = A_new, h_u = h_u_new) res <- solve(model) res$x## example, adapted from PIQP documentation library(piqp) library(Matrix) P <- Matrix(c(6., 0., 0., 4.), 2, 2, sparse = TRUE) c <- c(-1., -4.) A <- Matrix(c(1., -2.), 1, 2, sparse = TRUE) b <- c(1.) G <- Matrix(c(1., 2., -1., 0.), 2, 2, sparse = TRUE) h_u <- c(0.2, -1.) x_l <- c(-1., -Inf) x_u <- c(1., Inf) settings <- list(verbose = TRUE) model <- piqp(P, c, A, b, G, h_u = h_u, x_l = x_l, x_u = x_u, settings = settings) # Solve res <- solve(model) res$x # Define new data A_new <- Matrix(c(1., -3.), 1, 2, sparse = TRUE) h_u_new <- c(2., 1.) # Update model and solve again update(model, A = A_new, h_u = h_u_new) res <- solve(model) res$x
An S7 class wrapping the PIQP C++ Solver. Users will never
need to directly create instances of this class and should use the
more user-friendly functions piqp() and solve_piqp().
piqp_model(solver_ptr = NULL, dims = list(), dense_backend = logical(0))piqp_model(solver_ptr = NULL, dims = list(), dense_backend = logical(0))
solver_ptr |
external pointer to the C++ solver object |
dims |
a named list with elements |
dense_backend |
logical flag indicating if the dense solver is used |
Settings parameters with default values and types in parenthesis
piqp_settings( rho_init = 1e-06, delta_init = 1e-04, eps_abs = 1e-08, eps_rel = 1e-09, check_duality_gap = TRUE, eps_duality_gap_abs = 1e-08, eps_duality_gap_rel = 1e-09, infeasibility_threshold = 0.9, reg_lower_limit = 1e-10, reg_finetune_lower_limit = 1e-13, reg_finetune_primal_update_threshold = 7L, reg_finetune_dual_update_threshold = 7L, max_iter = 250L, max_factor_retires = 10L, preconditioner_scale_cost = FALSE, preconditioner_reuse_on_update = FALSE, preconditioner_iter = 10L, tau = 0.99, iterative_refinement_always_enabled = FALSE, iterative_refinement_eps_abs = 1e-12, iterative_refinement_eps_rel = 1e-12, iterative_refinement_max_iter = 10L, iterative_refinement_min_improvement_rate = 5, iterative_refinement_static_regularization_eps = 1e-08, iterative_refinement_static_regularization_rel = .Machine$double.eps^2, verbose = FALSE, compute_timings = FALSE )piqp_settings( rho_init = 1e-06, delta_init = 1e-04, eps_abs = 1e-08, eps_rel = 1e-09, check_duality_gap = TRUE, eps_duality_gap_abs = 1e-08, eps_duality_gap_rel = 1e-09, infeasibility_threshold = 0.9, reg_lower_limit = 1e-10, reg_finetune_lower_limit = 1e-13, reg_finetune_primal_update_threshold = 7L, reg_finetune_dual_update_threshold = 7L, max_iter = 250L, max_factor_retires = 10L, preconditioner_scale_cost = FALSE, preconditioner_reuse_on_update = FALSE, preconditioner_iter = 10L, tau = 0.99, iterative_refinement_always_enabled = FALSE, iterative_refinement_eps_abs = 1e-12, iterative_refinement_eps_rel = 1e-12, iterative_refinement_max_iter = 10L, iterative_refinement_min_improvement_rate = 5, iterative_refinement_static_regularization_eps = 1e-08, iterative_refinement_static_regularization_rel = .Machine$double.eps^2, verbose = FALSE, compute_timings = FALSE )
rho_init |
Initial value for the primal proximal penalty parameter rho (default = 1e-6) |
delta_init |
Initial value for the augmented lagrangian penalty parameter delta (default = 1e-4) |
eps_abs |
Absolute tolerance (default = 1e-8) |
eps_rel |
Relative tolerance (default = 1e-9) |
check_duality_gap |
Check terminal criterion on duality gap (default = TRUE) |
eps_duality_gap_abs |
Absolute tolerance on duality gap (default = 1e-8) |
eps_duality_gap_rel |
Relative tolerance on duality gap (default = 1e-9) |
infeasibility_threshold |
Threshold for infeasibility detection (default = 0.9) |
reg_lower_limit |
Lower limit for regularization (default = 1e-10) |
reg_finetune_lower_limit |
Fine tune lower limit regularization (default = 1e-13) |
reg_finetune_primal_update_threshold |
Threshold of number of no primal updates to transition to fine tune mode (default = 7) |
reg_finetune_dual_update_threshold |
Threshold of number of no dual updates to transition to fine tune mode (default = 7) |
max_iter |
Maximum number of iterations (default = 250) |
max_factor_retires |
Maximum number of factorization retires before failure (default = 10) |
preconditioner_scale_cost |
Scale cost in Ruiz preconditioner (default = FALSE) |
preconditioner_reuse_on_update |
Reuse preconditioner on problem update (default = FALSE) |
preconditioner_iter |
Maximum of preconditioner iterations (default = 10) |
tau |
Maximum interior point step length (default = 0.99) |
iterative_refinement_always_enabled |
Always run iterative refinement and not only on factorization failure (default = FALSE) |
iterative_refinement_eps_abs |
Iterative refinement absolute tolerance (default = 1e-12) |
iterative_refinement_eps_rel |
Iterative refinement relative tolerance (default = 1e-12) |
iterative_refinement_max_iter |
Maximum number of iterations for iterative refinement (default = 10) |
iterative_refinement_min_improvement_rate |
Minimum improvement rate for iterative refinement (default = 5.0) |
iterative_refinement_static_regularization_eps |
Static regularization for KKT system for iterative refinement (default = 1e-8) |
iterative_refinement_static_regularization_rel |
Static regularization w.r.t. the maximum abs diagonal term of KKT system. (default = .Machine$double.eps^2) |
verbose |
Verbose printing (default = FALSE) |
compute_timings |
Measure timing information internally (default = FALSE) |
a list containing the settings parameters.
Solves
s.t.
for real matrices P (nxn, positive semidefinite), A (pxn) with p number of equality constraints, and G (mxn) with m number of inequality constraints
solve_piqp( P = NULL, c = NULL, A = NULL, b = NULL, G = NULL, h_l = NULL, h_u = NULL, x_l = NULL, x_u = NULL, settings = list(), backend = c("auto", "sparse", "dense") )solve_piqp( P = NULL, c = NULL, A = NULL, b = NULL, G = NULL, h_l = NULL, h_u = NULL, x_l = NULL, x_u = NULL, settings = list(), backend = c("auto", "sparse", "dense") )
P |
dense or sparse matrix of class dgCMatrix or coercible into such, must be positive semidefinite |
c |
numeric vector |
A |
dense or sparse matrix of class dgCMatrix or coercible into such |
b |
numeric vector |
G |
dense or sparse matrix of class dgCMatrix or coercible into such |
h_l |
numeric vector of lower inequality bounds, default |
h_u |
numeric vector of upper inequality bounds, default |
x_l |
a numeric vector of lower variable bounds, default |
x_u |
a numeric vector of upper variable bounds, default |
settings |
list with optimization parameters, empty by default; see |
backend |
which backend to use, if auto and P, A or G are sparse then sparse backend is used ( |
A list with elements solution elements
Schwan, R., Jiang, Y., Kuhn, D., Jones, C.N. (2023). “PIQP: A Proximal Interior-Point Quadratic Programming Solver.” https://doi.org/10.48550/arXiv.2304.00290
piqp(), piqp_settings() and the underlying PIQP documentation: https://predict-epfl.github.io/piqp/
## example, adapted from PIQP documentation library(piqp) library(Matrix) P <- Matrix(c(6., 0., 0., 4.), 2, 2, sparse = TRUE) c <- c(-1., -4.) A <- Matrix(c(1., -2.), 1, 2, sparse = TRUE) b <- c(1.) G <- Matrix(c(1., 2., -1., 0.), 2, 2, sparse = TRUE) h_u <- c(0.2, -1.) x_l <- c(-1., -Inf) x_u <- c(1., Inf) settings <- list(verbose = TRUE) # Solve with PIQP res <- solve_piqp(P, c, A, b, G, h_u = h_u, x_l = x_l, x_u = x_u, settings = settings) res$x## example, adapted from PIQP documentation library(piqp) library(Matrix) P <- Matrix(c(6., 0., 0., 4.), 2, 2, sparse = TRUE) c <- c(-1., -4.) A <- Matrix(c(1., -2.), 1, 2, sparse = TRUE) b <- c(1.) G <- Matrix(c(1., 2., -1., 0.), 2, 2, sparse = TRUE) h_u <- c(0.2, -1.) x_l <- c(-1., -Inf) x_u <- c(1., Inf) settings <- list(verbose = TRUE) # Solve with PIQP res <- solve_piqp(P, c, A, b, G, h_u = h_u, x_l = x_l, x_u = x_u, settings = settings) res$x
Return the solver status description string
status_description(code)status_description(code)
code |
a valid solver return code |
a status description string
status_description(1) ## for solved problem status_description(-1) ## for max iterations limit reachedstatus_description(1) ## for solved problem status_description(-1) ## for max iterations limit reached
Update settings of a PIQP model
update_settings(model, ...)update_settings(model, ...)
model |
A |
... |
not used |
invisible NULL