From 98a9566fc4303e6306a7d576a897fc9e8275a384 Mon Sep 17 00:00:00 2001
From: Ray Speth <speth@mit.edu>
Date: Thu, 10 Jan 2019 20:10:01 -0500
Subject: [PATCH] [Input] Create StoichSubstance objects from YAML definitions

The (fixed) density is now read from the species definition, rather than the
phase definition.
---
 src/thermo/StoichSubstance.cpp | 23 ++++++++++++++++++++
 src/thermo/ThermoFactory.cpp   |  1 +
 test/data/thermo-models.yaml   | 39 ++++++++++++++++++++++++++++++++++
 test/thermo/thermoFromYaml.cpp | 23 ++++++++++++++++++++
 4 files changed, 86 insertions(+)
 create mode 100644 test/data/thermo-models.yaml

diff --git a/src/thermo/StoichSubstance.cpp b/src/thermo/StoichSubstance.cpp
index 59cd0f96f..051b9f39f 100644
--- a/src/thermo/StoichSubstance.cpp
+++ b/src/thermo/StoichSubstance.cpp
@@ -122,6 +122,29 @@ void StoichSubstance::initThermo()
                            "stoichiometric substances may only contain one species.");
     }
 
+    if (species(0)->input.hasKey("equation-of-state")) {
+        auto& eos = species(0)->input["equation-of-state"].as<AnyMap>();
+        if (eos.getString("model", "") != "constant-volume") {
+            throw CanteraError("StoichSubstance::initThermo",
+                "fixed-stoichiometry model requires constant-volume species "
+                "model for species '{}'", speciesName(0));
+        }
+        if (eos.hasKey("density")) {
+            setDensity(eos.convert("density", "kg/m^3"));
+        } else if (eos.hasKey("molar-density")) {
+            setMolarDensity(eos.convert("molar-density", "kmol/m^3"));
+        } else if (eos.hasKey("molar-volume")) {
+            setMolarDensity(1.0 / eos.convert("molar-volume", "m^3/kmol"));
+        } else {
+            throw CanteraError("StoichSubstance::initThermo",
+                "equation-of-state entry for species '{}' is missing 'density',"
+                " 'molar-volume' or 'molar-density' specification",
+                speciesName(0));
+        }
+    } else if (m_input.hasKey("density")) {
+        setDensity(m_input.convert("density", "kg/m^3"));
+    }
+
     // Store the reference pressure in the variables for the class.
     m_p0 = refPressure();
 
diff --git a/src/thermo/ThermoFactory.cpp b/src/thermo/ThermoFactory.cpp
index bdc7c24a4..7dd04c754 100644
--- a/src/thermo/ThermoFactory.cpp
+++ b/src/thermo/ThermoFactory.cpp
@@ -56,6 +56,7 @@ ThermoFactory::ThermoFactory()
     reg("Edge", []() { return new EdgePhase(); });
     reg("Metal", []() { return new MetalPhase(); });
     reg("StoichSubstance", []() { return new StoichSubstance(); });
+    m_synonyms["fixed-stoichiometry"] = "StoichSubstance";
     reg("PureFluid", []() { return new PureFluidPhase(); });
     reg("LatticeSolid", []() { return new LatticeSolidPhase(); });
     reg("Lattice", []() { return new LatticePhase(); });
diff --git a/test/data/thermo-models.yaml b/test/data/thermo-models.yaml
new file mode 100644
index 000000000..f21a9af83
--- /dev/null
+++ b/test/data/thermo-models.yaml
@@ -0,0 +1,39 @@
+phases:
+- name: NaCl(s)
+  thermo: fixed-stoichiometry
+  species: [NaCl(s)]
+  # specification of density was moved to the species object
+  state: {T: 300.0, P: 1 atm}
+
+- name: KCl(s)
+  thermo: fixed-stoichiometry
+  species: [KCl(s)]
+  state: {T: 300.0, P: 1 atm}
+
+
+species:
+- name: NaCl(s)
+  composition: {Na: 1, Cl: 1}
+  thermo:
+    model: Shomate
+    temperature-ranges: [250, 1075]
+    data:
+    - [50.72389, 6.672267, -2.517167, 10.15934, -0.200675, -427.2115, 130.3973]
+  equation-of-state:
+    model: constant-volume
+    # For constant density, specify either (mass) "density", "molar-density", or
+    # "molar-volume"
+    density: 2.165 g/cm^3
+
+- name: KCl(s)
+  composition: {K: 1, Cl: 1}
+  thermo:
+    model: NASA7
+    temperature-ranges: [300, 1000]
+    data:
+    - [5.3934311E+00, 2.6535242E-03, 9.6075655E-07, -5.0251843E-09,
+       4.0721228E-12, -5.4248389E+04, -2.1596814E+01]
+  equation-of-state:
+    model: constant-volume
+    molar-volume: 0.0376521717 L/mol
+
diff --git a/test/thermo/thermoFromYaml.cpp b/test/thermo/thermoFromYaml.cpp
index 42ce58cc5..93566584a 100644
--- a/test/thermo/thermoFromYaml.cpp
+++ b/test/thermo/thermoFromYaml.cpp
@@ -67,3 +67,26 @@ TEST(ThermoFromYaml, speciesAll)
     EXPECT_EQ(thermo->species(1)->name, "NO");
     EXPECT_EQ(thermo->species(2)->name, "N2");
 }
+
+TEST(ThermoFromYaml, StoichSubstance1)
+{
+    AnyMap infile = AnyMap::fromYamlFile("thermo-models.yaml");
+    auto phaseNodes = infile["phases"].asMap("name");
+    auto thermo = newPhase(*phaseNodes.at("NaCl(s)"), infile);
+    EXPECT_EQ(thermo->type(), "StoichSubstance");
+    EXPECT_EQ(thermo->nSpecies(), (size_t) 1);
+    EXPECT_EQ(thermo->nElements(), (size_t) 2);
+    EXPECT_DOUBLE_EQ(thermo->density(), 2165.0);
+    EXPECT_DOUBLE_EQ(thermo->cp_mass(), 864.8437519457644); // Regression test based on XML
+}
+
+TEST(ThermoFromYaml, StoichSubstance2)
+{
+    AnyMap infile = AnyMap::fromYamlFile("thermo-models.yaml");
+    auto phaseNodes = infile["phases"].asMap("name");
+    auto thermo = newPhase(*phaseNodes.at("KCl(s)"), infile);
+    EXPECT_EQ(thermo->type(), "StoichSubstance");
+    EXPECT_EQ(thermo->nSpecies(), (size_t) 1);
+    EXPECT_EQ(thermo->nElements(), (size_t) 2);
+    EXPECT_NEAR(thermo->density(), 1980, 0.1);
+}