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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "CollidingCloud.H"
29#include "CollisionModel.H"
30#include "NoCollision.H"
31
32// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
33
34template<class CloudType>
36{
37 collisionModel_.reset
38 (
40 (
41 this->subModelProperties(),
42 *this
43 ).ptr()
44 );
45}
46
47
48template<class CloudType>
49template<class TrackCloudType>
51(
52 TrackCloudType& cloud,
53 typename parcelType::trackingData& td,
54 const scalar deltaT
55)
56{
57 td.part() = parcelType::trackingData::tpVelocityHalfStep;
58 CloudType::move(cloud, td, deltaT);
59
60 td.part() = parcelType::trackingData::tpLinearTrack;
61 CloudType::move(cloud, td, deltaT);
62
63 // td.part() = parcelType::trackingData::tpRotationalTrack;
64 // CloudType::move(cloud, td, deltaT);
65
66 this->updateCellOccupancy();
67
68 this->collision().collide();
69
70 td.part() = parcelType::trackingData::tpVelocityHalfStep;
71 CloudType::move(cloud, td, deltaT);
72}
73
74
75
76template<class CloudType>
78{
79 CloudType::cloudReset(c);
80
81 collisionModel_.reset(c.collisionModel_.ptr());
82}
83
84
85// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
86
87template<class CloudType>
89(
90 const word& cloudName,
91 const volScalarField& rho,
92 const volVectorField& U,
93 const volScalarField& mu,
94 const dimensionedVector& g,
95 bool readFields
96)
97:
98 CloudType(cloudName, rho, U, mu, g, false),
99 constProps_(this->particleProperties()),
100 collisionModel_(nullptr)
101{
102 if (this->solution().active())
103 {
104 setModels();
105
106 if (readFields)
107 {
109 this->deleteLostParticles();
110 }
111
112 if
113 (
114 this->solution().steadyState()
116 )
117 {
119 << "Collision modelling not currently available "
120 << "for steady state calculations" << exit(FatalError);
121 }
122 }
123}
124
126template<class CloudType>
128(
130 const word& name
131)
132:
134 collisionModel_(c.collisionModel_->clone())
135{}
136
137
138template<class CloudType>
140(
141 const fvMesh& mesh,
142 const word& name,
144)
145:
146 CloudType(mesh, name, c),
147 collisionModel_(nullptr)
148{}
149
150
151// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
153template<class CloudType>
155{}
156
157
158// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
160template<class CloudType>
162{
163 cloudCopyPtr_.reset
164 (
165 static_cast<CollidingCloud<CloudType>*>
166 (
167 clone(this->name() + "Copy").ptr()
168 )
169 );
170}
171
172
173template<class CloudType>
175{
176 cloudReset(cloudCopyPtr_());
177 cloudCopyPtr_.clear();
178}
179
180
181template<class CloudType>
183{
184 if (this->solution().canEvolve())
185 {
186 typename parcelType::trackingData td(*this);
187
188 this->solve(*this, td);
189 }
190}
191
192
193template<class CloudType>
194template<class TrackCloudType>
196(
197 TrackCloudType& cloud,
198 typename parcelType::trackingData& td
199)
200{
201 // Sympletic leapfrog integration of particle forces:
202 // + apply half deltaV with stored force
203 // + move positions with new velocity
204 // + calculate forces in new position
205 // + apply half deltaV with new force
206
207 label nSubCycles = collision().nSubCycles();
208
209 if (nSubCycles > 1)
210 {
211 Info<< " " << nSubCycles << " move-collide subCycles" << endl;
212
213 subCycleTime moveCollideSubCycle
214 (
215 const_cast<Time&>(this->db().time()),
216 nSubCycles
217 );
218
219 while(!(++moveCollideSubCycle).end())
221 moveCollide(cloud, td, this->db().time().deltaTValue());
222 }
224 moveCollideSubCycle.endSubCycle();
225 }
226 else
227 {
228 moveCollide(cloud, td, this->db().time().deltaTValue());
229 }
231
232
233template<class CloudType>
235{
236 CloudType::info();
237
238 scalar rotationalKineticEnergy = rotationalKineticEnergyOfSystem();
239 reduce(rotationalKineticEnergy, sumOp<scalar>());
241 Info<< " Rotational kinetic energy = "
242 << rotationalKineticEnergy << nl;
243}
244
245
246// ************************************************************************* //
reduce(hasMovingMesh, orOp< bool >())
const uniformDimensionedVectorField & g
void deleteLostParticles()
Remove lost particles from cloud and delete.
Definition: Cloud.C:118
Adds coolisions to kinematic clouds.
void setModels()
Set cloud sub-models.
void storeState()
Store the current cloud state.
autoPtr< CollisionModel< CollidingCloud< CloudType > > > collisionModel_
Collision model.
void cloudReset(CollidingCloud< CloudType > &c)
Reset state of cloud.
void moveCollide(TrackCloudType &cloud, typename parcelType::trackingData &td, const scalar deltaT)
Move-collide particles.
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:75
Place holder for 'none' option.
Definition: NoCollision.H:57
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
A cloud is a registry collection of lagrangian particles.
Definition: cloud.H:60
const motionSolver & motion() const
Return the motionSolver.
Class used to pass data into container.
Reads fields from the time directories and adds them to the mesh database for further post-processing...
Definition: readFields.H:158
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:66
A class for managing sub-cycling times.
Definition: subCycleTime.H:51
void endSubCycle()
End the sub-cycling and reset the time-state.
Definition: subCycleTime.C:68
A class for handling words, derived from Foam::string.
Definition: word.H:68
U
Definition: pEqn.H:72
const volScalarField & mu
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
messageStream Info
Information stream (stdout output on master, null elsewhere)
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
bool isType(const Type &t)
Check is typeid is identical to the TargetType.
Definition: typeInfo.H:218
error FatalError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
CEqn solve()
const word cloudName(propsDict.get< word >("cloud"))