Supported by the State key 973-Project for fundamental Research and
Jiang Chun Bo , Chen Liqiu and Chen Yongcan
Dept. Hydraulic Eng., Tsinghua University, Beijing 100084
E-mail: jcb@mail.tsinghua.edu.cn
Supported by the State key 973-Project for fundamental Research and
Development Program(Grant No. G1999043607 ) and the National Science
Foundation of China, Grand No. 59979013
Abstract:
A fractional step finite volume method is proposed
to solve the weakly compressible flow equations. The stability analysis in
multi-dimensional case is carried out, which shows that the present method is
more accurate and stable compared with the conventional finite volume
formulations. It has third-order accuracy and uniform CFL condition.
Keywords: finite volume scheme, weakly compressible model, stability analysis.
The finite volume method (FVM) is becoming increasingly popular in applications to Euler equation, compressible Navier-Stokes equations and incompressible flow problems [1-3]. The most attractive feature of finite volume formulation, is that the resulting solution would imply that the integral conservation of quantities such as mass, momentum, and energy is exactly satisfied over any group of finite volumes, and of course, over the whole domain. The conventional FVM, as the Galerkin finite element method (FEM), is not suitable for convection-dominated flows. In order to overcome this difficulty, discretizations have to be either of the central-difference type with explicit artificial diffusion or of the upwind type which are naturally dissipative. The applications of FVM have been greatly extended by introducing unstructured grid, because this kind of grid can easily fit complex boundaries, especially when using the triangular grid. There are roughly two kinds of unstructured grids in FVM. In the cell-centered approach, the flow variables are stored at the centroids of the cells and the control volume is simply the cell itself. In this case, the surface integrations are not accurate and time consuming, because no values on the surfaces are defined. The vertex-centered scheme is more accurate on the non-uniform meshes than its cell-centered counterparts. The flow variables are stored at the nodes and the control volume is taken to be the summation of part of the neighboring cells sharing a common node.
The accuracy and stability properties are analyzed in details, the numerical accuracy is verified by a purely convection flow.
Considering the following compressible Navier-Stokes equations
(1)
(2)
where,
is the velocity component, p
is the pressure, c is the speed of the
sound,
is
coordinate, t is time,
. Noting that the sound speed is
given by
, and by using the continuity
equation (2) , the following formulation can be derived
(3)
that is
(4)
equation (4) is another
form of continuity equation. From here, we assume that the variations of ρand
c are small enough to be regarded as constants, which is the assumptions of the
weakly compressible flow. Equations (2), (4) are the basic equations for
imcompressible flows, their dimensionless forms can be written as
(5)
(6)
where
Re is Reynolds number defined as
, L and U
are the reference length and reference velocity respectively. The superscript
asterisk is omitted in the following discussion.
The finite volume formulation of cell-centered form of equations (5) and (6) can be obtained as
(7)
where
and
are average values on the cell e,
is the cells volume,
and
are defined as
(8)
In equation (8), the first order spatial derivatives have to be determined previously, this can be done by the following method
(9)
where
is the interpolation function as used by finite element method. The volume
integration term
can be approximated by taking the average value
on a cell, that is
(10)
By using the triangular grid, the finite volume formulations of vertex-centered for a node i can be obtained as
(11)
(12)
where i, j, k
is the three node number of the cell,
stands for the summation of the cells having a common node. In
finite element method, the variation of unknowns within an element can be
considered by using the shape function. For example, when the triangular element
is used for two dimensional computation, the mass coefficient matrix can be
expressed as
(13)
where
,
are shape functions. Compared with
the finite element method, equations (11) and (12) can be modified as
(14)
(15)
Equations (14) and (15) are more accurate, as it
consideres the variation of a quantity within a cell. While equations (11) and
(12) simply take the quantity in a
cell as constant, this is not in agreement with the practical fact.
The modified equations (14) and (15) make the
present finite volume scheme has third-order spatial accuracy, which will be
discussed in the following section. Thereafter it is refered as consistent
finite volume scheme, while
equations (11) and (12) are called as conventional finite volume scheme.
Equations(14) and (15) can be written as the following form
(16)
where
represents
or
;
represents
or
. As equation(16) possesses
diagonal dominated property, it can be solved by the Jacobian iteration
(17)
where
is the lumping matrix and k is
iteration step. In general, two times of iteration is enough in the present
computations. From above processes, we can find that it is unnecessary to form
the global coefficient matrix and only cell level matrix is needed in the
computations. When this formulation is combined with the 3-step time
integrations, much more stable and accurate solution can be obtained for
convection dominated flows. The three-step time integrations can be written as
(18)
The present method is proposed to compute multi-dimensional problems, its stability properties in multi-dimensional case have to be made. As it is very complicated in deriving the amplification factor in three-dimensional case, the diffusion term of equation (19) is omitted and only uniform grid, as shown in Fig.1, is used in the analysis. The simplified form of equation(19) can be written as
(19)
are the velocity components in x,
y, z directions respectively.
Node
(i, j, k)
Fig.
1 3
3
3 regular mesh
Using the uniform cubic vortex-centered grid as
shown in Fig.1, whose size is
. The stability analysis can be obtained by performing Fourier analysis.
Assuming
(20)
where
,
,
,
,
,
,
,
are wave numbers in x, y, z
directions. The spatial discretized formulation for node (i, j, k)
is
(21)
where M is defined as
(22)
The details explanations for equation (22) and the processes of deriving the amplification factor can be found in Appendix A. The final amplification factor in three-dimensional case can be obtained as
(23)
The parameter Y is defined as
(24)
where
,
For two-dimensional
case, all the terms containing z or
should be eliminated from the above formula. And for one-dimensional case,
all the terms containing y, z or
,
should be eliminated from the above formula. The amplification factor is
shown in Figs.2 (a) and (b). The result of two-step scheme (Lax-Wendroof) is
also plotted on the same figure. It can be seen that the present method has
enlarged stability domain both in 2-D and 3-D cases.
(a) Three-dimensional case
(b)
Two-dimensional case
Fig. 2 Stability domain in multi-dimensional cases
The proposed finite volume method has been used
to compute scalar transportation problems, the purely convection examples is
used to verify the accuracy and the stability properties of the present scheme.
Considering one-dimensional purely convection flow, the initial concentration is
given by the following cosine hillmain and third order accuracy.
Fig. 3 Simulations of 1-D Cosine Wave
, with
and
(25)
(b). It can be seen that the two-step (Lax-Wendroof) lumping scheme is unstable, while the three-step lumping scheme is better than the former. However, the oscillations are still intolerable. Much more stable and accurate results can be obtained by the present consistent finite volume formulation, as shown in Fig.3 (c).
A new finite volume formulation, which using
fractional step time integration and consistent volume integration matrix, has
been proposed to solve the convection diffusion problems. Its stability
properties and accuracy has been analyzed both in one dimensional and
multi-dimensional cases. Compared with the Lax-Wendroff scheme, the present
finite volume formulation has enlarged computational domain and has third order
accuracy.
APPENDIX
A : The stability analysis for
three dimensional cases
In stability analysis part, equation (22) has been derived as
(26)
where
(A-1)
(A-2)
(A-3)
(A-4)
(A-5)
(A-6)
(A-7)
(A-8)
(A-9)
(A-10)
(A-11)
(A-12)
(A-13)
and
,
,
are defined as
(A-14)
(A-15)
(A-16)
Acknowledgments
This research was supported by the State key 973-Project for fundamental Research and Development Program(Grant No. G1999043607 ) and the National Science Foundation (Grant No. 59979013).
References
[1] Charles C. S. Song, Yuan, M.S. (1988). A weakly compressible flow model and rapid convergence methods. ASME, J. Fluids Eng., 110, 441-445.
[2] Mingham C. G., Causon D. M. (1998). High-resolution finite-volume method for shallow water flows. J. Hydraulic Eng., 124(6), 605-612.
[3] Jiang,C.B. and Kawahara,M. (1993). The analysis of unsteady incompressible flows by a Three-step finite element method. Int. J. for Numerical methods in Fluids, 16, 793-811.
[4] Crumpton, P.I. and Shaw, G.J. (1992). Cell vertex finite volume discretizations in three dimensions. Int. J. for Numerical meth. Fluids, 14, 505-527.
[5] Donea, J. (1984). A Taylor-Galerkin method for convective transport problems., Int. J. Numer. Meth. Eng., 20, 101-119,.