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 complex PSD cone #194

Merged
merged 7 commits into from
Oct 17, 2024
Merged
Show file tree
Hide file tree
Changes from 5 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
1 change: 1 addition & 0 deletions src/COSMO.jl
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@ using Reexport

export assemble!, warm_start!, empty_model!, update!

const RealOrComplex{T <: Real} = Union{T, Complex{T}}
const DefaultFloat = Float64
const DefaultInt = LinearAlgebra.BlasInt

Expand Down
10 changes: 8 additions & 2 deletions src/MOI_wrapper.jl
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,8 @@ const IntervalConvertible = Union{GreaterThan, Interval}
const Zeros = MOI.Zeros
const SOC = MOI.SecondOrderCone
const PSDT = MOI.Scaled{MOI.PositiveSemidefiniteConeTriangle}
const SupportedVectorSets = Union{Zeros, MOI.Nonnegatives, SOC, PSDT, MOI.ExponentialCone, MOI.DualExponentialCone, MOI.PowerCone, MOI.DualPowerCone}
const ComplexPSDT = MOI.Scaled{MOI.HermitianPositiveSemidefiniteConeTriangle}
const SupportedVectorSets = Union{Zeros, MOI.Nonnegatives, SOC, PSDT, ComplexPSDT, MOI.ExponentialCone, MOI.DualExponentialCone, MOI.PowerCone, MOI.DualPowerCone}
const AggregatedSets = Union{Zeros, MOI.Nonnegatives, MOI.GreaterThan}
aggregate_equivalent(::Type{<: MOI.Zeros}) = COSMO.ZeroSet
aggregate_equivalent(::Type{<: Union{MOI.GreaterThan, MOI.Nonnegatives}}) = COSMO.Nonnegatives
Expand Down Expand Up @@ -449,7 +450,12 @@ end


function processSet!(b::Vector{T}, rows::UnitRange{Int}, cs, s::PSDT) where {T <: AbstractFloat}
push!(cs, COSMO.PsdConeTriangle{T}(length(rows)))
push!(cs, COSMO.PsdConeTriangle{T, T}(length(rows)))
nothing
end

function processSet!(b::Vector{T}, rows::UnitRange{Int}, cs, s::ComplexPSDT) where {T <: AbstractFloat}
push!(cs, COSMO.PsdConeTriangle{T, Complex{T}}(length(rows)))
nothing
end

Expand Down
12 changes: 6 additions & 6 deletions src/algebra.jl
Original file line number Diff line number Diff line change
Expand Up @@ -222,23 +222,23 @@ function symmetrize_full!(A::AbstractMatrix{T}) where {T <: AbstractFloat}
nothing
end

# this function assumes real symmetric X and only considers the upper triangular part
function is_pos_def!(X::AbstractMatrix{T}, tol::T=zero(T)) where T
# this function assumes real or complex Hermitian X and only considers the upper triangular part
function is_pos_def!(X::AbstractMatrix{<:RealOrComplex{T}}, tol::T=zero(T)) where {T}
# See https://math.stackexchange.com/a/13311
@inbounds for i = 1:size(X, 1)
X[i, i] += tol
end
F = cholesky!(Symmetric(X), check = false)
F = cholesky!(Hermitian(X), check = false)
return issuccess(F)
end

function is_neg_def!(X::AbstractMatrix{T}, tol::T=zero(T)) where T
function is_neg_def!(X::AbstractMatrix{<:RealOrComplex{T}}, tol::T=zero(T)) where {T}
@. X *= -one(T)
return is_pos_def!(X, tol)
end

is_pos_def(X::AbstractMatrix{T}, tol::T=zero(T)) where T = is_pos_def!(copy(X), tol)
is_neg_def(X::AbstractMatrix{T}, tol::T=zero(T)) where T = is_pos_def!(-X, tol)
is_pos_def(X::AbstractMatrix{<:RealOrComplex{T}}, tol::T=zero(T)) where {T} = is_pos_def!(copy(X), tol)
is_neg_def(X::AbstractMatrix{<:RealOrComplex{T}}, tol::T=zero(T)) where {T} = is_pos_def!(-X, tol)


"Round x to the closest multiple of N."
Expand Down
3 changes: 2 additions & 1 deletion src/chordal_decomposition/chordal_decomposition.jl
Original file line number Diff line number Diff line change
Expand Up @@ -75,7 +75,8 @@ function _analyse_sparsity_pattern(ci::ChordalInfo{T}, csp::Array{Int, 1}, sets:
end

DenseEquivalent(C::COSMO.PsdCone{T}, dim::Int) where {T <: AbstractFloat} = COSMO.DensePsdCone{T}(dim)
DenseEquivalent(C::COSMO.PsdConeTriangle{T}, dim::Int) where {T <: AbstractFloat} = COSMO.DensePsdConeTriangle{T}(dim)
DenseEquivalent(C::COSMO.PsdConeTriangle{T, T}, dim::Int) where {T <: AbstractFloat} = COSMO.DensePsdConeTriangle{T, T}(dim)
DenseEquivalent(C::COSMO.PsdConeTriangle{T, Complex{T}}, dim::Int) where {T <: AbstractFloat} = COSMO.DensePsdConeTriangle{T, Complex{T}}(dim)

function nz_rows(a::SparseMatrixCSC{T}, ind::UnitRange{Int}, DROP_ZEROS_FLAG::Bool) where {T <: AbstractFloat}
DROP_ZEROS_FLAG && dropzeros!(a)
Expand Down
12 changes: 11 additions & 1 deletion src/constraint.jl
Original file line number Diff line number Diff line change
Expand Up @@ -93,7 +93,17 @@ function Constraint{T}(
# this constructor doesnt work with cones that need special arguments like the power cone
set_type <: ArgumentCones && error("You can't create a constraint by passing the convex set as a type, if your convex set is a $(set_type). Please pass an object.")
# call the appropriate set constructor
convex_set = set_type{T}(size(A, 1))
n = size(A, 1)
# we can deduce whether the PsdCone must be real or complex from the dimension
if set_type <: PsdConeTriangles
if n == 1 || isqrt(n)^2 != n
convex_set = set_type{T, T}(n)
else
convex_set = set_type{T, Complex{T}}(n)
end
else
convex_set = set_type{T}(n)
end
Constraint{T}(A, b, convex_set, args...)
end

Expand Down
Loading