Application Center - Maplesoft

# The VectorCalculus Package

You can switch back to the summary page by clicking here.

The VectorCalculus Package in Maple 8

Maple 8 provides a new package called VectorCalculus  for computing with vectors, vector fields, multivariate functions and parametric curves.  Computations include vector arithmetic using basis vectors, "del"-operations, multiple integrals over regions and solids, line and surface integrals, differential-geometric properties of curves, and many others.

You can easily perform these computations in any coordinate system and convert results between coordinate systems. The package is fully compatible with the Maple LinearAlgebra  package.  You can also extend the VectorCalculus  package by defining your own coordinate systems.

 > restart: with(VectorCalculus);

```Warning, the assigned names <,> and <|> now have a global binding
```

```Warning, these protected names have been redefined and unprotected: *, +, ., Vector, diff, int, limit, series
```

 >

Creating Vectors and Vector Fields

It's usually convenient to specify a default coordinate system using the SetCoordinates  function.  (For a listing of the 37 coordinate systems predefined in Maple, see the help page ?Coordinates).

 > SetCoordinates( cartesian[x, y, z] );

You can define vectors in the same way you do with Maple's LinearAlgebra  package, using either the Vector  constructor, or <>  notation.
Note the output appears in terms of the cartesian basis vectors
,  and .  Vectors defined this way are recognized by the LinearAlgebra  pakcage.

 > v := Vector( [a, b, c] );

The MapToBasis  function converts vectors and vector fields to other coordinate systems.

 > MapToBasis( v, cylindrical[r, theta, z] );

A vector field is a data type in the VectorCalculus  package. The bars over the basis vectors distinguish vector fields from constant vectors.

 > F := VectorField( );

 > attributes( F );

Evaluate the vector field F  at a particular vector, say , using the evalVF  function

 > evalVF( F, <4, b, 1> );

Change the default coordinate system to polar.

 > SetCoordinates( polar[r, theta] );

Define a vector field in the polar coordinate system.

 > F := VectorField( <0, 1> );

Below we see that the value of the basis vector  is not constant over the plane.  This is a major difference from cartesian coordinates, in which the basis vectors  and  are constants over the plane (i.e.  and .)  You can use Maple to clarify this common point of confusion for students.

 > MapToBasis( F, cartesian[x, y]);

 >

Operations on Vectors and Vector Fields

Basic Arithmetic

The package is not limited to three dimensions.  The cartesian  option without variables specified allows you to compute with vectors of any dimension.

 > SetCoordinates( cartesian );

Vector addition and scalar multiplication over

 > v1 := a * <1, 2, 5, -1> + b * <3, 11, 0, 8>;

 > v2 := 2*v1 - v1/3;

 >

Dot Product, Divergence and Del

The .  operator and the DotProduct  function are equivalent ways to compute vector dot products.

 > . ;

 > DotProduct( , );

You can compute the divergence of a vector field using either the Divergence  function, or the dot product with the Del  or Nabla  operators.

 > F := VectorField( , cartesian[x,y,z] );

 > Divergence( F );

 > Del . F;

 > Nabla . F;

 >

Cross Product, Curl and Del

To compute the cross product of two vectors in , use the CrossProduct  function or the &x  operator.

 > CrossProduct( ,  );

 > &x ;

Compute the curl of a vector field F  using either the Curl  function or the expression Del &x F.

Note that Maple distinguishes the vector  0  from the number 0 by writing a basis-vector next to the 0.

 > Curl( VectorField( , cartesian[x,y,z] ) );

Here's a nonconservative (rotational) vector field.

 > F := VectorField( , cartesian[x,y,z] );

We know that it's nonconservative because its curl is nonzero.

 > Del &x F;

Routines from the plots  package can accept outputs from the VectorCalculus  package.

 > plots[fieldplot3d]( F, x=-1..1, y=-1..1, z=-1/2..1/2, scaling=constrained, shading=zhue );

Multivariate Differential Calculus

