0

I want to compute a symbolic gradient with sympy, e.g.,

import sympy as sym
x, y, z  = sym.symbols("x y z", real=True)

T = sym.cos(x**2+y**2)

gradT = sym.Matrix([sym.diff(T, x), sym.diff(T,y), sym.diff(T,z)])

Now I would like to create a lamddify function with this expression:

func = lambdify((x,y,z), gradT,'numpy')

To use the function I have:

gradT_exact = func(np.linspace(0,2,100), np.linspace(0,2,100), np.linspace(0,2,100))

and I receive the following error:

<lambdifygenerated-3>:2: VisibleDeprecationWarning: Creating an ndarray from ragged nested sequences (which is a list-or-tuple of lists-or-tuples-or ndarrays with different lengths or shapes) is deprecated. If you meant to do this, you must specify 'dtype=object' when creating the ndarray.
  return (array([[-2*x*sin(x**2 + y**2)], [-2*y*sin(x**2 + y**2)], [0]]))

If I change T to be a function of x,y,z it gives me no problems... Why is it giving warnings when T only depends on x and y and z is set to zero.

Thanks in advance!

Manuel Oliveira
  • 527
  • 5
  • 19
  • This is a bug in sympy: https://github.com/sympy/sympy/issues/21613 – Oscar Benjamin Jul 16 '21 at 00:23
  • @OscarBenjamin, I don't see how `lambdify` can adjust that `[0]` to accurately reflect possible shape(s) of the other terms (for a `x,y,z` in general). – hpaulj Jul 16 '21 at 01:13
  • In general it is impossible for lambdify to know what the user wants when they substitute symbols for arrays into a matrix-valued expression e.g. where should the extra array dimensions go? It would be possible to define something sensible but it probably needs to be a different function from lambdify that would have arguments for controlling how to make it work. – Oscar Benjamin Jul 16 '21 at 10:40

1 Answers1

0

The gradT expression:

In [84]: gradT
Out[84]: 
⎡        ⎛ 2    2⎞⎤
⎢-2⋅x⋅sin⎝x  + y ⎠⎥
⎢                 ⎥
⎢        ⎛ 2    2⎞⎥
⎢-2⋅y⋅sin⎝x  + y ⎠⎥
⎢                 ⎥
⎣        0        ⎦

and its conversion to numpy:

In [87]: print(func.__doc__)
Created with lambdify. Signature:

func(x, y, z)

Expression:

Matrix([[-2*x*sin(x**2 + y**2)], [-2*y*sin(x**2 + y**2)], [0]])

Source code:

def _lambdifygenerated(x, y, z):
    return (array([[-2*x*sin(x**2 + y**2)], [-2*y*sin(x**2 + y**2)], [0]]))

If x and y are arrays, then 2 terms will reflect their dimension(s), but the last is [0]. That's why you get the ragged warning.

lambdify does a rather simple lexical translation. It does not implement any deep understanding of numpy arrays. At some level it's your responsibility to check that the numpy code looks reasonable.

with scalar inputs:

In [88]: func(1,2,3)
Out[88]: 
array([[1.91784855],
       [3.8356971 ],
       [0.        ]])

but if one input is an array:

In [90]: func(np.array([1,2]),2,3)
<lambdifygenerated-1>:2: VisibleDeprecationWarning: Creating an ndarray from ragged nested sequences (which is a list-or-tuple of lists-or-tuples-or ndarrays with different lengths or shapes) is deprecated. If you meant to do this, you must specify 'dtype=object' when creating the ndarray.
  return (array([[-2*x*sin(x**2 + y**2)], [-2*y*sin(x**2 + y**2)], [0]]))
Out[90]: 
array([[array([ 1.91784855, -3.95743299])],
       [array([ 3.8356971 , -3.95743299])],
       [0]], dtype=object)

The result is object dtype containing 2 arrays, plus that [0] list.

To avoid this problem, the lambdify would have to produce a function like:

In [95]: def f(x,y,z):
    ...:     temp = 0*x*y
    ...:     return np.array([-2*x*np.sin(x**2 + y**2), -2*y*np.sin(x**2 + y**2)
    ...: , temp])

where temp is designed to give 0 value, but with a shape that reflects the broadcasted operations on x and y in the other terms. I think that's asking too much of lambdify.

In [96]: 

In [96]: f(np.array([1,2]),2,3)
Out[96]: 
array([[ 1.91784855, -3.95743299],
       [ 3.8356971 , -3.95743299],
       [ 0.        ,  0.        ]])
hpaulj
  • 221,503
  • 14
  • 230
  • 353