sci.math.num-analysis SUNA15, Non-uniform 2-D convection ================================== In the sequel, references [*] correspond with SUNA* articles. Skew Upwind Finite Element Difference schemes were developed for a _uniform_ velocity field [05,06,07]. This obviously imposes severe limitations upon the applicability of the theory. The present article overcomes this limitations. We will restrict ourselves in what follows, however, to incompressible fluids. According to [05], equation (5), the general SUFED scheme for 2-D pure steady uniform convection is (apart from a different offset for numbering): [ + (y3 - y1).u - (x3 - x1).v ].(T1 - T2) + [ - (y2 - y1).u + (x2 - x1).v ].(T1 - T3) = 0 According to [14] Finite Differences schemes for the Stream Function F are: F2 - F1 = (y2 - y1).u - (x2 - x1).v 3 F3 - F1 = (y3 - y1).u - (x3 - x1).v o / \ Consequently: / \ 1 o -------- o 2 (F3 - F1).(T1 - T2) - (F2 - F1).(T1 - T3) = 0 This expression can be written in several other ways: = F1.(T2 - T3) + F2.(T3 - T1) + F3.(T1 - T2) = T1.(F3 - F2) + T2.(F1 - F3) + T3.(F2 - F1) Note that the temperature at a node is associated with streamfunction values that are _opposite_ to this node. Which seems a sensible result, because the streamfunction has something to do with flow that comes in at a triangle side, and such a flow can only influence the opposite node. The above scheme was found independently, and earlier, in a more specialized form, by my ex-collegue Jan Tulp. We worked at Neratoom at that time, and he called it a "dupwind" scheme. (I don't know of any references. Neratoom does not exist anymore.) The form [ (F3 - F1).(T1 - T2) - (F2 - F1).(T1 - T3) ] has the following analytical equivalent: - dF/dk.dT/dh + dF/dh.dT/dk = 0 --> (dT/dh)/(dT/dk) = (dF/dh)/(dF/dk) With other words: the gradients of the isotherms are the same as the gradients of the streamlines, or the isotherms follow the streamlines. Which is of course true for pure steady convection. Another analytical analogue follows with the continuity equation (conservation of mass) and Green's theorem: u.dT/dx + v.dT/dy = d(u.T)/dx + d(v.T)/dy . Integrating: // / / || [d(T.u)/dx + d(T.v)/dy] dx.dy = O T.u.dy - T.v.dx = O T.dF // / / which resembles the form: T1.(F3 - F2) + T2.(F1 - F3) + T3.(F2 - F1) Concerning the other form: F1.(T2 - T3) + F2.(T3 - T1) + F3.(T1 - T2) reverse engineering is applied: / / // O F.dT = O F.[dT/dx.dx + dT/dy.dy] = || [d/dx(F.dT/dy) - d/dy(F.dT/dx)] dx.dy / / // // // = || [dF/dx.dT/dy - dF/dy.dT/dx] dx.dy = - || [ u.dT/dx + v.dT/dy ] dx.dy // // / / // Conclusion: O T.dF = - O F.dT = || [ u.dT/dx + v.dT/dy ] dx.dy , analytically. / / // The scheme will be applied now to a patch of triangles around a mesh point (0), according to the figure below: 5 o - o 4 Suppose for a moment there is NO upwinding. Then: / \ / \ 6 o - 0 - o 3 T0.(F2 - F1) + T1.(F0 - F2) + T2.(F1 - F0) + \ / \ / T0.(F3 - F2) + T2.(F0 - F3) + T3.(F2 - F0) + o - o T0.(F4 - F3) + T3.(F0 - F4) + T4.(F3 - F0) + ... + 1 2 T0.(F1 - F6) + T6.(F0 - F1) + T1.(F6 - F0) = 0 Simplified: T0.( 0 ) + T1.(F6 - F2) + T2.(F1 - F3) + T3.(F2 - F4) + .... = 0 | NO upwinding always results in a zero coefficient for the temperature ( T0 ) | to be calculated. Which is commonly recognized as weird. Rules for upwinding the non-uniform 2-D scheme will be derived in the sequel. For convenience, we first introduce FLUXES instead of stream-function values. The relationship between the former and the latter is established by means of: 3 dF1 = F3 - F2 o dF2 = F1 - F3 ; dF1 + dF2 + dF3 = 0 / \ dF3 = F2 - F1 or - dF2 - dF3 = dF1 2 * * 1 / \ ! Take a good look at o ---- * ---- o this way of numbering. 1 3 2 At this point the velocity field becomes non-uniform ! According to [14], the streamfunction is discretized along a meshline (1,2) as: F2 - F1 = (y2 - y1).u3 - (x2 - x1).v3 & same for (3,2) and (1,3) . ^ ^ It is emphasized that velocities are located at staggered, mid-side positions. Furthermore, it is advantageous to establish any numbering conforming the idea that temperatures Tk are influenced by opposite (uk,vk), according to [05], or the idea that temperatures and fluxes are kind of "conjugate" quantities. Upwinding is introduced next. It is remembered from above that: - dF2.(T1 - T2) - dF3.(T1 - T3) = 0 - ______________ | | | At the other hand, equation (3) in [05] reads: ==> | dF2 >= 0 | | | dF3 >= 0 | (U <= 0 ).(T1 - T2) + (V <= 0 ).(T1 - T3) = 0 - |______________| Upwinding Thus: T1.dF1 + T2.dF2 + T3.dF3 = 0 ^^ <0 >=0 >=0 : rule of coefficients [Patankar] OK Theorem: there is only ONE or NO upwind point in a triangle. Lemma: there can not be two or three upwind points. Proof: suppose there are two upwind points, say (1) and (2). Then: dF2 >= 0 & dF3 >= 0 for (1) --> F1 >= F3 & F2 >= F1 dF1 >= 0 & dF3 >= 0 for (2) --> F3 >= F2 & F2 >= F1 Conclusion: F1 = F2 = F3 <-- F1 >= F3 >= F2 & F1 <= F2 This means that 0.T1 + 0.T2 + 0.T3 = 0 so temperatures can be anything. The above is contradictory to the assumption that one of the <= must be < . Beware of the special cases! There are essentially _two_ of them: First case: F1 = F2 , F2 = F3 therefore F3 = F1 . 2______ 1 For example a CORNER in a heat exchanger. Then F1 = F2 , F2 = F3 , | <-- / therefore 0.T1 + 0.T2 + 0.T3 = 0 --> temperatures _undefined_. | / Because of other triangles in the mesh, this in fact only affects | / the corner temperature: T2 = undefined , while T3 = T1 . A more 3 severe situation occurs when a vortex is present. Then equations like T1 = T2 = T3 = T4 = T5 = ... = T1 occur, which render many temperatures _independent_ of any boundary conditions! Second case: dF1 = 0 , dF2 > 0 , dF3 < 0 (: apart from 1,2,3 permutations) For example if a mesh-line coincides with a streamline: F2 ------- F3 . Then dF1.T1 + dF2.T2 + dF3.T3 = 0 & dF2 + dF3 = 0 --> dF3.(T3 - T2) = 0 . Consequently: T2 = T3 , as it should be. This case especially occurs with flow BOUNDARIES. It may be considered as good practice to set up separate equations for the streamline(s) along a fixed boundary, because numerical inaccuracy with an upwinding procedure for the bulk-flow might too often "miss" points there. - * 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 * E-mail: Han.deBruijn@RC.TUDelft.NL --| Fax: +31 15 78 37 87 ----* Blood