Untitled -- Sage

Calculate the surface integral of flux type $\iint_S^+ P(x,y,z)\, dydz+Q(x,y,z) \,dxdz+R(x,y,z) \,dxdy $ where $S^+$ is a surface with given orientation.

This is equivalent to the problem: Find the flux of the vector field $\vec F(x,y,z)=\,\lt P,Q,R \gt\,$ over $S^+$.

We will solve the (relatively) simple problem of $\iint_S x dydz+y \,dxdz+z \,dxdy $ where S is top side of the triangle x+y+z=1, $x, \,y, \,z \ge 0$. Other examples follow.

Formula Sheet From the formula sheet, we see: $\iint_S^+ P dydz+Q \,dxdz+R \,dxdy = \iint_S \vec F \cdot d \vec S$, where $\vec F= \lt P,Q,R \gt$

Surface element is the vector product: $ d \vec S= \vec n \cdot dS = \pm \,\, d \vec S_u \times d \vec S_v \,\,\, du \, dv$

Remember that the $\pm $ depends on whether or not the parametrization is orientation preserving. We will find this by finding the sign of the dot product of the vectors: $\vec n$ and $d \vec S_u \times d \vec S_v$.

Integrand is the mixed product: $\vec F \,d \vec S = \pm \,\, \vec F \cdot d \vec S_u \times d \vec S_v \,\,\, du\, dv$

Preparation for SOLVER: We must parameterize the surface S with 2 parameters and find any vector $\vec n$ normal to S and in the positive direction of the given orientation.

Work particular to the given problem:

== We have an explicit z=f(x,y). Namely $z=1-x-y$ So we just let x=u, y=v and z=f(u,v). Other parametrizations will not be so easy.

