1.1.9.4. 2D Finite Volume Method: Non-Cartesian Grids

1.1.9.4.1. Summary

  • The cross product is used to compute the volume

  • The midpoint rule is used to compute the flow quantity and the source term

  • The FVM with the cell-centred or cell-vertex formulation can be used to compute the fluxes

  • If the FVM is applied on a Cartesian grid (every angle is a right angle), the FD formulas are recovered

1.1.9.4.2. Non-Cartesian Grid

Non-Cartesian, although shown as Cartesian, but we will not assume we know the volume of the cell ABCD

../_images/centred_scheme_2.png

For cell ABCD:

\[{\partial \over \partial t} \int_{\Omega_{ij}} U d \Omega + \oint_{ABCD} \mathbf{F} \cdot d \mathbf{S} = \int_{\Omega_{ij}} Q d \Omega\]
\[{\partial \over \partial t} U_{ij} d \Omega_{ij} + \sum_{ABCD} (\pm fdy \pm gdx) = Q_{ij} d \Omega_{ij}\]

where: \(f\) and \(g\) are the Cartesian components of the flux vector \(\mathbf{F}\)

The sign in the flux appears due to the surface vector, e.g for Side AB:

\[\mathbf{S}_{AB} = \Delta y_{AB} \mathbf{i} - \Delta x_{AB} \mathbf{j} = (y_B - y_A) \mathbf{i} - (x_B - x_A) \mathbf{j}\]

1.1.9.4.3. Evaluation of the Flow Quantity Terms

Can use midpoint rule which approximates a volume integral by the product of the centre value and the CV:

\[U_{i,j} = \int_{\Omega_{i,j}} U d \Omega = U_{i,j} \Delta \Omega_{i,j}\]

1.1.9.4.4. Evaluation of the Source Terms

Can use midpoint rule which approximates a volume integral by the product of the centre value and the CV:

\[Q_{i,j} = \int_{\Omega_{i,j}} Q d \Omega = Q_{i,j} \Delta \Omega_{i,j}\]

1.1.9.4.5. Evaluation of the Volumes

  • The vector product is the area of a parallelogram

  • Take the absolute value to ensure the value is positive (\(\Delta x\) and \(\Delta y\) can be negative)

\[\begin{split}\left| \mathbf{a} \times \mathbf{b} \right| = \left| \text{det} \begin{bmatrix} \Delta x_1 & \Delta x_2 \\ \Delta y_1 & \Delta y_2 \end{bmatrix} \right| = \left| \Delta x_1 \Delta y_2 - \Delta x_2 \Delta y_1 \right|\end{split}\]

  • The area we are looking for (the quadrilateral) is half the area of the parallelogram:

\[\Omega_{ABCD} = {1 \over 2} \left| \mathbf{AC} \times \mathbf{BD} \right|\]
\[\begin{split}\mathbf{AC} = \begin{bmatrix} \Delta x_{AC} \\ \Delta y_{AC} \end{bmatrix}\end{split}\]
\[\begin{split}\mathbf{BD} = \begin{bmatrix} \Delta x_{BD} \\ \Delta y_{BD} \end{bmatrix}\end{split}\]
\[\Omega_{ABCD} = {1 \over 2} \left| \Delta x_{AC} \Delta y_{BD} - \Delta x_{BD} \Delta y_{AC} \right|\]
../_images/parallelogram.png

1.1.9.4.6. Evaluation of the Fluxes

The evaluation of fluxes determines the difference between a variety of schemes, e.g. cell vertex or cell centered

It depends on the location of the flow variables w.r.t. the mesh and selected scheme

1.1.9.4.6.1. Central Scheme and Cell Centered FVM

We know the values at the black dots, we don’t know the values in red

../_images/centred_scheme_3.png

  • Option 1: Average of fluxes perpendicular to face - e.g. midpoint rule

\[f_{AB} = f_{i+1/2,j}\]

Take the average flux:

\[f_{i+1/2,j} = {1 \over 2} (f_{i,j} + f_{i+1,j})\]

where \(f_{i,j} = f(u_{i,j})\)

Flux of the average flow quantity:

\[f_{AB} = f \left( {u_{i,j} + u_{i+1,j} \over 2} \right)\]

(not the same, because f is generally non-linear)

Can also use Upwind, Linear Interpolation, QUICK and Higher Order Schemes for midpoint rule

  • Option 2: Average of fluxes parallel to face A and B - e.g. trapezoid rule

\[f_{AB} = {1 \over 2} (f_A + f_B)\]

Evaluate flow quantity in A and B (at the corner A)

\[u_A = {1 \over 4}(u_{i,j} + u_{i+1,j} + u_{i+1,j-1} + u_{i,j-1})\]

Or average the fluxes:

\[f_A = {1 \over 4}(f_{i,j} + f_{i+1,j} + f_{i+1,j-1} + f_{i,j-1})\]

(more flux evaluations - could be expensive computationally)

Can also use Simpson’s Rule, cubic polynomials for this type of rule

