SprayCloud.C
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-2023 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 "SprayCloud.H"
27 #include "AtomisationModel.H"
28 #include "BreakupModel.H"
29 #include "ConeInjection.H"
30 #include "parcelThermo.H"
31 
32 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
33 
34 template<class CloudType>
36 {
37  atomisationModel_.reset
38  (
40  (
41  this->subModelProperties(),
42  *this
43  ).ptr()
44  );
45 
46  breakupModel_.reset
47  (
49  (
50  this->subModelProperties(),
51  *this
52  ).ptr()
53  );
54 }
55 
56 
57 template<class CloudType>
59 (
61 )
62 {
63  CloudType::cloudReset(c);
64 
65  atomisationModel_.reset(c.atomisationModel_.ptr());
66  breakupModel_.reset(c.breakupModel_.ptr());
67 }
68 
69 
70 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
71 
72 template<class CloudType>
74 (
75  const word& cloudName,
76  const volScalarField& rho,
77  const volVectorField& U,
78  const dimensionedVector& g,
79  const fluidThermo& carrierThermo,
80  bool readFields
81 )
82 :
83  CloudType(cloudName, rho, U, g, carrierThermo, false),
84  cloudCopyPtr_(nullptr),
85  atomisationModel_(nullptr),
86  breakupModel_(nullptr)
87 {
88  setModels();
89 
90  if (readFields)
91  {
92  parcelType::readFields(*this, this->composition());
93  this->deleteLostParticles();
94  }
95 
96  if (this->solution().resetSourcesOnStartup())
97  {
98  CloudType::resetSourceTerms();
99  }
100 }
101 
102 
103 template<class CloudType>
105 (
107  const word& name
108 )
109 :
110  CloudType(c, name),
111  cloudCopyPtr_(nullptr),
112  atomisationModel_(c.atomisationModel_->clone()),
113  breakupModel_(c.breakupModel_->clone())
114 {}
115 
116 
117 template<class CloudType>
119 (
120  const fvMesh& mesh,
121  const word& name,
122  const SprayCloud<CloudType>& c
123 )
124 :
125  CloudType(mesh, name, c),
126  cloudCopyPtr_(nullptr),
127  atomisationModel_(nullptr),
128  breakupModel_(nullptr)
129 {}
130 
131 
132 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
133 
134 template<class CloudType>
136 {}
137 
138 
139 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
140 
141 template<class CloudType>
143 (
144  parcelType& parcel
145 )
146 {
147  CloudType::setParcelThermoProperties(parcel);
148 
149  const liquidMixtureProperties& liqMix = this->composition().liquids();
150 
151  const scalarField& Y(parcel.Y());
152  scalarField X(liqMix.X(Y));
153  const scalar pc = this->p()[parcel.cell()];
154 
155  // override rho and Cp from constantProperties
156  parcel.Cp() = liqMix.Cp(pc, parcel.T(), X);
157  parcel.rho() = liqMix.rho(pc, parcel.T(), X);
158  parcel.sigma() = liqMix.sigma(pc, parcel.T(), X);
159  parcel.mu() = liqMix.mu(pc, parcel.T(), X);
160 }
161 
162 
163 template<class CloudType>
165 (
166  parcelType& parcel,
167  const label injectori
168 )
169 {
170  CloudType::checkParcelProperties(parcel, injectori);
171 
172  // store the initial size and position
173  parcel.d0() = parcel.d();
174  parcel.mass0() = parcel.mass();
175  parcel.position0() = parcel.position(this->mesh());
176 
177  parcel.y() = breakup().y0();
178  parcel.yDot() = breakup().yDot0();
179 
180  parcel.liquidCore() = atomisation().initLiquidCore();
181 
182  parcel.injector() = injectori;
183 }
184 
185 
186 template<class CloudType>
188 {
189  cloudCopyPtr_.reset
190  (
191  static_cast<SprayCloud<CloudType>*>
192  (
193  clone(this->name() + "Copy").ptr()
194  )
195  );
196 }
197 
198 
199 template<class CloudType>
201 {
202  cloudReset(cloudCopyPtr_());
203  cloudCopyPtr_.clear();
204 }
205 
206 
207 template<class CloudType>
209 {
210  if (this->solution().canEvolve())
211  {
212  typename parcelType::trackingData td(*this);
213 
214  this->solve(*this, td);
215  }
216 }
217 
218 
219 template<class CloudType>
221 {
222  CloudType::info();
223  scalar d32 = 1.0e+6*this->Dij(3, 2);
224  scalar d10 = 1.0e+6*this->Dij(1, 0);
225  scalar dMax = 1.0e+6*max(scalar(0), this->Dmax());
226  scalar pen = this->penetration(0.95);
227 
228  Info<< " D10, D32, Dmax (mu) = " << d10 << ", " << d32
229  << ", " << dMax << nl
230  << " Liquid penetration 95% mass (m) = " << pen << endl;
231 }
232 
233 
234 // ************************************************************************* //
Templated atomisation model class.
Templated break-up model class.
Definition: BreakupModel.H:55
void deleteLostParticles()
Remove lost particles from cloud and delete.
Definition: Cloud.C:233
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:80
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:225
void info() const
Print cloud information.
Definition: DSMCCloud.C:996
Generic GeometricField class.
Templated base class for spray cloud.
Definition: SprayCloud.H:75
virtual ~SprayCloud()
Destructor.
Definition: SprayCloud.C:135
void setModels()
Set cloud sub-models.
Definition: SprayCloud.C:35
void storeState()
Store the current cloud state.
Definition: SprayCloud.C:187
void cloudReset(SprayCloud< CloudType > &c)
Reset state of cloud.
Definition: SprayCloud.C:59
void evolve()
Evolve the spray (inject, move)
Definition: SprayCloud.C:208
SprayCloud(const word &cloudName, const volScalarField &rho, const volVectorField &U, const dimensionedVector &g, const fluidThermo &carrierThermo, bool readFields=true)
Construct given carrier gas fields.
Definition: SprayCloud.C:74
void info()
Print cloud information.
Definition: SprayCloud.C:220
void restoreState()
Reset the current cloud to the previously stored state.
Definition: SprayCloud.C:200
void checkParcelProperties(parcelType &parcel, const label injectori)
Check parcel properties.
Definition: SprayCloud.C:165
void setParcelThermoProperties(parcelType &parcel)
Set parcel thermo properties.
Definition: SprayCloud.C:143
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:57
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
scalar mu(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture viscosity [Pa s].
scalar sigma(const scalar p, const scalar T, const scalarField &X) const
Estimate mixture surface tension [N/m].
scalarField X(const scalarField &Y) const
Returns the mole fractions corresponding to the given mass fractions.
scalar Cp(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture heat capacity [J/kg/K].
scalar rho(const scalar p, const scalar T, const scalarField &X) const
Calculate the mixture density [kg/m^3].
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:51
A class for handling words, derived from string.
Definition: word.H:62
U
Definition: pEqn.H:72
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
const dimensionedScalar c
Speed of light in a vacuum.
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
messageStream Info
T clone(const T &t)
Definition: List.H:55
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
static const char nl
Definition: Ostream.H:266
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
volScalarField & p
PtrList< volScalarField > & Y
const word cloudName(propsDict.lookup("cloudName"))