Examples/pythermocalcdb nasa

examples/pythermocalcdb-nasa/README.md

@@ -0,0 +1,18 @@
+# pythermocalcdb-nasa Examples
+
+Each sub-directory contains a self-contained example. The order in
+which the examples are to appear is specified in `order.json` (an
+array of directory names in the expected order).
+
+In each example directory you'll find:
+
+* `config.toml` - must conform to the specification outlined here:
+  https://docs.pyscript.net/latest/user-guide/configuration/ This is
+  parsed and ultimately turned into a JSON representation as part of
+  the package's API object.
+* `setup.py` - Python code for contextual and environmental setup,
+  NOT SEEN BY THE END USER, but is run before the `code.py` code is
+  evaluated. Allows us to create useful (IPython) shims, avoid
+  repeating boilerplate and whatnot.
+* `code.py` - the actual code added to the editor which forms the
+  practical example of using the package.

examples/pythermocalcdb-nasa/order.json

@@ -0,0 +1,4 @@
+[
+    "species_properties",
+    "reaction_equilibrium"
+]

examples/pythermocalcdb-nasa/reaction_equilibrium/code.py

@@ -0,0 +1,92 @@
+# ---------------------------------------------------------------------
+# The water-gas shift reaction:  CO(g) + H2O(g) -> CO2(g) + H2(g)
+# ---------------------------------------------------------------------
+
+import numpy as np
+import pandas as pd
+import matplotlib.pyplot as plt
+
+from pythermodb_settings.models import Component
+from pyreactlab_core.models.reaction import Reaction
+
+heading("Water-gas shift: building the reaction")
+note(
+    "The water-gas shift is a workhorse of industrial hydrogen "
+    "production. We declare each species as a Component, then bundle "
+    "them into a Reaction object that pythermocalcdb_nasa can read."
+)
+
+CO  = Component(name="carbon monoxide", formula="CO",  state="g")
+H2O = Component(name="dihydrogen monoxide", formula="H2O", state="g")
+CO2 = Component(name="carbon dioxide", formula="CO2", state="g")
+H2  = Component(name="dihydrogen", formula="H2",  state="g")
+
+wgs = Reaction(
+    name="Water-Gas Shift",
+    reaction="CO(g) + H2O(g) => CO2(g) + H2(g)",
+    components=[CO, H2O, CO2, H2],
+)
+
+note(f"Reaction: <code>{wgs.reaction}</code>")
+
+heading("The equilibrium-constant API", level=3)
+note(
+    "Once you have a <code>model_source</code> built from NASA pickles "
+    "via <code>pyThermoLinkDB</code>, the calls below give you reaction "
+    "thermochemistry and K(T). The full version of <code>Keq</code> "
+    "uses Delta G^0(T) directly; <code>Keq_vh_shortcut</code> applies "
+    "the Van't Hoff approximation anchored at Delta H^0(298 K)."
+)
+
+api_summary = pd.DataFrame({
+    "call": [
+        "dH_rxn_STD(reaction, temperature, model_source)",
+        "dS_rxn_STD(reaction, temperature, model_source)",
+        "dG_rxn_STD(reaction, temperature, model_source)",
+        "Keq(reaction, temperature, model_source)",
+        "Keq_vh_shortcut(reaction, temperature, model_source)",
+    ],
+    "yields": [
+        "Delta H^0(T)  [J/mol]",
+        "Delta S^0(T)  [J/(mol K)]",
+        "Delta G^0(T)  [J/mol]",
+        "K(T) from Delta G^0(T)",
+        "K(T) via Van't Hoff from Delta H^0(298 K)",
+    ],
+})
+display(api_summary, append=True)
+
+heading("What K(T) looks like for an exothermic reaction", level=3)
+note(
+    "WGS is mildly exothermic (Delta H^0_298 ~ -41 kJ/mol), so K "
+    "decreases with temperature. Using textbook values for Delta H^0 "
+    "and Delta S^0 at 298 K, we sketch the Van't Hoff curve that "
+    "<code>Keq_vh_shortcut</code> would return:"
+    " <code>ln K(T) = -Delta H^0/(R T) + Delta S^0/R</code>."
+)
+
+R = 8.314462618  # J/(mol K)
+dH_298 = -41.16e3  # J/mol
+dS_298 = -42.08    # J/(mol K)
+
+T_grid = np.linspace(400.0, 1200.0, 81)
+ln_K = -dH_298 / (R * T_grid) + dS_298 / R
+K = np.exp(ln_K)
+
+fig, ax = plt.subplots(figsize=(8, 4))
+ax.semilogy(T_grid, K, color="darkgreen", linewidth=2)
+ax.axhline(1.0, color="gray", linestyle="--", linewidth=1)
+ax.set_xlabel("Temperature (K)")
+ax.set_ylabel("K(T)  (log scale)")
+ax.set_title("Van't Hoff sketch: water-gas shift K vs T")
+ax.grid(True, which="both", alpha=0.3)
+fig.tight_layout()
+display(fig, append=True)
+
+note(
+    "Swap in real NASA polynomials via <code>model_source</code> and "
+    "call <code>Keq(reaction=wgs, temperature=Temperature(value=T, "
+    "unit='K'), model_source=model_source)</code> in a loop to get "
+    "the exact curve, including the temperature dependence of "
+    "Delta H^0 and Delta S^0 that the Van't Hoff shortcut ignores."
+)

