20.4 AVECTOR: A Vector Algebra and Calculus Package

This package provides REDUCE with the ability to perform vector algebra using the same notation as scalar algebra. The basic algebraic operations are supported, as are differentiation and integration of vectors with respect to scalar variables, cross product and dot product, component manipulation and application of scalar functions (e.g. cosine) to a vector to yield a vector result.

Author: David Harper.

20.4.1 Introduction

This package ([Har89]) provides REDUCE with the ability to perform vector algebra using the same notation as scalar algebra. The basic algebraic operations are supported, as are differentiation and integration of vectors with respect to scalar variables, cross product and dot product, component manipulation and application of scalar functions, e.g., cosine) to a vector to yield a vector result.

A set of vector calculus operators are provided for use with any orthogonal curvilinear coordinate system. These operators are gradient, divergence, curl and del-squared (Laplacian). The Laplacian operator can take scalar or vector arguments.

Several important coordinate systems are pre-defined and can be invoked by name. It is also possible to create new coordinate systems by specifying the names of the coordinates and the values of the scale factors.

20.4.2 Vector declaration and initialisation

Any name may be declared to be a vector, provided that it has not previously been declared as a matrix or an array. To declare a list of names to be vectors use the vec command:

   vec a,b,c;

declares the variables a, b and c to be vectors. If they have already been assigned (scalar) values, these will be lost.

When a vector is declared using the vec command, it does not have an initial value.

If a vector value is assigned to a scalar variable, then that variable will automatically be declared as a vector and the user will be notified that this has happened.

A vector may be initialised using the avec function which takes three scalar arguments and returns a vector made up from those scalars. For example

   a := avec(a1, a2, a3);

sets the components of the vector a to a1, a2 and a3.

20.4.3 Vector algebra

(In the examples which follow, v, v1, v2, etc. are assumed to be vectors while s, s1, s2, etc. are scalars.)

The scalar algebra operators +,-,* and / may be used with vector operands according to the rules of vector algebra. Thus multiplication and division of a vector by a scalar are both allowed, but it is an error to multiply or divide one vector by another.

Vector multiplication is carried out using the infix operators dot and cross. These are defined to have higher precedence than scalar multiplication and division.

The last two expressions are equivalent due to the precedence of the cross operator.

The modulus of a vector may be calculated using the VMOD operator.

   s := vmod v;

A unit vector may be generated from any vector using the vmod operator.

   v1 := v/(vmod v);

Components may be extracted from any vector using index notation in the same way as an array.

It is also possible to set values of individual components. Following from above:

   v(1) := b;

The vector v now has components ax, b, az.

Vectors may be used as arguments in the differentiation and integration routines in place of the dependent expression.

Vectors may be given as arguments to monomial functions such as sin, log and tan. The result is a vector obtained by applying the function component-wise to the argument vector.

20.4.4 Vector calculus

The vector calculus operators div, grad and curl are recognised. The Laplacian operator is also available and may be applied to scalar and vector arguments.

These operators may be used in any orthogonal curvilinear coordinate system. The user may alter the names of the coordinates and the values of the scale factors. Initially the coordinates are x, y and z and the scale factors are all unity.

There are two special vectors : coords contains the names of the coordinates in the current system and hfactors contains the values of the scale factors.

The coordinate names may be changed using the coordinates command.

   coordinates r,theta,phi;

This command changes the coordinate names to r, theta and phi.

The scale factors may be altered using the scalefactors operator.

   scalefactors(1,r,r*sin(theta));

This command changes the scale factors to 1, r and r sin(theta).

Note that the arguments of scalefactors must be enclosed in parentheses. This is not necessary with the coordinates command.

When vector differential operators are applied to an expression, the current set of coordinates are used as the independent variables and the scale factors are employed in the calculation. (See, for example, Batchelor G.K. ’An Introduction to Fluid Mechanics’, Appendix 2.)

Several coordinate systems are pre-defined and may be invoked by name. To see a list of valid names enter

   symbolic !*csystems;

and REDUCE will respond with something like

   (cartesian spherical cylindrical)

To choose a coordinate system by name, use the command getcsystem.

To choose the Cartesian coordinate system :

   getcsystem ’cartesian;

