ReactingHeterogeneousParcelI.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2018-2019 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class ParcelType>
33 :
34  ParcelType::constantProperties(),
35  hRetentionCoeff_(this->dict_, 0.0)
36 {}
37 
38 
39 template<class ParcelType>
42 (
43  const constantProperties& cp
44 )
45 :
46  ParcelType::constantProperties(cp),
47  hRetentionCoeff_(cp.hRetentionCoeff_)
48 {}
49 
50 
51 template<class ParcelType>
54 (
55  const dictionary& parentDict
56 )
57 :
58  ParcelType::constantProperties(parentDict),
59  hRetentionCoeff_(this->dict_, "hRetentionCoeff")
60 {}
61 
62 
63 template<class ParcelType>
65 (
66  const polyMesh& mesh,
67  const barycentric& coordinates,
68  const label celli,
69  const label tetFacei,
70  const label tetPti
71 )
72 :
73  ParcelType(mesh, coordinates, celli, tetFacei, tetPti),
74  F_(0),
75  canCombust_(0)
76 {}
77 
78 
79 template<class ParcelType>
81 (
82  const polyMesh& mesh,
83  const vector& position,
84  const label celli
85 )
86 :
87  ParcelType(mesh, position, celli),
88  F_(0),
89  canCombust_(0)
90 {}
91 
92 
93 template<class ParcelType>
95 (
96  const polyMesh& mesh,
97  const barycentric& coordinates,
98  const label celli,
99  const label tetFacei,
100  const label tetPti,
101  const label typeId,
102  const scalar nParticle0,
103  const scalar d0,
104  const scalar dTarget0,
105  const vector& U0,
106  const vector& f0,
107  const vector& angularMomentum0,
108  const vector& torque0,
109  const scalarField& Y,
110  const scalarField& F,
111  const constantProperties& constProps
112 )
113 :
114  ParcelType
115  (
116  mesh,
117  coordinates,
118  celli,
119  tetFacei,
120  tetPti,
121  typeId,
122  nParticle0,
123  d0,
124  dTarget0,
125  U0,
126  f0,
127  angularMomentum0,
128  torque0,
129  Y,
130  constProps
131  ),
132  F_(F),
133  canCombust_(0)
134 {}
135 
137 // * * * * * * * * * constantProperties Member Functions * * * * * * * * * * //
138 
139 
140 template<class ParcelType>
141 inline Foam::scalar
143 hRetentionCoeff() const
144 {
145  scalar value = hRetentionCoeff_.value();
146 
147  if ((value < 0) || (value > 1))
148  {
150  << "hRetentionCoeff must be in the range 0 to 1" << nl
151  << exit(FatalError) << endl;
152  }
153 
154  return value;
155 }
156 
157 
158 // * * * * * * * * * * ThermoParcel Member Functions * * * * * * * * * * * * //
159 
160 template<class ParcelType>
162 F() const
163 {
164  return F_;
165 }
166 
167 
168 template<class ParcelType>
170 F()
171 {
172  return F_;
173 }
174 
175 
176 template<class ParcelType>
177 inline Foam::label
179 {
180  return canCombust_;
181 }
182 
183 
184 template<class ParcelType>
186 {
187  return canCombust_;
188 }
189 
190 
191 // ************************************************************************* //
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:125
error FatalError
Error stream (stdout output on all processes), with additional &#39;FOAM FATAL ERROR&#39; header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Definition: dictionary.H:129
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:598
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:50
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:531
ReactingHeterogeneousParcel(const polyMesh &mesh, const barycentric &coordinates, const label celli, const label tetFacei, const label tetPti)
Construct from mesh, position and topology.
scalar hRetentionCoeff() const
Return const access to the fraction of enthalpy retained by.
const scalarField & F() const
Return const access to F.
dynamicFvMesh & mesh
Class to hold reacting particle constant properties.
volVectorField F(fluid.F())
label canCombust() const
Return const access to the canCombust flag.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
const volScalarField & cp
PtrList< coordinateSystem > coordinates(solidRegions.size())
PtrList< volScalarField > & Y
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74