How to compute a gradient, a divergence or a curl

This tutorial introduces some vector calculus capabilities of SageMath within the 3-dimensional Euclidean space. The corresponding tools have been developed via the SageManifolds project.

The tutorial is also available as a Jupyter notebook, either passive (nbviewer) or interactive (binder).

First stage: introduce the Euclidean 3-space

Before evaluating some vector-field operators, one needs to define the arena in which vector fields live, namely the 3-dimensional Euclidean space ${\mathbb{E}}^3$. In SageMath, we declare it, along with the standard Cartesian coordinates $(x,y,z)$, via EuclideanSpace:

Sage

sage: E.<x,y,z> = EuclideanSpace()
sage: E
Euclidean space E^3

Python

>>> from sage.all import *
>>> E = EuclideanSpace(names=('x', 'y', 'z',)); (x, y, z,) = E._first_ngens(3)
>>> E
Euclidean space E^3

Sage Live

E.<x,y,z> = EuclideanSpace()
E

Thanks to the notation <x,y,z> in the above declaration, the coordinates $(x,y,z)$ are immediately available as three symbolic variables x, y and z (there is no need to declare them via var(), i.e. to type x, y, z = var('x y z')):

Sage

sage: x is E.cartesian_coordinates()[1]
True
sage: y is E.cartesian_coordinates()[2]
True
sage: z is E.cartesian_coordinates()[3]
True
sage: type(z)
<class 'sage.symbolic.expression.Expression'>

Python

>>> from sage.all import *
>>> x is E.cartesian_coordinates()[Integer(1)]
True
>>> y is E.cartesian_coordinates()[Integer(2)]
True
>>> z is E.cartesian_coordinates()[Integer(3)]
True
>>> type(z)
<class 'sage.symbolic.expression.Expression'>

Sage Live

x is E.cartesian_coordinates()[1]
y is E.cartesian_coordinates()[2]
z is E.cartesian_coordinates()[3]
type(z)

Besides, ${\mathbb{E}}^3$ is endowed with the orthonormal vector frame $(e_x,e_y,e_z)$ associated with Cartesian coordinates:

Sage

sage: E.frames()
[Coordinate frame (E^3, (e_x,e_y,e_z))]

Python

>>> from sage.all import *
>>> E.frames()
[Coordinate frame (E^3, (e_x,e_y,e_z))]

Sage Live

E.frames()

At this stage, this is the default vector frame on ${\mathbb{E}}^3$ (being the only vector frame introduced so far):

Sage

sage: E.default_frame()
Coordinate frame (E^3, (e_x,e_y,e_z))

Python

>>> from sage.all import *
>>> E.default_frame()
Coordinate frame (E^3, (e_x,e_y,e_z))

Sage Live

E.default_frame()

Defining a vector field

We define a vector field on ${\mathbb{E}}^3$ from its components in the vector frame $(e_x,e_y,e_z)$:

Sage

sage: v = E.vector_field(-y, x, sin(x*y*z), name='v')
sage: v.display()
v = -y e_x + x e_y + sin(x*y*z) e_z

Python

>>> from sage.all import *
>>> v = E.vector_field(-y, x, sin(x*y*z), name='v')
>>> v.display()
v = -y e_x + x e_y + sin(x*y*z) e_z

Sage Live

v = E.vector_field(-y, x, sin(x*y*z), name='v')
v.display()

We can access to the components of $v$ via the square bracket operator:

Sage

sage: v[1]
-y
sage: v[:]
[-y, x, sin(x*y*z)]

Python

>>> from sage.all import *
>>> v[Integer(1)]
-y
>>> v[:]
[-y, x, sin(x*y*z)]

Sage Live

v[1]
v[:]

A vector field can evaluated at any point of ${\mathbb{E}}^3$:

Sage

sage: p = E((3,-2,1), name='p')
sage: p
Point p on the Euclidean space E^3
sage: p.coordinates()
(3, -2, 1)
sage: vp = v.at(p)
sage: vp
Vector v at Point p on the Euclidean space E^3
sage: vp.display()
v = 2 e_x + 3 e_y - sin(6) e_z

Python

>>> from sage.all import *
>>> p = E((Integer(3),-Integer(2),Integer(1)), name='p')
>>> p
Point p on the Euclidean space E^3
>>> p.coordinates()
(3, -2, 1)
>>> vp = v.at(p)
>>> vp
Vector v at Point p on the Euclidean space E^3
>>> vp.display()
v = 2 e_x + 3 e_y - sin(6) e_z

