JanevReactionRateI.H
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-2021 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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
27 
29 (
30  const scalar A,
31  const scalar beta,
32  const scalar Ta,
34 )
35 :
36  A_(A),
37  beta_(beta),
38  Ta_(Ta),
39  b_(b)
40 {}
41 
42 
44 (
45  const speciesTable&,
46  const dictionary& dict
47 )
48 :
49  A_(dict.lookup<scalar>("A")),
50  beta_(dict.lookup<scalar>("beta")),
51  Ta_(dict.lookup<scalar>("Ta")),
52  b_(dict.lookup("b"))
53 {}
54 
55 
56 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
57 
59 {}
60 
61 
63 {}
64 
65 
66 inline Foam::scalar Foam::JanevReactionRate::operator()
67 (
68  const scalar p,
69  const scalar T,
70  const scalarField&,
71  const label
72 ) const
73 {
74  scalar lta = A_;
75 
76  if (mag(beta_) > vSmall)
77  {
78  lta *= pow(T, beta_);
79  }
80 
81  scalar expArg = 0;
82 
83  if (mag(Ta_) > vSmall)
84  {
85  expArg -= Ta_/T;
86  }
87 
88  scalar lnT = log(T);
89 
90  for (int n=0; n<nb_; n++)
91  {
92  expArg += b_[n]*pow(lnT, n);
93  }
94 
95  lta *= exp(expArg);
96 
97  return lta;
98 }
99 
100 
101 inline Foam::scalar Foam::JanevReactionRate::ddT
102 (
103  const scalar p,
104  const scalar T,
105  const scalarField&,
106  const label
107 ) const
108 {
109  scalar lta = A_;
110 
111  if (mag(beta_) > vSmall)
112  {
113  lta *= pow(T, beta_);
114  }
115 
116  scalar expArg = 0;
117 
118  if (mag(Ta_) > vSmall)
119  {
120  expArg -= Ta_/T;
121  }
122 
123  scalar lnT = log(T);
124 
125  for (int n=0; n<nb_; n++)
126  {
127  expArg += b_[n]*pow(lnT, n);
128  }
129 
130  scalar deriv = b_[1];
131 
132  for (int n=2; n<nb_; n++)
133  {
134  deriv += n*b_[n]*pow(lnT, n-1);
135  }
136 
137  lta *= exp(expArg);
138 
139  return lta*(beta_+Ta_/T+deriv)/T;
140 }
141 
142 
145 {
146  return NullObjectRef<List<Tuple2<label, scalar>>>();
147 }
148 
149 
151 (
152  const scalar p,
153  const scalar T,
154  const scalarField& c,
155  const label li,
156  scalarField& dcidc
157 ) const
158 {}
159 
160 
161 inline Foam::scalar Foam::JanevReactionRate::dcidT
162 (
163  const scalar p,
164  const scalar T,
165  const scalarField& c,
166  const label li
167 ) const
168 {
169  return 0;
170 }
171 
172 
174 {
175  writeKeyword(os, "A") << A_ << nl;
176  writeKeyword(os, "beta") << beta_ << nl;
177  writeKeyword(os, "Ta") << Ta_ << nl;
178  writeKeyword(os, "b") << b_ << nl;
179 }
180 
181 
182 inline Foam::Ostream& Foam::operator<<
183 (
184  Ostream& os,
185  const JanevReactionRate& jrr
186 )
187 {
188  jrr.write(os);
189  return os;
190 }
191 
192 
193 // ************************************************************************* //
virtual Ostream & write(const char)=0
Write character.
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
Ostream & writeKeyword(Foam::Ostream &os, const keyType &kw)
Write the keyword to the Ostream with the current level of indentation.
Definition: keyType.C:155
dimensionedScalar log(const dimensionedScalar &ds)
void postEvaluate() const
Post-evaluation hook.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
void preEvaluate() const
Pre-evaluation hook.
const List< Tuple2< label, scalar > > & beta() const
Third-body efficiencies (beta = 1-alpha)
Janev, Langer, Evans and Post reaction rate.
dimensionedScalar exp(const dimensionedScalar &ds)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
static const char nl
Definition: Ostream.H:260
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar dcidT(const scalar p, const scalar T, const scalarField &c, const label li) const
Temperature derivative of the pressure dependent term.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void dcidc(const scalar p, const scalar T, const scalarField &c, const label li, scalarField &dcidc) const
Species concentration derivative of the pressure dependent term.
A wordList with hashed indices for faster lookup by name.
JanevReactionRate(const scalar A, const scalar beta, const scalar Ta, const FixedList< scalar, nb_ > b)
Construct from components.
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
volScalarField & p
scalar ddT(const scalar p, const scalar T, const scalarField &c, const label li) const
void write(Ostream &os) const
Write to stream.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844