From 242684dc3d396df725b2480193f71e04ae500129 Mon Sep 17 00:00:00 2001
From: Henry <Henry>
Date: Thu, 27 Oct 2011 17:19:26 +0100
Subject: [PATCH] Thermo: limit the temperature during energy->temperature
 conversion

---
 .../specie/thermo/eConst/eConstThermo.H       |  3 ++
 .../specie/thermo/eConst/eConstThermoI.H      | 10 +++++
 .../specie/thermo/hConst/hConstThermo.H       |  3 ++
 .../specie/thermo/hConst/hConstThermoI.H      | 10 +++++
 .../thermo/hPolynomial/hPolynomialThermo.H    |  3 ++
 .../thermo/hPolynomial/hPolynomialThermoI.H   | 10 +++++
 .../specie/thermo/janaf/janafThermo.H         |  6 +--
 .../specie/thermo/janaf/janafThermoI.H        | 40 ++++++++++---------
 .../specie/thermo/specieThermo/specieThermo.H |  5 ++-
 .../thermo/specieThermo/specieThermoI.H       | 33 ++++++++++++---
 10 files changed, 95 insertions(+), 28 deletions(-)

diff --git a/src/thermophysicalModels/specie/thermo/eConst/eConstThermo.H b/src/thermophysicalModels/specie/thermo/eConst/eConstThermo.H
index 9154126ada7..bd8bdb45738 100644
--- a/src/thermophysicalModels/specie/thermo/eConst/eConstThermo.H
+++ b/src/thermophysicalModels/specie/thermo/eConst/eConstThermo.H
@@ -135,6 +135,9 @@ public:
 
     // Member Functions
 
+        //- Limit the temperature to be in the range Tlow_ to Thigh_
+        inline scalar limit(const scalar T) const;
+
         // Fundamental properties
 
             //- Heat capacity at constant pressure [J/(kmol K)]