Sage Live

p = E((3,-2,1), name='p')
p
p.coordinates()
vp = v.at(p)
vp
vp.display()

Vector fields can be plotted:

Sage

sage: v.plot(max_range=1.5, scale=0.5)
Graphics3d Object

Python

>>> from sage.all import *
>>> v.plot(max_range=RealNumber('1.5'), scale=RealNumber('0.5'))
Graphics3d Object

Sage Live

v.plot(max_range=1.5, scale=0.5)

For customizing the plot, see the list of options in the documentation of plot(). For instance, to get a view of the orthogonal projection of $v$ in the plane $y=1$, do:

Sage

sage: v.plot(fixed_coords={y: 1}, ambient_coords=(x,z), max_range=1.5,
....:        scale=0.25)
Graphics object consisting of 81 graphics primitives

Python

>>> from sage.all import *
>>> v.plot(fixed_coords={y: Integer(1)}, ambient_coords=(x,z), max_range=RealNumber('1.5'),
...        scale=RealNumber('0.25'))
Graphics object consisting of 81 graphics primitives

Sage Live

v.plot(fixed_coords={y: 1}, ambient_coords=(x,z), max_range=1.5,
       scale=0.25)

We may define a vector field $u$ with generic components $(u_x,u_y,y_z)$:

Sage

sage: u = E.vector_field(function('u_x')(x,y,z),
....:                    function('u_y')(x,y,z),
....:                    function('u_z')(x,y,z),
....:                    name='u')
sage: u.display()
u = u_x(x, y, z) e_x + u_y(x, y, z) e_y + u_z(x, y, z) e_z
sage: u[:]
[u_x(x, y, z), u_y(x, y, z), u_z(x, y, z)]

Python

>>> from sage.all import *
>>> u = E.vector_field(function('u_x')(x,y,z),
...                    function('u_y')(x,y,z),
...                    function('u_z')(x,y,z),
...                    name='u')
>>> u.display()
u = u_x(x, y, z) e_x + u_y(x, y, z) e_y + u_z(x, y, z) e_z
>>> u[:]
[u_x(x, y, z), u_y(x, y, z), u_z(x, y, z)]

Sage Live

u = E.vector_field(function('u_x')(x,y,z),
                   function('u_y')(x,y,z),
                   function('u_z')(x,y,z),
                   name='u')
u.display()
u[:]

Its value at the point $p$ is then:

Sage

sage: up = u.at(p)
sage: up.display()
u = u_x(3, -2, 1) e_x + u_y(3, -2, 1) e_y + u_z(3, -2, 1) e_z

Python

>>> from sage.all import *
>>> up = u.at(p)
>>> up.display()
u = u_x(3, -2, 1) e_x + u_y(3, -2, 1) e_y + u_z(3, -2, 1) e_z

Sage Live

up = u.at(p)
up.display()

How to compute various vector products

Dot product

The dot (or scalar) product $u\cdot v$ of the vector fields $u$ and $v$ is obtained by the method dot_product(), which admits dot() as a shortcut alias:

Sage

sage: u.dot(v) == u[1]*v[1] + u[2]*v[2] + u[3]*v[3]
True

Python

>>> from sage.all import *
>>> u.dot(v) == u[Integer(1)]*v[Integer(1)] + u[Integer(2)]*v[Integer(2)] + u[Integer(3)]*v[Integer(3)]
True

Sage Live

u.dot(v) == u[1]*v[1] + u[2]*v[2] + u[3]*v[3]

$s=u\cdot v$ is a scalar field, i.e. a map ${\mathbb{E}}^3\to \mathbb{R}$:

Sage

sage: s = u.dot(v)
sage: s
Scalar field u.v on the Euclidean space E^3
sage: s.display()
u.v: E^3 → ℝ
   (x, y, z) ↦ -y*u_x(x, y, z) + x*u_y(x, y, z) + sin(x*y*z)*u_z(x, y, z)

Python

>>> from sage.all import *
>>> s = u.dot(v)
>>> s
Scalar field u.v on the Euclidean space E^3
>>> s.display()
u.v: E^3 → ℝ
   (x, y, z) ↦ -y*u_x(x, y, z) + x*u_y(x, y, z) + sin(x*y*z)*u_z(x, y, z)

