syringePressureFvPatchScalarField.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) 2011-2023 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 
27 #include "volMesh.H"
29 #include "fvPatchFieldMapper.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, false),
42  Ap_(dict.lookup<scalar>("Ap")),
43  Sp_(dict.lookup<scalar>("Sp")),
44  VsI_(dict.lookup<scalar>("VsI")),
45  tas_(dict.lookup<scalar>("tas")),
46  tae_(dict.lookup<scalar>("tae")),
47  tds_(dict.lookup<scalar>("tds")),
48  tde_(dict.lookup<scalar>("tde")),
49  psI_(dict.lookup<scalar>("psI")),
50  psi_(dict.lookup<scalar>("psi")),
51  ams_(dict.lookup<scalar>("ams")),
52  ams0_(ams_),
53  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
54  curTimeIndex_(-1)
55 {
56  scalar ps = (psI_*VsI_ + ams_/psi_)/Vs(db().time().value());
58 }
59 
60 
62 (
64  const fvPatch& p,
66  const fvPatchFieldMapper& mapper
67 )
68 :
69  fixedValueFvPatchScalarField(sppsf, p, iF, mapper),
70  Ap_(sppsf.Ap_),
71  Sp_(sppsf.Sp_),
72  VsI_(sppsf.VsI_),
73  tas_(sppsf.tas_),
74  tae_(sppsf.tae_),
75  tds_(sppsf.tds_),
76  tde_(sppsf.tde_),
77  psI_(sppsf.psI_),
78  psi_(sppsf.psi_),
79  ams_(sppsf.ams_),
80  ams0_(sppsf.ams0_),
81  phiName_(sppsf.phiName_),
82  curTimeIndex_(-1)
83 {}
84 
85 
87 (
90 )
91 :
92  fixedValueFvPatchScalarField(sppsf, iF),
93  Ap_(sppsf.Ap_),
94  Sp_(sppsf.Sp_),
95  VsI_(sppsf.VsI_),
96  tas_(sppsf.tas_),
97  tae_(sppsf.tae_),
98  tds_(sppsf.tds_),
99  tde_(sppsf.tde_),
100  psI_(sppsf.psI_),
101  psi_(sppsf.psi_),
102  ams_(sppsf.ams_),
103  ams0_(sppsf.ams0_),
104  phiName_(sppsf.phiName_),
105  curTimeIndex_(-1)
106 {}
107 
108 
109 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
110 
111 Foam::scalar Foam::syringePressureFvPatchScalarField::Vs(const scalar t) const
112 {
113  if (t < tas_)
114  {
115  return VsI_;
116  }
117  else if (t < tae_)
118  {
119  return
120  VsI_
121  + 0.5*Ap_*Sp_*sqr(t - tas_)/(tae_ - tas_);
122  }
123  else if (t < tds_)
124  {
125  return
126  VsI_
127  + 0.5*Ap_*Sp_*(tae_ - tas_)
128  + Ap_*Sp_*(t - tae_);
129  }
130  else if (t < tde_)
131  {
132  return
133  VsI_
134  + 0.5*Ap_*Sp_*(tae_ - tas_)
135  + Ap_*Sp_*(tds_ - tae_)
136  + Ap_*Sp_*(t - tds_)
137  - 0.5*Ap_*Sp_*sqr(t - tds_)/(tde_ - tds_);
138  }
139  else
140  {
141  return
142  VsI_
143  + 0.5*Ap_*Sp_*(tae_ - tas_)
144  + Ap_*Sp_*(tds_ - tae_)
145  + 0.5*Ap_*Sp_*(tde_ - tds_);
146  }
147 }
148 
149 
151 {
152  if (updated())
153  {
154  return;
155  }
156 
157  if (curTimeIndex_ != db().time().timeIndex())
158  {
159  ams0_ = ams_;
160  curTimeIndex_ = db().time().timeIndex();
161  }
162 
163  scalar t = db().time().value();
164  scalar deltaT = db().time().deltaTValue();
165 
166  const surfaceScalarField& phi =
167  db().lookupObject<surfaceScalarField>(phiName_);
168 
169  const fvsPatchField<scalar>& phip =
170  patch().patchField<surfaceScalarField, scalar>(phi);
171 
172  if (phi.dimensions() == dimFlux)
173  {
174  ams_ = ams0_ + deltaT*sum((*this*psi_)*phip);
175  }
176  else if (phi.dimensions() == dimMassFlux)
177  {
178  ams_ = ams0_ + deltaT*sum(phip);
179  }
180  else
181  {
183  << "dimensions of phi are not correct"
184  << "\n on patch " << this->patch().name()
185  << " of field " << this->internalField().name()
186  << " in file " << this->internalField().objectPath()
187  << exit(FatalError);
188  }
189 
190  scalar ps = (psI_*VsI_ + ams_/psi_)/Vs(t);
191 
192  operator==(ps);
193 
194  fixedValueFvPatchScalarField::updateCoeffs();
195 }
196 
197 
199 {
201 
202  writeEntry(os, "Ap", Ap_);
203  writeEntry(os, "Sp", Sp_);
204  writeEntry(os, "VsI", VsI_);
205  writeEntry(os, "tas", tas_);
206  writeEntry(os, "tae", tae_);
207  writeEntry(os, "tds", tds_);
208  writeEntry(os, "tde", tde_);
209  writeEntry(os, "psI", psI_);
210  writeEntry(os, "psi", psi_);
211  writeEntry(os, "ams", ams_);
212  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
213 
214  writeEntry(os, "value", *this);
215 }
216 
217 
218 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
219 
220 namespace Foam
221 {
223  (
226  );
227 }
228 
229 // ************************************************************************* //
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...
const dimensionSet & dimensions() const
Return dimensions.
Generic GeometricField class.
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:160
Foam::fvPatchFieldMapper.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:87
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:231
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:251
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:82
This boundary condition provides a pressure condition, obtained from a zero-D model of the cylinder o...
syringePressureFvPatchScalarField(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:306
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
const dimensionSet dimMassFlux
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
SurfaceField< scalar > surfaceScalarField
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
error FatalError
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
const dimensionSet dimFlux
label timeIndex
Definition: getTimeIndex.H:4
dictionary dict
volScalarField & p
Foam::surfaceFields.