Modern Fortran Edition of the quadprog solver
This is an updated version of the quadprog solver initially written by Berwin A. Turlach in Fortran 77. It can be used to solve strictly convex quadratic programs of the form
using an active set method. It is most efficient for small to moderate sized QP described using dense matrices.
Updates to the original code include:
- Sources have been translated from FORTRAN 77 fixed-form to Fortran 90 free-form.
- All obsolescent features (
goto, etc) have been removed. It is now 100% standard compliant (Fortran 2018). - It makes use of derived-type and easy to use interfaces. The
qp_problemclass is used to defined the quadratic program andsolveto compute its solution. - Calls to
blasfunctions and subroutines now replace some hand-crafted implementations for improved performances. - Calls to
lapacksubroutines now replace the original functionalities provided bylinpack. - Utility solvers for non-negative least-squares (
nnls) and bounded-variables least-squares (bvls) are provided.
The library can be built with the Fortran Package Manager fpm using the provided fpm.toml file like so:
fpm build --releaseOnly double precision (real64) is currently supported.
To use QuadProg within your fpm project, add the following to your fpm.toml file:
[dependencies]
QuadProg = { git="https://github.com/loiseaujc/QuadProg.git"}The library requires some blas and lapack routines which are not included. You thus need to have it available on your system, otherwise it'll automatically pull them from the fortran standard library stdlib.
The following program solves the QP
program example
use quadprog
implicit none
integer, parameter :: dp = selected_real_kind(15, 307)
! Size of the problem.
integer, parameter :: n = 3
! Quadratic cost.
real(dp) :: P(n, n), q(n)
! Inequality constraints.
real(dp) :: C(n, n), d(n)
! Convenience types.
type(qp_problem) :: prob
type(OptimizeResult) :: solution
! Miscellaneous
integer :: i
!> Setup the quadratic function..
P = 0.0_dp ; q = [0.0_dp, 5.0_dp, 0.0_dp]
do i = 1, n
P(i, i) = 1.0_dp
enddo
!> Setup the inequality constraints.
C(:, 1) = [-4.0_dp, 2.0_dp, 0.0_dp]
C(:, 2) = [-3.0_dp, 1.0_dp, -2.0_dp]
C(:, 3) = [0.0_dp, 0.0_dp, 1.0_dp]
d = [-8.0_dp, 2.0_dp, 0.0_dp]
!> Solve the inequality constrained QP.
prob = qp_problem(P, q, C=C, d=d)
solution = solve(prob)
if (solution%success) then
print *, "x =", solution%x ! Solution of the QP.
print *, "y =", solution%y ! Lagrange multipliers.
print *, "obj =", solution%obj ! Objective function.
endif
end program- The original source code was released under the GNU General Public Licence version 2 (GPL-2).
- Development continues on Github.
quadprog, an equivalent project written in JavaScript.quadprog, an equivalent project written in Rust.quadprog, some R bindings to the originalquadprogsolver.GoldfarbIdnaniSolver.jl, a port of the originalquadprogcode in Julia.
- D. Goldfarb and A. Idnani (1983). A numerically stable dual method for solving strictly convex quadratic programs. Mathematical Programming, 27, 1-33.
