1.4. Unstructured meshes

So far we have used a Cartesian coordinate system to define the algebraic expressions. However, complex geometries cannot be represented accurately in such coordinate system. One advantage of the Finite Volume Method is its suitability for complex geometries, such as structured curvilinear meshes and specially for unstructured meshes. It is important to mention that FDM can also be adapted to such meshes but it is more complicated to implement it.

The essence of the Finite Volume Method is the conservation of a quantity by interpreting the influence of physical mechanisms as fluxes through the surfaces of a defined volume. That quantity is defined at the centroid and is representative of the entire volume. Rectangular shapes were used to illustrate a 1D convection problem but any shape can be used, “two” dimensions shapes (it is still a volume by assuming a uniform depth of the 2D shapes) like rectangles, triangles, pentagons, hexagons or non-regular shapes. It is evident that the number of surfaces with possible fluxes is different in every shape.

Consider a rectangular metal piece that has a 60 degree cut at both ends. How can you represent this problem using Finite Differences in a structured grid? For the FVM it is simple to just use triangle volumes.

Orthogonal and non-orthogonal meshes

For a mesh composed only of equilateral triangles, just like for the metal with a cut of 60 degrees, the line connecting the centroids of two contiguous volumes aligns with the normal vector of the surfaces, this volumes are said to be orthogonal. However, if this is not the case, then the mesh is referred to as non-orthogonal. In practice, having a perfectly orthogonal mesh is too restrictive and a degree of non-orthogonality (non-uniform) is needed to represent complex geometries. A correction is needed!

Special Treatment of Diffusion in non-orthogonal meshes

The flux due to diffusion in orthogonal volumes can be approximated adequately with a central difference scheme. However, the flux between non-orthogonal volumes (the line connecting the centroids and the normal to the surface are not parallel) needs a correction. Consider the formulation of the diffusion flux at the surface

\[ \int_{S_i} \mathbf{n_i} \cdot (D \nabla \phi)dS = \mathbf{n_i} \cdot (D \nabla \phi) \Delta S_i \]

where the subindex \(i\) refers to the surface number. The divergence operation projects the rate of change of \(\phi\) in the normal direction of the surface \(\mathbf{n_i}\):

\[ \mathbf{n_i} \cdot (D \nabla \phi)=D\frac{\partial \phi}{\partial n} \]

The central differencing to define the flux at the surface using the values of \(\phi\) at the centroids is:

\[ \frac{\partial \phi}{\partial d_c} \approx \frac{\phi_{V2}-\phi_{V1}}{\Delta d_{c}} \]

This clearly shows that for orthogonal meshes the central differencing is adequate since \(\mathbf{n}=\mathbf{d}_c\) therefore \(\frac{\partial \phi}{\partial n} = \frac{\partial \phi}{\partial d_c}\). However, the more non-orthogonal the mesh is, the larger the deviation between \(\mathbf{n}\) and \(\mathbf{d}_c\). This is why cross-diffusion is required.

Cross-Diffusion

To reveal the component of \(\frac{\partial \phi}{\partial d_c}\) that contributes to \(\frac{\partial \phi}{\partial n}\) we can rewrite \(\nabla \phi\) in terms of the unit normal and tangential vectors of the surface:

\[ \nabla \phi = \frac{\partial \phi}{\partial x} \mathbf{i} + \frac{\partial \phi}{\partial y} \mathbf{j} = \frac{\partial \phi}{\partial n} \mathbf{n} + \frac{\partial \phi}{\partial s} \mathbf{s_t} \]

\(\mathbf{s}_t\) being the unit vector of the surface tangent. Then, \(\nabla \phi\) can be projected in the direction of the line connecting the centroids:

\[ \frac{\partial \phi}{\partial d_c} = \nabla \phi \cdot \mathbf{d}_c = \frac{\partial \phi}{\partial n} \mathbf{n} \cdot \mathbf{d}_c + \frac{\partial \phi}{\partial s} \mathbf{s_t} \cdot \mathbf{d}_c \]

which can be interpreted in geometric terms using the angle between \(\mathbf{n}\) and \(\mathbf{d}_c\):

\[ \frac{\partial \phi}{\partial d_c} = \frac{\partial \phi}{\partial n} \cos\theta - \frac{\partial \phi}{\partial s} \sin \theta \]

from where it can be deduced

\[ \frac{\partial \phi}{\partial n} = \frac{\partial \phi}{\partial d_c} \frac{1}{\cos\theta} - \frac{\partial \phi}{\partial s} \frac{\sin \theta}{\cos\theta} \]

The final expression is:

\[ \int_{S_i} \mathbf{n_i} \cdot (D \nabla \phi)dS = (D \cdot \Delta S_i) \left(\frac{\partial \phi}{\partial d_c} \frac{1}{\cos\theta} - \frac{\partial \phi}{\partial s} \frac{\sin \theta}{\cos\theta}\right) \]

Conclusions of the implications of non-orthogonal meshes

There is a dependence on \(\theta\) which can be computed directly as it required the information of the known location of four points \(1,2,a,b\). The derivative \(\frac{\partial \phi}{\partial d_c}\) is computed with the information of the location points \(1,2\) and the values \(\phi_{V1},\phi_{V2}\) and \(\frac{\partial \phi}{\partial s}\) requires the location points \(a,b\) and the values \(\phi_{a},\phi_{b}\). So we need to approximate the value of \(\phi\) at the vertices! This last computation is specially costly and adds a potential error. Also, note that the intersection between the line connecting the centroids and the surface is not the midpoint of the surface, so the accuracy is reduced as consequence. These reasons highlight the importance of keeping the non-orthogonality of the mesh limited as much as possible.

Attribution

This chapter is written by Robert Lanzafame and Jaime Arriaga Garcia. Find out more here.