Sage Live

s = u.dot(v)
s
s.display()

It maps points of ${\mathbb{E}}^3$ to real numbers:

Sage

sage: s(p)
-sin(6)*u_z(3, -2, 1) + 2*u_x(3, -2, 1) + 3*u_y(3, -2, 1)

Python

>>> from sage.all import *
>>> s(p)
-sin(6)*u_z(3, -2, 1) + 2*u_x(3, -2, 1) + 3*u_y(3, -2, 1)

Sage Live

s(p)

Its coordinate expression is:

Sage

sage: s.expr()
-y*u_x(x, y, z) + x*u_y(x, y, z) + sin(x*y*z)*u_z(x, y, z)

Python

>>> from sage.all import *
>>> s.expr()
-y*u_x(x, y, z) + x*u_y(x, y, z) + sin(x*y*z)*u_z(x, y, z)

Sage Live

s.expr()

Norm

The norm $\Vert u\Vert$ of the vector field $u$ is defined in terms of the dot product by $\Vert u\Vert =\sqrt{u\cdot u}$:

Sage

sage: norm(u) == sqrt(u.dot(u))
True

Python

>>> from sage.all import *
>>> norm(u) == sqrt(u.dot(u))
True

Sage Live

norm(u) == sqrt(u.dot(u))

It is a scalar field on ${\mathbb{E}}^3$:

Sage

sage: s = norm(u)
sage: s
Scalar field |u| on the Euclidean space E^3
sage: s.display()
|u|: E^3 → ℝ
   (x, y, z) ↦ sqrt(u_x(x, y, z)^2 + u_y(x, y, z)^2 + u_z(x, y, z)^2)
sage: s.expr()
sqrt(u_x(x, y, z)^2 + u_y(x, y, z)^2 + u_z(x, y, z)^2)

Python

>>> from sage.all import *
>>> s = norm(u)
>>> s
Scalar field |u| on the Euclidean space E^3
>>> s.display()
|u|: E^3 → ℝ
   (x, y, z) ↦ sqrt(u_x(x, y, z)^2 + u_y(x, y, z)^2 + u_z(x, y, z)^2)
>>> s.expr()
sqrt(u_x(x, y, z)^2 + u_y(x, y, z)^2 + u_z(x, y, z)^2)

Sage Live

s = norm(u)
s
s.display()
s.expr()

For $v$, we have:

Sage

sage: norm(v).expr()
sqrt(x^2 + y^2 + sin(x*y*z)^2)

Python

>>> from sage.all import *
>>> norm(v).expr()
sqrt(x^2 + y^2 + sin(x*y*z)^2)

Sage Live

norm(v).expr()

Cross product

The cross product $u\times v$ is obtained by the method cross_product(), which admits cross() as a shortcut alias:

Sage

sage: s = u.cross(v)
sage: s
Vector field u x v on the Euclidean space E^3
sage: s.display()
u x v = (sin(x*y*z)*u_y(x, y, z) - x*u_z(x, y, z)) e_x
 + (-sin(x*y*z)*u_x(x, y, z) - y*u_z(x, y, z)) e_y
 + (x*u_x(x, y, z) + y*u_y(x, y, z)) e_z

Python

>>> from sage.all import *
>>> s = u.cross(v)
>>> s
Vector field u x v on the Euclidean space E^3
>>> s.display()
u x v = (sin(x*y*z)*u_y(x, y, z) - x*u_z(x, y, z)) e_x
 + (-sin(x*y*z)*u_x(x, y, z) - y*u_z(x, y, z)) e_y
 + (x*u_x(x, y, z) + y*u_y(x, y, z)) e_z

Sage Live

s = u.cross(v)
s
s.display()

We can check the standard formulas expressing the cross product in terms of the components:

Sage

sage: all([s[1] == u[2]*v[3] - u[3]*v[2],
....:      s[2] == u[3]*v[1] - u[1]*v[3],
....:      s[3] == u[1]*v[2] - u[2]*v[1]])
True

Python

>>> from sage.all import *
>>> all([s[Integer(1)] == u[Integer(2)]*v[Integer(3)] - u[Integer(3)]*v[Integer(2)],
...      s[Integer(2)] == u[Integer(3)]*v[Integer(1)] - u[Integer(1)]*v[Integer(3)],
...      s[Integer(3)] == u[Integer(1)]*v[Integer(2)] - u[Integer(2)]*v[Integer(1)]])
True

