ProblemCompute the flux of the vector field |
> F := x^2*i+x*y*j+x^3*y^3*k;
2 3 3
F := x i + x y j + x y k
|
over the surface of the tetrahedron bounded by the plane |
|
|
Solution:
| Recall the expression of the Divergence Theorem: |
> ;
/ / / / /
| | | | |
| | F.N dS = | | | div F dV
| | | | |
/ / / / /
S V
| Let's start with the right hand side tripple integral. First we compute div F. |
> divF := diverge(evalm(F),[x,y,z]);
divF := 3 x
| and we then integrate over the solid tetrahedron, |
> I2 := Int(Int(Int(divF,z=0..1-x-y),y=0..1-x),x=0..1);
1 1 - x 1 - x - y
/ / /
| | |
I2 := | | | 3 x dz dy dx
| | |
/ / /
0 0 0
> I2 := int(int(int(divF,z=0..1-x-y),y=0..1-x),x=0..1);
I2 := 1/8
| Now the left hand side. Let us split the surface into the four triangles: S = Sxy + Sxz + Syz + Sxyz. Where Sab is the triangle on the ab plane and Sxyz is the triangle on the plane x+y+z=1. Let us denote the surface integrals over Sxy, etc... by J1,J2,J3 and J4. We have, that the outward unit normal for Sxy is (-k) and thus, |
> ;
/ /
| | 3 3
J1 := | | -x y dS
| |
/ /
Sxy
| which evaluates to: |
> J1 := int(int(dotprod(F,-k),y=0..1-x),x=0..1);
-1
J1 := ----
1120
| Now compute the two central surface integrals J2 and J3. The outward unit normals for Sxz and Syz are j and i respectively so, |
> ;
/ /
| |
J2 := | | x y dS = 0
| |
/ /
Sxz
| this integral is 0 since y=0 on the xz plane. For J3 we have, |
> J3 := Int(Int(dotprod(F,i),S),u);
/ /
| | 2
J3 := | | x dS = 0
| |
/ /
Syz
| again this integral is 0 since x=0 on the yz plane. For the last integral we see that the ouward unit normal is |
> N := evalm(vector([1,1,1])/sqrt(3));
[ 1/2 1/2 1/2]
N := [1/3 3 , 1/3 3 , 1/3 3 ]
| Notice that we need to divide by sqrt(3) (the norm of (1,1,1)) to get a UNIT normal. Hence, |
> ;
/ /
| | 2 1/2 1/2 3 3 1/2
J4 := | | 1/3 x 3 + 1/3 x y 3 + 1/3 x y 3 dS
| |
/ /
Sxyz
|
Computing the dS for this last integral we find
|
> J4 := Int(Int(dotprod(F,[1,1,1]),y=0..1-x),x=0..1);
1 1 - x
/ /
| | 2 3 3
J4 := | | x + x y + x y dy dx
| |
/ /
0 0
| which evaluates to: |
> J4 := int(int(dotprod(F,[1,1,1]),y=0..1-x),x=0..1);
141
J4 := ----
1120
| Finally combining the J's we get: |
> I1 := J1+0+0+J4;
I1 := 1/8
| The same as before!! |