The source files for all examples can be found in /examples.
Constrained Least Squares
This example setsup a constrained least squares problem using the quadratic program interface of GeNIOS.jl. It is from the OSQP docs.
Specifically, we want to solve the problem
\[\begin{array}{ll} \text{minimize} & \|Ax - b\|_2^2 \\ \text{subject to} & 0 \leq x \leq 1. \end{array}\]
using GeNIOS
using Random, LinearAlgebra, SparseArrays
using JuMPGenerating the problem data
Random.seed!(1)
m, n = 30, 20
Ad = sprandn(m, n, 0.7)
b = randn(m);Using JuMP
We can model this problem using JuMP and then pass it to GeNIOS.
model = JuMP.Model(GeNIOS.Optimizer)
@variable(model, x[1:n])
@objective(model, Min, sum((Ad*x - b).^2))
@constraint(model, 0 .<= x)
@constraint(model, x .<= 1)
optimize!(model)
# Show optimal objective value
println("Optimal value: $(round(objective_value(model), digits=4))")Starting setup...
Setup in 0.000s
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Iteration Obj Val r_primal r_dual ρ Time
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
0 0.000e+00 Inf Inf 1.000e+00 0.000
1 0.000e+00 0.000e+00 3.992e+01 1.000e+00 0.000
20 -1.486e+01 3.163e-01 3.174e-02 1.000e+00 0.032
40 -1.248e+01 1.417e-01 9.042e-03 1.000e+00 0.033
60 -1.116e+01 6.001e-02 1.152e-02 2.000e+00 0.034
80 -1.061e+01 2.496e-02 3.778e-03 2.000e+00 0.034
100 -1.043e+01 1.194e-02 1.446e-03 2.000e+00 0.035
120 -1.036e+01 6.127e-03 6.277e-04 2.000e+00 0.036
140 -1.033e+01 3.248e-03 3.045e-04 2.000e+00 0.037
160 -1.031e+01 1.011e-03 3.681e-04 4.000e+00 0.037
164 -1.031e+01 7.345e-04 2.664e-04 4.000e+00 0.038
SOLVED in 0.038s, 164 iterations
Total time: 0.038s
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Optimal value: 15.3832QPSolver interface
To make the formulations match, we will introdudce a new variable $y = Ax - b$. (GeNIOS does not support constants in the objective function.) The problem becomes
\[\begin{array}{ll} \text{minimize} & y^Ty \\ \text{subject to} & Ax - y = b \\ & 0 \leq x \leq 1. \end{array}\]
In OSQP form, this problem is
\[\begin{array}{ll} \text{minimize} & y^Ty \\ \text{subject to} & Ax - y = z_1 \\ & x = z_2 \\ & b \leq z_1 \leq b \\ & 0 \leq z_2 \leq 1 \end{array}\]
or more explicitly
\[\begin{array}{ll} \text{minimize} & \begin{bmatrix} x \\ y \end{bmatrix}^T \begin{bmatrix} 0 & \\ & I \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix} \\ \text{subject to} & \begin{bmatrix} A & -I \\ I & 0 \end{bmatrix} \begin{bmatrix} x \\ y \end{bmatrix} = z & \begin{bmatrix} b \\ 0 \end{bmatrix} \leq z \leq \begin{bmatrix} b \\ 1 \end{bmatrix}. \end{array}\]
P = blockdiag(spzeros(n, n), 2*sparse(I, m, m))
q = spzeros(n + m)
M = [
Ad -sparse(I, m, m);
sparse(I, n, n) spzeros(n, m)
]
l = [b; spzeros(n)]
u = [b; ones(n)];Building and solving the problem
Now we setup the solver and solve the problem
solver = GeNIOS.QPSolver(P, q, M, l, u)
res = solve!(solver)--- GeNIOSResult ---
Status: OPTIMAL
Obj. val: 15.38
num iters: 140
setup time: 0.013s
solve time: 0.22s
Examining the solution
We can check the solution and its optimality
xstar = solver.zk[m+1:end]
ls_residual = 1/√(m) * norm(Ad*xstar - b)
feas = all(xstar .>= 0) && all(xstar .<= 1)
println("Is feasible? $feas")
println("Least Squres RMSE = $(round(ls_residual, digits=8))")
println("Primal residual: $(round(solver.rp_norm, digits=8))")
println("Dual residual: $(round(solver.rd_norm, digits=8))")
println("Optimal value: $(round(solver.obj_val, digits=4))")Is feasible? true
Least Squres RMSE = 0.71620903
Primal residual: 0.00114969
Dual residual: 0.00031261
Optimal value: 15.3831This page was generated using Literate.jl.