Using NumPy in the Xpress Python interface

NumPy arrays can be used as usual when creating variables, functions (linear and quadratic) of variables, and constraints. All functions described in this manual that take lists or tuples as arguments can take array's, i.e., NumPy array objects, as well, as in the following example:

import numpy as np
import xpress as xp
N = 20
S = range(N)
x = np.array([xp.var()                  for i in S], dtype=xp.npvar)
y = np.array([xp.var(vartype=xp.binary) for i in S], dtype=xp.npvar)
constr1 = x <= y
p = xp.problem()
p.addVariable(x, y)
p.addConstraint(constr1)

The above script imports both NumPy and the Xpress Python interface, then declares two arrays of variables and creates the set of constraints x_i ≤ y_i for all i in the set S.

The NumPy arrays must have the attribute dtype equal to xpress.npvar (abbreviated to xp.npvar here) in order to use the matricial/vectorial form of the comparison (<=, =, >=), arithmetic (+, -, *, /, **), and logic (&, |) operators.

NumPy allows for multiarrays with one or more 0-based indices. Given that declaring a NumPy multiarray of variables can result in a long line of code, the xpress.vars function in its simplest usage returns a NumPy array of variables with one or more indices. Consider the following three array declarations:

import numpy as np
import xpress as xp
x = np.array([xp.var(name='v({0})'.format(i)) for i in range(20)], dtype=xp.npvar).reshape(5,4)
y = np.array([xp.var(vartype=xp.binary) for i in range(27)], dtype=xp.npvar).reshape(3,3,3)
z = np.array([xp.var(lb=-1, ub=1) for i in range(1000)], dtype=xp.npvar)

These can be written equivalently in the compact form as

import numpy as np
import xpress as xp
x = xp.vars(5, 4, name='v')
y = xp.vars(3, 3, 3, vartype=xp.binary)
z = xp.vars(1000, lb=-1, ub=1)

The only side effect is that the assigned names change. In order to preserve the naming convention of the Xpress library, one can specify the parameter setting name='' in the call to xp.vars. This also makes the creation of large arrays of variables much faster. We use this shorter notation in the remainder of this chapter.

The main advantage of using NumPy operations is the ability to replicate them on each element of an array, taking into account all broadcasting features. For example, the following script ``broadcasts'' the right-hand side 1 to all elements of the array, thus creating the set of constraints x_i + y_i ≤ 1 for all i in the set S.

constr2 = x + y <= 1

All these operations can be carried out on arrays of any number of dimensions, and can be aggregated at any level. The following example shows two three-dimensional array of variables involved in two systems of constraints: the first has two variables per each of the 200 constraints, while the second has 10 constraints and 20 variables in each constraint.

z = xp.vars(4, 5, 10)
t = xp.vars(4, 5, 10, vartype=xp.binary)
p.addVariable(z,t)
p.addConstraint(z**2 <= 1 + t)
p.addConstraint(xp.Sum(z[i, j, k] for i in range(4) for j in range(5)) <= 4
                 for k in range(10))

Finally, a note on sums of multi-dimensional NumPy arrays: in keeping with the way NumPy arrays are handled, the sum of a bi-dimensional array results in a one-dimensional array with the xpress.Sum operator. The result of such a sum is exemplified by the following code:

>>> a = np.array([[1, 2, 3], [4, 5, 6]])
>>> a
array([[1, 2, 3],
       [4, 5, 6]])
>>> sum(a)
array([5, 7, 9])

For the casual NumPy user, suffice it to say that the sum is done on the first dimension. Similarly, when creating a NumPy array of dimensions k of expressions, xpress.Sum returns a (k-1)-array resulting from the sum across the first dimension.

It is important to note the following: NumPy does not use the __iadd__ operator when computing these sums, but rather the __add__ operator. For reasons discussed above and in the entry regarding the xpress.Sum operator, this can have a huge impact on performance. Consider the following example:

m,n = 1000,10
a = np.random.random((m,n))
x = xp.vars(m, n)
sum_0d = xp.Sum(xp.Sum(a*x))
sum_1d = xp.Sum(a*x)

The above example has a poor performance, and it is advised to avoid using xpress.Sum as such on a multi-dimensional array. If a scalar sum of all elements of the array is sought, such as sum(sum(a)) for the numerical array above, we strongly advise to flatten the array first, and run instead xpress.Sum(b.flatten()) if b is a multiarray of expressions. The multiarray has dtype equal to xpress.npexpr in order to be used for array operations. If only one pass is required, then it is better to explicitly create a vector whose elements are defined with a call to xpress.Sum:

prod = a*x
sum_0d = xp.Sum(prod.flatten())
sum_1d = np.array([xp.Sum(prod[i,j] for i in range(m)) for i in range(n)], dtype=xp.npexpr)