Cones module
See the Cone interface and predefined cones section.
Hypatia.Cones
— ModuleProper cone definitions, oracles, and utilities.
Hypatia.Cones.Cone
— Typeabstract type Cone{T<:Real}
A proper cone.
Array utilities
Hypatia.Cones.vec_length
— Functionvec_length(_, len)
The dimension of the real vectorization of a real or complex vector of length len::Int
.
Hypatia.Cones.vec_copyto!
— Functionvec_copyto!(v1, v2)
Copy a vector in-place.
vec_copyto!(rvec, cvec)
Copy a complex vector to a real vector in-place.
vec_copyto!(cvec, rvec)
Copy a real vector to a complex vector in-place.
Hypatia.Cones.svec_length
— Functionsvec_length(side)
The dimension of the vectorized triangle of a real symmetric matrix with side dimension side::Int
.
svec_length(_, side)
The dimension of the real vectorized triangle of a real symmetric or complex Hermitian matrix with side dimension side::Int
.
Hypatia.Cones.svec_side
— Functionsvec_side(len)
The side dimension of a real symmetric matrix with vectorized triangle length len::Int
.
svec_side(_, len)
The side dimension of a real symmetric or complex Hermitian matrix with real vectorized triangle length len::Int
.
Hypatia.Cones.svec_idx
— Functionsvec_idx(row, col)
The index in the vectorized triangle of a symmetric matrix for element (row::Int
, col::Int
).
Hypatia.Cones.block_idxs
— Functionblock_idxs(incr, block)
The indices corresponding to block block::Int
in a vector of blocks with equal length incr::Int
.
Hypatia.Cones.scale_svec!
— Functionscale_svec!(arr, scal)
Rescale the elements corresponding to off-diagonals in arr::AbstractVecOrMat
, with scaling scal::Real
, for real symmetric matrices.
Hypatia.Cones.scale_svec_incr!
— Functionscale_svec_incr!(arr, scal, incr)
Rescale the elements corresponding to off-diagonals in arr::AbstractVecOrMat
, with scaling scal::Real
and block increment incr::Int
, for real symmetric matrices.
Hypatia.Cones.smat_to_svec!
— Functionsmat_to_svec!(vec, mat, rt2)
Copy a real symmetric matrix upper triangle to a svec-scaled vector in-place.
smat_to_svec!(vec, mat, rt2)
Copy a complex Hermitian matrix upper triangle to a svec-scaled real vector in-place.
Hypatia.Cones.svec_to_smat!
— Functionsvec_to_smat!(mat, vec, rt2)
Copy a svec-scaled vector to a real symmetric matrix upper triangle in-place.
svec_to_smat!(mat, vec, rt2)
Copy a svec-scaled real vector to a complex Hermitian matrix upper triangle in-place.
Cone oracles
Hypatia.Cones.dimension
— Functiondimension(cone)
The real vector dimension of the cone.
Hypatia.Cones.get_nu
— Functionget_nu(cone)
The barrier parameter $\nu$ of the cone.
Hypatia.Cones.set_initial_point!
— Functionset_initial_point!(arr, cone)
Set the array equal to the initial interior point for the cone.
Hypatia.Cones.is_feas
— Functionis_feas(cone)
Returns true if and only if the currently-loaded primal point is strictly feasible for the cone.
Hypatia.Cones.is_dual_feas
— Functionis_dual_feas(cone)
Returns false only if the currently-loaded dual point is outside the interior of the cone's dual cone.
Hypatia.Cones.grad
— Functiongrad(cone)
The gradient of the cone's barrier function at the currently-loaded primal point.
Hypatia.Cones.hess
— Functionhess(cone)
The Hessian (symmetric positive definite) of the cone's barrier function at the currently-loaded primal point.
Hypatia.Cones.inv_hess
— Functioninv_hess(cone)
The inverse Hessian (symmetric positive definite) of the cone's barrier function at the currently-loaded primal point.
Hypatia.Cones.hess_prod!
— Functionhess_prod!(prod, arr, cone)
Compute the product of the Hessian of the cone's barrier function at the currently-loaded primal point with a vector or array, in-place.
Hypatia.Cones.inv_hess_prod!
— Functioninv_hess_prod!(prod, arr, cone)
Compute the product of the inverse Hessian of the cone's barrier function at the currently-loaded primal point with a vector or array, in-place.
Hypatia.Cones.use_dder3
— Functionuse_dder3(_)
Returns true if and only if the oracle for the third-order directional derivative oracle dder3
can be computed.
Hypatia.Cones.dder3
— Functiondder3(cone, dir)
Compute the third-order directional derivative, in the direction dir
, the cone's barrier function at the currently-loaded primal point.
Predefined cone types
Hypatia.Cones.Nonnegative
— Typemutable struct Nonnegative{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Nonnegative cone of dimension dim
.
Nonnegative{T}(dim::Int)
Hypatia.Cones.PosSemidefTri
— Typemutable struct PosSemidefTri{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Real symmetric or complex Hermitian positive semidefinite cone of dimension dim
in svec (scaled upper triangle) format.
PosSemidefTri{T, R}(dim::Int)
Hypatia.Cones.DoublyNonnegativeTri
— Typemutable struct DoublyNonnegativeTri{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Real symmetric doubly nonnegative cone (intersection of nonnegative and positive semidefinite cones) of dimension dim
in svec format.
DoublyNonnegativeTri{T}(dim::Int, use_dual::Bool = false)
Hypatia.Cones.PosSemidefTriSparse
— Typemutable struct PosSemidefTriSparse{I<:Hypatia.Cones.PSDSparseImpl, T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Real symmetric or complex Hermitian sparse positive semidefinite cone of side dimension side
and sparse lower triangle row and column indices rows
, cols
in svec format. Note all diagonal elements must be present.
PosSemidefTriSparse{T, R}(side::Int, rows::Vector{Int}, cols::Vector{Int}, use_dual::Bool = false)
Hypatia.Cones.LinMatrixIneq
— Typemutable struct LinMatrixIneq{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Linear matrix inequality cone parametrized by list of real symmetric or complex Hermitian matrices mats
of equal dimension.
LinMatrixIneq{T}(mats::Vector, use_dual::Bool = false)
Hypatia.Cones.EpiNormInf
— Typemutable struct EpiNormInf{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Epigraph of real or complex infinity norm cone of dimension dim
.
EpiNormInf{T, R}(dim::Int, use_dual::Bool = false)
Hypatia.Cones.EpiNormEucl
— Typemutable struct EpiNormEucl{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Epigraph of Euclidean norm (AKA second-order) cone of dimension dim
.
EpiNormEucl{T}(dim::Int)
Hypatia.Cones.EpiPerSquare
— Typemutable struct EpiPerSquare{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Epigraph of perspective function of halved squared Euclidean norm (AKA rotated second-order) cone of dimension dim
.
EpiPerSquare{T}(dim::Int)
Hypatia.Cones.EpiNormSpectralTri
— Typemutable struct EpiNormSpectralTri{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Epigraph of real symmetric or complex Hermitian matrix spectral norm (i.e. maximum absolute value of eigenvalues) cone of dimension dim
in svec format.
EpiNormSpectralTri{T, R}(dim::Int, use_dual::Bool = false)
Hypatia.Cones.EpiNormSpectral
— Typemutable struct EpiNormSpectral{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Epigraph of real or complex matrix spectral norm (i.e. maximum singular value) for a matrix (stacked column-wise) of nrows
rows and ncols
columns with nrows ≤ ncols
.
EpiNormSpectral{T, R}(nrows::Int, ncols::Int, use_dual::Bool = false)
Hypatia.Cones.MatrixEpiPerSquare
— Typemutable struct MatrixEpiPerSquare{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Matrix epigraph of perspective function of real or complex matrix outer product for a matrix (stacked column-wise) of nrows
rows and ncols
columns with nrows ≤ ncols
.
MatrixEpiPerSquare{T, R}(nrows::Int, ncols::Int, use_dual::Bool = false)
Hypatia.Cones.GeneralizedPower
— Typemutable struct GeneralizedPower{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Generalized power cone parametrized by powers α
in the unit simplex and dimension d
of the normed variables.
GeneralizedPower{T}(α::Vector{T}, d::Int, use_dual::Bool = false)
Hypatia.Cones.HypoPowerMean
— Typemutable struct HypoPowerMean{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Hypograph of weighted power mean cone parametrized by powers α
in the unit simplex.
HypoPowerMean{T}(α::Vector{T}, use_dual::Bool = false)
Hypatia.Cones.HypoGeoMean
— Typemutable struct HypoGeoMean{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Hypograph of geometric mean cone of dimension dim
.
HypoGeoMean{T}(dim::Int, use_dual::Bool = false)
Hypatia.Cones.HypoRootdetTri
— Typemutable struct HypoRootdetTri{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Hypograph of real symmetric or complex Hermitian root-determinant cone of dimension dim
in svec format.
HypoRootdetTri{T, R}(dim::Int, use_dual::Bool = false)
Hypatia.Cones.HypoPerLog
— Typemutable struct HypoPerLog{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Hypograph of perspective function of sum-log cone of dimension dim
.
HypoPerLog{T}(dim::Int, use_dual::Bool = false)
Hypatia.Cones.HypoPerLogdetTri
— Typemutable struct HypoPerLogdetTri{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Hypograph of perspective function of real symmetric or complex Hermitian log-determinant cone of dimension dim
in svec format.
HypoPerLogdetTri{T, R}(dim::Int, use_dual::Bool = false)
Hypatia.Cones.EpiPerSepSpectral
— Typemutable struct EpiPerSepSpectral{Q<:Hypatia.Cones.ConeOfSquares, T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Epigraph of perspective function of a convex separable spectral function h
over a cone of squares Q
on a Jordan algebra with rank d
.
EpiPerSepSpectral{Q, T}(h::Hypatia.Cones.SepSpectralFun, d::Int, use_dual::Bool = false)
Hypatia.Cones.EpiRelEntropy
— Typemutable struct EpiRelEntropy{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Epigraph of vector relative entropy cone of dimension dim
.
EpiRelEntropy{T}(dim::Int, use_dual::Bool = false)
Hypatia.Cones.EpiTrRelEntropyTri
— Typemutable struct EpiTrRelEntropyTri{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Epigraph of matrix relative entropy cone of dimension dim
in svec format.
EpiTrRelEntropyTri{T}(dim::Int, use_dual::Bool = false)
Hypatia.Cones.WSOSInterpNonnegative
— Typemutable struct WSOSInterpNonnegative{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.Cone{T<:Real}
Interpolant-basis weighted sum-of-squares polynomial cone of dimension U
, for real or real-valued complex polynomials, parametrized by vector of matrices Ps
derived from interpolant basis and polynomial domain constraints.
WSOSInterpNonnegative{T, R}(U::Int, Ps::Vector{Matrix{R}}, use_dual::Bool = false)
Hypatia.Cones.WSOSInterpPosSemidefTri
— Typemutable struct WSOSInterpPosSemidefTri{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Interpolant-basis weighted sum-of-squares polynomial (of dimension U
) positive semidefinite matrix (of side dimension R
) cone, parametrized by vector of matrices Ps
derived from interpolant basis and polynomial domain constraints.
WSOSInterpPosSemidefTri{T}(R::Int, U::Int, Ps::Vector{Matrix{T}}, use_dual::Bool = false)
Hypatia.Cones.WSOSInterpEpiNormEucl
— Typemutable struct WSOSInterpEpiNormEucl{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Interpolant-basis weighted sum-of-squares polynomial (of dimension U
) epigraph of ℓ₂
norm (of dimension R
) cone, parametrized by vector of matrices Ps
derived from interpolant basis and polynomial domain constraints.
WSOSInterpEpiNormEucl{T}(R::Int, U::Int, Ps::Vector{Matrix{T}}, use_dual::Bool = false)
Hypatia.Cones.WSOSInterpEpiNormOne
— Typemutable struct WSOSInterpEpiNormOne{T<:Real} <: Hypatia.Cones.Cone{T<:Real}
Interpolant-basis weighted sum-of-squares polynomial (of dimension U
) epigraph of ℓ₁
norm (of dimension R
) cone, parametrized by vector of matrices Ps
derived from interpolant basis and polynomial domain constraints.
WSOSInterpEpiNormOne{T}(R::Int, U::Int, Ps::Vector{Matrix{T}}, use_dual::Bool = false)
Helpers for PosSemidefTriSparse
Hypatia.Cones.PSDSparseImpl
— Typeabstract type PSDSparseImpl
An implementation type for the sparse positive semidefinite cone PosSemidefTriSparse
.
Hypatia.Cones.PSDSparseDense
— Typestruct PSDSparseDense <: Hypatia.Cones.PSDSparseImpl
Dense Cholesky-based implementation for the sparse positive semidefinite cone PosSemidefTriSparse
.
Hypatia.Cones.PSDSparseCholmod
— Typestruct PSDSparseDense <: Hypatia.Cones.PSDSparseImpl
Dense Cholesky-based implementation for the sparse positive semidefinite cone PosSemidefTriSparse
.
EpiPerSepSpectral helpers
Hypatia.Cones.ConeOfSquares
— Typeabstract type ConeOfSquares{T<:Real}
A cone of squares on a Jordan algebra.
Hypatia.Cones.VectorCSqr
— Typestruct VectorCSqr{T<:Real} <: Hypatia.Cones.ConeOfSquares{T<:Real}
Real vector cone of squares.
Hypatia.Cones.MatrixCSqr
— Typestruct MatrixCSqr{T<:Real, R<:Union{Complex{T<:Real}, T<:Real}} <: Hypatia.Cones.ConeOfSquares{T<:Real}
Real symmetric or complex Hermitian positive semidefinite cone of squares.
Hypatia.Cones.vector_dim
— Functionvector_dim(
_::Type{var"#s68"} where var"#s68"<:Hypatia.Cones.VectorCSqr,
d::Int64
) -> Int64
The rank of the vector cone of squares, equal to the vector length.
vector_dim(
_::Type{var"#s67"} where var"#s67"<:(Hypatia.Cones.MatrixCSqr{var"#s66", R} where var"#s66"<:Real),
d::Int64
) -> Int64
The rank of the matrix cone of squares, equal to the side dimension of the matrix.
Hypatia.Cones.SepSpectralFun
— Typeabstract type SepSpectralFun
A univariate convex function defined on positive reals.
Hypatia.Cones.NegLogSSF
— Typestruct NegLogSSF <: Hypatia.Cones.SepSpectralFun
The negative logarithm function $x \to - \log(x)$.
Hypatia.Cones.NegEntropySSF
— Typestruct NegEntropySSF <: Hypatia.Cones.SepSpectralFun
The negative entropy function $x \to x \log(x)$.
Hypatia.Cones.NegSqrtSSF
— Typestruct NegSqrtSSF <: Hypatia.Cones.SepSpectralFun
The negative square root function $x \to -x^{1/2}$. Note this is a special case of the negative power: NegPower01SSF(0.5)
.
Hypatia.Cones.NegPower01SSF
— Typestruct NegPower01SSF <: Hypatia.Cones.SepSpectralFun
The negative power function $x \to -x^p$ parametrized by $p \in (0, 1)$
Hypatia.Cones.Power12SSF
— Typestruct Power12SSF <: Hypatia.Cones.SepSpectralFun
The power function $x \to x^p$ parametrized by $p \in (1, 2]$. Note for $p = 2$, it is more efficient to use EpiPerSquare
.