Quadratic Program (QP) Example in Julia

We will show how to solve the Basic QP example problem both natively in Clarabel.jl and also by solving with Clarabel.jl within either JuMP or Convex.jl.

Clarabel.jl native interface

using Clarabel, SparseArrays

P = sparse([3. 0.;0. 2.].*2)

q = [-1., -4.]

A = sparse([1. -2.;    #<-- LHS of equality constraint
            1.  0.;    #<-- LHS of inequality constraint (upper bound)
            0.  1.;    #<-- LHS of inequality constraint (upper bound)
           -1.  0.;    #<-- LHS of inequality constraint (lower bound)
            0. -1.;    #<-- LHS of inequality constraint (lower bound)
            ])

b = [0.;        #<-- RHS of equality constraint
     ones(4)]   #<-- RHS of inequality constraints

cones =
    [Clarabel.ZeroConeT(1),           #<--- for the equality constraint
     Clarabel.NonnegativeConeT(4)]    #<--- for the inequality constraints

settings = Clarabel.Settings()

solver   = Clarabel.Solver()

Clarabel.setup!(solver, P, q, A, b, cones, settings)

result = Clarabel.solve!(solver)

-------------------------------------------------------------
           Clarabel.jl v0.5.1  -  Clever Acronym
                   (c) Paul Goulart
                University of Oxford, 2022
-------------------------------------------------------------

problem:
  variables     = 2
  constraints   = 5
  nnz(P)        = 2
  nnz(A)        = 6
  cones (total) = 2
    : Zero        = 1,  numel = 1
    : Nonnegative = 1,  numel = 4

settings:
  linear algebra: direct / qdldl, precision: Float64
  max iter = 200, time limit = Inf,  max step = 0.990
  tol_feas = 1.0e-08, tol_gap_abs = 1.0e-08, tol_gap_rel = 1.0e-08,
  static reg : on, ϵ1 = 1.0e-08, ϵ2 = 4.9e-32
  dynamic reg: on, ϵ = 1.0e-13, δ = 2.0e-07
  iter refine: on, reltol = 1.0e-13, abstol = 1.0e-12,
               max iter = 10, stop ratio = 5.0
  equilibrate: on, min_scale = 1.0e-04, max_scale = 1.0e+04
               max iter = 10

iter    pcost        dcost       gap       pres      dres      k/t        μ       step
---------------------------------------------------------------------------------------------
  0  -3.2878e-01  -1.8973e+01  1.86e+01  9.76e-17  9.02e-17  1.00e+00  1.77e+00   ------
  1  -5.7453e-01  -1.8040e+00  1.23e+00  4.81e-17  3.23e-17  6.08e-02  1.15e-01  9.44e-01
  2  -6.4203e-01  -7.5345e-01  1.11e-01  1.63e-17  3.25e-17  5.21e-03  8.97e-03  9.90e-01
  3  -6.4286e-01  -6.4541e-01  2.55e-03  1.08e-16  3.26e-17  1.19e-04  2.00e-04  9.90e-01
  4  -6.4286e-01  -6.4288e-01  2.63e-05  4.89e-17  7.99e-18  1.23e-06  2.06e-06  9.90e-01
  5  -6.4286e-01  -6.4286e-01  2.63e-07  8.82e-15  1.81e-14  1.23e-08  2.06e-08  9.90e-01
  6  -6.4286e-01  -6.4286e-01  2.63e-09  2.30e-16  4.02e-16  1.23e-10  2.06e-10  9.90e-01
---------------------------------------------------------------------------------------------
Terminated with status = solved
solve time =  942μs

>>> Clarabel - Results
Status: SOLVED
Iterations: 6
Objective: -0.6429
Solve time:  942μs

result.x

2-element Vector{Float64}:
 0.428571428198843
 0.2142857140994215

Using JuMP

We can solve the same problem a little more easily by using Clarabel.jl as the backend solver within JuMP. Here is the same problem again:

using Clarabel, JuMP

model = JuMP.Model(Clarabel.Optimizer)
set_optimizer_attribute(model, "verbose", true)

@variable(model, x[1:2])
@constraint(model, x[1] == 2x[2])
@constraint(model,  -1 .<= x .<= 1)
@objective(model, Min, 3x[1]^2 + 2x[2]^2 - x[1] - 4x[2])

optimize!(model)

-------------------------------------------------------------
           Clarabel.jl v0.5.1  -  Clever Acronym
                   (c) Paul Goulart
                University of Oxford, 2022
-------------------------------------------------------------

problem:
  variables     = 2
  constraints   = 5
  nnz(P)        = 2
  nnz(A)        = 6
  cones (total) = 2
    : Zero        = 1,  numel = 1
    : Nonnegative = 1,  numel = 4

settings:
  linear algebra: direct / qdldl, precision: Float64
  max iter = 200, time limit = Inf,  max step = 0.990
  tol_feas = 1.0e-08, tol_gap_abs = 1.0e-08, tol_gap_rel = 1.0e-08,
  static reg : on, ϵ1 = 1.0e-08, ϵ2 = 4.9e-32
  dynamic reg: on, ϵ = 1.0e-13, δ = 2.0e-07
  iter refine: on, reltol = 1.0e-13, abstol = 1.0e-12,
               max iter = 10, stop ratio = 5.0
  equilibrate: on, min_scale = 1.0e-04, max_scale = 1.0e+04
               max iter = 10

iter    pcost        dcost       gap       pres      dres      k/t        μ       step
---------------------------------------------------------------------------------------------
  0  -3.2878e-01  -1.8973e+01  1.86e+01  9.48e-17  1.50e-16  1.00e+00  1.77e+00   ------
  1  -5.7453e-01  -1.8040e+00  1.23e+00  4.89e-17  6.32e-17  6.08e-02  1.15e-01  9.44e-01
  2  -6.4203e-01  -7.5345e-01  1.11e-01  5.91e-17  3.34e-17  5.21e-03  8.97e-03  9.90e-01
  3  -6.4286e-01  -6.4541e-01  2.55e-03  1.14e-16  0.00e+00  1.19e-04  2.00e-04  9.90e-01
  4  -6.4286e-01  -6.4288e-01  2.63e-05  4.71e-17  3.20e-17  1.23e-06  2.06e-06  9.90e-01
  5  -6.4286e-01  -6.4286e-01  2.63e-07  8.82e-15  1.81e-14  1.23e-08  2.06e-08  9.90e-01
  6  -6.4286e-01  -6.4286e-01  2.63e-09  1.96e-16  3.17e-16  1.23e-10  2.06e-10  9.90e-01
---------------------------------------------------------------------------------------------
Terminated with status = solved
solve time =  925μs

Here is the solution

JuMP.value.(x)

2-element Vector{Float64}:
 0.4285714281988431
 0.21428571409942154

and the solver termination status again

JuMP.termination_status(model)

OPTIMAL::TerminationStatusCode = 1

Note that in the Clarabel.jl output that follows the call to solve! using Convex.jl, the problem posed to the solver has been converted to a second-order cone program with a linear objective. You can see this because now nnz(P) == 0 (there is no quadratic term in the objective) and the solver reports two second order cone constraints.

Although the solution will be the same, the required number of iterations and solve time are slightly higher. When solving problems with quadratic objectives in Clarabel.jl, it is generally preferable to use either the native Clarabel.jl interface or JuMP, both of which handle quadratic terms in the objective directly.

This page was generated using Literate.jl.

Quadratic Program (QP) Example in Julia

Clarabel.jl native interface

Using JuMP

Using Convex.jl