{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# PDE python exercise 4"
   ]
  },
  {
   "attachments": {
    "image.png": {
     "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAEWCAYAAABG030jAAAAAXNSR0IArs4c6QAAAERlWElmTU0AKgAAAAgAAYdpAAQAAAABAAAAGgAAAAAAA6ABAAMAAAABAAEAAKACAAQAAAABAAABbqADAAQAAAABAAABFgAAAAC6i9fzAAAla0lEQVR4Ae2dCZQU1bnHvxFkERAImyg7grIpCIIgsgURgk8FBEEDgo8XlWfUHDVK8AkqScTlaGLcokEhJCCIEgMIqGwKIoIsEXTYF5FNAVmFGaZffYU9ds9091R1V01Xdf3uOT1dde+tu/xuzb9vfXWXrJDhBAcBCEAAAn4hMOYsv5SUckIAAhCAwBkCCDd3AgQgAAGfEUC4fdZgFBcCEIAAws09AAEIQMBnBBBunzUYxYUABCBgWbhPnz4tEydOlFGjRsny5ctjktu0aZPs3r07ZhieEIAABCDgDAHLwn3PPfdIgwYN5MYbb5SOHTvKhg0bokpw5MgRadmypYwcOTLKnxMIQAACEHCWgCXhnjZtmrRv394U7JMnT0pOTo4sXbo0qiRLliyRY8eOSatWraL8OYEABCAAAWcJWBLuqVOnysCBA82cZ86caX43a9YsqiSLFi0yz7t27RrlzwkEIAABCDhLwJJw33TTTVKiRAnRSZaTJk2Siy66SC6//PKokixYsECqVKkiLVq0iPLnBAIQgAAEnCVgSbjVrq1OzSPbt2+XwYMHR5Xi6NGjsnLlSunUqZNkZWVFhXECAQhAAALOErAk3OEsp0+fbh6GhTzsr/bt3NxcwUwSJsI3BCAAAfcI2BLu+fPnS506dUxTSWSRFi5caJ526dIl0ptjCEAAAhBwgYAt4c7OzjaH/BUshwp31apVpXnz5gWDOIcABCAAAYcJ2BLusmXLSvny5aOKsHHjRlmxYoV07twZ+3YUGU4gAAEIuEPAlnDrkMBly5bJDz/8YJZmy5YtMmDAANO+jZnEnQYiVQhAAAIFCWTZ2UhBJ9889thjsnXrVjnvvPPk3HPPlV27dslrr70m69atk6ZNmxZM39J5Xl6enHWWrd8QS+kSCQIQgEAGEhhjWbhVXHWoX+RwPx3X3ahRIylVqpSsX78+aT5jx46Vhx9+OOnruRACEIBAgAhY2wFHF4+qWbOmdO/ePYqNzqLcvHmzPPLII1H+dk50XPiYMWPM8eF2riMuBCAAgaASsGSfeO+992Tfvn1Su3btfE6HDx+We++9V4YMGZI/HT4/0MbBhAkTRFce1G8cBCAAAQgUTcCScLdt21Z02rvastXt2LFDevXqZX7Gjx9fdC5xYqipJSzY+q3nOAhAAAIQSEzAso1b1yJ5//33RW3d+hk0aFDKKwEuXrzYHEYYLqIuVKXT5nEQgAAEIBCXwJiScYMKBOh0dqentL/xxhtRueg5wh2FhBMIQAAChQhY7nEXujJFD127W4cU6gJVYaeTe/bs2SPlypULe/ENAQhAAALRBKyNKom+xpkzXbAqUrQ1VT0PL2TlTC6kAgEIQCDzCFh6OelGtcNmkvBIlfB32N+NPEkTAhCAQCYQSItw69htXZhK1z4JLxGr33qu/hqOgwAEIACB2ATSItzhoX99+/Y1p81r0XT6vJ5HDhGMXWR8IQABCASbQFqEe+LEiSb1oUOHRtEPn4fDowI5gQAEIAABk0BahHv27Nny+OOPS7du3cyetu5fqT1uPVd/DcdBAAIQgEBsAmkR7saNG5uLSumKgEeOHJHPPvvM/NZzXWxKw3EQgAAEIBCbQFqEO3ZR8IUABCAAASsEEG4rlIgDAQhAwEMEEG4PNQZFgQAEIGCFAMJthRJxIAABCHiIAMLtocagKBCAAASsEEC4rVAiDgQgAAEHCMyYMcOBVEQsL+vqSG4kAgEIQCBABHJycmTFihXy5ZdfyqRJk+TTTz8VXRk1VUePO1WCXA8BCEAgDoFTp07JunXrzLkp7du3jxPLvjfCbZ8ZV0AAAhCwRED3Fhg+fLh07NhRSpQoYekaK5EQbiuUiAMBCEDAQwQQbg81BkWBAAQKE9A9bp988kl56qmn5Nprr5Xc3Nz8SN9//725afm3336b7xeEA15OBqGVqSMEfExg/PjxpmA3aNBARo4cKWvWrJHWrVubNZozZ47oRuYVKlRIuYa6+5aum2TH6VLUbdu2tXOJI3ERbkcwkggEIOAWgezsbNNOrEPpSpYsKQ0bNszPav78+XLFFVdI6dKl8/30YNmyZeZLwTJlysgNN9xgaR/bli1bSuXKlaPSKeqkfv36RUVxJRzhdgUriUIAAk4ReOCBB8yk/v73v5tmkUqVKuUnrTtmDRo0KP9cD5544gnzReD9999v7qbVp08fefPNN4sUZf1BiPxRiErUYyfYuD3WIBQHAhCIJlC9enVz7LOu0z9w4MD8wN27d8uGDRuka9eu+X779u2TP/3pT3LfffdJVlaW1KtXz+yRP/vss/lxMuGAHncmtCJ1gECGE/joo49Ex0T37Nkzv6bhfWt1fPTBgwelVKlSoqaT888/X3Rt/7C78MIL5fnnn5fHHnss7BXz+y9/+Yt8+OGHMcNieeoPwx133CE9evSIFeyqH8LtKl4ShwAEnCCwbds2qV27tlSsWDE/ucWLF0ubNm1MwX7llVdk8ODBovEK2rt1E/IdO3bkXxfvYNiwYaJmFTtOnwbS4RDudFAnTwhAwBYBtT0fP37c3Exce7o6G3HevHny85//3ExHzSZq+9ZeuYZHOu19q39RTifL6Mctp0MWtRy661eqo2AQbrdaiXQhAAHHCFx99dUyYsQIufvuu0VNH+ecc47MnDlTRo0aJY8++qj079/fzEvFO3Kct3qqWEb21B0rlMWEtNy6Zsn+/fvNYY233367+aJUhxKGf3gsJpUfDeHOR8EBBCDgZQJjxowpVLy33347yk+H9OnLyUi3d+9eUf90uRdffNHxrH+y4DueNAlCAAIQKF4CV155pdkb37hxY37G+sJSe7mZ5NLe4966davJU2c/hUKhTGJLXSAAgTgEunTpIvpx2ulCTjpR55lnnpFWrVrJrl27pHfv3ub4b6fzSmd6aRdufQusbtGiRebHPOEPBCCQ8QTcEG6Fpi8y3TBPeKlB0i7cYRh169aVoUOHhk/5hgAEMpCAjr3WThouNQKeEW6d4RTr5UNq1eNqCEDASwT0fxzhTr1FeDmZOkNSgAAEMoWAsc6J3HSTyHvviRjLyXrVJS3cN998s+g0VBwEIACBjCFw771iLH4ixuBwkQsuEHnwQTE2jPRc9ZIS7lmzZsnkyZPNAeWeqxEFggAEIJAsAWMZWGMBEl0XVoyFvsWYhinGLBmRdu1EXnpJ5NChZFN29DrLNu6TJ0/Kb3/7W3OqKTYqR9uAxCDgLgHtPRq2ZS84LYV+zB6tlstrTocoGi9QCzljSr0sXy7G9E0xxi0XCi5uD8s9bh0f2alTJ7nzzjtlypQpjm58WdyVJj8IBIrA6NFnxEYFJ82fMUZZjD6s6He6yxIzf+1lhxmpicTorErNmj/1uI1VCL3gLPe4deeJfv365ZdZl1DEQQACEMgoAj/8IPL66yJvvCGyc6fIkCFirBUrcvHFnqqmZeH2VKkpDAQgAAE3CDz3nMjq1WdMOcbCVoZpwY1cUk4T4U4ZIQlAAAIZQ+Chh3xRFc8It05918H5OAhAIHMJ6MxJXOoEXBFuXR/Xrghv377dXFc39SqRAgQgAIHMJuCKcI823hjrx4rT9UkmTJggnTt3dmW1MCtlcDvO9OnTo17sup0f6UPA6wTcWmDK6/V2qnyuCLedwtWvX9+Mrjs1WxV7O+l7Ie7YsWPl4YcfFh2Zg4MABCCQKgHL47gLZqT7uunYblzRBHTPuzwPr3tQdA2IAQEIeIlA0l3AY8eOeakeni4Lwu3p5qFwEPAdgaR73L6raRoLrE8mp0+fTmMJyBoCEMgkAgh3MbQmPe5igEwWEAgQAYS7GBpbhRsHAQhAwCkCKIpTJBOkc+rUKUaUJOBDEAQgYI8Awm2PV1Kxc3Jy5Oyzz07qWi6CAAQgUJAAwl2QiAvn+mKSMdwugCVJCASUAMLtcsNrbxvRdhkyyUMgYAQQbpcbPDc3F+F2mTHJQyBoBBBul1tcX0yWLVvW5VxIHgIQCBIBhNvl1tYZpmV0A1IcBCAAAYcIINwOgYyXjAp3uXLl4gXjDwEIQMA2AYTbNjJ7FyDc9ngRGwIQKJoAwl00o5RiINwp4eNiCEAgBgGEOwYUJ70QbidpkhYEIKAEEG6X74Pjx49L+fLlXc6F5CEAgSARQLhdbu0DBw5IpUqVXM6F5CEAgSARQLhdbu1Dhw4h3C4zJnkIBI0Awu1yix88eFAqV67sci4kDwEIBIkAwu1yayPcLgMmeQgEkADC7XKjI9wuAyZ5CASQAMLtcqNj43YZMMlDIIAEEG6XG11HlVStWtXlXEgeAhAIEgGE2+XW3rNnj1SrVs3lXEgeAhAIEgGE2+XW3rt3r9SoUcPlXEgeAhAIEgGE28XW1unu6lgd0EXIJA2BABJAuF1s9H379kn16tVdzIGkIQCBIBJAuF1sdYTbRbgkDYEAE0C4XWz8/fv3S82aNV3MgaQhAIEgEkC4XWz1nTt3Itwu8iVpCASVAMLtYsvv2rVLLrjgAhdzIGkIQCCIBBBuF1sd4XYRLklDIMAEEG4XGx/hdhEuSUMgwAQQbhcb/+uvv5ZatWq5mANJQwACQSSAcLvY6vS4XYRL0hAIMAGE26XGP3r0qOTl5bH7jUt8SRYCQSaAcLvU+lu3bpU6deq4lDrJQgACQSaAcLvU+tu2bZN69eq5lDrJQgACQSaAcLvU+trjrl+/vkupkywEIBBkAgi3S61Pj9slsCQLAQgIwu3STYBwuwSWZCEAAYTbrXsA4XaLLOlCAAL0uF26B7Zv3y4NGzZ0KXWShQAEgkwA4Xah9XU5V3WVK1d2IXWShAAEgk4A4XbhDti4caM0atTIhZRJEgIQgIBg43bjJkC43aBKmhCAQJgAPe4wCQe/EW4HYZIUBCBQiADCXQhJ6h4Id+oMSQECEIhPAOGOzybpEIQ7aXRcCAEIWCCAcFuAZDfKpk2bpHHjxnYvIz4EIAABSwQQbkuYrEfSDYLLly8v5557rvWLiAkBCEDABgGE2wYsK1HXr18vTZs2tRKVOBCAAASSIoBwJ4Ut/kUId3w2hEAAAs4QQLid4ZifypdffilNmjTJP+cAAhCAgNMEEG6HidLjdhgoyUEAAoUIINyFkKTmgXCnxo+rIQCBogkg3EUzshxj7969UrJkSalWrZrla4gIAQhAwC4BhNsusQTxV69eLS1atEgQgyAIQAACqRNAuFNnmJ/C2rVr5ZJLLsk/5wACEICAGwQQbgepItwOwiQpCEAgLgGEOy4a+wEIt31mXAEBCNgngHDbZxbzipycHMnOzpbmzZvHDMcTAhCAgFMEEG6HSH711VdSv359KV26tEMpkgwEIACB2AQQ7thcbPuuWbNGWrVqZfs6LoAABCBglwDCbZdYnPgrVqxAuOOwwRsCEHCWAMLtEM+VK1dK69atHUqNZCAAAQjEJ4Bwx2djOSQvL09WrVqFcFsmRkQIQCAVAgh3KvR+vHbDhg1So0YNqVixogOpkQQEIACBxAQQ7sR8LIViJrGEiUgQgIBDBBBuB0Ai3A5AJAkIQMAyAYTbMqr4EXVESdu2beNHIAQCEICAgwQQ7hRh5ubmyueffy5t2rRJMSUuhwAEIGCNAMJtjVPcWLo+ic6YrFChQtw4BEAAAhBwkgDCnSLNTz/9VNq1a5diKlwOAQhAwDoBhNs6q5gxEe6YWPCEAARcJIBwpwgX4U4RIJdDAAK2CSDctpH9dMHhw4dl165dbFf2ExKOIACBYiCAcKcAecmSJeZokqysrBRS4VIIQAAC9ggg3PZ4RcX++OOPpWPHjlF+nEAAAhBwmwDCnQJhhDsFeFwKAQgkTQDhThLdqVOnRGdMdujQIckUuAwCEIBAcgQQ7uS4ia5PctFFF0n58uWTTIHLIAABCCRHAOFOjptgJkkSHJdBAAIpE0C4k0S4dOlSueqqq5K8mssgAAEIJE8A4U6Cne54s3DhQunUqVMSV3MJBCAAgdQIINxJ8Fu9erXUrFnT3PUmicu5BAIQgEBKBBDuJPAtWLBAunbtmsSVXAIBCEAgdQIIdxIM1UzSpUuXJK7kEghAAAKpE0C4bTI8ffq0LF68mB63TW5EhwAEnCOAcNtkqbvd1KlTR6pWrWrzSqJDAAIQcIYAwm2T46JFi6Rbt242ryI6BCAAAecIINw2Wc6ZM0e6d+9u8yqiQwACEHCOAMJtg+WJEyfkk08+ocdtgxlRIQAB5wkg3DaYqplEd3MvV66cjauICgEIQMBZAgi3DZ7z5s2THj162LiCqBCAAAScJ4Bw22CKcNuARVQIQMA1Agi3RbTffPON7N69W1q3bm3xCqJBAAIQcIcAwm2R69y5c+Waa66Rs84CmUVkRIMABFwigApZBPvuu+9K7969LcYmGgQgAAH3CCDcFtjqNmXz58+Xnj17WohNFAhAAALuEkC4LfDVRaVatGghVapUsRCbKBCAAATcJYBwW+A7a9YszCQWOBEFAhAoHgIItwXOM2fOlGuvvdZCTKJAAAIQcJ8Awl0E4+zsbNGlXNVUgoMABCDgBQIIdxGtMH36dLnuuuuKiEUwBCAAgeIjgHAXwfqdd96RPn36FBGLYAhAAALFRwDhTsB6586dsnXrVnZzT8CIIAhAoPgJINwJmGtvW80kJUqUSBCLIAhAAALFSwDhTsAbM0kCOARBAAJpI4Bwx0F/4MABWbVqlbk+SZwoeEMAAhBICwGEOw72t956y5ziXqpUqTgx8IYABCCQHgIIdxzuU6ZMkYEDB8YJxRsCEIBA+ggg3DHY79mzxzST9OrVK0YoXhCAAATSS6BkerP3Zu5Tp06V66+/XkqXLu3NAlIqCEAgJgF9N/XMM89ImTJlJC8vT77//nu5//775fzzz48Z36+eCHeMllMzySOPPBIjBC8IQMCrBI4cOSJXXXWVPP/889KtWzezmGvWrDHnYSxevDijxBtTSYG7cMeOHbJhwwbp3r17gRBOIQABLxMYN26cVK9ePV+0tayXXnqpdOjQQX73u995uei2y4ZwF0D2z3/+UwYMGCAlS/IwUgANpxDwLIFQKCSvv/66tG3btlAZ1U/Nn8eOHSsU5lcPhLtAy40fP15uvfXWAr6cQgACSkDNEV999ZWoUHrJbdmyRXRD72rVqhUq1nnnnScnTpyQlStXFgrzqwfCHdFyy5Ytk6ysLGnXrl2EL4cQgIASWLJkiTz33HPy5ptvyuDBg6Og6Lo++gJQBTQdLpxvhQoVCmV/zjnnmH6bN28uFOZXD+wBES03YcIEetsRPDiEQJiA7rs6bdo0U7jV7PCHP/xBJk2aFA4W3SXqu+++K/IF4A8//CAjRowQTc+Oa968uTz00ENxLzl06JAZFmtdobBfOE7cRHwUgHD/2FgnT540exJr1671UfNRVAgUDwHdBeqGG24wM/vXv/4lrVq1isp4wYIF5pOqDsNL5DRch+fp/5sdV7Vq1YTRdbMTdfrEXNCFzTrhOAXD/XiOcP/YanpjtmnTRmrVquXHdqTMEHCVQKdOnczNso8ePSozZsyQp59+Oio/3VD7V7/6VZRfvJOmTZvGC0raPzznQsduF3Rhwc6k5SsQ7h9b+a9//asMHTq0YJtzDgEIGATCPV7t4Ghv+cYbb8znsn79etm3b1/UMLz8wGI6CL+UVFNMQRf2C8cpGO7Hc4TbaDXdLOHzzz+Xd999149tSJkhUGwE5s6dK61bt44avbFo0SJzpqKOl87NzTU/8UwmKvojR460bSpp0qSJ3HXXXXHrqTZwNZMcPHiwUJy9e/eafs2aNSsU5lcPhNtoOe1t6xDA8OOWXxuTckPAbQI6FLCgfVtnJaqZUf9/dI/WBg0aFIoTLpeaK7S3btfGrUP6ErlKlSpJ+/btJdbIke3bt0vt2rUzasPvwAt3Tk6OOXBfbz4cBCCQmECVKlWiJqdt27ZN5syZI7fccot54fLly6Vv375xE9FesfbM3XD33XefOWJFTSPhHr++mJw9e7b85je/ifni0o1yFEeagRdufdGiL0saN25cHLzJAwK+JjBq1CgZPny4vPfee+ZknHXr1smTTz4pun69rvGj/0exRnYUR6X1B+M///mPDBs2TP785z+bTwA61V2fBu65557iKEKx5RF44X711Vfl9ttvLzbgZAQBPxNQc4T2qnXxposvvthcHkLro4s66exEtTWn040ePdqc2Tl58mTR0SS33XabXHbZZekskit5B1q49df5iy++kH79+rkCl0Qh4HsCTzwh8vbbIh99JEYX1qxOuXLlCpk7GjZs6Jmq6g+KfjLZBVq4dfquvqlmQalMvsWpW9IEDHOIYXMQWbEiX7STTosLHSUQWOHev3+/0ZF4WzZt2uQoUBKDQEYQ2LhRZMgQMcbIijGPPSOqlEmVCOwiUy+//LL079/fnA2WSQ1KXSCQMgFjBUC59loRNZMYNm2c9wgEssetQwBVuD/44APvtQglgkA6CehyrYMGiVx9tch//3c6S0LeCQgEssetC663bNlSdDYWDgIQiCCgW/YdPizGMoARnhx6jUDgetw6ROgJ4xEwcklKrzUK5YFAWgjo6BFjaWNjxwEx3tinpQhkao1A4FpHJwnUqVOn0HAma7iIBYEMJWAMizWW9xN5/30xFiLJ0EpmTrUCJdw6/VUXgNdhgDgIQOBHAgcOiFx3nRjbo4uxyAhYfEAgUDbuf//736KTB67WFy84CEBAjOmFYgyvOvPRl5I4XxAITI9bF1jXdQvGjRvni4ahkBAoFgIPPHDGnv3HPxZLdmTiDIHACPc//vEP0aUfe/fu7Qw5UoGA3wnonpHGNmTmy8izAvXw7feWk0AIt47bfsQY5jRx4kTfNxgVgIAjBHQau7HUqRibIBg9GkeSJJHiIxCIn1ndKEHHbF911VXFR5acIOBVArojTJ8+Iq+9Jsaaxl4tJeVKQCDje9zHjx83R5LoXnk4CASegPH0aYr2//yPyPXXBx6HXwFkfI9b1+ft2bNn3K2U/NpwlBsCSRH43/8V0W3A/u//krqci7xBIKN73Bs2bDC3JdNdqHEQCDwBY30eWbJE5NNPxdimJvA4/Awgo4VbtyvSrZaqV6/u5zai7BBIncDHH4vxhl7kk09EypdPPT1SSCuBjDWVqE1bd3f+9a9/nVbAZA6BtBP4+msxtlYXMYbEiod2qkk7Fx8XICN73EeM9YTvvPNOcyEpdrfx8d1J0VMnYOx4bk5n14k2zBhOnadHUsjIHrfOkOzVq5d07tzZI5gpBgTSREDX1NYhf/fdl6YCkK0bBDKux62bI/zLmA2mGwHjIBBoAs88I/LllyJLlwYaQyZWPqOE+7CxAPzw4cNl/PjxUrFixUxsL+oEAWsEdHenp54S+ewzkTJlrF1DLN8QyCjhHjZsmDG3oI90797dNw1AQSHgOIEtW0Ruvllk+nSR2rUdT54E008gY4T7pZdekq1bt8rkyZPTT5USQCBdBI4eFfmv/xJ57DEx1nhIVynI12UCGSHcnxhjU3URqeXLl0upUqVcRkbyEPAoAd3od/BgkY4dRe64w6OFpFhOEPC9cO/cudMYonqjufJf/fr1nWBCGhDwJ4HHHxfZt0/kzTf9WX5KbZmAr4VbX0bq+toPGGNUdfgfDgKBJfDuuyLGKpjm2to8dWb8beBb4T5q2PKuvPJKueaaa+Tee+/N+IaighCIS0CH/N12m8icOSI1asSNRkDmEPDlBBxdqrVChQryy1/+Up5++unMaQ1qAgG7BA4dOvMy8tlnRdq0sXs18X1KwFaPW1fbmzFjhixbtswUzh49esgtt9xSrFVX88iAAQNEh/49+OCDxZo3mUHAUwSMfVTlppvOTGnXl5K4wBCw3ONeuXKlXH755TJr1ixp27atlC1bVkaMGGEKaHHR2rx5s7Rr185cW/tvf/tbcWVLPhDwJoGRI0Vyc89MtPFmCSmVSwQs97hvvfVW06asq+6d9ePGojcZv/ZXGwvXDBw40LQ1u1RGM9kpU6aYtuyxY8easyPdzIu0IeB5Asb/g0ydeuZlZIkSni9ucRZQLQO/+MUvZMiQIaK6Vbdu3eLMvljystTjzs7OlnXr1pmmibBoa+m6du1q9oCnTZvmWmH3GvvjDRo0SB599FGzt69T2nEQCDSB1atF7rpLREeS/OxngUYRq/JTjR80fTrX3a90iHC3bt3M4cLHjh2LFd2XfpaE+4svvpCzzz5bOnToUKiSHY3B/m7sMLNr1y7zh6JZs2bSuHFjWbVqlbRu3bpQ/nhAIFAEvv32jE37lVdEWrQIVNWtVlZXB503b54x6/9mY5mWMrJgwQKz532esWXbbcbom8WLF0tIJyv52FkylRw4cMD4Yf+ZKd4F61q1alXR8FTd+++/L18bC77v3r1b9FFHJ9a0bNnSfNzRXv64ceNSzYLrIeBrAlmnT8utEybIdqMXucDoTIl+cIUI6IgzHcTQqFEjYwLpHfLWW2+ZeqJDiF9//XXz06BBA1PMfWtKMX55inQvvPBCqE6dOjHjGcPxQgagqLAxY8boz5nlT6VKlSzHtZMuca23Aay8z+ovxv/UTONj7BbJ/0sCBsYgCst8SpQoEXr88cdDp0+fjtIwj5+MttTjrlKlihzS8aIxnO42o73xSKe2Jf1YdQsXLpRFixZZir527Vqjo/GF+Rhk5QLiJ6YEH3/wafn559LR2Oj3HsNc2ch456RmACsuaO2rTLTHrS8n1WnPO9zjNj1+/BOIHrdhEzJ/wfbs2VPoh8hYRjWkn+Jymtcbb7xhOTviJ0YFH4/y6d8/FPrggzOF++STUKhq1VAoO9v8X+P+T9xm2ns2TK8hY45JyBi2nN/7Ll++fMiY/xEyOomhvLy8xIl4O3S0GumLdDk5OSGjVx168cUXo+J+9913ITVzvPrqq1H+bp0YNqqQ8WsaOnjwoKUsiJ8YE3w8ymf37pDxjxUKnTwZCu3aFQrVrBkKzZ4dor0St1c4VE0fRsfa/GRlZYWM0W+hCRMmmPzCcXz+PdrSqBLdcFeN/DpTUd/WGr9o5nAbnXJerlw56devX+RTiGvHOvlH1ycxfiws5UH8xJjg41E+ugGCrqmtrm9fkbvvFmMVNXM4LPf/GSyJ/urM6obGbvY6hFjX6J8/f745yEG1KlOcJRu3VlYnvqi9SEXa6IHLqVOnpHnz5qJ7PFauXLlYeKitSpdwteqIn5gUfDzKR+dF6Oa+t98uUqeOyEMPiTHkSvYb6/LcqH4WXdDaN4xFhw9v2rQpfJqZ33YfGXJzc0OrV68O7d+/3+6lKcU3FpYyzSRqnrHiiJ+YEnw8ykffI1WsGAo9+2wodMklodCRI6HQ738fyqtcOdTXsNdy/ydut4CEWhtVEvmTZQyfkUsvvTTSq1iOZ8+eLVdccUWhESzxMid+PDJn/OHjUT5qJtEtx/74R5EXXhDjn03kkktkzu9/L0feeYf7P3GzBSbUko3bCzSC9thHfRPfdRnLZ+JEEWM9IMNIK/LwwyIvvyxiCPZEY7YfZsLE90SgQv3waHHixAnj6bFiaN++fZaKS/zEmODjUT46gsQYDWF+jIltoVOnzILSXonbK4Ch9k0l6fhVmzt3rlx22WVSrVo1S9kTPzEm+HiUz5YtIv37izz/fNRONrRX4vYKYqjlUSXphKOrD/bXG9qiI35iUPDxKB/dnV0/BRztVQAIp/pc5m130piEoGaSWLM2Y5Wc+LGo/OQHn59YxDqCTywqP/l5jc9PJQvUkbUJOOn8gdMJPzqKpYbFTVCJn7i14AOfSAJ+vx8i6xKkY8+PKsnY0QNx7jLqGwfMj97wgU9iAgEJ9fIDhjE701wL5ZtvvrFUTOInxgQf+EQS8Pv9EFmXgB1721Ty4YcfStOmTaVmzZqWfkaJnxgTfOATScDv90NkXYJ27GlTCW/TE9+O8IFPJIGg3Q+RdQ/csVcfMXQpWWMDh9COHTssFZH4iTHBBz6RBPx+P0TWJYDH3jWV6AafF154odSuXdvSjynxE2OCD3wiCfj9foisSxCPPWsqYfRA4tsRPvCJJBC0+yGy7oE89uJjhi4da+weH9q2bZul4hE/MSb4wCeSgN/vh8i6BPTYm6YS3Ti4Xr16UrduXUs/psRPjAk+8Ikk4Pf7IbIuQT32pKkkaI991Dfxvx984JOYQABDvfaooTs0G9PbQ8bWQ5aKRvzEmOADn0gCfr8fIusS4GPvmUo+/vhjc8KNbvZpxRE/MSX4wCeSgN/vh8i6BPnYc6YSHosT347wgU8kgaDdD5F1D/Sxlx438vLyQsb09lB2dralYhE/MSb4wCeSgN/vh8i6BPzYW6aSpUuXijEMUBo3bmzpx5T4iTHBBz6RBPx+P0TWJejHnjKVBO2xj/om/veDD3wSEwhwqFceOfQxrlatWqH169dbKhLxE2OCD3wiCfj9foisC8ch75hKli9fLhUqVJAmTZpY+hklfmJM8IFPJAG/3w+RdeFYxDOmEh6LE9+O8IFPJIGg3Q+RdefYIOCVxw5jento7dq1lotD/MSo4AOfSAJ+vx8i68KxR0wlK1askNKlS0uLFi0s/ZgSPzEm+MAnkoDf74fIunB8hoAnTCVBe+yjvon//eADn8QECPWEqaRBgwahVatWWX4CIn5iVPCBTyQBv98PkXXh2CQwOku/+P2CAAQgAAHfEBjjCVOJb3BRUAhAAAIeIIBwe6ARKAIEIAABOwQQbju0iAsBCEDAAwQQbg80AkWAAAQgYIcAwm2HFnEhAAEIeIBASaMMYzxQDooAAQhAAALWCCz8f7uu8PwPZGkpAAAAAElFTkSuQmCC"
    }
   },
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this exercise we are looking for the solution of a confined, unsteady couette flow. This problem is controlled by a parabolic PDE and is illustrated below. We consider the unsteady Couette flow, confined by two walls, which has the following governing equation:\n",
    "\n",
    "$$ \\frac{\\partial u}{\\partial t} = \\frac{\\partial^2 u}{\\partial y^2}$$\n",
    "\n",
    "The system we are looking at is described in the following picture, we are looking a fluid constrained between two rigid walls, one fixed and one moving with a constant velocity $u_0$. We are interested by the distribution of the velocity along the y axis within the time frame [0,0.4].\n",
    "\n",
    "![image.png](attachment:image.png)\n",
    "\n",
    "where the velocity $u$: \n",
    "1. $u = 0$ when $y=0$\n",
    "2. $u = u_0$ when $y=1$\n",
    "3. driven by the PDE when $0<y<1$\n",
    "\n",
    "As initial conditions we have $u=0$ for $0<y<1$"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this exercise, you have to:\n",
    "\n",
    "1. Discretize the PDE using the FTCS scheme (i.e forward in time and centered in space)\n",
    "2. Implement a function which solve the PDE using the FTCS scheme. As input we want the initial and boundary conditions as well as the number of points along the y axis and the number of time step. The diffusion number $D$ will be computed inside this function. This function should return the complete history of the velocity $u$ (i.e. a matrix). In this function you can use `for` loops to fill out the \n",
    "columns of $u$. \n",
    "3. Run the problem with 10 grid points along $y$ and 1001 timesteps $\\Delta t$ and plot the distribution of the velocity along $y$ for the timesteps 100, 500 and 1000\n",
    "4. Run the problem with 10 grid points along $y$ and 11 timesteps $\\Delta t$ and plot the distribution of the velocity along $y$ for the timesteps 1, 5 and 10. What can you observe?\n",
    "5. Define a new function where you solve the PDE problem using array slicing instead of for loops and repeat question 3."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "base",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "version": "3.9.10"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "072f120aa02d08ea61d001fa93115a43c81ff98089754e42d3890c56566da5f7"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}