Merge branch 'master' into expr/unary

Author: Joao-Dionisio

CHANGELOG.md

@@ -3,11 +3,13 @@
 ## Unreleased
 ### Added
 - `Expr` and `GenExpr` support NumPy unary functions (`np.sin`, `np.cos`, `np.sqrt`, `np.exp`, `np.log`, `np.absolute`)
+- Added `getBase()` and `setBase()` methods to `LP` class for getting/setting basis status
 - Added `getMemUsed()`, `getMemTotal()`, and `getMemExternEstim()` methods
 ### Fixed
 - Removed `Py_INCREF`/`Py_DECREF` on `Model` in `catchEvent`/`dropEvent` that caused memory leak for imbalanced usage
-- Used getIndex() instead of ptr() for sorting nonlinear expression terms to avoid nondeterministic behavior 
+- Used `getIndex()` instead of `ptr()` for sorting nonlinear expression terms to avoid nondeterministic behavior 
 ### Changed
+- Speed up `constant * Expr` via C-level API
 ### Removed
 - Removed outdated warning about Make build system incompatibility
 

src/pyscipopt/__init__.py

@@ -28,6 +28,7 @@
 from pyscipopt.scip      import LP as LP
 from pyscipopt.scip      import IISfinder as IISfinder 
 from pyscipopt.scip      import PY_SCIP_LPPARAM as SCIP_LPPARAM
+from pyscipopt.scip      import PY_SCIP_BASESTAT as SCIP_BASESTAT
 from pyscipopt.scip      import readStatistics as readStatistics
 from pyscipopt.scip      import Expr as Expr
 from pyscipopt.scip      import MatrixExpr as MatrixExpr

src/pyscipopt/expr.pxi

@@ -337,22 +337,24 @@ cdef class Expr(ExprLike):
         cdef PyObject *v2_ptr = NULL
         cdef PyObject *old_v_ptr = NULL
         cdef Term child
-        cdef double prod_v
+        cdef double coef
 
         if _is_number(other):
-            f = float(other)
-            return Expr({v:f*c for v,c in self.terms.items()})
+            coef = float(other)
+            while PyDict_Next(self.terms, &pos1, &k1_ptr, &v1_ptr):
+                res[<Term>k1_ptr] = <double>(<object>v1_ptr) * coef
+            return Expr(res)
 
         elif isinstance(other, Expr):
             while PyDict_Next(self.terms, &pos1, &k1_ptr, &v1_ptr):
                 pos2 = <Py_ssize_t>0
                 while PyDict_Next(other.terms, &pos2, &k2_ptr, &v2_ptr):
                     child = (<Term>k1_ptr) * (<Term>k2_ptr)
-                    prod_v = (<double>(<object>v1_ptr)) * (<double>(<object>v2_ptr))
+                    coef = (<double>(<object>v1_ptr)) * (<double>(<object>v2_ptr))
                     if (old_v_ptr := PyDict_GetItem(res, child)) != NULL:
-                        res[child] = <double>(<object>old_v_ptr) + prod_v
+                        res[child] = <double>(<object>old_v_ptr) + coef
                     else:
-                        res[child] = prod_v
+                        res[child] = coef
             return Expr(res)
 
         elif isinstance(other, GenExpr):

src/pyscipopt/lp.pxi

@@ -531,6 +531,66 @@ cdef class LP:
 
         return binds
 
