Go to the documentation of this file.
37 template<
class CloudType>
44 this->subModelProperties(),
52 this->subModelProperties(),
60 this->subModelProperties(),
69 template<
class CloudType>
81 packingModel_(
nullptr),
82 dampingModel_(
nullptr),
83 isotropyModel_(
nullptr)
90 <<
"MPPIC modelling not available for steady state calculations"
99 this->deleteLostParticles();
105 template<
class CloudType>
113 packingModel_(
c.packingModel_->clone()),
114 dampingModel_(
c.dampingModel_->clone()),
115 isotropyModel_(
c.isotropyModel_->clone())
119 template<
class CloudType>
128 packingModel_(
nullptr),
129 dampingModel_(
nullptr),
130 isotropyModel_(
nullptr)
136 template<
class CloudType>
143 template<
class CloudType>
150 clone(this->
name() +
"Copy").ptr()
156 template<
class CloudType>
159 this->cloudReset(cloudCopyPtr_());
160 cloudCopyPtr_.clear();
164 template<
class CloudType>
169 typename parcelType::trackingData td(*
this);
171 this->
solve(*
this, td);
176 template<
class CloudType>
177 template<
class TrackCloudType>
180 TrackCloudType&
cloud,
181 typename parcelType::trackingData& td
188 td.part() = parcelType::trackingData::tpLinearTrack;
189 CloudType::move(
cloud, td, this->db().time().deltaTValue());
196 this->forces().setCalcNonCoupled(
false);
197 this->forces().setCalcCoupled(
false);
203 if (dampingModel_->active())
205 if (this->
mesh().moving())
208 <<
"MPPIC damping modelling does not support moving meshes."
213 td.updateAverages(
cloud);
216 dampingModel_->cacheFields(
true);
219 td.part() = parcelType::trackingData::tpDampingNoTrack;
220 CloudType::move(
cloud, td, this->db().time().deltaTValue());
223 td.part() = parcelType::trackingData::tpCorrectTrack;
224 CloudType::move(
cloud, td, this->db().time().deltaTValue());
227 dampingModel_->cacheFields(
false);
234 if (packingModel_->active())
236 if (this->
mesh().moving())
239 <<
"MPPIC packing modelling does not support moving meshes."
244 td.updateAverages(
cloud);
245 packingModel_->cacheFields(
true);
246 td.part() = parcelType::trackingData::tpPackingNoTrack;
247 CloudType::move(
cloud, td, this->db().time().deltaTValue());
248 td.part() = parcelType::trackingData::tpCorrectTrack;
249 CloudType::move(
cloud, td, this->db().time().deltaTValue());
250 packingModel_->cacheFields(
false);
257 if (isotropyModel_->active())
260 td.updateAverages(
cloud);
263 isotropyModel_->calculate();
271 this->updateCellOccupancy();
274 this->forces().setCalcNonCoupled(
true);
279 template<
class CloudType>
286 const scalar alphaMin =
gMin(
alpha().primitiveField());
289 Info<<
" Min cell volume fraction = " << alphaMin <<
endl;
316 Info<<
" Min dense number of parcels = " << nMin <<
endl;
Selector class for relaxation factors, solver type and solution.
void setModels()
Set cloud sub-models.
const word cloudName(propsDict.get< word >("cloud"))
A class for handling words, derived from Foam::string.
const dimensionedScalar mu
Atomic mass unit.
A class for managing temporary objects.
void storeState()
Store the current cloud state.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
void motion(TrackCloudType &cloud, typename parcelType::trackingData &td)
Particle motion.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Base class for collisional damping models.
void evolve()
Evolve the cloud.
#define forAll(list, i)
Loop across all elements in list.
const List< DynamicList< molecule * > > & cellOccupancy
messageStream Info
Information stream (stdout output on master, null elsewhere)
DSMCCloud< dsmcParcel > CloudType
dimensionedScalar alphaMax("alphaMax", dimless/dimTime, laminarTransport)
Base class for packing models.
reduce(hasMovingMesh, orOp< bool >())
Mesh data needed to do the Finite Volume discretisation.
const uniformDimensionedVectorField & g
errorManipArg< error, int > exit(error &err, const int errNo=1)
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const wordHashSet &selectedFields, LIFOStack< regIOobject * > &storedObjects)
Read the selected GeometricFields of the templated type.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
A cloud is a registry collection of lagrangian particles.
void restoreState()
Reset the current cloud to the previously stored state.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const dimensionedScalar c
Speed of light in a vacuum.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Base class for collisional return-to-isotropy models.
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
Type gMin(const FieldField< Field, Type > &f)
Adds MPPIC modelling to kinematic clouds.
virtual ~MPPICCloud()
Destructor.
Type gMax(const FieldField< Field, Type > &f)