OpenFOAM® v2012: New and improved numerics

23/12/2020

New run-time geometry calculation selection

This release adds user-selectable geometry calculation. The system/fvSchemes dictionary now has an optional geometry section which overrides the method to calculate the geometric properties:

geometry
{
    // Use specialised method for high-aspect ratio cells
    type        highAspectRatio;
    minAspect   10;     // when to start blending lower-order method
    maxAspect   100;    // when to use lower-order method only
}

The selected fvGeometryScheme is responsible for calculating

faceCentres, faceAreas
cellCentres, cellVolumes
linear interpolation weights
cell-centre distance
non-orthogonality corrections

The following schemes are available:

basic: This is the default scheme and provides the exact behaviour as previous versions. It calculates the mass centroid for faces and cells.
stabilised: Face centres are calculated using only positive triangle contributions and is supposedly more stable on concave/distorted faces. This is the default in openfoam.org.
highAspectRatio: This scheme blends between basic and an edge-length weighted, face-area weighted average. This avoids truncation errors on high aspect ratio cells, at the cost of being lower order. The blending is linear across the range minAspect, maxAspect (see example above).
averageNeighbour: Similar to highAspectRatio but followed by a pass to align cell centres on top of one another to minimise non-orthogonality.

Specifying the method to calculate the geometric properties instead of e.g. limiting the effect on discretisation means that it also affects the mesh generation. This is especially important for mesh generators that use the mesh quality e.g.snappyHexMesh.

High-aspect ratio meshes

The averageNeighbour scheme can be very effective for tetrahedral and high-aspect ratio prismatic boundary layer meshes. These might have a very high non-orthogonality angle in the prism boundary layer due to truncation errors and alignment.

The figure below shows a comparison of the distribution of the non-orthogonality angle when computed with the basic scheme (left) and averageNeighbour (right).

[Picture]

Note:

all solvers and applications using fvMesh or its derivatives, e.g.dynamicFvMesh, should work with the user-selectable geometry schemes. Any application using the underlying polyMesh level, e.g.blockMesh, will not and use the basic functionality.
any user-coded dynamicMesh might have to be adapted to be able to use the new schemes. If not it will use the basic functionality.
the highAspectRatio and averageNeighbour schemes are still under active investigation and development.
checkMesh -writeAllFields provides an easy way to visualise the non-orthogonality and estimated aspect ratio using the current scheme.

Source code: $FOAM_SRC/finiteVolume/fvMesh/fvGeometryScheme
Tutorial: $FOAM_TUTORIALS/mesh/snappyHexMesh/airfoilWithLayers

Lagrangian: Improved particle tracking on moving meshes

This release includes several improvements to particle tracking on moving mesh cases. These cases would often fail in earlier versions due to:

particles crossing AMI patches;
particle-wall interactions; and
fragility of the barycentric tracking (ported from openfoam.org)

Improved particle-cyclicAMI patch interaction

When a particle hits a cyclicAMI patch it is transported to its neighbour cyclicAMI patch. This could lead to a situation where the particle continually ’bounces’ between both sides. The behaviour can be avoided using new fraction keyword, specified in the cyclicAMI patch entry of the constant/polyMesh/boundary file, which advances the particle inside the receiving cell by a specified fraction. Note that the particle is tracked from fraction zero to one.

By defauly, the fraction is set to zero for backwards compatibility. If non-zero the fraction is added to the tracking fraction when hitting the face.

Tutorial: $FOAM_TUTORIALS/lagrangian/MPPICDyMFoam/denseRotor2DAMI/system/blockMeshDict

Improved particle-moving-wall patch interaction

Particle-moving-wall interactions can under certain conditions, lead the tracking to stall. The problem arises when a particle hits a moving wall and the resulting particle velocity is extremely close to the wall velocity, and the particle final location is very close to the wall.

To counter this condition a new check has been added, where, if the final particle velocity after hitting the wall is below a threshold value the particle is deleted.

The behaviour is controlled by the new UrMax, set in the patchInteractionModel. It corresponds to the relative velocity between a particle and wall below which the particle is deleted. The default value of zero signifies that the particle should not be removed.

Tutorial: $FOAM_TUTORIALS/lagrangian/MPPICDyMFoam/denseRotor2DAMI/constant/kinematicCloudProperties

More robust tracking

The barycentric tracking was made more robust to address issues when the cell barycentric coordinate decomposition $det(A) \relax \special {t4ht=$ becomes negative. This improvement was ported from openFoam.org.

Modifications on MPPIC cloud

The MPPIC sub-models have been transferred from the MPPIC to the kinematic cloud. The MPPIC sub-models, i.e. packing, damping, etc. provide a velocity correction to the particles near the packing zone. To avoid additional tracking of these particles the correction is now applied in the kinematic cloud.

This improves performance in highly packed beds where the velocity correction can be considerable.

This change is transparent to the user and the MPPIC solvers (MPPICDyMFoam and MPPICFoam) are still available.

Source code: $FOAM_SRC/lagrangian/intermediate/clouds/Templates/KinematicCloud
Tutorial: $FOAM_TUTORIALS/lagrangian/MPPICDyMFoam/column/
$FOAM_TUTORIALS/lagrangian/MPPICDyMFoam/denseRotor2DAMI/

Function1: new frequency or period input

The sine and square Function1 classes now include support for frequency or period specifications. For slow oscillations it can be more intuitive to specify the period, e.g.

{
    type        sine;

    t0          3888000;   // 45 days
    period      31536000;  // 365 days

    scale       1;
    level       0;
}

In addition, a Cosine Function1 is now available.

Better tailoring of square Function1, with optional separate mark and space keywords, e.g.

{
    type        square;

    t0          0;         // 0 days
    period      31536000;  // 365 days
    mark        100;       // 100 days on
    space       265;       // 265 days off

    scale       1;
    level       0;
}

Function1: new time limits

The new limitRange class provides a wrapper for Function1s that enables users to set minimum and maximum limits on input values. For example, applying limits to a polynomial entry:

limitedPolyTest        limitRange;
limitedPolyTestCoeffs
{
    min         0.4;
    max         1.4;

    value       polynomial
    (
        (5 1)
        (-2 2)
        (-2 3)
        (1 4)
    );
}

Here, the return value will be: