zhyang-dev/OpenFOAM-12
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

fluidSolvers: Stabilise rDeltaT calculation

Browse Source 
When the flow is stationary (e.g., at the beginning of a run) the
rDeltaT calculation now requires a maxDeltaT setting in the PIMPLE
sub-section of the fvSolution dictionary. This prevents floating point
errors associated with rDeltaT approaching zero.
This commit is contained in:
Will Bainbridge
2023-07-27 11:42:51 +01:00
parent 9cba18a073
commit 2f579ca041
6 changed files with 91 additions and 48 deletions
Download Patch File Download Diff File Expand all files Collapse all files

17
applications/modules/VoFSolver/setRDeltaT.C
View File

@ -42,20 +42,25 @@ void Foam::solvers::VoFSolver::setRDeltaT()
const volScalarField rDeltaT0("rDeltaT0", rDeltaT);
// Set the reciprocal time-step from the local Courant number
// and maximum and minimum time-steps
rDeltaT.ref() = fvc::surfaceSum(mag(phi))()()/((2*maxCo)*mesh.V());
if (pimpleDict.found("maxDeltaT"))
// Clip to user-defined maximum and minimum time-steps
scalar minRDeltaT = gMin(rDeltaT.primitiveField());
if (pimpleDict.found("maxDeltaT") || minRDeltaT < rootVSmall)
{
rDeltaT.max(1/pimpleDict.lookup<scalar>("maxDeltaT"));
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("maxDeltaT");
rDeltaT.max(clipRDeltaT);
minRDeltaT = max(minRDeltaT, clipRDeltaT);
}
if (pimpleDict.found("minDeltaT"))
{
rDeltaT.min(1/pimpleDict.lookup<scalar>("minDeltaT"));
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("minDeltaT");
rDeltaT.min(clipRDeltaT);
minRDeltaT = min(minRDeltaT, clipRDeltaT);
}
Info<< "Flow time scale min/max = "
<< gMin(1/rDeltaT.primitiveField())
<< ", " << gMax(1/rDeltaT.primitiveField()) << endl;
<< gMin(1/rDeltaT.primitiveField()) << ", " << 1/minRDeltaT << endl;
setInterfaceRDeltaT(rDeltaT);

26
applications/modules/incompressibleFluid/setRDeltaT.C
View File