Sage Live

all([s[1] == u[2]*v[3] - u[3]*v[2],
     s[2] == u[3]*v[1] - u[1]*v[3],
     s[3] == u[1]*v[2] - u[2]*v[1]])

Scalar triple product

Let us introduce a third vector field, $w$ say. As a example, we do not pass the components as arguments of vector_field(), as we did for $u$ and $v$; instead, we set them in a second stage, via the square bracket operator, any unset component being assumed to be zero:

Sage

sage: w = E.vector_field(name='w')
sage: w[1] = x*z
sage: w[2] = y*z
sage: w.display()
w = x*z e_x + y*z e_y

Python

>>> from sage.all import *
>>> w = E.vector_field(name='w')
>>> w[Integer(1)] = x*z
>>> w[Integer(2)] = y*z
>>> w.display()
w = x*z e_x + y*z e_y

Sage Live

w = E.vector_field(name='w')
w[1] = x*z
w[2] = y*z
w.display()

The scalar triple product of the vector fields $u$, $v$ and $w$ is obtained as follows:

Sage

sage: triple_product = E.scalar_triple_product()
sage: s = triple_product(u, v, w)
sage: s
Scalar field epsilon(u,v,w) on the Euclidean space E^3
sage: s.expr()
-(y*u_x(x, y, z) - x*u_y(x, y, z))*z*sin(x*y*z) - (x^2*u_z(x, y, z)
 + y^2*u_z(x, y, z))*z

Python

>>> from sage.all import *
>>> triple_product = E.scalar_triple_product()
>>> s = triple_product(u, v, w)
>>> s
Scalar field epsilon(u,v,w) on the Euclidean space E^3
>>> s.expr()
-(y*u_x(x, y, z) - x*u_y(x, y, z))*z*sin(x*y*z) - (x^2*u_z(x, y, z)
 + y^2*u_z(x, y, z))*z

Sage Live

triple_product = E.scalar_triple_product()
s = triple_product(u, v, w)
s
s.expr()

Let us check that the scalar triple product of $u$, $v$ and $w$ is $u\cdot (v\times w)$:

Sage

sage: s == u.dot(v.cross(w))
True

Python

>>> from sage.all import *
>>> s == u.dot(v.cross(w))
True

Sage Live

s == u.dot(v.cross(w))

How to evaluate the standard differential operators

The standard operators $\mathrm{grad}$, $\mathrm{div}$, $\mathrm{curl}$, etc. involved in vector calculus are accessible as methods on scalar fields and vector fields (e.g. v.div()). However, to allow for standard mathematical notations (e.g. div(v)), let us import the functions grad(), div(), curl() and laplacian():

Sage

sage: from sage.manifolds.operators import *

Python

>>> from sage.all import *
>>> from sage.manifolds.operators import *

Sage Live

from sage.manifolds.operators import *

Gradient

We first introduce a scalar field, via its expression in terms of Cartesian coordinates; in this example, we consider some unspecified function of $(x,y,z)$:

Sage

sage: F = E.scalar_field(function('f')(x,y,z), name='F')
sage: F.display()
F: E^3 → ℝ
   (x, y, z) ↦ f(x, y, z)

Python

>>> from sage.all import *
>>> F = E.scalar_field(function('f')(x,y,z), name='F')
>>> F.display()
F: E^3 → ℝ
   (x, y, z) ↦ f(x, y, z)

Sage Live

F = E.scalar_field(function('f')(x,y,z), name='F')
F.display()

The value of $F$ at a point:

Sage

sage: F(p)
f(3, -2, 1)

Python

>>> from sage.all import *
>>> F(p)
f(3, -2, 1)

Sage Live

F(p)

The gradient of $F$:

Sage

sage: grad(F)
Vector field grad(F) on the Euclidean space E^3
sage: grad(F).display()
grad(F) = d(f)/dx e_x + d(f)/dy e_y + d(f)/dz e_z
sage: norm(grad(F)).display()
|grad(F)|: E^3 → ℝ
   (x, y, z) ↦ sqrt((d(f)/dx)^2 + (d(f)/dy)^2 + (d(f)/dz)^2)

Python

