sci.math.numanalysis
SUNA33, Streamline Upwind & FLOTRAN
===================================
The last article I copied from the "International Journal for Numerical Methods
in Engineering" is written by Rita J. Schnipke and James G. Rice. It is called:
"A Finite Element Method for Free and Forced Convection Heat Transfer". I made
the copy because, even at a first glance, there are ideas of Unification in it.
Quoting without permission from vol. 24, page 119 (1987):
> This formulation is an iterative, segregated solution method _along the lines
> of most finite difference methods_. [ .. ] The method is essentially a finite
> element adaption of the SIMPLER algorithm proposed by Patankar.
"Reference 7" was underlined: J.G. Rice and R.J. Schnipke, "A Monotone Stream
line Upwind Finite Element Method for Convectiondominated Flows". After having
changed employers, I got on the right track again, but a couple of years later.
Several phonecalls were made in december 1991, with the vendor of a well known
finite element package: ANSYS. A module for calculating fluids, called FLOTRAN,
can be delivered with ANSYS. The folder about it made me curious. Documentation
about FLOTRAN was sent, ... with the abovementioned article being part of it!
BTW, it can be found in "Computer Methods in Applied Mechanics and Engineering"
vol. 48 (1985) page 313327.
The RiceSchnipke method of streamline upwinding will be compared in this post
with SUFED, the Skew Upwind Finite Element Differences which were invented by
myself. Moreover, It will be proven that the two methods are in fact (more or
less) _equivalent_ with each other. Thus encouraging me to order FLOTRAN ;)
F1  / First, there is a definition of "downwind nodes".
2 o  o  Quoting without permission:
 /  1 > The term 'downwind' defines a node for which
 /  > the negative of the velocity vector points back
F2  /  F4 > into the element. [ Node (1) in the figure ]
 /  > Once a downwind node has been identified in
 / Fn  > this manner, the advection term approximation
3 o / o 4 > requires the determination of the upstream
x',y' F3 > location illustrated in [ the figure ], denoted
> by x' and y'. Once x' and y' are located,
> the advection terms for the node is approximated [ in my own ASCII ] as:
Us/ds.(T1  T')
where: ds = sqrt{ (x1  x')^2 + (y1  y')^2 }
Us = sqrt{ u1^2 + v1^2 }
> The location of the x' and y' coordinates as well as the calculation of T'
Well, it is explained how these quantities have to be calculated. A significant
detail in the argument is that x', y' and T' :
> will be a function of either [ contributions from node 2 ] and [ node 3 ] or
> [ from node 3 ] and [ node 4 ] but never all three.
Huh? Doesn't this in fact mean that NO quadrilaterals, but only _triangles_ are
actually in effect? A little further inconsequence is found on the next page:
> .. the linear arc length calculation ... is used in this formulation.
^^^^^^
Other evidence is found in SUNA32, where it was shown that upwind differencing
in quadrilaterals (a la HdB) leads in a very natural way to skew upwind schemes
in triangles alone. So let's restrict our attention to triangles, and apply the
techniques of FLOTRAN, nevertheless. Reference is also made to my SUNA15 post.
At first, some properties of the "flow rates" F will be recapitulated here.
Note that these F do _not_ correspond to Stream Function values, as is the
case in SUNA15, but to differences of the stream function values instead. The
correspondence SUNA15 > RiceSchnipke is as follows. See ASCII figure below:
(F3  F2) > F1 ; (F1  F3) > F2 ; (F2  F1) > F3 .
Define instead, as has been done already in SUNA15:
dF1 = F3  F2 ; dF2 = F1  F3 ; dF3 = F2  F1
As a consequence: dF1 + dF2 + dF3 = 0 (: sum of Fluxes must be zero).
3 The upstream point in a parent triangle (displayed <)
 \ is determined a la RiceSchnipke as follows:
 \ T'
dF2   X dF1 k/h = dF3/dF2 ; h + k = 1
 /  \
 /  \ Thus: h.(1 + dF3/dF2) = 1 > h = dF2/(dF2 + dF3)
 1  2 k = dF3/(dF2 + dF3)
/  dF3
T  T1 = h.(T2  T1) + k.(T3  T1) >
Using dF1 + dF2 + dF3 = 0 : T' T1 =  dF2/dF1.(T2  T1)  dF3/dF1.(T3  T1)
Due to isoparametrics: x' x1 =  dF2/dF1.(x2  x1)  dF3/dF1.(x3  x1)
y' y1 =  dF2/dF1.(y2  y1)  dF3/dF1.(y3  y1)
Factors Fn and 1  Fn in the RiceSchnipke method can be equivalenced with
our  dF2/dF1 and  dF3/dF1 , which indeed sum up to 1 = ( dF2  dF3)/dF1 .
By definition, the relationship between global (u,v) and local velocities is:
^^^^^^^^^^
u = U.(x2  x1) + V.(x3  x1) =  dF2/J.(x2  x1)  dF3/J.(x3  x1)
v = U.(y2  y1) + V.(y3  y1) =  dF2/J.(y2  y1)  dF3/J.(y3  y1)
Where: J = (x2  x1).(y3  y1)  (y2  y1).(x3  x1) = 2 * triangle area
The upstreampoint coordinates are recognized herein:
u = dF1/J.(x' x1)
v = dF1/J.(y' y1) : velocity vector is directed along the streamline!
Let us calculate now the term Us/ds.(T1  T') :
Us = sqrt{ u^2 + v^2 } = dF1/J.sqrt{ (x' x1)^2 + (y' y1)^2 }
J.Us/ds =  dF1 , because dF1 < 0 .
We found: T' T1 =  1/dF1.{ dF2.(T2  T1) + dF3.(T3  T1) }
Hence: J. Us/ds.(T' T1) = dF2.(T2  T1) + dF3.(T3  T1)
= dF1.T1 + dF2.T2 + dF3.T3
Conclusion: RiceSchnipke = SUNA15 ... Whoa, not so fast!

It remains to be checked that no smuggling with signs has occurred in between.
Worse, the _real_ smuggling is in the silent replacement of (u1,v1) by (u,v).
The only possible way out is to _define_ velocities (u,v) _instead_ of (u1,v1).
With other words: it should be recognized that the assumption (u,v) = (u1,v1)
was indeed the (only) WEAK point in the RiceSchnipke theory. 
It is replaced by the above. The velocity at the streamline gets _defined_ by:
^^^^^^^^^^^^^^^^^^^^^^^^^^
 u   x' x1  But ... don't ask me where (u,v) is located.
  = dF1/J .   Certainly not in the middle of the streamline,
 v   y' y1  as was previously suggested (: SUNA32).
To be continued ...

* Han de Bruijn; Applications&Graphics  "A little bit of Physics * No
* TUD Computing Centre; P.O. Box 354  would be NO idleness in * Oil
* 2600 AJ Delft; The Netherlands.  Mathematics" (HdB). * for
* Email: Han.deBruijn@RC.TUDelft.NL  Fax: +31 15 78 37 87 * Blood