+    def getBase(self):
+        """Returns the basis status of columns and rows.
+
+        Status values are defined in SCIP_BASESTAT: LOWER, BASIC, UPPER, ZERO.
+
+        Returns
+        -------
+        tuple of (list of int, list of int)
+            Column basis statuses and row basis statuses.
+
+        """
+        cdef int ncols = self.ncols()
+        cdef int nrows = self.nrows()
+        cdef int* c_cstat = <int*> malloc(ncols * sizeof(int))
+        cdef int* c_rstat = <int*> malloc(nrows * sizeof(int))
+        cdef int i
+
+        PY_SCIP_CALL(SCIPlpiGetBase(self.lpi, c_cstat, c_rstat))
+
+        cstat = [c_cstat[i] for i in range(ncols)]
+        rstat = [c_rstat[i] for i in range(nrows)]
+
+        free(c_rstat)
+        free(c_cstat)
+
+        return cstat, rstat
+
+    def setBase(self, cstat, rstat):
+        """Sets the basis status of columns and rows.
+
+        Status values are defined in SCIP_BASESTAT: LOWER, BASIC, UPPER, ZERO.
+
+        Parameters
+        ----------
+        cstat : list of int
+            Column basis statuses (length must equal ncols).
+        rstat : list of int
+            Row basis statuses (length must equal nrows).
+
+        """
+        cdef int ncols = self.ncols()
+        cdef int nrows = self.nrows()
+        if len(cstat) != ncols:
+            raise ValueError(f"cstat has length {len(cstat)}, expected {ncols}")
+        if len(rstat) != nrows:
+            raise ValueError(f"rstat has length {len(rstat)}, expected {nrows}")
+        cdef int* c_cstat = <int*> malloc(ncols * sizeof(int))
+        cdef int* c_rstat = <int*> malloc(nrows * sizeof(int))
+        cdef int i
+
+        for i in range(ncols):
+            c_cstat[i] = cstat[i]
+        for i in range(nrows):
+            c_rstat[i] = rstat[i]
+
+        PY_SCIP_CALL(SCIPlpiSetBase(self.lpi, c_cstat, c_rstat))
+
+        free(c_rstat)
+        free(c_cstat)
+
     # Parameter Methods
 
     def setIntParam(self, param, value):

src/pyscipopt/scip.pxd

@@ -1532,6 +1532,8 @@ cdef extern from "scip/scip.h":
     SCIP_RETCODE SCIPlpiGetPrimalRay(SCIP_LPI* lpi, SCIP_Real* ray)
     SCIP_RETCODE SCIPlpiGetDualfarkas(SCIP_LPI* lpi, SCIP_Real* dualfarkas)
     SCIP_RETCODE SCIPlpiGetBasisInd(SCIP_LPI* lpi, int* bind)
+    SCIP_RETCODE SCIPlpiGetBase(SCIP_LPI* lpi, int* cstat, int* rstat)
+    SCIP_RETCODE SCIPlpiSetBase(SCIP_LPI* lpi, const int* cstat, const int* rstat)
     SCIP_RETCODE SCIPlpiGetRealSolQuality(SCIP_LPI* lpi, SCIP_LPSOLQUALITY qualityindicator, SCIP_Real* quality)
     SCIP_RETCODE SCIPlpiGetIntpar(SCIP_LPI* lpi, SCIP_LPPARAM type, int* ival)
     SCIP_RETCODE SCIPlpiGetRealpar(SCIP_LPI* lpi, SCIP_LPPARAM type, SCIP_Real* dval)

src/pyscipopt/scip.pxi

@@ -118,6 +118,12 @@ cdef class PY_SCIP_LPPARAM:
     POLISHING      = SCIP_LPPAR_POLISHING
     REFACTOR       = SCIP_LPPAR_REFACTOR
 
+cdef class PY_SCIP_BASESTAT:
+    LOWER = SCIP_BASESTAT_LOWER
+    BASIC = SCIP_BASESTAT_BASIC
+    UPPER = SCIP_BASESTAT_UPPER
+    ZERO  = SCIP_BASESTAT_ZERO
+
 cdef class PY_SCIP_PARAMEMPHASIS:
     DEFAULT      = SCIP_PARAMEMPHASIS_DEFAULT
     CPSOLVER     = SCIP_PARAMEMPHASIS_CPSOLVER

src/pyscipopt/scip.pyi

