MovingPhaseModel.H
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) 2015 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 Class
25  Foam::MovingPhaseModel
26 
27 Description
28  Class which represents a moving fluid phase. Holds the velocity, fluxes and
29  turbulence model. Provides access to the turbulent quantities.
30 
31  Possible future extensions include separating the turbulent fuctionality
32  into another layer. It should also be possible to replace this layer with a
33  stationary phase model, in order to model packed beds or simple porous
34  media. This would probably require extra functionality, such as returning
35  the inputs into the general pressure equation (A, HbyA, etc ...).
36 
37  Note that this class does not return the turbulence model, it just provides
38  indirect access to the turbulent data. This is so a layer without
39  turbulence modelling (such as a stationary model) could be substituted.
40 
41 SourceFiles
42  MovingPhaseModel.C
43 
44 \*---------------------------------------------------------------------------*/
45 
46 #ifndef MovingPhaseModel_H
47 #define MovingPhaseModel_H
48 
49 #include "phaseModel.H"
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 
57 /*---------------------------------------------------------------------------*\
58  Class phaseModel Declaration
59 \*---------------------------------------------------------------------------*/
60 
61 template<class BasePhaseModel>
62 class MovingPhaseModel
63 :
64  public BasePhaseModel
65 {
66  // Private data
67 
68  //- Velocity field
69  volVectorField U_;
70 
71  //- Flux
72  surfaceScalarField phi_;
73 
74  //- Volumetric flux
75  surfaceScalarField alphaPhi_;
76 
77  //- Mass flux
78  surfaceScalarField alphaRhoPhi_;
79 
80  //- Lagrangian acceleration field (needed for virtual-mass)
81  volVectorField DUDt_;
82 
83  //- Turbulence model
85 
86  //- Continuity error
87  volScalarField continuityError_;
88 
89  //- Phase diffusivity divided by the momentum coefficient.
90  // Used for implicit treatment of the phase pressure and dispersion
92 
93 
94  // Private static member functions
95 
96  //- Calculate and return the flux field
98 
99 
100 public:
101 
102  // Constructors
103 
105  (
106  const phaseSystem& fluid,
107  const word& phaseName,
108  const label index
109  );
110 
111 
112  //- Destructor
113  virtual ~MovingPhaseModel();
114 
115 
116  // Member Functions
117 
118  //- Correct the phase properties other than the thermo and turbulence
119  virtual void correct();
120 
121  //- Correct the kinematics
122  virtual void correctKinematics();
123 
124  //- Correct the turbulence
125  virtual void correctTurbulence();
126 
127  //- Correct the energy transport e.g. alphat
128  virtual void correctEnergyTransport();
129 
130  //- Return the momentum equation
131  virtual tmp<fvVectorMatrix> UEqn();
132 
133 
134  // Implicit phase pressure and dispersion support
135 
136  //- Return the phase diffusivity divided by the momentum coefficient
137  virtual const surfaceScalarField& DbyA() const;
138 
139  //- Set the phase diffusivity divided by the momentum coefficient
140  virtual void DbyA(const tmp<surfaceScalarField>& DbyA);
141 
142 
143  // Momentum
144 
145  //- Constant access the velocity
146  virtual tmp<volVectorField> U() const;
147 
148  //- Access the velocity
149  virtual volVectorField& U();
150 
151  //- Return the substantive acceleration
152  virtual tmp<volVectorField> DUDt() const;
153 
154  //- Constant access the continuity error
155  virtual tmp<volScalarField> continuityError() const;
156 
157  //- Constant access the volumetric flux
158  virtual tmp<surfaceScalarField> phi() const;
159 
160  //- Access the volumetric flux
161  virtual surfaceScalarField& phi();
162 
163  //- Constant access the volumetric flux of the phase
164  virtual tmp<surfaceScalarField> alphaPhi() const;
165 
166  //- Access the volumetric flux of the phase
167  virtual surfaceScalarField& alphaPhi();
168 
169  //- Constant access the mass flux of the phase
170  virtual tmp<surfaceScalarField> alphaRhoPhi() const;
171 
172  //- Access the mass flux of the phase
173  virtual surfaceScalarField& alphaRhoPhi();
174 
175 
176  // Turbulence
177 
178  //- Return the turbulence model
179  virtual const phaseCompressibleTurbulenceModel& turbulence() const;
180 };
181 
182 
183 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
184 
185 } // End namespace Foam
186 
187 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
188 
189 #ifdef NoRepository
190 # include "MovingPhaseModel.C"
191 #endif
192 
193 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
194 
195 #endif
196 
197 // ************************************************************************* //
virtual const phaseCompressibleTurbulenceModel & turbulence() const
Return the turbulence model.
virtual tmp< volVectorField > U() const
Constant access the velocity.
virtual tmp< fvVectorMatrix > UEqn()
Return the momentum equation.
virtual tmp< surfaceScalarField > alphaPhi() const
Constant access the volumetric flux of the phase.
multiphaseSystem & fluid
Definition: createFields.H:10
virtual ~MovingPhaseModel()
Destructor.
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:67
A class for handling words, derived from string.
Definition: word.H:59
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
virtual tmp< volVectorField > DUDt() const
Return the substantive acceleration.
virtual tmp< volScalarField > continuityError() const
Constant access the continuity error.
Namespace for OpenFOAM.
virtual void correctKinematics()
Correct the kinematics.
virtual tmp< surfaceScalarField > phi() const
Constant access the volumetric flux.
virtual void correctEnergyTransport()
Correct the energy transport e.g. alphat.
virtual const surfaceScalarField & DbyA() const
Return the phase diffusivity divided by the momentum coefficient.
Class which represents a moving fluid phase. Holds the velocity, fluxes and turbulence model...
MovingPhaseModel(const phaseSystem &fluid, const word &phaseName, const label index)
virtual tmp< surfaceScalarField > alphaRhoPhi() const
Constant access the mass flux of the phase.
virtual void correctTurbulence()
Correct the turbulence.
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
virtual void correct()
Correct the phase properties other than the thermo and turbulence.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
A class for managing temporary objects.
Definition: PtrList.H:118