jump-dev / Pajarito.jl

A solver for mixed-integer convex optimization
Mozilla Public License 2.0
129 stars 22 forks source link

Not able to solve an MISDP problem #423

Closed taranjitsinghb closed 2 years ago

taranjitsinghb commented 5 years ago

Hi, I am considering to solve a misdp problem. I found this example online and thought, I would try to see if this works. However, It always get stuck at creating a diagonal matrix of an array variable.


using JuMP, Pajarito, LinearAlgebra

# Set up QP JuMP model, solve, print solution
# xB is a bound on the absolute values of the estimate variables x_j
# Solver can be either MIQP/MINLP or MICP
# Model is a bit "low level" in order to be compatible with MathProgBase ConicToLPQPBridge

# Set up MISDP JuMP model, solve, print solution
# No bound xB is needed in this model
function misdp_cardls(m, d, A, b, k, rho, solver)
    mod = Model(solver=solver)
    @variable(mod, tau)
    @variable(mod, z[1:d],Int)
    @objective(mod, Min, tau)
    @constraint(mod, sum(z) <= k)
    Im = Matrix{Float64}(I,m,m);
    dz = diagm(z);
    @SDconstraint(mod, [(Im + 1/rho*A*dz*A') b ; b' tau] >= 0)
    solve(mod)
    println("  selected features (z) = \n$(getvalue(z))\n")
end

#=========================================================
Choose solvers and options
=========================================================#

mip_solver_drives = true
rel_gap = 1e-5

# using Cbc
# mip_solver = CbcSolver()

using CPLEX
mip_solver = CplexSolver(
    CPX_PARAM_SCRIND=(mip_solver_drives ? 1 : 0),
    # CPX_PARAM_SCRIND=1,
    CPX_PARAM_EPINT=1e-8,
    CPX_PARAM_EPRHS=1e-7,
    CPX_PARAM_EPGAP=(mip_solver_drives ? 1e-5 : 1e-9)
)

# using SCS
# conic_solver = SCSSolver(eps=1e-6, max_iters=1000000, verbose=0)

using Mosek
conic_solver = MosekSolver(LOG=0)

micp_solver = PajaritoSolver(
    mip_solver_drives=mip_solver_drives,
    log_level=3,
    rel_gap=rel_gap,
    mip_solver=mip_solver,
    cont_solver=conic_solver,
)

#=========================================================
Specify/generate data
=========================================================#

d = 6  # Dimension of feature space
m = 20  # Number of samples

#srand(100)       # Change or comment random seed to get different data
A = randn(m, d)  # Sample point matrix (rows are samples)
b = randn(m)     # Sample measurement vector
# @show A
# @show b

k = floor(Int, d/2)  # Number of features to select (||x||_0 <= k)

rho = 1.  # Ridge regularization multiplier

xB = 4  # Bound on absolute values of estimate variables (|x_j| <= xB)

#=========================================================
Solve JuMP models
=========================================================#

println("\n\n****MISDP model****\n")
misdp_cardls(m, d, A, b, k, rho, micp_solver)
mtanneau commented 5 years ago

It always get stuck at creating a diagonal matrix of an array variable.

Can you provide the error message you're getting? Are you using Julia v1 (or later)?

If you're using Julia v1.0 or later, then my best guess is that your error comes from this line (in the misdp_cardls function)

    dz = diagm(z);

The diagm syntac is deprecated in v1, instead you should use

    dz = Diagonal(z)
taranjitsinghb commented 5 years ago

@mtanneau Thank you for the tip. It was indeed getting stuck there. I changed it to:

dz = Diagonal(z);

However, it still gets stuck at the conversion of Diagonal as variable type for the expression

1/rho*A*dz*A'

it throuws the following error

ERROR: LoadError: MethodError: Cannot `convert` an object of type JuMP.GenericQuadExpr{Float64,Varia
ble} to an object of type JuMP.GenericAffExpr{Float64,Variable}
mtanneau commented 5 years ago

This looks more like a JuMP-related than Pajarito-related error... I would suggest you open an issue with JuMP directly

erikfilias commented 5 years ago

@mtanneau Thank you for the tip. It was indeed getting stuck there. I changed it to:

dz = Diagonal(z);

However, it still gets stuck at the conversion of Diagonal as variable type for the expression

1/rho*A*dz*A'

it throuws the following error

ERROR: LoadError: MethodError: Cannot `convert` an object of type JuMP.GenericQuadExpr{Float64,Varia
ble} to an object of type JuMP.GenericAffExpr{Float64,Variable}

You can make:

1/rho*A*(Matrix{Float64}(I, d, d).*z)*A'

It works for me.

odow commented 2 years ago

Closing as stale and because this doesn't look like a bug in Pajarito.