{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Mesh examples" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "this notebook illustrates the basic ways of interacting with the pyro2 mesh module. We create some data that lives on a grid and show how to fill the ghost cells. The pretty_print() function shows us that they work as expected." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import pyro.mesh.boundary as bnd\n", "import pyro.mesh.patch as patch\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline\n", "\n", "# for unit testing, we want to ensure the same random numbers\n", "np.random.seed(100)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setup a Grid with Variables" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "There are a few core classes that we deal with when creating a grid with associated variables:\n", "\n", "* `Grid2d` : this holds the size of the grid (in zones) and the physical coordinate information, including coordinates of cell edges and centers\n", "\n", "* `BC` : this is a container class that simply holds the type of boundary condition on each domain edge.\n", "\n", "* `ArrayIndexer` : this is an array of data along with methods that know how to access it with different offsets into the data that usually arise in stencils (like {i+1, j})\n", "\n", "* `CellCenterData2d` : this holds the data that lives on a grid. Each variable that is part of this class has its own boundary condition type." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We start by creating a `Grid2d` object with 4 x 6 cells and 2 ghost cells" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "2-d grid: nx = 4, ny = 6, ng = 2\n" ] } ], "source": [ "g = patch.Grid2d(4, 6, ng=2)\n", "print(g)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "tags": [ "nbval-ignore-output" ] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Help on Grid2d in module pyro.mesh.patch object:\n", "\n", "class Grid2d(builtins.object)\n", " | Grid2d(nx, ny, ng=1, xmin=0.0, xmax=1.0, ymin=0.0, ymax=1.0)\n", " | \n", " | the 2-d grid class. The grid object will contain the coordinate\n", " | information (at various centerings).\n", " | \n", " | A basic (1-d) representation of the layout is::\n", " | \n", " | | | | X | | | | X | | |\n", " | +--*--+- // -+--*--X--*--+--*--+- // -+--*--+--*--X--*--+- // -+--*--+\n", " | 0 ng-1 ng ng+1 ... ng+nx-1 ng+nx 2ng+nx-1\n", " | \n", " | ilo ihi\n", " | \n", " | |<- ng guardcells->|<---- nx interior zones ----->|<- ng guardcells->|\n", " | \n", " | The '*' marks the data locations.\n", " | \n", " | Methods defined here:\n", " | \n", " | __eq__(self, other)\n", " | are two grids equivalent?\n", " | \n", " | __init__(self, nx, ny, ng=1, xmin=0.0, xmax=1.0, ymin=0.0, ymax=1.0)\n", " | Create a Grid2d object.\n", " | \n", " | The only data that we require is the number of points that\n", " | make up the mesh in each direction. Optionally we take the\n", " | extrema of the domain (default is [0,1]x[0,1]) and number of\n", " | ghost cells (default is 1).\n", " | \n", " | Note that the Grid2d object only defines the discretization,\n", " | it does not know about the boundary conditions, as these can\n", " | vary depending on the variable.\n", " | \n", " | Parameters\n", " | ----------\n", " | nx : int\n", " | Number of zones in the x-direction\n", " | ny : int\n", " | Number of zones in the y-direction\n", " | ng : int, optional\n", " | Number of ghost cells\n", " | xmin : float, optional\n", " | Physical coordinate at the lower x boundary\n", " | xmax : float, optional\n", " | Physical coordinate at the upper x boundary\n", " | ymin : float, optional\n", " | Physical coordinate at the lower y boundary\n", " | ymax : float, optional\n", " | Physical coordinate at the upper y boundary\n", " | \n", " | __str__(self)\n", " | print out some basic information about the grid object\n", " | \n", " | coarse_like(self, N)\n", " | return a new grid object coarsened by a factor n, but with\n", " | all the other properties the same\n", " | \n", " | fine_like(self, N)\n", " | return a new grid object finer by a factor n, but with\n", " | all the other properties the same\n", " | \n", " | scratch_array(self, nvar=1)\n", " | return a standard numpy array dimensioned to have the size\n", " | and number of ghostcells as the parent grid\n", " | \n", " | ----------------------------------------------------------------------\n", " | Data descriptors defined here:\n", " | \n", " | __dict__\n", " | dictionary for instance variables (if defined)\n", " | \n", " | __weakref__\n", " | list of weak references to the object (if defined)\n", " | \n", " | ----------------------------------------------------------------------\n", " | Data and other attributes defined here:\n", " | \n", " | __hash__ = None\n", "\n" ] } ], "source": [ "help(g)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Then create a dataset that lives on this grid and add a variable name. For each variable that lives on the grid, we need to define the boundary conditions -- this is done through the BC object." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "BCs: -x: periodic +x: periodic -y: reflect-even +y: outflow\n" ] } ], "source": [ "bc = bnd.BC(xlb=\"periodic\", xrb=\"periodic\", ylb=\"reflect\", yrb=\"outflow\")\n", "print(bc)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "cc data: nx = 4, ny = 6, ng = 2\n", " nvars = 1\n", " variables:\n", " a: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: periodic +x: periodic -y: reflect-even +y: outflow \n", "\n" ] } ], "source": [ "d = patch.CellCenterData2d(g)\n", "d.register_var(\"a\", bc)\n", "d.create()\n", "print(d)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Working with the data" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we fill the grid with random data. `get_var()` returns an `ArrayIndexer` object that has methods for accessing views into the data. Here we use `a.v()` to get the \"valid\" region, i.e. excluding ghost cells." ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "a = d.get_var(\"a\")\n", "a.v()[:,:] = np.random.rand(g.nx, g.ny)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "when we pretty_print() the variable, we see the ghost cells colored red. Note that we just filled the interior above." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.12157 0.2092 0.17194 0.33611\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.0047189 0.89132 0.81168 0.81765\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.84478 0.57509 0.97862 0.94003\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.42452 0.13671 0.2197 0.4317\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.27837 0.82585 0.10838 0.27407\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.5434 0.67075 0.18533 0.81622\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\n", " ^ y\n", " |\n", " +---> x\n", " \n" ] } ], "source": [ "a.pretty_print()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`pretty_print()` can also take an argument, specifying the format string to be used for the output." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.122 0.209 0.172 0.336\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m0.00472 0.891 0.812 0.818\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.845 0.575 0.979 0.94\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.425 0.137 0.22 0.432\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.278 0.826 0.108 0.274\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.543 0.671 0.185 0.816\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\n", " ^ y\n", " |\n", " +---> x\n", " \n" ] } ], "source": [ "a.pretty_print(fmt=\"%7.3g\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "now fill the ghost cells -- notice that the left and right are periodic, the upper is outflow, and the lower is reflect, as specified when we registered the data above." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\u001b[31m 0.17194\u001b[0m\u001b[31m 0.33611\u001b[0m\u001b[31m 0.12157\u001b[0m\u001b[31m 0.2092\u001b[0m\u001b[31m 0.17194\u001b[0m\u001b[31m 0.33611\u001b[0m\u001b[31m 0.12157\u001b[0m\u001b[31m 0.2092\u001b[0m \n", "\u001b[31m 0.17194\u001b[0m\u001b[31m 0.33611\u001b[0m\u001b[31m 0.12157\u001b[0m\u001b[31m 0.2092\u001b[0m\u001b[31m 0.17194\u001b[0m\u001b[31m 0.33611\u001b[0m\u001b[31m 0.12157\u001b[0m\u001b[31m 0.2092\u001b[0m \n", "\u001b[31m 0.17194\u001b[0m\u001b[31m 0.33611\u001b[0m 0.12157 0.2092 0.17194 0.33611\u001b[31m 0.12157\u001b[0m\u001b[31m 0.2092\u001b[0m \n", "\u001b[31m 0.81168\u001b[0m\u001b[31m 0.81765\u001b[0m 0.0047189 0.89132 0.81168 0.81765\u001b[31m 0.0047189\u001b[0m\u001b[31m 0.89132\u001b[0m \n", "\u001b[31m 0.97862\u001b[0m\u001b[31m 0.94003\u001b[0m 0.84478 0.57509 0.97862 0.94003\u001b[31m 0.84478\u001b[0m\u001b[31m 0.57509\u001b[0m \n", "\u001b[31m 0.2197\u001b[0m\u001b[31m 0.4317\u001b[0m 0.42452 0.13671 0.2197 0.4317\u001b[31m 0.42452\u001b[0m\u001b[31m 0.13671\u001b[0m \n", "\u001b[31m 0.10838\u001b[0m\u001b[31m 0.27407\u001b[0m 0.27837 0.82585 0.10838 0.27407\u001b[31m 0.27837\u001b[0m\u001b[31m 0.82585\u001b[0m \n", "\u001b[31m 0.18533\u001b[0m\u001b[31m 0.81622\u001b[0m 0.5434 0.67075 0.18533 0.81622\u001b[31m 0.5434\u001b[0m\u001b[31m 0.67075\u001b[0m \n", "\u001b[31m 0.18533\u001b[0m\u001b[31m 0.81622\u001b[0m\u001b[31m 0.5434\u001b[0m\u001b[31m 0.67075\u001b[0m\u001b[31m 0.18533\u001b[0m\u001b[31m 0.81622\u001b[0m\u001b[31m 0.5434\u001b[0m\u001b[31m 0.67075\u001b[0m \n", "\u001b[31m 0.10838\u001b[0m\u001b[31m 0.27407\u001b[0m\u001b[31m 0.27837\u001b[0m\u001b[31m 0.82585\u001b[0m\u001b[31m 0.10838\u001b[0m\u001b[31m 0.27407\u001b[0m\u001b[31m 0.27837\u001b[0m\u001b[31m 0.82585\u001b[0m \n", "\n", " ^ y\n", " |\n", " +---> x\n", " \n" ] } ], "source": [ "d.fill_BC(\"a\")\n", "a.pretty_print()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can find the L2 norm of the data easily" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.5749769043407793" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "a.norm()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "and the min and max" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.004718856190972565 0.9786237847073697\n" ] } ], "source": [ "print(a.min(), a.max())" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## `ArrayIndexer`" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We we access the data, an `ArrayIndexer` object is returned. The `ArrayIndexer` sub-classes the NumPy `ndarray`, so it can do all of the methods that a NumPy array can, but in addition, we can use the `ip()`, `jp()`, or `ip_jp()` methods to the `ArrayIndexer` object shift our view in the x, y, or x & y directions.\n", "\n", "To make this clearer, we'll change our data set to be nicely ordered numbers. We index the `ArrayIndex` the same way we would a NumPy array. The index space includes ghost cells, so the `ilo` and `ihi` attributes from the grid object are useful to index just the valid region. The `.v()` method is a shortcut that also gives a view into just the valid data.\n", "\n", "Note: when we use one of the `ip()`, `jp()`, `ip_jp()`, or `v()` methods, the result is a regular NumPy `ndarray`, not an `ArrayIndexer` object. This is because it only spans part of the domain (e.g., no ghost cells), and therefore cannot be associated with the `Grid2d` object that the `ArrayIndexer` is built from." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "pyro.mesh.array_indexer.ArrayIndexer" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "type(a)" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "numpy.ndarray" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "type(a.v())" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [], "source": [ "a[:,:] = np.arange(g.qx*g.qy).reshape(g.qx, g.qy)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "tags": [ "nbval-ignore-output" ] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\u001b[31m 9\u001b[0m\u001b[31m 19\u001b[0m\u001b[31m 29\u001b[0m\u001b[31m 39\u001b[0m\u001b[31m 49\u001b[0m\u001b[31m 59\u001b[0m\u001b[31m 69\u001b[0m\u001b[31m 79\u001b[0m \n", "\u001b[31m 8\u001b[0m\u001b[31m 18\u001b[0m\u001b[31m 28\u001b[0m\u001b[31m 38\u001b[0m\u001b[31m 48\u001b[0m\u001b[31m 58\u001b[0m\u001b[31m 68\u001b[0m\u001b[31m 78\u001b[0m \n", "\u001b[31m 7\u001b[0m\u001b[31m 17\u001b[0m 27 37 47 57\u001b[31m 67\u001b[0m\u001b[31m 77\u001b[0m \n", "\u001b[31m 6\u001b[0m\u001b[31m 16\u001b[0m 26 36 46 56\u001b[31m 66\u001b[0m\u001b[31m 76\u001b[0m \n", "\u001b[31m 5\u001b[0m\u001b[31m 15\u001b[0m 25 35 45 55\u001b[31m 65\u001b[0m\u001b[31m 75\u001b[0m \n", "\u001b[31m 4\u001b[0m\u001b[31m 14\u001b[0m 24 34 44 54\u001b[31m 64\u001b[0m\u001b[31m 74\u001b[0m \n", "\u001b[31m 3\u001b[0m\u001b[31m 13\u001b[0m 23 33 43 53\u001b[31m 63\u001b[0m\u001b[31m 73\u001b[0m \n", "\u001b[31m 2\u001b[0m\u001b[31m 12\u001b[0m 22 32 42 52\u001b[31m 62\u001b[0m\u001b[31m 72\u001b[0m \n", "\u001b[31m 1\u001b[0m\u001b[31m 11\u001b[0m\u001b[31m 21\u001b[0m\u001b[31m 31\u001b[0m\u001b[31m 41\u001b[0m\u001b[31m 51\u001b[0m\u001b[31m 61\u001b[0m\u001b[31m 71\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 10\u001b[0m\u001b[31m 20\u001b[0m\u001b[31m 30\u001b[0m\u001b[31m 40\u001b[0m\u001b[31m 50\u001b[0m\u001b[31m 60\u001b[0m\u001b[31m 70\u001b[0m \n", "\n", " ^ y\n", " |\n", " +---> x\n", " \n" ] } ], "source": [ "a.pretty_print()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We index our arrays as {i,j}, so x (indexed by i) is the row and y (indexed by j) is the column in the NumPy array. Note that python arrays are stored in row-major order, which means that all of the entries in the same row are adjacent in memory. This means that when we simply print out the `ndarray`, we see constant-x horizontally, which is the transpose of what we are used to." ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "tags": [ "nbval-ignore-output" ] }, "outputs": [ { "data": { "text/plain": [ "array([[22., 23., 24., 25., 26., 27.],\n", " [32., 33., 34., 35., 36., 37.],\n", " [42., 43., 44., 45., 46., 47.],\n", " [52., 53., 54., 55., 56., 57.]])" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "a.v()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can offset our view into the array by one in x -- this would be like {i+1, j} when we loop over data. The `ip()` method is used here, and takes an argument which is the (positive) shift in the x (i) direction. So here's a shift by 1" ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "tags": [ "nbval-ignore-output" ] }, "outputs": [ { "data": { "text/plain": [ "array([[ 1., 2., 3., 4., 5., 6., 7., 8.],\n", " [11., 12., 13., 14., 15., 16., 17., 18.],\n", " [21., 22., 23., 24., 25., 26., 27., 28.],\n", " [31., 32., 33., 34., 35., 36., 37., 38.],\n", " [41., 42., 43., 44., 45., 46., 47., 48.],\n", " [51., 52., 53., 54., 55., 56., 57., 58.]])" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "a.ip(-1, buf=1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A shifted view is necessarily smaller than the original array, and relies on ghost cells to bring new data into view. Because of this, the underlying data is no longer the same size as the original data, so we return it as an `ndarray` (which is actually just a view into the data in the `ArrayIndexer` object, so no copy is made.\n", "\n", "To see that it is simply a view, lets shift and edit the data" ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "scrolled": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\u001b[31m 9\u001b[0m\u001b[31m 19\u001b[0m\u001b[31m 29\u001b[0m\u001b[31m 39\u001b[0m\u001b[31m 49\u001b[0m\u001b[31m 59\u001b[0m\u001b[31m 69\u001b[0m\u001b[31m 79\u001b[0m \n", "\u001b[31m 8\u001b[0m\u001b[31m 18\u001b[0m\u001b[31m 28\u001b[0m\u001b[31m 38\u001b[0m\u001b[31m 48\u001b[0m\u001b[31m 58\u001b[0m\u001b[31m 68\u001b[0m\u001b[31m 78\u001b[0m \n", "\u001b[31m 7\u001b[0m\u001b[31m 17\u001b[0m 27 37 47 57\u001b[31m 67\u001b[0m\u001b[31m 77\u001b[0m \n", "\u001b[31m 6\u001b[0m\u001b[31m 16\u001b[0m 26 36 46 56\u001b[31m 66\u001b[0m\u001b[31m 76\u001b[0m \n", "\u001b[31m 5\u001b[0m\u001b[31m 15\u001b[0m 25 35 45 55\u001b[31m 65\u001b[0m\u001b[31m 75\u001b[0m \n", "\u001b[31m 4\u001b[0m\u001b[31m 14\u001b[0m 24 34 44 54\u001b[31m 64\u001b[0m\u001b[31m 74\u001b[0m \n", "\u001b[31m 3\u001b[0m\u001b[31m 13\u001b[0m 23 33 0 53\u001b[31m 63\u001b[0m\u001b[31m 73\u001b[0m \n", "\u001b[31m 2\u001b[0m\u001b[31m 12\u001b[0m 22 32 42 52\u001b[31m 62\u001b[0m\u001b[31m 72\u001b[0m \n", "\u001b[31m 1\u001b[0m\u001b[31m 11\u001b[0m\u001b[31m 21\u001b[0m\u001b[31m 31\u001b[0m\u001b[31m 41\u001b[0m\u001b[31m 51\u001b[0m\u001b[31m 61\u001b[0m\u001b[31m 71\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 10\u001b[0m\u001b[31m 20\u001b[0m\u001b[31m 30\u001b[0m\u001b[31m 40\u001b[0m\u001b[31m 50\u001b[0m\u001b[31m 60\u001b[0m\u001b[31m 70\u001b[0m \n", "\n", " ^ y\n", " |\n", " +---> x\n", " \n" ] } ], "source": [ "d = a.ip(1)\n", "d[1,1] = 0.0\n", "a.pretty_print()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here, since d was really a view into $a_{i+1,j}$, and we accessed element (1,1) into that view (with 0,0 as the origin), we were really accessing the element (2,1) in the valid region" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Differencing" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "`ArrayIndexer` objects are easy to use to construct differences, like those that appear in a stencil for a finite-difference, without having to explicitly loop over the elements of the array.\n", "\n", "Here's we'll create a new dataset that is initialized with a sine function" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [], "source": [ "g = patch.Grid2d(8, 8, ng=2)\n", "d = patch.CellCenterData2d(g)\n", "bc = bnd.BC(xlb=\"periodic\", xrb=\"periodic\", ylb=\"periodic\", yrb=\"periodic\")\n", "d.register_var(\"a\", bc)\n", "d.create()\n", "\n", "a = d.get_var(\"a\")\n", "a[:,:] = np.sin(2.0*np.pi*a.g.x2d)\n", "d.fill_BC(\"a\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Our grid object can provide us with a scratch array (an `ArrayIndexer` object) define on the same grid" ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "pyro.mesh.array_indexer.ArrayIndexer" ] }, "execution_count": 20, "metadata": {}, "output_type": "execute_result" } ], "source": [ "b = g.scratch_array()\n", "type(b)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can then fill the data in this array with differenced data from our original array -- since `b` has a separate data region in memory, its elements are independent of `a`. We do need to make sure that we have the same number of elements on the left and right of the `=`. Since by default, `ip()` will return a view with the same size as the valid region, we can use `.v()` on the left to accept the differences.\n", "\n", "Here we compute a centered-difference approximation to the first derivative" ] }, { "cell_type": "code", "execution_count": 21, "metadata": {}, "outputs": [], "source": [ "b.v()[:,:] = (a.ip(1) - a.ip(-1))/(2.0*a.g.dx)\n", "# normalization was 2.0*pi\n", "b[:,:] /= 2.0*np.pi" ] }, { "cell_type": "code", "execution_count": 22, "metadata": { "tags": [ "nbval-ignore-output" ] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.125\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjgAAAGfCAYAAABShKg9AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAA9hAAAPYQGoP6dpAABtWklEQVR4nO3deVxU9f7H8dfMsAnKuKCAior7Bi6YCuaakua+lN26Vjez7N4W85pptqiVpv0q22zVrK5XvYmW5Zbmrqi5b7gvoIKIygCyw/n9cRBF3FBmvrN8no/HPDqOZ4b3OCYfzjnzfRs0TdMQQgghhHAiRtUBhBBCCCFKmww4QgghhHA6MuAIIYQQwunIgCOEEEIIpyMDjhBCCCGcjgw4QgghhHA6MuAIIYQQwunIgCOEEEIIpyMDjhBCCCGcjgw4QgghhHA6btZ88nXr1vHBBx+wfft24uPjWbhwIf369bvlY9auXcvIkSPZv38/VatWZfTo0QwfPrzIPlFRUbz55pscO3aMOnXq8N5779G/f/87zpWfn8/Zs2cpV64cBoPhbl6aEEIIIWxM0zRSU1OpWrUqRuNtjtFoVrRkyRJt3LhxWlRUlAZoCxcuvOX+x48f17y9vbWXX35ZO3DggPbtt99q7u7u2vz58wv32bRpk2YymbRJkyZpMTEx2qRJkzQ3Nzdt8+bNd5wrLi5OA+QmN7nJTW5yk5sD3uLi4m77vd6gabYp2zQYDLc9gvPaa6+xaNEiYmJiCu8bPnw4u3fvJjo6GoDBgweTkpLC0qVLC/fp3r07FSpUYM6cOXeUxWKxUL58eeLi4vD19b27FySEEEIIm0pJSSEoKIjk5GTMZvMt97XqKaqSio6OJjIyssh9Dz74IDNmzCAnJwd3d3eio6N55ZVXiu0zbdq0mz5vVlYWWVlZhb9OTU0FwNfXVwYcIYQQwsHcyeUldnWRcUJCAv7+/kXu8/f3Jzc3l6SkpFvuk5CQcNPnnTx5MmazufAWFBRU+uGFEEIIYTfsasCB4lPZlTNo195/o31uNc2NHTsWi8VSeIuLiyvFxEIIIYSwN3Z1iiogIKDYkZjExETc3NyoVKnSLfe5/qjOtTw9PfH09Cz9wEIIIYSwS3Z1BCc8PJwVK1YUue+PP/6gVatWuLu733KfiIgIm+UUQgghhH2z6hGctLQ0jh49WvjrEydOsGvXLipWrEiNGjUYO3YsZ86c4ccffwT0T0x9/vnnjBw5kmHDhhEdHc2MGTOKfDrq5ZdfpkOHDkyZMoW+ffvy66+/snLlSjZs2GDNlyKEEEIIB2LVIzjbtm2jRYsWtGjRAoCRI0fSokUL3nrrLQDi4+OJjY0t3D84OJglS5awZs0amjdvzjvvvMOnn37KwIEDC/eJiIhg7ty5fP/994SGhjJr1izmzZtHmzZtrPlShBBCCOFAbLYOjj1JSUnBbDZjsVjkY+JCCCGEgyjJ92+7ugZHCCGEEKI0yIAjhBBCCKcjA44QQgghnI4MOEIIIYRwOjLgCCGEEMLp2NVKxsLx7TtjYeHOM+S73ofzChkNBvo0q0qzoPKqowghhMuSAUeUGktGDkN/+ItzKVm339nJ/brrDCte6UgFHw/VUYQQwiXJgCNKzaTFMZxLyaJGRW96NwtUHUeZJXsTOJF0mYm/H+Djwc1VxxFCCJckA44oFeuPnGfetjgMBvjwkWbcV6ui6kjKdG3kz8AvN7Fw5xl6NwukS8ObF8EKIYSwDrnIWNyzy1m5jInaC8ATbWu69HAD0KJGBYbeHwzA6wv2kZKZoziREEK4HhlwxD37YPkhziRnUK18GUZ3b6g6jl0Y2a0BtSp5k5CSyeQlB1XHEUIIlyMDjrgnW09cZNamkwC8PzAEH0856wlQxsPElIGhAMzZGsvGo0mKEwkhhGuRAUfctcycPF6L2gPA4FZBtK9XWXEi+9KmdiWGtK0JwJgFe0jPzlWcSAghXIcMOOKufbzyMCeSLuPv68nrPRupjmOXXuvRkGrlyxB3MYMPlh9SHUcIIVyGDDjiruyOS+bbdccBeLdfCOYy7ooT2aeynm5MGhACwKxNJ9l28qLiREII4RpkwBEllp2bz+j5e8jXoE+zqnRrLB+DvpWO9SvzcFh1NA1GR+0hMydPdSQhhHB6MuCIEvti9VEOnUulko8H4/s0UR3HIbzRszGVy3ly/PxlPvnziOo4Qgjh9GTAESUSE5/CF6uPAjChbxMqShXBHTF7u/Nev6YAfLPuOHtPWxQnEkII5yYDjrhjuXn6qancfI3Ixv70DHHdOoa7EdkkgF6hgeTla7w6fzfZufmqIwkhhNOSAUfcsW/Xn2DvGQu+Xm68268pBoNBdSSHM6GPftTrYEIqX645pjqOEEI4LRlwxB05dj6Nj1ceBuDNXo2p4uulOJFjqlTWk7d7Nwbg89VHOJSQqjiREEI4JxlwSlteLmRfVp2iVOXna7w2fw/Zufl0qF+ZQWHVVUdyaH2aVaVrI39y8jRGz99Nbp6cqhJCOJlM9dcZyoBTmpJjYVZPWPQiaJrqNKXmx+iTbDt1CR8PE5MHhMipqXtkMBh4r39Tynm5sfu0hZkbT6iOJIQQpSf9IkwPh2WvQ26Wshgy4JSm1HNw+i/YFwW7ZqtOUyriLqYzZZm+Au+YhxpRrXwZxYmcg7+vF2/21E9VffjHYY6fT1OcSAghSoGmwa//gpQzcGQ55OUoiyIDTmkKug+6jNO3l7wKSY693ommaYxZsIeMnDzaBFfk8dY1VEdyKg+3qk77en5k5eYzJmov+fnOc9RPCOGi/voODi0BkwcMmgmeZZVFkQGntLUbAcEdICcd5v9D6eG5ezXvrzg2Hr2Al7uRKQNDMRrl1FRpMhgMTOofgreHia0nL/KfLadURxJCiLuXsA+WF/yQ33UCBDZTGkcGnNJmNEH/b8C7EiTshZXjVSe6K/GWDN5bHAPAqMgG1PLzUZzIOQVV9GZMj4YAvL/0IHEX0xUnEkKIu5CdDvOfhrwsqPcgtH1edSIZcKzCNxD6falvb54Oh5erzVNCmqbxxsJ9pGbl0jyoPP9oF6w6klP7e5uatK5VkfTsPF5fuBfNiS5QF0K4iOVjIekQlPWHftPBDj6MIgOOtdR/ENoUTLC/PA+pCWrzlMCvu87y58FEPExGpg4KxSSnpqzKaDTw/sAQPN2MrD+SxM/bTquOJIQQd27/L7B9FmCAAd+Aj5/iQDoZcKyp2wQICIH0C7DgWci3//VOzqdmMf63/QC82KUu9f3LKU7kGmpXLsvIbvUBeGfxAc6lZCpOJIQQdyA5Fn57Sd++fwTU7qQyTREy4FiTmycMnAnu3nBiLWycpjrRbY1ftJ/k9BwaB/oyvFMd1XFcytD7g2lW3UxqZi7jFu6TU1VCCPuWlwtRw/RF/aqFQedxqhMVIQOOtVWuDz2m6tur3oW4v9TmuYVl++JZvDcek9HA1EGhuJvkr4ctuZmMTB3UDHeTgZUx5/htT7zqSEIIcXNrp0DcZvAoBwNngMlddaIi5DuYLbT4OzQZAFoeRA21iyWsr5ecns0bv+inpoZ3rE3TambFiVxTg4ByvNC5HqAfTbuQ5rjLDAghnNiJ9bDuA3279zSoaH8fRpEBxxYMBv0vQPkakHwKfn/F7qocJv5+gKS0LOpWKcuLXeqpjuPSnu9Uh4YB5bh4OZu3F+1XHUcIIYpKv6hfV4oGzR+HkEGqE92QTQac6dOnExwcjJeXF2FhYaxfv/6m+z711FMYDIZityZNmhTuM2vWrBvuk5lpxxdmepn163EMpoIqh/+qTlRo9aFEFuw4g8EAUwaG4uVuUh3JpXm4GflgUDNMRgO/74ln+X7H+QSeEMLJaRr8+gKknoVKda9egmGHrD7gzJs3jxEjRjBu3Dh27txJ+/bt6dGjB7GxsTfc/5NPPiE+Pr7wFhcXR8WKFXn44YeL7Ofr61tkv/j4eLy8vKz9cu6NHVY5pGbm8PqCvQA83S6YsJoVFCcSACHVzTzboTYAb/yyD0u6uj4XIYQo9Nd3cGixXVQx3I7VB5yPPvqIoUOH8swzz9CoUSOmTZtGUFAQX3755Q33N5vNBAQEFN62bdvGpUuX+Mc//lFkP4PBUGS/gIAAa7+U0lFY5XBZX/VRcZXD+0sPEm/JpEZFb0ZFNlCaRRT18gP1qF3Zh/OpWby7+IDqOEIIV3duv11VMdyOVQec7Oxstm/fTmRkZJH7IyMj2bRp0x09x4wZM+jatSs1a9Yscn9aWho1a9akevXq9OrVi507d970ObKyskhJSSlyU+ZKlUOZipCwB1ZOUBYl+tgFZm/Rj6S9PzCEMh5yasqeeLmbmDowFIMBft5+mrWHz6uOJIRwVUWqGCLtoorhdqw64CQlJZGXl4e/v3+R+/39/UlIuP11BfHx8SxdupRnnnmmyP0NGzZk1qxZLFq0iDlz5uDl5UW7du04cuTGp3wmT56M2WwuvAUFBd39iyoNRaocvoDDf9g8Qnp2Lq9F7QHgsTY1iKhjHytPiqJa1arIk+G1AHh9wV7SsnLVBhJCuKblr8P5g3oVQ1/7qGK4HZtcZGy47g9C07Ri993IrFmzKF++PP369Styf9u2bfn73/9Os2bNaN++Pf/73/+oX78+n3322Q2fZ+zYsVgslsJbXFzcXb+WUtOgO7QZrm//MtzmVQ4f/nGY2IvpBJq9GFtQ9ijs0+juDQiqWIYzyRlMWXpQdRwhhKs58Cts/x4wQP+voWxl1YnuiFUHHD8/P0wmU7GjNYmJicWO6lxP0zRmzpzJkCFD8PDwuOW+RqOR++6776ZHcDw9PfH19S1yswtdJ4C/7asctp+6xMyNJwCYNCCEcl72tTiTKMrbw433B4QC8NPmU2w+fkFxIiGEy0iOhUUv6tvtXoY6ndXmKQGrDjgeHh6EhYWxYsWKIvevWLGCiIiIWz527dq1HD16lKFDh97262iaxq5duwgMDLynvDbn7qVfhX6lymHTJ1b/kpk5ebwWtQdNgwEtq9G5QRWrf01x79rV9eNvrfVTq2Oi9pCRnac4kRDC6V1fxdDlDdWJSsTqp6hGjhzJd999x8yZM4mJieGVV14hNjaW4cP10zNjx47liSeeKPa4GTNm0KZNG5o2bVrs9yZMmMDy5cs5fvw4u3btYujQoezatavwOR1K5frQY4q+vepdOL3Nql/us1VHOJqYhl9ZT97q1diqX0uUrrEPNSLA14uTF9L5aMUh1XGEEM5u3VS7rmK4HasPOIMHD2batGlMnDiR5s2bs27dOpYsWVL4qaj4+Phia+JYLBaioqJuevQmOTmZZ599lkaNGhEZGcmZM2dYt24drVu3tvbLsY4WQ/Qqh/xc/Sp1K1U57Dtj4au1xwF4t18Tynvf+tSfsC++Xu5MGqAP/DM2nGBn7CXFiYQQTuvkBruvYrgdg+aClcUpKSmYzWYsFov9XI+TkQxft9fPdzYdBAO/K9Wr1HPy8un7+UYOxKfQMySQLx5vWWrPLWzrlXm7WLjzDPWqlOX3l+7H000+3i+EKEXpF+Gr+yHljF7F0G+66kSFSvL9W7qo7EWZ8vohQIMJ9s2H3XNK9em/XnuMA/EplPd2Z3yfJrd/gLBbb/VqjF9ZD44kpvHFqqOq4wghnImm6RcVp5yBinXsuorhdmTAsSdBraHz6/r24lGQVDrfvI6cS+XTP/XnGt+7CZXLeZbK8wo1Kvh4MLGvfqpq+ppjHDircOFKIYRz2TYDDv4ORne7r2K4HRlw7M39r0Ct9gVVDv+45yqHvHyNV+fvITsvny4Nq9C3edVSCipUeigkkB5NA8jN13h1/m5y8myzxIAQwomd2w/LCn7I7jYBqjZXGudeyYBjb4wmGFB6VQ7fbzzBrrhkynm68V7/pne0wKJwDBP6NsFcxp39Z1P4Zt1x1XGEEI4sOx3mD9WrGOp2gzb2X8VwOzLg2CPfqlcv6rqHKoeTSZf5YLn+ceJxPRsRaC5TWgmFHahSzou3e+sf9f9k5RGOJqYqTiSEcFjLX4fzMeBTRa8SMjr+eOD4r8BZNegBrZ/Tt395vsRVDvn5Gq9F7SErN592dSsx+D7F/VvCKvq3qEanBpXJzstn9Pw95OW73IcihRD3qrCKARjgOFUMtyMDjj3rNrGgyiEJFj5XoiqH2Vtj2XLiImXcTbw/IFROTTkpg8HApP4hlPV0Y0dsMrM2nVQdSQjhSJLjrqliGAF1uiiNU5pkwLFn11Y5HF8Dmz69o4edSc7g/SUxwJWiRm8rhhSqVS1fhrEP6YWpHyw/yKkLlxUnEkI4hLxcWOC4VQy3IwOOvStS5fAOnN5+y901TWPsgr1czs6jVc0KPBley/oZhXJ/u68G4bUrkZmTz5iovbjg+p1CiJJa9wHERhdUMXzncFUMtyMDjiNoMQSa9NerHKKehsybr3sSteMM6w6fx8PNyJRBoRiNcmrKFRiNBt4fGEIZdxPRxy8wZ2uc6khCCHt2cqPeNQXQ62OoWFttHiuQAccRGAzQaxqYa8Clk/D7K/pqk9dJTMlk4m/7AXila33qVHbcBZpEydWs5MOoBxsAMGlJDGeTMxQnEkLYpfSL+qkpLR+aPQahD6tOZBUy4DiKMuUL+qluXOWgaRpv/LKPlMxcQqqZGdbe8YrRxL17KqIWLWuUJy0rl9cXyqkqIcR1rq9ieMhxqxhuRwYcR1KjDXQeq29fV+WweG88fxw4h5vRwNRBobiZ5K11RaaC99/DZGTNofMs3HlGdSQhhD3ZNvOaKoYZ4FlOdSKrke+Cjub+kcWqHC5ezubtX/VTU//qXJdGgXbSkC6UqFulHC93rQfAhN8OkJiaqTiREMIunDugL+gH0HU8VG2hNI61yYDjaK6vcvhzIhN+28+Fy9k08C/HvzrXVZ1Q2IFnO9SmSVVfLBk5hcOvEMKFZafD/KchNxPqdoW2/1SdyOpkwHFE11Y5RH+OZc8SjAb0UxNu8pYKcDcZ9VOVRgNL9yWwZG+86khCCJX+GHdNFcNXTlHFcDvO/wqdVYMeZIUNA+D/3L9iRBszzYLKq80k7EqTqmae71QHgLd+3cely9mKEwkhlDiwSL/2BpyqiuF2ZMBxYBMzBhOTXwM/Qwr/skwtUZWDcA0vdKlLvSplSUrLZuLvB1THEULYWnIcLHpB3273slNVMdyODDgOav2R88zekcgLOS+SZ/LCdGLtHVc5CNfh6WZi6qBQjAZYuPMMqw6eUx1JCGErebmw4Fm9iqFqS+jsXFUMtyMDjgO6nJXLmKi9ANzfNgLTlXUM7qDKQbieFjUqMPR+fV2k1xfsIyUzR3EiIYRNrP8/iN2kVzEMmgFuHqoT2ZQMOA5o6rKDnEnOoFr5Mozu3hBaPgGN+91RlYNwTSO7NaBWJW8SUjKZXFDEKoRwYqc2wdqCHsNeHzllFcPtyIDjYLaeuMgP0acAeH9gCD6ebnqVQ+9PrlY5LB55wyoH4brKeJiYMjAUgDlb49h4NElxIiGE1aRfhKgrVQx/g9BHVCdSQgYcB5KZk8drUXsAGNwqiPb1rrkS/toqh70/w+65akIKu9WmdiWGtK0JwJgFe0jPzlWcSAhR6gqrGE7rR20e+kB1ImVkwHEgH688zImky/j7evJ6z0bFdyhS5fDvIlUOQgC81qMh1cqXIe5iBh8sP6Q6jhCitBWpYpjp1FUMtyMDjoPYHZfMt+uOA/BuvxDMZdxvvOO1VQ5RT0OurH0irirr6cakASEAzNp0km0nLypOJIQoNS5WxXA7MuA4gOzcfEbP30O+Bn2aVaVbY/+b71xY5VAB4nfDnxNsF1Q4hI71K/NwWHU0DUZH7SEzJ091JCHEvcrJcLkqhtuRAccBfLH6KIfOpVLJx4PxfZrc/gG+VaHv1SoHjqywbkDhcN7o2ZjK5Tw5fv4yn/x5RHUcIcS9Wn5tFcOXLlHFcDvyJ2DnYuJT+GK1fi3N+D5NqOhzh+sYNHwIWj+rby8cDqmywJu4yuztznv9mgLwzbrj7D1tUZxICHHXYn6DbTP07f5fQdkqavPYCRlw7Fhunn5qKjdfI7KxP71CA0v2BN3egSpNID0JFj4nVQ6iiMgmAfQKDSQvX+PV+bvJzpW/H0I4HMtp+LWgiiHiJaj7gNo8dkQGHDv27foT7D1jwdfLjXf7NcVgMJTsCdy99Kvo3crA8dUQ/Zl1ggqHNaHgqODBhFS+XHNMdRwhREnk5+nr3WQm6xcUd3lTdSK7IgOOnTp2Po2PVx4G4M1ejani63V3T1SlIfR4X9/+cyKckSoHcVWlsp683bsxAJ+vPsKhhFTFiYQQd2zdBwVVDGVhoOtVMdyODDh2KC9fY/T8PWTn5tOhfmUGhVW/tyds+eTVKof5Q6XKQRTRp1lVujbyJydPY/T83eTmyakqIezetVUMPT+CSnXU5rFDMuDYoR+jT7L91CV8PExMHhBS8lNT1ytS5XAClowqnaDCKRgMBt7r35RyXm7sPm1hxoYTqiMJIW7l+iqGZoNVJ7JLMuDYmbiL6Uxdpq8wO+ahRlQrX6Z0nvjaKoc982DXnNJ5XuEU/H29eLOnfqrqoxWHOX4+TXEiIcQNaRr89pJUMdwBGXDsiKZpjFmwh4ycPFoHV+Tx1jVK9wvUaAOdpMpB3NjDrarTvp4fWbn5vBa1h/x8KWwVwu5s/17/WLhUMdyWTQac6dOnExwcjJeXF2FhYaxfv/6m+65ZswaDwVDsdvDgwSL7RUVF0bhxYzw9PWncuDELFy609suwunl/xbHx6AW83I1MHRiK0XiPp6ZupP1IqHm/VDmIYgwGA5P6h+DtYeKvk5f4z5ZTqiMJIa6VGAPLCn5I7fq2y1cx3I7VB5x58+YxYsQIxo0bx86dO2nfvj09evQgNjb2lo87dOgQ8fHxhbd69eoV/l50dDSDBw9myJAh7N69myFDhvDII4+wZcsWa78cq4m3ZPDe4hgA/t2tAbX8fKzzhaTKQdxCUEVvxvRoCMD7Sw8SdzFdcSIhBFC0iqHOA9D2X6oT2T2DpmlWPQ7dpk0bWrZsyZdffll4X6NGjejXrx+TJ08utv+aNWvo3Lkzly5donz58jd8zsGDB5OSksLSpUsL7+vevTsVKlRgzpzbX1uSkpKC2WzGYrHg6+tb8hdVyjRN45kftvHnwUSaBZVnwfMRmKxx9OZaBxfD3Mf07cejoF5X63494TDy8zUe/WYzW09epH09P358uvW9X+guhLg3i/8Nf30HPpXh+U0uu1pxSb5/W/UITnZ2Ntu3bycyMrLI/ZGRkWzatOmWj23RogWBgYE88MADrF69usjvRUdHF3vOBx988KbPmZWVRUpKSpGbPfl111n+PJiIh8nIB4NCrT/cADTsCfcN07d/kSoHcZXRaOD9gSF4uhlZfySJn7edVh1JCNcW87s+3IBUMZSAVQecpKQk8vLy8Pcv2n7t7+9PQkLCDR8TGBjIN998Q1RUFAsWLKBBgwY88MADrFu3rnCfhISEEj3n5MmTMZvNhbegoKB7fGWl53xqFuN/2w/Ai13qUt/fhheMRRZUOVw+rw85UuUgCtSuXJaR3eoD8M7iA5xLyVScSAgXZTkNvxacjop4UW8KF3fEJhcZX394W9O0mx7ybtCgAcOGDaNly5aEh4czffp0evbsyf/93//d9XOOHTsWi8VSeIuLi7uHV1O6xi/aT3J6Do0DfRneycYLNbmXuVrlcGyV3jwuRIGh9wfTrLqZ1Mxcxi3ci5XPZgshrlesiuEt1YkcilUHHD8/P0wmU7EjK4mJicWOwNxK27ZtOXLkSOGvAwICSvScnp6e+Pr6FrnZg2X74lm8Nx6T0cDUQaG4mxR8ar9IlcMEqXIQhdxMRqYOaoa7ycDKmEQW7T6rOpIQrmXd/0kVwz2w6ndUDw8PwsLCWLFiRZH7V6xYQURExB0/z86dOwkMvNqkHR4eXuw5//jjjxI9p2rJ6dm88Yt+amp4x9o0rWZWF6blk9C4r1Q5iGIaBJTjhc76Jxgn/HaAC2lZihMJ4SJORcPagh8+pYrhrrhZ+wuMHDmSIUOG0KpVK8LDw/nmm2+IjY1l+PDhgH766MyZM/z4448ATJs2jVq1atGkSROys7P5z3/+Q1RUFFFRUYXP+fLLL9OhQwemTJlC3759+fXXX1m5ciUbNmyw9sspNRN/P0BSWhZ1KvvwYpd6t3+ANV2pcjiz42qVw4Bv1GYSduP5TnVYui+egwmpvL1oP58/1lJ1JCGcW8YliHpGr2IIfVSqGO6S1c+JDB48mGnTpjFx4kSaN2/OunXrWLJkCTVr1gQgPj6+yJo42dnZjBo1itDQUNq3b8+GDRtYvHgxAwYMKNwnIiKCuXPn8v333xMaGsqsWbOYN28ebdq0sfbLKRWrDyWyYMcZDAaYOqgZXu4m1ZH0dXEGfgcGo17lsHuu6kTCTni4GflgUDNMRgO/74ln+f4bX8wvhCgFmgaLrqli6Pl/t3+MuCGrr4Njj1Sug5OamUPkx+uIt2TydLtg3urd2KZf/7bWToXV7+nnfJ9bJ4dFRaEpyw7y5ZpjVC7nycpXOmL2dlcdSQjns+17+H0EGN1g6AqoJkdMr2U36+CI4t5fepB4SyY1Knoz6sH6quMU1/7fULMdZKcVrJopVQ5C9/ID9ahd2YfzqVm8u/iA6jhCOJ/EGFg2Rt9+4G0Zbu6RDDg2tOlYErO36Kfj3h8YgreH1S+BKjmjCQZ8W1DlsAtWTVSdSNgJL3cTUweGYjDAz9tPs/bwedWRhHAeORn6hzxyM6FOFwh/QXUihycDjo2kZ+cyJmovAI+1qUFEHT/FiW7BXA36fqFvb/oMjq5Um0fYjVa1KvJURC0AXl+wl7SsXLWBhHAWf7wBifv1KoZ+X4FRvj3fK/kTtJEP/zhM7MV0As1ejC0oM7Rr11Y5LJQqB3HVqw82IKhiGc4kZzBl6UHVcYRwfNdXMZS783XixM3JgGMD209dYubGEwBMGhBCOS8HuThTqhzEDXh7uPH+gFAAftp8is3HLyhOJIQDkyoGq5EBx8oyc/J4LWoPmgYDWlajcwMHKkmTKgdxE+3q+vG31nqn25ioPWRk5ylOJIQDys+DBc/qVQyBzaWKoZTJgGNln606wtHENPzKevJWLzv7SPidqNIQuk/Wt6XKQVxj7EONCPD14uSFdD5acUh1HCEcz/oP4dRGfVmOQTOliqGUyYBjRfvOWPhq7XEA3u3XhPLeDvqXN+wpaNTnapVDVqrqRMIO+Hq5M2lAUwBmbDjBzthLihMJ4UBiN8Oagh8ee34oa45ZgQw4VpKTl8+r8/eQl6/RMySQ7k0Db/8ge2UwQJ9Pwbe6XuWweJTqRMJOdGnoT/8W1cjXYPT8PWTlyqkqIW6rSBXDYGj2qOpETkkGHCv5as0xYuJTKO/tzvg+TVTHuXdFqhzmSpWDKPRWr8b4lfXgSGIan686qjqOEPbtShWDJQ4qBMNDUsVgLTLgWMHhc6l8VvAP/fjeTahczlNxolJSMxw6FqyyufjfcOGY2jzCLlTw8WBiX/1U1ZdrjrH/rEVxIiHs2I4fIGaRXsUwaAZ42bYuyJXIgFPK8vI1Rs/fQ3ZePl0aVqFv86qqI5WuDqOuVjlEDZUqBwHAQyGB9GgaQG7B3/+cPFlSQIhiEg/C0murGMLU5nFyMuCUsu83nmBXXDLlPN14r39TDAaD6kil69oqh7M7pcpBFJrQtwnmMu7sP5vCN+uOq44jhH3JySjo98uQKgYbkQGnFJ1MuswHy/WPy77esxGB5jKKE1mJuRr0KVgTR6ocRIEq5bx4u7e+FMInK48QdzFdcSIh7Mgfb0oVg43Jn3ApSsvKpVr5MrSrW4lH7wtSHce6GvWC+57RtxcOh7REtXmEXejfohrhtSuRnZfPjA0nVMcRwj4cXAx/fatv95MqBluRAacUNa1mZsnL7Zk2uIXznZq6kch3oUpjvcphoVQ5CDAYDPyrc10A5v4Vy8XLco2WcHGWM1erGMJfgHpSxWArMuCUMi93k/N8aup2CqscvODYn7D5C9WJhB1oV7cSTar6kpmTz4/RJ1XHEUKdK1UMGZf0KoYH3ladyKXIgCPuTZVGV6scVk6AMzvU5hHKGQwGhnfUV2X9YdNJ6akSrmv9R3BqA7j7SBWDAjLgiHsX9g9o1Bvyc/SPjkuVg8vr0TSAoIpluJSew/+2xamOI4TtSRWDcjLgiHtnMEDvgiqHi8elykHgZjLybPvaAHy7/ji5si6OcCWFVQx5EPKIVDEoIgOOKB3eFYtWORxcojqRUGxQWBAVfTw4fSmDxXvjVccRwnb+ePNqFUPPD/UfAoXNyYAjSk/NcIh4Ud9e9a58qsrFlfEw8VRELQC+XnscTdPUBhLCFpKOwq7Z+nb/r6SKQSEZcETpajcCPH31Ba0O/KI6jVBsSNualHE3cSA+hfVHklTHEcL61r6vt4TX7w412qpO49JkwBGly7vi1SXIV0+CvFy1eYRSFXw8eLS1vujl1+uknFU4uXMHYO98fbvz62qzCBlwhBW0fV7vqrpwBPb+rDqNUGzo/cGYjAY2Hr3A3tPSNC6c2JpJgAaN+0JgM9VpXJ4MOKL0eflCu5f17bXvQ16O2jxCqeoVvOnTrCoAX8lRHOGszu6CmN8AA3QaqzqNQAYcYS2tn9VL5S6dvHrBnXBZz3bQPzK+dG88py5cVpxGCCtYPUn/b8jD+gKoQjkZcIR1ePjA/SP17bUfQG6W2jxCqUaBvnRqUJl8TV8XRwinEvcXHFkOBhN0GqM6jSggA46wnlZPQ7mqkHIatv+gOo1Q7LkO+kquP287TVKaDLzCiax+V/9v87/JisV2RAYcYT3uXtChYFXj9f8H2elq8wil2tauSLOg8mTl5vPDppOq4whROk5ugONrwOgOHUarTiOuIQOOsK4WQ6B8DUg7B9tmqE4jFDIYDAwvuBbnx+hTXM6SJQSEg9M0WPWevt3yCahQU20eUYQMOMK63Dyg42v69oaPpYjTxUU2CSDYzwdLRg5z/5ISTuHgjq2C2E1g8rx6tFrYDRlwhPWFPgoV60D6Bdjyteo0QiGT0cCwghLOGeuPkyMlnMJRaRqsLjh6c99Q8K2qNo8oRgYcYX0mt6vrQmz6FDKSlcYRag1oWQ2/sp6ctWTy2+6zquMIcXcOL4Mz28HdG+5/RXUacQM2GXCmT59OcHAwXl5ehIWFsX79+pvuu2DBArp160blypXx9fUlPDyc5cuXF9ln1qxZGAyGYrfMzExrvxRxt5oOgMqNINMC0V+oTiMU8nI38Y92tQAp4RQOKj//6rU3bZ6DslXU5hE3ZPUBZ968eYwYMYJx48axc+dO2rdvT48ePYiNjb3h/uvWraNbt24sWbKE7du307lzZ3r37s3OnTuL7Ofr60t8fHyRm5eXl7VfjrhbRhN0LjiKs/lLuHxBbR6h1N/b1MTHw8Shc6msOXRedRwhSiZmEZzbCx7lIOIl1WnETVh9wPnoo48YOnQozzzzDI0aNWLatGkEBQXx5Zdf3nD/adOmMXr0aO677z7q1avHpEmTqFevHr/99luR/QwGAwEBAUVuws417A0BoZCdCps+UZ1GKGT2duexNjUA+Gqt1DcIB5Kfd3XV4vB/6QXDwi5ZdcDJzs5m+/btREZGFrk/MjKSTZs23dFz5Ofnk5qaSsWKRf8SpaWlUbNmTapXr06vXr2KHeG5VlZWFikpKUVuQgGjETqP07e3fAOp59TmEUo9fX8w7iYDW05cZGfsJdVxhLgze+dD0iHwKg/h/1SdRtyCVQecpKQk8vLy8Pf3L3K/v78/CQkJd/QcH374IZcvX+aRRx4pvK9hw4bMmjWLRYsWMWfOHLy8vGjXrh1Hjhy54XNMnjwZs9lceAsKCrr7FyXuTf0HoVoryM3QPzYuXFaguQx9m1cD9GtxhLB7eTl6gTBAu5fAy6w2j7glm1xkbDAYivxa07Ri993InDlzGD9+PPPmzaNKlasXcbVt25a///3vNGvWjPbt2/O///2P+vXr89lnn93wecaOHYvFYim8xcXJ+hvKGAzQ5Q19e9sMsJxRm0co9VzBwn/LDyRw7Hya4jRC3MbuOXDxOHj7QevnVKcRt2HVAcfPzw+TyVTsaE1iYmKxozrXmzdvHkOHDuV///sfXbt2veW+RqOR++6776ZHcDw9PfH19S1yEwrV7gQ120Fetl7hIFxWPf9ydG1UBU2D76SEU9iz3CxYO1Xfvv8V8CyrNo+4LasOOB4eHoSFhbFixYoi969YsYKIiIibPm7OnDk89dRT/Pe//6Vnz563/TqaprFr1y4CAwPvObOwAYPh6rU4O36ESyeVxhFqPddRLyeM2n6GxBRZ6kHYqR0/giUOygboC/sJu2f1U1QjR47ku+++Y+bMmcTExPDKK68QGxvL8OHDAf300RNPPFG4/5w5c3jiiSf48MMPadu2LQkJCSQkJGCxWAr3mTBhAsuXL+f48ePs2rWLoUOHsmvXrsLnFA6gVjuo3Rnyc6/+VCRc0n21KhJWswLZefl8LyWcwh7lZMC6gqPNHUaBexm1ecQdsfqAM3jwYKZNm8bEiRNp3rw569atY8mSJdSsqZeSxcfHF1kT5+uvvyY3N5d//etfBAYGFt5efvnlwn2Sk5N59tlnadSoEZGRkZw5c4Z169bRunVra78cUZquXIuzew4kHVWbRSh15Vqc/2w+RWpmjuI0Qlxn20xISwBzkF6qKRyCQXPBZURTUlIwm81YLBa5Hke1/z4Kh5dC00EwSNrGXVV+vka3j9dy7PxlXn+oIc92qKM6khC6rDT4pBmkJ0Gfz2TAUawk37+li0qo1fl1/b/7ouDcAbVZhDJGo4HnCoaaGRtOkJWbpziREAW2fqMPNxWCodnfVKcRJSADjlArMBQa9wU0WDNJdRqhUN8WVfH39eRcSha/7pISTmEHMi2wsWDV9U5jweSuNo8oERlwhHqdxgIGiPkNzu5SnUYo4ulm4ul2wQB8vfYY+fkud/Zc2Jvo6ZCZDH4NIGSQ6jSihGTAEepVaQQhD+vbq+Uojit7rE0Nynm6cez8Zf48mKg6jnBl6Rdh83R9u/NYvTBYOBQZcIR96DQGDCY4shzitqpOIxQp5+XO4231T1hKCadQatOnkJUC/iHQqK/qNOIuyIAj7EOlOtC84AK+1e+pzSKUerpdLTxMRrafusS2kxdVxxGuKO08bPla3+78ul4ULByOvGvCfnQYDUZ3OL4GTqxXnUYoUsXXiwEt9RJOOYojlNjwMeSkQ9WW0KCH6jTiLsmAI+xHhZoQ9qS+vfo9cL0lmkSBYR1qYzDAyphEjpxLVR1HuJKUs/DXd/p2lzf0ahnhkGTAEfal/b/B5Amx0XBsleo0QpE6lcsS2Vgv5P16nZRwChta/yHkZUGNcKjTRXUacQ9kwBH2xbcq3PeMvr3qXTmK48KGF5Rw/rrrDPGWDMVphEtIjoXtP+jbcvTG4cmAI+zP/a+Auzec3QGHl6lOIxRpUaMCrYMrkpOnMXPDCdVxhCtYOxXycyC4I9S6X3UacY9kwBH2p2xlaPOcvr3qPcjPV5tHKPN8wVGc/26JxZIhJZzCii4cg13/1bevFAELhyYDjrBPES+Bpy+c2wsxv6pOIxTp1KAyDfzLcTk7j/9sPqU6jnBma94HLQ/qPQhBrVWnEaVABhxhn7wrQtt/6turJ0O+lC+6IoPBwHMdawPw/caTZObI3wNhBYkHYe/P+vaVAmDh8GTAEfYr/J/gVR6SDsHe+arTCEV6N6tKVbMXSWlZLNhxRnUc4YzWTAI0aNQbqjZXnUaUEhlwhP3yMkO7l/Ttte9DnlyD4YrcTUaGtteP4ny7/jh5UsIpSlP8HjjwK2CATnL0xpnIgCPsW+vnwNsPLh6H3XNUpxGKPHpfEOYy7pxIuswf+xNUxxHO5ErBb9OB4N9YbRZRqmTAEfbNsyy0H6lvr50KuVlq8wglfDzdeCL8agmnJusjidJwehscXgoGI3QaqzqNKGUy4Aj71+ppKBcIljjY8aPqNEKRJyNq4elmZPdpC5uPSwmnKAVXin2b/Q386qrNIkqdDDjC/rmX0SscANb9H+TIqrauyK+sJw+3qg7A1+ukhFPco1Ob9DoYoxt0HK06jbACGXCEY2j5BJiDIC0Bts1UnUYoMqx9bYwGWHPoPDHxKarjCEelaXoVDECLIVChltI4wjpkwBGOwc3z6k9Z6z+CrDS1eYQSNSv50CMkEIBvpIRT3K3ja+DURr3Yt8OrqtMIK5EBRziOZn+DirUhPQm2fq06jVBkeAe9vmHR7rOcvpSuOI1wONcevWn1NJirqc0jrEYGHOE4TO7QcYy+vfFTyLSozSOUCKlupl3dSuTla8yQEk5RUkf+gDPbwK2MXuwrnJYMOMKxhAwCvwaQmQzR01WnEYo8V3AUZ+7WOC5dzlacRjiM/PyrR2/aPAvl/NXmEVYlA45wLEYTdC5YryL6C0iXjwu7ovb1/Ggc6EtGTh4/SQmnuFMHf4OEPeBRFiJeVp1GWJkMOMLxNOoL/iGQnQqbPlWdRihwbQnnrE1SwinuQH6eXtwLepGvTyW1eYTVyYAjHI/RCF3G6dtbvoa0RLV5hBI9QwKpXqEMFy9n8/O2ONVxhL3btwDOx+gdd+H/Up1G2IAMOMIx1e8O1cIgJx02TFOdRijgZjIyrLCE8wS5efmKEwm7lZcLawqO3kS8CGXKK40jbEMGHOGYDAboXHAU56/vIOWs2jxCiYdbVaeCtzuxF9NZuk9KOMVN7JkLF4+BdyVoM1x1GmEjMuAIx1WnC9QIh7wsvcJBuBxvDzeejKgF6PUNUsIpisnNhjVT9O12I8CznNI4wnZkwBGOy2CALm/o2zt+hEvyaRpX9ER4Lbzcjew7k8LGoxdUxxH2ZudPYImFsv5w3zOq0wgbkgFHOLZa90PtTpCfA+umqk4jFKjo48Gj99UApIRTXCcn8+rR3fajwMNbbR5hUzLgCMfXueAozq45cEG+wbmiofcHYzIaWH8kiX1nZIVrUWD795B6FnyrQ9iTqtMIG5MBRzi+oPug3oOg5cGa91WnEQoEVfSmV6hewvm1lHAKgOzLsP5Dfbvjq3phr3ApNhlwpk+fTnBwMF5eXoSFhbF+/fpb7r927VrCwsLw8vKidu3afPXVV8X2iYqKonHjxnh6etK4cWMWLlxorfjCEXR+Xf/v3p8hMUZtFqHEsx30j4wv3nOW2AtSwunytn4Dl89DhVrQ/HHVaYQCVh9w5s2bx4gRIxg3bhw7d+6kffv29OjRg9jY2Bvuf+LECR566CHat2/Pzp07ef3113nppZeIiooq3Cc6OprBgwczZMgQdu/ezZAhQ3jkkUfYsmWLtV+OsFdVm0Oj3oB2db0L4VKaVDXToX5l8jX4boMcxXFpmSmw8RN9u+MYvahXuByDZuXPVbZp04aWLVvy5ZdfFt7XqFEj+vXrx+TJxb8RvfbaayxatIiYmKs/hQ8fPpzdu3cTHR0NwODBg0lJSWHp0qWF+3Tv3p0KFSowZ86cYs+ZlZVFVlZW4a9TUlIICgrCYrHg6+tbKq9T2IFzB+DLCECD59ZDYKjqRMLGNh1N4rHvtuDlbmTja12oVFZOS7ikNVNgzSTwqw//3Kx32AmnkJKSgtlsvqPv31Y9gpOdnc327duJjIwscn9kZCSbNm264WOio6OL7f/ggw+ybds2cnJybrnPzZ5z8uTJmM3mwltQUNDdviRhz/wbQ9OB+vbqSWqzCCXC61QitLqZzJx8foiWZQNcUvpFiP5c3+40RoYbF2bVAScpKYm8vDz8/YtW0vv7+5OQcONVRxMSEm64f25uLklJSbfc52bPOXbsWCwWS+EtLk56a5xWp7FgMMLhpXB6m+o0wsYMBgPPdagDwI/RJ0nPzlWcSNhc9OeQlQJVmkDj/qrTCIVscpGxwWAo8mtN04rdd7v9r7+/JM/p6emJr69vkZtwUn51odnf9O1V76rNIpTo3jSAmpW8SU7PYd5f8sOMS7mcBJsLPpTS+XW9mFe4LKu++35+fphMpmJHVhITE4sdgbkiICDghvu7ublRqVKlW+5zs+cULqbjaDC6wfHVcHKj6jTCxkxGQ2EJ53frT5AjJZyuY8PHkHMZAptDw56q0wjFrDrgeHh4EBYWxooVK4rcv2LFCiIiIm74mPDw8GL7//HHH7Rq1Qp3d/db7nOz5xQupkItaPmEvr36PZB+IpczKKw6fmU9OJOcweI98arjCFtIideLdwG6vKlXuQiXZvXjdyNHjuS7775j5syZxMTE8MorrxAbG8vw4Xqj69ixY3niiScK9x8+fDinTp1i5MiRxMTEMHPmTGbMmMGoUaMK93n55Zf5448/mDJlCgcPHmTKlCmsXLmSESNGWPvlCEfRfhSYPOHURji+RnUaYWNe7iaeKijh/GqtlHC6hA0fQW4mBLWBug+oTiPsgNUHnMGDBzNt2jQmTpxI8+bNWbduHUuWLKFmzZoAxMfHF1kTJzg4mCVLlrBmzRqaN2/OO++8w6effsrAgQML94mIiGDu3Ll8//33hIaGMmvWLObNm0ebNm2s/XKEozBXg1ZP69ur3pWjOC5oSNtaeHuYOJiQytrD51XHEdaUHAvbvte3u7whR28EYIN1cOxRST5HLxxY6jn4tDnkpMPf5kGD7qoTCRt75/cDzNhwgra1KzL32XDVcYS1LHoRdvwIwR3gyd9UpxFWZDfr4AihVDl/aD1M3179HuTLxaauZuj9wbgZDWw+fpFdccmq4whruHAMds7Wt68U7wqBDDjC2bUbAR7lIGEPHJSf7FxN1fJl6NO8KgBfr5Wmeae0dqpetFu3G9SQyxTEVTLgCOfmXRHaPq9vr54M+Xlq8wibu7Lw37L9CZxIuqw4jShV5w/B3v/p21cKd4UoIAOOcH7h/wIvM5yPgX0LVKcRNtYgoBxdGlZB0+CbdVLC6VTWTAYtHxr2gmotVacRdkYGHOH8ypSHiJf07TWTIU+W73c1wzvqR3GidpwmMTVTcRpRKhL2wv6FgEGO3ogbkgFHuIY2w8G7Elw8Bnvmqk4jbOy+WhVoUaM82bn5zNp4UnUcURpWT9b/26Q/+DdRm0XYJRlwhGvwLAv3v6Jvr5kCudlq8wibMhgMhUdxftp8irQsOYrn0M5sh0OL9WLdTmNVpxF2SgYc4TpaDYWy/mCJhZ0/qU4jbKxbI39qV/YhNTOXOVtib/8AYb9WT9L/GzoYKtdXm0XYLRlwhOvw8NYrHADWfQA5GWrzCJsyGg0810Ev4Zyx4QTZubIukkM6FQ1HV+qFuh1Hq04j7JgMOMK1hD0JvtUhNf7q0u7CZfRrUY0q5TxJSMnk111nVMcRd2P1e/p/W/wdKtZWm0XYNRlwhGtx84SOr+rbGz6CbFkXxZV4upl4+v5gQP/IeH6+yzXVOLbja+HkejB5QIdXVacRdk4GHOF6mj8OFWrB5fOw9RvVaYSNPdamBuU83TiSmMaqg4mq44g7pWl6cS5A2D/AXF1tHmH3ZMARrsfkDh3H6NsbP4HMFLV5hE35ernzWNsaAHy9TuobHMbRlXB6K7h5QfuRqtMIByADjnBNoY+AX33IuASbv1SdRtjY0+2CcTcZ+OvkJbafuqg6jrgdTYNV7+jbrYdBuQC1eYRDkAFHuCaj6er6GdGfQ7p8k3Ml/r5e9G9RDYCv1kp9g907+DvE7waPsnqBrhB3QAYc4boa9wP/ppCVog85wqU8W1DCueLAOY4mpilOI24qP//qujdthoOPn9o8wmHIgCNcl9F4tcNm81eQdl5tHmFTdauUpVtjfwC+kWtx7Nf+BZB4ADzNEPGC6jTCgciAI1xbg4egagvIuQwbp6lOI2zsSn3Dwp1nOJciJZx2Jy9XL8gFfbgpU0FtHuFQZMARrs1ggM5v6Nt/fQcp8WrzCJsKq1mB+2pVICdPY+aGE6rjiOvt/R9cOAplKuqnp4QoARlwhKj7AAS1hdxMWP+h6jTCxq4cxZm9JZaUzBzFaUShvBxY876+ff8I8PJVGkc4HhlwhDAYoMs4fXv7LEiWIkZX0rlBFepVKUtaVi6zN8t7bzd2/geST4FPFbhvmOo0wgHJgCMEQHAH/ZafoxdxCpdhNBp4ruAozsyNJ8jKzVOcSJCTefX/w/b/1otyhSghGXCEuOLKtTg7Z8MF+VSNK+nTrCqBZi/Op2axcIeUcCq3fRaknAHfahD2lOo0wkHJgCPEFTXaQN1uoOXB2qmq0wgb8nAzMlRKOO1DdvrVa+E6jAJ3L7V5hMOSAUeIa125FmfPPDh/SG0WYVOPtq6Br5cbx5Mu88eBc6rjuK6/voXLiVC+JjT/u+o0woHJgCPEtaq2gIa9AO3q+hvCJZT1dGNIeE0Avlp7DE2Tozg2l5UKG6bp2x1fAzcPpXGEY5MBR4jrdX4dMMD+hZCwV3UaYUNPRQTj4WZkV1wyW09IP5nNbf4KMi5CpboQOlh1GuHgZMAR4nr+TaBJf337SgeOcAmVy3kyKKw6AF+vkxJOm8q4BJs+07c7jQWTm9o8wuHJgCPEjXQaCwYjHFoCZ7arTiNsaFj72hgMsOpgIocSUlXHcR3RX0CWBao0hiYDVKcRTkAGHCFupHJ9CH1U3171ntoswqaC/Xzo0TQAgK+lhNM2Ll+AzV/q251f14twhbhH8rdIiJvpOBqMbnDsTzgVrTqNsKHnOugL/y3adZYzyRmK07iAjdMgOw0CmxVc5C/EvZMBR4ibqRgMLQo+prrqXZBP1biMZkHlCa9didx8KeG0utQE2Pqtvt35Db06RYhSIAOOELfS4VUwecCpDXBireo0woae61gbgDlbY0lOz1acxomt/whyM6D6fVCvm+o0wonIgCPErZirQ9g/9O1V78lRHBfSsX5lGgaUIz07j/9sPqU6jnOynIbt3+vbXeTojShdVh1wLl26xJAhQzCbzZjNZoYMGUJycvJN98/JyeG1114jJCQEHx8fqlatyhNPPMHZs2eL7NepUycMBkOR26OPPmrNlyJcWfuR4FYGTm+FIytUpxE2YjAYGF5Qwvn9xpNk5kgJZ6lb9wHkZUOt9hDcUXUa4WSsOuA89thj7Nq1i2XLlrFs2TJ27drFkCFDbrp/eno6O3bs4M0332THjh0sWLCAw4cP06dPn2L7Dhs2jPj4+MLb119/bc2XIlxZuQBo/Yy+vVquxXElPUMDqVa+DBcuZzN/+2nVcZzLxROw8z/6dudxcvRGlDqrraQUExPDsmXL2Lx5M23atAHg22+/JTw8nEOHDtGgQYNijzGbzaxYUfQn5M8++4zWrVsTGxtLjRo1Cu/39vYmICDAWvGFKKrdCNj2PcTvhoO/Q6PeqhMJG3A3GXmmfTATfjvAt+uP87fWNTAZ5RtxqVg7FfJzoc4DUDNcdRrhhKx2BCc6Ohqz2Vw43AC0bdsWs9nMpk2b7vh5LBYLBoOB8uXLF7l/9uzZ+Pn50aRJE0aNGkVq6s0X5MrKyiIlJaXITYgS8fGDts/r26snQX6+2jzCZgbfF0R5b3dOXUhn2b4E1XGcw/nDsGeuvn2l4FaIUma1ASchIYEqVaoUu79KlSokJNzZPxKZmZmMGTOGxx57DF9f38L7H3/8cebMmcOaNWt48803iYqKYsCAm698OXny5MLrgMxmM0FBQSV/QUKE/ws8zZB4APYvUJ1G2Ii3hxtPhNcCpISz1Kx9H7R8aPAQVAtTnUY4qRIPOOPHjy92ge/1t23btgH6RXrX0zTthvdfLycnh0cffZT8/HymT59e5PeGDRtG165dadq0KY8++ijz589n5cqV7Nix44bPNXbsWCwWS+EtLi6upC9bCChTASJe1LfXTIa8XLV5hM08GV4TL3cje89YiD52QXUcx3ZuP+yL0rc7v642i3BqJb4G54UXXrjtJ5Zq1arFnj17OHfuXLHfO3/+PP7+/rd8fE5ODo888ggnTpxg1apVRY7e3EjLli1xd3fnyJEjtGzZstjve3p64unpecvnEOKOtB0Om6fDhaOw93/Q/DHViYQNVCrrySOtgvgx+hRfrj1GRF0/1ZEc15UC28b9ICBEaRTh3Eo84Pj5+eHnd/v/ucPDw7FYLGzdupXWrVsDsGXLFiwWCxERETd93JXh5siRI6xevZpKlSrd9mvt37+fnJwcAgMD7/yFCHE3PMvB/SNgxVuw5n1oOgjcPFSnEjYwrH1t/rP5FOuPJLH/rIUmVc2qIzmeszv1i/QNRjl6I6zOatfgNGrUiO7duzNs2DA2b97M5s2bGTZsGL169SryCaqGDRuycOFCAHJzcxk0aBDbtm1j9uzZ5OXlkZCQQEJCAtnZ+kqix44dY+LEiWzbto2TJ0+yZMkSHn74YVq0aEG7du2s9XKEuOq+YeBTBZJPwa7/qE4jbCSoojc9Q6sC8PXa44rTOKgrxbUhj0Dl4p+kFaI0WXUdnNmzZxMSEkJkZCSRkZGEhoby008/Fdnn0KFDWCwWAE6fPs2iRYs4ffo0zZs3JzAwsPB25ZNXHh4e/Pnnnzz44IM0aNCAl156icjISFauXInJZLLmyxFC5+EN7f+tb6/7P8jJVJtH2MxzHfT6hsV744m7mK44jYOJ3QJHV4DBpBfZCmFlBs0FPxKQkpKC2WzGYrHc9voeIW4oJxM+awkpZ6D7FP3aHOEShszYwvojSTwZXpMJfZuqjuM4fugNJ9ZByyegz2eq0wgHVZLv39JFJcTdcPeCDqP07fUfQrb8NO8qrtQ3zNsWx8XLUsJ5R06s029Gd73AVggbkAFHiLvV/O9QviZcToS/vlWdRthIRJ1KNK3mS2ZOPj9sOqk6jv3TtKvX3oQ9BeVr3HJ3IUqLDDhC3C03D+g0Rt/eMA2ybr6atnAe15Zw/hh9kvRsWQ/plo7+CXGbwc3r6rVrQtiADDhC3IuQR6BSXci4CJu/Up1G2Ej3JgHUqOjNpfQcft4mJZw3pWl6QS3Afc+AryzlIWxHBhwh7oXJDTqN1bc3fQYZl9TmETbhZjIyrH0wAN+uP05unnST3dChJfraN+4+emGtEDYkA44Q96rJAKjSGLIsEP2F6jTCRh5uFUQlHw9OX8pg8d541XHsT37+1VWL2zwHZSurzSNcjgw4Qtwr4zWrsm7+Ei4nqc0jbMLL3cSTEbUA+GrtcSnhvN6BX+DcPvD0vdrhJoQNyYAjRGlo2AsCm0F2GmycpjqNsJEnwmtSxt1ETHwK64/IYFsoP08vpAUIfwG8K6rNI1ySDDhClAaDATq/oW9v/Q5SE9TmETZR3tuDR1sHAfDV2mOK09iRvT9D0mEoUwHaPq86jXBRMuAIUVrqdYPqrSE3A9ZOUZ1G2Mgz7WtjMhrYdOwCe04nq46jXk7G1aM37V4GL1ktXqghA44QpcVggAfe0re3zYQjK9TmETZRrXwZ+jSTEs5Cy8fBpZNQNgBaP6s6jXBhMuAIUZqC20Pr5/TthcPlVJWLeK6jXsK5dF88J5MuK06j0IFFsG2Gvt3/S/DwUZtHuDQZcIQobd0mgn8IpCfBwuf0j8sKp9YwwJdODSqTr+nr4rik5DhY9IK+3W4E1OmiNI4QMuAIUdrcvWDQTHD3huNrYNOnqhMJG7hS3/Dz9tOcT81SnMbG8nJhwbOQaYFqYdDlDdWJhJABRwirqFwfehRcaLzqHTi9XW0eYXVtgivSLKg82bkuWMK5/v8gdhN4lIOB34HJXXUiIWTAEcJqWgyBJv0hPxeinobMFNWJhBUZDAaeL7gW58fok1zOcpESzlObrn5qsNfHULG22jzCLhw7n6Z88UsZcISwFoMBek0Dcw39UyWLR+rlg8JpdWscQLCfDymZufy0+ZTqONaXfhGihoGWD80eg9CHVScSdiDuYjq9P9vAk9//RUpmjrIcMuAIYU1lyuuH7A0mffGz3XNVJxJWZDIaeL6Tfi3OJyuPEHcxXXEiK9I0+O0lSDkNFevAQ1NVJxJ2QNM0xi7YS3p2Hpk5eZT1cFOWRQYcIaytRhvoXNA4vvjfkHRUbR5hVYNaVqdNcEUycvJ4LWqP8sP0VrP9e4j5DYzuMGgGeJZTnUjYgZ+3nWbD0SQ83YxMGRiK0WhQlkUGHCFs4f6RUKs95FzWr8fJdbFP2bgQo9HAlIGheLkb2XTsAnP/ilMdqfSdOwDLCob2ruOhagulcYR9SLBk8s7iAwD8O7I+wX5q10GSAUcIWzCaYMA3UKYixO+GPyeqTiSsqJafD6MiGwAwaXEM8ZYMxYlKUU4GzH8acjOhbjdo+0/ViYQd0DSNN37ZS2pmLs2CyjP0fvUXm8uAI4St+FaFftP17ejPpcrByf2jXTDNg8qTmpXLuIX7nOdU1fJxcD4GfKpAvy/BKN9GBCzafZaVMYm4mwx8MCgUk8JTU1fI30whbKlBj+uqHM6pzSOsxmTU/6H3MBlZdTCRX3edVR3p3sX8drWKYcDXULay2jzCLiSlZTF+0X4AXuxSj/r+9nE9lgw4Qthat4ng31SqHFxAPf9yvPRAXQDG/7bfsVc4tpyGX69UMbwsVQyi0PhF+7mUnkOjQN/CTxHaAxlwhLC1K1UObmXg+GqI/kx1ImFFz3WsQ+NAX5LTcwp/ynU4+Xn6ejeZyVC1JXSWKgahW74/gd/3xBcesXQ32c9YYT9JhHAllRtcrXL4cyKckSoHZ+VuMjK14JqExXvjWbo3XnWkklt3TRXDoBng5qE6kbADlvQc3vhlHwDPdahN02pmxYmKkgFHCFVaPgGN++lVDvOHSpWDE2tazczwghqHN3/dT3J6tuJEJXAqGta+r2/3+kiqGEShdxYf4HxqFnUq+/DSA/VUxylGBhwhVDEYoPcnBVUOJ/RFAJ3lkzaimBe71KNulbIkpWUx8fcDquPcmfSLEPVMQRXD3yD0EdWJhJ1Ye/g887efxmCAqYNC8XI3qY5UjAw4QqhUpMrhf1Ll4MS83E1MHRSKwQALdpxh9cFE1ZFurUgVQ2146APViYSdSM3MYWzUHgD+ERFMWM2KihPdmAw4QqgmVQ4uo2WNCjzdLhiA1xfuJVVhEeFtFalimClVDKLQlGUHOWvJpEZFb0Y9WF91nJuSAUcIe1CsysGBrtEQJTIqsgE1K3kTb8lk8tKDquPcWGKMVDGIG9p8/AL/2RwLwPsDQvBWWKZ5OzLgCGEPCqscKhRUOUxQnUhYSRkPE+8PCAXgv1ti2XQsSXGi6xSpYugqVQyiUEa2XiAL8LfWNYio66c40a3JgCOEvfCtCn2vrXJYqTaPsJrwOpV4vE0NAMZE7SU9O1dxomv88QYkHpAqBlHMRysOcepCOoFmL8Y+1FB1nNuSv7lC2JOGD0HrZ/XtX6TKwZmN6dGQqmYvYi+m8+Efh1XH0cX8Dn99p2/3/wrKVlGbR9iNnbGXmLHhBACT+ofg6+WuONHtyYAjhL3p9g5UaQKXz0uVgxMr5+XOewNCAJi58QTbT11SG8hyGn79l74d8RLUfUBtHmE3snLzGD1/D/kaDGhRjc4NHWPwteqAc+nSJYYMGYLZbMZsNjNkyBCSk5Nv+ZinnnoKg8FQ5Na2bdsi+2RlZfHiiy/i5+eHj48Pffr04fTp01Z8JULYkFQ5uIzODaowsGV1NA1Gz99NZk6emiBFqhhaQJc31eQQdunzVUc5kpiGX1kP3uzVWHWcO2bVAeexxx5j165dLFu2jGXLlrFr1y6GDBly28d1796d+Pj4wtuSJUuK/P6IESNYuHAhc+fOZcOGDaSlpdGrVy/y8hT94yBEaavSEHoUrB4rVQ5O7c1ejfAr68mx85f5bNURNSEKqxjKwkCpYhBX7T9rYfqaYwC807cpFXwc5++G1QacmJgYli1bxnfffUd4eDjh4eF8++23/P777xw6dOiWj/X09CQgIKDwVrHi1UWELBYLM2bM4MMPP6Rr1660aNGC//znP+zdu5eVK298UWZWVhYpKSlFbkLYvZZPSpWDCyjv7cG7/ZoC8NXa4+w7Y7FtgGurGHp+BJXspw1aqJWTl8/o+XvIy9d4KCSAHiGBqiOViNUGnOjoaMxmM23atCm8r23btpjNZjZt2nTLx65Zs4YqVapQv359hg0bRmLi1RU/t2/fTk5ODpGRkYX3Va1alaZNm970eSdPnlx4msxsNhMUFHSPr04IG7i+ymHJKNWJhJV0bxpAz5BA8vI1Xp2/h5w8G113lXGpaBVDs8G2+brCIXyz7jj7z6ZQ3tudCX2aqo5TYlYbcBISEqhSpfiFSFWqVCEhIeGmj+vRowezZ89m1apVfPjhh/z111906dKFrKyswuf18PCgQoUKRR7n7+9/0+cdO3YsFoul8BYXF3cPr0wIG7q2ymHPPKlycGLj+zShgrc7MfEpfFVwSsCqNA0WSRWDuLGjial8slI/Zfp278ZULuepOFHJlXjAGT9+fLGLgK+/bdu2DQCDwVDs8Zqm3fD+KwYPHkzPnj1p2rQpvXv3ZunSpRw+fJjFixffMtetntfT0xNfX98iNyEcRo020OmaKocLNvjmJ2yucjlP3u7dBIDPVh3l8LlU637B7bMgZpFUMYhirhxJzM7Lp3ODyvRrXk11pLtS4gHnhRdeICYm5pa3pk2bEhAQwLlzxdfwOH/+PP7+/nf89QIDA6lZsyZHjuiTZEBAANnZ2Vy6VPQjlYmJiSV6XiEcSvuRUPN+yE4rWGVWqhycUd/mVXmgYRWy8/J5teDaB6tIjIFlY/Ttrm9LFYMoYtamk+yMTaacpxuTBoTc8qCEPSvxgOPn50fDhg1vefPy8iI8PByLxcLWrVsLH7tlyxYsFgsRERF3/PUuXLhAXFwcgYH6xU1hYWG4u7uzYsWKwn3i4+PZt29fiZ5XCIdSpMphl1Q5OCmDwcB7/UMo5+nG7rhkZhYsrFaqrq1iqPMAtP1X6X8N4bBOXbjMB8v1jrSxDzUi0FxGcaK7Z7VrcBo1akT37t0ZNmwYmzdvZvPmzQwbNoxevXrRoEGDwv0aNmzIwoULAUhLS2PUqFFER0dz8uRJ1qxZQ+/evfHz86N///4AmM1mhg4dyr///W/+/PNPdu7cyd///ndCQkLo2rWrtV6OEOqZq0HfL/RtqXJwWgFmL8b1bATA//1xiBNJl0v3CxRWMVTWVyuWKgZRID9fY0zUXjJz8omoU4m/tXbsD+RY9W/27NmzCQkJITIyksjISEJDQ/npp5+K7HPo0CEsFv1jkSaTib1799K3b1/q16/Pk08+Sf369YmOjqZcuavnhz/++GP69evHI488Qrt27fD29ua3337DZDJZ8+UIoV7DnnDfMH1bqhyc1uD7gmhXtxJZufm8FrWH/NI6VSVVDOIW5vwVS/TxC5Rx1wthHfXU1BUGTdOsdJLXfqWkpGA2m7FYLHLBsXA8ORnw7QOQuB/qdIHHo+SncCcUdzGdyI/XkZGTxzv9mjKkbc17e0LLafiynb5accSLEPluqeQUzuFscgaRH68jLSuXt3o15un7g1VHuqGSfP+WfxWFcDTuZa5WORxbpZ+uEk4nqKI3r3XXT+e/vySG05fS7/7J8vNgwbPXVDG8VTohhVPQNI3XF+4lLSuXsJoVeDKilupIpUIGHCEcUZEqhwlS5eCkngivRauaFbicncfrC/dx1wfc138IpzZKFYO4oQU7zrDm0Hk83IxMGRiKyejYp6aukAFHCEfV8klo3PdqlUOWlddNETZnNBqYMigUDzcj6w6fZ/72uygVjt0Maybr21LFIK6TmJrJxN8PADCiaz3qVimrOFHpkQFHCEdVWOUQpFc5LJYqB2dUp3JZXulaH4B3fj9AYkrmnT/42iqG0EelikEU89Yv+7Fk5NC0mi/Ptq+tOk6pkgFHCEdWpkJBlYMR9syVKgcnNax9MCHVzKRk5jLulzs8VXWlisESp1cx9Pw/6wcVDmXJ3niW7U/AzWhg6sBmuJmcayRwrlcjhCuq0VaqHJycm8nIBw+H4m4ysOLAOX7fE3/7BxVWMbjp191IFYO4xqXL2bz16z4A/tm5Lo2rOt8nimXAEcIZtP831GwnVQ5OrGGAL//sVBeA8Yv2cyEt6+Y7X1vF8MDbUK2lDRIKRzLx9wMkpWVT378sL3SuqzqOVciAI4QzMJpgwLdXqxxWTVSdSFjBvzrXpYF/OS5czmbCbwduvFNOhn7ReW6mvk5S+Au2DSns3p8x51i48wxGA0wd1AwPN+ccBZzzVQnhiq6tctj0GRyVKgdn4+FmZOqgUIwGWLT7LCsO3GAl6z/e1BeB9KkM/aSKQRSVkpnDuIX6qaln2temeVB5tYGsSP7mC+FMrq1yWDgc0hLV5hGlrllQeYZ10D/tMm7hXiwZOVd/8+Bi+Otbfbv/V1DOX0FCYc8mL4khISWTYD8fRnarrzqOVcmAI4SziXwHqjSBy+f1ISc/X3UiUcpe6VqfYD8fElOzmLQ4Rr/TcgZ+LWgGj3gR6kr5sChq49Ek5myNA+D9ASF4uTt3f6MMOEI4myJVDn/C5i9UJxKlzMvdxNRBoRgMMG9bHOsPJehVDBmXILC5VDGIYi5n5TJmwR4AngivSZvalRQnsj4ZcIRwRlUaQveC1WtXToAzO9TmEaXuvloVeaKggPPgz+Ph1Aa9imHQTKliEMV8sPwQcRczqFa+DKO7N1QdxyZkwBHCWYU9BY36QH6O/tFxqXJwOqO7N6S770mezilY4LHnh1LFIIrZdvIiP0SfBGDygBDKerqpDWQjMuAI4awMBujzKfhWlyoHJ+WTn8o0t88xGTQW5N3PVt9I1ZGEncnMyWP0/D1oGjzSqjod6ldWHclmZMARwplJlYPzKqhi8Eo/S5JHNd7KeYrXovaQmZOnOpmwI9NWHuF40mWqlPNkXM/GquPYlAw4Qji7muHQsWBVW6lycB47fiisYvB69Ht8fCtwIukyH684rDqZsBN7Tifz7frjALzXPwRzGXfFiWxLBhwhXEGHUVerHKKGSpWDo0s8CEuvVDG8RdnabXivXwgA364/zu64ZHXZhF3Izs1n9Pw95OVr9G5WlW6NXW9NJBlwhHAF11Y5nN0Jq95RnUjcrZzMgr6xjIIqhhcB6NrYn77Nq5Kvwavzd5OVK6eqXNmXa45xMCGVij4ejO/tWqemrpABRwhXYa4GfT7Xtzd9KlUOjuqPN25axfB27yZU8vHg8Lk0vlgtpyJd1aGEVD5ffQSA8X2aUKmsp+JEasiAI4QradQL7ntG35YqB8dzbRVDv+JVDBV9PJjQtwkA01cfJSY+xdYJhWK5efmMnr+bnDyNbo396R0aqDqSMjLgCOFqIt+FKo2lysHRXFvFEP4C1LtxFUPPkEAebOJPbr7G6Pl7yM2T99eVzNhwgt2nLfh6ufFuv6YYDAbVkZSRAUcIV1NY5eAlVQ6OIj+vaBXDA2/fdFeDwcA7fZvi6+XG3jMWvl1/wnY5hVLHz6fxUcGn6N7o1Rh/Xy/FidSSAUcIV1SlkVQ5OJL1H+lVDO4+d1TFUMXXi7d666eqPl55mGPn02yRUiiUn6/xWtQesnLzaV/Pj4fDqquOpJwMOEK4qrB/QKPeepVD1FCpcrBXsZthTcEwWoIqhoEtq9GxfuUiHxcWzuunzaf46+QlfDxMTB4Q4tKnpq6QAUcIV2UwQO+CKoeLx2HJq6oTietlJEPUM6DlQcgj0OzRO36owWBg0oAQfDxMbD91iR8LuoiE84m7mM6UZQcBGNOjIdUreCtOZB9kwBHClXlXhIHf6lUOu+fA7nmqE4krNA1+exkscVAhWD96U8KfyquVL8OYhxoBMHXZIeIuplsjqVBI0zTGLthLenYerYMr8nibmqoj2Q0ZcIRwdTUjrqlyGClVDvZixw9w4BcwusGgGeDle1dP83jrGrQJrkhGTh6vRe1B0+RUlTP5edtpNhxNwtPNyJSBoRiNcmrqChlwhBBS5WBvrqtioFrYXT+V0WhgysBQvNyNbDp2gbl/xZVSSKHauZRM3ll8AIBRkQ0I9vNRnMi+yIAjhCiocvgGvMpLlYNq11Yx1O5cWMVwL2r5+TAqsgEAkxbHEG/JuOfnFGppmsa4hXtJzcylWVB5nr4/WHUkuyMDjhBCZ64Ofa+tcvhTbR5XteJNvYrB2w/6f12kiuFe/KNdMM2DypOalcu4hfvkVJWDW7T7LCtjEnE3GfhgUCgmOTVVjAw4QoirGvWGVkP1balysL2DS2DrN/p2/+JVDPfCZNS/EXqYjKw6mMivu86W2nML27qQlsWE3/RTUy92qUd9/3KKE9knGXCEEEU9+B5UbgSXE+GX56XKwVZSzsKv/9S3w1+Aet1K/UvU8y/HSw/UBWD8b/s5n5pV6l9DWN/bi/Zz8XI2jQJ9eb7Tna2L5IpkwBFCFHVtlcPRlbB5uupEzq9IFUMz/cJiK3muYx0aB/qSnJ7D+EX7rfZ1hHUs35/A73viC4/IuZvk2/jNWPVP5tKlSwwZMgSz2YzZbGbIkCEkJyff8jEGg+GGtw8++KBwn06dOhX7/UcfvfMFsIQQt+HfGB6cpG+vHK9feCysZ8NHcHK9XsUwcCa4eVrtS7mbjEwtuGZj8d54lu6Nt9rXEqXLkp7DG7/sA+C5DrVpWs2sOJF9s+qA89hjj7Fr1y6WLVvGsmXL2LVrF0OGDLnlY+Lj44vcZs6cicFgYODAgUX2GzZsWJH9vv76a2u+FCFcT6unoWEvvcph/tNS5WAtsVtg9ZUqhv8Dv7pW/5JNq5kZ3rE2AG/+up/kdFkWwBG8s/gA51OzqFPZh5ceqKc6jt1zs9YTx8TEsGzZMjZv3kybNm0A+PbbbwkPD+fQoUM0aNDgho8LCAgo8utff/2Vzp07U7t27SL3e3t7F9tXCFGKDAbo8xmc3XW1yqH/V6pTOZeMZH3docIqhr/Z7Eu/2KUey/ef42hiGhN/P8BHjzS32dcWJbf28Hnmbz+NwQBTB4Xi5W5SHcnuWe0ITnR0NGazuXC4AWjbti1ms5lNmzbd0XOcO3eOxYsXM3To0GK/N3v2bPz8/GjSpAmjRo0iNfXmP11mZWWRkpJS5CaEuANS5WA9RaoYat1VFcO98HI3MXVQKAYDLNhxhtUH5RNz9iotK5fXF+wF4B8RwYTVrKg4kWOw2oCTkJBAlSpVit1fpUoVEhIS7ug5fvjhB8qVK8eAAQOK3P/4448zZ84c1qxZw5tvvklUVFSxfa41efLkwuuAzGYzQUFBJXsxQriymhHQ8TV9e/FI/WiOuHc7frxaxTBw5l1XMdyLljUq8HQ7fYG41xfuJTUzx+YZxO1NWXqQM8kZ1KjozagH66uO4zBKPOCMHz/+phcCX7lt27YN4IZ17Zqm3XGN+8yZM3n88cfx8vIqcv+wYcPo2rUrTZs25dFHH2X+/PmsXLmSHTt23PB5xo4di8ViKbzFxclS5UKUSPtRUCNCr3KYL1UO9+z8IVhaMDR2eROq330Vw70aFdmAmpW8ibdkMnnpQWU5xI1tPn6BnzafAuD9ASF4e1jtyhKnU+IB54UXXiAmJuaWt6ZNmxIQEMC5c+eKPf78+fP4+99+8ar169dz6NAhnnnmmdvu27JlS9zd3Tly5MgNf9/T0xNfX98iNyFECZjcrqly2AGr31WdyHEVqWLoBBEvKY1TxsPE+wNCAfjvllg2HUtSmkdclZGtF6QC/K11DSLq+ilO5FhKPAr6+fnh53f7P+Tw8HAsFgtbt26ldevWAGzZsgWLxUJERMRtHz9jxgzCwsJo1qzZbffdv38/OTk5BAYG3v4FCCHuTvkg/aLj/w2BjZ9AcEeo+4DqVI5nxVtwbl+pVzHci/A6lXi8TQ1mb4llTNRelo1oL0cK7MBHKw5x6kI6gWYvxj7UUHUch2O1/7MaNWpE9+7dGTZsGJs3b2bz5s0MGzaMXr16FfkEVcOGDVm4cGGRx6akpPDzzz/f8OjNsWPHmDhxItu2bePkyZMsWbKEhx9+mBYtWtCuXTtrvRwhBEDjPvrHx0GqHO7GwSWwtWBJi35fQjn7+STomB4NqWr2IvZiOv+3/LDqOC5vZ+wlZmw4AcCk/iH4erkrTuR4rPqjw+zZswkJCSEyMpLIyEhCQ0P56aefiuxz6NAhLBZLkfvmzp2Lpmn87W/FPzLp4eHBn3/+yYMPPkiDBg146aWXiIyMZOXKlZhM8rE5IazuwUlS5XA3rq1iaPsvqB+pNs91ynm5896AEAC+33SC7acuKU7kurJy8xg9fw/5GgxoUY3ODYt/YEfcnkFzwUrZlJQUzGYzFotFrscR4m6cOwDfdobcTIh8DyJeUJ3IvuXnwY999dWKA0LhmZVWXa34Xvz7f7uJ2nGaOpV9WPxSe1lvRYEP/zjEZ6uO4lfWk5UjO1De20N1JLtRku/f6k/+CiEcj1Q5lMy1VQyDvrfb4QbgzV6N8CvrybHzl/ls1Y0/uCGsZ/9ZC1+uOQbAO32byHBzD2TAEULcHalyuDMKqhjuRXlvD97t1xSAr9YeZ98Zy20eIUpLTl4+o+fvITdf46GQAHqEyAdn7oUMOEKIu3OlysG3WkGVw2jViexPRjJEPVNQxfCwTasY7kX3pgH0DAkkL1/j1fl7yMmT66xs4Zt1x9l/NoXy3u5M6NNUdRyHJwOOEOLueVeEAVeqHP4Le/6nOpH9KKxiiIXyNaHnRzatYrhX4/s0oYK3OzHxKXxVcMpEWM/RxFQ+WamfEny7d2Mql7Pf05iOQgYcIcS9qdUOOhQcvfldqhwK7fzpahXDIDVVDPeicjlP3u7dBIDPVh3l8Dk5BWktV46UZefl07lBZfo1r6Y6klOQAUcIce86vAo1wiE7VaocQK9iuHLKrssbUL2V2jx3qW/zqjzQsArZefm8On8Pefku96Fbm5i16SQ7Y5Mp5+nGpAEhd1xnJG5NlqoUQtw7k5t+quqrdnqVw5J/Q70HVadSZ83ka6oYXlad5q4ZDAbe6x/C1o/WsjsumZkbTjCsQ23VsZzKqQuX+WC53gE29qFGBJrLKE7kPGQdHFkHR4jSc2CRXuUg9CqG5zfa1WrFd2vu1ljGLNiLp5uRZSM6EOznozqSU8jP13j8uy1EH79ARJ1KzH6mjRy9uY2SfP+WIzhCiNLTuA90Ha9XErgykwd0GuMUww3A4PuC+G3PWTYevcBrUXuYO6wtRqN8I75Xc/6KJfr4Bcq464WnMtyULjmCI0dwhBDituIuphP58ToycvJ4p19ThrStqTqSQzubnEHkx+tIy8rlrV6Nefr+YNWRHIKsZCyEEKJUBVX05rXuelHy+0tiOH0pXXEix6VpGq8v3EtaVi5hNSvwZEQt1ZGckgw4Qggh7sgT4bVoVbMCl7PzeH3hPlzwBECpWLjzDGsOncfDzciUgaGY5HSfVciAI4QQ4o4YjQamDArFw83IusPnmb/9tOpIDicxNZMJvx0AYETXetStUlZxIuclA44QQog7VqdyWV7pWh+Ad34/QGJKpuJEjuWtX/ZjycihaTVfnm0vH7m3JhlwhBBClMiw9sGEVDOTkpnLuF/kVNWdWrI3nmX7E3AzGpg6sBluJvkWbE3ypyuEEKJE3ExGPng4FHeTgRUHzvH7nnjVkezepcvZvPXrPgD+2bkujavKJ3itTQYcIYQQJdYwwJd/dqoLwPhF+7mQlqU4kX2b+PsBktKyqe9flhc611UdxyXIgCOEEOKu/KtzXRr4l+PC5ezCC2dFcasOnmPhzjMYDTB1UDM83ORbry3In7IQQoi74uFmZOqgUIwGWLT7LCsOnFMdye6kZObw+gL91NQz7WvTPKi82kAuRAYcIYQQd61ZUPnCAs5xC/diychRnMi+TF4SQ0JKJsF+PozsVl91HJciA44QQoh78krX+gT7+ZCYmsWkxTGq49iNjUeTmLM1DoD3B4Tg5W5SnMi1yIAjhBDinni5m5g6KBSDAeZti2P9kfOqIyl3OSuXMQv2APBEeE3a1K6kOJHrkQFHCCHEPbuvVkWeDK8FwJiovVzOylUbSLEPlh8i7mIG1cqXYXT3hqrjuCQ31QGEEEI4h1cfbMDKmHOcvpTBsz9to75/OdWRlMjOzee/W2MBmDwghLKe8q1WBflTF0IIUSp8PN14f0Aof5+xhY1HL7Dx6AXVkZR6pFV1OtSvrDqGy5IBRwghRKm5v54f0x9vyf6zFtVRlPLxdCs8ZSfUkAFHCCFEqXooJJCHQgJVxxAuTi4yFkIIIYTTkQFHCCGEEE5HBhwhhBBCOB0ZcIQQQgjhdGTAEUIIIYTTkQFHCCGEEE5HBhwhhBBCOB0ZcIQQQgjhdKw64Lz33ntERETg7e1N+fLl7+gxmqYxfvx4qlatSpkyZejUqRP79+8vsk9WVhYvvvgifn5++Pj40KdPH06fPm2FVyCEEEIIR2TVASc7O5uHH36Y559//o4fM3XqVD766CM+//xz/vrrLwICAujWrRupqamF+4wYMYKFCxcyd+5cNmzYQFpaGr169SIvL88aL0MIIYQQDsagaZpm7S8ya9YsRowYQXJy8i330zSNqlWrMmLECF577TVAP1rj7+/PlClTeO6557BYLFSuXJmffvqJwYMHA3D27FmCgoJYsmQJDz74YLHnzcrKIisrq/DXKSkpBAUFYbFY8PX1Lb0XKoQQQgirSUlJwWw239H3b7u6BufEiRMkJCQQGRlZeJ+npycdO3Zk06ZNAGzfvp2cnJwi+1StWpWmTZsW7nO9yZMnYzabC29BQUHWfSFCCCGEUMquBpyEhAQA/P39i9zv7+9f+HsJCQl4eHhQoUKFm+5zvbFjx2KxWApvcXFxVkgvhBBCCHtR4jbx8ePHM2HChFvu89dff9GqVau7DmUwGIr8WtO0Yvdd71b7eHp64unpWWRf0A91CSGEEMIxXPm+fSdX15R4wHnhhRd49NFHb7lPrVq1Svq0AAQEBAD6UZrAwMDC+xMTEwuP6gQEBJCdnc2lS5eKHMVJTEwkIiLijr7OlQuW5VSVEEII4XhSU1Mxm8233KfEA46fnx9+fn53HepWgoODCQgIYMWKFbRo0QLQP4m1du1apkyZAkBYWBju7u6sWLGCRx55BID4+Hj27dvH1KlT7+jrVK1albi4OMqVK3fbI0Oi9Fy5uDsuLk4u7rYz8t7YL3lv7Je8N7anaRqpqalUrVr1tvuWeMApidjYWC5evEhsbCx5eXns2rULgLp161K2bFkAGjZsyOTJk+nfvz8Gg4ERI0YwadIk6tWrR7169Zg0aRLe3t489thjAJjNZoYOHcq///1vKlWqRMWKFRk1ahQhISF07dr1jnIZjUaqV69uldcsbs/X11f+MbBT8t7YL3lv7Je8N7Z1uyM3V1h1wHnrrbf44YcfCn995ajM6tWr6dSpEwCHDh3CYrEU7jN69GgyMjL45z//yaVLl2jTpg1//PEH5cqVK9zn448/xs3NjUceeYSMjAweeOABZs2ahclksubLEUIIIYSDsMk6OEJAydYvELYl7439kvfGfsl7Y9/s6mPiwrl5enry9ttvF/lEm7AP8t7YL3lv7Je8N/ZNjuAIIYQQwunIERwhhBBCOB0ZcIQQQgjhdGTAEUIIIYTTkQFHCCGEEE5HBhwhhBBCOB0ZcESpmj59OsHBwXh5eREWFsb69etvuu+CBQvo1q0blStXxtfXl/DwcJYvX27DtK6lJO/NtTZu3IibmxvNmze3bkAXVtL3Jisri3HjxlGzZk08PT2pU6cOM2fOtFFa11LS92b27Nk0a9YMb29vAgMD+cc//sGFCxdslFYUoQlRSubOnau5u7tr3377rXbgwAHt5Zdf1nx8fLRTp07dcP+XX35ZmzJlirZ161bt8OHD2tixYzV3d3dtx44dNk7u/Er63lyRnJys1a5dW4uMjNSaNWtmm7Au5m7emz59+mht2rTRVqxYoZ04cULbsmWLtnHjRhumdg0lfW/Wr1+vGY1G7ZNPPtGOHz+urV+/XmvSpInWr18/GycXmqZpMuCIUtO6dWtt+PDhRe5r2LChNmbMmDt+jsaNG2sTJkwo7Wgu727fm8GDB2tvvPGG9vbbb8uAYyUlfW+WLl2qmc1m7cKFC7aI59JK+t588MEHWu3atYvc9+mnn2rVq1e3WkZxc3KKSpSK7Oxstm/fTmRkZJH7IyMj2bRp0x09R35+PqmpqVSsWNEaEV3W3b4333//PceOHePtt9+2dkSXdTfvzaJFi2jVqhVTp06lWrVq1K9fn1GjRpGRkWGLyC7jbt6biIgITp8+zZIlS9A0jXPnzjF//nx69uxpi8jiOlYt2xSuIykpiby8PPz9/Yvc7+/vT0JCwh09x4cffsjly5d55JFHrBHRZd3Ne3PkyBHGjBnD+vXrcXOTfyas5W7em+PHj7Nhwwa8vLxYuHAhSUlJ/POf/+TixYtyHU4pupv3JiIigtmzZzN48GAyMzPJzc2lT58+fPbZZ7aILK4jR3BEqTIYDEV+rWlasftuZM6cOYwfP5558+ZRpUoVa8VzaXf63uTl5fHYY48xYcIE6tevb6t4Lq0k/9/k5+djMBiYPXs2rVu35qGHHuKjjz5i1qxZchTHCkry3hw4cICXXnqJt956i+3bt7Ns2TJOnDjB8OHDbRFVXEd+NBOlws/PD5PJVOwnm8TExGI/AV1v3rx5DB06lJ9//pmuXbtaM6ZLKul7k5qayrZt29i5cycvvPACoH9T1TQNNzc3/vjjD7p06WKT7M7ubv6/CQwMpFq1apjN5sL7GjVqhKZpnD59mnr16lk1s6u4m/dm8uTJtGvXjldffRWA0NBQfHx8aN++Pe+++y6BgYFWzy2ukiM4olR4eHgQFhbGihUrity/YsUKIiIibvq4OXPm8NRTT/Hf//5XzlNbSUnfG19fX/bu3cuuXbsKb8OHD6dBgwbs2rWLNm3a2Cq607ub/2/atWvH2bNnSUtLK7zv8OHDGI1GqlevbtW8ruRu3pv09HSMxqLfVk0mE6Af+RE2pu76ZuFsrnykcsaMGdqBAwe0ESNGaD4+PtrJkyc1TdO0MWPGaEOGDCnc/7///a/m5uamffHFF1p8fHzhLTk5WdVLcFolfW+uJ5+isp6Svjepqala9erVtUGDBmn79+/X1q5dq9WrV0975plnVL0Ep1XS9+b777/X3NzctOnTp2vHjh3TNmzYoLVq1Upr3bq1qpfg0mTAEaXqiy++0GrWrKl5eHhoLVu21NauXVv4e08++aTWsWPHwl937NhRA4rdnnzySdsHdwEleW+uJwOOdZX0vYmJidG6du2qlSlTRqtevbo2cuRILT093capXUNJ35tPP/1Ua9y4sVamTBktMDBQe/zxx7XTp0/bOLXQNE0zaJocNxNCCCGEc5FrcIQQQgjhdGTAEUIIIYTTkQFHCCGEEE5HBhwhhBBCOB0ZcIQQQgjhdGTAEUIIIYTTkQFHCCGEEE5HBhwhhBBCOB0ZcIQQQgjhdGTAEUIIIYTTkQFHCCGEEE7n/wGL376GZlB6KAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "plt.plot(g.x[g.ilo:g.ihi+1], a[g.ilo:g.ihi+1,a.g.jc])\n", "plt.plot(g.x[g.ilo:g.ihi+1], b[g.ilo:g.ihi+1,b.g.jc])\n", "print (a.g.dx)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Coarsening and prolonging" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "we can get a new `ArrayIndexer` object on a coarser grid for one of our variables" ] }, { "cell_type": "code", "execution_count": 23, "metadata": {}, "outputs": [], "source": [ "c = d.restrict(\"a\")" ] }, { "cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.65328 0.65328 -0.65328 -0.65328\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.65328 0.65328 -0.65328 -0.65328\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.65328 0.65328 -0.65328 -0.65328\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.65328 0.65328 -0.65328 -0.65328\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\n", " ^ y\n", " |\n", " +---> x\n", " \n" ] } ], "source": [ "c.pretty_print()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "or a finer grid" ] }, { "cell_type": "code", "execution_count": 25, "metadata": {}, "outputs": [], "source": [ "f = d.prolong(\"a\")" ] }, { "cell_type": "code", "execution_count": 26, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m 0.22 0.55 0.86 0.99 0.99 0.86 0.55 0.22 -0.22 -0.55 -0.86 -0.99 -0.99 -0.86 -0.55 -0.22\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m\u001b[31m 0\u001b[0m \n", "\n", " ^ y\n", " |\n", " +---> x\n", " \n" ] } ], "source": [ "f.pretty_print(fmt=\"%6.2g\")" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "celltoolbar": "Tags", "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.13" } }, "nbformat": 4, "nbformat_minor": 4 }