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-2018 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 )
39 :
40  fixedValueFvPatchScalarField(p, iF),
41  phiName_("phi"),
42  curTimeIndex_(-1)
43 {}
44 
45 
47 (
48  const fvPatch& p,
50  const dictionary& dict
51 )
52 :
53  fixedValueFvPatchScalarField(p, iF, dict, false),
54  Ap_(readScalar(dict.lookup("Ap"))),
55  Sp_(readScalar(dict.lookup("Sp"))),
56  VsI_(readScalar(dict.lookup("VsI"))),
57  tas_(readScalar(dict.lookup("tas"))),
58  tae_(readScalar(dict.lookup("tae"))),
59  tds_(readScalar(dict.lookup("tds"))),
60  tde_(readScalar(dict.lookup("tde"))),
61  psI_(readScalar(dict.lookup("psI"))),
62  psi_(readScalar(dict.lookup("psi"))),
63  ams_(readScalar(dict.lookup("ams"))),
64  ams0_(ams_),
65  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
66  curTimeIndex_(-1)
67 {
68  scalar ps = (psI_*VsI_ + ams_/psi_)/Vs(db().time().value());
70 }
71 
72 
74 (
76  const fvPatch& p,
78  const fvPatchFieldMapper& mapper
79 )
80 :
81  fixedValueFvPatchScalarField(sppsf, p, iF, mapper),
82  Ap_(sppsf.Ap_),
83  Sp_(sppsf.Sp_),
84  VsI_(sppsf.VsI_),
85  tas_(sppsf.tas_),
86  tae_(sppsf.tae_),
87  tds_(sppsf.tds_),
88  tde_(sppsf.tde_),
89  psI_(sppsf.psI_),
90  psi_(sppsf.psi_),
91  ams_(sppsf.ams_),
92  ams0_(sppsf.ams0_),
93  phiName_(sppsf.phiName_),
94  curTimeIndex_(-1)
95 {}
96 
97 
99 (
102 )
103 :
104  fixedValueFvPatchScalarField(sppsf, iF),
105  Ap_(sppsf.Ap_),
106  Sp_(sppsf.Sp_),
107  VsI_(sppsf.VsI_),
108  tas_(sppsf.tas_),
109  tae_(sppsf.tae_),
110  tds_(sppsf.tds_),
111  tde_(sppsf.tde_),
112  psI_(sppsf.psI_),
113  psi_(sppsf.psi_),
114  ams_(sppsf.ams_),
115  ams0_(sppsf.ams0_),
116  phiName_(sppsf.phiName_),
117  curTimeIndex_(-1)
118 {}
119 
120 
122 (
124 )
125 :
126  fixedValueFvPatchScalarField(sppsf),
127  Ap_(sppsf.Ap_),
128  Sp_(sppsf.Sp_),
129  VsI_(sppsf.VsI_),
130  tas_(sppsf.tas_),
131  tae_(sppsf.tae_),
132  tds_(sppsf.tds_),
133  tde_(sppsf.tde_),
134  psI_(sppsf.psI_),
135  psi_(sppsf.psi_),
136  ams_(sppsf.ams_),
137  ams0_(sppsf.ams0_),
138  phiName_(sppsf.phiName_),
139  curTimeIndex_(-1)
140 {}
141 
142 
143 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
144 
145 Foam::scalar Foam::syringePressureFvPatchScalarField::Vs(const scalar t) const
146 {
147  if (t < tas_)
148  {
149  return VsI_;
150  }
151  else if (t < tae_)
152  {
153  return
154  VsI_
155  + 0.5*Ap_*Sp_*sqr(t - tas_)/(tae_ - tas_);
156  }
157  else if (t < tds_)
158  {
159  return
160  VsI_
161  + 0.5*Ap_*Sp_*(tae_ - tas_)
162  + Ap_*Sp_*(t - tae_);
163  }
164  else if (t < tde_)
165  {
166  return
167  VsI_
168  + 0.5*Ap_*Sp_*(tae_ - tas_)
169  + Ap_*Sp_*(tds_ - tae_)
170  + Ap_*Sp_*(t - tds_)
171  - 0.5*Ap_*Sp_*sqr(t - tds_)/(tde_ - tds_);
172  }
173  else
174  {
175  return
176  VsI_
177  + 0.5*Ap_*Sp_*(tae_ - tas_)
178  + Ap_*Sp_*(tds_ - tae_)
179  + 0.5*Ap_*Sp_*(tde_ - tds_);
180  }
181 }
182 
183 
185 {
186  if (updated())
187  {
188  return;
189  }
190 
191  if (curTimeIndex_ != db().time().timeIndex())
192  {
193  ams0_ = ams_;
194  curTimeIndex_ = db().time().timeIndex();
195  }
196 
197  scalar t = db().time().value();
198  scalar deltaT = db().time().deltaTValue();
199 
200  const surfaceScalarField& phi =
201  db().lookupObject<surfaceScalarField>(phiName_);
202 
203  const fvsPatchField<scalar>& phip =
204  patch().patchField<surfaceScalarField, scalar>(phi);
205 
206  if (phi.dimensions() == dimVelocity*dimArea)
207  {
208  ams_ = ams0_ + deltaT*sum((*this*psi_)*phip);
209  }
210  else if (phi.dimensions() == dimDensity*dimVelocity*dimArea)
211  {
212  ams_ = ams0_ + deltaT*sum(phip);
213  }
214  else
215  {
217  << "dimensions of phi are not correct"
218  << "\n on patch " << this->patch().name()
219  << " of field " << this->internalField().name()
220  << " in file " << this->internalField().objectPath()
221  << exit(FatalError);
222  }
223 
224  scalar ps = (psI_*VsI_ + ams_/psi_)/Vs(t);
225 
226  operator==(ps);
227 
228  fixedValueFvPatchScalarField::updateCoeffs();
229 }
230 
231 
233 {
235 
236  os.writeKeyword("Ap") << Ap_ << token::END_STATEMENT << nl;
237  os.writeKeyword("Sp") << Sp_ << token::END_STATEMENT << nl;
238  os.writeKeyword("VsI") << VsI_ << token::END_STATEMENT << nl;
239  os.writeKeyword("tas") << tas_ << token::END_STATEMENT << nl;
240  os.writeKeyword("tae") << tae_ << token::END_STATEMENT << nl;
241  os.writeKeyword("tds") << tds_ << token::END_STATEMENT << nl;
242  os.writeKeyword("tde") << tde_ << token::END_STATEMENT << nl;
243  os.writeKeyword("psI") << psI_ << token::END_STATEMENT << nl;
244  os.writeKeyword("psi") << psi_ << token::END_STATEMENT << nl;
245  os.writeKeyword("ams") << ams_ << token::END_STATEMENT << nl;
246  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
247 
248  writeEntry("value", os);
249 }
250 
251 
252 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
253 
254 namespace Foam
255 {
257  (
260  );
261 }
262 
263 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
Foam::surfaceFields.
surfaceScalarField & phi
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
dimensionedSymmTensor sqr(const dimensionedVector &dv)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:362
Macros for easy insertion into run-time selection tables.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const dimensionSet & dimensions() const
Return dimensions.
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:265
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
const dimensionSet dimDensity
syringePressureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
label timeIndex
Definition: getTimeIndex.H:4
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
This boundary condition provides a pressure condition, obtained from a zero-D model of the cylinder o...
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576
const dimensionSet dimVelocity