kap_perf.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def show_book(figname):\n",
    "    print(r\"Siehe Abbildung §\\ref{fig_XXX}§\".replace(\"XXX\", figname))\n",
    "    plt.savefig(f\"{figname}_snip.pdf\")\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Wärmeleitung"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin heat_setup`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = 0.1  # K/(m²s)\n",
    "\n",
    "# Abtastzeit der Iteration ist 2 ms\n",
    "dt = 0.002\n",
    "# Gesamtzeit ist 1 s\n",
    "T = 1.0\n",
    "steps = int(T/dt) # Anzahl der Zeitschritte\n",
    "\n",
    "# Gitterdimension 1 x 1 m\n",
    "grid_dim = (1, 1)\n",
    "# 25 Gitterpunkte auf jeder Achse\n",
    "grid_samples = (25, 25)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Max dt: 0.003999999999999999\n",
      "Sim should converge\n"
     ]
    }
   ],
   "source": [
    "dx, dy = grid_dim[0] / grid_samples[0], grid_dim[1] / grid_samples[1]\n",
    "max_dt = 1/(2*a)*(dx*dy)**2/(dx**2+dy**2)\n",
    "print('Max dt:', max_dt)\n",
    "if dt >= max_dt:\n",
    "    print(\"Sim will freak out\")\n",
    "else:\n",
    "    print(\"Sim should converge\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "def plot_heat_grid(u, grid_dim):\n",
    "    plt.imshow(u.T, cmap='coolwarm', origin='lower', extent=(0, grid_dim[0], 0, grid_dim[1]))\n",
    "    plt.colorbar(label=\"Temperatur [°C]\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin heat_sim_numpy`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "def sim_heat_equation_numpy(a, grid_dim, grid_samples, dt, steps):\n",
    "    dx, dy = grid_dim[0] / grid_samples[0], grid_dim[1] / grid_samples[1]\n",
    "    # Anfangsbedingungen: alles 0°C\n",
    "    u = np.zeros(grid_samples)\n",
    "    \n",
    "    # Definiere Maske um Temperatur auf Scheibe in Schleife zu setzen\n",
    "    disk_mask = np.zeros_like(u, dtype=bool)\n",
    "    for xi in range(grid_samples[0]):\n",
    "        for yi in range(grid_samples[1]):\n",
    "            # Scheibe mit Radius 0.1 und Mittelpunkt (0.5, 0.5)\n",
    "            if (xi*dx - 0.5)**2 + (yi*dy - 0.5)**2 < 0.1**2:\n",
    "                disk_mask[xi,yi] = True\n",
    "\n",
    "    for i in range(steps):  # Simulationsschleife\n",
    "        # Setze Randwerte\n",
    "        u[0,:] = u[1,:]  # Neumann am linken Rand\n",
    "        u[-1,:] = u[-2,:] # Neumann am rechten Rand\n",
    "        u[:,0] = 0 # Dirichlet (0°C) am unteren Rand\n",
    "        u[:,-1] = 50 # Dirichlet (50°C) am oberen Rand\n",
    "        u[disk_mask] = 50  # Erhitzte Scheibe (50°C)\n",
    "        # Diskretisierte Wärmeleitungsgleichung\n",
    "        u[1:-1,1:-1] += a*dt*((u[2:,1:-1] + u[:-2,1:-1] - 2*u[1:-1,1:-1])/(dx**2) +\n",
    "                              (u[1:-1,2:] + u[1:-1,:-2] - 2*u[1:-1,1:-1])/(dy**2))\n",
    "    \n",
    "    return u"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin heat_sim_python_loops`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "def sim_heat_equation_python_loops(a, grid_dim, grid_samples, dt, steps):\n",
    "    dx, dy = grid_dim[0] / grid_samples[0], grid_dim[1] / grid_samples[1]\n",
    "    # Anfangsbedingungen: alles 0°C\n",
    "    u = np.zeros(grid_samples)\n",
    "\n",
    "    for i in range(steps):  # Simulationsschleife\n",
    "        for xi in range(grid_samples[0]):\n",
    "            for yi in range(grid_samples[1]):\n",
    "                if xi == 0: # Neumann am linken Rand\n",
    "                    u[0,yi] = u[1,yi]\n",
    "                elif xi == grid_samples[0] - 1: # Neumann am rechten Rand\n",
    "                    u[-1,yi] = u[-2,yi]\n",
    "                elif yi == 0: # Dirichlet am unteren Rand\n",
    "                    u[xi,0] = 0\n",
    "                elif yi == grid_samples[1] - 1: # Dirichlet am oberen Rand\n",
    "                    u[xi,-1] = 50\n",
    "                elif (xi*dx - 0.5)**2 + (yi*dy - 0.5)**2 < 0.1**2:\n",
    "                    u[xi,yi] = 50  # Erhitzte Scheibe (50°C)\n",
    "                else:\n",
    "                    # Diskretisierte Wärmeleitungsgleichung\n",
    "                    u[xi,yi] += a*dt*((u[xi+1,yi] + u[xi-1,yi] - 2*u[xi,yi])/(dx**2) +\n",
    "                                      (u[xi,yi+1] + u[xi,yi-1] - 2*u[xi,yi])/(dy**2))\n",
    "    \n",
    "    return u"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin heat_sim_python_loops_numba`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from numba import jit\n",
    "#!@jit\n",
    "#!def sim_heat_equation_numba(a, grid_dim, grid_samples, dt, steps):\n",
    "#!    # Code aus sim_heat_equation_python_loops 1:1 kopieren\n",
    "#!    ..."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "@jit\n",
    "def sim_heat_equation_numba(a, grid_dim, grid_samples, dt, steps):\n",
    "    dx, dy = grid_dim[0] / grid_samples[0], grid_dim[1] / grid_samples[1]\n",
    "    u = np.zeros(grid_samples)\n",
    "\n",
    "    for i in range(steps):  # Simulationsschleife\n",
    "        for xi in range(grid_samples[0]):\n",
    "            for yi in range(grid_samples[1]):\n",
    "                if xi == 0: # Neumann am linken Rand\n",
    "                    u[0,yi] = u[1,yi]\n",
    "                elif xi == grid_samples[0] - 1: # Neumann am rechten Rand\n",
    "                    u[-1,yi] = u[-2,yi]\n",
    "                elif yi == 0: # Dirichlet am unteren Rand\n",
    "                    u[xi,0] = 0\n",
    "                elif yi == grid_samples[1] - 1: # Dirichlet am oberen Rand\n",
    "                    u[xi,-1] = 50\n",
    "                elif (xi*dx - 0.5)**2 + (yi*dy - 0.5)**2 < 0.1**2:\n",
    "                    u[xi,yi] = 50\n",
    "                else:\n",
    "                    u[xi,yi] += a*dt*((u[xi+1,yi] + u[xi-1,yi] - 2*u[xi,yi])/(dx**2) +\n",
    "                                      (u[xi,yi+1] + u[xi,yi-1] - 2*u[xi,yi])/(dy**2))\n",
    "    \n",
    "    return u"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin timeit_python_loops`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.87 s ± 124 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
     ]
    }
   ],
   "source": [
    "#!%timeit sim_heat_equation_python_loops(a, grid_dim, grid_samples, dt, steps)\n",
    "print(\"1.87 s ± 124 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\")#!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin timeit_numpy`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "31.8 ms ± 11.2 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)\n"
     ]
    }
   ],
   "source": [
    "#!%timeit sim_heat_equation_numpy(a, grid_dim, grid_samples, dt, steps)\n",
    "print(\"31.8 ms ± 11.2 ms per loop (mean ± std. dev. of 7 runs, 10 loops each)\")#!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin timeit_numba`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "4.05 ms ± 86.6 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n"
     ]
    }
   ],
   "source": [
    "# Erster Aufruf ist Aufwärmung und kompiliert die Funktionenn\n",
    "#!sim_heat_equation_numba(a, grid_dim, grid_samples, dt, steps)\n",
    "# Jetzt (nach Kompilierung) ist Zeitmessung sinnvoll\n",
    "#!%timeit sim_heat_equation_numba(a, grid_dim, grid_samples, dt, steps)\n",
    "print(\"4.05 ms ± 86.6 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\")#!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Siehe Abbildung §\\ref{fig_heat_sim}§\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQ8AAADQCAYAAAAZMORwAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAeeUlEQVR4nO2df7RdZXnnP99zzr25+UmAQIgkCJZYG3EEzOLH6FREagN2CGvKsMiMHZxhmWkrLlus07R2QKNrqtNVdewwtXcKi+hSIXW0ZFHUMgjSYkFSUSBQNKYgwUiGHwk/k9xzzjN/7H3hcNjvPnufs8/e5yTPZ6297t7v773v2c9+3+d5n/eVmeE4jpOXWtUNcBxnPHHh4ThOX7jwcBynL1x4OI7TFy48HMfpCxcejuP0xdCFh6RrJO2WdH8gXpI+J2m7pHslnTrsNjnOwYCkhyXdJ+kHkrbGYUdIulnSj+O/hw+r/jJ6HtcCa1LizwVWxsd64M9LaJPjHCy8w8xONrPV8fUG4BYzWwncEl8PhaELDzO7HXgqJcla4AsWcSewWNKyYbfLcQ5S1gKb4vNNwAXDqqgxrIJzcCzwaMf1zjhsV3dCSeuJeidMobcsZ7KUBjpOJ7uZYa+1NEgZb6nNt2esFYzfzv5twL6OoGkzm+5KZsDfSjLgL+L4pWY2++78HFg6SDvTGAXhkZn44UwDrNSUfab+2opb5ByK/G7rkYHLeMZafLYR/v3+WvNH+zqGIiHeZmaPSToauFnSP3VGmpnFgmUojILweAxY0XG9PA5znIMW1UR9borW4NneZZjZY/Hf3ZK+DpwGPC5pmZntiof/uwtpcAKjYKrdAvyH2OpyBrC3o9vlOAcnNajPrQePXkiaL2nh7DnwLuB+ovfpkjjZJcANQ7qD4fc8JH0FOAtYImkncCUwAWBmnwduAs4DtgMvAP9x2G1ynMoRaGIgtclS4OuSIHqPv2xm35R0N7BZ0qXAI8BFA7c1wNCFh5mt6xFvwPuH3Q7HGSV6Dlt6YGY7gDcnhD8JvHOApmVmFHQefVGfqrHodfMqq7/dbufOYzP5dVftZr48/dRRJP18TWuNgb7AA9efl9rDBYz2Vex9V8HYCg/HGWckqE+Ogsqxf1x4OE4VSNQmeitGRxkXHo5TAQJqdR+2OI6TF4n6hA9bHMfJiYQPWxzH6QNF5tpxZmyFR22izsLXLBh6PRYwyVor2STabodNpRaIax0IO0iF6gmVlUY7UFY/hMbraS+EAnnSxv6q5+vap9VfK+hlre8cfLghH7Y4jtMvqrnwcBwnL97zcBynHySoNVxh6jhOXlxh6jhOf4haw4ctlVCfqLHg6GKsLWmWi1Bcu5lsIUmzaITytGbCTnahPEVae/oh9NVMs2iE5jWkfYFDeoFQnjQlZFEzOot46SVXmDqO0w+uMHUcpx9cYeo4Tt+4wtRxnPxI3vNwHKcf5ApTx3HyI+X32xk1xlZ41Bp15i1ZlCtP0MmtQFNt60AzWFY/eUJm3FCeNFNx6F5CzwXC5sTQeD3NgtCYk/xzS+u+h8yioTypjnk57yVvm3LhwxbHcfrFhy2O4+QmMtWOt/AopfWS1kh6SNJ2SRsS4o+TdKukeyTdK+m8MtrlONUhqNXCR9ZSpHr83twYX58g6a74Xbte0tB2gx+68JBUB64CzgVWAeskrepK9kfAZjM7BbgY+F/DbpfjVIqgVq8Hjxx8EHiw4/pTwGfM7ETgaeDSAlv9CsroeZwGbDezHWZ2ALgOWNuVxoBZ7edhwM9KaJfjVIeEGvXgka0ILQfeDfxlfC3gbOCrcZJNwAXFNz6iDJ3HscCjHdc7gdO70nwU+FtJHwDmA+ckFSRpPbAeYPlhC5g64tXWlmj3ynxYwAoCaQ5wgfAUy0lz34HE8DRrS3PfTGJ4yMmrNZN2LyGrSvgbEnSAC/zAQxYVgPpkclwoPC0urxUGUhSUOa0tytczCJeTXu8SSVs7rqfNbLorzWeB/wIsjK+PBPaY2ewPaifR+zcURkVhug641sz+VNKZwBclnWRmr/i1xw9vGuDk1xxV7b6KjjMAknoJoSfMbHVK/l8DdpvZP0o6q+DmZaIM4fEYsKLjenkc1smlwBoAM/sHSVPAEmB3Ce1znEoY0FT7VuD82LgwRTTs/x/AYkmNuPeR9K4VRhk6j7uBlbEWeJJIIbqlK81PiXf2lvRLRA/j/5XQNsephrjnETp6YWZ/YGbLzex4onfq22b274FbgQvjZJcANwzrFoYuPGIJeBnwLSKt8GYz2yZpo6Tz42QfAt4n6YfAV4D3Wj/KC8cZIwZVmAb4feBySduJdCBXF9LYBErReZjZTcBNXWFXdJw/QNQNc5xDAklExpHBMbPbgNvi8x1EFs6hMyoK09yoUWdicT7fFoIbOPVhbZkJ+JbsT7aoQNgSMPPC/mCeEEGfmxTfFtXyd+ZCbQ75sKRZThpTE7nCAepTyXOc6hMBK8xkuKzQcCCvb0tRDm0D9jAqZ2yFh+OMNdEiplW3YiBceDhORRQ1X6QqXHg4ThXEM0zHGRcejlMVBSlM81erUzMkmzGz+9ISuPBwnCqQoLphy3eI5l+lSa8TgOPTCnHh4ThVIFC9stfvbjM7Oy2BpG/3KmRshYdqNRqLFvZO2EFw3lmKqdZmkh3TQibZVto4to+l9kOm4uBCMilW3352jAstUahafie3kEl2Yv7ccJ65cxLDa3OSTbh9mWpzDh9yusyHWtPX76EIegmOrGnGVng4zjgTLYBczbBF0q8CC83sq13hFwJ7zezmLOWMt6HZccaWeJ5H6BguVxDpPbq5DdiYtRDveThOVVSnMJ1jZq9yPDWzJyTNz1qICw/HqYJqrS2LOtz2O5qkCSCsgOrChy2OUxXVDVu+Bvzvzl6GpAXA5+O4TIxtz0ONOrXDFufKE9zcqBleBpBmsrVF+/Ylh/fxNUlbBjHkmBaynKQuabg/Oc5SnOnqgWUFJwLfpzQns5CTW8iiAtAIWGJqU1PJ4QErTNSAYqwteVY3D1Jtz+OPgE8Aj0h6JA47jsh9/79mLWRshYfjjD0VzTCNhysbJH0MODEO3m5mL+Ypx4WH41SBBNVNEkPS6wH1moKehgsPx6mKWmXzPD5J9O63JcnMPtxPOS48HKcSVNmwhWi94LfF59/ttxAXHo5TASaw6hSmv0e0IZSAP+i3EBcejlMJ1a0kZmbfIXmGaS7GVnhYrY4tyLuGabJJUgFzLAAzyQ5wtUbAASvFjFcPOOa1UtY9DZk+QybZF58OK8wPPJ9ikg4wOT9kqk2+/6k0U21g3dH63GSzK0B9fvKEx9pUsnlXARMuAKH/Wd7hQ0E9BqtO53GMmf180DQ+ScxxqkCKFKahY7jc1DtJ7zRj2/NwnLGniMlm/fFmSc+kxAtIiwdK6nlIWiPpIUnbJW0IpLlI0gOStkn6chntcpzKkKKhd+DonV1Tkr4n6YfxO/OxOPwESXfF79r18S6Nr8DM6ma2KOVYaGY9N8geuvCQVAeuAs4FVgHrJK3qSrOSSOv7VjN7I/A7w26X41SNqRY8MrAfONvM3gycDKyRdAbwKeAzZnYi8DTRPtBDoYyex2lEU193mNkB4DpgbVea9wFXmdnTAGbmG1w7BzmD9Tws4rn4ciI+DDgbmF3kZxNwwRAaD5Sj8zgWeLTjeidwelea1wNIugOoAx81s292FyRpPbAeYMXSJbTnL87XEgssqdcKWyF0IOAA1wjsWJa2tFxgucP6i8l1QNgxLuTktnfH88Gyms+EHfBCNBYl1z/38IDDWsoyjCGrSn1e2Au8Nn9eYrjmJVthNBl2srOJgNNcXpNpEdaWWYVpmCWStnZcT5vZ9CuLUB34RyL/lKuAnwB7OlztdxK9f0NhVBSmDWAlcBawHLhd0pvMbE9novjhTQOc+oYTfSNsZ2wxwNJNxE+Y2erUMsxawMmSFgNfB96Qtf5Y8Gwzs8x5uilj2PIYsKLjenkc1slOYIuZzZjZPwM/IhImjnPQYqoHj1zlRB/ZW4EzgcWSZjsFSe/abJ4W8JCk4/ptfxnC425gZawFngQuBrZ0pflrol4HkpYQDWN2lNA2x6kIDaQwlXRU3ONA0lzgV4AHiYTIhXGyS4AbUoo5HNgm6RZJW2aPrHcw9GGLmTUlXQZ8i0ifcY2ZbZO0EdhqZlviuHdJegBoAR82syeH3TbHqQwNPMN0GbApHn7UgM1mdmP8Dl0n6RPAPUQL/ITIvPBPEqXoPMzsJrpmrJnZFR3nBlweH45z0GOxtaXv/Gb3AqckhO8gsnBmKWMg/5ZRUZjmxmp1ZuYdljNTso613kzxLZkIWAjqyX4SSuly1gLWlkaKtaUeWFYv6PPyYmCpRaA9k1/HHCovVH8jsNQghJcUrC9cEMxTW5C8sZfNC4RPhn1bLGAhy/0SFzR9PON8jqEh6Vki3S3AJJG593kzy+Q0NrbCw3HGG9HOqRgtGjN7SQIrWsh1LXBG1vzuGOc4VSDibeMCR8nEk87+GvjVrHm85+E4FWAj0POQ9G86LmvAaiA8hu7ChYfjVETVwgP41x3nTeBhXu06EsSFh+NUgipXmAJ/aWZ3dAZIeiuQybes8tY7zqGIAW3VgkdJ/FnGsETGtufRrtXZN3V4QkzYHKmAY1y9HXaMqzeTh4AW2HOjnuI0VQvVfyBsKp7Y+2xi+OT8ZAewWiOsbMvvFhcuL1T/xMJkRzaAxqJk82ptYdgyaAsXJ4a3p5Id49qTYSc7C5jX2zlNr4UsH6jqdB6SzgT+JXCUpM65VYuIJnJmYmyFh+OMOxUOWyaBBUTvf6dEf4aXp7b3pKfwkHQz8Htm9sO8LXQcJ5kqrS2zq6dLutbMHumZIUCWnsfvA5+V9DDwh2a2q9/KHMd5GaOyTZ9meUHSnwBvBF6ammtmZ2fJ3LPfZGbfN7N3ADcC35R0ZezF5zjOALSpBY+S+BLwT8AJwMeITLV3Z82cqZXx1NWHgD8HPgD8WNJv5G2p4zgRhmhTDx4lcaSZXQ3MmNl3zOw/ES1jmIksOo87iCTTNuBO4L1E0uqDkv6Vma3vq9kD0ladFybzbfqkgGNcox22djQayZ2syXqyA9hEQKMP0AhMO661w3aQOS8kb+I0b8lTyeHLw45hLz6+PxgXYu7SZKvKvCXJlpM5S44IllVbnGQdAzvsyGCe1rzk/3FzKtmZrtUI33+rVoy1pV0rxs4wAsOW2d3Odkl6N/AzIPwP7CLLU1gPPBC7zXfyAUkPZq3IcZxOVObwJMQnJB0GfIhofsci4HezZu4pPMxsW0r0u7NW5DjOKzGrrucRLyK00sxuBPYC78hbxkCiL154xHGcnBjVKkzjNUzXDVKGTxJznEoQbat82HKHpP8JXA+8tGeHmX0/S2YXHo5TAVHPo3KF6cnx340dYbMbR/XEhYfjVIFRec8jnr/VN2MrPNrUeK6dbC4MoYDTXKMWNpVO1JPNuDP1ZBPmZIqpcCrgTDeR4uPQCMQtDqwh+ospjmkzz2de5+UlJuYn38+8Fa9JDG8csyxYVvvwJYnhzfnJJlyAmalkU+2BiWTHuGbAHBvFJZvXe2y+9CqKmlZetalW0lLgvwGvMbNz4z2kz4znfvSk8kGX4xyKWKzzCB0lcS3RtiezX4IfkWOTeRcejlMRLVPwKIklZrYZaEO0xxI5Vm4oRXhIWiPpIUnbJW1ISffrkkxS6h6djnMwYKbg0QtJKyTdKukBSdskfTAOP0LSzZJ+HP8NjwnheUlHEi+CI+kMojkfmRi68Igno1wFnAusAtbFY6vudAuBDwJ3DbtNjlM1BQxbmsCHzGwV0XYJ74/fqw3ALWa2Erglvg5xOdHWr78Qu6F8gch3LRNl9DxOA7ab2Q4zOwBcR/Iiqx8HPkWO1ZsdZ5wxFDx65jXbNTsfw8yeJdqn9liid2tTnGwTcEFKGd8H3k60qth/Bt4Y70SXiTKsLccCj3Zc7wRO70wg6VRghZn9jaQPhwqStJ7I14aly1bwzMyrLQtp0lAK7BinFGtLLbDcXi3ZCjGVsmNZM2ChmZtioZkMLKvXmJ9shVh09DHBsmx/frmsOcltswUBh7UF4V7ygXnJcQcmwzvG7ZtIjjug5GfZtBRriyVbSfJOEy/C69WMXrqNJZK2dlxPm9l0UkJJxxNtPXkXsLRjzZ2fA0tDFUiaAn4beBvR0OXvJH3ezDL9UCo31Sran/HTRN66qcQPbxrgDSedmn/vRMcZIdrtVOHxhJn11P1JWgD8H+B3zOwZdZidzcwU+mJGfAF4lpcXPf53wBeBf9urXihHeDwGrOi4Xh6HzbIQOAm4Lb7xY4Atks43s07J6zgHERp4hqmkCSLB8SUz+1oc/LikZWa2S9Iy0rdROCnWmcxyq6QHstZfhs7jbmClpBMkTQIXEylpADCzvWa2xMyON7PjidYMccHhHNQYUc8jdPQiXqDrauBBM/t0R9QW4JL4/BLghpRivh9bWGbLPB3I/N4NvedhZk1JlxFNRqkD15jZNkkbga1mtiW9BMc5OBlwhulbgd8A7pP0gzjsD4FPApslXQo8AlyUUsZbgO9K+ml8fRzwkKT7iEY9/yKtAaXoPMzsJuCmrrArAmnPKqNNjlMlZtDK0MMI57e/h6D0eWfGYtb03QBGQGHaL622eHZ/ssY9RMiFoVFL3owJoB6ImxPwh3mxnuLbUku2nOybl+ynATA157Dk+ucn+4lMHPF8YjhArRlebjFEO7CsYmsi+V72p1hO9teT/W72E35m+y35f3ygldyuZjtsCQlZN/JOB2+1ixntlziTNBEzeySeRLaCDlngLvmOM+IEltQtDUkfJ7Jy/oSXt1p0l3zHGWXMNNCwpSAuAn4hnryZG3eMc5yKMAsfJXE/sLjfzN7zcJwKMKrXeQB/DNwj6X7gpX05zOz8LJldeDhORbTDevqy2ETkT3YfsVt+Hlx4OE4FmEG7+p7HC2b2uX4zj63waLXFnufzNT9kqq2naH4a9eQB6EQ9WVBPNsICfE7AMW4qxTFuTi3Z9DlnKtnJrDF3JjEcoE4zGBci5ATWtORnf8CSl/oD2N9Kjptphf+P+wOm15lWcniaEjIUl+YYl/TfbBb00o9Az+PvJP0x0azUzmGLm2odZ1QxBpskVhCnxH/P6AhzU63jjDTlWlWSmzDg6uluqnWcCoh6HuGjDCQtlXS1pG/E16tin5hMuPBwnIoYgXke1+KrpzvOmGHV9TwkzaorBlo9fWx1Hq02PP1cPoVTYJ8kaikitB7INNlIzhSyzgBMTiQ/7sl62EIx2Uj+X4asPWlOfjUlx6V9QUITmZoBZ7I0p7GZVnLcgUA4wEwzUH8rXziEv+jtnF/6VkodWTGglfk1LZzvAacy4OrpYys8HGfcqVBhOiv9uldPPwq4MGshLjwcpwKi9Twqkx5HSbo8Pv860Vo7IprrcQ6QaQV1Fx6OUxEV9jzqwAJevZhQeKPjBFx4OE4FmFWq89hlZhsHLcSFh+NURFnzORIoZGqrCw/HqQirTueRdY3TVMZWeDRbxp69r+73KeT9lkI9ZQOwWsBUG8rTqIfrn2gkxwUsuHF5yRU1gvWHf5AhU3Uaod93yKkrzV+jGeimh8LT4kJd/rSveVEKyrT2ZsWsup6HmT1VRDljKzwcZ9xpV9fzKAQXHo5TAWYj4ZI/EKVMT5e0RtJDkrZL2pAQf7mkByTdK+kWSa8to12OUx1GqxU+eiHpGkm74yUEZ8OOkHSzpB/Hf8O7jhfA0IWHpDpwFXAusApYJ2lVV7J7gNXxDlVfBf77sNvlOFVigJkFjwxcy6s3bdoA3GJmK4Fb4uuhUUbP4zRgu5ntiJd4vw5Y25nAzG41sxfiyzuJNsN2nIMXY6Ceh5ndDnQrPtcSrUtK/PeCQtvcRRk6j2OBRzuudwKnp6S/FPhGUoSk9cB6gEVHrGDPnvCSe3lIdYwLrFFYC1hVQo50UVnJcY2AFSYtLrR0olLupdaHJSpEO/B1TPvhhxSEzWY4T9CqEqgn7avdj34yyZya1t7M5fZuzxJJnZtOT5vZdI9il5rZrvj858DS/lvYm5FSmEp6D7AaeHtSfPzwpgGOee2p462qdg5tDNrpPYwnzGx138WbmaShviNlCI/HiPbCnGV5HPYKJJ0DfAR4u5nt7453nIMJM6NV/ESPxyUtM7NdkpYBu4uuoJMydB53AyslnSBpEriYyA34JSSdAvwFcL6ZDfWGHWdUsHb46JMtwCXx+SXADUW0M8TQex5m1pR0GdFyZ3XgGjPbJmkjsNXMtgB/QuTl91fxDNGfZt21ynHGldYAEz0kfQU4i0g3shO4EvgksDleh/QRor1oh0YpOg8zu4lozYDOsCs6zs8pox2OMypEk8T6V0mY2bpAVCF+K1kYKYVpHlqtNnueeqF3wgyE/Fcg7CujoM9LeCQYsrbUUvIEl04MWYFSBqL9+P1knHPwEmnWlpAjWJrPSShPyNqT1uUP3Uvel7goXUUPhenIM7bCw3HGmSEpTEvFhYfjVESFLvmF4MLDcSrALH2INw648HCcSsjswzKyuPBwnCqw4hSvVeHCw3EqwAzaTRceldBqttn75HOFlBUyu6bmCZg908y+IfNqWv0hZ7a+2pxmx82JBSY4pSkBQ+bV9HpCJtn8jnEh8ppqmzNFLHtufT2PUWJshYfjjDMGtH3Y4jhObizqPY8zLjwcpxLc2uI4Th+YQauIPRwqxIWH41SE9zwqotVssfeJPYWU1Y/DWLCsPiwaqY55BVpIynCMSy2rwL0GitzzJG+7iugxmJmbah3H6QODVoU7XReBCw/HqQDD3DHOcZw+8J6H4zj94j0Px3FyY2a0vefhOE5+zIctVdFqNnn2yT1VN8M5BGk1mwOXYebDFsdx+sGg7TNMHcfJj2ED7O40CpSxYxyS1kh6SNJ2SRsS4udIuj6Ov0vS8WW0y3GqYta3JXRkodd7NWyGLjwk1YGrgHOBVcA6Sau6kl0KPG1mJwKfAT417HY5TrVE1pbQ0YuM79VQKaPncRqw3cx2mNkB4DpgbVeatcCm+PyrwDtVpMOJ44wascI0dGQgy3s1VMrQeRwLPNpxvRM4PZQm3tt2L3Ak8ERnIknrgfXx5f7v3njW/UNpcTaW0NU+r/+Qqf8XBy3gub0Pfevvt/zykpQkU5K2dlxPm9l0x3WW92qojJXCNH540wCStprZ6qra4vUfuvV3vdR9YWZrimhLlZQxbHkMWNFxvTwOS0wjqQEcBjxZQtscZ1zJ8l4NlTKEx93ASkknSJoELga2dKXZAlwSn18IfNvGfaUUxxkuWd6roTL0YUusw7gM+BZQB64xs22SNgJbzWwLcDXwRUnbgaeIHkQvpnsnGSpe/6Fbf9X3HnyvymyD/APvOE4/lDJJzHGcgw8XHo7j9MXIC4+qp7ZnqP9ySQ9IulfSLZJeW2b9Hel+XZJJKsx8maVuSRfF979N0peLqjtL/ZKOk3SrpHvi539ewfVfI2m3pMT5RIr4XNy+eyWdWmT9I4+ZjexBpAj6CfA6YBL4IbCqK81vA5+Pzy8Gri+5/ncA8+Lz3yq7/jjdQuB24E5gdYn3vhK4Bzg8vj665Gc/DfxWfL4KeLjg398vA6cC9wfizwO+AQg4A7iryPpH/Rj1nkfVU9t71m9mt5rZC/HlnUT29qLIOgX540T+QPtKrvt9wFVm9jSAme0uuX4DFsXnhwE/K7B+zOx2IutfiLXAFyziTmCxpGVFtmGUGXXhkTQF99hQGjNrArNT28uqv5NLib5ERdGz/rirvMLM/qbAejPVDbweeL2kOyTdKanIWZNZ6v8o8B5JO4GbgA8UWH8W8v4+DirGanr6KCPpPcBq4O0l1lkDPg28t6w6u2gQDV3OIupx3S7pTWa2p6T61wHXmtmfSjqTaK7QSTbuC2WMCaPe86h6anumKcCSzgE+ApxvZvsLqjtL/QuBk4DbJD1MNO7eUpDSNMu97wS2mNmMmf0z8CMiYVIEWeq/FNgMYGb/AEwROcyVReVTxCulaqVL2kH0ZdsBnMDLSrM3dqV5P69UmG4uuf5TiBR7K6u4/670t1GcwjTLva8BNsXnS4i68EeWWP83gPfG579EpPNQwf+D4wkrTN/NKxWm3yv6NzDKR+UNyPDPO4/oi/YT4CNx2EairzxEX5u/ArYD3wNeV3L9/xd4HPhBfGwps/6utIUJj4z3LqJh0wPAfcDFJT/7VcAdsWD5AfCuguv/CrALmCHqZV0K/Cbwmx33f1XcvvuKfPbjcPj0dMdx+mLUdR6O44woLjwcx+kLFx6O4/SFCw/HcfrChYfjOH3hwsNxnL5w4eE4Tl+48DjEiNe/+JX4/BOS/qzqNjnjiTvGHXpcCWyUdDTR1PrzK26PM6b4DNNDEEnfARYAZ5nZs1W3xxlPfNhyiCHpTcAy4IALDmcQXHgcQsSrXH2JaAWs5wpevMc5xHDhcYggaR7wNeBDZvYg0dKFV1bbKmeccZ2H4zh94T0Px3H6woWH4zh94cLDcZy+cOHhOE5fuPBwHKcvXHg4jtMXLjwcx+mL/w/gqoslfoGymAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 360x216 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "u = sim_heat_equation_numpy(a, grid_dim, grid_samples, dt, steps)\n",
    "plt.figure(figsize=(5, 3))\n",
    "plot_heat_grid(u, grid_dim)\n",
    "plt.xlabel(\"$x$\")\n",
    "plt.ylabel(\"$y$\")\n",
    "plt.tight_layout()\n",
    "show_book('heat_sim')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Visualisierung Laplace-Operator"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Siehe Abbildung §\\ref{fig_laplace_sketch}§\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAUcAAADQCAYAAACUXaMkAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAx+UlEQVR4nO3dd1xX1f/A8df5sBFQcCsqbkEUFMS9NUeuNDWbaunXpmlp2vhp37JpmaXmN2dZWa7K1DQ1SVEpxRzkREXFiYQbZJ3fHxfIAcj4wGW8n48HDz7rnvv+wIc35957zvsorTVCCCFuZzE7ACGEKIwkOQohRAYkOQohRAYkOQohRAYkOQohRAYkOQohRAZszQ4gO8qVK6e9vLzMDkMIUcyEhYVd1FqXz+i5IpEcvby82Llzp9lhCCGKGaXUicyek8NqIYTIgCRHIYTIgCRHIYTIQJE45ygKp8TERKKiooiPjzc7FCGy5OjoiKenJ3Z2dtneJt+So1JqPtALuKC19k19zAP4HvACIoFBWuvY/IohS0kJEBMBV8/CjX8g7h+4EWPcTooDG3uw2IGNHVhswd4F3KpA6arg5ml8t3MyJfTCIioqCldXV7y8vFBKmR2OEBnSWhMTE0NUVBQ1a9bM9nb52XNcCMwAvrrlsQnARq31e0qpCan3X8nHGAwJN+Dkdji3F87vh/N/w8XDkJJ4xwsVOJUBWyfjueQESE769/adSpWHSo2gsj9UaQJV/KF0NSghiSI+Pl4Soyj0lFKULVuW6OjoHG2Xb8lRa71ZKeV1x8N9gQ6pt78EgsmP5Ki10Ss8sh4iNsCJrZCUeujnVhUqNoS6XY3vZaqDc1lw8jASo8Um4zYT4+HqGbh8Gq6chstREHsczu6BbZ9CSpLxOueyULMd1OkKdbqAa0Wrv73CRBKjKApy8zkt6HOOFbXWZ1NvnwOsnzlOhhK/5Ckcr50y7petCwHDjETlGQBO7rlr184RPGoZX3dKjDd6o2f/gqidcPQ3+PsH47lKjYxE6dPH6GFKMrGKmJgYOnfuDMC5c+ewsbGhfHljLO+ff/6Jvb29meHdJjg4GHt7e1q1amV2KAC4uLhw7dq1bL22X79+nDt3jtDQUKvHMXnyZFxcXHj55Zczfc3s2bNxdnbm8ccfz/D54OBgpk6dyqpVq6wen2kXZLTWWimVaaVdpdRIYCRA9erVs99waU+Oqep8m9iFUx6tGNG9E23qlstzvFmyczQSr2cANHsKUlLgfDhErIeIjbB1OoR8DOXqg99gaDQIylTL35iKubJly7J7924ge39k+S0pKQlb24z/nIKDg3FxcclRcsyqvYJy6dIlwsLCcHFx4dixY9SqlUHHIJ+NGjWqwPeZpqCH8pxXSlUGSP1+IbMXaq2/0FoHaq0D03oE2VLaE++xq2n10CscSy7Ho/P+YPjCHURcuJrn4LPNYoHKjaHtSzBsDYyLgF7TjF7rxv/CJ76wsBfs+d64MCSsIiwsjPbt2xMQEEC3bt04e9Y4SOnQoQNjxowhMDAQb29vduzYQf/+/albty6vv/46AJGRkTRo0IBHHnkEb29vHnzwQW7cuHHPdl988UUCAwOZPn06P//8M82bN6dJkyZ06dKF8+fPExkZyezZs5k2bRr+/v5s2bKFoUOHsmzZsvS4XVxcACOJtm3blj59+uDj43PX+3NxceG1117Dz8+PFi1acP78eYAM9wtw7do1hg0bRqNGjWjcuDHLly9Pbyujdu60YsUKevfuzUMPPcR3332X/nhm8Z89e5Z27drh7++Pr68vW7ZsAWDevHnUq1ePoKAgRowYwXPPPXfXvo4ePUr37t0JCAigbdu2HDx4EDD+8U2dOhWAiIgIunTpgp+fH02bNuXo0aPp7/PBBx9M//1ZbXUDrXW+fWFclQ6/5f6HwITU2xOAD7LTTkBAgM6N+MQk/b/fI7Tv/63VX247nqs2rC7mmNbB72s9vYnWk9y0/rCucf/qBbMjy7H9+/ffdn/Q7G13fX2V+nO/cTMpw+eX7DiptdY65trNu57LrkmTJukPPvhAt2zZUl+4YPwcv/vuOz1s2DCttdbt27fX48eP11pr/cknn+jKlSvrM2fO6Pj4eF21alV98eJFffz4cQ3okJAQrbXWw4YN0x9++KFOSEjIst2nn346PY5//vlHp6SkaK21njNnjh47dmx6fB9++GH665544gm9dOnS9PulSpXSWmu9adMm7ezsrI8dO5bh+wT0ypUrtdZajxs3Tr/11ltZ7nf8+PF69OjRt8WXVTt36tKli968ebM+dOiQ9vX1vWf8U6dO1W+//bbWWuukpCR95coVffr0aV2jRg0dExOjExISdJs2bfSzzz5718+lU6dO+vDhw1prrUNDQ3XHjh3vek1QUJBesWKF1lrruLg4ff36db1p0ybt5uamT506pZOTk3WLFi30li1bMnw/d35eU38WO3UmeSc/h/Isxrj4Uk4pFQVMAt4DliilngROAIPya/8ADrY2jGxXmwFNPXFzMsY3/fjXaS5cjeeJVl442GZy8SU/edSE9uOh7ctw7DcI/Rw2TYHNU6HRQGj5jHGhSOTIzZs3CQ8Pp2vXrgAkJydTuXLl9Of79OkDQKNGjWjYsGH6c7Vq1eLUqVOUKVOGatWq0bp1awAeffRRPv30U7p3755lu4MHD06/HRUVxeDBgzl79iwJCQk5GjaSJigoKNPt7O3t6dWrFwABAQGsX78+y/1u2LDhth6fu7t7lu3c6vz58xw5coQ2bdqglMLOzo7w8HB8fX0zjb1Zs2YMHz6cxMRE+vXrh7+/Pxs3bqR9+/Z4eHgAMHDgQA4fPnzbdteuXWPbtm0MHDgw/bGbN2/e9pqrV69y+vRpHnjgAcAYt3jrz8zT0xMAf39/IiMjadOmTaZxZld+Xq0ekslTnfNrn5kp6+KQfjsk4iLLwqJYFHqCiT286eFbyZwrrhaLcZGoTheIPgx/zIY9i2H31+DTDzpMgAreBR9XHnz/n5aZPudkb5Pl8x6l7LN8/l601jRs2JDt27dn+LyDg/EZsFgs6bfT7iclGSMN7vwcKKXu2W6pUqXSbz///POMHTuWPn36EBwczOTJkzPcxtbWlpSUFABSUlJISPj31Mqt7d3Jzs4uPUYbG5v0uLO733u1c6slS5YQGxubnmivXLnC4sWLmTJlSqbxt2vXjs2bN7N69WqGDh3K2LFjcXNzyzKWtDbKlCmTfg45p279fWb2fnKjxE0fnDrQj6+GB+FsZ8sz3+xi4OzthJ++bG5Q5etBr49hzN/Qbpwx/GhWS1g23Eic4p4cHByIjo5OT2KJiYn8/fffOWrj5MmT6dt/++23tGnThvr162e73cuXL1O1alUAvvzyy/THXV1duXr133PeXl5ehIWFAbBy5UoSE+8cb5szme23a9euzJw5M/1+bGz251ssXryYtWvXEhkZSWRkJGFhYem90MziP3HiBBUrVmTEiBE89dRT7Nq1i2bNmvH7778TGxtLUlLSbec907i5uVGzZk2WLl0KGP/o9uzZc9trXF1d8fT05McffwSMnmXaOeH8UuKSI0C7euVZM7ot7/ZvRGTMDaKv3rz3RgXB2QM6vQ4v7oM2Y+DQWpjVHFaMhEunzI6uULNYLCxbtoxXXnkFPz8//P392bZtW47aqF+/PjNnzsTb25vY2Fiefvpp7O3ts93u5MmTGThwIAEBAZQr9+8Iid69e/PDDz+kX5AZMWIEv//+O35+fmzfvj3L3mJ2ZLbf119/ndjYWHx9ffHz82PTpk3Zai8yMpITJ07QokWL9Mdq1qxJ6dKl+eOPPzKNPzg4GD8/P5o0acL333/P6NGjqVq1Kq+++ipBQUG0bt0aLy8vSpcufdc+v/nmG+bNm4efnx8NGzbkp59+uus1ixYt4tNPP6Vx48a0atWKc+fO5fRHlSNKF4F1qwMDA3V+1XO8kZCEk50NSik+3XiEhKQURnWojYtDIZh2fv2iMQzozy8ABa1fgNajwT5vf0zWcuDAAby9i9ahf2YiIyPp1asX4eHhZodS7Fy7dg0XFxeSkpJ44IEHGD58ePq5w4KU0edVKRWmtQ7M6PUlsud4K2d72/RzS6dj45ixKYKOU4P5fsdJklNM/sdRqhzc9xY8twMa9ITf34fPAo0hQKnnfIQo7CZPnpw+vKdmzZr069fP7JCypcT3HO/018lY3l59gLATsTSo5Mr7AxrjV61Mgez7nk6GwtoJcOYvqBpgjJ2s7GdaOMWp5yiKP+k55lGT6u4sG9WSmQ83JT4xGUc7Y7hPofgnUr0FPPUb9JttnIP8oiOsew0SrpsdmRDFTiE4sVb4KKW4v3FluvtWwsZiDHl4ackeXB1tGd2lHh6lTJy3a7GA/xCo3x02TIbtM2D/SqMXWbeLeXEJUcxIzzELaYkxJUXj4mjL13+cpP2Hm5iz+Rg3k5LNDc7JHXpPh2G/GHO7vxkAy540LuIIIfJMkmM2WCyK//b1Ze3otgTWcGfKmgN0/Xgze6MumR0a1GgFo0Kgw0Q4sBJmtYCDa8yOSogiT5JjDtSt6MqCYUF8NTyI8q4OVC5tVAI3vRdp62DMqBkZDC6V4Lsh8NOzEH/F3LjyUUxMDP7+/vj7+1OpUiWqVq2afv/WGSeFQXBwcI7HXFp7/6VLl8bf35/GjRvTpUsXLlzItObLPfXs2ZNLly5ZL8BCSpJjLrSrV57lT7eivKsDWmsenvMHY77fzZlLceYGVrEhjPjNqAa0+1uY3RoiQ8yNKZ+klSzbvXs3o0aNYsyYMen3zajlmNWUtdwkR2tNgUvTtm1bdu/ezd69e2nWrNltM2dyas2aNZQpU8Z6wRVSkhzzKDFZ06KWB6v3naXj1GCmrjvEtZvW/WDnxKGjx7no9wwMX2esfbOwF6z/P0jO2xS1oqA4lyxbuXJles+4fv366XOevby8uHjROM+8c+dOOnTokOXPSGvN1atX04tQXL9+neHDhxMUFESTJk3SZ6bcuHGDQYMG4ePjwwMPPEDz5s1JG06Xts/IyEi8vb0ZMWIEDRs25L777iMuLi79Z/PKK68QFBREvXr10suXFSVytTqP7G0tjOvWgCFB1flw3SFmbIrgux2n+Gp4ED5V7j3p3tq8vb3p2LEjGzduNM5FrnvVmGVzYhsMmAfuNfJnx79MgHP7rNtmpUbQ471svVRrzfPPP89PP/1E+fLl+f7773nttdeYP38+YFSi2blzJ9OnT6dv376EhYXh4eFB7dq1GTNmDACHDh1i3rx5tG7dmuHDhzNr1ixGjx6dZbsJCQnpSSM2NpbQ0FCUUsydO5cPPviAjz76iFGjRt1WjHfevHmZvo9du3YRHh5+V2WePn36pFcWGjRoEO3bt8/BDxK2bNmCv78/MTExlCpVinfeeQeAKVOm0KlTJ+bPn8+lS5cICgqiS5cufP7557i7u7N//37Cw8Px9/fPsN0jR46wePFi5syZw6BBg1i+fDmPPvooYPR+//zzT9asWcObb77Jhg0bchSz2SQ5WomnuzPTH2rC0FZezA05Tq3yxhS/2OsJuBfg0B+t9b9zTu1LGVe0a7aHn0fD/9pC35ng3bvA4ikoJaFkGcAHH3yAk5MTzz77bI7abdu2bfpSAu+//z7jx49n9uzZ/Prrr6xcuTK9oGx8fDwnT54kJCSE0aNHA+Dr60vjxo0zbLdmzZrpiTMgIIDIyMj05/r375/h40WFJEcra1LdnZkPG4cs8YnJ9PoshAaVXJnY05s6FVzMCcq3v7E64rJh8P2jEDQSur5lDAGylmz28PJLSShZtmHDBpYuXcrmzZszbCu764f36dOHAQMGAMbPbfny5dSvXz9b297pznJhaYfVtz5nzTJiBUnOOeYjpeDxljX48/g/dPtkM5N+Cuef6yZdSfWoCcN/hRbPGoUs5neDSyfNiSUfFPeSZSdOnODZZ59l6dKlODn9u176rW1lVA4sIyEhIdSuXRuAbt268dlnn6XPAPvrr78AaN26NUuWLAFg//797Ntn5VMmRYAkx3zkYGvDf9rXJnhcBx4Oqp4+iPxYdPZWfrM6W3vo/g48tBj+OQb/aw9Hs1fGqrAr7iXLFi5cSExMTHqF7Z49ewIwadIkRo8eTWBgIDY2mVe2Tzvn6Ofnx6JFi/joo48AeOONN0hMTKRx48Y0bNiQN954A4BnnnmG6OhofHx8eP3112nYsGGGpcaKtczWTyhMX7ldQ6awOXzuip6yen/6eh/Hoq+l37YWQPv4+Nz7hRcjtJ7RXOvJZbTe8rHWuYgjozU5iqrjx4/rhg0bmh1GoZGUlKTj4uK01lpHRERoLy8vffPmTZOjyptCs4aMuFvdiq682tOoChJ99Sa9Pt2CTxU3Xr/fp+Ar/5StDU9tgJXPG3O0T++CfrPAwbVg4xCF0o0bN+jYsSOJiYlorZk1a1ahWgu8IEhyNIm7sx2v9/Lho18P0XfmVh5oUpVx3epTpYzTvTe2FgcXeHA+eAbCr2/AnM4wZLGROEsYLy8vKXR7C1dXVwqqTGBhZco5R6XUGKXU30qpcKXUYqWUFS+bFg22NhaGBFVn08sdeKZDbVbvO0vnj37nwpXsXXG0GqWg5bPw+E9wPRrmdobjm++9nRDFXIEnR6VUVeAFIFBr7QvYAA8VdByFhaujHeO7N+C3l9rzas8GVHAz/k/siPynYCuR12xrTD10qQiLHoAdmQ9UvpUuDHUuhbiH3HxOzbpabQs4KaVsAWfgjElxFBqe7s481tILgKPR1xj8v+3c/+kWtkbcuwTZ9evXuXbtGteuGVfBb968mX7/1nFn9+RRE55cD7U7weqxsGYcJGc+Ps3R0ZGYmBhJkKJQ01oTExNz21rX2WHKMglKqdHAFCAO+FVr/UgGrxkJjASoXr16wIkTJwo2SBNprVm97yzv/XKQqNg4OjeokOkg8gMHDuDj44OtrbEWTtqYOTs7O7TW2NnZER0dnbMV7lKSjfnY22dArQ4w8EtwKnPXyxITE4mKisr24GMhzOLo6Iinpyd2dna3PZ7VMgkFnhyVUu7AcmAwcAlYCizTWn+d2TYFuYZMYRKfmMyCrZHM2hSBBrZP7ISr4+2/3MTEROrUqcPJk3cP6LaxsaFz586sW7cudwHsWgSrXoSydeGRJVCmeu7aEaKQKmxryHQBjmuto7XWicAKoJUJcRR6jnY2PN2hNpvGdWDaYH9cHY3e4Mo9Z9JrSNrZ2TF//nycnZ3v2t7e3p7//e9/uQ+g6WPw6Aq4cgbmdoEzu3PflhBFjBnJ8STQQinlrIzJrJ2BAybEUWSUc3Ggq09FAP46dYkXFv/FfdM2szb8LFprOnfuTPv27bG1/XdklrOzM2PHjsXLyytvO6/VHp5cBzb2sKAnHFqbt/aEKCLMOuf4JsZhdRLwF/CU1vpmZq8vqYfVmdl8OJopqw9w6PxVgrw8eL2XN2VSrtCgQYP0CzAVKlQgMjLytnm4eXL1PHw7CM7thR4fQNAI67QrhIkK22E1WutJWusGWmtfrfVjWSVGcbd29cqz+oU2vPNAI45dvMZTX+6kUhVPxo0bh5OTE6VKlWL27NnWS4wArhVh2Bqo2w3WvGzMqpGr1KIYM6XnmFPSc8zc1fhEjkZfx79aGa5ev4F3n6epkXSKkOCNd5XgsoqUZCM57pwPfg9Dn0/Bxu7e2wlRCGXVc5Tpg0Wcq6Md/qnzssOirmHbbBA3S9mzZOcpHgyolr68rNVYbOD+j42FvILfgRsXYeBCo7CuEMWIlCwrRjrUr8CKZ1pRvawzryzfR6/PQrI1iDzHlIIOr0CvTyBiA3zZB67HWH8/QphIkmMx07S6O8ufbsWMh5twNT6Rt1btJyW/piEGDoNBi+B8OMy/r1gVzxVCkmMxpJSiV+MqbBjbni8eC8RiUVy+kcg7aw4Qa+1K5N694LEfjaIV87pB9CHrti+ESSQ5FmOOdjZUL2sMDt929CJztxyj/YebmLP5WPogcquo0RKGroGUJJjf3agNKUQRJ8mxhOjRqDJrX2xH0xruTFlzgK4fb+aXfWett4NKvjB8rVEj8ss+cLzorVMsxK0kOZYg9Sq6snBYEF8OD8LRzsLyXVHW3UHZ2jB8HbhVga8HwKFfrNu+EAVIkmMJ1L5eeda80JYPH/QD4PjF67y0ZA9nLuWgvFlm3KrAsF+gYkP47hHYuyTvbQphAkmOJZStjQX3UsaaIHujLvHz3jN0+iiYj389xPWbeVxjuFRZeGIl1GgFK0ZC2MK8ByxEAZPkKOjrX5WNY9vT1acSn/4WQYepwSwPy+Mht4MrPLIU6naFn0fD9lnWCVaIAiLJUQBQzcOZz4Y0YfnTrfB0d2Lf6ct5b9TOCQZ/A969Yd1E2Pxh3tsUooBIchS3CajhzoqnWzGhRwMAQo/F8NSXOzgafS13Ddraw4MLofFg+O1t2PCmFKwQRYIkR3EXpRSOdjYAnL0cR+ixf+g2bTOTV/6du0HkNrbQbzYEDIWQj2HtBEmQotCTwhMiSw808aRt3fJMW3+Yr7ZHsmJXFBN6ePNw8xwumWCxGHOx7ZwhdBYkJ0DPj4zHhSiEJDmKeyrn4sCUBxrxRCsv3llzgLhEY3ZNSopGKbJfGk0p6PaOUVV86yeQnAi9P5UEKQolSY4i29IGkacVsli+K4qlO6N4vZc3jT3LZK8RpaDLZCNBbv7AmHLYd6ZRCk2IQkT+ZYscs6TWiLS3tXA0+hp9Zmxl7Pe7OXs5m4PIlYJOr0HH12DPYmMsZBbrYwthBuk5ilzr61+Vjg0qMGvTUeZvPc6a8LO81tObx1p6Za+B9uONHuSGSZCSCAPmSVVxUWiYkhyVUmWAuYAvoIHhWuvtZsQi8sbN0Y4JPRrwSPPqfLDuEBXdHAG4mZSMrcVy70rkbV40EuS6iaBTYMB8Y/iPECYz67B6OrBWa90A8EOWZi3y0gaRpy0hO3PTUXp/FsK27FQib/kMdH8fDvwMy4ZBkpVrTgqRCwWeHJVSpYF2wDwArXWC1vpSQcch8kfalWvvSq5cjkvk4bl/8NSXO+89iLzFKOjxIRxcBUufkAQpTGdGz7EmEA0sUEr9pZSaq5SS1ZmKmR6NKrPxpfa80r0Bocdi6DZtM1+Hnsh6o+YjoedUOLQGljwOSbJirzCPGcnRFmgKfK61bgJcBybc+SKl1Eil1E6l1M7o6OiCjlFYgaOdDU93qE3wuA4MalYtfZXEy3GJJCSlZLxR0Ai4/yM4/At8/5gkSGGaAl+3WilVCQjVWnul3m8LTNBa35/ZNrJudfEydsludp2IZUIPb7o1rJjxIPKd82HVGKjbDQYvAluHgg9UFHtZrVtd4D1HrfU54JRSqn7qQ52B/QUdhzBPH78q2NlYGPV1GIO/CGVfVAYVgAKHQ69pcGSdHGILU5h1tfp54Bul1F7AH3jHpDiECTrUr8Avo9vydj9fjl64Ru8ZIXy/I4NlXQOHw/0fw+G1sOQJSZCiQJkyzlFrvRvIsCsrSgZbGwuPtqhBH/8qzA4+Sof6FQCIir2Bu7M9pRxSP5rNngQ0rH7JSJCDvpJxkKJAyPRBYSo3RzvGd29ABVcHtNaMXbKHjlODWbLzFMmpc7hp9pRxFfvwLzLMRxQYSY6iUFBKoZTile4NqOruxPhle28fRB404t9hPkuHSoIU+U6SoyhU0iqRfzakSfog8pV7zhhPBo0wBoofWm3MpElONDdYUaxJ4QlR6Cil6O1Xha4+Ffnmj5N09TamJIafvkxV36G46xRY+wosf1KKVYh8I8lRFFqOdjY82aYmYBTWfX7xX8Rcu8kLne9jaNckbNe/BsoC/ecaSzEIYUVyWC2KBItFMfvRAPyqleHt1Qfosq0hBxu/An//AD9IPUhhfZIcRZFRv5Iri55szsJhzbCzsdD9Tz+O+o2D8OXw4yhISTY7RFGMyLGIKHI61K9Amzrl+CX8HDUb9YRyzrDxTeKSwGng/2TJBWEV9+w5phZ/mKeUelEp1VkpVb4gAhMiK7Y2Fnr7VcFiUdxsOZrPLUNwOrCU8NmPcz1ehvmIvMvOYXUfYClgD/wHiFRK3aP2lBAFx8HWhl7PfczqskPxvbCKje8PZsmOE/8OIhciF+6ZHLXWZ7TWa7XWH2itB2FM+5ub/6EJkX3VPJy5//npnPF7nj76N5J+Gk14VKzZYYkiLDuH1TVuva+1PgDUy7eIhMiDKv3eQrd5iYdtN+G39y3Qmu/+PMmxe1UiF+IO2bkgs1gpVR04DuwDLmEsjCVE4aMUqvMboJNh6yckpMC7u7pzPSGZR1vUYHTnuriXksIV4t7umRy11q2UUY20NtAI8AB653dgQuSaUtBlMugU7Ld9yvYmFt5OfoKvtkfyw1+neaFzXR5rUQN7WxnJJjKXrU+HNkRorX/QWs/TWkfld2BC5IlS0PW/0PI5nP+ayztO3/LLC21p7Fma9385yPkr8WZHKAo5Gecoii+l4L63jfWwQ2dR32LDouFvczzmBtU8nAGYuu4Q3X0r4Vu1tMnBisJGkqMo3pSCbu8YCXL7DFCKml3fAuDClXgW/3mSmcER9G/iybhu9alU2tHkgEVhIclRFH9KQff3jOmF2z4zilV0eZMKbo5sGteBWZuOMj/kOKv3nWFku9qMal8LZ3v50yjp5Iy0KBmUgp4fGlXFt06HDZNBa9wc7ZjQowEbX2pPF++KfPvHCZJk8LhAeo6iJFHKqCauNWz9xHisy2RQimoezsx4uCmx1xNwc7QjKTmFsUv28FBQNVrVLmdm1MIkpiVHpZQNsBM4rbXuZVYcooRJS5CkJkiloPMk4zukj4GMio0j7EQsK/ecoYt3RV7t2YBa5V3Mi1sUODMPq0cDB0zcvyipLBbo+REEDIOQabDxv0Zv8hZe5Uqx8aX2jO9en9BjMdw3bTNv/vw3NxKkbmRJYUrPUSnlCdwPTAHGmhGDKOEsFmNNbDSEpH6/pQcJRiXyZzrUYWBANT5ef5jtR2Owt5HT9CWFWYfVnwDjAVeT9i9EaoKcBiijB6lToMubtyVIgPKuDrzbvxE3k5KxtbFwJT6RYQt2MLJdLe7zqYi64/WieCjwf4NKqV7ABa112D1eNzK1luTO6OjoAopOlDhpPcjAJ42r2OvfuOsQO42DrVFE9/zleK7EJfKfRWE89EUo4acvF2TEooAonckHId92qNS7wGNAEuAIuAErtNaPZrZNYGCg3rlzZwFFKEokrWHNONgxB1o+Z8ysyaJHmJScwuIdp5i2/jCxNxLo38STd/r7pidQUTQopcK01oEZPVfgh9Va64nARAClVAfg5awSoxAFIm0cpLIYM2l0ijGzJpMEaWtj4bEWNejrX4WZmyI4Hn09PTEmp2hsLHKoXdTJOEch0igFPd43EmToLEhJgh4fZNmDdHO0Y2IPb9KOwE79c4OH54byQqe6DGjqiUWSZJFl6qU3rXWwjHEUhYpS0P1d49D6zy9g1RhIScnGZkYSjEtMpmwpB8Yt20vvGSFsO3oxvyMW+UTGJQhxp7RqPm3GQNgCWPl8tpd9rVfRlR+eacWnQ5pw6UYiD8/5g1GLwkiRKYlFjhxWC5GRtJkzNg7w+3uQkgh9Z4HNvf9klFL08avCfT4VmRdynCtxiemH1zcSkqSoRREhvyUhMqMUdJxoJMTf3obkBOg/B2zssrW5o50Nz3ask34/7MQ/DF+4UyqRFxHy2xHiXtqNM6qK//0DLB0KSTdz1UxpJzsae5bmrVX7uW/a76z7+xwFPZROZJ8kRyGyo/Vo48r1wVWweAgk3MhxE3UquPLV8CAWDGuGrY2F/ywKY+SiLOdCCBPJYbUQ2dX8P2DnbFyg+WYgPPwdOORsBqxSio71K9C2TjkW/3ky/VxkSorm4vWbVHCVSuSFhfQchciJpo/BgLlwcjt81Q/iYnPVjK2NhcdaevFIc2NZ+JV7ztD+g2CmrT8slX8KCUmOQuRUowdh8CI4txcW9oZreZ/7H1DDnU7eFZi+8QgdpwazLCxKhv+YTJKjELnR4H4Y8h3ERMCCHnA5b6sVV/NwZubDTVn+dEsql3bi5aV7ePbbXVYKVuRGgReeyA0pPCEKrRPb4dtB4OAGj/8E5erce5t70Frz896zuDrY0rFBBeISkjl3JZ6a5UpZIWBxq6wKT0jPUYi8qNEShq6CpHiY3w3O7slzk2mDyDs2qADA/K3H6frx7/z35/1cupGQ5/ZF9khyFCKvKvvB8HVg5wQLexm9SSsaGOjJwEBPFm47TvsPg5kfcpyEpHvP9xZ5I8lRCGsoVweGrwWXirDoATj8q9WaruDqyLv9G7NmdFsae5bmv6v288ryvVZrX2RMzjkKYU3XL8LX/eH839BnBvgPsWrzWmuCD0dTwdWBhlVKc+FqPBeu3MS3ammr7qekkHOOQhSUUuXgiVVQozX8OApCPsl02YXcSBtE3rCKkQxnbTpK7xkhvLRkD+cux1ttP0KSoxDW5+gGjywF3wGwYRKsezVbNSFzY0zXeoxsW4uf95yh41QZRG5NkhyFyA+2DtB/LjR/2qgqvuKpXBesyEppJzsm9vRm/dh2uMedZvrGI7z/y0Gr76ckkuQoRH6xWIyq4l0mQ/hyYz52fP6sVBiy9kd2TB+FTfB0RrStCcDh81cJPRaTL/srCSQ5CpGflDIqivf7HE5shXnd4NIpq+7i6tWrPP/88yQlJXE+fBs/frsAgFmbInjoi1BGfrWT4xevW3WfJYEkRyEKgv/D8OhyuHIa5naGM7ut1vTrr79OXFwcANevX2fixIlcvHiR9wY0Zly3+myNuMh902QQeU4VeHJUSlVTSm1SSu1XSv2tlBpd0DEIYYpaHYzB4jb2sKAnHF6X5yYPHjzIF198QXz8v1eqExISGDt2bHol8k3jOvBggDGI/KvtJ/K8z5KiwMc5KqUqA5W11ruUUq5AGNBPa70/s21knKMoVq6eM+Zjn9tnrJXd7KlcNaO1pk2bNoSGhpJyx9VwJycntmzZQkBAQPpjB85eobqHM6UcbNkacZHrN5Po6lMxfeXEkqhQjXPUWp/VWu9KvX0VOABULeg4hDCNayUYugbqdIXVL8Ga8ZCc8+E3q1atyjAxAsTFxfHkk0/e9ph3ZTdKORj1rRdsjWTkojCGzAkl/HT+XCQq6kytBK6U8gKaAH9k8NxIYCRA9erVCzYwIfKbgwsMWQy/vgGhMyHmCDy4AJzKZLsJHx+f2xLgnDlzqFixIn369AGgUaNGmW77+aNNWfznST7ZcITeM0IY0NSTl++rT6XSUok8jWnTB5VSLsDvwBSt9YqsXiuH1aJY2/UVrBoL7l5Gjchclj1TSvHuu+8yYcKEbG9zOS6RWZsiWLA1krcf8GVQYLVc7buoKlSH1QBKKTtgOfDNvRKjEMVe08eNWpBx/8DcTnB0U4HtOm0Q+W8vt2dAU08AloVFSSVyzLlarYB5wAGt9ccFvX8hCiWv1jDiN3CtAl8PgO2zrDon+1483Z2xSV3s6+c9Z3h56R76zAxh+9GSO4jcjJ5ja+AxoJNSanfqV08T4hCicHH3gid/hfo9YN1EWP4UJBT84O0FQ5sx/SF//rmWwJA5xiDyEzElbxB5gV+Q0VqHACV37IAQWXF0g0GLYOs02PgWXNgPg7+GsrULLASLRdHXvyrdGlZiXshxZm2K4MyleGqULVnLNMgMGSEKG4sF2r5kzKi5eha+6AiH1hZ4GGmDyLdN7EzL2mUBmLb+MAu2HicxufhXIpfkKERhVaczjPwd3GvA4sGw8b+5Gg+ZV6Wd7ABISdHsjbrEmz/vp9u0zazff56iUCw7tyQ5ClGYudcwzkM2eQy2fARf9srzMrC5ZbEo5g9txoKhzbBYFCO+2smQOaFEXLhmSjz5TZKjEIWdnRP0nWHUhzy3D2a3gUO/mBKKUoqODSrwy+i2vNW3IccvXk+/yl3cepGSHIUoKhoPhP9shtKesPghWDsRksypsmNnY+Gxll6EvNIpfT3t0d/t5pMNxacSuSRHIYqSsrXhyQ0QNNKoMD63E5w3ara0bt2aNm3aFGg4djZGCrmZlExSSgqfbDhCp6m/s7wYDCKX1QeFKKoOroaVL8DNK9D5/6DFM2CxMTWkHZH/8Paq/eyJuoxvVTemP9SE2uVdTI0pK4Vu+qAQwgoa3A/PhBrVfX59Hb7sDbHm1mts5uXBD8+05pPB/qSkQLlSDgAkFcGhP9JzFKKo0xp2fwu/vAJo6DYFmj5hLNFgalgapRTJKZr+s7YS6OXBC53qUtrZztS4biU9RyGKM6WgySPwzDao0gR+Hg0Le8HFIyaHZSTn+MRkvCu7MX/rcdpP3VRkBpFLz1GI4iQlBXZ/bRxmJ8ZBu3HQ+kWwtTc7MvafucKUNfvZGhFDrXKlWDCsmelTEqXnKERJYbEYJdCe3QHevWHTFGNc5IntZkeGTxU3vn6yOfOHBlKrvAtVyjgBcO1m4Rz6I8lRiOLItSI8OB8eXmr0IBd0h2VPWn1Z2JxSStGpQUXmPhGInY2FazeT6PxRMOOW7uH8lfh7N1CAJDkKUZzVuw+eDYV24+HgKpgRCL9NMaUUWmb6+lflp91n6Dg1mOkbjhSaQeRyzlGIkuLSKdgwCcKXg2tl6DwJGg8yfWwkwMmYG7y39gBr9p2jkpsjK59rTQW3/F/PRs45CiGgTDXjUHv4r0Zy/HEUzGoJ4SuMCzkmql7WmVmPBLB0VEt6NKpEeVdjfOSZS3GmxSTJUYiSpnpzeGqjsdqhUrBsGMxuDftXmp4km3l5MKl3Q5RSnLkUR8epwYz8aifHLxb8aQBJjkKURBYL+PaHp7fBgHmQnAhLHoMv2sHepaYVtLiVu7M9z3WsQ0jERe6b9jtvrdrP5RuJBbZ/OecohICUZNi3DDZ/aKyh7VoZmj0FAcOgVFlTQ7twJZ6P1x/m+52n8HC2Z9O4Drg5WmeWTVbnHE1Jjkqp7sB0wAaYq7V+L6vXS3IUooCkpMDRjUbFn6O/ga2jcdEmYChUaWrqlMT9Z67wx/EYhrWuCcC+1OIWKg8xFarkqJSyAQ4DXYEoYAcwRGu9P7NtJDkKYYILB+GP2bDnO0iKg7J1wW8wNBpkVCg30b6oy/SeEULLWmV57X5vfKuWzlU7hS05tgQma627pd6fCKC1fjezbSQ5CmGiuEuw/yfY+z2c2Go8Vr0VNOwHdboU6MqIaRKTU1j850mmrT/MpbhEHvIvz8tBDpSt6ZejdrJKjgW+NCtQFbh1mH4U0NyEOIQQ2eFUBgKeML5iT8C+JbB3Cfwy3njeo5aRJOt0Ba/WYJ//86XtbCw8Xi+ZgckRnN7xM1X378TmhAeMO2i1Q38zkmO2KKVGAiMBqlevbnI0QgjAOJxuN874ijkKERshYj3sWgR/fgHKAuUbQGV/o0JQFX+o2DDvCfN6DFz426h6fj4cIkMg9jhOQB2P2sQ1fRx7n+5G+TYrJUc5rBZC5F1ivHHIfTIUzvwFZ3fD9eh/n3fygNJVwc3T+O5SCWwdwMYebGzBYgcWW6Oq+Y2Y1K9/jO8Xj8C1c/+25VwWqgZC3a7G8rUetXIddmE7rN4B1FVK1QROAw8BD5sQhxDCWuwcjURVp7NxX2u4csZIlNEH4PJpuHIaLp2Ek9sh/lIWjSlwcjeSoLMH1O5o9D4r+EBFX3CpUCBXzQs8OWqtk5RSzwHrMIbyzNda/13QcQgh8pFSRg+xdFXw7nX380k3ja+UJGMAekqicdvBDRxLF4r53qacc9RarwHWmLFvIUQhYOtgfBViMn1QCCEyIMlRCCEyIMlRCCEyIMlRCCEyIMlRCCEyUCRKlimlooETOdysHHAxH8IpaMXlfYC8l8KouLwPyN17qaG1Lp/RE0UiOeaGUmpnZiPfi5Li8j5A3kthVFzeB1j/vchhtRBCZECSoxBCZKA4J8cvzA7ASorL+wB5L4VRcXkfYOX3UmzPOQohRF4U556jEELkWrFLjkqp7kqpQ0qpCKXUBLPjyS2lVDWl1Cal1H6l1N9KqdFmx5QXSikbpdRfSqlVZseSF0qpMkqpZUqpg0qpA6n1SYskpdSY1M9WuFJqsVLK0eyYskspNV8pdUEpFX7LYx5KqfVKqSOp393zso9ilRxTF++aCfQAfIAhSikfc6PKtSTgJa21D9ACeLYIvxeA0cABs4OwgunAWq11A8CPIvqelFJVgReAQK21L0b5wIfMjSpHFgLd73hsArBRa10X2Jh6P9eKVXIEgoAIrfUxrXUC8B3Q1+SYckVrfVZrvSv19lWMP8Kq5kaVO0opT+B+YK7ZseSFUqo00A6YB6C1TtBaXzI1qLyxBZyUUraAM3DG5HiyTWu9Gfjnjof7Al+m3v4S6JeXfRS35JjR4l1FMqHcSinlBTQB/jA5lNz6BBgPpJgcR17VBKKBBamnCOYqpfJ/Nal8oLU+DUwFTgJngcta61/NjSrPKmqtz6bePgdUzEtjxS05FjtKKRdgOfCi1vqK2fHklFKqF3BBax1mdixWYAs0BT7XWjcBrpPHQzezpJ6P64uR8KsApZRSj5oblfVoYxhOnobiFLfkeBqodst9z9THiiSllB1GYvxGa73C7HhyqTXQRykViXGao5NS6mtzQ8q1KCBKa53Wg1+GkSyLoi7Aca11tNY6EVgBtDI5prw6r5SqDJD6/UJeGituyTF98S6llD3GCeaVJseUK0ophXFu64DW+mOz48ktrfVErbWn1toL4/fxm9a6SPZQtNbngFNKqfqpD3UG9psYUl6cBFoopZxTP2udKaIXl26xEngi9fYTwE95aazQrludG8Vs8a7WwGPAPqXU7tTHXk1df0eY53ngm9R/vseAYSbHkyta6z+UUsuAXRgjI/6iCM2WUUotBjoA5ZRSUcAk4D1giVLqSYwqXoPytA+ZISOEEHcrbofVQghhFZIchRAiA5IchRAiA5IchRAiA5IchRAiA5IchRAiA5IchRAiA5IcRbGQWvuya+rtt5VSn5kdkyjaitUMGVGiTQL+q5SqgFHBqI/J8YgiTmbIiGJDKfU74AJ0SK2BKUSuyWG1KBaUUo2AykCCJEZhDZIcRZGXWp7qG4z6hNeUUneWzxcixyQ5iiJNKeWMUYvwJa31AeAtjPOPQuSJnHMUQogMSM9RCCEyIMlRCCEyIMlRCCEyIMlRCCEyIMlRCCEyIMlRCCEyIMlRCCEyIMlRCCEy8P/LAuh3s8QIygAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 360x216 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "from scipy.interpolate import interp1d\n",
    "xx = np.linspace(0, 10)\n",
    "straight_fun = interp1d([0, 10], [10, 0])\n",
    "yy_straight = straight_fun(xx)\n",
    "bendy_points = np.array([[0, 10], [2.5, 9], [5, 5], [7.5, 1], [10, 0]])\n",
    "bendy_fun = interp1d(bendy_points[:,0], bendy_points[:,1], 'cubic')\n",
    "yy_bendy = bendy_fun(xx)\n",
    "\n",
    "plt.figure(figsize=(5,3))\n",
    "plt.plot(xx, yy_straight, '--', label='Temperatur nach Ausgleich')\n",
    "plt.plot(xx, yy_bendy, label='Temperatur zu Beginn')\n",
    "plt.arrow(2.5, 8.8, 0, -0.8, head_width=0.3, fc=\"k\")\n",
    "plt.arrow(7.5, 1.2, 0, 0.8, head_width=0.3, fc=\"k\")\n",
    "plt.legend()\n",
    "plt.ylabel('$u$')\n",
    "plt.xlabel('$x$')\n",
    "#plt.tick_params(axis='x', which='both', bottom=False, top=False, labelleft=False, labelbottom=False)\n",
    "#plt.tick_params(axis='y', which='both', bottom=False, top=False, left=False, labelleft=False, labelbottom=False)\n",
    "show_book('laplace_sketch')\n",
    "#plt.plot(bendy_points[:,0], bendy_points[:,1], 'o')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Performance messen"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin perf_counter`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "31.8312628"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import time\n",
    "#!time.perf_counter()  # Zeit in s seit undefiniertem Punkt\n",
    "print(\"31.8312628\")#!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin perf_counter_sim`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Simulation braucht 4.069655700000002 s (Wanduhr-Zeit)\n"
     ]
    }
   ],
   "source": [
    "start_time = time.perf_counter()  # Startzeit merken\n",
    "sim_heat_equation_python_loops(a, grid_dim, grid_samples, dt, steps)\n",
    "time.sleep(2)  # \"Inaktive\" Zeit wird mit gemessen\n",
    "sim_time = time.perf_counter() - start_time  # Aktuelle Zeit - Startzeit\n",
    "#!print(f\"Simulation braucht {sim_time} s (Wanduhr-Zeit)\")\n",
    "print(\"Simulation braucht 4.069655700000002 s (Wanduhr-Zeit)\")#!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin timedelta`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'2:19:42.250000'"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import datetime\n",
    "long_time_seconds = 2*60*60 + 19*60 + 42 + 0.25 # 2h 19m 42s 250ms\n",
    "# Menschenlesbare Ausgabe in einzelnen Einheiten\n",
    "str(datetime.timedelta(seconds=long_time_seconds)) "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin process_time_sim`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Simulation braucht 1.84375 s (Prozessorzeit)\n"
     ]
    }
   ],
   "source": [
    "start_time = time.process_time()\n",
    "sim_heat_equation_python_loops(a, grid_dim, grid_samples, dt, steps)\n",
    "time.sleep(2)  # Inaktive Zeit wird nicht gemessen\n",
    "sim_time = time.process_time() - start_time\n",
    "#!print(f\"Simulation braucht {sim_time} s (Prozessorzeit)\")\n",
    "print(\"Simulation braucht 1.84375 s (Prozessorzeit)\")#!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## timeit"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin timeit_sim`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Simulation braucht 2.0994281249999993 s (Wanduhr-Zeit)\n"
     ]
    }
   ],
   "source": [
    "import timeit\n",
    "nr_runs = 4  # Wie oft soll die Funktion ausgeführt werden?\n",
    "# timeit.timeit nimmt eine ausführbare Funktion ohne Argumente\n",
    "#!time_per_run = timeit.timeit(lambda: sim_heat_equation_python_loops(a, grid_dim, grid_samples, dt, steps), number=nr_runs) / nr_runs\n",
    "#!print(f\"Simulation braucht {time_per_run} s (Wanduhr-Zeit)\")\n",
    "print(\"Simulation braucht 2.0994281249999993 s (Wanduhr-Zeit)\")#!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin timeit_magic`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.97 s ± 78.2 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
     ]
    }
   ],
   "source": [
    "#!%timeit sim_heat_equation_python_loops(a, grid_dim, grid_samples, dt, steps)\n",
    "print(\"1.97 s ± 78.2 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\")#!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin timeit_repeat`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([7.5460774, 7.0870385, 7.7462764])"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "run_times = np.array(timeit.repeat(lambda: sim_heat_equation_python_loops(a, grid_dim, grid_samples, dt, steps), number=nr_runs, repeat=3))\n",
    "run_times"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Zeit pro Aufruf: 1.8649493583333328 s ± 68.99028515097449 ms\n"
     ]
    }
   ],
   "source": [
    "print(f\"Zeit pro Aufruf: {np.mean(run_times/nr_runs)} s ± {np.std(run_times/nr_runs)*1000} ms\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# C-Code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "#!C:\\dev\\mingw-w64\\i686-8.1.0-posix-dwarf-rt_v6-rev0\\mingw32\\bin\\gcc.exe -c -fpic -lm heat_sim.c\n",
    "#!C:\\dev\\mingw-w64\\i686-8.1.0-posix-dwarf-rt_v6-rev0\\mingw32\\bin\\gcc.exe -shared -o libheatsim.so heat_sim.o"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin ct_loadlib`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "import ctypes\n",
    "lib = ctypes.cdll.LoadLibrary(r\"C:\\dev\\pfifd\\code\\kap_perf\\libheatsim.so\")#!\n",
    "#!lib = ctypes.cdll.LoadLibrary(\"libheatsim.so\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin ct_def_fun`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "# Hilfstyp für NumPy-Arrays und ctypes\n",
    "from numpy.ctypeslib import ndpointer\n",
    "c_heat_sim = lib.heatSim\n",
    "c_heat_sim.argtypes = [ndpointer(ctypes.c_float, ndim=2), ctypes.c_float,\n",
    "                       ctypes.c_uint, ctypes.c_uint,\n",
    "                       ctypes.c_float, ctypes.c_float,\n",
    "                       ctypes.c_float, ctypes.c_uint]\n",
    "c_heat_sim.restype = ctypes.c_uint"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin ct_wrapper`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "def sim_heat_equation_c(a, grid_dim, grid_samples, dt, steps):\n",
    "    dx, dy = grid_dim[0] / grid_samples[0], grid_dim[1] / grid_samples[1]\n",
    "    u = np.zeros(grid_samples, dtype=np.float32)  # Datentyp explizit auf C float setzen\n",
    "    # c_heat_sim als Python-Funktion aufrufen. Rückgabewert ist die Anzahl\n",
    "    # an Schleifendurchläufen\n",
    "    loop_counter = c_heat_sim(u, a, grid_samples[0], grid_samples[1], dx, dy, dt, steps)\n",
    "    return u"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin timeit_c`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "17.2 ms ± 377 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)\n"
     ]
    }
   ],
   "source": [
    "#!%timeit sim_heat_equation_c(a, grid_dim, grid_samples, dt, steps)\n",
    "print(\"17.2 ms ± 377 µs per loop (mean ± std. dev. of 7 runs, 10 loops each)\") #!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "libHandle = lib._handle\n",
    "del lib\n",
    "ctypes.windll.kernel32.FreeLibrary(libHandle)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Profiling"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin profiling_beispiel --input_prefix=\"\"`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "import time\n",
    "def schnarchnase():\n",
    "    time.sleep(0.5)\n",
    "\n",
    "def ruckizucki():\n",
    "    time.sleep(0.001)\n",
    "\n",
    "def profiling_beispiel(n):\n",
    "    schnarchnase()\n",
    "    for i in range(n):\n",
    "        ruckizucki()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin profile_fast`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "import cProfile"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "         26 function calls in 0.665 seconds\n",
      "\n",
      "   Ordered by: standard name\n",
      "\n",
      "   ncalls  tottime  percall  cumtime  percall filename:lineno(function)\n",
      "        1    0.000    0.000    0.508    0.508 <ipython-input-30-97c425bea2eb>:2(schnarchnase)\n",
      "       10    0.000    0.000    0.156    0.016 <ipython-input-30-97c425bea2eb>:5(ruckizucki)\n",
      "        1    0.001    0.001    0.665    0.665 <ipython-input-30-97c425bea2eb>:8(profiling_beispiel)\n",
      "        1    0.000    0.000    0.665    0.665 <string>:1(<module>)\n",
      "        1    0.000    0.000    0.665    0.665 {built-in method builtins.exec}\n",
      "       11    0.664    0.060    0.664    0.060 {built-in method time.sleep}\n",
      "        1    0.000    0.000    0.000    0.000 {method 'disable' of '_lsprof.Profiler' objects}\n",
      "\n",
      "\n"
     ]
    }
   ],
   "source": [
    "cProfile.run(\"profiling_beispiel(10)\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`begin profile_slow`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "         1006 function calls in 8.282 seconds\n",
      "\n",
      "   Ordered by: standard name\n",
      "\n",
      "   ncalls  tottime  percall  cumtime  percall filename:lineno(function)\n",
      "        1    0.000    0.000    0.508    0.508 <ipython-input-30-97c425bea2eb>:2(schnarchnase)\n",
      "      500    0.008    0.000    7.763    0.016 <ipython-input-30-97c425bea2eb>:5(ruckizucki)\n",
      "        1    0.010    0.010    8.281    8.281 <ipython-input-30-97c425bea2eb>:8(profiling_beispiel)\n",
      "        1    0.000    0.000    8.281    8.281 <string>:1(<module>)\n",
      "        1    0.000    0.000    8.282    8.282 {built-in method builtins.exec}\n",
      "      501    8.263    0.016    8.263    0.016 {built-in method time.sleep}\n",
      "        1    0.000    0.000    0.000    0.000 {method 'disable' of '_lsprof.Profiler' objects}\n",
      "\n",
      "\n"
     ]
    }
   ],
   "source": [
    "cProfile.run(\"profiling_beispiel(500)\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`end`"
   ]
  },
  {
   "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}