Skip to content

Instantly share code, notes, and snippets.

@ljmartin
Created March 14, 2025 05:38
Show Gist options
  • Select an option

  • Save ljmartin/3f10feae41923f71e0e0437a573a2fd6 to your computer and use it in GitHub Desktop.

Select an option

Save ljmartin/3f10feae41923f71e0e0437a573a2fd6 to your computer and use it in GitHub Desktop.
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"id": "865a8508",
"metadata": {},
"outputs": [],
"source": [
"from rdkit import Chem\n",
"from rdkit.Chem import rdDistGeom\n",
"\n",
"import os\n",
"\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt"
]
},
{
"cell_type": "markdown",
"id": "93d22223",
"metadata": {},
"source": [
"# use rdkit example where we can induce an intra-mol hbond\n",
"\n",
"example comes from: https://github.com/rdkit/rdkit/discussions/7658\n",
"\n",
"first up, random conformers. as in Greg's example, the =O and H are far apart in several conformers. "
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "e6b1a92a",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"3.5509453080751796\n",
"3.267209134144379\n",
"2.4030999716313475\n",
"4.4486523906077915\n",
"3.4748769955217416\n",
"3.7180844765319754\n",
"2.3374946203630773\n",
"4.845241756813917\n",
"3.790555985922051\n",
"2.507302462526262\n"
]
}
],
"source": [
"m1 = Chem.AddHs(Chem.MolFromSmiles('O=CCCO'))\n",
"ps = rdDistGeom.ETKDGv3()\n",
"\n",
"ps.randomSeed = 0xf00d\n",
"rdDistGeom.EmbedMultipleConfs(m1,10,ps)\n",
"\n",
"for conf in m1.GetConformers(): print((conf.GetAtomPosition(0)-conf.GetAtomPosition(10)).Length())"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "a213fa2e",
"metadata": {},
"outputs": [],
"source": [
"coords = []\n",
"# save the coordinates for later.\n",
"for confid in range(10):\n",
" coords.append(m1.GetConformer(confid).GetPositions())"
]
},
{
"cell_type": "markdown",
"id": "e1c7f7f7",
"metadata": {},
"source": [
"save these 10 conformers to directories `./inps/0` through `./inps/10`:"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "2dfe66c8",
"metadata": {},
"outputs": [],
"source": [
"os.mkdir('./inps')\n",
"\n",
"for i in range(m1.GetNumConformers()):\n",
" if not os.path.exists(f'./inps/{i}/',):\n",
" os.mkdir(f'./inps/{i}/')\n",
" Chem.MolToMolFile(m1, f'./inps/{i}/conf.sdf', confId=i)"
]
},
{
"cell_type": "markdown",
"id": "ed7cbbbd",
"metadata": {},
"source": [
"next up, restraint the =O H distance. we see that all of them show close interactions:"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "12954f3e",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"1.803187178185039\n",
"1.5405472628556918\n",
"1.7644039879418514\n",
"1.8168518320698919\n",
"1.625253266281199\n",
"1.7347552180856414\n",
"1.816385058315408\n",
"1.7517026470221735\n",
"1.8067746505159683\n",
"1.6663290159535473\n"
]
}
],
"source": [
"mbounds = rdDistGeom.GetMoleculeBoundsMatrix(m1)\n",
"\n",
"# lower bound on the O-H distance:\n",
"mbounds[10,0] = 1.6\n",
"\n",
"# upper bounds on the O-H distance:\n",
"mbounds[0,10] = 1.8\n",
"\n",
"ps.SetBoundsMat(mbounds)\n",
"\n",
"rdDistGeom.EmbedMultipleConfs(m1,10,ps)\n",
"\n",
"for conf in m1.GetConformers(): print((conf.GetAtomPosition(0)-conf.GetAtomPosition(10)).Length())"
]
},
{
"cell_type": "code",
"execution_count": 6,
"id": "83925f0b",
"metadata": {},
"outputs": [],
"source": [
"# also save the coordinates for later.\n",
"for confid in range(10):\n",
" coords.append(m1.GetConformer(confid).GetPositions())"
]
},
{
"cell_type": "markdown",
"id": "9b7c6a3e",
"metadata": {},
"source": [
"save them in directories `./inps/10` through `./inps/19`:"
]
},
{
"cell_type": "code",
"execution_count": 7,
"id": "a037325e",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"10\n",
"11\n",
"12\n",
"13\n",
"14\n",
"15\n",
"16\n",
"17\n",
"18\n",
"19\n"
]
}
],
"source": [
"for i in range(m1.GetNumConformers()):\n",
" j = i+10\n",
" if not os.path.exists(f'./inps/{j}/',):\n",
" os.mkdir(f'./inps/{j}/')\n",
"\n",
" Chem.MolToMolFile(m1, f'./inps/{j}/conf.sdf', confId=i)\n",
" print(j)"
]
},
{
"cell_type": "markdown",
"id": "becd9ac9",
"metadata": {},
"source": [
"now write a bash script that uses antechamber to calculate am1-bcc charges for each conformer, in each directory:"
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "fa977016",
"metadata": {},
"outputs": [],
"source": [
"formal_charge = Chem.GetFormalCharge(m1)\n",
"with open('./inputs.dat', 'w') as f:\n",
" for i in range(20):\n",
" f.write(f'cd ./inps/{i}')\n",
" f.write('\\n')\n",
" \n",
" s = f\"antechamber -i conf.sdf -fi sdf -o conf.mol2 -fo mol2 -pf yes -dr n -c bcc -nc {formal_charge}\"\n",
" f.write(s)\n",
" f.write('\\n')\n",
" \n",
" s = f\"antechamber -dr n -i conf.mol2 -fi mol2 -o charges.mol2 -fo mol2 -c wc -cf charges.txt -pf yes\"\n",
" f.write(s)\n",
" f.write('\\n')\n",
" \n",
" f.write('cd ../../')\n",
" f.write('\\n')"
]
},
{
"cell_type": "code",
"execution_count": 9,
"id": "c5525784",
"metadata": {},
"outputs": [],
"source": [
"!bash inputs.dat > log.dat"
]
},
{
"cell_type": "markdown",
"id": "87c0045e",
"metadata": {},
"source": [
"load up all the charges and plot:"
]
},
{
"cell_type": "code",
"execution_count": 10,
"id": "c0367a3d",
"metadata": {},
"outputs": [],
"source": [
"charges = []\n",
"for i in range(20):\n",
" fname = f'./inps/{i}/charges.txt'\n",
" chrgs = open(fname).read()\n",
" chrgs = np.array([float(i.removesuffix('\\n')) for i in chrgs.strip().split()])\n",
" charges.append(chrgs)"
]
},
{
"cell_type": "code",
"execution_count": 11,
"id": "a045ea0f",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"20\n"
]
}
],
"source": [
"print(len(charges))"
]
},
{
"cell_type": "code",
"execution_count": 12,
"id": "86e66570",
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAi8AAAGdCAYAAADaPpOnAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjAsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvlHJYcgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAbZ9JREFUeJzt3Xl4W+WZ///3kWRbXuV9l7fEWSAbSUjiJEDZQth3QtMGWihfGKBs007LMPMrw8yUlmlppkPZCxQIEChLoU0DKRASyL7vduI43mVLXiTvtqTz++PYhhAnsRNJR5Lv13X5urB8dM7HTrDuHD33cyuqqqoIIYQQQoQIg94BhBBCCCFGQooXIYQQQoQUKV6EEEIIEVKkeBFCCCFESJHiRQghhBAhRYoXIYQQQoQUKV6EEEIIEVKkeBFCCCFESDHpHcDXvF4vdXV1xMfHoyiK3nGEEEIIMQyqqtLW1kZ2djYGw4nvrYRd8VJXV4fVatU7hhBCCCFOQXV1Nbm5uSc8JuyKl/j4eED75hMSEnROI4QQQojhcLlcWK3WwdfxEwm74mXgraKEhAQpXoQQQogQM5wlHwFZsPv0009TWFiI2WxmxowZrF279oTH9/T08Mgjj5Cfn09UVBRjxozhpZdeCkRUIYQQQgQ5v995Wb58OQ888ABPP/008+bN47nnnuPSSy9l37595OXlDfmcm266iYaGBv74xz8yduxYGhsbcbvd/o4qhBBCiBCgqKqq+vMCs2fPZvr06TzzzDODj02cOJFrrrmGxx9//JjjV65cyc0338zhw4dJTk4e8fVcLhcWiwWn0ylvGwkhhBAhYiSv335926i3t5etW7eyYMGCox5fsGAB69atG/I5H374ITNnzuSJJ54gJyeHcePG8ZOf/ISuri5/RhVCCCFEiPDr20YOhwOPx0NGRsZRj2dkZGCz2YZ8zuHDh/nyyy8xm828//77OBwO7r77bpqbm4dc99LT00NPT8/g5y6Xy7ffhBBCCCGCSkAW7H575bCqqsddTez1elEUhWXLljFr1iwuu+wynnzySV555ZUh7748/vjjWCyWwQ/Z40UIIYQIb34tXlJTUzEajcfcZWlsbDzmbsyArKwscnJysFgsg49NnDgRVVWpqak55viHH34Yp9M5+FFdXe3bb0IIIYQQQcWvxUtkZCQzZsxg1apVRz2+atUq5s6dO+Rz5s2bR11dHe3t7YOPlZWVYTAYhtxxLyoqanBPF9nbRQghhAh/fn/b6KGHHuLFF1/kpZdeYv/+/Tz44INUVVVx1113Adqdk1tuuWXw+MWLF5OSksIPf/hD9u3bx5o1a/jpT3/KbbfdRnR0tL/jCiGEECLI+X2fl0WLFtHU1MRjjz1GfX09kyZNYsWKFeTn5wNQX19PVVXV4PFxcXGsWrWKH//4x8ycOZOUlBRuuukm/uu//svfUYUQQggRAvy+z0ugyT4vQgghROgJmn1ehBBCCCF8TYqXMPf/Pf0YD/7p13rHEEIIIXxGipcw1uZysWzC+byZdwm/+f2jescRQgghfEKKlzD23nsv06HEA1AXb9Y5jRBCCOEbUryEsSOu5sH/tiVZTnCkEEIIETqkeAljzoSv98WpjknVMYkQQgjhO1K8hLHm+LjB/64x5dDcZNcxjRBCCOEbUryEMXt0/OB/dykxvPrGszqmEUIIIXxDipcwZo9MPOrzOrP8cQshhAh98moWptpcLuyGNADG9h0CwJYsi3aFEEKEPilewtT7771MlxKDonqZ4qgEoCYuRedUQgghxOmT4iVMHXE2AZCsNpNhawGgJiKbNpdLz1hCCCHEaZPiJUy1JsQAkOZpIr/PiEH14FIsLH/7eZ2TCSGEEKdHipcw1RwfC0Bqt5OFi24h02sDoKKvS89YQgghxGmT4iVMOWK0NunUzg4yC8aQ29MAgC31xGPGhRBCiGAnxUuYskckAZDk7AAgx9UKQG1csl6RhBBCCJ+Q4iUMtblc2I3aOICE/uIl3e4EoDoyU7dcQgghhC9I8RKGPvzwVToVbc1LUVI2AOltPQA0GdJ4++2XdMsmhBBCnC4pXsJQef8Mo2RvE5ff+D0AFlxxE2neRgD2Nlbrlk0IIYQ4XVK8hKGBadJpnibiErU1LmMnTyO3V+s4apBFu0IIIUKYFC9haKBNOq3HOfiYYjCQ62oGoDY+UY9YQgghhE9I8RKGHNFxAKR0tB/1eHqztrtujTkj4JmEEEIIX5HiJQzZI/vbpF0dRz2e0qR9blMy+eLzFQHPJYQQQviCFC9hxuv1DrZJxzuP3k33kgsuJ1FtQVUMrNm1UY94QgghxGmT4iXMvPnqM3Qo2ttGxUlZR31twqy55PbWA9CQEh/wbEIIIYQvSPESZg63aGMAktRmrlh0y1FfM5pM5HZoi3brLImBjiaEEEL4hBQvYcZp6W+TdjcRE39sS3Rmk9aBVBOdHtBcQgghhK9I8RJmmhMG2qRbh/x6Uv+YgFpDNvv37QhQKiGEEMJ3pHgJM45obS3Lt9ukB8w5YxYxajsexcRfPn4vkNGEEEIIn5DiJczYIxMBSHZ1Dvn1OZdehtVdB0BDUmygYgkhhBA+I8VLGPF6vTQOtEl/a4+XAVHmGHI7mgCoS7IELJsQQgjhK1K8hJG3X3+GDkV722hsatZxj8tqaQWgJiYtELGEEEIIn5LiJYwcatIGLyaqzVxz84+Oe1yKvX9MgDGHmrqqgGQTQgghfEWKlzDiTIgBtDbpKHP0cY+blllEpNpDj2Jm+TsvBSqeEEII4RNSvISRr9uknSc8bv4115PrqQXAFhfl91xCCCGEL0nxEkbsMdpYgNTOodukB8RbksnpsgNQL4t2hRBChBgpXsKIo79N+tvTpIeS3dq/WV1ssj8jCSGEED4nxUuY0Nqkte6heOfQe7x8U1r/TrtVETm0uVx+zSaEEEL4khQvYeLt15+jvb9Nelxa7kmPHxthwai66VDiee31p/wdTwghhPAZKV7CRHmTtmuuRW3l2u/eftLjL1n0fbK99QDUGL1+zSaEEEL4khQvYaK1v0063e0gIjLypMcnZWST290IQH2KLNoVQggROqR4CRPN/cVLam/rsJ+T7dSOrY1L8kMiIYQQwj+keAkTjhhtvUvqcaZJDyXN0QpAdWS2PyIJIYQQfiHFS5gYSZv0gOwuBUX10qIk8/qrsmhXCCFEaJDiJQx4PR4ajSkAWFwnb5MecNXNt5Kuauteytqa/ZJNCCFE+Ghzubh8xZ+4Z9lvWfvlx7rlkOIlDLyz7HnaFG3R7dhU67Cfl1kwBmuPNszRlprgl2xCCCHCx9JXf8fW6Kl8mHUOtvpa3XJI8RIGDjq0v0AJaivXLz55m/Q3ZbtaAKiNl0W7QgghTmxTYR4AMzt3ceONt+mWIyDFy9NPP01hYSFms5kZM2awdu3aYT3vq6++wmQyMW3aNP8GDHGugWnSniZMEREjem6GQ9tptyYqw+e5hBBChI8/vfYU26InAzC7Qr+7LhCA4mX58uU88MADPPLII2zfvp1zzjmHSy+9lKqqqhM+z+l0csstt3DhhRf6O2LI+3qadOuIn5vq7AGgwZDJRx+95ctYQgghwsgXZi8excQYdzk///EvdM3i9+LlySef5Pbbb+dHP/oREydOZOnSpVitVp555pkTPu/OO+9k8eLFlJSU+DtiyHMMc5r0UK68ahFJ3iYAtlXs92kuIYQQ4aH0wE6+SpkCwLyqMp3T+Ll46e3tZevWrSxYsOCoxxcsWMC6deuO+7yXX36Z8vJyfvGLk1d2PT09uFyuoz5GG3tUIgDJI2iTHlB45hSsfdqYgMYUWbQrhBDiWC9+tQKnkkii2kLS53t47p5/obHSplsevxYvDocDj8dDRsbR6ykyMjKw2Yb+pg8ePMjPf/5zli1bhslkOuk1Hn/8cSwWy+CH1Tr8bptw4PV4sBtTAUgYxjTpb1MMBnLbtTbpWkuiL6MJIYQIE+vyxgEwr2k3fdlZNBk6aak7pFuegCzYVRTlqM9VVT3mMQCPx8PixYv5j//4D8aNGzescz/88MM4nc7Bj+rqap9kDhV/fv15XP1t0mPSck7pHJlN2t2qGnO6z3IJIYQID7/6v/+g3DQGo+pm7M5yvDFx9CUmk1UwRrdMJ7+1cRpSU1MxGo3H3GVpbGw85m4MQFtbG1u2bGH79u3ce++9AHi9XlRVxWQy8cknn3DBBRcc9ZyoqCiioqL8900EuYOOGsiDeNXJTYvvOKVzJDe1AVBnyGbrpjXMmHWuLyMKIYQIYRsLcwGY3rUbd68JTJDf0UpCVpZumfx65yUyMpIZM2awatWqox5ftWoVc+fOPeb4hIQEdu/ezY4dOwY/7rrrLsaPH8+OHTuYPXu2P+OGJGd/p1G6x4FxGG+zDeW8sy8gTm3DqxhZuW7VyZ8ghBBiVHhr+QtsidHao88uPYwnJh68Kpf86E5dc/n1zgvAQw89xJIlS5g5cyYlJSU8//zzVFVVcddddwHa2z61tbW8+uqrGAwGJk2adNTz09PTMZvNxzwuNC2D06Sdp3yO6edfiPXTP7M/cgKNyfG+iiaEECLEfUYHfUok+e5KPFVNKInpZDlbyDrjTF1z+b14WbRoEU1NTTz22GPU19czadIkVqxYQX5+PgD19fUn3fNFHJ8jdqBNuu2Uz2E0mcjtaGJ/JNQlJvoomRBCiFBWWVXOl6naXZf5NQdQEtNBVbnghut1ThagBbt33303R44coaenh61bt3LuuV+vqXjllVdYvXr1cZ/76KOPsmPHDv+HDFH2/mnSyc6Rt0l/U2Zz/0670WmnG0kIIUQYePbjt2g2pBCvuojfVwlAitNJ8Tnn6ZxMZhuFNE9vL3aT1iad2NZ1WudKdmgdR7XGbCqOHDztbEIIIULbuvxiAEpadkGM9g/b+ecFx8axUryEsPeXvYhTSQSgKPXU2qQHzBo/FbPaRa8SxTsfvOaDdEIIIULVb576T0ojxmFQPYzbcxgUhXhXO2ddu0jvaIAULyHtQJM2GCtedfHdJXed1rnmX3o1OZ46ABoTok87mxBCiNC1MU/bzmRqz17cHm1t5czxeXpGOooULyGszTIwTdpx2ueKMsdg7bQDUJeceNrnE0IIEZo+/HAZG+O0OUazDh4CgwFzRxfn3HGvzsm+JsVLCBucJt3b6pPzZbX0L9qNSfHJ+YQQQoSev3c00KuYyfHU4KnrBmByUjQGQ/CUDMGTRIxY02lMkx5Kir2/eDHl4Gxt9ck5hRBChI7mJjtfpWt7uMyv24fJHENEdy8Lf/qIzsmOJsVLCLNHaTONTmWa9FDOSMrGpPbRqcTy8qu/98k5hRBChI7fvvMsjYYMYtR2LAcbACg29GKMiNA52dGkeAlR32yTtrhGPk16KBffuJhsr7Zo1xZt9Mk5hRBChI4NhWMBKHHuxmCMx9jr5oqf/ZvOqY4lxUuI+stbf6RVSQKgMDXbJ+eMtyRj7WoEwJZs8ck5hRBChIalT/83eyMnoqheJuw/DEB+j4sYS6K+wYYgxUuI2t9YDUCs2sb3lvyTz86b7exf9xKb7LNzCiGECH4bsrVmjUm9+/H0xqK4PSy87wF9Qx2HFC8hqs3y9TRpX0q1twJQHembuzlCCCGC36effcjGBK09ena5tst6dlsr6flFesY6LileQlRzfP8eLz5qkx6QZzCjqB6cSiIvPP8bn55bCCFEcPqgrowuJYYMrw2lrhe8Xi68dbHesY5LipcQ5YiNB05vmvRQrrn5NjJVbYV5pcc3C4GFEEIErzaXiy+zzgBgvm0PxohIUludFE0PjjlGQ5HiJUQ5Btqk23zTJj0gKSOb3G6teLGlyKJdIYQId//z2lLqDdmY1S5SyxpAVTnn0u/oHeuEpHgJQe6uLuwmbWGVxXl606SHkuNqAaA2Psnn5xZCCBFcNhZp61pmt+1CUWKxONuYeum1Oqc6MSleQtBf3nyJFkXrBhqb5vuFtRkOreOoOjLT5+cWQggRPP7w/BPsipoIwKQDFQDMmBqci3S/SYqXEHSgSdtILlZt57vfP71p0kPJ6HAD4DCk8+5bL/r8/EIIIYLD+tQ4VMXIxN4DeLrMxLg6OOfWu/WOdVJSvISgNks04Jtp0kO57ubbSPVqE6Z3Omr8cg0hhBD6Wr/+U9Ynau3RJRVae/SZmXEoQTSA8XiCP6E4RsvgNOkWv5w/s2AMub02AOwpCX65hhBCCH29XbaVDiWOVK8dQ1UHUR3dXPLgw3rHGhYpXkJQU2z/NOku30yTHkpuWzMAtQmyaFcIIcJNm8vFlznaWpf59j2YIs0UR3kwRUTqnGx4pHgJQXY/tUl/U3qTC4CaqHS/XUMIIYQ+lr76O6qNViLVHjJLazF193BpEA5gPB4pXkJMT3s7jcaBNmn/bSKX2qLd1ak3ZLLmsxV+u44QQojA21BUAMDZHbtAjSOvt4PY+NDZ20uKlxDz17f+RItBK16KUv3XynzxeVeSoLaiKka+2L3Bb9cRQggRWH98+XdsN08CYErpYQy9fVz60EM6pxoZKV5CTGmz1v0To7Zz9VXf99t1JpXMxdpXD0BjcrzfriOEECKw1sab8CpGivsOorZHkdXuIi23QO9YIyLFS4gZmCad5nEQn+C/TiDFYCC3vQmA+sREv11HCCFE4Ozes411yf3t0ZVlKG43F9z+PZ1TjZwULyGmZaB46Wv1+7Uym/t32jXLol0hhAgHf9r8CS7FQpLaTFRVO6mtLsZMnaN3rBGT4iXEOGL626R9PE16KMkOreOo1phN2YGdfr+eEEII//oqdzwA8xy7MWKi5IrzdU50aqR4CTEOs/ZWkT/bpAdcMPMCYtQO3EoE7674s9+vJ4QQwn/++6n/oMJUiEntw1paQ6KzjekLg3sA4/FI8RJCup1O7P1t0oku/7VJD5h+/oXkurU5SvbkWL9fTwghhP9sKrQCMKNrN2pfNNPOKtY50amT4iWErFj+Ks3906QLk/0/8dloMpHbqc1PqpNFu0IIEbLeeOs5tkZPBmB66WFinR2ce4vvB/sGihQvIaS0tQ5VMRCtdnLN1f5rk/6mrBZt0W5NTGpArieEEML3PjN041YiKHRXoDqNTMyJx2A06h3rlEnxEkJcCTEApHntfm2T/qYUe3/xYsrB3lAXkGsKIYTwnfLDB/gqVbvrMq/6AOa2Di6+7190TnV6pHgJIQNt0ul+miY9lJmFE4hQe+lWonn1jecCdl0hhBC+8cJn79GiJJOgOokub6EoWiHKHKN3rNMixUsICcQ06W8778obyfXUAtAQHxWw6wohhPCNr/LGATC3eRfmHrjkJz/TOdHpk+IlhAxMk04JQJv0gChzDDlddgBsyaEztEsIIQQ88dRjHIwYi1F1U1haTa67G0tS6K9hlOIlRHQ2N2M3aZ1GlgC0SX9TTmsrADWxyQG9rhBCiNOzMT8LgGndezG4FC687590TuQbUryEiBXv/IlmRdvjJd8S2Ko5xa7ttFsdkYOzv5ARQggR3D74y6tsjtXmGM08WE5mRwe5RRN1TuUbUryEiLLWBlTFgFnt5Lprbw3otSfEpWBQPbQpCSxb9oeAXlsIIcSpWdnVRK8ShdVThVrTzvxbbtQ7ks9I8RIi2gc6jbz+nSY9lMu++wOyvPUA1Bi9Ab22EEKIkauz1fBl2iQA5tUeIKOjj4mzztM5le9I8RIiBtqkUwPYJj0gLiGJ3J5GAGwpsmhXCCGC3VN/eQWHIY1YtQ3LwQZmLJyrdySfkuIlRDhiteIlrcv/06SHkuPUiqbaOFm0K4QQwW59wVgA5rbuIqm5j9lXLNI5kW9J8RIiHP1t0oGYJj2UtP6ddqsjs3S5vhBCiOF58g//xf7ICSiqh3GlRzhjahGKIbxe7sPruwlTrfV1g23SgZgmPZS8Pm0GRrMhhWWvyKJdIYQIVhtz0wCY0rOf6No+Lrjtbp0T+Z4ULyHgH++9TZNObdIDrltyB+neBgBKOxy6ZBBCCHFin6x6l43xWnv0rPJDjM1KwBQRqXMq35PiJQSUtdlQFSNmtYvrA9wmPSApI5vcHhsADSmB7XYSQggxPH9prKBbiSbLW0fk/kYW3P/PekfyCyleQkCHJfDTpIeS09a/aDchSbcMQgghhtbcZOfLTK09en79PsbGmYmJC88O0YAUL08//TSFhYWYzWZmzJjB2rVrj3vse++9x8UXX0xaWhoJCQmUlJTw8ccfByJm0Bpok07rC3yb9DdlOPoX7UZl6ppDCCHEsX739rM0GDKJVjtJ213FRQ88oHckv/F78bJ8+XIeeOABHnnkEbZv384555zDpZdeSlVV1ZDHr1mzhosvvpgVK1awdetWzj//fK688kq2b9/u76hBS49p0kNJ6+gFoFFJ55OP3tI1ixBCiKNtKCwCYI5rJ9ZOL2lZeTon8h+/Fy9PPvkkt99+Oz/60Y+YOHEiS5cuxWq18swzzwx5/NKlS/mXf/kXzj77bIqLi/nlL39JcXExH330kb+jBi27WXurKMWlb/Gy6KbbSFSbURUDGw7v1zWLEEKIrz313K/YE6nNLTpz/xHOvfM2nRP5l1+Ll97eXrZu3cqCBQuOenzBggWsW7duWOfwer20tbWRnDz05mg9PT24XK6jPsJJc00VdpPWaZTY1qVrlsyCMVh7tTEBDanxumYRQgjxtfXpFlTFwJm9+0mramPMpOl6R/IrvxYvDocDj8dDRkbGUY9nZGRgs9mGdY7f/va3dHR0cNNNNw359ccffxyLxTL4YbVaTzt3MPn0veU4+tukc4Jg4VVOezMAdZZEfYMIIYQAYO2XH7PeMhmAOYcPcvYNl+ucyP8CsmBXUZSjPldV9ZjHhvLmm2/y6KOPsnz5ctLT04c85uGHH8bpdA5+VFdX+yRzsDjU2YSqGIlUu1l0/e16xyGrSVu0W2Me+s9DCCFEYP358E46lTjSvQ2kbKti+oVSvJyW1NRUjEbjMXdZGhsbj7kb823Lly/n9ttv5+233+aiiy467nFRUVEkJCQc9RFO2vvbpNN1bpMekNKsjSeoM2SzdcMXOqcRQojRrc3l4svsMwCY17iXs86ZrXOiwPBr8RIZGcmMGTNYtWrVUY+vWrWKuXOPP+HyzTff5Ac/+AFvvPEGl18e/hXkibQmBEeb9IBLzr2cWLUdj2Ji5fp/6B1HCCFGtd++9jtqjblEqt3kbz/Ced8L74W6A/z+ttFDDz3Eiy++yEsvvcT+/ft58MEHqaqq4q677gK0t31uueWWwePffPNNbrnlFn77298yZ84cbDYbNpsNp9Pp76hByRGn7zTpb5tUMhdrXy0A9mRZtCuEEHraWFQIwOz2XUzKzMJgNOqcKDD8XrwsWrSIpUuX8thjjzFt2jTWrFnDihUryM/PB6C+vv6oPV+ee+453G4399xzD1lZWYMf999/v7+jBiWHWd9p0t+mGAzkdDYBUJek/wJiIYQYrZ578TfsjDoTgKn7D3Px3aPnddIUiIvcfffd3H330FMtX3nllaM+X716tf8DhYjGw4e+nibdps806aFkNzshEWqj0/SOIoQQo9ZXSWa8ipHxfWUUt0OUOUbvSAEjs42C2Nr338WhaFOkc6OD5y5H8sCYAGM2VRWHdE4jhBCjz/bt61mfpE2PLjlSxoL7H9Q5UWBJ8RLEynqa8CpGItUerrvm+3rHGTRv6lwi1W56FTNvvfcnveMIIcSo8/rONbQpCaR4HeTvqyUp9cQdvOFGipcg1p6kzTRK89pJTgmet2jmXXQZVk8dAPbE0XObUgghgoGztZUvcyYAMM+xh8vv/rHOiQJPipcg1pqgFQbB0iY9wGgykdtpB6AuOXjezhJCiNHg93/6HZWmfCLUXsbtPELeuDP0jhRwUrwEsa+nSQdHm/Q3Zbdo615qY1N0TiKEEKPLprEFAMzs3MUll12tbxidSPESxBwD06Tb9Z0mPZTU/kW7VaYcWpuadE4jhBCjwwtP/ZJt0docoxkHDjN5/vk6J9KHFC9BqvbAPuwR/W3SruBpkx4wNWcsRtVNpxLHq68/pXccIYQYFdZnxOFRTIxxl3Ne/pl6x9GNFC9B6quP/jLYJp0dHadzmmNdfO0isr3aot266NGxo6MQQuhp/ccf8VWK1h49t7KU+dd/V+dE+pHiJUiV97XgUUxEqL3ceM2tesc5RpQ5BmtXIwA2WbTL7575L+59/Tc0N9n1jiKECFPvVe/GqSSSqLYwzdmLYhi9L+Gj9zsPch2W/plGQdYm/U3ZrlYAauOS9Q2iszaXi9eKz+LPORfxq/df0DuOECIMNVdUsC5/PABzm3az6L5/1TmRvqR4CVKtiVrxkt7XrHOS40uz9++0G5mF1+vVOY1+nnx1KXXGHACqU5N0TiOECEdPv/cC5aYxGFU3k0prMEVE6h1JV1K8BKmmWK14SekOvjbpAQWRCSiql1YlmZee/63ecXSzpdA6+N9VsaNrl0shhP91OBxsHV8EwPSu3dzx45/rnEh/UrwEIa/XO9gmnRok06SHcv33f0SG2gDAEW/wdUQFwjvvvMT26EmDn1cZcyk/fEDHREKIcPPH537Llhhtoe7ZZYeJt4zut+pBipegVLN7O40R2uZvwdgmPSAuIQlrtw2AhpTRuWj3U48TtxJBnqeSWLWdPiWSd/62XO9YQogw0dPWxp6xGfQpkeS7K7n5wuv1jhQUpHgJQl+t+Ptgm3SWOVbnNCeW3dYKQG3C6Fvr0eZysT5tIgAldaXk91UDUCfznoQQPrLiuT/wVap2d3d+7QHGTZmmb6AgIcVLEKp0OwfbpG+6NvjapL8pvX+n3ZrITJ2TBN5vXltKgyETs9qJpdxBXpsDgOqk0VfICSF8z9Pby8YEN02GVOJVFxfGZ+kdKWhI8RKE2pODv016QHaPCkCjIYP33/qjzmkCa0thPgBnt+8mor2XbEcrAJVm+QUjhDh9W996i/UF4wAoadnFZdffonOi4CHFSxByWvoHMrqDa5r0UG5e8k+keLU7Djsbq3VOEzivv/EMO8za1tyTyyqJdvVgtmlt7fWGLD797EM94wkhQpzX6+WT+v2URozDoHqYUSkbYH6TFC9BqCmu/85LEE6T/rak9Exye7VFu42p8TqnCZw1xh48iolCdwVqs5ei8ZMw272keB2oioG1+7frHVEIEcL2fvghOycWAjC1Zy/3P/QfOicKLlK8BBmP2409iKdJDyW3XbvjUGsZHWs9mpvsrEs9A4A5tWVEtrZy+T23E2+IJL+3FoD6lAQ9IwohQtzftqxhY5zWHj3rULnOaYKPFC9BpnL7JhwRWhEQzG3S35TR1L9oNypd5ySB8bu3n8VhSCdG7SCu3E5qRAJR0WasE8djdWqFXJUlReeUQohQdfDTT6k6I4dexUyOp4aHvv+A3pGCjhQvQWbzxx9jV7RFupkRZp3TDE9qq1Zk1RuyWPOPv+qcxv82FxQAMKttN1Ft3Vxy1/8DYO61V5PZ2ApAZWSOTumEEKFu5d9X8GVGf3t03T4s0sF4DClegkwFnbiVCExqH9dc9T294wzLVZctIl514lWMrNm9Se84fvXSq79nV5T2ltGZpZXEdajkTRwLQHJWGsaaZhTVQ7MhhT+99pSeUYUQIah2yxZsYyw0GjKIUdu5ynqG3pGCkhQvQaajfyBjmmonOzNX5zTDM3bKNKx9dQA0poT3ot11ZhWvYmSMuxy1qY+p3zn/qK/HtxvJ9tYDUNrZqkNCIUQo++yDD9hQVAxAiXM3F55/hc6JgpNJ7wDiaK2W/uIliKdJDyW3o5l9kVCXmKh3FL+ps9WwLkW7lVtSXYa5tZ3zvnvdUcckxSWQ32WjNi6XulE6MkF8TfV62fLh06TseoFItVvvOCOk4DDn05U5E8u4+RRNPQdzTJzeocKao6yMirg+9kZORFG9zKp36B0paEnxEmSa4/r3eAniadJDyWxyQhLURAf3pnqn4w8fvEzz+CuJVduIrnCQnZqDwXD0zcsJ585nW/NOiIOq+FSdkopgUF9Ziv3Nuzm7e4veUU5Zblc9VGyAiqfoXWmkNGIsLSlnEVlYQt60C0jNzNM7YlhZ/eZb7J9aAMCk3v3cf/e/6xsoiEnxEkQ8bjf2aO1tl5T24J0mPZSUZhcAtcZs9u3awhlTZuqcyPc2F2h7Lsx27SG6pYurHn3kmGNmLjifT36/GqxwJMJKm8tFfIK0TY8mXo+HzX/+Hybt+x1ZSjeH3Eksy7gFtynEft2qkNDZQUnXDiZ4DpCqtDLeXQoNpdDwFmyAOiWDuvgpeHJnk3bGueRPmIEx1L7PIOGqq2O3p40NCRcBMOewtEefiPwtCyKH1n+JPUIbdZ4UIm3SAxbMu4RnejvpVmL48B9/Drvi5YU/PsnuwvMAmHCgEos7goTUY8fSmyIjUCqdREzvpVOJ46VlT3H/P/1roOMKnVQf3Enb23czu28Pdb2x/KL456zKP4sGQ+jO/opRr8LqriW3w0F2s5O0lmYubv2EqcZKsmkg27UK9q2CfdCmRlMRfSYd6TOIL55HwdRziRuFQ1tPxdrXXqO5OJ0uJYYMr40HF/2T3pGCmhQvQWTrp//Afu7VAGQYInVOMzLT5p6H9R/vcTBiLI1J4fe++DpLBKpiZFxfGYq9m/mLbjrusYneKPI8NZSbiqiJUgKYUujF3dfLljf/k2nlzxDbBz8r/BdWFU6jzqi1zMeq7aR4Q2sdmwcDDYYMOpVYSiPGUZo4DhK1rz2tXofVXYu1005Wi5M0h4Obmt+hKKqVKd1boGoLVD2H5x8Kh0xFNCVNw1gwh9wp55OZV6zr9xWMulta2NTcxFcztH8gzbftITlloc6pgpsUL0GkytCNW4nAqLq54rIb9Y4zIorBQE6ng4OWsdiSwmuhamVVOeuTtTlGs6sOEePqYer58457fFZuLnkdDZRbiqhOln91hrvDezbief9uxnRU8ljhA6waM5VqoxXQipYL7FvIPlCHosbonHTkvO5uTIkGmtMSqUuyUB2dTo0xh24lmoMRYzloGQsWoACeVm/A6qklt9NOTksrKY1OZlRvYnJEBZO7/4rZ/i6GLdBACjVxk+nLmUXKhHMoOHM2EZFRen+ruvry1VfxWmOoM+ZgVru4Or1Q70hBT4qXINKZrN2xSFPt5OeF3tsu2a1OsEB1THgtVH1u5Zu0Fl9BvOoiutxO0fgzT3j82VddwcYNfwYLVMWOjl2HR6Oe7k62v/5vRO/9mHfGLWLVnClUmgoAiFY7Ob95C9Y9NcR2G5idX0hmUZG+gUdI9Xppb2mlraWZ9vo22g7V0tpbSbP6FWpqBC3pidQlWqiOTaPGmEOPYqbcNIbyhDGQAORDxMxryfXUktulFTTJjlZMtk6Sm3tIOryb1jXrqFW6cUdGo6bkED9+BuPOvZzkzNGzyaO7q4vNtgY2nzsRgNltu1hw9Z06pwp+UrwEkdaE0GyTHpDicEI+1JhyaKytIj0nPDoRNuVpLzpzWncT7ezmil/cccLj888cR+JfWqAIqo25lB7YyfgJUwMRVQTI7i9XUv7mH9k+4Qz+cdljlJvGABCldnNeyxaK9hwhpsvE7PxC5v3gB5iio3VO7Huevj7abTactbXs372bbS1HaEqNoy45kZqYVKpMOXQrMVSYCqmIL4R4IA+Mqpscby15XY1kt7aSYnfibu3D6ImGfXWs3Ps8UT09mPt6MHtVYqOisMRbiE+IJz4piYTUVCzZ2STkWokJgzubG19fRl+8m539m1/OrmnUOVFokOIliDTF9xcv3S6dk5ya2ePP4lm1ly4lhteX/5GHwmAK6h+ef4K9Y7XV/8VllaRFxBNhPvktbkNVO7Fz2ulQ4nhv1Qc8LMVLWDiwfgtfvPwi9RPT+OzqGymLGAdAhNrLuc4tjNtTgbnTxGxrEXN/8AMiY0LvraLhMkZEYLFasVit5M2ZwyXf+npzk51X33qOuiiDVtDEplAVkUOnEkeVMZ+quHyIA3K1gibbW09udwM5ra2k2r10tnoxKv0/v75eaGrSPg4dGryGye3G3OcmBpVYo5HYqCiMxtB6WTvY3s7BWfmoipGJvQd46B5pjx6O0PpTDmN9PT04+qdJJ4dYm/SA8y+9luzP/0qlqQBbfGgtOD6eTSkxqIqBib0HMNi6WHjP/cN6XorZQn5fNfsiJ1KXGL4vYKNBV1s76/68gr2rV9I6JZPV15/HvkjtFr9J7WN+21bO2F1OZIeJWdl5zH/oNiJjY3VOrb/klDQeuOffjnqszeXitbeeodLgoT4pgZq4FKoicmlX4qk2WqmOtUIskAOK6iHLa8Pa00COs4U0uxOl0YVBjUE1mlCNBtwmE+0mE+0DF1BVcPcF+ls9Ld4oL+sStenRJUcOneRoMUCKlyBR9tVq7P3TpEOtTXqA0WTC2m2nMq4AW3LoL9otP3yADUnajrqzqsqJ71DJHTd2WM8tOuss8toa2ZcCNUnHtlSL4KaqKhXb97Du3Q9pKN9C1/QJfPHdC9kdpa13MqpuStq3M2XPISLbDJydZWX+/T8kKj68x2OcrviEBO7+fz876rE2l4u33nmRCncXdSn9BU1kNi4lkTpjDnUxORADZIGieslQG7B228hxtZBhdxJZUYexrRsiYlBMMaCEVoefbXoGHUo8qV479179A73jhAwpXoLEji/WYT/vcgDSQ6xN+puyW1shDmriUvSOctpe+PRdnGMvx6K2EnWwkSnnf2fYz51z9WV8/ubvIAUqzaG7x8do09nWxvo/r2D/2lX0dNjonnYma2+9jh3myYB2N2BO5w7O2lNGhMvIzPRszr33h0TJRoSnLD4hgTtuf+iYx1969feUd7moT06gJj6Z6shsWpRkbEoWtpgsNscAmcBkSPc2YO2xkeNqJqYntO68rMnRCuL59j1kZ16sc5rQIcVLkKg2dNCnRGJU3VwVYm3S35TqcEIu1ERk09rURGJK6BYxG/O0RZhzWvYQ6+zjvJuvO8kzvhaXGE9MfStMgHpDFp+sepcFF1/vp6TidKiqyuHtu1n/7l9oKN8KqpueKeNZd9Z8tpqnoCoGFNXLzK6dnL27lAiXgRlpmZz7Lz/EHMazvPR22y33HfPY6288Q6nTQX1KArXxSVRFZtFkSKPRkEFjdAZbQ3RddKTaw8Lo0P1dqQcpXoJEZ7L2Hnmq6gjJNukBE5OzUFQPLsXCG8ue4e77/u3kTwpCv3/2l+wffxmK6mVsWSU5aZnHzDE6mViHlxSvnSZDGuvK9kjxEmS6XG2se/dv/XdZGgDoPaOYDTMnsyl6KqpiBGB6905Kdh/A2KowPSWNc39yW1h0uYSi7y8+dtfZd99/mZ311diS46lLSKLXEHova1Nrq7nm9of1jhFSQu9POUw5Q3Sa9LddffPtPL72U+qMOdRGePWOc8o2p2trds7sPYChrpMrfzHyXywpSSnk99TRFJ1GXYq8rRAMVFXl8Lad2l2Ww9tBdQPQO24MW2ZNYn3sdLz9RcuUnj3M37MPY4vC9MRkzvvn24gJ4TuJ4er6a3+I/LNg9JHiJUg0DU6TDs026QGRUWZyexqoi8nBlhqai3Z379nGhkRtjcPZlYdJdpuwDDHH6GSmXHQRXzVsYFs0VFvkRU9Pnc5W1r83sJbFPvh439hits2ewFdxM/Ao2q/DM3v3c+7evZiaVKZZEvnOg7cRmxpeGy8KEeqkeAkC3R3tOAamSXe0n+To4JfjaoUYqI0LzS6b1zZ+TFvRpSSpzUQerGfedxef0nkmnzOL9D/8HTLgSOTo2TE0WGh3WXZ84y6Lp/8rJtzFxeycWcTahBn0KdoC+Ql9pZy3bw8Rdi/T4hM4//4fEpsuOyQLEYykeAkCpWs+G2yTTnZ16Zzm9KXbnZAJ1ZGZeL3eEa8V0dsGa/9C3ea9WFwqU849/hyjEzGYjERWtaBM9tBiSOFPrz3FrUvu9WVUMYSO1hY2vP839q35B72djsHHDaZ0lAl5bD0zjS8SZ9CrmAEY4y7nwgO7ibT1MSU2lgt+fBtxmdIhJkQwk+IlCOz6ahP2C64EILX//fZQltf/L9kmQxpvvfIUi287tmsgWD35h/+i7IwrUFQPRWVVFE0447TOl9QTRba3nlpjLqWdrb4JKY6her2Ub9/Ohnc/7L/LMrDeKoKY5MlEnZHKmiyF1ckz6Va0lpQC9xEuKtuJsbqdSZFRXHLPPSRkZ+v3TQghhk2KlyBQF9lNrxKFUXVzzWWL9I5z2hbdejdLN2/GbkjnQGeL3nFGZEuOtjZlcu9+Imo6uHLp6Q1IsxaPJb/LRm1cLnWpiT5IKL6pvaWZjR+sYN+aT+jt/Hqxu8GUSc7E+ZiKI/ibwcHnKTPoVLRF8VZPFRcf2knEkTaKezu5/IF/ISkvPOZwCTFaSPESBDoTtV+qKWpTSLdJD4hLSCK314bdnE5jauh02ezcsZGNCdqOujMrKkiLjMMYGXFa55x99VV8uuZ1iIPKOFn06Quq18uhrVvZ+P5HNJTvAAbuskQSkzSFKRcuoDfVxZsN+/k0bSYdirYJWLanlosrdmAudzKurZ5ZS+6geM65en0bQojTIMVLEBhsk3aHdpv0N+W0tbDdDLXxiXpHGbbXtn9OR8FCUrwOIvfVs+DHx+76OVIZBdmkvtMKVqiMyKW5yU5yStrphx2lDm/fzd9+/wS937ijZzBlkT1hHiXXL8TWeZCXD27gU3UGrozzAcjw2lhwZDvRB1uY7DxC1qwZzLv1VyghthZLCPG1gPzf+/TTT1NYWIjZbGbGjBmsXbv2hMd/8cUXzJgxA7PZTFFREc8++2wgYurm62nSTp2T+E5Gk/a91JgzdE4yfBtytLlFJU17Se4xYB3mHKOTMR5xEqn20KnE8ery531yztFq9atvaYWLEkVM0tnMvv4R/um5pznj2in8z7bXucXVx/uZF+BSLKR6G1lc+TFXf/YZs7duY6H5AOf/4jfM/+HPpHARIsT5/f/g5cuX88ADD/DII4+wfft2zjnnHC699FKqqqqGPL6iooLLLruMc845h+3bt/Ov//qv3Hfffbz77rv+jqqLTlfr123SITpNeigZ7dp8EZuSyccfvqlzmpP7zVOPcShiLAbVQ35pFVMv/I7Pzp1hTsTqqQGgJjr0F2TrRfV6aW04CMCsa37MXc/8f2TNTOLn7/6OmxuaeSfnIlqVZJK9TSyq/oTrPv+Ukk3buNX9V4quupy5j35IaqasbREiHPi9eHnyySe5/fbb+dGPfsTEiRNZunQpVquVZ555Zsjjn332WfLy8li6dCkTJ07kRz/6Ebfddhu/+c1v/B1VFwdWf/qNadKh3yY9YMmt95CotqAqBjZW7Nc7zkltsWr7eUzp2UdUbRfnLbrBZ+eeMG8e+R2NANQkJ/rsvKNN1d4yVE8nEEHOjFx+8sqvuaaiirfyFtBkSCVRbeGG2n9ww2crmbthG/f1vklBsUrqv3zO9IU/0Du+EMKH/Fq89Pb2snXrVhYsWHDU4wsWLGDdunVDPmf9+vXHHH/JJZewZcsW+vpCa1rocOzfuJVGg7YGIo3QGuV+IknpmeT21gPQGORb469f/ymb4rUddWceLic33bd7fMxceCE5Ta0AVMaEzttowWbf2s30qX3YLj6b66oOsqxgIXZDOvGqk2ttn7Ho878zf91W7u19m/PSN1N76W85+8G3saTIz1yIcOPXBbsOhwOPx0NGxtG/PDIyMrDZbEM+x2azDXm82+3G4XCQlZV11Nd6enro6ekZ/NzlCq3t9WsieuhVzBhUDwvPv0LvOD6V29HMniiosyTqHeWE3indQmf+JaR5G4nYW8eV//1rn54/KiaKxIYWKIJqYy6lB3YyfsJUn15jNKjZv5vqqy/m3eyLAIhV27jQvpXMfbUU220siF1LRnoLG1Ov48xbniQzQYYnChGuArJqTVGOvqOgquoxj53s+KEeB3j88cexWCyDH1ar1QeJA6cz6es26XB7QcscWLQbHdxbrK/PGQdAiWMfGWo08Um+f9GLrOwgVm3DrUTw3j8+8Pn5w53H7abNUc729CIA5rdtYsnqvzDny73c3fFnlqR/SF+cmX2XvMXse18mTgoXIcKaX4uX1NRUjEbjMXdZGhsbj7m7MiAzM3PI400mEylDTHR9+OGHcTqdgx/V1dW++wYCwGXRBjKmuZt0TuJ7KS3anKZaQzbb1q3ROc3QnnjqMSpMhRhVN9YDVcxbfKNfrmPNKaCgT1u0W9e/r48Yvood++lNieOIMR+AmTtK+V7nSm5LeZu02BY2ZH6P1J9u4YySS3VOKoQIBL8WL5GRkcyYMYNVq1Yd9fiqVauYO3fukM8pKSk55vhPPvmEmTNnEhFx7IZhUVFRJCQkHPURSpoSBtqkQ+vtruG49tJFxKjteBQTf1+/6uRP0MHmPK2Inta9l9jaHiafM98v15l55WXktWlzdqoT5a7ASB34ajOuM8fgVYykexu4xfU+xXHVVBjyKb/6L8y562nMMXF6xxRCBIjf3zZ66KGHePHFF3nppZfYv38/Dz74IFVVVdx1112AdufklltuGTz+rrvuorKykoceeoj9+/fz0ksv8cc//pGf/OQn/o4acG1NDuxhNE3628ZOmYbVXQeAIyX4XlhWr17B5jhtoe6M8sOMmTjBb9cqmjqRLHsrAJXRWSc+WByjtnQPFVnawvZxHVUkx7SyPu//kfOzTYybfp7O6YQQgeb34mXRokUsXbqUxx57jGnTprFmzRpWrFhBfr52+7e+vv6oPV8KCwtZsWIFq1evZtq0afznf/4nv//977n++uv9HTXg9q9ehWOwTbpT5zT+kduhvR1Wl2TROcmxPqjaQ7cSQ4bXRsTeOi6/+/TmGJ2IwWAgrlbbQdmmZLJi5Tt+u1a4cff10t5cQVmctp6tqMHO3jn/Q8lt/0NklFnndEIIPQRkPMDdd9/N3XffPeTXXnnllWMeO++889i2bZufU+mvdMsOGi++FoBUb/i0SX9TVksrJEJNTPBtib8+azwAJfb95EYnYYqK9Ov1ktoiSPXacRjS2Fi+n8v8erXwcXDTbvryM6g15gJQVHmI8Yv+SedUQgg9yR7ZOqqLctOjmFFUD5ddcKXecfwipakNgBpjDocP7NE5zdcef+o/qDTlY1L7yNlTycV3/z+/X7PwzMnk92hvo9WH0MBKvZWu30LTOO1OrdVTzSzKpZtIiFFOihcddYVxm/SABbMvJFLtoUcx8+7f3tI7zqAt+dkAnNW1h8RGN7nFvpljdCJzrrkcq1N7G60q4djOOTG0+oN7OZyuTeQe11ZDz5g5OicSQuhNihcdheM06W+bPv98cj21ADQGSYvwJ6veZUustlD3rEMVTL7gnIBcNzE9mczGVgAqo3ICcs1Q19vdTWdrJaUx2p2XAlsjyZMv0TmVEEJvUrzoaHCadE/4TJP+NsVgIKfLDkB9kCza/ch2mB7FTJa3DvPeer5z86KAXTuqqgVF9dCiJPPSq78P2HVD1YF12+mbUIDDkI5RdXPG4e2MnXau3rGEEDqT4kUnLQ11X0+T7gifadJDyW7VirPa2GSdk0Cby8WGLK0luqRhPwUZuQG9flZECjlebd3LoTDc28fXDm7ciq1Iu0tV4KmkONGD0RSQPgMhRBCT4kUnBz5fhSNSW3SY7Azv4iXNrhUvVRE5tDQ06Jrl96/+jmqjlQi1l+w9lVzx4I8Dev0pF1xAfpe2g3RtSmJArx2KbIf3cShNGy9R7KxDHfcdfQMJIYKCFC86KduxZ3CadIpH5zB+NqNgIkbVTYcSz6tvPK1rli0F2p2W6V17SG2BhOTA3g2a8p0ScltaAKiKSw3otUNNV3s7bc4KysyFABTUNWKdebnOqYQQwUCKF53Yojx0K9EoqoeLzwvvX8gLrr6ZLG89AA0xx454CJSP/racrTGTADir7DDzvndTwDOYIiNItWnFS2WEleYme8AzhIp9a7fiOWsKLsVCpNrDtOrNZBf6bxdkIUTokOJFJ13JA23SzUyeNF3nNP5lNJmwdjcCUJ+i36Ldla3V9CpR5HhqiNnXwOT5/pljdDLRVe1Eqj10KrH86e3ndMkQCsq3bKc2X3vLaGzfYTJzgns6uRAicKR40YlrsE06/KZJDyXb2QpAbZw+m4u1uVysz5gIQIntABPPnKxLDoCi/HFYPdqE6RqzLD49nsaKfRxK1gZnFrfYiDpjgc6JhBDBQooXnTTHa4MKU3tGR8dJukNbtFsdmY2nry/g11/66u+oM+YQqXaTtauahXffFfAMA2ZffSX5HdqdqJrkRN1yBLOOVhft3TYORhYBYK1pYMysS3VOJYQIFlK86MBRdXhwmnRqGE6THsqY2BQU1UuLkswrL/wu4NffXKgN9ZvZsYecvkgiIv07x+hEcsYVkOPoX7Qbk6FbjmC254tN9M6cQpcSQ6zazuyGDSQkyq7EQgiNFC86OPDFp9gjE4HwnSb9bTctuZN0VbvbcETtCui1333/ZbZHawt1p5Ud4cJ/uiOg1x9KUv+i3WpjLrv3hP8Q0pEq37qdqlytG6+4p4LYcTN0TiSECCZSvOjg0K4D2PvbpJM9qs5pAiMyyoy1R9vfpCE1sIt2/9HdQp8SSZ6nirgDDVjHjQvo9YcS3+glTm3DrUTw0eoP9Y4TdBxV+zmYmAXA2KYGkqcu1DmRECKYSPGigwazly4lBkX1cvHc0bMIMdul3W2ojQ/cot02l4v16f0LdesOMP2iCwN27ROZUjKfgt5qAOoS43ROE1yc9iY6TR0cNhUAYK2uY+xZ5+kbSggRVKR40UFXsvZilaw2M3XabJ3TBE5Gk7Y4uSYqA6/XG5Br/va132EzZGFWu8jYVcV5NwZ+b5ehnH3FJeS1OwCoTtJ/bEIw2bN6Ix1TJ9KnRJKoNnN2+25MEfqtURJCBB8pXnTg6p+unOYZHW3SA7I9CgANhkzef+PFgFxzS6E2jXhm+27GxKaiGILjr3xMfAxZ9lYAqsyZ+oYJMhXbd1CZre0+PL7rCLFTLtI5kRAi2ATHb/JRRPV6aR6YJt0dvtOkh7LkBz8myasVbHuaa/x+vTfeeo7tZm2h7pSyI1x6f2DnGJ1MXK32519vyGbFynd0ThM8mmpKOZigDWMcY7eTMyO8d6AWQoycFC8BZqsowxHTP026M7wHMn5bXEIS1r7+MQEBWLT7BZ14FBOF7goyD7eTkBxcrbZWQyqpXq0Da+PhfTqnCQ5NtTa6kg0cMeYBkF9ZQU7RGTqnEkIEGyleAuzgms+xR4yOadJDyW1vBqAuIdGv12lusrM+TXvRm1NbRsniG/x6vVMxfeElFPTUAVCn49iEYLJn9UZcE8egKkYyvDbOjm0Lmrf6hBDBQ34rBFj5vjLsRu39/OS+wCxaDSaZA4t2zf6dU7N0+TM0GjKIVjvJ2FHJ5Pnn+vV6p2LCrKlYnVoxV50QXHeF9FK5cyeHM7X/P8Z1VGGedLHOiYQQwUiKlwBrjIJORVvzcu6cC3ROE3gpbdqmfHWGbNZ8/JHfrrO5sACAWW27mVSk/74uQzGYjGQ09k+YjsrROY3+VFWlua6UslhtN+SiBjtjZsl6FyHEsaR4CbDuFK1wSfY2MWumPlON9XTjNT8gTm3DqxhZu2+TX67x6utPsytKe8to0oEKFtyp3xyjk4mtasOgemhRkvnjy4EfmxBMGo5U0221UGfMQVG9FFUdxJKcpncsIUQQkuIlwFz9G5KNtjbpAXljx2HtqwWgMSXeL9dYG9GLRzExxl2OtcFDZJTZL9fxhYmFU8jxauteyntH3xqob9q7ehOOYu2ui9Vbw+TcRH0DCSGClhQvAXRUm3TP6GqT/qbcDq1wq0tM9Pm562w1rEs9E4DZNQe56J9+5PNr+NKcqy8jv0sbm1CTmqhvGJ1V7dnJ4bT+eUZtNSTKSAAhxHFI8RJAtWV7cURrd15SRsk06aFkNmuFW020798SePr9l2gypBGrtpO9rYrc4vE+v4YvpVozyW3unzAdN3rfIlFVFXvtXkpjtE0Fi+rsFE8/X+dUQohgJcVLAJWt+Qx7ZH+b9CiZJj2U1BatcKsx5rBv20afnntzQREAs127mTU/NNYUpfVPmK6MsNLcZNc5jT5qSw/jnphPkyENo+pmnG0fEZFRescSQgQpKV4CqK68ZrBNOqXXo3Ma/Vx1yTWY1S76lEj+8pnvOo5eeOl37I7ShjCeuf8I59x4s8/O7U8Wh0qk2k2XEsMrbz6rdxxd7P1iEw2F2hTpQvcRxk6YrHMiIUQwk+IlgGpNbjoU7W2jeTPm6ZxGPxOnzCLHoy1SbUqO9dl518cb8SpGivsOMrbbHDKbm80692Ly3Nq4hNrYCJ3T6KNm3y4OpWYAMM5VR87ZV+icSAgRzELjt3uYGGiTTlKbKSm5UOc0+rJ2am+P1CUn+uR8lVXlrE/WFurOqT7Egvvu9cl5A2HqhfPJ79B+HjXJSTqnCTzV68Vh20+ZuRCA/Fo71rFTdE4lhAhmUrwEUNvANGm3Q+ck+stq6V+0G+ObnWWfX/EGLYYU4lUXObtqSEoJncWvUTFR5DT1r3uJydA5TeAd2VVK3/QJtCkJRKrdTG49GDJ3zYQQ+pDfEAHicbtpTpA26QEpjv4xAaYcGquqTvt8GwvGADDLuYcLbwi+OUYnk1SvFS/Vxlx279mmc5rA2vflJmrztHERxX0V5M+8SOdEQohgJ8VLgFTt3Y4jun+a9Chukx5wzpRZmNQ+OpVYXn/7hdM619PP/5q9kRMAOHPfESafE3ottjl92l0jj2LiL599qHecgKo9sIdDSdodp+IWG0UyEkAIcRJSvARI+ZdrsUcmAl/P9xnNzr3oSrL7d5ZttJzeDribkmNQFQMTekuZGp/li3gBd+4N15Pfqy3ataXE6ZwmcDxuN00tVRyM1FrcrTU2ktJC889QCBE4UrwESH1lPY39bdJJPaO3TXqAYjCQ26UtUq1PtpzyecoPH2BD8iQAZlUd4uI7/8kn+QItb9IY8tq0tVDVick6pwmc8q176Z1ZTLcSTazaxlmMzrEZQoiRkeIlQOoj++hQtLeNSqbN1jlNcMhxtgJQG3vqL9YvfvourUoSCaqT/AONQT3H6EQMBgNZjlYAKs2j587DgXVbqMrVFleP66kgr+RqnRMJIUKBFC8B0pWsvRWQqDZzzvxLdE4THNLs2sLl6shs3F1dp3SOjXna2w1zWnZz1V3BOz16OBJrtUXMNkMWH370ls5pAqOubA8HLdkAjHU0UDxjdG8hIIQYHileAqTdEgNAmltuiw84M92KonpwKok8/9xvRvz8/3v6v9gXORFF9XLG3kqsxRP9kDJwpo05mzRvIwBbqkp1TuN/7t5eWjx2DpsKALDW1hNljtE3lBAiJEjxEgB9vT00J2h3XqRN+mvXLLqNTLUBgNoodcTP35ypbeg2sa+UeeNDf1Oz2VctJL9HW8Rcn3rq64BCRemGnXRMGYtbiSDJ28SM5NB8y08IEXhSvARA5a7N2GO04iW1U9qkBxhNJnK7teKlIWVkL9b79+1gQ2L/Qt3Kcubf+D2f5wu0hBQLVmczANUJvtm8L5iVbthKVba23mV81xHGnnO9zomEEKFCipcAKP9qHY7+NukkaZM+So5L25ytNn5k2+K/vG4FLsVCotrM2CpX2OzImtmg/TyOROXonMT/bIf2UhqvfZ9j7Hbyxp+lcyIhRKgIj9/4Qa6xuoFGo/YvzJRut85pgkuGY2DRbiaevr5hP29Tnraj7pzmvVx//z/7JZsekhvdGFQPrUoyL770pN5x/Ka3qwtXbBdVxjwArHU1YVOACiH8T35bBEBdpJv2/jbpWZNm6pwmuOSZtJEJDkM6b7781LCe8+T//gcHIsajqB4m7DtCUlr4zAM6b8E15HhrATjkDt+7dPu+3E7rxAJUxUCmt575UybrHUkIEUKkeAmA7v4dUy1qK9/5zmU6pwkuS267j1Svtlld2TAXM2/N0+5iTeo9wOXnh9fP84y508nv1NYB1aWG74Tpg5u2UJHRv79LRxWFMhJACDECUrwEQHuiVrykyzTpY0RGmcnttQFgT0046fFbNq1mo0X7V/qsinImnxNe+4KYIiOwNmvrXqriUnVO4z+NFfspi7MCUGRzkJKRq3MiIUQokeLFz3q6O2lO0PauSO1t1TdMkMpp1zpsahNOfqdh+e71tCvxJHubmOoKzzELqQOLdk1WmpvsOqfxva62dpwZBuoN2Siql4KmOr0jCSFCjBQvfnZk2wYcMdp6l1SZJj2kDIe2s2xNVDper/eEx26wjgWgpHkvV9//M79n00O+O4UotZtuJYaXlj2tdxyf2/vFFprHaHda8jzVzJu3QOdEQohQ49fipaWlhSVLlmCxWLBYLCxZsoTW1tbjHt/X18fPfvYzJk+eTGxsLNnZ2dxyyy3U1YXuv8wOb1g/2CadLG3SQ8rq0jqw6g2ZfPzuG8c97n+W/n8cjCjGoHqYsL8ybHdjvfDmm8lzaxOm6xPCb+O2Q1u3U56eDsC4thqKz75I50RCiFDj1+Jl8eLF7Nixg5UrV7Jy5Up27NjBkiVLjnt8Z2cn27Zt49///d/Ztm0b7733HmVlZVx11VX+jOlX9tomGo3ahmNJ3cNvBR5NbrjpNhLUVlTFyOba42+Lv61AG1g4pWcfNy36YaDiBVzWmFzyOrQxAdXJifqG8YOawzsojS4AoLDejjk6Vt9AQoiQY/LXiffv38/KlSvZsGEDs2drU5RfeOEFSkpKKC0tZfz48cc8x2KxsGrVqqMe+7//+z9mzZpFVVUVeXl5/orrN7YID22Ktnvs2WfIJlxDycrNx7p3HXsjE7GnDL1o98vP/sbGBG2h7szDh8m/9PhFcDjIcbRCIlRGh08bOEBbcyu9Y1NoNqRgVN2M75O3UoUQI+e3Oy/r16/HYrEMFi4Ac+bMwWKxsG7dumGfx+l0oigKiYmJQ369p6cHl8t11Ecw6UrTOo0S1FYuvCB07yD5W067NrCy/jh/zu9V7aFTiSPVa6ckzF7Qh5Lc0ApAjTGX7ds36BvGh/au3oyt/w5akbuCkotu1DmRECIU+a14sdlspPe/r/1N6enp2Gy2YZ2ju7ubn//85yxevJiEhKH/Rf74448PrqmxWCxYrdbTyu1rHRbtlniaR6ZJn0hWc/9Ou+Zj/854PR425hQDUOLYy6W33hHQbHqYaS0hXnXiUUysWLNC7zg+U759O4dS+9e7OOspPONsnRMJIULRiIuXRx99FEVRTvixZcsWABRFOeb5qqoO+fi39fX1cfPNN+P1enn66eN3XDz88MM4nc7Bj+rq6pF+S37T0eakOUErXtJ7WvUNE+RSnNpi5lpjNltX/+Oorz3xv7+g3FSEUXUzsawOg9GoR8SAmnvdJRT0aot264ax/02oqDuyk7KoIgDyaxtlJIAQ4pSMeM3Lvffey80333zCYwoKCti1axcNDQ3HfM1ut5ORceLb/n19fdx0001UVFTw2WefHfeuC0BUVBRRUVHDCx9ghzd/haN/mnSKTJM+oRuuuplnatroUmJYufVzZnzn6w6U7WO0u2lTe/bygzsf0ClhYMXEx5DX1sTuKKhOTNY7jk+02hx0T86jXYnHrHZxliU8u8WEEP434uIlNTWV1NST7/xZUlKC0+lk06ZNzJo1C4CNGzfidDqZO3fucZ83ULgcPHiQzz//nJSUlJFGDBrVW7ZiP+sMAJJdHTqnCW5F4yZjPfIOZRHFNKXEDz7+9/deZ1Ni/0Ld8gqSL83WK2LAZTW2QCpUmsPje969eiO1Vm0kwNjeCkou+Z7OiYQQocpv92wnTpzIwoULueOOO9iwYQMbNmzgjjvu4Iorrjiq02jChAm8//77ALjdbm644Qa2bNnCsmXL8Hg82Gw2bDYbvb29/orqN422JuxGrdBL7pI26ZPJ7dTGJ9R9Y9HuJ646upQY0r0NXFx4pk7J9JHi6AagwZDJ+++9om8YHziycweHkjIBKG6xkZqdr3MiIUSo8usbzsuWLWPy5MksWLCABQsWMGXKFF577bWjjiktLcXp1BZr1tTU8OGHH1JTU8O0adPIysoa/BhJh1KwaIz04upvk55SfIbOaYJfVkv/34MYreDr7exkfc44AErs+zjn8tHVmXL5ZYtJ92pvvW6vr9Q5zemrt5dxMFJb72KtbdQ5jRAilPltnxeA5ORkXn/99RMeo6rq4H8XFBQc9Xmo607V1rvEq04uWzi6XnhPRUpTGxRAjSmHyr07ef3jdzhy1g2Y1D4mHTl2/VS4G3vWGeR/vIXGmAzqUy16xzkt9qp6Oifl0aOYiVddnDNxot6RhBAhTJb6+1G7pX+atEemSQ/HwnMvIULtpVuJZvmKd9g+Vntb4ayuPdxx37/pnC7wDCYjVqc2tLIqIbQnTO9ZvYGqHG29S3F3BTO/c53OiYQQoUyKFz9xtjhoGZwm7dQ5TWiYPutccj21ANSlJ7AlbpL2ePkRzP1dW6NNRv9mdUeicmkLsg0YR6Jq9y4OWrSFx8WORqJj40/yDCGEOD4pXvykYtNaHLHaC25qZ5vOaUJHTpcdgI+tM+lWosny1rGoZPROHc73WjCoHpxKIstef0rvOKdEVVXquqo5bCoAINcm612EEKdHihc/qd6+E3v/NOkUl0yTHq7sVu0uVYui7W0yp/EAZ5w9T89Iurrqu7eS69XuRh0xheZ6sPpDVXRMtOJRTKR4HVwyiotRIYRvSPHiJ00NrdhN2jqFpK7Qa/PWS1rT12+xRai9TLWF7lslvpCcnUZepzZOoy41Ud8wp2jvmo0cydLWu4zrOsKkuRfonEgIEeqkePETe4QXp5IIwKTCCfqGCSGzx0/FoHoAmN61hzt+PPoW6n6btbkFgMq4NJ2TnJqqPbsoS8gBYEyjfVSMdxBC+JcUL37SnTbQJu3iyssX6ZwmdFy88HrGug9jUD2cXXoYo8mv3fwhIa2xFYBKUy5VleX6hhkhVVVpNDZRbdBGPOS3yOJ1IcTpk+LFT9oTB6ZJS5v0SCgGA4t3bubHW97ixz+4X+84QWHOmLlEqd10KzG89e7LescZkaq9h2gdn4uqGLTF1zffqXckIUQYkH/W+kGLvX5wmnRab6u+YULQXT99Qu8IQWXepQvI+/JDDkYUY0sOrZbxfWs3czhHW/s1rr2K9JzLdE4khAgHcufFDw5vWENTzECbtEyTFqcnKiaKvA6thbw6OUnnNCNTvW8nZbHaZoNFNrkLKYTwDSle/KBm527sUdp27jJNWvhCjkNbtFsVnalzkuHzejw4krppMGSiqB4mqIrekYQQYUKKFz9odrR9o026R+c0IhykObQiuMaYw7ovVuqcZngqdpTiKMwCIN9TzY3fv1fnREKIcCHFix/YTV5aFe32/qS8Yp3TiHBw1QWLSFCdeBQT/9j+ld5xhmX/V5soT08HYFxbLTFxoT1cUggRPKR48YOedG29S6zaxlVXfU/nNCIcFJ89mfzeGgAa0hL1DTNM5fu2UhZdAEBBvV3fMEKIsCLFi4+pXi/tSTJNWviWwWAgr037+1SdmKxzmpPz9Llps5ppUZIxqX2cnZ6rdyQhRBiR4sXH7LZqWuK1adLSJi18KWtgszpzFl6vV98wJ3Fw824a8jMAKHJXcNm1P9A3kBAirEjx4mNVG9fgiI0HZJq08K3cnkgAGgyZvPHK73VOc2Kl67dyKKV/vUtrveyULITwKSlefKx6934c/W3SKW3SJi18Z/Gtd5PubQCgtCu4t9k/eGg7ZVFFAOTVy9unQgjfkuLFx1qaO7CbUgBI7JA2aeE7CSkW8rvrAKgP4kW7vT29dBZZ6FDiMKudXDHnfL0jCSHCjBQvPuaI8NKiaAsqz8wr0jmNCDdWZzMAVQkpOic5vtL1O6i1ahOwi3srmD77Ap0TCSHCjRQvPtaTOtAm3c41V9+icxoRbjIbWgE4EplLU329vmGOo3T9Zg4maTsBFzfbdE4jhAhHUrz4kNfjoSNJpkkL/5mWMQGj6salJPL6G0/pHWdIR2r2cSiiEIC8hmad0wghwpEULz7UUFNOy+A06Rad04hwdNkNN5HjrQWgPiFG5zTH6u7ooq04hV7FTLzqZPF1t+odSQgRhqR48aHKjWtpiu2fJt0l06SF75kiI8jv1DqOalMS9Q0zhP1rt1KVq613Gd9dQV7+OJ0TCSHCkRQvPlS7p+zradLSJi38JLdZu6tXGZemc5JjlW3cTJklG4Axjkad0wghwpUULz7U6uyi0ah1gSRLm7Twk3SHtvlhlcnK3o3BNaSxovUQR4z5AOS3SAEvhPAPKV58yGFUaTFoxcu4DKvOaUS4uvyca4lSu+lWovlw7d/0jjOow9lG27h0PIqJVK+d++/5d70jCSHClBQvPtSbqa13iVHbWXjJDTqnEeFqyuwZ5LurAWhMteic5mt7v9jMkUztraxxnZUyEkAI4TdSvPiIx+2m0/J1m3R8QoLOiUQ4y2u3A1CTnKRzkq+VbdpCWbw2PXpMo13nNEKIcCbFi4/UVRygZaB46WvVN4wIezlN/Yt2YzLwuN06p9FUeeuoNuYBMM4boXMaIUQ4k+LFR6o3fYkjpr9NWqZJCz/LausDoNaQw1/ffFHnNOBytNBaqL1llO2p5Y47fqJzIiFEOJPixUfqD5TjMGtvFck0aeFv3110FwlqKx7FxA5Hjd5x2LN6E4cz+te7tFfrnEYIEe6kePGRFlc39v426SRpkxZ+lpGTSUGvVrQ0ZOi/7qV08ybKYrUW6aIGWe8ihPAvKV58pMmo0tw/TXpsWrbOacRokOdqAqDakqxzEqiLaaHRkIGiepidUah3HCFEmJPixUd6M+JQFQPRaieXX7pI7zhiFMiytwJQac6mrUm/QaDNdXaa8tMBKPBUcfU139ctixBidJDixQf6enu+bpP22qVNWgREsTkVgEZDBsteWqpbjl2fredQmla8jHPV6pZDCDF6SPHiA7XlewfbpNNlmrQIkO/fehcZXhsANTH6bQi3e/NXlEVrbxUV1Ot3B0gIMXpI8eIDNZu+kmnSIuAMJiP53fUA2NL0W7TrzIRWJYkItZdLZ5yrWw4hxOghxYsP1JceGZwmLW3SIpCsrc0AVCWk6HL9hopabHlai/QYdwVzSi7QJYcQYnSR4sUHWtt7sZu0jg9Le7fOacRokmF3AnAkMpf9X64J+PV3ffoVh5IzAChuqQ/49YUQo5MULz7QYoJmRfuXb1Fyus5pxGhy/pT5GFU3LiWRv639a8Cvv3PHRsqi+te7NMp6LyFEYEjx4gM9GTGoigGz2snVV3xP7zhiFDnnggXkerQOH3uAJ0yrqkpnYTSdShzRaiffu2pJQK8vhBi9pHg5TT3dnXRZtMW66V6ZJi0CL69L6ziqTQnsot3a0iPU5Gh3HIt7D1OQPzag1xdCjF5SvJym2tKdg23SqdImLXRgbdb+3lXFpdHXHbg1VztXreFgUhYAxU0NAbuuEEJI8XKaqjevwxHbv0Fdl0yTFoGX2dIFQKXJyhdvvxqw6+4r3055RP96l2bZIkAIEThSvJymhkM1OPrbpJOlTVro4Karb8GsdtGjmNlkqwjINb1eL+1FFnqVKBLUVu669cGAXFcIIcDPxUtLSwtLlizBYrFgsVhYsmQJra2tw37+nXfeiaIoLF261G8ZT5ezwz3YJp0obdJCBwVFY8hzVwNgT08MyDUrdx+iMlsbTzC+q0LWegkhAsqvxcvixYvZsWMHK1euZOXKlezYsYMlS4bXkfDBBx+wceNGsrODe0Jzi1Ghqb9NujAxTec0YrTKb7cDUJ0cmEW721d+ykFLDgBjHfaAXFMIIQb4rXjZv38/K1eu5MUXX6SkpISSkhJeeOEF/vrXv1JaWnrC59bW1nLvvfeybNkyIiIi/BXRJ3oyo1EVI2a1i2uvkmm6Qh85Ta0AVMVk0lx9xO/XK206QIUxH4Cx3arfryeEEN/kt+Jl/fr1WCwWZs+ePfjYnDlzsFgsrFu37rjP83q9LFmyhJ/+9KeceeaZJ71OT08PLpfrqI9A6epoo0umSYsgYO0zAlBjyOGvy17067W8Hg9thal4FSNp3kbuufNnfr2eEEJ8m9+KF5vNRnr6sbvNpqenY7PZjvu8X//615hMJu67775hXefxxx8fXFNjsViwWq2nnHmkqg/soHWgeOmTNmmhn+8vvguL2opXMVIRpfj1Wgc37+VIVv96l85Kv15LCCGGMuLi5dFHH0VRlBN+bNmyBQBFOfaXqKqqQz4OsHXrVv73f/+XV1555bjHfNvDDz+M0+kc/Kiurh7pt3TKamWatAgSlsRE8ntrAGjISPTrtbas/ISyuFwAihodfr2WEEIMxTTSJ9x7773cfPPNJzymoKCAXbt20dBw7MZVdrudjIyMIZ+3du1aGhsbycvLG3zM4/Hwz//8zyxdupQjR44c85yoqCiioqJG9k34SOOReuxF2vv+KW1SvAh95bma2JUG1ZZkPG43RtOI//celip3LTXGiwCYFp3sl2sIIcSJjPi3W2pqKqmpqSc9rqSkBKfTyaZNm5g1axYAGzduxOl0Mnfu3CGfs2TJEi666KKjHrvkkktYsmQJP/zhD0ca1e+cnV7sJq3TKLGtR+c0YrTLtrdCGlSasylf/wXjzrnQ59fo6+2jpUDrqsvx1LB48Z0+v4YQQpyM39a8TJw4kYULF3LHHXewYcMGNmzYwB133MEVV1zB+PHjB4+bMGEC77//PgApKSlMmjTpqI+IiAgyMzOPek6waDEqOPrbpK0JifqGEaPelAztLmCjIYOVn67wyzUOrN3K4QyteBnfHri3aIUQ4pv8us/LsmXLmDx5MgsWLGDBggVMmTKF11577ahjSktLcTqd/ozhN32ZZlTFSKTazY3X/EDvOGKUu+GGH5Dh1RbDN6b5p/Nt86pVlMZoRVJhQ5NfriGEECfjnzfF+yUnJ/P666+f8BhVPfEeEUOtcwkGbc7mb0yTljZpERzyu+toiMmk3k877dbFurAb0jGoHs4rmuiXawghxMnIbKNTVLNvG60J0iYtgou1tRmAqvhUOp2tPj13T3cPTVZtvVuhp5IFF1/v0/MLIcRwSfFyimq3bsQRJ9OkRXDJdmh/Fysjc1n3lm8nTO/4+xoOpWl7NxU7a316biGEGAkpXk6RvaoBh1mmSYvgsnD+JRhVNy7FwlZ7lU/PvWH1J5RFFwJQ0NDs03MLIcRISPFyilxd3sFp0kltXTqnEUIzY8Zccj3aZnVNGb4d0tiWZcKpJBKh9nLtvIU+PbcQQoyEFC+nqNkADkV7/z83VhbriuCR36VtDlmb6rvipau9k/r+9S5j+w4zddrskzxDCCH8R4qXU+TNiMGrGIlUe2SatAgq1uZWAKpi06jbv8cn51z39kccTNF2xi5uPf5sMiGECAQpXk5Bq8NGZ5LWJp3mtZOckqZzIiG+luXqBqDSZGX18jd9cs7tO9dzMLIIgAJ7q0/OKYQQp0qKl1NQvW8LrQkxgLRJi+Bz2/fuxqx20quYqYj2+OScHflxdCkxxKjt/OCGH/nknEIIcaqkeDkFDdu2fGOatLRJi+CSnJJGvltbtOvIOP3BiS57MzU52nqXcT0VZGfmnvY5hRDidEjxcgrs1U04zNoi3dR2mSYtgk9+ux2A6uQk3H19p3Wuz19dzsHELADGNB87KV4IIQJNipdT4OoFe4T2L9pEaZMWQSjH0QpAZUwGu1d+eFrnOlB7gPIIbX+XQqf8fRdC6E+Kl1PQYvy6TTrHHK9zGiGOVWCIAqDWkMNX69ec1rk6ChPpUyJJVFu48/v3+yKeEEKcFileRkj1enGnR+NRTESovVx/zRK9IwlxjMU3/T8S1Ra8ipGm9FPfh6ixoobKbK2bblzXERlAKoQIClK8jFBTYy3dlv6ZRtImLYJUfEIC+T3a/CHbaey0+48/vU5ZQg4AY+yNPskmhBCnS4qXEarfv4XWRK14Se+T+S4ieOW5HADUJCTRYju1QYoVPTaOGPMAmOA2+iybEEKcDileRsi2bRtNsVrxktItbdIieGU5nABUmrPZ/v7yUzpHR34KqmIk3dvAnT/6iS/jCSHEKZPiZYTstc1ft0nLNGkRxGbmjQeg0ZDB9trqET//0ObtVGRqC9PHd1T6NJsQQpwOKV5GyNVnoDEiBZA2aRHcrrryZjK99QC0ZCeO+Plr/vw+ZXFWAAobHb6MJoQQp0WKlxFq/UabdHZUjM5phDix/G6teKlPS0L1ekf03HpzJ7VGbTfd6fGyMF0IETykeBkBr8dDX6p5sE36xmtv1TuSECdkbdVmb1XFp1C++asRPbc1T7vDaPVUc/OiO3yeTQghTpUULyNgrz9CT5K0SYvQkdmkLSo/EpnL1g+Hv9Pulg//Rnl6//4ubTV+ySaEEKdKipcRqNuzGaelfyCjW6ZJi+B39QVXYVTdtCkWDpuH/7ytn6+lLCYfgEJbk5/SCSHEqZHiZQRsO3fRFNd/50WmSYsQMHnSdKwe7c5JS+bwN6uzp3hwGNIwqm4umjjdX/GEEOKUSPEyAs22Vuz9bdIpMk1ahIj8Tm0SdG1qIl3tJy+6+3p6cFi1t4wKPJV85zuX+TWfEEKMlBQvI+DqM+CI0P71muiSNmkRGnKb+xftxqaz86M/n/T4L155lUNp6QAUO+v8mk0IIU6FFC8j0KIYsCvav0izIkawgEAIHeV09AFQacplx7atJz1+1+7dlJkLASi0yQgMIUTwkeJlmDxuN57MKNxKBCa1j6uvWqx3JCGG5QffvYtotZNexUxTWvxJj+/Ii8WlWIhUe7juO1cGIKEQQoyMFC/D1Npko7d/IGOaaic7M1fnREIMT3JKGvlubTyA/SQTpl0OO/W52iaMY/sOM3mSLNYVQgQfKV6GyWw005rQX7zINGkRYvLa+ydMJydRf+jAcY/79IU/cjA5E4DiFltAsgkhxEhJ8TJMvb0emuP693iRadIixOQ4WgGojM5k14fvHve4Q7Y6DkZq610K7K5ARBNCiBGT4mWYIszR2KO19QIp7TJNWoSWMRHaHK5aQzb76uuPe1xHQQLdSgyxaju33nh7oOIJIcSISPEyTHHJCdgjkgFIcnXqnEaIkfnRDx8kUW3GqxhpzR163UvN3t3U5PSPBOg5LOu6hBBBS4qXYaqz1Qy2SWeaonROI8TIFfTUAtCYnoS7r++Yr6954y0OJmYBMKapMaDZhBBiJKR4GabVq/9GLB2Y1D4uX3iD3nGEGDGrS5tRVG1JZv/qlcd8/UhPO4dNBQAUObsDGU0IIUZEipdhWnzznZRe8B3eiXaRnzdG7zhCjFi2wwlAZVQ2pV+sOeprqtdLV0ESfUokiWozdyy5T4+IQggxLFK8jFBJyYV6RxDilMwecwYAdkM6h3Af9bXdH/+NymztbdHxXUeIT0gIeD4hhBguKV6EGCUuW3gjWV5tVlFbztGLdrf9YzVl8TkAjLXbA55NCCFGQooXIUaRvG5t47n61ERc31iUWxflpdJoBWCCN1KXbEIIMVxSvAgxilhbtN2hq+JT2fH+cgB6u7toy09BVYxkeG3ccduDekYUQoiTkuJFiFEku6UdgMrIXA7s2gfAxmUvUZGpzTMa11GlWzYhhBguKV6EGEWuPP8qTGofbUoCjWnauIsDuw9SFqu9ZVTU6NAznhBCDIsUL0KMIpMnTcfqqQHAlZWE6vXSkGykzpiDonqZlZKlc0IhhDg5KV6EGGXyOhsAqE1JYu/nf6M1V3vLyOqt4fprf6hnNCGEGBYpXoQYZXKbWwGojE1j8wcfUp6u7e9S3FajYyohhBg+KV6EGGVyu7UN6qpMVg4nJFMWXQBAka1Jx1RCCDF8UrwIMcrcetOdRKud9CpR9BYk0GRIxai6WTB5lt7RhBBiWPxavLS0tLBkyRIsFgsWi4UlS5bQ2tp60uft37+fq666CovFQnx8PHPmzKGqSlo4hfCF5JQ08vuqAVhdcCYARe4jnDP/Ej1jCSHEsPm1eFm8eDE7duxg5cqVrFy5kh07drBkyZITPqe8vJz58+czYcIEVq9ezc6dO/n3f/93zGazP6MKMarktWsjACr7p0gXu+p0TCOEECNj8teJ9+/fz8qVK9mwYQOzZ88G4IUXXqCkpITS0lLGjx8/5PMeeeQRLrvsMp544onBx4qKivwVU4hRKdfhhOSvP8+3tegXRgghRshvd17Wr1+PxWIZLFwA5syZg8ViYd26dUM+x+v18re//Y1x48ZxySWXkJ6ezuzZs/nggw/8FVOIUWlsdPzgf0eq3dx80TX6hRFCiBHyW/Fis9lIT08/5vH09HRsNtuQz2lsbKS9vZ1f/epXLFy4kE8++YRrr72W6667ji+++GLI5/T09OByuY76EEKc2G233EeSqs05Ku6rYPyEqTonEkKI4Rtx8fLoo4+iKMoJP7Zs2QKAoijHPF9V1SEfB+3OC8DVV1/Ngw8+yLRp0/j5z3/OFVdcwbPPPjvkcx5//PHBBcEWiwWr1TrSb0mIUWlMt7YIfnyTrHcRQoSWEa95uffee7n55ptPeExBQQG7du2ioaHhmK/Z7XYyMjKGfF5qaiomk4kzzjjjqMcnTpzIl19+OeRzHn74YR566KHBz10ulxQwQgzDdQ4XY1jFffOv1DuKEEKMyIiLl9TUVFJTU096XElJCU6nk02bNjFrlrZ/xMaNG3E6ncydO3fI50RGRnL22WdTWlp61ONlZWXk5+cP+ZyoqCiioqJG+F0IIW675T69IwghxCnx25qXiRMnsnDhQu644w42bNjAhg0buOOOO7jiiiuO6jSaMGEC77///uDnP/3pT1m+fDkvvPAChw4d4qmnnuKjjz7i7rvv9ldUIYQQQoQQv+7zsmzZMiZPnsyCBQtYsGABU6ZM4bXXXjvqmNLSUpxO5+Dn1157Lc8++yxPPPEEkydP5sUXX+Tdd99l/vz5/owqhBBCiBChqKqq6h3Cl1wuFxaLBafTSUJCgt5xhBBCCDEMI3n9ltlGQgghhAgpUrwIIYQQIqRI8SKEEEKIkCLFixBCCCFCihQvQgghhAgpUrwIIYQQIqRI8SKEEEKIkCLFixBCCCFCihQvQgghhAgpUrwIIYQQIqSMeKp0sBuYduByuXROIoQQQojhGnjdHs7UorArXtra2gCwWq06JxFCCCHESLW1tWGxWE54TNgNZvR6vdTV1REfH4+iKD49t8vlwmq1Ul1dLUMf/Uh+zoEhP+fAkZ91YMjPOTD89XNWVZW2tjays7MxGE68qiXs7rwYDAZyc3P9eo2EhAT5HyMA5OccGPJzDhz5WQeG/JwDwx8/55PdcRkgC3aFEEIIEVKkeBFCCCFESJHiZQSioqL4xS9+QVRUlN5Rwpr8nANDfs6BIz/rwJCfc2AEw8857BbsCiGEECK8yZ0XIYQQQoQUKV6EEEIIEVKkeBFCCCFESJHiRQghhBAhRYqXYXr66acpLCzEbDYzY8YM1q5dq3eksPP4449z9tlnEx8fT3p6Otdccw2lpaV6xwp7jz/+OIqi8MADD+gdJezU1tby/e9/n5SUFGJiYpg2bRpbt27VO1ZYcbvd/Nu//RuFhYVER0dTVFTEY489htfr1TtayFuzZg1XXnkl2dnZKIrCBx98cNTXVVXl0UcfJTs7m+joaL7zne+wd+/egGST4mUYli9fzgMPPMAjjzzC9u3bOeecc7j00kupqqrSO1pY+eKLL7jnnnvYsGEDq1atwu12s2DBAjo6OvSOFrY2b97M888/z5QpU/SOEnZaWlqYN28eERER/P3vf2ffvn389re/JTExUe9oYeXXv/41zz77LE899RT79+/niSee4H/+53/4v//7P72jhbyOjg6mTp3KU089NeTXn3jiCZ588kmeeuopNm/eTGZmJhdffPHgjEG/UsVJzZo1S73rrruOemzChAnqz3/+c50SjQ6NjY0qoH7xxRd6RwlLbW1tanFxsbpq1Sr1vPPOU++//369I4WVn/3sZ+r8+fP1jhH2Lr/8cvW222476rHrrrtO/f73v69TovAEqO+///7g516vV83MzFR/9atfDT7W3d2tWiwW9dlnn/V7HrnzchK9vb1s3bqVBQsWHPX4ggULWLdunU6pRgen0wlAcnKyzknC0z333MPll1/ORRddpHeUsPThhx8yc+ZMbrzxRtLT0znrrLN44YUX9I4VdubPn8+nn35KWVkZADt37uTLL7/ksssu0zlZeKuoqMBmsx312hgVFcV5550XkNfGsBvM6GsOhwOPx0NGRsZRj2dkZGCz2XRKFf5UVeWhhx5i/vz5TJo0Se84Yeett95i27ZtbN68We8oYevw4cM888wzPPTQQ/zrv/4rmzZt4r777iMqKopbbrlF73hh42c/+xlOp5MJEyZgNBrxeDz893//N9/97nf1jhbWBl7/hnptrKys9Pv1pXgZJkVRjvpcVdVjHhO+c++997Jr1y6+/PJLvaOEnerqau6//34++eQTzGaz3nHCltfrZebMmfzyl78E4KyzzmLv3r0888wzUrz40PLly3n99dd54403OPPMM9mxYwcPPPAA2dnZ3HrrrXrHC3t6vTZK8XISqampGI3GY+6yNDY2HlNxCt/48Y9/zIcffsiaNWvIzc3VO07Y2bp1K42NjcyYMWPwMY/Hw5o1a3jqqafo6enBaDTqmDA8ZGVlccYZZxz12MSJE3n33Xd1ShSefvrTn/Lzn/+cm2++GYDJkydTWVnJ448/LsWLH2VmZgLaHZisrKzBxwP12ihrXk4iMjKSGTNmsGrVqqMeX7VqFXPnztUpVXhSVZV7772X9957j88++4zCwkK9I4WlCy+8kN27d7Njx47Bj5kzZ/K9732PHTt2SOHiI/PmzTum1b+srIz8/HydEoWnzs5ODIajX8qMRqO0SvtZYWEhmZmZR7029vb28sUXXwTktVHuvAzDQw89xJIlS5g5cyYlJSU8//zzVFVVcdddd+kdLazcc889vPHGG/zlL38hPj5+8G6XxWIhOjpa53ThIz4+/ph1RLGxsaSkpMj6Ih968MEHmTt3Lr/85S+56aab2LRpE88//zzPP/+83tHCypVXXsl///d/k5eXx5lnnsn27dt58sknue222/SOFvLa29s5dOjQ4OcVFRXs2LGD5ORk8vLyeOCBB/jlL39JcXExxcXF/PKXvyQmJobFixf7P5zf+5nCxB/+8Ac1Pz9fjYyMVKdPny7tu34ADPnx8ssv6x0t7EmrtH989NFH6qRJk9SoqCh1woQJ6vPPP693pLDjcrnU+++/X83Ly1PNZrNaVFSkPvLII2pPT4/e0ULe559/PuTv5FtvvVVVVa1d+he/+IWamZmpRkVFqeeee666e/fugGRTVFVV/V8iCSGEEEL4hqx5EUIIIURIkeJFCCGEECFFihchhBBChBQpXoQQQggRUqR4EUIIIURIkeJFCCGEECFFihchhBBChBQpXoQQQggRUqR4EUIIIURIkeJFCCGEECFFihchhBBChBQpXoQQQggRUv5/dPyEzeL+UaYAAAAASUVORK5CYII=",
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"for chrg in charges:\n",
" plt.plot(chrg)"
]
},
{
"cell_type": "markdown",
"id": "be89b420",
"metadata": {},
"source": [
"further, we can plot the charge on each atom as a function of distance:"
]
},
{
"cell_type": "code",
"execution_count": 16,
"id": "1969d9ee",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"<matplotlib.legend.Legend at 0x135113ce0>"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAkMAAAGwCAYAAACq12GxAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjAsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvlHJYcgAAAAlwSFlzAAAPYQAAD2EBqD+naQAAQUdJREFUeJzt3XtclHXe//H3gAIeYDwgnhI0MQ+peSANrcxWPGSl1b15ytRcu63cTHPTfmZopVS7lrWt6XrIMlOrNXO7izTLdFU0D6xumHnGDDwHKIkK1+8PlomR08www8xwvZ6PxzyKa77XdX2uL9+LeXudxmIYhiEAAACTCvB2AQAAAN5EGAIAAKZGGAIAAKZGGAIAAKZGGAIAAKZGGAIAAKZGGAIAAKZWxdsF+Lq8vDz9/PPPCg0NlcVi8XY5AADAAYZhKCsrS40aNVJAQOnHfghDZfj555/VpEkTb5cBAABccPz4cV133XWltiEMlSE0NFRSfmeGhYV5uRoAAOCIzMxMNWnSxPY5XhrCUBkKTo2FhYURhgAA8DOOXOLCBdQAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUeAI1SpSbZ2j7kXM6lXVJEaEh6tKsjgIDyn6SpyvzubouR3l6+e5QVo3e3AZ/6D9X+ON2Fa45vGawZEhnLuaUu35/6YuS6vSX+vEbX/qdEYZQrMT/pGnGP1OUlnHJNq2hNUTx97RR37YN3Tqfq+vy9LZUpLJq9OY2+EP/ucIft6u4mgtztX5/6YuS6rz3poZa8+80n68fv/G1MWcxDMOo8LX6kczMTFmtVmVkZJjmu8kS/5Omx97fpWsHRkFef/uhTsUOVlfmc3VdjvL08t2hrBofvb2Z/r7xiFe2wR/6zxX+uF0l1VyYK/X7S184sv2F+Vr9+E1FjTlnPr+5Zgh2cvMMzfhnSrF/cAqmzfhninLz7Fu4Mp+r63KUp5fvDmXVaEhasKloECp4X/LcNvhD/7nCH7ertJoLc7Z+f+kLR7e/MF+qH7/x1TFHGIKd7UfOlXgIXsofrGkZl7T9yLlyz+fquhzl6eW7Q1k1SlJpfxM8uQ3+0H+u8MftcmScFHCmfn/pC2e2vzBfqR+/8dUxRxiCnVNZjv3BubadK/O5ui5HeXr57uCudXtiG/yh/1zhj9vlSi2OzOMvfVHe9Xu7fvzGV8ccYQh2IkJDXGrnynyurstRnl6+O7hr3Z7YBn/oP1f443a5Uosj8/hLX5R3/d6uH7/x1TFHGIKdLs3qqKE1RCXd3GhR/hX/XZrVKfd8Zc2jEtblKFe3pSI50gcBFnllG/yh/1zhj9vlyDgp4Ez9jiy3bo0gdY6q7WipHuHM9hfmi79Ls/PV/c/vwtDcuXPVrFkzhYSEqHPnztq0aZND823evFlVqlRRhw4dPFugg3LzDG09dFafJp/Q1kNnfeYCv8AAi+LvaSOp6Adwwc/x97Qp8iwIV+YrbZ4Cv17J1bqUdOc2ohw1VbSyarRIGnNbsxLfl9y7DYXH5fYj5zStv/P956tju4Az48JXtsWRfaXwe46OCUeWe/biZfX48zdK/E+a4wW7maPbX5iv7OOw56t/l/3q1vqVK1dq+PDhmjt3rrp376758+dr4cKFSklJUWRkZInzZWRkqFOnToqOjtbJkyeVnJzs8Do9cWu9rz1foTgV/ZyhKav26pfsK0Xec8etlpWhvytiG9zxDBd/6OsCvtDnzqrI5wwV5iu3qfOcocqjIvYvZz6//SoMde3aVZ06ddLbb79tm9a6dWsNHDhQCQkJJc43ePBgtWjRQoGBgVq9erVXw5C/PNNDqrgnUOfmGer+8nqlZ+YU+75FUgNriP41+c5K/XRdbz6Buqxx+behHVW7RnCp6/ansV2gpD715W3x1BOoL1/N0y0JX+ncxaL/KJHcsx+6A0+grjw8/Ttz5vPbb55AffnyZe3cuVNTpkyxm967d29t2bKlxPneeecdHTp0SO+//75eeumlMteTk5OjnJzfPpQzMzNdL/oaZT1fwaL85yvEtWngEztxYIBFsc3reny+7UfOlRiEJPtbLV2px5WavKGsGj21DY6Myxf/b1+pH4L+NrYLFNenvr4tnhoHO4+dLzEISe7ZD92hpO33h30c9nzpd+Y31wydOXNGubm5ql+/vt30+vXrKz29+GtKDhw4oClTpmjZsmWqUsWx3JeQkCCr1Wp7NWnSpNy1F/DV5yt4m6/eamkW7hiXlWlsV6ZtcQb7IczMb8JQAYvF/l9ihmEUmSZJubm5Gjp0qGbMmKEbbrjB4eU/++yzysjIsL2OHz9e7poL8MemeL56q6VZuGNcVqaxXZm2xRnshzAzvzlNFh4ersDAwCJHgU6dOlXkaJEkZWVlaceOHdq9e7fGjRsnScrLy5NhGKpSpYrWrl2rO++8s8h8wcHBCg4O9sg28MemeAW3WqZnXCr21ETBtQrcHusZ7hiXlWlsV6ZtcQb7IczMb44MBQUFqXPnzlq3bp3d9HXr1qlbt25F2oeFhWnv3r1KTk62vcaOHauWLVsqOTlZXbt2rajSbXz1+Qre5qu3WpqFO8ZlZRrblWlbnMF+CDPzmzAkSRMnTtTChQu1ePFi7du3TxMmTFBqaqrGjh0rKf8U18MPPyxJCggIUNu2be1eERERCgkJUdu2bVWjRo0Kr58/NiXr27ah3n6okxpY7f+13cAa4pN3IVUm7hiXlWlsV6ZtcRb7IczKb06TSdKgQYN09uxZvfDCC0pLS1Pbtm31+eefKyoqSpKUlpam1NRUL1dZuoI/Ntc+X6EBz8RQ37YNFdemAbfHeoE7xmVlGtuVaVucxX4IM/Kr5wx5gyceuij5x3NvYD7uGJeVaWxXpm0BzKbSPnTRGzwVhgAAgOc48/ntV9cMAQAAuBthCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmBphCAAAmJrfhaG5c+eqWbNmCgkJUefOnbVp06YS265atUpxcXGqV6+ewsLCFBsbqy+//LICqwUAAL7Or8LQypUr9dRTT2nq1KnavXu3brvtNvXr10+pqanFtt+4caPi4uL0+eefa+fOnerZs6fuuece7d69u4IrBwAAvspiGIbh7SIc1bVrV3Xq1Elvv/22bVrr1q01cOBAJSQkOLSMG2+8UYMGDdLzzz/vUPvMzExZrVZlZGQoLCzMpboBAEDFcubz22+ODF2+fFk7d+5U79697ab37t1bW7ZscWgZeXl5ysrKUp06dUpsk5OTo8zMTLsXAACovPwmDJ05c0a5ubmqX7++3fT69esrPT3doWXMnj1bFy9e1IMPPlhim4SEBFmtVturSZMm5aobAAD4Nr8JQwUsFovdz4ZhFJlWnOXLl2v69OlauXKlIiIiSmz37LPPKiMjw/Y6fvx4uWsGAAC+q4q3C3BUeHi4AgMDixwFOnXqVJGjRddauXKlRo8erY8++ki9evUqtW1wcLCCg4PLXS8AAPAPfnNkKCgoSJ07d9a6devspq9bt07dunUrcb7ly5dr5MiR+uCDD9S/f39PlwkAAPyM3xwZkqSJEydq+PDhiomJUWxsrP7+978rNTVVY8eOlZR/iuvEiRN67733JOUHoYcfflhvvPGGbrnlFttRpWrVqslqtXptOwAAgO/wqzA0aNAgnT17Vi+88ILS0tLUtm1bff7554qKipIkpaWl2T1zaP78+bp69aqeeOIJPfHEE7bpI0aM0JIlSyq6fAAA4IP86jlD3sBzhgAA8D+V8jlDAAAAnkAYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApkYYAgAApuZSGFq6dKm6d++uRo0a6dixY5KkOXPm6NNPP3VrcQAAAJ7mdBh6++23NXHiRN1111365ZdflJubK0mqVauW5syZ4+76AAAAPMrpMPTXv/5VCxYs0NSpUxUYGGibHhMTo71797q1OAAAAE9zOgwdOXJEHTt2LDI9ODhYFy9edEtRAAAAFcXpMNSsWTMlJycXmf7FF1+oTZs27qgJAACgwlRxdoY//elPeuKJJ3Tp0iUZhqHt27dr+fLlSkhI0MKFCz1RIwAAgMc4HYZGjRqlq1ev6plnnlF2draGDh2qxo0b64033tDgwYM9USMAAIDHWAzDMFyd+cyZM8rLy1NERIQ7a/IpmZmZslqtysjIUFhYmLfLAQAADnDm89vpI0OFhYeHl2d2AAAAr+MJ1AAAwNQIQwAAwNQIQwAAwNQIQwAAwNTcGobee+89HTp0yJ2LBAAA8Ci3hqGRI0eqTZs2+uMf/+jOxQIAAHiMW8NQXl6e9u/fr7Zt27pzsQAAAB5TrocumgEPXQQAwP848/nt0pGhQ4cO6bnnntOQIUN06tQpSVJiYqK+//57VxYHAADgNU6HoW+//Vbt2rXTtm3btGrVKl24cEGStGfPHsXHx7u9QAAAAE9yOgxNmTJFL730ktatW6egoCDb9J49e2rr1q1uLQ4AAMDTnA5De/fu1X333Vdker169XT27Fm3FAUAAFBRnA5DtWrVUlpaWpHpu3fvVuPGjd1SFAAAQEVxOgwNHTpUkydPVnp6uiwWi/Ly8rR582ZNmjRJDz/8sCdqBAAA8Binw9DMmTMVGRmpxo0b68KFC2rTpo1uv/12devWTc8995wnagQAAPAYl58zdOjQIe3evVt5eXnq2LGjWrRo4e7afALPGQIAwP848/ldxdWVNG/eXM2bN3d1dgAAAJ/gdBiaOHFisdMtFotCQkIUHR2tAQMGqE6dOuUuDgAAwNOcPk3Ws2dP7dq1S7m5uWrZsqUMw9CBAwcUGBioVq1aaf/+/bJYLPrXv/6lNm3aeKruCsNpMgAA/I9Hv45jwIAB6tWrl37++Wft3LlTu3bt0okTJxQXF6chQ4boxIkTuv322zVhwgSXNwAAAKCiOB2G/vznP+vFF1+0S1lhYWGaPn26Xn31VVWvXl3PP/+8du7c6dZCC8ydO1fNmjVTSEiIOnfurE2bNpXa/ttvv1Xnzp0VEhKi66+/XvPmzfNIXQAAwD85fc1QRkaGTp06VeQU2OnTp5WZmSkp/8GMly9fdk+FhaxcuVJPPfWU5s6dq+7du2v+/Pnq16+fUlJSFBkZWaT9kSNHdNddd2nMmDF6//33tXnzZj3++OOqV6+eHnjgAbfXBwDwT7m5ubpy5Yq3y4ATqlatqsDAQLcsy+lrhoYNG6atW7dq9uzZuvnmm2WxWLR9+3ZNmjRJ3bp109KlS7VixQr95S9/0Y4dO9xSZIGuXbuqU6dOevvtt23TWrdurYEDByohIaFI+8mTJ2vNmjXat2+fbdrYsWP173//2+HvUeOaIQCovAzDUHp6un755RdvlwIX1KpVSw0aNJDFYinynkdvrZ8/f74mTJigwYMH6+rVq/kLqVJFI0aM0Ouvvy5JatWqlRYuXOjsokt1+fJl7dy5U1OmTLGb3rt3b23ZsqXYebZu3arevXvbTevTp48WLVqkK1euqGrVqkXmycnJUU5Oju3ngqNdAIDKpyAIRUREqHr16sV+qML3GIah7OxsnTp1SpLUsGHDci3PqTCUm5urnTt36pVXXtHrr7+uw4cPyzAMNW/eXDVr1rS169ChQ7mKKs6ZM2eUm5ur+vXr202vX7++0tPTi50nPT292PZXr17VmTNniu28hIQEzZgxw32FAwB8Um5uri0I1a1b19vlwEnVqlWTJJ06dUoRERHlOmXm1AXUgYGB6tOnjzIyMlSzZk21b99eN910k10Q8rRrU7thGKUm+eLaFze9wLPPPquMjAzb6/jx4+WsGADgiwquEapevbqXK4GrCn535b3ey+nTZO3atdPhw4fVrFmzcq3YWeHh4QoMDCxyFOjUqVNFjv4UaNCgQbHtq1SpUuK/AoKDgxUcHOyeogEAPo9TY/7LXb87l76oddKkSfrss8+UlpamzMxMu5enBAUFqXPnzlq3bp3d9HXr1qlbt27FzhMbG1uk/dq1axUTE1Ps9UIAAMB8nA5Dffv21b///W/de++9uu6661S7dm3Vrl1btWrVUu3atT1Ro83EiRO1cOFCLV68WPv27dOECROUmpqqsWPHSso/xfXwww/b2o8dO1bHjh3TxIkTtW/fPi1evFiLFi3SpEmTPFonAAD+ZsmSJapVq5a3y/AKp0+TffPNN56owyGDBg3S2bNn9cILLygtLU1t27bV559/rqioKElSWlqaUlNTbe2bNWumzz//XBMmTNDf/vY3NWrUSG+++SbPGAIA+K2RI0fql19+0erVq+2mb9iwQT179tT58+dNG2pc5XQY6tGjhyfqcNjjjz+uxx9/vNj3lixZUmRajx49tGvXLg9XBQAwq9w8Q9uPnNOprEuKCA1Rl2Z1FBjAdUgFSnqUjS9x+jRZgezsbP3www/as2eP3QsAALNI/E+abn3law1ZkKTxK5I1ZEGSbn3layX+J83bpRXrl19+0aOPPqr69esrJCREbdu21WeffWbX5ssvv1Tr1q1Vs2ZN9e3bV2lpv23Ld999p7i4OIWHh8tqtRZ7wMFisWjevHkaMGCAatSooZdeekmS9NJLLykiIkKhoaH6wx/+oClTphR5FM8777yj1q1bKyQkRK1atdLcuXM90xHXcDoMnT59WnfffbdCQ0N14403qmPHjnYvAADMIPE/aXrs/V1Ky7hkNz0945Iee3+XzwWivLw89evXT1u2bNH777+vlJQUvfzyy3bP58nOztZf/vIXLV26VBs3blRqaqrddbZZWVkaMWKENm3apKSkJLVo0UJ33XWXsrKy7NYVHx+vAQMGaO/evXrkkUe0bNkyzZw5U6+88op27typyMhIu2+TkKQFCxZo6tSpmjlzpvbt26dZs2Zp2rRpevfddz3bMXLhNNlTTz2l8+fPKykpST179tQnn3yikydP6qWXXtLs2bM9USMAAD4lN8/QjH+mqLjvszIkWSTN+GeK4to08Mgps88++6zIM/5yc3NLneerr77S9u3btW/fPt1www2SpOuvv96uzZUrVzRv3jw1b95ckjRu3Di98MILtvfvvPNOu/bz589X7dq19e233+ruu++2TR86dKgeeeQR28+DBg3S6NGjNWrUKEnS888/r7Vr1+rChQu2Ni+++KJmz56t+++/X1L+db8pKSmaP3++RowYUXqHlJPTYejrr7/Wp59+qptvvlkBAQGKiopSXFycwsLClJCQoP79+3uiTgAAfMb2I+eKHBEqzJCUlnFJ24+cU2xz9z/dumfPnkWOrGzbtk0PPfRQifMkJyfruuuuswWh4lSvXt0WhKT8r7ko+MoLKf9Zfc8//7y+/vprnTx5Urm5ucrOzra7eUmSYmJi7H7ev39/ket9u3Tpoq+//lpS/lmn48ePa/To0RozZoytzdWrV2W1Wkus112cDkMXL15URESEJKlOnTo6ffq0brjhBrVr144LlQEApnAqq+Qg5Eo7Z9WoUUPR0dF203766adS5yn4+orSXHuhs8ViUeHvcx85cqROnz6tOXPmKCoqSsHBwYqNjdXly5eL1Hetkr4RQso/hSflnyrr2rWrXTt3fTN9aZy+Zqhly5bav3+/pPzvIJs/f75OnDihefPmlfuL0gAA8AcRoSFubVcR2rdvr59++kk//vijy8vYtGmTnnzySd1111268cYbFRwcrDNnzpQ5X8uWLbV9+3a7aTt27LD9f/369dW4cWMdPnxY0dHRdq+K+MYLl64ZKriyPD4+Xn369NGyZcsUFBRU7K3tAABUNl2a1VFDa4jSMy4Ve92QRVIDa/5t9r6iR48euv322/XAAw/otddeU3R0tH744QdZLBb17dvXoWVER0dr6dKliomJUWZmpv70pz85dMTpj3/8o8aMGaOYmBh169ZNK1eu1J49e+yuWZo+fbqefPJJhYWFqV+/fsrJydGOHTt0/vx5TZw40eXtdoTTR4aGDRumkSNHSpI6duyoo0eP6rvvvtPx48c1aNAgd9cHAIDPCQywKP6eNpLyg09hBT/H39PG55439I9//EM333yzhgwZojZt2uiZZ54p88LrwhYvXqzz58+rY8eOGj58uJ588knbpTOlGTZsmJ599llNmjRJnTp10pEjRzRy5EiFhPx25OwPf/iDFi5cqCVLlqhdu3bq0aOHlixZUiFHhixG4ZN2KCIzM1NWq1UZGRkKCwvzdjkAADe5dOmSjhw5ombNmtl9KDsj8T9pmvHPFLuLqRtaQxR/Txv1bculI6WJi4tTgwYNtHTpUpeXUdrv0JnPb6dPk+Xm5mrJkiVav369Tp06ZbvoqUDBleEAAFR2fds2VFybBjyBugzZ2dmaN2+e+vTpo8DAQC1fvlxfffVVkS9T9xanw9D48eO1ZMkS9e/fX23bti1ydTgAAGYSGGDxyO3zlYnFYtHnn3+ul156STk5OWrZsqX+8Y9/qFevXt4uTZILYWjFihX68MMPddddd3miHgAAUMlUq1ZNX331lbfLKJHTF1AHBQUVebYBAACAv3I6DD399NN64403xHXXAACgMnDoNFnB94QU+Prrr/XFF1/oxhtvLPK0ylWrVrmvOgAAAA9zKAxd+70g9913n0eKAQAAqGgOhaF33nnH03UAAAB4hdPXDB05ckQHDhwoMv3AgQM6evSoO2oCAACoME6HoZEjR2rLli1Fpm/bts32NR0AAMD/NG3aVHPmzPF2GRXO6TC0e/dude/evcj0W265RcnJye6oCQAAlOL48eMaPXq0GjVqpKCgIEVFRWn8+PE6e/ast0vzS06HIYvFoqysrCLTMzIynPqyNwAAKoW8XOnIJmnvx/n/zfPsZ+Hhw4cVExOjH3/8UcuXL9fBgwc1b948rV+/XrGxsTp37pxH1++sK1eueLuEMjkdhm677TYlJCTYBZ/c3FwlJCTo1ltvdWtxAAD4tJQ10py20rt3S/8Ynf/fOW3zp3vIE088oaCgIK1du1Y9evRQZGSk+vXrp6+++konTpzQ1KlTS51/zZo1iomJUUhIiMLDw4s8Pic7O1uPPPKIQkNDFRkZqb///e9270+ePFk33HCDqlevruuvv17Tpk2zCzzTp09Xhw4dtHjxYl1//fUKDg6WYRj64YcfdOuttyokJERt2rTRV199JYvFotWrV9vmPXHihAYNGqTatWurbt26GjBgQIVcj+x0GHr11Vf19ddfq2XLlho1apRGjRqlli1bauPGjfrzn//siRoBAPA9KWukDx+WMn+2n56Zlj/dA4Ho3Llz+vLLL/X444+rWrVqdu81aNBAw4YN08qVK0t8MPL//d//6f7771f//v21e/durV+/XjExMXZtZs+erZiYGO3evVuPP/64HnvsMf3www+290NDQ7VkyRKlpKTojTfe0IIFC/T666/bLePgwYP68MMP9Y9//EPJycnKy8vTwIEDVb16dW3btk1///vfi4S27Oxs9ezZUzVr1tTGjRv1r3/9SzVr1lTfvn11+fLl8nRbmZz+brI2bdpoz549euutt/Tvf/9b1apV08MPP6xx48apTp06nqgRAADfkpcrJU6WVFzoMCRZpMQpUqv+UkCg21Z74MABGYah1q1bF/t+69atdf78eZ0+fVoRERFF3p85c6YGDx6sGTNm2KbddNNNdm3uuusuPf7445LyjwK9/vrr2rBhg1q1aiVJeu6552xtmzZtqqefflorV67UM888Y5t++fJlLV26VPXq1ZMkJSYm6tChQ9qwYYMaNGhgqyUuLs42z4oVKxQQEKCFCxfavgT+nXfeUa1atbRhwwb17t3b8Y5yktNhSJIaNWqkWbNmubsWAAD8w7EtRY8I2TGkzBP57ZrdVmFlFRwRKggT10pOTtaYMWNKXUb79u1t/2+xWNSgQQOdOnXKNu3jjz/WnDlzdPDgQV24cEFXr15VWFiY3TKioqJsQUiS9u/fryZNmtiCkCR16dLFbp6dO3fq4MGDCg0NtZt+6dIlHTp0qNSay8ulMAQAgKldOOnedg6Kjo6WxWJRSkqKBg4cWOT9H374QbVr11Z4eHix8197aq04137NlsViUV5eniQpKSnJdmSpT58+slqtWrFihWbPnm03T40aNex+NgyjxIBWIC8vT507d9ayZcuKvFc4WHmC09cMAQBgejXru7edg+rWrau4uDjNnTtXv/76q9176enpWrZsmQYNGlRi8Gjfvr3Wr1/v8vo3b96sqKgoTZ06VTExMWrRooWOHTtW5nytWrVSamqqTp78LRx+9913dm06deqkAwcOKCIiQtHR0Xava78WzN0IQwAAOCuqmxTWSFJJRzssUljj/HZu9tZbbyknJ0d9+vTRxo0bdfz4cSUmJiouLk6NGzfWzJkzS5w3Pj5ey5cvV3x8vPbt26e9e/fq1VdfdXjd0dHRSk1N1YoVK3To0CG9+eab+uSTT8qcLy4uTs2bN9eIESO0Z88ebd682XYBdUFwGzZsmMLDwzVgwABt2rRJR44c0bfffqvx48frp59+crhGVxCGAABwVkCg1PeV//5wbSD67899X3brxdMFWrRooR07dqh58+YaNGiQmjdvrkcffVQ9e/bU1q1bS72Z6Y477tBHH32kNWvWqEOHDrrzzju1bds2h9c9YMAATZgwQePGjVOHDh20ZcsWTZs2rcz5AgMDtXr1al24cEE333yz/vCHP9guxA4JCZEkVa9eXRs3blRkZKTuv/9+tW7dWo888oh+/fXXItckuZvFKOn+O0iSMjMzZbValZGR4fFfBgCg4ly6dElHjhxRs2bNbB/ITktZk39XWeGLqcMa5wehNve6p9BKavPmzbr11lt18OBBNW/e3KVllPY7dObz26ELqDt27FjmhU8Fdu3a5VA7AAD8Xpt782+fP7Yl/2LpmvXzT4154IiQv/vkk09Us2ZNtWjRQgcPHtT48ePVvXt3l4OQOzkUhoq7Yh0AACg/+FTg7fP+KisrS88884yOHz+u8PBw9erVq8hdaN7CabIycJoMACont5wmg1e56zQZF1ADAABTc/qhi7m5uXr99df14YcfKjU1tcj3hfjat+UCAFAaTpD4L3f97pw+MjRjxgy99tprevDBB5WRkaGJEyfq/vvvV0BAgKZPn+6WogAA8LSCJy1nZ2d7uRK4quB3d+1Ts53l9JGhZcuWacGCBerfv79mzJihIUOGqHnz5mrfvr2SkpL05JNPlqsgAAAqQmBgoGrVqmX73q3q1as7fOc0vMswDGVnZ+vUqVOqVauWAgPLd/ee02EoPT1d7dq1kyTVrFlTGRkZkqS7777boQcvAQDgKwq+OLTwF5HCf9SqVcvuy19d5XQYuu6665SWlqbIyEhFR0dr7dq16tSpk7777jsFBweXuyAAACqKxWJRw4YNFRERoStXrni7HDihatWq5T4iVMDpMHTfffdp/fr16tq1q8aPH68hQ4Zo0aJFSk1N1YQJE9xSFAAAFSkwMNBtH6zwP+V+zlBSUpK2bNmi6Oho3Xtv5Xv0OM8ZAgDA/7j96zhKc8stt+iWW24p72IAAAC8wqEwtGbNGvXr109Vq1bVmjVrSm1bGY8OAQCAysuh02QBAQFKT09XRESEAgJKfjSRxWJRbm6uWwv0Nk6TAQDgf9x+miwvL6/Y/wcAAPB3Tj+B+r333lNOTk6R6ZcvX9Z7773nlqIAAAAqitN3kwUGBiotLU0RERF208+ePauIiAhOkwEAAK/z6LfWG4ZR7OPKf/rpJ1mtVmcXBwAA4FUO31rfsWNHWSwWWSwW/e53v1OVKr/NmpubqyNHjqhv374eKRIAAMBTHA5DAwcOlCQlJyerT58+qlmzpu29oKAgNW3aVA888IDbCwQAAPAkh8NQfHy8cnNzFRUVpT59+qhhw4aerAsAAKBCOHXNUGBgoMaOHatLly55qh4AAIAK5fQF1O3atdPhw4c9UQsAAECFczoMzZw5U5MmTdJnn32mtLQ0ZWZm2r0AAAD8idPPGSr8dRyFb7EvuOWe5wwBAABv8+i31n/zzTcuFwYAAOBrnA5DPXr08EQdAAAAXuF0GCqQnZ2t1NRUXb582W56+/bty10UAABARXE6DJ0+fVqjRo3SF198Uez7le2aIQAAULk5fTfZU089pfPnzyspKUnVqlVTYmKi3n33XbVo0UJr1qzxRI0AAAAe43QY+vrrr/X666/r5ptvVkBAgKKiovTQQw/p1VdfVUJCgidqlCSdP39ew4cPl9VqldVq1fDhw/XLL7+U2P7KlSuaPHmy2rVrpxo1aqhRo0Z6+OGH9fPPP3usRgAA4H+cDkMXL15URESEJKlOnTo6ffq0pPyHMe7atcu91RUydOhQJScnKzExUYmJiUpOTtbw4cNLbJ+dna1du3Zp2rRp2rVrl1atWqUff/xR9957r8dqBAAA/sfpa4Zatmyp/fv3q2nTpurQoYPmz5+vpk2bat68eR77vrJ9+/YpMTFRSUlJ6tq1qyRpwYIFio2N1f79+9WyZcsi81itVq1bt85u2l//+ld16dJFqampioyMLHZdOTk5ysnJsf3MgyQBAKjcXLpmqOBUU3x8vBITExUZGak333xTs2bNcnuBkrR161ZZrVZbEJKkW265RVarVVu2bHF4ORkZGbJYLKpVq1aJbRISEmyn4qxWq5o0aVKe0gEAgI9z+sjQsGHDbP/fsWNHHT16VD/88IMiIyMVHh7u1uIKpKen207NFRYREaH09HSHlnHp0iVNmTJFQ4cOLfVJlM8++6wmTpxo+zkzM5NABABAJebwkaHs7Gw98cQTaty4sSIiIjR06FCdOXNG1atXV6dOnVwKQtOnT5fFYin1tWPHDkn2X/1RoOArQMpy5coVDR48WHl5eZo7d26pbYODgxUWFmb3AgAAlZfDR4bi4+O1ZMkSDRs2TCEhIVq+fLkee+wxffTRRy6vfNy4cRo8eHCpbZo2bao9e/bo5MmTRd47ffq06tevX+r8V65c0YMPPqgjR47o66+/JtwAAAA7DoehVatWadGiRbbw8tBDD6l79+7Kzc1VYGCgSysPDw936IhSbGysMjIytH37dnXp0kWStG3bNmVkZKhbt24lzlcQhA4cOKBvvvlGdevWdalOAABQeTl8muz48eO67bbbbD936dJFVapUqZDn9rRu3Vp9+/bVmDFjlJSUpKSkJI0ZM0Z333233Z1krVq10ieffCJJunr1qv7nf/5HO3bs0LJly5Sbm6v09HSlp6cX+QoRAABgXg6HodzcXAUFBdlNq1Kliq5ever2ooqzbNkytWvXTr1791bv3r3Vvn17LV261K7N/v37lZGRIUn66aeftGbNGv3000/q0KGDGjZsaHs5cwcaAACo3CyGYRiONAwICFC/fv0UHBxsm/bPf/5Td955p2rUqGGbtmrVKvdX6UWZmZmyWq3KyMjgeiMAAPyEM5/fDl8zNGLEiCLTHnroIeerAwAA8CEOh6F33nnHk3UAAAB4hdNPoAYAAKhMCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDUCEMAAMDU/CYMnT9/XsOHD5fVapXVatXw4cP1yy+/ODz///7v/8pisWjOnDkeqxEAAPgfvwlDQ4cOVXJyshITE5WYmKjk5GQNHz7coXlXr16tbdu2qVGjRh6uEgAA+Jsq3i7AEfv27VNiYqKSkpLUtWtXSdKCBQsUGxur/fv3q2XLliXOe+LECY0bN05ffvml+vfvX+a6cnJylJOTY/s5MzOz/BsAAAB8ll8cGdq6dausVqstCEnSLbfcIqvVqi1btpQ4X15enoYPH64//elPuvHGGx1aV0JCgu1UnNVqVZMmTcpdPwAA8F1+EYbS09MVERFRZHpERITS09NLnO+VV15RlSpV9OSTTzq8rmeffVYZGRm21/Hjx12qGQAA+AevhqHp06fLYrGU+tqxY4ckyWKxFJnfMIxip0vSzp079cYbb2jJkiUltilOcHCwwsLC7F4AAKDy8uo1Q+PGjdPgwYNLbdO0aVPt2bNHJ0+eLPLe6dOnVb9+/WLn27Rpk06dOqXIyEjbtNzcXD399NOaM2eOjh49Wq7aAQBA5eDVMBQeHq7w8PAy28XGxiojI0Pbt29Xly5dJEnbtm1TRkaGunXrVuw8w4cPV69eveym9enTR8OHD9eoUaPKXzwAAKgU/OJustatW6tv374aM2aM5s+fL0l69NFHdffdd9vdSdaqVSslJCTovvvuU926dVW3bl275VStWlUNGjQo9e4zAABgLn5xAbUkLVu2TO3atVPv3r3Vu3dvtW/fXkuXLrVrs3//fmVkZHipQgAA4I8shmEY3i7Cl2VmZspqtSojI4OLqQEA8BPOfH77zZEhAAAATyAMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAUyMMAQAAU6vi7QJMKy9XOrZFunBSqllfiuomBQR6uyoAhbGfAqZAGPKGlDVS4mQp8+ffpoU1kvq+IrW513t1Ae748K8sAYL9FDANi2EYhreL8GWZmZmyWq3KyMhQWFhY+ReYskb68GFJ13a7Jf8/D77HH1p4hzs+/CtLgGA/BfyeM5/fXDNUkfJy8z8oivyB1W/TEqfktzOrvFzpyCZp78f5/zVzX1Skgg//wiFGkjLT8qenrKmYZfgC9lP2Q5gOp8kq0rEtRT8o7BhS5on8ds1uq7CyfEZlOargb8r88Lfkf/i36l/y6S53LMNXmH0/ZT+ECXFkqCJdOOnedpVJZTmq4I+c+fD35DJ8hZn3U/ZDmBRhqCLVrO/edpUFpyW8yx0f/pUpQJh1P2U/hIlxmqwiRXXLP9ycmabi/+BY8t+P6lbRlRXP1buCnJ2vIk5L+MMdTmXV6KltcPRD/cLJ/GtIilu3vwaI4vrU1/fTwjXXqCcZhpR9pvxj4sgm/zg9WNJ+4A/7OOz50O+MMFSRAgLzz7t/+LDy70op/If2v3ep9H3ZN3ZgV68bcGW+/Z87VpOrRxX84RqIsmr05DaU+eEvyRIgffn/Sl63rweI4pTWp766nxZXc2GujomUNdI/n3SsrTeP7pX0O2v7P9J/PvbtfRz2fOzvMrfWl8Htt9ZLJQyCxvl/YH1hx3X1tmJX5ktZI3043LG6Rnzm/L9I/eEW6bJq7PZHactfS37fHdtgq0HFrKc4xay7xGX4UF8XcGRcSL61n5ZYc2Eu9LVDyy3Elf3QHZyt0xfHHfJV0N9lZz6/CUNl8EgYknzq8GCRuua0LeVw+X//hf/U3qKncJydr8x5CglrXHSdntqWiuRIH1gCJCOvpDfdtw3FhXRn1+3rQV9yblxIvrGfOrOvODMmPLVcd3OqzsJ8YB+HvQr8u+zM5zenybwlINA3b8t19fodV+Yrc55CXDkt4Q+3SDvSByWGEcmt29Dm3vxb3ws+/C+ctD815si6r12GLwX9As6OC1/YT53ZV5wZE04tV947PehsnTY+sI/Dno/+XSYMwZ6rdwW5Mp+j89zyuGtHFfzhDid3rdtdyykc0vd+7Nq6fTXoF/CHcXEtV2pxZB5Hl1uttnTPm947ulfe34Uv/S7Nzkf3P8IQ7Ll6V5Ar8zk6T8u7HGtX2rrc0c4T3LVuT2yDP/SfK/xxu1ypxZF5HF3u79+Vru/hfA3uUt7fhS/9Ls3OR/c/njMEewV3BRVcyFaEJf8akGvvCnJlPlfX5ShPL98dyqxR+dfteGMb/KH/XOGP2+XIOLFxon5H+6LprY7X6glObX9hPvi7NDsf3f8IQ7BXcPu/pKKDtZTbil2Zz9V1OcrTy3eHMmu0SLHjSnlfntsGf+g/V/jjdpVac2FO1u8vfeHw9hfmQ/XjNz465ghDKKrNvfm3NoY1tJ8e1qj0Wx5dmc/VdTnK08t3h7Jq7P2i97bBH/rPFf64XSXVXJgr9ftLX5RYZ2Op25P/PdpQeLqP1Y/f+OCY49b6Mnjs1np/UFFPoC7Puhzlq48yKMxbT6B2R23+yh+3y1NPoPaXvuAJ1JWHh39nPGfIjUwdhgAA8FPOfH5zmgwAAJgaYQgAAJgaYQgAAJgaYQgAAJgaYQgAAJgaYQgAAJgaYQgAAJgaYQgAAJgaYQgAAJhaFW8X4OsKHtCdmZnp5UoAAICjCj63HfmiDcJQGbKysiRJTZo08XIlAADAWVlZWbJaraW24bvJypCXl6eff/5ZoaGhslgs3i7HJ2RmZqpJkyY6fvw439dWCvqpbPRR2egjx9BPZTNbHxmGoaysLDVq1EgBAaVfFcSRoTIEBATouuuu83YZPiksLMwUO1R50U9lo4/KRh85hn4qm5n6qKwjQgW4gBoAAJgaYQgAAJgaYQhOCw4OVnx8vIKDg71dik+jn8pGH5WNPnIM/VQ2+qhkXEANAABMjSNDAADA1AhDAADA1AhDAADA1AhDAADA1AhDJrNx40bdc889atSokSwWi1avXl1q+1WrVikuLk716tVTWFiYYmNj9eWXX9q1WbBggW677TbVrl1btWvXVq9evbR9+/Yiy5o7d66aNWumkJAQde7cWZs2bXLnprmNt/po+vTpslgsdq8GDRq4e/PcxhP9tGrVKsXExKhWrVqqUaOGOnTooKVLlxZZlpnHkiN95E9jyRN9VNiKFStksVg0cODAIu/5yziSvNdP/jSWyoMwZDIXL17UTTfdpLfeesuh9hs3blRcXJw+//xz7dy5Uz179tQ999yj3bt329ps2LBBQ4YM0TfffKOtW7cqMjJSvXv31okTJ2xtVq5cqaeeekpTp07V7t27ddttt6lfv35KTU11+zaWl7f6SJJuvPFGpaWl2V579+5167a5kyf6qU6dOpo6daq2bt2qPXv2aNSoURo1apTdH3GzjyVH+kjyn7HkiT4qcOzYMU2aNEm33XZbkff8aRxJ3usnyX/GUrkYMC1JxieffOL0fG3atDFmzJhR4vtXr141QkNDjXfffdc2rUuXLsbYsWPt2rVq1cqYMmWK0+uvSBXZR/Hx8cZNN93kQpXe56l+MgzD6Nixo/Hcc8/ZfmYsFXVtH/nrWHJnH129etXo3r27sXDhQmPEiBHGgAED7N7313FkGBXbT/46lpzFkSE4JS8vT1lZWapTp06JbbKzs3XlyhVbm8uXL2vnzp3q3bu3XbvevXtry5YtHq3XG1zpowIHDhxQo0aN1KxZMw0ePFiHDx/2dLleU1Y/GYah9evXa//+/br99tslMZauVVwfFTDLWCqpj1544QXVq1dPo0ePLjKP2caR5Fo/FTDDWOKLWuGU2bNn6+LFi3rwwQdLbDNlyhQ1btxYvXr1kiSdOXNGubm5ql+/vl27+vXrKz093aP1eoMrfSRJXbt21XvvvacbbrhBJ0+e1EsvvaRu3brp+++/V926dSui9ApVUj9lZGSocePGysnJUWBgoObOnau4uDhJjKUCpfWRZK6xVFwfbd68WYsWLVJycnKx85htHEmu9ZNknrFEGILDli9frunTp+vTTz9VREREsW1effVVLV++XBs2bFBISIjdexaLxe5nwzCKTPN35emjfv362f6/Xbt2io2NVfPmzfXuu+9q4sSJHq+9IpXWT6GhoUpOTtaFCxe0fv16TZw4Uddff73uuOMOWxuzj6Wy+sgsY6m4PsrKytJDDz2kBQsWKDw8vNT5zTCOpPL1k1nGEmEIDlm5cqVGjx6tjz76yO5oRmF/+ctfNGvWLH311Vdq3769bXp4eLgCAwOL/Ivr1KlTRf5l5s/K00fFqVGjhtq1a6cDBw54olyvKaufAgICFB0dLUnq0KGD9u3bp4SEBN1xxx2Mpf8qrY+KUxnHUkl9dOjQIR09elT33HOPbVpeXp4kqUqVKtq/f7+aNGliinEkla+fmjdvXmR5lXEsSdxNBgcsX75cI0eO1AcffKD+/fsX2+bPf/6zXnzxRSUmJiomJsbuvaCgIHXu3Fnr1q2zm75u3Tp169bNY3VXpPL2UXFycnK0b98+NWzY0N3leo0j/XQtwzCUk5MjibFUksJ9VJzKNpZK66NWrVpp7969Sk5Otr3uvfde9ezZU8nJyWrSpIkpxpFU/n4qTmUbSzZevHgbXpCVlWXs3r3b2L17tyHJeO2114zdu3cbx44dMwzDMKZMmWIMHz7c1v6DDz4wqlSpYvztb38z0tLSbK9ffvnF1uaVV14xgoKCjI8//tiuTVZWlq3NihUrjKpVqxqLFi0yUlJSjKeeesqoUaOGcfTo0YrbeAd5q4+efvppY8OGDcbhw4eNpKQk4+677zZCQ0N9so8MwzP9NGvWLGPt2rXGoUOHjH379hmzZ882qlSpYixYsMDWxuxjyZE+8qex5Ik+ulZxd0n50zgyDO/1kz+NpfIgDJnMN998Y0gq8hoxYoRhGPk7Q48ePWzte/ToUWp7wzCMqKioYtvEx8fbrftvf/ubERUVZQQFBRmdOnUyvv32W89vsAu81UeDBg0yGjZsaFStWtVo1KiRcf/99xvff/99xWy0CzzRT1OnTjWio6ONkJAQo3bt2kZsbKyxYsWKIus281hypI/8aSx5oo+uVdyHvGH4zzgyDO/1kz+NpfKwGIZhuHZMCQAAwP9xzRAAADA1whAAADA1whAAADA1whAAADA1whAAADA1whAAADA1whAAADA1whAAADA1whCASmP69Onq0KGD7eeRI0dq4MCBXqsHgH8gDAEo1fHjxzV69Gg1atRIQUFBioqK0vjx43X27NlyL/vo0aOyWCxKTk4u8t4dd9yhp556qlzLf+ONN7RkyRKH2laW4PTBBx8oMDBQY8eO9XYpgN8gDAEo0eHDhxUTE6Mff/xRy5cv18GDBzVv3jytX79esbGxOnfunLdLLJXValWtWrW8XUaFWrx4sZ555hmtWLFC2dnZ3i4H8AuEIQAleuKJJxQUFKS1a9eqR48eioyMVL9+/fTVV1/pxIkTmjp1qlfre/nll1W/fn2FhoZq9OjRunTpkt371x7t+fjjj9WuXTtVq1ZNdevWVa9evXTx4kVNnz5d7777rj799FNZLBZZLBZt2LBBkjR58mTdcMMNql69uq6//npNmzZNV65csS2z4NTc0qVL1bRpU1mtVg0ePFhZWVm2Nnl5eXrllVcUHR2t4OBgRUZGaubMmbb3T5w4oUGDBql27dqqW7euBgwYoKNHjzrdH0ePHtWWLVs0ZcoUtWrVSh9//LHTywDMiDAEoFjnzp3Tl19+qccff1zVqlWze69BgwYaNmyYVq5cKW991/OHH36o+Ph4zZw5Uzt27FDDhg01d+7cEtunpaVpyJAheuSRR7Rv3z5t2LBB999/vwzD0KRJk/Tggw+qb9++SktLU1pamrp16yZJCg0N1ZIlS5SSkqI33nhDCxYs0Ouvv2637EOHDmn16tX67LPP9Nlnn+nbb7/Vyy+/bHv/2Wef1SuvvKJp06YpJSVFH3zwgerXry9Jys7OVs+ePVWzZk1t3LhR//rXv1SzZk317dtXly9fdqpPFi9erP79+8tqteqhhx7SokWLnJofMK3yf/E9gMooKSnJkGR88sknxb7/2muvGZKMkydPuryOI0eOGJKMatWqGTVq1LB7BQQEGOPHjy9x3tjYWGPs2LF207p27WrcdNNNtp9HjBhhDBgwwDAMw9i5c6chyTh69GixyyvctjSvvvqq0blzZ9vP8fHxRvXq1Y3MzEzbtD/96U9G165dDcMwjMzMTCM4ONhYsGBBsctbtGiR0bJlSyMvL882LScnx6hWrZrx5ZdflllPgdzcXKNJkybG6tWrDcMwjNOnTxtVq1Y1Dhw44PAyALPiyBAAlxj/PSJksViKvJeamqqaNWvaXrNmzSp1WStXrlRycrLdKyYmptR59u3bp9jYWLtp1/5c2E033aTf/e53ateunX7/+99rwYIFOn/+fKnrkPJPrd16661q0KCBatasqWnTpik1NdWuTdOmTRUaGmr7uWHDhjp16pStzpycHP3ud78rdvk7d+7UwYMHFRoaauuvOnXq6NKlSzp06FCZ9RVYu3atLl68qH79+kmSwsPD1bt3by1evNjhZQBmVcXbBQDwTdHR0bJYLEpJSSn2LqsffvhBtWvXVnh4eJH3GjVqZHeHWJ06dUpdV5MmTRQdHW037dpTc+UVGBiodevWacuWLVq7dq3++te/aurUqdq2bZuaNWtW7DxJSUkaPHiwZsyYoT59+shqtWrFihWaPXu2XbuqVava/WyxWJSXl+fQduTl5alz585atmxZkffq1avn8PYtXrxY586dU/Xq1e2WvXv3br344osKDAx0eFmA2XBkCECx6tatq7i4OM2dO1e//vqr3Xvp6elatmyZBg0aVOyRoSpVqig6Otr2KisMuaJ169ZKSkqym3btz9eyWCzq3r27ZsyYod27dysoKEiffPKJJCkoKEi5ubl27Tdv3qyoqChNnTpVMTExatGihY4dO+ZUnS1atFC1atW0fv36Yt/v1KmTDhw4oIiICLs+i46OltVqdWgdZ8+e1aeffqoVK1YUOcJ24cIFffHFF07VDJgNYQhAid566y3l5OSoT58+2rhxo44fP67ExETFxcWpcePGdndEVbTx48dr8eLFWrx4sX788UfFx8fr+++/L7H9tm3bNGvWLO3YsUOpqalatWqVTp8+rdatW0vKP9W1Z88e7d+/X2fOnNGVK1cUHR2t1NRUrVixQocOHdKbb75pC0+OCgkJ0eTJk/XMM8/ovffe06FDh5SUlGS7uHnYsGEKDw/XgAEDtGnTJh05ckTffvutxo8fr59++smhdSxdulR169bV73//e7Vt29b2at++ve6++24upAbKQBgCUKIWLVpox44dat68uQYNGqTmzZvr0UcfVc+ePbV161aPHPFx1KBBg/T8889r8uTJ6ty5s44dO6bHHnusxPZhYWHauHGj7rrrLt1www167rnnNHv2bNs1NmPGjFHLli0VExOjevXqafPmzRowYIAmTJigcePGqUOHDtqyZYumTZvmdK3Tpk3T008/reeff16tW7fWoEGDbNcUVa9eXRs3blRkZKTuv/9+tW7dWo888oh+/fVXhYWFSZI2bNggi8VS4u32ixcv1n333aeAgKJ/0h944AF99tlnOnnypNN1A2ZhMQwv3RcLAHDIkiVLNHPmTKWkpBS5PglA+XFkCAB8XGJiombNmkUQAjyEI0MAAMDUODIEAABMjTAEAABMjTAEAABMjTAEAABMjTAEAABMjTAEAABMjTAEAABMjTAEAABMjTAEAABM7f8DHEWEcFsyFkkAAAAASUVORK5CYII=",
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"dists = []\n",
"for pos in coords:\n",
" dist = np.linalg.norm(pos[0] - pos[1])\n",
" dists.append(dist)\n",
" \n",
"plt.scatter(dists, np.array(charges)[:,10], label='H charge')\n",
"plt.scatter(dists, np.array(charges)[:,0],label='O charge')\n",
"plt.xlabel('O - H distance, A')\n",
"plt.ylabel('Partial charge, e')\n",
"plt.legend()"
]
},
{
"cell_type": "code",
"execution_count": 17,
"id": "a1030b06",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"20"
]
},
"execution_count": 17,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"len(coords)"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "e0f68bff",
"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.12.8"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment