sci.math.research

SUNA11, Theorem Preliminaries
=============================

This article is part of the "Series on Unified Numerical Approximations"


             UNIFICATION OF FINITE VOLUME AND GALERKIN METHODS

         FOR CONSERVATION EQUATIONS, IN THE TWO-DIMENSIONAL CASE

                          PART I : PRELIMINARIES


Consider the following partial differential equation:

         dQx/dx + dQy/dy = 0           ("d" denotes partial differentiation)

Here (Qx,Qy) can be interpreted as the heat-flow in a two-dimensional medium
or as the velocities (u,v) in a 2-D flow field.

According to the Finite Domain technique, a kind of Finite Difference Method,
this equation is first integrated over a certain control volume/area:

         //
        || (dQx/dx + dQy/dy) dx.dy = 0
       //

Partial integration (Green's theorem) results in an equivalent formulation,
which contains a line-integral:

       /
      O (Qx.dy - Qy.dx)           (counterclockwise)
     /

According to the (Galerkin) Finite Element Method, the conservation equation
must first be multiplied by a weighting function  F , and then integrated over
a finite element area, according to:

       //
      || F.(dQx/dx + dQy/dy) dx.dy = 0
     //

Partial integration results in an expression with line-integrals, all of which
become part of the boundary conditions, and a term for the bulk material:

      //
   - || (dF/dx.Qx + dF/dy.Qy) dx.dy           (: mind the minus sign)
    //

For the sake of generality, both the Finite Volume and the Galerkin formulation
are assumed to be applicable at a curvilinear grid, consisting of quadrilateral
shaped cells or elements.
For the sake of simplicity, we wish to restrict ourselves to regular molecules
and elements, such as the classical five point star, and square quadrilaterals.

The general curvilinear grid is mapped upon a regular grid by a transformation
of coordinates  x(h,k) , y(h,k) , in such a way that the _whole_ problem may be
assumed to be regular in the  (h,k)  plane. This is the Finite Volume point of
view, essentially a _global_ transformation.
From a Finite Element point of view, however, such coordinate transformations
are always carried out _locally_/element-wise: by isoparametrics. Restrictions
of orthogonality do not apply in a Finite Element context.

In either case, integrands of the Finite Volume and the Galerkin formulations
must be transformed first into such regular  (h,k)  coordinates.

The vector  (Qx,Qy)  is expressed in components at the (h,k) coordinate system
as follows:

  | Qx |      | dx/dh |      | dx/dk |
  |    | = Qh.|       | + Qk.|       |
  | Qy |      | dy/dh |      | dy/dk |

The inverse of this transformation is the one that's needed:

    J.Qh = + dy/dk.Qx - dx/dk.Qy
    J.Qk = - dy/dh.Qx + dx/dh.Qy

Where:  J = dx/dh.dy/dk - dx/dk.dy/dh   (assume that  J > 0 ).

Differentials of the global coordinates are transformed as follows:

     dx = dx/dh.dh + dx/dk.dk           (with proper interpretation of the
     dy = dy/dh.dh + dy/dk.dk            partial and non-partial "d's")

Last but not least:

       dT/dh = dT/dx.dx/dh + dT/dy.dy/dh
       dT/dk = dT/dx.dx/dk + dT/dy.dy/dk

Solving for  dT/d(x,y)  by Cramer's rule:

     J.dT/dx = dy/dk.dT/dh - dy/dh.dT/dk
     J.dT/dy = dx/dh.dT/dk - dx/dk.dT/dh

The integrand of the Finite Volume formulation can now be written as follows:

      Qx.dy - Qy.dx = Qx.[ dy/dh.dh + dy/dk.dk ] - Qy.[ dx/dh.dh + dx/dk.dk ]
                    = [ Qx.dy/dh - Qy.dx/dh ].dh + [ Qx.dy/dk - Qy.dx/dk ].dk
                    = J.[ Qh.dk - Qk.dh ]

Which gives a remarkebly simple result:

||     /                    /
||    O (Qx.dy - Qy.dx) =  O J.(Qh.dk - Qk.dh)
||   /                    /

Same kind of procedure for the Galerkin integral:

  { dF/dx.Qx + dF/dy.Qy }.dx.dy = { J.dF/dx.Qx + J.dF/dy.Qy }.dh.dk =
{ [ dy/dk.dF/dh - dy/dh.dF/dk ].Qx + [ dx/dh.dF/dk - dx/dk.dF/dh ].Qy }.dh.dk =
{ dF/dh.[ dy/dk.Qx - dx/dk.Qy ] + dF/dk.[ dx/dh.Qy - dy/dh.Qx ] }.dh.dk =
{ dF/dh.J.Qh + dF/dk.J.Qk }.dh.dk

Same kind of remarkebly simple result:

||      //                                 //
||   - || (dF/dx.Qx + dF/dy.Qy) dx.dy = - || J.(dF/dh.Qh + dF/dk.Qk) dh.dk
||    //                                 //

In order to shorten notation, we shall make everywhere the substitution:

    Qh := J.Qh
    Qk := J.Qk

Resulting in the two localized integrals:

    /                             //
   O (Qh.dk - Qk.dh)     ;     - || (dF/dh.Qh + dF/dk.Qk) dh.dk
  /                             //
           F.D.                              F.E.

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
* E-mail: Han.deBruijn@RC.TUDelft.NL --| Fax: +31 15 78 37 87 ----* Blood