VoidFraction.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-2020 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 "VoidFraction.H"
27 
28 // * * * * * * * * * * * * * Protectd Member Functions * * * * * * * * * * * //
29 
30 template<class CloudType>
32 {
33  if (thetaPtr_.valid())
34  {
35  thetaPtr_->write();
36  }
37  else
38  {
40  << "thetaPtr not valid" << abort(FatalError);
41  }
42 }
43 
44 
45 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46 
47 template<class CloudType>
49 (
50  const dictionary& dict,
51  CloudType& owner,
52  const word& modelName
53 )
54 :
55  CloudFunctionObject<CloudType>(dict, owner, modelName, typeName),
56  thetaPtr_(nullptr)
57 {}
58 
59 
60 template<class CloudType>
62 (
63  const VoidFraction<CloudType>& vf
64 )
65 :
67  thetaPtr_(nullptr)
68 {}
69 
70 
71 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
72 
73 template<class CloudType>
75 {}
76 
77 
78 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
79 
80 template<class CloudType>
82 {
83  if (thetaPtr_.valid())
84  {
85  thetaPtr_->primitiveFieldRef() = 0.0;
86  }
87  else
88  {
89  const fvMesh& mesh = this->owner().mesh();
90 
91  thetaPtr_.reset
92  (
93  new volScalarField
94  (
95  IOobject
96  (
97  this->owner().name() + "Theta",
98  mesh.time().timeName(),
99  mesh,
100  IOobject::NO_READ,
101  IOobject::NO_WRITE
102  ),
103  mesh,
105  )
106  );
107  }
108 }
109 
110 
111 template<class CloudType>
113 {
114  volScalarField& theta = thetaPtr_();
115 
116  const fvMesh& mesh = this->owner().mesh();
117 
118  theta.primitiveFieldRef() /= mesh.time().deltaTValue()*mesh.V();
119 
121 }
122 
123 
124 template<class CloudType>
126 (
127  parcelType& p,
128  const scalar dt,
129  const point&,
130  bool&
131 )
132 {
133  volScalarField& theta = thetaPtr_();
134 
135  theta[p.cell()] += dt*p.nParticle()*p.volume();
136 }
137 
138 
139 // ************************************************************************* //
virtual void postMove(parcelType &p, const scalar dt, const point &position0, bool &keepParticle)
Post-move hook.
Definition: VoidFraction.C:126
dictionary dict
Creates particle void fraction field on carrier phase.
Definition: VoidFraction.H:50
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
const dimensionSet dimless
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:372
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
fvMesh & mesh
static word timeName(const scalar, const int precision=curPrecision_)
Return time name of given scalar time.
Definition: Time.C:666
A class for handling words, derived from string.
Definition: word.H:59
virtual ~VoidFraction()
Destructor.
Definition: VoidFraction.C:74
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:34
errorManip< error > abort(error &err)
Definition: errorManip.H:131
virtual void preEvolve()
Pre-evolve hook.
Definition: VoidFraction.C:81
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
virtual void postEvolve()
Post-evolve hook.
Definition: VoidFraction.C:112
virtual void write()
Write post-processing info.
Definition: VoidFraction.C:31
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
VoidFraction(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
Definition: VoidFraction.C:49
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:75
Templated cloud function object base class.