setRDeltaT.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2013-2016 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29{
30 volScalarField& rDeltaT = trDeltaT.ref();
31
32 const dictionary& pimpleDict = pimple.dict();
33
34 // Maximum flow Courant number
35 scalar maxCo(pimpleDict.get<scalar>("maxCo"));
36
37 // Maximum time scale
38 scalar maxDeltaT(pimpleDict.getOrDefault<scalar>("maxDeltaT", GREAT));
39
40 // Smoothing parameter (0-1) when smoothing iterations > 0
42 (
43 pimpleDict.getOrDefault<scalar>("rDeltaTSmoothingCoeff", 0.1)
44 );
45
46 // Damping coefficient (1-0)
48 (
49 pimpleDict.getOrDefault<scalar>("rDeltaTDampingCoeff", 1.0)
50 );
51
52 // Maximum change in cell temperature per iteration
53 // (relative to previous value)
54 scalar alphaTemp(pimpleDict.getOrDefault<scalar>("alphaTemp", 0.05));
55
56 // Maximum change in cell concentration per iteration
57 // (relative to reference value)
58 scalar alphaY(pimpleDict.getOrDefault<scalar>("alphaY", 1.0));
59
60 Info<< "Time scales min/max:" << endl;
61
62 // Cache old reciprocal time scale field
63 volScalarField rDeltaT0("rDeltaT0", rDeltaT);
64
65 // Flow time scale
66 {
67 rDeltaT.ref() =
68 (
69 fvc::surfaceSum(mag(phi))()()
70 /((2*maxCo)*mesh.V()*rho())
71 );
72
73 // Limit the largest time scale
74 rDeltaT.max(1/maxDeltaT);
75
76 Info<< " Flow = "
77 << 1/gMax(rDeltaT.primitiveField()) << ", "
78 << 1/gMin(rDeltaT.primitiveField()) << endl;
79 }
80
81 // Heat release rate time scale
82 if (alphaTemp < 1)
83 {
84 volScalarField::Internal rDeltaTT
85 (
87 );
88
89 Info<< " Temperature = "
90 << 1/(gMax(rDeltaTT.field()) + VSMALL) << ", "
91 << 1/(gMin(rDeltaTT.field()) + VSMALL) << endl;
92
93 rDeltaT.ref() = max(rDeltaT(), rDeltaTT);
94 }
95
96 // Reaction rate time scale
97 if (alphaY < 1)
98 {
99 dictionary Yref(pimpleDict.subDict("Yref"));
100
101 volScalarField::Internal rDeltaTY
102 (
103 IOobject
104 (
105 "rDeltaTY",
106 runTime.timeName(),
107 mesh
108 ),
109 mesh,
110 dimensionedScalar(rDeltaT.dimensions(), Zero)
111 );
112
113 bool foundY = false;
114
115 forAll(Y, i)
116 {
117 if (i != inertIndex && composition.active(i))
118 {
119 volScalarField& Yi = Y[i];
120
121 if (Yref.found(Yi.name()))
122 {
123 foundY = true;
124 const scalar Yrefi = Yref.get<scalar>(Yi.name());
125
126 rDeltaTY.field() = max
127 (
128 mag
129 (
130 reaction->R(Yi)().source()
131 /((Yrefi*alphaY)*(rho*mesh.V()))
132 ),
133 rDeltaTY
134 );
135 }
136 }
137 }
138
139 if (foundY)
140 {
141 Info<< " Composition = "
142 << 1/(gMax(rDeltaTY.field()) + VSMALL) << ", "
143 << 1/(gMin(rDeltaTY.field()) + VSMALL) << endl;
144
145 rDeltaT.ref() = max(rDeltaT(), rDeltaTY);
146 }
147 else
148 {
149 IOWarningIn(args.executable().c_str(), Yref)
150 << "Cannot find any active species in Yref " << Yref
151 << endl;
152 }
153 }
154
155 // Update tho boundary values of the reciprocal time-step
156 rDeltaT.correctBoundaryConditions();
157
158 // Spatially smooth the time scale field
160 {
161 fvc::smooth(rDeltaT, rDeltaTSmoothingCoeff);
162 }
163
164 // Limit rate of change of time scale
165 // - reduce as much as required
166 // - only increase at a fraction of old time scale
167 if
168 (
170 && runTime.timeIndex() > runTime.startTimeIndex() + 1
171 )
172 {
173 rDeltaT = max
174 (
175 rDeltaT,
176 (scalar(1) - rDeltaTDampingCoeff)*rDeltaT0
177 );
178 }
179
180 // Update tho boundary values of the reciprocal time-step
181 rDeltaT.correctBoundaryConditions();
182
183 Info<< " Overall = "
184 << 1/gMax(rDeltaT.primitiveField())
185 << ", " << 1/gMin(rDeltaT.primitiveField()) << endl;
186}
187
188
189// ************************************************************************* //
Y[inertIndex] max(0.0)
surfaceScalarField & phi
pimpleControl & pimple
const word & executable() const noexcept
Name of executable without the path.
Definition: argListI.H:51
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
basicSpecieMixture & composition
PtrList< volScalarField > & Y
const volScalarField & T
scalar rDeltaTDampingCoeff(pimpleDict.getOrDefault< scalar >("rDeltaTDampingCoeff", 1.0))
scalar alphaY(pimpleDict.getOrDefault< scalar >("alphaY", 1.0))
const dictionary & pimpleDict
Definition: setRDeltaT.H:32
scalar rDeltaTSmoothingCoeff(pimpleDict.getOrDefault< scalar >("rDeltaTSmoothingCoeff", 0.1))
scalar alphaTemp(pimpleDict.getOrDefault< scalar >("alphaTemp", 0.05))
volScalarField rDeltaT0("rDeltaT0", rDeltaT)
dynamicFvMesh & mesh
engineTime & runTime
tmp< volScalarField > trDeltaT
Definition: createRDeltaT.H:3
scalar maxCo
label inertIndex
CombustionModel< rhoReactionThermo > & reaction
#define IOWarningIn(functionName, ios)
Report an IO warning using Foam::Warning.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Type gMin(const FieldField< Field, Type > &f)
Type gMax(const FieldField< Field, Type > &f)
Foam::argList args(argc, argv)
scalar Qdot
Definition: solveChemistry.H:2
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333