CollidingCloud.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 "CollidingCloud.H"
27 #include "CollisionModel.H"
28 #include "NoCollision.H"
29 
30 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
31 
32 template<class CloudType>
34 {
35  collisionModel_.reset
36  (
38  (
39  this->subModelProperties(),
40  *this
41  ).ptr()
42  );
43 }
44 
45 
46 template<class CloudType>
47 template<class TrackData>
49 (
50  TrackData& td,
51  const scalar deltaT
52 )
53 {
54  td.part() = TrackData::tpVelocityHalfStep;
55  CloudType::move(td, deltaT);
56 
57  td.part() = TrackData::tpLinearTrack;
58  CloudType::move(td, deltaT);
59 
60  // td.part() = TrackData::tpRotationalTrack;
61  // CloudType::move(td);
62 
63  this->updateCellOccupancy();
64 
65  this->collision().collide();
66 
67  td.part() = TrackData::tpVelocityHalfStep;
68  CloudType::move(td, deltaT);
69 }
70 
71 
72 
73 template<class CloudType>
75 {
76  CloudType::cloudReset(c);
77 
78  collisionModel_.reset(c.collisionModel_.ptr());
79 }
80 
81 
82 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
83 
84 template<class CloudType>
86 (
87  const word& cloudName,
88  const volScalarField& rho,
89  const volVectorField& U,
90  const volScalarField& mu,
91  const dimensionedVector& g,
92  bool readFields
93 )
94 :
95  CloudType(cloudName, rho, U, mu, g, false),
96  constProps_(this->particleProperties()),
97  collisionModel_(nullptr)
98 {
99  if (this->solution().active())
100  {
101  setModels();
102 
103  if (readFields)
104  {
105  parcelType::readFields(*this);
106  this->deleteLostParticles();
107  }
108 
109  if
110  (
111  this->solution().steadyState()
112  && !isType<NoCollision<CollidingCloud<CloudType>>>(collisionModel_())
113  )
114  {
116  << "Collision modelling not currently available "
117  << "for steady state calculations" << exit(FatalError);
118  }
119  }
120 }
121 
122 
123 template<class CloudType>
125 (
127  const word& name
128 )
129 :
130  CloudType(c, name),
131  collisionModel_(c.collisionModel_->clone())
132 {}
133 
134 
135 template<class CloudType>
137 (
138  const fvMesh& mesh,
139  const word& name,
141 )
142 :
143  CloudType(mesh, name, c),
144  collisionModel_(nullptr)
145 {}
146 
147 
148 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
149 
150 template<class CloudType>
152 {}
153 
154 
155 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
156 
157 template<class CloudType>
159 {
160  return !collision().controlsWallInteraction();
161 }
162 
163 
164 template<class CloudType>
166 {
167  cloudCopyPtr_.reset
168  (
169  static_cast<CollidingCloud<CloudType>*>
170  (
171  clone(this->name() + "Copy").ptr()
172  )
173  );
174 }
175 
176 
177 template<class CloudType>
179 {
180  cloudReset(cloudCopyPtr_());
181  cloudCopyPtr_.clear();
182 }
183 
184 
185 template<class CloudType>
187 {
188  if (this->solution().canEvolve())
189  {
190  typename parcelType::template
191  TrackingData<CollidingCloud<CloudType>> td(*this);
192 
193  this->solve(td);
194  }
195 }
196 
197 
198 template<class CloudType>
199 template<class TrackData>
201 {
202  // Sympletic leapfrog integration of particle forces:
203  // + apply half deltaV with stored force
204  // + move positions with new velocity
205  // + calculate forces in new position
206  // + apply half deltaV with new force
207 
208  label nSubCycles = collision().nSubCycles();
209 
210  if (nSubCycles > 1)
211  {
212  Info<< " " << nSubCycles << " move-collide subCycles" << endl;
213 
214  subCycleTime moveCollideSubCycle
215  (
216  const_cast<Time&>(this->db().time()),
217  nSubCycles
218  );
219 
220  while(!(++moveCollideSubCycle).end())
221  {
222  moveCollide(td, this->db().time().deltaTValue());
223  }
224 
225  moveCollideSubCycle.endSubCycle();
226  }
227  else
228  {
229  moveCollide(td, this->db().time().deltaTValue());
230  }
231 }
232 
233 
234 template<class CloudType>
236 {
237  CloudType::info();
238 
239  scalar rotationalKineticEnergy = rotationalKineticEnergyOfSystem();
240  reduce(rotationalKineticEnergy, sumOp<scalar>());
241 
242  Info<< " Rotational kinetic energy = "
243  << rotationalKineticEnergy << nl;
244 }
245 
246 
247 // ************************************************************************* //
void cloudReset(CollidingCloud< CloudType > &c)
Reset state of cloud.
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
DSMCCloud< dsmcParcel > CloudType
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject *> &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
bool hasWallImpactDistance() const
If the collision model controls the wall interaction,.
void evolve()
Evolve the cloud.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
void moveCollide(TrackData &td, const scalar deltaT)
Move-collide particles.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
A class for handling words, derived from string.
Definition: word.H:59
bool isType(const Type &t)
Check the typeid.
Definition: typeInfo.H:126
rhoEqn solve()
Templated collision model class.
void endSubCycle()
End the sub-cycling and reset the time-state.
Definition: subCycleTime.C:56
Adds coolisions to kinematic clouds.
void motion(TrackData &td)
Particle motion.
static const char nl
Definition: Ostream.H:262
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Place holder for &#39;none&#39; option.
Definition: NoCollision.H:49
void restoreState()
Reset the current cloud to the previously stored state.
void info()
Print cloud information.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
messageStream Info
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
void setModels()
Set cloud sub-models.
virtual ~CollidingCloud()
Destructor.
A class for managing sub-cycling times.
Definition: subCycleTime.H:48
autoPtr< CollisionModel< CollidingCloud< CloudType > > > collisionModel_
Collision model.
void storeState()
Store the current cloud state.