All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ode.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) 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 #include "ode.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class ChemistryModel>
32 (
33  const typename ChemistryModel::reactionThermo& thermo
34 )
35 :
37  coeffsDict_(this->subDict("odeCoeffs")),
38  odeSolver_(ODESolver::New(*this, coeffsDict_)),
39  cTp_(this->nEqns())
40 {}
41 
42 
43 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
44 
45 template<class ChemistryModel>
47 {}
48 
49 
50 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
51 
52 template<class ChemistryModel>
54 (
55  scalar& p,
56  scalar& T,
57  scalarField& c,
58  const label li,
59  scalar& deltaT,
60  scalar& subDeltaT
61 ) const
62 {
63  // Reset the size of the ODE system to the simplified size when mechanism
64  // reduction is active
65  if (odeSolver_->resize())
66  {
67  odeSolver_->resizeField(cTp_);
68  }
69 
70  const label nSpecie = this->nSpecie();
71 
72  // Copy the concentration, T and P to the total solve-vector
73  for (int i=0; i<nSpecie; i++)
74  {
75  cTp_[i] = c[i];
76  }
77  cTp_[nSpecie] = T;
78  cTp_[nSpecie+1] = p;
79 
80  odeSolver_->solve(0, deltaT, cTp_, li, subDeltaT);
81 
82  for (int i=0; i<nSpecie; i++)
83  {
84  c[i] = max(0.0, cTp_[i]);
85  }
86  T = cTp_[nSpecie];
87  p = cTp_[nSpecie+1];
88 }
89 
90 
91 // ************************************************************************* //
An abstract base class for solving chemistry.
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
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual void solve(scalar &p, scalar &T, scalarField &c, const label li, scalar &deltaT, scalar &subDeltaT) const
Update the concentrations and return the chemical time.
Definition: ode.C:54
rhoReactionThermo & thermo
Definition: createFields.H:28
const label nSpecie
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
ode(const typename ChemistryModel::reactionThermo &thermo)
Construct from thermo.
Definition: ode.C:32
const volScalarField & T
volScalarField & p
virtual ~ode()
Destructor.
Definition: ode.C:46