How to create a vector function in sympy with actual vectors, that can be transformed between coordinate systems?

[Andreas Schuldei]

I saw this https://stackoverflow.com/questions/46993819/how-to-create-a-vector-function-in-sympy which uses Matrix() as a workaround to create a vector. The author says, that it can not be transformed between coordinate systems, like real vectors, though.

I need to transform my input and output vector from one coordinate system to another (and back). How are vector functions done in that case? My function is simple:

def B_el(r_vec, I):
mu_0 = 4 * np.pi * 1e-7
a1 = -0.05
a2 = 0.0
C = mu_0 * I / np.pi
r1 = r_vec.i
r2 = r_vec.j
u = (2 * C * r2 * (a1 * r1 + a2 * r2)) / np.square(r1*r1 + r2*r2) - (C * a2) / (r1*r1 + r2*r2)
v = (C * a1) / (r1*r1 + r2*r2) - (2 * C * r1 * (a1 * r1 + a2 * r2)) / np.square(r1*r1 + r2*r2)
return Matrix([u, v, 0])

[Alan Bromborsky]

You might want to look at this link -

https://galgebra.readthedocs.io/en/latest/

Also if you could show me symbolically (not code) what you are doing perhaps I could give you an example of how to do it in galgebra.

--
You received this message because you are subscribed to the Google Groups "sympy" group.
To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/8db840a8-20f9-47ac-a1c3-11b6658f00bfn%40googlegroups.com.

[Andreas Schuldei]

I am putting together the components of a vector field (a magnetic field, caused by a current in several conductors) in cartesian coordinates. The field is derived from calculating the rotation of its magnetic vector potential, which can be expressed as

