Calculus IV with Maple
Copyright 2002, Dr. Jack Wagner
j.wagner@intelligentsearch.com

Lesson 12: The Integral Theorems of Green, Stokes and Gauss

Topics:   Orientation of volumes, surfaces and closed curves.  The classical theorems of Green, Stokes and Gauss  are presented and demonstrated.
Maple commands introduced: animatecurve

 =

 =  

 =  

The orientability of surfaces has previously been discussed in the context of surface integrals.  Because the integral theorems relate integrals over closed sets to integrals over their boundaries we must now be concerned with the orientation of the boundaries of closed sets relative to the sets.  

Recall that a closed curve, , is positively oriented if increasing values of t result in counterclockwise movement around the curve.  Alternatively, a closed curve, , which is the boundary of a subset, D, of   , is positively oriented if D is on the left as you move around the curve with increasing values of t.

Green's Theorem

Let D be a closed subset of R   with boundary   .  Let P and Q be functions,   R x R -> R , defined and continuous on D with continuous derivatives.  Then,      =  .   Alternatively, if X  is a vector field on D and n  the unit normal on , we have the geometric, coordinate free formulation,      =  

Example 11.1

 Find the area of an ellipse of semi-major axis a and semi-minor axis b.

As in most cases of closed curves, we will use polar coordinates to define the curve  that encloses the elliptic area in the plane.  Notice that this curve is positively oriented with respect to the enclosed area.

>	restart: with(LinearAlgebra):with(VectorCalculus):

Parametrization of the ellipse

>	e := <a*cos(theta), b*sin(theta)>;

Now we define a carefully selected vector field.

>	X := VectorField(<x/2, y/2>, 'cartesian'[x, y]);

>	Divergence(X, [x, y]);

Therefore, from the alternative form of Green's theorem,  =  = A

>	de := diff(e, theta);

>	nds := <de[2], -de[1]>;

Now reparametrize X in polar coordinates.

>	X := evalVF(X, e);

>	dA := X.nds;

>	int(dA, theta=0..2*Pi);

 In the same vein, setting P=-y and Q=x,   

 =  = 2 A

This formula is useful for working with parameterized curves, but working with curves  parameterized by trigonometric functions can be tricky.

Example 11.2

Find the area enclosed by the y axis and the curve, ,    = 0..  

>	restart: with(plots): with(LinearAlgebra): with(VectorCalculus):

>	C := <cos(theta)^2, sin(theta)^3>;

>	plot([C[1], C[2], theta=0..2*Pi]);

>	X := C[1]:Y := C[2]:

>	dX := diff(C[1], theta); dY := diff(C[2], theta);

>	int(X*dY-Y*dX, theta=0..2*Pi);

         So we try the top half.

>	plot([C[1], C[2], theta=0..Pi]);

>	int(X*dY-Y*dX, theta=0..Pi);

Now what is the problem?  The problem is that the curve we are looking at is really two curves superimposed on one another, each oriented oppositely.  The plots from 0 to   and from    to  are identical but in opposite directions.  To see this use the animatecurve function in the plots package.  This function shows the actual drawing of the curve for increasing values of the parameter.

animatecurve

Right click on the plot and select play  from the animate sub-menu.

>	animatecurve([C[1], C[2], theta=0..Pi/2], color=red);

>	animatecurve([C[1], C[2], theta=Pi/2..Pi], color=red);

This makes very explicit the fact that the plot from 0 to   traverses the curve twice; once in each direction.  The integral over this interval is therefore 0.  With this knowledge in hand we try again.

>	int(X*dY-Y*dX, theta=0..Pi/2);

What happens when you integrate from  to  ?   The critical role of orientation is clear!

>	int(X*dY-Y*dX, theta=Pi/2..Pi);

Stokes' Theorem

Stokes theorem    =     , is a generalization of Green's theorem to non-planar surfaces. To see this, consider the projection operator  onto the x-y plane. Let  be the unit tangent vector to , the projection of the boundary of the surface.  Then, =    =  =  =  =

         Let   be the angles between n  and the x, y, and z axes respectively.  Note that   dA  = dxdy   and

n=  . Considering only the projection onto the x-y plane,

 =  =

Stokes theorem is therefore the result of summing the results of Green's theorem over the projections onto each of the coordinate planes.

Example 11.3

Use Stokes' theorem to find the line integral of   around the intersection of the elliptic cylinder   and the plane  .  

In approaching any problem of this sort, a picture is invaluable.

>	restart: with(plots): with(LinearAlgebra):with(VectorCalculus):

>	S := <3*cos(theta), 2*sin(theta), z>;

>	P1 := plot3d(S, theta=0..2*Pi, z=-6..6, color=blue, style=wireframe, grid=[40,40]):

>	P := <x, y, 1-x/2-y/2>;

>	P2 := plot3d(P, x=-4..4, y=-4..4, color=pink, style=patchnogrid):

>	F := VectorField(<(x+y)^3, -(x-y)^3, z^3>, 'cartesian'[x, y, z]);

>	P3 := fieldplot3d(F, x=-4..4, y=-4..4, z=-6..6, color=black, arrows=SLIM, grid=[5, 5, 5]):

>	display(P1, P2, P3, axes=framed, labels=[x, y, z]);

To apply Stokes' theorem, we will integrate curl F  over the elliptic area cut out by the cylinder on the plane.

