Modeling min() in JuMP at scale: SOS1 stops converging, Big-M is numerically unstable

Specific Domains Optimization (Mathematical)
1

I’m implementing a formula engine in JuMP where users can define expressions like:

z = min(x, y)

and I translate these automatically into JuMP constraints.

I have tried two formulations:

1. SOS1 formulation

function minOf2(model, vars...)
    z = @variable(model)
    δ = @variable(model, [1:length(vars)])
    @constraint(model, δ .>= 0)
    for (i, v) in enumerate(vars)
        @constraint(model, z <= v)
        @constraint(model, z + δ[i] == v)
    end
    @constraint(model, δ in SOS1())
    return z
end

This works for small models, but once the number of min() expressions increases (and therefore SOS1 sets grow), SCIP fails to converge or stalls extremely early.

2. Big-M formulation

I also tried the classic binary + Big-M:

z ≤ v[i]
z ≥ v[i] - M * (1 - b[i])

This works only when M is tight.
In my application, user-provided values don’t always have good bounds, and large M values quickly introduce numerical instability in SCIP.

My question

Given that:

  • SOS1 scaling is poor in SCIP for many min/max operators
  • Big-M is unstable when M is not tight
  • I cannot reliably compute tight bounds because user formulas are arbitrary

What is the recommended way in JuMP/SCIP to model min() reliably in large models?

Are there alternative formulations that scale better?

Any practical advice for building a formula engine on top of JuMP would be greatly appreciated.

2

Welcome to Discourse,

but I think we lack motivation to implement this:

  • First, JuMP already has some high-level interfaces like SOS1 etc. And users aren’t well-motivated to add further wrappers.
  • Second, your constraint z = min(x, y) is nonconvex intrinsically.

You are mentioning large models, so the answer is simply:

  • There is no such a universally reliable way. You have to explicitly exploit specific structures of your problem. e.g. to decide tight value of big-M’s.

Users should try their best to build tight models themselves.

If so, the user should stop using a blind big-M method.

My overall opinion is: in real-life large scale MILP,s, we should only formulate our problems using standard MILP language (only Ax <= b constraints, excluding others, e.g. SOS1 constraint is considered vague and blind, z = min(x, y) is considered blind, etc.). And we should try our best to make the polyhedral tight. e.g. x <= 1/2, x Int should be tighten to x <= 0, x Int.

One of the most advisable criteria in mathematical programming is:

  • There is no foolproof method (e.g. you are asking for a foolproof interface here, which is hopeless).

BTW, Gurobi have some off-the-shelf foolproof constraints APIs

including the nonconvex z = min(x, y) you are mentioning here. But whether it is effective under large-scale problems—you can do some tests and give us feedback. (You can call these “general constraints” APIs from Gurobi.jl)