Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Support nonlinear constraints. #100

Merged
merged 39 commits into from
Mar 2, 2019
Merged

Support nonlinear constraints. #100

merged 39 commits into from
Mar 2, 2019

Conversation

rschwarz
Copy link
Collaborator

@rschwarz rschwarz commented Feb 22, 2019

Fixes #93.

This adds support for nonlinear constraints based on :ExprGraph. The Julia expression is converted to an intermediate array-based representation which is then converted to a SCIP-style expression graph.

The code is mostly translated from CSIP and the old MathProgBase wrapper. However, the handling of powers (with constant exponent) and unary minus is done earlier, to avoid memory leaks or worse.

I hope that the tests that I created use the same expression format that would also be generated by JuMP. I'm tempted to use JuMP for additional tests, but don't want cyclic dependencies.

EDIT: No longer using the intermediate representation, so the CSIP layer is basically removed now.
EDIT: Still waiting for input on how to fix the non-deterministic behavior related to constant expressions.

 - and scip/intervalarith.h for types
 - and nlpi/type_exprinterpret.h for types
 - use Vector{Ptr{SCIP_EXPR}} instead for Julia calls
 - will be converted automatically be ccall
 - don't use separate EXPR_CONST
 - this would lead to memory leak
 - see discussion at scipopt/CSIP#28
 - was resetting bound for binary variables
 - needed to change order (first change type, then reset to old bounds)
 - error only appeared when using SCIP in debug mode
 - don't recurse to create a subexpression for the exponent
 - instead, extract the numeric value and add directly
 - replace :(child) with :(0 - child) in the MOI wrapper
 - can not easily fix this in the ManagedSCIP or risk memory leak.
 - enable (and fix) test
@rschwarz
Copy link
Collaborator Author

Hm, the Travis tests are inconsistent. Locally, the results seem to be non-deterministic:

  • most of the time, all tests pass
  • sometimes, the first nonlinear test fails (status: unbounded)
  • sometimes, the second nonlinear test fails (status: infeasible)

Copy link
Collaborator

@mlubin mlubin left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We just had a discussion how to test MINLP solvers. Basically JuMP is the most practical format to write down small test problems at the moment. You can try running against https://github.com/lanl-ansi/MINLPTests.jl in a separate travis stage to avoid the dependency on JuMP, e.g., https://github.com/JuliaOpt/AmplNLWriter.jl/blob/0aaad316c7b99907f619b4b00337d87ab4a0facf/.travis.yml#L28.

@ccoffrin @odow

src/MOI_wrapper/nonlinear_constraints.jl Outdated Show resolved Hide resolved
src/MOI_wrapper/nonlinear_constraints.jl Outdated Show resolved Hide resolved
@rschwarz
Copy link
Collaborator Author

Output for the unbounded result:

feasible solution found by trivial heuristic after 0.0 seconds, objective value 0.000000e+00
presolving:
(round 1, fast)       0 del vars, 0 del conss, 0 add conss, 0 chg bounds, 0 chg sides, 0 chg coeffs, 1 upgd conss, 0 impls, 0 clqs
problem infeasible or unbounded: variable <t_> with objective -1 can be made infinitely large
presolving (2 rounds: 2 fast, 0 medium, 0 exhaustive):
 1 deleted vars, 0 deleted constraints, 0 added constraints, 0 tightened bounds, 0 added holes, 0 changed sides, 0 changed coefficients
 0 implications, 0 cliques
transformed 1/1 original solutions to the transformed problem space
presolving detected unboundedness
Presolving Time: 0.00

SCIP Status        : problem is solved [unbounded]
Solving Time (sec) : 0.00
Solving Nodes      : 0
Primal Bound       : -1.00000000000000e+20 (1 solutions)
Dual Bound         : -1.00000000000000e+20
Gap                : 0.00 %

@odow
Copy link
Contributor

odow commented Feb 22, 2019

Let me know if you need a hand interfacing MINLPTests. It is definitely the way to go. You should be able to pass all nlp and nlp-cvx tests. There are some issues with the mixed integed ones.

@rschwarz
Copy link
Collaborator Author

rschwarz commented Feb 23, 2019