@ -53,24 +53,32 @@ void Foam::solvers::incompressibleFluid::setRDeltaT()
const volScalarField rDeltaT0("rDeltaT0", rDeltaT);
// Set the reciprocal time-step from the local Courant number
// and maximum and minimum time-steps
rDeltaT.ref() = fvc::surfaceSum(mag(phi))()()/((2*maxCo)*mesh.V());
if (pimpleDict.found("maxDeltaT"))
// Clip to user-defined maximum and minimum time-steps
scalar minRDeltaT = gMin(rDeltaT.primitiveField());
scalar maxRDeltaT = gMax(rDeltaT.primitiveField());
if (pimpleDict.found("maxDeltaT") || minRDeltaT < rootVSmall)
{
rDeltaT.max(1/pimpleDict.lookup<scalar>("maxDeltaT"));
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("maxDeltaT");
rDeltaT.max(clipRDeltaT);
minRDeltaT = max(minRDeltaT, clipRDeltaT);
maxRDeltaT = max(maxRDeltaT, clipRDeltaT);
}
if (pimpleDict.found("minDeltaT"))
if (pimpleDict.found("minDeltaT") || maxRDeltaT > rootVGreat)
{
rDeltaT.min(1/pimpleDict.lookup<scalar>("minDeltaT"));
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("minDeltaT");
rDeltaT.min(clipRDeltaT);
minRDeltaT = min(minRDeltaT, clipRDeltaT);
maxRDeltaT = min(maxRDeltaT, clipRDeltaT);
}
Info<< "Flow time scale min/max = "
<< 1/maxRDeltaT << ", " << 1/minRDeltaT << endl;
// Update the boundary values of the reciprocal time-step
rDeltaT.correctBoundaryConditions();
Info<< "Flow time scale min/max = "
<< gMin(1/rDeltaT.primitiveField())
<< ", " << gMax(1/rDeltaT.primitiveField()) << endl;
if (rDeltaTSmoothingCoeff < 1.0)
{
fvc::smooth(rDeltaT, rDeltaTSmoothingCoeff);

38
applications/modules/isothermalFluid/setRDeltaT.C
View File

@ -51,24 +51,12 @@ void Foam::solvers::isothermalFluid::setRDeltaT()
const volScalarField rDeltaT0("rDeltaT0", rDeltaT);
// Set the reciprocal time-step from the local Courant number
// and maximum and minimum time-steps
rDeltaT.ref() = fvc::surfaceSum(mag(phi))()()/((2*maxCo)*mesh.V()*rho());
if (pimpleDict.found("maxDeltaT"))
{
rDeltaT.max(1/pimpleDict.lookup<scalar>("maxDeltaT"));
}
if (pimpleDict.found("minDeltaT"))
{
rDeltaT.min(1/pimpleDict.lookup<scalar>("minDeltaT"));
}
// Set the reciprocal time-step from the local acoustic Courant number
if (pimple.transonic())
{
surfaceScalarField phid
(
"phid",
fvc::interpolate(psi)*fvc::flux(U)
);
surfaceScalarField phid("phid", fvc::interpolate(psi)*fvc::flux(U));
rDeltaT.ref() = max
(
@ -77,12 +65,26 @@ void Foam::solvers::isothermalFluid::setRDeltaT()
);
}
// Update the boundary values of the reciprocal time-step
rDeltaT.correctBoundaryConditions();
// Clip to user-defined maximum and minimum time-steps
scalar minRDeltaT = gMin(rDeltaT.primitiveField());
if (pimpleDict.found("maxDeltaT") || minRDeltaT < rootVSmall)
{
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("maxDeltaT");
rDeltaT.max(clipRDeltaT);
minRDeltaT = max(minRDeltaT, clipRDeltaT);
}
if (pimpleDict.found("minDeltaT"))
{
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("minDeltaT");
rDeltaT.min(clipRDeltaT);
minRDeltaT = min(minRDeltaT, clipRDeltaT);
}
Info<< "Flow time scale min/max = "
<< gMin(1/rDeltaT.primitiveField())
<< ", " << gMax(1/rDeltaT.primitiveField()) << endl;
<< gMin(1/rDeltaT.primitiveField()) << ", " << 1/minRDeltaT << endl;
// Update the boundary values of the reciprocal time-step
rDeltaT.correctBoundaryConditions();
if (rDeltaTSmoothingCoeff < 1.0)
{

19
applications/modules/multicomponentFluid/setRDeltaT.C
View File

@ -45,21 +45,26 @@ void Foam::solvers::multicomponentFluid::setRDeltaT()
const scalar maxCo(pimpleDict.lookup<scalar>("maxCo"));
// Set the reciprocal time-step from the local Courant number
// and maximum and minimum time-steps
rDeltaT.ref() =
fvc::surfaceSum(mag(phi))()()/((2*maxCo)*mesh.V()*rho());
if (pimpleDict.found("maxDeltaT"))
// Clip to user-defined maximum and minimum time-steps
scalar minRDeltaT = gMin(rDeltaT.primitiveField());
if (pimpleDict.found("maxDeltaT") || minRDeltaT < rootVSmall)
{
rDeltaT.max(1/pimpleDict.lookup<scalar>("maxDeltaT"));
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("maxDeltaT");
rDeltaT.max(clipRDeltaT);
minRDeltaT = max(minRDeltaT, clipRDeltaT);
}
if (pimpleDict.found("minDeltaT"))
{
rDeltaT.min(1/pimpleDict.lookup<scalar>("minDeltaT"));
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("minDeltaT");
rDeltaT.min(clipRDeltaT);
minRDeltaT = min(minRDeltaT, clipRDeltaT);
}
Info<< " Flow = "
<< 1/gMax(rDeltaT.primitiveField()) << ", "
<< 1/gMin(rDeltaT.primitiveField()) << endl;
Info<< "Flow time scale min/max = "
<< gMin(1/rDeltaT.primitiveField()) << ", " << 1/minRDeltaT << endl;
}
// Maximum change in cell temperature per iteration

19
applications/modules/multiphaseEuler/multiphaseEuler/setRDeltaT.C
View File

@ -53,23 +53,32 @@ void Foam::solvers::multiphaseEuler::setRDeltaT()
}
// Set the reciprocal time-step from the local Courant number
// and maximum and minimum time-steps
rDeltaT.ref() = fvc::surfaceSum(maxPhi)()()/((2*maxCo)*mesh.V());
if (pimpleDict.found("maxDeltaT"))
// Clip to user-defined maximum and minimum time-steps
scalar minRDeltaT = gMin(rDeltaT.primitiveField());
if (pimpleDict.found("maxDeltaT") || minRDeltaT < rootVSmall)
{
rDeltaT.max(1/pimpleDict.lookup<scalar>("maxDeltaT"));
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("maxDeltaT");
rDeltaT.max(clipRDeltaT);
minRDeltaT = max(minRDeltaT, clipRDeltaT);
}
if (pimpleDict.found("minDeltaT"))
{
rDeltaT.min(1/pimpleDict.lookup<scalar>("minDeltaT"));
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("minDeltaT");
rDeltaT.min(clipRDeltaT);
minRDeltaT = min(minRDeltaT, clipRDeltaT);
}
Info<< "Flow time scale min/max = "
<< gMin(1/rDeltaT.primitiveField()) << ", " << 1/minRDeltaT << endl;
// Update the boundary values of the reciprocal time-step
rDeltaT.correctBoundaryConditions();
fvc::smooth(rDeltaT, rDeltaTSmoothingCoeff);
Info<< "Flow time scale min/max = "
Info<< "Smoothed flow time scale min/max = "
<< gMin(1/rDeltaT.primitiveField())
<< ", " << gMax(1/rDeltaT.primitiveField()) << endl;

