A Vector Algebra and Calculus Package for REDUCE

1 A Vector Algebra and Calculus Package forREDUCED avid HarperAstronomy UnitQueen Mary and Westfield CollegeUniversity of LondonMile End RoadLondon E1 4 NSEnglandElectronic IntroductionThis package1is written in RLISP (the LISP meta-language) and is intendedfor use with REDUCE It provides REDUCE with the ability to performvector Algebra using the same notation as scalar Algebra . The basic algebraicoperations are supported, as are differentiation and integration of vectorswith respect to scalar variables, cross product and dot product, componentmanipulation and application of scalar functions ( ) to a Vector toyield a Vector set of Vector Calculus operators are provided for use with any orthogonalcurvilinear coordinate system.

2 These operators are gradient, divergence,curl and del-squared (Laplacian). The Laplacian operator can take scalaror Vector important coordinate systems are pre-defined and can be invokedby name. It is also possible to create new coordinate systems by specifyingthe names of the coordinates and the values of the scale : Computer Physics Communications,54, 295-305 (1989)12 Vector DECLARATION AND INITIALISATION22 Vector declaration and initialisationAny name may be declared to be a Vector , provided that it has not previouslybeen declared as a matrix or an array. To declare a list of names to be vectorsuse the VEC command:VEC A,B,C;declares the variablesA,BandCto be vectors. If they have already beenassigned (scalar) values, these will be a Vector is declared using theVECcommand, it does not have aninitial a Vector value is assigned to a scalar variable, then that variable willautomatically be declared as a Vector and the user will be notified that thishas Vector may be initialised using theAVEC function which takes three scalararguments and returns a Vector made up from those scalars.

3 For exampleA := AVEC(A1, A2, A3);sets the components of the vectorAtoA1, Vector Algebra (In the examples which follow,V,V1,V2etcare assumed to be vectors whileS,S1,S2etc are scalars.)The scalar Algebra operators +,-,* and / may be used with Vector operandsaccording to the rules of Vector Algebra . Thus multiplication and division ofa Vector by a scalar are both allowed, but it is an error to multiply or divideone Vector by := V1 + V2 - V3;Addition and subtractionV := S1*3*V1;Scalar multiplicationV := V1/S;Scalar divisionV := -V1;NegationVector multiplication is carried out using the infix are defined to have higher precedence than scalar multiplication Vector ALGEBRA3V := V1 CROSS V2;Cross productS := V1 DOT V2;Dot productV := V1 CROSS V2 + V3;V := (V1 CROSS V2) + V3;The last two expressions are equivalent due to the precedence of modulus of a Vector may be calculated using := VMOD V;A unit Vector may be generated from any Vector using := V/(VMOD V).

4 Components may be extracted from any Vector using index notation in thesame way as an := AVEC(AX, AY, AZ);V(0);yields AXV(1);yields AYV(2);yields AZIt is also possible to set values of individual components. Following fromabove:V(1) := B;The vectorVnow has componentsAX,B, may be used as arguments in the differentiation and integrationroutines in place of the dependent := AVEC(X**2, SIN(X), Y);DF(V,X);yields (2*X, COS(X), 0)INT(V,X);yields (X**3/3, -COS(X), Y*X)Vectors may be given as arguments to monomial functions such asSIN,LOGandTAN. The result is a Vector obtained by applying the function component-wise to the argument := AVEC(A1, A2, A3);SIN(V);yields (SIN(A1), SIN(A2), SIN(A3))4 Vector CALCULUS44 Vector calculusThe Vector Calculus operators div, grad and curl are recognised.

5 The Lapla-cian operator is also available and may be applied to scalar and := GRAD S;Gradient of a scalar fieldS := DIV V;Divergence of a Vector fieldV := CURL V1;Curl of a Vector fieldS := DELSQ S1;Laplacian of a scalar fieldV := DELSQ V1;Laplacian of a Vector fieldThese operators may be used in any orthogonal curvilinear coordinate sys-tem. The user may alter the names of the coordinates and the values of thescale factors. Initially the coordinates areX,YandZand the scale factorsare all are two special vectors :COORDS contains the names of the coordinatesin the current system andHFACTORS contains the values of the scale coordinate names may be changed using R,THETA,PHI;This command changes the coordinate names toR, scale factors may be altered using (1,R,R*SIN(THETA));This command changes the scale factors to1,RandR SIN(THETA).

