{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Chemkin Input and Output\n", "This notebook describes pmutt's functionality to read and write Chemkin files. We will use the NH3 formation mechanism as a case study.\n", "\n", "## Topics Covered\n", "- Read species *ab-initio* data, reactions, and catalyst sites from a spreadsheet\n", "- Write the thermdat, gas.inp, surf.inp, T_flow.inp, EAg.inp, EAs.inp, tube_mole.inp files" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Input Spreadsheet\n", "All the data will be imported from the [`./inputs/NH3_Input_Data.xlsx`](https://github.com/VlachosGroup/pmutt/blob/master/docs/source/examples_jupyter/chemkin_io/inputs/NH3_Input_Data.xlsx) file. There are four sheets:\n", "1. `cat_sites` contains catalyst site properties for the adsorbed species\n", "2. `refs` contains *ab-initio* and experimental data for a handful of gas species to calculate references\n", "3. `species` contains *ab-initio* data for each specie\n", "4. `reactions` contains elementary steps" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "First, we change the working directory to the location of the Jupyter notebook." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import os\n", "from pathlib import Path\n", "\n", "# Find the location of Jupyter notebook\n", "# Note that normally Python scripts have a __file__ variable but Jupyter notebook doesn't.\n", "# Using pathlib can overcome this limiation\n", "try:\n", " notebook_path = os.path.dirname(__file__)\n", "except NameError:\n", " notebook_path = Path().resolve()\n", "os.chdir(notebook_path)\n", "excel_path = './inputs/NH3_Input_Data.xlsx'" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Below is a helper function to print tables easily." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "import pandas as pd\n", "from IPython.display import display\n", "\n", "def disp_data(io, sheet_name):\n", " data = pd.read_excel(io=io, sheet_name=sheet_name, skiprows=[1])\n", " data = data.fillna(' ')\n", " display(data)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Catalytic Sites**" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
namesite_densitydensitybulk_specie
0RU00012.167100e-0912.2RU(B)
\n", "
" ], "text/plain": [ " name site_density density bulk_specie\n", "0 RU0001 2.167100e-09 12.2 RU(B)" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "disp_data(io=excel_path, sheet_name='cat_sites')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**References**" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
nameelements.Nelements.Helements.RUT_refHoRT_refpotentialenergysymmetrynumberstatmech_modelatomsvib_wavenumbervib_wavenumber.1vib_wavenumber.2vib_wavenumber.3
0N2200298.150.000000-16.632IdealGas./N2/CONTCAR2744
1NH3130298.15-18.380253-19.543IdealGas./NH3/CONTCAR3534346417651139
2H2020298.150.000000-6.772IdealGas./H2/CONTCAR4342
3Ru001298.150.000000Placeholder
\n", "
" ], "text/plain": [ " name elements.N elements.H elements.RU T_ref HoRT_ref \\\n", "0 N2 2 0 0 298.15 0.000000 \n", "1 NH3 1 3 0 298.15 -18.380253 \n", "2 H2 0 2 0 298.15 0.000000 \n", "3 Ru 0 0 1 298.15 0.000000 \n", "\n", " potentialenergy symmetrynumber statmech_model atoms vib_wavenumber \\\n", "0 -16.63 2 IdealGas ./N2/CONTCAR 2744 \n", "1 -19.54 3 IdealGas ./NH3/CONTCAR 3534 \n", "2 -6.77 2 IdealGas ./H2/CONTCAR 4342 \n", "3 Placeholder \n", "\n", " vib_wavenumber.1 vib_wavenumber.2 vib_wavenumber.3 \n", "0 \n", "1 3464 1765 1139 \n", "2 \n", "3 " ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "disp_data(io=excel_path, sheet_name='refs')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Species**" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
nameelements.Nelements.Helements.RUphasestatmech_modelsymmetrynumberatomspotentialenergyvib_wavenumber...vib_wavenumber.2vib_wavenumber.3vib_wavenumber.4vib_wavenumber.5vib_wavenumber.6vib_wavenumber.7vib_wavenumber.8vib_wavenumber.9vib_wavenumber.10vib_wavenumber.11
0N2200GIdealGas2./N2/CONTCAR-16.632744...
1NH3130GIdealGas3./NH3/CONTCAR-19.543534...17651139
2H2020GIdealGas2./H2/CONTCAR-6.774342...
3N2(S)200SHarmonic-17.242197.19...347.343335.67462.07632.1794
4N(S)100SHarmonic-9.34549.105...504.323475.805459.081410.018
5H(S)010SHarmonic-41003.51...616.29
6NH3(S)130SHarmonic-20.433491.09...3364.521583.521582.071124.22570.212567.221333.09122.85983.828670.6251
7NH2(S)120SHarmonic-16.593469.3...1503.02698.869625.596615.94475.13298.12153.25
8NH(S)110SHarmonic-13.213403.13...710.581528.526415.196410.131
9TS1_NH3(S)130SHarmonic-19.243453.41...1723.851487.95959.151888.946594.089428.431227.032206.047142.136
10TS2_NH2(S)120SHarmonic-15.873426.44...922.831660.967525.595496.837330.674290.278
11TS3_NH(S)110SHarmonic-11.931201.6...462.016402.159242.139
12TS4_N2(S)200SHarmonic-14.67485.614...386.186280.943168.431
13RU(S)001SPlaceholder...
14RU(B)001SPlaceholder...
\n", "

15 rows × 21 columns

\n", "
" ], "text/plain": [ " name elements.N elements.H elements.RU phase statmech_model \\\n", "0 N2 2 0 0 G IdealGas \n", "1 NH3 1 3 0 G IdealGas \n", "2 H2 0 2 0 G IdealGas \n", "3 N2(S) 2 0 0 S Harmonic \n", "4 N(S) 1 0 0 S Harmonic \n", "5 H(S) 0 1 0 S Harmonic \n", "6 NH3(S) 1 3 0 S Harmonic \n", "7 NH2(S) 1 2 0 S Harmonic \n", "8 NH(S) 1 1 0 S Harmonic \n", "9 TS1_NH3(S) 1 3 0 S Harmonic \n", "10 TS2_NH2(S) 1 2 0 S Harmonic \n", "11 TS3_NH(S) 1 1 0 S Harmonic \n", "12 TS4_N2(S) 2 0 0 S Harmonic \n", "13 RU(S) 0 0 1 S Placeholder \n", "14 RU(B) 0 0 1 S Placeholder \n", "\n", " symmetrynumber atoms potentialenergy vib_wavenumber ... \\\n", "0 2 ./N2/CONTCAR -16.63 2744 ... \n", "1 3 ./NH3/CONTCAR -19.54 3534 ... \n", "2 2 ./H2/CONTCAR -6.77 4342 ... \n", "3 -17.24 2197.19 ... \n", "4 -9.34 549.105 ... \n", "5 -4 1003.51 ... \n", "6 -20.43 3491.09 ... \n", "7 -16.59 3469.3 ... \n", "8 -13.21 3403.13 ... \n", "9 -19.24 3453.41 ... \n", "10 -15.87 3426.44 ... \n", "11 -11.93 1201.6 ... \n", "12 -14.67 485.614 ... \n", "13 ... \n", "14 ... \n", "\n", " vib_wavenumber.2 vib_wavenumber.3 vib_wavenumber.4 vib_wavenumber.5 \\\n", "0 \n", "1 1765 1139 \n", "2 \n", "3 347.343 335.674 62.076 32.1794 \n", "4 504.323 475.805 459.081 410.018 \n", "5 616.29 \n", "6 3364.52 1583.52 1582.07 1124.22 \n", "7 1503.02 698.869 625.596 615.94 \n", "8 710.581 528.526 415.196 410.131 \n", "9 1723.85 1487.95 959.151 888.946 \n", "10 922.831 660.967 525.595 496.837 \n", "11 462.016 402.159 242.139 \n", "12 386.186 280.943 168.431 \n", "13 \n", "14 \n", "\n", " vib_wavenumber.6 vib_wavenumber.7 vib_wavenumber.8 vib_wavenumber.9 \\\n", "0 \n", "1 \n", "2 \n", "3 \n", "4 \n", "5 \n", "6 570.212 567.221 333.09 122.859 \n", "7 475.13 298.12 153.25 \n", "8 \n", "9 594.089 428.431 227.032 206.047 \n", "10 330.674 290.278 \n", "11 \n", "12 \n", "13 \n", "14 \n", "\n", " vib_wavenumber.10 vib_wavenumber.11 \n", "0 \n", "1 \n", "2 \n", "3 \n", "4 \n", "5 \n", "6 83.8286 70.6251 \n", "7 \n", "8 \n", "9 142.136 \n", "10 \n", "11 \n", "12 \n", "13 \n", "14 \n", "\n", "[15 rows x 21 columns]" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "disp_data(io=excel_path, sheet_name='species')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Reactions**" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
reaction_stris_adsorption
0H2 + 2RU(S) = 2H(S) + 2RU(B)True
1N2 + RU(S) = N2(S) + RU(B)True
2NH3 + RU(S) = NH3(S) + RU(B)True
3NH3(S) + RU(S)= TS1_NH3(S) = NH2(S) + H(S) + R...False
4NH2(S) + RU(S) = TS2_NH2(S) = NH(S) + H(S) + ...False
5NH(S) + RU(S) = TS3_NH(S) = N(S) + H(S) + R...False
62N(S) + RU(B) = TS4_N2(S) = N2(S) + RU(S)False
\n", "
" ], "text/plain": [ " reaction_str is_adsorption\n", "0 H2 + 2RU(S) = 2H(S) + 2RU(B) True\n", "1 N2 + RU(S) = N2(S) + RU(B) True\n", "2 NH3 + RU(S) = NH3(S) + RU(B) True\n", "3 NH3(S) + RU(S)= TS1_NH3(S) = NH2(S) + H(S) + R... False\n", "4 NH2(S) + RU(S) = TS2_NH2(S) = NH(S) + H(S) + ... False\n", "5 NH(S) + RU(S) = TS3_NH(S) = N(S) + H(S) + R... False\n", "6 2N(S) + RU(B) = TS4_N2(S) = N2(S) + RU(S) False" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "disp_data(io=excel_path, sheet_name='reactions')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Reading data\n", "Before we can initialize our species, we need the catalytic sites and the references.\n", "\n", "### Reading Catalytic Sites" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "{'bulk_specie': 'RU(B)',\n", " 'class': \"\",\n", " 'density': 12.2,\n", " 'name': 'RU0001',\n", " 'site_density': 2.1671e-09}\n" ] } ], "source": [ "import os\n", "from pprint import pprint\n", "from pathlib import Path\n", "from pmutt.io.excel import read_excel\n", "from pmutt.chemkin import CatSite\n", "\n", "cat_site_data = read_excel(io=excel_path, sheet_name='cat_sites')[0]\n", "cat_site = CatSite(**cat_site_data)\n", "\n", "# Print the properties of the catalyst site\n", "pprint(cat_site.to_dict())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Reading reference species" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "{'H': -129.34222830159834, 'N': -320.10077207763885, 'RU': 0.0}\n" ] } ], "source": [ "from pmutt.empirical.references import Reference, References\n", "\n", "references_data = read_excel(io=excel_path, sheet_name='refs')\n", "\n", "# Convert data to Reference objects and put them in a list\n", "refs_list = [Reference(**ref_data) for ref_data in references_data]\n", "\n", "# Assign the Reference objects to a References object so offsets can be calculated\n", "refs = References(references=refs_list)\n", "\n", "# Print out the offsets calculated\n", "print(refs.offset)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Reading species" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "from pmutt.empirical.nasa import Nasa\n", "\n", "# Range of data to fit the Nasa polynomials\n", "T_low = 298. # K\n", "T_high = 800. # K\n", "\n", "species_data = read_excel(io=excel_path, sheet_name='species')\n", "species = []\n", "for specie_data in species_data:\n", " specie = Nasa.from_model(T_low=T_low, T_high=T_high, references=refs, **specie_data)\n", " # If the species is a surface species, assign the catalyst site specified above\n", " if specie.phase.lower() == 's':\n", " specie.cat_site = cat_site\n", " specie.n_sites = 1\n", " species.append(specie)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The warning above is typical when empirical objects are fitting to `StatMech` objects with the `placeholder` preset." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Reading reactions" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "from pmutt import pmutt_list_to_dict\n", "from pmutt.reaction import ChemkinReaction, Reactions\n", "\n", "# Convert list of Nasa polynomials to dictionary of Nasa polynomials\n", "species_dict = pmutt_list_to_dict(species)\n", "\n", "reactions_data = read_excel(io=excel_path, sheet_name='reactions')\n", "reactions_list = []\n", "for reaction_data in reactions_data:\n", " reaction = ChemkinReaction.from_string(species=species_dict, **reaction_data)\n", " reactions_list.append(reaction)\n", "reactions = Reactions(reactions=reactions_list)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Writing Chemkin files\n", "Now that we have all the required objects, we can write the output files. All outputs can be found in the [./outputs folder](https://github.com/VlachosGroup/pmutt/blob/master/docs/source/examples_jupyter/chemkin_io/outputs)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Writing thermdat" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "from pmutt.io.thermdat import write_thermdat\n", "\n", "write_thermdat(filename='./outputs/thermdat', nasa_species=species)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The thermdat file can be return as a string by omitting ``filename``." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "THERMO ALL\n", " 100 500 1500\n", "N2 20200702N 2 G298.0 800.0 523.4 1\n", " 3.75207309E+00-1.30159027E-03 1.83848785E-06-1.17413430E-10-3.65391117E-13 2\n", "-1.67466046E+02 2.02487235E+00 3.41792195E+00 8.91688129E-04-3.46864753E-06 3\n", " 5.44691795E-09-2.46804017E-12-1.27218598E+02 3.46925446E+00 4\n", "NH3 20200702N 1H 3 G298.0 800.0 533.6 1\n", " 3.29427463E+00 2.57481517E-03 5.20883059E-07-1.00678313E-09 3.59356180E-13 2\n", "-8.39739729E+03 3.59518179E+00 4.57124680E+00-6.89724277E-03 2.70972941E-05 3\n", "-3.44331645E-08 1.62568430E-11-8.53630392E+03-1.78236482E+00 4\n", "H2 20200702H 2 G298.0 800.0 574.6 1\n", " 3.37903251E+00 7.92301907E-04-1.88662952E-06 1.85974220E-09-5.68447254E-13 2\n", " 1.70231785E+03-3.75384878E+00 3.50717925E+00-8.50397135E-05 3.79570576E-07 3\n", "-7.57943298E-10 5.72423558E-13 1.68725521E+03-4.30359616E+00 4\n", "N2(S) 20200702N 2 S298.0 800.0 482.4 1\n", " 3.42605836E+00 5.03710708E-03-6.92689694E-06 6.03246113E-09-2.18862675E-12 2\n", "-7.10908593E+03-1.39234158E+01 4.40553502E-01 2.97759950E-02-8.47698148E-05 3\n", " 1.16195177E-07-6.12938589E-11-6.81713036E+03-1.67535850E+00 4\n", "N(S) 20200702N 1 S298.0 800.0 502.9 1\n", "-8.29107669E-01 2.64469843E-02-4.38832016E-05 3.47600781E-08-1.07819310E-11 2\n", "-1.10025322E+04 5.84860624E-01-5.23174593E+00 6.11537184E-02-1.47762309E-04 3\n", " 1.74618161E-07-8.21997526E-11-1.05501788E+04 1.88655390E+01 4\n", "H(S) 20200702H 1 S298.0 800.0 461.9 1\n", "-2.11302416E+00 1.72572991E-02-2.60046478E-05 1.92388088E-08-5.68363784E-12 2\n", "-6.07714366E+03 8.35375655E+00-1.41307531E+00 1.00687718E-02 1.09758988E-06 3\n", "-2.55216914E-08 2.17566538E-11-6.12991563E+03 5.64663049E+00 4\n", "NH3(S) 20200702N 1H 3 S298.0 800.0 461.9 1\n", " 1.58312977E+00 1.57753291E-02-1.71587153E-05 1.15414255E-08-3.21958978E-12 2\n", "-1.42567384E+04-6.35076167E+00 1.16759533E+00 2.05732170E-02-3.66958772E-05 3\n", " 4.56112325E-08-2.49628744E-11-1.42311185E+04-4.80504278E+00 4\n", "NH2(S) 20200702N 1H 2 S298.0 800.0 523.4 1\n", "-6.86530719E-01 2.43298826E-02-3.68705157E-05 2.88483783E-08-8.86963660E-12 2\n", "-1.19475373E+04 1.03191178E+00-2.56305506E+00 3.84726874E-02-7.73503504E-05 3\n", " 8.09903199E-08-3.43595838E-11-1.17458655E+04 8.90806366E+00 4\n", "NH(S) 20200702N 1H 1 S298.0 800.0 543.9 1\n", "-1.33660645E+00 2.33106600E-02-3.77180253E-05 2.97486446E-08-9.12157176E-12 2\n", "-1.48903410E+04 3.58875263E+00-3.73687435E+00 4.05226889E-02-8.44226964E-05 3\n", " 8.65932632E-08-3.53036247E-11-1.46202196E+04 1.37781588E+01 4\n", "TS1_NH3(S) 20200702N 1H 3 S298.0 800.0 451.7 1\n", " 5.21221040E-02 2.11615812E-02-2.49542203E-05 1.66204937E-08-4.59897868E-12 2\n", "-2.50090675E+03-1.43098814E+00 4.96418515E-02 2.15674671E-02-2.75604959E-05 3\n", " 2.23225323E-08-8.79255935E-12-2.50469914E+03-1.46480767E+00 4\n", "TS2_NH2(S) 20200702N 1H 2 S298.0 800.0 554.1 1\n", "-1.48538501E+00 2.67287608E-02-3.86606417E-05 2.81558194E-08-8.15140597E-12 2\n", "-5.88676392E+03 4.23478570E+00-2.58881424E+00 3.43781302E-02-5.87419101E-05 3\n", " 5.18321984E-08-1.87340797E-11-5.75828326E+03 8.95651588E+00 4\n", "TS3_NH(S) 20200702N 1H 1 S298.0 800.0 482.4 1\n", "-1.24242592E-01 1.65859456E-02-2.43704704E-05 1.77594067E-08-5.19878424E-12 2\n", "-2.53350782E+03-1.18098922E+00-1.80701218E+00 3.07479132E-02-6.96816870E-05 3\n", " 8.30020912E-08-4.08170225E-11-2.37120889E+03 5.69776360E+00 4\n", "TS4_N2(S) 20200702N 2 S298.0 800.0 482.4 1\n", " 1.52431744E+00 1.40315985E-02-2.40821843E-05 1.96281604E-08-6.24114074E-12 2\n", " 2.15756709E+04-8.60900953E+00-1.75107618E+00 4.12189598E-02-1.09828935E-04 3\n", " 1.41320046E-07-7.17404084E-11 2.18956564E+04 4.82385770E+00 4\n", "RU(S) 20200702RU 1 S298.0 800.0 554.1 1\n", " 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 2\n", " 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 3\n", " 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 4\n", "RU(B) 20200702RU 1 S298.0 800.0 554.1 1\n", " 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 2\n", " 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 3\n", " 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 4\n", "END\n" ] } ], "source": [ "thermdat_str = write_thermdat(nasa_species=species)\n", "print(thermdat_str)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Writing gas.inp and surf.inp" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [], "source": [ "from pmutt.io import chemkin as ck_io\n", "\n", "ck_io.write_gas(filename='./outputs/gas.inp',\n", " nasa_species=species,\n", " reactions=reactions,\n", " act_method_name='get_G_act',\n", " act_unit='kcal/mol')\n", "ck_io.write_surf(filename='./outputs/surf.inp',\n", " reactions=reactions,\n", " act_method_name='get_G_act',\n", " ads_act_method='get_H_act',\n", " act_unit='kcal/mol')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "Note that `act_method_name` is 'get_G_act'. We use this formalism here since we do not include entropic effects in the preexponential factor.\n", "\n", "Similarly to ``write_thermdat``, the gas.inp and surf.inp file can written as a string by omitting the filename. Note there are no gas-phase reactions." ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "! File generated by pMuTT (v 1.2.21dev) on 2020-07-02 20:55:07.146448\n", "!Elements present in gas and surface species\n", "ELEMENTS\n", "H\n", "RU\n", "N\n", "END\n", "\n", "!Gas-phase species\n", "SPECIES\n", "N2\n", "NH3\n", "H2\n", "END\n", "\n", "!Gas-phase reactions. The rate constant expression is:\n", "!k = kb/h * (T)^beta * exp(-Ea/RT)\n", "!Each line has 4 columns:\n", "!- Reaction reactants and products separated by <=>\n", "!- Preexponential factor, kb/h\n", "!- Beta (power to raise T in rate constant expression)\n", "!- Ea (Activation Energy or Gibbs energy of activation in kcal/mol\n", "REACTIONS\n", "END\n" ] } ], "source": [ "gas_file = ck_io.write_gas(nasa_species=species,\n", " reactions=reactions,\n", " act_method_name='get_G_act',\n", " act_units='kcal/mol')\n", "print(gas_file)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "! File generated by pMuTT (v 1.2.21dev) on 2020-07-02 20:55:07.177427\n", "!Surface species\n", "!Each catalyst site has the following format:\n", "!SITE/[Site name]/ SDEN/[Site density in mol/cm2]/\n", "![Adsorbate Name]/[# of Sites occupied]/ (for every adsorbate)\n", "!BULK [Bulk name]/[Bulk density in g/cm3]\n", "SITE/RU0001/ SDEN/2.16710E-09/\n", "\n", " RU(S)/1/\n", " H(S)/1/\n", " N2(S)/1/\n", " NH3(S)/1/\n", " NH2(S)/1/\n", " NH(S)/1/\n", " N(S)/1/\n", "\n", "BULK RU(B)/12.2/\n", "END\n", "\n", "!Surface-phase reactions.\n", "!The reaction line has the following format:\n", "!REACTIONS MW[ON/OFF] [Ea units]\n", "!where MW stands for Motz-Wise corrections and if the Ea\n", "!units are left blank, then the activation energy should be in cal/mol\n", "!The rate constant expression is:\n", "!k = kb/h/site_den^(n-1) * (T)^beta * exp(-Ea/RT)\n", "!where site_den is the site density and is the number of surface species (including empty sites)\n", "!Each line has 4 columns:\n", "!- Reaction reactants and products separated by =\n", "!- Preexponential factor, kb/h/site_den^(n-1), or \n", "! sticking coefficient if adsorption reaction\n", "!- Beta (power to raise T in rate constant expression)\n", "!- Ea (Activation Energy or Gibbs energy of activation in specified units\n", "!Adsorption reactions can be represented using the STICK keyword\n", "REACTIONS MWON KCAL/MOL\n", "H2+2RU(S)=2H(S)+2RU(B) 5.000E-01 1.000E+00 0.000E+00\n", "STICK\n", "N2+RU(S)=N2(S)+RU(B) 5.000E-01 1.000E+00 0.000E+00\n", "STICK\n", "NH3+RU(S)=NH3(S)+RU(B) 5.000E-01 1.000E+00 0.000E+00\n", "STICK\n", "NH3(S)+RU(S)=NH2(S)+H(S)+RU(B) 9.615E+18 1.000E+00 2.429E+01\n", "NH2(S)+RU(S)=NH(S)+H(S)+RU(B) 9.615E+18 1.000E+00 1.217E+01\n", "NH(S)+RU(S)=N(S)+H(S)+RU(B) 9.615E+18 1.000E+00 2.450E+01\n", "2N(S)+RU(B)=N2(S)+RU(S) 9.615E+18 1.000E+00 8.647E+01\n", "END\n" ] } ], "source": [ "surf_file = ck_io.write_surf(reactions=reactions,\n", " act_method_name='get_G_act',\n", " ads_act_method='get_H_act',\n", " act_unit='kcal/mol')\n", "print(surf_file)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Writing T_flow.inp" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [], "source": [ "# Conditions used to write files\n", "T = [300., 400., 500.] # Temperature in K\n", "P = [1., 2., 3.] # Pressure in atm\n", "Q = [10., 20., 30.] # Standard volumetric flow rate in cm3\n", "abyv= [100., 50., 25.] # Catalyst surface area to reactor volume in 1/cm\n", "\n", "ck_io.write_T_flow(filename='./outputs/T_flow.inp', T=T, P=P, Q=Q, abyv=abyv)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "As shown before, we can return T_flow as a string by omitting the filename." ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "! File generated by pMuTT (v 1.2.21dev) on 2020-07-02 20:55:07.217922\n", "!Conditions for each reaction run\n", "!Only used when MultiInput in tube.inp is set to \"T\"\n", "!T[K] P[atm] Q[cm3/s] abyv[cm-1] Run #\n", "3.000E+02 1.000E+00 1.000E+01 1.000E+02 !1 \n", "4.000E+02 2.000E+00 2.000E+01 5.000E+01 !2 \n", "5.000E+02 3.000E+00 3.000E+01 2.500E+01 !3 \n", "EOF\n" ] } ], "source": [ "T_flow_str = ck_io.write_T_flow(T=T, P=P, Q=Q, abyv=abyv)\n", "print(T_flow_str)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Writing EAg.inp and EAs.inp" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [], "source": [ "# Convert T_flow inputs into list of dictionaries that can be used by write_EA.\n", "# In the future, this will be replaced by a function\n", "conditions = []\n", "for T_i, P_i, Q_i, abyv_i in zip(T, P, Q, abyv):\n", " condition = {\n", " 'T': T_i,\n", " 'P': P_i,\n", " 'Q': Q_i,\n", " 'abyv': abyv}\n", " conditions.append(condition)\n", "\n", "ck_io.write_EA(filename='./outputs/EAs.inp',\n", " reactions=reactions,\n", " write_gas_phase=False,\n", " act_method_name='get_GoRT_act',\n", " ads_act_method='get_HoRT_act',\n", " conditions=conditions)\n", "ck_io.write_EA(filename='./outputs/EAg.inp',\n", " reactions=reactions,\n", " write_gas_phase=True,\n", " act_method_name='get_GoRT_act',\n", " conditions=conditions)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Reminder that we use `act_method_name` as 'get_GoRT_act' for the [reason described above](#act_method_name_explanation)." ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "! File generated by pMuTT (v 1.2.21dev) on 2020-07-02 20:55:07.278920\n", "!The first line is the number of reactions. Subsequent lines follow the format\n", "!of rxn (from surf.out) followed by the EA/RT value at each run condition.\n", "!There may be one slight deviation from surf.out: any repeated species should\n", "!be included in the reaction string with a stoichiometric coefficient equal to\n", "!the number of times the species appears in the reaction. If not using\n", "!MultiInput, then only the first value is used.\n", " 7 !Number of reactions\n", "! 1 2 3\n", "H2+2RU(S)<=>2H(S)+2RU(B) 0.00E+00 0.00E+00 0.00E+00\n", "N2+RU(S)<=>N2(S)+RU(B) 0.00E+00 0.00E+00 0.00E+00\n", "NH3+RU(S)<=>NH3(S)+RU(B) 0.00E+00 0.00E+00 0.00E+00\n", "NH3(S)+RU(S)<=>NH2(S)+H(S)+RU(B) 4.08E+01 3.12E+01 2.55E+01\n", "NH2(S)+RU(S)<=>NH(S)+H(S)+RU(B) 2.04E+01 1.55E+01 1.26E+01\n", "NH(S)+RU(S)<=>N(S)+H(S)+RU(B) 4.11E+01 3.07E+01 2.44E+01\n", "2N(S)+RU(B)<=>N2(S)+RU(S) 1.45E+02 1.09E+02 8.79E+01\n", "EOF\n" ] } ], "source": [ "EAg_file = ck_io.write_EA(reactions=reactions,\n", " write_gas_phase=False,\n", " act_method_name='get_GoRT_act',\n", " conditions=conditions)\n", "print(EAg_file)" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "! File generated by pMuTT (v 1.2.21dev) on 2020-07-02 20:55:07.329254\n", "!The first line is the number of reactions. Subsequent lines follow the format\n", "!of rxn (from surf.out) followed by the EA/RT value at each run condition.\n", "!There may be one slight deviation from surf.out: any repeated species should\n", "!be included in the reaction string with a stoichiometric coefficient equal to\n", "!the number of times the species appears in the reaction. If not using\n", "!MultiInput, then only the first value is used.\n", " 0 !Number of reactions\n", "! 1 2 3\n", "EOF\n" ] } ], "source": [ "EAs_file = ck_io.write_EA(reactions=reactions,\n", " write_gas_phase=True,\n", " act_method_name='get_GoRT_act',\n", " conditions=conditions)\n", "print(EAs_file)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Writing tube_mole.inp" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "\n", "# Generating a list of conditions to input\n", "mole_frac_conditions = []\n", "for x_N2 in np.linspace(0., 0.25, 3):\n", " x_H2 = x_N2*3.\n", " x_NH3 = 1. - x_N2 - x_H2\n", " mole_fractions = {'N2': x_N2, 'H2': x_H2, 'NH3': x_NH3, 'RU(S)': 1.}\n", " mole_frac_conditions.append(mole_fractions)\n", " \n", "# Write the tube_mole.inp file\n", "ck_io.write_tube_mole(mole_frac_conditions=mole_frac_conditions, \n", " nasa_species=species, \n", " filename='./outputs/tube_mole.inp')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Below is the output of the function." ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "! File generated by pMuTT (v 1.2.21dev) on 2020-07-02 20:55:07.376300\n", "!Specify the 'species/phase/' pair /(in quotes!)/ and the associated\n", "!composition values. If the composition does not sum to 1 for each phase or\n", "!site type, it will be renormalized to 1. At the end of a calculation, a\n", "!file containing the complete composition and mass flux (the last entry) will\n", "!be generated. This file's format is completely compatible with the current\n", "!input file and can be used to restart that calculation.\n", "0 itube_restart -- will be >0 if a restart file is used or 0 for the first run\n", "4 Number of nonzero species\n", "! 1 2 3\n", "'N2/GAS/' 0.000 0.125 0.250\n", "'NH3/GAS/' 1.000 0.500 0.000\n", "'H2/GAS/' 0.000 0.375 0.750\n", "'RU(S)/RU0001/' 1.000 1.000 1.000\n", "EOF\n" ] } ], "source": [ "tube_mole_str = ck_io.write_tube_mole(mole_frac_conditions=mole_frac_conditions, \n", " nasa_species=species)\n", "print(tube_mole_str)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.6" } }, "nbformat": 4, "nbformat_minor": 2 }