Encyclopaedia Index



1. Introduction
2. Finite-Volume Equations
3. Basic Discretisation Schemes
3.1 Central-Differencing Scheme
3.2 Upwind-Differencing Scheme
3.3 Hybrid-Differencing Scheme
4. Higher Order Discretisation Schemes
4.1 Classification of Schemes
4.2 Flux-Limiter Formulation
4.3 Linear Higher-Order Schemes
4.4 Non-linear higher-order Schemes
5. Recommendations
6. Scheme Activation
7. Schemes in BFCs
8. Applications
9. Concluding Remarks


An important consideration in CFD is the discretisation of the
convection terms in the finite-volume equations.

The accuracy, stability and boundedness of the solution depends
on the numerical scheme used for these terms.

The default scheme used in PHOENICS for all variables is the
hybrid-differencing scheme (HDS), which employs:

The UDS is bounded and highly stable, but highly diffusive when
the flow direction is skewed relative to the grid lines.

The HDS is only marginally more accurate than the UDS, as the 2nd-
order CDS will be restricted to regions of low Peclet number.

The 2 approaches commonly used to remedy the problem of numerical
diffusion are:

For engineering problems, the necessary degree of grid refinement is
generally impractical as the UDS & HDS are sluggish to grid-
refinement tests.

Thus, schemes with higher-order truncation errors than UDS have
been proposed in an attempt to improve resolution.

PHOENICS V3.4 now has an extensive set of higher-order convection
schemes which are valid in all coordinate systems for the staggered-
grid option.

Two classes of higher-order scheme are now available:



These schemes may be used for both single- and two-phase flows,
although the facility has yet to be extended to include:

The other limitation is that there is no special treatment of the
domain and internal boundaries, at which all higher-order schemes
are reduced to the UDS.

The purposes of this lecture are:


The discretised equation for a general specific variable F is:

Jh - Jl + Jn - Js + Je - Jw

+ Dh - Dl + Dn - Ds + De - Dw = Sp (1)

where Sp is the source term for the control volume p, and Jf and
Df represent, respectively, the convective and diffusive fluxes
of F across the control-volume face f (f=h,l,n,s,e or w).

The convection fluxes through the cell faces are calculated as:

Jf = Cf*Ff (2)

where Cf is the mass flow rate across the cell face f.

The convected variable Ff is stored at the cell centres, and so
its value must be determined by interpolation, i.e. from a scheme.

In PHOENICS, the various higher-order schemes are implemented by
maintaining the UDS discretisation and adding extra higher-order
terms in the form of a source term.

Ff can be explicitly formulated in terms of its neighbouring nodal
values by a functional relationship of the form:

Ff = P(Fnb) (3)

where Fnb denotes the neighbouring-node F values.

From equations (1) through (3), the discretised equation becomes:

{ Dh + Ch*[P(Fnb)]h } - { Dl + Cl*[P(Fnb)]l } +

{ Dn + Cn*[P(Fnb)]n } - { Ds + Cs*[P(Fnb)]s } +

{ De + Ce*[P(Fnb)]e } - { Dw + Cw*[P(Fnb)]w } = Sp (4)

The higher-order (HO) schemes are introduced into (4) by using
the deferred-correction procedure, whereby:

Ff = Ff(U) + Ff' (5)
Ff' = Ff(H) - Ff(U) (6)

Here Ff' is a HO correction which represents the difference between
the UDS value Ff(U) and the HO value Ff(H).

If eqn (5) is substituted into eqn (4), the resulting discretised
equation is:

{ Dh + Ch*Fh(U) } - { Dl + Cl*Fl(U) } +

{ Dn + Cn*Fn(U) } - { Ds + Cs*Fs(U) } +

{ De + Ce*Fe(U) } - { Dw + Cw*Fw(U) } = Sp + Bp (7)

The deferred-correction source term, Bp, is given by:

Bp = Cl*Fl' - Ch*Fh' + Cs*Fs' - Cn*Fn' + Cw*Fw' - Ce*Fe' (8)

This treatment leads to a diagonally-dominant coefficient matrix
since it is formed using the UDS.

If eqn (8) is expanded in terms of nodal values, the final form of
the discretised equation is:

ap*Fp = sum (anb*Fnb) + Sp + Bp (9)

where: ap and anb are the convection-diffusion coefficients obtained
from the UDS: Fp is the cell-average value of F stored at the cell
centre; and the summation is over the immediate neighbouring nodes
nb (=L,H,S,N,W and E).

