surfintvec.mws
Surface Integrals of Vector Fields
Using Maple and the
vec_calc
Package
This worksheet shows how to compute surface integrals of vector fields using
Maple
and the
vec_calc
package. As examples we compute
* The Flux of a Magnetic Field
* The Expansion of a Gas
To start the
vec_calc
package, execute the following commands:
>
restart;
>
libname:="C:/mylib/vec_calc7", libname:
>
with(vec_calc): vc_aliases:
>
with(linalg):with(student):with(plots):
Warning, the protected names norm and trace have been redefined and unprotected
Warning, the name changecoords has been redefined
>
The Flux of a Magnetic Field
The magnetic field of a current
moving along an infinitely long straight wire has magnitude
where
is a constant and
is the distance from the wire. Its direction is counterclockwise around the wire as given by the right hand rule. Suppose two wires are parallel to the
-axis and located at
and
and each carries a current
but in opposite directions. We want to find the flux of the magnetic field through the rectangle
,
, and
.
The magnetic field of the wire at
is
>
B1:=MF( [x,y,z], evall(mu[0]*i/(4*Pi*((x-2*a)^2+(y)^2))*[-y,x-2*a,0]) );
The magnetic field of the wire at
is
>
B2:=MF( [x,y,z], evall(-mu[0]*i/(4*Pi*((x+2*a)^2+(y)^2))*[-y,x+2*a,0]) );
The total field is
>
B3:=MF( [x,y,z], evall(B1(x,y,z)+B2(x,y,z)) );
For the purpose of plotting the functions (and making numerical computations) we make the following assignments to the variables:
>
a:=1; b:=1; c:=1; i:=1; mu[0]:=1;
We first plot the two functions separately:
>
B1plot:=fieldplot3d(B1,-3..3,-1..1,-1..1, grid=[9,9,2], color=blue, arrows=SLIM, orientation=[-90,0], axes=normal, tickmarks=[7,3,3]):
>
B2plot:=fieldplot3d(B2,-3..3,-1..1,-1..1, grid=[9,9,2], color=red, arrows=SLIM, orientation=[-90,0], axes=normal, tickmarks=[7,3,3]):
>
display(array([B1plot,B2plot]));
Notice that for both fields, the arrows point in the negative
-direction between the two wires. So the fields add up. Now we plot the total field:
>
B3plot:=fieldplot3d(B3,-3..3,-1..1,-1..1, grid=[9,9,2], color=black, arrows=SLIM, orientation=[-90,0], axes=normal, tickmarks=[7,3,3]): B3plot;
Before we go on, we unassign the parameters:
>
a:='a'; b:='b'; c:='c'; i:='i'; mu[0]:='mu[0]';
The rectangle may be parametrized by
>
Rect:=[u,0,v];
where
, and
. The tangent vectors and the normal vector are
>
Ru:=diff(Rect,u); Rv:=diff(Rect,v); N:=Ru &x Rv;
Notice that the normal points in the negative
-direction, which is the direction we saw in the plot that the magnetic field points.
On the rectangle, the magnetic field is
>
B3(u,0,v); simplify(%);
So the flux is
>
interface(showassumed=0); assume(0<b,b<2*a,0<c);
>
Muint(B3(u,0,v) &. N,u=-b..b,v=0..c); Flux:=value(%);
>
Flux:=combine(%,ln);
Note: There is no guarantee that if you re-execute these commands that Maple will give the answer in the same form.
Another way to compute this surface integral is to use the
Surface_int_vector
command (or its alias
Siv
) from the
vec_calc
package which works directly with the parametrized surface and the vector field:
>
Ruv:=MF([u,v],Rect);
>
Siv(B3,Ruv,u=-b..b,v=0..c); Flux:=value(%);
There is a second way to compute this flux. The magnetic field is solenoidal because it is divergence-free:
>
DIV(B3);
So the magnetic field has a vector potential,
.
>
VEC_POT(B3,'A');
>
A(x,y,z);
By Stokes' Theorem, the flux is the line integral of the vector potential around the boundary curve. Here the boundary consists of four line segments:
>
r1:=MF([u],[u,0,0]); #for u=-b..b
>
r2:=MF([v],[b,0,v]); #for v=0..c
>
r3:=MF([u],[-u,0,c]); #for u=-b..b
>
r4:=MF([v],[-b,0,c-v]); #for v=0..c
So the flux is the total line integral:
>
Liv(A,r1,u=-b..b)+Liv(A,r2,v=0..c)+Liv(A,r3,u=-b..b)+Liv(A,r4,v=0..c); Flux:=value(%);
>
The Expansion of a Gas
The velocity field of a gas is
. We want to find the expansion of the fluid out through the sphere
which is defined as
We enter the fluid velocity as a function:
>
V:=MF([x,y,z],[x*y^2,y*z^2,z*x^2]);
The sphere may be parametrized as
>
Rsph:=[2*sin(phi)*cos(theta),2*sin(phi)*sin(theta),2*cos(phi)];
So the tangent vectors are:
>
R[theta]:=diff(Rsph,theta);
>
R[phi]:=diff(Rsph,phi);
and the normal vector is
>
N:= R[theta] &x R[phi];
We simplify the third component:
>
N[3]:=simplify(N[3]);
The components of
are negative in the first octant. So
points inward. We want the flux outward. So we reverse
:
>
N:=-N;
On the sphere, the velocity is
>
VR:=V(op(Rsph));
Note: We need the op command to strip the square brackets off
.
The dot product of the velocity field and the normal is:
>
VN:=VR &. N;
So the expansion is:
>
Muint(VN,theta=0..2*Pi,phi=0..Pi); Expansion:=value(%);
Another way to compute this surface integral is to use the
Surface_int_vector
command (or its alias
Siv
) from the
vec_calc
package which works directly with the parametrized surface and the vector field:
>
R:=MF([theta,phi],Rsph);
>
Siv(V,R,theta=0..2*Pi,phi=0..Pi); Expansion:=value(%);
Notice there is a minus sign difference in the answer because we did not reverse the normal.
There is a second way to compute this expansion. By Gauss' theorem, the expansion is the volume integral of the divergence of the velocity over the interior of the sphere. In this case, the divergence is:
>
DIV(V);
which in spherical coordinates is
>
divV:=rho^2;
The spherical volume Jacobian is
>
J:=rho^2*sin(phi);
So the expansion is
>
Muint(divV*J,rho=0..2,theta=0..2*Pi,phi=0..Pi); Expansion:=value(%);
>
>