plenumPressureFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2016-2024 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 "fieldMapper.H"
29 #include "volFields.H"
30 #include "surfaceFields.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
35 (
36  const fvPatch& p,
38  const dictionary& dict
39 )
40 :
41  fixedValueFvPatchScalarField(p, iF, dict),
42  gamma_(dict.lookup<scalar>("gamma", dimless)),
43  R_(dict.lookup<scalar>("R", dimGasConstant)),
44  supplyMassFlowRate_
45  (
46  dict.lookup<scalar>("supplyMassFlowRate", dimMass/dimTime)
47  ),
48  supplyTotalTemperature_
49  (
50  dict.lookup<scalar>("supplyTotalTemperature", dimTemperature)
51  ),
52  plenumVolume_(dict.lookup<scalar>("plenumVolume", dimVolume)),
53  plenumDensity_(dict.lookup<scalar>("plenumDensity", dimDensity)),
54  plenumTemperature_
55  (
56  dict.lookup<scalar>("plenumTemperature", dimTemperature)
57  ),
58  rho_(1.0),
59  hasRho_(false),
60  inletAreaRatio_(dict.lookup<scalar>("inletAreaRatio", unitFraction)),
61  inletDischargeCoefficient_
62  (
63  dict.lookup<scalar>("inletDischargeCoefficient", unitFraction)
64  ),
65  timeScale_(dict.lookupOrDefault<scalar>("timeScale", dimTime, 0.0)),
66  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
67  UName_(dict.lookupOrDefault<word>("U", "U"))
68 {
69  if (dict.found("rho"))
70  {
71  rho_ = dict.lookup<scalar>("rho", dimDensity);
72  hasRho_ = true;
73  }
74 }
75 
76 
78 (
80  const fvPatch& p,
82  const fieldMapper& mapper
83 )
84 :
85  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
86  gamma_(ptf.gamma_),
87  R_(ptf.R_),
88  supplyMassFlowRate_(ptf.supplyMassFlowRate_),
89  supplyTotalTemperature_(ptf.supplyTotalTemperature_),
90  plenumVolume_(ptf.plenumVolume_),
91  plenumDensity_(ptf.plenumDensity_),
92  plenumTemperature_(ptf.plenumTemperature_),
93  rho_(ptf.rho_),
94  hasRho_(ptf.hasRho_),
95  inletAreaRatio_(ptf.inletAreaRatio_),
96  inletDischargeCoefficient_(ptf.inletDischargeCoefficient_),
97  timeScale_(ptf.timeScale_),
98  phiName_(ptf.phiName_),
99  UName_(ptf.UName_)
100 {}
101 
102 
104 (
107 )
108 :
109  fixedValueFvPatchScalarField(tppsf, iF),
110  gamma_(tppsf.gamma_),
111  R_(tppsf.R_),
112  supplyMassFlowRate_(tppsf.supplyMassFlowRate_),
113  supplyTotalTemperature_(tppsf.supplyTotalTemperature_),
114  plenumVolume_(tppsf.plenumVolume_),
115  plenumDensity_(tppsf.plenumDensity_),
116  plenumTemperature_(tppsf.plenumTemperature_),
117  rho_(tppsf.rho_),
118  hasRho_(tppsf.hasRho_),
119  inletAreaRatio_(tppsf.inletAreaRatio_),
120  inletDischargeCoefficient_(tppsf.inletDischargeCoefficient_),
121  timeScale_(tppsf.timeScale_),
122  phiName_(tppsf.phiName_),
123  UName_(tppsf.UName_)
124 {}
125 
126 
127 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
128 
130 {
131  if (updated())
132  {
133  return;
134  }
135 
136  // Patch properties
137  const fvPatchField<scalar>& p = *this;
138  const fvPatchField<scalar>& p_old =
139  db().lookupObject<volScalarField>
140  (
141  internalField().name()
142  ).oldTime().boundaryField()[patch().index()];
143  const fvPatchField<vector>& U =
144  patch().lookupPatchField<volVectorField, vector>(UName_);
145  const fvsPatchField<scalar>& phi =
146  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
147 
148  // Get the timestep
149  const scalar dt = db().time().deltaTValue();
150 
151  // Check if operating at a new time index and update the old-time properties
152  // if so
153  if (timeIndex_ != db().time().timeIndex())
154  {
155  timeIndex_ = db().time().timeIndex();
156  plenumDensityOld_ = plenumDensity_;
157  plenumTemperatureOld_ = plenumTemperature_;
158  }
159 
160  // Calculate the current mass flow rate
161  scalar massFlowRate(1.0);
162  if (phi.internalField().dimensions() == dimVolumetricFlux)
163  {
164  if (hasRho_)
165  {
166  massFlowRate = - gSum(rho_*phi);
167  }
168  else
169  {
171  << "The density must be specified when using a volumetric flux."
172  << exit(FatalError);
173  }
174  }
175  else if
176  (
177  phi.internalField().dimensions()
178  == dimMassFlux
179  )
180  {
181  if (hasRho_)
182  {
184  << "The density must be not specified when using a mass flux."
185  << exit(FatalError);
186  }
187  else
188  {
189  massFlowRate = - gSum(phi);
190  }
191  }
192  else
193  {
195  << "dimensions of phi are not correct"
196  << "\n on patch " << patch().name()
197  << " of field " << internalField().name()
198  << " in file " << internalField().objectPath() << nl
199  << exit(FatalError);
200  }
201 
202  // Calcaulate the specific heats
203  const scalar cv = R_/(gamma_ - 1), cp = R_*gamma_/(gamma_ - 1);
204 
205  // Calculate the new plenum properties
206  plenumDensity_ =
207  plenumDensityOld_
208  + (dt/plenumVolume_)*(supplyMassFlowRate_ - massFlowRate);
209  plenumTemperature_ =
210  plenumTemperatureOld_
211  + (dt/(plenumDensity_*cv*plenumVolume_))
212  * (
213  supplyMassFlowRate_
214  *(cp*supplyTotalTemperature_ - cv*plenumTemperature_)
215  - massFlowRate*R_*plenumTemperature_
216  );
217  const scalar plenumPressure = plenumDensity_*R_*plenumTemperature_;
218 
219  // Squared velocity magnitude at exit of channels
220  const scalarField U_e(magSqr(U/inletAreaRatio_));
221 
222  // Exit temperature to plenum temperature ratio
223  const scalarField r
224  (
225  1.0 - (gamma_ - 1.0)*U_e/(2.0*gamma_*R_*plenumTemperature_)
226  );
227 
228  // Quadratic coefficient (others not needed as b = +1.0 and c = -1.0)
229  const scalarField a
230  (
231  (1.0 - r)/(r*r*inletDischargeCoefficient_*inletDischargeCoefficient_)
232  );
233 
234  // Isentropic exit temperature to plenum temperature ratio
235  const scalarField s(2.0/(1.0 + sqrt(1.0 + 4.0*a)));
236 
237  // Exit pressure to plenum pressure ratio
238  const scalarField t(pow(s, gamma_/(gamma_ - 1.0)));
239 
240  // Limit to prevent outflow
241  const scalarField p_new
242  (
243  neg(phi)*t*plenumPressure + pos0(phi)*max(p, plenumPressure)
244  );
245 
246  // Relaxation fraction
247  const scalar oneByFraction = timeScale_/dt;
248  const scalar fraction = oneByFraction < 1.0 ? 1.0 : 1.0/oneByFraction;
249 
250  // Set the new value
251  operator==((1.0 - fraction)*p_old + fraction*p_new);
252  fixedValueFvPatchScalarField::updateCoeffs();
253 }
254 
255 
257 {
259  writeEntry(os, "gamma", gamma_);
260  writeEntry(os, "R", R_);
261  writeEntry(os, "supplyMassFlowRate", supplyMassFlowRate_);
262  writeEntry(os, "supplyTotalTemperature", supplyTotalTemperature_);
263  writeEntry(os, "plenumVolume", plenumVolume_);
264  writeEntry(os, "plenumDensity", plenumDensity_);
265  writeEntry(os, "plenumTemperature", plenumTemperature_);
266  if (hasRho_)
267  {
268  writeEntry(os, "rho", rho_);
269  }
270  writeEntry(os, "inletAreaRatio", inletAreaRatio_);
271  writeEntry
272  (
273  os,
274  "inletDischargeCoefficient",
275  inletDischargeCoefficient_
276  );
277  writeEntryIfDifferent<scalar>(os, "timeScale", 0.0, timeScale_);
278  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
279  writeEntryIfDifferent<word>(os, "U", "U", UName_);
280  writeEntry(os, "value", *this);
281 }
282 
283 
284 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
285 
286 namespace Foam
287 {
289  (
292  );
293 }
294 
295 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic GeometricField class.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base class for field mapping.
Definition: fieldMapper.H:48
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:88
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:229
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:83
This boundary condition provides a plenum pressure inlet condition. This condition creates a zero-dim...
plenumPressureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
U
Definition: pEqn.H:72
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Type gSum(const FieldField< Field, Type > &f)
VolField< vector > volVectorField
Definition: volFieldsFwd.H:65
dimensionedScalar pos0(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
const dimensionSet dimMassFlux
const dimensionSet dimVolumetricFlux
const dimensionSet dimGasConstant
const dimensionSet dimless
SurfaceField< scalar > surfaceScalarField
const dimensionSet dimTemperature
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionSet dimTime
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionSet dimDensity
const dimensionSet dimVolume
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar neg(const dimensionedScalar &ds)
const dimensionSet dimMass
bool cp(const fileName &src, const fileName &dst, const bool followLink=true)
Copy, recursively if necessary, the source to the destination.
Definition: POSIX.C:753
error FatalError
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
static const char nl
Definition: Ostream.H:266
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const unitConversion unitFraction
label timeIndex
Definition: getTimeIndex.H:4
dictionary dict
volScalarField & p
Foam::surfaceFields.