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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2014-2019 OpenFOAM Foundation
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
30 #include "twoPhaseSystem.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
37  (
39  JohnsonJacksonParticleThetaFvPatchScalarField
40  );
41 }
42 
43 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
44 
47 (
48  const fvPatch& p,
50 )
51 :
52  mixedFvPatchScalarField(p, iF),
53  restitutionCoefficient_("restitutionCoefficient", dimless, Zero),
54  specularityCoefficient_("specularityCoefficient", dimless, Zero)
55 {}
56 
57 
60 (
62  const fvPatch& p,
64  const fvPatchFieldMapper& mapper
65 )
66 :
67  mixedFvPatchScalarField(ptf, p, iF, mapper),
68  restitutionCoefficient_(ptf.restitutionCoefficient_),
69  specularityCoefficient_(ptf.specularityCoefficient_)
70 {
71 }
72 
73 
76 (
77  const fvPatch& p,
79  const dictionary& dict
80 )
81 :
82  mixedFvPatchScalarField(p, iF),
83  restitutionCoefficient_("restitutionCoefficient", dimless, dict),
84  specularityCoefficient_("specularityCoefficient", dimless, dict)
85 {
86  if
87  (
88  (restitutionCoefficient_.value() < 0)
89  || (restitutionCoefficient_.value() > 1)
90  )
91  {
93  << "The restitution coefficient has to be between 0 and 1"
94  << abort(FatalError);
95  }
96 
97  if
98  (
99  (specularityCoefficient_.value() < 0)
100  || (specularityCoefficient_.value() > 1)
101  )
102  {
104  << "The specularity coefficient has to be between 0 and 1"
105  << abort(FatalError);
106  }
107 
108  fvPatchScalarField::operator=
109  (
110  scalarField("value", dict, p.size())
111  );
112 }
113 
114 
117 (
119 )
120 :
121  mixedFvPatchScalarField(ptf),
122  restitutionCoefficient_(ptf.restitutionCoefficient_),
123  specularityCoefficient_(ptf.specularityCoefficient_)
124 {}
125 
126 
129 (
132 )
133 :
134  mixedFvPatchScalarField(ptf, iF),
135  restitutionCoefficient_(ptf.restitutionCoefficient_),
136  specularityCoefficient_(ptf.specularityCoefficient_)
137 {}
138 
139 
140 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
141 
143 (
144  const fvPatchFieldMapper& m
145 )
146 {
147  mixedFvPatchScalarField::autoMap(m);
148 }
149 
150 
152 (
153  const fvPatchScalarField& ptf,
154  const labelList& addr
155 )
156 {
157  mixedFvPatchScalarField::rmap(ptf, addr);
158 }
159 
160 
162 {
163  if (updated())
164  {
165  return;
166  }
167 
168  // lookup the fluid model and the phase
169  const twoPhaseSystem& fluid = db().lookupObject<twoPhaseSystem>
170  (
171  "phaseProperties"
172  );
173 
174  const phaseModel& phased
175  (
176  fluid.phase1().name() == internalField().group()
177  ? fluid.phase1()
178  : fluid.phase2()
179  );
180 
181  // lookup all the fields on this patch
183  (
184  patch().lookupPatchField<volScalarField, scalar>
185  (
186  phased.volScalarField::name()
187  )
188  );
189 
190  const fvPatchVectorField& U
191  (
192  patch().lookupPatchField<volVectorField, vector>
193  (
194  IOobject::groupName("U", phased.name())
195  )
196  );
197 
198  const fvPatchScalarField& gs0
199  (
200  patch().lookupPatchField<volScalarField, scalar>
201  (
202  IOobject::groupName("gs0", phased.name())
203  )
204  );
205 
207  (
208  patch().lookupPatchField<volScalarField, scalar>
209  (
210  IOobject::groupName("kappa", phased.name())
211  )
212  );
213 
214  const scalarField Theta(patchInternalField());
215 
216  // lookup the packed volume fraction
218  (
219  "alphaMax",
220  dimless,
221  db()
222  .lookupObject<IOdictionary>
223  (
224  IOobject::groupName("turbulenceProperties", phased.name())
225  )
226  .subDict("RAS")
227  .subDict("kineticTheoryCoeffs")
228  );
229 
230  // calculate the reference value and the value fraction
231  if (restitutionCoefficient_.value() != 1.0)
232  {
233  this->refValue() =
234  (2.0/3.0)
235  *specularityCoefficient_.value()
236  *magSqr(U)
237  /(scalar(1) - sqr(restitutionCoefficient_.value()));
238 
239  this->refGrad() = 0.0;
240 
241  scalarField c
242  (
244  *alpha
245  *gs0
246  *(scalar(1) - sqr(restitutionCoefficient_.value()))
247  *sqrt(3*Theta)
248  /max(4*kappa*alphaMax.value(), SMALL)
249  );
250 
251  this->valueFraction() = c/(c + patch().deltaCoeffs());
252  }
253 
254  // for a restitution coefficient of 1, the boundary degenerates to a fixed
255  // gradient condition
256  else
257  {
258  this->refValue() = 0.0;
259 
260  this->refGrad() =
261  pos0(alpha - SMALL)
263  *specularityCoefficient_.value()
264  *alpha
265  *gs0
266  *sqrt(3*Theta)
267  *magSqr(U)
268  /max(6*kappa*alphaMax.value(), SMALL);
269 
270  this->valueFraction() = 0;
271  }
272 
273  mixedFvPatchScalarField::updateCoeffs();
274 }
275 
276 
278 (
279  Ostream& os
280 ) const
281 {
283  os.writeEntry("restitutionCoefficient", restitutionCoefficient_);
284  os.writeEntry("specularityCoefficient", specularityCoefficient_);
285  writeEntry("value", os);
286 }
287 
288 
289 // ************************************************************************* //
Foam::fvPatchField< scalar >
Foam::fvPatchScalarField
fvPatchField< scalar > fvPatchScalarField
Definition: fvPatchFieldsFwd.H:40
Foam::fvPatchField< scalar >::write
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:384
Foam::phaseModel
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:54
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::JohnsonJacksonParticleThetaFvPatchScalarField::rmap
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: JohnsonJacksonParticleThetaFvPatchScalarField.C:152
Foam::twoPhaseSystem
Class which solves the volume fraction equations for two phases.
Definition: twoPhaseSystem.H:53
JohnsonJacksonParticleThetaFvPatchScalarField.H
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::JohnsonJacksonParticleThetaFvPatchScalarField::JohnsonJacksonParticleThetaFvPatchScalarField
JohnsonJacksonParticleThetaFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Definition: JohnsonJacksonParticleThetaFvPatchScalarField.C:47
Foam::constant::atomic::alpha
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Definition: readThermalProperties.H:212
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
Foam::pos0
dimensionedScalar pos0(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:188
fluid
twoPhaseSystem & fluid
Definition: setRegionFluidFields.H:3
Foam::constant::electromagnetic::kappa
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Foam::JohnsonJacksonParticleThetaFvPatchScalarField
Robin condition for the particulate granular temperature.
Definition: JohnsonJacksonParticleThetaFvPatchScalarField.H:66
Foam::magSqr
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
Foam::Field< scalar >
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::JohnsonJacksonParticleThetaFvPatchScalarField::write
virtual void write(Ostream &) const
Write.
Definition: JohnsonJacksonParticleThetaFvPatchScalarField.C:278
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
alphaMax
dimensionedScalar alphaMax("alphaMax", dimless/dimTime, laminarTransport)
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
os
OBJstream os(runTime.globalPath()/outputName)
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::dimensioned< scalar >
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::JohnsonJacksonParticleThetaFvPatchScalarField::autoMap
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: JohnsonJacksonParticleThetaFvPatchScalarField.C:143
U
U
Definition: pEqn.H:72
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::constant::mathematical::pi
constexpr scalar pi(M_PI)
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::List< label >
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
Foam::constant::universal::c
const dimensionedScalar c
Speed of light in a vacuum.
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
Foam::IOobject::groupName
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::makePatchTypeField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Foam::JohnsonJacksonParticleThetaFvPatchScalarField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients.
Definition: JohnsonJacksonParticleThetaFvPatchScalarField.C:161
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
Foam::dimless
const dimensionSet dimless
Dimensionless.
Definition: dimensionSets.C:189