quadprog_ospq_comp_1.R

library(quadprog)
library(osqp)

# This gist compares the quadprog solver against OSPQ on a very simple problem.

##
##  The following example comes directly from `?? quadprog` documentation`
##
## Assume we want to minimize: -(0 5 0) %*% b + 1/2 b^T b
## under the constraints:      A^T b >= b0
## with b0 = (-8,2,0)^T
## and      (-4  2  0) 
##      A = (-3  1 -2)
##          ( 0  0  1)
## we can use solve.QP as follows:
##
Dmat       <- matrix(0,3,3)
diag(Dmat) <- 1
dvec       <- c(0,5,0)
Amat       <- matrix(c(-4,-3,0,2,1,0,0,-2,1),3,3)
bvec       <- c(-8,2,0)

# quadprog solution
quadprog_solution <- solve.QP(Dmat,dvec,Amat,bvec=bvec)

#
# We now compare this solution against OSPQ
#
# IMPORTANT! The default eps_abs, eps_rel tolerance settings for OSPQ
# do not give very high precision. So we change the defaults to get a 
# solution very close to the quadprog solution
settings <- osqpSettings(verbose = FALSE, eps_abs=1e-8, eps_rel = 1e-8)
osqp_solution <- solve_osqp(Dmat, -dvec, t(Amat), l=bvec, pars=settings)

# Compute the error
err <- sqrt(sum((quadprog_solution$solution-osqp_solution$x)**2))
cat(sprintf("quadprog vs. OSQP difference: %.4ef", err))