Skip to content

Instantly share code, notes, and snippets.

@benbovy
Last active October 29, 2018 14:05
Show Gist options
  • Save benbovy/2bb908c745657712e1da213c0a9badd9 to your computer and use it in GitHub Desktop.
Save benbovy/2bb908c745657712e1da213c0a9badd9 to your computer and use it in GitHub Desktop.
xarray-simlab: wrap existing, object-oriented code
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Wrap object-oriented code with xarray-simlab"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"import attr\n",
"import numpy as np\n",
"\n",
"import xarray as xr\n",
"import xsimlab as xs"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"'0.2.0+4.gb7054cd'"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"xs.__version__"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Example of existing classes"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"@attr.s\n",
"class Mesh:\n",
" vertices = attr.ib()\n",
" faces = attr.ib()\n",
" \n",
" def advect(self, velocity, dt):\n",
" self.vertices += velocity * dt\n",
"\n",
"\n",
"@attr.s\n",
"class Object3D:\n",
" velocity = attr.ib()\n",
" material = attr.ib()\n",
" \n",
" def compute_mesh(self):\n",
" # create some random mesh\n",
" self.mesh = Mesh(vertices=np.random.rand(3, 3),\n",
" faces=[(0, 1), (1, 2), (2, 0)])\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Example of wrapping the classes above as a xarray-simlab model"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"@xs.process\n",
"class CreateObjects3D:\n",
" initial_velocities = xs.variable(dims='objects_3d')\n",
" materials = xs.variable(dims='objects_3d')\n",
" objects_3d = xs.variable(dims='objects_3d', intent='out')\n",
"\n",
" def initialize(self):\n",
" self.objects_3d = [\n",
" Object3D(v, m)\n",
" for (v, m) in zip(self.initial_velocities, self.materials)\n",
" ]\n",
"\n",
"\n",
"@xs.process\n",
"class ComputeMeshes:\n",
" objects_3d = xs.foreign(CreateObjects3D, 'objects_3d')\n",
" meshes = xs.variable(dims='objects_3d', intent='out')\n",
" all_vertices = xs.on_demand(dims=('objects_3d', 'vertex', 'space_dims'))\n",
"\n",
" def initialize(self):\n",
" for obj in self.objects_3d:\n",
" obj.compute_mesh()\n",
" self.meshes = [obj.mesh for obj in self.objects_3d]\n",
" \n",
" @all_vertices.compute\n",
" def _get_all_vertices(self):\n",
" return np.array([m.vertices for m in self.meshes])\n",
"\n",
"\n",
"@xs.process\n",
"class AdvectObjects:\n",
" objects_3d = xs.foreign(CreateObjects3D, 'objects_3d')\n",
" # this is needed so that this process is run after ComputeMeshes\n",
" meshes = xs.foreign(ComputeMeshes, 'meshes')\n",
" \n",
" def run_step(self, dt):\n",
" for obj, mesh in zip(self.objects_3d, self.meshes):\n",
" mesh.advect(obj.velocity, dt)\n"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"model = xs.Model({'objects': CreateObjects3D,\n",
" 'meshes': ComputeMeshes,\n",
" 'advect': AdvectObjects})"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAkMAAACDCAYAAACOTh8GAAAAAXNSR0IArs4c6QAAQABJREFUeAHtXQd8FcX2/iC99wYJCYEEQu9NOgiIghQVfPhU7OU9FcEufxErFrA8O3YsiICABRDpvUMKCQkQSEjvvcL/nLnZcNPbDaScw2+zu7Mzs7Pfvcz99rRpc/ny5V8gIggIAoKAICAICAKCQCtFoA2Rocut9NnlsQUBQUAQEAQEAUFAEEBbwUAQEAQEAUFAEBAEBIHWjICQodb86cuzCwKCgCAgCAgCgoBohuQ7IAgIAoKAICAICAKtGwHj8o8fefInBO9aWr5YzgWBOiNgZeeF0XesqXM7aSAICAKCgCAgCFxNBCqQobzsJFhYOaJjt3FXcxxyrxaGQHpSJOKiAlvYU8njCAKCgCAgCLREBCqQIX5IYzMruLTr1hKfV57pKiFQXFQACBm6SmjLbQQBQUAQEAQagoA4UDcEPWkrCAgCgoAgIAgIAs0eASFDzf4jlAcQBAQBQUAQEAQEgYYgUKmZrCEdSltBQBAQBASBpoVAakYBEpLzkJ5ZgDTa0jML1T47twiFRZdQWEgb72ljMTZuC1PaTHgzaQtLc2PY25jAzsYU9ramsLM2gauTBRztTJvWg8poBIF6IiBkqJ7ASTNBQBAQBJoSAvkFl3DmQiZOn89Q++i4HFyMz0ZMYg74miZGbdvAytJEbRZmxjAyaktbGxiX7HG5DYovXUJR8SVcKr6s9nn5RcjOKUQWbcWXrixaYGbaFh4ulmjnaglPd0t07mADfx9bdPKygbmZkXZL2QsCTR4BIUNN/iOSAQoCgoAgUBGB8zHZOHYqBcdCUhAUnoqouGzw4kpmpkZETqzh7GgJ/07OGNLfAk525qTFMdcRINLyNETy8oqQlVuI1PQ8JKflISU1Fyl0fCgoFeu3RiEvvxht2gDtiSD18HdAvwBH9KGto6d1Q24rbQWBRkWgYf8rGnVo0rkgIAgIAoKAhkBGViH2HEvAzsMJOByUhDQyfZmSCaujlx26dXbB+OEd4elmAycHC0VGtHaG3psTmeLNme7jV0nnSUSOYuKzSCuVhXNRadh+ME4RJFsysw3o5oyRA90wvJ8rmdxMKmktRYLAtUFAyNC1wV3uKggIAoJAjQiwr8/GXTHYdiAOJ8KSieS0gV8He4we7E0mKTu097ABm72akjBJ4q1XVxc1LDarXYzPxNnzaTh1JgWvfnoCxWR+6+nviLGD3TFpRDs42Zs1pUeQsbRCBIQMtcIPXR5ZEBAEmi4CTBR2HonH79uisedognJg7tXFFXfP6IGunRzJF6d5TdtM1jp42Kpt9JAO5L9UjNCzyQg6nYzPVp3GBytCMKS3K6aO8VRaI3baFhEErjYCzet/1dVGR+4nCAgCgsBVQiA7pwhr/r6AH/44i5S0fHT1dcScm7uhN2lYOKKrpQj7NPXu6qq224q6IDAsCYdOxuL5945SxJopZt/oi5nXd4AtRayJCAJXCwEhQ1cLabmPICAICAKVIJBOvkDfrI3Ams3ncRltMLSvB0YN6gAHu5ZvOmItUL/urmpLz8jHjkNR+GZNBL5aHY7p4ztg7ozOcKBQfhFBoLEREDLU2AhL/4KAICAIVIJAAeX2+fnPc/iKfvyN2rTFxJG+GNqvHcxJc9Iaxc7WDFPHdSYfoo7YdzwGf+68gHX/XMBd0zpjzhRfipJrOdqx1vj5NvVnFjLU1D8hGZ8gIAi0OAT2HU/Eq5+dRAYlPxxDfjRjh3ZQIfEt7kHr8UCmRAZHDfLCsL7tsW3/BXzz2xms2hyJ5+7vhZEDXOvRozQRBGpGQMhQzRhJDUFAEBAEDIIAZ3xe9u0ppfHo190NM+/yg42VmIEqA5f9pCaM8MF1/dtj7eYIzF9yCJNHeeKpe7rD2lJ+uirDTMrqj4B8o+qPnbQUBAQBQaDWCERQduh5bx5Cbt4l3Hdbz9LQ81p30EorcrbsO6YFoG93F/zyZxhmzduBd58ZQA7mdq0UEXnsxkBAjLCNgar0KQgIAoKAHgI7DsXjnhf2kBbIHM88MFCIkB42tT3s7udM2A2Co4Ml7lu4F1v2xda2qdQTBGpEwCBkaNfeI5hx+2M4ePhkjTcsX+E/T76K+c+9Vb64wnlt62kNT0dE4ouvV+GuB57DV9+t0Yqv2j43Nw8b/96FFxd/gIlT7zf4feuCeV2xM/hgpUNBoBUjsGrTeSx46zAG9vLAI3f0gbWYxer9bbC0MMFDt/fBdQM88dzSo/hhw7l69yUNBQF9BAxiJouNS0RcfBJiYhP1+67x+BItBnj2XBSMTYxpTR0KKuUFbSqR2tbTb5qTk4/IyIsIjziPkdcN0L90VY6zs3NRUFCE/QdPICcn1+D3rC3m9cHO4IOVDgWBVorAWsob9NbyIEwd3wnjh3m3UhQM+9ht6RV+2vjOlJPIDO99F6IWmZ092cewN5HeWh0CBiFDt82YhFHDB8DN1blOALalb/XK75eC9/pE6PDRIBw/GYr77r5F9VdVvepu1qdXF+Tm5mLbroPVVWu0a87ODph64xis+2MrToWeMfh9KsOcCeV/57+GD955XmHKN60PdgYfrHQoCLRCBPYeS8QbXwTixtG+QoQa4fMfPdgL9D6Nd78OhpuzOcYMcm+Eu0iXrQUBg5jJGKy6EiENYBtrK1pJ2UI7VRqml179H33J6VuuJ+Xr6V2q8tDIyGCPV+U9arpg1LbxcoaUx/yTL37G0eMhauVq/XHVBzv99nIsCAgCdUMgPikPL75/FIN6e1D+IJ+6NZbatUZg7FAvDB/QHos+PI7ouJxat5OKgkB5BAyiGTp/IQab/9kDV9IM3UzakMysbOzdfwzbSSsz7z93Yeeuw9i17yhcXZzw8H2z4OioiwLIzy/AngPHsXnLbry+6Amw6efRea8iNS0DO3YfRmJSGmZMHQffjl5l6rG2gyUjMxtbdxzA8ROncCE6Ft27dsIjD9wOCwvz8s9Z6/OLMQn45PMfUURkjK12148ZhrGjh4Cfcfm3q1FYWIhpN43DkEG9VZ87aZw79xxGSmo6uvj5YM7sKeQTYFnt/QoKCvHr2k0IDolAXkEBOvt2wK3TJ4K1SfoScuoMNvy5FeejYuDr44VRIwZhYP8eqkp5zD9d/jNW/LxBXXv97c8Ia0c8eO8slMdYwy41PROrVv+F02RGtLQ0x+QJI0ufiTtJSEzG+j+2IYS0WpbmZnBzd8F/H5qj+pc/goAgUD0Ciz85SS95ZrhtcpfqK8rVBiMwc6I/zl/MwMsfncAXrwxtcH/SQetEoMGqEyZBLy5+H19/vxbx8TqfoTXrtuCtZV9hB5Ggx596A8eIrOTl5eOPjdvxxjtfKKT5x/itZV/ipVc+ABMKsvDAkkjMjJuvV9e7EbGZPHGEIjbl62kf1UuvfoifVv6OF595CAuffZgclnfj7fe+1i7Xa9++nSsm0X15TGGnIxUR4o68O7RTBIOPNSLEDtoriVBMvXEsRhNR+YHGcv8j/6cICNerTDIysnDvIy9iH/kS3XPXTNw1ZxoRyb2Yc8/TOBl0urTJn5t24vmXlmESkZT33noe585fJCxfV3Uqw3wE+UV5e7VT7W+6YTSGD+uPyjDmCjGxCXjg0f9T5PSuOTcjOTkNTz67BL//tV215z+LXvufMl0uffMZTKHn079WWkkOBAFBoAICOw/H4+DJRNw+JQDNedHR9PR0bNy4EYtfXow1a9dWeM6mUmBk1Aa339QVJ8JSsHmPRJg1lc+luY2jwWRowrjrlDZG/8H5B1ZzWv737VPx2stP4LMPF8HB3pZ8gU6pqg52NorAdPHvWNrUwcGOtCRe6tyDNBF9ewfAx7t9hXpcgTUebBJyc3OGsbEROlI9Ly8PBIWEl/ZX34PhQ/tjQL8eymTH2hlNjhwNxh23T1Gn7Jj97Q/r8NLzj6BXD3/lH9S/b3elxfmLSFlV8vEXP+HM2SgicA+iEz0rt1347ENKm/bakk9RVFiExMRUvPvBN7hp8mj07tkFpqYmmDPrJkVOws+cR2WYdw/oDAvS8LAwbnxeGcZ8/cNPV6BfnwBMmzIOPbv74+EHZnMxPvvqF7VnZ/jjJ8PIMVFn4hs8sJeqpy7KH0FAEKgWgS9WhaNvN1d09LSttl5Tv5iZmYWTJ0/i0OFDyM/La9LD9fSwwYCeHlj+65UXyiY9YBlck0PAIGYyU5OKqwubm+nKWKOiCRMb1hIV0g++CUWQsZiW7LU62l7fobqyemZmpvjx67dhYW6ObIrW2rh5F6Kj45BPZidDyC3TJoAduX/7/R90C+iEsPBzuHT5EnoE+KnuN5Jpj/2anlm4tPR2rP1q384N2WQmrEq2bNuvzGH6/j69e3al5zBD1MU4nImMwuEjQeT8nYfrhvQt7WYYHW/5/ctSE2BlmHPl8rhxmT7GHNnGGjv2I2ITGQs/B4+7LdkF2YTnTgSznYcrPidyxOay++feisUv/kfVlT+CgCBQNQKnzqYjlLYn7+lfdaVmcsXTsz0mT74Ru3btahYjZv+hJZ8dxPHQVPTpWtbloFk8gAzymiJgEDJU2yfQNA21rV9TPdY0Lf92DZGWQMy+ZbIyZXF+IUPI8GH9lBlpy9Z9eOyRO7D+922kSRlf2vWFC7HK8fvrT18rLavpIJX8ipiMWOs5jHMbdvTuRpqcI8eCERUVp0xiXG5jY8W7UqmNL1RlZKi0Azq4QISRZe6/p2H2rTeq48r+vPD0g1hIJszfNvyDv8mMt+CJuZg4fnhlVaVMEBAEShDYuj8OLo6W8PFsGdmRm0IQSm2/XO3drNHezQpb98cKGaotaFKvFIGrSoZK71qLg5p+1FmbwQkN48mk88n7L8HW1pr8WrbVoufaVWFH4+lTx+OzL1diHRGCXfuO4JEHby9tbGZuqjRS7Nys+epoF9knhzUr5YW1MTy5JCSlKBMca2A0cXF2VIfsXG5CZj8W9lv6F5nH9CUrO6daB+0qUjWVdmFuZqaOQ8LOlpZpB3n5+cjNyYMdmTDZzPb98iV4/6PvlXP84jc+gbOTA9gUKCIICAKVI3A0JAV+HQ2vlcjKysSOHTuV1nnsuLHYsGEDacKjMXzECIwaORLxCQl0fTtOh55Gt24BuOXWW8sMcN/+/dhPW1paGjp18sUtM2bC0kr3spWUlIRNmzbhdHi40rS7urrgnnvuLdOeT44eOUJa5V1ITUnB6DFjMJY2Tdi/aP2G9Th35pwy148bOw79+l/RjtX2Hlp/9d37dXTEkZDk+jaXdq0YgQb7DBkaO40EcdRWdbJ91yGV0HBg/56KCHHd3NwCkB+2wYTzBHFCyI8+/wlDKXpMPwWAf2cfdZ/Pl+v8bLSbBgafxsdUvzLhvvxK2mm+U1q90NNnwKa/zp280aljB1X8F2WwZnOZJmxi1Px6tDL9PWNXXHxJbfrl+see5CDOGqat2w+QmSyy9BLnKHrl9U9Kkmcm4MNPVigfr0UvPIrnn3pAJcXcQVFzIoKAIFA1ApEXM9HOtaxGt+ratbsSRaSHnZg//fRT/LXxL7zyyisoppdBJhhL330XS5ctxSsvv4zkpGSEhoXh2+++wx9//FHa+YoVK7B+3ToKDJmI64YNw5o1azB/wXwyietcCt55522Vj+zlRYswYcIE/P33ltK22sHmzZvxxZdfksndhEjTaSxbuhThp3X+OfHx8XhqwQK4OLvg1lm3ITU1DS9RX3///bfWHLW5R2nlBhx4uFpTZFlWA3qQpq0VAYOQIc62zBIXf4WRJ6dkqLJs0mRowo65LOkUSaZJVmlb3TXWPrDspZD7UNJe/LxK95+6fD2NmPy9bR94aQoOLWe/HnZAZp8YdlJOTEpVfaVRqH59hM1wY0cNVk05nF5fbr5pDIWkW6ikjo9RlNfaDVvwzntfgXMkPf7InaVVk5JT1DFHdrE89vC/1X7V2s3KP4dPGJfI8zG48183w5ZMYzeS4zSnIeBnuO+RhWpZkaeefweLXv8f7r9Ll4iyMsxdnOxV3xv+3KaivzgCjUUfOyZkt82YqPyEHnniFXxJ6QIY44ceexkeRJQCKIqPhX2bUlLS1bHmu9SbnK1FBAFBoGoEMrIKYWvg5Ta8PD2xaPFidVN2NXjttVdx37334tnnnlNlp0PD8BYRmkcffRTvLVumyiIiItT+7NmzWLVqFebPn4+AgABFdnr36k1apYvYum0rEhISERxMWZxL8qH169eP5oAA1Vb/Twdvb3z80UfqHnfdeZe6FBISovbLiST16NkTkyZNorZdcdeduvlvxfcr1PXa3kP/fvU9trEyIb/RS2qrbx/SrnUiYLSIRP/Rk6IPIjstEh7effWLqzxe/dtmfPndamRSzh/OE5RFzsNh4efx+8btpKEoRsTZC5R/pyO++n4NDh0JVP2Enj6ncvK8+/43KiKMC09R1Jafnw86+XhiJ5Ebjtbaf+gkbiYS8vFnP1aox1FW3E8ERVcdDzxFUVcBGEF+PgeoDZOAtmSO+vHn35UpK/pivPLV4QixuoozEQwOsWcnYn1hcxNHZHGIfMSZCyqvEudH4jD/zp06qIiwl1//WOXp4XanQiMoeqwDenTrjC7+PthIofP/7NivfIRY23MrEZS5/56uHKBNjI0xdEgfnAo7o3DgbNycdmDJq/MpuaUTKsOcNWSXSLvDGrM9lOOpfXt39KUs3G+89XkF7MaPGYqk5FRwpBxH5B04HIihg3vjUcrRxJF5GRRFsu73rco8FhUTR8RqBzjC7vbbqvYx0seGj7PS45CScAadB9xT/pKcCwItFgEtkszdxbDaIRPSyPz000/wJlJy4426/4c21tb45ZeV8PbxwQ033KAwtSJT/OrVv9J8YYmx48Zh7do1NLeeIsITrExhbA7LIJOWlaUl9eWDAQMGYOvWrdi7by9pdFLh7++PEaNGUoCLieqPtU+s4Rk6ZAj6ElFiYa39li1bwASpKxEsJmDxCfHYs2ePusfBQ4fA86MJRcFOnDiJNPc21d5DdWqgPwkpOTgSRIvizvSDMYXciwgCtUWgDZlHyliWQvd9iITIbeg34tr9iPGQmFhxeL1mNqvqgTiSTNMScZ3y51W1q0t5Gml17MmPpiqJpigwHidHZNVW2OeJtUE5uTlEkrxVNFllbZVGiepy2oHaCI+Vo8LYh6omYRMcE0XOrcRaLk0Y/+KiYqWKTyT/Jn6umj4Hra22j7twHBGBmzDpwV1akewFgRaPwLi5mzFhREeMHOhp8GedMmUK+QN1w5IlS0r7vpV8g/z8/PD666+Xlt1GZb7kF/Tmm0vwyuJXEBgchF9Wriy9Xv4gMDAQb7/9tiJDFkSSHn74YYwZPVpVYyL19NNPY8b06Zh7j+43ISw0FAueego3T5uG0aNHYd4T83Avaaqm0XlVUt09qmpTn/J9x2OxZmMYdv+gI4f16UPatE4EjJviY/MPb2UOyJWNVZ8I8fXy5/ptODEkOyVXJ5yN+aNlC8tUqY4IcUVP0sLUVdhB27djzRMm5wqqi9Q0Vv2+2HfIr7O3fpE6ZvzZnMZfjvo8W4UOpUAQaCUI+LS3QUxCdpN5WlPyQ8zNyVFmMQ6V1xf29XFxcUGXLl3xwQcfYPmXy7Fj+w7lD+To6IjevXrpV6/02NxUF5ARTn5E5SWfAzJofUhbW9sG3aN8v9WdxyZkwbt9zS+C1fUh11onAk2SDDXWRzGMTE/6eY8qu49ZSX6kyq5JmSAgCAgC1SHQr5sj/tgRU10Vg14rp9iv0HcnX1/s3r0b33//PZ577tnS66dIu8NO1XeSf8862j/00ENYMH8BEaDeihjt37e/VmTI3cNDRaDt3r0HM2acoUg1nc8hj2vp0mWYecsMRYgaco/SQdfiIPxcCkYOqBjJW4umUqWVI9CqyBD7LvEmIggIAoJAYyAwdog7vlkbgcjodIPmGuKQeJZs0vJowtFgrH1JSkzUimjZozyVeJZ9N1kmTpqoHKj37qVlkxYuVNFkkZGROExZpd9+621Vd+fOXZg1axaZ4h0wcOBA1Y7D81kSE5N0e/Id0iQpWRcok0xlxuTfeNPUKVj1yyo8Tw7dbDqzJFPbXvIfYodtfz9/xMbGorp7aP02dH8xPgsX47MxdohHQ7uS9q0QgQY7ULdCzOSRa4GAOFDXAiSp0uIQcHEwx64jCbSCerZaksMQDxgTE4P3yYwVHxennJ9TKM+PB2lk3nzjDRVen5WVhXNEcHzJkXrZe+8p8sF5fzj30NChQ9GnTx/KEXQU586dwyFybk5PT8Pjj8+DT0cfCnjJomCOv5R5LCYmFn+TY/SgwYMxnUgNR5l98cXnqk4C9cW+iByk8fnnV8rY13DG9JlITknGaQq1Z9+gY8eOKcfsu+fOVUv6VHcPQ+Cj9bFuC6UnIb/v/97RVSuSvSBQawSapAN1rUcvFZssAuJAbZiPhnznEZ+cS2+8OWqLTcpFakYB0mnjfWp6PrJzi1BQdAlFtBWWbHx3XiRUbSa6vbWlCextzeBoy3tTONqZwd3ZgrL2WqrN1dGcHOUNM+6W1Et0XA5++P0sRWgZw4Pw8nCxgLuLJdrR3sLcqMKj8kKt85ccxrx7BjSp9clYQ8MfsIf7FR9HFSxBUb8c0JFMGh93ulbXYAkNgFzSSsUScWPTmaWFhVascpRxZLEh7lHaabmD6NhMvL38EF59vB8mXCeaoXLw1Ok0KTUfUUTmo+h7H5uQg5SSeSY5vQBpNN/k5JXMN4WXSuady2qeMaV5xpg2U4rk5v8rTvZmtD6mKRxK5pr2rpbw9KCNsoQ7UnlTk1ZlJmtq4Mt4BAF9BJLT8hEWmYHT52ijPa9xFUOTUfElXcCnuZkRTSLmtK6cKSzNTShYwByuTjaURNNEhRFzdnNj2nRLKFBEYPFlFPKPEO+LLtMkVoisnEIkpRfifGy6Ok5OzUUBTWosxsZt4Oluha4d7eDvQ06vvNGxnY0uzFp/rK3p+OOfwvD33sr9gOxsTBVBcidi1M6VSJKzpSJL3f0d8MO6EDzz4CD1Q9EU8GJtUnlRwRJk6mKp7Hr5+tWd8zqRvuSjVF4MeY/yffM5f89/+j0Uvbs4ChGqDKAqyhi3iAuZ4PX01HYmHWejMkpzNDG5cXawUPONlYUpzQNEZNztYG5KwTX0omXEL1s037RtS8l+iUwXl7yIceLf3PxiNb+kEIG6QHNNZnYBklLz1AsbD8eSXiI6dbBDQCdbBPjSnraOnjbUVxWDvQrFQoauAshyC0GgMgTik/NwJDgZvITD4eAkXKQ3MRZHe3PS1Fiju58Lxg6zpgnJjN6yLGgZlsZ5m8rIIr+TtDyk0GTF0ThRsVnYfyIJaRn5ajwcnTOwhxP6dePNUb3xqQut5I+tddVkMD2TtHS0hZ7TJSgtD8nTb+7A9cN9MHm0+CqWx8ZQ5yuIdMYmZmHOTT3ALxSskRCpHAF+yToUmIyDJ5Nw9FQy8oi08EtWe3cb0tjYoE93D7jS2nrOZO61szEsjmRhRVomzTWUC4rzQV2My6KxpGLtlguUu+oSrdlpjP7dnTG4lzMG9nSCz1WOChQzWeXfGSltIAJsJgsP3Iix//4NZpZX1mBrYLfNvjlre7Yfise2g3E4S29l/Ibl3d6W3pLs0dnbHh3a2SitT1N40Cx6mzsfk4GI82k01nRcoGPWUnX1tcfYwW6US8cdnbxafhjzHzsuYtH/jtf7IzEzNcLbz46qd3tpWDUCW/dF4be/w8tUYLNvT9LM9aKtp789/EjDaUTai9YoTEBOhqViy75YbNkfQ0Qkn16qTODn7QB/WkOP5xxXJ6trah7nOSUuMRvhkamIoC38QhpFIBahHX2O1w/1wDjaWHPU2CJkqLERbqX9Mxk6feIPXL5UBAtrN9i794CDW0/a94SdS1fKEN44Wo6mCDf7+6zfGoXfd0QjgbRBrPnp7ueMXl2c4UskiP16moMUFBTj9PlUBIUlITg8mTQi+aQ2t8SU0V64aYwnvVGaN4fHqPMYz8dk45bHt9e5HTdgIpRPuE0d3wnjh3nXqw9pVDkC2w9EYc2mcPTws0dQuC7arrKarPno1sleEaNeXRyovkOT9FmpbOz1LYtNzMVvpHFZty0KyeQD5EEZ0Xt3dUHPABd4kRaoKQv7SV6ISceJU4k4EZZIBC4XHmSCnj6uA6aO9Wo0zZ+Qoab8rWjGY9M5UG/EwJveQ1pcIFLjg5BGW35OChEhY9g6+cPerUcpSbK0qzkBZXOCg992tuylbLh/XyAzWLJyXB7Q0w19u7s1+cmotjhz+PiRoAQy9cVRyHchhvR2xcwJHTBqYO0zsdf2Xle7Hjuih5HvFv/IBoWnKp8hnqTrItcPa4f5c7th64E4vLU8CCMHeWL6BL9Wq6WoC3bV1eXPYf3WCGzdewFP3NkNt9/YkbQJGQgkDUjg6TScPJ2qovmq68OTtUdEjDQNUmdvmxbxubB5+6c/zmHf8QTYWpthcB8PMjm5wc3ZsMvDVIetoa9FxWXi0Mk4teXlF9H84o45UzoqzZ8h7yVkyJBoSl+lCFQVTZabGYPUOCZGRJCIJGUknSYH30KYWtiT5ojIkdIe0d61O4xNm99/YF4kkrVA3647g0R6o+nd1VVNSF19na6pc2DpB9MIB+yIGRSeiP3HYpXGqEM7K9w1rTMmj2ivnLIb4ZYG75LfpPmHlIkP78Mi05UfQ31uxE7oz97fQ/k+aO13kGl04QfHSJNmi7kzuzea/5d2v5a6z8ktxDdrgnEuKg2L/tMH48mEUplwpCV/joFEjHgLjkhT/jGV1eUyTXvEmiM2rTFJ4igoQwr/PzFqpPXS9h5LxGe/nEYIPWdXXwcMH+CJ7v7OLYLgaZ8Bv6AcP5WA3Ycu4hy9iA3s4Yz7b/ND3wBHrUqD9kKGGgSfNK4KgarIUPn6ly4VIiMxrIQc6UhSTkasCu+1duiozGo6ktQDNk6c3bZp2v5ZE8RaoM9pQsohe/eQvu0wZnAHOJEjYmsStv1v3X8BhwPjKITfDI/+qysmjyq7DMS1xiM3r1hFz+h+KHUEiB1vqxN2ouYV6asTUxMj3HmzL+bO6AyOxCkvHLkz781DyM27hFtv8EcvMluI1B6B4PAk/PJnGEUwtcG7zwygH/3a+5Hw/8+I85lKa6Q+d9IiRZP5ujoxpPaIHZcffGm/isycP7c7RvQ3TJZsjgJb8kUwkb1U0nQ5YyI563cgH8SWLqfPpWLz7kiKvE2llw4XPHVvd3jTS1hDRMhQQ9CTtlUiUFsyVFkH+TnJyqTGZjXWHqUnhKCoMJc0RZakMeqmtEcO5IPEZjZTC8O8FVQ2jtqW7T2eiGXfhNDkmo0RA7xooU5vWFG4e2sW9if6a8c57D8eQ2H69lhA5qLeXR2uCSSRF7OUloC1Pmz2OhOVqcKxqxoMO5h270zaAfJFYQ1Bd9qforDj/756oKomGNTTGc+QNqiDR/UTMueEWvbtKaz75wL6kcl05iQ/2DRSlGCVg21mF9gEu3ZzBEVAxRKx9sRT95BmjSKPGio67VGqcjBmLVLImZq1R+p7oRyz2UHbXuXrqs043vvuFH7YcLa06pjB7lhAz1FfPzsm5h9Ryoe1f5+nIAZ7zJjoRzl8mrYvUOnDG/CAgzvWbg6nFCQUTTjFF/ff4qe0fPW5hZCh+qAmbWpEoCFkqELnly8hM+VMqXmNSVJWaqRK5mZp206PHLFzdhe0advwibLCGCopyMopwpIvg7Bx50WVbXgKOck6Uwi8yBUEeNHS9f9EIIQcrm+d5IPH/x1ATsUVtSZXWjTsiH8kgshUEFRiImHzSGZ21RodzpHi62WjiE8P/pEj4tPRs2KEHPcxbu7f6junP0IO4553Vzd6I2+nX1zj8T4i0K9+dhIZmYUYM6QDxg7toJyta2zYiipwuPU20jL+s+8C5dQywnP392rUdcdYexROGhwmRux3xBokDn6oTtgkqkWtsQ+SXwfbSs3h9764V5Eu/b44MeFDs7pg1g0+lbbRr6t/fCgoGQvfP0aRncDN4zqjP/kEtWZhH7K9Ry/ij+1nybRpgjfm9auT1lDDTsiQhoTsDYqAQclQJSMrKsgi7VEwOWYHkoO2zjm7IC+dbPKmsCVCxFojnfaoJyxsKvcrqKTbWhcdDEyicOsTlLDwMmbf1EVFh9W6cSuseCQwHr9uDFMZaV+jLMHdOtfexFEVXMr0QaYnRXxY60M/Yudjsqqqrso56zZHH2nEpxtpgDgBXG3k1id2gLVMLEyiZk7wxqO3d6Ef6vqRb052+fOf5/DVmggYtWmLcdd5Y2i/dpTUrnbjqc2Ym2MdjlrcRxrFLXs4/0yR8j/jt/7GJNFV4cRJA5kUKXJEpjU2S3FunqrEwtyYItfsKFL0inO2FZGeUXduormi8nac5PS5B3uqdlX1y+X8o//pyjC19l3vAFc173DyVREdAqyN/uG3EESQP9l/5wTgXzfVLbeXkCH5JjUKAo1NhiobdE56lDKraea1zOQImkCKKM+RUykxUhFsZGozMqm/BufnPyOx9JtgZea4dbJ/k8kLVBkmTamMkzv+uD6UsmunYOHDvTF5ZN18iXiZgMAS0sM/UDX9MJmQ3w5n0u5JodTsFMsh1Zwlur7yHTnFf7giVL11PvdAzT9etb1POmmzeHHXNZvP4zL5xA3t64FRgzoQcTRs0rvajuda1UunJJ87DkVh75EYtQba9PEdlP+VoR2ZG/J8TMDZ/4cj106WaB85S3x14uZsjvikvOqqKHJ9y0RvPMLkmshTeeHAjOeXHcW+E4m45YYuGEpRYiKVI/APRRluoGhDTvnBJLO2OaYqkKGw/R8h4sg3ld9FSgWBOiBgbGKFiQ9sr0MLw1a9VJxP/kahpdoj1iLlZSWQGc0INo6+Zcxr1g4+Nd6c38ze/ioIq+lHa+rYzvQm36HGNlKhLAKcBG7D1jP01n8e95J9/6FZ/mUrlJyx1oSzOmvmLvb1iaN12aoTXhJDn/h09bU1eA4n9jNprB/nbDK7shP+D3+cVWvOdenoiEG9PJSjNRO7ligcIRRIeasOkT/QqTOUgoKWN5k9uaPSulWX+bspYcHO91rkGmuQOLFqddqj6sbuQrm65t/dTSUa1Oqx6XfeG4dxNjoT983uBV/PhmtVtb5b6j6InO2/XR2sIs7enN+vVlrFCmSoIDcFKTHH6Q2FZi0RQaABCJhbuZBGplcDejB807zshFKzGpMjJkvFRXkwMbNW4fycFFKXHLI7lV2ZdJgI/d+HxylnTCyFRvcgM4tk1W7Ip3PwRBx+3BBCfkQdlXP1RXq71pm72N8nVa3RxgvPViWaOUKZvEo0Py1lGQYOwd55JJ5IYzT2HksAE6FeXVzpO+eErp0cyUG0ouagKpyaYjknoQw9m0yfczL50SQin3LHcI6qqZS4cyTlqGouSUirwpY/P15jUEWt0Xd5x6EEIkdFVVWvtPy6vq54+r4eKsHgwy/vJ9+lXDxyRx/KFm1ZaX0prIjAhdgMfPLDCfSj0Pt3nh5Qo19WBTJUsUspEQRaLgKXLxUjIzlchfaz7xEnh8xOu6Ae2Mreq4QY9cCa/W5YuSUL98/qRf4uTi0XkKv4ZAcpkRovZtqOVoBnMlSV8EKfHDbLxIeju3jfuQpH1ar6aK7lrInauCsG2yhx44mwZJVywo+ylnft7EwY2KE9RRDV1gxwrTDgF4no+AycpcifU2dSEE5ZzJkw9PR3pGVd3DFpRLtGyyp8rZ5Z/75THt5ao1ZTv752zNnLvSjDeyyZ2J64u1+zTpyoPdPV3p+lfEQff38MU8Z4qdxf1d1fyFB16Mi1VolAYX46aY805+xAJF4MxOrgGzBw5O3oQ46LIoZDYPfhiyp3jH6PbB7RQph1mh97Cj8XR1E2l+whTdHOwwk4HMQL6RaoCLSOlFfGh0wn7dyt1WKbTrTSOPHHayZJqbmIic8ibUaWSo7HmcrZbGRrY4IB3ZyV9md4P1eVc+eaDfIq3TiR1gKb/OCWBt1t7i09KVpVclLVF0Q2w36x8iReeKgXpo3zqrIbIUNVQiMXBAFQdFI2bp+/k9aV8sQNozsLJI2AwKc/n0D42VQ8RmH3Q/u4NDh5WiMMsUl2yd/NY6dScCwkRWXOjorLptB/3Xpo7Vyt4UyrjztS6L8jpXtwsjOn9bjMKT+PCcwp4qkhkpdXhCzKBJ2ankerxOchRW25SCQSFEv5Xpj4MBlr72qpovY4QzBvlaUsaMg4mkPbrfvj8My7R+o9VBtrUyz67zBlKq13J9IQq2kNuwMUobjqvVFwd648iELIkHxRBIFqEHjg//YjJZ0cGO+p2eZcTTdyqRoEOJ/Mks8PkvnLDkvI2VGkfghwxNEZSjVw+nyG2kfF5ZCGJhsxiTm0WOwV/ys2q1kRKeLNgvyPjIzaqmUijHnfVueoXUy2rSJKZHOJzFm8Z58XTn6YRRtHVGnC4e7uzpbgleJ50d7OHWxoRXRbteclLlq7/PJXJAVdBFcJAzuMs+M//0Brmwedf7EqHJRVAE/eO+CaavmqHHgzu8BzzFtfHER3SnuwZEHlc0zDXhGaGSAyXEGgLgjsOZpAb97JWEATUslvRF2a16tuSnIy9u3fj8OHD2PAgAG48cYbq+1n4cKFyo9k8eLF1dZryhfZQXjahM74/KeTlAU4vcZ8K035Wa7l2JiYcP6mynI4se9RQnIe0jILkE4bm9g4pJ8zYvMPBUd1aRs/g7FxW5jSxs7M/PlYkjbJnsxcdvTjbU9rdtmRKdPVyaLFr/7e0M+TM2afjSZtGTmN65Md7biy3Emc0T6CFp5dcP+1I0LHjh3D0aNHcYDmos+/+KKhMFzz9vwdnkoJKtlcxqkROOVGeREyVB4RORcEShD4bt1ZSqbodFXX+klKSlJEiMlQr169q/0sLpNNJDIyUtUpprd3fsNvrtLDz5lwtsGK9Wfx+ry+zfUxmuy4OR1AY6UEaLIP3QQGxsuG8KK9dZGVf0SqxVY7eFT8wa5LP/Wte4m0gomJidiy5R/k5GTXt5sm147XbmvnaoVfNp7Hiw/1rDC+5jt7VngUKRAEDIdAPL1FHw1JxnUD6pYYsKEj8O/SBddfP75W3XCU1Weffqo2QxKh4ydO4IcffqjVGAxZaXi/9hSGHKe0FYbsV/oSBJoLAqytO3AyEYP7XN15Rx+ftqQGnzBhArw825NG/NpQhMaagwb3aUfLu8RUujbhtXlSfeTlWBBoggjw2lG88niA79UPo+ekkLUVSysr8GYoSUhIxLvvvEMZgK/4mBiq75r66UmruHOyxSPBKTVVleuCQItE4Ehwssrx193v2i9A3YaIEC87c7WlMecgzg/Ha0qyOb68iJmsPCJyLggQAiG0wKdPezvlWGpoQAoKCvH7H78jLDQMBQX56Ojjgxun3AQnx7LEi9XVBw4cwI4dO5Cfl49Ro0dj5MgRpcPJz8/HwYMHcYi2J+fPLy1nn6P9tKWlpaFTJ1/cMmNmGcJUVFSEdevWIyQ4GEWUZ6l/v37U70jk5uXh+eeeU+327yPH8eQU3DD5Bvj7+YPNd5s2bcLp8HBYmJvD1dUF99xzb+k9DXFgZWFCzqSWavXwkQMkhYEhMJU+rg0CnJWa80N19LJG/25OtcqAzCPl7NUeLlYGT6yZlZ2F3bt3IygoCDHRMfDv4o+77r5b/V/WEAoMDMTmTZsRHx9H80YnZGZk0CUdGYqjsm++/oYc6otxmTTSo2i+GD58OKKio/HTjz+iqLAQkyZNQr/+/VV3Nc1B4adPY9PmzYiOjoK3tw+GDRuK3r37IDYurso5SBtnQ/YujhYqopKz2/MSPfoimiF9NORYEChBgFerdqL/OIaWrKxMIi5P4ij5BP3r9tm47bbbsJ3IzqOPPopTp06Vud2qX3/BF+S8mJGejqPHjuLtt9/C+vXrVR226X/44QdU9rYiTFrDFStWYP26dZg0cSKuGzYMa9aswfwF84l0FagquTk5ePqppxEXF4v5Ty3AtJtvVvd459131cR44+TJqp6fvx/Gjh0LVxcdKXnnnbeVyvzlRYuUCv3vvxuWO0Ubb/m9C4WD17TWU/k2ci4INDUE3vwiCO99F4LHXzuIcXM34/HXD2ElRZZxhF91EpuYSwkoDT/vvP3W21i7Zi3mPT4P856ch63btuHTjz8uHQqfv/rKKxg5eiRefe012NnbK6LDpngWdzd3NR8wyTl7JkIRIS738vSEs7PuJU4jQjXNQVv/+QevvfGG6m/x4lcRFRWFF19cqOY/ftGqag7i+xlCnCkPV2xCxaV9hAwZAl3po8UhwKrUxlg9/Gt6uzpPTs+PP/EEvEkjFBAQgCeemIfsrGy89957YK2NJn6dOmP58uVqcnr22WdV8Q/0FpaTmwsXFxcsWPAU2rdvp1XH2bNnsWrVKswnLRH3y3b/3uSEHR19kSa/rareCmp//kIk7rvvPlhaWKBv377w8/dH1PkLsLOzg3dHH1WPJ7+ePXvCniZFVlsHB4dQ2LXOfNePNEkBXQNUPUP/4ay7vEaXiCDQnBHIzC4sHT4vP8LLqrxDIfYz/rsN0/+7He98HQKOGtNPecANcii6z8zUsAYbfhFirY8rzRlGxpTV2ssLnu3a41RoqBpjZmYGPvvsMwwcNBADBwyEqakpZsyYAQvLskt/DBo0CH1691bzAWt2NDl5MhAzZ85UpzXNQRwt+yn5OV4/fjy6detG9zLB9BnTVUQst+X5prI5SLuXIfamNMfkUK6s8mJY1Mv3LueCQDNFgH+UeckAQ8vuXbuUOYzJjCbdu3enRHjmiImJwfnz57ViDBk6tPR48ODBcHJ2RjKZq6KoTpeuXdU1MzPz0jrbt28Dm9Zeozc7TfLJ9OXh7o6cbN0b6WYydXGUmpnZlRXR+a2xmNTf2lugaqt7IVSHbBJzc3PD9yu+R2JSIubMmYMFTz+l3cKg+yLC3NSk9j5TBr25dCYIGAgBTm5ZlURTcsyVf55TG88z/bs7YVhfF9pcVRqDzJziqprWq5zJzUekBTKn//OsGf6HtEAXY2NQWKIt3rlzJ80P2fTyc2UdSW7j4+2Nc+fOlbnn5JtuAjs3/7Vxo3qJOnPmjJpztPmopjmIteBsjmfipQkTsJW//FLGZKeu6c1BWl1D7IuLeI6pqAcSMmQIdKWPFoeAqxORk8Q8gz4X+/CwVseinMMzR4J1Ie3MiZMncfFiDEzobaky6Ubanl1EptKVLb9ijYvkC8Bvc8uWLat4kUoS4mnBSJqIrKzLOlzrku6VnRzalPgKaB09/vjjyiS3kSbBHTR5PvzwwxgzerR22WD72MRsFBWYYMWGs8q2z6HJ/MPCeysL3bEVHVuaC2EyGOjSkcERsLaq3U+rpjVizREQrL7rRsZtVO4nzo1jKGGtL/v2nCAiw6ZxTzJvnSNNDMuFC1Fq70YvTfrCL0fl+chg0g4500sZE6j77r1X+RGyr5AmNc1By0j7zWJjbaM1UXs2j5WXincvX6N+52mZ+XB2uPIyqPVSu09Mqy17QaCVIMCZdA+cTDLo01pbW6tcQKzdYdMTa1w0cXLS2d1ZTZxdRW4PCzJrsXi0u2Ia09rz3tTMVL35sVnMk8Ji9SU+Ph7GJrr/7kcOH0FxUbFSmWt1WJWuH0arPxGxtqlLl6744IMPsPzL5dixfQeWLV0KR0dHMsNdeZvU+qrvnjVxick5iE+8jCByYK9OeKJmcqQ2RY7omAkTJQe0tDCijYgTbby6PdexIPLEiQO5nIkUl/OesyTzsQXtS9wjqrutXGtFCHBkY25esTKp5NA+l0wrbF4pLcstpjQQnLiS6pB5ixNYavtwygReH+EM3ywHjsdi+MCy/4fr0x+34f+/S5YsQWJCAt56awmsiYhs/ntzaXeF5PzMEko+i2wGq054jph8w2R89/132Eha5oMHDmLu3XeXNqlpDjIx1s1B7Hs0Y/r00nZ8wNop/chY/TmoTMUGnOTQ55Wclgtfr7JkjLsUMtQAYKVpy0VgYE8nfPRjqFpssr2btUEe1JgmAt+OnRAeEY6g4CCMdR1T2u/piAhlq/ft1JHs+0Gl5foHgaQ58vDwQPsqyFAnX18VMfL999/jued0Pkbcnn0D2Kn6mWeegY2NLdhHYMeO7Rg7bpzqnifLtymcnt/0tAlImyC5AhOpddT+oYcewoL5C5QfEhMjjjgzJBmKoNXML+kt9aAGV8UfTjjJPxzqxyO5ikp1LGZipMgRLVGhO26rlqtQx2TOMDNrqxZGZdMGZw5mnzLes1nPjN7i+U1eHVOZOqbszaZ8nfYmVIdV87pjrkdlfJ32QsJq/qA4OzaTE86WXcDHtLwI79U5lxUWq+tcp3QjX518OmftC/vm8J4zQfMxr5+mNi4rOWayk0vHvOUR+eFlSK6VGFIrtG/fPhw9cgRTp05VRIifKY+iU/n/EEvHjj70F+Cs07Nnz1bH/IeXZLlUUqe0kA6un3A9fvzpR3z99deUE+36Mr5FNc1B3bt3U11t3boVN9xwA72I6DRCHOW2e/cemmMerHQO0r9/Q47DzqXSSx/Qy9+hQjdChipAIgWCANSq6bxm0IETsZgxwc9gkNxLjsvPPvsM/tjwO4ZfN1w5ELKWKJoiKubcMQfWVleIF2t4NGFtTjy92S1atKiMBke7zvuJkyYqB+q9e/fgxYULVTQZZ6g+fPgQ2C+I5V//mk3Okp/jfSIzgTQB2dra0iS0CxMnTFJ+Qdn0dsZy+MhhDBoyCMFBwWB/pZ07d2HWrFlwcHDAwIE6e3+3boZ1ot5/PA7+vrZ4/fF+JUSHFgQlR1R2ZueN37zVmzif87FeuXoj5zd3Kq8toVIPqvdH+1FMQ4FeaeMfGhm1KSVG2vIXas8kijZjMpvon5vQOa8jxuVqT+2NuB7lhOEyXl+M+1QblfFaZJwvRtvzjwFr1qgIbegPm0I0fzEmZnzOP5Pa7yD/aKpzIqrMVfmc1yfjVFS8540xZ38v3f6S8rfjc86Mzvsi8tNgcsGkho+1pT/0j7UyjfQoAlRCeBr/UzDMHRhn1kCiDZH17IpOutXdpRNpK26Z6I0f/4ykBaIzKPGiR3XVa31N0yizaasXaXJDw8LAvj4s+4ko+ft1UY7LISEhiuCw2Yu1RGElDtZr1/6GMWNGqzrchrXXPHdtpxeqSTdcMZHxtZrmINb8rF69RgWRzH/ySQqpH0bmunOIOHsGH/3vQ+5CaZx5rz8HTZs2jYsaLEcC4yjVgTNsaTmZ8iILtZZHRM4FgRIE2G/ls5Wn8RKtGs2LWhpKODcQa1bYNNaTJqeDlEtozJgx6q2Mf5Q4/P7jjz9R4fQc1cXOjOzr88ADD6JHj+5lhvEERaXFkuM1OyCyhNIE9uYbbyI5RacuYVPc4xRO26tXT3Wdf8h+Xb0aK3/+GZyniPu+9bZbMXuW7o2Qo9nmzZsHJlHs5P3SSy+pOhyOy0RtADk7xlJYPkeksCbJUJJEK56/9tF+/N+jvTF5ZMPMA0xqlMmCTRxsuighSZqpg80a+loAdVxiDtEIEe+VhoC1BCXag/qSLENhJP0YFgFef03T+rEZlTV9vGeTqc6MesW8Wmpm1UysbIol0qPzX2NfNja76urzKH/7JwqvfXqyxgEzmR09yA23TvJBv266RIs/Exn63w+hWPSYYeYd/j//yuJXFLlwINP2FFrv0MzcjFJqLKccP95YuPD/1LIbb5EpjfMGsSmsT58+yqxmamamfANZA6RvwuJ5hqPCOAK2vNQ0B/FL3tJlS6FFpHFOI46M1Uz7lc1BPM6GSmJKLl79aB/emNcP44ZWJJpChhqKsLRvsQjwD+K0/2xD187OmDW5i0Gfk01TnF8jlxyqfXx8VDRZZTdgp2cbG+syqmj9ekyGWJvD+Yj0JTY2ll/3VSSZfrl2zG/sCYkJcHXWhdtq5bznyTOBtFCurq5KY6A0ARRtxmNOptBYd3K01DQJ+u0acrz8l0BkZOZi5dJRSo3dkL4aqy1rMTRTi2ZuKSCixJoMNsewySafzTVswuFz0mqwSaaQ/LM0bYcy8ZTU1cw8SgNSUqbTjlzRnOjOK9em8OfS3ETTgmmaLU3jxeZCJies9TKl8G82E6njUnOirkzfvMjHVzbddc1cyXsOoea9ZtbU7UtMnER4WFPWWLJlXyyeW3q0yu5dHM0xfXwHzLi+A+UVMitTj+edqY9sRa8AN8ycaDitNEeS6YfLlz/nQfBLlAVFqDLx4YAP1gJVJRkUyMGa5aqkpjkonfKn8Xe4snuUn4Oqukddyr9aFUQLFOdUOceImawuaErdVoUAvzU+eXc3vPDeMfQJcEGXjro3N0OAwG9ftXnbcXXTJT2s6p755PjcwatDhcvsW1SdcAQZh9xXJkx0OJReEz5nfyeWmvrV2tRlfyQoHidDE/HxS0OaLBHi51FmKcKBtQFNQYibKtOTzgxFpik2VRFhY18PdkbXzFb6pizNtMU/NlS9jEmRy/T5FX3sZUgvm4DYxNaWLvDWho6ZUKiNzHLKFEd7Y6rEpMeYNyI4fF0z9TUF3K7GGKoKrecw+lvJFDZ6sHuVZIznnf/M6YpXSbPUm5ao6exdNSGpy7PoEyFuV/6cy/Sz4FdGUriOJtURIa5T01zBEW5VSfk5qKp6tS0/FpKA46cS8P4Lg6qcY5rG/+raPpHUEwSuGgJk8085iwDHIDw8chulqU+C6+03w8Gu7FvcVRtOyY04rf6LL7yIUaNGKzt+CkWm3aAX2nq1x9PQ+3Eo/co/QjF7ckcM7KGLqGton62lPRMTU/qjy8Rg1Foeu1k8p6ebZalWlaMZJ4/yVP5Avp5XfAKre5CpY72waXcMflh/Cs/cP5A0x/JTXR1e1V1Ly8jDL3+G4Ub6DIb1uRLBW76NmMnKIyLnrRKBgtxUpMUHIZW2tLhApCWEUL6bbApHp7VsnAOw8shQnEnugMfu6kfOd9eOELHvEGejZgdI9um57rrrcK8BfXeu5ofPfkIffHsUPu2s8AlphVh7ICIItBQEDlJqjqTUfKUFqk9eLF7f7M5nd8PO2gIP39FbOcu3FGyu1nOwPyDPMWambfD169epdBpV3VvIUFXISHmLReDypSJkJJ1GKpOe+EDaByEn46J6k7Oy94a9Ww84uNNSFLS3depMtoK2SEjJw30v7qU6bfHov/vCxsr0muLDpIizVjdXSUnLw4ffHYO9rQk+f3lIpdEdzfXZZNyCgKEQOBudhXtf2EOmMkfcNbN7laY1Q92vJfXDPnyf/HAc6Vm5+OaN4XAlP63q5JqTociTPyF419LqxijXDISAlZ0XRt+xxkC9NZ9ucjPjSkkPa38yEkPJp6IApuZ2RHi6l5IfO9fuMDGrmIxLe9K4pFw8svgARShdwgOzesLDQPmHtP5by/5cVDq+JIdpD1cL/O/FQUSIri2xbC24y3M2TwSOnUrBE68fhLenHebe0rNR1kxsnshUPerM7AJ8vvIkMshExlrnTpREtya55mQodN+HiAn7HR276RLA1TRguV4/BNKTIhEXFYhJD+6qXwfNpFVxUR7SycSl0/qQyYvIT152EjnNGcOGtDz27qT1cSOtD+2t7Co6Htf0mBlZhXjq7SMIOZOGf00NIMfq6h2ca+qvtV3fdzQGqzaG4bq+bnj18T4qtLm1YSDPKwjUFYGwcxl47LUD5Lxvhvtm9YKD7bUz1dd17Fe7fhz5IX7x80mKJmyDjxYORnvy36qNNAkylBC5Df1G3FOb8UqdeiIQd+E4IgI3tTgylJ12voy5KzPlDCWDK4a5lUupqYtNXnYuAWhrbJgJhCN13vk6GL9uOo9BvTwwc5KfyktSz4+mVTTLyMrHL3+EIfB0EuZO74yHZnchk2OreHR5SEHAIAhcTMjBk0uOIDS/cxMAABHcSURBVI72s2/qij7d5EWsPLB7jlzE2s3h6Oprj3ee7k+ksfZaZyFD5dFsoectgQwV5mcoTU8a+fikkq9PWnwwCvMzKfuuGZGdrqTt0fn5OJCvj7n1ldDwxvpI9x5LxOKPT1AoM1SWapmcKiLNodr7KYv3+i0R5AhqjMX/7UvhwhVT4VdsKSWCgCBQHgHOV7Xs2xDdi1hvD0wb34mSodb+B798fy3lPDU9H7+SxjmIXrbumemHB271rzKEvqpnFjJUFTItrLy5kaHLl4uRmRxxxdxFBCiLtEAsVvZeOj+fEnOXrZM/+Thfm9DidDKbffD9KWzYFoWOZNO/+Xo/2lediKyFfa2qfZywsyn47e8IsNqaM+w+cnsXleG32kZyURAQBGpEYOfhBCxZHqiW/LhhdEcMH+DZKp2rOSnpP/suYMvuSPDySS882Ks0k3eNIJarIMkLygEip9cGgfzsxFJtD/v7pCecopXV82Biag07cnL28BtPBKgn+fv0gAk5PjcVsaM1bhY+3Au30Y/90m9CsOyrw+jh74yxQzsYLFlaU3nW2o4jJCIZ/+w9j/DINAzv54b3nx8AbwqfFxEEBAHDIDBygCsG9RyNL1dHYMWGCOw5fBHXj/BB/+7uddaIGGZEV7cXXudu/7FYbKF5hleif3CWP/51o2+D0nOIZujqfobX7G5NSTN0qTifyE6oMnkpcxdpfXKz4lXYuo2jrzJ3Melhs5e1Q8drhll9brz3eCK+pgnqeGgKfEhTNHqwF3p1dW7xOUI4jPVocDy2H4jCxfgsDO3jirkzOqFvgOGydtfn85A2gkBLRyA6LodIUTj+2nkRjvYWGDvMi4iSh1rSpKU9O+cNOnA8Blv3RykSNI2WNGEfRGeHhvuDChlqad+WKp7nWpKhnPRondanxNcnMymclgEogpmlo07bQ5FdnNPHnkLbjSjJYUuQwNNp+Pa3M9h1JI5MQybo190NQ/q4w8ujZZnQOEx+//FYHAuJV0tDjBvaDnfd3Al+3jWHsraEz1meQRBoKghcjM/B12vPECmKpuVQ2mJAT3cM69eOoqlql/W6qTxHZeM4F52BveQczfMML/UybVwH3EnzjCFIkHY/IUMaEi18f7XIEGdtTksIVlmcdVqfYBTkpaGtkQk5OXdR5EeX1LAHLGzatXDUaeFDyiL7+/ZorNsahajYbLg6WZIZzYk2F/h62Tc7lfYliqQLP5+Kk2FJCA5PAidP7NTBliYnL9xAq82z2VBEEBAErh0C7Mf4+7Zo/Lo5Eqw1au9mRdppV/Shdc6aU2608zEZOEHriZ04lYTElBxyO7BV67pNGtG+2kzS9UVeyFB9kWtm7RqFDNGqj5m0fpdm6uKcPlmp59RKxJa27ZS2x0FpfXrClohQ27at+4cyKDwN2w/GYRttF2KywYs5+vk4EJmwQydvB7RztW5y4ea8GGh0fAYiItNx9kIqEaE0sKqak5iNHeyB0YPc4O/TsrRdzey/tgxXEKgSgeOhqdiyj3xr9sUgmZYG4ZcxnnO6dHRQeyuag5qKpGfmKz/D8MgUhJ1NRUp6Hs2JFhg/rB3GD/VAgG/j+ooKGWoq34RGHochyFBBbopauoJJDy9joVu/KwfGppYqj49uCQtyciYCZGohviLVfaTniQztPpKAI8HJOHYqGVk5ReAFHb08bNCO1Nrt3W3gSeTIzdmKtGpXJyEPR2bEJWUjJi5L+f1cjM8kbVYm8vKLYWtjgv4BzhjQwxHD+7upSaq655NrgoAg0LQQOBmWit1HE3DgZDJCiWzwiw5rjTzdbeHpYQ0v2rdzt74qGa5z8gpxkeaZqNgs0l7xPJOB+CT6LaH1CXt0dqD8bc4YOcCNSNvVe9ESMtS0vq+NNpq6kqFLlwpp2YownZMz+fow+cnJiCHNRRtyavbRmbuI9HA2Z3Z65vW7ROqHAOfiOX2eVcIpCKVMs2HnSAsTnYmiostKU+RgZw5ncox0tDdXG6+LZmVhovKLWNFq1uZ0bEyEyYT8BNhXwNi4LWnnQEuOXFJ+PBx5wX3l5hdSKG4hsin6IovS1WfScXJ6rjJ1paTlgnN1sJiatEUnLxtKXGantj7kBF3b1bbrh4C0EgQEgauJQDa9fB0JSSFzdwpOnc0gcpQGzq7PYmdtCmdHSzg66OYdPreyovmGsl9b097cxEjNMUY0z+jmnDa6eYZepjghbWFRsXqBysrRzTO8zyCtTxKZ1JNpcebE1Bw1D/G9HOxMSeNjj26d7NCri4MKuDA3uzZpUiS0nj8REeRmxirCw4uWMvFJJyJ0qbiQ1u+yp6iu7vAKmKpzciZHZ2NTCZM25FeG+CW6kKmJN014Ujl3MUuZ09gxMpqyzl4k+39QWALSMgtKJy6tfl32fD8b8u1xsDFDe1dL9Auwpz07WlqqEHjvdtbNzpepLs8vdQWB1o6AlaUxaV5c1aZhEZuYizNRrKWh+SY+W+2DTycgNSOfXpyKtGp12qu5hgiUo50ZPEkLNaS3I2mgvEgbZYnO5Gvo4tjwKLA6DaiaykKGqgGnpV4qLsxVJi5txXY2e+XnJNMPIK3f5eyncvl497hVLWdhSYu7ilx9BDhiojP55fBWmbCKOyOrQBGj7NwicGZaDm9nU1cBbW3onwmpnE1Iy2NCb3Cs7WEfJV4U1ZYmJ56kRAQBQUAQ0BDwoKSFvFUmxZcuIy2D5hvacvKKdfMMzTc877Dm2Zg00jzH8MbzDZMtXgqDtUptm4nRQMhQZZ98Cy27RFmdd638F2V2pvW7yPnZgpas4Mgu3753KHOXnWtX8k9pOky9hX4MBnksnmCY2MiK7waBUzoRBASBahAwatsGTvZmaqumWrO+JGSoWX98dR+8i9cQ+A28T5EfM1rMVEQQEAQEAUFAEGjtCAgZakXfgLZtjNB12GOt6InlUQUBQUAQEAQEgZoRaCbWvJofRGoIAoKAICAICAKCgCBQHwSEDNUHtVq22bX3CGbc/hgOHj5Zqxa5uXnY+PcuvLj4A0ycen+t2kglQUAQEAQEAUFAEGgYAmImaxh+1baOjUtEXHwSYmITq62nXczOzkVBQRH2HzyBnJxcrVj2goAgIAgIAoKAINCICIhmqBy4h48GYfk3v5Yrrd/pbTMmYe3PH2DalHG16sDZ2QFTbxwD7w4tf82uWgEilQQBQUAQEAQEgauAgJAhPZBZi/PSq/+jNOWUxMVA4ubqXOeejNpemwycdR6oNBAEBAFBQBAQBFoAAs3STHYy6DT+2rwTc2bdhJBTZ7Br31E4ONjitukT0c7DFdt3H8bOXYfU0gS3zZyELn4dSz+qwODT2EHXT4VGwNzcHNOnjsPwof1xMSYej857FalpGep6YlIaZtC1gK6daJmCTKxa/RdOR5yHpaU5Jk8YiSGDeqs+ecmD/YeOY8Of2/HEo3di1dpNOHsuCs899QByc/Kw+Z89cCVCdDNpfDSpagza9cr2CYnJWP/HNoSEnqEVe83g5u6C/z40p7KqUiYICAKCgCAgCAgCdUCg2WmGvv95AxY89xbW/b4VL7z8PrbtOgB7OxusXrsZ/yEy89QL72Djpl20FkoWEaZdePLZt9Qq6owJk6iHHnsZvj6eeP+t5xVxenbhMlyIjoWlhTlm3Hy9gq4bEaDJE0fA3c2Z/H0S8MCj/wdXFyfcNedmJCenUZ9L8Ptf21XdV978BE89/w52EsF67v+WYduOAzhw6CR++uVPcoR+H19/vxbx8Vd8hqobg+qwij+LXvufWhds6ZvPYMqNY0vvX0V1KRYEBAFBQBAQBASBWiLQ7MjQv2dPwU2TdVqWyRNH4o2Xn8T8x+7GhHHXISEpBUMH9sZbr80Hkwa+npqaDtaqsBw4dEIRo06+XjA2McZQ0u6wSSw09CxpluzQmcpZPEjr0rd3gCr78NMV6NcnQPn99Ozuj4cfmK3qfPbVL2q/6IVHVV0+YQK1+sf38dUnr+Kxh+fgkQduV3X0/1Q3Bv16+sdsvjt+Mow0XTrz2eCBvcBjEREEBAFBQBAQBASBhiPQLM1k1la69VM6dbyyblYHLw+FRqdOHUpR0RyRY+OSwb47c++YjhHD+iuzWfCpCKxd/4+qG19ClrSGvDI7C0d07dh1mBa1tFImMi5j8tS+nRvaUp2CAlrI1NQE7PjMMnL4AKW96drFV52bmpiovf6f2o5Bvw1rqNj89zkRMCZ298+9FYtf/I9+FTkWBAQBQUAQEAQEgXoi0CzJUGXPakakpLwYl2hSiot1K+6yNigzKwf3PPwi/Dt7Y8jg3tiz/2j5ZqXnF6Lj1PHcf0/D7FtvLC0vf6CRJ87wXJPUdQxafy88/SAWvvIBftvwD/7+Zy8WPDEXE8cP1y7LXhAQBAQBQUAQEATqiUCzM5PV8zlVs607D2LeM29iKpnZnp1/P3yqCGHXyI25mW7R0pCwsxVum5efr0xwFS7UUFDbMeh3w9qo7gGd8f3yJcocmE0aq8VvfIIjx4L1q8mxICAICAKCgCAgCNQDgRZDhi5fvlzj46/4ab0yc40bM0TVzc3LU3utqUaCCgsLVblnO1dYkGP11u0HyEwWqcr4D9/rldc/UQkVSwtreVDTGCrrhp24P/xkBRzsbcE+Ss9TpBqPYceew5VVlzJBQBAQBAQBQUAQqAMCzZIMcUQXSxZlbNYkPSNLHWrO0nySRiHxLBlZ2WpvXJK/5/OvVoE1NJ98/rMqD6Ew+517jsDZSef7s/fAcYSSNujX3zbhthkTFYF65IlX8OW3q/Hzqj9URJoHESUOu2dJS9Pd5yKRFn3hjNIscfE6B24+rmkMXCcpOYV3KqRfHdCfLdv2IyUlXZ1eN6Sv2vcWJ2oNHtkLAoKAICAICAL1RsBoEUm9WxugYVL0QWSnRcLDW/cDX1OXP//6J1b/thkFpL2JvHARHm4u2L7rIH5ZvVGVnaLIMBdyaD52/BS4bl5+ASIoPxCHxvftG4DDR4MpxD4M6elZyu/m4OFARJy5QCaz9hh5XX/spPXEwqn+fgqPv+fOmRgzcjCRk1SVz+jo8RAcoPpDydfoUYoUMzJqi5co5H0fkadLly7jZGAY5S4yIwdtHzXGL79bjczMbPCyHFlEyAb27wlHJ7sqx+Di5IhVazapXEKMA+dC6uTbQTlpcyoBzlkUFRNHYfU7VG6k22+r2o+pPI5Z6XFISTiDzgPuKX9JzgUBQUAQEAQEgVaNQBsyt9RsX2pEiEL3fYiEyG3oN+Lq/EgXFRWjiByqNX8gTppYVFQEMzNT9ZQMB5MXDq/XzGZ8gRdRjb4YT5FkrpR4URfNVl9YahpD+X55TMU07mLyHUqk9AEczaY/tvL1KzuPu3AcEYGbMOnBXZVdljJBQBAQBAQBQaDVItBioslq+wkaGxuBN01Yu2NkpCNCXMYkg8PYywv7DvlRBJohpKYxlL8Hj4mj0PjD8mzvXv6ynAsCgoAgIAgIAoJAAxBolj5DDXheaSoICAKCgCAgCAgCgkAZBIQMlYFDTgQBQUAQEAQEAUGgtSEgZKi1feLyvIKAICAICAKCgCBQBgEhQ2XgkBNBQBAQBAQBQUAQaG0ICBlqbZ+4PK8gIAgIAoKAICAIlEFAyFAZOOREEBAEBAFBQBAQBFobAkKGWtsnLs8rCAgCgoAgIAgIAmUQEDJUBg45EQQEAUFAEBAEBIHWhoCQodb2icvzCgKCgCAgCAgCgkAZBIQMlYFDTgQBQUAQEAQEAUGgtSEgZKi1feLyvIKAICAICAKCgCBQBgEhQ2XgkBNBQBAQBAQBQUAQaG0ICBlqbZ+4PK8gIAgIAoKAICAIlEFAyFAZOOREEBAEBAFBQBAQBFobAkKGWtsnLs8rCAgCgoAgIAgIAmUQMC5zdg1O2rRpi8yUKOxY9/I1uHvruqWxiVXremB5WkFAEBAEBAFBoBYItLlMUot6jValIDcFKTHHcZn+iTQuAuZWLnBw79W4N5HeBQFBQBAQBASBZobANSdDzQwvGa4gIAgIAoKAICAItDAExGeohX2g8jiCgCAgCAgCgoAgUDcEhAzVDS+pLQgIAoKAICAICAItDAEhQy3sA5XHEQQEAUFAEBAEBIG6IfD/RJoFtiZ37ScAAAAASUVORK5CYII=\n",
"text/plain": [
"<IPython.core.display.Image object>"
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"model.visualize(show_inputs=True)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Setup and run the model"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"# create a setup for two different Object3D to advect\n",
"\n",
"in_ds = xs.create_setup(\n",
" model=model,\n",
" clocks={'clock': [1, 2, 3, 4]},\n",
" input_vars={\n",
" 'objects': {'initial_velocities': ('objects_3d', [2.2, 0.8]),\n",
" 'materials': ('objects_3d', ['metal', 'stone'])}\n",
" },\n",
" output_vars={'clock': ('meshes', 'all_vertices')}\n",
")"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"<xarray.Dataset>\n",
"Dimensions: (clock: 4, objects_3d: 2)\n",
"Coordinates:\n",
" * clock (clock) int64 1 2 3 4\n",
"Dimensions without coordinates: objects_3d\n",
"Data variables:\n",
" objects__initial_velocities (objects_3d) float64 2.2 0.8\n",
" objects__materials (objects_3d) <U5 'metal' 'stone'"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"in_ds"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"out_ds = in_ds.xsimlab.run(model=model)"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"<xarray.Dataset>\n",
"Dimensions: (clock: 4, objects_3d: 2, space_dims: 3, vertex: 3)\n",
"Coordinates:\n",
" * clock (clock) int64 1 2 3 4\n",
"Dimensions without coordinates: objects_3d, space_dims, vertex\n",
"Data variables:\n",
" objects__initial_velocities (objects_3d) float64 2.2 0.8\n",
" objects__materials (objects_3d) <U5 'metal' 'stone'\n",
" meshes__all_vertices (clock, objects_3d, vertex, space_dims) float64 2.603 ..."
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"out_ds"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"<xarray.DataArray 'meshes__all_vertices' (clock: 4, objects_3d: 2, vertex: 3, space_dims: 3)>\n",
"array([[[[ 2.602522, 3.160774, 2.90911 ],\n",
" [ 2.597264, 2.202417, 3.093257],\n",
" [ 2.526336, 2.366329, 3.088938]],\n",
"\n",
" [[ 1.24707 , 1.08344 , 1.351336],\n",
" [ 1.383656, 1.041781, 0.992727],\n",
" [ 0.894755, 1.335864, 1.54705 ]]],\n",
"\n",
"\n",
" [[[ 4.802522, 5.360774, 5.10911 ],\n",
" [ 4.797264, 4.402417, 5.293257],\n",
" [ 4.726336, 4.566329, 5.288938]],\n",
"\n",
" [[ 2.04707 , 1.88344 , 2.151336],\n",
" [ 2.183656, 1.841781, 1.792727],\n",
" [ 1.694755, 2.135864, 2.34705 ]]],\n",
"\n",
"\n",
" [[[ 7.002522, 7.560774, 7.30911 ],\n",
" [ 6.997264, 6.602417, 7.493257],\n",
" [ 6.926336, 6.766329, 7.488938]],\n",
"\n",
" [[ 2.84707 , 2.68344 , 2.951336],\n",
" [ 2.983656, 2.641781, 2.592727],\n",
" [ 2.494755, 2.935864, 3.14705 ]]],\n",
"\n",
"\n",
" [[[ 7.002522, 7.560774, 7.30911 ],\n",
" [ 6.997264, 6.602417, 7.493257],\n",
" [ 6.926336, 6.766329, 7.488938]],\n",
"\n",
" [[ 2.84707 , 2.68344 , 2.951336],\n",
" [ 2.983656, 2.641781, 2.592727],\n",
" [ 2.494755, 2.935864, 3.14705 ]]]])\n",
"Coordinates:\n",
" * clock (clock) int64 1 2 3 4\n",
"Dimensions without coordinates: objects_3d, vertex, space_dims"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# show vertices (x,y,z) coordinates at each step of advection\n",
"\n",
"out_ds.meshes__all_vertices"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Python [conda env:fastscape_py36]",
"language": "python",
"name": "conda-env-fastscape_py36-py"
},
"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.6.2"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment