MichaelisMentenReactionRateI.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) 2018-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 speciesTable& st,
31  const dictionary& dict
32 )
33 :
34  species_(st),
35  Vmax_(dict.lookup<scalar>("Vmax")),
36  Km_(dict.lookup<scalar>("Km")),
37  s_(st[dict.lookup("S")])
38 {}
39 
40 
41 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
42 
44 {}
45 
46 
48 {}
49 
50 
51 inline Foam::scalar Foam::MichaelisMentenReactionRate::operator()
52 (
53  const scalar p,
54  const scalar T,
55  const scalarField& c,
56  const label li
57 ) const
58 {
59  return Vmax_/(Km_ + c[s_]);
60 }
61 
62 
64 (
65  const scalar p,
66  const scalar T,
67  const scalarField& c,
68  const label li
69 ) const
70 {
71  return 0;
72 }
73 
74 
76 {
77  return true;
78 }
79 
80 
82 (
83  const scalar p,
84  const scalar T,
85  const scalarField& c,
86  const label li,
87  scalarField& ddc
88 ) const
89 {
90  ddc = 0;
91  ddc[s_] = - Vmax_/sqr(Km_ + c[s_]);
92 }
93 
94 
96 {
97  writeEntry(os, "Vmax", Vmax_);
98  writeEntry(os, "Km", Km_);
99  writeEntry(os, "S", species_[s_]);
100 }
101 
102 
103 inline Foam::Ostream& Foam::operator<<
104 (
105  Ostream& os,
106  const MichaelisMentenReactionRate& mmrr
107 )
108 {
109  mmrr.write(os);
110  return os;
111 }
112 
113 
114 // ************************************************************************* //
Michaelis-Menten reaction rate for enzymatic reactions.
void preEvaluate() const
Pre-evaluation hook.
void postEvaluate() const
Post-evaluation hook.
void write(Ostream &os) const
Write to stream.
void ddc(const scalar p, const scalar T, const scalarField &c, const label li, scalarField &ddc) const
The derivative of the rate w.r.t. concentration.
MichaelisMentenReactionRate(const speciesTable &species, const dictionary &dict)
Construct from dictionary.
scalar ddT(const scalar p, const scalar T, const scalarField &c, const label li) const
The derivative of the rate w.r.t. temperature.
bool hasDdc() const
Is the rate a function of concentration?
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.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
A wordList with hashed indices for faster lookup by name.
const dimensionedScalar c
Speed of light in a vacuum.
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
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dictionary dict
volScalarField & p