@@ -481,6 +481,7 @@ class LP:
     def delCols(self, firstcol: Incomplete, lastcol: Incomplete) -> Incomplete: ...
     def delRows(self, firstrow: Incomplete, lastrow: Incomplete) -> Incomplete: ...
     def getActivity(self) -> Incomplete: ...
+    def getBase(self) -> Incomplete: ...
     def getBasisInds(self) -> Incomplete: ...
     def getBounds(
         self, firstcol: Incomplete = ..., lastcol: Incomplete = ...
@@ -505,6 +506,7 @@ class LP:
     def ncols(self) -> Incomplete: ...
     def nrows(self) -> Incomplete: ...
     def readLP(self, filename: Incomplete) -> Incomplete: ...
+    def setBase(self, cstat: Incomplete, rstat: Incomplete) -> Incomplete: ...
     def setIntParam(self, param: Incomplete, value: Incomplete) -> Incomplete: ...
     def setRealParam(self, param: Incomplete, value: Incomplete) -> Incomplete: ...
     def solve(self, dual: Incomplete = ...) -> Incomplete: ...
@@ -1815,6 +1817,13 @@ class PY_SCIP_LOCKTYPE:
     MODEL: ClassVar[int] = ...
     def __init__(self) -> None: ...
 
+class PY_SCIP_BASESTAT:
+    LOWER: ClassVar[int] = ...
+    BASIC: ClassVar[int] = ...
+    UPPER: ClassVar[int] = ...
+    ZERO: ClassVar[int] = ...
+    def __init__(self) -> None: ...
+
 class PY_SCIP_LPPARAM:
     BARRIERCONVTOL: ClassVar[int] = ...
     CONDITIONLIMIT: ClassVar[int] = ...

tests/test_expr.py

@@ -278,6 +278,11 @@ def test_mul():
     x = m.addVar(name="x")
     y = m.addVar(name="y")
 
+    # test Expr * number
+    assert str((x + y) * 2.0) == "Expr({Term(x): 2.0, Term(y): 2.0})"
+    assert str(2.0 * (x + y)) == "Expr({Term(x): 2.0, Term(y): 2.0})"
+
+    # test Expr * Expr
     assert str(Expr({CONST: 1.0}) * x) == "Expr({Term(x): 1.0})"
     assert str(y * Expr({CONST: -1.0})) == "Expr({Term(y): -1.0})"
     assert str((x - x) * y) == "Expr({Term(x, y): 0.0})"

tests/test_lp.py

@@ -1,5 +1,6 @@
 from pyscipopt import LP
 from pyscipopt import SCIP_LPPARAM
+from pyscipopt import SCIP_BASESTAT
 
 def test_lp():
     # create LP instance, minimizing by default
@@ -89,3 +90,25 @@ def test_lp():
 
     assert round(myLP.getObjVal() == solval)
     assert round(5.0 == solval)
+
+    # test basis get/set
+    binds = myLP.getBasisInds()
+    assert len(binds) == myLP.nrows()
+
+    cstat, rstat = myLP.getBase()
+    assert len(cstat) == myLP.ncols()
+    assert len(rstat) == myLP.nrows()
+    assert all(s in (SCIP_BASESTAT.LOWER, SCIP_BASESTAT.BASIC,
+                     SCIP_BASESTAT.UPPER, SCIP_BASESTAT.ZERO) for s in cstat)
+    assert all(s in (SCIP_BASESTAT.LOWER, SCIP_BASESTAT.BASIC,
+                     SCIP_BASESTAT.UPPER) for s in rstat)
+
+    # set the same basis back and re-solve
+    myLP.setBase(cstat, rstat)
+    solval2 = myLP.solve()
+    assert round(solval2, 10) == round(solval, 10)
+
+    # verify basis is preserved after set
+    cstat2, rstat2 = myLP.getBase()
+    assert cstat2 == cstat
+    assert rstat2 == rstat