I have a curl of a magnetic field that I calculated that I want to graph, using matplotlib but I end up with this AttributeError. Is there something I'm doing wrong?

I tried changing the variables from using the reference frame R[0] to just x, but also incorporating the vectors using the reference frame of cartesian coordinates.

import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.mplot3d import axes3d
import sympy 
from sympy import Symbol, diff, Array, sin, cos
from sympy import init_printing
from sympy.physics.vector import curl, ReferenceFrame
init_printing()


alpha = Symbol('\u03B1')
B0 = Symbol('B0')
R = ReferenceFrame('R')

# In order to get the curl, setting it as a vector
V = 0*R.x + B0*sin(alpha*R[0])*R.y + B0*cos(alpha*R[0])*R.z

# Calculates the curl of the vector, which results in a vector of 
# alpha*B, which is correct 
C = curl(V,R)
print('The curl of B is:',0*R.x + C,)

fig = plt.figure()
ax = fig.gca(projection = '3d')

x,y,z = np.meshgrid(np.arange(0.01, 1, 0.2), 
                    np.arange(0.01, 1, 0.2),
                    np.arange(0.01, 1, 0.2))
u = 0*R.x
v = B0*alpha*sin(alpha*x)*R.y
w = B0*alpha*cos(alpha*x)*R.z

ax.quiver(x, y, z, u, v, w, length = 0.1)
plt.show()

I'm expecting a 3D graph showing the vector field; below is the traceback of the error.

AttributeError                            Traceback (most recent 
call last)
<ipython-input-12-5974c739f1f4> in <module>
     10                     np.arange(0.01, 1, 0.2))
     11 u = 0*R.x
---> 12 v = B0*alpha*sin(alpha*x)*R.y
     13 w = B0*alpha*cos(alpha*x)*R.z
     14 

~/anaconda3/lib/python3.7/site-packages/sympy/core/function.py in 
__new__(cls, *args, **options)
    440 
    441         evaluate = options.get('evaluate', 
global_evaluate[0])
--> 442         result = super(Function, cls).__new__(cls, *args, 
**options)
    443         if evaluate and isinstance(result, cls) and 
 result.args:
    444             pr2 = min(cls._should_evalf(a) for a in 
result.args)

~/anaconda3/lib/python3.7/site-packages/sympy/core/function.py in 
__new__(cls, *args, **options)
    249 
    250         if evaluate:
--> 251             evaluated = cls.eval(*args)
    252             if evaluated is not None:
    253                 return evaluated

~/anaconda3/lib/python3.7/site- 
packages/sympy/functions/elementary/trigonometric.py in eval(cls, 
arg)
    293             return arg._eval_func(cls)
    294 
--> 295         if arg.could_extract_minus_sign():
    296             return -cls(-arg)
    297 

AttributeError: 'ImmutableDenseNDimArray' object has no attribute 
'could_extract_minus_sign'

1 Answers

0
Patol75 On

Having a look at this answer, the code below should provide a step forward.

import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.mplot3d import axes3d  # noqa
from sympy import cos, sin, Symbol
from sympy import init_printing
from sympy.physics.vector import curl, ReferenceFrame

init_printing()
alpha = Symbol('\u03B1')
B0 = Symbol('B0')
R = ReferenceFrame('R')
# In order to get the curl, setting it as a vector
V = 0 * R.x + B0 * sin(alpha * R[0]) * R.y + B0 * cos(alpha * R[0]) * R.z
# Calculates the curl of the vector, which results in a vector of alpha * B,
# which is correct
C = curl(V, R)
print('The curl of B is:', C)
fig = plt.figure()
ax = fig.gca(projection='3d')
x, y, z = np.meshgrid(np.arange(0.01, 1, 0.2),
                      np.arange(0.01, 1, 0.2),
                      np.arange(0.01, 1, 0.2))
B0 = 0.2
alpha = 5
u = 0
v = B0 * alpha * np.sin(alpha * x)
w = B0 * alpha * np.cos(alpha * x)
ax.quiver(x, y, z, u, v, w, length=0.1)
plt.show()

Additionally, I guess there should be a way to calculate u, v and w like this C.dot(R.y).subs([(B0, 1), (alpha, 1), (R[0], x)]), but I am not familiar enough with sympy to figure out why this returns a Sympy expression. By the way, why don't you use numpy ?