FreeFOAM The Cross-Platform CFD Toolkit
timeVaryingUniformTotalPressureFvPatchScalarField.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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
36 (
37  const fvPatch& p,
39 )
40 :
41  fixedValueFvPatchScalarField(p, iF),
42  UName_("U"),
43  phiName_("phi"),
44  rhoName_("none"),
45  psiName_("none"),
46  gamma_(0.0),
47  p0_(0.0),
48  timeSeries_()
49 {}
50 
51 
54 (
55  const fvPatch& p,
57  const dictionary& dict
58 )
59 :
60  fixedValueFvPatchScalarField(p, iF),
61  UName_(dict.lookupOrDefault<word>("U", "U")),
62  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
63  rhoName_(dict.lookupOrDefault<word>("rho", "none")),
64  psiName_(dict.lookupOrDefault<word>("psi", "none")),
65  gamma_(readScalar(dict.lookup("gamma"))),
66  p0_(readScalar(dict.lookup("p0"))),
67  timeSeries_(dict)
68 {
69  if (dict.found("value"))
70  {
72  (
73  scalarField("value", dict, p.size())
74  );
75  }
76  else
77  {
79  }
80 }
81 
82 
85 (
87  const fvPatch& p,
89  const fvPatchFieldMapper& mapper
90 )
91 :
92  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
93  UName_(ptf.UName_),
94  phiName_(ptf.phiName_),
95  rhoName_(ptf.rhoName_),
96  psiName_(ptf.psiName_),
97  gamma_(ptf.gamma_),
98  p0_(ptf.p0_),
99  timeSeries_(ptf.timeSeries_)
100 {}
101 
102 
105 (
107 )
108 :
109  fixedValueFvPatchScalarField(tppsf),
110  UName_(tppsf.UName_),
111  phiName_(tppsf.phiName_),
112  rhoName_(tppsf.rhoName_),
113  psiName_(tppsf.psiName_),
114  gamma_(tppsf.gamma_),
115  p0_(tppsf.p0_),
116  timeSeries_(tppsf.timeSeries_)
117 {}
118 
119 
122 (
125 )
126 :
127  fixedValueFvPatchScalarField(tppsf, iF),
128  UName_(tppsf.UName_),
129  phiName_(tppsf.phiName_),
130  rhoName_(tppsf.rhoName_),
131  psiName_(tppsf.psiName_),
132  gamma_(tppsf.gamma_),
133  p0_(tppsf.p0_),
134  timeSeries_(tppsf.timeSeries_)
135 {}
136 
137 
138 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
139 
141 (
142  const vectorField& Up
143 )
144 {
145  if (updated())
146  {
147  return;
148  }
149 
150  p0_ = timeSeries_(this->db().time().timeOutputValue());
151 
152  const fvsPatchField<scalar>& phip =
153  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
154 
155  if (psiName_ == "none" && rhoName_ == "none")
156  {
157  operator==(p0_ - 0.5*(1.0 - pos(phip))*magSqr(Up));
158  }
159  else if (rhoName_ == "none")
160  {
161  const fvPatchField<scalar>& psip =
162  patch().lookupPatchField<volScalarField, scalar>(psiName_);
163 
164  if (gamma_ > 1.0)
165  {
166  scalar gM1ByG = (gamma_ - 1.0)/gamma_;
167 
168  operator==
169  (
170  p0_
171  /pow
172  (
173  (1.0 + 0.5*psip*gM1ByG*(1.0 - pos(phip))*magSqr(Up)),
174  1.0/gM1ByG
175  )
176  );
177  }
178  else
179  {
180  operator==(p0_/(1.0 + 0.5*psip*(1.0 - pos(phip))*magSqr(Up)));
181  }
182  }
183  else if (psiName_ == "none")
184  {
185  const fvPatchField<scalar>& rho =
186  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
187 
188  operator==(p0_ - 0.5*rho*(1.0 - pos(phip))*magSqr(Up));
189  }
190  else
191  {
193  (
194  "timeVaryingUniformTotalPressureFvPatchScalarField::updateCoeffs()"
195  ) << " rho or psi set inconsitently, rho = " << rhoName_
196  << ", psi = " << psiName_ << ".\n"
197  << " Set either rho or psi or neither depending on the "
198  "definition of total pressure.\n"
199  << " Set the unused variables to 'none'.\n"
200  << " on patch " << this->patch().name()
201  << " of field " << this->dimensionedInternalField().name()
202  << " in file " << this->dimensionedInternalField().objectPath()
203  << exit(FatalError);
204  }
205 
206  fixedValueFvPatchScalarField::updateCoeffs();
207 }
208 
209 
211 {
212  updateCoeffs(patch().lookupPatchField<volVectorField, vector>(UName_));
213 }
214 
215 
217 write(Ostream& os) const
218 {
220  if (UName_ != "U")
221  {
222  os.writeKeyword("U") << UName_ << token::END_STATEMENT << nl;
223  }
224  if (phiName_ != "phi")
225  {
226  os.writeKeyword("phi") << phiName_ << token::END_STATEMENT << nl;
227  }
228  os.writeKeyword("rho") << rhoName_ << token::END_STATEMENT << nl;
229  os.writeKeyword("psi") << psiName_ << token::END_STATEMENT << nl;
230  os.writeKeyword("gamma") << gamma_ << token::END_STATEMENT << nl;
231  os.writeKeyword("p0") << p0_ << token::END_STATEMENT << nl;
232  timeSeries_.write(os);
233  writeEntry("value", os);
234 }
235 
236 
237 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
238 
239 namespace Foam
240 {
242  (
245  );
246 }
247 
248 // ************************ vim: set sw=4 sts=4 et: ************************ //