StationaryPhaseModel.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) 2015-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 
26 #include "StationaryPhaseModel.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class BasePhaseModel>
32 (
33  const phaseSystem& fluid,
34  const word& phaseName,
35  const bool referencePhase,
36  const label index
37 )
38 :
39  BasePhaseModel(fluid, phaseName, referencePhase, index)
40 {}
41 
42 
43 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
44 
45 template<class BasePhaseModel>
47 {}
48 
49 
50 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
51 
52 template<class BasePhaseModel>
54 {
55  return true;
56 }
57 
58 
59 template<class BasePhaseModel>
62 {
64  << "Cannot construct a momentum equation for a stationary phase"
65  << abort(FatalError);
66 
67  return tmp<fvVectorMatrix>();
68 }
69 
70 
71 template<class BasePhaseModel>
74 {
76  << "Cannot construct a momentum equation for a stationary phase"
77  << abort(FatalError);
78 
79  return tmp<fvVectorMatrix>();
80 }
81 
82 
83 template<class BasePhaseModel>
86 {
88  << "Cannot access the velocity of a stationary phase"
89  << abort(FatalError);
90 
91  return volVectorField::null();
92 }
93 
94 
95 template<class BasePhaseModel>
98 {
100  << "Cannot access the velocity of a stationary phase"
101  << abort(FatalError);
102 
103  return const_cast<volVectorField&>(volVectorField::null());
104 }
105 
106 
107 template<class BasePhaseModel>
110 {
112  << "Cannot access the velocity of a stationary phase"
113  << abort(FatalError);
114 
115  return volVectorField::null();
116 }
117 
118 
119 template<class BasePhaseModel>
122 {
124  << "Cannot access the flux of a stationary phase"
125  << abort(FatalError);
126 
127  return surfaceScalarField::null();
128 }
129 
130 
131 template<class BasePhaseModel>
134 {
136  << "Cannot access the flux of a stationary phase"
137  << abort(FatalError);
138 
139  return const_cast<surfaceScalarField&>(surfaceScalarField::null());
140 }
141 
142 
143 template<class BasePhaseModel>
146 {
148  << "Cannot access the flux of a stationary phase"
149  << abort(FatalError);
150 
151  return surfaceScalarField::null();
152 }
153 
154 
155 template<class BasePhaseModel>
158 {
160  << "Cannot access the face velocity of a stationary phase"
161  << abort(FatalError);
162 
164  return Uf_;
165 }
166 
167 
168 template<class BasePhaseModel>
171 {
173  << "Cannot access the face velocity of a stationary phase"
174  << abort(FatalError);
175 
176  return const_cast<surfaceVectorField&>(surfaceVectorField::null());
177 }
178 
179 
180 template<class BasePhaseModel>
183 {
185  << "Cannot access the face velocity of a stationary phase"
186  << abort(FatalError);
187 
188  return surfaceVectorField::null();
189 }
190 
191 
192 template<class BasePhaseModel>
195 {
197  << "Cannot access the flux of a stationary phase"
198  << abort(FatalError);
199 
200  return surfaceScalarField::null();
201 }
202 
203 
204 template<class BasePhaseModel>
207 {
209  << "Cannot access the volumetric flux of a stationary phase"
210  << abort(FatalError);
211 
212  return const_cast<surfaceScalarField&>(surfaceScalarField::null());
213 }
214 
215 
216 template<class BasePhaseModel>
219 {
221  << "Cannot access the flux of a stationary phase"
222  << abort(FatalError);
223 
224  return surfaceScalarField::null();
225 }
226 
227 
228 template<class BasePhaseModel>
231 {
233  << "Cannot access the flux of a stationary phase"
234  << abort(FatalError);
235 
236  return surfaceScalarField::null();
237 }
238 
239 
240 template<class BasePhaseModel>
243 {
245  << "Cannot access the flux of a stationary phase"
246  << abort(FatalError);
247 
248  return const_cast<surfaceScalarField&>(surfaceScalarField::null());
249 }
250 
251 
252 template<class BasePhaseModel>
255 {
257  << "Cannot access the flux of a stationary phase"
258  << abort(FatalError);
259 
260  return surfaceScalarField::null();
261 }
262 
263 
264 template<class BasePhaseModel>
267 {
269  << "Cannot calculate UgradU of a stationary phase"
270  << abort(FatalError);
271 
272  return tmp<fvVectorMatrix>(nullptr);
273 }
274 
275 
276 template<class BasePhaseModel>
279 {
281  << "Cannot calculate DUDt of a stationary phase"
282  << abort(FatalError);
283 
284  return tmp<fvVectorMatrix>(nullptr);
285 }
286 
287 
288 template<class BasePhaseModel>
291 {
293  << "Cannot access the continuityError of a stationary phase"
294  << abort(FatalError);
295 
296  return volScalarField::null();
297 }
298 
299 
300 template<class BasePhaseModel>
303 {
305  << "Cannot access the kinetic energy of a stationary phase"
306  << abort(FatalError);
307 
308  return volScalarField::null();
309 }
310 
311 
312 template<class BasePhaseModel>
315 {
317  << "Cannot access the dilatation rate of a stationary phase"
318  << abort(FatalError);
319 
320  static autoPtr<volScalarField> divU_;
321  return divU_;
322 }
323 
324 
325 template<class BasePhaseModel>
327 (
329 )
330 {
332  << "Cannot set the dilatation rate of a stationary phase"
333  << abort(FatalError);
334 }
335 
336 
337 template<class BasePhaseModel>
340 {
341  return volScalarField::New
342  (
343  IOobject::groupName("k", this->name()),
344  this->mesh(),
346  );
347 }
348 
349 
350 template<class BasePhaseModel>
353 {
355  << "Cannot access the pPrime of a stationary phase"
356  << abort(FatalError);
357 
358  return surfaceScalarField::null();
359 }
360 
361 
362 // ************************************************************************* //
Generic GeometricField class.
static const GeometricField< Type, PatchField, GeoMesh > & null()
Return a null geometric field.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
static word groupName(Name name, const word &group)
virtual tmp< fvVectorMatrix > UEqn()
Return the momentum equation.
virtual tmp< volScalarField > k() const
Return the turbulent kinetic energy.
virtual tmp< fvVectorMatrix > DUDt() const
Return the substantive acceleration matrix.
virtual tmp< volVectorField > U() const
Return the velocity.
virtual surfaceScalarField & alphaRhoPhiRef()
Access the mass flux of the phase.
virtual const autoPtr< volScalarField > & divU() const
Return the phase dilatation rate (d(alpha)/dt + div(alpha*phi))
virtual surfaceScalarField & phiRef()
Access the volumetric flux.
virtual tmp< surfaceScalarField > alphaPhi() const
Return the volumetric flux of the phase.
virtual tmp< volScalarField > K() const
Return the phase kinetic energy.
virtual tmp< surfaceScalarField > pPrimef() const
Return the face-phase-pressure'.
virtual const autoPtr< surfaceVectorField > & Uf() const
Return the face velocity.
virtual surfaceVectorField & UfRef()
Access the face velocity.
virtual ~StationaryPhaseModel()
Destructor.
StationaryPhaseModel(const phaseSystem &fluid, const word &phaseName, const bool referencePhase, const label index)
virtual surfaceScalarField & alphaPhiRef()
Access the volumetric flux of the phase.
virtual tmp< surfaceScalarField > alphaRhoPhi() const
Return the mass flux of the phase.
virtual tmp< surfaceScalarField > phi() const
Return the volumetric flux.
virtual tmp< fvVectorMatrix > UfEqn()
Return the momentum equation for the face-based algorithm.
virtual volVectorField & URef()
Access the velocity.
virtual bool stationary() const
Return whether the phase is stationary.
virtual tmp< fvVectorMatrix > UgradU() const
Return the velocity transport matrix.
virtual tmp< volScalarField > continuityError() const
Return the continuity error.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:73
A class for managing temporary objects.
Definition: tmp.H:55
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
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
dimensionedSymmTensor sqr(const dimensionedVector &dv)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const dimensionSet dimVelocity
error FatalError
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.