OK, so I feel there is still some work on the table, in several steps:

  • identify the source of the non-deterministic behavior in tests (current theory: a Julia array that is passed to SCIP in a ccall is GC'd before SCIP has created a copy, maybe here or here).
  • remove the intermediate array-based representation of the expression graph which was mostly useful for the CSIP layer in C.
  • add more comprehensive tests through JuMP with MINLPTests.

@rschwarz
Copy link
Collaborator Author

Unfortunately, adding a GC.preserve block with the arrays used around the SCIP calls did not change the behavior (see gist.)

@codecov
Copy link

codecov bot commented Feb 23, 2019

Codecov Report

Merging #100 into master will increase coverage by 1.28%.
The diff coverage is 34.89%.

Impacted file tree graph

@@            Coverage Diff            @@
##           master    #100      +/-   ##
=========================================
+ Coverage   11.62%   12.9%   +1.28%     
=========================================
  Files         130     134       +4     
  Lines        3967    4240     +273     
=========================================
+ Hits          461     547      +86     
- Misses       3506    3693     +187
Impacted Files Coverage Δ
src/SCIP.jl 100% <ø> (ø) ⬆️
src/wrapper.jl 100% <ø> (ø) ⬆️
src/MOI_wrapper.jl 83.33% <ø> (ø) ⬆️
src/wrapper/pub_expr.jl 1.62% <1.62%> (ø)
src/wrapper/expr_manual.jl 100% <100%> (ø)
src/managed_scip.jl 97.53% <100%> (ø) ⬆️
src/MOI_wrapper/variable.jl 76.52% <62.5%> (-0.76%) ⬇️
src/MOI_wrapper/nonlinear_constraints.jl 88.88% <88.88%> (ø)
src/nonlinear.jl 92.85% <92.85%> (ø)
src/wrapper/CEnum.jl 74.46% <0%> (-2.01%) ⬇️
... and 7 more

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 5507c00...574190f. Read the comment docs.

@rschwarz
Copy link
Collaborator Author

While having SCIP print the problem in its representation, I found that weird values are shown for the constants in the failing cases (rather than 2.0):

CONSTRAINTS
  [nonlinear] <c0>: ((4.94066e-324 * <x0>) + <x1>) <= 1.5;
END
CONSTRAINTS
  [nonlinear] <c0>: ((6.90208e-310 * <x0>) + <x1>) <= 1.5;
END

These look like access to unitialized memory. I checked the values just before the ccall where it was still OK. Also the value is passed as Cdouble, not via a Ptr or in an array.

Just to be sure, I added GC.enable(false) and GC.enable(true) around all the includes in run_tests.jl which did not make any difference, so I'm ruling that out for now...

@rschwarz
Copy link
Collaborator Author

In the other testset, it becomes more apparent that all constant values (inside the expressions) are corrupted:

CONSTRAINTS
  [nonlinear] <c0>: (<x0> + 4.94066e-324) == 2;
  [nonlinear] <c1>: (<x1> - <x2>) == 2;
  [nonlinear] <c2>: ((4.94066e-324 - <x3>) + 6.9532e-310) == 2;
  [nonlinear] <c3>: (<x4> + <x5> + <x6>) == 2;
  [nonlinear] <c4>: (<x7> * <x8> * <x9>) == 2;
  [nonlinear] <c5>: realpower((<x10> + <x11>), 6.91388e-310) == 2;
  [nonlinear] <c6>: (<x12> / <x13>) == 2;
  [nonlinear] <c7>: sqrt(<x14>) == 2;
  [nonlinear] <c8>: exp(<x15>) == 2;
  [nonlinear] <c9>: log(<x16>) == 2;
END

In particular, the second argument of realpower is affected, even though it is not a proper sub-expression of type SCIP_EXPR_CONST.

@rschwarz
Copy link
Collaborator Author

I also found that if I compile with cmake option -DNOBLKBUFMEM=ON, I can not reproduce the failing behavior.

@fserra, do you have an idea what this would imply?
I also ran valgrind on that test, but it does not look like it found anything relevant (gist).

@rschwarz
Copy link
Collaborator Author

I was able to narrow down the problem further:
After creating the SCIP expression for constant values, that value is sometimes corrupted directly afterwards. That is, after the call to SCIPexprCreate, I would get

value = 2.0
expr__ = Base.RefValue{Ptr{Nothing}}(Ptr{Nothing} @0x00000000027f7050)
SCIPexprGetOpReal(expr__[]) = 6.89937396681325e-310
SCIPexprGetOpIndex(expr__[]) = -1682932864
SCIPexprGetOpData(expr__[]) = Ptr{Nothing} @0x00007f019bb07b80

The three methods SCIPexprGet* all get the "same" value, which is a C union.
In SCIPexprCreate, that union opdata is actually stack allocated, but then assigned to a struct member of SCIP_EXPR*, so that shouldn't be a problem.

@fserra
Copy link
Collaborator

fserra commented Feb 25, 2019

does it also go away if you just disable block memory (-DNOBLKMEM)? The idea of block memory is that scip allocates memory itself and the re-uses it. I don't see how this can cause this though

@rschwarz
Copy link
Collaborator Author

rschwarz commented Feb 25, 2019

Yes, adding -DNOBLKMEM is enough to avoid the problem.

I tried to create minimal examples for both the C code and (equivalent) Julia code that are both self-contained: gist.

EDIT: I also asked on discourse as I don't know how to continue here...

@fserra
Copy link
Collaborator

fserra commented Feb 26, 2019

I guess the problem has to be the way julia handles variadic c functions. I implemented a SCIPexprCreate22 which receives a SCIP_Real instead of ... and only creates constant expressions.
Then I modified the julia script so that it calls that function instead and after running 400 times the test it never failed.

@mlubin
Copy link
Collaborator

mlubin commented Feb 26, 2019

The Julia signatures for the call to SCIPexprCreate don't seem to match what the documentation says to do for variadic functions: https://docs.julialang.org/en/v1/manual/calling-c-and-fortran-code/#man-bits-types-1. There should be ... somewhere. (I haven't done this myself.)

