Flux Limiters and Interpolation Schemes in Finite Volume Method

I. Interpolation Schemes in Finite Volume Method:

The approximation of surface and volume integrals may require values of the variable at locations other than the computational nodes of the CV. Values at these locations are obtained using interpolation formulae. Some of the possibilities are:  

1. Upwind Interpolation Scheme: (UDS)

The name upwind refers to the assumption that we are going to use the node that is upwind or upstream of the given node. it is generally used for convection dominated problems (i.e. High Peclet Number `Pe >> 1` flows). For the east face of a quantity it is given as:

`f_e = {(f_P if (v.n)_e > 0) , (f_E if (v.n)_e < 0):}`

The interface values are given by:  

for `v>0` we have: `u_(i-1/2) ~~ u_(i-1) ; u_(i+1/2)~~u_i`

`I_c ~~ v(u_i - u_(i-1))/(Deltax)` (Backward Differencing)

for `v<0` we have: `u_(i-1/2) ~~ u_i ; u_(i+1/2) ~~ u_(i+1)`

`I_c ~~ v(u_(i+1) - u_i)/(Deltax)` (Forward Differencing)

The resulting taylor series expansion is:

`vu_(i+1/2) = vu_i - v(Deltax)/2*((delu)/(delx))_(i+1/2)` where the second term is the leading truncation error which resembles a diffusive flux (Numerical Diffusion)

Upwind Differencing Scheme is the only approximation that satisfies boundedness criterion unconditionally. It is first order accurate in space.

2. Central Differencing Approximation: (CDS)

Approximate value of the variable at the Control Volume face centre by linear interpolation of the values at two nearest computational nodes:

`f_e = f_Elambda_e + f_P(1-lambda_e)` where, `lambda_e = (x_e - x_P)/(x_E - x_P)`

The linear interpolation is equivalent to the use of central difference formula of the first order derivative, and hence, this scheme is also termed as Central Difference Scheme (CDS).

Interpolation Polynomial:

`p_1(x) = u_L*(x_R - x)/(x_R - x_L) + u_R*(x - x_L)/(x_R - x_L)`

Average Interface Values are given by:

`u_(i-1/2) ~~ (u_(i-1)+u_i)/2` and `u_(i+1/2) ~~ (u_i + u_(i+1))/2`

`I_c ~~ v*(u_(i+1) - u_(i-1))/(2Deltax)` Central Differencing

Taylor Series Approximation:

`u_(i+1) = u_(i+1/2) + (Deltax)/2*((delu)/(delx))_(i+1/2) + (Deltax)^2/8*((del^2u)/(delx^2))_(i+1/2) + ....`

 `u_i = u_(i+1/2) - (Deltax)/2*((delu)/(delx))_(i+1/2) + (Deltax)^2/8*((del^u)/(delx^2))_(i+1/2) + ....`

 From the above expressions we have:

`u_(i+1/2) = (u_i + u_(i+1))/2 - (Deltax)^2/8*((del^2u)/(delx^2))_(i+1/2)+...`

This approximation is second order accurate. It is the simplest and most widely used interpolation forevaluation of gradients required for evaluation of diffusive fluxes.

3. Quadratic Upwind Differencing Scheme: (QUICK)

QUICK Stands for Quadratic Upwind Interpolation scheme for Convective Kinematics. It is a third order flux approximation. Although it has second order overall accuracy. It approximates the valuess of the variable at CV face centre by Quadratic Interpolation of the values at three nearest computational nodes (one Downstream node D and two upstream nodes U and UU)

`f_e = f_E + g_1(f_D- f_U) + g_2(f_(UU) - f_U)`

where `g_1` and `g_2` is given by:

`g_1 = ((x_e - x_U)(x_e - x_(UU)))/((x_(UU) - x_U)(x_(UU) - x_D))`

`g_2 = ((x_e - x_U)(x_e - x_D))/((x_(UU) - x_U)(x_(UU) - x_D))`

Other Interpolation schemes:

1. Hybrid Interpolation Schemes: (CDS and UDS)
2. Total Variation Diminishing (TVD) Scheme


