SRIFallOffFunctionI.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 b,
32  const scalar c,
33  const scalar d,
34  const scalar e
35 )
36 :
37  a_(a),
38  b_(b),
39  c_(c),
40  d_(d),
41  e_(e)
42 {}
43 
44 
46 :
47  a_(dict.lookup<scalar>("a")),
48  b_(dict.lookup<scalar>("b")),
49  c_(dict.lookup<scalar>("c")),
50  d_(dict.lookup<scalar>("d")),
51  e_(dict.lookup<scalar>("e"))
52 {}
53 
54 
55 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
56 
57 inline Foam::scalar Foam::SRIFallOffFunction::operator()
58 (
59  const scalar T,
60  const scalar Pr
61 ) const
62 {
63  const scalar logPr = log10(max(Pr, small));
64 
65  const scalar X = 1/(1 + sqr(logPr));
66 
67  const scalar psi = a_*exp(-b_/T) + exp(-T/c_);
68 
69  return d_*pow(psi, X)*pow(T, e_);
70 }
71 
72 
73 inline Foam::scalar Foam::SRIFallOffFunction::ddT
74 (
75  const scalar T,
76  const scalar Pr,
77  const scalar F
78 ) const
79 {
80  const scalar logPr = log10(max(Pr, small));
81 
82  const scalar X = 1/(1 + sqr(logPr));
83 
84  const scalar psi = a_*exp(-b_/T) + exp(-T/c_);
85  const scalar dpsidT = a_*b_/sqr(T)*exp(-b_/T) - 1/c_*exp(-T/c_);
86 
87  return F*(X/psi*dpsidT + e_/T);
88 }
89 
90 
91 inline Foam::scalar Foam::SRIFallOffFunction::ddPr
92 (
93  const scalar T,
94  const scalar Pr,
95  const scalar F
96 ) const
97 {
98  static const scalar logTen = log(scalar(10));
99 
100  const scalar logPr = log10(max(Pr, small));
101  const scalar dlogPrdPr = Pr >= small ? 1/(logTen*Pr) : 0;
102 
103  const scalar X = 1/(1 + sqr(logPr));
104  const scalar dXdPr = -sqr(X)*2*logPr*dlogPrdPr;
105 
106  const scalar psi = a_*exp(-b_/T) + exp(-T/c_);
107 
108  return F*log(psi)*dXdPr;
109 }
110 
111 
113 {
114  writeEntry(os, "a", a_);
115  writeEntry(os, "b", b_);
116  writeEntry(os, "c", c_);
117  writeEntry(os, "d", d_);
118  writeEntry(os, "e", e_);
119 }
120 
121 
122 // * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
123 
124 inline Foam::Ostream& Foam::operator<<
125 (
126  Foam::Ostream& os,
127  const Foam::SRIFallOffFunction& srifof
128 )
129 {
130  srifof.write(os);
131  return os;
132 }
133 
134 
135 // ************************************************************************* //
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
virtual Ostream & write(const char)=0
Write character.
The SRI fall-off function.
scalar ddT(const scalar T, const scalar Pr, const scalar F) const
void write(Ostream &os) const
Write to stream.
SRIFallOffFunction(const scalar a, const scalar b, const scalar c, const scalar d, const scalar e)
Construct from components.
scalar ddPr(const scalar T, const scalar Pr, const scalar F) const
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
const volScalarField & psi
volScalarField & b
Definition: createFields.H:27
const dimensionedScalar e
Elementary charge.
const dimensionedScalar F
Faraday constant: default SI units: [C/mol].
const dimensionedScalar c
Speed of light in a vacuum.
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar log10(const dimensionedScalar &ds)
dimensionedScalar log(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dictionary dict