CollidingCloud.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-2022 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 "subCycleTime.H"
28 
29 #include "CollisionModel.H"
30 #include "NoCollision.H"
31 
32 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
33 
34 template<class CloudType>
36 {
37  collisionModel_.reset
38  (
40  (
41  this->subModelProperties(),
42  *this
43  ).ptr()
44  );
45 }
46 
47 
48 template<class CloudType>
50 {
51  CloudType::cloudReset(c);
52 
53  collisionModel_.reset(c.collisionModel_.ptr());
54 }
55 
56 
57 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
58 
59 template<class CloudType>
61 (
62  const word& cloudName,
63  const volScalarField& rho,
64  const volVectorField& U,
65  const volScalarField& mu,
66  const dimensionedVector& g,
67  const bool readFields
68 )
69 :
70  CloudType(cloudName, rho, U, mu, g, false),
71  constProps_(this->particleProperties()),
72  collisionModel_(nullptr)
73 {
74  setModels();
75 
76  if (readFields)
77  {
79  this->deleteLostParticles();
80  }
81 
82  if
83  (
84  this->solution().steadyState()
86  )
87  {
89  << "Collision modelling not currently available "
90  << "for steady state calculations" << exit(FatalError);
91  }
92 }
93 
94 
95 template<class CloudType>
97 (
98  const word& cloudName,
99  const volScalarField& rho,
100  const volVectorField& U,
101  const dimensionedVector& g,
102  const fluidThermo& carrierThermo,
103  const bool readFields
104 )
105 :
106  CollidingCloud(cloudName, rho, U, carrierThermo.mu(), g, readFields)
107 {}
108 
109 
110 template<class CloudType>
112 (
114  const word& name
115 )
116 :
117  CloudType(c, name),
118  collisionModel_(c.collisionModel_->clone())
119 {}
120 
121 
122 template<class CloudType>
124 (
125  const fvMesh& mesh,
126  const word& name,
128 )
129 :
130  CloudType(mesh, name, c),
131  collisionModel_(nullptr)
132 {}
133 
134 
135 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
136 
137 template<class CloudType>
139 {}
140 
141 
142 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
143 
144 template<class CloudType>
146 {
147  cloudCopyPtr_.reset
148  (
149  static_cast<CollidingCloud<CloudType>*>
150  (
151  clone(this->name() + "Copy").ptr()
152  )
153  );
154 }
155 
156 
157 template<class CloudType>
159 {
160  cloudReset(cloudCopyPtr_());
161  cloudCopyPtr_.clear();
162 }
163 
164 
165 template<class CloudType>
167 {
168  if (this->solution().canEvolve())
169  {
170  typename parcelType::trackingData td(*this);
171 
172  this->solve(*this, td);
173  }
174 }
175 
176 
177 template<class CloudType>
178 template<class TrackCloudType>
180 (
181  TrackCloudType& cloud,
182  typename parcelType::trackingData& td
183 )
184 {
185  // Sympletic leapfrog integration of particle forces:
186  // + apply half deltaV with stored force
187  // + move positions with new velocity
188  // + calculate forces in new position
189  // + apply half deltaV with new force
190 
191  const label nSubCycles = collision().nSubCycles();
192 
193  if (nSubCycles > 1)
194  {
195  Info<< " " << nSubCycles << " move-collide subCycles" << endl;
196  }
197 
198  for (label subCyclei = 0; subCyclei < nSubCycles; ++ subCyclei)
199  {
200  td.stepFractionRange() =
202  (
203  scalar(subCyclei)/nSubCycles,
204  scalar(subCyclei + 1)/nSubCycles
205  );
206 
207  td.part() = parcelType::trackingData::tpVelocityHalfStep;
208  CloudType::move(cloud, td);
209 
210  td.part() = parcelType::trackingData::tpLinearTrack;
211  CloudType::move(cloud, td);
212 
213  // td.part() = parcelType::trackingData::tpRotationalTrack;
214  // CloudType::move(cloud, td);
215 
216  this->updateCellOccupancy();
217 
218  this->collision().collide();
219 
220  td.part() = parcelType::trackingData::tpVelocityHalfStep;
221  CloudType::move(cloud, td);
222  }
223 
224  td.stepFractionRange() = Pair<scalar>(0, 1);
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 deleteLostParticles()
Remove lost particles from cloud and delete.
Definition: Cloud.C:232
void move(TrackCloudType &cloud, typename ParticleType::trackingData &td)
Move the particles.
Definition: Cloud.C:275
Adds collisions to clouds.
void setModels()
Set cloud sub-models.
void motion(TrackCloudType &cloud, typename parcelType::trackingData &td)
Particle motion.
void storeState()
Store the current cloud state.
autoPtr< CollisionModel< CollidingCloud< CloudType > > > collisionModel_
Collision model.
void cloudReset(CollidingCloud< CloudType > &c)
Reset state of cloud.
CollidingCloud(const word &cloudName, const volScalarField &rho, const volVectorField &U, const volScalarField &mu, const dimensionedVector &g, const bool readFields=true)
Construct given carrier fields.
void evolve()
Evolve the cloud.
void info()
Print cloud information.
void restoreState()
Reset the current cloud to the previously stored state.
virtual ~CollidingCloud()
Destructor.
Templated collision model class.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:79
void info() const
Print cloud information.
Definition: DSMCCloud.C:970
Generic GeometricField class.
Place holder for 'none' option.
Definition: NoCollision.H:52
A cloud is a collection of lagrangian particles.
Definition: cloud.H:55
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:57
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:101
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:51
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
U
Definition: pEqn.H:72
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
const dimensionedScalar mu
Atomic mass unit.
const dimensionedScalar c
Speed of light in a vacuum.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
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
bool isType(const Type &t)
Check the typeid.
Definition: typeInfo.H:131
messageStream Info
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
T clone(const T &t)
Definition: List.H:55
error FatalError
static const char nl
Definition: Ostream.H:260
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
const word cloudName(propsDict.lookup("cloudName"))