msc_thesis/Archive/Untitled.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "9005a683-4dd4-49e8-b7fb-5bd7c695b8d3",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "#importing modules\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "from scipy.optimize import curve_fit"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "78801aba-2566-4d5d-81e9-22701bdf4dcd",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAisAAAGdCAYAAADT1TPdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAgO0lEQVR4nO3dcWyc910/8M85bex2ja9NTGJ3DY3XoQ7LYp3TZXXVwljVLAOZ7Q8mBkrVClGoWdi6Fe3XrJNSF0WhYoCQYKnGpm1SpF8rkRWwAKvRmrarlpClcTSMRYFiSLSdidKMcyiys9nP74/8bOL47MaO7+57d6+XdH/4ue85H3/1Ve/d5/k+nyeXZVkWAACJaqp2AQAASxFWAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJJ2TbULuFozMzPxgx/8INatWxe5XK7a5QAAVyDLsjh//nzcfPPN0dS09LmTmg8rP/jBD2Lz5s3VLgMAWIHTp0/HLbfcsuSYmg8r69ati4iLf2xra2uVqwEArsTExERs3rx57nt8KTUfVmYv/bS2tgorAFBjrmQLhw22AEDShBUAIGnCCgCQNGEFAEiasAIAJE1YAQCSJqwAAEkTVgCApNV8U7hqm57J4tjYuThzfjI2rmuJbZ3rY02TZxQBwGoRVq7C0EghBgZHo1CcnDvWkW+JPX1dsaO7o4qVAUD9cBlohYZGCtF/4MS8oBIRMV6cjP4DJ2JopFClygCgvggrKzA9k8XA4GhkJd6bPTYwOBrTM6VGAADLIayswLGxcwvOqFwqi4hCcTKOjZ2rXFEAUKeElRU4c37xoLKScQDA4oSVFdi4rmVVxwEAixNWVmBb5/royLfEYjco5+LiXUHbOtdXsiwAqEvCyiKmZ7I48sab8Vcnvx9H3nhz3mbZNU252NPXFRGxILDM/rynr0u/FQBYBfqslHAl/VN2dHfE/p09C8a167MCAKsql2VZTd9fOzExEfl8PorFYrS2tl7175vtn3L5pMyeI9m/s2deENHBFgCWbznf386sXOLt+qfk4mL/lPu72ucCyZqmXPTetqGSZQJAQ7Fn5RL6pwBAeoSVS+ifAgDpEVYuoX8KAKRHWLmE/ikAkB5h5RL6pwBAeoSVy8z2T2nPz7/U055vWXDbMgBQfm5dLmFHd0fc39WufwoAJEBYWYT+KQCQBpeBAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASRNWAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASRNWAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASRNWAICkCSsAQNKEFQAgacIKAJA0YQUASFoSYeVLX/pSdHZ2RktLS2zdujW+/e1vV7skACARVQ8rzz33XDz66KPxxBNPxPDwcNx7773xkY98JE6dOlXt0gCABOSyLMuqWcAHPvCB6Onpif37988d++mf/un42Mc+Fvv27Xvbz09MTEQ+n49isRitra3lLBUAWCXL+f6u6pmVCxcuxGuvvRbbt2+fd3z79u3xne98p0pVAQApuaaa//jZs2djeno6Nm3aNO/4pk2bYnx8vORnpqamYmpqau7niYmJstYIAFRX1fesRETkcrl5P2dZtuDYrH379kU+n597bd68uRIlAgBVUtWw0tbWFmvWrFlwFuXMmTMLzrbM2r17dxSLxbnX6dOnK1EqAFAlVQ0ra9euja1bt8ahQ4fmHT906FDcfffdJT/T3Nwcra2t814AQP2q6p6ViIjPfvaz8cADD8Sdd94Zvb298eUvfzlOnToVjzzySLVLAwASUPWw8iu/8ivx5ptvxlNPPRWFQiG6u7vjb//2b+PWW2+tdmkAQAKq3mflaumzAgC1p2b6rAAAvB1hBQBImrACACRNWAEAkiasAABJE1YAgKQJKwBA0oQVACBpwgoAkDRhBQBImrACACRNWAEAkiasAABJE1YAgKQJKwBA0oQVACBpwgoAkDRhBQBImrACACRNWAEAkiasAABJE1YAgKQJKwBA0oQVACBpwgoAkLRrql1AI5meyeLY2Lk4c34yNq5riW2d62NNU67aZQFA0oSVChkaKcTA4GgUipNzxzryLbGnryt2dHdUsTIASJvLQBUwNFKI/gMn5gWViIjx4mT0HzgRQyOFKlUGAOkTVspseiaLgcHRyEq8N3tsYHA0pmdKjQAAhJUyOzZ2bsEZlUtlEVEoTsaxsXOVKwoAaoiwUmZnzi8eVFYyDgAajbBSZhvXtazqOABoNMJKmW3rXB8d+ZZY7AblXFy8K2hb5/pKlgUANUNYKbM1TbnY09cVEbEgsMz+vKevS78VAFiEsFIBO7o7Yv/OnmjPz7/U055vif07e/RZAYAlaApXITu6O+L+rnYdbAFgmYSVClrTlIve2zZUuwwAqCkuAwEASRNWAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASRNWAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASRNWAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASRNWAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAk7ZpqF8DVm57J4tjYuThzfjI2rmuJbZ3rY01TrtplAcCqEFZq3NBIIQYGR6NQnJw71pFviT19XbGju6OKlQHA6nAZqIYNjRSi/8CJeUElImK8OBn9B07E0EihSpUBwOqpaljZsmVL5HK5ea/HH3+8miXVjOmZLAYGRyMr8d7ssYHB0ZieKTUCAGpH1S8DPfXUU/Hwww/P/XzDDTdUsZracWzs3IIzKpfKIqJQnIxjY+ei97YNlSsMAFZZ1cPKunXror29vdpl1Jwz5xcPKisZBwCpqvqelaeffjo2bNgQd9xxR+zduzcuXLiw5PipqamYmJiY92pEG9e1rOo4AEhVVc+sfPrTn46enp646aab4tixY7F79+4YGxuLr3zlK4t+Zt++fTEwMFDBKtO0rXN9dORbYrw4WXLfSi4i2vMXb2MGgFqWy7JsVXdgPvnkk28bJr773e/GnXfeueD4wYMH45d/+Zfj7NmzsWFD6X0WU1NTMTU1NffzxMREbN68OYrFYrS2tl5d8TVm9m6giJgXWGY7rOzf2eP2ZQCSNDExEfl8/oq+v1c9rJw9ezbOnj275JgtW7ZES8vCyxPf//7345ZbbomjR4/GBz7wgSv695bzx9YjfVYAqEXL+f5e9ctAbW1t0dbWtqLPDg8PR0RER4cv2Su1o7sj7u9q18EWgLpVtT0rR44ciaNHj8bP//zPRz6fj+9+97vxmc98Jn7pl34pfvInf7JaZdWkNU05tycDULeqFlaam5vjueeei4GBgZiamopbb701Hn744fjc5z5XrZIAgARVLaz09PTE0aNHq/XPAwA1oup9VgAAliKsAABJE1YAgKQJKwBA0oQVACBpwgoAkDRhBQBImrACACRNWAEAkiasAABJE1YAgKQJKwBA0oQVACBpwgoAkDRhBQBImrACACRNWAEAkiasAABJE1YAgKQJKwBA0oQVACBpwgoAkDRhBQBImrACACRNWAEAkiasAABJE1YAgKQJKwBA0q6pdgGkaXomi2Nj5+LM+cnYuK4ltnWujzVNuWqXBUADElZYYGikEAODo1EoTs4d68i3xJ6+rtjR3VHFygBoRC4DMc/QSCH6D5yYF1QiIsaLk9F/4EQMjRSqVBkAjUpYYc70TBYDg6ORlXhv9tjA4GhMz5QaAQDlIaww59jYuQVnVC6VRUShOBnHxs5VrigAGp6wwpwz5xcPKisZBwCrQVhhzsZ1Las6DgBWg7DCnG2d66Mj3xKL3aCci4t3BW3rXF/JsgBocMIKc9Y05WJPX1dExILAMvvznr4u/VYAqChhhXl2dHfE/p090Z6ff6mnPd8S+3f26LMCQMVpCscCO7o74v6udh1sAUiCsEJJa5py0XvbhmqXAQAuAwEAaRNWAICkCSsAQNKEFQAgaTbYUlHTM5m7jABYFmGFihkaKcTA4Oi8hyV25FtiT1+X/i0ALMplICpiaKQQ/QdOLHiq83hxMvoPnIihkUKVKgMgdcIKZTc9k8XA4GhkJd6bPTYwOBrTM6VGANDohBXK7tjYuQVnVC6VRUShOBnHxs5VrigAaoawQtmdOb94UFnJOAAai7BC2W1c1/L2g5YxDoDGIqxQdts610dHviUWu0E5FxfvCtrWub6SZQFQI4QVym5NUy729HVFRCwILLM/7+nr0m8FgJKEFSpiR3dH7N/ZE+35+Zd62vMtsX9njz4rACxKUzgqZkd3R9zf1a6DLQDLIqxQUWuactF724ZqlwFADXEZCABImrACACRNWAEAkiasAABJE1YAgKQJKwBA0oQVACBpwgoAkDRhBQBImrACACRNWAEAkiasAABJE1YAgKQJKwBA0q6pdgGwGqZnsjg2di7OnJ+MjetaYlvn+ljTlKt2WQCsgrKeWdm7d2/cfffdcf3118eNN95YcsypU6eir68v3vGOd0RbW1t86lOfigsXLpSzLOrM0Egh7nn6xfjVPz8an372ZPzqnx+Ne55+MYZGCtUuDYBVUNawcuHChfj4xz8e/f39Jd+fnp6OX/zFX4y33norXn311Xj22Wfj4MGD8dhjj5WzLOrI0Egh+g+ciEJxct7x8eJk9B84IbAA1IGyXgYaGBiIiIivf/3rJd9/4YUXYnR0NE6fPh0333xzRET84R/+YTz00EOxd+/eaG1tLWd51LjpmSwGBkcjK/FeFhG5iBgYHI37u9pdEgKoYVXdYHvkyJHo7u6eCyoRER/+8IdjamoqXnvttZKfmZqaiomJiXkvGtOxsXMLzqhcKouIQnEyjo2dq1xRAKy6qoaV8fHx2LRp07xjN910U6xduzbGx8dLfmbfvn2Rz+fnXps3b65EqSTozPnFg8pKxgGQpmWHlSeffDJyudySr+PHj1/x78vlFp6ez7Ks5PGIiN27d0exWJx7nT59erl/AnVi47qWVR0HQJqWvWdl165d8YlPfGLJMVu2bLmi39Xe3h5///d/P+/YD3/4w/jRj3604IzLrObm5mhubr6i309929a5PjryLTFenCy5byUXEe35i7cxA1C7lh1W2traoq2tbVX+8d7e3ti7d28UCoXo6OiIiIubbpubm2Pr1q2r8m9Qv9Y05WJPX1f0HzgRuYh5gWX2vNyevi6bawFqXFn3rJw6dSpOnjwZp06diunp6Th58mScPHky/vu//zsiIrZv3x5dXV3xwAMPxPDwcHzrW9+K3/3d342HH37YnUBckR3dHbF/Z0+05+df6mnPt8T+nT2xo7ujSpUBsFpyWZaVOoO+Kh566KH4xje+seD44cOH44Mf/GBEXAw0v/3bvx0vvvhiXHfddfFrv/Zr8cUvfvGKL/VMTExEPp+PYrEo4DQwHWwBastyvr/LGlYqQVgBgNqznO9vDzIEAJImrAAASRNWAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASVv2U5cBIFWeE1afhBUA6sLQSCEGBkejUJycO9aRb4k9fV2ewF7jXAYCoOYNjRSi/8CJeUElImK8OBn9B07E0EihSpWxGoQVAGra9EwWA4OjkZV4b/bYwOBoTM+UGkEtEFYAqGnHxs4tOKNyqSwiCsXJODZ2rnJFsaqEFQBq2pnziweVlYwjPTbYAlDTNq5rWdVx/K9U7q4SVgCoads610dHviXGi5Ml963kIqI9f/GLliuX0t1VLgMBUNPWNOViT19XRFwMJpea/XlPX5d+K8uQ2t1VwgoANW9Hd0fs39kT7fn5l3ra8y2xf2ePPivLkOLdVS4DAVAXdnR3xP1d7Unssahly7m7qve2DRWpSVgBoG6sacpV7Au0XqV4d5XLQADAnBTvrhJWAIA5s3dXLXbxLBcX7wqq5N1VwgoAMCfFu6uEFQBgntTurrLBFgBYIKW7q4QVAKCkVO6uchkIAEiasAIAJE1YAQCSJqwAAEkTVgCApAkrAEDShBUAIGn6rECFTc9kSTRZAqgVwgpU0NBIIQYGR6NQ/N9Hq3fkW2JPX1fF21cD1AqXgaBChkYK0X/gxLygEhExXpyM/gMnYmikUKXKANImrEAFTM9kMTA4GlmJ92aPDQyOxvRMqREAjU1YgQo4NnZuwRmVS2URUShOxrGxc5UrCqBGCCtQAWfOLx5UVjIOoJEIK1ABG9e1rOo4gEYirEAFbOtcHx35lljsBuVcXLwraFvn+kqWBVAThBWogDVNudjT1xURsSCwzP68p69LvxWAEoQVqJAd3R2xf2dPtOfnX+ppz7fE/p09DdVnZXomiyNvvBl/dfL7ceSNN90FBSxJUziooB3dHXF/V3tDd7DVGI8InZxZnlyWZTX9vzQTExORz+ejWCxGa2trtcsBljDbGO/y/+jMfkU12hmmRiWwErG872+XgYCK0BiPCJ2cWRlhBagIjfEQWFkpYQWoCI3xEFhZKRtsYRE2AK4ujfEQWFkpYQVKsAFw9c02xhsvTpa8DJCLi7dxa4xXvwRWVsplILiMDYDloTEeOjmzUsIKXMIGwPLSGK+xCayslMtAcInlbADsvW1D5QqrIxrjNbbZwHr5ZdZ2l1lZgrACl7ABsDLWNOVWPezZEF07BFaWS1iBS9gAWJtsiK495Qis1C97VuASNgDWHhuiof4JK3AJGwBriw3R0BiEFbiMO1Zqh46o0BjsWYESbACsDTZEQ2MQVmARNgCmz4ZoaAwuAwE1y4ZoaAzCClCzbIiGxiCsADXNhmiof/asADXPhmiob8IKUBdsiIb65TIQAJA0Z1aARXk4IJACYQUoycMBgVSU9TLQ3r174+67747rr78+brzxxpJjcrncgtczzzxTzrKAt+HhgEBKyhpWLly4EB//+Mejv79/yXFf+9rXolAozL0efPDBcpYFLMHDAYHUlPUy0MDAQEREfP3rX19y3I033hjt7e3lLAW4Qst5OKC7b4BKSOJuoF27dkVbW1u8//3vj2eeeSZmZmYWHTs1NRUTExPzXsDq8XBAIDVV32D7e7/3e3HffffFddddF9/61rfisccei7Nnz8YXvvCFkuP37ds3d8YGWH0eDgikZtlnVp588smSm2IvfR0/fvyKf98XvvCF6O3tjTvuuCMee+yxeOqpp+IP/uAPFh2/e/fuKBaLc6/Tp08v908AluDhgEBqln1mZdeuXfGJT3xiyTFbtmxZaT1x1113xcTERPznf/5nbNq0acH7zc3N0dzcvOLfD/VoNfuhzD4csP/AichFzNto6+GAQDUsO6y0tbVFW1tbOWqJiIjh4eFoaWlZ9FZnYL5y9EOZfTjg5b+3XZ8VoArKumfl1KlTce7cuTh16lRMT0/HyZMnIyLi3e9+d9xwww0xODgY4+Pj0dvbG9ddd10cPnw4nnjiifjN3/xNZ0/gCsz2Q7n8JuLZfihX89RhDwcEUpHLsqxszRIeeuih+MY3vrHg+OHDh+ODH/xgDA0Nxe7du+Nf//VfY2ZmJt71rnfFb/zGb8QnP/nJuOaaK8tRExMTkc/no1gsRmtr62r/CZCs6Zks7nn6xUVvM87FxTMhr/6fDwkYQHKW8/1d1rBSCcIKjerIG2/Gr/750bcd938fvks/FK6KZ0RRDsv5/q76rcvAyuiHQiV4RhQpSKIpHLB8+qFQbp4RRSqEFahR+qFQTp4RRUqEFahRs/1QImJBYNEPhau1nGdEQbkJK1DDZvuhtOfnX+ppz7dc1W3LYE8UKbHBFmqcfiiUgz1RpERYgTqwpinn9mRW1eyeqPHiZMl9K7N9fOyJohJcBgKosOmZLI688Wb81cnvx5E33ly1Taqr+XvtiSIlzqwAVFC5+pZ4RhT1TAdbgApZ7FlOs+cmVropuly/d5YOtpTDcr6/XQYCqIBy9S2pRD+U2T1RH73jndF72wZBhYoTVgAqoFx9S/RDoREIKwAVUK6+Jfqh0AiEFYAKKFffEv1QaATCCkAFlOtZTp4RRSMQVgAqoFx9S/RDoREIKwAVUq5nOXlGFPVOnxWACitX3xL9UKgly/n+1sEWoMLK9Swnz4iiXrkMBAAkTVgBAJImrAAASRNWAICkCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASRNWAICkeZAhwBI8yRiqT1gBWMTQSCEGBkejUJycO9aRb4k9fV2xo7ujipVBY3EZCKCEoZFC9B84MS+oRESMFyej/8CJGBopVKkyaDzCCsBlpmeyGBgcjazEe7PHBgZHY3qm1AhgtQkrAJc5NnZuwRmVS2URUShOxrGxc5UrChqYsAJwmTPnFw8qKxkHXB1hBeAyG9e1rOo44OoIKwCX2da5PjryLbHYDcq5uHhX0LbO9ZUsCxqWsAJwmTVNudjT1xURsSCwzP68p69LvxWoEGEFoIQd3R2xf2dPtOfnX+ppz7fE/p09+qxABWkKB7CIHd0dcX9Xuw62UGXCCsAS1jTlove2DdUuAxqay0AAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASRNWAICkCSsAQNJqvoNtlmURETExMVHlSgCAKzX7vT37Pb6Umg8r58+fj4iIzZs3V7kSAGC5zp8/H/l8fskxuexKIk3CZmZm4vXXX4+urq44ffp0tLa2VrukqpuYmIjNmzebj0uYk/nMx0LmZD7zMZ/5WOhq5yTLsjh//nzcfPPN0dS09K6Umj+z0tTUFO985zsjIqK1tdUiuoT5WMiczGc+FjIn85mP+czHQlczJ293RmWWDbYAQNKEFQAgaXURVpqbm2PPnj3R3Nxc7VKSYD4WMifzmY+FzMl85mM+87FQJeek5jfYAgD1rS7OrAAA9UtYAQCSJqwAAEkTVgCApNVMWPnSl74UnZ2d0dLSElu3bo1vf/vbS45/+eWXY+vWrdHS0hLvete74plnnqlQpZWxnPl46aWXIpfLLXj90z/9UwUrLp9XXnkl+vr64uabb45cLhd/+Zd/+bafqff1sdw5qfc1sm/fvnj/+98f69ati40bN8bHPvaxeP3119/2c/W6TlYyH/W8Rvbv3x8/8zM/M9fcrLe3N/7u7/5uyc/U69qYtdw5Kff6qImw8txzz8Wjjz4aTzzxRAwPD8e9994bH/nIR+LUqVMlx4+NjcUv/MIvxL333hvDw8Px+c9/Pj71qU/FwYMHK1x5eSx3Pma9/vrrUSgU5l4/9VM/VaGKy+utt96K9773vfGnf/qnVzS+3tdHxPLnZFa9rpGXX345PvnJT8bRo0fj0KFD8eMf/zi2b98eb7311qKfqed1spL5mFWPa+SWW26J3//934/jx4/H8ePH40Mf+lB89KMfjX/8x38sOb6e18as5c7JrLKtj6wGbNu2LXvkkUfmHXvPe96TPf744yXHf+5zn8ve8573zDv2W7/1W9ldd91Vthorabnzcfjw4Swish/+8IcVqK66IiJ7/vnnlxxT7+vjclcyJ420RrIsy86cOZNFRPbyyy8vOqaR1smVzEejrZGbbrop+8pXvlLyvUZaG5daak7KvT6SP7Ny4cKFeO2112L79u3zjm/fvj2+853vlPzMkSNHFoz/8Ic/HMePH48f/ehHZau1ElYyH7Pe9773RUdHR9x3331x+PDhcpaZtHpeH1erUdZIsViMiIj169cvOqaR1smVzMesel8j09PT8eyzz8Zbb70Vvb29Jcc00tqIuLI5mVWu9ZF8WDl79mxMT0/Hpk2b5h3ftGlTjI+Pl/zM+Ph4yfE//vGP4+zZs2WrtRJWMh8dHR3x5S9/OQ4ePBjf/OY34/bbb4/77rsvXnnllUqUnJx6Xh8r1UhrJMuy+OxnPxv33HNPdHd3LzquUdbJlc5Hva+Rf/iHf4gbbrghmpub45FHHonnn38+urq6So5tlLWxnDkp9/qomacu53K5eT9nWbbg2NuNL3W8Vi1nPm6//fa4/fbb537u7e2N06dPxxe/+MX42Z/92bLWmap6Xx/L1UhrZNeuXfG9730vXn311bcd2wjr5Erno97XyO233x4nT56M//qv/4qDBw/Ggw8+GC+//PKiX86NsDaWMyflXh/Jn1lpa2uLNWvWLDhrcObMmQXJdlZ7e3vJ8ddcc01s2LChbLVWwkrmo5S77ror/uVf/mW1y6sJ9bw+VlM9rpHf+Z3fib/+67+Ow4cPxy233LLk2EZYJ8uZj1LqaY2sXbs23v3ud8edd94Z+/bti/e+973xJ3/yJyXHNsLaiFjenJSymusj+bCydu3a2Lp1axw6dGje8UOHDsXdd99d8jO9vb0Lxr/wwgtx5513xrXXXlu2WithJfNRyvDwcHR0dKx2eTWhntfHaqqnNZJlWezatSu++c1vxosvvhidnZ1v+5l6XicrmY9S6mmNXC7Lspiamir5Xj2vjaUsNSelrOr6KMu23VX27LPPZtdee2321a9+NRsdHc0effTR7B3veEf27//+71mWZdnjjz+ePfDAA3Pj/+3f/i27/vrrs8985jPZ6Oho9tWvfjW79tprs7/4i7+o1p+wqpY7H3/8x3+cPf/889k///M/ZyMjI9njjz+eRUR28ODBav0Jq+r8+fPZ8PBwNjw8nEVE9kd/9EfZ8PBw9h//8R9ZljXe+siy5c9Jva+R/v7+LJ/PZy+99FJWKBTmXv/zP/8zN6aR1slK5qOe18ju3buzV155JRsbG8u+973vZZ///Oezpqam7IUXXsiyrLHWxqzlzkm510dNhJUsy7I/+7M/y2699dZs7dq1WU9Pz7xb7B588MHs537u5+aNf+mll7L3ve992dq1a7MtW7Zk+/fvr3DF5bWc+Xj66aez2267LWtpacluuumm7J577sn+5m/+pgpVl8fsLXOXvx588MEsyxpzfSx3Tup9jZSai4jIvva1r82NaaR1spL5qOc18uu//utz/z39iZ/4iey+++6b+1LOssZaG7OWOyflXh+5LPv/u4IAABKU/J4VAKCxCSsAQNKEFQAgacIKAJA0YQUASJqwAgAkTVgBAJImrAAASRNWAICkCSsAQNKEFQAgacIKAJC0/we7NwipjQ7Y6AAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "x_data = np.array([ 0.23547456, 0.15789474, 0.31578947, 0.47368421, 0.63157895, \n",
    "                   0.78947368, 0.94736842, 1.10526316, 1.26315789, 1.42105263, \n",
    "                   1.57894737, 1.73684211, 1.89473684, 2.05263158, 2.21052632, \n",
    "                   2.36842105, 2.52631579, 2.68421053, 2.84210526, 3.45454545 ])\n",
    "y_data = np.array([ 2.95258285, 2.49719803, -2.1984975, -4.88744346, -7.41326345, \n",
    "                   -8.44574157, -10.01878504, -13.83743553, -12.91548145, -15.41149046, \n",
    "                   -14.93516299, -13.42514157, -14.12110495, -17.6412464 , -16.1275509 , \n",
    "                   -16.11533771, -15.66076021, -13.48938865, -11.33918701, -11.70467566])\n",
    " \n",
    "plt.scatter(x_data , y_data)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "209c548e-0f5e-4695-9956-b2ecbeaee42e",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "def model_f(x,a,b,c):\n",
    "  return a*(x-b)**2+c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "b69882c8-40ba-4880-bb4f-85d69011da9c",
   "metadata": {
    "tags": []
   },
   "outputs": [],
   "source": [
    "popt, pcov = curve_fit(model_f, x_data, y_data, p0=[3,2,-16])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "af049b55-0340-47de-9568-08df509bfae4",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([  4.34571181,   2.16288856, -16.22482919])"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "popt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "28c0f5f3-8450-418e-8ae3-e10c68f4cb5c",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[ 0.19937578, -0.02405734, -0.1215353 ],\n",
       "       [-0.02405734,  0.00517302,  0.00226607],\n",
       "       [-0.1215353 ,  0.00226607,  0.29163784]])"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pcov"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "0b40b322-f464-4385-bd2b-1cde2c8f081f",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAisAAAGdCAYAAADT1TPdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAA6VElEQVR4nO3deVxVdf7H8fcFhOsC19QUTEvSRiOm3DLpp6NtZja2THvpaIspuaRppmYhpjlOtk2LbU6bM9mipbSQTrnUpJmKU0hZGY2mGKMUoIUo9/z++AaCLAJy7zn33tfz8biPzj33QB/PnPG++64uy7IsAQAAOFSY3QUAAADUhLACAAAcjbACAAAcjbACAAAcjbACAAAcjbACAAAcjbACAAAcjbACAAAcLcLuAo6V1+vVrl27FB0dLZfLZXc5AACgFizLUmFhodq2bauwsJrbTgI+rOzatUvt27e3uwwAAFAPO3bsULt27Wq8JuDDSnR0tCTzh42JibG5GgAAUBsFBQVq37592fd4TQI+rJR2/cTExBBWAAAIMLUZwsEAWwAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GgBvyic3Uq8ltZn5ym3sEito93qFd9C4WHsUQQAQEMhrByD9MwcpaZlKSe/qOxcnMetlMEJGpgYZ2NlAAAED7qB6ik9M0fJCzdVCCqStDu/SMkLNyk9M8emygAACC6ElXoo8VpKTcuSVcVnpedS07JU4q3qCgAAUBeElXpYn51XqUWlPEtSTn6R1mfn+a8oAACCFGGlHnILqw8q9bkOAABUj7BSD62j3Q16HQAAqB5hpR56xbdQnMet6iYou2RmBfWKb+HPsgAACEqElWqUeC2t3bZXSzfv1NpteysMlg0PcyllcIIkVQospe9TBiew3goAAA2AdVaqUJv1UwYmxmn+kO6VrotlnRUAABqUy7KsgJ5fW1BQII/Ho/z8fMXExBzz7ytdP+XIm1LaRjJ/SPcKQYQVbAEAqLu6fH/TslLO0dZPccmsn3JBQmxZIAkPcympY0t/lgkAQEhhzEo5rJ8CAIDzEFbKYf0UAACch7BSDuunAADgPISVclg/BQAA5yGslMP6KQAAOA9h5Qil66fEeip29cR63JWmLQMAAN9j6nIVBibG6YKEWNZPAQDAAQgr1WD9FAAAnIFuIAAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiEFQAA4GiElaOxLLsrAAAgpBFWqrNmjdS/v7Rggd2VAAAQ0ggr1dm4UVq9WpozRzp0yO5qAAAIWYSV6tx6q9SqlfTdd9KiRXZXAwBAyCKsVKdpU+mOO8zx7NmS12tvPQAAhCjCSk1Gj5aaN5e++kpassTuagAACEmElZrExEjjxpnjWbOYGQQAgA0IK0dz++1Ss2bSf/4jvf223dUAABByCCtH06KF6Q6SaF0BAMAGhJXauOMOqXFjaf16acUKu6sBACCkEFZqo3VrM5VZku67j9YVAAD8iLBSW3feKUVGSh9/bBaLAwAAfkFYqa0TTpBuucUcz5xpby0AAIQQwkpd3HWX1KiRtHKlaWEBAAA+R1ipixNPlIYPN8f33WdrKQAAhArCSl1NnSqFh0vLl0uffmp3NQAABD3CSl3Fx0t//rM5pnUFAACfI6zUx7RpUliY9M470saNdlcDAEBQc0RYefLJJxUfHy+3260ePXroo48+srukmnXqJF1/vTlmZhAAAD5le1h59dVXNX78eN19993KyMhQ3759ddFFF2n79u12l1az6dNN68qyZdKmTXZXAwBA0LI9rDz00EO6+eabdcstt+jUU0/VI488ovbt22v+/Pl2l1azzp0Pt66kptpbCwAAQczWsFJcXKyNGzdqwIABFc4PGDBAn3zyiU1V1UH51pWMDLurAQAgKNkaVvbs2aOSkhK1adOmwvk2bdpo9+7dVf7MgQMHVFBQUOFlm86dpeuuM8e0rgAA4BO2dwNJksvlqvDesqxK50rNmTNHHo+n7NW+fXt/lFi90taVpUtpXQEAwAdsDSutWrVSeHh4pVaU3NzcSq0tpaZOnar8/Pyy144dO/xRavW6dDncusLMIAAAGpytYSUyMlI9evTQihUrKpxfsWKFzj777Cp/JioqSjExMRVetittXXnrLVpXAABoYLZ3A91xxx167rnn9Pe//11ffvmlJkyYoO3bt2vUqFF2l1Z75VtXZsywtRQAAIJNhN0FXHPNNdq7d69mzpypnJwcJSYm6t1339VJJ51kd2l1c8890iuvmJlBGzdKPXrYXREAAEHBZVmWZXcRx6KgoEAej0f5+fn2dwn9+c/Syy9LF18svf22vbUAAOBgdfn+tr0bKKjcc4/Zkfmdd9iRGQCABkJYaUinnCINHWqOGbsCAECDIKw0tNLWlfR0ae1au6sBACDgEVYa2sknS8OHm+OUFFtLAQAgGBBWfGH6dCkiQlqxQvroI7urAQAgoBFWfKFDB+nmm83x9OlSYE+4AgDAVoQVX5k+XYqKktaskT74wO5qAAAIWIQVX2nXTho50hzfcw+tKwAA1BNhxZemTpUaN5bWrZPefdfuagAACEiEFV+KjZXGjDHHtK4AAFAvhBVfmzxZatbM7Mb85pt2VwMAQMAhrPhaq1bS+PHm+N57pZISW8sBACDQEFb8YeJEqXlzacsWadEiu6sBACCgEFb8oXlz6c47zXFKinTwoK3lAAAQSAgr/jJunNS6tbRtm/T883ZXAwBAwCCs+EuzZtK0aeb4vvukoiJ76wEAIEAQVvxp5EizWNwPP0hPPWV3NQAABATCij+53Yd3Yr7/fmnfPnvrAQAgABBW/G3YMKlTJ+l//5MeecTuagAAcDzCir81aiTNnGmO582T8vLsrQcAAIcjrNjhmmuk00+X8vOluXPtrgYAAEcjrNghLMyMWZGkv/1N2rnT3noAAHAwwopdBg2S+vQxU5hLu4UAAEAlhBW7uFzSnDnmeMEC6Ztv7K0HAACHIqzYqU8f6eKLzeaG99xjdzUAADgSYcWPSryW1m7bq6Wbd2rttr0q8VrS7Nnmw1dflTIy7C0QAAAHirC7gFCRnpmj1LQs5eQfXmY/zuNWyuAEDbz+eumf/5SmTpXS022sEgAA56FlxQ/SM3OUvHBThaAiSbvzi5S8cJNWDxkjRURI778vffihTVUCAOBMhBUfK/FaSk3LklXFZ6Xnpmz+Rd6RI82bu+6SrKquBgAgNBFWfGx9dl6lFpXyLEk5+UXaNHS01LSptGGD9MYb/isQAACHI6z4WG5h9UGlvJ1RMdKkSebN3XdLBw/6sCoAAAIHYcXHWke7a3/dxInS8cebNVcWLPBxZQAABAbCio/1im+hOI9brmo+d8nMCuoV30KKjpbuvdd8MGOGtG+fn6oEAMC5CCs+Fh7mUsrgBEmqFFhK36cMTlB42G/vbr1VOvlk6ccfpYcf9ludAAA4FWHFDwYmxmn+kO6K9VTsEor1uDV/SHcNTIw7fDIyUpo1yxz/9a9Sbq4fKwUAwHlclhXY82QLCgrk8XiUn5+vmJgYu8upUYnX0vrsPOUWFql1tOn6KWtRKc/rlXr1kjZulMaMkR57zP/FAgDgQ3X5/iasONXKldK555rF4rKypFNOsbsiAAAaTF2+v+kGcqpzzpEGDZIOHZKmTbO7GgAAbENYcbK5c6WwMLNI3Lp1dlcDAIAtCCtOlpgoDR9ujidPZhl+AEBIIqw4XWqq1Lix9NFH0rJldlcDAIDfEVacrl07acIEc3zXXSzDDwAIOYSVQHDXXWYZ/q1bpWeesbsaAAD8irASCGJiTHeQZJbhz8+3tRwAAPyJsBIobrlF6tJF2rNH+stf7K4GAAC/IawEikaNzPL7ktkzaPt2e+sBAMBPCCuB5I9/lPr3lw4cYKE4AEDIIKwEEpdLmjfPHP/jH9Jnn9lbDwAAfkBYCTQ9ekhDh5rjO+5goTgAQNAjrASi++83C8V9/LG0ZInd1QAA4FOElUDUrp10553mePJkM4YFAIAgRVgJVHfeKcXFSd99Jz32mN3VAADgM4SVQNWsmekOkqT77pP+9z976wEAwEcIK4Hsz3+WunWTCgrMyrYAAAQhwkogCwuTHnrIHD/9tLRli731AADgA4SVQNe/v3T55VJJCVOZAQBBibASDB54QIqMlJYvl955x+5qAABoUISVYNCxozR+vDm+4w6puNjWcgAAaEiElWBx991SmzbSN99Ijz9udzUAADQYwkqwiImRZs82x6mpUm6uvfUAANBACCvBZPjww1OZ77nH7moAAGgQhJVgEh4uPfqoOX72WWnzZlvLAQCgIRBWgk3fvtLVV5spzOPGMZUZABDwCCvB6IEHzK7MH30kvfaa3dUAAHBMCCvB6MQTpSlTzPGkSdL+/fbWAwDAMSCsBKs775ROOkn64Qdp7ly7qwEAoN4IK8GqcePD+wb99a9Sdra99QAAUE+ElWB2+eXSuedKBw6Y7iAAAAIQYSUIlHgtrd22V0s379TabXtV4v1tBpDLZaYyh4dLS5ZIK1bYWygAAPUQYXcBODbpmTlKTctSTn5R2bk4j1spgxM0MDFOSkyUxo6VHnnE/PPzz82mhwAABAhaVgJYemaOkhduqhBUJGl3fpGSF25SemaOOTFjhtS6tbR16+FF4wAACBC2hpUOHTrI5XJVeE0pnXKLGpV4LaWmZamqJd9Kz6WmZZkuIY/HDLKVpJkzpZ07/VUmAADHzPaWlZkzZyonJ6fsNX36dLtLCgjrs/MqtaiUZ0nKyS/S+uw8c2LoUCkpSdq3T5o82T9FAgDQAGwPK9HR0YqNjS17NWvWzO6SAkJuYfVBpcrrwsKkxx83g27/+U9pzRofVgcAQMOxPazMnTtXLVu2VNeuXTV79mwVFxfXeP2BAwdUUFBQ4RWKWke7635d9+7SyJHmePRo6eBBH1QGAEDDsjWs3H777Vq0aJFWrlypMWPG6JFHHtFtt91W48/MmTNHHo+n7NW+fXs/VessveJbKM7jlquaz10ys4J6xbeo+MGsWVLLllJmpmlpAQDA4VyW1bDb8s6YMUOpqak1XvPZZ5+pZ8+elc4vXrxYV155pfbs2aOWLVtW+bMHDhzQgQMHyt4XFBSoffv2ys/PV0xMzLEVH2BKZwNJqjDQtjTAzB/S3UxfPtJzz0kjRkjR0dJXX0lt2/q8VgAAyisoKJDH46nV93eDh5U9e/Zoz549NV7ToUMHud2VuzF27typdu3aad26dTrrrLNq9e+ryx82GB11nZWqeL3S//2ftG6ddO210iuv+KlaAACMunx/N/iicK1atVKrVq3q9bMZGRmSpLi4ar5kUcnAxDhdkBCr9dl5yi0sUuto0/UTHlZdB5HMYNsnnpDOPFNatMi0spx7rv+KBgCgDmxbwXbt2rVat26dzjnnHHk8Hn322WeaMGGCLrnkEp144ol2lRWQwsNcSupYdbdZtbp3l5KTTWgZPVr6z39Y2RYA4Ei2DbCNiorSq6++qv79+yshIUH33nuvRowYoVfokvCfWbPMyrZffSU9/LDd1QAAUKUGH7Pib6E+ZuWYvfSSNGyY1KSJlJUlnXSS3RUBAEJAXb6/bV9nBTYbOlTq10/65Rdp3Di7qwEAoBLCSqhzuaT586VGjaRly6SlS+2uCACACggrkE49VZo0yRyPHWv2DwIAwCEIKzCmT5c6dJB27DA7MwMA4BCEFRhNmhxefv+hh6QvvrC3HgAAfkNYwWEXXyz96U9SSYnZ8NDrtbsiAAAIKzjCo4+aPYPWrpWeecbuagAAIKzgCO3aSbNnm+MpU6ScHHvrAQCEPMIKKrvtNqlnTyk/Xxo/3u5qAAAhjrCCysLDTRdQeLj02mvSu+/aXREAIIQRVlC1bt0Ot6rcdpu0f7+t5QAAQhdhBdWbMUM68UTpv/+VUlLsrgYAEKIIK6hes2ZmKX7J7Mq8YYO99QAAQhJhBTUbNEi69lqz5sqIEdLBg3ZXBAAIMYQVHN2jj0otWkibN5sWFgAA/IiwgqNr3Vp68EFznJIiffutvfUAAEIKYQW1M2yYdN55UlGRWYrfsuyuCAAQIggrqB2XS3r6acntlj78UHr+ebsrAgCECMIKaq9jR2nmTHN8xx3Srl321gMACAmEFdTNhAmHl+K/7Ta6gwAAPkdYQd1EREh//7v559Kl0uuv210RACDIEVZQd7//vTRtmjkeM0bau9feegAAQY2wgvq5+27ptNOk//2PnZkBAD5FWEH9REaa7qCwMGnhQuntt+2uCAAQpAgrqL9evcysIMmsvfLTT/bWAwAISoQVHJuZM6Xf/c5MYy4NLgAANCDCCo5N48ZmgTiXS3rhBem99+yuCAAQZAgrOHZnn314kO2IEWYNFgAAGghhBVUq8Vpau22vlm7eqbXb9qrEe5TF32bNkjp1knbupDsIANCgIuwuAM6Tnpmj1LQs5eQXlZ2L87iVMjhBAxPjqv6hJk3M7KB+/cw/r7hCGjTITxUDAIIZLSuoID0zR8kLN1UIKpK0O79IyQs3KT0zp/of7ttXuv12czxiBLODAAANgrCCMiVeS6lpWaqqw6f0XGpaVs1dQvfff3h20LhxvigTABBiCCsosz47r1KLSnmWpJz8Iq3Pzqv+lzRuLL344uHF4t56q8HrBACEFsIKyuQWVh9U6nRd797S5MnmeORIac+eY6wMABDKCCso0zra3XDXzZghJSZKublScrJkHWU2EQAA1SCsoEyv+BaK87jlquZzl8ysoF7xLY7+y6KiTHdQRIT0xhvSP//ZkKUCAEIIYQVlwsNcShmcIEmVAkvp+5TBCQoPqy7OHKF7dyklxRyPHi3t2NEgdQIAQgthBRUMTIzT/CHdFeup2NUT63Fr/pDu1a+zUp0pU6SzzjKr2t54o+T1NmC1AIBQ4LKswB5MUFBQII/Ho/z8fMXExNhdTtAo8Vpan52n3MIitY42XT+1blE50tdfS127Sr/+Kv3tb9LYsQ1aKwAg8NTl+5uwAv948knTFeR2SxkZUpcudlcEALBRXb6/6QaCfyQnSwMGSEVF0tCh0sGDdlcEAAgQhBX4h8slPf+81KKFtGGDlJpqd0UAgABBWIH/tG0rPf20OZ4zR/r4Y3vrAQAEBMIK/OvKK6Xhw82soKFDzSwhAABqQFiBX5V4La2/PUX7TzhR+v57ecey2SEAoGaEFfhNemaO+sz9UFcv2qKh/ceqxBWmsJdf0uZ5T9ldGgDAwQgr8Iv0zBwlL9xUtqvzpnan6vGkqyVJJ0+fpFXLP7OzPACAgxFW4HMlXkupaVk6ckGfx86+VpvadlbMgf067tYbVVLMdGYAQGWEFfjc+uy8shaV8g6FR2jc4DtVGNlYZ/x3i3ZNvseG6gAATkdYgc/lFlYOKqV+aB6ruy8cLUlq99gD0r//7a+yAAABgrACn2sd7a7x82UJ/bX4tHPk8nqlG26Qfv7ZP4UBAAICYQU+1yu+heI8blW3DaJL0pNXjJd18snSf/8rjRwpBfaWVQCABkRYgc+Fh7mUMjhBkioFltL3d151plyvvCJFREivvSYtWODXGgEAzkVYgV8MTIzT/CHdFeup2CUU63Fr/pDuGpgYJ/XqJd1/v/lg3DhpyxYbKgUAOI3LsgK7vb0uW0zDfiVeS+uz85RbWKTW0W71im+h8LBy7S1erzRokPT++9Jpp0mffSY1bmxfwQAAn6jL93eEn2oCJJkuoaSOLau/ICxMeukl6YwzTMvKhAnSU6xwCwChjG4gOE/r1tLLL0sul9mledEiuysCANiIsAJnOv98ado0czxihPTNN/bWAwCwDWEFzjVjhvSHP0j79klXXSUVVb+4HAAgeBFW4FwREdIrr0jHHy/95z9m/AoAIOQQVuBsbdtKCxea8StPPWXCCwAgpBBW4HwDBhwev3LrrdLWrfbWAwDwK8IKAsOMGVK/fmb8yhVXSPv3210RAMBPCCsIDBERZgpzbKxZfyU5mf2DACBEEFYQOGJjTWAJDzfrsDz7rN0VAQD8gLCCwNKv3+H9g8aOlTZutLceAIDPEVYQeO68U7rkEqm4WLrySikvz+6KAAA+RFhB4HG5pBdflE4+Wfr+e+mGG6SSErurAgD4CGEFgal5c2nJErMjc3q6lJpqd0UAAB8hrCBwnXGG9Mwz5vi++6S0NHvrAQD4BGEFgW3IEGnMGEnSoRuG6F/LPtbabXtV4mVaMwAEC5+GldmzZ+vss89WkyZN1Lx58yqv2b59uwYPHqymTZuqVatWGjdunIqLi31ZFoLM+zdN0ucnnqaIwgK1v/kG3fzESvWZ+6HSM3PsLg0A0AB8GlaKi4t11VVXKTk5ucrPS0pKdPHFF2v//v36+OOPtWjRIi1evFgTJ070ZVkIIumZORr1aqZu+eNk5TY9Tp33bNe8dx7W7p9/VfLCTQQWAAgCPg0rqampmjBhgn7/+99X+fny5cuVlZWlhQsXqlu3bjr//PP14IMP6tlnn1VBQYEvS0MQKPFaSk3LkiUpN7qlRl02TcVhERr09SdKXve6JCk1LYsuIQAIcLaOWVm7dq0SExPVtm3bsnMXXnihDhw4oI3VLPZ14MABFRQUVHghNK3PzlNOflHZ+03tTtW9F4ySJE1a87L6bftMOflFWp/NOiwAEMhsDSu7d+9WmzZtKpw77rjjFBkZqd27d1f5M3PmzJHH4yl7tW/f3h+lwoFyC4sqnVvUdaAWdr1IYbL0t7R5is/bWeV1AIDAUeewMmPGDLlcrhpfGzZsqPXvc7lclc5ZllXleUmaOnWq8vPzy147duyo6x8BQaJ1tLvK86nn36rPTkhQzIH9enbxfYoTA7YBIJBF1PUHxowZo2uvvbbGazp06FCr3xUbG6tPP/20wrmffvpJBw8erNTiUioqKkpRUVG1+v0Ibr3iWyjO49bu/CKVH5VyMLyRki+fqmUvTlCnvB9kTR0tvZ1mNkAEAAScOrestGrVSl26dKnx5XZX/V+8R0pKSlJmZqZycg7P2Fi+fLmioqLUo0ePupaGEBMe5lLK4ARJ0pHtcHubHqdb/zRdJW63XOnvSVOm+L9AAECD8OmYle3bt2vz5s3avn27SkpKtHnzZm3evFn79u2TJA0YMEAJCQkaOnSoMjIy9MEHH2jSpEkaMWKEYmJifFkagsTAxDjNH9JdsZ6KATnW49aYSVcr/MUXzYl588x+QgCAgOOyLMtn8zqHDx+uF6v4gli5cqX69+8vyQSa2267TR9++KEaN26s66+/XvPmzat1V09BQYE8Ho/y8/MJOCGsxGtpfXaecguL1DrarV7xLRQe9lt7y733muX4IyOlVaukpCRbawUA1O3726dhxR8IKzgqr1e68krpzTel1q2lTz+VajmuCgDgG3X5/mZvIAS/sDDppZekrl2l3Fxp8GCJ9XkAIGAQVhAamjUzuzLHxUmZmdJ110klJXZXBQCoBcIKQke7dtKyZVLjxtK770qTJtldEQCgFggrCC09e5ouIUl65BHpySdtLQcAcHSEFYSeK6+UZs82x2PHSu+8Y289AIAaEVYQmqZOlW66ycwUuuYaKSPD7ooAANUgrCA0uVzSU09J550n7d8vXXyxxD5TAOBIhBWErkaNpDfekBISpJwcE1jy8+2uCgBwBMIKQlvz5mZmUGys9MUX0hVXSMXs0gwATkJYAU46SXr7balpU+mDDw6PZQEAOAJhBZCkHj1Ml1B4uPSPf0jTptldEYB6KPFaWrttr5Zu3qm12/aqxBvQO8rgNxF2FwA4xsCB0nPPSTfeKM2daxaRGzPG7qoA1FJ6Zo5S07KUk19Udi7O41bK4AQNTIyzsTIcK1pWgPKGD5dmzjTH48aZ1hYAjpeemaPkhZsqBBVJ2p1fpOSFm5SemWNTZWgIhBXgSNOnSyNHSpYl3XCD9OGHdlcEoAYlXkupaVmqqsOn9FxqWhZdQgGMsAIcyeWSnnhC+tOfzMygyy5j0TjAwdZn51VqUSnPkpSTX6T12Xn+KwoNirACVKV0oG3//lJhoRnP8u23dlcFoAq5hdUHlfpcB+chrADVcbult96SunaVcnOlAQOkXbvsrgrAEVpHuxv0OhzmlNlVzAYCauLxSO+9J/XpI23bZgLL6tVSy5Z2VwbgN73iWyjO49bu/KIqx624JMV63OoV38LfpQU0J82uomUFOJrYWGnFCqltW2nLFumii0zXEABHCA9zKWVwgiQTTMorfZ8yOEHhYUd+iuo4bXYVYQWojfh4E1hatpQ++0y69FKpiP5vwCkGJsZp/pDuivVU7OqJ9bg1f0h31lmpAyfOrqIbCKithAQpPV0691xp5Urp6qulxYvNhogAbDcwMU4XJMRqfXaecguL1DradP3QolI3dZldldTRP13itKwAddGzp5SWZgbfpqWZdVgOHbK7KgC/CQ9zKaljS13a9QQldWxJUKkHJ86uIqwAddWvn7RkiWlRef116eab2fgQQNBw4uwqwgpQHxddJL36qlmP5aWXpNtuMyveAkCAK51dVV2blEtmVpA/Z1cRVoD6uvxy6eWXzYq3Tz8tjR9PYAEQ8Jw4u4qwAhyL664zOzVL0t/+Jk2cSGABEPCcNruK2UDAsbrpJjPIduRI6eGHTdfQX/9qWlwAIEA5aXYVYQVoCLfeagbZJidL8+aZwDJnDoEFQEArnV1lN7qBgIYyapT0+OPmeO5caepUuoQAoAEQVoCGNHq0GbsimcBy550EFgA4RoQVoKGNHSs98YQ5fvBBZgkBwDEirAC+cNttZjqzZFpaRo9m4TgAqCfCCuArt94qLVhgBtnOny+NGCGVlNhdFQAEHMIK4Es33SS9+KIUFib9/e9mL6GDB+2uCgACCmEF8LWhQ83S/I0amX9ecYVU5L8NwAAg0BFWAH+48krprbfKdmv++fyBevuTr7V2216VeBl8CwA1IawA/jJokNY/8bJ+iXSr+b9X64QrL1Hy35arz9wPlZ6ZY3d1AOBYhBXAT9Izc3TN1411/TWz9ZM7Wt1ytuq1f0yR9cMPSl64icACANUgrAB+UOK1lJqWJUvS5radddUNc5XTrKV+t3e73lg4WfF5O5WalkWXEABUgbAC+MH67Dzl5B8eVPttqxN15ZAH9N1xbdWuIFev/WOyWn31hdZn59lYJQA4E2EF8IPcwsqzf3Z6WuuqG/6qL9p0VKtf8rXolalS+ns2VAcAzkZYAfygdbS7yvN7mzbXtdfN0ZoO3dT0YJHOGn+jWZcFAFCGsAL4Qa/4ForzuOWq4rP9UU10y5X36r2u5yvs0CFp+HDp/vvZTwgAfkNYAfwgPMyllMEJklQpsLgkHQxvJNeLL0h33WVO3n23WZ6f1W4BgLAC+MvAxDjNH9JdsZ6KXUKxHrfmD+mugaefIP3lL9Ljj5vl+RcskC6+WMrPt6li3ynxWlq7ba+Wbt7JwngAjsplWYHd1lxQUCCPx6P8/HzFxMTYXQ5wVCVeS+uz85RbWKTW0W71im+h8LAj2lveflu69lpp/34pMVF65x3pxBPtKbiBpWfmKDUtq8LsqDiPWymDEzQwMc7GyuBPtfr/AYJaXb6/CSuAU23aJP3xj1JOjtSmjbR0qXTWWXZXdUzSM3OUvHCTjvxLp/Qrav6Q7gSWEEBghVS372+6gQCn6t5d+vRT6fTTpR9/lPr1k155xe6q6q38wnhHKj3HwnjBrzSwlg8qkrQ7v4iVnFEtwgrgZO3bS//+t3TJJdKBA9L110v33CN5vXZXVmdHLox3JEtSTn4RC+MFsYAIrPv2BeU4sUBHWAGcrlkzackSafJk837WLLOLc2GhvXXVUVUL4x3LdQg8jg+s27ZJSUnSkCEB+R8EwYywAlTDUTNWwsOluXOl55+XIiOlN980f6l++619NdVRdQvj1fc6BB5HB9b335d69pQyM6UNG6Tt2/1fA6oVYXcBgBM5dgDg8OFS587SFVdIW7ZIZ54pLVokXXihfTXVUunCeLvzi6rsBnDJTOPuFd/C36XBTxwZWC1LeuABaepU05py1lmmJbNtW//VgKOiZQU4guMHACYlmf/y691b+vlnadAgsz6Lwyf2HW1hPElKGZzA9NUgVtNKzpJ5DuL8GVj37ZOuu84sxuj1SjfdJK1eTVBxIMIKUE5ADACUzF+mq1ZJN99s/pKdOlW6/HLHDww86sJ4TFsNao4KrF99ZVpRXn1VioiQnnhCeu45KSrK9/9u1BnrrADlrN22V9c9u+6o170yoreSOrb0Q0VHYVnmL9gxY6TiYqlTJ2nxYjPd2cFYECy02d7N+vrrphVl3z4pLk567TWpTx/f/3tRQV2+vxmzApTj6AGAVXG5zB5C3bqZGULffmu6h5580oxvcajwMFeDhz0CUOAYmBinCxJi/f+/V3Gx6fJ55BHzvn9/M+arTRvf/ntxzAgrQDmOHABYGz17Shs3mimX6enSjTdKK1eapu1mzeyuzuds/y911JkvAmuNsrPNFhbr15v3d91llgGI4GswEDBmBSjHcQMA66JlS7OH0OzZZiPEl14ys4W++MLuynzK8QOiYb8lS0zr4/r10nHHma0r/vIXgkoAIawA5ThqAGB9hIVJ06aZVpW2bc0gwl69pPnzHT9bqD4CZkA07PHrr2Y81xVXmMHnvXtLGRlmRWgEFMIKcISgmLHyhz9ImzdLAwdKRUXSbbdJl10m7dljd2UNyvErosI+n39uWhafeMK8v/NOac0a6aST7K0L9UIbGFAF2wYANqTjjzfdQo8+Kk2ZIi1bZmYJvfyydN55dlfXIAJuQDR8z7Kkxx4z21McOGAGz774YkAsnIjq0bICVKN0AOClXU9QUseWgRVUSoWFSRMmmN2bu3SRcnKk888353791e7qjlnADoiGb+zcKV10kXT77SaoXHyxaWEhqAQ8wgoQCrp2NbOFRo407x95ROre3ayEG8ACekA0Go5lSa+8IiUmmj1+3G7TupKWJrVubXd1aACEFSBUNGkiPfWU6RqKjTWDb3v3lmbMMOtPBKCAHxCNY7dnj5mSfP31ZvuJnj3NINoxY8w6RAgKhBUg1AwaZHaWveoqqaRESk01M4Y2bbK7snoJigHRqDvLMivRJiSYFWgjIsyz/MknpssTQYXl9oFQZVlm9c6xY6W9e6XwcDMQ9557AnJ/FFawDSG7d0ujR5v1UyTptNOkF14wrSoIGHX5/iasAKHuxx9NYHn9dfP+1FOlZ55hrxQ4j2WZUDJpkpSXZ1pTpk0zrwAM2KGuLt/fdAMBoa5NG9OM/sYbZjDil19Kffuawbg//WR3dYCxdat07rlmA8K8PLMi7Wefma4fgkrQI6wAMK64wgSVW24x7595Rtapp+rrR57R0owftHbbXlaChf8VFUkzZ5o1glatkho3lv76VzMdv2tXu6uDn9ANBKCyNWu0b9hNavb9NknSJyeernsvGKX9HX/H5oDwn3fflcaNk7aZ51AXXWR2FO/Qwday0DAc0w00e/ZsnX322WrSpImaN29e5TUul6vS66mnnvJlWQCOIr3FKepxxTw92OcGFUVE6uztn+u958fqxree1KQFH7M5IHzr++/N9hAXX2yCStu2ZjD4O+8QVEKUT8NKcXGxrrrqKiUnJ9d43fPPP6+cnJyy17Bhw3xZFoAalG4OeCCikR77v+t0/i3ztfyU3mrkLdGt65do5bO3amPKgyo5VGJ3qQg2hYXS3XebqcdLl5oBtJMmmTWBrrmGdVNCmE/3BkpNTZUkvfDCCzVe17x5c8XGxvqyFAC1dOTmgD942ujWP01X/22fKeWDZxT/U47uXvKg9p2xXM3mP242TQSOhddr9u+ZNs1MS5akc86RHn/crKOCkOeIAbZjxoxRq1atdOaZZ+qpp56S1+ut9toDBw6ooKCgwgtAw6lu079VHc/UgJuf1KxzblJBZBM1y/pC6tfPNNd/+aV/i0TwWL5c6tHDzPLZvVvq2FF6803pgw8IKihje1i577779Prrr+tf//qXrr32Wk2cOFH3339/tdfPmTNHHo+n7NW+fXs/VgsEv5o2/TsY3kjP9fqT+o98VruvG2Y2Sly61OzJcuut0q5dfqwUAW3jRumCC8wmg5s3SzEx0rx50pYtJgDT5YNy6hxWZsyYUeWg2PKvDXXYHG369OlKSkpS165dNXHiRM2cOVMPPPBAtddPnTpV+fn5Za8dO3bU9Y8AoAa12RwwKq6Njl/4vFm2/9JLTTP+s89KnTpJkyeb/VqAqnz5pRl/0rOn9K9/SY0aSePHm4G0EyeyZgqqVOcxK2PGjNG1115b4zUdjmG0du/evVVQUKAff/xRbdq0qfR5VFSUoniYgQoacqn50s0BkxdukktS+bUNKm0OeOqp0ltvSf/+twkpn3wiPfCANH+++QKaOFGqZiYgQsy335oF3P75TxNuXS7phhvMGirx8XZXB4erc1hp1aqVWrVq5YtaJEkZGRlyu93VTnUGUFF6Zo5S07IqDIqN87iPaT2U0s0Bj/y9sdX93v/7P+njj6X33jN7C23aJM2aJT32mFnK//bbJR/+vQEH27pVmjNHWrjQbJwpmW6e1FSz0BtQCz5dFG779u3Ky8vTsmXL9MADD+ijjz6SJHXq1EnNmjVTWlqadu/eraSkJDVu3FgrV67UxIkTNXz4cD366KO1+newKBxCWXpmjpIXbtKR/ycubQE51l2H69ViY1mmteXee003kSQ1bSqNGmVaWuJYUC4k/Oc/0v33mz2nSr9mLr7YhJQePeytDY7gmI0Mhw8frhdffLHS+ZUrV6p///5KT0/X1KlT9e2338rr9erkk0/WLbfcotGjRysionaNPoQVhKoSr6U+cz+s0PJRnkumJeTju861Z/dhr9cMvp01y7S0SFJkpDR0qAktp57q/5rgW5YlffihGSibnn74/KWXmvVTzjzTvtrgOI4JK/5AWEGoWrttr657dt1Rr3tlRG8ldWzph4qqYVnmi2vWLDOmpdQf/2hCS79+zPxwuKO2sBUXmxaUefPMzB7JzBS76ioTUn7/e1vqhrPV5fvbp4vCAfCd6tZDqe91PuNymT1dLrrIhJV580w30dtvm1diojRmjDRkiOkugqPUOCaqhSU9/bR5/fij+bBJE7NmyvjxZs0UoAHYvs4KgPqpaT2U+lznF2efLS1ZYpZPHzXKfLFlZprjE04wX3BbtthdJX5TOiaqfFBxWV6d/J91OnTVNfKedJKZzfPjj2b/nlmzpB07zMBqggoaEN1AQIAqHbOyO7+o0gBbyQFjVmrj55+l55+Xnnji8M66kpSUJI0YYboRmjWzrbxQduSYqOP35emqL/6lqz9foQ4/H97I0urbV64xY6TLLzdrpgC1xJgVIESU/pevVPV6KMc6G8hvvF7p/ffNwnLLlh2e4tq0qXTFFdKf/yz17y+Fh9taZihZu22vbnxylQZ8s05/ylypvt9nKNwyW6EURDbR0tP665UzBuqe6dfbOyYKAYuwAoQQX6yzYqvdu6UXXpAWLDALiZVq1066+mqz+umZZzIo11eKi6UVK7T96Rd13PtvK7r417KPNrbtole6DtQ7nfvo10jTvfjotV11adcT7KoWAYywAoSYhlzB1jEsS1q3TnrpJWnRItNlVOqkk0xwuewy6ayzaHE5Vr/+ajYOXLzYDH4ud6+3e9rozdPO1Zun9df3LSqHEttnmyFgEVYABJeiIrM67muvSWlp0v79hz9r3dpMg77kEuncc6XoaPvqDCS7d5t7umyZ2fn4l18OfxYbK++VVyq5uKOWN+8kq4pWrIAYEwVHI6wACF6//CK9+66ZVfTuu1J+/uHPGjUyS/9feKE0YIB0xhmObHXxVUtYjb/311+ltWvN2KD33zcrzJbXvr0JfFdfbe5heHjwjImCIxFWAISG4mLpo4/MSrnvvCN9913Fzz0e6Q9/MINz+/Y14SUy0pZSS/lqjNGRv7fpgV90bsH3ur3RLnX6cpPpUisurvhDPXpIgwebFWbPOKPKcUBBNyYKjkFYARCavv3WtBqkp0urV0uFhRU/d7vNF3Tv3lLPnlK3btIpp5jVVv3AV3s5vb/pez32t6VK+PE7dc3Zqm67tup3e7aXzd4pExcnXXCBaXm64ALp+ONr9fuDckwUbEdYAYBDh8zS76tWmdfatVJeXuXrmjY1rQqnnSYlJJg9i0491cw+asAQ0yB7OR08aNaj+fLLspf1+ecqydyiCG9Jpct/iDleG9qdpi87d9PkOaMU/rtTmEUFx2C5fQCIiDCtJz17SpMmmdlF33xjQsu6dWZzxc8/N4N1P/mk4r5Fkukuio+XTj7Z/POEE8wqrW3bmhaKli3NKyqqVuWsz86rNqjIsuQ+eEDa8T9tSVup0yN+lXbulHbtMivCfvedlJ1tjr0VW0tcMn+R/+xupi1tTtYXbTop44Qu2tS2i/7XrEXZdf0jWiqJoIIARVgBEBpcLul3vzOvYcPMuUOHpK+/Ni0wWVmmtSIry3QnFRdLW7eaV02aNpWaNzcr7TZtal5utwlLpS/LUnzePr2w4yeFe0vU+NABNTlYpMYHi9S0uEieon2KKjloft/8o/w5mjaVunQxr1NP1adN4jThG5d2RR9fY6uJ7XtEAceAsAIgdEVEmK6fhISK5w8dMq0Y2dmmVeP7700rx65dpsXjxx9Nl1JJiWmZKT+Vuhqxv71qUhwWIbVqqcj27Q634JxwgmndKW3hiY2tEEq82/ZqVy1233bUHlFAHRFWAOBIEREmGMTHm7VbquL1SgUF0t69ZhG1/fulffvMq7jYBJ7SV1iYvK4wzV7+jfKKSvRLRJR+bRSlXyLd+qVRY/3cuJny3dGKadVcH085T6rD4NVe8S0U53EfdY+oXvEtqvgUCAyEFQCoj7Aw0/3TvHntLpd0ZlLN65Y8eMlpdZ5lEx7mUsrgBCUv3CRXNb83ZXACs3cQ0PwzXw8AoIGJcZo/pLtiPRW7ZGI97mNaYM1XvxdwCqYuA4Cf2bKCLeAwTF0GAAcLD3P5ZPM/X/1ewG50AwEAAEcjrAAAAEcjrAAAAEcjrAAAAEcjrAAAAEcjrAAAAEcjrAAAAEcjrAAAAEcjrAAAAEcjrAAAAEcjrAAAAEcjrAAAAEdjI0MAqAE7GQP2I6wAQDXSM3OUmpalnPyisnNxHrdSBidoYGKcjZUBoYVuIACoQnpmjpIXbqoQVCRpd36RkhduUnpmjk2VAaGHsAIARyjxWkpNy5JVxWel51LTslTireoKAA2NsAIAR1ifnVepRaU8S1JOfpHWZ+f5rygghBFWAOAIuYXVB5X6XAfg2BBWAOAIraPdDXodgGNDWAGAI/SKb6E4j1vVTVB2ycwK6hXfwp9lASGLsAIARwgPcyllcIIkVQospe9TBiew3grgJ4QVAKjCwMQ4zR/SXbGeil09sR635g/pzjorgB+xKBwAVGNgYpwuSIhlBVvAZoQVAKhBeJhLSR1b2l0GENLoBgIAAI5GWAEAAI5GWAEAAI5GWAEAAI5GWAEAAI5GWAEAAI5GWAEAAI5GWAEAAI5GWAEAAI4W8CvYWpYlSSooKLC5EgAAUFul39ul3+M1CfiwUlhYKElq3769zZUAAIC6KiwslMfjqfEal1WbSONgXq9XW7duVUJCgnbs2KGYmBi7S7JdQUGB2rdvz/0oh3tSEfejMu5JRdyPirgflR3rPbEsS4WFhWrbtq3CwmoelRLwLSthYWE64YQTJEkxMTE8ROVwPyrjnlTE/aiMe1IR96Mi7kdlx3JPjtaiUooBtgAAwNEIKwAAwNGCIqxERUUpJSVFUVFRdpfiCNyPyrgnFXE/KuOeVMT9qIj7UZk/70nAD7AFAADBLShaVgAAQPAirAAAAEcjrAAAAEcjrAAAAEcLmLDy5JNPKj4+Xm63Wz169NBHH31U4/WrV69Wjx495Ha7dfLJJ+upp57yU6X+UZf7sWrVKrlcrkqvr776yo8V+86aNWs0ePBgtW3bVi6XS2+99dZRfybYn4+63pNgf0bmzJmjM888U9HR0WrdurUuu+wybd269ag/F6zPSX3uRzA/I/Pnz9fpp59etrhZUlKS3nvvvRp/JlifjVJ1vSe+fj4CIqy8+uqrGj9+vO6++25lZGSob9++uuiii7R9+/Yqr8/OztagQYPUt29fZWRkaNq0aRo3bpwWL17s58p9o673o9TWrVuVk5NT9jrllFP8VLFv7d+/X2eccYYef/zxWl0f7M+HVPd7UipYn5HVq1dr9OjRWrdunVasWKFDhw5pwIAB2r9/f7U/E8zPSX3uR6lgfEbatWunv/zlL9qwYYM2bNigc889V5deeqm2bNlS5fXB/GyUqus9KeWz58MKAL169bJGjRpV4VyXLl2sKVOmVHn95MmTrS5dulQ4N3LkSKt3794+q9Gf6no/Vq5caUmyfvrpJz9UZy9J1ptvvlnjNcH+fBypNvcklJ4Ry7Ks3NxcS5K1evXqaq8JpeekNvcj1J6R4447znruueeq/CyUno3yaronvn4+HN+yUlxcrI0bN2rAgAEVzg8YMECffPJJlT+zdu3aStdfeOGF2rBhgw4ePOizWv2hPvejVLdu3RQXF6fzzjtPK1eu9GWZjhbMz8exCpVnJD8/X5LUokWLaq8JpeekNvejVLA/IyUlJVq0aJH279+vpKSkKq8JpWdDqt09KeWr58PxYWXPnj0qKSlRmzZtKpxv06aNdu/eXeXP7N69u8rrDx06pD179visVn+oz/2Ii4vTM888o8WLF2vJkiXq3LmzzjvvPK1Zs8YfJTtOMD8f9RVKz4hlWbrjjjvUp08fJSYmVntdqDwntb0fwf6MfPHFF2rWrJmioqI0atQovfnmm0pISKjy2lB5NupyT3z9fATMrssul6vCe8uyKp072vVVnQ9UdbkfnTt3VufOncveJyUlaceOHZo3b57+8Ic/+LROpwr256OuQukZGTNmjD7//HN9/PHHR702FJ6T2t6PYH9GOnfurM2bN+vnn3/W4sWLNWzYMK1evbraL+dQeDbqck98/Xw4vmWlVatWCg8Pr9RqkJubWynZloqNja3y+oiICLVs2dJntfpDfe5HVXr37q1vvvmmocsLCMH8fDSkYHxGxo4dq2XLlmnlypVq165djdeGwnNSl/tRlWB6RiIjI9WpUyf17NlTc+bM0RlnnKFHH320ymtD4dmQ6nZPqtKQz4fjw0pkZKR69OihFStWVDi/YsUKnX322VX+TFJSUqXrly9frp49e6pRo0Y+q9Uf6nM/qpKRkaG4uLiGLi8gBPPz0ZCC6RmxLEtjxozRkiVL9OGHHyo+Pv6oPxPMz0l97kdVgukZOZJlWTpw4ECVnwXzs1GTmu5JVRr0+fDJsN0GtmjRIqtRo0bWggULrKysLGv8+PFW06ZNre+//96yLMuaMmWKNXTo0LLrv/vuO6tJkybWhAkTrKysLGvBggVWo0aNrDfeeMOuP0KDquv9ePjhh60333zT+vrrr63MzExrypQpliRr8eLFdv0RGlRhYaGVkZFhZWRkWJKshx56yMrIyLD++9//WpYVes+HZdX9ngT7M5KcnGx5PB5r1apVVk5OTtnrl19+KbsmlJ6T+tyPYH5Gpk6daq1Zs8bKzs62Pv/8c2vatGlWWFiYtXz5csuyQuvZKFXXe+Lr5yMgwoplWdYTTzxhnXTSSVZkZKTVvXv3ClPshg0bZvXr16/C9atWrbK6detmRUZGWh06dLDmz5/v54p9qy73Y+7cuVbHjh0tt9ttHXfccVafPn2sd955x4aqfaN0ytyRr2HDhlmWFZrPR13vSbA/I1XdC0nW888/X3ZNKD0n9bkfwfyM3HTTTWV/nx5//PHWeeedV/albFmh9WyUqus98fXz4bKs30YFAQAAOJDjx6wAAIDQRlgBAACORlgBAACORlgBAACORlgBAACORlgBAACORlgBAACORlgBAACORlgBAACORlgBAACORlgBAACORlgBAACO9v/Ae79vqwIziQAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "a_opt, b_opt, c_opt = popt\n",
    "x_model = np.linspace(min(x_data), max(y_data), 100)\n",
    "y_model = model_f(x_model, a_opt, b_opt, c_opt) \n",
    " \n",
    "plt.scatter(x_data, y_data)\n",
    "plt.plot(x_model, y_model, color='r')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "0a6d78fa-1a21-460c-92d8-2b62beb3a98f",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfoAAAGiCAYAAAAPyATTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAApEElEQVR4nO3dfXRV1Z3/8c/l6YYquQIxDwx5whkjgg+YKAkDIjIrEFatVurIqitCxzJDa7SSoQzRzghOO1mO1GFZFcoUYRQ7Q2eCth0oJWthwFkkhWD4aQdIfQCTASJCMUGUhOTu3x80t73NA0nOPTf3nv1+zdprOCd737t3z2q/97v3Puf4jDFGAADAk4YMdgcAAIB7CPQAAHgYgR4AAA8j0AMA4GEEegAAPIxADwCAhxHoAQDwMAI9AAAeRqAHAMDDCPQAAHiYq4H+7NmzKi4uViAQUCAQUHFxsT755JNe2yxatEg+ny+s5Ofnu9lNAAA8a5ibH/7Vr35V//d//6cdO3ZIkv76r/9axcXF+vnPf95ru7lz52rjxo2h4xEjRrjZTQAAPMu1QH/48GHt2LFDNTU1mjp1qiTpX//1X1VQUKD6+nrl5OT02Nbv9ys1NdWtrgEAYA3XAn11dbUCgUAoyEtSfn6+AoGA9u7d22ugr6qqUnJysq666irNnDlT3/ve95ScnNxt3dbWVrW2toaOg8Ggfvvb32rs2LHy+XyRGxAAICqMMTp37pzGjRunIUPcW2G+cOGC2traHH/OiBEjlJCQEIEeucO1QN/U1NRtcE5OTlZTU1OP7YqKinTfffcpMzNTR48e1d///d/rzjvv1IEDB+T3+7vULy8v16pVqyLadwDA4GtsbNT48eNd+ewLFy4oO/NKNZ3qcPxZqampOnr0aMwG+34H+pUrV142sO7fv1+Sus2ojTG9Ztr3339/6N+TJ09WXl6eMjMztW3bNt17771d6peVlam0tDR03NzcrIyMDP28epyuuJKbCrzuaz/9xmB3AVG08C/eGOwuIApaz7fr6dlVGjVqlGvf0dbWpqZTHTp6IFOJowYeK1rOBZWd+6Ha2tq8E+hLSkq0YMGCXutkZWXp7bff1kcffdTlbx9//LFSUlL6/H1paWnKzMzUu+++2+3f/X5/t5n+FVcO0ZUOLh7iw5AY/S8W3JFw5fDB7gKiKBrLr4mjhjgK9PGg34E+KSlJSUlJl61XUFCg5uZm7du3T7fddpsk6Ve/+pWam5s1bdq0Pn/fmTNn1NjYqLS0tP52FQCAXnWYoDqMs/axzrWfMRMnTtTcuXO1ePFi1dTUqKamRosXL9YXv/jFsI141113nV577TVJ0qeffqply5apurpax44dU1VVle666y4lJSXpy1/+sltdBQBYKijjuMQ6V+crXn31Vd1www0qLCxUYWGhbrzxRr3yyithderr69Xc3CxJGjp0qN555x3dfffduvbaa7Vw4UJde+21qq6udnWtBgBgp2AE/i/WufrAnDFjxmjz5s291jHm97+GRo4cqV/+8pdudgkAAKu4GugBAIhlHcaowwx8+t1J22gh0AMArOV0nd36NXoAADC4yOgBANYKyqjD4xk9gR4AYC2m7gEAQFwjowcAWItd9wAAeFjwd8VJ+1jH1D0AAB5GRg8AsFaHw133TtpGC4EeAGCtDiOHb6+LXF/cQqAHAFiLNXoAABDXyOgBANYKyqcO+Ry1j3UEegCAtYLmUnHSPtYxdQ8AgIeR0QMArNXhcOreSdtoIdADAKxlQ6Bn6h4AAA8jowcAWCtofAoaB7vuHbSNFgI9AMBaTN0DAIC4RkYPALBWh4aow0HO2xHBvriFQA8AsJZxuEZvWKMHACB2sUYPAADiGhk9AMBaHWaIOoyDNfo4eNY9gR4AYK2gfAo6mNwOKvYjPVP3AAB4GBk9AMBaNmzGI9ADAKzlfI2eqXsAADCIyOgBANa6tBnPwUttmLoHACB2BR0+Apdd9wAAYFCR0QMArGXDZjwCPQDAWkEN8fwDcwj0AABrdRifOhy8gc5J22hhjR4AAA8jowcAWKvD4a77DqbuAQCIXUEzREEHm/GCcbAZj6l7AAA8jIweAGAtpu4BAPCwoJztnA9GriuuYeoeAAAPI6MHAFjL+QNzYj9fJtADAKzl/BG4sR/oY7+HAABY6NixY3rooYeUnZ2tkSNH6pprrtGTTz6ptra2fn0OGT0AwFqx/D76I0eOKBgM6oc//KH+9E//VL/+9a+1ePFinT9/XqtXr+7z5xDoAQDWiuWp+7lz52ru3Lmh4wkTJqi+vl5r164l0AMA0BfO76O/1LalpSXsvN/vl9/vd9S37jQ3N2vMmDH9ahOVNfoXX3xR2dnZSkhIUG5urt58881e6+/evVu5ublKSEjQhAkTtG7dumh0EwCAAUlPT1cgEAiV8vLyiH/H+++/rx/84AdasmRJv9q5Hui3bNmixx57TE888YTq6uo0Y8YMFRUVqaGhodv6R48e1bx58zRjxgzV1dXp8ccf16OPPqqKigq3uwoAsEzQ+BwXSWpsbFRzc3OolJWV9fidK1eulM/n67XU1taGtTlx4oTmzp2r++67T1//+tf7NUbXp+6fffZZPfTQQ6GOrVmzRr/85S+1du3abn/xrFu3ThkZGVqzZo0kaeLEiaqtrdXq1as1f/78LvVbW1vV2toaOv7j6RMAAHoSdDh133kffWJiohITE/vUpqSkRAsWLOi1TlZWVujfJ06c0KxZs1RQUKD169f3u4+uBvq2tjYdOHBAK1asCDtfWFiovXv3dtumurpahYWFYefmzJmjDRs26OLFixo+fHjY38rLy7Vq1arIdhwAAJckJSUpKSmpT3WPHz+uWbNmKTc3Vxs3btSQIf3/UeLq1P3p06fV0dGhlJSUsPMpKSlqamrqtk1TU1O39dvb23X69Oku9cvKysKmSxobGyM3AACAp3W+ptZJccuJEyd0xx13KD09XatXr9bHH3+spqamHuNnT6Ky697nC7/P0BjT5dzl6nd3XnJvZyMAwPs65FOHg3vhnbS9nJ07d+q9997Te++9p/Hjx4f9rTMu9oWrGX1SUpKGDh3a5dfHqVOnumTtnVJTU7utP2zYMI0dO9a1vgIAEEsWLVokY0y3pT9cDfQjRoxQbm6uKisrw85XVlZq2rRp3bYpKCjoUn/nzp3Ky8vrsj4PAIATsTx1Hymu97C0tFQ/+tGP9NJLL+nw4cNaunSpGhoaQvcBlpWV6cEHHwzVX7JkiT788EOVlpbq8OHDeumll7RhwwYtW7bM7a4CACzTod9P3w+sxD7X1+jvv/9+nTlzRk899ZROnjypyZMna/v27crMzJQknTx5Muye+uzsbG3fvl1Lly7VCy+8oHHjxum5557r9tY6AADQu6hsxvvmN7+pb37zm93+bdOmTV3OzZw5U2+99ZbLvQIA2M7p9Hs8TN3zrHsAgLVi+aU2kUKgBwBYyzh8Ta1x8fa6SIn9nyIAAGDAyOgBANZi6h4AAA/7wzfQDbR9rIv9nyIAAGDAyOgBANbqcPiaWidto4VADwCwFlP3AAAgrpHRAwCsFdQQBR3kvE7aRguBHgBgrQ7jU4eD6XcnbaMl9n+KAACAASOjBwBYy4bNeAR6AIC1jMO31xmejAcAQOzqkE8dDl5M46RttMT+TxEAADBgZPQAAGsFjbN19qCJYGdcQqAHAFgr6HCN3knbaIn9HgIAgAEjowcAWCson4IONtQ5aRstBHoAgLV4Mh4AAIhrZPQAAGvZsBmPQA8AsFZQDh+BGwdr9LH/UwQAAAwYGT0AwFrG4a57EwcZPYEeAGAt3l4HAICH2bAZL/Z7CAAABoyMHgBgLabuAQDwMBsegcvUPQAAHkZGDwCwFlP3AAB4mA2Bnql7AAA8jIweAGAtGzJ6Aj0AwFo2BHqm7gEA8DAyegCAtYyc3QtvItcV1xDoAQDWsmHqnkAPALCWDYGeNXoAADyMjB4AYC0bMnoCPQDAWjYEeqbuAQDwMDJ6AIC1jPHJOMjKnbSNFgI9AMBavI8eAADENTJ6AIC1bNiMR6AHAFjLhjX6qEzdv/jii8rOzlZCQoJyc3P15ptv9li3qqpKPp+vSzly5Eg0ugoAgKe4ntFv2bJFjz32mF588UX9+Z//uX74wx+qqKhIhw4dUkZGRo/t6uvrlZiYGDq++uqr3e4qAMAyNkzdu57RP/vss3rooYf09a9/XRMnTtSaNWuUnp6utWvX9touOTlZqampoTJ06FC3uwoAsEzn1L2TEutczejb2tp04MABrVixIux8YWGh9u7d22vbKVOm6MKFC7r++uv1ne98R7Nmzeq2Xmtrq1pbW0PHLS0tkqSv/fQbGpKQ4HAEiHUTtn422F1AFO369hWD3QVEQbu5GLXvMg4z+ngI9K5m9KdPn1ZHR4dSUlLCzqekpKipqanbNmlpaVq/fr0qKiq0detW5eTkaPbs2dqzZ0+39cvLyxUIBEIlPT094uMAACBeRWXXvc8X/ovHGNPlXKecnBzl5OSEjgsKCtTY2KjVq1fr9ttv71K/rKxMpaWloeOWlhaCPQCgT4wkY5y1j3WuBvqkpCQNHTq0S/Z+6tSpLll+b/Lz87V58+Zu/+b3++X3+x31EwBgp6B88vFkvIEbMWKEcnNzVVlZGXa+srJS06ZN6/Pn1NXVKS0tLdLdAwDA81yfui8tLVVxcbHy8vJUUFCg9evXq6GhQUuWLJF0aer9+PHjevnllyVJa9asUVZWliZNmqS2tjZt3rxZFRUVqqiocLurAADL2PDAHNcD/f33368zZ87oqaee0smTJzV58mRt375dmZmZkqSTJ0+qoaEhVL+trU3Lli3T8ePHNXLkSE2aNEnbtm3TvHnz3O4qAMAyQeOTz+P30fuMcbINIfa0tLQoEAgo87vf4/Y6C3B7nV181f9vsLuAKGg3F1Wln6q5uTnswWmR1BkrJv/k2xr6hYHv8+r4rFW//stnXO2rUzzrHgBgLWMc7rqPg1SZ19QCAKwV60/G+9KXvqSMjAwlJCQoLS1NxcXFOnHiRL8+g0APAECMmjVrln7yk5+ovr5eFRUVev/99/WVr3ylX5/B1D0AwFqR2nXf+fj1TpF6xsvSpUtD/87MzNSKFSt0zz336OLFixo+fHifPoOMHgBgrc631zkpkpSenh72OPby8vKI9/W3v/2tXn31VU2bNq3PQV4iowcAWCxSm/EaGxvDdt1H8omtf/d3f6fnn39en332mfLz8/Xf//3f/WpPRg8AgEOJiYlhpbdAv3LlSvl8vl5LbW1tqP63v/1t1dXVaefOnRo6dKgefPBB9efOeDJ6AIC1LmX0Ttbo+9+mpKRECxYs6LVOVlZW6N9JSUlKSkrStddeq4kTJyo9PV01NTUqKCjo0/cR6AEA1hqMR+B2Bu6Bfd+lXxatra19bkOgBwAgBu3bt0/79u3T9OnTNXr0aH3wwQf6h3/4B11zzTV9zuYl1ugBABYzEShuGTlypLZu3arZs2crJydHf/VXf6XJkydr9+7d/drsR0YPALBWLL+97oYbbtCuXbscfw4ZPQAAHkZGDwCwl9P59zh4qQ2BHgBgL6cvpomD99ET6AEA1uI1tQAAIK6R0QMArBXLu+4jhUAPALCX8TlbZ4+DQM/UPQAAHkZGDwCwlg2b8Qj0AAB7WXAfPVP3AAB4GBk9AMBa7LoHAMDr4mD63Qmm7gEA8DAyegCAtZi6BwDAyyzYdU+gBwBYzPe74qR9bGONHgAADyOjBwDYi6l7AAA8zIJAz9Q9AAAeRkYPALCXBa+pJdADAKxlw9vrmLoHAMDDyOgBAPayYDMegR4AYC8L1uiZugcAwMPI6AEA1vKZS8VJ+1hHoAcA2Is1egAAPIw1egAAEM/I6AEA9mLqHgAAD7Mg0DN1DwCAh5HRAwDsZUFGT6AHANiLXfcAACCekdEDAKzFk/EAAPAyC9boXZ2637Nnj+666y6NGzdOPp9Pr7/++mXb7N69W7m5uUpISNCECRO0bt06N7sIAICnuRroz58/r5tuuknPP/98n+ofPXpU8+bN04wZM1RXV6fHH39cjz76qCoqKtzsJgAAnuXq1H1RUZGKior6XH/dunXKyMjQmjVrJEkTJ05UbW2tVq9erfnz53fbprW1Va2traHjlpYWR30GANjDJ4dr9BHriXtiatd9dXW1CgsLw87NmTNHtbW1unjxYrdtysvLFQgEQiU9PT0aXQUAeEHn7XVOSoyLqUDf1NSklJSUsHMpKSlqb2/X6dOnu21TVlam5ubmUGlsbIxGVwEAiAsxt+ve5wv/dWSM6fZ8J7/fL7/f73q/AAAeZMGu+5gK9KmpqWpqago7d+rUKQ0bNkxjx44dpF4BADzLgkAfU1P3BQUFqqysDDu3c+dO5eXlafjw4YPUKwAA4pergf7TTz/VwYMHdfDgQUmXbp87ePCgGhoaJF1aX3/wwQdD9ZcsWaIPP/xQpaWlOnz4sF566SVt2LBBy5Ytc7ObAABLdT4Zz0mJda5O3dfW1mrWrFmh49LSUknSwoULtWnTJp08eTIU9CUpOztb27dv19KlS/XCCy9o3Lhxeu6553q8tQ4AAEcsmLp3NdDfcccdoc103dm0aVOXczNnztRbb73lYq8AALBHTG3GAwAgqsjoAQDwLhveXhdTu+4BAEBkkdEDAOzl9DG2cfAIXAI9AMBerNEDAOBdrNEDAIC4RkYPALAXU/cAAHiY08fYxkGgZ+oeAAAPI6MHANiLqXsAADzMgkDP1D0AAB5GRg8AsBb30QMAgLhGoAcAIMa1trbq5ptvls/n08GDB/vVlkAPALCXiUCJguXLl2vcuHEDassaPQDAWpFao29paQk77/f75ff7HfTs937xi19o586dqqio0C9+8Yt+tyejBwDYLQLZfHp6ugKBQKiUl5dHpGsfffSRFi9erFdeeUVf+MIXBvQZZPQAADjU2NioxMTE0HEksnljjBYtWqQlS5YoLy9Px44dG9DnkNEDAOwVoTX6xMTEsNJboF+5cqV8Pl+vpba2Vj/4wQ/U0tKisrIyR0MkowcAWGsw7qMvKSnRggULeq2TlZWl7373u6qpqenyoyEvL08PPPCA/u3f/q1P30egBwAgipKSkpSUlHTZes8995y++93vho5PnDihOXPmaMuWLZo6dWqfv49ADwCwVww/6z4jIyPs+Morr5QkXXPNNRo/fnyfP4dADwCwlg2PwCXQAwAQB7KysmRM/39ZEOgBAPaK4an7SCHQAwDsZUGg5z56AAA8jIweAGAtNuMBAOBlFkzdE+gBAPayINCzRg8AgIeR0QMArMUaPQAAXsbUPQAAiGdk9AAAazF1DwCAlzF1DwAA4hkZPQDAXhZk9AR6AIC1fL8rTtrHOqbuAQDwMDJ6AIC9mLoHAMC7uL0OAAAvsyCjZ40eAAAPI6MHANgtDrJyJwj0AABr2bBGz9Q9AAAeRkYPALAXm/Gc2bNnj+666y6NGzdOPp9Pr7/+eq/1q6qq5PP5upQjR4642U0AgKU6p+6dlFjnakZ//vx53XTTTfra176m+fPn97ldfX29EhMTQ8dXX321G90DAMDzXA30RUVFKioq6ne75ORkXXXVVZHvEAAAf8iCqfuYXKOfMmWKLly4oOuvv17f+c53NGvWrB7rtra2qrW1NXTc0tIiSVr4F28o4crhrvcVg2vXt68Y7C4gin554uBgdwFR0HIuqNHXRue72HUfZWlpaVq/fr0qKiq0detW5eTkaPbs2dqzZ0+PbcrLyxUIBEIlPT09ij0GACC2xVRGn5OTo5ycnNBxQUGBGhsbtXr1at1+++3dtikrK1NpaWnouKWlhWAPAOgbC6buYyqj705+fr7efffdHv/u9/uVmJgYVgAA6BMTgRLjYiqj705dXZ3S0tIGuxsAAA+yYY3e1UD/6aef6r333gsdHz16VAcPHtSYMWOUkZGhsrIyHT9+XC+//LIkac2aNcrKytKkSZPU1tamzZs3q6KiQhUVFW52EwAAz3I10NfW1obtmO9cS1+4cKE2bdqkkydPqqGhIfT3trY2LVu2TMePH9fIkSM1adIkbdu2TfPmzXOzmwAAW1mwRu9qoL/jjjtkTM//KWzatCnsePny5Vq+fLmbXQIAIMRnjHy9xKm+tI91Mb8ZDwAADFzMb8YDAMA1TN0DAOBdNuy6Z+oeAAAPI6MHANiLqXsAALyLqXsAABDXyOgBAPZi6h4AAO+yYeqeQA8AsJcFGT1r9AAAeBgZPQDAavEw/e4EgR4AYC9jLhUn7WMcU/cAAHgYGT0AwFrsugcAwMvYdQ8AAOIZGT0AwFq+4KXipH2sI9ADAOzF1D0AAIhnZPQAAGux6x4AAC+z4IE5BHoAgLVsyOhZowcAwMPI6AEA9rJg1z2BHgBgLabuAQBAXCOjBwDYi133AAB4F1P3AAAgrhHoAQD2MhEoLsrKypLP5wsrK1as6NdnMHUPALBWPEzdP/XUU1q8eHHo+Morr+xXewI9AAAOtbS0hB37/X75/f6IfPaoUaOUmpo64PZM3QMA7BU0zouk9PR0BQKBUCkvL49YF59++mmNHTtWN998s773ve+pra2tX+3J6AEA9orQk/EaGxuVmJgYOh2pbP5b3/qWbrnlFo0ePVr79u1TWVmZjh49qh/96Ed9/gwCPQDAWj45XKP/3f9PTEwMC/S9WblypVatWtVrnf379ysvL09Lly4Nnbvxxhs1evRofeUrXwll+X1BoAcAIIpKSkq0YMGCXutkZWV1ez4/P1+S9N577xHoAQC4rEF4Ml5SUpKSkpIG9HV1dXWSpLS0tD63IdADAKwVy7fXVVdXq6amRrNmzVIgEND+/fu1dOlSfelLX1JGRkafP4dADwBADPL7/dqyZYtWrVql1tZWZWZmavHixVq+fHm/PodADwCwVwy/j/6WW25RTU2N488h0AMArOUzRj4Ha/RO2kYLD8wBAMDDyOgBAPYK/q44aR/jCPQAAGsxdQ8AAOIaGT0AwF4xvOs+Ugj0AAB7DcKT8aKNQA8AsFYsPxkvUlijBwDAw8joAQD2smDq3tWMvry8XLfeeqtGjRql5ORk3XPPPaqvr79su927dys3N1cJCQmaMGGC1q1b52Y3AQCW8gWdl1jnaqDfvXu3Hn74YdXU1KiyslLt7e0qLCzU+fPne2xz9OhRzZs3TzNmzFBdXZ0ef/xxPfroo6qoqHCzqwAAeJKrU/c7duwIO964caOSk5N14MAB3X777d22WbdunTIyMrRmzRpJ0sSJE1VbW6vVq1dr/vz5bnYXAGAbpu4jq7m5WZI0ZsyYHutUV1ersLAw7NycOXNUW1urixcvdqnf2tqqlpaWsAIAQJ+YCJQYF7VAb4xRaWmppk+frsmTJ/dYr6mpSSkpKWHnUlJS1N7ertOnT3epX15erkAgECrp6ekR7zsAAPEqaoG+pKREb7/9tv793//9snV9Pl/Ysfnd1Mgfn5eksrIyNTc3h0pjY2NkOgwA8LzOZ907KbEuKrfXPfLII/rZz36mPXv2aPz48b3WTU1NVVNTU9i5U6dOadiwYRo7dmyX+n6/X36/P6L9BQBYgjV6Z4wxKikp0datW7Vr1y5lZ2dftk1BQYEqKyvDzu3cuVN5eXkaPny4W10FAMCTXA30Dz/8sDZv3qwf//jHGjVqlJqamtTU1KTPP/88VKesrEwPPvhg6HjJkiX68MMPVVpaqsOHD+ull17Shg0btGzZMje7CgCwkdHv30k/kBL7Cb27gX7t2rVqbm7WHXfcobS0tFDZsmVLqM7JkyfV0NAQOs7Oztb27dtVVVWlm2++Wf/4j/+o5557jlvrAAARxxq9Q6YP/wFs2rSpy7mZM2fqrbfecqFHAAD8ASOHa/QR64lreKkNAAAexkttAAD2smDXPYEeAGCvoKSuj2jpX/sYx9Q9AAAeRkYPALCW053z1u+6BwAgplmwRs/UPQAAHkZGDwCwlwUZPYEeAGAvCwI9U/cAAHgYGT0AwF4W3EdPoAcAWIvb6wAA8DLW6AEAQDwjowcA2CtoJJ+DrDwY+xk9gR4AYC+m7gEAQDwjowcAWMxhRq/Yz+gJ9AAAezF1DwAA4hkZPQDAXkEjR9Pv7LoHACCGmeCl4qR9jGPqHgAADyOjBwDYy4LNeAR6AIC9WKMHAMDDLMjoWaMHAMDDyOgBAPYycpjRR6wnriHQAwDsxdQ9AACIZ2T0AAB7BYOSHDz0Jhj7D8wh0AMA7MXUPQAAiGdk9AAAe1mQ0RPoAQD2suDJeEzdAwDgYWT0AABrGROUcfCqWSdto4VADwCwlzHOpt9ZowcAIIYZh2v0cRDoWaMHAMDDyOgBAPYKBiWfg3V21ugBAIhhTN0DAIB4RkYPALCWCQZlHEzdx8PtdWT0AAB7dT4C10lx2bZt2zR16lSNHDlSSUlJuvfee/vVnoweAIAYVVFRocWLF+uf/umfdOedd8oYo3feeadfn0GgBwDYK2gkn/PNeC0tLWGn/X6//H6/k56pvb1d3/rWt/TMM8/ooYceCp3Pycnp1+cwdQ8AsJcxl26RG3C5FOjT09MVCARCpby83HHX3nrrLR0/flxDhgzRlClTlJaWpqKiIv3v//5vvz6HjB4AAIcaGxuVmJgYOnaazUvSBx98IElauXKlnn32WWVlZen73/++Zs6cqd/85jcaM2ZMnz6HjB4AYC0TNI6LJCUmJoaV3gL9ypUr5fP5ei21tbUKBi/t6H/iiSc0f/585ebmauPGjfL5fPrP//zPPo+RjB4AYC8TlBTdJ+OVlJRowYIFvdbJysrSuXPnJEnXX3996Lzf79eECRPU0NDQ5+9zNdCXl5dr69atOnLkiEaOHKlp06bp6aef7nUjQVVVlWbNmtXl/OHDh3Xddde52V0AgGVM0Mg42IxnBnB7XVJSkpKSki5bLzc3V36/X/X19Zo+fbok6eLFizp27JgyMzP7/H2uBvrdu3fr4Ycf1q233qr29nY98cQTKiws1KFDh3TFFVf02ra+vj5svePqq692s6sAAMSUxMRELVmyRE8++aTS09OVmZmpZ555RpJ033339flzXA30O3bsCDveuHGjkpOTdeDAAd1+++29tk1OTtZVV1112e9obW1Va2tr6Li5ufnS+fPt/e8w4k67uTjYXUAUtZyL/aeQwbmWTy9d54Fky/3VblodvZimXe7+b9AzzzyjYcOGqbi4WJ9//rmmTp2qXbt2afTo0X3/EBNF7777rpFk3nnnnR7rvPHGG0aSycrKMqmpqebOO+80u3bt6rH+k08+2flGAgqFQqF4qLz//vtuhCJjjDGff/65SU1NjUg/U1NTzeeff+5aX53yGROdV+8YY3T33Xfr7NmzevPNN3usV19frz179ig3N1etra165ZVXtG7dOlVVVXU7C/DHGf0nn3yizMxMNTQ0KBAIuDKWWNTS0qL09PQut3h4nY3jtnHMkp3jtnHM0qWZ2YyMDJ09e7ZPM7sDdeHCBbW1tTn+nBEjRighISECPXJH1Hbdl5SU6O2339b//M//9FovJycnbLNeQUGBGhsbtXr16m4DfU9PHwoEAlb9F6NT560dtrFx3DaOWbJz3DaOWZKGDHH3DvCEhISYDtCREpX76B955BH97Gc/0xtvvKHx48f3u31+fr7effddF3oGAIC3uZrRG2P0yCOP6LXXXlNVVZWys7MH9Dl1dXVKS0uLcO8AAPA+VwP9ww8/rB//+Mf66U9/qlGjRqmpqUnSpWn1kSNHSpLKysp0/Phxvfzyy5KkNWvWKCsrS5MmTVJbW5s2b96siooKVVRU9Ok7/X6/nnzyyYg8fjCeMG57xm3jmCU7x23jmCV7x+0WVzfj+Xy+bs9v3LhRixYtkiQtWrRIx44dU1VVlSTpn//5n7V+/XodP35cI0eO1KRJk1RWVqZ58+a51U0AADwrarvuAQBA9PFSGwAAPIxADwCAhxHoAQDwMAI9AAAe5olAf/bsWRUXFysQCCgQCKi4uFiffPJJr20WLVokn88XVvLz86PT4QF68cUXlZ2drYSEBOXm5vb6KGHp0tsDc3NzlZCQoAkTJmjdunVR6mnk9GfMVVVVXa6pz+fTkSNHothj5/bs2aO77rpL48aNk8/n0+uvv37ZNvF+rfs7Zi9c6/Lyct16660aNWqUkpOTdc8996i+vv6y7eL9Wg9k3F643oPJE4H+q1/9qg4ePKgdO3Zox44dOnjwoIqLiy/bbu7cuTp58mSobN++PQq9HZgtW7boscce0xNPPKG6ujrNmDFDRUVFamho6Lb+0aNHNW/ePM2YMUN1dXV6/PHH9eijj/b5eQSxoL9j7lRfXx92Xf/sz/4sSj2OjPPnz+umm27S888/36f6XrjW/R1zp3i+1p2v8a6pqVFlZaXa29tVWFio8+fP99jGC9d6IOPuFM/Xe1AN3vt0IuPQoUNGkqmpqQmdq66uNpLMkSNHemy3cOFCc/fdd0ehh5Fx2223mSVLloSdu+6668yKFSu6rb98+XJz3XXXhZ37m7/5G5Ofn+9aHyOtv2PufPPh2bNno9C76JBkXnvttV7reOFa/6G+jNmL1/rUqVNGktm9e3ePdbx2rY3p27i9eL2jKe4z+urqagUCAU2dOjV0Lj8/X4FAQHv37u21bVVVlZKTk3Xttddq8eLFOnXqlNvdHZC2tjYdOHBAhYWFYecLCwt7HGN1dXWX+nPmzFFtba0uXoz9d7gPZMydpkyZorS0NM2ePVtvvPGGm92MCfF+rZ3w0rVubm6WJI0ZM6bHOl681n0ZdycvXe9oivtA39TUpOTk5C7nk5OTQ4/c7U5RUZFeffVV7dq1S9///ve1f/9+3XnnnWGvvI0Vp0+fVkdHh1JSUsLOp6Sk9DjGpqambuu3t7fr9OnTrvU1UgYy5rS0NK1fv14VFRXaunWrcnJyNHv2bO3ZsycaXR408X6tB8Jr19oYo9LSUk2fPl2TJ0/usZ7XrnVfx+216x1tUXtNbX+tXLlSq1at6rXO/v37JXX/qF1jTI+P4JWk+++/P/TvyZMnKy8vT5mZmdq2bZvuvffeAfbaXX88nsuNsbv63Z2PZf0Zc39fcewlXrjW/eG1a93X13hL3rrWbr2+HOFiNtCXlJRowYIFvdbJysrS22+/rY8++qjL3z7++OMuv3x7k5aWpszMzJh8HW5SUpKGDh3aJZM9depUj2NMTU3ttv6wYcM0duxY1/oaKQMZc3fy8/O1efPmSHcvpsT7tY6UeL3Wna/x3rNnz2Vf4+2la92fcXcnXq/3YIjZQJ+UlKSkpKTL1isoKFBzc7P27dun2267TZL0q1/9Ss3NzZo2bVqfv+/MmTNqbGyMydfhjhgxQrm5uaqsrNSXv/zl0PnKykrdfffd3bYpKCjQz3/+87BzO3fuVF5enoYPH+5qfyNhIGPujg2vOI73ax0p8XatzQBe4+2Faz2QcXcn3q73oBqsXYCRNHfuXHPjjTea6upqU11dbW644QbzxS9+MaxOTk6O2bp1qzHGmHPnzpm//du/NXv37jVHjx41b7zxhikoKDB/8id/YlpaWgZjCJf1H//xH2b48OFmw4YN5tChQ+axxx4zV1xxhTl27JgxxpgVK1aY4uLiUP0PPvjAfOELXzBLly41hw4dMhs2bDDDhw83//Vf/zVYQ+i3/o75X/7lX8xrr71mfvOb35hf//rXZsWKFUaSqaioGKwhDMi5c+dMXV2dqaurM5LMs88+a+rq6syHH35ojPHmte7vmL1wrb/xjW+YQCBgqqqqzMmTJ0Pls88+C9Xx4rUeyLi9cL0HkycC/ZkzZ8wDDzxgRo0aZUaNGmUeeOCBLrdhSDIbN240xhjz2WefmcLCQnP11Veb4cOHm4yMDLNw4ULT0NAQ/c73wwsvvGAyMzPNiBEjzC233BJ2O8rChQvNzJkzw+pXVVWZKVOmmBEjRpisrCyzdu3aKPfYuf6M+emnnzbXXHONSUhIMKNHjzbTp08327ZtG4ReO9N5K9Efl4ULFxpjvHmt+ztmL1zr7sb7h/87ZYw3r/VAxu2F6z2YeE0tAAAeFve31wEAgJ4R6AEA8DACPQAAHkagBwDAwwj0AAB4GIEeAAAPI9ADAOBhBHoAADyMQA8AgIcR6AEA8DACPQAAHvb/AdkCXjsu/EZDAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 640x480 with 2 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(np.log(np.abs(pcov)))\n",
    "plt.colorbar()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e0ed3d63-7503-471f-a2a5-f3f9864d855f",
   "metadata": {},
   "source": [
    "Logistic curve model of cavitation erosion progress in metallic materials\n",
    "\n",
    "The authors previously found that the change in volume loss rate with the exposure time can be expressed by a logistic curve. In this study, the validity of this model is examined for various materials such as pure aluminum, carbon steels, stainless steels, cobalt alloys, and so on. The MDE (mean depth of erosion) d as a function of the exposure time can be expressed by three parameters α, β and c as in the following equation:\n",
    "\n",
    "The parameters α, β and c are derived from the relation between the nominal incubation period and the slope of the maximum rate stage, from the average thickness of the removed layer when the nominal incubation period is terminated, and from an arbitrary point (t0, d0) of the MDE in the maximum rate stage. We conclude that the calculated curve based on this model is in good agreement with the MDE data points for various materials, test conditions and test methods.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "f8f0c46b-74ce-489b-9788-f99ede808700",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x77d502175810>"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAigAAAGiCAYAAADNzj2mAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAB52ElEQVR4nO3deVxU9f7H8dewDTsICoigoqKC+5bhiimYWVm3X/u1xX3L0Mq1W2YFamV6c4v0lrflWvdWXm9lirvmhgtuuIO4gQgiIMsAM9/fH5NTk0ugwGH5PB+Peej3nO+ceXOAmQ/nfM/36JRSCiGEEEKIKsRG6wBCCCGEEH8kBYoQQgghqhwpUIQQQghR5UiBIoQQQogqRwoUIYQQQlQ5UqAIIYQQosqRAkUIIYQQVY4UKEIIIYSocqRAEUIIIUSVIwWKEEIIIaqcMhUojRs3RqfT3fAYO3YsAEopZsyYgb+/P05OToSHh3PkyJEKCS6EEEKImqtMBUp8fDypqamWR1xcHACPP/44AHPmzGHu3LksWLCA+Ph4/Pz8iIiIIDc3t/yTCyGEEKLG0t3NzQKjoqL44YcfOHnyJAD+/v5ERUUxefJkAAwGA76+vsyePZuRI0eWT2IhhBBC1Hh2d/rEoqIivvjiCyZOnIhOpyMpKYm0tDQiIyMtffR6Pb1792b79u23LFAMBgMGg8HSNplMXLlyBW9vb3Q63Z3GE0IIIUQlUkqRm5uLv78/NjZ3P8T1jguUlStXcvXqVV544QUA0tLSAPD19bXq5+vrS0pKyi23ExMTw1tvvXWnMYQQQghRhZw7d46AgIC73s4dFyjLli1jwIAB+Pv7Wy3/41EPpdRtj4RMnTqViRMnWtrZ2dk0bNiQc+fO4e7ufqfxhBBCCFGJcnJyCAwMxM3NrVy2d0cFSkpKCuvWreO7776zLPPz8wPMR1Lq169vWZ6enn7DUZXf0+v16PX6G5a7u7tLgSKEEEJUM+U1POOOThJ9+umn+Pj4MHDgQMuyoKAg/Pz8LFf2gHmcyubNm+nWrdvdJxVCCCFErVHmIygmk4lPP/2U559/Hju7356u0+mIiooiOjqa4OBggoODiY6OxtnZmWeeeaZcQwshhBCiZitzgbJu3TrOnj3LkCFDblg3adIkCgoKGDNmDFlZWXTt2pW1a9eW2/koIYQQQtQOdzUPSkXIycnBw8OD7OzsW45BUUpRUlKC0Wis5HRCaM/e3h5bW1utYwghhJXSfH6XxR1fxaOVoqIiUlNTyc/P1zqKEJrQ6XQEBATg6uqqdRQhhKgw1apAMZlMJCcnY2tri7+/Pw4ODjKZm6hVlFJcvnyZ8+fPExwcLEdShBA1VrUqUIqKijCZTAQGBuLs7Kx1HCE0Ua9ePc6cOUNxcbEUKEKIGuvu56LVQHlMoStEdSVHDYUQtYF80gshhBCiypECRQghhBBVjhQoQgghhKhypECpROnp6YwcOZKGDRui1+vx8/Ojf//+7NixA4Dw8HCioqJueN7KlSutxh382XYAFi9eTNu2bS33NAoLC2P16tU3bHvRokUEBQXh6OhIp06d2Lp16y3zx8TEWGYMLq39+/djb29Pz549S/0cLd3uayzLvhJCCHF3pECpRI899hgHDhxg+fLlnDhxglWrVhEeHs6VK1fKfTsBAQHMmjWLPXv2sGfPHu677z4GDRrEkSNHLH2+/vproqKimD59Ovv376dnz54MGDCAs2fP3vCa8fHxxMbG0rZt2zJlHT9+PK+++ioHDhygis0JeIPbfY1l2VdCCCHKgapisrOzFaCys7NvWFdQUKASExNVQUGBBsnuTlZWlgLUpk2bbtmnd+/e6uWXX75h+ffff6+uf6tKs51bqVOnjlq6dKmlfc8996hRo0ZZ9WnZsqWaMmWK1bLc3FwVHBys4uLibpnxZr788kv10EMPqXPnzilAnTp1qsyZK8uffY2l2Vdbt25VdnZ2qrCw0LIsKSlJAerMmTNKKaWSk5MVoL799lvVs2dP5ejoqDp27KiSk5PVxo0bVZcuXZSTk5MKDw9XmZmZN81anX8PhBA11+0+v+9EtT+CopQiv6hEk4cqwxEBV1dXXF1dWblyJQaD4Y6/3jvZjtFoZMWKFeTl5REWFgaY55TZu3cvkZGRVn0jIyPZvn271bKxY8cycOBA+vXrV+qceXl5TJs2jdmzZxMQEICHhwcJCQmlfv6diI6OtuyfWz1udVrmdl9jafdVQkICISEh6PV6q2Wenp40atTI0gbz6aLo6Gh27NhBZmYmgwcPZvbs2SxcuJBNmzZx6NAhli1bdre7RAghqq1qNVHbzRQUGwl9Y40mr504sz/ODqXbhXZ2dnz22WcMHz6cJUuW0LFjR3r37s1TTz1VptMmZdnOoUOHCAsLo7CwEFdXV77//ntCQ0MByMjIwGg04uvra/UcX19f0tLSLO0VK1awb98+4uPjS50R4N133+X+++8nJCQEgNDQUBISEnjsscfKtJ2yGDVqFE888cRt+zRo0OCGZX/2NZZ2Xx04cIAOHTpY9UlISKBdu3ZWferUqcOKFSuoW7cuAH369GHDhg0kJibi4uICQJcuXay2LYQQtU21P4JSnTz22GNcvHiRVatW0b9/fzZt2kTHjh357LPPKmQ7LVq0ICEhgZ07dzJ69Gief/55EhMTrfr8cdIvpZRl2blz53j55Zf54osvcHR0LHW+pKQkYmNjmTFjhmVZ69atK/wIipeXF82aNbvtw8nJyeo5Zfkab7evwFyMtG/f3qrP/v37rQqUhIQEHn74YUtxAnD27FmefvppS3FyfVlQUFCpv3YhhKhxyuVEUTkq6xgUk8mk8gzFmjxMJtNdf71Dhw5VDRs2VEop9dBDD6kXXnjhhj6ffvqpcnd3L/V2bqVv375qxIgRSimlDAaDsrW1Vd99951Vn/Hjx6tevXoppX4b+2Jra2t5AEqn0ylbW1tVUlJy09d5+OGHb3ieTqdTAQEBVv0GDRqknnjiCdW5c2fVrFkzlZCQcMvspen77rvvKhcXl9s+tmzZYvWc0nyNpdlXJSUlysnJSW3cuNGqT4MGDdQ//vEPSzsoKEgtXrzYqo+Xl5fVtgsKCpSdnZ3aunXrTfeFjEERQlRF5T0Gpdqf4tHpdKU+zVIVhYaGsnLlSgBatmx500uB4+PjadGiRam3cytKKcu4FQcHBzp16kRcXByPPvqopU9cXByDBg0CoG/fvhw6dMhqGy+++CItW7Zk8uTJN70PTFxcHL/88gv79+/Hzu6370t8fDxDhgwhMzMTb29vAA4ePMjo0aP5+uuvWb58OfPmzePTTz+9afbS9L2TUzyl+RptbW3/dF8dP36cgoIC/P39Let37NjBhQsXLEdQcnJyOHPmjNVpoJSUFK5cuWK17MiRIxiNRqsjL0IIUeuUS5lTjmrqVTwZGRmqT58+6vPPP1cHDhxQSUlJ6ptvvlG+vr5qyJAhSinzFR5OTk5qzJgxKiEhQR0/flwtWLBA6fV69c0335R6O0opNXXqVLVlyxaVnJysDh48qKZNm6ZsbGzU2rVrLX1WrFih7O3t1bJly1RiYqKKiopSLi4ulitObuZ2V/EUFxerkJAQNWvWrBvWpaSkKECtW7dOKWW+asbf399yFGbXrl1q0KBBN91uWfqWh5t9jX+2r7788ksFqHHjxqkTJ06on376STVr1kwBavv27UoppbZs2aJsbW1Vfn6+ZbsrV65Unp6eVq+1bNky1bRp01vmq86/B0KImkuOoFRTrq6udO3alQ8//JDTp09TXFxMYGAgw4cPZ9q0aQA0btyYrVu3Mn36dCIjIyksLKR58+Z89tlnPP7446XeDsClS5cYPHgwqampeHh40LZtW37++WciIiIsfZ588kkyMzOZOXMmqamptG7dmp9++slyxUlZffTRR2RmZjJu3Lgb1l2/A3VCQgJ9+/bl4MGDtGrVynIUZt++fbRp04Yvv/yS77//3vK80aNH4+TkdNO+lenP9lVCQgIREREkJyfTunVrQkNDmTVrFkOGDGHhwoWEhYVx4MABWrZsaTUO5o9jVMA8kPaPY1mEEKK20SlVtWbPysnJwcPDg+zsbNzd3a3WFRYWkpycbJnNU1RfS5Ys4cMPP+Tw4cNcvXqVfv36sWrVqpsWR2Xpq5X+/fvTsWNHYmJiKvy15PdACFEV3e7z+07IERShiYMHD/Loo4/SvXt3CgoK+OCDD25ZcJSlr1YOHDjACy+8oHUMIYSoMaRAEZo4ePAg//rXv5g1a1a59tVCWloaly5dKvNtAIQQQtyaFChCExcuXCAwMLDc+2rBz8+vyt9nSAghqhuZqE1oIjk5uUL6CiGEqBmkQBFCCCHuUGJm4p93EndEChQhhBDiDnx59Eue/OFJ/nH4H1pHqZFkDIoQQghRRsuPLOf9Pe8DcLXwqrZhaigpUIQQQogy+OTgJ/x9/98BGNF2BOPa3zg5pbh7UqAIIYQQpaCUYvGBxSw+sBiAse3HMqrdKI1T1VxSoAghhBB/QinF/H3zWXZ4GQATOk1gSOshGqeq2aRAEUIIIW5DKcV7e97j88TPAZjUZRKDQwdrnKrmkwJFCCGEuAWTMhGzK4YVx1cA8HrX13my5ZMap6odpEARQgghbsKkTMzcMZNvT36LDh0zus3gL8F/0TpWrSEFihBCCPEHRpORN7a/warTq7DR2fB297d5uOnDWseqVaRAEUIIIX6nxFTCtG3TWJ28GludLdE9onmgyQNax6p1ZCbZSpSens7IkSNp2LAher0ePz8/+vfvz44dOwAIDw8nKirqhuetXLkSnU5X6u0ALF68mLZt2+Lu7o67uzthYWGsXr36hm0vWrSIoKAgHB0d6dSpE1u3br1l/piYGHQ63U0z3sr+/fuxt7enZ8+epX5OZZsxYwY6nc7q4efnd0O/suwrIUT1VGwsZtKWSaxOXo2dzo73er8nxYlG5AhKJXrssccoLi5m+fLlNGnShEuXLrF+/XquXLlS7tsJCAhg1qxZNGvWDIDly5czaNAg9u/fT6tWrQD4+uuviYqKYtGiRXTv3p2PP/6YAQMGkJiYSMOGDa1eMz4+ntjYWNq2bVumrOPHj+fVV19l4cKFKKWsCq2qpFWrVqxbt87StrW1tVpfln0lhKieioxFvLLpFTad34S9jT1zw+cSHhiudazaS1Ux2dnZClDZ2dk3rCsoKFCJiYmqoKDgt4Umk1KGa9o8TKZSf11ZWVkKUJs2bbpln969e6uXX375huXff/+9uv6tKs12bqVOnTpq6dKllvY999yjRo0aZdWnZcuWasqUKVbLcnNzVXBwsIqLi7tlxpv58ssv1UMPPaTOnTunAHXq1KkyZ64Mb775pmrXrt1t+5RmX23dulXZ2dmpwsJCy7KkpCQFqDNnziillEpOTlaA+vbbb1XPnj2Vo6Oj6tixo0pOTlYbN25UXbp0UU5OTio8PFxlZmbeNMtNfw+EEHeloLhAjYwbqVp/1lp1+ryT2np+q9aRqp3bfX7fiep/BKU4H6L9tXntaRfBwaVUXV1dXXF1dWXlypXce++96PX6O3rJO9mO0Wjk3//+N3l5eYSFhQFQVFTE3r17mTJlilXfyMhItm/fbrVs7NixDBw4kH79+vHOO++UKmdeXh7Tpk1j9erVBAQE4OHhQUJCAk2bNi3lV1p20dHRREdH37bP6tWrb3q66eTJk/j7+6PX6+natSvR0dE0adIEKP2+SkhIICQkxOp7kpCQgKenJ40aNbK0wXy6KDo6GldXVx555BEGDx6Mq6ur5UjTAw88wLJly3jttdfuaF8IIUovvzif8RvHsyt1F052Tnx030d0rd9V61i1XvUvUKoJOzs7PvvsM4YPH86SJUvo2LEjvXv35qmnnirTaZOybOfQoUOEhYVRWFiIq6sr33//PaGhoQBkZGRgNBrx9fW1eo6vry9paWmW9ooVK9i3bx/x8fFl+nrfffdd7r//fkJCQgAIDQ0lISGBxx57rEzbKYtRo0bxxBNP3LZPgwYNbljWtWtX/vnPf9K8eXMuXbrEO++8Q7du3Thy5Aje3t6l3lcHDhygQ4cOVn0SEhJo166dVZ86deqwYsUK6tatC0CfPn3YsGEDiYmJuLiYC94uXbpYbVsIUTHyivMYu34sey/txdnOmUX9FtHJt5PWsQQ1oUCxdzYfydDqtcvgscceY+DAgWzdupUdO3bw888/M2fOHJYuXcoLL7xQ7ttp0aIFCQkJXL16lW+//Zbnn3+ezZs3W4oU4IYxIep340TOnTvHyy+/zNq1a3F0dCx1vqSkJGJjYzl8+LBlWevWrS1HDyqKl5cXXl5eZX7egAEDLP9v06YNYWFhNG3alOXLlzNx4kTLutvtKzAXI88884xVn/3791sVKAkJCTz88MOW4gTg7NmzPP3005bi5PqygQMHlvlrEUKUXm5RLqPXjebA5QO42ruyuN9i2vu01zqWuK6s54TOnz+vnn32WeXl5aWcnJxUu3bt1J49eyzrTSaTevPNN1X9+vWVo6Oj6t27tzp8+HCpt1/mMSjV3NChQ1XDhg2VUko99NBD6oUXXrihz6effqrc3d1LvZ1b6du3rxoxYoRSSimDwaBsbW3Vd999Z9Vn/PjxqlevXkqp38a+2NraWh6A0ul0ytbWVpWUlNz0dR5++OEbnqfT6VRAQIBVv0GDBqknnnhCde7cWTVr1kwlJCTcMntp+r777rvKxcXlto8tW7bcdh9d169fP8uYk9Lsq5KSEuXk5KQ2btxo1adBgwbqH//4h6UdFBSkFi9ebNXHy8vLatsFBQXKzs5Obd1683PgNfH3QIjKdrXwqnrqf0+p1p+1VmFfhalDlw9pHanaK+8xKGW6zDgrK4vu3btjb2/P6tWrSUxM5IMPPsDT09PSZ86cOcydO5cFCxYQHx+Pn58fERER5ObmlmNZVXOEhoaSl5cHQMuWLdmzZ88NfeLj42nRokWpt3MrSikMBgMADg4OdOrUibi4OKs+cXFxdOvWDYC+ffty6NAhEhISLI/OnTvz7LPPkpCQcMOVLtef/8svv7B//36r5y1btozz58+TmZlp6Xvw4EE6d+5MfHw8r7/+OvPmzbtl9tL0HTVqlNVr3uzRuXPn2+4jAIPBwNGjR6lfv36p99Xx48cpKCjA3/+38VA7duzgwoULliMoOTk5nDlzxuo0UEpKCleuXLFaduTIEYxGo9WRFyFE+ckqzGL42uEczjyMp96TZZHLaF23tdaxxB+VpZqZPHmy6tGjxy3Xm0wm5efnp2bNmmVZVlhYqDw8PNSSJUtK9Ro19QhKRkaG6tOnj/r888/VgQMHVFJSkvrmm2+Ur6+vGjJkiFLKfIWHk5OTGjNmjEpISFDHjx9XCxYsUHq9Xn3zzTel3o5SSk2dOlVt2bJFJScnq4MHD6pp06YpGxsbtXbtWkufFStWKHt7e7Vs2TKVmJiooqKilIuLi+WKk5u53VU8xcXFKiQkxOr7f11KSooC1Lp165RS5iuD/P39LUdhdu3apQYNGnTT7Zal75145ZVX1KZNm1RSUpLauXOnevDBB5Wbm5vVfvizffXll18qQI0bN06dOHFC/fTTT6pZs2YKUNu3b1dKKbVlyxZla2ur8vPzLdtduXKl8vT0tMqzbNky1bRp01vmrc6/B0JoLSM/Qz3630dV689aq14reqkTV05oHanG0PQqnlWrVtG/f38ef/xxNm/eTIMGDRgzZgzDhw8HIDk5mbS0NCIjIy3P0ev19O7dm+3btzNy5MgbtmkwGCx/1YP5r8yayNXVla5du/Lhhx9y+vRpiouLCQwMZPjw4UybNg2Axo0bs3XrVqZPn05kZCSFhYU0b96czz77jMcff7zU2wG4dOkSgwcPJjU1FQ8PD9q2bcvPP/9MRESEpc+TTz5JZmYmM2fOJDU1ldatW/PTTz9Zrjgpq48++ojMzEzGjRt3w7rAwECcnZ1JSEigb9++HDx4kFatWlmOwuzbt482bdrw5Zdf8v3331ueN3r0aJycnG7at7ycP3+ep59+moyMDOrVq8e9997Lzp07rfbDn+2rhIQEIiIiSE5OpnXr1oSGhjJr1iyGDBnCwoULCQsL48CBA7Rs2RInJyfLdv84RgXMA2nbt29fbl+fEMLscv5lhq0dRlJ2EvWc6rG0/1KaeDTROpa4lbJUM3q9Xun1ejV16lS1b98+tWTJEuXo6KiWL1+ulFLql19+UYC6cOGC1fOGDx+uIiMjb7rNN998UwE3PGraERRhbfHixap58+aqqKhIpaenq7Zt297yyE1Z+molMjLyhvljKor8HghRdqnXUtXA7waq1p+1Vn2/6avOZFet95CaQNMjKCaTic6dO1vmmujQoQNHjhxh8eLFPPfcc5Z+f3a1w+9NnTrV6kqJnJwcAgMDyxJLVEMHDx7k0UcfpXv37hQUFPDBBx/c8shNWfpq5cCBA2W6EksIUXkuXLvA0DVDuXDtAg1cG7A0cikBbgFaxxJ/okwFSv369a0uUQUICQnh22+/BbDcvyQtLc0ywBDM94754xwS1+n1+juetExUXwcPHuRf//oXs2bNKte+WkhLS+PSpUtlvg2AEKLincs5x9C1Q0nNSyXQLZBlkcuo71r/z58oNFemq3i6d+/O8ePHrZadOHHC8tdsUFAQfn5+Vlc7FBUVsXnzZsvVDkIAXLhwodRHysrSVwt+fn4opSz3OBJCVA3J2cm8sOYFUvNSaezemE/7fyrFSTVSpiMoEyZMoFu3bkRHR/PEE0+we/duYmNjiY2NBbDc6TY6Oprg4GCCg4OJjo7G2dn5hgmsRO2WnJxcIX2FEALg9NXTDFs7jIyCDJp6NGVp/6XUdar7508UVUaZCpQuXbrw/fffM3XqVGbOnElQUBDz5s3j2WeftfSZNGkSBQUFjBkzhqysLLp27cratWtxc3Mr9/BCCCHEHx2/cpwRcSO4UniF5nWa80nkJ3g5ln2WaaEtnVJKaR3i93JycvDw8CA7Oxt3d3erdYWFhSQnJxMUFFSmqdeFqEnk90CIW0vMTGRE3AiyDdmEeocSGxGLh95D61i1wu0+v+9EmcagCCGEEFXVwcsHGbZmGNmGbNrWa8snkZ9IcVKNSYEihBCi2tufvp8RcSPILc6lo09HYiNicXe4+7/ihXakQBFCCFGtxafFMzJuJHnFedzjdw+L+y3Gxd7lz58oqjQpUIQQQlRb2y9uZ8y6MRSUFNDNvxsL+i7A2d5Z61iiHEiBIoQQolracn4LL61/iUJjIb0CevH3+/6Ok53Tnz9RVAtSoAghhKh2NpzdwMsbX6bIVMR9gfcxL3weeluZlbwmkQJFCCFEtbLmzBpe2fQKJaYSIhtF8n74+9jb2msdS5QzKVCEEEJUGz8k/cCkLZMoUSU82ORBZveajb2NFCc1kRQolSg9PZ2RI0fSsGFD9Ho9fn5+9O/fnx07dgAQHh5OVFTUDc9buXKl1d2g/2w7AIsXL6Zt27a4u7vj7u5OWFgYq1evvmHbixYtskz41alTJ7Zu3XrL/DExMZbbGZTW/v37sbe3p2fPnqV+jhYuXLjAX//6V7y9vXF2dqZ9+/bs3bvXqk9Z9pUQovx9f/J7pm2dhkmZeKTZI7zT/R3sbMo0IbqoRuQ7W4kee+wxiouLWb58OU2aNOHSpUusX7+eK1eulPt2AgICmDVrFs2aNQNg+fLlDBo0iP3791tuavf1118TFRXFokWL6N69Ox9//DEDBgwgMTGRhg0bWr1mfHw8sbGxZb5j7/jx43n11VdZuHAhSimrQquqyMrKonv37vTp04fVq1fj4+PD6dOn8fT0tPQpy74SQpS/b45/w9s73wbgieZPMP3e6djo5G/sGk1VMdnZ2QpQ2dnZN6wrKChQiYmJqqCgQINkdycrK0sBatOmTbfs07t3b/Xyyy/fsPz7779X179VpdnOrdSpU0ctXbrU0r7nnnvUqFGjrPq0bNlSTZkyxWpZbm6uCg4OVnFxcbfMeDNffvmleuihh9S5c+cUoE6dOlXmzJVh8uTJqkePHrftU5p9tXXrVmVnZ6cKCwsty5KSkhSgzpw5o5RSKjk5WQHq22+/VT179lSOjo6qY8eOKjk5WW3cuFF16dJFOTk5qfDwcJWZmXnTLNX590CIO/FF4heq9WetVevPWquYXTHKZDJpHUncxO0+v+9EtS8/lVLkF+dr8lBluI2Rq6srrq6urFy5EoPBcMdf751sx2g0smLFCvLy8ggLCwOgqKiIvXv3EhkZadU3MjKS7du3Wy0bO3YsAwcOpF+/fqXOmZeXx7Rp05g9ezYBAQF4eHiQkJBQ6uffiejoaMv+udXjZqdlVq1aRefOnXn88cfx8fGhQ4cOfPLJJ5b1pd1XCQkJhISEoNfrrZZ5enrSqFEjSxvMp4uio6PZsWMHmZmZDB48mNmzZ7Nw4UI2bdrEoUOHWLZsWXnuHiGqpc8Of8as3bMAeLHVi0zuMrlKHokV5a/an+IpKCmg61ddNXntXc/sKvWEQHZ2dnz22WcMHz6cJUuW0LFjR3r37s1TTz1VptMmZdnOoUOHCAsLo7CwEFdXV77//ntCQ0MByMjIwGg04uvra/UcX19f0tLSLO0VK1awb98+4uPjS50R4N133+X+++8nJCQEgNDQUBISEnjsscfKtJ2yGDVqFE888cRt+zRo0OCGZUlJSSxevJiJEycybdo0du/ezfjx49Hr9Tz33HOl3lcHDhygQ4cOVn0SEhJo166dVZ86deqwYsUK6tY13/q9T58+bNiwgcTERFxczLNfdunSxWrbQtRGsQdj+Wj/RwCMaDuCce3HSXFSi1T7IyjVyWOPPcbFixdZtWoV/fv3Z9OmTXTs2JHPPvusQrbTokULEhIS2LlzJ6NHj+b5558nMTHRqs8ff9nV78aJnDt3jpdffpkvvviiTHfNTUpKIjY2lhkzZliWtW7dusKPoHh5edGsWbPbPpycbpzEyWQy0bFjR6Kjo+nQoQMjR45k+PDhLF682Krf7fYVmIuR9u3bW/XZv3+/VYGSkJDAww8/bClOAM6ePcvTTz9tKU6uLwsKCrqj/SBEdaeUYmHCQktxMrb9WF7q8JIUJ7VMtT+C4mTnxK5ndmn22mXl6OhIREQEERERvPHGGwwbNow333yTF154AXd3d7Kzs294ztWrV2+4dfXttnOdg4ODZZBs586diY+PZ/78+Xz88cfUrVsXW1vbG/5KT09Ptxwp2Lt3L+np6XTq1Mmy3mg0smXLFhYsWIDBYMDW1vaGvBMmTCAzM5OAgADLMpPJdMPRi0ceeQS9Xk9SUhJXr17lP//5j9WHeVn7RkdHEx0dfdPnX7d69eobriiqX7++5cjSdSEhIXz77bcApdpXRqORI0eO3HAEZd++fTz66KOW9oEDB5g0aZJVn4SEBMaNG2dpFxYWcuLEiRuKHSFqA6UU8/fNZ9lh8ynOqI5RDG0zVONUQgvVvkDR6XTV+r4LoaGhrFy5EoCWLVve9FLg+Ph4WrRoUert3IpSyjJuxcHBgU6dOhEXF2f1ARoXF8egQYMA6Nu3L4cOHbLaxosvvkjLli2ZPHnyTYuTuLg4fvnlF/bv34+d3W8/XvHx8QwZMoTMzEy8vb0BOHjwIKNHj+brr79m+fLlzJs3j08//fSm2UvT905P8XTv3p3jx49bLTtx4oRl3Ehp9tXx48cpKCjA39/fsn7Hjh1cuHDBUkjl5ORw5swZqyImJSWFK1euWC07cuQIRqPxlsWaEDWVUor39rzH54mfA/Ba59d4rtVzGqcSmimXobblqKZexZORkaH69OmjPv/8c3XgwAGVlJSkvvnmG+Xr66uGDBmilDJf4eHk5KTGjBmjEhIS1PHjx9WCBQuUXq9X33zzTam3o5RSU6dOVVu2bFHJycnq4MGDatq0acrGxkatXbvW0mfFihXK3t5eLVu2TCUmJqqoqCjl4uJiueLkZm53FU9xcbEKCQlRs2bNumFdSkqKAtS6deuUUuYrg/z9/VVJSYlSSqldu3apQYMG3XS7Zel7J3bv3q3s7OzUu+++q06ePKm+/PJL5ezsrL744gtLnz/bV19++aUC1Lhx49SJEyfUTz/9pJo1a6YAtX37dqWUUlu2bFG2trYqPz/fst2VK1cqT09PqzzLli1TTZs2vWXe6vx7IMStGE1G9c6OdyxX6/zr6L+0jiTKqLyv4qn2R1CqC1dXV7p27cqHH37I6dOnKS4uJjAwkOHDhzNt2jQAGjduzNatW5k+fTqRkZEUFhbSvHlzPvvsMx5//PFSbwfg0qVLDB48mNTUVDw8PGjbti0///wzERERlj5PPvkkmZmZzJw5k9TUVFq3bs1PP/1kOXJQVh999BGZmZlWpyuuCwwMxNnZmYSEBPr27cvBgwdp1aqV5SjMvn37aNOmDV9++SXff/+95XmjR4/Gycnppn3LS5cuXfj++++ZOnUqM2fOJCgoiHnz5vHss89a+vzZvkpISCAiIoLk5GRat25NaGgos2bNYsiQISxcuJCwsDAOHDhAy5YtrcbB/HGMCphPA8npHVGbmJSJmTtm8u3Jb9Gh482wN3msecUNqBfVg06pMlwrWwlycnLw8PAgOzv7hnEXhYWFJCcnW2bzFNXXkiVL+PDDDzl8+DBXr16lX79+rFq16qbFUVn6aqV///507NiRmJiYCn8t+T0QNYnRZOSN7W+w6vQqbHQ2vN39bR5u+rDWscQduN3n952QIyhCEwcPHuTRRx+le/fuFBQU8MEHH9yy4ChLX60cOHDAaoCyEOLPlZhKmLZtGquTV2OrsyWmZwwDggZoHUtUEVKgCE0cPHiQf/3rX8yaNatc+2ohLS2NS5culfk2AELUZsXGYiZvnUxcShx2Ojvm9J5DRKOIP3+iqDWkQBGauHDhAoGBgeXeVwt+fn5lmlVYiNquyFjEK5tfYdO5Tdjb2DM3fC7hgeFaxxJVjBQoQhPJyckV0lcIUbUVlhQStSmKXy78goONA/P6zKNnQNW+27nQhhQoQgghKkVBSQHjN4xnZ+pOHG0d+ajvR9xb/16tY4kqSgoUIYQQFS6/OJ+x68ey59IenOycWNh3IV38umgdS1RhUqAIIYSoULlFuYxZN4aEywm42ruyuN9i2vu01zqWqOKkQBFCCFFhsg3ZjIobxeHMw7g5uBEbEUvruq21jiWqgWpZoMgVE6I2k59/UV1kFWYxMm4kR68cxVPvSWxELCHeIVrHEtWEjdYBysLe3h6A/Px8jZMIoZ2ioiKAm96sUYiqIrMgk6Frh3L0ylG8HL1Y1n+ZFCeiTKrVERRbW1s8PT1JT08HwNnZGZ1Op3EqISqPyWTi8uXLODs7W90tWoiq5HL+ZYatHUZSdhL1nOqxNHIpTTybaB1LVDPV7h3Oz88PwFKkCFHb2NjY0LBhQynORZWUlpfGsLXDSMlJwdfZl2X9l9HIvWrdmkJUD9WuQNHpdNSvXx8fHx+Ki4u1jiNEpXNwcMDGplqdnRW1xIVrFxi6ZigXrl3A38WfZf2XEeAWoHUsUU1VuwLlOltbWzkHL4QQVcS5nHMMXTuU1LxUAt0CWRa5jPqu9bWOJaox+TNMCCHEXUnOTuaFNS+QmpdKY/fGfNr/UylOxF2rtkdQhBBCaO/01dMMXTOUzMJMmno0ZWn/pdR1qqt1LFEDyBEUIYQQd+T4leMMWTOEzMJMWtRpwT/u/4cUJ6LcyBEUIYQQZZaYmciIuBFkG7IJ9Q4lNiIWD72H1rFEDSJHUIQQQpTJwcsHGbZmGNmGbNrWbcsnkZ9IcSLKnRQoQgghSm3fpX2MiBtBbnEuHXw68HHEx7g7uGsdS9RAUqAIIYQolfi0eEatG0VecR73+N3Dkn5LcHVw1TqWqKHKVKDMmDEDnU5n9bg+syuYb2I2Y8YM/P39cXJyIjw8nCNHjpR7aCGEEJVr+8XtjFk3hoKSAsLqh7Gg7wKc7Z21jiVqsDIfQWnVqhWpqamWx6FDhyzr5syZw9y5c1mwYAHx8fH4+fkRERFBbm5uuYYWQghRebae38pL61+i0FhIr4BefNT3I5zsnLSOJWq4MhcodnZ2+Pn5WR716tUDzEdP5s2bx/Tp0/nLX/5C69atWb58Ofn5+Xz11VflHlwIIUTF23B2A+M3jqfIVMR9gfcxL3weelu91rFELVDmAuXkyZP4+/sTFBTEU089RVJSEgDJycmkpaURGRlp6avX6+nduzfbt2+/5fYMBgM5OTlWDyGEENpbc2YNr2x6hRJTCf0b9+f98Pext7XXOpaoJcpUoHTt2pV//vOfrFmzhk8++YS0tDS6detGZmYmaWlpAPj6+lo9x9fX17LuZmJiYvDw8LA8AgMD7+DLEEIIUZ5+TPqRSVsmUaJKeLDJg8zqOQt7GylOROUpU4EyYMAAHnvsMdq0aUO/fv348ccfAVi+fLmlzx9vAa+Uuu1t4adOnUp2drblce7cubJEEkIIUc5WnlrJ1K1TMSkTjzR7hHe6v4OdjczrKSrXXV1m7OLiQps2bTh58qTlap4/Hi1JT0+/4ajK7+n1etzd3a0eQgghtPHvE//mb7/8DYXi8eaP81a3t7C1kTvHi8p3VwWKwWDg6NGj1K9fn6CgIPz8/IiLi7OsLyoqYvPmzXTr1u2ugwohhKhYXx39ipk7ZgLwbMiz/O3ev2Gjk+myhDbKdMzu1Vdf5aGHHqJhw4akp6fzzjvvkJOTw/PPP49OpyMqKoro6GiCg4MJDg4mOjoaZ2dnnnnmmYrKL4QQohwsP7Kc9/e8D8ALrV5gYqeJtz09L0RFK1OBcv78eZ5++mkyMjKoV68e9957Lzt37qRRo0YATJo0iYKCAsaMGUNWVhZdu3Zl7dq1uLm5VUh4IYQQd2/poaXM3zcfgOFthvNSh5ekOBGa0ymllNYhfi8nJwcPDw+ys7NlPIoQQlQgpRSLDyxm8YHFAIxpP4bR7UZrnEpUV+X9+S3DsoUQohZSSjF/33yWHV4GQFTHKIa2GapxKiF+IwWKEELUMkop3tvzHp8nfg7Aa51f47lWz2mcSghrUqAIIUQtYlImYnbFsOL4CgCmd53OUy2f0jiVEDeSAkUIIWoJkzIxc8dMvj35LTp0vBn2Jo81f0zrWELclBQoQghRCxhNRt7Y/garTq/CRmfDO93f4aGmD2kdS4hbkgJFCCFquBJTCdO3Teen5J+w1dkS0zOGAUEDtI4lxG1JgSKEEDVYsamYyVsmE5cSh53Ojjm95xDRKELrWEL8KSlQhBCihioyFvHK5lfYdG4T9jb2zA2fS3hguNaxhCgVKVCEEKIGKiwpZMKmCWy7sA29rZ55febRo0EPrWMJUWpSoAghRA1TUFLA+A3j2Zm6E0dbRz7q+xH31r9X61hClIkUKEIIUYPkF+czdv1Y9lzag7OdMwv7LqSzX2etYwlRZlKgCCFEDZFblMuYdWNIuJyAq70ri/stpr1Pe61jCXFHpEARQogaINuQzeh1ozmUcQg3BzdiI2JpXbe11rGEuGNSoAghRDV3tfAqI+JGcPTKUTz1nsRGxBLiHaJ1LCHuihQoQghRjWUWZDIibgQnsk7g5ejFJ5Gf0LxOc61jCXHXpEARQohq6nL+ZYatHUZSdhJ1neqyLHIZTTybaB1LiHIhBYoQQlRDaXlpDFs7jJScFHycfVgWuYzGHo21jiVEuZECRQghqpmL1y4ydM1Qzl87j7+LP0v7LyXQLVDrWEKUKylQhBCiGjmXe46ha4aSmpdKgGsAy/ovw9/VX+tYQpQ7KVCEEKKaSM5OZtjaYaTnp9PYvTFLI5fi6+KrdSwhKoQUKEIIUQ2cvnqaoWuGklmYSVOPpiztv5S6TnW1jiVEhbHROoAQQojbO37lOEPWDCGzMJPmdZrzj/v/IcWJqPHkCIoQQlRhiZmJjIgbQbYhmxCvEGIjYvF09NQ6lhAVTo6gCCFEFXXw8kGGrRlGtiGbtnXbsrT/UilORK0hBYoQQlRB+y7tY0TcCHKLc+no05GPIz7G3cFd61hCVBopUIQQooqJT4tn1LpR5BXncY/fPSzutxhXB1etYwlRqaRAEUKIKmT7xe2MWTeGgpICwuqHsaDvApztnbWOJUSlkwJFCCGqiC3nt/DS+pcoNBbSK6AXH/X9CCc7J61jCaEJKVCEEKIK2HB2Ay9vfJkiUxH3Bd7HvPB56G31WscSQjNSoAghhMbWnFnDK5teocRUQv/G/Xk//H3sbe21jiWEpqRAEUIIDf2Y9COTt0ymRJUwsMlAZvWchb2NFCdCSIEihBAaWXlqJVO3TsWojAxqOoh3u7+LnY3MnykEyEyyQgihiX+f+Dczd8wE4PHmj/P6va9jo5O/GYW4TgoUIYSoZF8d/YqY3TEAPNPyGabcMwWdTqdxKiGqFilQhBCiEi0/spz397wPwAutXmBip4lSnAhxE1KgCCFEJVl6aCnz980HYHib4bzU4SUpToS4BSlQhBCigimlWHJgCYsOLAJgTPsxjGo7SooTIW5DChQhhKhASin+vv/vLD20FICojlEMbTNU41RCVH1SoAghRAVRSvHBng9YnrgcgNc6v8ZzrZ7TOJUQ1cNdXdMWExODTqcjKirKskwpxYwZM/D398fJyYnw8HCOHDlytzmFEKJaMSkTMbtjLMXJtK7TpDgRogzuuECJj48nNjaWtm3bWi2fM2cOc+fOZcGCBcTHx+Pn50dERAS5ubl3HVYIIaoDkzIxc8dM/nXsX+jQ8WbYmzzd8mmtYwlRrdxRgXLt2jWeffZZPvnkE+rUqWNZrpRi3rx5TJ8+nb/85S+0bt2a5cuXk5+fz1dffVVuoYUQoqoymoz87Ze/8e3Jb7HR2fB297f5v+b/p3UsIaqdOypQxo4dy8CBA+nXr5/V8uTkZNLS0oiMjLQs0+v19O7dm+3bt990WwaDgZycHKuHEEJURyWmEqZtm8aq06uw1dkS0yOGQc0GaR1LiGqpzINkV6xYwb59+4iPj79hXVpaGgC+vr5Wy319fUlJSbnp9mJiYnjrrbfKGkMIIaqUYlMxk7dMJi4lDjudHXN6zyGiUYTWsYSotsp0BOXcuXO8/PLLfPHFFzg6Ot6y3x+v7VdK3fJ6/6lTp5KdnW15nDt3riyRhBBCc0XGIiZumkhcShz2NvbMDZ8rxYkQd6lMR1D27t1Leno6nTp1siwzGo1s2bKFBQsWcPz4ccB8JKV+/fqWPunp6TccVblOr9ej1+vvJLsQQmiusKSQCZsmsO3CNhxsHJh/33x6NOihdSwhqr0yHUHp27cvhw4dIiEhwfLo3Lkzzz77LAkJCTRp0gQ/Pz/i4uIszykqKmLz5s1069at3MMLIYSWCkoKeGnDS2y7sA1HW0cW9lsoxYkQ5aRMR1Dc3Nxo3bq11TIXFxe8vb0ty6OiooiOjiY4OJjg4GCio6NxdnbmmWeeKb/UQgihsfzifMZtGEd8WjxOdk4s7LuQLn5dtI4lRI1R7jPJTpo0iYKCAsaMGUNWVhZdu3Zl7dq1uLm5lfdLCSGEJq4VXWPM+jHsT9+Pq70ri/stpr1Pe61jCVGj6JRSSusQv5eTk4OHhwfZ2dm4u7trHUcIIaxkG7IZvW40hzIO4ebgRmxELK3rtv7zJ4qaRyk49B9o+QA4uGidRnPl/fl9V1PdCyFEbXK18CrD1w7nUMYhPPWeLItcJsVJbXX5BCx/CL4bBptmaZ2mRpKbBQohRClkFmQyPG44J7NO4uXoxSeRn9C8TnOtY4nKVlwI2+bCtg/BWAR2TuB686tUxd2RAkUIIf7E5fzLDFs7jKTsJOo51WNp5FKaeDbROpaobEmb4YcJcOW0uR0cCQ+8D3UaaZurhpICRQghbiMtL41ha4eRkpOCr7Mvy/ovo5G7fCDVKnkZsGY6HFxhbrv6wYDZEDoIbjEJqbh7UqAIIcQtXLx2kaFrhnL+2nn8XfxZ2n8pgW6BWscSlcVkgoQvIO4NKMgCdNBlKPR9Axw9tE5X40mBIoQQN3Eu9xxD1wwlNS+VANcAlvVfhr+rv9axRGVJP2Y+nXP21xvd+raBh+ZDQKfbP0+UGylQhBDiD85kn2Ho2qGk56fT2L0xn0R+gp+Ln9axRGUoLoAt78Mv88FUDPbO0GcadB0NtvKRWZlkbwshxO+cvnqaYWuHkVGQQVOPpiztv5S6TnW1jiUqw+kN8MNEyEo2t5sPgAfmgGdDbXPVUlKgCCHEr45fOc6IuBFcKbxC8zrN+STyE7wcvbSOJSratXRYMw0O/dvcdvM3FyYtH5RBsBqSAkUIIYDEzERGxI0g25BNiFcIsRGxeDp6ah1LVCSTCfb/0zwItjAbdDZwzwjoMx0cZSZzrUmBIoSo9Q5dPsTIdSPJLcqlTd02LIlYgruDfEDVaOlH4X9RcG6nuV2/HTw4Dxp01DKV+B0pUIQQtdr+9P2MXjeavOI8Ovh0YFHfRbg6uGodS1SUonzY8h5s/zuYSsDeBe573XzkRAbBViny3RBC1FrxafGMXT+WgpICuvh1YcF9C3C2d9Y6lqgop9bBj69A1hlzu+WD5gnXPAI0jSVuTgoUIUSttOPiDsZvGE+hsZCw+mHMv28+TnZOWscSFSH3EqyZCoe/NbfdG8CAORDyoLa5xG1JgSKEqHW2nN/ChI0TKDIV0bNBTz7s8yF6W73WsUR5M5lg76ew7i0w/DoItuso87wmejet04k/IQWKEKJW2XB2A69sfoUSUwl9Avvwfu/3cbB10DqWKG+XjpgHwZ7fbW7Xb2+eCda/vYahRFlIgSKEqDXWnlnL5C2TKVElRDaKZFavWdjb2GsdS5SnojzYPBt2LDQPgnVwg75/gy7DwMZW63SiDKRAEULUCj8m/ci0bdMwKRMDmwzkne7vYGcjb4E1ysk4+HEiXD1rboc8BPfPBo8G2uYSd0R+O4UQNd7KUyt545c3UCgeafYIM8JmYCt/TdccOanw8xRIXGluewTCA+9BiwGaxhJ3RwoUIUSN9u8T/2bmjpkAPN78cV6/93VsdDYapxLlwmSEPf+A9TPBkAM6W7h3NIRPBb3MZVPdSYEihKixvjr6FTG7YwB4puUzTLlnCjq5t0rNkHoQfoiCC3vN7QadzDPB1m+rZSpRjqRAEULUSMuPLOf9Pe8D8EKrF5jYaaIUJzWB4RpsioGdi0EZQe8Ofd+AzkNkEGwNIwWKEKLGWXpoKfP3zQdgeJvhvNThJSlOaoLjq+HHVyHnvLkd+gjcPwvc62saS1QMKVCEEDWGUoolB5aw6MAiAMa0H8PodqM1TiXuWs5FWD0Jjv7P3PZsCA98AM0jtc0lKpQUKEKIGkEpxd/3/52lh5YC8HLHlxnWZpjGqcRdMRlh9yew4R0oygUbOwgbC72ngIPcM6mmkwJFCFHtKaV4f8/7/DPxnwC81vk1nmv1nMapxF25mGAeBHtxv7kdcA88NA98W2kYSlQmKVCEENWaSZmYtXsW/zr2LwCmdZ3G0y2f1jiVuGOGXNgYDbuWgDKB3gMiZkDHF8BGLg+vTaRAEUJUWyZlYuaOmXx78lt06Hgj7A3+r/n/aR1L3KljP8JPr0HOBXO71V/Mg2DdfLXNJTQhBYoQoloymoy8uf1N/nv6v9jobJjZbSaDmg3SOpa4E9nnYfVkOPaDue3ZCB6cC836aZtLaEoKFCFEtVNiKmH6tun8lPwTtjpb3u3xLgObDNQ6ligrYwnsjoWN70LRNfMg2G7joddrMghWSIEihKheik3FTNkyhbUpa7HT2TG712wiG8vlptXOhX3mQbCpB8ztwHvNg2B9QrRMJaoQKVCEENVGkbGIVze/ysZzG7GzseOD3h9wX8P7tI4lyqIwx3zEZHeseRCsowdEzIQOz8kgWGFFChQhRLVgMBqYsHECWy9sxcHGgXl95tEzoKfWsURpKWWeaG31JMhNNS9r8zj0jwZXH22ziSpJChQhRJVXUFLA+A3j2Zm6E0dbR/5+398J8w/TOpYoratn4adJcGK1uV0nyDwItqkc/RK3JgWKEKJKyy/OZ9yGccSnxeNk58TCvgvp4tdF61iiNIwlsGuxeV6T4nywsYceUdDzFbB30jqdqOKkQBFCVFm5RbmMWTeGhMsJuNi7sLjfYjr4dNA6liiN83vgf1Fw6ZC53bAbPPgh+LTUNJaoPqRAEUJUSdmGbEavG82hjEO4Objxcb+PaVOvjdaxxJ8pzIb1b0P8UkCBoydEvgPtn5VBsKJMpEARQlQ5VwuvMiJuBEevHMVD70FsRCyh3qFaxxK3oxQkroTVU+BamnlZ2ych8l1wradpNFE9SYEihKhSMgsyGR43nJNZJ/Fy9CI2IpYWXi20jiVuJysFfnoVTq41t72amgfBNgnXNJao3qRAEUJUGZfzLzNs7TCSspOo61SXpZFLaerZVOtY4laMxbBjIWyaBSUFYOsAPSZAj4lg76h1OlHNlemE4OLFi2nbti3u7u64u7sTFhbG6tWrLeuVUsyYMQN/f3+cnJwIDw/nyJEj5R5aCFHzpOWl8eKaF0nKTsLH2YdP+38qxUlVdm43fNwb1r1pLk4a9YBRv0CfaVKciHJRpgIlICCAWbNmsWfPHvbs2cN9993HoEGDLEXInDlzmDt3LgsWLCA+Ph4/Pz8iIiLIzc2tkPBCiJrh4rWLvPjzi6TkpODv4s9n939GY4/GWscSN1NwFX6YAMsiIf0IOHnBI4vhhR+gXnOt04kaRKeUUnezAS8vL9577z2GDBmCv78/UVFRTJ48GQCDwYCvry+zZ89m5MiRpdpeTk4OHh4eZGdn4+7ufjfRhBDVwLmccwxdO5TUvFQCXANY1n8Z/q7+WscSf6QUHP4Wfp4KeenmZe2fhYi3wcVb22yiSijvz+87HoNiNBr597//TV5eHmFhYSQnJ5OWlkZk5G837dLr9fTu3Zvt27ffskAxGAwYDAZLOycn504jCSGqmeTsZIatHUZ6fjqN3RuzNHIpvi6+WscSf3QlGX58BU6vN7e9g81zmgTJrQZExSlzgXLo0CHCwsIoLCzE1dWV77//ntDQULZv3w6Ar6/1m4uvry8pKSm33F5MTAxvvfVWWWMIIaq501dPM2ztMDIKMmjq0ZSl/ZdS16mu1rHE7xmLYftHsHk2lBSCrd48C2yPKLDTa51O1HBlLlBatGhBQkICV69e5dtvv+X5559n8+bNlvU6nc6qv1LqhmW/N3XqVCZOnGhp5+TkEBgYWNZYQohq5PiV44yIG8GVwis0r9OcTyI/wcvRS+tY4vfO7jTPBHv5qLkd1AsGfgh1m2kaS9QeZS5QHBwcaNbM/APauXNn4uPjmT9/vmXcSVpaGvXr17f0T09Pv+Goyu/p9Xr0eqnEhagtEjMTGRE3gmxDNiFeIcRGxOLp6Kl1LHFdQRasmwF7PzO3nb3Ndxxu+yTc5o9NIcrbXc87rJTCYDAQFBSEn58fcXFxlnVFRUVs3ryZbt263e3LCCFqgEOXDzFs7TCyDdm0rduWpf2XSnFSVSgFB/8NC7r8Vpx0GAzj9kC7p6Q4EZWuTEdQpk2bxoABAwgMDCQ3N5cVK1awadMmfv75Z3Q6HVFRUURHRxMcHExwcDDR0dE4OzvzzDPPVFR+IUQ1kZCewKh1o8grzqODTwcW9V2Eq4Or1rEEQOZp+HEiJG0yt+u2gIfmQSP541Jop0wFyqVLlxg8eDCpqal4eHjQtm1bfv75ZyIiIgCYNGkSBQUFjBkzhqysLLp27cratWtxc3OrkPBCiOohPi2esevHUlBSQBe/Liy4bwHO9s5axxIlBvjl77DlPTAazINge78G3V4GOwet04la7q7nQSlvMg+KEDXLjos7GL9hPIXGQsLqhzH/vvk42TlpHUuc+QV+iIKME+Z2k3AYOBe8ZfZecWeqzDwoQgjxZ7ae30rUxiiKTEX0bNCTD/t8iN5WBsVrKv8KxP0N9n9hbrvUg/4x0Ob/ZJyJqFKkQBFCVIiNZzfyyuZXKDYV0yewD+/3fh8HWzltoBml4MAKWDsd8jPNyzq9AP1mgFMdLZMJcVNSoAghyt3aM2uZvGUyJaqEyEaRzOo1C3sbe61j1V4Zp+DHCZC8xdz2CTXPBNvwXm1zCXEbUqAIIcrVj0k/Mn3bdIzKyANBD/Buj3exs5G3Gk2UGGDbh7D1AzAWgZ0ThE+GsHFgKwWjqNrkXUMIUW7+e+q//O2Xv6FQDGo6iLe6vYWtja3WsWqn5K3muw5nnjS3m/WDB94HryBtcwlRSlKgCCHKxX9O/IeZO2aiUDze/HFev/d1bHR3PRekKKu8TFj7Ohz4ytx29YX7Z0GrR2UQrKhWpEARQty1r45+RczuGACeafkMU+6Zctt7cIkKoBQkfAlr/wYFVwAddB4Cfd8AJ0+t0wlRZlKgCCHuyvIjy3l/z/sAPB/6PK90fkWKk8p2+YT5dE7KNnPbp5V5JtjAezSNJcTdkAJFCHHHPjn4CX/f/3cAhrcZzksdXpLipDIVF5oHwG77EEzF5kGwfabCvWNkEKyo9qRAEUKUmVKKxQcWs/jAYgDGth/LqHajNE5VyyRtgh8mwpXT5nZwpHkQbJ1GmsYSorxIgSKEKBOlFPP3zWfZ4WUARHWMYmiboRqnqkXyMmDNdDi4wtx29YMBsyF0kAyCFTWKFChCiFJTSvHenvf4PPFzAF7r/BrPtXpO41S1hMkECV+YB8EWXgV0cM9wuO91cPTQOp0Q5U4KFCFEqZiUiZhdMaw4bv7LfXrX6TzV8imNU9US6cfMN/Y7u8Pc9msDD86HgE6axhKiIkmBIoT4UyZl4u2db/OfE/9Bh443w97kseaPaR2r5isugC3vwS9/Nw+CtXeBPtOg6yiwlbdvUbPJT7gQ4raMJiNvbH+DVadXYaOz4e3ub/Nw04e1jlXznVoPP06ErDPmdvMB8MB74BmoaSwhKosUKEKIWyoxlTB923R+Sv4JW50t0T2ieaDJA1rHqtmupcPPU+Hwf8xtN394YA60fFAGwYpaRQoUIcRNFZuKmbxlMnEpcdjp7JjTew4RjSK0jlVzmUyw7zNYNwMKs0FnA/eMhPumg95N63RCVDopUIQQNygyFvHq5lfZeG4j9jb2zA2fS3hguNaxaq5LR+B/UXB+t7ldvx08OA8adNQylRCakgJFCGGlsKSQCZsmsO3CNvS2eub1mUePBj20jlUzFeXD5tmwYwGYSsDB1XzZcJfhMghW1HryGyCEsCgoKWD8hvHsTN2Jo60jH/X9iHvr36t1rJrpZJx5EOzVs+Z2ywdhwBzwaKBtLiGqCClQhBAA5BfnM3b9WPZc2oOznTML+y6ks19nrWPVPLlp5kGwR74zt90DzFfntJTBx0L8nhQoQghyi3IZs24MCZcTcLV3ZXG/xbT3aa91rJrFZIK9/4B1b4EhB3S2cO9oCJ8Kelet0wlR5UiBIkQtl23IZlTcKA5nHsbNwY3YiFha122tdayaJe2QeRDshT3mtn9HeGieeTCsEOKmpEARohbLKsxiZNxIjl45iqfek9iIWEK8Q7SOVXMU5cGmGNixCJQRHNyg35vQeQjY2GqdTogqTQoUIWqpzIJMhscN52TWSbwcvVgauZTgOsFax6o5jv8MP70K2efM7dBH4P5Z4F5f01hCVBdSoAhRC13Ov8ywtcNIyk6inlM9lvZfShOPJlrHqhlyLsLqyXB0lbnt0RAGvg/N+2ubS4hqRgoUIWqZtLw0hq0dRkpOCr7Ovizrv4xG7o20jlX9mYwQvxTWvw1FueZBsN3GQe/J4OCidTohqh0pUISoRS5cu8DQNUO5cO0CDVwbsDRyKQFuAVrHqv4uJsAPUXBxv7kd0MU8E6yfDDYW4k5JgSJELXE25yxD1w4lLS+NQLdAlkUuo76rjIe4K4Zc2BgNu5aAMoHewzwIttOLYGOjdTohqjUpUISoBZKzkxm2ZhjpBek0dm/Msv7L8HH20TpW9XbsR/jpNci5YG63fgz6x4Cbr7a5hKghpEARooY7mXWS4WuHk1mYSTPPZnwS+Ql1nepqHav6yj4PP02C4z+a23Uaw8APoFk/TWMJUdNIgSJEDXbsyjFGrB1BliGLFnVaEBsZi5ejl9axqidjCez+GDa8C8V5YGMH3V+GXq+BvZPW6YSocaRAEaKGOpJxhBFxI8gpyqGVdys+jvgYD72H1rGqpwt7zTPBph00twPvNc8E6yOT2glRUaRAEaIGSkhPYPS60Vwrvka7eu1Y3G8xbg5uWseqfgpzYMM7sDsWUODoCREzocNgGQQrRAWTAkWIGmZP2h7Grh9Lfkk+nXw7sbDvQlzsZR6OMlEKEv8LP0+B3FTzsrZPQuS74FpP22xC1BJSoAhRg+xM3cn4DeMpKCmga/2u/L3P33G2d9Y6VvWSlWK+OufkGnPbqwkMnAtN+2ibS4haRgoUIWqIbRe2EbUxCoPRQPcG3ZkXPg9HO0etY1UfxmLYuQg2zYLifLCxhx4ToOcrYC/7UYjKJgWKEDXApnObmLhpIsWmYsIDw/mg9wc42DpoHav6OLfbPAg2/Yi53ag7PPgh1GuhaSwhajMpUISo5uJS4pi0eRIlqoSIRhHM7jkbe1t7rWNVDwVXYf1M2PMPQIGTF0S+De2fBZ1O63RC1GplGoYeExNDly5dcHNzw8fHh0ceeYTjx49b9VFKMWPGDPz9/XFyciI8PJwjR46Ua2ghhNlPST/x2ubXKFElDGg8gDm95khxUhpKwaH/wMJ7YM8yQJmLknF7oMNfpTgRogooU4GyefNmxo4dy86dO4mLi6OkpITIyEjy8vIsfebMmcPcuXNZsGAB8fHx+Pn5ERERQW5ubrmHF6I2+++p/zJl6xSMysjDTR8mpmcMdjZyUPRPXUmCL/4C3w6Fa5fAuxk8/z94ZBG4eGudTgjxK51SSt3pky9fvoyPjw+bN2+mV69eKKXw9/cnKiqKyZMnA2AwGPD19WX27NmMHDnyT7eZk5ODh4cH2dnZuLu732k0IWq0/5z4DzN3zESheCz4Md4IewMbnczLcVslRbD977DlPSgpBFu9eQBsjyiw02udTohqr7w/v+/qz63s7GwAvLzMU2cnJyeTlpZGZGSkpY9er6d3795s3779pgWKwWDAYDBY2jk5OXcTSYga71/H/kX0rmgAnmrxFFO7TpXi5M+kbIcfJsDlY+Z2UG/zIFjvptrmEkLc0h0XKEopJk6cSI8ePWjdujUAaWlpAPj6Wt/N09fXl5SUlJtuJyYmhrfeeutOYwhRqyw/spz397wPwHOhz/Fq51fRyXiJW8u/AuvehH3/NLed60L/aGj7hIwzEaKKu+MCZdy4cRw8eJBt27bdsO6Pb5hKqVu+iU6dOpWJEyda2jk5OQQGBt5pLCFqrKWHljJ/33wAhrcZzksdXpLi5FaUgoNfw5rpkJ9hXtbxeeg3A5zlZolCVAd3VKC89NJLrFq1ii1bthAQEGBZ7ufnB5iPpNSvX9+yPD09/YajKtfp9Xr0ejn/K8StKKVYfGAxiw8sBmBM+zGMajtKipNbyTgFP06A5C3mdr0Q8+mcRmHa5hJClEmZTlwrpRg3bhzfffcdGzZsICgoyGp9UFAQfn5+xMXFWZYVFRWxefNmunXrVj6JhahFlFLM3zffUpxEdYxidLvRUpzcTInBPAvs4jBzcWLnBH3fhJFbpDgRohoq0xGUsWPH8tVXX/Hf//4XNzc3y5gTDw8PnJyc0Ol0REVFER0dTXBwMMHBwURHR+Ps7MwzzzxTIV+AEDWVUoo58XP44ugXAEzqMonBoYM1TlVFJW8xD4LNPGVuN+sHD7wPXkG3f54QosoqU4GyeLH5r7jw8HCr5Z9++ikvvPACAJMmTaKgoIAxY8aQlZVF165dWbt2LW5ucqt3IUrLpExE74rm6+NfA/B619d5suWTGqeqgvIyzONMDq4wt1194f4YaPUXGQQrRDV3V/OgVASZB0XUdkaTkZk7Z/Ldye/QoeOtbm/xaPCjWseqWkwmSPgC4t6AgixAB12GQt83wNFD63RC1EpVah4UIUT5KjGV8Ldf/sYPST9go7Phne7v8FDTh7SOVbWkH4MfouDsDnPbtw08NA8COmuZSghRzqRAEaKKKDYVM2XLFNamrMVOZ0dMrxjub3y/1rGqjuIC8yywv/wdTMVg7wx9pkHX0WArb2VC1DTyWy1EFVBkLOLVza+y8dxG7Gzs+KD3B9zX8D6tY1Udp9bBj69A1hlzu/kAeOA98JQ5k4SoqaRAEUJjBqOBCRsnsPXCVhxsHJjXZx49A3pqHatqyE2Dn6fCke/MbTd/eGAOtHxQBsEKUcNJgSKEhgpKChi/YTw7U3fiaOvIR30/4t7692odS3smI+z5B6x/GwzZoLOBrqPMp3T0ckWgELWBFChCaCSvOI+x68ey99JenO2cWdB3AV38umgdS3upB82DYC/sNbf9O8CD88C/vYahhBCVTQoUITSQW5TL6HWjOXD5AK72rizut5j2Pu21jqUtwzXYFAM7F4MygoMb9HsTOg8BG1ut0wkhKpkUKEJUsmxDNiPjRnIk8wjuDu7ERsTSqm4rrWNp69iP8NMkyDlvboc+AvfPAvf6t32aEKLmkgJFiEp0pfAKI9aO4HjWcero6/BJ5Ce08GqhdSztZJ83FybHfzS3PRvCAx9A80htcwlRCoYSIztOZ9LY24XGdV20jlPjSIEiRCW5nH+Z4WuHczr7NN6O3iyNXEqzOs20jqUNYwnsWgIbo6E4D2zsoNtL0GsSODhrnU6IW8rOL2bj8XTiEi+x+cRlrhlKGBPelEn3t9Q6Wo0jBYoQlSAtL41ha4eRkpOCj7MPyyKX0dijsdaxtHF+D/wvCi4dMrcDu5oHwfqGaplKiFu6cLWAuCNpxB29xK6kK5SYfrtDjI+bHhe9fJRWBNmrQlSwC9cuMHTNUC5cu4C/iz9L+y8l0K0WTjBWcBXWzzRfPowCR0+ImAkdBoONjcbhhPiNUoojF3OIS7xEXOIlElNzrNYH+7gSEepLZCs/2jbwwMZG5uSpCFKgCFGBUnJSGLZ2GGl5aQS6BbIschn1XWvZwE+l4PC35gnX8tLNy9o9DRFvg2s9bbMJ8atio4ndyVcsRcmFqwWWdTY66NzIi4hQX/qF+hIk400qhRQoQlSQ01dPM3ztcC4XXCbII4ilkUvxcfbROlblyjxtnqI+aaO57d0MHvwQgnppm0sI4JqhhM3HL7M2MY2Nx9LJKSyxrHO0t6FncD0iQn3p29IHb1e9hklrJylQhKgAx68cZ0TcCK4UXiG4TjCfRHyCt5O31rEqT4kBfpkPW94HowFs9dDrVej+MtjJG73QTlp2IXFHzUdJdp7OpMhosqzzdnGgb4gPEaF+9GhWFycHmX9HS1KgCFHOjmQcYUTcCHKKcgj1DuXjfh/j6eipdazKk7wVfpgAmSfN7SbhMHAueDfVNJaonZRSHL+US9yRS8QdvcTB89lW65vUdSEi1JeIUF86NKyDrYwnqTKkQBGiHCWkJzB63WiuFV+jXb12LOq3CHcHd61jVY68DFj7Ohz4l7nt4gP9o6HN/8mN/USlKjGaiD+TZR5PcjSNc1d+G0+i00GHQE8iQv2ICPWlmY+rhknF7UiBIkQ5iU+LZ+z6sRSUFNDJtxML+y7Exb4WDKYzmWD/5xD3BhReBXTm6en7vgFOnhqHE7VFnqGELScuE5d4iQ3H07maX2xZ52BnQ49mdc3jSUJ88HFz1DCpKC0pUIQoB79c+IWXN76MwWggrH4Y8++bj5Odk9axKt6lI+bTOed2mdu+beCheRDQWdNYonZIzylk3dF04hLT+OV0JkUlv40nqeNsz30tzaduejWvi7ODfNxVN/IdE+IubTq3iYmbJlJsKqZ3QG8+CP8AvW0NHwhalAebZsGOheYb+9m7wH3T4Z6RYCtvK6JiKKU4mX6NuMRLrE28xIFzV63WN/J2JiLEXJR0alQHO1uZX6c6k3cSIe7CmjNrmLJlCiWqhIhGEczuORt7W3utY1WsYz/B6kmQfc7cDnnIfGM/jwBtc4kaqcRoYm+KeTzJuqOXOJOZb7W+XaAnkb8Ocg32cUUn451qDClQhLhD/zv9P17/5XVMysQDQQ/wbo93sbOpwb9SV8/B6sm/3djPoyE8MAdaDNA2l6hx8gwlbD15mbWJl9h4LJ2s348nsbWhWzNv86RpIb74ust4kpqqBr+bClFx/n3i37y9420Uir8E/4U37n0DW5saOmeCsRh2Ljaf0rl+Y7+wcdB7EjjUgkHAolKk5xay/qj5JnzbTmVYjSfxdLbnvhY+9Av1pVfzerjKvW9qBfkuC1FGXx79klm7ZwHwdMunmXLPFGx0NfRc99ld5kGw6UfM7cB7zTPByo39xF36/XiSuMRLJPxhPElDL2fL/CSdZTxJrSQFihBlsPTQUubvmw/Ai61eZEKnCTXznHf+FVj3Juz7p7nt5GW+sV/7Z+XGfuKOlRhN7EnJYl2iedK0FBlPIm5DChQhSkEpxYKEBcQejAVgVLtRjGk3pua9gSoFCV9B3N8gP9O8rMNfod9McKlFU/WLcnPNUMLW28xP0r2pNxGhfvQL8cFHxpOI35ECRYg/oZTivT3v8Xni5wBM6DSBIa2HaJyqAqQfhR8mwtnt5na9EPPpnEZh2uYS1c6lnELW/Xq/m+2nrO934+lsz30tfYgM9aVncD1cZDyJuAX5yRDiNkzKxLs73+WbE98AMPWeqTwT8ozGqcpZUR5sngM7FoCpBOydIXwK3DsGavol06JcKKU4lpZruRT4j/e7kflJxJ2QAkWIWygxlfDm9jdZdXoVOnTM6DaDvwT/RetY5evYT+ZLh7PPmtstHoABs8Gzoba5RJVXbDSxO/mKZZDrhas33u+mX6gvkaG+NK0n40lE2UmBIsRNFBuLmbJ1CmtT1mKrs+XdHu8ysMlArWOVn6tnYfWU381pEggD5kDLB7TNJaq07IJiNp+4zLrES2w8nk5uYYllnaP9b/e7ua+lL/XcavhsyqLCSYEixB8YjAYmbprIlvNbsLex573e79G3YV+tY5WPkiLYudB8Sqc43zynSbeXoNdrMqeJuKlzV/JZf9R81c2upCuUmJRlnbeLA31DfIgI9aNHs7o4OdTQuYCEJqRAEeJ38ovzGb9xPLtSd6G31TO/z3y6N+iudazycWYb/PgKXD5mbjfqDgM/AJ8QbXOJKsVkUhy6kG0Z5HosLddqfTMfV/qF+BIR6kP7wDrY2sipG1ExpEAR4lc5RTmMXTeWhMsJONs5s6DvArr4ddE61t27lg5xb8CBf5nbznWh/7vQ9knzYAFR6xUWG9l+OoO4xHTWH71Eeq7Bss5GB10ae9EvxJd+ob4E1ZUjbaJySIEiBJBVmMXIuJEcvXIUNwc3lvRbQtt6bbWOdXdMRtj7KayfCYXZgA46vwh93wCnOlqnExrLuGZgw7F01iVeYuvJDAqKjZZ1Lg629Gpej4hQX/q08KGOi4OGSUVtJQWKqPXS89MZsXYEp7NP4+XoRWxELC28Wmgd6+5c2Ac/ToSL+83t+u1g4IcQ0EnbXEIzSilOpV8j7ugl1h9NZ9/ZLNRvw0mo7+FoOUpybxMv9HYynkRoSwoUUatduHaB4WuHcy73HD7OPiyNXEqQR5DWse5cQZb5iMmeTwEFeg/o+zfoPARq6s0MxS0VG03En7nC+qPprLvJ1PJtGnj8WpT4EFrfXS4FFlWKFCii1krOTmb42uFcyr9EgGsAn0R+QoBbgNax7oxlivo3ID/DvKzNExD5Drj5aptNVKrrlwKvP3qJjcfSyfndpcAOtjZ0a+ZNvxBf+ob4UN/DScOkQtyeFCiiVjp25Rgj40ZypfAKTTyaEBsRi69LNf0gTzsMP70KZ3eY2/VawgPvQ1BPbXOJSpOSmce6o+YBrruTrS8F9nJx4L6WPvQL8ZGp5UW1Ij+potZJSE9gzLox5BbnEuIVwscRH1PHsRoOGi3MgU0xsOtjUEaZor4WMZoUCeeyLFfdnEy/ZrW+mY8rfUPM97uRS4FFdSUFiqhVdlzcwcsbX6agpICOPh1Z0HcBbg5uWscqG6Xg8LewZjpcSzMvCx0E/aPBo5qeohJ/6vpdgdcdTWfT8XQy84os62xtdHRpXMdyV+BG3nIpsKj+ylygbNmyhffee4+9e/eSmprK999/zyOPPGJZr5TirbfeIjY2lqysLLp27crChQtp1apVeeYWoszWp6zntS2vUWwqppt/Nz4M/xBne2etY5VN+lH46TU4s9Xc9moKD8yBZv20zSUqxPmsfMsA111JV6zuCuzuaEd4Cx/6hvgQ3twHD2c5aiZqljIXKHl5ebRr144XX3yRxx577Ib1c+bMYe7cuXz22Wc0b96cd955h4iICI4fP46bWzX7S1XUGP899V/e2P4GJmUiolEEs3rOwsG2Gs3tYMiFzbNh52LzHYftnKDnK9B9PNjJPU9qCqNJceD8Vdb/einwH2dxbeztTN9fB7h2aeyFvdwVWNRgZS5QBgwYwIABA266TinFvHnzmD59On/5i/mur8uXL8fX15evvvqKkSNH3l1aIe7Al0e/ZNbuWQA82uxR3gh7AzubanJ28/rpnLWvQ26qeVnLB+H+GLnjcA1xzVDCtpPmUzcbj1mfurHRQedGXvQN8aFviC9N67nIpcCi1ijXd+nk5GTS0tKIjIy0LNPr9fTu3Zvt27fftEAxGAwYDL9Nq5yTk1OekUQtppRiyYElLDqwCIDBoYN5rfNr1ecN/lKi+XROyjZzu04QPPAeBEdom0vctes34Ft/LJ2dSZkUG3+76sbN0Y7ezevRL8SX8Bb18HSuRkf6hChH5VqgpKWZB+z5+lpfrunr60tKSspNnxMTE8Nbb71VnjGEwKRMzImfw5dHvwRgTPsxjGo7qnoUJ4XZsGnWb1fn2DlBr1cg7CWwd9Q6nbgDRpNi39ks1h9NZ8OxS5y4ZH3VjZy6EeJGFXKc+48fAkqpW34wTJ06lYkTJ1raOTk5BAYGVkQsUUsUm4p545c3+CHpBwCm3jOVZ0Ke0ThVKZhMcHAFxL0JeenmZSEPma/OkdM51U52QTFbTlxmw7F0Nh5P52p+sWWdrY2Ozo3q0C/El/tCfGhaz1XDpEJUTeVaoPj5+QHmIyn169e3LE9PT7/hqMp1er0evV4G+YnyUVhSyGubX2PT+U3Y6ex4p8c7DGwyUOtYf+5igvl0zvnd5rZ3MxgwW67OqUaUUpy+fI0Nx9JZfzSdPSlZGH83YZqHkz3hLepxX0u56kaI0ijXAiUoKAg/Pz/i4uLo0KEDAEVFRWzevJnZs2eX50sJcYOcohxeWv8S+9L3obfVMzd8Lr0Cemkd6/byMmHD27D3M0CBvQv0nmSebM1Oxh5UdYYSI7uSrrDhWDobjqVz9or1vW6CfVy5L8SH+1r40KlRHezk1I0QpVbmAuXatWucOnXK0k5OTiYhIQEvLy8aNmxIVFQU0dHRBAcHExwcTHR0NM7OzjzzTDU4xC6qrcv5lxm1bhQnsk7gZu/GR30/opNvFb5zr7EE9vwDNr5jHnMC0Pr/IPJtcPfXNpu4rfScQjYeNx8l2XYqg/wio2Wdg60NXZt4cV9LH/q29KWhdzWbZ0eIKqTMBcqePXvo06ePpX19/Mjzzz/PZ599xqRJkygoKGDMmDGWidrWrl0rc6CICnM25ywj4kZw4doF6jrVZUm/JbTwaqF1rFtL3gqrJ0P6EXPbt415srVG3bTNJW7q+twkm46ls+F4OocvWF9p6OOmp08LH+4L8aFHs7pyrxshyolOKaX+vFvlycnJwcPDg+zsbNzd3bWOI6q4o5lHGbVuFFcKrxDoFsjHER8T6FZFB1lfPWuezyTxv+a2Ux2473Xo9CLY2GqbTVjJzi9my8nLbDyWzuYTl63mJtHpoG2AJ31b+nBfSx9C67tjI/e6EaLcP7+l1BfV1o6LO5iwaQJ5xXm09GrJ4n6LqetUV+tYNyrKh1/mwS/zoaQQdDbQeQj0mQ7OXlqnE5gHuB6/lMuGY+lsOnaZvWetB7i66e3o1bwefVr6EN6iHnVdZWC/EBVNChRRLf2U9BPTf5lOiamEe/zuYV6feVXvpn9KwaH/wLo3IeeCeVnjnnD/LPBrrW02QZ6hhO2nM81FyfF0UrMLrdYH+7iar7hp4UPnxnVkbhIhKpkUKKLa+eeRf/LenvcA6N+4P9E9oqvefXUu7IXVU367bNijoXkAbOgg8zkCUenMlwHnsel4OpuOX2Z3svXN9xztbejWtC59WtQjvIUPgV4ywFUILUmBIqoNkzIxd89clicuB+DZkGeZ1GUSNroq9JdtzkVYPxMO/MvctneBnhMgbBzYO2mbrRbKLyphx+lMNh2/zMbj6ZzPKrBa39DL2VyQtPQhrIk3jvYyFkiIqkIKFFEtGIwGpm+bzpozawCY0GkCL7Z6sepMXV+UD9s/Mo81Kf51Lox2T0PfN8G9/m2fKsrP9cnSNh2/fNOjJNcvAw5v4UOfFvUIqis33xOiqpICRVR52YZsxm8Yz770fdjZ2DGz20weavqQ1rHMTCY49G9Y/9Zv40wCu0L/GAiowvOw1CDXDCVsP5XB5hPmouTCVeujJAF1nAhvUY8+LXwIa+qNs4O87QlRHchvqqjSLly7wOh1o0nOTsbV3pV5febRtX5XrWOZndkGa6ZDaoK57dEQIt6CVo/KOJMKpJTiWFoum09cZvPxy+xJuWJ1N+DrR0l6NzePJWlaT46SCFEdSYEiqqxDlw/x0oaXyCzMxNfZl0X9FtG8TnOtY0HGKYh7A47/aG47uJnHmdw7Vu42XEGy8orYeiqDLScus+XEZdJzDVbrG3o5E96iHuEt6nFvEzlKIkRNIL/FokqKS4lj2tZpFBoLaV6nOYv6LsLX5eY3nKw01y7D5lmw51NQRtDZQqcXIHwquNbTNlsNU2I0kXDuKltOXGbzyQwOnr/K76eUdLK3JaypN72C69K7hQ9BdV20CyuEqBBSoIgqRSnFp0c+5cO9HwLQs0FP3uv9Hi72Gn4AFeXDzoWwbT4U5ZqXBfc3XzZcrwpPqV/NnM3MZ8vJy2w9eZntpzLJNZRYrW/p50av5vXoFVyPLkF10NvJFTdC1GRSoIgqo9hYzLu73uXbk98C8HTLp5nUZRJ2Nhr9mBpLYP/nsHk25Kaal9Vvby5Mgqr4XZKrgZzCYnaczmTryctsPZlBSqb1nYA9ne3p3qwuvX8tSvw85PSZELWJFCiiSsgqzGLiponsubQHG50Nk7pM4tmQZ7UJoxQcXWWezyTz1zt3ezY0XzLc6i9gU4XmXalGio0mDpy7ytaTGWw7lUHCuatW08nb2ejo2LAOvZrXpWdwPVo38MBW7nEjRK0lBYrQ3Mmsk7y04SUuXLuAi70Lc3rNoVeARkcokjaZC5MLe81tZ2/oNQk6vwh2cv+Vsrg+J8m2XwuSnUlXuPaH0zZBdV3oGVyXXsH1uLepN65yJ2AhxK/k3UBoatO5TUzeMpn8knwC3QL56L6PaOrZtPKDnN9jnsskeYu5be8C3caZZ4B1lLtql1ZadiG/nMowP05ncCnH+mqb66dtejarS4/gugTUkenkhRA3JwWK0IRJmfjk4CcsTFiIQnGP3z180PsDPB09KzdI2mHYGP3bJcM29uY7Dfd8Bdw0vmqoGsjOL2ZHUgbbT2fyy6kMTl/Os1rvYGfDPY29zEVJcF1C67tjI6dthBClIAWKqHS5RblM2zaNTec2AfBkiyeZfM9k7G3sKy9E+lHYNAsSV5rbOhto9wyETzaPNxE3lWcoIf7MFXaczmT76UwOX8y2uvxXp4O2DTzo1qwuPZrVpVOjOnJ/GyHEHZECRVSq01dPE7UxijM5Z3CwceD1e1/n0eBHKy/A5ROwZQ4c+g+gAJ155tfwqVCvCkwCV8UUFhvZm5LFjtOZ7EjK5MC5q5T8bmArQNN6LnRvVpduTesS1sQbD+dKLDSFEDWWFCii0qxOXs2M7TPIL8nHz8WPeeHzaFW3VeW8ePpR2PIeHP4Oc2EChDxkLkx8KylDNVBYbGTf2Sx2Jl1hZ1ImCWevWt1sD6CBpxPdm3kT1tSbsCZ15fJfIUSFkAJFVDiD0cCc3XP45sQ3ANzjdw/v9X4PL0evin/xtEOw5X1I/C+WwqTlg9B7EtRvV/GvX8XlF5WwL+Uqu5Mz2Zl85aYFia+7nrAmvxUkDb1lYKsQouJJgSIq1Nmcs7yy+RWOXTmGDh3D2gxjTPsxFT/52tldsPUDOLnmt2UhD0Ov16B+24p97Sosp7CYvSlZ7E6+wq6kTA6ez77hlI2Pm56wpt50DfLm3iZeBNWVm+0JISqfFCiiwvyY9CNv73ybvOI86ujrENMzhu4NulfcCyoFp9bDtg8hZZt5mc4GQh+BXq/WylM56bmF7DljLkjiz1zhaGoOf6hHqO/hSNcgL7o28ebeJt409naWgkQIoTkpUES5yynK4Z2d77A6eTUAHX06MqfXnIq72V9JERz+FrZ/BOlHzMts7KH909A9Crw1mFdFA+aJ0fLYc+YKe1Ky2HPmCmf+MH08QGNvZ7o09uKeIC/ubeJNQB0nKUiEEFWOFCiiXMWnxTN923RS81Kx1dkysu1IhrcdXjGndAqyYN8/YdfHkHPBvMzeBTo9D2FjwSOg/F+zCiksNnLwfDZ7U7LYm3KFvSlZZOUXW/XR6aClnztdGtehc2MvugZ54esug1qFEFWfFCiiXBSUFLBw/0L+mfhPFIpAt0BiesbQrl4FDETNOAW7lkDCV1D868Rgrr7QdaR5kjWnOuX/mhpTSnHhagH7zl5lX0oW+89mceRizg3jR/R2NrQP9KRz4zp0buRFx0Z18HCSy36FENWPFCjirsWnxTNj+wzO5p4F4NFmjzL5nsm42LuU34uYjHBqHez+BE7F/bbcpxXcOxraPlGj7pWTZyjh4Pls9p/LIuHsVfafu8rlXMMN/eq56encqA6dGtWhY6M6tPb3wMFObmYohKj+pEARdyyvOI8P937I18e/BsDH2Yc3w94s3xv95WXC/s9hzz/gasqvC3XQ/H5zYRLUy3weoxorNpo4cSmXA+eyOXDuKgnnrnIyPfeGwax2NjpC/d3p2LAOHRp60rFhHRk/IoSosaRAEWWmlOLnMz/zfvz7pBekA/BY8GO80vkV3Bzc7v4FTCY4swX2LodjP4CxyLzc0QM6DDafxqmmA1+NJkXS5WscupDNwfPZHDh/lcSLORhKTDf0re/hSIeGnnQIrEP7hp609vfAyUGmjRdC1A5SoIgyOX31NNG7otmdthuAANcA3uz2JvfWv/fuN559Hg78C/Z/AVlnfltevz10GQatHwOH6jNJWInRxOnLeRy+kM2hC9kcvpBNYmoO+UXGG/q6OdrRLsCTdoEev/7rKYNZhRC1mhQoolSyCrOIPRjLimMrKFEl6G31DGszjBdbv4je9i7GfhiuwdH/wYGvIHkrltle9e7Q5nHo+Bz4ty+PL6FCFRQZOZaWw5GLOSSmmv89lnrzIyPODra08nenTQNzQdI2wJNGXs5yl18hhPgdKVDEbRWUFPBF4hf84/A/uFZ8DYA+gX2Y1GUSAW53eBlvicE84PXwt3B8NRT/bq6ORj2g/TPQ6hFwKMdBtuXk+tU0x1JzOZaWw9HUXI6m5XAmI++GMSMArno7Qv3dfy1IPGjTwIMm9VyxlWJECCFuSwoUcVMGo4H/nvovSw4s4XLBZQBaerVkQscJdGvQrewbLC6EpE1wdBUc/QEM2b+t82oC7Z4xX4lTp1H5fAHlIPOagROXrnEyPZdjabkcT8vlRFouuYaSm/av56YntL65GAmpby5IGsqRESGEuCNSoAgr+cX5/PvEv1l+ZLmlMGng2oBxHcbxQNAD2OjKcAlrYbb5SMnRH+DkWii69ts6N39o/Rfzw7+jZlfiKKXIuFbEqfRrnErP5WT6NU7+WpRkXCu66XPsbXU0redKSH13Wvq5mf+t74aPm4wZEUKI8iIFigAgoyCDfx//N18d+4qrhqsA+Dr78mLrF3m8+eM42Dr8+UaUgsxTcGINnPgZzu4A0++ONrj5Q8hDEDoIGoaBTeXN12EoMXLuSj5Jl/NIysjjdPo1Tl++xunLeWQXFN/yeYFeTrTwdSPY142Wfm609HMnqK6LzDUihBAVTAqUWkwpxf70/aw4toK4s3GU/FpMBLoFMrT1UB5u+jD2tn8yC2n+FfOpm6SNcHoTZJ+1Xu8dDC0Hmu8k7N+hQouSYqOJ81kFpGTmcSYjjzOZ+SRn5HEmM49zV/JvOkYEzAdvGno5E+zjSlMfV4J93Aj2caWZjysuevkVEUIILci7by2Unp/O6uTVrDq9ihNZJyzL29Vrx7MhzxLRKOLW987Jy4SUXyBlu/mOwWmHsVx5A+ab9DXubp5ILTiyXOcrUUqRXVDMuSsFnL2Sz7msfM5eyedsZj4pV/K4eLUQ462qEMDFwZYm9VwJqutC03quNPVxoZmPK429XXC0l/lFhBCiKpECpZa4WniVTec38WPSj+xK3YX6tahwtHXkgSYP8GSLJwn1DrV+kskEl4/B+d1wLt78b8aJGzfuEwpN+kDTPuZTN3rXO8poMikyrhm4cLWAi1cLuXA1n4tXCzmfVcD5rHzOZxVw7RYDVK9ztLehoZczjb1dCKrrQuO6LjT2dqFJPRd83PQy66oQQlQTUqDUUEopUnJS2Hx+MxvPbWR/+n5M6rc5OTr6dGRgk4H0b9wfD70HGIsh/SikHoSL+yE1AdIOWQ9sva5eS2jU3XykpFEPcPP90zzFRhOXcw1cyinkUk4hqdmFpOUUkpZdSOrVQi5mF3App5Bi462PgFxX11VPQy8nGno509DLmYBfC5JG3s5ShAghRA0hBUoNkpaXxu603exK3UV8WjypealW61vUaUG/hn0Z6NWWwIIcuHwcjkyC9CPm/xtvctWKvQs06AgBXSDwHvO/LnUBcxF0zVBCRkYel3MNZFwzcDnX/EjPLSQ910B6jvn/mXlFqD+vPbDRga+7Iw08nfD3dKJBHfO/gXWcCKjjTANPJ5nuXQghagEpUKqpq4VXOZZ1jMMZhzmccZhDGYdIz0+36mOns6WzcwPCbdwJLyymwbnzcOCNmxcigHJwxVgvlGvebchyDyHVJYTztg3IyDeSlVdE5oEirmxPIvPaMa7kFZFxzXDTmVJvxc5Gh4+bHj8PR/w8HPF1d8TP3ZH6nk74e5j/9XXTY2crV8gIIURtV2EFyqJFi3jvvfdITU2lVatWzJs3j549e1bUy9VIRpORtPw0zuac5VzuOZKzkzl15Tinrp4iw5B1Q39bBaElJu7Ju8Y9hYW0LzTgrJJv6Fesc+CSQyDnbANJIoAjxkASDA04luuJKccGTl/vmf3r4/ZcHGyp66annqueuq56fNzN//dx1+Pj5oiPux5fd0e8nB1k0jIhhBClUiEFytdff01UVBSLFi2ie/fufPzxxwwYMIDExEQaNmxYES9ZrZiUidyCLK7knOVS1lnSci5yOTeVy/mXuFSYweWiLDKM17isCinR3fq8SIPiElobDLQxFNHaUERIURHOv55HKVAOXFD1OaN8SVF+pCgfzipfTil/Lqq6mApufZTCzdEOLxcHPJ0dqONsj5eLA94uDni56PFysaeuqx5vVz3eLg54uzrg7CAH4oQQQpQvnVKlGRlQNl27dqVjx44sXrzYsiwkJIRHHnmEmJiY2z43JycHDw8PsrOzcXd3L+9oVpRSlJgUxSVGiouKMBQXUmQoxGAooLDoGoWGAgoNeRiK8ikqzsdQlIeh2Pz/4pJ8iooLKDLmU1xSSJGpgCJjIUUmA0XKQLEqwkARhRRj0JVQqDOSZ6O4ZqPItdFhKuVATnulCCguoWFJCQ2LiwkuKsa/SIebwYVcVYdLqg6pyotLyos0VYfzqh7nVT0MDp64OTrg4WSPu5Md7o72uDna4eFk/+sy8791nB3wdLbH89d/PZzssZdTLEIIIcqovD+/y/1P36KiIvbu3cuUKVOslkdGRrJ9+/Yb+hsMBgwGg6Wdk5NTrnk2713J2/tfR+kUJjA/dNzw/xKdDiOg7vQKEJtfH7el+/Vh5mIyUcdowsOow9Vog6vRHlejIy4mFxzxxFHnjd7enxK9N0ZXL4xO3lxw9SXD2R1XRztc9Ha46e2o72iHq978cHe0x9XRTm5GJ4QQolor9wIlIyMDo9GIr6/1pae+vr6kpaXd0D8mJoa33nqrvGNYFJcUc8neujAoK3ulsFVgr8AeHXbK/HDABntssMcWe2xx0Nmj1zmgt9HjYKPHydYFZ3sXXOzccNZ74O7ojYeLD3Xc/PFy98XV3QdnR0cc7WxkYKgQQgjxOxU2eOCPc1EopW46P8XUqVOZOHGipZ2Tk0NgYGC55ejYsjfRhijsbB2wt7fHztYBO1t7HBwccbBzwMHeEb2DIw4OehwdnLG3d8LWVo+tjR32tvbY6exkXg0hhBCikpV7gVK3bl1sbW1vOFqSnp5+w1EVAL1ej16vL+8YFl4ePjzUa2iFbV8IIYQQ5a/czys4ODjQqVMn4uLirJbHxcXRrVu38n45IYQQQtRAFXKKZ+LEiQwePJjOnTsTFhZGbGwsZ8+eZdSoURXxckIIIYSoYSqkQHnyySfJzMxk5syZpKam0rp1a3766ScaNWpUES8nhBBCiBqmQuZBuRuVOQ+KEEIIIcpHeX9+y7WtQgghhKhypEARQgghRJUjBYoQQgghqhwpUIQQQghR5UiBIoQQQogqRwoUIYQQQlQ5UqAIIYQQosqRAkUIIYQQVU6F3c34Tl2fNy4nJ0fjJEIIIYQoreuf2+U1/2uVK1Byc3MBCAwM1DiJEEIIIcoqMzMTDw+Pu95OlZvq3mQycfHiRdzc3NDpdOWyzZycHAIDAzl37pxMn4/sjz+S/WFN9oc12R/WZH9Yk/3xm+zsbBo2bEhWVhaenp53vb0qdwTFxsaGgICACtm2u7t7rf8B+j3ZH9Zkf1iT/WFN9oc12R/WZH/8xsamfIa3yiBZIYQQQlQ5UqAIIYQQosqpFQWKXq/nzTffRK/Xax2lSpD9YU32hzXZH9Zkf1iT/WFN9sdvyntfVLlBskIIIYQQteIIihBCCCGqFylQhBBCCFHlSIEihBBCiCpHChQhhBBCVDlSoAghhBCiyqnxBcqiRYsICgrC0dGRTp06sXXrVq0jaSImJoYuXbrg5uaGj48PjzzyCMePH9c6VpURExODTqcjKipK6yiauXDhAn/961/x9vbG2dmZ9u3bs3fvXq1jaaKkpITXX3+doKAgnJycaNKkCTNnzsRkMmkdrVJs2bKFhx56CH9/f3Q6HStXrrRar5RixowZ+Pv74+TkRHh4OEeOHNEmbCW43f4oLi5m8uTJtGnTBhcXF/z9/Xnuuee4ePGidoEr2J/9fPzeyJEj0el0zJs3r8yvU6MLlK+//pqoqCimT5/O/v376dmzJwMGDODs2bNaR6t0mzdvZuzYsezcuZO4uDhKSkqIjIwkLy9P62iai4+PJzY2lrZt22odRTNZWVl0794de3t7Vq9eTWJiIh988EG53E+jOpo9ezZLlixhwYIFHD16lDlz5vDee+/x0UcfaR2tUuTl5dGuXTsWLFhw0/Vz5sxh7ty5LFiwgPj4ePz8/IiIiLDc7LWmud3+yM/PZ9++ffztb39j3759fPfdd5w4cYKHH35Yg6SV489+Pq5buXIlu3btwt/f/85eSNVg99xzjxo1apTVspYtW6opU6ZolKjqSE9PV4DavHmz1lE0lZubq4KDg1VcXJzq3bu3evnll7WOpInJkyerHj16aB2jyhg4cKAaMmSI1bK//OUv6q9//atGibQDqO+//97SNplMys/PT82aNcuyrLCwUHl4eKglS5ZokLBy/XF/3Mzu3bsVoFJSUionlIZutT/Onz+vGjRooA4fPqwaNWqkPvzwwzJvu8YeQSkqKmLv3r1ERkZaLY+MjGT79u0apao6srOzAfDy8tI4ibbGjh3LwIED6devn9ZRNLVq1So6d+7M448/jo+PDx06dOCTTz7ROpZmevTowfr16zlx4gQABw4cYNu2bTzwwAMaJ9NecnIyaWlpVu+ter2e3r17y3vrr7Kzs9HpdLX2CKTJZGLw4MG89tprtGrV6o63U+XuZlxeMjIyMBqN+Pr6Wi339fUlLS1No1RVg1KKiRMn0qNHD1q3bq11HM2sWLGCffv2ER8fr3UUzSUlJbF48WImTpzItGnT2L17N+PHj0ev1/Pcc89pHa/STZ48mezsbFq2bImtrS1Go5F3332Xp59+Wutomrv+/nmz99aUlBQtIlUphYWFTJkyhWeeeabW3t149uzZ2NnZMX78+LvaTo0tUK7T6XRWbaXUDctqm3HjxnHw4EG2bdumdRTNnDt3jpdffpm1a9fi6OiodRzNmUwmOnfuTHR0NAAdOnTgyJEjLF68uFYWKF9//TVffPEFX331Fa1atSIhIYGoqCj8/f15/vnntY5XJch7642Ki4t56qmnMJlMLFq0SOs4mti7dy/z589n3759d/3zUGNP8dStWxdbW9sbjpakp6ffUPnXJi+99BKrVq1i48aNBAQEaB1HM3v37iU9PZ1OnTphZ2eHnZ0dmzdv5u9//zt2dnYYjUatI1aq+vXrExoaarUsJCSkVg4oB3jttdeYMmUKTz31FG3atGHw4MFMmDCBmJgYraNpzs/PD0DeW/+guLiYJ554guTkZOLi4mrt0ZOtW7eSnp5Ow4YNLe+tKSkpvPLKKzRu3LhM26qxBYqDgwOdOnUiLi7OanlcXBzdunXTKJV2lFKMGzeO7777jg0bNhAUFKR1JE317duXQ4cOkZCQYHl07tyZZ599loSEBGxtbbWOWKm6d+9+w2XnJ06coFGjRhol0lZ+fj42NtZvj7a2trXmMuPbCQoKws/Pz+q9taioiM2bN9fK91b4rTg5efIk69atw9vbW+tImhk8eDAHDx60em/19/fntddeY82aNWXaVo0+xTNx4kQGDx5M586dCQsLIzY2lrNnzzJq1Cito1W6sWPH8tVXX/Hf//4XNzc3y18/Hh4eODk5aZyu8rm5ud0w/sbFxQVvb+9aOS5nwoQJdOvWjejoaJ544gl2795NbGwssbGxWkfTxEMPPcS7775Lw4YNadWqFfv372fu3LkMGTJE62iV4tq1a5w6dcrSTk5OJiEhAS8vLxo2bEhUVBTR0dEEBwcTHBxMdHQ0zs7OPPPMMxqmrji32x/+/v783//9H/v27eOHH37AaDRa3l+9vLxwcHDQKnaF+bOfjz8WaPb29vj5+dGiRYuyvdDdXWBU9S1cuFA1atRIOTg4qI4dO9bay2qBmz4+/fRTraNVGbX5MmOllPrf//6nWrdurfR6vWrZsqWKjY3VOpJmcnJy1Msvv6waNmyoHB0dVZMmTdT06dOVwWDQOlql2Lhx403fL55//nmllPlS4zfffFP5+fkpvV6vevXqpQ4dOqRt6Ap0u/2RnJx8y/fXjRs3ah29QvzZz8cf3ellxjqllCpbSSOEEEIIUbFq7BgUIYQQQlRfUqAIIYQQosqRAkUIIYQQVY4UKEIIIYSocqRAEUIIIUSVIwWKEEIIIaocKVCEEEIIUeVIgSKEEEKIKkcKFCGEEEJUOVKgCCGEEKLKkQJFCCGEEFXO/wPOiHP91b9S7QAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "def Hattori2010855(a, b, c, t):\n",
    "    \"\"\"\n",
    "    Implements the Logistic curve model\n",
    "    Hattori2010855\n",
    "\n",
    "\n",
    "    a - (1/h)\n",
    "    b - (1/um)\n",
    "    c - (-)\n",
    "    t - (hr)\n",
    "    \"\"\"\n",
    "    \n",
    "    import numpy as np\n",
    "    cum_loss = (a/b)*t - (1/b)*np.log((1+c)/(1+c*np.exp(-a*t)))\n",
    "    \n",
    "    return cum_loss\n",
    "    \n",
    "    \n",
    "# From Table 6    \n",
    "fig, ax = plt.subplots()\n",
    "\n",
    "t = np.linspace(0, 12, 100) #hr\n",
    "mde = Hattori2010855(0.675, 0.323, 14.0, t)\n",
    "ax.plot(t, mde, label=\"SUS304 $A_{p–p} = 40 μm$\")\n",
    "\n",
    "mde = Hattori2010855(1.33, 0.317, 42.0, t)\n",
    "ax.plot(t, mde, label=\"SUS304 $A_{p–p} = 50 μm$\")\n",
    "\n",
    "mde = Hattori2010855(2.07, 0.316, 107, t)\n",
    "ax.plot(t, mde, label=\"SUS304 $A_{p–p} = 60 μm$\")\n",
    "ax.set_xlim(right=14)\n",
    "ax.set_ylim(top=70)\n",
    "ax.legend()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "024320a0-1946-4018-97db-514d154d287c",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "a - thickness of eroded layer\n",
    "K = constant in the stress-strain relationship\n",
    "l, L - depth of hardened layers\n",
    "n = exponent in the stress-strain relationship\n",
    "r - radius\n",
    "x = distance from the surface\n",
    "e = strain\n",
    "o - stress\n",
    "theta - metallurgical shape factor\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "rho - density\n",
    "S_o - exposed area\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "77f821fa-a1d5-42be-ac18-dd152c97a62c",
   "metadata": {},
   "source": [
    "\n",
    "Berchiche2002601\n",
    "\n",
    "\n",
    "\n",
    "material characterization\n",
    "- tensile test\n",
    "n = exponent in the stress-strain relationship\n",
    "K = constant in the stress-strain relationship\n",
    "sigma = stress\n",
    "sigma_e - elastic stress\n",
    "sigma_r - rupture stress\n",
    "sigma_s - surface stress\n",
    "\n",
    "microhardness measurements\n",
    "L - depth of hardened layers\n",
    "theta - metallurgical shape factor\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a3145e7e-725c-4a57-a05b-f7b864dccb5d",
   "metadata": {},
   "source": [
    "Berchiche2002601 - Cavitation Erosion Model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "id": "71641ee9-c98d-4158-9db7-89ac55cc1a97",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle o{\\left(0.47 \\right)} = 1017008914.03609$"
      ],
      "text/plain": [
       "Eq(o(0.47), 1017008914.03609)"
      ]
     },
     "execution_count": 44,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import sympy as sp\n",
    "\n",
    "a = sp.Symbol(r'a')\n",
    "l = sp.Symbol(r'l')\n",
    "L = sp.Symbol(r'L')\n",
    "theta = sp.Symbol(r'\\theta')\n",
    "e_s = sp.Symbol(r'\\epsilon_s')\n",
    "e_r = sp.Symbol(r'\\epsilon_r')\n",
    "\n",
    "x = sp.Symbol(r'x')\n",
    "e = sp.Function('e')(x)\n",
    "\n",
    "o = sp.Function('o')(e)\n",
    "o_e = sp.Symbol('\\sigma_e')\n",
    "K = sp.Symbol('K')\n",
    "n = sp.Symbol('n')\n",
    "#o_e = 400e6 # Pa sp.Symbol('\\sigma_e')\n",
    "#K = 900e6 # Pa sp.Symbol('K')\n",
    "#n = 0.5   # sp.Symbol('n')\n",
    "\n",
    "e_r = 0.47\n",
    "\n",
    "eq1 = sp.Eq(e, e_s*(1 - x/l)**theta)\n",
    "eq1\n",
    "\n",
    "eq2_RHS = L*( (e_s/e_r)**(1/theta) - 1 )\n",
    "eq2_LHS = a\n",
    "eq2 = sp.Eq(eq2_LHS, eq2_RHS)\n",
    "eq2\n",
    "\n",
    "eq3 = sp.Eq(e, e_r*(1 - x/L)**theta)\n",
    "eq3\n",
    "\n",
    "eq4 = sp.Eq(l, L*(e_s/e_r)**(1/theta))\n",
    "eq4\n",
    "\n",
    "eq6 = sp.Eq(o, o_e + K*e**n)\n",
    "#eq6.subs(e,e_r)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "id": "529f1f75-23b8-4393-a1bf-8ef42ddf31e9",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Eq(e(x), \\epsilon_s*(1 - x/l)**\\theta)\n",
      "Eq(e(x), \\epsilon_s*(1 - x/l)**5.0)\n",
      "\n",
      "Eq(a, L*((2.12765957446809*\\epsilon_s)**(1/\\theta) - 1))\n",
      "Eq(a, 0.000232600382236332*\\epsilon_s**0.2 - 0.0002)\n",
      "\n",
      "Eq(e(x), 0.47*(1 - x/L)**\\theta)\n",
      "Eq(e(x), 1.46875e+18*(0.0002 - x)**5.0)\n",
      "\n",
      "Eq(l, L*(2.12765957446809*\\epsilon_s)**(1/\\theta))\n",
      "Eq(l, 0.000232600382236332*\\epsilon_s**0.2)\n",
      "\n",
      "Eq(o(e(x)), 900000000.0*e(x)**0.5 + 400000000.0)\n",
      "Eq(o(e(x)), 900000000.0*e(x)**0.5 + 400000000.0)\n",
      "\n"
     ]
    }
   ],
   "source": [
    "# Substitute stuff\n",
    "for eq in [eq1, eq2, eq3, eq4, eq6]:\n",
    "    print(eq)\n",
    "    eq = eq.subs(o_e, 400e6) # Pa sp.Symbol('\\sigma_e')\n",
    "    eq = eq.subs(K, 900e6) # Pa sp.Symbol('K')\n",
    "    eq = eq.subs(n, 0.5)   # sp.Symbol('n')\n",
    "    eq = eq.subs(L, 200e-6) # m\n",
    "    eq = eq.subs(theta, 5.0)\n",
    "    eq = eq.subs(e_r, 0.47)\n",
    "    print(eq)\n",
    "    print()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "id": "d95a6b15-b6f4-4bea-8738-a06b54a03440",
   "metadata": {
    "tags": []
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle a = L \\left(\\left(2.12765957446809 \\epsilon_{s}\\right)^{\\frac{1}{\\theta}} - 1\\right)$"
      ],
      "text/plain": [
       "Eq(a, L*((2.12765957446809*\\epsilon_s)**(1/\\theta) - 1))"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "eq2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "723d48b3-da63-447f-a424-f8b6ee2d1f6d",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "40b92018-8698-47d7-862c-cfae9839d693",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "a3421df3-ee9f-4a04-bbe1-9b56f7722bc8",
   "metadata": {},
   "source": [
    "Micu2017894 - A New Model for the Equation Describing the Cavitation Mean Depth Erosion Rate Curve"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4b53e29d-0bf8-4199-a519-1bbdecc618f4",
   "metadata": {},
   "outputs": [],
   "source": [
    "def Micu2017894(A,B,C,t):\n",
    "    import numpy as np\n",
    "    cum_loss = A*t*(1-np.exp(-B*t))+(-C/2)*t*t    \n",
    "    ero_rate = A*(1-np.exp(-B*t))+A*B*t*np.exp(-B*t)-C*t\n",
    "    \n",
    "    return [cum_loss, ero_loss]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f4396d19-4365-4350-ae69-a949c26957a7",
   "metadata": {},
   "outputs": [],
   "source": [
    "Bordeacsu2006new - New Contributions to Cavitation\n",
    "Erosion Curves Modeling"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "882c5044-3447-4556-9a3a-13adfe26c403",
   "metadata": {},
   "outputs": [],
   "source": [
    "def Bordeacsu2006new(a, b, c, t):\n",
    "    \n",
    "    import numpy as np\n",
    "    \n",
    "    cum_loss = a*t - b*t*np.exp(-c*t)\n",
    "    ero_loss = (a - b*np.exp(-b*t) + b*c*t*np.exp(-b*t))\n",
    "\n",
    "    return [cum_loss, ero_loss]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "92ca0901-09fe-420f-964a-bc5ec766376e",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.11.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}