reduce-historical/r34/doc/orthovec.tex

\documentstyle[12pt]{article}
%\documentstyle{paper}
\setlength{\textheight}{9.0in}
\setlength{\textwidth}{6.5in}
\setlength{\oddsidemargin}{-0.0in}
\setlength{\topmargin}{0.0in}
\begin{document}
\begin{center}
{\Large {\bf ORTHOVEC} \\[3ex]
Version 2 of the REDUCE program \\
for 3-D vector analysis in \\ 
orthogonal curvilinear coordinates }\\
\vspace{1.25in}
{\large  James W.~EASTWOOD}\\[2ex]
{\large  AEA Technology}\\
{\large  Culham Laboratory}\\
{\large  Abingdon}\\
{\large  Oxon OX14 3DB}\\[2ex]
{\large  June 1990} \\
\normalsize
\vspace{0.50in}
Submitted to {\em Computer Physics Communications} \\
\end{center}
\vspace{0.75in}

\begin{abstract}
The revised version of ORTHOVEC is a collection of REDUCE 3.3 
procedures and operations which provide a simple to use environment
for the manipulation of scalars and vectors.
Operations include addition, 
subtraction, dot and cross products, division, modulus, div,
grad, curl, laplacian, differentiation, integration, ${\bf a \cdot 
\nabla}$ and Taylor expansion.  Version 2 is summarized in Comp.  Phys.
Commun. 64 (1991) 121-122.  It differs from the original (Comp.  Phys.
Commun. 47 (1987) 139-147) in revised notation and extended capabilities.
\end{abstract}

%\newpage
%\tableofcontents
\newpage
\begin{center}
{\Large{\bf New Version Summary}}
\end{center}
\begin{tabular}{ll}
\underline{Title of program}:&ORTHOVEC\\[2ex]
\underline{Catalogue number}:&AAXY\\[2ex]
\underline{Program obtainable from}: &CPC Program Library,\\
&Queen's University of Belfast, N.~Ireland\\[2ex]
\underline{Reference to original program}: &CPC 47 (1987) 139-147\\[2ex]
\underline{Operating system}:&UNIX, MS-DOS + ARM-OS\\[2ex]
\underline{Programming Language used}: &REDUCE 3.3\\[2ex]
\underline{High speed storage required}: &As for 
the underlying PSL/REDUCE \\
&system, typically $>$ 1 Megabyte\\[2ex]
\underline{No. of lines in combined programs and test deck}:&600 \\[2ex]
\underline{Keywords}: & Computer Algebra, Vector Analysis,\\
& series Expansion,  Plasma Physics, \\
&Hydrodynamics, Electromagnetics.\\[2ex]
\underline{Author of original program}: &James W. EASTWOOD\\[2ex]
\underline{Nature of Physical Problem}: 
&There is a wide range using vector\\
& calculus in orthogonal curvilinear coordinates\\
& and vector integration, differentiation\\
& and series expansion.\\[2ex]
\underline{Method of Solution}: & computer aided algebra using\\
&standard orthogonal curvilinear coordinates\\
&for differential and integral operators.\\[2ex]
\underline{Typical running time}: 
& This is strongly problem dependent:\\
&the test examples given took respectively\\
& 10,19 and 48 seconds on a SUN 4/310,\\
&SUN 4/110 and ACORN Springboard. \\[2ex]
\underline{Unusual Features of the Program}: 
&The REDUCE procedures use\\
&LISP vectors \cite{r2}
to provide a compact\\
&mathematical notation similar\\
& to that normally found in vector\\
& analysis textbooks.\\
\end{tabular}

\newpage
\section{Introduction}
The revised version of ORTHOVEC is, like the original \cite{r4}, a 
collection of REDUCE procedures and operators designed to simplify the
machine aided manipulation of vectors and vector expansions frequently
met in many areas of applied mathematics.  The revisions
have been introduced for two reasons: firstly, to add extra capabilities
missing from the original and secondly, to tidy up input and
output to make the package easier to use.

The changes from Version 1 include\\
i) merging of scalar and vector unary and binary operators, $+, - , *, /
$\\
ii) extensions of the definitions of division and exponentiation 
to vectors\\
iii) new vector dependency procedures\\
iv) application of L'H\^opital's rule in limits and Taylor expansions\\
v) a new component selector operator\\
vi) algebraic mode output of LISP vector components\\