>>> from sage.all import *
>>> grad(F)
Vector field grad(F) on the Euclidean space E^3
>>> grad(F).display()
grad(F) = d(f)/dx e_x + d(f)/dy e_y + d(f)/dz e_z
>>> norm(grad(F)).display()
|grad(F)|: E^3 → ℝ
   (x, y, z) ↦ sqrt((d(f)/dx)^2 + (d(f)/dy)^2 + (d(f)/dz)^2)

Sage Live

grad(F)
grad(F).display()
norm(grad(F)).display()

Divergence

The divergence of the vector field $u$:

Sage

sage: s = div(u)
sage: s.display()
div(u): E^3 → ℝ
   (x, y, z) ↦ d(u_x)/dx + d(u_y)/dy + d(u_z)/dz

Python

>>> from sage.all import *
>>> s = div(u)
>>> s.display()
div(u): E^3 → ℝ
   (x, y, z) ↦ d(u_x)/dx + d(u_y)/dy + d(u_z)/dz

Sage Live

s = div(u)
s.display()

For $v$ and $w$, we have:

Sage

sage: div(v).expr()
x*y*cos(x*y*z)
sage: div(w).expr()
2*z

Python

>>> from sage.all import *
>>> div(v).expr()
x*y*cos(x*y*z)
>>> div(w).expr()
2*z

Sage Live

div(v).expr()
div(w).expr()

An identity valid for any scalar field $F$ and any vector field $u$:

Sage

sage: div(F*u) == F*div(u) + u.dot(grad(F))
True

Python

>>> from sage.all import *
>>> div(F*u) == F*div(u) + u.dot(grad(F))
True

Sage Live

div(F*u) == F*div(u) + u.dot(grad(F))

Curl

The curl of the vector field $u$:

Sage

sage: s = curl(u)
sage: s
Vector field curl(u) on the Euclidean space E^3
sage: s.display()
curl(u) = (-d(u_y)/dz + d(u_z)/dy) e_x + (d(u_x)/dz - d(u_z)/dx) e_y
 + (-d(u_x)/dy + d(u_y)/dx) e_z

Python

>>> from sage.all import *
>>> s = curl(u)
>>> s
Vector field curl(u) on the Euclidean space E^3
>>> s.display()
curl(u) = (-d(u_y)/dz + d(u_z)/dy) e_x + (d(u_x)/dz - d(u_z)/dx) e_y
 + (-d(u_x)/dy + d(u_y)/dx) e_z

Sage Live

s = curl(u)
s
s.display()

To use the notation rot instead of curl, simply do:

Sage

sage: rot = curl

Python

>>> from sage.all import *
>>> rot = curl

Sage Live

rot = curl

An alternative is:

Sage

sage: from sage.manifolds.operators import curl as rot

Python

>>> from sage.all import *
>>> from sage.manifolds.operators import curl as rot

Sage Live

from sage.manifolds.operators import curl as rot

We have then:

Sage

sage: rot(u).display()
curl(u) = (-d(u_y)/dz + d(u_z)/dy) e_x + (d(u_x)/dz - d(u_z)/dx) e_y
 + (-d(u_x)/dy + d(u_y)/dx) e_z
sage: rot(u) == curl(u)
True

Python

>>> from sage.all import *
>>> rot(u).display()
curl(u) = (-d(u_y)/dz + d(u_z)/dy) e_x + (d(u_x)/dz - d(u_z)/dx) e_y
 + (-d(u_x)/dy + d(u_y)/dx) e_z
>>> rot(u) == curl(u)
True

Sage Live

rot(u).display()
rot(u) == curl(u)

For $v$ and $w$, we have:

Sage

sage: curl(v).display()
curl(v) = x*z*cos(x*y*z) e_x - y*z*cos(x*y*z) e_y + 2 e_z

Python

>>> from sage.all import *
>>> curl(v).display()
curl(v) = x*z*cos(x*y*z) e_x - y*z*cos(x*y*z) e_y + 2 e_z

Sage Live

curl(v).display()

Sage

sage: curl(w).display()
curl(w) = -y e_x + x e_y

Python

>>> from sage.all import *
>>> curl(w).display()
curl(w) = -y e_x + x e_y

Sage Live

curl(w).display()

The curl of a gradient is always zero:

Sage

sage: curl(grad(F)).display()
curl(grad(F)) = 0

Python

>>> from sage.all import *
>>> curl(grad(F)).display()
curl(grad(F)) = 0

Sage Live

curl(grad(F)).display()

The divergence of a curl is always zero:

Sage

sage: div(curl(u)).display()
div(curl(u)): E^3 → ℝ
   (x, y, z) ↦ 0

