From 4bc4fda90c01927b3cf7148b11ab62604abde783 Mon Sep 17 00:00:00 2001 From: Thomas Fiala Date: Tue, 17 Feb 2015 23:23:51 +0000 Subject: [PATCH] Add C++ and Cython support for computing elemental mass/mole fractions --- include/cantera/thermo/Phase.h | 35 +++++++++++++++++++ interfaces/cython/cantera/_cantera.pxd | 2 ++ interfaces/cython/cantera/onedim.py | 47 ++++++++++++++++++++++++++ interfaces/cython/cantera/thermo.pyx | 39 +++++++++++++++++++++ src/thermo/Phase.cpp | 25 ++++++++++++++ 5 files changed, 148 insertions(+) diff --git a/include/cantera/thermo/Phase.h b/include/cantera/thermo/Phase.h index 5a9fe873e..4fe685f4b 100644 --- a/include/cantera/thermo/Phase.h +++ b/include/cantera/thermo/Phase.h @@ -526,6 +526,41 @@ public: //! is the number of species in the phase. virtual void setConcentrations(const doublereal* const conc); + //! Elemental mass fraction of element m + /*! + * The elemental mass fraction \f$Z_{\mathrm{mass},m}\f$ of element \f$m\f$ + * is defined as + * \f[ + * Z_{\mathrm{mass},m} = \sum_k \frac{a_{m,k} M_m}{M_k} Y_k + * \f] + * with \f$a_{m,k}\f$ being the number of atoms of element \f$m\f$ in + * species \f$k\f$, \f$M_m\f$ the atomic weight of element \f$m\f$, + * \f$M_k\f$ the molecular weight of species \f$k\f$, and \f$Y_k\f$ + * the mass fraction of species \f$k\f$. + * + * @param[in] m Index of the element within the phase. If m is outside + * the valid range, an exception will be thrown. + * + * @return the elemental mass fraction of element m. + */ + doublereal elementalMassFraction(const size_t m) const; + + //! Elemental mole fraction of element m + /*! + * The elemental mole fraction \f$Z_{\mathrm{mole},m}\f$ of element \f$m\f$ + * is defined as + * \f[ + * Z_{\mathrm{mole},m} = \sum_k \frac{a_{m,k}}{\sum_j a_{j,k}} X_k + * \f] + * with \f$a_{m,k}\f$ being the number of atoms of element \f$m\f$ in + * species \f$k\f$and \f$X_k\f$ the mole fraction of species \f$k\f$. + * + * @param[in] m Index of the element within the phase. If m is outside the + * valid range, an exception will be thrown. + * @return the elemental mole fraction of element m. + */ + doublereal elementalMoleFraction(const size_t m) const; + //! Returns a const pointer to the start of the moleFraction/MW array. //! This array is the array of mole fractions, each divided by the mean //! molecular weight. diff --git a/interfaces/cython/cantera/_cantera.pxd b/interfaces/cython/cantera/_cantera.pxd index 29da5fd17..c318ec906 100644 --- a/interfaces/cython/cantera/_cantera.pxd +++ b/interfaces/cython/cantera/_cantera.pxd @@ -100,6 +100,8 @@ cdef extern from "cantera/thermo/ThermoPhase.h" namespace "Cantera": double moleFraction(string) except + double concentration(size_t) except + + double elementalMassFraction(size_t) except + + double elementalMoleFraction(size_t) except + # state setters void setState_TR(double, double) except + diff --git a/interfaces/cython/cantera/onedim.py b/interfaces/cython/cantera/onedim.py index b6c19fd31..543693153 100644 --- a/interfaces/cython/cantera/onedim.py +++ b/interfaces/cython/cantera/onedim.py @@ -158,6 +158,53 @@ class FlameBase(Sim1D): """ return self.profile(self.flame, 'lambda') + def elemental_mass_fraction(self, m): + r""" + Get the elemental mass fraction :math:`Z_{\mathrm{mass},m}` of element + :math:`m` at each grid point, which is defined as: + + .. math:: Z_{\mathrm{mass},m} = \sum_k \frac{a_{m,k} M_m}{M_k} Y_k + + with :math:`a_{m,k}` being the number of atoms of element :math:`m` in + species :math:`k`, :math:`M_m` the atomic weight of element :math:`m`, + :math:`M_k` the molecular weight of species :math:`k`, and :math:`Y_k` + the mass fraction of species :math:`k`. + + :param m: + Base element, may be specified by name or by index. + + >>> phase.elemental_mass_fraction('H') + [1.0, ..., 0.0] + """ + vals = np.empty(self.flame.n_points) + for i in range(self.flame.n_points): + self.set_gas_state(i) + vals[i] = self.gas.elemental_mass_fraction(m) + return vals + + def elemental_mole_fraction(self, m): + r""" + Get the elemental mole fraction :math:`Z_{\mathrm{mole},m}` of element + :math:`m` at each grid point, which is defined as: + + .. math:: Z_{\mathrm{mole},m} = \sum_k \frac{a_{m,k}}{\sum_j a_{j,k}} X_k + + with :math:`a_{m,k}` being the number of atoms of element :math:`m` in + species :math:`k` and :math:`X_k` the mole fraction of species + :math:`k`. + + :param m: + Base element, may be specified by name or by index. + + >>> phase.elemental_mole_fraction('H') + [1.0, ..., 0.0] + """ + vals = np.empty(self.flame.n_points) + for i in range(self.flame.n_points): + self.set_gas_state(i) + vals[i] = self.gas.elemental_mole_fraction(m) + return vals + def solution(self, component, point=None): """ Get the solution at one point or for the full flame domain (if diff --git a/interfaces/cython/cantera/thermo.pyx b/interfaces/cython/cantera/thermo.pyx index a4f4c3d87..c88be3fd8 100644 --- a/interfaces/cython/cantera/thermo.pyx +++ b/interfaces/cython/cantera/thermo.pyx @@ -328,6 +328,45 @@ cdef class ThermoPhase(_SolutionBase): def __set__(self, C): self._setArray1(thermo_setConcentrations, C) + def elemental_mass_fraction(self, m): + r""" + Get the elemental mass fraction :math:`Z_{\mathrm{mass},m}` of element + :math:`m` as defined by: + + .. math:: Z_{\mathrm{mass},m} = \sum_k \frac{a_{m,k} M_m}{M_k} Y_k + + with :math:`a_{m,k}` being the number of atoms of element :math:`m` in + species :math:`k`, :math:`M_m` the atomic weight of element :math:`m`, + :math:`M_k` the molecular weight of species :math:`k`, and :math:`Y_k` + the mass fraction of species :math:`k`. + + :param m: + Base element, may be specified by name or by index. + + >>> phase.elemental_mass_fraction('H') + 1.0 + """ + return self.thermo.elementalMassFraction(self.element_index(m)) + + def elemental_mole_fraction(self, m): + r""" + Get the elemental mole fraction :math:`Z_{\mathrm{mole},m}` of element + :math:`m` as defined by: + + .. math:: Z_{\mathrm{mole},m} = \sum_k \frac{a_{m,k}}{\sum_j a_{j,k}} X_k + + with :math:`a_{m,k}` being the number of atoms of element :math:`m` in + species :math:`k` and :math:`X_k` the mole fraction of species + :math:`k`. + + :param m: + Base element, may be specified by name or by index. + + >>> phase.elemental_mole_fraction('H') + 1.0 + """ + return self.thermo.elementalMoleFraction(self.element_index(m)) + def set_unnormalized_mass_fractions(self, Y): """ Set the mass fractions without normalizing to force sum(Y) == 1.0. diff --git a/src/thermo/Phase.cpp b/src/thermo/Phase.cpp index 2d6bc4583..3233dc30b 100644 --- a/src/thermo/Phase.cpp +++ b/src/thermo/Phase.cpp @@ -633,6 +633,31 @@ void Phase::setConcentrations(const doublereal* const conc) m_stateNum++; } +doublereal Phase::elementalMassFraction(const size_t m) const +{ + checkElementIndex(m); + doublereal Z_m = 0.0; + for (size_t k = 0; k != m_kk; ++k) { + Z_m += nAtoms(k, m) * atomicWeight(m) / molecularWeight(k) + * massFraction(k); + } + return Z_m; +} + +doublereal Phase::elementalMoleFraction(const size_t m) const +{ + checkElementIndex(m); + doublereal Z_n = 0.0; + for (size_t k = 0; k != m_kk; ++k) { + int nTotalAtoms = 0; + for (size_t l = 0; l != m_mm; ++l) { + nTotalAtoms += nAtoms(k, l); + } + Z_n += nAtoms(k, m) / nTotalAtoms * moleFraction(k); + } + return Z_n; +} + doublereal Phase::molarDensity() const { return density()/meanMolecularWeight();