1.1.9.4.6.2. Central scheme and Cell Vertex FVM

We know the values at the black dots, we don’t know the values in red - remember flux is not at a point, it’s through a face by integration.

../_images/vertex_scheme_4.png

  • Option 1 - fluxes perpendicular to face:

\[f_{AB} = f \left( {u_{i,j} + u_{i+1,j} \over 2} \right)\]

  • Option 2 - fluxes parallel to face:

\[f_{AB} = {1 \over 2} (f_A + f_B)\]

  • The last one corresponds to trapezoidal rule for the integral:

\[\int_{AB} f dy = {1 \over 2}(f_A + f_B)(y_B - y_A)\]

Summing over the sides ABCD:

../_images/vertex_scheme_6.png

Use trapezium rule: half the sum of the parallel sides times the distance between them:

\[\begin{split}\oint_{ABCD} \mathbf{F} \cdot d \mathbf{S} = {1 \over 2} \left( \begin{bmatrix} {f_A - f_C \\ g_C - g_A} \end{bmatrix} \cdot \begin{bmatrix} {\Delta y_{DB} \\ \Delta x_{DB}} \end{bmatrix} + \begin{bmatrix} {f_B - f_D \\ g_B - g_D} \end{bmatrix} \cdot \begin{bmatrix} {\Delta y_{AC} \\ \Delta x_{AC}} \end{bmatrix} \right) \\ = {1 \over 2} \left( (f_A-f_C) \Delta y_{DB} + (g_C-g_A) \Delta x_{DB} + (f_B-f_D) \Delta y_{AC} + (g_B-g_D) \Delta x_{AC} \right)\end{split}\]

1.1.9.4.7. 2D Example - FVM Reduces to FD

This is a generic example - we could use cell centered or cell vertex here (but just showing cell centered):

../_images/centred_scheme_6.png

For cell ABCD:

\[{\partial \over \partial t} \int_{\Omega_{i,j}} U d \Omega + \oint_{ABCD} \mathbf{F} \cdot d \mathbf{S} = \int_{\Omega_{i,j}} Q d \Omega\]
\[\begin{split}{\partial \over \partial t} u_{i,j} \Omega_{i,j} + \begin{bmatrix} {f_{AB} \\ g_{AB}} \end{bmatrix} \cdot \begin{bmatrix} {\Delta y_{AB} \\ \Delta x_{AB}} \end{bmatrix} + \begin{bmatrix} {f_{BC} \\ g_{BC}} \end{bmatrix} \cdot \begin{bmatrix} {\Delta y_{BC} \\ \Delta x_{BC}} \end{bmatrix} + \begin{bmatrix} {f_{CD} \\ g_{CD}} \end{bmatrix} \cdot \begin{bmatrix} {\Delta y_{CD} \\ \Delta x_{CD}} \end{bmatrix} + \begin{bmatrix} {f_{DA} \\ g_{DA}} \end{bmatrix} \cdot \begin{bmatrix} {\Delta y_{DA} \\ \Delta x_{DA}} \end{bmatrix} = q_{i,j} d \Omega_{i,j}\end{split}\]

Volume:

\(\Omega_{i,j} = \Delta x \Delta y\)

Surfaces:

\(\Delta y_{AB}\) = \(\Delta y\) \(\quad\) \(\Delta x_{AB}\) = \(0\) \(\quad\) \(\Delta y_{BC}\) = \(0\) \(\quad\) \(\Delta x_{BC}\) = \(\Delta x\) \(\quad\) \(\Delta y_{CD}\) = \(\Delta y\) \(\quad\) \(\Delta x_{CD}\) = \(0\) \(\quad\) \(\Delta y_{DA}\) = \(0\) \(\quad\) \(\Delta x_{DA}\) = \(\Delta x\)

Fluxes:

\(f_{AB}\) = \(f_{i+1/2,j}\) \(\quad\) \(g_{AB}\) = \(0\) \(\quad\) \(f_{BC}\) = \(0\) \(\quad\) \(g_{BC}\) = \(g_{i,j+1/2}\) \(\quad\) \(f_{CD}\) = \(-f_{i-1/2,j}\) \(\quad\) \(g_{CD}\) = \(0\) \(\quad\) \(f_{DA}\) = \(0\) \(\quad\) \(g_{DA}\) = \(-g_{i,j-1/2}\)

Hence:

\[{\partial \over \partial t} u_{i,j} \Delta x \Delta y + ({f_{i+1/2,j} - f_{i-1/2,j}}) \Delta y + ({g_{i,j+1/2} - g_{i,j-1/2}}) \Delta x = q_{i,j} \Delta x \Delta y\]

Dividing by the volume we have recovered the finite difference formula:

\[{\partial \over \partial t} u_{i,j} + {({f_{i+1/2,j} - f_{i-1/2,j}}) \over \Delta x} + {({g_{i,j+1/2} - g_{i,j-1/2}}) \over \Delta y} = q_{i,j}\]

Note: \(f_{i,j}\) and \(g_{i,j}\) do not appear, hence odd-even decoupling can occur