phaseScalarTransport.H
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) 2019-2020 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::functionObjects::phaseScalarTransport
26 
27 Description
28  Evolves a passive scalar transport equation within one phase of a
29  multiphase simulation. The scalar is considered to be a phase-intensive
30  property; i.e., its value represents an amount per-unit of the phase. In
31  addition to the scalar, the function also writes out the product of the
32  volume fraction and the scalar, as this provides a phase-extensive field
33  which is often more convenient to post-process.
34 
35  Most entries are the same as for the \c scalarTransport function. Refer to
36  its documentation for details. Entries specific to this function are
37  detailed below. Note that the phase-name will be determined by stripping
38  the extension from the supplied field name.
39 
40  If \c alphaPhi is specified and found in the database then it will be used
41  to transport the field. This is the likely mode of operation if this
42  function is used with Euler-Euler solvers, in which \c alphaPhi -like
43  fluxes are available. If \c alphaPhi is not found, then a pressure-like
44  equation will be solved in order to construct it so that it exactly matches
45  the time-derivative of the phase volume or mass. This is likely to be
46  necessary in volume-of-fluid solvers where \c alphaPhi is not part of the
47  solution procedure. The pressure field name will be required in this case.
48 
49 Usage
50  \table
51  Property | Description | Req'd? | Default
52  alpha | Name of the volume-fraction field | no\\
53  | alpha.<phase-name>
54  alphaPhi | Name of the phase-flux field | no\\
55  | alphaPhi.<phase-name>
56  p | Name of the pressure field | no | p
57  residualAlpha | Small volume fraction used to stabilise the solution\\
58  | no | rootSmall
59  writeAlphaField | Also write out alpha multiplied by the field\\
60  | no | true
61  \endtable
62 
63  Example specification for interFoam:
64  \verbatim
65  phaseScalarTransport1
66  {
67  type phaseScalarTransport;
68  libs ("libsolverFunctionObjects.so");
69 
70  field s.water;
71  p p_rgh;
72  }
73  \endverbatim
74 
75  Example specification for reactingMultiphaseEulerFoam:
76  \verbatim
77  phaseScalarTransport1
78  {
79  type phaseScalarTransport;
80  libs ("libsolverFunctionObjects.so");
81 
82  field s.water;
83  alphaPhi alphaRhoPhi.water;
84  rho thermo:rho.water;
85  }
86  \endverbatim
87 
88 See also
89  Foam::functionObjects::fvMeshFunctionObject
90  Foam::functionObjects::scalarTransport
91 
92 SourceFiles
93  phaseScalarTransport.C
94 
95 \*---------------------------------------------------------------------------*/
96 
97 #ifndef functionObjects_phaseScalarTransport_H
98 #define functionObjects_phaseScalarTransport_H
99 
100 #include "fvMeshFunctionObject.H"
101 #include "volFields.H"
102 #include "fvOptionList.H"
103 
104 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
105 
106 namespace Foam
107 {
108 namespace functionObjects
109 {
110 
111 /*---------------------------------------------------------------------------*\
112  Class phaseScalarTransport Declaration
113 \*---------------------------------------------------------------------------*/
114 
115 class phaseScalarTransport
116 :
117  public fvMeshFunctionObject
118 {
119  // Private Data
120 
121  //- Name of field to process
122  const word fieldName_;
123 
124  //- Name of the phase in which to solve
125  const word phaseName_;
126 
127  //- Name of phase volume-fraction field (optional)
128  word alphaName_;
129 
130  //- Name of the phase flux field (optional)
131  word alphaPhiName_;
132 
133  //- Name of the mixture flux field (optional)
134  word phiName_;
135 
136  //- Name of density field (optional)
137  word rhoName_;
138 
139  //- Name of the pressure field (optional)
140  word pName_;
141 
142  //- Diffusion coefficient (optional)
143  scalar D_;
144 
145  //- Flag to indicate whether a constant, uniform D_ is specified
146  bool constantD_;
147 
148  //- Laminar diffusion coefficient (optional)
149  scalar alphaD_;
150 
151  //- Turbulent diffusion coefficient (optional)
152  scalar alphaDt_;
153 
154  //- Number of corrector iterations (optional)
155  int nCorr_;
156 
157  //- Residual volume-fraction
158  scalar residualAlpha_;
159 
160  //- Name of field whose schemes are used (optional)
161  word schemesField_;
162 
163  //- Flag to indicate whether to write the field multiplied by the phase
164  // fraction
165  bool writeAlphaField_;
166 
167  //- Run-time selectable finite volume options, e.g. sources, constraints
168  fv::optionList fvOptions_;
169 
170  //- The field
171  volScalarField s_;
172 
173  //- The field multiplied by the phase fraction
174  autoPtr<volScalarField> alphaSPtr_;
175 
176  //- Potential field used to generate the phase flux
177  autoPtr<volScalarField> PhiPtr_;
178 
179 
180  // Private Member Functions
181 
182  //- Return the potential field used to generate the phase flux.
183  // Constructed on demand.
184  volScalarField& Phi();
185 
186  //- Return the phase flux. Tries to look it up, and generates it if the
187  // lookup fails. The generation creates a crude guess for alphaPhi
188  // then solves a Laplacian to correct the flux to match the time
189  // derivative of alpha.
190  tmp<surfaceScalarField> alphaPhi();
191 
192  //- Return the diffusivity field
193  tmp<volScalarField> D(const surfaceScalarField& alphaPhi) const;
194 
195 
196 public:
197 
198  //- Runtime type information
199  TypeName("phaseScalarTransport");
200 
201 
202  // Constructors
203 
204  //- Construct from Time and dictionary
206  (
207  const word& name,
208  const Time& runTime,
209  const dictionary& dict
210  );
211 
212 
213  //- Destructor
214  virtual ~phaseScalarTransport();
215 
216 
217  // Member Functions
218 
219  //- Read the settings from the given dictionary
220  virtual bool read(const dictionary&);
221 
222  //- Solve for the evolution of the field
223  virtual bool execute();
224 
225  //- Do nothing. The field is registered and written automatically.
226  virtual bool write();
227 };
228 
229 
230 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
231 
232 } // End namespace functionObjects
233 } // End namespace Foam
234 
235 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
236 
237 #endif
238 
239 // ************************************************************************* //
dictionary dict
const word & name() const
Return the name of this functionObject.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
engineTime & runTime
virtual bool read(const dictionary &)
Read the settings from the given dictionary.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
virtual bool write()
Do nothing. The field is registered and written automatically.
A class for handling words, derived from string.
Definition: word.H:59
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
List of finite volume options.
Definition: fvOptionList.H:64
A class for managing temporary objects.
Definition: PtrList.H:53
virtual bool execute()
Solve for the evolution of the field.
TypeName("phaseScalarTransport")
Runtime type information.
phaseScalarTransport(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Namespace for OpenFOAM.