The LISP vector primitives are again used to store vectors, although 
with the introduction of LIST types in algebraic mode in REDUCE
3.3 \cite{r1}, the implementation may have been more simply achieved 
using lists to store vector components.

The philosophy used in Version 2 follows that used in the original: 
namely, algebraic mode is used wherever possible.  The view 
is taken that some computational inefficiencies are 
acceptable if it allows coding to be 
intelligible to (and thence adaptable by)
users other than LISP experts familiar
with the internal workings of REDUCE.

Procedures and operators in ORTHOVEC fall into the five classes: 
initialisation, input-output, algebraic operations, differential
operations and integral operations.  Definitions are given in 
the following sections, and
a summary of the procedure names and their meanings are give in Table 1.
The final section discusses test examples.

\section{Initialisation}
The procedure VSTART initialises ORTHOVEC.  It may be
called after ORTHOVEC has been INputted (or LOADed if a fast load
version has been made) to reset coordinates.  VSTART provides a
menu of standard coordinate systems:-
\begin{enumerate}
\item cartesian $(x, y, z) = $ {\tt (x, y, z)}
\item cylindrical $(r, \theta, z) = $ {\tt (r, th, z)}
\item spherical $(r, \theta, \phi) = $ {\tt (r, th, ph) }
\item general $( u_1, u_2, u_3 ) = $ {\tt (u1, u2, u3) }
\item others
\end{enumerate}
which the user selects by number.  Selecting options (1)-(4) 
automatically sets up the coordinates and scale factors.  
Selection option
(5) shows the user how to select another coordinate system.
If VSTART is not called, then the default cartesian coordinates are
used.
ORTHOVEC may be re-initialised to a new coordinate system at any 
time during a given REDUCE session by typing 
\begin{verbatim}
VSTART $.
\end{verbatim}

\section{Input-Output}

ORTHOVEC assumes all quantities are either 
scalars or 3 component vectors.  
To define a vector $a$ with components $(c_1, c_2, c_3)$ use
the procedure SVEC as follows
\begin{verbatim}
a := svec(c1, c2, c3);
\end{verbatim}

The standard REDUCE output for vectors when using the terminator ``$;$''
is to list the three components inside square brackets
$[\cdots]$, with each component in prefix form.  A replacement for the
standard REDUCE procedure MAPRIN is included in
the package to change the 
output of LISP vector components to algebraic notation.  The procedure
VOUT (which returns the value of its argument)
can be used to give labelled output of components 
in algebraic form: e.g.,
\begin{verbatim}
b := svec (sin(x)**2, y**2, z)$
vout(b)$
\end{verbatim}

The operator {\tt \_} can be used to select a particular 
component (1, 2 or 3) for output e.g.
\begin{verbatim}
b_1 ;
\end{verbatim}

\section{Algebraic Operations}

Six infix operators, sum, difference, quotient, times, exponentiation
and cross product, and four prefix
operators, plus, minus, reciprocal
and  modulus are defined in ORTHOVEC.  These operators can take suitable
combinations of scalar and vector arguments,
and in the case of scalar arguments reduce to the usual definitions of
$ +, -, *, /, $ etc.

The operators are represented by symbols 
\begin{verbatim}  
+, -, /, *, ^, ><
\end{verbatim}

The composite {\tt ><} is an
attempt to represent the cross product symbol 
$\times$ in ASCII characters.
If we let ${\bf v}$ be a vector and $s$ be a scalar, then
valid combinations of arguments of the 
procedures and operators and the type of the result 
are as summarised below.  The notation used is\\
{\em result :=procedure(left argument, right argument) } or\\
{\em result :=(left operand) operator (right operand) } . \\

\underline{Vector Addition} \\
\begin{tabular}{rclcrcl}
{\bf v} &:=& VECTORPLUS({\bf v})  &{\rm or}& {\bf v} &:=&  +  {\bf v} \\
     s  &:=& VECTORPLUS(s)  &{\rm or} &      s  &:=&  +       s  \\
{\bf v} &:=& VECTORADD({\bf v},{\bf v})  &{\rm or }& {\bf v} &:=&  
{\bf v} +  {\bf v} \\
     s  &:=& VECTORADD(s,s)  &{\rm or }&      s  &:=&  s + s \\