All of the schemes provided in PHOENICS calculate the cell-face
values Ff using at most two upwind cell-centre values and one
downwind cell-centre value, as shown in this Figure .

This stencil involves the upstream, central and downstream grid
points, designated by u, c and d respectively.


3.1 Central-Differencing Scheme

The most natural interpolation assumption for the cell-face
value of the convected variable Ff would appear to be the CDS,
which calculates the cell-face value from:

Ff = 0.5*(Fc + Fd) (10)

This scheme is 2nd-order accurate, but is unbounded so that
unphysical oscillations appear in regions of strong convection
and also in the presence of discontinuities, such as shocks.

The CDS may be used directly in very low Reynolds-number flows
where diffusive effects dominate over convection.

3.2 Upwind-Differencing Scheme

The UDS (see Courant et al [1952]) assumes that the convected
variable at the cell face f is the same as the upwind cell-centre

Ff = Fc (11)

The UDS is unconditionally bounded and highly stable, but it is
only 1st-order accurate in terms of truncation error.

The scheme is therefore highly diffusive when the flow direction
is skewed relative to the grid lines.

3.3 Hybrid-Differencing Scheme

The HDS (Spalding [1972]) switches between CDS and UDS according
to the local cell Peclet number, as follows:

Ff = 0.5*(Fc + Fd) for Pe < 2  (12)

Ff = Fc for Pe > 2

The cell Peclet number is defined as:

Pe = r*abs(Uf)*Af/Df (13)

in which Af=cell-face area and Df=physical diffusion coefficient.

When Pe > 2, CDS calculations tend to become unstable so that the
HDS reverts to the UDS. Physical diffusion is ignored when Pe > 2.

The HDS scheme is marginally more accurate than the UDS, because
the 2nd-order CDS will be used in regions of low Peclet number.


4.1 Classification of Higher-Order Schemes

LINEAR SCHEMES are those whose coefficients are not direct
functions of the convected variable when applied to a linear
convection equation.

Linear convection schemes of 2nd-order accuracy or higher may
suffer from unboundedness and are not unconditionally stable.

NON-LINEAR SCHEMES analyse the solution within the stencil and
adapt the discretisation to avoid any unwanted behaviour, such
as unboundedness.

These two types of scheme may be presented in a unified way by
use of the FLUX-LIMITER formulation.

4.2 Flux-Limiter Formulation

The FLUX-LIMITER formulation calculates the cell-face value of
the convected variable from:

Ff = Fc + 0.5*B(r)*(Fc-Fu) (14)

where B(r) is termed a limiter function, and the gradient ratio r
is defined as:

r = (Fd-Fc)/(Fc-Fu) (15)

The generalisation of this approach to handle non-uniform meshes
has been given by Waterson [1994].

From equation (14) it can be seen that B(r)=0 gives the UDS, and
B(r)=r gives the CDS.

4.3 Linear Higher-Order Schemes

PHOENICS provides the following linear higher-order schemes:

These are unified for implementation purposes as members of the
Kappa class of schemes:

B(r) = 0.5*{(1+K)*r+(1-K)} (16)

where: K = 1 gives CDS, K = 0.5 gives QUICK, K = -1 gives LUS,
K = 0 gives Fromm and K = 1/3 gives CUS.

These schemes are plotted in the Flux-Limiter Diagram (FLD) in the
next panel, which takes the form of a plot of B(r) against r.

The two main regions of the flux-limiter diagram are given by r<0,
indicating an extremum, and r>0 indicating monotonic variation.

Linear flux limiter diagram

4.4 Non-Linear Higher-Order Schemes

The following QUICK-based non-linear schemes have been provided:

B(r) = max(0,min(2*r,0.75*r+0.25,4)) (17)

B(r) = 2*(r+|r|)/(r+3) (18)

B(r) = max(0,min(2*r,0.25+0.75*r,0.75+0.25*r,2)) (19)

B(r) = r*(3r+1)/(r+1)**2 for r > 0
B(r) = 0. for r <= 0 (20)

PHOENICS provides the following Fromm-based non-linear schemes:

B(r) = max(0,min(2*r,0.5+0.5*r,2)) (21)

B(r) = (r+|r|)/(r+1) (22)

B(r) = 3*(r2+r)/{2.*(r2+r+1)} (23)

B(r) = (r2+r)/(r2+1) (24)

The remaining non-linear non-linear schemes in PHOENICS are:

