solidParticle.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-2018 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 "solidParticleCloud.H"
27 
28 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32  defineTemplateTypeNameAndDebug(Cloud<solidParticle>, 0);
33 }
34 
35 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
36 
38 (
39  solidParticleCloud& cloud,
40  trackingData& td,
41  const scalar trackTime
42 )
43 {
44  td.switchProcessor = false;
45  td.keepParticle = true;
46 
47  while (td.keepParticle && !td.switchProcessor && stepFraction() < 1)
48  {
49  if (debug)
50  {
51  Info<< "Time = " << mesh().time().timeName()
52  << " trackTime = " << trackTime
53  << " steptFraction() = " << stepFraction() << endl;
54  }
55 
56 
57  const scalar sfrac = stepFraction();
58 
59  const scalar f = 1 - stepFraction();
60  trackToAndHitFace(f*trackTime*U_, f, cloud, td);
61 
62  const scalar dt = (stepFraction() - sfrac)*trackTime;
63 
64  const tetIndices tetIs = this->currentTetIndices();
65  scalar rhoc = td.rhoInterp().interpolate(this->coordinates(), tetIs);
66  vector Uc = td.UInterp().interpolate(this->coordinates(), tetIs);
67  scalar nuc = td.nuInterp().interpolate(this->coordinates(), tetIs);
68 
69  scalar rhop = cloud.rhop();
70  scalar magUr = mag(Uc - U_);
71 
72  scalar ReFunc = 1.0;
73  scalar Re = magUr*d_/nuc;
74 
75  if (Re > 0.01)
76  {
77  ReFunc += 0.15*pow(Re, 0.687);
78  }
79 
80  scalar Dc = (24.0*nuc/d_)*ReFunc*(3.0/4.0)*(rhoc/(d_*rhop));
81 
82  U_ = (U_ + dt*(Dc*Uc + (1.0 - rhoc/rhop)*td.g()))/(1.0 + dt*Dc);
83  }
84 
85  return td.keepParticle;
86 }
87 
88 
90 {
91  return false;
92 }
93 
94 
96 (
98  trackingData& td
99 )
100 {
101  td.switchProcessor = true;
102 }
103 
104 
106 {
107  const vector nw = normal();
108 
109  const scalar Un = U_ & nw;
110  const vector Ut = U_ - Un*nw;
111 
112  if (Un > 0)
113  {
114  U_ -= (1.0 + cloud.e())*Un*nw;
115  }
116 
117  U_ -= cloud.mu()*Ut;
118 }
119 
120 
122 {
124  U_ = transform(T, U_);
125 }
126 
127 
129 {
130  particle::transformProperties(separation);
131 }
132 
133 
134 // ************************************************************************* //
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
bool hitPatch(solidParticleCloud &cloud, trackingData &td)
Overridable function to handle the particle hitting a patch.
Definition: solidParticle.C:89
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
virtual void transformProperties(const tensor &T)
Transform the physical properties of the particle.
Definition: particle.C:972
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
dynamicFvMesh & mesh
bool switchProcessor
Flag to switch processor.
Definition: particle.H:116
Type interpolate(const cellPointWeight &cpw) const
Interpolate field for the given cellPointWeight.
A Cloud of solid particles.
const interpolationCellPoint< scalar > & rhoInterp() const
const interpolationCellPoint< scalar > & nuInterp() const
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
Definition: tetIndices.H:81
defineTemplateTypeNameAndDebug(IOPtrList< ensightPart >, 0)
const interpolationCellPoint< vector > & UInterp() const
Class used to pass tracking data to the trackToFace function.
Definition: solidParticle.H:86
labelList f(nPoints)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void hitProcessorPatch(solidParticleCloud &cloud, trackingData &td)
Overridable function to handle the particle hitting a.
Definition: solidParticle.C:96
PtrList< coordinateSystem > coordinates(solidRegions.size())
bool keepParticle
Flag to indicate whether to keep particle (false = delete)
Definition: particle.H:119
bool move(solidParticleCloud &, trackingData &, const scalar)
Move.
Definition: solidParticle.C:38
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void transformProperties(const tensor &T)
Transform the physical properties of the particle.
void hitWallPatch(solidParticleCloud &cloud, trackingData &td)
Overridable function to handle the particle hitting a wallPatch.
scalarField Re(const UList< complex > &cf)
Definition: complexFields.C:97
label nw
Definition: createFields.H:25
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477