== So: $\left\{ \begin{array}{l}x = u \\y = v\\z =1-u-v \end{array} \right.\,\,\,\,\,u \in [0,1 ] \,\,\,\,v \in [0,-x+1]$

== By examination, we see: $\vec n= <1,1,1>$ (That is, the vector pointing up from the triangle.)

SOLVER: Our program follows the hand solution method. Everything in red requires something particular to your problem!

Parameterize S and input as vector function. Requires parametrization - done above.
Find $\vec n$ and input as vector function. Requires normal - done above.
Input F.
Find the partials of S.
Determine $\pm$ by checking orientation preserving (+1 if orientation-preserving, else -1). Requires checking that you get $\pm 1$.
Substitute parameterization in F.
Find the mixed product $\vec F \,d \vec S = \pm \,\,\,\vec F \cdot d \vec S_u \times d \vec S_v \,\,\, du\, dv$
Find intervals of integration and integrate. Requires intervals of parametrization - done above.

Step 0: The program defines our variables. We are given variables (x,y,z). We need (u,v) as parameters.

var ('u v'); var('x y z')

\newcommand{\Bold}[1]{\mathbf{#1}}\left(u, v\right)

\newcommand{\Bold}[1]{\mathbf{#1}}\left(x, y, z\right)

\newcommand{\Bold}[1]{\mathbf{#1}}\left(u, v\right)

\newcommand{\Bold}[1]{\mathbf{#1}}\left(x, y, z\right)

Step 1: We define $\vec S$, $\vec n$ and $\vec F$. Changes with problem; you must input your parametrization in Sf and your vector in n.

Sf=vector((u,v,1-u-v)) n=vector((1,1,1)) F=vector((x,y,z))

Step 2: The program finds the partial derivatives.

Sprime_u=diff(Sf,u) view(Sprime_u)

\newcommand{\Bold}[1]{\mathbf{#1}}\left(1,\,0,\,-1\right)

Sprime_v=diff(Sf,v) view(Sprime_v)

\newcommand{\Bold}[1]{\mathbf{#1}}\left(0,\,1,\,-1\right)

Step 3: The program checks whether our parameterization is "orientation-preserving" or not. Check that we get 1 or -1.

If we get 0 here, we need to change the point, e.g. (u=2.0,v=2.0). Use "real" values with decimal points. (If we get something other than +1, -1 or 0, we have made an error in step 1.)

npar=Sprime_u.cross_product(Sprime_v) np=sign((n.dot_product(npar))(u=1.0,v=1.0)) view(np)

\newcommand{\Bold}[1]{\mathbf{#1}}1

Step 4: The program defines a standard sage function for "change of variables" for a 2-variable parametrization (surface).

The program changes the variables of $\vec F$.

def changevar(f, eqn, newvar1,newvar2): return f.substitute(eqn)

F=changevar(F,x==Sf[0],u,v) F=changevar(F,y==Sf[1],u,v) F=changevar(F,z==Sf[2],u,v) view(F)

\newcommand{\Bold}[1]{\mathbf{#1}}\left(u,\,v,\,-u - v + 1\right)

Step 5: Then program finds the mixed product and multiplies it by the orientation from step 3.

M=matrix(([F[0],F[1],F[2]],[Sprime_u[0],Sprime_u[1],Sprime_u[2]],[Sprime_v[0],Sprime_v[1],Sprime_v[2]])) Int=np*M.determinant() view(Int)

\newcommand{\Bold}[1]{\mathbf{#1}}1

Step 6: The program computes the integral (flux). Changes with problem - we must put in our intervals of integration.

integral(integral(Int,(v,0,-u+1)),(u,0,1))

\newcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2}

So our flux is: $\iint_S x dydz+y \,dxdz+z \,dxdy = \frac{1}{2}=0.5$

Notice that this last "makes sense" in our particular problem since the integrand is 1 and the area of the projected triangle is 1/2.

Scroll down to the bottom graph 2. It is a visualization of calculating the flux!

Graph 1: We must usually use the parameterization of S that we get to graph the surface. We can do this easily if the intervals for u and v are numbers.

Here: $v \in [0,-u+1] $ and $u \in [0,1] $ so not numbers. I don't know if one can graph with intervals with variables, so I just "drew" the triangle a couple of different ways.

S=polygon([(1,0,0),(0,1,0),(0,0,1)]) show(S,figsize=3)

T=sum([line3d([(0,0,1),(c,1-c,0)], color=hue((c+8)/8), width=0.8) for c in [0..1,step=0.01]]) show(T,figsize=3)

Graph 2: Now we want a visual interpretation of what we have done. Remember that with the change of variables from the parametrization, we have:

$\vec F= \lt u,v,1-u-v \gt$
$\vec n=\lt 1,1,1 \gt$
Integrand = 1

So the vectors that describe the normalized flux volume (and scaled by vector-scale=$vs$) are:

Starting point: $\lt u,v,1-u-v \gt$
End point: ${\lt u,v,1-u-v \gt }$ + $\lt \color{teal}{1} \cdot \color{purple}{1} \cdot vs,\color{teal}{1} \cdot \color{purple}{1}\cdot vs,\color{teal}{1} \cdot \color{purple}{1}\cdot vs \gt$= ${\lt u,v,1-u-v \gt }$ =${\lt u+1\cdot vs,v+1\cdot vs,1-u-v+1\cdot vs \gt }$

Remember that our surface with respect to u and v has $u \in [0,1]$ and $v \in [0,-u+1].

To do our programming interation, we substitute "c" for u and "d" for v and run d from [0,-c+1] and c from [0,1]. (Usually start tries with big steps and see how it goes and then reduce step size.)

The FLUX is the VOLUME between the surface S and the surface 'formed' by the endpoints of the arrows (multiplied by the orientation and divided by $\left\| {\vec n} \right\|$ )!

Here I have added a "blank" triangle so that the box is nice. Probably and easier way to do this...

vs=.5 vf=sum([sum([arrow3d((c,d,1-c-d),(c+1*vs,d+1*vs,1-c-d+1*vs), color=hue((c+8)/8), width=0.8) for d in [0..(-c+1),step=0.2]]) for c in [0..1,step=0.2]]) T_blank=polygon([(2,0,0),(0,2,0),(0,0,2)], color='white', opacity=0) show(T+vf+T_blank,figsize=3)

Untitled

446 days ago by pub