From 7c11e11982a35fe2c8528b7292ef420ecf8d8fac Mon Sep 17 00:00:00 2001
From: ignis <ignis@postech.ac.kr>
Date: Wed, 17 Nov 2021 22:34:31 +0900
Subject: [PATCH] added H2O for moisture and N2 for ash. reduced magic numbers

---
 coal.py | 53 +++++++++++++++++++++++++++++++++++++++--------------
 1 file changed, 39 insertions(+), 14 deletions(-)

diff --git a/coal.py b/coal.py
index 0bde57c..84f8c81 100644
--- a/coal.py
+++ b/coal.py
@@ -4,7 +4,9 @@ import cantera as ct
 
 # Air stream Temperatures and Mass flow rates
 t = list(map(float, ''' 348.15  315.65  600.79  308.54   318.03  306.11  339.45 '''.split()))
-m = list(map(float, ''' 114.4   6.94    362.92  7.25     7.25    7.25    7.25 '''.split())) 
+m = list(map(float, ''' 114.4   6.94    362.92  7.25     7.25    7.25    7.25 '''.split()))
+
+stdT = 298.15
 
 coalT = 348.15
 coalMfr = 56.8813
@@ -45,17 +47,18 @@ print("Total Air  flow rate    = ", airmix.mass)
 
 """#############################################################################
 Dummy gaseous coal object containing elementary composition
-    mass of moisture and ash contents is added to N2
+    mass of ash contents is added to N2
 #############################################################################"""
 
 coal = ct.Solution(thermo='IdealGas', species=system_species)
-coal.TPY = 348.15, ct.one_atm, '''\
-    C:  81.41,
-    H2:  5.47,
-    O2: 10.83,
-    N2: 36.9647930755,
-    S:   0.57
-    '''
+coal.TPY = 348.15, ct.one_atm, {
+    "C":  81.41,
+        "H2":  5.47,
+        "O2": 10.83,
+        "N2": 1.72 + 100 * ash / vm_fc,
+        "S":   0.57,
+        "H2O" : 100 * ash / vm_fc,
+}
 
 
 """#############################################################################
@@ -70,10 +73,32 @@ print("HV of Coal     , kJ/kg  = ", coalHV)
 print("HV of Coal(daf), kJ/kg  = ", coalHVdaf)
 
 # sum product of 1kg coal enthalpy of formation - kJ/kg
+def hf_product_coefs (product_name, element_name):
+    """ returns kJ/kg """
+    # J/kmol / kg/kmol / 1000
+    return (nasa_species[product_name].thermo.h(stdT)
+            / ct.Element(element_name).weight
+            / nasa_species[product_name].composition[element_name]
+            / 1000. )
+
+'''
+H(stdT) values from NASA polynomial
+-  32762.281048240053
+- 119952.68929829611
+-   9258.666766366705
+
+H_f values from google search data
+-  32762.45348
+- 141887.60121
+-   8919.38250
+
+Discrepency in Enthalpy of Formation for H2O is due to phase difference
+value above is for vapor and otherwise is for liquid water
+'''
 sum_product_hf = (
-        -  32762.45348 * coal.mass_fraction_dict()['C']
-        - 141887.60121 * coal.mass_fraction_dict()['H2']
-        -   8919.38250 * coal.mass_fraction_dict()['S'])
+          hf_product_coefs("CO2", "C") * coal.elemental_mass_fraction('C')
+        + hf_product_coefs("H2O", "H") * coal.elemental_mass_fraction('H')
+        + hf_product_coefs("SO2", "S") * coal.elemental_mass_fraction('S'))
 print("Sum(Hf_product), kJ/kg  = ", sum_product_hf)
 
 sum_coal_hf = - coalHV + sum_product_hf
@@ -98,8 +123,8 @@ print("Dummy Coal H    , kJ/kg = ", coal.enthalpy_mass/1000.)
 """#############################################################################
 Can't Set coal object's H to coal_enthalpy since it is composed of gaseous
     C S H2 O2 N2 and can't have H=coal_enthalpy with T >= 0 K.
-Therefore only difference between real coal enthalpy and dummy gas coal is 
-    caculated and enthalpy difference will be added later.
+Therefore only difference between real coal enthalpy and dummy gas coal is
+    caculated and it will be added later.
 #############################################################################"""
 
 enthalpy_added_after_mixing = (coal_enthalpy*1000 - coal.enthalpy_mass) * fuelMfr # J