SurfaceFilmModel.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 "SurfaceFilmModel.H"
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
32 template<class CloudType>
34 :
36  g_(owner.g()),
37  ejectedParcelType_(0),
38  massParcelPatch_(0),
39  diameterParcelPatch_(0),
40  deltaFilmPatch_(0),
41  nParcelsTransferred_(0),
42  nParcelsInjected_(0)
43 {}
44 
45 
46 template<class CloudType>
48 (
49  const dictionary& dict,
50  CloudType& owner,
51  const word& type
52 )
53 :
54  CloudSubModelBase<CloudType>(owner, dict, typeName, type),
55  g_(owner.g()),
56  ejectedParcelType_
57  (
58  this->coeffDict().lookupOrDefault("ejectedParcelType", -1)
59  ),
60  massParcelPatch_(0),
61  diameterParcelPatch_(0),
62  deltaFilmPatch_(),
63  nParcelsTransferred_(0),
64  nParcelsInjected_(0)
65 {}
66 
67 
68 template<class CloudType>
70 (
72 )
73 :
75  g_(sfm.g_),
76  ejectedParcelType_(sfm.ejectedParcelType_),
77  massParcelPatch_(sfm.massParcelPatch_),
78  diameterParcelPatch_(sfm.diameterParcelPatch_),
79  deltaFilmPatch_(sfm.deltaFilmPatch_),
80  nParcelsTransferred_(sfm.nParcelsTransferred_),
81  nParcelsInjected_(sfm.nParcelsInjected_)
82 {}
83 
84 
85 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
86 
87 template<class CloudType>
89 {}
90 
91 
92 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
93 
94 template<class CloudType>
95 template<class TrackCloudType>
97 {
98  const labelList& filmPatches = this->filmPatches();
99 
100  forAll(filmPatches, filmi)
101  {
102  const label filmPatchi = filmPatches[filmi];
103 
104  const fvMesh& mesh = this->owner().mesh();
105  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
106 
107  const labelList& injectorCellsPatch = pbm[filmPatchi].faceCells();
108 
109  // Get and cache the properties of the droplets ejected from the film
110  cacheFilmFields(filmi);
111 
112  const vectorField& Cf = mesh.C().boundaryField()[filmPatchi];
113  const vectorField& Sf = mesh.Sf().boundaryField()[filmPatchi];
114  const scalarField& magSf = mesh.magSf().boundaryField()[filmPatchi];
115 
116  if (massParcelPatch_.size())
117  {
118  forAll(injectorCellsPatch, j)
119  {
120  if (massParcelPatch_[j] > 0)
121  {
122  const label celli = injectorCellsPatch[j];
123 
124  const scalar offset = max
125  (
126  diameterParcelPatch_[j],
127  deltaFilmPatch_[j]
128  );
129 
130  const point pos = Cf[j] - 1.1*offset*Sf[j]/magSf[j];
131 
132  // Create a new parcel
133  parcelType* pPtr =
134  new parcelType(this->owner().pMesh(), pos, celli);
135 
136  // Check/set new parcel thermo properties
137  cloud.setParcelThermoProperties(*pPtr);
138 
139  setParcelProperties(*pPtr, j);
140 
141  if (pPtr->nParticle() > 0.001)
142  {
143  // Check new parcel properties
144  cloud.checkParcelProperties(*pPtr, -1);
145 
146  // Add the new parcel to the cloud
147  cloud.addParticle(pPtr);
148 
149  nParcelsInjected_++;
150  }
151  else
152  {
153  // TODO: cache mass and re-distribute?
154  delete pPtr;
155  }
156  }
157  }
158  }
159  }
160 }
161 
162 
163 template<class CloudType>
165 {
166  label nTrans0 =
167  this->template getModelProperty<label>("nParcelsTransferred");
168 
169  label nInject0 =
170  this->template getModelProperty<label>("nParcelsInjected");
171 
172  label nTransTotal =
173  nTrans0 + returnReduce(nParcelsTransferred_, sumOp<label>());
174 
175  label nInjectTotal =
176  nInject0 + returnReduce(nParcelsInjected_, sumOp<label>());
177 
178  os << " Parcels absorbed into film = " << nTransTotal << nl
179  << " New film detached parcels = " << nInjectTotal << endl;
180 
181  if (this->writeTime())
182  {
183  this->setModelProperty("nParcelsTransferred", nTransTotal);
184  this->setModelProperty("nParcelsInjected", nInjectTotal);
185  nParcelsTransferred_ = 0;
186  nParcelsInjected_ = 0;
187  }
188 }
189 
190 
191 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
192 
193 #include "SurfaceFilmModelNew.C"
194 
195 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Base class for cloud sub-models.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:79
const Boundary & boundaryField() const
Return const-reference to the boundary field.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
Templated wall surface film model class.
void inject(TrackCloudType &cloud)
Inject parcels into the cloud.
virtual void info(Ostream &os)
Write surface film info to stream.
virtual ~SurfaceFilmModel()
Destructor.
SurfaceFilmModel(CloudType &owner)
Construct null from owner.
CloudType::parcelType parcelType
Convenience typedef to the cloud's parcel type.
A cloud is a collection of lagrangian particles.
Definition: cloud.H:55
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:101
const volVectorField & C() const
Return cell centres.
const surfaceVectorField & Sf() const
Return cell face area vectors.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
Foam::polyBoundaryMesh.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:403
A class for handling words, derived from string.
Definition: word.H:62
dimensionedScalar pos(const dimensionedScalar &ds)
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:251
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void offset(label &lst, const label o)
static const char nl
Definition: Ostream.H:260
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dictionary dict
Foam::surfaceFields.