SurfaceFilmModel.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2017 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 "surfaceFilmModel.H"
28 #include "mathematicalConstants.H"
29 
30 using namespace Foam::constant;
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
34 template<class CloudType>
36 :
38  g_(owner.g()),
39  ejectedParcelType_(0),
40  massParcelPatch_(0),
41  diameterParcelPatch_(0),
42  UFilmPatch_(0),
43  rhoFilmPatch_(0),
44  deltaFilmPatch_(0),
45  nParcelsTransferred_(0),
46  nParcelsInjected_(0)
47 {}
48 
49 
50 template<class CloudType>
52 (
53  const dictionary& dict,
54  CloudType& owner,
55  const word& type
56 )
57 :
58  CloudSubModelBase<CloudType>(owner, dict, typeName, type),
59  g_(owner.g()),
60  ejectedParcelType_
61  (
62  this->coeffDict().lookupOrDefault("ejectedParcelType", -1)
63  ),
64  massParcelPatch_(0),
65  diameterParcelPatch_(0),
66  UFilmPatch_(0),
67  rhoFilmPatch_(0),
68  deltaFilmPatch_(owner.mesh().boundary().size()),
69  nParcelsTransferred_(0),
70  nParcelsInjected_(0)
71 {}
72 
73 
74 template<class CloudType>
76 (
78 )
79 :
81  g_(sfm.g_),
82  ejectedParcelType_(sfm.ejectedParcelType_),
83  massParcelPatch_(sfm.massParcelPatch_),
84  diameterParcelPatch_(sfm.diameterParcelPatch_),
85  UFilmPatch_(sfm.UFilmPatch_),
86  rhoFilmPatch_(sfm.rhoFilmPatch_),
87  deltaFilmPatch_(sfm.deltaFilmPatch_),
88  nParcelsTransferred_(sfm.nParcelsTransferred_),
89  nParcelsInjected_(sfm.nParcelsInjected_)
90 {}
91 
92 
93 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
94 
95 template<class CloudType>
97 {}
98 
99 
100 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
101 
102 template<class CloudType>
103 template<class TrackData>
105 {
106  if (!this->active())
107  {
108  return;
109  }
110 
111  // Retrieve the film model from the owner database
113  this->owner().mesh().time().objectRegistry::template lookupObject
115  (
116  "surfaceFilmProperties"
117  );
118 
119  if (!filmModel.active())
120  {
121  return;
122  }
123 
124  const labelList& filmPatches = filmModel.intCoupledPatchIDs();
125  const labelList& primaryPatches = filmModel.primaryPatchIDs();
126 
127  const fvMesh& mesh = this->owner().mesh();
128  const polyBoundaryMesh& pbm = mesh.boundaryMesh();
129 
130  forAll(filmPatches, i)
131  {
132  const label filmPatchi = filmPatches[i];
133  const label primaryPatchi = primaryPatches[i];
134 
135  const labelList& injectorCellsPatch = pbm[primaryPatchi].faceCells();
136 
137  cacheFilmFields(filmPatchi, primaryPatchi, filmModel);
138 
139  const vectorField& Cf = mesh.C().boundaryField()[primaryPatchi];
140  const vectorField& Sf = mesh.Sf().boundaryField()[primaryPatchi];
141  const scalarField& magSf = mesh.magSf().boundaryField()[primaryPatchi];
142 
143  forAll(injectorCellsPatch, j)
144  {
145  if (diameterParcelPatch_[j] > 0)
146  {
147  const label celli = injectorCellsPatch[j];
148 
149  const scalar offset =
150  max
151  (
152  diameterParcelPatch_[j],
153  deltaFilmPatch_[primaryPatchi][j]
154  );
155  const point pos = Cf[j] - 1.1*offset*Sf[j]/magSf[j];
156 
157  // Create a new parcel
158  parcelType* pPtr =
159  new parcelType(this->owner().pMesh(), pos, celli);
160 
161  // Check/set new parcel thermo properties
162  td.cloud().setParcelThermoProperties(*pPtr, 0.0);
163 
164  setParcelProperties(*pPtr, j);
165 
166  if (pPtr->nParticle() > 0.001)
167  {
168  // Check new parcel properties
169  // td.cloud().checkParcelProperties(*pPtr, 0.0, true);
170  td.cloud().checkParcelProperties(*pPtr, 0.0, false);
171 
172  // Add the new parcel to the cloud
173  td.cloud().addParticle(pPtr);
174 
175  nParcelsInjected_++;
176  }
177  else
178  {
179  // TODO: cache mass and re-distribute?
180  delete pPtr;
181  }
182  }
183  }
184  }
185 }
186 
187 
188 template<class CloudType>
190 (
191  const label filmPatchi,
192  const label primaryPatchi,
194 )
195 {
196  massParcelPatch_ = filmModel.cloudMassTrans().boundaryField()[filmPatchi];
197  filmModel.toPrimary(filmPatchi, massParcelPatch_);
198 
199  diameterParcelPatch_ =
200  filmModel.cloudDiameterTrans().boundaryField()[filmPatchi];
201  filmModel.toPrimary(filmPatchi, diameterParcelPatch_, maxEqOp<scalar>());
202 
203  UFilmPatch_ = filmModel.Us().boundaryField()[filmPatchi];
204  filmModel.toPrimary(filmPatchi, UFilmPatch_);
205 
206  rhoFilmPatch_ = filmModel.rho().boundaryField()[filmPatchi];
207  filmModel.toPrimary(filmPatchi, rhoFilmPatch_);
208 
209  deltaFilmPatch_[primaryPatchi] =
210  filmModel.delta().boundaryField()[filmPatchi];
211  filmModel.toPrimary(filmPatchi, deltaFilmPatch_[primaryPatchi]);
212 }
213 
214 
215 template<class CloudType>
217 (
218  parcelType& p,
219  const label filmFacei
220 ) const
221 {
222  // Set parcel properties
223  scalar vol = mathematical::pi/6.0*pow3(diameterParcelPatch_[filmFacei]);
224  p.d() = diameterParcelPatch_[filmFacei];
225  p.U() = UFilmPatch_[filmFacei];
226  p.rho() = rhoFilmPatch_[filmFacei];
227 
228  p.nParticle() = massParcelPatch_[filmFacei]/p.rho()/vol;
229 
230  if (ejectedParcelType_ >= 0)
231  {
232  p.typeId() = ejectedParcelType_;
233  }
234 }
235 
236 
237 template<class CloudType>
239 {
240  label nTrans0 =
241  this->template getModelProperty<label>("nParcelsTransferred");
242 
243  label nInject0 =
244  this->template getModelProperty<label>("nParcelsInjected");
245 
246  label nTransTotal =
247  nTrans0 + returnReduce(nParcelsTransferred_, sumOp<label>());
248 
249  label nInjectTotal =
250  nInject0 + returnReduce(nParcelsInjected_, sumOp<label>());
251 
252  os << " Parcels absorbed into film = " << nTransTotal << nl
253  << " New film detached parcels = " << nInjectTotal << endl;
254 
255  if (this->writeTime())
256  {
257  this->setModelProperty("nParcelsTransferred", nTransTotal);
258  this->setModelProperty("nParcelsInjected", nInjectTotal);
259  nParcelsTransferred_ = 0;
260  nParcelsInjected_ = 0;
261  }
262 }
263 
264 
265 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
266 
267 #include "SurfaceFilmModelNew.C"
268 
269 // ************************************************************************* //
Collection of constants.
virtual void info(Ostream &os)
Write surface film info to stream.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const surfaceVectorField & Sf() const
Return cell face area vectors.
type
Types of root.
Definition: Roots.H:52
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
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
label nParcelsInjected_
Number of parcels injected from the film model.
Base class for cloud sub-models.
virtual void cacheFilmFields(const label filmPatchi, const label primaryPatchi, const regionModels::surfaceFilmModels::surfaceFilmModel &filmModel)
Cache the film fields in preparation for injection.
dimensionedScalar pos(const dimensionedScalar &ds)
void toPrimary(const label regionPatchi, List< Type > &regionField) const
Convert a local region field to the primary region.
scalarListList deltaFilmPatch_
Film height of all film patches / patch face.
dynamicFvMesh & mesh
void inject(TrackData &td)
Inject parcels into the cloud.
A class for handling words, derived from string.
Definition: word.H:59
virtual const volScalarField & cloudMassTrans() const =0
Return the film mass available for transfer.
virtual const volScalarField & delta() const =0
Return the film thickness [m].
Foam::polyBoundaryMesh.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
const labelList & primaryPatchIDs() const
Return the list of patch IDs on the primary region coupled.
Definition: regionModelI.H:172
static const char nl
Definition: Ostream.H:262
List< vector > UFilmPatch_
Film velocity / patch face.
const dimensionedVector & g
virtual const volVectorField & Us() const =0
Return the film surface velocity [m/s].
Templated wall surface film model class.
virtual void setParcelProperties(parcelType &p, const label filmFacei) const
Set the individual parcel properties.
virtual ~SurfaceFilmModel()
Destructor.
dimensionedScalar pow3(const dimensionedScalar &ds)
Foam::KinematicCloud< Cloud< basicKinematicCollidingParcel > > ::parcelType parcelType
Convenience typedef to the cloud&#39;s parcel type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual const volScalarField & cloudDiameterTrans() const =0
Return the parcel diameters originating from film.
scalarList diameterParcelPatch_
Parcel diameter / patch face.
label nParcelsTransferred_
Number of parcels transferred to the film model.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
scalarList rhoFilmPatch_
Film density / patch face.
SurfaceFilmModel(CloudType &owner)
Construct null from owner.
const volVectorField & C() const
Return cell centres as volVectorField.
const labelList & intCoupledPatchIDs() const
Return the list of patch IDs internally coupled with the.
Definition: regionModelI.H:179
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:69
label ejectedParcelType_
Ejected parcel type label - id assigned to identify parcel for.
const dimensionedVector & g_
Gravitational acceleration constant.
const Switch & active() const
Return the active flag.
Definition: regionModelI.H:43
scalarList massParcelPatch_
Parcel mass / patch face.
virtual const volScalarField & rho() const =0
Return the film density [kg/m3].