Earlier versions of Maple could already perform multivariate differential calculus, but the VectorCalculus  package makes it especially simple and intuitive.  Working with gradients, laplacians and hessians of real-valued functions is straightforward, as the following worked examples show.

Laplacians, Gradients, and Hessians, Oh My!

 > SetCoordinates( cartesian ):

Example 1: Find the unit vector normal to the surface  at the point

We'll solve this problem as students are taught.  Introduce the 3-variable function  and treat the surface as the level surface

 > f := -z + x^2*y^2 + y + 1;

 > p := <0,0,1>;

The normal to the surface at p  is just (grad f)(p) , since the gradient of f  is normal to f' s level surfaces.  We compute the gradient with the function Gradient ...

 > Gf := Gradient( f, cartesian[x,y,z]);

...and evaluate the gradient at p  with the evalVF  function.

 > Gf0 := evalVF( Gf, p );

Finally, we normalize the gradient vector with the Normalize  function from the LinearAlgebra  package.

 > Gf0_unit := LinearAlgebra :- Normalize(Gf0, 2);

A plot of the result.

 > a := plots[arrow]( p, Gf0_unit, color=blue, thickness=3): S := plot3d( f+z, x=-1..1, y=-1..1, style=wireframe, axes=boxed): plots[display](S, a, scaling=constrained, orientation=[-20,65]);

 >

Example 2: Find the point on the paraboloid  closest to the point

 > SetCoordinates(cartesian[x,y,z]):

We'll solve this as a constrained minimization problem using the Lagrange Multiplier Theorem.
Interpret the paraboloid as a level surface of the function
.

 > g := (x^2+y^2)/5-z;

 > p := <1,1,0>;

Let the objective function f  be the squared distance from point p  to a point .  We wish to minimize f .

 > f := (x-p[1])^2 + (y-p[2])^2 + (z-p[3])^2;

Form the Lagrangian for the constrained minimization problem,

 > L := f - lambda * g;

Compute the gradient of the Lagrangian.

 > GL := Gradient( L );

Solve the four equations g = 0, GL = 0   for the four unknowns x, y, z   and

 > sol := fsolve( {g = 0, seq( GL[i] = 0, i=1..3 )},                {x,y,z,lambda} ) ;

Verify that the solution is a minimum using the generalized 2nd-derivative test, i.e.,

 > HL := Hessian( L, [x,y,z,lambda]);

 > LinearAlgebra :- Determinant(HL);

The determinant is negative at the solution found, so the solution is a minimum.

 > eval( %, sol );

A plot of the result.

 > P := plot3d( g+z, x=0..2, y=0..2, style=wireframe, axes=boxed ): L := plottools[line]( [1,1,0], eval([x,y,z], sol), thickness=3, color=red): plots[display]( P, L, scaling=constrained, orientation=[30,85]);

 >

Example 3: Solve the time-dependent heat equation for a cylindrical heat source with external heat generation

We find a general solution to the time-dependent heat equation

in cylindrical coordinates, where
is the temperature, q   is the thermal source strength, and k  is the thermal diffusivity.

The Laplacian  function is convenient for writing down the PDEs of classical physics, especially in non-cartesian coordinate systems.

 > Lap := Laplacian( u(r, theta, z, t), cylindrical[r,theta,z]);

 > HeatEqn := Lap + q/k = diff(u(r,theta,z, t), t);

Solve the partial differential equation.

 > pdsolve( HeatEqn, build );

 >

Multivariate Integral Calculus

The VectorCalculus  package greatly simplifies integration of both real- and vector-valued functions over geometric regions in  and , such as rectangles, parallelepipeds, circles, ellipses, spheres, triangles, tetrahedra, and regions bounded by curves and surfaces.

Definite Integrals of Multivariate Functions over Regions and Solids

VectorCalculus  overloads Maple's int  routine to include options for specifying the geometry of the region of integration.

Compute  over the region

 > int( sin(x)*cos(y)*tan(z), [x, y, z] = Parallelepiped( 0..Pi, 0..Pi/2, 0..Pi/4 ) );

