FreeFOAM The Cross-Platform CFD Toolkit
syringePressureFvPatchScalarField.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 
27 #include <finiteVolume/volMesh.H>
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
40 (
41  const fvPatch& p,
43 )
44 :
45  fixedValueFvPatchScalarField(p, iF),
46  curTimeIndex_(-1)
47 {}
48 
49 
51 (
53  const fvPatch& p,
55  const fvPatchFieldMapper& mapper
56 )
57 :
58  fixedValueFvPatchScalarField(sppsf, p, iF, mapper),
59  Ap_(sppsf.Ap_),
60  Sp_(sppsf.Sp_),
61  VsI_(sppsf.VsI_),
62  tas_(sppsf.tas_),
63  tae_(sppsf.tae_),
64  tds_(sppsf.tds_),
65  tde_(sppsf.tde_),
66  psI_(sppsf.psI_),
67  psi_(sppsf.psi_),
68  ams_(sppsf.ams_),
69  ams0_(sppsf.ams0_),
70  curTimeIndex_(-1)
71 {}
72 
73 
75 (
76  const fvPatch& p,
78  const dictionary& dict
79 )
80 :
81  fixedValueFvPatchScalarField(p, iF),
82  Ap_(readScalar(dict.lookup("Ap"))),
83  Sp_(readScalar(dict.lookup("Sp"))),
84  VsI_(readScalar(dict.lookup("VsI"))),
85  tas_(readScalar(dict.lookup("tas"))),
86  tae_(readScalar(dict.lookup("tae"))),
87  tds_(readScalar(dict.lookup("tds"))),
88  tde_(readScalar(dict.lookup("tde"))),
89  psI_(readScalar(dict.lookup("psI"))),
90  psi_(readScalar(dict.lookup("psi"))),
91  ams_(readScalar(dict.lookup("ams"))),
92  ams0_(ams_),
93  curTimeIndex_(-1)
94 {
95  scalar ps = (psI_*VsI_ + ams_/psi_)/Vs(db().time().value());
97 }
98 
99 
101 (
104 )
105 :
106  fixedValueFvPatchScalarField(sppsf, iF),
107  Ap_(sppsf.Ap_),
108  Sp_(sppsf.Sp_),
109  VsI_(sppsf.VsI_),
110  tas_(sppsf.tas_),
111  tae_(sppsf.tae_),
112  tds_(sppsf.tds_),
113  tde_(sppsf.tde_),
114  psI_(sppsf.psI_),
115  psi_(sppsf.psi_),
116  ams_(sppsf.ams_),
117  ams0_(sppsf.ams0_),
118  curTimeIndex_(-1)
119 {}
120 
121 
123 (
125 )
126 :
127  fixedValueFvPatchScalarField(sppsf),
128  Ap_(sppsf.Ap_),
129  Sp_(sppsf.Sp_),
130  VsI_(sppsf.VsI_),
131  tas_(sppsf.tas_),
132  tae_(sppsf.tae_),
133  tds_(sppsf.tds_),
134  tde_(sppsf.tde_),
135  psI_(sppsf.psI_),
136  psi_(sppsf.psi_),
137  ams_(sppsf.ams_),
138  ams0_(sppsf.ams0_),
139  curTimeIndex_(-1)
140 {}
141 
142 
143 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
144 
145 scalar 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().deltaT().value();
199 
200  const surfaceScalarField& phi =
201  db().lookupObject<surfaceScalarField>("phi");
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  {
216  FatalErrorIn("syringePressureFvPatchScalarField::updateCoeffs()")
217  << "dimensions of phi are not correct"
218  << "\n on patch " << this->patch().name()
219  << " of field " << this->dimensionedInternalField().name()
220  << " in file " << this->dimensionedInternalField().objectPath()
221  << exit(FatalError);
222  }
223 
224  scalar ps = (psI_*VsI_ + ams_/psi_)/Vs(t);
225 
226  operator==(ps);
227 
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 
247  writeEntry("value", os);
248 }
249 
250 
251 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
252 
254 (
257 );
258 
259 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
260 
261 } // End namespace Foam
262 
263 // ************************ vim: set sw=4 sts=4 et: ************************ //