Custom cone not working (not even tutorial example) without in_dual.

[Timeroot]

Hi, when I use a custom cone, I run into issues where it's trying to call in_pol_recc(::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::Nonpositives{Float64}, ::Float64) even though I haven't defined that method. The tutorial seemed to say that I wouldn't need to implement this.

Specifically, I'm running COSMO v0.8.1, and JuMP v0.21.9. If you run the snippet,

using COSMO, LinearAlgebra, SparseArrays, Test, JuMP

# define new cone type
struct Nonpositives{T} <: COSMO.AbstractConvexCone{T}
    dim::Int64
end

function COSMO.project!(x::AbstractVector{T}, C::Nonpositives{T}) where {T <: AbstractFloat}
    x .= min.(x, zero(T))
end
import Base.copy

struct NonPos <: MOI.AbstractVectorSet
	dimension::Int
end
Base.copy(set::NonPos) = set #this is needed for MOI

function COSMO.processSet!(b::Vector{T}, rows::UnitRange{Int}, cs, s::NonPos) where {T <: AbstractFloat}
    push!(cs, Nonpositives{T}(length(rows)))
    nothing
end
#tell MOI, that COSMO supports constraints with this new cone
MOI.supports_constraint(optimizer::COSMO.Optimizer, ::Type{<:Union{MOI.VectorOfVariables, MOI.VectorAffineFunction{Float64}}}, ::Type{NonPos}) = true

A1 = diagm(0 => ones(2))
b1 = [-3.; -2.]

model = JuMP.Model(COSMO.Optimizer);
@variable(model, x[1:20]);
@objective(model, Max, x[1] + x[2] + x[3]);
@constraint(model, A1 * x[1:2] .+ b1 in NonPos(2));
@constraint(model, x[1] == 1);
JuMP.optimize!(model)

You'll get the error

------------------------------------------------------------------
          COSMO v0.8.1 - A Quadratic Objective Conic Solver
                         Michael Garstka
                University of Oxford, 2017 - 2021
------------------------------------------------------------------

Problem:  x ∈ R^{20},
          constraints: A ∈ R^{3x20} (3 nnz),
          matrix size to factor: 23x23,
          Floating-point precision: Float64
Sets:     Nonpositives of dim: 2
          ZeroSet of dim: 1
Settings: ϵ_abs = 1.0e-05, ϵ_rel = 1.0e-05,
          ϵ_prim_inf = 1.0e-04, ϵ_dual_inf = 1.0e-04,
          ρ = 0.1, σ = 1e-06, α = 1.6,
          max_iter = 5000,
          scaling iter = 10 (on),
          check termination every 25 iter,
          check infeasibility every 40 iter,
          KKT system solver: QDLDL
Acc:      Anderson Type2{QRDecomp},
          Memory size = 15, RestartedMemory,	
          Safeguarded: true, tol: 2.0
Setup Time: 0.19ms

Iter:	Objective:	Primal Res:	Dual Res:	Rho:
1	-1.6000e+06	1.7200e+01	1.0000e+00	1.0000e-01
25	-7.3992e+07	4.6277e-06	1.0000e+00	1.0000e-01
ERROR: MethodError: no method matching in_pol_recc(::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::Nonpositives{Float64}, ::Float64)
Closest candidates are:
  in_pol_recc(::AbstractArray{T,1}, ::COSMO.ZeroSet{T}, ::T) where T at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/convexset.jl:33
  in_pol_recc(::AbstractArray{T,1}, ::COSMO.Nonnegatives{T}, ::T) where T at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/convexset.jl:79
  in_pol_recc(::AbstractArray{T,1}, ::COSMO.SecondOrderCone, ::T) where T at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/convexset.jl:119
  ...
Stacktrace:
 [1] in_pol_recc!(::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::Nonpositives{Float64}, ::Float64) at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/convexset.jl:874
 [2] (::COSMO.var"#35#36"{Float64})(::Tuple{SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true},Nonpositives{Float64}}) at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/convexset.jl:831
 [3] _all(::COSMO.var"#35#36"{Float64}, ::Base.Iterators.Zip{Tuple{Array{SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true},1},Array{COSMO.AbstractConvexSet,1}}}, ::Colon) at ./reduce.jl:820
 [4] all at ./reduce.jl:815 [inlined]
 [5] in_pol_recc! at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/convexset.jl:831 [inlined]
 [6] is_dual_infeasible!(::Array{Float64,1}, ::COSMO.Workspace{Float64}) at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/infeasibility.jl:61
 [7] check_termination!(::COSMO.Workspace{Float64}, ::COSMO.Settings{Float64}, ::Int64, ::Float64, ::Symbol, ::Float64, ::Float64, ::Float64, ::Int64) at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/solver.jl:345
 [8] optimize!(::COSMO.Workspace{Float64}) at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/solver.jl:161
 [9] optimize! at /home/timeroot/.julia/packages/COSMO/Sb5eV/src/MOI_wrapper.jl:155 [inlined]
 [10] optimize!(::MathOptInterface.Utilities.CachingOptimizer{Optimizer{Float64},MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.GenericModel{Float64,MathOptInterface.Utilities.ModelFunctionConstraints{Float64}}}}) at /home/timeroot/.julia/packages/MathOptInterface/YDdD3/src/Utilities/cachingoptimizer.jl:252
 [11] optimize!(::MathOptInterface.Bridges.LazyBridgeOptimizer{MathOptInterface.Utilities.CachingOptimizer{Optimizer{Float64},MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.GenericModel{Float64,MathOptInterface.Utilities.ModelFunctionConstraints{Float64}}}}}) at /home/timeroot/.julia/packages/MathOptInterface/YDdD3/src/Bridges/bridge_optimizer.jl:319
 [12] optimize!(::MathOptInterface.Utilities.CachingOptimizer{MathOptInterface.AbstractOptimizer,MathOptInterface.Utilities.UniversalFallback{MathOptInterface.Utilities.GenericModel{Float64,MathOptInterface.Utilities.ModelFunctionConstraints{Float64}}}}) at /home/timeroot/.julia/packages/MathOptInterface/YDdD3/src/Utilities/cachingoptimizer.jl:252
 [13] optimize!(::Model, ::Nothing; bridge_constraints::Bool, ignore_optimize_hook::Bool, kwargs::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at /home/timeroot/.julia/packages/JuMP/klrjG/src/optimizer_interface.jl:185
 [14] optimize! at /home/timeroot/.julia/packages/JuMP/klrjG/src/optimizer_interface.jl:157 [inlined] (repeats 2 times)
 [15] top-level scope at REPL[48]:1

If you play with the constants a bit, you'll instead get

ERROR: MethodError: no method matching in_dual(::SubArray{Float64,1,Array{Float64,1},Tuple{UnitRange{Int64}},true}, ::Nonpositives{Float64}, ::Float64)

It seems like something is happening where some intermediate processing tries to take a slice, in a way that isn't supported properly by custom cones. I'm guessing something wrong with how it choses whether to use in_dual or not.

[migarstka]

Thank you @Timeroot for the bug report. I will take a look at this.

[migarstka]

Hi @Timeroot, sorry for the delayed reply. I was on holiday.
I just ran the documentation script with COSMO v0.8.1, and I don't get any error (because in_dual and in_pol_recc are defined).

For your example a quick fix would be to add:

set_optimizer_attribute(model, "check_infeasibility", 50_000)

to turn the infeasibility checks off.