FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 
29 #include <finiteVolume/volFields.H>
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
40 (
41  const fvPatch& p,
43 )
44 :
45  mixedFvPatchScalarField(p, iF),
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 
56 (
58  const fvPatch& p,
60  const fvPatchFieldMapper& mapper
61 )
62 :
63  mixedFvPatchScalarField(ptf, p, iF, mapper),
64  accommodationCoeff_(ptf.accommodationCoeff_),
65  Twall_(ptf.Twall_),
66  gamma_(ptf.gamma_)
67 {}
68 
69 
71 (
72  const fvPatch& p,
74  const dictionary& dict
75 )
76 :
77  mixedFvPatchScalarField(p, iF),
78  accommodationCoeff_(readScalar(dict.lookup("accommodationCoeff"))),
79  Twall_("Twall", dict, p.size()),
80  gamma_(dict.lookupOrDefault<scalar>("gamma", 1.4))
81 {
82  if
83  (
84  mag(accommodationCoeff_) < SMALL
85  || mag(accommodationCoeff_) > 2.0
86  )
87  {
89  (
90  "smoluchowskiJumpTFvPatchScalarField::"
91  "smoluchowskiJumpTFvPatchScalarField"
92  "("
93  " const fvPatch&,"
94  " const DimensionedField<scalar, volMesh>&,"
95  " const dictionary&"
96  ")",
97  dict
98  ) << "unphysical accommodationCoeff specified"
99  << "(0 < accommodationCoeff <= 1)" << endl
100  << exit(FatalError);
101  }
102 
103  if (dict.found("value"))
104  {
106  (
107  scalarField("value", dict, p.size())
108  );
109  }
110  else
111  {
112  fvPatchField<scalar>::operator=(patchInternalField());
113  }
114 
115  refValue() = *this;
116  refGrad() = 0.0;
117  valueFraction() = 0.0;
118 }
119 
120 
122 (
125 )
126 :
127  mixedFvPatchScalarField(ptpsf, iF),
128  accommodationCoeff_(ptpsf.accommodationCoeff_),
129  Twall_(ptpsf.Twall_),
130  gamma_(ptpsf.gamma_)
131 {}
132 
133 
134 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
135 
136 // Map from self
138 (
139  const fvPatchFieldMapper& m
140 )
141 {
142  mixedFvPatchScalarField::autoMap(m);
143 }
144 
145 
146 // Reverse-map the given fvPatchField onto this fvPatchField
148 (
149  const fvPatchField<scalar>& ptf,
150  const labelList& addr
151 )
152 {
154 }
155 
156 
157 // Update the coefficients associated with the patch field
159 {
160  if (updated())
161  {
162  return;
163  }
164 
165  const fvPatchScalarField& pmu =
166  patch().lookupPatchField<volScalarField, scalar>("mu");
167  const fvPatchScalarField& prho =
168  patch().lookupPatchField<volScalarField, scalar>("rho");
169  const fvPatchField<scalar>& ppsi =
170  patch().lookupPatchField<volScalarField, scalar>("psi");
171  const fvPatchVectorField& pU =
172  patch().lookupPatchField<volVectorField, vector>("U");
173 
174  // Prandtl number reading consistent with rhoCentralFoam
175  const dictionary& thermophysicalProperties =
176  db().lookupObject<IOdictionary>("thermophysicalProperties");
178  if (thermophysicalProperties.found("Pr"))
179  {
180  Pr = thermophysicalProperties.lookup("Pr");
181  }
182 
183  Field<scalar> C2 = pmu/prho
184  *sqrt(ppsi*mathematicalConstant::pi/2.0)
185  *2.0*gamma_/Pr.value()/(gamma_ + 1.0)
186  *(2.0 - accommodationCoeff_)/accommodationCoeff_;
187 
188  Field<scalar> aCoeff = prho.snGrad() - prho/C2;
189  Field<scalar> KEbyRho = 0.5*magSqr(pU);
190 
191  valueFraction() = (1.0/(1.0 + patch().deltaCoeffs()*C2));
192  refValue() = Twall_;
193  refGrad() = 0.0;
194 
196 }
197 
198 
199 // Write
201 {
203  os.writeKeyword("accommodationCoeff")
204  << accommodationCoeff_ << token::END_STATEMENT << nl;
205  Twall_.writeEntry("Twall", os);
206  os.writeKeyword("gamma")
207  << gamma_ << token::END_STATEMENT << nl;
208  writeEntry("value", os);
209 }
210 
211 
212 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
213 
215 
216 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
217 
218 } // End namespace Foam
219 
220 // ************************ vim: set sw=4 sts=4 et: ************************ //