sci.math.research
SUNA12, Theorem Proving
=======================
This article is part of the "Series on Unified numerical Approximations"
UNIFICATION OF FINITE VOLUME AND GALERKIN METHODS
FOR CONSERVATION EQUATIONS, IN THE TWODIMENSIONAL CASE
PART II : THE PROOF
We are in the process of discretizing P.D.E.'s of the following form:
dQx/dx + dQy/dy = 0 ("d" denotes partial differentiation)
In part I, two formulations were presented, which serve as a starting point for
discretizing these equations according to the two mainstream numerical methods:
Finite Volume and Galerkin. We arrived at equivalent _regular_ formulations of
the corresponding integrals, using coordinate transformations. The integrals
are:
/ //
 [ Qh.dk  Qk.dh ] ;   [ dF/dh.Qh + dF/dk.Qk ] dh.dk
/ //
Assembling local Finite Element matrices results in a global system's matrix,
which can be easily looked upon as a collection of Finite Difference schemes:
each row of the global matrix corresponds to a F.D. method, so to speak.
Using such a F.D. method, the grid is seen to be built up from five node stars,
where the central node is enclosed by a control area. The lineintegral applies
to the boundaries of this control area.
Employing a F.E. viewpoint, on the contrary, the same grid is seen to be built
up from quadrilaterals. The Galerkin integral shall be evaluated for each of
these elements in the mesh.
Switching again to the F.D. viewpoint, all the quadrilateral elements contain
four pieces of a finite domain control area, each piece belonging to a quarter
of a nodal point:
3 4 3_________4
   
________ ________
 <    <  
________1________2 ________
   1 2
________

 F.D. view <=> F.E. view
In this way, Five Node Stars in the Finite Difference/Volume domain correspond
with Quadrilaterals in the Finite Element domain.
As a consequence, there are also four pieces wherein the Finite Domain integral
is splitted up:
/
 ( dk.Qh  dh.Qk )
/
These lineintegrals will successively correspond to the following vector of,
necessarely incomplete, equations (watch out for the signs!):
dk . Qh  dh . Qk
piece at node 1 : +1/2.1/2.(Qh1 + Qh2) + 1/2.1/2.(Qk1 + Qk3)
at node 2 : 1/2.1/2.(Qh1 + Qh2) + 1/2.1/2.(Qk2 + Qk4)
at node 3 : +1/2.1/2.(Qh3 + Qh4)  1/2.1/2.(Qk1 + Qk3)
at node 4 : 1/2.1/2.(Qh3 + Qh4)  1/2.1/2.(Qk2 + Qk4)
On the other hand, for the same quadrilateral element, there exists a Galerkin
formulation. The shape functions of a quadrilateral element can be chosen as
follows:
N1(h,k) = (1  h).(1  k) ; N2(h,k) = h.(1  k)
N3(h,k) = (1  h).k ; N4(h,k) = h.k
An arbitrary function F at the cluster will be interpolated by this polynomial
base:
F(h,k) = N1(h,k).F1 + N2(h,k).F2 + N3(h,k).F3 + N4(h,k).F4
The partial derivatives of F are calculated: o  o
 * * 
dF/dh = (1  k).(F2  F1) + k.(F4  F3)  * * 
dF/dk = (1  h).(F3  F1) + h.(F4  F2) o  o
Numerical integration will be employed. In order to accomplish this, four so
called integration points (*) have to be chosen. It is presumed, however, that
these integration points _don't_ have to be Gaussian per se. See figure above.
This means that there exist four points (*):
(1/2  a,1/2  a) , (1/2 + a,1/2  a) , (1/2 + a,1/2  a) , (1/2 + a,1/2 + a)
where 0 < a < 1/2 , also with weighting factor 1/4 , but NOT necessarily: a =
= 1/(2.sqrt(3)) , which would define numerical integration according to Gauss.
The Galerkin integral
//
  [ dF/dh.Qh + dF/dk.Qk ] dh.dk
//
can now be put in algebraic form:
  (1  k).Qh  (1  h).Qk 
  F1 F2 F3 F4  . 1/4 . Sigma  + (1  k).Qh  h .Qk 
over(*)   k .Qh + (1  h).Qk 
 + k .Qh + h .Qk 
Sigma
Here over(*) denotes summation over the integration points.
By partial differentiation to the nodal values of F , the Galerkin integral
reaches its stationary value. A vector of equations is then associated which
each quadrilateral, finally resulting in what is called the local "stiffness"
matrix and "load" vector of the element considered.
The precise choice a of the integration points was _deliberately_ left open.
We take the freedom to select a in such a way that the integration points are
precisely coincident with the nodal points: a = 1/2 (= 0.5 instead of 0.289).
Thus substitute h = 0  1 and k = 0  1 (with "" = "xor") and sum up:
h = 0 , k = 0 h = 1 , k = 0 h = 0 , k = 1 h = 1 , k = 1
   
  Qh1  Qk1    Qh2  0    0  Qk3    0  0 
 1/4.[  + Qh1  0  +  + Qh2  Qk2  +  + 0  0  +  + 0  Qk4  ]
  0 + Qk1    0 + 0    Qh3 + Qk3    Qh4 + 0 
 + 0 + 0   + 0 + Qk2   + Qh3 + 0   + Qh4 + Qk4 
This results in the following set of, necessarily incomplete, equations:
  (Qh1 + Qh2)  (Qk1 + Qk3) 
 1/4 .  + (Qh1 + Qh2)  (Qk2 + Qk4) 
  (Qh3 + Qh4) + (Qk1 + Qk3) 
 + (Qh3 + Qh4) + (Qk2 + Qk4) 
Comparing this with the vector for the Finite Difference (domain) Method shows
that these expressions are IDENTICAL.
^^^^^^^^^
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