4

I have the following 'issue' with sympy at the moment:

I have a symbolic expression like M = matrix([pi*a, sin(1)*b]) which I want to lambdify and pass to a numerical optimizer. The issue is that the optimizer needs the function to input/output numpy arrays of shape (n,) and specifically NOT (n,1).

Now I have been able to achieve this with the following code (MWE):

import numpy as np
import sympy as sp
a, b = sp.symbols('a, b')
M = sp.Matrix([2*a, b])
f_tmp = sp.lambdify([[a,b]], M, 'numpy')
fun   = lambda x: np.reshape( f_tmp(x), (2,))

Now, this is of course extremely ugly, since the reshape needs to be applied every time fun is evaluated (which might be LOTS of times). Is there a way to avoid this problem? The Matrix class is by definition always 2 dimensional. I tried using sympy's MutableDenseNDimArray-class, but they don't work in conjunction with lambdify. (symbolic variables don't get recognized)

Hyperplane
  • 1,422
  • 1
  • 14
  • 28

1 Answers1

2

One way is to convert a matrix to a nested list and take the first row:

fun = sp.lambdify([[a, b]], M.T.tolist()[0], 'numpy')

Now fun([2, 3]) is [4, 3]. This is a Python list, not a NumPy array, but optimizers (at least those in SciPy) should be okay with that.

One can also do

fun = sp.lambdify([[a, b]], np.squeeze(M), 'numpy')

which also returns a list.

In my test the above were equally fast, and faster than the version with a wrapping function (be it np.squeeze or np.reshape): about 6 µs vs 9 µs. It seems the gain is in eliminating one function call.

  • This is really helpful! There only seems to be a bug with the (1,1) case: if I use np.squeeze(M, axis=1) for a (1,1) Matrix, the result has shape () instead of the expected (1,). – Hyperplane Jul 03 '18 at 17:01
  • Also the second options seems to return lists, too and not np.array. (using sympy 1.1.1 and numpy 1.14.3) – Hyperplane Jul 03 '18 at 17:16
  • Right on both counts. The edge case of `(1, 1)` can be fixed with `np.squeeze(np.array(M), axis=1)` but the fact that we don't get a NumPy array anyway suggests the first version is better anyway. –  Jul 03 '18 at 22:00