Jump.jl: bikeshedding default solver behavior

[mlubin]

Ok, time to figure out the user-facing side of default solvers. As of https://github.com/JuliaOpt/JuMP.jl/pull/959, we have the following significant changes:

1. An informational message is printed whenever no solver is satisfied but one happens to be available.
2. On Julia 0.6, users have to explicitly load solver packages to make them available.
3. On Julia 0.5, there are no breaking changes yet, but we now print a warning when we have to explicitly load a solver in a way that will break on 0.6.

The question is, how do we present example code now?

using JuMP

m = Model()
@variable(m, 0 <= x <= 2)
@variable(m, 0 <= y <= 30)
@objective(m, Max, 5x + 3y)
@constraint(m, 1x + 5y <= 3.0)
status = solve(m)

will no longer work. Some options are:

Option 1:

using JuMP, Clp

m = Model()
@variable(m, 0 <= x <= 2)
@variable(m, 0 <= y <= 30)
@objective(m, Max, 5x + 3y)
@constraint(m, 1x + 5y <= 3.0)
status = solve(m)
# prints a message on solve that clp was chosen by default

Option 2:

using JuMP, Clp

m = Model(solver = ClpSolver())
@variable(m, 0 <= x <= 2)
@variable(m, 0 <= y <= 30)
@objective(m, Max, 5x + 3y)
@constraint(m, 1x + 5y <= 3.0)
status = solve(m)
# no messages printed by JuMP

This is the most similar to "production" code but feels a bit ugly to me. It's potentially annoying that you'd have to install Clp or modify the code (digging through it to find all uses of ClpSolver()) to run examples.

Option 3:

using JuMP

JuMP.loaddefaultsolvers() # Loads available solvers

m = Model()
@variable(m, 0 <= x <= 2)
@variable(m, 0 <= y <= 30)
@objective(m, Max, 5x + 3y)
@constraint(m, 1x + 5y <= 3.0)
status = solve(m)
# If an LP solver is installed, prints a message on solve that one was chosen by default. If not, prints an error with suggested solvers.

Should we even have a loaddefaultsolvers method?

Relevant prior discussion at https://github.com/JuliaOpt/MathProgBase.jl/pull/150.

CC @odow @chriscoey @blegat @chkwon

[mlubin]

My objections to option 2 could be alleviated by following a convention of putting

using JuMP, Clp
solver = ClpSolver()

at the top of each file so that it's easier to run the examples with a new solver.

[chriscoey]

I don't see the problem with Option 2. It's the most straightforward and transparent. Different solvers may have slightly different behaviors and I feel users should be aware of and actively make a decision about which solver to use. It's more verbose perhaps, but IMO less "ugly" than the other options.

The suggestion you just posted, Miles, is pretty much exactly what I was about to suggest.

[mlubin]

Should we even allow not specifying a solver then?

[chriscoey]

I think it creates headaches (as we have seen with Pajarito a little). Keep things simple, always specify a solver.

[chkwon]

I also vote for Option 2.

Disallowing a non-solver model creation could potentially affect JuMP extensions. In my Complementarity.jl package, I create a JuMP Model() without actually solving it. Can we do solver=nothing?

[mlubin]

@chkwon, good point. We can make it an error at solve time instead of when calling Model().

[joehuchette]

Personally I think the INFO message in option 1 will help obviate a lot of the potential confusion. Option 2 would be fine as well, but there's an elegance to Option 1 for teaching/workshops that I would miss.

[mlubin]

Nobody besides me and @joehuchette seems to favor the magical option of picking a default solver based on which packages are loaded, so seems like option 2 is the winner.

[odow]

What about

using JuMP, Clp
JuMP.adddefaultsolver!(ClpSolver())
# or
# JuMP.setdefaultsolvers!(ClpSolver())

m = Model()
@variable(m, x2)
@constraint(m, x == 1)
solve(m) # defaults to Clp

m2 = Model()
@variable(m2, x)
@nlconstraint(m2, exp(x) == 1)
solve(m2) # errors "No Non-Linear solver loaded as default"

using JuMP, Clp, Ipopt
JuMP.adddefaultsolver!(ClpSolver())
JuMP.adddefaultsolver!(IpoptSolver())
# or
# JuMP.setdefaultsolvers!(ClpSolver(), IpoptSolver())

m = Model()
@variable(m, x2)
@constraint(m, x == 1)
solve(m) # defaults to Clp

m2 = Model()
@variable(m2, x)
@nlconstraint(m2, exp(x) == 1)
solve(m2) # defaults to Ipopt

[mlubin]

@odow, that's problematic because we don't have a great way to query the capabilities of a solver given a solver object.

[mlubin]

We could have a global default but that doesn't seem much better than putting solver = ... at the top of your file and using Model(solver=solver) everywhere.

[odow]

Well why not have a way to query solver capability? Ref https://github.com/JuliaOpt/JuMP.jl/issues/83

MathProgBase.supports(solver::AbstractMathProgBaseSolver)
# returns a subset of [:linear, :integer, :quadratic, :conic, :smooth, :nonconvex, :callbacks]

[mlubin]

If you'd like to lead that effort, I'd support it :)

[mlubin]

Implementation at https://github.com/JuliaOpt/JuMP.jl/pull/970