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-2019 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  scalar X = 1.0/(1 + sqr(log10(max(Pr, small))));
64  return d_*pow(a_*exp(-b_/T) + exp(-T/c_), X)*pow(T, e_);
65 }
66 
67 
68 inline Foam::scalar Foam::SRIFallOffFunction::ddT
69 (
70  const scalar Pr,
71  const scalar F,
72  const scalar dPrdT,
73  const scalar T
74 ) const
75 {
76  scalar X = 1.0/(1 + sqr(log10(max(Pr, small))));
77  scalar dXdPr = -X*X*2*log10(Pr)/Pr/log(10.0);
78  return
79  F
80  *(
81  e_/T
82  + X
83  *(a_*b_*exp(-b_/T)/sqr(T) - exp(-T/c_)/c_)
84  /(a_*exp(-b_/T) + exp(-T/c_))
85  + dXdPr*dPrdT*log(a_*exp(-b_/T) + exp(-T/c_))
86  );
87 }
88 
89 
90 inline Foam::scalar Foam::SRIFallOffFunction::ddc
91 (
92  const scalar Pr,
93  const scalar F,
94  const scalar dPrdc,
95  const scalar T
96 ) const
97 {
98  scalar X = 1.0/(1 + sqr(log10(max(Pr, small))));
99  scalar dXdPr = -X*X*2*log10(Pr)/Pr/log(10.0);
100  scalar dXdc = dXdPr*dPrdc;
101  return F*dXdc*log(a_*exp(-b_/T) + exp(-T/c_));
102 }
103 
104 
106 {
107  writeEntry(os, "a", a_);
108  writeEntry(os, "b", b_);
109  writeEntry(os, "c", c_);
110  writeEntry(os, "d", d_);
111  writeEntry(os, "e", e_);
112 }
113 
114 
115 // * * * * * * * * * * * * * * * Ostream Operator * * * * * * * * * * * * * //
116 
117 inline Foam::Ostream& Foam::operator<<
118 (
119  Foam::Ostream& os,
120  const Foam::SRIFallOffFunction& srifof
121 )
122 {
123  srifof.write(os);
124  return os;
125 }
126 
127 
128 // ************************************************************************* //
dictionary dict
virtual Ostream & write(const char)=0
Write character.
The SRI fall-off function.
dimensionedScalar log(const dimensionedScalar &ds)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
scalar ddc(const scalar Pr, const scalar F, const scalar dPrdc, const scalar T) const
void write(Ostream &os) const
Write to stream.
const dimensionedScalar c
Speed of light in a vacuum.
stressControl lookup("compactNormalStress") >> compactNormalStress
SRIFallOffFunction(const scalar a, const scalar b, const scalar c, const scalar d, const scalar e)
Construct from components.
dimensionedScalar exp(const dimensionedScalar &ds)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
scalar ddT(const scalar Pr, const scalar F, const scalar dPrdT, const scalar T) const
dimensionedScalar log10(const dimensionedScalar &ds)
const dimensionedScalar e
Elementary charge.
Definition: doubleScalar.H:105