Note the quote which prefixes the name of the coordinate system. This is required because getcsystem (and its complement putcsystem) is a symbolic procedure which requires a literal argument.

REDUCE responds by typing a list of the coordinate names in that coordinate system. The example above would produce the response

   (x y z)

whilst

   getcsystem ’spherical;

would produce

   (r theta phi)

Note that any attempt to invoke a coordinate system is subject to the same restrictions as the implied calls to coordinates and SCALEFACTORS. In particular, getcsystem fails if any of the coordinate names has been assigned a value and the previous coordinate system remains in effect.

A user-defined coordinate system can be assigned a name using the command putcsystem. It may then be re-invoked at a later stage using getcsystem.

Example 1

We define a general coordinate system with coordinate names X,Y,Z and scale factors h1,h2,h3 :

   coordinates x,y,z;
   scalefactors(h1,h2,h3);
   putcsystem ’general;

This system may later be invoked by entering

   getcsystem ’general;

20.4.5 Volume and Line Integration

Several functions are provided to perform volume and line integrals. These operate in any orthogonal curvilinear coordinate system and make use of the scale factors described in the previous section.

Definite integrals of scalar and vector expressions may be calculated using the defint function.

Example 2

To calculate the definite integral of \(\sin (x)^2\) between 0 and 2\(\pi \) we enter

   defint(sin(x)**2,x,0,2*pi);

This function is a simple extension of the int operator taking two extra arguments, the lower and upper bounds of integration respectively.

Definite volume integrals may be calculated using the volintegral function whose syntax is as follows :

Example 3

In spherical polar coordinates we may calculate the volume of a sphere by integrating unity over the range \(r\)=0 to rr, \(\theta \)=0 to \(\pi \), \(\phi \)=0 to 2*\(\pi \) as follows :

Note the use of the special vector volintorder which controls the order in which the integrations are carried out. This vector should be set to contain the number 0, 1 and 2 in the required order. The first component of volintorder contains the index of the first integration variable, the second component is the index of the second integration variable and the third component is the index of the third integration variable.

Example 4

Suppose we wish to calculate the volume of a right circular cone. This is equivalent to integrating unity over a conical region with the bounds:

We evaluate the volume by integrating a series of infinitesimally thin circular disks of constant z-value. The integration is thus performed in the order : d(\(\phi \)) from 0 to \(2\pi \), dr from 0 to p*Z, dz from 0 to H. The order of the indices is thus 2, 0, 1.

   volintorder := avec(2,0,1);
   vlb := avec(0,0,0);
   vub := avec(p*z,h,2*pi);
   volintegral(1,vlb,vub);

(At this stage, we replace p*h by rr, the base radius of the cone, to obtain the result in its more familiar form.)

Line integrals may be calculated using the lineint and deflineint operators. Their general syntax is

where

\(\langle \)vector-function\(\rangle \)

is any vector-valued expression;

\(\langle \)vector-curve\(\rangle \)

is a vector expression which describes the path of integration in terms of the independent variable;

\(\langle \)variable\(\rangle \)

is the independent variable;

\(\langle \)lower-bound\(\rangle \)

\(\langle \)upper-bound\(\rangle \)

are the bounds of integration in terms of the independent variable.

Example 5

In spherical polar coordinates, we may integrate round a line of constant theta (‘latitude’) to find the length of such a line. The vector function is thus the tangent to the ‘line of latitude’, (0,0,1) and the path is (0,lat,phi) where phi is the independent variable. We show how to obtain the definite integral, i.e. from \(\phi =0\) to \(2 \pi \) :

deflineint(avec(0,0,1),avec(0,lat,phi),phi,0,2*pi);

20.4.6 Defining new functions and procedures

Most of the procedures in this package are defined in symbolic mode and are invoked by the REDUCE expression-evaluator when a vector expression is encountered. It is not generally possible to define procedures which accept or return vector values in algebraic mode. This is a consequence of the way in which the REDUCE interpreter operates and it affects other non-scalar data types as well : arrays cannot be passed as algebraic procedure arguments, for example.

20.4.7 Acknowledgements

This package was written whilst the author was the U.K. Computer Algebra Support Officer at the University of Liverpool Computer Laboratory.

Hosted by