Merge pull request #102 from SCIP-Interfaces/rs/more_nl_ops

Support more operators in nonlinear expressions.

Author: rschwarz

README.md

@@ -30,12 +30,12 @@ Only Linux is tested and officially supported. Contributions to supporting other
 operating systems are welcome.
 
 We recommend using one of the provided installers, e.g.,
-`SCIPOptSuite-6.0.0-Linux.deb` for systems based on Debian. Adding the SCIP.jl
+`SCIPOptSuite-6.0.1-Linux.deb` for systems based on Debian. Adding the SCIP.jl
 package should then work out of the box:
 
     pkg> add SCIP
 
-If you [build SCIP from source](https://scip.zib.de/doc-6.0.0/html/CMAKE.php)
+If you [build SCIP from source](https://scip.zib.de/doc-6.0.1/html/CMAKE.php)
 you should set the environment variable `SCIPOPTDIR` to point the the
 **installation path**. That is, `$SCIPOPTDIR/lib/libscip.so` should exist.
 
@@ -95,10 +95,18 @@ constraints by name (`SingleVariable`-set constraints are not stored as SCIP
 constraints explicitly).
 
 Support for more constraint types (quadratic/SOC, SOS1/2, nonlinear expression)
-is planned, but SCIP itself only supports affine objective functions, so we will
-stick with that. More general objective functions could be implented via a
+is implemented, but SCIP itself only supports affine objective functions, so we
+will stick with that. More general objective functions could be implented via a
 [bridge](https://github.com/JuliaOpt/MathOptInterface.jl/issues/529).
 
+Supported operators in nonlinear expressions are as follows:
+
+- unary: `-`, `sqrt`, `exp`, `log`, `abs`
+- binary: `-`, `/`, `^`, `min`, `max`
+- n-ary: `+`, `*`
+
+In particular, trigonometric functions are not supported.
+
 ## Old Interface Implementation
 
 A previous implementation of SCIP.jl supported

src/nonlinear.jl

@@ -1,14 +1,17 @@
 
 # Mapping from Julia (as given by MOI) to SCIP operators
 const OPMAP = Dict{Symbol, SCIP_ExprOp}(
-    :+ => SCIP_EXPR_SUM,
-    :* => SCIP_EXPR_PRODUCT,
-    :- => SCIP_EXPR_MINUS,
-    :/ => SCIP_EXPR_DIV,
-    :^ => SCIP_EXPR_REALPOWER,
-    :sqrt => SCIP_EXPR_SQRT,
-    :exp => SCIP_EXPR_EXP,
-    :log => SCIP_EXPR_LOG,
+    :+ => SCIP_EXPR_SUM,       # n-ary
+    :* => SCIP_EXPR_PRODUCT,   # n-ary
+    :- => SCIP_EXPR_MINUS,     # unary, binary
+    :/ => SCIP_EXPR_DIV,       # unary
+    :^ => SCIP_EXPR_REALPOWER, # binary (or INTPOWER)
+    :sqrt => SCIP_EXPR_SQRT,   # unary
+    :exp => SCIP_EXPR_EXP,     # unary
+    :log => SCIP_EXPR_LOG,     # unary
+    :abs => SCIP_EXPR_ABS,     # unary
+    :min => SCIP_EXPR_MIN,     # binary
+    :max => SCIP_EXPR_MAX,     # binary
 )
 
 """Subexpressions and variables referenced in an expression tree.
@@ -56,10 +59,13 @@ function push_expr!(nonlin::NonlinExpr, mscip::ManagedSCIP, expr::Expr)
 
             # Exponent (second child) is stored as value.
             @assert isa(expr.args[3], Number)
-            exponent = Cdouble(expr.args[3])
-
-            # Create SCIP expression
-            @SC SCIPexprCreate(SCIPblkmem(mscip), expr__, OPMAP[op], base, exponent)
+            if isa(expr.args[3], Integer)
+                exponent = Cint(expr.args[3])
+                @SC SCIPexprCreate(SCIPblkmem(mscip), expr__, SCIP_EXPR_INTPOWER, base, exponent)
+            else
+                exponent = Cdouble(expr.args[3])
+                @SC SCIPexprCreate(SCIPblkmem(mscip), expr__, SCIP_EXPR_REALPOWER, base, exponent)
+            end
 
         elseif op == :- && num_children == 1
             # Special case: unary version of minus. SCIP only supports binary
@@ -74,7 +80,7 @@ function push_expr!(nonlin::NonlinExpr, mscip::ManagedSCIP, expr::Expr)
             # Finally, add the (binary) minus:
             @SC SCIPexprCreate(SCIPblkmem(mscip), expr__, OPMAP[op], left, right)
 
-        elseif op in [:sqrt, :exp, :log]
+        elseif op in [:sqrt, :exp, :log, :abs]
             # Unary operators
             @assert num_children == 1
 
@@ -84,7 +90,7 @@ function push_expr!(nonlin::NonlinExpr, mscip::ManagedSCIP, expr::Expr)
             # Add this operator on top
             @SC SCIPexprCreate(SCIPblkmem(mscip), expr__, OPMAP[op], child)
 
-        elseif op in [:-, :/]
+        elseif op in [:-, :/, :min, :max]
             # Binary operators
             @assert num_children == 2
 

src/wrapper/expr_manual.jl

@@ -61,9 +61,20 @@ end
 function SCIPexprCreate(blkmem, expr, op, base::Ptr{SCIP_EXPR}, exponent::Cdouble)
     # WARNING: In the actual C function, the last two arguments given here are
     # part of the variadic arguments. But since Julia only supports variadic
-    # arguments of the same type, we specify the first one (nchildren::Cint)
+    # arguments of the same type, we specify the first one (base::Ptr{SCIP_EXPR})
     # directly, and pretend that only the last one is variadic.
     ccall((:SCIPexprCreate, libscip), SCIP_RETCODE,
           (Ptr{BMS_BLKMEM}, Ptr{Ptr{SCIP_EXPR}}, SCIP_EXPROP, Ptr{SCIP_EXPR}, Cdouble...),
           blkmem, expr, op, base, exponent)
 end
+
+# SCIP_EXPR_INTPOWER (mixed binary op)
+function SCIPexprCreate(blkmem, expr, op, base::Ptr{SCIP_EXPR}, exponent::Cint)
+    # WARNING: In the actual C function, the last two arguments given here are
+    # part of the variadic arguments. But since Julia only supports variadic
+    # arguments of the same type, we specify the first one (base::Ptr{SCIP_EXPR})
+    # directly, and pretend that only the last one is variadic.
+    ccall((:SCIPexprCreate, libscip), SCIP_RETCODE,
+          (Ptr{BMS_BLKMEM}, Ptr{Ptr{SCIP_EXPR}}, SCIP_EXPROP, Ptr{SCIP_EXPR}, Cint...),
+          blkmem, expr, op, base, exponent)
+end

test/MOI_nonlinear_exprs.jl

@@ -60,7 +60,7 @@ end
 
     @test MOI.supports(optimizer, MOI.NLPBlock()) == true
 
-    num_vars = 17
+    num_vars = 20
     x = MOI.add_variables(optimizer, num_vars)
     for i in 1:num_vars
         MOI.add_constraint(optimizer, x[i], MOI.Interval(0.1, 10.0))
@@ -74,11 +74,15 @@ end
         :(-x[$(x[4])] + 4.0                     == rhs), # MINUS (unary)
         :(x[$(x[5])] + x[$(x[6])] + x[$(x[7])]  == rhs), # SUM
         :(x[$(x[8])] * x[$(x[9])] * x[$(x[10])] == rhs), # PRODUCT
-        :((x[$(x[11])] + x[$(x[12])])^0.8       == rhs), # REALPOWER
+        :((x[$(x[11])])^3                       == rhs), # INTPOWER
+        :((x[$(x[12])])^0.8                     == rhs), # REALPOWER
         :(x[$(x[13])] / x[$(x[14])]             == rhs), # DIV
         :(sqrt(x[$(x[15])])                     == rhs), # SQRT
         :(exp(x[$(x[16])])                      == rhs), # EXP
         :(log(x[$(x[17])])                      == rhs), # LOG
+        :(abs(x[$(x[18])] - 11)                 == rhs), # ABS
+        :(min(x[$(x[19])], x[$(x[20])]) + 1     == rhs), # MIN
+        :(max(x[$(x[19])], x[$(x[20])]) - 1     == rhs), # MAX
     ]
 
     data = MOI.NLPBlockData(
@@ -101,9 +105,13 @@ end
     @test -sol[4] + 4.0             ≈ rhs  atol=atol rtol=rtol
     @test sol[5] + sol[6] + sol[7]  ≈ rhs  atol=atol rtol=rtol
     @test sol[8] * sol[9] * sol[10] ≈ rhs  atol=atol rtol=rtol
-    @test (sol[11] + sol[12])^0.8   ≈ rhs  atol=atol rtol=rtol
+    @test (sol[11])^3               ≈ rhs  atol=atol rtol=rtol
+    @test (sol[12])^0.8             ≈ rhs  atol=atol rtol=rtol
     @test sol[13] / sol[14]         ≈ rhs  atol=atol rtol=rtol
     @test sqrt(sol[15])             ≈ rhs  atol=atol rtol=rtol
     @test exp(sol[16])              ≈ rhs  atol=atol rtol=rtol
     @test log(sol[17])              ≈ rhs  atol=atol rtol=rtol
+    @test abs(sol[18] - 11)         ≈ rhs  atol=atol rtol=rtol
+    @test min(sol[19], sol[20])     ≈ 1.0  atol=atol rtol=rtol
+    @test max(sol[19], sol[20])     ≈ 3.0  atol=atol rtol=rtol
 end