Skip to content
Permalink
9a4273205b
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
387 lines (387 sloc) 68.5 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": 138,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"0.6111111111111112\n"
]
}
],
"source": [
"Ta = 65\n",
"\n",
"t0 = 0\n",
"T0 = 85\n",
"\n",
"deltat = 2\n",
"Tt = 74\n",
"\n",
"dT_dt = (Tt-T0)/deltat\n",
"\n",
"K = -(dT_dt/(Tt-Ta))\n",
"print(K)\n"
]
},
{
"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": 139,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"-5.5\n"
]
},
{
"data": {
"text/plain": [
"0.6111111111111112"
]
},
"execution_count": 139,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def measure_K(Temp_t1,Temp_t2,Temp_ambient,delta_t):\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",
" dT_dt = (Temp_t2 - Temp_t1) / delta_t\n",
" K = -(dT_dt/(Temp_t2 - Temp_ambient))\n",
" print(dT_dt)\n",
" return K\n",
"\n",
"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 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": 140,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"[0. 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1. 1.1 1.2 1.3 1.4 1.5 1.6 1.7\n",
" 1.8 1.9]\n",
"step is: 0.1\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3RU1fbA8e9OTyAQCEgLEBAUpIQuYgPhIagUlfpQERXs9fnsiPp+dkVFFARFsGEBUbFTpEnvXZAeek2AJKTt3x93wBBSJmVmUvZnrSxmztyyM2vYuXPuOfuIqmKMMab08PN1AMYYY7zLEr8xxpQylviNMaaUscRvjDGljCV+Y4wpZQJ8HYA7KlWqpNHR0b4OwxhjipVly5YdUtXKmduLReKPjo5m6dKlvg7DGGOKFRHZkVW7dfUYY0wpY4nfGGNKGUv8xhhTyhSLPn5jTMmRkpJCbGwsSUlJvg6lxAgJCSEqKorAwEC3trfEb4zxqtjYWMLDw4mOjkZEfB1OsaeqHD58mNjYWOrUqePWPiW3q2f11/BWY3guwvl39de+jsgYAyQlJREZGWlJv5CICJGRkXn6BlUyr/hXfw1TH4CUROd53C7nOUDTPr6LyxgDYEm/kOX1/SyZV/wzXvgn6Z+Wkui0G2NMKVcyE39cbN7ajTEmj8aPH899992X6zZ79uw58/yOO+5g/fr1eT7XrFmzuO666/K8X3ZKZuIvH5W3dmOM8YDMif/DDz/koosu8mFEjpKZ+Ds+C4GhZzWpwl+RV/koIGNMUdKzZ09atmxJo0aNGDNmDABly5bl6aefJiYmhrZt27J//34Apk6dysUXX0zz5s3p1KnTmfbTjh8/Tp06dUhJSQEgPj6e6OhovvnmG5YuXcqAAQNo1qwZiYmJtG/f/kz5mV9//ZUWLVoQExNDx44dAVi8eDHt2rWjefPmtGvXjr/++ssjv3/JvLl7+gbujBcgLpb0ctU5eCKN2lu+YMHMTlxyVXffxmeMAeD5qetYvye+UI95UfVyDOvWKMdtxo0bR8WKFUlMTKR169bceOONnDx5krZt2/Liiy/y2GOPMXbsWJ555hkuu+wyFi5ciIjw4Ycf8tprr/Hmm2+eOVZ4eDjt27fnp59+omfPnnz55ZfceOON9O7dm/fee4833niDVq1anXX+gwcPMnjwYObMmUOdOnU4cuQIAA0aNGDOnDkEBAQwffp0nnrqKSZPnlyo7w+U1MQPTvJ3/QHwA8KP7efQux1pPHsIi0PDaXNJB9/GZ4zxmREjRjBlyhQAdu3axebNmwkKCjrTj96yZUumTZsGOPMO+vbty969e0lOTs5yrPwdd9zBa6+9Rs+ePfn4448ZO3ZsjudfuHAhV1xxxZljVaxYEYC4uDgGDhzI5s2bEZEz3yIKW8lN/JmERVQh7c4fSRrdibq/DmR58GRatGjt67CMKdVyuzL3hFmzZjF9+nQWLFhAWFgY7du3JykpicDAwDPDIv39/UlNTQXg/vvv55FHHqF79+7MmjWL55577pxjXnrppWzfvp3Zs2eTlpZG48aNc4xBVbMcgjl06FA6dOjAlClT2L59O+3bty/w75uVktnHn43w86IJGvQDAX5Q5fv+rF6/ztchGWO8LC4ujgoVKhAWFsbGjRtZuHBhrtvXqFEDgAkTJmS73S233EL//v0ZNGjQmbbw8HCOHz9+zraXXHIJs2fPZtu2bQBnunoynmv8+PF5+r3yolQlfoDyNS8ifcAkIuQEYV/3YcOW7b4OyRjjRV26dCE1NZWmTZsydOhQ2rZtm+P2zz33HL179+byyy+nUqVK2W43YMAAjh49Sv/+/c+03Xrrrdx1111nbu6eVrlyZcaMGcMNN9xATEwMffv2BeCxxx7jySef5NJLLyUtLa2Av2n2RFU9dvDC0qpVKy3shVgOrZ1B+KR+bKIWIbf/SP2a1Qr1+MaYrG3YsIGGDRv6OoxCN2nSJL7//ns+/fRTn5w/q/dVRJapaqvM25a6K/7TKjXuSNx1Y2jIVo6O6832/Ud8HZIxppi6//77eeKJJxg6dKivQ3FLqbm5m5XzWl3PvoRjtJn5ELM+6E/gfZOoUTHc12EZY4qZd99919ch5IlHr/hF5GERWScia0VkooiEZHjtURFREcm+08wLql4xiL2XDKN9+kJWvn8rB+ISc9/JGGOKMY8lfhGpATwAtFLVxoA/0M/1Wk3gX8BOT50/L6pd/Qh7Yh7g2tTpzHrvbo6cTPZ1SMYY4zGe7uMPAEJFJAAIA04XrXgLeAwoMneWq/d8gX0X3kyf5Cl8/95/iUv0zMQJY4zxNY8lflXdDbyBc1W/F4hT1d9FpDuwW1VX5bS/iAwRkaUisvTgwYOeCjPjCanadwT7a13HoITxfPH+C5w8ler58xpjjJd5squnAtADqANUB8qIyC3A08Czue2vqmNUtZWqtqpcubKnwjybnx9VBo7nQJUruDN+BGmv1kVtBS9jTCH74YcfeOWVV/K1b3R0NIcOHSrQ+T3Z1dMJ2KaqB1U1BfgWGITzh2CViGwHooDlIlLVg3HkjX8g513cFxGhXHocgv6zgpclf2O8r4Qto5qamkr37t154oknfBaDJxP/TqCtiISJU5SiI/Ctqp6nqtGqGg3EAi1UdZ8H48i72a84CT8jW8HLGO87vYxq3C4oxIuw7du307BhQwYPHkyjRo3o3LnzOWWTDx06RHR0NOCUT+jZsyfdunWjTp06jBw5kuHDh9O8eXPatm17puTCli1b6NKlCy1btuTyyy9n48aNgDOD95FHHqFDhw48/vjjZy3isn//fq6//npiYmKIiYlh/vz5QNalowuLx8bxq+oiEZkELAdSgRVA4UbvKdms1KVxsdhKocYUol+egH1rsn89dgmknTq7LSURvr8PlmVTN6dqE+iaezfK5s2bmThxImPHjqVPnz65lj9eu3YtK1asICkpiXr16vHqq6+yYsUKHn74YT755BMeeughhgwZwujRo6lfvz6LFi3innvuYebMmQBs2rSJ6dOn4+/vf1YdngceeIArr7ySKVOmkJaWxokTJ4CsS0dHRkbm+nu5w6MTuFR1GDAsh9ejPXn+fCsf5brCOFucXwTByWmEBvn7IChjSqHMST+39jyoU6cOzZo1A5wyzNu3b89x+w4dOhAeHk54eDjly5enW7duADRp0oTVq1dz4sQJ5s+fT+/evc/sc+rUP3H27t0bf/9zc8fMmTP55JNPAKcqaPny5YGsS0cXi8RfbHV81vk6mWHBdkUITjvBi2M+5akhNxEWZG+dMQWW25X5W42zvAijfE0Y9FOBTh0cHHzmsb+/P4mJiQQEBJCeng5AUlJSttv7+fmdee7n50dqairp6elERESwcuXKLM9XpkwZt2PLrnR0YSm1tXpy1LQPdBvhfLgQKF8T6fIqaWWq8tjBJ/nfB5/bUE9jvCGLZVQJDHXaPSA6Opply5YBTtG1vChXrhx16tThm2++AZya+6tW5ThqHYCOHTsyatQoANLS0oiPj89z6ei8ssSfnaZ94OG18Nwx59+2d1J2yC/4l6nIE4ee5PkxX3DCkr8xnpXFRRjdRvyzvGohe/TRRxk1ahTt2rXL15DJzz//nI8++oiYmBgaNWrE999/n+s+77zzDn/88QdNmjShZcuWrFu3Ls+lo/Oq1JZlzrdjO0kYczXJJ+N4oeIrPH9nf8JDAn0dlTHFRkkty+xrVpbZkyJqETb4F4LDyvHMkSd55oOviE+y8g7GmOLDEn9+VIgmdPDPhIWVYdiRJ3h69FfEJVjyN8YUD5b486tiXULu+IWwsFCGHX2KJz74mmMJVtXTGHcUhy7m4iSv76cl/oKIPJ+Q238mPDSIF449yWOjJ1lJZ2NyERISwuHDhy35FxJV5fDhw4SEhOS+sYvd3C0MB/8i+aOuHEtM58lyL/PanTcQWTY49/2MKYVSUlKIjY0t1HHppV1ISAhRUVEEBp490CS7m7uW+AvLgQ0kf3QNR5OUx8Nf5vU7b6ByuCV/Y4zv2KgeTzuvIUG3/UjFEHj5+NM88sF3HDhuVzTGmKLHEn9hqtKIwEFTqRScxojjjyDDG1o9f2NMkWOJv7BVbULgZfcTISeorEesnr8xpsixxO8Jy8afW77Z6vkbY4oIS/yekEM9f2OM8TVL/J5QPirL5oNanvl/F2ytTGOMKShL/J6QRSlZRQiXRIaPn8hv64rWSpPGmNLFEr8nZFXP/+qXCIqoxqeBL/Lx55/x9dIsFpcwxhgvsGWkPKVpn3Nqhvs3up7gCd355MirDP42mbiEfgy+oq6PAjTGlFZ2xe9N5arhd9svBFRpwEdBb7L01wm8/ttGq1lijPEqjyZ+EXlYRNaJyFoRmSgiISLyuohsFJHVIjJFRCI8GUORU6YSfgOn4h/VgveD3mX37Ak8/d1a0tIt+RtjvMNjiV9EagAPAK1UtTHgD/QDpgGNVbUpsAl40lMxFFmhEcjNU/CrcylvBY2CpR/z4JcrSE5N93VkxphSwNNdPQFAqIgEAGHAHlX9XVVPL1a7EMh67GNJF1wW+fc3SP3OvBT4EVXXfcjgT5aSkGzr+BpjPMtjiV9VdwNvADuBvUCcqv6eabPbgF+y2l9EhojIUhFZevDgQU+F6VuBIdD3M7ioJ88Efk7zraO5+cNFtpqXMcajPNnVUwHoAdQBqgNlROSmDK8/DaQCn2e1v6qOUdVWqtqqcuXKngrT9wKCoNc4aDaAhwIm03Xv+/T9YD4H4q2ypzHGMzzZ1dMJ2KaqB1U1BfgWaAcgIgOB64ABakNawM8fuo+ENkO4w/9Hbj32Lr1H/cnOwwm+jswYUwJ5chz/TqCtiIQBiUBHYKmIdAEeB65UVctsp/n5QdfXIKgM/ea9Rf3EnQS+exjlMFI+ypkNnGlegDHG5IfHEr+qLhKRScBynC6dFcAYYB0QDEwTEYCFqnqXp+IoVkSg03NwZBst13/3T/vpss5gyd8YU2AenbmrqsOAYZma63nynCXC7mXntp0u62yJ3xhTQDZztyjKoayz3RIxxhSUJf6iKJuyzvs1gqemrCE1zSZ6GWPyzxJ/UZRlWWeI8E9mzZLZ3D5hKSdO2UQvY0z+WOIvirIq6/yvFwgJj2RK2EuwZSZ9Ri9gv431N8bkgxSHPuNWrVrp0qVLfR2G78Xvhc97k35gA0+l3cmc0I58PKgNF1YN93VkxpgiSESWqWqrzO12xV+clKsGg37CL7odr/i9x79Tv6XXqD/505ZzNMbkgSX+4iakPAyYDI17cV/aZ/xf8KcMGreQSctsIXdjjHtsBa7iKCAIbhgL4VXpsWAkVcvHccs3txN7NIEHO9bHNTHOGGOyZIm/uPLzg6tfhPBqXPz70/xcMZ7rp99L7NFEXrq+CUEB9mXOGJO1XLODiESKSDcRuVNEbhGRFmKXlEVHu/vgxo+om7SemRVfYd6yVdw2fgnxSVba2RiTtWwTv4hcLiI/46yYdT1OeeUWwP8Ba0VkqIiU9U6YJkdNeiE3TaZS6kFmRrzI4W0r6T1qAXuOJfo6MmNMEZTtcE4ReQt4V1W3ZvFaENAdQFUneTRCbDin2/atgc96kZqcwIenOtFD5lDVqnsaU2plN5zTxvGXNMd2wthO6Mn9nNUfFxjqTAqz5G9MqZHncfwi8lGGxzdlt50pYiJqgX8A59yESUlEZ7zgi4iMMUVMTjd3W2R4/IinAzGFKH5P1u1xsSQmp3k3FmNMkZNT4i/6fUAma9lU99yrFeg1er7d9DWmlMsp8UeJyHDXTd7Tj8/8eCtAkw9ZVPcEqBgaQPrhbXQfOY9lO474IDBjTFGQU+J/EmeZxLUZHmf8MUVVFtU9ufw/hJDMj6HPcknAX/Qbs5Cvl+zydaTGGB/IcVSPiIQAqqqnvBfSuWxUTyE5vAW+6IMe3cHYiAd5aU8Lbru0Dk9d04AAf5vpa0xJk59RPY8DPwBTReQxV9v7eTzpwyKyTkTWishEEQkRkYoiMk1ENrv+rZDXX8bkU+T5cMd0pHY7hhx5g89r/8zHf25h0PglxCXYTF9jSoucLvMaq2pnVe0MNHe1uZ2kRaQG8ADQSlUbA/5AP+AJYIaq1gdmuJ4bbwmtADdNhla3cen+z5hXexyrt+6mx3vz+PvAcV9HZ4zxgpwSf1lXXZ7WQH5LMwQAoSISAIQBe4AewATX6xOAnvk8tskv/0C4djh0eZUaB2ax4LzXKJO0j+vfm88fGw/4OjpjjIfllPj/A9yEc5X+oKttQvabn01VdwNvADuBvUCcqv4OVFHVva5t9gLnZbW/iAwRkaUisvTgwYPunta4SwTa3gX//oawhN38EPwsncrv4rYJS/hg9haKw4xuY0z+eKxkg6vvfjLQFzgGfANMAkaqakSG7Y6qao5dSHZz18MObHRu+p7Yz8eV/ssL2xtyffMavHxDE0IC/X0dnTEmn/Jzc/cPEblbRKpnag8QkStE5CMRGZTDOTsB21T1oKqmAN8C7YD9IlLNdaxqgPUt+Np5DWDwTKR6c27b9z++vuAPpqyIpe+YhbaguzElUE4LsVwL3AFMcd2oPQKEAsE4N2XfU9WcLsN3Am1FJAxIBDoCS4GTwEDgFde/3xf0lzCFoEwluOV7mPoQbVaNZUW1RSQe2E3l4Uc4VaY6wVc/ZwXejCkhsk38qpoAjABGiEgwTl98oqq6tbK3qi4SkUnAciAVWAGMwblR/LWI3I7zx6F3wX4FU2gCgqHn+5CaRIV131LBVekt+ORuUr67nwBALPkbU+y5tfSiawJXnqd5quowYFim5lM4V/+mKBKB2CXnNAemJ3H4h6cJuuB6wkMCfRCYMaaw2HRNc6642CybK6QcpMfIP/lrn433N6Y4s8RvzpVNdU8NDCUpMYGe7/3Jdyt2ezkoY0xhcSvxi0iUiHRwPQ4WkTKeDcv4VFbVPf0C8E9NYHalV2lfNZmHvlrJM9+t4VSq1fc3prjJNfGLyG04NXs+dDXVxkbilGxZVffsOQr6fk7g0S28f/IhXmx2hM8W7qTP6AXstvr+xhQruU7gEpGVQBtgkao2d7WtVtWmXogPsAlcRcqhzfDlADi8mY2N/0Ov1a0I9Pfj7X7NufKCyr6OzhiTQZ4ncGWQpKrJGQ7kD+cu6WpKiUr1YfAMaNiNBmteZ2G9T6kdrtz68WLenr6J9HQr9WBMUedO4v/TVZY5xNXP/xXwo2fDMkVacDj0ngD/eoGyW3/m28ChDGmYztvTN3Pr+CUcOZmc+zGMMT7jTuJ/DDgObMQp1jYDeNqTQZliQAQufRBunoJfwiGe2H0PE9odZOGWw3R7dx4rdx3zdYTGmGzkmPhd3TrjVHWUql6vqj1dj9O9FJ8p6uq2hyGzkcjzuXL5g8xt/Sd+mkbv0fP5dOEOq/JpTBGUY+JX1TSgmojYVE2TvYiaMOhXaH4TVVa+y8xq7/Fk5YV0+PkqeL4C6cMbweqvfR2lMcbFnZINW4G5IvI9ToE1AFR1hMeiMsVPYAh0Hwk1WhH44yMMYhbi51ztS3ws6d8/4FxlWK0fY3zOnT7+g8A0nBW0Kmf4MeZsItBqEJSthHB2F49fWiLHf3rWRv0YUwTkesWvqkO9EYgpQU5kvWJamaR93Dp+CW/2jqFyeLCXgzLGnObOzN1pIvJ75h9vBGeKqWxq/SSFVGLR1sN0fWcOszfZcprG+Io7XT3PAENdPy/iDOtc5cmgTDGXVa0fICz5GLM77qRiWCADxy3mpZ83kJxqA8SM8TZ3unoWZWqaLSKzPRSPKQlO38Cd8YJT4rl8lDPmf+OPVJ39GD9fdAMv1bqTMXO2smDLYUb0b06dSlb3zxhvcadWT7kMT/2AlsAoVb3Ak4FlZLV6Soj0dJg3HP54CSJqMb/FG9w9M53UtHT+17MxN7TIuovIGJM/BanVsw5Y6/p3Bc6s3cGFG54pFfz84IpH4dafIC2Zdn/0Y84Vf9GoWjke+XoVD3+1kuNJKb6O0pgSz50r/kBVTcnUFqCqqR6NLAO74i+BEo7Ad/fApl/QC6/hg4j/8Nqc/dSsGMaIfs2JqRnh6wiNKfYKcsWfuY8fYLEbJ7xQRFZm+IkXkYdEpJmILHS1LRWRNu78AqaECasI/SfC1S8jm6dx14aB/HR9EKlpyo2j5jN69hYb82+Mh2Sb+EXkPBGJAUJFpImINHX9XIYzmStHqvqXqjZT1WY49wUSgCnAa8DzrvZnXc9NaSQCl9wDt/8O/gE0/KUvM9oso3PDyrzyy0YGfryYA/FJvo7SmBInp1E91wK3AVHA+xnaj+MM7cyLjsAWVd0hIgqcvmFcHtiTx2OZkqZGC7hzDkx9iJA5/8d7df9kaYt2VF/3AZWGHyYxrBqhXZ63cg/GFBJ3+vj7qGqBKmyJyDhguaqOFJGGwG84i7n4Ae1UdUcW+wwBhgDUqlWr5Y4d52xiShpVWD4BfvwPZLqFlCzBpF77DmGt+vsoOGOKn+z6+HNN/K6drwYaASGn21T1JTdPHIRzVd9IVfeLyAhgtqpOFpE+wBBV7ZTTMezmbinzxgVwYv85zXupzN8DFnB5fSsVZYw78n1zV0TeBwYCjwChwE1AvTycuyvO1f7p/8kDgW9dj7/BWc/XmH+cOJBlcxUOcfNHixn63VoSkr02qMyYEsedUT2Xqeq/gcOugm0X4/T7u6s/MDHD8z3Ala7HVwGb83AsUxpkU+tHQspz+6XRfLZoB9e8M5dlO454OTBjSga3Fls//a+IVHU9j3bn4CISBvyLf67wwZn89aaIrAJewtWPb8wZWdX6ET8k6RhDjw1l0oA6pKYrvUcv4OVfNnAqNc03cRpTTLmzEMvPIhIBvAGsBNKACe4cXFUTgMhMbfNwhncak7Wsav10fBZOxcNvz9Byd1emd3mT57dewAeztzJr40GG942hUfXyvo3bmGIix5u7IuIHtD5dqE1EQoFQVfXqd2y7uWvOOPQ3TBkCu5dBkz7Mrf84/5m6nSMnk3mwY33ubn8+Af7ufJE1puTL181d16Lq72R4nujtpG/MWSrVg9t+h/ZPwdrJXD69OzNvELo2qcab0zZx4+gF/H3ghK+jNKZIc+fSaJqI9PB4JMa4yz8A2j8Od0yDwFDKfnUD71b4hvf7XsSOwye5dsRcxs3bZiUfjMmGOxO4juLMsD0FJOJMvFJVrej58BzW1WOylZwA04fB4jFQuQFHrh7Jo/Ng5sYDtK1bkdd7xVCzYq4VRowpkfI9gUtE/LNqV1WvDaWwxG9y9fd0+O5eSDiMtn+SRYdDqLVyOFU5TEJoVcK6voBfjJV8MKVLQWfu9gPqqupLIhIFVFHVZR6IM0uW+I1bEo7Ajw/D+u9A/ED/WdYxiWCOdXqDqpfd4sMAjfGugszcHQl0AG52NSUAows3PGMKQVhF6D0eQiuelfQBQjhF2rTnGTFjs63za0o9d27utlPVO3FN5HKN6gnyaFTG5JcIJB7N8qXqcpjh0zbR7d15rNiZ9TbGlAbuJP4U13h+BRCRSMAumUzRlUPJh49uaUF8Ugo3jJrP81PXcfKU1fwxpY87if89YDJQWUSeB+YBr3o0KmMKIpuSDyQdo+PCQUy7pRo3t63Nx39up/Nbc5i96aBv4jTGR3JN/Kr6CfAMTsmGI0BvVf3S04EZk29N+0C3EVC+JiDOv9d/AD1HwYENlB3XnhcipzFpSGtCAv0YOG4xj3y1kiMnk30duTFe4e6onqbAZTjdPX+q6mpPB5aRjeoxheb4fvj5P7BhKlSL4dS17/Le+hDen7WFcqGBDOt2Ed1jqiMivo7UmAIryKiep3HKKlfHKcf8hYg8WfghGuMF4VWg72fQ5xOI30vwuKt4JOAbfrq3NTUrhvHglyu5bfwSdh9L9HWkxniMOxO4NgAtXZU2T5daXqaqDb0QH2BX/MZDEo7Ab0/BqolQuQFp3d5lws7KvP7bX/gJPNalATe3rY2fn139m+Ip31f8wA7OLt8cAGwtrMCM8ZmwinD9aBgwCU6dwH9cZ247MZZp97WiRe0KDPthHb1Gz2f9nnhfR2pMoXLniv9boDXOAukKdMYZ2bMfQFUf8XCMdsVvPC8pHqY/B0s/ggrRaKMbSFg2kdDEfezVSBbWuZd/9bufciGBvo7UGLcVpFbP7Tm9rqofFTC2XFniN16zfR58cyucPHuIZ4IG8ZL/3bTqdic9mtnNX1M8ZJf4c12ByxuJ3ZgiI/oy8A8+pzlMknmAL2jz1cV8uWQn/+vRmPpVwn0QoDEF586oni4iskREDojIERE5KiK2GIspueJ3Z9lcOf0QL17fmA17j9P1nbm8/PMGm/lriiV3bu6OBO4EagCVgUquf3MkIheKyMoMP/Ei8pDrtftF5C8RWScirxXkFzCm0GVX8gFlQNJXzHywLTe2iOKDOVvp+OZsflq9F3fmwxhTVLiT+GOBlaqaoqppp39y20lV/1LVZqraDGdx9QRgioh0AHoATVW1Ec6MYGOKjqxKPgSEQI1W8MeLRH7agVdbHGHy3e2oWCaIe79Yzi3jFrP1oC35aIoHdxL/Y8BUEfmviDxw+ieP5+kIbFHVHcDdwCuqegpAVQ/k8VjGeFZWJR+6vwuDZ8BN3zolnz/pQcslj/LDwPN5vnsjVu46Rpe35/LGb3+RmOy1NYqMyRd3RvX8AqQAa8hQlVNVh7p9EpFxwHJVHSkiK4HvgS44pZ4fVdUlOe1vo3pMkZKSBH++DXOHQ0AwXPUMBxvczMu/buLbFbupERHKsG4X8a+LqtjoH+NTBRnOuUxVWxbgxEHAHqCRqu4XkbXATOBBnPkBX+Gs7qWZ9hsCDAGoVatWyx07duQ3BGM84/AW+PlR2DITqsXAdW+x6FQ0Q79fy6b9J7iqwXk8160RtSJtzV/jGwWZuTtDRK4qwLm74lzt73c9jwW+VcdinG8RlTLvpKpjVLWVqraqXDnXe8nGeF/k+U7XT6+PneJvYzty8foX+WlwE565tiGLth6m01uzee3XjZyw0T+mCHHniv8oUB7n5mwyIICqakW3TiDyJfCbqn7sen4XUF1VnxWRC4AZQK3MV/wZWVePKfKS4mHWyzAi7dgAABaySURBVLBoNIRFQsNupP71O/7Hd7M7PZJR/gNoes0d9GpZE3+r/WO8pCBdPf5ZtbszssdV0G0XTldOnKstCBgHNMP5Q/Koqs7M6TiW+E2xsXc1fH0LHN12VnMSwTyWfDt/V7mGodddxCXnR/ooQFOa5Lurx5XgewOPux5Xw0nauVLVBFWNPJ30XW3JqnqTqjZW1Ra5JX1jipVqTSE95ZzmEE7xSvnviEtMof/Yhdz56VJ2HD7pgwCNcW/m7kigA3CzqykBGO3JoIwp1uKynvkblriHGQ9fxqOdL2Du5kN0Gj6bl37eQHzSuX8ojPEkd27utlPVO3GGXqKqR4Agj0ZlTHGWzcxfgJCP2nNf9B5mPdqens1qMHbuVjq8PovPFu4gNS092/2MKUzuJP4UEfHDKcmMiESSYTy/MSaTrGb+BoZCmzsh+Th80p3zfr6d1zuWY+p9l3H+eWV55ru1XDtiHnM328LvxvOyTfwicrpy53vAZKCyiDyPU4v/VS/EZkzxlNXM324j4JrX4N4lcNVQ2PIHvNeGxuuH89XARowa0IKElFRu/mgxt49fwhYr/2A8KNtRPSKyXFVbuB43AjrhDOWcrqprvReijeoxJVD8XpjxAqz6AsqcBx2fJalxP8Yv2MnImX+TlJLGzZfU5sGO9YkIs55Vkz95Hs4pIitUtbnHI3ODJX5TYsUug1+fgNjFUK0ZdHmFgxVbMHzaJr5aspPwkEDu7XA+t1wSTUhgliOrjclWfhJ/LDA8uwOqaravFTZL/KZEU4U1k2D6MGctgEY3QM02JM97l8ATe9idHsnYwJto1OUObmhRgwB/d27NGZO/xL8XGIXTvXMOVX2+UCPMgSV+Uyokn4Q/R8CcN0HPHuJ5egLY+kpd+O/VF9LZCsAZN+Qn8Z/p4/c1S/ymVHmzARzfe05zYmh1rvUfxdZDJ2lRK4LHuzTg4ro2A9hkLz8zd+1ywhhfOL4vy+bQxL38/tDlvHxDE3YfS6TvmIUM+ngxG/bGezlAU9zllPg7ei0KY8w/sp0ApgR8ci39q+5l1qMdeLxLA5btOMo1I+by8Fcr2XUkwathmuIr28TvmqFrjPG27CaANRsAR7bCuM6EfnsLdzdKZe5jVzHkirr8vGYvV705i+d+WMfhE6d8E7cpNnKtzlkUWB+/KXVWf+2M84+Ldb4BdHzWmRiWfBIWvg/z3oGUk9D8Jmj/JHu1AiNmbObrpbGEBPgx+Iq63HF5XcoGB+R+LlNi5bssc1Fgid+YTE4egrlvwuKx4BcAbe+GSx/k7+MBvPn7X/yydh+RZYK4/6p69L+4FsEBNgegNLLEb0xJdHQ7zHwR1nwNoRXg8keh9R2s3JfEq79sZMHWw1QrH8I9HerRp1WU/QEoZSzxG1OS7V0F05+HLTOc2kD1OqF/T4O43Rz0q8z/JfViaXgn7r2qHr1b1iQowCaBlQaW+I0pDbbOgu/vh7idZzWn+YcwIux+3jnYnBoRodzboR69WkbZH4ASriCLrRtjiou67cmqarp/WhIP+U3kk9vacF65YJ6asoYOb8ziy8U7SbF1AEodS/zGlDTZrAAmcbFcEbiRb+9ux/hBrakUHswT3zp/AL5aYn8AShNL/MaUNNlNABM/mHAdMqEb7UP+5rt72vHxoNZElgni8clr6PjmbL5eusv+AJQCHkv8InKhiKzM8BMvIg9leP1REVERqeSpGIwplbKbANZ9JHR5FQ5tgo+7Ip/2pEPoNr6791LG3dqK8qGBPDZpNZ2Gz+abpbtsKcgSzCs3d0XEH9gNXKyqO0SkJvAh0ABoqaqHctrfbu4ak0fZTQADSE6ApePgz7fh5EE4vyN0eAqt0ZIZGw7w9oxNrN0dT3RkGPddVZ8ezaoTaKWgiyWfjuoRkc7AMFW91PV8EvA/4HuglSV+Y3wg+SQs+RD+fAcSDkP9q6HDk2i1ZkzfcIC3p29i3Z54oiqEcucVdendqqYtBlPM+DrxjwOWq+pIEekOdFTVB0VkO9kkfhEZAgwBqFWrVssdO3Z4PE5jSqVTJ2DxGJg/AhKPwoXXQI2W6LLxEBd7Zh7A/LCO3HZZNDe1rU25kEBfR23c4LPELyJBwB6gEXAc+APorKpxOSX+jOyK3xgvSIqHRR84pSBSE896Kc0/lDHlH+TVPU0JDwnglktqM+jSOlQqG+yjYI07fDmOvyvO1f5+4HygDrDKlfSjgOUiUtULcRhjchJSDq78L4RVPOcl/7RE7k77nKn3Xcbl9Svx/qwtXPrKTIZ9v5bYo1YOurjxRum+/sBEAFVdA5x3+gV3r/iNMV4Uvyfr9rhdNEleyfv/voIth07ywewtfLF4J58v2kn3ZtW5+8rzqV8l3LuxmnzxaFePiIQBu4C6qhqXxevbsa4eY4qWtxpD3K5z28UPNB1qtILLH4ELurIn/hQfzt3GxMU7SUxJo/NFVbinQz2a1YzwftzmHFarxxjjntVfw9QHICVDP39gKFzzJqSdgnlvw7EdULmh8weg0Q0cSUpn/J/bGD9/O/FJqVxaL5J72tej3fmRtii8D1niN8a4L6d5AGmpsG6KcxP44AaIqA2XPQQx/+ZEegBfLNrBh3O3ceD4KWKiyjP4irp0aVSVAJsL4HWW+I0xhSs9HTb9CnPfgN3LoGwVuOQ+aDWI5PU/kfLbc4Qm7WNPeiQfBt1E1JUD6dO6pg0F9SJL/MYYz1CFbXOcbwDbZkNAKKSnQHrqmU2SCOax5NuZGdievq1rcmu7aGpWDPNh0KWDJX5jjOfFLoPx154zDwAguWwN/hv1OT+t3ku6Kl0bV+OOy+vQvFYFHwRaOmSX+G0lZmNM4YlqCalJWb4UdGIP7/RrzuNdGjBhwXa+WLSTn9bspWXtCtxxWR06N6qKv5/dCPYGu9tijClc2ZWFRmFcV6rvm8mTV1/Awic78ly3izh4/BR3f76c9m/8wbh52zhxKjWb/U1hsa4eY0zhymo4aEAoNLwOdi5yloWsWBfa3gPN/k1aQBjT1u/no3lbWbL9KOHBAfS/uBYD20VTIyI0+/OYXFkfvzHGe7IbDpqWCht+gAUjnZFAIRHQ6jZoMwTKVWPlrmN8NG8bP6/ZC0DXxlW5tV00LWtXsPkA+WCJ3xhTdKjCrsWw4F3Y8CP4BUCTXnDJvVC1CUcWfob/zP8RnryfPemRfFF2ILXa30qPZjUIDbLS0O6yxG+MKZqObHWqgi7/FFJOQqUGcHQrpCWf2eT0cNBZQe3p3aomN7etTXSlMj4MuniwxG+MKdoSj8KyCU4Xkaad8/KpMjX4T43P+HXtPlLTlSsvqMwtl9Sm/YXn2WigbFjiN8YUD89FAFnlJYHnjnEgPomJi3fx+aIdHDh+iqgKodzUtjZ9W9WkQpkgb0dbpFniN8YUD9lVBwWo9y9oMxjqdSJFhd/X7eeTBdtZtO0IQQF+dI+pzi2X1KZplFUHBUv8xpjiIsvhoCFQv7NzQ/jEPqcwXOvbofnNEFaRv/Yd55MF25myYjcJyWnE1Ixg4CW1uaZJtVK9TrAlfmNM8ZHtcNAU2DDVWSR+x5/OH4TGvZxvAdWbEZ+UwrfLYvlk4Q62HjxJxTJB3NiiBv3a1OL8ymV9/Vt5nSV+Y0zJsn8dLB4Lq7+ClASIag2tB0Ojnuj67zn16zCCE/ayRyN5NaUP+2p1p1+bmqXqW4AlfmNMyZR4DFZNdL4FHP4bAstCWtJZ1UFT/EJ4JeBuPopvTbmQAK5v7nwLaFitnA8D9zxL/MaYki09HbbNgon9sywUp+VrsqD7LL5asotf1u4jOTWdmJoR9Gtdk24x1SkbXPJqVlriN8aUDtkOBwUeXg/la3D0ZDJTVuzmyyU72bT/BGWC/OkWU51+bWoRE1W+xJSH8HriF5ELga8yNNUFngVqAN2AZGALMEhVj+V0LEv8xhi35TQcVPygXidoMRAuuBr1C2D5zmN8uXgnP67eS2JKGg2qhtO/TS16NqtB+bDivVqYT6/4RcQf2A1cDFwIzFTVVBF5FUBVH89pf0v8xhi3ZbdY/FVDIeEIrPjMGRJatgo0GwAtboaKdTmelMIPq/bw5eJdrNkdR3CAH9c2qUavllG0rRuJXzGcHezrxN8ZGKaql2Zqvx7opaoDctrfEr8xJk9yWyx+8++wfILzr6ZDnSuh5UBocB0EBLNr9njKznuJ8skH2KPOmsFlW/+bG1tGUacY1QjydeIfByxX1ZGZ2qcCX6nqZ1nsMwQYAlCrVq2WO3bs8HicxphSJn4PrPgcVnwCx3ZCaEWo3hJ2zD3rBvEpCebx5Nv5Lu0yWtSK4MaWUVzXtDrlQ4t2V5DPEr+IBAF7gEaquj9D+9NAK+AGzSUIu+I3xnhUejps/QOWfwLrv8tyk7TwKD5s9QOTl8eyaf8JggL8+NdFVejVIorL61ciwL/oLWjoy8TfA7hXVTtnaBsI3AV0VNWE3I5hid8Y4zU5FYl79ggqwtrd8UxeHsv3K3dzNCGFSmWD6dmsOje2jCpScwN8mfi/BH5T1Y9dz7sAw4ErVfWgO8ewxG+M8ZqcRgWVrwlN+0JMf6hUj+TUdP746wDfLo9l5sYDpKQpF1Urx40to+jRrDqVygZ7N/ZMfJL4RSQM2AXUVdU4V9vfQDBw2LXZQlW9K6fjWOI3xnhNdmsGNx8AR7fDlpnODeGo1hDTDxrdAGEVOXIymamr9jB5eSyrY+MI8BPaX1iZHs1q0KlhFZ+sHGYTuIwxxl05jQo6vs95fdVEOLAe/IPgwq7Ot4B6ncA/kE37j7Ph9w9ps2UkVfQQe6nEjOp3EnXlQC6rV5mgAO/cD7DEb4wxhUkV9q2GlRNhzTeQcAjCKkGT3lCmEsx946xvDYkE8XjyHcwJ6UDXxtXoHlOdi+tU9Oj8AEv8xhjjKWkpsHma8y1g069nrRecUWJYdZ6o9QXT1u8nITmNKuWC6da0Ot2bVadJjcIvFWGJ3xhjvCHhCLxWJ5sXneUjE5JTmb7hAD+s3MPsTc5N4ejIMLrHOH8E6p0XXiihWOI3xhhvyW35yCa9oMG1EBxOXEIKv6zdyw+r9rBg62FU4aJq5ejerDrdYqpTIyI032FY4jfGGG/JcmRQMNRp79wQjtvlrB52wdXOCmL1O0NgCAfik/hxtfNHYOUup3bl6Jta0qVx1XyFYYnfGGO8KbuRQenpELsY1kxyZgmfPAhB4dDwOuePQN0rwT+Qg/M/JXj2/xF+aj+SeWSRmyzxG2NMUZOWCtvnwJrJzlrCp+IgLBKqNIadCyHt1D/bBoZCtxF5Sv6W+I0xpihLPeWMDFo7GdZNIcuyEeVrwsNr3T5kdom/6FUVMsaY0igg2Onu6f1x9tvExRbKqSzxG2NMUVM+Km/teWSJ3xhjipqOzzp9+hkFhjrthcASvzHGFDVN+zg3csvXBMT5N483dnMSUChHMcYYU7ia9im0RJ+ZXfEbY0wpY4nfGGNKGUv8xhhTyljiN8aYUsYSvzHGlDLFomSDiBwEduRz90rAoUIMp7BZfAVj8RWMxVdwRTnG2qpaOXNjsUj8BSEiS7OqVVFUWHwFY/EVjMVXcMUhxsysq8cYY0oZS/zGGFPKlIbEP8bXAeTC4isYi69gLL6CKw4xnqXE9/EbY4w5W2m44jfGGJOBJX5jjCllinXiF5EuIvKXiPwtIk9k8bqIyAjX66tFpIW7+3opvgGuuFaLyHwRicnw2nYRWSMiK0XEI+tOuhFfexGJc8WwUkSedXdfL8X33wyxrRWRNBGp6HrNo++fiIwTkQMikuU6eEXgs5dbfL7+7OUWn68/e7nF57PPXqFQ1WL5A/gDW4C6QBCwCrgo0zbXAL8AArQFFrm7r5fiawdUcD3uejo+1/PtQCUfv3/tgR/zs6834su0fTdgphffvyuAFsDabF732WfPzfh89tlzMz6fffbcic+Xn73C+CnOV/xtgL9VdauqJgNfAj0ybdMD+EQdC4EIEanm5r4ej09V56vqUdfThUDhrKtWSPF5aF9PxdcfmFjIMWRLVecAR3LYxJefvVzj8/Fnz533LztF4v3LxKufvcJQnBN/DWBXhuexrjZ3tnFnX2/El9HtOFeIpynwu4gsE5EhhRxbXuK7RERWicgvItIoj/t6Iz5EJAzoAkzO0Ozp9y83vvzs5ZW3P3vu8tVnz21F9LOXq+K8Apdk0ZZ5bGp227izb0G5fQ4R6YDzn++yDM2XquoeETkPmCYiG11XId6MbzlOrY8TInIN8B1Q3819vRHfad2AP1U14xWap9+/3Pjys+c2H3323OHLz15eFMXPXq6K8xV/LFAzw/MoYI+b27izrzfiQ0SaAh8CPVT18Ol2Vd3j+vcAMAXnK65X41PVeFU94Xr8MxAoIpXc2dcb8WXQj0xftb3w/uXGl589t/jws5crH3/28qIofvZy5+ubDPn9wfm2shWowz83eRpl2uZazr7Bttjdfb0UXy3gb6BdpvYyQHiGx/OBLj6Iryr/TPJrA+x0vZdF4v1zbVcepy+2jDffP9exo8n+5qTPPntuxuezz56b8fnss+dOfL7+7BX0p9h29ahqqojcB/yGc6d/nKquE5G7XK+PBn7GGV3xN5AADMppXx/E9ywQCbwvIgCp6lT5qwJMcbUFAF+o6q8+iK8XcLeIpAKJQD91PtFF5f0DuB74XVVPZtjd4++fiEzEGXlSSURigWFAYIbYfPbZczM+n3323IzPZ589N+MDH332CoOVbDDGmFKmOPfxG2OMyQdL/MYYU8pY4jfGmFLGEr8xxpQylviNMaaUscRvSgQRiRCRezI8ry4ikzx0rp6nq0WKyHgR6eWJ82Rz7utE5Hlvnc+UTJb4TUkRAZxJ/Kq6R1U9lZAfA9730LEBEBH/bF76CejuqhFjTL5Y4jclxSvA+a4a6K+LSPTpWuoicquIfCciU0Vkm4jcJyKPiMgKEVmYoY76+SLyq6u41lwRaZD5JCJyAXBKVQ9laL7CVdN+6+mrf3G87qrVvkZE+rra24vIjxmON1JEbnU93i4iz4rIPKC3iDwgIuvFqZn/JYBrEtMs4LrCfwtNaVFsZ+4ak8kTQGNVbQYgItGZXm8MNAdCcGbTPq6qzUXkLeAW4G2cRbPvUtXNInIxzlX9VZmOcylOAbGMquEUOWsA/ABMAm4AmgExQCVgiYi4U6grSVUvc/0Oe4A6qnpKRCIybLMUuBz42o3jGXMOS/ymtPhDVY8Dx0UkDpjqal8DNBWRsjiLk3zjmm4PEJzFcaoBBzO1faeq6cB6EaniarsMmKiqacB+EZkNtAbic4nzqwyPVwOfi8h3ONUpTzsAVM/lOMZkyxK/KS1OZXicnuF5Os7/Az/g2OlvDDlIxCnOld2xJdO/maVydhdrSKbXM9Z9uRZnJajuwFARaaSqqa59EnOJ05hsWR+/KSmOA+H53VlV44FtItIbzvTRx2Sx6QagnhuHnAP0FRF/EamMk8AXAzuAi0QkWETKAx2z2llE/ICaqvoHzs3kCKCs6+ULgCzXgjXGHZb4TYmgTj35P103U1/P52EGALeLyCpgHVkv6TcHaC4Z+oOyMQWnq2YVMBN4TFX3qeounL751cDnwIps9vcHPhORNa5t3lLVY67XOuCM7jEmX6w6pzF5JCLvAFNVdboPzl0Fp9Rvlt8UjHGHJX5j8siVfC9W1R98cO7WQIqqrvT2uU3JYYnfGGNKGevjN8aYUsYSvzHGlDKW+I0xppSxxG+MMaWMJX5jjCll/h8+ds3XDVmupwAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"import math\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"\n",
"\n",
"Temp_t1 = 85\n",
"#Temp_t2 = ????\n",
"Temp_ambient = 65\n",
"delta_t = 2\n",
"\n",
"n = .1\n",
"t = np.arange(0,delta_t,n)\n",
"print(t)\n",
"print('step is: ',n)\n",
"\n",
"\n",
"def T_analytical(Temp_t1,Temp_ambient,K,t):\n",
" Temp_t2 = Temp_ambient + (Temp_t1-Temp_ambient)*(math.e)**(-K*t)\n",
" return Temp_t2\n",
"\n",
"\n",
"def T_numericalfunc(Temp_t1,Temp_ambient,K,step): \n",
" T_numerical = np.zeros(len(t))\n",
" T_numerical[0] = Temp_t1\n",
" for i in range(1,len(t)):\n",
" T_numerical[i] = T_numerical[i-1] + step*(-K*(T_numerical[i-1]-Temp_ambient))\n",
" return T_numerical\n",
"\n",
"\n",
"\n",
"plt.plot(t,T_analytical(Temp_t1,Temp_ambient,K,t),'-',label='analytical');\n",
"plt.plot(t,T_numericalfunc(Temp_t1,Temp_ambient,K,n),'o-',label='numerical');\n",
"plt.legend();\n",
"plt.xlabel('time (hours)');\n",
"plt.ylabel('Temperature (ºF)');\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"b) As t approaches infinity, the temperature of the corpse will approach the ambient temperature as it reaches an equilibrium."
]
},
{
"cell_type": "code",
"execution_count": 141,
"metadata": {},
"outputs": [],
"source": [
"#c. #To solve for the time when the temperature was 98.6ºF, we can solve the analytical solution \n",
"# substituting 98.6 for T(t) shown below\n",
"#T(t) = Ta + (T(0)-Ta)e^(-Kt)\n",
"#98.6 = 65 +(85-65)e^(-Kt)\n",
"#solving for t yields -.8489 hours, which means that the time of death was 0.8489 hours before 11am, or\n",
"# at approximately 10:09am"
]
},
{
"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": 142,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEWCAYAAABhffzLAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3gVZfbA8e9JgUDoJDQDhEBo0puIQUFcZNGfXcGysuqKuop9Feuu7KrYG+66CHZAEQU7YpeywBIISCf0UJMgEEr6+f0xk2yIIbmE3Exy7/k8Tx7unXre4eZk7jsz5xVVxRhjTPAI8ToAY4wxlcsSvzHGBBlL/MYYE2Qs8RtjTJCxxG+MMUHGEr8xxgQZS/xVlIhsEZFzjjNvoIisq+yYTMUTkd0ikuB1HAAislBErvE6jgL+PDZVra2VzRJ/BRORH0XkVxGp6a99qOpcVe1QEdsq4w/M1SJyyP05KiL5Rd4fqoj9e0VEholIciXv830Rebgy91kV+HKsg/XYeMUSfwUSkVhgIKDABZ4GUwFUdYqq1lHVOsDvgZ0F791pVZKIhIiIXz/bIhLmz+1XlX2awGSJv2JdCywE3gJGFZ0hIm+JyD9F5Cv3jHm+iDQTkRfdbwhrRaRnse31FZHV7vw3RSTC3dYgEUkpsu0WIvKRiKSKyGYRub3IvL+JyHQReUdEMkRklYj0cee9C7QCPnNjuu9EGywiLUXkExFJE5FNInJzkXnjRWSKiHzgbj9JRNqIyF/d5beIyOAiyy8Ukb+LSKKIHHDbVL/I/IEiskhE9ovIUhE5o9i640RkEXAEaCEiN7nHNUNEkkXkenfZxsBMIK7IN5jGxc86i5+pul0P94rIKuBgWe0vdpxuBy4FHnH392GR2X1FZKXb5ikiUqPo/kXkERHZA/zLnX6riGwUkXQR+VhEmrrTO4pIbrH9FnZpiEiYiLzsrrdRRG4vvjzQ1l3noIh8KSINi2zrUvfzuF9EvhWReHd6hIioiMQUWfZ9EXn4eMf6ZI+Nu97FIrLCjWeuiHQu6di7y54nIhvcZZ8vNq+jON/U97m/Q2+LSF133iMiMqXY8q+LyPjj7ataUFX7qaAfIBn4M9AbyAGaFpn3FpDmzosAvgc24/yxCAX+AfxQZPktwEqgJdAImA/8w503CEhxX4cAicCjQA0gDtgEnOvO/xuQCQx39/MksLDYfs7xoW2F+ywyLRT4Bbjf3Xd7YBtwljt/PE4SHgyEAR+4bb7XfT8GWFNkewuBrUBHoA7wGTDJnRcLpAPnuG0eDqQCDYusuwnoAIS7278AaAOIu95R4FR3+WFAcrH2vA88XOT9McsAu4H/Ai2AWmW1v4RjeMz2i2xzPtAUiMb5DP2xyP5zgXHu9mu57d4NdMP5HE0EvnGX7wjkFtv+QuAa9/WdwHKgOdAY+Lno8u6y64C2QCSwAPibO68LkOF+DmoAjwBr3OMcgfMtN6aktpZ0rCvg2PQHduH8PoUCo4H1QFgJ224GHHY/D+HAA+5xvabIcTvbbVcz9ziMd+e1dttdx31fE/i14HNUXX/sjL+CiHMRqjUwXVUTgY3AVcUWm6mqiaqaiXMWlKmq76hqHk5SLH7GP0FVt6vqPuBx4MoSdt0XiFbVcaqaraqbgNeBkUWWmaeqX7r7eRfofpLNLZAARKjqU+6+1wNvFtv3d6r6g6rmAjOAesBz7vv3gY4iUqvI8m+q6lpVPQT8tUibRwEfq+q3qpqvql8Cq4GhRdadpKrrVDVHVXNV9VNV3ayOb4Gf3JhPxguqulNVj/rYfl+3uUdVU4EvgR5F5mUBf3e3fxS4Gpioqivcz9F9wBARaebDfq4AnlfVXaqaDjxdwjKvq+pGVT2M8/9VEMuVOJ/fH1U1G3gCiAL6nGBbT9Txjs1NOL8fiaqap6oTcZJy7xK2cQHwX/fzkIPT7n0FM93P2/fuMd4NvAic5c7bCiwBLnYX/z9gs6quqvimVh5L/BVnFDBHVdPc91Mp1t0D7Cny+mgJ74v3m28v8norzplmca1xujX2F/wAD+KcJRXYXeT1ESBCKqa/uDUQW2zfd+OcNRUo3sZUdU+d3PfgnF0WKN7m2m53T2vgmmL76sOxx6TouojIBSKy2P0Kvx/nrC6qfE0tcR++tN8Xxf9/in4OdrvJqkALnOMCgKrux+l2OsWH/bQoFv/2EpY5XizF95sH7PBxvyfjePG0Bh4sduyjjxPPMe0uEjtQ2FX6oYjsEJGDwCSO/Zy8DRTcAXQNzslTtWYXiyqAe8Z6BRAqIgUf1JpAAxHprqrLy7nplkVetwJ2lrDMdpwzkPhy7uNkyrNuB9aqateT2EZxxdt8RFUPiMh2nDP6MaWsW9gWEYkEPgQuA75S1VwRmY3T7XPMskUcBmoXeV9SAi+63om2vzzHuvg6O3GSHgDuH8V6OIksHOczWFNVs9xFirZhFxBT5H3RY12W4vsNxUmyO4BsnK7N4sdu5XHaUJITPTbbgS9U9Tkflt2F00UFOBf/OfYPxDM4//ddVPVXERmJ0/VaYAbwkoicivMNs8TrONWJnfFXjIuAPKAzzlfRHkAnYC5OH3553SoiMSLSCOcs/oMSllkMHBSR+0WkloiEikgXEenr4z724FwXKI95ACJyp3uBL0xEuolIr3JuD+CPItJeROrgXJ8oaPPbwOUiMsRtYy339fHOrmvhJMK9QL6IXECRX36cdjdx91MgCThfRBqIyCk41yBKc6LtP5ljXWAacKP7fxwBPAV873ZR7MS57nG1e4z+zLEJbjpwlzg3FTTGudbiqw+Ai0XkTBEJB8biXHNZoqr5ONc6Cvb7f8DpRdYt6VgXd6LHZiIwRkT6iKOO+w2vdgnLfopzkfh8N/a/4Fw3K1AXOITze9QK51tbIbfb8VOcY/+je6yrNUv8FWMUTt/0NlXdXfADTMD5ZSjvN6upwByci5abOPYsBCj82vp/OH9sNuNcQJ4E1C++7HE8CTzsfl0+kUSA2wUxHBiA0w2QinPnycnc6vkuzi/YDiAfuMfd1yacOz8ew2njVuAOjvMZdrvc7sW5QJyO88f5yyKLLMf5Zd7qtr0R8AbOBcRtwOduHMdVjvZPxElA+0Xk/dK2Xco+P8f5P/sUJ9E3A/7gzssD/oRzbSQN54w+scjqE3Au2K7GuUj9Oc41BF/2uwK4Afg3TjuHABe612oAbgNG4Fz4vNjddoGSjnVxJ3RsVHU+cLsbz36cC7tXUcI3B1XdhXPd5UU39qY4/fYFHsW5XnMA59rbRyXs8m2gKwHQzQMg/+tuNcZbIrIQ54Lde17HEgxE5GKcu1cq5GHAQCYi7XH+WDRT1SNex3Oy7IzfmCAhInVFZKjbHdMKeBjnDNeUwr2ecTfwXiAkfbCLu8YEkxCcZyva41zM/JQSug/N/7jdUttwn43xOJwKY109xhgTZKyrxxhjgky16OqJiorS2NhYr8MwxphqJTExMU1Vo4tPrxaJPzY2liVLlpS9oDHGmEIisrWk6dbVY4wxQcYSvzHGBBlL/MYYE2SqRR9/SXJyckhJSSEzM9PrUEw1FRERQUxMDOHh4V6HYkylqraJPyUlhbp16xIbG4uIlL2CMUWoKunp6aSkpNCmTRuvwzGmUlXbrp7MzEwaN25sSd+Ui4jQuHFj+8ZoglK1TfyAJX1zUuzzY4JVtU78xhgTiLJy81iQnMZTs9ey52DFfyu1xH+SZs6ciYiwdu3aE173j3/8IzNmzPjN9CVLlnD77beXO6YnnniixOmnnXYaPXr0oFWrVkRHR9OjRw969OjBli1byr0vf/v444/LdWyNqU5UlbW7DzJp7iaufWMx3R+bw1WTFvH6z5tYueNAhe+v2l7crSqmTZtGQkIC77//Pn/7298qZJt9+vShT5/yj2H9xBNP8OCDD/5m+qJFiwB46623WLJkCRMmTCj3PipSbm4uYWElfxQ//vhjQkJC6NixY4Vsz5iqYu/BTOYlpzFvQxpzk9NIzXDGxImLjmRk31YktIuif9vG1KlZ8Z9lO+M/CYcOHWL+/PlMnjyZ99//36BBP/74I2eddRZXXHEF7du3Z+zYsUyZMoV+/frRtWtXNm7cWLjst99+y8CBA2nfvj2ff/554frnn38+AIcPH+b666+nb9++9OzZk08++QRwkvcll1zCsGHDiI+P57777gNg7NixHD16lB49enD11Vf73JavvvqK008/nV69ejFixAgOHz4MQExMDA899BD9+/enb9++LF26lKFDh9K2bVtef/31wjYMHjyYiy66iM6dO3PrrbdSUPW1tO3+/e9/54wzzmDmzJm89tpr9O3bl+7du3P55Zdz9OhR5s6dy5dffsldd91V+M0kISGBpKQkAHbv3k27du0AmDRpEiNHjuT888/n97//PQDjx4+nX79+dOvWjXHjxp3If60xFe5odh4/rtvLPz5fzbkv/Ey/J77j7unL+WHdXvrHNebpS7uxYOzZfH/PIP52wamc07mpX5I+BMgZ/2OfrWL1zoMVus3OLerx1/87tdRlZs2axbBhw2jfvj2NGjVi6dKl9OrlDLe6fPly1qxZQ6NGjYiLi+NPf/oTixcv5qWXXuKVV17hxRdfBGDLli389NNPbNy4kcGDB5OcnHzMPh5//HHOPvts3njjDfbv30+/fv0455xzAEhKSmLZsmXUrFmTDh06MGbMGMaPH8+ECRMKk6Mv9u7dy/jx4/nuu++oXbs2jz/+OC+99FLht4bY2FgWLlzImDFjuOGGG5g3bx6HDh2ie/fu3HjjjYDzbWL16tW0bNmS3/3ud3zyyScMGDCg1O1GRkYyf/58ANLT07n5ZmcM67Fjx/LWW29xyy23MHz4cC677DIuuuiiMtvxn//8h6SkJBo2bMiXX37Jtm3bWLRoEarK8OHDWbBgAQMGDPD5uBhzMvLzlVU7DzI3OZV5G9JYsuVXsvPyqREaQt82DbmoZ0cGxkfRuXk9QkIq90aDgEj8Xpk2bRp33nknACNHjmTatGmFib9v3740b94cgLZt2zJ06FAAunbtyg8//FC4jSuuuIKQkBDi4+OJi4v7TX/2nDlz+PTTT3n22WcB5zbWbdu2ATBkyBDq13eG1u3cuTNbt26lZcuWJ9yOBQsWsHr16sKkmJ2dTUJCQuH8Cy64oDD23NxcIiMjiYyMJCQkhEOHDgHQv39/Ciqojhw5knnz5gGUut0RI0YUvl6xYgWPPvoo+/fvJyMjo/Abz4kYOnQoDRs2BJzj9tVXX9GzZ0/A+Xa2fv16S/zGr3bsP8q8DanM3ZDG/OQ0fj2SA0DHZnUZNaA1CfHR9IttRK0aoZ7GGRCJv6wzc39IT0/n+++/Z+XKlYgIeXl5iAhPP/00ADVr1ixcNiQkpPB9SEgIubm5hfOK31JY/L2q8tFHH9Ghw7HDoi5atOiYfYSGhh6z3ROhqgwbNox33y15HOmisRdvV8E+S2pHWduNjIwsfH3ttdfy1Vdf0aVLFyZNmsTChQtLXCcsLIz8/HyA39yDX3R7qsrDDz/MDTfcUOJ2jKkIGZk5LNy0rzDZb0pzujKj69ZkcMcmDIyP4ox2UTSpG+FxpMeyPv5ymjFjBtdeey1bt25ly5YtbN++nTZt2hSe6frqww8/JD8/n40bN7Jp06bfJPhzzz2XV155pbDPfNmyZWVuMzw8nJycHJ9jGDBgAD/99BObNm0CnOsKGzZsOIFWwMKFC9m2bRt5eXlMnz6dhISEE9ru4cOHadasGTk5OUydOrVwet26dcnIyCh8HxsbS2JiIkCJd0QVOPfcc5k8eXLhNYWUlBTS0tJOqE3GFJebl0/i1l956dsNXP7aAnqO+4Yb31nCB0u207JRbR4+rxNf33kmix8cwvNX9ODinjFVLulDgJzxe2HatGmMHTv2mGmXXnopU6dOPaYLoywdOnTgrLPOYs+ePbz22mtERBz7IXnkkUe488476datG6pKbGxs4UXg4xk9ejTdunWjV69eTJkypcwYmjZtyuTJkxkxYgTZ2dmAc2dQfHy8z+0YMGAA99xzD6tWrWLQoEFccMEFiIjP2x03bhz9+vWjVatWdOnSpfBs/sorr+Smm27iueeeY9asWfzlL39hxIgRvPnmmwwePPi48QwfPpy1a9fSv39/wPkDMnXqVKKionxukzGqytb0I8xNTmPehlQWbEwnIzMXEejSoj6jz4wjIT6K3q0bUjPM2+6bE1Etxtzt06ePFh+IZc2aNXTq1MmjiExR3377LRMmTGDWrFleh3LC7HNkitt/JJsFG9OZuyGNecmpbN93FIBTGtRiYHwUCfFRDGgbRaPIGh5HWjYRSVTV39wbbmf8xpiglp2bz9JtvxbeT/9Lyn7yFerUDKN/XGNuHBhHQrso2kRFBkyZD0v85qSdc845hbeYGlPVqSrJew+5Z/RpLNyUzpHsPEJDhO4x9RlzdjwD46Po3rIB4aGBeRm0Wid+VQ2Yv8Cm8lWHbk5TMdIOZTE/Oc1J9hvS2O3Wv4ltXJtLe8WQEB/F6W0bUy8iOMZmqLaJPyIigvT0dCvNbMqloB5/8YvpJjBk5uTx3y37nO6bDWms3uU84Fm/VjgJ7Zx++oR2UbRsVNvjSL1RbRN/TEwMKSkppKameh2KqaYKRuAy1V9+vrJm90Hmud03izfvIys3n/BQoXfrhvzl3A4ktIuiyyn1Ca3kp2SrIr8mfhFpAEwCugAKXK+q/xGRMcBtQC7whared6LbDg8Pt5GTjAliuw9kMndDKvOSnadk0w45twzHN6nD1ae1ZmB8FP3aNCLST/VuqjN/H5GXgNmqepmI1ABqi8hg4EKgm6pmiUgTP8dgjAkAh7NyWbQ5vbCffsNep1xIVJ0abvdNNAntomhW37rvyuK3xC8i9YAzgT8CqGo2kC0itwDjVTXLnb7XXzEYY6qvvHzllx0HCsshLN32Kzl5Ss2wEPq1acTlfWJIaBdNx2Z1K73IWXXnzzP+OCAVeFNEugOJwB1Ae2CgiDwOZAL3qup/i68sIqOB0QCtWrXyY5jGmKpi+74jhQ9OzU9O58BRp/RI5+b1uP6MNgyMj6ZPbEMiwqvPU7JVkT8TfxjQCxijqotE5CVgrDu9IdAf6AtMF5E4LXZvnapOBCaC8+SuH+M0xnjkwNEc/rMxnXlu6eIt6UcAaFYvgqGdm5LgFjmLqlOzjC2ZE+HPxJ8CpKjqIvf9DJzEnwJ87Cb6xSKSD0ThfDswxgSwnLx8krbvd/vpU1mecoC8fKV2jVD6xzVm1IBYBsZH0Ta6jt2m7Ud+S/yqultEtotIB1VdBwwBVgMbgbOBH0WkPVADsLKJxgQgVWVz2mHmuvfTL9yUzqGsXEIEusY04M+D2pLQLoqerRpSIywwn5Ktivx9V88YYIp7R88m4DrgMPCGiKwEsoFRxbt5jDHV177D2cx3x5Kdl5zGjv1OkbOWjWpxQY8WDGznFDmrXzs4npKtivya+FU1CShp1PBr/LlfY0zlycrNI3HLr27p4jRW7jyAKtSNCOOMtlHcMqgtA+OjaN04suyNmUphTzYYY06IqrJ+zyHmurdZLt68j6M5eYSFCL1aNeSuc9qTEB9Ft1PqExagRc6qO0v8xpgy7c3IdLpu3O6bvRlZAMRFRzKib0sS2kXRv21j6thTstWC/S8ZY37jaHYei7fsY+56pyTC2t3O8JcNa4eTEB/NQLfQWYsGtTyO1JSHJX5jDPn5yupdB/l5g3M//ZItv5Kdl0+N0BD6tmnI/cM6MjA+is7N69lTsgHAEr8xQWrH/qOF5RAWbExn32GnyFnHZnUZNaA1CfHR9IttRK0a9pRsoLHEb0yQyMjMYeGmfU6yT05jU+phAJrUrcmgDtEMdJ+SbVLXipwFOkv8xgSo3Lx8lqcccC/IprJs235y85WI8BD6xzXmqn6tGBgfTfum9pRssLHEb0wA2Zp+mJ/dcggLNqaTkZmLCHQ9pT6jz4wjIT6K3q0bUjPMum+CmSV+Y6qx/UeyWbAxvbCi5fZ9zlOypzSoxXldmztFztpG0TCyhseRmqrEEr8x1Yyq8taCLcxK2skvKfvJV6hTM4zT2zbmxoFxJLSLok1UpHXfmOOyxG9MNfNhYgqPfbaabjH1GXN2PAPjo+jesgHh9pSs8ZElfmOqkQ17MvjrJ6sY0LYx795wmg0cbsrFThGMqSYyc/K4beoyatcI5cURPSzpm3KzM35jqonHPlvNuj0ZvH19P5rUs3vtTfnZGb8x1cDnK3YybfE2bj6rLWe1j/Y6HFPNWeI3porbln6EBz76hV6tGnDP0PZeh2MCgCV+Y6qw7Nx8bpu2FBF4+cqedueOqRDWx29MFfbU7LWsSDnAa9f0IqZhba/DMQHCTh+MqaK+W7OHyfM2M+r01gzr0tzrcEwAscRvTBW068BR7vlwOZ2b1+OB4Z28DscEGEv8xlQxuXn53DEtiezcfCZc1ZOIcCuoZiqW9fEbU8W8/N0GFm/ZxwsjuhMXXcfrcEwAsjN+Y6qQ+clpvPJDMpf1juHinjFeh2MClCV+Y6qI1Iws7vwgibioSMZdeKrX4ZgA5tfELyINRGSGiKwVkTUicnqRefeKiIpIlD9jMKY6yM9X7p6exMGjObx6dS9q17BeWOM//v50vQTMVtXLRKQGUBtARFoCvwO2+Xn/xlQL//55E3M3pPH4xV3o2Kye1+GYAOe3M34RqQecCUwGUNVsVd3vzn4BuA9Qf+3fmOoices+np2zjvO6Nueqfq28DscEAX929cQBqcCbIrJMRCaJSKSIXADsUNXlpa0sIqNFZImILElNTfVjmMZ4Z/+RbG6flkSLBhE8eWlXGzXLVAp/Jv4woBfwL1XtCRwG/gY8BDxa1sqqOlFV+6hqn+hoq0ZoAo+qct+MFezNyGTClb2oFxHudUgmSPgz8acAKaq6yH0/A+cPQRtguYhsAWKApSLSzI9xGFMlvfOfrcxZvYf7h3Wke8sGXodjgojfEr+q7ga2i0gHd9IQYKmqNlHVWFWNxfnj0Mtd1pigsXLHAR7/Yg1nd2zCDQltvA7HBBl/39UzBpji3tGzCbjOz/szpso7lJXLmGnLaBRZg2cv7279+qbS+TXxq2oS0KeU+bH+3L8xVY2q8vDMX9iafphpN/anUWQNr0MyQcie3DWmEn2YmMKspJ3ceU57Totr7HU4JkhZ4jemkiTvzeCvn6xiQNvG3Dq4ndfhmCBWZlePiDQGBgAtgKPASmCZqtrDV8b4KDMnj1unLKN2jVBeHNGD0BDr1zfeOW7iF5GBwANAMyAJ2AtEACOB1iLyPvCCqh6qjECNqc7Gfb6adXsyeOu6vjSpF+F1OCbIlXbGfwlwm6puKj7DvUvnAmAYzv35xpjj+HzFTqYu2sZNZ8UxqEMTr8Mx5viJX1XvKmVeNpbwjSnTtvQjPPDRL/Rs1YB7h3YoewVjKsFxL+6KyOQir6+pnHCMCRzZufncNm0pIvDKlT0JD7V7KUzVUNonsVeR13f7OxBjAs3Ts9eyIuUAT1/WjZiGtb0Ox5hCpSV+u2vHmHL6bs0eJs3bzLWnt2ZYl+Zeh2PMMUq7uBsjIs8DUuR1IVW1bwHGlGDXgaPc8+FyOjevx4PDO3kdjjG/UVrif+A4r40xx5Gbl88d05LIzs1nwlU9iQgP9TokY36jtLt6JotIhPNSsyoxJmOqrZe/28DiLft4YUR34qLreB2OMSUq7a6e+4FPgc9E5D532j8rKzBjqpsFyWm88kMyl/WO4eKeMV6HY8xxldbV00VVhwKIyDR3WkP/h2RM9ZN2KIs7PkgiLiqScRee6nU4xpSqtMRfR0R6AaGAfWc15jjy85W7PkjiwNEc3rm+H7Vr+HuYC2NOTmmf0HuA23Bu67zDnfa23yMyppr598+bmLshjccv7kKn5vW8DseYMpV2cXcTxR7cUtXZfo/ImGokceuvPDtnHed1bc5V/Vp5HY4xPint4u4PInKLiLQoNj1MRM4UkckiYkMpmqB14EgOt09bRosGETx5aVcbQtFUG6V19ZwH/AmYKSKnAPuAWkBN4DvgVVVd4v8Qjal6VJX7PlrOnoOZzLhlAPUiwr0OyRifldbVcwR4GXhZRGoCTYCjqppWWcEZU1W985+tfL1qDw8N70SPlg28DseYE+LT7QfuA1zb/RyLMdXCyh0HePyLNZzdsQk3JLTxOhxjTpjViTXmBBzKymXMtGU0iqzBs5d3J8SGUDTVkN1wbIyPVJVHZq1ka/phpt3Yn0aRNbwOyZhy8emMX0RiRGSw+7qmiET6uF4DEZkhImtFZI2InC4iz7jvV4jITBGxDlJTLcxITGHmsh3cMaQ9p8U19jocY8qtzMQvItfj1OyZ5E5qDXzi4/ZfAmarakegO7AG+AanHEQ3YD1W+dNUA8l7M3j0k1WcHteY285u53U4xpwUX874bwf6AwcBVHU9zh0+pRKResCZwGR3vWxV3a+qc1Q1111sIWDVrEyVlpmTx61TllG7RigvjuxBqPXrm2rOl8Sf6Q6uDoCIhOIMzlKWOCAVeFNElonIpBK6iK4HvvI5WmM8MO7z1azbk8FzV3Snab0Ir8Mx5qT5kvjnu2WZI9x+/g+Az31YLwxn3N5/qWpP4DAwtmCmiDwE5AJTSlpZREaLyBIRWZKamurD7oypeF+s2MXURdu46aw4BnUo84uuMdWCL4n/PiADWItTrO074CEf1ksBUlR1kft+Bu4A7iIyCjgfuFpVSxzbV1UnqmofVe0THR3tw+6MqVjb0o8w9qMV9GzVgHuHdvA6HGMqTKm3c7rdOm+o6ijgXyeyYVXdLSLbRaSDqq4DhgCrRWQYcD9wlvt0sDFVTnZuPrdNW4oIvDyyJ+Gh9siLCRylJn5VzROR5iISrqo55dj+GGCKiNQANgHXAf/FqffzjVvUaqGq3lyObRvjN0/PXsuKlAO8dk0vWjaq7XU4xlQoXx7g2gTMFZFPcPrpAVDVl8taUVWTgD7FJtu9cKZK+37tHibN28y1p7dmWJfmXodjTIXzJfGn4tx7X9v9MSZg7TpwlHumL6dz83o8OLyT1+EY4xdlJn5VfaQyAjHGa7l5+dwxLU3cZrsAABdwSURBVIms3HwmXNWTiPBQr0Myxi/KTPwi8g3O8IvHKBiI3ZhA8fJ3G1i8ZR/PX9GduGgbZtoELl+6eh4u8joCuBTI8k84xnhjQXIar/yQzGW9Y7iklz1MbgKbL109i4pN+klEfvJTPMZUurRDWdzxQRJxUZGMu/BUr8Mxxu986eqpV+RtCNAbsFsdTEDIz1funr6cA0dzeOf6ftSuYZXKTeDz5VO+CqePX3BKLGwGbvRnUMZUlolzN/Hz+lT+cVEXOjWvV/YKxgQAXxJ/XPGHt0TETotMtZe49Vee+Xodw7s24+rTWnkdjjGVxpfn0Iv38QMsruhAjKlMB47kcPu0ZbRoEMGTl3TDfYrcmKBw3DN3EWmC05dfS0S68r9SzPWwB7lMNaaq3PfRcvYczGTGLQOoXyvc65CMqVSlddmch1MvPwb4Z5HpGYA91GWqrXcXbuXrVXt4aHgnerS0kT9N8Dlu4lfVN3EGUblCVadXYkzG+M2qnQf4x+drGNwhmhsS2ngdjjGe8OU+/ukici5wKs4DXAXTn/BnYMZUtENZudw2dRkNI8N57ooehNgQiiZI+XIf/z+BBjjj576J8+TuQj/HZUyFUlUembWSremHmXpjfxpF1vA6JGM848tdPQmqehWQ7hZsOw0bIN1UMzMSU5i5bAd3DGlP/7jGXodjjKd8Gmy94F8Raea+j/VbRMZUsOS9GTz6ySpOj2vMbWfbcBDG+PIg1pci0gB4FkgC8oC3/RqVMRUkMyeP26Yuo3aNUF4c2YNQ69c3pswxd0OAr1R1P/ChiHwO1FLVfZUSnTEnadznq1m7O4O3rutL03oRZa9gTBAotatHVfOBl4q8P2pJ31QXX6zYxdRF27jprDgGdWjidTjGVBm+9PF/IyIX+j0SYyrQtvQjjP1oBT1bNeDeoR28DseYKsWXPv7bgPoikgUcxSndoKrayK+RGVNO2bn5jJm2FBF4eWRPwkN9Ob8xJnj4kvij/B6FMRXoma/XsjzlAK9d04uWjayslDHFlXkqpKp5wOXA/e7r5kAPfwdmTHl8v3YPr8/dzB/6t2ZYFxsvyJiSlJn4RWQCMBj4gzvpCPCaP4Mypjx2HTjKPdOX06l5PR46r5PX4RhTZfnS+TlAVW/CfZDLvavHp+fdRaSBiMwQkbUiskZETheRRiLyjYhscP9teBLxGwNAbl4+d7yfRFZuPq9e1ZOI8FCvQzKmyvIl8ee49/MrgIg0BvJ93P5LwGxV7Qh0B9YAY4HvVDUe+M59b8xJefn7ZBZv3sc/LupCXHQdr8MxpkrzJfG/CnwERIvIY8A84KmyVnIHaT8TmAygqtnug2AX8r8nf98GLipH3MYUWrAxjVe+38ClvWK4pJeVkTKmLL6UZX5HRBKBc9xJl6vqSh+2HQek4tT07w4kAncATVV1l7vtXe5IX8acsANHcvjnj8m8uWALbaIiGXfhqV6HZEy14Oug6aFADk53j683RYcBvYAxqrpIRF7iBLp1RGQ0MBqgVSsbCNv8T2ZOHm8v2MKrPySTkZXLJT1j+Mu5HYis6evH2Zjg5ks9/oeAq4CZOA9vTRWRKar6ZBmrpgApqlowWPsMnMS/R0Sau2f7zYG9Ja2sqhOBiQB9+vRRn1pjAlpevjJz2Q6en7OOnQcyGdwhmvuGdaRT83peh2ZMteLLKdI1QG9VPQIgIo/jdNuUmvhVdbeIbBeRDqq6DhgCrHZ/RgHj3X8/OYn4TRBQVX5cl8pTs9eydncG3WPq89wVPTi9rdXVN6Y8fEn8W4stFwZs8nH7Y4ApIlLDXec6nK6i6SJyA7AN5+EwY0qUtH0/T365hkWb9xHbuDavXtWL4V2bIWLllY0pL18S/xFglYh8jdPHPxSYJyLPA6jq3cdbUVWTgD4lzBpSjlhNENmcdphnvl7Ll7/sJqpODf5+4amM7NfK6u4YUwF8SfxfuD8FbLxd4zepGVm89N163l+8nRphIdx5Tjx/GhhHHbtwa0yF8eV2zsmVEYgJboeycpn48yYmzd1Edm4+V/Zrxe1D4omuW9Pr0IwJOL7c1TMM+DvQ2l3eyjKbCpOdm8/7/93Gy99tIO1QNud1bc6953agTVSk16EZE7B8+f48AbgC+AXfSzUYUypV5YtfdvHM1+vYmn6E/nGNmDSqEz1aNvA6NGMCni+JPwVIcodhNOakLdiYxviv1rIi5QAdm9Xlzev6Mqh9tN2pY0wl8SXx3wd8JiI/AlkFE1X1ZX8FZQLTml0HGf/VWn5an0qL+hE8d3l3Lup5CqEhlvCNqUy+JP7HcMo1NMC6ekw5pPx6hOfnrGdm0g7qRYTz4PCOXHt6rJVONsYjviT+Jqra2++RmIDz6+Fs/vljMm8v2AoCo8+M489ntaN+7XCvQzMmqPmS+L8TkbNV9Xu/R2MCQmZOHm/O38I/f0zmcFYul/aK4a7ftadFg1peh2aMwbfEfyNwr4gcAbKx2znNceTlKx8lpvD8N+vZfTCTIR2bcN+wjnRoVtfr0IwxRfiS+KP8HoWp1lSV79bs5anZa9mw9xA9WjbgpZE9OC3OiqgZUxX58uRunoiMBOJU9QkRiQGa4lToNEEuceuvPPXVWhZv2UdcVCT/uroXw7pYETVjqjJfntydAITjDKP4BE7RtteAvv4NzVRlG1MP8czsdcxetZuoOjX5x0VdGNG3pRVRM6Ya8KWrZ4Cq9hKRZQCqus8ts2yC0N6Dmbz43QY++O92IsJCuPt37bkhoY2NfmVMNeLLb2uOiITglGRGRBpj9/MHnYzMHLeI2mZy8/P5Q//W3HZ2O6LqWBE1Y6qb4yZ+EQlT1VzgVeAjIFpEHsOp2/NYJcVnPJadm8+URVt55ftk9h3O5v+6t+Deoe1p3diKqBlTXZV2xr8Y6KWq74hIInAOzq2cl6vqykqJzngmP1/5bMVOnp2zju37jjKgbWPG/r4j3WKsiJox1V1pib/wtgxVXQWs8n84piqYtyGN8bPXsHLHQTo1r8fb13flzPgou1PHmABRWuKPFpHShlV83g/xGA+t3HGAp2avZe6GNE5pUIsXRnTnwu6nEGJF1IwJKKUl/lCgDkXO/E1g2r7vCM/NWcespJ00qB3Ow+d14g+nt6ZmmBVRMyYQlZb4d6nquEqLxFS6fYezmfB9Mu8t3EpICPx5UFtuHtSWehFWRM2YQOZTH78JLEez83hj/mZe+3Ejh7NzuaJPS+48pz3N6kd4HZoxphKUlviHVFoUplLk5uXzYWIKL3yznr0ZWfyuc1PuO7cD8U2tiJoxweS4iV9V91VmIMZ/VJU5q/fw9Oy1bEw9TO/WDXn16l70jbUCq8YEI78+Zy8iW4AMIA/IVdU+ItIDp9ZPBJAL/FlVF/szjmC2ZMs+nvxqLYlbf6VtdCT//kNvhnZuardmGhPEKqPAymBVTSvy/mngMVX9SkSGu+8HVUIcQSV5bwZPzV7HN6v30KRuTZ68pCuX944hzIqoGRP0vKispUA993V9YKcHMQSsvQczef6b9Uxfsp3aNcK4d2h7rk9oQ+0aVkTNGOPwdzZQYI6IKPBvVZ0I3Al8LSLPAiHAgJJWFJHRwGiAVq1a+TnMwLBwUzp/nrKUjMwcRg2IZczZ8TSKtEKqxphj+Tvxn6GqO0WkCfCNiKwFLgPuUtWPROQKYDJOHaBjuH8kJgL06dNH/RxntaaqvLdoG499uopWjWoz/ab+tGtid+oYY0rm18Svqjvdf/eKyEygHzAKuMNd5ENgkj9jCHTZufn89dOVTFu8ncEdonlxZE/q17IHsIwxx+e3K30iEikidQteA0OBlTh9+me5i50NbPBXDIEuNSOLq15fyLTF2/nzoLZMGtXXkr4xpkz+PONvCsx0bxsMA6aq6mwROQS8JCJhQCZuP745MStS9nPTu4n8eiSbl6/syQXdW3gdkjGmmvBb4lfVTUD3EqbPA3r7a7/BYNayHdz/0Qqi6tRkxs0D6HJKfa9DMsZUI3aPXzWSl688NXstE3/eRL/YRvzzml429KEx5oRZ4q8mDhzJYcz7y/h5fSrX9G/Fo+efSo0wexjLGHPiLPFXA8l7M/jT20vYsf8oT1zclatOs+cajDHlZ4m/ivt29R7u/CCJiPAQpt7Y3wqrGWNOmiX+KkpVefWHZJ77Zj2ntqjHxD/0oUWDWl6HZYwJAJb4q6Aj2bn85cMVfPHLLi7s0YLxl3SjVg0bBtEYUzEs8Vcx2/cd4cZ3lrBuTwYP/L4jo8+MsxLKxpgKZYm/ClmwMY1bpywlN1958499GdShidchGWMCkCX+KkBVeec/Wxn3+WpiG9fm9Wv7EBddx+uwjDEByhK/x7Jy83h01io+WLKdIR2b8OLIHtSNsHo7xhj/scTvob0HM7n5vUSWbtvPbYPbcffv2hMSYv35xhj/ssTvkeXbnSJrB47m8OpVvTivW3OvQzLGBAlL/B74eGkKYz/+heg6NfnolgF0blGv7JWMMaaCWOKvRLl5+Yz/ai2T5m2mf1wjXr2qF42tyJoxppJZ4q8k+49kM2baMuZuSGPU6a15+PzOhIdakTVjTOWzxF8J1u/J4MZ3lrBz/1HGX9KVkf2syJoxxjuW+P1szqrd3PVBErVqhPH+6P70bm1F1owx3rLE7yf5+cor3yfzwrfr6RZTn3//oTfN61uRNWOM9yzx+8HhrFzumb6c2at2c3HPU3jykq5EhFuRNWNM1WCJv4JtSz/C6HeXsH5PBg+f14kbEtpYkTVjTJViib8CzU9O49apS1GFt6/vx8D4aK9DMsaY37DEXwFUlbcWbOEfX6whLiqS16/tQ2xUpNdhGWNMiSzxn6Ss3DwenrmSDxNT+F3nprwwogd1atphNcZUXZahTsKeg5nc9G4iSdv3c/uQeO4cEm9F1owxVZ5fE7+IbAEygDwgV1X7uNPHALcBucAXqnqfP+Pwh2XbfuWmdxM5lJXLv67uxe+7WpE1Y0z1UBln/INVNa3gjYgMBi4EuqlqlohUu2GmZiSm8ODHv9C0fk3euWEAHZtZkTVjTPXhRVfPLcB4Vc0CUNW9HsRQLrl5+Tzx5VremL+ZAW0b8+pVvWgYWcPrsIwx5oT4u0qYAnNEJFFERrvT2gMDRWSRiPwkIn1LWlFERovIEhFZkpqa6ucwy/br4WxGvbmYN+Zv5rozYnnn+n6W9I0x1ZK/z/jPUNWdbnfONyKy1t1nQ6A/0BeYLiJxqqpFV1TVicBEgD59+igeWrfbKbK2+0AmT1/WjSv6tPQyHGOMOSl+TfyqutP9d6+IzAT6ASnAx26iXywi+UAU4P1pfQlmr9zF3dOXU6dmGO/f1J9erRp6HZIxxpwUv3X1iEikiNQteA0MBVYCs4Cz3entgRpA2vG245X8fOWFb9Zz83tLiW9al8/GJFjSN8YEBH+e8TcFZrp1asKAqao6W0RqAG+IyEogGxhVvJvHa4eycrn7gyTmrN7Dpb1iePziLlZkzRgTMPyW+FV1E9C9hOnZwDX+2u/J2pp+mBvfWcLG1MM8en5nrjsj1oqsGWMCij25W8S8DU6RNYC3r+tHQnyUxxEZY0zFs8SPU2Rt8rzNPPHlGuKb1GXitb1p3diKrBljAlPQJ/7MnDwenPkLHy/dwbmnNuX5K3oQaUXWjDEBLKgz3O4Dmdz0XiLLt+/nrnPaM+bsdlZkzRgT8II28Sdu/ZWb30vkSFYu//5Db849tZnXIRljTKUIysQ/fcl2Hp65kmb1I3jvhtPo0Kyu1yEZY0ylCarEn5OXz+NfrOGtBVtIaBfFhKt60qC21dsxxgSXoEn8+w5nc+uUpfxnUzo3JLThgd93JCzU3zXqjDGm6gmKxL9m10FufGcJezOyeO7y7lzaO8brkIwxxjMBn/i//GUX90xfTr1aYUy/6XR6tGzgdUjGGOOpgE78r/6QzDNfr6Nnqwb8+5reNKkX4XVIxhjjuYBO/LGNIxnRpyXjLjqVmmFWZM0YYyDAE/953ZpzXjcbBN0YY4qy21qMMSbIWOI3xpggY4nfGGOCjCV+Y4wJMpb4jTEmyFjiN8aYIGOJ3xhjgowlfmOMCTKiql7HUCYRSQW2lnP1KCCtAsOpDqzNwcHaHBxOps2tVTW6+MRqkfhPhogsUdU+XsdRmazNwcHaHBz80Wbr6jHGmCBjid8YY4JMMCT+iV4H4AFrc3CwNgeHCm9zwPfxG2OMOVYwnPEbY4wpwhK/McYEmaBI/CLydxFZISJJIjJHRFp4HZO/icgzIrLWbfdMEQn4wYZF5HIRWSUi+SISsLf8icgwEVknIskiMtbreCqDiLwhIntFZKXXsVQGEWkpIj+IyBr3M31HRW4/KBI/8IyqdlPVHsDnwKNeB1QJvgG6qGo3YD3wgMfxVIaVwCXAz14H4i8iEgq8Cvwe6AxcKSKdvY2qUrwFDPM6iEqUC9yjqp2A/sCtFfn/HBSJX1UPFnkbCQT8FW1VnaOque7bhUCMl/FUBlVdo6rrvI7Dz/oByaq6SVWzgfeBCz2Oye9U9Wdgn9dxVBZV3aWqS93XGcAa4JSK2n5Aj7lblIg8DlwLHAAGexxOZbse+MDrIEyFOAXYXuR9CnCaR7GYSiAisUBPYFFFbTNgEr+IfAs0K2HWQ6r6iao+BDwkIg8AtwF/rdQA/aCsNrvLPITztXFKZcbmL760OcBJCdMC/htssBKROsBHwJ3Fei5OSsAkflU9x8dFpwJfEACJv6w2i8go4HxgiAbIAxsn8P8cqFKAlkXexwA7PYrF+JGIhOMk/Smq+nFFbjso+vhFJL7I2wuAtV7FUllEZBhwP3CBqh7xOh5TYf4LxItIGxGpAYwEPvU4JlPBRESAycAaVX2+wrcfICeCpRKRj4AOQD5OeeebVXWHt1H5l4gkAzWBdHfSQlW92cOQ/E5ELgZeAaKB/UCSqp7rbVQVT0SGAy8CocAbqvq4xyH5nYhMAwbhlCjeA/xVVSd7GpQfiUgCMBf4BSdvATyoql9WyPaDIfEbY4z5n6Do6jHGGPM/lviNMSbIWOI3xpggY4nfGGOCjCV+Y4wJMpb4TUAQkQYi8uci71uIyAw/7esiEXnUff2WiFzmj/0cZ9/ni8hjlbU/E5gs8ZtA0QAoTPyqulNV/ZWQ7wP+6adtA4VVOEvyBXCBiNT25/5NYLPEbwLFeKCtO+bCMyISW1C7XUT+KCKzROQzEdksIreJyN0iskxEFopII3e5tiIyW0QSRWSuiHQsvhMRaQ9kqWpakclnisgCEdlUcPYvjmdEZKWI/CIiI9zpg0Tk8yLbmyAif3RfbxGRR0VkHnC5iNwuIqvdMRXeB3BLb/yIU4rDmHIJmFo9JuiNxRl/oAcUVjQsqgtOhcMIIBm4X1V7isgLOFVbX8QZ1PpmVd0gIqfhnNWfXWw7ZwBLi01rDiQAHXHKJ8zAGRegB9Ad52nT/4qIL+MEZKpqgtuGnUAbVc0qNpDOEmAgMN2H7RnzG5b4TbD4wa1rniEiB4DP3Om/AN3cKogDgA+dMimAU/KiuOZAarFps1Q1H1gtIk3daQnANFXNA/aIyE9AX6CsCotFy2evAKaIyCxgVpHpe4GAH0XO+I8lfhMssoq8zi/yPh/n9yAE2F/wjaEUR4H6pWxbiv1bXC7HdrFGFJt/uMjr84AzcQoLPiIip7qD60S4cRhTLtbHbwJFBlC3vCu7tc43i8jlUNhH372ERdcA7XzY5M/ACBEJFZFonAS+GKdIYGcRqSki9YEhJa0sIiFAS1X9AedicgOgjju7Pc4wk8aUiyV+ExBUNR2Y715Mfaacm7kauEFElgOrKHlIw5+BnlKkP+g4ZuJ01SwHvgfuU9Xdqrodp29+Bc7gOMuOs34o8J6I/OIu84Kq7nfnDca5u8eYcrHqnMacIBF5CfhMVb/1YN9NgamqWuI3BWN8YYnfmBPkJt/TVLXSB0ARkb5AjqomVfa+TeCwxG+MMUHG+viNMSbIWOI3xpggY4nfGGOCjCV+Y4wJMpb4jTEmyPw/GWdAVOOzC5wAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"import math\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"\n",
"Temp_t1 = 85\n",
"delta_t = 2\n",
"\n",
"time = [-3,-2,-1,0,1,2]\n",
"Temp_ambient = [55,58,60,65,66,67]\n",
"T_a = 65\n",
"\n",
"plt.plot(time,Temp_ambient,label='Ambient Temperature')\n",
"plt.title('Ambient Temperature throughout the day')\n",
"plt.legend();\n",
"plt.xlabel('time (hours)'); \n",
"plt.ylabel('Temperature (ºF)');\n"
]
},
{
"cell_type": "code",
"execution_count": 144,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"step is: 1\n",
"[169.5728738 126.11419753 100.27777778 85. 73.38888889\n",
" 69.4845679 ]\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEWCAYAAACJ0YulAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3hUZfbA8e9JgSQQAgmdBBK6lIQOUhTEgg1RpNgryNpXdxHs7i4r9p8NXVQWXRUFARGxIiJKFVR6lRo6CVKTkHJ+f9xJHMIkmZTJpJzP88yTmVvPTSAn7/ve+x5RVYwxxhiAAH8HYIwxpuywpGCMMSaHJQVjjDE5LCkYY4zJYUnBGGNMDksKxhhjclhSMKacEZG+IvK7iBwXkQH+jsdULJYUTIlw/YLKfmWJSIrb5+v8HV9xiMg+Eent7zjcjAOeVdXqqvqVv4MxFUuQvwMwFYOqVs9+LyLbgdtVda7/IvKOiASpakY5O0cTYG0JHu8MpfF9MWWTtRRMqRCRQBF5TES2isghEflARGq61rUWkQwRuU1EdotIkojcKiJni8gaEflDRF50O9YoEZknIv8RkaMisk5EznFbHyki77n+wt8lIk+ISECufV8XkcPAGNf554tIsogcFJF3RSTctf00oC7wjavVc6+IDBCRLbmuL6c1ISLjReRDEflYRI4Bw/O7/jy+X3e5uoiSRGSGiNRzLU8EGmbHk8e+sSIyy3WeQyLygtvP4CkR2Ski+0Vkktt1Zv8MRojILuALt2WjRGSviOwRkXvcztNLRH51/Qz2icjTbuv6iMhS18/uFxHp5cU/E1MWqKq97FWiL2A7cH6uZWOAH3F+oYUAk4H/uta1BhR4GagKDAROANOBKKAxcBjo7tp+FJAB3AkEAzcCyUAN1/ovgVeBMKAB8CtwU659RwCBQKjr/OcBVYD6wBJgvFvs+4Debp8HAFtyXV/ONsB4IA24BOcPr9D8rt/D9+8S1/HiXdtOBL7NK55c+wYD610xhLnO3dO17k7XuiZADeBz4K1cP4O33fbLXvau63NH1/c5+zp/BYa43oe7/XxigSTgfNf1XwIcBGr5+9+mvbz4/+vvAOxV8V55JIVtQC+3z3HASUDcfvlEua0/AVzh9nkOMMr1fhSwLdfxVwFDXL/wTgDBbutuAb5023dTAfEPBxa7fS5KUvjG2+v3cP4PgH+4fa4JZAH1PcWTa99+wG4gwMO6hcCtbp8TPPwMGrqtz14W67bsFeB11/tlwCPuPzfX8ieyk43bsh+AYf7+t2mvgl82pmB8TkQEiMHpknCfgTEApyUAkKmqSW7rUoD9uT5Xd/ucmOs0O3D+Cm+C89f1Qee0Oedx7+7ZlSu+hjitlJ44f/EGAHu9ubZ85JzDi+s/lGvfhsC87A+q+oeIHAUa4SSE/MTgJMwsD+sa4nyfsu3AaQFEuj5nqeqe/K7FtU/2oPtNwJPAJld32uOq+jXOz+AaERnitl+w6/ymjLOkYHxOVVVEdgNXqeqK3OtFpHYRDhud63NjYA/OL7DjOF0VeU0BnHv5cziti3aqelhEhgP/ymf7EzhdLACISDB//mI9Y5+Crt+DPTi/WLOPH4HT3bPbi313AbEiEuAhMZx2XJzvWQpOl1AdzrzObDE4rb/sffYAqOp6YJiIBOK0rmaISC1XDG+r6j0ejmXKOBtoNqXlTWC8iMQAiEhdEbm8GMeLcQ2ABonI9Ti/rL5R1W04YwLPiki4iASISAvJ/5bScJxEclREGgMP5Fq/H2jq9nk9ECki/V0J4SkK/r9UmOufAowQkXYiEgI8A8xT1YJaCQA/AceAf4pImIiEikhPt+P+TUQauwaY/wV8mE/yzPaE6zgJwA3Ax65ruFFEolQ1EziCk1SycMYghri+P4GuffuLSH0v4jd+ZknBlJZngbnAPNcdOYuATsU43gL+HPh8BLhSVY+41l2D0w+/wbX+Y6BePsd6HKdL5AgwE2eA2904YJzrTpq7VfUQcB9O338iTpdO7i6g3Ly+flX9HHga+Aznr/L6OL+MC6Sq6TgDuwmu2HYCV7lWvwHMcJ37d5zvTe4EmFsmsBRnTOQrnLGOBa51lwEbXdfzNDBUVTNUdSswGCdZHsLpcroP+31TLkjBfyQYU7aIyCjgalU939+xVGQi0hpYo6rWzVyJWOY2xhiTw5KCMcaYHNZ9ZIwxJoe1FIwxxuQo1wNItWvX1tjYWH+HYYwx5cqKFSsOqWodT+vKdVKIjY1l+fLl/g7DGGPKFRHZkdc66z4yxhiTw5KCMcaYHJYUjDHG5CjXYwrGmIojPT2dxMREUlNT/R1KhRESEkJ0dDTBwcFe72NJwRhTJiQmJhIeHk5sbCxu056bIlJVkpKSSExMJC4uzuv9rPvIGFMmpKamEhUVZQmhhIgIUVFRhW55WVIwxpQZlhBKVlG+n5UyKfxx8hRPzV7LkZR0f4dijDFlSqVMCjuTT/Luou08/cV6f4dijKmAJk+ezN13313gNnv2/Fn99Pbbb2fdunWFPtf8+fO57LLLCr1fXiplUoiPrsmIPk356OddLNxSUG0UY4wpebmTwttvv02bNm38GJGjUiYFgPvPb0lsVBhjZ6zm5KkMf4djjCkjBg0aROfOnWnbti0TJ04EoHr16jzyyCMkJCTQo0cP9u/fD8Ds2bPp3r07HTt25Pzzz89Znu3YsWPExcWRnu50VR89epTY2FimTZvG8uXLue666+jQoQMpKSn07ds3Z9qer776ik6dOpGQkED//v0BWLZsGT179qRjx4707NmTjRs3+uT6K+0tqaFVAhk/OJ7hE5fw4jebePQy/2doY4zjqdlrWbfnaIkes03DGjxxedsCt5s0aRKRkZGkpKTQtWtXBg8ezIkTJ+jRowfjxo1j9OjRvPXWWzz66KP07t2bJUuWICK8/fbbPPvss7zwwgs5xwoPD6dv377MmTOHQYMG8dFHHzF48GCGDBnC66+/zvPPP0+XLl1OO//BgwcZMWIECxYsIC4ujuTkZABat27NggULCAoKYu7cuTz88MNMn567cmzxVdqkANCjaRTXdW/MpIXbuDS+AR0b1/J3SMYYP3vllVeYOXMmALt27WLz5s1UqVIlp9++c+fOfPvtt4DzbMWwYcPYu3cvp06d8vg8wO23386zzz7LoEGD+O9//8tbb72V7/mXLFnCOeeck3OsyMhIAI4cOcJNN93E5s2bEZGc1kdJq9RJAWDMxa2Zt+EAD01fxex7elM1KNDfIRlT6XnzF70vzJ8/n7lz57J48WLCwsLo27cvqampBAcH59zeGRgYSEaG0+V8zz338MADDzBw4EDmz5/Pk08+ecYxe/Xqxfbt2/nhhx/IzMykXbt2+cagqh5vJX3sscfo168fM2fOZPv27fTt27fY1+tJpR1TyBYeEsy4K9uxaf9xJnz/u7/DMcb40ZEjR6hVqxZhYWFs2LCBJUuWFLh9o0aNAHj33Xfz3O7GG2/kmmuu4ZZbbslZFh4ezrFjx87Y9uyzz+aHH35g27ZtADndR+7nmjx5cqGuqzAqfVIAOK91Pa7o0JAJ87ewYV/J9mMaY8qPAQMGkJGRQXx8PI899hg9evTId/snn3ySIUOG0KdPH2rXrp3ndtdddx2HDx/mmmuuyVl28803M2rUqJyB5mx16tRh4sSJXHXVVSQkJDBs2DAARo8ezdixY+nVqxeZmZnFvNK8lesazV26dNGSKrKTdDyNC15aQEytUGbc2YvAAHuy0pjStH79es466yx/h+ETn3zyCbNmzeJ///tfqZ/b0/dVRFaoahdP21tLwSWqelWeHNiWlYlHmPTTNn+HY4ypIO655x7GjBnDY4895u9QvFLpB5rdXR7fgM9+280L327kgjb1iK1dzd8hGWPKuVdffdXfIRSKtRTciAj/HNSO4IAAxsxYRXnuWjPGmKKwpJBLg4hQxl5yFku2JvPRz7v8HY4xxpQqSwoeXNMthrObRvHvOevZd8SqQBljKg+fJQURmSQiB0RkjduyDiKyRER+E5HlItLNbd1YEdkiIhtF5CJfxeUNEeHpq9qTnpXFo5+utm4kY0yl4cuWwmRgQK5lzwJPqWoH4HHXZ0SkDTAcaOvaZ4KI+PXR4tja1XjwglbMXX+A2av2+jMUY0wl8dlnnzF+/Pgi7RsbG8uhQ8Wf9dlnSUFVFwDJuRcDNVzvI4DseWOvAD5S1TRV3QZsAbrhZ7f0iiUhOoInP1tL8olT/g7HGONu1VR4qR08WdP5umqqvyMqloyMDAYOHMiYMWP8GkdpjyncDzwnIruA54GxruWNAPdR3UTXsjOIyEhX19PygwcP+jTYoMAAnr06gWOp6fxj9lqfnssYUwirpsLse+HILkCdr7PvLXZi2L59O2eddRYjRoygbdu2XHjhhWdMa33o0CFiY2MBZ7qJQYMGcfnllxMXF8drr73Giy++SMeOHenRo0fOFBW///47AwYMoHPnzvTp04cNGzYAzlPNDzzwAP369eOhhx46rTjP/v37ufLKK0lISCAhIYFFixYBnqf2Lkml/ZzCX4C/qup0ERkKvAOcD3h6fNhjR76qTgQmgvNEs68Czdaqfjh39m3Oy99tZmCHhpzXup6vT2mM+XIM7Fud9/rEnyEz7fRl6Skw625YkcccRPXbw8UFd81s3ryZKVOm8NZbbzF06NACp6des2YNv/76K6mpqTRv3pxnnnmGX3/9lb/+9a+899573H///YwcOZI333yTFi1asHTpUu68807mzZsHwKZNm5g7dy6BgYGnzWl07733cu655zJz5kwyMzM5fvw44Hlq76ioqAKvy1ulnRRuAu5zvZ8GvO16nwjEuG0XzZ9dS353Z79mfLlmLw/PWMO3D0QSHhLs75CMqdxyJ4SClhdCXFwcHTp0AJxpsrdv357v9v369SM8PJzw8HAiIiK4/PLLAWjfvj2rVq3i+PHjLFq0iCFDhuTsk5b2Z5xDhgwhMPDMIdR58+bx3nvvAc7MrBEREYDnqb3Lc1LYA5wLzAfOAza7ln8GfCgiLwINgRbAslKOLU9VgwJ5ZnA8V72xiPFfbmDcle39HZIxFVtBf9G/1M7VdZRLRAzcMqdYp65atWrO+8DAQFJSUggKCiIrKwuA1NTUPLcPCAjI+RwQEEBGRgZZWVnUrFmT3377zeP5qlXzfuaEvKb2Lkm+vCV1CrAYaCUiiSJyGzACeEFEVgL/BkYCqOpaYCqwDvgKuEtVfTcNYBF0bFyLW3vF8cHSnSzZmuTvcIyp3Po/DsGhpy8LDnWW+0BsbCwrVqwAnMntCqNGjRrExcUxbdo0wKmXsHLlygL369+/P2+88QYAmZmZHD16tNBTexeFL+8+ukZVG6hqsKpGq+o7qvqTqnZW1QRV7a6qK9y2H6eqzVS1lap+6au4iuPBC1vSODKMMdNXkXKqTOUsYyqX+KFw+StOywBxvl7+irPcB/72t7/xxhtv0LNnzyLd9vnBBx/wzjvvkJCQQNu2bZk1a1aB+7z88st8//33tG/fns6dO7N27dpCT+1dFDZ1diEt2nKIa99eyh3nNGXsJRVzml9j/KEiT53tTzZ1to/1bF6b4V1jeOvHraxK/MPf4RhjTImypFAEYy85izrhVRn9ySpOZWT5OxxjjCkxlhSKICI0mH8Nas+Gfcf4zw9W19mYklKeu7PLoqJ8Py0pFNEFbepxWXwDXp23hc37zyy+bYwpnJCQEJKSkiwxlBBVJSkpiZCQkELtZ5XXiuHJgW35acshRk9fxSejelpdZ2OKITo6msTERHw9fU1lEhISQnR0dKH2saRQDLWrV+WJy9vw149X8u6i7dzaO87fIRlTbgUHBxMXZ/+H/M26j4ppUIdG9GtVh+e+3siu5JP+DscYY4qlciaFEpxyV0QYd2V7AgTGzrCCPMaY8q3yJQUfTLnbsGYoYy45i5+2HGLa8sSSi9UYY0pZ5UsK3/3DmWLXXXqKs7wYruvWmG6xkfxzzjr2H7W6zsaY8qnyJYUjefwln9dyLwUECOMHt+dURhaPfbrGupGMMeVS5UsKEXncnpXX8kJoWqc6f72gJd+s28+Xa/YV+3jGGFPaKl9S8DTlrgRAv4dL5PC3946jXaMaPD5rDYetrrMxppypfEkh95S7oZGgWbDtRyiBLp+gwACeHZzAHyfT+eecdcWP1xhjSlHlfHgtfujp867PHw/zn4aopnDO34t9+DYNazDq3Ga89v0WBiY0pG+rusU+pjHGlIbK11Lw5NyHIH4YzPsXrMm/SLe37unfnGZ1qvHIzDUcT8sokWMaY4yvWVIAEIGBr0Ljs2HmX2BX8ctDVw0K5Nmr49lzJIVnv9pQAkEaY4zvWVLIFlQVhn0ANRrClGsgeVuxD9m5SSQ3nR3Le4t3sGxbcgkEaYwxvmVJwV21KLjuE8jKgA+HQkrxK6v9/aJWRNcKZcz0VaSmW11nY0zZZkkht9rNYfgHTkth6o2QmV6sw1WrGsTTV7Vn66ETvPLd5hIK0hhjfMOSgiexvZ0xhm0/wJwHin2rap8WdRjSOZr/LNjKmt1HSihIY4wpeZYU8tLhGuf21F/eg4UvF/twj17ahshqVRj9ySrSM62uszGmbPJZUhCRSSJyQETW5Fp+j4hsFJG1IvKs2/KxIrLFte4iX8VVKH0fhrZXwdwnYN2sYh0qIiyYf17RlnV7jzJxwdYSCtAYY0qWL1sKk4EB7gtEpB9wBRCvqm2B513L2wDDgbaufSaISKAPY/NOQAAMegOiu8GMkZC4oliHG9CuARe3q8/L321my4HjJRSkMcaUHJ8lBVVdAOS+D/MvwHhVTXNtc8C1/ArgI1VNU9VtwBagm69iK5TgELhmClSvB1OGwx87i3W4p65oS2hwIGOmryIry2ZSNcaULaU9ptAS6CMiS0XkBxHp6lreCNjltl2ia1nZUK02XDcNMtLgw2GQWvTB4rrhITx2WRuW7zjM/5bsKMEgjTGm+Eo7KQQBtYAewN+BqSIigHjY1uOf0SIyUkSWi8jygwcP+i7S3Oq0gmHvwaFNMO0WyCz61BWDOzWiT4vaPPPVBhIPW11nY0zZUdpJIRGYoY5lQBZQ27U8xm27aGCPpwOo6kRV7aKqXerUqePzgE/TtC9c9hL8/h18+fci36oqIvz7yvYAPDzTCvIYY8qO0k4KnwLnAYhIS6AKcAj4DBguIlVFJA5oARR/AiJf6HQj9Loflk+Cxa8X+TAxkWE8NKA1CzYdZMYvu0swQGOMKTpf3pI6BVgMtBKRRBG5DZgENHXdpvoRcJOr1bAWmAqsA74C7lLVsjsnRP8noM0V8M2jsGFOkQ9zQ48mdGlSi398vo4Dx6yuszHG/6Q8d1106dJFly9f7p+Tp6fA5EvhwHq45Qto2LFIh9ly4DiXvPwj57epy4TrOpdwkMYYcyYRWaGqXTytsyeaiyo4FK75CMJqw4fD4UhikQ7TvG517ju/BV+s3sdXa/aWcJDGGFM4lhSKo3pduPZjSD/p3KqadqxIhxl5TlPaNKjBY7PWcuRk8SbgM8aY4rCkUFz12sCQyU430ie3FulW1eDAAJ69Op7kE6cY94XVdTbG+E+BSUFEokTkchG5Q0RuFJFOrmcLTLbm/eHS52HzN/D1w0U6RLtGEYw8pylTlyfy4+ZSfP7CGGPc5JkURKSPiHwBfAtcCcQBnYB/AWtE5DERqV46YZYDXW6Fs++GZf+Bpf8p0iHu69+CprWrMXbGak5YXWdjjB/k11K4CrhbVTup6q2qOkZV71fVS4COwHpyTXhX6V3wD2h9GXw1BjZ9XejdQ4IDGT84nsTDKTz/zUYfBGiMMfnLMymo6l9V1eMcz6p6SlU/UdVPfBdaORQQCFdNhPrxzlQYe1cV+hDd4iK58ewmTF60nRU7DvsgSGOMyVt+3UfvuL2/vnTCqQCqVHNuVQ2t6dyRdLTwt5mOHtCaBjVCeGj6KtIyyu4zfMaYiie/7qNObu8f8HUgFUqNBnDtVEg7ClOGwakThdq9etUgxl3Vni0HjvPavC0+CtIYY86UX1Iov486lwX128HV/4V9q2H67ZBVuL/4+7Wqy1UdG/HG/N9Zt+eoj4I0xpjT5ZcUokXkRRF5ye19zqu0AizXWl4IFz8LG7+Abx4r9O6PXdaGmmHBPDR9FRlW19kYUwrySwpjgbXAGrf37i/jjW4joPsoWPI6/Px2oXatVa0KTw1sx+rdR3jnp20+CtAYY/4UlNcKVX1HREKct075TFNEF/0bDm+HL0ZDzVhocb7Xu17Svj4XtqnHi99u4oI29Whaxx4NMcb4Tn53Hz2EU+dgtoiMdi2bUFqBVSgBgTD4HWdKjGk3w37vG1oiwj8HtaNKUABjZqy2us7GGJ/Kr/uonapeqKoX4jysBk4pTVMUVavDNR87Xz8cBsf2e71rvRohPHZpG5ZtS+bDZTt9GKQxprLLLylUd81z1BWwPouSENHImVX1ZBJMGQ6nvK/PPKRLNL2aRzH+yw3s+SPFh0EaYyqz/JLCg8D1wHDgPteyd30eUUXXIMHpStrzK8wcCVne3VUkIjx9ZTyZWcojM1dbXWdjjE/kN83FVlV9QFUfzJ7uQlW/Kr3QKrDWlziDz+tnw3dPer1b46gw/nZRK77feJBZv+3xXXzGmEorv4Hm70XkLyLSMNfyIBE5R0TeEZFbfB9iBdXjL9D1dlj4MqyY7PVuN/eMpWPjmjw1ey2HjttNYcaYkpVf99GlQDAwU0QSRWSViGwGtgK3AG+o6n9LI8gKSQQGPAPNL4DPH4Dfv/dqt8AA4dnB8ZxIy+Sp2VaQxxhTsvLrPjqpqq+oanegGU6SOFtVG6vqLaq6vNSirKgCg+DqSVCnNUy9CQ5s8Gq3FvXCufu85sxeuYdv13l/F5MxxhTEq3KcqpqmqrtU9ZCvA6p0Qmo4dyQFh8CHQ+C4d1XXRp3bjNb1w3n009UcSbG6zsaYkmE1msuCmjHOdNvHD8JH10B6wbecVgkK4JnB8Rw8lsb4L9eXQpDGmMrAkkJZ0agTDH4LEpfDp3/x6lbVhJiajOjTlCnLdrFoizXijDHF51VSEJFoEennel9VRKp5sc8kETkgIms8rPubiKiI1HZbNlZEtojIRhG5qDAXUWGcdblT0nPtTPj+X17tcv/5LYmNCmPMjNWknLKCPMaY4ikwKYjIrThzIGVP8dkEmOXFsSfjoYaziMQAFwA73Za1wXlIrq1rnwkiEujFOSqenvdAp5vgxxfg1/cL3Dy0SiBPXxXPzuSTvGB1nY0xxeRNS+FeoAdwFEBVNwF1C9pJVRcAyR5WvQSM5vQiPlcAH7kGtLcBW4BuXsRW8YjApS9A034w+z7YtqDAXc5uFsW13RszaeE2ft1pdZ2NMUXnTVJIVdVT2R9cf8FLUU4mIgOB3aq6MteqRsAut8+JrmWejjFSRJaLyPKDB727U6fcCQyGoe9CVHP4+AY4tLnAXcZe3Jp6rrrOpzKsII8xpmi8SQoLXVNnh7jGFT4GPi/siUQkDHgEeNzTag/LPE7uo6oTVbWLqnapU6dOYcMoP0IinDrPgcHwwRA4kZTv5uEhwYy7sh2b9h9nwnyr62yMKRpvksJo4BiwAWdivO9wfrkXVjMgDlgpItuBaOAXEamP0zKIcds2GrDJfWo1geFT4Nhe+OhaSE/Nd/PzWtfjig4Nef37LWzcd6yUgjTGVCT5JgVXV9EkVX1DVa9U1UGu94Xun1DV1apaV1VjVTUWJxF0UtV9OAPZw113NsUBLYBlhb+cCiimK1z5JuxaAp/dDQXMjvr4ZW0IDwlm9CcrybSCPMaYQso3KahqJtBARIILe2ARmQIsBlq55k66LZ/zrAWmAuuAr4C7XOc2AG2vhP6Pw+ppMH98vptGVa/KE5e3YWXiEf670Oo6G2MKJ88azW62Aj+KyCzgRPZCVX0lv51U9ZoC1sfm+jwOGOdFPJVT7wcgeSv8MB4im0LCsDw3HZjQkNkr9/D8Nxu5oE09mkQV+FiJMcYA3o0pHAS+BcKAOm4vU5pE4NKXILaP0420Y1E+mzp1nYMDAhgz3QryGGO8J+X5F0aXLl10+fJKNllrymF4+wI4eQhu/w6imuW56YdLd/LwzNU8fVV7runWuBSDNMaUZSKyQlW7eFrnzRPN34rIN7lfJR+m8UpoLbhuKkiAc6vqSU/PBzqGd42hR9NI/j1nPfuO5H/nkjHGgHfdR48Cj7le43BuTc398JkpTZFNYfiHcGSX83BbximPmwUECOOviic9K4tHP7VuJGNMwQpMCqq61O31g6reS2WdgqIsadwDBr0BO36C2ffmeatqbO1qPHhBK+auP8Dnq/aWcpDGmPLGm+6jGm6vmiLSH2hQCrGZgrS/Gvo9AiunwILn89zsll6xJERH8MRna9l68HgpBmiMKW+86T5aC6xxff0V52nmEb4MyhTCOX+H+OHOVNurP/G4SVBgAM8PSQBg0OsLWWi1F4wxefAmKTR11WWOUdU4VT0PWOjrwIyXRGDgK9CkF3x6J+xc6nGzFvXCmXVXLxpEhHLjpGX8b8mOUg7UGFMeeJMUPP2WsSkoypKgqjDsfYiIdsp5Jnt+kjkmMoxP/nI257asw2OfruGJWWvIyLQZVY0xf8ozKYhIXRFJAEJFpL2IxLtevXEeZDNlSVgkXDcNNAs+HOo8z+BBeEgwb93YhRF94nh38Q5umfwzR1LSSzlYY0xZlV9L4VLgNZwZSycAr7teD+PcnmrKmqhmMOwDp6Uw9cY8b1UNDBAeubQNzwxuz+Lfk7hywkK2HzrhcVtjTOVS4BPNIjJUVaeWUjyFUimfaPbGyo9g5h3Q8XoY+Joz7pCHJVuT+Mv7K8hSeOP6TvRsVjvPbY0xFUOxnmhW1akicpGIPCAiD2e/Sj5MU2IShsM5o50azwv/L99NezSN4tO7elEnvCo3vrOMKct25ru9MaZi8+Y5hQnATcADQChwPdDcx3GZ4ur3MLS7GuY+CWs/zXfTJlHVmHFnT3o2r83YGav5x+x1VovBmErKm+6jVaoaLyIrVTVBRMKB6ap6YemEmDfrPipAeiq8dwUkLodqkYlUBgMAAB2MSURBVHD8oHOHUv/HIX7oGZtnZGYx7ov1/Hfhdvq2qsMr13SkRkihS2kYY8q4YnUfAdkzqaW6SmemArElFJvxpeAQ55e/ZsLxA4A68yXNvhdWnTlMFBQYwBOXt2Xcle34afMhBk9YxM6kk6UftzHGb7xJCl+ISE3geeA3YDvg+dFZU/b89BKQqzWYngLf/SPPXa7r3oT3bu3GgWNpXPH6TyzdmuTbGI0xZUZBNZoDgC9V9Q9VnQbEAe1V1Qaay4sjiYVb7tKzeW0+vasXtapV4fp3ljL1510+CM4YU9YUVKM5C3jZ7XOKquY9gb8peyKiPS8PCsm3FgNAXO1qzLyzFz2aRjF6+ir+/cV6G4A2poLzpvvoWxG5wueRGN/o/zgEh56+LCAYMtJgwtmwZW6+u0eEBvPfm7ty49lNmLhgKyPfW87xtAwfBmyM8SdvksLdwEwRSRGRZBE5LCLWWigv4ofC5a9ARAwgztdBE+COH5wqbu8Phjl/g1N5DygHBQbwjyva8Y8r2jJ/00EGT1jErmQbgDamIvLmltRAT8tVNdMnERWC3ZJaTOmpMO+fsPg1iGoOV02ERp3z3eXHzQe564NfCA4M4D83dKZLbGQpBWuMKSnFfaI5ExgCPOR63wDoULIhGr8IDoGLxsGNnzl3JL19Acx/BjLz7h7q06IOM+/qRY3QYK59aynTV+Q/YG2MKV+8eaL5NaAfcINr0UngTV8GZUpZ03PhL4ug3WCY/2+YdBEk/Z7n5s3qVGfmnT3p3KQWD05byTNfbSDLBqCNqRC8GVPoqap34HqIzXX3UZWCdhKRSSJyQETWuC17TkQ2iMgqEZnpev4he91YEdkiIhtF5KIiXIspjtCaMPgtuHoSJG2BN3vD8kl51n6uGVaF927rxrXdG/PG/N8Z9f4KTtgAtDHlnjdJId31vIICiEgU4E1llsnAgFzLvgXaqWo8sAkY6zpmG2A40Na1z4S8xjKMj7UbDHcuhpju8PlfndoMx/Z73DQ4MIBxg9rx5OVtmLt+P1e/uZjdf6SUcsDGmJLkTVJ4HZgO1BGRp4CfgGcK2klVFwDJuZZ9o6rZf04uwanVAHAF8JGqpqnqNmAL0M27SzAlrkZDuH4GXPwcbFsAE3rAus88bioi3Nwrjkk3dyUx+SRXvLaQX3Z6LvBjjCn7vBlofg94FGeai2RgiKp+VALnvhX40vW+EeD+yGyia9kZRGSkiCwXkeUHDx4sgTCMRwEB0H0k3PEj1GwMU29wakCnHvW4ed9WdZl5V0/CqgQyfOISZv22u5QDNsaUBG9aCgCBQDpwqhD75ElEHgEygA+yF3nYzGNntqpOVNUuqtqlTp06xQ3FFKROS7h9rlOfYeUUeKMXbF/ocdPmdcOZdVcvOsbU5L6PfuP5rzfaALQx5Yw3dx89AkwBGuJ093woImOLekIRuQm4DLhO/3xIIhGIcdssGthT1HOYEhYYDOc9Ard+DQGBMPlS+PZx56noXGpVq8L/buvOsC4xvPb9Fu784BdOnrIBaGPKC2/+6r8e6Kqqj6rqIzh9/TcW5WQiMgB4CBioqu6PxH4GDBeRqiISB7QAlhXlHMaHYrrBqJ+g802w8GV46zzYv/aMzaoEBTB+cHsevfQsvl63j6H/WczeIzYAbUx54E1S2AEEuX0OArYWtJOITAEWA61EJFFEbgNeA8Jx5lP6TUTeBFDVtcBUYB3wFXBXWXhi2nhQtTpc/jJc8zEc3w8T+8KiVyHr9BvSRITb+zTlnZu6sP2QMwC9ctcf/onZGOM1b6a5mAF0Bb7G6ee/EOcOpP0AqvqAj2PMk01z4WcnDsHs+2DD5xDbBwa9ATVjzths475j3Pbuzxw8lsbzQxK4PKGhH4I1xmTLb5oLb5LCbfmtV9V3ihFbsVhSKANU4bcP4MuHQALgkuedSfjk9HsHko6nMer9Ffy8/TD39W/B/ee3QMTT/QXGGF8rVlIoyywplCGHt8PMUbBzMbQZBJe9BGGnT5aXlpHJIzPX8MmKRC6Nb8ALQxIICbZnFI0pbcWaEE9EBojIz64pK2zqbONZrVi4eQ6c/yRsmOOxVkPVoECeuzqesRe35ovVexn2n8XsP5rq6WjGGD/xZqD5NeAOnIfJ6gC1XV+NOV1AIPT+K4yYl2etBhHhjnObMfGGLmw+cJwrXlvImt1H/Bi0McadN0khEfhNVdNVNTP75evATDnWIB5Gzoez74af34L/9IHdK07b5II29Zj+l54EBghXv7mIL1bv9UuoxpjTeZMURgOzReTvInJv9svXgZlyzotaDWc1qMGnd/WiTYMa3PnBL7z63WbK8xiXMRWBN0nhKSATqInTbZT9MqZgBdRqqBNelQ9H9ODKjo144dtN3P/xb6SmW0PUGH8JKngT6qpq/jUajclPdq2GVgPg8wecWg0XjYPOt4AIIcGBvDg0geZ1q/Pc1xvZkXSSiTd2pm54iL8jN6bS8aal8J2InOfzSEzFl0+tBhHhrn7NefP6zmzcd4xBry1k7R4bgDamtHmTFEYAc0XkuN2SaoqtgFoNA9rVZ9qos1Hg6jcW8/Xaff6L1ZhKyJukUBsIBiKwW1JNSSigVkO7RhHMuqsXLeuHM+r9FUyYv8UGoI0pJd4U2ckEhgAPud43ADr4OjBTCeRTq6FujRA+HtmDy+Ib8uxXG3lw6krSMmwA2hhf8+aJ5teAfsANrkUngTd9GZSpRPKp1RASHMgrwzvw4AUtmfHrbq59aymHjp9Zw8EYU3K86T7qqap3AKkAqpoMVPFpVKbyyaNWg4hwT/8WTLiuE2v3HOGK1xayYZ/nkqDGmOLzJimki0gArvKYIhIFZOW/izFFkE+thkvaN2DqHWeTkZXF4AmLmLtuv7+jNaZCyjMpiEj2MwyvA9OBOiLyFE4thWdKITZTWbUaAHcugRYXwjePwnsD4Y9dxEfXZNZdvWlWtzoj/reciQt+twFoY0pYnlNni8gvqtrJ9b4tcD4gwFxVXVN6IebNps6u4PKo1ZCSnsXfpq1kzuq9DOkczbgr21MlyJtGrzEG8p86O78nmnMqoLjKZZ5ZjNcYXxKBjtdDbG+nVsPMkbDxC0Ive4lXr+lI87rVefm7zexIOskb13ciqnpVf0dsTLmXX0shEXgxrx1VNc91pcVaCpVIViYsegXmjYOwKBj0OjQ/n89W7uHv01ZSt0ZV3rmpKy3rhfs7UmPKvKIW2QkEqgPhebyMKT151GoYeFZNPr7jbFLTs7hqwiK+33jA35EaU655NaZQVllLoZJKT4V5/4TFr0FUc7hqInurt+H2d5ezfu9RHrm0Dbf2irUa0MbkoagtBfsfZcomD7UaGvz6CtNGduXCNvX55+freHjmak5l2J3TxhRWfkmhf6lFYUxR5KrVEPa/S5kwoAZ392vOlGW7uHHSUg6fOOXvKI0pV/LsPioPrPvI5Fgz3anVkHkKLvwXP+9Np+GK52ggSZwMqU+Vi56kSsfh/o7SmDKhqN1HxT3pJBE5ICJr3JZFisi3IrLZ9bWW27qxIrJFRDaKyEW+istUUO61GuY8QNdfx9JIDhGAUj11L5mz7uHLD18myeZOMiZfvnziZzIwINeyMcB3qtoC+M71GRFpAwwH2rr2mSAigT6MzVRE2bUaQmqCnj6eEMop2m94hZ7j5/HIzNVsO3TCT0EaU7b5LCmo6gIgdzGeK4B3Xe/fBQa5Lf9IVdNUdRuwBejmq9hMBRYQAKmeK7Y1Ckjiyo6NmLY8kfNemM8d/1vOih1WL8oYd6U9N0A9Vd0L4Ppa17W8EbDLbbtE17IziMhIEVkuIssPHjzo02BNORUR7XGxoIyX11lyY3Xu7tuMpduSGfzGYq6asJCv1uwlM6v8jq8ZU1LKyoQxnm5/9fg/VFUnqmoXVe1Sp44VgDMe9H8cgkNPXxYUAnH9YMMcIqdcxoPbbmfpRbv496VxHDp+ilHv/8J5L8znf4u3k3LKivmYyqu0k8J+EWkA4Pqa/fhpIhDjtl00sKeUYzMVRfxQuPwViIgBxPk68FW46VN4cANc9hIoVP3yAa796SLmt53De5eHUzOsCo/NWkvP8d/x4jcbraCPqZR8ekuqiMQCn6tqO9fn54AkVR0vImOASFUd7ZqF9UOccYSGOIPQLVzlP/Nkt6SaIlOFxJ/h57dh7UzIPIU26cnWJsN5bldLvt6QTHBgAIM7RXN7nzia1anu74iNKTH53ZLqs6QgIlOAvkBtYD/wBPApMBVoDOwEhrgquSEijwC3AhnA/ar6ZUHnsKRgSsSJQ/Dr+7B8EvyxA6rV4XDr4bx18lzeXp3BqYwszj+rHiPPaUrX2Fo2fYYp9/ySFEqDJQVTorKy4Pd5sPwd2PQVAGlNL2BO1Yv51/oGJKdkkhBTk5F9mjKgXX0CAyw5mPLJkoIxhfXHLlgxGX55D04cIKtmE36reyVPJXZkZXIwMZGh3N67KUO6RBNWJb+yJMaUPZYUjCmqjFOwYTb8PAl2/IQGVmFvo4t47eg5fLivIRGhVbihRxNu7NmEuuEh/o7WGK9YUjCmJBxY74w7rPwI0o5yslZrZgYO4Ond7TkVUI2rOjXi9j5xNK9r5UZM2WZJwZiSlHYc1nzi3Lm0bzVZwdVZFnEB/9rXkzUZjejfui4jzmlK97hIG5Q2ZZIlBWN8QRUSlzsD02tmQGYau2t0YMKxc5mW0onW0bUZ0acpF7erT1BgWXlO1BhLCsb43snkP29rPbyN1Cq1mK7n8cbxc6BmY27rHcfQLjFUq2qD0sb/LCkYU1qysmDrPPh5ErrpS1Dll6pdee3YufxapTPX9ojj5p6x1K1hg9LGfywpGOMPRxJhxbvwy7twfD+HguozKbUv07Uv53Row4hzmtKyng1Km9JnScEYf8pMhw2fw8/vwPYfyZAgvsrsxuT086neojcjz2nG2c2ibFDalBpLCsaUFQc3wvJJZP32IQFpR9lCYyan92dj3Yu5/tx2XNK+AcE2KG18zJKCMWXNqROw+hOyfn6HgH0rOUkIMzJ68U3YZZzTpy/DuzWmug1KGx+xpGBMWaUKu39Bf36LrDUzCMxMY3lWSz4JuIiorkO5oXdL6kfYoLQpWZYUjCkPTibDbx+StuQtqh7dTpKGMz2rLwdbXcfg83vRun4Nf0doKghLCsaUJ1lZsG0+JxdNJOT3b0Cz+CErnl/rDabbBcPp1bKuDUqbYrGkYEx5dWQ3qUsnkbl8MtVOHSJRa/Nd6MXUPXcE53drb4PSpkgsKRhT3mWmk77uc5Lnv0G9pKWc0kB+COxBasLN9L1wEOGhVWDVVPjuH87zERHRTq3q+KH+jtyUQZYUjKlAsg5sYvfc14ncPI1qeoItGs2J2vG0PzKPgIzUPzcMDnVqVVtiMLnklxSs7WlMORNQtyUx175MtbFbSOzzHMGh1UhI+uL0hACQnuK0HIwpBEsKxpRXVcKI7j+SJmOWoXgeeNYju/hl+wHSM7NKOThTXtnTMcZUABIRDUd2nbkcaPrfDnxBJ7bX7kfV1hfQuUU08dERVA0KLP1ATZlnYwrGVASrpsLse50uo2xBoZyMv4HkQ/uJ3D2PsMxjpGowP2a153u6cqjhebRp0ZTucVF0bFyTkGBLEpWFDTQbUxnkd/dRZjrsXEzq6s/QDXMIPbmHTAJYkdWSrzM7M59uRMW0onvTSHo0jaJT41qEVrEkUVFZUjDG/EkV9q6EDXPIXP85gQfXAbAtMI7ZpzrydUYXNgXEEh9di+5xkXRvGkXnJrVsLqYKxJKCMSZvyVthwxewYQ66awmiWRypUp8fA7vx0dF4Fme2goAg2jWKoEdcJN2bRtIlNpIaIcH+jtwUUZlLCiLyV+B2QIHVwC1AGPAxEAtsB4aq6uH8jmNJwZgSdvwgbPoKNsyB3+dBZhrpVWqyKaInX6R35v2DzTiSWYUAgTYNa9A9LorucZF0i4ukZlgVf0dvvFSmkoKINAJ+AtqoaoqITAW+ANoAyao6XkTGALVU9aH8jmVJwRgfSjvuJIYNc5xEkfoHGhTC4fq9+TnkbD451pYFuyEtIwsRaFUvnB5N/0wSUdWr+vsKTB7KYlJYAiQAR4FPgVeAV4G+qrpXRBoA81W1VX7HsqRgTCnJTIcdi5wEsWEOHE0ECSArpge76p7HfLrw7b4wVuw4TEp6JgAt6lane9NIpzXRNJK64TYFeFlRppICgIjcB4wDUoBvVPU6EflDVWu6bXNYVWt52HckMBKgcePGnXfs2FFaYRtj4LSBajbMgQNrneX12pHZ8hI21jqH+Ufqs3TbYZZvT+bEKSdJNK1d7bQk0SAi1I8XUbmVqaQgIrWA6cAw4A9gGvAJ8Jo3ScGdtRSMKQPcBqrZuRhQiIiB1peS2fIS1gS1ZemOIyzdmsyy7ckcS80AoHFkWM7dTd3jIomJDPPvdVQiZS0pDAEGqOptrs83Aj2A/lj3kTHlm4eBakJrQcsBTpKI68f6pEyWbktm6dYklm1P5o+T6QA0qhnqShJOa6JJVJjVjfCRspYUugOTgK443UeTgeVAYyDJbaA5UlVH53csSwrGlGGnDVR/CalHICgEmp0HrS+FlgPICo1i04FjLN2azNJtSSzdmkzSiVMA1K8RQje3JNGsTjVLEiWkTCUFABF5Cqf7KAP4Fef21OrAVJzksBMYoqrJ+R3HkoIx5UQeA9U0PttJEK0vhVqxqCq/HzzOkq3JOa2JA8fSAKhdveppLYkWdasTEGBJoijKXFIoKZYUjCmH8hmozkkQ9eNBBFVle9JJlm5NykkSe444U4RHVqtC19hadI+LIiGmJk2iwoiqVsVaE16wpGCMKbuSfoeN2QPVS3AfqKb1pdC4JwQ6U2yoKomHU1i8NSmnyynx8J+TAIZVCaRxZBgxkWE0iQyjcZTzvnFkGNG1Qm1mWBdLCsaY8uH4QWf8YcMc+P37MwaqaXYeVKl22i67/0hh476j7Ew6yc7kFHYmn2Bn8kl2Jp8kNf3POhIi0KBGSE6SaOxKGtnvIytRK8OSgjGm/Ek7Dr9/5/ZE9REICj1toJrfv8tzZlhV5eDxNFey+PO1K/kkO5JO5oxVZKtWJTAnYTRxJYvsz40qWCvDkoIxpnzLTIcdC90Gqnc7yyUA1K2qXCHqUqecyiTx8J/JYkeSkzCyP6dlnN7KaBgRSkxkaE7LIiYyjCZR1WgcGUatsOBy1cqwpGCMqThUYe9v8O5ASDt65noJhEadoUYDCG+Y62sDqNHQSR75yMpytTKST+a0NHYln2SHK2EczNXKqF41yNWqCKVJVLXTuqga1QylSlDZqnycX1KwCdKNMeWLCDTsCGnHPK/XTAgOgf3rYMt3cOr4mduE1HQlCA8JI7wBATUaUq96berVCKFrbOQZu588lUHi4RR2JjmJIruF8fvBE3y/8SCn3FoZAQINIkJPG8fIGQiPDKNmGWtlWFIwxpRPedSlJiIGbpr95+fUo3BsLxzdk+vrXji2Bw6sh+P7T++GAggIhvD6HpNHWHgDWtZoSMvmDaBKvdN2y8pSDhxLO20cY2eSM/j93YYDHDp+eisjPKeV4YxluLcyGvqhlWFJwRhTPvV//My61MGhznJ3ITWcV518Zs3JzIATB/5MFLm/5tvqiDitmyqgRgPqhzegfo2GdGvQAFo1hLDmEOD8cj95KoNdySnscCWK7FbG5gPHmLfxwBmtjIY1Tx/HyE4eTSKrERFW8oWOLCkYY8qn7MHkvOpSF0ZgkNN1VKMh0Dnv7bJbHe4tjaN7/2yB5NnqCILq9aGG08poVaMhrcIbQERDiM7utmpHVlAo+4+lehzHmLt+P4eOn2JgwE+MDppKjYCk4l1zHmyg2RhjSlK+rQ63hHLKw5jIaa2OBmd0XaVtX0bwvCcIyMjVOvLyjqtsNtBsjDGlxdtWR9qxvBNHHq0Oj7Xs0lOc1lIJtRYsKRhjjD9UDYc64VCnZd7bZGbAiYN/JoyPr/O83ZHEEgvLkoIxxpRVgUFO91GNBtAI584qj3dcRZfYKcvWExXGGGPy1v/xMx+883THVTFYUjDGmPIifqgzqBwRA4jztZCDzAWx7iNjjClP4oeWaBLIzVoKxhhjclhSMMYYk8OSgjHGmByWFIwxxuSwpGCMMSZHuZ77SEQOAjuKcYjawKESCqc8qGzXC3bNlYVdc+E0UdU6nlaU66RQXCKyPK9JoSqiyna9YNdcWdg1lxzrPjLGGJPDkoIxxpgclT0pTPR3AKWssl0v2DVXFnbNJaRSjykYY4w5XWVvKRhjjHFjScEYY0yOSp0UROSfIrJKRH4TkW9EpKG/Y/I1EXlORDa4rnumiNT0d0y+JiJDRGStiGSJSIW+bVFEBojIRhHZIiJj/B2Pr4nIJBE5ICJr/B1LaRCRGBH5XkTWu/5N31fS56jUSQF4TlXjVbUD8DlQcpUqyq5vgXaqGg9sAsb6OZ7SsAa4Cljg70B8SUQCgdeBi4E2wDUi0sa/UfncZGCAv4MoRRnAg6p6FtADuKukf8aVOimo6lG3j9WACj/qrqrfqGqG6+MSoOTq+JVRqrpeVTf6O45S0A3YoqpbVfUU8BFwhZ9j8ilVXQAk+zuO0qKqe1X1F9f7Y8B6nEKdJabSF9kRkXHAjcARoJ+fwylttwIf+zsIU2IaAe4FfBOB7n6KxfiYiMQCHYGlJXncCp8URGQuUN/DqkdUdZaqPgI8IiJjgbuBJ0o1QB8o6Jpd2zyC0xT9oDRj8xVvrrkSEA/LKnzrtzISkerAdOD+XD0exVbhk4Kqnu/lph8Cc6gASaGgaxaRm4DLgP5aQR5UKcTPuSJLBGLcPkcDe/wUi/EREQnGSQgfqOqMkj5+pR5TEJEWbh8HAhv8FUtpEZEBwEPAQFU96e94TIn6GWghInEiUgUYDnzm55hMCRIRAd4B1qvqiz45RwX5Q7FIRGQ60ArIwpmCe5Sq7vZvVL4lIluAqkCSa9ESVR3lx5B8TkSuBF4F6gB/AL+p6kX+jco3ROQS4P+AQGCSqo7zc0g+JSJTgL4400jvB55Q1Xf8GpQPiUhv4EdgNc7vLYCHVfWLEjtHZU4KxhhjTlepu4+MMcaczpKCMcaYHJYUjDHG5LCkYIwxJoclBWOMMTksKZgKTURqisidbp8bisgnPjrXIBF53PV+sohc7Yvz5HHuy0TkqdI6n6m4LCmYiq4mkJMUVHWPqvrql/VoYIKPjg3kzITqyRxgoIiE+fL8puKzpGAquvFAM1fNjOdEJDZ77n0RuVlEPhWR2SKyTUTuFpEHRORXEVkiIpGu7ZqJyFciskJEfhSR1rlPIiItgTRVPeS2+BwRWSQiW7NbDeJ4TkTWiMhqERnmWt5XRD53O95rInKz6/12EXlcRH4ChojIvSKyzlUT4yMA13Ql83GmLzGmyCr83Eem0huDUz+iA+TMLOmuHc5MkyHAFuAhVe0oIi/hzJ77fzgF0kep6mYR6Y7TGjgv13F6Ab/kWtYA6A20xplu4hOcug4dgAScp3B/FhFv6jykqmpv1zXsAeJUNS1XkaTlQB9gqhfHM8YjSwqmsvveNS/9MRE5Asx2LV8NxLtmo+wJTHOmnQGcaUJyawAczLXsU1XNAtaJSD3Xst7AFFXNBPaLyA9AV6CgmS7dpzhfBXwgIp8Cn7otPwBU+OqBxrcsKZjKLs3tfZbb5yyc/x8BwB/ZLY18pAAR+Rxbcn3NLYPTu3NDcq0/4fb+UuAcnEkcHxORtq7CSSGuOIwpMhtTMBXdMSC8qDu75qrfJiJDIGdMIMHDpuuB5l4ccgEwTEQCRaQOzi/3ZTgTMrYRkaoiEgH097SziAQAMar6Pc7Adk2gumt1S5zSo8YUmSUFU6GpahKw0DWw+1wRD3MdcJuIrATW4rnE5QKgo7j1MeVhJk73z0pgHjBaVfep6i6csYBVOIWPfs1j/0DgfRFZ7drmJVX9w7WuH85dSMYUmc2SakwJEZGXgdmqOtcP564HfKiqHlsYxnjLkoIxJcT1i7m7qpZ6YRsR6Qqkq+pvpX1uU7FYUjDGGJPDxhSMMcbksKRgjDEmhyUFY4wxOSwpGGOMyWFJwRhjTI7/B38uMCw8po//AAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"import math\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"\n",
"Temp_t1 = 85\n",
"delta_t = 2\n",
"#time = [-3,-2,-1,0,1,2]\n",
"Temp_ambient = [55,58,60,65,66,67]\n",
"T_a = 65\n",
"\n",
"n = 1\n",
"print('step is: ',n)\n",
"t = np.arange(-3,3,n)\n",
"\n",
"\n",
"def T_analytical(Temp_t1,T_a,K,t):\n",
" Temp_t2 = T_a + (Temp_t1-T_a)*(math.e)**(-K*t)\n",
" return Temp_t2\n",
"\n",
"\n",
"T_numerical = np.zeros(len(t))\n",
"T_numerical[3] = Temp_t1\n",
"for i in range(4,len(t)):\n",
" T_numerical[i] = T_numerical[i-1] + n*(-K*(T_numerical[i-1]-Temp_ambient[i]))\n",
"for i in range(2,-1,-1):\n",
" T_numerical[i] = T_numerical[i+1] - n*(-K*(T_numerical[i+1]-Temp_ambient[i]))\n",
"\n",
"print(T_numerical)\n",
"plt.plot(t,T_analytical(Temp_t1,T_a,K,t),'-',label='analytical');\n",
"plt.plot(t,T_numerical,'o-',label='numerical');\n",
"plt.title('Temperature of corpse');\n",
"plt.legend();\n",
"plt.xlabel('time (hours)');\n",
"plt.ylabel('Temperature (ºF)');\n"
]
},
{
"cell_type": "code",
"execution_count": 145,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"time of death was approximately 6.6 minutes after 10am, or at 10:06am\n"
]
}
],
"source": [
"#time of death is somewhere between hour -1 and hour 0 (based on values for T_numerical)\n",
"#do linear interpolation\n",
"Temp = 98.6\n",
"Temp0 = 85\n",
"Tempneg1 = 100.2777777\n",
"t0 = 0\n",
"tneg1 = -1\n",
"\n",
"#let times be y, temps are x\n",
"\n",
"tTemp = ((t0-tneg1)/(Temp0-Tempneg1))*(Temp-Tempneg1)+tneg1\n",
"minsbefore11 = -tTemp*60\n",
"minsafter10 = 60-minsbefore11\n",
"print('time of death was approximately',round(minsafter10,1), 'minutes after 10am, or at 10:06am')\n"
]
},
{
"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
}