1. Apart from QUICK all other interpolation schemes would produce inconsistent fluxes when quadratic interpolation formula is used. The results obtained from the schemes do not satisfy the overall equation. An example for the same is given below:


Piecewise Linear Interpolation functions must be used for evaluating the fluxes consistently as shown below:



For convective fluid flow it has been observed that low - order schemes are usually stable but quite dissipative in nature around points of discontinuity or shocks while the higher order schemes are unstable in nature and show oscillations in the vicinity of discontinuity. Highly accurate and Oscillation free schemes are known as High Resolution Schemes.

Flux limiter functions are used to tune high order and low order schemes in such a way that the resulting scheme gives a high order accuracy in the smooth region of the flow and maintains first order accuracy in the vicinity of shocks and discontinuities. For such a scheme TVD (Total Variation Diminishing) criteria is maintained. Their basic forms are shown below:

`F(u_(i+1/2)) = f^(low)(i+1/2) - phi(r_i)(f^(low)(i+1/2) - f^(high)(i+1/2))`

`F(u_(i-1/2)) = f^(low)(i-1/2) -phi(r_(i-1))(f^(low)(i-1/2) - f^(high)(i-1/2))`

where `f^(low)` = low precision flux (1st order accurate)

and `f^(high)` = high precision flux (higher order accurate)

`phi(r)` = Flux limiter function where `r = (u_i-u_(i-1))/(u_(i+1) - u_i)`

Common Forms of Flux Limiter Functions: CHARM, HCUS, HQUICK, Koren, minmod etc.

Projects by Priyotosh Bairagya

LAMINAR INCOMPRESSIBLE FLOW SIMULATION THROUGH A PIPE IN OPENFOAM: [PART 2/2]: SYMMETRY BOUNDARY CONDITION: OBJECTIVES:1. Creating the Mesh-Script (blockMeshDict file) with symmetry boundary condition.2. Simulation Results and post processing the velocity profile at Read more

LAMINAR INCOMPRESSIBLE FLOW SIMULATION THROUGH A PIPE IN OPENFOAM: [PART 1/2]: WEDGE BOUNDARY CONDITION: OBJECTIVES:1. Calculation of Pre-eliminary quantinites related to the flow.2. Creating the Mesh script (blockMeshDict file) for specifying the geometry and bounda Read more

ANALYSIS OF NUMERICAL STABILITY OF VARIOUS ITERATIVE SOLVERS FOR TRANSIENT 2D HEAT CONDUCTION: [Part: 3/3] INTRODUCTION: The criterion of stability of a numerical scheme is determined by the way the errors propagate while the solution moves from one time-step to th Read more

MESH GENERATION AND ANALYSIS USING BLOCKMESH FOR FLOW OVER A BACKWARD FACING STEP: The purpose of the following project is to generate the geometry for a variation of the incompressible cavity flow problem in OpenFOAM. For this purpose we have modified the lid-driven c Read more

ANALYSIS OF VARIOUS ITERATIVE SCHEMES FOR THE SOLUTION OF A 2D HEAT CONDUCTION PROBLEM: [PART: 2/3] In the previous part we had explained the problem statement and the MATLAB Program in detail. In this Part we are going to explain the outputs from the 2D Heat Conduct Read more

NUMERICAL SOLUTION OF 1D SUPERSONIC NOZZLE FLOW SIMULATION BY MAC-CORMACK METHOD PROJECT OBJECTIVES: i. Numerical solution of the governing equations in both conservative and non-conservative forms. ii. Creating user defined functions for calculating the flow quan Read more

Simulation of a 2D Heat Conduction problem in steady and unsteady/transient forms using iterative methods. Project Objectives: 1. Solving the 2 Dimensional Heat conduction equation in the generalized form using various iterative techniques: i. Explicit Solver (for Read more


WEEK 4: (Effect of Grid-Size on output for the solution of 1D linear wave equation) 1. Problem Setup: Given Partial Differential Equation: `(delu)/(delt) + c(delu)/(delx) = 0` Numerical Discretization:  `u_(i,n+1) = u_(i,n) + (cDeltat)/(Deltax)*(u_(i-1,n+1) Read more


The End