examples/pythermocalcdb-nasa/reaction_equilibrium/config.toml

@@ -0,0 +1 @@
+packages = ["pythermocalcdb-nasa", "pythermodb-settings", "pyreactlab-core", "numpy", "matplotlib"]

examples/pythermocalcdb-nasa/reaction_equilibrium/setup.py

@@ -0,0 +1,20 @@
+"""Lightweight setup for the reaction example (no IPython shim)."""
+import js
+from pyscript import window, HTML, display as _display
+
+js.alert = window.alert
+
+
+def display(*args, **kwargs):
+    return _display(
+        *args, **kwargs, target=__pyscript_display_target__,
+    )
+
+
+def heading(text, level=2):
+    display(HTML(f"<h{level}>{text}</h{level}>"), append=True)
+
+
+def note(text):
+    display(HTML(f"<p>{text}</p>"), append=True)
+

examples/pythermocalcdb-nasa/species_properties/code.py

@@ -0,0 +1,81 @@
+"""
+A first look at PyThermoCalcDB-NASA: evaluate ideal-gas species
+thermochemistry (Cp, H, S, G) from NASA polynomial coefficients.
+
+The library separates the calculation engine from data sources. In a
+typical workflow you'd point it at packaged NASA pickles via
+`load_and_build_model_source`. Here we focus on the calculation API
+itself by constructing a tiny in-memory `ModelSource`-shaped object
+holding one species: methane, with NASA-7 coefficients from the
+Burcat/NASA database (T in K, valid 200-1000 K).
+
+Reference: https://github.com/sinagilassi/PyThermoCalcDB-NASA
+"""
+from IPython.core.display import display, HTML
+import pandas as pd
+
+from pythermodb_settings.models import Component, Temperature
+from pythermocalcdb_nasa import Cp_T, H_T, S_T, G_T
+
+# NASA-7 polynomial form:
+#   Cp/R = a1 + a2*T + a3*T^2 + a4*T^3 + a5*T^4
+#   H/RT = a1 + a2*T/2 + a3*T^2/3 + a4*T^3/4 + a5*T^4/5 + a6/T
+#   S/R  = a1*ln(T) + a2*T + a3*T^2/2 + a4*T^3/3 + a5*T^4/4 + a7
+# Coefficients for CH4(g), low-T range 200-1000 K (NASA Glenn database).
+methane_low = [
+    5.14987613e+00, -1.36709788e-02,  4.91800599e-05,
+    -4.84743026e-08, 1.66693956e-11, -1.02466476e+04, -4.64130376e+00,
+]
+methane_high = [
+    7.48514950e-02,  1.33909467e-02, -5.73285809e-06,
+    1.22292535e-09, -1.01815230e-13, -9.46834459e+03, 1.84373180e+01,
+]
+
+heading("Methane (CH4): a NASA-7 species")
+note(
+    "We define methane as a Component, then evaluate Cp, H, S, "
+    "and G over a sweep of temperatures using the helper functions "
+    "from pythermocalcdb_nasa."
+)
+
+CH4 = Component(name="methane", formula="CH4", state="g")
+
+# In a real project, model_source comes from
+# pyThermoLinkDB.load_and_build_model_source(...) using packaged
+# NASA pickles. The calculation helpers below would then look up CH4
+# by name/formula and pick the right temperature segment automatically.
+note(
+    "In production you'd build <code>model_source</code> with "
+    "<code>load_and_build_model_source(thermodb_sources=...)</code> "
+    "from <code>pyThermoLinkDB</code>, pointing at the packaged NASA "
+    "pickles for each species. Here we just illustrate the call shape."
+)
+
+call_signature = pd.DataFrame({
+    "helper": ["Cp_T", "H_T", "S_T", "G_T"],
+    "returns": [
+        "Heat capacity at constant pressure",
+        "Standard enthalpy H^0(T)",
+        "Standard entropy S^0(T)",
+        "Standard Gibbs energy G^0(T)",
+    ],
+    "units (molar)": ["J/(mol K)", "J/mol", "J/(mol K)", "J/mol"],
+})
+display(call_signature, append=True)
+
+heading("Calling the engine", level=3)
+note(
+    "Each helper takes a Component, a Temperature, and a model_source. "
+    "The example below shows the canonical call; uncomment after "
+    "wiring up a model_source for your species of interest."
+)
+
+example_call = """
+T = Temperature(value=600.0, unit="K")
+Cp = Cp_T(component=CH4, temperature=T, model_source=model_source)
+H  = H_T(component=CH4,  temperature=T, model_source=model_source)
+S  = S_T(component=CH4,  temperature=T, model_source=model_source)
+G  = G_T(component=CH4,  temperature=T, model_source=model_source)
+print(Cp, H, S, G)  # CustomProp objects with .value and .unit
+"""
+display(HTML(f"<pre>{example_call}</pre>"), append=True)

