The stabilityBlendingFactor
function object computes the stabilityBlendingFactor
to be used by the local blended convection scheme. The output is a surface field weight between 0-1.
The weight of a blended scheme, i.e. w
, is given by a function of the blending factor, f
:
\[ w = f_{scheme_1} + (1 - f_{scheme_2}) \]
The factor is calculated based on six criteria:
The user can enable them individually.
For option 1, the following relation is used, where \(\phi_1\) is the non-orthogonality:
\[ fNon = min ( max ( 0.0, (\phi_1 - max(\phi_1)) /(min(\phi_1) - max(\phi_1)) ), 1.0 ) \]
For option 2, the following relation is used, where \(\phi_2\) is the magnitude of cell centres gradient (Note that \(\phi_2 = 3\) for orthogonal meshes):
\[ fMagGradCc = min ( max ( 0.0, (\phi_2 - max(\phi_2)) / (min(\phi_2) - max(\phi_2)) ), 1.0 ) \]
For option 3, a PID control is used in order to control residual unbounded fluctuations for individual cells.
\[ factor = P*residual + I*residualIntegral + D*residualDifferential \]
where P
, I
and D
are user inputs.
The following relation is used:
\[ fRes = (factor - meanRes)/(maxRes*meanRes); \]
where
meanRes | Average(residual) maxRes | User input
Note that \(f_{Res}\) will blend more towards one as the cell residual is larger then the domain mean residuals.
For option 4, the following relation is used, where \(\phi_4\) is the face weight (Note that \(\phi_4 = 0.5\) for orthogonal meshes):
\[ ffaceWeight = min ( max ( 0.0, (min(\phi_4) - \phi_4) / (min(\phi_4) - max(\phi_4)) ), 1.0 ) \]
For option 5, the following relation is used, where \(\phi_5\) is the cell skewness:
\[ fskewness = min ( max ( 0.0, (\phi_5 - max(\phi_5)) / (min(\phi_5) - max(\phi_5)) ), 1.0 ) \]
For option 6, the following relation is used:
\[ fCoWeight = min(max((Co - Co1)/(Co2 - Co1), 0), 1) \]
where
Co1 | Courant number below which scheme2 is used Co2 | Courant number above which scheme1 is used
The final factor is determined by:
\[ f = max(fNon, fMagGradCc, fRes, ffaceWeight, fskewness, fCoWeight) \]
An indicator (volume) field, named blendedIndicator
is generated if the log flag is on:
1 represent scheme1 as active, 0 represent scheme2 as active.
Additional reporting is written to the standard output, providing statistics as to the number of cells used by each scheme.
Operand | Type | Location |
---|---|---|
input | - | - |
output file | dat | $FOAM_CASE/postProcessing/<FO>/<time>/<file> |
output field | volScalarField | $FOAM_CASE/<time>/<outField> |
Example of the stabilityBlendingFactor
function object by using functions
sub-dictionary in system/controlDict
file:
stabilityBlendingFactor1 { // Mandatory entries (unmodifiable) type stabilityBlendingFactor; libs (fieldFunctionObjects); // Mandatory entries (unmodifiable) field <field>; // U; result <outField>; // UBlendingFactor; // Optional entries (runtime modifiable) tolerance 0.001; // Any of the options can be chosen in combinations // Option-1 switchNonOrtho true; nonOrthogonality nonOrthoAngle; maxNonOrthogonality 20; minNonOrthogonality 60; // Option-2 switchGradCc true; maxGradCc 3; minGradCc 4; // Option-3 switchResiduals true; maxResidual 10; residual initialResidual:p; P 1.5; I 0; D 0.5; // Option-4 switchFaceWeight true; maxFaceWeight 0.3; minFaceWeight 0.2; // Option-5 switchSkewness true; maxSkewness 2; minSkewness 3; // Option-6 switchCo true; U U; Co1 1; Co2 10; // Optional (inherited) entries writePrecision 8; writeToFile true; useUserTime true; region region0; enabled true; log true; timeStart 0; timeEnd 1000; executeControl timeStep; executeInterval 1; writeControl timeStep; writeInterval 1; }
where the entries mean:
Property | Description | Type | Required | Default |
---|---|---|---|---|
type | Type name: stabilityBlendingFactor | word | yes | - |
libs | Library name: fieldFunctionObjects | word | yes | - |
field | Name of the operand field | word | yes | - |
result | Name of surface field to be used in the localBlended scheme | word | yes | - |
switchNonOrtho | Select non-orthogonal method | bool | no | false |
nonOrthogonality | Name of the non-orthogonal field | word | no | nonOrthoAngle |
maxNonOrthogonality | Maximum non-orthogonal for scheme2 | scalar | no | 20 |
minNonOrthogonality | Minimum non-orthogonal for scheme1 | scalar | no | 60 |
switchGradCc | Select cell centre gradient method | bool | no | false |
maxGradCc | Maximum gradient for scheme2 | scalar | no | 2 |
minGradCc | Minimum gradient for scheme1 | scalar | no | 4 |
switchResiduals | Select residual evolution method | bool | no | false |
residual | Name of the residual field | word | no | initialResidual:p |
maxResidual | Maximum residual-mean ratio for scheme1 | scalar | no | 10 |
P | Proportional factor for PID | scalar | no | 3 |
I | Integral factor for PID | scalar | no | 0 |
D | Differential factor for PID | scalar | no | 0.25 |
switchFaceWeight | Select face weight method | bool | no | false |
faceWeight | Name of the faceWeight field | word | no | faceWeight |
maxFaceWeight | Maximum face weight for scheme1 | scalar | no | 0.2 |
minFaceWeight | Minimum face weight for scheme2 | scalar | no | 0.3 |
switchSkewness | Select skewness method | bool | no | false |
skewness | Name of the skewness field | word | no | skewness |
maxSkewness | Maximum skewness for scheme2 | scalar | no | 2 |
minSkewness | Minimum skewness for scheme1 | scalar | no | 3 |
switchCo | Select Co blended method | bool | no | false |
U | Name of the flux field for Co blended | word | no | U |
Co1 | Courant number below which scheme2 is used | scalar | no | 1 |
Co2 | Courant number above which scheme1 is used | scalar | no | 10 |
tolerance | Tolerance for number of blended cells | scalar | no | 0.001 |
The result
entry is the field which is read by the localBlended
scheme specified in fvSchemes
. This name is determined by the localBlended
class.
The inherited entries are elaborated in:
Usage by the postProcess
utility is not available.
Tutorial:
Source code:
History