JohnsonJacksonParticleSlipFvPatchVectorField.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 "JohnsonJacksonParticleSlipFvPatchVectorField.H"
28 #include "twoPhaseSystem.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
35  (
37  JohnsonJacksonParticleSlipFvPatchVectorField
38  );
39 }
40 
41 
42 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
43 
46 (
47  const fvPatch& p,
48  const DimensionedField<vector, volMesh>& iF
49 )
50 :
51  partialSlipFvPatchVectorField(p, iF),
52  specularityCoefficient_("specularityCoefficient", dimless, 0)
53 {}
54 
55 
58 (
59  const JohnsonJacksonParticleSlipFvPatchVectorField& ptf,
60  const fvPatch& p,
61  const DimensionedField<vector, volMesh>& iF,
62  const fvPatchFieldMapper& mapper
63 )
64 :
65  partialSlipFvPatchVectorField(ptf, p, iF, mapper),
66  specularityCoefficient_(ptf.specularityCoefficient_)
67 {}
68 
69 
72 (
73  const fvPatch& p,
74  const DimensionedField<vector, volMesh>& iF,
75  const dictionary& dict
76 )
77 :
78  partialSlipFvPatchVectorField(p, iF),
79  specularityCoefficient_
80  (
81  "specularityCoefficient",
82  dimless,
83  dict.lookup("specularityCoefficient")
84  )
85 {
86  if
87  (
88  (specularityCoefficient_.value() < 0)
89  || (specularityCoefficient_.value() > 1)
90  )
91  {
93  << "The specularity coefficient has to be between 0 and 1"
94  << abort(FatalError);
95  }
96 
98  (
99  vectorField("value", dict, p.size())
100  );
101 }
102 
103 
106 (
107  const JohnsonJacksonParticleSlipFvPatchVectorField& ptf
108 )
109 :
110  partialSlipFvPatchVectorField(ptf),
111  specularityCoefficient_(ptf.specularityCoefficient_)
112 {}
113 
114 
117 (
118  const JohnsonJacksonParticleSlipFvPatchVectorField& ptf,
119  const DimensionedField<vector, volMesh>& iF
120 )
121 :
122  partialSlipFvPatchVectorField(ptf, iF),
123  specularityCoefficient_(ptf.specularityCoefficient_)
124 {}
125 
126 
127 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
128 
130 (
131  const fvPatchFieldMapper& m
132 )
133 {
134  partialSlipFvPatchVectorField::autoMap(m);
135 }
136 
137 
139 (
140  const fvPatchVectorField& ptf,
141  const labelList& addr
142 )
143 {
144  partialSlipFvPatchVectorField::rmap(ptf, addr);
145 }
146 
147 
149 {
150  if (updated())
151  {
152  return;
153  }
154 
155  // lookup the fluid model and the phase
156  const twoPhaseSystem& fluid = db().lookupObject<twoPhaseSystem>
157  (
158  "phaseProperties"
159  );
160 
161  const phaseModel& phased
162  (
163  fluid.phase1().name() == internalField().group()
164  ? fluid.phase1()
165  : fluid.phase2()
166  );
167 
168  // lookup all the fields on this patch
170  (
171  patch().lookupPatchField<volScalarField, scalar>
172  (
173  phased.volScalarField::name()
174  )
175  );
176 
177  const fvPatchScalarField& gs0
178  (
179  patch().lookupPatchField<volScalarField, scalar>
180  (
181  IOobject::groupName("gs0", phased.name())
182  )
183  );
184 
185  const scalarField nu
186  (
187  patch().lookupPatchField<volScalarField, scalar>
188  (
189  IOobject::groupName("nut", phased.name())
190  )
191  );
192 
193  const scalarField nuFric
194  (
195  patch().lookupPatchField<volScalarField, scalar>
196  (
197  IOobject::groupName("nuFric", phased.name())
198  )
199  );
200 
201  word ThetaName(IOobject::groupName("Theta", phased.name()));
202 
203  const fvPatchScalarField& Theta
204  (
205  db().foundObject<volScalarField>(ThetaName)
206  ? patch().lookupPatchField<volScalarField, scalar>(ThetaName)
207  : alpha
208  );
209 
210  // lookup the packed volume fraction
212  (
213  "alphaMax",
214  dimless,
215  db()
216  .lookupObject<IOdictionary>
217  (
218  IOobject::groupName("turbulenceProperties", phased.name())
219  )
220  .subDict("RAS")
221  .subDict("kineticTheoryCoeffs")
222  .lookup("alphaMax")
223  );
224 
225  // calculate the slip value fraction
226  scalarField c
227  (
229  *alpha
230  *gs0
231  *specularityCoefficient_.value()
232  *sqrt(3*Theta)
233  /max(6*(nu - nuFric)*alphaMax.value(), small)
234  );
235 
236  this->valueFraction() = c/(c + patch().deltaCoeffs());
237 
238  partialSlipFvPatchVectorField::updateCoeffs();
239 }
240 
241 
243 (
244  Ostream& os
245 ) const
246 {
248  os.writeKeyword("specularityCoefficient")
249  << specularityCoefficient_ << token::END_STATEMENT << nl;
250  writeEntry("value", os);
251 }
252 
253 
254 // ************************************************************************* //
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void write(Ostream &) const
Write.
fvPatchField< vector > fvPatchVectorField
JohnsonJacksonParticleSlipFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
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 > &)
virtual void rmap(const fvPatchVectorField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:362
Macros for easy insertion into run-time selection tables.
stressControl lookup("compactNormalStress") >> compactNormalStress
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
static const char nl
Definition: Ostream.H:265
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.
Field< vector > vectorField
Specialisation of Field<T> for vector.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
volScalarField & nu
Namespace for OpenFOAM.