FreeFOAM The Cross-Platform CFD Toolkit
turbulentMixingLengthFrequencyInletFvPatchScalarField.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) 2006-2011 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 
30 #include <finiteVolume/volFields.H>
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace incompressible
38 {
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
42 turbulentMixingLengthFrequencyInletFvPatchScalarField::
43 turbulentMixingLengthFrequencyInletFvPatchScalarField
44 (
45  const fvPatch& p,
47 )
48 :
49  inletOutletFvPatchScalarField(p, iF),
50  mixingLength_(0.0),
51  phiName_("undefined-phi"),
52  kName_("undefined-k")
53 {
54  this->refValue() = 0.0;
55  this->refGrad() = 0.0;
56  this->valueFraction() = 0.0;
57 }
58 
59 turbulentMixingLengthFrequencyInletFvPatchScalarField::
60 turbulentMixingLengthFrequencyInletFvPatchScalarField
61 (
63  const fvPatch& p,
65  const fvPatchFieldMapper& mapper
66 )
67 :
68  inletOutletFvPatchScalarField(ptf, p, iF, mapper),
69  mixingLength_(ptf.mixingLength_),
70  phiName_(ptf.phiName_),
71  kName_(ptf.kName_)
72 {}
73 
74 turbulentMixingLengthFrequencyInletFvPatchScalarField::
75 turbulentMixingLengthFrequencyInletFvPatchScalarField
76 (
77  const fvPatch& p,
79  const dictionary& dict
80 )
81 :
82  inletOutletFvPatchScalarField(p, iF),
83  mixingLength_(readScalar(dict.lookup("mixingLength"))),
84  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
85  kName_(dict.lookupOrDefault<word>("k", "k"))
86 {
87  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
88 
89  this->refValue() = 0.0;
90  this->refGrad() = 0.0;
91  this->valueFraction() = 0.0;
92 }
93 
94 turbulentMixingLengthFrequencyInletFvPatchScalarField::
95 turbulentMixingLengthFrequencyInletFvPatchScalarField
96 (
98 )
99 :
100  inletOutletFvPatchScalarField(ptf),
101  mixingLength_(ptf.mixingLength_),
102  phiName_(ptf.phiName_),
103  kName_(ptf.kName_)
104 {}
105 
106 turbulentMixingLengthFrequencyInletFvPatchScalarField::
107 turbulentMixingLengthFrequencyInletFvPatchScalarField
108 (
111 )
112 :
113  inletOutletFvPatchScalarField(ptf, iF),
114  mixingLength_(ptf.mixingLength_),
115  phiName_(ptf.phiName_),
116  kName_(ptf.kName_)
117 {}
118 
119 
120 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
121 
122 void turbulentMixingLengthFrequencyInletFvPatchScalarField::updateCoeffs()
123 {
124  if (updated())
125  {
126  return;
127  }
128 
129  // Lookup Cmu corresponding to the turbulence model selected
130  const RASModel& rasModel = db().lookupObject<RASModel>("RASProperties");
131 
132  const scalar Cmu =
133  rasModel.coeffDict().lookupOrDefault<scalar>("Cmu", 0.09);
134 
135  const scalar Cmu25 = pow(Cmu, 0.25);
136 
137  const fvPatchScalarField& kp =
138  patch().lookupPatchField<volScalarField, scalar>(kName_);
139 
140  const fvsPatchScalarField& phip =
141  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
142 
143  this->refValue() = sqrt(kp)/(Cmu25*mixingLength_);
144  this->valueFraction() = 1.0 - pos(phip);
145 
146  inletOutletFvPatchScalarField::updateCoeffs();
147 }
148 
149 
150 void turbulentMixingLengthFrequencyInletFvPatchScalarField::write
151 (
152  Ostream& os
153 ) const
154 {
156  os.writeKeyword("mixingLength")
157  << mixingLength_ << token::END_STATEMENT << nl;
158  os.writeKeyword("phi") << phiName_ << token::END_STATEMENT << nl;
159  os.writeKeyword("k") << kName_ << token::END_STATEMENT << nl;
160  writeEntry("value", os);
161 }
162 
163 
164 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
165 
167 (
170 );
171 
172 
173 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
174 
175 } // End namespace incompressible
176 } // End namespace Foam
177 
178 // ************************ vim: set sw=4 sts=4 et: ************************ //