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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2015 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
31#include "fvPatchFieldMapper.H"
32#include "volFields.H"
33#include "basicThermo.H"
35
36// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
37
39(
40 const fvPatch& p,
41 const DimensionedField<scalar, volMesh>& iF
42)
43:
44 mixedFvPatchScalarField(p, iF),
45 UName_("U"),
46 rhoName_("rho"),
47 psiName_("thermo:psi"),
48 muName_("thermo:mu"),
49 accommodationCoeff_(1.0),
50 Twall_(p.size(), Zero),
51 gamma_(1.4)
52{
53 refValue() = 0.0;
54 refGrad() = 0.0;
55 valueFraction() = 0.0;
56}
57
58
60(
61 const smoluchowskiJumpTFvPatchScalarField& ptf,
62 const fvPatch& p,
63 const DimensionedField<scalar, volMesh>& iF,
64 const fvPatchFieldMapper& mapper
65)
66:
67 mixedFvPatchScalarField(ptf, p, iF, mapper),
68 UName_(ptf.UName_),
69 rhoName_(ptf.rhoName_),
70 psiName_(ptf.psiName_),
71 muName_(ptf.muName_),
72 accommodationCoeff_(ptf.accommodationCoeff_),
73 Twall_(ptf.Twall_),
74 gamma_(ptf.gamma_)
75{}
76
77
79(
80 const fvPatch& p,
81 const DimensionedField<scalar, volMesh>& iF,
82 const dictionary& dict
83)
84:
85 mixedFvPatchScalarField(p, iF),
86 UName_(dict.getOrDefault<word>("U", "U")),
87 rhoName_(dict.getOrDefault<word>("rho", "rho")),
88 psiName_(dict.getOrDefault<word>("psi", "thermo:psi")),
89 muName_(dict.getOrDefault<word>("mu", "thermo:mu")),
90 accommodationCoeff_(dict.get<scalar>("accommodationCoeff")),
91 Twall_("Twall", dict, p.size()),
92 gamma_(dict.getOrDefault<scalar>("gamma", 1.4))
93{
94 if
95 (
96 mag(accommodationCoeff_) < SMALL
97 || mag(accommodationCoeff_) > 2.0
98 )
99 {
101 << "unphysical accommodationCoeff specified"
102 << "(0 < accommodationCoeff <= 1)" << endl
103 << exit(FatalIOError);
104 }
105
106 if (dict.found("value"))
107 {
108 fvPatchField<scalar>::operator=
109 (
110 scalarField("value", dict, p.size())
111 );
112 }
113 else
114 {
115 fvPatchField<scalar>::operator=(patchInternalField());
116 }
117
118 refValue() = *this;
119 refGrad() = 0.0;
120 valueFraction() = 0.0;
121}
122
123
125(
126 const smoluchowskiJumpTFvPatchScalarField& ptpsf,
127 const DimensionedField<scalar, volMesh>& iF
128)
129:
130 mixedFvPatchScalarField(ptpsf, iF),
131 accommodationCoeff_(ptpsf.accommodationCoeff_),
132 Twall_(ptpsf.Twall_),
133 gamma_(ptpsf.gamma_)
134{}
135
136
137// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
138
139// Map from self
141(
142 const fvPatchFieldMapper& m
143)
144{
145 mixedFvPatchScalarField::autoMap(m);
146}
147
148
149// Reverse-map the given fvPatchField onto this fvPatchField
151(
152 const fvPatchField<scalar>& ptf,
153 const labelList& addr
154)
155{
157}
158
159
160// Update the coefficients associated with the patch field
162{
163 if (updated())
164 {
165 return;
166 }
167
168 const fvPatchScalarField& pmu =
169 patch().lookupPatchField<volScalarField, scalar>(muName_);
170 const fvPatchScalarField& prho =
171 patch().lookupPatchField<volScalarField, scalar>(rhoName_);
172 const fvPatchField<scalar>& ppsi =
173 patch().lookupPatchField<volScalarField, scalar>(psiName_);
174 const fvPatchVectorField& pU =
175 patch().lookupPatchField<volVectorField, vector>(UName_);
176
177 // Prandtl number reading consistent with rhoCentralFoam
178 const dictionary& thermophysicalProperties =
179 db().lookupObject<IOdictionary>(basicThermo::dictName);
180
182 (
183 "Pr",
184 dimless,
185 thermophysicalProperties.subDict("mixture").subDict("transport")
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 os.writeEntryIfDifferent<word>("U", "U", UName_);
213 os.writeEntryIfDifferent<word>("rho", "rho", rhoName_);
214 os.writeEntryIfDifferent<word>("psi", "thermo:psi", psiName_);
215 os.writeEntryIfDifferent<word>("mu", "thermo:mu", muName_);
216
217 os.writeEntry("accommodationCoeff", accommodationCoeff_);
218 Twall_.writeEntry("Twall", os);
219 os.writeEntry("gamma", gamma_);
220 writeEntry("value", os);
221}
222
223
224// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
225
226namespace Foam
227{
229 (
231 smoluchowskiJumpTFvPatchScalarField
232 );
233}
234
235
236// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:239
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Definition: Ostream.H:251
static const word dictName
Definition: basicThermo.H:256
const Type & value() const
Return const reference to value.
virtual bool write()
Write the output fields.
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Smoluchowski temperature jump boundary condition.
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.
virtual void rmap(const fvPatchField< scalar > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
volScalarField & p
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
Definition: fvPatchField.H:676
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
constexpr scalar piByTwo(0.5 *M_PI)
const std::string patch
OpenFOAM patch number as a std::string.
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:83
const dimensionSet dimless
Dimensionless.
List< label > labelList
A List of labels.
Definition: List.H:66
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
fvPatchField< vector > fvPatchVectorField
fvPatchField< scalar > fvPatchScalarField
dictionary dict
dimensionedScalar Pr("Pr", dimless, laminarTransport)