\end{tabular} \\

\underline{Vector Subtraction} \\
\begin{tabular}{rclcrcl}
{\bf v} &:=& VECTORMINUS({\bf v})  &{\rm or}&
 {\bf v} &:=&  +  {\bf v} \\
 s  &:=& VECTORMINUS(s)  &{\rm or} &      s  &:=&  +       s  \\
{\bf v} &:=& VECTORDIFFERENCE({\bf v},{\bf v})  &{\rm or }& {\bf v} &:=&
  {\bf v} +  {\bf v} \\
 s  &:=& VECTORDIFFERENCE(s,s)  &{\rm or }&      s  &:=&  s + s \\
\end{tabular} \\

\underline{Vector Division}\\
\begin{tabular}{rclcrcl}
{\bf v} &:=& VECTORRECIP({\bf v})  &{\rm or}& {\bf v} &:=&  /  
{\bf v} \\
 s  &:=& VECTORRECIP(s)  &{\rm or} &      s  &:=&  /       s  \\
{\bf v} &:=& VECTORQUOTIENT({\bf v},{\bf v})  &{\rm or }& {\bf v} &:=&  
{\bf v} /  {\bf v} \\
{\bf v} &:=& VECTORQUOTIENT({\bf v},    s  )  &{\rm or }& {\bf v} &:=&  
{\bf v} /     s    \\
{\bf v} &:=& VECTORQUOTIENT(   s   ,{\bf v})  &{\rm or }& {\bf v} &:=&  
   s    /  {\bf v} \\
     s  &:=& VECTORQUOTIENT(s,s)  &{\rm or }&      s  &:=&  s / s       
       \\
\end{tabular} \\

\underline{Vector Multiplication}\\ 
\begin{tabular}{rclcrcl}
{\bf v} &:=& VECTORTIMES(   s   ,{\bf v})  &{\rm or }& {\bf v} &:=&     
s    *  {\bf v} \\
{\bf v} &:=& VECTORTIMES({\bf v},   s   )  &{\rm or }& {\bf v} &:=& {\bf
 v}  *     s    \\
   s    &:=& VECTORTIMES({\bf v},{\bf v})  &{\rm or }&    s    &:=& {\bf
 v}  *  {\bf v} \\
   s    &:=& VECTORTIMES(   s   ,   s   )  &{\rm or }&    s    &:=&     
s    *     s    \\
\end{tabular} \\
 
\underline{Vector Cross Product} \\
\begin{tabular}{rclcrcl}
{\bf v} &:=& VECTORCROSS({\bf v},{\bf v})  &{\rm or }& {\bf v} &:=& {\bf
 v} $\times$   {\bf v} \\
\end{tabular} \\
 
\underline{Vector Exponentiation}\\
\begin{tabular}{rclcrcl}
   s    &:=& VECTOREXPT ({\bf v},   s   )  &{\rm or }&    s    &:=& {\bf
 v}  \^~ s   \\
   s    &:=& VECTOREXPT (   s   ,   s   )  &{\rm or }&    s    &:=&    s
     \^~ s   \\
\end{tabular} \\

\underline{Vector Modulus}\\
\begin{tabular}{rcl}
   s    &:=& VMOD (s)\\
   s    &:=& VMOD ({\bf v}) \\
\end{tabular} \\

All other combinations of operands for these operators lead to error 
messages being issued .  The first two instances of vector
multiplication are scalar multiplication of vectors, the third is the
product of two scalars and the last is the inner (dot) product.  The
prefix operators  {\tt +, -, /} can take either scalar or vector
arguments and return results of the same type as their arguments. 
VMOD returns a scalar.

In compound expressions, parenthesis may be used to specify the order of 
combination.  If parentheses are omitted the ordering of the
operators, in increasing order of precedence is
\begin{verbatim}
+ | - | dotgrad | * | >< | ^ | _
\end{verbatim}
and these are placed in the precedence list defined in REDUCE \cite{r1}
 after $<$.
The differential operator 
DOTGRAD is defined in the
following section, and the component selector {\tt \_} was introduced in
section 3.

Vector divisions are defined as follows:  If ${\bf a}$ and ${\bf b}$ are
 vectors and $c$ is a scalar, then
