Correct way to express SOSII constraints

I have an application in which I would like to use SOS II constraints to model a piecewise linear approximation to non-linear function. In particular I have a log(z) term in my objective which I would like to approximate with a piecewise linear function of z, over the interval [zl,zu].

I first specify a set of k breakpoints over [zl,zu]

k=10
zl=10
zu=1000
gridpoints=collect(linspace(log(zl),log(zu),k))
breakpoints=map(exp,gridpoints)
fpoints=map(log,breakpoints)

fpoints is the value of the log function evaluated at the breakpoints. Suppose I would like to maximize log(z). I can set up a model in jump as follows

model = Model(solver=GurobiSolver())
@variable(model, z)
@constraint(model, z<=zu)
@constraint(model, zl<=z)
@variable(model, x[1:k])
addSOS2(model,breakpoints.*x)
#addSOS2(model,x)
@constraint(model,sum(x[i] for i=1:k) == 1) #Convexity Row
@constraint(model, sum(breakpoints[i]*x[i] for i=1:k)==z) #Reference Row
@objective(model,Max,sum(fpoints[i]*x[i] for i=1:k)) #Function Row
solve(model)
getvalue(z)

My question is whether to use addSOS2(model,x) or addSOS2(model,breakpoints.*x)? It is true that x are the variables which need to specified as a specially ordered set of type II, but the documentation for jump specifies that it should be specified with unique weights. In the original paper by Beale and Tomlin the weights are used to determine the branching, and in the paper it says that the weights should be taken as the breakpoints in the reference row.

So which of addSOS2(model,x) or addSOS2(model,breakpoints.*x) is better?

Use addSOS2(model,breakpoints.*x), since the weights do have to be unique.

Alternatively, you could use my PiecewiseLinearOpt.jl package, and compare your solver’s native SOS2 support against a number of MIP formulations.

2 Likes

For the case of maximizing a concave function like log, there’s no need for a MIP formulation at all. You could try using a mixed-integer convex solver like Pajarito or adding tangent hyperplanes manually.

1 Like

Thanks for your responses @joehuchette and @miles.lubin. The logarithmic term is part of a bigger problem I’m working on. In particular I am trying to minimize a MIQCP program, which has a log term in the objective. My approach is to build a piecewise linear approximation to the logarithm so that I can still solve it as a MIQCP instead of reaching for more general MINLP methods. Is there a different approach @miles.lubin?

Yes, you can approximate \log(x) as a piecewise linear concave function f(x). The constraint t \le f(x) can be modeled directly with linear programming, no need for SOS. In your objective you’ll be maximizing t so equality will hold in the optimal solution.

@miles.lubin this is a great suggestion for the problem I presented. Unfortunately when I looked at the more complicated problem that I am actually working on, its a maximization problem containing a -log(x) term in the objective. If I understand your suggestion correctly, t < \log(x) can be very well approximated using a set of linear constraints. For my problem of Maximizing t s.t. t\leq -\log(x), I think I cannot build a good approximation to the constraint using linear constraints as -\log(x) is convex. This is what led me to using a piecewise linear approximation using SOS or binary variables. Do you think it is the correct place to end up or are you able to suggest any other tricks?

In that case it’s pretty reasonable to use SOS or binary variables.

Thanks for confirming :slight_smile:

You could also try COUENNE. I’ve used it through https://github.com/JuliaOpt/CoinOptServices.jl, and my limited experience with the solver was actually pretty good. You pass it the original problem (nonlinearities and all), and internally it automatically sets up a linear relaxation. The cool thing is that it adaptively tightens the relaxation inside a spatial branch-and-bound framework, and can actually solve the original problem to global optimality (insert computational complexity caveat here).