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-2016 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 
29 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
30 
31 template<class CloudType>
33 {
34  collisionModel_.reset
35  (
37  (
38  this->subModelProperties(),
39  *this
40  ).ptr()
41  );
42 }
43 
44 
45 template<class CloudType>
46 template<class TrackData>
48 (
49  TrackData& td,
50  const scalar deltaT
51 )
52 {
53  td.part() = TrackData::tpVelocityHalfStep;
54  CloudType::move(td, deltaT);
55 
56  td.part() = TrackData::tpLinearTrack;
57  CloudType::move(td, deltaT);
58 
59  // td.part() = TrackData::tpRotationalTrack;
60  // CloudType::move(td);
61 
62  this->updateCellOccupancy();
63 
64  this->collision().collide();
65 
66  td.part() = TrackData::tpVelocityHalfStep;
67  CloudType::move(td, deltaT);
68 }
69 
70 
71 
72 template<class CloudType>
74 {
75  CloudType::cloudReset(c);
76 
77  collisionModel_.reset(c.collisionModel_.ptr());
78 }
79 
80 
81 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
82 
83 template<class CloudType>
85 (
86  const word& cloudName,
87  const volScalarField& rho,
88  const volVectorField& U,
89  const volScalarField& mu,
90  const dimensionedVector& g,
91  bool readFields
92 )
93 :
94  CloudType(cloudName, rho, U, mu, g, false),
95  constProps_(this->particleProperties()),
96  collisionModel_(NULL)
97 {
98  if (this->solution().steadyState())
99  {
101  << "Collision modelling not currently available for steady state "
102  << "calculations" << exit(FatalError);
103  }
104 
105  if (this->solution().active())
106  {
107  setModels();
108 
109  if (readFields)
110  {
111  parcelType::readFields(*this);
112  }
113  }
114 }
115 
116 
117 template<class CloudType>
119 (
121  const word& name
122 )
123 :
124  CloudType(c, name),
125  collisionModel_(c.collisionModel_->clone())
126 {}
127 
128 
129 template<class CloudType>
131 (
132  const fvMesh& mesh,
133  const word& name,
135 )
136 :
137  CloudType(mesh, name, c),
138  collisionModel_(NULL)
139 {}
140 
141 
142 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
143 
144 template<class CloudType>
146 {}
147 
148 
149 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
150 
151 template<class CloudType>
153 {
154  return !collision().controlsWallInteraction();
155 }
156 
157 
158 template<class CloudType>
160 {
161  cloudCopyPtr_.reset
162  (
163  static_cast<CollidingCloud<CloudType>*>
164  (
165  clone(this->name() + "Copy").ptr()
166  )
167  );
168 }
169 
170 
171 template<class CloudType>
173 {
174  cloudReset(cloudCopyPtr_());
175  cloudCopyPtr_.clear();
176 }
177 
178 
179 template<class CloudType>
181 {
182  if (this->solution().canEvolve())
183  {
184  typename parcelType::template
185  TrackingData<CollidingCloud<CloudType>> td(*this);
186 
187  this->solve(td);
188  }
189 }
190 
191 
192 template<class CloudType>
193 template<class TrackData>
195 {
196  // Sympletic leapfrog integration of particle forces:
197  // + apply half deltaV with stored force
198  // + move positions with new velocity
199  // + calculate forces in new position
200  // + apply half deltaV with new force
201 
202  label nSubCycles = collision().nSubCycles();
203 
204  if (nSubCycles > 1)
205  {
206  Info<< " " << nSubCycles << " move-collide subCycles" << endl;
207 
208  subCycleTime moveCollideSubCycle
209  (
210  const_cast<Time&>(this->db().time()),
211  nSubCycles
212  );
213 
214  while(!(++moveCollideSubCycle).end())
215  {
216  moveCollide(td, this->db().time().deltaTValue());
217  }
218 
219  moveCollideSubCycle.endSubCycle();
220  }
221  else
222  {
223  moveCollide(td, this->db().time().deltaTValue());
224  }
225 }
226 
227 
228 template<class CloudType>
230 {
231  CloudType::info();
232 
233  scalar rotationalKineticEnergy = rotationalKineticEnergyOfSystem();
234  reduce(rotationalKineticEnergy, sumOp<scalar>());
235 
236  Info<< " Rotational kinetic energy = "
237  << rotationalKineticEnergy << nl;
238 }
239 
240 
241 // ************************************************************************* //
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
error FatalError
virtual bool hasWallImpactDistance() const
If the collision model controls the wall interaction,.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
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.
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
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
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
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.