JohnsonJacksonParticleThetaFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2014-2016 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 
73 (
74  const fvPatch& p,
75  const DimensionedField<scalar, volMesh>& iF,
76  const dictionary& dict
77 )
78 :
79  mixedFvPatchScalarField(p, iF),
80  restitutionCoefficient_
81  (
82  "restitutionCoefficient",
83  dimless,
84  dict.lookup("restitutionCoefficient")
85  ),
86  specularityCoefficient_
87  (
88  "specularityCoefficient",
89  dimless,
90  dict.lookup("specularityCoefficient")
91  )
92 {
93  if
94  (
95  (restitutionCoefficient_.value() < 0)
96  || (restitutionCoefficient_.value() > 1)
97  )
98  {
100  << "The restitution coefficient has to be between 0 and 1"
101  << abort(FatalError);
102  }
103 
104  if
105  (
106  (specularityCoefficient_.value() < 0)
107  || (specularityCoefficient_.value() > 1)
108  )
109  {
111  << "The specularity coefficient has to be between 0 and 1"
112  << abort(FatalError);
113  }
114 
116  (
117  scalarField("value", dict, p.size())
118  );
119 }
120 
121 
124 (
125  const JohnsonJacksonParticleThetaFvPatchScalarField& ptf
126 )
127 :
128  mixedFvPatchScalarField(ptf),
129  restitutionCoefficient_(ptf.restitutionCoefficient_),
130  specularityCoefficient_(ptf.specularityCoefficient_)
131 {}
132 
133 
136 (
137  const JohnsonJacksonParticleThetaFvPatchScalarField& ptf,
138  const DimensionedField<scalar, volMesh>& iF
139 )
140 :
141  mixedFvPatchScalarField(ptf, iF),
142  restitutionCoefficient_(ptf.restitutionCoefficient_),
143  specularityCoefficient_(ptf.specularityCoefficient_)
144 {}
145 
146 
147 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
148 
150 (
151  const fvPatchFieldMapper& m
152 )
153 {
154  mixedFvPatchScalarField::autoMap(m);
155 }
156 
157 
159 (
160  const fvPatchScalarField& ptf,
161  const labelList& addr
162 )
163 {
164  mixedFvPatchScalarField::rmap(ptf, addr);
165 }
166 
167 
169 {
170  if (updated())
171  {
172  return;
173  }
174 
175  // lookup the fluid model and the phase
176  const twoPhaseSystem& fluid = db().lookupObject<twoPhaseSystem>
177  (
178  "phaseProperties"
179  );
180 
181  const phaseModel& phased
182  (
183  fluid.phase1().name() == internalField().group()
184  ? fluid.phase1()
185  : fluid.phase2()
186  );
187 
188  // lookup all the fields on this patch
190  (
191  patch().lookupPatchField<volScalarField, scalar>
192  (
193  phased.volScalarField::name()
194  )
195  );
196 
197  const fvPatchVectorField& U
198  (
199  patch().lookupPatchField<volVectorField, vector>
200  (
201  IOobject::groupName("U", phased.name())
202  )
203  );
204 
205  const fvPatchScalarField& gs0
206  (
207  patch().lookupPatchField<volScalarField, scalar>
208  (
209  IOobject::groupName("gs0", phased.name())
210  )
211  );
212 
214  (
215  patch().lookupPatchField<volScalarField, scalar>
216  (
217  IOobject::groupName("kappa", phased.name())
218  )
219  );
220 
221  const scalarField Theta(patchInternalField());
222 
223  // lookup the packed volume fraction
225  (
226  "alphaMax",
227  dimless,
228  db()
229  .lookupObject<IOdictionary>
230  (
231  IOobject::groupName("turbulenceProperties", phased.name())
232  )
233  .subDict("RAS")
234  .subDict("kineticTheoryCoeffs")
235  .lookup("alphaMax")
236  );
237 
238  // calculate the reference value and the value fraction
239  if (restitutionCoefficient_.value() != 1.0)
240  {
241  this->refValue() =
242  (2.0/3.0)
243  *specularityCoefficient_.value()
244  *magSqr(U)
245  /(scalar(1) - sqr(restitutionCoefficient_.value()));
246 
247  this->refGrad() = 0.0;
248 
249  scalarField c
250  (
252  *alpha
253  *gs0
254  *(scalar(1) - sqr(restitutionCoefficient_.value()))
255  *sqrt(3.0*Theta)
256  /max(4.0*kappa*alphaMax.value(), SMALL)
257  );
258 
259  this->valueFraction() = c/(c + patch().deltaCoeffs());
260  }
261 
262  // for a restitution coefficient of 1, the boundary degenerates to a fixed
263  // gradient condition
264  else
265  {
266  this->refValue() = 0.0;
267 
268  this->refGrad() =
269  pos(alpha - SMALL)
271  *specularityCoefficient_.value()
272  *alpha
273  *gs0
274  *sqrt(3.0*Theta)
275  *magSqr(U)
276  /max(6.0*kappa*alphaMax.value(), SMALL);
277 
278  this->valueFraction() = 0.0;
279  }
280 
281  mixedFvPatchScalarField::updateCoeffs();
282 }
283 
284 
286 (
287  Ostream& os
288 ) const
289 {
291  os.writeKeyword("restitutionCoefficient")
292  << restitutionCoefficient_ << token::END_STATEMENT << nl;
293  os.writeKeyword("specularityCoefficient")
294  << specularityCoefficient_ << token::END_STATEMENT << nl;
295  writeEntry("value", os);
296 }
297 
298 
299 // ************************************************************************* //
U
Definition: pEqn.H:83
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)
const Type & value() const
Return const reference to value.
virtual void write(Ostream &) const
Write.
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Macros for easy insertion into run-time selection tables.
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
dimensionedScalar pos(const dimensionedScalar &ds)
JohnsonJacksonParticleThetaFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
stressControl lookup("compactNormalStress") >> compactNormalStress
fvPatchField< scalar > fvPatchScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static word groupName(Name name, const word &group)
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 > &)
static const char nl
Definition: Ostream.H:262
dimensionedScalar alphaMax(laminarTransport.lookup("alphaMax"))
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.
Namespace for OpenFOAM.
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:363