Python

>>> from sage.all import *
>>> div(curl(u)).display()
div(curl(u)): E^3 → ℝ
   (x, y, z) ↦ 0

Sage Live

div(curl(u)).display()

An identity valid for any scalar field $F$ and any vector field $u$ is

$$\mathrm{curl}(Fu)=\mathrm{grad}F\times u+F\mathrm{curl}u,$$

as we can check:

Sage

sage: curl(F*u) == grad(F).cross(u) + F*curl(u)
True

Python

>>> from sage.all import *
>>> curl(F*u) == grad(F).cross(u) + F*curl(u)
True

Sage Live

curl(F*u) == grad(F).cross(u) + F*curl(u)

Laplacian

The Laplacian $\mathrm{\Delta }F$ of a scalar field $F$ is another scalar field:

Sage

sage: s = laplacian(F)
sage: s.display()
Delta(F): E^3 → ℝ
   (x, y, z) ↦ d^2(f)/dx^2 + d^2(f)/dy^2 + d^2(f)/dz^2

Python

>>> from sage.all import *
>>> s = laplacian(F)
>>> s.display()
Delta(F): E^3 → ℝ
   (x, y, z) ↦ d^2(f)/dx^2 + d^2(f)/dy^2 + d^2(f)/dz^2

Sage Live

s = laplacian(F)
s.display()

The following identity holds:

$$\mathrm{\Delta }F=\mathrm{div}\left(\mathrm{grad}F\right),$$

as we can check:

Sage

sage: laplacian(F) == div(grad(F))
True

Python

>>> from sage.all import *
>>> laplacian(F) == div(grad(F))
True

Sage Live

laplacian(F) == div(grad(F))

The Laplacian $\mathrm{\Delta }u$ of a vector field $u$ is another vector field:

Sage

sage: Du = laplacian(u)
sage: Du
Vector field Delta(u) on the Euclidean space E^3

Python

>>> from sage.all import *
>>> Du = laplacian(u)
>>> Du
Vector field Delta(u) on the Euclidean space E^3

Sage Live

Du = laplacian(u)
Du

whose components are:

Sage

sage: Du.display()
Delta(u) = (d^2(u_x)/dx^2 + d^2(u_x)/dy^2 + d^2(u_x)/dz^2) e_x
 + (d^2(u_y)/dx^2 + d^2(u_y)/dy^2 + d^2(u_y)/dz^2) e_y
 + (d^2(u_z)/dx^2 + d^2(u_z)/dy^2 + d^2(u_z)/dz^2) e_z

Python

>>> from sage.all import *
>>> Du.display()
Delta(u) = (d^2(u_x)/dx^2 + d^2(u_x)/dy^2 + d^2(u_x)/dz^2) e_x
 + (d^2(u_y)/dx^2 + d^2(u_y)/dy^2 + d^2(u_y)/dz^2) e_y
 + (d^2(u_z)/dx^2 + d^2(u_z)/dy^2 + d^2(u_z)/dz^2) e_z

Sage Live

Du.display()

In the Cartesian frame, the components of $\mathrm{\Delta }u$ are nothing but the (scalar) Laplacians of the components of $u$, as we can check:

Sage

sage: e = E.cartesian_frame()
sage: Du == sum(laplacian(u[[i]])*e[i] for i in E.irange())
True

Python

>>> from sage.all import *
>>> e = E.cartesian_frame()
>>> Du == sum(laplacian(u[[i]])*e[i] for i in E.irange())
True

Sage Live

e = E.cartesian_frame()
Du == sum(laplacian(u[[i]])*e[i] for i in E.irange())

In the above formula, u[[i]] return the $i$-th component of $u$ as a scalar field, while u[i] would have returned the coordinate expression of this scalar field; besides, e is the Cartesian frame:

Sage

sage: e[:]
(Vector field e_x on the Euclidean space E^3,
 Vector field e_y on the Euclidean space E^3,
 Vector field e_z on the Euclidean space E^3)

Python

>>> from sage.all import *
>>> e[:]
(Vector field e_x on the Euclidean space E^3,
 Vector field e_y on the Euclidean space E^3,
 Vector field e_z on the Euclidean space E^3)

Sage Live

e[:]

For the vector fields $v$ and $w$, we have:

Sage

