ParticleErosion.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 "ParticleErosion.H"
27 #include "mathematicalConstants.H"
28 
29 using namespace Foam::constant;
30 
31 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
32 
33 template<class CloudType>
35 (
36  const label globalPatchi
37 ) const
38 {
39  forAll(patchIndices_, i)
40  {
41  if (patchIndices_[i] == globalPatchi)
42  {
43  return i;
44  }
45  }
46 
47  return -1;
48 }
49 
50 
51 template<class CloudType>
53 {
54  if (QPtr_.valid())
55  {
56  QPtr_->write();
57  }
58  else
59  {
61  << "QPtr not valid" << abort(FatalError);
62  }
63 }
64 
65 
66 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
67 
68 template<class CloudType>
70 (
71  const dictionary& dict,
72  CloudType& owner,
73  const word& modelName
74 )
75 :
76  CloudFunctionObject<CloudType>(dict, owner, modelName, typeName),
77  QPtr_(nullptr),
78  patchIndices_(),
79  p_(this->coeffDict().template lookup<scalar>("p")),
80  psi_(this->coeffDict().template lookupOrDefault<scalar>("psi", 2.0)),
81  K_(this->coeffDict().template lookupOrDefault<scalar>("K", 2.0))
82 {
83  const wordList allPatchNames = owner.mesh().boundaryMesh().names();
84  wordList patchName(this->coeffDict().lookup("patches"));
85 
86  labelHashSet uniquePatchIDs;
87  forAllReverse(patchName, i)
88  {
89  labelList patchIDs = findStrings(patchName[i], allPatchNames);
90 
91  if (patchIDs.empty())
92  {
94  << "Cannot find any patch names matching " << patchName[i]
95  << endl;
96  }
97 
98  uniquePatchIDs.insert(patchIDs);
99  }
100 
101  patchIndices_ = uniquePatchIDs.toc();
102 
103  // Trigger creation of the Q field
104  preEvolve();
105 }
106 
107 
108 template<class CloudType>
110 (
112 )
113 :
115  QPtr_(nullptr),
116  patchIndices_(pe.patchIndices_),
117  p_(pe.p_),
118  psi_(pe.psi_),
119  K_(pe.K_)
120 {}
121 
122 
123 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
124 
125 template<class CloudType>
127 {}
128 
129 
130 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
131 
132 template<class CloudType>
134 {
135  if (QPtr_.valid())
136  {
137  QPtr_->primitiveFieldRef() = 0.0;
138  }
139  else
140  {
141  const fvMesh& mesh = this->owner().mesh();
142 
143  QPtr_.reset
144  (
145  new volScalarField
146  (
147  IOobject
148  (
149  this->owner().name() + ":Q",
150  mesh.time().name(),
151  mesh,
154  ),
155  mesh,
157  )
158  );
159  }
160 }
161 
162 
163 template<class CloudType>
165 {
166  const fvMesh& mesh = this->owner().mesh();
167  if (!p.onBoundaryFace(mesh)) return;
168 
169  const polyPatch& pp = mesh.boundaryMesh()[p.patch(mesh)];
170  const label patchi = pp.index();
171  const label localPatchi = applyToPatch(patchi);
172 
173  if (localPatchi != -1)
174  {
175  // Get patch data
176  vector nw, Up;
177  this->owner().patchData(p, pp, nw, Up);
178 
179  // Particle velocity relative to patch
180  const vector& U = p.U() - Up;
181 
182  // Quick rejection if the particle is travelling away from the patch
183  if ((nw & U) < 0)
184  {
185  return;
186  }
187 
188  const scalar magU = mag(U);
189  const vector UHat = U/magU;
190 
191  // Impact angle
192  const scalar alpha = mathematical::pi/2 - acos(nw & UHat);
193 
194  // Get the face value to accumulate into
195  const label patchFacei = pp.whichFace(p.face());
196  scalar& Q = QPtr_->boundaryFieldRef()[patchi][patchFacei];
197 
198  // Finnie's model
199  const scalar coeff = p.nParticle()*p.mass()*sqr(magU)/(p_*psi_*K_);
200  if (tan(alpha) < K_/6)
201  {
202  Q += coeff*(sin(2*alpha) - 6/K_*sqr(sin(alpha)));
203  }
204  else
205  {
206  Q += coeff*(K_*sqr(cos(alpha))/6);
207  }
208  }
209 }
210 
211 
212 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllReverse(list, i)
Reverse loop across all elements in list.
Definition: UList.H:446
Templated cloud function object base class.
const CloudType & owner() const
Return const access to the owner cloud.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:80
const fvMesh & mesh() const
Return references to the mesh.
Definition: DSMCCloudI.H:41
Generic GeometricField class.
bool insert(const Key &key)
Insert a new entry.
Definition: HashSet.H:109
List< Key > toc() const
Return the table of contents.
Definition: HashTable.C:227
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
Function object to create a field of eroded volume, Q, on a specified list of patches....
virtual ~ParticleErosion()
Destructor.
virtual void preFace(const parcelType &p)
Pre-face hook.
label applyToPatch(const label globalPatchi) const
Returns local patchi if patch is in patchIds_ list.
virtual void write()
Write post-processing info.
virtual void preEvolve()
Pre-evolve hook.
ParticleErosion(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
bool empty() const
Return true if the UList is empty (ie, size() is zero)
Definition: UListI.H:325
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
const word & name() const
Return const reference to name.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:418
const polyMesh & mesh() const
Return reference to polyMesh.
Definition: fvMesh.H:441
label index() const
Return the index of this patch in the boundaryMesh.
wordList names() const
Return the list of patch names.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:404
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
label whichFace(const label l) const
Return label of face in patch from global face label.
Definition: polyPatch.H:360
void reset(const label nPoints, const label nInternalFaces, const label nFaces, const label nCells)
Reset this primitiveMesh given the primitive array sizes.
const dictionary & coeffDict() const
Return const access to the coefficients dictionary.
Definition: subModelBase.C:128
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label patchi
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
#define WarningInFunction
Report a warning using Foam::Warning.
Collection of constants.
dimensionedScalar tan(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
bool findStrings(const wordReListMatcher &matcher, const std::string &str)
Return true if string matches one of the regular expressions.
Definition: stringListOps.H:52
dimensionedSymmTensor sqr(const dimensionedVector &dv)
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
errorManip< error > abort(error &err)
Definition: errorManip.H:131
dimensionedScalar sin(const dimensionedScalar &ds)
const dimensionSet dimVolume
dimensioned< scalar > mag(const dimensioned< Type > &)
error FatalError
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar acos(const dimensionedScalar &ds)
dictionary dict
volScalarField & p