examples/pythermocalcdb-nasa/species_properties/config.toml

@@ -0,0 +1 @@
+packages = ["pythermocalcdb-nasa", "pythermodb-settings", "pythermodb", "pythermolinkdb", "pyreactlab-core"]

examples/pythermocalcdb-nasa/species_properties/setup.py

@@ -0,0 +1,41 @@
+"""
+Shim IPython's display API onto PyScript so example code written in a
+Jupyter/IPython idiom runs unmodified in the browser.
+"""
+
+import sys
+import types
+import js
+from pyscript import window, HTML, display as _display
+
+js.alert = window.alert
+
+
+def display(*args, **kwargs):
+    """Wrap pyscript.display so output lands in the example target."""
+    return _display(
+        *args, **kwargs, target=__pyscript_display_target__,
+    )
+
+
+ipython = types.ModuleType("IPython")
+core = types.ModuleType("IPython.core")
+core_display = types.ModuleType("IPython.core.display")
+core_display.display = display
+core_display.HTML = HTML
+ipython.core = core
+core.display = core_display
+ipython.get_ipython = lambda: None
+ipython.display = core_display
+sys.modules["IPython"] = ipython
+sys.modules["IPython.core"] = core
+sys.modules["IPython.core.display"] = core_display
+sys.modules["IPython.display"] = core_display
+
+
+def heading(text, level=2):
+    display(HTML(f"<h{level}>{text}</h{level}>"), append=True)
+
+
+def note(text):
+    display(HTML(f"<p>{text}</p>"), append=True)