SuppressionCollision.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) 2013-2016 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 "kinematicCloud.H"
28 
29 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
30 
31 template<class CloudType>
33 {
34  const kinematicCloud& sc =
35  this->owner().mesh().template
36  lookupObject<kinematicCloud>(suppressionCloud_);
37 
38  volScalarField vDotSweep(sc.vDotSweep());
39 
40  dimensionedScalar Dt("dt", dimTime, dt);
41  volScalarField P(type() + ":p", 1.0 - exp(-vDotSweep*Dt));
42 
43  forAllIter(typename CloudType, this->owner(), iter)
44  {
45  typename CloudType::parcelType& p = iter();
46  label celli = p.cell();
47 
48  scalar xx = this->owner().rndGen().template sample01<scalar>();
49 
50  if (xx < P[celli])
51  {
52  p.canCombust() = -1;
53  p.typeId() = max(p.typeId(), suppressedParcelType_);
54  }
55  }
56 }
57 
58 
59 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
60 
61 template<class CloudType>
63 (
64  const dictionary& dict,
65  CloudType& owner
66 )
67 :
68  StochasticCollisionModel<CloudType>(dict, owner, typeName),
69  suppressionCloud_(this->coeffDict().lookup("suppressionCloud")),
70  suppressedParcelType_
71  (
72  this->coeffDict().lookupOrDefault("suppressedParcelType", -1)
73  )
74 {}
75 
76 
77 template<class CloudType>
79 (
81 )
82 :
84  suppressionCloud_(cm.suppressionCloud_),
85  suppressedParcelType_(cm.suppressedParcelType_)
86 {}
87 
88 
89 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
90 
91 template<class CloudType>
93 {}
94 
95 
96 // ************************************************************************* //
dictionary dict
type
Types of root.
Definition: Roots.H:52
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
virtual void collide(const scalar dt)
Update the model.
const word suppressionCloud_
Name of cloud used for suppression.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Inter-cloud collision model, whereby the canReact flag can be used to inhibit devolatilisation and su...
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:453
const label suppressedParcelType_
Suppressed parcel type - optional.
dimensionedScalar exp(const dimensionedScalar &ds)
virtual const tmp< volScalarField > vDotSweep() const =0
Volume swept rate of parcels per cell.
SuppressionCollision(const dictionary &dict, CloudType &owner)
Construct from dictionary.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:218
virtual ~SuppressionCollision()
Destructor.
Virtual abstract base class for templated KinematicCloud.
Templated stochastic collision model class.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
volScalarField & p
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:69