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 "multiphaseEuler.H"
27 #include "fvcSmooth.H"
28 #include "fvcSurfaceIntegrate.H"
29 
30 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
31 
32 void Foam::solvers::multiphaseEuler::setRDeltaT()
33 {
34  volScalarField& rDeltaT = trDeltaT.ref();
35 
36  const dictionary& pimpleDict = pimple.dict();
37 
38  const scalar maxCo
39  (
40  pimpleDict.lookupOrDefault<scalar>("maxCo", 0.2)
41  );
42 
43  const scalar rDeltaTSmoothingCoeff
44  (
45  pimpleDict.lookupOrDefault<scalar>("rDeltaTSmoothingCoeff", 0.02)
46  );
47 
48  surfaceScalarField maxPhi("maxPhi", phi);
49 
50  forAll(movingPhases, movingPhasei)
51  {
52  maxPhi = max(maxPhi, mag(movingPhases[movingPhasei].phi()));
53  }
54 
55  // Set the reciprocal time-step from the local Courant number
56  rDeltaT.internalFieldRef() =
57  fvc::surfaceSum(maxPhi)()()/((2*maxCo)*mesh.V());
58 
59  // Clip to user-defined maximum and minimum time-steps
60  scalar minRDeltaT = gMin(rDeltaT.primitiveField());
61  if (pimpleDict.found("maxDeltaT") || minRDeltaT < rootVSmall)
62  {
63  const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("maxDeltaT");
64  rDeltaT.max(clipRDeltaT);
65  minRDeltaT = max(minRDeltaT, clipRDeltaT);
66  }
67  if (pimpleDict.found("minDeltaT"))
68  {
69  const scalar clipRDeltaT = 1/pimpleDict.lookup<scalar>("minDeltaT");
70  rDeltaT.min(clipRDeltaT);
71  minRDeltaT = min(minRDeltaT, clipRDeltaT);
72  }
73 
74  Info<< "Flow time scale min/max = "
75  << gMin(1/rDeltaT.primitiveField()) << ", " << 1/minRDeltaT << endl;
76 
77  // Update the boundary values of the reciprocal time-step
78  rDeltaT.correctBoundaryConditions();
79 
80  fvc::smooth(rDeltaT, rDeltaTSmoothingCoeff);
81 
82  Info<< "Smoothed flow time scale min/max = "
83  << gMin(1/rDeltaT.primitiveField())
84  << ", " << gMax(1/rDeltaT.primitiveField()) << endl;
85 
86  if (faceMomentum)
87  {
88  trDeltaTf.ref() = fvc::interpolate(rDeltaT);
89  }
90 }
91 
92 
93 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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 fvMesh & mesh
Region mesh.
Definition: solver.H:101
const surfaceScalarField & phi
Reference to the mass-flux field.
tmp< volScalarField > trDeltaT
Optional LTS reciprocal time-step field.
const phaseSystem::phaseModelPartialList & movingPhases
Reference to the moving phases.
Switch faceMomentum
Cell/face momentum equation switch.
tmp< surfaceScalarField > trDeltaTf
Optional LTS reciprocal face time-step field.
scalar maxCo
Provides functions smooth spread and sweep which use the FvFaceCellWave algorithm to smooth and redis...
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField.
tmp< VolField< Type > > surfaceSum(const SurfaceField< Type > &ssf)
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
void smooth(volScalarField &field, const scalar coeff)
Definition: fvcSmooth.C:35
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
SurfaceField< scalar > surfaceScalarField
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)