@rschwarz
Copy link
Collaborator Author

Yes, that must be it. Changing the signature for ccall does seem to resolve issue for me.
Unfortunately, several operators (SUM, PRODUCT, REALPOWER) require variadic arguments of different types, so I guess we can't get around introducing another layer of C inbetween.

rschwarz added 8 commits March 1, 2019 21:49
Turns out that the expression that JuMP produces is actually this:

  Expr
    head: Symbol ref
    args: Array{Any}((2,))
      1: Symbol x
      2: MathOptInterface.VariableIndex
        value: Int64 2

While I just tried `dump(:(x[MathOptInterface.VariableIndex(2)]))`, which gives
something completely different:

  Expr
    head: Symbol ref
    args: Array{Any}((2,))
      1: Symbol x
      2: Expr
        head: Symbol call
        args: Array{Any}((2,))
          1: Expr
            head: Symbol .
            args: Array{Any}((2,))
              1: Symbol MathOptInterface
              2: QuoteNode
                value: Symbol VariableIndex
          2: Int64 2
Copy link
Contributor

@odow odow left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Any reason to only test a couple of MINLPTests instead of the full set (excluding the ones with trig functions etc)?

test/MINLPTests/run_minlptests.jl Outdated Show resolved Hide resolved
@rschwarz
Copy link
Collaborator Author

rschwarz commented Mar 2, 2019

Any reason to only test a couple of MINLPTests instead of the full set (excluding the ones with trig functions etc)?

Yes: I had to change the target_termination to OPTIMAL or INFEASIBLE, rather than the LOCAL_ variants.

Also, I handpicked the few that are supported, which is a smaller list than the excluded ones:

  • SCIP itself does not support trigonometric functions.
  • SCIP.jl only supports affine objective functions ( I guess it can also do Min -x but not Min x). So with a suitable objective bridge in MOI, I could add many more.

@odow
Copy link
Contributor

odow commented Mar 2, 2019

So with a suitable objective bridge in MOI

It's not very had to add support for SingleVariable objective functions as a temporary stop-gap? It can be removed once there is a bridge.

@rschwarz
Copy link
Collaborator Author

rschwarz commented Mar 2, 2019

It's not very had to add support for SingleVariable objective functions as a temporary stop-gap? It can be removed once there is a bridge.

Well, I meant that with an objective bridge, we could support quadratic and nonlinear objectives. The objective bridge itself will probably require SingleVariabe objectives, so that it makes sense to support it manually.

rschwarz added 6 commits March 2, 2019 21:02
 - operator `\`, is not even supported through JuMP/MOI?
 - previously only allowed one type of constraint per variable
 - now also allow simultaneous LessThan and GreaterThan constraints
 - allow JuMP's sum(...) for singleton sets
@rschwarz
Copy link
Collaborator Author

rschwarz commented Mar 2, 2019

So, the major thing left to do is the C wrapper for SCIPexprCreate.

Actually, I'm not so sure how to do that after all. For SCIP_EXPR_SUM, I will need to call SCIPexprCreate with an unknown number of arguments. The first one will be int nchildren, but then there is a set of SCIP_EXPR*. It's not possible to simply forward variadic arguments in C function calls, without providing another version of the called function (with a different signature).

The only workaround I can figure out now is to pass an array of SCIP_EXPR* and then prepare many different calls, based on the size of the array (say, from 1 to 20). But this would introduce an arbitrary limit on the size of expressions used 😞

UPDATE: I'd rather merge this PR now, and deal with that issue later (#101).

@rschwarz rschwarz merged commit 9ec6257 into master Mar 2, 2019
@rschwarz rschwarz deleted the rs/nonlin_cons branch March 3, 2019 18:02
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

Support constraints with general nonlinear expression.
4 participants