From 7211f088912ec34678963a2ddbf5071983b17ed9 Mon Sep 17 00:00:00 2001
From: Ray Speth <speth@mit.edu>
Date: Fri, 14 Dec 2018 18:31:46 -0500
Subject: [PATCH] [Input] Parse YAML entries for NASA 9-coefficient polynomials

---
 src/thermo/SpeciesThermoFactory.cpp     | 25 +++++++++++++++++++++++++
 test/thermo/thermoParameterizations.cpp | 25 +++++++++++++++++++++++++
 2 files changed, 50 insertions(+)

diff --git a/src/thermo/SpeciesThermoFactory.cpp b/src/thermo/SpeciesThermoFactory.cpp
index 229f23ce0..1ca8003f8 100644
--- a/src/thermo/SpeciesThermoFactory.cpp
+++ b/src/thermo/SpeciesThermoFactory.cpp
@@ -343,6 +343,27 @@ static SpeciesThermoInterpType* newNasa9ThermoFromXML(
 }
 
 
+void setupNasa9Poly(Nasa9PolyMultiTempRegion& thermo, const AnyMap& node,
+                    const UnitSystem& units)
+{
+    setupSpeciesThermo(thermo, node, units);
+    vector_fp Tranges = units.convert(
+        node.at("temperature-ranges").asVector<AnyValue>(2, 999), "K");
+    const auto& data = node.at("data").asVector<vector_fp>(Tranges.size()-1);
+    map<double, vector_fp> regions;
+    for (size_t i = 0; i < data.size(); i++) {
+        if (data[i].size() != 9) {
+            throw CanteraError("setupNasa9Poly", "Wrong number of coefficients "
+                "for NASA9 polynomial. Expected 9, but got {}", data[i].size());
+        }
+        regions[Tranges[i]] = data[i];
+    }
+    thermo.setMinTemp(Tranges.front());
+    thermo.setMaxTemp(Tranges.back());
+    thermo.setParameters(regions);
+}
+
+
 SpeciesThermoInterpType* newSpeciesThermoInterpType(const XML_Node& thermo)
 {
     std::string model = toLowerCopy(thermo["model"]);
@@ -410,6 +431,10 @@ unique_ptr<SpeciesThermoInterpType> newSpeciesThermo(
         unique_ptr<ShomatePoly2> thermo(new ShomatePoly2());
         setupShomatePoly(*thermo, node, units);
         return unique_ptr<SpeciesThermoInterpType>(move(thermo));
+    } else if (model == "NASA9") {
+        unique_ptr<Nasa9PolyMultiTempRegion> thermo(new Nasa9PolyMultiTempRegion());
+        setupNasa9Poly(*thermo, node, units);
+        return unique_ptr<SpeciesThermoInterpType>(move(thermo));
     } else {
         throw CanteraError("newSpeciesThermo",
             "Uknown thermo model '{}'", model);
diff --git a/test/thermo/thermoParameterizations.cpp b/test/thermo/thermoParameterizations.cpp
index 082c68abe..4740375d8 100644
--- a/test/thermo/thermoParameterizations.cpp
+++ b/test/thermo/thermoParameterizations.cpp
@@ -181,3 +181,28 @@ TEST(SpeciesThermo, Shomate2FromYaml1) {
     EXPECT_DOUBLE_EQ(h_RT, -5.747001162488512);
     EXPECT_DOUBLE_EQ(s_R, 29.878976075583324);
 }
+
+TEST(SpeciesThermo, Nasa9PolyFromYaml) {
+    AnyMap data = AnyMap::fromYamlString(
+        "model: NASA9\n"
+        "temperature-ranges: [200.00, 1000.00, 6000.0, 20000]\n"
+        "reference-pressure: 1 bar\n"
+        "data:\n"
+        "- [2.210371497E+04, -3.818461820E+02, 6.082738360E+00, -8.530914410E-03,\n"
+        "   1.384646189E-05, -9.625793620E-09, 2.519705809E-12, 7.108460860E+02,\n"
+        "   -1.076003744E+01]\n"
+        "- [5.877124060E+05, -2.239249073E+03, 6.066949220E+00, -6.139685500E-04,\n"
+        "   1.491806679E-07,  -1.923105485E-11, 1.061954386E-15, 1.283210415E+04,\n"
+        "   -1.586640027E+01]\n"
+        "- [8.310139160E+08, -6.420733540E+05, 2.020264635E+02, -3.065092046E-02,\n"
+        "   2.486903333E-06, -9.705954110E-11, 1.437538881E-15, 4.938707040E+06,\n"
+        "   -1.672099740E+03]");
+    UnitSystem U;
+    double cp_R, h_RT, s_R;
+    auto st = newSpeciesThermo(data, U);
+    EXPECT_DOUBLE_EQ(st->refPressure(), 1e5);
+    st->updatePropertiesTemp(2000, &cp_R, &h_RT, &s_R);
+    EXPECT_DOUBLE_EQ(cp_R, 4.326181187976);
+    EXPECT_DOUBLE_EQ(h_RT, 3.3757914517886856);
+    EXPECT_DOUBLE_EQ(s_R, 30.31743870437559);
+}