JohnsonJacksonParticleThetaFvPatchScalarField.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) 2014-2018 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 "JohnsonJacksonParticleThetaFvPatchScalarField.H"
28 #include "twoPhaseSystem.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
35  (
37  JohnsonJacksonParticleThetaFvPatchScalarField
38  );
39 }
40 
41 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
42 
45 (
46  const fvPatch& p,
47  const DimensionedField<scalar, volMesh>& iF
48 )
49 :
50  mixedFvPatchScalarField(p, iF),
51  restitutionCoefficient_("restitutionCoefficient", dimless, 0),
52  specularityCoefficient_("specularityCoefficient", dimless, 0)
53 {}
54 
55 
58 (
59  const JohnsonJacksonParticleThetaFvPatchScalarField& ptf,
60  const fvPatch& p,
61  const DimensionedField<scalar, volMesh>& iF,
62  const fvPatchFieldMapper& mapper
63 )
64 :
65  mixedFvPatchScalarField(ptf, p, iF, mapper),
66  restitutionCoefficient_(ptf.restitutionCoefficient_),
67  specularityCoefficient_(ptf.specularityCoefficient_)
68 {
69 }
70 
71 
74 (
75  const fvPatch& p,
76  const DimensionedField<scalar, volMesh>& iF,
77  const dictionary& dict
78 )
79 :
80  mixedFvPatchScalarField(p, iF),
81  restitutionCoefficient_
82  (
83  "restitutionCoefficient",
84  dimless,
85  dict.lookup("restitutionCoefficient")
86  ),
87  specularityCoefficient_
88  (
89  "specularityCoefficient",
90  dimless,
91  dict.lookup("specularityCoefficient")
92  )
93 {
94  if
95  (
96  (restitutionCoefficient_.value() < 0)
97  || (restitutionCoefficient_.value() > 1)
98  )
99  {
101  << "The restitution coefficient has to be between 0 and 1"
102  << abort(FatalError);
103  }
104 
105  if
106  (
107  (specularityCoefficient_.value() < 0)
108  || (specularityCoefficient_.value() > 1)
109  )
110  {
112  << "The specularity coefficient has to be between 0 and 1"
113  << abort(FatalError);
114  }
115 
117  (
118  scalarField("value", dict, p.size())
119  );
120 }
121 
122 
125 (
126  const JohnsonJacksonParticleThetaFvPatchScalarField& ptf
127 )
128 :
129  mixedFvPatchScalarField(ptf),
130  restitutionCoefficient_(ptf.restitutionCoefficient_),
131  specularityCoefficient_(ptf.specularityCoefficient_)
132 {}
133 
134 
137 (
138  const JohnsonJacksonParticleThetaFvPatchScalarField& ptf,
139  const DimensionedField<scalar, volMesh>& iF
140 )
141 :
142  mixedFvPatchScalarField(ptf, iF),
143  restitutionCoefficient_(ptf.restitutionCoefficient_),
144  specularityCoefficient_(ptf.specularityCoefficient_)
145 {}
146 
147 
148 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
149 
151 (
152  const fvPatchFieldMapper& m
153 )
154 {
155  mixedFvPatchScalarField::autoMap(m);
156 }
157 
158 
160 (
161  const fvPatchScalarField& ptf,
162  const labelList& addr
163 )
164 {
165  mixedFvPatchScalarField::rmap(ptf, addr);
166 }
167 
168 
170 {
171  if (updated())
172  {
173  return;
174  }
175 
176  // lookup the fluid model and the phase
177  const twoPhaseSystem& fluid = db().lookupObject<twoPhaseSystem>
178  (
179  "phaseProperties"
180  );
181 
182  const phaseModel& phased
183  (
184  fluid.phase1().name() == internalField().group()
185  ? fluid.phase1()
186  : fluid.phase2()
187  );
188 
189  // lookup all the fields on this patch
191  (
192  patch().lookupPatchField<volScalarField, scalar>
193  (
194  phased.volScalarField::name()
195  )
196  );
197 
198  const fvPatchVectorField& U
199  (
200  patch().lookupPatchField<volVectorField, vector>
201  (
202  IOobject::groupName("U", phased.name())
203  )
204  );
205 
206  const fvPatchScalarField& gs0
207  (
208  patch().lookupPatchField<volScalarField, scalar>
209  (
210  IOobject::groupName("gs0", phased.name())
211  )
212  );
213 
215  (
216  patch().lookupPatchField<volScalarField, scalar>
217  (
218  IOobject::groupName("kappa", phased.name())
219  )
220  );
221 
222  const scalarField Theta(patchInternalField());
223 
224  // lookup the packed volume fraction
226  (
227  "alphaMax",
228  dimless,
229  db()
230  .lookupObject<IOdictionary>
231  (
232  IOobject::groupName("turbulenceProperties", phased.name())
233  )
234  .subDict("RAS")
235  .subDict("kineticTheoryCoeffs")
236  .lookup("alphaMax")
237  );
238 
239  // calculate the reference value and the value fraction
240  if (restitutionCoefficient_.value() != 1.0)
241  {
242  this->refValue() =
243  (2.0/3.0)
244  *specularityCoefficient_.value()
245  *magSqr(U)
246  /(scalar(1) - sqr(restitutionCoefficient_.value()));
247 
248  this->refGrad() = 0.0;
249 
250  scalarField c
251  (
253  *alpha
254  *gs0
255  *(scalar(1) - sqr(restitutionCoefficient_.value()))
256  *sqrt(3*Theta)
257  /max(4*kappa*alphaMax.value(), small)
258  );
259 
260  this->valueFraction() = c/(c + patch().deltaCoeffs());
261  }
262 
263  // for a restitution coefficient of 1, the boundary degenerates to a fixed
264  // gradient condition
265  else
266  {
267  this->refValue() = 0.0;
268 
269  this->refGrad() =
270  pos0(alpha - small)
272  *specularityCoefficient_.value()
273  *alpha
274  *gs0
275  *sqrt(3*Theta)
276  *magSqr(U)
277  /max(6*kappa*alphaMax.value(), small);
278 
279  this->valueFraction() = 0;
280  }
281 
282  mixedFvPatchScalarField::updateCoeffs();
283 }
284 
285 
287 (
288  Ostream& os
289 ) const
290 {
292  os.writeKeyword("restitutionCoefficient")
293  << restitutionCoefficient_ << token::END_STATEMENT << nl;
294  os.writeKeyword("specularityCoefficient")
295  << specularityCoefficient_ << token::END_STATEMENT << nl;
296  writeEntry("value", os);
297 }
298 
299 
300 // ************************************************************************* //
fvPatchField< vector > fvPatchVectorField
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:362
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Macros for easy insertion into run-time selection tables.
JohnsonJacksonParticleThetaFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
stressControl lookup("compactNormalStress") >> compactNormalStress
virtual void write(Ostream &) const
Write.
fvPatchField< scalar > fvPatchScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static word groupName(Name name, const word &group)
const Type & value() const
Return const reference to value.
List< label > labelList
A List of labels.
Definition: labelList.H:56
errorManip< error > abort(error &err)
Definition: errorManip.H:131
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensionedScalar pos0(const dimensionedScalar &ds)
static const char nl
Definition: Ostream.H:265
dimensionedScalar alphaMax(laminarTransport.lookup("alphaMax"))
U
Definition: pEqn.H:72
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionedScalar c
Speed of light in a vacuum.
virtual void updateCoeffs()
Update the coefficients.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Namespace for OpenFOAM.