{ "cells": [ { "cell_type": "markdown", "metadata": { "collapsed": true }, "source": [ "# 06 - FWI with total variation (TV) minimization as constraints" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Introduction\n", "\n", "There is a lot of research on regularization to improve the quality of the final result beyond the simple box constraints implemented in the previous tutorials. In this tutorial we look at how more advanced FWI techniques such as [total variation denoising]( https://en.wikipedia.org/wiki/Total_variation_denoising) applied as a constraint can be implemented using [Devito](http://www.opesci.org/devito-public) and [scipy.optimize.minimize](https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.minimize.html) using [scikit-image](http://scikit-image.org/). This is a variant of the FWI with TV constrains algorithm described in [Peters and Herrmann 2017](https://doi.org/10.1190/tle36010094.1).\n", "\n", "[Dask](https://dask.pydata.org/en/latest/#dask) is also used here to speed up the examples.\n", "\n", "This tutorial uses the same synthetic datasets and model setup as the previous two tutorials, so check back if you get lost on parts of the code specific to Devito, SciPy.optimize or Dask. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setting up (synthetic) data\n", "We are going to set up the same synthetic test case as for the previous tutorial (refer back for details)." ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": true }, "outputs": [], "source": [ "#NBVAL_IGNORE_OUTPUT\n", "from examples.seismic import Model, demo_model\n", "\n", "# Define the grid parameters\n", "def get_grid():\n", " shape = (101, 101) # Number of grid point (nx, nz)\n", " spacing = (10., 10.) # Grid spacing in m. The domain size is now 1km by 1km\n", " origin = (0., 0.) # Need origin to define relative source and receiver locations\n", "\n", " return shape, spacing, origin\n", "\n", "# Define the test phantom; in this case we are using a simple circle\n", "# so we can easily see what is going on.\n", "def get_true_model():\n", " shape, spacing, origin = get_grid()\n", " return demo_model('circle-isotropic', vp=3.0, vp_background=2.5, \n", " origin=origin, shape=shape, spacing=spacing, nbpml=40)\n", "\n", "# The initial guess for the subsurface model.\n", "def get_initial_model():\n", " shape, spacing, origin = get_grid()\n", "\n", " return demo_model('circle-isotropic', vp=2.5, vp_background=2.5, \n", " origin=origin, shape=shape, spacing=spacing, nbpml=40)\n", "\n", "\n", "from examples.seismic.acoustic import AcousticWaveSolver\n", "from examples.seismic import RickerSource, Receiver\n", "\n", "# This is used by the worker to get the current model. \n", "def get_current_model(param):\n", " \"\"\" Returns the current model.\n", " \"\"\"\n", " model = get_initial_model()\n", " model.m.data[:] = np.reshape(np.load(param['model']), model.m.data.shape)\n", " return model\n", "\n", "# Inversion crime alert! Here the worker is creating the 'observed' data\n", "# using the real model. For a real case the worker would be reading\n", "# seismic data from disk.\n", "def get_data(param):\n", " \"\"\" Returns source and receiver data for a single shot labeled 'shot_id'.\n", " \"\"\"\n", " true_model = get_true_model()\n", " dt = true_model.critical_dt # Time step from model grid spacing\n", "\n", " # Set up source data and geometry.\n", " nt = int(1 + (param['tn']-param['t0']) / dt) # Discrete time axis length\n", "\n", " src = RickerSource(name='src', grid=true_model.grid, f0=param['f0'],\n", " time=np.linspace(param['t0'], param['tn'], nt))\n", " src.coordinates.data[0, :] = [30, param['shot_id']*1000./(param['nshots']-1)]\n", "\n", " # Set up receiver data and geometry.\n", " nreceivers = 101 # Number of receiver locations per shot\n", " rec = Receiver(name='rec', grid=true_model.grid, npoint=nreceivers, ntime=nt)\n", " rec.coordinates.data[:, 1] = np.linspace(0, true_model.domain_size[0], num=nreceivers)\n", " rec.coordinates.data[:, 0] = 980. # 20m from the right end\n", "\n", " # Set up solver - using model_in so that we have the same dt,\n", " # otherwise we should use pandas to resample the time series data. \n", " solver = AcousticWaveSolver(true_model, src, rec, space_order=4)\n", "\n", " # Generate synthetic receiver data from true model\n", " true_d, _, _ = solver.forward(src=src, m=true_model.m)\n", "\n", " return src, true_d, nt, solver\n", "\n", "# Define a type to store the functional and gradient.\n", "class fg_pair:\n", " def __init__(self, f, g):\n", " self.f = f\n", " self.g = g\n", " \n", " def __add__(self, other):\n", " f = self.f + other.f\n", " g = self.g + other.g\n", " \n", " return fg_pair(f, g)\n", " \n", " def __radd__(self, other):\n", " if other == 0:\n", " return self\n", " else:\n", " return self.__add__(other)\n", "\n", "from devito import Function, clear_cache\n", "\n", "# Create FWI gradient kernel for a single shot\n", "def fwi_gradient_i(param):\n", " # Need to clear the workers cache.\n", " clear_cache()\n", "\n", " # Communicating the model via a file.\n", " model0 = get_current_model(param)\n", " src, rec, nt, solver = get_data(param)\n", " \n", " # Create symbols to hold the gradient and the misfit between\n", " # the 'measured' and simulated data.\n", " grad = Function(name=\"grad\", grid=model0.grid)\n", " residual = Receiver(name='rec', grid=model0.grid, ntime=nt, coordinates=rec.coordinates.data)\n", " \n", " # Compute simulated data and full forward wavefield u0\n", " d, u0, _ = solver.forward(src=src, m=model0.m, save=True)\n", " \n", " # Compute the data misfit (residual) and objective function \n", " residual.data[:] = d.data[:] - rec.data[:]\n", " f = .5*np.linalg.norm(residual.data.flatten())**2\n", " \n", " # Compute gradient using the adjoint-state method. Note, this\n", " # backpropagates the data misfit through the model.\n", " solver.gradient(rec=residual, u=u0, m=model0.m, grad=grad)\n", " \n", " # Copying here to avoid a (probably overzealous) destructor deleting\n", " # the gradient before Dask has had a chance to communicate it.\n", " g = np.array(grad.data[:])\n", " \n", " # return the objective functional and gradient.\n", " return fg_pair(f, g)\n", "\n", "import numpy as np\n", "from distributed import Client\n", "\n", "# Dumps the model to disk; workers will pick this up when they need it.\n", "def dump_model(param, model):\n", " np.save(param['model'], model.astype(np.float32))\n", "\n", "def fwi_gradient(model, param):\n", " # Dump a copy of the current model for the workers\n", " # to pick up when they are ready.\n", " param['model'] = \"model_0.npy\"\n", " dump_model(param, model)\n", "\n", " # Define work list\n", " work = [dict(param) for i in range(param['nshots'])]\n", " for i in range(param['nshots']):\n", " work[i]['shot_id'] = i\n", " \n", " # Start dask cluster\n", " client = Client()\n", " \n", " # Distribute worklist to workers.\n", " fgi = client.map(fwi_gradient_i, work)\n", " \n", " # Perform data reduction.\n", " fg = client.submit(sum, fgi).result()\n", " \n", " # Shutdown cluster\n", " client.close()\n", " \n", " # L-BFGS in scipy expects a flat array in 64-bit floats.\n", " return fg.f, fg.g.flatten().astype(np.float64)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## FWI with L-BFGS-B\n", "Equipped with a function to calculate the functional and gradient, we are finally ready to define the optimization function." ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": true }, "outputs": [], "source": [ "from scipy import optimize\n", "from skimage.restoration import denoise_tv_chambolle as denoise\n", "\n", "# Define bounding box constraints on the solution.\n", "def apply_box_constraint(m):\n", " # Maximum possible 'realistic' velocity is 3.5 km/sec\n", " # Minimum possible 'realistic' velocity is 2 km/sec\n", " return np.clip(m, 1/3.5**2, 1/2**2)\n", "\n", "# Many optimization methods in scipy.optimize.minimize accept a callback\n", "# function that can operate on the solution after every iteration. Here\n", "# we use this to apply constraints and to monitor the true relative\n", "# solution error.\n", "relative_error = []\n", "\n", "def fwi_tv_callbacks(x):\n", " # Apply boundary constraint\n", " x.data[:] = denoise(x.reshape(181, 181), weight=5.0e-3).flatten()\n", " x.data[:] = apply_box_constraint(x)\n", " \n", " # Calculate true relative error\n", " true_x = get_true_model().m.data.flatten()\n", " relative_error.append(np.linalg.norm((x-true_x)/true_x))\n", "\n", "def fwi(model, param, ftol=1e-6, maxiter=20):\n", " result = optimize.minimize(fwi_gradient,\n", " model.m.data.flatten().astype(np.float64),\n", " args=(param, ), method='L-BFGS-B', jac=True,\n", " callback=fwi_tv_callbacks,\n", " options={'ftol':ftol,\n", " 'maxiter':maxiter,\n", " 'disp':True})\n", "\n", " return result" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We now apply our FWI function and have a look at the result." ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " fun: 225.78709341897513\n", " hess_inv: <32761x32761 LbfgsInvHessProduct with dtype=float64>\n", " jac: array([ 0., 0., 0., ..., 0., 0., 0.])\n", " message: b'CONVERGENCE: REL_REDUCTION_OF_F_<=_FACTR*EPSMCH'\n", " nfev: 18\n", " nit: 9\n", " status: 0\n", " success: True\n", " x: array([ 0.16000003, 0.16000003, 0.16000003, ..., 0.16000003,\n", " 0.16000002, 0.16000002])\n" ] } ], "source": [ "#NBVAL_SKIP\n", "\n", "# Change to the WARNING log level to reduce log output\n", "# as compared to the default DEBUG\n", "from devito import configuration\n", "configuration['log_level'] = 'WARNING'\n", "\n", "# Set up inversion parameters.\n", "param = {'t0': 0.,\n", " 'tn': 1000., # Simulation lasts 1 second (1000 ms)\n", " 'f0': 0.010, # Source peak frequency is 10Hz (0.010 kHz)\n", " 'nshots': 16} # Number of shots to create gradient from\n", "\n", "model0 = get_initial_model()\n", "\n", "# Apply FWI with TV.\n", "result = fwi(model0, param, ftol=1e-6, maxiter=100)\n", "\n", "# Print out results of optimizer.\n", "print(result)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAdkAAAGDCAYAAABnUmqTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzsvXn8XFV9//98k4SsSAyRLUUWoQjUlmpqqUtZ+uML9ssi\nX3GrILVuXVyoX2qluGBApUpxeWircakotEVQWZSKIqCt/aaKK4JGwKQQQDCEQHaS8P79cWYyZ27m\n3Dlz7535fJJ5PR+PecyZe88598ydO3Pmvt7v836buyOEEEKI5tllogcghBBC7KxokhVCCCGGhCZZ\nIYQQYkhokhVCCCGGhCZZIYQQYkhokhVCCCGGxMgnWTPbz8yuMrNHzewxM/uSmT01s+0MM/uAmT1g\nZhvM7P+Z2R8Oe8xCCCFEFUY6yZrZLOAm4OnAWcCZwCHAzWY2O6OLTwOvBd4JnAQ8ANxgZkcOZ8RC\nCCFEdWyUwSjM7M3AJcCh7n5Xa9uBwJ3AW939kpK2vwP8CPgzd//n1rapwO3AUnc/ZdjjF0IIIQZh\n1HLxKcCS9gQL4O7LgO8Ap2a03QxcEbXdAvwbcIKZTW9+uEIIIUR1Rj3JHgH8tMf224HDM9ouc/f1\nPdruChxcf3hCCCFEc4x6kp0HPNJj+yrgyTXatvcLIYQQk4apEz2AYWNmrwNeB8CU2c/iSU8PO+K/\nF1ZslOqswgByTd6eUc6tl9tm0DqD1MtpP8i+JtuIZqnyvchpX6XfsjaWUc7tu6x9Tr0q48wl9zv2\nROt53XJ800oDONhsO6lwUB6AG9z9xJrd7DSMepJ9hN53rKm71GLb/RNtoXNH24W7LwYWA9i8hc4J\nt4YdM6JKxbMwNaOcy5bMfTllgI0V2gx6nCpjLquXs73fvibbiGap+yuSal/Wb902ud/rnDZlvx+p\n35kqx8wl9zvW/i25YeG2TRuAv6pwyJi3w/yaXexUjHqSvZ1gWy1yOHBHRtvTzGxWwS57OPA4cFfv\nZhFG73c8zEk2btPEhND0l64pmp5wq9QbVvtxYljXV1m/cfuy70tqX5O/YlV+C4bZpmxsMbnnsFd/\ndZUIUcqobbLXAkeZ2UHtDWZ2APDc1r4yrgOmAS+O2k4FXgp83d03NT1YIYQYJ4zwI1vnIboZ9ST7\nSWA5cI2ZnWpmpwDXAPcCn2hXMrP9zWyLmb2zvc3df0hYvvMhM3uNmf0RYfnOgcC7RvgehBBip6Qt\n9tV5iG5Gek7cfZ2ZHQd8EPg84TP9JnC2u6+Nqhowhe3/BLwKeA9wITAX+DFworv/YODB1LWJTNTV\nlHPcMpmtChMht45KShbd1D2fqWstR7Ys1iv2lSMRN33tp/oalVycO/4qZqkecnH7TlY0x8inCne/\nB3hRnzrL6WEpcPcNwFtaDyGEEGJSo7t7IYQQQNo3VFRH51MIIQQguXgYaJKFanaPyWyTzR3bZLO1\n5oxHa2knntylJLltmqTuMrsyW2lum5y1sXGdfv01SS+7dsEmq0mhWXQ+hRBCALqTHQYjT9ouhBBC\njAu6k+1HlVBoMU1Hpollpo2JOlUiw6TqTAa0nGdykbscp6xNzrKdJpbjDLo0punlOLlhFUclF/dZ\n6iO5uHl0PoUQQgCSi4eBJlkhhBCAJtlhoEm2F0168I5K+sk95qhk7VEF9R+WR7JIU9ccMUi9YVHF\nDJSSfnM9hXPbNEkFyV2TQrPI8UkIIYQYEvrTIoQQApBcPAzGa5J1OvLJMCXEiV54X0UurhJgoC7D\nDCzR5JjHVW4eZrD/UTGoR3Gu129KOi7bV9ZmWEk8ypIs9EDexc2j8ymEEALQnewwkE1WCCGEGBK6\nkxVCCAFILh4GOp+wvZ2iii2pbjSaKqSiP9W1yea6/Y/KVqnlPBNDlWD/w4zelGpfZTlO3cD9uW2q\nRHxqkgGv91HIxWZ2OnAG8CxgPnAP8CXgve6+pk/bpwIXAMcCTwHuBb4AvM/d10X1dgH+Fng9sDew\nFFjk7l9s/A31QZOsEEIIYGR3sucA9wHnAiuAI4HzgWPN7Dnu/kTPsZnNBm4k/A94B2Fy/j3g3cAh\nwEuj6he0jnMe8H3gZcCVZnaSu18/hPeURJOsEEIIYGSOTye7+6+j17eY2SrgUuAY4KZEu+cSJtMT\n3f2G1rabzWwecI6ZzXL39Wa2J2GCvcjdL47qHQxcBGiSHQmjkg1HFfWmTL7KkZJzJbcqY25oecFA\nfSvvbLMMU/rNOWbuceoG+C9bjpMr/ebUa0Iu3kGv18IE2+Z7recFJU13bT2vLmxfTXDibWfGPaFV\n97JCvcuAz5jZge6+LH/E9ZB3sRBCCKAjF9d5VOTo1vPPSurcCNwJvN/MDjezOWZ2HPBm4OORTfYI\nYBNwV6H97a3nw6sPc3DG905WCCFEFw3JxfPN7Nbo9WJ3X5w8ptkCYBFwo7vfmqrn7hvN7HnAF+lM\nmACfAt4QvZ4HrHZ3L3SxKto/MsZvkm1LLGXy15ZEuWlP4ybPfhWvykHrlO1rIsdoXYblXbyDynKV\nqJJAIsejeJi/NFWu/bp5XlMyclmb3O/YMOn12xZNRQ1NsivdfWFORTObA1zTGs2r+tSdAVwB7AWc\nSXB8ejbwzlb7v6gx5qExfpOsEEKIJKOaFMxsJnAdcBBwtLuv6NPk1QTHqEPcvS0Ff9vMHgUWm9nH\n3f3HwCPAXDOzwt1s+w52FSNENlkhhBAjxcymAVcBC4E/dvfbMpo9gyADF22t3209H9Z6vh2YDjyt\nUK9ti71j8BFXZ7zuZFMJAiabHNh0YIuUzJXyOs7tK1dCbJKm5d5R5b2drOR+bhPhXVw2hipeu1Wk\n3yrexcOSi5u49nPk4rqfZ78kBCFQxOXAccBJ7r4ks+dfEe5QDy5MtL/fer6v9fw1YDPwCsIa2jZn\nAD8dpWcxjNskK4QQIokZTB3yJAt8DHgx8B5gnZkdFe1b4e4rzGx/4G5ClKZFrX2fBd4CXG9m7yHY\nZBcSAlN8H/gOgLs/ZGaXAOea2RrgB4RAFccBp9R8dwOjSVYIIQQQJtlpU4Z+mBe0ns9rPWLeTYj+\nZMAUIpOmuy9vTcjnAxcSQjLeCywG3lOIFHUesJawvKcdVvEl7v6Vht9LXzTJCiGEGBnufkBGneV0\ngkvE2+8AXpLRfithIr5w8BE2y/hOsrl2urrLeWJyIynltsmtN+iynbrLgcr2TXRSgWHaZ3dEmk7A\n3uR3pMp12HTEp0GXAA3S97AY9Lcttsk2IReLLnQ6hRBCAA05PokudDqFEEIE2pZQ0RjjNcnGS3hy\nIz6lGGZ+zNyx5EpRgy7byU0wkLuso4rUmEvdKE/jKBHnMhEScW693Gs3tQQnp1x8XVdiTtUp0uSy\nskF/25S1vXEUjEIIIYQYEvrPIoQQIqA72cYZ39OZK6nkehfX9Q5O0bRMlyOzxfJwWZu6UaKq1ksx\nrAQBopsmTCWDeryX7RuVd3HdiE+5VPke5HrT50S6G99ZYSjodAohhAjI8alxNMkKIYQISC5unPE6\nnbF3cUwVuXgyBEmv2z7H67j4OldmS0nZVeTvutLxUCXiYl7oHOoOqMpFsF3wnMFo0jRRVq9se5Ny\nb653cd1jDpOc36ai6adPMArRPOM1yQohhEijO9nG0ekUQgjRQTbZRtEkK4QQIqA72cYZr9OZG/GJ\nxL7cepNt2U5O38NcOpHaXiVK1IQswSkzWqUGtHkYA+nR97TMNvFJrGmfbYJhLcepG7GpeMwqEZ/I\n2FclClmu78jGxPZUG9lkh8p4TbJCCCHS6E62cXQ6hRBCdJBNtlHGd5LNjd6UK91UkTqrRIapEng/\nJmfZTlFyq0tOggHIO2+5kbZqE2toxYOkpODcwdSVkmOJuMpasiFKx1WW46S2l0m3qXLdJTxVJOoi\ng34EUC3/careoNHs4ktdd7KNo9MphBAioEm2cZSFRwghhBgS4/WfpUo+2SoyTF1JN5e6x8mJ/lSs\nN6zA7MXXTebUzCYlERfl3bJ9gw6gTDpOeQ6XfdiDeh4X21eQj6tIv6l6VaIvNZkPNvc6TtUZJk3/\nTvVCd7KNo9MphBCigxyfGkWTrBBCiIDuZBtnfE9nlQAUVXLQjuoM53oxV9k+EXJxqty4N3GORFzm\nXZxq00QwilQf00rq5AanSFGzfV1TRe71kdtmUI/iKl7QVajrNVy2b9AAOvIuHipyfBJCCCGGxMgn\nWTPbz8yuMrNHzewxM/uSmT01o93vmdmnzexOM1tvZveY2eVmduAoxi2EEDs97aTtdR6ii5EKA2Y2\nC7gJ2AScRRAqLgRuNrPfdvd1Jc1fChwBfAS4DdgXeAdwq5kd6e73DnXwQgixsyO5uHFGfTpfCxwE\nHOrudwGY2U+AO4HXA5eUtH2/u58TbzCz7wDLWv2+s7FRNukaP6yA/oMwaAL1JiI+5RxzUiz7ybHD\nbii0SdleywbQZMKAJo3+RRts0UBXg7qB98uW4wwr6XrZtV8lQluTEZuaXsKTShAw5FnBzE4HzgCe\nBcwH7gG+BLzX3deUtDsfeFdi9yZ3nxHVXQ7s36Peae5+dbWRV2PUk+wpwJL2BAvg7stak+WplEyy\n7v5Qj23/Y2a/BhYMY7BCCCEa5xzgPuBcYAVwJHA+cKyZPcfdn0i0+xTwtcK22a1t1/aof0Or35il\n1YZcnVFPskcA1/TYfjvw4kE7M7PDgD2Bn9UclxBCiLZNdric7O6/jl7fYmargEuBYwgmxe1w9xWE\nSXkbZnYmYR67tEeTle6+pJER12DUk+w84JEe21cBTx6kIzObCnwc+DXw6eyGOfJN6qzUzS1bJalA\nFXLbl0V5SvWXK+OmZLey9zmSq7GojeVIxEWpd3OiXhW5uMoHn5sxIT7mrERfFT6EKjJ/lehNuRJz\nbsSnYeWGLaOKOSNX+t2YUZ6EEZ8KE2yb77WeB1UlzwIeJNy1Tkp25CU8HwWeA5zh7r0mbgDM7HVm\ndquZ3crmXp+tEEIIoDPJ1nlU4+jWc7YqaWb7AccCl7t7r78PJ7dWomwysyVm9sLKo6vBqCfZR+h9\nx5q6w+2JmV0EvA74M3f/elldd1/s7gvdfSHTnjLQYIUQYuwY8SRrZguARcCN7n7rAE3PIMxhvaTi\n64A3AicAryDc33/ZzM4YfIT1GLVcfDvBLlvkcOCOnA7M7Dzgb4E3uvvnGxxbh2F5+pUpe01KxFVk\n6TLpuK63cXzMGYntxfFUkfNqk/IaLvMuTsnFZd7Eubphjt0i94TkJjKIX1eI/lTFKzzXUzhH+q3i\nXZzrUTxMci6JYXkXN898M4sny8XuvrhXRTObQ/DT2QK8asDjvBL4obv/pLjD3d9YOM6XgSXAe4HL\nBjxOLUZ9J3stcJSZHdTeYGYHAM+lt3dYF2b2JsK62vPc/aNDGqMQQownzQSjWNlWD1uP1AQ7k3DH\neRBwQsuxKW+YZs8Gnk7vu9jtcPetwJXAfma2T+5xmmDUk+wngeXANWZ2qpmdQvgXcy/wiXYlM9vf\nzLaY2TujbS8DPkRw177JzI6KHoeP9F0IIcTOyIhssmY2DbgKWAj8sbvfNuBIzyLIM/8yYLuRM1K5\n2N3XmdlxwAeBzxM+0m8CZ7v72qhq+/9U/CfgxNb2E1uPmG8RXL/zaSJAd069svbD8iKu229RPovl\n41yZLUciLrbPlZVjBn6vufpZWbD/DYlybg7aKjTpil43iUCBup7CdeXiKscZpkdxzm9G3d+S3HqD\nXnYj8C42s12Ay4HjgJMGXWZjZrsCLwP+PeGp3KvNVELUwHvc/YEBh1yLUdtkcfd7gBf1qbOcQrgZ\nd/9T4E+HNS4hhBCMYp3sxwhxEd4DrDOzo6J9K9x9hZntD9wNLHL3RYX2JxGcZXtKxWb28lad6wlB\nL/YG/gp4JvDyJt9IDiOfZIUQQow1L2g9n9d6xLybEKWpl5rZ5ixCbIWvJPpfBuxDiCA4D1gH3Aqc\n6O4jX087XpOs05FPqqz7b9prr4rqN6pAFTEpWTglI5e1qSKTNepdXBZYYn1UTknCZfvK5OK6wShi\nibeuXByPpdh+Zv/mTQSjGJZcXNZmQjzWI5o2N9X9zZmgfLLufkBGneUkgme7+6l92i4hSNGTgvGa\nZIUQQqRRFp7G0ekUQggR0CTbODtyWEUhhBBiUqP/LFBtiUiVfLLDTBAwqk+yyehPTUR8qm1bSw2o\nik22LEpU6phlpC6Q+Di5y3HiemXtM/rekWyyg9phq/hqVCHXvpob7L+pJTwwCu/isUKTrBBCiIDk\n4sbR6RRCCBHQJNs443s6c+XZYUZfaTJBQBlNxpmPqSIdVzmHZUuFssYdr1EoHjQV2SklHRdf50jH\nvY47KPEbjZfZlPWbkojjcabyzJYcvmxfk9LvnAptcpeSDTMhR4q6S3iKiTuqSMQ5Y5Vc3ChyfBJC\nCCGGxPjeyQohhOhGcnHj6HT2IiWj1pWL6x6/SF0Pydz2OfJV3WQBuX0XJbOBKYvElCMD59Ybplyc\nW29DYntKRoasceZ6F0+GBAFNRnaq+13O7WtUCQImKOLTuKHTKYQQItCOGCwaQzZZIYQQYkjoThaq\nBYkYptxTdpwqEvGwPuUyGTcnH2xxe45EXEkOLPsQUgH+y6TfKnJxajxlpN5QWaCLmFRSgbJcuRlD\nqSIXDzMYxTCDlVT5/uf01fTqhLrWiDaSixtHp1MIIUQHzQqNotMphBAioDvZxtHpFEIIEZDjU+OM\n1yRbJWl7TNnSmiZtKrnHmaz2Wchb3lM8Hyk7bq7NrTYpW2XRBpqql4oYVXacupRdVFXGlkHZZ5Cy\nlZbZV+cktjdhkyVzX5tc/4zcPnJ/C1KJAHJ/P3KOX9ZGDI3xmmSFEEKkkVzcODqdQgghOmhWaJTx\nPZ25kmxO++Lr3CUro5Jucj7lKgGFqiwhKDvvqWU7ZRF9yNy3jbIIRymptYm1E1Uk4pw2ZdptlYsy\n4zBln9uolvAMmie2bF/TvwU59XKjnU3EEh7ZZBtlfCdZIYQQ3UgubhxFfBJCCCGGhP6zQH5UpbI2\nw/IurhKNqowqn3jdNimv4bLzkeulqit49OR6F1eRfqsE+68b1aluJKcidU1HShCwU6HTKYQQooNs\nso2iSVYIIURAd7KNo9PZj5TcUnbm6kg3/eoNi6avhBxP4Vy5uEpe0i4JMJVLtVeHvertDF+TXH01\n8V4nIjdsrsRMxvYyhvl9G6b0q8ASOwQ7w6+HEEKIJtCdbOPIu1gIIUSgPcnWefQ7hNnpZna1md1r\nZhvMbKmZvc/MduvT7nwz88RjY6HuLmZ2rpktN7ONZvZjM3vRYCejGfSfRQghRIfhOz6dA9wHnAus\nAI4EzgeONbPnuPsTiXafAr5W2Da7te3awvYLWsc5D/g+8DLgSjM7yd2vb+JN5KJJthc5S2PK7CO5\nQXfquvrX/fRG1T4nqTekk7PXXQqy1qIXMwsdTEvs25CoUzxQXC6L0BT3kRv9KWUXnpYoQ/d7mJmo\nF2+flW6/MyRgH3R5Tm7Ep1ElYG/ap6PfORiNXHyyu/86en2Lma0CLgWOAW7q1cjdVxAm5W2Y2ZmE\nEV8abduTMMFe5O4XtzbfbGYHAxcBI51kJRcLIYQYGYUJts33Ws8LBuzuLOBB4IZo2wnArsBlhbqX\nAc8wswMHPEYtdCcrhBAiMHGOT0e3nn+W28DM9gOOBT7k7vE9+hHAJuCuQpPbW8+HA8sSfU4FXgic\nCBwF7EuQd1YCS4FvAVe4+y9yx6lJtikmQiJK9TWqqE65VMknG+cY3ZjYXrZvY6K8pehbEcvCOXli\ne73uRfHNVQn2H5OSe4vyd/z6SVF5t4ztEH5lW8yNNs9JlKHZ3LBlbeou28mJljbMiE+pwP9lbepI\nv7l1itS3yc43s1uj14vdfXGqspktABYBN7r7ral6PTiDoMZeWtg+D1jt7l7YviraXxzDTOAtwJuA\n+cAvgB8C3yT8SMwDDgT+L3C+md0C/J27/3e/QWqSFUIIEWjmTnaluy/MOpzZHOAawt+BVw14nFcC\nP3T3nwzYrhd3AQ8DFxLuVB/qVcnMDPhDwgR/k5md7e6fLOtYk6wQQojACOXi1t3jdcBBwNEtx6bc\nts8Gng6c3WP3I8BcM7PC3Wz7DnZVjzZvcvcv9jtuq79vAd8ys3cB+/drM16TrNNbPsmN3lRGlTyr\nVVJ8pvoe1SdZN0FBWcQnEvvmJrZDtwSXqheXV8eexsCWPaMXuR688evHEttjGbrYdy45EnEs/UK3\n/BurYnskthfe29wBy8XXKVm5TGKu4l3cJHU9eKtIv12rOjPb7EQRnsxsGnAVsBA43t1vG7CLswg2\nmH/pse92YDrwNLrtsoe3nu8oNsiZYHu0uR+4v189eRcLIYQIjCYYxS7A5cBxwAvdfclAQzTblbDu\n9d8TnspfI0zAryhsPwP4qbv3dHrqc8zfNLNTzWzvQduO152sEEKIUnz4wSg+BrwYeA+wzsyOivat\ncPcVZrY/cDewyN0XFdqfRJBjig5PALj7Q2Z2CXCuma0BfgC8lDCpn9JvcGb2YWCau/9l6/WpwJWE\n+fJRM/sjd/9B7pvVJAvV8smW9VElgMWoctBWoUmvzjJP4bpe0angCfFxVhfax/Lx2khS3bgHjVI3\n52ld2TB1PoqfQY5EPL9Cm+JxUl7Iw5SLq5y3yRyMYgiyshtsHf6s8ILW83mtR8y7CdGfjODn3Ett\nPYtgV/1KyTHOA9YCbwb2Jiy/eYm7l7Vp878J3s5tLiDcHb8LuLg1xpMz+gE0yQohhBgh7n5ARp3l\ndK0n69p3akb7rQRP4QsHHB7APsBy2La86LeA17r7D83sQ0CpN3ERTbJCCCECo7mTnexsJMREhhAk\nYw2diFRr2N7bsBSdTiGEEECQi7dMqesPm4rvv8PwA+AvzWwZ8JfAN6KkBQcADwzSmSbZXlSJBpOz\nhGeiEjUPutSmrE7d6DplS3hS9ri4XGbbi22Aa6NyWaSd1PFTNszi61Q5N5FB2fWVWvIRl+P3WXyd\nKheXj6RILbNp+nzkJGMv2zeqpS1VknjU/S0YJu3jRCtJ3YytU+tOC4/XbD/hvIOQROB2whq9N0T7\nXkjnrjYLTbJCCCG2sXXK8N2LJzPuvsTMDgAOA5a6e+wu+RlCyMVsNMkKIYQYa8zseuDLwLXu/qC7\nPwZsF5fY3Yt5a/sy0CTbWogbZyVY5u47vDbQCDkS86iiP5XJsE1GiaoiHceUSWYpebKKfB73VYxW\nFC1H2WX+um3lPfZ6OGrySKFJvK/zJ3c31mwrz2R9V5vpkYS2K5syBg2PM31beX2U93VNdELWFAL8\nr47e4GqevK388KbOkqRHfxUtT1pZiPgU/2ePJeZcmT03d3DdZWFNBvivIuPWDfafex03LTH3qecY\nW0eQtX0Scj9hqc4/mdl3gauBqwfJtpOi7yVtZguB1xBy9D21sPtxM/se8K/AZe6+ptheCCHEjoFj\nbBnDSdbdX9MK/v9c4FTCnPc+M1tKmHC/7O4D2WLbJCfZ1uR6MSHjwG2Ehb8/BH5Nd+qf3ydkm7/I\nzN4P/IO757pWCCGEmERsHVMrYiv4/3+2Hn9jZr9FcHQ6FXibmd0PXEuQlW8u5LBNUnY2v0VYdPsX\n7l6aSNfMZrQG8lZChI4Lcg6+Q1BFXs1tMxGRYapIxznRipr2SM6tM6fjGjltbkdImbtHR/fciwe3\nlfctxPPeJ3q9H/f2rLcX3Vmv9oz6S0rHm7pFnVnrOssaLDPJw+aOWsya3Tt6aywRP1wIv/RwlAjg\nQfbaVr5/+j7byg/sv29ne1QGuJ/O6wfpJE94eGvnOKtXdmvuT6ztSNls7Bk/oNloS2X1itfHoFGR\nJioS07CSAgzYfozl4u1w958CPwUubAWlOI0wz30VWAeRPaaEsp+1p7n7rzIHsxG4ArjCzPbqV18I\nIYTYUXD3+4CPAh81s7mE0ItZJCfZ3Am2R7sH+9cSQggx2dCdbDdmtifbu+3h7pfn9lFJfG+lKioe\nNCvMh5ntB3wQOJ4Qm/JG4Gx3v2fAMbwNeB/wHXd/3iBtgXxvybrtR+VRWKVNXa/j3GDuOWWAqdGq\n+BkdD9wZczqeunN3747wv29C7o3LB3P3tvIBISRp9LqT9erArZ19T7pzc6dS8cqMX8fqc5wK+tFC\nm3VROVMunja7U543u+PmMG+PTnn/PQuZvmL1N3ZTPLBTvG/PTj7ZZRzQ1Xx5VHF5tG/5lE753r32\n62pz/16dg3Z5N6/rlNfHkjLwxJboh3xjpItvieTmKp7kdb9juQFBUttz+x5mwImafY/7JGtm84CP\nAC8Cdk1Uyz5JWT+nrQz27yKkJ/qNHu08py8zmwXcBGwiZFJwQgDnm83st919XVn7qJ+DgLdDwVgm\nhBCiMuPqXVzgU8D/BywGfk7NEFa59yz/SEiAex3wbzUO+lrgIOBQd78LwMx+AtwJvB64JLOffyIk\n/T0UBdQQQgjRHMcBb3b3f26is9wJ6hTgHHf/SM3jnQIsaU+wAO6+zMy+Q/Da6jvJmtmfAM8EXg58\naaCjO4NLTlUk1boeyaOSi5uUiMuCDcyIZOCpnYNOm7Ohq8msWBae3ts7uCj3Pi2Sgg9m22XFb7J0\nW/nQKAraU35eCPYb+83fGZV/GZWX0U1CIl4fScRrCppMHJoidQkWP4JYYJ0ZKapP2j3aMY9uYrfD\nWDqO5OIFB3YGveCwWOOG5x7WyUV917zf2FZeym9uK9/NwV1t7uJp28qx3Hz/7I5H88Ozu72gV2/q\nLSVv3hipc7GMDPU9l6uYZ1IScdNt6no+N0SwyY79fctqoJJPUi9y0y1sovvnqCpHEFyii9wOHN6v\nsZk9mWDPfau7r+pXXwghxGBsZUqtx07AxwiqayPk/mX5LPAy4Bs1jzcPCjHqAqvIW3P0AUJw5s/m\nHtDMXge8DoBpxYBVQggh2si7GNz9A2b2D2Z2G8Extzhnubtnx4LInWTfQYjp+HXghh4Hxd0/k3vQ\nKpjZ84FXAs9sRebIwt0XEwzY2MyF2e2EEGLcCBa18Z5kzexE4M8JMfqP6FHFGSDgUu4k+yyCPXVP\ngtdVr4N4NzB0AAAgAElEQVTmTLKP0PuONXWHG/MJ4NPAitZiYAjjn9J6vcHd8yKvj5KJsK+OwHYD\nlC/Hieywu0S21tjuOmd2d1SkOLJS9xKcjq01tsECHM4d28qHRnbYBcsia8JtUYOfF8aZssNGdlcv\n+LA/FHX9GL3L3dbm7o9gM3nEoftnRlf2zGg88wpjm9c5bcyK7bXxe4sFnfj9AxYZhQ55xopt5QOe\n3in/YvflXW3iZVT7Rvms745stcWlQg9N7xiPV07vRKlau64TzWr92q1dbZ4gEVmqSoKA1PelbDlO\nKqdvE0t4cu2wYhR8EPgRIY/syLyLPw48TNCp6xz0dnr/Mzgcol/L3hzWevx5j32PAH8NfKjiuIQQ\nQsjxCWB/gnfxD5voLPdsPh043d2vr3m8a4GLzewgd/8lQCs57nOBt/Vpe2yPbR8iLAp+I0S3O0II\nIQZGNlkg3MXu07dWJrmT7FJgdt9a/fkk4Rb8GjN7Ox1t+16CHAyAme0P3A0scvdFAO5+S7EzM1sN\nTO21rzGqLHOp0neV5QFVIktVOWa8PCe5hKfb3B1LxLvFgfsTS3Oge3lOaglOLAkDHLq18/pJt0VC\nbOzDHvvFF5fjJKI3rY9k2AcLy3FikTslERfl4lgizlUA49PbJR2XHOexaKx7ROW9YkPKukQZuiNV\nReVpnTwIHPGMWHuGOftFn29mft1Z0cin0JGFp87ulohj1sZRorZsF+muPzkSce5ynNx8sjnLeYpU\nkYgblJU1yXI28Bkz+7m7b5e4fVByp423Ae83s++6+/9UPZi7rzOz4wia9+cJYRW/SQirGC9gNMId\nau4SIyGEEDXRnSwQkt08GfgvM3uM3t7FT9u+WW9yJ9m3E5yefmFmv0gc9Oicjloxil/Up85ywkTb\nr69jco4phBBCZPIdgsraCLmT7Fa298scb6pofsU2OXJvrhRVdpwmg5SnojzN6HbsTknEsSdqMXpT\nLAUn5eJN3XLx7NuivBTxFRp7zd6TKENSIl4VyajFqCexXBzLtesT2yF9qmMZeVqiDnSf9jLpOem5\nHL23pnNR7r+lk6Rg6oEduTeWgYtMLdnXZsvs7juqLZFcvLErwUDJmRvUPJLrXVw34lORUXkR9zF/\nKXYxuPsZTfaXNcnqjlEIIcaDcfcuNrMj3P32kv2nu/tVuf1l2TzN7Df67M+SioUQQkxe2jbZYYZV\nNLPTzexqM7vXzDaY2VIze5+Z7da3cWh/mJldaWYro/ZvLtRZbmbe4/HCjEN8LTXnmdmLCMlpssn9\ny3KDmT3X3VcXd7QiMX0FyDpBY0eZPJOTiGCYASxSdYokPYo74mSc8xVgt+kdUTX2Io6DTGyf27Xz\n+sCoHOd8nb2skLY4Cr6QzPMaB2wo5Hn1RFD/WCIuSr/x65R3cVG2zQlGkfsR5DItUZ4aeQrvsV0Q\nkag8PbG96NgbvV4wvXPmNu3b+dweL6TljH+MN0UHircXf7C37t4Z7ONR8oAnpkbvrumEHDnexbnB\nKFJ1qo4th8kgUW/POcB9wLnACuBI4HzgWDN7TllucjNbSEiXegvwGsI3+hBgTo/qN7T6jVnao16R\n24Cvt+a8bf5HZnYa8K+EXLPZ5F6Sa4Gvmtkfufu2S8rMngdcT1j/KoQQYgdmRN7FJ7v7r6PXt5jZ\nKuBS4BjCJLodZrYL8Dngm+5+WrTr5sRxVrr7kgrjO73VZ3vO22BmpxLSvP6ju58zSGe5S2T+N7AH\ncGXrjWJmzyFMsF8FGjUUCyGEmBi2MKXWox+FCbbN91rPC0qaHkOI+pebd7wS7r6e7jnvNOALwGJ3\nP3vQ/rImWXdfCZxIyOP6aTP7A+DfCbfjrxgkYL8QQojJSTufbJ1HRdp+PWUpVZ/Xep5hZkvMbLOZ\nPWRmHzGzmT3qn2xm681sU6t+jj0W2Dbn/S/gd4GrgE+5+xtz28dknxF3X25mLwC+BbwCuA54mbv3\n98XfWagb/Sl3CU+V5Tip5TxlbWou4dllRieE9W67dydD34OO4W/PyCgaL+HZr7CeJt63T1Sed3/0\n5mJbK3TbW2NDaiJyUdEmuyFaeZRaglO0yabqbU5sh7xTXXZ5xBT7TpEypc+MvrGbCxGfpsXnZ3ai\nHCeNL76O6u27eyfv9YbZcQp6WB8F++8uz+y5vfh6w5xOvbVR0nfWliyvr+vPkLs0p64dNpehJW0f\n7RIeM1sALAJudPdbS6ru23q+AvgoIUjSwlbb/YBYQr6OcHe8jLBq7Q3Al83sTHe/rMcY3pk45n8R\n/gA8FNVpJtWdmf1ZYte1wAuArwNnmVn7qENNdSeEEGKHYL6ZxZPl4lbK0e0wsznANYS/DK/q029b\neb3M3dsT3i1mNgW4yMwOc/efARTvOs3sy8AS4L3AdpMs2ztIFXlXVG4s1d2n+rT9p8JBNckKIcQO\nTgN3sivdfWG/Si2J9zrgIOBod1/Rp0lbHvtGYfvXgYsIXso95WZ332pmVwJ/b2b7uPsDhSplsWBq\nUTbJHjisg04Kml4vUXcMOdJxcV+TOWhjisdMJAKIc8PuRndu2DhQfLyEJ5aO50eSctjXqTd/U7Qv\nloS7m3S/jhXrWAbdlCgDGyIJMFf6zamXHYmphFRkp5mJ7cXXqXKXFF5YfjItfp1KJFCQ3Ls+gyiH\n7Yzoc9vjwO4PLjYnzGfltvJq5vYsb7dvRrQvjjY2NTNxQJWIT01+34a5fKZG36OK+GRm0wi2zoXA\n8e5+W58mEFKlDoVhmj2T00mdRABCCCF2PHwE+WRbK1QuB44DThpgmc2/E/4mn0C4A25zYuv5e9u1\n6BxzKvBS4J4ed7GY2XR337R9y3Jy2o13/CwhhBBdjMDx6WPAi4H3AOvM7Kho3wp3X5FId/qwmb0P\neEcrO85NhDvhdwKXuvtdAGb2cuAkwhLT+4C9gb8irI55eWJMy1p9f9bd1yTqbMPMnk1wvPohfeyz\nZY5PPwLeDVyds0SnFYbqrYST9P5+9SeMQd08m/Aizukr17s4ZzxNezF2ycWdP20zZ3eEx+3l4kei\n8uq+24uvZz8aBX2J5cli/tNUntRERB4vvP8NkUiUUvPKpN6cwP9l9cpIJQVIeQ2XtUl5S28o/Ad/\nUvw6Rzouvo4/q0ghnrtnt8Y8f/bKqNoe28ormd9pU7g+4mtstymd8uqpHen4idzvWJXA/XXbizYv\naD2f13rEvJvghJRKd7qIkKPjLwmRox4APkD3RLeMkHT9EoIBYx1wK3Ciu9+QGNPZhEn/IjO7HvgP\n4MfArwm/Mk8m2I6fDZzcKn+O/r5LpT/ZnyMkWf+omX0h46BHE3LDfrTfQYUQQkw+RrGEx90PyKiz\nnB7pTls3fJdQEpCiJT8fN+CYvmBmXyKkYX018H62d4Yywp3xFcAn3P1OMiizyV5iZp8mxId8NfBm\nts+xZ4QJ9xrgj9z9WzkHFUIIMfkY56Tt7r6FMIFeYWYzCPLyvoTo3A8DP3f3ZSVd9KRUfHT3R4F/\nAP7BzJ4KHFU8KPDdKgbjnY4ciTlX+h2md3FOX0XiABRTO/rq9EirnUUhQUDk6tsl8yW2hz4iITPl\nHVyUKmNJM/YP3NR7+5aCD+HmjHKRKrEGUsfMpUr+iJzjFOvEcrrleuBmfFYzCp/bzNm9PdNT5eLr\nXaMD7RqZMDZOnc3AVPnu5PaXW2dUnscZjHs+WYBWjP7/aqKvQSI+3cP2Ka+FEEIIkUDexUIIIYDR\nLOEZN3Q2hRBCAONtkx0WmmQnkibtuCMitn/tSidBwMyCTTa218b1UtsBpm+NDHqppSTF87EpsS9R\n3lxi78o1cefYcZuO+BRTJf5btr05sllPi+3XqTJ0v9mELbwYaWt69NnH10Fs259ViLUV22Hj9tOj\nRBVFc3H296rf9uK+KjbUSZIEoB+aZJtFk6wQQghgdGEVx4ncpO1CCCHETo+Z/VkiP20ldCc7WagS\nZapu1JoUZWGEIqZEnU0taIhTotdxvdR2gClbEuGX4q6rRNfZiaPuFKXfsuQBOe27yJVHU59PSZtY\nIk5fH93X1PSuNol47k1GaOu3r27fkxA5PgEhCNM/mNnnCEEn7qjTWfbZNLODgJcATyWsk41xd391\nnYEIIYSYeGST5enA64CzgDeY2X8BHweudPfHS1v2IGuSNbMXAl8gyMsPsZ0bw3aRoIQQQuxgyLsY\nWuES/8bM/g44Hfhz4PPAh8zss4Qk9FkhFSH/TvYC4BbgFe7+64FGLHYMMq+EqVOHlnZRZJJKFtA4\nk1g17JKYdU02yrhPsm3cfTPwr8C/mtnTCXezbwHeYma3AO8vSTiwjVzHp4OAizXBCiGEGBfMbLaZ\nvY6Q//YPgduAdwFzgOvN7F39+sidZH8OUT4qIYQQOx3tJTx1HjsDZnakmf0TcD/wEcIc+Hx3P9Ld\nL3T33ycovG/s11euIPRWgh793+7+y6oDFyVUkebKEotWqZfBli2dL1HshVj8csWSU1wvtR1g69So\nj6lRPtm46zLP59R7m8SyZxVyg1HkvO3swBZl5zn1+ZS0eZxdt5XT10fmNbUl84c9Z2yZnvXZlOWG\nnoTIuxjM7LvAs4B7gYuATyVU3K8REsaXUpa0/duFTXsAPzOzO4FVhX3u7kf3O5gQQojJjWyyrARe\nCHyllb82xQ+AQ/p1VvaX5Qm6vYaXZg1PCCGE2HG5EPhxrwnWzGYDv+Pu/9VaznN3v87KkrYfU2eU\nQgghdiy0hAeA/wD+APhuj31Pb+3PPkm562RfCXzV3R/usW8ecJK7fy73oGIA6tpqm6gX2ZK2bult\nP3uc6V1NNkWvY/tbajvApilRH9OjhSpx6JNiGJQB7YHTSt5zruk6tmNuSGxvwvwW95caT5VkAWVt\nYrM4OWXo/kymZ5SBTdFnH18H65nVc3vZvi6b7DDtnrl23NQYqthnR2zTVexiAKxk33S2T5FRSu7P\n7D8TZvbtJlngwNZ+TbJCCLGDM46OT2b2VOCAaNPvmlnxL/1M4NUEh6hscs9m2cw+mx3Cb04IIUQZ\nYywXv4qw/tVbj3/sUccId7F9l+3ElHkXHwk8M9p0spn9VqHaTOBlQHaIqZ2WnL8rucsDqizHGdES\nnlia29Ql83UnrdgQvV7DblF5Ts/txX1Pmb22s2N2VKkgO3ZJlVMS9eKVQYXfj2mJckyudJyb46HJ\nPLFF4uPmvLfidsu5vor/71OycPS5bYw/Q2BDJP12Xx+9y9AtF8flTRuLF0UGVZbwVGmTIx2XoduX\nUfE54D8JE+nXgTcBPyvU2QQsHTQoU9lHfSphZocws5+XqPcw4RZaCCHEDsy43sm6+zJgGYCZHQ98\n193XNNF32ST7IeCzhJn9l8D/AX5YqLMJeLDPWiIhhBA7COM4yca4+zeb7K9sCc+jwKMAZnYg8ECV\nND+TjvY7LpNhmozyUtfTt2n5qorcHJ2rJzZ2JOINmzqS3Ybps+IWPMLcbeW5UXk1T47Kc7vaxPLg\nut0f2laevXsU/akgO3a9jssJj2QrvM+ZsZQc+QyWSa1VTmEVUtJv2TFT447rxcL+zKLSGr+Oz2Hq\nPBdf7967vHp2vANWMr+zr+v66F2GbnNCLBdvjq7J7O913e9OGak2O4D0O67exWb2C+B0d/9JK+BS\n2Y2ju/uhuX1nXTbu/j+tgRxL8DJeANwH/D93vzn3YEIIIcQk5L+BNVG5MXU2d53sPOBK4FhCJKhH\ngCeHXXYz8BJ3L4ZaFEIIsQMxrrGL3f3MqHxGk33nns2PAL8HnEHIDr/ZzKYBLyG4On8YOLOk/eRj\nMlxHORJxE57GTeqbUTCKTZFMt2Z6tyfozEjqWxklcJrL6qj8SFebuezVKU/v1Js9L3LmK+aCiv/a\n7Z7YnvB4BZi5rlOeFpVjSTUOOAFpj+Km1cCU3FtFyp6ZKhc9hXMk4m7lt/sz2bNT3BiVHy58cPHr\nB6PPPd4eS8rQbWrYsC56F1WCUaS8pTcW6sX7tiS2F4n7qGCSqU1NiXrYNlkzO50wlzwLmA/cA3wJ\neG+Os5GZHQYsItz0zW61/0d3/3BUZxfgb4HXA3sTwgIvcvcvNvtu+pP7M3sycK67/0t7Qyuh7eWt\nu9wLhzE4IYQQo2NE3sXnEMyN5wIrgCOB84Fjzew57v5EqqGZLQRuAm4BXkPwGzoEIoN94ILWcc4D\nvk9YanqlmZ3k7teXDc7MLgb2dPdX9th3KcHZ963932Ygd5LdSnot7FIGDDMlhBBi8jEix6eTC2tN\nbzGzVcClwDGESXQ7WnennwO+6e6nRbtuLtTbkzDBXuTuF7frmNnBhNR1pZMscBrw7sS+Gwnp7bIn\n2dyk7dcAL03sexlwde4BhRBCjC+JYA7faz0vKGl6DHAYcEmfQ5wA7ApcVth+GfCM1mqZMhYQJOhe\n3NNnjNuReyd7HfBBM/sqwQHqQWAvgk32CODNZnZcu7K79/wnssPT5NKesn1l9WJb0AQs4WFLJ8Lm\n41GknTUzum2yu07vrPaaFVk14/JM1ne1ifftRsc0s9u+nfK8xwpGs0cHLEd2V+hewvKkaF9sh+2O\nZdUdsSkuD/PyiG2vKfsqwKyMevGypWnF5TixvXVeVN4rKke2VgD27V2+f/be28r3sl9Xk/j1/ewT\nlTsdFO248RKvrihPWzLPfJVoVik7bFl4r5zjFG2lkyh5wAQ5PrXzkRejLMU8r/U8w8yWEGy6jwD/\nBvytu7e/tkcQYjjcVWh/e+v5cFqBJxKsBg4iSNJFDma7X5Fycs/mVa3n/YAX9NjfNiYbwfV5/BZa\nCSHEDk5DNtn5ZnZr9Hqxuy9OVTazBQRHphvd/dZUPTp/4a4APgq8DVjYarsfQeaF8PdwdY8gSaui\n/WV8EzjPzK6L77rNbD7Bjnxjn/Zd5E6yxw7SqRBCiB2PhibZle6+MKeimc0hmCO3EIL0l9E2b17m\n7u9slW8xsynARWZ2mLuX3Qnn8g6CfH2nmV1LcM5aQAg1vJl0iOGe5Aaj+NaAgxRtciXi3DpVJK8m\nwxJFam0c/Wn92u6IT7Fc/GCkL06JfOSmFPzlpib27Uqnr+kHLu1qM3tj5IgYK8nrEuVNXc2xaN9u\n0b710fZiQP9UgP+YptXAXLk4fv2kqBwLr/NSkjB0S8FxeZ+oXLRoRa/v27fT4fJox7KuLGKwPHod\nS8fxtVKM+BQv2+mO8lSSJGxQU0mZXJwjHVc5ZrGPJpMHVDjOqCI+mdlMginyIOBod1/Rp0k71eo3\nCtu/TnBoOpIgNz8CzDUzK9zNti/O0pgO7v5LM/s9wqqZF7TaPQx8BXhHK85xNgP95LZul48ifGev\nc/dVrZx7j5e5XQshhBBtWnEWriLIvce7+20ZzW7vX2VbvenA0+i2yx7eer6jXwfu/kvgTzKPV0qW\nd7EFPkC4bb4W+AydBLfXMODtsxBCiMlHO+JTnUc/WktxLgeOA17o7ksyh/fvBC3qhML2E1vPbQ/l\nrxEEp1cU6p0B/HSQO1EzO9TM/sDMfjO3TZHcO9lzgTcQDMzfIMR2bHMdIdrTBTkdmdl+wAeB4wmO\nUjcCZ7t7ymW62L5vtI+RUEUGruJdnCv9pryOy9rUjUazsSNibp7aLf2ueTRaG16MENSiaPuJX8eS\nVbx96/TuNoc+oyMfPylHyC2u6I7eWyx473F/VKUgMeeodMUoUanTW6bepSI+xePs9unuloi7HIIj\nWdhiGTj2DAa6nIAPicqHReVndDf5n/2esq38Czpx05d2lbt/o2IpOZaLY4/iNZsK+WRjk0SXdzG9\ny0VyvjvF9ql9ZdGfUh9q7gqAYSYS6POLP6JgFB8DXgy8B1hnZkdF+1a4+woz2x+4mxClaRGAuz9s\nZu8D3mFmjxHW0y4krFu91N3vatV7yMwuAc41szXADwhLUI8DTskZoJn9aWt8e0fbfgX8nbtfOsib\nzZ1kX0N4s+9rGZlj7iLclvfFzGYRTswm4CyCJ/KFhIXCv+3upa7RA0T7EEIIUYERTLLtFSrnsb0K\n+m5C9CcjrFIpqq2LCIH8/5IQcOIB4ANsf5N3HrAWeDOdsIovcfev9Bucmb2coNZ+izCB/6rVxyuA\nz5jZRne/ol8/bXIn2QVA6pb+cbZPfpXitQQj96Htfx1m9hNCNKnXU7LIODfahxBCiMmLux+QUWc5\nYaItbnfCPFEakMLdtxJu4KqE/P1b4F/dvSg3f9rMLicsHWp8kr0P+C16T2i/Q/nC3phTgCXtCRZC\nRnoz+w7BPbrsxB1DEK1en3msZqiSQ7JK31VkpZzAFIPUyyElmW3sTkwaO/pujYIFPD6n4xX6+PTu\nNusj39i4vDYSRdcUBNLVUzoeqIcf2fFnWLB75ECYmws1DkwRnZuZD9HFzCi4RezN+1hULsrFxdc5\n5HgUF52D50U3IXvEenEsC6ck4eLrSBbe+Fud8tLZ3Y1SEvFdkcBVDEbxQDSgOPfw2nWdz7fosf5E\nyqO4irzaZDCKXIm5ym/JiANTjEgunuwcSphoe/F54MuDdJYbVvFK4J1m9txom7eMwf+XEHEjhyOA\nn/bYfjsdz68UXdE+zGyzmT1kZh9puYILIYSogRP8Ieo8dgLWkg6duG9rfza5k+z5wM+Bb9NJFHAl\ncFvr9UWZ/cyDQn6zwCqI8lj1Jo728XWC49T7CbbZf0k1MrPXmdmtZnYrW3uFzBRCCBEYvnfxDsAN\nwHvN7A/ija21sxcQvJyzyQ1GscHMjiGsGzqB4Oz0cOuAl7v7MP3h2lSK9tEK57UYwGYubCzbfTZN\nxxQeVptcknJxwXyypaOZbY7yfT4aSX6bIukYYM2cjlS4ekrnP9dDkZ/s/QV32FiGjMu/GQWtOPTA\nX2wrP+Wwwp/QWFe5s3d52i+7m/xG5Hm8IJKSH4oU6ofpJpaLU1+W4seRCiwxL5K4ZxX14vj0PDUq\nxwEkYk/hgly8OZKF7959/23lWAa+u+DneBcHbyvHQSZK4xBv6i0Lb42ulS55GLpjFFeRUat8d3Ly\nyRbH0qR38ai8jltILgZChp1vA/9pZv9DcK7am7Bs9ZekpeSeZP/MtgzJn289qvIIve9YU3e4MbnR\nPoQQQohKuPv9ZnYkQSV9PmF++hHwYeAz7j6QXJw1ybaiOi0kBFdzwsz+fXffWNpwe24n2GWLHE7/\nKBy50T6EEEJURHey0JpIP9R61KLUJmtm083swwSb6bcI9tAvEG6lHzazi81s17I+ClwLHGVmB0XH\nOAB4bmtfGbnRPoQQQlSgnbR9zB2fGqXfnexXCFEyriFkk7+HsHZpP+Ak4K8Jd6F/nHm8TxIiR11j\nZm8n3BVfANwLfKJdqU60j1KMzjvOjQyTs71qvVSbuvbVKgkCciM+pWxEpRHxo8UoUzvljWu7B7px\nasdkvmpGx4b30NxOiKJlsw/oahPbCu+IHNT3495t5YOj8KUHHLK8q338+oBoJdq8OyORprhALYpN\nZpFNdq+4XDTKxmFW4nMVR6Aq/j7FIVZS0f7LcrtGdlg/qFO+Z14nQtPyQrT/2N4aB/VPRWiC7nyw\nD2+dv628fm3HqryxsByHyPaazAdbDPyfut6GZZ8tvq6ST7ZuDtrcJUAN2WvbYRXHDTO7kzAf5eDu\nfmj/aoHk2TSzFxNCF57u7r3WBX3KzF4EXGFm/8fdv5QxsnWt5O4fJNh2jZC77+yCzl032ocQQogK\njKlc/N/kT7IDUfaX5eXAFxITLADu/kUzu5IQbqrvJNtqcw/woj51llMj2ocQQgiRi7ufMay+yybZ\n3wXentHHV6gWumpiqauIlMk9qe1NLMfJid5UZWy5VJGlctt0yYOdN7p2RVSe+hRifjWjo4P+eO6R\n28rz9u7otXtNeXBbeR+i9TfAvjwQlTv79j2kdxlgTzq68HxWbivPZfW28m5b13S1mbm2k7xgWiHh\nQBsvfB6boghUa2Z3tOM1kY78MPPjJl1LZeIlNA92LYPqyLtlkZjiNg8+2tGlN64upCVYG5kDRrGY\nr4nj1DXP5C7hyVnO0+t1r+0jiP6kJTzNU+b49BS6rE9J7mF7y5AQQogdDMfY+sSUWo+dATP7bTP7\ngpn9ysweN7NntrZfaGb/a5C+yibZWQRv3n48zvauNkIIIXY0HLZsmVLrsaNjZs8h2Gh/h2AGjd/U\nLsCfD9JfP8FwQbzcJsFvDHLAHY6J9hRuOnpTk46DZVJUal9um40Z24v7oihTq6Z2Qo+umtEp/2zO\nM7vbx2rr/I6ku3skN+8xvdtVeG4UN6VLIqYjEe82pXu9+vTdO/9Xp0QuxXG5KNNtIkqmEGUviJMk\nrI6C6xdfx+WVj3Zk5I2/isJEddTudqMO8Vso+9xSJoyyHMdl+/ptLzIs6RjyvIOb+F6n2pd5F+d4\n/Ze1ESn+nuCQewrbT6q3sn0y+FL6XcZXZfRhDMkrSwghxOhwt66sWWPKs4AXufsTZlZ0wF0J7NWj\nTZKys/mqQUcmhBBixyVMsju+5FuTTXSHDo/ZG3g0sa8nyUnW3S8dpKOxoIp3YG5/VY5T17u4CnW9\ni8uCCGzMKBejhsavU1JnSnouvo4CZTw6Y+9Oec7edDFnwDLkfVZlY6tyPlLlsmCoKbm37L1VOR8p\nSbPMa7fpFQE5deoG+2+yTWnAF/Lo18bRJAv/CbzJzK6OtrXV2j+jd171JGOvCwghhAi4G1s2j/0k\n+07CRPtDQkpXB84ws/cDRwHPHqSz3HyyQgghxE6Pu/8QOIagjZ1P8Ds6m6CzHNsrpWoZupMVQgjR\nwnhiq6YFd/8ecLSZzSKsQXjE3df0adYTnc1eNGn7GZWrf902ZeQEKa+7NAfSdsPViXLZvpUZdYqv\nc+2WOyIzEuW5hXpzM8rdQabS9XKvidQq+4n4dWo6KtvGRL0mjpM6v7m22na92H/W6U7eMCaY2WeA\nz7r7t+Pt7r6evKBMSTTJCiGECLiN5SQLvBQ4y8zuAT4HfD4rs1sGsskKIYQIOCGGeJ3HjslewGuA\n5Ue/oCAAACAASURBVISY/UvN7Dtm9loz271Ox7qThWbOwrDc9oe5VKgKVSI55UZvypGIi3Lvyhrl\n4niIIzutSpQBHkuUN6Q6BjaX7GtT/KCmJfbFS/jiRLPF11Eg/41REtqNUcSn1fExSMvCVeTzKqaJ\nmLJgrVX6q7L8rO53NLUkqXg+m4zepqhOA9NKtfrPwD+b2X7AmcAZhDznHzaza4FLgRvc/YlB+tad\nrBBCiA5baj52cNz9Xnd/r7sfTliy8xngOELGufvM7OJB+tMkK4QQIuCM/SQb4+7fdfc3AAuADxIy\nzv31IH2Ml1xsDP6OmwzI30Sw/7rjyaFKZJncnJqpyEVl+1KRiyBPVo4l4i3FMNsPReVcuTiWhWO5\neEuiDnTLxbnkyMXF48Tj2SNRLx5LJB0DrJ7V+5ApCTS3Xu61XxYRLIcdxSO5eA6r5KDNOT+DtvGM\nOjUxs9MJUuyzCEaJewjZbt7bb5mMmaXi5P+uu/8oqrcc2L9HvdPc/eoe21PHOxh4ZWu8BxC+YF/I\nbQ/jNskKIYSYaM4B7gPOBVYARxKCPhxrZs/JsHl+lmArjflFj3o3tPqNWdpvcGb2ZOBlhMn12YS/\nHt8A/g642t0H8lDQJCuEECLgVBNdBuNkd/919PoWM1tFcCw6BripT/v73H1JxnFWZtbDzKYBJxEm\n1hcAuwJ3AG8DLnP3B3L66YUm2arkSkRV2tSVgZsMMlGsl9qeKxfnehfnBsTPkZW7xlZUox5L7IvL\nZZJsyqO4+EtVN7NCFVfSxxLbY7l5Wnrf2mg5RixvFj+DVCKBXNNA7vWaE7RimFLnMM07w1ppMOj5\ncIhSHA+FwgTb5nut5wU99o2CB4HdCbahxcCl7v79JjqW45MQQogOE+P4dHTrOScu8F+Y2SYzW29m\nN5nZ8xP1Tm7V2WRmS8zshSV9fgt4EbCvu7+pqQkWdCcrhBCizQgcn4qY2QJgEXCju9/ap/plhKU0\n9xMcm/4GuMnMjnf3W6J61xHujpcRAk28AfiymZ3p7pcVO3X302q/kQSaZIUQQjTJfDOLJ8vF7r64\nV0UzmwNcQ5jaX9WvY3c/M3r5H2Z2DfBT4ALg+VG9NxaO82VgCfBewkQ9MjTJ9mJYtpem7TB1o97U\ntWXl2mRzbK259aq06aJoK908YBnSb7wsqlOT3iRFG3FqX2xvTS37WZ9uvzFazpP7GcT22lRSgmKb\nnOU8vV5PVibit6DBpO0NnOeV7r6wXyUzm0m44zwIONrdVwx6IHdfY2ZfJSRTL6u31cyuBP7ezPap\n48g0KJpkhRBCBEYkF7e8ea8CFgLHu/ttwz/qxKBJVgghRGA0wSh2AS4nhCo8KXeZTaKvJxGW3ny3\nT72phEw794zyLhY0yQaaWCYzLOm36UDoKbmoyrKfKkt4hikXJ8fjqR0ljXIC+u9IVPngEpuLn0Hd\nzzr3+xK3KUseEDOMZS79jlHltyCVgzaVYKBffyl6feeLiXOGf7l/DHgx8B5gnZkdFe1b4e4rzGx/\n4G5gkbsvAjCzc4CDgZsJS272JwS22Bt4RbsDM3s5YeK9nhD0Ym/gr4BnAi8f7lvbHk2yQgghRskL\nWs/ntR4x7yZEaTJgCt3LTJcCpwGnE9a0PgZ8B3i1u8d3ssuAfYBLCHFD1wG3Aie6+w1NvpEcNMkK\nIYQIjCDik7sfkFFnOYV7bHe/juAo1a/tEoIUPSnQJDsIE+FdPMgY6pArEeWqjimpsErEp7qRpSrp\nX6ng/GPKMD/rupJqk4kymvLS7demLEFASiKukmRhUEYQ8Wnc0K+HEEKIwAQEo9jZUVhFIYQQYkiM\n751sFem3bvsm5OK61P2X2mSw/1zZsaxN1vup+2HntmnamJWSrKclymVtal5swww8UmVoud7FdamS\n33ZUwShyjj9obmjdyTbO+E6yQgghutEk2ziaZIUQQnTQJNsommSFEEIEdCfbOOM1yRr13nEVm8qM\nxPYqx6xK6kuT6/afs0ymaTtd3chSSYonNGXfjAPq5yYIaPrXKWWHzb0Qy2y3NSizpQ8z4lOK3KUx\nTf7aNWmfza1XtqQpdRkW22gCHTnjNckKIYRIozvZxtEkK4QQIjCCiE/jhiZZGG6w/4lYwlPmtp+q\nV0Uurrscp4llP7WJT04q52qUVzV7ELk6XfyLVibppsY5s1Avfp2SmFN1ivUimjYNxHln6177qUhQ\nuf2V1Wnyo879LciRgYuvq0R/6pUgQBGfGkeTrBBCiA6SixtFEZ+EEEKIIaE72X7kRnyqEiVqMkV5\nypWL6wb7X5vYXtamSi7TLmI9LDdCUpl3cVGi7UWuhpjTV7G/Mrl4Vka9ml7HTeQOjt/OWvKoIv2O\nKjJU2Rja5OZ5HlWUKEV8GgmaZIUQQgQ0yTaOJlkhhBABeRc3zvhOsrkx41PbhxnguwrDCvxftq+u\np3Cu9Js7tiyaCEaR03dxYJtL9uX0l5K1i57POV7EDV+IqWsi1+u3ye9LUR4eVt7ZmLqBKcr2Nd1G\nd6kjZ3wnWSGEEN1oCU/jaJIVQgjRQXe7jaJJVgghRECOT42jSRby7bNNR2/KXerTJKkv0DCD/ddd\n9tNoxKeyCEcpO2zZMpu4vzhKVNGOW/dDzbEdF4+TKpfZZxNLeura7Jv2YWgiecCg7ZukyvnItTcP\nclzYPuKTHJ8aRcEohBBCiCGhO1khhBABOT41zsgnWTPbD/ggcDxBqLgRONvd78lo+1TgAuBY4CnA\nvcAXgPe5+7r+B2fwqDGjir5S5ZMY1rKdUQX7r7KEp0jtdK6xPBp3kBuJKRWlveklPDG5cvEkivg0\nquQYRVLLiCYiElQZORJx8byn3ludJTyyyTbOSCdZM5sF3ARsAs4ifKQXAjeb2W+XTZRmNpswIU8D\n3gHcA/we8G7gEOClwx29EEKMAZpkG2XUd7KvBQ4CDnX3uwDM7CfAncDrgUtK2j6XMJme6O43tLbd\nbGbzgHPMbJa7rx/e0IUQYidHjk+NM+pJ9hRgSXuCBXD3ZWb2HeBUyifZXVvPqwvbVxMcuIxBmIjo\nK2VjGCaDehSXSb912+TmoB1ZxKccD9wycuXiJr2Ly8aZ2pdqX3Y+MqkSqatKJKa63vg50nGZvFr3\nI8wN1l/lNyc3B628cEbOqL2LjwB+2mP77cDhfdreSLjjfb+ZHW5mc8zsOODNwMezbLJCCCHStB2f\n6jz6YGanm9nVZnavmW0ws6Vm9j4z2y2jrSceRxbq7WJm55rZcjPbaGY/NrMX5Z+I5hj1/5p5wCM9\ntq8CnlzW0N03mtnzgC8SJuU2nwLe0NgIhRBiXBmN49M5wH3AucAK4EjgfOBYM3uOuz/Rp/1ngU8U\ntv2i8PqC1nHOA74PvAy40sxOcvfra41+QHYY8cDMZgBXAHsBZxIcn54NvJNwWfxFot3rgNcBMP2p\nvTuvspA7V3Wscpwq5H4xcuTeKsEocj2Sq+SGre1BXEbKuziXiZaLi+R4EVfwKM4l15wwNbFvVMFb\nmvQ6biIIf9mlkzpOXbPWxOWTPdndfx29vsXMVgGXAscQnGPLuM/dl6R2mtmehAn2Ine/uLX5ZjM7\nGLgI2Kkn2UfofceausONeTXhAzgksul+28weBRab2cfd/cfFRu6+GFgMYLst9KoDF0IIUZ/CBNvm\ne63nBQ0c4gSCD89lhe2XAZ8xswPdfVkDx8li1DbZ2wl22SKHA3f0afsMYHXsNNXiu63nw2qOTQgh\nxpu2d3GdRzWObj3/LKPuX5jZJjNbb2Y3mdnzC/uPICwTLc4VbTNjP/+fRhn1ney1wMVmdpC7/xLA\nzA4gLM95W5+2vwLmmtnBhYn291vP92WNYGrhuaxObrmsv6bPcI6Ukxu8ocmYwmXSb67EPFRZuE3R\nCb1KJIRYbo1/Vcq0uCr5aVPHjCmLPZwjMQ8xGEUV7+Iy6TinTZHUR5K6vop9pbx2c+TdquT+5lQJ\nRtHrOMWvRP2IT/PN7Nbo9eKWotgTM1sALAJudPdbU/VaXAZ8Bbgf2B/4G+AmMzve3W9p1ZlHuCEr\nKperov0jY9ST7CcJTkrXmNnbCf+bLiBEbtpmyDaz/YG7gUXuvqi1+bPAW4Drzew9BJvsQkJgiu8D\n3xnRexBCiJ2TZmyyK919YU5FM5sDXNM66qv61Xf3M6OX/2Fm1xBWrFwAFO9oJwUjlYtby2yOI3iC\nfR64HFgGHOfua6OqBkyJx+fuy4GjgB8RokRdTwhusRg4PsMjTQghRBntSbbOIxMzmwlcRwhQdIK7\nrxh4uO5rgK8Sov+1eYSgehbv0dt3sKsYISP3Lm7FKC5dr9SaULcLLuHudwAvGc7IhBBCjAIzmwZc\nRVAjj3f32xrs/nZgOvA0uu2ybVtsP/+fRtlhlvAMlVFFfGoguE6SKkH0q0RVajJKVBO2vUbtYXWX\ntqRstUWqfPCpsZWNORX9qeyiHCxw2nZU+dzq2lrLtpf10aZs2U6OHXaYv6JVIj7lJmZILeEZclhF\nM9uFoGIeB5xUthwno68nASfRcYAF+BrhXbyCENu+zRnAT0fpWQyaZIUQQrQZTaq7jwEvBt4DrDOz\no6J9K9x9RS+/HDM7BzgYuBl4kOD4dA6wN2FCDW/B/SEzuwQ418zWAD8gJJA5jhDad6RokhVCCNFh\n+MEoXtB6Pq/1iHk3IfrTdn45wFLgNOB0YHfgMYLD66vd/bvd3XAesJYQdnfvVtuXuPtXGnsXmYzX\nJBvnkx1VgoBhkhP4v/g6Z9lOlSU8TUu/dZfzZC+xiOXRujJuajlPcRB1ZemycQ5p2U4ZOdJp7med\nKxGn6tX97hWl45zoS00E98ppk/ubk5KOU21qWggGxd0PyKiznMLI3P06gqNUzjG2EhxkLxx8hM0y\nXpOsEEKINEra3jiaZIUQQgSUT7ZxxneSzY34lMo1OVGewjlt6kq3o5J+q0R8Guq/7FzpuErEp2EF\n5c8dZ6rNELXCMhm5bpSomCajP5UdM+WF3ESCgByqmLWKY+7leRxfAqNxfBorxneSFUII0Y3k4sYZ\ndYIAIYQQYmwYrzvZ2Ls4pooMU8aoPI3ryrB1c8NW8UgeWWCJuuRKxykNskwebjJBQBmpcU+ARFzm\ngdvkdVwkV0pOkUoQUMaoflWreBf3ej/Fy2FSfQ93fMZrkhVCCJFGjk+No0lWCCFEQI5PjaNJVggh\nRECOT40zvpPsZA72X0bdpQ9Vgv03Gbi/SiSnCf/S59owqxjth8mIQ/kUyV3Ck9rexLKwnGU/VRIZ\nlCUVSNF0JKjU+R30t22CL5OdnfGdZIUQQnSjO9nG0SQrhBAiIMenxhmvSXaYCQKapMncsFXaTFSw\n/7pJAWJyA9U3yk6su9W93svOe13ptoxU3xOxnKeKfJ7bd52IT0Xk+NQoCkYhhBBCDInxupMVQghR\njk/0AHYuxmuSrSIX52yvSo5EVMUDN1dmy43eVEUuTrXP2V6VnByyowrmvjMzzO9B2WdYxSN5UFk4\nN5FBFU/j7BzHDSLv4glHcrEQQggxJDTJCiGEEENivOTiFLnexblt6lLXu7iu3Fs3mMUg9UZBbqD6\nmHGVkZsOmDAodc0OdT3jc+XlKp7GZdS9DhWMYtKiSVYIIUQLLZRtGk2yQgghWijkU9NokhVCCNFC\nd7JNM16TbJMRn5omd5lLyn7UpE216ehNwyQnok6VSDsTFdFrWNR9P5PhlyLn820iSlRdcpb0TNRS\nsr6fo+5km0bexUIIIcSQmAz/T4UQQkwKJBc3zfhNsjlycZNUUV6GGax/WMtx6i47GiZVlvDkUjew\n+0TT9JKdicirXLZ90GU/TVz7w8o72/RSoZ5LeDTJNo3kYiGEEBFbaj7KMbPTzexqM7vXzDaY2VIz\ne5+Z7TbIKM3sbWbmZvafPfYtb+0rPl44yDGaYDL9pxZCCLHzcw5wH3AusAI4EjgfONbMnuPuT/Tr\nwMwOAt4OPFRS7YZWvzFLK4y3FuM7yeZKXhMhf1WRr1IycJU2TQT7H5YUXOaVWSVBQKqvquPZmch9\nb6OK+JRzzLpycRMeyLFEXFc6LqOu6WXi5OKT3f3X0etbzGwVcClwDHBTRh//BFwOHEr6aljp7kvq\nDLQJduafCCGEEAMx/CU8hQm2zfdazwv6tTezPwGeCbwc+FKDQxsKsskKIYRo0b6TrfOoxNGt55+V\nVTKzJwMfBN7q7qv69Hmyma03s01mtmQi7LEwbneycTCKmCbOQlnAg151yvbVDSzRdIKAifYOrkKV\nz2C8vg3VmQgZuUiuOSBHLi6TdOsmFUhJx8NMMFKL0QejMLMFwCLgRne/tU/1DwC/AD7bp951hLvj\nZcBewBuAL5vZme5+Wb0RD8ak+WiFEELsFMw3s3iyXOzui3tVNLM5wDWEmf1VZZ2a2fOBVwLPdHcv\nq+vubyy0/TKwBHgvoElWCCHERNCI49NKd1/Yr5KZzSTccR4EHO3uK/o0+QTwaWCFmc1tbZsKTGm9\n3uDum3o1dPetZnYl8Pdmto+7P5D7ZuqiSVYIIUSL0cjFZjYNuApYCBzv7rdlNDus9fjzHvseAf4a\n+FBjg2yI8Z1km470k7N8pKyPYUagadImOxHkRmzKGefOcMU3+Xk0fT5S/U10JKjiGP7/9s492o6q\nvuOfbx4QHraB0hAgkBvlYQIK1agIbQgggsKiq2vxUAgSbXxQH1gFQYMhiVEQ6ouihRjKK60BFRUf\naVKqNxFbkKwWxIAhaJSCoUkI0QghQfj1j70POZmcc+65Z2bOzJz7+6w169zZsx+/+d2Z+c3es/fv\n18lStGbfXlt9x81zCU+zNlNdH/kv4ZE0jLD85gTgtEEsszm+QdoXgeHAB4FHW7Q5AjgbeKybvVjo\njUeO4ziOkwld6cl+GTgT+DTwjKSj6449bmaPSxoP/BKYa2ZzAcysP1mRpE3AiPpjkt4OnAb8gOD0\nYizwfrYv++kqbmQdx3GcbvKW+DszbvXMIXhpEqGH2sky0zXAfsDngb2BZ4AVwClmtqSD+lLhRnYg\n0i4VaLfuTpbj5BVPthVlXs6Tp2P1MlH0XZvnEp6ig3W06+w/6yU8zepuVx+ZffrJf7jYzPrayPNr\nEr6omuSb2iDtHsJQdCko+nZ1HMdxSkWvvZUWixtZx3EcJ+Kh7rJmaBnZZh6fOqHVLNdW+ZodKyJA\nQLt50g6FN6MTZ/+dOPgfWld5sRTtDSrrAAPNjrWaNdxs+Li+zVZlOtFNJ59xasdaunVw0uKPH8dx\nHCfiPdmscSPrOI7jRLrvu7jXcSML+Tpf6GT4Kc+Zwt2aEdyJc460dbVzNfvzo32yfDq0G785r/aT\ndDKs3K4Di3aOJYeUm8WwznLouJU8L+E92axxI+s4juNEvCebNR5P1nEcx3FywnuyjuM4TsSHi7Om\n60ZW0jjgEkL0hSOB3YAJ0cPHQGVHAZ8CpgGjgfuBS8xseVuN14+EtFoWknb5SSfOx7sVIKBo0jr7\nT+tpy18rs6dsnp2ypBMvUfU0W5rTbrCAtMvX2pnTscMSHh8uzpoihosPBs4ihCb68SDL3gC8G5hF\ncAC9Flgi6ahMJXQcxxmS1HqyaTanniLeJ5eb2b4AkmYAb26nkKQjgXOAd5nZjTFtGbASmAucno+4\njuM4QwXvyWZN142smb3YYdHTCa9Jt9XV9UdJi4BLJe1qZlvbri3tUpYsyqQNEJDXsHQZ7rEsl+20\nogzn2g5VGV6tJ+th5HbqK0JPnVxDnSzhSRvEo9kx9/iUK1W6dQ8H1pjZs4n0lcAuhGHolV2XynEc\np2fwiU9ZUyUjuzfhO26SjXXHHcdxnI7x4eKsqZKR7QhJ7wHeE3e3slA/L1KeErAPsKFoIQrGdeA6\nqOF6gMO2/7l2CczeJ2V9Q12fO1AlI/s0ML5Beq0Hu7HBMcxsPjAfQNIKM5ucj3jVwHXgOgDXQQ3X\nQ9BB7W8zO6VIWXqRKnl8WglMkLR7In0SsA14tPsiOY7jOE5zqmRkvwuMBM6sJUgaAZwNLB3UzGLH\ncRzH6QKFDBdLOiP++dr4+xZJ64H1ZrZM0njgl8BcM5sLYGb/I+k24IuSRgJrgAuACcC5bTY9P7OT\nqC6uA9cBuA5quB5cB7kis+4vkpLUrNFlZjZVUh/BiM4xs9l15XYDPk1wSjEaeIDgVrE/T3kdx3Ec\npxMKMbKO4ziOMxSo0jfZhkg6UNI3JP1O0u8l3SHpoDbLjpJ0taS1krZI+i9JU/KWOWs61YGk10m6\nQdJqSc9KekzSv0ia0A25syTNdZCo51JJJunuPOTMm7R6kDRR0tclbYj3xCpJF+Ypc9akfCYcJOnm\neC9skfSIpHmS9shb7iyRNE7SP8Zn2rPxmu5rs2xPPBfLQqWNbJxp/EPglcD5wHnAIcCP2rwpKh9w\nIKUOziZ40roGeCtwKfAaYIWkA3MTOmMyuA5q9bwcuAxYl4eceZNWD5ImA/cCuwIzCNfE54Dhecmc\nNWl0EI/fBUwBPkk4/wXAR4F/zlHsPPBALGXBzCq7ARcCLwAH16VNILgs+cgAZY8kuDd5Z13aCGAV\ncGfR59YlHYxpkDYeeJEw6azw88tbB4l6lgDXA/3A3UWfV5evhWHAQ8C3ij6PAnXw5vhMODmRfmUs\nv3vR5zcIPQyr+3tGPK++Nsr1xHOxTFule7KEoAH3mNlLa2TNbA3wE+Cv2yi7U8ABYBFwsqRdsxc3\nFzrWgZnt1GMzs98A64EDMpYzT9JcBwBIOofQi/94LhJ2hzR6mApMBD6fm3TdIY0Odom/mxLpmwgv\nIcpKyLyxjAOxUL3nYmmoupE9HGjkJnElwUnFQGXXWOuAA1UgjQ52QtJEYAzwcEq5ukkqHUjaC/gC\n8DEza+g5rCKk0cNfxt9Rku6R9LykdZKuibP6q0IaHdwFrAaukjRJ0p6STiD0jq8zs2eyFbWU9Mpz\nsTRU3ci2ChqwV4qyteNVII0OdiA697iO0JO9Ib1oXSOtDq4GHgFuylCmIkijh/3j723AUuAk4CrC\nUOO/ZiVgF+hYB2b2HOFlYxjBqGwG/gP4HvCBbMUsLb3yXCwNVfJd7OTPtcAxwKlm1uhG6zkk/RXw\nDuA1Fj9ADVFqL9wLzWxW/Ltf0nDgSkkTzaxKoxuDRtIowkvGvoQJU48BrydMAPojwfmN4wyKqhvZ\np2n8dtrsbSxZdtABB0pIGh28hKQrCdGKzjezpRnJ1i3S6OB6Qq/9cUmjY9oIYHjc32LVcdmZRg9P\nxd9/T6QvJUz8OYpqfEJIo4O/JXybPqTum+5ySb8D5ku6zsweyEzSctIrz8XSUPXh4pWEbwhJJhFm\nSg5UthcCDqTRAQCSZgKXAB8ys1szlK1bpNHBROB9hIdLbTsWODr+XaXeS9r7oRdIo4NXAZvqJ01F\nfhp/J6aUrQr0ynOxNFTdyN4JHB3XNwIQF1wfG4+1olcCDqTRAZI+BMwDZprZtTnJmDdpdHB8g+0B\nwuSZ44FvZC9ubqTRw2JgK3ByIr0W+uy+bETMnTQ6eBIYLSk5uecN8feJjGQsM73yXCwPRa8hSrMB\nexDerB4kTM8/nfCA/BWwZ12+8YRvKrMS5RcReiszgBMJD9TnCN/nCj+/vHUAvI2wJnYxoedWv00q\n+ty6dR00qK+faq6TTXs/XB7TPwO8ieCcZAtwU9Hn1g0dAH3A7wmT4M4nvGRdHNNWULf2tAobcEbc\n/omw9vWCuH/cANdB5Z+LZdoKFyD1CcBBwDfjjbAZ+DaJRdfx5jFgdiJ9N8K6wCfjRXQvMLXoc+qW\nDgizaa3J1l/0eXXrOmhQVyWNbFo9ENaBfiQaqW3Ab4C5wMiiz6uLOpgE3A78L+EF4xHgH4C9ij6v\nDvTQ8t7u9ediWTYPEOA4juM4OVH1b7KO4ziOU1rcyDqO4zhOTriRdRzHcZyccCPrOI7jODnhRtZx\nHMdxcsKNrOM4juPkhBtZpxRIuk3SRkljE+nDJd0naXWZQq5J6pNkkqbXpU2X9K4GeafHvH1dFLHW\n9jBJ90u6qC5tdpQnN9/lkj4s6UFJ/oxxhjR+Azhl4YOEhfFfSaRfBLwWmGFmW7ouVXPWAm8Evl+X\nNh3YycjGPG+MZbrNNGA/dtZr3lwP/DnBc5LjDFncyDqlwMzWAX8P/I2kMwEkHQrMBq43s2UFircT\nZrbVzO4xs/Vt5F0f8xbh9/Ui4BbbOQh3rsQXolti+44zZHEj65QGM7sF+DfgWkn7EELQrQc+NlDZ\nuiHZKZK+LekPkp6S9OXkMLOk/STdImmDpK2SfiZpWiLPWEk3S/ptzLNW0vckjYnHdxgultQPHAcc\nG9MtpjUcLpY0UtI8Sb+WtC3+zpM0si5PrY33SpobZdgk6buSxrWhkzcQIssMGHRd0ilRZ9fGIeZa\n2++TdIWkJyVtlrRQ0u6SDpa0JJZ5VFKjHusiYJKkYwZq33F6larHk3V6j/cSwm3dC7ycEEB+8yDK\nLyT4nv0K2wNu70EYykXSHsAyQszRTxB81E4DbpW0u5nNj/XcSnCgfnHMsy/BWXoyBFiNv4ttD4/n\nAMF3bjNuBs4iOOO/GzgGmBnP+ZxE3o8D/0kYih4DfC62NbVF/RAi6GwmOMhviqR3AAuAuWY2L6bV\nt91PGPadBFxFCCrxF8BXCX59LwBulLTCzOpD5t0f2z8lyu84Q4+inSf75ltyA64gfJ/95iDKTI9l\nrkukzwReAA6N+x+I+aYm8t0FrAOGx/0/EOLrNmuvL9YzvS6tnwaBBepk64v7R9DYMftlMf3ViTb6\nE/kuiun7D6CTxcBPGqTPjuVHEEYJnid88250fj9MpN8R06fVpe1FiOZyeYO2fkwIkVb4deWbb0Vs\nPlzslApJfwKcR3iQv07SywZZxe2J/UWEzyKvj/tTgCfMrD+RbyFhos6kuH8fcLGkCyW9SnVduwyY\nUtdmUgYIw871/CCx/2D8PWiAdvYnDLc34wvAHOAMM1vQJM/ixP4v4u+SWoKZPU14QTmwQfn1nKVv\nnwAAAuFJREFUUQ7HGZK4kXXKxtWEntGphKHRKwZZ/v+a7B8Qf/em8SzfJ+uOQwhSfSehp/cz4AlJ\nszJaklJrIylHUoYaGxP7tQlUowZoZ1Rd3ka8nRCc/q4WeZ5O7G9rkd5Ini2E0GmOMyRxI+uUBklT\ngXcDl5nZYmAecMEgJ87s22T/ifi7ERjLzoytO46ZrTOz95vZAcArCbF357D9e2saakYzKcfYxPG0\nPEV4YWnGiYTe8GJJe2bUZpK9gQ051e04pceNrFMK4gzgrxKGab8Ukz9LmAS1QNIubVZ1VmL/bYSJ\nOvfG/WXAOEnHJvKdQxjyfChZoZmtMrNPEHpvR7Roeyvt9dqW18lWz7nxt7+NOtrhF4SJVM1YSZg8\ndQj5GdoJwKoc6nWcSuBG1ikLcwmzeWeY2YsAZvY8MAM4jDCBqR3eKulqSSdJmglcTlgnujoevwlY\nDdwhaUZcunIrcBLwSTN7QdKfRi9TH47HT5R0DaFXuLRF2w8BR0g6W9JkSYc1ymRmPwe+BsyWdHmU\ndRZhQtLXzOzBRuU6YDnwCkl/1iyDmT1MMLSvAJZ08A28KZJGA4ey/aXCcYYcvoTHKRxJkwmOKD6T\nNDBm9lNJXwIulXS77bhEpBHTgI8SlpVsI/SOX3KIYGbPSDqOsBTlSuBlhJ7WeWZWm3j0HPDfhKHr\n8YSe8CrgXDP7Tou2P0t4IVgA7EnoNU9tknc68CvCspzLgN/G8nMGOL/B8B3CuZxGWDLUEDNbFXXy\nI2CppJMzav9Uwv/gWxnV5ziVQ2ZWtAyOk5roFOJG4BAze7RgcUqDpJuAcWb2pgLaXgxsMLPzut22\n45QF78k6Tm8zB3hY0mQzW9GtRiUdBZwAHN6tNh2njPg3WcfpYcxsDWFoekyXmx5LcNThowrOkMaH\nix3HcRwnJ7wn6ziO4zg54UbWcRzHcXLCjazjOI7j5IQbWcdxHMfJCTeyjuM4jpMTbmQdx3EcJyf+\nHziTDoWdVDhqAAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#NBVAL_SKIP\n", "\n", "# Show what the update does to the model\n", "from examples.seismic import plot_image, plot_velocity\n", "\n", "model0.m.data[:] = result.x.astype(np.float32).reshape(model0.m.data.shape)\n", "model0.vp = np.sqrt(1. / model0.m.data[40:-40, 40:-40])\n", "plot_velocity(model0)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAakAAAFuCAYAAAA7wedXAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztvXuwLUd1p/ktwb16IISuJAyNhCwrwJ6QhkcbTdttMSMk\n2iBAD9wIN3Rj2lKAANPQAwG4GRw92AH0GDTG2JjHbeO2Y8TDAQYjdQzmKYSj+2IjPBgk2zwlC4Gh\nJesiJIG4AuX8Ufvo5sm7V551smrvXbvO74s4cXJXZVZlZVXt3PnLlWtZSgkhhBBijBy26goIIYQQ\nHuqkhBBCjBZ1UkIIIUaLOikhhBCjRZ2UEEKI0aJOSgghxGhZSSdlZg81s/eZ2W1m9l0ze7+ZnbyK\nugghhBgvtux1UmZ2FPDXwA+AXwMS8BrgKOCRKaU7l1ohIYQQo+W+Kzjnc4FTgZ9KKX0FwMw+D3wZ\neB7wWyuokxBCiBGyipHUx4EjUkpnFtuvBkgpnbXUCgkhhBgtqxhJnQ58cM7264CnRw5wmFnaqLhl\n221e5sC+7VLr1r190TItx/YYoj0i+YZs2yFYdX12qqOxluuOlonkazl/9B3x8o3pO+eHwD0pGcC5\n556bbrnlll7n/OxnP/vhlNK5vQ4yAKvopI4D9s/ZfiuwxytkZpcClwLcB3jgbPuuLM+uosx9nHRO\n1HLkniz9o2LfjwL5yjJ5vrsbyuT7otfmtUetDSL5vPPXaCkTZdUmq/dsnWWtKJ/Dvvlq79K8PLV8\nkXevxHumy2cy/5x/t0TLtLwvLe2xkf52tu2WW27mmmv+onKmrTHbdUKvAwzEKjqpJlJKe4G9ALvN\n0sYDULvhi/qCrX35R1/U2vG2e6xIhxVl6M5jiPZdN9a1/l7nGv1SbXkPomX6vmM5XidTXudhTjpa\nJvIjuaTPda5aQVgUq+ik9jN/xOSNsIQQQmybH666AoOwik7qOrp5qZLTgL9Zcl2EEGKCJNRJtXMF\ncJmZnZpS+hqAmZ0CnAn8hz4Hjkp/OVFpJs839JxDi8yS01f+aKGlrXNaJLFFzmNNmejzEbkn5bPq\nyVM12WpR71JULvS2l9e/y8nnba8dr+VdrpU/4JQ5yHQ6qVXI5/8ZuAH4oJldaGYX0Fn7fR14+wrq\nI4QQYqQsvZOaeZQ4B/gS8P8A7wSuB85JKd2x7PoIIcT02BhJ9fkbByux7ksp3Qg8bRXnFkKI6TMd\nuW9tTNBb8IaJLVpxX2p6eXRtRGSfp2mPgWh9NPc0LC3t6T2T5T3Mn8PI/FTt2NH6eGujanjvfMs6\nqZppeblWc955akTn6/Lzzs+nTkoIIcSomUYnNbYf20IIIcS97JiRVFTiW5TUtLv4fMDZl2+v/YLY\nrpuYGn1/qdTabF1cKe1EWky2l2VO3peadOelo27Vaiboke+ZIeRPz63aQVLjkcfHjumkhBBi56A5\nKSGEEKNFndRkGFp2anEk6W3vZ91Tl2ZyoufpS9+2aUGTrv0dx9bKRPINLTpFogHUnEBHPZp7ERai\nEmFO3+ew9o5uHPtQB7PT6KT0DgshhNg2ZnaSmf2ume0zs++ZWZq5uKuVedss3+XR86iTEkKIybEU\njxMPA36RLnrFn2+V2czOBJ4FfHc7V7Lj5b6SIXvtvlY70XwRS7+yzNiCDi7KKa1oC+zZItcNaelX\nCywaOX80NlRuWVtKdy0OZvs6CvDarbyHW9/TpcxJfSql9CAAM3sO8AQvo5ntovPN+lrgeds5id57\nIYSYHIsfSaWUtvNb5OV0/fZl2ygDaCQlhBBigZjZw4BfA56SUrrbbHsxhNVJCSHE5BhE7jvBzK7J\nPu9NKe1tOM5bgfenlK5qqcTadlIb48xcq/acOw5BSxA3z5NESXTMHA2INi9/rTyB7UPT4qWi9Xji\nILW5He85XqQ5edRDSsTsvGYann/25qFaHMzWTN2j72Xk2W1r996d1C0ppTP6HMDMngX8L8BPtR5j\nbTspIYQQHqtfzGtmRwO/Bfwm8AMzO3a26zBg1+zznSml+Z6dssxCCCEmxSiCHp4APBB4HZ2Z+sbf\nQzlouv6UrQ4yqZFUOSTuK/9FHUR6RE2+vXrWZIGIRFjm8aSV6HFX/YumRdJbdZ3HQPQ59OTgqEeS\nqCQVeXajJuie7F2T7rwyQ5igRxjawexGOjUcZ8F8Czh7zvb3AF+gM0e/dquDTKqTEkIIAcuS+8zs\nolnyMbP/TzKzm4GbU0pXA5+cU+Yu4NsppUP2zUOdlBBCTI6lzUm9t/j8ltn/q4HHDXGCteyk8kgp\n+dC9HK7PGxLDYi3B+h57yBXrNclku5aC5bGHdha73bhTQ8h4U7YIjFhylrLbkFZ83rNWq0PU0i9i\n0ddi3VdzFhuV+CLPZa3da/k8/PZdfCeVUtreoqeuzCnbyb+WnZQQQogaq7fuGwrNKwshhBgtazuS\nuqf4X6Zh2IW+UWnIW2gYpa+D2Zp0F5FWoiHFo4tCW2iR+PpKd+su/UUXiNZkwBansFGJzivjba9d\njyfR1Z59z7ovGk+qtoDXO09OtN29tql9t81nOiOpte2khBBCeOQz9+uNOikhhJgc0xlJaU5KCCHE\naJnUSKrmcaLF5LqFvsdrMdf1PEZEdeyaXu7NQ9V+3XhzE0N6jKgda6d4pmjxJLFIZ7Eetbmq6JxW\nxLNENIChZ45ee/ZbnsO+1Obo5u071OPENEZSk+qkhBBCwJTkPnVSQggxOdRJrRTP40Q5JI+Yg9di\nQw1JVBaoxaA6zMnXQn7smuQScTq6LPPtFokvej/X3QS9ZFnxoDxqUpUnW9fe174m6BFz8vKdGnJa\noHYPIjJ8+V7WZP2O6XRS6yjFCyGE2CGs5UhKCCFEjemMpNa2k9oYCtdiIUU8LAwt87TE4WnJk4ey\n9NqgJmX2LVOTL4aIsbPd/MuyAoycs4UWbw99jx21CIwStdSL5KtZ53lOpWvWfX1DwedE73vLPfXK\nyLpPCCHEhNBISgghxGhRJzUaagtMPQez0cWrHlGpqe+C15pc4Fn31CykImWillhRq8oWSSwiMUad\nfPaN/bMqy8WIo9G+kusQln59naN672jNOm+3ky8aT2pRcaJqRC1o8/Tdzvatjjc11r6TEkIIUaKR\nlBBCiNGiTkoIIcSoUSe1MhIHddjaHI6n70bLLHKls3fs6Pn7BnSMOKit1adWxtu3yPmdvnNXOYsM\nrpgTDRzZ4j1iUZ4lWrwllHm8d9GbK6rtawl6OOSc5dB4c03lOzZvvupQE/RpsJadlBBCiBqS+4QQ\nQowWdVKjIWpyHTF7XRUtsoLniLYmeeTt4V131NS1Jkt6dYg4/K0RlXO8fC0y3tByUEsMqCFpMVeO\n1iU/9t1OusST7lpMyKMeJ6JLRhZlxl9b5uG9Y6Wz6a1N/9VJCSGEGDXL8nm/WOQFXQghxGhZ25HU\nvN8INeeohzn5og5VCWxvzdfXorBFOvNkvFqcqly2icqKLdaKLXF8IhJfTSKMHKtGS/h2L0/t2BFr\nulYiUmDUgrb2jrV4f/CkvOi9bnk+hqQmz3ltFfVSsbFvs3Wf5D4hhBCjRZ2UEEKI0aJOapTUnGx6\nFjTRGFSeXFiyLLkuQu3aotZ9kTaMyqw12Wq77Ra1uqst3OzryLalPjk1S7+aXLbVsUpanuO+Dohr\n19ZinefFkPLK1yz1FinxRdotIt1B3MHsfKbTSclwQgghxGiZ1EhKCCEEaCQlhBBi5Pyw518dMzvJ\nzH7XzPaZ2ffMLJnZKUWef2Fm7zKz683s+2b2VTN7q5n9WPQqBhtJmdlFwLOAxwAnADcC7wdel1K6\nPcu3B3gD8FTgSGAf8JKU0hdazht1sukFPazNjeRlluUoNaflPC0eBWpzVZ7G7pmjw+LapzZXFHEg\nGvVC0NdLRYn37PUNVhcN+Nn32LX5x4jZ+RFFmYiD2KjHib4m6AS2w7DLYmtteMDJF51zPMhSRlIP\nA34R+Czw58AT5uR5HnAs8Frgy8DDgV8Hnmhmj0wp3bHVSYaU+14GfAN4JXAT8Gjg1cDZZvZzKaV7\nzMyAK4FTgBcB+2f5rzKzR6eUbhqwPkIIIRbHp1JKDwIws+cwv5P6lZTSzdnnq83sS8DVdB3cH2x1\nkiE7qfOLynzSzG4F/gh4HPAJ4ALgTOCclNJVAGa2D7geeAXw4gHrI4QQO5TFj6RSSluKAUWfsMFn\nZv9PjJxnsE4qWJkLgG9udFCzcreZ2ZXAhSywk4qa4Q7pELVGZDKwFqepBa8Nak52vafQM5UF32y9\nJk/1dfTryTk1LwaRWEQ12agF7/mqeWWImLDX9kXNySPeLMp77e2LtnvUBN3bF/E+Ue5bpFwfkUmj\nDmZrMuvWUvGoDSfOmv3/20jmRRtOlJU5Hbh2Tr7rgJPN7OgF10cIIXYAG51UL8OJE8zsmuzv0r61\nMrP7A79N1yf8aaTMwkzQzexE4DeAj6WUrpltPg64YU72W2f/9wBzJ9JmDXQpgA1aUyGEmCK9R1K3\npJTOGKImAGZ2X+DddMramSmlUAUX0knNRkQfpGuli4c4ZkppL7AX4D5maWPI3uIstuZxoiV+0qLk\ng1WHrwbfMs3bXn72rABrVoQeNckmYuUVtRhriVtVw3uO8vPXZLTIsVpiPrUcuxbGPG+bqHTn5Ssd\nHW9X4os6E+5L1Oou+p0Tcdg77/PYMbPD6OwT/gXwlJTS56NlB++kzOxIOgu+U4GzCou9/XSjpZLj\nsv1CCCF6Mbo5qbcB/wq4KKX08e0UHLSTMrNdwPuAM4Cfn7P26TrmmymeBtwYsZkXQgixFePppMzs\n/waeA/zblFJoHipnyMW8hwHvBM4BzkspfXpOtiuAi83srJTS1bNyxwDnA+8aqi4btDiL9WSsmpQQ\niZk0xKLQvnGn8vJe3KgW+aKUczwJqCZ5bJeo9VbfkORDW/dtHfa77Tw1CamljGe9WYYx9xbt7nK2\nl/v6ho8f0ioz+kz2dcwblceji3nnf38sp5OaOXGAzokDwJPM7Gbg5pTS1Wb2q8BL6dZDfdnMfjYr\nfnNK6atbnWPIkdTvAU+nW1l8Z1GZm2ay3xV0HiYuN7OXc3AxrwGvH7AuQgixg1naSOq9xee3zP5f\nTbc+9kmzz5fM/nL+CPjlrU4wZCe1UZlXzf5yfh149czrxHnAZXQXcwRdp3V2SunrA9ZFCCHEgkkp\nVY2tU0qP63uOIRfznhLMdyvze9UmIqG1y321Mn0X8C7Sx1/keLmsUVqM5RJfLtvULOByogsNI/7c\nau3eIgV6loM1i0LPsmxZ1n2185T3bqvyQ1uC5s9HTaZtke4ismDNz2J0Ae92KctHnsMW6a+8ty3v\nS4x1swGcj0J1CCHE5BiP4URf1EkJIcTkmE4npXhSQgghRsukR1IRk+1Sd/Z05NocjDe30GIyvkjn\nl0PGL2oxqd3t5IkSnSvyPB9EvRjUTNAjcyBRM/6aRwGvTH4Nfec/a3OJ+bG9GEfgzyNFzcmj7d7i\nZSLC0LM2fb1HRN/ReW2w2YJhOiOpSXdSQgixc1EnJYQQYpRoJLVSjIND3Jrk0eIsNuKZohYLyctX\nDuNbJgOHnEBsid8Ulfsi3grK83uyi3fNNRmuRXbKZauaSX7LPcjbw5OnSk8OOV5bt3gxqUlLnjzl\nxYkqP+dS5BFOuizjyYU1uS+nJU5UdDmKR4uz2FqZyBTD9plOJyXDCSGEEKNlLUdSQgghakxnJLX2\nnVRNUstpcRYbdY66KC8T0WFui+Vgi1VTJMZR+XmXs70mEXrnqVnWRZzFRj0fRK37WvA8OdSIOH4d\nwrIt4uGgZiEZac9avvzYiwz/3vIetEjd0dhfkWM3OTpO8jghhBBirPRdczISNCclhBBitOzIkVSL\n7OSF/S73RePb4OwbIu6Ul39Rg/9ae+QSzp1ZuiX+UU5N/ui7qDQqO0Xx2qPFIbL3fEXrGXWI6ln0\nHV6UiTiLrUmEXlt78c7KfC0METfKyxNZzFs6mO0bZ23uvU4s7oVfMjuykxJCiEmjTkoIIcSomcic\nlDopIYSYGhpJjZOaVt1iaup5nKg5mG3x5OARncca+tgRoqvuPXPuvg5mS7z5v9q8oLcv6si2hehc\nkefFJJ/PGNq5qveMR+eXvHzlPJYXbLJlji1K3+9r73kvn31v7qnFBF10TKqTEkIIMUNynxBCiFEi\nuW/1bAz/+3pOqMkK0SF+ZJV5eZ6o81qvbt6x+hI9Vk2iyOWcu7J0Lu2UZrhRZ5wbROM81eS+iFxW\nk5la2j3q+cTzYtJX4otKVVGvHRFz/5pEWHPmm7Pdtm4ZRESdxfb1HhFdftFb/lQnJYQQYpQkJiP3\nyeOEEEKI0bL2I6maxZjXA0eH3i3WfV75Ia3+ymPnDB1+3pMva9KdJ/X8IEvnMmDtPF6eRTJ0G64a\n79ktY1h597cm93mWerUynsTX19tKVDbvKwW2OJiNerkZFMl9QgghRokMJ4QQQoyaicxJ7fhOKrow\ntzbE9yRHL5ZS+XlZC3OHlLG8mFEQk43KukQWPkatqiLSYY1amUgbDvHdELEMq7XHAWdf1LovGhsq\nIvHVQsFHrRUji6D7LhYv8e5jrQ0POOnos9tiybmxzyp51pkd30kJIcTkkNwnhBBi1EjuE0IIMUo0\nkhontSB/EXP0kqjHiYjmX3MgOiSLXPjmzTXVNPaoWXKu33tOOmuB47Y7pwVx7yAeLcHzonNnXgDC\nlvbI0147g28aXnMw6+3zTNOhLTCoR995qL7eI6IBDO92ti+MJXRSZnYS8KvAGcCjgCOBn0gp3VDk\n2wO8AXjqLM8+4CUppS9EzqPFvEIIIVp4GPCLwH7gz+dlMDMDrgTOBV4EPI3ud8tVs05uSyY1khJC\nCDFj8UO2T6WUHgRgZs8BnjAnzwXAmcA5KaWrZnn3AdcDrwBevNVJJt1JtUgBEQko6nHCk71qx47K\nH0NKI9HhdCn1zDsWbG7DXA6qyVPfydK5N4o7s/QuJ0+5zzN/7uuQFXxz4RqePJRvz71xwObr89Lf\ny9J5O9Xy1crcP0t7Juj59nKfZ1pec45ac6Lqsajv3tp3ROQelp+949XO09ck/16WIPellCK34gLg\nmxsd1KzcbWZ2JXAhgU5Kcp8QQkyRH/X8G4bTgWvnbL8OONnMjt7qAJMeSQkhxI5kGC/oJ5jZNdnn\nvSmlvds8xnHADXO23zr7vwe4o3YAdVIFERmu5nEi6hFgGU5Mh5AL8zbIZY5c+qt50/DkvrI98n23\nZ+moXDdk7J+8zjVvCaVk6eFdtyfjgS/L5dtvd9Ll5+8620vp7n6BdBkKvsXjRIvEl5OXGdLCsvYc\nt4SCj3o+ibyLK3J6fEtK6YzVnPog6qSEEGKKjGOd1H660VLJcdn+KuqkhBBiaown6OF1zLf6Ow24\nMaVUlfpgTTspY/7wt++QujxGxMFktHytTHRRqSc59rX0a6Hv81+TefJ9uTzlWf2V+Y7K0rVFpR7R\neGMtP1T7Whvmi3Fzi8CyPTy5MJfucmu+8nPehkc46fKzF0esb5wo2H5bt8iAJS1Oer13OVqfQWW9\ncYykrgAuNrOzUkpXA5jZMcD5wLsiB1jLTkoIIcTqMbOLZsnHzP4/ycxuBm6edUpX0HmYuNzMXk4n\n772Sbqzx+sg51EkJIcTUWJ7vvvcWn98y+3818LiU0j1mdh5w2WzfEXSd1tkppa9HTqBOSgghpsgS\n5qRSSluGsUop3QpcMvvbNpPqpGoOZvtSmyuKeKaomap6gRaHJqJ3l+f3gjPW2sO7npqD2XxuI58b\nyedd8jmpsj3zfV7wvtJ8Ot+Xz9V8NCWWwfl28P2ueYzI03l75HMm5b3N2/C4LF2bX4rMQ9U8TngO\naksiAQxbiDqL7usstq9T2vL6vTm7XvNT8oIuhBBi1Eykk5JbJCGEEKNl0iOpFpNtb4hec/zoOZKt\nrTj36lOrZ4tTWI+oxOC1VU2izM2P8zrnclIpl3veG7zz3K/Il8tbxzrpf5f+tCh1i5P+j1m6XMaR\nf/4hMXKR7KCrsivTq7PteU0BTrg39Qf2rHvTuSPe27J07lUCfNP9mpcMzxNEzZx8SGexUYfQEU8O\nUU8Q0Xe0JUbZItmyfcezTqo3k+6khBBixzIRuU+dlBBCTA2NpNaPqJTQNwZVdMV5JF8tblVfzwU5\n0WNFrRBbrAg9KS+3THugkwZ4cJZ+UJbek/5Z9ulXNhdK3zyY9ry1lp5fPRPD/MaVWowXkKnm/uHI\ng+LoJeknD57SvnRv+uYs+7eK4vnnPF8uEebSIfiyYDTOkkf0u7Lvd2rNE0RfB8TesWvn8a4nKpnW\nyoSYyEhKhhNCCCFGy44ZSQkhxI5B66RimNmfAU8EXptS+rVs+x7gDcBTgSPp3GS8JKX0hUXWpw9R\nGS4acj6SL+pgthaqG2df1CrL21eTKLxFu7nSVdqy5fLdQ7L0Q7P0iYc5O2Bz7M9vZ+lP/+XBdK57\nAfxjls51sFzuq62yjeo5EbnvAUWZYzMXsccflPju89cHNz84a7QHP6Uo//cHk3feejD9jSxLJnYC\nm5snK1JVPyNxlkqLwpbvzogVX9Tq7u5AuvzsOZuNWgSuJB7UROakFib3mdkzgUfN2W7AlcC5wIuA\np9G9xleZ2UmLqo8QQuwYNkZSqw8f35uFdFKzkdIbgZfO2X0BcCbwSymld6eU/my27TDgFYuojxBC\niPVkUXLfbwLXppTebWZlzJALgG+mlK7a2JBSus3MrgQuBF4cOUGkd42MdqMLe2vbdzn7vLDy0Xy1\nMp705i2OLPFkvNLqzvPHli/wrPnEy2W9XNIr1br88/G59HWqk35ZcYAPZenctC2X/nJ5D3xNK5f4\nahqQ94C06J+14E7HZ+m8EXMzxtcV5T9wMHm/rx1M/+RXs3QmCQLcmD08uRSYN2feZBBrtrI5IjJa\njZaFuRG5riYRenLhIhW1XiMIzUn5mNljgWczR+qbcTqbZxA2uA54tpkdHYnWKIQQosJE5qQG7aTM\nbDfwduCylNIXnWzHATfM2b7xI20Ph/qiwcwuBS4F2c0LIUQVjaRcXkFnrffagY9LSmkvsBdgl9ly\nYikIIcS6opHUZszsZOBVwHOAw80sn6o43MyOpZOw99ONlko2nAvsH6pO+YirxZFlNH5SS9wZT7PP\ntfNonKa+5OffXezzplBq0ym5l4jcnPyULH1iOZGVzzc93Nn+2Cz9maJ8PveU21J7ZuaweRLFC+BU\nvugH2D6Rib3S/UNunp5PBHluJkp78n+epfMJv9x+NpurAjg5+/zgmw6mb8jyZJsPqU5+CTUr/pzI\n/FKtTNTsPWJ2XjOVH9LBbNSJs+gYUjk7le6763K6jmbjD7pp7v3AI+jmnk6fU/404EbNRwkhRE8m\nZII+pNz3OeDsOduvouu43gF8BbgCuNjMzkopXQ1gZscA5wOlJaAQQogWRtTR9GGwTiql9B3gk+X2\nbu0uf59S+uTs8xV0HiYuN7OX042wXgkY8Pqh6lPiSX8tRE3Dax4nIjGkak5pPVkx6vg1EuIdNst6\neTwoz7QcNqtLP56l95yQfchlPPAlvjydy1O5vAcx3SlPgy/x1czMI8GMSiINX+qseZx4r55euvz8\n01k612aPKcpk+3ZncuNPZu1+RKGJ5aptfjk1iSaX/2rSWY73znqm4S0eJ2oxzqJOoL1riEpWg0lb\n8oLeTkrpHjM7D7gMeAvd9+A+4OyU0teXXR8hhJgkGknFSCnZnG23ApfM/oQQQoi5yAv6wERkCfAt\nkTxLP9i+I9uo9OcZn4Ev8eXxm0rvEbnEd0zuFcGT8Wr7bnbSpfcIT9bzLPjKz55FX3kTWiSUFvci\nQ/JXWTqX/mrmo45EeXJhEXifH7AlLc9+9HhR674DTvpHzvba8SIxo2CxVntbWihL7hNCCDFqJPcJ\nIYQYJfI4MR2iQ/KaFVLE2Wu0TE3+8CyUWmQFz8jsfkW+3BjMi/lUxlc5Jl/N+9BAujzgbYF0Fm4J\n8GU9z0oOfDOvoYMhRR+KCN6DU+q0nuVg7jXzp4oynoVjxQTuxC9n2e6Zn61lYW3Nyi0i8dXOGXUW\n6+1bxXd/xFn0IZP/E2HHd1JCCDFJNCclhBBilEjuE0IIMVrUSa03Ld4naiafnl7teZ+oHa/mpSJi\nXh41ic2nKTyvEuB7lshN0PeUzmIf5KTzA+TzVuWJcjPniBPY8nM0ql5kHirqQdQzMy/L1PJ5eBOI\n+bWVpuD5nF1+g/N7Vd6DvK2/52wvvcVmnx+UeZ+9a36WQ/bV8kWIOnSO3N6aCXo0AGqON324NIey\nE5H7FJpJCCHEaFEnJYQQU2NJXtDN7Ewz+4iZ/Q8zu93M/srMBvUktPZyX23ovApJNuos1pP4arJi\nNKaNhyc/5mpQaYKeq3DHOelDPMwe7xzAS8NmW/dcaoo4gYWYB9GaObmnG9Vu3JBySvkge+fxrqfU\nqvLPnkl+bb1Bfn/y+1nG5Mo+H5k59j0uC7hT+vXNP+fPXi5B197ryK2qOWf2pL9amch26L8cpO+x\nNrFguc/MHgl8DPg08Fy6N/ci4B1mdnhK6a1DnGftOykhhBAFyzGceAZdX3p+Fgfwo7PO69mAOikh\nhBAOi++kdtON2cul9bcxP/p6E2vZSRmxoXA0THxfItZ5UV+eNTWprz/SFus+L/zQ7qOcTOVnL52f\nFOBvsnQuVUU9QXhWfFHpLhowKNLwfZ3Qgv8gRWTA8rMnC/7Xosxjs/RRTrq81w/I0plEeGwm991a\nFMkPkYf+isajyolK4H0t9SZizT00fwi8APgdM3sdXWf1dODxwC8NdZK17KSEEEJUGMYL+glmdk32\neW9Kae+9p0jpWjN7HPAB4IWzzXcDz08pvaf32WeokxJCiCnSf/h3S0rpDG+nmT0c+BPgOuD5wPeB\nC4G3mdldKaV39q4B6qRWRiT8e1/rvlIS9daE5us7S4Ov/PMmKfAoJ1O5L5f1djnbIebNMyrDRV/O\nSJnar9GW8PHRY+VtFVmxGm2PmkWgd3/ch6DYl6V3Z8c6qpAi8yLe4xG12o0u5vWk8haJby2kv+UY\nTryO7onihjQhAAAgAElEQVQ6P6W08TR93MyOB95kZu9OKfUez2mdlBBCiBYeAXw+66A2+Eu6hQs/\nNsRJNJISQogpsni3SN8CHmlmu4uO6mfoVueVNjNNqJMSQoipsRy5783Ae4ErzewtdHNSFwDPBN44\nZ4TVhDqpCjWPER55vlp5T3OPOJFtxYuR55mjw+b5qk378kJlRb19NRvjlnmoCC3lh7ZFjpSptUff\nukS9aXj3yntYyn2OI9sjijmpI5x0dE4qQm21gZev1uRr6at1wZVOKb3PzJ4M/Crw+3S386t0ln5v\nH+o86qSEEGJqLClUR0rpQ8CHFnkOGU4IIYQYLRpJFbTEeVoUUcmiL9G4V5skGE8OKg/iHSDqUNVj\nZHbAuaGtreKnX8vDUuaJ3Kvy2jwJd/f8LGW2aJylIZt0kY+Od20rYWTvSCvqpIQQYmoM43FiFKiT\nEkKIKaKRlFg0fSWOWvmVSxFjYyK/OqeCJst7siTDiWWgZ0EIIcRo0UhKCCGmyETUAXVSBS1+RpfF\nkBJdNOy257TzPrWw7C0OUWsWZPMYWq9chfnmkNTMMj1qFpbRmFyeY+AD8zeXn6Onyel76/PyefXL\nJuy7dtw751KYkNynTkoIIabGhKz7NCclhBBitGgkJYQQU0Ry3/Rpucfe0DRqDh7VrlvK5HhzAaXb\nYm/qafcBZ0d5EO8ApRThLdXvu4R/SVpBk5eJISdXakTaFvx7dVeWrj0gd83fnm+uZGtSp6KXFpmH\nqr3vazdlqTkpIYQQo2Yic1LqpIQQYmpoJCX64qlDNRkvou7kx62Zk+f77nLSAN/L0ndm6ft5O8pC\n3sHLE3mOSqPSX0QWLLeXMuVW5UuiXwIR8/ra9bQ47PXK1GJDedJdns7vLWy+91n6wA/mZ6kdOlcS\nW5o22hxe00RN5T0n1GIxqJMSQogpMpEeVJ2UEEJMjQmtk1rLTiqXW2vKzLJ+SETUpZp0F7UIjEgW\nNfJn1jPkKtWc27P0d7L0sVnG3bdVCuXpmkR4Wpb+iyztyVO1kPWZ1FRt6Kju49HXOq92E7cb26ls\nj0j49/OKMt/I0vm98u4nQH7vvzM/W+3x8Cz9opZ2UTXYa0JP0qsdb22++ycyktJiXiGEEKNlLUdS\nQgghKkjuGw9jG9F68kM0LPsuZ3utTMtwOLekyiWXUoX7bpbOZZtbs/SD/7Eo9AAn/R1nO2zWgI5w\n0rkEVYtJ7kldpTVfVPcZkuhNjDxINUu9/HPehodn6VK687TdPF3e6+zzgTvmZ8uLl6fJn70WlTVi\nxAgx37nl++Y5Xq75Sc4/R9/L/NjeouEmZXlsX46NrH0nJYQQokDrpIQQQoyaich9MpwQQggxWnbk\nSKrlB0aLSWyLQ4CaVXFEi4+a5HseJ0oT9HweypsqOjw3+Qb2fDv7cD8nfVRxonzf/Z0K1TxW5BXK\n556iExUepUPVmksPj5YHxJtvOsLZns81web2zdv2mCxdThbd7KS/7aSLz99ysuXzl7B53tPzOFFr\n2sg8VG36MU/nTVjz0BI1R48GbvTqNhiS+4QQQowWdVJCCCFGzUTmpHZ8J1UzIfXucdR7RHQFvGd2\nHlWA+saTyqWR0gQ9IjGWE5u7M33nfhHTcNgsY52SpY/L0p7dPPguNGpuDFq0mb52wd7NrnmM8HRW\nT9KDzZJpbu7/P2fpvy/KfDNLfz1L3+CkgW9k7ZsXz9XCUlXMFdxoPKntenVpMUGvlYmYo8/7PBTl\ncec5uU15hgmNpGQ4IYQQYrTs+JGUEEJMEsl98zGzJwP/Afhpumb6EvCKlNInZvv3AG8AngocCewD\nXpJS+sLQdVk1NUnMC59U8zhR82ARyePJFLnkUhjqHSL/zaNUFfLznHTLwfSemk7ivVCnOGWiwX9a\nJA8vzhS0aQ8RTaomf3qy3rFZOpdFAR6YpX86S+cS39eKMl/N0l+Zn+/Gom1yVTC36Ms9TtTiSUWd\nyuZErGFrMpyXLh8p73jRR6ClTJQt22qJct9W3/l9GbTtzOx5wAeBzwK/ADwdeC+z18zMDLgSOBd4\nEfA0umfhKjM7aci6CCHEjuZHPf8CbPWdPwSDjaTM7BTgt4GXp5R+O9v14Sx9AXAmcE5K6apZuX3A\n9cArgBcPVR8hhBCLI/id35sh5b5L6IZ6b6vkuQD45kYHBZBSus3MrgQuZIGd1JDybNRZbDSelGep\nV1vMG42j4+GpY+Xa1RxPJinVMc9h7UMzq78T/7oolHsdvdPZ/tgsncth5eeouaTXcPkFlDehJgV6\nRCS+8nq8xc25rHd8ln5QUf4hWfraLJ1LeqXcl0l8389M9XKFMJf3wLfi89ZgQ0ziq/2Q95rzHmd7\neTzPUq+2MNdbwFtzSuuds/aOes5mt81yvKBHvvN7M6Tc91jg74BnmNlXzeyHZvYVM3thlud0Nr8y\nG1wHnGxmRw9YHyGE2LksXu6LfOf3ZsiR1ENmf28A/g+632xPB95sZvdNKb2J7nfgDXPKbvzG3gPc\nMWc/ZnYpcCnIbl4IIaoMYzhxgpldk33em1Lam32OfOf3ZshO6jA6YeKXU0rvn237xEy3fCXQq8Kz\nxtkLsMssbZEd8Ee70XvnDfejVnc1P3yRfDWLwKi/Po++1kq5ZFPKOV7o8NyH27cKM8JT/vZg+vhc\nU/K0ppcVJ81Dn+erSnOTs1ybgs0maHlFPedy5WfvATtkdXOWzn3seZIe+JZ7D87SucR3YlH+j7N0\nLuvdcDCZrt9cxGvq3Cdf1A9fbZGu56MvKvF572XNus/z0ReNDRWNQRUxLI1eZ2/6y323pJTOqOxf\n6Hd+fpKh2HjlP1ps/wjwIDP7J8B+utFSycZruH/A+gghhFgcke/83gzZSV0XzHP6nO2nATemlOZK\nfUIIIbbBhty32DmpyHd+b4bspD4w+//EYvu5wE0ppX8ArgBONLOzNnaa2THA+bN9QgghhuCenn9b\nE/nO782Qc1L/L3AV8HYzO4FOCX868ATg4lmeK+g8TFxuZi+nk/deCRjw+u2cbCjryqjjgxotJuie\n2XnUbD3ifaLEm4eqXbOnv9fmpPJ8uVlyPiVUmjLn0yYPycbTuan6j2dp+4viAF/K0l5go+icVF7p\n0l2C57w2p+YZODc195zAwmbz8nxOKvckkc9PPb4onzXwd2+bu3nT1B1sbrboFF1kfqdm2u0dK4o3\nV1UzQffmrsp6eWHJap5kvHekZQ64lzn6cjxORL7zezNYJ5VSSmb2VOA/Ab9ON/f0d8C/SSm9a5bn\nHjM7D7gMeAvdK7sPODulVH5vCSGEaGXBnVTkO38IBvXdl1L6LvDC2Z+X51a6RWCXDHluIYQQyyXy\nnd8XeUEPEvXqEJUFImbnUY8TXp4aUfN8b6V9zfGC9wMuV9vKenr+VHNL7FzpeuBNm8s/JCuUq2DH\np5/LPhWF8s93ZFeX29DnrhNgs/10X7kvv9A8rDvA7nzDQbeWB+xv7k17qma5L2/33IS8jPOUS3kt\nP8K9Z2/oH/QR7w/l8+2ZoEfjt3nLL2oeJzwpsqzbwsLHywu6EEKI0TKRoIfqpIQQYmpMKDLvjumk\nWrxMeNSku5bw8X29R/T1OFEjIgvWLCQ9maRUyiIr+nP5pRYtPZcIH2D//d70S9N/LUplwa6O/s78\ndDUaUmnrtkH5WuWfc/eUXnAogBPuTb3dnnFvOq9Zns6t8QC+m6Vz5x61++75vs2lspozDc8CLfp8\ntnyneucsjxXx8BK17svL1Bwy5+epXZu3ryZf7iTXcDumkxJCiB2F5qSEEEKMlYmofdPupCIOHms/\nNloW2fZdmNviyDaKJ6O1SKE1x6DeQt+aUJZ/9mIRRa3PPGnmv9h5m/Llhna58PbJFPJf3JufN7s3\n3dIeNck1b4MjAmnYbGwYLZPfa88B8iqkqprVXcTSD3xZz4sFB75EWHtfWtpjq0HShKakpt1JCSHE\nTmUiat+Omn8TQgixZmgkJYQQE0Ny30gZ+qYM6WWiZrYePeeQK9NbTNC9OZDSnNybR8pNoXMTadg8\n73K7k87nCPL5E9jssCGPK5jPoUTb83/N5opKovOZOS2m2Xm98+vOry1vz9JQPm/P3ONELc7i/Z30\nmL/sol4uIu9bzVlsPveUP++1wKRRDy2R56jl3Z+K3DepTkoIIYRGUkIIIUaMOqkV03IDoubT3rA6\nKt1Fzcn7epmI0GL6W4v9E/UekX/2JKnSb2vuPSGX+DxPElGpKpe3Sokw3+fJazUnv1571trQa7ea\nCbrXbvn2vM1gsweKw53tY2PI96A8lhfPqWXJR9+Yb1PpPJbFWnZSQggh6mhOSgghxCiR3LdDWKTH\niUj5oRexRX5ZlQ92RKqKWvd5UlW5z7NA89KwWQr0wjSVcl9E4is9CvT1DuC1Z+nJIW83rz59Lb7K\neFKe49Wa/FmLrRShxfNJX7wQ79H3Om+D8tmPvP81OThKpK2n0klpMa8QQojRopGUEEJMjAkF5p12\nJxW16MvxrHiii2wjsaXmfR4Kz4lsLV9t+z1OOs9XxtSJOJjNrf5KvMW4nkPY2r78WKVU1RI/KSLh\n1ixG+zoGjiyuBv8e1JzFeveqxdqxrzPkRcag8s5TWwAckf5g87vgla9J6lEi7/lU5L5Jd1JCCLET\n0UhKCCHEqJnKSEqGE0IIIUbLpEZSQ5hyRuahSh06uho9cp4WFhXMEDbPbRxwttc8Tnjly3ksb07I\nmw+pOfb07kFpTu7t6+vkdwjPCfk9iZiD1+bbIu0Jm+ehossNvH3RNhjSy8TQ8lbLfJsXULE2J+7t\ni34vbJTJQ3VqnZQQQohRozkpIYQQo0QjqQnT4j0iIjUN7UQ2/5XU13llLe6NJwXWJCBP1quVyc2h\nI/JprT2jccAickrf8tD2izbihSTqBSUqEXqeGKL3OvIewOb2OFDJ14foe1B7Vu5x9nmSHvixpmrn\n8WJN1Z6brdpqSp2UDCeEEEKMFo2khBBigmhOaiTUhrQRLxPR2FDRMNMtZaJ4VnwtD2PUG8eQ1n2e\nfFJ+XkasranjtVtUto5ItuXniNRVUkqOEVbhYSH6XnvOZ2vWvIuINTUluW/tOykhhBCHMpWRlOak\nhBBC9MbM/szMkpm9Zsjjru1Iqs9QdujYUIuSp2rh3/s6mIzEOAJf6sm3l6HPI4t5a5JH3zaMOmGN\nLKKsOYvt+0u1JrN69Y46Ss5pcZTsLewt7/UuZ1/L/czL1J79CENIXd41RGXrvG1q75h37Jb7u8Gy\n5T4zeybwqEUcWyMpIYSYID/q+RfFzPYAbwReOlTdc9RJCSHExNjwgt7nbxv8JnBtSundg1S+YG3l\nPiGEED7LkPvM7LHAs1mQ1AcT66TK3j9idl7quZEV/dH5lKHnqrbrZaLWHt48Ry2AYYvTUU+LLwPu\nbbd9aqby3jlLc2fPZHpoInNKNc8F3rVF57Q8ok56vXlF8OehWtqzdADs4b0H0TIe5TMYcdzcMqdd\nPofee9niIHpgTjCza7LPe1NKezc+mNlu4O3AZSmlLy6qEpPqpIQQQgxmOHFLSumMyv5XAEcCr+1/\nKh91UkIIMUEWuU7KzE4GXgU8BzjczA7Pdh9uZscCt6eUeveVO6aTikp3Ee8RQzgdjdDX/LkmAXkS\nX026u8tJlxJQ7XgbRON45dTkMU/i82SrkhbZyCtTuzdePWuOfb0ytfbwqEnQ+fV4nhNqz8cPgueJ\nUHP2OmT8tWgdvPatSaaRZ7L83GKOPq89lmCCfiqdan/5nH0vm/39U+BzfU+0YzopIYTYSSzY48Tn\ngLPnbL+KruN6B/CVIU6kTkoIIcS2SCl9B/hkud3MAP4+pXTIvlbWvpNqWaGdUxt673a2R0OKD20x\nFrm+Wp6IBVzUgWitzA+cfbVYRn3x5EtPtooeq2aJFSVi0Re1quxrCZYTlZ087xPlvr7Wffk1LNLS\nzytfI+J9ovzsWfSV99o7dvR5nYcczAohhBg1q+ikUko29DHVSQkhxMTY8DgxBdayk2q5Ad6Qeuh4\nUn0X7fb99RNd4NkSHjxaxpOhcsmjxcFs9Nq8MPct1Cwsc2rOUb3jeTIe+BKbF59riMW8kZDzpYNZ\nTxZcZEywiHPlZdHy/VFKmZE4XJH7m+ZVcAKsZSclhBCijuakhBBCjBIZTgghhBg1mpNaAyJeJoYI\neriKeCfeA1gzS47MKZXmsT9w9kW9VHjzHFEzfo/ofFv03njt1uIZoyQShDFq+l8rE6H2HHtm0vkc\nSvnceYEsy7mrSH2iRM3Tt3ueFk8UtTlLbzlLed+8ZRI1U/et6jqlkZTiSQkhhBgtkx5JCSHETkQm\n6A5mdibwfwKPpnPh/mXgzSmlP8jy7AHeADx1lmcf8JKU0hdazhl1vBiV7qKyYIRomZaHybvumllz\nROJrMUEvJcL8vHncqBavHR61lf59JbGaF4Mh5b6aCXlE4qs5mO37BeVJfDXPGC3vWF886a8WG2pI\nao5wPVmwFtfMk/5q8rZ3bZL7CszskcDH6Nr5ucC/BD4DvMPMXjDLY8CVwLnAi4CnzfJfZWYnDVUX\nIYTYyWzMSfX5GwtDjqSeQffj4fyU0h2zbR+ddV7PBt4KXACcCZyTUroKwMz2AdfTBdB68YD1EUKI\nHYvkvkPZTacGfK/YfhuwZ5a+APjmRgcFkFK6zcyuBC5kG53UvJ4+Guep5j0i4pliWfGkSiLh36Mx\nlxbpYNaz4os6II3Ek6o5fu3rsLNm3df3xY/cw3JfxMtEyy/fmpTptWcpr3nePVpk0b7vTs0arq90\nHyUSg65mVelJf+X97etYe50Y8jv1DwEDfsfMHmJmx5rZc4HHA2+c5TkduHZO2euAk83s6AHrI4QQ\nOxLJfXNIKV1rZo8DPgC8cLb5buD5KaX3zD4fB9wwp/its/97gDvm7MfMLgUuha4nFEII4TOmjqYP\ng3VSZvZw4E/oRkXPB75PJ+G9zczuSim9s8/xU0p7gb0A9zEL+VJsWZgbzeeVidLyAG1X4ovGKOor\n95UsQyYt2887j7coFXwZrRaDqiV2WMvCa09Gi1oHekRlVu/Zr8ms+QLe8tmLnLOlbkPOu0QtAlus\nOj2rv/Kzt6C6LLOVdC0T9Pm8jq7tzk8pbbTvx83seOBNZvZuYD8H56dyjpv93z9gfYQQQqw5Q85J\nPQL4fNZBbfCXwPHAj9GNsk6fU/Y04MbMKlAIIUQPNCd1KN8CHmlmu4uO6mfolIBbgSuAi83srJTS\n1QBmdgxwPvCuvhWILrys+ZCLyA/L8iUVtejxJL4Wf3C1hbkHnO2lBBTx19fXwqq2eNWT+Mpr82Sj\nmnXfkPG+aseNWPHV6hLJV7NcjCxKLffVZNKcqNWbdx6P8jmMUIvZ1JeoRXHEd1/NCtmLLSW571De\nDLwXuNLM3kI3J3UB8EzgjSmlA2Z2BZ2HicvN7OV08t4r6WwhXj9gXYQQYkczptFQH4a07nufmT0Z\n+FXg9+k84nyVztLv7bM895jZecBlwFtmefYBZ6eUvj5UXYQQYiczJS/og/ruSyl9CPjQFnluBS6Z\n/QkhhBAua+8FPRrnKWpS62nffc3Ma+W9Xzw1LwQtTku9+SXPFLu2bxXxtVrmioY09Yf+1xadk9qu\nA+FFzj/UlhF4XhG8+S2IOaUt2e7713d+ami82FLgL3+ILoXw7r3mpIQQQowSyX1CCCFGizqpEVEb\noke8R9TM1hcl/dWoyXWRmElRE3TPo0FZxjNljoaC7xufy3N02lfGK+vgxZCqlelL7RoiJujRZyWy\nvSTq4cF7JrznsyTqiDbyLnrPyhBEj+fVs9bunnl6LeT8PNl3qu7i1r6TEkIIcSiakxJCCDFKJPeN\ngHnD6qj3iJrsFIknVTKkxVdLPKgWqyrPe0RLePBoTK6+1ls1J50e0Rc1aok55K/Tlrq1OJVtcT7r\ntUdUMvXCuoP/HLa8Y0PGoIpeWws1jxMt8nhEOtdISgghxCiZ0khqWW7ohBBCiG2zliMp42Dv2mKp\nF40nFcWzDIuW8agt8IzEeao5i/Ws+6JWSEO34XZZ5K/E2v1c1q/T6KJfL/+y6hmxros6svWk5do+\nz3KwJWR9bYH4kM90rT0ioeTLzxvlS+u+RT8DZnYR8CzgMcAJwI3A+4HXpZRuH+o8a9lJCSGE8FmS\nF/SXAd+gcxJ+E/Bo4NXA2Wb2cymlQaqgTkoIISbIEkbT56eUbs4+f9LMbgX+CHgc8IkhTqJOSggh\nJsYyDCeKDmqDz8z+nzjUeda+k4rOjUTNyYckOtZt8SKwXXN0iAU9rJnhem24CtPhIXSEyFxi7UXv\nO0/R8iXS10FtFM/8ujaf4pWpzZd5c0otyxUic1Wt+RY1PwX+vFjUwezIrN/Omv3/26EOuPadlBBC\niEMZ4IfcCWZ2TfZ5b0ppr5fZzE4EfgP4WErpGi/fdlEnJYQQE2Mgue+WlNIZkYxmdjTwQeCHwMX9\nT32Qte2kNoa7Ndkp4mWiZrbe4hB1SImvPJbnPLbFBD0ai2hR3iOieC9aTXbqe56h7/WQtMTNavFS\nkVNzuOvti3pY8OpWPruHVfbNO26U3DNG7X7WPFNsl+iyhu0urylN0Jf1fJrZkcCVwKnAWSmlm4Y8\n/tp2UkIIIVaLme0C3gecAfx8SukLQ59DnZQQQkyMZVj3mdlhwDuBc4DzUkqfXsR51rKTMg4OcWsS\n1G5nX1S2GpKW2FC1MtuNN1Qe2ztPzUIy6qTXo5bPe6GiVnctsaEidVnW81FjSCvAvsdt8cpQ1qXv\ns+vVx5MBwX9ea20QCUE/xPMRkT9rHie8617COqnfA54OvBa408x+Ntt301Cy38isF4UQQvRlw+NE\nn78AT5r9fxWwr/h7zhDXAWs6khJCCFFnCYt5T1nwKYAJdFJR2anvkDEqAbU4APXS0VDuEUu/WvkW\nJ701WqTACLV7sKgXsrxvq5YeWqxHhz52TiTmUjTeWE757HploouBa4t2I+TyWouE3JeIowKFjxdC\nCLEWTCmelDopIYSYGEvygr4U1EkJIcQE0UhqxcybG4hq30POT9WorfSPeJwoy2zXQWzNhD2nNkcQ\naau+gSJr52kxJx/SO0DJqs3TW65nWWVa7o/3TEZN0L3z952Dqh27b77oCKc297bVvN6URlKrngcW\nQgghXNZ2JCWEEMJHct8KiXqcWIUcE5Xu8nwHnO1Rs/WIJ4rycyQmT8miTMtLvHpG4j+V5fvSstxg\nFSzS5LyvaXrteBEHteW+iLeWKDVpeLsSdCuRexd5L3MTdFn3CSGEGDWakxJCCCEWzNqOpDZ615rc\n11ee8ob1tV8oUYeZ3r5aPKiIZ4m+caJKluUxYrtlWn4ltvwiW1fJpKV9Itc6RHt4El9UIqzJgt6x\nvDJDWwF6DG1xutX7I7lPCCHEqFEnJYQQYpRMaZ3U2ndSy5pU6yvx1RbmtkiEEWeztZhL3vZVTVJG\nLKaiC0Qj59iptCxEXyTefYuGqffkumXJeEMQsQhuefanMpKS4YQQQojRsvYjKSGEEJuR3CeEEGLU\nTEXuW8tOKvc4kbOs1d+1MlGPE5GghaUJesSzRK3+y/IYkeOZGNfmy6LB87zziDa8dq+17aLavWZC\n7j1HtTngRc1r1OaNW5ZZ1I69HWSCLoQQYtRMRe6T4YQQQojRsmNGUtGhb2SIXnMWGzENj5apxdSJ\neJaIyg1940kNQUQm6fvLMHotq/gF2redh67zIh3WthAxs454oiipPftDymV9j1WbLph3bMl9Qggh\nRos6KSGEEKNmKnNSk+6k8ps0pIQU9QSRy3NlmQNOvqhFYNSiL6fmjHdIIvGCWuI09a3zmF/aZdVt\nWaHkl0XE0q9WZmjPFN4z3mLpV7ueMT/LQzPpTkoIIXYikvuEEEKMmqmMtibVSZU3JRLDpRyG913Q\n6El8tYW5UYlwuxZ9Y1hfEJH+SqbuMHOsDGnR1/de9bXEbFmAXJPko8f24s71bdua9L+RL2XbNJIS\nQggxaqbSSY3hx7YQQggxl1AnZWYnmdnvmtk+M/uemSUzO2VOvj1m9vtmdouZ3WlmHzOzR8zJd4SZ\nvcHM/sHMvj877v/W/3KEEEJseEHv8zcWonLfw4BfBD4L/DnwhDKDmRlwJXAK8CJgP/BK4Coze3RK\n6aYs+zuApwAvB74GvBD4sJn985TS59ou5VCGXuW9QXSu6ICzvbYvqmMvayjfEowwYu7f11nsshzk\nrooxSTXRLyyvzov0XtESINPbV3sm8/d1t7MdfEfJUefOLd85Q7b7GIl2Up9KKT0IwMyew5xOCrgA\nOBM4J6V01SzvPuB64BXAi2fbHgX8a+CSlNJ/mW27GrgO+I3ZcYQQQjQyJcOJkNyXUor8sLkA+OZG\nBzUrdxvd6OrCIt/dwB9n+X4IvAd4opkdHqmTEEIIn50m90U4Hbh2zvbrgGeb2dEppTtm+a5PKX1v\nTr7ddNLidbUT5b8SImbmNfpKGeUxPBkv6pS2ds4xPTg17xF9ncVGfjlN5Vfisun7DEXbPfq+LIvI\nEpRovTzpL3r+lu+cyPdHYpoM2UkdB9wwZ/uts/97gDtm+fZX8h037+BmdilwKUx/PkIIIfowJblv\nbdZJpZT2AnsBdptN9UeDEEIMwpiUlz4M2UntpxstlRyX7d/4/+OVfLfO2beJu+GWm+BO4JbtVnJi\nnMDOboOdfv2gNgC1wcb13/u9eg98+M5uex9G0aZDdlLXMd/q7zTgxtl81Ea+XzCzo4p5qdPo5N6v\nbHWilNIDzeyalNIZfSu9zuz0Ntjp1w9qA1AbzLv+lNK5q6rP0AzpceIK4EQzO2tjg5kdA5w/27fB\nlcAu4OlZvvsC/wr4SErpBwPWSQghxBoTHkmZ2UWz5GNm/59kZjcDN6eUrqbriPYBl5vZyzm4mNeA\n128cJ6X0/5nZHwO/bWa76NZRvQD4CeDf9LweIYQQE2I7ct97i89vmf2/GnhcSukeMzsPuGy27wi6\nTpE3PbwAAAWGSURBVOvslNLXi7IXA68FXgMcC/w1cG5K6a+2UZ+928g7VXZ6G+z06we1AagNJn39\nlpIM5YQQQowTeUEXQggxWtRJCSGEGC1r1UmZ2UPN7H1mdpuZfdfM3m9mJ6+6XovAzC4ysz81s6/P\nwpl80cz+k5ndv8gXCo8yBczsz2ZhYl5TbJ90G5jZk83sU2Z2x+y5v8bMzsn2T/36zzSzj5jZ/zCz\n283sr8zskiLPJNpAYZEOZW06KTM7CvgE8D8B/xb4JeDhdKFA7rfKui2Il9F5Nnkl8CTgrXRWkB81\ns8NgU3iUc+nCozyNzrz/KjM7aRWVXhRm9kzgUXO2T7oNzOx5wAfpwuT8At3SjfcCR832T/36Hwl8\njO6angv8S+AzwDvM7AWzPFNqg42wSPvpwiIdwjav9x107fYfgfOAf6ALi/TohdR+EaSU1uIP+Pd0\nX9oPy7b9BPBD4KWrrt8CrveBc7Y9m84t1zmzzxfOPp+d5XkAndeO31n1NQzYFnuAbwHPnF3va7J9\nk20Duths3wf+90qeyV7/7FpeR7fI/+hi+z5g39TaADgsSz9ndl2ntNxzuh91Cbg423Zf4IvAFau+\n1ujf2oyk6EJ8fDqldK9HipTS9cB/Y3MokEmQUrp5zubPzP6fOPsfDY+y7vwmcG1K6d1z9k25DS6h\nc8H2tkqeKV8/dI7GDwBl1ITbOKgETaYNksIiHcI6dVK1UCCnLbkuq2LDm8ffzv7X2uRkMzt6KbVa\nIGb2WLoR5AudLFNug8cCfwc8w8y+amY/NLOvmFneFlO+foA/pHMI8Dtm9hAzO9bMngs8HnjjLM/U\n26Aker2RsEijZ506qVqIj3mObSeFmZ1IF7n4Yymla2abtwp7stbtYma7gbcDl6WUvuhkm3IbPIRu\n3vUNwP9F5xvzo8Cbzezfz/JM+fpJKV0LPA54KvANumv9PeD5KaX3zLJNug3mEL3eprBIY2NtQnXs\nZGa/jD5IN/928Yqrs0xeARxJ551kJ3IYcH/gl1NK759t+8TM2uuVwJtWVK+lYWYPB/6E7tf/8+nm\n6C4E3mZmd6WU3rnK+onFs06dVC0UyLxfC5PAzI6k05pPBc5KKd2U7Y6GR1k7ZksLXkU3eXx4oZ8f\nbmbHArcz4TYA/pFuJPXRYvtHgHPN7J8w7euHznDibuD8lNJGUNyPm9nxwJvM7N1Mvw1KlhYWaQys\nk9x3HZ3GWnIa8DdLrstSmDngfR9wBvDklNIXiiy1NsnDo6wjp9L5f7yc7mXb+IPOPH8/8Aim3QbX\nBfNM9fqhu8efzzqoDf4SOB74MabfBiXR670O+InZ8p0yXygs0hhYp07qCuBnzezUjQ0z2eNMNocC\nmQSztVDvBM4BnppS+vScbNHwKOvI54Cz5/xB13GdTfeSTbkNPjD7/8Ri+7nATSmlf2Da1w/d0oNH\nzuYnc34GuItuNDD1NijZWWGRVm0DH/0D7kf3pfQFOk36Ajrv6V+jWEMxhT+6xbuJzlP8zxZ/J83y\nHAb8d+DrwDPovsw+SffiPnTV17CgdinXSU22Deis2j5BJ/s9n85w4j/P2uCXp379s+u7aHa9H569\n908A3jzb9ltTbIPZNV+UfQe8YPb5rO1eL525+X462fzxdMrMXcBPr/o6w+2x6gps8+adTDeJ+l26\n+Yg/pVjoNpU/4IbZAzrv79VZvuOAP5g9oN8DPg48atX1X2C7bOqkpt4GwDF01mzfppNoPg/8651y\n/bPre9LsS/jm2Xv/OeBXgPtMsQ0q7/0nt3u9dIZHv0U3Ir0L+Au60Eorv87on0J1CCGEGC3rNCcl\nhBBih6FOSgghxGhRJyWEEGK0qJMSQggxWtRJCSGEGC3qpIQQQowWdVJCCCFGizopIYQQo+X/B/qk\nGBwNrfWtAAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#NBVAL_SKIP\n", "\n", "# Plot percentage error\n", "plot_image(100*np.abs(model0.vp-get_true_model().vp.data)/get_true_model().vp.data, vmax=15, cmap=\"hot\")" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGWCAYAAACqxYPqAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xd81dX9x/HXJ4uEvcMMYe9pWAEEBFFR0TpwoIK1buvo\nsFprW6tWa61atVq1/QlWxK2IAwUURAOy9x5h7xE2ZJzfH99LG0PGDdzke5P7fj4e30eS73xfRfPh\nnPM9x5xziIiIiPglyu8AIiIiEtlUjIiIiIivVIyIiIiIr1SMiIiIiK9UjIiIiIivVIyIiIiIr1SM\niIiIiK9UjIhEIDPraWajzWytmR0JbCvN7FUz6+V3PhGJLKZJz0Qih5lFA88AdwMngCnAMiAbaAEM\nBqoCNzjn/uNXThGJLDF+BxCRUvVnvEJkLnCFcy4990EzqwL8Bqhe+tFEJFKpm0YkQphZS+BXwG7g\ngryFCIBz7qBz7nfAq7mu62xmH5jZLjM7bmarzexxM6uc5/7JZuYC3T8tzOwjM9tnZofNbLKZdc5z\n/loz22tmcQXkXRi4tnKufRXM7NeBY0fMLCNw7/75XD81kCfezP5sZuvNLNPM7s11ziAzSwvca6eZ\nvW5mtcws3cxO+edjZlXN7DEzW2Fmx8xsj5l9nPezBc5ND2yVzezvZrY18M9vkZldUcBnrmZmj5rZ\nUjM7Gvjn84OZ/Sqfcwea2eeBDMfMbJmZPWBm+kumlDkqRkQixyi8/+Zfcc7tKuxE59xxgMAv+RnA\nRcBE4FlgL/BbYKqZJeRzeTIwE6gJ/B8wCRgEfGNmibnOewuoAQzNewMz6wB0AsY75w4F9sUH7vUU\ncByvYHo3cN4UM7usgI/zEXA9MBl4HtgcuN9Q4EugMzAOeB04K/CM2Hwy1Q58roeArcA/gAmBz5Zm\nZr3zeXYs8BUwBPgAeBNoDrxrZkPy3D8RmA38DjgMvBg4/zDwYJ5z78LrYksBxgeyHASeCPwzESlb\nnHPatGmLgA34BnDAOUGeHw2sBXKA/rn2GzA6cK8/5NqfHNjngN/kudejgf0P5NrXOrDv/Xye/UTg\n2NB89j2Q59w6QDqwC0jItX9q4Py5QPV8PttGIBPolmf/V4Hr0vNcMy6w/5o8+1sAGcDiPPvTA+d/\nDMTl2j8osH9invM/DOz/XT7/PBrm+r59IPdMoFqefy8vBu5xhd9/3rRpK87mewBt2rSVzgYsD/yi\nah3k+f0D54/P51gDvNaJdbn2nSxG1gFRec4/eeyDPPvnAMfy+aW6AdgJxAT2RQH7gKUFZL0rcP+L\ncu2bmndfrmMDAsfeyedYz7zFCFAbb5Dv5wU8/+nANR1y7TtZjDTN5/x0YE+un+vhFX0rgOgi/r08\nH7hvj3yOVQ3c55QCT5u2cN7UtygiBTk5DmJa3gPOua1mthpob2ZVnHMHcx1e4JzLyXPJ5sDXvANj\n38Tr+rkcr0sHoC+QBLzonMsK7GsduHajmf0xn6wtA1/bAJ/mOTYnn/NPfrYZ+RybDWTl2dcdryCq\nVMDz2+Z6/pJc+/c759bnc/5mIHe3TgpeEfa1cy47n/NzO1ksXRToasrraCCHSJmhYkQkcmzH+yXV\nEFgZxPlVA193FHK/9oHzchcjB/Ke6JzLMjPwukFyG4fXqnAd/ytGRgS+vpnrvJqBr50CW0Eq5bMv\nv/xVAl9PGTvjnMsxs915dp98/tmBLdjnZxRwXhY/HrNXLfB1ayH3zp3FgIeLkUMkrGkAq0jkSAt8\nPSfI808WFYkFHE/Mc16xOed24A3E7G9mDc0sFrgSWOOc+yGfLO8456yQ7ZF8npHfZEoni6c6eQ+Y\nWRRQK8/uk8//SxHPH1OMj5/b/sDXBkGcewCvZaRSITmanmYOEV+oGBGJHKPxxhPcEngzpEBmVgFY\nEPjxlJYAM6sPtMIbM3Iw7/FiehPv/0XXAhfg/c1/bJ5zluMVECmBidvO1MLA1/xmm03h1LdpZuMV\nACU1O+3cwP3PCeLzzcJrGelRQllESp2KEZEI4ZxbjdclUgf4zMya5D3HzKqY2ePALcB3eINRh5lZ\n3zynPg7EAW+EINpHwBG87pn8umgIjB35J95rsU/k9ws7MMV9xSCfOR1v3MbluecICdz3T3lPds5t\nB97Ha8H5eT7PtvzmOglW4P4f4Y2NeTDvcTNrmOvHl/AG075gZqe0pJhZopm1zbtfJJxpzIhIZPkt\nkAD8HFhlZpPxpoPPwftFfy7eGJDrA2Mnfoo3v8hkM3sX2IL3JkovvL/NP3WmgZxzh8xsPHAN0A6Y\n5Zxbk8+pv8drtfg1XoE0HW/Ok0Z484O0BurjFTZFPTPbzG7He+32ezMbB+zBa5nJxhu7kXcQ7u14\nY26eN7MbgR+AQ3iDbXvhdVvFF+Oj53UH3niYR81sGN7bQHF443K6Eeg6cs4tDhREL+L9O/wM7+2c\nGnivGffDG0+y/AyyiJQqtYyIRBDnXLZz7m68NznG4f1yvRO4B+8Nk/eBns65NwPnTwuc+zlwIfBL\nvJaVJ4ABzrmjIYp2siUkljytIrmyH8ObPOwuvDEWV+FNbd8b75XYkXizywbFOfcpXvGxGG8A7U3A\nPP5XkB3Mc/6ewLMewusmuQGvQOmGNx7n2mCfXUCeHXhdL0/ivTl0D95kbVWAx/Kc+zLeW0df4HWj\n3QcMAyrizeny1plkESltWihPRCQXM2sOrAHec84N9zuPSCRQy4iIRKTAmjF519epAPwt8OP40k8l\nEpnUMiIiEcnMUvCmyJ8IrMcbc3EO0AxvordBQUxAJiIhoGLkNAWacscAdfEWsrrZOZffTI8iEoYC\nC9M9hTfmIhGvpXgd3kJzfwnheBgRKYKKkdNkZpOAd51zr5nZuXgj29sUMMGSiIiIFMCXMSNm1sfM\nvjKznWZ20MzmBV4hLKnnNTKzF8xshpkdMTNnZsmFnN/YzN43swwzO2BmH5pZUq7jdfBe5RsN4Jyb\nhDe6/qyS+gwiIiLlVanPM2JmnYDJeMtf34w3J8AVwL/NrELglbVQawEMx5sXYTre64EF5asIfI23\nIulIvFkRHwO+MbNOzrnDePMKbHPOZea6ND2wP9+umtq1a7vk5OQz/iAiIiJlxdy5c3c7505ZdiEv\nPyY9uxpvsayLnXOHAvsmBYqUG4B8i5HA9NNNnXNp+RyrCPR3zn1RwDO/dc4lBs79GYUUI3gFUjO8\nZdbXBK5ZBKwGbgWeKeLz5Ss5OZk5czSkREREIoeZbQjmPD+6aeKAE5w6S2IGhee5F/gy77TUZpYA\nTABez/ua3kn5LGdemGHAzNwzQAaWAP8euCSwayNQP7Co10nJgf0iIiJSDH4UI6Pxxlc8b2YNzKy6\nmd0MDAKeLeS6h/BewfvCzPrAjwqRzsCQXC0tZ6I9sCSf/UvxpqrGObcLb7GqUYEc5wY+09wQPF9E\nRCSilHox4pxbgre2xaV461zsA/4B3Oace7uQ67Lw1q74Eq8gGYRXiHTBmw9gUYgi1gxkymsv3jwE\nJ90G3Ghmq4C/AiPye5PGzC42s1czMjJCFE9ERKR88WMAa0vgA7yWhtuAo3jdH/80s2POubxLh/+X\ncy7LzK7GmwdgMl7RcI5zbmFB15SUwAqoqUGcNwGYkJKScnPJpxIRESl7/BjA+mcgE28A64nAvilm\nVgv4u5mNK2KMRwxQGW9FzVi8haFCaR8/bgE5qaAWExERETkDfowZ6QgsylWInDQLb4nsugVdaGbx\neEt+d8Ob52MSMNHMeocw31K8cSN5tcNbal1ERERCyI9iZDvQyczi8uzvCRzDG5txisACVh8B3YHB\nzrnZeEuIT8YrSHqFKN8nQC8za5br2clAn8AxERERCSE/ipEX8ebxmGBml5jZEDN7EW9w6sv5tJic\n9DjQA68QWQAQmHTsKmBK4H75vtoLYGZXmNkV/G+W1AsC+/rnOfU1vAnMxgfyDcNbvXMT8MppfF4R\nEREphC9r05jZBcBv8LpD4oG1wKvAKwWtkmlmNYBGzrnF+RyLBVKcczMKeWZBH3Sac25AnnOT8F4z\nPvnK7hTgXudceuGfrGApKSlOk56JiEgkMbO5zrmUos7zYwArgZlSC5ottaBr9lHAANJAC0mBhUjg\nHCvGszYClxcnn4iIiJweXxbKExERETlJxYiIiIj4SsWIiIiI+ErFSBn1zcqdZGYXZ/0/ERGR8OTL\nAFY5M0u2ZHDj67OpW6UCI3o24dqeSdSpUsHvWCIiIqdFLSNlULv6VXl9VHfa1q/Ks5NXkfrkFO59\nez7zN2q2ehERKXvUMlIGRUUZA9vUZWCbuqzbdYg3Zmzg/bmb+XjBVjo3rs6o1CYM7VifCjHRfkcV\nEREpki+TnkWikp707NDxLD6ct5kxaems3XWY2pXjuLZHEiN6NSGxanyJPVdERKQgwU56pmKklJTW\nDKzOOb5bs5sxaelMWbGTaDPO71CPUanJnNWkBmZBz/0mIiJyRsJ6BlYpOWZGv5Z16NeyDhv3HOGN\nGem8M2cTny7aRvsGVRmZmsywzg2Ij1UXjoiIhAe1jJQSP9emOXIii4/mb2FMWjqrdhyiRsVYru6R\nxHW9mtCweoIvmUREpPxTN02YCYeF8pxzzFi3hzFp6UxatgOAIe3qMTI1mV7NaqoLR0REQkrdNHIK\nMyO1eW1Sm9dm874jvDlzI2/P3sjEpdtpU68KN/RO5iddG5IQpy4cEREpPWoZKSXh0DKSn2OZ2Xyy\nYCuj09JZtu0A1RJiuap7Y67v1YTGNSv6HU9ERMowddOEmXAtRk5yzjFnwz5Gf5/OxKXbyXGOQW0S\nGZWaTJ8WtdSFIyIixaZuGikWM6N7ck26J9dkW8ZRxs7cyLhZG5m8fAct6lZmZO8mXNatEZUq6I+M\niIiEllpGSkm4t4zk51hmNp8t2saYGeks2pxBlQoxXJHSiJG9k0muXcnveCIiEubUTRNmymIxcpJz\njvmb9jMmLZ3PF28jM9sxoHUdRqYm079lHaKi1IUjIiKnUjESZspyMZLbzgPHeGvWRsb+sJFdB4/T\ntHYlbujdhCvOakSV+Fi/44mISBhRMRJmyksxctKJrBy+WLKNMWnpzNu4n0px0Vx+ViNu6J1Mi7qV\n/Y4nIiJhQMVImClvxUhuizbvZ3RaOp8u3MaJ7Bz6tazNyN7JDGxTl2h14YiIRCwVI2GmPBcjJ+0+\ndJy3Z23kzZkb2X7gGI1rJnBDr2SGpzSmWkV14YiIRBoVI2EmEoqRkzKzc/hq6Q7GpKUzK30vCbHR\nXNq1IaNSk2ldr4rf8UREpJSoGAkzkVSM5LZs6wHGpKXz8YItHM/KoVezmoxKbcrgtnWJiY7yO56I\niJQgFSNhJlKLkZP2HT7BO3M28Z8ZG9iy/ygNqydwXa8mXN29MTUqxfkdT0RESoCKkTAT6cXISdk5\njsnLvS6ctLV7qBATxSVdGjAyNZn2Dar5HU9EREJIxUiYUTFyqpXbDzJmRjofzdvC0cxsuifXYGRq\nMue1r0esunBERMo8FSNhRsVIwTKOZPLe3E28MWMDG/ceoV7VeEb0TOKanknUrlzB73giInKaVIyE\nGRUjRcvOcUxduZPRaelMX72buOgoLupUn5GpyXRuXN3veCIiUkxatVfKnOgoY1DbRAa1TWTtrkO8\nkZbO+3M38+H8LXRpXJ0b+yRzQYf6xMWoC0dEpDxRy0gpUcvI6Tl4LJMP5m7mjRkbWLf7MHWqVODa\nHkmM6JlE3arxfscTEZFCqJsmzKgYOTM5OY7pa3YzJi2db1buJNqMoR29LpxuSdUx07TzIiLhRt00\nUq5ERRn9W9Whf6s6pO8+zH9mbuDdOZv4ZOFWOjasxsjUZC7qVJ/42Gi/o4qISDGpZaSUqGUk9A4f\nz+Kj+VsYk5bO6p2HqFkpjmt6NOa6Xk2oXy3B73giIhFP3TRhRsVIyXHOkbZ2D6PT0pmyfAdmxnnt\nExnZO5keTWuqC0dExCfqppGIYWb0aVGbPi1qs2nvEd6cuYG3Z2/i88XbaVu/KiN7N+GSLg1JiFMX\njohIOFLLSClRy0jpOnoim/ELtjA6LZ0V2w9SvWIsV6V4XTiNa1b0O56ISERQN02YUTHiD+ccs9bv\nZcyMdL5cugPnHIPbJjIqNZnezWupC0dEpASpm0YErwunZ7Na9GxWi637jzL2hw2Mm7WJr5btoFVi\nZW7oncxl3RpSMU7/KYiI+EUtI6VELSPh41hmNp8u2saYtHQWb8mgSnwMw1Mac0PvJjSpVcnveCIi\n5Ya6acKMipHw45xj3sb9jE5L54vF28h2joGt6zIyNZl+LWoTFaUuHBGRM6FiJMyoGAlvOw4cY+wP\nG3nrh43sPnScZrUrcUPvJlx+ViOqxMf6HU9EpExSMRJmVIyUDcezsvli8XZGp6WzYNN+KleI4fJu\nDbkhNZnmdSr7HU9EpExRMRJmVIyUPQs37WdMWjqfLtrGiewc+rWszajUZAa2rqsuHBGRIKgYCTMq\nRsquXQeP8/asjbz5wwZ2HDhOk1oVub5XE65MaUy1BHXhiIgURMVImFExUvZlZufw5dLtjElLZ3b6\nPhJio7msW0NGpSbTMrGK3/FERMKOipEwo2KkfFmyJYMxaemMX7iVE1k5pDavxcjUZAa3TSRaXTgi\nIoCKkbCjYqR82nv4BG/P3sibMzawNeMYDasncH3vJlzdvTHVK8b5HU9ExFcqRsKMipHyLSs7h8nL\ndzA6LZ2Z6/ZSISaKS7s0ZGRqMu0aVPU7noiIL1SMhBkVI5FjxfYDjEnbwEfzN3MsM4ceyTUZ1SeZ\nIe0SiYmO8jueiEipUTESZlSMRJ6MI5m8O2cTb8xMZ9Peo9SvFs91vbwunFqVK/gdT0SkxKkYCTMq\nRiJXdo7j6xU7GZOWzndrdhMXE8XFnRowKjWZjo2q+R1PRKTEaNVekTARHWWc2y6Rc9slsmbnQcak\nbeCDeZv5YN5muiVVZ2RqMhd0qE9cjLpwRCQyqWWklKhlRHI7cCyT9+ds5o0Z6aTvOUKdKhUY0TOJ\na3smUbdKvN/xRERCQt00YUbFiOQnJ8cxbfUuxqSlM3XlLmKjjaEd6zMqNZmuSTX8jicickbUTSNS\nBkRFGQNb12Vg67qs332YN2ak8/6czYxfsJXOjaoxMjWZCzvVp0JMtN9RRURKjFpGSolaRiRYh45n\n8dG8zYxOS2ftrsPUrhzHNT2SGNGzCfWqqQtHRMoOddOEGRUjUlzOOb5bs5sxaelMWbGTaDPO61CP\nUanJpDSpgZmmnReR8KZuGpEyzszo17IO/VrWYeOeI/xnZjrvzN7EZ4u20a5+VUalJjOsSwPiY9WF\nIyJlm1pGSolaRiQUjpzI4uP5WxmTls7KHQepUTGWq7oncX3vJjSsnuB3PBGRH1E3TZhRMSKh5Jxj\n5rq9jElL56tl2wE4t10iI1OT6d2slrpwRCQsqJtGpBwzM3o3r0Xv5rXYsv8ob87cwNuzNvLl0h20\nTqzC7y9uR58Wtf2OKSISFE35KFLGNayewG/Ob8OMBwfx1BWdOJGdw63/mUv67sN+RxMRCYqKEZFy\nIj42muEpjXnzZz2JiTbuGDuPY5nZfscSESmSihGRcqZh9QSeGd6ZZdsO8MiEZX7HEREpkooRkXLo\nnDaJ3D6gOeNmbeSj+Zv9jiMiUigVIyLl1C/PbUWP5Jr89sMlrN5x0O84IiIFUjEiUk7FREfxwrVd\nqRgXzR1j53HkRJbfkURE8qViRKQcS6waz9+v7sqaXYf43UdL0LxCIhKOVIyIlHN9W9bmnkEt+XD+\nFt6ZvcnvOCIip1AxIhIBfn5OS/q1rM3vP1nK0q0ZfscREfkRFSMiESA6ynj2qi7UqBjLnWPncfBY\npt+RRET+S8WISISoXbkCL1zTjU37jvLAB4s1fkREwoaKEZEI0qNpTX59Xms+W7yNN2Zs8DuOiAig\nYkQk4tzSrxmD2tTlsc+WsXDTfr/jiIioGBGJNFFRxt+Gd6ZulXjuGDuPjCMaPyIi/lIxIhKBqleM\n48Vru7Lz4DF++d4CcnI0fkRE/KNiRCRCdU2qwW+HtmXy8p28Nn2d33FEJIIVWYyYWZyZbTWzi0oj\nkIiUnlGpyQztWI+nvlzJ7PS9fscRkQhVZDHinDsBxAHHSj6OiJQmM+PJyzvRuEYCd701jz2Hjvsd\nSUQiULDdNBOAy0oyiIj4o2p8LP8Y0Y19RzK5950FZGv8iIiUsmCLkQ+Ai83sTTO7wsz6mFlq7q0k\nQ4pIyWrfoBqPDGvP9NW7+cc3a/yOIyIRJibI8z4JfL02sOX+q5MFfo4OYS4RKWVXd2/MrPV7eXby\nKs5qUoM+LWr7HUlEIkSwxcgFJZpCRHxnZjx2aQcWb8ngnrfn89nd/UisGu93LBGJAKb1KUpHSkqK\nmzNnjt8xRIq0esdBhr34PR0bVeOtn/UkJlozAIjI6TGzuc65lKLOK9b/ZcysipkNMrMrA1+rnH5E\nEQlHLROr8PhPOjBr/V6embTK7zgiEgGCLkbM7HfANuAr4B1gErDNzB4qoWwi4pPLujXimh6NeWnq\nWr5ZsdPvOCJSzgVVjJjZncCfgI+AoUBXvHEkHwF/MrPbSyyhiPjiDxe3p239qtz37gK27D/qdxwR\nKceCbRm5C3jJOXe9c+5L59zCwNfrgZeBn5dcRBHxQ3xsNC+N6EZWtuOut+ZxIivH70giUk4FW4w0\nA8YXcGx84LiIlDNNa1fiL5d3Yv7G/fxl4gq/44hIORVsMbIXaF3AsdaB4yJSDl3YqT6jUpP593fr\nmbhku99xRKQcCrYY+Rh4PPAWjZ3caWY/AR4NHBeRcurBoW3o3Kgav35vIRv2HPY7joiUM8EWIw8A\nK/DeojliZhvM7AjwPrAycFxEyqkKMdG8eG03zOCOsfM4lpntdyQRKUeCKkaccxlAKjAceBX4FngN\nuBLo45w7UGIJRSQsNK5ZkWeGd2Hp1gM8+ukyv+OISDlS5HTwZhYH3AhMd869j9caIiIRaHC7RG7t\n34xXpq2jR9OaXNKlod+RRKQcKLJlxDl3AngOqFPycUQk3P1qSGu6J9fgwQ8Xs2bnIb/jiEg5EOyY\nkZVAUkkGEZGyITY6iheu6UZ8bDR3jJ3L0RMaPyIiZybYYuQR4Pdm1qokw4hI2VCvWjzPXdWF1TsP\n8fD4JX7HEZEyrsgxIwF3AJWBpWa2Em+NmtzL/Trn3HmhDici4evsVnX4+TkteX7Kano0rcnwlMZ+\nRxKRMirYYqQKkB7YwCtMRCTC3TOoJXPS9/Lwx0vo2LAabetX9TuSiJRB5pwr+iw5YykpKW7OnDl+\nxxAJuV0HjzP0+elUrhDDJ3f1oUp8rN+RRCRMmNlc51xKUecVOWbEzOLMbJyZ9Q1NNBEpT+pUqcAL\n13Rlw57DPPjhYvQXHBEprmBf7b0IiC75OCJSFvVqVotfndeaTxdt482ZG/yOIyJlTLBv0/wA9CjJ\nICJStt12dnMGtq7Do58uZ9Hm/X7HEZEyJNhi5B7gFjP7mZnVLslAIlI2RUUZzwzvQu3Kcdwxdh4Z\nRzL9jiQiZUSwxcgCoCnwCrDDzDLN7ESu7XjJRRSRsqJGpTheHNGN7RnH+NX7CzV+RESCEuyrvX/j\nx/OKiIjkq1tSDR4c2pZHP13Gv6av5+azm/kdSUTCXFDFiHPugZIOIiLlx0/7JDNr/R6enLiCrknV\nSUmu6XckEQljwXbT/FfgVd9EM9PbNSKSLzPjqSs607B6Ane9NZ89h9STKyIFC7oYMbNzzSwNOAxs\nAToH9r9kZsNLKJ+IlFHVEmJ5aUQ39h45wX3vLiQnRz29IpK/oIoRMxsKfAFk4i2al/u6bcCNoY8m\nImVdh4bV+MPF7fh21S5emrrG7zgiEqaCbRn5EzDWOdcfeDLPsUVAx5CmEpFy49oeSVzSpQHPTFpF\n2trdfscRkTAUbDHSHngz8H3ettZ9gOYeEZF8mRl//klHmtauxN3jFrDz4DG/I4lImAm2GDkEFDQc\nPgnQX3dEpECVKsTw0oizOHQ8k7vHzSdb40dEJJdgi5EpwG/MrHKufc7MYoE7gC9DnkxEypXW9arw\n2KUdmbluL89OWuV3HBEJI8FOevYQ3vo0K4AJeF01v8B7oyYRuKpE0olIuXLFWY2YtX4PL36zhrOS\nazCwdV2/I4lIGAiqZcQ5txboDkwDrgQMbyXfpUAv59ymEksYpsysuZl9Z2arzGy+maX4nUmkLHhk\nWAfa1KvCL95ZwNb9R/2OIyJhIOh5Rpxz651zI5xztYEY51wN59w1zrn1JZgvnP0TGOOcawXcD4w1\nM/M5k0jYS4iL5h8junEiK4e73ppHZnaO35FExGfFnoEVwDl32v/3MLOpZuYK2Cae7n2LeGYjM3vB\nzGaY2ZHAs5ILOb+xmb1vZhlmdsDMPjSzpFzH6wC9gNEAzrlJeK1FZ5VEfpHypnmdyjx5eSfmbdzP\nUxNX+B1HRHx2WsXIGboD6J1n+0Xg2Ccl9MwWwHC815CnF3aimVUEvgbaACOB64GWwDdmVilwWhKw\nzTmXe4309MB+EQnCxZ0bcEPvJrw2fT1fLd3udxwR8VGwA1hDxjm3LO8+M7sZOAG8XdB1ZlYfaOqc\nS8vnWEWgv3PuiwIu/9Y5lxg492fAkEIi3gw0A1o759YErlkErAZuBZ4p5FoRKYaHLmzL/I37+eV7\nC/msXlWSalX0O5KI+MCPlpEfCRQSVwITnHN7Czn1XuBLM+ub5/oEvDd8Xs/z6vF/FbNbaRgw82Qh\nErh+PfA9cElg10agfuDV5pOSA/tFJEgVYqJ5aUQ3AO58ax7Hs7J9TiQifvC9GAF+AlQBxhRx3kPA\nROALM+sDPypEOgNDnHOHQpCnPbAkn/1LgXYAzrldwCxgVCDHuXhjRuaG4PkiEaVxzYr87crOLN6S\nwWOfLvc7joj4oNjFiJnFmVmimUWHKMMNwE68hfgK5JzLAq7Bm2DtCzMbhFeIdAEGOecWhShPTbyx\nJXntBWrk+vk24EYzWwX8FRjhnDtlWkkzu9jMXs3IyAhRPJHyZ0j7etzcryn/mbmBTxZu9TuOiJSy\noIsRMzvAgD+zAAAgAElEQVTXzNKAw8AWvNYIzOwlMxt+Og83swbAYLxF+LKKOj9wztXA5MDWDRjs\nnFt4Os8/E8651c65VOdcK+dcF+fcrALOm+Ccu6VatWqlHVGkTLn//DZ0S6rOgx8sYu2uUDRyikhZ\nEVQxYmZD8VouMoFH8ly3DbjxNJ9/XeBeRXXR5BYDVAZygFgg1CPe9vHjFpCTCmoxEZEQiI2O4sVr\nuxEXE8WdY+dx9ITGj4hEimBbRv6E13rRH3gyz7FFQMfTfP5IYGGwLRtmFg98jNci0guYBEw0s96n\n+fz8LMUbN5JXO+CUN4FEJHQaVE/g2au6sHLHQf7wSX5Dt0SkPAq2GGkPvBn4Pu+4iH1A7eI+ODB9\nejuCbBUxswrAR3jT0g92zs3GWxNnMl5B0qu4GQrwCdDLzJrlenYy0IeSmwdFRAIGtK7LXQNb8O6c\nzbw3J+JWmhCJSMEWI4fwuinykwTsPo1n3wBkAWODPP9xoAdeIbIAIDDp2FV4qwpPKOjVXgAzu8LM\nruB/s6ReENjXP8+pr+FNYDbezC4xs2HAeGAT8EqQWUXkDNw7uBW9m9Xi4fFLWLn9oN9xRKSEWT4v\ngJx6ktnbQCvgbOAo3tiRs/BegZ0GLHfO3RT0Q735ObbizedxcZDX1AAaOecWF3C/FOfcjEKuL+iD\nTnPODchzbhLwLHDyld0pwL3OufRgsuYnJSXFzZkz53QvF4k4Ow8eY+jfv6NqQgyf3NWXyhVKfY5G\nETlDZjbXOVfkQrLBFiPNgR+AY3iv094CvIX3Rk0iXiGg9tRCqBgRKb60tbu57l8/cGGnBjx/dRe0\nFqVI2RJsMRJUN41zbi3eWI1peLOlGnAR3mDPXipERKQkpDavzS/ObcWEhVsZ+4MmOBYpr4Ju9wxM\niT4CwMyizmTlXhGRYN0xoAWz0vfxpwnL6NK4Oh0aas4ekfIm2HlGbjWz6id/ViEiIqUlKsp47qou\n1Kocxx1j55FxNLPoi0SkTAn2bZqXgG1m9r6ZDTMzjSQTkVJTs1IcL17bla37j3L/+wsJZqybiJQd\nwRYjTYHH8OYb+RivMHnBzHqUWDIRkVzOalKTBy5ow5dLd/B/36f7HUdEQijYAawbnXOPO+faAj2B\ncXgDWWeY2Uoze6gkQ4qIANzUtynntkvkic+XM2+jVmcQKS+KvWqvc262c+5uoCFwOd7aMH8KdTAR\nkbzMjKev6Ez96vHcNXYe+w6f8DuSiIRAsYsRADPrCTwHvIpXlCwNZSgRkYJUqxjLP67txu5DJ7jv\n3QXk5Gj8iEhZF3QxYmZNzexhM1sJpOG1irwBdHXOdSqpgCIieXVqVJ2HL2rL1JW7eHnaWr/jiMgZ\nCuqtGDP7Hm+V3KN4A1jvBibpFV8R8ct1vZrww/q9/O2rlZzVpAa9mtXyO5KInKZgW0aOAj8F6jnn\nrnPOfalCRET8ZGY8eXknkmtV4ufj5rPr4HG/I4nIaQr2bZrBzrkxzrlDJR1IRCRYlSvE8NJ13Thw\nNJN73p5PtsaPiJRJBRYjZlbXzKJzfV/oVnqRRUT+p029qjx6aQfS1u7h71NW+x1HRE5DYWNGtgG9\ngVnAdqCov3JEhyqUiEhxDE9pzKz1e3nh69WkNKnB2a3q+B1JRIqhsGLkDmBdru/V/ikiYevRSzqw\neHMG976zgM/u7kv9agl+RxKRIJnWeCgdKSkpbs6cOX7HECnX1uw8xLAXv6Nd/aqMu6UXsdGnNZWS\niISImc11zqUUdV6wq/Z+bmatCjjWwsw+L25AEZFQa1G3Mk9c1pE5G/bx9Jcr/Y4jIkEK9q8N5wPV\nCzhWDTgvNHFERM7MJV0aMqJnEq98u45Jy3b4HUdEglCcNsyC+nOaAIdDkEVEJCQevqgdHRpW5Zfv\nLmDT3iN+xxGRIhT2au/1ZvaVmX0V2PXCyZ9zbdPxpoRPK5W0IiJBiI+N5qVrz8IBd701j+NZ2X5H\nEpFCFNYyEgdUCWwAlXL9fHLLxitGbirBjCIixZZUqyJ/vaIzCzdn8MTnK/yOIyKFKPDVXufcv4F/\nA5jZDOCnzrnlpRVMRORMnd+hHjf1bcq/v1tP9+SaXNipvt+RRCQfwU4H31uFiIiURb85vw1dk6rz\nmw8WsW6XVrQQCUfFegnfzFqb2TAzG553K6mAIiJnIi4mihev7UZMtHHH2Hkcy9T4EZFwE+w8I1XM\n7GtgGfAR8HZgG5drExEJSw2rJ/Ds8C6s2H6QP36y1O84IpJHsC0jjwNJwBDAgKuBocAHwFqgT4mk\nExEJkYFt6nLHgOa8PXsTH8zd7HccEckl2GJkKF5BMjXw81rn3ETn3HBgOnBLCWQTEQmpX5zbih5N\na/K7j5ewasdBv+OISECwxUgDYI1zLhs4DlTOdewdYFiog4mIhFpMdBQvXtOVShWiuWPsPA4fz/I7\nkogQfDGyA6gZ+H4D0DPXseRi3EdExFd1q8bz/NVdWbvrEA99tBgtFiriv2CLiO+B7oHvxwKPmNnf\nzexvwF+BSSURTkSkJKS2qM19g1vx8YKtjJu1ye84IhGvwEnP8ngUaBT4/imgHnAtkAB8BdwV+mgi\nIiXnroEtmJ2+lz9OWEqnRtXo0LCa35FEIlawk56tdM5NCXx/3Dl3p3OujnOusnPuMufcrpKNKSIS\nWlFRxnNXdaFmxTjufGseB45l+h1JJGJprIeIRKxalSvwwrVd2bzvKL95f5HGj4j4pMBuGjO7vxj3\ncc65v4Ygj4hIqeqeXJP7z2vNE1+sYHRaOjf2aep3JJGIU9iYkSeLcR+HN5BVRKTMublfM2an7+XP\nny+nU6PqnNWkht+RRCJKYd00CcXYKpZsTBGRkhMVZTx9ZWcSq8Zz1Ssz+MU7C1i5XZOiiZSWAltG\nnHPHSzOIiIifqleM44PbU3ll2jrGzdrIh/O3MLhtXW4f0JyzmtQs+gYictqsOAO2zGwIcDZQC/iz\nc26TmfUC1jvndpRQxnIhJSXFzZkzx+8YIhKEfYdPMGZGOqPT0tl/JJMeyTW5bUAzBraui5n5HU+k\nzDCzuc65lCLPC6YYMbOqwASgH9508HFAd+fcPDN7C9jpnLv3DDOXaypGRMqeIyeyeHvWJv41fR1b\nM47Rpl4VbuvfnIs61ScmWi8jihQl2GIk2P+angJaAYOAqngr9540KbBfRKRcqRgXw0/7NmXa/QP5\n25Wdyc5x3PvOAgY8PZU3ZqRz9ES23xFFyoVgi5GfAL91zn0D5OQ5tgFICmkqEZEwEhsdxeVnNeLL\ne8/mtRtSqFulAr8fv5S+f/maF79eTcYRTZgmciaCnQ6+KlDQAg4VgOjQxBERCV9RUca57RIZ3LYu\ns9P38fLUNTz91SpenrqWa3smcVPfZtSrFu93TJEyJ9hiZDVwDjA5n2P9gKUhSyQiEubMjB5Na9Kj\naQ+WbT3AK9+u5d/frWd0WjqXdW3ELf2b0bxOZb9jipQZwQ5gvRP4G/AQ8BawBe+tmsbAq8DPnXOj\nSy5m2acBrCLl26a9R3ht+jremb2JE9k5nN++Hrf1b07nxtX9jibim5C+TRO44XPAz0/+iDfrqgP+\n7pz75ekGjRQqRkQiw+5Dxxn9fTpvzEjnwLEsUpvX4o4BLejTopZeC5aIE/JiJHDTlsB5QF1gD/CV\nc275aaeMICpGRCLLwWOZjJu1kX9/t54dB47ToWFVbu/fgvM71CM6SkWJRIaQFSNmFgc8ArzvnJsb\nonwRR8WISGQ6npXNx/O38Mq0dazbfZjkWhW5tX9zLuvWkAoxGvsv5VuoJz07ApzvnPs2FOEikYoR\nkciWneP4aul2Xp62lkWbM6hbpQI/7duUET2TqBIf63c8kRIR6knPFgLtziySiEjkio4yLuhYn/F3\n9mHsz3rSul4VnvxiBalPfs1TE1ew66CWA5PIFWzLSD/gDeBm51x+r/dKEdQyIiJ5Ld6cwT+nreXz\nJduIjY5ieEojbunXnKRaWghdyodQd9OsBmrjTX52BNiO9ybNSc451/o0s0YEFSMiUpD1uw/z6rdr\n+WDuFrJycriwUwNu69+M9g2q+R1N5IyEuhh5mx8XH6dwzl0TfLzIo2JERIqy88Ax/v39esbO3Mih\n41n0b1WH2wc0p2fTmnotWMqkEnm1V06fihERCVbG0UzenLmB179fz+5DJ+iaVJ3b+zdncNtEovRa\nsJQhKkbCjIoRESmuY5nZvDd3M69+u5ZNe4/Som5lbj27GZd0aUhcTLDvH4j4R8VImFExIiKnKys7\nh88Wb+Of09axfNsB6leL52f9mnF198ZUqhDsEmMipU/FSJhRMSIiZ8o5x7RVu3h56lp+WL+X6hVj\nGdk7mZGpydSsFOd3PJFTqBgJMypGRCSU5m7Yxz+nrWXSsh0kxEZzVffG3Hx2MxpWT/A7msh/qRgJ\nMypGRKQkrN5xkFe+XcfH87cAMKxLA27r35xWiVV8TiaiYiTsqBgRkZK0df9R/jV9PeNmbeRoZjaD\n2yZy+4DmnNWkht/RJIKFvBgxs0TgHuBsoCZwhXNumZndAcxyzuk3bSFUjIhIadh3+ARvzNjA6LT1\n7DuSSY+mNbm9f3MGtK6juUqk1IV0bRozawMsBm7Hm4G1NRAfONwauPc0c4qISAjVqBTHPYNb8v0D\n5/CHi9uxee8Rbhw9mwv+Pp3xC7aQlZ3jd0SRUwT7ovrTwHqgKTAUyF1efw/0DnEuERE5AxXjYrix\nT1Om3T+Qv13Zmewcxz1vL+Dyl9PIVEEiYSbYYqQ/8Gfn3H5OnRZ+O1A/pKlERCQkYqOjuPysRnx5\n79k8/pMOLNycwduzNvodS+RHijOFX3YB+2sBR0OQRURESkhUlHFtjyR6N6vFs5NXk3E00+9IIv8V\nbDEyB7i+gGOXAzNDE0dEREqKmfHQhW3Zd+QE//hmjd9xRP4r2GLkceByM5sAXInXVXO2mb0CDAf+\nXEL5REQkhDo0rMYV3Rox+vt0Nuw57HccESDIYsQ5Nxmv6OgMvIU3gPUZ4EJguHPu+xJLKCIiIfWr\n81oTE208+cUKv6OIAMUYM+Kc+xBoAnQCBgNdgSTn3McllE1EREpAYtV4buvfnC+WbGfW+r1+xxEp\n1gBWnGeJc+5r59xC55zeDxMRKYNu7teM+tXieeyzZeTkaCZu8VdQa0+b2fCiznHOvXvmcUREpDQk\nxEVz//mtue+dhXy8YAuXdWvkdySJYEEVI8DbBezPXU6rGBERKUMu6dyQ179P56mJK7mgQ30S4qL9\njiQRKthumrb5bH2BvwBrA9+LiEgZEhVl/O7Cdmw/cIzXpq/zO45EsKBaRpxzKws4lGZm2Xhr1swI\nWSoRESkVPZrWZGjHerw8dS1XdW9MYtX4oi8SCbFiDWAtwDfAsBDcR0REfPDA+W3JznE8/WVBf+8U\nKVmhKEZS8FbyFRGRMiipVkVu7JPM+/M2s2RLht9xJAIF+zbN/fnsjgM6AD8BXgtlKBERKV13DGzB\ne3M389hnyxh3cy/MrOiLREIk2LdpnsxnXzawBXgWeCRkiUREpNRVS4jlvsEteXj8UiYt28GQ9vX8\njiQRJNhumoR8tgrOuWTn3APOOa3aKyJSxl3TI4kWdSvz58+XcyJLc1pK6SmyGDGzOOCPQAfn3PFc\nm/6kioiUIzHRUTx0YVvS9xzhPzM3+B1HIkiRxYhz7gRwD1Cp5OOIiIifBrSqQ7+WtXl+ymr2Hznh\ndxyJEMF20ywE2pVkEBER8Z+ZNxHawWOZPDd5td9xJEIEW4zcD/zGzAaXZBgREfFf63pVuLpHEm/O\n3MDaXYf8jiMRINhi5P+A6sCXZnbQzFab2apcm2bKEREpR35xbiviY6N54vMVfkeRCBDsq71z+fGi\neCIiUo7VrlyBOwe24C8TV5C2ZjepLWr7HUnKsWDXprm6pIOIiEh4ubFPMmN/2MCjny3n05/3JTpK\nE6FJySiwm8bM1plZ59IMIyIi4SM+NpoHLmjD8m0HeH/uJr/jSDlW2JiRZKBCKeUQEZEwdGHH+nRL\nqs7TX63i0PEsv+NIORWKhfJERKScMjMevqgduw4e55Vpa/2OI+VUUcWIBq2KiES4rkk1uKRLA179\ndh1b9mv1Dwm9ogawPmJmu4O4j3POjQxFIBERCT/3n9+GiUu289eJK3ju6q5+x5FypqhipAtwPIj7\nqAVFRKQca1g9gZ/1a8o/vlnLqD5N6dK4ut+RpBwpqpvmUudc0yC2ZqWSVkREfHP7gBbUrlyBxz5d\nhnP6O6iEjgawiohIUCpXiOFXQ1oxZ8M+Pl+83e84Uo6oGBERkaBdmdKYNvWq8MQXyzmWme13HCkn\nVIyIiEjQoqO8VX037zvK6LR0v+NIOVFgMeKci3LOzSrNMCIiEv76tqzNoDZ1+cfXa9h9KJh3HEQK\np5YREREptt9e2Jajmdk8O2mV31GkHFAxIiIixda8TmWu69WEcbM2smrHQb/jSBmnYkRERE7LPYNa\nUrlCDI9/ttzvKFLGqRgREZHTUqNSHHcPasm0VbuYunKn33GkDFMxIiIip+2G3skk16rI458tJys7\nx+84UkapGBERkdMWFxPFg0PbsnrnIcbN3uR3HCmjVIyIiMgZGdIukZ5Na/LspFUcOJbpdxwpg1SM\niIjIGTEzHr6oHfuOnOAf36zxO46UQSpGRETkjHVoWI3LuzXi9e/S2bjniN9xpIxRMSIiIiHxqyGt\niY4y/jJxhd9RpIxRMSIiIiFRr1o8t/ZvxmeLtzEnfa/fcaQMUTEiIiIhc8vZzahXNZ5HP11GTo7z\nO46UESpGREQkZCrGxfDr81qzcHMGnyzc6nccKSNUjIiISEj9pGtDOjasxl8mruDoiWy/40gZoGJE\nRERCKirK+N2FbdmWcYx/TV/ndxwpA1SMiIhIyPVsVovz29fj5Wlr2XngmN9xJMypGBERkRLx4NA2\nZGbn8PRXK/2OImFOxYiIiJSIJrUqMSo1mffmbmbp1gy/40gYUzEiIiIl5q5zWlI9IZbHP1uOc3rV\nV/KnYkREREpMtYRY7ju3FWlr9zB5+U6/40iYUjEiIiIl6poeSTSvU4k/f76cE1k5fseRMKRiRERE\nSlRsdBQPXdiW9bsPM/aHDX7HkTCkYkRERErcwNZ16deyNs9NXs3+Iyf8jlMsOw8c4/73F3LxC9+R\ncTTT7zjlkooREREpcWbGQxe25eCxTJ6fssbvOEE5lpnNi1+vZsDTU/lo/haWbs3gGb2mXCJUjIiI\nSKloU68qV3VP4o0Z6azbdcjvOAVyzjF+wRYG/W0aT3+1in4tazPpvv5c36sJ/5m5gSVb9JpyqKkY\nERGRUvOLc1tRISaKJ75Y4XeUfM3buI/LXk7jnrcXUL1iLONu7sUr16eQXLsSvxjSmpqVKvDQx0u0\nInGIqRgREZFSU6dKBe4Y2IJJy3aQtna333H+a8v+o9zz9nwueymNzfuO8tQVnfjkrr70bl7rv+dU\nS4jldxe2ZeGm/bw9e5OPacsfFSMiIlKqburblIbVE3js0+Vk+9zCcPh4Fn/7aiXnPD2ViUu2c9fA\nFnzzqwEMT2lMdJSdcv4lXRrQq1lN/jJxBXsOHfchcfmkYkREREpVfGw0v7mgDcu2HeCDeZt9yZCT\n43hvziYGPj2VF75ew5D29Zjyy/786rzWVK4QU+B1Zsajl3Tg8PEsngzTrqaySMWIiIiUuos71adr\nUnWe/nIlh49nleqzf1i3h2H/+I5fv7+IBtUT+OD2VF64piuNalQM6vqWiVX4Wb9mvDd3M3PS95Zw\n2sigYkREREqdmfHwRe3YefA4r0xbWyrP3LjnCLe/OZerXp3JnkMneO6qLnx4eypnNalR7HvdPagF\nDarF87uPl5CVrVllz5SKERER8UW3pBpc3LkBr05fx9b9R0vsOQeOZfLE58sZ/Mw0pq7cxS/ObcXX\nvxzApV0bEpXPuJBgVIyL4fcXt2fF9oOMTksPbeAIpGJERER885vzW5Pj4K9fhn4ysazsHMb+sIGB\nf53KK9+uY1iXBkz99QDuHtSShLjoM77/ee0TGdi6Ds9OWsX2jGMhSBy5VIyIiIhvGtWoyM/6NuWj\n+VtYuGl/yO47ffUuLnz+Ox76aAnN61Rmwl19efrKziRWjQ/ZM8yMPw5rT2aO47HPloXsvpFIxYiI\niPjq9gHNqV05jsc+W4ZzZ/aq79pdh7hp9Gyu//csjmRm8fKIbrxzay86NqoWorQ/1qRWJe4c0IJP\nF23ju9XhM29KWaNiREREfFUlPpZfDmnN7PR9fLFk+2ndY/+REzwyYSnnPfstP6zfywMXtGHSff25\noGN9zE5vXEiwbu3fjORaFfn9+CUcz8ou0WeVVypGRETEd8NTGtOmXhWe+GJ5sX6hZ2bn8Pr36+n/\n16mMSUvnypTGfPOrAdzWvznxsWc+LiQY8bHR/OmSDqzbfZjXvl1XKs8sb1SMiIiI76KjvFV9N+09\nypgg3k5xzjFl+Q7Oe+5bHpmwjI4Nq/H5Pf144rKO1KlSoeQD53F2qzoM7ViPF75ew6a9R0r9+WWd\nihEREQkL/VrW4Zw2dXlhyppCp1pfuf0gN/zfLG4aMwcc/HtkCv+5qQdt6lUtxbSneviidkRHGY9M\nWOprjrJIxYiIiISN3w5tw5HMbJ6bvPqUY7sPHee3Hy3mgr9/y6LNGfz+onZ8ed/ZDGqbWOLjQoJR\nv1oC9w1uxeTlO5m0bIffccqUgifgFxERKWUt6lZhRM8kxv6wkRt6N6FlYhWOZ2Uz+vt0Xvx6DUcz\ns7mhdzL3Dm5J9Ypxfsc9xag+ybw3dxN//GQpfVrUomKcfs0GQy0jp8nMmpvZd2a2yszmm1mK35lE\nRMqDewe3omJcNI99tpyJS7Zx7jPf8sQXK+jetCYT7z2bPw5rH5aFCEBsdBSPXdqRLfuP8uLXa/yO\nU2aoGDl9/wTGOOdaAfcDYy0c2glFRMq4mpXiuPuclkxbtYvb3pxHQmw0/7mpB/83qjst6lb2O16R\nejStyeXdGvHa9HWs2XnI7zhlgm/FiJkNNbNvzeyQmR0wszlmdk4JPauRmb1gZjPM7IiZOTNLLuT8\nxmb2vpllBLJ9aGZJuY7XAXoBowGcc5MAA84qifwiIpHmhtQmDE9pxGOXduCzu/vSr2UdvyMVy4ND\n25AQG83vxy8544ncIoEvxYiZ3QqMB+YCPwGuBN4Dglu/ufhaAMOBfcD0IrJVBL4G2gAjgeuBlsA3\nZlYpcFoSsM05l5nr0vTAfhEROUMVYqJ56orOXNerCTHRZa8Rv3blCvz6/Dakrd3DJwu3+h0n7JX6\nyJpAi8RzwK+dc8/lOvRlEdfVB5o659LyOVYR6O+c+6KAy791ziUGzv0ZMKSQR90MNANaO+fWBK5Z\nBKwGbgWeKSyniIgIwLU9knhvziYe+2w5A9vUpWp8rN+RwpYf5eZPgRy8MRfFcS/wpZn1zb3TzBKA\nCcDrZpZvZ6JzLqcYzxkGzDxZiASuXw98D1wS2LURqG9muf9kJQf2i4iIEB1lPHZpB3YfOs6zk1b5\nHSes+VGM9AVWAFeb2VozyzKzNWZ2ZxHXPQRMBL4wsz7wo0KkMzDEOReKkULtgSX57F8KtANwzu0C\nZgGjAjnOxRszMjcEzxcRkXKiU6PqjOiZxJi0dJZuzfA7TtjyoxhpgDcG46/Ak3hdJpOAF83snoIu\ncs5lAdfgded8YWaD8AqRLsAg59yiEOWriTe2JK+9QI1cP98G3Ghmq/A+ywiXzyglM7vYzF7NyNAf\nQhGRSPTrIW2oUTGOhz9eQk6OBrPmx49iJAqoAtzqnHvNOfe1c+52vFaPBwu7MFCQXA1MDmzdgMHO\nuYUlnDm/LKudc6nOuVbOuS7OuVkFnDfBOXdLtWols3y1iIiEt2oVY/nt0LbM27if9+Zu8jtOWPKj\nGNkT+Dopz/6vgMTAQNXCxACV8cadxBL6N3D28eMWkJMKajEREREp1GXdGtIjuSZPfrGCfYdP+B0n\n7PhRjJz2CkJmFg98jNci0guvoJloZr1DlA28fO3z2d8OWBbC54iISIQwMx69tAMHjmXx1Jcr/I4T\ndvwoRj4KfD0vz/7zgc3OuW35XWRmFQLXdsfrmpkNXIXXXTPRzHqFKN8nQC8za5br2clAn8AxERGR\nYmtdrwo39W3KuFmbmLdRDe25+VGMfA58A7xiZreZ2RAzew1vIOvDhVz3ONADrxBZABCYdOwqYAow\noaBXewHM7Aozu4L/zZJ6QWBf/zynvoY3gdl4M7vEzIbhTdC2CXilmJ9VRETkv+4Z1JJ6VeP53UdL\nyMouzqwT5Zv5MU2tmVUFngCuwBufsQJ40jn3ViHX1AAaOecW53MsFkhxzs0o5PqCPug059yAPOcm\nAc8CJ1/ZnQLc65xLL+RjFSolJcXNmTPndC8XEZFy4vPF27hj7Dz+cHE7buzT1O84JcrM5jrnilxI\n1pdiJBKpGBEREQDnHCNfn838DfuY8sv+1K0a73ekEhNsMVL2JvwXEREpw8yMPw1rz/HsHB7/fLnf\nccKCihEREZFSlly7Erf1b874BVtJW7Pb7zi+UzEiIiLigzsGNCepZkUeHr+EE1mRPZhVxYiIiIgP\n4mOjeeSS9qzddZh/fbfO7zi+UjEiIiLik4Gt63Je+0Sen7KazfuO+B3HNypGREREfPT7i9tjGH+a\nELmTfKsYERER8VHD6gncM7glXy3bwZTlO/yO4wsVIyIiIj77aZ+mtKxbmT9OWPr/7d19tFTVecfx\n7w+4KKAIGOI73iCK9f11LTG+VG19wdS32kSbqugyNSZZS2qt1aUxVO0yrbZpGlOTaBu1GhNfMcYG\nAQ0YY4ga31FBVFJNQVEEVG7ggk//2PuayTAz94Iz98zc+/usddYw++x99nP2HO48c84+M/yuc23R\n4fQ6JyNmZmYFGzxoAJcfvxuvL+3gP362oOhwep2TETMzsyYwYYfNOXHvbfjO7Fd5dcn7RYfTq5yM\nmGLrFLkAAA53SURBVJmZNYmLJ+7MRm0DuOzeufSnn2txMmJmZtYkPrnpxlxw5HgeWfA29z+3qOhw\neo2TETMzsybyVwdsz65bD+eKn7zA+6vWFB1Or3AyYmZm1kQGDhBXnrAbb723in+bMb9h/axas5a3\nVvyO+W++x2OvLWX63MV8UFDyM6iQXs3MzKyqvceM5JT9x/D9Rxdy8n7bsvOWwyvWiwhWrl7Lso5O\nlq1czfKVnfnfnSzrWM3yjs5Ulp8vW9nJ8ry+o8ItxNMmH1y1r0ZyMmJmZtaELjxqPA/MXczf/OgZ\nJozdPCUXXclETiiWd6ymc231ia6DBw1gxJA2RgxtY8SQwWw3aii7dz0fOpjNStaNGNpG++bDenEP\nf8/JiJmZWRMaOWwwl31mF86//WleX7ry94nD0DbGD9+U4R8lEulxs5xQlCYXG7cNLHo3esTJiJmZ\nWZM6Ye9tOG7PrRkwQEWH0lCewGpmZtbE+noiAk5GzMzMrGBORszMzKxQTkbMzMysUE5GzMzMrFBO\nRszMzKxQTkbMzMysUE5GzMzMrFBORszMzKxQTkbMzMysUE5GzMzMrFBORszMzKxQTkbMzMysUIqI\nomPoFyQtB16u82Y/Abxd521a37AZsLzoIFpUXx+7Vtu/Zoq3qFh6q99G9LN9RIzurtKgOndq1f0o\nIv66nhuU9ERE7FfPbVrfIOl79T7e+ou+Pnattn/NFG9RsfRWv0WOtS/T9J77ig7A+hUfbxuur49d\nq+1fM8VbVCy91W9hY+3LNC3MZ0bMzKwv8JmR1va9ogMwMzP7uHxmxMzMzArlMyNmtl4k7SDpEUnz\nJT0lyZcKe6ivj11f379G6u9j52TEzNbXd4CbImIn4ELgVkkqOKZW0dfHrq/vXyP167FzMtJH9fcs\nuxlJOlnSVEmvS+qQNE/SVZI2bWCf20r6lqRfSlopKSS116i/naQ7JS2XtELS3ZLGlKwfDRwA3AgQ\nETMAAfs2ah9yv0dJekjSYkmrJL0h6XZJuzSwz8LGTtK03N+Vdd6t0j76xLFR0v9ESQ9Lej/H94Sk\nwxvUV58au2bgZKTv6tdZdpO6AFgLXAwcA1wHnAvMkNSo/4vjgM8C7wI/r1VR0lDgIWBn4AzgNGBH\n4GeShuVqY4BFEdFZ0nRhLm+kUcCvga8AR5LGcFdgjqTtG9RnIWMn6VRgzzrE352+cmwg6RzgXtIx\nciLwF8AdwNAGddlnxq5pRISXJliAbYFvAb8EVgIBtFepux1wJ+mb8lYAdwNjStaPBt4D2krK5gP7\nFb2f/XkBRlcoOz2/1ofXaLcVcGCVdUOBY2q0HVDy77O7Oa7OIyVL40rKPgWsAc7Pz/cF5pe1mw6c\nVMB4js/787d9ZeyAkcBi4NTc35XdjEFL7V+DjoN2oAOYvJ7t+v3YNdPiMyPNo0eZtrPs1hURSyoU\nP54ft6nRdDLwgKSDSgslDSF9SdH3JW1Spc8P1yPE44A5EbGgpP1rwC+A43PR/wJbSWoradeey3vb\nO/lxTY06rTZ2/wQ8HxG39bDfVtu/RjgL+JB0Nnh9eOyaiJOR5vFwRGwRERNJpxer+QIwFjghIqZG\nxL2kg3174JxeiNPq69D8+GKNOpcA04CfSvo0/MEfzD2BIyPi/TrEsivwfIXyucAu8FFC9RgwKcfx\np6Rr27+uQ//dkjRQ0mBJOwLfJZ1FqPXG3TJjl98UTwe+vB79tsz+NdBBwEvAKZJekbRG0gJJ3Y2j\nx66JOBlpEuuRaTvL7iMkbQNcDsyMiCeq1YuINaTT9g+Q/nAeQfqDuRdwREQ8W6eQRpHOzJVbSrp8\n0OWLwJmS5gNXA5+PfF65F/wKWEW67LgH6fLWW9Uqt8rYSRpMSq6uiYh5Pe20VfavTjFUszXp7PDV\nwNdJc4pmANdKOq9aI49dc/EP5bWeXUkTtcrNJU3aIiKWSOrKsq/vj1l2s8unf+8lXWI4s7v6EbFG\n0inA7cBM0h+3wyPimYYGWjmWl4EDe7vf7DRgOOns4AWkyb8HRcTCag1aZOwuBIYA/7gB22yF/Wuk\nAcCmwKSIuDuXPZTvbrkY+Ga1hh675uEzI63HWXaLKzkVPBY4KiLe6GHTQcAmpOvjbdT/ToF3+cNj\nqEu1Y67XRcSLEfGrPKfiCNJ4XNSDpk07dvkWz0uArwIbSRohaURe3fV8YDf9N+3+9YKuuUMzysqn\nA1tI2qqb9v157JqGk5E+KiJejogDI2KniNgrIh4rOiaDfOnsTmA/YGJEPNfDdhsDU4F9SN9HMAOY\nJmlCHcObSzrzVm4X4IU69lMXEbEMWECa/F1VC4zdWGBj4BbSm1PXAunsz7vA7tUat8D+NdrcDW3o\nsWseTkZaj7PsFpW/S+RW4HDSBOQ5PWy3EXAPsD/wJxHxOPA50mnlaZIOqFOIPwYOkDS2pO924NN5\nXVORtAXprrJXatRphbF7GjiswgIpQTmMlHSto0X2r9HuyY9HlZUfDbwREYsqNfLYNZmi7y32su5C\njfvWSbf1PlKhfBYwu+jYvdR8Xa/Lr+uVpE9hpcu2NdpdQzoVvXdZeRvpO2aWAJvUaH9yXrr6Pzc/\nP7Ss3jDSm95zpMnQxwHPAK/W2n4vjd09pMsYx5PenM8h3UGxDNipL45d17HSTZ2W3b86HhvKfxff\nIV2ePhK4Pu/PJI9dayyFB+ClwotSOxmZTJr0OLakrB3opMaXP3kpfiF910tUWabUaDcS2L3KujZg\nQjf9VutzVoW6Y4C7SF+m9x7pFPY6x2EBY/f3pAnYy0hfCjiPdPdJzdhaeezoWTLSsvtX5+NjOPBt\n4E1gNfAs8Jceu9ZZlAfKmoCkk/M/jyBl+F8iZedLImJ2rjOMlFV3AJeSDv4rSLPJ94j63BdvZmbW\na5yMNBFJ1V6M2RHxxyX1xgDfALpu2X2Q9FXICxsdo5mZWb05GTEzM7NC+W4aMzMzK5STETMzMyuU\nkxEzMzMrlJMRMzMzK5STETMzMyuUkxEzMzMrlJMRMzMzK5STEbM+SNIkSSFpXEnZZEknFRjTCElT\nJO1TYd0sSbMKCKup5NfsyqLjMOttg4oOwMx6zWTgEdKPgBVhBPA14A3gybJ1X+r9cMysWTgZMbMN\nJmmjiFj1cbcTES/UIx6rTZKAtohYXXQsZqV8mcasH5C0ENge+Hy+FBCSbixZv6ekH0t6V1KHpF9I\nOrhsGzdKekPSBEmPSuoA/jmvO0XSQ5KWSHpf0lOSzihp2w68lp9eXxLDpLx+ncs0ksZLukfSshzT\nHElHl9WZkrezo6T7c9+/kXSZpJp/3yS157bnSLpc0qLc132Sti2rG5KmVGk/qcIY7dc1RpLmSTo2\nrz9f0kJJKyTdK2l05dB0Sd5Oh6SHJe1VodJJeUxW5rjvyL9bVVpnoaRbJJ0l6SXSL9oeW2tczIrg\nZMSsfzgRWAw8AEzIyxUAeQ7Ho8Ao4AvAnwPvADMl7Vu2nc2AHwK3AccAP8jlO5B+Fv004ATgPuAG\nSV/M6xcBXfNVriqJ4f5KwUramnRJaU/gK8BngWXA/ZKOqdDkHuCh3PdU4B+AMyrUq+RiYBxwFnBe\njuuWHratZDhwM3ADadzfAu6S9C/AYcCXSZfMDiP97H2504GJpP2eBGwBPChpVFeFPK53AS8AJwPn\nALsBsyVtWra9w4DzSWNyNPDsx9g3s8aICC9evPSxhfQmFsC4krKFwC0V6j4IvAgMLikbmMumlpTd\nmLd5fDd9DyBdAr4eeKakvD23P7tCm1nArJLn1wBryuIfCMwDniwpm5K3eWbZ9p4DpncTZ1c8s8rK\nL8jlW5eUBTClSvtJFcbokJKyPXLZPGBgSfm/Ap1lZQG8DQwr66cTuCI/3wRYDvxXWTyfIp35mFz2\nmq8Etiz6mPTipdbiMyNm/ZikIcChwB3Ah5IGSRoECJgJHFLWpBP4SYXt7CjpNkm/zXU6gbOB8RsY\n2iHAnIhY0FUQEWtJZ2T2kjS8rH75GZbngTH0zP+UPX8uP/a0fbkPIuLhkucv5ceZeR9KywcBW5XH\nExEfdD2JiIXAHNIZG/LjcODWrtcrv2av522Wv2ZzImLxBu6LWa/wBFaz/m0U6YzDV/OyDkkDIuLD\n/HRJ2RsqkjYBZpA+gV8EvEL6hH4u6dLHhsb1VIXyxaREaSSwoqR8aVm9VcDGPeyrUlvWo325ZaVP\nImJ1mjfKu2X1uiaRlvfzZoVtvgnsmv/9yfw4s0r/5f0sqhqpWZNwMmLWvy0DPiTNXbi5UoWSRATS\nZYRyE0iTYw+OiEe6CvOn9Q21FNiyQvmWOYbyN9xGWwUMLivbvEF9bVGl7Lf53+/kx0nA3Ap13yt7\nXuk1M2sqTkbM+o9VwJDSgoj4QNLPSRNFnyxLPHpqaH7s7CqQNBI4vkL/lMdQxWxgsqT2fJkCSQOB\nzwFPRcSKWo0b4DekCaKlGnVXykRJw7ou1eQ7kQ4Avp7XP0pKOMZFxE0NisGsVzkZMes/XgAOlvQZ\n0uWOt/Mb/fnAw8ADkv6TdFr/E8A+pMmVF3Wz3UdJl0y+LelrwDDgUtJEzM1K6r1J+lR/iqRngQ+A\n1yLiHdb1DdIn/xl5mytIX4y2E8XcmvpD4FJJl5DmbxwMnNqgvjqA6ZKuBjYi3QWzgjQmRMQKSX9H\nGu/RwE9JE1q3Ic3/mRURP6i4ZbMm5QmsZv3HxaQ7Om4HHifdiUJEPAnsT0oU/h2YDnwT2J2UpNQU\nEUtIt7AOBO4k3bp7A2W3x+azLmeT5nvMzDH8WZVt/h9wEOkyxHV5u6OAYyNiWo/3uH6uAq4l3W47\nFfgj0m3MjXAzaULutcBNwBLgiIj4aG5LRHwXOI40Qfi/SZNwp5A+YD7doLjMGkYRvpxoZmZmxfGZ\nETMzMyuUkxEzMzMrlJMRMzMzK5STETMzMyuUkxEzMzMrlJMRMzMzK5STETMzMyuUkxEzMzMrlJMR\nMzMzK9T/A8Scaw5ddXhCAAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "#NBVAL_SKIP\n", "import matplotlib.pyplot as plt\n", "\n", "# Plot objective function decrease\n", "plt.figure()\n", "plt.loglog(relative_error)\n", "plt.xlabel('Iteration number')\n", "plt.ylabel('True relative error')\n", "plt.title('Convergence')\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This notebook is part of the tutorial \"Optimised Symbolic Finite Difference Computation with Devito\" presented at the IntelĀ® HPC Developer Conference 2017." ] } ], "metadata": { "anaconda-cloud": {}, "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.4.3" }, "widgets": { "state": {}, "version": "1.1.2" } }, "nbformat": 4, "nbformat_minor": 1 }