{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "
\n", " \n", " \"QuantEcon\"\n", " \n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# General Purpose Packages" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Contents\n", "\n", "- [General Purpose Packages](#General-Purpose-Packages) \n", " - [Overview](#Overview) \n", " - [Numerical Integration](#Numerical-Integration) \n", " - [Interpolation](#Interpolation) \n", " - [Linear Algebra](#Linear-Algebra) \n", " - [General Tools](#General-Tools) " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Overview\n", "\n", "Julia has both a large number of useful, well written libraries and many incomplete poorly maintained proofs of concept\n", "\n", "A major advantage of Julia libraries is that, because Julia itself is sufficiently fast, there is less need to mix in low level languages like C and Fortran\n", "\n", "As a result, most Julia libraries are written exclusively in Julia\n", "\n", "Not only does this make the libraries more portable, it makes them much easier to dive into, read, learn from and modify\n", "\n", "In this lecture we introduce a few of the Julia libraries that we’ve found particularly useful for quantitative work in economics\n", "\n", "Also see [data and statistical packages](https://lectures.quantecon.org/data_statistical_packages.html) and [optimization, solver, and related packages](https://lectures.quantecon.org/optimization_solver_packages.html) for more domain specific packages" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Setup" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "hide-output": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\u001b[32m\u001b[1mActivated\u001b[0m /home/qebuild/repos/lecture-source-jl/_build/website/jupyter/Project.toml\u001b[39m\n", "\u001b[36m\u001b[1mInfo\u001b[0m quantecon-notebooks-julia 0.1.0 activated, 0.2.0 requested\u001b[39m\n" ] } ], "source": [ "using InstantiateFromURL\n", "github_project(\"QuantEcon/quantecon-notebooks-julia\", version = \"0.2.0\")" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "hide-output": true }, "outputs": [], "source": [ "using LinearAlgebra, Statistics\n", "using QuantEcon, QuadGK, FastGaussQuadrature, Distributions, Expectations\n", "using Interpolations, Plots, LaTeXStrings, ProgressMeter" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Numerical Integration\n", "\n", "Many applications require directly calculating a numerical derivative and calculating expectations" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Adaptive Quadrature\n", "\n", "A high accuracy solution for calculating numerical integrals is [QuadGK](https://github.com/JuliaMath/QuadGK.jl)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "hide-output": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(value, tol) = quadgk(cos, -2π, 2π) = (-1.5474478810961125e-14, 5.7846097329025695e-24)\n" ] } ], "source": [ "using QuadGK\n", "@show value, tol = quadgk(cos, -2π, 2π);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This is an adaptive Gauss-Kronrod integration technique that’s relatively accurate for smooth functions\n", "\n", "However, its adaptive implementation makes it slow and not well suited to inner loops" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Gaussian Quadrature\n", "\n", "Alternatively, many integrals can be done efficiently with (non-adaptive) [Gaussian quadrature](https://en.wikipedia.org/wiki/Gaussian_quadrature)\n", "\n", "For example, using [FastGaussQuadrature.jl](https://github.com/ajt60gaibb/FastGaussQuadrature.jl)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "hide-output": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "w ⋅ f.(x) = 0.6666666666666667\n" ] } ], "source": [ "using FastGaussQuadrature\n", "x, w = gausslegendre( 100_000 ); # i.e. find 100,000 nodes\n", "\n", "# integrates f(x) = x^2 from -1 to 1\n", "f(x) = x^2\n", "@show w ⋅ f.(x); # calculate integral" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The only problem with the `FastGaussQuadrature` package is that you will need to deal with affine transformations to the non-default domains yourself\n", "\n", "Alternatively, `QuantEcon.jl` has routines for Gaussian quadrature that translate the domains" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "hide-output": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "w ⋅ cos.(x) = -3.0064051806277455e-15\n" ] } ], "source": [ "using QuantEcon\n", "\n", "x, w = qnwlege(65, -2π, 2π);\n", "@show w ⋅ cos.(x); # i.e. on [-2π, 2π] domain" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Expectations\n", "\n", "If the calculations of the numerical integral is simply for calculating mathematical expectations of a particular distribution, then [Expectations.jl](https://github.com/QuantEcon/Expectations.jl) provides a convenient interface\n", "\n", "Under the hood, it is finding the appropriate Gaussian quadrature scheme for the distribution using `FastGaussQuadrature`" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "hide-output": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "E(f) = 2.3614304974639855e-17\n" ] }, { "data": { "text/plain": [ "true" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "using Distributions, Expectations\n", "dist = Normal()\n", "E = expectation(dist)\n", "f(x) = x\n", "@show E(f) #i.e. identity\n", "\n", "# Or using as a linear operator\n", "f(x) = x^2\n", "x = nodes(E)\n", "w = weights(E)\n", "E * f.(x) == f.(x) ⋅ w" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Interpolation\n", "\n", "In economics we often wish to interpolate discrete data (i.e., build continuous functions that join discrete sequences of points)\n", "\n", "The package we usually turn to for this purpose is [Interpolations.jl](https://github.com/JuliaMath/Interpolations.jl)\n", "\n", "There are a variety of options, but we will only demonstrate the convenience notation" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Univariate with a Regular Grid\n", "\n", "Let’s start with the univariate case\n", "\n", "We begin by creating some data points, using a sine function" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "hide-output": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdZ2BUVfow8OfcmWTSe2Yy6QklEDqBNBIIXToowqqI7CrWrFhYO5ZVUezu6rrrq9j/igWQLr2mUELovaZMMpM+qVPueT+MxpAESJmZe+69z+8T3CRzz8yce57TD6GUAkIIISRXnNAJQAghhISEgRAhhJCsYSBECCEkaxgIEUIIyRoGQoQQQrKGgRAhhJCsYSBECCEkaxgIEUIIyRoGQoQQQrKGgRAhhJCssRIIjxw5InQSmGC1Wq1Wq9CpYIvZbBY6CWyxWq08zwudCoZQSi0Wi9CpYIvFYsHtMzuOlUA4ePBg/NoAoL6+vr6+XuhUMITn+YqKCqFTwRaj0djU1CR0KhhiNpurq6uFTgVbqqqqsHLQcawEQoQQQkgQGAgRQgjJGgZChBBCsoaBECGEkKxhIEQIISRrdg6EVqu1T58+La9UVlZOmzYtICBg+vTplZWV9r0dQggh1E32DIQffvhhamrqmTNnWl5ctmxZVFSUTqeLjIx866237Hg7hBBCqPvsGQgHDhy4ZMmSVhdXrVqVmZmpUqkyMzNXrlxpx9shhBBC3UfsvoydkGte08vLy2AwuLu7NzQ0aDSampqa6/3VP/7xD0JIq+sDBgyYMWOGfVPIstraWgDw8vJq+6NKEzleRU7VcN5KOi6ED3aTxf4DPM+XlZWp1WqhE8KQ6upqV1dXd3d3oRNyI3UWUlAPpQ2kpJFYKdwSyge4OirHmkym2tragIAAB72+GJWXl/v4+Li4uAidEOF5eHgoFIob/47S0YmglNrCG6X0xpuHubq6tg2ESqWy7UUJs73Zlm/ZSuF/ZxUfnFbWW2i8H433pWVN3FN5yp7e9JYw+kic1Vsp5YhI/iB0QhjC+GdSb4GPzig/PsOp3UDtRrUe0GSFp/OU6Wp6Z4x1Yijvau8peox/IILAz6RZRz4EhwfC0NDQgoKCXr16FRUVhYWF3eA3X331VfzabI3p5hbh0Qq6cI/VQwmbJyv6+//54Zh52FtKvz7Hj9qiXDFGMSRQsp8bz/MNDQ3tNpFly2KxqFQqBluEPIVvz/MvHORTNOTgLC7G+89saTTDz5f4T84q3j4Fv05QRHjaM8eaTCZKKWaSlhobG728vLBF2EEOXD6xc+dOAJg2bdry5csppcuXL5dVJ2f3vZxnHb/Rcn8fbvsUZcsoCAAuHIzWki9GKl5L4G7ZZPnkFG7BjARm5uH2bdb/nOJXjFWsGKNoGQUBwNsF/tqb2z1VeUcPLvlXa7Zeyt0YSHQcOEZo+0dVVdVdd9115MiRoUOHfvPNN76+vtf7K57n5dwizMvLe/bZZ7OysgAgNTU1av7rhz0Hr52oDLlZvf9cNZ273TowgHwxSiG9j4/neYPBoNFohE4IQ6qqqlhrEVp4uHOHtdEKP49T3LTnc30B/esuy7vJirt72qcibjKZjEZjYGCgXV5NGsrKynx9fbFF2EH2D4RdI/NAmJWVNXXCuOdSYqf2CgGAdedKXtl3ccWGrZMzUjvy541WmLDRMkpLXk24yZiw6GAgbIu1QGilcPdOa5WJrhqnVHUsA56opNM2W19J4OwSCzEQtoWBsFMwEDJh7NixE0jpXf3Dm698d7xwM9Vs27atg69gaISkXy1Lh3F/6SGp3YIwELbFVCDkKSzYZS1toL9OULp1php2vJKOWW/ZMlk5KKC7Tz0GwrYwEHaKpApN8crNzZ0Ye80KgYmx6tzc3I6/QrAb/DpesSjHesDARM0GycS/TvDna+jq8Z2LggDQ3598mKKYvdVaZXJMyhDqMAyETOA4znLtmeNmnue4zn07AwLIZ+mKW7dai+owFiJnOF1Fl+Zbv8lQuHdp+vkdPbgpkeTunRYeMywSFAZCJqSnp687V9ryyvpzpenp6Z19nWmR3P19uAf33Wi9JkJ2YeHhnl3WfyYoevh0vW/z7URFlQmW5uO0ZyQkh68jRB2xdOnSjJEjrZRO6aUBgPXnSj85qtu264cuvNQzg7jBK/lfr/AzorCWgxzojSO8vwoe6NutbObCwYoxisErLbNjSB8/mU4RQILDQMiE3v0G+b52cM/2l7/euhcA0tPT93/5z9jY2C68lAsH/01T3L3TOi6M88SvFznG4XL67xPWQ7OU3Y9doR7k2cGKv2dbt0xiLr++//77X3/9tdCp6Aqr1cpxnMSmH44ePfq9995zxCszl/Pk6cVD1hH9Y7/7+3dGoxEAvL29u/Nq6SEkPYS8fti6dLjUVlMgFvAU7t1tfS/ZbhvEZMZz/++U9cXlayzncwghaWlpkyZNsssrd1NRUdHo0aPvvvtuoROCYOfOnbt27XLQi2MgFN6xCvrtef7Ybfac6Px2kmLgL+b5vTjsbkJ2t+Ii76qAu+y0HB4AGmprXN6fsfpC3uSeGgqw+NN/vz8keeXKlSzsmhYaGjpkyBChU4HgypUrjguEOIwkvJfy+KcGKYLc7PmaIe6wZIjiwX1WnI6H7MvCwyt5/NJh9tzG6Pnnn48sO73uL8mZw2P/Pjx2w50pcOHIK6+8Yr87IHQjGAgFdric7jfQB/vY/4t4OJ4rb4RNBRgKkT19fpaP9IIxofbsafj555+fSO7J/TGgpSDkieQeP/30kx1vgdANYCAU2EuH+KcHcl1bhnVjCgIvDOFez8elFMhuGizw2mH+9WH2HHumlBoMhgifa7pEInzcDQaDHe+C0A1gIBTSoTKaV07vi3PUtzA7hitrhJ06bBQi+/j3ST5ZTYYH27M5SAiJiYk5YTC2vHjCYOzarGmEugADoZCWHLK+MNghzUEbBYGnB2GjENlHtQnePWb9Z4L9C40HHnjgxV2nyht+32ytrN708q7TDzzwgN1vhFC7cNaoYLL19GQlrB7v2LrIvJ7cPw/zuXqapMbpo6hb/n2CnxzB9XXAPOQnnniioqJi7EcfDlN78RRySuoeeeLJRx55xO43QqhdGAgF80qe9YUh3E0Pb+smFw6eGsi9ns+vmYBrClHXmXj45BS/eZJDchHHcUuXLs3MzDx48CAhxNs0NLC3VmKLwRHLMBAK41QVPVpB10xwRtf033pzr+db8sro0CAsWVAX/XSR7+cP/fwdmIVCQ0OnT58OAFEV9JZNlsx4rrMnWiDUNThGKIyPT/L3xTm8OWijUsAT/bm3juK+xqjr/nWC/3s/JxUXAwPI4EDy7XnMsTfRkUbzkiVLQkJCnHAjUcNAKACjGb6/wD/ggLWD13N/H25LEV9cj9NHUafV1dVl62lFE0yJcF6OfXqg4u2jPB7PdGNPP/30TX/nk08+OX78eJdvMW7cuA7eSNQwEArg63P82FAuzE77NHaElwvMjuG+OIvlCuqoysrKzMzMgICAYD/fsfFhcfuWmU1NTrv7KC3xdYVNhZhjb+TNN9+86e+Ul5cHBQV1+Rbbtm3r4I1EDQOhAD45xWc6q5ep2UN9uU9P87jlGuoIs9k8duxYw47V22cPPvvIuN+m9qxZ/f6CBQucmYaFcdznZ7B39BrLly/XarVBQUH/+te/4I8eS0LIDz/8MGjQoMDAwA8++KDl78+cORMABg8e3LJv83p/VVZWNnv27ODg4Li4uF9++aXtn1dUVMyfP1+r1YaGht5zzz0VFRU3vruI4GQZZ9taRAEgPcTZfe6DA4naDTYX0kkREu/uR933yy+/EN2lt+ck2rY9i/Hz+HTK4JFfrTp+/Hj//v2dk4Y7enDPHDDr6hVaD+fcsEMarfCX7U5amOuuhO9HXzNf6Mknn9y9e7dKpXrkkUceffTR5utXr17Nz8/fsWPH1KlTH3vssebrq1evJoTk5+e3O8jX6q8ef/zx2NjYH374ITs7e+bMmVOnTm3154899pirq+vFixcB4NFHH33yySe/+OKLG9xdRDAQOtt/TvF/78cJEose7Mv97zQ/KQKn4qGbOHToUEZ0ENei9PRyVSaG+R06dMhpgdDLBW6N5r46xz8ziKGOKxcOFvR20uOrbBO9Ro4c+dxzz82fP3/jxo0trz/00EOEkDFjxjQ0NHT89Vv91caNG8+cOaNUKtPT08+dO6dUto4OGzZsOHnypLu7OwC89tprAwcO7M7dmYKB0Kmu1tLdOv6bDHueuNRxf+nBPbXfXFjHhTtxeBKJkYuLS42ldbun0cK3LRwd6r44bt5O69ODhKk4tktBYGaUYIF59erVmzdv/vLLLz/99NMtW7Y0X+/4Caa2E0/b/SvbWb62f5eUlKhUKk9Pz1Z/3tyyJIRYrdbO3p1ZDFW15OCzM/y8XoIdHO+phDt6cJ+fwXFCdBMZGRmbLugbW8TCYmNjnt6Ynp7uzGQkqYmHEjfL/VNMTExMTMyLL7546NChTv2hSqXasWMHpfSTTz653u+MGzfu7bfftlqtubm5iYmJJtPvO96ZzWbbPyZNmvT88883NjY2NDQ8//zzkydP7vIbYQ0GQuehAN+dp3/tLeRnfn8f7rMzOGUG3cSECRMGjZ541+pDewvKC2oaNl0ovXPVwcwn/hEZGenklNwXx312GqfM/O6pp55KSUnJyMh4++23O/WHr7322m233TZw4ECNRnO93/n3v/997NgxtVp9xx13fPHFF/7+/gAwefLkHj162H7hgw8+aGhoiI6Ojo2NNZlM4p0a0xahlIlCkRDC87y0l23uLaEP7bMeu+1G7UFbx4VDuxpGrLU8PZCbLlz3TqfwPG8wGG7w9MpQVVWVSqWyDdU4Ds/zf/nnZ3t+/T+PmoLevXtnZmZOmTLFoXdsV2UTxK4wX5jrEqC67u+YTCaj0RgYGGj3uy9evDgkJGTx4sV2f2XUWatXr/7yyy9Xr17tiBfHMULn+fY8f3cv4cPPgl7ct+fp9Cih04HYxnHc2YF/++GBhaO0QlZP/VUwNZL75hy/qL/wzw6SKsxbTtJkhV8u83f2EL7JOzuG21LEG81CpwOxLb+cVpthpKBR0Oa+OO6Ls9g7ihwIA6GTrC/gBwYQFqZr+qtgpJZbdRlLFnQjX5/j5/dkYqxipJZUmeBEJRODOEiSMBA6ybfn6d09Wfm074gl31/AQIiuy8LDDxf5u9jIsQRgdgz58aIUcmxTU1PzbEzEDiYyuuRVNsH2Yn5WNCuf9vQoLkdPS8W6+BU53KZC2sOH9PZloUEIADA3lltxUcQtwsrKyueffz42NtbPy9Pbw71v375vvPFGXV2d0OlCv2OlaJa2Hy/xkyI4X1eh0/EHDyVMjeR+viSFKjZyhK/P8fPZaA7aDAsmFgqHy0UZC0+cOJEwZHD15u+/GRV+7pFxFzLHfzws8MyK/yYnJV29elXAhHV8lv6Nf1MCs/0ZyusS9s05fh5LxQoAZCgvvfnSc3Pnzl28eHFnF+ciaau3wOYi/rYYhnIsAZgbS1aIsHe0qqpqxtQpbyRonhnRO8r3911Tewd6vZbR95Eo5Yzp06XdU2o7xYl9DOV1qbpspOdq6IQwhipN33777TNThsy88lOG8Ywie83kjLTXX39d6EQhVmwo4JOCyQ3W7QnC1jsquibhu+++OylYMSIioO2PpvYK6WUu+/zzz52fKqexneLEPgyEDvfLZTozinNh5pMuKyv7+yMPfzV5wNOpvab1Dvn78Nh1c5PeW/rq0aNHhU4aYsLPl+hslpqDNgMDiKcScvUiC4UrVqy4vW/o9X46Nz5sxYoVHX+1VscwrVmzZvDgwX5+flqt9p133gEAQshzzz2n1Wr/+c9/vvzyy7179/b19V26dKntR2+99ZZarU5PT798+XLza5aUlNx+++1qtTo2Nnb+/Pk6nQ4AKioq7rrrrsDAwJ49e9ru1VK7P22bmOZTnNr+iDXMZXfpWXWZoWkyALBt27ZB/qqBGp/mK1ovt6k9NWvXrhUwVYgRjVbYXMTPYHLjoTmxnLh6Ry0Wy5VLF6P9rnuOVA9/z9OnT3f8BZ988snNmzdnZWXZntYlS5bMmzevvLx8w4YNzz//vO13BgwYsHnz5pdeesnPz+/48eM///zzq6++avtRbW1tcXFxenr6448/3vyaCxcunD179pUrV/Ly8nr06LFw4UIAWLRoEQBcunTpyJEjeXl5rZLR7k/bJsa2C0x+fn676WQK7izjWCUNcLKKjgllqF+0urpa7dF63o7G07W6ulqQ9CCmbC7kBwUQtWO3b+uiubFkzHr+3SRg6DSKG+I4jnAKnl43wRaed3HpxFk0rY5hOnz48IEDB7766qtdu3Y1jzXedtttttfMzMxUKpXjxo1rbGy0/WjBggVKpfKJJ56Ii4trfs0dO3asW7eu+b/BwcEAsGHDhhMnTvj4+ADAm2+++dVXX7VMRrs/bTcxNjf4ESNYrPdJyarL/JQIzpWljzkuLu5waTV/7WhLXkl1y2cDyRab/aI2cb5E4w67S0TTO8pxXJ++fc+U117vF06V1Q4YMKDjL7h69erMzMyVK1dOmjQJAObMmfPhhx8GBwe/8cYbzb/j6upqm8ZpOzOr7ZROjuOaT1ACAH9//wsXLlBKKaW1tbUHDx60/U7zL7T7Cm1/2m5ibvojRjCa4yVj1WV+VjRb1de0tDS/2L6v7T1rtvIAQAG+OHL1tMn19ttvFzppSGAmHtYXMJdjW7o9llspqk2R5s+f/9XR666R+Pro1Xvuuafjr9bqGKYtW7Y8//zzU6dO3bRpEwBYLJYb//kXX3xhsVjef//9kSNHNl+89dZb33jjjfr6er1eP2PGDFusmjx58uLFi41GY11d3bPPPtvqddr96fUSYzabO5tO58NA6EBVJthvoBPD2fqQFQrFmjVrDJGDk5bvnvPLgdTlu9c1BW7ZssXW0YHkbGsRjfcnoR7sBsIZUeTXK2KaO/rwww+fcQlecaKo7Y8+OnCR9Bw8Z86cjr9aq2OYli5dmpGRMWDAgPLy8okTJy5YsODGf242m0NCQrZv3/7hhx82X3zttdesVmtMTEy/fv2io6Ntk1nef/99nuejo6MHDhw4atSoVq/T7k/bTYztFKfOptP58BgmB/rqHL/mCv1lnKLjf+KEY5iaFRUVnT9//oPi0ImDYx/s24lEOhMew9SW445hune3dUAAeYztcx7ifrKsGKMYHPhnWcH4MUw6nW7mzJkB5Vfu6B8WH+RtpXBMX/310QL3PgkrVqzw8/OzV2pvjBBWSvuuwWOYxGrVZTo7ht3QHhYWFhYWZrjELz/LP9hX6NQgoVl4WFfAvziU9TJhWiRZc5W2DISM02q1+/bt+/bbb//vhx9Obj6pUCgGDhz46LvP33bbbRKr+osX65levOossFPHfzmqE1PCBHFLOHffHmuNGXxYTylyrF0lNMqLRHmxXjRPj+KeyLG+OITpZmsrSqVywYIFwnYJfv/99wLenXFiykzisrGAT1YTP2b2F70eLxdI05BNBWKagIAcYc0Vtha8Xs8IDblSSwvrRNzLJ4i//OUvQieBXSLI9yK18jK9VQzFCgDMiOJWX8FiRe7WF9ApEaw3BwFAQWByBLcGcyyyH3GU1KJj4eG3Qn5apDg+3ulR3KZC3oRtQhk7XUVNVhgYIIJACADTI8maq5hfkd2Io6QWnSw9jfEm2uvurMQWjTv09SM7dVjFlq91BXRqpDiiIABMDOdy9LTGLHQ6kFRgIHSIjQX8JDH0MjWbGcX9egWr2PK14So/OUI0pYGXC4zQkN8KMcci+xBN1heXDQVURMUKAEyLJOuvYotQpqpNcKiMrR1xb2p6FA4TIrsRU2EtFoV1tLieJgaLqVjp40cIgZNVWLLI0W+FfFoI8RDVWqppkWRjAW/BNqGcmM3ml19+eebMmXfffffXX3/N83b7+kWV90ViQwGdGM4pxBQHAQAmhZNNBTTeT2zpRt22oYBOEVUHBgCEepAYb5Ktp+khjs2xu3fvxmXvLDh+/PiB/fuDCo9PiwysNVg+enbdF198sXHjRjc3t+6/OAZC+9tQQOfEiu/JmRhOPj7JPzFAZAUi6iaewqZC/pUE8RUFE8PJ5iI+PcSBuwO+9NJLL7/8cnFxseNu4SCNjY2urq4tj4kQu61bt2aoVa+P/n0TrFlx2jtXHfzoo4+6uQGeDSu7z0lmr1ETD+pvzRfmugSquvLnztxrtJU6C4R+Zy6608WLpS1mcK/Rtuy712iOnt6/13r0VvEFwj0l9Ikc64GZSsftNSpeZWVlvr6+nTrskGWlpaUx4WEnHxqrbHG048bzpd/X+2/fvr37ry+d+gIjduloP3/StSgoLE8lDA/GRRSys6GAF8U6+rZS1OSCkeobhE4Hcjyj0ejn5qK89oDjIA/Xmpoau7w+BkI721AgpmnorUwM5zbhlHSZWXuVThXJzg+tKDkYreW2FGGOlb6IiAgjKAtrrqn1HCiuio+Pt8vri/IBYNnGAjpZnPVrAJgUQTYWYItQRnT1cLWWJolqhnNLt4STTYWYY6VPpVItWrQoc9NRXW2j7cruq+Wf5Bc89thjdnl98Q0MsOyikdaYxXRATCv9/YmFwtlq2ttXrG8BdcqWIn5MKKcUbX34lgjywiErj6FQBl555RU3N7dRby7r4cHVmMxumoiff107dOhQu7w4BkJ7Wn+VTgrnRB1DJoaTjQUYCOViSxEdHybi7zrCkwSqyJEK2pP5Y15QNykUihdeeOGr2MdfDD83LMI3OjrajpMrRVsVZNKWIjohXMTFCgDcEk5+w0EXeaAA24p5UQdCAJgYTrboxP0WUAddNtJa6jp75OCYmBj7LjHAQGg3Fh72lPBjQsX9kY4L4/aV0HqL0OlAjnesgnoqSYy3uKPIxHBui/iW+aGu+K2ITnRMl5u4S22m5BpojDcJtsMuB0LycYEhQWR3CY66SJ/Y+0VtRoWQ/HIwWkT/RtBNOS7HYiC0m61FdJz4ixUAGB+GU9JlYUuR6PtFAcBdCUnBdI/BgfvLIBZYeNhRzI8Lc0jMwkBoN1uL+fGO+ZKcbGwo2VaELUKJa7RCVinN0Eohx04II9tLcd6fxO030CgvorHPfkqtSeExYIHRDEfKaZpG9PVrABgWRK7UUkOj0OlAjrSvlPb3J/4i3AKprXFaugMDodQ5tCcfA6F97NTxicHEXRIPo5KDkVpuRzH2jkrZliJe7DOcm/X1g0YeLhuxG0PKNhfxE8MdFbAwENrH1iLqoM5rQYwNJduKsViRss2FVBo9+QBAANKCrNtxm1zpqjLB8Qo6wmFdbhJ5EgS3rVgiM2VsxoaSrThMKF2GRrhkFNnZ0TeWHmzBgW0J21HMjwghKofNiMJAaAfF9bSkng4V7c5qbcX7kwYrvYR9TRK1tYgfpeVcJPT0j1JbthXjVmuStb2YOnSJtoQeBeFsLaKjQzlxb612LQIwJpTD3lGpksxSn2YRHtRDSU5VYY6Vph06OlrrwByLgdAOpFesAA4TStp2HR0XKrUcOyaUbMccK0WlDVBc79jDDDAQ2sF2nRR26GhlbCjZjn1NUnS1ljZZaR8/qeVYDIRStaOYHxXCKRyZYTEQdtepKurKQazIN2xsK9KL+LqS4xVYskjN9mKJrKNvZUwot0vHWzHDSs5OHc1wZL8oYCDsvh3Fju28FtAYLfaOSpCjh1uEEuIOIe4kvxxzrNRs19ExDu7Jx0DYXbtK6CgpFisAMDaMbMNl9ZKzy/H1a6GMDcPeUakprqeVTbSfPwZChlGAXTpeqsXKaC23p4RaMBRKyIUaaqHQS6IHL48Jxaqb1Gwrphlah8/Jx0DYLaeqqJcLifSSZrES5AbRXiQP+5okRKr9ojYZWi67lDZZhU4Hsp8dxQ7vFwU7BsLKyspp06YFBARMnz69srKy5Y9SU1PJHx588EF73ZEFO4sl28tkM0pLduHOVRLihHkHAvJzhd6+ZL8Bc6x0OKfqZrdAuGzZsqioKJ1OFxkZ+dZbbzVfp5SeOXNGp9MZjUaj0fjBBx/Y644skPAAoc0oLdmlw74m6ZB2IASAMaFkB1bdpOKS0UlLfewWCFetWpWZmalSqTIzM1euXNl8XafTWSyWKVOmaLXaefPm1dTU2OuOgrMNEI4KkXKxMjKE21dKcUq6NJytphxADx8p59hRWm5PCVbdJGJ7MR3tlKU+djs3qKioKCoqCgBs7cLm6zqdLiEh4b333ouMjHz88ccXLVr0/ffft/sK9913HyGtH9GEhIQ777zTXom0r9M1xEPh5mutqa6222vW1tYCAM+z8iS7AIS4ue27YhzkL0ySeJ6vqalxc3MT5O5sqqmpUalUJpOps3+48aIyLYirrq5zRKoEZDKZamtrlUolAAz0IDmlbobKald5z3+wNTlcXFyETki3bL6qSgu2VldbuvMiXl5eCsVNtuu2WyCklNrCGKXUav1ztDohIWH79u22f7/55pv9+vW73isMGDCg7cWoqChmv8vsci5dTe2bPNurMfWWR4VAVrnLMLVggdDFxYWpD0RwLi4uSqWyC59JVrlyjJaX3odJKW3OJAEu0MsHjhtdk4Jk3Y/h8gehE9Itew2KJQO7m2Pbtq/aslsgDA0NLSgo6NWrV1FRUVhYWPP1Q4cONTU1paamAoCrq6tKdd0jsRctWtSRFLMju8I6NZJ4eNjzkG9bHcLDw8OOr9lNYyP4b8/TZzwcdgLKDfE87+7uztQHIjiTyaRSqdzd3Tv7h3v05mVJKg8PMT1lHaFUKq1Wa3MmyQi15lYqR0fKuklYX1/v4eEh6kB4oYYqiLW/xhnPvt3yyrRp05YvX04pXb58+YwZMwBg586dAFBXVzdr1qxTp06ZTKZXX3115syZ9rqjsOQwQGgzKoTbW8LjrqNid7qKunAkWnJ7AbY1UktwmFACdpXQkc6a2GW3QPjiiy8ePXo0IiLixIkTL7zwAgCMHj0aANLT05csWTJt2rSwsLDKysply5bZ647CkvYKwpbU7qB2J0dx01GR21XijPVYLEgP4bJwhpf47XN8An4AACAASURBVC2h6c5qadita9TPz2/9+vUtr1BKAYAQkpmZmZmZaa8bMWKH1FcQtjQqhOwucewxKMjRduskeFhYuwJVEO5J8stpQpAs3q9U7S6h/xjopP5tWXejd8cuncRXELaEy+olYG+p8+rXgsMcK3bF9bTG5LzDwjAQdtHuEn6knIqV3SU4SihiV2qpmac9Jb2CsKWRIWR3CWZYEduto2khjt5h9E8YCLvibDV14UiUDAYIbUI9iJ8rOVmJJYtY7S6h6SEyethHark9OMNLzPaUUme2NGT0bNjRnhKnfkkswL4mUdvjxHkHLAhxhyA3cgKrbqK1W4eBkHl7S2manIoVABilJTsxEIqWMyfgMWKUluzC3lFxKm+Cgjo6yImz8zAQdoUMW4RpGlybJVZljaCrp/0dfLQpa0aGkN1YdROnvSV8qoYonJhhMRB2WkkDVDTRvs6azsSIGG+i5MiFGixZxGd3CT/CucUKCzK0ZBfO8BKnPU4f0sZA2Gl7Svg0jcNPTGZQWgjZW4oFi/jsldlMGZtwT+KhJGeqMMeKz+4Smq5xagkru8ej++Q276BZqprsw0AoQruduFUVU9I0mGPFp9YMp6rosGAMhGyTbSBMDyF7cfaB2NSa4Wy1TPdYGYGBUISy9HRYEHFz7ib/GAg7p9oE52voUFkWKwMCSHE9NTQKnQ7UGVl6mhBE5Hk43wgNduaLz54S3vktDVk+H92QrafD5FqsKAikqMm+Upw7KiaCFCuM6OdPyhqpvkHodKDO2FNC05w+pC3LEr0b5FysAMAIDZeFVWxR2a2jI+U3U8aGI5CiJnux6iYeZh4OldEk5w4QAgbCztorRG2FHWkhZA8OE4pHkxXyymmyWtZVNxwmFJG8MtrDm/i6Ovu+8i3Tu8BWrKTIuFhJCibHKmi9Reh0oI45WEbjfImXiE8p7660EJwvIyb7BNq0CwNhJ9iKFW8ZFyvuShgQQA4YsGQRh6xSOsK567FYMzyIHMeqm3jsEyjHYiDshL0lsttitK00DfaOikZWqaw7MADAXQkDA8h+rLqJRFYpn4qBkHFZepom7/o1AIwIwYmjopGt50fIvuo2QoPrX8XhfA1VCHS8HQbCjqIA2aW8zOvXAJCm4bL11IoFC/POVVOVgkR4yj3HYtVNLPaVCtbSwEDYUWerqacLCZN9sRLkBloPcqwCIyHrhBpuYQ1W3cRCwByLgbCjsmU/3NIMt3AUhWw95liAP6pux/GQXuZll1JBBggBA2HHYf26WaqGZOuxWGEdzu1qhsOE7Ktsgiu1dFAABkK2Yf26WaqG4P4yjKtsgoI6OkBmh/FeD+6+zb5sPU0MJkqBIhIGwg6pbIKCWjpQoNoKa3r7khoT1dULnQ50fdl6Oly4YoU1I7Dqxrx9pbyAXW74oHRIlp4mqrFY+R0BSNGQLJyJx7B9pTwu9WnWy5fUW2hRHcZCdgm7eyUW7R2SXcqnYr9oCylqDocJWbavlKZq8On+HQFIVnM5mGNZZeYhr5wmOn2v7Wb4qHRIVilNwWKlBRwmZJmZh7wyWe+13VaKhmRhIGTV4XJh9tpuhoX7zVl4OIjFyrUSg8nRCtpoFTodqD355TRa0GKFQSlqko1VN1ZlCbdwwkYp4L3Zd/HixZ9++invfKGXNcY0aR6o1UKniBUeSujrRw6V4ZISFgm4Qwezmqtubgqhk4LayCql06OEzLHYIryu//73v0lDBl1Y8cmA89tGH3i/X98+a9euFTpRDEnB1YSsytbTFAyE1/JQQh8/kleGOZZFOXqBu9ywRdi+M2fOPPPEY7/ePryHv6ftyvSCir/ec8+FS5d8fX2FTRsjUtXkx0sUBgidDtRGVil9YzjWcVuzbQQhbBccaqu4njbxtKcPtgjZs3r16ik9gpujIACMiAiI8yRbt24VMFVMScUVFEwqrqcmnsZ6Y3HfWooa+zBYtK+UpqgFjkQYCNtXVlYW7u3e6mKYj3tZWZkg6WFQpBdRcuSiEUsWtmSV0mShixU24XwZNgneLwoYCK8nNjb2uKGm1cUT+prY2FhB0sMm3LCDQVmlFNe8tivamxACl7HqxhjBp4wCBsLrmTNnzv5Ky5qzOtt/KcBHBy7yQWEZGRmCpost2NfEIJwpcwPJmGMZ02SFYxU0IQgDIZMCAwPXb9jw8RXL5O+zH95wZMzXe/cS9dq1a11cXIROGkNS1NgiZAsjxQqzsOrGmkNltI8f8RR61qbQ92dYYmLi/rz84Nf3PxVf3C+u5+DBgwnB8uUaQ4LI+RpqNIM3Vg/YcLicxjFQrDArRU2+v4AzvBiSw8apPvjE3MixamXfoUnzZuKn1D5XDgYHkgMGOiZU+KyMwLYXIAPFCrMSgsjpKlprBi+surEhS09vixY+x2LX6I1kMzCdiXHJwQT3MmYHnpp5YyoFDAwgB3FZPTMYaRFiILyRXAyEN5OsxkDIEMyxN5WiwRzLiqu11MrTaAbWvGIgvBHsaLqpFA3J1vNYrrDAtkNHD0F36GAfVt3Ywc6pPkwkgk3F9bTBisXKTYR6EHcFuViDJYvwWNihg33JapKtx/kyTGCnJx8fm+vKLqUpao6Jb4ltuDaLEdnYgdEBEZ7ElSOXcFk9AzAQikCOniax8SUxDvuaGIE5toOSMMcyoNEKJytZWfOKgfC6cKP6DsJAyIImKxytoMODMcfeHOZYFhwqo/H+xJ2NtWkYCNtn5uFoBR3ORm2FcQlB5Ew1rbcInQ55yytnYocOUcA1PyxgYa/tZhgI25dfTmO9Ca667QiVAvr5k0O4NktQudgv2mEJQeREJW20Cp0OedtvoInMdGBgIGxfVinuXNwJ2NckOEYWJouCuxL64mn1QmNqbhcGwvblGhhqtrMP+5oEl8tS/Zp9WHUTVnE9bbTSWGYWp2EgbB8783pFAddmCUtXD0YT7eWLObajMBAKK7uUJrO0OA0DYTv0DVCDxUpnRHsTAuRqLZYswsg18MlqPBulEzAQCivXwNaQNgbCduTo+SQsVjopMRiX1QsmV0+TcE+ZzujhQ0xWPv+izmrFOTMCYGrKKGAgbBcOt3QB7mUsIJwy2il1dXVPPfVU+f2B4wdE+/r43H///eXl5UInSkYsPOSXs7XmFZcdtSNXT58coBA6FSKTFEyePYDDhAKwUjhUhlW3Tpg9ezY5fzj3rmFqT1V1k3np3o3jxx/Mzc11ccH1Us5wtIJGeREflj5sbBG2xlM4WEYTsX7dScODydEK2oT9TE53vJKGexI/V6HTIRI5OTnHs3b9+5aBak8VAPiqXN4YE091F3/99VehkyYX2XrmFqdhIGztZBXVuJNAldDpEBsPJfTyJUcqsHfU2XL1WG/rhPz8/OQwf5Xiz6KPIyQ9IjA/P1/AVMkKU0vpbTAQtoYLk7ssKZjk4jCh07E274BxKpWqzty646LWZFGpsPLrJAwuTsNA2BrOO+iyJDXJNWAgdDYMhJ2SkZGRU1xVbGxsvmJssmy+qB87dqyAqZKP8ibQN9C+fmzlWAyEreXoaRJjzXaxwLVZzldtgsI62o+xYoVlMTExj/7jmbkrD6w+oztbXvvbBf2cXw7cMueu1NRUoZMmCzl6mhhMGFpLDwAYCFsxmuFyLR0QwNi3JBJxvqS8kZY13vw3kb3sN9ChQUSJz3FnvPzyy5/+uHq7e48HDlUvKta++PFnn332mdCJkov9ep61AULA5ROtHDDQwYHEBYuVLuEIDA8muQY6JYK5jC5V2C/aNePHjx8/fjxPIeAb89jpLE3kl7pcA82MZ25xGhb512BwFFdcktUkFzcddaIcPY+BsMs4AglBZD8ObDsLZXLKKGAgbIXNL0lEktQcThx1GluxkhSMT3HXJatxqrPznK2m/q5E7S50OtrAR+gauVi/7p5kNdlvoDwWLE5xvpp6KonWQ+h0iFlSMMnBPgxnYbYnHwPhny4ZqYKQcE8WvyexCFRBoBs5U42R0BlY28JfjJLV3H4DxfzqHMwuTsNA+KdcVmsr4oLL6p0Ge/K7T+0OPq7kHFbdnILZ8wwwEP4pB0+ltwdcVu80zHY0iUtSMOZYZ2iwwNlqOiSQxRyLgfBPzDbbxQWX1TtHkxVOVdGhTBYr4pKE82Wc4lAZ7edPVMwtnQDAQNjMxMOxCpoQhMVKdw0OJOeqaa1Z6HRIXV45jfMl7rgSuNtw4qhz5BjY3bQLA+HvjpTTnj7EE4uVbnPlYEAAySvHksWx9usZHW4RnSGB5FQVbbAInQ6pYznHYiD8HQ632BH2NTkBDmnbi5sC4v3JYay6ORjLZSwGwt/hTHQ7wtkHToCTnO0oKRgHth1LVw/1Fhrrw2iOxUD4O5yJbkfYInQ0QyOpMtFevphj7QOnOjtaroFPUhNm8ysGQgCACiaPyBKvWG9i5mlRHZYsjnKgHBKD2S1WRAcXvzpajp4mq9kNN+ymzJly9XQ4e0dkidpw7B11pANlBPtF7ainL6mzUF290OmQLsYXp2EgBLA127Ff1K6SgnH3bUfR6/XZBUbca9uOCMDwYHLAgJuOOgRP4XA5Hc7w4jR7PkuVlZXTpk0LCAiYPn16ZWXlTa+zAwcI7Q4HXRzh008/1Wq1aUMG5s3Tvjg7df/+/UKnSDqSgjnMsQ5yopKGuBN/ldDpuD57BsJly5ZFRUXpdLrIyMi33nrrptcZQQH262kSw/3XYjQ8mOSVUSsWLPbz7rvvvvPsk1+Pjc1bmHHyobHz/esnjB1z/PhxodMlETjDy3HYn5Nvz9J/1apVmZmZKpUqMzNz5cqVN73OiHPV1NuVaNg7IkvU/FwhzJOcqMSSxT4sFssbb7zx8aSBfYO8AYAAzOqjva+/9p133hE6aRKRGEwOGLDq5hAsryC0sedOKkVFRVFRUQBga//d9Hort956a9uLKSkpCxcutGMi29pxVTnEV1lZaXToXTqotrYWACwWKexyMdTPbcflxgjSrc3WeJ6vrq52dXW1V6pEqqCggDQYbVGwWVpk4EuHD7M53OBMJpOptraW47pVrScAwSqP/Ver+/hIYaSwqqqK53kXFxehEwIAkF3icVdYY2WlMB+sj4+PQnGTHU7tGQgppYQQ2z+sVutNr7dy6623kjYTwqOiotzdHdtYy6/hUkPA0XfpIFsIZCQx3ZSsJnkVige6txsmz/Nubm7S+EC6w9/fv85ktfBU2WJyc02T2cPDAz8chUJhtVq7/zkkBZOjtW5DNFJoFbq7u7u7u7MQCI1muFrPDQtRuQg0+tQ2rLRlz0AYGhpaUFDQq1evoqKisLCwm15vZd68eR1Jsd0dLLfMj1O4uTHRcjebzQDg5uYmdELsIC2MfnrO6ubWrcYcz/MqlUoaH0h3RERE9O7Xf81Z3a19Qpsv/niyaOKdD+KHw3Gc2Wzu/ueQquXzKugDbkyej9BJtqeGhUCYXUkHB1q9PZjOpfaM0dOmTVu+fDmldPny5TNmzACAnTt3tnudHY1WOFVFBwcwEQUlZoA/uVxLjXgMhZ385z//+WfOlXdzzh/SVe25Wv7A+vxCj5AnnnhC6HRJB86XcYRchvfabmbPQPjiiy8ePXo0IiLixIkTL7zwAgCMHj263evsyCujff3wLBuHUHIwOJAcwCnpdpKampp//MSB2EkPnuSWG/3GPPiP3NxcHx8fodMlHYMCyAUjniBmZ+zPlAH7do36+fmtX7++5RVKabvX2ZGLW/g7km337TGh+AnbR3h4eMi9H/xN3XBHTwUODdqdCwcD/MmhMjpKiznWbg4Y6L9SWF+cxnr6HI3xjX/EDvua7C5HTxMCrjvjDHUTbgRhX1drKU9ppBfrZazcAyHuKeNQScEkVy+FyeiMKK6nJiuN8sSS2lFw92372m8Qx14lIkii4+gboAbPsnGkSC/CEXKlFksW+2B8C38JSFLjwYT2JIoBQpB5IMzW84kMH5ElDYlYxbYf7Ml3tBhvYqV4gpjdsL+5mo2sAyH2izoBDrrYEeZYJ0jE3bftxMJDfjlNYPjQiWayDoTY0eQEyThfxk6sFA6V0eEYCB0MZ3jZy9EKGu1FfIRf039z8g0DPBYrTjE8mBypoCacMdNtxytpuCfxk/uuqw6XFIzDhPYhln5RkHMgPFlFNe4kkOEjsqTBUwmx3uRYBZYs3ZUrknkHYpekJofLqQWrbt2Wq6diOfBcvoEQh1ucBmfi2QXmWOfwdoFIL3IcTxDrtv0GmiiSqpt8AyFOwHMaXJtlF5hjnQZ7R7uvygRFdbS/vzhyrHwDYY6epmCx4hQ4cbT7jGa4UksHiKRYETvMsd2Xq6cJQUQhkgwr00BYa4aLRjoAD51wing/UtpAK5qEToeYHTDQwYFEKdPn1dlw4mj3iahfFGQbCA+W0UEBxFWm797ZOALDgsh+rGJ3QzbOlHGi/v6ksI5WYtWtG3L0vIi63GQaCnJwuMW5ktQkBzcd7QYRTcCTAAWBhCByoAyrbl1EAfYbxLQ4TaaBECfgOVlSMNlzTn/p0iXbyVyos/YbeKy6ORP2jnbHhRrqqSShHqLJsTINhLgky5l27dr1zJTB+X+NHjW4nzo4+OOPP8Zw2CkXjdSFI+GemGOdBw9O6Q7R9eTL8Wj2K7WUAo1i/ogsaTh48ODMKZPeHdNnwrjRAHCuovbBl55paGhYvHix0EkTDbFs4S8lyWpu4R4rBcDPvQtEt9RHji1C3GLUmd59991HhkRMiFXb/tsrwOvftwx86623rFY8XbajRFesSIDWAzxdyMUa7LroCtF1uckxHuw34LwD5zl27FhqeGDLK/FB3mZjVWlpqVBJEp1czLFCwGX1XdNohVNVdEigmHKsHANhdqnIaiui5unpWdVkbnnFbOXrzFYPDw+hkiQuTVY4USmOs2wkBpfVd82hMtrPn7gphE5HZ8guEJp4OFqBxYrzTJw48btjBS2v/HiqeOiw4X5+fkIlSVzyymmcL/GQ42i+wPAEsa4R4+I02T1eR8ppL1/iJYYjsqRh8eLFGevW/W3t4Tv7h3u4KHZcLvvlcvXG3zYLnS7R2C/CYkUahgaSE5W00QriatwIbr+BTosUWY6VXYswG7cYdS4fH5+cnJwpjz63koYtvuJ1asBtx06cTEhIEDpdoiG6meiS4a6Evn4kD5fVd1J2qfjKWNkFQpyA53yurq6LFi1as2bNWz9sDZr3pkajETpFYpJroMk4U0YgyXiCWCfp6qHeQmN9RJZjZRcIcUmWgJLVJLsUi5VO0DeA0UR7+mKOFQYGws7K0fMpGiK6/CqvQKhvgEoT7Y3FikCivQlPaUEdliwdla3nE9XiK1YkAwNhZ+UaaGKw+MKK+FLcHTl6PikYixUhJas5bBR2HG7+IKwePsTE00KsunWYGAcIQW6BMNeAA4QCw7VZnSK6HTqkJzGYw0ZhB1l4yC8X06ETzWQWCPU0SYTNdilJwb6mDrNSOFRGh+OaV0Fh72jHHa2gkV7E11XodHSejKICT+FgGbYIBTYsmBwpp024z2gHHKugYZ7EXyV0OuQNA2HHZelpikaUBayMAuHxSqr1IAFYrAjKUwm9fcmRCixZbg7XvLIgSU2OVFATnsjUAeKdky+jQIhnfDMiCavYHYMDhCzwVEIPb3KkHHPszYk3x8ooEIq32S4x2NfUQVl6moo5lgGYYztC3wBljbSPOBenySgQinRer/RgsdIR5U2gb6B9/TDHCg9zbEfk6PlkNeHEmWHlEggrmkBXT/v5i/NbkpbevsRopiUNQqeDbdmlNClYrMWKxGAg7IhsMa95FWu6OytHT4cHEwUWKwwgAMODSa4epx/ciK1+LXQqEABAnB+pMtFSrLrdkHgHCEE+gTBbz+NwCztS1Fw2VrFvKKuUpmrk8ngyjgAkBpMcrLpdn5XCwTKaKNrZiHJ50rJKaYpom+3Sk6ImWbjR2vXZltKLt1iRnmQ17i9zI0craISY17zK4mBeW7GCS+nZkaQmh8upiQdXrJy053glLqVnS4qGLM3HFmE7Nm/evHbt2tzL5e7a+PLRDwUGBgqdoq6QRTl0vJKG4lJ6lni7QA9vko9rs64jC2c4MyZZTQ6VUTOGwhZ4nr/nnnse+sutQce23m0+0Tfn0/g+cfv27RM6XV0hixYhLpxgUKqGZJVi71/7skvpSC1+MgzxcYEYb3Kkgg7DrV//8NNPP2WvW7npzhR3pQIAZvXRJp3RLViw4OzZs0RsZ/zIokWYjUvp2ZOiIThf5npwczUGpeLA9rXWrFmzYFCELQrazIzT1pYUnjhxQsBUdY0sAqF4d8CTsBQ8rf46yhrB0IhL6ZmDVbdWqqqq1J6tB5zUnq5VVVWCpKc7pB8IyxqhtIHGY7HCmJ4+xIxHnrYnR49L6VmELcJW4uLi8nTVLa/UmiwXKht69eolVJK6TPqBEIsVZiWpOSxZ2sKl9Gzq6UsaLLQIq25/eOCBB/7vjH77ZYPtv8Ymy+Obj8+8fY5GoxE2YV0g/UCYjcUKq1Kxr6k9WaU0BZfSs4cApGi4LMyxf4iLi1u9dt3rp+vGfrtv9s/7R3y9N3r8rP/9739Cp6srpD9rNKuUPjNIcfPfQ06XoiaLc3FC+jUsPBwqw5kyjLINbN8eI3Q6mJGRkfHljmN//en0W4Mr4+Pjg4KChE5RF0k8ENqKFVxKz6ZhQeREJW2wgLvEs2En5FfQKG/i6yp0OlB7RmjIP/Zj1e0a+8sVYxLiR44Qd2ND4j0wRypolBfxw2KFSe5K6OdPDpZhX9Of9pbQEbjUh1XDgsjxClpvETodLNlXKoVTMyUeCPeVYrHCNBwmbCULcyzDbFW3Q1h1ayFbL4UcK4NAGCL6L0nCcDVhK9l4Kj3bsOrW0mUjNfM0xlv0OVbigTCrlKbiACHDUjQkG0+3+cOVWmrhaaz4ixUJw4NTWsrW0xGSmOEshfdwPVdqqZmnPXywWGFXhCdRKci5aixZAH7vwJDyIykBaSFYdfuTZMaepPzU7S2haVisMG+EhuzFKjYA2OYdYAcG20I9iIoj52swxwL8fny0FHKslONEth6LFREYoSH7MBACAEC2VIoVabMdnCJ0KoRXa4ZzNXRIoBRyrJQD4Z4SmoYzZZiXFoKBEADAaIZzNXQonvLDPKy62eQa6JBAohL3AsLfSTYQ1pjhopEOlkRtRdoG+JOSeqpvEDodQsvR06GBxFWyT6R0pIWQvSUYCCWygtBGso9ddilNCMJiRQQ4AslqnIAA+0p5acw7kLyBAaS4nhoahU6H0LJLecmMPUk2UGRJ6EuSvBEaDvuaskvx+GhxUBBIUpOsUllX3XhqO+dVIhFEIm+jrb2lOGVUNEbIvq/JSiHXQFOkUqxIHlbdjlZQjTtRuwudDjuR5oNn4eGgAU+lF42kYHK0gjbIeAvH45VU60GC3IROB+qYdNlX3faW0nQJTUWUZiA8WkEjvEiASuh0oI7xUEI/f3JAxls47iuhadgvKh6JweRYpayrbnulNSdfmoFwbykWKyIj85l4u0voSC3mWNGwVd32G+SbY/eUYIuQeZLZ70A+RmhkPfsAq26ik64he+RadbtQQwFASpviSjMQ7imhIyVUW5GDNA2Xpae8LAsWW7EigS38ZWVECMmS65ofiTUHQZKB8HwN5QhEY7EiKmp3CHIjJ6vkGAmx3iZGaRouq5Ra5ZhhpTZTBiQZCPeU0FE43CJCaRqZDhNKr34tB0FuoPUgxyowx0qBNAOhxL4kmZDtpqPSq1/LRHqIHA9O0TeAoZH295dUjpVgINyNHU3ilK4hu3SyK1ZKGqCskcb7YY4Vn1S1HKtue0r4VDXhpJVhpRYIi+tpjYn2wWJFhHr5EiuFS0Z5lSx7S/i0EE5ixYpMpIfIceKoJDftktr72a2jaSFYqojVSC3ZLbOSBXvyxauHDwGAi7KruklwqY/UAuGeUuwXFbGRIWS3zHpHsSdf1OSWY2vNcLqaDguWWo6VWiDcrcNiRcRGhsirRVhjhgs1eGqmiI3Skl1yyrFZejo0kLhJ4jDeliQVCMuboKCODsJiRbTi/Um1iRbVyaVk2VdCE4Px1EwRGxUirxlee0t4KW0x2kxSj+DeEj5VQxQS/JrkggCkh3DyaRTuLuHTJTfvQFb6+JFGK71SK5ccu0tHM7QSzLGSekt7SigWK2Inq97RPbjXtvilaeRSdWu0wuFyaW7jbLewUVlZOW3atICAgOnTp1dWVrb8UWpqKvnDgw8+aK87trW7hKZL8UuSlVFaufQ1NVjgaAVNkty8A7nJkE2OzdbTAQHEUyl0OhzAboFw2bJlUVFROp0uMjLyrbfear5OKT1z5oxOpzMajUaj8YMPPrDXHVupNcOpKglOZ5KbgQGkpIHqG4ROh+PtK6WDAomHFIsVWRkpm0C4s5gfLdEODLsFwlWrVmVmZqpUqszMzJUrVzZf1+l0FotlypQpWq123rx5NTU19rpjK1l6mhAkwelMcsMRSNNwe0qkv6//Th2fIdFiRVb6yWaG104dHSnFAUIAsFt1tKioKCoqCgBs7cLm6zqdLiEh4b333ouMjHz88ccXLVr0/ffft/sKGRkZbS+OGjXqscce60gCNl5UJfnRigpTV1LPjNraWgAwm81CJ0RICT6umy+T0b5NAMDzfFVVlYuLi9CJsr8tBR7PxZsqKjp9zHl1dbVKpWpokEGruWNMJlNtbS0hgtUqkgLcN5xvuC2CoRPrKysrrVarHR+cJp7klXnGuxorxLbPuK+vr0JxkxZStwJhnz59zpw5AwCUUkqpLSNSSq1Wa/PvJCQkbN++3fbvN998s1+/ftd7taeffrptVg4JCfH09OxIYvaVk6VDwdNT3CUmz/MA0MG3LFVjIyAzl/P0VAIAz/MNDQ3S+0DqLXCqhhsVrnJXqjr7t2azWaVSubu7OyJhYuTi4kIpFTCTjA4j+6sU8/swFCEaGxs9m3QIggAAGa5JREFUPT3tGAgPlJAB/hDs42GvF3Qajrt5K7ZbgfD06dPN/w4NDS0oKOjVq1dRUVFYWFjz9UOHDjU1NaWmpgKAq6urSnXdx37SpEldrtMZzXC62pwW5qISedeoyWQCgBt8SnKQrIXLteY6UAWogOf5G2cbkdppoEODrH6eXXlfqj/YPVUiRQgxmUwCfiBjI+j/O2dVqRga77U9NXYMhFll1jFhoBJ7CXsdduvwnTZt2vLlyymly5cvnzFjBgDs3LkTAOrq6mbNmnXq1CmTyfTqq6/OnDnTXndsaXcJTQzGAUKJUHKQrCbSHibEAUIpGeBPDA20RNJ91RIeIAQ7BsIXX3zx6NGjERERJ06ceOGFFwBg9OjRAJCenr5kyZJp06aFhYVVVlYuW7bMXndsaUcxPzpUsl+SDI3WcjskPRNvh46Olm6xIjccgREh3G6dZKtuthWEI6S7OM1ubXk/P7/169e3vEIpBQBCSGZmZmZmpr1u1K4dOvrvFCxWpGNMKPnbbskWK0YznKykSWrJFisylKElu0ronFih0+EYEl5BaCOF4FFlggs1dDiuIJSQoUGkuF6yfU17S+hw7MmXlgwt2V4s2T4MCa8gtJFCINxRzKeoiYsU3gr6nYJAegi3s1iajcIdOl6SGzbK2eBAUt5ICyW6mnCnjo6SdI6VwnvboaM4QCg9Y0LJdokOE+4optKuX8sQAcjQcpJsFEp4i9FmUogfO3VYrEjQ2FCyrUiCxUq1Cc5U00QcIJScsWFkmxQD4b5SiQ8QggQCob4BCmrp0CAsVqQm3p/UW+hlo9RKlj0lNFmNZxBK0NhQslWKVbetRfy4UIkXsKJ/HHfq+PQQDs8glB4CMDqU21EidDrsDQcIpaqnD3Hh4HSV1GLh1iI6LkziOVb0b2+Hjo6Rem1FtsaEkh2S62vaUkTHYo6VqLGhUps7Wt4E52postR78sUfCIvpaCxWJGpsKNmuo1IqV0oaoKgODwuTrLFhZKu0AuH2Yj5NI/05+eJ+f4V1tKKJDvDHYkWaYryJm4Kcr5XOgrutRXyGFnvyJWtMKLdLx1slFArl0C8KYg+EvxXS8WEch8WKdI3Wwt4ycZ8o0tLmQjohHPOrZIW4g9aDHC6XTiTcWkTHh0k/x4o7EG4pwmJF4saEEskEQgqwtZiXQ7EiZ1Ja9nPZSOstNF4GXW4iDoQ8hW1YrEjdaC3ZV6aQRl/TiUrqriCx3phjpWxcGNkqlR2RNhfRcWGy6HETcSA8VEa1HiTUQw5fk3yFuEOYOz1gkEIkxH5ROcjQcrl62mi9+W+yb4s8+kVB1IFws2y+JJkbrTZtKpRCFXtrsfQXJiNvF+jvT7JKRV914ynsKOZlstRHzIGwkJ8YLuL0ow4arbb8Vij6YsXEw74S3BRXFiaEk81Foq+65ZVTjTsJ88RAyDCjGfLLabqk94FFNokBllNVtLxJ6HR0z94SGu9PAlRCpwM53i3h3MYC0VfdthbR8bLpyRdrINxRzCepibuk94FFNi6EjtJyW0Vexd5axI/Dnnx5GB5MiutpkciPZNpaJJd+UbDjCfXOwfP82rVr8/LytpS5J6aPAUgUOkXIGSaGkd8K6VwxH/+9pYi+myydnQHQDSgIjA/jNhXSe+PEGkjqLHDAQOWzKa6Y3mdJSUlycvKrDy2o3fx/g3M//en+jAceeMBqlcT0LHRDkyLIpkJevBVs24aNKVLfsBE1uyWcbBLzwPbWIj5RTbwlsoL35sTUInz44Yf7mUr/OSfJVpwsSupx64/ff56QcP/99wucMuRgMd7Ey4Ucq6ADA0QZS7YU8iNDOMlv2IiaTQznHssxW3iFUpxf+sYCOklOUxFF81br6urWr1v3dGqv5oLQ00WxKKnHDz/8IGSykLPYekeFTkUXrSugUyNFGcJR12jcIcab5OjFmmM3FdLJETLKsaIJhOXl5T4uxMv1miZspK97aWmpUElCzjQxnPtNnKsJrRR+K+RlVawgAJgUTjaKM8ceq6AKAn38ZJRjRRMINRpNPVWU1ZtaXjxdVhsdHS1QipBTjQ4lBwy01ix0OjovR0/DPUm4PNZjoWa3hHMiHSbcUCCv5iCIKBCqVKp5d9/97PaTjZbfZ8cUGRvfyT533333CZsw5ByeShgWTHaViK9kWX+VnyKzYgUBQIqGXDFSXb3Q6ei8DQX85AjRhAa7ENNkmffff3/hwtqMb1Ylhfo3WKwHDA3Pvvz6rFmzhE4XcpJbwrmNBfyUCJEtQthQQP8zQmRpRt2nIDAmlNtcxN/TS0xBpcoE+eU0QyuvqpuYAqGHh8d333138uTJw4cPe3h4fJGaqtFohE4Ucp5pkWTcRv7fqSCiZ7SgjhbX0yRcOCFLt0SQTYX0nl5Cp6MzNhfy6SGy26tEfG83Pj4+Pj5e6FQgAfTxI55KyC+nQwJFE1fWX6W3hOOR9DJ1Szh5KtcqrkUUGwvpJJn1i4KIxggRAoDpkeTXK2KaibfuKj8FF07IVagHifEme8VzEgVPYVMBP0l+Q9oYCJGYTI/i1lwRTbHSYIE9JRTPSJGzGVGciKpuh8pogEqOZ0fjI4rEJFVDCuvo1VpxxMIdOjokiPi5Cp0OJJyZ0WS1eKpuGwqoDJuDgIEQiYuCwJRIbu1VcZQs667yUyPxEZO1/v7ElYP8cnHk2FWX+ZlRcsyxcnzPSNSmR5I1YuhrorJcmIzamhZJfhVDo/CikZY20FRZHvKKgRCJzMRwLtdAq0w3/01hHTRQlQLi5bRPFWrXTJEME/5yic6M5jhZZlgMhEhkPJSQHkLY33f0l8v87BhZFiroWqkaoqunl4ysNwpXXeZnRcs0IohvHSFC0yO5NVdYP6d31WX6/WjcUAYBR2Bo+e4H/7FjgFfj0KFD58yZo1QyV/CWNMCZatltKNNMpvEfidr0KG5jIW9muE2YX05NPAwJkmmxgpqZzea5c+eefnnGoKMrPPav+88zf09ISCgqKhI6Xa39comfGsm5yjUgyPV9IzHTuEOcL9mpY7evaeVlfnY0wTCIPvroowt7Nm+9K/XRxB73D43+4dbhiYrqzMxModPV2uor/Mwo+WZYDIRIlGbHcD9eZLdJ+PMlOjsGHy4EP//889+Hx7oo/swMTyT13LB+fUNDg4CpaqW8CfbrZb3zg3zfORK1ubFk1WXexGQoPFVFjWZIxI22EYDBYAj3cW95xVul9FBAdXW1UElqa91Vfnw458HcwKXzYCBEohTuSeL9yWYmDz796RKdHYP9oggAIDY29rihpuWVQmMD7+oeFBQkVJLa+uUSvTVa1hkWAyESq7mx3Aome0dXXuJvles0dNTKgw8+uGzfuSvVv5/Pa2yyPLX1xL333svOxNFaM+wu4afI78SJllj5MhDqrNkx3JJD5gaLgqmz087X0NIGOkKW23OgtmbOnFn05jszl7wQ763wdFXu19Xcfvc9b775ptDp+tOvV/j0EOIr7x1xWSpCEOoMjTskBJENBfxtLE1L+eUSnSXX7TlQux555JE77rjj4MGDVca67IJBS+6PcXVlKH98e55f0JuhJ0gQcn//SNTmxnIrLrI1TPjdBf6OHvhYoWsEBARMmDBhzm2zbh0axVSOLWmA/QY6XfZbw8v9/SNRuy2a21zE15iFTscf8sqo0QxpIQzV9xFT7ujBfX+BoYHt7y/wM6M4pgYXBIGBEImYvwrSNGTdVVZKlq/O8Qt6Ybcouq4xoaS4Hs5Ws9Io/OYcP68XRgEMhEjk5vZgpXfUwsOPF/l5PTEOouviCMyOIT+wkWNPVVF9I4zEDgwMhEjsZkZxe0p4PQPbdGwo4Hv5kh4+WKygG7kjlpXe0a/P8fN6EgVmWAyESOy8XWBWFPfVOSFLlqKioqNHj35xov4e7GVCN5OoJmYeDpUJ3CjkKfzfBTqvJ+ZYAAyESALujeM+O8MLUq4cPHhw2LBhQ+N6/GVs2sZZQUf+93RjY6MQCUGiQQDm9+I+PyNwo3B3CQ1UQX9/bA8CYCBEEpCqIS4c7Clxdii8cuXK+LFj5wWaDt6Xsfmu1Jx7Uo+t/PKhhx5ycjKQ6CyM4364yBsFne38zTkem4PN8INAUnBfHPf/Tju7iv3f//53Rozf7L6htkp1kIfrR5MG/vR/35WWljo5JUhctB6QoRVypLDGDKuv4ILXP+EHgaTg7l7cuqt8ZZNTb3ry5MmkUP+WV3xVLr0DPE6fPu3UdCAReqAP9z+nV92afX2OHxfGaT2Euj9zMBAiKQhUwS0R3HfOrWJ7e3uXN5haXSyrN3l7ezszGUiMJoSTGhMcMAgwtE0BPj7J/z0eC/8/4WeBJGJhHPepc6vYU6ZM+f54ocn65023XzaAT+CgQYOcmQwkRgTg3jhhGoVbi6iKw/2ProGBEEnE6FDSYHFqFXvu3LnxGRNnrMj96VTx9suGt7LOPbnz/Geff65QKJyWBiRe98ZxKy/zVa37FBzu45P8I9gcvBZ+HEgiCMB9cdx/Tjmvis1x3A8//qS/+7PtAUNXmEPcR88+fPzEhAkTnJYAJGrBbjAhjPvuvFMbhVdq6d4S/k6cL3ot2W+2iiTk/j5crx/NBXVchKeTun3WXuV7ZsxaO+1259wOScwDfblHs6wPxztvf9r/nuLn9+I8seC/FtYLkHT4q2BBb+79Y86rYn94nH+0Hz5EqIsytERBYEOBk/rzG62w/Cz/MPaLtoGfCJKUx/tzX53jK5yyjuJIBb1QA7dF40OEuogAPD+Ee/2w1Tm3++ECnxBEeuJ2uG3gM4wkJcyT3OKjn73opblz5z7yyCNbtmxx3L3+dZx/KJ5T4jOEuuG2aK7SBNuLHd4otFJ44wi/eADO5GoHPsRIUrZt27ZxQf/Y3K8yas+EnNzx4NxZf/vb3yi1fylz2Uh/vcLf3wefINQtHIFnB3Gv5zu8UfjteT7UA8aEYnOwHThmiqTDYrH8dcGCd9JjJsSqbVfmxIdN+f7H1dOmzZo1y773ejmPz+zHBars+6pIju7owb2cx+8rpSM0jopSVgpL8/lP07A52D6szyLpyM/PV9RWNEdBAPB0UdzVP3zt2rX2vdGZarqpkH8Ce5mQPbhw8NRA7s0jDmwUfnGWj/KCUVpsDrYPAyGSjurq6mCP1m20YE9VdXW1fW/03AF+8UCFj4t9XxXJ1197c4fL4XC5Q0YKTTy8ns+/PBTrbdeFgRBJR+/evU+X19abr6lZ5+mq+vTpY8e77DfQ/Qb6SF98dpDdqBTwjwHckoMOaRR+cZbv4wupDut3lQB8mJF0RERETJox64ktx2pNFtuV9edKV12qXLhwoR3v8vxB65IhnDsOryO7ejieu2SE1VfsvAq2yQpL8/lXErA5eCP4NCNJ+fzzzxcvXpzy1Ze9/T30dU2F3rFf/rQpOjraXq+/rZherYW/9cYaJLIzFw4+HqFYsMs6IYzzsF/B/Fq+NTGYJAZjc/BGiCNmlncBIYTneULk/m0ZjUYAwHN8mvE8bzAYNBpNp/6qoqLi5MmTISEhmxqjV1+FrZPtU67UW2DoKsubidzMKCEDYVVVlUqlcnd3FzANTDGZTEajMTAwUOiE2MHdO60RnrB0eLcacAaDITc3t0ePHo1BvW7ZTI/c6hKCmeWGsGKLJCggICAtLa1nz54P9VMaGmHVZft0Nz1zwDosmAgbBZG0vZOk+Pwsf6qqi+2Tmpqae++9t3dU+KsPzp+UMmzE8CH3u+zHKHhT2DWKpExB4P1kxd92W0dpuYDurfnbXkxXXaZHb8VHBjmQxh2WDFFkZlm3dakbY/78+fyp3OwF6V6uSgBYe7bkuYWTFiYf///t3X9Mk/kdB/DvQ1sKCLRUUGxPmMcoeJLMo2GbXM6oBHf+eKgm29yyi9G7LOqSxVvMgnOGeNFoNNNo3Kab2Hl/6A03WcRhxgEp42b4kXMjDNaCu6FC7fGzAm2xwPN890djg6WAQNvnoc/79YepfUr55OvH593v8+PbtLS0YFcaUfDZFiLcFi3zvTXMj8yT/CJOAoxMkA8/5669K0vCHfQQYofWRg2Pk1//Z96HMbq6uur+dv98YY43BQkhrD6V/Zq6tLQ02DVGGgQhRL4zebIXHPn4nwu/Nv2jBu47Oua9N6R+DhvCQMaQPxXITrdwNbb5fXbr7Ox8KyUhVv7K+cXcVHVHR0dQC4xACEKIfPIo8sct8j900soFfd/N76z851/RX30LF6BDmKxJYMq2yN+vm+wcnkfHqlSqPpf/F6/0uTwqlSqo1UUgBCFIwspYUlYg+7B+8suR+WWhqZM/3cJXbZPFYx0ZCKN3U5nTebKizzjHa3+nmMFg8MSpa7v6fc84xyc/be8xGo0hKTGCIAhBKjasYE4aZJsqueb+183CTx7xJx7yNdtkbybgoCiE2wf6qB1pzHdrJ4fHX+v1CoXinV+YDpm7Pq633nv0lanlybZPG977wd4dO3aEuNIlL8hByHGc33JWDoeDZVmNRlNUVORwOIL76wDm5cfZUb/Jj9pZNXnry7mvRPjkEf/LL/jq7bJMFVIQhHHum7KcJCb3L5NfDMzx6c09SX7WyLWlbvlXu2X1rg8+i17Tl1NgKv/rlStXwlPqkhbMG+ovXbp069at5ubmqe959OhRp9N5/vz5I0eOJCQknDlzJnAduKGeEIIb6qdZ2A31s/v3EDVWcz/MYErelikDnfhrd9CfN3P/HSF3C2Vr1aLrSdxQ7yeSbqgPqPwxf+gBd3y97KfrAk9dqm304D+4b69gLm6QpcQQQsjAwIBKpVIocED/tQQzCM1ms8vlYll26ntmZWXdvXs3OzvbarUajcaZrl9CEHohCP2EIggJIf0vyP6/Tzb10++/GfX+16M2rGB4SvpeELub/t7Klz/mj31D9pO3ohSiPHWAIPQT8UFICPnfKN1TyzEMMaZHbXuDeTuZoZR0jdLWIVr+mD7opb99Rzb1qmYE4bwEf4k1hnnlPePj4/v7+2NjY8fGxlauXDkyMjLTT+Xk5Ex/vrCw8NixY8GtUMycTichJD4+XuhCxILn+cHBwZSUlFC8eY+b+XOP8na3ou8F4+aYJAWfrKSbV0x+pPeoo0Wx9GBAw8PDSqUyJiZG6ELEYnx83Ol0ajQaoQsJrQmePBiQ1/bKa3oVQ+PEwzPLo+naRN6g4Q5keOJkr3Ts0NBQQkICgpAQolar5fI5VidYbBBmZ2d7J3m+9/ELwmXLlg0ODsbExLjd7pSUFJfLFbgOhmlqapo+I9RoNJJaEwEzQj8hmhH66RsjGiWRi3L+Nx1mhH6kMCP00+MmiQoyyzdiYkboI5fL5zzWuNj1oqxW6+wv0Gq13d3dmZmZNptNp9PN8sq8vDwcGvU2LtrXh+d5hUIR6gHRLanxVrwkdCFiQSmV2oCsmevOQDTJvITwM3BdXR0hhGVZk8lEKTWZTLidBQAAxCaEQbh582ZCSElJSWtr6+rVq9vb248fPx66XwcAALAAwQ9C3wlC7wO1Wl1ZWdnT01NRUYGVfuZkNpvNZrPQVYiI2+2+cOGC0FWIy507d9ra2oSuQkSePn1648YNoasQl9LSUrvdLnQVS8YSuTxAMhobGxsbG4WuQkTcbve1a9eErkJc7t+/b7FYhK5CRHp6em7fvi10FeJy8+bN3t5eoatYMhCEAAAgaQhCAACQNAQhAABIWvBXllmY1NTUrVu3Cl2F8Nrb2wkh69atE7oQsfB4PNXV1Tt37hS6EBFpampatWqVpBaamN3AwEBbW9umTZuELkREamtrDQaDWq0WuhDhnTp1as7/LGIJwtbW1paWFqGrAACAiMKybFJS0uyvEUsQAgAACALnCAEAQNIQhAAAIGkIQgAAkDQEocDy8/OZlw4ePOh73uFwsCyr0WiKioocDoeAFYZfeXl5Tk6OWq3euHFjZ2fn1E0zDVekmqkN0B5oDx/sQxYPQSgkSmlHR4fdbh8dHR0dHb148aJv09mzZ9PT0+12e1pa2rlz5wQsMsy6urr27dtnMpnsdjvLsvv37/dtmmW4ItVMbYD2QHt4YR8SHBSEY7PZEhMTc3Nz4+PjjUZjb2+vb5Ner7dYLJRSi8Wi1+uFqzHcampqDhw44H3c19e3fPly36ZZhitSzdQGaA+K9qCUYh8SJJgRCslutxsMhuvXrz958kSlUh0+fNi3yWazpaenE0K8n+mEqzHcCgoKrl69SgjhOK6kpGTPnj2+TbMMV6SaqQ3QHmgPL+xDgkPoJJacrKysgCP/7NmzpKQk31/j4uLGxsYopS6XKy4uLqwlht30Mamqqlq/fn1xcfHExETAH/Ebrkg1UxtIqj2mQ3sEJOV9yCLJhYtgibJarb7HDx8+9Hg8+fn5hJDo6GilUunbpNVqu7u7MzMzbTabTqcToNAwmjomlNLi4uKGhoaysjK9Xj/1ZbMMV6SaqQ0k1R5ToT38YB8SFDg0KiSXy7V7926LxTI+Pn7y5Mldu3YRQurq6gghLMuaTCZKqclkMhqNAhcaRvX19RUVFffu3dNqtU6n0+l0kpdjEnC4Itv0NkB7oD2mwj4kOIScjkoez/OXL1/OyMhITk7eu3fv8PAwpdT7j+JwOLZv367T6ViWff78udCVhs+JEyemt6j3z4DDFdmmtwHaA+0xFfYhQYG1RgEAQNJwaBQAACQNQQgAAJKGIAQAAElDEAIAgKQhCAEAQNIQhAAAIGkIQgAAkDQEIQAASBqCEAAAJA1BCAAAkoYgBAAASUMQAgCApCEIAQBA0hCEAAAgaf8HiusJ0cDZI+YAAAAASUVORK5CYII=" }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "using Interpolations\n", "using Plots\n", "gr(fmt=:png);\n", "\n", "x = -7:7 # x points, coase grid\n", "y = sin.(x) # corresponding y points\n", "\n", "xf = -7:0.1:7 # fine grid\n", "plot(xf, sin.(xf), label = \"sin function\")\n", "scatter!(x, y, label = \"sampled data\", markersize = 4)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "To implement linear and cubic spline interpolation" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "hide-output": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "li(0.3) = 0.25244129544236954\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdd3hb1d0H8HO1riTvIe8hS5YtWZLt2M7eJDQFsgij0KaUAgFaUiAhZSdAaYFAWqAF2pdCKFBKKTSBhrDJgMQ4iaeWZcvyXrIky0O25r33/UOpcRMn8ZB0rqTzefr0CRrnfu04/umeiVEUBRAEQRAkUjFgB0AQBEEQmFAhRBAEQSIaKoQIgiBIREOFEEEQBIloqBAiCIIgEQ0VQgRBECSioUKIIAiCRDRUCBEEQZCIhgohgiAIEtFQIUQQBEEiGl0KYUNDA+wItEAQBEEQsFPQi8fjgR2BXgiCIEkSdgoaoSjK6/XCTkEvXq8XbZ85fXQphKWlpeivDQAwPj4+Pj4OOwWNkCQ5ODgIOwW9jI6Oulwu2CloxOPxDA8Pw05BL0NDQ+jDwfTRpRAiCIIgCBSoECIIgiARDRVCBEEQJKKhQoggCIJENFQIEQRBkIiGCiGCIAgS0fxcCAmCkEqlkx+x2WwbNmxITEzcuHGjzWbz7+UQBEEQZI78WQhffPHFJUuWNDU1TX5w7969ubm5fX19OTk5zz77rB8vhyAIgiBzh/lxGfvRo0fHxsY2bNgwuc3CwsKPPvpIKpXq9fpNmzadUya/z4FhDz/8MIZh5zyuUCjWr1/vr4T0Z7fbAQDR0dHnPE5S5IDD0jPa123vi8WjF6SWRbH5MAIGG0mSFoslJSUFdhAaGR4e5nA4PB4PdhC6cLvddrs9MTERdhAasVqtsbGxbDYbdhD4+Hw+g3GJWz6WH6+3evXq8x/s6enJzc0FAPjuCy/ydq/Xe34hjLSNgnxf7OQv2U243zN89KHxk1hObGZ0WlZ0RoPF8krDG4UJ+UszFqzJXsFm+PMvkW6o/4IdhEbo/z1xEe4vOo91jHRZnYMWx6CH9M5PLV2euUgSL8LAuf/G547+35DgQ9+TCdP5JvjzjvBsi9j/tBkVFWW1Wrlc7vj4uEAgGBsbu9C7SJI8vxBGmtHRUQBATEyM7z9r+1W/P/1KQaL4V+W3JfISJl7m9DqremsOt3w55hl/Yvn9An4ynLiBR5Kk2WxOTU2FHYRGhoaGcByn5x0hSVFftB19veHvsuSCirTSZH5SCj+ZpMgT3aeOdpxwE+4bi7ZsLrjSvxd1u92jo6NJSUn+bTakWSyWuLg4dEc4TQG/mcjIyOjq6pJIJD09PZmZmYG+XNggKer5M38+3Vu7Y8GdizIqznmWy+Kuylm6MmfJP3UH7/jsvkeW7CxPK4GSE0Em6K2G5069xGNxn1j+QFFy4eSnChLFtxT/2GBr3fvdHzVm/a6Fd3FZOKycCHKOAC6fOHbsGABgw4YN+/fvpyhq//79mzZtCtzlQp1KpVq/fn1OTk5OTs769esf+/zp7tHev61/6fwqOAED2I1FW/Ys+/VT373wru5AMNMiyDnU5sYHj/3mx/JrX/rB3nOq4ARJgujldc8yGYxffv7rntGLDZQgSDAFsBD6hgz37NmjUqmys7O1Wu2jjz4auMuFtKqqqoUrLjuccMXQ4/qhx/WqMtlnqq838S7nsbiXfG9piuLVH/7hUMvnX7QdDUJUBDmfakC7+5unHl1y35rc5Rd/Jc7kPLT43k0FV9z1xQMGW2tw4iHIxfl/jHB2InyM8PLLL/8qdRNYdQcAIIX1jZD1j4ajsuWdn3/55ZfTbKFjuOuerx5+auWjF/owHqLQGOH56DZGWD+gefzbvbuX7ppR//zxzsqXa/e/esXv4/G4OQZAY4TnQ2OEM4J2lqGFqqoqMG8TAGCd871C5msa9yMu+Q1VVVXTbyE3LvvBxffs/vYZ87glYDER5FytQx2Pf7v38WX3z3SUemXOknV5qx/7Zq+XRCdRI5ChQkgLGIYBimJgblb0RztbxpYN9QKSmOn98aKMiuukGx85/pTTi05tRYLBQ3iePLnvjnk3l6YqZ/H2nxf/mM/mv1Tzmt+DIciMoEJIC0uXLgW1B0TM96Sj7t8n7Hmm7+XrdS8tXbp0pu3cILs6Ny7rL3V/C0BGBDnXqw1v58RmXSFaM7u3MzDs0aU7a/obDhunOwSAIIGACiEt/O53v0us/X0m9bHUlHfCm3J9r+w+Zs3zW2ez3Oruitu/6apssrb4PSSCTFbT33C048R9C345l0ai2PzfrXzk/+reHEBd+gg8aLIMLVCAuvPQfRm1BstX7f82WpctW/bUA/dFf/Yqr3hp3IZbZtraZ61HDjR9/Jcf7mNgIf9BB02WOR8dJsvY3WO3fHL3fQt+uTCjfO6t/U39bsdw92PLfj27twdusszzzz//1ltv+b3ZICAIgsFghNlv1NWrV//hD38IRMvhvDtXCPmy7RjFom718GL++fFTsYm+nWWIggLLXx6hPK74q+8EM/mBXida/Ynxy0Mtn2+SXBGwyEhEe/7MX1ZkL/ZLFQQA3Fh0zc8+vuv/Dr3eW93JYDCWLFly+eWX+6XlOerp6Vm9evVPf/pT2EEQcOzYsePHjweocVQI4bO7x/5S97cnhNdyko8wYr/fOJgZkyC461nLq7tt7/8p4bpfTb8WYgDbseAXO756ZEX2kgTuXOemI8g5Gq3NDQPadzb+xV8NusddlkNdfy3qrlXNpwgCvHrvlaXC999/n8+Hv7N8RkbGvHnzYKdAQEdHR+AKYch3nYWBA80fL8gozzK284rPnR3D4EcLfvk0Ye0ffPsZivBOv828uJwfitb8ue4NvyZFEAAA+Evdmz8vvhFncvzV4O7du7+u47mTKzI2lYGNu8Fj1Z90g9/85jf+ah9BLg4VQsicXueBpsM/kW1xaqp4ysXnvwDjcJNue5x0uaz7n6Q87um3/DPlDbX9qqZBNGsG8aeq3uoh59CsZ4pO6b333gObnjB6b81hfZBI9gIGC2x67L333vPjJRDkIlAhhOw/LZ/PS1Wm2kYZ/BiWYOpNyTE2J+mWRxncKMv/PUq5HNNsmcfi3li05W3Nv/wXFol0JEW9WvfWttKb/DgPi6KosaHBq3j9z3T884fmsbXE8wAAkCwcGBjw1yUQ5OJQIYTJQ3jea/zwxqItDnUlr3jJRV6JMVmJW3/NTs0xv/wgOTYyzfY35P+g0WowDrX7ISuCAPB52xGchS/NWuCX1sixkbFTX1hfe6L61pU39n90Okpehd0+HNXOxoZBt1ooFPrlKghySagQwvRp69cFCaKCRLFDVclTXmr5PIbFX3sXXlBq/tOviZHB6bTPYXKuk258R/O+H7IiEc9NuN9QvfvLslvmeLguMThgP/6h+aX7+397i7PxDL981cdZK7cePPN3fMGX0avLbVi+9x/gvV233367v5IjyMWhWaPQEBTxru7Ao0t3evraAeFlZ4kv/R4Mi1v/cwYvyvzH+5J/8TQrKe2S79gsueIG3e2dIz05segwSGRODrV8IUnIUwpks3u7p6/doap0qCuJIQtPsSjmsmtxSSnG5gAA7i1d0Wux/fHxEk/+ssR5ntT5n2/46bq7777br/ER5IJQIYTm6/Zv0qJT5MnSkc//wVVerF/0HDFrrmdwo8x/2pV85+/YabkXfzGXxb2mcP072vcfWnzv3PIiEY2kqA/0/9m9dNfM3kZRrjadU13pUFcCCvCUi+O33IkLiwDjf/qiGAzGvn377r777tOnT6eYmgc4+oqtS8JsMTgAwOVydXR0sFis7OxsdC4EraBCCM27ugPby28DADjUJ+OvvnNG741aehXGi7K88mDSbU9wcgou/uIthetv/Oj2Xnt/RvSl7yARZEqVPacTuPFFyZf4YfOhvB5Xc71DXenUVDFiE3nKxUk/383OFF38Xb5TqYkhC/OlX7wY99HVBVexGEx/ZIevtbV1z549Bw8dHuenAoqMcw/ecP21jz/+eFoa+idJC6gQwqE2N3pJoiyt2DtoIoaseJ58pi3wy1YxcJ7l1T1JNz+C519s7/8oNn9zwZX/0P5718K75hAZiWjv6/9znXTjxV9DOsedjWecqkqnvoadkcdTLo5Z+6PpdOBPxoxPzmcnZXKiv2o//kPRZXOITBfffPPNxh/9dHjt/eCplwAeBQAYHh/6v29f/2jBki8Of6RUzubgDr/AsOlusXnxV06/HdpCk2XgONTy+Yb8H2AAc6pO8hSLzukpmiaufGHSzx+x/u23DnXlxV95TeH6o50nRt32WYVFIprL5WqxtfbZ+1fkTLHOFQBAjo2Mn/nK8tfH+h/fOvbdZxyhLPWhVwW/ei561ZaZVkEfrrR8C5n5ru7fFAjtX68AgL6+vs033jR854dg5e2+KggAAPx4sO6+/htfX795i90ezv8q165dCzvCtKBCCMGo236y+9QPRWsAAA5V5fkbykwfLlYm3/Hboff/NH76YgfZxOGxizPmf9Z6ZNYXQiLN8PDwjh07BAIBNzruR7+9KbEvmvSSk1/gHTTZj39ofuWh/id/Pl7/Lb90edrjfxf88unolZuZk3YKnAVcWl5o7GYz2DV9DXP7IuDbu3evbfHtIFMxxXOSZZ15V7zyyitBDxU8X3/9NewI04IKIQSftx5dnDE/Fo8h7UOevna8oHQurXGyJcm/fGb4kzfHKj+5yMs2Stb9x/BZGHzERoLA6/WuXbv2hVNmy/1n2H/p5pZkvf5X1S233AIA8HQbRz592/TsL8wv7PCYOmMuuzb9d+8lb3uCP38tg+ufrUFxkcLb17ZeuCoMzik8ePAgWHDDBZ9e+KODBw9Ov7X9+/enp6cnJyf/8Y9/BAD85z//KS0tjY+PT09P37dvHwAAw7CHH344PT39N7/5zeOPP15QUBAXF/fUU0/5nnr22WdTUlKWL1/e3t4+0WZ/f/91112XkpIiEoluuummvr4+AMDg4OBPfvKTpKSk/Px837Umm/LZ88Ns3rwZAFBaWnr+U3SDCiEEHxu/WC/5AQDAoaniSssx1lznj7HTcgXbnxv9+v3RIx9c6DXFKXImg9lg0szxWkgkOHDgQLUVAze/BhIyM1ifmakVFTc8Ke7Xde3+sfXNpyivO+H6u9OfeCfh+ru50nKM6eepBhibwxHJl7tizvTVDbmG/dt4MHk8nq7efpCQdcFXJOcZjcbpN3jfffd98cUXlZWVhw4dAgDs3r1769atVqv1k08+eeSRR3yvUSqVX3zxxWOPPRYfH6/RaD744IMnn3zS95Tdbu/t7V2+fPmOHTsm2ty2bdu1117b0dFRW1srFou3bdsGALjnnnsAAG1tbQ0NDbW1tefEmPLZ88N8+OGHAID6+vopc9IKmiwTbGpzI0ESJSlyAIBDVcmf7589G1nJ6YK791n+/BDlGIu96mdTvmZj/g8/MnxWmgptcB4JFWfOnAGKHwIGk0uNF2Af7mpiEKTu0+iMU1mLrr0jGMv7uNIKT5NqSe6CL9uOX3KSDm2xWCwOi+EiPIB1gQ3KPc4ZnbCxYsWKhx9++Kabbvr0008BAHV1dWfOnHnzzTePHz/udp/diPiaa67xrc3Yvn07i8Vau3at0+n0PXXzzTezWKydO3cWFhZOtHn06NGPP/544j8FAgEA4JNPPtFqtbGxsQCAZ5555s0335wcY8pnpwzjc5GnaALdEQbbxDQZyuVwt2q4svn+apkZlyT41T6nvnro36+AqSZxrROtPtNXZ3MO+euKSLhisViAcAMA5hGHMhzYLVkvXCl+/k+tlCva/4ffTokrq3Dqa9aLLw/p3lEMw+aVloKOmgu+or26oqJi+g1++OGH27dvP3DgwBVXXAEAuP7661988UWBQPD0009PvIbD4fiWYLJYLF+GcxphMBgEQUz8Z0JCgtFopCiKoii73V5dXe17zeSv4vwWzn92yjCXfIomUCEMKrt77GT3qR+IVgMAHLrTHJHcX8MqPoyo2OS79np62wb/sQ+QxDnPRrH5K7IXf2L8yo9XRMLSypUrQd1HwOMUsL5lOeU9bAEY6mW2VS1btiw4AViCTIzJlpHRJEXqLM3BuWgg3HbbbeCLF6Z+jiLBl3+89dZbp99aXl5eXl7enj17ampqAABffvnlI488sn79+s8++wwA4PVe4qS2N954w+v1Pv/88ytWrJh4cMuWLU8//fT4+PjAwMCmTZt8terKK6/ctWvX6Ojo2NjYQw89dE47Uz57oTAej2emOYMPFcKg+rzt6KKMing8DgDgnNt80QthcPnJd/6OHBuxvvE7yus559lNkisOtXxOhviiHyTQ1q1bt3F+Afv1TQ62yWAvAw2HwfNX3L/99tzcS+xk5EdcabmzsfoK0ZrDxi+CdlG/u/nmm9ckjID/PHluJw3hAf+45ycLhb57u2m6//77Fy9evGrVqueeew4A8NRTT61atUqpVFqt1nXr1t18880Xf7vH40lLSzty5MiLL7448eBvf/tbgiDy8vLkcrlQKPRNZnn++edJkhQKhcXFxStXrjynnSmfnTLMlVdeKRaLZ5oTAooeAAAkScJOEXB3frbrVG8NRVGkx93z4DXEqO2cF4yMjIyMjMz9QqTXY/3bUwMvP0g6x895atunO0711Mz9EsFBEER/fz/sFPRis9nGx8/9a/U7r9d776u7Hnr+yoK83DVr1hw8eDDQVzzHuPo78ysPDTps6//143GP4yKvdLlcFoslEBnuu+++5557bo6NjIyMbN68GQjLwdaXwAPHwP1HwI/2gQzZrbfe6nQ6/ZJzOujz2352Dh48uGnTpgA1ju4Ig6fX3t9nN5WnlQAAXM117Iw8RnR8gK6FMVmJP32AlZhq/vPD5Pj/rNi9Snz5521HA3RdJGwwmcyh5JHLcWFTa/tXX33lmwofTFxJiau9MZ7BK0lVHOn4NshX96OYmJiDBw9+9dre27j1y04+uvy7x7YL2qsOvPHaa6/hOA47HQIAmjUaTEfav12Vs5SJMQEADtXJQPSL/g8GI+FH9wwfet38x53Jv3iaGXd2msOqnKWv1r/l9Dq5LG5gAyChzDjUPuocnpe3HlYADOdxsiWuFtVV4sv/rv3XVeLLYSXxizVr1qxZ458p4rPz7rvvQrw6zaE7wuD5quObNcIVAABAkg7NKa5iUcAviWFxG2/jz19r/uMur6XP91gcHqtIlp3oPhXwqyOh7LPWIyvtHF5+CcQMXGm5U1+9IGNez2h//xg6sH5Obrjhwuv6Ix4qhEHSOtQx7hlXCKQAAFeblhmfNLttGGchZs310auuNr/8gNfU5XtkjXDF1+3fBOfqSCgiKOLrtmOLukY5QinEGFzZfGdjNRNjLs9edKzzJMQkSHhDhTBIvm4/viZ3he9ob4fqJK84SNPQfaKXb4xb/3PzS/e7O5sAACuyF6nNjcOukWBmQELI6d66VGa0MK1w7tsezQU7I49yO73WvtU5y0J6mBChOVQIg4EC1JGOE5cJl/v+06mp4hXP4CRev+CXr0644R7ra4+7WrVcFndBRhn6iI1cyOdtR1YTibgEZr8oAABgGF5Y5tTXlKYqzOPWXns/5Dyw+Vavh9+RxdChQhgMjZZmFoMlSRABADzdLYDBvOTJ8oHAlS9K3Hq/df+TTn2Nkifdf/ztrVu3Pvjggw0NIb/HP+JHTq/rdG9tWccgF+oAoQ9XWu5qrGZgjJXZS450nIAdhxYeeOAB2BHCDSqEwfBV+zdrhWeXnTpUJ4N/OzgBL5iXfNtjfa//9q1tO80Mx78Tlu3Vc8tX/XDv3r2wIiF0c6q3Rpog4g/0s3OmdR59QHELy1wtaorwXpa77CjqHQUAAPDMM8/AjhBuUCEMOJKijndVXpZ7dlAwGAsnLmo0JmXL+6cfWlEaTxQKypPA+keIh04+/Nvn1Go1xFQIfRzvqlzKzuTkyf1+rMQsMKJiWanZ7jadQlA05BrpGO6CnQi+iQ7Sf/7znyUlJUlJSS+88AK4wIFKUx7V9M477/g210Z84P+ghz21WZfIjc+OzQQAeM09pGOMk1N4yXcFztdff10fLb9W8oen+h58O6evC2wGidlkxbWHDh1SKtHBFJHOQ3hO99behM2HP0D4X74NuPH84tU5S491nvyZEuYyAIrwDn/01+BcC2Ox4zbedpEXdHZ21tfXHz16dP369ffee++2bdu2bt361ltvuVyuF198cdu2bR9//PHu3bt/+tOf7tixQ6VSLVq0aNeuXQCAkydPHj2KdtX4HiqEAXeiq2pZ1tklgw5VJU+5GEAd67bZbCAurYuT+oDg2XLwizxC3cZUgvj0oSF0KgUCqvsb8uJzeTV6/CdXwc5yFldabnv/T3Hrf746d9lzVS/BLYQAAFZyRnAudMk78l/84hcYhl122WUOhwNc4EClKY9AeuKJJ9Ad4WSoEAbcie5TT608exalQ3Uy7qqbocYBhYWFoP1FQJG9bMEVYykyxhdtTCVoO1NwVbD30EJo6HjnyZWp84iRf3Ey82FnOYuTU0gMWYhha1FyocPrbBvuzIvLgRUGY7KiV2yCdfVzxMTETP7PhIQElUolEokAAGNjY1arFQBw/fXXczicG2+88emnn37rrbd8r0RV8BxojDCwjEPtAIC8+FwAADFs9Zp7OGIF3EjLly+fl8oFBx4FhGfUO4/N0YBj/5c+UHvdddfBDYZAR1BEZc+Z+U4eLlYCBm1+OTAYeEGpq6kWA9iqnKXHOtCyn6lNeaAS/Y9AogPa/KyHqW86v1uRvdj3Z4e6kitfCH0CApPJPHTo0HqWHjyYr/7nRy7myOK+d7/88su4uDi4wRDo6vrVmTFpse1tOA0WTkzGlZY79TUAgGXZCyt7TsOOQ1NTHqgUAkcg0QDqGg2sE92n7p1/u+/PDtXJ6OUb4ebxyczMPHToUEdHR4vB8GXzS0te3CWXy2GHQuA73lW5InuJq/owTX5QJ3Bl84c/eg2QpDxZOjBmMY0NpEalwA4FAfXf05Qm/n/y4zExMfv37z/nLXfdddddd93l+/Ovf/1rAMDf//734KQNIeiOMIBMYwNWh1WeLAUAkA67p9PALSyDHep7ubm5a9auXczN+bb1OOwsCHwkRZ3sPrUsoYgcH4Wy4cNFMGMTmXFJ7q5mBsZYlFn+XU817ERIWEGFMIBOdJ9akrmAgTEAAE5NFV5QgnFod/LRopxFGnunw+uEHQSBTG3WJXATkrq7cUkJ3InNU/ItogAALM1aeBL1jiJ+hQphAH3bVbUse6Hvzw5VJU8Jcx39hSTklxU4Wad6a2AHQSCr7D69PHuhy9BAnxWEk3GlZwvh/PR5WrN+3OOAnQgJH6gQBsqIa9Rga61IKwUAUG6Xy9DAlS+AHWoK7Oz8MovrBJqJF/Eqe84szpzvalFxJaWws0yBI1Z4+zrI8VEei6sQyE731cJOhIQPVAgD5WTP6Yq0Ug6TAwBw6ms4uYUMfswl3xV8GJO1KEbyXU+NlyRgZ0Gg6bOb7J4xEcWnvB5WShbsOFPAmCyOWOFqqgMALM1aUNlzBnYiJHygQhgoJ7tPL836b7+o+iQX3kbbl5SeV5qJcRsGNLCDINB811O9MKPcbVDjtLwd9JlYRLEkc35VTzVJkbATIWECFcKA8BCe2v6GRRnlAACK8Dq1p3nyRbBDXRAuUpSPME90n4IdBIGmqrd6UUY5bQcIfbiyCqe+GlCUgJ+cFpWiMTfCTkQ76MDC2UGFMCBUZl1efE4sHgMAcBvVLEEmMz4ZdqgL4ghlpd2DVd2orylCOb0ujblxfvo8V0sDLimGHeeCWMkZGBv39HcAAJZkzT/ZjeaOTg0dWDhTqBAGRFVvzcKMct+f4R5AOB0YB89LEHq9ru7RXthZEAhqTQ2FSfn4oBUwmaykdNhxLoYrLXc2VgMAlmYt/La7CnYcmgrXAwu9Xu+TTz55zTXX/PznP3/nnXcmbykwR2hnmYA41VvzyJIdAABAUQ71d4Ltz8JOdAm4SDHPqz3VW5NVGKSd9RH6qOqpWZRR4WxpoMOR9BfHlVWMHj8Yc9m1+Ql5HsLTNdLjO+AsaLwk8e+mQ8G5FovBuqZw/eRH9u/f/8gjj3g8nj179tx9990Yhu3du3ffvn2FhYVvv/22UCj0vQzDMIqiMAx79913n3766e7u7t27d9977739/f2/+tWvjh8/Hh0dvWzZsr1796an0/pzzzlOnTp92J4Bim4AzpG/7fnTG2+88fHHH3O5flicje4I/a/Pbhp12wsSxQAAd0cTgxfNEgT13+os4GKFctBzurcOdhAEgqremsWZ812GBjrPlPHBJSWezmbK5cAAtiCj7HRfsH9iKUBZHIPB+d+gY/Ccq993331ffPFFZWXloUNni7Hdbu/t7V2+fPmOHTvOT+s7sPD9999/+OGHAQDbtm279tprOzo6amtrxWLxtm3bAv3t8i9PYi7Y+jIo2wyW3AQePP51D/HKK6/4pWV0R+h/Vb3VC9PLMIABABxquveL+nDyFLJ3ev8ca3ERbpzJgR0HCR7jUDsDw3JiMvqMmvjNd8COcwkYh8vOLnC1qLjyhfPT533W+vU590yBxmaw7iq7JZhXnGzFihUPP/zwTTfd9Omnn/oeufnmm1ks1s6dOwsLpzjuezoHFoaSuEn3rxgDrP7lxx//eefOnXNvGN0R+t/p3tpJA4SVXGUIFEIGPzomLkXET0OLKCJNVU/14sz5nv4OjMun85SuCVxZhW+YsCKtVDWg8xAe2ImC58MPP9y+ffuBAweuuOKKyY8zGAyCmGId8PkHFhqNRoqiKIqy2+3V1aG2ZSuGAQC2eu6WeaoBACA2ZWRkxC8No0LoZx7C0zCgLU8vAQB4+jsowsPJossBpxeHixXzqPgzqHc0wnzXU704s8LVXE/nhROTcaXlzsYzAIBoTpQwLkcdSYso8vLy8vLy9uzZU1Nzdk/EN954w+v1Pv/88ytWrLjk26c8sDCUeJxs4DFH9ZiYWQAA0FIpk8n80jAqhH5WN6AWJ+TFcmLAxP6iIbKmhyNSFFucp9DOVZFk1BmPjK4AACAASURBVG1vHWovTVG4Wmi9gnAydkYeRRBeSx8AYGFGWUTttXb//fcvXrx41apVzz33nO8Rj8eTlpZ25MiRF1988ZJvn/LAwlDS16gYPhztYQ0y0kDjEfbnz91zzz1+aRgVQj87FVILJybD84szWtvsbnuf3QQ7CxIktf0qhUDGwVguowbPp+8KwnNwC8uc+moAwPz0stO9EVQIf/nLX1qtVpPJdOutt/oeeeaZZywWy8mTJ/Py8sD0Diw0mUxms/m1116LioqC8DXMgTA1UdL0B26fAzxSlP+f7Yc+eLeiosIvLaNC6GdVPTW+DWWIwQHCNoDnFcFONF3M2EQmhz8vXnIm6DPxEFiq++sr0kvdPS3MuCRmTALsONPFlZY7G2sAANKkfLPDanXaYCdCgkEoFEqWSVflzzec/FSv169bt85fLaNC6E+99n6H1yFOEALffFHlYsBgwg41A7hYUeaNCf6UdASW6r76+WmlLoMqVPpFffDCMrdRRXk9DIxRkVZa3V8POxEc7777LuwIwdZE2tYsuTo/P5/J9OevVlQI/em7nupFGRX/XTjxHS8U5otOxhEr5CZ7nUmFTqKIBH12k5NwCeNzXIYGnPZL6Sdj8KNZ6UJ3qxYAsCB9Xk2kFsIbbrgBdoSg8rgd4wxKmOP/1a6oEPpTdV99RXopAIAcG/H0tOIF82AnmhlcpMCNjdkxmRqzDnYWJOCq++sr0kowknS36XCxAnacmZk0TDivdkBN+m+3LYS27ONDBVic707Dv1Ah9BsvSTQMaMrTSgAADvV3uLQMY4fYynSWIBNQVEVCwZlI/YgdUar76ivSSt0dTczkdEZULOw4MzOxmjCZn5TIjTeOtMFOhATcmNuuiBcFomVUCP1GZ9FnxWTE4bEAAIe6MuT6RX04IrnCxa3tV8EOggQWSZG1JlV5WomzpYEbOvNFJ3ByComRQWLYCgCoSC2tM6ONIMLfGOkuzpkfiJZRIfSb6v4G3+0g5XK4jRpu0QLYiWYDz5OL+gbbhzvHPOOwsyAB1DxoTOImJPOTXIYGvIDuW4xOAcO4hfMmekcbrFrYgZCAc2GkXLoyEC2jQug3Z/rq5qfPAwA4G6s5eUUMLh92otnAxUqyVSdLKkB7rYU333g25fW4O5o4ohAbIPThSs/2jsqTCttHO8c9DtiJkMDiUAw+Hh2IltGm2/5hd4+1DXXIBVIAgEN1gleyFHaiWWJn5JGjg2WJK2r7VUsyQ/KmFpmOM/31NxZd7W5vZKfnhuqHNmnF0IevApLkMDmiWKHarJvYy8KP1Gr1Bx984PdmkZk6c+ZMFBMPUOOoEPpHnUmlFBThTA7l9Tj1NXFX3wk70WxhGEdYpHRzX7Ccgh0FCRSn19U82FKSonB9+X5oLZyYjBmbwExIcXfoQWa+MlFW26/yeyG8++67d+7c+c9//tO/zQaB2+1ms9lYiOzvOB2GXs38eYEazEaF0D+q+xt8CydczfXsNGEIbdJxPlysyO6zWjyDgw5bIi+EvxDkQhoGNAUJYh6LazfUx6z7Cew4s8eVljubarmZ+cVJRX8z+L9c5eTkhOjtoMViiYuLY7PZsIP4CUlufvvq31y1N0DNozFC/5hYQRha+4tOCRcpvEZNcUpRHRomDFNnBwjdLndPawjtAng+rrTcN0woic/rGe0bdvnnUB6EbrrbakgMy0jOC1D7qBD6gWlsYMwzLorPBRTl1J7iKhfDTjQn7JwCz0D3vCRZbX8D7CxIQNSa1OVpJa5WDSdTjHG4sOPMHkck95q6yLERJsZUCorqTeijW3hqMFYWsZIC1z4qhH5wpq++Iq0UA5irTcuITWQlpV/6PTSGMVmcbInSy0erCcPSiHu0z95fkJjvalGF5MKJSTAmC89Xug31AICytOI6kxp2IiQgNAM6eVJh4NpHhdAPqvvrfSfxOlSVod4v6oOLlek9JjfhRkcyhZ8Gk0YpKGIxmCF0GO9F4NIKT3MdAKAsrbgG9WGEKb3bXCwKYE8bKoRzRVJUXb/at5Teqf6OVxyqCycmw8UKV6umNFVZa0I3heGmzqQpTVVQLod3oIuTK4UdZ664sgpXUw2gKHF83ohr1DJuhZ0I8bOxwd4eNiHLC+BqLlQI58o41BbHjU3hJ3t6jADD2OlC2In8gCOUubtbygRy9BE7/NSb1PNSla6WBnZOIcYK+VmFrKQ0DOeTpk4GhpWmKtBHt/CjbjwqBFEcVqAWEQJUCOeutl81L1UJwqhfFACAcbjs1JxiIrq2X0UBtK9/+Bh2jZjGzZJEkdOg4oZ+v6gPp7CMbFUBAMrSitHAdvjR9NTLYnICeglUCOeqzqQqSy0GZxdOhEO/qA8uVsR3d/NY3I7hbthZEL+pN2mUAhkTY7oMDbgktGfKTMALywijCgAwL7UY3RGGH91oZ3FWWUAv4bdCaLPZNmzYkJiYuHHjRpvNNvmpJUuWYP91550hu+XKVAiKUJsbS1MVXksvOTYcBiMuEzgihatVU5KqQJuOhpP6AXVpqpIcHyWs/ezsfNhx/IMtUpK9raRzPCc2k6SoXns/7ESI35BuRzPLUVywIqBX8Vsh3Lt3b25ubl9fX05OzrPPPjvxOEVRTU1NfX19o6Ojo6OjL7zwgr+uSAdN1pa0qNQ4PNahquQql4Aw2tAIFyvc7fpSgbwBrc0KI3UmzbwUpcvQwBHJMWaYbCyFcXBmlsTV4rspVKL1r+GkVX+Si7EEsWkBvYrfCuHBgwe3b9+O4/j27dsPHDgw8XhfX5/X673qqqvS09O3bt06MhJWWz/U9qvK0v7bLxqaBxBeCIMfw4xPllMx9eiOMFwMuYbN4xZJoshlaAiDhROTMcTFLn01AGBeqrLehI5kCh+qtio5N+Ars/32kbCnpyc3NxcA4LsvnHi8r6+vvLz8D3/4Q05Ozo4dO+6555533313yhbuuOOO87eILSsru+GGG/wV0u9O99RuyFtn6+7wmLqcKULX8PAcG7Tb7QAAkiT9kW6usKwCntHIAAx9b3N6VCqUDCRJjoyMcLkhvPuJ342MjOA47na7Z/rGyr7TsviC0ZFRR3Mdd8v24Tn/uNKE2+12pImJj/+CXT4s5ue+1v/3sPnSZs13yxEGe42qB5slmaVz+QuNjo5mMpkXf43fCiFFUb4yRlEUQRATj5eXlx85csT352eeeUYul1+ohcLCKTYOyMzMpO3fpYf0Ng8ZS1MVoPZbdkEZh8ube5u+L5YmXzIpknu0VcoCWeOwISc+C04GkmSz2TT5htAEm81msViz+J5obU2lKXKma4wcteFZ+YARJnPlKIpipeYAimIMm7OTMxkYZnXb0qJSYOeCif1fsIPMDUU1k8NX5y+fyxcynSM4/FYIMzIyurq6JBJJT09PZmbmxOM1NTUul2vJkiUAAA6Hg+MXXAuyc+fO0Do0pH5AkxefI4hLtjRVxyz+IY/vh0PdfJ8h+P5oau7wogrT4TfKV97cYG68WnYVlAwkSfJ4PJp8Q2jC7XbjOM7jzfiDl8bauFl6JaOtmZtfwo8OyAGnULBYLIIgGNJyrF3Lz8kvTpU3jRpFAiHsXDCNj4/z+fxQL4S2rkYrDuR5ZUzsErd0c+S3j4QbNmzYv38/RVH79+/ftGkTAODYsWMAgLGxsauvvrqxsdHtdj/55JObN2/21xWhq+tXl6UWkw67u13PlVbAjuN/zLgkBs5TMJPr0RaOoc/mHLY6bOJ4YfgNEPpwZeVOfTUAoDRFgWZ4hQd103EJIy7QVRD4sRDu2bNHpVJlZ2drtdpHH30UALB69WoAwPLly3fv3r1hw4bMzEybzbZ3b6AOlAq+WpOqNFXp1J7GJcUY7od+URriiOQCk5mgiP6xAdhZkDmpN6lLUuQMjBGuhRAvmOdu1VIed2mqAs3wCg/qfo0sXhyEC/mtazQ+Pv7w4cOTH6EoCgCAYdj27du3b9/urwvRhNPrMgwalQLZ2KfPhdl80clwscJlVJfkKhpMmjTRZbDjILPXMKAtSZETQ2bSOc5Oy4Udx/8YvGhWep67TZtdUOr0ukxj5tQoAexQyJw0OnpvFAZjUCZMRsuDT23WSRLFOMVwNddz5QthxwkUXKx0GzXoI3YYUJt1xSlyV3M9V1ISTgteJ+NKy536GgxgJSnyhgG0iCK0eYatLbhXKQrGdl2oEM5SnUldllrs1NewsyWMqFjYcQKFJcikvB4FNwMdeRrSxjzjfXaTJFHkbFGFZb+oD1dW4TuwviRFoUKFMMQ1648mMbix3GD8dkWFcJbqTKp5aUqHOnw22r4QTp48bcDq8DrM4xbYWZBZUg3oZMkFTIzpagnPAUIfTraEGB0ibAOlaGvA0KfurJXxMy/9On9AhXA2HF5n21BnUYLYqT3FUwbwuEg6wMUKt1GjFBShvqbQpRrQFgvkXksfIAiWIEi/XCDAMG7hPGdTbV58zpBrxOIYhB0ImT3tSLsyI0gf2lAhnA21WVeYJKHamlhJ6cz4MB+Qx0UKd6umNFWJekdDl8qsLUkpchnq8YIwOXHiQiYPE6Le0dBFedx6zF4sCexe2xNQIZyNBpOmNEXhUJ0M+35RAAA7U0SMDCqic1RmHewsyGy4CLfR1i5NKnAZGvD8sO0X9cGl5a7mekASxWi+TCgbaKkeZzOFSXnBuRwqhLNRZ9KUpsgdmqowXjjxPQzj5MqyBketjsEhV6Rv4RiKGi1N4gQhl4W7jOowHiD0YcYkMBNT3B16NHE0pNUbT8g4SRgI0vRmVAhnzOl1tg615zsZDC6PlZoNO04w4GKF16gtSirUmPWwsyAz1jCgUwqKPKZOjMVmJcLZPD2YuNIKp75GkiAyj1uGXWF13E3k0FmaZUlT7D4dIKgQzpjGrJckiknNGW4k3A4CAM4uq9coBDKNuRF2FmTGVGZtSYrCZVCFfb+oj28RBQNjFCUXqgZQf34IoqhGj7lUHIwVhD6oEM5Y/YCmNEURCQsnJrCzCzwDXcoEMfq1EnIIimi0NMsFheG6s9r5OEKZ19xL2oeV6KNbaBrva23ngaLssqBdERXCGas3aRTsFMrl4GRJYGcJEozF5mTni0dJ41C70+uCHQeZAcNga1pUSiw72mUM56X0k2FMFp5f7GyqLRYUoRleoUjX9E0mFsVjBe8UUlQIZ8ZFuFtsraLufl7x0nDdp2pKuEiBtetF8TlNgy2wsyAzoBrQlqTKPb1tDH4MMy4Jdpwg8S2ikCUXtg51oI9uIUfd2yCPEwbziqgQzozWrBcnCCn1aV5x8Pqv6YAjUriMaqWgSI16R0OKyqwrFshdhnquJMxXEE6Gyypc+hqcwRbF5+qtzbDjIDPTON6tzCoP5hVRIZyZhgFtSazYO2jCRXLYWYIKzytyd7UoEiVqNOgSOihAqQcalSlFTkOk9Iv6sBJTMV6Up6cV9Y6GHGLU1sx2FwdxpgxAhXCm6gc0hSNenmIRYAT8rEhawXAeOyVb6sa1Fj1JUbDjINPSOdzNY3OT8Xh3m4YjVsCOE1RcWYVTX60QyNCan9DS1fwdYLHSY9KCeVFUCGfAQ3iarIY8Q1vkzBedDBcreJ1t8dy4tuEO2FmQaVGbG5WCInd3CzMumRmTADtOUHGl5U59dXFKkcbcSFIk7DjIdDW0n5JzM4J8UVQIZ0Bnbc6LzWR1teEFwZvXSx8cscJl1CgFMjRMGCo0Fr1CII2chROT4fklnm5jDMVK4iW2DqGPbiGjcchYlBrsgSdUCGeg3qQpouJwaRnG5sDOAgEuUrjbdMpkGRomDBVac6NCIHMZ6vFIminjg7E5HKHMZagvTilSo2HCEEF53HpqpCR/WZCviwrhDDQMaAr6RyJtvugERlQsMz5JRkWjLRxDwrBrxOqwCaMz3O2NuFgJOw4EuLTcqa9RCNBHt5Ax0q7p5TIKU6RBvi4qhNNFUITe2iw2tnNlFbCzQMMRKQR9AwRFmMbMsLMgl6A2N8oFUqLLwBJkMvjRsONAwJVWuBqri9FRmqFDZfhWxIxnM9lBvi4qhNPVPGhMZcbGZ8sYvEj8neKDi5Uuo0aRLEU7V9GfxtyoSJZF2sKJydhpORRFpY57SYo0jQ3AjoNcmnZAp0jMD/51USGcrnqTRjbOiNh+UR9crHS3ahTJUo0FTUmnO7VZpzw7QBihhRD8d4sZtJowNFBUo8tUnLco+FdGhXC6VCatpMfMVSyGHQQmZnwyxuZIWUlatDaL3jyEp8XWLo0TujubcVFkrSCcbGI1oXoA9WHQnWegq4UHlDkQxp5QIZwWkqLUJo0iKocZG1mLsc7HEStyrfaOkW6n1wk7C3JBTYPGnNhMZncbO12I4TzYcaDBC+a523TKBAnqzKe/lqbKaCYngRsf/EujQjgtbcMdsSQjTbEcdhD4cJGCbGvMTxA2Wg2wsyAXpDbrlIIiV0tDRG0xej4Gl8/OEGXb7P1jA2OecdhxkItR9dQoonOgXBoVwmlp6NdIbW5exJzEexG4SOE2qhXJ6KQ3WtOYGxUCqbM5EpfSnwOXlXv1dZIEkc7SBDsLcjG6kQ5F5jwol0aFcFrqO76TUTGs5HTYQeBjpWSRLqeMn6FF82XoigKU1qKXx4s9va0coQx2HMh8e60pBFK06SidkWMjTRxXiQjOzQYqhNOithrKRBE9TeZ7GIaLFQWjhMaMdt+mqe6RXg6TE9fXx8nKxzg47DiQcbIkhH1ExktHH93ozGyoHmZjwoRcKFdnQblqqGhra/vggw+MljYgdcbkzocdhy5wkZzZ3hrDie4a6c6Ny4YdBzmXxqJXCGSRucXoFDCMKy3LN9t1liaSIhkY+vRPR6q276ScVFh/O+hn4oJeffXVonkL7v+86wyDJRympKuuOHz4MOxQtICLla5WjTxZqkWDLrTkW0of4SsIJ+NKK9gGbRIvoW24E3YWZGpaq0Ee9J3VJqBCOLXm5ua7dj7gvO8IuOH3uXmUHSu23fTWT2/62fDwMOxo8LEzRMSwRR4nRH1N9KQ16+XxQo+pi5ML7TcLrXCl5a7menlyIZrhRU+U19NI2kqDexjvZKgQTu3DDz/0ll0D0qUAAILdrWJdBmSX2VKUX3/9NexoNMBgcHKlBeMM9GuFhuzuMdO4OcsyjAtlGCvYezbSEyM6jpWcLgVxqA+DnhydTe18IEuDtvMDKoRTM5vNICkHAJBNtBAMohpfBgAAiTlmM9psGgAAcLEyo2/APG4dcY/CzoL8j0Zrc0GC2GPQoH7RybjSCrHFjj660ZPecCKdGRPF5sMKgArh1EQiEehqAADIiaMJznjS943qahCJRJCT0QNHpPAYtdIkCVqbRTdai14ukLoM9Xg+KoTfw2XlgmbDiGvU5kSjG7Sj6lMXxeVBDIAK4dSuv/76xPZvQPUHiZjW5REDigKfPSeNcq1atQp2NFrg5BR4TJ2KxHytGRVCetGY9fKYHMJqYmdD2MWftnChjLT2yhJE6KaQdihK5+wrEcKclo8K4dSSkpI+Ofyx7JsnXYyO9m8bwGMlS/oOHzp0iM1Ggy4AAICx2Jys/AIvT2NBv1ZohAKU3moQD7s5IjnGRIujJmEw8fySQiIKzfCiG6+lt5lPFucshJgBFcILWrhw4Tcnv7LzGDs3b6n5+J0TJ07k56OP2N/DxQqxeVRvNRAUATsLclb7cFcsHsNtNaABwvNxZRUSsx3tL0M33YYqwGSmR6dCzIAK4cVoDN9KSP4NN24tKyvDMAx2HHrhiBSs1qYUfnLbUAfsLMhZWrNeIZC5Whq4qBCehystz2luMdhaPYQHdhbkew2d1XJ+FtwMqBBeTEN3jSIazpY/9IfnFbm7DEWJEjQlnT60Fr0sOpsYGWRnoFld52ImpPB4cVl4ksHWCjsL8j3tUJsivRhuBlQIL0Zj7yzNKoOdgqYwnMdOyS5kxKNCSB86S5NklMDziwED/dOeAldWISWjNGiYkDbIcXsTywFrr+0J6F/LBTk9zg7MoZReBjsIfXHEivwhF1pBQRN299jAuCWtqxctnLgQrrRcZB5FP7H0MWys7+dikmQJ3BioEF6Qqul4rofFjxfADkJfuEiR0tk96LANu0ZgZ0GA1qKXJkm8aK/tC8PFSnGvWTuApjrThdpYKeYksRmQZzijQnhBDR2n5Nw02CloDRcrvG2N0sR8naUZdhYEaC16eXQu6XSw09DA9tQwNicro8jtcZrHLbCzIAD4TkpJKoSdAhXCC1PZDMVpkIdwaY4RFcuITSzipaG+JjrQmPUSB+BKSgCa4XxhXGm5hOCjgW06oAhvo9danLcIdhBUCC+AoIgWaqRUugp2ELrDxQrxGECLlKEjKUpvNeT2WlG/6MVxpRViCxompAVXl8EYjSky4f/EokI4tabOuiQPSEwvgB2E7nCxIq/PqrcaSIqEnSWitQ13JPES8BYdKoQXx0rNLnBzNH0NsIMgoLmlMoHBi8fjYAdBhfAC6gzfyJnJqIvpknBxMd6iT+QltA93wc4S0RotTbKYbEAQLEEm7Cx0V5RTbhzu8pBe2EEinbanXhYrhJ0CAFQIL0Rl0RcnodvBS2PGJ2Nsjiw6W4d6R6HSWpryXRx0Ozgd8bKFaV52C1pWD5tuvLc4mxYLtVEhnAIFqEbPQGn+cthBQgNHJJd4cDT7AC6dpUk8MIxLSmEHCQHcgnn5Q05NvxZ2kIjmtfQ18QilEOZe2xNQIZxCp7Wd4yUzxBWwg4QGXKQQW+xo9gFEY57xgXFLaksLuiOcDgznybipmo7TsINEtIGWM2NsRm5sNuwgAKBCOKW6piMyMgZjc2AHCQ24WJlubDePW0fddthZIpTWoi+IzmKxcFYizC38Q4gyd4F2CHWNwlTfcVrGS2fQYx4GKoRTUPWpFFCPSw4trJQs4HJIYnMarWhZPRw6S1MByUc7q01fnnyV0+tq7NCTJJrtDId2qFWZpoSd4ixUCKegGe8pES6AnSJ0YBgukheCWNQ7CovO0iyyjuMFqBBOy9jY2J4//TVn2L3ljuuiY2Jvv/12q9UKO1RkIR12PcMOfa/tCagQnss6PjhGufMLlsEOEkpwsSJ/xIPmy0BBAarR2pzb2onuCKfpuuuu2/tlE+XNkVy9yfF061/1xOWXX+71otUUwTPaqu7iMaTJ8DdX80GF8Fz1hm+kLg4rNhF2kFDCESmEHX2NlmYKULCzRJyukV4+g52IxzHjkmBnCQFVVVWfntaBbW93MRfijFbATwA/faXOxvzoo49gR4sg2pbKXHY8l4XDDnIWKoTnaug6UwT7uOSQw8kUR1ut0Wxe10gv7CwRp9HSVAji0HzRaaqvrwcFywGbewr/wTBuZ1MugDFA0dq6ujrY0SKI1tIoT4R89NJkqBCeSz3cVpyBFmPNEIPBEUoLOQI0TBh8WkuTeNjLRSsIpwfHceCyAwAszJQ4N7vEcxIAAJyjOE6Xu5OwRxFendtcSpsBQoAK4TkcXmcvOSYvXAE7SOjhiBWScSYqhMGnszTldfdzxArYQULDqlWrmIZvga0HAMD2pIrI08AxDBoOrVmzBna0SOHuMRqiMSWd7jdQIfwfqs4zIgfGT0drJ2YMFyvz+gfRCoogc3pdXSPdYk4yMyYBdpbQkJeX98iOu8C+y8GZf1nH0qMIDdh3+S2b1i5ZQqMblPDW2lzJZ+JJPBr9xKJC+D8ajCdlHAHaa3sWODmFWd39nSM9Tq8TdpYIoh80CBlxUQXzYAcJJU888cTnb7+8eeCD8S8/Go0e/9dzD73++uuwQ0WQhp56eSy9zo5GhfB/qKxNxQIZ7BQhCWOx+RliITe5edAIO0sE0VmaJOMALZyYqR/84AcHDx785vMzYxzm6kVy2HEii268uzirHHaK/4EK4fe8JNHitSnFaK/tWeKIFQUED60mDCadWS/qG8TFdNmhI7RgAMtnxqmbjsMOEkG81v5mHlEihH8q/WSoEH6vaUCf4qISRcWwg4QqXKQQWcbRMGEwaQd0Ul4agx8NO0ioKkrI15k0sFNEEEtLzTALE8bRYq/tCagQfq++5VsZGYNx0CzqWeLkFQm7+tEdYdCYxszA68kU0auXKbQo8hbpHf2wU0QQVftpGTeNgdGr9NArDVzqfrUiIR92ihDG4PIz4jK9Hqd53AI7S0TQWZryXWy0lH4ulHkLW7iE24I2gggStc2gSKPdUh9UCM+iAKVz9qG9tucIFysLmQnopjA4dAM6sXUMF9Hu10oIicNjYxl4q/5b2EEiAukcb2KM0Wev7QmoEJ7VNdzD9hKZBfQawg05uFiRP0I2WtAwYTBo+hoKeWkYzoMdJLTJorI1nTWwU0SEsVZ1Ox8rSqHdzHxUCM+qbz0pdXKYaK/tucFFCmGvWYfmywSelyRa7X3y3Pmwg4S8osySxpEO2CkigtZYmc2O47G4sIOcCxXCs1TdtfIYek1kCkWM6DgJM8FgNXhJAnaWMNdia03zMuMkaKbMXBXnLTKwXcTIIOwg4U8zoFPSaa/tCagQnqUd7SjJQNtz+EF8njKVEdU61A47SJjTmLTiYTdHSLteppCTnyg24WC4pR52kHBHEjr3AN1WEPqgQggAADbn8BDhyC9Eh/H6AS5SSNxstPt2oOk6zhRy09Bqn7ljMZh5eJKutRJ2kDDn6jEaojBlJh3vN1AhBAAAdVd1wRiGpwthBwkHuFiZZxpGw4SBphtqVWSghRP+UZRSpDWjj26B1dZSxWVykvl0PD4aFUIAAGhoqyrCU9Be237BTEzJd+M6kxZ2kHA26raPEA6RFG0H6B/K3PkGMEo67LCDhDN1d508Rgg7xdRQIQQAAPVgkyIFbbzrN+KsYqvDNuIahR0kbOlMGtEYxc2Vwg4SJooEspYYhrutEXaQcKa2d5Vkl8FOMTVUCIGLcLd7RxTipbCDhA+eWJFP8fVWA+wgYauxrUrKEWAsNuwgYSI1EWUPSQAAIABJREFUSsBgsrtbTsMOErYI20Az11uStxh2kKmhQgh0Jm3OOBWbh+4I/YYjUohtLp0VDboExMDAgNaiL0J9GH4ljRNqehpgpwhb5paaYTbt9tqe4M9CaLPZNmzYkJiYuHHjRpvNdsnHaaLeeEIGYtD2HH7ETs0W2UltPxom9LNXX301PT29cOFKg2fw9b8cOH0a3cH4TVFWWZN7gHK7YAcJT/VtVVIe7fbanuDPWHv37s3Nze3r68vJyXn22Wcv+ThNqExaOdpr278wTCGQNVqbKUDBjhI+fv/739/xxAv9dxxO2/tFDMn4V9qtq9au02jQEUL+oUgtao3juDtRN0ZAaG0txan07cPwZyE8ePDg9u3bcRzfvn37gQMHLvk4HZAU1eyk6RrPkCYQlfIIrHsEbervH16v9+mnnwbb/g6ylMXeb2Jc8d6FP3Gsvnvfvn2wo4WJwiRJO9sz1qKGHSQMUS6HnjFanEe7vbYnsPzYVk9PT25uLgDAd/93ycfPsWXLlvMfXLx48bZt2/wY8hwdI10xHhJPkdCkz9ZutwMAvF4v7CBzRQhyxAaqurMuOpM/l3ZIkhweHuZwOP4KFqK6urqsbibIUgIABJjWTuYDAEDRmrpPfkWTH12I3G633W5nMOb6sT6Vk6AzVsltV/olFVxDQ0MkSbLZtJhO5TTWd/AYmXgGlJ/V2NhYJpN58df4sxBSFIVhmO8PBEFc8vFzbNmyBTtvJV9ubi6PF8DRO0NLQ4GTxU/NDNwlZsRXAgP6JQdJniz/U49xyHhF/tq5NEOSJJfLDYdvyNwkJCQAlx2QXsBgedimFs/VAAAwPsTn89E3h8lkEgQx9++DLLWopf1EBc4BjEv83qQ/Ho/H4/FoUgh1vbVZzJj46HgoVz+/rJzPn4UwIyOjq6tLIpH09PRkZmZe8vFzbN26dTqJ/UvT16CMyeVy6bIbusfjAQDQJ89cyOKEb5h1c/xaSJLEcTw8viFzkZ2dXSKTNJz+V1rFahvfowZLAQDg5Jvr1q9D3xwGg+HxeOb+fSjOLKmKr2ZYejg5BX4JBpHvXw1NCqHWrC/JLqTzD6o/xwg3bNiwf/9+iqL279+/adMmAMCxY8emfJw+tKOdykyarvEMddKc8k6H2UW4YQcJE6+88krUv3dV6J6IdeKE7jj4849KvMadO3fCzhU+ipILW6IwlxENE/oVSTbSda/tCf4shHv27FGpVNnZ2Vqt9tFHHwUArF69esrHacLiGHSSblHBQthBwlO0uCTLw2weNMIOEiaWLFmiV9eVZdp4PSNr1X/Y9+Olp06dio2NhZ0rfOTGZg9hHmsrOobCn9x9rc3RQJlN6/PC/Nk1Gh8ff/jw4cmPUBQ15eM0oeqqKbRTnHQR7CDhiSOU5n/i1pm0SgE6Ksg/srKyxhIY1626+fIHN6OhQb9jYJg0UaJr0wkpCu087C+G5u9imNwkXgLsIBdD0+WNwdHQViXjpIA5TzZDpoSx2AWcVG13Lewg4YOwDbTgHll2BewgYUueJjfGcTymTthBwoeqp1Yekws7xSVEdA3Q2AwKGq/xDAOKjGKdrQ12ivDRpTtBslgCvgB2kLBVlFxgTMDdRrRNgd/oxrqKs2jdLwoiuRA6vM4u74giH+21HUC5kgVOwml1RPpCN3/RtJ+S8TNgpwhnsqTCFmzc2Yrmy/gHMWxt4nhou9f2hMgthNoBnXCcispTwA4SzvA8eb6d1A6gTUf9o9FmlNPygO+wkcCNi+JGt3epYAcJEybDmXEOIzeepnttT4jcQqgyVhZS0Wiv7YBicKMkIFbTjvaG9gOvucfAJRU5aIAwsOQCWQvu9Q6aYAcJB3Vt38nwdAzQfeZR5BbCBpOmOFECO0X4kycV6Ab0sFOEg7Hmug4eVZiENogPrKJkaVtqvButJvQHjc1YkqaEneLSIrQQkhTZ7DQphWgFYcAp85e0uAcI6oJb6yHT1Gz8Lp2TEMWe096tyCUVJRcYcK8LzZeZM8rt1GOjJeJlsINcWoQWwtah9gQPJZCgXqaAS5BUJLmoNls77CAhjqK01qaiVLQiM+AkCaJuYnSkDd0RztVwm7qHh0lTpLCDXFqEFsL61qpCB5MZj6ahBxwjOk7iwdWtVbCDhDaPqbOFj8kzS2EHCX9sJluUkGskh4mRQdhZQpu65aSYlcBm0mK/04uL0EKo7qmVR+fAThEpZLFCTXcd7BShzWVoaI1hFiWF/GbQIaEoubAtU+BuQ7Od50Q9oFUmhsZPbIQWQq29szgLTUMPEnl2mX6sG3aK0GZpqbFh3tw4uk9DDw9FSYXGWBYaJpwTitK6B0pEdF9B6BOJhdA8bnF53bkSWu+GHk4KZCsHSYfdPQY7SMiiqMZ+bWGimIFF4j/Y4CtKLmwih12tqBDOnrOvtZUPlDl031PGJxL/Xal66grHKDxTDDtIpOAkpec5GZo2NEw4S56e1pY4jiIUpqGHh/ToVAIDA0M9pBN9epslfdO3yQx+LCcGdpBpicRC2NBWJeUIwuAQ6hAi5aWrUSGcLaehvjWBW5RcCDtIBClKLmjPyXC36WAHCVUNvfXK2DzYKaYrEguhdrBFmYJ2VgsqWWqR1tzc2tpKkiTsLKHHaWhoxuwyNFMmiIqSCtuSo9Ew4azpxnpofgbhZBFXCJ1eZ5d3WCYOjSHc8HDs2LEXnn5TTwyKF68VpKS+9NJLvoMqkWkhifZebRw3LoEbBztKBJElFzSxnG40TDgrxKhNj7vK8kNgKb2PPw/mDQnagUbhOBUtQneEQVJdXb1u4xbPz17bRL6Z/ZsPu0zgV8/82Ol07tq1C3a00ODuMrQKouUpaCl9UMmSCoyOfkevh/K4MTYHdpwQ064/wWCyU6NTYQeZroi7I2xorSwkoxjcKNhBIsW+ffvc6x6gSjbEueOLiBMgQwZufXPv3r0EgTZdmxaXocGYHCNLQgOEQcVn89KiU3uz0tydTbCzhJ769tMKXibsFDMQcYVQa9IWJaD5osGj0WhA4UoAwBghTsZ0AACQXWxxkCYT2t1/WlwtKgPbVZSMBgiDTZ5c2J6RjIYJZ0Ez0locUpOcI6sQkhTV6OwrEaIVhMETFRUFxgYBAAbmIoLdCwAAXjdw2fl8tHn0pVGEd7Szsc8zJE4ImQl4YUOWVGDgkWiYcKYoj1uPjc4rWAk7yAxEViFsH+6M9VApaK/tIFq3bh345jUAgIa51IITQk8HOPnmooqy+Ph42NFCgLtD35aRnJ8gZjMibjgfOnlyYbPb4m5vBCTqxp8Bs7F2mMMQJYfSeWGRVQhV7acKx5msxJAZwg0Du3btmke1g5euJlRfRDnw+don4j974uWXX4adKzS4muvb0hJRvygUuXE5VtewIynZ3WOEnSWU1BtPyFjJobULUihlnbuGrhpZdBbsFJElNjb21KlTL97yg42df+P32ktyx5u06rKyMti5QoOrpcHA9aKl9FAwMEyaJGnPzXCjYcKZUJn1iuQQOHppssgqhLrRjpIMdJBNsLHZ7Lvvvvujjz667cc7bNGEICkRdqLQQHnc7q4WvdMkR4UQEnmytCWeg+bLzABF6byWeeKlsHPMTAQVQptzaJRw5kmWwA4SueSZpc3RAG1bNU3uNt3g/7d33/Ft1ef+wL9H62h4W7K8ty3ZGl6ZzoAQZsFJKLRQCim0BUJJy6XQH5RVbhkt3DJKyi0tJaVAy7rQS7jsQIDEjrNsa1my5b1kyZLlJUvnaJzfHwpumh1b9tE553n/0VewZZ1P1BM91nc837w8Po+vkMrpzsJRlXJVB5omeswIWkCcnemRriExVcmQXttzOFQIjSNt5T4Kz4W9E7RRyhQCAd7fvpfuIMxA2A39uRkaBcNGmdikUl5um+hFElnQNUh3FmYwdXxVjCWK+AxrQcChQmjo2V8hkGN8WH1Hp8rUEvPAYbpTMEPA3mZP4MEEIY2S8aRUccpYcRFME54lw6hBm1xMd4pzxqFCaPF06hTQp4pm2vxlNmoiPDVOd5B4RxH+oKOvM+SBCUJ6VcpVXYpEmCY8S5bZEX3BSrpTnDOuFEIyTPYGvZWlMEFIM42iojtVTNiO0B0k3hE9Ziy/tGdysDyNSfux2EerUNuFJNFtpDsIAwRnvN14sKZ8Pd1BzhlXCmGHuzPXTyUVManrDyup0kqG+OSk7RDdQeIdYTcMFOQVJufhTJtuYZlKuco2M4goKjzuojtLvLO0fy5H4kScGYfxHosrhbCtp1EdkfKkzPt/iGVEfFFBUq5tsAWW4Z0eYTfYU4Qapu3HYp+SlELn7BhZrCKg19qZGAaOaGV5dKeYD64UQpPDpGHgFC4raZXarlQxOdhJd5D4FfHPhMaGO8JeWDJKOx7GK08t6c+WE90murPEO/N0X1VODd0p5oMThZBClHXWoSuAFqNxQSNXdafLAlZYO3pKRJdJVFhh9XRqYaVMHNAo1F1SDD4Rnl4kRHbwfDWq8+kOMh+cKIRDUyOiUDi7DA6diAsahboDmw7AeplTI7oMk8WlFKKUsgy6swBUKVfZSFdkZjIyM0F3lvjVbWvEMYEyhZE9LDlRCI39B1V+nkCeRXcQgBBCmbIMnkA06u6LzE7TnSVOEZ1tXWlSLez2iQ8aubrd3SEoUBM9FrqzxK+2viYNnkl3inniRCE0DRyukDLpuGTWq5Creotyic42uoPEo4hvKjzh7qAmYSt9nEgVJyeIZJ7CPNhNeBqmcXtVRiXdKeaJE4XQPNmjy66iOwX4F41c1S2XBWwwTXgShN0gKtaaPTatHD4RxguNXN2VIoZDek+JoiwRb3XZ+XTnmCf2F8JpcsYT8atKV9MdBPyLVlFh5c0ErIdhE8WJCLsBlWr6J4fK02Cdc7yolKs6qcmgaygSmKU7SzwaGbIEeagwV0N3kHlifyG0OExFs5Q4H0aZ4ogqrWTA5yRFwuBoP91Z4k7A3jaQkVyUnM+4zsUsppWrLe5OUX452WelO0t8+eSTT7Zv3/7UHx/M8Qk9bg/dceaJ/YWwrbuxkpcGvbbjiogvKkkpHCwrgU0UxwlPjUd8UzY0pYUdhPGkNK1o1OcKF6mg+/acSCRyww03XHrjz54fL3XLpaP9WHmFZt++fXTnmg/2F0Kz26ZJh4+DcUerUHfBNOEJCHsbXqq3uDugp0xc4WP8stTiHkUibKuf89Zbb732RQt66BC64PbZZH9PwW3ebz994403RiIRuqOdM5YXwjAV7gq6dSUwQRh3NHJ1B5oODnRShJ/uLHGEsBvwsqp2TycsGY03WoW6Qxggh7upIEl3lriwa9cutGEbEkrk4ZFpYahduAot/263x9/ezryTt1leCLs9PXKCSithZNcfdtMqKszjnYL8cqIL+vr/C2E3enNyeRhPKVPQnQX8G41CbRnvEmYWQHfAqMnJSZSchRCqC+1OCySGkBAhhFKyJiaY13aA5YWwrWefKijhyZLoDgKOly5JlQmlntJSaDEzJ+x1RQi/DU1oYVw0/mjlFe3uDmGJBqYJo1QqFepuRggVUIdmwiqEEPJP8UY7ysvLaU527lheCM0jxsqkArpTgJPTyiu65QmwXmZOwG4Ql1e1ezrhMN44lIQnpktSR7MzYZowatu2bZKml4TG9wncYRJcgGYn0F9/+P3vXJmRwby+gCwvhO2+IX1eLd0pwMlpFCpr0E0FiZB7hO4scYGwt+GlVUaXVatgaocOdtMoKjolYbLPhhi4HiTmysvLP/6/975remxIFBl45pf8+9U/WZ37wgsv0J1rPthcCN3+8UCYLFStoTsIODmtosI8ZhWr62B0NIqwG8LF6uHpkTLYSh+XtHJV+1QfP1URHOmhO0tcWL9+/Q9vubxImPbVn3/j7Ot8/vnnpVIp3aHmg82F0NB3QD3LE8qhy2icKkkpHJv1BMs0MDqKEAqNDSOEdVJTZWklQh5se41HGrnaPGbDS7QwOnoURZnHLCuK6teuXZuenk53mvljcyE09R2skGQhDKM7CDg5HsZTp5d1p0vJbhMVCtIdh2aE3YCXV5vGrDo4dCJeFabkTxJTswXF0H07iuxrtybwqgtX0B1kodhcCM0TXdpMLd0pwOloFer2qQFBZgHZy/UDbgi7QVyqN7utcPpS3MIQViEv70oWET1maJOLEJoy7usTR1jQ/IG1hZAIkwPh6cpSmCCMaxq52uyGaUKEKIroNgpLdFY3LBmNa1q52uob5uHSkGuI7iz0s3TuzUvMkgoldAdZKNYWQuuoJc8fSciHX67jmkau7vB0CdW1gfZDdGehU9A5gInE/Tx/uiQtGYdtr/FLq6gwu20imCZEKDjS2y4OVeWwYVk+awthW/c+NZaKCaF/f1xLEMmUMsWATBCeGg9PMrV1/cIRdgNeWmUaa4cJwjhXkV5uH+/mFVcSnD+b0G9s7FTIqjKYevTSsVhbCC3Odk1aCd0pwJnpFJVGj02squHy6Chhb8PLqsxjNiiEcU4qlOQkZg3Jk6AQ+kxNdmxGz4o9r+wshBSirIRTXwS9thlAl1FhHrPi6mUEZwshRRHdZrxUbx6zaqAQxj2toqI95EahUNjrojsLbULjzl7SI0/ISMIT6c4SA+wshIOTw5JQOLtsOd1BwJnpFZUGl0WsXhboaOFmw47gcA8vIcUrwvyhQF5SNt1xwBnoFBXmMauoWMPlD4V+wz57SR47xkURWwuhsXtfOYHzElLoDgLOTCnLEPIEo1iAn6ogB7jY1z9gbxOXVZlc7TpFBYZg22u8q8rQGF3teLGWy7sJ/cYmWyJfD4UwnpmGWisS8uhOAc6WLqPSONYuVi/j5jm9RJcRL9Ob3VYNnErPBAqpXMQXurMzOVsIw9PeoKu/3T8CnwjjmmW6X58HZxAyhm6u6SgHe61FwmSPBS+tMo1ZdaxYd8AFOkWlFU1FpscjM5N0Z6FBwLR/TFWZIJLJJWl0Z4kNFhbCGdLnpvzl5bCVnjF0ikqjq11UrAk5ByO+KbrjLCly0M5PyyBxUf/kkArWOTOETlFhcttEhZVEDxc7IvmNjR258mole/p2sbAQmgeOFM0isbKQ7iDgbBWnFI4HvJMhH16qC3S00B1nSRF2A15WbfPYS1IKRXzY9soMuoxK01g7XqLl4HqZSMBH9lnb+TNVGVAI45ixt1EjUkKvbQbhYZhGrja5uLiJgrAbxGV6g6tdnwHjooxRlFwwEZjy5xVx8LT6gOWgqFRncNtYM0GIWFkIzR47tC1mnKPThBV1AdsR7rQzpsIhos8qKtYaxyysWYDHBTwMq5CX23Ai6BqkCD/dcZaU39Q0pqqUCMRKmYLuLDHDtkIYoSL20LimDCYIGSa6cFSQnoXhkuBIL91xlgjZbxMq8yixxOru1MKSUUbRKSrMnk5RXhnR2053lqVDBUmio7UjTcSmCULEvkLY5banERF5YRXdQcC5qUgv75noD4QIcQWHNlEQnW14WXXXeG+mLCNJxIYOHdzxzTShjlOjo4GOFmFuicHbVc2iCULEvkJosH+tQkmYCKc7CDg3OF9UklJo83SK1cu403SU6DLgpXoDjIsyUEV6ec/EAFWo4tQxFAFTk0RXb3RZ4BNhXDM7zJrkIrpTgPnQZ1Qax9rxUj050MmFeRcqSJKDXXixxggrZRgo+qtbdxKfHO6mQkG64yyJSNhvOTBWXMTn8ZWyDLrTxBLbCqF1dkSbv4zuFGA+dIpKk8uKiXBRgTpgN9AdZ9GRvRZhdjHCxSZXOzta+HONPqPSNN4lVOZzpDUg0W0WpCnNAUe1Ukd3lhhjVSH0+r0+iihRr6U7CJgPfUalxW0LU2Gxuo7gwDRhwG4Ql1cNTA5JhRK5NJ3uOOCcaRUV5jErXqIluTE66jc1SfRr2pzmGnZNECKWFUJD195yUihIgvcURkoUJWQlKDvHu8UVy7jQay16GK9xDMZFmUqnqLC4bQKOHENBUX7Tfomuvs1lrmLXBCFiWyEcOFwpyaE7BZg/fYbG4LIIswqpcCg0Nkx3nEVEEf6go09UWGF0WdjUoYNTkvGkDKl8MC2B7LWy/gQxcsiOCUVOqQBDWHZCJt1xYoxVhdAy0aPL0tOdAsxfdYbW4DQjhFjfgJvoMYvyyzGhCHrKMFqVUmuc7uWnKoIjPXRnWVx+Y5NEv6bVaaph3cdBxKZCGAwH+yPTGtV6uoOA+atSakxj1ggVYf0miui4qNM3FowEcxPhMF6miv7qxoWzCf3GRol+jcFlZuUABnsKoW3EkE1gidlldAcB85eCJ6dLUrsn+nBVDdljZvGqdMJuEJdXGccssF6U0aoyNMaxdtZPE4bGhqnArCivrGXUWJvJwlE39hRCg32fWpAGvbaZripDa3BaeJIEQVYhydI3l4h/JugaEuarYAch06VJUlPwJIcyleg2sbhHrt+wT6JfMzA9wucJ2DdBiNhUCC1jVk16Od0pwEJVKbVtrm+mCVk6Okp0GfFiDcYXGF3QU4bx9Bkao2+Ih0tYvLzLb2yU6OtbncZa1u0gjIpZIfR6vQ0NDWlpaZs2bfJ6vcd+q76+HvvGtm3bYnXF41hJl76kfpGeHCyZqgyN0WWhEMXiTRTRCUJvYNLj95akFNIdByxItVJrdFlExVq29loLT3pCnlG8RMfWcVEUw0L4xBNPFBQUOByO/Pz8J598cu7rFEV1dHQ4HI7p6enp6elnn302Vlc81sjEEAqH80pXLMaTg6Ukl6QlihL6JgdFeeXhaW94wk13otgj7Aa8rKrNadJnVPIw9ozKcFN1hrbNaRax95Bev6lJrFlJ8XgGl6UGPhGe3j//+c/t27fjOL59+/Z333137usOhyMUCl1++eVZWVnXX3/91NRUrK54rLbOL9VhKSYSL8aTgyVWlaExOM0Iw8TlNQTrDqyPzEyEJ9yi3NJWp4l9rao4SCGVS4Ril1LO1mMo/IZGiX5Nz0R/gkimkMrpjrMoBLF6ouHh4YKCAoRQ9HPh3NcdDkddXd3TTz+dn59/55133nHHHa+//vpJn2HdunUnfnHDhg133nnnGa/e1n+kTJIzPj4+3/jxYmZmBiEUDLJ2teTZKEsoOjTUul6+KpynnjI2TacVCYVCukPFTNjSjOWVj09MHHEYzstYNY+bdnJyEsdxv5/9fcnPEkmSMzMzGH0L5SqSyw5M9q8nCU9vJ5YcF6XC6/WGw+GF/8OhZqfJgQ6/onBfz5eVKeVMfI9NTk7m8/mnf8yCCqFare7o6EAIURRFUVT0RqQoKhwOzz2mrq7uiy++iP75t7/9rUZzyqUBDz744Im3slKplMlkZ0zS6R+5qHDT2TwyzkUiEYQQC/4iC7E8t/ZvHW/JZLKIbvX4p69JxTibXpDp4U6JqoYQhKfI6crMCt65v30Hg0EcxyUSyWLEYyKhUEhRFI03SW121RGn4cIijcDZK84uoCvGsQKBgEwmW3ghJGwHImXVspRUa3vnhry1TPyXyOOdeeBzQYXQZrPN/Tk7O3twcLCsrGx4eDgn5199zo4cOUIQRH19PUJIJBLh+ClPCrz44ovn9zudP+h3YH6dZuNpnpwpSJJECLHgL7IQuXg2zhc5CXe+PEeQrqTcQ3hBId2hYsbbbUpet6nFa9NnVErE8xnMx78R82wMhWEYSZI0viDLsqtfNr8uKbsyONiJr7qErhjHir7ZLrwQzrQfkNWsF4qEZrftntV3sPWui9kcYUNDw86dOymK2rlz5+bNmxFCX375JULI5/NdeeWVVquVJMlHHnlky5YtsbriHEt3UwEpEKezcHcLZ9Uoda1OI0IIV9VRveyZeglPeiKz08LsojanuUbJzgV4HJSVoOTz+O4sJcuOoaBIgugyiStW2Md75NL0VHEy3YkWS8wK4UMPPWQ0GvPy8iwWywMPPIAQ2rBhA0Jo3bp1Dz74YENDQ05OjtfrfeKJJ2J1xTltvc2VOFRBVqlR6lpGo4WwNtJtpDtOzEQ3TiAMa3UaazJhpQx7VGVoLNhUeGo8MjNJd5aYCVgPiQrVPGlCi9NYy+rf22K2WCYlJeWDDz449isURSGEMAzbvn379u3bY3WhE1m9XZcqaxfv+cHSq8user5lZ4SiREWVyOuMzEzwElLoDhUD0Y0THr93gpgqSo6LySQQE1UZ2jaXZWVBBdFrkehYsqHZb2yK/l1anMZNpZfSHWcRMX4PE4UoW8irL4Ve26wil6Yn4Yk9E30YX4DlqwOdbXQnig3C3oaX6VudxuoM7TyWyYC4VZupbxk14CVa1myioMKhgO2wWLc6TIUtY7ZqNh46MYfxhbDf2SELU8oCGGVim1qlvmXUgBDCirTs6LUWHndRoaAwI6/NaYYdhCyTnZDJ5wlc2UrWbKsn7AZBRi4/Ka3d3ZGbmJ0oSqA70SJifCFs7fhKjSWjs1ggC5ilNlPf6jQhhLASPWE7zIKOxgF7K14WnSA0sbVnI5fVZeqNvJmgc5Ai2LDF0286Oi56yNFWl1lFd5zFxfj6YRk1a1JL6E4BYq9WqTe4LKFIGEuWY2IZCw4+jU4Quv3jM6SvMCWf7jggxmoz9S0ukyivlOiz0p1lwSgqYNov0a1GCB12tC3LqqY70OJifCFsDzh0BdBilIWS8MSshMyOcTtCiB0NuKNLRltGjdVKLYZggpBtlmVWtzpNwmINC44PI/usPFmSQJHjC872TvZrFRV0J1pczC6EU/4JLyJValgpw05134yOitV1ARuzC2FobBjx+AJ5VqvTyO51B5yVKk6RS9IGs+QsOIbCb2qS6OsRQq1Ok0auxvkiuhMtLmYXQqNtT0lIzBczr+sPOBs1Sn1LdFt9aVVwqDsS8NGdaP6IzjZxWRVCqGXUWJfJ8oEmzqrJ1Jv4PnKwiwoxu12w37Rfol+DuDEuipheCE0DLZWynDM/DjBTVYbG5rGTkSAmFIkKKwi7ge5E8xfoMuBlVQOGnOFdAAAfM0lEQVRTw2Eqkp8ENy071Sr1rW6rUJkXHOykO8v8BR19KBwS5pQghA6Pti3jwO9tzC6ElqkebTbLlzNxmVQoKU4p7JzoRgjh6jrCxtgjmSiK7DLhJfrDjrblHPj9mrNqlDqz28orqiSYvJvQb2yMfhx0+samiZliDpwdzeBCGKEi3WhGpz6P7iBgEdUq9WavDSEkVi8LWA/RHWeegqP9GC7hp2UcGW1bllVDdxywWBJEsvyk3J6sVEbvJvQbm8TRcdHRtmVZ1Vzo/MDgQtjZezA9xEtVFNEdBCyi2ky9edyGEBJm5lMUFXIN0Z1oPgi7AS+vDlPhNpeZ3T0bQV1mlVnoJ3vaUSRCd5b5CI07I1PjeGEFik4Qsn0HYRSDC6Gxq0klTKc7BVhcmnTVkG9khvSh6NpRZm6iIOwGvFRvddszZUoWt/AHCKG6zKoWj42fkh4c6aU7y3wEjI1izUrE40UoqmXUWMeNkXwGF0Kz26ZLV9GdAiwuIV+oSi6Nrh0VVyxjZK81iiK6TXhp1eHRVpggZD2doqJ7ojdcXMnQ0dG59aJd3p4kPDFDKqc70VJgcCG0Bcd0xWvoTgEWXbVcc2DkCEIIL68hey1UkKQ70bkhh7r4SWn8pNTDDgPrW1UBEV+kSiu1ZyYycb1MZGYiONyDl1WhbyYI6U60RJhaCN3jgz4UKS6DnjLsV5WuPeRoRQjxxFJhdjHjditHx0Vng/7uiV59hobuOGDRLc+qaeP7yB4z4xrk+s3N4oplmFCEEOLUCmemFkKjbY8KyTBezM5TBHErW6rEENY/OYgQwivqGDc6SnQZ8PLqNpepIr2c9R06AEJoZXbtoXErJhSF3CN0Zzk3fmOjWF+PEAqEAjaPvTqDK63hGVsIhw2aRFgvyhUrsmsPOloRQmJ1HcGsXmuRMNnTjpfoDjnalsPGCW4oSS3yB/3jJaXMGh2lCD/Z0y6uWI4QanEa1ellUqGE7lBLhGGFMBKJvPPOO/fdd1/bRHcixYlZXIAQWpFVc9DRghAS5ZZFfNNhr4vuRGeLHOjkyzN5siSOtKoCCCEMYcuzasxyMbPWy/gtB0QlWp5YihBqHj6yMruO7kRLh0mF0OFwLF++/Op7n36qUziCh2+9+w8//vGPw+Ew3bnAoqvNrLKM2cgwiTAMV9UwaHSUsBvEpVVjs+5JYqokBcYwuGJFdm0LmiAZNZ8912gbIXRgBAphvLrttttaUuvRPV/qN65KJfGJX7S89OnhF198ke5cYNHJhNKS1CKDy4KiLWaYUwgDdgNeVnVgpIUjHTpA1LKsauNkN0kGwhNjdGc5K1QoSHS0SDSrEELR+fjC5Dy6Qy0dxhRCn8/3wYcfoS3/iRAqoA4KSSUSSdEV97355pt0RwNLYXlW9aG5acLONiocojvRmVHhENlvExVrmkcOr85ZRnccsHSSRIkFyXk9xXlEj4XuLGeF6GgRZhfxEpIRQvtHDnPq4yBiUCH0eDwhSSoSJyKEBlFdO7YZIYTkRS4XY6aLwEKszKo7ONKCEOIlJAvkWWS/je5EZ0b2WYXK/LAIbxk1rsiqpTsOWFKrsuuM6SKSIetl/Kb9Et3RcdGDIy2rcqAQxiWlUikNTaPJUYTQQdElrcINCCE0ZCosLKQ3GFgaZWklE8Ska9aNEBKrlxFWBoyORluMGlyWopT8ZDyJ7jhgSS3Pqm2JuJmxXiYS8ZubxdrVCCF/KGDz2Gs41hGXMYUQx/GtW7ei125H5OzRL3n60Xu/uvnmm2nNBZYID8NqM6sOO9oQc3YTRrfSN48cXpUN46Kco04v84ZmXT53xDdFd5YzIHos/BS5ID0TIXTY0VYpV0kEYrpDLSnGFEKE0FNPPXWDLgV7UIv+shU9f1XCE/XPPnT3li1b6M4Flsiq7LrmkcMIIbywIuQZicxM0J3odCiSIIe68GJN88iRVTBByD08DFuWVd1emBH/04R+U9PcuCjX1otGMakQSqXSV155xdq0+x93bNr165t7Oix33HEH3aHA0lmVs+zIqIEMk4jHx0urAvF9Ti/Z1y7KKXaQE7PB2dJU2DjBRSuyaoxJPDLuR0cDx2ycOOhoWQWFMP6pVKrvfe97DQ0NCoWC7ixgSSWJEktSCludJhQ9iSK+j2QKdLbhZVWNQwdX5yzHEGyc4KIV2bXGsNsX3+tlgkPdiC8QZhUihHon+vkYPy8ph+5QS415hRBw2ZrclY1DBxFC4orlgY4j8dzUmOgy4GXVMEHIZanilNzkXJOvnyL8dGc5Jb9x39zHweaRI1xbLxoFhRAwydrclY1DByhE8VPkPFkyOdRFd6KTowh/0NEfySm0ujvh6CUuW5u3qjU7KZ53+/iNTRLd0fPs9g01r8lZSW8eWkAhBEySk5glFUo7x7vR0XN643R0lOg2ifJVh93tFfJy7nQuBidam7vyoJQIxGuvtZB7JOKfERWoEELewMTA5HCVUkt3KBpAIQQMsyZ3RdPQIRRtMROvuwkJuwEvq2oeObw6ezndWQCdCpLzpCKprT9Ob1S/YZ9EtxphGEJo72Dzyuw6ISfPtoNCCBhmTe6KxqEDCCG8RBcc6YkEfHQnOomA3SAq1TcPH+bmjAs41tqCNc3+vvhsCug3Nkn0R8dFG4cOrsnl6FHnUAgBw2jkFW7/uNPnwoQiUVEl0dlGd6LjRfwzYbfDLkOJeGJuYjbdcQDN1hasPZIqCA7a6Q5yvPCkJzQ2JCrRIYT8oYBprJ2DOwijoBAChuFh2Krsuv3DhxFCYnU8tpgh7EZRUcXXwwfOy6unOwugX6VcNSXi9Xc20R3keH7TfrFmJcYXIIT2Dx/SZ1Rydj4bCiFgnvrcFY3DB1H0SKb2g3THOR7RZcBLq/YNNp+XD4UQIB6GlUUU77Z+cOedd7766qvBYJDuREf5jY3HrBc9sDaXi+tFo6AQAuZZkVVjGbPNBv0CZR7GFwadA3Qn+jeE3dCfnc7DeMUpBXRnATQjSfKqq656+7XDXUmhHY6crU+8WlNTMzQ0RHcuFPHPBAc6xepahFAoEj400lqfw9EJQgSFEDCRWCDWZ1Tuj/YdVdUQ8dRiJjIzEZ5w7yeG4OMgQAjt2LHjXaNz8MpPB2RYxbqN6D8+sORecvvtt9OdCwXMB/DyKkwkRgi1Oo35yTlpklS6Q9EGCiFgpPPz1+7p34eO7iaMo2lCwm7ES7RfD+6HQggQQu+88w667BcRvizVr1BRnyGE0OX3ffjRx34/zb1m/KYmiXZ19M97B5vX5q6iNw+9oBACRlqXt6pl1OALzuJl1WSflQqSdCc6KtBlGC3MJ8LB8rQSurMA+nk8HpRegBAai6xI5x3MI51IkhQSJUxOTtKYiiIJorNNrF2FEKIQ1Th8kMsThAgKIWAomVBardQ2Dh3kiaXC3BKiy0h3oqOIzrZmaeC8vNXQaBsghIqLi1F/C0LoK+G1HYnU3wd+oR7dn8QPyeVyGlMFbEeE+eU8aSJCyOruTBTKONho+1hQCAFTXVCwbk//XhRdOxofo6PhSU9kdnqv13pe/hq6s4C4cNttt6F/Poic9jCSuKgVj2fXvD7y+EM3Xy8Q0NnAxW9qnNtH/1nfVxsK1tEYJh5AIQRMtSZ3pcFlmSKnxeplAeshuuMghBBhN4yXqSaJqUq5iu4sIC5s2rTphd8+lPLsBvT0pa5de0fHP/xKkH0t7ibo6z5KhUOB9kMS3WqEUISivh5ouqBgLV1h4gQUQsBUEoF4WVZ14+ABYU4xFfCFxp10J0KE3XBIga/PW83DYFwUHHXrrbf2dXXu3nHfzptuza7Muu63T8pv/KXnr48F2un57Y3sNgnk2fzkdIRQq9OYLk3j+LgogkIIGO2CgnVfDOxDGIar6og4GB0l7G1fhoYuLFxPdxAQX5KTkzdu3Ljp8obzC9Z83r8XL6uW3/pr7xvPzLZ+tfRh/MZ/jYt+3vf1hQVwu0IhBEy2KntZu7tjgpgUq+toP7A+5BntFgQIFNYqKuhNAuLWhYXnfd73NUJIlFcuv+3xyf/9s2//R0uagKL8pv3Rk3iDkdC+oQMXcH6CEEEhBIwmFuArs2r3DTaL1XVEl5HeBv9El2FffsplxRthvSg4laoMrTcw2Tc5iBASZhUqfvq76d1vTn/xP0sWgBzo4ElkAkUOQujAyJGilAK5NH3Jrh63oBACZttQsPaL/n08WZJAkU32WmlM4utsaxJMXFy0gcYMIM7xMOyCgrVf9O+N/qdAnqW44+nZQ7sn39+5NAGOHRf9ou/rjfBxECEEhRAw3crsui5vr9M3RvuB9QdGWwqS8rISlDRmAPFvY+H63X3/mhrkJ6UpfvpfRJdx4u0/IIpa7Kv7TfujhTAQChxwtED/oygohIDZRHzRxsL1H/XsFqvraCmEFEX19fW1ffHxl0mhS1WXLn0AwCyqtFI+xmt3d859hSdNlN/2eHC0f/wfv0OR8OJdOjjaT4VIYU4JQqhx6KBWXpGMJy3e5RgECiFgvCtKL/6w+3NBvio87gpPeZfy0s3NzTU1NUV163/35OMGnPjov3fNzs4uZQDARJeVXPh+18fHfoUnlsq3PRbxTXn++hgVWqxzmvzGJomuHmEYQmh331cbYXnzN6AQAsYrSSlMFScfdhnxsmqiY+k2UfT29l50yWWGFXejJ7ozLqniE2V/2N1z6623LlkAwFCXFV+4d7B5hvQd+0VMKEr/0a+QQOB+4X6KWJSW3H7jvui4qNs/bh6zrcvjdKPtY0EhBGzwrZILP+zevcQnUfzpT3+aWfZ9tOIajKJcCYPdWAP68cv/eOe90dHRJcsAmChVnLwiq+aT3j3HfR3jC9JvuFcgzx77432R2ZnYXjQ87gpPuPEiDULow+7PNhSslQjEsb0Ec0EhBGxwYeF5hx1tgeLyQEfLEqw4iLJarVhZfY2/467xZyeFWA9vDZIkRzIrOjo6liYAYK5N5Zftsp9sByGPl3rNHXhR5dgffhGejuU4v9/UJNGsQjxehKI+6PrsitKLY/jkTAeFELCBTChdm7vyc4+Bn5hKDtoX+3JUOBToaNmaxTskeuu/hnf0JQ+NhC6mov+appxJSbAAAZxBdYYWQ5jBZTnJ9zAsefPNkur1Y8/dHR53xeqKflNjdB/9IUdLijgZjgk7FhRCwBKXl170ftenuLp28VrMUGTAb9g7/uqTjge/N/XRq8pyzXc+7b+08OGRBI8ZuxohhAwf5Ekjer1+kQIANrmi9JJd9o9P9d2ki7+XsG6Ta8fdIdfQwq8VmZkMDvfi5TUIofe7PmmAj4P/js6jQACIIZ2iEsOw3jxlzt6vky65LobPHPFN+c37/cb9ZLdJVKiW6Nckb/4xPylNHomsPNQdbvyuW6snhxqRfV/ioVf++tbrfD4/hlcHbHVp8QUvm173BiZTxcknfUDC+s08sXTs+Xvkt/w6uudh3vyWZlxdiwlF435vm9N83+o7F/Js7AOFELDHFSUXfThuv8nRF5mdjh46uhBhrytgPey3HCC7TaISnbR6Xdr1d/MkCXMP4PF4r/79tc1vXJ9+aFbVv1NXp7t9Z1t2dvYCrws4IkEkW5e36qOe3ddVXnWqx0hXXMSTJbr/9ED6TQ+IijTzvpbf0ChdvhEh9EH3ZxsK1kqFknk/FStBIQTscXnpxde+9/a3S0pT7W2Sqnn2jgp5HAHzgdm2vSHXoLhyuWz5hek/+CUmOvn6uq8GGsuVJc888+gCUgPu2lx22cP7nry24koedspZKrFmVdoNEs/OR1Ovu0tcsWweV6EIP9ljSdt6b4Si/q/r00fW/3IBkdkJCiFgD5lQ2lB6yYf91h9aD59bIaQosq/db9rvNzYiihLr6pMbbsILKxHvDJPob9t2bdV9d0GhAYep08sypPIv+vdeWHjeaR6Gl1Wl//BBz85HUq6+XVJ1zofoBtoPiYoqeWLpwZGWZHESLJM5ERRCwCpXqRpu7PyowT6QehYPpkJBwt7mN+0PmPfzElMluvr0mx44+8kY85htmpxZlb18IYEBx23VXrPjyIsXFKw//WHOoqJK+bZH3X9+KEL4ZSsuOqdLzK0Xfcv23pXlly8oLktBIQSski5JXa6o/d/+fX3XbAmnKDdt2nT55cf/y48EZgPWQwFjU8B2RJhVKNbXK+54WpCeda7Xetv23lWqBjiMHizEsqxqqVCyd3D/GftfC3NKFNufHPvjLym/L+G8LWf5/FQoGLAdSb5ym81jH5gavqjw/IUmZiPYPgFY5dNPP33utud2y1EoU/fnWf0VN999ww03UBSFEApPe31NH7r/9ODowzfMHv4CV9Vm3v+S4me/Szz/2/Oogt0Tfcax9stKNi7CXwJwyw3aa14xv0WhMzeCEChyMn721Ezj/019/NpZPjlhbxNmFvITU1+zvP29yisFPFjSfBLwiRCwRygU+uGPfuS5ckcw0ijPdyHdQ2j19XufXbf/uYdL0UzI2S+uWC5beXH6jfdh+EJXzf2p9W/Xa74DTarAwtXnLH/Z+Pr+4cP1OWceZuenKDJ+9ruxF+6P+H0pW25BZxqQ8BubxPr6/snBdnfHg2vujlFktoFPhIA9Wltbh0kx0l/eTm1tT3Pc5Xrls/7/9+4lRc4OU9Il12U98kbaDfdIqtctvAoaXJaBqaFNpZfEJDYA39dc/YrpzbN8MC8hRXH7E+RAp/eNZ1AkcrqHRiIBc7NEV/+q+e3vqDfjfFEMsrIRFELAHlNTUygpAyHkQmphKNWT0HNvzvblgatfc2JidR3Gj834B4WoP7b89eaqG4R8YUyeEIB1easDocAhR+tZPp4nSVDc9lh40uP522+ocOhUDyP62nlJaS4cHXK0bi67LEZhWQgKIWAPtVqNjVhQYBoh9CHv4Y60AaMsK9K9v6KiIoZX+WqgKUyFLyic5z5FAE7Ew7Af6K79c9srEeq0n/COgYnE6Tf/J0LI/cIDpzq2yW9skujX/MPyzpbyb8Em+tOAQgjYIycn59rN30I7f4hmvbNUrjN8XtHEo0mtb95yyy2xukQoEv6L4bVtNTdiCBaLgljaULA2QST7Z+eHZ/8jGF+Q/oNfCtKVY8/fG5mdPvEBAdP+qXL1V4NNV6mviF1SFoJCCFjlxRdf3L4mT3C/Gj2+tv/+x5UJ9pd3/TU/Pz9Wz7/L/rFSpqjLrIrVEwIw587l214xvTnuP5fTl3i81Gv+Ay/WjP33vZGZiWO/ExruQRj2fP+H16i3JIkW2nGQ3aAQAlaRyWQ7duzwDPcdeGNH174v7rvork9n957NwvSzMTQ98jfTGz+tuzkmzwbAcfKTcq8ovej5lp3n9mMYlrzlFmnNea4dvwhPjCGERkZG9u7dO7jnvbaKouEZxzWVVy5KXBaBQghYKCkpacWKFYWFhZcUbxDwBB/3fLHw5wxT4ceanvmB7trC5LyFPxsAJ7VVd63FbWsZNZ7rDyZu/G7C2gbnsz+/88brctVV3/7l73v2ffKka/+F4jVCHmyTOwMohIDNMIRtr/vRn9tecfoWesDpK6a3EkSyK1XfikkwAE4K54t+WnfzM4deCEZOuRb0VBLWbdrZMf59qbfi4d2FP3vlgDbFQah+8u2f9PX1LUJSVoFCCFhOlVb6/cqr7//q8UCImPeTtLs7dnV9fO+qO2CNDFhsa3JXFKXk7zj84rn+YE9Pz4PvN/9n4Z1/H/ntjyb+si8d606837fsupdeemkxcrIJFELAflerG0pSC59sfm5+k4W+4OyjTU//fPm2dMnZtPIGYKHuXXWH0WV5z/7ROf1UZ2cnytX9X9qGn+f9xC5vC85eGKSSUdGKzs7ORcrJGlAIASfcteInwzOjb7T/81x/0BecvfuLX63JWbEub/ViBAPgRFKh5PHzH3jZ9Ear03T2P5WamoomRjAUcaZ+3S2oauJvRQihieGUlJTFCsoWUAgBJ4j4osfW3/c/He/vHWw++5+aIX0///zBSnn5T+p+uHjZADhRdkLmr9bc/evG343MjJ7lj9TV1eXhRMnkQ3zkb4rcF8SEyD+J9r60ZcvZHlXBWVAIAVfIpemPr7//94f//DfTm2czRjpFTv/88wf1GZqf1t0MU4Ng6VUrdTfqrr1nz6+Hpx1n83iBQHDzH25Pixxp/7ubOvA2+vQZ7JEVP732issug+ZqZxDjQhgOh9Vq9bFf8Xq9DQ0NaWlpmzZt8nrPZaMoALGmSi/906VPHXK03vflozOk7zSPbBk1/vTTe2sz9bfXwmdBQJvNZZd9t2Lz7Z/+v72D+0//yAgVeb39XXO489XrXnjk/IKrJ96/M3f0q3f+9txzzy1NVEbDoke1xcTvf//7f/zjHwcPHjz2Oe+9996ZmZmnnnrqrrvuSkxM/M1vfnPyHBgWiUQwzp9xOj09jRBKTIQ2EEdFIpGxsTGlUhnD5wxFwv/dsrN55PB1ld9en1efhP/bq907OfBCy8uD08O3VG89P39NDK8bKxMTEziOSyTQOvIokiSnp6fT09PpDrJYOjxdv9r3xHl59bfUbOVjJzlQsMvb8+SBPyQIZb9YuT0rQYkQcrvdycnJQiH0hT8rsSyEe/bs8fl8DQ0Nxz6nSqV677331Gq1zWbbvHlzR0fHyXNAIUQIQSE8wWIUwqiDIy0f9Xx+YOSIRq6uydRNBqbc/nGnzzU07bhe853N5ZfF7TZkKITHYX0hRAhNkdOPNT7TPzW4Pm/1urzVGrmah2Eev3dwaqhx+NDuvq9urf7BJcUb5sbwoRCek1gWwqPPiP3bcyYkJIyNjUkkEr/fr1Qqp6amTvVTKpXqxEJ4ySWXPPDAA7FNGM9mZmYQQgkJCXQHiReRSMTj8SgUikV6fiJMtriNHd6uFDw5FU/JkKQXJOZJBXFdYyYnJ3EcF4vhTOCjSJKcmZlJS0ujO8ii65sePOBsOehq8RDecCQswATZssyS5IKrSxqShP/22/P4+HhiYiIUQoRQSkqKQHCGX2oXWgjVanX0Q97c8xxXCGUymcfjEYvFs7OzCoXC5zv5xAyGYS0tLScWwpSUlJycnIUkZBb4RHicxftEyFzwifA4XPhEeBy334Pz8UTRKX9jhk+EcwQCwRnHGhc6+GOz2U7/gOzs7MHBwbKysuHh4dOXtOrqahgajd64cPvOiUQiQqEQXpBjCb9Bd5B4QVEU116QLGHm6R8AN8k5WcTtE19++SVCqKGhYefOnRRF7dy5c/PmzYt3OQAAAGAeFrEQbtiwASH00EMPGY3GvLw8i8XCqdk+AAAAjBD7Qjg3QRj9Q0pKygcffDA0NLRr167k5OSYX45l9uzZs2fPHrpTxJHZ2dmnn36a7hTx5Z133jGbzXSniCMDAwMvv/wy3Sniy1/+8heH46y24QMEnWXiTXNzc3PzOfQAY73Z2dkXXzznNvzs9uGHH1qtVrpTxJGhoaG33nqL7hTx5e9//7vT6aQ7BWNAIQQAAMBpUAgBAABwGhRCAAAAnBb7zjLzk5mZefHFF9Odgn4WiwUhpNFo6A4SLwiC+Oyzz6644gq6g8SRAwcOZGVl5efn0x0kXrjdbrPZfP7559MdJI58/vnndXV1cBIhQujRRx894z+WeCmERqOxra2N7hQAAABYpaGhITU19fSPiZdCCAAAANAC5ggBAABwGhRCAAAAnAaFEAAAAKdBIaRZfX099o1t27bNfd3r9TY0NKSlpW3atMnr9dKYcOm9++67Wq02JSVl/fr1nZ2dx37rVC8XW53qNoDbA26POfAesnBQCOlEUVRHR4fD4Zienp6enn722WfnvvXEE08UFBQ4HI78/Pwnn3ySxpBLrLe398Ybb9y5c6fD4WhoaLjpppvmvnWal4utTnUbwO0Bt0cUvIfEBgXoMzw8nJSUVFtbm5CQsHnzZqfTOfet8vJyq9VKUZTVai0vL6cv41LbvXv3rbfeGv2zy+VKT0+f+9ZpXi62OtVtALcHBbcHRVHwHhIj8ImQTg6Ho66u7qWXXurv709OTr7jjjvmvjU8PFxQUIAQiv5OR1/GpbZx48YXXngBIRQOhx966KFrrrlm7lunebnY6lS3AdwecHtEwXtIbNBdiTlHpVKd9JUfGRlJTU2d+0+pVOr3+ymK8vl8Uql0SSMuuRNfk08++aS6uvqee+4JBoMn/ZHjXi62OtVtwKnb40Rwe5wUl99DFkhAXwnmKJvNNvfnI0eOEARRX1+PEBKJRDiOz30rOzt7cHCwrKxseHg4JyeHhqBL6NjXhKKoe+65Z//+/W+++WZ5efmxDzvNy8VWp7oNOHV7HAtuj+PAe0hMwNAonXw+35VXXmm1WkmSfOSRR7Zs2YIQ+vLLLxFCDQ0NO3fupChq586dmzdvpjnoEvr666937dr1/vvvZ2dnz8zMzMzMoG9ek5O+XOx24m0AtwfcHseC95DYoPPjKOdFIpEdO3aUlJTI5fKtW7dOTk5SFBX9P8Xr9X7rW9/KyclpaGiYmJigO+nSefjhh0+8RaP/e9KXi91OvA3g9oDb41jwHhIT0GsUAAAAp8HQKAAAAE6DQggAAIDToBACAADgNCiEAAAAOA0KIQAAAE6DQggAAIDToBACAADgNCiEAAAAOA0KIQAAAE6DQggAAIDToBACAADgNCiEAAAAOA0KIQAAAE6DQggAAIDT/j+TveA4VW5DFAAAAABJRU5ErkJggg==" }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "li = LinearInterpolation(x, y)\n", "li_spline = CubicSplineInterpolation(x, y)\n", "\n", "@show li(0.3) # evaluate at a single point\n", "\n", "scatter(x, y, label = \"sampled data\", markersize = 4)\n", "plot!(xf, li.(xf), label = \"linear\")\n", "plot!(xf, li_spline.(xf), label = \"spline\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Univariate with Irregular Grid\n", "\n", "In the above, the `LinearInterpolation` function uses a specialized function\n", "for regular grids since `x` is a `Range` type\n", "\n", "For an arbitrary, irregular grid" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "hide-output": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAyAAAAGQCAIAAADZR5NjAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3deVxWZf7/8euc+wZBUVDE2ETct9zLjNRQK01FQatvlpm5TJtT03cqZ2pss1KbcnKs+fV1JqdBzRwTXFvMUsGlNDVRU0LFDW8kBJTl5oZzzvX74zYid+C+ubfX848edO7DOR/rDt69z3Wfo0gpBQAAABxHdfUAAAAA3oaABQAA4GAELAAAAAcjYAEAADgYAQsAAMDBCFgAAAAORsACAABwMAIWAACAgxGwAAAAHIyABQAA4GD1HbCklBkZGXU5gqZpPN4H1VVWVrp6BLgRKaWmaa6eAm7EMAxd1109BdyIruuGYTj7LPUdsCoqKm6++ea6HKGoqIifnqguPz/f1SPAjei6XlhY6Oop4EZsNltxcbGrp4AbKSsrKysrc/ZZuEQIAADgYAQsAAAAByNgAQAAOBgBCwAAwMEIWAAAAA5GwAIAAD4hIyMjISEhNja2TZs2Y8eOzczMdN65CFgAAMD7bdmyJT6uX1zZ0Y339FifeOON+T/G9b35hx9+cNLpzE46LgAAgPt44YUXXri1zf1do+1/O7VXKynlSy+9tHr1amecjgYLAAB4OSnlju++G9b2huob7253w/bt2510RhosAADgYWw2W4MGDa5nT0OKjAK52WJolzwdx6kP3qPBAgAAnsFms73xxhstW7Zs0qhhRETE9OnTS0pKLt2t0hDbzsg5e42RX2qhiyof3KgfOifa9+r3xZEz1Xf74siZuLg4J41KgwUAADzDAw88cHZX2uLBndo27Xqq2PpmavLd27Zt3rxZVdUyTXz3s0yzyLRcY+fPskOwMiBcmdxR/fftaliAEEJse3dWwtA7Kw05ot0NhhSpmZb3957emP5fJ41KwAIAAB7g+++///arzzY+1L+hn0kIEd048O/Dug9fuv1/5q093W54RoHs0UwZGK482810W7jSxO/ib4+Li9u07dsZM2a8s2KLyWQaOHDg9p2rOnTo4KRpCVgAAMAD7Ny5My66mT1d2ZlV5fZWzQ8d+v7N+xP6himB1wo13bp1W7lyZXFxsaIoQUFBTp2WgAUAANxanlWk5RrLj6tB+sUr1St0Y3C03+0RiksGuwoCFgAAqA+GYWRkZJw4caJNmzZdu3ZVlKulolyrSLMYm3PlptPydJkcEK72jRvw4fzC/LKK5g397fuUVuobsvMW3X57vYxfMwQsAADgdHv37p00adL5Y5mtQxr9dLYkpnufhQsXtmvXrvo+p8vkZovcbJGbc2WeVfYPVwdFKJM7qD1CFZMihOhSPvWxB5cs/NNtHTo3DzpaVPb29qxbh40aOHCgi/5MV6NIp94F4hI2m61JkyY2m63WR8jPzw8ODvbzu2T1GnyVxWKJiIhw9RRwF5qmFRYWhoWFuXoQuAur1Wqz2UJCQlw9iE87f/58x44dn+7SbHy3lkIIQ8r3d2an5ov9+/fnVfptssi0XLnJIgttckC4Gh+h3B6hdGuqqJdruP773/9+8MEHWVlZsbGxDz/88KRJk1S1ZvecYg0WAADwBikpKZ0DNHu6EkKoivL7vm0+X7o99uUv9RuHDQxXB4YrT3VVuza96lVDIYQQ991333333efsgeuOgAUAAJzr8OHD3VsEX7Sxxw3BLUOOvvGgn9stUHcEAhYAAHCKM1ax8bSx0SJTTzWNLym/6NXTxeXDopt7ZboSBCwAAOBAZ21is8XYeFputMjTZXJguDooQhn1+9HjB/7lSGFp26aN7Lv9kHtud0H5J3fd5dppnYeABQAALkNKuXbt2l27dgUGBsbHx99yyy1X2vNchUjLNb45LTdZ5NHzsn+4MihS/U8HteeFT/8JIdrNnvvu2Of+eF+HFm2aNjyYX7Iqu3DBvz704s+j8ClCeDw+RYjq+BQhLsKnCGvn559/Hj16dNnRHwfFNrdpxpqs3DvH3v/Pf/7TZLpwI/WSSrHljLRfATxUJPu1UAZFqoMilJuaK+YrfKrv8OHDn3zyybFjx9q1a/fQQw9FRUXV35+nmvr5FCEBCx6PgIXqCFi4CAGrdsaNG6fs3zJnSFd7A2XV9LHLd0x84fUe907beNrYZJF7C2Sf5sqgCHVwpHJLC8W/ZrdKcCVu0wAAAFygoqJiZWrqjom3Va1ADzSbnrml3dT3lvZr8/igSOXVPqa4Ftd+9p8v458NAAD4jbMFhf5Cbxrwm4tFsSGBMUbulgSSw3XhHxMAABBCiKPF8usc+fVp+c3JJiVq4Oni8sjGAVWvHsovaR3byoXjeRbPuWQKAAAc7edy8d+jxu+26G2Xaf3XaGm5cli0svuegCcmP/zixh+tmm7f7Uypbc62rClTprh2Wg9S4wZL1/WuXbseOnSo+sa4uLjt27fbv3700Uc/+OADx0wHAAAcrVQT6bny6xxjw2mZXSwHhqt3RClPd1W7Nv31rp9vvfXWo+fOxSd/2i+qqU3XvztT+sc/zxg3bpwLx/YsNQtY8+bN+/jjjzMzM6tvlFJmZmZaLBb7gnyzmcuOAADUk1OnTmVlZUVFRbVr1+4qjz3WDLEzX27IkV+fNnbly96hyh1R6vtxat+wy99VITAwMDk5+eDBP+/evTsgIOCfcXF8XrtGahaGunfv3rZt24SEhOobLRaLpmkjRoz46aefhgwZsmDBgoCAgCsdQQghpZw/f/5VdoiOjh46dOiVXrVarX5+ftymAVWsVmtZWZmrp4C70DSNtwSqs1qtFRUV/v7+tT7C/Pnzly5d6sCRHEXX9bNnz9qsVrOq6FKqZr/mzZtf9CetMES5Lqy6sOnCrIgAkwg0idYmcU4RK4RY4arRXWTgwIGzZ8+2Wq1CiKuE0WsKCAi45rfXLGANGjTo0o0Wi6VPnz5z586NiYl55plnnn766au/EaWUu3fvvsoOVqt18ODBV3q1srKysrLy+meG1+Mtgeo0TeMtgerq/pbIyckZPHjwQw895MCpUP82bdr09ddf298MiqLU5S3RoEGDa+7jgMt5ffr0+eabb+xfz549u2vXrlffX1XVf//737U+XWVlJTcaRXVlZWXBwRc/pB0+S9M0wzB4S6CKv7+/zWary1vC398/PDy8V69eDpwK9e/48eNmszk4OFhV1Xq40agDPkW4a9eubdu22b/29/e/nlgHAADgxeoUsDZt2iSEKC0tTUpKOnjwYEVFxcyZMxMTEx0zGgAAgGeqU8CyL8kaMGDAjBkzEhISoqKiCgsL58yZ46DZAAAAPFJt1mBVPR/a/oWiKNOmTZs2bZoj5wIAAPBY3MkdAADAwQhYAAAADkbAAgAAcDACFgAAgIMRsAAA8GaKolT9FfWGgAUAgPebPn26q0fwLQQsAAC83+zZs109gm8hYAEA4P2qLhR+8sknPXr0CA0Nfffdd4UQubm59957b4sWLdq0aTNhwgSLxSKEWL16dc+ePUNCQiIiIt5++237Ny5ZsiQsLMy1fwoP4oCHPQMA4IPKdXH/N7qrp/hVoFksHWS65m4nTpz44YcfNm7cOHLkyD/84Q9Tp04dP358cnKyzWabN2/e1KlT165dO2PGjIceeuiZZ57JyMjo16/fs88+K4TYunXrxo0bnf/n8BIELAAAasNPFRM7uNHKcfP1LWN//PHHFUUZPHiw1WoVQmzcuHHt2rVVr9o7qj179uzcufM///nP5s2bKyoq7C+9+uqrNFjXj4AFAEBtmBSR2MrzVto0bty4+t82bdo0IyOjTZs2QojS0tKzZ88KIe677z5/f/9x48bNmjUrOTnZvifpqkY8750BAED9KNPE6uPGxtPS1YM40ZgxY2bNmlVWVpaXlzd69OhZs2YJIb766qsXX3xx5MiRX3zxhRBC0zRXj+l5CFgAAPxGdrF870fj7i+08CWVfz9ghDRw9UDO9Prrr+u63rp1665du8bGxtqXtL/55pvx8fHdunU7e/bs0KFDJ06c6OoxPQ+XCAEAEJohtpyRn5001p6QBTY5vKU6paO6bIjaxE88+40bLbSqBSnlRX+tvr1x48YLFy686FuefPLJJ5980v71c889J4RYvHhx/UzrNQhYAADflWcVn58y1p2QG04b7ZooI1qqyfFq71BF9exMBdcjYAEAfIsUYne+XHdSrjth/HRO3hGlDm+pzI/zuyHQ1ZPBixCwAAA+obhSfJVjfHZSfnbSCPZXRrRUZt1sGhCu+LEaGU5AwAIAeKotW7akp6fbbLa+ffvefffdl32ecdZ5sfaE38Z87bs8eWsLZWSM+uce5rZN3PQSoM1mUxTF39/f1YOgrghYAADPU1FR8fDDD2//fHVCh3B/k/qn//e3OZ16rly5smnTpkIImy7Sci+sWLdq5jvDtWld1JV3qo3c9ZdeYWHh22+/vXTpUsvJE4aUbdp3mDBhwlNPPdWoUSNXj4Zactf3GgAAVzZ37twjaV9+89Bt/iZVCPF03zb/u37/40//8Y4X//XZSfnNaaNrU2VES3X5ELVjw3KbzRYS4r4LrA4cOJAwYvjwMNOi26NbBXcQQvx0tiR52Qf9lixZ99lnMTExrh7wNxRFqf5RxFrvef3H8VAELACA51m2bNmf+rWzpyshhKoo029rf8uS5WLi/42NNS3o79c84MKeVqvLhrweRUVFo0eOmNXnhttaNqva2CE06PX4zmuzckePGvXdjh2+dsXwjjvu2LBhg6unqCuW9gEAPM+ZM2daBf+mlAoPCvCrLPvnzeUPtlOr0pX7e+edd+4OM1VPV1VGtg9vX5n/4Ycf1v9UrvX111+7egQHIGABADxGmSZSjhkPbdJ/btTqYH5J9ZeOFJY2Dml60YP23N+yZcvu7Rx5pVf/p0vUsmXL6niKhQsXRkRENG/e/O9//7sQYvXq1T179gwJCYmIiLDft11RlBdeeCEiIuK111575ZVXOnToEBwc/Oabb9pfeuutt1q0aDFgwIBjx45VHTM3N/fee+9t0aJFmzZtJkyYYLFYhBAFBQUPPvhgaGhou3bt7Oeq7rKvXjpMYmKiEKJnz56XvuRZCFgAAHdXaBOLDhtjNugRSyo/OGjE3aC89cepb27J/LnMZt+hpEKbsfHglClTXDtnTWmadjz7aGxIwyvt0LZpo0OHDtXxLH/84x/Xr1+/bdu2NWvWCCFmzJgxfvz4s2fPfvbZZy+++KJ9n27duq1fv/7ll18OCQnZv3//p59+OnPmTPtLJSUlp0+fHjBgwDPPPFN1zKlTp95zzz3Hjx/fvXt327Ztp06dKoR4+umnhRDZ2dl79+7dvXv3RWNc9tVLh1m5cqUQ4ocffrjsnB6kvpeY2Wy2Jk2a2Gy2Wh8hPz8/ODjYz8/PgVPBo1ksloiICFdPAXehaVphYWFYWJirB4EDWMrEquNGyjHjuzw5OFJNilUSYtSmvzwZ8NVXX53/zl/7hQeZVfXbnKJR48a/9957ly5XslqtNpstJCSk1mM8++yz4eHhzz77bK2PcCWGYTQMDDz0aLz5CneOzy0pT/ryyMmTJ+tyltGjRxuGMWHChKSkJLPZbBjGzp07Dxw4sHnz5uTkZCmloig2m83Pz09V1crKSrPZLKVUVdX+0pEjR9q0aZOfn9+xY8ezZ8/aF6cHBQWVlpZWnSIsLCwvLy80NPTAgQPh4eFCiNzc3IiIiOoZ47KvXjqM+GX9+2VfqouVK1d+9NFHK1euLC4uVhQlKCiojge8Oha5AwDcy9FimXJMph4zDhXJu1uqj3ZSV96pNrzk99XLL788ZcqU7du3V1ZW/vWmm9q3b++KYetEVdVOnTtnni3pGnb5K5sH80u6detWx7OsXLly/fr1H3300YIFC7766qv77rvP399/3Lhxs2bNSk5Otu9TFUzNZrMQ4tI7iqmqqut61d82bdo0IyOjTZs2QojS0tKzZ8/a96na4bJHuPTVyw5zzZc8ApcIAQBuYV+BfHW30TNFu221dvicfKmXyfKg3+J409jWl0lXdlFRUffcc8+4ceM8MV3ZTZgw4T8ZJ670anLGiYcffriOp2jdunXr1q1feumlXbt2CSG++uqrF198ceTIkV988YUQQtO0q3/7v//9b03T/va3vw0cOLBq45gxY2bNmlVWVpaXlzd69OhZs2YJIYYPH/7ss88WFxeXlpb++c9/vug4l331SsNUVlbWdE53Q8ACALiMFOLbPDl9h97+v9ror/TzlfK9OFPOA34f9DcNjVb8feB31BNPPJHpF7bsQM6lL72386jSrud9991Xx1M8//zzt956a3x8/F//+lchxJtvvhkfH9+tW7ezZ88OHTp04sSJV//2ysrK8PDwb775Zt68eVUbX3/9dV3XW7du3bVr19jYWPsi9L/97W+GYcTGxnbv3v3222+/6DiXffWywwwfPrxt27Y1ndPdsAYLHo81WKiONVgeQTPE5lyZcsxYdVyG+IukWGVMrNor1CmPr3HnNVh2FoslMTGx2dnj426M6tK8sS7FvrxzyRknAzv1WbZsWV0mrztvuh0oa7AAAN7Jqon1OUbqMbn2hNG2iTImVv1muNoh2E0fC1hvIiIitm7dunjx4o8/+eTH9T+aTKbu3bs/9c6LY8eOvezTFeERCFgAAOc6VyE+O2mkHJNf5Rh9mitJserrN5mjGxEdfmU2mydOnOiGV8GWLl3q6hE8FQELAOAUeVax6riRetzYmitvj1CTYpUP+vuFNrj2N8J93H///a4ewVMRsAAAjnS8RKYek6nHjIwCOaylOrG9+t/BahDrZuFjCFgAAAf4sehCrjpRIke1Up/vYbojUmlgcvVYgIsQsAAAtSSF2JUvU7KN1OOyTBOJrZS3bzENCFdMLK+CzyNgAQBqRpciPVemHjNWHpcNzSKplbIo3tSnOR94A35FwAIAXBebLjaclqnHjNXHjZggJSlW/XyY2iWEWAVcBgELAHA1xZXi85NGyjH55SmjR6iS1Eqd0cvcKohcBVwNAQsAfNHJkyeFEC1btrzSDvnlYvUJI/WYkWaR/cOVpFh1fpxfWEA9juhO0tLSuOenp9u/f399no6ABQC+ZcWKFf/7v/9bWZCnKMIUEvb2229Xf9rdqVL52Um55oSRZpEDI5R7W6uL49VgfxfO63ovv/zyK6+8cvr0aVcPgjpp1qzZpU9IdB4CFgD4kJUrVz45cfz/u7vHzZGdhRC7LEWPTX5YVdXud45NOSZTjxlHi+XIGHVqR3X5EDWAmywIIYRo3LjxO++84+op4GEIWADgQ2bNmvXqwE43R154fnCfiJDX4zs/PP3NRuWjh0Wrf+mlDotW/VTXzgh4AwIWAPiQ/fv3xz3Ur/qWuOhmleu3nBnvxwojwIH4/xQA8BU/lwu1QaOi8srqGwvLK4MaNSRdAY5FwAIAL2fTxZoTxn1f6+2WVYb2uWvJvpPVX/14/8mhQ4e6ajbAW3GJEAC81vY8uSjL+O9Ro2eoMqG9unCgX9GQ2QP69y/4av/oDuFCiDVZudvPq+kr5rh6UsDbELAAwNscL5GLsuSiw4aqiAnt1T1jzC0bXbgGGBQdvW///nnz5i3eskVK2f/+xH8+/XTjxo1dOzDgfQhYAOAlrJpYe9JIzjK2n5FjW6sfDjDdFn6Zm2MGBQW9+OKLLpgP8CUELADwbLoUX5+Wi7KMtSeM+Ah1Ugd1xR2qPytsAZciYAGApzpYJJcdNf6TJRuaxIT26ju3+LUIdPVMAIQQBCwA8DiFNrE820jOMo6ViLGxyso7TT2acZsFwL0QsADAM9h0sfaEkXxYpucaCTHqq31MgyIUlWQFuCUCFgC4O/vdFpZnGz2aKRPaq0vi/YL8XD0TgKsiYAGAmzpVKpcclgt/MhQh7m+r7Bhtbt2YwgrwDAQsAHAvxZViRbbxnyxjf6H8nzbqonhT3zByFeBhCFgA4BYMKb4+LZN/udvCU13VETHcbQHwVAQsAHCxQ0XyE+62AHgXAhYAuEZ+uVh6xEjOMnKtYnw75bOhps4hXAoEvAQBCwDqlU0X604ayVlys8VIiFFn3WwaHMndFgBvQ8ACgHqyK18mZxnLjhrtmigT2quL4v0ac7cFwEsRsADAuXJK5afZcuFPhk0X97dVto0yt+FuC4C3I2ABgFPY77aQnGXsK5T/00Zd0N90SwtyFeArCFgA4EjV77Zwe4T6e+62APgkAhYAOIb9bgvJWTKQuy0APo+ABQB1UmgTy7ON5CzjWIkYG6uk3GHqGcqlQMDXEbAAoDYuutvCy71NQ7jbAoBfELAAoGYOFMpFh42PfuJuCwCuiIAFANfFfreFf/9klHO3BQDXQsACgKuxamLtSSM5y9h+Rg5vqb59i2lIlEKwAnB1BCwAuAxDim9Oy+QsY80JY2CE+kgH9dMhagOTq8cC4CEIWADwGxfdbeFt7rYAoOYIWAAgBHdbAOBQBCwAPq3CEF+eMhZlyS9PGUOj1ek91LujVTM3XgdQNwQsAD7quzy56LCx7KjRranycAf1w4HcbQGAwxCwAPiWEyVy8WGZnGUIIR5qr+5KNMcEcSkQgIMRsAD4hIvutvBeHHdbAOBEBCwA3syQYtuZC5cCbw5THmqnLhusNuQnHwAn48cMAO90sEgmZxmLD8vwQDGhvfr6fX5hAa6eCYDPIGAB8Cr55eKj7AYp27TTZWJ8O+XLu01dQrgSCKC+EbAAeIyNGze+++67mZmZkZGR99577+9+9zuT6cK91avdbUHGh5ln9OJuCwBciR8/ADzD3LlzHxo9fLAt+4NbwqaEFCe/8ZeEhATDMA4Uyj/t1GOWVs7Za9wRpRy7T/nXzaUJMaQrAK5EgwXAA+Tn57/20oyVY29q17SREKJd00b9opqO+GRL2z+tkL2T7m+rbBtlbtNYEUJomlHo6mkBgIAFwANs3769c0iAPV3Z+ZvUhPY3ZBekLb//XhcOBgCXRYcOwANUVFQE+l388yrQz9zMVOGSeQDg6ghYADxAr969fzhz/pytsvrGzcfz+/Tp46qRAOAqCFgA3N2xYvm7gy3N/R/83dofjhWVCSFKKrTX0g7lBoSOHz/e1dMBwGUQsAC4LynEgkNG31XanVHqibX/N/r30+/9MuvGD765ZfGOyh6DNm/e3LBhQ1fPCACXoUgp6/N8NputSZMmNput1kfIz88PDg728+Op97jAYrFERES4ego4XnaxnJKul2li4UBT52o3Cy0qKgoODlau8CBBTdMKCwvDwsLqa0y4O6vVarPZQkJCXD0I3EVxcbGiKEFBQU49Cw0WALdjL65uWaXdFaVuSTB3/u2t2ENCQq6UrgDATXCbBgDuJbtYTk7Ty3WRNtLciafcAPBMNFgA3EVVcTU0Wk1PIF0B8GA0WADcwtFiOTlNrzBEeoK5YzDRCoBno8EC4GKGFH/bb9yySkuKVdNHkq4AeAMaLACulHlOTkrT/VTx7Shz2yZEKwBeggYLgGsYUiw4ZPRfo42KUb8ZTroC4FVosAC4wJHzclKabgixLcHcnmuCALxObRosXdc7depUfUthYWFCQkKzZs1GjRpVWFjooNkAeCFdinn7jX6rteEt1c0jSFcAvFONA9a8efPi4uIyMzOrb5wzZ06rVq0sFktMTMxbb73luPEAeJUfi2Tcam3FMWP7KPP0HqpKuALgpWocsLp37z5jxoyLNqampk6bNq1BgwbTpk1LSUlx0GwAvIdmiDl7jfi12gNt1U0jzO1YcQXAq9XyWYSK8ptvDAoK+vnnnwMDA61W6w033HD+/PkrfaPNZmvUqNHVnxwXFxc3f/78K71aUFDQuHFjnkWIKmfOnLnhhhtcPQWuJrPY9NSehiF+cm7PsqhAw6nn0jTt3LlzoaGhTj0LPEh5ebnNZgsODnb1IHAXJSUliqI0atSo1kdo1qyZ2XyNVeyOWeQupbQ/GkxKqev61Xc2mUybNm26yg6BgYFX+eEopeRhz6iuoqKC36ZuSzPE3P3ybwfka72VKR0VRTRw+hk1TVVV3hKowsOecRF/f/86PuxZVa99AdAxASsyMvLkyZPt27fPycmJioq65v5t27at9blMv6j1EeBleD+4rf2F8pHNemiA+D7RHBNUT9cEpZS8JVAdvzVwEZPJpCiKs98Sdb0Plr2LSkhIWLhwoZRy4cKFo0ePdsBcADxZpSFe22MMXqc93kX9Ylj9pSsAcBN1bbAGDRokpXzppZcefPDBli1b9u7de9GiRQ6ZDICH+uGsfCRNj2wo9iSZoxoRrQD4oloGrKoV7vYvQkJC1q1b57ChAHgmzRDv7DPm7tdn9jH9rhMPigDgu7iTOwDHyCiQj6TpNwSKXYnmaIorAL6NgAWgrioNMZfiCgCqIWABqJOMAjlxsx7RUOxOZMUVAFxAwAJQSxRXAHAlBCwAtbG3QE7crEdRXAHA5RCwANQMxRUAXBMBC0ANfJcnJ6XpbZqIPUnmyIYUVwBweQQsANelXBev7NY/zDTeuIniCgCugYAF4Nq+zZOPpOndmio/3uMXFuDqaQDA7RGwAFyNVRMzdukfHzH+fqvpntYUVwBwXQhYAK5oe56clKZ3a6rsHUNxBQA1QMACcBlWTby6R0/OMubfahpLcQUANUTAAnCxbWfkpDS9ezMlY4xfc4orAKg5AhaAX9mLq0VZcn6cOiaW4goAaomABeCCrWfkpDS9RzMlY6w5tIGrpwEAT0bAAiDKNPHaHn1RlnwvTk2iuAKAOiNgAb5uS66clK73pLgCAMchYAG+q6q4+sdt6uhWFFcA4DAELMBHpefKSWl6r1Bl31hzM4orAHAoAhbgc+zF1eLD8h9x6iiKKwBwAgIW4FvW58ip6fotYcq+MeamFFcA4BwELMBXnK8Uf/xW/ypH/nOA6a4oxdXjAIA34+oA4BO+OCW7rdBMisgYayZdAYCz0WABXu5chXh+h74+R344wHQH0QoA6gUNFuDNvjglu6doQoiMMWbSFQDUGxoswDvZi6uvcvYnftoAABdESURBVOTCgaYhkUQrAKhXNFiAF/r8pOy2QhNC7B1jJl0BQP2jwQK8SlGFmL5D/ypHfnS7aTDRCgBchAYL8B6fnZTdV1xYcUW6AgAXosECvIG9uNqQI/8TbxoUQbQCABejwQI83rpqK65IVwDgDmiwAA9WVVwtijfFE60AwG3QYAGeau2JC8VVxlgz6QoA3AoNFuB5Cm3iTzv1r0/LxfGm24lWAOB+aLAAD7PmhNEt5cKKK9IVALgnGizAY+RZxZPb9AOF8tMhpn4tiFYA4L5osADPsOyo0SOlsm0TsTvJTLoCADdHgwW4uzyreGKbfrBQrrrL3DeMaAUAHoAGC3Bry7ONHimV7ZqI3UmkKwDwGDRYgJs6YxVPbNUzz8nVd5lvJloBgEehwQLc0fJso2dKZftgsSuRdAUAnocGC3AvuVbxxFY965xcM9R8U3OiFQB4JBoswI3Yi6sOweL7RNIVAHgwGizALeRaxeNb9CPFct1Qcx+iFQB4OBoswPXsxVXHEPF9IukKALwBDRbgSpYy8fhWPbtYfjbU3JtoBQDeggYLcJnl2UbP1MpOIWJnIukKALwKDRbgApYy8egW/XiJ/GKYuVco0QoAvA0NFlCvpBDJWUbP1MouTcXORNIVAHgnGiyg/hwvkVPT9Tyr+HKYuSfRCgC8Fw0WUB+kEAsOGb1Ttd7NlZ2JpCsA8HI0WIDTZRfLKel6qSbSE8xdQohWAOD9aLAAJ5JC/ONHo+8qbVi0upV0BQA+gwYLcJZjxXJyul6mibSR5s5EKwDwJTRYgOPZV1z1XaXdFaVuSSBdAYDPocECHCy7WE5O08t1kTbS3IloBQA+iQYLcBh7cXXLKm1otJqeQLoCAN9FgwU4xtFiOTlNrzBEeoK5YzDRCgB8Gg0WUFdVxdWwaDVtJOkKAECDBdTNkfNycrquGWJrgrkD0QoAIISgwQJqzZBiwSGj32rt7mg1bSTpCgDwKxosoDaOnJeT0nRDiG0J5vZEKwDAb9FgATVTVVwNb6luHkG6AgBcBg0WUAM/FslHNusNTGL7KHO7JkQrAMDl0WAB10UzxJy9xoA12phYddMI0hUA4GposIBrO1AoH0nTQ/zFniRzTBDRCgBwDTRYwNXYi6tB67QpHdUv7yZdAQCuCw0WcEX7C+Ujm/XQAPF9ItEKAFADNFjAZdiLq8HrtKmd1M+Hka4AADVDgwVcbF+BnJSmhwaIXUnmlo2IVgCAGiNgAb+qNMSbPxjv/6i/1dc0sQP9LgCglghYwAUZBfKRNP2GQLE7yRxNcQUAqAMCFiAqDTF3nzF3vz6zj+l3nSiuAAB1RcCCr7MXV+GBYlcixRUAwDEIWPBdFFcAACchYMFH7S2Qj2zWIxuK3YnmKIorAIBDEbDgcyiuAADORsCCb/kuT05K09s0EXuSzJENKa4AAE5BwIKvKNfFK7v1DzONN26iuAIAOBcBCz7h2zw5KU1v20TsHUNxBQBwOgIWvJy9uErOMv5+q+me1hRXAID6QMCCN9ueJyel6d2aKnvH+IUFuHoaAIDPIGDBO1k18eoePTnLmH+raSzFFQCgfhGw4IW2nZGT0vTuzZSMMX7NKa4AAPWOgAWvYi+uFmXJ+XHqmFiKKwCAaxCw4D22npGT0vQezZSMsebQBq6eBgDgwwhY8AZVxdV7cWoSxRUAwNUIWPB4Owr8ntui9aS4AgC4DQIWPFiZJl7bo390qPEHA9XEVhRXAAB3we8keKr0XNkzVTt6Xmy8vYh0BQBwKzRY8Dz24mrxYfmPOHVUK9ViMVw9EQAAv0HAgodJy5WT0/ReoUrGGHMzVlwBANwSAQse43ylePY7ff0puWCA6a4oHtgMAHBfrFyBZ/jylOy+QivXxJ4kM+kKAODmaLDg7s5Xiue+09fnyH8NMN1BtAIAeAIaLLi1L0/Jbis0IUTGGDPpCgDgKWiw4KbOVYjnd+jrc+TCgaYhkUQrAIAnqVmDVVhYmJCQ0KxZs1GjRhUWFlZ/KS4uTvnFY4895tAh4XM+Pym7p1workhXAACPU7OANWfOnFatWlkslpiYmLfeeqtqu5QyMzPTYrEUFxcXFxe/++67jp4TvqKoQjy6RX9ym/7vgab/629q7OfqgQAAqDlFSnn9e3fs2HHVqlWdOnU6dOjQ6NGjMzMz7dtPnz7duXPndu3a/fTTT0OGDFmwYEGLFi0uewSbzda4ceN//etfVzlLRETEgAEDrvTq2bNnmzRp4ufHL14v9EWO8vvv1KFRcnZvI+i6/w3n5uaGh4c7cy54Ek3TioqKmjdv7upB4C6sVmtFRUVwcLCrB4G7KC4uVhQlKCio1kfw9/dX1WtUVDVbg5WTk9OqVSshhL3HqtpusVj69Okzd+7cmJiYZ5555umnn166dOmVDiKlXLVq1VXOcuONN950001XetVqtfr5+WmaVqPJ4ebOVSqv7PPfnGd+/ybrgDBdaMJ63f+Gy8vLrVarM6eDJ9E0zWq18pZAlfLycpvN5u/v7+pB4C7Ky8sVRTGZTLU+wvW0PDVrsBo1anT27NmAgICysrKwsLDS0tJL97FYLF27di0oKLjsEWw2W5MmTWw22/Wf9CL5+fnBwcE0WN5k3Un52BZ9eEvlnVtM119cVbFYLBEREU6YCx5J07TCwsKwsDBXDwJ3YbVabTZbSEiIqweBu6h7g3U9atZgRUZGnjx5sn379jk5OVFRUVXbd+3aZbPZ4uLihBD+/v4NGvAEE1yXogoxfYe+IUcuijfFR7CYHQDgJWq2yD0hIWHhwoVSyoULF44ePVoIsWnTJiFEaWlpUlLSwYMHKyoqZs6cmZiY6IxZ4WXWnvjlHldjzaQrAIA3qVnAeumllzIyMlq2bHngwIG//OUvQohBgwYJIQYMGDBjxoyEhISoqKjCwsI5c+Y4ZVh4i0KbeHSL/odv9cXxpv/rb2rE7dgAAN6lZmuw6o41WFhzwnh8qzGipTK3n2OiFWuwUB1rsHAR1mDhIu64BguoizyrmLZN331WfjzINDCca4IAAK/FswhRT5ZnGz1SKiMaiowxZtIVAMC70WDB6fKs4olt+sFCueouc98wohUAwPvRYMG57MVVuyZidxLpCgDgK2iw4CxnrOKJrXrmObn6LvPNRCsAgC+hwYJTLM82eqZUtg8WuxJJVwAAn0ODBQfLtYontupZ5+SaoeabmhOtAAC+iAYLjmQvrjoEi+8TSVcAAN9FgwXHyLWKx7foR4rluqHmPkQrAIBvo8GCA9iLq44h4vtE0hUAADRYqBtLmXh8q55dLD8bau5NtAIAQAhBg4W6WJ5t9Eyt7BQidiaSrgAA+BUNFmrDUiYe3aIfL5FfDDP3CiVaAQDwGzRYqBkpxD8PGT1SKvu1UHYlkq4AALgMGizUwPESOTVdz7OK9XebexKtAAC4AhosXBcpxIJDRu9UrXdzZWci6QoAgKuhwcK1HSuWU7foxZViS4K5cwjRCgCAa6DBwtXYi6u+q7Q7ItWtpCsAAK4PDRau6FixnJyul2li80iiFQAANUCDhcuoKq7uilK5LAgAQE3RYOFi2cVycpperou0keZORCsAAGqOBgu/shdXt6zShkar6QmkKwAAaokGCxccLZaT0/QKQ6QnmDsGE60AAKg9Giz8WlwNi1bTRpKuAACoKxosX3fkvJycrmuG2Jpg7kC0AgDAEWiwfJchxYJDRr/V2t3RatpI0hUAAA5Dg+WjjpyXk9J0Q4htCeb2RCsAAByKBsvnVBVXw1uqm0eQrgAAcDwaLN/yY5F8ZLPewCS2jzK3a0K0AgDAKWiwfIVmiDl7jQFrtDGx6qYRpCsAAJyIBssnHCiUj6TpIf5iT5I5JohoBQCAc9FgeTl7cTVonTalo/rl3aQrAADqAw2WN9tfKB/ZrIcGiO8TiVYAANQfGizvZC+uBq/TpnZSPx9GugIAoF7RYHmhfQVyUpoeGiB2JZlbNiJaAQBQ32iwvIq9uBry2YXiinQFAIBL0GB5j4wC+UiafkOg2J1kjiZaAQDgOgQsb1BpiLn7jLn79Zl9TL/rRCsJAICLEbA8nr24Cg8UuxIprgAAcAsELA9GcQUAgHsiYHmqvQXykc16ZEOxO9EcRXEFAIA7IWB5HoorAADcHAHLw/xwVj6Spkc3EnuSzJENKa4AAHBHBCyPUa6LV3brH2Yab9xEcQUAgFsjYHmGb/PkpDS9bROxdwzFFQAA7o6A5e7sxVVyljHvVtO9rSmuAADwAAQst7Y9T05K07s1VfaO8QsLcPU0AADg+hCw3JRVE6/u0ZOzjPm3msZSXAEA4FEIWO5oS66clK7f3FzZN9YvtIGrpwEAADVEwHIvZZp44Xt9+VH5/m1qYiuKKwAAPBIBy41sOyMfSdN7NFP2jjE3Z8UVAAAei4DlFuwrrhZlyffi1KRYiisAADwbAcv1tuTKyel6j2ZKxlgzK64AAPACBCxXKtPEa3v0RVmsuAIAwKsQsFwmPVdOTtd7NlP2jTU3o7gCAMCLELBcwF5cLT4s349TR1NcAQDgdQhY9S09V05K03uFKhljKK4AAPBOBKz6c75SPPedvvaE/KC/mhBDcQUAgNfi13w9WZ8ju6/QrJrYP9ZMugIAwLvRYDmdvbhanyP/OcB0Z5Ti6nEAAIDTUaU415enZLcVmhAiY4yZdAUAgI+gwXKWcxXi+R36+hz54QDTHUQrAAB8CQ2WU3xxSnZPuVBcka4AAPA1NFgOZi+uvsqRCweahkQSrQAA8EU0WI70+clfV1yRrgAA8Fk0WI5RVCGm79A35Mj/xJsGRRCtAADwaTRYDvDZSdl9hSaE2DvGTLoCAAA0WHVSVVwlx5viiVYAAEAIQYNVF2tP/LLiaqyZdAUAAKrQYNVGoU38aaf+9Wm5ON50O9EKAAD8Fg1Wja05YXRLubDiinQFAAAuRYNVA3lW8ex3+rY8+fEg08BwohUAALg8GqzrtTzb6J5SGWgWGWPMpCsAAHAVNFjXlmcVT27TfyyUq+8y9w0jWgEAgGugwbqG5dlGj5TKtk3E7iTSFQAAuC40WFd0xiqe3KYfKpKr7zLfTLQCAADXjQbr8pZnGz1TKts1EbsSSVcAAKBmaLAudsYqHt+qZ52juAIAALVEg/Ub9hVXHYLF9xRXAACgtmiwLsi1ise36EeK5dqh5puaE60AAEDt0WAJ8cuKq44h4vtE0hUAAKgrX2+wLGXi8a360WK5bqi5D9EKAAA4gk83WMuzjV6plZ1CxPeJpCsAAOAwPtpgWcrEY1v1Y8Xys6Hm3kQrAADgUL7YYC3PNnqmVnYOETsTSVcAAMDxfKvBOl4ip6breVbxxTBzr1CiFQAAcApfabCkEAsOGb1Ttd7NlZ2JpCsAAOBEPtFgHS+RU9L1/HLx9XBzT6IVAABwMi9vsOzF1c0rtTsi1e8TSVcAAKA+eHODdaxYTknXSzWxeaS5cwjRCgAA1BPvbLDsxVXfVdqdUeqWBNIVAACoV17YYGUXyynpulUTaSPNnYhWAACg3nlVg2Uvrm5Zpd0VpaYnkK4AAIBreE+DdbRYTknTbYZITzB3DCZaAQAAl/GGBsteXPVbpQ2NVtNGkq4AAICLeViDdeLEie+++y42NrZnz55+fn5CiKPFcnKaXkFxBQAA3EaNG6zCwsKEhIRmzZqNGjWqsLDwmtsdpbCwcNy4cX06d5j/7OMP3T24a9euGzdttq+4GhatplNcAQAAt1HjgDVnzpxWrVpZLJaYmJi33nrrmtsdZfz48dretG2PDPg46aYN42/7c6egYcNHLEg/sjXBPL2HqhKuAACA26hxwEpNTZ02bVqDBg2mTZuWkpJyze0OceTIke3ffDV7SNdAs8m+5c42LR7oEHrXsY86UFwBAAA3U+M1WDk5Oa1atRJC2Puqa26/lK7rXbp0ucoO/fr1u6gD27VrV+fmjRuYfhMHe4YHbzxwID8/v6Z/BHiZgoIC+4I8QAihadq5c+cUhf/1wgXl5eU2m03TNFcPAndRUlKiKEp5eXmtjxASEmI2XyNB1ThgSSntP7mklLquX3P7pVRVXbx48VV2aNy4cXBwcPUt0dHRuaW2i3azFJe3iGlx0Z7wQWVlZbwNUEXTNMMweEugir+/v81m4y2BKqqqKooSFBRU6yOYTKZr7lPjgBUZGXny5Mn27dvn5ORERUVdc/ulFEXp3bt3jU7at29frXHo54fP3N3uBvuW8zZtyf5T/3ppLNUF/Pz8eBugiqIovCVQnT1z85ZAFT8/P/sPCqeepcZrsBISEhYuXCilXLhw4ejRo4UQmzZtuux2BzKbzYsWLZqx49Sfv/lxxaHT//g+e+jH2+6d/NjQoUMdeyIAAIC6U6SUNfqGoqKiBx98cO/evb179160aFFwcLCiKFLKS7df9tttNluTJk1stouv912PgoKCDz/8cPfu3dHR0WPGjLn11ltrcRB4H4vFEhER4eop4C40TSssLAwLC3P1IHAXVqvVZrOFhIS4ehC4i+Li4jpeIrweNQ5YdVSXgGU3e/bsBx54ICYmxoFTwaM99dRT8+bNY1Ez7E6fPv3hhx/OmDHD1YPAXezYsWP//v2TJk1y9SBwF6mpqYGBgcOGDXPqWTzvUTlLliw5c+aMq6eAG3n//fcNw3D1FHAX+fn5S5YscfUUcCMHDx5ct26dq6eAG0lPT9+xY4ezz+J5AQsAAMDNEbAAAAAcjIAFAADgYDW+D1YdKYoSFhY2YcKEWh+hoKBg5syZfB4EVRo2bDhx4kQWucPu/Pnz586dq8sPGXiZEydOWCwW3hKosm/fPj8/v8OHD9f6CK+//vo1P2xX358iFELs2bNn37599XxSAAAAh0hISGjatOnV93FBwAIAAPBurMECAABwMAIWAACAgxGwAAAAHIyABQAA4GCeFLB0Xe/UqZOrp4C7SElJufHGG0NCQgYOHPjTTz+5ehy43ueff96lS5eQkJAuXbqsX7/e1ePAXezfv79Ro0aungKuFxcXp/zisccec/bpPOZThPPmzfv444937NjhKQPDqbKzs3v06LFhw4Zu3bq99957K1eu3Lp1q6uHgivput68efNPP/100KBBKSkpTz/9dE5OjquHgusVFRUNGTJk9+7d/O7wcVLK5s2bHzhwICgoSAhhNpsDAgKcekaPabC6d+8+Y8YMV08Bd3H06NEHHnigb9++gYGBEydOzMzMdPVEcDFN05YsWTJ48ODS0tIGDRpwL2IIIQzDePjhh1944QVXDwLXs1gsmqaNGDEiIiJi/Pjx58+fd/YZPabBslMUDxsYzqbr+rRp01RVff/99109C1yvpKSkcePGiqJs2bIlLi7O1ePAxd54442ioqK//vWv/O7Arl27nnvuublz58bExDzzzDMVFRVLly516hk97D3HfySobv369dOnTx86dOjrr79uNtf3c5/gnkpLS+fNm5eamrpz505XzwJX2rBhw5tvvrl+/Xqz2czvDlRnsVi6du1aUFDg1LN4zCVCoDop5fPPPz9z5sxly5bNnj2bdIWjR49Onz5dCNGoUaPJkycfPHjQ1RPBxTZs2LBx40Y/Pz/7g0rtvaarh4LL7Nq1a9u2bfav/f39GzRo4OwzErDgkdLS0lavXr1mzZrIyMiSkpKSkhJXTwQXi4yMXLBgQXp6upRy2bJlvXr1cvVEcLHZs2fLXwghpJT9+/d39VBwmdLS0qSkpIMHD1ZUVMycOTMxMdHZZ+T/++GRNm3alJmZWf1Zm/T/Pi4gICA1NfUPf/hDdnZ2p06dFi5c6OqJALiRAQMGzJgxIyEh4dy5c8OHD58/f76zz8hlaQAAAAfjEiEAAICDEbAAAAAcjIAFAADgYAQsAAAAByNgAQAAOBgBCwAAwMEIWAAAAA5GwAIAAHAwAhYAAICDEbAAAAAcjIAFAADgYAQsAAAAByNgAQAAONj/B8TONZNBn/9iAAAAAElFTkSuQmCC" }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x = log.(range(1, exp(4), length = 10)) .+ 1 # uneven grid\n", "y = log.(x) # corresponding y points\n", "\n", "interp = LinearInterpolation(x, y)\n", "\n", "xf = log.(range(1, exp(4), length = 100)) .+ 1 # finer grid\n", "\n", "plot(xf, interp.(xf), label = \"linear\")\n", "scatter!(x, y, label = \"sampled data\", markersize = 4, size = (800, 400))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "At this point, `Interpolations.jl` does not have support for cubic splines with irregular grids, but there are plenty of other packages that do (e.g. [Dierckx.jl](https://github.com/kbarbary/Dierckx.jl) and [GridInterpolations.jl](https://github.com/sisl/GridInterpolations.jl))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Multivariate Interpolation\n", "\n", "Interpolating a regular multivariate function uses the same function" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "hide-output": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "interp_linear(3, 2) = 1.6094379124341003\n", "interp_linear(3.1, 2.1) = 1.6484736801441782\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "interp_cubic(3, 2) = 1.6094379124341\n", "interp_cubic(3.1, 2.1) = 1.6486586594237707\n" ] } ], "source": [ "f(x,y) = log(x+y)\n", "xs = 1:0.2:5\n", "ys = 2:0.1:5\n", "A = [f(x,y) for x in xs, y in ys]\n", "\n", "# linear interpolation\n", "interp_linear = LinearInterpolation((xs, ys), A)\n", "@show interp_linear(3, 2) # exactly log(3 + 2)\n", "@show interp_linear(3.1, 2.1) # approximately log(3.1 + 2.1)\n", "\n", "# cubic spline interpolation\n", "interp_cubic = CubicSplineInterpolation((xs, ys), A)\n", "@show interp_cubic(3, 2) # exactly log(3 + 2)\n", "@show interp_cubic(3.1, 2.1) # approximately log(3.1 + 2.1);" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "See [Interpolations.jl documentation](https://github.com/JuliaMath/Interpolations.jl#convenience-notation) for more details on options and settings" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Linear Algebra" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Standard Library\n", "\n", "The standard library contains many useful routines for linear algebra, in\n", "addition to standard functions such as `det()`, `inv()`, `factorize()`, etc.\n", "\n", "Routines are available for\n", "\n", "- Cholesky factorization \n", "- LU decomposition \n", "- Singular value decomposition, \n", "- Schur factorization, etc. \n", "\n", "\n", "See [here](https://docs.julialang.org/en/v1/stdlib/LinearAlgebra/) for further details" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## General Tools" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### LaTeXStrings.jl\n", "\n", "When you need to properly escape latex code (e.g. for equation labels), use [LaTeXStrings.jl](https://github.com/stevengj/LaTeXStrings.jl)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "hide-output": false }, "outputs": [ { "data": { "text/latex": [ "an equation: $1 + \\alpha^2$" ], "text/plain": [ "L\"an equation: $1 + \\alpha^2$\"" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "using LaTeXStrings\n", "L\"an equation: $1 + \\alpha^2$\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### ProgressMeter.jl\n", "\n", "For long-running operations, you can use the [ProgressMeter.jl](https://github.com/timholy/ProgressMeter.jl) package\n", "\n", "To use the package, you simply put a macro in front of `for` loops, etc.\n", "\n", "From the documentation" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "hide-output": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "\r", "\u001b[32mComputing... 14%|█████▌ | ETA: 0:00:07\u001b[39m" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\r", "\u001b[32mComputing... 38%|██████████████▉ | ETA: 0:00:04\u001b[39m" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\r", "\u001b[32mComputing... 58%|██████████████████████▋ | ETA: 0:00:02\u001b[39m" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\r", "\u001b[32mComputing... 78%|██████████████████████████████▍ | ETA: 0:00:01\u001b[39m" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\r", "\u001b[32mComputing... 98%|██████████████████████████████████████▎| ETA: 0:00:00\u001b[39m" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\r", "\u001b[32mComputing...100%|███████████████████████████████████████| Time: 0:00:05\u001b[39m\n" ] } ], "source": [ "using ProgressMeter\n", "\n", "@showprogress 1 \"Computing...\" for i in 1:50\n", " sleep(0.1) # some computation....\n", "end" ] } ], "metadata": { "filename": "general_packages.rst", "kernelspec": { "display_name": "Julia 1.2", "language": "julia", "name": "julia-1.2" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.2.0" }, "title": "General Purpose Packages" }, "nbformat": 4, "nbformat_minor": 2 }