Idiomatic way to evaluate a scalar field at a given vector value?

[Ryan D Johnson]

I've been using sympy in jupyter notebooks as a simple free symbolic algebra system while helping my kids review calculus and physics in anticipation of starting as freshman at an engineering school in the fall.

We made it through all of single variable calculus without too much trouble, but now we're working through multivariable calculus and a specific set of exercises from Khan Academy is giving me trouble.

I have code working that arrives at the correct solution, but it involves doing an awkward python list comprehension.

I am seeking advice on whether there's a more idiomatic way to evaluate a scalar field at a vector, whether that means expressing things differently or using a different set of functions, etc.

I have isolated an example problem into a stand-alone python file.

I couldn't quite figure out how to best include a code snipped to this list. I will link to a gist and then also include the code below between triple backticks.

Thank you.

/rdj

https://gist.github.com/rdj/af9a5da6a10bab1d821c82ae8d5e13ff

```python
#!/usr/bin/env python3

# The Khan Academy exercises related to the multivariable chain rule
# ask you to find the derivative of the composition of a scalar field
# and a vector-valued function. Their article is available here:
#
# https://www.khanacademy.org/math/multivariable-calculus/multivariable-derivatives/differentiating-vector-valued-functions/a/multivariable-chain-rule-simple-version?modal=1
#
# Here is one of the exercises from this section:
#
# Let f(x,y,z) = z*cos(y) + z^2*x and g(t) = (t, -t^2, -t).
#
# h(t) = f(g(t))
#
# Find h'(t).

import sympy as sp
import sympy.vector as spv

R3 = spv.CoordSys3D('')
x,y,z = R3.base_scalars()

# Let f(x,y,z) be a scalar field mapping R3 to R:
f = z*sp.cos(y) + z**2*x
print(F"f(x,y,z) = {f}")

# Let g(t) be a vector valued function mapping R to R3:
t = sp.Symbol("t", real=True)
g = t * R3.i + -t**2 * R3.j + -t * R3.k
print(F"g(t) = {g}")

# Find the derivative with respect to t of the composition:
# d/dt[f(g(t))] = ∇_{g'(t)}f(g(t))

# First, get g'(t), a new vector-valued function.
gp = sp.diff(g, t)
print(F"g'(t) = {gp}")

# Next, get the directional derivative of f(x,y,z) in the direction
# g'(t), a new scalar field.
del_gp_f = spv.directional_derivative(f, gp)
print(F"∇_{{g'(t)}} f(x,y,z) = {del_gp_f}")

# Now, evaluate the directional derivative at g(t), yielding the result
# d/dt[f(g(t))], which should be a real function mapping R to R.

# HERE - Is there a less awkward way to evaluate the scalar field at a
# given vector?
#
# Note: Of course we could calculate the desired result a different
# way by evaluating gradient(f) at g then dotting with g', but again
# mechanically that's evaluating a scalar field at a given vector so
# we'd have to jump through the same hoops.

dfog_dt = del_gp_f.subs([(s, g.components[v]) for v,s in zip(R3.base_vectors(), R3.base_scalars())])
print(F"d/dt[f(g(t))] = {dfog_dt}")
```

[Alan Bromborsky]

Look at the following link.  Galgebra is built on top of sympy -

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

You can define a set of basis vectors say e1, e2, and e3, and coefficients a1, a2, and a3, then a general vector a1*e1+a2*e2+a3*e3.  The coefficients can be functions of the coordinates.  You can also define a scalar function of the coordinates say f(x1,x2,x3) then the gradient operator is grad and grad*f is the gradient of f, etc..

--
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/c9e45103-82a5-49b1-acb0-260adeb579bfn%40googlegroups.com.

[Alan Bromborsky]

Here is an example of galgebra output -

On 5/12/21 6:18 PM, Ryan D Johnson wrote:

