#!/usr/bin/env python
# coding: utf-8

# # Vector calculus with SageMath
# 
# This worksheet illustrates the operators regarding scalar and vector fields on pseudo-Riemannian manifolds introduced in Trac ticket [#24622](https://trac.sagemath.org/ticket/24622).
# 
# Since SageMath 8.3, it is rather obsolete regarding calculus in Euclidean spaces. See these 
# [vector calculus examples](https://sagemanifolds.obspm.fr/vector_calculus.html) instead. 

# In[1]:


version()


# In[2]:


get_ipython().run_line_magic('display', 'latex')


# In[3]:


from sage.manifolds.operators import *   # to get the operators grad, div, curl, etc.


# ### The Euclidean space as a 3-dimensional Riemannian manifold

# In[4]:


M = Manifold(3, 'M', structure='Riemannian', start_index=1)
X.<x,y,z> = M.chart()


# In[5]:


g = M.metric()
g[1,1], g[2,2], g[3,3] = 1, 1, 1
g.display()


# ### Gradient of a scalar field

# In[6]:


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


# In[7]:


grad(F).display()


# In[8]:


curl(grad(F)).display()


# In[9]:


norm(grad(F)).display()


# ### Laplacien of a scalar field

# In[10]:


laplacian(F).display()


# In[11]:


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


# ### Vector field

# In[12]:


v = M.vector_field(name='v')
v[1] = function('v_x')(x,y,z)
v[2] = function('v_y')(x,y,z)
v[3] = function('v_z')(x,y,z)
v.display()


# In order not to clutter the outputs, we omit the coordinate arguments in the display of chart functions:

# In[13]:


M.options.omit_function_arguments=True


# In[14]:


v.display()


# Norm of a vector field:

# In[15]:


s = norm(v)
print(s)


# In[16]:


s.display()


# Scalar product of two vector fields:

# In[17]:


v.dot(grad(F)).display()


# Cross product of two vector fields:

# In[18]:


v.cross(grad(F)).display()


# Divergence :

# In[19]:


s = div(v)
s.display()


# Curl:

# In[20]:


w = curl(v)
print(w)


# In[21]:


w.display()


# To use the notation `rot` instead of  `curl`, simply do

# In[22]:


rot = curl


# An alternative is

# In[23]:


from sage.manifolds.operators import curl as rot


# We have then

# In[24]:


rot(v).display()


# In[25]:


rot(v) == curl(v)


# The divergence of a curl is always zero:

# In[26]:


div(curl(v)).display()


# Laplacian of a vector field:

# In[27]:


laplacian(v).display()


# In[28]:


curl(curl(v)).display()


# In[29]:


grad(div(v)).display()


# Check of a famous identity:

# In[30]:


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


# Two other identities regarding any scalar field `F` and any vector field `v`:

# In[31]:


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


# In[32]:


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


# The left-hand side is

# In[33]:


curl(F*v).display()


# In[ ]: