Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

add to docs and edit readme #709

Merged
merged 1 commit into from
May 18, 2021
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
47 changes: 11 additions & 36 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand Down
8 changes: 6 additions & 2 deletions docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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(
Expand Down
29 changes: 29 additions & 0 deletions docs/src/api.md
Original file line number Diff line number Diff line change
@@ -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]
```
170 changes: 166 additions & 4 deletions docs/src/index.md
Original file line number Diff line number Diff line change
@@ -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
```
2 changes: 1 addition & 1 deletion src/Models/Models.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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}()
Expand Down
6 changes: 3 additions & 3 deletions src/Solvers/steppers/combined.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand Down