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-2013 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  (
102  "Foam::CollidingCloud<CloudType>::CollidingCloud"
103  "("
104  "const word&, "
105  "const volScalarField&, "
106  "const volVectorField&, "
107  "const volScalarField&, "
108  "const dimensionedVector&, "
109  "bool"
110  ")"
111  ) << "Collision modelling not currently available for steady state "
112  << "calculations" << exit(FatalError);
113  }
114 
115  if (this->solution().active())
116  {
117  setModels();
118 
119  if (readFields)
120  {
121  parcelType::readFields(*this);
122  }
123  }
124 }
125 
126 
127 template<class CloudType>
129 (
131  const word& name
132 )
133 :
134  CloudType(c, name),
135  collisionModel_(c.collisionModel_->clone())
136 {}
137 
138 
139 template<class CloudType>
141 (
142  const fvMesh& mesh,
143  const word& name,
145 )
146 :
147  CloudType(mesh, name, c),
148  collisionModel_(NULL)
149 {}
150 
151 
152 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
153 
154 template<class CloudType>
156 {}
157 
158 
159 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
160 
161 template<class CloudType>
163 {
164  return !collision().controlsWallInteraction();
165 }
166 
167 
168 template<class CloudType>
170 {
171  cloudCopyPtr_.reset
172  (
173  static_cast<CollidingCloud<CloudType>*>
174  (
175  clone(this->name() + "Copy").ptr()
176  )
177  );
178 }
179 
180 
181 template<class CloudType>
183 {
184  cloudReset(cloudCopyPtr_());
185  cloudCopyPtr_.clear();
186 }
187 
188 
189 template<class CloudType>
191 {
192  if (this->solution().canEvolve())
193  {
194  typename parcelType::template
195  TrackingData<CollidingCloud<CloudType> > td(*this);
196 
197  this->solve(td);
198  }
199 }
200 
201 
202 template<class CloudType>
203 template<class TrackData>
205 {
206  // Sympletic leapfrog integration of particle forces:
207  // + apply half deltaV with stored force
208  // + move positions with new velocity
209  // + calculate forces in new position
210  // + apply half deltaV with new force
211 
212  label nSubCycles = collision().nSubCycles();
213 
214  if (nSubCycles > 1)
215  {
216  Info<< " " << nSubCycles << " move-collide subCycles" << endl;
217 
218  subCycleTime moveCollideSubCycle
219  (
220  const_cast<Time&>(this->db().time()),
221  nSubCycles
222  );
223 
224  while(!(++moveCollideSubCycle).end())
225  {
226  moveCollide(td, this->db().time().deltaTValue());
227  }
228 
229  moveCollideSubCycle.endSubCycle();
230  }
231  else
232  {
233  moveCollide(td, this->db().time().deltaTValue());
234  }
235 }
236 
237 
238 template<class CloudType>
240 {
241  CloudType::info();
242 
243  scalar rotationalKineticEnergy = rotationalKineticEnergyOfSystem();
244  reduce(rotationalKineticEnergy, sumOp<scalar>());
245 
246  Info<< " Rotational kinetic energy = "
247  << rotationalKineticEnergy << nl;
248 }
249 
250 
251 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
Adds coolisions to kinematic clouds.
void restoreState()
Reset the current cloud to the previously stored state.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
void setModels()
Set cloud sub-models.
virtual bool hasWallImpactDistance() const
If the collision model controls the wall interaction,.
rhoEqn solve()
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
void storeState()
Store the current cloud state.
autoPtr< CollisionModel< CollidingCloud< CloudType > > > collisionModel_
Collision model.
A class for handling words, derived from string.
Definition: word.H:59
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void info()
Print cloud information.
DSMCCloud< dsmcParcel > CloudType
messageStream Info
void evolve()
Evolve the cloud.
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
void readFields(const Mesh &mesh, const IOobjectList &objects, PtrList< GeoField > &fields)
Templated collision model class.
void moveCollide(TrackData &td, const scalar deltaT)
Move-collide particles.
void cloudReset(CollidingCloud< CloudType > &c)
Reset state of cloud.
void motion(TrackData &td)
Particle motion.
virtual ~CollidingCloud()
Destructor.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
error FatalError
A class for managing sub-cycling times.
Definition: subCycleTime.H:48
void endSubCycle()
End the sub-cycling and reset the time-state.
Definition: subCycleTime.C:56