advectiveFvPatchField.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) 2011-2019 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 "advectiveFvPatchField.H"
28 #include "fvPatchFieldMapper.H"
29 #include "volFields.H"
30 #include "EulerDdtScheme.H"
31 #include "CrankNicolsonDdtScheme.H"
32 #include "backwardDdtScheme.H"
33 #include "localEulerDdtScheme.H"
34 
35 
36 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
37 
38 template<class Type>
40 (
41  const fvPatch& p,
43 )
44 :
46  phiName_("phi"),
47  rhoName_("rho"),
48  fieldInf_(Zero),
49  lInf_(-great)
50 {
51  this->refValue() = Zero;
52  this->refGrad() = Zero;
53  this->valueFraction() = 0.0;
54 }
55 
56 
57 template<class Type>
59 (
60  const advectiveFvPatchField& ptf,
61  const fvPatch& p,
63  const fvPatchFieldMapper& mapper
64 )
65 :
66  mixedFvPatchField<Type>(ptf, p, iF, mapper),
67  phiName_(ptf.phiName_),
68  rhoName_(ptf.rhoName_),
69  fieldInf_(ptf.fieldInf_),
70  lInf_(ptf.lInf_)
71 {}
72 
73 
74 template<class Type>
76 (
77  const fvPatch& p,
79  const dictionary& dict
80 )
81 :
83  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
84  rhoName_(dict.lookupOrDefault<word>("rho", "rho")),
85  fieldInf_(Zero),
86  lInf_(-great)
87 {
88  if (dict.found("value"))
89  {
91  (
92  Field<Type>("value", dict, p.size())
93  );
94  }
95  else
96  {
97  fvPatchField<Type>::operator=(this->patchInternalField());
98  }
99 
100  this->refValue() = *this;
101  this->refGrad() = Zero;
102  this->valueFraction() = 0.0;
103 
104  if (dict.readIfPresent("lInf", lInf_))
105  {
106  dict.lookup("fieldInf") >> fieldInf_;
107 
108  if (lInf_ < 0.0)
109  {
111  (
112  dict
113  ) << "unphysical lInf specified (lInf < 0)" << nl
114  << " on patch " << this->patch().name()
115  << " of field " << this->internalField().name()
116  << " in file " << this->internalField().objectPath()
117  << exit(FatalIOError);
118  }
119  }
120 }
121 
122 
123 template<class Type>
125 (
126  const advectiveFvPatchField& ptpsf
127 )
128 :
130  phiName_(ptpsf.phiName_),
131  rhoName_(ptpsf.rhoName_),
132  fieldInf_(ptpsf.fieldInf_),
133  lInf_(ptpsf.lInf_)
134 {}
135 
136 
137 template<class Type>
139 (
140  const advectiveFvPatchField& ptpsf,
142 )
143 :
144  mixedFvPatchField<Type>(ptpsf, iF),
145  phiName_(ptpsf.phiName_),
146  rhoName_(ptpsf.rhoName_),
147  fieldInf_(ptpsf.fieldInf_),
148  lInf_(ptpsf.lInf_)
149 {}
150 
151 
152 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
153 
154 template<class Type>
157 {
158  const surfaceScalarField& phi =
159  this->db().objectRegistry::template lookupObject<surfaceScalarField>
160  (phiName_);
161 
162  fvsPatchField<scalar> phip =
163  this->patch().template lookupPatchField<surfaceScalarField, scalar>
164  (
165  phiName_
166  );
167 
169  {
170  const fvPatchScalarField& rhop =
171  this->patch().template lookupPatchField<volScalarField, scalar>
172  (
173  rhoName_
174  );
175 
176  return phip/(rhop*this->patch().magSf());
177  }
178  else
179  {
180  return phip/this->patch().magSf();
181  }
182 }
183 
184 
185 template<class Type>
187 {
188  if (this->updated())
189  {
190  return;
191  }
192 
193  const fvMesh& mesh = this->internalField().mesh();
194 
195  word ddtScheme
196  (
197  mesh.ddtScheme(this->internalField().name())
198  );
199  scalar deltaT = this->db().time().deltaTValue();
200 
202  this->db().objectRegistry::template
203  lookupObject<GeometricField<Type, fvPatchField, volMesh>>
204  (
205  this->internalField().name()
206  );
207 
208  // Calculate the advection speed of the field wave
209  // If the wave is incoming set the speed to 0.
210  const scalarField w(Foam::max(advectionSpeed(), scalar(0)));
211 
212  // Calculate the field wave coefficient alpha (See notes)
213  const scalarField alpha(w*deltaT*this->patch().deltaCoeffs());
214 
215  label patchi = this->patch().index();
216 
217  // Non-reflecting outflow boundary
218  // If lInf_ defined setup relaxation to the value fieldInf_.
219  if (lInf_ > 0)
220  {
221  // Calculate the field relaxation coefficient k (See notes)
222  const scalarField k(w*deltaT/lInf_);
223 
224  if
225  (
228  )
229  {
230  this->refValue() =
231  (
232  field.oldTime().boundaryField()[patchi] + k*fieldInf_
233  )/(1.0 + k);
234 
235  this->valueFraction() = (1.0 + k)/(1.0 + alpha + k);
236  }
237  else if (ddtScheme == fv::backwardDdtScheme<scalar>::typeName)
238  {
239  this->refValue() =
240  (
241  2.0*field.oldTime().boundaryField()[patchi]
242  - 0.5*field.oldTime().oldTime().boundaryField()[patchi]
243  + k*fieldInf_
244  )/(1.5 + k);
245 
246  this->valueFraction() = (1.5 + k)/(1.5 + alpha + k);
247  }
248  else if
249  (
251  )
252  {
253  const volScalarField& rDeltaT =
254  fv::localEulerDdt::localRDeltaT(mesh);
255 
256  // Calculate the field wave coefficient alpha (See notes)
257  const scalarField alpha
258  (
259  w*this->patch().deltaCoeffs()/rDeltaT.boundaryField()[patchi]
260  );
261 
262  // Calculate the field relaxation coefficient k (See notes)
263  const scalarField k(w/(rDeltaT.boundaryField()[patchi]*lInf_));
264 
265  this->refValue() =
266  (
267  field.oldTime().boundaryField()[patchi] + k*fieldInf_
268  )/(1.0 + k);
269 
270  this->valueFraction() = (1.0 + k)/(1.0 + alpha + k);
271  }
272  else
273  {
275  << ddtScheme << nl
276  << " on patch " << this->patch().name()
277  << " of field " << this->internalField().name()
278  << " in file " << this->internalField().objectPath()
279  << exit(FatalError);
280  }
281  }
282  else
283  {
284  if
285  (
288  )
289  {
290  this->refValue() = field.oldTime().boundaryField()[patchi];
291 
292  this->valueFraction() = 1.0/(1.0 + alpha);
293  }
294  else if (ddtScheme == fv::backwardDdtScheme<scalar>::typeName)
295  {
296  this->refValue() =
297  (
298  2.0*field.oldTime().boundaryField()[patchi]
299  - 0.5*field.oldTime().oldTime().boundaryField()[patchi]
300  )/1.5;
301 
302  this->valueFraction() = 1.5/(1.5 + alpha);
303  }
304  else if
305  (
307  )
308  {
309  const volScalarField& rDeltaT =
310  fv::localEulerDdt::localRDeltaT(mesh);
311 
312  // Calculate the field wave coefficient alpha (See notes)
313  const scalarField alpha
314  (
315  w*this->patch().deltaCoeffs()/rDeltaT.boundaryField()[patchi]
316  );
317 
318  this->refValue() = field.oldTime().boundaryField()[patchi];
319 
320  this->valueFraction() = 1.0/(1.0 + alpha);
321  }
322  else
323  {
325  << ddtScheme
326  << "\n on patch " << this->patch().name()
327  << " of field " << this->internalField().name()
328  << " in file " << this->internalField().objectPath()
329  << exit(FatalError);
330  }
331  }
332 
334 }
335 
336 
337 template<class Type>
339 {
341 
342  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
343  writeEntryIfDifferent<word>(os, "rho", "rho", rhoName_);
344 
345  if (lInf_ > 0)
346  {
347  writeEntry(os, "fieldInf", fieldInf_);
348  writeEntry(os, "lInf", lInf_);
349  }
350 
351  writeEntry(os, "value", *this);
352 }
353 
354 
355 // ************************************************************************* //
Local time-step first-order Euler implicit/explicit ddt.
dictionary dict
Second-oder Crank-Nicolson implicit ddt using the current and previous time-step fields as well as th...
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:438
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
surfaceScalarField & phi
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
word phiName_
Name of the flux transporting the field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
scalar lInf_
Relaxation length-scale.
word rhoName_
Name of the density field used to normalise the mass flux.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Type fieldInf_
Field value of the far-field.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
This boundary condition provides a base class for &#39;mixed&#39; type boundary conditions, i.e. conditions that mix fixed value and patch-normal gradient conditions.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
ITstream & ddtScheme(const word &name) const
Definition: fvSchemes.C:360
label k
Boltzmann constant.
Macros for easy insertion into run-time selection tables.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
void write(Ostream &, const label, const dictionary &)
Write with dictionary lookup.
const dimensionSet & dimensions() const
Return dimensions.
dynamicFvMesh & mesh
Pre-declare SubField and related Field type.
Definition: Field.H:56
A class for handling words, derived from string.
Definition: word.H:59
advectiveFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
Foam::fvPatchFieldMapper.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
static const zero Zero
Definition: zero.H:97
Second-order backward-differencing ddt using the current and two previous time-step values...
virtual tmp< scalarField > advectionSpeed() const
Calculate and return the advection speed at the boundary.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:265
virtual void write(Ostream &) const
Write.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
const dimensionSet dimDensity
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
This boundary condition provides an advective outflow condition, based on solving DDt(W...
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
Basic first-order Euler implicit/explicit ddt using only the current and previous time-step values...
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:583
IOerror FatalIOError
const dimensionSet dimVelocity