MAIN4_TP2.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy\n",
    "import os\n",
    "\n",
    "import scipy.sparse as sparse\n",
    "import scipy.sparse.linalg\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "#import matplotlib.animation as animation\n",
    "from matplotlib import animation, rc\n",
    "from IPython.display import HTML"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 197,
   "metadata": {},
   "outputs": [],
   "source": [
    "Nx = 200\n",
    "Nt = 2*Nx\n",
    "\n",
    "# domaine spatial et maillage\n",
    "X_min = 0\n",
    "X_max = 1\n",
    "h = 1/Nx\n",
    "X = numpy.zeros(Nx)\n",
    "for i in range(0,Nx):\n",
    "    X[i] = X_min + i*h*(X_max-X_min)\n",
    "\n",
    "# domaine temporel\n",
    "Temps_final = 1\n",
    "dt = Temps_final/Nt\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 198,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = 0.5\n",
    "eps = 0.5\n",
    "\n",
    "# fonction initiale\n",
    "def u_ini(x):\n",
    "    return 1+ eps*numpy.sin(2*numpy.pi*x)\n",
    "\n",
    "# champ vitesse\n",
    "def a_func(x):\n",
    "    return 1+A*numpy.cos(2*numpy.pi*x)\n",
    "\n",
    "# solution approchee \n",
    "u_tx = numpy.zeros((Nt+1, Nx))\n",
    "\n",
    "# tableau des masses\n",
    "mass = numpy.zeros((Nt+1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 199,
   "metadata": {},
   "outputs": [],
   "source": [
    "# sol exacte  --  lorsqu'elle est connue !\n",
    "def u_ex(t,x):\n",
    "    return 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 200,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_times_array_for(array):\n",
    "    Nta = len(array)-1\n",
    "    dta = Temps_final/Nta\n",
    "    times = numpy.zeros(Nta+1)\n",
    "    for n in range(0,Nta+1):\n",
    "        times[n] = n*dta\n",
    "    return times"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_sol(nt, plot_u_ex=False, Y_min=0, Y_max=2, save_fn=None):\n",
    "    if False: #plot_u_ex:\n",
    "        plt.plot(X, u_ex(nt*dt,X), '-', color='r')\n",
    "        \n",
    "    plt.xlim(X_min, X_max)\n",
    "    plt.ylim(Y_min, Y_max)\n",
    "    plt.xlabel('x')\n",
    "    plt.plot(X, u_tx[nt,:], '-', color='k')\n",
    "    if save_fn is not None:\n",
    "        plt.savefig(save_fn)\n",
    "    plt.show()\n",
    "\n",
    "def plot_x_curves(\n",
    "    array1=None, label1=None, \n",
    "    array2=None, label2=None, \n",
    "    array3=None, label3=None,\n",
    "    name=\"some values along x\", save_fn=None, logscale=False\n",
    "):\n",
    "    \"\"\"\n",
    "    to do\n",
    "    \"\"\"\n",
    "    \n",
    "    \n",
    "def plot_t_curves(\n",
    "    array1=None, label1=None, \n",
    "    array2=None, label2=None, \n",
    "    array3=None, label3=None,\n",
    "    name=\"some values along time\", save_fn=None, logscale=False\n",
    "):\n",
    "    plt.xlabel('t')\n",
    "    plt.ylabel(name)\n",
    "    if logscale:\n",
    "        plt.yscale(\"log\")\n",
    "    assert array1 is not None\n",
    "    times = get_times_array_for(array1)\n",
    "    plt.plot(times, array1, '-', label=label1)\n",
    "    if array2 is not None:\n",
    "        times = get_times_array_for(array2)\n",
    "        plt.plot(times, array2, '-', label=label2)\n",
    "    if array3 is not None:\n",
    "        times = get_times_array_for(array3)        \n",
    "        plt.plot(times, array3, '-', label=label3)\n",
    "    if save_fn is not None:\n",
    "        plt.savefig(save_fn)\n",
    "    if label1 is not None:\n",
    "        plt.legend()\n",
    "    plt.show()\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 202,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEKCAYAAAAfGVI8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3Xd4lYX9/vH3JwkkaFFZVn8ISi1WlqzIqAhJEEoQga+IBUcRB1VAu76odRREVBDrYBS0QhEHiiCCBcoIIjIlCCirFikiqwoi5QthJHx+f+SgKQVyICd5kjz367pykfOsc5/nCrnzbHN3REQkvOKCDiAiIsFSEYiIhJyKQEQk5FQEIiIhpyIQEQk5FYGISMjlWwRmVs3M3jezdWa21sx+dYJpzMyGmdlGM/vEzBrlGdfDzP4R+eoR6w8gIiIFY/ldR2BmFwIXuvvHZlYeWAF0dvd1eaZpD9wLtAeaAi+4e1MzqwhkAsmAR+Zt7O57CuXTiIjIact3i8Ddd7j7x5Hv9wHrgarHTdYJGO+5lgLnRQrkZ8Acd/8m8st/DtAupp9AREQKJOF0JjazS4CGwLLjRlUFvszzemtk2MmGn2jZvYBeAGeffXbjyy+//HSiiYiE2ooVK3a5e5UzmTfqIjCzHwCTgV+7+7/P5M1Oxd1fAl4CSE5O9szMzFi/hYhIqWVmX5zpvFGdNWRmZcgtgdfd/Z0TTLINqJbn9UWRYScbLiIixUQ0Zw0ZMAZY7+7PnmSyacAvImcPNQP2uvsOYBbQ1swqmFkFoG1kmIiIFBPR7Bq6CrgV+NTMVkWGPQRUB3D30cAMcs8Y2ggcAHpGxn1jZo8DyyPzDXT3b2IXX0RECirfInD3hYDlM40DfU4ybiww9ozSiYhIodOVxSIiIaciEBEJORWBiEjIqQhEREJORSAiEnIqAhGRkFMRiIiEnIpARCTkVAQiIiGnIhARCTkVgYhIyKkIRERCTkUgIhJyKgIRkZBTEYiIhJyKQEQk5FQEIiIhpyIQEQm5fB9VaWZjgQ7AV+5e9wTj+wE351leLaBK5HnFm4F9QA6Q7e7JsQouIiKxEc0WwTig3clGuvtQd2/g7g2A3wMfHPeA+tTIeJWAiEgxlG8RuPsC4Jv8povoDkwoUCIRESlSMTtGYGZnkbvlMDnPYAdmm9kKM+sVq/cSEZHYyfcYwWm4Dlh03G6hFu6+zczOB+aY2YbIFsZ/iRRFL4Dq1avHMJaIiJxKLM8a6sZxu4XcfVvk36+AKUCTk83s7i+5e7K7J1epUiWGsURE5FRiUgRmdi7QCpiaZ9jZZlb+2PdAW2BNLN5PRERiJ5rTRycAKUBlM9sK9AfKALj76Mhk/wPMdvf9eWb9ITDFzI69zxvu/rfYRRcRkVjItwjcvXsU04wj9zTTvMM2AfXPNJiIiBQNXVksIhJyKgIRkZBTEYiIhJyKQEQk5FQEIiIhpyIQEQk5FYGISMipCEREQk5FICIScioCEZGQUxGIiIScikBEJORUBCIiIaciEBEJORWBiEjIqQhEREJORSAiEnIqAhGRkFMRiIiEXL5FYGZjzewrM1tzkvEpZrbXzFZFvv6QZ1w7M/u7mW00swdjGVxERGIjmi2CcUC7fKb50N0bRL4GAphZPDASSAdqA93NrHZBwsrpycnJ4fDhw7h70FFEpBhLyG8Cd19gZpecwbKbABvdfROAmb0JdALWncGyJA9354svvmDDhg2sX7+ezz77jG3btrFjxw527tzJ3r17OXjwIEeOHAHAzEhKSqJcuXJUqVKFCy64gAsvvJBLLrmEyy+/nFq1alGrVi3Kly8f8CcTkSDkWwRRam5mq4HtwP+6+1qgKvBlnmm2Ak1PtgAz6wX0AqhevXqMYpUOR44cYcmSJSxYsIClS5eydOlSdu/e/d34ihUrUq1aNS644ALq1q1LhQoVSEpKIikpifj4eA4dOsTBgwc5cOAAX331FTt27GD58uVMmjSJ7OxsILcsateuTfPmzWnevDmtW7fm4osvDuoji0gRikURfAxc7O7/Z2btgXeBmqe7EHd/CXgJIDk5OfT7Mr7++mv++te/MmPGDGbPns2///1vAGrVqkWnTp1o0qQJtWvXplatWlSuXPmM3uPIkSNs2rSJDRs2sHr1apYuXcqkSZN4+eWXv3uv9PR0OnbsyNVXX01cnM4tECmNLJr9x5FdQ39197pRTLsZSCa3DAa4+88iw38P4O5P5beM5ORkz8zMzDdXabN//36mTp3K66+/zqxZs8jJyeH//b//R/v27UlPTyc1NZUKFSoUaoajR4+yfv16Zs2axcyZM1mwYAGHDx+matWqdOvWje7du9OoUSPMrFBziMjpMbMV7p58RjO7e75fwCXAmpOMu4DvC6UJsAUwcrc2NgE1gLLAaqBONO/XuHFjD5PVq1d7r169/Oyzz3bAq1Wr5g8++KCvXLnSjx49Gmi2ffv2+YQJE/y6667zMmXKOOBXXHGFjx492vft2xdoNhH5HpDpUfx+PdFXNCUwAdgBHCF3P/8dwN3A3ZHxfYG1kV/0S4Gf5pm3PfAZ8DnwcLShwlAE2dnZ/vbbb3vLli0d8KSkJO/Zs6d/8MEHnpOTE3S8E9q9e7ePHj3a69ev74CXL1/e+/bt659//nnQ0URCr1CLIIiv0lwER44c8XHjxvlll13mgF9yySU+dOhQ3717d9DRonb06FFfvHix33zzzV6mTBmPj4/3W2+91detWxd0NJHQKkgR6OhfEcnOzmbMmDFcdtll3HbbbZQrV46JEyeyceNG/vd//5eKFSsGHTFqZkbz5s157bXX+Oc//8l9993H5MmTqVOnDl27dmXDhg1BRxSR06AiKGTuzowZM2jQoAF33nknlSpVYurUqaxcuZKuXbsSHx8fdMQCqVq1Ks8++yybN2/moYceYtasWdStW5e7776bHTt2BB1PRKKgIihEK1eupE2bNlx77bUcOnSISZMm8dFHH9GxY8dSd9ZNlSpVGDRoEBs3bqR3796MGTOGH//4xwwYMICsrKyg44nIKagICsHevXvp27cvjRs3ZtWqVbzwwgusXbuWLl26lLoCON7555/PsGHD2LBhAx06dOCxxx6jTp06TJ8+PehoInISKoIYcnfefPNNLr/8ckaNGkXfvn3ZuHEj9913H2XLlg06XpG69NJLeeutt5g3bx5JSUl06NCBzp0788UXXwQdTUSOoyKIkc2bN5Oenk737t2pWrUqy5YtY9iwYZx33nlBRwtUamoqq1atYvDgwcyZM4fatWszYsQIjh49GnQ0EYlQERSQu/Pyyy9Tr149Fi9ezLBhw1i2bBnJyWd2gV9pVLZsWR544AHWr19Py5Ytuffee2ndujWbNm0KOpqIoCIokO3bt3Pttddy1113ceWVV/Lpp59y7733lvgzgQpL9erVmTFjBmPGjOHjjz/miiuuYOTIkbpNtkjAVARn6O2336Zu3brMnz+f4cOHM3fuXN2tMwpmxu23386aNWto0aIFffv2pUOHDnz99ddBRxMJLRXBacrKyuKee+7hxhtv5LLLLmP16tX07dtXd+Y8TdWqVWPmzJmMGDGCjIwM6tevz7x584KOJRJK+u11GjZs2ECzZs0YPXo0999/Px9++CE1a572Hbclwszo06cPy5Yt49xzz+Waa67hoYce+u6BOiJSNFQEURo/fjyNGzdm+/btzJgxgyFDhlCmTJmgY5UK9evXJzMzkzvuuIOnnnqKlJQUtm/fHnQskdBQEeTj8OHD9O7dmx49epCcnMyqVatIT08POlapc/bZZ/PnP/+ZCRMmsHr1aho1asSHH34YdCyRUFARnMLOnTtp3bo1o0aNol+/fmRkZFC1atWgY5Vq3bp1Y9myZZxzzjmkpaUxfPhwnVUkUshUBCfx0UcfkZyczIoVK5gwYQJPP/00CQmxesSznEqdOnVYvnw57du357777uPWW2/lwIEDQccSKbVUBCfwl7/8hauvvpoyZcqwZMkSunXrFnSk0Dn33HOZMmUKjz/+OG+88QatWrXScQORQqIiyOPo0aP069eP22+/nZYtW5KZmUn9+vWDjhVacXFxPPLII0ydOpX169fTpEkTVq5cGXQskVJHRRBx4MABunbtyjPPPEPv3r2ZOXMmlSpVCjqWANdddx2LFi0iLi6OFi1a8O677wYdSaRUybcIzGysmX1lZmtOMv5mM/vEzD41s8VmVj/PuM2R4avMLDOWwWNp586dpKSkMGXKFJ577jlGjBih4wHFTP369Vm2bBl16tTh+uuvZ+jQoTqILBIj0WwRjAPanWL8P4FW7l4PeBx46bjxqe7ewN2L5V3Y1q5dS7NmzVi7di1Tpkzh17/+dal/ZkBJdeGFF/LBBx/QtWtX7r//fu6++26ys7ODjiVS4uX7Z6+7LzCzS04xfnGel0uBiwoeq2hkZGRw/fXXc9ZZZ7FgwQIaN24cdCTJR7ly5ZgwYQKXXnopTz31FDt37mTChAmcddZZQUcTKbFifYzgDmBmntcOzDazFWbW61QzmlkvM8s0s8yiuAHZxIkTSU9Pp3r16ixbtkwlUILExcXx5JNPMnz4cN577z2uueYadu/eHXQskRIrZkVgZqnkFsEDeQa3cPdGQDrQx8xanmx+d3/J3ZPdPblKlSqxinVCo0aNolu3bjRt2pQPP/yQ6tWrF+r7SeHo27cvEydO5OOPP6ZFixZ6+pnIGYpJEZjZFcDLQCd3/+5PM3ffFvn3K2AK0CQW73em3J3HHnuM3r1706FDB2bPnh36J4iVdDfccAOzZ89mx44dNG/enNWrVwcdSaTEKXARmFl14B3gVnf/LM/ws82s/LHvgbbACc88Kgo5OTn07duXAQMGcNttt/HOO+9Qrly5oOJIDLVs2ZKFCxcSFxdHq1atWLRoUdCRREqUaE4fnQAsAX5iZlvN7A4zu9vM7o5M8gegEvCn404T/SGw0MxWAx8B0939b4XwGfJ16NAhbrrpJv70pz/Rr18/xo4dq9NDS5m6deuyePFizj//fNq2bcucOXOCjiRSYlhxPBc7OTnZMzNjc9nB/v37+Z//+R/mzJnD008/Tb9+/WKyXCme/vWvf9G2bVs2bNjAW2+9RefOnYOOJFIkzGzFmZ6mX6qvLN63bx/p6elkZGQwduxYlUAI/PCHP2T+/Pk0bNiQG264gddeey3oSCLFXqktgm+//ZY2bdqwePFi3njjDXr27Bl0JCkiFSpUYM6cObRq1Ypf/OIXjB49OuhIIsVaqSyCXbt2kZaWxscff8ykSZP4+c9/HnQkKWLly5dn+vTpdOjQgXvuuYchQ4YEHUmk2Cp1RbBz505SU1NZv349U6dO1T7iEEtKSmLy5Ml0796dBx98kMceeyzoSCLFUqk6dWbr1q20bt2arVu3Mn36dNLS0oKOJAErU6YMr776KomJiQwYMICcnBwee+wx3U9KJI9SUwSbN28mLS2NXbt2MWvWLFq0aBF0JCkm4uPjGTNmDAkJCTz++ONkZ2fzxBNPqAxEIkpFEWzcuJG0tDT27dvH3LlzadIk0AuYpRiKi4vjxRdfJD4+nqeeeoqcnBwGDx6sMhChFBTB559/TkpKCgcPHuT999+nQYMGQUeSYiouLo5Ro0aRkJDA008/TXZ2Ns8884zKQEKvRBfBpk2bSE1N5eDBg2RkZOixkpIvM2P48OHEx8fz7LPPkp2dzfPPP68ykFArsUXwz3/+k9TUVPbv368SkNNiZjz//PPEx8fz3HPPkZ2dzfDhw4mLK3Un0YlEpUQWwebNm0lNTWXfvn1kZGRod5CcNjPjj3/8IwkJCQwdOpSjR4/ypz/9SVsGEkolrgi++OILUlNT2bt3L3PnzqVhw4ZBR5ISyswYMmQI8fHxDB48mPj4eIYPH64ykNApUUXw5Zdfkpqayp49e5g7d66eKiYFZmY8+eST3x04TkhI4LnnnlMZSKiUmCLYunUrKSkp7N69m7lz55KcfEY32RP5L2b23VlEx44d6GwiCZMSUQTbtm0jJSWFXbt2MXv2bK688sqgI0kpY2Y8++yz5OTk8OyzzxIfH8+QIUNUBhIKxb4Itm3bRmpqKl999RWzZ8+madOmQUeSUsrMeOGFF8jOzmbo0KEkJCToCmQJhWJdBDt27CAtLY0dO3Ywe/ZsmjVrFnQkKeXMjBEjRpCTk8NTTz1FmTJldLM6KfWKbREcu4votm3bmDVrFs2bNw86koTEsSuQs7OzGThwIPHx8fzhD38IOpZIoYnqChozG2tmX5nZCR8+b7mGmdlGM/vEzBrlGdfDzP4R+eoRzfsdOXKE1NRUtm7dysyZM7nqqqui+zQiMRIXF8ef//xnevToQf/+/XnyySeDjiRSaKLdIhgHjADGn2R8OlAz8tUUGAU0NbOKQH8gGXBghZlNc/c9p3qzzz77jLi4OGbOnMnVV18dZUSR2IqLi2PMmDHk5OTw8MMPk5CQwP333x90LJGYi6oI3H2BmV1yikk6AePd3YGlZnaemV0IpABz3P0bADObA7QDJpzq/Q4fPsy8efNo2bJlNPFECk18fDzjxo0jJyeHBx54gISEBH77298GHUvkP3zyyScFmj9WxwiqAl/meb01Muxkw/+LmfUCekHuA8hbtWoVo2giBRMfH8/48ePJzs7md7/7HfHx8fzqV78KOpYIAJ9++mmBH8JVbA4Wu/tLwEsAycnJHnAckf+QkJDA66+/Tk5ODr/+9a9JSEigT58+QceSkFu7di2tW7cmMTGxQMuJ1e0WtwHV8ry+KDLsZMNFSpwyZcowYcIEOnXqRN++fRk1alTQkSTE1q9fT1paGgkJCbz//vsFWlasimAa8IvI2UPNgL3uvgOYBbQ1swpmVgFoGxkmUiKVLVuWiRMnct1119G7d29efPHFoCNJCP39738nLS0NM2PevHlcdtllBVpeVLuGzGwCuQd+K5vZVnLPBCoD4O6jgRlAe2AjcADoGRn3jZk9DiyPLGrgsQPHIiVV2bJlefvtt+nSpQt33303ZkavXr2CjiUh8Y9//IPU1FSOHj3K/Pnzufzyywu8zGjPGuqez3gHTrjD1N3HAmNPP5pI8ZWYmMjkyZPp0qULv/zlLzEz7rrrrqBjSSm3ceNGUlNTOXLkCPPnz6dWrVoxWa4eySRyho6VQfv27enVqxdjxowJOpKUYnkfzTtv3jzq1KkTs2WrCEQK4FgZtGvXjrvuuouxY7XxK7F37KmMBw4cICMjg3r16sV0+SoCkQJKSkpiypQp/OxnP+POO+/kL3/5S9CRpBTZsmULqamp/Pvf/2bOnDmF8nx2FYFIDBwrgzZt2nDHHXcwbty4oCNJKZD3qYxz5syhUaNG+c90BlQEIjGSlJTEu+++yzXXXMPtt9/OK6+8EnQkKcG2bNnyHw/kKsynMqoIRGKoXLlyTJ06lWuuuYaePXvy6quvBh1JSqDNmzfTqlUrdu/ezZw5c2jSpEmhvp+KQCTGjpVB69at6dGjB6+99lrQkaQE2bRpE61atWLv3r1kZGQUegmAikCkUBwrg9TUVHr06MHrr78edCQpATZu3EirVq34v//7PzIyMmjcuHGRvK+KQKSQnHXWWbz33nukpKRw66236gCynNJnn31Gq1atyMrKYt68eTRs2LDI3ltFIFKIjpVBmzZt6Nmzp25UJye0YcMGUlJSOHLkCO+//36hnCJ6KioCkUJ21llnMW3atO9uVPfcc88FHUmKkXXr1pGSkvLdvYNifbFYNFQEIkUgMTGRSZMm0bVrV37729/yxBNPBB1JioFVq1aRkpJCXFwc8+fPp3bt2oHkKDYPphEp7cqWLcsbb7xBUlISjzzyCAcOHGDQoEGYWdDRJACLFy+mffv2nHPOOcydO7fAt5IuCBWBSBFKSEhg3LhxJCUl8eSTT5KVlcUf//hHlUHIzJ07l06dOlG1alXmzp1L9erVA82jIhApYnFxcbz44ouUK1eO5557jqysLEaOHElcnPbUhsHUqVO58cYb+clPfsKcOXP44Q9/GHQkFYFIEMyM559/nnLlyjFkyBC+/fZbXnnlFcqWLRt0NClEr732GrfddhvJycnMmDGDihUrBh0JUBGIBMbMGDx4MJUqVeL+++/nm2++YfLkyfzgBz8IOpoUglGjRtGnTx9SUlKYOnUq5cuXDzrSd7QtKhKwfv36MXbsWObOnUvr1q3ZtWtX0JEkhtydQYMG0bt3b6699lpmzJhRrEoAoiwCM2tnZn83s41m9uAJxj9nZqsiX5+Z2bd5xuXkGTctluFFSouePXsyZcoUPvnkE1q0aMGWLVuCjiQxkJOTQ+/evXn00Ue59dZbeeedd0hKSgo61n/JtwjMLB4YCaQDtYHuZvYfJ7u6+2/cvYG7NwCGA+/kGZ11bJy7d4xhdpFSpWPHjsyaNYudO3dy1VVXsW7duqAjSQFkZWVxww03MHr0aB544AFeeeUVypQpE3SsE4pmi6AJsNHdN7n7YeBNoNMppu8OTIhFOJGwadmyJQsWLCA7O5sWLVqwaNGioCPJGfjmm29o06YNU6dOZdiwYQwePLhYnyIcTRFUBb7M83prZNh/MbOLgRrAvDyDk8ws08yWmlnnk72JmfWKTJf59ddfRxFLpHS64oorWLx4MZUrVyYtLY0JE/R3VUmyZcsWWrRowfLly3nrrbe49957g46Ur1gfLO4GTHL3nDzDLnb3ZOAm4Hkzu/REM7r7S+6e7O7JVapUiXEskZKlRo0aLFmyhGbNmnHTTTfxxBNP4O5Bx5J8LF++nGbNmrF9+3ZmzZpF165dg44UlWiKYBtQLc/riyLDTqQbx+0WcvdtkX83AfOBoru3qkgJVqlSJWbPns0tt9zCI488wu23387hw4eDjiUn8fbbb9OyZUsSExNZtGgRKSkpQUeKWjRFsByoaWY1zKwsub/s/+vsHzO7HKgALMkzrIKZJUa+rwxcBegImEiUEhMTGT9+PAMGDGDcuHG0a9eOPXv2BB1L8jh2euiNN95Io0aN+Oijj6hTp07QsU5LvkXg7tlAX2AWsB6Y6O5rzWygmeU9C6gb8Kb/5/ZrLSDTzFYD7wOD3V1FIHIazIz+/fvz6quvsnDhQpo1a8b69euDjiXAwYMHufXWW3n00Ue55ZZbyMjIoCTu2rbiuN8xOTnZMzMzg44hUuwsXLiQLl26kJWVxauvvkqnTqc6gU8K07Zt27jhhhtYunQpgwYN4qGHHgr0zCAzWxE5HnvadGWxSAnSokULMjMzufzyy+ncuTP9+/fn6NGjQccKnQ8++IBGjRqxZs0aJk2axMMPP1ysTw/Nj4pApISpVq0aCxYs4LbbbmPgwIF07tyZvXv3Bh0rFNyd559/ntatW3PeeeexbNkyunTpEnSsAlMRiJRASUlJjB07lhEjRjBz5kyuvPJKVq5cGXSsUm3//v3ccsst/OY3v6FDhw4sX748sCeKxZqKQKSEMjP69OnDvHnzOHDgAM2aNWPkyJG63qAQrFq1isaNGzNhwgSeeOIJ3nnnHc4555ygY8WMikCkhLv66qtZtWoV11xzDX379qVLly46xTRG3J1hw4bRtGlT9u3bR0ZGBg899FCpe4hQ6fo0IiFVuXJl3nvvPZ555hnee+89GjZsyMKFC4OOVaLt2rWLjh078qtf/Yq2bduyevVqUlNTg45VKFQEIqVEXFwcv/vd71i4cCFxcXG0bNmSfv36cfDgwaCjlTh//etfueKKK5g9ezYvvPAC06ZNo3LlykHHKjQqApFSpmnTpqxevZq77rqLZ555hsaNG7NixYqgY5UIe/bsoUePHlx33XVUrlyZZcuWcd9995XoU0OjoSIQKYXKly/Piy++yMyZM/n2229p2rQpjzzyCFlZWUFHK7amT59O3bp1ef3113n00UfJzMykQYMGQccqEioCkVKsXbt2rFmzhptvvpknnniCunXrMmvWrKBjFStffvklN9xwAx06dKBixYosW7aMgQMHUrZs2aCjFRkVgUgpV6FCBV555RUyMjJISEigXbt2dOvWje3btwcdLVBHjhxh6NCh1KpVi+nTpzNo0CAyMzNp3Lhx0NGKnIpAJCTS0tL45JNPGDhwIO+++y6XXXYZAwcOZP/+/UFHK1Luzt/+9jcaNGjA/fffT1paGuvWrePhhx8mMTEx6HiBUBGIhEhiYiKPPvooa9asoV27dvTv35+aNWsyZswYcnJy8l9ACZeZmUnr1q1JT0/n4MGDTJs2jWnTplGjRo2gowVKRSASQj/+8Y+ZNGkSixYt4uKLL+bOO++kfv36vPnmm6WyENasWcPPf/5zrrzySj799FOGDRvG+vXrue6664KOViyoCERC7Kc//SmLFy9m4sSJHD16lO7du1OnTh3Gjx9PdnZ20PEKbPny5XTu3Jl69eoxY8YMHnnkET7//HPuvffeUB0Mzo+KQCTkzIyuXbuyZs0aJk6cSGJiIj169ODSSy9l8ODB7Nq1K+iIpyU7O5spU6bQunVrmjRpwoIFC+jfvz9ffPEFjz/+eKm6R1CsqAhEBMi9Mrlr166sXLmSqVOncumll/L73/+eiy66iJ49e7Jo0aJifUO7bdu2MWjQIGrUqMH111/PP/7xD4YMGcIXX3zBgAEDqFixYtARiy09oUxETmrt2rWMHDmS8ePHs3//fmrUqMFNN93EzTffTK1atYKOx+7du5k8eTJvvPEGCxYswN1p06YNffr04dprryUhISHoiEWmIE8oi6oIzKwd8AIQD7zs7oOPG38bMBTYFhk0wt1fjozrATwSGT7I3V/J7/1UBCLFy759+3j33Xd57bXXmDt3LkePHqVmzZqkp6fTvn17WrVqRVJSUqHncHfWrFnDzJkzmTlzJgsXLiQ7O5uf/OQndO/enZtuuomaNWsWeo7iqFCLwMzigc+ANsBWYDnQPe9D6CNFkOzufY+btyKQCSQDDqwAGrv7Ke+RqyIQKb527tzJ5MmTmTFjBvPmzePgwYMkJibSuHFjmjdvTvPmzalXrx4/+tGPCvQXubvzr3/9i3Xr1vHRRx+xZMkSlixZwtdffw1AvXr1aN++PTfeeCMNGza78M2tAAAJnklEQVQs9fcDyk9hF0FzYIC7/yzy+vcA7v5Unmlu48RF0B1IcfdfRl6/CMx39wmnek8VgUjJkJWVxfz585k7dy5LlixhxYoVHD58GICyZctSs2ZNfvSjH3HBBRdw4YUXUqVKFc466yySkpJISkriyJEjHDx4kKysLPbs2cPOnTvZsWMHX375JRs2bODbb7/97r1q1qxJ8+bNufrqq2nXrh0XXXRRUB+7WCpIEURT11WBL/O83go0PcF0XcysJblbD79x9y9PMm/VE72JmfUCegFUr149ilgiErRy5cqRnp5Oeno6AIcOHWL16tWsW7eO9evXs379erZs2cKyZcv4+uuv8z3YfM4553DBBRdQtWpVunXrRq1atahVqxaNGjWiUqVKRfGRQilWR1LeAya4+yEz+yXwCpB2Ogtw95eAlyB3iyBGuUSkCCUmJtKkSROaNGnyX+OOHDnCnj17yMrK4uDBgxw8eJAyZcqQlJREuXLlOOecczj77LMDSC3RFME2oFqe1xfx/UFhANx9d56XLwNP55k35bh5559uSBEp+cqUKcP5558fdAw5gWiuI1gO1DSzGmZWFugGTMs7gZldmOdlR2B95PtZQFszq2BmFYC2kWEiIlJM5LtF4O7ZZtaX3F/g8cBYd19rZgOBTHefBtxnZh2BbOAb4LbIvN+Y2ePklgnAQHf/phA+h4iInCFdUCYiUgoU5Kwh3WJCRCTkVAQiIiGnIhARCTkVgYhIyKkIRERCTkUgIhJyKgIRkZBTEYiIhJyKQEQk5FQEIiIhpyIQEQk5FYGISMipCEREQk5FICIScioCEZGQUxGIiIScikBEJORUBCIiIRdVEZhZOzP7u5ltNLMHTzD+t2a2zsw+MbMMM7s4z7gcM1sV+Zp2/LwiIhKsfB9eb2bxwEigDbAVWG5m09x9XZ7JVgLJ7n7AzO4BngZ+HhmX5e4NYpxbRERiJJotgibARnff5O6HgTeBTnkncPf33f1A5OVS4KLYxhQRkcISTRFUBb7M83prZNjJ3AHMzPM6ycwyzWypmXU+g4wiIlKI8t01dDrM7BYgGWiVZ/DF7r7NzH4EzDOzT9398xPM2wvoBVC9evVYxhIRkVOIZotgG1Atz+uLIsP+g5ldAzwMdHT3Q8eGu/u2yL+bgPlAwxO9ibu/5O7J7p5cpUqVqD+AiIgUTDRFsByoaWY1zKws0A34j7N/zKwh8CK5JfBVnuEVzCwx8n1l4Cog70FmEREJWL67htw928z6ArOAeGCsu681s4FAprtPA4YCPwDeNjOALe7eEagFvGhmR8ktncHHnW0kIiIBM3cPOsN/SU5O9szMzKBjiIiUGGa2wt2Tz2ReXVksIhJyKgIRkZBTEYiIhJyKQEQk5FQEIiIhpyIQEQk5FYGISMipCEREQk5FICIScioCEZGQUxGIiIScikBEJORUBCIiIaciEBEJORWBiEjIqQhEREJORSAiEnIqAhGRkFMRiIiEXFRFYGbtzOzvZrbRzB48wfhEM3srMn6ZmV2SZ9zvI8P/bmY/i110ERGJhXyLwMzigZFAOlAb6G5mtY+b7A5gj7v/GHgOGBKZtzbQDagDtAP+FFmeiIgUE9FsETQBNrr7Jnc/DLwJdDpumk7AK5HvJwGtzcwiw99090Pu/k9gY2R5IiJSTCREMU1V4Ms8r7cCTU82jbtnm9leoFJk+NLj5q16ojcxs15Ar8jLQ2a2JopsYVAZ2BV0iGJA6+F7Whff07r43k/OdMZoiqBIuPtLwEsAZpbp7skBRyoWtC5yaT18T+vie1oX3zOzzDOdN5pdQ9uAanleXxQZdsJpzCwBOBfYHeW8IiISoGiKYDlQ08xqmFlZcg/+TjtummlAj8j3NwDz3N0jw7tFziqqAdQEPopNdBERiYV8dw1F9vn3BWYB8cBYd19rZgOBTHefBowBXjWzjcA35JYFkekmAuuAbKCPu+dEkeulM/s4pZLWRS6th+9pXXxP6+J7Z7wuLPcPdxERCStdWSwiEnIqAhGRkAusCApy24rSJop18VszW2dmn5hZhpldHETOopDfusgzXRczczMrtacORrMuzOzGyM/GWjN7o6gzFpUo/o9UN7P3zWxl5P9J+yByFgUzG2tmX53sWivLNSyyrj4xs0b5LtTdi/yL3IPOnwM/AsoCq4Hax03TGxgd+b4b8FYQWYvJukgFzop8f0+Y10VkuvLAAnIvVkwOOneAPxc1gZVAhcjr84POHeC6eAm4J/J9bWBz0LkLcX20BBoBa04yvj0wEzCgGbAsv2UGtUVQkNtWlDb5rgt3f9/dD0ReLiX3eozSKJqfC4DHyb2f1cGiDFfEolkXdwEj3X0PgLt/VcQZi0o068KBcyLfnwtsL8J8RcrdF5B7dubJdALGe66lwHlmduGplhlUEZzothXH33riP25bARy7bUVpE826yOsOctu+NMp3XUQ2c6u5+/SiDBaAaH4uLgMuM7NFZrbUzNoVWbqiFc26GADcYmZbgRnAvUUTrVg63d8pxecWE5I/M7sFSAZaBZ0lCGYWBzwL3BZwlOIigdzdQynkbiUuMLN67v5toKmC0R0Y5+5/NLPm5F7XVNfdjwYdrCQIaougILetKG2iug2HmV0DPAx0dPdDRZStqOW3LsoDdYH5ZraZ3P2f00rpAeNofi62AtPc/Yjn3t33M3KLobSJZl3cAUwEcPclQBK5N6QLo9O+tU9QRVCQ21aUNvmuCzNrCLxIbgmU1v3AkM+6cPe97l7Z3S9x90vIPV7S0d3P+GZbxVg0/0feJXdrADOrTO6uok1FGbKIRLMutgCtAcysFrlF8HWRpiw+pgG/iJw91AzY6+47TjVDILuGvAC3rShtolwXQ4EfAG9HjpdvcfeOgYUuJFGui1CIcl3MAtqa2TogB+jn7qVuqznKdfE74M9m9htyDxzfVkr/cMTMJpD7B0DlyDGR/kAZAHcfTe4xkvbkPv/lANAz32WW0nUlIiJR0pXFIiIhpyIQEQk5FYGISMipCEREQk5FICIScioCEZGQUxGIiIScikAkH2Z2ZeS+7klmdnbk3v91g84lEiu6oEwkCmY2iNzbFpQDtrr7UwFHEokZFYFIFCL3uFlO7jMQfuruOQFHEokZ7RoSiU4lcu/3VJ7cLQORUkNbBCJRMLNp5D4ZqwZwobv3DTiSSMzowTQi+TCzXwBH3P0NM4sHFptZmrvPCzqbSCxoi0BEJOR0jEBEJORUBCIiIaciEBEJORWBiEjIqQhEREJORSAiEnIqAhGRkPv/w3R/ErSg4HMAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "mass =  1.0\n"
     ]
    }
   ],
   "source": [
    "## initialisation de la solution avec u_ini\n",
    "u_tx[0,:] = u_ini(X)\n",
    "plot_sol(0, plot_u_ex=True, save_fn=\"u_ini.png\")\n",
    "\n",
    "mass[0] = h*numpy.sum(u_tx[0,:])\n",
    "print(\"mass = \",mass[0])\n",
    "\n",
    "## vitesse\n",
    "a = numpy.zeros((Nx,))\n",
    "a[:] = a_func(X)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 203,
   "metadata": {},
   "outputs": [],
   "source": [
    "# choix du type de schema: implicite ou explicite\n",
    "implicit_scheme = False\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 204,
   "metadata": {},
   "outputs": [],
   "source": [
    "# assemblage d'une matrice creuse M (pour un schema implicite)\n",
    "if implicit_scheme:\n",
    "    row = list()\n",
    "    col = list()\n",
    "    data = list()\n",
    "    for i in range(0,Nx):\n",
    "\n",
    "        # M_i,i = qque chose (exemple)\n",
    "        row.append((i))\n",
    "        col.append((i))\n",
    "        data.append( 1+0.5*dt/h )  \n",
    "    \n",
    "        # M_i,i-1 = autre chose (exemple)\n",
    "        row.append((i))\n",
    "        col.append((numpy.mod(i-1,Nx)))  # modulo pour conditions periodiques\n",
    "        data.append( -0.5*dt/h )\n",
    "\n",
    "    row = numpy.array(row)\n",
    "    col = numpy.array(col)\n",
    "    data = numpy.array(data)      \n",
    "    M = (sparse.coo_matrix((data, (row, col)), shape=(Nx, Nx))).tocsr()\n",
    "else:\n",
    "    M = None"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 205,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD8CAYAAAB0IB+mAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAADMVJREFUeJzt3W2M5fVZh/HrWxZqTLGY7hjJ7tLBuDViNSlOkMZEMagBYtgXRQNJpRB0Eyw+NiZEk2LwVWOsCSkprrIBGgWUmmaM25DGYlDjEgZKkYdgRqxlkYQp1G0bbHH19sU5cSfTWc5h5jwMe1+fZMJ5+M353/ll9poz539mSFUhSTr9vW3eA0iSZsPgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqYte8Drx79+5aXFyc1+El6S3pscce+0pVLWzlc+cW/MXFRVZWVuZ1eEl6S0ry71v9XF/SkaQmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJamJk8JPsS/JQkmeSPJ3k1zdZkyS3JVlN8mSSC6czriRpq8b588gngI9U1eNJzgYeS/K5qnpm3ZrLgf3Djx8DPjn8ryRphxj5DL+qXqqqx4eXvw48C+zZsOwAcE8NHAXOSXLuxKeVJG3Zm3oNP8ki8D7gkQ137QFeWHf9GN/+TYEkB5OsJFlZW1t7c5NKkrZl7OAneQfwaeA3quprWzlYVR2qqqWqWlpY2NL/oUuStEVjBT/JmQxi/2dV9VebLHkR2Lfu+t7hbZKkHWKcd+kEuBN4tqo+foply8C1w3frXAwcr6qXJjinJGmbxnmXzo8Dvwj8c5Inhrf9DnAeQFXdARwBrgBWgdeA6yc/qiRpO0YGv6r+AciINQV8eFJDSZImz9+0laQmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNjAx+ksNJXk7y1CnuvyTJ8SRPDD8+OvkxJUnbtWuMNXcBnwDueYM1f19VPzeRiSRJUzHyGX5VPQy8OoNZJElTNKnX8N+f5ItJPpvkhyb0mJKkCRrnJZ1RHgfeXVXfSHIF8Blg/2YLkxwEDgKcd955Ezi0JGlc236GX1Vfq6pvDC8fAc5MsvsUaw9V1VJVLS0sLGz30JKkN2HbwU/yvUkyvHzR8DFf2e7jSpIma+RLOknuBS4Bdic5BtwCnAlQVXcAVwE3JjkB/BdwdVXV1CaWJG3JyOBX1TUj7v8Eg7dtSpJ2MH/TVpKaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1MTL4SQ4neTnJU6e4P0luS7Ka5MkkF05+TEnSdo3zDP8u4LI3uP9yYP/w4yDwye2PJUmatJHBr6qHgVffYMkB4J4aOAqck+TcSQ0oSZqMSbyGvwd4Yd31Y8PbJEk7yExP2iY5mGQlycra2tosDy1J7U0i+C8C+9Zd3zu87dtU1aGqWqqqpYWFhQkcWpI0rkkEfxm4dvhunYuB41X10gQeV5I0QbtGLUhyL3AJsDvJMeAW4EyAqroDOAJcAawCrwHXT2tYSdLWjQx+VV0z4v4CPjyxiSRJU+Fv2kpSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDUxVvCTXJbkuSSrSW7e5P7rkqwleWL48UuTH1WStB27Ri1IcgZwO/AzwDHg0STLVfXMhqX3V9VNU5hRkjQB4zzDvwhYrarnq+p14D7gwHTHkiRN2jjB3wO8sO76seFtG30gyZNJHkiyb7MHSnIwyUqSlbW1tS2MK0naqkmdtP1rYLGqfgT4HHD3Zouq6lBVLVXV0sLCwoQOLUkaxzjBfxFY/4x97/C2/1dVr1TVt4ZX/xT40cmMJ0malHGC/yiwP8n5Sc4CrgaW1y9Icu66q1cCz05uREnSJIx8l05VnUhyE/AgcAZwuKqeTnIrsFJVy8CvJbkSOAG8Clw3xZklSVuQqprLgZeWlmplZWUux5akt6okj1XV0lY+19+0laQmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNGHxJasLgS1ITBl+SmjD4ktSEwZekJgy+JDVh8CWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqQmDL0lNjBX8JJcleS7JapKbN7n/7UnuH97/SJLFSQ8qSdqekcFPcgZwO3A5cAFwTZILNiy7AfhqVX0/8EfAxyY9qCRpe8Z5hn8RsFpVz1fV68B9wIENaw4Adw8vPwBcmiSTG1OStF3jBH8P8MK668eGt226pqpOAMeBd01iQEnSZMz0pG2Sg0lWkqysra3N8tCS1N44wX8R2Lfu+t7hbZuuSbILeCfwysYHqqpDVbVUVUsLCwtbm1iStCXjBP9RYH+S85OcBVwNLG9Yswx8aHj5KuDzVVWTG1OStF27Ri2oqhNJbgIeBM4ADlfV00luBVaqahm4E/hUklXgVQbfFCRJO8jI4ANU1RHgyIbbPrru8jeBn5/saJKkSfI3bSWpCYMvSU0YfElqwuBLUhMGX5KaMPiS1ITBl6QmDL4kNWHwJakJgy9JTRh8SWrC4EtSEwZfkprIvP5sfZKvA8/N5eA7z27gK/MeYodwL05yL05yL076gao6eyufONafR56S56pqaY7H3zGSrLgXA+7FSe7FSe7FSUlWtvq5vqQjSU0YfElqYp7BPzTHY+807sVJ7sVJ7sVJ7sVJW96LuZ20lSTNli/pSFITUw9+ksuSPJdkNcnNm9z/9iT3D+9/JMnitGealzH24reSPJPkySR/m+Td85hzFkbtxbp1H0hSSU7bd2iMsxdJfmH4tfF0kj+f9YyzMsa/kfOSPJTkC8N/J1fMY85pS3I4yctJnjrF/Uly23Cfnkxy4VgPXFVT+wDOAP4V+D7gLOCLwAUb1vwKcMfw8tXA/dOcaV4fY+7FTwHfObx8Y+e9GK47G3gYOAoszXvuOX5d7Ae+AHz38Pr3zHvuOe7FIeDG4eULgC/Ne+4p7cVPABcCT53i/iuAzwIBLgYeGedxp/0M/yJgtaqer6rXgfuAAxvWHADuHl5+ALg0SaY81zyM3IuqeqiqXhtePQrsnfGMszLO1wXA7wMfA745y+FmbJy9+GXg9qr6KkBVvTzjGWdlnL0o4LuGl98J/McM55uZqnoYePUNlhwA7qmBo8A5Sc4d9bjTDv4e4IV1148Nb9t0TVWdAI4D75ryXPMwzl6sdwOD7+Cno5F7MfwRdV9V/c0sB5uDcb4u3gO8J8k/Jjma5LKZTTdb4+zF7wEfTHIMOAL86mxG23HebE+A+f6mrU4hyQeBJeAn5z3LPCR5G/Bx4Lo5j7JT7GLwss4lDH7qezjJD1fVf851qvm4Brirqv4wyfuBTyV5b1X977wHeyuY9jP8F4F9667vHd626Zokuxj8mPbKlOeah3H2giQ/DfwucGVVfWtGs83aqL04G3gv8HdJvsTgNcrl0/TE7ThfF8eA5ar676r6N+BfGHwDON2Msxc3AH8BUFX/BHwHg7+z081YPdlo2sF/FNif5PwkZzE4Kbu8Yc0y8KHh5auAz9fwrMRpZuReJHkf8McMYn+6vk4LI/aiqo5X1e6qWqyqRQbnM66sqi3/DZEdbJx/I59h8OyeJLsZvMTz/CyHnJFx9uLLwKUASX6QQfDXZjrlzrAMXDt8t87FwPGqemnUJ031JZ2qOpHkJuBBBmfgD1fV00luBVaqahm4k8GPZasMTlJcPc2Z5mXMvfgD4B3AXw7PW3+5qq6c29BTMuZetDDmXjwI/GySZ4D/AX67qk67n4LH3IuPAH+S5DcZnMC97nR8gpjkXgbf5HcPz1fcApwJUFV3MDh/cQWwCrwGXD/W456GeyVJ2oS/aStJTRh8SWrC4EtSEwZfkpow+JLUhMGXpCYMviQ1YfAlqYn/A0NQEvP+I7r7AAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# First set up the figure, the axis, and the plot element we want to animate\n",
    "fig, ax = plt.subplots()\n",
    "\n",
    "# on garde les bornes pour afficher les solutions\n",
    "u_min = 0 #min(u_tx[0,:])\n",
    "u_max = 2 # max(u_tx[0,:])\n",
    "Y_min = u_min - 0.1*(u_max-u_min)\n",
    "Y_max = u_max + 0.1*(u_max-u_min)    \n",
    "\n",
    "ax.set_xlim(( 0, 1))\n",
    "ax.set_ylim((Y_min, Y_max))\n",
    "\n",
    "line, = ax.plot([], [], lw=2)\n",
    "line2, = ax.plot([], [], lw=2, color='r')\n",
    "\n",
    "# initialization function: plot the background of each frame\n",
    "def init():\n",
    "    line.set_data([], [])\n",
    "    line2.set_data([], [])\n",
    "    return (line,line2)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'Nt' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-2-3a529cf2e56f>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m## schema numerique\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      2\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0;32mfor\u001b[0m \u001b[0mnt\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mNt\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      4\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0mimplicit_scheme\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m         \u001b[0mu_tx\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mnt\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0msparse\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mlinalg\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mspsolve\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mM\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mu_tx\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mnt\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'Nt' is not defined"
     ]
    }
   ],
   "source": [
    "## schema numerique  \n",
    "\n",
    "for nt in range(0,Nt):               \n",
    "    if implicit_scheme:\n",
    "        u_tx[nt+1,:] = sparse.linalg.spsolve(M, u_tx[nt,:])\n",
    "\n",
    "    else:\n",
    "        for i in range(0,Nx):        \n",
    "            u_tx[nt+1,i] = 0.9*u_tx[nt,i]   # todo\n",
    "            \n",
    "            \n",
    "    mass[nt+1] = h*numpy.sum(u_tx[nt+1,:])\n",
    "\n",
    "plot_sol(Nt, save_fn=\"u_final_\"+repr(Nt))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 208,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEKCAYAAADjDHn2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAGjBJREFUeJzt3Xu0FeWd5vHvg4LEEUMExlFRQRuiCCrmeEHbFqXbVkzrNIm3aBAnq+3R2NrRJF5iYmuWSWcSWT2ZZILE0Iq2xnsGBUPUjoNJJHowykUWikj0oNPQGvEWlctv/qg6uqXh1Os5u/auw34+a+11dl127V95wIe33rfeUkRgZmbWlT7NLsDMzKrPYWFmZoUcFmZmVshhYWZmhRwWZmZWyGFhZmaFHBZmZlbIYWFmZoUcFmZmVmjbZhdQL4MHD45hw4Y1uwwzs15lwYIF/x4RQ4r222rCYtiwYbS3tze7DDOzXkXS71P282UoMzMr5LAwM7NCDgszMyu01fRZmJmlWrduHR0dHbzzzjvNLqVh+vfvz9ChQ+nbt2+3Pu+wMLOW09HRwYABAxg2bBiSml1O6SKCV155hY6ODoYPH96tY/gylJm1nHfeeYdBgwa1RFAASGLQoEE9akk5LMysJbVKUHTq6fk6LMzMrJDDwszMCjkszMyskMPCzKwJVq5cyT777MOUKVMYOXIkZ5xxBg8++CBHHHEEI0aM4LHHHuOxxx5j3LhxjB07lsMPP5xly5YBsGTJEg455BAOPPBA9t9/f5599lneeustTjjhBA444ABGjx7NbbfdVtd6PXTWzFraVfcu4emXXq/rMUftuiNX/tV+hfstX76cO+64gxkzZnDwwQdzyy238Ktf/YpZs2bxrW99i5kzZ/LII4+w7bbb8uCDD3L55Zdz1113MW3aNC688ELOOOMM3nvvPTZs2MCcOXPYddddmT17NgBr166t6zk5LMzMmmT48OGMGTMGgP32248JEyYgiTFjxrBy5UrWrl3LWWedxbPPPosk1q1bB8C4ceO45ppr6OjoYNKkSYwYMYIxY8Zw8cUXc8kll/DpT3+aI488sq61OizMrKWltADKst12273/vk+fPu8v9+nTh/Xr1/P1r3+do48+mnvuuYeVK1cyfvx4AD73uc9x6KGHMnv2bCZOnMh1113HMcccwxNPPMGcOXO44oormDBhAt/4xjfqVqvDwsysotauXctuu+0GwA033PD++hUrVrDXXntxwQUX8MILL7Bw4UL22WcfdtppJ84880wGDhzI9ddfX9da3MFtZlZRX/3qV7nssssYO3Ys69evf3/97bffzujRoznwwANZvHgxkydPZtGiRe93el911VVcccUVda1FEVHXAzZLW1tb+OFHZpZi6dKl7Lvvvs0uo+E2d96SFkREW9Fn3bIwM7NCDgszMyvksDCzlrS1XIJP1dPzdViYWcvp378/r7zySssERufzLPr379/tY3jorJm1nKFDh9LR0cGaNWuaXUrDdD4pr7scFmbWcvr27dvtJ8a1Kl+GMjOzQqWFhaQZklZLWryF7ZL0fUnLJS2UdNAm23eU1CHpB2XVaGZmacpsWdwAHNfF9uOBEfnrHOBHm2z/JjCvlMrMzOwjKS0sImIe8GoXu5wEzIzMfGCgpF0AJH0K2Bn4RVn1mZlZumb2WewGvFiz3AHsJqkPcC3w5aZUZWZm/0EVO7jPA+ZEREfRjpLOkdQuqb2VhsCZmTVaM4fOrgJ2r1kemq8bBxwp6TxgB6CfpDcj4tJNDxAR04HpkE0kWH7JZmatqZlhMQs4X9JPgUOBtRHxMnBG5w6SpgBtmwsKMzNrnNLCQtKtwHhgsKQO4EqgL0BETAPmABOB5cDbwNll1WJmZj1TGBaSLgT+GXgDuB4YC1waEV2OVIqI0wu2B/DFgn1uIBuCa2ZmTZTSwf3fIuJ14FjgE8DngX8stSozM6uUlLBQ/nMicFNELKlZZ2ZmLSAlLBZI+gVZWMyVNADYWG5ZZmZWJSkd3F8ADgRWRMTbkgbhzmgzs5ZSGBYRsVHSvwGjJHlKczOzFpQyGuo7wKnA08CGfHXgSf7MzFpGSkvhvwKfjIh3yy7GzMyqKaWDewX5zXRmZtaaUloWbwNPSnoIeL91EREXlFaVmZlVSkpYzMpfZmbWolJGQ90oqR8wMl+1LCLWlVuWmZlVScpoqPHAjcBKsju3d5d0Vv4kPDMzawEpl6GuBY6NiGUAkkYCtwKfKrMwMzOrjpTRUH07gwIgIp7Bo6PMzFpKSsuiXdL1wM358hlAe3klmZlZ1aSExblkz53oHCr7CPC/S6vIzMwqJ2U01LvA1PxlZmYtaIthIen2iDhF0iKyuaA+JCL2L7UyMzOrjK5aFhfmPz/diELMzKy6tjgaKiJezt+eFxG/r30B5zWmPDMzq4KUobN/sZl1x9e7EDMzq66u+izOJWtB7C1pYc2mAcCvyy7MzMyqo6s+i1uA+4FvA5fWrH8jIl4ttSozM6uULYZFRKwF1kq6ZJNNO0jaISJeKLc0MzOripSb8maTDZ0V0B8YDiwD9iuxLjMzq5CUm/LG1C5LOgiPhjIzaykpo6E+JCKeAA4toRYzM6uolOdZXFSz2Ac4CHiptIrMzKxyUvosBtS8X0/Wh3FXOeWYmVkVpfRZXNWIQszMrLq6uinvXjYzgWCniDixlIrMzKxyumpZfK8nB5Y0g2wSwtURMXoz2wX8T2Ai8DYwJSKekHQg8CNgR2ADcE1E3NaTWszMrGe6uinv/3a+l9QPGJkvLouIdQnHvgH4ATBzC9uPB0bkr0PJAuJQsuCYHBHPStoVWCBpbkS8lvCdZmZWgpTRUOOBG4GVZDfm7S7prIiY19XnImKepGFd7HISMDMiApgvaaCkXfJnfHce4yVJq4EhgMPCzKxJUkZDXQscGxHLACSNBG4FPtXD794NeLFmuSNf1zk1OpIOAfoBz/Xwu8zMrAdSbsrr2xkUAPm//PuWV1JG0i7ATcDZEbFxC/ucI6ldUvuaNWvKLsnMrGWlhEW7pOsljc9f1wPtdfjuVcDuNctD83VI2pHsfo6vRcT8LR0gIqZHRFtEtA0ZMqQOJZmZ2eakhMW5wNPABflrSb6up2YBk5U5DFgbES/nnen3kPVn3FmH7zEzsx5KuSnvXWAqMFXSTsDQfF2XJN0KjAcGS+oAriS/fBUR04A5ZMNml5ONgDo7/+gpwJ8BgyRNyddNiYgn00/LzMzqKWU01MPAifm+C4DVkn4TEV/q6nMRcXrB9gC+uJn1NwM3F9VlZmaNk3IZ6uMR8TowiezS0KHAhHLLMjOzKkkJi23zkUmnAPeVXI+ZmVVQSlhcDcwFlkfE45L2Ap4ttywzM6uSlA7uO4A7apZXAJ8psygzM6uWj/ykPDMzaz0OCzMzK+SwMDOzQoVhIWlnST+RdH++PErSF8ovzczMqiKlZXED2WioXfPlZ4C/L6sgMzOrnpSwGBwRtwMbASJiPdkT7MzMrEWkhMVbkgaRP4+7c9K/UqsyM7NKSXn40UVkM8TuLenXZE+tO7nUqszMrFJSwmIJcBTwSbLHqi7Do6jMzFpKyv/0H42I9RGxJCIWR8Q64NGyCzMzs+rYYstC0n8heyb2xySNJWtVAOwIbN+A2szMrCK6ugz1l8AUssedTq1Z/wZweYk1mZlZxWwxLCLiRuBGSZ+JiLsaWJOZmVVMyqyzd0k6AdgP6F+z/uoyCzMzs+pIme5jGnAq8Hdk/RYnA3uWXJeZmVVIymiowyNiMvCHiLgKGAeMLLcsMzOrkpSw+GP+821JuwLrgF3KK8nMzKom5aa8+yQNBL4LPEE27cf1pVZlZmaVktLB/c387V2S7gP6R4TnhjIzayGFYSFpG+AEYFjn/pKIiKldfc7MzLYeKZeh7gXeARaRT1NuZmatJSUshkbE/qVXYmZmlZUyGup+SceWXomZmVVWSstiPnCPpD5kw2YFRETsWGplZmZWGSlhMZXsRrxFEREl12NmZhWUchnqRWCxg8LMrHWlhMUK4GFJl0m6qPNV9CFJMyStlrR4C9sl6fuSlktaKOmgmm1nSXo2f52VfjpmZlaGlLB4HngI6AcMqHkVuQE4rovtxwMj8tc5wI8AJO0EXAkcChwCXCnpEwnfZ2ZmJUm5g/uq7hw4IuZJGtbFLicBM/PLW/MlDZS0CzAeeCAiXgWQ9ABZ6NzanTrMzKznunqs6j9FxN9LupdsPqgPiYgTe/jdu5H1h3TqyNdtaX0pXn9nHZfcubCsw5uZlW7Y4P/EJcftU+p3dNWyuCn/+b1SK+gBSeeQXcJijz326NYxNm4MnlvzZj3LMjNrqL7bpPQo9ExXj1VdkL9dGhGra7dJ+mQdvnsVsHvN8tB83SqyS1G16x/eQo3TgekAbW1t3RqtNXD7fvziS0d156NmZi0jJY4ekXRK54Kki4F76vDds4DJ+aiow4C1EfEyMBc4VtIn8o7tY/N1ZmbWJCk35Y0Hpks6GdgZWEo2SqlLkm7NPztYUgfZCKe+ABExDZgDTASWA28DZ+fbXpX0TeDx/FBXd3Z2m5lZc6SMhnpZ0s+By8hmnb00Igov8kfE6QXbA/jiFrbNAGYUfYeZmTVGyvMsHgReAkaT9TH8RNK8iPhy2cWZmVk1pPRZ/CAiJkfEaxGxCDgc8JPyzMxaSMplqJ9J2hk4OF/1WM2jVs3MrAUUtizykVCPAScDpwC/lfTZsgszM7PqSBkN9TXg4M57LSQNAR4E7iyzMDMzq46UPos+m9yU90ri58zMbCuR0rL4uaS5fDCR36lk90iYmVmLSOng/oqkScCf5qumR0Q97uA2M7NeosuwkLQN8GBEHA3c3ZiSzMysarrse4iIDcBGSR9vUD1mZlZBKX0WbwKL8ocQvdW5MiIuKK0qMzOrlJSwuBtfgjIza2kpHdw3SuoH7EP2xLxlEfFe6ZWZmVllpEwkOBG4DngOEDBc0t9GxP1lF2dmZtWQchlqKnB0RCwHkLQ3MBtwWJiZtYiUO7Hf6AyK3ArgjZLqMTOzCkppWbRLmgPcTtZncTLweH6jHhHhzm8zs61cSlj0B/4NOCpfXgN8DPgrsvBwWJiZbeVSRkOd3YhCzMysulKeZ7GXpHslrZG0WtL/kbRXI4ozM7NqSOngvoWsv2IXYFfgDj6YgdbMzFpASlhsHxE3RcT6/HUzWT+GmZm1iJQO7vslXQr8lKxD+1RgjqSdACLi1RLrMzOzCkgJi1Pyn3+7yfrTyMLD/RdmZlu5lNFQwxtRiJmZVVdKywJJo4FR1PRVRMTMsooyM7NqSZlI8EpgPFlYzAGOB34FOCzMzFpEymiozwITgP+X36B3AOAn55mZtZCUsPhjRGwE1kvaEVgN7F5uWWZmViWpEwkOBH4MLCB7zOqjpVZlZmaVkjIa6rz87TRJPwd2jIiF5ZZlZmZVknIZ6n0RsfKjBIWk4yQtk7Q8v7Fv0+17SnpI0kJJD0saWrPtf0haImmppO9L0kep1czM6ucjhcVHIWkb4Idko6dGAadLGrXJbt8DZkbE/sDVwLfzzx4OHAHsD4wGDuaDKdLNzKzBthgWknp6M94hwPKIWBER75FNF3LSJvuMAv41f//Lmu1Bdk9HP2A7oC/ZMzXMzKwJumpZ3Akg6aFuHns34MWa5Y58Xa2ngEn5+78GBkgaFBGPkoXHy/lrbkQs7WYdZmbWQ111cPeRdDkwUtJFm26MiKl1+P4vAz+QNAWYB6wCNkj6E2BfoLMP4wFJR0bEI7UflnQOcA7AHnvsUYdyzMxsc7pqWZwGbCALlAGbeRVZxYfvxxiar3tfRLwUEZMiYizwtXzda2StjPkR8WZEvAncD4zb9AsiYnpEtEVE25AhQxJKMjOz7thiyyIilgHfkbQwIu7vxrEfB0bkfR+ryMLnc7U7SBoMvJrf9HcZMCPf9ALwN5K+DYisc/ufulGDmZnVQcpoqN9ImiqpPX9dK6lwuo+IWA+cD8wFlgK3R8QSSVdLOjHfbTywTNIzwM7ANfn6O4HngEVk/RpPRcS9H+nMzMysbhQRXe8g3QUsBm7MV30eOCAiJm35U43X1tYW7e3tzS7DzKxXkbQgItqK9kuZ7mPviPhMzfJVkp7sfmlmZtbbJE0kKOlPOxckHQH8sbySzMysalJaFv8dmFnTT/EH4KzySjIzs6pJmUjwKeCAfHpyIuL10qsyM7NKSXqsKjgkzMxaWWkTCZqZ2dbDYWFmZoUKw0LSyZIG5O+vkHS3pIPKL83MzKoipWXx9Yh4Ix8+++fAT4AflVuWmZlVSUpYbMh/ngBMj4jZZM+ZMDOzFpESFqskXQecCsyRtF3i58zMbCuR8j/9U8gmA/zLfPrwnYCvlFqVmZlVSsp9FrsAsyPiXUnjyZ6LPbPUqszMrFJSWhZ38cHT66aTPdDollKrMjOzSkkJi435sykmAf8rIr5C1towM7MWkRIW6ySdDkwG7svX9S2vJDMzq5qUsDib7PnX10TE8/ljUm8qtywzM6uSlFlnnwYuqFl+HvhOmUWZmVm1FIaFpBHAt4FRQP/O9RGxV4l1mZlZhaRchvpnsuk91gNHkw2bvbnMoszMrFpSwuJjEfEQoIj4fUT8A9nUH2Zm1iJSbsp7V1If4FlJ5wOrgB3KLcvMzKokpWVxIbA9WSf3p4DP42dwm5m1lJTRUI/nb98kG0ZrZmYtJmU0VBvwNWDP2v0jYv8S6zIzswpJ6bP4F7JZZhcBG8stx8zMqiglLNZExKzSKzEzs8pKCYsrJV0PPAS827kyIu4urSozM6uUlLA4G9iHbPLAzstQATgszMxaREpYHBwRnyy9EjMzq6yU+yx+I2lUdw4u6ThJyyQtl3TpZrbvKekhSQslPSxpaM22PST9QtJSSU9LGtadGszMrOdSWhaHAU9Kep6sz0JAFA2dlbQN8EPgL4AO4HFJs/JZbDt9D5gZETdKOoZswsLP59tmkk2L/oCkHfBILDOzpkkJi+O6eexDgOURsQJA0k+Bk4DasBgFXJS//yXws3zfUcC2EfEAQES82c0azMysDgovQ+WTB/6HV8KxdwNerFnuyNfVeorsca0Afw0MkDQIGAm8JuluSb+T9N28pWJmZk2Q0mdRpi8DR0n6HXAU2SSFG8haPEfm2w8G9gKmbPphSedIapfUvmbNmoYVbWbWasoMi1XA7jXLQ/N174uIlyJiUkSMJZtShIh4jawV8mRErIiI9WSXpw7a9AsiYnpEtEVE25AhQ8o6DzOzlldmWDwOjJA0XFI/4DTgQ3eCSxqcT38OcBkwo+azAyV1JsAxfLivw8zMGqi0sMhbBOcDc4GlwO0RsUTS1ZJOzHcbDyyT9AywM3BN/tkNZJegHpK0iGwE1o/LqtXMzLqmiGh2DXXR1tYW7e3tzS7DzKxXkbQgItqK9mt2B7eZmfUCDgszMyvksDAzs0IOCzMzK+SwMDOzQg4LMzMr5LAwM7NCDgszMyvksDAzs0IOCzMzK+SwMDOzQg4LMzMr5LAwM7NCDgszMyvksDAzs0IOCzMzK+SwMDOzQg4LMzMr5LAwM7NCDgszMyvksDAzs0IOCzMzK+SwMDOzQg4LMzMr5LAwM7NCiohm11AXktYAv+/BIQYD/16ncnqLVjvnVjtf8Dm3ip6c854RMaRop60mLHpKUntEtDW7jkZqtXNutfMFn3OraMQ5+zKUmZkVcliYmVkhh8UHpje7gCZotXNutfMFn3OrKP2c3WdhZmaF3LIwM7NCLRUWko6TtEzSckmXbmb7dpJuy7f/VtKwxldZXwnnfJGkpyUtlPSQpD2bUWc9FZ1zzX6fkRSSev3ImZRzlnRK/rteIumWRtdYbwl/tveQ9EtJv8v/fE9sRp31ImmGpNWSFm9huyR9P//vsVDSQXUtICJa4gVsAzwH7AX0A54CRm2yz3nAtPz9acBtza67Aed8NLB9/v7cVjjnfL8BwDxgPtDW7Lob8HseAfwO+ES+/J+bXXcDznk6cG7+fhSwstl19/Cc/ww4CFi8he0TgfsBAYcBv63n97dSy+IQYHlErIiI94CfAidtss9JwI35+zuBCZLUwBrrrfCcI+KXEfF2vjgfGNrgGust5fcM8E3gO8A7jSyuJCnn/DfADyPiDwARsbrBNdZbyjkHsGP+/uPASw2sr+4iYh7wahe7nATMjMx8YKCkXer1/a0UFrsBL9Ysd+TrNrtPRKwH1gKDGlJdOVLOudYXyP5l0psVnnPePN89ImY3srASpfyeRwIjJf1a0nxJxzWsunKknPM/AGdK6gDmAH/XmNKa5qP+ff9Itq3Xgax3k3Qm0AYc1exayiSpDzAVmNLkUhptW7JLUePJWo/zJI2JiNeaWlW5TgduiIhrJY0DbpI0OiI2Nruw3qiVWhargN1rlofm6za7j6RtyZqurzSkunKknDOS/hz4GnBiRLzboNrKUnTOA4DRwMOSVpJd253Vyzu5U37PHcCsiFgXEc8Dz5CFR2+Vcs5fAG4HiIhHgf5kcyhtrZL+vndXK4XF48AIScMl9SPrwJ61yT6zgLPy958F/jXynqNeqvCcJY0FriMLit5+HRsKzjki1kbE4IgYFhHDyPppToyI9uaUWxcpf7Z/RtaqQNJgsstSKxpZZJ2lnPMLwAQASfuShcWahlbZWLOAyfmoqMOAtRHxcr0O3jKXoSJivaTzgblkIylmRMQSSVcD7RExC/gJWVN1OVlH0mnNq7jnEs/5u8AOwB15X/4LEXFi04ruocRz3qoknvNc4FhJTwMbgK9ERK9tNSee88XAjyV9iayze0pv/sefpFvJAn9w3g9zJdAXICKmkfXLTASWA28DZ9f1+3vxfzszM2uQVroMZWZm3eSwMDOzQg4LMzMr5LAwM7NCDgszMyvksDArkaSBks5rdh1mPeWwMCvXQLLZjM16NYeFWbn+Edhb0pOSvtvsYsy6yzflmZUof4DWfRExusmlmPWIWxZmZlbIYWFmZoUcFmbleoNsWnSzXs1hYVaifGbXX0ta7A5u683cwW1mZoXcsjAzs0IOCzMzK+SwMDOzQg4LMzMr5LAwM7NCDgszMyvksDAzs0IOCzMzK/T/Aa/l7hoHh6CUAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot_t_curves(array1=mass, label1=\"mass\", name=\"mass of approximate solutions\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'numpy' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-3-a12facf460f9>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m# estimations d'erreurs locales:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m \u001b[0merr_loc\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnumpy\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mzeros\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mNt\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mNx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      3\u001b[0m \u001b[0merr_loc_norm\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnumpy\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mzeros\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mNt\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0;32mfor\u001b[0m \u001b[0mnt\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mNt\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'numpy' is not defined"
     ]
    }
   ],
   "source": [
    "# estimations d'erreurs locales:\n",
    "err_loc = numpy.zeros((Nt+1,Nx))\n",
    "err_loc_norm = numpy.zeros((Nt+1))\n",
    "\n",
    "for nt in range(0,Nt):\n",
    "    for i in range(0,Nx):                    \n",
    "        ip = numpy.mod(i+1,Nx)\n",
    "        im = numpy.mod(i-1,Nx)\n",
    "        err_loc[nt,i] = 0        # todo\n",
    "\n",
    "    # norme (L2) des erreurs locales, pour l'instant t = nt * dt\n",
    "    tmp_sum = 0\n",
    "    for i in range(0,Nx):\n",
    "        tmp_sum += err_loc[nt,i]**2\n",
    "    err_loc_norm[nt] = numpy.sqrt(h*tmp_sum)\n",
    "        \n",
    "    # err_loc_norm[nt+1] = numpy.sqrt((h**2)*numpy.sum(u_tx[nt+1,:]**2))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### stockage de quelques erreurs pour comparaison"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 162,
   "metadata": {},
   "outputs": [],
   "source": [
    "err_loc_norm_1 = err_loc_norm\n",
    "label_1 = repr(Nt)+\"_\"+repr(Nx)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 175,
   "metadata": {},
   "outputs": [],
   "source": [
    "err_loc_norm_2 = err_loc_norm\n",
    "label_2 = repr(Nt)+\"_\"+repr(Nx)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 210,
   "metadata": {},
   "outputs": [],
   "source": [
    "err_loc_norm_3 = err_loc_norm\n",
    "label_3 = repr(Nt)+\"_\"+repr(Nx)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'err_loc_norm_1' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-4-d645cfb7c8b7>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m plot_t_curves(\n\u001b[0;32m----> 2\u001b[0;31m     \u001b[0marray1\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0merr_loc_norm_1\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlabel1\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mlabel_1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      3\u001b[0m     \u001b[0marray2\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0merr_loc_norm_2\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlabel2\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mlabel_2\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m     \u001b[0marray3\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0merr_loc_norm_3\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mlabel3\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mlabel_3\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m     name=\"local error norms\", logscale=True)\n",
      "\u001b[0;31mNameError\u001b[0m: name 'err_loc_norm_1' is not defined"
     ]
    }
   ],
   "source": [
    "plot_t_curves(\n",
    "    array1=err_loc_norm_1, label1=label_1, \n",
    "    array2=err_loc_norm_2, label2=label_2, \n",
    "    array3=err_loc_norm_3, label3=label_3, \n",
    "    name=\"local error norms\", logscale=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### animation de la dernière solution calculée"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {},
   "outputs": [],
   "source": [
    "# animation function\n",
    "def animate(nt):\n",
    "    line.set_data(X, u_tx[nt,:])\n",
    "    #line2.set_data(X, u_ex(nt*dt,X))\n",
    "    return (line,line2)\n",
    "#plt.show()    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'animation' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-5-6bab59656043>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      4\u001b[0m \u001b[0;31m# see https://matplotlib.org/api/_as_gen/matplotlib.animation.FuncAnimation.html\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 6\u001b[0;31m anim = animation.FuncAnimation(fig, animate, init_func=init,\n\u001b[0m\u001b[1;32m      7\u001b[0m                                frames=Nt+1, interval=20, blit=True)\n\u001b[1;32m      8\u001b[0m \u001b[0mHTML\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0manim\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mto_html5_video\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'animation' is not defined"
     ]
    }
   ],
   "source": [
    "# call the animator. \n",
    "# interval: Delay between frames in milliseconds. Defaults to 200.\n",
    "# blit=True means only re-draw the parts that have changed.\n",
    "# see https://matplotlib.org/api/_as_gen/matplotlib.animation.FuncAnimation.html\n",
    "\n",
    "anim = animation.FuncAnimation(fig, animate, init_func=init,\n",
    "                               frames=Nt+1, interval=20, blit=True)\n",
    "HTML(anim.to_html5_video())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "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.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}