7_2_Mixture_fractions.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Mixture Fractions\n",
    "\n",
    "This Jupyter Notebook was written by Dan Haworth. It is intended as a tutorial to be used in parallel with Chapter 7 of the book \"An Introduction to Combustion: Concepts and Applications\" by Stephen R. Turns and Daniel C. Haworth. That book is referred to as \"Turns 4th ed.\" throughout this Notebook.\n",
    "\n",
    "This Notebook was last updated 19 February 2020.\n",
    "\n",
    "The objective of this tutorial is to review the basic concept and definition of a mixture fraction, to explain the relationship between mixture fraction and equivalence ratio, and to show how one can define a mixture fraction for a reacting ideal-gas mixture using Cantera's `Solution()` object. We will use the element and species data that are included in gri30.cti.\n",
    "\n",
    "The reader should be familiar with the material in the previous tutorial Notebook (7_1_Molecular_transport_properties) before working through this Notebook."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Concept and definitions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As discussed in Chapter 7 of Turns 4th ed., a *mixture fraction* is an example of a *conserved scalar*: that is, the equation governing the spatial and temporal evolution of a mixture fraction in a combustion system contains no chemical source term. A mixture fraction is a particular conserved scalar that is defined such that its local value at any point in a chemically reacting mixture corresponds to the local mass fraction of mixture that originated from the fuel (versus from the oxidizer). As we will see later in Chapters 9 and 13, mixture fractions are useful in describing nonpremixed (or diffusion) flames. \n",
    "\n",
    "A mixture fraction can be written as a linear combination of species mass fractions. The particular linear combination that is needed is one such that: its value is equal to one in pure fuel; its value is equal to zero in pure oxidizer; and its net chemical source term is equal to zero. \n",
    "\n",
    "For the idealized case of a system that contains three generic species (fuel, oxidizer, and product), an appropriate definition for the local value of the mixture fraction $f$ in terms of the local values of the species mass fractions ($Y_{F}$, $Y_{Ox}$, and $Y_{Pr}$ for fuel, oxidizer, and product, respectively) is:\n",
    "\n",
    "$$ f \\equiv Y_F + \\frac{Y_{Pr}}{\\nu + 1} $$\n",
    "\n",
    "where $\\nu$ is the mass of oxidizer required per unit mass of fuel to form $ 1 + \\nu $ mass units of product. This definition satisfies the three key requirements, since the local fuel, oxidizer, and product mass fractions must sum to one everywhere: $Y_F$ and $Y_{Pr}$ are both equal to zero in pure oxidizer (where $Y_{Ox} = 1$), while $Y_{Ox}$ and $Y_{Pr}$ are both equal to zero in pure fuel (where $Y_{F} = 1$).\n",
    "\n",
    "A more general definition of a mixture fraction for a system that contains an arbitrary number of species $N_{spec}$ is in terms of appropriately normalized elemental mass fractions. As defined in tutorial Notebook 2_1_Mass_moles_molecular_weights_mixtures, the mass fraction of element $ \\alpha $ in a mixture containing $ N_{elem} $ elements and $N_{spec}$ species, denoted as $z_{\\alpha} $, is a linear combination of species mass fractions. Specifically:\n",
    "\n",
    "$$ z_{\\alpha} = W_{\\alpha} \\sum _ {i=1}^{N_{spec}} \\bigl( n_{i \\alpha} Y_i / MW_i \\bigr) $$\n",
    "\n",
    "where the notation is the same as that used in Notebook 2_1_Mass_moles_molecular_weights_mixtures. The equation for $z_{\\alpha}$ contains no chemical source term; element mass fractions are conserved in chemical reactions. It remains, then, to choose an appropriate element (or elements) and normalization such that the value is equal to one in pure fuel and to zero in pure oxidizer. For the case of a hydrocarbon fuel (a fuel molecule that contains only carbon and hydrogen atoms) and our standard approximation for air oxidizer (so that the oxidizer contains only oxygen and nitrogen atoms), it is evident that all carbon and all hydrogen in the system come from the fuel, while all oxygen and nitrogen in the system come from the oxidizer. For example, if reactions occur that eventually form $CO_2$ and $H_2O$, the carbon in the $CO_2$ and the hydrogen in the $H_2O$ originated from the fuel, while the oxygen in both came from the oxidizer. In that case, three possible definitions for a mixture fraction are as follows:\n",
    "\n",
    "$$ f = f_C \\equiv \\frac{z_{C}}{z_{C,fuel}} $$\n",
    "\n",
    "where $z_{C,fuel} $ is the mass fraction of carbon in the fuel molecule;\n",
    "\n",
    "$$ f = f_H \\equiv \\frac{z_{H}}{z_{H,fuel}} $$\n",
    "\n",
    "where $z_{H,fuel} $ is the mass fraction of hydrogen in the fuel molecule; or:\n",
    "\n",
    "$$ f = f _ {CH} \\equiv z_{C} + z_{H} $$\n",
    "\n",
    "where no further normalization is required, since $z_{C+H,fuel} = 1$. A somewhat more general definition is needed in cases where an element is present in both the fuel and in the oxidizer. In that case, a mixture fraction based on the elemental mass fraction for element $ \\alpha $ can be defined as:\n",
    "\n",
    "$$ f = f _ {\\alpha} \\equiv \\frac{z_{\\alpha} - z_{\\alpha,oxid}}{z_{\\alpha,fuel} - z_{\\alpha,oxid}} $$\n",
    "\n",
    "where $ z_{\\alpha,fuel} $ and $ z_{\\alpha,oxid} $ are the elemental mass fractions for element $ \\alpha $ in the fuel and in the oxidizer, respectively.\n",
    "\n",
    "In the case where all species diffuse at equal rates in the gas mixture (no differential diffusion), the local values of all elemental-mass-fraction-based mixture mixtures will be the same. In the more general case where different species diffuse at different rates, the values will be different for different elements. As we saw in tutorial Notebook 7_1_Molecular_transport_properties, different species do, in fact, diffuse at different rates. We will look at this, and at other aspects of mixture fractions, when we get to Chapter 9 (Laminar Diffusion Flames) and Chapter 13 (Turbulent Nonpremixed Flames) of Turns 4th ed. The advantages of introducing mixture fractions should become clearer when we get to that point."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Mixture fractions in Cantera"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's go ahead and define a gas mixture using gri30.cti and Cantera's `Solution()` object:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Running Cantera version: 2.4.0\n"
     ]
    }
   ],
   "source": [
    "# access modules\n",
    "import cantera as ct\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# report Cantera version\n",
    "print(\"Running Cantera version: {}\".format(ct.__version__))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# define the gas mixture\n",
    "gas = ct.Solution('gri30.cti')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As we saw in Notebook 2_1_Mass_moles_molecular_weights_mixtures, the [`elemental_mass_fraction()`](https://cantera.org/documentation/docs-2.4/sphinx/html/cython/thermo.html#cantera.ThermoPhase.elemental_mass_fraction) attribute can be used to access the elemental mass fractions of each element in the gas mixture:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.0, 1.0, 0.0, 0.0)"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# print the elemental mass fractions of carbon, hydrogen, oxygen, and nitrogen in the gas mixture\n",
    "gas.elemental_mass_fraction('C'), gas.elemental_mass_fraction('H'), \\\n",
    "gas.elemental_mass_fraction('O'),gas.elemental_mass_fraction('N')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Recall that the composition of \"gas\" is currently the default composition of 100% H2. So the elemental mass fraction of hydrogen in the gas mixture is equal to one, and the elemental mass fractions of all other elements are equal to zero.\n",
    "\n",
    "For the remainder of this tutorial, we will consider ethane (C2H6) fuel and use our standard approximation for air as our oxidizer. The temperature and pressure of the mixture are not relevant, for our purposes in this tutorial.\n",
    "\n",
    "We reset the gas mixture to correspond to pure fuel, for purposes of finding the elemental mass fractions of carbon and hydrogen in the fuel:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "  gri30:\n",
      "\n",
      "       temperature             300  K\n",
      "          pressure          101325  Pa\n",
      "           density         1.22149  kg/m^3\n",
      "  mean mol. weight         30.0696  amu\n",
      "\n",
      "                          1 kg            1 kmol\n",
      "                       -----------      ------------\n",
      "          enthalpy     -2.7853e+06       -8.375e+07     J\n",
      "   internal energy     -2.8683e+06       -8.625e+07     J\n",
      "           entropy          7633.8        2.295e+05     J/K\n",
      "    Gibbs function     -5.0755e+06       -1.526e+08     J\n",
      " heat capacity c_p          1753.2        5.272e+04     J/K\n",
      " heat capacity c_v          1476.7         4.44e+04     J/K\n",
      "\n",
      "                           X                 Y          Chem. Pot. / RT\n",
      "                     -------------     ------------     ------------\n",
      "              C2H6              1                1         -61.1855\n",
      "     [  +52 minor]              0                0\n",
      "\n"
     ]
    }
   ],
   "source": [
    "gas.TPX = None , None , {'C2H6':1.0}\n",
    "\n",
    "gas()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.7988788691849985, 0.20112113081500146)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# store and print the elemental mass fractions of carbon and hydrogen in the C2H6 fuel\n",
    "z_Cfuel = gas.elemental_mass_fraction('C')\n",
    "z_Hfuel = gas.elemental_mass_fraction('H')\n",
    "\n",
    "z_Cfuel , z_Hfuel"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These are very close to the values that one would get using approximate atomic weights for the elements: $ z_{C,fuel} \\approx 2 \\times 12 / ( 2 \\times 12 + 6 \\times 1 ) = 0.80 $ and $ z_{H,fuel} \\approx 6 \\times 1 / ( 2 \\times 12 + 6 \\times 1 ) = 0.20 $."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's reset the composition to correspond to a stoichiometric mixture of ethane (C2H6) and air, using our standard approximation for air:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "  gri30:\n",
      "\n",
      "       temperature             300  K\n",
      "          pressure          101325  Pa\n",
      "           density         1.17478  kg/m^3\n",
      "  mean mol. weight         28.9198  amu\n",
      "\n",
      "                          1 kg            1 kmol\n",
      "                       -----------      ------------\n",
      "          enthalpy      -1.622e+05       -4.691e+06     J\n",
      "   internal energy     -2.4845e+05       -7.185e+06     J\n",
      "           entropy            6998        2.024e+05     J/K\n",
      "    Gibbs function     -2.2616e+06        -6.54e+07     J\n",
      " heat capacity c_p          1053.8        3.048e+04     J/K\n",
      " heat capacity c_v          766.32        2.216e+04     J/K\n",
      "\n",
      "                           X                 Y          Chem. Pot. / RT\n",
      "                     -------------     ------------     ------------\n",
      "                O2       0.198188         0.219289         -26.2923\n",
      "              C2H6      0.0566251        0.0588766         -64.0568\n",
      "                N2       0.745187         0.721835         -23.3272\n",
      "     [  +50 minor]              0                0\n",
      "\n"
     ]
    }
   ],
   "source": [
    "# set the gas mixture to a stoichiometric mixture of ethane and air\n",
    "fuel = 'C2H6'\n",
    "air  = {'O2':1.0, 'N2':3.76}\n",
    "phi  = 1.0\n",
    "gas.set_equivalence_ratio(phi,fuel,air)\n",
    "\n",
    "# print the gas mixture thermochemical state\n",
    "gas()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the mass fraction of fuel in the stoichiometric fuel-air mixture is approximately 0.0589. By definition, then, this is the value of the mixture fraction in a stoichiometric mixture of C2H6 and air: $ f_{stoich} = 0.0589 $.\n",
    "\n",
    "Now compute three different values for the mixture fraction, as defined in Section 1 above:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.05887663098913928, 0.05887663098913927, 0.05887663098913928)"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# compute and print three different elemental-mass-fraction-based mixture fractions\n",
    "f_C  = gas.elemental_mass_fraction('C') / z_Cfuel\n",
    "f_H  = gas.elemental_mass_fraction('H') / z_Hfuel\n",
    "f_CH = gas.elemental_mass_fraction('C') + gas.elemental_mass_fraction('H')\n",
    "\n",
    "f_C , f_H , f_CH"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "All three values are the same, as they must be.\n",
    "\n",
    "Now reset the gas mixture composition to correspond to equilibrium products of constant-pressure-adiabatic combustion, and recompute the mixture-fraction values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "  gri30:\n",
      "\n",
      "       temperature          2259.6  K\n",
      "          pressure          101325  Pa\n",
      "           density        0.150333  kg/m^3\n",
      "  mean mol. weight         27.8743  amu\n",
      "\n",
      "                          1 kg            1 kmol\n",
      "                       -----------      ------------\n",
      "          enthalpy      -1.622e+05       -4.521e+06     J\n",
      "   internal energy      -8.362e+05       -2.331e+07     J\n",
      "           entropy          9758.2         2.72e+05     J/K\n",
      "    Gibbs function     -2.2212e+07       -6.191e+08     J\n",
      " heat capacity c_p          1488.2        4.148e+04     J/K\n",
      " heat capacity c_v          1189.9        3.317e+04     J/K\n",
      "\n",
      "                           X                 Y          Chem. Pot. / RT\n",
      "                     -------------     ------------     ------------\n",
      "                H2     0.00346831       0.00025083         -25.5386\n",
      "                 H    0.000460984      1.66693e-05         -12.7693\n",
      "                 O    0.000292645      0.000167974         -17.1951\n",
      "                O2     0.00561035       0.00644051         -34.3901\n",
      "                OH     0.00320194       0.00195365         -29.9643\n",
      "               H2O       0.158434         0.102397         -42.7336\n",
      "               HO2    5.98268e-07      7.08427e-07         -47.1594\n",
      "              H2O2    4.76267e-08      5.81184e-08         -59.9287\n",
      "                CO       0.011672         0.011729         -38.4612\n",
      "               CO2      0.0974841         0.153915         -55.6563\n",
      "               HCO    1.14553e-09      1.19255e-09         -51.2305\n",
      "              CH2O    1.64711e-11      1.77427e-11         -63.9998\n",
      "                 N    2.12024e-08      1.06541e-08         -13.8413\n",
      "                NH     3.1152e-09      1.67803e-09         -26.6106\n",
      "               NH2    1.06135e-09       6.1008e-10         -39.3799\n",
      "               NH3    2.44115e-09      1.49149e-09         -52.1492\n",
      "               NNH    9.18595e-10      9.56399e-10          -40.452\n",
      "                NO     0.00225453       0.00242696         -31.0364\n",
      "               NO2    4.35517e-07      7.18806e-07         -48.2315\n",
      "               N2O    1.20246e-07      1.89866e-07         -44.8778\n",
      "               HNO    3.99463e-08       4.4446e-08         -43.8057\n",
      "                CN    1.25965e-13      1.17576e-13         -35.1075\n",
      "               HCN    2.74668e-11      2.66306e-11         -47.8768\n",
      "              HOCN    1.57343e-12      2.42866e-12         -65.0718\n",
      "              HNCO    4.95119e-10      7.64237e-10         -65.0718\n",
      "               NCO    2.47794e-11      3.73519e-11         -52.3025\n",
      "                N2        0.71712         0.720701         -27.6827\n",
      "     [  +26 minor]    4.29325e-16      4.30846e-16\n",
      "\n"
     ]
    }
   ],
   "source": [
    "# reset the gas mixture to correspond to equilibrium products of constant-pressure adiabatic combustion\n",
    "gas.equilibrate('HP')\n",
    "\n",
    "# print the gas mixture thermochemical state\n",
    "gas()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(0.05887663099016989, 0.058876630989412246, 0.05887663099001751)"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# compute and print three different elemental-mass-fraction-based mixture fractions\n",
    "f_C  = gas.elemental_mass_fraction('C') / z_Cfuel\n",
    "f_H  = gas.elemental_mass_fraction('H') / z_Hfuel\n",
    "f_CH = gas.elemental_mass_fraction('C') + gas.elemental_mass_fraction('H')\n",
    "\n",
    "f_C , f_H , f_CH"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Again, the three values are the same, and they are equal to the values in the initial fuel-air mixture. Elemental mass fractions and mixture fractions do not change as a result of chemical reactions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we note that Cantera does have a built-in function that returns the value of an element-mass-fraction based mixture fraction: the [`mixture_fraction()`](https://cantera.org/documentation/docs-2.4/sphinx/html/cython/onedim.html#cantera.CounterflowDiffusionFlame.mixture_fraction) attribute. However, to use this function, it is necessary to tell Cantera what the fuel is and what the oxidizer is, so that it can do the appropriate normalization. The attribute `mixture_fraction()` is not even defined at this point in our Notebook, as you can verify by trying to access it. It is defined only when we define an appropriate flame object, such as a [`CounterflowDiffusionFlame()`](https://cantera.org/documentation/docs-2.4/sphinx/html/cython/onedim.html#counterflowdiffusionflame). We will do that when we get to Chapter 9."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Relationship to equivalence ratio"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The equivalence ratio $ \\Phi $ was defined in Chapter 2 of Turns 4th ed. (see also tutorial Notebook 2_2_Stoichiometry) in terms of the fuel-to-oxidizer mass ratio (fuel-to-air ratio $(F/A)$, for the case of air oxidizer) as:\n",
    "\n",
    "$$ \\Phi \\equiv \\frac{(F/A)}{(F/A)_{stoich}} $$\n",
    "\n",
    "where $(F/A)_{stoich}$ is the stoichiometric fuel-to-air ratio. The local value of the fuel-to-air ratio $(F/A)$ and the local value of the mixture fraction $f$ are related by:\n",
    "\n",
    "$$ (F/A) = \\frac{f}{1-f} $$\n",
    "\n",
    "Note that $(F/A)=0$ for $f$ = 0 (pure air), $(F/A)\\rightarrow \\infty$ for $f$ = 1 (pure fuel), and $ (F/A)_{stoich} = \\frac{f_{stoic}}{1-f_{stoich}}$, where $ f_{stoich} $ is the value of the fuel mass fraction in a stoichiometric fuel-air mixture: $ f_{stoich} = \\bigl( F/(F+A) \\bigr) _ {stoich}$. It then follows that:\n",
    "\n",
    "$$ \\Phi = \\frac{f ( 1 - f_{stoich} )}{f_{stoich} ( 1-f) } $$\n",
    "\n",
    "This is consistent with the limiting cases: $ \\Phi = 0 $ for $ f = 0 $, $ \\Phi \\rightarrow \\infty $ for $ f = 1 $, and $ \\Phi = 1 $ for $ f = f_{stoich} $.\n",
    "\n",
    "Let's plot the relationship between $ \\Phi $ and $f$ for ethane-air reactant mixtures:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# compute equivalence ratio as a function of mixture fraction\n",
    "fmin = 0.\n",
    "fmax = 1.\n",
    "finc = 0.01\n",
    "farr = np.arange(fmin,fmax,finc)\n",
    "nf   = len(farr)\n",
    "phiarr = np.zeros(nf)\n",
    "\n",
    "# f_st is the stoichiometric mixture-fraction value, which we computed earlier\n",
    "f_st = f_C\n",
    "\n",
    "for i, f in enumerate(farr):\n",
    "    phiarr[i] = f*( 1. - f_st ) / f_st / ( 1. - f )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEWCAYAAACXGLsWAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deZxcVZ3//9e7l3T2PXRCEghL2GVtFsetAVlVwvwGHRiXgGhQER3G+aLofAWZcb6oM8PguE2AiKCAyDgSFQVEehDHsO8JSxuykZCFrJ2kl+r+/P44pzo3naquSndXVVfX5/l41KPuPXc753Z1feqcc++5MjOcc8653lSVOgPOOecGPw8WzjnncvJg4ZxzLicPFs4553LyYOGccy4nDxbOOedy8mAxhEnaT1KLpOp+7meWJJNUM1B5Gwwk/UbS3FLnY7AZzOdF0ghJv5S0RdLPinC8H0j6v4U+TjmQ32dRepKWAfVAZyL5VjP7bGlytDtJs4DXgVozS5U2N30j6VrgYDP7SKnzMlRIagR+bGYzinjMjwJXAH8x0J9FSRcDnzCzdw7kfoeKIfVLscx9wMx+V+pMlCNJNeUYxCSJ8IOtq9R5KYU+/t32B17Ntl25fhbKgpn5q8QvYBnw3izLqoF/ATYAS4HLAQNqMm0LXEv4tQcwK70ucCHwZI99XwksjNPvA54BtgIrgWsT683qccxxwC3AGuAN4J+A6rjsYuDRmOdNhBrJOYl9TQR+CKyOy3+RWPZ+4FlgM/C/wNG9nDOL5+I14PWYdmPM+1bgKeBdMf1soB3oAFqA52J6E+GXJIQm2X8AlgPrgNuAcVmOvQR4f2K+Jv59jo/zp8T8bwaeAxoT6zYBXwf+COwEDo7nbCmwLZ6vD/f8W2b5O2TcLkN+rwV+Bvw4rvsCcAhwdSzrSuDMHnlMn5fvA/ckln0DeAgYFfPfFc9pC7AvcCvwT4n1G4FVPT7rXwSeB9riudsX+C9gfSzH57KU42s9/o6XxnPwR+AGYCPhs3gQ8Hvgrfh3+QkwPrGfmcDP4/HeAr4DHA60Emr3LcDmuG7P8nwSaI7HWgjs2+Mz+SnCZ3IT8F1i681QeJU8A/7KGSw+BbwcP+ATgYfpW7AYGb8oZifWfQK4ME43Am8jfGkeDawFzu+5nzj/C+A/4xfGPsDjwGVx2cXxn/mThED3aUJgSDd5/hr4KTABqAXeE9OPJ3xxnRy3mxvLVpflvBjwYDwnI2LaR4BJsbxfAN4Ehvc8L4l9NLHrS/Hj8UvgQGB0/DK5Pcuxvwr8JDH/PuDlOD2d8AV0bjyXZ8T5KYljrgCOjPkcRwhuh8bl04AjM+W5x99zVLbtMuT3WsIX4Vlx29sIX8pfiX+DTxIDbobzMhJ4Nf5d30X48p2R+Mys6nGsW8kdLJ4lfJ5HxHP0VDynw+L5Xwqc1UtZkufkYiBFaJqqifs8OJ73OmAK8Ajw73H9akIAvyGew+HAOxP7ejRbeYDTYvmPj/v+D+CRHp/JXwHjgf0IwejsUn+/DNj3VKkz4K/uf6AWwi/R9OuTcdnvgU8l1j2TPgSLOP9j4KtxejYheIzMkqd/B27ouR9C30ob8Qs6Lr8IeDhOXww0J5aNjNtOJXyhdQETMhzv+8A/9kh7hRhMMqxvwGk5zusm4Jie5yWxvIldX4oPAZ9JLDuUEPRqMuz34OS5I/xyTZ/XL9IjyAD3A3MTx7wusWxU/Hv/VfKcZsozewaLjNtlyO+1wIOJ+Q/Ez1u6Njgm7nd8z/MS508i/JJeDlyUSG+kb8Hi44n5k4EVPfZxNfDDXsrSM1isyLRuYp3zgWfi9NsJX+KZ/q4X03uwuAX4ZmLZ6PgZmZX4TL4zsfxu4Eu95a2cXn411OBxvpmNT7xuiun7EpoJ0pb34xh3EL7YAf6G0AS0A0DSyZIelrRe0hZCjWZyhn3sT/g1ukbSZkmbCbWMfRLrvJmeSO+f8I81E9hoZpuy7PcL6X3G/c4klD+b5HlB0hckLYlXymwm/GrPVIZM9mX3c7ucXcFxN2bWTGiK+oCkkcB5hHObLscHe5TjnYRAuUe+zWw78NeE871G0q8lHZYrs33Ybm1ieiewwcw6E/MQ/kaZjvU44de+CF+A/ZX8u+0P7NvjfH2ZDOc9z/0haR9Jd0l6Q9JWwo+k9OdgJrDc+tavsdtnxMxaCLXG6Yl13kxM7yDLOS1HHiwGvzWED3jafj2Wbyf8ek+b2su+HgAmSzqWEDTuSCy7g9AGO9PMxgE/IHw59LSSULOYnAhsY83syDzKshKYKGl8lmVf7xEwR5rZnb3sz9ITkt5F+FX/IULNZTywJVEG23Pz3awmfHGl7Udo3libeXXuJJzDOcDiGEDS5bi9RzlGmdn1mfINYGb3m9kZhIDyMpD+odDr37aX7QaUpMsJzS6rgauSWciwej6fx+R2KwlNYMnzNcbMzt2LLPbMx/+LaUeb2VhC82T6c7AS2C/LZeB79RmRNIrQ7PnGXuS1bHmwGPzuBj4naYakCcCXeix/FrhQUq2kBuCCbDuKv6buAb5FaOt/MLF4DOFXf6ukkwg1j0z7WEMIOv8qaaykKkkHSXpProLEbX8DfE/ShJjnd8fFNwGfijUcSRol6X2SxuTabyL/KWITg6SvAmMTy9cCsyRl+8zfCVwp6QBJo4F/Bn7ayy/QuwhNgp9m96D7Y0KN4yxJ1ZKGS2qUlPHyUkn1ks6LXzxthOah9C/+Z4F3x/tlxhGaZ/LZbsBIOoTQafwR4KPAVfHHBoRzOinmLe1Z4FxJEyVNBf42xyEeB7ZK+mK8h6Ja0lGSTuxHtscQm3UlTQf+T4/jrQGuj5+x4ZLekSjPDEnDsuz3DuASScdKqiN8Rh4zs2X9yGvZ8GAxePwy3kCXfv13TL+J0Ob9HPA0oeM16f8Srv7YRLha5A56dwfwXuBnPb4IPwNcJ2kbobOxt+aGjxE6IxfH497D7s0svfkooZ33ZUKH9t8CmNmThI7W78R9NhPakPN1PyEQvUpoKmhl9+aJ9A1cb0l6OsP2C4DbCZ2hr8ftr8h2sBj4/gT8BaHDPp2+klDb+DIhcK0kfFll+1+rInTGryb0C7yH8LfAzB6M+36e0An8q3y2Gyjx1/ePgW+Y2XNm9los1+2S6szsZUKQXRqbkPYlnMPnCH0TD5A4N5nEprAPAMcSzvsG4GZCE2JffY3QCb2FcEFF9/9M4ngHEy40WEVozoPQP/gS8KakDRny+hDh/+2/CAHnIMJVhhXBb8orM0PhBjnnXPnxmoVzzrmcPFg455zLyZuhnHPO5eQ1C+ecczkNyYEEJ0+ebLNmzerz9tu3b2fUqFEDl6EyUIllhsosdyWWGSqz3Htb5qeeemqDmU3JtGxIBotZs2bx5JNP9nn7pqYmGhsbBy5DZaASywyVWe5KLDNUZrn3tsySso4Q4c1QzjnncvJg4ZxzLicPFs4553IqWLCQtEDSOkkv9ki/QtIrkl6S9M1E+tWSmuOysxLpZ8e0Zkk9x0VyzjlXBIXs4L6VMM7PbekESacSxs052szaJO0T048gjLFyJGEY4N/FAcwgPG3qDMIYLk9IWmhmiwuYb+eccz0ULFiY2SNxHKOkTwPXm1lbXGddTJ8D3BXTX5fUTHjgCoQH6SwFkHRXXNeDhXPOFVGxL509BHiXpK8TRvX8ezN7gvDwkEWJ9Vax64EiK3ukn5xpx5LmAfMA6uvraWpq6nMmW1pa+rV9OarEMkNllrsSywyVWe6BLHOxg0UN4dnLpwAnAndLOpDMD9kxMvepZByfxMzmA/MBGhoarD/XU/v12JWjEstdiWWGyij3h29eRGeX8cOLT2LEsOoBLXOxg8Uq4OcWBqR6XFIX4XGHq9j9aXAzCOP000u6c865hMeWbiTVZVRXZfr93T/FvnT2F8Bp0P0ErmGEh50sJDztrU7SAcBswhOtngBmx6eXDSN0gi8scp6dc27Q6+wyUl2GBLXVAx8sClazkHQn0Eh45vMq4BrC08gWxMtp24G5sZbxkqS7CR3XKeDy9MPkJX2W8BS0amCBmb1UqDw751y5ak91AVBXU4VURsHCzC7KsugjWdb/OvD1DOn3AfcNYNacc27IaUuFR7DX1VQXZP9+B7dzzg0BbYmaRSF4sHDOuSGgrSMGi1oPFs4557JIN0MNq/Zg4ZxzLotdzVDeZ+Gccy6L7g5ub4ZyzjmXTXefhXdwO+ecy8aboZxzzuXkl84655zLaVefhdcsnHPOZeE1C+ecczmlg8UwDxbOOeeyaetIjw3lwcI551wWfjWUc865nLzPwjnnXE5+B7dzzrmcdt3BXWbNUJIWSFoXn4rXc9nfSzJJk+O8JH1bUrOk5yUdn1h3rqTX4mtuofLrnHPlrL2zfJuhbgXO7pkoaSZwBrAikXwO4bnbs4F5wPfjuhMJj2M9GTgJuEbShALm2TnnylLZjg1lZo8AGzMsugG4CrBE2hzgNgsWAeMlTQPOAh40s41mtgl4kAwByDnnKl2h7+Au2DO4M5F0HvCGmT3X44Hi04GViflVMS1beqZ9zyPUSqivr6epqanP+WxpaenX9uWoEssMlVnuSiwzDP1yv/FmKwCvvbyYpk2vAgNb5qIFC0kjga8AZ2ZanCHNeknfM9FsPjAfoKGhwRobG/uWUaCpqYn+bF+OKrHMUJnlrsQyw9Av94Klj8O69Rx/7NE0HroPMLBlLubVUAcBBwDPSVoGzACeljSVUGOYmVh3BrC6l3TnnHMJQ+YObjN7wcz2MbNZZjaLEAiON7M3gYXAx+JVUacAW8xsDXA/cKakCbFj+8yY5pxzLqFs7+CWdCfwJ+BQSaskXdrL6vcBS4Fm4CbgMwBmthH4R+CJ+LoupjnnnEso9B3cBeuzMLOLciyflZg24PIs6y0AFgxo5pxzbohJXw013O/gds45l017uTZDOeecKx4fSNA551xO6auh/OFHzjnnsirbq6Gcc84Vh5n5Y1Wdc871Lj3ibG21qK7KNPBF/3mwcM65MlfoJijwYOGcc2Wv0MOTgwcL55wre93Dk3uwcM45l033DXkFepYFeLBwzrmy130lVLXXLJxzzmXR3cFdoHGhwIOFc86VvUI/ywI8WDjnXNnzS2edc87lVOhBBMGDhXPOlb3uS2fLsc9C0gJJ6yS9mEj7lqSXJT0v6b8ljU8su1pSs6RXJJ2VSD87pjVL+lKh8uucc+Vq10155dkMdStwdo+0B4GjzOxo4FXgagBJRwAXAkfGbb4nqVpSNfBd4BzgCOCiuK5zzrkoPTZUWTZDmdkjwMYeaQ+YWSrOLgJmxOk5wF1m1mZmrxOexX1SfDWb2VIzawfuius655yLinE1VMGewZ2HjwM/jdPTCcEjbVVMA1jZI/3kTDuTNA+YB1BfX09TU1OfM9bS0tKv7ctRJZYZKrPclVhmGNrlXvx6OwBr17xBU9OG7vSBLHNJgoWkrwAp4CfppAyrGZlrPpZpn2Y2H5gP0NDQYI2NjX3OX1NTE/3ZvhxVYpmhMstdiWWGoV3u5ztfg1de5eADZtHYeGh3+kCWuejBQtJc4P3A6WaW/uJfBcxMrDYDWB2ns6U755yjOAMJ9hosJA0nfLG/C9gX2Am8CPzazF7a24NJOhv4IvAeM9uRWLQQuEPSv8XjzAYeJ9Q4Zks6AHiD0An+N3t7XOecG8q6r4Yq4KWzWYOFpGuBDwBNwGPAOmA4cAhwfQwkXzCz57NsfyfQCEyWtAq4hnD1Ux3woCSARWb2KTN7SdLdwGJC89TlZtYZ9/NZ4H6gGljQlyDlnHNDWTHu4O6tZvGEmV2bZdm/SdoH2C/bxmZ2UYbkW3pZ/+vA1zOk3wfc10s+nXOuopW0GcrMft0zTdJUM3szLl9HqG0455wrocE46qz/wnfOuUGmfRAOJJjpElfnnHMlNBgffnRTQXLhnHOuz0o6kKCk0T3TzOx7udZxzjlXXKUeSPBeSf8q6d2SRqUTJR0o6VJJ97PnQIHOOeeKrBjPs+jtaqjTJZ0LXAa8Q9IEwj0QrwC/Buamr4xyzjlXOsVohur1Dm6/x8E55wY/f6yqc865nHb1WQyeq6Gcc84NMsW4g9uDhXPOlbnum/JqSzM2VDdJxxBGngX4g5k9V7AcOeec2yuD4qY8SZ8nPKRon/j6saQrCpYj55xzeUt1dpHqMiSorS7cIBv51CwuBU42s+0Akr4B/An4j4LlyjnnXF7aO3d1bsdHPxREPnUWAZ2J+U58jCjnnBsUinH3NuRXs/gh8Jik/47z59PLcymcc84VTzHu3oY8ahZm9m/AJcBGYBNwiZn9e67tJC2QtE7Si4m0iZIelPRafJ8Q0yXp25KaJT0v6fjENnPj+q/F53c755yLinH3NvQ+kODY+D4RWAb8GLgdWB7TcrmVPceO+hLwkJnNBh6K8wDnEJ67PRuYB3w/cexrgJOBk4Br0gHGOedcce7eht5rFnfE96eAJxOv9HyvzOwRQm0kaQ7wozj9I0KTVjr9NgsWAeMlTQPOAh40s41mtgl4EB+80DnnuhXj7m3ofSDB98f3AwbwePVmtibud018jjfAdGBlYr1VMS1b+h4kzSPUSqivr6epqanPmWxpaenX9uWoEssMlVnuSiwzDN1yv7YpNEO1bt+zfANZ5pwd3JIeMrPTc6X1U6arq6yX9D0TzeYD8wEaGhqssbGxz5lpamqiP9uXo0osM1RmuSuxzDB0yz2seQM89hhTJo2nsfHtuy0byDL31mcxPPYZTJY0IXZOT5Q0C9i3j8dbG5uXiO/rYvoqYGZivRnA6l7SnXPOMTj6LC4j9E8cFt/Tr3uB7/bxeAuB9BVNc+O+0ukfi1dFnQJsic1V9wNnxmA1ATgzpjnnnKM4gwhC730WNwI3SrrCzPb6bm1JdwKNhJrJKsJVTdcDd0u6FFgBfDCufh9wLtAM7CBcqouZbZT0j8ATcb3rzKxnp7lzzlWstiIMIgh59FmY2X9IOgo4AhieSL8tx3YXZVm0R1+HmRlweZb9LAAW5Mqnc85VopJfDZUm6RpCDeEIQg3gHOBRoNdg4ZxzrvCK1QyVz94vINQG3jSzS4BjgLqC5so551xeBkMHd9pOM+sCUvGu7nXAgQXNlXPOubzs6rMocTMU8KSk8cBNhKuhWoDHC5or55xzeSnGg48gR7BQGBz9/5nZZuAHkn4LjDWz5wuaK+ecc3kp+UCC0H2V0i8S88s8UDjn3OBRrOdZ5BOKFkk6saC5cM451yfFep5FPn0WpwKXSVoObCeM12RmdnRBc+accy6nkt/BnXBOQXPgnHOuzwbTHdzLC5oD55xzfVasO7gLu3fnnHMFNZju4HbOOTdIdd9nMRiChaT9Jb03To+QNKaguXLOOZeX9sEy3IekTwL3AP8Zk2aQuPfCOedc6RTr0tl89n458A5gK4CZvQbs0+sWzjnniiLdZzG8lHdwp/NiZu3pGUk1ZHkOtnPOueIaTHdw/4+kLwMjJJ0B/Az4ZX8OKulKSS9JelHSnfF53wdIekzSa5J+KmlYXLcuzjfH5bP6c2znnBtKBlMz1JeA9cALhOdy3wf8Q18PKGk68DmgwcyOAqqBC4FvADeY2WxgE3Bp3ORSYJOZHQzcENdzzjlH8tLZ0tcsRgALzOyDZnYB4RGnI/p53BpCTaUGGAmsAU4jdKQD/Ag4P07PifPE5afH0XCdc67iDabnWTwEvJfwHAsIgeIB4C/6ckAze0PSvwArgJ1xX08Bm80sFVdbBUyP09OBlXHblKQtwCRgQ3K/kuYB8wDq6+tpamrqS/YAaGlp6df25agSywyVWe5KLDMMzXKbWfels//76CNU9fgdPZBlzidYDDezdKDAzFokjezrASVNINQWDgA2E/pAMo0/le5Ez1SL2KOD3czmA/MBGhoarLGxsa9ZpKmpif5sX44qscxQmeWuxDLD0Cx3W6oT7v8ttdXitFNP3WP5QJY5n3rLdknHp2cknUCoEfTVe4HXzWy9mXUAPyfUUsbHZikI93KsjtOrgJnx2DXAOGBjP47vnHNDwtadoTFmVF0+v/v7J58j/C3wM0npL+9pwF/345grgFNi7WQncDrwJPAwcAFwFzAXuDeuvzDO/yku/318KJNzzlW09dvaAJgyuq7gx8pn1NknJB0GHEpoEno51gj6xMwek3QP8DSQAp4hNB/9GrhL0j/FtFviJrcAt0tqJtQoLuzrsZ1zbijZ0BKDxZhBECyiE4FZcf3jJGFmt/X1oGZ2DXBNj+SlwEkZ1m0FPtjXYznn3FDVXbMYDMFC0u3AQcCzQGdMNqDPwcI551z/rY81i8mDoRkKaACO8H4C55wbXDYUsWaRz9VQLwJTC50R55xzeyddsxgUHdzAZGCxpMeBtnSimZ1XsFw555zLKd1nMXkw9FkA1xY6E8455/behsFUszCz/5G0PzDbzH4X748o7IhVzjnncirm1VB9eVLedPxJec45V1IdnV1s2tFBlWDiqGEFP54/Kc8558rQWy3hmXQTR9VRXVX4gbj9SXnOOVeGitkEBSV6Up5zzrn+2dB9Q17hm6CgBE/Kc84513/FrlnkczVUF3BTfDnnnBsE1hdxEEHoJVhIeoFe+ibM7OiC5Mg551xOxRyeHHqvWby/KDlwzjm31wZNzcLMlhclB8455/ZasWsW+dyUd4qkJyS1SGqX1ClpazEy55xzLrPuq6EG0aWz3wEuAl4DRgCfAP6jPweVNF7SPZJelrRE0tslTZT0oKTX4vuEuK4kfVtSs6Tnk88Dd865SjXoahYAZtYMVJtZp5n9EDi1n8e9EfitmR0GHAMsIVyi+5CZzQYeivMA5wCz42se8P1+Hts558paa0cn21pT1FaLcSNqi3LMfEad3SFpGPCspG8Ca4BRfT2gpLHAu4GLAeLd4e2S5gCNcbUfAU3AF4E5wG3x4UuLYq1kmpmt6WsenHOunKWboCaNqqOqCEN9QH7B4qOEGshngSuBmcBf9eOYBxJu8vuhpGOAp4DPA/XpAGBmaySlx5+aDqxMbL8qpu0WLCTNI9Q8qK+vp6mpqc8ZbGlp6df25agSywyVWe5KLDMMrXL/eXN4wvVw2nst00CWOZ9gcTxwn5ltBb42QMc8HrjCzB6TdCO7mpwyyRQ297j/w8zmA/MBGhoarLGxsc8ZbGpqoj/bl6NKLDNUZrkrscwwtMrd/tKbsOgpDtx3Mo2NJ2ZdbyDLnE+fxXnAq5Jul/S+OJBgf6wCVpnZY3H+HkLwWCtpGkB8X5dYf2Zi+xnA6n7mwTnnytaGOOJsscaFgjyChZldAhxMGEDwb4A/S7q5rwc0szeBlZIOjUmnA4uBhcDcmDYXuDdOLwQ+Fq+KOgXY4v0VzrlKVuxxoSC/ZijMrEPSbwjNPyMInc6f6MdxrwB+EjvOlwKXEALX3ZIuBVYAH4zr3gecCzQDO+K6zjlXsYr5ONW0nMFC0tnAhYTLZZuAm4EP9eegZvYs0JBh0ekZ1jXCA5icc86xq2ZRrBvyIL+axcXAXcBlZtZW2Ow455zLZf1grFmY2YXFyIhzzrn8bCjyIILQSwe3pEfj+zZJW3u+Fy2HzjnndjOomqHM7J3xfUzRcuOcc65X29tS7GjvpK6mijF1/b2TIX/5dHDvlyndzFYMfHacc871ZvXmnQDUjx2OVJyhPiC/Du5fJ6aHAwcArwBHFiRHzjnnslq8JvQCHDq1uI0++XRwvy05H4cIv6xgOXLOOZdVOlgcPm1sUY+b1xDlSWb2NJB9MBLnnHMFs2TNNgCOmDbIahaS/i4xW0UYx2l9wXLknHMuqyWxZnHEtHFFPW4+fRbJ8JUi9GH8V2Gy45xzLpv129pYv62N0XU1zJgwoqjHzqfPYiCGJXfOOddP6VrFYVPHFO2hR2n5NEMt7G25mZ03cNlxzjmXzZISdW5Dfs1QrwNTgR/H+YuAZcD9BcqTc865DLr7K/YdnMHiODN7d2L+l5IeMbMvFypTzjnn9lSqy2Yhv0tnp0g6MD0j6QBgSuGy5JxzrqfWjk7+vH47VYJD64s/ClM+NYsrgSZJS+P8LPymPOecK6rmdS10dhkHTRnFiGHVRT9+PldD/VbSbOCwmPTyQDzXQlI18CTwhpm9P9ZY7gImAk8DHzWzdkl1wG3ACcBbwF+b2bL+Ht8558rJ4tWla4KC3ocovyoxe56ZPRdfbZL+eQCO/XlgSWL+G8ANZjYb2ARcGtMvBTaZ2cHADXE955yrKKXsr4De+yySDz26useys/tzUEkzgPcRHtGKwtCJpwH3xFV+BJwfp+fEeeLy01XMoRadc24QKOWVUNB7M5SyTGea31v/DlzFrrvDJwGbzSwV51cB0+P0dGAlgJmlJG2J62/YLUPSPGAeQH19PU1NTX3OXEtLS7+2L0eVWGaozHJXYpmhvMttZrywcgcAm19/kaY1+Q3rN5Bl7i1YWJbpTPN5k/R+YJ2ZPSWpMZ3cy/F7W7YrwWw+MB+goaHBGhsbe66St6amJvqzfTmqxDJDZZa7EssM5V3ulRt3sOP+h5k4ahjnn3Vq3s+xGMgy9xYsjomPTxUwIvEoVRGea9FX7wDOk3Ru3M9YQk1jvKSaWLuYAayO668CZgKrJNUA44CN/Ti+c86VlUebQ0PKcTPHF/WBR0lZ6zJmVm1mY81sjJnVxOn0fG1fD2hmV5vZDDObRegX+b2ZfRh4GLggrjYXuDdOL4zzxOW/N7M+12ycc67cPLh4LQBnHFFfsjzs9fMsCuiLwN9Jaib0SdwS028BJsX0vwO+VKL8Oedc0W1vS/Fo8wYkOO3wfUqWj+I97TsDM2sCmuL0UuCkDOu0Ah8sasacc26Q+MNr62lPdXHcfuPZZ0x/egD6ZzDVLJxzzvXwwCBoggIPFs45N2ilOrv4/cvrADjTg4VzzrlMnly+ic07Ojhg8igOmjK6pHnxYOGcc4NU+iqo9x6+T8kumU3zYOGcc4OQmSUumZ1a4tx4sHDOuUHp5Te3sWLjDiaOGsYJ+08odXY8WDjn3GB06x+XAXDu26ZSXVX6sVM9WDjn3CCzdmsr//3MG0jwiXcemHuDIvBg4Zxzg8yCP75Oe2cX52Z/j0wAABHMSURBVBw1lVmTR5U6O4AHC+ecG1S2tnZwx6IVAFz27oNKnJtdPFg459wgcudjK9jWluKUAydyzMzxpc5ONw8Wzjk3SLSlOlnwx9cB+NR7Bk+tAjxYOOfcoHHb/y5n7dY2Dps6hvccMqXU2dmNBwvnnBsElr+1nX998BUArjr70JLfsd2TBwvnnCsxM+Pqn79Aa0cXc47dl9MOK+2ggZl4sHDOuRK7+8mV/O+f32LCyFq++v4jSp2djIoeLCTNlPSwpCWSXpL0+Zg+UdKDkl6L7xNiuiR9W1KzpOclHV/sPDvnXKGs3drKP/16CQDXnnckk0bXlThHmZWiZpECvmBmhwOnAJdLOoLwuNSHzGw28BC7Hp96DjA7vuYB3y9+lp1zbuDtaE/xydueZFtrilMPncJ5x+xb6ixlVfRgYWZrzOzpOL0NWAJMB+YAP4qr/Qg4P07PAW6zYBEwXtK0ImfbOecGVGeX8fm7nuX5VVuYOXEE3/rgMYOuUztJZla6g0uzgEeAo4AVZjY+sWyTmU2Q9CvgejN7NKY/BHzRzJ7ssa95hJoH9fX1J9x11119zldLSwujR5f2QSPFVollhsosdyWWGQZfue9Y0sYDy1OMqoWvnDyCfUcP/G/3vS3zqaee+pSZNWRaVjNgudpLkkYD/wX8rZlt7SWiZlqwR4Qzs/nAfICGhgZrbGzsc96amproz/blqBLLDJVZ7kosMwyecpsZ336omQeWv0pttbjlkpM55cBJBTnWQJa5JMFCUi0hUPzEzH4ek9dKmmZma2Iz07qYvgqYmdh8BrC6eLl1zrmB0dllXLvwJW5ftBwJvnXBMQULFAOtFFdDCbgFWGJm/5ZYtBCYG6fnAvcm0j8Wr4o6BdhiZmuKlmHnnBsArR2dXHHn09y+aDnDaqr4/oeP5/zjppc6W3krRc3iHcBHgRckPRvTvgxcD9wt6VJgBfDBuOw+4FygGdgBXFLc7DrnXP80r2vhc3c+w+I1WxkzvIabP9bAyWVSo0grerCIHdXZOihOz7C+AZcXNFPOOVcAZsZPn1jJ1365mJ0dnew3cST/+dETOHza2FJnba+VrIPbOeeGsqXrW7juV4tpemU9AH953HSum3MkY4bXljhnfePBwjnnBlBLW4rv/L6ZWx5dSkenMaauhuvOP5K/PG5GqbPWLx4snHNuAGxr7eC2Py3n5j8sZdOODgA+1DCDq84+jMmDdAiPveHBwjnn+mH15p3c8dgKbvvTMra2pgA4Yf8J/MP7Due4/SaUNnMDyIOFc87tpVRnF39o3sCdj63gd0vW0hVvEz75gIl8/vTZvP2gSYN66I6+8GDhnHN56Ooynlm5mV8+t5pfPb+aDS3tANRUifcdPY25b9+fhlkTS5zLwvFg4ZxzWWxr7WDR0o08tGQtv1uyjg0tbd3LDpwyiv/vuOl86MSZ7DNmeAlzWRweLJxzLtrRnuKZFZt5/PWN/LF5A8+s3Exn166h6GZMGMFZR07l/GOnc9T0sUOuqak3HiyccxWpq8t4/a3tPL9qM8+u2MyzKzfz0uqtpBLBobpKnLD/BN49ewpnHlnPYVPHVFSASPJg4Zwb8jbvaOeVjZ2sXLScV9/cxuI1W1myZis72jt3W69KcNT0sZw4ayKnHDiJtx80ibFlehPdQPNg4ZwbElraUqzcuIPlb21n2Vs7WLZhO0vXb+fP61t4a3t7XOvF3baZOnY4b5sxjmNnjufYmeM5esa4sr3DutA8WDjnBr1UZxfrW9p4c0srb25pZfWWVlZv3snqzTt5Y/NOVm7c0X0jXCYjaquZOsI4Yfa+HFI/miOmjePwaWMG7fOuByMPFs65kmhPdbFpRztvtbSzcXs7b21v462Wdja0tLF+WxvrW9pYt7WNddvaeGt7G7ke6llXU8WMCSOYNWkU+08axf6TRnLglFEcNGU0U8cO55FH/ofGxmOKU7ghyIOFc67POjq72NaaYltrB9taU2xt7WDrzhRbd3awJfHanH7f0c6mHe1s2t5BS1sq7+NIMHl0HdPGDWfquOFMGzec6eNHsG98zZwwgsmj66iqqszO52LwYOFchejo7GJnRyebW7tYtmE7O9o72dmRYkd7Z5iO7zvaU2xv62R7e4qWthTb4ytMd7K9LcW2thAgWju6+pyf6ioxYeQwJo0axoRRtUweXcekUcOYNLqOKWPqmBLf68cOZ/LoYdRUF/1ZbS7Bg4VzRdLVZbR3doVXatero7OLtlR4tafC8raOzvi++3xrRxdtqU7aOrpo7X7vorWjk9aOXelhPgSH1vZOWlOddHQm2nGamgakTFWCMcNrGTO8hrHxfczwWsaN2PUaPzK8xo6oZcLIYUwYWcv4kcMYU1fjNYEyUjbBQtLZwI1ANXCzmV1f4iy5AjEzUl1GZ3yluoyurl1pqa6uLOlGqrNrt/nOri46OvdcnuoMy5Ys6+DVR/5MR2dIS3WvH95TXV2kOo2OTqOjs6t7eUdnV3zF9Pjenk5PhW3TX/7pPJRSlWDksBqqLMW40SMYWVvDiGHVjBxWzchhNYwcVs2ouuR0DaPistHDaxhdV8OoupoQEOL0yGHVFXvfQaUpi2AhqRr4LnAGsAp4QtJCM1s8kMdpS3WyrTXF5rYu1m5txQy6zOgyw4zwwuiK6WbJaXZbL71dl7Hbesl1d9tPF7ut37nbtruWmxGXJdbt2nM6vd/Orl3HTa7X1WPZylVt3L/xBbq64rG7duWhqyu97a79dCbSd1s/vnd2kSFt1yu9n/QXfmdiPldH5oB7+eWiHEaCYdVV4VUTXrXp6URaXXyF6eru6eG11d3rDa8N83U1u97raqsZXlPdvWx4bQgEYTocQxJNTU00NjYWpcxu6CiLYAGcBDSb2VIASXcBc4ABDRYPLl7LZ+94Jsw8/NBA7ro8rFxR6hx0q64S1RLVVaKmSlRVidpq7UqvFjVVVVQJaqurQnpct6YqzNdUh/nqqqpd28bl6fm1a1Yza7+Z1FSHtJqqqu7teqbVVova6qr4Su+nimE1YTp8+Senw3w6GFRXyX+Fu7JVLsFiOrAyMb8KODm5gqR5wDyA+vp6mvrQJvvauhRjhgFmVFVVIULVHcJDw6Vd71U95sN7WLkqLqfHelXqMU2YqUpvn3GdOC+oQuE9uf/Eut3HkhLTIb1qj2ntdryO9jZG1NXt2r+gKrHObsfpZZ3kS0C1duW5OrlMYVnPbap7nMv89L2TtaW2g9Gj1+25wIBUfPVy1Lb4KictLS19+v8od5VY7oEsc7kEi0zfHLs1VpjZfGA+QENDg/Wlmt0IXPkhKrKaXollhsosdyWWGSqz3ANZ5nK5Fm0VMDMxPwNYXaK8OOdcxSmXYPEEMFvSAZKGARcCC0ucJ+ecqxhl0QxlZilJnwXuJ1w6u8DMXipxtpxzrmKURbAAMLP7gPtKnQ/nnKtE5dIM5ZxzroQ8WDjnnMvJg4VzzrmcPFg455zLSVb0gXgKT9J6YHk/djEZ2DBA2SkXlVhmqMxyV2KZoTLLvbdl3t/MpmRaMCSDRX9JetLMGkqdj2KqxDJDZZa7EssMlVnugSyzN0M555zLyYOFc865nDxYZDa/1BkogUosM1RmuSuxzFCZ5R6wMnufhXPOuZy8ZuGccy4nDxbOOedyqthgIelsSa9Iapb0pQzL6yT9NC5/TNKs4udy4OVR7r+TtFjS85IekrR/KfI5kHKVObHeBZJM0pC4vDKfckv6UPx7vyTpjmLncaDl8fneT9LDkp6Jn/FzS5HPgSRpgaR1kl7MslySvh3PyfOSju/Tgcys4l6EYc7/DBwIDAOeA47osc5ngB/E6QuBn5Y630Uq96nAyDj96XIvdz5ljuuNAR4BFgENpc53kf7Ws4FngAlxfp9S57sIZZ4PfDpOHwEsK3W+B6Dc7waOB17Msvxc4DeEJ46eAjzWl+NUas3iJKDZzJaaWTtwFzCnxzpzgB/F6XuA07V3D4YejHKW28weNrMdcXYR4amE5SyfvzXAPwLfBFqLmbkCyqfcnwS+a2abAMwsw8PIy0o+ZTZgbJwexxB44qaZPQJs7GWVOcBtFiwCxkuatrfHqdRgMR1YmZhfFdMyrmNmKWALMKkouSucfMqddCnhF0k5y1lmSccBM83sV8XMWIHl87c+BDhE0h8lLZJ0dtFyVxj5lPla4COSVhGej3NFcbJWUnv7f59R2Tz8aIBlqiH0vIY4n3XKTd5lkvQRoAF4T0FzVHi9lllSFXADcHGxMlQk+fytawhNUY2EGuQfJB1lZpsLnLdCyafMFwG3mtm/Sno7cHssc1fhs1cyA/JdVqk1i1XAzMT8DPasjnavI6mGUGXtrapXDvIpN5LeC3wFOM/M2oqUt0LJVeYxwFFAk6RlhDbdhUOgkzvfz/i9ZtZhZq8DrxCCR7nKp8yXAncDmNmfgOGEwfaGsrz+73Op1GDxBDBb0gGShhE6sBf2WGchMDdOXwD83mJvURnLWe7YJPOfhEBR7m3YkKPMZrbFzCab2Swzm0XopznPzJ4sTXYHTD6f8V8QLmhA0mRCs9TSouZyYOVT5hXA6QCSDicEi/VFzWXxLQQ+Fq+KOgXYYmZr9nYnFdkMZWYpSZ8F7idcQbHAzF6SdB3wpJktBG4hVFGbCTWKC0uX44GRZ7m/BYwGfhb781eY2Xkly3Q/5VnmISfPct8PnClpMdAJ/B8ze6t0ue6fPMv8BeAmSVcSmmIuLvcfgZLuJDQlTo59MdcAtQBm9gNC38y5QDOwA7ikT8cp8/PknHOuCCq1Gco559xe8GDhnHMuJw8WzjnncvJg4ZxzLicPFs4553LyYOHKVhwh9vbEfI2k9ZJ+FefP622U2bjOxZL2LWAePydpiaSfDMC+dsurpJslHdHf/cZ9nS/pq1mW/U7ShIE4jitfHixcOdsOHCVpRJw/A3gjvdDMFprZ9Tn2cTGwV8Ei3tGfr88A55rZh/uxj7SLSeTVzD5hZov7sJ9MrgK+l2XZ7YRyuArmwcKVu98A74vTFwF3phfEX+LfidP3SvpYnL5M0k8kXUAY/+onkp6VNELSsng3M5IaJDXF6WslzZf0AHCbpGpJ35L0RHxGwGU9MybpB4ThshdKujLDPmZJ+oOkp+PrLxLbXiXpBUnPSbo+S16b0sOSSLoorv+ipG8k9tMi6etxP4sk1WfI5yFAm5ltyHKOF8Zz6yqYBwtX7u4CLpQ0HDgaeCzLevOAr0p6F+Eu3ivM7B7gSeDDZnasme3McawTgDlm9jeEMYa2mNmJwInAJyUdkFzZzD5FGIPnVDO7IcM+1gFnmNnxwF8D3waQdA5wPnCymR0DfLO3vMamqW8ApwHHAidKOj8uHgUsivt5hDAseU/vAJ7OVug4hHmdpHIfddn1Q0UO9+GGDjN7XuEphhcRhjXItt7a2Cb/MPCXZtaXQSEXJr6kzwSOjr/4IQw0ORt4fS/2UQt8R9KxhOE2Donp7wV+mH6uSB55PRFoMrP1ALF/5N2EsZ/agfTQ608Rmup6mkbu8ZHWEZrAynY4ENc/HizcULAQ+BfC+Di9/fp9G+HLrrc+ihS7atzDeyzbnpgWoXZy/17ldPd9XAmsBY6Jx0w/eEns3RDSvT2UqyMx9lEnmf/ndxKCHZKqCUEFQmBLd3oPj+u5CuXNUG4oWABcZ2YvZFtB0knAOcBxwN8nmoy2EYYpT1tGaCoC+Ktejnk/8GlJtXH/h0gatZf5Hgesic9S+Chh8DuAB4CPSxoZ9z0xS17THgPeI2ly/LK/CPifvcjHEuBgADPrjM1cx6YDhcKIklMJ58ZVKA8WruyZ2SozuzHbckl1wE3Ax81sNaHPYkH8ErwV+EG60xj4GnCjpD8QfolnczOwGHha0ouEYd33tqb+PWCupEWEJqjtsTy/JdSWnpT0LPD3cf2eeU2Xfw1wNaGJ7TngaTO7dy/y8QhwXDwfmZxA6PdI7cU+3RDjo84655B0I/BLM/tdlmULzeyh4ufMDRZes3DOAfwzMDLLshc9UDivWTjnnMvJaxbOOedy8mDhnHMuJw8WzjnncvJg4ZxzLicPFs4553L6/wGirf4X2LZudQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot equivalence ratio as a function of mixture fraction\n",
    "fig, ax = plt.subplots()\n",
    "ax.plot(farr, phiarr, label='Phi',   lw=2)\n",
    "ax.set(xlabel='Mixture fraction (-)',\n",
    "       ylabel='Equivalence ratio (-)',\n",
    "       title ='Equivalence ratio versus mixture fraction')\n",
    "#ax.legend(handlelength=4,loc='best')\n",
    "ax.grid(True)\n",
    "fig.savefig(\"Phi_vs_f.pdf\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This scaling obscures the behavior for stoichiometric and fuel-lean reactant mixtures, but it serves to make one important point regarding mixture fractions versus equivalence ratios: mixture fraction is limited to the range 0-1, while equivalence ratio ranges from 0 to $\\infty$. \n",
    "\n",
    "Replot for a narrower range of equivalence ratios around $ \\Phi = 1 $:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZMAAAEWCAYAAACjYXoKAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3hUVfrA8e+bQkKvoZeAhNBsgGA3IF1F3fW3Yq9rW1372rEr9nUtq6isa191V6UpYInYKSo9CaGHFjoJJW3e3x/nxgwxZWAymUzyfp4nT+b2955M5p17z7nniKpijDHGBCMq3AEYY4yJfJZMjDHGBM2SiTHGmKBZMjHGGBM0SybGGGOCZsnEGGNM0CyZ1DIi0llEckUkOsj9JIqIikhMVcVWE4jIpyJyUbjjqGlqcrmISH0RmSwiO0Xkg2o43ksick+oj1PbiD1nEh4isgpoAxT5zX5dVa8NT0T7E5FEYCUQq6qF4Y3m4IjIfUB3VT0/3LHUFiKSArylqh2r8ZgXANcBx1b1e1FELgYuV9Xjq3K/dVGt+tYZgU5T1c/DHUQkEpGYSExyIiK4L3G+cMcSDgf5d+sCZJS3XaS+F2odVbWfMPwAq4Ch5SyLBp4EtgArgL8ACsSUtS1wH+7bIkBi8brAWGBuqX3fCEzyXp8C/ALsAtYC9/mtl1jqmE2B14ANwDrgISDaW3Yx8K0X83bcFc0ov321AP4FrPeWf+y37FTgV2AH8D1wWAVlpl5ZLANWevOe9WLfBcwDTvDmjwTygQIgF5jvzU/FfRMFd5v3bmA1kA28ATQt59hLgVP9pmO8v08/b/poL/4dwHwgxW/dVOBh4DtgL9DdK7MVQI5XXueV/luW83coc7sy4r0P+AB4y1t3IdADuMM717XA8FIxFpfLP4EP/ZY9BnwBNPTi93llmgu0B14HHvJbPwXIKvVevw1YAOR5Zdce+C+w2TuPv5ZzHveX+jte5pXBd8AzwDbce/EQ4Etgq/d3eRto5refTsD/vONtBZ4HegH7cHcHcoEd3rqlz+fPQKZ3rElA+1Lvyatw78ntwAt4d3zq2k/YA6irP1ScTK4C0rx/gBbAVxxcMmngfZAk+a07BxjrvU4BDsV9qB4GbALOKL0fb/pj4GXvA6U1MBu40lt2sffP/mdcIrwalziKb6NOBf4DNAdigZO8+f1wH2yDvO0u8s4trpxyUWCmVyb1vXnnAy29870Z2AjEly4Xv32kUvKhean3IdENaOR92LxZzrHHAW/7TZ8CpHmvO+A+oEZ7ZTnMm07wO+YaoI8XZ1Nc8kv2lrcD+pQVc6m/Z8Pytisj3vtwH5QjvG3fwH1o3+X9Df6Ml5DLKJcGQIb3dz0B9+Hc0e89k1XqWK9TeTL5Ffd+ru+V0TyvTOt55b8CGFHBufiXycVAIe7WV4y3z+5euccBCcAs4O/e+tG4BP+MV4bxwPF++/q2vPMBhnjn38/b93PArFLvySlAM6AzLlmNDPfnSzh+wh5AXf3x/sFycd9ki3/+7C37ErjKb93hHEQy8abfAsZ5r5NwyaVBOTH9HXim9H5wdTt5eB/g3vJzgK+81xcDmX7LGnjbtsV94PmA5mUc75/Ag6XmpeMlmzLWV2BIJeW6HTi8dLn4LU+l5EPzC+Aav2XJuKQYU8Z+u/uXHe6bb3G53kapJARMBy7yO+YDfssaen/vP/qXaVkx8/tkUuZ2ZcR7HzDTb/o07/1WfDXZ2Ntvs9Ll4k0PxH0TXw2c4zc/hYNLJpf6TQ8C1pTaxx3Avyo4l9LJZE1Z6/qtcwbwi/f6GNyHfFl/14upOJm8Bjzut6yR9x5J9HtPHu+3/H3g9opiq60/1porvM5Q1WZ+P69489vjbkMUWx3EMd7BffADnIu7xbQHQEQGichXIrJZRHbirohalbGPLrhvsxtEZIeI7MBdpbT2W2dj8Yvi/eP+8ToB21R1ezn7vbl4n95+O+HOvzz+5YKI3CwiS72WPjtw3/rLOoeytGf/sl1NSfLcj6pm4m51nSYiDYAxuLItPo//K3Uex+MS6e/iVtXdwNm48t4gIlNFpGdlwR7Edpv8Xu8Ftqhqkd80uL9RWceajbtaENwHZLD8/25dgPalyutOyij3APeHiLQWkfdEZJ2I7MJ9iSp+H3QCVuvB1avs9x5R1VzcVWcHv3U2+r3eQzllWttZMqmZNuD+AYp1LrV8N+7bf7G2FexrBtBKRI7AJZV3/Ja9g7sH3ElVmwIv4T48SluLuzJp5Zf4mqhqnwDOZS3QQkSalbPs4VIJtYGqvlvB/rT4hYicgLsq+BPuyqcZsNPvHPT3m+9nPe6DrVhn3O2TTWWvzru4MjwdWOIlmOLzeLPUeTRU1fFlxQ2gqtNVdRgu4aQBxV8kKvzbVrBdlRKRv+Bu66wH/uYfQhmrB/J+9N9uLe4Wm395NVbV0QcQYuk4HvXmHaaqTXC3P4vfB2uBzuU0cz+g94iINMTdVl13ALHWCZZMaqb3gb+KSEcRaQ7cXmr5r8BYEYkVkQHAWeXtyPs29iHwBK6uYabf4sa4q4Z9IjIQd+VS1j424JLSUyLSRESiROQQETmpshPxtv0UeFFEmnsxn+gtfgW4yrtCEhFpKCKniEjjyvbrF38h3i0MERkHNPFbvglIFJHy3ufvAjeKSFcRaQQ8Avyngm+w7+FuOV7N/kn5LdwVywgRiRaReBFJEZEym8+KSBsRGeN9MOXhbj8VXzH8CpzoPS/UFHf7J5DtqoyI9MBVap8PXAD8zfsyAq5MW3qxFfsVGC0iLUSkLXBDJYeYDewSkdu8Z0iiRaSviBwVRNiN8W4bi0gH4NZSx9sAjPfeY/Eicpzf+XQUkXrl7Pcd4BIROUJE4nDvkZ9UdVUQsdZKlkzCa7L3gGHxz0fe/Fdw99znAz/jKob93YNrvbId19rlHSr2DjAU+KDUB+U1wAMikoOrDK3odsaFuMrSJd5xP2T/2zgVuQB3nzkNV+F+A4CqzsVVBD/v7TMTdw87UNNxiSoDdytiH/vf/ih+wG2riPxcxvYTgTdxlbUrve2vK+9gXmL8ATgW16CgeP5a3NXKnbjEthb3YVbe/1cUrrHAely9xEm4vwWqOtPb9wJcJfWUQLarKt6397eAx1R1vqou887rTRGJU9U0XBJe4d2iao8rw/m4upEZ+JVNWbxbbacBR+DKfQvwKu4W5cG6H1dJvhPX4OO3/xm/43XHNYTIwt0uBFc/uRjYKCJbyoj1C9z/239xCekQXCtJU4o9tBgBasMDhMaY2s2uTIwxxgQtZMlERDp5LYWWishiEbm+jHVERP4hIpkiskBE+vktu0hElnk/F4UqTmOMMcEL2W0uEWkHtFPVn70K1Xm4prBL/NYZjbtHPRrX9vxZVR0kIi2AucAAXGuLeUD/cpqXGmOMCbOQXZmo6gZV/dl7nYNro9+h1GqnA2+o8yPQzEtCI3APXBU/nzAT1z2GMcaYGqhaOnr0KpCPBH4qtagD+7e+yfLmlTe/rH1fAVwBEB8f379z59KPZNRNPp+PqCirErNyKGFlUcLKokRGRsYWVU0Idj8hTyZe+/3/Ajeo6q7Si8vYRCuY//uZqhOACQDJycmanp4eRLS1R2pqKikpKeEOI+ysHEpYWZSwsighIsH0sPGbkKZmEYnFJZK3VbX0sxLgrjj8n/TuiGtDX958Y4wxNVAoW3MJrpO0par6dDmrTQIu9Fp1HQ3s9B4Mmw4M956Ybo576nh6qGI1xhgTnFDe5joO9+TzQhH51Zt3J14/U6r6EjAN15IrE9dB2iXesm0i8iCuu3RwPa5uC2GsxhhjghCyZKKq31J23Yf/Ooob7KisZRNx3V0YY4yp4aw5gzHGmKBZMjHGGBM0SybGGGOCZsnEGGNM0CyZGGOMCZolE2OMMUGzZGKMMSZolkyMMcYEzZKJMcaYoFkyMcYYEzRLJsYYY4JmycQYY0zQLJkYY4wJmiUTY4wxQbNkYowxJmiWTIwxxgTNkokxxpigWTIxxpg6andeYZXtK2TD9orIROBUIFtV+5ax/FbgPL84egEJ3vjvq4AcoAgoVNUBoYrTGGPqoi/TNnH3R4uqbH8hSybA68DzwBtlLVTVJ4AnAETkNOBGVd3mt8pgVd0SwviMMabOyc7Zx/2TlzB1wYYq3W/IbnOp6ixgW6UrOucA74YqFmOMqet8PuXd2WsY+tTXTF2wgQb1orn7lF5Vtv9QXpkEREQaACOBa/1mKzBDRBR4WVUnhCU4Y4ypBTKzc7nzo4XMXum+36ckJ/DQGX3p2LwBf66iY4iqVtGuyti5SCIwpaw6E791zgbOV9XT/Oa1V9X1ItIamAlc513plLX9FcAVAAkJCf3ff//9KjyDyJWbm0ujRo3CHUbYWTmUsLIoUVfKotCnTF1RwOTlBRQqNKkH5/WKY2DbaEQEgMGDB8+rinrpsF+ZAGMpdYtLVdd7v7NF5CNgIFBmMvGuWiYAJCcna0pKSkiDjRSpqalYWVg5+LOyKFEXymLuqm3c8b+FLMsuAODsAZ24Y3RPmjWoF5LjhTWZiEhT4CTgfL95DYEoVc3xXg8HHghTiMYYE1F27SvgsU/TePunNQB0a9WQh888lGMOaRnS44ayafC7QArQSkSygHuBWABVfclb7Uxghqru9tu0DfCRdwkWA7yjqp+FKk5jjKkNVJXpizcy7pPFZOfkERMlXJ1yCH8Z3J342OiQHz9kyURVzwlgnddxTYj9560ADg9NVMYYU/ts2LmXcZ8sZuaSTQD069yMR/9wGMltG1dbDDWhzsQYY8xBKPIpb/24miemp5ObV0ijuBhuG9WT8wZ2JipKqjUWSybGGBOB0jbu4o7/LeSXNTsAGNGnDfeP6UvbpvFhiceSiTHGRJB9BUU89+UyXv56BYU+pU2TOO4f05eRfduGNS5LJsYYEyG+z9zCnR8tZNXWPYjABUd34daRyTSJjw13aJZMjDGmptu+O59Hpi3lg3lZAPRo04hH/3Ao/bu0CHNkJSyZGGNMDaWqTJq/ngcmL2Hr7nzqRUdx3ZDuXHnSIdSLqVkjiFgyMcaYGmjttj3c/fEivs7YDMCgri145A+HckhCzewGxpKJMcbUIIVFPv713SqenpnB3oIimtaP5a7Rvfi/AR1/60+rJrJkYowxNcSidTu5/X8LWLRuFwCnHd6ecaf2JqFxXJgjq5wlE2OMCbM9+YU8MzOD175diU+hQ7P6PHRGXwb3bB3u0AJmycQYY8IoNT2buz5axLode4kSuOz4rtw0rAcN4yLr4zmyojXGmFpiS24eD05Zwie/rgegd7smjP/joRzWsVmYIzs4lkyMMaYa+XzK+3PX8uinaezcW0B8bBQ3Du3Bpcd3JTa6ZjX3PRCWTIwxppqkb8zhro8WMnf1dgBOSGrFw2ccSueWDcIcWfAsmRhjTIjtzS/i2S+W8eo3rj+tVo3iGHdab047rF2Nbu57ICyZGGNMCH2Ztolxnywma/teROD8oztz64ieNK0f/v60qpIlE2OMCYGNO/dx/+TFfLpoIwC92jXhkTP7cmTn5mGOLDQsmRhjTBUq8in//n4VT81IZ3d+EQ3qRXPTsB5cfGwiMRFcwV4ZSybGGFNFFmTt4M6PFv72BPvw3m24d0wfOjSrH+bIQi9kyUREJgKnAtmq2reM5SnAJ8BKb9b/VPUBb9lI4FkgGnhVVceHKk5jjAnWrn0FPDU9nTd+XI0qtG8az/2n92VY7zbhDq3ahPLK5HXgeeCNCtb5RlVP9Z8hItHAC8AwIAuYIyKTVHVJqAI1xpiDoapMW7iR+ycvJjsnj+go4bITunL9yUkR9wR7sEJ2tqo6S0QSD2LTgUCmqq4AEJH3gNMBSybGmBpjzdY9jJu0iNR010X8EZ2a8ciZh9K7fZMwRxYe4U6dx4jIfGA9cIuqLgY6AGv91skCBpW3AxG5ArgCICEhgdTU1NBFG0Fyc3OtLLBy8GdlUSKYsij0KZ+tLOCT5QUU+KB+DPxfj3qkdMonO+NnsjOqNtZIEc5k8jPQRVVzRWQ08DGQBJT1BI+WtxNVnQBMAEhOTtaUlJQQhBp5UlNTsbKwcvBnZVHiYMti9spt3PXRQpZlFwBw+hHtueuUXrRuHF/FEUaesCUTVd3l93qaiLwoIq1wVyKd/FbtiLtyMcaYsNi+O59HP13K+3PdGOyJLRvw4Bl9OSEpIcyR1RxhSyYi0hbYpKoqIgOBKGArsANIEpGuwDpgLHBuuOI0xtRdqsqH87J4ZNpStu8poF50FFelHMI1KYcQHxsd7vBqlFA2DX4XSAFaiUgWcC8QC6CqLwFnAVeLSCGwFxirqgoUisi1wHRc0+CJXl2KMcZUm8zsHO76aBE/rdwGwNHdWvDQGYfSvXXNHIM93ELZmuucSpY/j2s6XNayacC0UMRljDEV2VdQxPNfZvLyrOUUFCktGtbj7lN6ceaRHWpNp4yhUGEyEZF43IOHJwDtcVcQi4CpdrVgjKltZmVs5p5PFrF66x4Axh7VidtH9aRZg3phjqzmKzeZiMh9wGlAKvATkA3EAz2A8V6iuVlVF4Q+TGOMCZ3sXft4YMoSpizYAEBym8Y8fGZfBiS2CHNkkaOiK5M5qnpfOcueFpHWQOeqD8kYY6pHkU9556fVPP5ZOjl5hcTHRnH9yT24/ITIHvUwHMpNJqo6tfQ8EWmrqhu95dm4qxVjjIk4i9bt5K6PFzF/7Q4ABicn8MDpfenUIvJHPQyHA62Anwb0C0UgxhhTHXLzCnl3aR4zp3+LT6FNkzjuO60PI/u2tQr2IBxoMrGSNsZEpOJOGR+csoSNuwqJErjkuERuGtaDxvG1a9TDcDjQZPJKSKIwxpgQWrE5l3snLeabZVsA6Nokin9ceCyHdmwa5shqj4paczVS1Vz/ear6YmXrGGNMTbE3v4gXUzN5+esV5Bf5aFo/lltHJNN+7wpLJFWsoiuTT0TkV9wAVvNUdTeAiHQDBgN/wl2pfBjyKI0x5gB9vmQT901eTNb2vQD8X/+O3D6qJy0bxZGaurKSrc2Bqqg118leb75XAseJSHOgEEgHpgIXFbfsMsaYmmLttj3cP3kJny/dBEDPto156Ax7ZiTUKqwzsW5NjDGRIq+wiFe/WclzXy5jX4GPRnEx3DisBxcd04UYe2Yk5MI9OJYxxgTt22VbGPfJIlZs2Q3AaYe35+5TetGmiY0zUl0smRhjItbGnft4aGpJNyjdEhry4Ol9Oa57qzBHVvdYMjHGRJyCIh///n4Vz8zMYHd+EfGxUVw3JInLT+hKXIyNMxIOASUTETkc13MwwDeqOj90IRljTPnmrNrGPR8vIm1jDgDDerdh3Km9rRuUMKs0mYjI9cCfgf95s94SkQmq+lxIIzPGGD9bcvN4dFoa//3ZDZ3bqUV97h/ThyE924Q5MgOBXZlcBgzye87kMeAHwJKJMSbkinzKO7PX8MRnaezaV+iGzj2pG9cM7m5D59YggSQTAYr8pouwPrqMMdVg/tod3PPJIhZk7QTgxB4J3D+mD11bNQxzZKa0QJLJv4CfROQjb/oM4LXKNhKRibhRGrNVtW8Zy88DbvMmc4Gri+tiRGQVkINLXIWqOiCAOI0xtcTOPQU8Pj2Nd2avQRXaNoln3Gm9GWU9+9ZYlSYTVX1aRFKB43FXJJeo6i8B7Pt13Bjvb5SzfCVwkqpuF5FRwARgkN/ywaq6JYDjGGNqCZ9P+e/PWYz/NI2tu/OJiRIuO6Erfz05iYZx1vi0Jquoo8cmqrpLRFoAq7yf4mUtVHVbRTtW1VkikljB8u/9Jn8EOgYWsjGmNlq6YRf3fLyIuau3AzCwawseOqMvPdo0DnNkJhCiqmUvEJmiqqeKyErAfyUBVFW7Vbpzl0ymlHWbq9R6twA9VfVyb3olsN077suqOqGCba8ArgBISEjo//7771cWVp2Qm5tLo0aNwh1G2Fk5lKipZbG3UPl4WT4z1xTiU2hSTxjbsx7HtIsO2S2tmloW4TB48OB5VVGVUG4yqQqBJBMRGQy8CByvqlu9ee1Vdb03zvxM4DpVnVXZ8ZKTkzU9Pb1KYo90qamppKSkhDuMsLNyKFHTykJVmbxgAw9NWUJ2Th5RAhcc3YWbhifTtH5oB6uqaWURTiJSJckkkOdMvlDVkyubdzBE5DDgVWBUcSIBUNX13u9sr+J/IFBpMjHGRIblm3MZ98kivst0//ZHdGrGQ2f0pW8HG2MkUlVUZxIPNABaed3PF19vNgHaB3tgEemMexDyAlXN8JvfEIhS1Rzv9XDggWCPZ4wJv735RTz35TJe+WYFBUVKswax3DayJ2cP6ERUlLXSimQVXZlcCdyASxzzKEkmu4AXKtuxiLwLpOCSURZwLxALoKovAeOAlsCL3n3R4ibAbYCPvHkxwDuq+tmBnpgxpuZQVaYv3siDU5aybocbrGrsUZ3428ietGhYL8zRmapQ0eBYzwLPish1B9N1iqqeU8nyy4HLy5i/Ajj8QI9njKmZMrNzuX9yyfjrvds14cEz+tK/S/MwR2aqUiDPmTwnIn2B3kC83/zynh8xxhhy8wp57otlvPbtSgp9StP6sdwyIplzB3Ym2m5p1TqBVMDfi7td1Rs36uIo4FvKfxjRGFOHqSqT5q/n4alLyc7JQwTOGdiZW0ck2y2tWiyQR0rPwt12+kVVLxGRNrgWWMYYs5+lG3Zx76TFzF7pnmk+olMzHji9D4d1bBbmyEyoBZJM9qqqT0QKRaQJkA1U+sCiMabu2Lm3gGdmZvDmj6sp8iktG9bjtpE9Oat/R2ulVUcEkkzmikgz4BVcq65cYHZIozLGRASfT/lwXhaPfeb60ooSuPjYRG4c2oOmDUL74KGpWSpMJuLa5z6qqjuAl0TkM6CJqi6oluiMMTXWgqwd3PPJYuav3QG4vrTuH9OHXu2ahDkyEw4VJhNVVRH5GOjvTa+qjqCMMTXXtt35PDE9jffmrEUV2jSJ487RvRhzeHvrHr4OC+Q2148icpSqzgl5NMaYGqvIp7zz02qenJHBzr0Frnv4E7ty3ZAkGln38HVeIO+AwcCVIrIa2E1Jr8GHhTQyY0yNMXfVNsZ9spglG3YBcEJSK+49rQ/dW1vPu8YJJJmMCnkUxpgaKTtnH+M/TeN/P68DoEOz+txzam9G9Gljt7TMfgJ5An51dQRijKk5Cop8/Pv7Vfz982Xk5hVSLyaKq07sxtUp3alfLzrc4ZkayG50GmP2813mFu6btJhl2bkADO3VhnGn9qZzywZhjszUZJZMjDEArN+xl4enLmXqwg0AJLZswL2n9WFwz9ZhjsxEgoCSiYh0AZJU9XMRqQ/EqGpOaEMzxlSHvMIiXv1mJc9/mcnegiLqx0Zz7ZDuXH5CV+Ji7JaWCUwgHT3+GTfGegvgEKAj8BIQ9EiLxpjw+iotm/snL2bV1j0AnHJYO+4a3Yv2zeqHOTITaQK5MvkLbtjcnwBUdZk3NrsxJkKt2bqHB6Ys5vOl2QAktW7E/WP6cGz3VmGOzESqQJJJnqrmFzcDFJEYQEMalTEmJPbmF/HPr5fz0tfLyS/00SguhhuGJnHRsYnERkeFOzwTwQJJJl+LyJ1AfREZBlwDTA5tWMaYqqSqzNtUyF1Pf/3bsLl/6NeB20f1pHXj+Eq2NqZygXwVuR3YDCzEjQs/Dbg7kJ2LyEQRyRaRReUsFxH5h4hkisgCEennt+wiEVnm/VwUyPGMMb+XsSmHC16bzXO/5LFux156t2vCh1cdw9N/OsISiakygVyZ1AcmquorACIS7c3bE8C2rwPPU/6ojKOAJO9nEPBPYJCItADuBQbgbqnNE5FJqro9gGMaY4Cdewp45vOSMUYaxsLto/tw7qAuNmyuqXKBJJMvgKG4cUzAJZIZwLGVbaiqs0QksYJVTgfeUFXFdSjZTETa4YYJnqmq2wBEZCYwEng3gHiNqdOKfMq7s9fw1Ix0tu8pIErggqO7MKjBZk49JjHc4ZlaKpBkEq+qxYkEVc0Vkap6FLYDsNZvOsubV9783xGRK3BNl0lISCA1NbWKQotsubm5VhbUvXJI21bE20vzWZvjA6BniyjO6xVHp8ZbyM3dXafKoiJ17X1RHQJJJrtFpJ+q/gwgIv2BvVV0/LKutbWC+b+fqToBmACQnJysKSkpVRRaZEtNTcXKou6UQ9b2PTw6Le23p9c7NKvP3af0YmTftr91yFhXyiIQVhZVL5BkcgPwgYis96bbAWdX0fGzgE5+0x2B9d78lFLzU6vomMbUGnvzi3jJa+qbV+gjPjaKa1K6c8WJ3YiPtafXTfUJpNfgOSLSE0jGXTGkqWpBFR1/EnCtiLyHq4DfqaobRGQ68IiINPfWGw7cUUXHNCbiqSpTF27gkalLWb9zHwBjDm/P7aN62tPrJiwC7ejxKCDRW/9IEUFVy2uh9RsReRd3hdFKRLJwLbRiAVT1JVwz49FAJq512CXesm0i8iBQPLrjA8WV8cbUdYvX7+T+yUuYvdL9S/Rp34T7xvThqMQWYY7M1GWB9M31Jq5Prl+BIm+2Un5z39+o6jmVLFdcdy1lLZsITKzsGMbUFdt25/PkjHTem70Gn0KLhvW4dUQyfxrQyZr6mrAL5MpkANDb++A3xlSzgiIfb/24mmdmZrBrXyExUcLFxyZy/dAkmtaPDXd4xgCBJZNFQFtgQ4hjMcaU8s2yzTwweclvA1Wd2COBcaf2onvrxmGOzJj9BZJMWgFLRGQ2kFc8U1XHhCwqY+q41Vt389DUpcxcsgmALi0bcM8pvTm5V2sbe93USIEkk/tCHYQxxtmdV8gLX2Xy6jcryS/y0bBeNNcOSeLS4xNtoCpTowXSNPjrUiMtNgDsXW1MFfL5lI9/Xcf4T9PIznE3AP7YryO3jUymdRPrjNHUfAcz0mIHbKRFY6rM/LU7uG/yYn5ZswOAwzs1477TenNk5+aVbGlMzWEjLRoTJtk5+3jis3Q+mJcFQELjOG4f2ZMzj+xAlDX1NRHGRlo0pprlF/r413cree7LTHLzCqkXHcWlx3fl2iHdaRQX6HPExtQsNtKiMdXoy7RNPDhlKSu37AZgaK/W3D8hiQcAABtwSURBVH1KbxJbNQxzZMYEJ5BkcjtwGfuPtPhqKIMyprZZvjmXB6csITV9MwCHJDRk3Gl9OKlHQpgjM6ZqBNKaywe84v0YYw7Azr0FPPfFMl7/fhWFPqVxfAw3DO3Bhcd0ITY6kFGzjYkM5SYTEVlIBXUjqnpYSCIyphYoLPLx7py1PDMzg2278xGBcwZ24ubhybRqFBfu8IypchVdmZxabVEYU4vMytjMQ1OXkLHJdYEyqGsL7jm1N307NA1zZMaETrnJRFVXV2cgxkS6zOxcHpm2lC/TsgHo1KI+d43uxYg+ba0LFFPrBfLQ4tHAc0AvoB7u6ffdqtokxLEZExF27Mnn758v460fV1PoUxrFxXDdkO5cfJx1gWLqjkBacz0PjAU+wHVHfyHQPZRBGRMJCop8vP3jap75fBk79xYQJXDOwM7cNKwHCY2tXsTULQE9IaWqmSISrapFwL9E5PsQx2VMjfZVejYPTVnC8s3ueZHjurfk7lN606udXbCbuimQZLJHROoBv4rI47hxTewJK1MnZWzK4aGpS5mV4Z4XSWzZgLtO6c1Q6xre1HGBJJMLgCjgWuBGoBPwx0B2LiIjgWdx9Syvqur4UsufAQZ7kw2A1qrazFtWhHtQEmCNjZ9iwmnb7nyemZnBO7PXUOQ9L3L9yUlceEwi9WLseRFjAkkm/YBpqroLuD/QHYtINPACMAzIAuaIyCRVXVK8jqre6Lf+dcCRfrvYq6pHBHo8Y0Ihv9DHGz+s4tkvlpGzr5DoKOGCo7tw47AetGhYL9zhGVNjBJJMxgB/F5FZwHvAdFUtDGC7gUCmqq4AEJH3gNOBJeWsfw5wbwD7NSbkVJUvlmbz8LSSfrROSGrFPaf2pkcbGzLXmNJEtfIOgEUkFhgFnA0cD8xU1csr2eYsYGTxeiJyATBIVa8tY90uwI9AR6+SHxEpBH4FCoHxqvpxOce5AjfeCgkJCf3ff//9Ss+nLsjNzaVRo0bhDiPsDqYc1ub4eDctjyVbfQC0bSic07Meh7WKjuh6EXtPlLCyKDF48OB5qjog2P0E2pqrQEQ+xXWvUh93hVFhMgHK+q8rL3ONBT4sTiSezqq6XkS6AV+KyEJVXV5GbBOACQDJycmakpJSSVh1Q2pqKlYWB1YOW3LzeGpGBv+ZswafQtP6sdw4NInzjq4d/WjZe6KElUXVC+ShxZG4D/vBQCqux+A/BbDvLFxlfbGOwPpy1h2LG4TrN6q63vu9QkRScfUpv0smxgQrr7CI179bxfNfZpKTV0hMlHDhMV24YWgSzRpYvYgxgQjkyuRiXF3JlaqadwD7ngMkiUhXYB0uYZxbeiURSQaaAz/4zWsO7FHVPBFpBRwHPH4AxzamUqrK9MUbeWRaGmu27QFgSM/W3Dm6F91b2y0QYw5EIF3Qjz2YHatqoYhcC0zHNQ2eqKqLReQBYK6qTvJWPQd4T/evvOkFvCwiPlyz5PH+rcCMCdaidTt5cMoSflq5DYCk1o24+9TeNr6IMQepoi7ov1XV40UkB1fXIf6/A+mbS1Wn4QbT8p83rtT0fWVs9z1waCAnYMyByM7Zx5PT3bjrqtC8QSw3DU/mnKM6EVML6kWMCZeKeg0+3vtt7SBNxNtXUMRr367kxa8y2Z1fRGy0cNExiVx3chJN68eGOzxjIl4gFfCdy5qvqmuqPhxjqpaqMmXBeh6dlsa6HXsBGNa7DXeO7kVXG3fdmCoTSAX8VL/X8UBXIB3oE5KIjKki89fu4NHZ+8jY/gsAPds25p5Te3Nc91ZhjsyY2ieQCvj96i5EpB9wZcgiMiZIWdv38MT0dD751bVEb9mwHjcPT+bsozoRHRW5Dx0aU5MF9NCiP1X9WUSOCkUwxgRj174CXvxqORO/W0l+oY96MVEM7RTN+ItSaBJv9SLGhFIgdSY3+U1G4Tp+3ByyiIw5QAVFPt6dvYa/f76MbbvzATj9iPbcOiKZzPmzLZEYUw0CuTLxb81ViKtD+W9owjEmcMWdMT7y6VJWeINUHZXYnLtO6c0RnZoBkBnOAI2pQwKpMwm423ljqsuidTt5eOpSflixFXCDVN0+qhcj+rSJ6M4YjYlUgdzmmlTRchu0ylSn9Tv28uT0dP73yzoAmjWI5fqTkzhvUBcbpMqYMArkNtdKoC3wljd9DrAK102KMdUiN6+Ql1KX88o3K8gr9FEvOoqLju3CtYOTaNrA6kSMCbdAksmRqnqi3/RkEZmlqneGKihjihUW+fjP3LU8MzODLbmucv2Uw9px24iedG7ZIMzRGWOKBZJMEkSkm9+IiV0B6w3PhJSqkpqxmUemLmVZdi4A/To3465TetO/S/MwR2eMKS2QZHIjkCoiK7zpROyhRRNCS9bv4pFpS/k2cwsAnVrU5/aRvRh9aFurXDemhgqkNddnIpIE9PRmpR3guCbGBGTTrn08NaOkR98m8TFcNySJC4/tQlxMdLjDM8ZUoKIu6P+mqsUDUo1R1Q/8lj1idSamquzJL+Tlr1cwYdYK9hYUERMlXHBsF/46JInmDW2kQ2MiQUVXJmMpGd3wDuADv2UjAUsmJihFPuW/87J4ckY62TnuYndkn7bcNqqn9ehrTISpKJlIOa/LmjbmgHyzbDMPT11K2sYcAA7v2JS7TunNwK4twhyZMeZgVJRMtJzXZU0bE5CMTTk8PHUpX2e47t06NKvP30Ymc9ph7YmyHn2NiVgVJZPDRWQX7iqkvvcabzo+kJ2LyEjgWdwY8K+q6vhSyy8GngDWebOeV9VXvWUXAXd78x9S1X8HckxTM23OyePpmRn8Z84afAqN42K4ZnB3LjkukfhYq1w3JtJVNGxvUP/hIhINvAAMA7KAOSIySVWXlFr1P6p6baltWwD3AgNwV0HzvG23BxOTqX5784t49ZsVvPT1cnbnFxEdJVx4dGeuPzmJlo3iwh2eMaaKHPB4JgdgIJDp97Dje8DpQOlkUpYRwExV3eZtOxNX6f9uiGI1VcznUz76ZR1Pzkhnw859AAzt1ZrbR/Wie+tGYY7OGFPVQplMOgBr/aazgEFlrPdHETkRyABuVNW15WzboayDiMgVwBUACQkJpKamBh95LZCbmxu2sli8pYj3M/JZvcsHQJcmUYxNrkevlrvJWjKXrEC+TlSRcJZDTWNlUcLKouqFMpmUVZtauuJ+MvCuquaJyFXAv4EhAW7rZqpOACYAJCcna0pKykEHXJukpqZS3WWxZP0uxn+Wxiyvcr1tk3huHZHMmUd2CFvlejjKoaaysihhZVH1QplMsoBOftMdgfX+K6jqVr/JV4DH/LZNKbVtapVHaKrEuh17eWpGOh/9sg71KtevHnwIlxzblfr1rHLdmLoglMlkDpDkdQy5DvcQ5Ln+K4hIO1Xd4E2OAZZ6r6cDj4hIcY9+w3EPTpoaZOfeAl5MzeRf360iv9BHbLRw/tFduG5IEi3syXVj6pSQJRNVLRSRa3GJIRqYqKqLReQBYK6qTgL+KiJjcMMBbwMu9rbdJiIP4hISwAPFlfEm/PIKi3jzh9U8/1UmO/YUAHDa4e25dXiydQtvTB0VyisTVHUaMK3UvHF+r++gnCsOVZ0ITAxlfObA+HzK5AXreWJ6Olnb9wIwqGsL7hzdi8O9MdeNMXVTSJOJqT2+z9zCI58uZdE69+xqjzaNuH1UTwYnt7Zu4Y0xlkxMxdI27mL8p2mkprsWWm2axHHzsGT+2L8j0db9iTHGY8nElGn9jr08PTOD//7sxhZpFBfD1SmHcOlx1kLLGPN7lkzMfnbtK+CfqcuZ+O1K8gp9xEQJ5x/TheuGdLfuT4wx5bJkYgDXQuvtH9fw3JfL2O610DrlsHbcOjyZRBtbxBhTCUsmdZzPp0xZuIEnpqexdptroTWwawvuGNWTIzs3r2RrY4xxLJnUYd8v38L4T9NYkLUTgO6tG3H7yJ6c3MtaaBljDowlkzoofWMO4z9dyldeC63WjeO4aVgPzurfkZjoqDBHZ4yJRJZM6pANO/fyzMwMPpyXhU+hYb1orjrpEC47oSsN6tlbwRhz8OwTpA7Yta+Al79ezmvfrmRfgddC6+jO/PXkJFpZCy1jTBWwZFKL5Rf6ePun1fzji5IWWqMPbcutI3rS1VpoGWOqkCWTWkhVmb2hkHFPf82abXsAOCqxOXeM7kU/a6FljAkBSya1zA/LtzL+szTmr80DoFtCQ24f2ZNhvdtYCy1jTMhYMqklFq3byePT038b5bBpnHDb6L78aYC10DLGhJ4lkwi3eutunpqRwaT5bhDLRnExXHFiN3poFiMHdQ5zdMaYusKSSYTanJPHc18u452f1lDoU+pFR3H+0V34y+BDaNkojtTUdeEO0RhTh1gyiTA5+wp4ZdYKXv12JXvyixCBP/TrwE3DetCxuY1yaIwJD0smESKvsIi3flzDC19lsm13PgBDe7XmlhHJ9GzbJMzRGWPqupAmExEZCTyLGwP+VVUdX2r5TcDluDHgNwOXqupqb1kRsNBbdY2qjgllrDVVkU/5+Jd1PD0zg3U7XEeM/bs05/ZRPTkqsUWYozPGGCdkyUREooEXgGFAFjBHRCap6hK/1X4BBqjqHhG5GngcONtbtldVjwhVfDWdqvJlWjaPf5ZO+qYcwA2V+7cR1hGjMabmCeWVyUAgU1VXAIjIe8DpwG/JRFW/8lv/R+D8EMYTMeat3sb4T9OYs2o7AB2a1efGYT0488gONlSuMaZGCmUy6QCs9ZvOAgZVsP5lwKd+0/EiMhd3C2y8qn5c9SHWLBmbcnj8s3Q+X7oJgOYNYvnL4O6cf3QX4mNtqFxjTM0VymRS1ldoLXNFkfOBAcBJfrM7q+p6EekGfCkiC1V1eRnbXgFcAZCQkEBqamrQgVe3rXt9fJRZwHfrClGgXjSMSIxlVGIsDYrW8ON3aw54n7m5uRFZFlXNyqGElUUJK4uqF8pkkgV08pvuCKwvvZKIDAXuAk5S1bzi+aq63vu9QkRSgSOB3yUTVZ0ATABITk7WlJSUqjuDENu+O58XvsrkjR9Xk++Nt37OwM5cd3J3WjeOD2rfqampRFJZhIqVQwkrixJWFlUvlMlkDpAkIl2BdcBY4Fz/FUTkSOBlYKSqZvvNbw7sUdU8EWkFHIernK8V9uQXMvHblbz89Qpy8goBOO3w9tw8rIeNt26MiUghSyaqWigi1wLTcU2DJ6rqYhF5AJirqpOAJ4BGwAde66TiJsC9gJdFxAdE4epMlpR5oAhSUOTjvTlr+ccXy9ic4y7CTkhqxW0je9K3Q9MwR2eMMQcvpM+ZqOo0YFqpeeP8Xg8tZ7vvgUNDGVt18vmUqQs38NSMdFZtdV3CH96xKbeN7Mmx3VuFOTpjjAmePQEfYt8s28xjn6WxaN0uALq1asitI5IZ2betPStijKk1LJmEyIKsHTz2WRrfZW4FoE2TOG4Y2oP/629dwhtjah9LJlVsxeZcnpqRwdSFGwBoEh/D1SndufjYROrXs2dFjDG1kyWTKrJx5z7+8eUy/jNnLUU+JS4miouPS+Sak7rTtEFsuMMzxpiQsmQSpG278/lnaiZv/LCavEIfUQJjj+rE9UOTaNe0frjDM8aYamHJ5CDl7Cvg1W9W8tq3K8n1nhU55dB23DisB91bNwpzdMYYU70smRygfQVFvPnDal5MzWT7ngIAUpITuGV4sj0rYoypsyyZBKigyMcHc7P4xxfL2LhrHwBHJTbn1hE9GdjVxhUxxtRtlkwq4fMpkxes5+mZGaz2Hjjs3a4Jt45MJqVHgj0rYowxWDIpl6ry+dJsnpqRTtpGNzhVt4SG3DwsmVF92xJl44oYY8xvLJmU4fvlW3hiejq/rNkBQPum8dwwtAd/6NfBHjg0xpgyWDLx8+vaHTw5PZ1vM7cA0LJhPa4d0p1zB3UmLsYeODTGmPJYMsGNcPjUjHSmL3YjHDaOj+HKE7txyXFdaRhnRWSMMZWp05+Ua7bu4e+fZ/DRr+tQhfjYKC45ritXntiNZg3qhTs8Y4yJGHUymWTvcl2fvDd7LYU+JTbajXB47eDutG4S3AiHxhhTF9WpZLJ9dz4vzVrOv79fxb4C1/XJH/t15IahSXRq0SDc4RljTMSqE8kkN6+Qf327kgmzSobJHdmnLTcP70FSm8Zhjs4YYyJfrU4m+wqKePunNbz4VSZbd+cDbpjcW4Ync3inZmGOzhhjao9amUwKi3x8OC+LZ79YxoadruuTfp2bceuInhxzSMswR2eMMbVPSJOJiIwEngWigVdVdXyp5XHAG0B/YCtwtqqu8pbdAVwGFAF/VdXpgRxz8nzX9cnKLbsB6Nm2MbeOSGZIz9bW9YkxxoRIyJKJiEQDLwDDgCxgjohMUtUlfqtdBmxX1e4iMhZ4DDhbRHoDY4E+QHvgcxHpoapFFR1zfa6P6979BYDElg24aXgypx7azro+McaYEAvllclAIFNVVwCIyHvA6YB/MjkduM97/SHwvLjLh9OB91Q1D1gpIpne/n6o6ID5PmjbJJ7rhyZxVv+OxFrXJ8YYUy1CmUw6AGv9prOAQeWto6qFIrITaOnN/7HUth3KOoiIXAFc4U3m/XTX0EXn3gXnBh9/pGsFbAl3EDWAlUMJK4sSVhYlkqtiJ6FMJmXdW9IA1wlkWzdTdQIwAUBE5qrqgAMJsraysnCsHEpYWZSwsighInOrYj+hvA+UBXTym+4IrC9vHRGJAZoC2wLc1hhjTA0RymQyB0gSka4iUg9XoT6p1DqTgIu812cBX6qqevPHikiciHQFkoDZIYzVGGNMEEJ2m8urA7kWmI5rGjxRVReLyAPAXFWdBLwGvOlVsG/DJRy89d7HVdYXAn+prCWXZ0IoziVCWVk4Vg4lrCxKWFmUqJKyEHchYIwxxhw8aztrjDEmaJZMjDHGBC0ikomIjBSRdBHJFJHby1geJyL/8Zb/JCKJfsvu8Oani8iI6ow7FA62LERkmIjME5GF3u8h1R17VQvmfeEt7ywiuSJyS3XFHCpB/o8cJiI/iMhi7/0R0YP6BPE/Eisi//bKYKnXpVNEC6AsThSRn0WkUETOKrXsIhFZ5v1cVHrb31HVGv2Dq7xfDnQD6gHzgd6l1rkGeMl7PRb4j/e6t7d+HNDV2090uM8pTGVxJNDee90XWBfu8wlXWfgt/y/wAXBLuM8njO+LGGABcLg33bIO/4+ci+t5A6ABsApIDPc5hbgsEoHDcH0knuU3vwWwwvvd3HvdvKLjRcKVyW/dsqhqPlDcLYu/04F/e68/BE4u3S2Lqq4EirtliVQHXRaq+ouqFj+rsxiI9zrajFTBvC8QkTNw/yCLqyneUAqmLIYDC1R1PoCqbtXAWk7WVMGUhQINvWfe6gP5wK7qCTskKi0LVV2lqgsAX6ltRwAzVXWbqm4HZgIjKzpYJCSTsrplKd21yn7dsgD+3bJUtm0kCaYs/P0R+EVd32eR6qDLQkQaArcB91dDnNUhmPdFD0BFZLp3u+Nv1RBvKAVTFh8Cu4ENwBrgSVXdFuqAQyiYz78D3jYSxjOplm5ZIkQwZeEWivTB9c48vArjCodgyuJ+4BlVza0lwxIEUxYxwPHAUcAe4AsRmaeqX1RtiNUmmLIYiBvyoj3u1s43IvK5ep3VRqBgPv8OeNtIuDKxbllKBFMWiEhH4CPgQlVdHvJoQyuYshgEPC4iq4AbgDu9B2wjVbD/I1+r6hZV3QNMA/qFPOLQCaYszgU+U9UCVc0GvgMiuf+uYD7/DnjbSEgm1i1LiYMuCxFpBkwF7lDV76ot4tA56LJQ1RNUNVFVE4G/A4+o6vPVFXgIBPM/Mh04TEQaeB+sJ7H/MBGRJpiyWAMMEachcDSQVk1xh0IgZVGe6cBwEWkuIs1xdzIqHqAw3C0OAmyVMBrIwLVMuMub9wAwxnsdj2uVk4lLFt38tr3L2y4dGBXucwlXWQB34+4H/+r30zrc5xOu94XfPu4jwltzBVsWwPm4hgiLgMfDfS7hKgugkTd/MS6h3hruc6mGsjgKdxWyGzfa7WK/bS/1yigTuKSyY1l3KsYYY4IWCbe5jDHG1HCWTIwxxgTNkokxxpigWTIxxhgTNEsmxhhjgmbJxEQUEVERedNvOkZENovIFG96TFm9o5bax8Ui0j6EMf7V63X27SrY136xisirItI72P16+zpDRMaVs+xz7/kCYwJiycREmt1AXxGp700PA9YVL1TVSao6vpJ9XIzrMiNg3gN9gboGGK2q5wWxj2IX4xerql6uqlX1UOHfgBfLWfYm7jyMCYglExOJPgVO8V6fA7xbvMD7Jv+89/oTEbnQe32liLztjdkwAHhbRH4VkfoiskpEWnnrDRCRVO/1fSIyQURmAG+ISLSIPCEic0RkgYhcWTowEXkJ1+X3JBG5sYx9JIrIN16nij+LyLF+2/5N3Fga80VkfDmxporIAG/9c7z1F4nIY377yRWRh739/CgibcqIsweQp6pbyinjSV7ZGhMQSyYmEr2H6yYnHjcWw0/lrHcFME5ETgBuBq5T1Q+BucB5qnqEqu6t5Fj9gdNV9VzgMmCnqh6Fe3L4z143Pb9R1atwfRgNVtVnythHNjBMVfsBZwP/ABCRUcAZwCBVPRz3JHq5sXq3vh4DhgBHAEeJ61YfoCHwo7efWcCfyziv44Cfyztpdd2Ox4lI6R6njSlTJPQabMx+VHWBuNHxzsF1TFjeepu8OoGvgDP14LoTn+T3IT4c149V8Yh0TXH9va08gH3EAs+LyBG4Hmp7ePOHAv9S19kiAcR6FJCqqpsBvPqZE4GPceNwTPHWm4e7FVhaO2BzJcfIxt1i21rJesZYMjERaxLwJJDC78dr8Xco7sOwojqSQkqu0ksPWbvb77Xgrm4q7vDu9/z3cSOwCTjcO+Y+v30fSN9GFfWdX6Al/SQVUfb/+V5cMkREonFJB1ziK66Uj/fWM6ZSdpvLRKqJwAOqurC8FURkIDAKN2TxLX63pHKAxn6rrsLdigI3cFh5pgNXi0ist/8eXu+yB6IpsEFVfcAFuKFVAWYAl4pIA2/fLcqJtdhPwEki0spLBucAXx9AHEuB7gCqWuTdRjuiOJGIiABtcWVjTKUsmZiIpKpZqvpsecvFDUn8CnCpuuGKbwYmeh+SrwMvFVdq4wbLelZEvsF9ky/Pq7jeZH8WkUXAyxz41f2LwEUi8iPuFtdu73w+w11tzRWRX4FbvPVLx1p8/huAO3C38OYDP6vqJwcQxyzgSK88ytIfV+9SeAD7NHWY9RpsTB0lIs8Ck1X183KWTdLIHXHRVDO7MjGm7noEaFDOskWWSMyBsCsTY4wxQbMrE2OMMUGzZGKMMSZolkyMMcYEzZKJMcaYoFkyMcYYE7T/B/X5bqDlRN3qAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot equivalence ratio as a function of mixture fraction\n",
    "fig, ax = plt.subplots()\n",
    "ax.plot(farr, phiarr, label='Phi',   lw=2)\n",
    "ax.set(xlim=(0.,0.1),\n",
    "       ylim=(0.,2.0),\n",
    "       xlabel='Mixture fraction (-)',\n",
    "       ylabel='Equivalence ratio (-)',\n",
    "       title ='Equivalence ratio versus mixture fraction')\n",
    "#ax.legend(handlelength=4,loc='best')\n",
    "ax.grid(True)\n",
    "fig.savefig(\"Phi_vs_f_zoomed.pdf\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Over this range, the relationship between $\\Phi$ and $f$ is close to (but not exactly) linear. Note that $\\Phi=1$ where $f=f_{stoich}$. As we have seen, the value of $ f_{stoich} $ varies depending on the fuel and oxidizer that are considered, but for most hydrocarbon-fuel/air-oxidizer mixtures, the value of $ f_{stoich} $ remains in a relatively narrow range around 0.06.\n",
    "\n",
    "The equivalence ratio $ \\Phi $ normally is defined for a particular unburned fuel-oxidizer reactant mixture. From this development, it can be seen that the value of $\\Phi$ (or of $f$) does not change as a result of chemical reaction: its value in the equilibrium products of combustion is equal to its value in the unburned reactants, and to its value at any intermediate partially burned state. The local value of $f$ tells us something about the state of the mixing between fuel an oxidizer, then, rather than something about the degree of completeness of chemical reactions. \n",
    "\n",
    "The mixture fraction provides no new information in a premixed system, compared to quantities that we have introduced earlier. Mixture fractions are useful principally in nonpremixed systems, as we will see when we get to Chapters 9 (Laminar Diffusion Flames) and 13 (Turbulent Nonpremixed Flames)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}