SuppressionCollision.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) 2013-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 "SuppressionCollision.H"
27 #include "parcelCloud.H"
28 
29 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
30 
31 template<class CloudType>
33 (
34  typename CloudType::parcelType::trackingData& td
35 )
36 {
37  const parcelCloud& sc =
38  this->owner().mesh().template
39  lookupObject<parcelCloud>(suppressionCloud_);
40 
41  volScalarField vDotSweep(sc.vDotSweep());
42 
43  dimensionedScalar Dt("dt", dimTime, td.trackTime());
44  volScalarField P(typedName("p"), 1.0 - exp(-vDotSweep*Dt));
45 
46  forAllIter(typename CloudType, this->owner(), iter)
47  {
48  typename CloudType::parcelType& p = iter();
49  label celli = p.cell();
50 
51  scalar xx = this->owner().rndGen().template sample01<scalar>();
52 
53  if (xx < P[celli])
54  {
55  p.canCombust() = -1;
56  p.typeId() = max(p.typeId(), suppressedParcelType_);
57  }
58  }
59 }
60 
61 
62 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
63 
64 template<class CloudType>
66 (
67  const dictionary& dict,
68  CloudType& owner
69 )
70 :
71  StochasticCollisionModel<CloudType>(dict, owner, typeName),
72  suppressionCloud_(this->coeffDict().lookup("suppressionCloud")),
73  suppressedParcelType_
74  (
75  this->coeffDict().lookupOrDefault("suppressedParcelType", -1)
76  )
77 {}
78 
79 
80 template<class CloudType>
82 (
84 )
85 :
87  suppressionCloud_(cm.suppressionCloud_),
88  suppressedParcelType_(cm.suppressedParcelType_)
89 {}
90 
91 
92 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
93 
94 template<class CloudType>
96 {}
97 
98 
99 // ************************************************************************* //
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:459
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:80
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:225
Generic GeometricField class.
Templated stochastic collision model class.
Inter-cloud collision model, whereby the canReact flag can be used to inhibit devolatilisation and su...
SuppressionCollision(const dictionary &dict, CloudType &owner)
Construct from dictionary.
virtual void collide(typename CloudType::parcelType::trackingData &td)
Update the model.
virtual ~SuppressionCollision()
Destructor.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
virtual const tmp< volScalarField > vDotSweep() const =0
Volume swept rate of parcels per cell.
Virtual abstract base class for parcel clouds. As parcelCloudBase but with additional virtualisation ...
Definition: parcelCloud.H:57
dimensionedScalar exp(const dimensionedScalar &ds)
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
const dimensionSet dimTime
word typedName(Name name)
Return the name of the object within the given type.
Definition: typeInfo.H:176
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dictionary dict
volScalarField & p