[Species] Add molecular weight attribute

The molecular weight is computed for each species when the molecular
weight is accessed on the species or added to a phase. The Phase
private attribute m_mwt is not changed.

include/cantera/thermo/Species.h

@@ -51,21 +51,42 @@ class Species
     //! species, where it represents the number of sites occupied.
     double size;
 
+    //! The molecular weight [amu] of the species.
+    const double molecularWeight();
+
+    //! Set the molecular weight of the species.
+    /*!
+     * Since phases can have custom element weights, the phase will always call this
+     * method when a species is added to that phase. The species may also call this
+     * method the first time the molecularWeight method is called if the species has
+     * not been added to a phase.
+     *
+     * @param weight: The weight of this species to assign, unless `compute` is `true`
+     * @param compute: If `true`, the weight will be computed from standard element
+     * weights from Elements.cpp, and the value of `weight` is ignored.
+     */
+    void setMolecularWeight(double weight, bool compute=false);
+
     shared_ptr<TransportData> transport;
 
     //! Thermodynamic data for the species
     shared_ptr<SpeciesThermoInterpType> thermo;
 
     //! Input parameters used to define a species, for example from a YAML input file.
     AnyMap input;
+
+protected:
+
+    //! The molecular weight of the species, in atomic mass units. Includes
+    //! electron mass for charged species.
+    double m_molecularWeight = Undef;
 };
 
 //! Create a new Species object from an AnyMap specification
 unique_ptr<Species> newSpecies(const AnyMap& node);
 
 //! Generate Species objects for each item (an AnyMap) in `items`.
 std::vector<shared_ptr<Species>> getSpecies(const AnyValue& items);
-
 }
 
 #endif

interfaces/cython/cantera/_cantera.pxd

@@ -233,6 +233,7 @@ cdef extern from "cantera/thermo/Species.h" namespace "Cantera":
         Composition composition
         double charge
         double size
+        double molecularWeight() except +translate_exception
         CxxAnyMap parameters(CxxThermoPhase*) except +translate_exception
         CxxAnyMap input
 

interfaces/cython/cantera/test/test_thermo.py

@@ -1232,6 +1232,19 @@ def test_standalone(self):
         self.assertEqual(len(c), 2)
         self.assertEqual(c['C'], 1)
         self.assertEqual(c['H'], 4)
+        self.assertNear(s.molecular_weight, 16.043)
+
+    def test_molecular_weight_with_electron(self):
+        yaml = """
+            name: Li+[elyt]
+            composition: {Li: 1, E: -1}
+            thermo:
+              model: constant-cp
+              h0: -278.49 kJ/mol
+              s0: 13.4 J/mol/K
+        """
+        s = ct.Species.from_yaml(yaml)
+        assert np.isclose(s.molecular_weight, 6.939451420091127)
 
     def test_defaults(self):
         s = ct.Species('H2')

interfaces/cython/cantera/thermo.pyx

@@ -176,6 +176,11 @@ cdef class Species:
         def __get__(self):
             return self.species.size
 
+    property molecular_weight:
+        """ The molecular weight [amu] of the species. """
+        def __get__(self):
+            return self.species.molecularWeight()
+
     property thermo:
         """
         Get/Set the species reference-state thermodynamic data, as an instance
@@ -2055,6 +2060,7 @@ def _element_names():
     names = elementNames()
     return tuple(pystr(n) for n in names)
 
+
 class Element:
     """
     An element or a named isotope defined in Cantera.

src/thermo/Phase.cpp

@@ -846,6 +846,7 @@ bool Phase::addSpecies(shared_ptr<Species> spec) {
     wt = std::max(wt, Tiny);
     m_molwts.push_back(wt);
     m_rmolwts.push_back(1.0/wt);
+    spec->setMolecularWeight(wt);
     m_kk++;
 
     // Ensure that the Phase has a valid mass fraction vector that sums to

src/thermo/Species.cpp

@@ -1,6 +1,7 @@
 // This file is part of Cantera. See License.txt in the top-level directory or
 // at https://cantera.org/license.txt for license and copyright information.
 
+#include "cantera/thermo/Elements.h"
 #include "cantera/thermo/Species.h"
 #include "cantera/thermo/SpeciesThermoInterpType.h"
 #include "cantera/thermo/SpeciesThermoFactory.h"
@@ -13,6 +14,8 @@
 #include <limits>
 #include <set>
 
+using namespace std;
+
 namespace Cantera {
 
 Species::Species()
@@ -34,6 +37,47 @@ Species::~Species()
 {
 }
 
+const double Species::molecularWeight() {
+    if (m_molecularWeight == Undef) {
+        setMolecularWeight(-1.0, true);
+    }
+    return m_molecularWeight;
+}
+
+void Species::setMolecularWeight(double weight, bool compute) {
+    if (compute) {
+        weight = 0.0;
+        const auto& elements = elementSymbolToWeight();
+        const auto& isotopes = isotopeSymbolToWeight();
+        for (const auto& comp : composition) {
+            auto search = elements.find(comp.first);
+            if (search != elements.end()) {
+                if (search->second < 0) {
+                    throw CanteraError("setMolecularWeight",
+                        "element '{}' has no stable isotopes", comp.first);
+                }
+                weight += search->second * comp.second;
+            } else {
+                search = isotopes.find(comp.first);
+                if (search != isotopes.end() && search->second > 0) {
+                    weight += search->second * comp.second;
+                }
+            }
+        }
+    }
+    if (m_molecularWeight != Undef) {
+        double weight_cmp = fabs(weight - m_molecularWeight) / max(weight, m_molecularWeight);
+        if (weight_cmp > 1.0e-9) {
+            throw CanteraError(
+                "Species::setMolecularWeight",
+                "Molecular weight of a species cannot be changed."
+            );
+        }
+    }
+
+    m_molecularWeight = weight;
+}
+
 AnyMap Species::parameters(const ThermoPhase* phase, bool withInput) const
 {
     AnyMap speciesNode;