Add tests for atoms/second_order_cone (#561)

src/atoms/second_order_cone/EuclideanNormAtom.jl

@@ -0,0 +1,32 @@
+mutable struct EuclideanNormAtom <: AbstractExpr
+    children::Tuple{AbstractExpr}
+    size::Tuple{Int,Int}
+
+    EuclideanNormAtom(x::AbstractExpr) = new((x,), (1, 1))
+end
+
+head(io::IO, ::EuclideanNormAtom) = print(io, "norm2")
+
+Base.sign(::EuclideanNormAtom) = Positive()
+
+monotonicity(x::EuclideanNormAtom) = (sign(x.children[1]) * Nondecreasing(),)
+
+curvature(::EuclideanNormAtom) = ConvexVexity()
+
+evaluate(x::EuclideanNormAtom) = norm(vec(evaluate(x.children[1])))
+
+function new_conic_form!(context::Context{T}, A::EuclideanNormAtom) where {T}
+    x = only(AbstractTrees.children(A))
+    t = conic_form!(context, Variable())
+    obj = conic_form!(context, x)
+    f = operate(vcat, T, sign(A), t, obj)
+    MOI_add_constraint(context.model, f, MOI.SecondOrderCone(length(x) + 1))
+    return t
+end
+
+function LinearAlgebra.norm2(x::AbstractExpr)
+    if sign(x) == ComplexSign()
+        return EuclideanNormAtom([real(x); imag(x)])
+    end
+    return EuclideanNormAtom(x)
+end

src/atoms/second_order_cone/GeoMeanAtom.jl

@@ -1,16 +1,19 @@
+_to_expr(x::AbstractExpr) = x
+_to_expr(x::Value) = constant(x)
+
 mutable struct GeoMeanAtom <: AbstractExpr
-    children::NTuple{N,AbstractExpr} where {N}
+    children::Vector{AbstractExpr}
     size::Tuple{Int,Int}
 
     function GeoMeanAtom(args::AbstractExprOrValue...)
-        args = Tuple(arg isa Value ? constant(arg) : arg for arg in args)
-        sz = size(first(args))
-        if any(!=(sz), size.(args))
-            error("geo mean must take arguments of the same size")
-        elseif any(x -> sign(x) == ComplexSign(), args)
-            error("The arguments should be real, not complex")
+        new_args = AbstractExpr[_to_expr(arg) for arg in args]
+        sz = size(first(new_args))
+        if any(!=(sz), size.(new_args))
+            error("[GeoMeanAtom] geomean must take arguments of the same size")
+        elseif any(x -> sign(x) == ComplexSign(), new_args)
+            error("[GeoMeanAtom] the arguments must be real, not complex")
         end
-        return new(args, sz)
+        return new(new_args, sz)
     end
 end
 
@@ -19,19 +22,14 @@ head(io::IO, ::GeoMeanAtom) = print(io, "geomean")
 Base.sign(::GeoMeanAtom) = Positive()
 
 function monotonicity(q::GeoMeanAtom)
-    return fill(Nondecreasing(), length(q.children))
+    return ntuple(_ -> Nondecreasing(), length(q.children))
 end
 
 curvature(::GeoMeanAtom) = ConcaveVexity()
 
 _geomean(scalar_args...) = prod(scalar_args)^(1 / length(scalar_args))
 
-function evaluate(q::GeoMeanAtom)
-    n = length(q.children)
-    children = evaluate.(q.children)
-    c1 = first(children)
-    return [_geomean((children[i][I] for i in 1:n)...) for I in eachindex(c1)]
-end
+evaluate(q::GeoMeanAtom) = _geomean.(evaluate.(q.children)...)
 
 function new_conic_form!(context::Context{T}, q::GeoMeanAtom) where {T}
     n = length(q.children)

src/atoms/second_order_cone/HuberAtom.jl

@@ -5,9 +5,9 @@ mutable struct HuberAtom <: AbstractExpr
 
     function HuberAtom(x::AbstractExpr, M::Real)
         if sign(x) == ComplexSign()
-            error("Argument must be real")
+            error("[HuberAtom] argument must be real")
         elseif M <= 0
-            error("Huber parameter must by a positive scalar")
+            error("[HuberAtom] parameter must by a positive scalar. Got `M=$M`")
         end
         return new((x,), x.size, M)
     end
@@ -38,7 +38,6 @@ function new_conic_form!(context::Context, x::HuberAtom)
     s = Variable(c.size)
     n = Variable(c.size)
     add_constraint!(context, c == s + n)
-    # objective given by s.^2 + 2 * M * |n|
     return conic_form!(context, square(s) + 2 * x.M * abs(n))
 end
 

src/atoms/second_order_cone/QolElemAtom.jl

@@ -5,7 +5,7 @@ mutable struct QolElemAtom <: AbstractExpr
     function QolElemAtom(x::AbstractExpr, y::AbstractExpr)
         if x.size != y.size
             error(
-                "elementwise quad over lin must take two arguments of the same size",
+                "[QolElemAtom] elementwise quad over lin must take two arguments of the same size",
             )
         end
         return new((x, y), x.size)
@@ -27,39 +27,32 @@ function evaluate(q::QolElemAtom)
 end
 
 function new_conic_form!(context::Context{T}, q::QolElemAtom) where {T}
-    sz = q.children[1].size
-    t = Variable(sz[1], sz[2])
-    t_obj = conic_form!(context, t)
     x, y = q.children
-    for i in 1:length(q.children[1])
+    t = Variable(x.size)
+    for i in 1:length(x)
         add_constraint!(
             context,
             SecondOrderConeConstraint(y[i] + t[i], y[i] - t[i], 2 * x[i]),
         )
     end
     add_constraint!(context, y >= 0)
-    return t_obj
+    return conic_form!(context, t)
 end
 
 qol_elementwise(x::AbstractExpr, y::AbstractExpr) = QolElemAtom(x, y)
 
-function Base.Broadcast.broadcasted(::typeof(^), x::AbstractExpr, k::Int)
-    if k != 2
-        error("raising variables to powers other than 2 is not implemented")
+function square(x::AbstractExpr)
+    if sign(x) == ComplexSign()
+        error(
+            "square of a complex number is not DCP. Did you mean square_modulus?",
+        )
     end
-    return QolElemAtom(x, constant(ones(x.size[1], x.size[2])))
+    return QolElemAtom(x, constant(ones(x.size)))
 end
 
-function Base.Broadcast.broadcasted(
-    ::typeof(Base.literal_pow),
-    ::typeof(^),
-    x::AbstractExpr,
-    ::Val{k},
-) where {k}
-    return Base.Broadcast.broadcasted(^, x, k)
-end
+sumsquares(x::AbstractExpr) = square(norm2(x))
 
-invpos(x::AbstractExpr) = QolElemAtom(constant(ones(x.size[1], x.size[2])), x)
+invpos(x::AbstractExpr) = QolElemAtom(constant(ones(x.size)), x)
 
 function Base.Broadcast.broadcasted(::typeof(/), x::Value, y::AbstractExpr)
     return dotmultiply(constant(x), invpos(y))
@@ -72,13 +65,18 @@ function Base.:/(x::Value, y::AbstractExpr)
     return MultiplyAtom(constant(x), invpos(y))
 end
 
-sumsquares(x::AbstractExpr) = square(norm2(x))
-
-function square(x::AbstractExpr)
-    if sign(x) == ComplexSign()
-        error(
-            "Square of complex number is not DCP. Did you mean square_modulus?",
-        )
+function Base.Broadcast.broadcasted(::typeof(^), x::AbstractExpr, k::Int)
+    if k != 2
+        error("raising variables to powers other than 2 is not implemented")
     end
-    return QolElemAtom(x, constant(ones(x.size[1], x.size[2])))
+    return square(x)
+end
+
+function Base.Broadcast.broadcasted(
+    ::typeof(Base.literal_pow),
+    ::typeof(^),
+    x::AbstractExpr,
+    ::Val{k},
+) where {k}
+    return Base.Broadcast.broadcasted(^, x, k)
 end

src/atoms/second_order_cone/QuadOverLinAtom.jl

@@ -3,8 +3,10 @@ mutable struct QuadOverLinAtom <: AbstractExpr
     size::Tuple{Int,Int}
 
     function QuadOverLinAtom(x::AbstractExpr, y::AbstractExpr)
-        if x.size[2] != 1 && y.size != (1, 1)
-            error("quad over lin arguments must be a vector and a scalar")
+        if x.size[2] != 1 || y.size != (1, 1)
+            error(
+                "[QuadOverLinAtom] quadoverlin arguments must be a vector and a scalar",
+            )
         end
         return new((x, y), (1, 1))
     end

src/variable.jl

@@ -102,9 +102,10 @@ sign!(x::AbstractVariable, s::Sign) = x.sign = s
 Returns the current value of `x` if assigned; errors otherwise.
 """
 function evaluate(x::AbstractVariable)
-    return _value(x) === nothing ?
-           error("Value of the variable is yet to be calculated") :
-           output(_value(x))
+    if _value(x) === nothing
+        error("Value of the variable is yet to be calculated")
+    end
+    return output(copy(_value(x)))
 end
 
 """