KochFriedlander.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) 2019-2023 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 "KochFriedlander.H"
28 #include "fvmSup.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace diameterModels
35 {
36 namespace shapeModels
37 {
38 namespace sinteringModels
39 {
42 }
43 }
44 }
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
52 (
53  const dictionary& dict,
54  const fractal& fractalShape
55 )
56 :
57  sinteringModel(dict, fractalShape),
58  dict_(dict.subDict(type() + "Coeffs")),
59  Cs_(dict_.lookup<scalar>("Cs")),
60  n_(dict_.lookup<scalar>("n")),
61  m_(dict_.lookup<scalar>("m")),
62  Ta_(dict_.lookup<scalar>("Ta")),
63  dpMin_(dict_.lookupOrDefault<scalar>("dpMin", 0))
64 {}
65 
66 
67 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
68 
71 {}
72 
73 
74 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
75 
78 {
80  (
82  (
83  "tau",
84  fractal_.group().mesh(),
86  )
87  );
88 
89  volScalarField::Internal& tau = tTau.ref();
90 
91  const sizeGroup& fi = fractal_.group();
92  const volScalarField& kappai = fractal_.fld();
93  const volScalarField& T = fractal_.group().phase().thermo().T();
94 
95  forAll(tau, celli)
96  {
97  const scalar dp = 6/max(6/fi.dSph().value(), kappai[celli]);
98 
99  tau[celli] =
100  Cs_*pow(dp, n_)*pow(T[celli], m_)*exp(Ta_/T[celli]*(1 - dpMin_/dp));
101  }
102 
103  return tTau;
104 }
105 
106 
109 {
110  const sizeGroup& fi = fractal_.group();
111  const volScalarField& kappai = fractal_.fld();
112  const volScalarField& alpha = fi.phase();
113 
115  (
116  IOobject
117  (
118  typedName("R"),
119  fi.time().name(),
120  fi.mesh()
121  ),
122  fi.mesh(),
124  );
125 
126  volScalarField::Internal tau(this->tau());
127 
128  forAll(R, celli)
129  {
130  R[celli] = fi[celli]*alpha[celli]/tau[celli];
131  }
132 
133  return fvm::Sp(R, kappai) - 6/fi.dSph()*R;
134 }
135 
136 
137 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Mesh & mesh() const
Return mesh.
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &, const Field< Type > &)
Return a temporary field constructed from name, mesh,.
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const Time & time() const
Return time.
Definition: IOobject.C:318
Class for modelling the shape of particle aggregates using the concept of fractal geometry....
Definition: fractal.H:106
Abstract base class for modelling sintering of primary particles in fractal aggregates.
virtual tmp< fvScalarMatrix > R() const
Sintering source term.
virtual tmp< volScalarField::Internal > tau() const
Characteristic time for sintering.
KochFriedlander(const dictionary &dict, const fractal &fractalShape)
Construct from a dictionary and a fractal shape model.
Single size class fraction field representing a fixed particle volume as defined by the user through ...
Definition: sizeGroup.H:102
const dimensionedScalar & dSph() const
Return representative spherical diameter of the sizeGroup.
Definition: sizeGroupI.H:51
const velocityGroup & group() const
Return const-reference to the velocityGroup.
Definition: sizeGroupI.H:44
const phaseModel & phase() const
Return const-reference to the phase.
Definition: sizeGroupI.H:37
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
const Type & value() const
Return const reference to value.
const word & name() const
Return const reference to name.
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
Calculate the matrix for implicit and explicit sources.
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
addToRunTimeSelectionTable(sinteringModel, KochFriedlander, dictionary)
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
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)
static scalar R(const scalar a, const scalar x)
Definition: invIncGamma.C:102
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dictionary dict