ReactingCloudI.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) 2011-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
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// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
29
30template<class CloudType>
33{
34 return *cloudCopyPtr_;
35}
36
37
38template<class CloudType>
41{
42 return constProps_;
43}
44
45
46template<class CloudType>
49{
50 return constProps_;
51}
52
53
54template<class CloudType>
57{
58 return *compositionModel_;
59}
60
61
62template<class CloudType>
65{
66 return *phaseChangeModel_;
67}
68
69
70template<class CloudType>
73{
74 return *phaseChangeModel_;
75}
76
77
78template<class CloudType>
81{
82 return rhoTrans_[i];
83}
84
85
86template<class CloudType>
87inline
90{
91 return rhoTrans_;
92}
93
94
95template<class CloudType>
98{
99 return rhoTrans_;
100}
101
102
103template<class CloudType>
105(
106 const label i,
108) const
109{
110 if (this->solution().coupled())
111 {
112 if (this->solution().semiImplicit("Yi"))
113 {
114 tmp<volScalarField> trhoTrans
115 (
117 (
119 (
120 this->name() + ":rhoTrans",
121 this->db().time().timeName(),
122 this->db(),
123 IOobject::NO_READ,
124 IOobject::NO_WRITE,
125 false
126 ),
127 this->mesh(),
129 )
130 );
131
132 volScalarField& sourceField = trhoTrans.ref();
133
134 sourceField.primitiveFieldRef() =
135 rhoTrans_[i]/(this->db().time().deltaTValue()*this->mesh().V());
136
137 const dimensionedScalar YiSMALL("YiSMALL", dimless, SMALL);
138
139 return
140 fvm::Sp(neg(sourceField)*sourceField/(Yi + YiSMALL), Yi)
141 + pos0(sourceField)*sourceField;
142 }
143 else
144 {
146 fvScalarMatrix& fvm = tfvm.ref();
147
148 fvm.source() = -rhoTrans_[i]/this->db().time().deltaTValue();
149
150 return tfvm;
151 }
152 }
153
155}
156
157
158template<class CloudType>
161{
163 (
165 (
167 (
168 this->name() + ":rhoTrans",
169 this->db().time().timeName(),
170 this->db(),
171 IOobject::NO_READ,
172 IOobject::NO_WRITE,
173 false
174 ),
175 this->mesh(),
177 (
178 rhoTrans_[0].dimensions()/dimTime/dimVolume, Zero
179 )
180 )
181 );
182
183 if (this->solution().coupled())
184 {
185 scalarField& rhoi = tRhoi.ref();
186 rhoi = rhoTrans_[i]/(this->db().time().deltaTValue()*this->mesh().V());
187 }
188
189 return tRhoi;
190}
191
192
193template<class CloudType>
196{
198 (
200 (
202 (
203 this->name() + ":rhoTrans",
204 this->db().time().timeName(),
205 this->db(),
206 IOobject::NO_READ,
207 IOobject::NO_WRITE,
208 false
209 ),
210 this->mesh(),
213 rhoTrans_[0].dimensions()/dimTime/dimVolume, Zero
214 )
216 );
217
218 if (this->solution().coupled())
219 {
220 scalarField& sourceField = trhoTrans.ref();
221 forAll(rhoTrans_, i)
223 sourceField += rhoTrans_[i];
224 }
225
226 sourceField /= this->db().time().deltaTValue()*this->mesh().V();
227 }
228
229 return trhoTrans;
231
232
233template<class CloudType>
236{
237 if (this->solution().coupled())
238 {
240 (
242 (
244 (
245 this->name() + ":rhoTrans",
246 this->db().time().timeName(),
247 this->db(),
248 IOobject::NO_READ,
249 IOobject::NO_WRITE,
250 false
251 ),
252 this->mesh(),
253 dimensionedScalar(dimMass/dimTime/dimVolume, Zero)
254 )
255 );
256
257 scalarField& sourceField = trhoTrans.ref();
259 if (this->solution().semiImplicit("rho"))
260 {
261
262 forAll(rhoTrans_, i)
263 {
264 sourceField += rhoTrans_[i];
266 sourceField /= this->db().time().deltaTValue()*this->mesh().V();
267
268 return fvm::SuSp(trhoTrans()/rho, rho);
269 }
270 else
271 {
272 tmp<fvScalarMatrix> tfvm(new fvScalarMatrix(rho, dimMass/dimTime));
273 fvScalarMatrix& fvm = tfvm.ref();
274
275 forAll(rhoTrans_, i)
276 {
277 sourceField += rhoTrans_[i];
278 }
279
280 fvm.source() = -trhoTrans()/this->db().time().deltaT();
281
282 return tfvm;
283 }
284 }
285
286 return tmp<fvScalarMatrix>::New(rho, dimMass/dimTime);
287}
288
289
290// ************************************************************************* //
Templated reacting parcel composition model class Consists of carrier species (via thermo package),...
Class to hold DSMC particle constant properties.
Definition: DSMCParcel.H:82
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
Templated phase change model class.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
Templated base class for reacting cloud.
Definition: ReactingCloud.H:72
const parcelType::constantProperties & constProps() const
Return the constant properties.
const PtrList< volScalarField::Internal > & rhoTrans() const
Return const access to mass source fields.
const CompositionModel< ReactingCloud< CloudType > > & composition() const
Return const access to reacting composition model.
const ReactingCloud & cloudCopy() const
Return a reference to the cloud copy.
tmp< fvScalarMatrix > SYi(const label i, volScalarField &Yi) const
Return mass source term for specie i - specie eqn.
tmp< volScalarField::Internal > Srho() const
Return tmp total mass source for carrier phase.
const PhaseChangeModel< ReactingCloud< CloudType > > & phaseChange() const
Return const access to reacting phase change model.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
Field< Type > & source() noexcept
Definition: fvMatrix.H:458
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:66
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
dynamicFvMesh & mesh
word timeName
Definition: getTimeIndex.H:3
const dimensionSet dimless
Dimensionless.
dimensionedScalar pos0(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
dimensionedScalar neg(const dimensionedScalar &ds)
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333