sciann_poisson_L_shape.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "sciann_poisson_L_shape.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "authorship_tag": "ABX9TyNV9nlhfe3MJC+ZEaoB+4GS",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    },
    "accelerator": "GPU"
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/swampxx/sciann-examples/blob/main/sciann_poisson_L_shape.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "FuJ_7U_8dAhh"
      },
      "source": [
        "# Install sciann\n",
        "%%capture sciann_installation\n",
        "!pip install sciann"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "n1qa4cyNUhKG",
        "outputId": "dd6c46cd-c0d8-41ee-c17a-0d934f6a7e17"
      },
      "source": [
        "import numpy as np\n",
        "from matplotlib import cm\n",
        "import matplotlib.pyplot as plt\n",
        "\n",
        "import sciann as sn\n",
        "\n",
        "%matplotlib inline"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "---------------------- SCIANN 0.6.3.1 ---------------------- \n",
            "For details, check out our review paper and the documentation at: \n",
            " +  \"https://www.sciencedirect.com/science/article/pii/S0045782520307374\", \n",
            " +  \"https://arxiv.org/abs/2005.08803\", \n",
            " +  \"https://www.sciann.com\". \n",
            "\n",
            " Need support or would like to contribute, please join sciann`s slack group: \n",
            " +  \"https://join.slack.com/t/sciann/shared_invite/zt-ne1f5jlx-k_dY8RGo3ZreDXwz0f~CeA\" \n",
            " \n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "JKNvJw5KDVC4"
      },
      "source": [
        "In this notebook, we solve the Poisson equation over an L-shaped domain $\\Omega = [-1, 1]^2 \\backslash [0, 1]^2$:\n",
        "\n",
        "$-\\Delta u(x,y) = 1, (x,y) \\in \\Omega$\n",
        "\n",
        "The boundary conditions:\n",
        "\n",
        "$u(x,y) = 0, (x,y) \\in \\partial \\Omega$\n",
        "\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "ibIeNKDQUmvb"
      },
      "source": [
        "# Define the L-shape domain\n",
        "\n",
        "x_data_1, y_data_1 = np.meshgrid(\n",
        "    np.linspace(-1, 1, 100), \n",
        "    np.linspace(0, -1, 50)\n",
        ")\n",
        "\n",
        "x_data_2, y_data_2 = np.meshgrid(\n",
        "    np.linspace(-1, 0, 50), \n",
        "    np.linspace(0,  1, 50)\n",
        ")"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "WF7RZQs2XCsf"
      },
      "source": [
        "# Prepare training data\n",
        "\n",
        "x_data = np.concatenate([x_data_1.flatten(), x_data_2.flatten()])\n",
        "y_data = np.concatenate([y_data_1.flatten(), y_data_2.flatten()])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "XKpE2_N8Xjtl"
      },
      "source": [
        "# Define Input variables x, y\n",
        "\n",
        "x = sn.Variable('x')\n",
        "y = sn.Variable('y')\n",
        "\n",
        "# Define the underlying network\n",
        "u = sn.Functional('u', [x, y], 8*[20], 'tanh')"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "j9Wh7gfYccR9"
      },
      "source": [
        "# The governing equation\n",
        "L1 = -sn.diff(u, x, order=2)- sn.diff(u, y, order=2) -1\n",
        "\n",
        "# Impose boundary conditions\n",
        "TOL = 0.001\n",
        "\n",
        "C1 = (1 - sn.sign(x+1 - TOL))*u\n",
        "C2 = (1 + sn.sign(x-1 + TOL))*u\n",
        "C3 = (1 + sn.sign(y-1 + TOL))*u\n",
        "C4 = (1 - sn.sign(y+1 - TOL))*u\n",
        "C5 = (1 + sn.sign(x + TOL))*(1 + sn.sign(y + TOL))*u"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "MWXfCyGTvnKf"
      },
      "source": [
        "# Construct SciANN Model\n",
        "model = sn.SciModel([x, y], [L1, C1, C2, C3, C4, C5])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "XxedXkImwhtd"
      },
      "source": [
        "# Training\n",
        "history = model.train([x_data, y_data], 6*['zeros'], learning_rate=0.002, epochs=5000, batch_size=256, verbose=2)"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "rTukSGBQwtey"
      },
      "source": [
        "# Genereate predictions\n",
        "u_pred = u.eval(model, [x_data, y_data])"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 295
        },
        "id": "y-4OGzULWDlM",
        "outputId": "7bc52898-36b2-4535-8827-2c9593a1b3ae"
      },
      "source": [
        "# Plotting\n",
        "\n",
        "# fill empty places with 0s\n",
        "empty_points_x, empty_points_y = np.meshgrid(\n",
        "    np.linspace(0, 1, 50), \n",
        "    np.linspace(0, 1, 50)\n",
        ")\n",
        "\n",
        "empty_points_x = empty_points_x.flatten()\n",
        "empty_points_y = empty_points_y.flatten()\n",
        "\n",
        "empty_u = np.zeros(50*50)\n",
        "\n",
        "# prepare pre\n",
        "plot_u_pred = u_pred.copy()\n",
        "plot_u_pred = np.concatenate([np.where(plot_u_pred>=0, plot_u_pred, 0), empty_u])\n",
        "\n",
        "plot_x = np.concatenate([x_data, empty_points_x])\n",
        "plot_y = np.concatenate([y_data, empty_points_y])\n",
        "\n",
        "plt.title('PINN Solution')\n",
        "plt.xlabel('x')\n",
        "plt.ylabel('y')\n",
        "plt.tricontourf(plot_x, plot_y, plot_u_pred)\n",
        "plt.colorbar()\n",
        "plt.show()"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYYAAAEWCAYAAABi5jCmAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3de7xdZX3n8c+XhCQqGgLhmkASMDYGcQBDaHVEyjU6I6E1aOiowYFJvdBOS6mEFx3AqBVoR9qOMJoR5GJroLFo+iIQw622VWwCIhiuJ+F2IuGWEFQgIeQ3f6xn4zo7+7L2Xte99+/9ep1X9l7XZ69z8nz38zzrIjPDOeecq9ml7AI455yrFg8G55xzI3gwOOecG8GDwTnn3AgeDM4550bwYHDOOTeCB4NzLUh6XNLxXa57oKRfSRqVdbmcy5MHg8tdqFxfCZXkM5KulrRbmHenpDPD62MkmaQr6tb/N0mnh9enh2U+X7fMsKRjmux/sqTvSnpe0hZJP69tL4fP+UaImNmTZrabmb2e9b6cy5MHgyvKh81sN+AIYBbwF02W+zXwCUlTW2xrE/B5SW9NuO/rgKeAKcCewCeAZxKu69zA8WBwhTKzDcDNwLuaLPIicDVwYYvNPAj8GDg74W6PBK42s1+b2XYz+6mZ3VybKelkSWslvRhaMO9stJHQ0vlS7P0xkobD6+uAA4F/Di2jz0uaGlo3o8My+0taLmmTpCFJ/yO2rYsk3SDpWkm/DOWZlfDzOZcpDwZXKEkHAB8CftpisS8DH5H0Wy2W+V/An0jaI8Fu7wIulzRf0oF15XkH8B3gT4C9gBVElfuYBNt9g5l9AniS0DIys0sbLLYUGAb2B+YBfynp2Nj8k8MyuwPLga91UgbnsuLB4IryPUkvAv8G/Avwl80WNLONwNeBxS2WuRdYBZybYN+nAv9KFCaPSbpX0pFh3seAm8xslZm9Bvw18CbgvQm2m1gIxPcB55rZq6H83wQ+GVvs38xsRRiTuA74T1mWwbmkPBhcUU4xs93NbIqZfdbMXmmz/CXASZJaVY4XAJ+RtE+rDZnZZjNbZGaHAPsA9xIFlYi+vT8RW3YH0XjEpASfqRP7A5vM7JexaU/U7Wdj7PXLwLhaN5RzRfJgcJVkZi8AfwN8scUyDwH/BJzfwXafJ2oV7A/sAfyCaFAagBAWBwAbGqz+a+DNsff71m++xa5/AexRN2B+YJP9OFcqDwZXZV8l6tJpOBgcfAH4FFG/fEOSLpH0LkmjQ8X8GWAohM8NwH+RdJykXYE/A7YCP2qwqXuBD0naQ9K+ROMScc8ABzUqg5k9Fbb5FUnjJL0bOAP4dovP5lwpPBhcZZnZS8ClRN/smy3zGFF//FtabOrNwI1EZzytJ2ohnBzWfxj4OPB/gOeBDxMNIG9rsJ3rgJ8BjwM/AK6vm/8V4C/C2U3nNFj/NGAqUevhRuBCM7u1RbmdK4X8QT3OOefivMXgnHNuhFKDQdJVkp6V9PMm8yXp78LFQPdJOiI2b4GkR8PPguJK7ZxzxZM0R9LDoT5c1GD+0ZLukbRd0ry6eQdK+oGkByU90ObOAqW3GK4G5rSY/0FgevhZCPxfgHBR04XAUcBs4EJJE3ItqXPOlSTciPFyojpxJnCapJl1iz0JnA78Q4NNXAv8lZm9k6jOfLbV/koNBjP7IdF9b5qZC1xrkbuA3SXtB5wErDKzTWa2mehCp1YB45xzvWw20Zl068OJEUuJ6sc3mNnjZnYfsCM+PQTIaDNbFZb7lZm93GpnVb94ZhLRxUY1w2Fas+k7kbSQqLXBKEa/5y2jvWHhnGvvpe3PPW9me3W7/vuPGWebN+1ovyCw9v7X1gKvxiYtMbMlsfeN6ryjEhblHcCLkv4JmAbcCixqddffqgdDauHgLgEYv+ve9t6Jp5ZcIudcL7hl4xVPtF+quc2bdvDdmyYmWnbGgU+/amZ53TRxNPB+4HCi7qbribqcrmy2QtljDO1sILoKtWZymNZsunPO9aM0dd4wcG/ohtoOfI/o9vdNVb3FsBw4S9JSombTFjN7WtJKojtT1vqFTgTOK6uQWds2Y3LZRUhtzEPDZRfBuX6yGpguaRpRIMwH/qCDdXeXtJeZPQccC6xptUKpwSDpO8AxwMRwX/sLgV0BzOzrRLdA/hAwRHRTsU+FeZskfZHoAwMsNrNWg9iV0g8VfztJPqOHh3PJmNl2SWcBK4FRwFVmtlbSYmCNmS0Pdwy+EZgAfFjSF8zsEDN7PVyJf1u4F9jdwP9rtb9Sg8HMTmsz34DPNZl3FXBVHuXK2iAEQTeaHRcPDOd2ZmYriL4sx6ddEHu9mqiLqdG6q4B3J91X1buSepIHQTr1x8+DwrlieTBkxMMgP/Fj6yHhXP48GFLyQChW7Xh7QDiXHw+GLnkglMsDwrn8eDB0qCqBsOXgsaXuf/y6raXuv8YDwrnseTB0oOhQKLvyb6Vd2YoODg8I57LjwZBQ3qFQ5RDoRqPPU0RYbJsx2cPBuZQ8GBLIKxT6LQzaiX/ePEPCw8G5dDwY2sg6FAYtDJrJOyQ8HJzrngdDC1mGggdCc3mFhIeDc93xYGgiq1DwQOhM7XhV5awn5wZR1W+7XYosQmHLwWM9FFLI6thV5fRi53qJtxhyUOVA+OUUjXj/1iespJK0560H58rhwVAnzTfMsgKhvrLPc90ygmTLwWNThYOPNTjXGQ+GjBQVCmlCIK/9FxEWacPBOZecB0NMt62FPEOh7CBIoqjuqTTh4K0G55LzYEgpj1DohTBoJV7+rEPCWw7O5c/PSgq6aS1kHQq/nKKeD4V6eXymKg/uO9cPSg0GSXMkPSxpSNKiBvMvk3Rv+HlE0ouxea/H5i0vtuTZ6sdAqJf1Z/RwcIMmQX15tKR7JG2XNK/B/LdJGpb0tXb7Kq0rSdIo4HLgBGAYWC1puZk9UFvGzP40tvwfAYfHNvGKmR1WVHnrZVEx9XsYNFL7zGWc3eTjDK5XJakvgSeB04Fzmmzmi8APk+yvzBbDbGDIzNab2TZgKTC3xfKnAd/JoyCddiN5KKSXxef3VoMbIG3rSzN73MzuA3bUryzpPcA+wA+S7KzMYJgEPBV7Pxym7UTSFGAacHts8jhJayTdJemU/IqZvUEPhRo/Ds4llri+rCdpF+B/07wlsZNeOStpPrDMzF6PTZtiZhskHQTcLul+M1tXv6KkhcBCgHG77Ja6IGm+pXpFuLO0XUt+lpKrqs2vv5llLx2RcOmbJkpaE5uwxMyWZFSUzwIrzGxYSlYHlRkMG4ADYu8nh2mNzAc+F59gZhvCv+sl3Uk0/rBTMISDuwRg/K57l3b/Bw+F1n45RZW+PYdzOXvezGa1mN9JfVnvd4D3S/ossBswRtKvzGynAeyaMoNhNTBd0jSiDzgf+IP6hSTNACYAP45NmwC8bGZbJU0E3gdc2k0hOhlf6La1UHQobD1wW2bbGvvkmMy21U634eCtBjcAEtWXjZjZf6u9lnQ6MKtVKECJwWBm2yWdBawERgFXmdlaSYuBNWZWOwV1PrDUzOI1xjuBb0jaQTROcnHd6HxlFBEKWQZBu20XGRTOuUiS+lLSkcCNRF+kPyzpC2Z2SDf7K3WMwcxWACvqpl1Q9/6iBuv9CDg018LVqdoZMHmGQdL95hESeXcp+Smrrle1qy/NbDVRF1OrbVwNXN1uX37lc47yaC1sPXBbaaFQL6+ydHPcqhbczvWyXjkrqedkHQpVCYNGamXzbibn+oO3GBIo89tolVoI7WRZVj+Ly7nyeDDkIKtKrVcCoZ6Hg3O9zYOhgnqpldBMr5ffuUHmwZCxtN9y+6lCzSLgOjmePgDtXDYGOhjSPN85D/0UCnH9+rmc61cDHQxJdPItNE1rod8rz37/fM71Ew+GCvBKszUfhHauWB4MGfHKqz0PQOd6gwdDyQatshy0z+tcL/JgKNGgVpLdfG5vkTlXHA8G55xzI/i9kjLQzbfZIloLUyc/1/W6jw/vlWFJdrb1wG253Fup3bMZ/M6qzrXnwdBCXhdM5RkKacKg2XbyDgnnXLV4V1KfmDr5ucxCoahtdxqQPs7gXDE8GAqWR2shr0BotJ+i9uWcK48HQ0plfostq6LOcp+DemaWc1XmwdCjyv7mXvb+nXP5KTUYJM2R9LCkIUmLGsw/XdJzku4NP2fG5i2Q9Gj4WVBsybuT1bfjqlTKVSmHc4MgQX15tKR7JG2XNC82/TBJP5a0VtJ9kj7Wbl+lnZUkaRRwOXACMAyslrTczB6oW/R6Mzurbt09gAuBWYABd4d1NxdQ9FJVrTKulSfNmUt5nbrqXL9IWF8+CZwOnFO3+svAJ83sUUn7E9WXK83sxWb7K7PFMBsYMrP1ZrYNWArMTbjuScAqM9sUwmAVMCenclZG1UKhDH5mkhtQbetLM3vczO4DdtRNf8TMHg2vfwE8C7T8JlfmdQyTgKdi74eBoxos9xFJRwOPAH9qZk81WXdSo51IWggsBBi3y24j5o15aLiwZzKk7UaqeihMnfyctxqci3lp+zhWbZyRcOmbJkpaE5uwxMyWxN4nrS9bkjQbGAOsa7Vc1Qef/xmYambvJmoVXNPpBsxsiZnNMrNZY3Z5U6aFK+rba9VDoaZXyulcBT1fq6fCz5L2q3RG0n7AdcCnzGxHq2XLDIYNwAGx95PDtDeY2QtmVru/wTeB9yRd15XDw8G5XKSq8yS9DbgJON/M7mq3fJnBsBqYLmmapDHAfGB5fIGQcDUnAw+G1yuBEyVNkDQBODFMq6Q03Ui9WNF2W2a/psG5ptrWl82E5W8ErjWzZUnWKS0YzGw7cBZRhf4gcIOZrZW0WNLJYbE/DqdY/Qz4Y6IRd8xsE/BFooO1GlgcpvWVXgwF51z2ktSXko6UNAycCnxD0tqw+keBo4HTY6f+H9Zqf6XeRM/MVgAr6qZdEHt9HnBek3WvAq7KtYCua2kHo51zIyWoL1cTdTHVr/dt4Nud7Kvqg88Dy1sLzfkpq87ly4MhZ4Pcb95NuA3y8XKuKjwYKqifWgv99FmcGxQeDM4550bwYHC5q0qrwR/r6VwyHgxdymsAtCqVaJl8nMG5cnkwuEJ44DnXOzwYctTpN1+vPJ1zVeDB4HqSX8vgXH5KvfLZFeuEfR9qu0zy2wR3rpOrof023M6Vx4OhIvLsRkoSCI2WzTMk8rTl4LGMX7e1/YLOuYY8GPpYJ4HQav0sA8LvoeRc9fkYQ59KGwp5bcs5V30eDH0oj4o8y20m7Tbz6xmcK4cHQwVkOb6Q57f7E/Z9qFKtBz8zybl8eDDkpIxvu0VV2lUKB+dc9gY+GPrl/jlFV9Zp9+cX8zlXXQMfDP2grG/wRezXxxmcK16pwSBpjqSHJQ1JWtRg/tmSHpB0n6TbJE2JzXs99vzSRA/FdtVSdKth24ydnnroXM9IUF8eLekeSdslzaubt0DSo+FnQbt9lRYMkkYBlwMfBGYCp0maWbfYT4FZZvZuYBlwaWzeK2Z2WPg5uZBC5yBt5Vh2f3/Z+/cBaDcIEtaXTwKnA/9Qt+4ewIXAUcBs4EJJE1rtr8wWw2xgyMzWm9k2YCkwN76Amd1hZi+Ht3fR4EHXrnxpwiFJMHbTnbTl4LHdFMe5qkpSXz5uZvcBO+rWPQlYZWabzGwzsAqY02pnZV75PAl4KvZ+mCjRmjkDuDn2fpykNcB24GIz+16jlSQtBBYCjNtlt1QFrpqyv63HnbDvQz17Cw3n8rBt2+hOrvKfGOqzmiVmtiT2vtP6Mq7RupNardATt8SQ9HFgFvCB2OQpZrZB0kHA7ZLuN7N19euGg7sEYPyue1shBR5QHg7Ode15M5tVdiFqyuxK2gAcEHs/OUwbQdLxwPnAyWb2xp3RzGxD+Hc9cCdweJ6FrZoqtRbSStud5OMMbgAkqi+zWrfMYFgNTJc0TdIYYD4w4uwiSYcD3yAKhWdj0ydIGhteTwTeBzxQWMkz0o/n8vdTYDlXIW3ryxZWAieGenMCcGKY1lRpwWBm24GziAr4IHCDma2VtFhS7SyjvwJ2A/6x7rTUdwJrJP0MuINojKHngqFbVa98uylfHiHpA9CuXySpLyUdKWkYOBX4hqS1Yd1NwBeJwmU1sDhMa6rUMQYzWwGsqJt2Qez18U3W+xFwaL6l655flJXPeEMWD+/ZNmNy31zt7gZLgvpyNU3O3DSzq4Crku7Lr3zuMVVvLaSRptXg4wzOZceDweUmjxDrtDXm3UnOdc6DweWq03DwVoNz5fNgKEk3FWA/dyMVye+Z5FxrHgxd8G+mncm61eDdSc7ly4PBFaKo1k7S0PZWg3PNeTD0iEHrRvJWg3Pl8WBwhfFWg3O9wYPBFaqTcPBWg3Pl8GBwhSsiHLzV4Fz3PBh6wKCNL2SlUTh4q8G59jwYStCPd1XtVJW6lLzV4NxIHgyuNFXqUnLO/UZPPMHNFWPe2+7Zadqyl47IdZ9FPPXtl1PEW58Y+fC+LQePZfy6rU3WcG6weTBkLOtbbuc9vtAoDJrNzzsk2pk6+bmWz9Dt9LbcHg7ONeZdSQNq3tvuaRsKWayTRJldSrXxBh9ncO43PBgGUNrKPY+A8PEG56rDg2HAZFmhZx0QRZyW66ewOtdeqcEgaY6khyUNSVrUYP5YSdeH+T+RNDU277ww/WFJJxVZ7l6VRzdQ1ttNGg5pTmH1loPrRd3Wl5J2lXSNpPslPSjpvHb7Ki0YJI0CLgc+CMwETpM0s26xM4DNZvZ24DLgkrDuTGA+cAgwB7gibM81kVco5LH9MsIhPtbg4w2uatLUl8CpwFgzOxR4D/CH8S/ZjbQNBkl/JGlCJx8iodnAkJmtN7NtwFJgbt0yc4FrwutlwHGSFKYvNbOtZvYYMBS211ey6lrJOxTi+ylqXzVZtxzqAyL+E5/uXMHS1JcGvEXSaOBNwDbgpVY7S9Ji2AdYLemG0JTJqh0+CXgq9n44TGu4jJltB7YAeyZcFwBJCyWtkbRm245XMip69wbhqucswqGoK6N/OUVvBEStwq+FQ7uxBw8Il6GJtXoq/Cysm5+mvlwG/Bp4GngS+Gsz29SqMG2vYzCzv5D0v4ATgU8BX5N0A3Clma1rt37ZzGwJsARg/K57W5vF+07R3+Dj+0173UMnF7+lvcYhCoffBEGjUIgHQbwFUTPmoeFEZXWDQdvUyXU1z5vZrJyKMht4HdgfmAD8q6RbzWx9sxUSjTGYmQEbw8/2sPFlki5NUdgNwAGx95PDtIbLhGbQeOCFhOu6klWx5dBN11L8p5VGXU/OZSRNffkHwC1m9pqZPQv8O9AyhJKMMfxPSXcDl4YNHmpmnyEaxPhIkk/UxGpguqRpksYQDSYvr1tmObAgvJ4H3B5CajkwP4zCTwOmA/+Roix9qazWQn0Z0pYjy3CA5F1LjSQNCWg+TuFcF9LUl08CxwJIegvw20DL/1RJbomxB/D7ZvZEfKKZ7ZD0XxOs35CZbZd0FrASGAVcZWZrJS0G1pjZcuBK4DpJQ8AmooNBWO4G4AGiFsznzOz1bsvi8pe2aynLbiVI2rXETvdYiquFQye31egkHPLqmqpyQHl3XGNp6kuis5m+JWktIOBbZnZfq/0pCpTBMH7Xve29E0/daXqz/yjNvhW2+kbZ7l5JnQw+pz0rqQothnppxx06ueFeu3AAEvUBtwqHen7vpXyUERi3bLzi7jT9/uMmHWBTPn12omUfueDsVPvKml/53ETZV8P2YyhA+nLl0a3ULszbdS/FdTIm4ZLz7rhieTC4wlUtHCDZXXFrAdHJldP1QeFhkZ6HRP78ttt9qKqthbiixxygfddSLRySdC8lGYNoplU4eFdUZ2rh4GMT2fJgcKUpMhwg2aA0dBcQ0F1I1Ou0ReFBEtk2Y7KHQ4Y8GFypaq2bbgOi1q2UdesBOgsIyD4kkiira6qKgeSth+x4MHTI78yZj6q2HqDzgICd/06KCoqiJA2kMgLEWw/peTBkKOvHeg6aMsIBkrUeoLuAqGn0haLfwqKRZgGSd2B4OKTjZyX1mV4YeG6lyDOWajq9sWHtFNe0XwTiZzl1erZTryviTC0/a6l73mJwlZNFywE6uxiu09ZDTZpWRCPtwqFfWxndXEWehLccujPwweDfKqqp6Luz1qQNCMguJBpJ0qro5fCItx6yCgkPh84NfDAUaRCexZCltGcsQXfhAN0HBBQXEs0k7ZKqeoDk1Ypw7XkwuMoro2upJk1AwM4nJJQRFM30Sutjy8FjU4eDtxo644PPFZTVIz37SdHPdqg3dfJzmbT44gPXvXAWW6MB8jIGyf12IsXyYHA9I6twqEJA1NQHRS+EBZR3RlWacPDxxOQ8GPpIr5+qmkRWnzFtq6wWEHmMG/VDWOTJWw758zEG13OyGJSGdGMPcWnHIZJoFA5VGq+ol/ftQbodd/CxhmQ8GFzPyuKUVsg+ICDfkKhp1pKoWmDkFRJZDEq7xjwYXE/LKhwgu4CA4kMirsqtizS3K2/EwyEfpYwxSNpD0ipJj4Z/JzRY5jBJP5a0VtJ9kj4Wm3e1pMck3Rt+Div2E7gqyXpsJeuzwvIcj0iqauMWZd7+o1cHoSXNkfSwpCFJixrMHyvp+jD/J5Kmxua9O1af3i9pXKt9lTX4vAi4zcymA7eF9/VeBj5pZocAc4C/kbR7bP6fm9lh4efe/Ivsqmze2+7JNCDSnr3UTBVCoqbsoMhqoHoQBqMljQIuBz4IzAROkzSzbrEzgM1m9nbgMuCSsO5o4NvAp0N9egzwWqv9lRUMc4FrwutrgFPqFzCzR8zs0fD6F8CzQLFtctdz8mg95HVdSZVCAsoLCg+HRGYDQ2a23sy2AUuJ6tG4eL26DDhOkoATgfvM7GcAZvaCmb3eamdlBcM+ZvZ0eL0R2KfVwpJmA2OAdbHJXw5dTJdJavpXIWmhpDWS1mzb8Urqgrvqy+O03TwDAkaGRBWDIm9Fdy1VsDtpYq2eCj8L6+ZPAp6KvR8O0xouY2bbgS3AnsA7AJO0UtI9kj7frjC5DT5LuhXYt8Gs8+NvzMwkNR2JkrQfcB2wwMx2hMnnEQXKGGAJcC6wuNH6ZrYkLMP4XfdONOJV9rePVRtndFUJLXvpiIG4liGJLAel47IcoG6lzMHrRoq4/1PagemqDUSP2tbRZ3nezGblVJTRwH8GjiTqor9N0t1mdlurFXJhZsc3myfpGUn7mdnToeJ/tslybwNuAs43s7ti2661NrZK+hZwToZFd30iq+sdGikqIGDnmy+WHRRZ32rcJbIBOCD2fnKY1miZ4TCuMB54gah18UMzex5A0grgCKLx3YbK6kpaDiwIrxcA369fQNIY4EbgWjNbVjdvv/CviMYnfp5raV1Py7MVlXcXUyNV6XLKq5spTbdS2a39HK0GpkuaFurG+UT1aFy8Xp0H3G5mBqwEDpX05hAYHwAeaLWzsoLhYuAESY8Cx4f3SJol6ZthmY8CRwOnNzgt9e8l3Q/cD0wEvlRs8V2vybuLrYyAgGqERB4BMUhPs0sijBmcRVTJPwjcYGZrJS2WdHJY7EpgT0lDwNmEsz3NbDPwVaJwuRe4x8xuarU/RYEyGMbvure9d+Kpb7xvNgDV6ltHqz/YJP85kv4H7raS8TGG9vLoWqpXRBdTK2V1N2XdvdTNeEMn4wytbo9xy8Yr7k7T7/+WvQ6wd87900TL3n3ln6XaV9b8Jnpu4BQRnmW1IGrKaklk3XropuXQx91JhfFgcAOpqJZV2QEB2d8qPImyr6x26Xgw9Jkiukn6RZHdboMYEFmFg483FM+DwQ20osdkBi0gymo5eHdSOh4MFVX24OUgKWPAfpACIotw8FZDsTwY+pB3J3WurLO5BiUgqjrmUMFbY1SCB4NzQZmn+lYlIPJU1XBwO/NgcC6m7OtAqhAOZV9R3Yx3JxXHg6FPeXdS96oQDlUIiDwU2WrwAejueTBUmA9Al6fscIBqtB7cYPJgyFCS2wEUeasCbzX0vrJbD3mEQ5pWg3cnFcODwbkmqtBqqCm79eAGiwdDn/NWQ/8oKxyq1mpw+fNgqLgsxhk8HPpHP4WDq66BDQa/sMUlUaXuJBfxcYb8DWwwlKnTAWhvNbi4fmk1FNGdlOSUVf+SuDMPhgHi4dA/fDDa5cmDYcB4OHTGj5cbRKUEg6Q9JK2S9Gj4d0KT5V6PPe95eWz6NEk/kTQk6frwcOy+luXFbl7ZuW75IHR5JM2R9HCo9xY1mD821IdDoX6cWjf/QEm/knROu32V1WJYBNxmZtOB28L7Rl4xs8PCz8mx6ZcAl5nZ24HNwBn5Fje5pM+8LeuZvDUeDu35MXJVIWkUcDnwQWAmcJqkmXWLnQFsDvXiZUT1ZNxXgZuT7K+sYJgLXBNeXwOcknRFSQKOBZZ1s34vy/oWGV7xOdczZgNDZrbezLYBS4nq0bh4vboMOC7Ul0g6BXgMWJtkZ2UFwz5m9nR4vRHYp8ly4yStkXRX+GAAewIvmtn28H4YmNRsR5IWhm2s2bbjlUwKX6Y8wsEDYiQ/Jq4EE2v1VPhZWDd/EvBU7H2jeu+NZUL9uAXYU9JuwLnAF5IWZnSnpU9K0q3Avg1mnR9/Y2YmyZpsZoqZbZB0EHC7pPuJPmxiZrYEWAIwfte9m+0nsbc+YZmdR/348F6V6bNd9tIRfs4+3opy2Rn1qjF+3dakiz9vZrNyKspFRF3vvwoNiLZyCwYzO77ZPEnPSNrPzJ6WtB/wbJNtbAj/rpd0J3A48F1gd0mjQypOBjZk/gEqbNXGGbmcrlirFAcxIHotEPzOuwNnA3BA7H2jeq+2zLCk0cB44AXgKGCepEuB3YEdkl41s68121lZXUnLgQXh9QLg+/ULSJogaWx4PRF4H/CAmRlwBzCv1fplSjoADd0PQudZMQxSV8ogfdYslH3SxABbDUwPZ2SOAeYT1aNx8Xp1HnC7Rd5vZlPNbCrwN8BftgoFyLHF0MbFwA2SzgCeAD4KIGkW8GkzOxN4J/ANSTuIAuxiM3sgrH8usFTSl4CfAlcW/QGqIK+WQ028wuy3VkQvh4G3FgaPmW2XdBawEhgFXGVmayUtBtaY2XKievA6SUPAJhiOR3EAAA+WSURBVKLw6EopwWBmLwDHNZi+BjgzvP4RcGiT9dcTjdL3hTRjDXmHQ00/hEQvh0GNh8LgMrMVwIq6aRfEXr8KnNpmGxcl2VdZLYa+N/bJMYXdWriocKjplZDohyCI81BwRfFgqIi0ZygVHQ41zSrfIgOj3wKgkSqEQtbjC52MxblieTB0IctTVrNUVjg0MgiVdVH6MRRctflN9HLU6TeiLP7zVaEScdlYtXGG/z4beOsTqS9Hcm14MDTQwUUpmcsqHLxC6V1V+/15a2HweFdSzroZhM7qiugqdS251qoUBHF5hYKPL1SbB0OXqjrOUK9W4XhAVFNVAwG8pTDIPBgKUGarocYDohqqHARxeYZCka2FMruFe5kHQ4XlcZO9eMXkIVGMXgkDyL+VkDYUfOC5GB4MKXTSndTtBW+1/6h53IXVWxHZ66UQqOddR67Gg6GJ8eu2suXgsZluM83V0Hneoru+MvOgSKaXQyCuqEDwAefe4cHQQ4p6foMHxW/0S+XfSJEthCxCwbuRiuPBkFKnZyelvYdSnl1LzTSrHPshMPq54m+m6C6jsloKPvDcvYENhjEPDbNtxuRS9p3FDfaq8PS3pJVqUQEyiJV8UmWNH2QVCt5aKNbABkMSSccZurmmIatwgGJbD93wCrscZQ8m98qYwpiHhssuQuV4MGSkrHCA3gkIl6+ygyAuy1Dw1kLxPBhKluVzG+IVg4dE/6tSENRk3UroNhR8fCEdD4Y2OjlttdvbZOTxUB8Pif5RxQBopFe6jlx7pQSDpD2A64GpwOPAR81sc90yvwtcFps0A5hvZt+TdDXwAWBLmHe6md2bc7ETSRMOQC5PffOQqL5eqfwbySsQvAtpJElzgL8leubzN83s4rr5Y4FrgfcALwAfM7PHJZ0AXAyMAbYBf25mt7faV1kthkXAbWZ2saRF4f258QXM7A7gMHgjSIaAH8QW+XMzW5amEGWemdRM3o8EbVQBeVjkq5cr/VbybCGkCYV+7EaSNAq4HDgBGAZWS1puZg/EFjsD2Gxmb5c0H7gE+BjwPPBhM/uFpHcBK4FJrfZXVjDMBY4Jr68B7qQuGOrMA242s5fzLVZjnV4FnfbOq3m2HhppVnF5YDTXr5V9O0V0FxXZUuihM5JmA0Nmth5A0lKiejQeDHOBi8LrZcDXJMnMfhpbZi3wJkljzaxpgpYVDPuY2dPh9UZgnzbLzwe+Wjfty5IuAG4DFrX6kFkoOhyg+ICol6Ty67XwGNQKPY0ixw7ShkI/thaCScBTsffDwFHNljGz7ZK2AHsStRhqPgLc066+zC0YJN0K7Ntg1vnxN2Zmkpr+NUjaDziUqPlTcx5RoIwBlhC1NhY3WX8hsBBg3C67dfAJ0svqmQ1lB0QrXtH2pzIGkvttTEGvbuukRTJR0prY+yVmtiTT8kiHEHUvndhu2dyCwcyObzZP0jOS9jOzp0PF/2yLTX0UuNHMXottu9ba2CrpW8A5LcqxhCg8GL/r3jv95XUyztDNjfVqf+xZBgRUMyRcbyr7bKKsAqHHWwvPm9msFvM3AAfE3k8O0xotMyxpNDCeaBAaSZOBG4FPmtm6doUp65nPy4EF4fUC4Pstlj0N+E58QggTJAk4Bfh5DmXMVNbfhsY+OeaNH+eSiv/dVOHvp8xQ6KHxBYDVwHRJ0ySNIepeX163TLxenQfcHnpkdgduIupy//ckOytrjOFi4AZJZwBPELUKkDQL+LSZnRneTyVKwH+pW//vJe0FCLgX+HSawuTdaqjJ63Gg9f+5vTUxuMqu6JPK8otSj7cUEgljBmcRdamPAq4ys7WSFgNrzGw5cCVwnaQhYBNReACcBbwduCCMywKcaGZNe2pKCQYzewE4rsH0NcCZsfeP0+C0KjM7Ns/ytZM2HCCbrqVmPCj6R69U9En12zhCkcxsBbCibtoFsdevAqc2WO9LwJc62Zdf+Rx0ek1D2gf5FBEQNa0qFw+N7PVbZZ6FvAKh29ZCj3UjFc6DIYUsnvJWZEA00k0l1ith4hV0ufJuHQxCF1JZPBhiurkSOqtHgJYdEJ3wCtc1UlQ3UdpA8NZCex4MGcjy+dDx/1y9EBJucJUxXuCthGJ4MNTp9v5JtT/YrAICdv6P50HhylCVAeMsQsFbC8l4MDSQ5uZ6WbYe6jX6D+ph4bpVlQq/naxaCR4KyXkwNJE2HCDb1kMzzf5ze2D0tl6ptPOUZbeRh0JnPBhaSHtb7iIDol43FUvVw8Qry/7nYwjV4MHQRhbPbIj/sZcREkl5xevKkGcYeEuhOx4MCWT5QJ9eCQnn8lREy8BDoXseDAnl8bQ3Dwk3CMroHvJQSMeDoQO1P7Y8Hgfa6D+Ph4XrBVUbF/BQSM+DoQtFPSu63X+4QQ+OqlVIrlweCNnxYOhSnq2HpLxidM4DIQ8eDCnF/yjLDAnnBomHQb48GDLkIeFcfjwMiuPBkJP6P2IPCuc640FQHg+GgjT6I/ewiHgF0Nig/H347796PBhKlPQ/RBkVhP9nLV9Vfget/v6qUkaXrVKCQdKpwEXAO4HZ4VnPjZabA/wt0cOvv2lmF4fp04ClwJ7A3cAnzKw3HivWBf/P58rkf3/V0Kw+jM0fC1wLvAd4AfiYmT0e5p0HnAG8Dvyxma1sta9dMi99Mj8Hfh/4YbMFJI0CLgc+CMwETpM0M8y+BLjMzN4ObCb6wM4515fa1Ic1ZwCbQ714GVE9SVhuPnAIMAe4ImyvqVKCwcweNLOH2yw2Gxgys/WhNbAUmCtJwLHAsrDcNcAp+ZXWOedK17A+rFtmLlF9CFH9eFyoL+cCS81sq5k9BgyF7TVV5TGGScBTsffDwFFE3Ucvmtn22PRJzTYiaSGwMLzdesvGK36eQ1k7MRF4vuQyQDXKUYUyQDXKUYUyQDXKUYUyAPxWmpVf2v7cyls2XjEx4eLjJMW71JeY2ZLY+2b1IY2WMbPtkrYQ1ZeTgLvq1m1aZ0KOwSDpVmDfBrPON7Pv57XfeuHgLgllWmNms4radyNVKENVylGFMlSlHFUoQ1XKUYUy1MqRZn0zm5NVWYqWWzCY2fEpN7EBOCD2fnKY9gKwu6TRodVQm+6cc/2qWX3YaJlhSaOB8UT1ZZJ1Ryhr8DmJ1cB0SdMkjSEaPFluZgbcAcwLyy0ACmuBOOdcCRrWh3XLLCeqDyGqH28P9eVyYL6kseGMzunAf7TaWSnBIOn3JA0DvwPcJGllmL6/pBUQ9ZEBZwErgQeBG8xsbdjEucDZkoaI+tCuTLjrJe0XyV0VygDVKEcVygDVKEcVygDVKEcVygDVKUfT+lDSYkknh8WuBPYM9eLZwKKw7lrgBuAB4Bbgc2b2eqv9KQoU55xzLlLlriTnnHMl8GBwzjk3Qt8Fg6RTJa2VtENS01PeJM2R9LCkIUmLYtOnSfpJmH59GOjptAx7SFol6dHw74QGy/yupHtjP69KOiXMu1rSY7F5h3VahqTlCMu9HtvX8tj0oo7FYZJ+HH5v90n6WGxe18ei2e84Nn9s+FxD4XNOjc07L0x/WNJJnX7uDstxtqQHwme/TdKU2LyGv5scynC6pOdi+zozNm9B+P09KmlB/boZl+OyWBkekfRibF5Wx+IqSc9KanhNkyJ/F8p4n6QjYvMyOxaVZmZ99UN0/6XfAu4EZjVZZhSwDjgIGAP8DJgZ5t0AzA+vvw58posyXAosCq8XAZe0WX4PYBPw5vD+amBeBsciUTmAXzWZXsixAN4BTA+v9weeBnZPcyxa/Y5jy3wW+Hp4PR+4PryeGZYfC0wL2xnV5e8gSTl+N/a7/0ytHK1+NzmU4XTga03+NteHfyeE1xPyKkfd8n8EXJXlsQjbORo4Avh5k/kfAm4GBPw28JOsj0XVf/quxWDVuN1G/NL0JNuYB9xsZi93sa8sy/GGIo+FmT1iZo+G178AngX26mJfcYXeQiBNOczsjtjv/i6i88yzlORYNHMSsMrMNpnZZmAV0f12iijHacB3utxXU2b2Q6IvYs3MBa61yF1E103tR7bHotL6LhgSanR5+SQ6vN1GC/uY2dPh9UZgnzbLz2fn/wBfDs3YyxTdNbEbScsxTtIaSXfVurMo6VhImk30bXJdbHI3x6LZ77jhMuFzxm8h0G7dpDrd1hlE31ZrGv1u8irDR8JxXiapdkFUKccidKdNA26PTc7iWCTRrJxZHotKq/K9kppSBW630aoM8TdmZpKanhMcvokcSnR+cs15RJXoGKJzqc8FFudYjilmtkHSQcDtku4nqiQTyfhYXAcsMLMdYXLiY9HrJH0cmAV8IDZ5p9+Nma1rvIVU/hn4jpltlfSHRC2pY3PYT1LzgWU28nz7oo7FwOvJYLAK3G6jVRkkPSNpPzN7OlR2z7Yoy0eBG83stdi2a9+wt0r6FnBOs5WzKIeZbQj/rpd0J3A48F0KPBaS3gbcRBTub9zwq5NjUafQWwikLAeSjicK0g+Y2dba9Ca/m04rw7ZlMLMXYm+/STQ2VFv3mLp17+xw/4nLETMf+FxdGbM4Fkk0K2eWx6LSBrUrKe/bbcQvTW+3jZ36UUMFWuvnP4Xo+RXdaFsOSRNq3TOSJgLvAx4o8liE38GNRP26y+rmdXssCr2FQJpySDoc+AZwspk9G5ve8HeTUxn2i709mejqWohasieGskwATmRk6zbTcoSyzCAa3P1xbFpWxyKJ5cAnw9lJvw1sCV9QsjwW1Vb26HfWP8DvEfX9bQWeAVaG6fsDK2LLfQh4hOgbx/mx6QcRVQJDwD8CY7sow57AbcCjwK3AHmH6LKInL9WWm0r0LWSXuvVvB+4nqgS/DezW5bFoWw7gvWFfPwv/nlH0sQA+DrwG3Bv7OSztsWj0Oybqhjo5vB4XPtdQ+JwHxdY9P6z3MPDBlH+T7cpxa/hbrX325e1+NzmU4SvA2rCvO4AZsXX/ezhGQ8Cn8jwW4f1FwMV162V5LL5DdObba0R1xRnAp4FPh/kieijOurCvWbF1MzsWVf7xW2I455wbYVC7kpxzzjXhweCcc24EDwbnnHMjeDA455wbwYPBOefcCB4MzjnnRvBgcM45N4IHgxtIko4MN4wbJ+ktip4F8a6yy+VcFfgFbm5gSfoS0dXPbwKGzewrJRfJuUrwYHADK9yvZzXwKvBeG3knT+cGlncluUG2J7Ab8FailoNzDm8xuAGm6LnBS4keCLOfmZ1VcpGcq4SefB6Dc2lJ+iTwmpn9g6RRwI8kHWtmt7db17l+5y0G55xzI/gYg3POuRE8GJxzzo3gweCcc24EDwbnnHMjeDA455wbwYPBOefcCB4MzjnnRvj//m9XhvsQ4MwAAAAASUVORK5CYII=\n",
            "text/plain": [
              "<Figure size 432x288 with 2 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "PqkiDt5IRTsH"
      },
      "source": [
        ""
      ],
      "execution_count": null,
      "outputs": []
    }
  ]
}