Compute  over the circle with center  and radius r , and  over all of

 > int( exp(-x^2-y^2), [x, y] = Circle( <0, 0>, r ) ), int( exp(-x^2-y^2-z^2), [x, y, z] = Sphere( <0, 0, 0>, infinity ) );

Find the volume of the solid bounded by the paraboloid , the plane , and the surface .  We use the Region  option

 > int( 1, [x,y,z] = Region( 0..1, x^2..x, 0..x^2+y^2 ) );

Compute  where S  is the tetrahedron with vertices

 > int( x+y+z, [x,y,z] = Tetrahedron( <0,0,0>, <1,0,0>, <0,1,0>, <0,0,1> ) );

Find the area of the section of the ellipse  from the angle  through the angle

 > int( x, [x,y] = Sector( Ellipse( x^2/4 + y^2/9 - 1 ), alpha, beta ) );

 >

Line and Surface Integrals of Vector-Valued Functions

VectorCalculus  has top-level functions for path, line and surface integrals.  You can specify the paths of integration with options such as LineSegments , Path , and Arc .

 > SetCoordinates( cartesian[x,y] );

Compute the line integral  for  over the unit square, demonstrating that F  is a nonconservative vector field.

 > LineInt( VectorField( <-x-y, x-y> ), LineSegments( <0,0>, <1,0>, <1,1>, <0,1>, <0,0> ) );

Compute the line integral  for  over the path  as t  ranges from 0 to 2 .

 > LineInt( VectorField( ), Path( , t=0..2 ) );

Compute the line integral  for  over the ellipse   in the first quadrant .

 > LineInt( VectorField( ), Arc( Ellipse( x^2/4+y^2/9-1 ), 0, Pi/2 ) );

 > restart: with(VectorCalculus):

```Warning, the assigned names <,> and <|> now have a global binding
```

```Warning, these protected names have been redefined and unprotected: *, +, ., Vector, diff, int, limit, series
```

Compute the surface area of a sphere

 > SurfaceInt( 1, [x,y,z] = Surface( , s=0..Pi, t=0..2*Pi, coords=spherical ) ) assuming r>0;

Compute the surface area of the plane  over the triangle in the xy-plane that has vertices

 > SurfaceInt( 1, [x,y,z] = Surface( , [s,t] = Triangle(<0,0>,<1,0>,<0,1>) ) );

 >

Operations on Parametric Curves

The VectorCalculus  package provides tools to analyse the differential-geometric properties of curves, such as torsion, normal vectors, radius of curvature, TNB frames, and others.

Example: Draw the evolute of an ellipse

The evolute  of a curve is the locus of its radius of curvature.  To compute the evolute of an ellipse, we first express the ellipse as a parametric curve.

 > c := <3*cos(t)/2, sin(t)>;

 > assume( t::real ): interface( showassumed=0 ):

Compute the normal vector of the ellipse as a function of t .

 > n := simplify( PrincipalNormal( c, t ) );

Scale the normal vector to a unit vector by dividing it by its length.

 > n_unit := simplify( n / LinearAlgebra :- Norm( n, 2) );

Compute the radius of curvature of the ellipse

 > r := RadiusOfCurvature( c, t );

The evolute is then the normal vector scaled by the radius of curvature at each value of t .

 > evolute := c + r * n_unit;

 > plot([[c[1], c[2], t=0..2*Pi], [evolute[1], evolute[2], t=0..2*Pi]], scaling=constrained);

Extending the VectorCalculus Package

You can add your own coordinate systems to the database of predefined systems in VectorCalculus using the AddCoordinates function .

 > restart: with(VectorCalculus): interface(showassumed=0):

```Warning, the assigned names <,> and <|> now have a global binding
```

```Warning, these protected names have been redefined and unprotected: *, +, ., Vector, diff, int, limit, series
```

Let's define the polarbear coordinate system   as

 > assume( R >= 0, delta >=0, delta <=2*Pi ):