[Davide Sandona']

Hello Ryan,

# HERE - Is there a less awkward way to evaluate the scalar field at a
# given vector?

Not that I'm aware of. I would have done the same intricate loop that you did.

Note that you were "lucky": if your gp vector would have had at least one zero-component, then gp.components would

not contain the zero component, resulting in a KeyError which would have complicated things even further.

Vector Module definitely needs to be polished a little further in order to become user-friendly. I fear this discussion

might get lost with time. I would suggest opening either an issue or a discussion suggesting improvements to this

module, providing the aforementioned example. https://github.com/sympy/sympy/discussions

Davide.

To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/dd7d4528-1cf6-9ca2-160d-0467efcd2f5a%40gmail.com.

[Oscar Benjamin]

On Thu, 13 May 2021 at 10:02, Davide Sandona' <sandona...@gmail.com> wrote:

Hello Ryan,

# HERE - Is there a less awkward way to evaluate the scalar field at a
# given vector?

Not that I'm aware of. I would have done the same intricate loop that you did.

Note that you were "lucky": if your gp vector would have had at least one zero-component, then gp.components would

not contain the zero component, resulting in a KeyError which would have complicated things even further.

A better way to do this might be:

In [95]: dict(zip(R3.base_scalars(), g.to_matrix(R3)))

Out[95]: 

⎧             2        ⎫

⎨x_: t, y_: -t , z_: -t⎬

⎩                      ⎭

 

Maybe there should be some easier way to do that. Conceptually this is a bit muddled though as it conflates vectors with points. Formally what is wanted here is to evaluate a scalar field at a point represented by a particular position vector wrt some reference frame. I don't know if the vector module has the concept of a point as opposed to a vector.

Vector Module definitely needs to be polished a little further in order to become user-friendly. I fear this discussion

might get lost with time. I would suggest opening either an issue or a discussion suggesting improvements to this

module, providing the aforementioned example. https://github.com/sympy/sympy/discussions

It certainly should be easier to use the vector module. I think some thought needs to be given to how you can use it without needing to explicitly create or use CoordSys3D. Most end user applications like this have no need for multiple reference frames so needing to operate with them is a distraction.

Oscar

[Ryan D Johnson]

Oscar makes a good point here:

Conceptually this is a bit muddled though as it conflates vectors with points. Formally what is wanted here is to evaluate a scalar field at a point represented by a particular position vector wrt some reference frame. I don't know if the vector module has the concept of a point as opposed to a vector.

It turns out sympy.vector does indeed have a Point class. It occurs to me that I could also dot with each basis vector. Both of these handle the case of a zero-valued component without introducing an if clause to the list comprehension (a problem with my original example, as Davide pointed out). 

But as far as I can tell, there’s no way to get the coordinates in a format compatible with subs, so however you get the components you still end up having to zip them with the base vectors from CoordSys3D. Here are three conceptually more pleasing ways to get a tuple of components:

In [1]: import sympy

In [2]: import sympy.vector

In [3]: R3 = sympy.vector.CoordSys3D("R3")

In [4]: t = sympy.Symbol("t", real=True)

In [5]: g = t * R3.i - t**2 * R3.j - 0 * R3.k

In [6]: # Three different ways to get the components as a tuple:

In [7]: R3.origin.locate_new(‘P', g).express_coordinates(R3) # Using Point
Out[7]: (t, -t**2, 0)

In [8]: tuple(b.dot(g) for b in R3.base_vectors()) # dot product with each basis
Out[8]: (t, -t**2, 0)

In [9]: tuple(g.to_matrix(R3))
Out[9]: (t, -t**2, 0)

If I were just making stuff up, it would be great if there were operations to (a) get the point at the end of a vector given a starting point and (b) convert a point to a dictionary or tuple in a format that subs is happy to consume. So something like:

g.to_point(R3.origin).to_dict(R3)

But even if we just had (b) we could write 

R3.origin.locate_new(‘P', g).to_dict(R3)

Seems a shame to keep passing the coordinate system, but Point doesn’t look like it remembers it. 

/rdj

[Alan Bromborsky]

Here is galgebra/sympy code for vector derivatives in rectangular and spherical coordinate with output annotated using latex -

def derivatives_in_rectangular_coordinates():
    #Print_Function()
    X = (x,y,z) = symbols('x y z') #coordinates
    o3d = Ga('e_x e_y e_z',g=[1,1,1],coords=X) #Rectangular coordinate basis
    (ex,ey,ez) = o3d.mv() #basis vectors
    grad = o3d.grad #gradient operator

    f = o3d.mv('f','scalar',f=True) #scalar function
    A = o3d.mv('A','vector',f=True) #vector function
    B = o3d.mv('B','bivector',f=True) #bivector function
    C = o3d.mv('C','mv') #general multivector function
    gprint('f =',f)
    gprint('A =',A)
    gprint('B =',B)
    gprint('C =',C)

    gprint('\\nabla f =',grad*f) #gradient of scalar function
    gprint('\\nabla\\cdot A =',grad|A) #divergence of vector function
    gprint('\\nabla A =',grad*A) #geometric derivative of vector function

    gprint('\\nabla\\times A = -I*(\\nabla\\wedge A) =',-o3d.I()*(grad^A)) #curl of vector function
    gprint('\\nabla B =',grad*B) # geometric derivative of bivector function
    gprint('\\nabla\\wedge B =',grad^B) # exterior derivative of bivector function
    gprint('\\nabla\\cdot B =',grad|B) # divergence of bivector function
    return

def derivatives_in_spherical_coordinates():
    #Print_Function()
    X = (r,th,phi) = symbols('r theta phi')
    s3d = Ga('e_r e_theta e_phi',g=[1,r**2,r**2*sin(th)**2],coords=X,norm=True)
    (er,eth,ephi) = s3d.mv()
    grad = s3d.grad

    f = s3d.mv('f','scalar',f=True)
    A = s3d.mv('A','vector',f=True)
    B = s3d.mv('B','bivector',f=True)

    gprint('f =',f)
    gprint('A =',A)
    gprint('B =',B)

    gprint('\\nabla f =',grad*f)
    gprint('\\nabla\\cdot A =',grad|A)
    gprint('\\nabla\\times A = -I*(\\nabla\\wedge A) =',(-s3d.E()*(grad^A)).simplify())
    gprint('\\nabla\\wedge B =',grad^B)

    return

--
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/73A0AD28-CA29-4D87-923D-213355CDD1C6%40gmail.com.

[Ryan D Johnson]

Hi, Alan.

> Here is galgebra/sympy code for vector derivatives in rectangular and spherical coordinate with output annotated using latex -

I did have a look at the galgebra package, and I was able to define my scalar field and vectors and compute the gradient or directional derivative, just like with sympy.vector. Though it did require significantly more scaffolding that is likely to be opaque to my students. As I said, I’m just working through an intro multivariable calculus course with my kids, so I only need R3.

What I didn’t find in galgebra was anything that addressed my core problem: idiomatically evaluating a scalar field at a point specified by a position vector from the origin.

/rdj