sage: laplacian(v).display()
Delta(v) = -(x^2*y^2 + (x^2 + y^2)*z^2)*sin(x*y*z) e_z
sage: laplacian(w).display()
Delta(w) = 0

Python

>>> from sage.all import *
>>> laplacian(v).display()
Delta(v) = -(x^2*y^2 + (x^2 + y^2)*z^2)*sin(x*y*z) e_z
>>> laplacian(w).display()
Delta(w) = 0

Sage Live

laplacian(v).display()
laplacian(w).display()

We have:

Sage

sage: curl(curl(u)).display()
curl(curl(u)) = (-d^2(u_x)/dy^2 - d^2(u_x)/dz^2 + d^2(u_y)/dxdy
 + d^2(u_z)/dxdz) e_x + (d^2(u_x)/dxdy - d^2(u_y)/dx^2 - d^2(u_y)/dz^2
 + d^2(u_z)/dydz) e_y + (d^2(u_x)/dxdz + d^2(u_y)/dydz - d^2(u_z)/dx^2
 - d^2(u_z)/dy^2) e_z
sage: grad(div(u)).display()
grad(div(u)) = (d^2(u_x)/dx^2 + d^2(u_y)/dxdy + d^2(u_z)/dxdz) e_x
 + (d^2(u_x)/dxdy + d^2(u_y)/dy^2 + d^2(u_z)/dydz) e_y
 + (d^2(u_x)/dxdz + d^2(u_y)/dydz + d^2(u_z)/dz^2) e_z

Python

>>> from sage.all import *
>>> curl(curl(u)).display()
curl(curl(u)) = (-d^2(u_x)/dy^2 - d^2(u_x)/dz^2 + d^2(u_y)/dxdy
 + d^2(u_z)/dxdz) e_x + (d^2(u_x)/dxdy - d^2(u_y)/dx^2 - d^2(u_y)/dz^2
 + d^2(u_z)/dydz) e_y + (d^2(u_x)/dxdz + d^2(u_y)/dydz - d^2(u_z)/dx^2
 - d^2(u_z)/dy^2) e_z
>>> grad(div(u)).display()
grad(div(u)) = (d^2(u_x)/dx^2 + d^2(u_y)/dxdy + d^2(u_z)/dxdz) e_x
 + (d^2(u_x)/dxdy + d^2(u_y)/dy^2 + d^2(u_z)/dydz) e_y
 + (d^2(u_x)/dxdz + d^2(u_y)/dydz + d^2(u_z)/dz^2) e_z

Sage Live

curl(curl(u)).display()
grad(div(u)).display()

A famous identity is

$$\mathrm{curl}\left(\mathrm{curl}u\right)=\mathrm{grad}\left(\mathrm{div}u\right)-\mathrm{\Delta }u.$$

Let us check it:

Sage

sage: curl(curl(u)) == grad(div(u)) - laplacian(u)
True

Python

>>> from sage.all import *
>>> curl(curl(u)) == grad(div(u)) - laplacian(u)
True

Sage Live

curl(curl(u)) == grad(div(u)) - laplacian(u)

How to customize various symbols

Customizing the symbols of the orthonormal frame vectors

By default, the vectors of the orthonormal frame associated with Cartesian coordinates are denoted $(e_x,e_y,e_z)$:

Sage

sage: frame = E.cartesian_frame()
sage: frame
Coordinate frame (E^3, (e_x,e_y,e_z))

Python

>>> from sage.all import *
>>> frame = E.cartesian_frame()
>>> frame
Coordinate frame (E^3, (e_x,e_y,e_z))

Sage Live

frame = E.cartesian_frame()
frame

But this can be changed, thanks to the method set_name():

Sage

sage: frame.set_name('a', indices=('x', 'y', 'z'))
sage: frame
Coordinate frame (E^3, (a_x,a_y,a_z))
sage: v.display()
v = -y a_x + x a_y + sin(x*y*z) a_z

Python

>>> from sage.all import *
>>> frame.set_name('a', indices=('x', 'y', 'z'))
>>> frame
Coordinate frame (E^3, (a_x,a_y,a_z))
>>> v.display()
v = -y a_x + x a_y + sin(x*y*z) a_z

Sage Live

frame.set_name('a', indices=('x', 'y', 'z'))
frame
v.display()

Sage

