patchEjection.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) 2015-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 "patchEjection.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace regionModels
34 {
35 namespace surfaceFilmModels
36 {
37 
38 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 
40 defineTypeNameAndDebug(patchEjection, 0);
41 addToRunTimeSelectionTable(ejectionModel, patchEjection, dictionary);
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
46 (
48  const dictionary& dict
49 )
50 :
51  ejectionModel(type(), film, dict),
52  deltaStable_(coeffDict_.lookupOrDefault<scalar>("deltaStable", 0.0))
53 {
54  const polyBoundaryMesh& pbm = film.regionMesh().boundaryMesh();
55  patchIDs_.setSize
56  (
57  pbm.size() - film.regionMesh().globalData().processorPatches().size()
58  );
59 
60  if (coeffDict_.found("patches"))
61  {
62  wordReList patchNames(coeffDict_.lookup("patches"));
63  const labelHashSet patchSet = pbm.patchSet(patchNames);
64 
65  Info<< " applying to patches:" << nl;
66 
67  label pidi = 0;
68  forAllConstIter(labelHashSet, patchSet, iter)
69  {
70  label patchi = iter.key();
71  patchIDs_[pidi++] = patchi;
72  Info<< " " << pbm[patchi].name() << endl;
73  }
74  patchIDs_.setSize(pidi);
75  patchEjectedMasses_.setSize(pidi, 0);
76  }
77  else
78  {
79  Info<< " applying to all patches" << endl;
80 
81  forAll(patchIDs_, patchi)
82  {
83  patchIDs_[patchi] = patchi;
84  }
85 
86  patchEjectedMasses_.setSize(patchIDs_.size(), 0);
87  }
88 
89  if (!patchIDs_.size())
90  {
92  << "No patches selected"
93  << exit(FatalError);
94  }
95 }
96 
97 
98 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
99 
101 {}
102 
103 
104 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
105 
107 (
108  scalarField& availableMass,
109  scalarField& massToEject,
110  scalarField& diameterToEject
111 )
112 {
113  // Do not correct if no patches selected
114  if (!patchIDs_.size()) return;
115 
116  const scalarField& delta = film().delta();
117  const scalarField& rho = film().rho();
118  const scalarField& magSf = film().magSf();
119 
120  const polyBoundaryMesh& pbm = film().regionMesh().boundaryMesh();
121 
122  forAll(patchIDs_, pidi)
123  {
124  label patchi = patchIDs_[pidi];
125  const polyPatch& pp = pbm[patchi];
126  const labelList& faceCells = pp.faceCells();
127 
128  // Accumulate the total mass removed from patch
129  scalar dMassPatch = 0;
130 
131  forAll(faceCells, fci)
132  {
133  label celli = faceCells[fci];
134 
135  scalar ddelta = max(0.0, delta[celli] - deltaStable_);
136  scalar dMass = ddelta*rho[celli]*magSf[celli];
137  massToEject[celli] += dMass;
138  availableMass[celli] -= dMass;
139  dMassPatch += dMass;
140  }
141 
142  patchEjectedMasses_[pidi] += dMassPatch;
143  addToEjectedMass(dMassPatch);
144  }
145 
147 
148  if (writeTime())
149  {
150  scalarField patchEjectedMasses0
151  (
152  getModelProperty<scalarField>
153  (
154  "patchEjectedMasses",
156  )
157  );
158 
160  Pstream::listCombineGather(patchEjectedMassTotals, plusEqOp<scalar>());
161  patchEjectedMasses0 += patchEjectedMassTotals;
162 
163  setModelProperty<scalarField>
164  (
165  "patchEjectedMasses",
166  patchEjectedMasses0
167  );
168 
170  }
171 }
172 
173 
175 {
176  // Do not correct if no patches selected
177  if (!patchIDs_.size()) return;
178 
179  scalarField patchEjectedMasses
180  (
181  getModelProperty<scalarField>
182  (
183  "patchEjectedMasses",
185  )
186  );
187 
189  Pstream::listCombineGather(patchEjectedMassTotals, plusEqOp<scalar>());
190 
191  forAll(patchIDs_, pidi)
192  {
193  label patchi = patchIDs_[pidi];
194  patchMasses[patchi] +=
195  patchEjectedMasses[pidi] + patchEjectedMassTotals[pidi];
196  }
197 }
198 
199 
200 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
201 
202 } // End namespace surfaceFilmModels
203 } // End namespace regionModels
204 } // End namespace Foam
205 
206 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
scalar delta
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const labelList & processorPatches() const
Return list of processor patch labels.
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
const word & name() const
Return name.
Definition: IOobject.H:303
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
Base class for film ejection models, handling mass transfer from the film.
Definition: ejectionModel.H:56
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch IDs corresponding to the given names.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
virtual const volScalarField & rho() const =0
Return the film density [kg/m^3].
Macros for easy insertion into run-time selection tables.
const surfaceFilmRegionModel & film() const
Return const access to the film surface film model.
virtual bool writeTime() const
Flag to indicate when to write a property.
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:334
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
wordList patchNames(nPatches)
virtual const volScalarField & delta() const =0
Return the film thickness [m].
scalarField patchEjectedMasses_
Ejected mass for each patch at which the film is removed.
Definition: patchEjection.H:67
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:55
patchEjection(surfaceFilmRegionModel &film, const dictionary &dict)
Construct from surface film model.
Definition: patchEjection.C:46
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1394
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
Foam::polyBoundaryMesh.
void addToEjectedMass(const scalar dMass)
Add to ejected mass.
Definition: ejectionModel.C:44
static const char nl
Definition: Ostream.H:260
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label patchi
const volScalarField::Internal & magSf() const
Return the face area magnitudes [m^2].
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
messageStream Info
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
virtual void patchEjectedMassTotals(scalarField &patchMasses) const
Accumulate the total mass ejected for the patches into the.
scalar deltaStable_
Stable film thickness - mass only removed if thickness exceeds.
Definition: patchEjection.H:61
defineTypeNameAndDebug(kinematicSingleLayer, 0)
addToRunTimeSelectionTable(surfaceFilmModel, noFilm, mesh)
Namespace for OpenFOAM.
labelList patchIDs_
List of patch IDs at which the film is removed.
Definition: patchEjection.H:64