06_01_heat_flow.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Heat Flow\n",
    "=========\n",
    "\n",
    "The problems consiste in determining the temperature variation along a\n",
    "bar of length $l$ at any instant of time, given the initial gradient of\n",
    "temperature.\n",
    "\n",
    "![bar](figures/bar.png)\n",
    "#### A bar with the ends in contact with heat reservoirs.\n",
    "\n",
    "The problem is described by the “Heat Equation” that can be derived as\n",
    "follows:\n",
    "\n",
    "At any instant of time, the heat flow through the bar equals the\n",
    "variatioon of energy inside the bar:\n",
    "$$\\mathrm{Heat\\,\\,Flow=Variation\\,\\,of\\,\\,internal\\,\\,energy}$$ or\n",
    "$$\\nabla \\cdot Q = -\\frac{dE}{dt}.$$ The variation of the internal\n",
    "energy is given by the body’s ability to store heat by raising its\n",
    "temperature: $$\\frac{dE}{dt} = \\rho c \\frac{dT}{dt},$$ where $\\rho$ is\n",
    "the density, and $c$ is the specific heat of the material. Fourier’s Law\n",
    "of heat conduction states that $$Q=-K\\nabla T,$$ where $K$ is the\n",
    "thermal conductivity. Hence, we obtain\n",
    "$$\\nabla \\cdot \\left( K \\nabla T\\right) = \\rho c \\frac{dT}{dt},$$ or\n",
    "$$\\frac{dT}{dt}=\\alpha \\nabla^2T,$$ with $\\alpha = K/(c\\rho)$. In 1d\n",
    "this equation is written\n",
    "$$\\frac{\\partial T(x,t)}{\\partial t} = \\alpha \\frac{\\partial ^2T}{\\partial x^2}.\n",
    "$$ We must solve this equation given the initial\n",
    "condition $$T(x,t=0) = f(x)\\,\\,\\,\\mathrm{(initial\\,\\,condition)}\n",
    "$$ and the boundary condition $$\\begin{eqnarray}\n",
    "T(x=0,t) &=& T_1, \\\\\n",
    "T(x=l,t) &=& T_2.\n",
    "\\end{eqnarray}$$\n",
    "\n",
    "Finite differences solution\n",
    "---------------------------\n",
    "\n",
    "![heat](figures/heat.png)\n",
    "#### The finite differences algorithm for the heat equation.\n",
    "\n",
    "\n",
    "The numerical solution is based on converting the differential equation\n",
    "into an approximate finite-diffence one. Following a derivation similar\n",
    "to the one we used for the wave equation we approximate the derivatives\n",
    "by finite differences: $$\\begin{eqnarray}\n",
    "\\frac{\\partial T}{\\partial t} &=& \\frac{T(x,t+\\Delta t)-T(x,t)}{\\Delta t}, \\\\\n",
    "\\frac{\\partial^2 T}{\\partial x^2} &=& \\frac{T(x+\\Delta x,t)+T(x-\\Delta x,t)-\n",
    "2T(x,t)}{(\\Delta x)^2}.\\end{eqnarray}$$\n",
    "\n",
    "Replacing in \\[heat\\_1d\\] we obtain the discrete equivalent:\n",
    "$$T(x,t+\\Delta t)=T(x,t)+\\frac{\\alpha}{C}\\left[T(x+\\Delta x,t)+T(x-\\Delta x,t)\n",
    "-2T(x,t) \\right], \n",
    "$$ with the constant $C=(\\Delta x)^2/\\Delta t$.\n",
    "We see in Fig. \\[heat\\_dif\\] that the temperature at the point\n",
    "$(x, t+\\Delta t)$ is determined by the temperatures at three points of\n",
    "the previous time step. The boundary conditions impose fixed values\n",
    "along the perimeter. The initial condition \\[heat\\_ini\\] is used to\n",
    "generate the temperature gradient at time $t=\\Delta t$, and the equation\n",
    "\\[heat1d\\_discrete\\] is used for the time evolution.\n",
    "\n",
    "The stability condition for a numerical solution is iven by\n",
    "$$\\alpha\\frac{\\Delta t}{(\\Delta x)^2} \\leq \\frac{1}{4}.$$ This means\n",
    "that if we make the time step smaller we improve convergence, but if we\n",
    "decrease the space step without a simultaneous quadratic increase of the\n",
    "time step, we worsen it.\n",
    "\n",
    "### Exercise 6.1: Finite differences program \n",
    "\n",
    "Solve the temperature distribution within an iron bar of length $l=50$cm\n",
    "with the boundary conditions $$T(x=0,t)=T(x=l,t)=0,$$ and initial\n",
    "conditions $$T(x,t=0)=100^{\\circ}\\mathrm{C}.$$ The corresponding\n",
    "constants for iron are:\n",
    "$$c=0.113 \\mathrm{cal/(g^{\\circ} C)},\\,\\, K=0.12 \\mathrm{cal/(sg^{\\circ}C)},\n",
    "\\,\\,\\rho=7.8 \\mathrm{g/cc}.$$\n",
    "\n",
    "1.  Write the program to solve the heat flow equation using the finite\n",
    "    differences method\n",
    "\n",
    "2.  Plot the temperature gradient along the bar at different instants\n",
    "    of time. Use 100 space divissions for the calculation. Choose an\n",
    "    appropiate time step such that the numerical solution is stable.\n",
    "    Check that the temperature diverges with time is the constant $c$ is\n",
    "    made larger that $0.5$.\n",
    "\n",
    "3.  Repeat the calculation for aluminum, $c=0.217cal/(g^{\\circ}C)$,\n",
    "    $K=0.49\n",
    "    cal/(g^{\\circ}C)$, $\\rho=2.7g/cc$. Note that the stability condition\n",
    "    requires you to change the size of the time step.\n",
    "\n",
    "4.  Analize and compare the results for both materials. The shape of the\n",
    "    curves may be the same but not the scale. Which bar cools faster and\n",
    "    why?\n",
    "\n",
    "5.  Pick a sinusoidal initial gradient: $$T(x,t=0)=\\sin{(\\pi x/l)}.$$\n",
    "    Compare with the analytic solution\n",
    "    $$T(x,t)=\\sin{(\\frac{\\pi x}{l})}e^{\\pi^2\\alpha t/l^2}.$$\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0871341048332\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYcAAAEPCAYAAACp/QjLAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAFCpJREFUeJzt3X+MZeV93/H3h/3Bb9issXcJUC3GQTUWaZwW12laeew4\nFrUSzF9gS6a0tfKPHcd1K6tLKpmNKtnGUWs3qvJHm2BtUEy0dVyK00hhTTwyVVyIE4wJC8VuIAWc\nHRx+rgM2P/bbP+Yu97Bn7lxm7t6598x5v6QRzz1zzplnHi33M9/nOefcVBWSJDWdNOsOSJLmj+Eg\nSWoxHCRJLYaDJKnFcJAktRgOkqSWqYVDkhuTLCW5t7Ht15Pcn+SeJF9Ocnbje9cl+U6SB5K8Z1r9\nkiSNN83K4QvA5cdtuw14S1X9PeBB4DqAJJcAVwOXDI75zSRWNZI0I1N7A66qO4Cnjtt2sKqODl7e\nCZw/aL8PuLmqXqyqh4HvAm+bVt8kSaub5V/n/xL4w0H7x4FHG997FDhvw3skSQJmFA5J/h3wQlV9\ncZXdfK6HJM3I1o3+gUn+OfBe4Ocamx8DLmi8Pn+w7fhjDQxJWoeqylr239DKIcnlwCeA91XVDxvf\nuhV4f5LtSS4EfgK4a6VzVJVfVVx//fUz78O8fDkWjoVjsfrXekytckhyM/AO4JwkjwDXs3x10nbg\nYBKAb1TVh6vqUJIDwCHgJeDDtd7fSJI0samFQ1V9YIXNN66y/6eAT02rP5Kk1857CTpqYWFh1l2Y\nG47FkGMx5FhMJl2avUnibJMkrVESap4XpCVJ3WA4SJJaDAdJUovhIElqMRwkSS2GgySpxXCQJLUY\nDpKkFsNBktRiOEiSWgwHSVKL4SBJajEcJEkthoMkqcVwkCS1GA6SpBbDQZLUYjhIkloMB0lSi+Eg\nSWoxHCRJLYaDJKnFcJAktRgOkqQWw0GS1GI4SJJaphYOSW5MspTk3sa2nUkOJnkwyW1JdjS+d12S\n7yR5IMl7ptUvSdJ406wcvgBcfty2vcDBqroYuH3wmiSXAFcDlwyO+c0kVjWSNCNTewOuqjuAp47b\nfAWwf9DeD1w5aL8PuLmqXqyqh4HvAm+bVt8kSavbusE/b1dVLQ3aS8CuQfvHgf/d2O9R4LxxJ3vy\nSfjIR+Cuu+D55+Hll09sZ/so2ZhjNsvPsF/TP2YjfsZm+V1OOQWuugo++cnJf/ZGh8MrqqqS1Gq7\nrHb8c8/Bnj1w5MiJ7Zckddn118PXvw5f/epk59nocFhKsruqDic5F3h8sP0x4ILGfucPtrXs27cP\ngAMH4MiRBWBhWn2VpE66/fZFrrlmkYsuWv85UrXqH+gTSbIH+EpVXTp4/Vngiaq6IcleYEdV7R0s\nSH+R5XWG84CvAm+q4zqX5JVNJ58ML7ywvH33bvjgB5dLKq3f0aNrP2at/3zW88/Nfq2N/Zr+MfP4\nu7z8Mtx003B6/bLLlqfcAZJQVWuapJpa5ZDkZuAdwDlJHgE+CXwGOJDkQ8DDwFUAVXUoyQHgEPAS\n8OHjg6HpySeHwQBw771wzjlT+kUkqSN+7Mfgc59bbv/lX052rqlWDifascrhwAG4+urlbdu2vToo\nJKmv/uAP4Bd/cbm9ZQu89NJyez2VQyfvJbjjjmF7587Z9UOS5sm73jVsv/wyfO976z9XJ8PhnnuG\n7T17ZtYNSZorp50G27cPXx88uP5zdTIcHnpo2L700tn1Q5LmzetfP2x/4xvrP08nw+EHPxi2L754\ndv2QpHnTDIfDh9d/nk6Gw49+NGy/4Q2z64ckzZszzxy2n356/efpZDi8+OKwvXv37PohSfPm7LOH\n7WefXf95OhkOzWconTf2CUyS1B/NKzibU/Br1blwOHr01XcOnn/+7PoiSfPmda8btp97bv3n6Vw4\nHD+HtmPHyvtJUh81w+H559d/ns6Fw6OPDtsb8chcSeqSXbuG7ebFO2vVuXB4rPGs1i1bZtcPSZpH\nzYt0Jnm0UOfC4fHHh+1t22bXD0maR81wmOQD0DoXDktLw/bJJ8+uH5I0j5oX6azn0eLHdC4cvv/9\nYdvPb5CkVzv+xuD1XrHUuXD4m78Ztk87bXb9kKR5dNJJr75Y55FH1nmeE9OdjfPUU8P2GWfMrh+S\nNK9OaryzP7biBy6/hnOcmK5snGeeGbbPOmt2/ZCkedW8WKe5TrsWnQuHI0eGbSsHSWprhkNzKn4t\nOhcOzet2Tz11dv2QpHm1deuw/bd/u75zdC4cmk9k9WolSWprhsMPf7i+cxgOkrTJNMOhN5eyGg6S\ntLpmOKz34XudC4eXXhq2XXOQpLbmgnRvppWa4WDlIEltzXDoZeXgHdKS1NYMh/U+trtz4dB8yuDp\np8+uH5I0r7ZvH7Z7M63UDAfXHCSpzcrBykGSWpofZ7DeD/zpXDg0n09u5SBJbc1w6FTlkOS6JPcl\nuTfJF5OcnGRnkoNJHkxyW5IdKx3bDAefrSRJbc01h+a9YWux4eGQZA/wS8BPV9WlwBbg/cBe4GBV\nXQzcPnjd0gwHp5Ukqa15mX+XppWeBV4ETkuyFTgN+B5wBbB/sM9+4MqVDq4atq0cJKmtOa3Umcqh\nqp4E/gPw/1gOhaer6iCwq6qOPXl8Cdi18vHD9plnTrWrktRJzfXY9VYOW8fvcmIluQj4V8Ae4Bng\nvyX5YHOfqqoktcLhVO17pX3PPQtcdNHCtLoqSZ30xBOLwCKw/s9z2PBwAP4B8CdV9QRAki8DPwMc\nTrK7qg4nORd4fOXD973SuvzyaXdVkrrnTW9aABaA5bXZZ5/9tTWfYxZrDg8Ab09yapIA7wYOAV8B\nrh3scy1wy7gT+WwlSWprTis1Hzm0FhteOVTVPUl+B/gmcBT4c+C/AGcCB5J8CHgYuGrcuU7q3F0a\nkjR9zefONW8cXotZTCtRVZ8FPnvc5idZriIkSRM4EZVDZ//2TmbdA0maT817wNZbORgOkrTJNCuH\n5o3Da9HZcHC9QZJW1rxBuHeVg+EgSStrhoOVgyQJePXVSrXi7cTjdfYtdsuWWfdAkuZT89FCvasc\nDAdJWlnzaiUrB0kScGIeStrZcNg6k9v3JGn+nXXW5OfobDhYOUjSypqfBLdenQ2Hbdtm3QNJ2rw6\nGw5OK0nSaJM+RcJwkCS1dDYcXHOQpNF6WzkYDpI0muEgSWrpbTi45iBJo/U2HKwcJGk0w0GS1NLb\ncHBaSZJGMxwkSS2TfuaN4SBJm5CVgySppbeVgwvSkjRab8PBp7JK0mi9DQenlSRptN6Gg5WDJI1m\nOEiSWjoZDkl2JPlSkvuTHEryD5PsTHIwyYNJbkuyY7VzuCAtSaNNLRySTPNv8/8E/GFVvRn4SeAB\nYC9wsKouBm4fvB7pRHxGqiRtVpP+Ab1attw52alXluRs4J9U1Y0AVfVSVT0DXAHsH+y2H7hytfM4\nrSRJo00zHCa8v26kC4HvJ/lCkj9P8l+TnA7sqqqlwT5LwK7VTmI4SNJok04rrXZB6OuT/GtWDomq\nqv84wc/8aeCXq+pPk3ye46aQqqqS1MqH7wPgzjthcXGBhYWFdXZDkjanxcVFnn56caJzrBYOW4Az\nJzr7yh4FHq2qPx28/hJwHXA4ye6qOpzkXODxlQ/fB8DCwvKXJOnVFhYWeMMbFnjqqWNbfm3N51gt\nHA5X1drPOMbgzf+RJBdX1YPAu4H7Bl/XAjcM/nvLaudxQVqSRpt0zWFW9xl/FPjdJNuB/wv8C5Yr\nlQNJPgQ8DFy12gkMB0kabZrh8O7JTj1aVd0DXDbJzzQcJGm0qV2tVFVPTHbq6TIcJGm0aV7KOtd8\n8J4kjTbpe2Rnw8HKQZJG6204nHzyrHsgSfOrt9NKVg6SNJqVgySppbfh4LOVJGm03oaDlYMkjTbp\nH9CGgyRtQr1dkDYcJGk0KwdJUktv1xxckJak0XpbOZxyyqx7IEnzq7eVg9NKkjSalYMkqcVwkCS1\nTPqIoc6Gg89WkqTRels5+HkOkjRabxekJUmj9XZaSZI0muEgSWrpZTgks+6BJM031xwkSS1WDpKk\nlkmfImE4SNIm1MvKQZK0ul6Gg5WDJK2ul3dIn9TJXkvSxunsmkOSLUnuTvKVweudSQ4meTDJbUl2\njD524/opSV3U2XAAPgYcAmrwei9wsKouBm4fvF6R4SBJq+tkOCQ5H3gv8FvAsbf6K4D9g/Z+4MrR\nx0+1e5LUeV1dc/gc8AngaGPbrqpaGrSXgF2jDnbNQZJWN+ln3mz4g6+T/ALweFXdnWRhpX2qqpLU\nSt+DfbzwAuzbBwsLCywsrHgKSeqtxcVFvvSlxYnOkaoR78FTkuRTwDXAS8ApwFnAl4HLgIWqOpzk\nXOBrVfV3jzu2oDjrLHjmmQ3ttiR1ykMPwRvfeOxVqKo1Tchv+ARNVf1qVV1QVRcC7wf+uKquAW4F\nrh3sdi1wy6hzuOYgSavr6ppD07HS5TPAzyd5EHjX4PWKDAdJWt2kT2Xd8GmlSRybVtq5E554Yta9\nkaT59eST8LrXHXvVgWmlE8HKQZJW18tnK3kpqyStznCQJLX08pPgDAdJWt2k75OdfJt1zUGSpquT\n4WDlIEnT1cm3WcNBkqark2+zhoMkTVcn32YNB0kab5L12U6+zW7ZMuseSNLm1slw8GolSZquToaD\n00qSNJ7TSpKkE6qT4WDlIEnjWTlIkloMB0lSS+/CwauVJGm83oWDlYMkjWc4SJJOKMNBkjYpKwdJ\nUovhIElqmeSesE6GgzfBSdJ4Vg6SpBbDQZLU0rtw2Lp11j2QpPnXuzUHKwdJGq93lYPhIEnjdapy\nSHJBkq8luS/JXyT5lcH2nUkOJnkwyW1Jdow6h+EgSeN1KhyAF4GPV9VbgLcDH0nyZmAvcLCqLgZu\nH7xekeEgSeN1alqpqg5X1bcG7R8A9wPnAVcA+we77QeuHHUOF6QlabyuVQ6vSLIHeCtwJ7CrqpYG\n31oCdo06znCQpPE6VTkck+QM4PeBj1XVkeb3qqqAGnWs4SBJ401SOczkbTbJNpaD4aaqumWweSnJ\n7qo6nORc4PGVj97Hn/0Z7NsHCwsLLCwsbESXJakzFhcXWVxc5MiR8fuOkuU/0jdOkrC8pvBEVX28\nsf2zg203JNkL7KiqvccdW1B89KPwG7+xod2WpM7Zswf+6q8AQlWtaZJpFpXDzwIfBL6d5O7BtuuA\nzwAHknwIeBi4atQJnFaSpPEmubJzw99mq+p/MXqt492v5RyGgySN19mrldbLcJCk8Tp5tdIkDAdJ\nGm+SaSXDQZI2qd5NK/n4DEkar3eVw7Zts+6BJM2/3lUOTitJ0nhWDpKkFisHSVJL78LBykGSxuvd\ntJJXK0nSeL0LBysHSRrPcJAktfQuHFyQlqTxehcOVg6SNJ7hIElq6d2lrE4rSdJ4k7xXdjIcrBwk\naTynlSRJLb0Lh+3bZ90DSZp/TitJklp6VzkYDpI0Xu/CwauVJGm83oWDaw6SNJ5rDpKklt6Fg5WD\nJI1nOEiSWgwHSVJL78LBq5UkabxNc7VSksuTPJDkO0n+7aj9DAdJGm9TVA5JtgD/GbgcuAT4QJI3\nr7TvJI+h3SwWFxdn3YW54VgMORZDjsUmCQfgbcB3q+rhqnoR+D3gfTPu09zyH/6QYzHkWAw5Fpsn\nHM4DHmm8fnSwTZK0DpPcEzZP4VCz7oAkbSaTVA6pmo/35CRvB/ZV1eWD19cBR6vqhsY+89FZSeqY\nqspa9p+ncNgK/B/g54DvAXcBH6iq+2faMUnqobm5KLSqXkryy8AfAVuA3zYYJGk25qZykCTNj3la\nkF7Va71BbjNKcmOSpST3NrbtTHIwyYNJbkuyY5Z93ChJLkjytST3JfmLJL8y2N678UhySpI7k3wr\nyaEknx5s791YHJNkS5K7k3xl8LqXY5Hk4STfHozFXYNtaxqLToTDWm6Q26S+wPLv3rQXOFhVFwO3\nD173wYvAx6vqLcDbgY8M/i30bjyq6ofAO6vqp4CfBN6Z5B/Tw7Fo+BhwiOHVj30diwIWquqtVfW2\nwbY1jUUnwoGe3yBXVXcATx23+Qpg/6C9H7hyQzs1I1V1uKq+NWj/ALif5fth+joezw2a21leq3uK\nno5FkvOB9wK/BRy7MqeXYzFw/NVJaxqLroSDN8i17aqqpUF7Cdg1y87MQpI9wFuBO+npeCQ5Kcm3\nWP6dv1ZV99HTsQA+B3wCONrY1texKOCrSb6Z5JcG29Y0FnNztdIYrpqvoqqqb/eAJDkD+H3gY1V1\nJBn+kdSn8aiqo8BPJTkb+KMk7zzu+70YiyS/ADxeVXcnWVhpn76MxcDPVtVfJ3k9cDDJA81vvpax\n6Erl8BhwQeP1BSxXD322lGQ3QJJzgcdn3J8Nk2Qby8FwU1XdMtjc2/EAqKpngP8J/H36ORb/CLgi\nyUPAzcC7ktxEP8eCqvrrwX+/D/x3lqfm1zQWXQmHbwI/kWRPku3A1cCtM+7TrN0KXDtoXwvcssq+\nm0aWS4TfBg5V1ecb3+rdeCQ559gVJ0lOBX4euJsejkVV/WpVXVBVFwLvB/64qq6hh2OR5LQkZw7a\npwPvAe5ljWPRmfsckvxT4PMMb5D79Iy7tGGS3Ay8AziH5bnCTwL/AzgA/B3gYeCqqnp6Vn3cKIOr\ncb4OfJvhdON1LN9R36vxSHIpywuLJw2+bqqqX0+yk56NRVOSdwD/pqqu6ONYJLmQ5WoBlpcOfreq\nPr3WsehMOEiSNk5XppUkSRvIcJAktRgOkqQWw0GS1GI4SJJaDAdJUovhIElqMRwkSS2GgzSBJJcl\nuSfJyUlOH3wA0SWz7pc0Ke+QliaU5N8DpwCnAo9U1Q0z7pI0McNBmtDgKbHfBJ4Hfqb8n0qbgNNK\n0uTOAU4HzmC5epA6z8pBmlCSW4EvAm8Ezq2qj864S9LEuvJJcNJcSvLPgB9V1e8lOQn4kyQLVbU4\n465JE7FykCS1uOYgSWoxHCRJLYaDJKnFcJAktRgOkqQWw0GS1GI4SJJaDAdJUsv/B5JZDtMJQxZa\nAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1062c5f10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "import numpy as np\n",
    "from matplotlib import pyplot\n",
    "import matplotlib.animation as animation\n",
    "from JSAnimation.IPython_display import display_animation\n",
    "\n",
    "K = 0.12\n",
    "c = 0.113\n",
    "rho = 7.8\n",
    "alpha = K/c/rho\n",
    "\n",
    "l = 50    #length of the string\n",
    "dx = 0.125  # space step\n",
    "nx = int(l/dx)  # number of points in space\n",
    "x = np.arange(0,l+dx,dx) # the +1 is necessary to store the value at l\n",
    "dt = 0.01\n",
    "C = dx**2/dt\n",
    "r = alpha/C\n",
    "\n",
    "print r\n",
    "\n",
    "t0 = np.zeros(nx+1)\n",
    "t1 = np.zeros(nx+1) # these arrays will contain the new displacements at t, and t+delta\n",
    "\n",
    "#Initial conditions\n",
    "t0[:] = 100.\n",
    "t0[0] = 0.\n",
    "t0[nx] = 0.\n",
    "\n",
    "\n",
    "fig = pyplot.figure()\n",
    "ax = pyplot.axes(xlim=(0, l), ylim=(0, 120), xlabel='x', ylabel='T')\n",
    "points, = ax.plot([], [], marker='', linestyle='-', lw=3)\n",
    "\n",
    "def evolve(i):\n",
    "    global t0, t1\n",
    "\n",
    "    for ix in range(1,nx):\n",
    "        t1[ix] = t0[ix] + r*(t0[ix+1]+t0[ix-1]-2*t0[ix])  \n",
    "\n",
    "    points.set_data(x, t1)\n",
    "\n",
    "    for ix in range(nx):\n",
    "        t0[ix] = t1[ix]\n",
    "\n",
    "    return points\n",
    "\n",
    "\n",
    "#anim = animation.FuncAnimation(fig, evolve, frames = 2000, interval=10)\n",
    "\n",
    "#display_animation(anim, default_mode='once')\n",
    "\n",
    "for i in range(20):\n",
    "    evolve(i)   \n",
    "pyplot.plot(x, t0, color='blue', ls='-', lw=3);\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "#### Challenge 6.1:\n",
    "\n",
    "Solve Exercise 6.2, part 5)\n",
    "\n",
    "### Exercise 6.2: Two bars in contact \n",
    "\n",
    "Two identical bars, 25cm long each, are incontact. One bar is initially\n",
    "at $50^{\\circ}C$, and the other at $100^{\\circ}C$. The free ends are\n",
    "kept at $0^{\\circ}C$. Calculate and plot the temperature distribution as\n",
    "a function of time.\n",
    "\n",
    "\n",
    "#### Challenge 6.2:\n",
    "\n",
    "Modify the initial conditions in the program to solve Exercise 6.2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}