sage: frame.set_name(('hx', 'hy', 'hz'),
....:                latex_symbol=(r'\mathbf{\hat{x}}', r'\mathbf{\hat{y}}',
....:                              r'\mathbf{\hat{z}}'))
sage: frame
Coordinate frame (E^3, (hx,hy,hz))
sage: v.display()
v = -y hx + x hy + sin(x*y*z) hz

Python

>>> from sage.all import *
>>> frame.set_name(('hx', 'hy', 'hz'),
...                latex_symbol=(r'\mathbf{\hat{x}}', r'\mathbf{\hat{y}}',
...                              r'\mathbf{\hat{z}}'))
>>> frame
Coordinate frame (E^3, (hx,hy,hz))
>>> v.display()
v = -y hx + x hy + sin(x*y*z) hz

Sage Live

frame.set_name(('hx', 'hy', 'hz'),
               latex_symbol=(r'\mathbf{\hat{x}}', r'\mathbf{\hat{y}}',
                             r'\mathbf{\hat{z}}'))
frame
v.display()

Customizing the coordinate symbols

The coordinates symbols are defined within the angle brackets <...> at the construction of the Euclidean space. Above we did:

Sage

sage: E.<x,y,z> = EuclideanSpace()

Python

>>> from sage.all import *
>>> E = EuclideanSpace(names=('x', 'y', 'z',)); (x, y, z,) = E._first_ngens(3)

Sage Live

E.<x,y,z> = EuclideanSpace()

which resulted in the coordinate symbols $(x,y,z)$ and in the corresponding Python variables x, y and z (SageMath symbolic expressions). To use other symbols, for instance $(X,Y,Z)$, it suffices to create E as:

Sage

sage: E.<X,Y,Z> = EuclideanSpace()

Python

>>> from sage.all import *
>>> E = EuclideanSpace(names=('X', 'Y', 'Z',)); (X, Y, Z,) = E._first_ngens(3)

Sage Live

E.<X,Y,Z> = EuclideanSpace()

We have then:

Sage

sage: E.atlas()
[Chart (E^3, (X, Y, Z))]
sage: E.cartesian_frame()
Coordinate frame (E^3, (e_X,e_Y,e_Z))
sage: v = E.vector_field(-Y, X, sin(X*Y*Z), name='v')
sage: v.display()
v = -Y e_X + X e_Y + sin(X*Y*Z) e_Z

Python

>>> from sage.all import *
>>> E.atlas()
[Chart (E^3, (X, Y, Z))]
>>> E.cartesian_frame()
Coordinate frame (E^3, (e_X,e_Y,e_Z))
>>> v = E.vector_field(-Y, X, sin(X*Y*Z), name='v')
>>> v.display()
v = -Y e_X + X e_Y + sin(X*Y*Z) e_Z

Sage Live

E.atlas()
E.cartesian_frame()
v = E.vector_field(-Y, X, sin(X*Y*Z), name='v')
v.display()

By default the LaTeX symbols of the coordinate coincide with the letters given within the angle brackets. But this can be adjusted through the optional argument symbols of the function EuclideanSpace, which has to be a string, usually prefixed by r (for raw string, in order to allow for the backslash character of LaTeX expressions). This string contains the coordinate fields separated by a blank space; each field contains the coordinate’s text symbol and possibly the coordinate’s LaTeX symbol (when the latter is different from the text symbol), both symbols being separated by a colon (:):

Sage

sage: E.<xi,et,ze> = EuclideanSpace(symbols=r"xi:\xi et:\eta ze:\zeta")
sage: E.atlas()
[Chart (E^3, (xi, et, ze))]
sage: v = E.vector_field(-et, xi, sin(xi*et*ze), name='v')
sage: v.display()
v = -et e_xi + xi e_et + sin(et*xi*ze) e_ze

Python

>>> from sage.all import *
>>> E = EuclideanSpace(symbols=r"xi:\xi et:\eta ze:\zeta", names=('xi', 'et', 'ze',)); (xi, et, ze,) = E._first_ngens(3)
>>> E.atlas()
[Chart (E^3, (xi, et, ze))]
>>> v = E.vector_field(-et, xi, sin(xi*et*ze), name='v')
>>> v.display()
v = -et e_xi + xi e_et + sin(et*xi*ze) e_ze

Sage Live

E.<xi,et,ze> = EuclideanSpace(symbols=r"xi:\xi et:\eta ze:\zeta")
E.atlas()
v = E.vector_field(-et, xi, sin(xi*et*ze), name='v')
v.display()

Make live