adiabatic-flame.ipynb

{
 "metadata": {
  "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.4"
  },
  "orig_nbformat": 2,
  "kernelspec": {
   "name": "python374jvsc74a57bd0379683f2fa419c33f6058c07b3f6f7ca16fbfcba7a3f6633d97c52aad41ac1dc",
   "display_name": "Python 3.7.4 64-bit ('base': conda)"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2,
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "source": [
    "# Adiabatic flame temperature notebook.\n",
    "\n",
    "Today, we are going to demonstrate the use of Cantera to calculate adiabatic flame temperatures for a range of fuels and conditions, and compare several methods in Cantera to how one would do it without Cantera.\n",
    "\n",
    "## Adiabatic Flame Temperatures\n",
    "\n",
    "As a bit of background, the adiabatic flame temperature is the temperature achieved when a given fuel-oxidizer combination is combustion in a perfectly insulated reactor. If no heat is lost to the surroundings and the reaction occurs at constant pressure, conservation of energy dictates that the reaction proceeds at constant enthalpy.  Namely, the enthalpy at the end of combustion (_state 2_) must equal that at the beginning (_state 1_):\n",
    "\n",
    "$$ h_{\\rm mix, 2} = h_{\\rm mix, 1}$$\n",
    "\n",
    "Wehn we consider the enthalpy of the mixture, there are two components to consider: the enthalpy of formation $\\Delta^\\circ_{\\rm f}h$ and the sensible enthalpy $h_i$.\n",
    "\n",
    "To make this more concrete, let's pick a specific reaction: stoichiometric combustion of methane, CH$_4$ with air::\n",
    "$${\\rm CH_4 + 2\\left(O_2 + 3.76\\,N_2\\right) \\leftrightharpoons 2\\,H_2O + CO_2 + 7.52\\,N_2 }$$\n",
    "\n",
    "From the reaction stoichiometry (note that the net stoichiometric coefficient for N$_2$ is zero), the heat of combustion released by this reaction is:\n",
    "$$ q_{\\rm rxn} = -\\sum_k \\nu_k \\Delta^\\circ_{{\\rm f}}h_k = \\Delta^\\circ_{\\rm f}h_{\\rm CH_4} + 2\\Delta^\\circ_{\\rm f}h_{\\rm O_2} - \\Delta^\\circ_{\\rm f}h_{\\rm CO_2}-2 \\Delta^\\circ_{\\rm f}h_{\\rm H_2O}$$\n",
    "\n",
    "Since the reaction must occur at constant enthalpy, the heat released must be stored in the products. In essence, the heat is used to heat the products to a final temperature:\n",
    "$$q_{\\rm heating} = \\sum_k\\nu_{{\\rm prod},k}\\left[h_k(T_2) - h_k(T_1)\\right] = h_{\\rm CO_2}(T_2) - h_{\\rm CO_2}(T_1) + 2\\left[h_{\\rm H_2O}(T_2) - h_{\\rm H_2O}(T_1)\\right] + 7.52\\left[h_{\\rm N_2}(T_2) - h_{\\rm N_2}(T_1)\\right] $$\n",
    "where, again, states 1 and 2 represent the final and initial states, respectively.  $\\nu_{\\rm prod}$ is the reverse stoichiometric coefficient (i.e. that for all product species). Note that this assumes complete combution - there are no reactant species left to heat.\n",
    "\n",
    "For a purely adiabatic reaction, we will have $q_{\\rm rxn} = q_{\\rm heating}$, and hence:\n",
    "$$\\sum_k\\nu_{{\\rm prod},k}\\left[h_k(T_2) - h_k(T_1)\\right] + \\sum_k \\nu_k \\Delta^\\circ_{{\\rm f}}h_k = 0$$\n",
    "or:\n",
    "$$h_{\\rm CO_2}(T_2) - h_{\\rm CO_2}(T_1) + 2\\left[h_{\\rm H_2O}(T_2) - h_{\\rm H_2O}(T_1)\\right] + 7.52\\left[h_{\\rm N_2}(T_2) - h_{\\rm N_2}(T_1)\\right] + \\Delta^\\circ_{\\rm f}h_{\\rm CO_2}-2 \\Delta^\\circ_{\\rm f}h_{\\rm H_2O} - \\Delta^\\circ_{\\rm f}h_{\\rm CH_4} + 2\\Delta^\\circ_{\\rm f}h_{\\rm O_2} = 0$$\n",
    "\n",
    "Note that only a few variables are dependent on the final temperature $T_2$.  We might be tempted to write the above equation in terms of those variables, but let's hold off for now. Our next step really depends on what method we use to solve the problem. In short, we need a means of calculating the enthalpy values, so that we may determine which value of $T_2$ satisfies our adiabatic condition.  \n",
    "\n",
    "## Approach 1 - Complete combustion, calculating \"by hand\" and using tabulated data.\n",
    "\n",
    "One popular source for these values are the [JANAF tables](https://janaf.nist.gov/), a thermodynamic database maintained by NIST.  These provide a range of thermodynamic properties, including enthalpis of formation, total sensible enthalpy, specific, and others, for a range of temperatures and for a wide variety of species.  For the curious, a collection of these tables are stored in the `data` folder that accompanies this notebook.\n",
    "\n",
    "We will take the enthalpies of formation at the initial condition of 298.15 K:\n",
    "\n",
    "| Species $k\\hspace{3em}$ | $\\Delta^\\circ_{\\rm f}h_k(298.15 K)\\,\\frac{\\rm kJ}{\\rm mol}$|\n",
    "|---|---|\n",
    "|CH$_4$| -74.873 |\n",
    "| O$_2$| 0.0 |\n",
    "| N$_2$| 0.0 |\n",
    "| CO$_2$| -393.522 |\n",
    "| H$_2$O | -241.826|\n",
    "\n",
    "From these, let's solve for the heat released by the reaction.  Give it the variable name `q_rxn`:"
   ],
   "cell_type": "markdown",
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Entahlpies of formation, kJ/mol\n",
    "h_f_CH4 = -74.873\n",
    "h_f_CO2 = -393.522\n",
    "h_f_H2O = -241.826\n",
    "h_f_O2 = 0\n",
    "\n",
    "q_rxn = h_f_CH4 + 2*h_f_O2 - h_f_CO2 - 2*h_f_H2O\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "The heat of reaction is 802.30 kJ/mol.\n"
     ]
    }
   ],
   "source": [
    "print('The heat of reaction is {:.2f} kJ/mol.'.format(q_rxn))"
   ]
  },
  {
   "source": [
    "This heat must be stored in the sensible energy of the products, which reach an elevated temperature $T_2$, after combustion.\n",
    "\n",
    "To calculate the final temeprature, we need to estimate how the species enthalpies change as a function of temperature. The JANAF tables actually give us $h_k(T) - h_k(298.15\\,{\\rm K})$ as a function of $T$, allowing for a high-precision solution.  This is demosntrated in the `adiabatic-flame-detailed.ipynb` notebook in this repo.\n",
    "\n",
    "However, a simpler and nearly as acurrate solution uses the definition of the specific heat to calculate the change in enthalpy: \n",
    "$$ h_k(T_2) - h_k(T_1) = \\int_{T_1}^{T_2}C_p dT$$\n",
    "If $C_p$ is constant, this reduces to:\n",
    "$$ h_k(T_2) - h_k(T_1) = C_p\\left(T_2 - T_1\\right)$$\n",
    "While $C_p$ is not, in fact constant, we can take the value at some intermediate temperature and recognize that the variation in $C_p$ from this average value is not _too_ significant.\n",
    "\n",
    "Finally, to simplify a little bit further, let us recognize that, on a molar basis, our products are $\\frac{7.52}{10.52} = 71.5%$ nitrogen.  So let us consult the JANAF tables once more, and choose the $C_p$ of N$_2$ at $T_{\\rm avg} = 1000$ K (a guess value) as our $C_{p,{\\rm avg}}$ for the mixture.\n",
    "$$C_{p,{\\rm avg}} = 32.697\\,\\frac{\\rm J}{\\rm mol-K}$$\n",
    "\n",
    "Thefefore, $q_{\\rm heating} = C_{p,{\\rm avg}}\\,n_{\\rm tot}\\left(T_2 - T_1\\right)$, where $n_{\\rm tot} = 10.52$, and we can therefore solve directly for $T_2$:\n",
    "\n",
    "$$T_2 = T_1 + \\frac{q_{\\rm rxn}}{n_{\\rm tot}C_{p,{\\rm avg}}}$$"
   ],
   "cell_type": "markdown",
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "Initial T_ad = 2630.61 K.\nT_avg = 1464.38\n"
     ]
    }
   ],
   "source": [
    "T_1 = 298.15 # Initial temperature, K\n",
    "C_p_avg = 32.697*10**(-3) #Convert to kJ/mol-K\n",
    "n_tot = 10.52\n",
    "\n",
    "T_2 = T_1 + q_rxn/n_tot/C_p_avg\n",
    "print('Initial T_ad = {:.2f} K.'.format(T_2))\n",
    "print('T_avg = {:.2f}'.format(0.5*(T_1+T_2)))"
   ]
  },
  {
   "source": [
    "Not bad!  If this is correct, our average temperature would $T_{\\rm avg} = \\frac{298.15 + 2630.61\\,{\\rm K}}{2}= 1464.4$ K.  So our average $C_p$ was taken at a temperature that is perhaps a bit low.  What if we re-did it, using $C_{p,\\,{\\rm N_2}}(1400\\,{\\rm K})=34.518\\,\\frac{\\rm J}{\\rm mol-K}$? \n",
    "\n",
    "Right away, we can see that $C_p$ hasn't changed all _that_ much, in going from 1000. to 1400 K.  But let's check, anyway:"
   ],
   "cell_type": "markdown",
   "metadata": {}
  },
  {
   "source": [
    "C_p_avg = 34.518*10**(-3) #Convert to kJ/mol-K\n",
    "\n",
    "T_2 = T_1 + q_rxn/n_tot/C_p_avg\n",
    "print('Initial T_ad = {:.2f} K.'.format(T_2))\n",
    "print('T_avg = {:.2f}'.format(0.5*(T_1+T_2)))\n"
   ],
   "cell_type": "code",
   "metadata": {},
   "execution_count": 20,
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "Initial T_ad = 2507.56 K.\nT_avg = 1402.85\n"
     ]
    }
   ]
  },
  {
   "source": [
    "We see that the adiabatic flame temperature estimate _did_ drop by 123 K, and now our actual average temperature is quite close to our guess value.\n",
    "\n",
    "Again a still more accurate approach would be to use the species enthalpy values directly, and search the tabulated data to find a $T_2$ that satisfies our adiabatic condition.  This involves a _slightly_ more sophisticated coding approach, and is left as an extension code (`adiabatic-flame-detailed.ipynb`) for you to explore on your own.\n",
    "\n",
    "We will note here that the detailed approach returns an adiabatic flame temperature of $T_{\\rm ad} = 2327$ K, which is not _too_ far off from our rather simple estimate, here!\n",
    "\n",
    "## Approach 2: complete combustion, using Cantera."
   ],
   "cell_type": "markdown",
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ]
}