6 Note that the arguments ofSCALEFACTORS must be enclosed in is not necessary Vector differential operators are applied to an expression, the currentset of coordinates are used as the independent variables and the scale factorsare employed in the calculation. (See, for example, Batchelor AnIntroduction to Fluid Mechanics , Appendix 2.)Several coordinate systems are pre-defined and may be invoked by see a list of valid names enterSYMBOLIC !*CSYSTEMS;and REDUCE will respond with something like4 Vector CALCULUS5(CARTESIAN SPHERICAL CYLINDRICAL)To choose a coordinate system by name, use the choose the Cartesian coordinate system :GETCSYSTEM CARTESIAN;Note the quote which prefixes the name of the coordinate system.

7 This is re-quired becauseGETCSYSTEM(and its complementPUTCSYSTEM) is aSYMBOLIC procedure which requires a literal responds by typing a list of the coordinate names in that coordi-nate system. The example above would produce the response(X Y Z)whilstGETCSYSTEM SPHERICAL;would produce(R THETA PHI)Note that any attempt to invoke a coordinate system is subject to the samerestrictions as the implied calls toCOORDINATESandSCALEFACTORS. In par-ticular,GETCSYSTEM fails if any of the coordinate names has been assigneda value and the previous coordinate system remains in user-defined coordinate system can be assigned a name using the commandPUTCSYSTEM. It may then be re-invoked at a later stage 1We define a general coordinate system with coordinate namesX,Y,Zand scalefactorsH1,H2,H3:COORDINATES X,Y,Z;SCALEFACTORS(H1,H2,H3);PUTCSYSTEM GENERAL;This system may later be invoked by enteringGETCSYSTEM GENERAL;5 VOLUME AND LINE INTEGRATION65 Volume and Line IntegrationSeveral functions are provided to perform volume and line integrals.

8 Theseoperate in any orthogonal curvilinear coordinate system and make use ofthe scale factors described in the previous integrals of scalar and Vector expressions may be calculated 2To calculate the definite integral of sin(x)2between 0 and 2 we enterDEFINT(SIN(X)**2,X,0,2*PI);This function is a simple extension of theINTfunction taking two extraarguments, the lower and upper bounds of integration volume integrals may be calculated using theVOLINTEGRAL functionwhose syntax is as follows :VOLINTEGRAL(integrand, vectorlower-bound, vectorupper-bound);Example 3In spherical polar coordinates we may calculate the volume of a sphere byintegrating unity over the ranger=0 toRR, =0 toPI, =0 to 2* as follows:VLB := AVEC(0,0,0);Lower boundVUB := AVEC(RR,PI,2*PI);Upper bound inr, , respectivelyVOLINTORDER := (0,1,2);The order of integrationVOLINTEGRAL(1,VLB,VUB);Note the use of the special vectorVOLINTORDER which controls the orderin which the integrations are carried out.

9 This Vector should be set tocontain the number 0, 1 and 2 in the required order. The first component ofVOLINTORDER contains the index of the first integration variable, the secondcomponent is the index of the second integration variable and the thirdcomponent is the index of the third integration 4 Suppose we wish to calculate the volume of a right circular cone. This isequivalent to integrating unity over a conical region with the bounds:5 VOLUME AND LINE INTEGRATION7z = 0 to H(H = the height of the cone)r = 0 to pZ(p = ratio of base diameter to height)phi = 0 to 2*PIWe evaluate the volume by integrating a series of infinitesimally thin circulardisks of constant z-value. The integration is thus performed in the order :d( ) from 0 to 2 , dr from 0 to p*Z, dz from 0 to H.

10 The order of the indicesis thus 2, 0, := AVEC(2,0,1);VLB := AVEC(0,0,0);VUB := AVEC(P*Z,H,2*PI);VOLINTEGRAL(1,VLB,VUB); (At this stage, we replaceP*HbyRR, the base radius of the cone, to obtainthe result in its more familiar form.)Line integrals may be calculated using theLINEINTandDEFLINEINT func-tions. Their general syntax isLINEINT( Vector -function, Vector -curve,v ariable);DEFLINENINT( Vector -function,vec tor-curve,variable,lower-bound,upper-bou nd);wherevector-functionis any Vector -valued expression; Vector -curveis a Vector expression which describes the path of integrationin terms of the independent variable;variableis the independent variable;lower-boundupper-boundare the bounds of integration in terms of the 5In spherical polar coordinates, we may integrate round a line of constanttheta ( latitude ) to find the length of such a line.