{
"cells": [
{
"cell_type": "markdown",
"metadata": {
"id": "o4JZ84moBKMr"
},
"source": [
"# Navier-Stokes Forward Simulation\n",
"\n",
"Now let's target a somewhat more complex example: a fluid simulation based on the Navier-Stokes equations. This is still very simple with ΦFlow (phiflow), as differentiable operators for all steps exist there. The Navier-Stokes equations (in their incompressible form) introduce an additional pressure field $p$, and a constraint for conservation of mass, as introduced in equation {eq}`model-boussinesq2d`. We're also moving a marker field, denoted by $d$ here, with the flow. It indicates regions of higher temperature, and exerts a force via a buouyancy factor $\\xi$:\n",
"\n",
"$$\\begin{aligned}\n",
" \\frac{\\partial \\mathbf{u}}{\\partial{t}} + \\mathbf{u} \\cdot \\nabla \\mathbf{u} &= - \\frac{1}{\\rho} \\nabla p + \\nu \\nabla\\cdot \\nabla \\mathbf{u} + (0,1)^T \\xi d\n",
" \\quad \\text{s.t.} \\quad \\nabla \\cdot \\mathbf{u} = 0,\n",
" \\\\\n",
" \\frac{\\partial d}{\\partial{t}} + \\mathbf{u} \\cdot \\nabla d &= 0 \n",
"\\end{aligned}$$\n",
"\n",
"\n",
"Here, we're aiming for an incompressible flow (i.e., $\\rho = \\text{const}$), and use a simple buoyancy model (the Boussinesq approximation) via the term $(0,1)^T \\xi d$. This approximates changes in density for incompressible solvers, without explicitly calculating $\\rho$. We assume a gravity force that acts along the y direction via the vector $(0,1)^T$. \n",
"We'll solve this PDE on a closed domain with Dirichlet boundary conditions $\\mathbf{u}=0$ for the velocity, and Neumann boundaries $\\frac{\\partial p}{\\partial x}=0$ for pressure, on a domain $\\Omega$ with a physical size of $100 \\times 80$ units. \n",
"[[run in colab]](https://colab.research.google.com/github/tum-pbs/pbdl-book/blob/main/overview-ns-forw.ipynb)\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Implementation\n",
"\n",
"As in the previous section, the first command with a \"!\" prefix installs the [phiflow python package from GitHub](https://github.com/tum-pbs/PhiFlow) via `pip` in your python environment. (Skip or modify this command if necessary.)"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "da1uZcDXdVcF",
"outputId": "1082dc87-796c-4b57-e72e-5790fc1444c9"
},
"outputs": [],
"source": [
"!pip install --upgrade --quiet phiflow==2.2\n",
"#!pip install --upgrade --quiet git+https://github.com/tum-pbs/PhiFlow@develop\n",
"\n",
"from phi.flow import * # The Dash GUI is not supported on Google colab, ignore the warning\n",
"import pylab"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "BVV1IKVqDfLl"
},
"source": [
"## Setting up the simulation\n",
"\n",
"The following code sets up a few constants, which are denoted by upper case names. We'll use $40 \\times 32$ cells to discretize our domain, introduce a slight viscosity via $\\nu$, and define the time step to be $\\Delta t=1.5$. \n",
"\n",
"We're creating a first `CenteredGrid` here, which is initialized by a `Sphere` geometry object. This will represent the inflow region `INFLOW` where hot smoke is generated."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"id": "WrA3IXDxv31P"
},
"outputs": [],
"source": [
"DT = 1.5\n",
"NU = 0.01\n",
"\n",
"INFLOW = CenteredGrid(Sphere(center=tensor([30,15], channel(vector='x,y')), radius=10), extrapolation.BOUNDARY, x=32, y=40, bounds=Box(x=(0,80),y=(0,100))) * 0.2\n"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "ExA0Pi2sFVka"
},
"source": [
"The inflow will be used to inject smoke into a second centered grid `smoke` that represents the marker field $d$ from above. Note that we've defined a `Box` of size $100 \\times 80$ above. This is the physical scale in terms of spatial units in our simulation, i.e., a velocity of magnitude $1$ will move the smoke density by 1 unit per 1 time unit, which may be larger or smaller than a cell in the discretized grid, depending on the settings for `x,y`. You could parametrize your simulation grid to directly resemble real-world units, or keep appropriate conversion factors in mind. \n",
"\n",
"The inflow sphere above is already using the \"world\" coordinates: it is located at $x=30$ along the first axis, and $y=15$ (within the $100 \\times 80$ domain box).\n",
"\n",
"Next, we create grids for the quantities we want to simulate. For this example, we require a velocity field and a smoke density field."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"smoke = CenteredGrid(0, extrapolation.BOUNDARY, x=32, y=40, bounds=Box(x=(0,80),y=(0,100))) # sampled at cell centers\n",
"velocity = StaggeredGrid(0, extrapolation.ZERO, x=32, y=40, bounds=Box(x=(0,80),y=(0,100))) # sampled in staggered form at face centers "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We sample the smoke field at the cell centers and the velocity in [staggered form](https://tum-pbs.github.io/PhiFlow/Staggered_Grids.html). The staggered grid internally contains 2 centered grids with different dimensions, and can be converted into centered grids (or simply numpy arrays) via the `unstack` function, as explained in the link above.\n",
"\n",
"Next we define the update step of the simulation, which calls the necessary functions to advance the state of our fluid system by `dt`. The next cell computes one such step, and plots the marker density after one simulation frame."
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/",
"height": 282
},
"id": "WmGZdOwswOva",
"outputId": "3ae4d68d-b586-4bbe-eca9-a223d7720949"
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Max. velocity and mean marker density: [0.15530995, 0.008125]\n"
]
},
{
"data": {
"text/plain": [
""
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAANAAAAD4CAYAAACdW2gvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAALCElEQVR4nO3df+hdd33H8edrMVWZQu3MQrDd7GZR8scaoZTK+ofr7Mj6TyvIsIMRWCEO7FCQQXGw6TZhwmb+GoOIXfOHqzrb0jK6HzELusFITWusaeOW1rWsIU3otNj+05n0vT/u+ZbY5fZ7v/d97zf3fr/PB1zuueece877kLw453y+h/dNVSFpOj9zqQuQlpkBkhoMkNRggKQGAyQ1vGk9d5bEIT8toxeqatvFFqxrgEa2rP8upZbzz45b4iWc1GCApAYDJDUYIKnBAEkNBkhqMEBSgwGSGgyQ1GCApAYDJDUYIKnBAEkNBkhqMEBSgwGSGgyQ1GCApAYDJDUYIKlh1QAleUuSR5J8N8kTST47zL8nyX8lOTa8ds29WmnBTNKV5xXgpqp6OclW4N+S/MOw7A+q6uvzK09abKsGqEY/3/Dy8HHr8LK/m8SE90BJtiQ5BpwFDlbVkWHR55I8nmRfkjeP+e7eJEeTHJ1NydLiyFp+HyjJ5cADwO8D/wM8D1wG7Aeerqo/WeX7ZWNFLZ/zj1bVdRdbsqZRuKp6ETgM7K6q0zXyCvA3wPXtOqUlM8ko3LbhzEOStwI3A99PsmOYF+A24Pj8ypQW0ySjcDuAA0m2MArc16rq75P8S5JtQIBjwO/Nr0xpMa3pHqi9M++BtJRmdA8k6acZIKnBAEkNBkhqMEBSgwGSGgyQ1GCApAYDJDUYIKnBAEkNBkhqMEBSgwGSGgyQ1GCApAYDJDUYIKmh09r36iRHkjyV5KtJLpt/udJimeQMtNLa91pgF7A7yQ3A54F9VfUe4EfAHXOrUlpQqwZo6P12sda+NwErfbEPMGptJW0qU7X2BZ4GXqyqc8MqzwHvGvNdW/tqw5ooQFV1vqp2AVcy6kD6vkl3UFX7q+q6cW2BpGU2bWvfDwCXJ1lpzHglcGq2pUmLb9rWvicYBekjw2p7gAfnVKO0sDqtfZ8EvpLkz4DvAF+aY53SQrK1r7QqW/tKc2GApAYDJDUYIKnBAEkNBkhqMEBSgwGSGgyQ1GCApAYDJDUYIKnBAEkNBkhqMEBSgwGSGgyQ1GCApIZJmopcleRwkieH1r6fGOZ/JsmpJMeG1y3zL1daLJM0FTkHfKqqHkvyduDRJAeHZfuq6i/mV5602FYNUFWdBk4P0y8lOcGYLqTSZrOme6Ak7wbeDxwZZt2Z5PEkdyd5x5jv2NpXG9bEba2SvA34JvC5qro/yXbgBUaN5v8U2FFVv7vKNmxrpSXUbGuVZCtwH/DlqrofoKrODD2zXwW+yKhntrSpTDIKF0ZdR09U1RcumL/jgtU+DByffXnSYptkFO5Xgd8Bvjf8xAnAp4Hbk+xidAn3DPCxOdQnLTRb+0qrsrWvNBcGSGowQFKDAZIaDJDUYICkBgMkNRggqcEASQ0GSGowQFKDAZIaDJDUYICkBgMkNRggqcEASQ0GSGrotPa9IsnBJCeH94v2hZM2sknOQCutfXcCNwAfT7ITuAs4VFXXAIeGz9KmsmqAqup0VT02TL8ErLT2vRU4MKx2ALhtTjVKC2uStlaveV1r3+1D32yA54HtY76zF9jbqFFaWBMPIgytfe8DPllVP75wWY16Y120P1ZV7a+q68a1BZKW2dStfYEzK91Jh/ez8ylRWlxTt/YFHgL2DNN7gAdnX5602FbtTJrkRuBfge8Brw6zP83oPuhrwC8AzwK/VVU/XGVbdibVEhrfmdTWvtKqbO0rzYUBkhoMkNRggKQGAyQ1GCCpwQBJDQZIajBAUoMBkhoMkNRggKQGAyQ1GCCpwQBJDQZIajBAUoMBkhomaSpyd5KzSY5fMO8zSU4lOTa8bplvmdJimuQMdA+w+yLz91XVruH18GzLkpbDJK19vwW8YbcdabPq3APdmeTx4RLPX2bQpjRtgP4a+GVgF3Aa+MtxKybZm+RokqNT7ktaWFMFqKrOVNX5qnoV+CJw/Rusa29sbVhTBWilJ/bgw8DxcetKG9mqP2+S5F7gg8A7kzwH/DHwwSS7GP0iwzPAx+ZXorS4bO0rrcrWvtJcGCCpwQBJDQZIajBAUoMBkhoMkNRggKQGAyQ1GCCpwQBJDQZIajBAUoMBkhoMkNRggKQGAyQ1GCCpwQBJDdP2xr4iycEkJ4d3GytqU5q2N/ZdwKGqugY4NHyWNp1pe2PfChwYpg8At822LGk5rNoXboztVXV6mH4e2D5uxSR7gb1T7kdaaNMG6DVVVaN+b2OX7wf2w0pfOGnjmHYU7sxKe9/h/ezsSpKWx7RnoIeAPcCfD+8PzqyiJfSTc9+Yy3a3vulDc9muZmeSYex7gX8H3pvkuSR3MArOzUlOAh8aPkubzqpnoKq6fcyiX59xLdLS8UkEqcEASQ3tYezNZl4DBmvZl4MLi8MzkNRggKQGAyQ1GCCpwUGEMdZzsGCtLlabAwuXhmcgqcEASQ0GSGowQFKDAZIaDJDUYICkBgMkNRggqcEASQ2tR3mSPAO8BJwHzlXVdbMoSloWs3gW7teq6oUZbEdaOl7CSQ3dABXwz0keHVr4/j9J9iY5muRoc1/Swulewt1YVaeS/DxwMMn3h2b0r7G1rzay1hmoqk4N72eBB4DrZ1GUtCymDlCSn03y9pVp4DeA42/8LWlj6VzCbQceSLKynb+tqn+cSVXSkpg6QFX1A+DaGdYiLR2HsaUGAyQ12JVnjHFdbhahW48deBaHZyCpwQBJDQZIajBAUoMBkhochVuji42A+Svdm5dnIKnBAEkNBkhqMEBSg4MIM+DN/ublGUhqMEBSgwGSGgyQ1NAKUJLdSf4jyVNJ7ppVUdKy6HTl2QL8FfCbwE7g9iQ7Z1WYtAw6Z6Drgaeq6gdV9b/AV4BbZ1OWtBw6AXoX8N8XfH5umPdTbO2rjWzuf0i1ta82ss4Z6BRw1QWfrxzmSZtG5wz0beCaJFczCs5Hgd9e5TsvwPlnh+l3jj5vOB7X8lnt2H5x3IJOZ9JzSe4E/gnYAtxdVU+s8p1tK9NJjm7EX7TzuJZP59ha90BV9TDwcGcb0jLzSQSp4VIGaP8l3Pc8eVzLZ+pjS5Ujy9K0vISTGgyQ1LDuAdpIT3AnuTvJ2STHL5h3RZKDSU4O7++4lDVOI8lVSQ4neTLJE0k+Mcxf6mNL8pYkjyT57nBcnx3mX53kyPB/8qtJLpt0m+saoA34BPc9wO7XzbsLOFRV1wCHhs/L5hzwqaraCdwAfHz4d1r2Y3sFuKmqrgV2AbuT3AB8HthXVe8BfgTcMekG1/sMtKGe4K6qbwE/fN3sW4EDw/QB4Lb1rGkWqup0VT02TL8EnGD0oPBSH1uNvDx83Dq8CrgJ+Powf03Htd4BmugJ7iW3vapOD9PPM/ox5qWV5N3A+4EjbIBjS7IlyTHgLHAQeBp4sarODaus6f+kgwhzVKO/ESzt3wmSvA24D/hkVf34wmXLemxVdb6qdjF6+Pl64H2d7a13gDbDE9xnkuwAGN7PXuJ6ppJkK6PwfLmq7h9mb4hjA6iqF4HDwAeAy5OsPNa2pv+T6x2g157gHkY6Pgo8tM41zNtDwJ5heg/w4CWsZSpJAnwJOFFVX7hg0VIfW5JtSS4fpt8K3Mzo/u4w8JFhtbUdV1Wt6wu4BfhPRteef7je+5/xsdwLnAZ+wuja+Q7g5xiNUJ0EvgFccanrnOK4bmR0efY4cGx43bLsxwb8CvCd4biOA380zP8l4BHgKeDvgDdPuk0f5ZEaHESQGgyQ1GCApAYDJDUYIKnBAEkNBkhq+D8TRzqZ7oeclQAAAABJRU5ErkJggg==",
"text/plain": [
""
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"def step(velocity, smoke, pressure, dt=1.0, buoyancy_factor=1.0):\n",
" smoke = advect.semi_lagrangian(smoke, velocity, dt) + INFLOW\n",
" buoyancy_force = (smoke * (0, buoyancy_factor)).at(velocity) # resamples smoke to velocity sample points\n",
" velocity = advect.semi_lagrangian(velocity, velocity, dt) + dt * buoyancy_force\n",
" velocity = diffuse.explicit(velocity, NU, dt)\n",
" velocity, pressure = fluid.make_incompressible(velocity)\n",
" return velocity, smoke, pressure\n",
"\n",
"velocity, smoke, pressure = step(velocity, smoke, None, dt=DT)\n",
"\n",
"print(\"Max. velocity and mean marker density: \" + format( [ math.max(velocity.values) , math.mean(smoke.values) ] ))\n",
"\n",
"pylab.imshow(np.asarray(smoke.values.numpy('y,x')), origin='lower', cmap='magma')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"A lot has happened in this `step()` call: we've advected the smoke field, added an upwards force via a Boussinesq model, advected the velocity field, and finally made it divergence free via a pressure solve.\n",
"\n",
"The Boussinesq model uses a multiplication by a tuple `(0, buoyancy_factor)` to turn the smoke field into a staggered, 2 component force field, sampled at the locations of the velocity components via the `at()` function. This function makes sure the individual force components are correctly interpolated for the velocity components of the staggered velocity. Note that this also directly ensure the boundary conditions of the original grid are kept. It internally also does `StaggeredGrid(..., extrapolation.ZERO,...)` for the resulting force grid. \n",
"\n",
"The pressure projection step in `make_incompressible` is typically the computationally most expensive step in the sequence above. It solves a Poisson equation for the boundary conditions of the domain, and updates the velocity field with the gradient of the computed pressure.\n",
"\n",
"Just for testing, we've also printed the mean value of the velocities, and the max density after the update. As you can see in the resulting image, we have a first round region of smoke, with a slight upwards motion (which does not show here yet). \n",
"\n",
"## Datatypes and dimensions\n",
"\n",
"The variables we created for the fields of the simulation here are instances of the class `Grid`.\n",
"Like tensors, grids also have the `shape` attribute which lists all batch, spatial and channel dimensions.\n",
"[Shapes in phiflow](https://tum-pbs.github.io/PhiFlow/Math.html#shapes) store not only the sizes of the dimensions but also their names and types."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Smoke: (xˢ=32, yˢ=40)\n",
"Velocity: (xˢ=32, yˢ=40, vectorᵛ=2)\n",
"Inflow: (xˢ=32, yˢ=40), spatial only: (xˢ=32, yˢ=40)\n"
]
}
],
"source": [
"print(f\"Smoke: {smoke.shape}\")\n",
"print(f\"Velocity: {velocity.shape}\")\n",
"print(f\"Inflow: {INFLOW.shape}, spatial only: {INFLOW.shape.spatial}\")\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Note that the phiflow output here indicates the type of a dimension, e.g., $^S$ for a spatial, and $^V$ for a vector dimension. Later on for learning, we'll also introduce batch dimensions.\n",
"\n",
"The actual content of a shape object can be obtained via `.sizes`, or alternatively we can query the size of a specific dimension `dim` via `.get_size('dim')`. Here are two examples:"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Shape content: (32, 40, 2)\n",
"Vector dimension: 2\n"
]
}
],
"source": [
"print(f\"Shape content: {velocity.shape.sizes}\")\n",
"print(f\"Vector dimension: {velocity.shape.get_size('vector')}\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The grid values can be accessed using the `values` property. This is an important difference to a phiflow tensor object, which does not have `values`, as illustrated in the code example below."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Statistics of the different simulation grids:\n",
"(xˢ=32, yˢ=40) float32 0.0 < ... < 0.20000000298023224\n",
"(xˢ=(31, 32), yˢ=(40, 39), vectorᵛ=2) float32 -0.12352858483791351 < ... < 0.15530994534492493\n",
"Reordered test tensor shape: (3, 5, 2)\n"
]
}
],
"source": [
"print(\"Statistics of the different simulation grids:\")\n",
"print(smoke.values)\n",
"print(velocity.values)\n",
"\n",
"# in contrast to a simple tensor:\n",
"test_tensor = math.tensor(numpy.zeros([3, 5, 2]), spatial('x,y'), channel(vector=\"x,y\"))\n",
"print(\"Reordered test tensor shape: \" + format(test_tensor.numpy('vector,y,x').shape) ) # reorder to vector,y,x \n",
"#print(test_tensor.values.numpy('y,x')) # error! tensors don't return their content via \".values\""
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Grids have many more properties which are documented [here](https://tum-pbs.github.io/PhiFlow/phi/field/#phi.field.Grid).\n",
"Also note that the staggered grid has a [non-uniform shape](https://tum-pbs.github.io/PhiFlow/Math.html#non-uniform-tensors) because the number of faces is not equal to the number of cells (in this example the x component has $31 \\times 40$ cells, while y has $32 \\times 39$). The `INFLOW` grid naturally has the same dimensions as the `smoke` grid.\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Time evolution\n",
"\n",
"With this setup, we can easily advance the simulation forward in time a bit more by repeatedly calling the `step` function."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "0hZk5HX3w4Or",
"outputId": "f7811af7-4b58-4ff6-a8b6-6e7bedefaa6e"
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Computed frame 0, max velocity 0.4609803557395935\n",
"Computed frame 1, max velocity 0.8926814794540405\n",
"Computed frame 2, max velocity 1.4052708148956299\n",
"Computed frame 3, max velocity 2.036500930786133\n",
"Computed frame 4, max velocity 2.921010971069336\n",
"Computed frame 5, max velocity 3.828195810317993\n",
"Computed frame 6, max velocity 4.516851425170898\n",
"Computed frame 7, max velocity 4.861286640167236\n",
"Computed frame 8, max velocity 5.125314235687256\n",
"Computed frame 9, max velocity 5.476282119750977\n"
]
}
],
"source": [
"for time_step in range(10):\n",
" velocity, smoke, pressure = step(velocity, smoke, pressure, dt=DT)\n",
" print('Computed frame {}, max velocity {}'.format(time_step , np.asarray(math.max(velocity.values)) ))\n"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "GMKKWQBLHIwP"
},
"source": [
"Now the hot plume is starting to rise:"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/",
"height": 282
},
"id": "Mfl80CjZxZcL",
"outputId": "92f3a9ba-d403-4799-a543-132ee8ed234c"
},
"outputs": [
{
"data": {
"text/plain": [
""
]
},
"execution_count": 9,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAANAAAAD4CAYAAACdW2gvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAAOd0lEQVR4nO3de4xU93nG8e8DWS4lKRcbE+RLcWqnEaoKVh3kyG7l0tJSt5FxFUVxo5a2lkikuEqkKJKVVq2TNFKjpvFfVVUiU1MptZ3aRqAoTkIIqu3UgoBNMJi4GMcXKBdzv3rZy9s/5qy14ZxhZ+edWWZ2n4802pl3zsz8js2zZ+a3v3mPIgIza86kKz0As27mAJklOEBmCQ6QWYIDZJbwnrF8MUme8rNudDQi5lbdMaYBqpk89i9pljLwRr17/BbOLMEBMktwgMwSHCCzBAfILMEBMktwgMwSHCCzBAfILMEBMktwgMwSHCCzBAfILMEBMktwgMwSHCCzBAfILMEBMktwgMwSHCCzhBEDJGmapK2Sfippt6QvFfVHJP1c0o7isrjtozXrMI105ekFlkbEWUk9wHOSni7u+0JEPNG+4Zl1thEDFLXTN5wtbvYUF/d3M6PBz0CSJkvaARwBNkbEluKur0raKekhSVPrPHaVpG2StrVmyGadQ6M5P5CkWcA64K+BY8AhYAqwGtgXEV8e4fHhxorWfQa2R8StVfeMahYuIk4Cm4HlEXEwanqBfweWpMdp1mUamYWbWxx5kDQdWAb8TNL8oiZgBbCrfcM060yNzMLNB9ZKmkwtcN+OiO9I+pGkuYCAHcCn2zdMs840qs9A6RfzZyDrSi36DGRmv8gBMktwgMwSHCCzBAfILMEBMktwgMwSHCCzBAfILMEBMktwgMwSHCCzBAfILMEBMktwgMwSHCCzBAfILMEBMkvItPa9UdIWSa9KelzSlPYP16yzNHIEGmrtuwhYDCyXdBvwNeChiLgJOAHc17ZRmnWoEQNU9H6rau27FBjqi72WWmsrswmlqda+wD7gZET0F5vsB66t81i39rVxq6EARcRARCwGrqPWgfRDjb5ARKyOiFvrtQUy62bNtvb9CDBL0lBjxuuAA60dmlnna7a17x5qQfpYsdlKYH2bxmjWsTKtfV8GHpP0D8CLwMNtHKdZR3JrX7MRubWvWVs4QGYJDpBZggNkluAAmSU4QGYJDpBZggNkluAAmSU4QGYJDpBZggNkluAAmSU4QGYJDpBZggNkluAAmSU4QGYJjTQVuV7SZkkvF619P1vUH5R0QNKO4nJX+4dr1lkaaSrSD3w+Il6Q9D5gu6SNxX0PRcTX2zc8s842YoAi4iBwsLh+RtIe6nQhNZtoRvUZSNIC4BZgS1G6X9JOSWskza7zGLf2tXGr4bZWkt4L/Dfw1Yh4StI84Ci1RvNfAeZHxF+N8Bxua2VdKNnWSlIP8CTwrYh4CiAiDhc9sweBb1LrmW02oTQyCydqXUf3RMQ3htXnD9vsHmBX64dn1tkamYW7Hfgz4KXiFCcAXwTulbSY2lu414FPtWF8Zh3NrX3NRuTWvmZt4QCZJThAZgkOkFmCA2SW4ACZJThAZgkOkFmCA2SW4ACZJThAZgkOkFmCA2SW4ACZJThAZgkOkFmCA2SW4ACZJWRa+86RtFHS3uJnZV84s/GskSPQUGvfhcBtwGckLQQeADZFxM3ApuK22YQyYoAi4mBEvFBcPwMMtfa9G1hbbLYWWNGmMZp1rEbaWr3rkta+84q+2QCHgHl1HrMKWJUYo1nHangSoWjt+yTwuYg4Pfy+qPXGquyPFRGrI+LWem2BzLpZ0619gcND3UmLn0faM0SzztV0a19gA7CyuL4SWN/64Zl1thE7k0q6A3gWeAkYLMpfpPY56NvADcAbwMcj4vgIz+XOpNaF6ncmdWtfsxG5ta9ZWzhAZgkOkFmCA2SW4ACZJThAZgkOkFmCA2SWMKrV2AZVfwiWqv84PLXn6oaftbfvaLkYg+UaEAxUVm3s+QhkluAAmSU4QGYJDpBZgicRAFX8Z/jA7D+q3PbPr15Uqn302hOV215/zcnya02q/rB/6O1fLtV+dOiqym3XHDhQqu08/VipNjh4vvLx1jo+ApklOEBmCQ6QWYIDZJbQSFORNZKOSNo1rPagpAOSdhSXu9o7TLPO1EhTkd8GzgL/ERG/XtQeBM5GxNdH9WId0ROh/PofnH1PqfbSp3sqH/2ev7yzVIv58yu3jem/1PCodKE8Y6Yj1Z3CYt1zpdrSr8ws1Z49/W/Vj493Gh6XQaonQkQ8A1y2247ZRJX5DHS/pJ3FWzyfmcEmpGYD9K/ArwKLgYPAP9fbUNIqSdskbWvytcw6VlMBiojDETEQEYPAN4Ell9nWvbFt3GpqKY+k+cPOzHAPsOty23eSmTM+WKr9eNn0Um3yF1ZUP8Fbb5VKk7Zsr962r6/xgU2dUirF7PLEAAB/UV5m9L03nyzVblx7S+XDD596vvFx2WWNGCBJjwJ3AldL2g/8PXCnpMXUvsX1OvCp9g3RrHONGKCIuLei/HAbxmLWdbwSwSzBATJLcIDMEibcF+pu6CnPps/52w+XN1y3sfLx57eeLNVOHSzP4gGc761eDlRlxvSLpdrM91cvuZn+4VdKtamr7ijVblpb3RXoMJ6FaxUfgcwSHCCzBAfILMEBMksYt5MIVZ12AG6fsaBcfLH8ofztp6s72ry0/7pS7eJg9e8h1R9eycDx8tbTDlW18IXFxw+XanMm7S7VfnPm3MrH/8/pimVDUZ7EsJH5CGSW4ACZJThAZgkOkFmCA2SWMG5n4eoZqGhCdOqH5d7Wz/38hsrHn+kvd/WZ2VM9W1ZlsE4TpKiYszv4TvVSoHOvXVuq3f50uV/27lPVS4wiGh+vXZ6PQGYJDpBZggNkluAAmSU00lRkDfDHwJFhrX3nAI8DC6g1Ffl4RFSfZeoKCfor6+vP/qBU++S+20u1rcfLy10ArplWngU401/9e2haxcm0qiYxAM4NlCcRzvdXLwZ6/u3y/7ZTfeVJj5/0PVH9YpVn+bZmNHIEegRYfkntAWBTRNwMbCpum004zfbGvhtYW1xfC6xo7bDMukOzfweaN6yx4iFgXr0NJa0CVjX5OmYdLf2H1IiI2mlL6t6/GlgNQ6c3MRs/mp2FOyxpPtTa/ALVJ7IxG+eaPQJtAFYC/1j8XN+yEbXZmd7ykpd/erliVitOVj5+0fvK/aqXvb+6B/asKeUvqfUOVJ9g7NjF8rKdDfurZ+Ge6Xu2VFv/2ulS7cz5fZWPt9Zp5BSPjwLPA78mab+k+6gFZ5mkvcDvFbfNJpxme2MD/G6Lx2LWdbwSwSzBATJLmHDfB7rQu79U+37/6lJt+tTqjjbLppff0R6/WO/M4+XlQBcGqn9nHa14jqumVj/r/7394zqvdykv2Wk3H4HMEhwgswQHyCzBATJLmHCTCFX6B06VanOn/lbltoculJfzTZ9cPYnQU/Hr6ULF934ATlws14/1Vk8CLJj9B6Xa6ye+W7mttZePQGYJDpBZggNkluAAmSU4QGYJnoUDoDyzNlhnGczZvsFS7Vhv9SzctMnl30/v1JmFO1Ixu9c3WH4tgKC6bmPPRyCzBAfILMEBMktwgMwSUpMIkl4HzlD74kl/RNzaikF1gnof1M/0lScXztVp7XuyYnnOxTqf/8/3lycRwk3AOl4rZuF+JyKOtuB5zLqO38KZJWQDFMAPJG0vWviWSFolaZukbcnXMus42bdwd0TEAUnXABsl/axoRv8ut/a18Sx1BIqIA8XPI8A6YEkrBmXWLZo+AkmaAUyKiDPF9d8HvtyykXWo1wfLbcBn9M6v3LZnUnmJT1+d03Sfulie3Xsl3hzl6GysZd7CzQPWSRp6nv+MiO+1ZFRmXaLpAEXEa8CiFo7FrOt4GtsswQEyS/D3gep440T1x7kps/+kVDvUO6ty24GY3vDrvTlYXsxxTtUnPq83Nht7PgKZJThAZgkOkFmCA2SW4ACZJXgWbpT2nniqVJsy+5OV2x7rK58hK+p0+7mgc6Wa+113Ph+BzBIcILMEB8gswQEyS1CMYeuX2jdS653R2qxTDWyv13HKRyCzBAfILMEBMktwgMwSUgGStFzSK5JelfRAqwZl1i2aDpCkycC/AH8ILATulbSwVQMz6waZI9AS4NWIeC0iLgKPAXe3Zlhm3SEToGuBt4bd3l/UfoFb+9p41vbV2G7ta+NZ5gh0ALh+2O3riprZhJE5Av0EuFnSjdSC8wngT0d4zFEYeKO4fnXt9rjj/eo+I+3br9S7I9OZtF/S/cD3qS1wWxMRu0d4zNyh65K2jacz2g3xfnWfzL6lPgNFxHcBf23SJiyvRDBLuJIBWn0FX7udvF/dp+l9G9PvA5mNN34LZ5bgAJkljHmAxtMKbklrJB2RtGtYbY6kjZL2Fj9nX8kxNkPS9ZI2S3pZ0m5Jny3qXb1vkqZJ2irpp8V+famo3yhpS/Fv8nFJUxp9zjEN0Dhcwf0IsPyS2gPApoi4GdhU3O42/cDnI2IhcBvwmeL/U7fvWy+wNCIWAYuB5ZJuA74GPBQRNwEngPsafcKxPgKNqxXcEfEMcPyS8t3A2uL6WmDFWI6pFSLiYES8UFw/A+yhtlC4q/ctas4WN3uKSwBLgSeK+qj2a6wD1NAK7i43LyIOFtcPUTsZc9eStAC4BdjCONg3SZMl7QCOABuBfcDJiOgvNhnVv0lPIrRR1P5G0LV/J5D0XuBJ4HMRcXr4fd26bxExEBGLqS1+XgJ8KPN8Yx2gibCC+7Ck+QDFzyNXeDxNkdRDLTzfioihjvrjYt8AIuIksBn4CDBL0tCytlH9mxzrAL27gruY6fgEsGGMx9BuG4CVxfWVwPorOJamSBLwMLAnIr4x7K6u3jdJcyXNKq5PB5ZR+3y3GfhYsdno9isixvQC3AX8L7X3nn8z1q/f4n15FDgI9FF773wfcBW1Gaq9wA+BOVd6nE3s1x3U3p7tBHYUl7u6fd+A3wBeLPZrF/B3Rf0DwFbgVeC/gKmNPqeX8pgleBLBLMEBMktwgMwSHCCzBAfILMEBMktwgMwS/h+wk0eEFYj+fgAAAABJRU5ErkJggg==",
"text/plain": [
""
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"pylab.imshow(smoke.values.numpy('y,x'), origin='lower', cmap='magma')"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "wnbQJvA-HPSL"
},
"source": [
"Let's compute and show a few more steps of the simulation. Because of the inflow being located off-center to the left (with x position 30), the plume will curve towards the right when it hits the top wall of the domain."
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/",
"height": 489
},
"id": "tkhCOzc0ITsj",
"outputId": "f6366c12-1eb5-4ff6-e0d7-94b806bfd8e4"
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Computing time step 0\n",
"Computing time step 1\n",
"Computing time step 2\n",
"Computing time step 10\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA54AAADvCAYAAACaGQseAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAAA1ZUlEQVR4nO3deZwcd3nn8e/Tc2s0ntFlWYcv5FOY2AbhOBgIMTjY8CI4WcJyOXjjrDfXBrI5OAIBErKwCYFkE5bFWTs2iSGEw7HDgokDBh8YY9mW5QtbhyVbx2g0mvvu47d/TGuZmqekac1MdXfVfN6vl15SPfPr6l+Pvl3Tv6mupy2EIAAAAAAAkpKr9QQAAAAAANnGwhMAAAAAkCgWngAAAACARLHwBAAAAAAkioUnAAAAACBRLDwBAAAAAIli4bnIzOxmM/t4recBJIWMI2vINLKMfCPryHh6sPCsI2a2x8xed5yvv8bM9i1g/79rZt1mNmRmN5lZy3z3BcxHkhk3s2vNrGhmIzP+vGa+cwUqkXCmLzCzb5tZr5m5D902s5VmdpuZjZrZXjN7x3zuBziWGuf7e2Y2MeN4/sx87gc4noQz/m4ze7j8unufmf25mTXO+PqSO4az8FwizOz1kt4v6bWSTpf0Ikkfq+mkgMX3QAhh+Yw/36v1hIAFyEv6Z0nXHePrn5U0JWmtpHdK+pyZvbhKcwMWaq58S9Jvzzien1uleQGLZZmk90paLemnNf0a/PdnfH3JHcNZeC6QmV1sZo+Y2bCZfVlS63HGbjKz75rZkfJv+G41s67y1/5B0mmS/rX8m70/nHXbdknfkrR+xm//1p/AVN8t6cYQwpMhhH5Jfyrp2hN6sFiSUpRxoCJpyXQI4ZkQwo2SnoyZV7uk/yDpwyGEkRDCfZLukHRNpftHNmUh38DxpCjjnwsh3BtCmAoh7Jd0q6TLZux7yR3DWXgugJk1S/oXSf8gaaWkr2g6RMe8iaRPSFov6XxJp0r6qCSFEK6R9LykN5V/s/fnM28YQhiVdJWkAzN++3fAzN5hZgPH+XNaeRcvlvTYjF0+Jmmtma1a0DcBmZayjEvSxeUfLM+a2YdnvqUFkFKZ6WM5R1IhhPDsjNpjmj7WY4nKUL6P+kT5mH6/cekElPqMv1o/+UXLkjyG86JsYS6V1CTpr0IIQdJXzey/HWtwCGGnpJ3lzcNm9mlJH1nIBEIIX5T0xQqGLpc0OGP76L87JB1ZyByQaWnK+D2SLpC0V9MH7i9LKmj6Bw5wVJoyfTzLJQ3Nqg1q+piOpSsr+Zak90l6StNvRXybps9KXRRC2LUI+0Z6pTLjZvarkrZI+rVyaUkew1l4Lsx6SfvLwT9q77EGm9laSX8t6VWaDlZOUn+iM/yJEUknzdg++u/hKt0/0ik1GQ8h7J6x+biZ/YmkPxALT0SlJtNzmH1MV3mbY/rSlpV8K4Tw4IzNW8zs7ZLeIOlvajQl1IfUZdzMrtb0a5HXhRB6y+UleQznrbYLc1DSBjOzGbXjvYXkv0sKkl4SQjhJ0rs0/RaAo1xXt1niur6906JdPGf/OTqfJyVdOOOmF0o6FELgbCeOJ00Zj9uXHeNrWLrSnOmZnpXUaGZnz6hdKK6XW+qyku9j3RfHdKQq42Z2paS/0/TbeR+fsZsleQxn4bkwD2j6rXy/Y2ZNZvZLki45zvgOTf+GY9DMNmj6bMxMhzTdbfZYDklaZWadRwshhFtndfGc/ef58tAvSLrOzDbb9EXVH5J08wk8VixNqcm4mV1V/s2mzOw8SR+WdPsJPl5kX5oybWbWKqm5vN1q5Y/BKl979HVJf2Jm7WZ2maQ3a/q6Jyxdmci3mXWZ2evLtUYze6emr4+780S+GcikNGX8ck03FPoPIYQfzdzpUj2Gs/BcgBDClKRf0nR32D5J/1HTITqWj0l6qabfw/1/Y8Z+QtKHbPrC5N+ffeMQwo8lfUnS7vKYE+msdaekP5d0t6YvpN6rBb7HHdmXpoxruk35djMblfTN8n3/9xO4PZaAlGX6dEnj+slvwMclzfwsw9+U1Capp3wfvxFCyPRvy3F8Gcp3k6SPSzosqVfSf5V09axGLFiCUpbxD0vqlPTNGWdDvzXj60vuGG7Rt0gDAAAAALC4OOMJAAAAAEgUC08AAAAAQKJYeAIAAAAAEsXCEwAAAACQqMZq3pmZ0ckI1dAbQlhTizsm46iSmmScfKNKOIYj68g4si4241VdeE5rqP5dYokp7q3t/ZNxJK2WGSffSBrHcGQdGUfWxWect9oCAAAAABLFwhMAAAAAkCgWngAAAACARLHwBAAAAAAkioUnAAAAACBRLDwBAAAAAIli4QkAAAAASBQLTwAAAABAolh4AgAAAAASxcITAAAAAJAoFp4AAAAAgESx8AQAAAAAJKqx1hMAAAAAsJSZqzQ1rnK1TR2Xu9rPtZ/rauva/P4aZpWG8n4WzwxOudqT4VlX2zf2I1ebmOyObAcV/B0scZzxBAAAAAAkioUnAAAAACBRLDwBAAAAAIli4QkAAAAASNSczYXMrFXSPZJayuO/GkL4iJndLOlnJQ2Wh14bQtiW0DyBxJBxZBn5RtaRcWRdFjPe2rw+sv1b637NjfmjV+90tY7zG1zNlu/3tYaYc2u5Wd2FSsENCfmiq5X621zth3e91dX++pmWyPZt/X/j91Ua8fNaQirpajsp6fIQwoiZNUm6z8y+Vf7aH4QQvprc9ICqIOPIMvKNrCPjyDoyjkyYc+EZQgiSji7Pm8p//K8IgJQi48gy8o2sI+PIOjKOrKjoGk8zazCzbZJ6JN0VQniw/KU/M7PtZvYZM2s5xm2vN7OtZrZ1caYMLD4yjiwj38g6Mo6sI+PIgooWniGEYgjhIkkbJV1iZhdI+oCk8yS9XNJKSe87xm1vCCFsCSFsWZwpA4uPjCPLyDeyjowj68g4sqCSazz/vxDCgJndLenKEMKnyuVJM/t7Sb+/6LMDqoyMI8vIN7KOjCPr6j/j5iqrOy52tb/adGVk+5ff6BsJNaxo9ruPO2Vm/j5dIyFJap7VmKgxplFRzL4aVra72mUbR13tkp0HI9tf+fZ73Jg/eO77rtY9eL+rZfWd1HOe8TSzNWbWVf53m6QrJP3YzNaVaybpaklPJDdNIDlkHFlGvpF1ZBxZR8aRFZWc8Vwn6RYza9D0QvWfQwjfMLPvmtkaTf9qY5ukX09umkCiyDiyjHwj68g4so6MIxMq6Wq7XZI7Rx5CuDyRGQFVRsaRZeQbWUfGkXVkHFlRUXMhAAAAAADm64SaCwEAAABAnNbmda72x6e90dV++apoMyFr8Y1+SsN5V7Mmf87MWmMm0hhzbm12w6GGmH3FNByKbVTU2uRKTS+JNkN6+7oDbszkTa9xtffuOOhqI+O7/H1mAGc8AQAAAACJYuEJAAAAAEgUC08AAAAAQKJYeAIAAAAAEkVzIQAAAAAnyDfd2dL2Flf7z1fudLVQim4XugtzjpGkxk5fbGj084hpByQ1RZc9FtMgSKUQd8uKuMZEq5e7Mde80TcNuv/v3+5qN0/8uauFMDXvudULzngCAAAAABLFwhMAAAAAkCgWngAAAACARLHwBAAAAAAkiuZCAAAAGWVxL/XMn3fIQuMSVFdz0xpX+8j5zX5gTMOesQPRDI4Pt7gxrW2+4ZDka7m2oh+2Iqa90GQ+shmKvlGRtfn5h9EJP67FNyYK+VnzaPDPs8ZzulztE6983tW+ddclrtY9eJ+rpQ1nPAEAAAAAiWLhCQAAAABIFAtPAAAAAECi5rzG08xaJd0jqaU8/qshhI+Y2ZmS/knSKkkPS7omcIEAUoiMI+vIOLJsKeTbzF93tmr5Ba72y51vdLXXn5J3tTgHJ6IvCR/s9WP+beIeVzs0vM3VisXhmHvw1/mhMvWa8dXt57naxWd0u9p4T4OrvXCwK7I9MuUzfkrHiKut7fTZyi3311ta5zJX06zrMkOf39fUw37+cWIuk1auM7r/3Kq2mEH+2tMV5/vrVl/x3Ve52m36YWQ7xFzvWu8qOeM5KenyEMKFki6SdKWZXSrpf0j6TAjhLEn9kq5LbJZAssg4so6MI8vIN7KOjCMT5lx4hmlHf+XQVP4TJF0u6avl+i2Srk5igkDSyDiyjowjy8g3so6MIysqusbTzBrMbJukHkl3SdolaSCEcPQc7z5JG45x2+vNbKuZbV2E+QKJIOPIuvlmnHwjDTiGI+vIOLKgooVnCKEYQrhI0kZJl0jyb+o+9m1vCCFsCSFsmd8UgeSRcWTdfDNOvpEGHMORdWQcWTBnc6GZQggDZna3pJ+R1GVmjeXftGyUtD+JCQLVRMaRdWQcWZbGfMc1DnpZ57sj2584r9ONecWWfa7WfMoev/9lvrFLXIMTlaLNf35touiGFAf9CbXv3n2pq/3ZM2Oudv/gZ2dVaDY0H/WU8de3vtrVGpp3u1pPd4erbTvSFdl+Ydzn9HVWcrWOIxOu1rzeN9mxQZ/BMB5ttFXY7/d1770+472T/jl6ervf/zkbox25Tto86sY0rGr1tTW+9vI1vmHS7YPRhknF4pAbU+/mPONpZmvMrKv87zZJV0h6WtLdkt5SHvZuSbcnNEcgUWQcWUfGkWXkG1lHxpEVlZzxXCfpFjNr0PRC9Z9DCN8ws6ck/ZOZfVzSo5JuTHCeQJLIOLKOjCPLyDeyjowjE+ZceIYQtku6OKa+W9PvMQdSjYwj68g4sox8I+vIOLKiouZCAAAAAADM1wk1FwIAAMD8WMzLrp876bdc7YuvOxLZXn1prxujnG9IYm2+IYlaY17q5WLOO5SijVys4Bu75Dp9w6Er3nLI1S55xI/74D0fimzffPjzbsxUvsfPC3Uhl1vmaud3+RyV8r62f3i5qz3QGx1379guN+b0ZWe42jmn+edC4zlrXC1M5l2t/zvjke0btp3txtx4+FG/f7W42rna5GrvGjs1sv3zzXvdmOXNU65m7f55+9KucVfraDstsj0w8oQbU+844wkAAAAASBQLTwAAAABAolh4AgAAAAASxcITAAAAAJAomgsBAAAsOnOVC7t+xdW+9gu+Oc/yzXO/PIttJNQR03CoOa65kJ+bSuH424pv2JJr9Ocwurb4Biqfat8Z2Z74xvVuzD8e/rSfVmnM1VB9TQ0nudr61oKrTU02uFrflM/qntGJyPYLU1vdmMH8i1ytJabBldqafW3W/iXp77efGdn+2HP/241ZuewcV+vKrXe17tKAqz05tDayvXlflxtz9rp+V2ts9d+zTSsGXa1YmnS1tOGMJwAAAAAgUSw8AQAAAACJYuEJAAAAAEgUC08AAAAAQKJoLgQAALDIGhtXuNqNF7a7WvuZvmFImCrNfQcFPyYX10ioJaYJUS7mvEMh2rQlDI/7MaWYecU1KmrxzVLazo02gPnDPX1uzPcfebWr7e2/y+9fMQ1mkKjOttNdbU2LbzY1PNpS0f42tEXHrbcL3ZjBvM+Wxa1c8r7JUfGFYVf7X4een3WzITdmky5ytZcs88/lmJ5aWt0Sbcg1mvfPvalef7uGDj9/i3mgrY1dkW3/COsfZzwBAAAAAIli4QkAAAAASBQLTwAAAABAouZceJrZqWZ2t5k9ZWZPmtl7yvWPmtl+M9tW/vOG5KcLLD4yjiwj38g6Mo6sI+PIikqaCxUk/V4I4REz65D0sJkdvdL7MyGETyU3PaAqyDiyjHwj6+oy4xcsv9rVXnyp7yxSmvC3DbN6jeSafTOdXExTnzAV06QkrgtK0e8vjMxqchTTsCUMxkw2rrlQKbiSzRp39qUDbsx/O/hyV3vf6HZXm5g64O8z22qe8abcMldb2eIbY03m/dKiJefzsLE9mssVY2vcmFEfQTW0V5a3u7+13tX2Df3fyHYIvjnSwdxeV3t1y0pXe/FJfnIdTdHnVVuTH1OY8vMvxTzQwyP++91gsxs3xXwv5L8X9WTOhWcI4aCkg+V/D5vZ05I2JD0xoFrIOLKMfCPryDiyjowjK07oGk8zO0PSxZIeLJd+28y2m9lNZuZ7DU/f5noz22pmWxc2VSB5ZBxZRr6RdWQcWUfGkWYVLzzNbLmkr0l6bwhhSNLnJG2SdJGmfwvzl3G3CyHcEELYEkLYsvDpAskh48gy8o2sI+PIOjKOtKvkGk+ZWZOmg35rCOHrkhRCODTj638n6RuJzBCoAjKOLMt+vv11LictO8fV3tb1Vld71xmjrjZc8D8aH+xrjWx/6cgTbsyuwTtdrVQaczUsvtpnvMFV/vDMU1wtFJ5ztZjPsNfEYDSDDS3+msyO9pKr2Sn+WjQ1x7zUi7kWVP3R58LIfYNuSLHgn2uNzf6asuaYaTR0RueRa/Pfs9euO+JqKw+c5WoHlt41njXIePT/p9VOciPirmGcKvr/1zizU9OfO+zGrGn1QWo8q9Pva8Jfq3nH/uWuVigOzDmv54fvd7XvlU5ztXM6ulxtRVN0Hq2NMc+zkn8OFYb9c/mFMX+N51h+9veovq/njFNJV1uTdKOkp0MIn55RXzdj2C9K8j+FgRQg48gy8o2sI+PIOjKOrKjkjOdlkq6R9LiZbSvXPijp7WZ2kaaX23sk/ZcE5gdUAxlHlpFvZB0ZR9aRcWRCJV1t71N8v95vLv50gOoj48gy8o2sI+PIOjKOrDihrrYAAAAAAJyoipoLAQBQG75RxWWdvxHZ/sab+9yYzmt9cyEVfYMW9cc0Z2jytas6o7+n/ciqV7sx1r3Z1f76d/1Jig/t/kdXG53wTWeQHo0NHa522YZDrhbXSKhnn2+Ccmi4PbK9om3CjWk72Tf/aejp93fQ2uRK4aC/7fBD0UZY9z91qhuza3T2B9hLp7X5xi4/dXKvq51y1kh0rp3+3MdZL/XzX/Pwma52QPe4GhaXzTrBuq50uhsTQlxzIf//enDCZ/CJvmhuDuV/7Mac3+Hv007ucrXSLv9c2zocc7yvQKHgnxvbJ/2J5bsOvMPVNm6KPvaGnP9ZEmKaC+VH/Pfs2WH/PRse3+tqacMZTwAAAABAolh4AgAAAAASxcITAAAAAJAoFp4AAAAAgETRXAgAUCd8I6FXdv6Wq91z28rIdjj9Ir+rh/znqJe6fdOIMFVyNcvFfGpB7nB0TIdvsqIzTnal3/n6eld75fVvcrXXPxptltI3/Jjfv2IaIaEuLG/b6GqNjb65yXBvs6s91bvS1Z4diY47bVmbG7O+z+e55WzfEEgFP4/i9m5X+/LDL4psf+H5ATdmrz3saueGC13tnVNrXe31zdFGNCs3TboxIe8z/lMt/nu73fz3MYQpV8MCWPTcVF6+iVQp5vxV36Q/Nt7b4/9f7xz9YmT7pFaf3Q3to35euU5fK/j9d+fm27DN/0yYyB9xtR8UHnW1lw2+LLK9sWPEjYkzcNg/vx/u9d/vLOCMJwAAAAAgUSw8AQAAAACJYuEJAAAAAEgUC08AAAAAQKJoLgQAqAubV7zV1e75ygo/cGA4shkeuc8NiW3D0+x/15rr8E1KVPK3Dvlog5ZS/4S/3ZG9fv/PHnS1l/7NBa729d+I/jh+/UO++ctk3tdQH05vfJmrNTT6JiWH+5a72qMDvhnLw73RRjmbV/icXtkck9Ocz7iN+wYtIztcSbe9EG3288DwTW7M8pgGMLta/Pz3jG1xtUMDHZHtk4Z8c6Fck59Xo/lmX2b+5SvNhRZXCNFjXpP8f04IPuMP9be62h0jX/G3VXT/a22TG7PypDE/sdJJrmRtPg/FMN/mPP55VSr5bA3kn3e1J/pfGtm+ZKVvGtQ+4vf1jT2+Cd338v/qZzbvx1Q/OOMJAAAAAEgUC08AAAAAQKJYeAIAAAAAEjXnwtPMTjWzu83sKTN70szeU66vNLO7zGxH+e+YC3GA+kfGkXVkHFlGvpF1ZBxZUUlzoYKk3wshPGJmHZIeNrO7JF0r6TshhE+a2fslvV/S+5KbKpAYMo6sq8OMN7jK5zef4mqhu8/VCs8ciWxbzjcfyXX4RhjW5O9TzTG1mOZCs+8hjBXcmDDpG20URn1zjIbvbHO1V33hFZHtl2x5kxuzdeD/uNox2igtNTXIdzQRG7TWT2rqsKsdmfCNV54bjmnQEn4Q3VffpW7M0CG/r/buXj/VQd9c6Ic7NrjaXcP/O7Id1zhmfcuFrvYLXee52sY2f9swK6r5sZhzHzmf570xzZFCTLOXjKv5MbxdvsHV8FTR1b7TM+hqE1P9rtbS1BnZ7gidbkzHKv8cUszxPv+8z8hUacTfdr6CP95P5Ydd7bnx6H0+0t/hxvzgiK997tBWV+sbfjxmIv77nTZznvEMIRwMITxS/vewpKclbZD0Zkm3lIfdIunqhOYIJIqMI+vIOLKMfCPryDiy4oQ+TsXMzpB0saQHJa0NIRztE98txfy6b/o210u6fgFzBKqGjCPrTjTj5BtpwjEcWUfGkWYVNxcys+WSvibpvSGEoZlfCyEEHeP9PiGEG0IIW0II/kOdgDpCxpF188k4+UZacAxH1pFxpF1FC08za9J00G8NIXy9XD5kZuvKX18nqSeZKQLJI+PIOjKOLCPfyDoyjiyY8622ZmaSbpT0dAjh0zO+dIekd0v6ZPnv2xOZIZAwMo6sq8eMr+28xNUu+5R/l9jkVx/zN57V/CfX6ptNhJimF9bZ4msrfaOHMOQbApWOjEe3x/3+izFNYsxPTaWxSVdrvvfRyPbHzvHNKX/hYT/XYnHI1Zaaesj3XtvvalOFNlcbK/iXXY0xzVKOjD4b2X6s3TfLGp/Y5CfS7u8z9PqM3NvrGxOVij73sw2UDriayTcXOn2Zz/iy5ny0ENNIaKjPz39X7mlXC/LNXrKsFhm3WQ20BuQb+IwWlrtaT26vqxVLEzG16P91c0zDuZb1MQfQk1e6Uv8On8vFFGJOJE8VfBOlbcV/j2x3H9zsxvQV/fenb/SZmHv1P0+yoJJrPC+TdI2kx81sW7n2QU2H/J/N7DpJeyW9NZEZAskj48g6Mo4sI9/IOjKOTJhz4RlCuE++k/xRr13c6QDVR8aRdWQcWUa+kXVkHFlRcXMhAAAAAADm44Q+TgVAWvlflFrM9RTLl53pamc3/2xF97C78MPI9uCsa5QkKYS8qx2jCR8y7nz9tKvlv7Hd1QoxlzBODkd/dDW2+mth2v3lkLIOf/2YihV+IHchmtPR3f45NTTQ7mq5mOvYuk4edzU9MRDZvOIN/oGXHpqaY5Konuj/64T8h9XnY65/izvarYm5Rtksel6gd8xf59jafJqrlc72x/DcsL92s3fCz8Rfx+afG72jT7naneavNT23Y6OrbVgefUwTI/661Vue3eBqfZP3uFr8yT9+liym2XkYzB1xY4rBH2jz8tf3loK/JneqMBzZnmzxrw+syZ8fK53iewEUCz6XTbllrjZ//rlQKA672uDozui2droxivlehJj9ZzXPnPEEAAAAACSKhScAAAAAIFEsPAEAAAAAiWLhCQAAAABIFM2Fjsk3XjGL+XDbptUV7W0y3xstBN8MYyldXIxkNTVGc/mu1b/uxnz80n2udsov+uYotvkMfwelmFzuil7IP3zXS9yQv7jnbFf7qwO3utrI+C6/f2TK5Sef5Gr5I92udni/b9BycCja0GJjl2/E09zlP+zcOnwziNxpq1ytdMg30Rh+Jnp83vnCGjfm/l7/mFY2+WP9ZfnDrrauZdZjKPmfBy1N/oPTJ6aS/eB0VGZ56HK1fDHmdUTMbXcN+WYjk/m+6O3MnyewmMZVNuobCanR37Z5nqcdCoVBV3t2/Duu9oW9b3G14UL0uXboOf/66W+7v+xvN7bjRKaIRRM9BvWXXnAjDk34n+mt8sfsnPnlRqEYbcg1kOt3YywXk/thfxxfsdE3bGuQb161uPwxOsS8tveW9ut6zngCAAAAABLFwhMAAAAAkCgWngAAAACARLHwBAAAAAAkakk2F7JZD/tFK97oxvzK6gtd7U0b/IXPp5484Pcfc8F/9+Fo04nvdvuGFjft3+9q24f+ydVKpZjmAViyGhu6XO0DZ/xmZPvDT77CjbGpKVcLQ75xhMb9RftxwqkbItvtb+pwYz7StszVXvxy31jjuqd8w6Gxyb0VzQP1x6zZ1ZY1+uPk1LBvxvLkYX+sfHSgJbJ9wUSLG3PFBt9cKHfeej+5Lt8QKDfgb3vg29H5fn6Hz/e9E0+72sUN57ja+jZ/n2uL0eZC/37f6W5MseSfs6gPB0pPudrTA69zNYtpLPJ48Tm/wzDruBjzHOob8sfTdUf86xSNTrhSW6NvAGOzWh/Ft0DxDVXGJw+62oP6qqvtOBBtRDMwsceNGYvZF+rDkRF/fLv70OWutqrkm0ZZTHOhEPKR7f7gmxeVJje6WsOAf53Sclarqy2zLldL3tJuHFQJzngCAAAAABLFwhMAAAAAkCgWngAAAACARM258DSzm8ysx8yemFH7qJntN7Nt5T9vSHaaQHLIOLKOjCPLyDeyjowjKyppLnSzpL+V9IVZ9c+EED616DNadL5Zxdkrro5sP/7rTW5M43/a4Gph3RZfi2mWEqdrPNoQ6PyeHjfmN2/zF1Zf/qf/2dXuHfq8n0fwzQNQsZuVmoybq1y6/Fdc7SNfWxm91QM/cmMK33vG1Yaf8fsfH/XPj7gGWm3t0UYBy8/1+2p8xSZXe+u/+OYr37rkGlf7h0OfjGwH+aZEOKabVcOMh+Abkjw94DM0MeZ/JO0e9U1Vvt8zEtneO+KPwy87pd3V2g/HNM8a8cfO0sFhV9t2JNrk4tYjN7kxU4UBV9uy6lxX62r2TYKaV0afL/920D/ufKHX1SCpDo7hfaP+eHrXwStd7YIV/rZHSntcLcxuUhJKbkxjg6+pz2c89PtmWU0531xo/iprOHRgsjuyzTH8hNysGme8UPTZ+vrgF11tReuLXC0X01xotp6RJ1zNGk71A/v98VlN/rX+GcE3aPuxew1FM6Bqm/OMZwjhHkl9VZgLUBNkHFlHxpFl5BtZR8aRFQu5xvO3zWx7+fR/zO/wppnZ9Wa21cy2LuC+gFog48i6OTNOvpFiHMORdWQcqTLfhefnJG2SdJGkg5L+8lgDQwg3hBC2hBD8+1SB+kXGkXUVZZx8I6U4hiPryDhSZ14LzxDCoRBCMYRQkvR3ki5Z3GkBtUXGkXVkHFlGvpF1ZBxpVElzIcfM1oUQjl45/ouS/BXBdaKz3Tcuuf+Ktsh2wx9c7W/4gm/0k3vwYT8un/e1OC3RRhFhRacfc+0bXenO57/mamfecrGrHRp8oLJ5oCL1mnGLaZb1to0+S3a4P7Ld+/ndbszOg2tcrVCq7HdRMS0tnOa9ftRZj/t5rP5V/xy65kzffOUfD896DpVoTLEQ1c24bz7ywMQuV2td5vPX1uCbP2wP90a2fzzuG6Vc23+Rq63f5y+Ryq3xjYkmdvlM7h6N/ricmDrg95Xz+3qo8LSrPXjkxa725S91Rbb/z6Eb3BhUrtrH8GLRN/D5ytDtrvb46KtcbWRW0504IUy62rIWn9Ow/mRXs0E/txW+d5Ust5jHWP+cp43L4qr+6xT/PzgW00QqH/NcKIW5s1QsDPnbjfvXEeGU1b722D5XO6O91dVsIPoaigZX1TfnwtPMviTpNZJWm9k+SR+R9Bozu0jTKdwj6b8kN0UgWWQcWUfGkWXkG1lHxpEVcy48QwhvjynfmMBcgJog48g6Mo4sI9/IOjKOrFhIV1sAAAAAAOY0r2s80+S0Jt/Aa+WHXh4t3HaXGzP2owFXGzzY5sdNNlU0j/a26DVrnaf4Dy1ve7n/AOqW61/pamfd4t/ffkhc47kUNDSe5Gq/+todrjZ8a/RDkh95fqMbM1b014tWqlTBxTq52Z/TLGn0hbWudum/+Os+X/Uqf+PmR7oi2xNTY3NPAnVr7/gPXe0vHrjG1dpjfkqNzLquaHxyvxuze+QVrvaKK05ztXCyP54uG/u+qz3yNX/d8Wylkr8Ob0f/ba72O4N3uprN+oD1YtFf74R65q9pHBh91tUeaexxtUIp7lg2e3/+eB1ijsPhrE2+tnqVq/38uh+52gd3RK/xLIljLI4vBP9aNp+PuW7S4s5zzb5+02e8OOpD3hiT8YYpf73z5lv9taa57uh1+Bxnq48zngAAAACARLHwBAAAAAAkioUnAAAAACBRLDwBAAAAAInKVHMhi3k4l7Wf4Qc+Gm3ic/hb/gL6x/f5ZixTJb9Oj+mfEqvYFx3Z2u0bEVzUd8jVVuaedLWXda5xtR8Mzfrg5zB3IwykT4P5T/3e/eRKVxuajI57dsR/kHJrzl+0f9oy3xxlZYvP0mTRPxcOT7ZEtl8Y8423Jkt+/l1P+w887x73jbyKJTKdJaMTz7naXzz/t662qv1cV5uY6p1z/y+M+0YVdmTAD4zpghWK/kPLJ0r+mO1VMkYqxTaTQdaE4I+nk/nDi7b/nuHlrnZqj29epFGftxdfNexqyx6MNn8bGhuY99ywdAXFNBcKFbxajmlANNbjX0e09vX52zb64/27LvKNCz+4O/r6eWSc5kLVxhlPAAAAAECiWHgCAAAAABLFwhMAAAAAkCgWngAAAACARGWquVCcou+fosF/749s3/fcaW7McMFfqNzZVGHjiJj7DLPaEB2c8BdMj+7e4GqXfWu/qz056BuvhFDZ3JBuk/luV/va86tdbUNbtDnKo0d8KLes9hf7T8Y00Bqc8lmNa7Q1u+HQ4Um//0eOTLjaxjbfIOOOff62+cIRV0O2FAqDrnZo8IcV3NLnuynm16qh2zelsJjmQsUB38jqwhXR4+6dAxVMC0tczIuBChtQzWYxrQxXtfumQdYb03hl0jc5CmO+AczalhdHtofGnnFjgPmJey7M5g/ay1bnXc16Y14LDI64UvuL/H2ub704sv3s+K4K5oXFxBlPAAAAAECiWHgCAAAAABLFwhMAAAAAkKg5F55mdpOZ9ZjZEzNqK83sLjPbUf57RbLTBJJDxpF1ZBxZR8aRZeQbWVFJc6GbJf2tpC/MqL1f0ndCCJ80s/eXt9+3+NM7MUH+YvnbR/7N1d6567LI9o/6mt2Yk1v9RcnDBb9Ob835cXENjUaL0cYAYwXfKOCBw/6/YzDvGx89lP+qv4N5NiyApBRlXDENJm47vM/Vzm44JbL9eGmnG7Np6jxXazTfSKgrpkvLRNHXDk9Faz3j/onw/YmvuNr2517kat3DD7saFuRmpSLj8z2O+edFc8yxOXQP+1u2+syHcT+Pl6+c1XDouROYHqrhZqUi4/NjOf86pb3dNw1ST0xzobxv0BImS6527ZoLItt/1H9bzEwqaRKDBNysDOdbkhpyra5mPvZSb7+vDfpGW6Uxn9WrOjZHtp/t941EeT2drDnPeIYQ7pE0+0j2Zkm3lP99i6SrF3daQPWQcWQdGUfWkXFkGflGVsz341TWhhAOlv/dLWntsQaa2fWSrp/n/QC1QsaRdRVlnHwjxcg4sozXKUidBX+OZwghmNkx33sRQrhB0g2SdLxxQL0i48i642WcfCMLyDiyjNcpSIv5drU9ZGbrJKn8d8/iTQmoC2QcWUfGkXVkHFlGvpE68z3jeYekd0v6ZPnv2xdtRotseHK/q/3FU9GHPRgG3JgLOzpd7YpT/AX6Xc1TrjZZ9BcrH5mKNrC4Y59vhnFP/l5Xu333kKsNj+1yNSy6Os24/0XltoG/d7UdrWdEthtiGlM8N3yOq/3MqglXW7vMX7Q/MuUbsrSOLots33HA325s4gVXG5143tW4uL8q6jTj8+F/h9rZ6Jun2NWvcLXiGb65VeOmx13tqsEHo2MeW+nGFAoxjV1QS5nJeHNjl6sVYxoeFq96XUX7azzTZ/x3+qIZ/+juVW5MvtBb0f5RFZnJtyS1Nvu8Fcf9uOIVl7uajfjXyi2nP+tqHz28PbJ9w5fXuTHjk75hIxZPJR+n8iVJD0g618z2mdl1mg75FWa2Q9LryttAKpFxZB0ZR9aRcWQZ+UZWzHnGM4Tw9mN86bWLPBegJsg4so6MI+vIOLKMfCMr5nuNJwAAAAAAFWHhCQAAAABI1II/TqXexV0k/O3CDZHttpY1bswVbf5dDX1TvmmQ5Ju2jBf9er531m1Xtfg9HTh8f8z+49B4BTP5hkOjE89FtnO55W7MxvU+p/1T/pDQ1uAzPpT3zYWGC9GMn72sw425Z7DgasBCtbX4BhFvOs83rbIDrb52km8kZ71zNwn6X+f/pqtd//jH57wdMB/r21/mamtf61+T2H7/mic0xow7MuDvpBRtevhnZ/mM/+GP/+Q4swTm77TWS1yt/ef8R5OWDh5wNRsd9TvsG/TjZr3s+cBpv+rG/PGOP42ZHZ9As1g44wkAAAAASBQLTwAAAABAolh4AgAAAAASlflrPOMUitH3fa9peZUb0z3u38/d1uCvk2iKWbqPF83V+qeitSOT/jrNM1a83tX29H/T3wFwgs7s9B3X+6f8uL6YazyXN8Zdz+nHDReiGR8vlNyYTSve7Gq7+lP9mdeoA+9Yea2rFQov+IHd/tpN6/DXC4VDA742K845f5jX5Z2/62rfHfyMHwicoOvWXOxqxe6drtb0wn5Xs5jXLmF/r99fPhrqziZ/DH9153tc7Z7Bv3Y14ES9Z+MmVysd8Ncs5/bGHNsn/AuasO+IqxUnoy/a17b6jF/e+V5X4zi+eDjjCQAAAABIFAtPAAAAAECiWHgCAAAAABLFwhMAAAAAkKgl2Vxo9gfBluQb/Yzk/QXHRyb9BfqtDX7tPhHTXKhnVrOifMnvP8jXgKRMFHwDrf68z3NnPq6RkH8u9EzEdFsBqqDBfPbGJ3xTrOIe31yoocnnu9Q97Gojvc2R7RwfKI4qmox5eTAW02PlpF2+WZa1xDwX9g662sRI9LlQIuKoopGCP45P7p5wtdbVB/2N8wVXKuzxx/GJ8ZbozUr+Pqdi1gRYPJzxBAAAAAAkioUnAAAAACBRLDwBAAAAAIla0DWeZrZH0rCkoqRCCGHLYkwKqBdkHFlHxpF1ZBxZRr6RJovRXOjnQgi9i7Cfmolr6jOc9xcXjxb8CeKBqZgLk+OaAMxq5BK4aD9NUp/xOOPFuNz7pkEDMbW4JgDD+WioY/pu0UCrfmUu4z3D7a62ZkePq7XlDrvaxE7f0GL/4TWR7ZJoppUyqc74REy/k4MHOl2t7YmYBlqd/qXe6E5/LN5/eEVkOy7jJZpq1atU51uSRmNeVxzavdzVNnQccbUQ89Ji+Dmf+/190edMPibOvE5JFm+1BQAAAAAkaqELzyDp38zsYTO7Pm6AmV1vZlvNbOsC7wuoBTKOrDtuxsk3MoCMI8t4nYLUWOhbbV8ZQthvZidLusvMfhxCuGfmgBDCDZJukCQz4z0aSBsyjqw7bsbJNzKAjCPLeJ2C1FjQGc8Qwv7y3z2SbpN0yWJMCqgXZBxZR8aRdWQcWUa+kSbzPuNpZu2SciGE4fK/f17SnyzazGpsT8k3oWifXOdqTTnfeCVf8r9MGpyKdgZ4Jjy/gNmhGrKU8Zx8Tp+bGHa15U2+WUVzzv9+Kq6B1tCs4pOTh/wg+rHUlaxkfLTgA3l/70mu1rln0tXWjvnnQX9fh6s9NxJtVrStnzCnQVYy3j/pX1fce2iVqy17asrXlvlaT1+Xqz0z1DFrm4zXu6zkW5IGfUx17/61rvaaXLerNeT8z4CDA/44/vRgtFnRLn/4p4FWwhbyVtu1km4zs6P7+WII4c5FmRVQH8g4so6MI+vIOLKMfCNV5r3wDCHslnThIs4FqCtkHFlHxpF1ZBxZRr6RNnycCgAAAAAgUSw8AQAAAACJWujHqWTC3n7/dvjmFb/kat2TXa5WDG0V3cfzpd7I9qj1VzQPYDHs6P+6qzWveKerHRxb5mrFUpOrxV16v2dWs6LhXJ8bs7v/X48zS2B+bu35M1fL2R+52tqWla62cdRnfqzgfzQ+0Nsc2d4xNObGfHfwM8edJzBff3fgT10tX/qQq61qXuNqnU0FVxuOyfgPeqNN6J4a8p1X7hv8n8edJzBff7nX90S6bspnvKPJZ7wlprnQYN5n/IdHoufbHh8edGMeGPzsceeJheGMJwAAAAAgUSw8AQAAAACJYuEJAAAAAEgUC08AAAAAQKJoLnQMlTZjOZJvcbWgoquN22hke0//NxcwO2Dhnuy/1dVaun7N1Q6N+oyX5C/kH8kNRbbjnkNAtfzDId9wqLXhw67W1rjc1Yol3z5r/9hUZPvbg59awOyAhbu5++Ou1pTzGW/O+QZxeX8I176xycg2zbJQazce9BmXfMOhnJmrFWKO4/vGo03haJZVfZzxBAAAAAAkioUnAAAAACBRLDwBAAAAAImyEOI+Cj6hOzMLUsPcA4EFKT4cQthSi3sm46iO2mScfKM6OIYj68g4si4+45zxBAAAAAAkioUnAAAAACBRLDwBAAAAAIla0MLTzK40s2fMbKeZvX+xJgXUCzKOrCPjyDLyjawj40iTeS88zaxB0mclXSVps6S3m9nmxZoYUGtkHFlHxpFl5BtZR8aRNgs543mJpJ0hhN0hhClJ/yTpzYszLaAukHFkHRlHlpFvZB0ZR6osZOG5QdILM7b3lWsRZna9mW01s60LuC+gFsg4sm7OjJNvpBjHcGQdGUeqNCZ9ByGEGyTdIB397CAgW8g4sox8I+vIOLKOjKNeLOSM535Jp87Y3liuAVlBxpF1ZBxZRr6RdWQcqbKQM54PSTrbzM7UdMjfJukdc9ymVyrulbR6+t+plfb5S+l/DMeb/+mLdB9kPL2yPv9aZfxovqXsf4/rXZbnXw/HcCnb3+M0yPL86yHjWf7+pkHa5y/NI+PzXniGEApm9tuSvi2pQdJNIYQn57jNGkkys60hhC3zve9aS/v8pfQ/hmrMn4wz/1qp1vxPNONH813NOSaF+ddWvR/DqzXHJDH/2qr3jPP9ra20z1+a32NY0DWeIYRvSvrmQvYB1DMyjqwj48gy8o2sI+NIk4Vc4wkAAAAAwJxqtfC8oUb3u1jSPn8p/Y+h3udf7/ObC/OvrTTMPw1zPB7mX1tpmH8a5ng8zL+26n3+9T6/uTD/2jvhx2Ah0FUZAAAAAJAc3moLAAAAAEgUC08AAAAAQKKqvvA0syvN7Bkz22lm76/2/Z8oM7vJzHrM7IkZtZVmdpeZ7Sj/vaKWczweMzvVzO42s6fM7Ekze0+5norHYGatZvYjM3usPP+PletnmtmD5Rx92cyaaz3Xo8h4dZHx6kpbvqV0Zzzt+ZbIeNLSnG8p/RlPW74lMl5tZPwnqrrwNLMGSZ+VdJWkzZLebmabqzmHebhZ0pWzau+X9J0QwtmSvlPerlcFSb8XQtgs6VJJv1X+nqflMUxKujyEcKGkiyRdaWaXSvofkj4TQjhLUr+k62o3xZ8g4zVBxqskpfmW0p3xtOdbIuNJu1npzbeU/oynJt8SGa8RMn5UCKFqfyT9jKRvz9j+gKQPVHMO85z3GZKemLH9jKR15X+vk/RMred4Ao/ldklXpPExSFom6RFJPy2pV1JjuR7JVY3nSMZr/1jIeHLzS2W+y3PNRMbTnO/yXMl4MvPORL7L801txus933FzIeM1eSxLNuPVfqvtBkkvzNjeV66lzdoQwsHyv7slra3lZCplZmdIuljSg0rRYzCzBjPbJqlH0l2SdkkaCCEUykPqKUdkvIbIeOKykm8pRfk4Kq35lsh4DaQqH0elNeMpyrdExmtqqWec5kILFKaX+XX/mTRmtlzS1yS9N4QwNPNr9f4YQgjFEMJFkjZKukTSebWd0dJS7/k4ioxjvuo9H1K68y2R8VpKQz6kdGecfNdWvefjKDJe/YXnfkmnztjeWK6lzSEzWydJ5b97ajyf4zKzJk0H/dYQwtfL5VQ9BkkKIQxIulvTp/O7zKyx/KV6yhEZrwEyXjVZybeUonxkJd8SGa+iVOUjKxlPQb4lMl4TZHxatReeD0k6u9wFqVnS2yTdUeU5LIY7JL27/O93a/q92nXJzEzSjZKeDiF8esaXUvEYzGyNmXWV/92m6ffEP63p0L+lPKye5k/Gq4yMV1VW8i2lJx+pzrdExmskTflIdcZTlm+JjFcdGZ+hBhelvkHSs5p+b/AfVfv+5zHfL0k6KCmv6fcvXydplaa7T+2Q9O+SVtZ6nseZ/ys1fep+u6Rt5T9vSMtjkPRTkh4tz/8JSX9crr9I0o8k7ZT0FUkttZ7rjDmT8erOn4xXd76pynd5zqnNeNrzXX4MZDzZ+aY23+X5pzrjact3eW5kvLrzJ+PlP1a+IQAAAAAAiaC5EAAAAAAgUSw8AQAAAACJYuEJAAAAAEgUC08AAAAAQKJYeAIAAAAAEsXCEwAAAACQKBaeAAAAAIBE/T/UNKy+rr7o1AAAAABJRU5ErkJggg==",
"text/plain": [
""
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"steps = [[ smoke.values, velocity.values.vector[0], velocity.values.vector[1] ]]\n",
"for time_step in range(20):\n",
" if time_step<3 or time_step%10==0: \n",
" print('Computing time step %d' % time_step)\n",
" velocity, smoke, pressure = step(velocity, smoke, pressure, dt=DT)\n",
" if time_step%5==0:\n",
" steps.append( [smoke.values, velocity.values.vector[0], velocity.values.vector[1]] )\n",
"\n",
"fig, axes = pylab.subplots(1, len(steps), figsize=(16, 5))\n",
"for i in range(len(steps)):\n",
" axes[i].imshow(steps[i][0].numpy('y,x'), origin='lower', cmap='magma')\n",
" axes[i].set_title(f\"d at t={i*5}\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We can also take a look at the velocities. The `steps` list above already stores `vector[0]` and `vector[1]` components of the velocities as numpy arrays, which we can show next."
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA6IAAAD2CAYAAADMDzr7AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABcCUlEQVR4nO29e7Bs6Vne97zr0t37dvY+9zlz0YykEZJGgEbJIKQC2UgESoiiwIQkQEyUsohwYlJQUDEKTmIgYEOwUcoFJhZGHpHiGiMHGdvEihARsrGkkTS6DtKM5qKZ0Zz72fe+95c/dg/sXu9z5nynu/fq7nWeX9Wpc/Z71uVb33rWt9bavb/fthAChBBCCCGEEEKIskhm3QAhhBBCCCGEELcWehEVQgghhBBCCFEqehEVQgghhBBCCFEqehEVQgghhBBCCFEqehEVQgghhBBCCFEqehEVQgghhBBCCFEqehEVQgghhBBCCFEqehEVQgghhBBCCFEqehFdAMzsT8zsB1/g/+8xs2Bm2Zjb/34ze8rM9szs/zazE+O3Voib5ygzbmbfZGYDM9s99Oetk7VYiBfmiDN9zszeZ2ZfGW7jnsL/183s3Wa2bWbnzezHxjgEIV6QGWf8QTPrFMb1dIzDEIJyxPn+djP7sJltDsfof2pma4f+/5YZw/UieotjZq8C8E8A/ACAswD2AfzjmTZKiOnzlRDC6qE/75l1g4SYgAGAPwLwn17n/38KwMsA3A3gjQD+tpm9uZymCTEVbpRxAPjfCuN6v6S2CTEp6wB+FsDtAF4J4A4Av3jo/38Kt8gYrhfRCIbf8bj30NcPmtnP3mCdnzCzjzz/nRIz+2/N7HNm1iDLHjezPzSzS2Z2bfjvO4f/93MA3gDgl4ff8ftlsrsPDf/eHC7z+ps4vP8SwL8MIXwohLAL4H8G8N2HvzMjqk/FMy5uQaqc6RDChRDCPwbwsess8lYA/2sI4VoI4REAvwbgv47dvlgMbvGMi4pT8Xz/Vgjhj0II+yGEazgYo7/h0CK3zBiuF9Gj4xcBtAH8T2b2MgB/D8BfDyG0yLIJgH+Gg+98vAhAE8AvA0AI4e8A+FMAPzz8jt8Pk/X/yvDvjeEyf2Zm3zj8yP96f75xuM6rAHzq+Q2FEL4EoAPgqybtAFF5FiXjAHDGzC6Y2RNm9k4zW5lKD4iqsUiZppjZcQDncGhcH/77VTEdICrPwmf8EP+dmV01s4+b2Qt9cipuHRY1338FwOeAW28MH2tOobgxIYSBmf1XAD4B4L/AwY+QfPI6y14B8PvPfz38TswHJ9z/hwFsRCy6CmCrUNsCoE9ExQuyQBn/cwD3D/++G8B7APwSgB+aZP+ieixQpl+I1eHfh8d1jekCQGUyDgD/CMCP4yDb3wrgd83sfAjh301h22JBWcR8m9m34OAT0K8flm6pMVyfiB4hIYQncRDqewD8yvWWM7NlM/sndiAM2sbBx/0bVs7E+10Axwq1YwB2Sti3WHAWIeMhhPMhhM+HEAYhhCcA/G288LwjcQuzCJm+AbvDvw+P6xrTxV9QgYwjhPCJEMKVEEIvhPCvAfwmgO+edbvE7FmkfJvZ6wD8FoDvCSF8cVi+pcZwvYjGsQ9g+dDXt8WsZGbfDuD1AD6A0UnIRX4cwMsBfH0I4Rj+8uN+G/4dbrAr9/9m9gYbtckV/7xhuOjnALz60HovAVAH8MXiNkWlqXLG2bY09lWfWynTf7nRg/lGz+HQuD789+dutK5YOG7JjL/AvuyGS4lFotL5NrPXAHgfgL8RQvjAX2z0FhvD9TAWx8MAvt/MUjuwVv3VG61gZqcA/FMAP4iDj9y/w8zecp3F13Dws+mbdvCrU/5u4f8vAHjJC+zuEg4Mc3+xTAjhTws2ueKfPx0u+pvDtr3BDubN/QyA94YQKvmdF3FdHkZFM25mbzSzu+2AuwD8PIA/uNHxiYXnYVQ008O2NnDwTUMAqNuojOM3cDBH6riZvQLAfwPgwRc+erGAPIxbNONm9j1mtmpmiZl9K4C/joOHelEdHkZF821mX40DK/R/H0L4l2Tbt84YHkLQnxv8AfAADr4TsQPg/wTw2wB+9gbrvBfA/3Ho628D8BUAJ8mytwP4Exx8HP9FHMxdCwCy4f+/fli/BuAfXWd/P4ODi2ITwOtu8vi+H8CXAezh4AH9xKz7XH/K/VPljAP4MQDP4uC7q0/jYG7R2qz7XH+O9k+VMz1cNxT/HPq/OoB3A9jGwcPUj836fOjP9P/c4hn/UxzMm9vGgcjle2d9PvRnun+qnG8cSJIGw30//+dzh/7/lhnDbXjAQgghhBBCCCFEKehHc4UQQgghhBBClIpeRCfAzH7yOpOR/82s2ybENFDGRdVQpkXVUcZFlVG+q4V+NFcIIYQQQgghRKnoE1EhhBBCCCGEEKWSlbmz9Ww5nKmvl7nLm8Lm+DdQldG0WX02Pu0P5R/bP385hHB6uluN49gcZbyMzMzzNTNtpp3TSTb3pRllnOV72hGIzRRfzPcqWy5+H3FnadzrIHb7sYQJzgbLd+z2+Lpx22PrznIMX8uWw+naRun7ZRli+WBnJDG/XEIW5MuRWkL2S5azwnL0Opj2AMGyNiC1wLLma4PIdQdse+Tg2HKDOcv4SroSNvKNkVqftLFHGs6W68N34uA6I0ARI32YkM/IMhKuGgl5LfX7aJA819Keq+V1fxxJvdCWWuqWQUZqCfmcj10LLBx9Espen9T8cqFLzoU/VPR7vn39Aav5RvdJxlntieZzURkv9UX0TH0d73zV3xhrXTZYxhI7DsbuY5JxddzjYDeVacOuh/h1Ix9YSlj3Oz72956K2uARcKa+jl+8721jrRv74wnsgWCS7U2S+zKumVlQRk7JrYY/AJHlvvuhn5tJxs/U1/FLhTF8kpzFZootl5GHi9R8b6Vke6zG9kvXZftNyAMMbvyQPsl9jUEfltlzToh94PDL9SIfTLrkoYbV2Lpv+ejfn9kYfrq2gZ97xQ+O1Ni9cZKXfvaCmZNc1dkDNMlaI/UPqau5f/pcyrq+ViPLNTp+vzW/j6wxWktzcl2RmkUOGuwFc9D1/d7rkEy2/ctBp+Nr7Y5/DG53fa3ZzV2t1fPba/b9uvtkuVk+p2zkG/ibL/qhkdqWP+W42vYnYKvjc7DT9yvvo+1qPfh1c/h+XbW6q52u11zt9mXfry9Z9W2+d7Xlai8+vulqZ+/ZdbX6vY2Rr5M7j7tlcIp8+LCy7GvsQb5NOn7LtyNc2fG1i3657nO+31uX/H53rjVcbWvP1zbb/lxsdvy52CTXzA88/DNRGdeP5gohhBBCCCGEKBW9iAohhBBCCCGEKBW9iAohhBBCCCGEKJVS54jeakwy/6eMOaEx+5xk3mj8fv1OYufjiWpTxhy6WJTTajOR3IVJOArbo3Nfp5xv9p3lAZnLSCU4ZBKrkQmmrJ/YPM/YvmPrzhtcSjP+9lJ6r2Xzw/1Opn1LphIidu6IBKZYo/NByVMmq9G2kTmiRUHS9fHzEf3MNiCQOc8DMpe51/fb6xEhjRGbz7TvY9Og2KIukxWR67898P3Qhp9n3La4OaJ9sm5G5qbv93xoWn0/R7RLzieb105hA2hWKNb9nFYs+bmVgc0RZVKjJd8nRubD0lrN90le23a1pNb0tdTXGMwT0Bn449gl86Bj0SeiQgghhBBCCCFKRS+iQgghhBBCCCFKRS+iQgghhBBCCCFKRS+iQgghhBBCCCFKRbIiMRUkcimPQPp12tKTWTBPQgeW3Plp3XxBfCLR3+FkfcrGjTQyG2xd2pYSslbMMxXATFvGxYq0T4gEh7SFDeFcTMSWY42pDqyvidvlOuvGjTCsr1nGJ5EpcVkRWS7zR5xkxWX8eglxu0TbGCOvDyoI6pM+IaIjJj+K7hNyzubpPnY9Anw+BiS8fRLyXvDFJlqutm+7rtaCrzHa2HC1pHfK1Y51vCBnp8fkOr7GhFRGbhaWF/axXHfLUDERqYXMXyDW87ImliBrd3yRiJOsTvZR9weW1r0kKc99LUv8+U6JRWwiOevYawohhBBCCCGEEGOgF1EhhBBCCCGEEKWiF1EhhBBCCCGEEKWiF1EhhBBCCCGEEKWyMLIiKqFYgEnhMcTO258FrG2xEgSGJDDXZxLpy7S3N0kkY9edxfU7K6kWPVayX68KWEyipTmRGeBSI19Lp3wqqaAkVmRSWI4d67QlY+z6ZhIiJsth4hV2Nd8KYqJitnrk+JhIKLYfIocDSkrC1h2wtsSN9iyXaUpkJN4LgyQvZLxO8lKjO41qW/AeF6BDBEFEiBTISUtSUouUFVXlefN6sOx2ySDbhj8p+4mXEG2Hi365/hVXG5CT3M5Ou9pJnHS1YzWf8dsa/i56dsnLlI6tNV0tWSK5zAr7YEajBhEYra/75VJ/EYVu19VsQO6ge/u+ViOvcMX2AjBWS3w/MXFXbO4neY/RJ6JCCCGEEEIIIUpFL6JCCCGEEEIIIUpFL6JCCCGEEEIIIUrlhi+iZtYws4+a2afM7HNm9tPD+oNm9oSZPTz8c/+Rt1aII0AZF1VG+RZVRxkXVUcZF1UlRlbUBvCmEMKumeUAPmxm/2b4f/9DCOGfH13zpkOcgmG+hUhcJDFdQqRmZtoCozlgahkP8H0xzzKqeWfafcdyGiswYk2JXTf28qDymsllSnMxhk8iMIo/H3QnpDK+Ni0hy7E2F2vs3KZEDhFLIH0yIG1LE9LzA/896ECaQiVMZEF2f0omUp7dNFPNeDEfTEzUJd3KlqO5JzAJETufxDGCnGSrTwRGVFJFThOT+qR1IjAqyF2SZbL9nIlSYmVFfp+B5Zld40RWZF1SI9cguy7nQGA05YwXvmayouD7ug0v19nDNVfb6X7Fr9vdcrUTy/e62gP29a72NSe99eqeFd++Fy23Xe3Mmpcprd7hjyN/0aqr2e0nRgvrfplYwgpZt0+MXJ2OrzX88cdeR/PODV9EQwgBwPNnMR/+WezXDiEOoYyLKqN8i6qjjIuqo4yLqhI1R9TMUjN7GMBFAO8PIXxk+F8/Z2afNrN3mpn3FwuxICjjosoo36LqKOOi6ijjoopEvYiGEPohhPsB3AngtWb21QD+RwCvAPB1AE4A+Am2rpm93cweMrOHtnrk9+AIMQdMK+PbyriYQ5RvUXWUcVF1ppXx/f5eWU0W4obclDU3hLAJ4IMA3hxCeC4c0AbwzwC89jrrvCuE8EAI4YH1bHniBgtxlEya8WPKuJhjlG9RdZRxUXUmzfhyulJia4V4YW44R9TMTgPohhA2zWwJwLcA+AUzOxdCeM7MDMB3Afjs0TbVM89yoesRO7e4DDnRuPuMlxrFSUdmzVxnnNTYd4+YuISKb2YEl7mMt61Jchor2pqn7E76O7amnW+WtXEh3hWa+YycD5YCdi6JswXpjRo2hOWAyV2YEKgoPEnIMpPcr6jAiYmy6C58WwbB9wprX38CCdFRuTWmmXEmnOuTPuyQYBEXDhcYkeVS0jc9cvEPgi/mRK6zTJaLvU+kNZLnJd+WdG10H7bsHymtEXm1sU7pkNGA2JoSMpIEctKS9tGLiY7qdyIe9XMK634qMDIv0ukELwPq9VuudnL55a52X/iPXK2W+l68Rs7dbQ0msfPL1WpeCJSukDO1Sn6qeWm0Fmq5XyYlGW/6n6yw3K8blsg3vpbZxUZe15K4tAVycgfk5huY4Czyfj+JsDTGmnsOwHvMLMXBNfZ7IYQ/NLM/Hl4YBuBhAH9z/GYIMVOUcVFllG9RdZRxUXWUcVFJYqy5nwbwGlJ/05G0SIiSUcZFlVG+RdVRxkXVUcZFVTmqnyIQQgghhBBCCCEoehEVQgghhBBCCFEqMXNExRETI1+ZlYSJCjFIexNi8Jhk8nKsLGbe5VRAfD/EijxiBUbThopbIpeLYdqCrtjtsexyxs8kE5ZUBZbHia79CWQ4bE0uW2BZHl9aQpcryImmLUVhx8VkNH0qoPAjBhMuMVnJgAxUfbIglaVFisFmSYCXE3VJH7ZJrUUEID1yeD3m4CExrZEi215m/nwuE+ELy0ya+sZkS0RWtOa3l2zURr62NSJ7ySMfM/u+86zlJTMh67ragAw4RkRHSY1cb63xr/HYW8eijv4DcsH2iBhqEHxfrzfucrWXh/tdrWZe9LOU+aydWfK9eLbhs3Bmuelqa2farpbdvupqdvsJVwtnTo1+vX7MLcNkRdb1bcPFi762vu5rGZEfRT8gUmOaK4UeuS+QMY3L1qabaH0iKoQQQgghhBCiVPQiKoQQQgghhBCiVPQiKoQQQgghhBCiVPQiKoQQQgghhBCiVEqVFQVwsUUMsW/MsfIQpkeYZPrttIULMdub9ncRqAQnUhoUKzCKFb7Ml76iHCaRGtFzxyaZU4HI0UugWJuLmZmdkMvXwpT1ErFSp2lLAI6CYtZi+y8232wxug8yvhgR7mQTtCUhrWFSn4wJXwq1LCMyFnZdEKlRIBIJshiVFaFHxBeRd+L+gAhqqGApanPXEYjNX+YHhTa1SXcxMdGed+vQ5TqRAcxJx7b67H7pz9NS6h/vTpPzySQ86ZJvS3q85mp2alT4YmtkxVrkY2aLyF2aXjLDwmbE4JQQWUzSJsKhjEiNyPXMMh4rNVpUmDgugx9P6okX/7x08EpXW7Lc1RpE9LOc+f2uZL5fl1N/cS3XO66Wb5CTt06yurriayvLo1+T9lIGZNDY3XUl297xy5047mtMftQhA07H90kg4q5+l4jQyL2C3QOmLTDSJ6JCCCGEEEIIIUpFL6JCCCGEEEIIIUpFL6JCCCGEEEIIIUpFL6JCCCGEEEIIIUqlVFnRJMRKjtibdazAaNrECFpi28KOa9qT4pnchrKgkpWqwDwXsbKQMhg397HfFYvNPRW30A2yPPvFJpFvVYUA3zdMXMBlBuPvNyV9H4itx8h+UyL6yem59LDFUrJfJivK81FpREakKExMxPYZmJiIHFef1ChEShHITpiYiZ1bfn3HNWXeCAD6hcPpkH5lYqLtru+HPVJr9li/kv4nYVjN/UjZC77WSPw5voOcd5rxNSI1OuWFNHbbidHCul8m1LygBj0i7mq2/HI7e345vxSMCFqMWKK4mIhsj93DKiQhAnw/klihlvjicm/Z1W4PL3W1c3W/HINlvJYSWREZY9cyfxGurXrBVXrcZ9COk6wWxUSAkw5Zmwi0WG1v3++T1LBDaltEYMQg6w72vNSo3/TZ7Xd88Ht9IiZi95kpC0b1iagQQgghhBBCiFLRi6gQQgghhBBCiFK54YuomTXM7KNm9ikz+5yZ/fSw/mIz+4iZPWZmv2tm/pdMCbEAKOOi6ijjosoo36LqKOOiqsR8ItoG8KYQwqsB3A/gzWb2OgC/AOCdIYR7AVwD8LYja6UQR4syLqqOMi6qjPItqo4yLirJDWVF4cBgsDv8Mh/+CQDeBOD7h/X3APgpAL86/SbeHExqNInAiC2XliE6KnzNJspPWwYxiJxuTKVGUxa+sMnQR8VRZzxQvcJ0GTCbCdlt7LmLhV0zsUKuI894ZCaJq4ISm2e+XFzG2VjllRs3z1FmnB1HlwoO/LqDyGuDjbk5ja3vwYzlIGqvXNaTpf6MFMVEAFCrj9aSnIlSxhdvBXIQSdcff6yYjw0hKbkA00HcPsYV843DVPMdgF6hvzukr5vkwtwmC17reHnIdvBinn3z4pF6qLvai2zd1daJaebeVb+PO9a9BGX1ZMfVkiX/aGgZGZ1Wlka+DGdOu0XCyopfr+P3id1dV7KUjM5E9GRt38dGpC1JTsZccg0yiRijTIHRNDNu8PdWdq03Ml9b6/tMnmmsudqxGhuLSVtIbYW8maxmvq8bmb8IswaRVNUiPyQmg6pdvjpaoFItIhxqeoFR6BDDGckuBlddydaWXC1cI7KiLb+93o7v5E7bX1vtru/4Lhnvpy0jjJojamapmT0M4CKA9wP4EoDNEMLzvfoMgDvGb4YQs0UZF1VHGRdVRvkWVUcZF1Uk6kU0hNAPIdwP4E4ArwXwitgdmNnbzewhM3tou0e+ayDEHKCMi6ozbsaVb7EITGsM3+0r42I+mVbG9/r+1+IIMStuypobQtgE8EEArwewYWbPf457J4Bnr7POu0IID4QQHjiWxf1eISFmhTIuqs7NZlz5FovEpGP4aqqMi/lm0oyvpOTHpYWYETHW3NNmtjH89xKAbwHwCA4ugu8ZLvZWAH9wRG0U4khRxkXVUcZFlVG+RdVRxkVVuaGsCMA5AO8xsxQHL66/F0L4QzP7PIDfMbOfBfBJAL8+TgOYhGHaDCJEKTdD7JzccaUtAJBGTJZn25+ElJwKJtphUqNYCU6shKg8VRGAqWbcouREk0zsZnBpTpzAiLV22r9gOCbjZchNjBwtP19x++ViIracLy5mxs1dw1RWRAVGfmvR4wHJd06yEch5Y6KjRohLOJMVxYiJACBfHq2ldSJFibkDAwg9IlwispwBE690/LGyPAYmmIoUE7H7VRI5yE3pOpjaGB7gxVpt0td7JNCbREx0Eddc7Wpy3tUaWHW1V9XPutoDJ/05uXvZy39uW2662hoRGOUnXAnJSS9GwVrEJ8U9L2OhsqKN477WaERtz4gsxogsxogsh92IjC1WooToJpjqs3jxeY91VyP1xY1a7mqnl/xyJ73TCCtkfKqRsYOJ7Y7lPgvrDZ+FxkkiHFon2WoTYdbDj/q2PL098nXvGsk4abCRB+qEOZOYwYmQrPljDXu+Ld0r/vjbO/6cNVu+Ma2evyF1yD2gKHMDJpNyxlhzPw3gNaT+OA5+Rl2IhUYZF1VHGRdVRvkWVUcZF1Vl2h98CCGEEEIIIYQQL4heRIUQQgghhBBClIpeRIUQQgghhBBClEqkKqFciNNiIphIh7gH+LpTnrQeKyYqiluo5GjKthPmluAylvEFRnQ5Ktph+x1ffjRLeL+yyd7T3QfPB5NlxYmmGGwfWaTMJLXRK51ta9pCLi4SYvtgV2pcW6jgay7dFzdPCD67TELUJuKbFhl0mfQgNstMrsG6OSU2kmVmZmPrktzWakRWtOpr2WphDG/QiyWqHWCyohapERkIu6Pm5PwM+r6W9v26TOBkVGrkW7II3/lmsqIOyXiTnJOd0Ha1S8nTfh/Bb/AVyVe72qlG6mpXiWNlLfPL3dH3vZ01fE7Tdf8YaMeJYOj4miuFlYLAqOtlTXb5sl/v7G2+tuq3jzX/+y5tc8svR44fmT9+JiZisPszowzZ5lFgiJMVLZHxaTX3fX2m4a+Fs3Wf8TN1H9613GfmStubjlYyL+apE4ERI+z7ffSf9sKwJ/69F4Y9tnlu5Gt2z1onx3BqycuFNla8QGxpxa+bkr7LV/3Y0m+S++y2FxPt7vn+3G17WVGTyIpafX++u+T+QYbDaBbhviCEEEIIIYQQokLoRVQIIYQQQgghRKnoRVQIIYQQQgghRKnoRVQIIYQQQgghRKmULisqTu5mYqJYsQeTgjC4IMcvx7Y2yZs6k+sYq5H2ZQVJBltmEpESm4zP/B2sn3rkpDHhUFH4cLBcnHBoUcVEDNbXrG9i88xg+YjdHutrno+47LKMF8VEbHtMDDPtjAciJjLaT2xk8qMB8xzFjjfsyFht3nJfzC4TEzWJmGi/55djEhh2bTBy0qldItdgvbdCRCaMLCXSiGV/cPmGXzddH5VG2DIRqnB7lisF0lFGpEnYJdcZ2d6ASIgycn56xEDBru9kylKxWVPsMpbJ9sD34b7tu1o/eBnJq8JrXG2JZHKbnPe13C+3RCRVS0Tkki+THK02XA1rS67kxEQAsFxYl9iAbI8Ih778pN/+6TN++ym5ZlLy2Bp7HZHxZkB8N4GMaUxMxO8x0xURHgnmu4zECsSVRaVGxzJ/hCdqRODT8AKfVSIwWibZ7RBpTkKeGbrb5Dx9wWfwM58662qfvOaFWdu9G98rVjMv/jlLZEDniITpRNMLjJbrvu9q275Pul3fJ822lxXtdXxb9rr+OtqPlBV1iJyOOh8j0SeiQgghhBBCCCFKRS+iQgghhBBCCCFKRS+iQgghhBBCCCFKRS+iQgghhBBCCCFKpXRZUREmJmKSFbYcmyjOiBWF0HVJje2ViVxi21IUEx1sb3RWPXEBTSRySUkPUKkOmdzPPB9MYETbzCb807MxX4KWm6GYVS6lIWIQdi1EnmIjnR17cZP4XUdMRNalYqK4Wl4QwTChEctQLFwkxDJOAs1cRaQ4IHIONmaw8x0r5WFymFkR4OVJHSL2YGKiXSIF2SehbxMHD4ONQ+2cic98ba3nBQzsHNVqvtG1Y0QIdNLLIJKTo8IXW/GiCiRE+tDx+7SmF3qEmq8BXnKBge/QQZfI4DrkuiWyJnYPYzCR2aJQbDl9TiG1Lvw5OYk7XW0j93npkAFhmYScSQUZbGxOiHzGclLMxhRrDUj+2m2/z6te6oQdL5TBxrqv9clAQq4ZJvgKJON9IqNhtW6fLEfGPnrrmEOKLWe5Iu43KitazfxR1xIyTrKHDUJG1qWSNFLr7vvsPvXEMVd7dGfV1ZgYsPjswpZh0r1tIhLKE3/dM3oka91dLwtj+WNyoVjhUJvVWFvIfnsTSBX1iagQQgghhBBCiFLRi6gQQgghhBBCiFK54Yuomd1lZh80s8+b2efM7EeG9Z8ys2fN7OHhn7ccfXOFmD7KuKgyyreoOsq4qDrKuKgqMdPIegB+PITwCTNbA/BxM3v/8P/eGUL4B0fXPCFKQRkXVUb5FlVHGRdVRxkXleSGL6IhhOcAPDf8946ZPQLgjnF3WJzTS6UtZCIsm+LMBAIMPgGZLEg+H2YfGUeLiei6ZAI3qxUmdZchcmHCG3YUsQIjNmufeUjYquycHRXTzHgIXvrEjqVL+oZN9o7NPfNI9EmNuFyohCiPDBc7d0y+Vcv8mc8LQoI0UlAQC5MQUTER2S87skDGJTYWBCb9InsoK+FTzTe8ZKlN+qVJDnibiEI2ieWh2fMrd4NfLieiqGad3dL8cmtkwGLXVa3u25Kd8OKH5Nyaq9ltx0cLK0tuGYa1iIRoe9fXMi93MSJ/SkgfJ00msiGCMiYtozXfvDKZ9nNK8XDY+MqOuRYarvbS2ilX65OwpZEWqNjnnuixM0ZCBMD2iGComNUukRU1W762TwRG3ct+uTMk9zu+HWHH7yPseYFRjxxCv+2vZyaLYTV2b+8HNrZMfoFMM+OG4EQ/KQl0nTircnILjX0mbhFJXIeIdNg7wBqRs7F+3dv1QqBHt72saJcI9daIdGm9sFiTPVgR2PjABEG7XX/PYtf4Tjd3tRZ7xiF9wvqTyRJ5nv1y7Fk1Vr7IuKk5omZ2D4DXAPjIsPTDZvZpM3u3mR2//ppCLAbKuKgyyreoOsq4qDrKuKgS0S+iZrYK4PcB/GgIYRvArwJ4KYD7cfBdmn94nfXebmYPmdlD2+zbUULMCdPI+I4yLuaUaeR7V/kWc8w0Mr6njIs5RhkXVSPqRdTMchwE/zdDCO8FgBDChRBCP4QwAPBrAF7L1g0hvCuE8EAI4YFjmf89OELMA9PK+JoyLuaQaeV7VfkWc8q0Mr6ijIs5RRkXVSTGmmsAfh3AIyGEXzpUP3dosb8G4LPTb54QR48yLqqM8i2qjjIuqo4yLqpKjDX3GwD8AIDPmNnDw9pPAvg+M7sfB/6KJwH80DgNYJNyYyfRMsUIl/CQ5di6ZINGxCt+yjCHyQLYpO4sJSKOKYpcmDyFTfKeROTCBEZUoBM5aT/SnzAt4ctUM15sJvGH0Mne3cgJ4CzjDCbkoqeYLRe3CyfVAoCMZKaWeulLUWDEMx7XDtYnTAjSH/h99IgsIXYfAyLMoTIXetL8wR2R82Wq+S5ml7hwsNP1x3ut7TNwqdt0tau25Wr7tu1qS2HF1bLu7a527zF/ju5c8uKLY3UiCSIEcmyW+wyF20YlNeH0abIxkvktf6xW87dqI/lOWl7QEvZ9zYg8DORaXiCmmvHimMgELfXEF4/BZ/IYWblFBnY2Dtcix+uieAa4zn2fnPaw53NvV3wGcc2Lg0KnkK22lxUFkkk6mJKBxHb8+ICeX65/wf+oafcKkZ7t+uu02fRPdM2OrzGpTpfcUJkYJvZ+egOmmvEoIRc5TXnkOMGeMZlwZ6fnx7YGeRZokDFrr+XFRDvtuqtdbse86vA+2MhH83uMvADsEwlTLExg1Or447pCjmGr59flktS4tkzyvDUJMdbcD4M/H/3r6TZFiNmgjIsqo3yLqqOMi6qjjIuqclPWXCGEEEIIIYQQYlL0IiqEEEIIIYQQolT0IiqEEEIIIYQQolTiZvBOiQAvnGHyGlaLFbkwoQ2DTd6tk9dyI+/qbLK2n1rMxUQ5ERMxkUtxuYwsQyVMtE/I5HkyQZoJDxJmtyGE4JcLZDoDmwxORTNRe51PYjLeIUKDNjEaMKkRcTXQfqWyIl/ikgzSPpD5+KnF5bmee2FFvTZaSzO/LZZJmuc+6SciELB+XJ6pzItKQmJrTE5CZEqRMq9ZEYIfd1vkYt0jQp/NrpeinE/Ou9rl/uOulif+1w08kL/c1b7hLBMTtV3teM23pVHzopWswUYikr/zXu6S1L48Wlhd9evdfY/f/IoX3rAbmzX9cWHLS1ssI0KtJO5+Egsbw9n4vwgY/NjJpEHHSPF43Z+7GhmIl7K4vmGLRa5KRTrdfV+rXfI56l+65Ne95k9yZ290e4HJU3I/rmd1Mh6SWnrBt81q5Pnwit9Ha9M/3u7v+Sc1JrxpEqlOq++3x/qYjeFMHDVLDD5H7P7GHgPYkbBjbpG+udrxfXil45e7Z9mPxXtdv26TiI6utv35bJEDYU8CXSY4LBzbGnuWIXKldmQ2WD9d6/hnlyf22PZciY5VTLZGxxayHHuPYRlgzz2x6BNRIYQQQgghhBClohdRIYQQQgghhBClohdRIYQQQgghhBClohdRIYQQQgghhBClUqqsCDAnMIgVE7XJZGMmbelFyorYRF028ZdRj5x4nhGpUUYmNTORSy0blUkwkUssA9Z3Uxa5BHKsA9KhBtbJfnuxy80jxYx3yaljYqIWsQVQgRGpsd5iE9S7kZP2G5Ezz3OS5yUifVlq+Fp9aTT3SU7ywmRFpJ8GXSLJIMKchEgAmEgiJceVMoER6c9Fzu6NCPBjNs03ufY3sedqTEzU7nvxzxvyN7na3Wv+XD7ld4Htbt3V7jvml7sNu65W2yBj+IvWXM3OrPsN9kfHcPvUI26RJPO34MEdd/ptrfu22eWrfrnI6zawMYRcV+zaYJkfMIETE+exttAWzg4DkBfGnSUiHFrNfa1BhG6x8hDWX+x5hglF2Li+1/HCnZ1rS67W2vXPHxc3fcaf3vUipuKz2jIR1W0QMdjxpZarrTT8cgmRO+ZEftRu+Wt8v+2Pn0mI2kR40ybPQp2BP7lsOfb8Ov7T29FRjHSdPMPFyvjYc/xXmr6vnyTCnRXyFtIhAp9tcu6YLOpa15+nDjkBfKT020sLjy5M3rOc+dwPiMGN9dMOeXb5wravPbnrn6EaKetPX2NyNDZW1en45ddl41c6wUiuT0SFEEIIIYQQQpSKXkSFEEIIIYQQQpSKXkSFEEIIIYQQQpSKXkSFEEIIIYQQQpRKqbKiEIDeIEbk4ifHtvxcYHTIcmx7DOIeoJN32cTsBtkv34dvTIOIiRp1InJpFEQXsbIi0rY+mQzNJDAsDkxWMSCmCzYJOyH7YDXalkgh1LwR4AUTTF7AxER7PhpostwTEQwTbbGMrxDBRmK+tkxWpokhMoklIp1YXve1fG10i0baxgjkIu83Sa724ybP94nwgAlZEiKmYHm+1WDZaw/8OWrZvqv1Bm1X++rsm11tmQgYPnXVC0/ONGqkRkQNREbFshxt61htuFI4sTHytV3zEib7yKd87RuISI70J1iNmMzY9TLwtxx6n2BSO369+BqT2rHavGHmZRzLmQ95jxwLExOtkXVXyP2cfSqw0/NVJrBj9989IuG5vLPsalsdf818YcdLjfbJPau43xVyrKe6vh17PS+e2Wj7saBOJDDbbd9edvxUhMmyG7kuyy5bjgmmigLDWXOQ8dGGpqTdxsSUZLktMnY8sctEpH7dDX86sU/FmUQqRe7TVzt+3f1e3PNMk5ymViEzrYHP7rHMj5Ps2WCLiJSYmOgjm5uutmNbrnayd8LVjvf9vWiNiLZ6RGBEc0/eixgWKcpj6BNRIYQQQgghhBClohdRIYQQQgghhBClcsMXUTO7y8w+aGafN7PPmdmPDOsnzOz9Zvbo8O/jR99cIaaPMi6qjjIuqozyLaqOMi6qSswnoj0APx5CuA/A6wD8LTO7D8A7AHwghPAyAB8Yfi3EIqKMi6qjjIsqo3yLqqOMi0pyQ1lRCOE5AM8N/71jZo8AuAPAdwL4puFi7wHwJwB+4obbK3zNJvx3yGT8Jpkov0/kLmzdWJHLMpm8y1hO436iuZb6ifZc5OLNEfny6IEYOVPEGUQlFEmbiFxacZIVNjGdilzo1O/FYJoZD/B9xoVcvrZL8rxDAr1PhCRMDpMSCdEakY+w70etETkMYyn3gVta87XGbWSvx+sjXxsxfQRy/IFYnYxca4H0yYCMIz0irkmMCGPmyzdxU0x7HD8M9flEyjmW05Ou9pL6uqs93vKihnXz4pUGGdiZZKRDpCVdcm0EYtcILZK/XS9OsnRntLC24pbB5o4rJR/9tN/nS+7w627t+eV2vPClv+2vje6+P/5Oxx8/65NOn4gvIiUwTNoyDYHRNPOdAKgXRC5M2NFIfDbYMa8RMdGJGpEWkjHsMhHzbBIxDGOLSIJa5Nx9uen3cbkdJ2JqMMNNASaeudrxwpculU/6dS+14xybReEUABh5TmHj1zwO9dPOeI3kt0hORXS+d6514qRNTAjKno92iFyoH8iYHdkWJjtl7wUD2iWj28vIg0CWENEbOa7zxD75Z91PutrV/hOutpqd9S0jzykpEZwlqLuakec+9ozDZJbs/YnVYrmpOaJmdg+A1wD4CICzwwsDAM4D8L0kxIKhjIuqo4yLKqN8i6qjjIsqEf0iamarAH4fwI+GEEY89CGEAG5Dhpm93cweMrOHdnpe3y/EvDCNjO8q42KOGSfjI/nuK99ifpnGGL6tMVzMMXoWF1Uj6kXUzHIcBP83QwjvHZYvmNm54f+fA3CRrRtCeFcI4YEQwgNrmf9RKiHmgWllfFUZF3PKuBkfyXeqfIv5ZFpj+DGN4WJO0bO4qCIx1lwD8OsAHgkh/NKh/3ofgLcO//1WAH8w/eYJcfQo46LqKOOiyijfouoo46KqxMz6/gYAPwDgM2b28LD2kwB+HsDvmdnbADwF4D+P2WFxMnCXTO5vEaEI8ZNgh4gkmmSGdIfMQDYyHX0lZ7NtfW01UmpUz32jl9a8rKhx2rcvXS9M5if7ZCKXwZ6fDJ3sEOELmUidkEnObPKyGelPupyvzSlTzXhRwNUhE+rZ5HkmJrra9hm6Nmi62rZ56QnjZG/D1Qz+u6PHiZmCSUUaTLpB8pzds+b3e6YgpcmIyWDH/wjR4IqXtABeFjPokmuB9HGS+hrL+IIzlYwbgLTQN0zesERkV+ske7fhtKudb/lzmZNbFRteWsRAwe4n7L7TI1KVgR+uAZIrNL0kyA2yV7f9MsSU0v3Es66W9/0+B89c8+s+69vRuuKPq7njBTXNFql1fL93iMCISWVYrUfGwildaVMbw82AWkFgxjLOZCdMzJMTKUxOBGmstpr5887EjUwg06WyIr/c+RYTSLkSlZEUH5mWyVjKPu1okX5q973A6NFdfwzE90IlOEtMMEXkSkzClJPxPyM11idMflQcM8dkehmHF3KxnDK2uuz69wfNZFHs0ZntleWU1ZicdJsIOzvkYm1H1lr90Z3s9f0OrsI/f11MnnK153Y/4Wrd7lVXS7NjrpYt+fF5P/O/qWc/+Oe5GpGU5Ym/3nIS3tjzyK6jWGKsuR/G9SVi3zz+roWYD5RxUXWUcVFllG9RdZRxUVUmeIcVQgghhBBCCCFuHr2ICiGEEEIIIYQolbjfDDwlAvwvfGa/9JXNn2NzRLfID4hvdvyCOwM/b2ZAfjJ9vdPwO4H/uex1Mpc0dv7c0ikyN+SeVVezk75WJOz4uYJ2xdcw8D/TnpH5tX3Snz0ybyUhv0SXti9yWkSYwi80nxeYO53NJWIZ3ycXw+XBrqs9mzzuatc6/hcgDwK5Fur3udqJ7qt8Y8j3qBIyz4XV0jUy1/hFp1wtfNVLyH4L61287PeZnPcLtsl80G1fm2TeMsszu+4D+ekpthxfd74x83NGlsldZI2Mkae6fu7KWbLyU7t+vE5JHllfsbltbDk2p6nZ9XNm2lu+ffklMh8UZP5nsW275CZGBofBvq91P+kz37nsx4vmVX8Mu3v+l5m3On65Vo/MKSTzQXtsPij5ZfNsPmjx/g/w62CWJAh0bmaRogsA4PMDWSZ3Sb9ukTmdxXl8AB/DOnTeqF9urxfn32BtJvoNd+76JAc1cgxsjvaz+7720FV//6ub76dTdf+cdrzu28LGJT6X1NfqdP4cmf/L5oOyYMyQxIKbz7uS+SCwecuAH0+eNTY32q/Jaqxn2NxPcoun2WWemOI8T4A/b+32/LPyJkadFFeTC26ZKz3/TLa9+6Sr9fqbrkZ7gEhcesELCzrwDo22+ftTJ/hz1h74kNfZXHNyLmokz5M8u+gTUSGEEEIIIYQQpaIXUSGEEEIIIYQQpaIXUSGEEEIIIYQQpaIXUSGEEEIIIYQQpVKqrAgAivNe2QR4NlF5j8yYvdr2E4ufhZebXDb/i2U7RAJzCl6e0mi+2NVetBLXbXnuZ1dnZ/2kevuqc64W7rp9tND2E5DtKxddLcEVvy1ixkn2xhe5UCFIpMglVkzBhC8LQfD9w2RF7Bcs7/b8zPtriT+fl5qPuNpey8uKlut3u9p9eIWvHfd5vmvJt6WR+szst730ZP9Jn9XVjedczV5818jXg3tf5pdp+En2yba/dnHJ/0LpWEEEk2X1yaT9AckkFxORfdB145gnmVcCL1A5lhMZSYO12eeMXRtnG36MZL9ovE8GHSbDYPvY7/sF94isaGfb5w8k3+EJf7102qMyiEHfb9+I7KTe8Nt65sKGq7GxtEOOq0vkQkwaxEQzTDjEJVHjL8d0KLMkMaBRELnUEj/2MZHLDsnQ0/s+QztEDMXGl+PkGYJ9esAy3iQCI+bLYvci9gy21fHL+ScQD5U1dX1DPh4+6WrNcM3VTif3ulq/5Z+h+kTQwkRb/VqcVCtkRCzGnlOYkCfW3FgSBqBWyO96zT9Pr9X9WMekRo/vbUTtt0b6ht2md3xTqKxol7wXMLkOy/he3x9HUUwEeDnRZv9pt8x+218J/YEXCTEMXhqUJOQ9gQr7iNQIvqP6ZDl2/2TjCKvFigJj0SeiQgghhBBCCCFKRS+iQgghhBBCCCFKRS+iQgghhBBCCCFKRS+iQgghhBBCCCFKpXRZURE2KZxNNm6RGbNboeVq5/FFV7u0+2lXYxN/X7T2ta52OxETnaz5BjJxQbPlxQXd8160km75CdLhPz4z+vXSklsmMX8MtuUnSFvGJk37Cc1sPn2fTe6n8gs2kTpOAsBkFVz4EldbVPrkaPax5WrNziVXM/OT29/Q+M9c7eXrXuCwTEaBWFkUEw0MukSi8njT1Wr/9mMjX9vGht9YveFrqZ/cj75PUeiSyfgdv26362sdIhPp9Mi65FpgNSZ9mUTmNSvMApYK0qrUfD+vpv44zjR8v1xsxQmg+D3B17pk4GBjSYecj8ttIojYWvPt2zrmahda/roq7ncp9QdR7EsAWCYykM+SdrTJMSyTfeTk/OREksSuZYscYSO9YAvxnW+z4MRs60TacmzZP390un7c+ErztKu1yLnLSR92qTQtDiYcYtcMe7baIxbJ7a43yGyF0WeLS8l5t8y1vpdFXtn1wr0w8P253LjL1ZaWjrvaJvw1WSPjep74Z7KUhJfFOSE2R3Z9pGS5fN5kReafWZmY6OSGfzY929h2ta80/fPpBfL8y8Q3XXLP2+v5GhObsu31yIDfI/3fDf5iaJvvg3YYfWbvDvyzzCD4Mdss7vUqSfwzTspkReR5nxHIs/2A5JQ+75MikyVOIiZiLMJ9QQghhBBCCCFEhdCLqBBCCCGEEEKIUtGLqBBCCCGEEEKIUrnhi6iZvdvMLprZZw/VfsrMnjWzh4d/3nK0zRTi6FDGRdVRxkWVUb5F1VHGRVWJmU37IIBfBvAbhfo7Qwj/YNIGJESkwKaKs4m1bfiJxXvdC67W7/vJ1S85/h2u9pql21ztpHdQUNFDRo6DHVt3yx9b9hnf5mz94ZGv+2/8q74ha6u+Ribeh563FjChTJeJXIigpdcnghZWo1IjJjDyy5UscnkQU8q4GVD0tGSk2Rn5FlDdfF834M/xUs3LL+5a+jpXu5sIrrY7LKe+gaukgew8ZUS2ki/7vCVLJG9PjU76r/0/H3bL4Gvv9bXLXuDUv0IkA5u+vft7XgKwTyQ1bZL7zsDX2iT3POPji7umxIOYQsZTA1azoqzIixqYvI191/OLybKrMQEFcR+h1ffFqx1fY+sy0c+ltpdrXO34W+Rm1x/Jlbbf3lIhLqfqvk9q5F7ybNPn7DNXfR+faPjljtd8bd3HG0up32+dtIWLjvz2mLCK1ZjUiC03Bg9iSmN4goBGOtrfTEy0fspLS/J1fyz3d3yuPn9t3dVic98kuWciFyr4Igt2iHlkr+fzdi14cc3F5JmRr690v+SW2Wl+2dUGAy9tBHx2+wMyrpN1W+aFjK2BH1uaZFyvkVDmTDhEcp+TdZnwZRAp/7sBD2Jazynw12I99+e8tuLv70t3+O1984oXUj356AlXe2Z3xdUuk3F3j8gC96hMzcP6nxGYwIeKfojhq0BCxERMQoTgL8qULJel/tktM/8yksL3nZHrqAxihXV03RstEEL4EICr4+9CiPlGGRdVRxkXVUb5FlVHGRdVZZI5oj9sZp8e/riA92kLsfgo46LqKOOiyijfouoo42KhGfdF9FcBvBTA/QCeA/APr7egmb3dzB4ys4d2euz3WQoxlyjjoupEZfxwvre6yrdYGMYaw691/Y/cCjGnjJXxTY3jYo4Y60U0hHAhhNAPIQwA/BqA177Asu8KITwQQnhgLfM/sy/EPKKMi6oTm/HD+V7PlW+xGIw7hh/P/fwsIeaRcTO+oXFczBExsiKHmZ0LITw3/PKvAfjsCy3/F+vBi1vYZPwaeT2ukQXr8JN3G9mGqyUrr3S11+df62pMIMMlAONLRhLS46Hj1+792RMjX6cnNtwy1vYT+cMVP5G/t+knW7d3fEOaLT/xudX1y7XIRHImcmH9xIQ3PbIcm27OakclMBo344AXVzWIGGQt9+0+Uff9elvzdldrLfk8f13tZa52ueXlA10yWR7wNpMzZJ49O5+MpE4kJSsk+AU7R+tjV9wijZ7PbvdJLx/be9Jvfvuqf6DcafoxY6/rc9/pEzFRtHwrTrTFzkQoUWA0TsYTBKxko7laIZKLlVrH1VaX/Xi13vC1L24dczUmtGEinUbqz+UmE7ORTt0i+W4SGdWXd5nc5cZCC+bQYOv9h8Efu1pu/sHxq3b8PeyOBhG0kCwfI+PPMrGqUakRORCyOSpyYVI/I5KkaTD2c4oBWToakEaj65ZjYqLaPX7MeeXGpqsd/7j/1PXClhfTXWn78apP5C5szGEZZ2IiBrtP7Jt/tmiGa6Pr9bzQKNB7js9kkvhjTcwfKxPKdM2PN+3gx6X2gIz1ZFzvRvYn8UCCPOIcGeNmPIA/d/ntx8lEa7f7+/vdwU9nzZ8gMqBNL+7a6flO3OwQWSIZY1jEE9LmlHwOx0Q/aSGDmfnnJSYXMosLQkpynyVEVkSWY7KihBwXO/554oYvomb22wC+CcApM3sGwN8F8E1mdj8O8vwkgB86uiYKcbQo46LqKOOiyijfouoo46Kq3PBFNITwfaT860fQFiFmgjIuqo4yLqqM8i2qjjIuqkqJP0QghBBCCCGEEELoRVQIIYQQQgghRMmMJSuahKJ0okFEBUzkcrzu35nPdf3k/v3+q/265uUXNTLJuU8mnrMJ6n0yGbpLJrz3iCSi1yYCh5bfiRUkFvb+T/llVsjE+y95WcDus3657R0/GXq37Sdht3r+GJjIpTuBhIiJXMoUE02TBECtkOm1LE4QkSe+D5ezNVe7rf0aVzMyGX2r70UwDfOXPLkUaF+zjDc7JIPbfrl0xUsKrCBH6e769S6910tCdlsbvh1EqtUkUi0mHGLHRYVDrjKZXCjSGzJXJBawXJAVrS+13HLHT/hfD7B6r9/e2TU/Xq39sc/t5r4fr5j8J0+8ZasfvORhp+fP2x6pXfFNwYWmz+TlgRe5nE+eHvn6Yvvzvh37j7paknjh0MbKV7napfSyq9Xat7maEanfgIzX7L7GxgEy/KNBv6XNxER+qWwORRrF6zoh0iZjpkUymKYv8b/a8Vxty9VWPuPlLrULXuTSHRCpUYfIf2y6/TogI9uAiogK7Uj8cwUjJddumvp1mVCGMSD565NBg9bI4MzuCfQ5hdXmbKwPAegVjqdLnvUGfXJ/axEZYeKPOvXDDhXWbex70dTVjr93N1LyXBEZcSZOy80fbx7880ytIIprk/E5T70YKyHPWgwmK8qZrAjk+iCvcGnwx8WeD9nwwJ4Fo5fzpWj0iagQQgghhBBCiFLRi6gQQgghhBBCiFLRi6gQQgghhBBCiFLRi6gQQgghhBBCiFIpVVZkBuQFkcsqEbkk5ic+M7nQWu6bf65zO1nXt+UyEQR1yIzyDpms3Y0UuXT7vjboklqbSAD2R9ty+Qt+AvJeiwhauqdcjQmHWkw4FCltoSKXSJFQrMhlUTED6gWxRZ54Uc9q5vvrLJncn3pXEZWKfPQqEWMN/IR3Mt8fK6QtRakYALRIPq42yT4u+LO8tuslBd3uaJsv76y4Zd5/3os+WNLWMr/POhGhRfpFYERDwZZjsMUS0p+M2H3MCjMgS0f7eqnhZRP1dZ/55KSXPCQvPuNq97yMjM3//mlX23nCj3+1y/6C2SdjXavvpRT73jeBKy1/HLsDf7yXkvN+3e6XRr7eaz3nd0AIRADT7Xup0162Sdrmr5dlIvLKmTSImIRYHpkEh11XxO9DawNmnZohIQD9gmgvEPEeM9AEkhdre7mVLftzsnS66WrHd31tm0gF19rkWkh87tl5Ik1GQkaxFP6aSW20lqVeqDII/vgDExMRaUvGpC3ml7MJPlNhIiEqpotdd77iTAnw4r5m15/fbpMIffbI+SQPJWSYRELuyY2MPR/5M7BC3lZqJNApGbTY+czJeFcjsqK8IHurJ14WFuCPoT+JrMiIEMn8tZATEV9GXusycn2wa5xJjdiVNe3HFH0iKoQQQgghhBCiVPQiKoQQQgghhBCiVPQiKoQQQgghhBCiVPQiKoQQQgghhBCiVMqVFSGgloxOGy5+DQBrkSKXnKzLtvf5bT/JNwT/Dt4ms8xX/Nxl+vbORC7Xmn5CPi770tqeF7nsNEcP+AvXNtwyn9z00oKNmj+GNSKEyok8hQlquMjFw2QssRIYBpVkkOUsUgJTFgkCGhEZr6d+cvtK5m0px5dbrtaoeVlAntzmak/t+/DuEiHLkvcRoE46e7fni8/s+4xfavmLdXDFb+9CQbb1EFnm0/sXXO2e7KSrnSMmgw1/eVDhQYMIFDJy/Oyayeh15NflghcieCBxzkj7ZklSuIYzIpYw0veUHjGlrHtpVe3Vp11tDZdcrdPZd7VTLZ/RYvYALmFr9uL0aj34MXwwGL1OjcgxjHRUSkQuTMbShx8H2qQdrYHfXouI9FpEOJKRgZjJ/8gtGznz+/hStOiuLALMifu6VDJI1iXGq5D5BUPP97WRcbhW99tr0HuH79lGymRFvq87fb89JnLJB/5+UhS3dFIvcmEwIVeS+O3XUj8WMFlRFvz1zGQs4oBBMLQLY8AekRXt7pG+rvlzl3V9rd8l0iAm/SKwe+MyMZ01yM2WiU33yH2GjW11kvFaGH1/6Bq5nolgtZeSlwdCQiRgNSIrqoPV/Niek9e6lFzPbCygz+yRgsdJRIv6RFQIIYQQQgghRKnoRVQIIYQQQgghRKnoRVQIIYQQQgghRKnc8EXUzN5tZhfN7LOHaifM7P1m9ujwb/9btIVYEJRxUXWUcVF1lHFRZZRvUVViZEUPAvhlAL9xqPYOAB8IIfy8mb1j+PVP3GhDicGJXJhwiE7GJ4KWDSJyWT/VdLVXXvUTrv/06XOudqnjJ/cXxRwAUCPyECZ/OE8kGRdJbf+yX/ex3dFT868uejHHCpmofPeSn9x/eslvf53Mo14mUiMmpmDSlpwJX8jkZTYJnc1xZoKWATkX2XSEBA9iShk38/llEqJj9Y6rHV/zopWNO4ms6NXrrvZdjz/pao9+1N+TvrS95mqbHT8MtIlUoNX3tc2uX3ez45d7bNtf05/uPDXy9SO7/8otkxFxy/7KG33bdu50tduWvAhmvebbtkKCukxGRnYtMDFCrNSIze6PFRiNwYOYUsYH415zHaKqaRLjS6T1IFnx43VjyV9rS+T6q5Pxhe2WCh3I8WcgUo90VHKRE5FLz3zQWObThNj6CH3z11mfiGGIlwjENxK9HLsO2LohHKl460FMIeMhwMmKOmSM7DX9GJ7ukc4J/tmF7pd4uwIZh/PU74ONG1ya5mv0vkoWXCLSqwYKsqLkmFuGibYC0VYl5Fqom79f1Yi0JYcf6zOyX3btTiJZiWVK+3gQUxrDAaBXkITtknv5NhFuZiR/jY7POJOQtcl11CcZZ/JLJs2JzXifjDtUyAV/T6mH0bG3a15+OiDjbhLpgk3IPmvw+8iDzzitke2x+1hKnz98++i6ZcuKQggfAnC1UP5OAO8Z/vs9AL5r/CYIMVuUcVF1lHFRdZRxUWWUb1FVxp0jejaE8Nzw3+cBnJ1Se4SYF5RxUXWUcVF1lHFRZZRvsfBMLCsKBz9rc92ftzGzt5vZQ2b20FbX/+ihEPOOMi6qzgtl/HC+r3X81AchFoHYjG/2NIaLxeNmnlO2lXExR4z7InrBzM4BwPDvi9dbMITwrhDCAyGEB9Zz/7P9QswpyrioOlEZP5zv4zU/d0WIOeamM76RaQwXC8NYzynHlHExR8TNpvW8D8BbAfz88O8/iFnJEFCPELmsEpHLxpr/Tvyx273oovE1frL8sTNe7vLtf/iEqz362ClXu9D0D147PT8ZeK/n3+kvtb0R6FLbz+h9dMv3wX/of3Tk66d3PuyWOb32Na7Wa77W1ToDL8loNfwxHMuZyMWVQFZFjUxM51Ijvy5bDmQCe0bWZaKFKTF2xosCruXcn9+1JZ/d1eO+VrvdZ8hu8xKi/KVevvWKk0+62pmP7LjaF77ic//Mvs/9xbY/8c/u+/P0+LYXFzwcPutqT+38fyNf9/u+bWnuf9JoK5x3tcthw9Xylh8L2JAXSHZZqgZkhv6ASBACm7VPsmtkXTbf36Yj5GLcdMZDAPoFkUufiNqCjzxCywsdwh6RFTE7ExHDsOWMCNeWyPXXSOPkauxUtuG3txx81hrpaK2f+2PNBv46SxMvoKilXkLHBElMxsLok75jmafCIbIcvV5ojbXvSAVGN59xGDqD0bGu2fHj8MqeH0vSuh/70h6xEBH6TTK+kPtgar5nmXyLSRXZPZRJW/b7RPBFxs5aGM1v3fyzBhv7mKwohe/jnIhh6iFO5JKSHTPxChtf2Sc0saNwGfKjQ4z1nDKAl47t9f39fbtN+pXkr0fuAUw41CZCpOK1BlxHYEmFXOx8xtEa+OsyJddCFkbblwc/7vaNXPckz4yEpI3toxa8OCon1ySTdLFrPFZMRB9nSG2S2Mf8+pbfBvBnAF5uZs+Y2dtwEPpvMbNHAfwnw6+FWEiUcVF1lHFRdZRxUWWUb1FVbviJaAjh+67zX9885bYIMROUcVF1lHFRdZRxUWWUb1FVJpYVCSGEEEIIIYQQN4NeRIUQQgghhBBClMq4sqKxMAPywgT6BpEVrTS8mGJ53dfyM775dpJMlr/zjCutfpufqHzfv/uKq60/4iUUj13z8qPNrp9c/BUiH/jzTX8cHw8fc7XzOx8f+drIZOPuwAucrmSXXW25wyac+/Za5PclBkTuMiDyDyptYdDd+u0R9wKdcD1LDEBamFRfS33WajWf+2yFyFeYGYrJL2p+Ynzy8ttcbWP3aVe7u7Xlajtdv71nW74t19peXPDk4JKrXek+5mpFDORYCd3g1fP7tutqe30vtWh0iYwgIRP5SXZ51NjkfnIeyTXD8swuhWKeZk2/IFDpE1HbgEjZBk2fW9vyAh/rEKlR19cGTSI8YXI1IiuqJX5dJlIzIm9ggrSMZLdho/eJfuaFFv2BP/7E/LVXtzVXy4msKA3knhgpoGDiLTrWM7fcfEV0IgKAbiHjTTIeNpu+lu6QXIU4WVGPXDM9IkasZX57DXKPYXlmMhJGi9nGCPWCVKVjxMZKdhkrK6qByIoQJ23JzfddrIyFdVOstIWN4cnRCrnGole4tltEOLTTI0KuDpH1BHZfJVKjgV+uTfbL7nlMvsWe/zJyUsjQhtbAZ5zJ3tLC2M7GeiYXGhCpEyMjY3YOIqxjMi/SloyKiZj0jNTomBFXY889segTUSGEEEIIIYQQpaIXUSGEEEIIIYQQpaIXUSGEEEIIIYQQpaIXUSGEEEIIIYQQpVKurAiAFSZtZ6mf0FurTyByYXS9JAKnNvx+79tztZNXdlztwu6Kq/WCn0C/1fFtfjZccbWdrpckJcnoZOUQfD+xWhdeftEK/vibPSLEIDOQs8jJ/QldbnxpS0ZlGqw2ZxIA8/2Tkkn7WeZrFns1tkme97y4CgMihDjt5Q+rG17+s3yFySr8hPxW3/f/vm37phD5RZqMXjNMlpWlvr1MatQzv/0ukYS0SZ+0ej5sOQl57KR9Lgbw/US6jm5vnhIeYOgWhBMdIlTptkgf7JHM5z7L1iOShx6RoZGIEk8DpSjNA3jfM+FLjeRvF0RIVhCt9IhwqJ/6+waTttTNS/jqwV8bOVmXCTgYbAyfhHkbmmMJAegX7lNNIm3Za3mhSEZEQowkJ3nuM3kIGcNJdtk9mYpcov2BfsEmyH2nuE8iXqkzgRGBjf8s40VBEsClLTmTtkSP63GyOlqLlB/NkoOMj9baA9/I/T4RPnWZ1Mzvgwm0QqTcj9b8LmjuY2HSuTb8TaVPxvYibMxOyPMHe3bh8iM/tsQKuXIiX2TPM1xmNn5tkojrE1EhhBBCCCGEEKWiF1EhhBBCCCGEEKWiF1EhhBBCCCGEEKWiF1EhhBBCCCGEEKVSqqyIETsZP1ZCgQ4xWOx7gQ9lyU8Qrh3z7VvKmMjF0yY2ki46rpaYn+icFmRFbJks8ZP2GT34Y+gQaUu77zu5TXxQGZEqMAnCJNIW4iahE6TnjsAn7rvFyDJMvhJaZML/XsvVLCMXCGlIILWs5ms1cl2y7mfHkZDJ9/XUi1oGg9EDDgmRr2ReDJaZz30S/PEzGUGfNJjlr0/6qSgwue66kbWYnMwbTOTS7vrbSKdJxAp7RExEZHUJ6ywGzXekmIfVIscwBhNOFK+D3Hy+MyIAY+KLRvDXQY0I8phIiY/DsQKjqMUqRQCckIvdG5tdf55qZLxmXZ3346RGKZPaEVlMauQ6IuOfRcoH2dgZyHNEvyCJ49KWuIc3KumiYqI4aQvLOBe0xAoZfY3ddhfh050AP453iCNuv+ePhgl3YoeJnDxXsGeNfmRmJqHP8kzERL0IWVFCzjp7DjJyXBnLM62RjJOep+N9ZO5jJY3sfE9yr1iEa0YIIYQQQgghRIXQi6gQQgghhBBCiFLRi6gQQgghhBBCiFKZaI6omT0JYAdAH0AvhPDANBolxLygjIuqo4yLqqOMiyqjfItFZhqyojeGEC6Pu3IgApBYkQubXR32vJjItnf9uimx8PSJGIB8ZpwnfvIyEwOwybtMMJETScsgG5W7pERMVEtWXS0jE/4ZsSKXHpF/xMqFJpG2sAzMUO5yUxkvZro/8CHqE/nFwLtcMGgS+cWeF15FzxTvEvkRkRRkRBhGZS7k+lgOPpcsv0URkZEf0KilRNJiy377IS73LFdcYERyT/qJnEaaU1ojU/4Da2A5RGWciVw6fT+Wdtr+1lLbJ1lO46QtsbK6JPKOlhBpGruCWC0lOWXjaTHPORETMaFFGvxB1OHvG3WyPSq0IFIKJtCZtpco0odUJtEZL97POmQMb/ZIxjtMPEgGDkJGxEQpkXllpMYkgLFyPyYwYowrbWE/eJcGn9NJpC05GSC4eCVW5BJXo2MGEx2R5Y6Am3pOKaaoS5716PhHbno5eV7IBkTgxsRd5PrIiUyTZXySIYaN2R34h7BARGBFWJ65wIj0XXTGSS1h54IIkaLFXb7GhFzTlnTpR3OFEEIIIYQQQpTKpC+iAcC/NbOPm9nb2QJm9nYze8jMHrrW3Z9wd0KUzk1lfLPbLLl5QkzMC2b8cL63lG+xmERnfLun5xSxcNzUc8puXxkX88OkP5r7jSGEZ83sDID3m9mfhxA+dHiBEMK7ALwLAF65dtsC/uY8cYtzUxl/xaoyLhaOF8z44Xy/XPkWi0l0xl+6ck4ZF4vGTT2nvGjpdmVczA0TfSIaQnh2+PdFAP8CwGun0Sgh5gVlXFQdZVxUHWVcVBnlWywyY38iamYrAJIQws7w398K4GdeaJ0AIBSmFw+YFISKXPxy0SKXzP8YgtW9VCBst3yNSZIIVORCJggzwcRSevyG28+I7KVuXgqTh5pvW+SU7lhBEJUVTSJyIe3j08P9ysU8TZNxM95zsiKS8R4RGLVJ7nd9TyR1ZjUiE/lzP7k9EGHMoB/Xh+y7VizjLIONZN23pWBJSsxfk0xMVIOXe2Vkcn9s7hk8p3HLsW81z+u3n2824wHmRE6dHpEVdf2tpdP2g2mSs54hmfeuHhgxKzDxlkVKLmJ9XyxXMVnLA5MVMXmFv37YNVWjYiIihmGCFlJjx0/VMxMMuUzqd9Tc/Dhu7r7SJrKiGrnBNUnumayI+4H89VGrk2ccKiYishhSiz3HjAG5LovCF3YdMPnW1KUt0bl3pega7TvWx5Hiymkx7nNK8TmOqaiMPJ93yPNMiwmMSN/kRPzTIPIto+fTty9WIMVgsiImJiouFysmYs8kKalldN3xx3EuHGKSrvGvj2nHeZIfzT0L4F8MjWsZgN8KIfzRVFolxHygjIuqo4yLqqOMiyqjfIuFZuwX0RDC4wBePcW2CDFXKOOi6ijjouoo46LKKN9i0dGvbxFCCCGEEEIIUSp6ERVCCCGEEEIIUSqT/vqWmyN4cUufSAAGTO7CZEX7ZJLzthe50LftuhcDhD2/LpMkMWIn7zLBxBKO+QXT4pde5NKAlxXVgrd6sInPsVCVSKSgJVrkMq8ml7HwoovuwJ/zHpnc32uRSetNkvEdIn1hsqIaERMRwVeIlBVxwYaHSSdq8NKhfjJ6vWXmZS45fJ7rwcuKmMzFJpB0TTv3sftl4q55Ux05WRHJd4tIW/KWH8PSJO7YMiowYusyKYNfl8ldmFAktuctMPnKaJ6Z0IKJLzJy/dTJ+M9FLmQMiRQTses79pqfRILDhC+zJASgV3gGYWNJhzy7tPrknPRiBUZEQkRELqy/crK9nJwAJh5hEjYu5PIbDIW1mZhoEmkLExNFS1vo8Y8v7oqVthylmOgoYfctJqHskv7qkmd2dn10yTjJhKVsfI7N+CTSuZjlpi0mSqmYKE5WlCfkmplAQhQv7iLPmxNJ7IQQQgghhBBCiBLRi6gQQgghhBBCiFLRi6gQQgghhBBCiFLRi6gQQgghhBBCiFIpVVYU4CfGF8UXANDr+Um+/TaZ5NwmE/73vIyFzcK2zAtfBm0id/GLUaEIndxOXvOZYKIRvHQotVE5RRIhwwCAGpFasMnQsRO1GZOJXIh0iizJtzf/IpcAL7ooCroAoEekFr1enMDIMp/70PPZtYyIiTqkr/txwwD7rhWbyM4m5DfCCll3dL8Jm/BP5Bc5vJgoVgwwiZBl2rmft+zGEII/FiqqIPlud8k56viaEYFRINKMpE2WI/2ckO1lpMYEGVRQQgJTC37czUh2XduYXImKL+KEFjzfvkglK7yJUevGMonQokyKSWDPKW0iY8mZuIsYX4pj3/Vg4kbWhzmRGmVU+OLXpTIScpZzJswq5J5LjuJq7NmIXh9UTMQyHvecRsVdvhQtJlqQiLt7F8s4u0fR8Z6ckzZZLieSRnYtsLG4TsbsnI7trkSvGZpBIj10bYvMKauxTMaKiahUi7SPnUV+jZPaBNfCJOgTUSGEEEIIIYQQpaIXUSGEEEIIIYQQpaIXUSGEEEIIIYQQpaIXUSGEEEIIIYQQpVKqrAjwogsml2Byl0GfLNcmk4EzIgAZENOF0yYB/aZft9fyXdQhMiU2uXqFzJpeSf2E/43+MVdrhvbI1wPS3oycvtiJ1LHSiEFgQpXpio6qRozooksm7feJrGjQJblv+lroxfUsu2Y6bSKWoVIBknFiv1hJ/IT/zuC4q7XQHPmaSauorIJIYLiYKE4gMO8wIdIs8fkmki0icmGCri4RGBmTrJB7AoMJv7pkvGZ92iDflqX57pDbZn/JlTph1HTH8s1gYiKab7Iul2aQ5SKlXWw5xuJdVS9MMR9MMkLHdZLTLrkWmMyL0Sfrshq7nzAh1xLZ7RrJeLPnn1OMHVsgcsgI+DNJrKAlLuNcNMbaQmqR10K84GvOBnH4cZy1kNX6pNiLvBaY4As9P56y8YRdbzlZcDnzxZXMB78zINJDci/rk2fvIrFiIjaOxwq56HKkOzMm7oqUqcbmftroE1EhhBBCCCGEEKWiF1EhhBBCCCGEEKUy0Yuomb3ZzL5gZo+Z2Tum1Sgh5gVlXFQdZVxUGeVbVB1lXCwyY7+ImlkK4FcAfBuA+wB8n5ndN62GCTFrlHFRdZRxUWWUb1F1lHGx6EwiK3otgMdCCI8DgJn9DoDvBPD5m9kIFRPRCfpEFkBELtb2M6mZyIVJYFpbfoL+tS0vodgjk6vrRAxwpuFK6Az8PmpNf7yb3dHlWqHrN0ZgoovU4iZIc1nFdGcqL5jUYioZLwq6AC6c6BGpSr9LZBAkawmRebGM7+/4CfpXdpddjWW8kfr9niUZ7xMJQKPpj22zIMRowmc8EFGAsYzTWmzGfY3llMphypjIf7T7uOmMF/PMxD+9SGkLExglZLymUjsiaGm2/fi61fIh7ZC2rGQ+a3csMzEduU8Q4ddOd/QaavX9tdyPFhgRUUfijyEnYeFSCiLDILVJMj/lW8e4TOc5hZwmJk/hAiN/AtqRsiImNWLb2+n6TLLr8ljui+eWiVQl8eP/dse3Za83mukOEUNy4aGH5Y9JW+Iz7mvRAqPY+8SM5C4FppJxBnt2GRCZHLs+umRdJkFkOQ1k5GmRdRnH/KWAPrF01ciJ2u363LcGoxnvkzzHiujYc3dOns9Z7lOW+8gxO3a52OhOIu6i25tg3TsAPH3o62eGNSGqgjIuqo4yLqqM8i2qjjIuFpojlxWZ2dvN7CEze2ir27zxCkIsGKMZ3591c4SYKiP57infonoczvi2Mi4qyOGM7ynjYo6Y5EX0WQB3Hfr6zmFthBDCu0IID4QQHljP/Y+5CjHHjJFx/2OuQswxN8z4SL4z5VssFDc9hh9TxsVicdMZX1HGxRwxyYvoxwC8zMxebGY1AN8L4H3TaZYQc4EyLqqOMi6qjPItqo4yLhYaC5ETyenKZm8B8L8DSAG8O4TwczdY/hKApwCcAnB57B3PDzqO+eHwMdwdQjg9jY0q45U4jiocAzAHGT+U72J7FpkqHEcVjgH4y+OYhzH8cHsWmSocA1C945iHjFetTxedKhzHTT+nTPQiOi5m9lAI4YHSdzxldBzzw7wdw7y1Z1yqcBxVOAZg/o5j3tozLlU4jiocAzB/xzFv7RmHKhwDoOM4CuapLZOg45gfxjmGI5cVCSGEEEIIIYQQh9GLqBBCCCGEEEKIUpnVi+i7ZrTfaaPjmB/m7RjmrT3jUoXjqMIxAPN3HPPWnnGpwnFU4RiA+TuOeWvPOFThGAAdx1EwT22ZBB3H/HDTxzCTOaJCCCGEEEIIIW5d9KO5QgghhBBCCCFKpfQXUTN7s5l9wcweM7N3lL3/cTGzd5vZRTP77KHaCTN7v5k9Ovz7+CzbeCPM7C4z+6CZfd7MPmdmPzKsL9pxNMzso2b2qeFx/PSw/mIz+8gwW787/J1as2ifMj4jlPFS2qZ8z5AqZHye8z1shzI+I6qQb0AZPyqU8flhWhkv9UXUzFIAvwLg2wDcB+D7zOy+MtswAQ8CeHOh9g4AHwghvAzAB4ZfzzM9AD8eQrgPwOsA/K1h/y/acbQBvCmE8GoA9wN4s5m9DsAvAHhnCOFeANcAvK3shinjM0cZP0KU77mgChmfy3wDyvgcUIV8A8r4UfEglPF5YSoZL/sT0dcCeCyE8HgIoQPgdwB8Z8ltGIsQwocAXC2UvxPAe4b/fg+A7yqzTTdLCOG5EMInhv/eAfAIgDuweMcRQgi7wy/z4Z8A4E0A/vmwPqvjUMZniDJ+5CjfM6YKGZ/jfAPK+EypQr4BZfyoUMbnh2llvOwX0TsAPH3o62eGtUXlbAjhueG/zwM4O8vG3Axmdg+A1wD4CBbwOMwsNbOHAVwE8H4AXwKwGULoDReZVbaU8TlBGT8SlO85YpEzPqf5BpTxuWGR8w0o4yWycNl4HmVcsqKpEQ70wwuhIDazVQC/D+BHQwjbh/9vUY4jhNAPIdwP4E4cfHfvFbNtUfVZlGwAyri4eRYlF8+z6BlXvstnEXLxPIueb0AZnwWLkg1AGX+esl9EnwVw16Gv7xzWFpULZnYOAIZ/X5xxe26ImeU4CP5vhhDeOywv3HE8TwhhE8AHAbwewIaZZcP/mlW2lPEZo4wfKcr3HFCljM9ZvgFlfOZUKd+AMl4CC5cNZfwvKftF9GMAXjY0KtUAfC+A95XchmnyPgBvHf77rQD+YIZtuSFmZgB+HcAjIYRfOvRfi3Ycp81sY/jvJQDfgoOfsf8ggO8ZLjar41DGZ4gyfuQo3zOmChmf43wDyvhMqUK+AWW8ZBYtG8r4YUIIpf4B8BYAX8TBzxH/nbL3P0G7fxvAcwC6OPiZ57cBOIkDs9WjAP5fACdm3c4bHMM34uCj/k8DeHj45y0LeBxfC+CTw+P4LID/ZVh/CYCPAngMwP8FoD6j9injszsGZfzo26Z8z/Y4Fj7j85zvYTuU8dkdw8Lne3gcyvjRtFsZn5M/08q4DVcSQgghhBBCCCFKQbIiIYQQQgghhBClohdRIYQQQgghhBClohdRIYQQQgghhBClohdRIYQQQgghhBClohdRIYQQQgghhBClohdRIYQQQgghhBClohdRIYQQQgghhBClohdRIYQQQgghhBCl8v8DgNTGVKNgnJwAAAAASUVORK5CYII=",
"text/plain": [
""
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAA54AAADrCAYAAAABiEkIAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Il7ecAAAACXBIWXMAAAsTAAALEwEAmpwYAABWnElEQVR4nO29e7At6Vne97zd67Lv5zZnzoxmBCOQsCycMKqaCFGAg3FwCaVcxDZxWcSU4igIl00CFVXMxbGBBKq42Mh2YVMWkTxywjUBg0wBiaKIyBAiGGAQuoDQZQbN9cy57LNv69r95Y+9xuz+3mef9e112WutPs+v6tQ5+zt9+br76a+791rvry2EACGEEEIIIYQQYl5ki+6AEEIIIYQQQoh6owdPIYQQQgghhBBzRQ+eQgghhBBCCCHmih48hRBCCCGEEELMFT14CiGEEEIIIYSYK3rwFEIIIYQQQggxV/TgKYQQQgghhBBirujBcwUws18zs//6Lv//iJkFM2tMuPxvMLOnzezQzH7BzC5P3lshzs48M25mX2VmpZkdnPjz1ul6LMTdmXOmHzSz95nZc6NlPBL9f9vM3mNme2b2gpn9dxNsghB3ZcEZf9zM+tG4nk+wGUJQ5pzv/9TMft3Mdkdj9P9sZtsn/r+2Y7gePO9xzOyLAfxLAN8I4BqAIwD/YqGdEmL2PBdC2Drx572L7pAQU1AC+FUAf+2U//8eAK8B8PkA/gKAv2dmbzqfrgkxE8ZlHAB+KBrXi3PqmxDTcgHA9wF4BYA/C+AhAD984v+/BzUdw/XgmcDoNxqvPvHz42b2fWPm+aiZ/eUTPzfN7IaZvZ5Me8nMfsnMXjKz26N/Pzz6v+8H8JUAfnT0G70fJav70Ojv3dE0X3aGzfsvAPzbEMKHQggHAP4BgL968jcvov7UPOPiHqTOmQ4hvBhC+BcAfvuUSd4K4H8KIdwOIXwCwI8D+C9Tly9Wg3s846Lm1DzfPxlC+NUQwlEI4TaOx+gvPzFJbcdwPXjOj38N4G+e+PnNAJ4PIfwemTYD8K9w/JuNzwPQAfCjABBC+PsA/h2Abxn9Ru9byPx/fvT3xdE0v2lmXzH6CP+0P18xmueLAfz+ywsKIXwaQB/AF0285eJeYVUyDgD3m9mLZvZZM3unmW1Ot+mipqxSpilmdgnAgzgxro/+/cXj5hX3BCuf8RP8HTO7ZWa/Y2Z3+2RU3Dusar7/PICPAfUfwyeqCRRJ/K8A/oGZ7YQQ9nD8Vdb/hU0YQrgJ4Ode/nn0m5YPTrPyEMKvA7iYMOkWgDtR2x0A+sRTjGNVMv6HAB4d/f35AN4L4EcAfPM06xe1ZFUyfTe2Rn+fHNc1pouXqUPGAeCfAXgHjrP9lwD8jJm9EEL4jRksW6wuK5dvM/saHH/C+aWjplqP4frEc06EEJ4D8BsA/pqZXQTwtQB+gk1rZhtm9i/tWPCzh+OP7y/a+RTKHwDYidp2AOyfw7rFCrMqGQ8hvBBC+HgIoQwhfBbA38Pd64bEPcqqZHoMB6O/T47rGtMFgNpkHCGE3w0h3AwhDEMIv4zjbfiri+6XWCyrlm8zeyOAnwTw9SGET46aaz2G68EzjSMAGyd+fiBxvvfi+CP//xzAb4YQnj1luncA+DMAvjSEsIM//fjeRn+HMetx/29mX2lV21v85ytHk34MwJecmO8LALQBfDJepqg1dc44W5bGvvpzL2X6Txd6XC/0PE6M66N/f2zcvGLluCczfpd12dipxCpR63zbcd3p+wD8VyGED/z7hdZ8DNdXbdN4EsA3mNnHAHwNgP8YwBMJ8/0Cjg2x1wD80F2m28bxd8t37fhVJt8d/f+LAL7gLvO/hGMD3Bdg9MAYQvh3+NOP6+/GTwD4zdHJ8LsA/kcAPx9CqMVvVkQyT6KmGTezvwDgMwD+BMDDAH4AwC+Om0+sPE+ippkGADNbA/Dyb+bbZrYWQuiOfv7XAP4HM3titB3fBOBvpSxXrBRP4h7NuJl9PY6tt0cA/hMcP2j8Zbogsao8iZrm28z+HI7z+9+EEP4tmaS2Y7h+65/Gt+J4QNvFsQX2F1JmCiF0cPz98VcB+Pm7TPpPAKwDuAHg/8NxGE/yTwF8vR1bt/4ZWc8RgO8H8Buj4uU3pvRvNO/HAPxtHD+AXsfxifh3UucXtaG2GQfwegD/L4DD0d9/AOC/PcP8YjWpc6aB4xuml7+S9Yejn1/muwF8GsDTAP4fAD8cQoj7J1afeznj3wrgWRxv+w8D+KYQwq+dcfliualzvt8B4CqAd5/4NPTkJ5q1HcMthHGfJItpMLN/COCLQgh/c+zEQqwgyrioG8q0qDvKuKgzyvfyoq/azpHRR/dvw7FVS4jaoYyLuqFMi7qjjIs6o3wvN/qq7RSY2XedUjz8K2b2TQA+B+BXQggfGrcsIZYRZVzUDWVa1B1lXNQZ5Xu10VdthRBCCCGEEELMFX3iKYQQQgghhBBirujBUwghhBBCCCHEXDlXudBGvhkuNC6e5yrFijDLL3y/2H/uRgjh6gwXmYwyvhqseoHBojKufIvz4IUFj+E7UcaNT7pwbAEdW0R11qqP14xF3qesk4zPktRYxtOxPE+6rDPNHJGa8ZTJ2LLYfHXM+PVTMn6uD54XGhfx1oe++TxXKVaEWZ50P/TZ7356hos7E8r4arDqg/yiMn6hcRF/6+HZ5HtZb+bPg0Xkr1yh0P/gAsfwHTKGzzKrWeLCUiZLXdYsSc1RymSzXNZZlrcM/KOnFpvxtzwwm3GcZZB9lTInjXk0bzNhGgBosHUmtqUwLH0bi9aATFeEs09zWhvL87wjPssh5Z/+Cc+4vmorhBBCCCGEEGKu6MFTCCGEEEIIIcRc0YOnEEIIIYQQQoi5cq41nmJ1WaHSCXEPo5yeP4uo1VxEbVsKqTVmrPuTZnfSurbzOFdWqebubqRsBjumLKep08VtqZGf97mRWnfGpmN1bJMui0lbSDld0rJmCVv8sg1XZQg4GlZ7OmlNMavBzMmERnZMXKvJajc3Gn7GFvnIrJ356VjfYlgeBsHPyOo+u6WfrleMXz5b1oBMSGtBo3n5+ZIW8izRTjbr/OoTTyGEEEIIIYQQc0UPnkIIIYQQQggh5ooePIUQQgghhBBCzBU9eAohhBBCCCGEmCsrKxdaVrnENCyLiGFJuuFYlv2zbKzSubDqx3DV+z8L5hm3SbO8SucAwHMUb0KqZIWRMtk0WU50V6wsKfsmztw0uyRFQsREKdPIixgp212QttRtn1iElSgSipefKl6ZJs8pQqNlowzAETPcRHjBlQ9Sk3x81aTjm583lgSxZa3nvm2n4fu+2fDJzJnRKCIQkVCPSIOOCt+5fOjbymh5fbKb2Z5nh6NH7ELDKKyp4i1O2oQpDqKzXIP1iacQQgghhBBCiLmiB08hhBBCCCGEEHNl7IOnma2Z2W+Z2e+b2cfM7HtH7Y+b2WfN7MnRn0fn3lsh5oAyLuqM8i3qjjIu6o4yLupCSo1nD8BXhxAOzKwJ4NfN7FdG//ffhxD+92k6sGq1OfNklnUYqcy7VGdFauLmmvEU7pXzgG3nsmRk7i8WX9x2LjzfjNTMp0y3LKdPar5neR5MWl+XmsdZ1rDN8RybWcYD/Ivb2W/o42lyMhHbXlaryfIbT8fma5B1srq21PMjnjOuVwOAjL3UntVgJuQ+dT6WQTZvEc3M5qP1otPUO08+61mZWcbLEHA4rNZEZgkpyUmx35AMZgWpy2T5DdE613K/Ny83fe3mfe2+a9tsDP0KCPE6B6Ses0s2oJn5Nj4uVFs7Bdk/JDSsnpO1DaKwDoJP+TT3Gkn1nFNeccc+eIYQAoCD0Y/N0Z8luVUUYnqUcVFnlG9Rd5RxUXeUcVEXkmo8zSw3sycBXAfw/hDCh0f/9f1m9hEze6eZtefVSSHmjTIu6ozyLeqOMi7qjjIu6kDSg2cIoQghPArgYQBvMLM/B+A7AbwWwH8E4DKAb2fzmtnbzewJM3viqDicTa+FmDHKuKgzM8t3qXyL5WRWGe9oDBdLyqwy3gtH59VlIRxnstqGEHYBfBDAm0IIz4djegD+FYA3nDLPu0IIj4UQHtvIN6fusBDzRBkXdWbqfGfKt1hups34usZwseRMm/G2bZxjb4WoMrbG08yuAhiEEHbNbB3A1wD4QTN7MITwvJkZgP8MwEdTVnivSFTqziKEMPOSs8w64yksw3kw6y5MengWIRya5YvMl51Z5zveBbPMEctCyvJTz6dp+ppy6OctEkol5aXiqdKg1L4u8tSYdcZjpwfbV/FxZe+qZ0IVlgfWFs/bJB8TNInph67TN1GJSJyRggx4A3agiaCFNFGpyrg+nDZfLBJi03F5kW+81zJeIOCwGIxdZyyRyc0nqVn6tgFpY3KeregJ5OF1Lwi62u65NiYS2mz57ckzf+bGwqze0D8GHQ6aro3RK/x2xoKkFjm5WQaZSOio8GKlflndpkHw05RTJDVFHDR3uRCABwG818xyHI9fPxtC+CUz+79HJ4IBeBLA356qJ0IsDmVc1BnlW9QdZVzUHWVc1IIUq+1HALyetH/1XHokxDmjjIs6o3yLuqOMi7qjjIu6cKYaTyGEEEIIIYQQ4qzowVMIIYQQQgghxFxJqfGsHUvgdZmqSH1SWcU067yXZCyzZt4ioWXIM5DWj9Q4LELIEjNNdlPFLXWA7ab48M36HHBil4RpZs00IqFZ5jtFJAT4TNL5kteZNuUqjv8hBCeuoSKeFAnHFMe5EX0swERCLSYcIpajBjMfEcpom5icx5g1iCSnCH66eH+wXg0TRUIDMsjGMiQuF/JtgfSkzvc8JQIO0B07XSyRseAD1yBtm6Hl2l7R8I8bV9vVg7hGZEAbDS/PubLlXwezvkbkQrlfXhEJgTpdLxIycr4MSO7buRcm5dHy6TlKBoYBCRwTQHVQbSvA5EKT34Fk5PPIFJmQneFzTH3iKYQQQgghhBBirujBUwghhBBCCCHEXNGDpxBCCCGEEEKIuaIHTyGEEEIIIYQQc6X2cqFlEa/EpPZrEXXqy1JUfy8JWlKZZZ5nLV9JyU2almK2pOY5Jb/TZPK85UjLRKpgJ1USlCIvmrdwiC2KyUzYOlOmY9Mky39YWzQzWxaTBqWIilJZBfFKgBfX0HzFAhIyEXGloOF9JDRLcVuTTLROVrCW+52cE1kK26ZhJFAZEkEQHcPJdIOEcy1VcMXOBSZjiaeLZUNAep4DmbAuY3iJAge2d+b5mHgmh5fzbOKya7vS9suLhVlMJHRxvePbLvm25pafNyPnWhkFs3HHH/2y9Nt52PfbyT65izPIIsOEQwwmgDqyg8rPQxu6aQJJdKr8JyOyqPi4G/yOZdk4dR3JUwohhBBCCCGEEBOgB08hhBBCCCGEEHNFD55CCCGEEEIIIeaKHjyFEEIIIYQQQsyVWsmFJvVGzFqyksKqF6mnSiIkCFoMi8h0yjpZ7hchHJqlSGiac3kVh4GU7Y2zME0cUyRBOZkoXTjEJCLje0zzTkLDxCizPD2pLCVB2lKQibhkJa0fq35Ne5kQgH60Mex4NVwA/A4wEpJG4n6KhUCtzM/IREKbRNCyxixHhCLKfZdIVkDEIgNyvrBVOsFVojSoR06iIZk5lkLxjKcKh9IO1CoIs2IKFDi03UpbCZ+bmIwc+xbWXdvF5jXX1iain0YUm/Xci3Iu7niR0Pr9vq/5Bb8CIxaf8jCaNxu4aXpd/2jUIOcfIxZyxZkE+HV/Ld4ZAHpFz7Xt263KzwV8/xlMLsSOp9n46fiyJBcSQgghhBBCCLEk6MFTCCGEEEIIIcRcGfvgaWZrZvZbZvb7ZvYxM/veUfurzOzDZvYpM/sZM2vNv7tCzB5lXNQdZVzUGeVb1B1lXNSFlE88ewC+OoTwJQAeBfAmM3sjgB8E8M4QwqsB3Abwtrn1Uoj5ooyLuqOMizqjfIu6o4yLWjBWLhRCCAAORj82R38CgK8G8A2j9vcC+B4APzb7LnKmkTPMUrzCljWJfOO0+VLFK/HyUkUPqdMti4xlHsw742x7553Bea4vldTtXpY8pOQ3+XyZtA9zMlXMKuMBadsW76dZ5y+WCTG5UCxnOUs/2LxOjMJmpAH3TWxeJiEa14dTFs+XH83M5isS+8qYV3ZTmOUYXgLoRgeD/YY+ltu0EoRUANAkGWHHPp6sTeQmG7k/YBeaXjay00oTkPSL6pYeDptuGpbBXum3qUN3RyReIcsakgzGsqfj6ZhcqNpWkEwOQ9qdysRjOJkzm4FObLYZH+KovD1+nZFwyIiMJs98Ri6106r4LNpXa0SMtXHFZ7fxYNu1ZZfW/Apa/hHH7lRlRc2+lxc1bvp+xH0FgGFCfvskbj3icWolXpzi41YELyBi8h8mDaLH0/zxjJeXKiU6jaQpzSw3sycBXAfwfgCfBrAbQnhZQfUMgIeS1yrEkqGMi7qjjIs6o3yLuqOMizqQ9OAZQihCCI8CeBjAGwC8NnUFZvZ2M3vCzJ44Kg8n66UQc2ZmGS+UcbGcTJrxk/nuKN9iSZnVGN4rj+bVRSGmYlYZL8u0T8CFmAdnstqGEHYBfBDAlwG4aGYvf479MIBnT5nnXSGEx0IIj21km9P0VYi5M3XGc2VcLDdnzfjJfK8r32LJmXYMb2cb59NRISZk2oxn5OuxQpwXY2s8zewqgEEIYdfM1gF8DY6LmT8I4OsB/DSAtwL4xZQVTlr/OCmLqH+b5fKmqfuclEXUcy6yzG/WGU8hpQZu3nXM8y77nOY8SKkPnWV9cmof5l3vNq8yuVlmnNUBxuTRrzTZvmR1mQyWo7iN1WSS93HTdVpiP+Jjw2sySW0oWQHprp+P5Y+0sX7E9ZxsOrZ8VhM3yzF8XvXus8x3GQI6w6i2jfSxWVZDXpBf47PtZTWeawnnx3bTF4bdv9Z1bRfXfM3XGqnxLEvf4e6gekvIjs2Q1HP2ybL2h74OLM5gXJMJAF0ywLC2funb4vwWpJ6T1WCytmVjlhkPoUCv2Bu7zhDd6bHawQfsi1xbmwSHnebxZK3cZ7x1zecou7bll3//BbICclJGF6dst+8mKQvf2wHJOMt9XPc5zTjbDP4XBP3yoPLzsPA1qgz2y4bM/CNgSo0nrRcN/jidxtgHTwAPAnivmeU4/oT0Z0MIv2RmHwfw02b2fQB+D8C7k9cqxHKhjIu6o4yLOqN8i7qjjItakGK1/QiA15P2z+D4O+ZCrDTKuKg7yrioM8q3qDvKuKgLZ6rxFEIIIYQQQgghzooePIUQQgghhBBCzJWUGs/asQziFVbKnipZSJU9zJK4RP88pEGLfCH5pAT4bZ53llJIKew/y7yMeLsXktPEdTI/TjwvW1RqJicWGk0227kRwF+a7Yg2JJYNAenCIZa/WCbUJMtvEoMPW36DTFeStcb9HQSyMHIAiYsFxF3hoPljbSRsVDgUtQ3JAeCZZ31LC3jK+cglTYujRMBR9LqJjOShHwk2YtnQ8bK8cKNBTqBLLb/8a+1h5ecHiEhop+3FKFvrXi7UipYFAIEEs9mvyl3YONYjFqVO4bdzg5z0eWRp6pODfzT0Ke8UXjrTC36biugMWVWR0LwJocRggldjMRnNxbDj2tbJk0UzY+NxLIPymcwuENkNEQmFhx7wKx14qZYdVmU8oUfy1l13bV2ScUZ8X8XleGmK0DZarq0oq+d3f7jvpmHyn6wkcqHMH6iCHOMsEg7R5Z/hc0x94imEEEIIIYQQYq7owVMIIYQQQgghxFzRg6cQQgghhBBCiLmiB08hhBBCCCGEEHPlXOVCqeKVWEiQKkGZlGnEK5OyCPFKKpMKTlI3ad6ClmWDbUYcL5YHlkEmR2HTxU0p05w23aSklc+n52aWEhIuTImnIfKVKeRFjFXLeAiBSmncdHGQyA4hfoLkC1KcLSYSahOZRYu08cwzYU91wgbZDx3yu9whOciWoOxiuzlFGnRaW3zc+Hzzz/yyU6LEHo4qbUyckYeqbKQVvHxkOPBCj+1m27U9sOb33oVmVZ7Tzr1gZ5uIhHaueAlRgyy/HPgMNvar0w2ZSGjotymWxAD8nFyPdhGTrDCR0EHw29mDFysV5oVDs8TCZJ/VnEW8ch4ElBgWkWQn+IzEEpl200t9tjOf500ykG/m48fo3pBZiYjU575Lrincf9W12aEXKIVPPFX5+eAz/tjcOvRyoUDG7Fbm99lmJNXqkpu2I3btIxeiJpGTFWU198PiyE3DLq5MDJUHLy+KRUIAkGXV88pInrPMz3cay3U2CCGEEEIIIYSoHXrwFEIIIYQQQggxV/TgKYQQQgghhBBirujBUwghhBBCCCHEXDlXuRBjluKVVFnKIsQr8TYwScy8xStUqJIqjoimY7PNWhq0qrKKFHFUnKVZu6xSlp8qKkolRQrG9s2kq0yWWSXOG+eX5XTW0qBllowxAoB+tBOoIC7ewyQMxD2SvD8a0eKa5Feoa0RmsZb7I8gkKIxBJBfqEvGK224Aw3KyhLMMFWQHDUkoB2S6WBxUkPnouUI6knqc2P6IyWY++k1HgQL72W6ljco0ItFME16yUoYt13apvebaNkguY2HPVtvLdLYvdVzb2jUi0FrzkpKyS2QykUSk3x+4aRqdtFFwGPxxjYeBNvHGDAKRC5mXxPTMb3uY8K4h5fgCQMaMaG5ZfqOWTi4UAoall1DFxDIlJp652PKPEdsNn8HtBhl7IznPoCT7iQxU4ZKXHIULvg0Dn9/yU7cqP9+8temm6RYkmIR1ct5uRdvZJeN/t/Bth8SL1TRy3pbVCcvgxwUwCZalnRuBTBeLp2LpFHC2c2+5zgYhhBBCCCGEELVDD55CCCGEEEIIIebK2AdPM3ulmX3QzD5uZh8zs28dtX+PmT1rZk+O/rx5/t0VYvYo46LOKN+i7ijjou4o46IupNR4DgG8I4Twu2a2DeB3zOz9o/97ZwjhH82ve0KcC8q4qDPKt6g7yrioO8q4qAVjHzxDCM8DeH70730z+wSAhyZd4aqIV2YtHIon43KTtHXOElYOTCVEbho/0TQyFkaqoGVaZpnxAC96Yl8riPfLNFIfNmvcxkRCDVZ/TpaVEflKScQR8TqY8IotP1WqlSLQohmkbePzy0UracunfUubbObMMt9lALqR7MHIUY2dC4G4Gth8TBLEiLPWzvze3Wx4SckWaVvLfRs7pv2yuhE5Ez+Q+XoZkY2QEyHO94Asa0BC2SMnxjBBLjRMHMMZKdKgVIoZLGumGccQB+Hm2Okyq94+5fDilR0iF7riHUSUOJeba14isnbFB6Lx0IZrs5yIcjrEZlL2Kj+27vhzI8/SrupsvI7z2yQnQoucVwe269p6OBjbByoNIvIfduyYQIVJgtjyUvpxVmZ7Lx5QlkRKE+Eybv6R4RIxRG2QMXWDjL2N6N6CXktJY9jwQiA0W67Jbt1ybb1nqyHsD/02FeT+hvUtFoABQCs6PzbIudeJ7XgA1shNWhGI6AfV/RjYcSTZpWctkxAlQPN8hpv9M63VzB4B8HoAHx41fYuZfcTM3mNml86yLCGWEWVc1BnlW9QdZVzUHWVcrDLJD55mtgXg5wB8WwhhD8CPAfhCAI/i+Lcw//iU+d5uZk+Y2RNHhVdiC7EszCLjHWVcLCmzyHevPDqv7gpxZmaR8SLhkyAhFsUsMr66L6sTdSDpwdPMmjgO+k+EEH4eAEIIL4YQinD8gpcfB/AGNm8I4V0hhMdCCI9t5OTjcSGWgFllfF0ZF0vIrPLdzvzX+IRYBmaV8TzzX9kTYhmYVcb1QguxSMbWeJqZAXg3gE+EEH7kRPuDo++cA8BfAfDRccuaZf0b/771uB4ck1L/xtaZWvfJiOedpi4s5WvwdL4p1hnXxKXWcybXDC2qAA6zzTjgM872S8qwn1K7CaTll9VzsvoEuizaOVL3Gc9Ha3xIPSBZ/KT1zjSDrC0hv+S91cnVaKyGlE+XuMApmGW+ywAcDas7JqO1mtW21M1kNV9s3ji7a/RF3r6m6GKr59q213xbQV5k3hlU68AM/gGF1QYdFewMYvU81Z8HJBxdUjjXJ9MNSr8/4hrP+OfTSM3ypGQ2vclgthkf4qi4kbDOuObXF29ezP+saxuScYXdWzSjWrG1jYGf5j5fm5g9cMEvjF0A9ruuKe9Wz5m87c8hBss9b7v7zwCwlvsrwrD052inuD22X5n5/ZOzNlLjGdc3njpddAXjdaXTP+jN9j4lIARS4xvh78X99l8iv6dZz/2BbZPa4NgfwTJjbX8cwsa6XynBnnrOtfUPovM2sWZ5QK4JzH8Rl282aR2ob2PXvh78MYq/kRHIFdJIbSiM3dCwGtK0tkmmeZkUq+2XA/hGAH9gZk+O2r4LwFvM7FEc3xc8BeCbk9cqxHKhjIs6o3yLuqOMi7qjjItakGK1/XXwDxt+efbdEeL8UcZFnVG+Rd1RxkXdUcZFXdAXvYUQQgghhBBCzBU9eAohhBBCCCGEmCspNZ4zJfYgTCpemVQkBPhCflbYz5bPZCxsXka83UWilmje8hG2fCaTiJuYFCB1+Qv0CM2dEAKKaKOZOyMkHP/UbDG8QIsUs5PlMwdFlnjEYjEAPY9JUX1ZpG1oPCfPLmvzjSy/KYIutqx7KeMlAg6LqvCAyoUiEQMTM5Tk5dUsf1RIFTVdbnnxyrUN/+qX7XUvKWm1iMBh6DvSjARGTEoRvxAdAI7Isl4i296LPC6doV9Wh5hp+kQkNCDSiDJKZZizNCgVS1b0nQ9lKNAd3jnzfI3cC0921tjL6dOWF8uFWpte9JNd3fYzPnDZt2U+b9bc9223q+eMmV9nIAIYdh+XInpj4+lanvZy+l6x59pi4RMTBDWIBCpkvi0LXqQTiKAllgnFsiE2zeIJQIJcCFbtd4Psp52mP4brRPbGzvJ+dF3YIEI4XCAioaY3GlnPy7LCM7d8W9TdFlnnWu7bjoY+S+y6FpMqcWT3EUfmr2FFWd1OOn7asuWtynL3TgghhBBCCCHEyqMHTyGEEEIIIYQQc0UPnkIIIYQQQggh5ooePIUQQgghhBBCzJVzlQuFEDCckXiFMalwiM3XJBIUJmNh/WfEBc0DJnYgEg3SlCQnoJIVMl2RKF6JJ0sVqjB5QCrzFivNgwBgEO1oJgmyuOCcTETcJci9t4AuP5a0sIy3WMbJr6JYcTwjzk28HwDQQBdkncSr4mDRSs0zlwuF6GeyTtKP1JyGGiiHCpTYD51KGxNnNCJ5TrtkkhUv78jNB/z+NZ+Zz9+oSoJeueNFI+2mF0S017yEqL3hJRvlgMhYovOlYFlu+fmuFn6bdge+7bnII3FITqD9wvd/EPx2FmS0j+VCDCaKmma6FHFQRmQ1iySgRH9YFe8EImuySODRzDfdNGw8TSWP8pav++NnRC4U7r8vafk2JCKXiKLvN6BP8szGQLbp8XUoJzdQrI2NMcPCy2TcsjIvoSHuH3pzlBGRUCAT5qiOY0xAFE+zDPjrEdkJIRbJ+e3YbrD7CL+sDsnNIDr3P6914Kax7S3fLd9ToNNxTeVtL5OL49UgIqT1pr8mtId+2/cH/roW37ukPs90yU3Jbbvu2txYRERCRqRarI3PO37QYtNkWfrjpD7xFEIIIYQQQggxV/TgKYQQQgghhBBirujBUwghhBBCCCHEXNGDpxBCCCGEEEKIuXK+ciHMULxyyvJTiBfHRELEEYEmkaw0EsUrw6iI2qhQgYlR0gqT4zm5eCWxjdWYRz8zaVC6ZKW+hAD0Y4EWma4RhTCwA0aMQDmV4vi2eE6W5zWysDVyLjRIG2MYyVa6RBJRkr2RUwEVmTeajAmIWJ6HJM+x5IzNy2RAsxZepbi3UgVm50GJAnt2Z+x0sUyjSURCvWLdtW0GL2h59ZYXPdy3VhWLBDJOtlp+vq37+q6tcYFkrefnzW5VgzQcelnGXrft5yPXiKttv/zPNav76LMHfpo74ci1DcxvUwoWiCCCiaKIjYXJXlKEQ2yaSWWC8yKEEsMi3s/MlFbdBzmRdzBRDrvnYcTyqoz5abY3XFO4fNlPN/RZshu3XFu5W83S0aHPMxsDG2Sb6PUluqYx+RLbP0zOU5ReHJMCF6OkSVYCyX0ZtbHLZWFeCrZ8MJteNfcZGQs2iJxnUPr9dJsI1S5HY/SFTSKM2rji20pyPjK50KHPfWyXyokIqZV78dYaaesSM2Ivui4wyWK38CF/qeP7ul+84NqySJgVyDjOxmcmEsqYhIiN7ZE4yEgO2LJOQ594CiGEEEIIIYSYK3rwFEIIIYQQQggxV8Y+eJrZK83sg2b2cTP7mJl966j9spm938z+ePT3pfl3V4jZo4yLuqOMizqjfIu6o4yLupDyiecQwDtCCK8D8EYAf9fMXgfgOwB8IITwGgAfGP0sxCqijIu6o4yLOqN8i7qjjItaMLYaNITwPIDnR//eN7NPAHgIwNcB+KrRZO8F8GsAvv3uywIGCXaOWLxCi55JBToTrzCcXIgUs6+zgmNSNc7ERIxBJArol0SyQIqEB6RvrPg+SbzChCqk8LlIEAexrU6XCy2XXmiWGS8BdCNLDfvtzjDaxw0ioWBSmSZpDL7O2wlN2iSnTAqw1fAF9C1yLpRE5hILtBqF7xiTZbFCe0Ys/2GCICYSYmMOy3h8LjCdyEKYwekyq4wXGOAOXqy0MQlHLCjIzctB+uGia3v92o5ru9j0woluUb103dfy0p3t+718pPUKf8nLdnzfwhGTUlQFIRtdL/VZO/Iylm7HS5SYlOKB9eqB/v3b/sDfybwQZgC/nUwQ0YgELTm5/FOJCxFJ5IHsxwThUCxiYdNMwizHcKBEKCN5FTkJLep3RjLO5DmszcjyY0FL8EMzsOmzFba3/fJ7RNpC8ju4UV3JkOSUjf3sGsGuL4fRzUssGwKAJrnB6cOf3yGMH6HZNKytJKN9ljidWRlN48lnUNU224z7/LLLTAALXRV2P3qrT25KCNuN6ji7vk4kTE1i1SqILOvw0LWVxLsWX/qN3BtR4VDG7o3YdlbbeuRe/ybxYj0z2HNtg9ILk1p59fxOlWzR6zS9dhOZXDQdmy/P/LXvNM50NpjZIwBeD+DDAK6NTgQAeAHAtbMsS4hlRBkXdUcZF3VG+RZ1RxkXq0zyg6eZbQH4OQDfFkKoPJqH43dB0N/Lm9nbzewJM3uiW/rfWgmxLMwi4z1lXCwxk2T8ZL4L9itkIZaEWYzh9X7hl1h1lHGx6iQ9eJpZE8dB/4kQws+Pml80swdH//8ggOts3hDCu0IIj4UQHlvL/DunhFgGZpXxtjIulpRJM34y33n0DjEhloVZjeH8zctCLB5lXNSBFKutAXg3gE+EEH7kxH+9D8BbR/9+K4BfnH33hJg/yrioO8q4qDPKt6g7yrioC2PlQgC+HMA3AvgDM3ty1PZdAH4AwM+a2dsAPA3gr49bUADQiwwhTKASS0RSpTU5Wdgaqf2Np2KCoHViKtog4pU2KUJm9CJRwNHQdywWzgBAlyyeFffH+4zJUwbJ4hU/XSxXSD0mK8LsMh4CusX4gvw8KjhvZ2nfeicuBrRJ8XrcwuRC2yTPF1q+uH+9wUQrnm6U6XzgpQCxZAsADsm5wDIY57dP8hyPLwCXEKUItMoZfx1pFhKVKZhJxsswxFFxo9rGxByIZQT+GGe5b3tk0+/z3YG/TF1dq8pStra9YKF5lQgQrvlvJNgF3xYOvYylMTyo/Ny648+L1q4/p5rkGnFApC1xy3buz5+9wn+YMQx+2xvmRQ9xW8v8djfJsQxEJFSaP045sZzFOcjIORVmc17MbAw/vtrFx3F8vwPZdy0yYK8xmQnZBUfD6n4vOsw2R8QrTfKthCEZw/cO/GSH1XXkVKjotzMneWD3VY0o5C1y2WPjdae87SckpEhQWNuKMMOMGyz69koox1/nC/j7gyMyljFh4HbD56YdCajyNrl3apHHlKGfzvZ9nsOAnLfRPUj882nEwkaA5z6+P98f+v3zzKHv/4Htu7ZWvuna4nGmKP35niLeOo0UWSCbppl50dlppFhtfx2nfy7/F5PXJMSSooyLuqOMizqjfIu6o4yLurCyv/oRQgghhBBCCLEa6MFTCCGEEEIIIcRcSanxnBllALrxW9oJca3mkLyhltUYsqfoASmeiL+rzWoudpr+++47Tf/99rV8fD0fAHSLav1Lk7xsehh82/5wfD0nQOrfSLf6ZKexGk9S/uaqVwKbiGCsiHeGLNtvTgoEHBbVnLDavjz6jjyvQyR1YOS79ewl3DGbpJ7zctvXhl3e8C8s3tpIe0HxUadaQ9Y48t/5j2udAeAGedl0j+S3M6zuo6OhH0v65OQY0JeI++VPSmrEJ13lnE+hM1GGAp1Btd6KvqQ9qpFjL6W+uvFqv3yyzoK9sD7Kc3PDBya/suba7IGLfgWXtv10d/zLyPNu9ZrQvOFrilhNXCo3+9XtvEIEBd09X+s2DP4VN6zeJmRb1YbUXLHp2GYmTMfqgRvnexsyN4bkRe6shnGN+CNY3eQgGisHR2RhqbKFwt/PhH1fxxxCXONJ6lbJtaTJBBIEX+PrQ3OH3LyUwfc/y3xusui+ilm442mO+5VWC8quv3H9W0bGuniaRWOWo5FXx4N+6a/9MQPyurg7g7SB5CKtkazmJiNOFgqrWT70/SexQRHVdBbknmRY+I7E5yPg7+sB4DCqeb1N3kDGnoNawV+v1rNLri2+lg7NH5OCbDirQU8lpcazFV9f7sJynQ1CCCGEEEIIIWqHHjyFEEIIIYQQQswVPXgKIYQQQgghhJgrevAUQgghhBBCCDFXzlcuhICjhJfU5lHBeYMU9TLxCnuOZuKV+AWva7kvur3U8hXBV4h4ZXPNTxcvHwA6vWpBe7PjC4nZi3jZNvVJjXAsXjkk4hUmsBkSy8qk4hUqQSELY7IdBnFKOdLUTudHiYADeGlDTBZlpEle0t4tfRsTUDXJjgqRW+QiEWM9eMHLUS7c5zPe3CKBI03rB9V1rO36dbIC/aePvADicOhzsxcJLA6JYICLhGYnwmLZJc4EZGRZKb3g51DCjOdEGQr0hneqjVQuNF5ksFnuuLbdvt8B9/l4OBqbfifZfUR28MB9rils+5d0W4MIQvarEods3QuIGCzzLJKxn4U5w/qFX+ewIOONHya8ICLz/coTxStUxkKOeXweGLluM+HQYjFYtB8CETjFJ2ZB5ELs1GXHNZBxpRGJqvo9Yl4Z+DGWiYTQJYK4I7ZNVTJyb9QiQsU2mS4nA2MZbecB6erzxZ5rY3KedvOinzlhvtz8dTXP2mQ6di6wtupxYeKVnJ2QC8Qsx1qzKq4ZDHfddCFUj3WXTHO967N7xe9Omvv41oVcSoAhuctjue+Re3EyazGMpF1Df171iTSIi4RIWyQEZaJPdh/Rht9pm/ByoTjSAyoqZeKwya7TjIycQ22TXEgIIYQQQgghxJKgB08hhBBCCCGEEHNFD55CCCGEEEIIIeaKHjyFEEIIIYQQQsyVc5YLlTgIpMg9Ii68bRAhQY/IWIakjYlXLrerbdsNX/37wI4Xr1y678i1tXbSxCv9g+o2rN3yVfWHQ18k/JlDb9bYH/ii/TtRBXO3IHIhWlycRrwXmTwlFuYAXJbC1smmIy4kv85EUdF5UaDAXrZ75vkaJLtMSHA43HBtzf62a9tpVPfLax+44ae55s/FxgVyXNd9AX0oiDSqXc1Xo+0L3B/Jbru2pw/XXdvvDvw5f71fPf86IDIBJjhhggn4bYrHmZzM12BSFSYSItktSfLdecRERUuV8RJlWT2uIfixLJZSMIL5Y0ViBeKZQojGmqxNxp4LXhoUrlz0bTsX/AqIeMdejLJLxvmiJDI7Pxkd75rRKneJRY7JIAJRrJWll2+ErDovX5ZvY+dUSY5dnpD5jMwHKgpcJDnyvCrKGA79uBVTlH48YtdqlnF2rWtl1eMaSLbQ9eu0PmnrEfFRz5+3sYwlPs/OAtumONIvHPk8DMz3fz3zkpWMSFXYuRBjZOxnIiEuJho/XSwbAoAmEccsktya2Gpdq7R1+jfddEU01jOR2dMHfp/vNMn1NRt/Q1eS6z76PqfW8RJElGScIrK6waDat+7A33sdsbYha/P97UbnKbt+5eS5ZJ3c73WDv9+Lr5tdI3kmyyrNH6fU60kMk2xtBC8LPH1+IYQQQgghhBBijujBUwghhBBCCCHEXBn74Glm7zGz62b20RNt32Nmz5rZk6M/b55vN4WYH8q4qDvKuKgzyreoO8q4qAspn3g+DuBNpP2dIYRHR39+ebbdEuJceRzKuKg3j0MZF/XlcSjfot48DmVc1ICxcqEQwofM7JFZrKxAiQOrSntYISsr/o5pkGkOhl5SUpDi3FdtV4U9r76066a59oiXCzWu+Of0rE0K3ElVfaNT3c7mxqGbZlj45X/kjt+mGx2/z14aVsUrQ6q08DBxCROoxKIVJl7JybLY8qlcKLFvMaWl6pFOZ5YZLzHEIaoiClrAHYmeyoSCbgBoZT7PD4X/wLU9dmWv8vOFB71cIt/2+ze/ROQJm74NRF5lzarQhElnNgdeHPHFe3dc2wde8DKJ57JnKz/34QUDDSJxaAd/DjHZQzu0o2n88BhIenMmRzEyVpA8l8xC5JY1fpJxzC7jwR1XLhIan+cD88cdeNCvkeyiQTRWEm8asO6PcbjocxU2vYQIRMYSH4bi0HesJDKWBhmjeFt13hukD2XpzykmWcmy8RIUJk9hMDlXNqEQyOYkEprlGJ5ZjrVmNSdHpZeqxGN4Ufrj9dyhlzy9dscfm7zl89DOIxkUGyo6RNZ46O8tQGQsgcirYrlLLGIBgH7h27qkrUPa9gbVjLNr4wa5Z9uyK64tJxkvEF+D0u6DjIzX7Pxg4qAkuVCYXi40y4znaOGSvbLStt963k03GFazlGX+mvjpvhcXfsHwmmsbMvFalOluhzySsIwzudCQ5Lnw6+xF4qDeMC3PR+T+PBYJAUAvFnS5Kfgnfk0iCVoj9y5OUkjuDwbkmhCfG8fLShMOxTB50Vrpz9vTmOYq8C1m9pHRx//+ai7E6qOMi7qjjIs6o3yLuqOMi5Vi0gfPHwPwhQAeBfA8gH982oRm9nYze8LMnhiS3xoKsaRMlPEi+E/0hFhSkjJ+Mt/n2DchpmWiMTzldQJCLAkT3qfoXlwsjokePEMIL4YQinD8PYYfB/CGu0z7rhDCYyGExxrZ2qT9FOJcmTTjufn3rgqxjKRm/GS+z7eHQkzOpGN4SqmPEMvA5PcpuhcXi2OiB08zO1mI81cAfPS0aYVYRZRxUXeUcVFnlG9Rd5RxsYqMlQuZ2U8B+CoA95nZMwC+G8BXmdmjOK6bfQrAN6esrESBQ9t1bTFxcSsriuV9JcXgxWtc26u3qtW4D33hnpumeb//rWd22ReI2zoRrxC5UNiPvoKZ+a9k3l94odGrb/jpfnnoBROfyz5b+ZkV/3Khii9ejiUrx/NWo9IMfv+EBCkRAFoMnREpR4o4iM13Vmad8U5ZlQuVRLJThGqmmZiCzdfI/PG6b/NR1/bK+3fH9BTINknGr/oCcbtE5Ct9IgVrRuKwDpFckG5d2Tpybfeve5nErcOnKz8PSi8YaOc7rm0j82Uva9jyHYmh8Rs7ZJ4+b4JwiKU5SUA0htll3GAW7QMi8AhuS/w27JbPubYyvHZ8FwAMyuq+LDtkH2Vkb7aJ5KNJvqVA5g2H1a+nDY788Wzmfl+s5/5caZKxbS+6zL1gL7lpWo1t18ZoZH4747acCCjYdYPJUqhwiLaNH5/ZfGdllmN4ljWx2bq/0tYf7o+dr5H7sfnZcNO13ex7gdYr1v1+suicKYncBD1/f2BULkS+WklkLP1e9VjHIhYA6JC2AyJoORiOF68YyccG+dZQr/TjekbG04FVr6OF+WtoKpPKt9j50sb0nzDOMuNNNHCtrGb8dsvn8ii7VZ2P3H90zV+/X+r68W1/0+/Po2E1SwdHZHw+9Nm1I5Jnci6wS2dnUB3j4j4AXIzVJedfj8iFioTLdU6i1SR5bgfftyJBVJWT+xR2LgzJs1WKXIg9S6yfIeMpVtu3kOZ3J69BiCVHGRd1RxkXdUb5FnVHGRd1YT5ucyGEEEIIIYQQYoQePIUQQgghhBBCzBU9eAohhBBCCCGEmCuJpozZUKJAJ9yJ2nxxayxeYRKRYeHbmIzlwdarXNvrr1YL/vNNXyCc7RDJwn1ESHKBiFeG5D1ga9WC/4zIWbKm3xdfdGnXtRXwApgb/U9Wfm7lvq/rRLJS2EXXVjLxSiQhSpFGAFweQJ1BCcIhm94jNHdCKNAvqseaiYOKsloIPyi8EGJY+KL9Ru6P/auIb2Q4rP5OifhfkF0isqxXXPYTXr7o2wY+q9aqFuRnsVALQPmMP0eZwOLBdR+So8GNys9M+BHafkOZhIKLVRrRz0yqQsYKlnEi08mI6SB2E/BlLQ9mORrR2DKEl6LRwEV0i13Xdqvn95GRfDey6nT9A3+M1zv+vEvF+kS6cLt6PpaDyY8Mk1I8f1S9JjDRz07zFa5tGPx2MjFKw6rne5PIIJiUIieCC/bKEZbdeF4mEmLn2SLJ0cSF7KFKW7e166aL7zfWGhfdND14Ccp1IsLa3/D7pRNJT7p9n4ew7++DcOivG0bamEOk26uu47DvRT+HQ9+PIyIX6ieIV3Lm/zK/rHUmVCH3Eb3onEkVU5Y2frwCgCyM/6yGZbwZlus1azkyXGpUz/9L5SvddNaoHgs2JrG2w4Hfny/1/HG93KrOu931xznc9NcXu8/nPhylvUO9G2WVy4X8MeyWvq1PYhNnnMkBc3Iz2yBCOybxbIdof1Ovnl/WMPhnDjbep5wLTEq6TnJwGvrEUwghhBBCCCHEXNGDpxBCCCGEEEKIuaIHTyGEEEIIIYQQc+VcazxDKNAvq9/XLoP/Dn5c0zkY+vqE/vCOayuD/473Rst/Z3xne6/aQB6/7QJ5Ger9vkYyXPIvNjZS4xl/V972yAtwM78v1lq+Jq5JamIOey9Wfh7kpKaDlElkOanpIfVvzeh75an1b4F8vz3MsFizZF9wXyABwdV0xvWcgK9RZvWcZelrG4akMKdJduf+YfVg7zR8HVh21RfPhQeu+rb7rrg265MXl5fVvmUv7vlpyEuM90ldxx6pnxsMq3WwJdmvtJ6WjDGF+ba4toG+SHm54nbumOVYa1bHwQ6p57SorjYzP86wl5E/RV4W/iWkFrmMjsOg48ejcJvUnpLaZDRIbRCpiStvVvvW6/htOiR1eHsD0jYkdWBZdT/uBF9rX+Caa+uYrw9nxLVnrDarQeo5aV0XqT1itW1xG7tGNJbs998NNHG1fKDSdtC87qYbRmPNWnbBTdMmdbSHQ3++3Br4fXA7ys0WrX/zde7Znm/Dvs9IiE8i+DrSQ1J/f0DqObuknpPVMcdrZHVtrczvi7W4ro0tDL4Gc0iyG0gNG6vDmxR2brTPUP92HjQy4GLkZLh65K/9RXRPWoDcj5Z+23ql38d7A79fdqN8bfd8xovrN11b40EytnfIPQm5N+oW1X50Ct8v1sbyXATfRk4r3y/S1iCdbRI/RSBjdEzsRwGAHP65pCD11OxciMftNqm5XsvSHyeXa8QXQgghhBBCCFE79OAphBBCCCGEEGKu6MFTCCGEEEIIIcRc0YOnEEIIIYQQQoi5cq5yISA4mdCQCIFiGcug8IXxReGLiwMpfGbigk6n+jLfHTKfXdjwy7/fS1bCFS8cApMLRcXWdn3Xz1f67dzv+CLeHulvWVbFF4X5FxYzyUr8EmyAF9/HbYG+IH65XgS+KGIpDdtXXlyT9gJrtqynDn0h+UFUpE/cLgDJOC6TPF++7PvR8xIfJ2TZJFKYgX/x8z4Rsnzyzvj90ch9/5nEhr3sft4wiUrKdGwu9jLoRZFbE1utquCGSpgimvmma9vMiciKmBle6PoxfDd6sf3Boc/a5RvkxeO3b/vOMVHWrpfXFfvV7Tzq+DF2r+fbDohIqOcvEcijw7xJxvBh8AKbJsh0RCQRw2RAk0qDgHThXAyT5S2SBjJctOrYsmMPuOk6eVWe1oaXQbVLIkshpwuTqe32q2PZhaYfJ4cv+XG4dZvIhYi0Kwx8R7rD6joPh348PWJyoYLI4BIuaSwdTSIcapdEHpYgQRmSW9wi+HNjlpJClvk1ImxcJJkZtprVc/hi7kVYnaJ6z9s1L11rkH08JOP44dC37Ua5v9D0yxq85I9XfpvIsvaJJI4wKKvbfVSQ8ZmIhPokzyzjRUKUSMSRJ8qFSquuoKRjtl+WkTG7QeZl54KTCxFZ1lqe/jmmPvEUQgghhBBCCDFX9OAphBBCCCGEEGKujH3wNLP3mNl1M/voibbLZvZ+M/vj0d/k+3lCrAbKuKg7yrioO8q4qDPKt6gLKZ94Pg7gTVHbdwD4QAjhNQA+MPpZiFXlcSjjot48DmVc1JvHoYyL+vI4lG9RA8bKhUIIHzKzR6LmrwPwVaN/vxfArwH49ll1KhaoGCmwpbYUUmz+TPa8a7t9WP2l0P1DXzCNTV9ojQs7vu2i/wVTGJAi5/2q6MLI8ocdXxD8/KEXqDyVfdy1NfNq3xr5upsmJwXudN8mwAqVF0GqxOVuzDLjBi+4CZmvQM9CdZosI3ljyye5/8yez9vtSC5UeKcPQGQVYd3nJrR93+heb0SF6uR87HZ8/58+9Ot8vus7vLlWFXwwMVYr8xKblvlzqAEv/WhExyQPTFTERCt+b9Di/gklQcuU8RxNJ1oZNohoKmIt81Kc7eBlba3Myw46xJNzo1fN7lWSoYdevOHamtd9G7a9eAU3vVyo7FXzfKfjz4v9gT+nmLyCCShiqcpGTsQoBcutX/6QCJ9SJFDT5Jvh5Vl+vnzCa1DMzDJu5kQ+l/pehJVn1WmawR+bJrnFYkehSzK+H0mpYqEWAAxu+fmaN3x20SdCnYPxeThgIiEqXvFtLOOxc4ZKVkgcWhmRrBAplYXqAnMi4ytIxlPlQinTMaElk8SclVnep2QANqJo7hCxT6eoXk8PiIyGjRlsPzGh2n4kF7ozIHm77fvVuunvD8IhkW4Oxu93Jsbisqw0gVZKkpgwMM/8nDlZZzM+F0gfhoHcf5Cxt0iUWsbjOJNltWI73l2XNxnXQggvP9G9AODa3SYWYgVRxkXdUcZF3VHGRZ1RvsXKMfXrVEIIwcxOfcg3s7cDeDsAZEQPL8Syc5aM25Jp04VI4W4ZP5nvRkZewSPECpCa8Q3y6bwQy85Z7lN2Gsq4WByTfuL5opk9CACjv6+fNmEI4V0hhMdCCI+xd+wJsaRMmPHleiedEHchKeMn892wtK+EC7EknDnja/rlilgdJrpP2SAlKUKcF5M+eL4PwFtH/34rgF+cTXeEWBqUcVF3lHFRd5RxUWeUb7FyjP0I0sx+CsfFy/eZ2TMAvhvADwD4WTN7G4CnAfz1tNUZsuiriBkpbs2z6ldyY9kQAIQ8rSj2ELdd2x/d+bzKz1+4d9NN0ybF7IHJWFLFK1EVfegO3CR3bnhBxkfu+LY1bLm2WLzCREKtzM/Xgv/tbjP4r0THohUml8jJ7zGYUIXtHy5jOfs0kzDLjBtytPLqfh6YL4SPpU556cUUgP86TKux7dpull6O9amDqmzqi2/6Y7p5SIxDTIxVknONTGeH1X4Uz+65aV64fdG1PX3kh6Hc/Plxf+OLKj934JfPRAcs4+3gz6tWqJ7LDSKvYOKI1NznTNISTbfsGc/RwKWyKlrp5z5HJaomCTZmbZVe1ubECeAChzuR9ORWz58/PfK7/+bzL7k263i5ULjhs1VE8re9vh9j9xNlLOz7cLGbgcka1ojcI2PSi0BkMokClRQmlguRMM9CvDJa9kwynhmw2az2aafvx4syGhczJpUhYwiDy4Us+pmIt3b92LnOxCtDfxIND8Yfw07h18lkWT2SwSGJW+ybY2NbzjJCIlIyEVbUjyEVqpB7OyLCm/R8YbIsJkc6KzO9TzGgFXVps+n31dawOsaV5KCm7qcesU11IonP4dD34fDA37ts3fYZLzp++X1yXcijbyN3Sn9suuSa0ydtKQItBjvzmPyH5j4+11i2SF/Z8pmYjhGP2yzPjTPcqKRYbd9yyn/9xeS1CLHEKOOi7ijjou4o46LOKN+iLizHOzGEEEIIIYQQQtQWPXgKIYQQQgghhJgrevAUQgghhBBCCDFXzvX9Jma5E9wYKW6Ni/TzzBcINxtEipN5AcAmLrm2P9yvFkx/6XN+WZu3D1ybdUjRPpOxkOnszn7l5+I5L6/43O0rru25I1+we1/5gGsrG9Vq4gJezpLDyzDWg9dqx5IVAGhG8zJpAiuqZ0KVBpmO1SV7McX4aRZNZjnWs0tRmz/NSpJVvyw/31Z2v2trl366F7rVffzCTS9yufrMLddmr/IyLjCB1p07froXqpKu3p/4DD535M81VrS/mfttujKsbvuheWFNgLd0sNw3A2mLpmuxjDO5UGLuY5EQkJZxJttYFDkyXLDqMeyEq266ePxpw2doHURiRraVCRxiCcXuwB/Pozu+bePZXdeWEdFb+ZIXdvX2q5kcECnF4dC3sXwT14s79s3M74u1nMiFmOgt+LYyEqjMTjV0TEpK2TkwC/HKLMnMsN6o9nOr4bM0HFbH8ED2KBsbGANiJInFK0eF33cHh/7eaOelQ78CJks58FkqYjkPCQkTIaVmPIbtHeLUQk4mbDK5UHQM2FstG0wIQyaMzxeAnzNxd9n9DZNVLhKDF9esMZlZtOP75F5jQOSfjCHJ+FEUMCapOuj4jF++7cfnoufn7fV8f0N0xMhhRo+ca0xyx65NbHkxZGinuWfuH5d7Jn9k8h/SryEbGAhxplme2fXqNJbrbBBCCCGEEEIIUTv04CmEEEIIIYQQYq7owVMIIYQQQgghxFw51xrPDDk2rFr/1jNfO1FEL48P5HvI7EXx63bBtV0qfd3kblSW+bldX//20FM3XVv+8IuuLWuQXbjv60PxbPVt5t1P+7qi58l32RmXSH1gUV6r/DwwUntKaAZfZ9UksYhrOpvG6t/S6jlZHVdK/Sb7Cjl7IfkiydDAllUzx2o140yzPLfg6yF3gq9Z3shIfduw+vNzh35Zf+ZTN1xb+/Oec20M2/U1nuVTL1V+3nvR5/mAvASd1U5sNPz+uFBUt6EZ/H5lNQusBo7WxUXHoJFaz5lQu3nadDFskrO8mHne5DBs59W8HQ233XQDVAPIxpQ1MvazOpEipf6N1FYeHPi60ovP7Lu2Rs8XrQ1f9OPn0VG1Hp7VeLKXkbN6IVYGFNc/speHF7SAiGSS1ErF86a++H2WsPOiuYT1b62oS2w8GpbVayerdUutzy7IGBiP4QNyTA+7/vrdu+5rPFmpaefQn3/dwo/PMX3SD1bjmZKu1Jr2BtuNmV+Du2cgRXJDVrvJaqJZgWgCtP9nqH87Dwy+ppCNN3HdZ5/kIyN5YMeVVRPGQy+rmTzo+4x3d8m1n1wDOn2f8V5RnY71ldU2p9Z4plVNelhGWPmsyyUZPzNyzWR5zsP48x3wmW6RglSWn9NYrhFfCCGEEEIIIUTt0IOnEEIIIYQQQoi5ogdPIYQQQgghhBBzRQ+eQgghhBBCCCHmyrnKhXLk2AoXqx0ggonSyBuKI5rBi0s2Sv9C+QvmBRNx3e1LPb+szh97+c/mg168YkPS133ycttI5HL7BS8I6hb+9wCsUHmTCI2KQVV80SP7h73g2pgQKEGqwgroU9uo7CVBHMR+S7JE3hUAQI4GLpSXK22WIM9oEFFOK/jsbhPh0EaDvAg8OtR3Bn75+88QucuniUCLZDzs+YwPnqq27e5fdNN0iYiAKRxaJBCbeXUbUgQqAJeosAzGLexl96l5TpFlHa/DT5eyzkWRGbAeiVa2Cz/W9KM8MylTi5wX7EXxjPgl9uzF43td36/ui1781u53XVvnhu/IfiR/6xKR0JAM2Ey8QtwPLh3sheLsvGDpID4jFJFAJaS86Ry8r5PC8n4WKcV5YPB9WiPB7OfVHcPGIz4OeNi4FUt8mMhkl9y7XNn1MpYs9zMzMVEsF2JCI5YHlvtJxStsDAgs5WRh8TDD/EA5EQlxaVfatseZZuP1smUc5vvNBE6xRIZJZVii2TWXjTe9KCQlOTb7A3+fcrDvc1+Q8fiIyIUGbhx0k/CMJ04Xw+4jjASTZpXu7liglXZ/U5K21GtAnGkmPMzPINBattNBCCGEEEIIIUTN0IOnEEIIIYQQQoi5MtVXbc3sKQD7AAoAwxDCY7PolBDLgjIu6o4yLuqOMi7qjPItVolZ1Hj+hRCCfxO9EPVBGRd1RxkXdUcZF3VG+RYrwbnKhTLk2AlVAVAevBjFz+e/EdyELxrehC+W3yIinpj9ge/Djec3XVv7D2+5tkbXS4jCYd+1dT9dne7WwY6fhhRHs+LiddIYIplHiyyLF9B7WGFyXFzMpAms4DgWBB0v30MlDAm1ykyOtEhyZNhBJI5KMCrk8Blsk9NzK/MZZwX/8aE+JOKqG7e9jGvzk7f98vvXXVu573N/8Llqf3eJ3CUWZgC8QJ9tUyy1yYjXa0gynhh7xzSCoNTcJ61ziSKeGbARWSg6Q5/TRjT+sO1qM7lQokgpFpcwadVuz58reze9sGuz58frwwM/70EkquiRc2pAssYkKwwn+6CBYfuHyCWIpCP2y7BuzVIklAq7biwSM3/dbZFj0c7jjCeKOhILnIroYFCBFhGv3Dnw4sI88xchKt+K5EIsu6wtNTbxpjPxChOepEhWALgAs+WzjLO70EnPBS4XWq6MM1gX/Xkw/l4D4CKblNgzgc8+ESPe6fhxnHFIzo9Y4pk6Ps+SREcTQM55f7aR+24q6CLiRXKd4DmoNrI8nyXh09Z4BgD/p5n9jpm9fcplCbGMKOOi7ijjou4o46LOKN9iZZj2E8+vCCE8a2b3A3i/mf1hCOFDJycYnQRvB4C2bU+5OiHOnTNlfM38J9lCLDl3zfjJfG/lFxbVRyGmITnjOw1lXKwcZ7pPuaiMiwUy1SeeIYRnR39fB/BvALyBTPOuEMJjIYTHmpl/B6EQy8xZM95SxsWKMS7jJ/O9rnyLFeQsGd/MfZmNEMvMWe9TNhvKuFgcEz94mtmm2fFHmGa2CeAvAfjorDomxKJRxkXdUcZF3VHGRZ1RvsWqMc1Xba8B+DcjiUYDwE+GEH71bjPkMGxn1YL2nEhw4iLYnMmFzJeDr2VE0EIq+eO62KPCl8W+sOfFKzuf7bq2jc6eayv9ZLgVyYpudn1xNJMHsELftQYp8I72UU62iTluWCE/I5alTCMIYkKgNJFQ2vJnyJkznsGwmUcF7USCU8YZJ9XgbZLx9YZvYwX/MV2Sh+tHXkJx4bmOa9vu9VzbsO+Xd+N2NeO3+17QkppxVrzedpvu91ljhnIhxjTynxQRFj9f5sqZMm5mTvwUS58AIKdShCpNMogw6QLb57H4Y0AGt10ilrh14D+x7RN5xWGfSFuiPDOhUUFkDamkiFdyJrAhocmY8CM6EZhYgoVtlu4NtnfOQbxy9nHciZ5SxqM0WMZTpFos43dIdm8T8UqDBOKAzHswrG7UgGRkllIqttV5gjTo1OVF8/K+pklWJj2V+fGdbFmJnDnfCGm71Etl6KIcTGSZct/HM+5PtNtEjMXuSfeJ+C6WCzGh0TRjXtwPdt+dKtVKEQ7x0yUt46k3F/H4xwRp5LHk9OWlT1olhPAZAF8y6fxCLDvKuKg7yrioO8q4qDPKt1g1prXaCiGEEEIIIYQQd0UPnkIIIYQQQggh5ooePIUQQgghhBBCzJVp3+N5JjIzbOTVQuEsobqV1d03iYyFSXc2SFtc6M0kES/1fPHy9nX/jsYLh94kVJLlXd+vyopuEXlFl0hoWFE6kzGEqP6azVdQ4cRkVe/TSFBSXRKsAHvSdZ4XmTHZij/WXi7kt6SZkYyTAxvLXgBf/N0lFe43el7+s33bv2u30/XTFeycOayKW5j4ok/6QQUTpLEVbROTFfCM+7ZJ5QHT5G1Sh0qKeOS8MADNqDsskynbmjq2sXEgnozJs5iUgkndjmhO/bm3F4kqjgo/zZCZJAh0/IzaqLuGzGhM5MEynzDWp0piZnn+MFHFssHEGd6plXaepu6DWLLDxs47Az/jTXLv0sz8zN3CJ6wTnUdM9sIykuoDijPObkCZ7CVVOBQvn/la2HWDjUWpGY9nZV2l/V8gAWnnepx7JhdiiWbXYUZ8baNyoaFf6S4RFxoZ9GKREODv94dkXKTZ9U089/E0ZD62/GSpVrTAaTKeSjw+pUoAT2MFhnwhhBBCCCGEEKuMHjyFEEIIIYQQQswVPXgKIYQQQgghhJgrevAUQgghhBBCCDFXzlcuBC9eSRFnULkQKcRlbUy8EhdIMzHFzZ7fNWvZpmu7TWQVjDuDZvSzL+xnRc5UwMFqkGPxCtlpRaL4IqUufNYiIbqOhHmX7TcnmXnZSkZEWHHNONtPDbID2jTP43PTIxm/1fcZX++su7bdvpdVsJr3/Ui+wsQXAyIlojIGKvOIG30ncpJMVnyfKlGZlFm6JJbILQTAiwZaxIITC4G4TIG0JY538XRMSLJH8nej50VfzcyfBzTf0fKYDI5JNVKdEZPKH5hIgq9zfOgTJjmVlMsLG6+XTbzCSLkHych4xPZJ6j6Ic8/EK/tEvMLG9ZyIVwpyvxELs/pkvGZMej/A9s80wqF4KGKZT8144u2SO54sK9OIXeZFigDICc/ouZp2HWbE+4V16XDI7l38RYftd5bx+H6fnVeM1HHczccaU68JKTI5Jsaa8f2Ny0GC8O9uLNt9uxBCCCGEEEKImqEHTyGEEEIIIYQQc0UPnkIIIYQQQggh5sr51ngasB69kZbVp6XAa8B8G6sPipvYd93ZS2ut619ay57c2SbF3yPfJ99bT335OKuNakZbldEaorQvfqd8P3zWdWfsJfGriMHXeLLvw8e7OPml4ok1cCkZ3yMZzGjG03IT1ygfkrrS9BdL+zZXrkFONLb4lNq2ZWaZzg2DzxYbwyct22PHnZ0HKYs/Ivl7idTus93LMh/Xu6XWv01K6j6kLygnbSk+hXnXPjOWscQzpQY/HnctcTxim9tIzH0MqzO+RWqbUz9h6EWZZvVvqRFJqYmj05BlsRvVlNslVv+cWpCaup0pi1u2T3gC/LWYnfvx8WHj8zS5j2dl06RmPHUYie9TUu9JZlnHTLQIE4/jdMxm909kMkbKcZo2z8t2PgghhBBCCCGEqBl68BRCCCGEEEIIMVemevA0szeZ2R+Z2afM7Dtm1SkhlgVlXNQdZVzUGeVb1B1lXKwSEz94mlkO4J8D+FoArwPwFjN73aw6JsSiUcZF3VHGRZ1RvkXdUcbFqjGNXOgNAD4VQvgMAJjZTwP4OgAfP20GM3PilWLCR99ZFuizAnr20tpuwUqC04iLc5lIiL0kmUFFNCnFv0QukfqS5ElZ8e9yT5TxViwXStgJqYXrk2acFdAfkYz3iJBlUlEAW2fqy5oZ8XYyVxb1CC0g95OygPPlzBmPj3OTSbCi48COC5X6TOgCYVnrk4O8OyBZoFKNyUQ8qaIKhpNSUEGcbyPvSKfyipjU5c+bOcuFzpzvEPy+SZG/pY5Hs8w9u2fY7c9O9sUyMksBVcp15LQJU1yJqXmepX8u9fjOkIkyHl+LU67zVP7GxqQpch/DxtR9Mo4zUgRXqYd+UjlP6rFPlRC5aRL3/zTEx27S8erfzz9FXx4C8LkTPz8zahOiLijjou4o46LOKN+i7ijjYqWY+y/YzeztZvaEmT3RKQ7nvTohzh1lXNSZk/k+KpVvUT+UcVF3lHGxLEzz4PksgFee+PnhUVuFEMK7QgiPhRAeW883p1idEOeOMi7qztiMn8z3RqZ8i5XizGO4Mi5WDGVcrBTTPHj+NoDXmNmrzKwF4G8AeN9suiXEUqCMi7qjjIs6o3yLuqOMi5XCwhRVqGb2ZgD/BMc1sO8JIXz/mOlfAvA0gPsA3Jh4xYtn1fsPrP423K3/nx9CuDqLlSjjK0vd+7+QjJ/Id0oflx31f7Es+xg+ro+rgPq/WJY943Xev6vAqvcfmCDjUz14ToqZPRFCeOzcVzwjVr3/wOpvw7L3f9n7Nw71f7GsQv9XoY93Q/1fLKvQ/1Xo491Q/xfLsvd/2fs3DvV/8UyyDSv+tgshhBBCCCGEEMuOHjyFEEIIIYQQQsyVRT14vmtB650Vq95/YPW3Ydn7v+z9G4f6v1hWof+r0Me7of4vllXo/yr08W6o/4tl2fu/7P0bh/q/eM68DQup8RRCCCGEEEIIce+gr9oKIYQQQgghhJgr5/7gaWZvMrM/MrNPmdl3nPf6z4qZvcfMrpvZR0+0XTaz95vZH4/+vrTIPt4NM3ulmX3QzD5uZh8zs28dta/ENpjZmpn9lpn9/qj/3ztqf5WZfXiUo58Zvb9qKVDGzxdl/HxZtXwDq53xVc83oIzPm1XON7D6GV+1fAPK+HmjjP8p5/rgaWY5gH8O4GsBvA7AW8zsdefZhwl4HMCborbvAPCBEMJrAHxg9POyMgTwjhDC6wC8EcDfHe3zVdmGHoCvDiF8CYBHAbzJzN4I4AcBvDOE8GoAtwG8bXFd/FOU8YWgjJ8TK5pvYLUzvur5BpTxefM4VjffwOpnfGXyDSjjC0IZf5kQwrn9AfBlAP6PEz9/J4DvPM8+TNjvRwB89MTPfwTgwdG/HwTwR4vu4xm25RcBfM0qbgOADQC/C+BLcfzC2saovZKrBfdRGV/8tijj8+vfSuZ71NdaZHyV8z3qqzI+n37XIt+j/q5sxpc936wvyvhCtuWezfh5f9X2IQCfO/HzM6O2VeNaCOH50b9fAHBtkZ1JxcweAfB6AB/GCm2DmeVm9iSA6wDeD+DTAHZDCMPRJMuUI2V8gSjjc6cu+QZWKB8vs6r5BpTxBbBS+XiZVc34CuUbUMYXyr2eccmFpiQcP+YvvRrYzLYA/ByAbwsh7J38v2XfhhBCEUJ4FMDDAN4A4LWL7dG9xbLn42WUcTEpy54PYLXzDSjji2QV8gGsdsaV78Wy7Pl4GWX8/B88nwXwyhM/PzxqWzVeNLMHAWD09/UF9+eumFkTx0H/iRDCz4+aV2obACCEsAvggzj+OP+imTVG/7VMOVLGF4Ayfm7UJd/ACuWjLvkGlPFzZKXyUZeMr0C+AWV8ISjjx5z3g+dvA3jNyILUAvA3ALzvnPswC94H4K2jf78Vx9/VXkrMzAC8G8AnQgg/cuK/VmIbzOyqmV0c/Xsdx9+J/wSOQ//1o8mWqf/K+DmjjJ8rdck3sDr5WOl8A8r4glilfKx0xlcs34Ayfu4o4ydYQFHqmwF8EsffDf77573+Cfr7UwCeBzDA8feX3wbgCo7tU38M4P8CcHnR/bxL/78Cxx/dfwTAk6M/b16VbQDwHwL4vVH/PwrgH47avwDAbwH4FID/DUB70X090Wdl/Hz7r4yfb39XKt+jPq9sxlc936NtUMbn29+Vzfeo/yud8VXL96hvyvj59l8ZH/2x0YxCCCGEEEIIIcRckFxICCGEEEIIIcRc0YOnEEIIIYQQQoi5ogdPIYQQQgghhBBzRQ+eQgghhBBCCCHmih48hRBCCCGEEELMFT14CiGEEEIIIYSYK3rwFEIIIYQQQggxV/TgKYQQQgghhBBirvz/LHMYK2UF5/IAAAAASUVORK5CYII=",
"text/plain": [
""
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"fig, axes = pylab.subplots(1, len(steps), figsize=(16, 5))\n",
"for i in range(len(steps)):\n",
" axes[i].imshow(steps[i][1].numpy('y,x'), origin='lower', cmap='magma')\n",
" axes[i].set_title(f\"u_x at t={i*5}\")\n",
" \n",
"fig, axes = pylab.subplots(1, len(steps), figsize=(16, 5))\n",
"for i in range(len(steps)):\n",
" axes[i].imshow(steps[i][2].numpy('y,x'), origin='lower', cmap='magma')\n",
" axes[i].set_title(f\"u_y at t={i*5}\")\n",
" "
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "ooqVxCPM8PXl"
},
"source": [
"It looks simple here, but this simulation setup is a powerful tool. The simulation could easily be extended to more complex cases or 3D, and it is already fully compatible with backpropagation pipelines of deep learning frameworks. \n",
"\n",
"In the next chapters we'll show how to use these simulations for training NNs, and how to steer and modify them via trained NNs. This will illustrate how much we can improve the training process by having a solver in the loop, especially when the solver is _differentiable_. Before moving to these more complex training processes, we will cover a simpler supervised approach in the next chapter. This is very fundamental: even when aiming for advanced physics-based learning setups, a working supervised training is always the first step."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Next steps\n",
"\n",
"You could create a variety of nice fluid simulations based on this setup. E.g., try changing the spatial resolution, the buoyancy factors, and the overall length of the simulation run."
]
}
],
"metadata": {
"colab": {
"collapsed_sections": [],
"name": "overview-ns-forw.ipynb",
"provenance": []
},
"kernelspec": {
"display_name": "Python 3",
"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.5"
}
},
"nbformat": 4,
"nbformat_minor": 1
}