modules/solid: Improved diffusion number calculation

The calculation of the diffusion number has been put into a form
consistent with finite-volume, rather than the finite-difference form
that was used previously.

This difference in formulations is analogus to that of the Courant
number in the fluid solvers. Whilst a textbook will typically define the
courant number as equal to 'U*deltaT/deltaX', in a finite-volume context
it is more appropriate to define it as 'Sum(phi)/V*deltaT' (where 'Sum'
is a sum over the cell's faces). Similarly, the finite-difference
Fourier number, 'kappa/rho/Cp*deltaT/deltaX^2', is more consistently
expressed for finite-volume as 'Sum(Sf*kappa*deltaX)/(V*rho*Cp)*deltaT'.

This makes the calculation of the diffusion number less sensitive to the
presence of small, poor quality faces, and therefore makes time-step
adjustment more robust on arbitrary polyhedral meshes.
This commit is contained in:
Will Bainbridge
2023-05-09 14:56:56 +01:00
parent 472451a0ca
commit b5142dca74

View File

@ -24,6 +24,7 @@ License
\*---------------------------------------------------------------------------*/ \*---------------------------------------------------------------------------*/
#include "solid.H" #include "solid.H"
#include "fvcSurfaceIntegrate.H"
#include "fvMeshMover.H" #include "fvMeshMover.H"
#include "localEulerDdtScheme.H" #include "localEulerDdtScheme.H"
#include "addToRunTimeSelectionTable.H" #include "addToRunTimeSelectionTable.H"
@ -61,19 +62,25 @@ void Foam::solvers::solid::correctDiNum()
: mag(thermo.Kappa())() : mag(thermo.Kappa())()
); );
const surfaceScalarField kapparhoCpbyDelta const volScalarField::Internal DiNumvf
( (
sqr(mesh.surfaceInterpolation::deltaCoeffs()) fvc::surfaceSum
(
mesh.magSf()
*fvc::interpolate(kappa) *fvc::interpolate(kappa)
/fvc::interpolate(thermo.rho()*thermo.Cp()) *mesh.surfaceInterpolation::deltaCoeffs()
)()()
/(mesh.V()*thermo.rho()()*thermo.Cp()())
*runTime.deltaT()
); );
DiNum = max(kapparhoCpbyDelta).value()*runTime.deltaTValue(); const scalar meanDiNum = gAverage(DiNumvf);
const scalar meanDiNum = const scalar maxDiNum = gMax(DiNumvf);
average(kapparhoCpbyDelta).value()*runTime.deltaTValue();
Info<< "Diffusion Number mean: " << meanDiNum Info<< "Diffusion Number mean: " << meanDiNum
<< " max: " << DiNum << endl; << " max: " << maxDiNum << endl;
DiNum = maxDiNum;
} }