{ "cells": [ { "cell_type": "markdown", "id": "heading", "metadata": {}, "source": [ "## Convert RKF Trajectory to DCD Format\n", "\n", "This example converts `ams.rkf` to DCD and PSF, mapping atoms back into molecules in each frame before writing.\n", "A short MD run is included first to produce the RKF input file.\n" ] }, { "cell_type": "markdown", "id": "imports-md", "metadata": {}, "source": [ "### Initial Imports" ] }, { "cell_type": "code", "execution_count": 1, "id": "imports-code", "metadata": { "tags": [] }, "outputs": [], "source": [ "from pathlib import Path\n", "\n", "from scm.flexmd import PSFTopology, pdb_from_plamsmol\n", "from scm.plams import (\n", " DCDTrajectoryFile,\n", " RKFTrajectoryFile,\n", " AMSJob,\n", " Settings,\n", " from_smiles,\n", " packmol,\n", " view,\n", ")" ] }, { "cell_type": "markdown", "id": "ac7af88b-b15d-4b2e-90b8-fc3154b55668", "metadata": {}, "source": [ "### Run MD Simulation" ] }, { "cell_type": "markdown", "id": "toy-md-md", "metadata": {}, "source": [ "Run a short toy MD simulation to generate an `ams.rkf` trajectory.\n" ] }, { "cell_type": "code", "execution_count": 2, "id": "5c10507d-d5d7-4207-8118-b83d7da89704", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAASwAAAEsCAYAAAB5fY51AABpuklEQVR4Ae29CZwdVZn+X1X39t6ddHfW7qxkhSSQsAqoJOyMIgpMooI6juMoAg6K+ndhxiQgm/MDZ1RUdFwYRRyigoqyConsu5AVErJ3Onsn6e50p/veqv/3OX1v02m6k16h+/Z7Pp9zq24tp6qeU+ep933Pe97jeZYMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMgYMR8A/+a//aQMCfP39+j+G0YMGCqI1rdHmT7/sqT/fX3XJ7oowuP4edaAh0BIEea4gduVh/OyaKIvjAEUJ/u/Xu3m+PvRc9Sfbdfaj2zrePSHvI9L3tPfZi9r1H694dpcnq73//e8Hdd999VFFRUWNHSwzDUNJO2Or4APLbT7nVDQ0NfnZ2dqckotbHqwxS4sCBA9WtrtPlv+Xl5dFnP/vZDj9nly9kJ3YFgc601U69W125mXfqnM6A8E7d49t+XaSCYOHCheHtt99etnTp0st27tw5aebMmc9BRCKddl+GNMmxjHFcbhs3HrIvmT6ujf1tbgqCwOeceMudqTKSbKvrbHkty0mvqwyuo/vbd6hnTB/fxlLnHoQN5ej9qmJzI9ip/IP2t1HGQZs4T8/XnPgf8GdvLBara97YwRX7iHQQqD5+mBFWqwpKk9Utt9wyKpFIzNu4ceP43Nzcv0+ZMuVXy5YtC0pLSw/b6AYNGhSrr68f1JZU1OpyHfqrhhqPx0uRqgIaa/P1WReRqRF3JZVwbja5ZXkeDTvWmcLSRAQ3FVFUPue2JVl2+j1LJpNxnk9lNt9fZ+4rfSynizdVTEZ/RBobGxs2b95cNmHChIe+9rWvLU0/dxqHTFke9NXOlIfq6nOkyermm28uh3D+iUpfyguwEuIa1gVVqdNSwGHuu/Iw+zNuN/WR3VMPdbiPSJMw2PGrQahtfkRS0rAk8Y4X9uaRb/mIqByRLfnNo1Jr+lgo5+fnH9i1a9epGzZsOI/1enYv7eL133KNvrbBCCtVIy3JCoL6BF/3tXyp/nTdddedw3oaJ72Fb31z3qzV9H590t/c2gNrGIa71AJ64NIdKkKG6x5uJBFqeUOHLt7xgzLyI3LDDTecguRfe8wxx9xfXFy8q+Nw9L8j0w2x/915D95xa7Ki6LVIVr/TJURWNMSOsk/6uJ5uvLqVdNla73MJcumNe0p/ALpTdrqMjPqITJ8+3Z83b17yxhtvPAsp6zQA+s3w4cNn8LFty3baHfz61LkDnrDaIitsRb+dO3dumiBa22T6VAVm+M2k66A7j5kuI2M+InpnRVbXX3/96ZDVe/ioLkIbWIE2cCLrGf2+dtVg250XqM+c2x5ZqYdw0aJFAxqbPlNJdiMHIZB+Z0VWSP6nQVC/SxvZM52sBMSAbZTpipeBXTYrsFgryUpkxbq/fPny9JdZOFkyBN5xBNLvbJqssJM6stL2Tpgt3vHn6M4NDEjCSld8e2QFoEZW3Xmr7NweRyD9zrYkq2uuuWZZavuAeV8HHGGlK97IqsfblBXYOwhoLKtzZG6HrDLaZtUa0gFFWOmKN7Jq/RrY/z6KgHNlkZnCyKqphgZML2GarG666abR9KxcyuMfZLPi/4ARq/to47TbOhgB546h9xZn0DNwVH2PbFYt1MABJVmloRkQElaarPhKzaDSP8nDG1ml3wBb9kUEHFkx6D7G8K6LjKzerKKMl7BakhU9KRcjXS0nAsM9uC243kCgMMnqzffB1t55BJrJavXq1Rfyzk43yerNSslowkqTFWrg0VS6yGoF463+MG3atAjH0JiWQOHsBG9C0rwWcL4GAyvKgDN8rlixwk+d03xQyxWLq9QSDVvvAgJGVocBLWMJC4Li4+SHt91221HV1dXv5/8L6P/3HQaPlrudjeDWW289wMaGb3zjG4e1GfTS8JSekABbltEeQbd89g6tQ+g9VlaHLtiFg/rLRyT1vroPqUlW7Vd0n3/h2r/19vekK/9//ud/xj3zzDNfzcnJqRk3btwjhOBoYjGkJs4+ZFwl7AZ+Xl5ect26de9FyhoyceLEexXehXU1/tbkJcc9C87XfpXYnkMgkH5fJfUfd9xxF/F/GrkzBna9zxGDoD+uy/Bx/SULt03/MyllpIQFebg6Wrly5Vk7duwYMmPGjJUQzZFsVKynCBXP7eelcMu2fiSekRp27959TF1dXdH48eP3cLzObzeuEvuSWVlZnXpRpG5CpK3rQWUkIdoeC87Hc4S80H0iOJ9IH2wPIn2wVaUc8iPC/jZTf/+IgAWPH/m421zM8gge8k60gTVpk0abDz1AN7ZuKBkFQ0lJSeOpp5563xFHHPF/jz/+uN+R4HstAAiGDh2aPPLII89l28jBgwffsWXLFr+/x1WS5Jh+RsiyfcbmIIjARQmlQRXx96DgfGphbCtln88YtnSRHVpySg65gIObSYv/7tymYjtUjOJEcZqLpNGvPyIiXOKtlY8dO3YPgRrv/PrXv77ZyKrtdyCjCUsNkgbqMbK9y3GVMNgrxnljqwB+GRlXqe1XpHe20iCzFUNeH4HuXqE/f0R4P+GoeEjAyBORgpdhM90sdwZFY+guLpl4fkYTVqsKU8M4pETR8vj0F05ShlJqn1t2RgpoWWZ76xiG0+W3d0h3tmtYR+JQBbzwwgtZf/rTn9ptIDJcvwnBoUrq8D4LztcKqu985zulTCiiDh4beN8Km5Z/BxJhtXzu7qz3dOPVvXSYSDtz4zLi4m+WOOmkkwa9733v+yh2uAuQaoag3nrYx/bzdX+ioqLif0444YSNLdSrt1yil3o/O/UBectNNW1Il9Fvg/OlXWWwkyrwXvPztPO8A36zEVaGvgIpskp+6EMfmk3o3J/MmjVrMrY8b+TIkR62PY/OAT356fR+/hudE/ORoP4b0pItSHalXiHQVlD3xDXSZfSbj0grDPTXkRSRQ4E//ThtHGWbHAIDYmjOQKtrqbNIVslTTjnljNGjR9/PjD+T6UBoJN53ErJST6bUP+UE0tZgOhb+C9K6BbLSFGSZ9k44p19hwrOJHCz1YwQy7eXsx1XRY7cum1WED1k5xuh7R40alcfg2WRhYWEWZBWDrJxrBlfTUhK2JKrGSZMmXY0t68IUaXWu248C+loSQcl4zX3JXhYq82ySYtLb+9ot2/10AAFTCTsAUn86ZPbs2cGSJUsSZWVlX0KqKmLwbCO9UFnYrGS3autR9NGKOCaC4K77zGc+cx//D2mkb6uQvrRN6jAE5ToRIOKcz33uc0Mh6mymwtoBadXoXlOSpFRJ08P6UuUd5l6MsA4DUD/b7YusuOc8XDo+jH+UGmPznHY0zjYfh8Yr21WEdjj9/PPPP5r1l/pr13radnfJJZccje3uSkY4nD1ixIjhkHG8qKhoF9tf3Ldv3/d4xocFBs/OqpFWmy9GH9xohNUHK6UbtyRpKYmtagoNcSRDiTz8ewKyRy+Up/8pY/tBl6DB6r886zXb8nGsvzRs2LC22e2gM/vWnzRZXXzxxZcxSP07dDTkohJ7TH/lgYmHtFnOHSuG/wfWr19/14MPPvhp/vfIaIK+hUTm3o3ZsDKwbiGlAhplrKamxlNm8Le3Z88eD8dEea+3+cQc79XW1nptDBNq8/i+tjFNVnQ0fAb3jR/SI5qDVJXAbhcOGTIkgqwkbSrLUzOB5PXRs88++945c+bI9UPtoEsETYH+3Uixj2EP1JL/gdbT//saTv39fkzC6u81ePD9015oLfH4bqSpRmw28aqqKkdUNF6voKDASVjYtnRM85mQlLdz505/27Zt3vbt2/c37+g/K65XlKFXR0FQtyEdhnQ66HnjkqxaDR1yKjJSZQOkdvZ3v/vdK2bOnPnfjz32WPz000/vjO1ORBXAcuke15ZoNX8VqBB/EcfNazqu5TG23gUE3nxru3CyndLnEFBD8XFjWLNq1ao3IKGpW7duDekhDHJzc13D1VAlSV0iMJGWJCv9h9iCNWvWJNeuXfuMnmrx4sXNja7PPWWrG0K68pGSpPb9O4QVR8JMQFJxpKpDqcB690Nw+OrPf/7zO+YUFdU+Nn9+fA7GekjIEX+ryzT/nYskBQklRUIfLSsbOjMI5gyNolMHed747CgqLfb9VXRJbqyJosf8LVue5sSkpC8jrWYIu7xihNVl6PrsiTEZ3vFo//3+/fu/QcDCBNIGNvjAqYOyZ0k9hMRcY5aKiH0rgToYf/311x9n7ORquQTIDaDPPuHBNyaykpRTiNR0tgzokHJMzyVyPkRy9j7UxbKCQYNO9E84wRnhdXyEO8Ti5cv9xRBaaxzSZHXisGEj35ud/eUxnveJI3x/WBl2wGLOLSDne95sFb6P/NqoUQ9vj6Jvv3fLlkdgQUlkIsNDEiL7LbWDgBFWO8D0481qpQHS0/9Devow4XUmst5IA86S0V32LKlJ+GY5wsI4n2QImw+xeQSO+zrnRgwXUXvrL0n3qo6GiTyjxuOpo8FXJwOE7TobJGm1k3yIOhw7YsR3lz7yyC83r1nz2E133fWyP29effp49SIuxs41Z8mScA64Ilkl5owcOXt6LHbHVN8fVwYxIl0lkaqiwUHgF/FhwHnExSEqDEPFsTi7LIrOXlpWdqtfWfmlFGmpeCOtNMidWBphdQKsfnKoGkKwcePGKrzcP0ojfqiysrKYBhzRne/Txd9MWJJGsF/FIDVv06ZNlz/yyCPPpI3X/eRZm28TAs6HoGM8YyRSZukkSdTiZnJuPji1ArkFevbpM2YcWRgE149CRf7fT31q/Z4zz3x2T0XFAxu2bFkCRus4PG3bCo8bNeqMMVF03wRcRwYnk40QVLw4CGIlUFMJ5zvbIKSFv4SHmIdLbmOyCBvhjHj86pVlZbmQ1hVUUCxl+2p9S/b/MAgYYR0GoB7e7c/HxjSdPIyssnfwaqNmKASqiKankpOyNm/e/Dyq4cmQ1g9omKfv2LkzKsjP92V8x4UhIvuoTasxzn8R6erP/ZSsnOoKUaxHgqzlWQrIERKXj9+VJ+O7jO4iatnxICCnGksSo4PBSWE5WVkhPh9hyciRMUhn/Ni8vPGc+OGtnle/5oILXt67bdtjycrKP3+usbF+cBTdO8r38/KRqvKDIAui8wohqkE45WI08wBWvR5N9Yh90DtwIIaYF7FsPDIev/y5srIqSOvfjbS69qobYXUNt86eJaKKLeRLTW6TmPQCU6gGu/WU7Sg677zzch5++OHXsF99ZcIRR7w4ZfLkKJsGhf0qSeMNampr/3DnnXd+nOtW91OyUj0IT0mU2/DuX4nrxnFIlBHPGJNnv2x36gWV1KWOBqmH6mjYjwvH9h07on1793qbVq0KBicSQQEklr1/fxgkkxGE5Y0cPjwXdjuFE045UFPzjX/atSvxkO/HC9AToaUYUQi9HMgwnzJzIENEOeSuPI+L8DmiJkVYqKXcBH8YBlVfnxzj+19/cPTo33ibNy9PGeJ7qr6FRcYnI6xeruK0kRaikloRu37UqOnFnjdpqOcNicJwf24UvfYX31/BV5c327W+7qoLGuyroSmJBx54QPGVZHw+chCN6bTjj08eN2tWfOyYMRqK4z/z4osbfv3rX1f/15VX5lz1ve+5Y3V8P0yyY6ln8DtITnfiFNqIZBUTWcnwLnvWXogp7daBxCk7V4Qk5m9Yt87bsHRpsiw3Nxri+0EJEVlLGxuDIjon8vbvj7KiKISAIpy6YmFNTTy7oSHKQkyLQ0gMD/DwtG3qiZRkJcKiM8NJWVILIUrEOyoVTk0kfHI0ki2bk8mvwWAfa/PL1Q/Bfztv2QirF9GezVdVRtpzhw4tm5mdfeU4z/vgsCiaXs6LXkLO5mUmPII3MYo2XFFe/sDmILjN37x5KS9yV3uT1HDVsyVyzLn22munTZ8+fdr3v//9f9i3dWvUiOtCbUWFt1cNievv27TpNRmVjyktlQrZn5OeNyBUzt30dP4z9quzWDYgSWXJAC/pimFHzbYsnjkBscVx/Wj8y5//XFNSWFgybtgwbzQq5CgkpTJUuuGQ2hA6I0rq62MlSEm5kN4LSF75kJWIRhnJyWViRHuwVpNklSYu7ZOEJbKiswPWRMGsj4Vsy4ui9/9k9OhS6nr3fErgY+aKU5GWDo2AEdah8eny3pRklSBi3j9P87zrpwdB2Uhe3mJyCULPIF7SvKYvdcBXfBxf688OTyb/6e9lZQuRtm7iDRb5yNDV0ZfZkRWn5Fx99dWfxxfrsxjdJxGcz3vttde83/zwh15OY6OfhSoUxyCdQPKIamrU0Dtavm6nTyfcOZKohR+FkP6EX9nJ9HbKaTZJDDDZtBxhIXXFIbI4aiPcXfGh9Rs2vLIeN4S/+/5ZgwsLZw/Jy5s0jIgWo5hMZALkNRYyGg/h7MZXbQ3VMQ2ih+6bYvOwTPBf2amAKQJzRCYSk9FdyzezD+BhUSxWXBpFp3D6nz0Gq3v0QLJuqQMIGGF1AKTOHpJWA+eUlV03yff//QgKGByGiUGoHKWxmHKMl9ZDZ2l60TV1WCKRRF3IJd+4vLz8aBwOLxVpiU06QFqaZixkEoPxZ5111q+ZJegUpA0Pj28Zm5P0APrxZDKIIykEqEZqXAka4f4dO1xDWdzZB+ybx+tZAohoJ6R1Jnhci0Psv6AeFsuDX+qg7FeoiElsWo9y3FfoiHiF4/Q0v4XkfrunujpGPi5/5swH1vl+yStbtnghIwDGQTwhQQ/rUPkkESMvedKflev5X8d+6tcLJFEpS4IVWWldvmBa5zhl1kJqPaA+xnK6V15T425A65YOj4AR1uEx6tQRabKaMXLkZeNisX8fAVEV8oIWYfYQSQ1GdRhEDmSYVW+SGozmSdRUXwcOqDcpMS0Wu+TVpt6kK3nN+US7D3pzSIEFTGA6nVmodWOPlJQENFAfNTA2derUu5GsTmSISgMG5jiEJalLNjEvm4YjwopBWHhWxuplDN6x4+/s8ihLvJgJCT7wEFCdPfDLDHz+L4joXaiEx2PDKsB29Qr2rGfxvVquh011NOicgMHSMSSyhhNPPPESQvSUHnvssclC7FZ7kUafe+op7/f33ushFXkyNKLcucxElF4tZFQNtkV8AIow2rs6FTnpYySykioolRDy4l400tqdi1QmQc1SJxEwwuokYIc5XI6FYenIkUeVBMF3SyErTLExqX65vMB5EFSBeq7Uo6TeJNk79GLrC6yXvbbWJ9hJll9f3zg1Frvi8ZEjH/W3bv39C1hITuC73fwpJkBfi/tw9ickqqtxYTgRd6QGesey5cmurBTScLJpMDHUQUlYPl34Af5XI2fOnOzdf/8Lw6ZNay66Rbn9dVXY6HkCVL7NLJV/R25OEAc14qc95LVdzrINqI3H0Hv6BSQxecsH5bhCXHD++d7HL7nEe+TRR72aXbu8Wj42DLlxuZp6KwDbfCSqHOpP7hO4OzQRVMt6TZFWLXWwj4tVk3dHUZYubKlzCBhhdQ6vQx6NdKWBriGq3zeH4KOTA2EhR7keJfUmZfOyE0nvzd4kkZakLH2J6ZUSeRETU/8DeqNCyriOC/7xBMiKZe7njj46L7u0NCoqLR1HyJQcwojmZ5WWTq1BzaS37/+DqELZaNQ7JgfGLDItz3t95UqvAEkuwPgccB+4d3uNbB/2vvd9jca7CCnPkd4hH65/7RRp6ZkkYYq8WhKyIo9KqtIxLqWCHoYQ1UfktwV+OjfeCMHI8VT+XGNGj/Ze3LrV24tkPJg6KuIAkVMupJXFcQFSlAosZT0f3KkEz2cfnrlegrrdz/69YF5FxNNt1MWuWGyjLr6lsLD5PvTf0qERMMI6ND6d2SvpKplbVjYWkvpQFnYpWotTx6T2MSOp+wIHkJaTrERW2FUcYUnCgkgccfGljuRsyIsOeU375ZgxLxeUltYSXW9kvKRkEATl+YzwzSXGU2zIEK8UW9VybC1qXJAPp9FAWG+kHBqmpx4sdb/n8j+O4Rhjl5dkyAgxlKPyyZOn33///VPe53krONfZwTrzwP3gWBHTYRPGenccRHWc3CBwQPXlWKoslwhFtziRzovHn3jC2wN5FVFf+ZSaA95x1mW7EuskWK8D/0IwzwVngsh76hWku9Jt38c7sQ+cN4bhntWx2NPuxszgftj6aXmAEVZLNLq3rq95WBRFs5Fucvm6JkVDfM6dJ6iWSUgr0pdX5CTJSlm2LBEWDcV1jUNobj/HQW3RrPz8Gdk0kpBeLhl+mb4anbPUI55xRBiGKA+pSsZcGlfMjaHDNqUxdA1cby8EtU89ZZWVXj7lRzK4q3Gxvh+1lNjJseTJJ2fR2ILFCxa4++8eBP32bCflQNrDsW+56BWKYCF3CA0UVxDA8z/wAe8H9LRCON52HhOlXvXjBao7kj4SB1ivpd5wffCyycSxcXWrHtl69tXwTuxDAq4Iwz//FJeGu3k/5nGEK8B+OoSAEVaHYOr4QXylj9fbD/1E0uOUG8gHRCrkBJkQJE0klSaq9H+Oc4mXXQly8+tRXw7g/8OX2k/y5Q/54kfYukK63THc+jU5OcEuAvMdQO3YS+OqpsFV0dCe+uUvvbqHHvIGrVnjMaDQSXUlNMJCJABdfyjDV2JXXaXereO42isRk6UyW4Ou7Bqvu4GB8yPA9dw75GCqQIeKI4ZjqbMDrsO5lMCAHgH/vL/cf78nKVc2QTG8TpJk1QCmjpRAEJulx2fIi7NNUq/qnN7EqBbVdCOOqOuC4CZ2y3xgqZMIGGF1ErAOHJ4vksIi1ZR5WetYV++QepRkeMU+1WyzksQjCckZ3VO9SZK42OrtYt/6ujqsYVlRI40A47lCK0QJJKjk3r0RcWM8n0ZVQU9WRE8jQ01iJZs3e3sfftgbvmqVNxlpbbQ83Dl3EFlSARY019Ag1iDx6qthXXb2z1aUlb1nQXn55xcykIQ74cgBR1pOugSSR5GsziGGWLPfFrG1qI7QqddXQfDPPP20V83HgRCmwtwZyvRBkntDLVm2QuEswkKOdoSV+mg1VgVB9uuJxA2/2bp12Vx2Q1iylYksLXUQASOsDgLV0cN4+wK9oNXk/WS+qu5FRh3w9kFA+ZBSDi88A2ebpCyphSIskZVcDVgmOU4NYBuqxJOQURbDRQ4gWTXs3es30FAaaEQNqJUHKKOBrC+7jOnyNxqBPWssUtYExs0N59qDydkcU8B11EvpeiXZJhWUDoCAmDPJQVH0qUsTiXHTR4/+wKLNmxu4m54ejK0r9uXk1DIE1ntQ7W5k0LhPD2sEWTkySY9FxBHXm79ggcfQJ68a+1YyOztEavapNV/1vI9cSL21JCw+PCHvQ7gPslpL5IwHt279jxRZmSrYhTfCCKsLoB3qFGhou6QjvbwiLWUGy3p5kEoORCTpSoZvYii57nARjSMsVDVJWepR2se6pLEk++6vrk7U19VVO49qz9vI8gAqBrJWtEJOkBjQ6/ldPnjo0BmTo+jzYxEHhkJE+ZwvUlQkgRJIDltXk71MpCUZCsJyUh2hUeihbJgcj59Zl0jcNtPzPgVhcVPu689iQCSRR2zLli2v4xryPQzt/0b0VQlOaH58PFQnfDgUikYzZ1922WUeUUqT++vqYjvAuFrxsECVHOTrJDIYOuMUYw6DGvKuKPrB8/n5V89n+8KmXTrEUicRMMLqJGCHOzxKJJbKkL6XA6vIheRcyMn1KEES+D+5F7mB9UEQFLFA3WymUg31JVcvUw376AJPFGRlxc8aNOh7P9qx4zqKUV3tIB+UnAjAFibeu3I0RISDaqRrKZJALoSXjzSmsCrO70uDc7k3J2WJsOQfJMJUGyPq6MRY7J//Nnr03Yxxe+Bu9swbeKQV4Iv1FexY0zG4n6khTTK6b6HTQj2FhObxXnn1Va+Uzg8iYcQeevDBdUS8GFIfjw/aBojbwJ1OFyfVllEPeKp6Sz1v5U6sg5AhkDqyUpUZWQmMLiQEAks9gQBqghPxB5WUXJIFCWCz8sUuu1N5Fy+v8g6IaRtSViVf7QpUwM2ob8oVkMcWCGwrhLUdMtkOrazj5R+Rm/tr3vAqKmpH6k2Xx7tWRWBZp7HMJ1wvEQRuUgNhY6Ad8gOSI6N8v4iz8qbvF4OAPXocGQ3sEW+lKbqACA31B/tLhBH+Gk73UFsGWqPS88qBVJEe3o96+N+Q1YE6IpKecMwx3uX/8i/eFVdc4V14/vnhpvXr9z/zzDM/Om327Jn0zs7AfeGy2v37n59aXR19CZn3LvD9PUOifllY6N/Q2HiTyAq2kuBlZKWXqxtJ77al7iMQMDRGhDWKITFnDaEX6fXVq4M9EMHWJhJpelNZV4+Regv3IwGp+1v2DlWCHEZR87wGjsF4m0CNiL8ehotvrah4qYXNQ/Gy0kQizTO+hE6q0qysc5HUCvAJ0gQKeIb6BNZiDdJqdqEQaaUlLamFkrCUZDuDPMkxBtlFRZ538r1jx07zN27UeB1FjdBzDZQkbP2UOiiP98VDp0+/59gJE6KjmKcxb+LExJQxY+I5VVX/995f/vJzL774ov+zn/2smh7F2zmvIb+o6GfDiDU2uLQ0SH8MkHCHo0vGhs2eHTLIOV13AwXPHn9Ok7B6BlKHI2P4TsXmkUvI3STr/n4kJtk4tnANqQxblSElkVglhFGJtLWF5RaWFVqSN9P9zfHBUkhrXSz2eU4JpzWRVFsvu2tgWZ53JsRCybhYcYKYTEu3niIuZ7fiXpwKKDUwnUVeqXWkiuRg/IRyGxpmc7q3mPvQcoAlYRpDYsZ5vWF9A64gtXRk7N+40fPeeMMLcROJV1c3aP+Pf/xjfWuyAEvL6iqk2aqCAr9WEiwDz33sXbhATGKfqsJSDyAwEF/IHoCt7SJQwYZisMWMlYjks5ODvagWFW8bhLCJU0RcLkMtW8gV6QypbCFvwueK48KV9K/jen7FPZs3u+7vhW+Vcvy7586NzZ89G55CmIqiqUhnGk3rIy81RxNAbnLuFPITchIVhCgPeudKAUm6bdonqiOlSQ5j/9Fuw8D9ia699lpJljkhHx2i/3mM0XE5YOkTV0wS9We2bBFw4d+avhHra/izPz8/YOQ5RsWhvldW5mUNHTpWMM6ZM0flWeomAvoyWOohBJBQGuVtTrQAj1H/3qWXXuoxxs85cjYyDKcOUpKndDHXk18UYw01ENkN4YBwIgazBfv4Su8Iw88/WVn5Y325FzU1Bt2hSMp9YOYxrZWyNpLyIaxcmpVzo3AkpXWus5/yZcCvgaSKpfrJyC5yklQlopIqqO0p8tK5ytwLIpglEIhEWCGG90iEpQ4LMCOeWBp7gZQmoi2o8jWMICg8wCgECIsaw4G0pqacYwIPx1xvIZ8eS91CwAirW/A1n+xEFHrjXtL4M7q/Y8Rh8oj26V155ZXegw884C1dtsz5Ye2CkDTWLBsSGYxNqQzfHgYtR4xajjbu2LHmmdWr/4VSH797/vzseQsXNkqSIg6KN2/evJYkVfC1iRPfOzSR+NCYmprzfhKGY5fGYtF+VJkaTqbhNEcTUCQBjSNUR0CBblcSllRAEZbISr5f7Mc3wlP0gWoO2UNT06GWIG+w0ZCmCJcGR/Ry2t2+vdBhM326q3enbkfRPjpa9tVmZxfWM5QqgR2LGV29rD17Rp130kmFfMz2cTAU1myDNHi7gIARVhdAa+MUfXF94k8tpUdoOSP8jyK0SQiBxRRI718+/WlNo+U9jZf0RpZ7aACTZszwzjzjDG/a0Ud7w7F3FBYV+aiQBfm5uR/5+/PP7zvzvPNe4SX3nCTFrMakgmvOOee9wxlYXZJInDeyoWHcML72ozGmb6SX8WlIah8HiXC0dJEEICE3pk3ERCrlP64SzZEERF5qkPXs34uUtYcydnFNFKDXdfxi/QzUJElUEhbYSMJyKqEmkwDrmu3bpzpY5s6VdIXDOyNFGSXFEKotNfF4eWNhYZgkTlkcqbqgtjb39AkTCh947jmGI6Qoy51sP11BwAirK6i1fU5Ar1EjA2W/iQ3rd5oIgUHQjI0NXUxxeUl/8pOf1GSl3hbiqp9y6qnehAkT3PRTw7B3yFeKiA5SHy4//uSTP7Vs2bKrzjzzzDs/M3fuacOzsz9YnEyeN5KJO0dCOvkYglFLwlh1tcLJBO8iHE0RX/6dkM0gCsgnyxcrixzj+IhGh4Opk6LSkQTcODf2NUBaad8vCCvYxLkVvv84C2+FSQPxEHxCPjChJCwlCCuLYVFNf5p+F+jLAlYMft5B725TZNKSEqaCLo7yDhwoGjlhgup1y6J586TSt1Qnmwqw3w4jYITVYagOe6BexBiS1e8hre+jGl4JOYVETkgSAcDHGTF4+eWXtfTe/e538y4Xu2mnFMZYMZhSSQ0hyYQJ+HzGbr/tllv+Hb+oMWXYTvJoKPROhVkQVUT8K7zbg7r9+wPsXV5WY2M4GRPYct+P0az8XAqh57A5koBzl4C89ovs6CnEwcF1a8noLoM8koE865P78MheH4ZLn54y5cWoAtoauI0rTUhVSKEHUAmzCZUMJTGYmY9FTBJXC7V5YVNvaljf2LihhmPq8/OjhAzvw4eH6IKxwVu2jL/77rtfpq59r0la1vmWuoCAEVYXQDvEKVIRFOny85AWvqH11zIFfFyj/wlh5QbRyhivgHCaKl7z5LUgKxWrL3WcAzWvnj/7jDPG1GzaFObivBjt3h0Qhz2gQQS1SEzVNTXhPghsV01NsJcB0tizAp8yiWSpUJZ+TKWRkPAIVeo3Dc5l6UKfsN1VPI1Lkpf8wmo5FAkttp6Y8sSHSixoOiTBoQMxMbeXFyv18DxpaKhBJRySRRRSdU74EFYWkS4AJcIipdA8zN/laXxhcN11171RhYpdh5quUEBeeXlEmGUvWLGiSDZIVQ1Z6mOaEAcitt16ZiOsbsH3lpObXuQm0rqBcWl/QRL6CI6FszHEH0nM9WJNhoAWh90bw7t6ndpKKRUjr6goqmVcYPWePUEtIWPwUAz3MnHE7m3bAmYtDip27gzWQoabams37QzDB3Ly8s5N+P5YfL1ENHGJfI1kDaSWIZ6vvSdZTld1hCbC0n4OqyEmxMpkcsmvpkxZNL+yMkBqGKhkpYk/5KAraJLTwVEdKevBWbGyslhvrKlxAizHCF5v9rhxcdwc1MFaGeOj0UDeCLFVgyNRX6OxxxzzGeyXm6nWh3W8kZZQ6Foywuoaboc7S5KWBtP+naWyBs1+G4L6ChIPJq4Ek6Yk0cjUJtpIUsaQkhgY7W/cti3as3FjctemTbEdmzYFRBII1mED27Rz54bttbUPbKmtvZdjnyDXlEXRezjzoSRTqVdq7j3uAUcGRRN4cxA26yIsGVPUMLmDBEHnslAF1y2PxS6NmCprQZOkxxEDK80Fr0WQFBj6Xxg9+tLyMLxyhOcVDUciws7n7yQOPmM0o6F5eUevzMt7Fan0F3/Kz//xwjVr9p13+ukTj5w16/PHzpoV0fMb5EBaSM+BJOiRZWUnM070Ieya9z744IOXQlx4nZik1ZW3ywirK6h17ByxkXhBX2NG3IRj9KXWyH+F3VVUUEW1HEKYY0lbbSVJYuu3bPGfXLIkvoXexc0VFcu37dr11y1VVX/heBnGxUVOjzwKHlpRWfkExPgPlHcHYtw4pC6peskqiGkwhxVxrCQsbkgShAIAxhgilLUjip5ek0xeum3Llgq2655FuAMqzU2R1aljxkw8Jgx/MdHz3jOOD0cpUugQkBhMHRVRH5C9z9jBLBrOdED6z2G1tZ8dPm7cDX8bM+brR0+bNnlUeXlIHQRjxoxpGnTehCLxFrOicePGfYj0R2Y5uoDNaOKRuhcHFM7dfVgjrO4ieOjz1fClWulzmiuCktFdoUoUBUATHMiOpUgArZO85UVqGMWrF/35zz/YXVX1R455juzEMr3o/4iDJxKBJKVwRZODe4zgc0sw7p6MNPcf7LqkPgiKt7GCbcsZ22VwV4x33RhRMjfgXPr99fn53/PWrEEgG9hkNY1AhhPC8LeTfX/EMCRhCD4ooYODGZC8EpxtFVfMRboAS8TjkJ6PaGwiMWnFsGE/kz9dXm5uWFBYqIgPLckKWOUbLJOX38DkIWcyFdutrF/GNn2pVBWWOoiAEVYHgerGYe4TirS0DqLC/3CHk6pkeJe6oIiWUg31kmtdiTFsIrYQ21eQl5//MmT1dUlbeulPO+20OEZxVqNQ6kur+9L/GNfYyvIKDP83ctL7wyh6NwaWKTg2ityStJylHPjoztzcB721a/fqBtk+ICUrHj2YxuOXlJePGeV5vxnreSMKISsiV8QlURVBVMXUS7YGj6t+JA2DI0EWgxj1dB/18BwzGX0QWyORYZluMrs1WXEJxLImSUoVnOSD8q8PP/zwbWxbSu9hTAZ5d5D9HBYBI6zDQtTtA+ACveexB1EJv0gPop/+AmsqLvXiSZJSL6IITC+21EbiMUVIYYoc8ATTUMWPycuLfe+BBxrUg3eYO9LLLw5Kz8t3O+vKbaaUKqSv/ID80vP8/kJIBNL677IgGJUfho25UJNmvVE8MTna0pnRFE9MSw1rIvkaMYBqf199fZCX+tDog6L6PERy3wYGxmObz/8gxy3FF08fCiOsQ4DWcpcRVks0emfdEQiq32NITEuZ2OBoJjVIQFZEg8GlALuWVEWph855lIaCOshcE3X+qlWr6sg/EkktaSIUR34duE0dp+uqMSjrf7olOTJLb2tDSmPXgEnOyF48evRsQhtfSFTYJA0iK04dKJ6YhjNpADsV0zQlm6Zlk6TFfo0BDSGvtdu3O/Kqb2GflI1SEnR7SRI1dTxd+48/3s1Z0t6htr0VAkZYrQDppb+44qxowNv907ysjzFMh5EzkQZKx+VUCpG5F5wxhXJ5YHKcBs3YEoOsvvrb3/52c2pK9a58hduSnFqSVy89br8pVuStyVA/RjyxiGm51BnhkiRdxRILJFGlSYueP5EXBknPZ5sccXfjm6UZtfGJc9OCyT6p2XZw/tUY0VRpBy2c8zDHjUttTX9IDjrI/rSNgBFW27j09FaRTQyXhOewK30A0vojpFWglxt7U4jRPeSL7GPDUpuJydEUl4jrn3zyye91g6x6+hkysbwEPXe5iQMHzkEMlTdn4BgeQsI1pCkjaTlDu2xXyhCYJsWV+lcAeQ3l/+vUVxVZ9aZOFPnaiazoDWyWmgUeUpVsk/gA75a6/5q2EVPLDO8CooPJCKuDQPXAYY60sGE9iqR1Murg13hxPwRpFci3SsZavup0CoZSHb8NYT0i72kcErsiWfXA7WZ8Ec5BdHoyOZohScNlGMQq5WY8amD9gDKkpIy6iEINlclvTpn/TLnmnG9nIHU9BVFtI6SQvNql1re2TYq8pP7LVqleYj5SPu/BcxmPcC88oBFWL4B6iCJFPnL+XMbyYxhcx2LDms7XdjRf312ogiuwZ63S+ZKsjKyExKGTtLhFixbJTifMpO7KhUTL9pLPccHIl1+O06MhXsorjcVyGbgc0ZPqk13QQ/WounhiSEVFGNizsVGlDe7OhsV2xcb6BLaqn2J/3Mh8kDlIViIrSV9S6zVztGxZIixto67lqhKg6u999tln79UN8mGyD5KA6GAywuogUD14mLQO18DWrl27kXXl5qQGSApohPYiN6Py1hVJn8QbE1bC6SCswNCpWexrJi4dPwfcz7j22kQK2+SUwsKhH0kkznkJkluFmicv3NpU1jySig+2B2JSPLGh7Edfl88JpQQePlheEuP6MRzzWYjqFiSnbMIna3JVkZWG8cg2yZhQNwSLeyHwQ2MC+1b20qVL599zzz2V6v21jxKAdyIZYXUCrB48VKSlJOKSnVc5LR1o30ENkP+WWiDQUvokquugk046aTKe5UTSSdb/9Kc/fQ1ykOTkycdJS0leIoyFTT2l/n+ceurJozzvn4Y3NFx8ZEPD0MtraqKPEyl2K8eWkhUVtjmeGNKVeg2lBhZDUtmyYek/xETQdzeP5Dewba1g+8M4A9esWhXuJaKG7FdyVRFhpXztiAi0Lxvb5W1/+ctfvqtngDiliaruLXUQASOsDgLVS4eliauXis+8YlMNPcly0jHHHPMfqNXn0JExUmF6pHr98Ic/XIfq9fvnn3/+FhwyK9MI/POZZ06YNWbM3LJYbN7oKDpuIjsGYVMKiHiRXV0dfI5Zs79MJ8goyEh9e7kQVBakFEPC0pekkfX9rBciXWnOR5//jAr1CCnjHYC4vsU24pV5D6PyrV23LtiBu0MRNq20bx33xHjomhuJc/ajlG3S6j5dOZ1YGmF1Aiw79J1FIE1WZ5111j8T+eLWmTNn4kI1WkEQI5wxRVgIP/4R3OWX6Hm95JFHHvns4gcfjE0aMeJfy+LxOVMLC/MJw+Pl1NZG+JPgk7A7Rs9fkMQWdQLqI15WfgUnZylDWERfdGTlpmZDuqrlXMUTwzOraR/7G9ku9wbsXckLmUU7P5n80/+G4f10Ax5XtWePU1cxuD+Hof0eTttnZAUK3UhGWN0Az059+xBIk9WMGTPef8QRR/yMLHUrwaiBGOSEW5RoxCVJLiGDyss49o8lqGSTNLYPaSqb2a09QvXgpRtEjY3xOghoH7HGdqHC7aypiQ9jfOBqhtcUQFZOT2OpeGGNLOtZakKPPKQoR1hskyGdbl2F74lqUBrXMUyH8ZlfraqsXFmVupmWCz2D2axaItL5dSOszmNmZ7z9CKgnUEQ0GIL6H6IhRLgPhPg7xeWg2YKsdGeyCyo0dYTkpV67cJA81BmX6W3dGm+gR68aX6gq/N+2bdvmb96xI3iD4IjL6NGrZNYiTOrhZp1PIRKPIKDm4IeU4ia+dYTFulMVtR+NcV8QZK9KJm//9datKyfhD7/mrXbIpHWkAFo3kxFWNwG003sfAXrTYhqehCH7wxDUSIzYCcZmxhWWJz1gvPVdMGKgyZhNCJ1qpCufgee7mRhk54YNUcXatZrVKLaCmY1e37Zt7aZ9++5eXVt759Dy8s8gQX1+D6MQsIbHIKJAbg77yUUQUz7LXLJURhUOqSlSYoLgh9lIV3992vOunst5i+A5dovPLPUwAkZYPQyoFdfzCEBWEng0gPxiiEiz1LgIF/QKHm6wsZfNwOTVK1dGVRs3+ptWrfJfW7HCW7F69Z4NlZUPrt2x45fVDQ2PUbQ4ydu5Zcu/QYoBTHiFIrRuIghiNdxUQh5MloQlgzyEpfuJGvHMUqRWZvJ+8IkouogJKfdDVuIyIytA6I1khNUbqFqZPY2AIwCIqpTeNl8Dw+U1rix/J6mFh0obiJDxf3fdFW1av/7pDZs3/2/F9u1/5vjNOkfschrEs6SJZCImwb2S0NbLMN5/PRkEY3dBXHs5Jo9lHkv1EOIrEdPQHSbu2FOXTN68prLy2+wKuUmpo45cWVrqBQSMsHoBVCuydxCAsEKNv5QHuZaMEPCIbe98nTSGr3XCUVOzFjFvRM2K//v97z/B/pd0jPyo/vEf/1F+UBKFQshK/lDpFOB9/iPcJe6CHD+CPf0C9LujcOwajKoogpOW+HfsW49s9/1FuKpvZpvK0cLIKo1iLy31RbBkCPR1BHw5gUI0CQ0wJqqqC8eDq4AjLE1Sq+0QlOu5k8qoEC8aWK5j+L+Yc1+KXnghS97lkJCM+M7lgAeHaw5KIp0YoxD2Mt7vdsjr/QX5+UfhMDoJ9XJSorFxaiXbtldUfMerqNgsm1WqgNblHFSo/ekZBEzC6hkcrZTeQ8AnZpSicjbihHkXvlanEjssUix8jdsTSWkojNRDSVkafAwhuf9IYb56AgmC+De2BQv+9Cfm2DhsAEQ9ichMEpM85UPOl+e8855nqeR6Ilm2FfXVHWA/vYOAEVbv4Gql9gwCTtvSjNqf//znL500adJcfK4i4qLHRE4yukNIHoPJnYQlElOvoVQ+SCwJicXwLF95//333wNhaXhOS9XvcHcoiSl9vO6jZZIUZupfS0TepnUjrLcJaLtMpxGAX/zouOOOy0fCuguv9guOPPJIebW7nJ7fER5CM6vwHn30UQ9ykpQlosF5vT4uInvjjTcuF+EhoUlakuTUlWTqXldQ64VzzIbVC6Bakd1GwIWAgYxyGSf4R4jqApxFGwnVomCHzsiOe4O7iKQpDc+5+OKL3X/GEPqvvvpq/IknnqgnX/aHP/xhcdpLvtt3ZQW84wiYhPWOV4HdQGsEIBgXXmfatGnXQlhnMgSnAdtUttwXpAqKpFomGdklcX34wx+OPvGJT+yh5/BpjPDXYHD/u5FVS6T6/7pJWP2/DjPtCRxZYWAfxcwyV5Fxh9KM73EXjSEtWbV86NS2xPjx4/1vfetbv3jttdfeD2EZWbUEKUPWjbAypCL7y2NgDFJUvebcxn27dxLCOhMDegFdfgqK58tWJUnqEMlNVorKeN4Xv/jFvLvuukt+Vl21WR3iMrbrnUTACOudRH+AXFsENR9vcpYa4KeZaZqztj3GPvYf9C4iUR1PL2BUXVcX7WUaNM1KI9cFTeTQTnKkhrPnKNTGPE1Oqo7Bdo61zf0UAbNh9dOK6y+3LcdKWEOSTmIhP3OHDSs8Li+P4Xiet6murtHfsaOGVSc6QV6xEyC3FwhnXHzbbcTVq/ZL8LOaqklmCR+zFRcGRfGU+0JbSfM7QmomVbUFToZsM8LKkIrsg4/hAqqLrN4/evSo46No3qgoOqvY82YVNzbmEvnAxzBV96ny8tcSvv+nVb5/t795cwXPkfQXLvROmzx5z/unTfP+YcYMb+z48V7R9OneMmamwQPdPaoM8LJrKckfS0N1NB8gXvCNeL4zLtlSJiJghJWJtfrOP5NUMaePfWzUqC9OjqJrjvT9IWVsLCYX0ctXSCZUSwmDicvR3U4vD8P5L5eV3XhtY+P/TiwtXTBn8uRPzpo8ORxSVBTPwX7F7KPeWEjqiaVLHUFJysIgr4ln3ezZSGNJptDS3I8v/uxnP6tGHQzoTTyk0YtbsdTPEDDC6mcV1g9u10lW5WVl+acEwV0Q1QXjMJaXEqplEAwyGIYpgWSYaTk9uUNI7PRoDLvKo+imLxCxM2fQoJKpuCkUcl6cMYEhQff2EsZY4td2pKw9kJfGDmoYjghLdi2NHWTcn4bh/D9hhA1LpGkpwxAwwsqwCn2HH8fHZuUM6+/1/Xsn+v7ZpWHYWIRRfVAQxAcTcK8UNa6Q4TOMoeFI7OyakkuG9IaGKNbYGJ0ahiWVu3YlC8vKYnEM7QcY6Ly7pia5EYns8YqK2GYecNjw4d4uSEyERYysCJWwEftV9vLly//v5z//+aPme/UOvwW9eHkjrF4Ed6AVLbIigF1ySlnZV0cHwdmDw7Ahz/OyNXFDHmQlosqHZGAaAksRCk82KNQ9RjBrUlIfi7nPlFpRAWMAt69b54ngtqDqvU7U0PsYZvNkZWVUPn58UtJVIROUQljqGYxhaM9GFfwrQ3D+xSZ5yOy3zggrs+v37Xw6R1Z5o0aNJpzwfxThNMXLlZWNZKScA/nk0tMX5EFhirGO/clJWbpDoi0Qia/pXvG5KsWIvgED+isQ3bLc3Nhje/fe9/DOnX52Ts77q1etim/iXNmvFCIZdXAnQfx+zGzK11LAASZ5cPazpsLsN9MQMMLKtBp9h55nNtLVEtwTBkfR5UVBUJAdhglGG8flXOVDPCKXuNRASVZp0tI6ZOYIS0uIyklbqIYjsWs9g8vDz+rrb3i2qupmPda4ceNOYVDzTAhqFrarXYSWeQlp62/s2qHhOim/K0Q2S5mKgN4nS4ZAtxGArJIewfEgqQuZ0l2KXiDmIKaLyyGkFUFaTg0UcWlaLhFWWkXUUtukJmKXyoKRJhFmXWQV4Zc1n/I2bNjwNG4LP2KM4GWElbkGsvodl9ghm5WRVbersF8UYITVL6qpz9+k3qNoxLp1YxBzxod4mCMrBfJJl0NUy9zuk0jCgtRSOThAeaiSY38xZsx0/LLC6chinAvjedIK0ln/09FDTbICjExPqnhLhkB3ERCZePFkcgSEk8tcfhGE45NdmE4mGvXqyMyO7BXSIxjIyH5AR6RIiv9MbOr56i3kGDGPzkU6i+eGIZ2MzQketDSQETDCGsi13/PPLlXQzZlVl1pquixlzZq8F1LKh6gGSe0TWYmgpCqy9PGz0pTxzoaFEV4m+P0cwxbn/Fl1/PHB3V/9qjds2DBHjnPmzAnNMbTnK7Cvl2iE1ddrqH/dn2ZtT0I2QTX3XUPWUlO/50FYORBWFoQlA3kh/33Zslj3WV+Gv1Ulx2YToG84232iiO6F2LYWFAR3Y6Oat2hRI96gHPFmkt2KiAyBBjq/udXWMhkBI6xMrt2379mcFEQc9TcY37erzveHV6HRDYaOHFmh5sm1IS61j3tK8p/xhE2+WZDSQ5DVM0RimFBe7hURoK9m6tSoYPJkf8vzz+/82po1z3lr1oRfuvzycVNmzpxGiOTRTE2/H3eGFyG+VSpOM+oYab19lf1OXskI651EP3OuLU0wRu9dDZOQPkHEvQt3hmFyEGYt+v68LCSsGCTlyIr/DazXQFhyDH0DdfA2fK6GpyKJjmB/qdwgjjvOG3XGGbk3VVb+695k8uQZ06bNY67AfGZmdiGScW9IbNy48Xmii15/wgkn/Blpy8YOgl2mJyOsHqxheVlTXDB9+vSIrnZJHf7ixYsDGnI0UCQApJ7b6Cm8CHXQ3wYAOCq4rj0fkqL30Gsk17OeD1E1QFjLsVs9wjjA6aiBxL7y9pFHMrvzaPYTA6sQ+9ePjps1Sz5YHmQVMmNOiIQl7tO7ewoS131EGL2V634J0tLchaYeAkymJiOsHqpZ+QLhZa3G4tSjVLGSPNx/2VsWLFhAD/3Clvt76Op9ohg9e8AMNo+Wjxz5oB+Pn7sD51FesLjYRQ8tsjoAWam3sAAVcS+kxGBobyTE9TyzOO/imN2oiPuJyrDh9de9HUhhJ59ySrKY/4STCYgmGkBWzhUHPGULS/I/mjJlytVLly7dxv9vG2kBYgYnI6weqNzU+LXkRRddNO6UU075CNLACaWlpaMZPpIg1O9LjHX7A43pUS5Fe8p41cWn5+9fcU94iW68oZVEaYDJ4vLFkqvCfoimiKVsW7iNenkQ17+WlnqLiGeF23q4C6P8rmTSzyos9N91yine4KKiWA5OpcyY4wY7c4pL4KklmmYY4WOaYLKKa//yl7/cx/qKAYBxEwgD8NcIq3uV7kNWTmq69NJLv3bsscd+bdasWYPHjBmj7nc3wwtDUt5D6JN/I6jcw6+88splNLS1GdygJEgFhHnZBIHMg55/jy9V8VZICzeHoAY1sYQswmIkoeJhMc97FGUx4Pnc4uKoJC8v9gSS1nKC9J144omM4MnjaPFf4CZI5fC3JPY55iISaQ7q4aco7stSwzkwUyXZt2AwkDY48XogPXBPPitqYCAV75xzzvk+U1LdOHHixMEElkvQeJKoMEnICuHCk5SVZO68s9/1rnc98vjjj0+gkYUpe1dP3k5fKUtEEUM1fAzyOAXSeoLnje/1/WA99LOaXj1y4g2y1tewbRn7niEiQ11u7ksY388lCsPTcdREwsaE9Aa6eFdaHiKJtCRpzWbpz5kzR7hbykAEjLC6WKmyWWlWlsmTJ/8rZHQFRNXA3HgREyDEIawYjUfDRpTjNFwtG7HBHIH0dc83v/nNbOxZbMrYJMKIVVZWrqrYsuU0nv/jMMrTqIXJPeBSAYG9QX6N9WUQ1wr2rUomr3pq8+ZTNu3e/RAE/0I9qiGSaaiJJxjs7EIgM9i5PcACXCp8jj3iqquuGowUyyVtAor2wOrP200l7FrtyWFRkkQxk3xeizQV4n8UV++V1Jh0rPF00Sl7C16SXiP2rWPOP//8j7Ptp4899lj89NNPT0jaQipo/niw7oz1anjpMvrhUqSlZwqRtn7F8lf08h0JlUznoUrd0Bts7FDQsl0VFa/p+eZCcovI2KUKNKEEIY9drqqq8nYSsE92LNQ+pyLq+HSS9KVjON7fg0ppKXMRMMLqWt1KYkpgWD8d9WUk0pSmlNI258V9iCLVgCOI7VMsf6bjJAmQwrZ6D1Um+/qzeiNSl7rmiEsSF+vKrZOOEVlp2QieO4jR7uFj5WbIEVGlPwSacALcCewghwkNSTwgogohLLmPrLrjjjv26gMgTN0B9pNRCBhhdaM6aVgnQSqM222M1HBwZnRZE36yr62SNdmn9JWpV1999RCkq5066Ctf+cpxRx111CyktZmokzuwfa2km/5vNDoX54nyWO230pakxDTpChSRkrKS9imLXGSkcqSPhPo7VLyvIpkFxcXFbsZnxdMSrlIRtU0B/MDEERbSWARh+evWrXs0VZ7KMcIChExLRljdqFFIaSR2Fl/qi7KkAqkmmtFFUkFbSZMl0Oiyb7nlll0rVqw48dRTT70BO9hZU6dOdeqOpq+C0DwIbOe55557y5FHHnmzyEpE149JKw3F4UjEqZGEO34Bj/nFSE5zmNarEbJinorQkZOm88JO6Fwc9FHgYxGCaYDzaA29sLdzIddrm76gLTMLASOs7tXnVn3dCSTniEpLNSapLxji3RCSlsVLOtDceRiRd55xxhmfPP7447+HV3wBLgARtpmkXCGYnl2n+EhbQ0eOHHkjQevmLFmy5GK2EaHFOUtKIsn4BDl/BpJ6kYB9RaiBjUivMbAO9EEAG0ldxAOMJZHEYqiOmi3nynvvvXdzujMk4wEaoA9ohNWNiseO8hwNydf0UrKrpBqSU1XUo5WeO4/G56QDNTZJCE8++WQRktQPjjjiiFw5l6LixDmfCMKOrNwd0VjVTd84duzYc08++eSfUsZHICypOmn1qht33qdPlRQmj/nVkP4HwOF/IKRJ6ilM2bQiSa/q4GBfHEzr2f7FZ5555g4jqz5drz1yc0ZYXYNRpOHzhX+ML/xWGsxwHEM1xk3GXuc3JGlKPVuytaRUF6mCEQ3Rh7CGMpBXjS7CaB9PE13LW+Ec6YEEOfAT48eP/zAN8sesPzpAIhM40kI1XAJhn4SEdQ3q9kVIWEcIP3ATpoz8Ce/H0P7/IK2lRlYt357MXTfC6lrdSi2Loa7tQZ37JtLUj7G1NEAoWTQuYtHVe9VIUgqVAiGlCSvEMB+g3nkaDH3MMcfoypziuwkatGyd0tuQvCKI7wr2a8691odl6n+RVoyIDFUsvzxp0qRrsFVNJedKSkX9Xst212lhZJWpr8Bbn8sI662YdHSLpCx5dP8E0poFSV2OHcU5OcqWJRUx3R2PJBYiiQVIYn+HsP7G8f+G71AI0flSHdnnDO1pgmp1A1IDxWbvoTexiGOqM8QA3+ox2/zrJFn2BGAr161XWx6VIqpIDrwtt9t65iJghNW9unWqC6R1BSS0CZXvKxs3bSrG3uLvwYlRBngR0769ewOkq6fodr+Q/18XSXGsvLgD9S5KYpBUJmN9WylFalkY4d80crV1YGZukzSbJq6WLhGhEVVmVvihnsoI61DoHH6fGpNyFqR1E/5D1WfMnv39S+bOTbzrhBPi5WPGJHds3Rr85Be/uPfm73znYqSjkN7AOCQVySNbWYZ4VERnsIf00r2EzVduxFHS9ULu2hVg02nePgBX0lgPwEe3R04jYISVRqJ7SzWmGNJTQw6uB8X8GYZNahC6WxFhUk6bMWP1tyErks/wlOeQrK6ku95P+2vJzqXhPJzvehYZk+hUxAZ6xqp27oyq1BO5deuOlStXEs8uI/yxuoe2nT1gETDC6rmqb1JbsEcxTsRLMKwEk7pXx3g4RKnVIhpdCmJ6FGKqhrDy8WpXLyETIwfNXtwaPxdqggZUxBH0hmXv2RNWM+xkdUXFLx584IEDDJpWnR0ydEHPPZKVZAj0LQSMsHq4PmLYpwJsUgGGd7r/RFYeI3JFVtEn58zJRa2rQPX7MVLWl9SziEE+m55F7wAExTZv44YN3qbVq72ZhATO0Qwy6JFTNAh4y5adEuPmkBeSLRkCAxEBI6yerHUkJU0S6kNYvgiLhLHKa0ytnzp1avKOJUti9AouwJdoFvarM1etWtWooSX4GQUKqfLUk096IUN8Bp18sp+PV3w8Pz8YnkxGM2KxWx688MInT1+4cNntxx+fteXFF6M57gpv/vBfnBbCjlpaMgQyDgEjrB6sUo3sJTCTk7AcYSE5RaiHUhGVppSVOULByF6DLesCVMG7UAEvYByct23bNmd8r2S5BSlrVEFBmDd6dJA3aJBPRM5wRENDQXl9/e9vOv74Ez/74ot7VV57khYXkSuEEZdAspRRCBhh9WB1iiViSEnNEhbrkrCSKQnLW7xYVxNp+YRa2c/ygww/uRSXh48zvOd4/LRKsWntr8N4/+SyZaWDDhyIMHT5OXjQxzF8Tdu/fzKi2A83jhhx/eRY7LSSKJpRzHTuBZSJVLUfc/yTa8NwCYU7h8q7uZ15mT+UR5haGiAIGGH1ZEXHYpEkLF+GdtQ5jFJOumrU/4OTIy1twqZ1J4s7GaozGBtWOWrivnxIaENV1eJn1q4dXzRiRJifkxMMjqL4U7hDLK+t/ej0WOyjU1A/h3HiYLLio8tBC0v8F8f7/s5XRo2649F4/Fvz8MSXtAWZmWMl+Fjq/wgYYfVcHfq4NGRJwgqQqvAC9Qgp4IWshziRtnEZkZaSBLMIA7xYrSWzzVu+e/eTxWEYK8rNjd7A+L47kfAncOwk1ocwhZamyBpMT2Q+S+dRSk8k7hRDKedLwxobP3j8yJGf9rduXWKkBSKWMgIBmV0sdQ8B/7wU6dQkk1OTSFX0CnrZUgNxCGU8jnegqiop0hi2Y0dbeEv6CckiNe3XUnP5Pc8MyVe9RA/j4wcOhBuxh02EnMZAVhBVNBj/iNLs7PhgBk/n5OXFg/x85RiDFyMGLzaWxWKTpgXB/U+OGnUGZen6bV2bS1kyBPoPAvYSd6Ou5jeRQPQAkXqPLyvLP7WwcPR4pB0iyvmEcPDwtfL2YHDPqqsrEGnMWLGiQcTVziXZ5YhLy8RpkBbOWT/MLSx8tDgnJzaO8wvZgbroFRLJYRChaAo1sJrYUAThIrp8MfrhYI8BjL5XUJCFy3xySCyWNz6K7n1sxIjxnIpLq5EWOFjqxwiYStjFypsL8SyERCaVlmpeqq8i/Vw6dtCgcdiVonziVsmvCrEprjn3Ts3Pv255fv57KsLwZlS0F0QcEBiLdt0PgiWUnTdq1OiJnnf0cGJjEcE8yIGsciDEXLzisWt5gdROwtcgVWHEQikkKicjqT1GYCOn+TEoKlEOhe2Koh9yvX/ggiwsGQL9FwGTsLpQdyIrJkxIvnvUqJPPzst7dkYQfGOS540biQwzHFIYg8PoERDIRDzVx8Xj/hGxWNG0WOwfZ8Ziz7xUVvYNWCNMkUebBDI7JbmVhuEXCM8wHCpK8mWBgbgwhKVhPFmahEGExTAeJ2FJupKUJWlL4ZlFZKiMEFeyPAjOe2j06Nm67t3tS3hdQMJOMQTeXgSMsDqJ93zIRGQ1fdSo01G3HkUCOhLDeAKKiEogqlKIqhTpJx8yiSl4H9mXXSkeTwxHlTs2Hr/+1bKy7zvymDs3SwH5NGwnlQNN/bU4ihQuuRCi+ijBzOWyIIlMUhNrgcdYHs+HtCCkJulKElZa2hJRaZ17cFKX7F1IeQVh+DEVMcykLMFgqZ8iYCph5youWICU8t3y8jFIU3djAM8rQO0qwA0Bu5JXBIkMhkSyRBYiE5GKLEcE9SNGcpwYMhG5cVpW1hWvjhjx2jGLFn3PWwT9vZk4uGm2l7FTpsz0a2vR6FwJmkfMDSCUhT4BcSXJcRGYsoYAKSvpmun/WuKdiruDprGY/flJk3LOaIor1XRshv7qozKdZ25NznPAFrSEsbKlfoiAEVYnKm327Nm+v2SJN9Xz/pvJCIfmI1kh28SdbQlyyIeksiXtpKUckZYIheE6TGMs8tAkg3FmUwgJ4H7D9ddcUzJu+vR3jRw2bMjg4mLNV7ibuFiP4gX/A0IDj544dqxXV1+fJGJfvIayDkCKhGvwKMk7gL0qh55DyLCpfPa5a2lKd9myRJTK9Coe4HjWhpfW1xewPOC2ZmCjnQ/CC3hWEAeA9hPPL980HcOqpf6EgBFWJ2pryZIliSLsVgWed2FeGCaRX+LMdOoknTiElY06iPrXZFeSbUnracJC8lELYSCh/1hZmf/CuecWznr3uxeOZ5CzZjPW9F6yTZHeRwDAT1/1hS+svPOOO7wihudMZv9UpLbxOKBm0fNYCyntI+djYI9j3Jea6IhLZ4sctU3Gd8gM1whHcPpXTORTHZKJaS7QLuQZyd7lI0eOmxiLnVDieccS4mdoFlJwbhC8SnifF/6hsvIlyAqR1/O+CXIcn7GYZGI9G2F1slZRAS/l5Y/cjKipc5nepsmuJBVMvXWSrNJ2JW2DRNQqCDvq/RGj+J/f+17vPTNnRsOGDEkSYsbXxKBEbWi+E+JhHXnDjTceeRxzFU5IJIIxkE9hRYWXs2mTl42f13omt9jBdXIgplKuHZOUpesqScKSNMe1kqxXQ1hIaN4+GjS+XJlos/SRrHyIJzm3vHzWkcR/H+X7Hxrn+wUjgGMQHwyses6XpBasVpWXv7TP979zUkXFnSIrRCzZB4209O70g2SE1blKwmzkvwd7kKZhdnPIQxVu3EvaruSkHUk8IqBU5svuBRDKCxDXD0aO9M4pK/NK8/P9gsLCeGuySt1OmA/hXfixjwX7X3rJKyLSKPGUvQQ+XfuR3BIM+6mCkJhSx0tARsUQWg7SmWMjSIoA8l4j+6u53h4iPUBWfhVurCtzc2tpnFKZuP2MSGmyCi8eNeobkz3v2iNx5xgG3kPpuCjhw+IIC5xE56W4m/DCH9cYRb98Y9SoTy4Pw08x7nKjkVb/eReMsDpQV+nge9iVioivPl7khIoVkD3Zh1yGOJjp1CtC2gkk8UjSUWZb2s70I/yjBk+a5GxdcaQwRRrVlFVtJAlwXg77dzOZRRwn1FrIr5oD93LtTRBiQNk+pNTAtWohqAK2uZJorJKs6iErqY5IWOFebGfbPe/lBzC43w2Nzmvi2DYu2782zW1S6ZJn0euKH9wV43hWVEANWYoVk0vBrFBErg8IuDnbXjIZZjFb9IQoOhNVcfHDY8ee5W/cuHZ+U1kmafXxV8AIq5MVBP2EWIicXUjGb+X9kEutVC8IAkbzCiCSdA+h1EXYx6tguM7L7DseguILD49REgmJzS1b/2h7No0tgYvCaiaq0LRhm5CwliFdabziZMqVdFcPYVWznsfxkiIC7kMSnYhsv+6JS2wlbwyCn+sat7GuZX9PkJXzhZuOX9voWOwKxlY2YDXMUo9tEUTF0CWvCLsfkxim3TuaCIvZi1CXA/XWjvG8Iw40Nv7h5+PGnfjJDRsaFqg6Mkf67O9V3Ob9G2G1CUvbGyECNfbYfn5qyLWpXAP57IMc8iCsHMhE9qhcERKSTwiZyDC/DMLaCWExt5ebBVqz5BC4D/NWo6Zdp6S3JpFaQFnPrV0bvvTss8Ha9es9n+B+s2iE29hXlyIlGqqXmyYstodkkSI9islaVKQVyeRzj5aU/HV+RYWMzIh9/T45sioZM2Z6SRhehx9cEgSzwMDPFfZkjQTw0721sieyzUlY+pg0jQTIIppG46R4fMa+Awf+g4q9hhrjIIts0ZffDiOsDtQO7cCJQ6hwdODVbqtn/PFe2ES2oSKIgR5DLxfyyIZ8RE5KJZBXPl95SVhKlTQSqW61SEua2osIo27WHM2co17CthJxshyp3btkSfDiK6+EZUQfHU9D3J4iJCfZcaLGF4ryJGHp6iIs1NVQxLoRgnojiq5YwTjGFU0Nkq39OyFdeYvIRYnEv+OsG2SjBvIi+/JL00gASaZxfQTkQKuxliIutjn1vMm9pIm8EpxWXx+iPl71u7Fjv49qWDnfVMM+/XIYYXWiejSZZ3l5+TOh70/eHUVhMS83TcHTTDl4pDvbVQQpJVjfD3kV8lVnLKE3mP2EOfa0bycq3XAiOewuKXETrWoqe0lkmnjV2VpS9yMJTNN/7d61K9x/4MDWvNzc8v00urWUNYyclvJojh5jFx1hqTJFWLKgQarx7Vx7fTL52b8xfjGtQqWK788LN9JA4yzB/EPxppEA8qtySaq0RgEEIiipg5KuRFxaBydAbrIpyvWD2YiQeMMR9CgO9f0PMurgJxNuvjlY+OKLZsvqo2+IEVYnK4YG8X+c8glJV9tZSUs2slNFkJJ67eTUKYN3Ok5VlQiFbTiaehWEQB46bJhHD6FTBUVWUgslcWmmaP3XTNBMrhqxzX/t9dd3PPvss+MZqnMjJPalGtSfxry8mLroRYSS8KQSynlV+oxaGsb/+K4w3IMkdvWyysqfz0a+QCLJBFWQp3OdoWFhMnkSfmu5qHlJnjmmh1NnSLq3Vqq4+wCAp9TBSEuw8umtjSCziKXUbWLwR7l1dVFQV/eBefPm/Uh0hvDsoJRkTfKZsJXNnrd8+fJo4cKFRmYC4x1KRlidAz6Gg+cjhIx5jtZwEoSQ5LutpkBbiPSmO9tRPcRUQ+ORIVyEpq+/pKDjaUSLMZyXErNdPYBqUJoxR7YsqYaStuQ8qvkJpQ5u2rw5oicw94UlS676r5/+dPQQpLKnn37aX0UM+K2yfXHlPDIyBK01IUu+GuSeWF7efduC4LrqiorVKckqU8hKT+oS4ymnizl4sAj112ska4mFyvXWymEW6Uu4eBFZRCUbVpKPgfzW1MNaD17bx48P9owY4dcXFs64b9680/74ox+9AFFJgPU0zpN19W0ouyQCYyVIbU9ttcXbhYARVieQ1tAcvN0bmDL+CkjqaRqIv4VmIEVDjKBGo0aAsdujt8oRCYTmRAK60L2TWH8BVXI9TqCiOchKU9g76Uqe7pqqXhKWiKeKKKWbILYFX/nK4PsXL/7PYUz5deYZZ3gfwzdLcbZeffVVb9myZV4VRKc0bPjw6NRTTw2mTJmy/pxzzvlESCPlfuOL8M53B2TYD4xxhB4Mxc7lOpZyK1HPaC2EpB7bXMhJ8zuqp3AdzraPbdni5eIHdwS+bFWMMNhx0knekJkzg+FjxniTiovH4hqyZObMmWu/+J//+cNvfvOb30Xiavjyl788curUqWM16S11hZbtr+dSksKcFsp/qt/S24WAEVYnkIas9KWNQRgvMN7vCsjldklTm3mB+bLH1GhqyUVkPOIdYfFdd4RF4wqlV5yEWPUUDWrTxo2ahzDJrDkR6p4/CA94CAvNJebX7N8fbK+s9K/48Ie9OlSWe+67LznrmGM0dXSsEFXy/PPPd5niWiYVr/s79qXnnvvKzGOP/c8FCxZ4p59+estjMmYdZW1PIx+Fap5ImCurtxaV2dsHvnmQFYETPflirWBI0zdWrvQYWeAVsX0xYzSPfNe7vJmMJDiivNwbMXw4UXoYZklMfoqZQP7Pb33rWxfxAdg5bty4OcwjWSQHX6qunlDWL9Np8gOI6lccJ8maVSMtYfF2JCOszqMsUogz682PMcA38qX/HpETCrZi/KXxJIuRugZzQCEZUy9hQ51XuY9tJbaXBrY3mfxfpK7fNyYSC5naa6ZmehZZyX4l9wapgwmkgtkzZyZ82OsHP/mJ30CkBzmZSvqSyigbl9RJZRoLV5G2E7rGw7aocNCgf/v+97//wzlz5iDsZWaDIibPCrHLXvK+VJYtLw/SygHDLLCRyjgIzH7Cx2EJmJ4AkWWB3ynjx3uTkFiHgvsQlgUsSQJSWadFSKqnZKFGqkcXwopUP/yX9n2K8saNG+fef//9H2UdwS4zMebZ+lwywupalUgbiTE1189xSXgeErkR5ng/Pk9xuRrsZKfebKcOQihiOAzxa2GPazdVVt7BX2/u3Ln3PfXUUxdivzoLV4mZECBCmvc6xvZJxQUFp2NEj29lOM5S7FVTkARkmHdkRmMUSYm4WiaRGUlSVoS0MHrMmDEzOO4Z2WHYplvIlCRCUU/g06jUjTV8E3bz330gwD4b4tYD+5C67FgamrOZ9T2o3q/rQzBqlFcMQWWjfqvjQ7lVEoZKyeFIXmAYaLwnSdvEkbp+CL4XnHvuufcuWLDgfWR3T6n9LCz1FgIHv/W9dZXMLFckEENKWsbyA6iIxyHdnIcP1LtRDYeIuEhJXCBeoOdw8e5Y7AHm9Kpjq179gJ4nnf/bVNaxzQkSO65y587PDF6xYi6NpiQBWaE++mT3xddSNq92krsvyG0G+5+RutnOcf11s8gh2LRp0/LysrJXkkFw/A46P5CuYvpAIHl56M6OWWRPVKfHdCSl37LOpLReAUSFwb6Z9CWltpPou+BYakvSa+qDICy5RCRjfAPq4tkXXXTRFaz/twIvnn766Yl2yrLNPYSAEVb3gBQ5OKkGCekl1pXbTXN52XnjdY6yXn6dq6X7j5E8kJ2MBvASjeIy8r0Yeu9HbQyRvJjucK+3Gx8u+WepMbU1DlEqjHoc5SaRwUm4qbPjVpjj10SjiCrZoJdZYPLRcL5wB1jKxngiqvbHsV89A3GJfNKSqqRW/W+PtLRdxyIBO7WdyzgCE4mlLhdiU/wy6vdP52Sw+q2H7SvJCKv7NZFWB9SIXENi6cQrlmlS0v9wURMxseqStomo0imCrFxZNAzZ6kOkozVIU40QVHwXDqeSqlD3mofySGWRK4QaltwjRFY7d+6MYV9hroutkvxEbul7SV8nE5aSZAI+Er9ByvoY/lTv2xmGjYCdJQDVW9sAuUjK1SgADVMqhrQOQD71GOLlRiIpVaSuddmn2koiM/XiirRaJ+pIrg0RhJXJ6nfrx37H/xth9VwVqK2kyStdqiOq9J9OLNXmgjfeeGMtquYrkNXxFRUVScgpJkOwkgzvkrhSxmBn40ISCEVaq1evrnjxxReX0aicwOFOyNAfpKlPIiY9i5p3xA5IC+DiuDr46rGtIcv/TarhVnIDNixJn8JNy/RHQJ0draVVSVNItR4zcXt8ON6CXkrKkjQsk9l0ciaq32957nd6gxHWO10D7V9f0loC28mtfOV/vX79+pCGpcbhyEkSghqTJCyRWEoaIGRWVfZrr732X9jIahYsWBDHM/ut4kH71+xPe/RxkE/aDjo+zmQc4Z1ImqdUIRXtB7c9kFYxByA7BVVIvDALjNXgbduxIzYcaVU9s8JORCVyEinpv4hIEpWIX6o3NkpvwoQJlESBTaqgW9eP6kDkh6TW9BVp3mMrvYWAEVZvIdv9cp3as3nz5t/QrX4JDeh8iKgRySqGGhPo6592NoWwEDSiBA0oe+XKlc//9a9//cH8+fMDyKqlytn9O+p7JTjSglTWTZo06XSe/5uIlJc3BkGxhk1VQUS5yjCNXnQZ5Tk2ZNwmUzrmOaIXWUntk8QlHzcZ19N2q+3EIdNHQcQkO5ZcS9Lqt7ahfgfUDzEAK18XNBmqfuvR+kwywuozVdHmjTj7Ew3oI7gx/BrSugDScoZ3JIIkzowRjYhw8jGGwtVl02AWL1269OOoj3WQlSS0TLRftQbKkZYGprPjGsj9h4hM74eIzmLDONTBifzfyb6NGJ0akbL+Yd26dQlIPi6ykmotMlJnhiQsuYvIHqjB5ziJus4NeiSdeiipTNKsCA2pKpIURn3sYdTBK7opXFVamwS02VIPImCE1YNg9kJRIhwfqaCW5QdxVL2MhnQlDWg6X/9Y2pmUfWvYfjs2r++xroabNv6zOiCSiEL2ugCy3szy9lT2iBJb8q6NG/ctSnVwQGgPQTZnS1pFsgqQlALUuuaRBhAW2nXoOi80BGrixIkeuDZ3eIjQUlJZI+dlv/766z/95S9/uavFuEMubam3EDDC6i1ke65cR1oqDkfVH7H4MY3uaAhqCpJBPl3zf0edWcX2A7Kx0JgGGlkJGiXhJBVYz5+WLkN6TKs2soEk+18EGX0Ite4uSOsCyMYZ3ocOHRql1Gsg9GPs8zZs2OAmCIGgAklSafWb/7pOoyRaSG/1c8899y2p34w7NOlKKPdyMsLqZYB7qHg1EiU1uiRShFQQp4Zoo5IGOqfGOg70hqPnb4mBJC8lkZnsTdjkvQ9C+h+D9C/HBnWSpFUZ3yEjBiSEz/IhINJ0+B6pidgKQwa7hyIsJFo8KGIasJ6N9LUKwns/PbJ7XnnllTRB6jqWehEBI6xeBLcXinaNjnKd+pMqX41TPlyZ2hvYXRjTZK9ytO4IDNL/Feu/wlv9SIhrkKQobFfbUfPWs13jB7+IunglEtoECK3ZSM8xRLpuvGPVqlXf4rA9ZJFVS4Lkr6XeQsAIq7eQ7b1y1eiUrZF0DeM0gUlaDVH9pE43J1RqqYURhPad0aNHSwU/ETI7khxBVJvpNXyebTLiKxlZNeHwtv0aYb1tUNuF+hgCaZcPkU5abRRbpT8EMdRF+Z/+LZWbb5/ewBh+bjoufWzzPlvpXQSMsHoXXyu97yPQHumk1e+DCI3HCVMD1/v+k2XgHRphZWCl2iP1GAJSH9OSWI8VagV1HQF9PSwZAoaAIdAvEDDC6hfVZDdpCBgCQsBUwg68B+o5knMgk5H606ZNS/cydeDMrh/CwOUevQ7GZN2MytRKd8vuiTJ0P5YMgU4hYIR1aLiIWhJEREmo/8Y3vtGecfbQJXRxL2MBu3jmYU/rLlnpAukyHAse9oodOEDDXSDpHiuvA5fs9CG99BHp9H0M5BOMsNqp/enTp6vxhDgN5jAq/z033XTTbga9MqFwPKHBsekkpx2SpsSrSm/rzJJzQ8pLE4CcF5mVKptZuhK7ta8zZaWP5fzm8rSNAb4+Xto4be9788bTB3dyyXjGiCFCPoSqZz7oOp0s6qDDeVb977HyDiq8h/704kdEoWvkUuGen+t0tt77NG49BL8rxgirDTRTY8OSkNRoiGM0ToMVLMuQtjQV8CBOcdEm06dKOuhskpqpF5SUx7mKEe4KgRBdbCuGgDDjvaIAdyyly6OcGA6OB9km2aboA/VIim4+vfS1OlbywUcp6B3lBNdffz2TTwcaaN2ppAfmhNYNUngSIDSq5pmdVNuZQrkPeP5NLqYspX7zEdFHgAHuiqmVdfvtt2fxQXBOrR3FgCFE/nvf+96IKBSKldanpdSOPlN7xxlhtUJGZKUv3A033DCJBnQxBLKOl+iHNABHHqkXqsdeChp/EbcQJ2zJQeQkqajVrR3yr17UlGRVxH3nc3BbpKCGfMhyWu+ECESogznvoPK4RkDucGFpIqIcNcwSzm35vCqbibL9AsWj6kzSPZAHcU6//YgwHCj58ssvn0KMrmpico2jHrPBoiU+7UKi52cYUZLxjBET676HuGAPtXtwBuwwwmpRiWmyQnqYweaLySsINfx7vTx6MbT87Gc/SxTeHk27erQ0YtT1cHl9vrj+/hGBpENCYRcjKe7iPXucHGO9XcJKf0RSx0QsfSJ2nEWMtF1HHXXUI6owynhLhNQ+X5EduEEjrBRILckKYrqYCl/OV+8eeTWz7siqBZ6SLNp9oVocd6hVV4bKPtRBnd3X1w3Xeh4k2O5i1xKWjPiI3HrrrcfwUJVXX3310pYP15F1tIHzmSdxN/H9f/7pT396Qxvva0eK6RfHGGFRTW2R1eTJk++55pprZBhpTVaq2J5ocK4MSW0qsAdTT5fXg7fWq0X1y49I2lUGEyNzwIY5ehdBSfkgFbwlculzCJUTI9rEhbxD09j/uyuuuGKlzud/u+e2LKc/rg94wmqPrAjI5siKSh2oBNDf3ueeqCdXxtv8EXFEe+ONN0q1i2Q/5Z2UFNoe6ej4UAOwZ82a9UHu9SiI7rd8XJel3+X+VnGdud8BTVjpCpbNiop3aqAkKyOrzrxCduzbiIAjN723uKlcxHWnQla/GyhkJZwHLGG1JisMn08QsO2xVKjbnlAv3sb32C41ABBoJis6Ef+R5z2Cj+wvIKvN6Xd5AGAwMAkrXcEtJKu/ffWrX32MCtdLYWQ1EN78/vWMrclqAurj//LObkm/y/3rcbp+twc5GHa9mP5zZrqCW5IVX6nHtD31FD1hC+k/gNid9nUE3kJW+AYOSLJSRQ0olfBQZHUII2dff6Ht/jIXASOrVnWblipabc68v5CVL1JqS7Iyssq8+s6AJzKyaqMSB4SElXKkCxkbeDS9Khfx/29pNdDIqo23wja9owikB95LI0gZ2CcMZDWwZWVkPGExjMGNMRNZsa7hNkuMrFq+Arbe1xDQlPdGVm3XSkYTFt2+CqtShxo4DbKaS8/K377+9a8/pmnFly9fHumlaBuW3t3aS3GV1Fng1Ihu3n1PlNHNWxiwp7v3keFgkqz0cTXJqtWrkNGEhfrXQDyrKSUlJaN57vshq6f1/CnH0FZQvH1/UUN762I90cOZLkPE1SMJFdzr62Mc3+mPiMwWfFQTAnzt2rVzWYw3NfCtr1+PvZRvLfod3eKkhC984QtXrFy58pKjjz76AaYbfwrH0BwNf1Bu6+7YzjujETlNSRIaqd/EVUrfd1eXrYLzdbUYO68LCNx8880fh6B21NXVZfHOHc36Lwain9XhoMtIwkoZ2SPZrfhancMMvpXEm8o+FBipcwYhlR0UV+lQ57S3L3199uex3hycL308L2SXg/NRxlvUWK5RT5ltRZVIX7IzS8W56rPB+bg3RX2tZVndmYdKHwtOB0V4ZYSDJBuFDKpKH9PZJerbQR/AzkR4VfC9GTNmhK+//vrHGAA9a+zYset5vru/9KUvVchkYZ1CB9dGRhLWwY/Y8X+9EVeJhhVv64Xu+F01hThWGTS2doPzIR12ui5pqO0G5+tMeZC8ixLK/bUVnM+Fu+B530K0h8MA7CjS7/WPCNcR4Th17HD3lN6fujfdYIxtb3k29nfoI6JywDokAN9JQ4cOXYM28CNirlUaWaWRPnjZ6Zf84NP7/L/OPN9BX8k+/2QD5AYHykekuro6C9vqSiIwSHoUUdv72MY73pkG3cbpGblJmHT3ZXFl6MXrSYT6uuFaz4oK013sWkLWk2W1LLfPrhtZHbpqerRBHfpSttcQ6DICA+IjkiL7AUfSXX4r7ERDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDwBAwBAwBQ8AQMAQMAUPAEDAEDAFDAAT+fy6wNVEc8WdBAAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "water_box = packmol(from_smiles(\"O\"), n_molecules=8, density=1.0)\n", "view(water_box, direction=\"tilt_z\", width=300, height=300)" ] }, { "cell_type": "code", "execution_count": 3, "id": "toy-md-code", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[01.04|15:07:41] JOB md_for_dcd STARTED\n", "[01.04|15:07:41] JOB md_for_dcd RUNNING\n", "[01.04|15:07:42] JOB md_for_dcd FINISHED\n", "[01.04|15:07:42] JOB md_for_dcd SUCCESSFUL\n" ] } ], "source": [ "md_settings = Settings()\n", "md_settings.input.ams.Task = \"MolecularDynamics\"\n", "md_settings.input.ReaxFF.ForceField = \"Water2017.ff\"\n", "md_settings.input.ams.MolecularDynamics.NSteps = 30\n", "md_settings.input.ams.MolecularDynamics.TimeStep = 0.5\n", "md_settings.input.ams.MolecularDynamics.Trajectory.SamplingFreq = 1\n", "md_settings.input.ams.MolecularDynamics.InitialVelocities.Temperature = 300\n", "md_settings.runscript.nproc = 1\n", "md_settings.runscript.preamble_lines = [\"export OMP_NUM_THREADS=1\"]\n", "\n", "md_job = AMSJob(name=\"md_for_dcd\", molecule=water_box, settings=md_settings)\n", "md_results = md_job.run()\n", "\n", "rkf_path = Path(md_results.rkfpath())" ] }, { "cell_type": "markdown", "id": "e8dfc9b7-61d2-4a64-ab42-2f78aeed74a8", "metadata": {}, "source": [ "### Convert Frames" ] }, { "cell_type": "markdown", "id": "convert-md", "metadata": {}, "source": [ "Convert RKF frames to DCD while remapping atoms into molecular entities.\n" ] }, { "cell_type": "code", "execution_count": 4, "id": "convert-code", "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "NSteps: 31\n", "Created files 'ams.dcd' and 'ams.psf'\n" ] } ], "source": [ "rkf = RKFTrajectoryFile(str(rkf_path))\n", "mol = rkf.get_plamsmol()\n", "print(f\"NSteps: {len(rkf)}\")\n", "\n", "pdb = pdb_from_plamsmol(mol)\n", "psf = PSFTopology(pdb=pdb)\n", "psf.write_psf(\"ams.psf\")\n", "\n", "dcd = DCDTrajectoryFile(\"ams.dcd\", mode=\"wb\")\n", "for i in range(len(rkf)):\n", " coords, cell = rkf.read_frame(i)\n", " mol.from_array(coords)\n", " mol.map_atoms_to_bonds()\n", " dcd.write_next(coords=mol.as_array(), cell=cell)\n", "\n", "print(\"Created files 'ams.dcd' and 'ams.psf'\")" ] }, { "cell_type": "code", "execution_count": null, "id": "6d3266d0-3fed-4708-a4cc-6899718bf20e", "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.8.16" } }, "nbformat": 4, "nbformat_minor": 5 }