\begin{eqnarray*}
{\bf a} /  {\bf b} & = &  \frac{{\bf a} \cdot {\bf b}}{  \mid {\bf b} 
\mid^2}\\
c / {\bf a}   & = &  \frac{c {\bf a}  }{ \mid {\bf a} \mid^2}
\end{eqnarray*}

Both scalar multiplication and dot products are given by the same symbol
, so braces are advisable to ensure the correct
precedences in expressions such as $({\bf a} \cdot {\bf b}) 
({\bf c} \cdot {\bf d})$.

Vector exponentiation is defined as the power of the modulus:\\
${\bf a}^n \equiv  {\rm VMOD}(a)^n =   \mid {\bf a} \mid^n$
\section{Differential Operations}
Differential operators provided are div, grad, curl, delsq, and dotgrad.
All but the last of these are prefix operators having a single
vector or scalar argument as appropriate.  Valid combinations of 
operator and argument, and the type of the result are as follows
\begin{tabular}{rcl}
s & := & div ({\bf v})  \\
{\bf v} & := & grad(s) \\
{\bf v} & := & curl({\bf v})  \\
{\bf v} & := & delsq({\bf v}) \\
 s  & := & delsq(s) \\
{\bf v} & := & {\bf v}  dotgrad {\bf v}  \\
 s & := & {\bf v}  dotgrad  s 
\end{tabular}

All other combinations of operator and argument type cause error
messages to be issued.  The differential operators have their usual
meanings \cite{r3}.  The coordinate system used by these operators is
set by invoking  VSTART (cf. Sec.1).  The names {\tt h1}, 
{\tt h2}  and {\tt h3 } are
reserved for the scale factors, and {\tt u1}, {\tt u2} and {\tt u3} are 
used for the coordinates.

A vector extension, VDF , of the REDUCE procedure DF allows the 
differentiation of a vector (scalar) with respect to a scalar to be
performed.  Allowed forms are
VDF({\bf v}, s)  $\rightarrow$  {\bf v}   and
VDF(s, s)  $\rightarrow$   s , 
where, for example\\
\begin{eqnarray*}
{\tt vdf( B,x)} \equiv \frac{\partial {\bf B}}{\partial x}
\end{eqnarray*}

The standard REDUCE procedures DEPEND and NODEPEND have been redefined
to allow dependences of vectors to be compactly
defined.  For example
\begin{verbatim}
a := svec(a1,a2,a3)$;
depend a,x,y;
\end{verbatim}
causes all three components {\tt a1},{\tt a2} and {\tt a3} of {\tt a} 
to be treated as functions of {\tt x} and {\tt y}.
Individual component dependences can still be defined if desired.
\begin{verbatim}
depend a3,z;
\end{verbatim}

The procedure VTAYLOR gives truncated Taylor series expansions of scalar
or vector functions:-
\begin{verbatim}
vtaylor(vex,vx,vpt,vorder);
\end{verbatim}
returns the series expansion of the expression 
VEX  with respect to variable VX
about point VPT  to order VORDER.  Valid
combinations of argument types are as follows:\\

\begin{tabular}{cccc}
VEX & VX & VPT &  VORDER \\[2ex]
{\bf v} & {\bf v} &  {\bf v} &  {\bf v}\\
{\bf v} &  {\bf v} & {\bf v} & s\\
{\bf v} & s & s & s \\
s & {\bf v} &  {\bf v} & {\bf v}   \\
s & {\bf v} & {\bf v} & s\\
s & s & s & s\\
\end{tabular}\\

Any other combinations cause error messages to be issued.  Elements of
VORDER must be non-negative integers, otherwise error
messages are issued.  If scalar VORDER is given for a vector expansion, 
expansions in each component are truncated at the same
order, VORDER.

The new version of Taylor expansion applies 
L'H\^opital's rule in evaluating coefficients, 
so handle cases such as $\sin(x) /  (x) $ , etc.
which the original version of ORTHOVEC could not. The procedure used for
this is LIMIT,
which can be used directly to find the limit of a scalar 
function {\tt ex} of variable {\tt x} at point {\tt pt}:-

\begin{verbatim}
ans := limit(ex,x,pt);
\end{verbatim}

