{ "cells": [ { "cell_type": "markdown", "id": "968e63dc", "metadata": {}, "source": [ "# Exercise sheet\n", "\n", "Some general remarks about the exercises:\n", "* For your convenience functions from the lecture are included below. Feel free to reuse them without copying to the exercise solution box.\n", "* For each part of the exercise a solution box has been added, but you may insert additional boxes. Do not hesitate to add Markdown boxes for textual or LaTeX answers (via `Cell > Cell Type > Markdown`). But make sure to replace any part that says `YOUR CODE HERE` or `YOUR ANSWER HERE` and remove the `raise NotImplementedError()`.\n", "* Please make your code readable by humans (and not just by the Python interpreter): choose informative function and variable names and use consistent formatting. Feel free to check the [PEP 8 Style Guide for Python](https://www.python.org/dev/peps/pep-0008/) for the widely adopted coding conventions or [this guide for explanation](https://realpython.com/python-pep8/).\n", "* Make sure that the full notebook runs without errors before submitting your work. This you can do by selecting `Kernel > Restart & Run All` in the jupyter menu.\n", "* For some exercises test cases have been provided in a separate cell in the form of `assert` statements. When run, a successful test will give no output, whereas a failed test will display an error message.\n", "* Each sheet has 100 points worth of exercises. Note that only the grades of sheets number 2, 4, 6, 8 count towards the course examination. Submitting sheets 1, 3, 5, 7 & 9 is voluntary and their grades are just for feedback.\n", "\n", "Please fill in your name here:" ] }, { "cell_type": "code", "execution_count": 1, "id": "b388fbcf", "metadata": {}, "outputs": [], "source": [ "NAME = \"Kees van Kempen\"\n", "NAMES_OF_COLLABORATORS = \"Bart Steeman\"" ] }, { "cell_type": "markdown", "id": "8a01a93a", "metadata": {}, "source": [ "---" ] }, { "cell_type": "markdown", "id": "41d26cde", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "markdown", "checksum": "3239a87542ad0d212698947a7b60cbb3", "grade": false, "grade_id": "cell-34662207849ac9d6", "locked": true, "schema_version": 3, "solution": false, "task": false } }, "source": [ "**Exercise sheet 6**\n", "\n", "Code from the lecture" ] }, { "cell_type": "code", "execution_count": 2, "id": "cb41d2a1", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "34c958a728224ac04d230ccc781dba94", "grade": false, "grade_id": "cell-be52f3b5932e9c30", "locked": true, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pylab as plt\n", "\n", "rng = np.random.default_rng()\n", "%matplotlib inline\n", "\n", "def aligned_init_config(width):\n", " '''Produce an all +1 configuration.'''\n", " return np.ones((width,width),dtype=int)\n", "\n", "def plot_ising(config,ax,title):\n", " '''Plot the configuration.'''\n", " ax.matshow(config, vmin=-1, vmax=1, cmap=plt.cm.binary)\n", " ax.title.set_text(title)\n", " ax.set_yticklabels([])\n", " ax.set_xticklabels([])\n", " ax.set_yticks([])\n", " ax.set_xticks([])\n", "\n", "from collections import deque\n", "\n", "def neighboring_sites(s,w):\n", " '''Return the coordinates of the 4 sites adjacent to s on an w*w lattice.'''\n", " return [((s[0]+1)%w,s[1]),((s[0]-1)%w,s[1]),(s[0],(s[1]+1)%w),(s[0],(s[1]-1)%w)]\n", "\n", "def cluster_flip(state,seed,p_add):\n", " '''Perform a single Wolff cluster move with specified seed on the state with parameter p_add.'''\n", " w = len(state)\n", " spin = state[seed]\n", " state[seed] = -spin \n", " cluster_size = 1\n", " unvisited = deque([seed]) # use a deque to efficiently track the unvisited cluster sites\n", " while unvisited: # while unvisited sites remain\n", " site = unvisited.pop() # take one and remove from the unvisited list\n", " for nbr in neighboring_sites(site,w):\n", " if state[nbr] == spin and rng.uniform() < p_add:\n", " state[nbr] = -spin\n", " unvisited.appendleft(nbr)\n", " cluster_size += 1\n", " return cluster_size\n", "\n", "def wolff_cluster_move(state,p_add):\n", " '''Perform a single Wolff cluster move on the state with addition probability p_add.'''\n", " seed = tuple(rng.integers(0,len(state),2))\n", " return cluster_flip(state,seed,p_add)\n", "\n", "def compute_magnetization(config):\n", " '''Compute the magnetization M(s) of the state config.'''\n", " return np.sum(config)\n", "\n", "def run_ising_wolff_mcmc(state,p_add,n):\n", " '''Run n Wolff moves on state and return total number of spins flipped.'''\n", " total = 0\n", " for _ in range(n):\n", " total += wolff_cluster_move(state,p_add)\n", " return total\n", "\n", "def sample_autocovariance(x,tmax):\n", " '''Compute the autocorrelation of the time series x for t = 0,1,...,tmax-1.'''\n", " x_shifted = x - np.mean(x)\n", " return np.array([np.dot(x_shifted[:len(x)-t],x_shifted[t:])/len(x) \n", " for t in range(tmax)])\n", "\n", "def find_correlation_time(autocov):\n", " '''Return the index of the first entry that is smaller than \n", " autocov[0]/e or the length of autocov if none are smaller.'''\n", " smaller = np.where(autocov < np.exp(-1)*autocov[0])[0]\n", " return smaller[0] if len(smaller) > 0 else len(autocov)" ] }, { "cell_type": "markdown", "id": "3317e002", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "markdown", "checksum": "5fc7aaa40ff9892f4cefca6bedcc5fd1", "grade": false, "grade_id": "cell-03072464d602c105", "locked": true, "schema_version": 3, "solution": false, "task": false } }, "source": [ "## MCMC simulation of the XY model\n", "\n", "**(100 Points)**\n", "\n", "_Goal of this exercise_: Practice implementing MCMC simulation of the XY spin model using an appropriate cluster algorithm and analyzing the numerical effectiveness.\n", "\n", "The **XY model** is a relative of the Ising model in which the discrete $\\pm 1$ spin at each lattice site is replaced by a continuous $2$-dimensional spin on the unit circle \n", "\n", "$$\n", "S_1 = \\{(x,y)\\in\\mathbb{R}^2: x^2+y^2=1\\}.\n", "$$\n", "\n", "To be precise, we consider a $w\\times w$ lattice with periodic boundary conditions and a XY configuration $s = (s_1,\\ldots,s_N) \\in \\Gamma = S_1^N$, $N=w^2$, with Hamiltonian that is very similar to the Ising model,\n", "\n", "$$\n", "H_{XY}(s) = - J \\sum_{i\\sim j} s_i \\cdot s_j.\n", "$$\n", "\n", "Here, as in the Ising model, the sum runs over nearest neighbor pairs $i$ and $j$ and $s_i \\cdot s_j$ is the usual Euclidean inner product of the vectors $s_i,s_j \\in S_1$. We will only consider the ferromagnetic XY model and set $J=1$ in the remainder of the exercise. Note that nowhere in the definition the $x$- and $y$-components of the spins are related to the two directions of the lattice (one could also have studied the XY model on a one-dimensional or three-dimensional lattice and the spins would still have two components). As usual we are interested in sampling configurations $s\\in \\Gamma$ with distribution $\\pi(s)$ given by the Boltzmann distribution\n", "\n", "$$ \n", "\\pi(s) = \\frac{1}{Z_{XY}} e^{-\\beta H_{XY}(s)}, \\qquad \\beta = 1/T.\n", "$$\n", "\n", "The XY model admits a (local) **cluster algorithm** that is very similar to the Wolff algorithm of the Ising model. It amounts to the following recipe:\n", "1. Sample a uniform seed site $i_{\\text{seed}}$ in $1,\\ldots,N$ and an independent uniform unit vector $\\hat{n} \\in S_1$.\n", "2. Grow a cluster $C$ starting from the seed $i_{\\text{seed}}$ consisting only of sites $j$ whose spin $s_j$ is \"aligned\" with the seed, in the sense that $s_j\\cdot\\hat{n}$ has the same sign as $s_{i_{\\text{seed}}}\\cdot \\hat{n}$, or $(s_j\\cdot\\hat{n})(s_{i_{\\text{seed}}}\\cdot \\hat{n})>0$. Like in the Ising model this is done iteratively by examining the neighbors of sites that are already in the cluster, and adding those that are aligned with appropriate probability. The difference with the Ising model is that this probability depends on the spins $s_i$ and $s_j$ that are linked (meaning that $s_j$ is an aligned neighbor of $s_i$) via the formula\n", "$$ p_{\\text{add}}(s_i,s_j) = 1 - \\exp\\big( -2\\beta\\,(s_i\\cdot\\hat{n})(s_j\\cdot\\hat{n})\\big).$$\n", "3. Once the cluster $C$ is constructed, all of its spins are \"flipped\" in the sense that they are reflected in the plane perpendicular to $\\hat{n}$, i.e. $s_j \\to s_j - 2(s_j\\cdot\\hat{n})\\hat{n}$.\n", "\n", "__(a)__ Verify by a calculation, to be included using markdown and LaTeX below, that the probabilities $p_{\\text{add}}(s_i,s_j)$ are the appropriate ones to ensure detailed balance for the Boltzmann distribution $\\pi(s)$. _Hint_: Follow the same reasoning [as for the Ising model](https://hef.ru.nl/~tbudd/mct/lectures/cluster_algorithms.html#cluster-algorithm-for-the-ising-model-the-wolff-algorithm). Compare the probabilities involved in producing the cluster $C$ in state $s$ and state $s'$. Why do the probabilities only differ at the boundary edges in the cluster $C$? **(25 pts)**" ] }, { "cell_type": "markdown", "id": "c1561680", "metadata": { "deletable": false, "nbgrader": { "cell_type": "markdown", "checksum": "1be75e4b28f44fd0e2a702c3f586003d", "grade": true, "grade_id": "cell-55323d84e19389a2", "locked": false, "points": 25, "schema_version": 3, "solution": true, "task": false } }, "source": [ "The condition for detailed balance is $$\\frac{P(s\\rightarrow s')}{P(s'\\rightarrow s)}=\\frac{\\pi(s')}{\\pi(s)}.$$\n", "Let's work to that expression by following the idea of the proof for the Ising model proof.\n", "\n", "Starting from the left-hand side, consider the same cluster from the same seed and unit vector $\\hat{n}$ but between different states $s$ and $s'$. Flipping between $s$ and $s'$ is done by the same cluster $C$, but from opposite sign. The probability of this happening is given by the boundaries of the clusters, namely that they are not added to the cluster, with probability $1 - p_{add}$ per boundary site. This gives a ratio $$\\frac{P(s\\rightarrow s')}{P(s'\\rightarrow s)} = \\frac{\\prod_{^+} (1-p_{add}(s_i, s_j))}{\\prod_{^-} (1-p_{add}(s_i, s_j))}$$ for pairs $^\\pm$ with $s_i \\in C$ and $s_j$ neighbouring $s_i$ just outside the cluster with $\\pm(s_i \\cdot \\hat{n})(s_j \\cdot \\hat{n}) > 0$. Substituting the formula for $p_{add}$ in, we find $$\\frac{P(s\\rightarrow s')}{P(s'\\rightarrow s)} = \\frac{\\prod_{^+} \\exp(-2\\beta(s_i \\cdot \\hat{n})(s_j \\cdot \\hat{n}))}{\\prod_{^-} \\exp(2\\beta(s_i \\cdot \\hat{n})(s_j \\cdot \\hat{n}))} = \\prod_{} \\exp(-2\\beta(s_i \\cdot \\hat{n})(s_j \\cdot \\hat{n})) = \\exp(-2\\beta\\sum_{} (s_i \\cdot \\hat{n})(s_j \\cdot \\hat{n}))$$ with $$ describing both types of pairs.\n", "\n", "For the right-hand side, $$\\frac{\\pi(s')}{\\pi(s)} = e^{\\beta \\left[H_{XY}(s)-H_{XY}(s')\\right]}$$ just from the definition of $\\pi$.\n", "This difference in energy we can write as\n", "$$\n", "H_{XY}(s)-H_{XY}(s') = \\sum_{} \\left[ s_i \\cdot s_j - s'_i \\cdot s_j \\right]\n", " = 2\\sum_{} (s_i \\cdot \\hat{n})(s_j \\cdot \\hat{n})\n", "$$\n", "by summing over all pairs $$ as above, as pairs not across the boundary do not give rise to an energy change. The last step is due to the flips from both signs to be accounted for.\n", "This shows $$\\frac{P(s\\rightarrow s')}{P(s'\\rightarrow s)}=\\frac{\\pi(s')}{\\pi(s)},$$ thus that detailed balance holds." ] }, { "cell_type": "markdown", "id": "66a9278e", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "markdown", "checksum": "696630fba38d4eedeebe2daa397e02e1", "grade": false, "grade_id": "cell-bf935a4be9ad6e06", "locked": true, "schema_version": 3, "solution": false, "task": false } }, "source": [ "__(b)__ In order to implement the cluster update described above, we take the state to be described by a Numpy array of dimension $(w,w,2)$, for which we have already provided a function `xy_aligned_init_config` to generate an all-aligned initial state. Write the function `xy_cluster_flip`, that grows and flips a cluster starting from the given seed site and $\\hat{n}$ and returns the cluster size, and `xy_cluster_move`, that performs the previous function to a random seed and direction $\\hat{n}$ and also returns the cluster size. **(20 pts)**" ] }, { "cell_type": "code", "execution_count": 3, "id": "6e4a6d63", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "0732e9388e50a79370df1d621764a8a8", "grade": false, "grade_id": "cell-4d68efbe5b0730b1", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [], "source": [ "def xy_aligned_init_config(width):\n", " '''Return an array of dimension (width,width,2) representing aligned spins in x-direction.'''\n", " return np.dstack((np.ones((width,width)),np.zeros((width,width))))\n", "\n", "def xy_cluster_flip(state,seed,nhat,beta):\n", " '''Perform a cluster move with specified seed and vector nhat on the state at temperature beta.'''\n", " w = len(state)\n", " \n", " # Let's flip along the way, starting with the seed.\n", " state[seed] -= 2*state[seed]@nhat*nhat\n", " # Use a double-ended queue to keep track of what elements to try to expand from.\n", " unvisited = deque([seed])\n", " cluster_size = 1\n", " while unvisited:\n", " site = unvisited.pop()\n", " # As the site is already added to the queue, it has already been flipped\n", " # so we need to take a minus into account.\n", " s_i_n = -state[site]@nhat\n", " for nbr in neighboring_sites(site,w):\n", " s_j_n = state[nbr]@nhat\n", " one_minus_p_add = np.exp(-2*beta*s_i_n*s_j_n)\n", " if np.sign(s_i_n) == np.sign(s_j_n) and rng.uniform() > one_minus_p_add:\n", " state[nbr] -= 2*s_j_n*nhat\n", " unvisited.appendleft(nbr)\n", " cluster_size += 1\n", " return cluster_size\n", "\n", "def xy_cluster_move(state,beta):\n", " '''Perform a single Wolff cluster move on the state with addition probability p_add.'''\n", " seed = tuple(rng.integers(0,len(state),2))\n", " phi = 2*np.pi*rng.uniform()\n", " nhat = np.array([np.cos(phi), np.sin(phi)])\n", " return xy_cluster_flip(state,seed,nhat,beta)" ] }, { "cell_type": "code", "execution_count": 4, "id": "3f810041", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "8abc44b0a19bcdb3c66ac52abd25b747", "grade": true, "grade_id": "cell-2e1ebc198d10ea5b", "locked": true, "points": 15, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "from nose.tools import assert_almost_equal\n", "assert 1 <= xy_cluster_flip(xy_aligned_init_config(4),(0,0),np.array([np.cos(0.5),np.sin(0.5)]),0.5) <= 16\n", "assert_almost_equal(np.mean([xy_cluster_flip(xy_aligned_init_config(3),(0,0),\n", " np.array([np.cos(0.5),np.sin(0.5)]),0.3) \n", " for _ in range(200)]),5.3,delta=0.7)\n", "assert_almost_equal(np.mean([xy_cluster_flip(xy_aligned_init_config(3),(1,2),\n", " np.array([np.cos(0.2),np.sin(0.2)]),0.2) \n", " for _ in range(200)]),4.3,delta=0.6)" ] }, { "cell_type": "code", "execution_count": 5, "id": "bde560e4", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "6524b7571176499732d477ffb6b7e40d", "grade": true, "grade_id": "cell-9445d8ff6dc2a16f", "locked": true, "points": 5, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "assert 1 <= xy_cluster_move(xy_aligned_init_config(4),0.5) <= 16\n", "assert_almost_equal(np.mean([xy_cluster_move(xy_aligned_init_config(3),0.3) \n", " for _ in range(200)]),3.6,delta=0.75)\n", "assert_almost_equal(np.mean([xy_cluster_move(xy_aligned_init_config(3),0.9) \n", " for _ in range(200)]),6.3,delta=0.75)" ] }, { "cell_type": "markdown", "id": "daa0b50e", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "markdown", "checksum": "4289e943aa451475a6d9402abbe6bcb1", "grade": false, "grade_id": "cell-aaed2919abf002dc", "locked": true, "schema_version": 3, "solution": false, "task": false } }, "source": [ "__(c)__ Estimate and plot the average cluster size in equilibrium for a 25x25 lattice ($w=25$) for the range of temperatures $T = 0.5,0.6,\\ldots,1.5$. It is not necessary first to estimate the equilibration time: you may start in a fully aligned state and use 400 moves for equilibration and 1000 for estimating the average cluster size. It is not necessary to estimate errors for this average. Store your averages in the data set `\"cluster-size\"` (an array of size 11) in the HDF5-file `xy_data.hdf5`, just like you did in Exercise sheet 5. Then read the data from file and produce a plot. **(20 pts)**" ] }, { "cell_type": "code", "execution_count": 6, "id": "1f91dd67", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "4d6b29807b163b30a89dbe9668117083", "grade": false, "grade_id": "cell-6e4765df49afeb66", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "T = 0.5\n", "\tequilibrating took 5.3386430740356445 seconds\n", "\tmeasuring took 13.652144193649292 seconds\n", "\n", "T = 0.6\n", "\tequilibrating took 5.194417953491211 seconds\n", "\tmeasuring took 12.401273727416992 seconds\n", "\n", "T = 0.7\n", "\tequilibrating took 4.591171503067017 seconds\n", "\tmeasuring took 11.19827151298523 seconds\n", "\n", "T = 0.8\n", "\tequilibrating took 3.9237494468688965 seconds\n", "\tmeasuring took 10.341700792312622 seconds\n", "\n", "T = 0.9\n", "\tequilibrating took 3.6344528198242188 seconds\n", "\tmeasuring took 9.163604021072388 seconds\n", "\n", "T = 1.0\n", "\tequilibrating took 2.604698896408081 seconds\n", "\tmeasuring took 6.472257137298584 seconds\n", "\n", "T = 1.1\n", "\tequilibrating took 1.5426356792449951 seconds\n", "\tmeasuring took 3.485330104827881 seconds\n", "\n", "T = 1.2000000000000002\n", "\tequilibrating took 0.4475827217102051 seconds\n", "\tmeasuring took 1.3380379676818848 seconds\n", "\n", "T = 1.3\n", "\tequilibrating took 0.294419527053833 seconds\n", "\tmeasuring took 0.7950901985168457 seconds\n", "\n", "T = 1.4\n", "\tequilibrating took 0.2504415512084961 seconds\n", "\tmeasuring took 0.46738409996032715 seconds\n", "\n", "T = 1.5\n", "\tequilibrating took 0.14740800857543945 seconds\n", "\tmeasuring took 0.39183998107910156 seconds\n", "\n" ] } ], "source": [ "temperatures = np.linspace(0.5,1.5,11)\n", "width = 25\n", "equilibration_moves = 400\n", "measurement_moves = 1000\n", "\n", "import h5py\n", "import time\n", "\n", "t_0 = time.time()\n", "def toc():\n", " \"\"\"Return how many seconds the function was last called.\"\"\"\n", " \n", " if not 't_0' in globals():\n", " return -np.inf\n", " \n", " global t_0\n", " dt = time.time() - t_0\n", " t_0 = time.time()\n", " return dt\n", "\n", "with h5py.File(\"xy_data.hdf5\", \"a\") as f:\n", " if not \"cluster-size\" in f:\n", " state = xy_aligned_init_config(width)\n", " cluster_sizes = np.zeros(len(temperatures))\n", " \n", " for idx, T in enumerate(temperatures):\n", " toc()\n", " print(\"T =\", T)\n", " beta = 1/T\n", " \n", " # Equilibrate\n", " for _ in range(equilibration_moves):\n", " xy_cluster_move(state,beta)\n", " print(\"\\tequilibrating took\", toc(), \" seconds\")\n", " \n", " # Measure\n", " for _ in range(measurement_moves):\n", " cluster_sizes[idx] += np.norm(xy_cluster_move(state,beta))\n", " cluster_sizes[idx] /= measurement_moves\n", " print(\"\\tmeasuring took\", toc(), \" seconds\")\n", " print()\n", " \n", " f.create_dataset(\"cluster-size\",data=cluster_sizes)" ] }, { "cell_type": "code", "execution_count": 7, "id": "f5ea3143", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "77f0ca1931f96503b4141e3789bb1d68", "grade": true, "grade_id": "cell-69f351f86b801f63", "locked": true, "points": 10, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "with h5py.File('xy_data.hdf5','r') as f:\n", " assert f[\"cluster-size\"][()].shape == (11,)\n", " assert_almost_equal(f[\"cluster-size\"][4],225,delta=40)\n", " assert_almost_equal(f[\"cluster-size\"][10],8,delta=8)" ] }, { "cell_type": "code", "execution_count": 8, "id": "e84152c5", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "dda3e8db7a6a4b7a7d6c58a3a03e487d", "grade": true, "grade_id": "cell-88fb1bb0cd597195", "locked": false, "points": 10, "schema_version": 3, "solution": true, "task": false } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEtCAYAAAAbeVcBAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABAkElEQVR4nO3dd3wc1bn/8c9XxZLcJBdZ7r1hDBhsWugYgyEQ05IAISEkgZDAhUAS2g/S6KSSm0BCC3BDINyEUBya4WJ6s4lxNxj3glxwL7IsPb8/zgjW8koayVqtVnrer9e+tDv1mdXMPnNmzpwjM8M555yrS1a6A3DOOZcZPGE455yLxROGc865WDxhOOeci8UThnPOuVg8YTjnnIvFE4ZzzrlYPGE455yLxROG2yOSFkk6Lt1xxCVplqSjG3mZwyT9R9ImSZc25rIT1pFR33NLJukBSTfGnLZF/d88YaSIpMmS1knKS3cszV1THlRmtreZTW7kxV4JTDazDmb2+z1dWEv7kYlDUvtou89JGNZB0hJJZ0p6WNL91eY5StJaST2aPuLWyRNGCkjqDxwBGPClFCw/p7GXmamayXfRD5jVkBmbSfy1UrB/kuH7ScpujHWY2WbgQuAOScXR4NuBKWb2D+BS4CRJ46J15wP3AD80s5WNEYOLwcz81cgv4CfAG8BvgIkJw68G/lFt2juA30fvewL/BFYDC4FLE6ZbBFwFTAfKgJxoeR8Dm4DZwGkJ0x8A/Cca97/A34EbE8bXuK4k29MHeDyadi3wh2pxHRe9N2BwwrgHqq3zKmB5FNM8YCzwP0AlsA3YDFzZkO8iScy7rStJvF+N1ln1KiOUFGJ/P8D/ARXA9mgZQ4G9gMnAekIi+VK1eWqMv5bvYxHwo2ieDdH/M7+B/8/NQJ/o/UXR/60k+vxj4N5q0w+Iljs+YdgRwBpg5J4uv9q8DwCPAEcT9rUeCeO+HG1bO+AW4NlalrMoWtd0YAtwH1ACPBvtEy8CnRKmr/F/BuwPvB/N93fgUWIeSyTsby3hlfYAWuILmA98HxgNlCccLP2ArUDH6HM2sBI4hFDam0pINm2AgcAC4IRo2kXANMKPd0E07MvRzppF+PHbAvSI5l8MXAbkAqcDO6p28rrWVW1bsoEPgN9GB2o+cHjC+M8OCGpJGMAwYCnQM/rcHxhUfRlx4kv2XVSLOfa6EubpCMwBvluf7yeadzLwneh9bvT/vzaa91jCD82wat9ZbfHvFmM07N3o/905ivWi+v4/o+mXAiMAATOAj4Dh0ef5wH5J5jmS8AN+DHAQ1RLIni4/Yd5OhGNiDXB+kvH/AJ6KYulby3IWAW8TkkQvYBXhR39/II+Q6H9a1/+Mz4+ly6PpziQc07GOpZr2t0x9pT2AlvYCDo92qK7R57nA5QnjXwe+Eb0fB3wcvT8YWFJtWdcAf4neLwK+Vce6pwETooN7OaBq670xzrqqDT80+nHY7Sw+Ia44CWNwdNAeB+TWtIzG+C7qs65oWBYwEbirvt9PNG4ynyeMI4BPgKyE8Y8AP6sWQ23xJ4txEXBuwufbgT81MN4ZhB/9E4B/Ra9DgBOBV2qJ6wTCD/UqEkqzjbX8hPlfJJxYFSYZV0IowVxWxzIWAV9L+PzPqv9v9Pm/gCfq+p8RjqUV7HosvUnMYynZ/zKTX34Po/GdB7xgZmuiz3+LhpHw+ezo/TnRZwilj56S1le9CGc8JQnzLk1ckaRvSJqWMP1IoCvhLHS5RXtsknnjrKtKH2Cxme2se9NrZmbzgR8QDsJVkh6V1LOGyev9XezBugBuAjoQrpPHXX9NegJLzawyYdhiwlluohrjr8UnCe+3Au0bGO+6aN4fEC6JbiSc2X8fqO2m/RJgJ6GksKiW6Rq6fCSdSygRvgjcVn28mZUSSh9x7hmVJrzfluRz1fdX2/8s2bG0OOH9nuwrGafZ33DLJJIKgK8A2ZKqDu48oEjSfmb2AeF+wq8l9QZOI5zBQ/gBWWhmQ2pZxWc7raR+hJt+Y4G3zKxC0jTCwbwS6CVJCTt6H8L9jrjrqrIU6CspJ0bS2Aq0TfjcHVj2WfBmfwP+Jqkj8GfCD8LXE7erHvFVn2fXkTWvaxeSziIk8APNrLwe66/JCqCPpKyEH6C+wIf1iL/WbUuivvGuB8YQ7g9MlnQGsB/hhOOJZDNIGgRMItx72QQ8I+k4M0v2w13v5Ufr6Ea49PkVQsl8lqS/mdmrMberoWr7nyU7lvrSsGMp43kJo3GdSrgBOgIYFb32Al4DvgFgZqsJlzD+QtjR5kTzvgtslHSVpAJJ2ZJGSjqwhnW1I/ywrAaQdD7hgAR4K4rjEkk5kiYQLhFUqc+63iUcNLdKaicpX9JhNcQ0DTgnWt544KiqEdGzCsdG1Yy3E87wKqLRpYRrvw2Jbzd1rCtxuv2B/wZOjf4vjbH+dwj3kq6UlKvwzMcphBulcVX/PupS33jXEa7JV53tbyTc7/qzmSX7nnoCLwE3mdkDZvZPwg34FyQli7Ney0/wB8Jlopct1Hy6ErinCaqm1/Y/e4tQqro0OpZOp+HHUsbzhNG4ziNcu1xiZp9UvQgHwtcSqlD+jXB9vepyFNGBdAohySwkFLvvBQqTrcjMZgO/JuzQpcA+hJpZmNkOwo3ubxPO9s4lXKMvq++6EqYdTLgksYxwgz2Zy6Jp1wNfY9ezyTzg1mhdnwDdCEV3CDVerouK9D+q73eRRG3rSjSBcKnkdUmbo9eze7L+6Lv/EuF6/RrgTsI9q7kxY4dq30eMddY33nWEqwtV+99GoIhQYk1mLaH66l0J63yYcIlpVSMsH0mnEu7//ThhHfcS9ref1DRfY6jtf5ZwLH2TsF1fJdQYrJp3T/fVjKJdL825lkrSO4SbpH9JdyzOuczkJYwWSuEp2O5RMfo8YF/guXTH5ZzLXH7Tu+UaBjxGqAnyMXCm+ROxzrk94JeknHPOxeKXpJxzzsXiCcM551wsnjCcc87FkjEJQ9Itkn6Q7jhaGkl9o+cPsqPPkyV9J3r/NUkvJExrkgY3VSzOudST9CVJsR4szYiEodA+/jcITTykY/0/ljRToUe1hZJ+XG38IknbEh7+eqGmZdWxnsckjZOUl9C0SEpFDxm2T/YErpk9bGbHp2rdqtZRUG2xNFeSOktaLen1asNHSZoqaWv0d1S18ZdL+kTSBkn3N+RpZkkjJL0Rvf+FUtTbX2uUeOLUhOu8RNIUSWWSHkgyfqykudE+9bJC80BV4yTpNoUOpdZKul2SEsb3j+bZGi0j8bh7Chgpad+6YsyIhEF4yvIZM9uWpvWLkLA6AeMJTW6cVW2aU6Ifu/Z78CM7mtAE877AzAZH2wwoAzoGaiS3EZoa/4ykNsCTwF8J+8yDwJPRcCSdQOjLZCyhob2BwM8bsO7RhKa1q96/34BlNBtesmQFcCNwf/URkroSnjC/ntC8/RRC3xxVLiQ0TbQf4ffjZEJT/VUeIfSP0wX4f8A/9HlHVVXjL6wzwnQ3lxvnRWi7PrFp51eAM6L3hxPaVDop+nwcMC3F8fwe+O+Ez4vYwyaMCT8sC6L33wNurza+kNAJzEpC0+U3AtnRuGzgV4RmCRYAF0ffSU6y+AituP41et+/2rST+byp7m8CryfMZ4QWXRdE6/olUZPQ0bRvEBqP+zSKb1D0v1sbTf8wUBRNv1tHQUli6Uno++BTQn8FF1TbhseAhwiN4c0CxiSMT9qBUiPvB4cSmmY5v9r3dDy7Ny+/hKj/CEKTGTcnjBsLfNKA9d8BnBe9XwG0Txh3NKFZjSsJzXesJPygnERoVO9T4NqE6bP4vEOutdF32zlh/P8SmlnZALwK7J0w7iRCB16bou3+UbL9J2EfGhy9fwC4C3iG0JbTcdTeGdHPojj+Gq1rBqHDqmuibVwKHB/zmPkmocn/XxGa/FgInBiNu4ldO8X6A+Gk8bfRejYQOmbarfOoRtqvbgQeqDbsQuDNhM/tCMfO8Ojzm8CFCeO/DbwdvR9KaBaoQ8L414j6U4k+H0Zo267W2DKlhLEP4aCv8grhgIDQXv0CPm/o7sho/G4knaOEZoiTvPrWFUhUzDuC3ZtXfji6NPGCpP3iblhUzFxP2Nl7R+/vAC6OYqrargcJjaANJnQCczxQVWS+gHBGsT+hldAz466/AU6L1nEAoS2mbyWMO5jwv+hGOOhEaBepJ6ERxj6Egx4z+zrhR7SqZHZ7knU9QvjR60nYppsljU0Y/yVCA3FFhMTyBwiNDwKXEFqg7UDol2FRso2RdHVt+0RNX0J0NvzHaD3VH2baG5hu0ZEYmR4Nrxr/QcK4D4ASSV1qWl+1dU+KYrsY+G9JGwnNaS+T9GzCpN0JHV71IrTHdA+hXbHRhH34J/q88cBLCQnlKML3vS7avirPAkMI/9v3Ccm/yn3Ad6PveiThJCGuc/i8efk3gacJ30cvQiL9QVQiq3IK4WSjE+GM+XlCsusF/IJdL1vXdsxA2F/nEboEuB24T5LM7P8RflAvifbNS6J5jyT8+BYR2pRam2yDJN1Zyz41vR7fTaJd9hkz20JI7rXtU4njFpjZphrGQygl91do3blmqciQKci45USZ1D4/I5sevX+OsBNUZdNXgNNTGMvPoy87r1p2LiA07X0N4UysqJ7L/SvhB7gd4QywbcK4EsIZQkHCsLOBl6P3/8euZwvHk7oSRmI3nd8HXkqYdkkd23gq8J+Ez9Xj+iwWQnKpYNezoluIzryibXgxYdwIYFv0vsYOlBpxP7iczztcqv49XQ88Wm36h4k6USIc6InfY2603f3rsf6hhP6uITSs+ONq448mnIFWnVF3iNZxcMI0Uwkt9UL4wRibMK4H4bhL1v1tUbSswujzEsLlj47Vptvle0nYhxJLGA8ljKurM6KfAZMSxp1CKAFU38Yi6j5mvgnMTxjXNpq3e/XjIPp8LOG4PISEjpZStG8lK2HcB9xabdgbwDej9xXs+hs5JNoeEZr1f7vavDclriNhH6yxF0OzzClhrCPsDFXeAoZKKiG0EvkQoT37roSmh1PSfr6kSwj3Mr5oZmVVw83sDTPbZmZbzewWQmutR8Rc5rLobPFswhnRKkKnLCsk/SaarB/hH7oy4cz3z4SzPYg6gElYbGIHL42t+np61jAOSd0UOi9aHp0F/5VwNhdHT+BT2/WsqHpHRNU7FMpX6Lejvh0o1Uu0rEsJ14KT2Uzo8jVRR8JllGTjq95vog7RjdH1RGeI0fsb+Lx1224Jk6+1zysQVN3/q6kToX7AvxL2rzmEH6EShSa7b5X0cfR/XBTNU/W/PINwWWqxpFckHUp89e3Yq3r8a5JsY3vqPmYgYf8xs60J8+7GzP6PUIL9I1Aq6e46z8YbV0P2qc0WskFd88Lnv6/rawsiUxLGdMIZFfDZP3cqoTntmRaaIH4TuILQ5emaZAtRqCa6uZZXjZekJH2L6EalmS2rabqqEAmZvU5m1ptwI/1FMysC7gYuNrMiM7simmwp4WypazS8yMw6mllVkXIl4Yy8SvXt2MLuHRs1VPX1rEjcnGrT3hIN29fMOhIuh6iW6ROtADpLSjxR6Eu4Fl0nM/ubmR1O+OEwkvTeBiDp2tr2iRoWfxDhDHy2Qm22O4CDFGo9ZRMuV+6bWEuFcCOy6jLmLMLNySr7AaVmlvQSR7Xt+kO0n7xCOOvtR+gRrjDaL5I1Nx7HUsI1/KKEV76ZLSdcNppAKLEVEkqCEP0vzew9M5tA+DF+gnD/A6rtd5KS7XeJ+0BVZ0SJMXQws5MauD21HTN12W3fNLPfm9lowqWcoSQ0xZ5I0p9q2afi9BSYzC77jKR2hHuEte1TieMGVjuWEsdDuGS8yMw21hZEpiSMZ0jojCfyCuH6cdX9isnVPu/GQjXR9rW8liSbT9LXgJuBcWa2oNq4vpIOk9RGoXOhHxPOvKqqO/ZXeH6hfy3bl1jD5QBCDYjEuFcCLxB66usoKUvSoIT7G48ROnjpLakTIbElmgacpdA5zJ7e4/ixpE6S+hAS9t9rmbYD4exmvaRe7H6A1dhRkJktJZwE3BJ9r/sSbuQ9nGz6RIrZgVK0nptr2ydqWMWzhB/NUdHrJ4Tr6aOis93J0fouVagifUk0X9W1/YeAbytUi+0EXEe4PFMV/wNKUq2ymv0IpYwDaJzaUX8CblJUVVNSsULHWxD+j2WEa/ZtCcdCVaxtohOxQgs9Fm7k8++6qhQ0SlI+0f2rWjRaZ0Qxjpm67LJvSjpQ0sGScgmJcDs171MX1bJP1ZiwFFqWzidUYsmO9vuq2ob/IlR9PSOa5ieEy/JV/aw8BFwhqVdUAv4h0T5lZh8SfgN+Gi3zNMIJzD8TVn8UYb+uVaYkjIeAkxS6QK3yCmFHfrWGz43pRkJ1tPcSzhT+FI3rQKjpsY5w9juecKZWdbbYh3AppbYz49HA+9EZ6XCS91f8DaANoTbKOuAfhLNcCDcznyccoO+T0MFL5HrC2cg6wj2Yv9FwTxJKd9OAfxOurdbk54QftA3RtNXjqqujoLMJP8wrCAfMT81sUowY43ag1CBmVma7dpC1ASiP3hOVeE8l/M/WEyoGnBoNx8yeI9xkfZmwbywGfpqwij5EJxzJRCXhT6OS9gF8XrV2T9xBqDjwgqRNwNuEewoQjr+qfXh2NC7R14FF0eWqiwglyaofql8Q+uf+iFArqUbW+J0R1XbM1OUO4ExJ6yT9nnAJ555oOYsJyfNXDYyrJtcRTm6uJnyH26JhWOgR8gzCvYd1hP9NYtX+PxMqDMwgVMn/N7tWADiLUFllHeHYONN27WXybGI855YxrdVKuhlYZWa/S3cs9SHpOmC1mTXZQ4dRaWYh4YZvXf1wu2ZE4VmNDwiX8crrmt65PSXpFODrZvaVOqfNlITh4vOE4ZxLhUy5JOWccy7NvIThnHMuFi9hOOeciyWjG4jr2rWr9e/fP91hOOdcRpk6deoaMyuue8pdZXTC6N+/P1OmTKl7Quecc5+R1KDWIFrNJam7Jn/MkrVb657QOedcUq0mYRjGdU/OxG/yO+dcw7SahHHBEQMp3bCdidNXpjsU55zLSK0mYeRmZ3Hz6SO5YeJsNmzzB2idc66+Wk3CABjdrzPHjSjh9ufm1j2xc865XbSqhAFw1QnDmTS7lKmL16U7FOecyyitLmEUts3lupNHcO3jMyivqEx3OM45lzFaXcIAOGXfHpQU5nPvawvTHYpzzmWMjH5wr6EkcciAztz23Fxui3E/47KxQ7h83NA6p3POuZasVSYMM+PthZ9y3qH9OPvgvrRrk0Ofzm3rntE551qxVnlJ6unpKyndsJ3rTh5B/y7teHLacj4s3eQP9TnnXC1aXcLYsLWcGyfO5ubT9yE3O4v83GwuOXYIvTsV8IuJs1m+flu6Q3TOuWap1SWM256fy7gRJYzu12mX4W3b5HDV+OFs2l7Oy3NXeWnDOeeqaVUJY+riT3lxdilXjh+edHx+bjbDu3ekfX4OD7+zxJOGc84laDU3vcsrKrn28Zlcf/IICgtya532wP6dObB/Zx55dwl5OVmctn8vJDVRpM451zy1mhLGPa8toHthPifv2yP2PGcf1Je9exYyZfE6Vm3cnsLonHOu+Ws1CUOIGyaMrHdJYVj3Dozo0ZG/vr2YDdvK/TKVc67VajUJ43tHD6Jvl4Y9a9EuL4crjh/G2s1l3PTvOazZXNbI0TnnXPOXsoQhKV/Su5I+kDRL0s+j4T+TtFzStOh1UsI810iaL2mepBNSFVtDDSxuzw/GDWXbjgqeneH9ajjnWpdU3vQuA441s82ScoHXJT0bjfutmf0qcWJJI4CzgL2BnsCLkoaaWUUKY6y39nk5tM/LYcGaLbw8dxWHDe5Km5xWU1BzzrViKfuls2Bz9DE3etV2A2AC8KiZlZnZQmA+cFCq4ttTRw0t5pjh3fjjy/OZNLs03eE451zKpfTUWFK2pGnAKmCSmb0TjbpE0nRJ90uqeoKuF7A0YfZl0bDqy7xQ0hRJU1avXp3K8GO5fNxQ9ulVyMTpK1i3Zcdu4++a/DFL1m5NQ2TOOde4UpowzKzCzEYBvYGDJI0E7gIGAaOAlcCvo8mTVV/arURiZneb2RgzG1NcXJySuOure2E+Rwwu5rEpS9levusVNMO47smZXrvKOZfxmuTiu5mtByYD482sNEoklcA9fH7ZaRnQJ2G23sCKpoivMRS2zeW7Rw1i2tL13PrsXLbtCInjgiMGUrphOxOn+01y51xmS2UtqWJJRdH7AuA4YK6kxCfnTgNmRu+fAs6SlCdpADAEeDdV8aXKIQO7cNFRA1m6bisvz1tFbnYWN58+khsmzmbDtvJ0h+eccw2WyhJGD+BlSdOB9wj3MCYCt0uaEQ0/BrgcwMxmAY8Bs4HngIubWw2puIratmFoSQfKyiv5sHQTQ0s6cNyIEm6P0VmTc841V8rka+tjxoyxKVOmpDuMWm3dsZPfTvqQwwZ35cp/TOeuc0fv1lKuc841JUlTzWxMfefzBwhSrG2bHP7fF0cwul8nDh/clav/OZ3yisp0h+Wcc/XmCaOJdMjP5eoTh5OdJX734ofpDsc55+qt1TRv3hx065jPdV8cwbn3vcMfX/64zukvGzuEy8cNbYLInHOubp4wmpCZcfdrC7hq/HAuOGIAUxev4+PVWzjrwD5kZXl/G8655s0TRhN6evpKSjds5ztHDCAnO4uDB3ahXV4OH63azM7KSvbuWZjuEJ1zrkZ+D6OJbNhazo0TZ3Pz6fuQm/351z6yVyGDu7XnP0vW89bHa6mozNxaa865ls1LGE3ktufnMm5ESdIqtdlZ4txD+mFm/PL5eYzp34ljh5ekIUrnnKuZlzCawNTFn/Li7FKuHD+81ukkceX44YzoUciT05azfP22JorQOefq5gkjxcorKrn28Zlcf/IICgtyY83TvTCf40d0Z9KsT1i1aTs7/bkN51wz4Akjxe55bQHdC/M5ed8edU+coKBNNt88bACbtu/kxn/PYdXG7SmK0Dnn4vGEkWJC3DBhJFLDqs0OKm7PT04eAYJfvzDPGzB0zqWNtyWVQdZuLmPlhu2s3bKDI4d0bXAScs61bt6WVCvQpX0eI3sVkpstHnl3qd/bcM41KU8YGegLg7pyzsF9eXLaCu55dUHS3vy8a1jnXGPzhJHBzhjdmwmjejJ18TreXrB2l3HeNaxzrrF5wshw3TrmM7pfJ0o3bmfRmi18umUH4F3DOucanyeMFkASE0b1omuHPP7yxkLeXrDWu4Z1zjU6TxgtSPu8HH54/DAO7N+Z37wwj8KCXO8a1jnXaFKWMCTlS3pX0geSZkn6eTS8s6RJkj6K/nZKmOcaSfMlzZN0Qqpia+mys8R/jR3C2s07OHpoMS/MLmXq4nXpDss5l+FSWcIoA441s/2AUcB4SYcAVwMvmdkQ4KXoM5JGAGcBewPjgTslZacwvhYtN2o+fUTPjhzQt4hrHveuYZ1zeyZlCcOCzdHH3OhlwATgwWj4g8Cp0fsJwKNmVmZmC4H5wEGpiq+16N2pLX86dzTt8nI45563+WSDNzHinGuYlDZvHpUQpgKDgT+a2TuSSsxsJYCZrZTULZq8F/B2wuzLomHVl3khcCFA3759Uxl+iyGJH44bxrn3vcMht7xU5/TeNaxzLpmUJgwzqwBGSSoC/iVpZC2TJ2vnYreHCMzsbuBuCE2DNEacLV1i17DfO3oQZsbqzWW8t3Ada7eU8dUD+5CX41f/nHO1a5JaUma2HphMuDdRKqkHQPR3VTTZMqBPwmy9gRVNEV9Ll9g1LIQSR7cO+Xxx3x4cNbSY1ZvKuOPFj1i5wfvfcM7VLJW1pIqjkgWSCoDjgLnAU8B50WTnAU9G758CzpKUJ2kAMAR4N1XxtRY1dQ1bpV+XdvTu1JZvHtafjdt28vA7i5m+bH3TB+qca/ZSeUmqB/BgdB8jC3jMzCZKegt4TNK3gSXAlwHMbJakx4DZwE7g4uiSltsDtXUNm6iwIJfCglwGFrdj8rzVzFqxgaWfbmPciBKys7xVXOdcChOGmU0H9k8yfC0wtoZ5bgJuSlVMrU1V17CTrjgq9jy52VmMG1GCmbF1RwXzPtnE7JUbOXFkd9rlxdtd7pr8MV/cpwd9u7RtaOjOuWbIn/RuoRrSNWwiSRzYvzMjenbkoP6deWbGSqYtXc+KGP2Me8OHzrVMnjBaqIZ2DZtM3y5t+fKYPgzo2o5nZ37C4rVb+GDp+hqn94YPnWuZPGG0UHvaNWwyhQW5fPvwAfQqKmDJp1t5+oMVLFyzhYrKXUsS3vChcy2TJ4wW6ntHD0rZPYSc7CxO2a8np+zXk/Vbd/DL5+exbssONpft/Gya0f06e8OHzrUwsROGpHapDMRlpv37duLqE4dTacY9ry7gzY/XsCVKHFedMJxJ3vChcy1GnQlD0hckzQbmRJ/3k3RnyiNzGaVL+zwuHzeULwzqypPTVvDfL31E27xsrjt5BNc+PsMbPnSuBYhTwvgtcAKwFsDMPgCOTGVQLrOdc3BfLjl2MHNXbuKDJevokJ/Dva8tTHdYzrk9FKtivZktrXbz1B+oc7WSxD69CxnZqyOT563i/AemcFuM+xne8KFzzVechLFU0hcAk9QGuJTo8pRzcfzlzcX8+IShjN2rhEVrtnDU0G4UtPHGDp3LNHEuSV0EXExoanwZoTOki1MYk2tBqho+vPDIQQzv3pH9+hTx4FuL2LpjZ90zO+ealTgJo8DMvmZmJWbWzczOJXSG5FytkjV82KOwgIuOGsSMZRu46d+zd6mK65xr3uIkjIWSHolanK3yTKoCci1HbQ0fHjywC5ePG8qaTWU8/M7i3R7+c841P3ESxgzgNeB1SYOiYd58qatVVcOHV44fXuM0bdvk0L9rO0b1KWL2io18VLqpCSN0ztVXnIRhZnYn4Wb305JOIUlPeM5VqW/Dh3v3LGTvnh15b9E6npy2vAkidM41RJxaUgIwszckjQX+DtR82uhavYY0fJiVJc45uC9mxv2vL6R3pwKO37t7CqN0ztVXnBLGSVVvzGwlcCyhq1XnktqThg8l8a3DB9CzqIC3Pl4bqzl151zTqLGEIelcM/srcHYNB/6rKYvKZbTvHT2o7onqMLJXIdt2VHD/Gws568A+FBbkkpOki1nnXNOp7QisamywQw0v51KqoE02Fx8zmB0VldwwcTbzV21Od0jOtWpKVa9okvoADwHdgUrgbjO7Q9LPgAuA1dGk15rZM9E81wDfJjQ9cqmZPV/bOsaMGWNTpkxJSfyueamsNDaV7eR/3lrEeV/oT4f83W+me9ewzsUjaaqZjanvfHFaq71dUkdJuZJekrRG0rkxlr0T+KGZ7QUcAlwsaUQ07rdmNip6VSWLEcBZwN6EeyR3SvL2IxwQbooXFuTy1QP78t6iT/mwdNNuXcB617DOpVaci8LHm9lG4GRC0yBDgR/XNZOZrTSz96P3mwjtT/WqZZYJwKNmVmZmC4H5wEEx4nOtSHGHPI4dXsKqjWXc8dJHu4zzrmGdS604CaOq7H8S8IiZfVrflUjqD+wPvBMNukTSdEn3S6p6DLgXsDRhtmXUnmBcK3b4kK5cNnYIL89bxf2vL6S8otK7hnUuxeIkjKclzQXGAC9JKga2x12BpPbAP4EfRCWVu4BBhEYMVwK/rpo0yey7XVuQdKGkKZKmrF69OsksrrWQxDHDunHk0GIWrN7C1MXrvGtY51KozoRhZlcDhwJjzKwc2Eq4fFQnSbmEZPGwmT0eLa/UzCrMrBK4h88vOy0D+iTM3htYkSSeu81sjJmNKS4ujhOGa+EGd2vP0JL2LFyzhenL1nPpsUO8a1jnUiBWxXYzW2dmFdH7LWb2SV3zKDy8cR8wx8x+kzA88fHf04CZ0fungLMk5UkaAAwB3o23Ga61k8SZo3szokdHHnxrEafu38u7hnWukcXqca+BDgO+DsyQNC0adi3hQcBRhMtNi4DvApjZLEmPAbMJNawurkpSzsWVk53FVeOHs3HbDv5v7ip+9fw8rjlpr3SH5VyLUOtzGFEpobeZLa1xojTy5zBcbV6eW8r5D8TbP7xrWNeaNPQ5jFpLGGZmkp4ARjc0MOfSwcz4y5uLuWr8MMp2VnLpsUPIyvJW+Z3bE3HuYbwt6cCUR+JcI6rqGvY7Rwxk3IgS3ltU79rgzrlq4iSMYwhJ4+Po2YkZkqanOjDnGqp617B79yykqG0bXv9oTbpDcy6jxbnpfWLKo3CuESXrGnZY9w68OKeUQwZ29lZvnWugOM9hLCY8H3Fs9H5rnPmcS4fauoa9+JjBTJy+kh07vaqtcw0Rp/HBnwJXAddEg3KBv6YyKOcaIk7XsAf07cSj7y1p4sicaxnilBROA74EbAEwsxV4fxiuGYrTNWzfLm35ypg+TJ63qgkjc65liJMwdlh4WMMAJLWrY3rn0iJu17D5udnMWrGRpZ9ubaLInGsZ4iSMxyT9GSiSdAHwInBvasNyrv6+d/Sg2J0nXXTUINZu2cGWsp0pjsq5liPOTe9fAf8gNCI4DPiJmf0+1YE5l0rZWaJXUQG/r9anhnOuZnFuet9mZpPM7Mdm9iMzmyTptqYIzrlUKu6QxyXHDmbm8g3pDsW5jBDnktS4JMP82QzXInTIz2Xq4nW8v8SbQneuLjUmDEnfkzQDGBY94V31Wgj4k96uxfjGof3omJ/Dui070h2Kc81abSWMvwGnEPqpOCXhNdrMzm2C2JxrEpIo6ZjP7178kIrKmltvdq61qzFhmNkGM1sEXAd8Ej3lPQA4V1JR04TnXNPokJ/Ltw4fwJrNZekOxblmK849jH8CFZIGE3rQG0AofTjXovTr0o63F6z1h/qcq0GchFFpZjuB04HfmdnlQM2P0jqXwSaM6kW7vBzKdnpnj85VFydhlEs6G/gGMDEalryhHudagNF9O3Hzv+d40nCumjgJ43zgUOAmM1soaQAxGh+U1EfSy5LmSJol6bJoeGdJkyR9FP3tlDDPNZLmS5on6YSGbpRzeyIrS1xw5EAWrN6S7lCca1biPOk928wuNbNHos8LzezWGMveCfzQzPYCDgEuljQCuBp4ycyGAC9Fn4nGnQXsDYwH7pSU3ZCNcm5P9e7Ulm3lFTw5bXm6Q3Gu2YjzpPdCSQuqv+qaz8xWmtn70ftNwBygFzABeDCa7EHg1Oj9BOBRMyszs4XAfOCgem+Rc43kgL6dMAv9gzvn4vW4NybhfT7wZaBzfVYiqT+wP/AOUGJmKyEkFUndosl6AW8nzLYsGuZc2py6fy/+/MrHnHNwXzrk+60717rFuSS1NuG13Mx+BxwbdwWS2hOq5v7AzDbWNmmy1SdZ3oWSpkiasnr16rhhONdgZ47uzeR5vq85F+eS1AEJrzGSLiJmB0qScgnJ4mEzezwaXCqpRzS+B1BV6X0ZoSvYKr2BFdWXaWZ3m9kYMxtTXFwcJwzn9kiX9nkcNayYR9/1nvpc6xbnktSvE97vBBYBX6lrJoVebO4D5pjZbxJGPQWcB9wa/X0yYfjfJP0G6AkMAd6NEZ9zKdcxPxcDVm3aTrcO+ekOx7m0qDNhmNkxDVz2YcDXgRmSpkXDriUkisckfRtYQrgngpnNkvQYMJuQmC42M68I75qNsw/qy0tzStm3tyjukJfucJxrcqqpBoikK2qbsVqpIS3GjBljU6ZMSXcYrhXZUraTP7+6gCvGDU13KM41mKSpZjam7il3Vds9jA51vJxrddrl5XD5cUN4bubKdIfiXJOr8ZKUmf28KQNxLlNIorzCeGP+Gg4b3DXd4TjXZOLUknowsTlzSZ0k3Z/SqJxr5k7Zrye9igpYsnZrukNxrsnEaUtqXzNbX/XBzNYRHsJzrlXrWVTAva8vYMfOynSH4lyTiJMwsqo1ENiZeNVxnWvR2uRkcdX44XxYuom7Jn/spQ3X4sV9DuNNSf8gPHn9FeCmlEblXIZol5fDuq07mLF8PW8tWMuD5x9IeATJuZYnTtMgDwFnAKXAauB0M/ufVAfmXKY4YkgxV40fzrJ1W5k43WtPuZYr1qUlM5tNeKDOOZdEn05tGV7SgV88PZsjhxZTWOANFbqWJ849DOdcHbKyxM2n78PBAztz+3Nz0x2OcynhCcO5RlLUtg1fO6gfT05bwdTF69IdjnONLlbCkNRP0nHR+wJJ/qS3c0kcOrgLFxwxgKv/OZ3yCq9u61qWOA/uXQD8A/hzNKg38EQKY3Iuo106dghlOyu448WP0h2Kc40qzk3viwldpb4DYGYfJfSS55yrRhLXn7w3Fzw0hT+8PL/O6S8bO4TLvTFDlwHiJIwyM9tRVbdcUg5JesJzzgVmxv+8vZjvHjmQvl3a8rWD+6U7JOcaRZx7GK9IuhYokDQO+F/g6dSG5Vzmenr6Sko3bOdHJwwjNzvL72W4FiNOwria8MDeDOC7wDPAdakMyrlMtWFrOTdOnM3Np+9DbnYWXxnTh0ffW0rZTu8LzGW+OD3uVQL3RC/nXC1ue34u40aUMLrfZ82vcfTQYp7+YCVnju6dxsic23N1JgxJM9j9nsUGYApwo5mtTUVgzmWaqYs/5cXZpUy64qhdhvfp3JaCNtm8OX8NX/D+M1wGi3NJ6lng38DXotfTwKvAJ8ADNc0k6X5JqyTNTBj2M0nLJU2LXicljLtG0nxJ8ySd0MDtcS4tyisqufbxmVx/8oikzYJ0bZ/HKx+tZuuOnWmIzrnGEaeW1GFmdljC5xmS3jCzwySdW8t8DwB/AB6qNvy3ZvarxAGSRgBnAXsDPYEXJQ01M7/w6zLCPa8toHthPifv26PGaa48YThTF6/jwP6dvEVbl5HilDDaSzq46oOkg4D20ccaT5fM7FXg05hxTAAeNbMyM1sIzCc8++FcRhDihgkja00E2VkiS/DYlKVNGJlzjSdOCeM7wP2S2gMCNgLfkdQOuKUB67xE0jcI90B+GPXg1wt4O2GaZdEw5zLC944eFGu6Mf0707tTW1as30bPooIUR+Vc44rTH8Z7ZrYPMAoYZWb7mtm7ZrbFzB6r5/ruAgZFy1pJ6JwJQiLabdXJFiDpQklTJE1ZvXp1PVfvXPp1bd+GP73ysXft6jJOrP4wJH2RcH8hv6rIbWa/qO/KzKw0YZn3ABOjj8uAPgmT9gZW1LCMu4G7AcaMGeNPnLuMk5OdxUVHDWLtljJ6FHopw2WOOI0P/gn4KvBfhJLAl4EGtXUgKfGO4GlAVQ2qp4CzJOVJGgAMAd5tyDqcywQ9iwp4Z8GnvDl/TbpDcS62ODe9v2Bm3wDWmdnPgUPZtTSQlKRHgLeAYZKWSfo2cLukGZKmA8cAlwOY2SzgMUKvfs8BF3sNKdfSTRjVkwozzLyg7DJDnEtS26O/WyX1BNYCA+qayczOTjL4vlqmvwm4KUY8zrUIkjhkYBdufXYuV5843KvaumYvTgnjaUlFwC+B94FFwCMpjMm5ViM3O4uxe5WwbN22dIfiXJ1qLWFIygJeMrP1wD8lTQTyzWxDUwTnXGtw0IDOPDtjJTsrjQFd26U7HOdqVGsJI2p48NcJn8s8WTjX+I4bUcJLc0rrntC5NIpzSeoFSWfIL7A6lzK52Vl8+/ABPP7+snSH4lyN4tz0vgJoB1RI2kaoWmtm1jGlkTnXykiiotKYtWIDe/csTHc4zu0mzpPeHcwsy8xyzaxj9NmThXMpcObo3pRXGBu2lac7FOd2E+fBPUk6V9L10ec+UQOEzrlGJon+Xdpy58vz0x2Kc7uJcw/jTsLDeudEnzcDf0xZRM61ckVt2/BfY4cw75NN6Q7FuV3ESRgHm9nFRA/wRa3LtklpVM61cu3zcnhh1icsWbs13aE495k4CaNcUjZR67GSigFvZtO5FPvuUYPYsmMnOyv8cHPNQ5yE8XvgX0A3STcBrwM3pzQq5xxtcrLIzc7i3tcXpjsU54AY1WrN7GFJU4GxhCq1p5rZnJRH5pxjcLf2rN1cxs6KSnKy45zfOZc6cWpJ3QF0NrM/mtkfPFk417QOHtiF30z6kE3bvaqtS684pyzvA9dJmi/pl5LGpDoo59yuvnPEQGYs81Z5XHrFeXDvQTM7CTgI+BC4TdJHKY/MOfeZzu3aUFKYz9MfJO2I0rkmUZ+LooOB4UB/YG5KonHO1WhQcXs+3bKDikrvcMmlR5x7GFUlil8As4DRZnZKyiNzzu3mvC/0559Tl3nScGkRp/HBhcChZuadDzvXDIzqW8QT/1nOGaN7pzsU18rEqVb7J0mdovaj8hOGv5rSyJxzSQ0t6UBRQS7Tl61n395F6Q7HtSJxLkl9B3gVeB74efT3ZzHmu1/SKkkzE4Z1ljRJ0kfR304J466JamLNk3RCQzbGudaiuEMej7+/nP9+6SNvPsQ1mTg3vS8DDgQWm9kxwP7A6hjzPQCMrzbsakKXr0OAl6LPSBoBnAXsHc1zZ9QciXMuCUlcfeJwVm8u47onZ2Lm9zRc6sVJGNvNbDuApDwzmwsMq2um6JLVp9UGTwAejN4/CJyaMPzRqAvYhcB8QjVe51wN8nOzOWJIV+aXbmLi9JXpDse1AnFuei+TVAQ8AUyStA5oaGXwEjNbCWBmKyV1i4b3At5OXGc0bDeSLgQuBOjbt28Dw3CuZRg3ojs7dlby06dmceTQYgoLctMdkmvB4jy4d5qZrTeznwHXA/fxecmgsSTrLzxpGdvM7jazMWY2pri4uJHDcC7zjN2rhM7t2nDbs/54lEuterVmZmavmNlTZrajgesrldQDIPq7Khq+DOiTMF1vGl6Kca5Vyc/N5rYz9uWF2Z8wdfG6dIfjWrCmbv7yKeC86P15wJMJw8+SlCdpADAEeLeJY3MuY+3ftxNj9yrhir9Po9z7z3ApkrKEIekR4C1gmKRlkr4N3AqMi54cHxd9xsxmAY8Bs4HngIvNrCJVsTnXEt186kja5+dw72vef4ZLDWVydbwxY8bYlClT0h2Gc83Gy3NLOf+BeMfEZWOHcPm4oSmOyDVHkqaaWb1bHo9TS8o5lwHMjL+8uZiLjhrICXt3Z2SvQnK90yXXiHxvcq6FeHr6Sko3bOeHxw+jX5d23PzMHHbs9PsZrvF4wnCuBdiwtZwbJ87m5tP3ITc7i87t2vCtwwawcM0Wynb67UDXODxhONcC3Pb8XMaNKGF0v8+aZ6NP57Z0yM/hlmfmetJwjcIThnMZburiT3lxdilXjh++27ieRQV896iBLFyzhe3lnjTcnvGE4VwGK6+o5NrHZ3L9ySNqbBakR2EBHfNzucXvabg95AnDuQx2z2sL6F6Yz8n79qh1up5FBVw+bigzlm/wkoZrME8YzmUwIW6YMBIpWXNsuypq24aeRfn88vl53hy6axBPGM5lsO8dPYi+XdrGnr5HYQFXnzicF2aXeknD1ZsnDOdamdzsLPbrXcQ9ry5Idyguw3jCcK4V6l6Yz3+NHcLf3lnCth1e0nDxeMJwrhU7bq9uPDFtebrDcBnCE4ZzrVi3jvmcfVBf7njxIy9puDp5wnDOcfbBfXjto9Vee8rVyhOGc45uHfI5dng3fv70bLbu2JnucFwz5QnDOQdATnYWFx8zmDkrN3mVW5eUJwzn3GeKO+QxpKQ9Nz8zhy1lXtJwu/KE4ZzbRcf8XC4bO4TFa7d60nC7SEvCkLRI0gxJ0yRNiYZ1ljRJ0kfR3051Lcc5lxpd2udR0jGPW5+dy2ZPGi6SzhLGMWY2KqFf2auBl8xsCPBS9Nk5lyZd2udx+bihLF67xZOGA5rXJakJwIPR+weBU9MXinMOoHO7NvQoLOC2Z+f6jXCXtoRhwAuSpkq6MBpWYmYrAaK/3ZLNKOlCSVMkTVm9enUThetc69W5XRuuPnE4M5ZvYNP28nSH49IoXQnjMDM7ADgRuFjSkXFnNLO7zWyMmY0pLi5OXYTOuc+0y8thcHF7fjvpo88e7rtr8scsWbs1zZG5ppSWhGFmK6K/q4B/AQcBpZJ6AER/V6UjNudccp3ateG6L+7FczM/YeP2cgzjuidn+tPhrUiTJwxJ7SR1qHoPHA/MBJ4CzosmOw94sqljc87VLitLHDqoCw+9uYgLjhhI6YbtTJy+Mt1huSaSjhJGCfC6pA+Ad4F/m9lzwK3AOEkfAeOiz865ZqaobRsuOXYIj767hGu/uBc3TJzNhm1+b6M1UCYXJ8eMGWNTpkxJdxjOtUobtpbz7MyVTF++AQE3nbZPukNyMUmamvBIQ2zNqVqtcy6DFLbN5ayD+tIhL4cXZpUydfG6dIfkUswThnNuj3z/mMGctn8vrnl8OuUVlekOx6WQJwzn3B4pLMjlyvHD2FpWwQ8f+4ANW8upqEzdpW6vzps+OekOwDmX+XKys7j1jH059753eOqDFXVOf9nYIVw+bmiD1lVVnffB8w9EUoOW4RrGE4Zzbo+ZGXe/toCrxg/ne0cP+mz4ms1lfFi6iXZtcpg0u5Tj9y6hS/s8ehbmN3hdFxwxkCf/s4KJ01dyyn49GyN8F5MnDOfcHnt6+kpKN2znO0cM2GV41/Z5dG2fB8B+fYqorDQmzSnlf6cs5bT9e/Hqh6s5elg3+nRuG3tdudlZ3Hz6SL731/c5cmgxhQW5jbotrmZ+D8M5t0c2bC3nxomzufn0fcjNrv0nJStLnLB3d35w3FD6dWnHcSNKKNtZyRP/Wc6vnp/H2s1lLFi9uc6nx0f368xxI0q4/bm5jbkprg5ewnDO7ZHbnp/LuBEljO5X/y5sehQWADC4W3sAtpdX8PK81cxYvoGO+bms2rSd8SN7JC1FXHXCcMb99hVOP6B3g9ZdX3dN/pgv7tODvl3il4ZaGi9hOOcabOriT3lxdilXjh/eKMvLz83mzNG9mTCqF0cPK+bQgV0xM34z6UMefHMRy9dvY9m6UEOqsG0u1508gmsfn9Ek1XnT2XZWc6kZ5gnDOdcg5RWVXPv4TK4/eURK7iNIom+XthS1bcMV44Zy3hf6k5slXv1wDTOWbeD+1xeys6KCko553PvawkZff3XpbDuruTT06JeknHMNcs9rC+hemM/J+/ZosnV265jPOQf3BWBkr44s+XQrlQY/+t/p3BbjfsaeVOdN58325lIzzBOGc65BhLhhwsi0PQshib6d23L9k7N2qc67YVs5bbKzeG7WSsp3GkNK2vPy3FUcu1cJqzeV8dbHazmwfyeys1Tv2BNvtjdl21nNpWaYNz7onMtYT32wgj/+33wmXnp4nTW0AEo3bufjVZsZXNKev761mKHdO9A+L4cZyzZw+ujebNxWTu9OBXTIr/kHecPWcsb99hXuOnd0k9xsT3Ttv2Y0SkOPDW180BOGcy4jNeYP946dlRjGax+uoXTTdg4e0JmnPljJwQM6s728AjM4bHBXcrJFbnZWvRNVY2msbW5owvBLUs65jLQn1Xmra5MTfvSPG1Hy2bArxnUAYOP2chau3kLZzgrueW0xJR3zGNi1HdvKK7jlmbl8YVAXynZWcsTQrkxfuoHuhXl0zM9l4/ZyehQWsGNnJbk5WbRrk73Hl+8Sa4Y1dbICTxjOuQxUVZ130hVHpXxdHfNz2a9PEQCXjh3y2fCbThvJ1+97l/vfqLuG1pFDunLLGfvy9/eWsnfPjuRmi2lLN3DGAb3494yVZCk80PjEf5Yzqm8RZeUVzFm5ibMP6suj7y2hqCCXI4YW8/pHaxjVp5C2ednc+9rCXZphaQp+Sco5l1HKKyo5+fevc8mxg9NWY8jMOO8v73HowC5N8qNdUWlUVBqbtpfTJieLVz5czSV/+0+seZPVDPNLUs65ViEd1Xmrq6ntrFTJzhLZWaJL+zzMjMemLNutocem0Owe3JM0XtI8SfMlXZ3ueJxzzUu6q/PWp+2sVGjqZJWoWZUwJGUDfwTGAcuA9yQ9ZWaz0xuZc665aOqz6uoa82Z7fVUlq7vOHZ2WZNXcShgHAfPNbIGZ7QAeBSakOSbnnAMav+2s+kpnsoLmlzB6AUsTPi+Lhn1G0oWSpkiasnr16iYNzjnXeqW67ay6pDtZQfNLGMkuSu5SjcvM7jazMWY2pri4uInCcs61dum82Z7uZFWluSWMZUCfhM+9gbo7CHbOuRRL58325lAzDJrZTW/gPWCIpAHAcuAs4Jz0huScc+m92Z7ummFVmlXCMLOdki4BngeygfvNbFaaw3LOubRKd82wKhn9pLek1cDidMfRAF2BNekOoon5NrcOrW2bM3V7+5lZvW8CZ3TCyFSSpjTksfxM5tvcOrS2bW5t29vcbno755xrpjxhOOeci8UTRnrcne4A0sC3uXVobdvcqrbX72E455yLxUsYzjnnYvGE4ZxzLhZPGCkUp28PSUdLmiZplqRXmjrGxlTX9koqlPS0pA+i7T0/HXE2Jkn3S1olaWYN4yXp99F3Ml3SAU0dY2OLsc1fi7Z1uqQ3Je3X1DE2trq2OWG6AyVVSDqzqWJrSp4wUiShb48TgRHA2ZJGVJumCLgT+JKZ7Q18uanjbCxxthe4GJhtZvsBRwO/ltSmSQNtfA8A42sZfyIwJHpdCNzVBDGl2gPUvs0LgaPMbF/gBlrGjeEHqH2bq46B2wgtVbRInjBSJ07fHucAj5vZEgAzW9XEMTamONtrQAeFBnHaA58CO5s2zMZlZq8StqMmE4CHLHgbKJKU3hbk9lBd22xmb5rZuujj24RGRDNajP8zwH8B/wQy+TiulSeM1Kmzbw9gKNBJ0mRJUyV9o8mia3xxtvcPwF6EFohnAJeZWWXThJc2cb6XluzbwLPpDiLVJPUCTgP+lO5YUqlZNT7YwtTZtwfh+x8NjAUKgLckvW1mH6Y6uBSIs70nANOAY4FBwCRJr5nZxhTHlk5xvpcWSdIxhIRxeLpjaQK/A64ys4p0tyibSp4wUidO3x7LgDVmtgXYIulVYD8gExNGnO09H7jVwsM/8yUtBIYD7zZNiGnRKvt4kbQvcC9wopmtTXc8TWAM8GiULLoCJ0naaWZPpDWqRuaXpFLns749ohu7ZwFPVZvmSeAISTmS2gIHA3OaOM7GEmd7lxBKU0gqAYYBC5o0yqb3FPCNqLbUIcAGM1uZ7qBSSVJf4HHg6xlaWq43MxtgZv3NrD/wD+D7LS1ZgJcwUqamvj0kXRSN/5OZzZH0HDAdqATuNbNaq+01V3G2l1Bj5gFJMwiXaq4ys0xsGvozkh4h1PjqKmkZ8FMgFz7b5meAk4D5wFZCKSujxdjmnwBdgDujM+6dmd6ia4xtbhW8aRDnnHOx+CUp55xzsXjCcM45F4snDOecc7F4wnDOOReLJwznnHOxeMJwzjkXiycM55xzsXjCcC2CpCJJ3093HHVJdZySLoj6V5kmqTLh/W9StU7XeviDe65FkNQfmGhmI5tBLCIcW7u1xNuQOGtbXi3z9ALeNLN+cedxri5ewnAtxa3AoOhs+peSzpX0bvT5z1HnNkjqL2mupHslzZT0sKTjJL0h6SNJByVM82DUa9w/ora+iJax27KjeeZIuhN4H+gj6Ymo2fpZki6sIc7+ib24SfqRpJ/VsLyk21SDkYQm5J1rNJ4wXEtxNfCxmY0C7ge+ChwWfa4AvpYw7WDgDmBfQmu55xCa4P4RcG00zTDg7qjXuI3A9wEk7VXLsocROkva38wWA98ys9GElkwvldQlMU4z+3Ed2/TZ8oC2dWxTdfsAGdkumWu+vPFB1xKNJfQz8l7U+F0Bu/aCttDMZgBImgW8ZGYWNYrYP5pmqZm9Eb3/K3Ap8Ktalv0qsDjqVa/KpZJOi973IXTT+kk9tiNxeXVtU3UjgUn1WJdzdfKE4VoiAQ+a2TU1jC9LeF+Z8LmSz4+J6jf3qj4nXXZ0b2JLwuejgeOAQ81sq6TJQH6SWHaya0k/cZotCe/r2qbq9gF+G3Na52LxS1KupdgEdIjevwScKakbgKTOkup787evpEOj92cDr9dz2YXAuihZDAcOSRInQCnQTVIXSXnAyTXEE3ubJGURSjNz42yoc3F5wnAtQtSr2xvRDeTzgeuAFyRNJ1ya6VHPRc4Bzovm7wzcFa1ndsxlPwfkRNPcALxdPU5JvzSzcuAXwDvARGr4ka/HeiHco1lmZmU1jHeuQbxarXPVNKcqus41J17CcM45F4uXMJxzzsXiJQznnHOxeMJwzjkXiycM55xzsXjCcM45F4snDOecc7F4wnDOOReLJwznnHOx/H8Dct9Lorn/yQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Plotting\n", "with h5py.File(\"xy_data.hdf5\", \"r\") as f:\n", " plt.figure()\n", " # My girlfriend told me to use this marker.\n", " plt.plot(temperatures, np.array(f[\"cluster-size\"]), marker=\"4\", markersize=20, linestyle=\"dashed\", linewidth=.5)\n", " plt.xlabel(\"temperature $T$\")\n", " plt.ylabel(\"average cluster size\")\n", " plt.title(\"Average cluster size for the $w \\\\times w$ XY model \\n\"\n", " \"(w = {}, #equilibrations = {}, #measurements = {})\"\n", " .format(width, equilibration_moves, measurement_moves))\n", " plt.show()" ] }, { "cell_type": "markdown", "id": "bd5fa10f", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "markdown", "checksum": "bcf91e9b587d185ef00c6b25028e18e5", "grade": false, "grade_id": "cell-101fca705096fe8e", "locked": true, "schema_version": 3, "solution": false, "task": false } }, "source": [ "__(d)__ Make an MCMC estimate (and plot!) of the **mean square magnetization per spin** $\\mathbb{E}[ m^2(s) ]$ for the same set of temperatures, where \n", "\n", "$$\n", "m^2(s) = \\left(\\frac{1}{N}\\sum_{i=1}^N s_i\\right)\\cdot\\left(\\frac{1}{N}\\sum_{i=1}^N s_i\\right).\n", "$$\n", "\n", "To choose the equilibration time and time between measurement, use the average cluster size from (c) to estimate how many moves correspond to 1 _sweep_, i.e. roughly $N = w^2$ updates to spins. Then use 100 equilibration _sweeps_ and 200 measurements of $m^2(s)$, with 2 _sweeps_ between each measurement. Store the measured values of $m^2(s)$ in the data set `\"square-magn\"` of dimension $(11,200)$ in `xy_data.hdf5`.\n", "Then read the data and plot estimates for $\\mathbb{E}[ m^2(s) ]$ including errors (based on batching or jackknife). If the errors are too small to see, you may multiply them by some number and indicate this in the title of the plot. **(20 pts)**" ] }, { "cell_type": "code", "execution_count": 9, "id": "71323e81", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "f347e7079633e91d931730dfbcc231d4", "grade": false, "grade_id": "cell-7421b7a140318565", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "T = 0.5\n", "\tequilibrating took 2.6317245960235596 seconds\n", "\tmeasuring took 13.053837060928345 seconds\n", "\n", "T = 0.6\n", "\tequilibrating took 2.8574910163879395 seconds\n", "\tmeasuring took 12.967190027236938 seconds\n", "\n", "T = 0.7\n", "\tequilibrating took 3.4188122749328613 seconds\n", "\tmeasuring took 16.535584449768066 seconds\n", "\n", "T = 0.8\n", "\tequilibrating took 3.129852294921875 seconds\n", "\tmeasuring took 14.666396617889404 seconds\n", "\n", "T = 0.9\n", "\tequilibrating took 2.749675750732422 seconds\n", "\tmeasuring took 12.164067506790161 seconds\n", "\n", "T = 1.0\n", "\tequilibrating took 2.770905017852783 seconds\n", "\tmeasuring took 11.673817873001099 seconds\n", "\n", "T = 1.1\n", "\tequilibrating took 2.59846830368042 seconds\n", "\tmeasuring took 9.312909126281738 seconds\n", "\n", "T = 1.2000000000000002\n", "\tequilibrating took 2.9586286544799805 seconds\n", "\tmeasuring took200] 10.83272385597229 seconds\n", "\n", "T = 1.3\n", "\tequilibrating took 2.4720122814178467 seconds\n", "\tmeasuring took 9.212384223937988 seconds\n", "\n", "T = 1.4\n", "\tequilibrating took 3.3535523414611816 seconds\n", "\tmeasuring took 10.767858028411865 seconds\n", "\n", "T = 1.5\n", "\tequilibrating took 3.0648155212402344 seconds\n", "\tmeasuring took 9.897510051727295 seconds\n", "\n" ] } ], "source": [ "measurements = 200\n", "equil_sweeps = 100\n", "measure_sweeps = 2\n", "\n", "with h5py.File(\"xy_data.hdf5\", \"a\") as f:\n", " if not \"square-magn\" in f:\n", " state = xy_aligned_init_config(width)\n", " square_magn = np.zeros((len(temperatures), measurements))\n", " \n", " for idx, T in enumerate(temperatures):\n", " toc()\n", " print(\"T =\", T)\n", " beta = 1/T\n", " # We let the sweep size depend on the width and the average cluster size\n", " # from before:\n", " sweep_size = np.ceil(width**2/f[\"cluster-size\"][idx]).astype(int)\n", " \n", " # Equilibrate\n", " for j in range(equil_sweeps*sweep_size):\n", " xy_cluster_move(state,beta)\n", " print(\"\\tequilibrating... [{}/{}]\".format(j, equil_sweeps*sweep_size), end='\\r')\n", " print(\"\\tequilibrating took\", toc(), \" seconds\")\n", " \n", " # Measure\n", " for j in range(measurements):\n", " xy_cluster_move(state,beta)\n", " square_magn[idx][j] = np.linalg.norm(np.mean(state, axis=(0,1)))**2\n", " print(\"\\tmeasuring... [{}/{}]\".format(j, measurements), end='\\r')\n", " \n", " # Sweeps between measurements\n", " # TODO: Skip sweeps between measurements for last measurement in set.\n", " for _ in range(measure_sweeps*sweep_size):\n", " xy_cluster_move(state,beta)\n", " print(\"\\tmeasuring took\", toc(), \" seconds\")\n", " print()\n", " \n", " f.create_dataset(\"square-magn\",data=square_magn)" ] }, { "cell_type": "code", "execution_count": 10, "id": "35782ff3", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "e565b97bdc2fe455715db76f3e958a84", "grade": true, "grade_id": "cell-5d3a2f49615613b2", "locked": true, "points": 10, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "with h5py.File('xy_data.hdf5','r') as f:\n", " assert f[\"square-magn\"][()].shape == (11, 200)\n", " assert_almost_equal(np.mean(f[\"square-magn\"][4]),0.456,delta=0.02)\n", " assert_almost_equal(np.mean(f[\"square-magn\"][9]),0.023,delta=0.01)" ] }, { "cell_type": "code", "execution_count": 11, "id": "1015f826", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "4e39eed6235d63377caca1b6ff0f28da", "grade": true, "grade_id": "cell-1f68f5bac39efe6c", "locked": false, "points": 10, "schema_version": 3, "solution": true, "task": false } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAEtCAYAAAAWZydGAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABANUlEQVR4nO3deXxU1fnH8c83ISEgS5B9CauIrKKgYOuCRSpordraumstFmnr0lq3/n6ttbXW2vZXq1WLuFRrXWvdanErFlFZBBRRQGQTCCggO4SQkDy/P+5NHMJkkrlMMlme9+s1r8xd5tznTu7MM/ece8+RmeGcc85VJiPdATjnnKvbPFE455xLyBOFc865hDxROOecS8gThXPOuYQ8UTjnnEvIE4VzzrmEPFE455xLyBOFa/QkLZQ0qh6V20/Se5J2SLoy1eWH2/hE0kk1UbZLnqSHJP26muum/H/niSIUvrlFktpVmD9fkknqmabQXIpV/CCZ2UAzm5bKMlNVbiWuA6aZWUszu/NAC2uMSUFSi3C/z4uZ11LSaklnSXpU0oMVXnOCpE2SOtd+xOnliWJfK4FzyyYkDQaapS+c9JHUJN0xuEr1ABZGeWF9+L8qcESc+YdLykzFNsxsJzABuENS+3D274C5ZvY0cCVwiqQx4bZzgPuAn5jZp6mIoV4xM38E/V19AvwMmBMz7w/A/wIG9AzndQH+CWwkSCxXxqx/A7Ac2AEsAs6Ms41rgAXANuBJIKeSeK4H1oZlLQFGh/OPAN4N5z8JPAH8OuZ1BhwSM/1Q2fJqxnd9GN8eoEmi/a3kPbw2fP0u4AGgI/BSuM3/AG2q834BRwLvhcv+Ee7rr6vzXlbxP3oEKAV2AzsJfp1/ApwULj87nF/22EPw673SeOOVGRNjWbn9gWnAVoIv+a9HOTaA14ESoDDc1qHVLHuf/2ui96M68SR5XOwE8sLnEwmO0Y7h9LXA/RXW7xWWOzZm3nHA58CgAy2/wmsfAh4HRgGbgM4xy74V7ttBwK3ASyk89qv6n1X6Oa/qvSfmuEvZ92MqC6vPj7I3l+BLuT+QCawh+PVmQE+CM7B5wI1ANtAbWAGcHHNgdQnXOzs8YDpX2MY74ToHA4uBiXFi6Rduu0s43RPoE25zFfBjIAs4Cyim+omiOvHNB/IIzqQS7m8l7+Gs8APSFdgQHuxHAE0JvuR+UeGDuF88Mft5Vbif3wCK2D9R7PdeVifmih+kyj5YQKuw3Muqev/ilcEXx1QWsAz4nzCmrxB8AfRL9tgI150GXBo+r27Z5f/Xyo79OPPixhPhuFgDDAAEfAAsBQ4Lp5cBh8d5zfEEX9wnAkdTIXEcaPkxr20DfEqQhC6Js/xp4IUwlu5VfH9U69iv6n9Ggs95dd77eP/PA3141dP+HgEuAsYAHxH8qi9zFNDezH5lZkVmtoLgdPQcADP7h5mtM7NSM3uS4IA9ukL5d4brbAb+BQyNE0MJwcE1QFKWmX1iZsuBkQQHzp/MrNiCU+Q51d2xJOJbY2a7q9rfSvzZzNab2VrgTWC2mb1nZnuAZwk+OFXFM5LgbObOcD+fIfjSqijeexkl5v1IygAeIzibuLeKeKsyEmgB/DaM6XXgRWKqORPsTyrLLvu/Vldl8ST7Hm8NY/wqwRfkh0AuMBZYa2bvV3yBmU0HziP4on4RmGBmL6eq/JjtbCH4Rd8ceCbOKj8k+CL/lZmtrqycUHWP/ar+Z4k+5yk5vpNV5+sr0+ARYDrB6e/fKizrAXSRtDVmXibBQYGki4CrCc4AIDgY9mkcBz6LeV5A8IttH2a2TNKPgJuAgZJeCcvtQnDgx/YNv6p6u1Xt+NbEPE+4v5VYH/N8d5zpFtWIpxn772dsXGXivZdRYo7nFqAlQV11VfFWpQuwxsxKY+atIvjlGavKY+MAyo73/lWlsniSfY+3ELxXPwJuAy4h+CX/AyBRY/xqYC/Br+hPEqwXtXwkXUDw//xP+NqJscvNbL2kz6lem1B1j/2q/meJPuepOr6T4omiAjNbJWklcAowvsLiNcBKM+tb8XWSehBk9tHATDMrkTSf4PQ3ShyPAY9JagXcS3AQ3w90laSYg6g7Qb15mQKCX0dlOgH5ScRX8cs57v4eqCri+ZT99zOPffezMtWJOeEgLJLOIfh1d5SZFVcj3qrKXAfkScqI+XLoDnxc1c5UQ3XLThRfsoPSJHtcbAWGE1TTTZP0TeBwYBDwXLwXSOoDvEbQtrIDmCLpJDOL94WddPnhNjoAtwPfJqg9WCjpsfBspiZV9T+Ld/yXfc5r7DOZiFc9xTce+IqZ7aow/x1gu6TrJTWTlClpkKSjCBq8jKAuFUmXEByoSQuvk/+KpKYEjZa7CaqjZhL8wrpSUhNJ32D/qo/5wHlhbGOBE8L5UeJLtL8HKlE8Mwn29/JwP0+nelU81Y15PUHd7n7Cq23+DJxhZhurGW/CMoHZBO0Z10nKCu+tOI2ggfJApaLsRLHHk+xxsYWgvr3s1/12gvane82spOLKkroAU4FbzOwhM/snQcP6q5LixZlU+THuAp4zs/9acCXTdcB94eeuJlX1P0v0Oa/Jz2SlPFHEYWbLzWxunPklBP/QoQRXG3xO8Cu/tZktAv6P4J+8HhgMvB0xhKbAb8PyPwM6AP9jZkUEDbvfIfhwnM3+9apXhTFuBc4n/EUVJb5E+xtxv2LLrjSemP0cH+7HBQR1uHuqUW51Yr4V+JmkrZKuqVDE6QTVFm9J2hk+XqrG+1dpmeH+fB0YF8ZzD3CRmX1U1f5UY39TUXai9yPeNpM9LrYQ1F48Fk5vJ2hDuK+S9TcRXIb6l5htPkpQlbQhBeUj6QzgWIIrlcq2cT+QT9BQXGOq+p8l+pzX5GcyEe1bDebqG0kPAflm9rN0x1KTJM0GJpnZX9Mdi3ONjZ9RuDopvAu2U3jqfTEwBKjsqhfnXA3yxmxXV/UDniK4UmQ5cJY1xjtinasDvOrJOedcQl715JxzLiFPFM455xLyROGccy4hTxSNgKRbwy5BnGvUJF0p6bfpjqO+8UTRwIV97V9E0A1IOrZ/raQPFYzGtlLStRWWfyJpd8zNba9G3M5TksZIairps6pf4apD0ihJ+bW8zVMlvRXeAPiZpPsktYxZ3lTSg5K2h8uvrvD6oZLmSSoI/w6NWTwZuCDsvsNVkyeKhu87wJQkew1NJREkqjYEvXleHvalFOs0M2sRPr4acTvDCLp1HkLQe2i9pXowuFANa03QpXYXgi7/uwG/j1l+E9CXoIO8Ewm6whgLICkbeB74O8Ex9zDwfDgfMyskGCPiotrYkYbCE0XDNw54o2xC0hthp2lIOlbBMK+nhNMnKejoLmXM7Hdm9q6Z7TWzJQQf4i+nchuS2hBc6r2JoHO4dyssnybp15JmhGct/5LUVsFwl9slzVHMULeSDpP0mqTNkpZI+nbMslMVjFe9XdIaSTfFLMuR9HcFw2VuDcvtGC7bZ7hRSTdJ+nv4vGf4fxgvaTXB2AVI+q6kxZK2SHpFQceEZa83ST+QtDQ8W7tZUh9JM8PYnir7cgzX/5qCYX23hu/DkJhln0i6RtICSdskPRnuy0EEX6pdYs74ukg6WtLccDvrJf3xQP+HsczsMTN72cwKwm7A72PfY+Yi4GYz22Jmi8Pl3wmXjSK4P+xPZrbHgqFiRdBVeJlpwKmpjLmh80TR8A0mGIypzBsEHyYIBodZwRcdBx5PTFKJJem88Eumskf3qgKRJIKRyir2APqopI2SXpV0eHV3TNJoBd0trwG6hc/vAH4YxnRCzOrnABcSdOXch6DPpr/yxaA8vwjLPIig19LHCPrYOhe4R9LAsJxdBF9UuQRfNt9X0G8QwMUEv4bzgLYEXVYncyZ3AsEv6JPDMv+HoM+f9gTdSD9eYf2xBGdSIwk6tJtM0L9XHkGHheeG+3Qk8CBwWRjXvcAL2rfzu2+H5fUiOCv7Ttgp5jhgXcwZ3zqC9/gOM2sVvpdPxdsZSd2rOGbOi/e6OI4nPGbCHwVdgNgxJt4Hyv4/A4EFFbroXhCzHIL/d7WPM+eJojHIJeimucwb7JsYbo2ZPoFKEkX4Ky83waOqQV0gqDLIIPiCLnM+wXgAPYD/Aq9Iyq3OjpnZVDPLJej48FsESeATgoFdcs0sdl/+Gnb2uI3gV/JyM/uPme0lGGq1bFCZrwGfmNlfw7OgdwmGnTwr3OY0M/vAgsGLFhB8eZe9f8UEX8SHmFmJmc0zs+3V2ZfQTWa2K6wmvAy41cwWhzH+Bhgae1YB3GZm28Outz8EXjWzFTH7WLZP3yPoSXV2GNfDBB0sjowpK5lBk4qBQyS1M7OdZjYr3kpmtrqKY+axeK+LpWDM6ov5oqO+sjEdtsWsto1g7JCy5bHLKi6H4PNQo53oNTSeKBq+Lez7IZkJHBpWiQwlGJwpT1I7gq6Ma6QvfkmXE/wSP9WCEb8AMLO3zWx3WM1wK0FvscdVs8z88CziXIK66A0ECWddnOqQ6g4q0wMYEfvLlyCZdQq3OULSf8MzoG0EZw1lgxc9ArwCPCFpnaTfScqqzr6EKg4adUdMDJsJqlBiByRKZp9+UmGf8th3YKSKgxS1oHLjCcbq/iisXvtaVTsWhaSRBGd2Z5lZ2VgNO8O/rWJWbcUXP4Z2VlhWcTkEn4eKycQl4Imi4VtA8KEGwMwKCMbcvQr4MOzSeAbByG3LzezzeIVIOj+mnjreo9KqJ0nfBW4ARptZVVfQGNUc7MnMuhFUl/wnPLOYDPww/LV6dcIXV24N8EaFX74tzOz74fLHCMZQzjOz1sCksngtGLbyl2Y2APgSwdlJWaPpLvYfUGq/XaoQx2UV4mhmZjMi7tMtFcpqbmYVq7Li2a+PHzNbambnElTN3QY8HVbZ7SOsekp0zJxf2UYVjAvyAvBdM5sas+0tBAP7xFYdHc4X1ZkLgSFhNWeZIexb3dmffauuXBU8UTR8U/iiaqTMG8DlfFHNNK3C9H7M7NGYeup4j7hVT+GXwW+AMRaM7xu7rLukL0vKDhtPryX4df52uLyskbdngv0ru9oJ4Ehgv3FEkvQiwRnXhQoGlcmSdJSk/uHylsBmMyuUdDTBuM5l+3OipMGSMgnGRCgmGIAJggGlzgnLG05YlZXAJOCnZW0jklpL+lbEfboPmBieDUnSQQoa5VtW+crgLKWtpPKqGkkXSGpvwehsW8PZ+w0QFFY9JTpmHo23QUmDCHoKvsLM/hVnlb8RjJ/RRtJhBFVrD4XLpoWxXKngMtrLw/mvx7z+BIKqOVdNnigavr8Bp0hqFjPvDYIvvOmVTKfSrwnq7efE/JKcFC5rCfyFoHpsLcHZwTgLrl6CoHpkVbisMsOAd8NfkIdRvbGNK2VmO4CvEjR+ryOokrmNYDApCAbP+ZWkHQT15rENuZ2ApwmSxGKC9/Xv4bKfEzT8bgF+yReD7FQWx7Phdp+QtJ2gDWJcxH2aS/Blele4/WV8cZVQVa/9iKAdZkVYbdWF4P+0UNJOgobtcyy47DRVfkLQgP9AzDET+3/9BUGPwqsI3uPfm9nLYbxFwBkEZ3Jbge8SjFZYBMGVaQTDHD+cwngbPO89thGQ9Btgg5n9Kd2xJEPSz4CNZpaWmwVdwyPpCoJqw+vSHUt94onCOedcQl715JxzLiFPFM455xLyROGccy6hBtf5WLt27axnz57pDsM55+qVefPmfW5m7eMta3CJomfPnsyde6CX0jvnXOMiaVVly7zqyTnnXEKeKJxzziXkicI551xCniicc84l5InCOedcQp4onHPOJeSJwjnnXEKeKJxzziXU4G64S9btr33MHVOXVrr8qtF9+fGYQytd7pxzDV2D62Z8+PDhFvXO7LPvnQnAk5cdk8qQnHOuzpM0z8yGx1vmVU/OOecS8kThnHMuIU8UzjnnEvJEESopNbYUFLF2y26mLl5PSWnDartxzrmoPFEQJIkLH5jNsg07yd+6mysef48LH5jtycI55/BEAcC0JRuYv2YrZXmhoKiE+Wu2Mm3JhvQG5pxzdYAnCmDhuu3sLirZZ97uohIWrduepoicc67u8EQBDOzSimbZmfvMMyB/624Ki0viv8g55xoJTxTAqH4dGJqXS4aC6aZNMmjdrAlPzlnDuffNSm9wzjmXZo2+Cw+AzAzxyPgRjLtjOgV7Svjl6QMZ1a8Dcz/ZTEFYJVVYXMI/5q7hW8PzyMnKrKJE55xrODxRhDIzRJvm2bRpDqP7dwRgRO+25ctf/2gDP39+IXe+voyJJ/Th/BHdPWE45xqFtFY9SRoraYmkZZJuiLP8Wknzw8eHkkokHZyOWE8Z3JknJ4ykb4cW3PziIo773X954K2Vfgmtc67BS9sZhaRM4G5gDJAPzJH0gpktKlvHzH4P/D5c/zTgx2a2OR3xQnCG8VjvtsxasYk7/rOUf72/ju9+uScApaVGRlkjRzV4r7XOufoinVVPRwPLzGwFgKQngNOBRZWsfy7weC3FltDI3m0ZOaEtOwqLkcSmnXv45l9mcOExPatdJfXjMYeWJwLvtdY5V5els+qpK7AmZjo/nLcfSc2BscA/K1k+QdJcSXM3btyY8kAr0zInC4DthXvp3LrZPlVSflmtc66hSOcZRbx6msoq/E8D3q6s2snMJgOTIRiPIpkg4lUB9bzh3+XPq1MF1KvdQTw+YWR5ldTNLy5i8vTlvHb1CbQKk4lzztVX6UwU+UBezHQ3YF0l655DDVU7xVYBHaiyKqnZKzbxzsrN5Univ0s2cEzvtn6VlHOuXkqYKCSVmFlNfbvNAfpK6gWsJUgG58WJoTVwAnBBDcWRciN6ty2/tHbN5gK++9Ac2rVoyvdP6MN5flmtc66eqaqNovqX8STJzPYClwOvAIuBp8xsoaSJkibGrHom8KqZ7aqpWGpS3sHNefx7IzmkfQt+5W0Yzrl6qKqqpxq9ScDMpgBTKsybVGH6IeChmoyjppVVSZW1Yfzu5Y84bUhncrIyy8fBKNhTwtTF6xnVrwOZSVxm65xzNU1mleeC2KonSWOAbwN3m9l8SRPCRuQ6Zfjw4TZ37tx0h5HQms0F5B3cnL0lpYz4zVQ27yrCgObZmQzNy+WR8SM8WTjnapWkeWY2PN6yZC6P/QFwLXCBpK8AQ1MQW6OUd3BzAF7+8DO2FBSVn7b5OBjOuboomUSx0cy2mtk1wFeBo2oopkZjxee7qHhCV1BUwozln6cnIOeciyOZRFF+c4GZ3QD8LfXhNC7xxsHIEAzr0QaA4pLSdITlnHP7qHaiMLPnK8yaFHdFV20Vx8Fonp3JyN5tOXlgZ4r2lnLKHW/y6xcXsa2gOL2BOucatUhdeEi6H1gvaY2k2ZLuk3RFimNr8MrGwTikQwu65Tbjz+ceUd6QXbi3hGE92vDA2ys54Q//5aG3V/oZhnMuLaL29XQc0NHM8oBvAM8CB6UsqkakbByMrm2aMbp/x/KrnVrlZPHbbw7h31ccx4DOrbjpX4s4+fbpfLptd5ojds41NlG78JgFtAE2mNlagjurpyR+iYtiQJdWPHrpCF7/aAPPvreWji1zANheWOz9SDnnakXUM4rJwBuSrpF0XNjNhqshkhjdvyN3nXckGRliy64iTvjdf7n+6QVs2FGY7vCccw1c1DOKvxNc9dSE4P6KIZJyzKxPyiJzlcrMFN88shsPz/yEfy1Yxw9G9eHS43pXuw8pHzTJOZeMat+ZXWH+m2Z2XIV5Tc1sTw3EmJT6cGd2RVEHLlr5+S5++9JiXlm4ni6tc5hy1XHkNs+ulW075xqWRHdmRz2jmC/pKjO7o2xGXUgS9UmqxsG498LhzFqxiTeXbixPEvlbCujWpnnqg3bONUpRE0VH4CRJ1wPvAu8D883sHymLrIFL+TgYYbfmyzbsZOyfpnPywE5cP/Ywurf1hOGcOzCRuhk3s2+bWX+gF3Aj8DHBGNguzbrk5nD5Vw7h9Y82cNIf3+DWKYvZXug37Dnnokt4RmFmCRNJWN30bvhwdUDz7Cb86KRDOeeo7vzh1SVMfnMFz89fx7RrR/mASc65SNI5FKqrQZ1a5/CHbx3Od77UkwX528qTxPtrtnJ4Xm56g3PO1StR76Nw9cSgrq05b0R3AGYu38Tpd7/NRQ++w5LPdpQPmrR2y26mLl5PSWmNjlPlnKunEl4eW+Mbl8YCdwCZwP1m9ts464wC/gRkAZ+b2QmJyqyPl8fWlqK9pfxt5ifcOXUpOwr30qFVUzbu2EOp+aBJzjV2Kb88VlJT4JtAz9gyzOxXSZSRCdwNjAHygTmSXjCzRTHr5AL3AGPNbLWkDlHidYHsJhlcelxvvnlkN675x/tM/eiLAZJiB00a3b9jGqN0ztU1UauengdOB/YCu2IeyTgaWGZmK8ysCHgiLDPWecAzZrYawMx86LcUaHNQNofn5e53SVtBUQkfrt2Wlpicc3VX1MbsbmY29gC33RVYEzOdD4yosM6hQJakaUBL4A4z22/AJEkTgAkA3bt3P8CwGoeyQZMKikr2mf/EnDUM6tqarxzWAcmroJxz0c8oZkgafIDbjvctVLHBpAkwDDgVOBn4uaT97lIzs8lmNtzMhrdv3/4Aw2oc4g2a1K9jS5o2yWD8w3M5e/Is3lu9Jb1BOufqhKiJ4lhgnqQlkhZI+kDSgiTLyAfyYqa7AevirPOyme0ys8+B6cDhEWN2MeINmjTlquN47eoTuPmMQazYuJPrnl5AqV8J5VyjF7XqaVwKtj0H6CupF8F4FucQtEnEeh64S1ITIJugaur2FGzb8cWgSW2aU96AnYm4cGQPzjyiK59t201Ghti5Zy9/nrqU7x3fm3YtmqY5audcbYuUKMxslaTDCUa6A3jTzN5Psoy9ki4HXiG4PPZBM1soaWK4fJKZLZb0MrAAKCW4hPbDKDG75LRo2oRDOrQEYNbyTdz/1kr+PmsVE0/ow/jjetE8O/lDx7s3d65+inQfhaSrgO8Bz4SzzgQmm9mfUxhbJH4fRXKq28348o07+d3LH/HKwvV0aNmUH510KOcenRe5wdu7N3eubkl0H0XUNorxwAgzu9HMbgRGEiQO10D1ad+Cey8cztMTjyHv4Oa89OGnflWUc41E1DYKAbHXVZZQSU+zrmEZ3vNgnp54DDv27AWCsS+u/+cCrh7Tj2E92qQ5OudcTYiaKP4KzJb0bDh9BvBASiJydZ4kWuVkAbB6UwEfr9/JN/8yg3GDOnHtyf3o3b5FmiN0zqVS1MbsP0p6A/gywZnEJWb2XkojczUmFaPrlfnSIe2Yds0o7n9zJZOnL+fVReu5cGQPfnHaAK+acq6BiNzNuJnNA+alMBZXS1I5uh7AQU2bcNVJfTlvRHfunLqUUrPyJLFnbwlNm/g4GM7VZ0klCklvmdmxknaw713UAszMWqU0OlevtG/ZlJvPGETZlXTzVm1h4t/ncdXovpx9VB5Zmd6rvXP1UVKfXDM7Nvzb0sxaxTxaepJwZcrOJpplZdKr7UH87LkPOfn26bz84WeYmY+D4Vw9E/U+itvM7Pqq5qWD30dRt5gZ/1m8gdte/ohlG3ZywqHtKC4xZq3Y5ONgOFeH1MR9FGPizEtFtx6ugZHEmAEdefmq4/jtNwaTd3Bz5q/ZStlJROw4GM65uimpRCHp+5I+APqFnQGWPVYCH9RMiK4haJKZwTlHd6dDyxx2V+jafHdRCYvWbU9TZM65qiR71dNjwEvArcANMfN3mNnmlEXlGqx442BI0KGVdzboXF2VbGP2NjP7xMzOBbYDHYEewCBJx9dEgK5hqTgORnZmBgJ+9a9FvPB+xV7mnXN1QaQ2CkmXEowN8Qrwy/DvTakLyzVUFcfB+MsFRzL1J6MY3K01OU388lnn6qKoN9xdBRwFzDKzEyUdRpAwnKtSvHEwHv/eyPLLah+e8QkdWjZl3ODO6QzTOReK+hOu0MwKASQ1NbOPgH6pC8s1NmVJoqTUeOH9dXz/0Xf50RPvsa2gOM2ROeeiJop8SbnAc8Brkp5n/2FMnUtaZoZ4YsJIfnRSX/614FNO/tN0pn+8Md1hOdeoRUoUZnammW01s5uAnxP0HHt6KgNzjVdWZgY/OulQnv3Bl2iR04RLHppD/paCdIflXKMVqY1C0lTg/8xsipm9Ec6bDExIspyxwB0EQ6Heb2a/rbB8FMG42SvDWc+Y2a+ixOzqnyHdcnnximOZuWIT3do0B2DN5gLyDm6eVDk+BKtzByZqY3Yv4HpJR5lZWSN23Fu/KyMpE7ib4C7vfGCOpBfMbFGFVd80s69FjNPVczlZmZzYrwMAM5Z9zoUPvsNlx/fmqpP6VrtX2tjecn0IVueSF7WNYiswGugo6V+SWkco42hgmZmtMLMi4Am8+solMLhba846shv3TFvO6Xe9zeJP/W5u52pD5KFQzWwv8ANJ3wHeApIdB7MrsCZmOh8YEWe9YyS9T9BYfo2ZLdwvGGkCYbVX9+7dkwzD1ZYDHTCpZU4Wt501hK8O7Mj1//yAr9/1FjeM68/4Y3vVWMzOueiJYlLZEzN7KOz/6YdJlhGvq9CKXdm+C/Qws52STiG4yqrvfi8ymwxMhqD32CTjcLUkVQMmje7fkVd/3IafPfcB2X6TnnM1LupQqPdWmJ4HfDfJYvKBvJjpblS4xNbMtsc8nyLpHkntzOzzJLflGpiDD8rm7vOOLJ9+fv5athfu5YIR3X0IVudSLNneY98K/+6QtD3msUNSshXGc4C+knpJygbOAV6osL1OCj/1ko4O492U5HZcAyWpPCm8umg9P3/uQy568B0+3bY7zZE517CkbYS7sI3jcoJ+ohYDT5nZQkkTJU0MVzsL+DBso7gTOMeijLTkGry7zj2Cm88YxNxPtnDy7dN57r21+KHiXGpE7RTwturMq0p4H8ahZtbHzG4J500ys0nh87vMbKCZHW5mI81sRpR4XcMniQtH9uClq46jb8eW/OjJ+cxfszXdYTnXIPgId65B6dnuIJ667BgeuHg4R3QPLsRbvcnv6nbuQPgId67BycxQea+0H6/fweg/TuO6p99na0ERWwqKWLtlN1MXr6ek1KumnKsOJVOPG95Y14Y6PMLd8OHDbe7cuekOw9URe/aWcMd/lvKXacvJysyguKQUA5pnZzI0L5dHxo8gM8OvknJO0jwzi9vDRqQR7oDzgeOAi81sFdAivCrJuTqlaZNMrht7GDeMO6w8SQAUFJUwf81Wpi3ZkNb4nKsPorZR3A0cA5wbTu8I5zlXJ+3ZW7rfvN1FJSxa592AOFeVqIlihJn9ECgEMLMtQHbKonIuxQZ2aUWz7H07EczJyqR/55Zpisi5+iNqoigOe381AEntgf1/sjlXR4zq14GhebmUNUc0y8qgqKSURd6xoHNVipoo7gSeJeg99haCTgF/k7KonEuxzAzxyPgRHNKhBd1ym3HnuUfw9SGd+eNrS3lk1qp0h+dcnRa1r6dHJc0j6Goc4AwzW5y6sJxLvcwM0aZ5Nm2aw5gBnRjVrwM79uzlxuc/pHWzLL5+eJd0h+hcnRT1zuymwJFAa6At8C1JN6YyMOdqWlZmBneddyRH9TyYq5+cz7xVdeIKb+fqnKjdjD8PbAPmAXtSF45ztSsnK5P7Lx7O3f9dxsAuUcbfcq7hi5ooupnZ2JRG4lyatMrJ4qfj+gOwraCYjTv3cEiHFmmOyrm6I2pj9gxJg1MaiXN1wJVPvMd5983y/qGcixE1URwLzJO0JOzr6QNJC1IZmHPp8LNT+1NUUsqFD85mw47CdIfjXJ0QNVGMIxiS9KvAacDXwr/O1Wt9O7bkoUuOZuOOPVz0wDts212c7pCcS7tIicLMVsV7pDo459JhaF4uky8czvKNO/n5cx+mOxzn0i5SY7akq+PM3gbMM7P5BxSRc3XAsX3bMfmi4fTvlNTAjc41SFGrnoYDE4Gu4WMCMAq4T9J11S1E0tiwnWOZpBsSrHeUpBJJZ0WM17mkndivA51a51BSajw5ZzWlPn6Fa6SiXh7bFjjSzHYCSPoF8DRwPMG9Fb+rqoCwr6i7CUbLywfmSHrBzBbFWe82grG1nUva7a99zB1Tl+4zr+cN/y5/ftXovvx4zKGVvv61RZ9x/T8/YOG67fzy6wORfPwK17hETRTdgaKY6WKgh5ntllTdG/COBpaZ2QoASU8ApwOLKqx3BfBP4KiIsbpG7sdjDk2YCKpy8sBOTDi+N5OnryC3eTZXH0BZztVHURPFY8AsSc+H06cBj0s6iP2/6CvTFVgTM50PjIhdQVJX4EzgKyRIFJImEFR/0b1792pu3rnqkcRPxx3G1oIi7py6lDbNs7jky72q9dp4ZzOxqjqbca4uiNop4M2SphDcTyFgopmVjT96fjWLiXf+XrES+E/A9WZWkuh038wmA5MhGAq1mtt3rtok8ZszB7O1oJjfv7KEU4d0pkPLnCpfF3s2c/a9MwF48rJjajRW51It6hkFwAogE8gBmks63symJ/H6fCAvZrobsK7COsOBJ8Ik0Q44RdJeM3suctTORdQkM4M7zz2C5Rt3VitJONdQRL089lLgKoIv9/nASGAmQRVRdc0B+krqBawFzgHOi13BzMrP7yU9BLzoScKlU05WZnnngY/NXk2f9gcxonfbNEflXM2KennsVQRtBqvM7ETgCGBjMgWY2V7gcoKrmRYDT5nZQkkTJU2MGJdztaKwuIQH3lrBpQ/P5cO129IdjnM1KmqiKDSzQgjGpjCzj4B+yRZiZlPM7FAz62Nmt4TzJpnZpDjrfsfMno4Yr3MplZOVySPjR9Aypwnf+es7rPx8V7pDcq7GRE0U+ZJygeeA18Krnyq2LzjXoHXJbcYjl46g1OCC+2fz2TbvRNA1TFH7ejrTzLaa2U3Az4H7Ce6BcK5R6dO+BQ9fcjTbdhcz/eOkal+dqzeiNmYPB/4X6BGWIeAWYEjqQnOufhjcrTX/vWYU7Vs2TXcoztWIqJfHPgpcC3wAlKYuHOfqp7Ik8c7KzTw0YyW3nz2Upk0y0xyVc6kRNVFsNLMXUhqJcw3Aqk27mPLBZ4j3ufPcI8jM8H6hXP0XNVH8QtL9wFSgvG8nM3smJVE5V099a3geWwuKuWXKYlo1y+I3Zw7yTgRdvRc1UVwCHAZk8UXVkwGeKFyj973je7OloIh7pi2nTfMsrht7WLpDcu6ARE0Uh5vZ4JRG4lwDcu3J/di6u5jVmwsoKTWvgnL1WtREMUvSgIpjRzjnApK4+fRBAGRmiD3FJWwpKKJgTwlTF69nVL8OnjxcvRE1URwLXCxpJUEbhQAzM7881rlQWSJYu2U3J/3xDQqLSzDgisffY2heLo+MH+HJwtULURPF2JRG4VwDNn/NFvbsLSnvQ7+gqIT5a7YybckGRvfvmNbYnKuOqONRrEp1IM41VMs37sIqjJKyu6iEReu2e6Jw9ULUvp6cc9U0sEsrmmXve/Nds+xMBnRplaaInEuOJwrnatiofh0YmpdLWXNEhmBIt9aM6tchvYE5V02eKJyrYZkZ4pHxIzikQwu65uYw6YJhPHrpSG/IdvVGpEShwAWSbgynu0s6OrWhOddwZGaINs2z6damOV8d2Iltu4t5as6adIflXLVEPaO4BzgGODec3gHcnZKInGsEHp7xCdf9cwFvLf083aE4V6WoiWKEmf0QKAQwsy1AdrKFSBoraYmkZZJuiLP8dEkLJM2XNFfSsRHjda5O+f6oPvRufxDX/3MBO/fsTXc4ziUUNVEUS8ok6N8JSe1Jsrvx8PV3A+OAAcC5kgZUWG0qQXchQ4HvEgyQ5Fy9l5OVye/POpx123Zz65TF6Q7HuYSiJoo7gWeBDpJuAd4CfpNkGUcDy8xshZkVAU9QYZQ8M9tpVn4F+kFAhavRnau/hvVow/gv9+LR2auZscyroFzdlfQNdwr6TJ4OzANGE3TfcYaZJfuzqCsQ25qXD4yIs70zgVuBDsCplcQ0AZgA0L179yTDcC59fvLVfhSVlHJIhxbpDsW5SiWdKMzMJD1nZsOAjw5g2/GuDdzvjMHMngWelXQ8cDNwUpx1JgOTAYYPH+5nHa7eaJadya/CzgPNzMeucHVS1KqnWZKOOsBt5wN5MdPdgHWVrWxm04E+ktod4Hadq3PWby/knMmzmLViU7pDcW4/URPFicBMScvDq5I+kLQgyTLmAH0l9ZKUDZwD7DO8qqRDwqouJB1JcGWVf5Jcg9MypwmfbivkuqcXUFDkV0G5uiVq77HjDnTDZrZX0uXAK0Am8KCZLZQ0MVw+CfgmcJGkYmA3cHZM47ZzDUbz7Cb87qwhnDN5Fr9/ZQm/OG1gukNyrlzk3mMltQH6Ajkxi5LqVdbMpgBTKsybFPP8NuC2KDE6V9+M7N2Wi47pwUMzPmHcoM4c3evgdIfkHBC9C49LCa58egX4Zfj3ptSF5VzjdP3Yw+ia24w/v7403aE4Vy5q1dNVwFHALDM7UdJhBAnDOXcADmrahL9+5yi65DZLdyjOlYvamF1oZoUAkpqa2UdAv9SF5Vzj1bdjSw5q2oTC4hJWbdqV7nCci3xGkS8pF3gOeE3SFhJc2upcY3X7ax9zx9R9q5F63vDv8udXje7Lj8ccGve1Ex6ZR/6WAqZceRw5WZlx13GuNuhALyKSdALQGnjJzIpTEtUBGD58uM2dOzfdYTh3wN5cupELH3iHy07ozU/H9U93OK6BkzTPzIbHWxbpjKJsHIoKhgK/ilKec25/x/Vtz7lH53Hf9BWMHdiJI7q3SXdIrpGK2kaxK+ZRQnBfRc8UxeScC/30lP50bJXDdU8voLC4JN3huEYq6n0U/xc7LekPVLir2jl34FrlZHHrNwbzp/8sZWtBMZ1ae1uFq31RG7Mrag70TlFZzrkYo/p14Pi+7cnwMbZdmkS94e6DsI+nBZIWAkuAO1IbmnOuTEaG2LyriFtfWsyevV4F5WpX1DOKr8U83wusNzPvycy5GvT+mq3c+8YKmmZmcPVX/bYlV3sinVGY2aqYx1pPEs7VvBMP68A3juzK3dOW8+HabekOxzUiUS+PvTrRcjP7Y7RwnHOJ3Pi1Aby59HOu+cf7vHD5sWQ3iXrhonPVF7XqaThBX09lVzqdRtBJ4JpKX+GcO2C5zbP5zZmD+d7f5nLfmyv44YmHxF0v3h3hsRLdEe5cRZHuzJb0KvBNM9sRTrcE/mFmY1McX9L8zmzXGPz17ZWcOqQzHVrmVLnu2ffOBODJy46p6bBcPZbyO7OB7kBRzHQRfsOdc7Xmki/3AqC01Cg1o0mmV0G5mhM1UTwCvCPp2XD6TOBvqQnJOVcdu4tKuPjBdziubzuuGN033eG4BizqVU+3AJcAW4DNwMVm9ptky5E0VtISScsk3RBn+fkx92vMkHR4lHida4iaZWfSoVVT7nx9KUs+25HucFwDFvWGu28BS83sDoKeY2+UdESSZWQCdxP0EzUAOFfSgAqrrQROMLMhwM3A5CjxOtdQ/fLrA2mVk8U1/3ifvSWl6Q7HNVBRKzZ/bmY7JB0LjAEeBiZV8ZqKjgaWmdkKMysCngBOj13BzGaY2ZZwchbQLWK8zjVIbVs05eYzBvHB2m1MfnNFusNxDVTURFHWh8CpwCQzex7ITrKMrux7OW1+OK8y44GX4i2QNEHSXElzN27cmGQYztVvpwzuzCmDO/HMu2sp9rMKVwOiNmavlXQvcBJwm6SmJJ904vVwFvdaXUknEiSKY+MtN7PJhNVSw4cPP7CRmJyrh35z5mCym2SQ5Vc/uRoQ9aj6NvAKMNbMtgIHA9cmWUY+kBcz3Y04w6lKGgLcD5xuZpsiRetcA5fbPJvm2cE42zOX+8fEpVbUq54KzOwZM1saTn9qZq8mWcwcoK+kXpKygXOoMKaFpO7AM8CFZvZxlFida0x++9JHXPzXd1i2YWe6Q3ENSNrOU8OOBC8nODNZDDxlZgslTZQ0MVztRqAtcI+k+ZL8lmvnEvjBiX1onp3JdU+/T0mp18K61EhrhaaZTTGzQ82sT3hvBmY2ycwmhc8vNbM2ZjY0fMS9vdw5F+jQModfnDaAd1dv5a9vr0x3OK6BiNp7bFPgmwTddpSXYWa/Sk1YzrmozhjalX8v+JTfv7KE0f07pjsc1wBEPaN4nuCeh73ArpiHcy7NJHHLmYMZ0bstUTr9dK6iqJfHdqsLPcU65+Lr2CqHv333aEpKjS0FRRTsKWHq4vWM6teBTB972yUp6hnFDEmDUxqJcy6lSkqNcyfP5OP1O8nfupsrHn+PCx+Y7Y3cLmlRE8WxwLywQ78Fkj6QtCCVgTnnDsy0JRv4YO328umCohLmr9nKtCUb0hiVq4+iVj2NS2kUzrmUW7huO4XFJfvM211UwqJ1272R2yUlUqIws1WpDsQ5l1oDu7SiWXYmBUVfJAsDmmZ5Nx8uOVHPKJDUBugLlI/FaGbTUxGUc+7AjerXgaF5ucxasYlSg5ysoC+osYM6pzs0V89EHY/iUmA6wV3Vvwz/3pS6sJxzByozQzwyfgSHdGhBt9xm3H3ekbz38zF0P7g5paXGQ2+v3K9qyrl4op6DXgUcBawysxOBIwDv39u5OiYzQ7Rpnk3XNs0Y3b9j+djas1du5qZ/LeLb985k3dbdaY7S1XVRE0WhmRVCcJe2mX0E9EtdWM65mnRMn7bce+EwVmzcxdfveovZK7zHWVe5qIkiX1Iu8BzwmqTnidNFuHOu7jp5YCee++GXaJWTxfn3z+apuWuqfpFrlKJe9XRm+PQmSf8lGDf75ZRF5ZyrFYd0aMlzl3+Z6/6xgF7tDkp3OK6OitqYLUkXSLrRzN4A5gNDUxmYc652tMrJYtKFwziq58EAPDZ7tbdbuH1ErXq6BzgGODec3gHcnZKInHNp8/nOPdw6ZbG3W7h9RE0UI8zsh0AhgJltAbJTFpVzLi3atWjKszHtFg/P+MR7oHWRE0WxpEyCGz2R1B4oTVlUzrm0KWu3OOHQ9vzihYX873Mfpjskl2ZRE8WdwLNAR0m3AG8DtyZbiKSxYceCyyTdEGf5YZJmStoj6ZqIsTrnktQqJ4v7LhrOlaP7Mrhr63SH49Is6lVPj0qaB4wOZ309vJei2sIzkruBMUA+MEfSC2a2KGa1zcCVwBlR4nTORZeRIa4ec2j59MsffsbBB2VzdK+D0xiVS4eoQ6EOB/6XL4ZCvUwSZjYkiWKOBpaZ2YqwzCcIRs0rTxRmtgHYIOnUKHE651KjtNT48+tLWfLZDm48bQAXjuyBFH8ApNtf+5g7pi6ttKyrRvflxzEJyNV9UTsFfBS4FviA6G0TXYHYO3zygRFRCpI0AZgA0L1794jhOOcqk5EhHp8wkh8/MZ8bn1/IB/nbuPmMQeRkZe637o/HHFqeCM6+dyYAT152TK3G61IrahvFRjN7wcxWmtmqskeSZcT7ORLp8gozm2xmw81sePv27aMU4ZyrQnm7xVcO4R/z8jn73pnsLvJOBRuDqGcUv5B0PzAV2FM208yeSaKMfCAvZrob3g2Ic3VaRoa4+qv9GNi1Ne+t3kqz7P3PKFzDEzVRXAIcBmTxRdWTAckkijlAX0m9gLXAOcB5EeNxzsWI107Q84Z/lz8/0HaCkwd24uSBnQD4cO023lu9hQsStFu4+i1qojjczAYfyIbNbK+kywnGssgEHjSzhZImhssnSeoEzAVaAaWSfgQMMLPtlZXrnNu3naCmPf7Oah6dvZoFCdotXP0WNVHMkjSgwqWsSTOzKcCUCvMmxTz/jKBKyjlXR918+iDaHpTNna8v4+MNO5l0wZF0bt0s3WG5FIramH0sMD+8WW6BpA8kLUhlYM65+qGs3eLeC4exbP0OTvvzW6z8fFe6w3IpFPWMYmxKo3DO1XsnD+zE85d/mfvfXElem2aUlBpbCooo2FPC1MXrGdWvA5kZ3oZRH6mhdfg1fPhwmzt3brrDcK5RKyk1zpk8k7mfbMGA5tmZDM3L5ZHxIzxZ1FGS5pnZ8HjLolY9OedcpaYt2cD7a7aV3xhVUFTCOys3c9frS9lb4v2H1jeeKJxzKbdw3XaKKySEvaXGn19fRmmYPaYt2cCMZZ9TWOw37dV1UdsonHOuUgO7tKJZdiYFMXduN8/O5CdjDiW7SfD79I+vfcyC/G1kN8ngyO65HNO7Hccf2o4jureJtE3vY6rmeBuFcy7lSkqNCx+YzawVmyi1+G0U2wuLmbNyMzOXb2Lmik0s+nQ7Ywd24i8XDAPgbzM/YVDX1gzp2pommclVfngfU8lL1EbhZxTOuZTLzBCPjB/BuDumU7CnhF+ePnC/q55a5WQxun9HRvfvCMDWgiJ2FO4FYMOOQm58fiEAB2VnclSvgzmmd1vGDupEj7YH1f4ONXKeKJxzNSIzQ7Rpnk2b5pQng0Rym2eT2zwYUblDyxzm/uwkZq/YzIzlnzNzxSamLdlIuxZN6dH2INZsLuDVRev5Up+29OvYkoyYBOSX5aaeJwrnXJ3UrkVTTh3SmVOHdAZg/fZCmoedEM5euZmbXww6hmjTPIuRvdtyTJ+2nHZ4F3746Lss27CTUoMrHn/PL8tNAU8Uzrl6oWOrnPLnZw3rxjF92gbtG8s3MWvFJl5e+Bm5zbKYv2Zr+ZVVBUUlvLt6C8/NX8uZQ7vuc+bhqs8ThXOuXuqa24yzhnXjrGHdMDM+3VbI0/Py9xsjo7C4lJ889T63vfQR7/zvSQDc/+YK1m8vpGOrHDq3bkan1jl0zQ3+RpGuK65qa7ueKJxz9Z4kuuQ2i3tZbtMmGXxrWDf6dGhRPu/d1VuYungDe/Z+ca/HwC6t+PeVxwFw3dPvs333Xjq1zqFz6xw6tc6hT/sWDOraOu72y3rrLSm1hA34qVZbowl6onDONRij+nVgaF7ufpfl/vL0Qft8Yd9z/jDMjG27i/l0WyGfbSskdiiNPXtLWbZxJ28u3ciuMOmc1L8j918cXD067o43ycyATq2alSeSod1yuXvasi/aRx57jyF5rfn7+BFJX96brJpuwPf7KJxzKZXuG99S/at+R2FxeSI5pENLzIyfPfcha7fu5rNthXy2vZCtBcWc2K89s1du3udspkx2kwwmntCHq8ccyvbCYr5xzwyaNskgJyuTpk0yaNokg7OG5XHqkM5sLSji9tc+/mJZ+PdLfdoxoEsrthcWM3P5pvLXN8kQv/73Ihbkb6v0npXq8PsonHO1pjYHTYon2ctyq9IyJ4uWOVnl05K45cx9x23bXVTCn19fyrQlG/d7/TG923J4Xi5HdM8FwAwO7diCPcWl7NlbSmFxCTsK97KrKLiHZEfhXp5/fx2FxSXs2VtK2W/5m88YxIAurVizuYDLHplXabwFRSXMX7OVaUs2pGT/wROFc84dsGbZmQzr0SZutyWXHtdrny/s1s2yuOf8YZWWlXdwc+bf+FUAzIyikiChZIfVV33at+DFK45lz95S9hSX8I95+Tz73tp9ythdVMKiddtTlijS2imgpLHh4EfLJN0QZ7kk3RkuXyDpyHTE6ZxzVSlrHymr7SmrAhrVr0PkMiXRtEkmrXKyyoeYzcnKZFDX1gzr0YYvHdKOrw3pXH5/SZlm2ZkM6NIq8nYrStsZhaRM4G5gDJAPzJH0QoXhVccBfcPHCOAv4V/nnNtHvLaRnjf8u/x5TbaNxNt2QVEJM5Zvos//TKnRbVfWgH8gCaqitDVmSzoGuMnMTg6nfwpgZrfGrHMvMM3MHg+nlwCjzOzTysr1xmznXGOTigb8utqY3RVYEzOdz/5nC/HW6QrskygkTQAmAHTv3j3lgTrnXF0U70xm/MNf/FBuCDfcxUt3FU9vqrMOZjYZmAzBGcWBh+acc3VfbV1hls7G7HwgL2a6G7AuwjrOOedqUDoTxRygr6RekrKBc4AXKqzzAnBRePXTSGBbovYJ55xzqZe2qicz2yvpcuAVIBN40MwWSpoYLp8ETAFOAZYBBcAl6YrXOecaq7TecGdmUwiSQey8STHPDfhhbcflnHPuC2m94c4551zd54nCOedcQp4onHPOJeSJwjnnXEINbjwKSRuBVemOI4J2wOfpDqKW+T43Dr7P9UMPM2sfb0GDSxT1laS5lfWz0lD5PjcOvs/1n1c9OeecS8gThXPOuYQ8UdQdk9MdQBr4PjcOvs/1nLdROOecS8jPKJxzziXkicI551xCnihqmaSxkpZIWibphkrWGSVpvqSFkt6o7RhTrap9ltRa0r8kvR/uc73uJVjSg5I2SPqwkuWSdGf4fiyQdGRtx5hq1djn88N9XSBphqTDazvGVKtqn2PWO0pSiaSzaiu2VPNEUYskZQJ3A+OAAcC5kgZUWCcXuAf4upkNBL5V23GmUnX2maCH4EVmdjgwCvi/cIyS+uohYGyC5eOAvuFjAvCXWoippj1E4n1eCZxgZkOAm2kYjb0PkXify47/2wiGU6i3PFHUrqOBZWa2wsyKgCeA0yuscx7wjJmtBjCzDbUcY6pVZ58NaClJQAtgM7C3dsNMHTObTrAPlTkd+JsFZgG5kjrXTnQ1o6p9NrMZZrYlnJxFMFplvVaN/zPAFcA/gXr9OfZEUbu6AmtipvPDebEOBdpImiZpnqSLai26mlGdfb4L6E8wzO0HwFVmVlo74aVFdd6Thmw88FK6g6hpkroCZwKTqlq3rkvrwEWNkOLMq3h9chNgGDAaaAbMlDTLzD6u6eBqSHX2+WRgPvAVoA/wmqQ3zWx7DceWLtV5TxokSScSJIpj0x1LLfgTcL2ZlQQny/WXJ4ralQ/kxUx3I/gVXXGdz81sF7BL0nTgcKC+Jorq7PMlwG/DEQ2XSVoJHAa8Uzsh1rrqvCcNjqQhwP3AODPblO54asFw4IkwSbQDTpG018yeS2tUEXjVU+2aA/SV1CtsrD0HeKHCOs8Dx0lqIqk5MAJYXMtxplJ19nk1wRkUkjoC/YAVtRpl7XoBuCi8+mkksM3MPk13UDVJUnfgGeDCenx2nBQz62VmPc2sJ/A08IP6mCTAzyhqlZntlXQ5wRUQmcCDZrZQ0sRw+SQzWyzpZWABUArcb2YJL7+ry6qzzwRXwTwk6QOCapnrzay+ddFcTtLjBFdvtZOUD/wCyILy/Z0CnAIsAwoIzqjqtWrs841AW+Ce8Bf23vreu2o19rnB8C48nHPOJeRVT8455xLyROGccy4hTxTOOecS8kThnHMuIU8UzjnnEvJE4ZxzLiFPFM455xLyROEaLEm5kn6Q7jiqUtNxSvpeOL7JfEmlMc//WFPbdA2L33DnGixJPYEXzWxQHYhFBJ+3/XrFjRJnovISvKYrMMPMelT3Nc6Bn1G4hu23QJ/w1/PvJV0g6Z1w+t5wUBkk9ZT0kaT7JX0o6VFJJ0l6W9JSSUfHrPNwOErb02FfXIRl7Fd2+JrFku4B3gXyJD0Xdh+/UNKESuLsGTtqmqRrJN1USXlx96kSgwi6cXcuKZ4oXEN2A7DczIYCDwJnA18Op0uA82PWPQS4AxhC0HPteQRdYV8D/E+4Tj9gcjhK23bgBwCS+icoux/BIEVHmNkq4LtmNoygZ9ErJbWNjdPMrq1in8rLA5pXsU8VDQbqbb9hLn28U0DXWIwmGOdjTtgpXTP2HXVspZl9ACBpITDVzCzsqLBnuM4aM3s7fP534ErgDwnKng6sCkexK3OlpDPD53kEw6F+lsR+xJZX1T5VNAh4LYltOQd4onCNh4CHzeynlSzfE/O8NGa6lC8+JxUb9Mqm45Ydtj3sipkeBZwEHGNmBZKmATlxYtnLvmf7sevsinle1T5VNBi4vZrrOlfOq55cQ7YDaBk+nwqcJakDgKSDJSXbqNtd0jHh83OBt5IsuzWwJUwShwEj48QJsB7oIKmtpKbA1yqJp9r7JCmD4Ozlo+rsqHOxPFG4BiscRe3tsGH4EuBnwKuSFhBUwXROssjFwMXh6w8G/hJuZ1E1y34ZaBKuczMwq2Kckn5vZsXAr4DZwItU8uWexHYhaIPJN7M9lSx3rlJ+eaxz1VCXLrV1rrb5GYVzzrmE/IzCOedcQn5G4ZxzLiFPFM455xLyROGccy4hTxTOOecS8kThnHMuIU8UzjnnEvJE4ZxzLqH/B2SlrH3TNIpzAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Plotting\n", "\n", "# As copied from last week:\n", "def jackknife_batch_estimate(data,observable,k):\n", " '''Divide data into k batches and apply the function observable to each \n", " collection of all but one batches. Returns the mean and corrected \n", " standard error.'''\n", " batches = np.reshape(data,(k,-1))\n", " values = [observable(np.delete(batches,i,0).flatten()) for i in range(k)]\n", " return np.mean(values), np.sqrt(k-1)*np.std(values)\n", "\n", "with h5py.File('xy_data.hdf5','r') as f:\n", " k = 100\n", " mean_sq, std_sq = np.array([jackknife_batch_estimate(data, np.mean, k) for data in f[\"square-magn\"]]).T\n", " plt.figure()\n", " # I do plot a line to visualize the trend.\n", " plt.errorbar(temperatures, mean_sq, yerr=10*std_sq,\n", " fmt='.--', markersize=10, capsize=4)\n", " plt.xlabel(\"temperature $T$\")\n", " plt.ylabel(\"mean square magnetization $\\\\overline{m^2}$\")\n", " plt.title(\"Mean square magnetization for the $w \\\\times w$ XY model \\n\"\n", " \"(w = {}, #measurements = {})\"\n", " .format(width, measurements))\n", " plt.show()" ] }, { "cell_type": "markdown", "id": "326cc9c7", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "markdown", "checksum": "9f5b8ab1f6e5969e99680e87dd4ae2d4", "grade": false, "grade_id": "cell-9fc8c563421f36e4", "locked": true, "schema_version": 3, "solution": false, "task": false } }, "source": [ "__(e)__ Produce a single equilibrated state for each temperature and store them in the data set `\"states\"` of dimension $(11,25,25,2)$ in `xy_data.hdf5`. Then read them and produce a table of plots using the provided function `plot_xy`, which shows colors based on the angle of the spin, each with a title to indicate the temperature. Can you observe the [**Kosterlitz–Thouless transition** of the XY model](https://en.wikipedia.org/wiki/Classical_XY_model)? **(15 pts)**" ] }, { "cell_type": "code", "execution_count": 12, "id": "3d5a419b", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "fa2024da2b5516d1dc942aeabe1a6c1d", "grade": false, "grade_id": "cell-ef87967a80cd8dbe", "locked": false, "schema_version": 3, "solution": true, "task": false } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "T = 0.5\n", "\tequilibrating took 5.539074659347534 seconds\n", "\n", "T = 0.6\n", "\tequilibrating took 5.275792598724365 seconds\n", "\n", "T = 0.7\n", "\tequilibrating took 6.957661867141724 seconds\n", "\n", "T = 0.8\n", "\tequilibrating took 6.382277250289917 seconds\n", "\n", "T = 0.9\n", "\tequilibrating took 5.454374074935913 seconds\n", "\n", "T = 1.0\n", "\tequilibrating took 5.399614572525024 seconds\n", "\n", "T = 1.1\n", "\tequilibrating took 4.563677787780762 seconds\n", "\n", "T = 1.2000000000000002\n", "\tequilibrating took 5.712708473205566 seconds\n", "\n", "T = 1.3\n", "\tequilibrating took 5.183827877044678 seconds\n", "\n", "T = 1.4\n", "\tequilibrating took 6.481664180755615 seconds\n", "\n", "T = 1.5\n", "\tequilibrating took 6.137166976928711 seconds\n", "\n" ] } ], "source": [ "width = 25\n", "state = xy_aligned_init_config(width)\n", "equil_sweeps = 200\n", "\n", "with h5py.File(\"xy_data.hdf5\", \"a\") as f:\n", " if not \"states\" in f:\n", " states = np.zeros(temperatures.shape + state.shape)\n", " for idx, T in enumerate(temperatures):\n", " toc()\n", " print(\"T =\", T)\n", " beta = 1/T\n", " sweep_size = np.ceil(width**2/f[\"cluster-size\"][idx]).astype(int)\n", " \n", " # Equilibrate\n", " for j in range(equil_sweeps*sweep_size):\n", " xy_cluster_move(state,beta)\n", " print(\"\\tequilibrating... [{}/{}]\".format(j, equil_sweeps*sweep_size), end='\\r')\n", " print(\"\\tequilibrating took\", toc(), \" seconds\")\n", " states[idx] = state.copy()\n", " print()\n", " \n", " f.create_dataset(\"states\", data=states)" ] }, { "cell_type": "code", "execution_count": 13, "id": "2f322a5f", "metadata": { "deletable": false, "editable": false, "nbgrader": { "cell_type": "code", "checksum": "42eb725fdf7ff4b35fbc1be664deb579", "grade": true, "grade_id": "cell-8968d7e317e522f7", "locked": true, "points": 5, "schema_version": 3, "solution": false, "task": false } }, "outputs": [], "source": [ "with h5py.File('xy_data.hdf5','r') as f:\n", " assert f[\"states\"][()].shape == (11, 25, 25, 2)" ] }, { "cell_type": "code", "execution_count": 14, "id": "69d891e8", "metadata": { "deletable": false, "nbgrader": { "cell_type": "code", "checksum": "e9db0fa47469fbe95861db34bd66f432", "grade": true, "grade_id": "cell-cfcff21eb3fb562f", "locked": false, "points": 10, "schema_version": 3, "solution": true, "task": false } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsgAAAK0CAYAAAAeS8vMAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAACN60lEQVR4nO39d7hV1bn+/9+DrqAIiEAUUbEhGjBsNAaxF2KPvaPGXtEolohRNEbUIEbsnViILYoNjV2JUeAgMUgMGgViEBERBZU6fn9s8v3Nk89zD91brOf9ui6vc/KM/cwx91pjzjVYOm9SzlkAAAAAajX4pk8AAAAA+DZhgwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABABRtkAAAAoIINspFSmlP5Z3FK6dPK/z5gKc7TNKV0U0rpo5TSuymlUwo/u8WSc6meW9+ldS74fvk2ruElP982pXRHSunDlNKslNLtS+tc8P3zbVzHKaWz/uu8Pl1ybisurfPB98e3cQ0v+fkTUkpvLfn5MSmlTZfWuXwfNPqmT+DbKufc4j//f0rpbUmH55yf+AqmOlfSWpI6SWov6emU0ms555Hm5/+dc17lKzgPfM98i9fwfZJGL/n5TySt/xWcE74nvo3rOOd8oaQLK+d1rqTNcs7vfwXnhe+4b+MaTiltLOkiSZtJ+h9JR0v6Y0qpfc550Vdwbt85fIP8zTtY0vk551k554mSrpd0yDd7SkCdfOE1nFLaTlJHSaflnGfnnBfknMd9facKWPW6F6eUkqSDJN361Z4e8LnqsoZXkzQh5zw21/6VysMkrShppa/jRL8L2CAvJSmlM5b8K+PwH9PTStIPJI2vlMdL6lqYaqWU0vQl/1rkspRS86X3W+D/sq9pDf9Y0uuSbk0pzUwpjU4pbb40fw/83/Y13ov/o7ekdpLu/bLnDkhf2xp+VFLDlNLGKaWGkg6T9Iqkd5faL/IdxwZ5Kck5X5RzXsH9Y9r+869dZldqsyUtZ37+75K6S+ogaStJPSQN/vJnD3xta3gVSdtJelq1/wrwt5Ie4L/dxNLyNa3jqr6S7sk5z/kSpw38f76mNfyxav9Q94KkeZJ+JenIJd8mQ2yQv2n/uaEuX6ktr9qF+//IOb+bc34t57w45/yWpP6S9vyKzxEoqdMalvSppLdzzjcu+c8rhkuaKqnXV3iOwOep6zqWJKWUlpG0l/jPK/DNq+saPly13xp3ldRE0oGSHkop/eArO8PvGDbIS0nwVPP/+ifqyTnPkjRNUrdKuZukCV9w2iwpfbkzB2p9TWv4r6pdt8BX4mu+F+8u6QNJzyyNcwekr20Nd5P0YM75H0u+dBu5pP8nS/N3+S5jg7yU5JwvzDm3cP8UWodJOjul1CqltK6kIyTdEv1gqo15WzXV6qjaJ1AfWOq/DP5P+jrWsKQ/SmqVUuqbUmqYUtpT0sqSRi3VXwb/Z31N6/g/+koaxr+WxtL0Na3h0ZJ2TCmtsWRPsa2ktSX9ban+Mt9hbJC/eb+S9KakyZKelXRJNZJlyZ8Yey/5nz+S9KKkuZL+rNqFfOLXe7rA/+MLr+Gc8weSdpF0qmr/+7gzJO1KPBa+BepyL1ZKaWXVPgsy7Os+UcCoyxoeJmm4av/tx0eSfifpqJzz37/WM/4WS/zBFwAAAPj/4xtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACoa1eWH04oNslZrGI4118Kw/mnheB1NfXlTl6SZpr6o0DPf1Fcv9Lhk7Xlqb3vamr/CfEbx7/KIx1bUYtvhfh/fIcXvWm22i7OsqS8o9LixDoWedqbe4D3fM3aq3s85ty0cNtQkrZiX1Wrh2Gx3NTT1x2s6N66Xft/PTL10Mbr3b0X9s9Dl/gKlroWeT8Lq2Ab+qmxtFt4HpXfH/cXSpYAhtyAbF3rMpddwlm9Z1MQMLON7Gs6O66X3dJ7G1msNpyYrZjVbLRxbybzl7/nblrTyX83Aqr5n4gphuUm8fGrHTH1O4UVK8ceKco/CRHYm9+khSSuZujkBSbV/p0LgX4XXzV38a87zPa+aG1DhdetgXp5ppWvF/KrLFIKuPq3nGpak1GrFrB+sFg9ONU2dCgc0a7/BZN+y2PwlzI0Lf3/iAncttfY97r62VuFD2w1NL0wzx9zveuQPCl3x+pqn5rbD7Wj+XZjFfc6/VliTy5kNxcelzVt9/uq0f8bruE4bZK3WUBrTKhxaXzPC+muFw51t6tsXetzf51n4vNMUU7+tsKvewOxIJulQ27OPfhPWr7E3bGmxeQt2l9l5yf8+vkOK3zXpUbv1ktY3f+yI/xhQHjum0HOyqTcf6nvSCSrc9rxltZo215hwbIR7kdb0x1v1xbju1rbkr4nSHsYlwx+u/QpdT5v6k4We0WG1wbJ9bEef8O91ku7YozDN4XG54aa+ZdEPzYC780p2r7TCXb5l5spmYD3fs8LDcb30mTlJqV5rWM1Wk34cr+ED/hS3XHZY4Xi/XsUMXOp7anYLyx3G+hb3hcgL7rqT1Mjc2BeMKUxkN/alvw26n6mXtiMXxOXTrvYtE039obd8zxpmN1BYXIebl+d89+cASQ3Nnx/WdZt6SePqu4Yl6QerScPjdaxfmJ5rC8czt7vmP/ctH/84rq9oriNJmuaupX19TzJ/N90V5t4p+c/zy3yLXjR/lhrz2R2FrvgD7k1tZDvcjuacwiz9TL17YU32fCeuPzWwMJHfbnn7xOuY/8QCAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgIq6pViolaT40fSXNMpM8A97tGNNrsyPCqFtL2mzsN5Jz9meLU192YY+xqSBjfd50/ZcZWJRFmsN2+NC24bOf912DDJPabpkHMnnGSwoPKk6RXFEwzTzHkhSV/M+lJJJXCrHWqUn708ojBU0ViEtwkRFdDVJFZLUy9T7nlU4CRdjUXhgfoQ9B5+XsVA3hvVG+q3tOUEDwvrQwtPWO5v6lGt8zwvmeItKESm943LDcb7lh+YJ/G6FaW45wgwM8j0zTYzTzJ8WJioEFxTNkRo+Hw895HpKT31/+q+4vkzhBG+Oy/1d0oik4/5gBlyyg6QFLtPyrR6+6QZT73Wq73E3oWEu0kTS7nFaxS6F8I8R2bymj/jcqrVMy6RC7I19kN/F4UhawyQGjLvc9+ikwtjnWWaWtIGJk2m3d1xfvfDirh6/IB83OdC2vHxQXHfLTpJ6XRjHfUy8sI3t2cvUSxkpfY82AyYBSJLUI47fSE/sb1u6bhvXHyhM425rEzb3PT2fjeuDzLqTfFDLeuZ9k6Sh7rO3lHZk8A0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABUp5/yFf7htTcq7m7863T18flzheO7v7S49bDhCbcP6Lpphe26p4/ySNFQbmJHS861dwuohut923GJ+n96F38cFCpT+CvLOpn63GtqeZiZN5GOtY3s6KU7f6G87vMJD2tojaWzOuaaux2yQanITxYt4nvu1Sm+5i0h50Ld0Mm/t5MKD+QvGfhLWr9Oytscll4z309gHpOOMmlou12WCu4Qk6WBTX6/Q4xb+B77l9GPiugs6kKSZ5gH6kvRIXO9bSP+4Ralea7hRqskrmDU8Mw4ukQ6b7Q/4RMu4fkrhJFY19dKFvpl7XP163/PEuXH9gsI8Jqlh04d9ywvZxcQU8k5Wiq+9kf7WrT7uCftSeou7/xRu+OffFNcH+LAM9TFpGaUAlI3quYYlKdU0yhpj3qwZH8b1wo2o6c/i+rxCssJWJlnB3dYlaYipl3ruNb9m68L9wV1Kp79XmMi9Pru96ntWim/UWxXWsbvnTzavpyR1MO9DT9+iEbvF9ab3+555/4zrXVdPtmeC2U/wDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgIpGdfnh99VU12k1c6C3w3p/zbPH+x8TMXaCiReTpDtN/Nlo2yGtpKZhfaF+UOhyISulmeJsr1sKsWj3mli0OwuzvGLqpZgZN9az8FpfbI8Vn7Pkk4eO/si2qMF9cX3FQ3xPfbWWtLMZu8WdvInokeQzatwkkj4wcUQfjPU9C9UsrF/lW7STqd+tk23PEF0W1vcozNP3RDNQiLo78rS4fmP8a0qSFq1lBqb4HveWvq+NbM8Fd70c1guBX9rlj3F9inmvv4xFy0gz1zSDZ5h6axPlJklHxuX9CzFPd7h1v9ltvmnGgXG9bWvfs7+pT/AtLgLuhfwn3/PItnG9lJ3ZPS73KcVCuuMVXgIbY1i4x6xpYt5KRl4Y10e5aLovraHcL35O2ziSa+BuC+zR5u1mtjN9/Rk8ZT7mp+7oe1xk6j1xYqskad6MeFEcWlhgQ+o8IOkNUx9VyNw0m4M5d/mWyUebgc38TWNajk9uxM838RPdGH8ozptdyERtGX+4FNJALb5BBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKAi5Zy/+A/XLJ81ZmMzur2p31A4Ypew2kn3247DTf08k1QhSQv1UzPiHx0+X/FjwL9Sc9tzheaG9SM/sS1qvOzKZqSUG+Aek/bPaTbQw2G9mzlnSepo6u4BekkarPhR3gdN+ogktTP1vQrzTE4am3OuKfxIqFmqyatqTDjmUh8GFI7X+vd1PQPpsoPi+puFnis004y8a3s+M9fXjMKfi6809fa2w4dVjCr0zHMPLhee5m8TL2F9UEgZafRZXN96oe9x6zt+NWsdt7oZKF3Gl6Z6reEGqSY3MWvYmefOT9Iub8X1ESZ0QpJ/c+8u9PzC1A8u9Lhb3axCz7GmXgogeszU3QeO5GN+SnEnJnamjQ8GsjkH/QrTTDT1WwphJm1mx/XSSzBI9VvDktSmJuU+Zhm7a22ACukFi8zBdvUtTc09pZ9v0SCT9mE/yCQ1/Xlcn1e6jM3Wac9rfIsLlylxl/LZhR73ebDqjYWm9Uz9x4W4nEdMNMgOhRfB/EYNZKKGJC1Oc8N1zDfIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq6hTztnZNykNNLMn2utx0DSsccV9TL+XxuAydLX3LbJMp1fLZwjyrmvoU27G/tgjrpeiV/moY1hepme1ZrDgnrIN89ImLChpoO6S+irOhOsjkQkmaps3Celc9Z3smKH5/ttaLtufJesa8pfY1WQeZRWwS0za9zR9vP1M/rhSRZZb3nYWop311vhkpZJx1PiUsb/2mv+anm3rh6tLVZqkuOqPQ5GK1CjFvNtHO5RFK0s6mXjo3k1jZ8B7fsuj4uH7Ipb7nlnpGZC2fanJPE/P2fKO454BCrN0tu5mBPxbil1Yx8UulDEjzPhVf1z5xfa37fc+kvc3AH97xTeeZuM1SZJtbw08XevrF5XzMBNsyR13DegttVJgozsdbW6/YjmtNfctC0GXSBfWOeWtUk/IK5lY8y8Sp7l6IJb1HG8QD7/zV9hy/Sly/03ZIH5lrbEEpD8+slcv8R6lOdpF8bmsi2e3J8SbCT/IxfqVoz1tNfU6hx6Uizi98hOVfm4HCyTXYOz5gl8LamWD2E3yDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEBFnVIsSk+dztQ2pss/othB94f1aTqzcBbuEc7SM5cuR+KNQs8Npv6m7WhqUhcuKsziHoSeVeh5Se+Zkf1tTxs9EdZnqhS30MvUS7kcLoHkwUKPy9hwdUnprno9Pb1hqslPmQQAlxNxmXuaWNKt5ungvoWepqbnE31ke2Zr+bDeSsf5ibRHWJ2qrWyHC4oohT48ZUJibEiNVI7FcNxtobFvafxQXF9wsO85/6a4PiAXUh3uMqkOcZBIrXfql2KRGtbktGy8hnuax8hfvrpwQJMcstU+vuUpExjQeKLvOcokaQw90fd0/V1cn3C077G3jfhyqPWaqReesNeBn4TlpmlZ2+JSg0rTHOuu8eZX+aa5cRTLCRppW67QoLC+tk63PZPqmcQiSW1qUu5j9hMuPGR44Xjj9GHdT+In5kb950Lsy20mxuLAQirWTpuH5U4P+xZ3+yyFpJxs6tcXekxoj8lBqeWyyQq3VQ0wL/VGhYSNl/PYeGCFHr7pw/g+3cCkcknS4jSXFAsAAADg87BBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCiTjFvqaZZ1pjVzKgLqvlN4Ygut6kUI3avqT/gW141OVQbFKJcrL0KYx+Yuo+Ge0PvhHUXtyVJvzD1l7SZ7dlQz4X1cXq8MNNGph7HAUmSHrkwru+wRmGefqbuwmQkpbH1ihdqlmryqibmbSfTc0TheF1sbJKPHTxD3cP6RXrZ9nym+Fe9sPBnXBf5U4pse8XU5xZ6JpioJl1QaJpq6ncXes6Jy40L2U/rm0t8XCGGz7w95VRIl4joUg8l6a36RWQ1STV5RbOGXTTT4YXjuUjJ0q/r3r5JJv5Nko5/Na6XAiAnu3i6o83BJOk8cxK7FyYycYA608RMSbKfRQeZe6CkDrfF9X/rnsI88eIarg62w0V3tdI023OOOZ7/9JLu+BIxb2m9mqzfx+v4hR4prG8qlykp2Tf4oFN9y+/jqL7ib32SWV+lXDT3hhQusrVOi+uTdvQ9G5rYuC18i73+zy70bGjqH9lPEP+55y49SZrgIiAvf903PbFOXHf3aElaNl7HfIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQEWjuv14A0lNzJh7HPOEwvHiRzsb6EbbsVi/MyOFxIMNfmoGSvkEnU29te1YTveH9blqbnt6m/o0mScxJUmDTf0Q2zFO7kneBYV55sTlv/intLWDefRWOxfmedrUryj0/KQw5s1bVprUJR67zISQ3PeWP96UFn3Cek/z0knSaBP4MmgZlxoitZkR19/XBD+R7gur7TTAdoxtFtdX+MzPkraI67mwhJcz4QDtC2Enk3IcAbDg7ANtzzg3cKifR1NM/V+zfc81cSzGWnf5lkmFUyhpK+kYMzbE1CcWjufyh0pXbDc3TyFcwiXzvF1Yw+nirvHA0T4lRub61gY32ZYGG8SPyy/Wtn4erRmX3QsqadrecX3IXXvanqtM/R+FxIDRJjFgQ/n4lj1MvZxi8SUsO1XqcVI4dLFtKmWr9I3LbssiSdcsG9ePLryJ5uNqq8Laf2o/M1BYxpPys2bE5U5I4/RYWG+eXByM9IIJ5ZhS+H1cZthEGwHkd2jTTPKGJGmHF+P6WJ9msvU2cQLKkzO+eGLbf/ANMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACAijrGvJVcZ+pXFnrivJTF6uFbXj0srrc3dUlq284M3OB7bP7KTrbjY11uRnwE3aqK866mFeLkGmuXsL6glC+k/nF5Rlvf0tZkcfVcuTCPi+EpZe2MN/XphZ76afCJ1Nz8Wh9vHtfnFGLeskkkvGcf37O/iYArRcP5QMKpvknXh9WBOtx2dP+sQ1jfw8RTSZJGxeWNzOssSWe7+QvTrOqCx0o5VG7ZDSz0nGLqO/mILLWKy6Ursr7+naRzm8Zji5aL6yPiRM1a5iTvuPoT27JniuOxtixMEwdQSe1lotwkNZwW1xcV5tF8N9DLtizWufHAX1w8pqQfrxvXr3arW+qU4jVUWifPm/r4wtXibiUvycSaFexxVmGwkPj5+T6RTFSd+bRSE5m8S0n3uBvRjS6UUJLujcvP+fd9QxN/5j7FJPlt0B8u9T2zzTm0XKkwUXwFvlB6n8xL8FRhS/WUyY3cyKSyFd1ZGBtt4twK24kP3PaxrYsI9vgGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAICKOqVYNNanWlHxI5xbKn5s/w79unBEl+6wnm/ZYLY5tza2ZYH2MyP9/DzqaeqlNAb3nLb3kk4M633kn7h8QuYR9kLyha4xaRVHF6IGTMqIXiukjGzgnp8uRQ24x86bF3rqp7OkK8zYOc/G9dIT5s/2jeurloI+1ozL95j5JamDOYkH5vSxPTtrcFj/TC7VxT7QrOXusi02G+DJwt1lD1P/YKHv6WASAKat7ntc8EXxyWmX+OCCbSQb3rJ7oWVcYayopbRoCzPmLrNCQsobt8X1Ljf4xAOTF1DIiZBeNiEgL5cusAluwEZVFH7Xwv3xfpMYsFvpyfdr4/KnPu1ksknKOaLw/ujduDziHd+yszYN6/fpBduznamfX0hAGOCHvoDOku4JRzY16RalRCj72VxTuEE8at73QvDFOHO45QppRzrO1EspKeatanzme7almbl/flwI4LE3otKavC8u31NIsVjVbEHWMvcfSZq0jhm4yPeM01/jgbs28E06KazyDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgIo6xbwt0FqapivDsTv0kOkyUWGSpLNN/apCT3y8BRpa6HHRbDcUetzYA4WeVU3d5PRIklqF1ZH6le3opPPCen/5nJnjjj7IjBxre6SLw+o5G5xmOwZqm7DeWFNszwKbxfW47amv2fJhfC4EygX+SdKcz+L6RYUIpq3Mcjy+MM9AE7mzczFoad+w+lbhz8V3/MkMdPazjF8jri/o73vuNdFRv/AtmrZ3XG9TiKCbad68NoX0rpku3ugD33Pv/XF972a+R2btfJ4OH0qHm/lGm56R7hKTtKbL3CskM/773JfD+gnayDcNNPVSfJ67WPcuZFAdfYwZKGRDPWjqu5kcx1LThr4jvz45rN+gTrbHRepN9NNoZ/PZ2rHQ4+ZxqYfSl415WySfJTbetPzeH85tNcZ84ns+NVGG5/iWy8zH7MmlaE+3BSit/QVx+ahCFObQMWYg/iiv5bY6hSTFOy+N610L01w2I66fXEiNtR++xfjFeE9zepxELEkaZOp8gwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABARZ1SLLRoeWn2tvFYy7lh+TL5pw2f1hNhfYRWL5zEBaZ+d6HnjbC6v160He4B/GFa1/a4J8gLD4OqlUmkGKkTbY97mP4kNS3MZFIkXt3Et2ywZVi+UoXYAPPbLtAE3zJ7nbje8rbCPPXzvvxDu73qWJekk019l0LP/uYp6DvMpSVJHV26hFwEgdSgRfxM8dmFAAB7eRUef59mnkLuapIqJB+K4a4hSfYx+49KdzGTgjC40NLXPFneuBB64xJIFvlwAun1wljBtA7S+UeZwThkSGsV3ot/mqSNRW18zw9MWsW0lr7HxsQUEgPsG7WFuWdI0sFXx/X6BCqVEoh+eGBcLwUqbRtH2ByuV2zLa+oe1vtrv8JE8X2hEMSiPi6FoV+hyQcafa6VxzbViSn+rD+r2alhfdFlhQO6xJNtXD6HpKfNTbdwvzs5u4iidr5Ja5t6IRXr5Q3C8tBCWkbjH8f1DQvJFy5IY8Rffc9+7loq3FhPdtf5MybeQpKeaGsGChE7d8UfCBfIp5mQYgEAAAB8AWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKioW8xbw/FSy5XMYIuwerJ87k8HzQ7rDfSe7VlsA9iutz236jdh3SVaSVJ3bRPWrzTRdJJ0r6lP1ma2Zzk9Z0ZutT2rmvrZmmd7bjDzPLnB5rankV4K6zNViFhSHA1nY+YkabQ53jalUKL66SbZV3xrU99ek+3x7lOc4+USrSRpZ1O/4zXfc5lJPhwyo7vtucfEue3h0w0lE6O0fyFx73adH9a31gDb093dFtbz8zgL3AUhqYM571Kc3KbXxHV355GkM0192u2FpprCWMl8+cspvg1rUpwuVsv8Ysud5FumxQlUxVQ0PVgYM67cMa4fd0mhyby5jYf7lgUmCWzTt/y97oW94/rxhbhGd3LDCsGQg/VUWF9fd9qeLqY+3nbIJ0a6W/qX1E6LdYqJ3jp9w2WX3rm87t+QxrvG9QaFWLR5w03OWiEPdMNL/xnWSxGwLmnOXOKSpHtujOsvF1LRNnIRkIWe/cfG9ds395+VW5vPyqcecVFukvY39QmX+56JcXl7mTVVwDfIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSknPMX/+GadbLGXGdGp5v6GYUjxs/arqVHbMckxY8ON7YZEtJ9WhTW/TPA/oHU+wo97gHuCSYRo5bLOzARBJJuNa9PKfPBJXYUHrzVCG1iRt4odLln/a/yLb+Jn/DVsYVpVkhjc851zgFomWryTzQmHLvb9LTQRoUjPmDq5lFaSbO0VVgvrcdjdVxYf1NX2h43cpl5uSX5RXRCocclKri4Dqk2TiRyQ6HHnNtab/kW9yR4ad27wIdSCIO7Kz5d6DlZ9VvDqX1N1kHxGrbvn0sokKQd3gnLHyfztL6ktU19mkl2kCRNNfXSY/kDTf3H9xeazDt1kHnEv6QUM1BKBjGON6kqQwvBQC+8Htc3LbzWne6K6z19i+55OK67pAdJWrCwfmtYkpqmmtzB3Ivd7aHE3XHbFXpe6BHX1zIpDZLUytRfLrwf15v3w92fJH+5nFPocR+ZpdfzIVO/uZBxttiMNfrM97idzr6+RZeZW9Cm8S1Lkk+XKd6MZ8TrmG+QAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKhggwwAAABU1DHmrU3WmB3MaJzVs5yuscf72MaIlYJZ4qyOXTTbdjyipmH9Rs2zPX3VMqyvVZjHJQJNlnvNpN4msq2ULvSSTjUjpddtvKmfXeg5NKyeohdtx2BtYEZcCJ4kHW7qhVC99Hq94oVWSjV5HxMttLvpGVU4njvDaws9PTU5rF+lTrbHxfS46BxJWs/FuQ0rNP3Kre9CiODrq8f1WYV5Gpt6D3+/SMsdHdaziWqSpD7PxvXS6+be7ysKPY8Xxpyh9Y1561STdaaJeXPpkG8WDnh1/Jp3TfHrLUl/bxbXF21cmMddFKULzKWHlnL6nOsLY783635s/DkgSW3MOzfzk8I8JnvxstN8y2hTL91RXeTZe+Z9k6RFvc1AKdXzrfrHvKUNarJGmHW8+jFhecN0tT1eZ1MvLf1xJtquuCbNG9LhT76ln6mX7kPuE9vdOiXpIrNzGFXIJFzP1NvpfNvTTAPC+pq2w6cilpbXR5oW1o9SB9tznbnnd9rczzPZ3Iv5BhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACAikZ1+/H5kqaEI1vrubBeyki40qQhPFjo8c9iegs1PKwfqgMLXceG1UmFx1sbaKwZ8WftnuB0qQWStIsuDesjZNIEJEkXxeVfruNbfh2X3YPyknS9Xq1zz8kmVqGD3rE98bOtn6+d/BPF7hz7FzJFppr3thQU8aBJq+hf6HFPNHfRTNtz/BptwvprhXmeamKe2n/aP82vnqb+6xm+Z5+2cb27T04YZ96gKeapZUkaYurmapDkH2AvPaEdZ/iUe4YWxooWyy/WU++K6z/bu3DAOEVmwm6+o+FIM+CiBCTJpSS8V3j8X1Pj8qWH+ZaH4vKehXVyTwuzvq9+3fbMfNjcO/fy87gIghsKLduZ+h2Fj6+1bovrHT/zPTubt+Eq31Lv+7AkqekCaXVzj38iTqsYZ24bkrSFud3EuVdL3GrqjxV6PoyvsWk/9NdY3/hjsXhuO5l6KZXjBvN5dJwLl5LUyZzb24VPir6mXtq7uVvDJ4Vkld6fxWkVpYAkHRKXJy8s9JidMN8gAwAAABVskAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoqFPMW0PN0Qomzq2j6SnFmNyn5mF9sdawPfNMWMgIDS7MFFusEwujcbhN40Jg2QIb1OWDcla1x/JGauWwvrHesj0vufCqg0vxT63CainKxUVklWK1Bpk4t9NViKCTj18qaSSptRlz5+6icyQXkFV+jUab+nKF2KYOJrZpD8VRbpKP1SlFLx5yVlwf+c9Ck8ugM+tHks0JSvv4FvealmKP1jP1jppgey5T17DevTDPtXohrN+gTQtd9dPsX1Ln0+KxCW+a63lu4YDXmI+BQqZmIxMXtqjwZjT82PT8ZlvfdOZJcd1dxJL0QVy+pxARZn/XUgzmQFM/vjDPfnF53/t9i/tV+5h7giRNNPVS7KD7lNqj0FPvqEJJ+ldj6bT4s+z0OMlUgwpxZXNNzNvMeApJ0kiTiliKtuuZ4mtswI6+Z2sTpXZFYR53DnfrFduzmrlLuSg3yad0SsfZke1N/VrdY3se055h/aFC9OC15ta04N++x246GxbeIINvkAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoSDnnL/zDNRum/PKz8ViDqXH9ofiBcEnSz9Q0rC/U8oWzeNLU5xd63GP2PmvgSN0f1kupHJM0xoxsXeiKnwfNI8zjtZLSLoeF9Ua63fa01byw/r4a2p4Fam9GCo+3G3MKCRsuF+TOwvFOThqbc66p63m0TTV5d/M+XatbwvooHWKP51ISfIf/fUtpDM+b+qqb+J6uL8b111r4ntzDDLjHliXptbi8VuEp+ymFwznz4mVfjAwZYJ5sNyEakvzqvr7Qs6Wpb6LFtmcZNazXGu6RfpRfMqkZjfSXsH6CtrLHG5rNveaXhYSbg039KN9io47cRSRJ7Ux9iG9pY57Yn/kH3zPSpKdM9y261dSfKoQjbfi7uN64MI8JvtDJm/ueAeZz2nxMS5Ju1vlhPW0ywDe9mOq1hiUp1WyYNeaZePCXLcNy44v98VZcGNenXVI4CXMvbHOMb+li6i/EpyxJmjI7rrct3B9uMt9fXuCnsbfpUqqS+5QvcZ8hW/qQL11n6vsX5nn5ajNQmKehWa6L+hQmuj9ex3yDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACgok4xbxvWpPyUSTJrdUdcn2ficyRpoUkYW8HEv0nSQu1qRh7zE731YVjecPVkW8bpcjNSyJnRzmF1OV1jO1x8V7db/Cw/OCSuT9MGtqeP4uyjiX4aDTX1Kws9IxVn3Rwvk3Mj6bJP4vrJy/p5htYz5q0mdc5jNMiMugCdR+3xpqtDWG+/oz+HZCKYJheia94w9eG+Re+aeilOboJbQqWYt2Pjcp81fEtPU59VmMZd4f0LPS788d5Cj4v2alXoaW3q5qWRJPVR/SKyOqWafLqJKuzXKO5Z0K1wwIGmPtq37H9uXL+jEDuoP5s4uZMKcXJOKW/T3Li6FmLR/qqPwvoKhcjROSbqKvvkTGlkXO5QeN3+rQlhfWP5DNURpt5Ohcg2c1X+VLvYjpH1XMOSlFrXZG1rNhQue8wltkrqZO6rhduqZppYyw5jfc+0deL6gNd9z0DNNCOl/USci5i2Pd63mJt+fjVeQ5I03ayjUnCuW3mlnnk3xvXGhWjIjU10Xymazo0VdoiaZNYx3yADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUGGeeY69IdkMibv3j+tNCsdrZZ7o/UOfebZnf8VPQs/TYbZnOZNWMU5X257Gip8UbaBFtmeeeer0Y821PT9V87De+hDboummp41JqpCkI0x9X/lHro83v2spNWCuSaswD3xLkhrdH9c7mzX1ZfxbrXSO9gzHBqpXWB9vkiokqZt5Ojk928b25JXj+taFp6DjlVV+QvslTQ7re6mT7Xntrbi+/aV+npEfmPofCj03xPXz/+R73O86yrfY8A2XOiFJF5l6T5MmUOs+Uy/FLdRPa0n7mbEPzFPfAwo34rzjK2H9FHW3PZeZJ/mLN/xjTFpF6dF3EzdyiA/F0S3mvuFSXSRpXZNWsWWhZ8TvzYBZ25Ls77pHoeUHJjPAJcFI0mWm3kTn256ztLjO85iP8C9mjQXSH96Jx35ubpLX+sNNnhrXm27re3YxaRWr+hb1M/fpzoXrxb2Kw3S97ThYa4b1eX860vY00SQz4j+BXSrVQJloMkm3K77IzFsgSTrOLKSfmnuWJD2gjcL6j/Sy7XH3xtK+xb3ffIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCiTjFvq0j6rRlzUUulkKMpfeL6eoWenU19S91ke05SUzPio08W6FEzcoHtkSaa+ma2w8VNlWKoepnYuL8UIttGmci2BYWAockmtGW+XrQ9/Ux9d5c+I2m6iWUqRcPV1w80VwNtREycUfVgIebtSMVxbmsW8tcmmZf83Qt9jzvcv/WCb1KrsNql0HGEmei61X1PQ5MG1MNfknZ9D/MtmmaSn7qYpCipNpoyUoqGc/esnsUMMzdTqad+xktqZ+7cq5nIpPyiCz+SdO3wsDz4qPdty2UfxOu+TSGqcKZZW3uaqC1JusdEqd1SWihmnpmFG4qL7nJ3dEk6/mdxfWi+3zf9bLew3Lgwz7Sz4nqXwv3iIhP0trZOtj3n32i+LzvsWT9RnKD6xcxuLD1iLurxpmeddf3xfvr3sDyvcO8aYdZK1xm+5zVTf7K5O2lJc+OQMRflVqtffKjCPeVNEwnYxUR+Sn6NDzdRbpL0oKmX9i26Ki6P2Ny3HPVs/Fm9U2ke46f16OEbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAqUs75C/9wi5qU1x9Ttwm6Fcau/CSu913W9/Qy9fsK87jkjdKD0KNNvXmh50mTIrFAWxa6upv6vbbjXr0V1s8ozOJet1LKyGRtENY76VXbc5Gpdy7MM8TU71Bb35RmjM051xQOG1or1eTBihfxzjourP9IV9rjuadpzQO7kqS7Tf3QQs/ftTisN5CJLZA0zDzt3LEwz6qmvl7hrZBLKDmg0LNFXB5wmm9xIQSl58B3MK/bnwrfDezyr7i+1ip+ngdMvb1vUWuleq3h1KQmayVzIz7CNP1qW3/A3/0pLC93km/Z3tRLCSkXmeSNtf0S1namftmNhYncInY3dcnfCOMgmFrutn504bVeI36t14pv6ZKkSYfF9T0LKTGPmXppPdr53ysMrlS/NSxJadmarDXNOh5smrbxMSn7p3XC+h2bFE6in6kfX+jpHpfzn2YWmt6My403th0PLoj3ZjtrgO2ZpfPDeiuZqCFJ55i0ilKKlNvXnVPo2cjUS/eMG0x9fqFnd1M/3wSmSJLeidcx3yADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKgwwTuxxfLxGs+b+mulyU1a2JY++cRGsz2vHrZntMaG9TgEq6wU2DZLi8L6S5pbj5n82bm0op6Fo42v8yxSnhm/QY+38T07qmlY30DzbM/fTDzeIZphe27xp1DUSD7qaLSJc3NRU5JPh5pZiBaa+GJc71+Yp5kmxsdSV9tzuCabkeGFmeJcq64zjrQdE9zi2qIwjcllHFJo+ThOcVJ+fYLtST3i12e5+JZQa9e47KLcJGl1Eyd301fxHURrSfvVtalw5zoxztv8eJbP27znsrje0F/mOuYzP+a4mCetV2hyEV1nF3pGmXp82dVy2ZmPxFFukmwk4qSG7/ieg+J8qnviFM5a5oP6bJ+SptPdPeunhXm+jPbyN71t9gnLL+sue7iNWprI2kJeWdd4Gj3qW9TxTy6w7HDbs1j7hvUGC963PTtrKzNyrO3Z09Tbmyg3yUe2lfYGJrSueFlOMfVS7Oi4lnG9zWzfc/7qpl6IUnTBeXyDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEBFnVIsPlUrjdO24djW5unSdoXjnW3SKh4r9Dw3Iq7/ahf/WPoR2tuMuEeXpVMUP1XsUgsk/5z4fJnYAkm7m7HSE6TdTb10btNN3T3BWjrgY+5RWUkLtVpY7yL/+HRrk/6xqp+m3mZLetCMmWAFTTiscEAXHXKObznF1EtPTksPhdXHCikWXWzmzMGFeeLnk/sVOo7cMa73nON7jjN1f0VK15njrV94DQaY28KQwjzuIfH1LvAtaUb8XcMehdegvhpOl1a4NB6beaBp+s1Z/oDu4ftf3eR7+scXxaIOvmWoC0J51/fYyBm3tCVtZd7zp3Z71jddtXlcLyVfuBvUjwtREZ+aKJYP4qQKSdJFcXmrVXyLu9+f7sOepN1NPb71fHmzJN1rxib+ISxv9Ovr/PFcaI+74Ut61wSOlPYghys+h3Rud9vT59y4/qi9E0pv6qmw3klx6ozkQ1dKCVfuFjVQPh3IzXSfzdHwr+kZhZSRQ0xaxS3mcpWkjcxlPsDdGyXptrjMN8gAAABABRtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACpSzvmL/3BKM+TDVICvU6ecc9u6NrGG8S3CGsZ3Xb3WsMQ6xrdKuI7rtEEGAAAAvu/4TywAAACACjbIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCCDbKRUppT+WdxSunTyv8+YCnO0zSldFNK6aOU0rsppVMKP5tSSr9MKU1Z8vPDU0rLL61zwffL17iG904p/Tml9ElK6Zkv8PP7p5Qmp5TmppTuTym1Xlrngu+fb+M6Til1SCmNSCn9O6WUU0qrLa3zwPfPt3QN75hSeiGl9OGSvcf1KaXllta5fB+wQTZyzi3+84+kKZJ2rtRuX4pTnStpLUmdJG0pqX9KqY/52YMlHSSpl6QfSFpG0hVL8VzwPfI1ruEPJA2RdNHn/WBKqauka1W7jttJ+kTSVUvxXPA9821cx5IWSxopaY+lOD++p76la7ilpAtUu5foImkVSZcsxXP5zmv0TZ8AdLCkQ3POsyTNSildL+kQ1d58/9vOkm7MOU+VpJTSIElPpZSOyTl/8nWdMFCVc35CklJKh3+BHz9A0oM55+eW9AyQNDGltFzO+eOv8DSBorqs45zzdElXpZT4DMW3Rh3X8B2V//nJkr3HeV/VuX0X8Q3yUpJSOmPJv6oI/zE9rVT7p7fxlfJ4SV3dNEv+qf7vpqr9Bhr4UuqzhuuhqyrrPef8pqT5ktZeSsfH/3Ff0zoGvjLf0BreTNKEr+jY30n86XcpyTlfpC/2rzWqWiz5v7MrtdmS3H8H9Khq/xOMuyTNknT6kvqydZwX+H/Ucw3XVQv97/Uuldc8UCdf0zoGvjJf9xpOKW0rqa+kjb+uOb8L+Ab5mzVnyf+tPmi3vCT3r5pvknSnpGdU+ye9p5fU//VVnBzwFZij/73epfKaBwB8RVJKP5Z0h6Q9c87/+KbP59uEDfJSklI667+eVP1f/0Q9S/6742mSulXK3WT+NUfOeXHO+Vc559Vyzqss+bl3lvwDfCn1WcP1MEGV9Z5SWkO1/5kQN2YsFV/TOga+Ml/XGk4pbShphKTDcs5PLq3jfl+wQV5Kcs4XVp9U/e9/Cq3DJJ2dUmqVUlpX0hGSbol+MKXUOqXUOdVaT9JgSQNzzouX+i+E/3Pqu4ZTSg1TSs1U+59sNUgpNUspNTY/fruknVNKvVNKzSUNlHQfD+hhafma1rGW/GzTJf+z6ZL/DXxpX8caTimtr9owgBNyzg9+Nb/Jdxsb5G/eryS9KWmypGclXZJz/v8SLJb8ibH3kv+5oqRHJM1V7X+PfFPO+bqv+XyB/3aQpE8lXS2p95L///r/DFbXcM55gqSjVbtRfk+1/+3xsV/3CQOBL7yOl/hU////TO7vS/438E2qyxr+haS2km6sfDvNQ3oVKef8TZ8DAAAA8K3BN8gAAABABRtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVNTpr5peMTXPq6m1GW1q6g0LR/zE1P/7L9qqcn9pXPNCz4ywOud//f0c/9vrTeL6uvP9LC7jp/TbOKW/Gq+DqS9jX09JMie+9pu+5R/x6zOusGw2bDY2HvhsRT/PKu/H9Q99y9g5ej/n3Nb/RKzRiik3Xi0ec39fd/vC8V4z18OPFn9ge9I4cwKN3rY9H5mluvwU2yKtHJf/Xbgk3d/33KKQtP0388fs0hVpLi9N10qFrmXC6nKabDs+nt4jHmg3qTDPKnH5zXh+SVJndx0t8D1j59ZrDa+Y2uTV1NGMut+rsz+NxubdWOTPoYFZDxv+0Fz/kvTXlnF9+f/+m8crzMfKZzPM+yppQo+JYb1T4f4419Tf149sT7OxKT639WyL/5Awt0BJamMSwtsUpllO75oR/3dMjG21ZlhfcZaf532NrdcalqTUcMWsRqvFgxvEn9nSvwtHjD8Zm2iq7Zhv7inSuoV5Xjd1c8OV5N7gHmP9p8vYteJPpB8s76+xf6tdWF9G023Pp+YzrIX8Z9gcrRrWe3xS+EAyf53Z22v5lpmKrzF94K9LtX41rs/YwPdMiddxnWLealLHPEanmFF3Ay79HRnjTX27Qs8vTL1noef6sDrK3kSkTVeP6y+/5Wdxv832vsVyr7IknW3q3fRKocss3Cd39S1bxxf18oVb80ddzIKeeLif57c3xPVCdHl6RmNzzjX+J2LL1KTceUw85v641L9wvO46MKwv+OQ229Oo+c3xQLtDbc/jZqlud5xt0eLfxPXzCn9i29LUNy382WtN8yeL0hUZ316lwTqx0BW/Q1vr57bjyUvN/e3UHQvzXBSXf1a4wf7xZ2bAfzApvVivNVyTuucxcn/p1U9N/W5/Git3iuuFfWsLs8/66N/m+pekH+wc17crXOjxnk0Tr/KfW+uZl/R6+Y3FaFO/rvAHnK4p/qJgwl9ti/+QGOZbDvlTXD+4MM2WGmRGRtmetPeIsH7kXX6e65TqtYYlKTWtyVrZ3Iz/eY3pOqdwxPiTsZNOsh2T5a7pvxTm2czUzQ1XkhR/xuV0pu1If+oe1s/fxl9jA3RqWN9Ql9qeceYzrLf8Z9jzujqs53HH2B6dEZcPfcy33Gb+hLzwts9804FrxPVr/ul7jonXMf+JBQAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVdUqx6Jhq8i8UP3XaTyNNk3uqWj744hn3LL0kHWvqhUchrVKAl3vc2GVISP4JW5/gMFF7hvUuOsJP86RJfdjaRT9J0sCwOlWH2A6X8XFUYZb/2Td+wnb4cL/O9svx4/I5+bC7pPXr9fR0Wrcm67p4DXfaLD73CYXjNU9mET9XiM8zftHbj7lLpfD8v33+/kmZuC1Juyh+Lx4opaKZX/XxPr5le+0W1s/R/bZnoEm4uFK/sz3HaYewngc8Ynsanh8H1HWxYWCyoWulO8wtqX5JLA1TTW5u7sMfaSvTVcoD7BVW0+a3+hYT+bKRfyv0kotsMvPXMnEZh7v7s6QbTG7QA4XPCHMKUwtZalun+J62u2/RRRoQ1hfKpxk0bhnHxJxfSBk52yYa+Zi3tM6m8YBLdpWkF79EikW3mqxHTIrFyvEh1yokkbxtImUXlD5LTSxan3Su7XjUfCKk/JLtecEk7Wyqo/2pqUtYPaWQyjFYvw7ruckvbU+a79aeWQ+S9Fx8X73TfIZKfj/Rz4TbSJJuicvjC9dld5dMMrsQL7MCKRYAAADA52KDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgIr4L5M32ulT9bPP9F8Vl/crHNA9WL37075najx232ifkrCdqbfQ+bZnvOK/l76b7rE9PuEifhpVku409e11ve3pdYhJsdAdtsfpWMho6KgLwvr/FJ7Nn2jSKtzvKUkbpThV4eJC2kK9zZYNPGmxWVwfVTjc5gvjCIcb4weqJUlHzI/r+xbmGW3qI4tPQbsn/f2T7C30alg/dC0/y3wz5vJmasWPLg/U9GJX5ATFqROSdKvitIoFcZiAJGnRfnFaxQmFRdzf1FfVH3yT9imMeeuqdKWbeAl94A+43bCw/O5jcV2S2iWTFGFSAST5eIf1ClfYBXfH9RuG+B6X+VL4WLlvl7i+exz+I0n6xxXxE/vp+D/bnkHrxp85T//dP/2fXVJOMVGpian7hKiur8f1vzfzsywqnMHnWbnxWJ24cvx7n67fh/VJ8skqG+uJsP5S8XU6I6yOzM/ajuHawozsZnvM7ki6/2rb02G3+LUZrANtz8aK0yrSI4W0svNM/Q3fot//LCzve1GhZ6O4nB7cu9Dk7jM+xSbfYqJ0uvhrzI3wDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgIo6xbxJCyUbw3RwXL7Yx8ro5lPius1EkTT6+bC8u14oNLmYHJczJ3XT6WF9vAbZnuNN/YVNbIsOeTGub1roaTgujmwptOhZfRTWG+ifha6Jpu7fU5OgVggW87P4MLn667zyWF386zjUpaPp6RkvOUnSrr3j+smFc1irafz+HZt9DM3pNvKuFA7ncnr8undX95O6xPZ01WlhfVMdZnskE1W4nY/I0uO3heX+MrE+8vGCfZv4eKU2dx4T1mfK5ABKGqo1w/qCT3yUW2M7UraMZqibrjOj74bVvbKPRepp6vFvVGsX3RfWG+kfvmmOiY27wMdgSheber9Cz71x2f2iknZ38VS/Klz8JgYznzDWt/ztJ3H9pcI0t8dRknr6UN9zQ3xfuED725YxWhzWm332C9vj71ifb54k85vpAx0U1ltrB3u8l/RePPB6W38SLeI4OZP+Jkkab1r21P22p10yMWt5oe2Z9ru4p/eJ/lV/Xr+OB7a5y/boThOz9nsfdedi1tIZ5sWRpNkmnu4vvuWQH8e/6y263DeZ7cmsTr7F4RtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrqmGLRVC4R4hRtFdYHFI7Wyj2E29ikW0jyTzX/1nZMVIew3kVNCvPET/o/WOh4XtPC+hkvxvNL0iDzAOmVhYdO3Svg0iAk6RotH9Y7q7vt2V6PmhH37LHUz7wGq5r3QJJ216SwPkdr2Z6+dqRskXyiRs8JZiB+YF+SdJ1JsWjnjiXp7fPiJ3P9qyqdrv5mxL/rjfV0WB+tRbYn7pCe1HDbM0HbmBGXHiN11U3xsR7f0fa4d+43hQSAvTaO62spTqqQpHGm3knP2Z6ZJv+j8bL+qXvpkcJYQfMpUvej4rFRcfnuWwrHc7fbWbsXmrYz9W6+5XEztoVP2Ej7mKf/C0lH+VVzJZ1RuJCnHh6WR2lT29KrmM1jrG8esX/fXXmSZsXlBVf6liYrxPkt11/qe9Y135e9vd0Q3/S4H/o8nca20rVp+3Bs4Vxzv1nW3yWX00ph/dh1zRqSNGj1uN70Lduii34fJz/ldKHtSTYsqvCp3WuDsPx8KR3o0rPiemvfooFuoIXvOcb8Qv5S1pUmseu4/LrtueUt896tHqecSNJDJq3Cfxp5fIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCijjFv78nl67TXoLAeh6XVauVC0xYcW+jqFVabNPKRYPMXXhfW35TPyemsuKcQxqMBm8dRZiOf9T0X3XWHGfERdMdqzbD+dCGybctucbTYDeN9BM7TNprNR7ZNNfUD9ZHt+cwEsCzX1rZIMwpjBc0kdTFjv+oa17td5o9nQpuUuh7tm7o+Zgb28D023M+dgbTAZPv01AeFWeIIuAc11va4q+jRgU/Ynnm/jOubNPTRZ+MUxx791ES5ST5EqRTw2GR+XH8/TiOUJA3rGkcVNZGPMNqvcA5FzeQzi26Oy1uv7a/zJzvG9wbNiu+BkqT+K8b1i/16VCuTAfXMtbYlP/OKGXnXz6N74/JUm2clKY6027SQOnjIQ3Gm3s1NTNaeJC14Pq5fbPIiJaldXG5ciHnLD8Vxm5JZ3JKelsnHevxsP5EuKIx9nkVyGXbnLxt35Dv89bRgz7jeuItZ35IGLY6vi/mFhLN3k8mGK2nvBm6wLVv3+F1Yb5z8tTwym/v0cz1szy4rx6/PiEf8PB2uieuFl03H2Y/EA3zT6vFnVQO9Z1teMfXt/dK3+AYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgIo6pli0knvS3j27vGrxeC5PoPAktN4Mq3MWlh5RXC+slp5kH68jw/rThcc0jzdpFfFzuv/hEg26F3riZ/PXLPWMfzQsH64RtmWWdgnrrTTZz6NTTP0+29FMO4f1j2f4c1uucAYln8knG5xnfq3h5uFuSWqtE8N6J8VPIEv+mnhev/ET2dgC/zT/Wro/rB9cmOVkbWZGSgkb8bU365zTbIfJGdD/mIf8JWlW71fD+l6+RZPN7/OYnrM9r5kbQ7chfp6J18f1OHPnS5q5jDRs7XDoglvjpIgnS5kZU+P77fzCVdbE/L42ckGSZsVJEaVPidHmntbTpCbVcmkVpTt+fGPv+rDvuDnF81ylAbbn2OYmUaG5n0cmFGPqQz5l4AxTd3cRSbpdJ5uRUqpU/VMsxnZZS+n3I+NBc+s479Jb7PGafHZIPPCaiVwozJMv2dy29Db3jkn6l59H5nhP+IirLtvEnyFDPylM466lzXwci3vXR+/g0z+mHRivvctusy26yrwNk6726TIdFJ/3NO1te1rrrrDepMmB/uQUnzjfIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqEg5+6iY/1aTVshjbAxUT1MvxP7YeJ84Yq2WCyWpT6DSy3YkrXN8WH/tdX+0c0y9FJC1vam3MtF0kqRecVDZwlFzbUuj7eIcoVGP+/d/uKmXwsiamMivJo2b+qYJpr62DyRLGjY251xTOJXQ6jUpnzcmHnPv392F47mQoxH6Qx3O6j/cNSS5+Lx83f2+pZupX1WYZljctFmO48Mk6XltY0YeKEz0C1P38+SRL4b16X38LG1MJFKjwnU8fcO43l6/tj15xC/D+lFxUqIk6bqkeq3hmrRsHqM45k07mdfvIb+KT9CeYd2HL0mt7EW7wDcdbl7YNwoTPTPYDBSa+pkFPiSOuqzl4jZL0XAm/uzw1XyLSdsbtpVviYNNpfMKkYja7Ekz4GMhXbzq4kIGXUO1rNcalqRUU5P1krkZv2aaHvLHa3NWXJ+Zb/JN5x0W1y/zLTrO1H+9QqHpCFP3exBpI1Pv51tqVo7rYwpRd25NrHCub/kwjlLTcz5+zS6jHgf5ntm/D8tNW/oIunm60YwU1n76ZbiO+QYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgIo6plisncfod/Fgt5/G9fHv2+MtVuOwfqGWtz1xfoN/YF+S+it+lDy1HGF78mw35pMG1leHsN7PdkhdTL2XZtqex9QmrBdyL3S9qQ+Ufw38E9ztCz2tTL0UnXCxqbsn2KWkU+r19PSyNSmvax6cjp/hLq8tl0LifiNJGmxSHx7UE7bHPdR9ZuEJ80Xj4lSTA0yYgCT1NfXtXy3cJza4LSw31iG2ZYHODutX6jzbc6w7uf62RVq/Y1gelafaFrdS75B5QlySS9HZU+Zpb0n31DPFokePlP8SB3qoSZOWYT0fMdser/cN8Xsb59HUcs/kH15KBsqnhfW883N+oofc54f5vJFk71vXjvItR+1uBlxkgTQqbx3WC0EL+s1kM1C6YUwx9aG+Zf1Ocf1vT/meXU2SxojSuk/vfIkUi+5ZY/4Ujr2c2ob1Xo388eYvjD8z03vx56Ukqe3CsNxwGT/Rok8vNSMuCUXS6/HvWQjtkfZ2yRNDbMta6e9hfdKphXkuiVMkzlB8X5ekixaZz4OGb/l5frd6XD9xW9vSwXwmTvu08Hm0zKtx/Zcb+J4LEykWAAAAwOdhgwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABARZ1i3nqkH+WX9EI4NknLhvVyIJiLBHrQ9szS+WH9zsI8x5pzHqJNbU8/TQrr52gt2xMHPUnb2/AsyQfXzfEtD5ieXccV5onj1+6TyQOStLt5fxZqfdsz0KyDs7TY9rxl/pwWB3TVWk5xLMvnSWvUZA2Mc97yNSluGlg4oEtZKyT+zDonrrcqzLO26YnD0urPreF3Cz2b3hbfQ7oeaF5PSRP057B+q35iew52qYMLHrU9LrQt/bMQb7i6iWb7+d625Ywb49/1Il3u50kn1WsN13RMeczJZvAXLnjwosIR4yttqomTlKRVTTLTSJNmJUnbP2nWw9Z3257ltWdYd2lpkrScPgnrM829SZLadY7PLb1xie15UHFs3U6lvE13wX5Q6HG5g4X8ydHmA/GCwjR3mPrThZ6d6xlVKEmpY03WSfG9OJm0x8VzbikcMV77qWUcvypJdgvgcjUl6Z/uHlWIK7zNxMa5e5okjTb1Sx7xPU/sEJYbb+Nj6xaoezywj8lDlfxn4jrxtVcrjh3toJVsxzTtZkbcfU6S4g/L49N7tmOo2U/wDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFXVKsWidavLWip9sdGkVlxTSC35m9ueDC+fQRdPMiHvkU5LuDau9davtcNkSnQuzuAfL97XPB0vqdUBcL4RYyD20/4N2haYJpu5ft0PVJ6zfqI9szydaPqyPtx1SL91SGI0lHVq/FIuaFbPG7BqO7ZJuDOsPPOfTGE7oHddLDycPtskG/hnzPpoR1u8pzNNCPcL6KRpre9z7tG9hniPM08mNzVqQpAXayIz4x/m31uth/cm7CvewFqb+hm/RifFrXT63dcO6CzmRpBH1TACoSY3zGJMwMcvkjbTSoMIR4ziE+8z1L0nzTX1fvVKYJ74/DZdPGdjvrLg+4EI/i3u+vZd8csni2fE9YbZfwmrVOq5PLCRSuM/J632L+rtf6HGf85NynNixi7mGJH/tv93dtiiN/xIpFivUZG1qkhIeMm/8O4U3/k1T32y277mmZVjOxwy1Lanl8fHAEX4a7RGX8yY+56u3uZZfiE9ZkpQWmXkKPTLruNOrvmXVHH8mPq8TbU8jXRvWF2prP5G2jMt/OdW3uE2a29RJ0uakWAAAAACfiw0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABWN6vLDa2i67tZl4VhqGYecDf3wJX/A2ZuE5UdXcFFukvTTsLp2IV7ot9olrJeiuG439e01oNAVh2Ftrf1tx5OjpsYDV5zhp3En3nO67xkdR4hdZd5PSbpZ+8UDrYbbnhazDg/rvdTXn5t6hdXxWqvQc2hhrKShXPbXiPxIWE8HFWLEel8alq/XabZlkE4K66fLZ/GYlCDtaTsk6YGwOlgb2o49TZzckOI8j4fVBZ8u9C3mgBue6SP1nqwx70O85Grtbc7hudKtr2tYbWpeG0l6UjuYkVLA4TuFsYK1Fkq/i6/1Vj+9zjSVAirjqMfdCxGQD7r7YDe/tuaPnxfX/YkpX3i6GZlS6HLn7a4iqUHLLmG9lYnAk6RddWdYf2APv4a1Xlzu7zukWXE55T/YliP1k7AeB4fVmqwN4oGLCnlf8cfxF7OMXMKgdJqJc+tXOJ5bSOf5+2qHc+P61jJRbpKNRVNP39L0x25N+OvyhRzH052e/O8zxNRL11h7cxuaHKeESpLeTv+IBw5e2zeZyy/t4nIMJT1hYuP8lkp6z8T6vVnKuovxDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFSnnwhP6/6UmLZvHyDyleLN5WvvpwgGH/c0MzC00zQmrs7SV7Wi1b/wE6a7D/e/ungVvpzitQ5JSW5MIMdS2SPHD08oNCk9Crx8/Lvt49k+db5feNyP3+nncm7edT7HQKFOfe3Bhnp1NPX6vJSnp0LE555rCQeO+mtWyxpxjRt0TxYUkgmPMU7buV5KkHdYwA3EKiiTpU/NU9zI/sy0ddH9YP9bPYvI9pH4jfE9631xHh8VpHZKkuy6P63vv6HtkEl+Kj7bHCSnlFIQ3w2oDnWo7Ft9m1uqBhSWaxtZrDXerSfnxMfFYuy3NVGsU7vPmEfd820jbktbpEw8UooEGmTCEUrLCYE02I4UUm33NfWv42YWZ3HP+E23HfYovCnd3Kbm2MDbM1NsVes4/K673MbcRSXpUL8QDv+lte9KZqtcalqRU0yNrjEu5MvuMR/7pD+gu6aMLKRxvmeSO1Qv3LhubE983asWf2RtrFdvx0l3mmjXXuCTllSaF9dSjkAg12NQ3e9HPs16ckqKJ/oMv5U/NyD22p6tJ7NjOdkiX3RjXGx3WzPYsTPPCdcw3yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKuoU87Z2qslXKM4X2l6D4qaDz7DHW3FYPPf7uqNwFt3i8r7r+xaT0tVkkv/d91oY12/UYtuzhvnzxrT8uu3ZM60T1u9Wd9vjY8duL/QMictn+2g4XXC36dnLtnx2waKw/lbhz2Jd9HJYf1Mb2Z41leoZ8+ajhXqrcVh//hi/TtJtcT1/vJLt6aMZYX2kFtge6ai4PNZk2khSj/jlOUVjbctgHW1GHvTzjP2Xmd9H0Om2P4blBge6oDmps4l/LCSLaYIONCMuj1CSngyrG8vF80kvKY77+zSZ6EdJy6hhvdbwMjUpdzYxb3eneK120RGFIw6Myzv9wHY89lA8z/aKY6ZquWg2k3UpSZpu6jcXesxa7VaIhnPpVB8Uprkqvj9O1562pf1hcb33jT7Wc455T/+nWyEKdLz5hW4vXMcH9A/Lz+SLbcuW6UvEvLWuydrWLOS+pmmHQgzkzx+O6zcWYt5+GMe8bfxX/9q+tE/8frx2l5+mq7mtlT4npHFx+byVfcuVpv7eO77ndXO8dX7ue8znTq6ZaVtS/p+w3kD+c2Kx+6w6KN43SfLRwi67V5J2jPcTfIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQEWdUixWSDV5c5Nicbbp6al7CkfcMqzepza2Y3dNMCOz/DRdesf1ifFTu7X6hdUH1cF27PIHM9DZz3Kyef53U9+i3TUyrB+qPrbnZu0RDxx8n59o2D/MQOFpcLUOqw1ajLAd/efEdfcgsyStV88UiwapJjcxa3jehabpzLf8Ae9aPa73KpzEypfG9ZtO9T2HzQ7L52sF2+LSHfoXnua94ci43t63aGfFiRRtCk8nm7e8uO53MvWT1bLQ5RRW1+zL43rLnxSOd3hcnmFiCyRppfqt4VVSTT7RrOEhpufawvHiu7DUQv6alaaE1d463nY8b+/dw/003S4IyxeM959bZ2vVsL5Xnmp77ja3R923r+3xXPKGNEtPhfVWeqFwPJOE0q1w7x7/gBnwKRaj8w1hvWfh4k/Tv0SKxQ9qso4yKRbuYn+zcEC3kH0Ih008uNUH/ejwRnF9QeGef9mzcf3kQsKVdoqTGk43YR2Sv6utd5bv0a/jJI3GhQiX+UfEaVWjb/DX5UbbxvWt/uRP7SkTWrLg4U9sT+O2y8YDz/t5tC4pFgAAAMDnYoMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQEWdYt5Ss5qs1eJYlqdMWsmWOt8eb7oGhPU42KfWfFMfVehx/lIYc4FAJgWr2OPitiSpham3qkePq0s2hEpbbeJ7Or0Y199Wd9szUa+E9V39NOpm6vesXGh6p34RWalm7awxV8aDj5gcmh18pIw00dS7FHrWDqt99I7tcOlGTxZmudXUC+FQetTU22uzQtfcsHqlfFbScephRjoW5hlv6j6HqrfiRfz8I4X7nomFetdENUlSe7dW/3Wbb0oH1WsN16QWeYy7Bv9i7oQ/9hlUN5iIsZ0L59BG8TUxQiZiSVJzU99ePncw7RbfcfP9k22Pj6F0OWCSy2U8QVvZjt1NPQ5Lq3W7CeJKO7qrVcoPH2dGSvcYd030LPSMNnWThSYp6ar6x7x1qsk6M95PbHpM3PNCIeLseBMJ1q5wDueYD839XA5lwR2FKLVBJkK0RU6258gU31dvKVxjTil1dL28eVg/Q8/Znt+kOEZwuHaxPW7vtr+5l0jS8+Z3dceSpO1NlGLKf/NN6SRi3gAAAIDPwwYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFY3q9NPLykYOzDIpFj8ySRWSNK5lXH95tj8F9zztG75FT2paWE9tO/gm83RrYxdVIWmBe6zZPwSsO2fE9f0LkRRHmCds+/kWmwwy70X/POgwm7/R3/a4jrNth9T3T3F9KxMqIck8p/oFzF1e+kt84P3NU9CPFZ4a7mce4B6gDQon0Tqs7lxIsXDL7sHCLMNMfZKOtj176Roz4p6kl9xNoZ+6Fnp2MnW/Hl/Q/WF9U31ge951Ay4QQ1Ibk1ZRTMoZaOrrHlhoOqh0RGue1tWbJnni0I3jp+Kfe9Kf/eFbjwzrZ6iP7blQC8P67nrZ9tjV+pzPLcqbxckKU3W97Zli7g69NNOfmkm+uMJm7HgT1caODXfZMr/xx0sd49SdXdylKukiU++iy3yT+RR9Uyb1R5J0VWGsrHXbsepzdLxe77gvTpk539yjJeksfRTWz9Xytqev+Sy9tnAfmmU+5UabpApJ6q9V44FUyti4KaweXMjSOsGc2wWFtAzpkrD6YCHFomeOc6n2TW/7aXquFtdHu5un9JjZP15USL7xn3xtCz0xvkEGAAAAKtggAwAAABVskAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFCRco7jVMIf7lqTNXxMPLiByZz5iY+UchlMbd7yLSYsRVv4Fg3W5LD+oDrZnummXkhf037uVy1kcfUxyV6FZDjNM/FnI0xcmiQdYOqTCvO4wKb2hZ6e2i+sj9KdtudiUx9hYgAlSbPT2JxNxlpBx5qUf2GW8OgUXwt35J8VjjgxrL4sk3so/96eU5hlnomo6S2TEyjpWlPvUkj8aZIbhvXxWmR7hpv6vX4a3W3q62lv29NUd4X1efJRavvrtrB+l5ranoWXfhYPnPqI7dkw7RDWx31iW6Rl67eGm9SkvKJZw3uYnivSYHu8Rs1ODusLP4vjMSVpjuKIzBZa0/Zo3zfD8vrD/WeQC/Dr7GexkYi9Cj2365V44LwNbc+sX8Xn3UpxbF6tOFcsnbin7eh6eXzBHmzuV5LUX7uYkd/aHm2xdlh+8xk/z5qq3xqWpDVqUr7QrOODzPW5IPlr8AJtFdYLyazm3ZD2LfT8rFlcP9fcNiRpgP0McWcg7Z/iD/Tb1d1PpJ5h9eJCLOLpeZ+wfqu530rSFFPfubAmu9mIwca2R+5+0vynvmWu+eR7+yjbklZTuI75BhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACAijqlWHRINbmv4sdOXbKBe9pR8k+/l5423t7UD9cE23Ofuob1PXbz8xxyf1z3z5y6PANpQiHII1/zSljfuPCkahNTLzw/rltWNgN9fU/6XVyfUHgRTjH1ob5F7nnUSS8Wmjap39PTqeYHWWOODMfa6LywPrOQktDI5DEseGqePwnzQPEoH/ShTXViWO8j8yZJGqnNwnp+6jk/0a2mHj8cLUna6/i4/phv0ce60YzML3TFKQhSu0LPBWH1Bc22HZv+PL4n5psG2J6U344HZv/e9miF+q3hmtQ1j3FPmN++flx/o3DAX7k77sG2Zbria6h9XqEwUbe4/MtnfYtZj8uZ9B9JetLUexbvQu6TJU52kGSX3cTp/jPVfUaUUouu0OlhPZ06yDf1N/W2PplAq8QJMvmdOAVKkpJWq3eKRYOalBuNiVNzVjOpOZMOKuxXRpn6eoWT2N3Uryr0uM3Oa77l3n/GSSS731KY59ArwvLTMjdc+UQolxokSS06mlijUlTMQFMvRTE9c2xY/nH2L/ZOpn528h9IZ+bRYd3cfSRJ+yVSLAAAAIDPxQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgIpGdfnhFvIRbL8w9bcLM2y9MK7fXYzjmWvqLq/Fx3tsdb+fxSWctPYt+mMzM1DI8El7d48HChlZA0xC1UO+xWcfxSlYkqRP5iwO65sX/lz18gZxvfOrPobvWhPD12QTf26b+qHP8aGk+8KRmbokrJ+i0+zRBpsIuB9tdZvtabdVXB+ptrZHejesLih05KEmzm3LQpNb+IXsRRcveI962J4++nlYn+Wn0Wg1D+vt7D1BmqbHw/qmNhBJ0gemfoW/WLJZ3qlr6Y5RP/O1jKaaa6bjAS5/qRSFF8ciTTRRbpJNKlRO/i6ULjdXbTy9JGmEiXOb6lts2tdGF/p4LBfRdeVtPlas5/S4fq+fRRcpjmbbI79ve4Zq33jgpsJEbg23HWJbNvzXPmH9hPTFo2DrooGk5V2cm7l3XPl7E0kmqb05zz1cLqrk4yvjRLJaZ8TlDjN8y+6D4/qbLhdV0gGHnBDW+yuuS9JI/Susn5Nczqs0eOrz8UCL3rZnbTP0j8KFOW9hHOf2kt6zPS+ZENhNTZSbJF2kX4f15fRLf3IG3yADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUJFy/uJPqKaONVknjYkH3zRNpceNH4jL7xaSL8wDpGpVmGawRoT1q7SL7TnOpDEs96qf5+PNzcBOvudeE5BQehL6RsXpEgsLf95ZziVCmAdYJSk3Ot+MmEe+JUnmcd2CH6lDWB+3TqHp9TQ251xT17lSTeOsMfFq6aP4MeQn1dAeb4HNcJhve/bUW2H9bNshdZd7A6+wPdcrfnmGF+Z5UpuF9TzTJGJIWrtNXH/TpE5I0mIdE9Y31qW25y/T4nrqcLXt0adHx/VldvQ95p14TT+xHeupZZ2OJUlKp9VrDTdLNXlVxffhSeZhdZcGIflknkLOh1xWRjeNtD2j1Cesb2rutZI0wtxvd9bLtqdJjt+nBceY2CRJckkHt/oWF7b06UHx/VmSmjWP7yVpjls/kj79MK538i1az9RbFHoeOjeu/9DUJenV+t2HJSnVdMgac2g8uOjCuN4wvndKksauHtcLSU173h/X78mFteJyUt5yGwApnxWnb/z4Tj/LX+6J62nPE33Tq5fH9aN8S74sPre0sU++kPqH1Vt1ku1wl1ghxEYXxiEnauyDi9Rw+R3Cegc9YnumJYXrmG+QAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKhggwwAAABUFALVAu9JusqMbRmXmz7sDzfvyrjetXAKLjyrZ6Gnt4lze6nw299q4oX6PluY6DJTd6+ZpD1MMkx+1mRaSRpiYtH6aSvbkw95Oqw/3agU89crLnc5x7e0MyFihViWfqPjczjydd8zzw8V9Zi+UGMGx3Fu2iMuX9zJZM1I2kPxSc4pnEN3xTE9TfQ72/OKXgzrG5soN8mnG01WIRbNBCkON1Fukk+OWnSLf9PTIW+E9b+Y+4gkaee4/NlJcWScJK28TDw285eFdf/rOCZsPe1tW+borrB+r0yOo6S+/gyK2kg62IwNOC6u73yWy/2ShpnYxhHb+nN44E8uArKz7XG3wZcL0ZkPmvouh21ke5rKXK/+1NTh1Djqatrlf/ZNah9WTzvIxI1J6j83XndN4+klSfOyuV9Nbmt70kpxPf/ez9NJ58XT/LWQDVc47881d2XpL3Gc27s/jg/cXoVMwIf+GpZvvd+3uH1D6+Q3B9edFX9ojzDJdJI0ytxuCktS6uYGCrmxG5wb1580dUlpGbfG7QkoPx9HeA7rbVt0vKl39y061ySsXrS8vxfPnx/fi5s0ObUwUxwvyjfIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSknEspBv/1w41qslqMCcdemR33dD/LH+8189TnqMI5dDT1QlCErjP1ONui1kt6Iayncze1PRueG9dd8oYk/c1FJ2hWocsc8djCK3fV4WF5uq63LW30SVhvpHdsz8VaK6z31+m25wwNCuvulZGkjZTG5px9hIPRviblA+IlrFtNz0z5p9KX01thfUThHLZT/GjuOe7pe0lTTf3aSX6ex+O3QoUMEl1k6v55Zsk9uFxK8piso8P6LrrG9pxs6qXf57mZcb1ZIZVjnjYJ651MkogkvW1CZ4bHgTOSpP2S6rWGU/OarC7xIs5jR5quJvZ4s0z6jQkZkiStaer3FnruiV9WKb41SZLyzwfEA1u4jBYpnRd/pp2xmY9cuMhe40Ntj9O48Q52bMEC80E5u6U/oLndb2iSKiRpC1MfLJ+OlF40izUOnKl1UP3uw5KUalbKGuPu8vFnWW/5yJPnzfuuOHChlgvoKG1Cfm82NU/4GIu8bfewPjGPtz3TTX1L9bA9WhTfF5o29GvfBN9o8OzC3rBlvAfYU6vYlrsnxPUflGLLDH83k942t8DU55++Ka0RrmO+QQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUNGoLj+85iJpsEmpGW168oUudkg6QX3C+hXqXuh5Jawfazu8lwvJOs1mx3FuW53re1qber/iWXQPq9Nl4o0ktdM9ZqCQTdPxhrD82lQf8/aYlg3rF+lN29PXxLyNMlFutceLg/je1JG2p77myq/VmdrGjLgsIOkzTQnrWwz1kW2PHx+PPW07fOTdZibKTZLambqLcpN8ItKKhaykxbrUjPgII6lnWB2hVrZjP/0mrHcpzDLaxLmZ9DdJUgs9GtYn2yBAKXWIs7gWfXSa7dmvcA5F8yWz7DTL3FNLkW0ujq8UNDnAJU25HE5JJ98f1yf69Dxp3zjObdjwQgTVZjMKB3TcNd6r0PNyWF2w4Oe+5Ykb4/o2j/ien8WxceN29C2/fTiuN5PPHbzSxPAd51Muv5wPVpVuuzoe2z4uP3+GP9zLN8X1vQqn8K9mcX3RxoWm2SbObRt/rf/CxLn5O4o0U783I4XQzYZx7Og8mddZ0uC34shNrb7Qz2PiYUv3jNQ13nCdIbOplHSZqe9cmMfFEnbSGrZlsqnzDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFSnnwtPA/6Umtc1jtJsZHWjqHxSOeHFYHVJ4tjPukP6tpwrzxPkAK+p82zHJ1FsXnhzuYJ4cLp9b/ETqVJnH7yV11DQzMqTO84zW/rajp5lnSOFJ6H46zo5Yd68d1/eKkw4kKWn02JxzjT+o6atpkjVmpXBsQ70T1gvPDNvVvWah54UUX3MvFbIVUv4onj/5R8xbqX1YH5+n257uahvWO8gnA0wzPdLZtkeaaOo+y+MMvR7WDy/M0sTU7y30nKwP44FXC7E3G8RPjy9XeHL646T6reEWNVk/HBOO9TGJEO7uLEkbbWAGSvEgJlYlj02+5+K3zcBVhYnid2qh/mo7GpvEAGl4YR6XnlLKGjk0rHZNf7Yd7n4xbZ3CNEPj8qbb+pY7TH3XwjTjPozraRXfk+ekeq1hSUo162WNMUkN28UxKV0f9+trws/NXubGa2zPrTomrPe1+xxJJ/0xrl9u0i0kHaJfhvWbXXyCpNQpvrefo/heI0kD5eJlSvfi+BprIPN7SrpZc8P6wcnH2KRsPqv+8qw/tR+/GpbzHT/082wbr4Pl2vq14+7FfIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCiTjFvKaUZkgrBJMDXplPO2eWLWaxhfIuwhvFdV681LLGO8a0SruM6bZABAACA7zv+EwsAAACggg0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKhgg2yklOZU/lmcUvq08r8PWIrz7J1S+nNK6ZOU0jOf87NbppReTSl9mFKamVL6Y0pp5aV1Lvh++Tau4f/quzmllFNKay6tc8H3z7dxHaeUtlhyLtVz67u0zgXAN6/RN30C31Y55xb/+f9TSm9LOjzn/MRXMNUHkoZIWlfSVp/zs69J2j7n/O+UUlNJ50u6WtIuX8F54TvuW7qG/3M+m0rq/BWcC75nvsXr+N8551W+gvMA8C3ABvkb9p8bfUrp8C/ws9P/q7RIEt++4RtVlzW85OcaSbpCUl9J47/CUwO+sLquYwDfb2yQl5KU0hmSznDjOecVltI8q0r6q6TlVbtBPmJpHBf4utawpJMlPZdz/mtKaSkdEqj1Na7jlVJK0yV9Iul+SWfnnOcupWMD+IaxQV5Kcs4XSbroa5hniqQVUkqtVbs5/vtXPSf+b/g61nBKqaOkoyT1+Crnwf9dX9O9+O+Sui/5v50k3SppsGrXNoDvAR7S+47KOX+g2pvyA0v+lTXwXTBE0sCc8+xv+kSA+so5v5tzfi3nvDjn/Jak/pL2/KbPC8DSwwZ5KUkpnfVfTzT/r3++omkbSVpJtf+5BfClfE1reGtJl6SU3k0pvbuk9mJKaf+ldHz8H/cN3YuzJP57IeB7hG8el5Kc84WSLqxrX0qpoaTGqn0vGqSUmklalHNeEPzs7pImSJokqY1q/5XeuCXfJgNfytexhiWtrf/9B/NpknYWD+thKfma7sVbSPqnpKmSVlHtf9LxQP3PGsC3Dd8gf/MOkvSpauPaei/5/6//z+CSbz16L/mfK0saKeljSa9KWizpZ1/r2QL/ry+8hnPO7y3519Pv5pz/8w3y+znnT7/ukwb+S13uxT+S9KKkuZL+LOlvkk78Ws8WwFcq5Zy/6XMAAAAAvjX4BhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACAijrlIK+YmubVtGw4NlOd6zzBu6Y+p8dfC11d4vLsxr6l5T/NwAqFeZqF1dU00Xa0cUmuC1ezPWNXbhMPtJ9ue3r8819h/d+z/N/ea2ZRU02zPXPUIazHK6DWDFNvV5hn7IbxPGowz0809m/v55zbFk4llFZslbXaD8KxxmOXCevrFo7XpPHYeKDwV7d8ulpcX2bs2oWZzN9B0Ph137KgpRmYb1te6RGnrXX1s6jxwrg+u3DxtxwbX6//0A9tz9prxK/1olZ+noYz4/ob7oKQtKqpl/5Od/MSKJfugGMX1msNr5ha5tXU3ozOCquT7W9Vm1MWKeUbdegRvxelQPbZ9vuY5rZnZX0c1t8x96Za8dpadax/qduaM5/d4y3b02hsfL+N7yK1XjH1XPp1VjD1zwo9rV1q4j9sS4/x8SqesNB/rnymsfVaw8B3QZ1i3mpSqzxGW4Vjw3RvWC98duliU38hr1LoeikuP7Kyb9lhHzOwc2GeeCN+q2psx8HuM2v6zbYnXXhIPHDmpbYn73daWD9nuH8v+5p6Z51ve0ZpQFjvZjukG0y9X2Ge9Ek8j5bxH05Ka4zNOfs3w7XVdM0aMzwc65A2COtmxUmSOnY0G9ctfc/4W+N6t/RkYaYm5gR6x3VJmurW9xTbsmKO/5Q3wc+idmYT+lBhE7pT6hjWty6c25N3x6/1rMJf8Nvqlri+6yG+Z6ipb+xb9L4ahvUFau2b0ox6reGatE4eo6vNaHwfPkpX2uONNnX/xyjp7By/F3cWeh4yG+HF2sT2DNITYf10/aowU3wjvjIdbTuO1R1h/aF8gO1pl+L7rfkKR5LU2vx5acHZhaadTN1/VyMd+KoZ2Nq25Pbx1xvrT/efKxOU6rWGge8C/hMLAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrqlGKRWtVkbTEmHlzPNBUehX7FBDV0z/f7pjV2i+v/LPS4QLnn/FPNW28WP6X9xKaFaUZdVBh04nSC6fkU29Eu7R7WHzNPsEvS9hoUD3Q5w5+ae7L6gNJz2rE05jU79rF5BtrFT0lS+3o+PZ1qNswa80w4tmmKY9HuMU/s157HYWG9g26yPdNOiq+5Tpf7ed6+Ja6POsS2qFdy67GXb9ouTsU443F/nxiU4yC8W+Uj6A4eaQZ+uq/t0b5x+ki6s5Bgo4NNfaNCzwWm7m9mhyhODbhF7/lp0kr1WsPLpZrcXfF9+B7T06GFP14+MK63ucb39DfXxJmFyLaGJpfsYi2yPSfrcjPiQkIl6biwmtM7tqO3WQ+FDBI9oNPDenrY3GtLB/zxGr7nuTim9JXNfYtLJjncxgNKDVrEr+nCOR/ZnoZqSYoFvrf4BhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACAirqlWKxSk3WcSbE4809x/blt/QFfjsu3nuZbxpunpwe/6n+PDX8Y18fFoQW1PvxJWD45/dm2DNaaYX1jvWF7XtKIeKDbrv7cxnczA64uqf+wuH64b9FEU9/1pUKT8Y+N/djat5sBn7aQtFq9np5eO9XkK0wCgHvye1YhxaKnqe+nbQpnMTSsPqg4DUKSdnIpMRP9a/QjvRDWxz1rW6TNfh7Xf3ijbdkzDnDQ3bv71033xWtoSCFdop953Z7JJ9ieESm+Lwxt5E9twQJzL9OqvumgdeL67z/xPal5vdZwzTopj7k6Hrthq7h+VeF4454w985tVrA9D2p2WN9ZlxRmejOs5jt8XMYz+8f1KYVZ+s6If5+8Upw6Uat/WJ2qNrZj1fxWWH9BPpGiV+oX1t/MQ2zPmiYJZUHyiSGNL1k2rOfTLrM98036RxM1tT1JIsUC31t8gwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABABRtkAAAAoKJuMW8/qMk6ysS8HWuatvbHa2rioT7SfNvTRCPD+g3axfbcaeo+IEsaqL5hfXndans+1EdhfYaWtz3tNCke2Gltf3IPdTQDz9uWx9QprG9voulqxbE/w3Sy7TjYxdZplO2ZrkFhvf3ltkU6KdUrXijVrJE1ZqAZPcPUm9jjnaI46mmw/mB7GuiwsL5owlzb466jp6f769ddRQf4WbSzqd+yQaFpe1O/qNBjlkN+rxAN187ULyjM83gcxLeay5iUNDnHN6acfByZzouD1NKOhftrTf3W8Co1KZ9obsP97zBNQ/zxHjIvxU7TCidxSlxe291sJT1p6h2f8j1DTGzdEN+it83l+uYC/150Nvc6/8EmDVHXsN7v7sIa7mzqPqlQ6YVTw3obXWp7ZimOgFs0uXCPWc2EVvZ3AZhSupiYN3x/8Q0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABV1S7Go+UHWmCPN6ONh9WW9aI+30V3x3Ffu48/hQVMf2cP3HD/WjzlXaHJYH59Xsz1TUvz7+AwEqU9+JKx30o62xz2lrQVb2p69FD8qPth2SB3VPayn3NL2DEjPhvWztNj2rGH+nDYtv+NPLq1SzxSLdllj9g3HjtTvwvp19/trZN7P4qyIpuZ9lSS9s1tYvnNl//T7vvY9v932nKH9w/qgv9sW5eXjc0gvFe4Tu5lrfMYmvqdtnP6xodawLf+TXIyFSyWRUn4gHnj1YdujDWbE9fPa2pZXzo3r3Vf30+it+qVY1DRPeUwcoCCN7mIG7rbHm5vXD+vNXyqchAtDON+3XPxMvIZONek/ktTwQpMAdOZCP5GOMnWT0iDpeh0T1kv37r46MKzn5BOIVsznhPX3h/p51j8+rs/xLZqs+APxDPkPw26mvq+LtpGUHiLFAt9ffIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCiTjFvK6WavI/GhGNDXfJXKSfnvVfD8gdpA9sSB1dJTxemmbd5XN80TiSTJPUy9UF5Bd+03YdxfYpv0d9N7M6lhdy67eNy/uFM23Kf2oT1l/0sGrRyXL/+Xz6O7IhV4vV0eiGxbXCjuL5hIcnpZdUvIivVLJ81ZmMzavKM7j/RHi//LH4FU36ucBYmoGnsubZjuR7xa/7Ren6WxRNnh/WG55roLEl6xdRLeYCrxxdYB/nXYJrixTVFfqF0NDFvK+bptmemzowHtrvQ9uhYU98tjgGUpOV0kj+e8XGqX0RW+5qUD4hvw/ptiuPcRmlPe7wbcry2bpHPqMt3xDF9OsDcnCT11siwXlpaG7nb4EWFpndN/cB1bcuCT14P632X9dPctCiuNzX3M0lK+VdhPa96nm8yEXB77eJbHlDDsD5E5qQlHevuJYXPrzSXmDd8f/ENMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVdUqxSC1qsn5oHp/+c/z0+YhkohAk7XJWXM8XDrI9o3V6WO+pp2xPenireMA97Vw7Uezqn/ueT28My20KT0K3N/UJuRCxofFhNaddbcdn6hjWl+nh/4x0vAnYuKKfT7HQa3H5wcf9OtvL1Odt66fRn+qZYtGkJmsls4b/9cO4vspf7fGuNKELp5hkAEmak+aG9UflF8pEUz/dB2zoXRO60P5c3+MSUi77ceH3MfUBctE20r2KEzZ2n2RbpK5xefH7vqXBNmZg9KO2Z1j+aVgvJeXYxIdr/umbjqnfGl67JuWhZglvZdbWYvlohaYHmqihIwonsdm5YTmnT23LVMX39dWb+WkWfjY5rKf8gG/6SXxRDPqzX8MtTH24n0W/NfWeJgxHkpo9GNc/M0kVkjT++LjezYcW6cw4tEgXaRPbk9O/wvroQozFRvVMEwK+C/gGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABAhc/+iXSQdLYbHBVWd7lxb388k111n4lyk6QL7IiJcpMkk2Kkx3yLt4cfMsebWUhsmxm/bBqRNrc93RWPlSKJTPqa8taFyLax18b1IXYRyOXj7awRtuMi7RLWT37cZKhJUuG0Sxr8cKyaj4mbXapVq3/5yY5d1dTT32zPEBPn5sOUpMdN3UW5SVI7HRfW87kuNE5K3eK4xH77+XlcXOLZG8VRbpKkzqZ+VDvf0296WG5wqG9Rczfwhm0xSXfFVMhOKY5zm7xBoameJi3uoe0/NTlvf4/LI9f1x2tzW1yfWci12/+dc8N62s336CJTX+ca39P9mLA8L82zLU3Na97az6KrTH2c/mx7RusnYb3nLD/PZ8+bgXt9z3AT89a9r4/OzJfF96zfrH2K7UkH7hkf67Y1/ckB32N8gwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABARd1SLD6WZJ5svnPHOK1iv3P94QZdGtf3LZzV/IUjw/ry6uOb5vghp+tdcX3CRTv4piPj8gczfMv+pr7LOr5nw9fj+rg/+B5tGZcHrnRwoamXqU8t9Lh0ABdbIPVTHJFw8qV3Fuapn+aSNjJjvx0Y1/c6xx8vnRs/Sd70577HvBV61OZoSC495YbCuj9FV4b1j1sWpnExMRcXerbeN64/YyJaJKm5WUO3x0kVkqQhpt7Nt6Sa+P3Jz1xme1x4Q+G30bGm3uVV3xNnt3wBb0s6xIyZJeRzS6SZD8f1fLRPb1leJkHhj2MLM50RVp/WE76lv5l//6a2ZUNTP0I9bM8His97L5NUIUnHxoEvxbiTX/WO6+cN9T3u4+uNh/z7MyzF78/Bus72LLrto3ig3Zu2R4XLFfiu4xtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVKWcT1xNokGpyo0ZjwrEPFsY9zfSJPd5gLRvWSyFiQw+L62/c5HvWPDCub3RbYR5T38gnBSmZ+LWzCzFzLjqqFMvkTDuxMHh5/D7k1s19zyxX9mvGRWS94mexMUZrFnqOUxqbc64p/EioZu2Ux/zODNYndytOX9Paw/1rNCnPDuvz0jK2p4nuCOsX27wv6fRLzMCpj9iec7RjWJ9iO+xLoJ1SF9uzcG68whuZa0iSvViuOt63HHebeR9u8D13PhPHZ5kUR0nSBBOrVXrdNq3nGk5Na7JWju/D+qe5qa1iboKSDnknrg/MPkbsXvP79tP5tudiDQjrcchjrQvMOZSSO28/Pa5fNcj3HGt+1TMLH4+/cS/PFr5n1jPxAVu9719rdTL1uVf4nlYnhOWJH/gWp0t63o4l9a7XGga+C/gGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAICKOqVYpAY1uWHT+Onp330W9xyrwtO5eimep+VGvuW6uJz3cc/SS3vp3rB+tp/Fpkjst7nvmffs/LDeW01sjwmK0KRCWobujst5Df88+BzdGdZbqPBauxdbpxR6Bpr6K4We4WG1t16wHS/UN8VixZTHxEEN0mOmPr1z4Yi9wmo6+lbfsr2pF2ISTp4R1y9b3fck0zP4Y39NuhXUwk9jVpZ0caHHpTvM+4u/H+VNTDzB3YWZ9to3LA/P8ZqTpD3iy1iNL/LT/OKcuP7bdLDtSRpWvxSL7jVZfzIpFuY81LpwQBc982cTTSRJmh6Xf7myb/n12LCc0122JeVHw/qVetX2rGrqOz1lW7T8VnH9Y91oe+bP/3lYb9z0fdszPa8Y1tvLxDNJymmvsD5VfWxPx57mGh/d3/aktnHMRxtzH5GkmfW8DwPfBXyDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACgolFdfnjtLF1j4txSNrEy6cnCEY+Py9e9bDvyPjPNyGDb4wLguvtkHWmqqR/rW3Z9No5zG1qYpqfWDOu9x75he15Y/bawvr4N3JJuN/VuhQg6aX9Tj6OXarkIuLm+peeouGN0YZr6+kTS+Hho+Ltxfd/n3/THOyEeO+XqYbZlsHst5i9re0bvYwZM5J8kLa6JV15jNbU9wzQvrPurS1KKo9n6u3uCJJOkpjs3KU1kwuFeK7S0i+Pc9tOpvueaS8Ly1uf432cLN3CsXwe6yg8VvS3pUDP20Flm4Bl/vEv/bAbcuyQ9qFXC+pwLfUzffr+O4/jS1X+wPdIGYfU4m00nSfG9M3d+znZ8skx83oeYzztJavxcHPOmdnGUmyS11x/NSCHj0UTaddT1vmX0tXH97qNsy/57xTFvPqhQhaA54LuPb5ABAACACjbIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqEg5+6eO/9saNSlfOCYea296tphUOOC9pn7m9oWm35r6Ob5l9/vCcrq08Lu7cIfHfMsU81BzR3X3Te3iSIUHp/tzc8/FFwI2tMAkCswxCQSS1MXUhxTm2cnUS4EU518e16880acGHJc0NudcUzhsaMVUk3dWvIhdUkMr7VI4YvwbN8v+afHTTX3AJ36WRs3juIoztKftGZQXmpEzbc+GKU5w+KCQSDFZcU/pYrlST4T17QvrcUNT/0hH2B53Z0qXn287rjwprh+rk/00XYbE9Qd8S1q7fms4tarJ2sLciP8Yx5100F32eNPeiV/zd1f273l7/TMeWGl123PyjLh+WX7H9khxwo3U0becFkeh5K7+9zn00Pg1ONzPol66JR441kWMSGlcPM8Lf/bn1iu5O9MHtkeKP0PTWZv6FvcB4hKdJGmTVK81DHwX8A0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAICKOsW8pTY1WTuYeKHfrxCWc+vZ/oD3mPqDhZMYEmfRjNKVtqXX3XGEzoN7+d99lzw2rN8qn2hzsIvjef8U26M5cTl1MjFKkvRzE6XUwre0uTx+Dd53mWOS0qDN4oEtnrU9856dH9an2tw8ac38aljPM3/oz23FekZkNazJadl4DS/+h4lautMfb5Z5a1sVXtejBsX1s32LOr5kBn7su+bmC8J68+P8PJuZy+i5ob5nrxPi6+iWQjScS3jsVYh566ytwvqbesr2mMur6GlT73ez/30mmpiwJoXXYM16RhWmrjVZw819eIP4fvtxammPt5yLF1xmxzqembRLetiO3WLqrUvTPHRMXH/natvSeLW4Pn/heoWJdjZ1FzMn6R9m7E3fstDcUhsve6DtySn+LHpMr9me7c2rnS48xPbozHge6Srfk24i5g3fW3yDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEBFnVIs2tekfIB5ePq3aft4oNVj/oCzDo/rnW/wPSYoYuE2vqVR845h/QRNsT1D1/HHcwa8HtcH6gXbs7Y2Des9C/O8a+pPbeB7Xv5r/DR9z1SaaY+w2qSRj2jovTCuP6mXbc9ReeM6zF5r+3omAKyTavI1ihfxlnrFdLlcA0kuPaVj4VH2qeaJ+YN9fMsNt8b1I/Shn+e0OLlgwXk+WaHxsieG9Sn6ne3pONkMrNbZ9szSG2G9g+2QFjaL6/t85ntGm/qkP/mevG17M7Kqb9LdYfUodbId1ynVL8ViuZqsHuZG/MxKYfl4zbDHG6pfhfVTdJ7tibOEpM4pfh0k6T7tGda72Q4faHSyD9KRNrswLF+fzrIth+t8M+JzUG5QHEdzuDayPTo8XpHp+ktsy8fpVDO/N8ukp3QupMT03cQM7FuY6KT6rWHgu4BvkAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVNQp5q3Hj1J+ySSWvbBsXI+DcGo92toMzPc9vzCpOxcu8j0LGsW/44q+RWua+naFHvfrnK2Rtme4+oT1OJiu1qYPx/Xjd/DxXVeYFD4d4ec5Za/4dRvvW9Tf1M8p9LzsIvX+fpdvSvvULyKrXU3WvnFEVv6de9e9Jo3iuLLjTdydJA1uHr9PJ8z112JrE9v0SiG2aaLp+UdqZ3uG5elhva/+YHucKWlvO9ZRu4T1dOoIf0CTSHjI3n7d36L34oGd2tqe/PB6YX2IXrM9/TQ0rO+q423PiPrGvNVskDXGvE5PrB6WG5pkQUlatFxcH+ST4dTC1I970feocVzONaVAxzguMeVzCz0XmPruvuUdEwG38km+57zLw3I+199H9jLxhvcUYuvmbR5/II7OTW1Pr6fi+rCtCzFv+WfxwFt/tD1ag5g3fH/xDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFXVKsahZJ+UxV5vBVqZ+qD/eqFfieq/uvme66bnVt6ibedJ/W31kexr0axnWHxziX69fmPo/CikWC7VZWG/UqrntueGDuF7KYOhi6u1e8j1zN47rLcw5S5LGmsexZ/kWHRKXn/6XTyfYMqleT083qUl5xTjEQtP0eFjPg312SbonXg8fFJ7mb6WTw/ri2UNszwot43kKb5/Wy+4J/C0LXWeYeiEB4K4L4/qQwjQD4/Iu2/qWEfmtsH6k1rA915pomW3MNSRJvzX1bgN8j4udefNI37JmPddw6lKTdXO8iBuat3bRnYUD7nZbWN4/HWhbblccZ5Ty7bankW4O6wvGzfPn5uJv4qCR2nMYdElY763TbM/z+rMZiVM0ao93UFi/uZAss+apZuCSQlqGiW8538wv+USlIYVZJumwsH56utH2DKpnEgvwXcA3yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKhrV5Yf/tZz0i63iscHaJKzv/4rPu/q9S1nb15/Du6Z+ujbwTXl2XP9JHOUmSZddFtf7Xemjx3ZeEGc9naEptmerHB9vu0J+l4tzu9K36O72ZmB6O9vTXHFT1gG256rCOTgugm6Lp+pxsM/xw7Gra0y6IBxLOc6OSqf4dZJ/0TesP1YIHtzehDA1aPmo7ZmbW4T1LimOjJOkp9LlYX1LTbY9o/LrYX1X/cb2zNw7Pjft/bTtkY4Iq621t+24M60e1veN08MkSemDHmbkDdvTQuZ+UbgvXdU1rh/rEx7r7xNJr8RDCy+J7ycLdvCHa1ITx5LdrmmFk3gtrOZp19iOGzrE9bShyw6VdO7Rcb3HJ75Hw8LqnEJHIxN9uCD9zDc9F5dTXsH33PahGXB3QeldE+fWbpKfxjluLZczJ0nxZ4EPugO+3/gGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAICKOqVYvC/5Z/PX/XNYvv1in/qgXePUB/WbalvcyDl61fYMXMGkEMQP0kuSTtYfwvrxH+5jexqtFJ/dcXP8a/Bgip8g3z6evtZQU9/Dt6R34yehm6bf257uJmHjL4W31DluTGGwx45h2YRKfEmfSZoYjvQ26RzPpZ72aMPNFeE7JPU8J66P3tm2XKG58cD7cSKHJG25ootd2NX29Er/iKfptbbtSS+YC+mRs2yPdvh5WB5fSLG4VvPjgUN9+kc+9JCw/lN9aHs6vx0v8GdMUoUkjXcDj/me+kpTpCbHmLG8W1jP6Uf2eNndUN4/wZ+EeS3uM0kVktTXvH2t3S8jaY8ed7oRP5FuCKtPFjpaX/NZWD+j0HPuZreH9Sx/s5to7qnrPRJ/DkhSu8Pi1+chF+kkaefX3fH856Rzz26FwfvrfDjgO4NvkAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVNQp5m2xpDlmLK8Yx9c8tIs/3ms5jkUbXjiHccfE8TX7X+2zx/LsVcN62nmK7ZmiOM5tvWX9uf3jorjeMcVxQJJ03CX+eJZLOCrE/uitOM5tXp5hW/rbkdZ25EE3cKhtkf5q6tcWejYrjBXNkTQqHPFxbr3s0Q42V9DGC/0ZvPC2WcO+RU8rXt/H7llo2mL9uP5MO9syXmuF9W6jBhYmihde0x3b2o71dGNYH5c/sT3DUnzxHd7xp7Zn+anxa/3RsfXIKiy4doQZKF2T9bSCpK3N2Cibu9WscMQ4TzG12cF2XD/94bB+xA/9LOf8NX7N/X1GOkTPhfVbdJ3tOUUnhfU5JlJTklKLuL6RPzU1a9UwrC+eNdv2dDFJjns+5OdJm8Tnfau5J0iS1tncDMQRl5LUJr0X1mfmu/w8S/cyAr5V+AYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgIo6pVi0ktTHjC14Kq7fWzje06Y+WX/0TVe/GJbvUPxUtSTdcbR5Gru5n2bVd8wTzyvfZnsuPuHAsH76Jn4e3WfqN/uWEevG9ccL0zhDH/ZJAy/vGL8Gu2uS7XlU8eP8qWchzkRxesQBvR8p9NTPR1pbjylerNvriLjp4CH2ePOHxfUzdJnteWHNuD74Pf9IeBsX7uASTSTJPDGv1drblm67m3P4wE/TKQ0I65Pdg/SSxh1uBq7xMTGH65V4YOrZtifOr5F01du2R+NXC8sdX/Atmh+XV7uz0FOK6ylo2GOsWo+J36d/n2yazi5MdsGjcX2n+H2VpAO1OKyPf9V/5zLwtPh+MvANf2r6o4uD2dW2XKE4XeK3Pf311Xd0fG4XZN/TxdXTCbYn/c58rjz4ju1Rzcph+c5CKof9ONzB31NnnutGbvXzAN9jfIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKAi5VyIivkvzWpSXnVM3SYYXBjbWSZ+TfsVujrH5Z8VstRcntyHr9qW3vphWH9usp/GZdp99otFtqWJ5oT1m7S87dnS1Dub+LJaF4fVtMFI25FfnWBGXDad1ExxNNS8dWyL2rwe19/fwkcspWc0Nudc448aWyvV5MGKF/EJJtLppcLx2qUr4vPb+3jfZOLXphzkWw4x9Seb+9dIz8flizf0LTeY+iTFUVOSNCD9K6yfn33O2/Xp2bB+ROEyvjdOeNTu8q/BEMX3t36lfLzGZn0v8C3a19SHuwEpaXi91vAqNSmfaO7Db5qeUlDXvEXmM6DhNbYnp7fMSHw/k6TU8sqw/u5s2yIXDvlSv8K6dxGChUsyTTCvwXs7+iaNMvW7fcuibcNyp4b+95n8aXxue/pERN3jPqev8j260b2nTXxPWqVeaxj4LuAbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq6pRiUdMo5TEmXKGRedJ20ad32eP11j5hfafCOZz+jjnf1oWmZU6K65deblsanNoirLfSXNvz/tC4ft8J/jXevUv89HJ6tvC+nBGXGw/zLQseMAMf+B7FL4E2/ZlvMQENercwjcvEeFs+MiRptXo9PV2TWuUx2sqMml9Ys/wBt3swrj/+tm05QZ3C+kA/i1r/Ka7nbWfanpTj/I3H5J/M30kNw/r8oT6JZY5Z3wfYDmlESzPw4RqFroNNfYjtOERxRILPWpDumB/XGzftVugyF/+1vW1HOqp+SSxpjZqsgXFMQYMD4zXcv3Dfes3U2yd/D7pWl5kRn9rhVsT0QvrOBaZ+hTYqzHO2qU+xHdNNxEX7UojFQ+eagfa+Z6ejzbF8YohNkXj9MN+yTrys9tRY23LPFub9fiZO3pAkpSdIscD3Ft8gAwAAABVskAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoqFvMW+qSx+jWcGyiid05p3C89jmOODMpS5Kk6/RwPPDWDrZnrdXjef7+kZ/nJBNnN1Qn2p6nUhwbVwgJUxdTXy/H8VSStH+KM7Ju1zTbc7E6hPX+esH2zNGmYb0U3/WA7gjrvbW/7VnT1LcszNNXqX4xb01SHmNSmM6ZGl8LA7vF66eosIhHmVytXnEqmyQp7RefWz7En1s6Z5OwfrxetD1D1SMeqIljxSRJnU39Dz7iMaf4OppeCAR0I916Fd6fXnH5B4N8y7+TadpplO1p8qCJx5vm4/HSD5Z+zJsOdO/teoUjxmFqOa1Y6HFxbkcVeuJMyfv0su3Y3dxPtFPhLmQS97Tav23L0+b+uNUf/DR3ximl2i8f45tOuzos50v9+3NCnhjWh6qtn8cEbs7RTbbjelM/WQf6adJtxLzhe4tvkAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoqFuKRZeUx9wSj3324/hp7WZ6yB/w7F3jkzo/fvpeknTen+P6rxb6ntcbxfV1bvM9K5gnd7f3LfrD/WG5TdrNthxu6oN80IByip/aH/Zj/14e7B7tbnWC7blvVny8bv7UtOYf4/qVP/M9T5v6H+RjRhqqZb2enl69JuXzTADAwbeYpgcLB7wvflpc2t33ND40rpvwBEnSM+PMwPm25b58X1jffQ8/TbrLpGU0usf2TNeeYf1IP40eUHczMrjQ5TJxbrYdD+W1w7rPypAON8EJJuyh1t1xeeHqvqVx83qmWDSvyepiFvGYa+L6X472BzTxBX184IFGxoE9yicN8E2aElZ/lIfZjnH6dVgfpF/anrNlEkWST2J506zHNXe0Lf5G+Ot1bUve7/V4YLSfJr1hojRm7217tm4Zf0bs5KdRPxOC9NM44EOSNDLVbw0D3wV8gwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABABRtkAAAAoKJOMW8ppRmSJn91pwN8YZ1yzm3r2sQaxrcIaxjfdfVaw8B3QZ02yAAAAMD3Hf+JBQAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKhggwwAAABUsEEGAAAAKtggAwAAABVskAEAAIAKNsgAAABABRtkAAAAoIINMgAAAFDBBhkAAACoYIMMAAAAVLBBBgAAACrYIAMAAAAVbJABAACACjbIAAAAQAUbZAAAAKCCDTIAAABQwQYZAAAAqGCDDAAAAFSwQQYAAAAq2CADAAAAFWyQAQAAgAo2yAAAAEAFG2QAAACggg0yAAAAUMEGGQAAAKj4/wHctV5lJahM4gAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "def plot_xy(state,ax,title=\"\"):\n", " '''Plot the XY configuration given by state. Takes an Axes object ax from \n", " matplotlib to draw to, and adds the specified title.'''\n", " angles = np.arctan2(*np.transpose(state,axes=(2,0,1)))\n", " ax.matshow(angles, vmin=-np.pi, vmax=np.pi, cmap=plt.cm.hsv)\n", " ax.title.set_text(title)\n", " ax.set_yticklabels([])\n", " ax.set_xticklabels([])\n", " ax.set_yticks([])\n", " ax.set_xticks([])\n", "\n", "# Make a table of plots\n", "with h5py.File(\"xy_data.hdf5\", \"a\") as f:\n", " fig, axs = plt.subplots(3, 4, figsize=(10, 10))\n", " axs = axs.ravel()\n", " for idx, T in enumerate(temperatures):\n", " plot_xy(f[\"states\"][idx], axs[idx],\"T = {:.1f}\".format(T))\n", " \n", " # Hide the axes of unused plots.\n", " [ax.axis('off') for ax in axs[len(temperatures):]]\n", " plt.tight_layout()\n", " plt.show()" ] }, { "cell_type": "code", "execution_count": null, "id": "4b31ffe5", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.12" }, "vscode": { "interpreter": { "hash": "31f2aee4e71d21fbe5cf8b01ff0e069b9275f58929596ceb00d14d90e3e16cd6" } } }, "nbformat": 4, "nbformat_minor": 5 }