FVM Interpolation schemes and flux limiters

Interpolation Schemes

In the finite volume formulation, evaluation of surface and volume integrals sometimes requires the values of unknown variables at locations other than the computational nodes of the control volume. These values are obtained by employing certain interpolation schemes. As an example let us consider the simple case of steady one-dimensional convection and diffusion.

`d/dx(\\rhou\\phi)=d/dx(\\Gamma(d\\phi)/dx)`

where `\\phi` is our solution variable (a scalar), `\\rho`is the density, `u` is the velocity and `\\Gamma` is the diffusion coefficient.

The above figure shows the three neighbouring control volumes, which we are going to use for deriving the discretization equation.

Integration of the governing differential equation over the control volume around grid point P yields:

`(\\rhou\\phi)_e-(\\rhou\\phi)_w=(\\Gamma(d\\phi)/dx)_e-(\\Gamma(d\\phi)/dx)_w`

Evaluation of left hand side of the above equation requires the values of `\\phi` over faces `e` and `w`. These values are estimated by employing certain interpolation schemes.

Some interpolation schemes used in CFD are as follows:

1. Upwind interpolation scheme

2. Linear interpolation scheme

3. Quadratic Upwind Interpolation for Convective Kinematics (QUICK) scheme

4. Hybrid Scheme

5. Total Variation Diminishing (TVD) scheme

For demonstration purposes, we will explain the Upwind interpolation scheme.

The Upwind interpolation scheme takes into account the direction of fluid flow and based on it, the value of solution variable at the face is approximated by its value at the upstream computational node. For instance, if the fluid is moving from west to east (refer figure), then:

`\\phi_e\\approx\\phi_P`

And if the fluid is moving from east to west (refer figure), then:

`\\phi_e\\approx\\phi_E`

Flux Limiters

Numerical approximation using higher order spatial discretization schemes results in better accuracy but exhibits spurious oscillations near discontinuities such as shocks. On the other hand, lower order spatial discretization schemes remain stable near discontinuities but dampen out the distributions because of their dissipative nature. Hence in order to capture discontinuities, high resolution schemes are desired - which are highly accurate and stable throughout the domain.

Flux limiters operate functionally as a switch and tune the numerical fluxes obtained from lower order scheme and higher order scheme such that near the discontinuity lower order (stable) scheme is employed and away from it higher order (accurate) scheme is employed.

Suppose `F_e` represents flux at face `e` calculated using high resolution scheme, `L_e` represents flux at `e` calculated using lower order scheme and `H_e` represents flux at `e` calculated using higher order scheme. Then `\\psi(r)` represents the flux limiter function such that:

`F_e=L_e+\\psi(r)(H_e-L_e)`

Flux limiter function depends upon the smoothness parameter `r` which is defined as:

`r_P=(\\phi_P-\\phi_W)/(\\phi_E-\\phi_P)`

The flux limiter function is constrained to be greater than or equal to zero. It assumes a value of 0 near discontinuities and a value of 1 away from them. Hence it switches between higher order scheme and lower order scheme for numerical representation of flux, thereby keeping the scheme stable and accurate throughout.