>	curlF := Del &x F;

Find the unit normal to the plane.

>	dP[x] := diff(P, x);

>	dP[y] := diff(P, y);

>	N := dP[x] &x dP[y];

Normalizing...note n is consistent with the orientation of  the x-y plane

>	n := Normalize(N, 2, inplace);

The dot product (curl F) . n

>	curlFn := evalVF(curlF, <x, y, z>).n;

Find dA for the cross section of the cylinder

>	C := <3*r*cos(theta), 2*r*sin(theta)>;

>	J := Jacobian(C, [r, theta]);

How does this compare to the Jacobian for the circular transformation ?

>	dA := simplify(Determinant(J));

Reparameterize curlFn .

>	curlFn := subs({x=S[1], y=S[2], z=S[3]}, curlFn);

Put it all together!

>	Int(Int(curlFn*dA, r=0..1), theta=0..2*Pi);

>

value(%);

Example 11.4

Evaluate the line integral of     around the curve formed by the intersection of the elliptic paraboloid,  , and the plane  .
It is not always necessary to reparameterize in terms of alternative parameterizations.  Sometimes Cartesian coordinates do just as nicely.

>	restart: with(plots): with(LinearAlgebra):with(VectorCalculus):

>	S := <x, y, 2*x^2+3*y^2>;

>	P := <x, y, 3*y+5*x+2>;

Now we need a picture of the intersection and a computation of the limits of integration.

>	P1 := plot3d(S, x=-2..4, y=-Pi..Pi, color=blue, style=wireframe):

>	P2 := plot3d(P, x=-2..4, y=-Pi..Pi, color=pink, style=patchnogrid):

>	eq := S[3]=P[3];

Intersection of Plane and Surface

>	e1 := solve(eq, y);

>	eq1 := 33-24*x^2+60*x=0;

x limits of integration

>	X := fsolve(eq1, x);

>	e2 := solve(eq, x);

>	eq2 := 41-24*y^2+24*y;

y limits of integration

>	Y := fsolve(eq2, y)

The equation of the curve of intersection comes in two parts.

>	E1 := subs(y=e1[1], P);

>	E2 := subs(y=e1[2], P);

>	P3 := spacecurve({evalm(E1), evalm(E2)}, x=-2..4, color=black, thickness=3, numpoints=1000):

>	display({P1, P2, P3}, labels=[x, y, z], axes=framed, view=[-2..4, -2..4, 0..25]);

Now we compute the integral.

>	F := VectorField(<y^2, x^2, z^2>, 'cartesian'[x, y, z]);

>	curlF := Del &x F;

>	dP[x] := diff(P, x);

>	dP[y] := diff(P, y);

>	N := dP[x] &x dP[y];

>	n := Normalize(N, 2, inplace);

This is our unit normal to the surface enclosed by the curve.

>	J := Jacobian(S, [x, y, z]);

>	J := DeleteColumn(J, [3]);

>	dA := simplify(sqrt(Determinant(Transpose(J).J)));

>	curlFn := evalVF(curlF, <x, y, z>).n;

>	Int(Int(curlFn*dA, x=X[1]..X[2]), y=Y[1]..Y[2]);

>

evalf(%);

 If you try this with some trigonometric parameterization you will find that it is not only more complicated, but Maple will not produce a numerical answer even after a much longer processing time.

Divergence Theorem

Green's theorem also generalizes to volumes. Let V be a closed subset of   with a boundary consisting of surfaces oriented by outward pointing normals.  Then,      The idea is to slice the volume into thin slices.  On each slice, Green's theorem holds in the form,  .  By summing over the slices and taking limits we obtain the divergence theorem.

Example 11.5

Find the flux, through the surface of a sphere of radius R , of the vector field F=[x, y, z] .  

>	restart: with(plots): with(LinearAlgebra): with(VectorCalculus):

The sphere is parameterized with geographic coordinates.

>	V := <r*cos(phi)*cos(theta), r*cos(phi)*sin(theta), r*sin(phi)>;

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

Here is a picture of what we are talking about.

>	S := subs(r=2, V);

>	P1 := fieldplot3d(F, x=-2..2, y=0..2, z=0..2, color=red, grid=[4, 4, 4], arrows=SLIM):

>	P2 := plot3d(S, theta=0..Pi, phi=0..Pi/2, color=blue, style=wireframe):

>	display(P1, P2, scaling=constrained, axes=framed);

>	J := Jacobian(V, [r, phi, theta]);

>	dV := simplify(sqrt(Determinant(Transpose(J).J))) assuming real;

>	divF := Divergence(F, [x, y, z]);

Integrate over the positive quadrant.

>	Int(Int(Int(divF*dV, r=0..R), phi=0..Pi/2), theta=0..Pi);

>	simplify(4*value(%), symbolic);

Practice

1. Given the ellipsoid,    , find the flux through the surface of the vector fields defined by:

     a)    

     b)      

2. Integrate   counterclockwise around the square formed by the lines,  y=1, y=3, x=1, x=3.

3. Consider the curve formed by the intersection of the surface   and the cylinder of radius 1 centered at the point (2, 2). Find the integral of     around this curve.

5. Integrate   around the curve formed by the intersection of
the cylinder of radius 2 whose long axis is the x axis, and the surface,  ,  .