ReactingCloudI.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) 2011-2021 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 "fvmSup.H"
27 
28 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
29 
30 template<class CloudType>
33 {
34  return cloudCopyPtr_();
35 }
36 
37 
38 template<class CloudType>
39 inline const typename CloudType::particleType::constantProperties&
41 {
42  return constProps_;
43 }
44 
45 
46 template<class CloudType>
47 inline typename CloudType::particleType::constantProperties&
49 {
50  return constProps_;
51 }
52 
53 
54 template<class CloudType>
57 {
58  return phaseChangeModel_;
59 }
60 
61 
62 template<class CloudType>
65 {
66  return phaseChangeModel_();
67 }
68 
69 
70 template<class CloudType>
73 {
74  return rhoTrans_[i];
75 }
76 
77 
78 template<class CloudType>
79 inline
82 {
83  return rhoTrans_;
84 }
85 
86 
87 template<class CloudType>
90 {
91  return rhoTrans_;
92 }
93 
94 
95 template<class CloudType>
97 (
98  const label i,
99  const volScalarField& Yi
100 ) const
101 {
102  if (this->solution().coupled())
103  {
104  if (this->solution().semiImplicit("Yi"))
105  {
106  tmp<volScalarField> trhoTrans
107  (
109  (
110  this->name() + ":rhoTrans",
111  this->mesh(),
113  )
114  );
115 
116  volScalarField& sourceField = trhoTrans.ref();
117 
118  sourceField.primitiveFieldRef() =
119  rhoTrans_[i]/(this->db().time().deltaTValue()*this->mesh().V());
120 
121  const dimensionedScalar Yismall("Yismall", dimless, small);
122 
123  return
124  fvm::Sp(neg(sourceField)*sourceField/(Yi + Yismall), Yi)
125  + pos0(sourceField)*sourceField;
126  }
127  else
128  {
130  fvScalarMatrix& fvm = tfvm.ref();
131 
132  fvm.source() = -rhoTrans_[i]/this->db().time().deltaTValue();
133 
134  return tfvm;
135  }
136  }
137 
139 }
140 
141 
142 template<class CloudType>
145 {
147  (
149  (
150  this->name() + ":rhoTrans",
151  this->mesh(),
153  (
154  rhoTrans_[0].dimensions()/dimTime/dimVolume,
155  0
156  )
157  )
158  );
159 
160  if (this->solution().coupled())
161  {
162  scalarField& sourceField = trhoTrans.ref();
163  forAll(rhoTrans_, i)
164  {
165  sourceField += rhoTrans_[i];
166  }
167 
168  sourceField /= this->db().time().deltaTValue()*this->mesh().V();
169  }
170 
171  return trhoTrans;
172 }
173 
174 
175 template<class CloudType>
178 {
179  if (this->solution().coupled())
180  {
181  tmp<volScalarField> trhoTrans
182  (
184  (
185  this->name() + ":rhoTrans",
186  this->mesh(),
188  )
189  );
190 
191  scalarField& sourceField = trhoTrans.ref();
192 
193  if (this->solution().semiImplicit("rho"))
194  {
195 
196  forAll(rhoTrans_, i)
197  {
198  sourceField += rhoTrans_[i];
199  }
200  sourceField /= this->db().time().deltaTValue()*this->mesh().V();
201 
202  return fvm::SuSp(trhoTrans()/rho, rho);
203  }
204  else
205  {
207  fvScalarMatrix& fvm = tfvm.ref();
208 
209  forAll(rhoTrans_, i)
210  {
211  sourceField += rhoTrans_[i];
212  }
213 
214  fvm.source() = -trhoTrans()/this->db().time().deltaT();
215 
216  return tfvm;
217  }
218  }
219 
221 }
222 
223 
224 // ************************************************************************* //
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
tmp< GeometricField< Type, fvPatchField, volMesh > > SuSp(const volScalarField &sp, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:102
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
Templated phase change model class.
Definition: ReactingCloud.H:57
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
const PtrList< volScalarField::Internal > & rhoTrans() const
Return const access to mass source fields.
const dimensionSet dimless
tmp< fvScalarMatrix > SYi(const label i, const 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.
dimensionedScalar neg(const dimensionedScalar &ds)
const dimensionSet dimTime
const ReactingCloud & cloudCopy() const
Return a reference to the cloud copy.
dynamicFvMesh & mesh
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
const parcelType::constantProperties & constProps() const
Return the constant properties.
dimensionedScalar pos0(const dimensionedScalar &ds)
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
Field< Type > & source()
Definition: fvMatrix.H:292
const dimensionSet dimMass
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const tmp< volScalarField::Internal > & Sp
Definition: alphaSuSp.H:7
Templated base class for reacting cloud.
Definition: ReactingCloud.H:72
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet dimVolume
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
A class for managing temporary objects.
Definition: PtrList.H:53
const PhaseChangeModel< ReactingCloud< CloudType > > & phaseChange() const
Return const access to reacting phase change model.
Calculate the matrix for implicit and explicit sources.