and the following CUS-based flux limiters:

B(r) = 1.5*(r+|r|)/(r+2) (27)

B(r) = max(0,min(2*r,2*r/3+1/3,2)) (28)

Classical and CUS-based non-linear schemes

Most limiters fall into the following 2 categories:

Limiter functions are designed to fulfil particular boundedness
criteria, usually either the Total Variation Diminishing (TVD)
or Positivity conditions.

All flux-limited schemes provided in PHOENICS are positive, and
the following schemes are TVD: Koren, MUSCL, van Leer harmonic,
Minmod, Superbee and UMIST.

The LINEAR HO schemes offer good resolution, but are unbounded
and may produce unphysical oscillations, which can lead to severe
convergence problems. Therefore, it is recommended that NON-LINEAR
(bounded) schemes always be applied to:

For incompressible flows, it is recommended that a LINEAR HO scheme
is applied to the momentum equations, unless severe convergence
difficulties are encountered.

LINEAR HO schemes:

NON-LINEAR schemes:


It should always be possible to obtain convergence when using
the HO schemes from the very start of the calculation.

However, typically, the false time-step requirements are between
0.01 and 0.1 of the value required by the UDS or HDS.

If convergence proves particularly problematic, then it is
suggested that the user try restarting the calculation from a UDS
or HDS solution.


The default scheme for all variables is the HDS, which is
activated by the setting DIFCUT=0.5.

The UDS is activated for all variables by setting DIFCUT=0.

The HO schemes can be activated from the MENU or Q1; the MENU
has a default HO setting of van-Leer harmonic for all variables.

The PIL command SCHEME is used to select the required HO scheme
for particular variables.

The syntax of the SCHEME command is:

SCHEME(NAME,variable name 1,variable name 2,...etc.)

The 1st argument NAME identifies the required scheme, and the
2nd argument permits the user to specify those SOLVEd variables
which will use the selected scheme.

The 1st argument NAME identifies the required scheme, as follows:

(a) Linear schemes

(b) Non-linear schemes
or HCUS.

For example, the PIL commands:


select QUICK for U1 and V1, and SMART for H1,C1 and C2, and UDS
for any SOLVEd variables which do not appear in a SCHEME command.

If ALL is entered as the 2nd argument, then the selected scheme is
applied to all SOLVEd-for variables.


The schemes described thus far are also available for use with the default
staggered-mesh BFC option.

The SCHEME command and the associated FORTRAN coding in gxhocs.for are
implemented only for the staggered-grid option in PHOENICS.

The CCM BFC option uses its own set of higher-order convection schemes for
single and multi-block meshes. These are activated by using the READ Q1
facility, i.e. by inserting the following commands in the Q1 file:


CCM supports the following schemes:

QUICK - quadratic upwind scheme
SUPERB - superbee scheme

Please note that the conventional SCHEME COMMAND is not compatible with
the CCM option.

The GCV option has its own set of higher-order convection schemes for both
single-block and multi-block meshes. They are activated by inserting the
following type of command in the Q1 file:


etc. The 3rd argument is the variable name and the 5th argument can be
one of:

CDS - central differencing scheme
SOUP - second-order upwind scheme? i.e linear upwind
QUICK - quadratic upwind scheme
MUSCL - van Leer scheme, VANL1 in staggered-mesh schemes
SUPERB - superbee scheme

Please note that the conventional SCHEME COMMAND is not compatible with
the GCV option.


The 'Numerical-Algorithms' library contains a number of Q1 files
exemplifying the use of the HO schemes, including:

2d Diagonal Scalar Convection N101
2d Laminar Wall-Driven Cavity N102
2d Laminar Backward Facing Step N103
2d Scalar Convection with Recirculation N104
2d Turbulent Backward Facing Step N105
2d Laminar Flow Over a Thin Fence N108
2d Turbulent Flow through an Orifice Plate N110
1d Shocked transonic flow in a laval nozzle N131
2d Transonic underexpanded free jet

Any of these cases may be loaded from the LIBRARY MENU, or from
the SATELLITE by using,for example, LOAD(N110).

The results of some of these cases may be found in the POLIS
Applications Album, and some will now be presented here.


Some recent PHOENICS developments have been presented for HO
convection-discretisation schemes on staggered meshes.

The schemes comprised 5 linear schemes and 12 non-linear schemes.

Some successful applications of these schemes were reported.

More development work needs to be done, including:

Sources of further information are: