From f6754d948a9e5e172e460f4cc7ed1db838631cfb Mon Sep 17 00:00:00 2001 From: Chris Coey Date: Tue, 18 May 2021 23:41:44 +1000 Subject: [PATCH] add to docs and readme --- README.md | 47 ++------- docs/make.jl | 8 +- docs/src/api.md | 29 ++++++ docs/src/index.md | 170 ++++++++++++++++++++++++++++++- src/Models/Models.jl | 2 +- src/Solvers/steppers/combined.jl | 6 +- 6 files changed, 216 insertions(+), 46 deletions(-) create mode 100644 docs/src/api.md diff --git a/README.md b/README.md index f179e6d8b..8e3934fc0 100644 --- a/README.md +++ b/README.md @@ -5,60 +5,35 @@ [![](https://img.shields.io/badge/docs-stable-blue.svg)](https://chriscoey.github.io/Hypatia.jl/stable) [![](https://img.shields.io/badge/docs-dev-blue.svg)](https://chriscoey.github.io/Hypatia.jl/dev) -Hypatia is a highly-customizable open source interior point solver for generic conic optimization problems, written in Julia. -For more information, please see our [working paper](https://arxiv.org/abs/2005.01136) and our [cones reference](https://github.com/chriscoey/Hypatia.jl/wiki/files/coneref.pdf). +# Hypatia.jl -primal (over x,s): -``` - min c'x : duals - b - Ax == 0 (y) - h - Gx == s in K (z) -``` -dual (over z,y): -``` - max -b'y - h'z : duals - c + A'y + G'z == 0 (x) - z in K* (s) -``` -where K is a convex cone defined as a Cartesian product of recognized proper -cones, and K* is its dual cone. -An objective offset can be provided as the keyword arg `obj_offset` (default 0). - -The primal-dual optimality conditions are: -``` - b - Ax == 0 - h - Gx == s - c + A'y + G'z == 0 - s'z == 0 - s in K - z in K* -``` +Hypatia is a highly-customizable interior point solver for generic conic optimization problems, written in Julia. +It is licensed under the MIT License; see [LICENSE](https://github.com/chriscoey/Hypatia.jl/blob/master/LICENSE.md). +To learn how to model and solve conic problems with Hypatia, see the many applied examples in the [examples folder](https://github.com/chriscoey/Hypatia.jl/tree/master/examples). +For more information about conic optimization, Hypatia's algorithms, and proper cones, please see our [working paper](https://arxiv.org/abs/2005.01136) and our [cones reference](https://github.com/chriscoey/Hypatia.jl/wiki/files/coneref.pdf). Hypatia is an experimental solver and a work in progress, and may not run with older releases of Julia. If you have trouble using Hypatia or wish to make improvements, please submit an issue or contribute a PR. Default options/parameters are not well-tuned, so we encourage you to experiment with these. -To learn how to model using exotic cones in Hypatia, look through the examples folder. -Our examples are set up using either [JuMP](https://github.com/jump-dev/JuMP.jl) or Hypatia's native interface. -Modeling with JuMP is generally more user-friendly, though it may make sense to try the more-expressive native interface for large dense or structured models. - -Here is a simple example (from D-optimal experiment design) that sets up a JuMP model and calls Hypatia: +Here is a simple example from D-optimal experiment design that sets up a JuMP model and calls Hypatia: ```julia using LinearAlgebra using JuMP using Hypatia -# setup model -V = rand(2, 3) +# setup JuMP model model = Model(Hypatia.Optimizer) @variable(model, x[1:3] >= 0) @constraint(model, sum(x) == 5) @variable(model, hypo) @objective(model, Max, hypo) +V = rand(2, 3) Q = V * diagm(x) * V' -@constraint(model, vcat(hypo, [Q[i, j] for i in 1:2 for j in 1:i]...) in MOI.RootDetConeTriangle(2)) +aff = vcat(hypo, [Q[i, j] for i in 1:2 for j in 1:i]...) +@constraint(model, aff in MOI.RootDetConeTriangle(2)) -# solve +# solve and query solution optimize!(model) termination_status(model) objective_value(model) diff --git a/docs/make.jl b/docs/make.jl index 68f29c1d9..c52c9f536 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -4,8 +4,12 @@ using Hypatia makedocs( sitename = "Hypatia", - format = Documenter.HTML(), - modules = [Hypatia] + format = Documenter.HTML(prettyurls = get(ENV, "CI", nothing) == "true"), + modules = [Hypatia], + pages = [ + "Home" => "index.md", + "API" => "api.md" + ], ) deploydocs( diff --git a/docs/src/api.md b/docs/src/api.md new file mode 100644 index 000000000..f1d7aee30 --- /dev/null +++ b/docs/src/api.md @@ -0,0 +1,29 @@ +# Hypatia API + +```@autodocs +Modules = [Hypatia] +``` + +## Cones module + +```@autodocs +Modules = [Hypatia.Cones] +``` + +## Models module + +```@autodocs +Modules = [Hypatia.Models] +``` + +## Solvers module + +```@autodocs +Modules = [Hypatia.Solvers] +``` + +## PolyUtils module + +```@autodocs +Modules = [Hypatia.PolyUtils] +``` diff --git a/docs/src/index.md b/docs/src/index.md index ef6368dbd..ab3b9f8ee 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -1,7 +1,169 @@ -# Hypatia.jl +# Hypatia -Documentation for Hypatia.jl +Hypatia is a highly-customizable interior point solver for generic conic optimization problems, written in Julia. +It is licensed under the MIT License; see [LICENSE](https://github.com/chriscoey/Hypatia.jl/blob/master/LICENSE.md). -```@autodocs -Modules = [Hypatia, Hypatia.Cones, Hypatia.Models, Hypatia.Solvers, Hypatia.PolyUtils] +To learn how to model and solve conic problems with Hypatia, see the many applied examples in the [examples folder](https://github.com/chriscoey/Hypatia.jl/tree/master/examples). +For more information about conic optimization, Hypatia's algorithms, and proper cones, please see our [working paper](https://arxiv.org/abs/2005.01136) and our [cones reference](https://github.com/chriscoey/Hypatia.jl/wiki/files/coneref.pdf). +Hypatia is an experimental solver and a work in progress, and may not run with older releases of Julia. +If you have trouble using Hypatia or wish to make improvements, please submit an issue or contribute a PR. +Default options/parameters are not well-tuned, so we encourage you to experiment with these. + +## Conic form + +The primal conic form over variable ``x \in \mathbb{R}^n`` is: +```math +\begin{aligned} +\min \quad c'x &: +\\ +b - Ax &= 0 +\\ +h - Gx &\in \mathcal{K} +\end{aligned} +``` +where ``\mathcal{K}`` is a proper cone. + +The corresponding conic dual form over variables ``y \in \mathbb{R}^p`` and ``z \in \mathbb{R}^q`` is: +```math +\begin{aligned} +\max \quad -b'y - h'z &: +\\ +c + A'y + G'z &= 0 +\\ +z &\in \mathcal{K}^* +\end{aligned} +``` +where ``\mathcal{K}^*`` is the dual cone of ``\mathcal{K}``. + +## Generic cone interface + +The cone interface allows specifying proper cones. +Hypatia predefines many proper cones that are practically useful; see the [Cones folder](https://github.com/chriscoey/Hypatia.jl/tree/master/src/Cones). +Defining a proper cone requires implementing several key oracles: +- an initial interior point inside the cone, +- a feasibility test, which checks whether a given point is in the interior of the cone, +- gradient and Hessian evaluations for a logarithmically homogeneous self-concordant barrier function for the cone. +Additional optional oracle can be specified to improve speed and numerical performance. +Defining a new cone automatically defines its dual cone (through the `use_dual` option) also. + +## Model interface + +A model specifies the data in the primal conic form above: +- ``\mathcal{K}`` is a vector of Hypatia cones, +- ``c \in \mathbb{R}^n``, ``b \in \mathbb{R}^p``, ``h \in \mathbb{R}^q`` are vectors, +- ``A \in \mathbb{R}^{p \times n}`` and ``G \in \mathbb{R}^{p \times n}`` are linear operators. +An objective offset can be specified using the keyword arg `obj_offset` (the default is 0). + +## Solver interface + +Hypatia's native solver interface provides low-level functions for solving models and querying solve information and conic certificates. +Below is a simple example of a spectral norm optimization problem: +```julia +using LinearAlgebra +import Hypatia +import Hypatia.Cones +import Hypatia.Solvers + +T = BigFloat +(Xn, Xm) = (3, 4) +dim = Xn * Xm +c = vcat(one(T), zeros(T, dim)) +A = hcat(zeros(T, dim, 1), Matrix{T}(I, dim, dim)) +b = rand(T, dim) +G = -one(T) * I +h = vcat(zero(T), rand(T, dim)) +cones = [Cones.EpiNormSpectral{T, T}(Xn, Xm),] +model = Hypatia.Models.Model{T}(c, A, b, G, h, cones); +``` +Now we call optimize and query the solution: +```julia +julia> solver = Solvers.Solver{T}(verbose = true); + +julia> Solvers.load(solver, model); + +julia> Solvers.solve(solver); + + iter p_obj d_obj | abs_gap x_feas z_feas | tau kap mu | dir_res step alpha + 0 2.0000e+00 0.0000e+00 | 4.00e+00 5.00e-01 6.14e-01 | 1.00e+00 1.00e+00 1.00e+00 | + 1 2.5147e+00 2.2824e+00 | 1.07e+00 2.24e-01 2.74e-01 | 6.71e-01 7.44e-01 3.14e-01 | 3.45e-77 co-a 7.00e-01 + 2 3.0958e+00 3.0966e+00 | 3.39e-01 7.40e-02 9.08e-02 | 6.08e-01 2.70e-01 1.01e-01 | 1.73e-77 co-a 7.00e-01 +... + 33 3.2962e+00 3.2962e+00 | 1.33e-30 1.86e-30 2.29e-30 | 2.21e-01 1.88e-30 3.50e-31 | 4.85e-50 co-a 5.00e-01 + 34 3.2962e+00 3.2962e+00 | 2.29e-31 2.77e-31 3.40e-31 | 2.23e-01 1.56e-31 5.28e-32 | 2.70e-48 co-a 8.50e-01 + 35 3.2962e+00 3.2962e+00 | 7.32e-32 1.21e-31 1.49e-31 | 2.04e-01 1.15e-31 1.93e-32 | 2.52e-49 co-a 6.00e-01 +optimal solution found; terminating + +status is Optimal after 35 iterations and 1.417 seconds + +julia> Solvers.get_status(solver) +Optimal::Status = 3 + +julia> Solvers.get_primal_obj(solver) +3.296219213377718379486912616497183695150874915748434424139285907044225666610375 + +julia> Solvers.get_dual_obj(solver) +3.29621921337771837948691261649702242580041815859477664034990154727149325280482 + +julia> Solvers.get_x(solver) +13-element Vector{BigFloat}: + 3.296219213377718379486912616497183695150874915748434424139285907044225666610375 + 0.2014951884389319019635492500560556305705409904059068283971229702024509946542355 + 0.2974304558864173403751380894865665103907350340298129119333915978313596204101965 +... + 0.5076818262444516526146124557277599761994226912376378984924714206880970861448005 + 0.4719189060586783692091058031965586401783925037880741775651825386863356120491953 + 0.6377366379371803537625028513405673487028387050401566067448589228065218202592247 +``` + +## MathOptInterface/JuMP + +[JuMP](https://github.com/jump-dev/JuMP.jl) is generally more user-friendly than Hypatia's native interface, though it may make sense to try the more-expressive native interface for large dense or structured models. +Below is a simple example from D-optimal experiment design: + +```julia +using LinearAlgebra +using JuMP +using Hypatia + +opt = Hypatia.Optimizer(verbose = false) +model = Model(() -> opt) +@variable(model, x[1:3] >= 0) +@constraint(model, sum(x) == 5) +@variable(model, hypo) +@objective(model, Max, hypo) +V = rand(2, 3) +Q = V * diagm(x) * V' +aff = vcat(hypo, [Q[i, j] for i in 1:2 for j in 1:i]...) +@constraint(model, aff in MOI.RootDetConeTriangle(2)) +``` +The model is now: +```julia +julia> model +A JuMP Model +Maximization problem with: +Variables: 4 +Objective function type: VariableRef +`AffExpr`-in-`MathOptInterface.EqualTo{Float64}`: 1 constraint +`Vector{AffExpr}`-in-`MathOptInterface.RootDetConeTriangle`: 1 constraint +`VariableRef`-in-`MathOptInterface.GreaterThan{Float64}`: 3 constraints +Model mode: AUTOMATIC +CachingOptimizer state: EMPTY_OPTIMIZER +Solver name: Hypatia +Names registered in the model: hypo, x +``` +Now we call optimize and query the solution: +```julia +julia> optimize!(model) + +julia> termination_status(model) +OPTIMAL::TerminationStatusCode = 1 + +julia> objective_value(model) +1.4303650845824805 + +julia> value.(x) +3-element Vector{Float64}: + 2.499999987496876 + 2.499999992300735 + 2.0202389761081463e-8 ``` diff --git a/src/Models/Models.jl b/src/Models/Models.jl index deea31edc..11fc0b298 100644 --- a/src/Models/Models.jl +++ b/src/Models/Models.jl @@ -31,7 +31,7 @@ mutable struct Model{T <: Real} b::Vector{T}, G, h::Vector{T}, - cones::Vector{Cones.Cone{T}}; + cones::Vector{<:Cones.Cone{T}}; obj_offset::T = zero(T), ) where {T <: Real} model = new{T}() diff --git a/src/Solvers/steppers/combined.jl b/src/Solvers/steppers/combined.jl index f6daac0b4..00b04f2af 100644 --- a/src/Solvers/steppers/combined.jl +++ b/src/Solvers/steppers/combined.jl @@ -83,19 +83,19 @@ function step(stepper::CombinedStepper{T}, solver::Solver{T}) where {T <: Real} if iszero(alpha) # recover - println("trying combined without adjustment") + solver.verbose && println("trying combined without adjustment") stepper.unadj_only = true solver.time_search += @elapsed alpha = search_alpha(point, model, stepper) if iszero(alpha) - println("trying centering with adjustment") + solver.verbose && println("trying centering with adjustment") stepper.cent_only = true stepper.unadj_only = false solver.time_search += @elapsed alpha = search_alpha(point, model, stepper) if iszero(alpha) - println("trying centering without adjustment") + solver.verbose && println("trying centering without adjustment") stepper.unadj_only = true solver.time_search += @elapsed alpha = search_alpha(point, model, stepper)