\section{Integral Operations}
Definite and indefinite vector, volume and scalar line 
integration procedures are included in ORTHOVEC.  
They are defined as follows:
\begin{eqnarray*}
{\rm VINT} ({\bf v},x)  & = & \int {\bf v}(x)dx\\
%
{\rm DVINT} ({\bf v},x, a, b) & = & \int^b_a {\bf v} (x) dx\\
%
{\rm VOLINT} ({\bf v}) & = & \int {\bf v} h_1 h_2 h_3  du_1 du_2 du_3\\
%
{\rm DVOLINT}({\bf v},{\bf l},{\bf u},n) & = & \int^{\bf u}_{\bf l} {\bf
 v} h_1 h_2 h_3 du_1 du_2 du_3\\
%
{\rm LINEINT} ({\bf v, \omega}, t) & = & \int {\bf v} \cdot {\bf dr} 
\equiv \int v_i h_i \frac{\partial \omega_i}{\partial t} dt\\
%
{\rm DLINEINT} ({\bf v, \omega} t, a, b) & = & \int^b_a v_i h_i 
\frac{\partial \omega_i}{\partial t} dt\\
\end{eqnarray*}

In the vector and volume integrals, ${\bf v}$ are vector or scalar,  $a,
b,x$ and $n$ are scalar.  Vectors ${\bf l}$ and ${\bf u}$ contain
expressions for lower and upper bounds to the integrals.  The integer 
index $n$ defines the order in which the integrals over $u_1, u_2$
and $u_3$ are performed in order to allow for functional dependencies in
 the integral bounds:
 
\begin{tabular}{ll}
n & order\\
1 & $u_1~u_2~u_3$\\
%
2 & $u_3~u_1~u_2$\\
%
3 & $u_2~u_3~u_1$\\
%
4 & $u_1~u_3~u_2$\\
%
5 & $u_2~u_1~u_3$\\
otherwise & $u_3~u_2~u_1$\\
\end{tabular}

The vector ${\bf \omega}$ in the line 
integral's arguments contain explicit paramterisation 
of the coordinates $u_1, u_2, u_3$ of the
line ${\bf u}(t)$ along which the integral is taken.

\begin{table}
\centering
\caption{Procedures names and operators used in ORTHOVEC version 2}
\vspace*{5mm}
\begin{tabular}{|l c l|} \hline 
Procedures & & Description \\ \hline 
VSTART & & select coordinate system \\ & & \\ 
SVEC   & & set up a vector \\ 
VOUT & & output a vector \\
VECTORCOMPONENT & \_ & extract a vector component (1-3) \\
 & & \\
VECTORADD & + & add two vectors or scalars \\ 
VECTORPLUS & + & unary vector or scalar plus\\ 
VECTORMINUS & - & unary vector or scalar minus\\ 
VECTORDIFFERENCE & - & subtract two vectors or scalars \\ 
VECTORQUOTIENT & / & vector divided by scalar \\ 
VECTORRECIP & / & unary vector or scalar division (reciprocal)\\ 
VECTORTIMES & * & multiply vector or scalar by vector or scalar \\ 
VECTORCROSS & $><$ & cross product of two vectors \\ 
VECTOREXPT & \^~ & exponentiate vector modulus or scalar \\
VMOD & & length of vector or scalar \\ & & \\ 
DIV & & divergence of vector \\ 
GRAD & & gradient of scalar \\ 
CURL & & curl of vector \\ 
DELSQ & & laplacian of scalar or vector \\ 
DOTGRAD & & (vector).grad(scalar or vector) \\ & &  \\ 
VTAYLOR & & vector or scalar Taylor series of vector or scalar \\ 
VPTAYLOR & & vector or scalar Taylor series of scalar \\ 
TAYLOR & & scalar Taylor series of scalar \\
LIMIT & & limit of quotient using L'H\^opital's rule \\ & &  \\
VINT & & vector integral \\ 
DVINT & & definite vector integral \\ 
VOLINT & & volume integral \\ 
DVOLINT & & definite volume integral \\ 
LINEINT & & line integral \\ 
DLINEINT & & definite line integral \\  & & \\
MAPRIN & & vector extension of REDUCE MAPRIN \\
DEPEND & & vector extension of REDUCE DEPEND \\
NODEPEND & & vector extension of REDUCE NODEPEND \\ \hline 
\end {tabular}
\end{table}

