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-2015 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_(p.size()),
52  specularityCoefficient_(p.size())
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  (
102  "Foam::JohnsonJacksonParticleThetaFvPatchScalarField::"
103  "JohnsonJacksonParticleThetaFvPatchScalarField"
104  "("
105  "const fvPatch& p,"
106  "const DimensionedField<scalar, volMesh>& iF,"
107  "const dictionary& dict"
108  ")"
109  ) << "The restitution coefficient has to be between 0 and 1"
110  << abort(FatalError);
111  }
112 
113  if
114  (
115  (specularityCoefficient_.value() < 0)
116  || (specularityCoefficient_.value() > 1)
117  )
118  {
120  (
121  "Foam::JohnsonJacksonParticleThetaFvPatchScalarField::"
122  "JohnsonJacksonParticleThetaFvPatchScalarField"
123  "("
124  "const fvPatch& p,"
125  "const DimensionedField<scalar, volMesh>& iF,"
126  "const dictionary& dict"
127  ")"
128  ) << "The specularity coefficient has to be between 0 and 1"
129  << abort(FatalError);
130  }
131 
133  (
134  scalarField("value", dict, p.size())
135  );
136 }
137 
138 
141 (
142  const JohnsonJacksonParticleThetaFvPatchScalarField& ptf
143 )
144 :
145  mixedFvPatchScalarField(ptf),
146  restitutionCoefficient_(ptf.restitutionCoefficient_),
147  specularityCoefficient_(ptf.specularityCoefficient_)
148 {}
149 
150 
153 (
154  const JohnsonJacksonParticleThetaFvPatchScalarField& ptf,
155  const DimensionedField<scalar, volMesh>& iF
156 )
157 :
158  mixedFvPatchScalarField(ptf, iF),
159  restitutionCoefficient_(ptf.restitutionCoefficient_),
160  specularityCoefficient_(ptf.specularityCoefficient_)
161 {}
162 
163 
164 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
165 
167 (
168  const fvPatchFieldMapper& m
169 )
170 {
171  mixedFvPatchScalarField::autoMap(m);
172 }
173 
174 
176 (
177  const fvPatchScalarField& ptf,
178  const labelList& addr
179 )
180 {
181  mixedFvPatchScalarField::rmap(ptf, addr);
182 }
183 
184 
186 {
187  if (updated())
188  {
189  return;
190  }
191 
192  // lookup the fluid model and the phase
193  const twoPhaseSystem& fluid = db().lookupObject<twoPhaseSystem>
194  (
195  "phaseProperties"
196  );
197 
198  const phaseModel& phased
199  (
200  fluid.phase1().name() == dimensionedInternalField().group()
201  ? fluid.phase1()
202  : fluid.phase2()
203  );
204 
205  // lookup all the fields on this patch
207  (
208  patch().lookupPatchField<volScalarField, scalar>
209  (
210  phased.volScalarField::name()
211  )
212  );
213 
214  const fvPatchVectorField& U
215  (
216  patch().lookupPatchField<volVectorField, vector>
217  (
218  IOobject::groupName("U", phased.name())
219  )
220  );
221 
222  const fvPatchScalarField& gs0
223  (
224  patch().lookupPatchField<volScalarField, scalar>
225  (
226  IOobject::groupName("gs0", phased.name())
227  )
228  );
229 
231  (
232  patch().lookupPatchField<volScalarField, scalar>
233  (
234  IOobject::groupName("kappa", phased.name())
235  )
236  );
237 
238  const scalarField Theta(patchInternalField());
239 
240  // lookup the packed volume fraction
242  (
243  "alphaMax",
244  dimless,
245  db()
246  .lookupObject<IOdictionary>
247  (
248  IOobject::groupName("turbulenceProperties", phased.name())
249  )
250  .subDict("RAS")
251  .subDict("kineticTheoryCoeffs")
252  .lookup("alphaMax")
253  );
254 
255  // calculate the reference value and the value fraction
256  if (restitutionCoefficient_.value() != 1.0)
257  {
258  this->refValue() =
259  (2.0/3.0)
260  *specularityCoefficient_.value()
261  *magSqr(U)
262  /(scalar(1) - sqr(restitutionCoefficient_.value()));
263 
264  this->refGrad() = 0.0;
265 
266  scalarField c
267  (
269  *alpha
270  *gs0
271  *(scalar(1) - sqr(restitutionCoefficient_.value()))
272  *sqrt(3.0*Theta)
273  /max(4.0*kappa*alphaMax.value(), SMALL)
274  );
275 
276  this->valueFraction() = c/(c + patch().deltaCoeffs());
277  }
278 
279  // for a restitution coefficient of 1, the boundary degenerates to a fixed
280  // gradient condition
281  else
282  {
283  this->refValue() = 0.0;
284 
285  this->refGrad() =
286  pos(alpha - SMALL)
288  *specularityCoefficient_.value()
289  *alpha
290  *gs0
291  *sqrt(3.0*Theta)
292  *magSqr(U)
293  /max(6.0*kappa*alphaMax.value(), SMALL);
294 
295  this->valueFraction() = 0.0;
296  }
297 
298  mixedFvPatchScalarField::updateCoeffs();
299 }
300 
301 
303 (
304  Ostream& os
305 ) const
306 {
308  os.writeKeyword("restitutionCoefficient")
309  << restitutionCoefficient_ << token::END_STATEMENT << nl;
310  os.writeKeyword("specularityCoefficient")
311  << specularityCoefficient_ << token::END_STATEMENT << nl;
312  writeEntry("value", os);
313 }
314 
315 
316 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
JohnsonJacksonParticleThetaFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:354
Namespace for OpenFOAM.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static word groupName(Name name, const word &group)
static const char nl
Definition: Ostream.H:260
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
stressControl lookup("compactNormalStress") >> compactNormalStress
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Macros for easy insertion into run-time selection tables.
fvPatchField< scalar > fvPatchScalarField
errorManip< error > abort(error &err)
Definition: errorManip.H:131
dimensionedScalar alphaMax(laminarTransport.lookup("alphaMax"))
virtual void write(Ostream &) const
Write.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
virtual void updateCoeffs()
Update the coefficients.
rDeltaT dimensionedInternalField()
error FatalError
List< label > labelList
A List of labels.
Definition: labelList.H:56
const dimensionedScalar c
Speed of light in a vacuum.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
dimensionedScalar pos(const dimensionedScalar &ds)
fvPatchField< vector > fvPatchVectorField
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)
U
Definition: pEqn.H:82
const Type & value() const
Return const reference to value.