20
applications/modules/shockFluid/setRDeltaT.C
View File

@ -50,17 +50,31 @@ void Foam::solvers::shockFluid::setRDeltaT(const surfaceScalarField& amaxSf)
// Set the reciprocal time-step from the local Courant number
rDeltaT.ref() = fvc::surfaceSum(amaxSf)()()/((2*maxCo)*mesh.V());
if (pimpleDict.found("maxDeltaT"))
// Clip to user-defined maximum and minimum time-steps
scalar minRDeltaT = gMin(rDeltaT.primitiveField());
if (pimpleDict.found("maxDeltaT") || minRDeltaT < rootVSmall)
{
rDeltaT.max(1/pimpleDict.lookup<scalar>("maxDeltaT"));
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("maxDeltaT");
rDeltaT.max(clipRDeltaT);
minRDeltaT = max(minRDeltaT, clipRDeltaT);
}
if (pimpleDict.found("minDeltaT"))
{
const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("minDeltaT");
rDeltaT.min(clipRDeltaT);
minRDeltaT = min(minRDeltaT, clipRDeltaT);
}
Info<< "Flow time scale min/max = "
<< gMin(1/rDeltaT.primitiveField()) << ", " << 1/minRDeltaT << endl;
// Update the boundary values of the reciprocal time-step
rDeltaT.correctBoundaryConditions();
fvc::smooth(rDeltaT, rDeltaTSmoothingCoeff);
Info<< "Flow time scale min/max = "
Info<< "Smoothed flow time scale min/max = "
<< gMin(1/rDeltaT.primitiveField())
<< ", " << gMax(1/rDeltaT.primitiveField()) << endl;
}