\section{Test Cases}

To use the REDUCE source version of ORTHOVEC, initiate a REDUCE session 
and then IN the file {\em orthovec.red}
containing ORTHOVEC.  However, it is recommended that for 
efficiency a compiled fast loading version be 
made and LOADed when required (see 
Sec.~18 of the REDUCE manual \cite{r1}).  
If coordinate dependent differential and
integral operators other than cartesian are needed, then VSTART must
be used to reset coordinates and scale factors.

Six  simple examples are given in the Test Run Output file 
{\em orthovectest.log} to illustrate the working of ORTHOVEC.  
The input lines were taken from the file
{\em orthovectest.red} (the Test Run Input), but could 
equally well be typed in at the Terminal.

\subsection{Example 1} 

Show that
\begin{eqnarray*}
({\bf a}  \times {\bf b}) \cdot ({\bf c} \times {\bf d}) - ({\bf a} 
\cdot {\bf c})({\bf b} \cdot {\bf d})
 + ({\bf a} \cdot {\bf d})({\bf b} \cdot {\bf c}) \equiv 0
\end{eqnarray*}
 
\subsection{Example 2}

Write the equation of motion
\begin{eqnarray*}
\frac{\partial {\bf v}}{\partial t} + {\bf v} \cdot {\bf \nabla v} 
+ {\bf \nabla} p - curl ({\bf B}) \times {\bf B}
\end{eqnarray*}
in cylindrical coordinates.

\subsection{Example 3}

Taylor expand
\begin{itemize}
\item $\sin(x) \cos(y) +e^z$
about the point $(0,0,0)$ to third order in $x$, fourth order in $y$ and
fifth order in $z$.

\item $\sin(x)/x$ about $x$ to fifth order.

\item ${\bf v}$ about ${\bf x}=(x,y,z)$ to fifth order, where
${\bf v} = (x/ \sin(x),(e^y-1)/y,(1+z)^{10})$.
\end{itemize}

\subsection{Example 4}
Obtain the second component of the equation of motion in example 2, and 
the first 
component of the final vector Taylor series in example 3.

\subsection{Example 5}

Evaluate the line integral 
\begin{eqnarray*}
\int^{{\bf r}_2}_{{\bf r}_1} {\bf A} \cdot d{\bf r}
\end{eqnarray*}
from point ${\bf r}_1 = (1,1,1)$ to point
${\bf r}_2 = (2,4,8)$ along the path $(x,y,z) = (s, s^2, s^3)$ where\\
\begin{eqnarray*}
{\bf A} = (3x^2 + 5y) {\bf i} - 12xy{\bf j} + 2xyz^2{\bf k}
\end{eqnarray*}
and $({\bf i, j, k})$ are unit vectors in the ($x,y,z$) directions.

\subsection{Example 6}

Find the volume $V$ common to the intersecting cylinders $x^2 + y^2 
= r^2$ and $x^2 + z^2 = r^2$ i.e. evaluate
\begin{eqnarray*}
V = 8 \int^r_0 dx \int^{ub}_0 dy \int^{ub}_0 dz
\end{eqnarray*}
where $ub = \overline{\sqrt { r^2 - x^2}}$

\begin{thebibliography}{9}

\bibitem{r4} J.~W.~Eastwood,{\em Computer Phys.~Commun.~}\underline{47} 
(1987) 139-147.

\bibitem{r1} A. C. Hearn,{\em REDUCE User's Manual}, 
(Rand Publication CP78 (Rev. 7/87), 1987.

\bibitem{r2} J. B. Marti, A. S. Hearn, M. L. Griss and C. Griss,\\
{\em Standard LISP Report}, 
ACM SIGPLAN Notices \underline{14} (1979) 48.

\bibitem{r3} M . Speigel,{\em Vector Analysis}, (Scheum Publishing Co., 
1959).

\end{thebibliography}

\newpage
\appendix
\begin{center}
{\Large{\bf Appendices}}
\end{center}
\section {Test Run Input}
$\bullet$ Insert file orthovectest.red

\section {Test Run Output}
$\bullet$ Insert file orthovectest.log

\section {Source Listing}
$\bullet$ Insert file orthovec.red

\end{document}