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  dict
100  ) << "unphysical accommodationCoeff specified"
101  << "(0 < accommodationCoeff <= 1)" << endl
102  << exit(FatalIOError);
103  }
104 
105  if (dict.found("value"))
106  {
108  (
109  scalarField("value", dict, p.size())
110  );
111  }
112  else
113  {
114  fvPatchField<scalar>::operator=(patchInternalField());
115  }
116 
117  refValue() = *this;
118  refGrad() = 0.0;
119  valueFraction() = 0.0;
120 }
121 
122 
124 (
125  const smoluchowskiJumpTFvPatchScalarField& ptpsf,
126  const DimensionedField<scalar, volMesh>& iF
127 )
128 :
129  mixedFvPatchScalarField(ptpsf, iF),
130  accommodationCoeff_(ptpsf.accommodationCoeff_),
131  Twall_(ptpsf.Twall_),
132  gamma_(ptpsf.gamma_)
133 {}
134 
135 
136 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
137 
138 // Map from self
140 (
141  const fvPatchFieldMapper& m
142 )
143 {
144  mixedFvPatchScalarField::autoMap(m);
145 }
146 
147 
148 // Reverse-map the given fvPatchField onto this fvPatchField
150 (
151  const fvPatchField<scalar>& ptf,
152  const labelList& addr
153 )
154 {
156 }
157 
158 
159 // Update the coefficients associated with the patch field
161 {
162  if (updated())
163  {
164  return;
165  }
166 
167  const fvPatchScalarField& pmu =
168  patch().lookupPatchField<volScalarField, scalar>(muName_);
169  const fvPatchScalarField& prho =
170  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
171  const fvPatchField<scalar>& ppsi =
172  patch().lookupPatchField<volScalarField, scalar>(psiName_);
173  const fvPatchVectorField& pU =
174  patch().lookupPatchField<volVectorField, vector>(UName_);
175 
176  // Prandtl number reading consistent with rhoCentralFoam
177  const dictionary& thermophysicalProperties =
178  db().lookupObject<IOdictionary>(basicThermo::dictName);
179 
181  (
182  "Pr",
183  dimless,
184  thermophysicalProperties.subDict("mixture").subDict("transport")
185  .lookup("Pr")
186  );
187 
188  Field<scalar> C2
189  (
190  pmu/prho
192  *2.0*gamma_/Pr.value()/(gamma_ + 1.0)
193  *(2.0 - accommodationCoeff_)/accommodationCoeff_
194  );
195 
196  Field<scalar> aCoeff(prho.snGrad() - prho/C2);
197  Field<scalar> KEbyRho(0.5*magSqr(pU));
198 
199  valueFraction() = (1.0/(1.0 + patch().deltaCoeffs()*C2));
200  refValue() = Twall_;
201  refGrad() = 0.0;
202 
203  mixedFvPatchScalarField::updateCoeffs();
204 }
205 
206 
207 // Write
209 {
211 
212  writeEntryIfDifferent<word>(os, "U", "U", UName_);
213  writeEntryIfDifferent<word>(os, "rho", "rho", rhoName_);
214  writeEntryIfDifferent<word>(os, "psi", "thermo:psi", psiName_);
215  writeEntryIfDifferent<word>(os, "mu", "thermo:mu", muName_);
216 
217  os.writeKeyword("accommodationCoeff")
218  << accommodationCoeff_ << token::END_STATEMENT << nl;
219  Twall_.writeEntry("Twall", os);
220  os.writeKeyword("gamma")
221  << gamma_ << token::END_STATEMENT << nl;
222  writeEntry("value", os);
223 }
224 
225 
226 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
227 
228 namespace Foam
229 {
231  (
233  smoluchowskiJumpTFvPatchScalarField
234  );
235 }
236 
237 
238 // ************************************************************************* //
dimensionedScalar Pr("Pr", dimless, laminarTransport)
fvPatchField< vector > fvPatchVectorField
virtual void write(Ostream &) const
Write.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
Macros for easy insertion into run-time selection tables.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
stressControl lookup("compactNormalStress") >> compactNormalStress
fvPatchField< scalar > fvPatchScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
List< label > labelList
A List of labels.
Definition: labelList.H:56
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
dimensioned< scalar > magSqr(const dimensioned< Type > &)
static const char nl
Definition: Ostream.H:262
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
smoluchowskiJumpTFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:396
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
const scalar piByTwo(0.5 *pi)
Namespace for OpenFOAM.
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:363
IOerror FatalIOError