A Study on Finite Volume Method and its Interpolation schemes
OBJECTIVE: To briefly discuss finite volume method and the need for interpolation schemes and flux limiters in finite volume method(FVM).
Theory:
In computational fluid dynamics, there are 3 commonly used numerical method based solvers to solve the governing equations .i.e Navier strokes equations.These three methods are:
1. Finite Difference method:
Finite difference methods are discritizations used for solving differential equations by approximating them with difference equations that finite difference approximate the derivatives.The properties and governing equations are considered to be concentrated at a single point. For FDMs the grid is linear and should be alligned in such a way that all the grid points lie in a straight line.
Drawbacks of Finite difference method:
1.1. We cannot use FDM directly to a body whose grid points are not alligned correctly along the x,y or z directions.These situations may occur due to complex designs,folds and cuts in a geometry.
1.2. Since grid points are not alligned, we cannot calculate the properties at the next adjacent point exactly,but we can use some correction factor in order to approximate the difference with that adjacent point.This method is known as Variable transformation method.This approach seeks to create a grid using other space variables such as eta and epsilon.These grids looks similar to cartesian grid.A suitable mapping function is defined and this function can be easily transformed through cartesian to new space variables and vise versa.The varaible transformation method gives an approximate solution of the property and also time consuming and complex.For complex profiles and moving parts this method is less preffered.Due to this drawback FVM is used widely.
2.. Finite element Method:
Finite element method subdivides a large system into smaller, simpler parts that are called finite elements.The finite element method finally results in a system of algebraic equations.This method approximates the unknown function over the domian.The simple equations that model these finite elements are then assembled into a larger system equations that models the entire problem.
3.. Finite volume method:
Finite difference method is applied over a volume so that the governing equations is valid over the entire volume .The properties are assumed to be concentrated at the geometric centre of the volume. The volume integrals in a partial differential equation containing a divergence term are converted into the surface integrals using the divergence theorem.The obtained terms are considered as fluxes at the surfaces of each finite volume because when we look through these discretised equations the flux entering a given volume is identical to that leaving the adjacent volume, which is known as conservative.FVM is widely used as it can be applied for complex geometry and moving parts.
Another advantage for the finite volume method is that it is easily formulated to allow for unstructured meshes.Fininte volume refers too the small volume surrounding each node point on a mesh.
In FVM instead of writing equation for a single point, we write equations for an entire volume.Also for a considered volume no assumptions are made over the shape of the volume.Inside this particular volume the flow is conserved,flow quantities are conserved,energy,mass momemtum is conserved.This is called as finite volume representation.
Solution process of fvm:
1.Discretize the solution domian by a grid(i.e. a set of finite volumes) and define the computational node at which variables are to be evaluated.
2. Apply the integral form of conservation law to each finite volume.
3. Approximate the surface and volume integrals using appropriate quadrature formulae in terms of the function values at computational nodes.
4. Collect the algebraic equations for all the finite volumes to obtain a system of algebraic equations in terms of unknown values of the variable at the computational nodes.
5. Solve the resulting system of algebraic equations to obtain values of the variable at each computational node.
Consider a steady state heat conduction equation with a source term(heat generated):
`frac{del}{delx}(alphafrac{delT}{delx})+S=0`
where, S is the source term
`alpha` is thermal diffusivity
`frac{delT}{delx}` is temperature gradient along x direction
The above governing equations is valid at every point in the domian, i.e. entire control volume.
If we consider a small control volume dV and it is multiplied with the governing equation and intergrated for entire control volume then that should also be equal to zero.i.e the integral form of governing equations.Taking the conservation equation validated at single point and if it is validated at single point then its is validated at entire volume and in this entire volume the properties are conserved.
Multiplying the equation by dv and intergrating it , we obtain
`int(frac{del}{delx}(alphafrac{delT}{delx})+S)dV =0`
`int(frac{del}{delx}(alphafrac{delT}{delx}))dV + intSdV=0`
since dV= A.dx
`intfrac{del}{delx}(alphafrac{delT}{delx})A.dx + intSdV=0`
`[alphaA(frac{delT}{delx})]_w^e + bar(S) dV=0`
`bar(S)` is the average value of source term(S) which when multiplied by dv will be equal to `intSdV`
`alpha_eAfrac{delT}{delx_e}-alpha_wAfrac{delT}{delx_w} +bar(S)dV=0`
Assuming `alpha` as constant and not a function as temperature
`alphaA[frac{T_e-T_p}{Delta x}]-alphaA[frac{T_p-T_w}{Deltax}]+bar(S) dV =0`
The source term which is to be added in a cell,if it is a function of temperature we have to linearize it.If we consider `alpha` as a function of temperature ,we have to linearize it as well.
(`bar(S)` linearization = `S_u+S_pT_p`)
`alphaA[frac{T_e-T_p}{Delta x}]-alphaA[frac{T_p-T_w}{Deltax}]+bar(S) dV =0`
[Heat flux]out - [Heat flux]in + source term=0
Therefore the amount of heat removed from a particular box is conserved.
Different schemes are used to solve the equation:
Gradient scheme:
From the above figure when we calculate temperature gradient on the east face,we basically use temperature from neighbouring point trying to evaluate or numerically represent gradient of temperature which is called gradient scheme.
Interpolation scheme:
We have seen that we need evaluation of surface and volume integrals in our finite volume formulation, but the surface as well as volume integrals they require the value of the unknown variable at locations other than the computational notes of the CV. So somehow we have to cover with some sort of schemes or mathematical approximation to compute the values at these locations, and that is why we need various types of interpolation formulae.
Various types of interpolation schemes are:
1. Upwind interpolation
2. Linear interpolation
3. Quadratic upwind interpolation(QUICK)
4. Hybrid, TVD(total variation dimension) and ENO(essentially non-oscillatory) interpolation scheme.
1. Upwind interpolation:
upwind interpolation scheme is essentially used for what we called convection dominated problems, and here the velocity at a given point direction of that velocity that determines the variable values to be used in the interpolation process, the name upwind refers to that we are going to make use of the node which is upwind or upstream of the given node.
The upwind interpolation(UDS) for approximating the value of variable f at the east face of the control volume is given by,
UDS is equivalent using a backward or forward finite difference approximation depending on the flow direction.
Advantages:
Upwind interpolation is the only approximation that satisfies boundedness criterion unconditionally, unconditionally in the sense that respect irrespective of what the velocity values are and irrespective of what is the grid spacing the solution will always remain bounded that is this use of upwind interpolation will never yield oscillatory solutions.
Disadvantages:
1. UDS is numerically diffusive.i.e. if we have got sharp peaks or oscillatory values which are actually present in the solution those would be damped out because of the presence of numerical diffusion, so the numerical diffusion is magnified in multi-dimensional flows especially of the flow is oblique to the grid space.
2. Peaks or rapid variations in the variables will be smeared out because of the introduction of additional diffusive coefficient what we called numerical diffusive coefficient by the use of upwind interpolation schemes.
3. First order accurate, hence very accurate grids are required to obtain accurate solutions.
2. Linear interpolation(CDS):
Approximate the value of the variable at control volume face centre by linear interpolation of the values at two nearest computational nodes.
`f_e= f_Elamda_e+f_p(1-lamda_e)`
where, `lamda_e= frac{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).
Characteristics:
1. CDS scheme is second order accurate.
2. It may produce oscillatory solutions therefore small grids is preferred.
3. Simplest and most widely used interpolation.
4. leads to a very simple approximation for evaluation of gradients, and gradients are the ones which we require in the evaluation of the diffusive fluxes.
5. Aproximations are second order accurate for uniform grids and it is formally first order accurate for non-uniform grids.
3. Quadratic upwind interpolation(QUICK)
Aproximate the value of the variable at control volume face centre by quadratic interpolation of the values at the three nearest computational nodes(one downstream node D and two upstream nodes U,UU).
`f_e= f_U+g_1(f_D-f_U)+g_2(f_(UU)-f_U)`
`g_1= frac{(x_e-x_U)(x_e-x_(UU))}{(x_D-x_U)(x_D-x_(UU))}`
`g_2= frac{(x_e-x_U)(x_e-x_D)}{(x_(UU)-x_U)(x_(UU)-x_D)}`
Quadratic of Upwind Interpolation for Convective Kinematics, convective kinematics because this particular interpolation would be used for convection dominated problem, and quadratic because we have obtain this interpolation based on what quadratic interpolation formula.
this particular scheme it is third order accurate in both uniform and non-uniform grids
4. Hybrid interpolation scheme
Blend of CDS and UDS scheme based on local peclet number
`f_e= gammaf_(CDS)+(1-gamma)f_(UDS)`
Flux Limiters:
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 unnstable in nature and show oscillations in the vincinity of discontinuity. Highly accurate and oscillation free schemes are known as high resolution schemes.
Flux limiter functions are used to tune higher order and low order schemes in such a way that the resulting scheme gives a higher order accuracy in the smooth region of the flow and maintains first order accuracy in the viscinity of shocks and discontinuities.For such a scheme TVD( total variation diminishing)criteria is maintained.Their basic forms are shown below,
where, `f^(low)`= low precision flux(1st order accurate)
`f^(high)`= high precision flux(higher order accurate)
`phi(r)`= flux limiter function where `r=frac{u_i-u_(i-1)}{u_(i+1)-u_i`
common forms of flux limiter functions:CHARM,HCUS,HQUICK,koren,minmod etc.