A_z  = -Const * dot(r,a)(dot(r,r)

A =(0,0,A_z

and then the rot(A):
B = curl(A)

gives me (after some simplifications)

u = (2 * C * r2 * (a1 * r1 + a2 * r2)) / np.square(r1*r1 + r2*r2) - (C * a2) / (r1*r1 + r2*r2)
v = (C * a1) / (r1*r1 + r2*r2) - (2 * C * r1 * (a1 * r1 + a2 * r2)) / np.square(r1*r1 + r2*r2)

w = 0

with B(u,v,w) being the vector field I am interested in.

My initial question ("How to create a vector") is mostly a sympy syntactical one. I look for a function like

Sys = CoordSys3D("Sys")

O = Sys.origin

u = (2 * C * r2 * (a1 * r1 + a2 * r2)) / np.square(r1*r1 + r2*r2) - (C * a2) / (r1*r1 + r2*r2)
v = (C * a1) / (r1*r1 + r2*r2) - (2 * C * r1 * (a1 * r1 + a2 * r2)) / np.square(r1*r1 + r2*r2)

w = 0
V = O.vector(u,v,w)

^^^^^^^^^^^^^^^^^^^

where I can specify a vector, relative to a coordinate system, by its components. The O.vector() function would return a vector that can then safely be transformed into other (resting) coordinate systems.

[Oscar Benjamin]

A SymPy Vector is constructed algebraically from the unit vectors i, j
and k of the coordinate system. For a vector field you also use the
coordinate system base scalars x, y and z.

In [11]: from sympy.vector import CoordSys3D, dot

In [12]: O = CoordSys3D('O')

In [13]: r = O.x*O.i + O.y*O.j + O.z*O.k

In [14]: r
Out[14]: (x_O) i_O + (y_O) j_O + (z_O) k_O

In [15]: dot(r, r)
Out[15]:
2 2 2
x_O + y_O + z_O

I often see confusion about this so clearly this is not very intuitive
and maybe it is not very clearly documented:
https://docs.sympy.org/latest/modules/vector/index.html

I think there should be an easier way where you can just use a
"standard" coordinate system and just do:

from sympy.vector.stdcoords import x, y, z, i, j, k
r = x*i + y*j + z*k
print(dot(r, r))

The vector docs should just begin by showing how to do that and then
how to do simple calculations in a single coordinate system.

--
Oscar

> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/9797a914-06a6-4f8e-a99f-5403504862bfn%40googlegroups.com.

[Alan Bromborsky]

I don't know if this would help but you problem cried out for cylindrical coordinates.  I assumed the vectors a and R were from the origin and had no theta component (if not let me know and I will run that case) then this code (snippet) -

def vector_potential_in_cylindrical_coordinates():
    #Print_Function()
    X = (r,th,z) = symbols('r theta z')
    (k,ar,az) = symbols('k a_r a_z')
    gprint('k =',k)
   #g=[1,r**2,1] is the diagonal of the metric tensor for cylindrical coordinates
    c3d = Ga('e_r e_theta e_z',g=[1,r**2,1],coords=X,norm=True)
    (er,eth,ez) = c3d.mv()
    grad = c3d.grad
    R = r*er+z*ez
    gprint(r'\bs{R} =',R)
    a = ar*er+az*ez
    gprint(r'\bs{a} =',a)

    A = (-k*(a|R)*(R|R)*ez).simplify()

    gprint(r'\bs{A} = -k(\bs{a}\cdot \bs{R})\bs{R}^{2}\bs{e}_{z} =',A)

    B = (-c3d.E()*(grad^A)).simplify()

    gprint(r'\bs{B} = \nabla\times \bs{A} = -I*(\nabla\wedge \bs{A}) =',B)

produces this output -

The actual output was a pdf file generated from a LaTeX generated by the galgebra print statement 'gprint.'  I am attaching the actual code and output files.

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/9797a914-06a6-4f8e-a99f-5403504862bfn%40googlegroups.com.

[Alan Bromborsky]

I realized you want 'a' to be a constant vector which my definition of 'a' in cylindrical coordinates is not.  I will develop a solution for a constant 'a'.

[Alan Bromborsky]

Ran the problem in three coordinate systems, rectangular, cylindrical, and non-orthogonal.  In the non-orthogonal the normalized basis vectors are e_a, e_p, and e_z.  e_a is in the direction of the 'a' vector, e_z is in the direction of the vector potential, and e_p is perpendicular to e_a and e_z -

[Alan Bromborsky]

Attached is the code and pdf output for all three cases.

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/9797a914-06a6-4f8e-a99f-5403504862bfn%40googlegroups.com.

[Chris Smith]

This is beyond what I am involved with regularly, but I wonder if this would be good in a tutorial for a pertinent module. Thanks @brombo

/c

[Andreas Schuldei]

I managed to install GAlgebra, and after a bit of tweaking, I managed to get some lovely output on windows. Thank you for doing this!

You are correct, I also did the cylindrical form first, just because this is a cable geometry, and it is kind of obvious to do it that way. But the geometry of my measurement setup pushes me towards doing the actual implementation and math in cartesian coordinates. I do like that i get LaTeX output, which i can use for my paper.

Thank you for your code; it gives me an example of how to deal with vectors.

[Andreas Schuldei]

I am trying to figure out the basics here.

I understand O.x, O.y and O.z have special meaning, first because they magically appear, just by instanciating O, and also because they seem to distinguish themselves from other generic symbols/skalars:

O = CoordSys3D('O')

r = O.x*O.i + O.y*O.j + O.z*O.k

dot(r,r)  # this works
O.x**2 + O.y**2 + O.z**2
from sympy import symbols
a, b, c = symbols("a, b, c")
p = a*O.i +b*O.j + c*O.k
p
a*O.i + b*O.j + c*O.k
dot(p,p)   # this does not work
Traceback (most recent call last):
  File "<input>", line 1, in <module>
  File "C:\Users\Andreas Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\physics\vector\functions.py", line 31, in dot
    raise TypeError('Dot product is between two vectors')
TypeError: Dot product is between two vectors
p = a*O.x*O.i + b*O.y*O.j + c*O.z*O.k
dot(p,p)  #  this does not work, either!
Traceback (most recent call last):
  File "<input>", line 1, in <module>
  File "C:\Users\Andreas Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\physics\vector\functions.py", line 31, in dot
    raise TypeError('Dot product is between two vectors')
TypeError: Dot product is between two vectors

Why is this? What factors are allowed with vector components and unit vectors? 

I feel a little restricted if my only vector component variables allowed are O.x, O.y, and O.z.

[Oscar Benjamin]

I don't see the exception that you showed:

In [9]: O = CoordSys3D('O')
...: r = O.x*O.i + O.y*O.j + O.z*O.k
...: dot(r,r) # this works
...: O.x**2 + O.y**2 + O.z**2
...: from sympy import symbols
...: a, b, c = symbols("a, b, c")
...: p = a*O.i +b*O.j + c*O.k
...: p
...: a*O.i + b*O.j + c*O.k
...: dot(p,p) # this does not work
Out[9]:
2 2 2
a + b + c

--
Oscar

> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/49cefd33-34db-4c97-9010-cc9be4f1f131n%40googlegroups.com.

[Andreas Schuldei]

p.dot(p)
a**2 + b**2 + c**2
dot(p,p)

Traceback (most recent call last):
  File "<input>", line 1, in <module>
  File "C:\Users\Andreas Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\physics\vector\functions.py", line 31, in dot
    raise TypeError('Dot product is between two vectors')
TypeError: Dot product is between two vectors

i am confused. why does this happen?

[Chris Smith]

You might have copied `a**2 + b**2 + c**2` and called them `p`? You can check by typing `type(p)`. If it is an Add then that is what happened somehow. It's always best to show the actual self-verifying code snippet. Here is what I used and everything work for me:

```

>>> from sympy.vector import dot

>>> from sympy.vector.coordsysrect import CoordSys3D

>>> O = CoordSys3D('O')

>>> r = O.x*O.i + O.y*O.j + O.z*O.k

>>> dot(r,r)  # this works

O.x**2 + O.y**2 + O.z**2

>>> from sympy import symbols

>>> a, b, c = symbols("a, b, c")

>>> p = a*O.i +b*O.j + c*O.k

>>> dot(p,p)   # this does not work

a**2 + b**2 + c**2

>>> p.dot(p)

a**2 + b**2 + c**2

```

/c

[Chris Smith]

Not sure if the following might be of any help, but it allows for work in 4 coordinate systems: https://github.com/DocNan/SymFields

/c