smoluchowskiJumpTFvPatchScalarField.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) 2011-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 
28 #include "fvPatchFieldMapper.H"
29 #include "volFields.H"
30 #include "basicThermo.H"
31 #include "mathematicalConstants.H"
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
36 (
37  const fvPatch& p,
38  const DimensionedField<scalar, volMesh>& iF
39 )
40 :
41  mixedFvPatchScalarField(p, iF),
42  UName_("U"),
43  rhoName_("rho"),
44  psiName_("thermo:psi"),
45  muName_("thermo:mu"),
46  accommodationCoeff_(1.0),
47  Twall_(p.size(), 0.0),
48  gamma_(1.4)
49 {
50  refValue() = 0.0;
51  refGrad() = 0.0;
52  valueFraction() = 0.0;
53 }
54 
55 
57 (
58  const smoluchowskiJumpTFvPatchScalarField& ptf,
59  const fvPatch& p,
60  const DimensionedField<scalar, volMesh>& iF,
61  const fvPatchFieldMapper& mapper
62 )
63 :
64  mixedFvPatchScalarField(ptf, p, iF, mapper),
65  UName_(ptf.UName_),
66  rhoName_(ptf.rhoName_),
67  psiName_(ptf.psiName_),
68  muName_(ptf.muName_),
69  accommodationCoeff_(ptf.accommodationCoeff_),
70  Twall_(ptf.Twall_),
71  gamma_(ptf.gamma_)
72 {}
73 
74 
76 (
77  const fvPatch& p,
78  const DimensionedField<scalar, volMesh>& iF,
79  const dictionary& dict
80 )
81 :
82  mixedFvPatchScalarField(p, iF),
83  UName_(dict.lookupOrDefault<word>("U", "U")),
84  rhoName_(dict.lookupOrDefault<word>("rho", "rho")),
85  psiName_(dict.lookupOrDefault<word>("psi", "thermo:psi")),
86  muName_(dict.lookupOrDefault<word>("mu", "thermo:mu")),
87  accommodationCoeff_(readScalar(dict.lookup("accommodationCoeff"))),
88  Twall_("Twall", dict, p.size()),
89  gamma_(dict.lookupOrDefault<scalar>("gamma", 1.4))
90 {
91  if
92  (
93  mag(accommodationCoeff_) < SMALL
94  || mag(accommodationCoeff_) > 2.0
95  )
96  {
98  (
99  "smoluchowskiJumpTFvPatchScalarField::"
100  "smoluchowskiJumpTFvPatchScalarField"
101  "("
102  " const fvPatch&,"
103  " const DimensionedField<scalar, volMesh>&,"
104  " const dictionary&"
105  ")",
106  dict
107  ) << "unphysical accommodationCoeff specified"
108  << "(0 < accommodationCoeff <= 1)" << endl
109  << exit(FatalIOError);
110  }
111 
112  if (dict.found("value"))
113  {
115  (
116  scalarField("value", dict, p.size())
117  );
118  }
119  else
120  {
121  fvPatchField<scalar>::operator=(patchInternalField());
122  }
123 
124  refValue() = *this;
125  refGrad() = 0.0;
126  valueFraction() = 0.0;
127 }
128 
129 
131 (
132  const smoluchowskiJumpTFvPatchScalarField& ptpsf,
133  const DimensionedField<scalar, volMesh>& iF
134 )
135 :
136  mixedFvPatchScalarField(ptpsf, iF),
137  accommodationCoeff_(ptpsf.accommodationCoeff_),
138  Twall_(ptpsf.Twall_),
139  gamma_(ptpsf.gamma_)
140 {}
141 
142 
143 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
144 
145 // Map from self
147 (
148  const fvPatchFieldMapper& m
149 )
150 {
151  mixedFvPatchScalarField::autoMap(m);
152 }
153 
154 
155 // Reverse-map the given fvPatchField onto this fvPatchField
157 (
158  const fvPatchField<scalar>& ptf,
159  const labelList& addr
160 )
161 {
163 }
164 
165 
166 // Update the coefficients associated with the patch field
168 {
169  if (updated())
170  {
171  return;
172  }
173 
174  const fvPatchScalarField& pmu =
175  patch().lookupPatchField<volScalarField, scalar>(muName_);
176  const fvPatchScalarField& prho =
177  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
178  const fvPatchField<scalar>& ppsi =
179  patch().lookupPatchField<volScalarField, scalar>(psiName_);
180  const fvPatchVectorField& pU =
181  patch().lookupPatchField<volVectorField, vector>(UName_);
182 
183  // Prandtl number reading consistent with rhoCentralFoam
184  const dictionary& thermophysicalProperties =
185  db().lookupObject<IOdictionary>(basicThermo::dictName);
186 
188  (
189  "Pr",
190  dimless,
191  thermophysicalProperties.subDict("mixture").subDict("transport")
192  .lookup("Pr")
193  );
194 
195  Field<scalar> C2
196  (
197  pmu/prho
199  *2.0*gamma_/Pr.value()/(gamma_ + 1.0)
200  *(2.0 - accommodationCoeff_)/accommodationCoeff_
201  );
202 
203  Field<scalar> aCoeff(prho.snGrad() - prho/C2);
204  Field<scalar> KEbyRho(0.5*magSqr(pU));
205 
206  valueFraction() = (1.0/(1.0 + patch().deltaCoeffs()*C2));
207  refValue() = Twall_;
208  refGrad() = 0.0;
209 
210  mixedFvPatchScalarField::updateCoeffs();
211 }
212 
213 
214 // Write
216 {
218 
219  writeEntryIfDifferent<word>(os, "U", "U", UName_);
220  writeEntryIfDifferent<word>(os, "rho", "rho", rhoName_);
221  writeEntryIfDifferent<word>(os, "psi", "thermo:psi", psiName_);
222  writeEntryIfDifferent<word>(os, "mu", "thermo:mu", muName_);
223 
224  os.writeKeyword("accommodationCoeff")
225  << accommodationCoeff_ << token::END_STATEMENT << nl;
226  Twall_.writeEntry("Twall", os);
227  os.writeKeyword("gamma")
228  << gamma_ << token::END_STATEMENT << nl;
229  writeEntry("value", os);
230 }
231 
232 
233 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
234 
235 namespace Foam
236 {
238  (
240  smoluchowskiJumpTFvPatchScalarField
241  );
242 }
243 
244 
245 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:387
dimensionedScalar Pr("Pr", dimless, laminarTransport)
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void write(Ostream &) const
Write.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:354
Namespace for OpenFOAM.
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
stressControl lookup("compactNormalStress") >> compactNormalStress
IOerror FatalIOError
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
virtual void rmap(const fvPatchField< scalar > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Macros for easy insertion into run-time selection tables.
fvPatchField< scalar > fvPatchScalarField
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
List< label > labelList
A List of labels.
Definition: labelList.H:56
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
fvPatchField< vector > fvPatchVectorField
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
smoluchowskiJumpTFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
const scalar piByTwo(0.5 *pi)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)