diff --git a/src/thermophysicalModels/specie/thermo/eConst/eConstThermoI.H b/src/thermophysicalModels/specie/thermo/eConst/eConstThermoI.H
index f031aa63f6a..31bd392b5dd 100644
--- a/src/thermophysicalModels/specie/thermo/eConst/eConstThermoI.H
+++ b/src/thermophysicalModels/specie/thermo/eConst/eConstThermoI.H
@@ -89,6 +89,16 @@ Foam::eConstThermo<EquationOfState>::New(const dictionary& dict)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class EquationOfState>
+inline Foam::scalar Foam::eConstThermo<EquationOfState>::limit
+(
+    const scalar T
+) const
+{
+    return T;
+}
+
+
 template<class EquationOfState>
 inline Foam::scalar Foam::eConstThermo<EquationOfState>::cp
 (
diff --git a/src/thermophysicalModels/specie/thermo/hConst/hConstThermo.H b/src/thermophysicalModels/specie/thermo/hConst/hConstThermo.H
index 014de5eef25..9a8352a064d 100644
--- a/src/thermophysicalModels/specie/thermo/hConst/hConstThermo.H
+++ b/src/thermophysicalModels/specie/thermo/hConst/hConstThermo.H
@@ -133,6 +133,9 @@ public:
 
     // Member Functions
 
+        //- Limit the temperature to be in the range Tlow_ to Thigh_
+        inline scalar limit(const scalar T) const;
+
         // Fundamental properties
 
             //- Heat capacity at constant pressure [J/(kmol K)]
diff --git a/src/thermophysicalModels/specie/thermo/hConst/hConstThermoI.H b/src/thermophysicalModels/specie/thermo/hConst/hConstThermoI.H
index 76ec41b620c..0a4a35358cd 100644
--- a/src/thermophysicalModels/specie/thermo/hConst/hConstThermoI.H
+++ b/src/thermophysicalModels/specie/thermo/hConst/hConstThermoI.H
@@ -89,6 +89,16 @@ Foam::hConstThermo<equationOfState>::New(const dictionary& dict)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class EquationOfState>
+inline Foam::scalar Foam::hConstThermo<EquationOfState>::limit
+(
+    const scalar T
+) const
+{
+    return T;
+}
+
+
 template<class equationOfState>
 inline Foam::scalar Foam::hConstThermo<equationOfState>::cp
 (
diff --git a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.H b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.H
index 7a4dcb3867a..31954a4a567 100644
--- a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.H
+++ b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermo.H
@@ -151,6 +151,9 @@ public:
 
     // Member Functions
 
+        //- Limit the temperature to be in the range Tlow_ to Thigh_
+        inline scalar limit(const scalar T) const;
+
         // Fundamental properties
 
             //- Heat capacity at constant pressure [J/(kmol K)]
diff --git a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermoI.H b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermoI.H
index 620019e6d13..a164392c617 100644
--- a/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermoI.H
+++ b/src/thermophysicalModels/specie/thermo/hPolynomial/hPolynomialThermoI.H
@@ -82,6 +82,16 @@ inline Foam::hPolynomialThermo<EquationOfState, PolySize>::hPolynomialThermo
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class EquationOfState, int PolySize>
+inline Foam::scalar Foam::hPolynomialThermo<EquationOfState, PolySize>::limit
+(
+    const scalar T
+) const
+{
+    return T;
+}
+
+
 template<class EquationOfState, int PolySize>
 inline Foam::scalar Foam::hPolynomialThermo<EquationOfState, PolySize>::cp
 (
diff --git a/src/thermophysicalModels/specie/thermo/janaf/janafThermo.H b/src/thermophysicalModels/specie/thermo/janaf/janafThermo.H
index 9c414afaac3..b552fad0945 100644
--- a/src/thermophysicalModels/specie/thermo/janaf/janafThermo.H
+++ b/src/thermophysicalModels/specie/thermo/janaf/janafThermo.H
@@ -119,9 +119,6 @@ private:
         //- Check that input data is valid
         void checkInputData() const;
 
-        //- Check given temperature is within the range of the fitted coeffs
-        inline void checkT(const scalar T) const;
-
         //- Return the coefficients corresponding to the given temperature
         inline const coeffArray& coeffs(const scalar T) const;
 
@@ -153,6 +150,9 @@ public:
 
     // Member Functions
 
+        //- Limit the temperature to be in the range Tlow_ to Thigh_
+        inline scalar limit(const scalar T) const;
+
         // Fundamental properties
 
             //- Heat capacity at constant pressure [J/(kmol K)]
diff --git a/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H b/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H
index 05f0c0c0ccf..5740816c239 100644
--- a/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H
+++ b/src/thermophysicalModels/specie/thermo/janaf/janafThermoI.H
@@ -52,22 +52,6 @@ inline Foam::janafThermo<EquationOfState>::janafThermo
 }
 
 
-template<class EquationOfState>
-inline void Foam::janafThermo<EquationOfState>::checkT(const scalar T) const
-{
-    if (T < Tlow_ || T > Thigh_)
-    {
-        FatalErrorIn
-        (
-            "janafThermo<EquationOfState>::checkT(const scalar T) const"
-        )   << "attempt to use janafThermo<EquationOfState>"
-               " out of temperature range "
-            << Tlow_ << " -> " << Thigh_ << ";  T = " << T
-            << abort(FatalError);
-    }
-}
-
-
 template<class EquationOfState>
 inline const typename Foam::janafThermo<EquationOfState>::coeffArray&
 Foam::janafThermo<EquationOfState>::coeffs
@@ -75,8 +59,6 @@ Foam::janafThermo<EquationOfState>::coeffs
     const scalar T
 ) const
 {
-    checkT(T);
-
     if (T < Tcommon_)
     {
         return lowCpCoeffs_;
@@ -112,6 +94,28 @@ inline Foam::janafThermo<EquationOfState>::janafThermo
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class EquationOfState>
+inline Foam::scalar Foam::janafThermo<EquationOfState>::limit(const scalar T) const
+{
+    if (T < Tlow_ || T > Thigh_)
+    {
+        WarningIn
+        (
+            "janafThermo<EquationOfState>::limit(const scalar T) const"
+        )   << "attempt to use janafThermo<EquationOfState>"
+               " out of temperature range "
+            << Tlow_ << " -> " << Thigh_ << ";  T = " << T
+            << endl;
+
+        return min(max(T, Tlow_), Thigh_);
+    }
+    else
+    {
+        return T;
+    }
+}
+
+
 template<class EquationOfState>
 inline Foam::scalar Foam::janafThermo<EquationOfState>::cp
 (
diff --git a/src/thermophysicalModels/specie/thermo/specieThermo/specieThermo.H b/src/thermophysicalModels/specie/thermo/specieThermo/specieThermo.H
index c6f0f7a3168..93268b9e6b4 100644
--- a/src/thermophysicalModels/specie/thermo/specieThermo/specieThermo.H
+++ b/src/thermophysicalModels/specie/thermo/specieThermo/specieThermo.H
@@ -103,14 +103,15 @@ class specieThermo
 
     // Private Member Functions
 
-        //- return the temperature corresponding to the value of the
+        //- Return the temperature corresponding to the value of the
         //  thermodynamic property f, given the function f = F(T) and dF(T)/dT
         inline scalar T
         (
             scalar f,
             scalar T0,
             scalar (specieThermo::*F)(const scalar) const,
-            scalar (specieThermo::*dFdT)(const scalar) const
+            scalar (specieThermo::*dFdT)(const scalar) const,
+            scalar (specieThermo::*limit)(const scalar) const
         ) const;
 
 
diff --git a/src/thermophysicalModels/specie/thermo/specieThermo/specieThermoI.H b/src/thermophysicalModels/specie/thermo/specieThermo/specieThermoI.H
index f8632376c20..b2d2edd5c14 100644
--- a/src/thermophysicalModels/specie/thermo/specieThermo/specieThermoI.H
+++ b/src/thermophysicalModels/specie/thermo/specieThermo/specieThermoI.H
@@ -43,7 +43,8 @@ inline Foam::scalar Foam::specieThermo<Thermo>::T
     scalar f,
     scalar T0,
     scalar (specieThermo<Thermo>::*F)(const scalar) const,
-    scalar (specieThermo<Thermo>::*dFdT)(const scalar) const
+    scalar (specieThermo<Thermo>::*dFdT)(const scalar) const,
+    scalar (specieThermo<Thermo>::*limit)(const scalar) const
 ) const
 {
     scalar Test = T0;
@@ -54,7 +55,8 @@ inline Foam::scalar Foam::specieThermo<Thermo>::T
     do
     {
         Test = Tnew;
-        Tnew = Test - ((this->*F)(Test) - f)/(this->*dFdT)(Test);
+        Tnew =
+            (this->*limit)(Test - ((this->*F)(Test) - f)/(this->*dFdT)(Test));
 
         if (iter++ > maxIter_)
         {
@@ -276,7 +278,14 @@ inline Foam::scalar Foam::specieThermo<Thermo>::TH
     const scalar T0
 ) const
 {
-    return T(h, T0, &specieThermo<Thermo>::H, &specieThermo<Thermo>::Cp);
+    return T
+    (
+        h,
+        T0,
+        &specieThermo<Thermo>::H,
+        &specieThermo<Thermo>::Cp,
+        &specieThermo<Thermo>::limit
+    );
 }
 
 
@@ -287,7 +296,14 @@ inline Foam::scalar Foam::specieThermo<Thermo>::THs
     const scalar T0
 ) const
 {
-    return T(hs, T0, &specieThermo<Thermo>::Hs, &specieThermo<Thermo>::Cp);
+    return T
+    (
+        hs,
+        T0,
+        &specieThermo<Thermo>::Hs,
+        &specieThermo<Thermo>::Cp,
+        &specieThermo<Thermo>::limit
+    );
 }
 
 
@@ -298,7 +314,14 @@ inline Foam::scalar Foam::specieThermo<Thermo>::TE
     const scalar T0
 ) const
 {
-    return T(e, T0, &specieThermo<Thermo>::E, &specieThermo<Thermo>::Cv);
+    return T
+    (
+        e,
+        T0,
+        &specieThermo<Thermo>::E,
+        &specieThermo<Thermo>::Cv,
+        &specieThermo<Thermo>::limit
+    );
 }
 
 
-- 
GitLab