{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "d217c6ec",
   "metadata": {},
   "source": [
    "# Tutorial: simple uniaxial test\n",
    "\n",
    "## Material creation\n",
    "We start by importing the library and creating a neo-Hookean material model"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "11249daa",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pymecht as pmt\n",
    "mat = pmt.MatModel('nh')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "21e19c1a",
   "metadata": {},
   "source": [
    "We can access the parameters of the material model and print them. It shows that the parameters are of type `ParamDict` which includes their current value, whether they are fixed or not, and their lower/upper bounds (if applicable)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "90334cd9",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Material model with 1 component:\n",
      "Component1: NH\n",
      " ------------------------------------------------------------------\n",
      "Keys              Value       Fixed?      Lower bound Upper bound \n",
      "------------------------------------------------------------------\n",
      "mu_0              1.00        No          1.00e-04    1.00e+02    \n",
      "------------------------------------------------------------------\n",
      "\n"
     ]
    }
   ],
   "source": [
    "mat_params = mat.parameters\n",
    "print(mat,mat_params)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a86c62b1",
   "metadata": {},
   "source": [
    "## Sample  creation\n",
    "Next we can create a uniaxial extension sample from the material and print its parameters, which will include the material parameters as well as the geometric parameters (length and cross-sectional area)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "e6abb551",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "An object of type UniaxialExtensionwith stretch as input, force as output, and the following material\n",
      "Material model with 1 component:\n",
      "Component1: NH\n",
      "\n",
      "------------------------------------------------------------------\n",
      "Keys              Value       Fixed?      Lower bound Upper bound \n",
      "------------------------------------------------------------------\n",
      "L0                1.00        No          1.00e-04    1.00e+03    \n",
      "A0                1.00        No          1.00e-04    1.00e+03    \n",
      "mu_0              1.00        No          1.00e-04    1.00e+02    \n",
      "------------------------------------------------------------------\n",
      "\n"
     ]
    }
   ],
   "source": [
    "sample = pmt.UniaxialExtension(mat)\n",
    "sample_params = sample.parameters\n",
    "print(sample)\n",
    "print(sample_params)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7fc1e435",
   "metadata": {},
   "source": [
    "## Simulation by setting parameters\n",
    "\n",
    "We can set the parameters values to what we would like and then apply a stretch to the sample to calculate the force. We import `numpy` and `matplotlib` libraries for the stretch vector and plotting the force-strech plot."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "e39f3239",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAAsTAAALEwEAmpwYAAApJUlEQVR4nO3deXQV9f3/8eebsEVAFgGRhJ2wirJcwX1XsH4r1Gql1qXVFttKtbXSitpWsa0LrbULVflZqrVV1IoYV0RRq1YrQZaQsAVkSSJr2BOyvn9/3Ate4w0Ekpu75PU4J4c785m59z1nkvti5jMzH3N3REREqmsS6wJERCQ+KSBERCQiBYSIiESkgBARkYgUECIiElHTWBdQXzp27Og9e/aMdRkiIgllwYIFW929U6S2pAmInj17kpWVFesyREQSipmtq6lNp5hERCQiBYSIiESkgBARkYgUECIiEpECQkREIkqaq5hERBqb2QsLmDpnBYU7SujaLpVJo/szblhavb2/AkJEJAHNXljA5FnZlJRXAlCwo4TJs7IB6i0konqKyczGmNkKM8szs9sitH/fzLLNbJGZvW9mg0Lze5pZSWj+IjN7JJp1iogkmqlzVhwIh/1KyiuZOmdFvX1G1I4gzCwFmAZcAOQD880s091zwxZ7yt0fCS1/CfAgMCbUttrdh0arPhGRRFawoyTi/MIa5h+JaB5BjATy3H2Nu5cBM4Gx4Qu4+66wyVaARi8SETmIpQU7ufpv/6uxvWu71Hr7rGj2QaQBG8Km84FR1RcysxuBW4DmwLlhTb3MbCGwC7jT3d+LYq0iInFt3ba9/P6NlWQuLqTdUc0Ye2JX5uRuZF951YFlUpulMGl0/3r7zJh3Urv7NGCamV0J3AlcC3wGdHf3bWY2AphtZoOrHXFgZhOACQDdu3dv4MpFRKJvy+5S/jxvFU/9bz1NU4wbz+nDhDP70Da1WUJfxVQAdAubTg/Nq8lM4GEAdy8FSkOvF5jZaqAf8IWn8bn7dGA6QCAQ0OkpEUkau/eV8//+s4bH3v+U0ooqrjipGzefl8GxR7c8sMy4YWn1GgjVRTMg5gMZZtaLYDCMB64MX8DMMtx9VWjyYmBVaH4noMjdK82sN5ABrIlirSIicaG0opJ/frSeaW/nUbS3jIuHHMdPL+xH706tG7yWqAWEu1eY2URgDpACzHD3HDObAmS5eyYw0czOB8qB7QRPLwGcCUwxs3KgCvi+uxdFq1YRkVirrHJmLyzgwbkrKdhRwml9j+HnYwZwQnq7mNVk7slxZiYQCLjGgxCRROPuvL1iMw+8voLlG3dzfNrR/HzMAM7IiDiGT70zswXuHojUFvNOahGRxmrBuu3c/9pyPl5bRM9jjuLP3xzGxUOOo0kTi3VpgAJCRKTBrdq0mwfmrGBu7iY6tm7BPeOOZ/xJ3WiWEl/PT1VAiIg0kMIdJfxh7kqe/ySfVs2bcuuF/bju9F4c1Tw+v4rjsyoRkSSyfW8Zf30njyc+XAcO153Wix+e05cOrZrHurSDUkCIiERJSVklMz74lEfeXc2e0gouHZbOTy7IIL39UbEurVYUECIi9ay8sopnszbwxzdXsXl3KecP7Myk0QPo36VNrEs7LAoIEZF64u68mr2R372xgk+37iXQoz3TvjWck3p2iHVpR0QBISJymCI9A6lTmxbc//pyluTvpN+xrXnsmgDnDeyMWXxcsnokFBAiIoch0khutzy7iCqHtHap/O7yE/nasDRS4uRehrpQQIiIHIZII7lVObRNbcpbPz2Lls1SYlRZ/YuvuzJEROJcTSO27SqpSKpwAB1BiIjUSllFFf/4cC0YEce+rM+R3OKFAkJE5CDcndeWbuT+15ezblsx/Y9tzdptxZRWRG8kt3ihgBARqcEn67fzm1eWsWDddvof24YnrhvJWf06RX0kt3ihgBARqWZDUTH3v76cl5d8Rqc2Lbjv0iFcNiKdpqGH6UV7JLd4oYAQEQnZWVLOtLfzePyDtTRpAjedl8ENZ/amVYvG+VXZOLdaRCRMWUUV//rfOv741ip2lpRz2fB0fnphf7q0bXnolZOYAkJEGi13Z07OJu5/fTmfbt3LaX2P4favDGRw17axLi0uKCBEpFFavGEHv3llGR+vLaJv59b8/dsncXb/Tgn9aIz6FtWAMLMxwB+BFOAxd7+vWvv3gRuBSmAPMMHdc0Ntk4HrQ203ufucaNYqIo1D/vZips5ZwYuLCunYujm/+drxXBHodqADWj4XtYAwsxRgGnABkA/MN7PM/QEQ8pS7PxJa/hLgQWCMmQ0CxgODga7Am2bWz92/eH+7iEgt7doX7ID++wdrMWDiOX254azetGnZLNalxa1oHkGMBPLcfQ2Amc0ExgIHAsLdd4Ut34rP708cC8x091LgUzPLC73fh1GsV0SSUHllFU9/vJ6H3lxF0d4yLh2exq0X9k/KO5/rWzQDIg3YEDadD4yqvpCZ3QjcAjQHzg1b96Nq637pomMzmwBMAOjevXu9FC0iycHdeXPZZu59bRlrtuzl5N4duPPiQRyfpg7o2op5J7W7TwOmmdmVwJ3AtYex7nRgOkAgEIjwdBQRaYyy83fym1dz+WhNEb07tUqKsRliIZoBUQB0C5tOD82ryUzg4SNcV0SEwh0lTJ2zghcWFtChVXPuGTuY8SO700wd0EckmgExH8gws14Ev9zHA1eGL2BmGe6+KjR5MbD/dSbwlJk9SLCTOgP4OIq1ikgC272vnIffWc3f3v8UB35wdh9+cHYfjlYHdJ1ELSDcvcLMJgJzCF7mOsPdc8xsCpDl7pnARDM7HygHthM6vRRa7lmCHdoVwI26gklEqj8k76cXZLC3vIqH5q5k294yxg3tyq2j+5Pe/qhYl5oUzD05Tt0HAgHPysqKdRkiEiXVh/qEz4dmGNmzA3dcPJATu7WLVXkJy8wWuHsgUlvMO6lFRGoj0lCfDnRo1ZxnbjhZHdBRoJ4bEUkINQ31uX1vmcIhSnQEISJxrbyyin99tK7Gdt3wFj0KCBGJW++v2srdL+WwavMe+nVuzbqixjHUZ7xQQIhI3Fm/rZh7Xsllbu4munVI5dGrR3DhoGN5cVFhoxjqM14oIEQkbuwtrWDa23k89t6nNE0xJo3uz/Wn96JlsxSg8Qz1GS8UECISc1VVzuxFBdz32nI27y7l0mFp/GzMgEY/olusKSBEJKYWb9jBXS/lsHD9Dk5Mb8vDV41gRI/2sS5LUECISIxs3r2PB15fwb8X5NOxdQseuOwELhueTpMmumQ1XiggRKRBlVZU8vgHa/nzvDxKKyq54czeTDy3rwbuiUMKCBFpEO7OvOWbueflXNZuK+a8AZ254+KB9O7UOtalSQ0UECISdXmb93DPy7m8u3ILvTu14vHvnMTZ/TvHuiw5BAWEiETNzpJy/vTWKp7471pSm6Vw58UDufbUnhqfIUEoIESk3lVWOc9mbeB3c1ZQVFzGFYFu3Dq6Px1bt4h1aXIYFBAiUq/mry3irswccgp3EejRnicuGalxoBOUAkJE6kXhjhLufW05Ly0u5Li2Lfnj+KFccmJXPWk1gSkgRKRO9pVX8ui7a3j43Tzc4aZz+/L9s/twVHN9vSQ67UEROSLuzmtLN/KbV5ZRsKOErwzpwuSLBtKtg4b7TBYKCBGplfDxoDu1aUGblk1ZvWUvA7q04anvjeLUPh1jXaLUs6hea2ZmY8xshZnlmdltEdpvMbNcM1tiZm+ZWY+wtkozWxT6yYxmnSJycPvHgy7YUYIDm3eXsnrLXi4bkc7LPzpd4ZCkonYEYWYpwDTgAiAfmG9mme6eG7bYQiDg7sVm9gPgAeCKUFuJuw+NVn0iUnsPzFn+pfGgAT5cvY2muqchaUVzz44E8tx9jbuXATOBseELuPvb7l4cmvwISI9iPSJyBLLzd1K4Y1/EtprGiZbkEM2ASAM2hE3nh+bV5HrgtbDplmaWZWYfmdm4SCuY2YTQMllbtmypc8Ei8rmdxeX8YvZSLpn2PjU9YFXjQSe3uOikNrOrgABwVtjsHu5eYGa9gXlmlu3uq8PXc/fpwHSAQCDgDVawSBJzd2Z9UsBvX13G9uIyrj2lJ/27tGbKS8u+cJpJ40Env2gGRAHQLWw6PTTvC8zsfOAO4Cx3L90/390LQv+uMbN3gGHA6urri0j9WbFxN7+YvZSP1xYxrHs7nrju87ugU5s11XjQjUw0A2I+kGFmvQgGw3jgyvAFzGwY8Cgwxt03h81vDxS7e6mZdQROI9iBLSJRsKe0gj++uZIZH6ylTcum3HfpEL4R6PaFwXs0HnTjE7WAcPcKM5sIzAFSgBnunmNmU4Asd88EpgKtgedCt+Ovd/dLgIHAo2ZWRbCf5L5qVz+JSD1wd17N3siUl3PYtKuU8Sd142djBtChVfNYlyZxwNyT49R9IBDwrKysWJchkjDWbNnDrzJzeG/VVgYddzT3jDteY0E3Qma2wN0DkdriopNaRBrOvvJKpr2dx6PvrqFF0ybc9dVBXHVyD93PIF+igBBpRN5atom7XsphQ1EJ44Z25faLB9K5TctYlyVxSgEh0gjkby/m7pdymZu7ib6dW+vZSVIrCgiRJFZaUclj733Kn+etwjBuu2gA153Wi+ZNdTpJDk0BIZKkPsjbyi9eXMqaLXsZM7gLv/jqINJ057McBgWESJLZtGsfv35lGS8tLqR7h6P4+3dO4pz+nWNdliQgBYRIkqiorOKJD9fxh7krKaus4ubzMvjB2X1o2Swl1qVJglJAiCSBrLVF3Dl7Kcs37uasfp24+5LB9OzYKtZlSYJTQIgksG17SrnvteU8tyCf49q25JGrhjN6cBdCTyYQqRMFhEgCqqpynp6/ngdeX8He0gpuOKs3N52bQasW+pOW+qPfJpE4Fz4WdNd2qYw/qRtvLt/M4g07GNWrA78edzwZx7aJdZmShBQQInFs/1jQ+8dhKNhRwu/nrqR1ixQeumIoY4d21ekkiRoFhEgcmzpnRcSxoNu0bKZHb0vU6XZKkThW05jPG3dGHiNapD7pCEIkDlVUVvHY+59S08P4NRa0NAQFhEicyc7fyW2zlpBTuIvjux5N3uY97KuoOtCusaCloSggROJEcVkFf5i7kr+9/ynHtG7Bw98azpjju/DiokKNBS0xoYAQiQP/WbmF21/IJn97Cd8c2Z3bLhpA29RmgMaClthRQIjEUNHeMu55OZcXFhbQu1MrnplwMqN6HxPrskSAKF/FZGZjzGyFmeWZ2W0R2m8xs1wzW2Jmb5lZj7C2a81sVejn2mjWKdLQ3J1Zn+Rz3u/f4eUlhdx0bl9evekMhYPElagdQZhZCjANuADIB+abWaa754YtthAIuHuxmf0AeAC4wsw6AL8CAoADC0Lrbo9WvSINZUNRMbe/kM17q7YyvHs77vv6CfTTndASh6J5imkkkOfuawDMbCYwFjgQEO7+dtjyHwFXhV6PBua6e1Fo3bnAGODpKNYrElUVlVXM+OBTHpy7khQzpowdzFWjetCkie6ElvgUzYBIAzaETecDow6y/PXAawdZ90u9dGY2AZgA0L1797rUKhJVSwuCl64uLdjF+QM7M2Xs8bqXQeJeXHRSm9lVBE8nnXU467n7dGA6QCAQqOmeIpGYKSmr5A9vBi9dbX9Uc6ZdOZyvDNHjuCUx1DogQh3IGe7+ppmlAk3dffdBVikAuoVNp4fmVX/f84E7gLPcvTRs3bOrrftObWsViQfvrQpeurqhqITxJ3Vj8kUDaXtUs1iXJVJrtQoIM/sewVM5HYA+BL+wHwHOO8hq84EMM+tF8At/PHBltfcdBjwKjHH3zWFNc4Dfmln70PSFwOTa1CoSa0V7y/j1K7nM+qSA3h1bMXPCyZysq5MkAdX2COJGgp3O/wNw91VmdtBR0N29wswmEvyyTwFmuHuOmU0Bstw9E5gKtAaeCx1yr3f3S9y9yMzuIRgyAFP2d1iLxCt358VFhUx5OZddJeVMPKcvE8/tqzGhJWHVNiBK3b1s/3lTM2sKNT5H7AB3fxV4tdq8X4a9Pv8g684AZtSyPpGY2lBUzJ2zl/Luyi0M7daO+74+hAFdjo51WSJ1UtuAeNfMbgdSzewC4IfAS9ErSyQxVFRW8fh/1/L7N1bSxOCurw7i6lN6kqJLVyUJ1DYgbiN4GWo2cAPBo4LHolWUSCLIKdzJbc9nk12wk/MGdOaecbp0VZJLbQMilWAfwv+DA3dJpwLF0SpMJF6VlFXy0Fsreey94KWrf7lyGBcPOU6XrkrSqW1AvAWcD+wJTacCbwCnRqMokXj1Qd5Wbn8hm3Xbirki0I3bv6JLVyV51TYgWrr7/nDA3feY2VFRqkkkLsxeWHBgHIYubVuS3i6V+eu206tjK57+3smc0keXrkpyq21A7DWz4e7+CYCZjQAiD5YrkgRmLyxg8qxsSsorAfhs5z4+27mPCwZ25s9XDtelq9Io1DYgbiZ4r0IhYEAX4IqoVSUSY1PnrDgQDuFyP9utcJBG45ABEeqQPgMYAOwfCHeFu5dHszCRWHF3CnZEPkAurGG+SDI65IBB7l4JfNPdy919aehH4SBJadOufXz3iawa23UZqzQmtT3F9IGZ/QV4Bti7f+b+PgmRROfuzF5UwF2ZuZRWVDJuaFfm5GykpLzqwDKpzVKYNLr/Qd5FJLnUNiCGhv6dEjbPgXPrtRqRGNi8ex93vLCUubmbGNGjPVMvO4HenVp/4Sqmru1SmTS6P+OGfWlYEpGkVauAcPdzol2ISENzdzIXF/KrzBxKyiq58+KBfOe0XgcekzFuWJoCQRq12j7uuy3BMaLPDM16l+ATVndGqzCRaNqyu5Q7Z2czJ2cTw7q343eXn0ifTq1jXZZIXKntKaYZwFLgG6Hpq4G/A5dGoyiRaHp5SSG/mL2UvWWVTL5oAN89o7cericSQW0Doo+7fz1s+m4zWxSFekSiZtueUn75Yg6vZH/Gielt+d3lJ5JxbJtYlyUSt2obECVmdrq7vw9gZqehO6klgbyW/Rl3zl7K7n0V/GxMfyac0ZumKYe8ylukUattQHwf+EeoLwJgO3BtdEoSqT9Fe8v45YtLeXnJZwxJCx419O+iowaR2jhoQJhZd3df7+6LgRPN7GgAd9/VINWJ1MHrSzdy5+xsdpaUc+uF/bjhrD4001GDSK0d6ghiNjAcwMyer9YPIRKXtu8t466XcnhxUSGDux7Nk9ePYuBxGv5T5HAd6r9T4Zd29D7cNzezMWa2wszyzOy2CO1nmtknZlZhZpdVa6s0s0Whn8zD/WxpnObmbuLCh/7DK0s+4yfn92P2jacpHESO0KGOILyG14cUesjfNOACIB+Yb2aZ7p4btth64NvArRHeosTdhx7OZ0rjtbO4nLtfymHWwgIGHnc0j3/nJAZ3bXvoFUWkRocKiBPNbBfBI4nU0GtC0+7uB/uv2Uggz93XAJjZTGAscCAg3H1tqK0q0huI1MZbyzYxeVY2RXvLuPm8DG48py/Nm6qvQaSuDhoQ7l6XB9+nARvCpvOBUYexfkszywIqgPvcfXb1BcxsAjABoHv37kdeqSSknSXlTHkpl+c/yWdAlzbM+PZJHJ+mowaR+lLby1xjoYe7F5hZb2CemWW7++rwBdx9OjAdIBAIHNYpMElsb6/YzOTns9myp5QfnduXH52boaMGkXoWzYAoALqFTaeH5tWKuxeE/l1jZu8Aw4DVB11Jkt6ufeX8+uVcns3KJ6Nza6ZfM4IT0tvFuiyRpBTNgJgPZJhZL4LBMB64sjYrmll7oNjdS82sI3Aa8EDUKpWE8J+VW/j580vYtGsfPzy7Dzefn0GLphr+UyRaohYQ7l5hZhOBOUAKMMPdc8xsCpDl7plmdhLwAtAe+KqZ3e3ug4GBwKOhzusmBPsgcmv4KElyu/eV89tXl/H0xxvo06kVs354GkO7tYt1WSJJz9yT49R9IBDwrKyah4qUxBE+UE+HVs2prHJ27Svne2f25ifn96NlMx01iNQXM1vg7oFIbfHcSS2N0OyFBUyelU1JeSUA2/aWYcDN52fw4/P7xbY4kUZGl31IXJk6Z8WBcNjPgeey8mNTkEgjpoCQuFFRWUXBjshPkS+sYb6IRI8CQuLC+m3FfOPRD2ts79outQGrERFQQEiMuTvPZW3goj/+h1Wb93D1KT1IrdYJndoshUmj+8eoQpHGS53UEjM7isu444WlvJL9GaN6deDBK4aS1i6VEd3bH7iKqWu7VCaN7s+4YWmxLlek0VFASEz8N28rtzy7mK17SvnZmP7ccGYfUpoEny4/bliaAkEkDiggpEGVVlTy4Bsrmf7eGnp1bMUL15zGkHQ9YE8kHikgpMHkbd7NTU8vIvezXXxrVHfuuHggRzXXr6BIvNJfp0Sdu/PPj9bx61eW0apFUx67JsD5g46NdVkicggKCImqLbtL+fnzS5i3fDNn9evE1MtPoHOblrEuS0RqQQEhUTNv+SYmPbeE3aUV3H3JYK45pQdmdugVRSQuKCCk3pWUVfLbV5fx5EfrGNClDU9POJl+x7aJdVkicpgUEFKvlhbs5OaZC1m9ZS/fO6MXt47urzEbRBKUAkLqRVWVM/29Nfz+jRV0aNWcf14/itMzOsa6LBGpAwWE1FnhjhJ++uxiPlyzjYuO78JvvzaE9q2ax7osEakjBYTUyctLCrl9VjYVVc4DXz+BywPp6ogWSRIKCDkiu/eVc1dmLs9/ks/Qbu146Iqh9OzYKtZliUg9UkDIYVuwrogfP7OIgu0l3HReBj86ty/NUvRgYJFkE9W/ajMbY2YrzCzPzG6L0H6mmX1iZhVmdlm1tmvNbFXo59po1im1U1FZxR/mruTyRz7EHZ694RRuuaCfwkEkSUXtCMLMUoBpwAVAPjDfzDLdPTdssfXAt4Fbq63bAfgVECA44uSC0Lrbo1WvHNy6bXv58TOLWLh+B5cOT+PuSwbTpmWzWJclIlEUzVNMI4E8d18DYGYzgbHAgYBw97Whtqpq644G5rp7Uah9LjAGeDqK9UoE7s6/F+RzV2YOKU2Mv1w5jP87oWusyxKRBhDNgEgDNoRN5wOj6rDulwYIMLMJwASA7t27H1mVUqMdxWXc/kI2r2Zv5OTeHXjwG0M19KdII5LQndTuPh2YDhAIBDzG5SS82QsLDozkdkzr5pRVVFFSXsltFw3ge2f0PjCgj4g0DtHsXSwAuoVNp4fmRXtdOQKzFxYweVY2BTtKcGDrnjJ276vgpvMy+P5ZfRQOIo1QNANiPpBhZr3MrDkwHsis5bpzgAvNrL2ZtQcuDM2TKJk6ZwUl5ZVfmOfAzI83RF5BRJJe1ALC3SuAiQS/2JcBz7p7jplNMbNLAMzsJDPLBy4HHjWznNC6RcA9BENmPjBlf4e1REfBjpKI8wtrmC8iyS+qfRDu/irwarV5vwx7PZ/g6aNI684AZkSzPgmOEX3vq8trbFentEjjpTucGrH124q57OEPefy/azmrX0daNvvir0NqsxQmje4fo+pEJNYS+iomOXKvZX/Gz/69BDN49OoRjB7c5QtXMXVtl8qk0f0ZN+xLVxeLSCOhgGhk9p9Sevy/azkxvS1/uXI43TocBcC4YWkKBBE5QAHRiKzfVszEpz9hSf5OrjutF7ddNIDmTXWWUUQiU0A0Eq8v/YxJ/16C8fkpJRGRg1FAJLmDnVISETkYBUQS21BUzMSnPmGxTimJyBFQQCSp15duZNK/FwM6pSQiR0YBkWTKKqq497Vl/P0DnVISkbpRQCSR8FNK3zmtJ5MvGqhTSiJyxBQQSSL8lNIjV41gzPE6pSQidaOASHDhp5ROSG/LNJ1SEpF6ooBIYNVPKd120QBaNE2JdVkikiQUEAlqTs5GJj23GAceuWo4Y44/LtYliUiSUUAkmLKKKu57bTkzPviUE9Lb8pdvDqf7MTqlJCL1TwGRQMJPKX371J5M/opOKYlI9CggEoROKYlIQ1NAxLmyiiruf305f3v/U4akBa9S0iklEWkICog4tqGomIlPL2Txhh06pSQiDS6qt9ma2RgzW2FmeWZ2W4T2Fmb2TKj9f2bWMzS/p5mVmNmi0M8j0awzHr2Rs5GL//Qeazbv4eFvDeeuSwYrHESkQUXtCMLMUoBpwAVAPjDfzDLdPTdsseuB7e7e18zGA/cDV4TaVrv70GjVF0/Ch/o8rm1LMo5tzbsrtzIkrS1/uXIYPY5pFesSRaQRiuYpppFAnruvATCzmcBYIDwgxgJ3hV7/G/iLmVkUa4o7sxcWMHlWNiXllQAU7txH4c59nJHRkceuDeioQURiJpqnmNKADWHT+aF5EZdx9wpgJ3BMqK2XmS00s3fN7IxIH2BmE8wsy8yytmzZUr/VN5Cpc1YcCIdwa7bsVTiISEzF66M+PwO6u/sw4BbgKTM7uvpC7j7d3QPuHujUqVODF1kfCneUHNZ8EZGGEs2AKAC6hU2nh+ZFXMbMmgJtgW3uXuru2wDcfQGwGugXxVpjYk9pBS2aRd4FXdulNnA1IiJfFM2AmA9kmFkvM2sOjAcyqy2TCVwben0ZMM/d3cw6hTq5MbPeQAawJoq1NrhPt+7la9M+YF95Fc1SvtjtktoshUmj+8eoMhGRoKh1Urt7hZlNBOYAKcAMd88xsylAlrtnAn8DnjSzPKCIYIgAnAlMMbNyoAr4vrsXRavWhjZv+SZunrmIpk2Mp747is27Sw9cxdS1XSqTRvdn3LDq3TUiIg3L3D3WNdSLQCDgWVlZsS7joKqqnGlv5/HgmysZdNzRPHr1CNLb665oEYkdM1vg7oFIbbqTuoHsKa3glmcW8UbuJsYN7cq9l55AanNdpSQi8UsB0QDWbNnDhCcX8OnWvfzi/wZx3Wk9aWS3e4hIAlJARNlbyzbx45mLaNa0CU9eP5JT+3SMdUkiIrWigIiSqirnz/Py+MObKzk+7WgeuUr9DSKSWBQQUbB7Xzm3PLuYubmb+NqwNO69dAgtm6m/QUQSiwKinq3esocJ/8hi7bZifvXVQXz7VPU3iEhiUkDUo7m5m7jlmWB/wz+vH8UpfY459EoiInFKAVEPqqqcP81bxUNvrmJIWlseuXoEaXpUhogkOAVEHe3aV84tzyzmzWWbuHR4Gr/9mvobRCQ5KCDqIG/zHiY8mcW6bcXc9dVBXKv+BhFJIgqII/RGzkZueXYxLZo24V/fHcXJvdXfICLJRQFxmKqqnIfeWsWf3lrFCelteeSqEXo0t4gkJQXEYdi1r5yfzFzEW8s3c9mIdH497nj1N4hI0lJA1FLe5t1M+McC1hcVc/clg7nmlB7qbxCRpKaAqIU5ORv56bOLadks2N8wSv0NItIIKCAOoqrKeejNlfxpXh4nprflYfU3iEgjooCowc6Scn7yzCLmLd/M5SPSuUf9DSLSyCggIli1aTcTnlzAhqJi7hk7mKtOVn+DiDQ+CohqXl+6kZ8+u4jU5ik89b2TGdmrQ6xLEhGJiagGhJmNAf4IpACPuft91dpbAP8ARgDbgCvcfW2obTJwPVAJ3OTuc6JR4+yFBUyds4LCHSW0btGU3aUVnNitHY9cNZzj2qq/QUQarybRemMzSwGmARcBg4BvmtmgaotdD2x3977AH4D7Q+sOAsYDg4ExwF9D71evZi8sYPKsbAp2lODA7tIKUpoYV43qrnAQkUYvagEBjATy3H2Nu5cBM4Gx1ZYZCzwRev1v4DwLnuwfC8x091J3/xTIC71fvZo6ZwUl5ZVfmFdZ5Tz05qr6/igRkYQTzYBIAzaETeeH5kVcxt0rgJ3AMbVcFzObYGZZZpa1ZcuWwy6wcEfJYc0XEWlMohkQUefu09094O6BTp06Hfb6Nd3ToHsdRESiGxAFQLew6fTQvIjLmFlToC3BzurarFtnk0b3J7XavQ2pzVKYNLp/fX+UiEjCiWZAzAcyzKyXmTUn2OmcWW2ZTODa0OvLgHnu7qH5482shZn1AjKAj+u7wHHD0rj30iGktUvFgLR2qdx76RDGDfvS2SwRkUYnape5unuFmU0E5hC8zHWGu+eY2RQgy90zgb8BT5pZHlBEMEQILfcskAtUADe6e2XED6qjccPSFAgiIhFY8D/siS8QCHhWVlasyxARSShmtsDdA5HaErqTWkREokcBISIiESkgREQkIgWEiIhElDSd1Ga2BVhXh7foCGytp3JiKVm2A7Qt8SpZtiVZtgPqti093D3incZJExB1ZWZZNfXkJ5Jk2Q7QtsSrZNmWZNkOiN626BSTiIhEpIAQEZGIFBCfmx7rAupJsmwHaFviVbJsS7JsB0RpW9QHISIiEekIQkREIlJAiIhIREkfEGY2w8w2m9nSGtrNzP5kZnlmtsTMhoe1XWtmq0I/10Zav6HUcTsqzWxR6Kf6I9cbXC22ZYCZfWhmpWZ2a7W2MWa2IrSdtzVMxTWr47asNbPs0H6J+ZMma7Et3wr9bmWb2X/N7MSwtrjZL3XcjkTbJ2ND27IoNLrm6WFtdf/+cvek/gHOBIYDS2to/wrwGmDAycD/QvM7AGtC/7YPvW6faNsRatsT6/1wmNvSGTgJ+A1wa9j8FGA10BtoDiwGBiXitoTa1gIdY70/DmNbTt3/NwBcFPa3Elf75Ui3I0H3SWs+70s+AVgeel0v319JfwTh7v8hONZETcYC//Cgj4B2ZnYcMBqY6+5F7r4dmAuMiX7FkdVhO+LOobbF3Te7+3ygvFrTSCDP3de4exkwk+B2x0wdtiXu1GJb/hv6WwD4iOBIjxBn+6UO2xF3arEtezyUCEArYP/revn+SvqAqIU0YEPYdH5oXk3z49XB6m0ZOvz8yMzGNXhl9SfR9smhOPCGmS0wswmxLuYwXU/wiBUSe7+Ebwck4D4xs6+Z2XLgFeC60Ox62SdRG1FO4koPdy8ws97APDPLdvfVsS5KOD20XzoDc81seeh/jHHNzM4h+MV6+qGWjWc1bEfC7RN3fwF4wczOBO4Bzq+v99YRBBQA3cKm00Pzapofr2qs1933/7sGeAcY1tDF1ZNE2ycHFbZfNgMvEDxVE9fM7ATgMWCsu28LzU64/VLDdiTkPtkvFGS9zawj9bRPFBCQCVwTugroZGCnu39GcCztC82svZm1By4MzYtXEbcjVH8LgNAvzmkEx/pORPOBDDPrZWbNCY5hHvOrso6EmbUyszb7XxP8/Yp4pUq8MLPuwCzgandfGdaUUPulpu1I0H3S18ws9Ho40ALYRj19fyX9KSYzexo4G+hoZvnAr4BmAO7+CPAqwSuA8oBi4DuhtiIzu4fgLz/AFHc/WCdxVB3pdgADgUfNrIrgfwjuc/eYBsShtsXMugBZwNFAlZn9mOBVMbvMbCLBX/QUYIa758RgEw440m0h+HjmF0J/202Bp9z99QbfgDC1+B37JXAM8NdQ3RXuHnD3injaL0e6HcCxJN4++TrB/xiWAyXAFaFO63r5/tKjNkREJCKdYhIRkYgUECIiEpECQkREIlJAiIhIRAoIERGJSAEhUgtmdoeZ5YQ9OXOUmf3YzI46gve6vRbLPG5mlx1ZtSL1QwEhcghmdgrwf8Bwdz+B4KMMNgA/BiIGhJmlHOQtDxkQIvFAASFyaMcBW929FMDdtwKXAV2Bt83sbQAz22NmvzezxcApZnaVmX0cOuJ41MxSzOw+IDU071+h9a4JHZksNrMnwz73TAuOV7BGRxMSC7pRTuQQzKw18D7Bo4U3gWfc/V0zWwsEQoGBmTnBO1mfNbOBwAPApe5ebmZ/BT5y93+Y2R53bx1aZzDBZ/6c6u5bzaxD6C7+xwk+vvkKYACQ6e59G3TDpdFL+kdtiNSVu+8xsxHAGcA5wDMWedS0SuD50OvzgBHA/NCjG1KBzRHWORd4bn/IVHscwmx3rwJyzezYetkYkcOggBCpBXevJPgk3HfMLBuINITjvtByEBzZ7wl3n1yHjy0Ne211eB+RI6I+CJFDMLP+ZpYRNmsosA7YDbSpYbW3gMtC4wpgZh3MrEeordzMmoVezwMuN7Nj9i9X3/WLHCkdQYgcWmvgz2bWDqgg+MTcCcA3gdfNrNDdzwlfwd1zzexOgqOTNSE45OiNBINlOrDEzD5x92+Z2W+Ad82sElgIfLuBtkvkoNRJLSIiEekUk4iIRKSAEBGRiBQQIiISkQJCREQiUkCIiEhECggREYlIASEiIhH9f0bb0Q+3GMkOAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "sample_params.set('mu_0',5)\n",
    "sample_params.set('A0',0.1)\n",
    "sample_params.set('L0',10)\n",
    "applied_stretch = np.linspace(1,1.3,10)\n",
    "force = sample.disp_controlled(applied_stretch,sample_params)\n",
    "plt.plot(applied_stretch,force,'-o')\n",
    "plt.xlabel('Stretch')\n",
    "plt.ylabel('Force')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e8ce4660",
   "metadata": {},
   "source": [
    "## Parameter fitting\n",
    "\n",
    "If we have some force measurements (say from experiments), we can also do parameter fitting. We should fix the values of `A0` and `L0` to ensure a unique solution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "0c3e5bc3",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Parameter fitting instance created with the following settings\n",
      "------------------------------------------------------------------\n",
      "Keys              Value       Fixed?      Lower bound Upper bound \n",
      "------------------------------------------------------------------\n",
      "L0                10.00       Yes         -           -           \n",
      "A0                0.10        Yes         -           -           \n",
      "mu_0              5.00        No          1.00e-04    1.00e+02    \n",
      "------------------------------------------------------------------\n",
      "\n",
      "1 parameters will be fitted.\n",
      "   Iteration     Total nfev        Cost      Cost reduction    Step norm     Optimality   \n",
      "       0              1         4.1049e-01                                    1.20e+01    \n",
      "       1              2         2.2983e-02      3.88e-01       5.00e+00       2.62e+00    \n",
      "       2              3         1.0599e-03      2.19e-02       1.48e+00       4.24e-02    \n",
      "       3              4         1.0540e-03      5.94e-06       2.48e-02       1.19e-05    \n",
      "       4              5         1.0540e-03      4.65e-13       6.93e-06       2.01e-10    \n",
      "`gtol` termination condition is satisfied.\n",
      "Function evaluations 5, initial cost 4.1049e-01, final cost 1.0540e-03, first-order optimality 2.01e-10.\n",
      "Fitting completed, with the following results\n",
      "------------------------------------------------------------------\n",
      "Keys              Value       Fixed?      Lower bound Upper bound \n",
      "------------------------------------------------------------------\n",
      "L0                10.00       Yes         -           -           \n",
      "A0                0.10        Yes         -           -           \n",
      "mu_0              11.51       No          1.00e-04    1.00e+02    \n",
      "------------------------------------------------------------------\n",
      "\n",
      "Results after fitting\n",
      "------------------------------------------------------------------\n",
      "Keys              Value       Fixed?      Lower bound Upper bound \n",
      "------------------------------------------------------------------\n",
      "L0                10.00       Yes         -           -           \n",
      "A0                0.10        Yes         -           -           \n",
      "mu_0              11.51       No          1.00e-04    1.00e+02    \n",
      "------------------------------------------------------------------\n",
      "\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAAsTAAALEwEAmpwYAAAqmUlEQVR4nO3dd3xUVf7/8deHAIYm0sQSEAxFQUKLlCgdpa6IqCii6PpdfuhiWV0VBbHs2gurgousIquIoiCKggpKUQwgoRiaIkGBgEgVDDUh5/fHHdgACQQykzuTeT8fjzyYW2bmc5gwb+49955jzjlERCR6FfO7ABER8ZeCQEQkyikIRESinIJARCTKKQhERKJccb8LOFmVK1d2NWrU8LsMEZGIsnDhwq3OuSq5bYu4IKhRowYpKSl+lyEiElHMbG1e23RqSEQkyikIRESinIJARCTKRVwfQW4yMzNJT09n3759fpcS8WJjY4mLi6NEiRJ+lyIihaRIBEF6ejrlypWjRo0amJnf5UQs5xzbtm0jPT2dmjVr+l2OiBSSInFqaN++fVSqVEkhUEBmRqVKlXRkJRJGRs5OIzlt6xHrktO2MnJ2WtDeo0gEAaAQCBL9PYqEl4S48gwct/hwGCSnbWXguMUkxJUP2nsUiVNDIiJFVVJ8ZYb3aczAcYvp27w6Y+evY3ifxiTFVw7aexSZI4JwkJ6eTo8ePahduzbx8fHcddddHDhwgDFjxjBw4EC/yztG2bJl/S5BRPIhKb4yfZtX5+UZq+nbvHpQQwCiMAhCdb7NOcdVV13FlVdeyU8//cSqVavIyMhg8ODBBXrdvGRlZYXkdUUk/CSnbWXs/HXc2b4WY+evO+Y7rKCiLghCdb5txowZxMbGcssttwAQExPDsGHDGD16NHv27GH9+vW0bduW2rVr89hjjwGwe/duunXrRsOGDbnooosYP348AAsXLqRNmzY0bdqUTp068euvvwLQtm1b7r77bhITE3niiSc477zzyM7OPvxa1apVIzMzk7S0NDp37kzTpk1p1aoVP/zwAwA///wzLVu2pEGDBgwZMqRA7RWRwnHoO2p4n8bcc3ndw6eJghkGUddHEKrzbcuXL6dp06ZHrDv99NOpXr06WVlZfPfddyxbtozSpUtz8cUX061bN9auXcs555zDlClTANi5cyeZmZnccccdfPzxx1SpUoXx48czePBgRo8eDcCBAwcOj7W0aNEiZs+eTbt27fj000/p1KkTJUqUoH///owcOZLatWszf/58br/9dmbMmMFdd93Fbbfdxk033cSIESMK1F4RKRyp6TuP+I469B2Wmr4zaKeIoi4I4MjzbXe2rxX08225ueyyy6hUqRIAV111FXPmzKFr167ce++9PPDAA3Tv3p1WrVqxbNkyli1bxmWXXQbAwYMHOfvssw+/Tu/evY94PH78eNq1a8d7773H7bffTkZGBsnJyVxzzTWH99u/fz8A3377LRMnTgTgxhtv5IEHHgh5u0WkYAa0iT9mXVJ85aB+b0VlEBx9vq1FfKUC/6XWq1ePCRMmHLFu165drFu3juLFix9zWaaZUadOHRYtWsTUqVMZMmQIHTp0oGfPntSvX5+5c+fm+j5lypQ5/PiKK67goYceYvv27SxcuJD27duze/duzjjjDJYsWZLr83V5qEgE2vELzH4WWv4VqtYP+suHtI/AzDqb2Y9mttrMBuWyvbqZzTSzxWaWamZdQ1kPhO58W4cOHdizZw9vvfUW4P1P/t577+Xmm2+mdOnSTJ8+ne3bt7N3714++ugjLrnkEjZu3Ejp0qXp27cv9913H4sWLaJu3bps2bLlcBBkZmayfPnyXN+zbNmyXHzxxdx11110796dmJgYTj/9dGrWrMkHH3wAeJ3Y33//PQCXXHIJ7733HgDvvPNOgdorIoVg16/w6T3wSiIsnQC/fh+StwlZEJhZDDAC6ALUA643s3pH7TYEeN851xi4Dng1VPUccrzzbQVhZkyaNIkPPviA2rVrU6dOHWJjY3nyyScBaNasGb169SIhIYFevXqRmJjI0qVLadasGY0aNeKxxx5jyJAhlCxZkgkTJvDAAw/QsGFDGjVqRHJycp7v27t3b8aOHXvEKaN33nmHN954g4YNG1K/fn0+/vhjAF566SVGjBhBgwYN2LBhQ4HaKyIhtHsbfDEYXm4Ei/4LTW6Eu5ZAoz4heTtzzoXmhc1aAo865zoFlh8EcM49lWOf14A1zrlnAvu/4JxLOt7rJiYmuqMnplm5ciUXXnhhsJsQtfT3KeKTfTsheTjMexUy90BCb2jzAFQs+NhfZrbQOZeY27ZQ9hGcC6zPsZwOND9qn0eBaWZ2B1AG6BjCekREwtOB3TD/Nfj2Jdj3O9TrAe0GQ5W6hfL2fncWXw+Mcc69EDgieNvMLnLOZefcycz6A/0Bqlev7kOZIiIhkLUfUt6Eb16A3Zuh9uXQfgic3bBQywhlEGwAquVYjgusy+lWoDOAc26umcUClYHNOXdyzo0CRoF3aihUBYuIFIqDWbDkHe9KoF3pUKMV9H4bqrfwpZxQBsECoLaZ1cQLgOuAo3s61gEdgDFmdiEQC2wJYU0iIv7JzoZlE2HWk7B9DZzbFHoMh/Pbgo+XdocsCJxzWWY2EPgCiAFGO+eWm9njQIpzbjJwL/AfM/sb4ICbXah6r0VE/OIc/DAFZj4Bm1dA1YvgunehbhdfA+CQkPYROOemAlOPWjc0x+MVwCWhrEFExDfOQdoMmPFP2LgIKtWCXm9A/augWPgM9RY+lUS4mJgYGjVqdPjnl19+ISnJuxL2l19+Ydy4cYf3XbJkCVOnTs3rpfLUtm1bjr50VkTC1Nq5MKYbjL0Kdm+BK4bD7fOhwdVhFQLg/1VDRUapUqWOGdbh0I1gh4KgTx+vi2TJkiWkpKTQtWvIb6QWkcK2cbF3BLD6SyhbFbo8B037QfHT/K4sTwqCECpbtiwZGRkMGjSIlStX0qhRI66//npGjBjB3r17mTNnDg8++CDdu3fnjjvuYNmyZWRmZvLoo4/So0cP9u7dyy233ML333/PBRdcwN69e/1ukojkZfNKrw9g5SdQqgJ0fAya9YeSpf2u7ISKXhB8Ngg2LQ3ua57VALo8fdxd9u7dS6NGjQCoWbMmkyZNOrzt6aef5vnnn+fTTz8FoGrVqqSkpDB8+HAAHnroIdq3b8/o0aP5/fffadasGR07duS1116jdOnSrFy5ktTUVJo0aRLcdolInkbOTiMhrvwRA1Imp20lNX3nkSOCbl8Ds56G1PehZFloMwha3g6xwZtTONSKXhD4JLdTQ/k1bdo0Jk+ezPPPPw/Avn37WLduHV9//TV33nknAAkJCSQkJASrXBE5gUOTWB0amyzngJUA7NwAXz8Li8dCsRKQdAdc+jcoXdHfwk9B0QuCE/zPPRw555g4cSJ16xbO7eQicmJ5TmJV1cHnD8KCN8BlQ9NboPXfodxZfpd8ysKr67qIKleuHH/88Ueey506deKVV17h0C0UixcvBqB169aHrzZatmwZqamphVi1iOScxOrWpmeQ9Mur8FJDmD8SGlwDdyyEbs9HdAiAgqBQJCQkEBMTQ8OGDRk2bBjt2rVjxYoVNGrUiPHjx/Pwww+TmZlJQkIC9evX5+GHHwbgtttuIyMjgwsvvJChQ4ceMxWmiIRWctpWJs77kbF1vubG73p4YwLV6QR//Q6uHAEVzvO7xKAI2TDUoaJhqENPf58iMHfVRma/+zz3nPYRJfdtY3tcB27b2IW7buhZKNPbBptfw1CLiESe7GxYMYl6U4bS0qVD1Uuh46NUrHYxdwWuGorEIDgeBYGIyCFrZsOXj8DGxZQ/sz5c9S+o1fHweEDBnjQ+XBSZIHDOaWL2IIi0U4UiQbFpKUx/BNK+gvLV4MqRkHAtFIvxu7JCUSSCIDY2lm3btlGpUiWFQQE459i2bRuxsbF+lyJSOHas9e4GTn3fuwHs8n/CxX+BEtH1b6BIBEFcXBzp6els2aKpDAoqNjaWuLg4v8sQCa3d2+Cb52HB62DF4NK74ZK7odQZPhfmjyIRBCVKlKBmzYJP7iwiRdyBPd7E8N++BAcyoNEN0PZBKH+u35X5qkgEgYjIcR3MgiVjYeZTkLEJ6naDDkPhzAv8riwsKAhEpOhyDn74FL58DLb9BNWawzVj4LyWflcWVhQEIlI0rZ0L04dC+ndQuQ5cNw7qdg2LqSHDjYJARIqWzSu9I4BVn0G5s+FPL3t9ATH6usuL/mZEpGjYuQFmPQlLxnnzAnQYCs1vi4iJYfymIBCRyLZ3B8wZBvNf84aFbnE7tLo3IucF8IuCQEQiU+Y++G6UNyLovp2Q0BvaD4YzqvtdWcRREIhIZMk+CKnjYcYTsCsdal0GHR/xppSVU6IgEJHI4Bz8NA2+fBQ2r4BzGkPPf0PN1n5XFvEUBCIS/tJTvEHh1s6BCjXh6jehfk9dChokCgIRCV9bf4KvHoeVk6FMFej6PDS9GWJK+F1ZkaIgEJHwk7EFZj0FC8dAiVLeeEAtB8JpZf2urEhSEIhI+Mjcx7x3nyBx/WiKZ+2FxD9Dm/tJ/q0YqfN+Y0AbBUEoaPJ6EfGfc7B8EoxoRos1L5OcWZfFf/oMuj1P8m/FGDhuMQlx5f2ussjSEYGI+Ct9IXzxIKyfD1Uvghs/ojgNuHXcYvpu/ZGx89cxvE/jIjlFZLhQEIiIP35fD189Bks/gDJnwhWveGMCFYshCejbvDovz1jNne1rKQRCTEEgIoVr/x8w518wd7i33Orv3gxhp5U7vEty2lbGzl/Hne1rMXb+OlrEV1IYhJCCQEQKR/ZBWDwWZvwTdm+GBtd6A8OdUe2I3ZLTtjJw3OLDp4NaxFc6YlmCT0EgIqGXNhOmDYHflnmTw1z/LsQl5rpravrOI770k+IrM7xPY1LTdyoIQkRBICKhs2UVTH8YVn3uDQZ3zRiod+Vx7wge0Cb+mHVJ8ZUVAiGkIBCR4Nuz3bshbMEbULIMdHwMmg+AErF+Vya5UBCISPBk7Yfv/gNfP+t1Cje9xbsruGwVvyuT41AQiEjBOQcrP/HmCN7xszc09OX/gDMv9LsyyQcFgYgUzMbF8MVgWPstVLkQ+k6EWh39rkpOQkiHmDCzzmb2o5mtNrNBeexzrZmtMLPlZjYulPWISBDt3ACTBsCotrDlR+g+DAbMUQhEoJAdEZhZDDACuAxIBxaY2WTn3Ioc+9QGHgQucc7tMLMzQ1WPiATJgd3w7Uvw7cveHMGX/g0uvQdiT/e7MjlFoTw11AxY7ZxbA2Bm7wE9gBU59vkLMMI5twPAObc5hPWISEFkZ8P373rzA2RsgvpXQcdHocJ5flcmBRTKIDgXWJ9jOR1oftQ+dQDM7FsgBnjUOff50S9kZv2B/gDVq2tiapFC9/M38MVDsCkVzk2E3m9DtWZ+VyVB4ndncXGgNtAWiAO+NrMGzrnfc+7knBsFjAJITEx0hVyjSFQYOTuNhLjyR9y4tWjxAip8+wQ1t86E8tWg1xtwUS9NEVnEhLKzeAOQcxCRuMC6nNKByc65TOfcz8AqvGAQkUKWEFeegeMWk5y2FfZsZ+N7d9Pg405U//07b0yggQugwdUKgSIolEcEC4DaZlYTLwCuA/octc9HwPXAm2ZWGe9U0ZoQ1iQieUiKr8zw6xrwzdgnaRLzAVWzMthS+1rO6vEPKKvrOIqykAWBcy7LzAYCX+Cd/x/tnFtuZo8DKc65yYFtl5vZCuAgcJ9zbluoahKR4/j5G5KmP0CSW86c/fVZ02QwN/Xs5ndVUghC2kfgnJsKTD1q3dAcjx1wT+BHRPywM90bGXT5JPaVOZfBxe7j3JZXM/a79dRK2KrB3qKA353FIuKXzH0w9xX45kVw2axLuIvey5rxwg0tvHkAalXWPABRQkEgEm2c84aF/nwQ7PgFLvwTXP4EU1OzeOGG8poHIAopCESiydbVXgCsng6V68KNH0F8OwAGtDl2d80DEB0UBCLRYP8f8PVzMPdVKFEKOj0JzfpDTAm/K5MwoCAQKcqcg6UfeMND//ErNLoBOjwC5ar6XZmEEQWBSFH16/cw9X5YPw/OaQzXvg3VLva7KglDCgKRombPdpjxD1g4BkpVhCtegUZ9oVhIR52XCKYgECkqsg/Cwjdhxj9h3y6vD6Dtg1DqDL8rkzCnIBApCtbOhc/ug01LoUYr6PIMVK3vd1USIRQEIpFs10aY/ggsfR9OPxeufhPq99TAcHJSFAQikShrP8x7FWY/B9lZ0Po+b6awkmX8rkwikIJAJNKsmubdFLY9Dep29e4JqFjT76okgikIRCLFtjRvlrBVn0OlWnDDRKitieKl4BQEIuHuwG745gVIfgViSsJlj0Pz26B4Sb8rkyJCQSASrpyD5R/CtIdh1wZI6A0dH4PTz/a7MiliFAQi4WjTMvjsAVg7B85qAFePhuot/K5KiigFgUg42bsDZj4JC16H2PLQfRg06QfFYvyuTIowBYGIz0bOTiPh3HIk/THNGxxu7w5+rd2Hz8+8lVsSm/hdnkQBBYGIz5qX/Q039gZgJVRrwZKEwfz58/0Mb17d79IkSigIRPxyYA98/SyNk18hs2RZHsu8jXJx/Rj7ebqmh5RCpSAQ8cOqL2Dq3+H3ddDoBkpc9jjlvt3GyzNWc2f7WgoBKVQal1akMO3cAOP7wrhroXgpuHkKXPkqyZtg7Px13Nm+FmPnryM5bavflUoU0RGBSGE4mAXfjYKZT3hjA3UYCi3vgOIlSU7bysBxiw+fDmoRX+mIZZFQUxCIhFr6Qvj0btiUCrU6QtfnjxgbKDV95xFf+knxlRnepzGp6TsVBFIoFAQiobL3d2+msAVvQLmz4Jr/Qr0exwwRPaBN/DFPTYqvrBCQQpPvIDCz84DazrkvzawUUNw590foShOJUM7Bsonw+YOwZys0HwDtHoLY0/2uTCRX+QoCM/sL0B+oCMQDccBIoEPoShOJQNvSYMq9sGamN2H8De97f4qEsfweEfwVaAbMB3DO/WRmZ4asKpFIk7Uf5vzLGyW0+GleP0DinzU0hESE/AbBfufcAQuc2zSz4oALWVUikWTNLO8oYNtquKiXN1FMubP8rkok3/IbBLPN7CGglJldBtwOfBK6skQiQMZm+GKwN19whZrQd6J3VZBIhMlvEAwCbgWWAv8PmAq8HqqiRMJadjYsGgNfPuoNE9H6fmh1D5Qo5XdlIqckv0FQChjtnPsPgJnFBNbtCVVhImFp01L49G+QvgBqtIJuL0KVOn5XJVIg+Q2Cr4COQEZguRQwDUgKRVEiYWd/Bsx6Cub9G0pVgJ6jIOHaY+4JEIlE+Q2CWOfcoRDAOZdhZqVDVJNI+HAOfpgCn93vTRfZ9Gbo8AiUruh3ZSJBk98g2G1mTZxziwDMrCmwN3RliYSB39fB1Pth1WdwZn24ZgxUa+Z3VSJBl98guAv4wMw2AgacBfQOWVUifjqYCXNHwOxnvOXL/+ndHRxTwt+6RELkhEEQ6BhuBVwA1A2s/tE5lxnKwkR8sW6e1xm8eQVc0B06Pw1nVPO7KpGQOmEQOOcOmtn1zrlhwLJCqEmk8O3ZDl8+AovegtPj4Lp34YKuflclUijyOzHNt2Y23MxamVmTQz8nepKZdTazH81stZkNOs5+vczMmVlivisXCYKRs1az6svRMDwRFr8DSXcyr+tURv5W98RPFiki8ttH0Cjw5+M51jmgfV5PCJxSGgFcBqQDC8xssnNuxVH7lcPrg5ifz1pEgmPnBq796T4qbpjBH5UaUu6mySTvPuvwpDAi0SJfQeCca3cKr90MWO2cWwNgZu8BPYAVR+33D+AZ4L5TeA+Rk5edDYv+C9OHUvFgJj83Hcw1ixPok1qCsfM1M5hEn3ydGjKz8mb2opmlBH5eMLPyJ3jaucD6HMvpgXU5X7cJUM05N+UE79//0Htv2bIlPyWL5G77GnjrCm/GsHMawe1zqfmn++nToiYvz1hN3+bVFQISdfLbRzAa+AO4NvCzC3izIG9sZsWAF4F7T7Svc26Ucy7ROZdYpUqVgrytRKvsg5A8HF5Ngl+/hz+9DDdNhoo1SU7bqonjJarlt48g3jnXK8fyY2a25ATP2QDkvO4uLrDukHLARcCswPDWZwGTzewK51xKPusSObHfVsDkgbBhIdTpAt1fhNPPAdDE8SLk/4hgr5ldemjBzC7hxHcWLwBqm1lNMysJXAdMPrTRObfTOVfZOVfDOVcDmAcoBCR4sg7ArKfhtdaw4xfo9QZc/+7hEIDjTxwvEi3ye0QwAHgrR7/ADqDf8Z7gnMsys4HAF0AM3uily83scSDFOTf5eM8XKZANC+Hjgd6NYQ2ugc7PQJlKx+ymieNFwJzLe6IxM6vunFuXY/l0AOfcrkKoLVeJiYkuJUUHDZKHA3tg5hMw71UoexZ0HwZ1O/tdlYjvzGyhcy7Xe7VOdETwEdAk8CITj+onEAkvP38Dn9zpXRnU9Ba47DGIPdHFbSJyoiDIOdj6+aEsROSU7dsJ0x+BhW96U0b2+wRqtva7KpGIcaIgcHk8FgkPq76AT+6GjE3QciC0GwwlNVWGyMk4URA0NLNdeEcGpQKPCSw759zpIa1OJC+7t8Hng7yJ48+sB73HQlxTv6sSiUjHDQLnXExhFSKSL87BsonejGH7dkHbB+HSe6B4Sb8rE4lY+b18VMR/u36FKffAj1PhnCbQYwRUred3VSIRT0Eg4c85b56AaQ/DwQPejGEtbodiOmAVCQYFgYS37T97l4T+/DXUaAV/egkqHXsTmIicOgWBhKfsgzD/NZjxD7AY6P4vaNIPiuV3VBQRyS8FgYSfzT94g8SlL4Danby7g8ufe+LnicgpURBI+Mg6AN/+C75+Dk4r5w0Sd1EvMDvhU0Xk1CkIJDxsWAST74Dflnlf/l2ehTIa+E2kMCgIxF+Ze2HmkzB3OJStCte9Cxd09bsqkaiiIBD/rF8Ak/4fbE/zOoIvexxKneF3VSJRR0EghS/rAMx+GuYMg9Pj4KaP4fy2flclErUUBFK4flvuHQVsWgqN+0KnpyBWQ1aJ+ElBIIUj+6DXDzDjn94cAeoLEAkbujtHQm/7zzCmG0wfCnU6we3zGPlbXZLTth6xW3LaVkbOTvOpSJHopSCQ0HEOUt6Ef18Cv62AnqPg2rehTGUS4sozcNziw2GQnLaVgeMWkxCnGcVECptODUlo/LHJmzx+9XSo2QaufBXKxx3enBRfmeF9GjNw3GL6Nq/O2PnrGN6nsSaNF/GBgkCCb9mH3nDRmfugy3Nw8f/lOkZQUnxl+javzsszVnNn+1oKARGf6NSQBM+e7TDhVphwC1SMhwFzoHn/PAeKS07bytj567izfS3Gzl93TJ+BiBQOHRFIcKz+0jsVtHsLtB8Cl/wNYvL+9TrUJ3DodFCL+EpHLItI4dERgRTM/gz49G8wthfEngH/9xW0vu+4IQCQmr7ziC/9Q30Gqek7C6FoEcnJnHN+13BSEhMTXUpKit9lCMC6+d7NYTt+gaSB0G4IlIj1uyoRyYWZLXTOJea2TaeG5ORl7YdZT8G3L3lXAt08BWpc4ndVInKKFARycjYt844CflsGTW6CTk96cweISMRSEEj+ZB/0jgBmPgmlKkCf9727hEUk4ikI5MS2pcFHt8H6+VCvB3QbBmUq+V2ViASJgkDy5hykjIZpQyCmBFz1OjS4WlNHihQxCgLJ3a6N3n0BaV/B+e2gxwhNIC9SRCkI5FhLJ8CUe72rg7o+7w0RoaMAkSJLQSD/s2e7N0bQ8kkQdzH0fA0qxftdlYiEmIJAPKumweSBXhh0GApJd53w7mARKRr0Lz3a7c+AaYNh4Rg4sx7cMAHOTvC7KhEpRBprKMqMnJ32v1E+186FkZfgFv6XxdX7Qf9ZCgGRKKQgiDIJceW5+53vSP/gfnizC/sys7m12OPsbTMUip/md3ki4gOdGooySRV28VWFpyi3PJXUqldy25ZePHdDkoZ+FoliIT0iMLPOZvajma02s0G5bL/HzFaYWaqZfWVm54WynqiX+j6MbE253WuZXPcZrlh7Lb1a1FUIiES5kAWBmcUAI4AuQD3gejOrd9Rui4FE51wCMAF4NlT1RLX9GTDpNvjwL3DWRaR0/ZRHV8drZjARAUJ7RNAMWO2cW+OcOwC8B/TIuYNzbqZzbk9gcR4QhwTXxiXwWmtIfQ/aDCK51Rj6f/wbw/s05p7L6x6eQF5hIBK9QhkE5wLrcyynB9bl5Vbgs9w2mFl/M0sxs5QtW7YEscQizDmYOwJe7wiZe6HfJ9DuQVI37tbMYCJyhLDoLDazvkAi0Ca37c65UcAo8GYoK8TSIlPGFvj4dvhpGtTtBj2GQ+mKAAxoc+ydwknxldVPIBLFQhkEG4BqOZbjAuuOYGYdgcFAG+fc/hDWEx3SZnoTx+z9XeMEiUi+hDIIFgC1zawmXgBcB/TJuYOZNQZeAzo75zaHsJai72AmzHwC5vwLKteBvh/CWRf5XZWIRICQBYFzLsvMBgJfADHAaOfccjN7HEhxzk0GngPKAh+Y97/Wdc65K0JVU5G1/WeYeCtsWAhN+kHnp6Fkab+rEpEIEdI+AufcVGDqUeuG5njcMZTvHxWWToBP/wYYXDMG6vf0uyIRiTBh0Vksp+DAbph6PywZC9WaQ6/X4YzqflclIhFIQRCJfk2FCX+Gbauh9X3QZpCGjBaRU6Zvj0jiHMx/DaY/DKUrQb/JULO131WJSIRTEESK3Vvho9vhpy+gThdvDuEylfyuSkSKAAVBJFgzGz7sD3u3Q5dnoVl/3RsgIkGjIAhnBzNh1lPwzYtQqRb0nQBnNfC7KhEpYhQE4WrHWu/egPQF0PhG6PIMlCzjd1UiUgQpCMLRsg/hk7sBB1ePhot6+V2RiBRhCoJwcmA3fD4IFr0FcRd79wZUqOF3VSJSxCkIwsWmpd69AVt/glb3QtsHIaaE31WJSBRQEPjNOfjuPzBtCJSqADd9BOe39bsqEYkiCgI/7dkOH/8VfpwKtS+HK/8NZTQvgIgULgWBX37+xrs3YM9Wb7TQ5gN0b4CI+EJBUNgOZsHsp+Hr56FSPFw/Hc5p5HdVIhLFFASF6Y9N8MEtsC4ZGt3g3SV8Wlm/qxKRKKcgKCy/fAsTboH9f0DPUdCwt98ViYgACoLQcw7mDofpj0DFmnDjR1C1nt9ViYgcVszvAoq0fbvg/Ztg2hDWVG7D/I4TjgiB5LStjJyd5mOBIiIKgtDZvBL+0w5+mAKX/5NNnUZx24TVJKdtBbwQGDhuMQlx5X0uVESinU4NhcLSCTD5DihZ1ps8psalJAHD+zRm4LjF9G1enbHz1zG8T2OS4nXfgIj4S0EQTFkHYNpg+G4UVG/pTSZf7qzDm5PiK9O3eXVenrGaO9vXUgiISFjQqaFg2bkBxnT1QqDlQOj3yREhAN7poLHz13Fn+1qMnb/u8GkiERE/6YggGNbM8gaMy9rvHQXU73nMLof6BA6dDmoRX+mIZRERv+iIoCCys+GbF+DtnlCmCvxlZq4hAJCavvOIL/2k+MoM79OY1PSdhVmxiMgxzDnndw0nJTEx0aWkpPhdBuz9HSYNgFWfeRPH/Oll3SUsImHLzBY65xJz26ZTQ6fi11R4/0bYma7J5EUk4ikITtbid2DKPd7cATdPherN/a5IRKRAFAT5lbkPPrsfFv0XarSCq9+EslX8rkpEpMAUBPmxY603VMSvS+DSv0G7IRCjvzoRKRr0bXYiP30JH/4fZB+E68bBBd38rkhEJKgUBHnJzoavn4VZT0PV+nDtW95EMiIiRYyCIDd7tsOHf4HVX0LD66Hbi1CytN9ViYiEhILgaBsWwfv9IGMTdB8GTW/RpaEiUqQpCA5xDhaO8a4MKlsV/vw5nNvU76pEREJOQQBwYA9MuRe+Hwfx7eGq16FMJb+rEhEpFAqC7Wtg/E3w21Jo84D3UyzG76pERApNdAfBD1O98YLMoM8HUOdyvysSESl00RkEB7Ng5hMw50U4u5F3aWiF8/yuSkTEFyEdhtrMOpvZj2a22swG5bL9NDMbH9g+38xqBLuGkbPTjpwAJmMLv/+nuxcCTfrBn79QCIhIVAtZEJhZDDAC6ALUA643s3pH7XYrsMM5VwsYBjwT7DoS4sozcNxiLwzWL2D/q5cSu2khq1s+A1e8DCVig/2WIiIRJZRHBM2A1c65Nc65A8B7QI+j9ukB/DfweALQwSy4F+0fmgDm87HDODi6M5v3OH7sNoFanQYE821ERCJWKIPgXGB9juX0wLpc93HOZQE7gWOu2zSz/maWYmYpW7ZsOelCkuIrc2H9RkzPaszkZu/Q8OI2J/0aIiJFVURMVemcG+WcS3TOJVapcvJDPyenbeW5lWewotUI3lj4uyaNFxHJIZRBsAGolmM5LrAu133MrDhQHtgWzCJyThp/z+V1Gd6n8f/6DEREJKRBsACobWY1zawkcB0w+ah9JgP9Ao+vBma4IE+irEnjRUSOL2T3ETjnssxsIPAFEAOMds4tN7PHgRTn3GTgDeBtM1sNbMcLi6Aa0ObYoaOT4isfDgYRkWgX0hvKnHNTgalHrRua4/E+4JpQ1iAiIscXEZ3FIiISOgoCEZEopyAQEYlyCgIRkShnQb5aM+TMbAuw9hSfXhkoKjcQqC3hp6i0A9SWcFWQtpznnMv1jtyIC4KCMLMU51yi33UEg9oSfopKO0BtCVehaotODYmIRDkFgYhIlIu2IBjldwFBpLaEn6LSDlBbwlVI2hJVfQQiInKsaDsiEBGRoygIRESiXJEIAjMbbWabzWxZHtvNzF42s9VmlmpmTXJs62dmPwV++uX2/MJUwLYcNLMlgZ+jh/wudPloywVmNtfM9pvZ34/a1tnMfgy0c1DhVJy7ArbjFzNbGvhMUgqn4rzloy03BH6vlppZspk1zLEtbD6TQD0FaUukfS49Am1ZEpit8dIc2wr+Heaci/gfoDXQBFiWx/auwGeAAS2A+YH1FYE1gT8rBB5XiMS2BLZl+P1ZnGRbzgQuBp4A/p5jfQyQBpwPlAS+B+pFWjsC234BKvv9WZxEW5IO/RsAuuT4txJWn0lB2hKhn0tZ/tenmwD8EHgclO+wInFE4Jz7Gm8+g7z0AN5ynnnAGWZ2NtAJmO6c2+6c2wFMBzqHvuK8FaAtYedEbXHObXbOLQAyj9rUDFjtnFvjnDsAvIfXbl8UoB1hJx9tSQ78WwCYhzezIITZZwIFakvYyUdbMlzgmx8oAxx6HJTvsCIRBPlwLrA+x3J6YF1e68PZ8WqODRw2zjOzKwu9suCJxM8lLw6YZmYLzay/38WcpFvxjj4h8j+TnG2BCPxczKynmf0ATAH+HFgdlM8lpBPTSKE7zzm3wczOB2aY2VLnXJrfRUW5SwOfyZnAdDP7IfC/v7BmZu3wvjwvPdG+4S6PtkTc5+KcmwRMMrPWwD+AjsF67Wg5ItgAVMuxHBdYl9f6cJZnzc65Q3+uAWYBjQu7uCCJxM8lVzk+k83AJLxTLGHNzBKA14EezrltgdUR+Znk0ZaI/FwOCQTW+WZWmSB9LtESBJOBmwJX3LQAdjrnfsWbT/lyM6tgZhWAywPrwlmubQm04TSAwC/IJcAKPwstgAVAbTOraWYl8eay9v0qqJNlZmXMrNyhx3i/X7leFRIuzKw68CFwo3NuVY5NEfeZ5NWWCP1capmZBR43AU4DthGk77AicWrIzN4F2gKVzSwdeAQoAeCcG4k3b3JXYDWwB7glsG27mf0D75cc4HHn3PE6akPuVNsCXAi8ZmbZeAH/tHPO1yA4UVvM7CwgBTgdyDazu/GuRNllZgPxfqFjgNHOueU+NAE49XbgDRk8KfDvtzgwzjn3eaE3IId8/H4NBSoBrwbqznLOJTrnssLpM4FTbwtQlcj7XHrh/QcwE9gL9A50HgflO0xDTIiIRLloOTUkIiJ5UBCIiEQ5BYGISJRTEIiIRDkFgYhIlFMQiASY2WAzW55jlMfmZna3mZU+hdd6KB/7jDGzq0+tWpHgURCIAGbWEugONHHOJeDdvr8euBvINQjMLOY4L3nCIBAJFwoCEc/ZwFbn3H4A59xW4GrgHGCmmc0EMLMMM3vBzL4HWppZXzP7LnAE8ZqZxZjZ00CpwLp3As+7KXCk8b2ZvZ3jfVubN1b+Gh0diF90Q5kIYGZlgTl4//v/EhjvnJttZr8AiYFgwMwc3l2d75vZhcCzwFXOuUwzexWY55x7y8wynHNlA8+pjzeeTZJzbquZVQzc1T4Gb0jh3sAFwGTnXK1CbbgIRWSICZGCcs5lmFlToBXQDhhvuc/CdRCYGHjcAWgKLAgMV1AK2JzLc9oDHxwKk6OGAPjIOZcNrDCzqkFpjMhJUhCIBDjnDuKN2jrLzJYCuU37ty+wH3izxP3XOfdgAd52f47HVoDXETll6iMQAcysrpnVzrGqEbAW+AMol8fTvgKuDoxpj5lVNLPzAtsyzaxE4PEM4Bozq3Rov2DXL1IQOiIQ8ZQFXjGzM4AsvNFd+wPXA5+b2UbnXLucT3DOrTCzIXgzXRXDm6ryr3gBMgpINbNFzrkbzOwJYLaZHQQWAzcXUrtETkidxSIiUU6nhkREopyCQEQkyikIRESinIJARCTKKQhERKKcgkBEJMopCEREotz/B2BxLZS9OD+IAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "force_measured = np.array([0,0.09,0.21,0.29,0.4,0.49,0.58,0.65,0.75,0.84])\n",
    "def sim_func(params):\n",
    "    force_modeled = sample.disp_controlled(applied_stretch,params)\n",
    "    return force_modeled\n",
    "\n",
    "sample_params.fix('A0')\n",
    "sample_params.fix('L0')\n",
    "param_fitter = pmt.ParamFitter(sim_func,force_measured, sample_params)\n",
    "param_fitter.fit()\n",
    "print(\"Results after fitting\")\n",
    "print(sample_params)\n",
    "plt.plot(applied_stretch,force_measured,'x',label='Observed')\n",
    "plt.plot(applied_stretch,sim_func(sample_params),'-',label='Fitted')\n",
    "plt.xlabel('Stretch')\n",
    "plt.ylabel('Force')\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "78297270",
   "metadata": {},
   "source": [
    "## Random Parameters: Monte Carlo simulation\n",
    "We can also approach the problem probabilistically by treating the parameteres as random variables. Either we can generate samples `RandomParameters` and propagate them through the model (Monte Carlo simulation). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "d3215990",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Keys              Type              Lower/mean        Upper/std         \n",
      "------------------------------------------------------------------------\n",
      "L0                fixed             10.00             10.00             \n",
      "A0                fixed             0.10              0.10              \n",
      "mu_0              uniform           1.00e-04          1.00e+02          \n",
      "------------------------------------------------------------------------\n",
      "\n",
      "random_samples are a list of parameters dictionaries: <class 'list'> 300\n",
      "results are of shape (300, 10)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXgAAAEGCAYAAABvtY4XAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAAsTAAALEwEAmpwYAAAqi0lEQVR4nO3deXxkZZ3v8c+vKpU9vWbtLJ1eoEGWhiaytIiKgoIgqIygMorjTOtcHcfrzMtR79zrdud19S4zzow60ioqi4CitOyCgriAjI2y9JJ0ls6+VPbaUvtz/6hKk2660+mkTp1TVb/368WrK1VJzu9Q3d88+Z3nPI8YY1BKKZV/XHYXoJRSyhoa8Eoplac04JVSKk9pwCulVJ7SgFdKqTxVZHcBC1VXV5vW1la7y1BKqZzx/PPPTxhjao73mqMCvrW1lb1799pdhlJK5QwR6TvRa5a2aERkjYjcJyLtInJQRC6x8nhKKaVeYfUI/l+Bx4wxN4hIMVBu8fGUUkqlWRbwIrIauAy4BcAYEwWiVh1PKaXU0axs0WwCxoHvicifROQ7IlJx7CeJyC4R2Ssie8fHxy0sRymlCouVAV8E7AD+wxhzPhAEPnPsJxljdhtj2owxbTU1x70QrJRSahmsDPhBYNAY81z64/tIBb5SSqkssCzgjTGjwICIbEs/9WbggFXHU0opdTSrZ9H8DXBXegZND/Ahi4+nlFIqzdKAN8a8ALRZeQyllMpVxhiGZuaoriyh1OPO+Pd31J2sSilVKMZ8Ybq9AULRBJdsKbbkGBrwSimVRZOBCF3eAP5w3PJjacArpVQW+MIxurwBpgLZu99TA14ppSwUisbp9gYZ84WzfmwNeKWUskAknqBnPMjwzBzG2FODBrxSSmVQPJGkdzLEwFSIRNKmZE/TgFdKqQxIJg2D03McngwSiyftLgfQgFdKqRUxxjAyG6ZnPEg4lrC7nKNowCul1DJ5/WG6vUGCEeunPC6HBrxSSp2imVCULm+AmVDM7lIWpQGvlFJLFIjE6fIGmPBH7C5lSTTglVLqJMKxBN3jAUZnw7ZNeVwODXillDqBaDxJ72SQwekQSWdMjDklGvBKKXWMRNLQPxWidzJIIpFDQ/ZjaMArpVSaMem57BNBog6Zy74SGvBKKcXRy/fmCw14pVRBy+byvdmmAa+UKkj+cIzOLC/fm20a8EqpgpKrUx6XQwNeKVUQnLTKY7ZowCul8tr8xtY94/kxM+ZUaMArpfLWRCBC51jAsYuBWU0DXimVdwrhAupSaMArpfJGOJbaJm9k1r5t8pxEA14plfMSSUPvZJD+ycK5gLoUlga8iPQCfiABxI0xbVYeTylVWIwxDM+G6RkPEIkV1gXUpcjGCP5NxpiJLBxHKVVAJgMROr0BAnl4B2qmaItGKZVTApE4nWN+Jgv8AupSWB3wBnhcRAxwqzFm97GfICK7gF0ALS0tFpejlMpVkXiCbq9eQD0VVgf8pcaYIRGpBZ4QkXZjzK8XfkI69HcDtLW16dumlDpKImnomwzSNxXK6bXZ7WBpwBtjhtJ/ekXkfuBC4NeLf5VSSqUuoI7MhunWC6jLZlnAi0gF4DLG+NOPrwS+ZNXxlFL5YyoY5dCYXy+grpCVI/g64H4RmT/OD40xj1l4PKVUjtMLqJllWcAbY3qA7VZ9f6VU/ojEU3egDs/oBdRM0mmSSinb5Mvm1k6lAa+UssXI7Bzd3iDhWP7sgeo0GvBKqayaCUU5NBbANxezu5S8pwGvlMqKcCxBlze1VZ7KDg14pZSlkklD31SI3omgrvSYZRrwSinLeH1hOr0B5qLaZ7eDBrxSKuMCkTgdo36mgzqf3U4a8EqpjIklknSPBxia1vnsTqABr5RaMWMMg9NzdI8HiOt8dsfQgFdKrYiuG+NcGvBKqWWZiybo9Prx+iJ2l6JOQANeKXVKEknD4Ykg/VNBkrqKr6NpwCullmx0Nkyn16/rs+cIDXil1En5wjEOjfqZCenyArlEA14pdUK6D2pu04BXSr1KMmkYmA7RM6HL+OYyDXil1FEmAhEOjfkJRXR5gVynAa+UAiAUTS0voNvl5Q8NeKUKXDyR5PBEkIHpkE57zDMa8EoVsKGZObq9AaJxTfZ8pAGvVAGaDcXoGPPrrkp5TgNeqQKiuyoVFg14pQqAMYb+KZ32WGg04JXKczOhKO2jutpjIdKAVypPReNJOr1+Rma0HVOoLA94EXEDe4EhY8w1Vh9PqUKnm2+oedkYwf8tcBBYlYVjKVXQZudidIzq7BiV4rLym4tIE/B24DtWHkepQhdLJDk44mNv75SGuzrC6hH814BPA1Un+gQR2QXsAmhpabG4HKXyz/DMHJ3eADG9WUkdw7IRvIhcA3iNMc8v9nnGmN3GmDZjTFtNTY1V5SiVd/zhGHt7pzgw7NNwV8dl5Qj+dcA7RORqoBRYJSJ3GmNutvCYSuW9eCJJz0SQgamQrtGuFmXZCN4Y81ljTJMxphW4CXhSw12plRnzhXm2Z5L+SQ13dXI6D16pHBCKxmkf9TOlS/mqU5CVgDfG/Ar4VTaOpVQ+SSQNhyeC9E8FdSlfdcp0BK+UQ3n9YTrHAsxFdWcltTwa8Eo5zFw0QceYnwl/xO5SVI7TgFfKIZJJQ99UiN6JIImkXkEtKLEQlGT+Zn8NeKUcYDIQoWPUT0jbMQWjeG6MuoHHqBt4hDIC8PE/gEhGj6EBr5SNwrEEnWMBxny64mMhKA5PUDv4c+oGHmHN+F4Eg3/NGcS3v4fiZBzcnoweTwNeKRsYYxiYmqN7IqAbcOQ5T2SK2sHHqRt4hLXj/4mYJIFVW+k5628Ya76K0KotXLJlPcXuzMexBrxSWaYbcOS/ougstYNPpELd+ywukyBYtYnDZ/41Y81XEVx9enbqyMpRlFK6AUeec0f91Az/krqBR1g/9jtcyRihimb6tn2YsearCaw5M+M99pPRgFcqC4Zm5ugc8+sGHHnGHQtSPfIUdf0Ps370N7iTUebKN9B/2gcYa74a/9qzsx7qC2nAK2WhYCRO+6iP6aCu0Z4vXPE5qkd+Rd3AI1SP/Ap3IkK4rI6hLe9jtOVqfOu22xrqC2nAK2WBZNLQOxmkd1KXGMgHrkSE9SNPUzfwKNUjT1EUDxEprWZ4058x1nwVM9UXgFi6f9KyaMArlWGzoRgHRnwEI3oRNZdJIsr6sd9RN/AwNUO/pCgeJFqyltGWdzDWcjXT1a8Fl9vuMhelAa9UhsQTSbrGAwxNz+lSvjlKkjHWjT1L3cAj1Aw9gSfmJ1a8mrHmqxhrfjvTtRdhXLkTm7lTqVIO5vWH6Rj1E4lpPybXSDLO2vHnqOtPhXpxdIa4pxJv4xWMNV/FVO1OjLvY7jKXRQNeqRWIxBN0jPrx+nRhsJxikqyZeJ66/oepG3yM4sgU8aIKxjdczljz1UzWvz5nQ30hDXillkmnPuYYY6ic7aC+70HqBx6iNDRCwl3K+IY3pUP9DSSLSu2uMqM04JU6RaFonIMjOvUxV5QGBqjvf4j6/oeo9HWSlCIm6y+l85y/Z2LD5SQ8FXaXaBkNeKWWaH4538MTAZ366HCe8CR1A49S3/8gayb/BMB09QUc3PEFvM1vI1ayzt4Cs0QDXqklmA3FODjq0/VjHMwdC1Az9Avq+x9k3dgzuEwC/+ptdJ7z94y1vJ1wRaPdJWbdkgNeRDYCpxljfiEiZUCRMcZvXWlK2S+eSNI9HmRwOqRTHx1IElHWj/6G+v4HqRl+EncizFx5I33b/pLRjddmbVEvp1pSwIvIXwG7gHXAFqAJ+BbwZutKU8pe4/7UJhzhmG7C4SgmyZrxvdT3P0Dd4M/xRGeJlqxluPVdjG58B7Prz3fMUgF2W+oI/mPAhcBzAMaYThGptawqpWwUiSc4NKqbcDiKMVTNHKSu/0Hq+x+idG6MeFE5441vYbTlWqbqdmJcmd0sIx8sNeAjxpiopH8qikgRoL+wqryjUx+dpSzQn5oB0/cAFf6e1AyYhsvo3P4PjG+4nGRRud0lOtpSA/5pEfkcUCYiVwD/BXjQurKUyi6d+ugcxXPjR2bArJ56EYDpmtfSf/oteJveSqxkrc0V5o6lBvxngA8DLwMfAR4BvrPYF4hIKfBroCR9nPuMMZ9ffqlKZZ5OfXQGdyxA7dDj1Pc9xDrvM4hJ4l9zJp3nfprRlrcTKW+wu8SctNSALwNuM8Z8G0BE3OnnQot8TQS43BgTEBEP8FsRedQY8/sVVaxUhszOxTg4olMf7SKJKNWjT1Pf9yDVI0/hTkQIVTTRe8ZHGG25luDqrXaXmPOWGvC/BN4CBNIflwGPAztP9AXGGLPg8z3p/7SxqWynUx9tlJ4B09C3h9rBn+OJ+YmUrGdo03sY3XitozbLyAdLDfhSY8x8WJMelZ/06kZ6pP88sBX4hjHmueWVqVRmTAQitI/o1MdsK/cfpr73ZzT0/Yyy0BDxogq8jVcwuvFapmsvyakleHPJUv+vBkVkhzHmjwAicgEwd7IvMsYkgPNEZA1wv4icbYzZt/BzRGQXqTn2tLS0nErtSi1ZJJ6gcyzA6KxOfcyWosgMdQOP0tC3hzWTf8KIi8m619F1zqcYb3wLyaIyu0vMe0sN+L8Ffiwiw4AA9cCNSz2IMWZGRJ4C3gbsO+a13cBugLa2Nv2FWWXc8Mwch3TqY1bM31na0LeHmuEncSVjBFadRue5n2Zk47VEy+rsLrGgnDTg022W1wNnANvST3cYYxadTyYiNUAsHe5lwBXAV1dYr1JLFo4lODDiYyoQtbuU/GYMVdP7aOjdQ/3AQxRHpomUrGdwy/sZab0e/5ozta9uk5MGvDEmISLvNcb8C8eMvk+iAfhB+geEC/iRMeahZdap1CnRG5asVxIapb7vZzT07aHS103CVcx445sZ2Xg9U/WX6p2lDrDUFs3vROTrwL1AcP7J+Z788RhjXgLOX1l5Sp2acCzBwREfkzpqt4Q7FqRm6Aka+vawbuxZBMN09QUcuODLeJuvIl68yu4S1QJLDfjz0n9+acFzBrg8o9UotQLaa7dIMsHa8edo6N1D7dDjFMVDhCqa6Tnr44xuvI65Sp0c4VRLCnhjzJusLkSp5QrHErSP+pnw676omVTu66Khdw8NfQ9QOjdKzFPFaMs1jGy8ntnqC7SvngOWulzwauDzwGXpp54GvmSMmbWqMKWWYmR2jo5RHbVniicyRV3/wzT03s/q6X0kxc1k/es5dN5nmGi4PO/2LM13S23R3EbqAut70h//OfA94F1WFKXUyUTiCdpH/IzrqH3FJBGleuQpGnr3UD3yNC4Tx7fmNXSc9znGWq4hWlptd4lqmZYa8FuMMe9e8PEXReQFC+pR6qRGZ8N0jPmJxXV1sGUzhlVTL6anNj6MJzpLpLSW/tNvYWTjdQTXbDv591COt9SAnxORS40xvwUQkdexhDtZlcqkSDxBx6gfr09H7ctVGhyioe9n1PfuoSLQS8JdirfxCkZar2eqdie43HaXqDJoqQH/UeD2dC8eYBr4oDUlKfVqY74w7aM6al8OVyJCzdATbOj5Meu9zwIwVXMRvWd+BG/TW0l4Km2uUFll0YAXkRZjTL8x5kVgu4isAjDG+LJSnSp40XiS9lGfjtqXoWp6PxsO30d934N4Yj7mKproPusTjLS+k3BFo93lqSw42Qh+D7ADQER+ckwfXilLeX1hDuqo/ZR4ItPU9z/IhsM/oWrmIAlXMd6mtzK86Qamay8CcdldosqikwX8womum60sRKl50XiSjlG/bnq9VMkE67zPsuHwfdQOPYErGcO39mwO7vgCYy1vJ168+qTfQuWnkwW8OcFjpSzh9YdpH/ET1VH7SZUGBtjQ+1Maeu+nLDRMtHgNg1vey/CmdxNYc6bd5SkHOFnAbxcRH6mRfFn6MemPjTFGF55QGRFLpEbtul774lzxcOqC6eH7WO99FoMwWX8pndv/gfENb8a4i+0uUTnIogFvjNE5U8pyOmo/CWNeuWDa/xCemI9QRRPdZ/0tw5vepRtSqxPSfbKUbXTUvjhPZJr6vgfZ0HsfVTPtJNwleBvfyvDmG5iuuVAvmKqT0oBXthj3R2gf9RGJ6aj9KMkE67zPpC+Y/gJXMsbskQum1+hyvOqUaMCrrIolkhwa8zMyo6P2heYvmG7o/SmloZH0BdP3pS+YnmF3eSpHacCrrJkIRDg4oqP2ea54mNqhx9lw+D7WeX9/5ILpoe2f0QumKiM04JXl4okkh8YCDM/o8kVHXzB9EE/Mn7pgevYnGW59p14wVRmlAa8sNRmIcHDETziWsLsUW6UumD7AhsP3UTXbkbpgOn+HqV4wVRbRgFeWSCQNh8b8DE0X8KjdGNZM7KWx+x7qBh9bcMH0i+k7TPWCqbKWBrzKuNm5GPuHZglFC3PUXhSdpaF3D40991Dp6ybmqWJo840MbX6PXjBVWaUBrzLGGEPvZIie8QCm0Ba2MIbVky/Q2H03dYOP4k5EmF23nf2v/V+MNV9NsqjM7gpVAdKAVxkxF02wf3iWmVDM7lKyyh3109D3Mxp77qVqtoN4UQUjre9icPONBNa+xu7yVIHTgFcrNjI7R/uon0ShbHxtDKumXqKx517q+x/GnZjDt/ZsDlzwZcZariHhqbC7QqUADXi1ArFEkvaRwlnW1x0LUN//EE3dd1M1c5B4UTmjLdcwuOUm/OvOsbs8pV7FsoAXkWbgdqCO1FLDu40x/2rV8VR2TQej7B/2FcT0x6rp/TR230N9/0MUxYP415zBwR1fYHTjO3S7O+VoVo7g48DfGWP+KCJVwPMi8oQx5oCFx1QWSyYNPRMB+iZDeX0h1RUPUd//MI0997B66mUS7lLGmq9mcMtN+NZtB5GTfxOlbGZZwBtjRoCR9GO/iBwEGgEN+BwVjMTZNzSLPxy3uxTLVM6009h9Dw39D1AUCxBYdRod5/8jIxuv13nrKudkpQcvIq3A+cBzx3ltF7ALoKWlJRvlqGUYnA7RORYgkcy/YbsrHqZu8BEau+9lzeSfUvuYNr+Nwc03MVt9gY7WVc6yPOBFpBL4CfBJY4zv2NeNMbuB3QBtbW35lx45LhpPcmDEx4Q/YncpGVcx20Vjzz009O7BE/MRrNrEoe2fZaT1emIla+0uT6kVszTgRcRDKtzvMsb81MpjqcybCEQ4MOzLq52WXIkItYOP0dh9D2snnifp8uBtvJKhLTel14TR0brKH1bOohHgu8BBY8w/W3UclXnJpKHTG2BgKmR3KRlT7uuhsedeGnrvpzg6Q6hyI53nfprh1ncRK11nd3lKWcLKEfzrgD8HXhaRF9LPfc4Y84iFx1Qr5A/H2DfkIxjJ/QupkoxRM/QETV13s278OZJSxHjjWxjaciNTtZfoCo4q71k5i+a3gP6+m0P6JoN0jwdI5nhHpjg8QWP3vTT23E3pnJe5iia6zvkUw63vJlpWY3d5SmWN3smqCMcS7B/2MR2M2l3KiqyafJHmrjuoG3gUVzLGZN2ltF/wZSbqLwOX2+7ylMo6DfgC5/WFOTDiI56j68hIIkrdwCM0d93B6qmXiRdVMLj5Jga3vp/Qqs12l6eUrTTgC1Q8kaQjhze/LgmN0tR9N40991IcmSJYtZn28/8HI63X6/IBSqVpwBeg2VCM/cM5uCFHeoek5s47qBl6AjFJJjZczsBpNzNVu1OnOCp1DA34AmKM4fBEkMMTwZxaR8YVn6O+7wGau+6karaDWPFq+k//EINb3ku4stnu8pRyLA34ApGLG3KUBfpp6vohGw7fhyfmw796Gwfa/iejLdfqDklKLYEGfAHIqQ05TJJ1Y8/Q3HUn1cNPYcSFt+lKBrfezEx1m7ZhlDoFGvB5LJc25HDHAmzo/SlNXXdR4T9MpGQ9h8/8a4a23ESkvN7u8pTKSRrweWo6GGXf8CyRmLPvWir3ddPcdRcNvfdTFA8yu+5c9l30fxhrugrjLra7PKVymgZ8njHG0D3u8A05kgmqR5+mufMO1o/9jqTLw1jz1QxsvRnf+u12V6eUpUSgpMhNWbGb8vR/JUXWLJuhAZ9HwrEE+4aceyG1KDLDhsM/oan7LsqDg4TL6ug6+78ytPk9xErX212eUhkjAqWeVIiXedxUFBdRWuyioriIMo8blys715I04PPEuD/CgREfMQcu7Vs5005z5x3U9z+IOxFmuua1dJ37acYb34xxeewuT6llWRji5cVuyj1FRx5nM8QXowGf44wxdHlTLRknkWScmqEnaO68g7UTe0m4SxnZ+A4Gt76fwJoz7S5PqSXJhRBfjAZ8DnNiS8Yd9dN4+Ec0d95BWWiYUEUTh7b/A8Ot7yZessbu8pR6lROFeEWJm9Ii54f4YjTgc9REIML+Yee0ZEqDQzR33k5jz48oigeZrrmQjvP/OxMNb9SVHJUjlHhclBcXUV6c6onnS4gvRgM+x8zPkumdcEZLZtXUS7R03Ebt4M8BGGu+mv7TP4R/3dk2V6YKkafIdWRmynyYzz9252mIL0YDPoc4piWTTFAz8iQtHbexduJ54p5K+k+/hYHTPkCkvMHe2lTec7skNfpeMAov9xRRXuLG49ZduhbSgM8RTmjJuOIhNhz+KS2dP6A80MdceSMd532O4U036BK9KqNcLigtclNeUvSqEXmpR1t+S6UB73BOaMkUz3lp7rqTxu57KI7OMLvuXF665GuMN16JcelfIbV88xc350fh84/LPG5E1x1aMf3X6WB2t2QqZ9ppOfQ96vsfQpJxxhuvoG/bh5hdv0MX/VJL5nYL5R43FenReMWRUXlh9sWzSQPeoWxryRjD+tHf0HLoe6wf+x0JdxlDm2+k//RbmKtsyW4tKqeUetyUl6R64wuDXFsq9tGAdxi7WjKSiNLQ9wAth75Hpa+TSGktnef8HUObb9T56+oIt0uOCm8djTubBryD2NGS8USmaOq6m6buuygJT+BfvY39F36V0ea362qOBazE43olxBeMyHU0nls04B0i2y2Zcl8PLYe+T0Pf/bgTESbqL2Pftr9guvYS7a8XiIXTDY+0VtJ/6mg8P2jA2yyrLRljWDP+BzYeuo2a4SdJuIoZ3Xgd/ad/iODqrdYfX9nCU+SisiTVRqkoLqKiREfjhcKygBeR24BrAK8xRm9rPI5stWQkGaN24DE2Hvoeq6b3ES1ZS89rPsbg1vcTLa229Ngqe+bbKpXp1krqzyKKLVprXDmflSP47wNfB2638Bg5KxstGXfUT1PPvTR33k7p3CjBqk0cvOBLjGy8nmRRqWXHVdYRIbW+eElqJF5eXJR6XOymSO/iVMewLOCNMb8WkVarvn+uykZLpjQ4SMuh29lw+McUxYNM1V5M+wVfZKLhDSAaArnA5YIyT3o0XuI+MiqvKC7K24WxVObZ3oMXkV3ALoCWlvyeZ211S6ZipoPW9t3UDzyMwZVa+Gvbh/CvPcuS46mVc7vlyCyVhWGud3KqTLA94I0xu4HdAG1tbU7dRXTFJgMR9lnUklk1+QKbDn6LmuEniReV03/aLfSffguR8vqMH0stT5FbqCyZb6fohU6VHbYHfL5LtWSC9E4EM/2NWTf2DK0Hv8W68eeIFq+h+6xPMLD1Zr0xyUbzQV5enGqvVJamwrykSINcZZ8GvIUsacmYJDVDT7Dp4K2smt5HuKyWQ9s/y9Dm95DwVGTuOGpRbpccudBZmZ65oiNydSq+9XQ35zatZueWV2ayPdM9wUuDs3z0DVsycgwrp0neDbwRqBaRQeDzxpjvWnU8p8l0S0aSMer7H2Rj+7ep9HUTqmzhwAVfZqT1nXrHqYVcLtItlVdCvLIkteqhUitxbtNqPv7DP/H1953Pzi3VPNM9ceTjTLFyFs17rfreTpbplowrHmbD4R+zseO7lIWG8a/exssX/zPeprfpUr0Z5HJxpK2ycGSuFzuVVXZuqebr7zufj//wT9x8UQt3Ptd/JOwzRRMig8KxBPuHZ5kOrrwl4476aer+IS2Hvk9JZJKZ9Tto3/F5JhveqEsJrIAIlBW/0laZD/TyYg3yQpGN1shS7dxSzc0XtfBvT3bxicu3ZjTcQQM+Y2ZCUV4anCW6wpaMJzxJS+cPaOq6C0/Mz0T963n5zI8yU92mwX6K5oP8lfaKziNX2WmNLNUz3RPc+Vw/n7h8K3c+18/FW9brCN5pBqZCdHr9JFeQ7SXBYTZ2fJfGwz/GlYjgbXorvWfs0s2rl2B+5kplaSrIq0o8VJTonZ3q+LLRGlmKhT9Ydm6p5uIt64/6OBM04FcgmTS0j/oZnplb9vco9/XQ2v5t6vt+BsDoxuvoPeMvCa3K7q+KuWD+7s6q0ldaK1WlOnMllzilPWJ1a2QpXhqcPSrM53/wvDQ4qwFvt3AswUuDs/jmltdvr5reT+vBW6kd/DlJdwmDW95L37YPE6nYkOFKc9P8wllVC0bm2l7JfU5pj1jdGlmK4/1A27mlWls0dlt2v90Y1kzspfXgt6ge/Q1xTyW9Z36E/tM+SKx0vTXFOtz8fPLK9Gh8/rGugJhZTho5290eyUZrxCk04E/RsvrtxrB+5Fdsar+VNRN/JFqyjs5z/o7BLe8jUVxlWa1OsnAVxMrSV0bmOg0xO5wycgb72yPZaI04hRjjnOVf2trazN69e+0u47iSScPBUR8jM+FT+KIEdYOP0tq+m6qZdubKN9C37cMMb7qBZFGZdcXazO0WVpUWUVniOdJeqSzRXYLsNh/qdl5YdFId+UJEnjfGtB3vNR3BL8Gp9ttTG1jvobX925QH+ghWbWb/hV9htOVajMtjcbXZVepxp0bk6VF5ValH7/JcwCmtEbB/5AyF1R5xAg34k5gORnl5aGn9dncsSGPPj2g59F1K57z41p7Nizv/nfHGK3J+HXaR1J2eVfNhXurRXvkSOKk14oQLi4XUHnECbdEsYmAqxKExPyf7X+SOBWnqupONHd+lODrDVM1F9J75UabqdubkzUlutxzpkc8HeVWJzmBZLie0JI4dOR/7scpd2qI5RUvtt7viIZq7fsjGjm9THJlmov4yDr/mY8xWZ390tlwlHld6BovnyOi8vDj3/1poa+RoOnIuTLn/LznDltJvd8XnaOq+m43t36YkMslE/evpOetv8K0/L3uFnqL5NVhWzY/IS1Mj9Hxdp1xbI0fLxpxr5Twa8AtMB6O8NDR7wiV+XfEwjT330tq+m5LwOJN1O3nprE8wW70jy5Uuzu2SV27bL03dul9ZWlizWJww3xr0oqKylwZ82mL9dlciwoaeH7Hp4K2UhL1M1V7My5d8jZma12a/0GO4XXLkoueqstSfFboyIqCtEaUKPuCTScOBER+js6/ut0siyobD97Hp4H9QOjfGdHUb+y7+v0zXXmxDpa+E+aoyz5FQd1qYO6n3ra0RVegKOuDDsQQvDszgD8ePel4SUTb0/jQV7KERZqp3sP/CrzJde0nWZsUcG+arSj05sWa5U3rf2hpRqoCnSR6v3y7JGA29e9h04JuUhYaYWX8ePWd9gqm611ka7PPTEnMtzE/ECdMCnfSbhFJW0mmSx+ifTK0nM/+zTZJx6vseYNOBb1AeHGB23Tm0X/AFJusvy3iwWxXmTgo0J/S+tTWiVIEFfCJpOLig3y7JOHX9D7H5wDcoD/ThW3sWL5x/KxMZ2hZvfk2WqlIPq9LzzK0amTulNQLO6H0rpQoo4I/qtycT1A88zKYD36DCfxj/mjN48XXfZHzDm5cd7At75laH+fHotECl1LEKIuCn0uvJxGJx6gYeYfP+r1Ph78G/etuy1ooRgYqSVHtldbknvXJike09cye0RnRaoFLOkfcB3z8ZonNslpqBn7N5/9ep9HUSWHUaL13yb3ibrlxSsJd63KwqK2J1enS+qsxz1E1DTul/O6E1or1vpZwjbwM+kTQcHJ4hceBBLtz/dapmOwhWbebli/+FsearThjsRW450maZD/WT3c7vhP63tkaUUsfKy2mS4Wic/mfuo+GFr1E1c5Bg1SZ6XvNxxpqvBtcrYe1yQWWJ56jReUXJ8n7m2T010Cm/RSilsmuxaZKWBryIvA34V8ANfMcY85XFPv9UA/5VoWYMf/zFPdQ+/y80hTsIVbakgr3lGowrddFzVZnnSJhXlWZ2Cdx/frzjSP/7U1duy9j3VUqpE7FlHryIuIFvAFcAg8AfROQBY8yBTB3jSGvkveex07zA9MNfZMfMy8yWNtJ+8VcJn3kDq8pLOa8s1Tf3uK3bnMIJ/W+llFrIyh78hUCXMaYHQETuAa4DMhbwO7dU8x83bKXszquBQ4SoYWDHlzntil2cUVaaqcOclPa/lVJOZOV+a43AwIKPB9PPHUVEdonIXhHZOz4+fsoHueiMVjzVm/hs7MP8+OI9nPuOT1CWxXCHxacGKqWUXWzfUNMYs9sY02aMaaupqTnlr3+mZ5L3T/8VNW/4CLf/YYRnuicsqHJxH33DlleN1HduqdaLm0opW1kZ8ENA84KPm9LPZczC1sinrtx25E5OO0JeKaWcxsqA/wNwmohsEpFi4CbggUweQFsjSil1YpZdZDXGxEXk48DPSU2TvM0Ysz+Tx9C7JpVS6sQsvZPVGPMI8IiVx1BKKXV8tl9kVUopZQ0NeKWUylMa8Eoplac04JVSKk85ajVJERkH+pb55dVAvkyAz5dzyZfzAD0XJ8qX84CVnctGY8xx7xJ1VMCvhIjsPdGKarkmX84lX84D9FycKF/OA6w7F23RKKVUntKAV0qpPJVPAb/b7gIyKF/OJV/OA/RcnChfzgMsOpe86cErpZQ6Wj6N4JVSSi2gAa+UUnnK8QEvIreJiFdE9p3gdRGRfxORLhF5SUR2LHjtgyLSmf7vg9mr+vhWeC4JEXkh/V9Gl10+VUs4jzNE5FkRiYjI3x/z2ttEpCN9jp/JTsUntsJz6RWRl9PvydJ3i7fIEs7l/em/Vy+LyDMisn3Ba455X1Z4Hrn2nlyXPpcX0jvbXbrgtZXnlzHG0f8BlwE7gH0neP1q4FFAgIuB59LPrwN60n+uTT9em4vnkn4tYPd7cQrnUQu8Fvgn4O8XPO8GuoHNQDHwIvCaXDyX9Gu9QLXd78cpnMvO+X8DwFUL/q046n1Z7nnk6HtSySvXQs8F2tOPM5Jfjh/BG2N+DUwt8inXAbeblN8Da0SkAXgr8IQxZsoYMw08AbzN+opPbAXn4ignOw9jjNcY8wcgdsxLRzZiN8ZEgfmN2G2zgnNxnCWcyzPpfwsAvye1yxo47H1ZwXk4zhLOJWDSiQ5UAPOPM5Jfjg/4JTjR5t5L2vTbYRaruTT9K9zvReT6rFeWGbn4nizGAI+LyPMissvuYk7Rh0n9tgi5/b4sPA/IwfdERN4pIu3Aw8BfpJ/OyHti6YYfKqM2GmOGRGQz8KSIvGyM6ba7qAJ3afo9qQWeEJH29IjN0UTkTaSC8dKTfa6TneA8cu49McbcD9wvIpcBXwbekqnvnQ8j+BNt7m35pt8WOGHNxpj5P3uAXwHnZ7u4DMjF9+SEFrwnXuB+Uq0ORxORc4HvANcZYybTT+fc+3KC88jJ92Re+gfRZhGpJkPvST4E/APAB9IzUC4GZo0xI6T2gr1SRNaKyFrgyvRzTnbcc0mfQwlA+s1/HXDAzkKXyfKN2LNFRCpEpGr+Mam/X8edKeEUItIC/BT4c2PMoQUv5dT7cqLzyNH3ZKuISPrxDqAEmCRD+eX4Fo2I3A28EagWkUHg84AHwBjzLVJ7vl4NdAEh4EPp16ZE5Muk/vICfMkYs9gFTsst91yAM4FbRSRJ6ofyV4wxtgX8yc5DROqBvcAqICkinyQ1K8MnFm/EfqqWey6klne9P/1vswj4oTHmsayfwAJL+Pv1P4D1wDfTdceNMW3GmLiT3pflngdQR+69J+8mNaiLAXPAjemLrhnJL12qQCml8lQ+tGiUUkodhwa8UkrlKQ14pZTKUxrwSimVpzTglVIqT2nAq4IgIv9NRPYvWLnvIhH5pIiUL+N7fW4Jn/N9EblhedUqlRka8CrvicglwDXADmPMuaRuBR8APgkcN+BFxL3ItzxpwCvlBBrwqhA0ABPGmAiAMWYCuAHYADwlIk8BiEhARP6fiLwIXCIiN4vIf6ZH/LeKiFtEvgKUpZ+7K/11H0j/ZvCiiNyx4LiXSWq98h4dzSs76I1OKu+JSCXwW1Kj9V8A9xpjnhaRXqAtHfiIiCF1J+GPRORM4H8D7zLGxETkm8DvjTG3i0jAGFOZ/pqzSK15stMYMyEi69J3UX+f1PKvNwJnAA8YY7Zm9cRVwXP8UgVKrZQxJiAiFwCvB94E3CvH37UoAfwk/fjNwAXAH9K3vpcB3uN8zeXAj+d/SBxzO/keY0wSOCAidRk5GaVOgQa8KgjGmASpVTh/JSIvA8fbAi2c/jxI7ar1A2PMZ1dw2MiCx7KC76PUsmgPXuU9EdkmIqcteOo8oA/wA1Un+LJfAjek1xVHRNaJyMb0azER8aQfPwn8mYisn/+8TNev1HLpCF4Vgkrg30VkDRAntVrnLuC9wGMiMmyMedPCLzDGHBCRfyS1O5CL1JZ9HyP1g2E38JKI/NEY834R+SfgaRFJAH8CbsnSeSm1KL3IqpRSeUpbNEoplac04JVSKk9pwCulVJ7SgFdKqTylAa+UUnlKA14ppfKUBrxSSuWp/w8d6a8ryeeVswAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEHCAYAAACp9y31AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAAsTAAALEwEAmpwYAAATEklEQVR4nO3dcZBlZX3m8e8joICaAKHDTmaYNAZKykrCQNqJKZOUQdklkFWsdXdlTZZNUZlkVyq6MYlIUglWJVVYFUWTylqZCGF01YCgwqJJFhFjuZUFZ3ACA6OCOiaMI9NGCWClwIFf/jhnKp2Z7uk7Q5977ft+P1W35pz3ntPndzz49On3vvc9qSokSe141qQLkCSNl8EvSY0x+CWpMQa/JDXG4Jekxhw96QJGcfLJJ9fs7Oyky5CkVWXbtm3fqKqZA9tXRfDPzs6ydevWSZchSatKkq8u1m5XjyQ1xuCXpMYY/JLUGINfkhpj8EtSYwx+SWqMwS9JjTH4JakxBr8kNWZVfHNXkiZp9vKPTezYu666cMV/pnf8ktSYwYM/yVFJPpfk1n79tCR3JnkwyfVJnj10DZKkfzGOO/43ADsXrL8NuLqqTge+BVw6hhokSb1Bgz/JOuBC4D39eoBzgRv7TbYAFw1ZgyTpXxv6jv+dwG8CT/fr3wc8UlX7+vWHgLWL7ZhkU5KtSbbOz88PXKYktWOw4E/yc8Deqtp2JPtX1eaqmququZmZg54jIEk6QkMO53wp8MokFwDHAt8DvAs4IcnR/V3/OmD3gDVIkg4w2B1/Vb2lqtZV1SzwWuCTVfU64A7gNf1mlwA3D1WDJOlgkxjH/2bg15I8SNfnf80EapCkZo3lm7tV9SngU/3yl4GN4ziuJOlgfnNXkhpj8EtSYwx+SWqMwS9JjTH4JakxBr8kNcbgl6TGGPyS1BiDX5IaY/BLUmMMfklqjMEvSY0x+CWpMQa/JDVmLNMyS9JKmL38Y5MuYSp4xy9JjRnyYevHJrkryd8muS/JW/v265J8Jcn2/rVhqBokSQcbsqvnCeDcqno8yTHAZ5L8Rf/eb1TVjQMeW5K0hMGCv6oKeLxfPaZ/1VDHkySNZtA+/iRHJdkO7AVuq6o7+7d+P8k9Sa5O8pwl9t2UZGuSrfPz80OWKUlNGTT4q+qpqtoArAM2Jvlh4C3AmcCLgZOANy+x7+aqmququZmZmSHLlKSmjGVUT1U9AtwBnF9Ve6rzBPBnwMZx1CBJ6gw5qmcmyQn98nHAecDnk6zp2wJcBOwYqgZJ0sGGHNWzBtiS5Ci6XzA3VNWtST6ZZAYIsB34lQFrkCQdYMhRPfcAZy/Sfu5Qx1zMJL/pt+uqCyd2bElait/claTGGPyS1BiDX5IaY/BLUmMMfklqjMEvSY0x+CWpMQa/JDXG4Jekxhj8ktQYg1+SGmPwS1JjDH5JaozBL0mNMfglqTEGvyQ1ZshHLx6b5K4kf5vkviRv7dtPS3JnkgeTXJ/k2UPVIEk62JB3/E8A51bVWcAG4PwkLwHeBlxdVacD3wIuHbAGSdIBBgv+6jzerx7Tvwo4F7ixb99C98B1SdKYDPmwdfoHrW8DTgf+GPgS8EhV7es3eQhYu8S+m4BNAOvXrx+yTGlVmuTzpLW6DfrhblU9VVUbgHXARuDMw9h3c1XNVdXczMzMUCVKUnPGMqqnqh4B7gB+Ajghyf6/NNYBu8dRgySpM+SonpkkJ/TLxwHnATvpfgG8pt/sEuDmoWqQJB1syD7+NcCWvp//WcANVXVrkvuBP0/ye8DngGsGrEGSdIDBgr+q7gHOXqT9y3T9/RpIix/67brqwokct8X/rbX6+c1dSWqMwS9JjTH4JakxBr8kNcbgl6TGGPyS1BiDX5IaY/BLUmMMfklqjMEvSY0x+CWpMQa/JDXG4Jekxhj8ktSYkYI/yY8MXYgkaTxGveP/X0nuSvI/knzvoBVJkgY1UvBX1U8BrwNOBbYl+UCS8w61T5JTk9yR5P4k9yV5Q99+ZZLdSbb3rwue8VlIkkY28hO4quqBJL8NbAX+EDg7SYArqurDi+yyD3hTVd2d5Pl0vzBu69+7uqr+4JkWL0k6fCMFf5IfBX4RuBC4Dfj3faD/APA3wEHBX1V7gD398mNJdgJrV6pwSdKRGbWP/4+Au4Gzqur1VXU3QFV9Dfjt5XZOMkv3/N07+6bLktyT5NokJx5+2ZKkIzVq8F8IfKCq/gkgybOSHA9QVe871I5JngfcBLyxqh4F3g38ELCB7i+Cty+x36YkW5NsnZ+fH7FMSdJyRg3+TwDHLVg/vm87pCTH0IX++/d/DlBVD1fVU1X1NPCnwMbF9q2qzVU1V1VzMzMzI5YpSVrOqMF/bFU9vn+lXz7+UDv0H/xeA+ysqncsaF+zYLNXAztGL1eS9EyNOqrn20nO2d+3n+THgH9aZp+XAr8A3Jtke992BXBxkg1AAbuAXz7MmiVJz8Cowf9G4ENJvgYE+DfAfz7UDlX1mX7bA338cAqUJK2skYK/qj6b5EzghX3TF6rqO8OVJUkayshf4AJeDMz2+5yThKp67yBVSZIGM+oXuN5HNwRzO/BU31yAwS9Jq8yod/xzwIuqqoYsRpI0vFGHc+6g+0BXkrTKjXrHfzJwf5K7gCf2N1bVKwepakrMXv6xSZcgSQcZNfivHLIISdL4jDqc86+T/CBwRlV9op+n56hhS5MkDWHURy/+EnAj8Cd901rgowPVJEka0Kgf7r6ebgqGR6F7KAvw/UMVJUkazqjB/0RVPbl/JcnRdOP4JUmrzKjB/9dJrgCO65+1+yHg/wxXliRpKKMG/+XAPHAv3WyaH2eEJ29Jkr77jDqqZ/9DU/502HIkSUMbda6er7BIn35VvWDFK5IkDepw5urZ71jgPwInrXw5kqShjdTHX1X/sOC1u6reSfcAdknSKjNqV885C1afRfcXwCH3TXIq3bTNp9B1E22uqnclOQm4nm5u/13Af6qqbx125ZKkIzJqV8/bFyzvow/sZfbZB7ypqu5O8nxgW5LbgP8G3F5VVyW5nG7E0JsPq2pJ0hEbdVTPzxzuD66qPcCefvmxJDvppnp4FfCyfrMtwKcw+CVpbEbt6vm1Q71fVe9YZv9Z4GzgTuCU/pcCwNfpuoIW22cTsAlg/fr1o5QpSRrBqF/gmgP+O90d+1rgV4BzgOf3ryUleR5wE/DGqnp04Xv9E70WnfqhqjZX1VxVzc3MzIxYpiRpOaP28a8DzqmqxwCSXAl8rKp+/lA7JTmGLvTfX1Uf7psfTrKmqvYkWQPsPbLSJUlHYtQ7/lOAJxesP8kSXTT7JQlwDbDzgK6gW4BL+uVLgJtHrEGStAJGveN/L3BXko/06xfRfTB7KC8FfgG4N8n2vu0K4CrghiSXAl9l+dFBkqQVNOqont9P8hfAT/VNv1hVn1tmn88AWeLtl49eorQ8n28sjW7Urh6A44FHq+pdwENJThuoJknSgEZ99OLv0o21f0vfdAzwv4cqSpI0nFHv+F8NvBL4NkBVfY1lhnFKkr47jRr8Ty4cc5/kucOVJEka0qjBf0OSPwFOSPJLwCfwoSyStCotO6qnH49/PXAm8CjwQuB3quq2gWuTJA1g2eCvqkry8ar6EcCwl6RVbtSunruTvHjQSiRJYzHqN3d/HPj5JLvoRvaE7o+BHx2qMEnSMJZ7itb6qvo74N+NqR5J0sCWu+P/KN2snF9NclNV/Ycx1CRJGtByffwL59p5wZCFSJLGY7ngryWWJUmr1HJdPWcleZTuzv+4fhn+5cPd7xm0OknSijtk8FfVUeMqRJI0HoczLbMkaQoMFvxJrk2yN8mOBW1XJtmdZHv/umCo40uSFjfkHf91wPmLtF9dVRv618cHPL4kaRGDBX9VfRr45lA/X5J0ZCbRx39Zknv6rqATl9ooyaYkW5NsnZ+fH2d9kjTVxh387wZ+CNgA7AHevtSGVbW5quaqam5mZmZM5UnS9Btr8FfVw1X1VFU9Tfcgl43jPL4kaczBn2TNgtVXAzuW2laSNIxRp2U+bEk+CLwMODnJQ8DvAi9LsoFu+oddwC8PdXxJ0uIGC/6quniR5muGOp4kaTR+c1eSGmPwS1JjDH5JaozBL0mNMfglqTEGvyQ1xuCXpMYY/JLUGINfkhpj8EtSYwx+SWqMwS9JjTH4JakxBr8kNcbgl6TGGPyS1JjBgj/JtUn2JtmxoO2kJLcleaD/98Shji9JWtyQd/zXAecf0HY5cHtVnQHc3q9LksZosOCvqk8D3zyg+VXAln55C3DRUMeXJC1u3H38p1TVnn7568ApS22YZFOSrUm2zs/Pj6c6SWrAxD7craoC6hDvb66quaqam5mZGWNlkjTdxh38DydZA9D/u3fMx5ek5o07+G8BLumXLwFuHvPxJal5Qw7n/CDwN8ALkzyU5FLgKuC8JA8Ar+jXJUljdPRQP7iqLl7irZcPdUxJ0vL85q4kNcbgl6TGGPyS1BiDX5IaY/BLUmMMfklqjMEvSY0x+CWpMQa/JDXG4Jekxhj8ktQYg1+SGmPwS1JjDH5JaozBL0mNMfglqTGDPYjlUJLsAh4DngL2VdXcJOqQpBZNJPh7P1NV35jg8SWpSXb1SFJjJhX8BfzfJNuSbJpQDZLUpEl19fxkVe1O8v3AbUk+X1WfXrhB/wthE8D69esnUaMkTaWJ3PFX1e7+373AR4CNi2yzuarmqmpuZmZm3CVK0tQae/AneW6S5+9fBv4tsGPcdUhSqybR1XMK8JEk+4//gar6ywnUIUlNGnvwV9WXgbPGfVxJUsfhnJLUGINfkhpj8EtSYwx+SWqMwS9JjTH4JakxBr8kNcbgl6TGGPyS1BiDX5IaY/BLUmMMfklqjMEvSY0x+CWpMQa/JDXG4Jekxhj8ktSYiQR/kvOTfCHJg0kun0QNktSqSTxs/Sjgj4GfBV4EXJzkReOuQ5JaNYk7/o3Ag1X15ap6Evhz4FUTqEOSmjT2h60Da4G/X7D+EPDjB26UZBOwqV99PMkXjvB4JwPfOMJ9VyvPuQ2ecwPyNuDIz/sHF2ucRPCPpKo2A5uf6c9JsrWq5lagpFXDc26D59yOlT7vSXT17AZOXbC+rm+TJI3BJIL/s8AZSU5L8mzgtcAtE6hDkpo09q6eqtqX5DLgr4CjgGur6r4BD/mMu4tWIc+5DZ5zO1b0vFNVK/nzJEnf5fzmriQ1xuCXpMZMbfC3MC1EklOT3JHk/iT3JXlD335SktuSPND/e+Kka11pSY5K8rkkt/brpyW5s7/e1/cDB6ZKkhOS3Jjk80l2JvmJab/WSf5n/9/2jiQfTHLstF3rJNcm2Ztkx4K2Ra9rOn/Yn/s9Sc45kmNOZfA3NC3EPuBNVfUi4CXA6/vzvBy4varOAG7v16fNG4CdC9bfBlxdVacD3wIunUhVw3oX8JdVdSZwFt35T+21TrIW+FVgrqp+mG4wyGuZvmt9HXD+AW1LXdefBc7oX5uAdx/JAacy+GlkWoiq2lNVd/fLj9EFwVq6c93Sb7YFuGgiBQ4kyTrgQuA9/XqAc4Eb+02m8Zy/F/hp4BqAqnqyqh5hyq813cjD45IcDRwP7GHKrnVVfRr45gHNS13XVwHvrc7/B05IsuZwjzmtwb/YtBBrJ1TLWCSZBc4G7gROqao9/VtfB06ZVF0DeSfwm8DT/fr3AY9U1b5+fRqv92nAPPBnfRfXe5I8lym+1lW1G/gD4O/oAv8fgW1M/7WGpa/rimTbtAZ/U5I8D7gJeGNVPbrwverG607NmN0kPwfsraptk65lzI4GzgHeXVVnA9/mgG6dKbzWJ9Ld4Z4G/ADwXA7uEpl6Q1zXaQ3+ZqaFSHIMXei/v6o+3Dc/vP/Pv/7fvZOqbwAvBV6ZZBddF965dH3fJ/TdATCd1/sh4KGqurNfv5HuF8E0X+tXAF+pqvmq+g7wYbrrP+3XGpa+riuSbdMa/E1MC9H3bV8D7Kyqdyx46xbgkn75EuDmcdc2lKp6S1Wtq6pZuuv6yap6HXAH8Jp+s6k6Z4Cq+jrw90le2De9HLifKb7WdF08L0lyfP/f+v5znupr3Vvqut4C/Nd+dM9LgH9c0CU0uqqayhdwAfBF4EvAb026noHO8Sfp/gS8B9jevy6g6/O+HXgA+ARw0qRrHej8Xwbc2i+/ALgLeBD4EPCcSdc3wPluALb21/ujwInTfq2BtwKfB3YA7wOeM23XGvgg3WcY36H7y+7Spa4rELoRi18C7qUb8XTYx3TKBklqzLR29UiSlmDwS1JjDH5JaozBL0mNMfglqTEGvyQ1xuCXBtbCFOFaXRzHLw2onyL8i8B5dF/O+SxwcVXdP9HC1DTv+KUFksz2Dzq5LskXk7w/ySuS/L/+oRgbk1yZ5NcX7LOjnx11MU1MEa7VxeCXDnY68HbgzP71X+imx/h14IrD/FnNTRGu734Gv3Swr1TVvVX1NHAf3ZOQim5ulNmJViatAINfOtgTC5afXrD+NN28+Pv41//fOfYQP6uZKcK1ehj80uHbRTcXPv3Drk87xLZNTBGu1cXglw7fTcBJSe4DLqMbtbOo6h4ReBnwV3TPRL6hqu4bS5XSEhzOKUmN8Y5fkhpz9PKbSFpOkv1PTDrQy6vqH8Zdj3QodvVIUmPs6pGkxhj8ktQYg1+SGmPwS1Jj/hl5tynsUxk7jgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "random_params = pmt.RandomParameters(sample_params)\n",
    "print(random_params)\n",
    "random_samples = random_params.sample(300)\n",
    "print('random_samples are a list of parameters dictionaries:', type(random_samples),len(random_samples))\n",
    "results = np.array([sim_func(s) for s in random_samples])\n",
    "print('results are of shape', results.shape)\n",
    "\n",
    "#We can plot the mean and variation of these\n",
    "mean = np.mean(results,axis=0)\n",
    "sd = np.std(results,axis=0)\n",
    "plt.plot(applied_stretch,force_measured,'x',label='Observed')\n",
    "plt.plot(applied_stretch,mean,'-')\n",
    "plt.fill_between(applied_stretch,mean-sd,mean+sd,alpha=0.3)\n",
    "plt.xlabel('Stretch')\n",
    "plt.ylabel('Force')\n",
    "plt.show()\n",
    "\n",
    "mu_samples = [s['mu_0'] for s in random_samples]\n",
    "plt.hist(mu_samples,bins=10)\n",
    "plt.xlabel('mu_0')\n",
    "plt.ylabel('Frequency')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ee207834",
   "metadata": {},
   "source": [
    "## Random parameters: Markov Chain Monte Carlo simulation\n",
    "Or we can perform a Markov Chain Monte Carlo (MCMC) simulation through some likelihood function. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "4d9a831e",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "MCMC instance created with the following settings\n",
      "------------------------------------------------------------------\n",
      "Keys              Value       Fixed?      Lower bound Upper bound \n",
      "------------------------------------------------------------------\n",
      "L0                10.00       Yes         -           -           \n",
      "A0                0.10        Yes         -           -           \n",
      "mu_0              11.51       No          1.00e-04    1.00e+02    \n",
      "------------------------------------------------------------------\n",
      "\n",
      "1 parameters will be varied.\n"
     ]
    }
   ],
   "source": [
    "def prob(params):\n",
    "    force_modeled = sample.disp_controlled(applied_stretch,params)\n",
    "    diff = force_modeled - force_measured\n",
    "    diff_norm2 = np.dot(diff,diff)\n",
    "    var = 2.\n",
    "    prob = np.exp(-diff_norm2/2./var) #in the case of MCMC we can skip the proportionality factor\n",
    "    return prob,force_modeled\n",
    "\n",
    "mcmc = pmt.MCMC(prob,sample_params)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "012a657d",
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|████████████████████████████████████████████████████████████████████████████████████████████| 1000/1000 [00:00<00:00, 2575.02it/s]"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "MCMC sampling completed. Acceptance rate: 0.812\n",
      "Number of samples: 812\n",
      "To access the samples, use get_samples()\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "\n"
     ]
    }
   ],
   "source": [
    "mcmc.run(1000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "55fcf9e0",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(812, 10)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAAsTAAALEwEAmpwYAAAu80lEQVR4nO3deZhcdZno8e9bVb1vWTrp7vSSdBYSFlkbQUQWB7iBUUHFEXAXjdwRx3HuzNU741UvzjyPc5e5V8UZjNzIgEoYRJggQUQY4CrCEFSWdJJO7/u+VlV31/beP6o6dmJ6r9O1vZ/n6SdV55w69Z5U93nr/H6/8/5EVTHGGJO5XIkOwBhjTGJZIjDGmAxnicAYYzKcJQJjjMlwlgiMMSbDeRIdwFKVlpbqli1bEh2GMcaklFdffXVQVTecbl3KJYItW7Zw6NChRIdhjDEpRUTa5lpnTUPGGJPhLBEYY0yGs0RgjDEZzhKBMcZkOEsExhiT4SwRGGNMhrNEYIwxGc6xRCAi+0SkX0TeXGC7i0UkJCI3OxWLMcaYuTl5RXAfsHu+DUTEDfw98HMH4zDGmJQ2FQzTMex3bP+OJQJVfQEYXmCzzwGPAP1OxWGMMakqGI5wvG+CF5sG6R2fcux9ElZiQkQqgfcCVwMXL7DtHmAPQE1NjfPBGWNMAoUjSsewn9YhH6Gw87NIJrLW0P8BvqiqERGZd0NV3QvsBairq7O5NY0xaUlV6R6bonnAy3Qwsmrvm8hEUAfsjyWBUuAGEQmp6mMJjMkYYxKif3yKxgEv/unwqr93whKBqtbOPBaR+4CfWhIwxmSaYV+Axn4v45PBhMXgWCIQkQeBq4BSEekEvgpkAajqPU69rzHGpIKJqSDH+70MewOJDsW5RKCqty5h2487FYcxxiSTyUCYpgEvvWPOjQJaqpSbmMYYY1LRdChM66CfrlE/kdXrB14USwTGGOOgUDhC27Cf9mE/4VUYCroclgiMMcYBkYjSOTJJy5CPYCjJLgFOYYnAGGPiSFXpHZ+iecDHZGD1h4IuhyUCY4yJk0HvNI39XrxToUSHsiSWCIwxZoXG/EEaByYY8SXuXoCVsERgjDHL5JsO0djvZWBiOtGhrIglAmOMWaKpYJjmAR89Y5Nocg4EWhJLBMYYs0jBcITWQR8dI8l3L8BKWCIwxpgFhCNK+7CftlUqC73aLBEYY8wcElUWerVZIjDGmNMY8k5zPAWHgi6HJQJjjJnFOx2ioW8iKaqCrhZLBMYYQ/qNBFoKSwTGmIwWjihtQz7ahvyEIxmWAWIsERhjMpKq0jM2RVOadwQvhiUCY0zGyaSO4MWwRGCMyRje6RDH+yYYyqCO4MWwRGCMSXvToWhHcPdo5nUEL4bLqR2LyD4R6ReRN+dY/yEReV1E3hCRF0XkPKdiMcZkpnBEaRn08WLTEF0jlgTm4lgiAO4Dds+zvgW4UlXfAnwd2OtgLMaYDNMzNsmLTYM09XuTdorIZOFY05CqviAiW+ZZ/+Kspy8BVU7FYozJHMO+AMf7JpiwjuBFS5Y+gtuBJ+daKSJ7gD0ANTU1qxWTMSaF+KZDHO/3MpjicwMkQsITgYhcTTQRXD7XNqq6l1jTUV1dnV3jGWNOCIQiNA96rQ9gBRKaCETkXOBe4HpVHUpkLMaY1BKJlYZuGfJZH8AKJSwRiEgN8BPgI6rakKg4jDGpp2dskqZ+H1PBcKJDSQuOJQIReRC4CigVkU7gq0AWgKreA3wFWA/8o4gAhFS1zql4jDGpb8QX4Hi/l/HJ1JwkPlk5OWro1gXWfwr4lFPvb4xJH/5AiON9qT9JfLJKeGexMcbMJRCK0DLoo2s0veYITjaWCIwxSScSUTpHJmke9KblHMHJxhKBMSapDHqnaeibwD9tHcGrxRKBMSYp+AMhGvrshrBEsERgjEmoUDjaD9AxYv0AiWKJwBiTMN2jkzT2ewmELAMkkiUCY8yqG/MHOdY3YfcDJAlLBMaYVTMVDNPY76V3bCrRoZhZLBEYYxwXiShtw35arS5QUrJEYIxxVP/EFMf7vEwGbDhosrJEYIxxhHc6REPfBMM2UXzSs0RgjImrYDhC84CPzhG/zQ+QIiwRGGPiQnWmLISPoA0HTSmWCIwxKzbiC3CsbwKvzROckiwRGGOWbSoY5nifl75xGw6ayiwRGGOWLBxRWod8tA/5CUesIyDVWSIwxixJ79gUjf1emyYyjVgiMMYsyvhUkIbeCUb9VhYi3VgiMMbMKxCK0NjvpWds0oaDpimXUzsWkX0i0i8ib86xXkTkWyLSKCKvi8iFTsVijFk6VaV9yM+LTYN0j1oSSGeOJQLgPmD3POuvB3bEfvYA/+RgLMaYJRjyTvPr5iEa+iZsqsgM4FjTkKq+ICJb5tnkRuB+VVXgJRFZIyIVqtrjVEzGmPlNBcM09E3QP26zhGWSRPYRVAIds553xpb9QSIQkT1ErxqoqalZleCMySSRiNI+7Kdl0GfDQTOQk01DcaOqe1W1TlXrNmzYkOhwjEkrI74AL7cM09jvtSSQoRJ5RdAFVM96XhVbZoxZBdOh6F3BNkmMSeQVwQHgo7HRQ5cCY9Y/YIzzVJWOYT+/bhqyJGAAB68IRORB4CqgVEQ6ga8CWQCqeg9wELgBaAT8wCecisUYEzXmD3K0d5wJKw5nZnFy1NCtC6xX4LNOvb8x5veC4QjH+7x0j04mOhSThOzOYmPSXNfoJI39XpsjwMzJEoExaWpiKsjR3gnGrDaQWYAlAmPSTCgcocmmijRLYInAmDTSOzZFQ98EAWsGMktgicCYNOCbDnG0d4IRXyDRoZgUZInAmBQWjigtg17ah/1E7CLALJMlAmNSVP/EFA29NlOYWTlLBMakmMlAmKO94wx5rRkoU2RND7Oh4xlwXQw1l8R9/5YIjEkRkdiE8a1DPmsGygCewBgbO5+mrOMga/t/jUvDoJ+1RGBMphr0TtPQO4E/YM1A6cwd9LKh+xnK2g+yvu+XuCJB/AVVtO38FL4z3sM5F7zdkfe1RGBMErOJYtKfK+RnQ/e/UdZxkPU9z+OOBJjKr6Bj+0foq7mB8bVvARFK8rNAxJEYLBEYk4Rsopj05gpNsb73eco6DrKh+znc4UmmczfSte0W+qpvYGz9+SCrVxzaEoExSWbEF+Bo7wS+aasQmk4kHGB9368o63iCDV3P4An5COSspWfLTfRV38BIaR243AmJzRKBMUkiEIrQ0DdhcwSkEYkEWdv/EuXtB9nQ9TRZwXGC2SX0VV9PX/UfM7LxEtQ1/2n4yTd7qC0t4JKt608se7FpkNc7x7jjym1xidMSgTFJoGdskoY+qxCaFiJh1g78O2UdB9nY+RTZgVFCWYX0b7qGvpobGN54GerOXvTuaksLuOf5ZgpyPFy8ZR0vNg1y549+y923XRC3kC0RGJNAk4EwR3rHGbZ7AlKbRigZ/A3lsZN/ztQAIU8+gxVX01dzA0PlVxBx5yxr17vKi7njyq38w9MNjPgC/ODldu6+7QIu21Yat/AtERiTAKrRzuDmAesMTlmqFA+/TlnHE5R1/IzcyV7C7hwGK66ir/oGBiuuIuLJi8tb7Sov5rqzy/jWs4382Tu3xzUJgCUCY1bd+FSQI902XWRKUqVotJ6y9oOUdT5Jnq+TiCuLofJ3cPzcv2Jw09WEswrj/rZHe8f5+eE+/uyd2/nBy+1cum29XREYk4rCEaV5IFogzuYJSC0FYw3Rk3/HQQq8rUTEw3DZ22g+604GKq8hlF3s2Hsf7R3nnueb+U/XncEn3l7LpdvWn+gjiFcycDQRiMhu4JuAG7hXVb9xyvoa4J+BNbFtvqSqB52MyZhEGPJOc7R3gkm7Mzhl5Po6KW//KeVtj1M4fhwVF8MbLqFt5+0MVF1HMGftqsTRMujjjiu3ck5lCQCXbSvl7tsu4PXOsbglAlGHvpqIiBtoAK4FOoFXgFtVtX7WNnuB36rqP4nIWcBBVd0y337r6ur00KFDjsRsTLzZkNDUkjU9TFnHk5S3P86awd8AMFp6Ib3V76K/ejeB3Pi2zS9FSX4WF29Zt+zXi8irqlp3unVOXhG8FWhU1eZYEPuBG4H6WdsoMHNNVQJ0OxiPMauqd2yKY30TNiQ0yblCfjZ0PUN5++Os7/0lLg3hLd5B41v+gt6adzFVUJXoEB3nZCKoBDpmPe8ETi2b9zXg5yLyOaAAuOZ0OxKRPcAegJqamrgHakw8WZno5CeRIOv6fkV52+Ns7PoF7vAkU/kVtJ/xCXo3vxtvyU7H6voko0R3Ft8K3Keq/0tE3gY8ICLnqOpJX6FUdS+wF6JNQwmI05gFqSodw5M0DXhtSGgyUqVk6LeUtx2grPNJsqdHCGaX0LP5PfRufg+jpRetan2fZLLoRCAim4EdqvoLEckDPKo6Mc9LuoDqWc+rYstmux3YDaCqvxaRXKAU6F9sXMYkg4mpIEd6JhifDCY6FHOKgrFGytsPUN7+U/J8nYTdOQxseie9Ne9hqPwdS7rLN10tKhGIyKeJNs2sA7YRPanfA/zRPC97BdghIrVEE8AtwG2nbNMe28d9InImkAsMLOUAjEmkmTmD24ZsSGgyyfH3UN7+BOXtBygaPYqKi6Gyt9N89ufor7zWkbH+qWyxVwSfJdr5+zKAqh4XkY3zvUBVQyJyJ/AU0aGh+1T1sIjcBRxS1QPAfwK+JyJfINpx/HF1ahiTMXE27AtwtGfcJotJEp7pUTZ2PkVF++OsGXgFQRlbdx7Hzv8yfTU3LHnEz0yxt13lv79H4GjvOC2DPq4/pyLe4SfUYhPBtKoGJNZ5IiIeoifuecXuCTh4yrKvzHpcDzgz5Y4xDgmGo0NCe0ZtSGiiuUJTlPY8S3nb45T2voArEsRXVEvz2Z+jt+bdTBZtXva+Z4q93XHlVnaVF5+4seuOK7fG8QiSw2ITwfMi8tdAnohcC/wp8LhzYRmTnGxIaOJJJBQt7dx2gI1dT+MJ+ZjO3UjH9g/TW/NuJtaeHZcRPzPF3u55vpmrdm7guWMDJ5JCullsIvgS0Y7dN4DPEP2Wf69TQRmTbKaCYY702JDQhIkVeCtvf5yyjoPkTA0Syiqkr/p6emvezciGtzoyqcuu8mKu2rmBn77ew7vOrUjLJACLTwR5RNv4vwcn7hrOA/xOBWZMMrAhoYmVP9FCeVt0xE++t42IK4vBiqvp2fxuhiquWnZp58U62jvOc8cGeNe5FTx3bIBd5UVpmQwWmwieIXqzlzf2PA/4OXCZE0EZkwxsSGhiZE2PUNb+BBVtj1Iy/AaKMLLxElp3fYb+quscLfA22+w+gV3lxewqLzrpeTpZbCLIVdWZJICqekUk36GYjEmoSERptiGhq0rCAUp7n6ei9VFKe57HFQkysWYXDed9kd6adxHIK1v1mGaKvc2c9Gf6DFoGfRmbCHwicqGq/gZARC4CJp0Ly5jEsCGhq0iV4pE3qGh9jLL2n5IdGGU6t5SO7R+mZ8tNeNecmdDwTjdENHplkF5JABafCD4PPCwi3YAA5cAHHYvKmFVmQ0JXT46/l/K2f6Wi7TEKx5sIu7IZqLyGni3vZbjs7QtO5m7ib8H/8VjH8DuAXcDO2OJjqmoNpyYt9E9McbRngoANCXWMK+RnY+fPqWh7jHV9v0ZQRkovov6ir9Nfff2qtfub01swEahqWERuVdX/Dby5CjEZsyqC4QjHem2uAMdohLUD/05F66Ns7HwKT8iPv6CKlrM+S8+Wm5gstErCyWKx12C/EpG7gYcA38zCmT4DY1LNwMQ0R3vHmQ7aVUC85Y83U9H2r5S3/St5/u7YeP8/pmfLTRld4TOZLTYRnB/7965ZyxR4Z1yjMcZh1hfgDM/0KGUdB9nU+iglw6+dKPLWeO5fMrDpGiKe3ESHaOaxqESgqlc7HYgxThv0TnOkx64C4kUiQdb3vEBF22Ns6H4WVySIt+QMGs79Ir2b300gb966lCaJLLYMdQnwVeCK2KLngbtUdcypwIyJl1A4QkOfl+5RG/G8YqoUjRymou0xytsfJ3t6hEDOOjq3fYieLTcxsebMFdf5yaSqn8lisU1D+4h2FP9J7PlHgO8D73MiKGPiZcg7zZGeCaaCdl/ASmRP9lHRdoCK1scoHD9OxJXFwKY/omfLexkqvxx1ZcXtvTKp6meyWGwi2Kaq75/1/L+JyO8ciMeYuLCrgJVzhSbZ0PULKtoeY33frxCNMLr+Ao5cdBd91dcTyi5x5H0zqepnslhsIpgUkctV9ZcAIvJ27M5ik6SGfQHqu8ftKmA5VCkZfJVNrT+hrONJPCEfk/mVtOy6Izrks2jLqoSRKVU/k8ViE8EdwP2xvgKAEeBjzoRkzPKEwhGO93vpGrHvKEuVPTlAReujbGp9hIKJFkKeAvqrdtOz5SZGNly86kM+M6XqZ7KYNxGISI2qtqvqa8B5IlIMoKrjqxKdMYs07AtwpGecSasRtGgSCVLa8zybWh5hfc9zuDTMSOlFtO76DH3Vu4l4ElNXMpOqfiaLha4IHgMuBBCRR07pJzAm4cIR5Xj/BJ3DdhWwWPnjzWxqeYSKtsfImRpgOncD7Ts/SXftzfiLahMdXkZV/UwWCyWC2ePAltxlLyK7gW8Snbz+XlX9xmm2+RPga0RvUHtNVW9b6vuYzDTiC1BvVwGL4g762Nj5MypbHmbN4G+IiJvBiqvorr2ZoYor4jrqZ6UyqepnslgoEegcjxcUK1b3HeBaoBN4RUQOxCasn9lmB/BfgLer6oiI2B0oZkHhiNLY76Vj2CbIm5cqxcOvUdn8MGUdB/GEfPiKajl+7l/Rs/kmAnkbEh2hSRILJYLzRGSc6JVBXuwxseeqqvOl6LcCjaraDCAi+4EbgfpZ23wa+I6qjhDdYf8yjsFkkFF/dESQzRcwt6ypISraHmNTyyMUjjcSdufRV309XbU3M1Z6UVwmdjfpZd5EoKormQ26EuiY9bwTuOSUbc4AEJFfEW0++pqq/uzUHYnIHmAPQE2NVSzMROGI0jQQvQqwWcNOIxJmfd//o7L5x5R2P4tLQ4yuP5/6ur+lr/oGwlmFC+7C7ujNXImeAcID7ACuAqqAF0TkLao6OnsjVd0L7AWoq6uz00CGGfMHOdw9ZlcBp5HnbWdTy4+paH2U3Mk+Ajnr6NjxUbprb8ZXsn1J+7I7ejOXk4mgC6ie9bwqtmy2TuDl2CQ3LSLSQDQxvOJgXCZFRGJXAe12FXASV2iSjV0/Z1Pzj1k38HK00mf5Ozh2wZcZrLgadWcva792R2/mcjIRvALsEJFaogngFuDUEUGPAbcC3xeRUqJNRc0OxmRSxJg/yOGeMfzTdhUAxIq9vUlly48pa/8pWcEJ/IU1NJ7zBXq2vJfp/PK4vI3d0ZuZHEsEqhoSkTuBp4i2/+9T1cMichdwSFUPxNZdJyL1QBj4K1Udciomk/wiEaV50EvbkF0FAGRNj1De9jibWh6maOwYYXcO/VX/ge7aDzhyx6/d0ZuZRFPsr62urk4PHTqU6DCMA8Ymg9R3j+ObDiU6lMTSCOv6XmRTy8Ns7PoFrkiQsbXn0F37Afpq/tix+X1PvaP31OcmsUrys7h4y7plv15EXlXVutOtS3RnsTGxqwAfbUO+jL4KyPH3UNn8MBWtPyHP300gew2d226lu/ZmvGt2Of7+dkdv5rJEYBJqfCrI4a4MvgqIDfusatpPac9zoMpw2ds5ft5/ZmDTNcvu+F0Ou6M3c1kiMAmhqrQN+Wke9BLJwJkjsyf72dTyCJXND5Hn72Y6t5TWXXvo2vonTBVUJTo8k2EsEZhVNxUMc7h7nBFfINGhrC6NsK7/11Q27WdD1zO4NMTQxrdx/LwvMVD5R0lV78dkFksEZlX1jU9xpGecUDhzOgOypobZ1PoTKpv3k+9tJ5C9hvYzPkbX1g+u2kQvxszHEoFZFaFwhKO9E/SOTSU6lNWhyprBQ1Q1PcjGzqdwRYKMlNbRdPbnGai6jog7J9ERGnOCJQLjuFF/gMPdmVEu2hMYo6L1MSqb91M43kQwq4jObbfStfWWeUs+WJ0fk0iWCIxjVKPDQlsH03xYaKzcc1XTg5R1HMQdnmZs3Xkcvvgb9FVfT8STt+AurM6PSSRLBMYR/kCIw93jjPmDiQ7FMe6gl/K2A1Q176do9CghTwE9m99L17ZbmFh71pL2ZXV+TCJZIjBx1z06ybG+CcJp2iFcNHKYyqYHKW//KZ6Qn4k1Z3LkorvorXnXoso9z8Xq/JhEsURg4iYYjnCkZ5z+8elEhxJ3rpCf8vYnqGzeT8nwG4TdufRV/zGd225hfN25cZnsxer8mESxRGDiYtgX4HD3GNPB9Lo7rGCsgaqm/VS0PYYn6MVbvJ1jF3yZns03xbXmz6l1fXaVF1mdH7NqLBGYFZmZM6BtKH3mD3aFp9nY8SRVzfujE727suirup7Obbc4NtWj1fkxiWSJwCybdzrEm11jeKfSo05Q/ngzlc0PUdH6KNmBUXyFW2g474v0bHkvwZzlV31cDKvzYxLJEoFZlo5hP439XsKR1O4QlkiI0u5nqG78Iev6XyIiHgYqr6Vz2y2MbLwk7vX+jUlGlgjMkkyHwtR3jzPkTe06QVnTw2xqfpiqpgfJ83czlV9B4zlfoLv2ZgJ5GxIdnjGryhKBWbSBiWnqe8YJhlK3Q7hwpJ7qxh9Q3v447vA0wxsvpeH8v2Zw0ztRl/05mMxkv/lmQeGI0tA3QdfIZKJDWRaJBNnQ9TTVxx9g7eCrhN159Gx+Lx07Poyv5IxEh2dMwlkiMPManwryZldqTiKfNTVEVfNDVDY9SO5kH/6CahrO+xLdte8nlF2S6PCMSRqOJgIR2Q18k+jk9feq6jfm2O79wI+Bi1XVJiROAqk8cUzx8OtUH3+Aso6DuCJBhsou5+hFdzFYfgW43Cdta8XejHEwEYiIG/gOcC3QCbwiIgdUtf6U7YqAzwMvOxWLWZroxDFjjPhSp06QhAOUdf6M6uMPUDL8GiFPAV1bP0jH9g/hL9425+us2Jsxzl4RvBVoVNVmABHZD9wI1J+y3deBvwf+ysFYzCKl2sQx2ZP9VDU9SGXzQ+RMDeIrquXYBV+me8v7FlX3x4q9GeNsIqgEOmY97wQumb2BiFwIVKvqEyIyZyIQkT3AHoCamhoHQjUpNXGMKiVDv6X6+ANs7HwK0TBDFVdQv/0jDJVfvuSx/1bszWS6hHUWi4gL+Afg4wttq6p7gb0AdXV1qfFVNYWkysQxrvA0ZR1PUH38BxSPvEkoq5DO7R+mY/uHmCzavOz9WrE3k+mcTARdQPWs51WxZTOKgHOA5yRau6UcOCAi77EO49WhqjQN+GgbSu6JY3L8vVQ1/YjK5ofInh7BW7ydIxd+jd7NNxLOKljRvq3Ym0lWLhfkZXnIz3ZTkOOmOC/LsfdyMhG8AuwQkVqiCeAW4LaZlao6BpTOPBeR54C/tCSwOiYDYd7sHkveiWNic/5WH3+ADV1PIygDFVfTseMjjGx8W9wKv1mxN5NILhfkZrnJz46e8KM/0ce5We6FdxAnjiUCVQ2JyJ3AU0SHj+5T1cMichdwSFUPOPXeZn7941PUJ2mHsCs0SXn741Q3/oCi0aMEs0toP+MTdG6/jamCqri/nxV7M04TgbwsN3nZf3jCz81yIQ5Us10qR/sIVPUgcPCUZV+ZY9urnIzFREtGH+/30jGcfCWjc32dVDX+iE0tPyY7MMpEyU7q6/6W3pp3L2rOX2MSSQRyPG7yc2In+SwPebEmnVyPG5cr8Sf7+didxRliMhDmja4xxieTqClIlbUDL1N9/H42dD+LIgxUXkvHjo8wWlrnSN1/Y1YiJ8t1UvPNiW/5Wcl/sp+PJYIMkGxNQRIJUtZxkM3H9lE0eoRAzlpad+2hc9utTOfb3bwmsTxuoSDHQ16Wm4IcDwWzTvjuFD7Zz8cSQRpLtqYgT2CcyuZ/ofr4P5M72Ye3eBv1dX9H7+b3EHHnJDo8k0FmOmkLsj0U5LjJy46e8POzPWR7Mm8OCksEaSqZmoJyfZ3UNNzPppaH8YR8DG+8lCN1X2eo/Aqb+MU4KtqUExuCme050Yafl+VOik7aZGGJIA31T0xR3534pqDiodeoadhHWedTKC76am6g/YxPMrH2rITGZdKL2y3kx5pxTrTf57jJz3LjcdsXjcWwRJBGkqIpSCNs6H6WmmP7WDt4iGBWEW1nfJKOHR9lOr88cXGZlDYzBDM/5/fDLwuyoyNzVnO8fbqyRJAmEt0U5ApNUtH6KDUN91HgbWUyv5Jj5/813bU3n7b4m5V/Nqcz01E7uymnIDvacZvKo3KSnSWCNJDIpqDsqUGqGn9IVeMPyQ6MMrbuLbz+tv/DQOV18079aOWfM1tuVnTMfeEpJ/0cj327TwRLBCksElEaB7y0D61+U1DBWCM1DfsobzuAKxJkYNM7ad/5yUWP/7fyz+nP5YL8bM9J3+xn/k3XYZipyhJBikpIU5Aqa/t/zeZj+yjtfYGwO5fu2pvpOONj+Itql7w7K/+cHrI8rhNDLwtyfv+vjcxJHZYIUtBqNwVJOBC9Aazh+xSNHmE6Zz1NZ3+ezu23EsxZt+z9Wvnn1CEyUxwtdpNVTmaPu083lghSyGo3BUVvAHuI6uP3x/0GMCv/nJxmmnNOtN3PGpJpzTnpyxJBiljNpqDoDWD/zKaWHzt2A5iVf04st0tO/nYf67i15pzMZIkgBaxWU1Dx0GtsPraPjV3O3wBm5Z9Xh9v1++GYhbGTfmFO8pQ/NsnBEkESW5WmoEiYDT0zN4C9ajeApSi3W04055x8wrfhmGZhlgiS1FQwzOudzjUFuUKTbGr9CTUN95HvbVvwBjCTHDwnTvie2Ak/2rxjJ3yzEpYIkpCTTUGewDhVjT+k5vh9ZE+PLPoGMLO6sjwuCmNDMWe+4a/29IUmc9hffhJxsikoa2qImob7qG76IZ6gl8HyK2g98zM2AUyCzXzDn2nKmXlsQzLNarJEkCSmgtFRQfGeTD7H183mY/+XypZ/wRUO8FrxlRx9y6co3f7WE9tYjR/nuVzEat97KMq1NnyTXBxNBCKyG/gm0cnr71XVb5yy/i+ATwEhYAD4pKq2ORlTMnKiKSh/vJktR79Hedu/AtC7+UZad32a3/hLo+P1C8etxo8DRCAv201RTtaJIZmFuTYs0yQ3xxKBiLiB7wDXAp3AKyJyQFXrZ232W6BOVf0i8h+B/w580KmYko0TTUGFI/XUHvkuGzt/RsSdTee2W2nbeTvTBZsA2FWM1fiJk7zYOPzCHDeFsRN/QbbHqmSalOPkFcFbgUZVbQYQkf3AjcCJRKCq/zZr+5eADzsYT1KJd1NQyeCr1B65h9Ke5wllFdK6aw8dZ3yMQG7pH2xrNX6WJtvjOrlJJ1ZLxyY9MenCyURQCXTMet4JXDLP9rcDT55uhYjsAfYA1NTUxCu+hBmYmOZw99jKm4JUWdf3S2qP3MPagVcI5Kyl8Zwv0Ln9Q4Sy5z65W42f03O7haJZHbcz/1rHrUl3SdFZLCIfBuqAK0+3XlX3AnsB6urqEjv/4gqoKk0DPloHfSvcUYSNXU+z5cg9FI8cZiqvjGPn/w1dWz9AxJM/70utxs/J7fiFudGTfVGuddyazOVkIugCqmc9r4otO4mIXAP8DXClqk47GE9CBcMRDnePMzix/EOUSJDy9p+y5ch3KZhoxl+4mfq6v6Nn842oO3tR+8i0Gj8et1CU66Fw1km/MMcKqBkzm5OJ4BVgh4jUEk0AtwC3zd5ARC4AvgvsVtV+B2NJKO90iNc7RvEHwst6vSs0xabWR9h89F7y/F1MlOzkjUv/N31Vu8G1tG+x6Vrjx77lG7N8jiUCVQ2JyJ3AU0SHj+5T1cMichdwSFUPAP8DKAQejg2ta1fV9zgVUyL0j09xuGec8DL6A9xBL1VND1LT8H1ypgYZXX8BRy/8CkMVV2X0TWD2Ld+Y+HK0j0BVDwIHT1n2lVmPr3Hy/RNpJf0BWdMjVB+/n+rjD5AVHGeo7HLeuPQORjdcnFEJwL7lG7M6kqKzON0EwxHe7BpjyBtY0uty/L3UNHyfyuaH8IT89FdeS8uZdzCx7i0ORZo87Fu+MYljiSDOltMfkOdtZ/PR77Gp9SegEfpq3kXrrs/gK9nuYKSJk5vlpijXE/vJsm/5xiSYJYI4Wmp/QMFYA1uOfJfyjieIiIeu2g/QtvN2pgqrF35xChCJTnt46kk/y27EMiapWCKIg6X2BxQPvcaWI/ewsfsZQp4C2s74JO1nfJxA3kaHI3WOy0W0WSfWjl+cG23isaYdY5KfJYIVWkp/QPHQ79j25rdY3/dLAtlraDr7z+jY/mFCOWucDzSOZtrzZ77hz7TnW1E1Y1KTJYIVWGx/QNHwG2w7/C1Ke54nkLOOhnO/SNe2WwhnFaxSpMuX7XGdOOkX50YraeZn26+NMenE/qKXqW88Wjo6HJm7P6BopJ6th7/Nhu5n8LmLebH2TqbPv/1EAki2eQDyst0nfdMvyvWQ47FOXGPSnSWCJYr2B3hpHZy7dHTB6DG2Hf42G7t+TjCrmMZzvsAzxTfy7V/1c0d1mF3lJHQegNmduMUzzTvWiWtMxrJEsATBcIQ3usYYnqM/oGCskdr6b1Pe8SShrEKazv4c7Ts+Tji7iG3AHVcWJmQegPwcN8W5WdGfvOg3fuvENcbMsESwSBNTQV7vHGPyNP0B+ePNbK2/m7L2Jwh78mk+809p3/kJQtklJ223GvMA5GW7TzrhF+d6rG6+MWZelggWYa7+gLyJNrbW3015++OE3bm07vo07Ts/STBn3Wn3E+95AHKz3Ced8Ivzsqx5xxizZJYI5jFXf0Cut4PaI/9EReujqCuLtjM+QdvOTxHMXT/nvlY6D0BOlutEe35xXrSZxyZMMcbEgyWCOZyuPyDX18WWI//EppafoOKic/tHaN31aQJ5Gxbc31LmAZgZsjlzwrcSDMYYJ1kiOI1T+wNy/L1sOXIPlS0PA9C17RZazvwMgbyyRe9zrnkAzqksOXHCL86LjuKxk74xZjVZIjjF7P6A7Ml+thz5LlXN+wHorr2ZljPvYDp/eeP+3S45cbKfOfnnZdtJ3xiTWJYIYlSVxn4vbUN+sqcG2Xp0L1VNDyKRED1b3kfLWf+RqYKqJe2zIMdDcZ6HkrwsSvKyrAyDMSYpWSLg9/0BE4O9bD92L9WNP8AVCdCz+SZazvpTJgtrFtxHlsd14oRvI3iMMakk4xPBxFSQ+qZWNr7xPc5tfAB3aJLeze+m5azP4i+qPe1rXC5iQzazTpz8rYnHGJOqMjoR9PX1MvnCt7jo2H24Q376qm+g+ezP4i8+eUKYvGx37Jt+9KRflOvBZXfmGmPShKOJQER2A98kOnn9var6jVPW5wD3AxcBQ8AHVbU1njHc83wT51aVcNm20hPLXq5vYeqXd3PZwENkBSfoq9pN89l34is5A7dbWJd38rd9G69vjElnjiUCEXED3wGuBTqBV0TkgKrWz9rsdmBEVbeLyC3A3wMfjGcc51aVcOePfsvdt13AZVU5tBz8B3a99l1K8NFfeQ29F3yBrMpzqYmd9Auy3daha4zJKE5eEbwVaFTVZgAR2Q/cCMxOBDcCX4s9/jFwt4iIqi5ursdFuGxbKXffdgH7f3gv57v+kdrQGD1lVyLX/VfW19ax0Zp4jDEZzslEUAl0zHreCVwy1zaqGhKRMWA9MDh7IxHZA+wBqKlZeATPqS7bVsrR8y7ipUO19J7/eW57//uWvA9jjElXKdH4rap7VbVOVes2bFi4nMOpXmwa5O7X4Xfv2Mv/rC/ixabBhV9kjDEZwslE0AVUz3peFVt22m1ExAOUEO00jpsXmwZP9BH8xXU7ufu2C7jzR7+1ZGCMMTFOJoJXgB0iUisi2cAtwIFTtjkAfCz2+Gbg2Xj2DwC83jkW7SiOjRqa6TN4vXMsnm9jjDEpy7E+glib/53AU0SHj+5T1cMichdwSFUPAP8XeEBEGoFhoskiru64ctsfLLtsW+lJw0mNMSaTOXofgaoeBA6esuwrsx5PAR9wMgZjjDHzS4nOYmOMMc6xRGCMMRnOEoExxmQ4SwTGGJPhJM6jNR0nIgNA2zJfXsopdy2nMDuW5JQux5IuxwF2LDM2q+pp78hNuUSwEiJySFXrEh1HPNixJKd0OZZ0OQ6wY1kMaxoyxpgMZ4nAGGMyXKYlgr2JDiCO7FiSU7ocS7ocB9ixLCij+giMMcb8oUy7IjDGGHMKSwTGGJPh0iIRiMg+EekXkTfnWC8i8i0RaRSR10XkwlnrPiYix2M/Hzvd61fTCo8lLCK/i/2cWvJ71S3iWHaJyK9FZFpE/vKUdbtF5FjsOL+0OhHPbYXH0ioib8Q+l0OrE/HpLeI4PhT7vXpDRF4UkfNmrUu1z2S+Y0mazyQWz0LHcmPsWH4nIodE5PJZ61Z+DlPVlP8BrgAuBN6cY/0NwJOAAJcCL8eWrwOaY/+ujT1em4rHElvnTfRnscRj2QhcDPwd8JezlruBJmArkA28BpyViscSW9cKlCb681jkcVw28zcAXD/rbyUVP5PTHkuyfSaLPJZCft+ney5wNPY4LuewtLgiUNUXiM5nMJcbgfs16iVgjYhUAP8BeFpVh1V1BHga2O18xHNbwbEknYWORVX7VfUVIHjKqrcCjararKoBYD/R406YFRxLUlnEcbwY+1sAeInozIKQmp/JXMeSdBZxLF6NnfmBAmDmcVzOYWmRCBahEuiY9bwztmyu5clsvphzY5eNL4nITaseWfyk4ucyHwV+LiKvisieRAezBLcTvfqE1P9MZh8LpOBnIiLvFZGjwBPAJ2OL4/K5ODoxjVl1m1W1S0S2As+KyBuq2pTooAyXxz6XjcDTInI09g0waYnI1URPnpcvtG2ym+NYUu4zUdVHgUdF5Arg68A18dp3plwRdAHVs55XxZbNtTyZzRmzqs782ww8B1yw2sHFSSp+LnOa9bn0A48SbWZJWiJyLnAvcKOqDsUWp+RnMsexpNxnMlssYW0VkVLi9LlkSiI4AHw0NuLmUmBMVXuIzqd8nYisFZG1wHWxZcnstMcSO4YcgNgvyNuB+kQGugKvADtEpFZEsonOZZ3wUVDLISIFIlI085jo79hpR4YkAxGpAX4CfERVG2atSrnPZK5jSbXPBEBEtouIxB5fCOQAQ8TpHJYWTUMi8iBwFVAqIp3AV4EsAFW9h+i8yTcAjYAf+ERs3bCIfJ3oLznAXao6X0et45Z7LMCZwHdFJEI0wX9DVROaCBY6FhEpBw4BxUBERP6c6EiUcRG5k+gvtBvYp6qHE3AIJyz3WIiWDX409jfsAX6kqj9b9QOIWcTv11eA9cA/xmIOqWqdqoZS7TNhjmMBykiizwQWdSzvJ/oFMAhMAh+MdR7H5RxmJSaMMSbDZUrTkDHGmDlYIjDGmAxnicAYYzKcJQJjjMlwlgiMMSbDWSIwJkZE/kZEDs+q8niJiPy5iOQvY19/vYht7hORm5cXrTHxY4nAGEBE3ga8C7hQVc8levt+B/DnwGkTgYi459nlgonAmGRhicCYqApgUFWnAVR1ELgZ2AT8m4j8G4CIeEXkf4nIa8DbROTDIvLvsSuI74qIW0S+AeTFlv0w9rqPxq40XhORB2a97xUSrZXfbFcHJlHshjJjABEpBH5J9Nv/L4CHVPV5EWkF6mKJARFRond1/ouInAn8d+B9qhoUkX8EXlLV+0XEq6qFsdecTbSezWWqOigi62J3td9HtKTwB4FdwAFV3b6qB24MaVJiwpiVUlWviFwEvAO4GnhITj8LVxh4JPb4j4CLgFdi5QrygP7TvOadwMMzyeSUEgCPqWoEqBeRsrgcjDFLZInAmBhVDROt2vqciLwBnG7av6nYdhCdJe6fVfW/rOBtp2c9lhXsx5hlsz4CYwAR2SkiO2YtOh9oAyaAojle9gxwc6ymPSKyTkQ2x9YFRSQr9vhZ4AMisn5mu3jHb8xK2BWBMVGFwLdFZA0QIlrddQ9wK/AzEelW1atnv0BV60Xky0RnunIRnabys0QTyF7gdRH5jap+SET+DnheRMLAb4GPr9JxGbMg6yw2xpgMZ01DxhiT4SwRGGNMhrNEYIwxGc4SgTHGZDhLBMYYk+EsERhjTIazRGCMMRnu/wNy1bj4PsAbOwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEHCAYAAABBW1qbAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAAsTAAALEwEAmpwYAAAVKklEQVR4nO3dfbRldX3f8fdHBnmIJoAzUspA7pCglCQ+TC+ELB+KECOCcUhiLUSTqaWZNsFUa6wMJEvMH6yFaSJq25CMQgBDQMQHqJgmQIisdlVweBCGJ5nAIDMCcw1FNLrAkW//OHu2x/HemXMv95x9mPN+rXXW3fu39zn7y29x72f2b++zf6kqJEkCeF7XBUiSxoehIElqGQqSpJahIElqGQqSpNaSrgt4NpYuXVpTU1NdlyFJzym33HLLN6pq2WzbntOhMDU1xfr167suQ5KeU5I8NNc2h48kSS1DQZLUMhQkSS1DQZLUMhQkSS1DQZLUMhQkSS1DQZLUMhQkSa3n9DeaNT9Ta6/p7Nibzj2ps2NLGpxnCpKklqEgSWoZCpKklqEgSWoZCpKklqEgSWoZCpKklqEgSWoNLRSSXJhka5INO7T/bpJ7k9yV5I/62s9MsjHJfUneMKy6JElzG+Y3mi8C/jtwyfaGJK8DVgEvr6qnkry4aT8SOAX4GeCfA9cleUlVfX+I9UmSdjC0M4WquhF4fIfm3wbOraqnmn22Nu2rgMur6qmqehDYCBw9rNokSbMb9TWFlwCvSXJTki8mOappPxh4uG+/zU3bj0iyJsn6JOtnZmaGXK4kTZZRh8IS4ADgGOC/AFckyXw+oKrWVdV0VU0vW7ZsGDVK0sQadShsBj5TPTcDzwBLgS3AIX37LW/aJEkjNOpQ+BzwOoAkLwGeD3wDuBo4JcleSVYAhwM3j7g2SZp4Q7v7KMllwLHA0iSbgbOBC4ELm9tUnwZWV1UBdyW5Argb2Aac7p1HkjR6QwuFqjp1jk1vn2P/c4BzhlWPJGnXnHlNI9HVrG/O+CbNj4+5kCS1DAVJUstQkCS1DAVJUstQkCS1DAVJUstQkCS1DAVJUstQkCS1DAVJUstQkCS1DAVJUstQkCS1DAVJUmtooZDkwiRbmwl1dtz2e0kqydJmPUk+mmRjkjuSrBxWXZKkuQ3zTOEi4IQdG5McAvwS8LW+5jfSm4LzcGANcP4Q65IkzWFooVBVNwKPz7LpPOB9QPW1rQIuqZ4vAfslOWhYtUmSZjfSmdeSrAK2VNVXkvRvOhh4uG99c9P2yCyfsYbe2QSHHnro8Iodoq5mIZOkXRnZheYk+wJnAe9/Np9TVeuqarqqppctW7Y4xUmSgNGeKfwUsALYfpawHLg1ydHAFuCQvn2XN22SpBEa2ZlCVd1ZVS+uqqmqmqI3RLSyqh4FrgZ+s7kL6Rjgm1X1I0NHkqThGuYtqZcB/xd4aZLNSU7bye5fAB4ANgIfA35nWHVJkuY2tOGjqjp1F9un+pYLOH1YtUiSBuM3miVJLUNBktQyFCRJLUNBktQyFCRJLUNBktQyFCRJLUNBktQyFCRJLUNBktQyFCRJLUNBktQyFCRJLUNBktQyFCRJrWFOsnNhkq1JNvS1/dck9ya5I8lnk+zXt+3MJBuT3JfkDcOqS5I0t2GeKVwEnLBD27XAz1bVy4CvAmcCJDkSOAX4meY9f5pkjyHWJkmaxdBCoapuBB7foe1vq2pbs/olYHmzvAq4vKqeqqoH6U3LefSwapMkza7Lawr/DvjrZvlg4OG+bZubNknSCHUSCkl+H9gGXLqA965Jsj7J+pmZmcUvTpIm2MhDIcm/Bd4EvK2qqmneAhzSt9vypu1HVNW6qpquqully5YNtVZJmjQjDYUkJwDvA95cVd/p23Q1cEqSvZKsAA4Hbh5lbZIkWDKsD05yGXAssDTJZuBsencb7QVcmwTgS1X1H6vqriRXAHfTG1Y6vaq+P6zaJEmzG1ooVNWpszRfsJP9zwHOGVY9kqRd8xvNkqSWoSBJahkKkqTW0K4pSONgau01nR1707kndXZsaaE8U5AktQwFSVLLUJAktQwFSVLLUJAktQYKhSQ/N+xCJEndG/RM4U+T3Jzkd5L8xFArkiR1ZqBQqKrXAG+j93jrW5L8VZLXD7UySdLIDXxNoaruB/4AOAP4V8BHk9yb5FeHVZwkabQGvabwsiTnAfcAxwG/XFX/olk+b4j1SZJGaNDHXPw34OPAWVX13e2NVfX1JH8wlMokSSM3aCicBHx3+8Q3SZ4H7F1V36mqTwytOknSSA16TeE6YJ++9X2btjkluTDJ1iQb+toOSHJtkvubn/s37Uny0SQbk9yRZOV8/0MkSc/eoKGwd1V9e/tKs7zvLt5zEXDCDm1rgeur6nDg+mYd4I305mU+HFgDnD9gXZKkRTRoKPxT/7/ek/xL4Ls72Z+quhF4fIfmVcDFzfLFwMl97ZdUz5eA/ZIcNGBtkqRFMug1hXcDn0rydSDAPwP+zQKOd2BVPdIsPwoc2CwfDDzct9/mpu0RdpBkDb2zCQ499NAFlCBJmstAoVBVX05yBPDSpum+qvreszlwVVWSWsD71gHrAKanp+f9fknS3OYz89pRwFTznpVJqKpL5nm8x5IcVFWPNMNDW5v2LfS+Lb3d8qZNkjRCg3557RPAHwOvphcORwHTCzje1cDqZnk1cFVf+282dyEdA3yzb5hJkjQig54pTANHVtXAwzVJLgOOBZYm2QycDZwLXJHkNOAh4K3N7l8ATgQ2At8B3jHocSRJi2fQUNhA7+LywP96r6pT59h0/Cz7FnD6oJ8tSRqOQUNhKXB3kpuBp7Y3VtWbh1LVCEytvabrEiRp7AwaCh8YZhGSpPEw6C2pX0zyk8DhVXVdkn2BPYZbmiRp1Aa9++i3gCuBP2+aDgY+N6SaJEkdGfQxF6cDrwKehHbCnRcPqyhJUjcGDYWnqurp7StJlgB+m1iSdjODhsIXk5wF7NPMzfwp4H8OryxJUhcGDYW1wAxwJ/Af6H3ZzBnXJGk3M+jdR88AH2tekqTd1EChkORBZrmGUFWHLXpFkqTOzOfZR9vtDfxr4IDFL0eS1KWBrilU1T/2vbZU1YeBk4ZbmiRp1AYdPlrZt/o8emcO85mLQZL0HDDoH/Y/6VveBmziB4+9liTtJga9++h1wy5EktS9QYeP3rOz7VX1ofkcNMl/Bv49vTua7qQ3qc5BwOXAi4BbgN/o/xa1JGn45nP30VH0ps0E+GXgZuD++R4wycHAf6I3k9t3k1wBnEJv5rXzquryJH8GnAacP9/Pl8ZFV3N2bDrXe0C0cIOGwnJgZVV9CyDJB4Brqurtz+K4+yT5HrAvvRndjgN+vdl+Mb05HAwFSRqhQR9zcSDQP5TzdNM2b1W1Bfhj4Gv0wuCb9IaLnqiqbc1um+k9nvtHJFmTZH2S9TMzMwspQZI0h0HPFC4Bbk7y2Wb9ZHr/mp+3JPsDq4AVwBP0Hq53wqDvr6p1wDqA6elpn9QqSYto0LuPzkny18BrmqZ3VNVtCzzmLwIPVtUMQJLP0JurYb8kS5qzheXAlgV+viRpgQYdPoLe2P+TVfURYHOSFQs85teAY5LsmyTA8cDdwA3AW5p9VgNXLfDzJUkLNOh0nGcDZwBnNk17An+5kANW1U30pva8ld7tqM+jNxx0BvCeJBvp3ZZ6wUI+X5K0cINeU/gV4JX0/pBTVV9P8sKFHrSqzgbO3qH5AeDohX6mJOnZG3T46OmqKprHZyf5seGVJEnqyqChcEWSP6d3Mfi3gOtwwh1J2u3scviouRj8SeAI4EngpcD7q+raIdcmSRqxXYZCVVWSL1TVzwEGgSTtxgYdPro1yVFDrUSS1LlB7z76eeDtSTYB/wSE3knEy4ZVmCRp9HYaCkkOraqvAW8YUT2SpA7t6kzhc/SejvpQkk9X1a+NoCZJUkd2dU0hfcuHDbMQSVL3dhUKNceyJGk3tKvho5cneZLeGcM+zTL84ELzjw+1OknSSO00FKpqj1EVIknq3nwenS1J2s0ZCpKklqEgSWoZCpKkViehkGS/JFcmuTfJPUl+IckBSa5Ncn/zc/8uapOkSdbVmcJHgP9VVUcALwfuAdYC11fV4cD1zbokaYRGHgpJfgJ4Lc0czFX1dFU9AawCLm52uxg4edS1SdKk6+JMYQUwA/xFktuSfLyZ3vPAqnqk2edR4MDZ3pxkTZL1SdbPzMyMqGRJmgxdhMISYCVwflW9kt6juH9oqKh/PugdVdW6qpquqully5YNvVhJmiRdhMJmYHNV3dSsX0kvJB5LchBA83NrB7VJ0kQbeShU1aPAw0le2jQdD9wNXA2sbtpWA1eNujZJmnSDzry22H4XuDTJ84EHgHfQC6grkpwGPAS8taPaJGlidRIKVXU7MD3LpuNHXIokqY/faJYktQwFSVLLUJAktQwFSVLLUJAktQwFSVLLUJAktQwFSVLLUJAktQwFSVLLUJAktQwFSVLLUJAktQwFSVLLUJAktToLhSR7JLktyeeb9RVJbkqyMcknmwl4JEkj1OWZwruAe/rWPwicV1U/Dfw/4LROqpKkCdZJKCRZDpwEfLxZD3AccGWzy8XAyV3UJkmTrKszhQ8D7wOeadZfBDxRVdua9c3AwbO9McmaJOuTrJ+ZmRl6oZI0SUYeCkneBGytqlsW8v6qWldV01U1vWzZskWuTpIm25IOjvkq4M1JTgT2Bn4c+AiwX5IlzdnCcmBLB7VJ0kQbeShU1ZnAmQBJjgXeW1VvS/Ip4C3A5cBq4KpR1ybtDqbWXtPJcTede1Inx9XiGqfvKZwBvCfJRnrXGC7ouB5JmjhdDB+1qurvgb9vlh8Aju6yHkmadON0piBJ6pihIElqGQqSpJahIElqGQqSpJahIElqGQqSpJahIElqGQqSpJahIElqGQqSpJahIElqGQqSpJahIElqGQqSpFYXczQfkuSGJHcnuSvJu5r2A5Jcm+T+5uf+o65NkiZdF2cK24Dfq6ojgWOA05McCawFrq+qw4Hrm3VJ0giNPBSq6pGqurVZ/hZwD3AwsAq4uNntYuDkUdcmSZOu02sKSaaAVwI3AQdW1SPNpkeBA7uqS5ImVWehkOQFwKeBd1fVk/3bqqqAmuN9a5KsT7J+ZmZmBJVK0uToJBSS7EkvEC6tqs80zY8lOajZfhCwdbb3VtW6qpquqully5aNpmBJmhBd3H0U4ALgnqr6UN+mq4HVzfJq4KpR1yZJk25JB8d8FfAbwJ1Jbm/azgLOBa5IchrwEPDWDmqTpIk28lCoqv8NZI7Nx4+yFknSD/MbzZKklqEgSWoZCpKklqEgSWoZCpKklqEgSWoZCpKklqEgSWoZCpKkVhePuZC0G5pae01nx9507kmdHXt345mCJKllKEiSWg4fSXrO62roancctjIUJGmBdsfrKA4fSZJahoIkqTV2oZDkhCT3JdmYZG3X9UjSJBmrUEiyB/A/gDcCRwKnJjmy26okaXKMVSgARwMbq+qBqnoauBxY1XFNkjQxxu3uo4OBh/vWNwM/379DkjXAmmb120nuW+CxlgLfWOB7h21ca7Ou+RvX2qxr/saqtnywXVxIXT8514ZxC4Vdqqp1wLpn+zlJ1lfV9CKUtOjGtTbrmr9xrc265m9ca1vsusZt+GgLcEjf+vKmTZI0AuMWCl8GDk+yIsnzgVOAqzuuSZImxlgNH1XVtiTvBP4G2AO4sKruGtLhnvUQ1BCNa23WNX/jWpt1zd+41raodaWqFvPzJEnPYeM2fCRJ6pChIElqTWQojOujNJJsSnJnktuTrO+4lguTbE2yoa/tgCTXJrm/+bn/mNT1gSRbmn67PcmJHdR1SJIbktyd5K4k72raO+2zndQ1Dn22d5Kbk3ylqe0Pm/YVSW5qfj8/2dx0Mg51XZTkwb4+e8Uo6+qrb48ktyX5fLO+uP1VVRP1oncB+x+Aw4DnA18Bjuy6rqa2TcDSrutoanktsBLY0Nf2R8DaZnkt8MExqesDwHs77q+DgJXN8guBr9J7VEunfbaTusahzwK8oFneE7gJOAa4Ajilaf8z4LfHpK6LgLd02WdNTe8B/gr4fLO+qP01iWcKPkpjAFV1I/D4Ds2rgIub5YuBk0dZE8xZV+eq6pGqurVZ/hZwD71v6HfaZzupq3PV8+1mdc/mVcBxwJVNexd9NlddnUuyHDgJ+HizHha5vyYxFGZ7lMZY/JLQ+x/vb5Pc0jzOY9wcWFWPNMuPAgd2WcwO3pnkjmZ4aeTDWv2STAGvpPcvzLHpsx3qgjHos2Yo5HZgK3AtvbP4J6pqW7NLJ7+fO9ZVVdv77Jymz85Lsteo6wI+DLwPeKZZfxGL3F+TGArj7NVVtZLeU2JPT/LarguaS/XOVcfiX0/A+cBPAa8AHgH+pKtCkrwA+DTw7qp6sn9bl302S11j0WdV9f2qegW9pxccDRzRRR072rGuJD8LnEmvvqOAA4AzRllTkjcBW6vqlmEeZxJDYWwfpVFVW5qfW4HP0vslGSePJTkIoPm5teN6AKiqx5pf4meAj9FRvyXZk94f3kur6jNNc+d9Nltd49Jn21XVE8ANwC8A+yXZ/sXaTn8/++o6oRmKq6p6CvgLRt9nrwLenGQTvWHv44CPsMj9NYmhMJaP0kjyY0leuH0Z+CVgw87fNXJXA6ub5dXAVR3W0tr+R7fxK3TQb83Y7gXAPVX1ob5NnfbZXHWNSZ8tS7Jfs7wP8Hp61zxuAN7S7NZFn81W17194R564/Yj7bOqOrOqllfVFL2/W39XVW9jsfur6yvpXbyAE+ndhfEPwO93XU9T02H07oT6CnBX13UBl9EbVvgevXHK0+iNX14P3A9cBxwwJnV9ArgTuIPeH+GDOqjr1fSGhu4Abm9eJ3bdZzupaxz67GXAbU0NG4D3N+2HATcDG4FPAXuNSV1/1/TZBuAvae5Q6uIFHMsP7j5a1P7yMReSpNYkDh9JkuZgKEiSWoaCJKllKEiSWoaCJKllKEiSWoaC1JGM6SPcNdn8noLUgSR70PsC5evpfQnvy8CpVXV3p4Vp4nmmIA0gyVSSe5uJVr6a5NIkv5jk/zQT6BzdTFzz3r73bGieTDobH+GusWQoSIP7aXpPEz2ief06vcdIvBc4a56fNc6PcNcEMxSkwT1YVXdW78midwHXV2/89U5gqtPKpEViKEiDe6pv+Zm+9WeAJcA2fvh3au+dfNbYPsJdk81QkBbPJnrzR5NkJbBiJ/uO5SPcJUNBWjyfBg5IchfwTnp3F82qetMnvhP4G3pzCFxRVXeNpEppJ7wlVZLU8kxBktRasutdJC1Uku0zr+3o+Kr6x1HXI+2Kw0eSpJbDR5KklqEgSWoZCpKklqEgSWr9f5NDpj0BJ4lzAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "mcmc_samples = np.array(mcmc.get_samples()) #convert the list into an array\n",
    "mcmc_probs = np.array(mcmc.get_probs()) #convert the list into an array\n",
    "mcmc_force = np.array(mcmc.get_values()) #convert the list into an array\n",
    "print(np.array(mcmc_force).shape)\n",
    "#We can plot the mean and variation of these\n",
    "mean = np.mean(mcmc_force,axis=0)\n",
    "sd = np.std(mcmc_force,axis=0)\n",
    "plt.plot(applied_stretch,force_measured,'x')\n",
    "plt.plot(applied_stretch,mean,'-')\n",
    "plt.fill_between(applied_stretch,mean-sd,mean+sd,alpha=0.3)\n",
    "plt.xlabel('Stretch')\n",
    "plt.ylabel('Force')\n",
    "plt.show()\n",
    "\n",
    "plt.hist(mcmc_samples)\n",
    "plt.xlabel('mu_0')\n",
    "plt.ylabel('Frequency')\n",
    "plt.show()"
   ]
  }
 ],
 "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.9.16"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}