HW3-blank.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Homework 3\n",
    "Due Friday February 5. Hopefully short. Based on question 5 of the exam.\n",
    "Rename it like `HW3-west.r` so the filename is `HW3-west.r.ipynb` where `west.r` is *your* husky ID (the bit before @husky in your email address). Then upload the `ipynb` file to BlackBoard. Help available via Slack.\n",
    "\n",
    "First set up the kernel and import some libraries."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "## The following allows you to do plt.savefig('levenspiel.pdf') before you plt.show()\n",
    "## and makes it draw nice vector SVG figures instead of PNG figures. Try it out.\n",
    "#from IPython.display import set_matplotlib_formats\n",
    "#set_matplotlib_formats('pdf', 'svg')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First we define volumetric flow rate, initial feed concentration, and a rate expression that gives $r_A$ as a function of concentration $C_A$.  For now we can just use it as a \"black box\" function. We give it $C_A$ and it gives us $r_A$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "v = 1.0 # volumetric flow rate (m3/s)\n",
    "Ca0 = 1.0 # concentration (mol/m3)\n",
    "\n",
    "# Rates are reported at a range of concentrations.\n",
    "# Then we'll make a function to intepolate.\n",
    "concentrations = np.array([ 1.00     ,  0.9890101 ,  0.9790202 ,  0.9690303 ,  0.9590404 ,\n",
    "        0.94905051,  0.93906061,  0.92907071,  0.91908081,  0.90909091,\n",
    "        0.89910101,  0.88911111,  0.87912121,  0.86913131,  0.85914141,\n",
    "        0.84915152,  0.83916162,  0.82917172,  0.81918182,  0.80919192,\n",
    "        0.79920202,  0.78921212,  0.77922222,  0.76923232,  0.75924242,\n",
    "        0.74925253,  0.73926263,  0.72927273,  0.71928283,  0.70929293,\n",
    "        0.69930303,  0.68931313,  0.67932323,  0.66933333,  0.65934343,\n",
    "        0.64935354,  0.63936364,  0.62937374,  0.61938384,  0.60939394,\n",
    "        0.59940404,  0.58941414,  0.57942424,  0.56943434,  0.55944444,\n",
    "        0.54945455,  0.53946465,  0.52947475,  0.51948485,  0.50949495,\n",
    "        0.49950505,  0.48951515,  0.47952525,  0.46953535,  0.45954545,\n",
    "        0.44955556,  0.43956566,  0.42957576,  0.41958586,  0.40959596,\n",
    "        0.39960606,  0.38961616,  0.37962626,  0.36963636,  0.35964646,\n",
    "        0.34965657,  0.33966667,  0.32967677,  0.31968687,  0.30969697,\n",
    "        0.29970707,  0.28971717,  0.27972727,  0.26973737,  0.25974747,\n",
    "        0.24975758,  0.23976768,  0.22977778,  0.21978788,  0.20979798,\n",
    "        0.19980808,  0.18981818,  0.17982828,  0.16983838,  0.15984848,\n",
    "        0.14985859,  0.13986869,  0.12987879,  0.11988889,  0.10989899,\n",
    "        0.09990909,  0.08991919,  0.07992929,  0.06993939,  0.05994949,\n",
    "        0.0499596 ,  0.0399697 ,  0.0299798 ,  0.0199899 ,  0.01      , 0.])\n",
    "rates = np.array([-0.06048653, -0.0731699 , -0.08604612, -0.09886635, -0.11142732,\n",
    "       -0.12356894, -0.13516981, -0.14614191, -0.15642518, -0.16598255,\n",
    "       -0.17479537, -0.18285957, -0.19018233, -0.19677939, -0.20267278,\n",
    "       -0.20788904, -0.21245778, -0.2164105 , -0.21977972, -0.22259831,\n",
    "       -0.22489892, -0.22671361, -0.22807356, -0.22900885, -0.22954833,\n",
    "       -0.22971954, -0.22954864, -0.22906041, -0.22827824, -0.22722418,\n",
    "       -0.22591895, -0.22438198, -0.2226315 , -0.22068456, -0.2185571 ,\n",
    "       -0.21626401, -0.21381918, -0.21123557, -0.20852525, -0.20569948,\n",
    "       -0.20276872, -0.19974272, -0.19663055, -0.19344064, -0.19018084,\n",
    "       -0.18685842, -0.18348016, -0.18005235, -0.17658083, -0.17307102,\n",
    "       -0.16952795, -0.1659563 , -0.16236039, -0.15874425, -0.15511159,\n",
    "       -0.15146588, -0.1478103 , -0.14414782, -0.14048117, -0.13681289,\n",
    "       -0.13314533, -0.12948065, -0.12582084, -0.12216776, -0.11852311,\n",
    "       -0.11488844, -0.11126521, -0.10765475, -0.10405826, -0.10047686,\n",
    "       -0.09691158, -0.09336335, -0.08983303, -0.0863214 , -0.08282915,\n",
    "       -0.07935693, -0.07590533, -0.07247485, -0.06906598, -0.06567912,\n",
    "       -0.06231465, -0.05897291, -0.05565417, -0.0523587 , -0.04908671,\n",
    "       -0.04583838, -0.04261387, -0.03941331, -0.03623679, -0.03308441,\n",
    "       -0.02995619, -0.0268522 , -0.02377243, -0.02071689, -0.01768555,\n",
    "       -0.01467838, -0.01169533, -0.00873634, -0.00580133, -0.00289022, -0.])\n",
    "import scipy.interpolate\n",
    "# Check the documentation for this. There's a lot of handy stuff in scipy!!\n",
    "fitted_rate_function = scipy.interpolate.interp1d(concentrations, rates, 'quadratic')\n",
    "# We could just use 'fitted_rate_function' directly,\n",
    "# but instead we'll wrap it in a function of our own called 'rate',\n",
    "# that handles a few common errors making it easier to debug things later.\n",
    "def rate(Ca):\n",
    "    \"\"\"\n",
    "    Returns rate $r_A$ as a function of concentration $C_A$.\n",
    "    If Ca is provided in mol/m3 then rate is returned in mol/m3/s.\n",
    "    \"\"\"\n",
    "    if not isinstance(Ca, float):\n",
    "        raise TypeError('Function only works on one number at a time! (Give it a float, not a {0})'.format(type(Ca)))\n",
    "    if Ca == 0: return 0.\n",
    "    if Ca < 0: raise ValueError('Negative concentration requested!')\n",
    "    if Ca > 1: raise ValueError('Rate unknown for concentrations above 1 mol/m3')\n",
    "    return float(fitted_rate_function(Ca))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First let's check that  it gives a sensible number for the initial rate. It should be negative, because A is consumed"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-0.060486529999999074"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rate(Ca0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "OK. And when there's no A left ($C_A = 0$) the reaction should have stopped ($-r_A = 0$)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rate(0.0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 1.\n",
    "\n",
    "Using the  variables `v` and `Ca0` defined above along with the `rate(Ca)` function that gives the rate law, reproduce the Levenspiel plot from the exam, as close as possible.\n",
    "Here is some code that may help with some of the formatting:\n",
    "```python\n",
    "plt.xlim(0,1)\n",
    "plt.xlabel(r'$X$',fontsize=20)\n",
    "plt.ylabel(r'$\\frac{F_{A0}}{-r_A}}$',fontsize=24, rotation='horizontal', labelpad=15)\n",
    "plt.grid()\n",
    "plt.tight_layout()\n",
    "```\n",
    "\n",
    "and a picture of what you're aiming for:\n",
    "![levenspiel](https://raw.githubusercontent.com/rwest/CHME4510/master/images/HW3-levenspiel.png)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Your graph goes here"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Question 2\n",
    "Again using `v`, `Ca0`, and the function `rate(Ca)` defined above, model the network of two 1 m<sup>3</sup> CSTRs in parallel followed by a PFR of volume $V_P$. \n",
    "![network](https://raw.githubusercontent.com/rwest/CHME4510/master/images/HW3-network.png)\n",
    "\n",
    "Define a function `Vprf(Xfinal)` that returns volume of the PFR $V_{PFR}$ required to reach the final conversion $X_{final}$.\n",
    "\n",
    "Evaluate `Vpfr(0.8)` to see if it's close to your estimate from the exam. (My estimate was 3.3 m<sup>3</sup>)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3.3"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# These could be helpful. Check their documentation.\n",
    "from scipy.optimize import fsolve\n",
    "from scipy.integrate import quad\n",
    "\n",
    "def Vpfr(Xfinal):\n",
    "    \"\"\"\n",
    "    Return the Volume of the PFR required to reach a final conversion of Xfinal,\n",
    "    when it is placed in series after two 1m3 CSTRs operating in parallel.\n",
    "    \"\"\"\n",
    "    volume = 3.3 # my guess from the exam. Replace this with some calculations!\n",
    "    return volume\n",
    "\n",
    "Vpfr(0.8)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now check that your notebook runs correctly if you restart the kernel, clear the output, and run all the cells from the top."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}