setRDeltaT.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2022-2024 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "VoFSolver.H"
27 #include "fvcSurfaceIntegrate.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 void Foam::solvers::VoFSolver::setRDeltaT()
32 {
33  volScalarField& rDeltaT = trDeltaT.ref();
34 
35  const dictionary& pimpleDict = pimple.dict();
36 
37  const scalar maxCo
38  (
39  pimpleDict.lookupOrDefault<scalar>("maxCo", 0.9)
40  );
41 
42  const volScalarField rDeltaT0("rDeltaT0", rDeltaT);
43 
44  // Set the reciprocal time-step from the local Courant number
45  rDeltaT.internalFieldRef() =
46  fvc::surfaceSum(mag(phi))()()/((2*maxCo)*mesh.V());
47 
48  // Clip to user-defined maximum and minimum time-steps
49  scalar minRDeltaT = gMin(rDeltaT.primitiveField());
50  if (pimpleDict.found("maxDeltaT") || minRDeltaT < rootVSmall)
51  {
52  const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("maxDeltaT");
53  rDeltaT.max(clipRDeltaT);
54  minRDeltaT = max(minRDeltaT, clipRDeltaT);
55  }
56  if (pimpleDict.found("minDeltaT"))
57  {
58  const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("minDeltaT");
59  rDeltaT.min(clipRDeltaT);
60  minRDeltaT = min(minRDeltaT, clipRDeltaT);
61  }
62 
63  Info<< "Flow time scale min/max = "
64  << gMin(1/rDeltaT.primitiveField()) << ", " << 1/minRDeltaT << endl;
65 
66  setInterfaceRDeltaT(rDeltaT);
67 
68  // Limit rate of change of time scale
69  // - reduce as much as required
70  // - only increase at a fraction of old time scale
71  if
72  (
73  pimpleDict.found("rDeltaTDampingCoeff")
75  )
76  {
77  // Damping coefficient (1-0)
78  const scalar rDeltaTDampingCoeff
79  (
80  pimpleDict.lookup<scalar>("rDeltaTDampingCoeff")
81  );
82 
83  rDeltaT = max
84  (
85  rDeltaT,
86  (scalar(1) - rDeltaTDampingCoeff)*rDeltaT0
87  );
88 
89  Info<< "Damped flow time scale min/max = "
90  << gMin(1/rDeltaT.primitiveField())
91  << ", " << gMax(1/rDeltaT.primitiveField()) << endl;
92  }
93 }
94 
95 
96 // ************************************************************************* //
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:28
virtual label startTimeIndex() const
Return start time index.
Definition: Time.C:785
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
virtual const dictionary & dict() const
Return the solution dictionary.
pimpleNoLoopControl pimple
PIMPLE inner-loop controls.
Definition: solver.H:107
const Time & runTime
Time.
Definition: solver.H:104
const fvMesh & mesh
Region mesh.
Definition: solver.H:101
const surfaceScalarField & phi
Reference to the mass-flux field.
Definition: VoFSolver.H:212
tmp< volScalarField > trDeltaT
Optional LTS reciprocal time-step field.
Definition: VoFSolver.H:141
virtual void setInterfaceRDeltaT(volScalarField &rDeltaT)=0
Adjust the rDeltaT in the vicinity of the interface.
scalar maxCo
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField.
tmp< VolField< Type > > surfaceSum(const SurfaceField< Type > &ssf)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
messageStream Info
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
Type gMin(const FieldField< Field, Type > &f)
Type gMax(const FieldField< Field, Type > &f)