Skip to content
Permalink
d3fa6b55f1
Switch branches/tags

Name already in use

A tag already exists with the provided branch name. Many Git commands accept both tag and branch names, so creating this branch may cause unexpected behavior. Are you sure you want to create this branch?
Go to file
 
 
Cannot retrieve contributors at this time
368 lines (368 sloc) 73.4 KB
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Computational Mechanics Project #01 - Heat Transfer in Forensic Science\n",
"\n",
"We can use our current skillset for a macabre application. We can predict the time of death based upon the current temperature and change in temperature of a corpse. \n",
"\n",
"Forensic scientists use Newton's law of cooling to determine the time elapsed since the loss of life, \n",
"\n",
"$\\frac{dT}{dt} = -K(T-T_a)$,\n",
"\n",
"where $T$ is the current temperature, $T_a$ is the ambient temperature, $t$ is the elapsed time in hours, and $K$ is an empirical constant. \n",
"\n",
"Suppose the temperature of the corpse is 85$^o$F at 11:00 am. Then, 2 hours later the temperature is 74$^{o}$F. \n",
"\n",
"Assume ambient temperature is a constant 65$^{o}$F.\n",
"\n",
"1. Use Python to calculate $K$ using a finite difference approximation, $\\frac{dT}{dt} \\approx \\frac{T(t+\\Delta t)-T(t)}{\\Delta t}$. "
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"K = 0.6111111111111112\n"
]
}
],
"source": [
"T0 = 85; T = 74; Ta = 65\n",
"dt = 2\n",
"dTdt = (T - T0) / dt\n",
"K = -dTdt/(T-Ta)\n",
"print('K =',K)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"2. Change your work from problem 1 to create a function that accepts the temperature at two times, ambient temperature, and the time elapsed to return $K$. "
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"def measure_K(Temp_t1,Temp_t2,Temp_ambient,delta_t):\n",
" ''' \n",
" \n",
" Determine the value of K based upon temperature of corpse \n",
" when discovered, Temp_t1\n",
" after time, delta_t, Temp_t2\n",
" with ambient temperature, Temp_ambient\n",
" \n",
" Arguments\n",
" ---------\n",
" Temp_t1: starting body temp\n",
" Temp_t2: current body temp\n",
" Temp_ambient: ambienet temperature\n",
" delta_t: time elapsed between starting and current temp\n",
" \n",
" Returns\n",
" -------\n",
" K: emperical constant for Newton's Law of Cooling\n",
" \n",
" '''\n",
" dTdt = (Temp_t2 - Temp_t1) / delta_t\n",
" K = -dTdt/(Temp_t2 - Temp_ambient)\n",
" return K "
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.6111111111111112"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"measure_K(85,74,65,2)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"3. A first-order thermal system has the following analytical solution, \n",
"\n",
" $T(t) =T_a+(T(0)-T_a)e^{-Kt}$\n",
"\n",
" where $T(0)$ is the temperature of the corpse at t=0 hours i.e. at the time of discovery and $T_a$ is a constant \n",
" ambient temperature. \n",
"\n",
" a. Show that an Euler integration converges to the analytical solution as the time step is decreased. Use the constant $K$ derived above and the initial temperature, T(0) = 85$^o$F. \n",
"\n",
" b. What is the final temperature as t$\\rightarrow\\infty$?\n",
" \n",
" c. At what time was the corpse 98.6$^{o}$F? i.e. what was the time of death?"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"a)\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOzdd3yV9dn48c+Vkx0SICEBkoDsDbIEZCOCDMFCoQWtorjax1Gtvypqrfr4PFatttJHrHXiBEGhWIbgwC2ylQ1hJCRhb8hOrt8f90lM4CQ5gRxOxvV+vc7rnHt87/s6jFy5v1NUFWOMMeZsAf4OwBhjTNVkCcIYY4xHliCMMcZ4ZAnCGGOMR5YgjDHGeBTo7wAqU4MGDbRZs2b+DsMYY6qNNWvWHFbVWE/HalSCaNasGatXr/Z3GMYYU22ISHJpx6yKyRhjjEeWIIwxxnhkCcIYY4xHliCMMcZ4ZAnCGGOMRz5NECJyr4hsEpGNIjJLREJF5DERSROR9e7XqFLKjhCRbSKSJCLTfBmnMcaYc/msm6uIJAB3Ax1UNVNE5gCT3If/rqrPllHWBcwAhgGpwCoR+UhVN/sqXmMATp48ycGDB8nNzfV3KMZcsKCgIOLi4oiKijqv8r4eBxEIhIlILhAOpAPNvCjXC0hS1V0AIjIbuAao9ASRm5PNmlmPE9GsJ50Hja/sy5tq5OTJkxw4cICEhATCwsIQEX+HZMx5U1UyMzNJS0sDOK8k4bMqJlVNA54FUoB9wAlVXeY+fKeI/CQir4tIfQ/FE4C9xbZT3fvOISK3ichqEVl96NChCscZGBhE+90zyfhxfoXLmprl4MGDJCQkEB4ebsnBVHsiQnh4OAkJCRw8ePC8ruGzBOH+wX8N0ByIByJE5DfAP4GWQFecxPGcp+Ie9nlc2UhVX1bVnqraMzbW42jxsuMMCGBfcDMiTyZVuKypWXJzcwkLC/N3GMZUqrCwsPOuMvVlI/WVwG5VPaSqucA8oK+qHlDVfFUtAF7BqU46WyrQpNh2Ik71lE+cjmpNQm4y+fkFvrqFqSbsycHUNBfyb9qXCSIF6CMi4eJEOBTYIiKNi50zDtjooewqoLWINBeRYJzG7Y98FairYXvqyhlSUnb56hbGGFPt+LIN4gfgA2AtsMF9r5eBZ0Rkg4j8BAwB7gUQkXgRWewumwfcCSwFtgBzVHWTr2Ktd0kXAPYnrfPVLYwxptrx6TgIVX1UVdupaidVvV5Vs93vnVW1i6qOVdV97nPTVXVUsbKLVbWNqrZU1f/1ZZzxbboBkJHq6WHGmOpHVWnevDkiQlKS79rXHnvsMRo0aFChMitXruSxxx6rlGt5Y/DgwUyYMKHSr1sb2EhqIKRuI05IJIFHtvk7FGMqxffff8+ePXsAmD17tn+DOcvKlSt5/PHHz9l/yy23sHTpUj9EZEpjCQJAhENhLah/xtogTM0wa9YsIiIi6N27N7NmzfJ3OF5JTEykR48e/g7DFGMJwi27fhuaFaRw7HS2v0Mx5oLk5+czd+5cxo4dy9SpU9m8eTM//fRT0fGZM2ciImzYsIFhw4YRERFBu3btmDdvXonrLFq0iGHDhhWNxO3Tpw/Lli07+3ZF8vLyiI+P9/h0MGjQIMaPH8/MmTO56667AKd3jYgwePBgwHMV05EjR7j99ttp3LgxoaGhtG3blueff77o+HPPPcdll11G3bp1adiwIWPGjPFplVptYwnCLaRxB6Ikk507t/s7FGMuyOeff86BAweYNGkSEyZMICgoyONTxLXXXsvYsWOZP38+rVu3ZtKkSaSmphYd3717N2PGjOHtt9/mww8/pG/fvowcOZJvv/3W430DAwOZMmUKM2fORPXnYUu7du3i66+/5qabbmL06NHcd999gFMN9v333/Piiy96vF5mZiaDBw/m3//+N4888giLFy/mvvvuIz395x7vqamp3HnnnSxYsIBXXnmF/Px8+vXrx4kTJ87rz86UVKOWHL0QsS26wmo4tPtHuLSzv8MxVcTj/9nE5vSTfrl3h/goHh3TscLlZs2aRb169RgxYgTBwcEMGzaM2bNn8+STT5boE3/vvfcydepUAHr06EHDhg1ZuHAhv/3tbwG48847i84tKChgyJAhbNq0iddee41+/fp5vPfUqVN56qmn+OKLLxgyZAjgPLHExcUxcuRIAgMDKVw3vk+fPmV+j7feeotNmzaxdu1aunbtCsAVV1xR4py///3vRZ/z8/OLnngWLFjADTfc4M0flymDPUG41XV3dc3Z57PetMb4XHZ2NvPnz2fcuHEEBwcDMHnyZPbs2cOKFStKnDt8+PCizzExMcTFxZV4gkhNTWXKlCkkJCQQGBhIUFAQy5YtY/v20p+yW7duzcCBA5k5cybg9KZ66623uP766wkMrNjvo59//jndunUrSg6erFixgmHDhhETE0NgYCDh4eGcPn26zBiN9+wJolBEA04G1CP02A5/R2KqkPP5Dd6flixZwvHjxxk1ahTHjx8HnG6eISEhzJo1i8svv7zo3Hr16pUoGxwcTFZWFuA8MYwdO5ZTp07x3//937Rq1YqIiAj+/Oc/lzuvz80338zvfvc7XnjhBVauXElycjI33XRThb/LkSNHaNy4canHU1JSGD58OL169eJf//oX8fHxBAcHM3r06KLvYS6MJYhijkW0oOGJ3eTkFRAcaA9XpvopbGuYOHHiOcfmzJlTokqmLElJSaxbt44lS5YwYsSIov2ZmZnllp04cSJ33303c+fOZfny5fTu3ZsOHTp4+Q1+FhMTU2aD88cff0xGRgYLFiwgIiICcBrKjx49WuF7Gc/sp2Ax+Q3a0VLS2HnwlL9DMabCTp8+zcKFC5k8eTLLly8v8frb3/7GgQMHWL58uVfXKkwEISEhRfuSk5NLbaAuLiwsjMmTJzNjxgzmzZt3ztNDYdVXeb/lDx06lHXr1pXogXV2jAEBASWqrubMmUNeXl65MRrvWIIopk5iJyIlkz27bMCcqX4WLFhARkYGv//97xk8eHCJ15133klMTIzXYyLatWtHYmIi9913H4sWLWL27NkMHz6chASPs+6f4+abb2bt2rWoKpMmTSpxrF27dgBMnz6dVatWsW2b5/9vN9xwA507d2b48OG89NJLLF++nNdff51p05wFJq+44gry8/O56aab+Oyzz/jHP/7BtGnTzqk6M+fPEkQx0c2dhuoTyRv8HIkxFTdr1ixat25N7969zzkWFBTEr371K+bNm0d2dvljfUJCQpg3bx6BgYFMmDCBRx55hAcffJBBgwZ5FUvPnj1JSEhg/Pjx1K1bt8SxAQMG8Mc//pHp06fTu3dvbr/9do/XCA0N5fPPP2fMmDH8+c9/ZuTIkTzzzDPEx8cD0LlzZ9544w1++OEHrr76at577z3mzp17zv3M+ZPi/ZWru549e+rq1avP/wIZR+GZ5rxX91auvbfUFVFNDbVlyxbat2/v7zBqhM2bN9OxY0c+/fRThg4d6u9war2y/m2LyBpV7enpmDVSFxcezanAaOqcTEJVbW0AYyroyJEjbNu2jUceeYROnTqdM27BVC9WxXSWU1GtaJqfzMFTNuWGMRX1n//8h/79+7Nv376iKT1M9WUJ4iwBce1pLWlsTj/u71CMqXZuvPFGCgoK2Lx5s028VwNYgjhL3Us6EyHZpO62kZjGmNrNpwlCRO4VkU0islFEZolIqIj8VUS2ishPIjJfRDz2SRORPe6V59aLyAW0PFdMWHwnAM7Y4kHGmFrOZwlCRBKAu4GeqtoJcOGsLf0J0ElVuwDbgQfLuMwQVe1aWgu7T8Q5fbQDDm+9aLc0xpiqyNdVTIFAmIgEAuFAuqouc685DbACSPRxDBUTVp/TQQ2IydhFVm6+v6Mxxhi/8VmCUNU04FkgBdgHnFDVs1cbmQosKe0SwDIRWSMit5V2HxG5TURWi8jqQ4cOVUboZNVvQ2tJZdt+m3LDGFN7+bKKqT5wDdAciAciROQ3xY4/DOQB75ZyiX6q2h0YCdwhIgM9naSqL6tqT1XtGRsbWymxBzfuQCtJZ4v1ZDLG1GK+rGK6EtitqodUNReYB/QFEJEpwNXAdVrKUG5VTXe/HwTmA718GGsJdZp0Jlyy2Z9sczKZ6sXTsp3FffHFF4gIGzdWj04YgwcPRkT4y1/+cs6xBg0a8Nhjj138oM6TiPDCCy9U6jV9/ffpywSRAvQRkXBxRssMBbaIyAjgAWCsqmZ4KigiESISWfgZGA5ctH/RAXHOkPTs9M0X65bGXBTdu3fn+++/p2XLlv4OpUL+/ve/k5Hh8cdFtfH99997nIa9KvNlG8QPwAfAWmCD+14vAy8AkcAn7i6sLwGISLyILHYXbwh8IyI/AiuBRar6sa9iPUdsWwBCjm+nJs1VZUxUVBR9+vQhLCzM36EUKW+Nicsvv5yTJ0/y8ssvX6SIKlfh9+vTpw8NGzb0czQV49NeTKr6qKq2U9VOqnq9qmaraitVbeLuvtpVVX/rPjddVUe5P+9S1Uvdr46q+r++jPMcYfXICImjaX4KqcfKXyDFmOrCU5WEiDB9+nQeeughYmNjiYuL44477jhn1teUlBQmTZpEdHQ04eHhXHXVVedM1T1t2jQ6d+5MnTp1SExM5LrrrmP//v0lzmnWrBn33XcfTzzxBImJiURFRZUZc3x8PDfddBPPPvtsmTPRDh48mAkTJpT5fffs2YOIMHv2bG666SaioqJITEzknXfeASiaLTY2NpYHHniAgoKCEtfbuHEjo0ePJjIyksjISCZOnFji+xXeb+nSpYwdO5Y6deoUre3tqYpp/vz59OrVi7CwMGJiYhg1ahTJyckAbN26lUmTJtGkSRPCw8Pp2LEjzz///Dkx+ZKNpC5FXoO2tJFUtuzzz4L1xlxMzz33HOnp6bzzzjv88Y9/5F//+hfTp08vOn706FH69+/Ptm3beOmll5gzZw5nzpzhyiuvLPEEcPDgQR566CEWLVrE888/z65du4rWbSjuvffe48svv+TFF1/k/fffLze+Bx54gAMHDvDGG29Uyvd94IEHaNy4MR9++CEDBgxgypQp3HfffaxcuZLXX3+de+65h2eeeYY5c+YUlUlKSqJfv35kZWXx9ttvM3PmTDZt2sSYMWPOqWm4+eabufTSS/noo4+4+eabPcbw9ttvM378eFq2bMmcOXN44403aNOmDYW9MdPS0mjbti0vvvgiixcv5tZbb+XRRx/l6aefrpQ/A2/YbK6lCEvoRKvUFXyWfozhHRv5OxzjL0umwX4/rQ/SqDOMfOqi3KpZs2bMnDkTgKuuuopvv/2WefPmcf/99wNOG8CZM2dYv3490dHRAPTr149mzZrx+uuvc8cddwDw+uuvF10zPz+fyy+/nMTERL799lsGDizZEXHhwoWEhoZ6Hd91113H008/zS233FJiFbnzccUVV/Dkk08C0Lt3bz744AM++ugjtm7disvlYsSIESxYsID58+cXLXj0+OOP06hRI5YsWVK0Kl6XLl1o164dixcvZvTo0UXXnzhxIk888USp9y8oKGDatGmMGzeuxCJOY8eOLfo8dOjQoqnSVZX+/fuTkZHBK6+8woMPljW+uPLYE0Qpghp1IFRyObTX5mQyNd/w4cNLbHfo0IHU1NSi7U8//ZRhw4YRFRVFXl4eeXl5REZG0qNHD4qvwbJkyRL69u1L3bp1CQwMJDHRGQe7fXvJ/0dDhw71OjkUeuihh0hJSeHdd0vrGe+94mtUREVFERsby6BBg3C5XEX7W7VqRVpaWtH2p59+yrhx4wgICCj6M2jevDnNmjXj7HVoiicLT7Zt20Z6evo5y7EWl5WVxaOPPkqrVq0ICQkhKCiIhx9+mN27d1+0ZVXtCaI0sU5Pppz0zTjDOUytdJF+g/e3s5fpDA4OLrFm9OHDh1mxYoXH6qDCH7arVq1i7NixjBs3jmnTphEXF4eI0KdPn3PWnz6fxto2bdowYcIE/vKXv3D99ddXuHxxnr6vN38GTz/9tMcqnr1795bYLu/7HTlyBIDGjRuXes4DDzzAq6++yqOPPkr37t2pV68eCxYs4H/+53/IysqiTp06Zd6jMliCKI27J1ODzF3sO5FJ47pVp9eHMRdbdHQ0Y8eO5ZFHHjnnWGRkJOA0uMbGxvL+++8XrQNR2OB6tvNdJ+Lhhx+ma9eufPDBB+ccCw0NJScnp8S+o0ePntd9PImOjmbcuHHccsst5xw7e+xJed8vJiYGgH379pV6zty5c7nrrruKqvkAFi1aVJGQL5gliNKERpETEU+bk6msTT7O6C6WIEztNXToUObMmUPHjh1L7SKbmZlJUFBQiR+OlVEdVFyXLl0YM2YMTz755DkNw4mJiXz11Vcl9n3yySeVdu+hQ4eyceNGevToccELIbVt25aEhATefPNNxowZ4/GczMxMQkJCirbz8/OZPXv2Bd23oixBlCGwcUc6nN7C+ynHGN2l9EdBY6qSnJwcj79hDxo06Lyv+Yc//IF33nmHK664grvuuouEhAQOHDjAl19+Sf/+/Zk8eTLDhg3j+eef55577mHMmDF89913Rd1HK9PDDz9M7969z9k/btw4XnvtNe69915Gjx7N8uXLWbp0aaXd97HHHqNXr16MHj2aqVOn0qBBA9LS0vjkk0+48cYbGTx4sNfXCggI4JlnnuG6667juuuuY/LkyYgIn3/+OZMnT6Znz54MGzaMGTNm0KpVK6Kjo5kxY0aZ3Xx9wRqpyxCQ0J1Wks7mPen+DsUYr506dYqJEyee89q0adN5X7NBgwasWLGCdu3ace+99zJ8+HDuv/9+Tpw4QZcuXQAYNWoUTz/9NB9++CFjx47lyy+/ZOHChZX1tYr06tWLYcOGnbN/9OjRPPnkk3zwwQeMGzeO5ORknn/++Uq7b5s2bVixYgXh4eHcdtttjBw5kkcffZSQkBBatWpV4etde+21fPjhh2zdupUJEyZwww03sHXrVgrnlPu///s/BgwYwB133MHUqVPp1KnTReu9VEhq0kjhnj176tm9CS7Ito9h1q+5Nu9R3njs94QEusovY6qtLVu20L59e3+HYUylK+vftoisKW3NHXuCKEtCdwDaaxIb02zAnDGmdrEEUZY6ceRHJnBpwC7WpRzzdzTGGHNRWYIohyuxO90Dd7MuxdaGMMbULpYgyhPfnUTdT1Jyir8jMcaYi8oSRHnc7RBxp7ew74TN7FrT1aROG8bAhf2btgRRnsZdAegiu1ibbNVMNVlQUFC5axMYU90UDmA8H5YgyhNWD41uSTfXLtZaQ3WNFhcXR1paGhkZGfYkYao9VSUjI4O0tDTi4uLO6xo+HUktIvcCtwCKs6rcTUA48D7QDNgD/EpVz/nJ616adDrgAl5VVb/NmiYJ3el+/HNmWIKo0QoXrklPTyc3N9fP0Rhz4YKCgmjYsGG5izKVxmcJQkQSgLuBDqqaKSJzgElAB+AzVX1KRKYB03DWqC5e1gXMAIYBqcAqEflIVf2zSHR8d6I3zOVQ2h6y8/rYgLkaLCoq6rz/MxlT0/i6iikQCBORQJwnh3ScubPfdB9/E/iFh3K9gCT30qM5wGz8Oee2u6G6ne5kU7oNmDPG1A4+SxCqmgY8C6QA+4ATqroMaKiq+9zn7AM8VY4lAMUnWE917zuHiNwmIqtFZHXhUn2VrlEXVFx0CdjJ2mSrZjLG1A7lVjGJSFdgABAPZAIbcaqITpRTrj7Ob/3NgePAXBH5jZdxeZpL12Oroaq+DLwMzlxMXl6/YoLDkbj29DqUzNs2YM4YU0uU+gQhIr8RkTXA40B9IBk4CVwJfCEir4lIYhnXvhLYraqHVDUXmAf0BQ6ISGP3PRoDBz2UTQWaFNtOxKme8p/4bnRmJ2uTK28BEmOMqcrKeoKIAQaq6hlPB0WkJ9Ae54e5JylAHxEJx3nyGAqsBs4AU4Cn3O8LPJRdBbQWkeZAGk7j9rXlfhtfSuhOxLq3CTyVYivMGWNqhVKfIFR1emnJwX18taqWulyTqv4AfACsxeniGoBTFfQUMExEduD0UnoKQETiRWSxu2wecCewFNgCzFHV85/MvjLEOw3Vl9qAOWNMLVHqE4SILFHVke7P96vqMxW9uKo+Cjx61u5snKeJs89NB0YV214MLK7oPX0mrgPqCqFboDNgzlaYM8bUdGX1YmpU7PMkXwdS5QUGI4060Sckxab+NsbUCmUlCJtr4Gzx3WmVv5PNacfJzsv3dzTGGONTZTVStxCReThdTgs/F1HV8T6NrCpK6E7IqldILEhlU/pJujet7++IjDHGZ8pKEL8s9vkFXwdSLZRoqD5mCcIYU6OVmiBU9bOLGUi10KA1BNehryuZz23AnDGmhvNmJPU6zm2POIEzpuEvqlp7Ro4FuKBxV3oc2MNfk4+hqoh4GvRtjDHVnzdzMX0CfAbc7H59AnwPHANm+iyyqiqhG02ykzhy8jTJRzL8HY0xxviMN9N991XV/sW214nIN6raX0Q2+CqwKiu+Oy7NpY3s5eukwzRrEOHviIwxxie8eYKIFJEehRsi0h0onDA/zydRVWXuqb8HRezl2x2H/RyMMcb4jjdPELcDb4tIEE6X1xzgZhGJACo8urraq3cJhEUzOHQv7+w8TH6B4gqwdghjTM1TboJQ1RVABxGJAURVi//aPMtnkVVVIhDfjbaHkjiZlceGtBN0bVLP31EZY0ylK7eKSURiReRfwJuqelhEOojIjb4PrQpL7EnUqR1EksE3O3y0SJExxviZN20QM4Ev+Xl9hh3Afb4KqFpo1h/RAn4Zk8I3SdYOYYypmbxJEHGq+h5QAOBe/Kd2T0SU2AtcIYyM3MGa5GNk5NS+tnpjTM3nTYI4IyLRuAfLichlwCmfRlXVBYVCk150yP6R3Hxl5e7aM1bQGFN7eNOL6f8B/8GZsO9LIAGYUF4hEWkLvF9sVwvgz8DlQFv3vnrAcVXt6qH8HpxElA/kqWpPL2K9eJoNoM4XfyE2MINvdhxmcNs4f0dkjDGVypteTKtFZAjO8qICbFbVHC/KbQO6AoiIC2fp0Pmq+nzhOSLyHM60HaUZclavqaqj+QDkiyeZHLeXZUkN/R2NMcZUurJWlBtbyqGmIoKqflSB+wwFdqpqcrHrC/Ar4IoKXKfqSOgBgWFcGbaNf+xsy8FTWcRFhvo7KmOMqTRlPUFMdL83APoCy3GeIAbh9GqqSIKYxLljJgYAB1R1RyllFFgmIgr8S1Vf9nSSiNwG3AbQtGnTCoR0gQJDoGlvWh9bD4zlu6Qj/KJbwsW7vzHG+FipjdSqer2qXg/kAh1U9Reqeg3QkQpMsSEiwcBYYO5ZhyZT9kC7fqraHRgJ3CEiA0uJ82VV7amqPWNjY70Nq3I0H0jYsa00D8uw7q7GmBrHm15MLVQ1rdh2Oj83MntjJLBWVQ8U7hCRQGA8JRuxS1DVdPf7QWA+0KsC97w4mjk569qGe/lmx2FUbZVWY0zN4U2C+EpEFonIb0TkOpyqpa8qcA9PTwpXAltVNdVTARGJEJHIws/AcGBjBe55ccR3heA6DAreyv6TWew8dMbfERljTKXxJkHcgTOaujfQB3jLva9cIhIODAPmnXXonDYJEYkXkcXuzYbANyLyI7ASWKSqH3tzz4vKFQRNL6f5qbUANu2GMaZG8aabq+K0H5zdhlAuVc0AYjzsv9HDvnRglPvzLuDSit7PL5oPICjpE7rVz+KbpCPc2K+5vyMyxphKUeoThIgsF5HfiUj8WfsDRWSgiLwmIjf5PsQqrtkAAH7dIJkVu46Qm1/g54CMMaZylFXFNBoIAuaLSKqI/CQiO4BdwE3AP1X1jYsRZJXW+FIIqUtf12ZOZ+fxU+pxf0dkjDGVotQqJnf10D+Af4hICBAHZFbZkc3+EuCCS/qScHA1IuP4esdhelwS7e+ojDHmgnnTSI2qZqvqXksOpWg+ANfx3QxulMs3tgypMaaG8CpBmHK42yF+Gb2LdXuPcyor188BGWPMhbMEURkadoKw+vRiE/kFylfb7SnCGFP9eZUgRCTRPaMrIhLiHrxmCgUEQLP+xB5ZSXREMB9v2u/viIwx5oJ5syb1VJzR06+6d10CLPBlUNVSs4HI8RR+3aqAz7ccICu3di+6Z4yp/rx5grgbZwT1SQBV3Y7To8kU19xph7im3k7O5ORbY7UxptrzJkFkFV8gyL34j/gupGoqth1ExNI6Yz2RoYFWzWSMqfa8SRDfisj9QKi7HeJ9YKFvw6qGRKBZf1zJ33Bluzg+2XzARlUbY6o1bxLE/ThrQ28Ffg98Bjzsy6CqrWYD4GQa45tlcSIzlx92HfV3RMYYc97KnKzPXZ30uqpOAf55cUKqxloPA6BPzg+EBbVjycZ99G/dwM9BGWPM+SnzCUJV84HGIhJ0keKp3uo1hUZdCNqxhCHtYlm66QD5BbaIkDGmevKmimkX8LWIPCgidxe+fB1YtdXuatj7A9e0CuLw6WzWphzzd0TGGHNevEkQh4BPgHAgttjLeNJuNKAMYg3BrgCWbLDeTMaY6smbBYMeOZ8Li0hbSq453QL4M1APuBUn8QA8pKqLzyqOiIwApgMu4FVVfep84rjoGnaEepcQmrSEAa3/yNJN+3nk6vaIWM9gY0z14s1I6k9EZNnZr/LKqeo2Ve2qql2BHkAGMN99+O+Fx0pJDi5gBjAS6ABMFpEOFflifiPiVDPt+oKr20WSdjyTDWkn/B2VMcZUWLlPEMCfin0OBX4JZFfwPkOBnaqa7OVv0r2AJPfSo4jIbOAaYHMF7+sf7UbDihkMD9mIKyCcJRv30yWxnr+jMsaYCin3CUJVfyj2+lJV78b5AV4Rk4BZxbbvdK9Q97qI1PdwfgKwt9h2qnvfOUTkNhFZLSKrDx065OmUi69JbwiPIWLXUi5vEcPHG/fjLO1tjDHVhzdVTFHFXvVEZCjQ2NsbiEgwMBaY6971T6Al0BXYBzznqZiHfR5/wqrqy6raU1V7xsZWkbZzVyC0GQnblzKyQwy7D59h+4HT/o7KGGMqxJteTJuAje73dTijqG+twD1GAmtV9QCAqh5Q1XxVLferG1YAACAASURBVABewfPTSCrQpNh2IpBegXv6X7vRkH2CUVE7EYElG/f5OyJjjKkQbxJEC1VtqqpNVLW5ql4BfFuBe0ymWPWSiBR/+hiHk3zOtgpoLSLN3U8gk3CmHK8+Wg6BoHDqJ39Cz0vq8/FG6+5qjKlevEkQP3jYt9Kbi4tIODAMmFds9zMiskFEfgKGAPe6z40XkcUAqpoH3AksBbYAc1R1kzf3rDKCwqDlFbB1ESM6NmLr/lPsPnzG31EZY4zXSu3FJCJxOG0NYSLSmZ/bBaJwBs2VS1UzgJiz9l1fyrnpwKhi24uBc7rAVivtroatCxkTd4AngP/8mM7dQ1v7OypjjPFKWd1cRwNTcer/Xyy2/xRwXoPnap02V4G4iEv9lH6thjFn9V7uHNKKgAAbNGeMqfpKrWJS1TdUdQBws6oOKPYapapzSytnigmPhkv6wtZF/PqypqQey+S7nUf8HZUxxnjFm6k25ojIVUBHnIFyhfuf9GVgNUa7q+HjB7iq0WnqhgUxe1WKTQFujKkWvBkH8SIwBfgDEAb8Bmjl47hqjnZOs0pI0hLGdUtg2aYDHDuTU04hY4zxP296MfVX1WuBI+6J+3rjtEsYb7jXiHCqmZqQk1/A/HVp/o7KGGPK5U2CyCp8F5FG7u1mPouoJmp3NexdSfs6mVyaWJf3V+21qTeMMVWeNwlisYjUA54F1gN7gA98GVSN414jgi0f8evLmrLtwCl+TLUZXo0xVVuZCUJEAoAlqnrc3XOpOdBZVR+6KNHVFA07QlxHWP8eYy5tTFiQi/dXpfg7KmOMKVN5a1IX4CzaU7idqapHfR5VTSMC3a+H9LVEntjO6C6N+Wh9Omey8/wdmTHGlMqbKqZPROQan0dS03X+FQQEwbp3mHRZE87k5LNog03gZ4ypurxJEHcC80UkU0SOisgxEbGniIqKiHHaIn6cTY/ECFrERjBn1d7yyxljjJ94kyAaAEFAHSDWvV1FFl6oZrpdD5lHkW1LmHRZE1YnHyPp4Cl/R2WMMR55s6JcPjAReMD9uTHOYj+moloOgagEWPcO47snEhggvG9PEcaYKsqbkdQv4EzLXTgLawbwki+DqrECXND1Wtj5GQ3yD3Nl+4bMW5tGTl6BvyMzxphzeFPF1FdVb8c9YM7diynYp1HVZF2vAy2AH9/j172acORMDp9tOeDvqIwx5hzeJIhc93gIBRCRGMB+5T1f0c2h2QBY9w4DW8UQXzeUmd/t8XdUxhhzDm8SxAzgQyBWRB4HvgGeLq+QiLQVkfXFXidF5B4R+auIbBWRn0RkvnuUtqfye9wrz60XkdUV+lZVXbfr4dgeXCnfcfOAFvyw+yhrko/5OypjjCnBm0bqt4A/4Uy1cRSYqKqzvSi3TVW7qmpXoAdO28V84BOgk6p2AbYDD5ZxmSHua/Qs/6tUIx3GQkhdWPc2ky5rQr3wIF76cqe/ozLGmBK8eYIAcAG5QE4FyhQ3FNipqsmqusy95jTACmrjzLBBYdD5l7B5ARF6himXN+OTzQfYccC6vBpjqg5vejE9DMwC4nF+mL8nImX91u/JJPc1zjYVWFJKGQWWicgaEbmtjPhuE5HVIrL60KFDFQzLj7pdD3lZsPFDpvRtRliQi5e+3OXvqIwxpog3TwO/AS5T1T+p6sNAL+AGb28gIsHAWGDuWfsfBvKAd0sp2k9VuwMjgTtEZKCnk1T1ZVXtqao9Y2Or0fi9+G7OBH5r3yY6IphJvZqwYH0aaccz/R2ZMcYA3iWIZEouTRoIVORX3ZHAWlUt6sspIlOAq4HrtJSFEVQ13f1+EKftolcF7ln1FZvAjwObuGVACwBe+cqeIowxVYM3CSID2CQir4rIK8AG4LiI/E1E/uZF+ckUq14SkRHAA8BYVc3wVEBEIkQksvAzMBzY6MW9qpfCCfzWzCShXhjXdE1g9qoUjtqSpMaYKsCbBLEIeAz4HqdR+b+Bz4FN7lepRCQcGAbMK7b7BSASZ5bY9SLykvvceBFZ7D6nIfCNiPwIrAQWqerH3n6paiMiBjpPhHXvwJkj/G5wC7JyC2xchDGmSggs7wRVfe18L+5+Qog5a1+rUs5NB0a5P+8CLj3f+1Yr/e+BH2fBDy/R6oqHGd6hIW9+t4fbB7YgIqTcvx5jjPEZb3oxjRCRVSJy0Kb79oHYts404Cv/BVkn+e3glpzIzGXWSltxzhjjX95UMb0A3A4kYNN9+8aAP0DWCVjzBt2b1qdPi2he/Xo32Xn5/o7MGFOLeZMgUoH1qpqrqvmFL18HVqsk9IDmg+D7GZCbxe8Gt2L/ySwWrEv3d2TGmFrMmwRxP/AfEfmjiNxd+PJ1YLXOgPvg9AH48T0Gtm5Ax/goZnyRZFOBG2P8xpsE8TiQD9TDqVoqfJnK1Hyg8yTx7XSkIJ//d1Vbko9k8Nb3e/wdmTGmlvKmm0ycqvbweSS1nQj0/wO8fx1s/jdDOk9gYJtYpn+2g/HdE4mOsCU4jDEXlzdPEJ+JyBU+j8RA21EQ2w6+/huo8qfR7cnIyWf6p9v9HZkxphbyJkHcCnwqIqetm6uPBQRAv3vg4CbYvpQ2DSOZ3KsJ7/yQQtJBm+nVGHNxeZMgGgBBQF2sm6vvdZ4AdZvCN85TxL1XtiE82MX/Ltri78iMMbWMNwsG5QMTgQfcnxsDXX0dWK3lCoK+d8HeHyD5O2LqhHDXFa1Yvu0QX22vRtOZG2OqPW9GUr8ADAGud+/KAF7yZVC1XvfrISIWvvorAFP6NqNpdDj/s2gzefnW7dUYc3F4U8XUV1VvB7IAVPUoYF1qfCkoDPrfC7uWw/alhAS6eHBkO7YfOM37q/f6OzpjTC3hTYLIFZEAnBXeEJEYwH6N9bXLboWY1vDxg5CXw4hOjejVPJq/LdvOyaxcf0dnjKkFSk0QIlI4RmIG8CEQKyKPA98AT1+E2Gq3wGAY8RQc3Qk//BMR4ZHRHTiakcOM5Un+js4YUwuU9QSxEkBV3wL+BDwLHAMmqursixCbaX0ltBkBX/4VTh2gc2JdxndL5I1v9rBtv3V7Ncb4VlkJQgo/qOomVZ2uqs+rqlcru4lIW/eCQIWvkyJyj4hEi8gnIrLD/V6/lPIjRGSbiCSJyLQKfq+a46onIS8LPnscgIdGtSMyNJA/zFlPrjVYG2N8qKwEESsifyjtVd6FVXWbqnZV1a5AD5zeT/OBacBnqtoa+My9XYKIuHCqtkYCHYDJItKh4l+vBohpCZf/F6x/F1LXEFMnhCfHd2ZT+kn+73OrajLG+E5ZCcIF1MFZHtTTqyKGAjtVNRm4BnjTvf9N4Bcezu8FJKnqLlXNAWa7y9VOA/8IdRrCkvuhoICrOjZifLcEZixP4se9x/0dnTGmhiprsr59qvrflXSfScAs9+eGqroPQFX3iUich/MTgOL9OVOB3pUUS/UTEglXPgb//h389D50ncyjYzvy/a4j3Df3Rxbe1Z/QIJe/ozTG1DBetUFcCBEJBsYCcytSzMM+LeX6t4nIahFZfehQDR5p3GWSMx34p49C9inqhgXx9C+7kHTwNM8u3ebv6IwxNVBZCWJoJd1jJLBWVQ+4tw+ISGMA9/tBD2VSgSbFthMBj8urqerLqtpTVXvGxtbgKaICAmDkM86iQl89C8DANrH8pk9TXvt2Nyt2HfFzgMaYmqbUBOEeMV0ZJvNz9RLAR8AU9+cpwAIPZVYBrUWkufsJZJK7XO2W2BMuvdZZmnTfjwA8NKo9TaPD+eMHP3I6O8/PARpjahJvRlKfNxEJB4YB84rtfgoYJiI73Meecp8bLyKLAVQ1D7gTWApsAeao6iZfxlptXPW/ENEAPrwFcjIIDw7kuYmXknos02Z8NcZUKp8mCFXNUNUYVT1RbN8RVR2qqq3d70fd+9NVdVSx8xarahtVbamq/+vLOKuV8GgY9xIc3gFLHwKgZ7NobhvQglkrU1i8YZ+fAzTG1BQ+TRDGR1oMhn53w5o3YMtCAP4wvA3dm9bjvjk/sin9RJnFjTHGG5Ygqqshf4LGXeGjO+FkOiGBLl66vgf1woO47a01HD6d7e8IjTHVnCWI6iowGH75KuRlw/zboaCAuMhQXrmhJ0fOZPPbt9eQk2dTcRhjzp8liOqsQWtnxtfdX8F3/wCgU0Jdnp14KauTj/HIvzei6nH4iDHGlMsSRHXX/QZoPwY+fwLS1wFwdZd47rqiFe+v3svM7/b4Nz5jTLVlCaK6E4Ex/4CIOPjgZsh2pgG/98o2DO/QkCcWbubrHTV4hLkxxmcsQdQE4dEw/mU4tgfm3gj5uQQECH//dVfaNIzkjnfXsvPQaX9HaYypZixB1BTNB8DVf4OkT2HhvaBKREggr9zQkyBXANe+soLdh8/4O0pjTDViCaIm6XGjMzX4urfhy2cAaBIdzru39iY3X5n88gr2WJIwxnjJEkRNM+RhuHQyfPEkrHsXgHaNonj3lt5k5+Uz+ZUVJB+xJGGMKZ8liJqmsNG6xWD4z92Q9BkA7RtH8e4tfcjMzWfyyytIOZLh1zCNMVWfJYiaKDAYfvU2xLaDOTfAvp8A6BDvPEmcyXGeJPYetSRhjCmdJYiaKjQKrpsLoXXh3YlwLBmAjvF1efeW3pzOzmPSy5YkjDGlswRRk0XFw3UfQF4mvDESDjkrz3VKqMs7N/fmVFYuv/znd7autTHGI0sQNV3DDnDjYijIg9dHQOoaADon1mXub/sSHBjAr/71PQt/8rhgnzGmFrMEURs06gRTP3aqnd4cAzuXA9C2UST/vqMfnRPqcud763j+0+02d5MxpoivV5SrJyIfiMhWEdkiIpeLyPsist792iMi60spu0dENrjPW+3LOGuF6BYwdSnUbwbv/Qo2/RuABnVCePfW3ozvnsDzn+7grlnryMrN92+sxpgqIdDH158OfKyqE9xrS4er6q8LD4rIc0BZq9sMUdXDPo6x9ohsBDctgvd+7UzJkfU89LiRkEAXz028lDYNI3n6463sPZrBKzf0JC4q1N8RG2P8yGdPECISBQwEXgNQ1RxVPV7suAC/Amb5KgbjQVh9uH4+tLoS/vN7WP4kFOQjIvx2UEte+k0Pdhw8zej/+4bl2w76O1pjjB/5soqpBXAIeENE1onIqyISUez4AOCAqu4opbwCy0RkjYjcVtpNROQ2EVktIqsPHbJZS70SHAGTZ8Gl18KXT8M74+G0kwyu6tiIef/Vl+jwYG56YxUPz99ARk6enwM2xviDLxNEINAd+KeqdgPOANOKHZ9M2U8P/VS1OzASuENEBno6SVVfVtWeqtozNja2kkKvBVxB8IsXnVHXKSvgpf6w+2vAmZpjwZ39uG1gC95bmcLof3zD2pRjfg7YGHOx+TJBpAKpqvqDe/sDnISBiAQC44H3Syusqunu94PAfKCXD2OtnUSgxxS45TMIiYS3xsKXf4WCfEKDXDw0qj2zbu1DTl4BE/75HX9bto3cfFvG1JjawmcJQlX3A3tFpK1711Bgs/vzlcBWVU31VFZEIkQksvAzMBzY6KtYa71GneC2L6DjeFj+PyWqnPq0iGHJPQMY1y2Rf3yexLgXv7WnCWNqCV+Pg7gLeFdEfgK6Ak+690/irOolEYkXkcXuzYbANyLyI7ASWKSqH/s41totJBJ++SqMme5UOb14Oax/D1SJCg3iuV9dyku/6c7Bk9mMf/E77pm9jn0nMv0dtTHGh6QmDYzq2bOnrl5tQyYu2IHNTg+n1JXQtC+Mfs4ZkQ2cyc7jxS+SeOXr3bhE+N3gltw2sAWhQS4/B22MOR8iskZVe3o8ZgnCeFRQAOvfgU/+7Kxz3ee/YNADEFIHgL1HM3hy8RaWbNxPQr0wHhzVjtGdG+P0XjbGVBeWIMz5O3MEPn3UWaUuKgFG/AXaj3UauIEVu47w+H82s2XfSdo1iuS/hrRidOfGuAIsURhTHViCMBcuZQUs/AMc3ASNuzpPE21Hggj5BcqC9Wm8+MVOkg6epllMOL8b3JJx3RIJDrTpvoypyixBmMqRnwc/zYav/grH9kCjLk6iaDcaRCgoUJZt3s+M5TvZkHaCxnVDuXVAC359WRMiQnw9q4sx5nxYgjCVKz8XNsx1EsXRXdCwMwz6I7S7GgJcqCpf7zjMjOVJ/LD7KBHBLsZ2TWByryZ0Tqhr7RTGVCGWIIxv5OfBxg+cRHEkCeo2hR43QLfrnYkBgbUpx3jvhxQW/pROVm4BHRpHMblXE67plkBUaJCfv4AxxhKE8a2CfNjyH1j9Ouz+EgICnfaJHjdBiyEQEMCJzFw+Wp/GrJV72bzvJKFBAVzVsRGjOjdmUJtY6yZrjJ9YgjAXz5GdsGYmrH8XMo446090+TV0uAbiOqDAhrQTzFq5lyUb93E8I5fwYBdXtItjVOfGDG4bS3iwtVcYc7FYgjAXX16281Sx9k3Y8w1oAcS0chJFh2ugURdyC5QVu46weMN+lm3az5EzOYQGBTCoTSwDWscysHUsTWPC/f1NjKnRLEEY/zp9ELYuhM0LnBljNR/qXeL0fmoxBC7pS15gOCv3HGXxhn18vuUg6SeyALgkJpz+rRowoHUsl7eMoW6YtVsYU5ksQZiq48wR2LbYSRZ7voa8LKfNIrEXtBgMLQajCd3ZdTSbr7cf4pukw3y/8whncvIJEGjTMJIel9QvejWNDrdeUcZcAEsQpmrKzYK9K2DXF84rfT2gEBQB8d0gsSck9iSnUXfWHQvlu51HWJtyjHUpxzmd7Sxi1KBOMN2a1qdzQl06NI6iQ3wUjeuGWtIwxkuWIEz1kHEUdn8Fyd9B2mrY9xMU5DrHohIhoTs06kx+bHv2uJqy4lgUa1JOsC7lOLsPnym6TL3wIDo0jqJ94yjaNoykRWwELWPrUD8i2E9fzJiqyxKEqZ5ys2D/Bkhd5SSMtDXOCO5CQeEQ2w4adiC7XktSacTm7FhWnazLj/tz2Lr/FNl5Py9wFB0RTMvYCFo0qEOzBhE0iQ6jSf1wmkSHUz88yJ46TK1kCcLUHNmn4dBWOLAJDm52v2+BjMMlz4tKQOs350x4AocCYtlbEMP2rHpsOBXJqqNhpGeUTAYRwS6aRIeTWD+MhlGhNK4b6n4Po1HdEBpGhVInJNCSiKlxykoQ1uHcVC8hdYraJkrIPO5M+1HsJUd2Uif1K+qc2k9zlOKLmmu9+uSExXI6KJpjAdEcLKhLam4kyQci2LMnjA1ZYRyjDsc0ktOEAUJIYAAN6oTQIDKEBhHBNKgTQkydYKIjgqkXHkz98CDqhQe5PwcTFRpIoMsmKzTVl08ThIjUA14FOgEKTAWuAm4FDrlPe0hVF3soOwKYDriAV1X1KV/Gaqq5sHpOG0VC93OP5eXAqXQ4kQon0uBECnJqPyGnDxBy+iAxpzfT6tQByCu2Ql7Izx8LJJDMoHpkBkRwWiI4cSqcY8fDOZIfysHcUI4WhJJCGGc0jNOEcoYwzmgoZwhFgsIIDI0gMDSCkNA61Al3nkQiQgKpE+IiIiSQiGBnOyLERXhwIOHBLkKDXIQHO68w93ZooIsgl9hTjLlofP0EMR34WFUniEgwEI6TIP6uqs+WVkhEXMAMYBiQCqwSkY9UdXNpZYwpVWCwM6K7frPSz1GFnNPOmI2Mo84o8EznPSDjCBEZR4jIOkmDrBOQdQKy0iDrBMoJpLAhvTQ57tdJyCGILELIJohMDSZTg8gmiCyCydYgcgjiJIEcJogcDSSHIHIIJIdA8nCRRyAFAcHgCgT3u7g/iysIcQUS4H4v/BzgCiTA5SLAFYTLFYi4XAS4AnEFuN/d24XnuQJcBAS6CAgoPObCJQHOPnHhcgUQEOC8iwQ45wQEEBDgfHYFBCAiBAi4AoQAEQICnO0AEcT9HuA+R0rsc96h5LmCswSJJceLy2cJQkSigIHAjQCqmgPkePkX3AtIUtVd7mvNBq4BLEEY3xBx1uUOiYSYlt4XU3VGjeecdlbeyznttJPknIacM5CbCblnICcDcjMJzj1DcG6mM/4jNwvNyyIvJ4OCnCwKcjIhLxvNPwN52Uh+DpKfQ0BBDlKQS4Dm4tJ858YK5LtfVVS+CopQgPMOggLq3j77cwFCvofjFO0rJMX2O8ekxLnOOYWK9pfys0cpvv/ccoXXlLP3e7ycN/fw5vyKORNYl44Pf1fBUuXz5RNEC5xqpDdE5FJgDfB797E7ReQGYDVwn6oeO6tsArC32HYq0NvTTUTkNuA2gKZNm1Ze9MZ4QwSCQp1XRIOKFwcqNDa8oAAK8iA/x3kV5Lm3cz18zndGrRfuL8gv2qcFeeTn55Ofl0d+QR4FefkUFORRkJ/vvArcr/w8VJWCgny0oAAtKKCgoAAtcParKlqQ734vcK7t3o8qqgWoFhR9dn6iF6AUHndeUmwbfi7vcB8D9z51/wT9uYeac1wRLfxhXkDJ/jda7KduQeEeRAuPnvUjudi9pfhmWT+6S+3wU8r+ip5fhrygyAqX8YYvE0Qg0B24S1V/EJHpwDTgBeAJnD+FJ4DncNomivOUVj3+qanqy8DL4PRiqpzQjamiAgKcqqXACxvTITj/Qa2XiimLL7tYpAKpqvqDe/sDoLuqHlDVfFUtAF7BqU7yVLZJse1EIN2HsRpjjDmLzxKEqu4H9opIW/euocBmEWlc7LRxwEYPxVcBrUWkubtxexLwka9iNcYYcy5fP2HeBbzr/iG/C7gJ+IeIdMWpMtoD3A4gIvE43VlHqWqeiNwJLMXp5vq6qm7ycazGGGOKsZHUxhhTi5U1ktqGeRpjjPHIEoQxxhiPLEEYY4zxyBKEMcYYj2pUI7WIHAKSz7N4A+BwuWf5T1WPDyzGylDV44OqH2NVjw+qVoyXqGqspwM1KkFcCBFZXVpLflVQ1eMDi7EyVPX4oOrHWNXjg+oRI1gVkzHGmFJYgjDGGOORJYifvezvAMpR1eMDi7EyVPX4oOrHWNXjg+oRo7VBGGOM8cyeIIwxxnhkCcIYY4xHtT5BiMgIEdkmIkkiMs3f8ZxNRJqIyHIR2SIim0Tk9+WXuvhExCUi60Rkob9j8URE6onIByKy1f1nebm/YzqbiNzr/jveKCKzRCTUz/G8LiIHRWRjsX3RIvKJiOxwv9evgjH+1f33/JOIzBeRelUtxmLH/p+IqIhUfDnCi6BWJwgRcQEzgJFAB2CyiHTwb1TnyMNZlrU90Ae4owrGCM5yslv8HUQZpgMfq2o74FKqWKwikgDcDfRU1U4409xP8m9UzARGnLVvGvCZqrYGPnNv+9NMzo3xE6CTqnYBtgMPXuygzjKTc2NERJoAw4CUix2Qt2p1gsBZzS5JVXepag4wG7jGzzGVoKr7VHWt+/MpnB9sCf6NqiQRSQRGA6/6OxZPRCQKGAi8BqCqOap63L9ReRQIhIlIIBCOn1dRVNWvgKNn7b4GeNP9+U3gFxc1qLN4ilFVl6lqnntzBc6KlH5Typ8jwN+B+zmfRagvktqeIBKAvcW2U6liP3yLE5FmQDfgh7LPvOiex/mHXlDeiX7SAjgEvOGuBntVRCL8HVRxqpoGPIvz2+Q+4ISqLvNvVB41VNV94PzyAsT5OZ7yTAWW+DuIs4nIWCBNVX/0dyxlqe0JQjzsq5LZXETqAB8C96jqSX/HU0hErgYOquoaf8dShkCgO/BPVe0GnMH/VSMluOvyrwGaA/FAhIj8xr9RVW8i8jBOFe27/o6lOBEJBx4G/uzvWMpT2xNEKtCk2HYifn6s90REgnCSw7uqOs/f8ZylHzBWRPbgVNFdISLv+Dekc6QCqapa+OT1AU7CqEquBHar6iFVzQXmAX39HJMnBwrXlXe/H/RzPB6JyBTgauA6rXqDvVri/CLwo/v/TSKwVkQa+TUqD2p7glgFtBaR5u51sycBH/k5phJERHDqzreo6t/8Hc/ZVPVBVU1U1WY4f36fq2qV+s1XVfcDe0WkrXvXUGCzH0PyJAXoIyLh7r/zoVSxhnS3j4Ap7s9TgAV+jMUjERkBPACMVdUMf8dzNlXdoKpx/7+9ewmt6orCOP7/VLADX2BFHChKxUc7UAQVSRQFnVgREcWBCEInHVTpqDOlTiSgIkIEBw4ctD5CxIBO+oCqSAQjvoIKlUoLtrTVYSSGapaDvS893hw1xiT3xnw/uNx9T+45ewdyzso+j7UiYnbebx4BS/LfaV0Z1QEiX8j6CviBtDO2RMTd2o6qjwZgB+k/81v5tb7WgxqBdgHfS7oDLAb213g8r8izm1bgBtBJ2jdrmo5B0ingKjBf0iNJXwBNwDpJD0h34DTV4RibgYnAT3l/OVaHYxwRnGrDzMxKjeoZhJmZvZ4DhJmZlXKAMDOzUg4QZmZWygHCzMxKOUDYqCJpauF24b8l/Vn43D6I/WyStDe3T0jaMljb7kffGyTtG67+7MPl21xt1JL0LdAVEQeHYNvtpAe1nkg6AVyIiNZB7mNsRLwoWS7S8xQN9figmI0cnkGYZZK68vtqSZcktUj6VVKTpO2SrknqlPRJ/t40SWcldeRXQ14+D+iJiCeFza+S1C7pYWU2oeRArv/QKWlbof8LhXE1S9qZ279L2ivpCrBV0m5J93Ltg9MAObXERVKqCbMBG1frAZjVqUXAQlKa5ofA8YhYplSwaRfwNanGxOGIuCJpFumJ/IWkp99vVG1vBtAILCClq2gFNpOe6l4EfAx0SLrcj7E9i4hGAEl/AXMioqeqMM51YCXQ8s6/uVnmAGFWrqOS1lrSb0Al9XYnsCa31wKfpjM6AEySNJEUDB5Xba8tInqBe5Km52WNwKl8mugfSZeApcDbsvWeKbTvkFKItAFtheX/krLCmg2YA4RZuZ5Cu7fwuZf/95sxwIqI6C6uKKkbmPyG7anqvdpzXj39W1169Gmh/TmpGNJGYI+kz3KOsY+Abszeg69BmA3cj6RkjwBIWpyb94G5/Vj/uYPHPAAAAMlJREFUMrBNqZ73NNKB/hrwB2lmMl7SZFJm1z4kjQFmRsQvpIJNU4AJ+cfzgD41kM3ehWcQZgO3GziaM8SOIx3wv8zvhyTpLbUIzgErgNukQlXfVFI+S2ohnT56ANx8zfpjge9yEBHpekillOoaal+L2UY43+ZqNgQkHQHOR8TPNeh7OnAyIkpnHmb95QBhNgTyQXp5RAx7ASpJS4H/IuLWcPdtHxYHCDMzK+WL1GZmVsoBwszMSjlAmJlZKQcIMzMr5QBhZmalXgIX4B3ab6ba1gAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"import numpy as np\n",
"from numpy import exp\n",
"import matplotlib.pyplot as plt\n",
"\n",
"T0 = 85; T = 74; Ta = 65\n",
"\n",
"t = np.linspace(0,15,50)\n",
"temp_analytical = Ta + (T0-Ta)*exp(-K*t)\n",
"\n",
"temp_numerical=np.zeros(len(t))\n",
"temp_numerical[0]= T0\n",
"for i in range(1,len(t)):\n",
" temp_numerical[i] = temp_numerical[i-1] + -K*(temp_numerical[i-1]-Ta)*(t[i]-t[i-1])\n",
" \n",
"plt.plot(t,temp_analytical, label='Analytical')\n",
"plt.plot(t,temp_numerical, label='Linear Numerical')\n",
"plt.xlabel('Time(hours)')\n",
"plt.ylabel('Temperature (degF)')\n",
"plt.legend(loc='best', prop={'size': 15});\n",
"print('a)')"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"b) The final temperature as time goes to infinity is 65 degF, or the ambient temperature.\n",
"c) Time: 0 hr 50 min 56 s\n",
" Time of death: 10:09:04 am\n"
]
}
],
"source": [
"import math\n",
"temp_death = 98.6\n",
"time = np.log((temp_death - Ta)/(T0-Ta))/K*3600\n",
"\n",
"hours = math.floor(time/3600)\n",
"minutes = math.floor((time/3600-hours)*60)\n",
"seconds = int(round(((time/3600-hours)*60-minutes)*60,0))\n",
"print('b) The final temperature as time goes to infinity is 65 degF, or the ambient temperature.')\n",
"print('c) Time:',hours,'hr',minutes,'min',seconds,'s')\n",
"print(\" Time of death: 10:09:04 am\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"4. Now that we have a working numerical model, we can look at the results if the\n",
"ambient temperature is not constant i.e. T_a=f(t). We can use the weather to improve our estimate for time of death. Consider the following Temperature for the day in question. \n",
"\n",
"|time| Temp ($^o$F)|\n",
"|---|---|\n",
" |8am|55|\n",
" |9am|58|\n",
" |10am|60|\n",
" |11am|65|\n",
" |noon|66|\n",
" |1pm|67|\n",
"\n",
" a. Create a function that returns the current temperature based upon the time (0 hours=11am, 65$^{o}$F) \n",
" *Plot the function $T_a$ vs time. Does it look correct? Is there a better way to get $T_a(t)$?\n",
"\n",
" b. Modify the Euler approximation solution to account for changes in temperature at each hour. \n",
" Compare the new nonlinear Euler approximation to the linear analytical model. \n",
" At what time was the corpse 98.6$^{o}$F? i.e. what was the time of death? \n",
" \n",
" "
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"def current_temp(time):\n",
" t = np.array([-3,-2,-1,0,1,2])\n",
" temp = np.array([55,58,60,65,66,67])\n",
" index = np.where(t==time)[0][0]\n",
" return temp[index]"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAESCAYAAAAWtRmOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deXhV1fX/8fdiCPMcRhkDAbTMIiIGZ61VK3W2WsWhoq1SbbW2amut1Tpb0fL9VautM87YVhS1jiAKZVYEJDLPhCGEIYQk6/fHOWgMCeSEm5zk5vN6nvvc3DPddRnuyj5777XN3RERESmqVtwBiIhI1aPkICIie1FyEBGRvSg5iIjIXpQcRERkL0oOIiKyFyUHERHZi5KDiIjspU55TjKzjsDBQGN3H5/YkMonNTXVu3btGncYIiLVyowZM7LcvXXx7ZGSg5n1Bh4FMsJNvucaZnY08DJwsbtPPLBwo+vatSvTp0+v7LcVEanWzGxZSdvLfFvJzLoBk4EhwBPAe4AVOWQyUACcU/4wRUSkKojS5/AHoDEwzN1HAZ8U3enuBeG2oYkLT0RE4hAlOZwEjHf3Wfs4ZgXQ/sBCEhGRuEVJDq2AEu9NFeFA/fKHIyIiVUGU5LAB6LafYw4GVpU/HBERqQqiJIePgNPDjum9mFk/4ETg3UQEJiIi8YmSHO4KnyeZ2UVAKoCZdQlfTwB2Ag8kNkQREalsZZ7n4O5fmNkFwNPAk+FmAxaHP+8Aznf3rxMaoYiI7GVXfgEzlm5mUmYWlwzrStumie3ujTQJzt3Hm1kP4KcEQ1ZbAdnAZ8Bj7r4modGJiAgA7s7CdTlMXpTFx4uymLZkI7m7C6lTyxjcpUW8ySEMcB1wZ0KjEBGRvazfmsvkzCwmL8piUmYWG3J2AZDWuhHnH9aZjB6pDO3eisb1ylUJaZ/KfEUzywZec/dLEx6FiIiwM6+AqUs2BslgURYL1+UA0KJhXTLSWzO8RyoZ6al0aN6gwmOJmm7WVkgUIiI1UGGhM2/1ViZlbmDyoiymL91MXkEhKbVrcVi3FvxoYG+Gp6dySPum1Kpl+79gAkVJDnOAXhUViIhITbBqy04mL9rApEVZfJKZxeYduwHo3a4JI4d1ISO9NUO6tqRBSu1Y44ySHB4AXjazDHefXFEBiYgkk5zc3Xy2eNM3CWFx1nYAWjepx7G92zA8PZUje6TSpknVKi4RJTnUBiYC75vZ88D/CG4zefED3f21xIQnIlK95BcUMmdlNpMXZTE5cwOzlm8hv9CpX7cWh3drxQWHd2Z4emt6tm2MWeXeKooiSnJ4hSARGHBx+CieGCzcFm97SESkkrg7yzbuYFJmFpMXbWDK1xvJyc3HDPp0aMaoo9LISE/l0C4tqFen+nw1RkkOoyssChGRamTLjjymfL2RSWHrYMWmnQAc1LwBp/ZtT0Z6KsO6p9KyUUrMkZZflBnSYysyEBGRqiovv5CZyzd/M9/g85VbKHRoXK8OQ9NaccXwNDJ6pNIttVGVvlUUReJnToiIVHPuTub6bWHLIIvPFm9kR14BtWsZ/Ts2Y/Rx6QxPT6V/p+bUrR2lRF31oeQgIgJkbdvFJ5nB5LPJi7JYuzUXgK6tGnLWoI5kpKdyRPdWNK1fN+ZIK0eUGdJzy3iou3v/csYjIlIpcncX8L+lm76Zjfzlmq0ANGtQl4xwJnJGj1Q6tWwYc6TxiNJy6EAJw1YJ1pWuF/68BSiIEoCZNQceB/qE17/M3T81s9HANUA+MMHdb4xyXRGRogoLnflrt4ZDTLOYtmQTu/ILqVvbOLRLC379/V5k9Eilz0HNqF3Js5Groigd0qklbTezWkA/4C/ALuCHEWMYA0x097PNLAVoaGbHAiOAfu6+y8zaRLymiAhrs3OZtGgDkzOD2chZ2/IASG/TmAsP78Lw9FSGdGtJowooXFfdHfCfiLsXArPN7IfAF8AtwG1lOdfMmgJHAZeE18oD8szsZ8Dd7r4r3L7+QOMUkeS3fVc+U5ds/KbfYNH6bQCkNk4JbxW1JqNHKu2aVa3ZyFVRwtKlu28zs7cJJsfdVsbT0gjWpv6nmfUHZgDXAj2B4WZ2J5AL3ODu/0tUrCKSHAoKnc9XZX9TmmLm8s3sLnDq1anFkG4tOWdwRzJ6tKZ3uyaVXriuukt0WyqXoG8iyvsPAka7+1QzGwP8NtzegmBBocOAl8wszd2/0+dhZqOAUQCdO3dOQPgiUtWt2LTjm8lnn2RuJHtnULjukPZNuezIbgxPb83gri2oX7f6zEauihKWHMysGUE/weoIp60EVrr71PD1KwTJYSXB2hEOTDOzQoI1qzcUPdndHwMeAxg8eHBJneUiUs1l79zNp19vZHJY1nrpxh0AtGtan5MOaUtGWLgutXG9/VxJoogylPVX+7hGJ+BsoA1lv6WEu681sxVm1svdFwLHA18CXwPHAR+aWU8gBcgq63VFpPraXVDI7BVbwn6DDcxZmU1BodMwpTZD01oxclhXhqen0r111S5cV91FaTncz7eF90qSCzwE3BExhtHAc+FIpcXApcB24B9m9gWQB4wsfktJRJKDu7MkazuTwvkGny3eyLZd+dQy6NuxOT8/pjsZPVIZ2LkFKXWSczZyVRQlOZQ2RLUQ2Ax87u7bowbg7rOBwSXs+knUa4lI9bBpex6fhGsjT87MYtWWoHBdp5YNOH1AB4b3CArXNWtYM2YjV0VR5jlMqMhARCR57covYMbSzWFZ6yy+WJ2NOzSpX4cju6fys2O6Mzw9lS6tGsUdqoSi9jlMcffP9nHM4cCR7v5gIoITkerJ3flq3TYmhUNMpy3ZxM7dBdSpZQzq3IJfntCTjPRU+h3UjDpJWriuuova53AbUGpyAE4AbgeUHERqmPU5ucFtovBW0fqcXQCktW7EeYd1IqNHKkO7t6KxZiNXC4n+W6pD0AchIkluZ14B05ZuYtJXQXmKBWtzAGjRsC4Z6a0ZHhav69C8QcyRSnkkOjn0BTYl+JoiUgUUFjpfrtnKx4uC+QbTl24mr6CQlNq1OKxbC35zcm+Gp6dySPummo2cBPaZHMzs38U2XWBmJY0sqg10Bg4hmMgmIklg1Zad35SmmPL1RjZtDwrX9W7XhJHDupCR3pohXVvSIEWzkZPN/loOpxX52QlqHvUs5didwATgugTEJSIxyMndzWeLNwUJITOLxRuC0eltmtTjmF6tGR7ORm7TRIXrkt3+kkOT8NmArcCfgbtKOK7A3XMTGZiIVLz8gkLmrMwOO5E3MGv5FvILnfp1azE0rRUXDOnM8PTW9Gyr2cg1zT6TQ9FJbeHiO1PLM9FNRKqOZRu383FYmmLK1xvJyc3HDPoe1IxRR6WRkZ7KoV1aUK+ObhXVZFEmwY2tyEBEpGJs2ZHHlK83flPJdMWmYDbyQc0bcGrf9kHhuu6ptGiUEnOkUpWUa7SSmbUgKM1dYhlEd595IEGJyIFxd56cspTXZ6/m85VbKHRoXK8OR3RvxRXD08jokUq31Ea6VSSlipQczCwDeICSayEVpfaoSIxenrGSP/7nS/p1bMbo49IZnp5K/07NqavZyFJGUcpnDAT+C+QATxJUT/0UWAIMA7oCbwLzEx2kiJTdonU5/OFf8xjWvRXPXH44tTXnQMohyq8RtwAFwBB3vzzc9ra7/4RgeOuDwJGEi++ISOXL3V3ANc/PomFKbR46b4ASg5RblOSQAfzb3ZcU2WYA7p4P/BpYBvwpceGJSBR//M+XLFyXw4PnDaBNU81FkPKLkhxaENxC2mM38E193XAxno+AYxMTmohE8cbc1Yybtpyrju7O0T1bxx2OVHNRkkMW0KzI6/VAtxKup4LsIpVs+cYd3PTq5wzq3JzrTyqtiIFI2UVJDouAtCKv/wecaGZdAMysFXAmwfrPIlJJ8vILuWbcTMzg4R8P1IgkSYgo/4omAseY2Z7WwyME5TVmm9kHBKOU2gF/TWyIIrIv90xcwNyV2dx7dj86tmgYdziSJKIkh8cICvHt6YT+ABgJZANHA7uAX7v73xMdpIiU7L3563hi8hJGHtGFk/u0jzscSSJRymdsAt4rtu1Z4Fkzq+3uBYkOTkRKtyZ7J9e/PIdD2jflplMOjjscSTJlbjmY2cNm9rOS9ikxiFSu/IJCrh03m7z8Qv56wUDq11VRAkmsKLeVrgS6VFQgIlJ2D7+3iGlLN3HnGX1Ia9047nAkCUVJDsuBVhUViIiUzSeZWTzyQSZnH9qRMwZ2jDscSVJRksOLwPfNrMl+j4zAzJqb2StmtsDM5pvZEUX23WBmbmapiXxPkepqQ84urntxNmmpjbh9xPfiDkeSWJTkcAfwFfCumR1jZoma7DYGmOjuvYH+hIX7zKwTcCJBi0WkxissdH710my27tzN2AsH0TClXBX3Rcokyr+u9QTJpCHhqCUz20GwtnRR7u7NKAMzawocBVwSnpgH5IW7/wLcCPwrQowiSevRjxczaVEWd57Rh97tmsYdjiS5KMnhK/ZOBAcqDdgA/NPM+gMzgGuB44FV7j5Hi5GIwIxlm7j/nYWc2rc9FwzpHHc4UgNEmeewvwV+yvv+g4DR7j7VzMYAtxG0Jk7a38lmNgoYBdC5s/7DSHLasiOPX4ybTYfm9bnrrL5avU0qRdxFWFYCK919avj6FYJk0Q2YY2ZLgY7ATDNrV/xkd3/M3Qe7++DWrVWFUpKPu3PjK3NZn5PLX388iKb168YdktQQ5UoOZlbXzNLD1eHKzd3XAivMrFe46Xhgpru3cfeu7t6VIIEMCo8VqVGe/nQZ73y5jt+c3Jv+nZrHHY7UIJGSg5m1NrOngC3AAoLKrHv2DTOzaWY2LGIMo4HnzGwuMAD4c8TzRZLSF6uyuXPCfI7r3YbLM4pXxxepWFHWkG4NTCVYK/q/QHPg0CKHzADSgfOBKWW9rrvPBkrtzwhbDyI1yrZd+YweN4uWjVK4/5z+6meQShel5XArQfmMH7r7ScCEojvdfRfwMUGFVhEpJ3fnd+M/Z9nG7Yw5fwAtG6XEHZLUQFGSww+B/7j7hH0csxzocGAhidRsL89YyeuzV3PdCT05PE0VayQeUZJDe4J+hn3ZSbAAkIiUQ+b6HP7wr3kM696Kq4/tEXc4UoNFSQ6b2X+rIB1YV/5wRGqu3N0FXP3cLBqm1Oah8wZQu5b6GSQ+UZLDZ8Dp4VrRezGzrsAPgI8OPCyRmuf2N75k4bocHji3P22a1o87HKnhoiSH+wluGb1vZsMJRzpZYDhBDaRawIMJj1Ikyb0xdzXPT13OlUencUyvNnGHIxKpfMZkM7uWoCDeh0V27SZYV7oQ+Hk4NFVEymj5xh3c9OrnDOzcnBtO6rX/E0QqQaSav+7+VzObBFwDDCVY/Ceb4JbTGCUGkWjy8gu5ZtxMzOCRHw+kbu24K9qIBCIXhHf3OcAVFRCLSI1z78QFzF2Zzd9+MoiOLRrGHY7IN/RrikhM3pu/jscnL+HiI7pwcp/2cYcj8h2RWw5hsb1LgIFAM4LbSrOAJ919VkKjE0lSa7J3cv3LczikfVNuPuXguMMR2Uuk5GBmdwG/Zu8WRwZwtZnd6+43Jyo4kWSUX1DIteNmk5dfyF8vGEj9urXjDklkL2W+rRQurPMbYA1BJdW+QOvw+RfAWuA3ZvbTCohTJGk8/N4ipi3dxJ1n9CGtdeO4wxEpUZSWwzUE60gPcvcNRbZvBOaZ2cvAXIJE8XjiQhRJHlMys3jkg0zOPrQjZwzsGHc4IqWK0iGdDrxcLDF8w93XAS+Hx4lIMVnbdnHti7NJS23E7SO+F3c4IvsUpeWwGdixn2N2AJvKH45IciosdH754myyd+7m6cuG0DAl8lgQkUoVpeXwBnCamZXYe2ZmdYDTwuNEpIhHP17MpEVZ/OGHh3Bw+6ZxhyOyX1GSw2+BAuANM+tXdIeZ9SdICvnhcSISmrFsM/e/s5BT+7bngiGd4w5HpEyitG0/BBoAJwEnmVkOQXnutny7hsNi4KNiSxq6u/c/8FBFqp/sHbv5xbhZdGhen7vO6qvlPqXaiJIcOgDOd/sUWhIU3tuzrXn4EKnx3J0bX53Duq25vPKzYTStXzfukETKLEpV1tSKDEQk2Tz96TLenreOW045mAGd9DuTVC+qrSRSAb5Ylc2dE+ZzXO82XJ7RLe5wRCJTchBJsG278hk9bhYtG6Vw/zn9qaXlPqUaKk/hve5AP6AjUOJNVHfXanBSI7k7v3/9C5Zt3M64K4bSslFK3CGJlEuZk4OZNSEoi3EWwcpvpf065ERYKtTMmofX7ROeexlwJvBDIA/4GrjU3beU9ZoicXllxkrGz1rFL0/oyeFpJS63LlItRGk5PAScA0wBngdWEcxrOFBjgInufraZpQANgXeBm9w938zuAW4iKPonUmVlrs/h1n/N44i0VlxzXI+4wxE5IFGSwwjgf8Bwd/dEvLmZNQWOIlgfAnfPI2gtvFPksM+AsxPxfiIVJXd3AVc/N4uGKbV56PwB1FY/g1RzUTqk6wEfJSoxhNKADcA/zWyWmT1uZo2KHXMZ8FYC31Mk4W5/40sWrsvhgXP707Zp/bjDETlgUZLDbKBLgt+/DjAI+H/uPhDYTpHyG2Z2C8Gtq+dKOtnMRpnZdDObvmFDicViRSrchLlreH7qcq48Oo1jerWJOxyRhIiSHP4MnG5mhybw/VcCK919avj6FYJkgZmNJCjkd2FprRV3f8zdB7v74NatWycwLJGyWb5xB799dS4DOzfnhpN6xR2OSMJEmSH9VrjK21tm9iIwk2D96JKOfa2M11xrZivMrJe7LwSOB740s5MJOqCPdvf9lQkXiUVefiHXjJuJGTx8/kDq1ta0IUkeUYay1gKGEtROuppg2Oleh4XboyyKOxp4LhyptBi4lKDjux7wblio7DN3vyrCNUUq3L0TFzB3ZTZ/+8kgOrVsGHc4IgkVZbTSnwiSwlKCFd9Wk4ChrO4+GxhcbLPGAUqV9v6CdTw+eQkXH9GFk/u0jzsckYSLkhwuBuYDg919ZwXFI1LlrcneyfUvzeGQ9k25+ZSD4w5HpEJEuUnaCnhLiUFqsvyCQq4dN5td+YX89YKB1K8b5Q6qSPURpeWwENA4PanRHn5vEdOWbuLBc/uT1rpx3OGIVJgoLYcxwAgz61oxoYhUbVMys3jkg0zOPrQjZw7qGHc4IhUqSsthLkHNo6lmdh8wg9KHss5MQGwiVUbWtl1c++Js0lIbcfuI78UdjkiFi5IcphMMUzXgnv0cqxuxkjQKC51fvTSH7J27efqyITRMiVzpXqTaifKv/EFKntsgktQem7SYj7/awB0/6sPB7ZvGHY5IpYgyQ/qGigxEpCqasWwz9729kFP6tuPCwzvHHY5IpdF8f5FSZO/YzS/GzaJD8/rcdWY/wtn6IjVCeZYJrQUMAw4GGrv7X8LttYF6qoUkycDdufHVOazbmssrPxtGswYlrogrkrQitRzM7ChgEfAR8Chwf5HdRwA5ZnZO4sITicczny3j7Xnr+M3JvRnQqXnc4YhUujInBzPrD0wEmgC/Jyiv/Q13nwwsB85IZIAilW3e6mzueGM+x/ZqzeUZ3eIORyQWUVoOtwK7CWor/RmYV8Ix0wjXYxCpjrbtyuea52fRolFdHjh3ALW03KfUUFGSw1HAeHdfvo9jlgMqUSnVkrvz+9e/YNnG7Yw5fyAtG6XEHZJIbKIkhybAuv0cU49ydHKLVAWvzFjJ+FmruPb4ngxNaxV3OCKxipIc1hCMUNqXfgTrPYhUK5nrc7j1X/M4Iq0V1xyn5UREoiSHd4AflLaGtJkdS3DraUIiAhOpLLm7C7jm+Vk0TKnNQ+cPoLb6GUQiJYc/A9uBD8zs94SrtZnZ0eHr8UAW8EDCoxSpQLe/8SUL1ubwwLn9adu0ftzhiFQJUcpnLDOzU4CXgD/ybRG+98PnVcCZ7r6/fgmRKmPC3DU8P3U5Vx6dxjG9tFyJyB6ROo/dfYqZdQfOBoYSrA6XDXwGvOjuuYkPUaRiLN+4g9++OpeBnZtzw0m94g5HpErZZ3IIZ0QvLTp81d13Ac+FD5FqKS+/kNHjZmIGD58/kLq1VWZMpKj9/Y/4ALikEuIQqVT3vb2AOSuzuffsfnRq2TDucESqnP0lBw3bkKTz/oJ1/H3SEi4a2oWT+2jOpkhJ1JaWGmVN9k6uf2kOB7dvyi2n7m/ajkjNFXtyMLPmZvaKmS0ws/lmdoSZtTSzd81sUfjcIu44pfrLLyjk2hdmsyu/kLEXDKR+Xa1mK1KasiSHil4adAww0d17A/2B+cBvgffcPR14L3wtckAefj+TaUs2cceP+pDWunHc4YhUaWVJDreZWUGER35Z39zMmhLMqn4CwN3z3H0LMAJ4KjzsKeBHET+XyHdM+TqLR95fxFmDOnLmoI5xhyNS5ZVlnkPUTukox6cBG4B/hutFzACuBdq6+xoAd19jZpqdJOWSvWM3//dhJv+cspRuqY24fcT34g5JpFooS3K4zd1vr8D3HwSMdvepZjaGCLeQzGwUMAqgc2ct/i7fyt1dwFNTljL2g0xyduVz5sCO/Pr7vWhUT0WDRcoi7v8pK4GV7j41fP0KQXJYZ2btw1ZDe2B9SSe7+2PAYwCDBw+u6L4RqQYKCp3xs1bx4DsLWZ2dy7G9WnPjyb05uH3TuEMTqVZiTQ7uvtbMVphZL3dfCBwPfBk+RgJ3h8//ijFMqQbcnQ8XbuCeiQtYsDaH/h2b8cC5Aziiu9ZlECmPuFsOAKOB58wsBVgMXErQUf6SmV1OsLrcOTHGJ1Xc7BVbuOvN+UxdsomurRoy9oJBnNK3HWaawylSXrEnB3efDQwuYdfxlR2LVC9LsrZz39sLePPztaQ2TuFPI77H+UM6q06SSALsLzl0A7ZURiAiZbUhZxdj3vuKF6atIKVOLa47IZ2fDk+jsTqbRRJmn/+b3H1ZZQUisj/bduXz2MeLeXzSYvLyC/nxkM784vh0WjepF3doIklHv2pJlZeXX8gL/1vOw+8tImtbHqf2bc8N3+9Ft9RGcYcmkrSUHKTKcncmfL6G+95eyLKNOxia1pLHRx7MgE7N4w5NJOkpOUiVNOXrLO5+awFzV2bTu10T/nnpYRzTs7VGIIlUEiUHqVLmr9nK3W8t4KOvNtChWX0eOKc/Pxp4ELVrKSmIVCYlB6kSVm7ewYPvfMX42atoWr8uN5/Sm4uP6Kqy2iIxKXNyMLMCgjpLf9rHMbcAf3R3JR0pk83b8/i/DzN5asoyMBh1VBo/P7oHzRrWjTs0kRotype4UbaKq2r/y37l7i7gn58s5f8+zGT7rnzOGtSRX57Ykw7NG8QdmoiQ+NtKLYDcBF9TkkhBofPqjJU8+O5XrN2ay/G923Djyb3p1a5J3KGJSBH7TA5mdlSxTV1L2AZQG+gMXAgsTFBskkTcnffmr+eeiQtYtH4bAzo1Z8z5Azg8TYXxRKqi/bUcPuTbZUKdoELqyFKONaAQuD4hkUnSmLFsM/e8tYBpSzeRltqI/3fhIE7uo8J4IlXZ/pLD7QRJwYBbCZLFRyUcVwBsBD5w9wWJDFCqr683bOO+iQuZOG8tqY3rcceP+nDeYZ1UGE+kGthfbaXb9vxsZiOB19394YoOSqq39Vtzeei9Rbz4vxXUr1OLX53Yk8szumkVNpFqpMz/W929W0UGItVfTu7usDDeEvILC7loaBeuOa4HqY1VGE+kutGvcnLA8vILeW7qMh55P5NN2/P4Yf8O3HBST7q0UmE8keoqUnIws3TgWmAIwbDVkqavurt3T0BsUsUVFjr/mbua+99ZyIpNOxnWvRW//UFv+nVUYTyR6i7KDOkjgP8CDYB8YF34vNehiQlNqrLJi7K4e+J8vli1lYPbN+Wpy/pyVHqqRiCJJIkoLYe7gHrAVcA/3L2kxCBJ7otV2dwzcQGTFmVxUPMG/OW8/ozofxC1VBhPJKlESQ6HAa+4+2MVFYxUXSs27eCBdxby+uzVNG9Yl9+dejAXHdGFenVUGE8kGUVJDnnA8ooKRKqmTdvz+Ov7mTz72TJq1YKfH9Odq47pTtP6KownksyiJIcpwMCKCkSqlp15BfzjkyX87cOv2Z6Xz7mDO3HdCT1p16x+3KGJSCWIkhxuBqaY2UXu/kxFBSTxyi8o5OUZK/nLu1+xPmcXJx7Slhu/34v0tiqMJ1KTREkOI4D3gSfN7KfADGBLCcf5vtZ8kKrJ3Xnny3XcO3EBX2/YzqFdWjD2wkEc1rVl3KGJSAyiJIfbivw8PHyUxIEyJwczWwrkENRnynf3wWY2APgbUJ9guOzP3X1ahFglgulLN3HXWwuYsWwz3Vs34tGLDuWkQ9pqWKpIDRYlORxbYVHAse6eVeT1vQQryr1lZqeEr4+pwPevkTLX53DPxIW8++U62jSpx11n9uWcQztSR4XxRGq8KLWVSqrGWlEcaBr+3AxYXYnvnfTWb83lwXe/4qXpK2iYUocbTurJZRndaJiiaioiEqgK3wYOvGNmDjwazqO4DnjbzO4HagHDSjrRzEYBowA6d+5cSeFWb58t3sjPn5tJTu5uRg7ryujj0mnZKCXusESkiomcHMysH3ABcDDQyN1PCLd3Jai59K67b45wySPdfbWZtQHeNbMFwNnAL939VTM7F3gCOKH4iWEieQxg8ODBXny/fMvdeXbqcv7473l0btmQl64cSo82GoEkIiWLWnjvdoIhrXtuShf9Qq4FjCP4rf+Rsl7T3VeHz+vNbDxBghlJUOAP4GXg8Shxynfl5Rfyh39/wbhpKzi2V2seOn8gzRpoEpuIlK7MPY9mdj7wO+BdYABBraVvuPtiYDpweoRrNjKzJnt+Bk4CviDoYzg6POw4YFFZrynftSFnFxf8/TPGTVvBz4/pzuMjD1NiEJH9itJy+AWQCYxw9zwzO6OEY+YTbVRRW2B8OGSyDvC8u080s23AGDOrA+QS9itINHNXbuHKZ2aweUceD/94IKf37xB3SCJSTURJDn2BJ909bx/HrCb4wi+TsLXRv4Ttk4FDI8Qmxbw+axW/eXUuqY3r8cpVw+hzULO4QxKRaiRKcjCgcD/HtK0yNN0AABJGSURBVCX4TV9iUlDo3DNxAY99vJghXVvyfz8ZpGU6RSSyKMlhEaUMKQUws9pABjDvQIOS8snesZvRL8zi46828JOhnbn1tO+RUkcT2kQkuijfHC8Bg8zs+lL23wT0AJ4/4Kgkssz1OYwYO5lPv87iz2f05Y4f9VViEJFyi9JyeAg4B7g3nHvgAOFEteHAYOAzwnkHUnn+++U6rntxNvXr1uL5K4aqWJ6IHLAo5TN2mtmxwBjgQmDPEmC/IuiLeBa4RsuHVh53Z+wHmTzw7ld8r0NTHrtoMB2aN4g7LBFJApEmwbl7NnCJmf2KYNnQVkA2MM3dN1RAfFKKHXn5/PrluUz4fA0jBnTg7jP70SBFS3aKSGKUq7aSu28C3k5wLFJGKzbt4Iqnp7NwXQ43/aA3o45KU3ltEUmoqlB4TyKY8nUWVz83k/xC55+XHMYxvdrEHZKIJKFSk4OZ/YOg0/lmd18Xvi4Ld/fLExKdfMPdefrTZdz+xpd0bdWQv188mLTWjeMOS0SS1L5aDpcQJId7gHXh67JwQMkhgXblF3Dr6/N4cfoKju/dhofOH0CT+qqPJCIVZ1/JoVv4vKrYa6lE67fmctWzM5i5fAvXHNuDX53Yk1q11L8gIhWr1OTg7sv29Voq3pwVQeG87J27GXvBIE7t1z7ukESkhlCHdBX12syV/Pa1z2nduB6v/mwYh3Rouv+TREQSpDwrwWUAlwIDCdZ3zgZmElRsnZzY8Gqe/IJC7n5rAY9PXsLQtJaMvWAQrVQ4T0QqWdSV4B4Bfk5QobWoAcClZjbW3X+RqOBqmi078hg9bhaTFmUx8ogu/O60Q6hbW/WRRKTyRVkJbjRwNbCEoOXQDWgQPl8Wbr/azK6ugDiT3lfrchgx9hM+W7yRu8/syx9H9FFiEJHYRGk5XEWwmM9gd99SZPsy4Ekz+zfwOUHLYmziQkx+78xbyy9fnE2DlDq8MGooh3ZR4TwRiVeUX03TgFeLJYZvhCU1Xg2PkzIoLHTG/HcRo56ZQfc2jfnP6COVGESkSojSctgI7GuJUML9WeUPp+bYviuf61+aw8R5azlj4EHcdWZf6tdV4TwRqRqiJIfXgdPN7GZ33118p5mlAKeHx8k+LN+4g1HPTOerdTn87tSDuTyjmwrniUiVEuW20s0Ew1b/a2bDLPw2s8CRwH+BzeFxUopPMrM4fexk1mTn8tRlQ/jpcFVUFZGqZ1+F9xaXsDkFaA9MAvLNLAtILXKdNcBsoHuC46z23J0npyzljgnzSUttxN8vHkzX1EZxhyUiUqJ93VaqRbgUaBG7geXFtq0u4TwpYld+Ab8b/wUvz1jJiYe05S/nDaBxPU1OF5Gqa1+1lbpWYhxJa93WXK58ZgazV2zhF8enc93x6SqcJyJVXuy/vprZUiAHKADy3X1wuH00cA2QD0xw9xtjC7KcZi3fzJXPzGDbrnz+34WD+EFfFc4Tkeoh9uQQOtbdvxkCa2bHAiOAfu6+y8yq3XJnr8xYyc2vfU7bZvV4+vJh9G6nwnkiUn2Up/BeP6A/0BEoacUZd/c/HWBcPwPudvdd4QXXH+D1Kk1+QSF/fnMB//hkCcO6t2LsBYNo0Sgl7rBERCIpc3Iws5bAM8DJezaVcqgDUZKDA++YmQOPuvtjQE9guJndCeQCN7j7/0qIaRQwCqBz584R3rJibN6exzXjZvJJ5kYuPbIrt5xyMHVUH0lEqqEoLYeHgB8QzGd4lmCFuPwExHCku68Obx29a2YLwrhaAEOBw4CXzCzN3b8zeipMJI8BDB48uPjIqkq1cG0OVzw9nbXZudx7dj/OHdwpznBERA5IlORwGjDF3U9KZADuvjp8Xm9m44EhwErgtTAZTDOzQoL5FBsS+d6JMvGLNfzqpTk0rleHF64cyqDOLeIOSUTkgES551EbmJLINzezRmbWZM/PwEnAFwQlOI4Lt/ckmHxX5Wo2FRY6f3n3K656dibpbZvwn9EZSgwikhSitBxmkviKq22B8WH5iDrA8+4+MazT9A8z+4KgmN/I4reU4rZtVz6/enE273y5jrMGdeTOM/qocJ6IJI0oyeFPwJtmlpGo5UDdfTHByKfi2/OAnyTiPSrCso3bueLp6Xy9YTu3nnYIlx7ZVfWRRCSplDk5uPv7ZnY+wW/6bxC0JLJLOfbpBMVX5UxelMXVz88E4KlLh5CRnhpzRCIiiRdlKGsKwcS0FsDI8FH8Vo+F25IuObg7T0xewp/fnE96myY8dvGhdGmlwnkikpyi3Fa6iyAhfAm8SFBwLxFDWau83N0F3Dz+c16buYrvf68tD547gEYqnCciSSzKN9z5BGtEHxb2CdQIa7NzufLZGcxZsYVfntCT0cf1UOE8EUl6UZJDc4LRRDUmMcxYtpmrnp3Bjl35PHrRoXz/e+3iDklEpFJESQ7zCRb6qRFemr6C343/gnbN6vPs5YfTq12TuEMSEak0UZLDA8Dfzaynu39VUQHFbXdBIXdOmM+TU5aS0SOVv14wkOYNVThPRGqWKMlhFTARmGpmY4AZlD6U9eMExFbpNm3P4+rnZvLp4o1cntGNm37QW4XzRKRGipIcPiQYpmrArew9jLWoajdVeP6arVzx9HTW5+zigXP6c9ahHeMOSUQkNlGSw+3sOyFUW29+vobrX5pD0wZ1eOnKIxjQqXncIYmIxCrKDOnbKjCO2Iz9IJP73l7IwM7NefQnh9Kmaf24QxIRiV1Cb6ibWS0zG5HIa1a0rq0acd7gTrwwaqgSg4hIKCHTfM2sC/BT4FKC4a7Vps/h1H7tObVfjRmhKyJSJuVODmZWm6DW0ijgBIJWiBOsFCciItVY5ORgZmkErYRLCNZjgGAhnkeBJ9x9WcKiExGRWJQpOZhZHeAMglbCsQSthDzgNeAs4F/ufmtFBSkiIpVrn8nBzNKBKwiqsaYSzHGYCTxJUGdpU7i+s4iIJJH9tRwWEvQjrAf+AvzT3edVeFQiIhKrsgxldeBN4BUlBhGRmmF/yeH3wDKCIaqfmNmXZnajmWnsp4hIEttncnD3O929O/ADYDzQHbgbWG5mE8zs3EqIUUREKpm5l71ckpm1AS4jGMqaxre1lmYAP3P3GQmPsOyxbSBo5ZRHKsFw3JpEn7lm0GeuGQ7kM3dx99bFN0ZKDt850ex4gqGtI4AUgkQxF3jc3ceWM8hYmNl0dx8cdxyVSZ+5ZtBnrhkq4jOXu7aSu7/n7ucBHYEbga+A/sDDCYpNRERicsCF99w9y93vd/eDgeOAcQceloiIxCkhhff2cPcPCRYFqm4eizuAGOgz1wz6zDVDwj9zufscREQkeWmBZBER2YuSQ8jM/mRmc81stpm9Y2Yd4o6popnZfWa2IPzc480s6ddHNbNzzGyemRWaWdKOaDGzk81soZllmtlv446nMpjZP8xsvZl9EXcslcHMOpnZB2Y2P/w3fW0ir6/k8K373L2fuw8A3gBqQpXZd4E+7t6PYLTZTTHHUxm+AM4EPo47kIoSrrUylmDy6iHAj83skHijqhRPAifHHUQlygeuDwcDDQWuTuTfs5JDyN23FnnZiG8n+CUtd3/H3fPDl58RDEtOau4+390Xxh1HBRsCZLr7YnfPA14gmI+U1Nz9Y2BT3HFUFndf4+4zw59zgPnAQYm6fkJHK1V3ZnYncDGQTbBuRU1yGfBi3EFIQhwErCjyeiVweEyxSCUws67AQGBqoq5Zo5KDmf0XaFfCrlvc/V/ufgtwi5ndBFwD/KFSA6wA+/vM4TG3EDRRn6vM2CpKWT5zkrMStiV9S7imMrPGwKvAdcXugByQGpUc3P2EMh76PDCBJEgO+/vMZjYSOA043pNkXHOEv+dktRLoVOR1R2B1TLFIBTKzugSJ4Tl3fy2R11afQyhc9W6P04EFccVSWczsZOA3wOnuviPueCRh/gekm1k3M0sBzgf+HXNMkmBmZsATwHx3fzDh10+SXxYPmJm9CvQCCgmqu17l7qvijapimVkmUA/YGG76zN2vijGkCmdmZwCPAK2BLcBsd/9+vFElnpmdAjwE1Ab+4e53xhxShTOzccAxBBVK1wF/cPcnYg2qAplZBjAJ+JzgewvgZnd/MyHXV3IQEZHidFtJRET2ouQgIiJ7UXIQEZG9KDmIiMhelBxERGQvSg4iIrIXJQepFszsQzOrsuOuw/g+N7NaRbZdYmZuZpfEGFqFMbOGZrbWzJ6JOxZJPCUHqVThl2WUxyVxx7w/ZnY2cDTBpKvC/R2fLMJZ9XcBF5rZkLjjkcSqUbWVpEr4YwnbrgOaAWMIZi0XNTt8vhhoWIFxlUtYwuAOgvUwxsccThweJahBdgdwUsyxSAIpOUilcvfbim8LWwfNgIfcfWkp5y2v0MDK7wSCsiu3JEvhwijcPdfMXgSuNLN0d18Ud0ySGLqtJNVCSX0OZnZMeOvpNjMbbGYTzSzbzDab2atm1ik8Ls3MXjCzDWa2M1xasX8p79PQzG4Kl4vdbmbbzOxTM/txKaFdHj7vcy0MMzs2/Aw5ZrbVzCaY2cGlHNvezMaa2VIzywvjfs3MDi3h2NvCP4NjStjXNdz3ZLHtT4bb08xsdLhM7E4z+zDcb2Y20symhO+da2YrzOxtMzuvhJBfICgTftm+/gykelFykGRwGEEBMoC/A9MIlgJ9z8x6h687Ak8TlGI/Gng3rIP/DQvW0J4M/BkoAP4BPEVQpO95M7uj2PEGHAesdfev9xHfacA7wFbgb2GspwAfmVlqsWt2A6YDPwe+Bh4A3gZOBaaY2Wll+yMpkzHAnwgKt40BPgm330mw5GY74CXgQeC/BIsInVPCdaYBu4ETExibxM3d9dAj1gewlGAxmq77OObD4J/rd7YdE57nwIXF9j0Rbt9EcMun6L7fh/uuLbb9yXD7jcW21wcmElS+HFBke+/w+P+UEvMl4f58gvUyiu67q5T3ejvcXjzmYeF1NgKNi2y/LTz+mBLev2u478lSPucqoFsJ520kWBOiYQn7Ukv5rLMIEmqTuP896ZGYh1oOkgwmu3vxVeyeCp+zgbuL7Xs6fB6wZ4OZtQJ+Akx393uLHuzuuQTrXhhwQZFdncPnNfuJ7wV3f6/YtsfC529G+ZhZR4JO3eVA8RimAOOAlgStokS4192XlLJvN8GX/Xe4e1Ypx68luBORsDWMJV7qkJZkML2EbXtWPpvt7sW/5Pas09GxyLbDCNY+cDO7rYTr1Q2fi/YTtAqfN5cjvj1rPLcosm1g+DzJ3XeXcM77BAlsIN8muAMxrZTtzwGjgXlm9jLwEfCpu2fv41qbwufUfRwj1YiSgySDkr608kvb5+75QXfBN1/48O0X/WHhozRF+yl2hs/19xNf8eG5RWOoXWRzs/C5tJbInu3N9/N+ZbW2lO2/JOjvuAz4bfjIN7M3gevdPbOEcxqEzztL2CfVkG4riQT2JJG/uLvt43FskXPWh8+tSIw9MbQrZX/7YsfBtyuAlfSL3v6SSIlDb929wN3HuHt/oC1wFsEcjtOBiWZWr4TT9vwZrC9hn1RDSg4igWkEX7TDI5wzj+C+fO8ExTArfM4ws5K+7PckpplFtu25pdWphOMHH2hA7r7e3V9z93MJbmt1B/qUcGgvvu3IliSg5CBC8CVIcK99sJn9vqQvZzPrHg413XNONsEM7n5m1qD48eWIYSXwLsEoo+uKvffhBJ3hm/nuTOw9/QaXFo05nONxa9QYzKyemR0fDtMtur0uQWc4wI5i+7oRtDA+dPcaNxEwWanPQeRb1wDpwO3ARWY2mWCh+g4EHdGHAT8Gio7weRU4lGC+w4QExHAVwXyD+8zsJILO7E4E8wsKgUvdPWfPwe4+1cw+Bo4CppnZ+wRf1D8kGBZbUotiXxoQzGlYamZTgWUEfSonEvwZ/Nvd5xc7Z0/ZjFcjvpdUYWo5iITcfSvBBLnRQBbBvfZfEdzOySHoqH232GlPAHkEtZ8SEcNigttBfyO4VXMD8AOCeRZHuvu/SjhtBPA4weir0QSjmW4kGH4b1fbwvAUEcyuuJWixbAV+RsmT4EYCG1BySCqmVqDIgTGzRwm+ILu6e2kjgJKSmfUD5gC/d/c79ne8VB9KDiIHyMzaAouAp9x9dNzxVCYze53gtlpPd9cw1iSi20oiB8jd1xFMTlttRRb7SXZm1pBghNVFSgzJRy0HERHZS435LUdERMpOyUFERPai5CAiIntRchARkb0oOYiIyF6UHEREZC//H67I+TOCPjTmAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"n = np.arange(-3,3,1)\n",
"temperature = np.zeros(len(n))\n",
"for i in range(0,len(n)):\n",
" temperature[i] = current_temp(n[i])\n",
" \n",
"plt.plot(n,temperature)\n",
"plt.xlabel('Time(hours)',size=20)\n",
"plt.ylabel('Ambient Temperature',size=20);"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Yes the graph looks correct. Adding more data would smooth out the straight lines into curves as the temperature change should be gradual and continuous. "
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAEdCAYAAADuCAshAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nOzdd3hUZfbA8e+ZTDqkkAJJ6KFIkxYFpYgiKCLY0IVdBQEXdNe67rroTxTLuqjrWlZFscKuoiAoCoqoWCkiTar03gk1pJDMvL8/7k2YJJMwqZNyPs9zn5l731vOBJ0z975NjDEopZSq3Rz+DkAppZT/aTJQSimlyUAppZQmA6WUUmgyUEophSYDpZRSaDJQXoiI8WHp4+84K4uIJIrIBBFp6O9YSkJE+ojIKhHJFJHMIvY5z/73vLwC41grIq9V1PlV+XD6OwBVJV3k8T4UWAA8Ccz12L6+UiPyr0TgUWAesMfPsZTEm8BW4G7AazJQKpcmA1WIMWZJ7nsRqWO/3eq5vboTkQDAYYzJ9sO1g4AcY4y7Aq8RCDQHnjHG/FBR11E1hz4mUmUiIs1EZIaIHBeR0yIyV0SSPcpzH0PcICL/E5FTIrJLRH5nl/+fiOwXkUMi8oSIiMexE0VkT4HHHctFpLuXOO4QkQ0ikiUi20Xk3gLlH4jITyJyk4hsALKATiLSSESm2MdkiMhGEXnU/jJFRM4DfrFPs9j+LJl22e32urPAtQ6IyJMe60vsz36niGwHMoAYu6yTiMyz/y4nRGSaiMT58HfvLyK/2H+TAyLykoiE2mVXAmcAAV63YzzXY5po+2+UZp/vwZJc02OfTvbnzRSRdSIyoED5DSKSIyJJBbbn/nfS/1yfXVUMTQaq1EQkHlgINAVuA4YBscB8+9evp+ewHllcj/Xl+j8R+TfQHrgVeBV4GLimwHERwDvAf4AbsR53zBORGI84xgMvANOBgcBbwDMicluBc7UCHgeeAK4CdgPxwAHgXuBK4HngDuBf9jE7gJH2+9uwHqH1Pucfp7C+wHDgfvszpotIG+BHu/wP9vm7Ah8XdyIR6Yz1yG4v1t/zCTvGafYui4FL7Pf/tGN++hzxvQAcAW4ApgBPicjoElwTEakLzAcCgaHARKx/1wSP63wKpAK3FLj+SKxHcF+fI05VUYwxuuhS5ALUAQxwq5eyZ7G+SCM9tsUBacBoe/08+/hJHvvEAG5gLSAe21cDUzzWJ9rHXu+xLQo4BUyw1+th/dL+e4HYngF2eax/YF+zTTGfVbAenY6yrxFgb0+x4+heYP/b7e3OAtsPAE96rC+x/yYxBfabAazxPB5oZ8fZt5g4PwHWYT3myt023I6ls70eYq/fdo5/39x/n08LbP8vsL2E1/wLVrKu77FPX3uf1zy2/Qv4zWM9ANjn+TfTpfIXvTNQZXE5VqXqaRFx2o9LjgG/Yn2Bevom940xJhU4Dnxn7G8D2xYgqcBxLqxfk7nHHseq0L7Q3tQL64tvRm4MdhzfAI1EpL7HubYZYzZ4nlxEHCLyNxH5DSupZGPdWdQh/y/aslpif25PlwOz7Dhy494I7Kfw38/ThcBMk7/OYTrWl27PUsZX8G5kFtDUvvvz9ZoXYn3Og7k7GGO+AU4WOPdbQGsRyW2ocAXW3/rdUsauyoEmA1UWscAIrC9Qz+VioFGBfY8XWD9TxLaQAtuOGWNyCmw7xNkv6lj7dWuBGObZ2z3jOEhhfweeAj4EBmF9od1nlxWMpSzyXduuwI4CHqHw3y+Rwn+/3OMEqF/wfMaYTKwv3XqljO9QEesJJbhmAy/nKXRuOyEv5uzjt5HAD8aYLaWMXZUDbU2kyuIo1iMQb8+jT5TTNaJFxFkgIcRj/XrOjQGgP9ZdSUGedwLexmu/EXjPGPNo7gYR6eJjbLnNNYOAHPtYASK97Jvv2sYYl4icxPo1/F8v+3v7UsUYY0TkINbfII+IhGDVrxz1dpwP4otY31+Cax7ASgjnOjdYzV6fF5HHsJLw7aWMW5UTvTNQZfENVgXwamPMsgLL5nK6RgAwOHdFRCKBy4Cl9qafsO4oGniJYZkx5vQ5zh+K1bLI0x8KrJ+xXwveKeT2OWjjsa23l/2K8g3Qroi4dxVz3M/ADZ4tr7CSmmD9PUrjugLr1wM7jTG5ScmXa/4CdPd8NCcifbESRkHTsb5/PsC6G/qolHGrcqJ3BqosnsFqNfKNiLyC9Wu9AdAH+NoYM7McrnEK+LeIRGP9Wh6HVcH6MoAx5rCI/AOYJCItsL6YnEBr4GJjzE3nOP9XwGgRWQHsxHrsVbCn8XashDBSRLKALGPMCvtah4BX7F+48Vithc6VgHKNB5aIyKdYdwhH7WtfgVXhvqiI4x7H+uKdKSJvAM2wKttnG2NW+njtgrqKyH+Az7AqfW8GxpTwmpOx/n0+F5EngLr2cYXuVowxaSLyITAaeMcYk1bKuFV58XcNti5Ve6GY1kR2eSNgKtaXYibWF+cUoLVdntta5fICx+VrcWNv+wD4yWN9Itav78uwWt1kASuxvuQLxjHSLsvE+vJZDNxV1Lk9tkdgPaY5htXk8TWsX8UGaFHg/FuwfsVmemy/GFgBpAPLsOocvLUm+l8Rf792WJW3x7AqsDdjN8c8x7/LFfb1srCe5b8EhHmUl7Q10Y1YrZvS7PM9XNJr2vt0wbqLyMJ6RDcQq9XYa17Od7V97Z7+/u9cF2M161OqKhKRicDNxphqNSaQ8o2IvARcaYxp5e9YlD4mUkpVMrtXd3usTnYP+DkcZdNkoJSqbO8CHYGZWI/lVBWgj4mUUkpp01KllFLV+DFRbGysadq0qb/DUEqpamX58uVHjDGFRsattsmgadOmLFu2zN9hKKVUtSIiO71t18dESimlNBkopZTSZKCUUgpNBkoppdBkoJRSCj8kAxG5R0TW2pNl32tvqyciX4nIZvs1urLjUkqp2qxSm5aKSHvgj1gjO57Bmth8rr3tG2PMRBEZhzUM7t8rMzblHydPnuTQoUNkZ2f7OxSlqrXAwEDi4+OJiPA2fcS5VXY/gzZYc6SmA4jI91iTalyDNQY+WMMff0cFJYMFvx3ktwOn+FOfFhVxelUCJ0+e5ODBgyQlJREaGkr+eVOUUr4yxpCRkcHevXsBSpUQKvsx0Vqgt4jEiEgYcBXWePj1jTH7AexXb9PkISJjRGSZiCw7fPhwqQJYtCWVF7/eTI7Lfe6dVYU6dOgQSUlJhIWFaSJQqgxEhLCwMJKSkjh0yOuMqedUqcnAWBNhP401u9Q84FfsuWN9PH6yMSbFGJMSF1eoN7VP2iREkJXjZkeqr5NRqYqSnZ1NaGiov8NQqsYIDQ0t9SPXSq9ANsa8ZYzpYozpjTUj1WbgoIgkANivpUttPmiTYN0+rd9/qqIuoUpA7wiUKj9l+f/JH62J4u3XxljTC04DPsWaexb7dXZFXb9FfB0CA4T1+05W1CWUUqra8cdAdTNFJAZrLtk/G2OO2dMbTheR0cAurPlYK0SQ00FyXB027NdkoJRSufzxmKiXMaatMaajMeYbe1uqMaavMaal/Xq0ImNomxChyUCVK2MMzZo1Q0TYsmVLhV1nwoQJxMbGluiYpUuXMmHChHI5ly/69OnDkCFDyv28qmLVyh7IbRMjOHQqi9S0LH+HomqIxYsXs2PHDgA++OAD/wZTwNKlS3nssccKbb/tttv48ssv/RCRqopqZTLIrUTeoJXIqpxMmzaN8PBwunXrxrRp0/wdjk8aNmxI165d/R2GqiJqeTLQR0Wq7FwuFzNmzGDw4MGMGjWK9evXs3r16rzyd999FxFhzZo19OvXj/DwcM477zxmzZqV7zxz586lX79+eb1Iu3fvzvz584u8bk5ODomJiV5/9V9yySVcf/31vPvuu9x1112A1dJEROjTpw/g/TFRamoqY8eOJSEhgZCQEFq3bs0LL7yQV/7cc89xwQUXEBkZSf369Rk0aFCFPhZTladWJoN64UHUjwjWZKDKxYIFCzh48CBDhw5lyJAhBAYGer07+P3vf8/gwYP5+OOPadmyJUOHDmXPnj155du3b2fQoEH897//ZebMmVx88cUMGDCAhQsXer2u0+lkxIgRvPvuuxhj8rZv27aNH3/8kZEjRzJw4EDuv/9+wHqUtXjxYl599VWv58vIyKBPnz588sknjB8/ns8//5z777+fffv25e2zZ88e7rzzTmbPns0bb7yBy+WiR48enDhxolR/O1V1VNtpL8uqTUIE6zUZVDmPfbbOb81+2yZG8OigdiU+btq0aURFRXHllVcSFBREv379+OCDD3jqqafytfu+7777GDVqFABdu3alfv36zJkzh9tvvx2AO++8M29ft9vNpZdeyrp163jrrbfo0aOH12uPGjWKiRMn8t1333HppZcC1p1IfHw8AwYMwOl0kjtXePfu3Yv9HFOnTmXdunWsWLGCTp06AXDZZZfl2+f555/Pe+9yufLuZGbPns3w4cN9+XOpKqpW3hmAlQy2HEojK8fl71BUNZaVlcXHH3/MddddR1BQEADDhg1jx44dLFmyJN++/fv3z3sfExNDfHx8vjuDPXv2MGLECJKSknA6nQQGBjJ//nw2bdpU5PVbtmxJ7969effddwGrVdPUqVO55ZZbcDpL9ltvwYIFdO7cOS8ReLNkyRL69etHTEwMTqeTsLAw0tLSio1RVQ+1+s4gx23YciiNdomR/g5H2Urzy9yfvvjiC44fP85VV13F8ePHAatpZXBwMNOmTeOiiy7K2zcqKirfsUFBQWRmZgLWncDgwYM5deoUjz/+OC1atCA8PJxHHnnknGPNjB49mjvuuIOXX36ZpUuXsnPnTkaOHFniz5KamkpCQkKR5bt27aJ///5ceOGFvP766yQmJhIUFMTAgQPzPoeqvmptMmjr0aJIk4Eqrdy6gRtvLNxPcvr06fkeqxRny5YtrFy5ki+++IIrr7wyb3tGRsY5j73xxhu5++67mTFjBt9++y3dunWjbdu2Pn6Cs2JiYoqtDJ43bx7p6enMnj2b8PBwwKrEPnq0QrsFqUpSa5NBs9hwQgIdWomsSi0tLY05c+YwbNgwxowZk69s5cqV/OUvf+Hbb7/16Vy5X/rBwcF523bu3MnChQs5//zziz02NDSUYcOG8corr/Dbb7/x73//O1957uOrzMxMQkJCijxP3759mTFjBqtXr/Z6zYyMDBwOR77HT9OnTycnx+exJlUVVmvrDAIcQuv6dTUZqFKbPXs26enp3HPPPfTp0yffcueddxITE+Nzn4PzzjuPhg0bcv/99zN37lw++OAD+vfvT1JSkk/Hjx49mhUrVmCMYejQoYXODfDiiy/yyy+/sHHjRq/nGD58OB06dKB///689tprfPvtt7z99tuMGzcOsCqTXS4XI0eO5JtvvuGll15i3LhxhR5/qeqp1iYDsOoNNuw/ma9ZnlK+mjZtGi1btqRbt26FygIDA7npppuYNWsWWVnn7ukeHBzMrFmzcDqdDBkyhPHjx/Pggw9yySWX+BRLSkoKSUlJXH/99URG5n/s2atXL/72t7/x4osv0q1bN8aOHev1HCEhISxYsIBBgwbxyCOPMGDAAJ555hkSExMB6NChA++88w4///wzV199Ne+//z4zZswodD1VPUl1/SJMSUkxy5YtK9M5pizawaOfrmPJg31pEFn07bOqGBs2bKBNmzb+DqNGWL9+Pe3atePrr7+mb9++/g5H+dG5/r8SkeXGmJSC22v9nQHA+v3aYUZVT6mpqSxatIi77rqL9u3bF+oXoJSvanUyOC+hLqBjFKnq67PPPqNnz57s378/b9gLpUqjVieDiJBAGtUL1Z7Iqtq69dZbcbvdrF+/XgedU2VSq5MBQJsGOreBUkppMkiIYMeR02Sc0WEplFK1V4mSgYgkiUhXEekhIm1EpNo3wWmTEIHbwMaDWm+glKq9zpkMRKSviEwRkb1Y8xMvBX4E1gInRGSJiPxVRBpUcKwVoq3ObaCUUkUnAxG5QUTWAp8BkcCLwFXAhcD5wCXAaGAJcBuwQ0ReqW5JoWF0KHWDnX4bNlkppaqC4sYmegp4FphmjDldxD4/Af8DEJGOwL3ArcDEcoyxQjkcwnkJOiyFUqp2K+4x0XnGmDeLSQT5GGN+NcaMBJ4un9AqT5uECH47cAq3u3r2xlb+M2HCBESEK664olDZkCFD8qaYrCwiwssvv5y33qdPH4YMGVKpMZTWd999h4gQGxtLWlpavrKXX365WvWh8DalaHmoyH/PIpOBKeU4FaU9zp/aJESQlpXDnmPnHi5YKW/mz5/PL7/84u8wCnn11Vf55z//6e8wSiQ1NZVJkyb5O4wyue222/jyyy/9HUaJFFdnsFpE2hfY9nsRqXFDFJ4dlkIfFamSq1evHueffz7/+Mc//B1KIW3btqVly5b+DiOfc02E06dPH5577rlqOWFOdnY2LpeLhg0bVrtOgMU9JmoPhOWuiEgA8F+geUUHVdla16+LQ7RFkSodEeGhhx7i008/Zc2aNcXuu2rVKvr27UtYWBjR0dH84Q9/4ODBg3nlO3bsQESYPn06Y8eOJTIykoYNG/Loo4/idrtLHFvBxwq5jy9WrlxJ9+7dCQsLo3Pnzvz444+Fjn3zzTdp164dwcHBNGnShGeeeSZf+eLFixk8eDCJiYmEh4fTqVMn3nvvvXz75A6RsXTpUvr06UNoaCjPPvtssTE/8MADHDt2jDfffLPIfXIfKa1du7bYz3vrrbeSkpLC3Llzadu2LWFhYQwcOJCjR4+yZcsWLr30UsLDw0lJSWH16tX5zuV2u5k4cSItWrQgODiYVq1aMWXKFK/Xmzx5MsnJyYSEhLBv3z6vj4lSU1MZO3YsCQkJhISE0Lp1a1544YW88ueee44LLriAyMhI6tevz6BBg4qdbKi8lbTTWfV5aFcCoUEBNI0N12SgSu3GG2+kVatWxd4dHD58mD59+pCens7777/Pf/7zH77//nv69evHmTNn8u37wAMPUKdOHT766CNuvvlmHn/8cT766KNyiTU9PZ0RI0YwduxYZs6cSXBwMNdddx3p6el5+zz77LPccccdXHvttcyZM4c77riD8ePH56uP2LlzJz169ODNN9/ks88+44YbbmDkyJFe53AYNmwYV199NZ9//jlXX311sfE1atSI4cOH88wzz5CdnV3mz7tr1y4eeeQRnnzySSZPnsyiRYsYM2YMQ4cOZejQoXz00Ufk5OQwdOjQfMPZ33XXXTz55JOMGTOGuXPnct111zFq1CjmzJmT7/wLFy5k0qRJPP3003z22Wdeh/TOyMigT58+fPLJJ4wfP57PP/+c+++/n3379uXts2fPHu68805mz57NG2+8gcvlokePHpw4UTkDadbamc4KapsQwardx/0dhvpiHBwo/td1hWnQAQaUriGcw+Fg3LhxjB49mscff5xWrVoV2ue5554D4MsvvyQiwno02apVK7p168bMmTMZNmxY3r69e/fO279fv37MmzePWbNmcdNNN5UqPk8ZGRm88MILeSOcJiQk0LlzZ3744QeuvPJKTp48yWOPPcbDDz/Mo48+mhdDeno6Tz75JHfccQcBAQH5JtExxtC7d2/27NnDG2+8ke+zANx9993cc889Psc4btw43nnnHaZOncro0aPL9HmPHj3K4sWLSU5OBmD16tU8++yzTJkyheHDh+fFP3DgQH777TfatGnDli1bmDRpEu+88w4jRowA4PLLL2f//v089thj+RLa8ePHWblyJQ0aFN2qfurUqaxbt44VK1bQqVMngEIjzHpOkepyuejXrx/x8fHMnj07L86KdK47A4eIOETEAQQU3Oa5VHCcFa5NQgR7jmVwMrPsv0RU7XTzzTfTuHHjIitsly5dSv/+/fMSAcCFF15I06ZN+emnn/Lt279//3zrbdu2Zc+ePXnrOTk5eYvLVbKhVAIDA/O1csqdLzn3/IsXL+b06dPceOON+a5z2WWXcfDgwbz9jh07xt13302TJk0IDAwkMDCQyZMns2nTpkLXHDhwYIliTE5OZujQoUycOLHEn6+gpk2b5iUCgBYtWgD5v4xzt+3duxeAb775BofDwXXXXZfvb9C3b19WrVqVL6auXbsWmwgAFixYQOfOnfMSgTdLliyhX79+xMTE4HQ6CQsLIy0tzevfsyKc685goZdtPxexb0AR26uFDknWrd2aPSfo0aL8m4QpH5Xyl3lV4HQ6eeCBB7j77ruZMGFCofL9+/fTrl27Qtvr169faFL5glNJBgUF5VWo7tixg2bNmuWVNWnShB07dvgcZ0REBA7H2d9vnnMkAxw5cgTAa6wAu3fvpkmTJtx6660sWbKE8ePH07ZtWyIiIpg0aRKzZ8/2+hlL6qGHHqJ9+/Z8+OGHJT7Wk7e/ZcHt3v4GLperyFnc9u/fT8OGDQHfPltqaioJCQlFlu/atYv+/ftz4YUX8vrrr5OYmEhQUBADBw6stIr04pLBHZUSQRXRsZH1H8aKncc0GahSGzVqFE8++SRPP124u01CQgKHDh0qtP3gwYMlanmSmJiYrxlrcHBw6YItQr169QCYM2eO1y+61q1bk5mZydy5c3n55Ze5/fbb88qKquQuTR+Btm3bct111/HUU08VmqozJMQaFq1gXcvRo0fLpX1/vXr1cDqdLFy4MF/izBUfH5/33pfPFhMTU2xl8Lx580hPT2f27NmEh4cD1t1fwR8JFanIZGCMeb3SoqgCIkMDaRlfhxW7jvk7FFWNBQcH89e//pUHH3yQrl27EhgYmFfWrVs3Jk2axKlTp6hb15pY6ZdffmHHjh307NnT52sEBQWRklJo1sJyc9FFFxEaGsq+ffuKfLxz4sQJXC5XvkR06tQpPv3003LtHPbwww/TpUsXPv7443zbc3+Vb9iwgS5dugDWHcvGjRu91teU1GWXXYbL5eLEiRP069evzOfr27cvM2bMYPXq1Zx//vmFyjMyMnA4HDidZ7+Sp0+fTk5OTpmv7atKr0AWkfuwxjIywBpgJDAO+CNw2N7tIWPM55UdW5fG0Xy5/gDGmGrV21FVLWPHjuWpp55i0aJF+Sa0/8tf/sKkSZO44oor+Pvf/05aWhrjxo2jQ4cO3HDDDX6MOL+oqCgmTJjAPffcw86dO+nduzdut5tNmzbx7bff8vHHHxMZGckFF1zA448/nvfYaeLEiURGRnLyZPm1yuvcuTMDBgzgiy++yLe9YcOGXHDBBYwfP56wsDDcbjdPPfVU3l1NWbVu3Zrbb7+doUOH8sADD5CSkkJmZibr1q1j06ZNxTZ79Wb48OG88sor9O/fnwkTJtC6dWu2b9/Opk2bmDhxYl7yGTlyJKNHj2bdunX861//KvSIqyJVasWviCQBdwMpxpj2WPUMuU0SnjfGdLKXik0E2d6fwXVpEsXx9Gy2HfFpBA6lvAoLC+O+++4rtD0uLo5vv/2WkJAQhg0bxp///Gd69erFV199lffMuqp44IEHmDx5Ml988QXXXHMNw4YN47333qNXr155+7z//vs0a9aM4cOHc88993DDDTdUSKuXhx9+2Ov2999/n8aNG3PzzTfz0EMP8cgjj9C6detyu+4rr7zC+PHjmTp1KldddRW33norc+fOpXfv3iU+V0hICAsWLGDQoEE88sgjDBgwgGeeeYbExEQAOnTowDvvvMPPP//M1Vdfzfvvv8+MGTOKrLOoCFKZo0fYyWAJ0BE4CXwCvARcDKQZY/7l67lSUlLMsmXLSh7EF+Ng41y4t3Dzxc0HT9Hv+R94dsj53JjSqOTnViWyYcMG2rRp4+8wlKpRzvX/lYgsN8YUes5YqXcGxpi9wL+w5kXYD5wwxsy3i++0h8B4W0SivR0vImNEZJmILDt8+LC3Xc4tMgmO74K0whV5yXF1iAhxsmKX9jdQStUulf2YKBq4BmgGJALhInIzMAlIBjphJYnnvB1vjJlsjEkxxqTExcWVLogku9XG3hWFihwOoVPjaFZqJbJSqpap7M5ilwPbjTGHjTHZwCzgYmPMQWOMyxjjBt7AmkCnYiR0BHHAvsLJAKBzoyg2HTxFWlbl1eIrpZS/+ZQMRCRDRNKLWNJEZJ+IfCEiA85xql1AdxEJE6u5Tl9gg4h49sa4DmtKzYoRFA5xbWDvcq/FXZpE4zbwqw5NoZSqRXy9M5gAHAH2AK8DTwKTgX3AUeBdoC4wR0SKHDzFGPMz8BGwAqtZqcM+zzMiskZEVgOXAoWbYpSnpC5WMvBSed7Jo/OZqnjVcPoLpaqssvz/5Gs/gwhgKXCj5+Q1dp+Bj4AgY0xPEZkGPAhMLybYR4FHC2y+pURRl1VSV1j5Xzi2HerlH5FbO59VnsDAQDIyMggLCzv3zkqpc8rIyMjX0bEkfL0z+CPwRsFZzOz1yVjzHgO8B5xXqkgqUzGVyGB1Plu5+7j+aq1g8fHx7N27l/T0dP1bK1UGxhjS09PZu3dvvqEySsLXO4NgoAXgbR63lh7nyQSyShVJZYpvA84QKxl0KDyfaJcmUXy4bDfbjpwmOa6OHwKsHXJH79y3b1+5jFuvVG0WGBhI/fr1842KWxK+JoPpwEQRMcBnWMNGxGE1E50IvG/v1xnYWKpIKlNAoNWqqIhK5M6NrW4OK3cd12RQwSIiIkr9H69Sqvz4+pjoLqwv/OeBHcBp+/U5rEdDd9v7rQLuLdcIK0pSV9j/K7gKNyFtEVeHuiFOrTdQStUaPiUDY0ymMWYs0BgYCIwCrgKaGGPGGmMy7f2+MsYsrrBoy1NSV8jJgMMbChU5HEKnRlHaokgpVWuUqNOZ3TlsnjFmijHmS2PMgYoKrMIldrZei+pv0DhaO58ppWoNn5OBiLQWkSkisk5ETohIJ3v7BBEp+4Dfla1ecwiJ0s5nSimF7z2QLwdWA62xRhqt63GsAL7PdF1ViFiPivau9Fqsnc+UUrWJr3cGTwMfGGO6U7jD2AqsVkTVT1JXOLQezhSevyAyNJAW2vlMKVVL+JoM2gH/td8X7B10HIgpt4gqU1IXMC7Yv9prcZfGUdr5TClVK/iaDI4ATYooa4M1ZlH1k2jNnYdqKbUAACAASURBVFpcJfLx9Gy268xnSqkaztdkMB14QkQ8Z8cxItIM+DvwQblHVhnq1ofIRkUOZ92lidX5TCe7UUrVdL4mg/8D1gE/A5vsbTOADcAW4PHyD62S5I5g6oV2PlNK1RY+DUdhjMkA+onIQKw5CGKxhq7+BphTcAC7aiWxC6yfDadTITx/1Yd2PlNK1Ra+jk0EgDFmLjC3gmLxj9wRTPetgJaFu0t0bhzNyws2k5aVQ53gEv25lFKq2ijy201ESjQOqjGm8Azz1UFiJ0CsEUy9JIMujaPyOp/1aBFb+fEppVQlKO6n7gEKNyMtTkAZY/GP4LoQd16xPZEDHMLiramaDJRSNVZxyeBGj/fhwD+AbcDHwCGgPtZ8xc2wKpirr6QusOlLaxpMkXxFESGBdGwYyU9bjvDXK1r7KUCllKpYRbYmMsbMzF2AXsCXxphLjDEvGGPeN8Y8b4zpDczHmre4+krqAulH4Pgur8U9W8axes9xTqTrBCxKqZrJ16alNwEfFlH2IdYdQvXlWYnsRc8WsbgNLN6WWolBKaVU5fE1GWQB3Yoo6w5U75/M8e0gILiYmc+iCA8K4Kcthys5MKWUqhy+tpV8E3hERKKAT7HqDOKxpr28G3i2YsKrJM4gaNDBalHkRWCAg27NY/hp85FKDkwppSqHr8ng/4ATwP3AX7BaGQmQCozHGtW0ekvqCiv/B24XOAo3jOrZIpYFvx1i99F0GtUL80OASilVcXyd9tIYY54GkrAGprvMfk00xkys1j2QczVMgezTcGCN1+KeLa1mpQu36N2BUqrmKem0l9nGmI3GmO/t1+pdV+Cp2SXW69YFXotbxtchvm4wP2kyUErVQEUmAxH5m4iEl+RkInKxiAwoe1h+ULc+1O9QZDIQEXq2iGXR1lTc7up/I6SUUp6KuzPoB+wRkTdEpL+IRBTcQSxtReR+EfkFmI1Vl1A9JV8Ku5ZAVprX4p4tYzl6+gzr95+s5MCUUqpiFdfprD9Wa6ForBZEx0Rkp4gsF5FFIrIeOAmsAcYCM4HmxpjPKyHuitGiL7izYedCr8W5w1HooyKlVE1TbGsiY8wPwA8iEg30BroADYAQrCGsNwILjTHea12rm0bdwRkKW76BVlcUKq4fEUKr+nVYuOUIt1+S7IcAlVKqYvg6n8ExrEdAsys2HD8LDIGmPYqsNwDr7uD9n3eRme0iJLB6js2nlFIFlag1UXkQkftEZJ2IrBWRaSISIiL1ROQrEdlsv0ZXdlx5kvtC6uYixynq1TKWrBw3y3XCG6VUDVKpyUBEkrB6LKcYY9pjDXs9FBgHfGOMaYk1e9q4yowrn+TLrNci7g4ubBaD0yH8qL2RlVI1SKXfGWA9mgoVEScQBuzDqqieYpdPAa71Q1yWuNYQkVRkMqgT7KRL42gdp0gpVaNUajIwxuwF/gXsAvYDJ4wx84H6xpj99j77scY9KkRExojIMhFZdvhwBX0Zi1hNTLd9B64cr7v0aBHLun0nOXr6TMXEoJRSlayyHxNFY90FNAMSgXARudnX440xk40xKcaYlLi4uIoK03pUlHkC9q30WtyzZSzGwKKt+qhIKVUzlCgZiMilds/kl0Skob2tu4jU9/EUlwPbjTGH7aEsZgEXAwdFJME+XwLWqKj+0/xSQGDrN16LOzaMpG6wU8cpUkrVGD4lAxGJFZEfgK+Be4A/c/ZRzp+wRi71xS6gu4iEiYgAfYENWJ3aRtj7jMDfTVjD6kFi5yLrDZwBDronx/Dj5iPUhDH6lFLK1zuDl7DmPO4ANCX/kBPzsYauOCdjzM/AR8AKrJ7LDmAyMBHoJyKb7XNN9DGuitOiL+xZBhnHvRb3ahnLnmMZ7DqaXsmBKaVU+fM1GVwFPGSMWY81l4Gn3UBDXy9ojHnUGHOeMaa9MeYWY0yWMSbVGNPXGNPSfj3q6/kqTPJlYFyw/QevxblDU2gTU6VUTeBrMnBgTX3pTT0gs3zCqUIaXgBBdYt8VNQ8NpyG0aEs+M2/1RtKKVUefE0GC4E77Of8uXLvEG4FvivHmKqGgEBo1tuqRPZSLyAiXNGuAT9tPsKpzJozrYNSqnbyNRk8CPQCVgEPYyWC4SIyH+hjb6t5ki+1hqU4us1r8YD2DTjjcuvdgVKq2vN12stVQDdgE1ZrIgFGAqeAi4wxGyosQn9q0dd63eK9iWmXxtHE1Q1m3toDlRiUUkqVv3MmA3sCmxhghzHmRmNMPcBpjIkyxtxgVyrXTPWaQ3TTIusNHA7hinb1+W7jYTLOuCo3NqWUKke+3Bk4gQPApbkbjDHuCouoqknuCzt+hBzvQ08MaJ9ARraL7zfpWEVKqerrnMnA7im8Gwiq+HCqoOTL4Ewa7P7Za/GFzeoRFRbIvLX7KzkwpZQqP75WID8HjBORqIoMpkpq3gcCw2DtTK/FgQEO+rWpzzcbDpGVo4+KlFLVk6/JoDvQGNglIvNF5L8iMtVjmXKuE1RbwXWgzSBYNwuyvXenGNChAaeycli0NbWSg1NKqfLhazJoAewE1gF17fWWBZaaq+NQaxTTTV94Le7RIpY6wU7mrdFWRUqp6snXOZAvquhAqrRml0DdBPj1A2h3XaHiYGcAl50Xz/z1B/iHqz3OAH/MGaSUUqWn31q+cATA+TfB5q8gzXuroQHtG3AsPZulO/w/rJJSSpWUT3cGIvL4ufYxxjxS9nCqsI7DYOGLsPYj6H5HoeJLWscREuhg3toDXJwc64cAlVKq9HxKBsAfvWyLBEKANHup2ckgvg0kdIJV73tNBmFBTi5pFceX6w4wYVA7HA7xchKllKqafB2OIsHLEobVEW07MLhCo6wqOg6DA6vh4DqvxQPaJ3DwZBYrd3ufA0EppaqqMtUZGGO+B/4NvFo+4VRxHYaAw2lVJHtxWZt4AgNEO6Appaqd8qhAPgi0LYfzVH3hsdCyP6yeDu7CHcwiQgLp0SKWeesO6HSYSqlqxdc5kB1elhAR6YxVV/BbxYZZhXQcCmkHYNu3XosHtG/A7qMZrNt3spIDU0qp0vP1ziAHyC6wnAaWAc2BOyskuqqo1ZUQElXko6J+bRvgEHRYa6VUteJra6I/UXju40xgD/CTMaaoKTFrHmcwtL/BalWUeRJCIvIV1wsPonvzGD5bvY+/9GulrYqUUtWCrz2QX6voQKqVjsNg2Vuw4VPofHOh4htTGnLfh7/y8/ajXJQc44cAlVKqZHytM0gXkZQiyjqLSHr5hlXFNUyBesmwaprX4gHtE6gb4mT6st2VHJhSSpWOr3UGIcXsGwwElE841YQIdBoGO3+CYzsLFYcEBnBtpyQ+X7OfExnZfghQKaVKpshkICKJInKhiFxob2qbu+6x9AbGYo1oWruc/zvr9Vfvdwe/u6ARWTluPl21txKDUkqp0imuzuCPwKNYFccGeMvLPgKcwUoItUtUY2jRD355C3rcC4Eh+YrbJ0XSLjGCD37ZzS0XNfVPjEop5aPiHhNNBi4AumF96Y8CLiywdARijTE1d3Kb4lx8F5w+BGumey0eekEj1u07ydq9Jyo5MKWUKpkik4ExZr8xZrkxZhnQBvjAXvdc1hhj0iov3CqmWW9ocD4sehnc7kLFgzslEex08OEvWpGslKrafB2obqMxJkssiSLSvOBS0YFWSSJw8d1wZCNs+apQcWRoIFd1SOCTVXvJzNb5kZVSVZevTUudIvI8cBzYDWz2stRO7a6FiIaw6D9ei29KacSpzBy+0MHrlFJVmK9NSx8Cfgfci1V/8BesXskLgR3ADRURXLUQEGjNb7DjR9i7olBx9+b1aBoTxgdL9VGRUqrq8jUZ/B6YAEy1138yxrxujOkN/Az0q4DYqo8uwyE4Aha/XKhIRLjpgkb8vP0o24+c9kNwSil1br4mg8bABmOMC8gCojzKpgA3+XISEWktIqs8lpMicq+ITBCRvR7bryrZx/CzkAjoeius+8RrJ7QhXRoS4BDtkayUqrJ8TQYHsKa5BOuxUA+Psia+nseuiO5kjOkEdAXSgY/t4udzy4wxn/sYV9XR7XarQvnnwsM4xUeEcGnreD5avoccV+FWR0op5W++JoMfOJsA3gYeFpG3RWQS8BwwpxTX7gtsNcbUjN7LkUnQfggsnwIZxwoV/+6CRhw+lcW3Gw/7ITillCqer8ngYWCG/f5fWD2TOwKXYCWH0sxnMBTwHMvhThFZbSeZaG8HiMgYEVkmIssOH66CX6oX3wnZp2H5u4WKLm0dR3zdYD78ZVflx6WUUudwzmQgIk6gPtbcBRjLP40xXY0xbY0x9xhjTpXkoiISBAzmbIKZBCQDnYD9WHcbhRhjJhtjUowxKXFxcSW5ZOVo0AGa94Elr0HOmXxFzgAHN6Y0ZMFvh9iZqhXJSqmqxZc7AzewGDi/HK87AFhhjDkIYIw5aIxxGWPcwBtYQ11UTxffZU2LufajQkUjLmqK0+Fg8g/b/BCYUkoV7ZzJwP6C3gLEluN1h+HxiEhEEjzKrgPWluO1KldyX4hvBz/+G1z5h6+Ojwjhhq4NmbF8D4dOZfopQKWUKszXOoNHgUdEpFVZLygiYVj9EmZ5bH5GRNaIyGrgUuC+sl7Hb0TgsochdbPXuoOxvZuT43LzzsIdlR6aUkoVxdc5kO8GYoD1IrINOEiBOZHtDmjnZIxJt8/lue0WH+OoHloPgKa94Lt/wvk3QUhkXlHT2HAGdEjgf4t3ckefZCJCAv0YqFJKWXy9M9gDLMCq8F1ur+8tsKhcItD/CUhPhZ+eL1R8xyXJnMrK4X9LakarWqVU9efTnYExZlhFB1LjJHaG84fC4lchZZQ1GY6tfVIkvVrG8vZPOxjVoxkhgbVr1lClVNXj651BHhGJE5EuIhJaEQHVKH3HW3cJ3zxRqOiOPskcScvio+V7/BCYUkrl53MyEJFRIrIDa2iKX7AmvEFEPhSRP1dMeNVcZEO46M/WTGh7l+cruqh5DB0bRTH5h206RIVSyu98nc/gHqyOYR8CV2ENY51rEdaopsqbnvdBeBx8+TCYs3XuIsKf+iSz62g6n6894McAlVLK9zuDe4DHjDF/B74uULYRaF2uUdUkwXWhz4OwaxH8NjdfUb829UmOC2fSd1sxxhRxAqWUqni+JoNEYEkRZTlAWPmEU0N1GQFx58FXj+QbpsLhEG6/JJkN+0/y3aYqONaSUqrW8DUZbAN6FlHWE9hQPuHUUAFO6PcEHN0Ky9/JV3RNpyQSIkOY9N1WPwWnlFK+J4P/AA+JyF+BRva2KBH5A3A/8FJFBFejtOwHzS6xOqKdOpi3Ocjp4LZezVm6/SiLth7xY4BKqdrM10lpJgH/AB4Hcn/CfgW8CTxjjJlSMeHVICJw1b8gOwM+uztfZfIfujUmMTKEJ+dswOXWugOlVOXzuWmpMeYJIAm4HrgNuAFoZG9XvohrBX0fhU3zYOX/8jaHBAbw9wHnsX7/SWau0H4HSqnKV6JOZ8aYY8aY2caYt40xnxhj9LlGSXW73Rq3aN64fPMlD+6YSKdGUTz75UZOZ+X4MUClVG1Ukk5n0SLyiIjMEZHl9ut4EalXkQHWOA4HXPsqIPDJn8BtdTgTEcZf3ZbDp7J47XutTFZKVS5fO511w6or+CtWU9Ll9uvfgK12ufJVVGMYMBF2/gQ/T8rb3LVJNIM6JjL5h23sPZ7hxwCVUrWNr3cGrwDrsOoIrjXGjDHGXIvVsmi9Xa5KotMfoNUA+PoxOPRb3ua/X2n133t23m9FHamUUuXO12TQDphojDnhudFen2iXq5IQgcEvQXAd+Hhs3qxoDaPDuK1XMz5ZtY9Vu4/7OUilVG3hazL4DShqBvo4YHP5hFPL1ImHq5+H/avgx+fyNt/RpwVxdYN5Ys56HaZCKVUpfE0G9wLjReQaEXEAiIhDRK4FHsaaCU2VRttr4PzfwffPwI6FANQJdvLX/q1YvvMYc1bv93OASqnawNdk8D8gFmve4kwROQhkAjOxprCcKiK7cpeKCbUGG/AM1GsOH/4BUq2WREO6NqJtQgQTv/iNzGyXnwNUStV0vs6B/B4F5jxW5Sg0Cn7/IbzZF6YNhdFfERAaxcNXt+H3b/zMK99u4f7+OjCsUqri+Drt5biKDqTWi0mG3/0Ppl4LM26FP8zg4uRYbujSkFe/28rlberTsVGUv6NUStVQJZ72UlWgpj1h0Auw7Vv44gEwhkcGtSW+bjB/mb5KHxcppSpMSXogDxaRd0Rkvoj8UHCpyCBrlc43Q497Ydnb8PPrRIYG8uyQjmw9fJpnv9zo7+iUUjWUrz2QHwc+AboAp4C9XhZVXvo+CuddDV8+CJu/omfLWIZf1IS3F25nybZUf0enlKqBxJd27HbroVeMMY9XfEi+SUlJMcuWLfN3GBXnzGl4+0o4uh1Gfk56TFuuevFHctyGeff2pk6wr3X/Sil1logsN8akFNzu62MiF7C4fENSxQoKh2EfQEgkTBlE2OFfee6mjuw7nsE/5q73d3RKqRqmJGMT3VqBcShvIpNg5OdW09Mp19CVjYzpncy0pbv5duMhf0enlKpBfHpMBCAizwH9gO+AgoPmGGPMo+UbWvFq/GMiTyf3wZTBcHIvZ256n0FzAjiWfob59/UmKizI39EppaqRoh4T+VpncBNWL+QArERwpsAuxhiTWB6B+qpWJQOAtEMw9Ro4uo0dl7/O5Z8Gcdl58bx2c1ccDvF3dEqpaqKsdQbPArOBOGNMjDEmocBSqYmgVqoTDyPmQGwrms6/jddS9jN//UFe+HqTvyNTStUAviaDaOA1Y8zRslxMRFqLyCqP5aSI3Csi9UTkKxHZbL9Gl+U6NVZ4DIz4DBI70XfN33g6eQ0vLdjCnNX7/B2ZUqqa8zUZzAZ6lvVixpiNxphOxphOQFcgHfgYGAd8Y4xpCXxjrytvQqPglo+Rpj343d5/8mr0NB6csZy1e0+c81CllCqKr8ngY2CEiLwsIteLyGUFl1Jcuy+w1RizE7gGmGJvnwJcW4rz1R7BdeHmj+GiO7kq4zPecz7Bg1Pmc+hUpr8jU0pVU75WILvPsYsxxgSU6MIibwMrjDEvi8hxY0yUR9kxY0yhR0UiMgYYA9C4ceOuO3fuLMkla6a1s3B98meOZgfyYvRDjL/zjwQ7S/RPoZSqRcramuic4ycbY3weOEdEgoB9QDtjzEFfk4GnWteaqDiHNpA2dSghp3YxN+FPDB7zBOLQMQiVUoWVqTWR/ay/2KWE8QzAuis4aK8fFJEEO9AEQHtUlUR8G+rc+QM7YnpzzYGX2T5piNUUVSmlfFSSUUudIjJSRF4RkU9FJNnefp2ItCzhdYcB0zzWPwVG2O9HYFVYq5IIiaT5n2fxccwYkg59T9aLXWHl/0DnUFZK+cDXUUubAxuAl4COwEAg0i7uBzzk6wVFJMw+ZpbH5olAPxHZbJdN9PV86ixHQABX3TGRJ5Je59esBJj9Z5g6OG8qTaWUKoqvdwYvAalAM6AP4Nnl9Tugt68XNMak2x3XTnhsSzXG9DXGtLRfy9SfoTYLdgYwftR1vNb0JR7KHs2Z3Stg0sXw0/PgyvZ3eEqpKsrXZNAHeNIYc4TCcyEfABLKMyhVNsHOAF695QL2Jg+lZ9rT7IntCV9PgNcvgU1f6qMjpVQhviaDbCCwiLIE4GT5hKPKS0hgAK/f0pXWLVvSa+coFnZ9AbJPw/s3wVv9Ydv3/g5RKVWF+JoMvgbGiUgdj21GRJzAn4F55R6ZKrOQwADeGJ5Czxax3LwonpkXfwJXvwAn91p1CVMGwe6l/g5TKVUF+JoM/gY0AjYBb2A9KhoH/Ao0B/6vQqJTZZabEHokx3L/zPW8cqoX5q7lcOVEOLQB3uoH790E23/Ux0dK1WK+9jPYgdWK6D2gE9acx62x7gi6GmN0DuQqLCQwgDdHpDC4YyLPfrmR+2dtJCtlDNy9yppvec9SmHI1vHoR/PImZKX5O2SlVCUrsgeyiPTG6hhWJb8ZtAdyyRlj+M+CLfz7q010bRLN67d0JbZOMGRnwNpZsPR12P8rBEdAp9/DBbdBbEm7kCilqrISD0chIi7gImNMlXyorMmg9Oau3s/9M1YREx7MW7emcF6DCKvAGNizDJZOhnUfgzsbGnWH9tdD22uhbn3/Bq6UKrPSJAM30F2TQc20es9x/jh1GWmZObw0rDN92xT4ok87BCv/C2tmwqF1IA5o0gPa3wBtBltzKyilqh1NBqqQAycyuW3qL6zbd5KxvZO5r19L7yOeHvoN1s2CtTMhdQtIADTtAcl9oUVfqN8eRKfeVKo6KG0yeBzY5ssFjDFTyxRhCWkyKB8ZZ1w89tk6PvhlN20SInjhd51o3aCu952NgQNrrMSwab51xwBQpwEkX2YlhuZ9IDy2ssJXSpVQaZOBr0o8n0FZaTIoX1+vP8i4Was5mZHD365ozeiezXA4zvFr/+Q+2LoAtnwD276FjGPW9nrJ0KgbNLrQeo07D3RIbaWqhNImg0sBn75xjTGnyxRhCWkyKH+paVk8OGsN89cfpFuzejx3U0caRof5drDbBftWwY4frY5su3+G9CNWWXAkJHWBhPOhfgdo0B5iWkKAs+I+jFLKK60zUD4xxjBj+R4e/2w9AA9c2ZrfX9gYZ0AJf9kbA0e3nU0Me5fB4Y3gOmOVBwRD/HlWcohtCTEtrNfopuAMLt8PpZTKo8lAlcjuo+n8feZqFm1NpVX9OjxydTt6tixjXYArG45sggNr4eAa+3UdnPaYiEccENXYSg7RTSGyEUQ1gsjG1mt4vD5yUqoMNBmoEjPG8OW6gzz1+QZ2HU3n8jb1eXhgG5rGhpfvhTJPWK2UUrfar1vgyGY4vgsyj+ffNyAI6ibYS32r8rquvdSJh7BYqwI7LAYCQ8s3TqVqgDLNgVwVaTKoPJnZLt5euJ1XFmzhjMvNyB7NuPOyFkSEFDWQbTnKOgXHd8OJ3VZyOL4LTu2HUwesJe0gZBUxaG5guNUfIiwGQqIgJBJC7deQKOt9cAQE14WgOtZr7hIUDs4QbTKrahxNBqrMDp3M5NkvN/LRij3UDXYy4uKmjOzRjHrhQf4N7MxpOzEcgvRUq+I6PRVOp9rrqdYdRsbxs69uHyb6EQcEhllLUJiVXAJDrcUZAoEh4Aw9++oMsupCnCEe7+3XgCAICLQX+73DXnc4z746Aq2KdYe9SAA4Auz1gPzbNFGpUtBkoMrN2r0neHnBFuatO0BoYAC/79aYP/ZqToPIEH+H5htjrPGYMo9bg/JlnYIzp6zXvPU0yE6HM+nWPBBn0q317HTIzoScDOs1Ox1yMu3ljPVaaP6niiJnk4Q4rCQhDqtOJe99bpnDSh557+2Fgtvw2CYe+xR4X+gVj+TkpazYbQXf2+v5zulZ5GW/0mwv9rz5Cip4/1K4+G6rVV4pFJUMtG2fKrH2SZG8dktXNh88xaTvt/Luoh1MXbyDIV0bMqZ3Ms3Ku06hvIlYv/SDfGw2WxLGgDsHcrKsllO5r65s+9XjvTsbXDnWqzvH2u7OKbC47MVeN66z24xHmXFb187dZtzWe+P2KLPfu12A8V5mjF1mb8vbzxR4z9ltbns/z225r+faVvA9nm89ywsVlnF7oX+4IjYXdUw57V9aBevSyoHeGagy2300ndd/2Mr0ZXs4k+Pm4uQYhl3YmP7t6nsf3kIp5Tf6mEhVuEOnMvlw6W4++GU3e49nEB0WyA1dGjL0wsa0iK9z7hMopSqcJgNVadxuw49bjvDB0l18tf4gOW5DSpNorj4/gSvbJ1SfugWlaiBNBsovDp/K4qPle5i1Yg+bD1nzJKU0iWZAhwQGtG9AYpT2BVCqMmkyUH63+eApvlh7gM/X7Oe3A6cA6Ngoiktbx9GrZRwdG0aWfNgLpVSJaDJQVcr2I6f5Yu1+vlx7gNV7T2AM1A1xcnFyDL1axtGrZSxNYqp4qySlqiFNBqrKOnb6DAu3HuGnzUf4cfMR9h7PACAhMoQuTaLp2jiark2iaZsYQaDeOShVJtrPQFVZ0eFBXH1+Ilefn4gxhm1HTvPT5iMs23mMFTuPMXf1fgCCnQ46NoyiY6NI2iZG0DYhkuZx4ZoglCoHemegqrwDJzJZsesYy3day/r9JzmTY829FOR00Lp+XdomRNC6QV1axNeheVw4iZGh556cR6laSB8TqRojx+Vm25HTrN93kvX7T7J+30nW7TvBsfSz4w2FBDpoHmslhuZxdWhcL4xG0aE0qhdG/YgQAjRRqFpKk4Gq0YwxHE7LYtvh02w9nMbWQ6fZdiSNrYfT2HMsI98oAYEBQlKUlRgSIkNoEBFCg8hQGkQGUz/CWo8OC9I7C1UjaZ2BqtFEhPi6IcTXDaF785h8ZVk5LvYey2D3sQx2H01nz7EMdh9LZ8/RdDYeOMXhtKxCQ8o4HUJMnSBiwoOJrRtMbJ0g4uoEUy88iKiwQKLCgogOCyI6LJDIsECiQoMIcmrdhaq+Kj0ZiEgU8CbQHmv0plHAFcAfgcP2bg8ZYz6v7NhUzRTsDKB5XB2ax3kfEiPH5eZwWhYHTmRay8lMDp/KIjXtDEfSsjiSlsXWQ2kcScsiy66r8CYk0EFESCB1Q5xEhAbmva8T7CTcXuoEB1jvg5yEBQUQFuQkNMhBaGDuegDBgQGEBDoICnAgOky1qiT+uDN4EZhnjBkiIkFAGFYyeN4Y8y8/xKNqOWeAg4TIUBIii+8NbYwh/YyLY+lnOJ6ezfH0bPu9tX4qK4eTGdmczMzmVGYOx9PPsOtoOmlZOZzOyiH9jKtEcYlAiNNKDCGBAQQ7HQQ5HQQ7A+xXaz0o4OxrYICDQKcQGGCtOwMEp8NBYIC1zRlgvXc6HDgdQoBDcAbYrw4HAQ4hwAEBDgcBIva6tc1hrzs8tjvE2p67TYS8fcQuC7DLxT6HQF65eBwv9mfWjO3KcAAACDFJREFUBOgflZoMRCQC6A3cCmCMOQOc0X98VR2ISN4v/IbRJT/e7TakZ7s4nZVDWlYOGWdcpJ9xkX7Gep+Rba1nZltLVo7bfu/OWz+T4yYrx8UZl5usbDdpWTlkZbvJdrvJdrnJzjFku9yccVn75rgNLnf1qxcUwU4OVsIRrA35EoldnjsFQ+6LeJTjsY947Ji/PPdY8XjvPSl5bvI8Lv+65z5n14r9lhOvb4tMjE9d14ELm9Ur7owlVtl3Bs2xHgW9IyIdgeXAPXbZnSIyHFgG3G+MOVbwYBEZA4wBaNy4ceVErFQ5cTiEOsHWY6P6lXhdt9uQ47aSRI7LkO1247K3uVyGHLc7r9ztBpcxuNxuXG5w2cnEZQzuAu/dxtrXGGu722Bvt/cx1t1U7r5uY5UZYz0fzntvDC43GM6WkXs8ueex3pN7rNtg8NiO5zQLJm/2gNzyvPPm7VfwmPznyb+v/d7LnAtn9zP51gsfWzTPRjz59ivmoPDg8h8avlJbE4lICrAE/r+9u42RqyygOP4/aTVawBAEau02Fg0pkEZLg4jWGKU0VCUthi+Imkb8YIwiGo2CJGogMfgSFAPRKOqSsKkxK0ZiorSpEjDRSq2FtlReAkinlm3xBYwNgcrxw73Lzs7Obtlpmedu5vySzcy9M7N7dmZnz9w7d56HVba3SroReAa4CXiK6te/Dlhk+/KZvleOJoqImL3pjibq9+EPLaBle2u9PAqstD1m+3+2XwB+CJzb51wREQOtr2Vg+0lgr6Rl9arVwAOSFrVd7QPArn7miogYdCWOJroCGKmPJHoU+CjwXUkrqHYTPQ58vECuiIiB1fcysL0D6Nxf9ZF+54iIiAn5yGRERKQMIiIiZRAREaQMIiKCOTyEtaSDwN96vPnJVB9ya7KmZ2x6Pmh+xqbng2Q8FpqW7w22T+lcOWfL4GhI2tbtE3hN0vSMTc8Hzc/Y9HyQjMdC0/ONy26iiIhIGURExOCWwQ9KB3gJmp6x6fmg+Rmbng+S8Vhoej5gQN8ziIiIyQZ1yyAiItqkDCIiYvDKQNJaSQ9KekTSVaXztJO0RNLvJO2RtFvSlUe+VRmS5kn6i6Rflc7SSdKJkkYl/bW+L99eOlMnSZ+tH+NdkjZKelUDMv1Y0gFJu9rWnSRps6SH69MeJvx8WfN9s36c75f0C0knlso3Xca2yz4vyZJOLpHtSAaqDCTNA24G3gucBXxQ0lllU01ymGrKzzOB84BPNixfuyuBPaVDTONG4De2zwDeQsNySloMfBo4x/ZyYB5wadlUAAwDazvWXQVssX06sKVeLmWYqfk2A8ttvxl4CLi636E6DDM1I5KWAGuAJ/od6KUaqDKgmkHtEduP2n4O+CmwvnCmF9neb3t7ff4/VP/EFpdNNZWkIeD9wC2ls3SS9BrgXcCPAGw/Z/vfZVN1NR94taT5wALg74XzYPtu4J8dq9cDt9bnbwUu7muoNt3y2d5k+3C9+EdgqO/BJufpdh8CfBv4AjNPh1zUoJXBYmBv23KLBv6zBZC0FDgb2DrzNYv4DtUf9gulg3TxRuAg8JN6N9Ytko4rHaqd7X3At6heJe4Hnra9qWyqaS20vR+qFyvAqYXzzORy4NelQ3SStA7YZ/u+0llmMmhloC7rGtfUko4Hfg58xvYzpfO0k3QRcMD2n0tnmcZ8YCXwPdtnA/+l7K6NKer97uuB04DXA8dJ+nDZVHObpGuodrOOlM7STtIC4Brgy6WzHMmglUELWNK2PEQDNs/bSXoFVRGM2L69dJ4uVgHrJD1OtZvtfEm3lY00SQto2R7fohqlKocmuQB4zPZB288DtwPvKJxpOmPjc5TXpwcK55lC0gbgIuBDbt4Hp95EVfr31c+ZIWC7pNcVTdXFoJXBvcDpkk6r52C+FLijcKYXSRLVvu49tm8onacb21fbHrK9lOr++63txryqtf0ksFfSsnrVauCBgpG6eQI4T9KC+jFfTcPe5G5zB7ChPr8B+GXBLFNIWgt8EVhn+1DpPJ1s77R9qu2l9XOmBays/04bZaDKoH6j6VPAnVRPvp/Z3l021SSrqOaDPl/SjvrrfaVDzUFXACOS7gdWAF8rnGeSeqtlFNgO7KR6HhYfskDSRuAPwDJJLUkfA64H1kh6mOpomOsblu8m4ARgc/18+X6pfDNknBMyHEVERAzWlkFERHSXMoiIiJRBRESkDCIigpRBRESQMoiYQtJXJT01zWXDkrb1O1PEyy1lEBERKYOIJqrni3hl6RwxOFIGEUdB0gpJWyQdkvQvSSOSFrZd/u56QpPlHbe7S9Jo2/KwpG2SLpa0G3gWeFv/fpMYdPNLB4hoqnqugSmr2y4/BbiLamiTy4DjqYZr2CzpnHrOjNlYCnwDuBYYAx6bfeqI3qQMIrp7LfD8NJeND9/9ufr0wvGhxiU9RDUHxSXAxh5+5gW2d8zydhFHLWUQ0d3TVENNd/oKsKg+fy6wqX3OCdt/qocqfiezL4N9KYIoJWUQ0d1h21MOIZX0DybKYBHQbdTbMeCkHn7mWA+3iTgm8gZyRO+mmwZyIRPz4D5bn3YeGdStLDKEcBSTMojo3VbgQkknjK+Q9FaqN4J/X69q1adntl1nCTA++U5EI2Q3UUTvbgA+Adwp6etMHE20k2rqUmy3JN0LXCfpENULsC8xseUQ0QjZMojoke2DwHuodgVtBG4G7gHWdBxWehnVVJe3Uc26di3wYH/TRswsM51FRES2DCIiImUQERGkDCIigpRBRESQMoiICFIGERFByiAiIkgZREQE8H8e8EGcDVvF/gAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"N = 50\n",
"t = np.linspace(0,15,N)\n",
"\n",
"temp_analytical = Ta + (T0-60)*exp(-K*t)\n",
"temp = np.array([55,58,60,65,66,67])\n",
"temp_numerical=np.zeros(len(t))\n",
"temp_numerical[0]= T0\n",
"\n",
"for i in range(0,len(t)):\n",
" if i == 1:\n",
" temp_numerical[i] = temp_numerical[i-1] + -K*(temp_numerical[i-1]-temp[3])*(t[i]-t[i-1])\n",
" if i == 2: \n",
" temp_numerical[i] = temp_numerical[i-1] + -K*(temp_numerical[i-1]-temp[4])*(t[i]-t[i-1])\n",
" if i >= 3:\n",
" temp_numerical[i] = temp_numerical[i-1] + -K*(temp_numerical[i-1]-temp[5])*(t[i]-t[i-1])\n",
"\n",
"plt.plot(t,temp_analytical, label='Analytical')\n",
"plt.plot(t,temp_numerical, label='Non-linear Numerical')\n",
"plt.legend(loc='best', prop={'size': 15})\n",
"plt.xlabel('Hour',size = 15)\n",
"plt.ylabel('Temperature (deg F)', size=15)\n",
"plt.title('Temperature of body', size=15);"
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"b) Time: 0 hr 42 min 39 s\n",
" Time of death: 10:17:21 am\n"
]
}
],
"source": [
"time = np.log((temp_death - 60)/(T0-60))/K*3600\n",
"\n",
"hours = math.floor(time/3600)\n",
"minutes = math.floor((time/3600-hours)*60)\n",
"seconds = int(round(((time/3600-hours)*60-minutes)*60,0))\n",
"print('b) Time:',hours,'hr',minutes,'min',seconds,'s')\n",
"print(\" Time of death: 10:17:21 am\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"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.3"
}
},
"nbformat": 4,
"nbformat_minor": 4
}