> 
#
# Initialize: load VectorCalculus
#
restart;
with(VectorCalculus):

> 
#
# Define OP's "toroidal" coordinate system
#
X := a*sinh(tau)*cos(phi)/(cosh(tau)  cos(sigma)):
Y := a*sinh(tau)*sin(phi)/(cosh(tau)  cos(sigma)):
Z := a*sin(sigma)/(cosh(tau)  cos(sigma)):
#
# Add OP's toroidal coordinate system to those known
# by the VectorCalculus package.
#
# NB this system which the OP calls "toroidal" is
# slightly different from the Maple's inbuilt toroidal
# system  so it has to be given a different "name"
#
AddCoordinates
( OPToroid[ sigma, tau, phi ],
[ X, Y, Z ]
):

> 
#
# Compute the Laplacian for an arbitrary function in
# OP's toroidal system
#
LP:= simplify
( expand
( Laplacian
( f(sigma, tau, phi),
OPToroid[sigma, tau, phi]
)
)
);
#
# Painfully, type in the Laplacian quoted on the
# Wikipedia page
#
# https://en.wikipedia.org/wiki/Toroidal_coordinates
#
wiki:= simplify
( expand
( ( cosh(tau)cos(sigma) )^3/(a^2*sinh(tau))
*
( sinh(tau)*diff(diff(f(sigma,tau,phi), sigma)/(cosh(tau)cos(sigma)),sigma)
+
diff( sinh(tau)*diff(f(sigma,tau,phi),tau)/(cosh(tau)cos(sigma)),tau)
+
diff(f(sigma,tau,phi),phi$2)/(sinh(tau)*(cosh(tau)cos(sigma)))
)
)
);
#
# Check whether the Laplacian derived by Maple from the
# definition of the coordinate system above (ie LP) is
# the same as the Laplacian quaoted by Wikipedia (ie wiki)
#
# They *look* the same, but are they???
#
is( wiki=LP );
#
# Laplacian of a specific (but arbitrary) function in the
# "new" coordinate system
#
simplify
( expand
( Laplacian
( exp(sigma)*sin(tau)*cos(phi),
OPToroid[sigma, tau, phi]
)
)
);








(1) 
> 
#########################################################
# Gradient for a general function in the "new" coordinate
# system
#
Gradient
( f(sigma, tau, phi),
OPToroid[sigma, tau, phi]
);
#
# Gradient of a specific (but arbitrary) function in the
# "new" coordinate system
#
Gradient
( exp(sigma)*sin(2*tau)*cos(phi+Pi/4),
OPToroid[sigma, tau, phi]
);


(2) 
> 
###################################################
# Divergence of a general vector field in the "new"
# coordinate system
#
simplify
( expand
( Divergence
( VectorField
( < f(sigma, tau, phi),
g(sigma, tau, phi),
h(sigma, tau, phi)
>,
OPToroid[ sigma, tau, phi ]
)
)
)
);
#
# Divergence of a specific (but arbitrary) vector
# field in the "new" coordinate system
#
simplify
( expand
( Divergence
( VectorField
( < exp(tau)*sin(sigma),
exp(tau)*cos(sigma),
sin(sigma)*cos(phi)
>,
OPToroid[ sigma, tau, phi ]
)
)
)
);


(3) 
> 
#######################################################
# The Divergence of the Gradient of a general function
# in the "new" coordinate system.
#
# Note that the Gradient of a (scalar) function returns
# a vector function and the Divergence of a vector
# function returns a scalar function, so
#
# Divergence(Gradient())
#
# accepts a scalar function and returns a scalar functon
#
d:= simplify
( expand
( Divergence
( VectorField
( Gradient
( f(sigma, tau, phi),
OPToroid[ sigma, tau, phi ]
),
OPToroid[ sigma, tau, phi ]
)
)
)
);
#
# Check that the Divergence of the Gradient of a
# general function is equal to the Laplacian
# computed earlier
#
is( d=LP );


(4) 
> 
########################################################
# The Gradient of the Divergence of a general function
# in the "new" coordinate system.
#
# Not that the Divergence of a (vector) function returns
# a scalar function and the Gradient of a scalar function
# returns a vector function, so
#
# Gradient(Divergence())
#
# accepts a vector function and returns a evctor functon
#
simplify
( expand
( Gradient
( Divergence
( VectorField
( < f(sigma, tau, phi),
g(sigma, tau, phi),
h(sigma, tau, phi)
>,
OPToroid[ sigma, tau, phi ]
)
),
OPToroid[ sigma, tau, phi ]
)
)
);


(5) 
