{ "cells": [ { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Intro to Generalized Stability Theory\n", "\n", "\n", "#### Navid C. Constantinou\n", "#### RSES, ANU, 2018" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Linear Systems\n", "\n", "Consider the linear system that describes perturbations about a basic state,\n", "\n", "\\begin{align}\n", "\\frac{\\mathrm{d} \\boldsymbol{\\phi}}{\\mathrm{d}t} = \\mathbb{A}\\,\\boldsymbol{\\phi},\\quad \\boldsymbol{\\phi}(t_0)=\\boldsymbol{\\phi}_0. \\tag{1}\n", "\\end{align}\n", "\n", "Above, $\\boldsymbol{\\phi}$ is the state vector that collectively describes the perturbations about the basic state and $\\mathbb{A}$ is a linear operator (which implicitly depends on the basic state).\n", "\n", "We are interested on how the energy of the perturbations evolves with time.\n", "\n", "Energy is a quadratic quantity and is, usually, defined through an inner product. Let's begin for simplicity with the case in which the perturbation energy is given through the inner product:\n", "\n", "\\begin{align}\n", "E(t) = \\big(\\boldsymbol{\\phi}(t), \\boldsymbol{\\phi}(t)\\big). \\tag{2}\n", "\\end{align}\n", "\n", "For example, for a finite-dimensional state space, e.g., $\\mathbb{C}^{n\\times 1}$, the inner product can be simply the Euclidean inner product \n", "\n", "\\begin{align}\n", "(\\boldsymbol{\\psi}, \\boldsymbol{\\phi}) \\equiv \\boldsymbol{\\psi}^\\dagger \\boldsymbol{\\phi} = \\sum_{j=1}^n \\psi_j^* \\phi_j. \\tag{3}\n", "\\end{align}\n", "\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Motivation\n", "\n", "Some questions that we will try to answer are:\n", "\n", "\n", "1. Can perturbations grow if the linear operator has no modal instabilities?\n", "\n", "2. What is the best way to initiate linear system (1) so that we get maximum energy for $t\\to\\infty$? Specifically, what's the initial condition of unit norm $\\boldsymbol{\\phi}_0$, with $\\|\\boldsymbol{\\phi}_0\\|=1$ and how much can energy grow?\n", "\n", "3. If we are interested in maximum energy growth for very short times $t\\ll 1$ do we have to initiate the linear system differently?\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Linear Systems\n", "\n", "The linear system (1) has solution:\n", "\n", "\\begin{align}\n", "\\boldsymbol{\\phi}(t) = \\mathbb{P}(t,t_0)\\,\\boldsymbol{\\phi}_0 , \\tag{4}\n", "\\end{align}\n", "\n", "where $\\mathbb{P}(t,t_0)$ is a *linear* operator that maps the state at time $t_0$ to the state at time $t$. Linear map $\\mathbb{P}$ is called the *propagator*. \n", "\n", "When the linear operator in (1) is time-independent then the propagator is nothing else than\n", "\n", "\\begin{align}\n", "\\mathbb{P}(t,t_0) = \\exp{(\\mathbb{A}(t-t_0))} , \\tag{5}\n", "\\end{align}\n", "\n", "defined as\n", "\n", "\\begin{align}\n", "\\exp{\\big(\\mathbb{A}(t-t_0)\\big)} = \\sum_{n=0}^{\\infty} \\frac{\\mathbb{A}^n(t-t_0)^n}{n!} . \\tag{6}\n", "\\end{align}\n", "\n", "\n", "**Exercise**: Verify that expression (4) with $\\mathbb{P}$ given by (5) solves (1).\n", "\n", "From hereafter let's take $t_0=0$ for simplicity.\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Adjoint and self-adjoint operators\n", "\n", "**Definition 1 (Adjoint operator)**: The adjoint operator of $\\mathbb{A}$ is denoted as $\\mathbb{A}^\\dagger$ and is defined through:\n", "\n", "\\begin{align}\n", "(\\boldsymbol{\\psi}, \\mathbb{A}\\boldsymbol{\\phi}) = (\\mathbb{A}^\\dagger\\boldsymbol{\\psi}, \\boldsymbol{\\phi})\\ \\ \\mathit{for}\\ \\mathit{every}\\ \\boldsymbol{\\phi},\\boldsymbol{\\psi} \\text{ in the vector space}. \\tag{8}\n", "\\end{align}\n", "\n", "Note that the adjoint operator depends on the inner product!\n", "\n", "**Definition 2 (Self-adjoint operator)**: An operator $\\mathbb{A}$ is called *self-adjoint* or *Hermitian* if and only if $\\mathbb{A}=\\mathbb{A}^\\dagger$.\n", "\n", "**Theorem I**: Self-adjoint operators have real eigenvalues.\n", "\n", "*Proof*: If $\\lambda$ is an eigenvalue of $\\mathbb{A}$ and $\\boldsymbol{u}$ the corresponding eigenvector, $\\mathbb{A}\\,\\boldsymbol{u} = \\lambda\\,\\boldsymbol{u}$. By forming:\n", "\n", "\\begin{align*}\n", "(\\boldsymbol{u}, \\mathbb{A}\\boldsymbol{u}) = (\\boldsymbol{u}, \\lambda \\boldsymbol{u}) = \\lambda (\\boldsymbol{u}, \\boldsymbol{u}).\n", "\\end{align*}\n", "\n", "On the other hand, using that $\\mathbb{A} = \\mathbb{A}^\\dagger$ we have\n", "\n", "\\begin{align*}\n", "(\\boldsymbol{u}, \\mathbb{A}\\boldsymbol{u}) = (\\mathbb{A}^\\dagger \\boldsymbol{u}, \\boldsymbol{u}) = (\\mathbb{A} \\boldsymbol{u}, \\boldsymbol{u}) = (\\lambda \\boldsymbol{u}, \\boldsymbol{u}) =\\lambda^* (\\boldsymbol{u}, \\boldsymbol{u}).\n", "\\end{align*}\n", "\n", "From these we get that $(\\lambda-\\lambda^*)(\\boldsymbol{u}, \\boldsymbol{u})=0$ and since $\\boldsymbol{u}\\ne 0$ (why is that?) we get that $\\lambda=\\lambda^*$." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Normal operators\n", "\n", "**Definition 3 (Normal operator)**: An operator $\\mathbb{A}$ is called *normal* if and only if it commutes with its adjoint: $\\mathbb{A}\\mathbb{A}^\\dagger=\\mathbb{A}^\\dagger\\mathbb{A}$.\n", "\n", "**Corrolary**: Self-adjoint operators are necessarily normal.\n", "\n", "**Theorem II**: Normal operators have orthonormal eigenvectors.\n", "\n", "The proof for this one is a bit cumbersome. We will prove it here for the case of self-adjoint operators, which are a subset of normal operators.\n", "\n", "*Proof*: Consider a non-degenerate operators. Assume $\\boldsymbol{u}_1$ and $\\boldsymbol{u}_2$ are two eigenvectors of $\\mathbb{A}$ that correspond to eigenvalues $\\lambda_1,\\lambda_2$ ($\\lambda_1\\ne \\lambda_2$). Then we have\n", "\n", "\\begin{align*}\n", "\\big(\\boldsymbol{u}_1,\\mathbb{A}\\boldsymbol{u}_2\\big) = \\lambda_2 \\big(\\boldsymbol{u}_1,\\boldsymbol{u}_2\\big).\n", "\\end{align*}\n", "\n", "But also we have\n", "\n", "\\begin{align*}\n", "\\big(\\boldsymbol{u}_1,\\mathbb{A}\\boldsymbol{u}_2\\big) = \\big(\\mathbb{A}^\\dagger\\boldsymbol{u}_1,\\boldsymbol{u}_2\\big)= \\big(\\mathbb{A}\\boldsymbol{u}_1,\\boldsymbol{u}_2\\big) = \\lambda_1 \\big(\\boldsymbol{u}_1,\\boldsymbol{u}_2\\big).\n", "\\end{align*}\n", "\n", "Combining the two we have $(\\lambda_1-\\lambda_2) \\big(\\boldsymbol{u}_1,\\boldsymbol{u}_2\\big) = 0$ from which it follows that $\\big(\\boldsymbol{u}_1,\\boldsymbol{u}_2\\big) = 0$." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "The following code snippet demonstrates the properties of self-adjoint operators." ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "slideshow": { "slide_type": "skip" } }, "outputs": [], "source": [ "using LinearAlgebra, Printf, PyPlot\n", "Base.show(io::IO, f::Float64) = @printf(io, \"%1.2f\", f)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A=\n", "2×2 Array{Float64,2}:\n", " -1.50 -0.50\n", " -0.50 -1.50 \n", " \n", "norm(A-A†)= 0.00, eigvalues of A: [-2.00, -1.00], and (u₁, u₂) = 0.00\n", " \n", "A=\n", "2×2 Array{Int64,2}:\n", " 0 1\n", " 1 0 \n", " \n", "norm(A-A†)= 0.00, eigvalues of A: [-1.00, 1.00], and (u₁, u₂) = 0.00\n", " \n", "A=\n", "2×2 Array{Complex{Float64},2}:\n", " 0.97+0.00im 0.52-0.32im\n", " 0.52+0.32im 0.33+0.00im \n", " \n", "norm(A-A†)= 0.00, eigvalues of A: [-0.04, 1.35], and (u₁, u₂) = 0.00 + 0.00im\n", " \n" ] } ], "source": [ "function analyze_operator(A)\n", " S, U = eigen(A); u1, u2 = U[:, 1], U[:, 2]\n", " println(\"A=\"); show(stdout, \"text/plain\", A);\n", " println(\" \"); println(\" \")\n", " println(\"norm(A-A†)= \", norm(A-A'), \", eigvalues of A: \", S, \", and (u₁, u₂) = \", dot(u1, u2));\n", " println(\" \");\n", "end\n", "\n", "A = [-3/2 -1/2; -1/2 -3/2]; analyze_operator(A)\n", "A = [0 1; 1 0]; analyze_operator(A)\n", "A = rand(2, 2) + im*rand(2,2); A = (A+A')/2; # this is to convert the random A to a self-adjoint one\n", "analyze_operator(A)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Adjoint operators\n", "\n", "**Theorem III (Properties of the adjoint)**: The adjoint operator $\\mathbb{A}^\\dagger$ has *(i)* eigenvalues which are complex conjugates of the eigenvalues of $\\mathbb{A}$. Furthermore, *(ii)* the eigenvectors of $\\mathbb{A}$ and $\\mathbb{A}^\\dagger$ that *do not* correspond to a complex conjugate eigenvalue pair are orthogonal to each other.\n", "\n", "*Proof*:\n", "\n", "*(i)* The eigenvalues $\\lambda$ of $\\mathbb{A}$ solve the characteristic polynomial $\\mathrm{det}(\\mathbb{A}-\\lambda\\mathbb{I})=0$. Using that $\\mathrm{det}{(\\mathbb{A}^\\dagger)} = \\mathrm{det}{(\\mathbb{A})}^*$ we can show that\n", "\n", "\\begin{align}\n", "0 = \\mathrm{det}(\\mathbb{A}-\\lambda\\mathbb{I})^*=\\mathrm{det}\\big[(\\mathbb{A}-\\lambda\\mathbb{I})^\\dagger\\big]=\\mathrm{det}(\\mathbb{A}^\\dagger-\\lambda^*\\mathbb{I}),\n", "\\end{align}\n", "\n", "which implies that the eigenvalues of $\\mathbb{A}^\\dagger$ are $\\lambda^*$, i.e., the complex conjugates of the eigenvalues of $\\mathbb{A}$." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "*(ii)* Consider now the eigenvector/eigenvalue relations for $\\mathbb{A}$ and $\\mathbb{A}^\\dagger$:\n", "\n", "\\begin{align*}\n", "\\mathbb{A}\\,\\boldsymbol{u}_j = \\lambda_j\\,\\boldsymbol{u}_j\\quad\\text{and}\\quad \n", "\\mathbb{A}^\\dagger\\boldsymbol{\\upsilon}_j = \\lambda_j^*\\,\\boldsymbol{\\upsilon}_j. \\tag{10a,b}\n", "\\end{align*}\n", "\n", "\n", "Take for simplicity that there is no degeneracy. By taking the inner product of (10a) with the eigenvector $\\boldsymbol{\\upsilon}_i$ we have\n", "\n", "\\begin{align*}\n", "(\\boldsymbol{\\upsilon}_i, \\mathbb{A}\\boldsymbol{u}_j) &= (\\boldsymbol{\\upsilon}_i, \\lambda_j\\boldsymbol{u}_j) = \\lambda_j\\,(\\boldsymbol{\\upsilon}_i,\\boldsymbol{u}_j).\n", "\\end{align*}\n", "\n", "But on the other hand, using the definition of $\\mathbb{A}^\\dagger$, \n", "\n", "\\begin{align*}\n", "(\\boldsymbol{\\upsilon}_i, \\mathbb{A}\\boldsymbol{u}_j) &= (\\mathbb{A}^\\dagger\\boldsymbol{\\upsilon}_i, \\boldsymbol{u}_j) = (\\lambda_i^*\\boldsymbol{\\upsilon}_i, \\boldsymbol{u}_j) = \\lambda_i\\,(\\boldsymbol{\\upsilon}_i,\\boldsymbol{u}_j).\n", "\\end{align*}\n", "\n", "These imply that\n", "\n", "\\begin{align*}\n", "(\\lambda_i-\\lambda_j)\\,(\\boldsymbol{\\upsilon}_i,\\boldsymbol{u}_j) = 0,\n", "\\end{align*}\n", "\n", "and since there is no degeneracy (i.e., $\\lambda_i\\ne\\lambda_j$) it follows that\n", "\n", "\\begin{align*}\n", "(\\boldsymbol{\\upsilon}_i,\\boldsymbol{u}_j) = 0\\text{ for }i\\ne j.\\tag{11}\n", "\\end{align*}\n", "\n", "\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "The following code snippet demonstrates the properties of the adjoint operators." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "A=\n", "2×2 Array{Complex{Float64},2}:\n", " 0.67+0.48im 0.11+0.82im\n", " 0.56+0.20im 0.39+0.59im \n", " \n", "eigvalues of A: λ₁=0.09 + 0.00im , λ₂=0.97 + 1.06im\n", " \n", "A†=\n", "2×2 Adjoint{Complex{Float64},Array{Complex{Float64},2}}:\n", " 0.67-0.48im 0.56-0.20im\n", " 0.11-0.82im 0.39-0.59im \n", " \n", "eigvalues of A†: μ₁0.09 - 0.00im , μ₂=0.97 - 1.06im\n", " \n", "inner products between eigenvectors of A and of A†\n", "(u₁, v₁) = -0.82 - 0.52im, (u₁, v₂) = 0.00 + 0.00im\n", "(u₂, v₁) = -0.00 + 0.00im, (u₂, v₂) = 0.82 + 0.52im\n" ] } ], "source": [ "A = rand(2, 2) + im*rand(2,2); \n", "\n", "S, U = eigen(A); u1 = U[:, 1]; u2 =U[:, 2]\n", "Sadj, V = eigen(copy(A')); v1 = V[:, 1]; v2 =V[:, 2]\n", "\n", "println(\"A=\"); show(stdout, \"text/plain\", A);\n", "println(\" \"); println(\" \")\n", "println(\"eigvalues of A: λ₁=\", S[1], \" , λ₂=\", S[2]); println(\" \"); \n", "\n", "println(\"A†=\"); show(stdout, \"text/plain\", A');\n", "println(\" \"); println(\" \")\n", "println(\"eigvalues of A†: μ₁\", Sadj[1], \" , μ₂=\", Sadj[2]); println(\" \");\n", "\n", "println(\"inner products between eigenvectors of A and of A†\");\n", "println(\"(u₁, v₁) = \", dot(u1, v1), \", (u₁, v₂) = \", dot(u1, v2))\n", "println(\"(u₂, v₁) = \", dot(u2, v1), \", (u₂, v₂) = \", dot(u2, v2))" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Rayleigh's quotient\n", "\n", "**Definition 4 (Rayleigh's quotient)**: For any self-adjoint operator $\\mathbb{M}$ the quantity:\n", "\n", "\\begin{align}\n", "R(\\mathbb{M}, \\boldsymbol{\\phi}) = \\frac{(\\boldsymbol{\\phi}, \\mathbb{M} \\boldsymbol{\\phi})}{(\\boldsymbol{\\phi}, \\boldsymbol{\\phi})}, \\tag{9}\n", "\\end{align}\n", "\n", "is called *Rayleigh's quotient*.\n", "\n", "**Properties of Rayleigh's quotient**:\n", "\n", "Rayleigh's quotient is a function that maps vectors $\\boldsymbol{\\phi}$ to scalars. It satisfies the following properties:\n", "\n", "1. $R(\\mathbb{M}, \\boldsymbol{\\phi})$ is real for *any* $\\boldsymbol{\\phi}$ in the vector space.\n", "2. $R(\\mathbb{M}, c\\boldsymbol{\\phi}) = R(\\boldsymbol{\\phi})$.\n", "3. $\\lambda_{\\min}\\le R(\\mathbb{M}, \\boldsymbol{\\phi})\\le \\lambda_{\\max}$, where is $\\lambda_{\\min},\\lambda_{\\max}$ are the minimum and maximum eigenvalues of $\\mathbb{M}$.\n", "4. The stationary values of $R$ are the eigenvalues of $\\mathbb{M}$ and they are obtained when $\\boldsymbol{\\phi}$ is the corresponding eigenvector of $\\mathbb{M}$." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "**Exercise**: Show property 1 of Rayleigh's quotient.\n", "\n", "**Exercise** (bit harder): Show property 4 of Rayleigh's quotient. That is, assume that $R(\\mathbb{M}, \\boldsymbol{\\phi})$ is a stationary point and show that $\\boldsymbol{\\phi}$ is necessarily an eigenvector of $\\mathbb{M}$. \n", "\n", "Hint: To do so, write consider the Rayleigh's quotient evaluated at the slightly perturbed state $R(\\mathbb{M}, \\boldsymbol{\\phi}+\\delta\\boldsymbol{\\phi})$ with $\\|\\delta\\boldsymbol{\\phi}\\|\\ll\\|\\boldsymbol{\\phi}\\|$. Since $R(\\mathbb{M}, \\boldsymbol{\\phi})$ is a stationary point then it must be that\n", "\n", "\\begin{align*}\n", "R(\\mathbb{M}, \\boldsymbol{\\phi}+\\delta\\boldsymbol{\\phi}) - R(\\mathbb{M}, \\boldsymbol{\\phi}) = \\mathcal{O}\\big[(\\delta\\boldsymbol{\\phi})^2\\big].\n", "\\end{align*}\n", "\n", "From properties 3 and 4 it is evident that the maximum (or minimum) value of Rayleigh's quotient $R$ is obtained for the eigenvector of $\\mathbb{M}$ that corresponds to $\\lambda_{\\max}$ (or $\\lambda_{\\min}$)." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "Now we have all tools available to attack the question of maximum energy growth for linear system (1).\n", "\n", "### Perturbation energy growth\n", "\n", "The energy growth for linear system (1) is $E(t)\\big/E(0)$. Using solution (4)-(5) we get\n", "\n", "\\begin{align*}\n", "\\frac{E(t)}{E(0)} &= \\frac{\\big(e^{\\mathbb{A}t}\\boldsymbol{\\phi}_0\\,,\\,e^{\\mathbb{A}t}\\boldsymbol{\\phi}_0\\big)}{\\big(\\boldsymbol{\\phi}_0\\,,\\,\\boldsymbol{\\phi}_0\\big)}\\\\\n", "&= \\frac{\\big(\\boldsymbol{\\phi}_0\\,,\\,(e^{\\mathbb{A}t})^\\dagger e^{\\mathbb{A}t}\\boldsymbol{\\phi}_0\\big)}{\\big(\\boldsymbol{\\phi}_0\\,,\\,\\boldsymbol{\\phi}_0\\big)}\\\\\n", "&= \\frac{\\big(\\boldsymbol{\\phi}_0\\,,\\,e^{\\mathbb{A^\\dagger}t} e^{\\mathbb{A}t}\\boldsymbol{\\phi}_0\\big)}{\\big(\\boldsymbol{\\phi}_0\\,,\\,\\boldsymbol{\\phi}_0\\big)}.\\tag{12}\n", "\\end{align*}\n", "\n", "But that's nothing else than the Rayleigh quotient of operator $e^{\\mathbb{A^\\dagger}t} e^{\\mathbb{A}t}$!\n", "\n", "**Exercise**: Show that operator $e^{\\mathbb{A^\\dagger}t} e^{\\mathbb{A}t}$ is self-adjoint and, furthermore, that its eigenvalues are positive.\n", "\n", "\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Maximum energy growth\n", "\n", "**Theorem IV (Bounds on perturbation energy growth)**\n", "\n", "The energy growth $E(t)\\big/E(0)$ of linear system (1) is bounded from above by the maximum eigenvalue of $e^{\\mathbb{A^\\dagger}t} e^{\\mathbb{A}t}$. The inital state $\\boldsymbol{\\phi}_0$ that produces maximum energy growth at time $t$ is the eigenvector of $e^{\\mathbb{A^\\dagger}t} e^{\\mathbb{A}t}$ that corresponds to its maximum eigenvalue.\n", "\n", "Similarly, $E(t)\\big/E(0)$ is bounded from below by the minimum eigenvalue of $e^{\\mathbb{A^\\dagger}t} e^{\\mathbb{A}t}$ with the corresponding eigenvector being the initial condition that produces this minimum energy growth.\n", "\n", "*Proof*: This follows after the remark that $E(t)\\big/E(0) = R\\big(e^{\\mathbb{A^\\dagger}t} e^{\\mathbb{A}t}\\,,\\, \\boldsymbol{\\phi}_0\\big)$ and then employing the properties of the Rayleigh quotient." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "**Corrolary 1**\n", "\n", "For $t\\ll1$, all eigenvectors of $e^{\\mathbb{A^\\dagger}t} e^{\\mathbb{A}t}$ coincide with the eigenvectors of $\\mathbb{A^\\dagger} + \\mathbb{A}$.\n", "\n", "**Corrolary 2**\n", "\n", "For $t \\to \\infty$, the eigenvector of $e^{\\mathbb{A^\\dagger}t} e^{\\mathbb{A}t}$ that corresponds to the largest eigenvalue coincides with the eigenvector of $\\mathbb{A^\\dagger}$ that corresponds to the eigenvalue with maximum real part.\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Back to the initial motivation\n", "\n", "We are now in a position to answer some of the motivational questions we've posed.\n", " \n", "What is the best way to initiate linear system (1) so that we get maximum energy for $t\\to\\infty$?\n", "\n", "- The answer to this comes from Corrolary 2 and it's the eigenvector of **the adjoint operator** $\\mathbb{A}^\\dagger$ that corresponds to the eigenvalue with maximum real part.\n", "\n", "If we are interested in maximum energy growth for very short times $t\\ll 1$ do we have to initiate the linear system differently?\n", "\n", "- Yes! From Corrolary 1, the most energy growth for short times comes about if we initiate our linear system with the eigenvector of operator $\\mathbb{A} + \\mathbb{A}^\\dagger$ that corresponds to the maximum eigenvalue." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "### Remarks\n", "\n", "The theorem on bounds on perturbation energy growth is general and holds for both normal and non-normal operators.\n", "\n", "For normal operators we have that operators: $\\mathbb{A}$, $\\mathbb{A}$, $(\\mathbb{A}+\\mathbb{A}^\\dagger)$, and $e^{\\mathbb{A}^\\dagger t}e^{\\mathbb{A}t}$ *all* commute with each other and thus share the same eigenvectors. Therefore one needs not need to distinguish between them; only computing the eigenvectors of $\\mathbb{A}$ is enough to give you everything you need. Thus, for normal operators the notion of perturbation growth is *equivalent* with the notion of modal instability. However, this is not the case at all for non-normal operators. The theorem above generalizes the notion of modal instability to take into account non-normality and transient energy growth effects. It lies in the heart of the, so called, Generalized Stability Theory.\n", "\n", "Before proving Collolaries 1 & 2 let's see all these in action.\n" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "**Example**\n", "\n", "Let's take the non-normal operator\n", "\n", "\\begin{align*}\n", "\\mathbb{A} = \\begin{pmatrix} -1 & 5 \\\\ 0 & -2\\end{pmatrix}\n", "\\end{align*}\n", "\n", "and compute $E(t)/E(0)$ with the following initial conditions:\n", "\n", "1. eigenvector of $\\mathbb{A}$ that has eigenvalue with maximum real part,\n", "2. eigenvector of $\\mathbb{A}^\\dagger$ that has eigenvalue with maximum real part, $\\mathbb{A}^\\dagger$,\n", "3. eigenvector of $\\mathbb{A}+\\mathbb{A}^\\dagger$ that has the maximum eigenvalue, and \n", "4. eigenvector of $e^{\\mathbb{A}^\\dagger t_*}\\,e^{\\mathbb{A} t_*}$ that has the maximum eigenvalue for some finite value of $t_*$." ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "**Example**\n", "\n", "1. The eigenvector of $\\mathbb{A}$ with eigenvalue with maximum real part is: $\\begin{pmatrix}1\\\\0\\end{pmatrix}$.\n", "2. eigenvector of $\\mathbb{A}^\\dagger$ with eigenvalue with maximum real part is: $\\frac1{\\sqrt{26}}\\begin{pmatrix}1\\\\5\\end{pmatrix}\\approx \\begin{pmatrix}0.196\\\\0.981\\end{pmatrix}$.\n", "3. eigenvector of $\\mathbb{A}+\\mathbb{A}^\\dagger$ that has the maximum eigenvalue is: $\\approx \\begin{pmatrix}0.773\\\\0.634\\end{pmatrix}$.\n", "\n", "We see that all of the above states are thoroughly different to one another. Let's see what's the evolution of the energy if we initiate our linear system with the above." ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [], "source": [ "# A function that integrates (1) and computes E(t)\n", "\n", "function integrateforward(x0; tfinal=1)\n", " t = range(0, stop=tfinal, length=1000)\n", " energy = zeros(length(t))\n", " upperbound = zeros(length(t))\n", " dt = t[2]-t[1]\n", " xt = x0\n", " P = Matrix(I, 2, 2);\n", " expAdt = exp(A*dt)\n", " for j=1:length(t)\n", " energy[j] = dot(xt, xt)\n", " upperbound[j] = maximum(eigvals(copy(P')*P))\n", " P = expAdt*P\n", " xt = expAdt*xt\n", " end\n", " return t, energy, upperbound\n", "end;" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "slideshow": { "slide_type": "slide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAA1sAAAFMCAYAAAA9cKYPAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAA9hAAAPYQGoP6dpAACJx0lEQVR4nOzdeXxU9b3/8dfsk8m+T1YSdkgUZBVBEBUSrAtqtbW9XK17q22tt1i1t3VpXW6tS6+i1S4/bW+1trWotSYRlVUUAQOasAYSsk72PbOf8/tjkkmGBEggIdvn+Xicx8ycc+bMd47jkPd8v+f70aiqqiKEEEIIIYQQYlBph7sBQgghhBBCCDEWSdgSQgghhBBCiCEgYUsIIYQQQgghhoCELSGEEEIIIYQYAhK2hBBCCCGEEGIISNgSQgghhBBCiCEgYUsIIYQQQgghhoB+uBswWiiKQmVlJaGhoWg0muFujhBCCCGEEGKYqKpKa2sriYmJaLUn7r+SsNVPlZWVpKSkDHczhBBCCCGEECNEWVkZycnJJ9wuYaufQkNDAd8JDQsLG+bWCCGEEEIIIYZLS0sLKSkp/oxwIhK2+qlr6GBYWJiELSGEEEIIIcQpLy+SCTKEEEIIIYQQYghI2BJCCCGEEEKIISBhSwghhBBCCCGGgIQtIYQQQgghhBgCEraEEEIIIYQQYghI2BJCCCGEEEKIISBhSwghhBBCCCGGgIQtIYQQQgghhBgCUtT4FNatW8e6devwer3D3RQhRqyGhgZaWlqw2+24XC7cbrd/8Xg8XHTRRWi1vt92PvvsM44cORKwvef+99xzD2azGYB//etf7Nq1C0VRUFXVv3Q9fvDBBwkPDwfg7bffZtOmTSiKAoBWqw1YfvzjHxMXFwfAxo0b2bx5c8B2jUbjv79mzRqsVisAe/fuZdeuXRgMhj6XefPmERkZCUB9fT02m63XPiaTCbPZjMlkOmXxQyGEEEKMHRpVVdXhbsRo0NLSQnh4OM3NzYSFhQ13c4Q4JUVRsNvtBAcHo7rdKB0dfPn559RXVGBvbsHR0oKzrRVXWxvutnY0Hg/XXXUlqsOJ6nKS969/YSstQ+Nxo/N40Xg9aLwKGq8XHfC1rCxUjwfV46Hwyy9paWhAr9Gg14AeDQaNxv/YGhePBkBVaW1pwe1yAb6udQ2g0YAGDRogyGxGo6qogOJ20/UVpQCqquIF3zZVJSQsDK1ej0ajoa2jg/aODryq6tu3c5+u+xPS0zEFB6PR6ym32SgpL8ejqrhRcasqHhXfY1Xla1deSYzVisZgYGf+F2zcstW3DyouRcWpqjhVBaeq8uDDD5Nx3nloTSbe+te/ePKZZ3Cpvn0ciuK/71RV3nvvPb72ta8B8Le//Y21a9diNpsxm80EBQX575vNZu677z4uuOACwBf43njjDf82i8VCcHAwwcHBWCwW5s6dS1JSEgB2u53m5mb/Np1Od1Y+b0IIIcR40t9sIGGrnyRsibPN4/HQWFtLU3k5LVVVtNuq6aitRWu3M2/GdJS2dpT2dj758ENaa2vROp3o3C4MHi9Gr5cgIESvJ9JsRu0MN2J4qXo9uqAgNGYT7U4X5bW1dCgKHYqCXVWxK4pvURUuu/oaps2ehTYoiM+//JIXfvd77KrSvb+iYld9+//vK69ww003odFqeeedd1i9erX/Nc1msz94BQcH84tf/IKvf/3rABQUFPDrX/+a4OBgwsLCCA0NJSwszH9/9uzZpKenA+B2u3E4HAQHB/t7KYUQQojxSsLWIJOwJc6EqqrUFpdQf/QoTaXHaKusxF5djaeuHm9TI0FehVmTJuFtaUFpaaaq6AhmrxfLIP9R69FocKgqHq0Wj1aLotej6PWoBj0YjEzJzEBjMqE1mbE1NtDh9qALCupczBjMXUsQ0fFxaA0G0OvR6A1oDHo0et9y/Dq0OtCARqvt7MbyLZqu+53rNT22oe3s9+oadacCihc6hxCiKL77igpq1/2etz3Xq6geN3T2xKkeD6rL3Xnfjeru3Ob2oLpPsd7pRHE5UZ0uVIfjJPd9t3QOazwbtMHBOHRaymrraFcU2jqXdsXru+9VuPz665i/dBnakBB2FRbwk0ceCdzX68Xdebznn3+eu+++G4AtW7awbNkyNBoNISEh/mDWdXvnnXf6Q1xlZSV/+tOfCAsLIyIiotcSFRXlHyoqhBBCjEYStgaZhC1xPFVRqC4qoubAAZqOFtNWVoqzyoanvh5amglVVKZarXga6vE2NoHHc1qvo6gq7UCHRoNDr8drMjFr8WK0oSFog4M5XFZGu9eLMTwCc2QElqgoLNHRBMfE+JaoKLTBwWgtFjRG46Ceg3FHVTsXxbdAZzjsCoadi393FTweFKcL1enwBTWnE9Xp9AUyhwOlw45i70C121E6Ojof91jX3tH5uHNdR+d+netUu33Q36YbsGs0mKOiCIuPRxseRk17OzlbttKqeGn1KjR33rYoXlq8Xr7/wAPc+N3vog0P55MdO1i6dOkJj//oo4/ys5/9DIBDhw7xH//xH70CWWRkJBERESxcuJA5c+b42uV2U1dXR0REBGazWa5/E0IIMWwkbA0yCVvjh6qq1JSWUrV3L/UHD9FSegxHZRVKXS1hXoXMpEQ8NbV46urA7T71AXtoVxRaNBo69HqcZjMeiwXCwzBGRXHR5ZejCwtDFxZGvdOJMSqKyORkjBERaMb6dTeqCl43uDvAbQeP3XcbsHSAx9G9j7sD3A7fOq8bvK4et/257wbF0x2cVG/3fUXpe706gF6qngEs4FYbeF+rA60edAbfbdeiM3RuM5z0sYoeVdWheHUobh1ejwbFrUFxgeJW8TpVFKfXtzg8KE4PXrsLxeFC6XCi2B14Oxwo7R2oDufg/Oc0GmlXVTq0GtpUaPF6aXK7qXc6qbfbuejyr7HsssvQRURQcKyU//zed2nyemnu0avWpWcwKygo4JxzzgHAaDQSFRVFVFQU0dHRREVF8c1vfpNvfvObALS1tZGbm+vf1nVrsVgG5T0KIYQY3/qbDWQ2QjGueL1eqioqqCwspG7/AVqLj2JqaeH8SZNwV1bhrqqi/uBBQgATkNi59OSoqwt43KgoNOt0dJhMuENCIDICY3QMloQEFn9tFd6wYJxhZpzBJrQ6D3qPnSCPHZfXhUtx4fQ6cXvd5HmduLw1uJRyXBYXLo8L11GXb7/OfV1eF26vG4/qwat48aqdi+JFUZVe6xVVwaN4fPcV3/bB+H1Fq9Gi0+jQaXW++4BW9X2haFUVnaqixXerUxW0ioJO9aJTFLSKF73qxej1oFc8GD1uDF4XRkXBgIpBVTGqKgYVDGrnY7ofGzvXGVQwoGJWVMxq16JgUn3rDHSPQBwWXcFsiH/O6oxtJ6/joQHMnctJqAooHg2KW4vXa0RRg/AqZhSvEa9Hj+Ix4HVr8bq0KC7wOhS8DgXF4cFrd6PYXaCCxuUiBAgB4roOrtVCUJBv2f4pVds/BSACeDd9or8Nbr0eh8FAu1ZLGxD31VfYHn0UXUQEzsYmrgoPp9Hj8YWz+noqamrYryiowNy5c/3HOXbsGNddd12v92g2m4mKiuIHP/gBP/nJTwBobGzkySef9AeyrnAWHR1NTEwMMTEx6PXyz6UQQoiBk56tfpKerdHB6/VSXl5OQ0kJ0yMicJWW4SorZeMbb6CpqSHS4yVOr8fQj+FHLapCo1FHR7AJR7gFV1Qw7thQYmZNoSlEQ0MI1AV5aMOB3WPH7rbT4enw3e+xiOGhRYNJqydIa8CkNWDSGgnSGTHpTJg7l573zXoTFp0Fiz4IiyGIIL0Fi8FCsMGCRR+MxRCMRW/BYgzBYrBg1Jk6X0ntHl7Icbeq0nsdXcMQu64p8/p62JTOnjav58SPT7TN6wSPq/O2a3H4evA8Dt+2vh57nJ3P6bHuDNKhqoLi1vjCWOet161FcXWt0+J1afB6DHjdBl9PnFOD1+ELbqf9uhoNHpMJXUQEwfHx6CIiaNXABzt2UONwUNnWRnlLC7UuF41eLw1eLz956CF+9vOfA4E9Zn354Q9/yHPPPQdAXV0dt9xyC7GxscTExBAbG+tfYmJimDBhgr/EgBBCiLFLerbEmKZ6PPz7T3+m9PPPcZUeQ1Ndg6WlhWiPh2S9nlCdjpIe+58DviFbBt/v/x6NSoNZS2OYjqZwPa6UMKpCvFSFuCmz2KkJVXCYuv738AAtnUsVcAgc+JYBCNIHBSxmnRmjzti9aLvvm3QmDFqD/75RZ8SgNfjX67V6tBoteq3e18PU2cvU875/u9eDzt6ArqMJrb0RraMJjb0ROhqg5629AVxtpz73aPACiga8aPBqfNOyewwWFFMoXmMwXmMwisGM12DGq/ctit6MV2/Cqzfi0Rvw6Ay4NHrcWh0urRa3BjxocKke3IoHl+LCrbhxe93++x5v5/rOdS6vr2fQ6XXi8DhweB0onT1KCip2xY1dGdhQz/7Sa/W+8GWw+G677hu67wfrgwkxhhBqDCXYEEyoIZQQY4hvnSGUYEMYocZQjLoRci2dqnYO1bSDq73HcM3OxdU1hLNzm6sjYLvG1YHObUfn394Ozlbf58rZ6gt7J3pppTuoeZ3agFuPU4vi6l7ncWlQOrcpHi0aVcXgcIDNht1mA3ydeVldBw8O8S09/e3vHP7oI3RR0RiCgnh3+XIaFYV6t9sX0NrbKW1q4lhjI4lhYaiKgkarpaqqinffffeE7+Oee+7h2WefBaC6uprLLrssIJj1vM3IyGDatGln8B9MCCHESCc9W/0kPVtnj9vt5tixYxw5coSjBw5Q/9VXOI4cIby5hRuWLsVRdBj3sdJTXi/VEKrBFqFSHanBFqGhJgJqwzXUhUFjCKjak/du6TQ6Qo2hhBh8fywff7/rscVgIUgfhEVvIcjQHaYsekt3sNKb0WoGebpsRYH2WmithFYbtFRCaxW0VPluW6t86xxNAzuu1gDBsRAc41uCoiAosscScdzjSDBHgH74A4OqqngUD3avHafHicPrwOFx4PQ6sXvsAaGs5/auW7vHToe7gw5PB+3udjo8HXS4OwLWO72Dc11TT0at0R/KQgwhvsUY4v+8dd3vWh9qDCXcFE6YMYwwoy+wDfrnayh4XN3Byx/C2sDV+djZdoLtbeBs8d13tvrue3y/diheOoc1dgUxLUrnrdfRuc7ZGdwcvvuq9zTOlVaLLjwUwsJo8Hpp0eppVlXqPR6qHQ7K29sobW7mqhtv5M6f/ARtcDBffvkls2bNOuEhf/SjH/HMM88AYLPZmDNnDvHx8cTFxREfHx+wzJo1i3PPPfe0TrsQQojBJxNkDDIJW4PL4XBQVFREWVkZq1atwtvWjuPIYX7z07W0HNxHutFAusZAokeH7gRX3rh1UBMO1ZEaqiOgOkJDdSTYIjXUhIPb4HueTqMj3BROhCmCCFMEkeZI//0IUwQR5ggiTZGEm8IDAlWQPmh4ZzvzOKG5HJpKobnMd9tU1n2/tco3hKw/DBYItUJIfGeQiu0OVCFxgY/NEQEz6olAbsUdEL7sbntAMOu6bfe00+HuoNXVSpu7jTZXG63uVtpd7bS6W2lztdHh6RiUNmnQEGoMJcwY1h3CTGGEG8MJM4X5Q1lf2yx6y+ic1c/t8IUuRzM4Wnw/KjiafYt/fc+le53S1oy3zeELYY4eYcypxePQ4nXqAkKa4h54ONMY9ejCQ3CYDDhMJlr1BurQUe32YrPbKW9r45Jrr+Xbd92FNiSEL7/8ktmzZ5/weD2DWWVlJeeee+4Jg9mcOXM477zzTvPECiGE6A8JW4NMwtbp2717Nzt27KCgqICysn1o6iqIcdhJ1xlI9+qZogQR1nziXqo2M1REQ3mMhopoDRXRUBNnhPhYooNjiQ6KJiYoxr9EB0UTbY4myhzlD1Aj8ld/rwdayqHhKDSWdIeprnDVauPU189ofAEqLAFCO5ewBAhNDFxnDpcANQJ5FS/tnnZfEOsRynqGM//jrnWuVlrdrTQ7m2l1tZ7xdYF6jT4gkIWZwog0RRJhjgj4UeL4HykMOsMgnYVh4nX7esl6hjR7U/eQWntj5zDbRpTWBryNDXibmvE0t+K1Kz2CWeetQ4fH4bs/0J4zjV6LLjQIxWLCbTHTYTTSpNVRo2gptbspa21j6ZVXcf1tt6IND2fPnj3+6fD7cu+99/L0008DUFFRwYwZM3oFM6vVSkJCAvPmzTvpsYQQQvRNwtYgk7B1Yu3t7Rw6dIjCA4XkF+dzqOog3752FY6jRXiLS3DuO0Jco5uURgg7yd+FjcFQEaOhPs5MR3I0yoRE9BPTiExMJy443h+kYoJiCDWEjo5f471uX3hqONp7aTzmm+jgZAwWCE+BiBSISO28n+pbwpJ8QUsnl16OZy6vixZXCy3OFlpcLTQ7m32Pe953ttDsau61j/sMrmcLNgT7QpgpknBzuC+g9Qhm4aZwf49xV1AbMdemnQlV9V2j5g9jPYJZZ0hTmmrx1NXjqa/H09iCt7kNT4sdjx1fT5pDi8ehw+vwXXM2EBqdBl2oEYJNeILNOMxGX6+ZqqXM5eVwi4MlV1/PN+64HY1OxxdffBEwS+Px/uu//otf//rXgC+YnXfeeSQkJPjDWNdtQkIC55xzDjNmzDij0yeEEGOFhK1BNt7DVtfHpN3dzrHWY/zfv/+PzfkbMbQ3kuD2kqbomdChJalOJakeLCe+Dp6GcD0tiWF4JljRpqVimTKVyOnnEB8/iRhLDAbtKPvFXFWhpQJqD0DtIWg40h2omsp8dZpORGeCqHSITIOICb5Q1TNQWaKlR0oMCVVVcXgdAYGsK4g1OZtodDbS7Gym0dFIk7PJtziaaHY1+ychGSiL3kKkOZJIUySR5kiizFFEBUURbY7uftxjGRPhrIuq+iYM6aiDjnpor4eOepTGKjw1VXhrbHjqG/A0NOJpasPb0oGn3RPQazag4Ywa0Fm06DuDmTfYjCPITIveQK2q5ZhTYX+TnQu//m2+/Z3vALBr1y7mz59/wkP2DGaVlZUsW7asVyDruj99+nQmTJhwRqdMCCFGMglbfXjvvff4r//6LxRF4Sc/+Qm33nprv587XsKWqqoUlxWz9aut5Jfkc7juMDUdlYQ52jknNIKI6naS61SS6lUS68F4ghyhaKE9PgxvagL6iemETZ1B3Iy5hE+dgXa0FhVVvL5eqtqDncHqINQd9N2ebBY/gwWiJvpCVdTEwCU00Vd/SIhRQlEVWl2tvkB2XBDz3z9uW7OzGe/JfnQ4gRBDiD949Qxj0UHRRJoiiQrqDmYRpgj02jHWy+v1+HrLOuqhvQ6lyYbXVobHVom7utrXc9bQjKepHU+LE0+bB49dg9epBbX/P9LoTCr6ED26UBNKsBF3kIkOo5FGnY4qNxxpd/FVfTuX/cet3HSz79/NUwWztWvX8qtf/QqAqqoqrrjiChITE0lMTCQpKSngdsKECYSHh5/ZuRJCiLNMwtZxPB4PM2fOZOPGjYSFhTFnzhx27NhBVFRUv54/1sKWoijsK9lHo7aRCkcFR5uPsjN/C56SY6S4NSQ3QFI9JNWpWBtBd4JPiVevxZUUh2nyJMKmziBsWgbmSZMwTpiAxjhKf5X2uqGh2BeousJU7QGoO+yfAa0XrR6iJkHsVIiZGhioQuKld0qMa4qq0OZuo8nh6zFrdPiWekc9DY4GGhwNNDoaffftvscetZ+Tv3TSoCHcFN4rnPmHH5sDr+scU71mXVQVnK2ozTY8lUfxVBTjqSrDU12Np64WT31TZzBz4Gnz4LVrUJX+fzdpDQr6YDCEGtCGm/GYDTgsJlr0RqrRUmL3sr+xg/0Vzdxw693cfsedAOzcuZMFCxac8Lg//vGPeeqppwDfdPnf//73+wxlSUlJBAcHn9k5EkKIQSJh6zjbt2/nqaeeYv369YCvSOX555/PDTfc0K/nj9awpSgKe4/sZeu+rewt30tJ81HU9nriPW5S2rUk1qskNUBivUrUSTpnPGYD2rRUQqbNIGTKNIyTJmGaNAlDUhIane7svaHB5HFCfVF3L1XXMMD6ohNfS6Uz+cJU7LQey3RfqBrtEwYIMUKoqkqLq6U7gJ1g6dre6GhEHWAx5jBjWPfkOubu60GPn3An0hSJTjtKv+NOQXU78FYdxVNWhKe8BI+tAk+NzXe9Wc8es3YvqqefoUyjojcr6C0K+hA9hggLmogQ2gw6WgxGatBS4lA4WN/OgfJGviqu4cf3/zc//OEPAfj8889ZuHDhCQ//05/+lF/+8peAr8D0008/3SuUWa1WDAb5PhZCDK0xV9R4y5YtPPXUU+zevZuqqirWr1/P6tWrA/Z58cUXeeqpp6iqqiIjI4PnnnuOCy+8EPCNL09KSvLvm5ycTEVFxdl8C0PK6/Wy+9Buth3Yhj5cwdVYQWtJEc5jpUQ1uYlvgq83qSQ0gMn/g3Hv6y7c4SEETZpMyJSpmCZNxjhpIqbJk9HHxY2OCSn64mqHukO+INUzWDUW+yqp9sUQ7Oulip3uC1QxncEqMg3G6B9eQowUGo2vlyrcFE5aeNop9/cqXv/Qxa4g1tVrVm+vp95eT529jjpHHXX2OjyKxz+JSHFz8UmPrdVoiTRFBoSwnj1lPdeHGcNG1fekxmBGnzoTferMk+6nqipKSwue0kN4yo7gLi/BU1WBp7oad10DnoYWPE12PO0eUDV47Do8dh3UA8fsgG9mJDOQ2rlcpFfQBynop4Siz3uBmh1/QB8VRnp4GF/86CKqVS0l7R4O17VzoKyBL4trKK9vJzY21t+uoqIinnzyyd7vS6MhLi6OBx54wB/iGhsbWb9+PcnJyf5lNP1wKoQYvUZN2Gpvb2fWrFl85zvf4dprr+21/c033+See+7hxRdfZPHixbz88susWrWKffv2kZqaSl8deKPpH8UuboedD//yK8prammoLMbTUI3ZbifSrRDVAfNbIKoVTnYVkFerwWuNJWzaDIKnTMM4MR3TxIkY09PRhYaetfcy6BzNPQLVgc6AdcB3ndWJmMIDe6i6wlVYklxLJcQoodPqfCUfgqJPuW9Xr5k/gHUtjrrAYGavo8HRgKIq1DvqqXfUc6jx0EmPrdfqiQmKITYo1rdYfLdxljj//VhLLBGmiJFZjuIENBoNuvBwdOfMx3TOia/TUj0e33VkVZV4Sg/jLi/GU1mOp9qGp64Bd30znmY7ilNB8WhxtWpxteqhBsCNL53VYwYmABNQudisYAjSoJ8Vg94Sg37Xb2n+0V/Rx8UyKTSMf/3XMoqbXBy2tbGvtIE9xTXUt7mprq4O+Df+4MGD3HLLLQHtDQ0N9Qev2267jeuuuw6Ajo4OioqKSElJISIiYlT+rSCEGDlG5TBCjUbTq2dr4cKFzJkzh5deesm/bsaMGaxevZonnniiz2GECxcu5Fvf+lafr+F0OnE6nf7HLS0tpKSkDPswwrw/PUbq4/93yv08JgNKYizGlBRCJ0zGMiEdY0oyhtRUjCkpaPSjJmcHUlVf/am6Q51h6mD3/daqEz/PEtMZpHr0VsVOl+uphBAn5FE8NDmbegcze2cwc3Svb3G19Pu4XaEsLigwhHUFs5igGOIscUSYxuYf+kp7O+6aGt+1ZOXFeMpLcFd2XVtWj7uhBU+zA5T+/XmiM3nRW7wYLF4MQQr6EA2EB6OPi8aSkoo+KZXKDvjzOx9xoKKJvcW1HChvwtljzpb//d//5fvf/z4An332GYsWLQLAYrEE9IYlJydzxRVXcP755/vei6Kg0WjG5H8nIcTJjblhhCfjcrnYvXs3999/f8D6lStXsn37dgAWLFhAQUEBFRUVhIWF8f777/Pzn//8hMd84okneOSRR4a03acj2ByG3QhtwSqeYBVtkBeL2YvebMZknULcnAuJvOBydKnTR/eXv8fl65GqO+SbpKLucGewOgzO5hM/LzShRy9Vj+F/wTFnr+1CiDGhKxTFBMUwjWkn3dfldflDWa29ltqOWmrsNdR21Pof19prfZN/KB5s7TZs7baTHtOgNRAbFEuMpXcw6/k43BQ+qr7vtcHBmNLTMaWnA+f3uY+qKHgbGnyzLlbX4Kkqx116FE9lGe4qG57aetwNrahuL16nDq9Th7Mx4AhAXefyBTqTl29aOkNZZij6BRZ0YUa8YcHYg4IIsXwKHzdDaALmkmounhHFvrJGato7OHToEIcOdfdqWq1Wf9jasWMHy5cv9wexlJSUgGA2Z84cUlJShuhMCiFGgzERturq6vB6vcTHxwesj4+Px2bz/WOm1+t5+umnWb58OYqicN999xEdfeIhJw888AD33nuv/3FXz9Zwu+Drd6Fc+11KDxdStuOfhJd9TIa7EIOmAaiEg5vh4C9pDkrBkDwHS9pcSJgNCedCUORwNz+Qs9U3619jce/b5vITX0+l0UJkeudEFZ2z/8VMg5gpEBRxVt+CEEIAGHVGEkISSAhJOOl+bq+bekc9NR013aGs635XKOuopdHZiFtxU9leSWV75UmPadAaAnrE4i3xviU43v84zhI3qmZg1Gi16GNi0MfEQEZGn/uoqoq3qcl37VhVle+2ogxPWQluWyWe6lrc9c2oLo8/kNF4/FG8QBtt5FNr2oXB4iXa4uUvk4LRZ5rRB6sQGoQzxEKTwUS1Q0NaVBHseQPCEmk68iVar5MjR45w5MiRXm3s2WOWn5/PD3/4Q38oS01N9S8pKSlERkaOqtAshOifMRG2uhz/JaWqasC6K6+8kiuvvLJfxzKZTJhMpkFt32DQarVo0TJx+iwmTp8FPEJpRSWFn+XiKtrCxPZ8ZmqOEW4vg8NlcPid7ieHWCF6MsRM9t2Gp0BYIoRafdv0g/APcefUw9gbwN4IbbXQWgktVb7Cv61V0FLpWxxNJz+WwdLZ3s7Z/2I6g1X0JNCPvP82QghxKgadAWuwFWuw9aT7ub1u6ux11NhrqOuo69VL1rW+K5RVtFVQ0XbySZ+izFEB4auvYBZiCBk1f/BrNBr0kZHoIyMxT5/e5z4Bgcxmw2Oz+XrGKo7hrqzwPa5rPEUgA/AAHqwmL96t71FueQeDxcuCYC+2/4xDE2rCGRZMkz6IaruO0kYPh2s6WBTn8I3KCEukqKiIrVu3nvD9vPDCC9x1110AHDlyhL/85S8BgSw5ORmz2XzmJ04IcVaNibAVExODTqfz92J1qamp6dXbNVDr1q1j3bp1eL0DL8h5tqQmJZJ67c3AzVQ123kz/xAlX27BUPMVGZpiztEUk6KthTabbzm2re8DmcLBHAam0O5Fb/Zd06TRgkbnu1XcvmnT3XZf3Sm33bc4mnwBSxlAfRxLtK+XKiq983Zi9/2QOLmeSggxLhl0hn71lLm8Ll8o6+wdq+moobqjmur2av/9mo4anF6nf6bGAw0HTng8i97iC2HB8f5Q5r/tXBdljho1E3z0N5Apzc24bbbuQGaz4amq8vWU2Wy4a+pQXW5/IHM0nOD1tB3EW7wkW7wsC/Zi+NtvaHrvGfTBXi4PM9HwyHSaNEHYOrSUNnk4ZGunoLSRwrJG0qyRvh8sNRp2797NQw891Ov4cXFxpKam8vDDD/O1r30NgPr6eoqKikhNTSU+Ph6tTO4kxIgypibImDt3Li+++KJ/3cyZM7nqqqt44oknzvg1R2Odrfo2Jxv2VZNbaOPLomOkKJWka2yka6vIMNUyzdxMnKYRo70Gjdc1uC+uD/INWwyO9s3sF5rg60ULS+xxP8kX7oQQQgwZVVVpdjb7Qlhn+PKHsh7BrL+TfOg1emItsQE9ZNZga3dPWWcvmUE7dmpd+QNZVRXuKhvuqkpfGKv0BTJ3ZSWe+oZ+TeqhM3rRB3dO6GHxYuhxXx9uQh9vpZkQCsuaKGl0c8jWxlfHGjhSa6e8RaXBrvL2229z1VVXAfC3v/2Nb3zjGwAYDAaSk5P9QxNTU1O54YYbyMzM9L+P0dJzKcRIN+aKGre1tVFUVATAeeedxzPPPMPy5cuJiooiNTWVN998kzVr1vDb3/6WRYsW8corr/C73/2OwsJCJkyYcMavPxrDVk8tDjcf768ht8DGpkM1ONzd10PFhhi5erqZrHQTs+K06N1t4GgBZwt4XaB4fb+2qQqoXtAafMP4DEG+nq+u26BI32KJ8q0TQggxatg9dl8Ia68OCGZdYayr96w/BaQ1aIgJiiEhOIH44Hjf0EmL1T+E0hpsJSYoZtT0kPWH6nbjqanpDGSdQayq0nc9WUUF7spKlA77qQ+kVTEEdQcx/XGhTBusQROdiD4yFcKTKShr4k/vbuSrY40ca/JS2qzQ7u4+3PHB7Pbbbw+4VqzrdsKECZx77rmEh4cP0RkSYmwZc2Fr06ZNLF++vNf6G2+8kVdffRXwFTX+1a9+RVVVFZmZmTz77LMsXbr0jF635zDCQ4cOjdqw1ZPd5WXzIV/w+mh/Da3O7mF/ERYDl86IJzvDypIpMZgNUsBXCCGEj0fxUGev6xXEbB22gKDmVtynPJZeo/f3igUsPULZWJv+3tva6g9h/p4xfzirwFNTA94TTA7Vg9aoHNcz5vHddi6eIDMthFDtMJByzgVEpMyE8BTeyNnK/U+8QGWriqePl+k5amjz5s388Y9/ZMKECQFLamrqiLymXYizbcyFreE22nu2TsTlUfjkSB15BTY+2FdNQ3v3cMJgo47l0+PIzrSyfFocwaYxcYmfEEKIIaSoCg2OBqrbq31T3HfY/FPd29ptVLVXUWuvRTnRjLM9mHVmX8/Ycb1iPUNZiDHkLLyrs0P1ePDU1gb0jAWEsspKlNbWUx5Ho1MDA5jF2yOMedCZVToMYTR4LFR1aClucHOgqp0bv/8g6bMvhPAUfr3uD6y9774+j5+QkMDrr7/ORRddBPgm9Dhw4IA/kIWGhg7maRFiRJKwNcjGatjqyeNV2FnSSF6hjdwCG7YWh3+bUa9l6ZRYsjOtrJgRT7hl7IzFF0IIcXZ19ZD1DGFdoayqvQpbu42GE81CcZwQQwjWYKs/lCUEJ/QKZibd2OmJ8ba1+QJYVwir9A1VdHcOVfRUV/uG/p9M11DFHgEsIJRZvCh6E01qMNUOI2XNCodr7BSUNXG41kFZs8rf8j5h9jxfvbFf//rXrF271n/4yMhI0tLS/OHr+9//PpMmTQLA4/Gg0+nGVI+lGJ8kbA2y8RC2elIUlb3lTeR2Bq9j9R3+bXqthkWTosnKsLIyI564UJmKVgghxOByeV2+3rE+esa61rW6Tt3LAxBtjiYh2De7Y2Jwom+mx2DfkhiSSJgxbMz88a+6XLirq3FXVPoDWMCtzQanmmFZo6IP6qtXzIsx2IPe4kWrAywxEJ5ESaOX7YWl7K9sZX9VO8eaFUqaVOo6fH9i7t69mzlz5gDw9NNP89BDDwUMTewZzGbNmoXFYhnisyTEmZOwNUjG4jVbA6WqKgdsreQW2MgrtHHA1v2Pm0YD8yZEkpVhJTvTSnKkfEEKIYQ4OzrcHb16xnqGsuqOauyeU09KYdFb/GGsK4D1DGMxQTHotWNjKL3q9fom8ugriFVU4q6sRHX345o78wl6xToXrV7FjYEGJZio9FkYYtIhPIXX3tnIK3/P41iTQlWb2msCx57B7B//+Afvv/++P5Clp6eTlpZGUlISOp1cUy6Gl4StQTbeerZO5mhtG3mFvinl95Y1BWw7Jymc7Exf8JoUO3bG0QshhBh9uqa9r2yvpKq9iqq2KirbK7G126hs863rz3BFnUZHvCU+oEfs+F6yIP3YmIVXVRQ8dXV4Kitx9QpknWHMfuoAqzP1DmCGYA/GkO4wpqClSQ3B5jBwrEnhoK2DW3/8KCFJMyAihR89/DTPPf9ir2Pr9XpSU1PJyclh6tSpABw4cICGhgbS09Ol3pg4KyRsDTIJW32rbLL7r/HaWdIQ8AvVlLgQf/CamTB2hmgIIYQYOxwehy+IdYaxrvtdYay6vRqP6jnlcSJNkViDrb16xbqCWaQpckz8O6iqKt7GRn/w6quHTGlrO+VxdGYFg6UrfHVP5mEM8WCweNHoQEVDK8HUOE2UtqgcrrbzZVkzRxs8HGtS2PpVCdHWFAB+8IMf8PzzzwNgMplIS0sL6A27/fbbiYyMHNJzI8YXCVuDTMLWqdV1FVEusLH9SB1ub/dHKyUqiOwMK9mZCZyXEoFWO/r/wRFCCDH2eRUvtfbagN6wnmGsqr2Kdnf7KY9j1pn7DGNJIUkkhiQSGxSLTjs2hsZ5W1oCe8MqKnBVlPvul5f3a0ZFX32x7p6wrkBmDPGiD/LiL9EWFAURKRRWtPLp/goKylsoafRyrFnhWJNCY+dcX3V1dURHRwPwwAMP8O677/qD2PG3kZFjIxiLoSVha5DINVunp9nu5uMDvuC1+VBtQBHluFATWRlWVmVaWZAehV4nXf1CCCFGJ1VVaXW3+oYo9ghgPYct1tnrTnkcvVYfGMCCE8duGGtu9gWw8nLc5RW4y8u7A1l5BarDcfIDaMFgUTBY3BhCPBi7himG+AKZ3qzQlZWcqoE6TxCJMxehiUiFiFSefPmv/P3DHZQ0qTTYe/8ZXFtbS0xMDOCrPXbs2LGAQCZ/BwqQsDXopGfr9HW4PGw+WEtuoa+IcluPIsqRnUWUV51jZfHkGEz6sfEPiRBCCNHF5XX5J+3oCmRdtxVtFdjabXjVk88QOF7CmKqqeOvrTxzGKqvgFBN4aHRgCFExWJz+GRQNPXrIdEbVH8YcqoEal5mSZpWD1XaK6tw8+dvX0USmQUQqq795I++8807A8SMjI/3h689//rN/9sSamhqCg4MJDg4eilMjRhgJW4NMwtbgcHq8bC+qJ6egig37qmns6P7CDDHpWT49jlWZVpZNjZUiykIIIcYFj+KhtqOWirYKKtsrfbdtlVS2VUoYO47q9foKP5eXd4exCl8gc1WU47FVg3LygtlaowZDiIIhyNGrV8wY7EVr6P7T2KEJorJDT0mTygFbB/uq2ilpUiluUqh1malubPMPOVy9ejXvvPMOsbGxTJw4kYkTJ5Kenu6/v2zZMpm4YwyRsDXIJGwNPo9X4fOSBv+U8tUtTv82k17L0qmxrMq0csl0KaIshBBi/JIw1n+q243bZusOYxUV/t4xV0U53tpTD+nUmTUYQrwYLU7ftWIhnRN3dE5xr+mZlyzREDEBIlL5y/vb2FZYRnGjr87YsWYFR+dgnqCgINrb2/3B7NFHH6WioqJXKJPrxUYPCVuDTMLW0FIUlT3lTeQW+GY2LG3oXUR5VWYCK2bGExtqGsaWCiGEECPLqcJYf2ZUPFUYi7PEodWM/l4ZxeHwTdzRRxhzV1TgbWo6+QG0YAjVYrS4MVgcviDWI5DpjIF/Vjd7zVR0+IYqXrT6PyEiFSImcOWN3yf300Lcx3XChYeHk5GRwSeffOJft2fPHiwWCxMmTMBkkr+BRgoJW4NMwtbZo6oq+6tayS2oIrfQxqHq7ilkNRqYnxZFdoaVrEwrSRFjo66JEEIIMVQGI4wZtAaSQpL8S3Josu9+aBLJIcmEGcdGiRdvW1v3sMSyMtxl5bjKfbfu8nJUl+ukz9eadRhDwWBxYAzyhTF/IAsO7BVTVGj0milr1XK41sn+qg6KmxTUsBRefftjCEsCnZ45c+aQn5+PRqMhKSnJ3xs2ceJEpk2bxvXXXz/EZ0X0RcLWIJHZCIffkdo2/1DDL8ubA7bNSg4nK9NKdoaViVJEWQghhBiwruntu0JYzzBW3lber2GKoYZQkkJ7h7HkkGQSQxIx681n6d0MHVVRfNeLlZXhKiv33ZZ3B7JTDlHUgCHcgCFYwWDu8A1TDOmc3j7Ei87YPYsiAFo9hCezs6iOr8pbOFTr4mijwtFGheLOmRSnT5/O/v37/U/5xje+QUdHR8C1Yl3DFGXijsElYWuQSc/WyFDe2EFeYTV5BTZ2Hmug56d3WnyoP3jNSAgdE7+wCSGEEMPNo3io7qimorWC8rZyylvLqWiroKKtgvLWcuod9ac8RmxQbEBPWM9AFm+JHxPXiykdHb5ZE8vKcZcfF8j6MaW91qTDEK7HGOzGYG7DYHEFFH0+/hR1KHqaiCAxY5HvurHINK674z6+LGuhpEnBdVw+XrRoEdu3b/c/fueddwgLC2PixIkkJyej043+/wZnk4StQSZha+SpbXXywT7fNV6fHqnHo3R/lCdEW/xDDWcnSxFlIYQQYqjYPXZ/j1hZa5kviHUGs4q2ilMWfe66Xqxnr1hXIEsKTSLSNPonjVBVFW9dXY8gFjhE0VNdffIDaEAfZsIYhm8WRVOrfxZFY4gXnSmwV0xRocFj4liLloPVDg5U2wlJOYf7Hn8BItMgJJ7omBgaGhoAMBgMTJgwwd8TNn/+fG6++eahOyFjgIStQSZha2Rr7nDz4f5qcgttbDlUi9PTfcWpNcxMVkY8WZlWFqRJEWUhhBDibFFVlWZnsy+ItZVR0VoR0CtW2V6JRzn59WIWvcU/RDE5JLn7erHOxWKwnKV3M3QUp9N3rVjAEMXuW7Wj46TP1xh1GCMMvunrTW2dk3d4/D1jx/eKqfogSlvgaKOXwsoOiuo9FDd1DlFsVFhycRa5ubn+/WfOnElYWBiTJk1i4sSJTJo0yb9YrdZxOaW9hK1BJmFr9Gh3etjUWUT54/3VtPfoR48KNrJiRjzZmVYumBwtRZSFEEKIYdR1vZi/R6xnr1hrBTX2mlMeI8oc1WtoYtdwRWuwFb12dNftVFUVb0NDdxArDwxkHpsNTvbnfFevWLjGd52YqcU3VPEEMygC2LUhBCVMh8g0nMEJfO+nv/JfL1beotJjMBGXXXYZ//73v/2PX3nlFZKSkpg0aRLp6eljdgZFCVuDTMLW6ORwe/mkqI7cAhsb9lfT1KOIcqhJz8Uz4sjOsLJsWiwW4+j+MhZCCCHGGqfX6Zuoo8d1Yl29YuVt5bS6Wk/6fL1GT0JIAimhKSSHJPtuQ329YymhKQQbRv+kEYrL1XsGxbLSztuyU/aK6YKNGCIMvlpi5jaMpjZfj1ioB735uEk7AK+qweYwcrRRYV9VB3FT53P1zfdAZBqthhjC4lL9+3bNoNjVC3bppZdyww03DMFZOPskbA0SmY1w7PB4FXYUdxdRrmntLqJsNmhZNjWW7EwrF0+PJzxIiigLIYQQI12LqyWgJ6y8rdx/v6KtArfiPunzI02RpISmkBSa5A9kXUFsLNQWU1UVb309rtIy3GWluErLfEHsWCmusjK8nddsnYjGqMcYZcYQhm8qe0MjxmCXv8hzX6en3auntFXLgWonh+pcHGnw9YgdaVTIvv4WXnr5d7792ttJTk4OGJbY8/5In7RDwtYgk56tsUVRVPLLGsktsJFTYKO80e7fZtBpuGBSDNmZVlbMjCcmZGx2fwshhBBjmaIq1HTUUNZaRnlrue+2czbF8tZyGp2NJ32+UWv0D0f094h13k8KTSJIP/prfXrb2v0hLCCMlZbhrqoCRTnxk7UaDJHBGCP0GELcvkk7DE0YO68V0xp6RwxFo0MbkQqRadSr4Tzx8l/9Qexoo0JbjzJmt912G6+88goAdrudl156iZUrV5KZmTnYp+G0SNgaZBK2xi5VVSmsbCGv0Dez4eGa7iLK2q4iyplWsjKsJEoRZSGEEGJMaHO1+cPX8YGsqq3qlIWeY4Ni/b1gPXvEkkOTiTZHj/4ZFF0u3JWVuMrKcJX6ApirtNQ/RFF1Ok/6fF1YEMYoM8ZQFUNQB0Z9va/Ycx+zJ3ZpchsoaYb9Njups5ay+PL/gMg0Dta5mTF/OS+sW8f3vve9IXrHAyNha5BJ2Bo/imra/MHrq4rjiiinRJCdYWVVppW0mNE/zlsIIYQQvXkUD7Z2G+Vt3UGsZyhrdZ/8WrEgfRBJIUm9esS6JvAw6oxn6Z0MDX+B59LA3jBXWRnu0lK8zc0nfb7WbMQQbcEYrsVgcWI0NGE0NWMI8WKw9D080aVo8YQmY7n0fjjv20P0zvpPwtYgk7A1PpU1dJBX6LvGa9exxoDJfqZbQ8nKsLLqHCvT4qWIshBCCDFeNDubew1N7Apitg4binri4XcaNMQHx/fuEesMZOGm8FH/N4W3pSVwaGLpMX8Y81RXn3z2RJ0WY0xo56QdbgymNoy6eowhLt809lf/BubedNbey4lI2BpkErZETYuDD/ZVk1doY/uRerw95j1Ni7aQlWllVWYCs5JH/5ekEEIIIU6P2+umsr2yz+GJZa1l2D32kz4/xBDSa3hialgqKaEpWC1WdMcXzRplFKfTN3Niaalv+vqe14mVl6O6Tz6pSdzdtxF9971nqbUnJmFrkEnYEj01dbj4cH8NuQU2thyuxdWjiHJCuJmsDCvZmVbmp0Wh00rwEkIIIYTvOvEGR0PA8ET/MMW2cmo6Tl5XTK/V+3vAUkJT/CEsJTRlbAxP9HrxVFf3CGCBwxSVtjYSnnyCiNWrh7upErYGm4QtcSJtTg+bDvqC18YDNQFFlKODjazMiCcrw8oFk2Iw6kf3FLJCCCGEGDoOj4PKtsqAnrCupby1/KRT2WvQkBDsqymWEtYZxkK7w5jFYDmL72TwqaqKt6kJjcGILmT4r5uXsDVIpM6WGAiH28u2w3XkFtrYsK+aZnuPIspmPZfO8AWvZVNjCTKO7mEAQgghhDh7vIqXmo4aSltLKWsto7S1lPLWckpbfI87PCcvXhxtjg7oCesZxsbCdWJnm4StQSY9W2Kg3F6FHUcbyCmoIq+wmrq27ilSgww6LprmK6K8fHocYWYpoiyEEEKI06OqKvWOel/46gpjLaX+x03OppM+P9QY2qsnrGuJtcSO+uLOQ0HC1iCTsCXOhFdRyS/tLqJc0dR9caxRp2Xx5GiyM61cOiOeaCmiLIQQQohB1OJq6R6S2FLm7xkray075XViZp3ZP2FHQCALSyEhOAG9Vn+W3sXIImFrkEnYEoOlq4iyL3hVcaS23b9Nq4GF6dH+IsrWcPMwtlQIIYQQY53D4wjoEetaSltKqWqvwqt6T/hcvUZPYkhi4NDEzqGKyaHJmHRj9wdkCVuDTMKWGCpFNa3+Hq/CypaAbeel+oooZ2damRA9/BeDCiGEEGL8cCtuqtqqugNYVyBr8U3g4fQ6T/hcDRriLHEnvE4sxBhyFt/J4JOwNcgkbImzoauIck6Bjd3HGgO2zUgI8wevqfEhciGrEEIIIYaNoirUdNT06g3rut/mbjvp86PMUd3hK8x3mxqaSmpYKuGm8LP0Lk6fhK1BJmFLnG3VnUWUcwuq+OxoQ0AR5YkxwWRlWsnOsHKuFFEWQgghxAiiqipNzqaAnrCePWMNjoaTPj/MGBYYwsJS/T1iUeaoEfF3j4StQSZhSwynxnYXH+6vJrfAxtbDdbi83UWUE8PN/uA1T4ooCyGEEGKEa3e3+3vBes6eWNpaesoJO+6afRd3zrrzLLX0xCRsDTIJW2KkaHW42XSw1ldE+WANHT2KKMeEGFkx0zfUcNHEaCmiLIQQQohRxe6x++uHlbb6lrIWXyiztdt4bMljXDHpiuFupoStwSZhS4xEDreXrYfryCmo4sN91bQ4PP5tYV1FlDN9RZTNBimiLIQQQojRq2tCjpEwy6GErUEmYUuMdG6vwmdH68kpsPFBoY26Npd/W5BBx/LpsWRlWLl4ehyhUkRZCCGEEOK0SdgaJOvWrWPdunV4vV4OHTokYUuMCl5F5YvSRnK+spFX2LuI8pIpMWRnWFkxM57IYOMwtlQIIYQQYvSRsDXIpGdLjFaqqlJQ0UJOQRW5BTaO1nUXUdZpNSxMj2JVppWVGVbiw6SIshBCCCHEqUjYGmQStsRYoKoqRTVt5BTYyC2wsa8qsIjynNQIVmUmkJVhJTXaMkytFEIIIYQY2SRsDTIJW2IsKq3vILfQ1+P1RWlTwLaZCWGsyvTNbDg5ToooCyGEEEJ0kbA1yCRsibHO1uzgg32+Hq8dxccVUY4N9gWvjAQyk8IkeAkhhBBiXJOwNcgkbInxpKHdxYf7qskttLHtuCLKSRFBZHf2eM1JjZQiykIIIYQYdyRsDTIJW2K8anW4+fhADXmFNjYeqMXu7llE2URWRjzZmVbOnxiNQSdFlIUQQggx9knYGmQStoTwFVHefKiWvAIbG/ZX09qjiHJ4kIFLZsSxKjOBC6fESBFlIYQQQoxZErYGmYQtIQK5PAqfHq0nt8DGhn2BRZQtRh3Lp8WRnWll+fQ4Qkz6YWypEEIIIcTgkrA1yCRsCXFiXkVlV0kDuYU28gpsVDY7/NuMei0XTo4hO9PKpTOkiLIQQgghRj8JW4NMwpYQ/aOqKl+WN5Nb6JvZsPi4IsrnT4wiOzOBrJnxxEkRZSGEEEKMQhK2BpmELSEGTlVVDlW3kVtgI6egigO2Vv82jQbmpEayKtNKVoaVlCgpoiyEEEKI0UHC1iCTsCXEmSupayev0EZOgY09ZU0B2zISexZRDh2eBgohhBBC9IOErUEmYUuIwVXVbOeDwmpyCqr4vLiBHjWUmRQbzKrMBLIzrWQkShFlIYQQQowsErYGmYQtIYZOfZuTD/dXk1Ng45OiOtze7q+l5MggsjO6iyhrpYiyEEIIIYaZhK0+XH311WzatIlLLrmEf/zjHwN6roQtIc6OFoebjQdqyPnKxqZDNTjcin9bbGhnEeWMBBZOjJIiykIIIYQYFhK2+rBx40ba2tp47bXXJGwJMQrYXb4iyrkFVXy0v4ZWZ3cR5QiLgUtnxJOdYWWJFFEWQgghxFnU32wwriqNLl++nE2bNg13M4QQ/RRk1JHdOWmGy6Ow/UgduQU2PthXTUO7i3/sLucfu8sJNupYPr2ziPK0OIKliLIQQgghRoARMwZny5YtXHHFFSQmJqLRaHj77bd77fPiiy+Snp6O2Wxm7ty5bN269ew3VAgxLIx6LRdNi+PJa8/l8wcv4a+3n89NF6RhDTPT7vLy3pdV3P16Puf9YgO3vraLf+wup7nDPdzNFkIIIcQ4NmJ+/m1vb2fWrFl85zvf4dprr+21/c033+See+7hxRdfZPHixbz88susWrWKffv2kZqaCsDcuXNxOp29nvvBBx+QmJg4oPY4nc6AY7W0tAzwHQkhhopep+X8idGcPzGan18+ky8rmskpqCK3wMax+g4+3F/Nh/ur0Ws1LJoUTVaGlZUZ8cSFShFlIYQQQpw9I/KaLY1Gw/r161m9erV/3cKFC5kzZw4vvfSSf92MGTNYvXo1TzzxRL+PvWnTJl544YVTXrP18MMP88gjj/RaL9dsCTFyqarKwepWcr6ykVdo61VEed6ESLI6ZzZMjpQiykIIIYQ4PWPqmi2Xy8Xu3bu5//77A9avXLmS7du3D8lrPvDAA9x7773+xy0tLaSkpAzJawkhBodGo2G6NYzp1jB+tGIqxXXt5BbYyC20sbesiZ0ljewsaeSX/97POUnh/uvBJsWGDHfThRBCCDEGjYqwVVdXh9frJT4+PmB9fHw8Nput38fJysriiy++oL29neTkZNavX8/8+fP73NdkMmEymc6o3UKI4ZUeE8x3L5rEdy+aRGWTnbxCG7kFNnaWNPBVRTNfVTTzVN5BpsSF+IPXzAQpoiyEEEKIwTEqwlaX4/8AUlV1QH8U5eXlDfg1161bx7p16/B6vQN+rhBi5EiMCOI7i9P5zuJ06tqcbNhXTW6Bje1H6jhc08bhj4t4/uMiUqK6iigncF5KhBRRFkIIIcRpGxVhKyYmBp1O16sXq6ampldv12C76667uOuuu/zjMoUQo19MiIkbFqRyw4JUmu1uPj7gC16bD9VS1mDnd1uL+d3WYuJCTWRlWFmVaWVBehR6KaIshBBCiAEYFWHLaDQyd+5cNmzYwNVXX+1fv2HDBq666qphbJkQYrQLDzJw9XnJXH1eMh0uD5sP1pJbaOOj/TXUtDr582fH+PNnx4jsLKK86hwriyfHYNJLEWUhhBBCnNyICVttbW0UFRX5HxcXF7Nnzx6ioqJITU3l3nvvZc2aNcybN49FixbxyiuvUFpayp133jmk7ZJhhEKMHxajnlXnJLDqnAScHi/bi+o7iyjbaOxw8/fd5fx9dzkhJj3Lp8exKtPKsqmxUkRZCCGEEH0aMVO/b9q0ieXLl/daf+ONN/Lqq68CvqLGv/rVr6iqqiIzM5Nnn32WpUuXnpX29Xd6RyHE2OPxKnxe0kBe58yG1S3dNfhMei1Lp8ayKtPKJdPjCbcYhrGlQgghhDgb+psNRkzYGukkbAkhABRFZU95E3kFNnIKbJQ2dPi36bUaLpgcQ3aGlRUz44kNlRlNhRBCiLFIwtYgk7AlhDieqqrsr2olt9BGbkEVh6rb/Ns0GpifFkV2hpWsTCtJEUHD2FIhhBBCDCYJW4Ok5zVbhw4dkrAlhDihI7Vt/lpeX5Y3B2yblRxOVqaV7AwrE6WIshBCCDGqSdgaZNKzJYQYiIomu+8arwIbO4810PObdlp8qD94zUgIlSLKQgghxChzVsKW2+3GZrPR0dFBbGwsUVFRp3uoEU/ClhDidNW2+ooo5xRU8emRejxK99fuhGiLf6jh7GQpoiyEEEKMBkMWttra2vjLX/7CG2+8weeff47T2T0rV3JyMitXruT2229n/vz5p9/6EUjClhBiMDR3uPnoQDU5BTa2HKrF6VH826xhZrIy4snKtLIgTYooCyGEECPVkIStZ599lscee4y0tDSuvPJKFixYQFJSEkFBQTQ0NFBQUMDWrVtZv349559/Ps8//zxTpkwZlDc03CRsCSEGW7vTw+ZDteQU2Ph4fzXtru56flHBRlbMiCc708oFk6OliLIQQggxggxJ2Lruuuv4+c9/zjnnnHPS/ZxOJ3/4wx8wGo3ceuut/W/1CCQTZAghzgaH28v2I3XkfGVjw/5qmjrc/m2hJj0Xz4gjO8PKsmmxWIxSRFkIIYQYTjJBxiCTni0hxNni8Sp8XtxAToGNvEIbNa3dw7XNBi3LpsaSnWnl4unxhAdJEWUhhBDibJOwNcgkbAkhhoOiqOSXNZFbUEVuoY2yBrt/m0Gn4YJJMWRn+ooox4RIEWUhhBDibBjSsFVeXs5LL73E9u3bsdlsaDQa4uPjueCCC7jzzjtJSUk5o8aPRBK2hBDDTVVV9lW1kNs5pfzhmu4iytquIsqZVrIyrCRKEWUhhBBiyAxZ2Nq2bRurVq0iJSWFlStXEh8fj6qq1NTUsGHDBsrKysjJyWHx4sVn/CZGEglbQoiRpqimu4jyVxXHFVFOiSA7w8qqTCtpMcHD1EIhhBBibBqysDV//nyWLFnCs88+2+f2H/3oR2zbto2dO3cOrMUjlEyQIYQYDcobO8jtvMZr17HGgCLK062hZGVYWXWOlWnxUkRZCCGEOFNDFraCgoLYs2cP06ZN63P7gQMHOO+887Db7X1uH62kZ0sIMVrUtDr4oLCavEJbryLKadEWsjMTyM60Mis5XIKXEEIIcRr6mw0GPH9wQkIC27dvP2HY+vTTT0lISBjoYYUQQgySuFAz/3H+BP7j/Ak0dbj4cH8NuQU2thyupaS+g99uPsJvNx8hIdxMVoaV7Ewr89Oi0GkleAkhhBCDacBh68c//jF33nknu3fvZsWKFcTHx6PRaLDZbGzYsIHf//73PPfcc0PQVCGEEAMVYTHy9bnJfH1uMm1OD5sO+oLXxgM1VDU7eHV7Ca9uLyE62MjKjHiyMqxcMCkGo1473E0XQgghRr3Tmo3wzTff5Nlnn2X37t14vV4AdDodc+fO5d577+X6668f9IYONxlGKIQYSxxuL9sO15FbaGPDvmqa7T2KKJv1XDrDF7yWTY0lyKgbxpYKIYQQI89ZqbPldrupq6sDICYmBoNh7BbXlLAlhBir3F6FHUcbyC2sIq+wmtoeRZSDDDoumuYrorx8ehxh5rH7PS+EEEL015CErdLSUlJTU/vdiIqKCpKSkvq9/0gksxEKIcYTRVH5orSR3AIbOQU2Kpq6Jzsy6rQsnhxNdqaVS2fEEy1FlIUQQoxTQxK24uPjufLKK7nttttYsGBBn/s0Nzfzt7/9jd/85jfccccdfP/73x9460cg6dkSQow3qqpSWNnSGbyqOFLb7t+m1cDC9Gh/EWVruHkYWyqEEEKcXUMSthoaGnj88cf54x//iMFgYN68eSQmJmI2m2lsbGTfvn0UFhYyb948/vu//5tVq1YNypsZCSRsCSHGu6KaVnILbOQW2iioaAnYdl6qr4hydqaVCdFSRFkIIcTYNqTXbDkcDt5//322bt1KSUkJdrudmJgYzjvvPLKyssjMzDyjxo9EEraEEKJbWUMHeYU2cgts7C4NLKI8IyHMH7ymxodILS8hhBBjzlmZIGM8kbAlhBB9q2lxkLevmrwCG58ercfbo4jyxJhgsjKtZGdYOVeKKAshhBgjhixsrVmzhpdffhmLxXLGjRxNJGwJIcSpNba7+HB/NXmFNrYcrsPlUfzbEsPN/uA1T4ooCyGEGMWGLGzpdDqqqqqIi4sD4I477uDJJ58kMjLSv4/b7R5z08BL2BJCiIFpc3rYeKCG3EJfEeUOl9e/LSbEyIqZvqGGiyZGSxFlIYQQo8qQhS2tVovNZvOHrbCwMPbs2cPEiRMBqK6uZsKECTgcjjNo/sgjYUsIIU6fw+1l6+E6cgtsfLg/sIhyWFcR5UxfEWWzQYooCyGEGNn6mw30Z/pCfWU1l8t1pocdMXrW2RJCCHF6zAYdK2bGs2JmPG6vwmdH68ktsJFXWE1dm5N/5lfwz/wKggw6lk+PJSvDysXT4wiVIspCCCFGsTPu2QoNDWXv3r0BPVuJiYljLpxIz5YQQgw+b48iyrl9FFFeMiWG7EwrK2bEExlsHMaWCiGEEN2GtGfr9ddfZ+nSpZxzzjkAMruUEEKI06LTapifFsX8tCj++2szKKhoIbewipwCG0dr2/n4QA0fH6hBp9WwMD2KVZlWVmZYiQ+TIspCCCFGvgH3bC1dupS9e/fS2tqKwWDA4/Fw/fXXs2TJEubMmUNsbCzTpk2Tni0hhBCnTVVVimrayC2wkVNgY19VYBHlOakRrMpMICvDSmr0+JodVwghxPAb8jpbhw4dYvfu3eTn5/tvm5qa/L1cEraEEEIMltJ6XxHlnIIqvihtCtg2MyGMVZm+mQ0nx0kRZSGEEENvyMLWgw8+yOrVq1mwYEGvbcXFxezatYv8/Hwef/zxgbd6BJOwJYQQI0N1i4MPCn09XjuKGwKLKMcG+4JXRgKZSWESvIQQQgyJIQtb3/nOd/j3v/+NTqfjiiuuYPXq1VxyySWYTKYzbvRIJmFLCCFGnoZ2Fx/uqya30Ma2w3W4vN1FlJMigsju7PGakxopRZSFEEIMmiEdRqiqKtu2beNf//oX7777LhUVFaxYsYIrr7ySyy+/nJiYmDNq/EgkYUsIIUa2Voebjw/UkFdoY+OBWuzunkWUTWRlxJOdaeX8idEYdFJEWQghxOkb8mu2etq/fz//+te/eOedd9i1axcLFy7kyiuv5IYbbiApKelMDz8iSNgSQojRw+H2svlQLXkFNjbsr6bV4fFvCw8ycMmMOFZlJnDhlBgpoiyEEGLAzmrY6qm2tpZ3332Xd999lwsvvJAf//jHg3n4YSNhSwghRieXR+HTziLKG/bZqGtz+bdZjDqWT4sjO9PK8ulxhJhOqyKKEEKIcWbYwtZYJWFLCCFGP6+isqukgdxCG3kFNiqbHf5tRr2WCyf7iihfKkWUhRBCnMRZCVtutxubzUZHRwexsbFERUWd7qFGvP6cUEVRcLlcfW4TYiwyGAzodDIES4xOqqryZXkzuYU2cgtsFNe1+7fptBrOnxhFdmYCWTPjiZMiykIIIXoYsrDV1tbGX/7yF9544w0+//xznE6nf1tycjIrV67k9ttvZ/78+aff+hFk3bp1rFu3Dq/Xy6FDh054Ql0uF8XFxSiK0sdRhBi7IiIisFqtMsW2GNVUVeVQta+Icm6hjf09iihrNDAnNZJVmVayMqykREkRZSGEGO+GJGw9++yzPPbYY6SlpXHllVeyYMECkpKSCAoKoqGhgYKCArZu3cr69es5//zzef7555kyZcqgvKHhdrITqqoqpaWluN1uEhMT0Wpllisx9qmqSkdHBzU1NURERJCQkDDcTRJi0JTUtZNX6Ate+ccVUc5I7FlEOXR4GiiEEGJYDUnYuu666/j5z3/OOeecc9L9HA4Hf/zjHzEajdx66639b/UIdrIT6na7KSoqIjExkfDw8GFqoRDDo76+npqaGqZOnSpDCsWYVNVs54PCanILbOworqdHDWUmx4WQneELXhmJUkRZCCHGiyG/ZmvRokXk5eWNm8kiTnZCHQ4HxcXFpKWlERQUNEwtFGJ42O12SkpKSE9Px2yW61rE2Fbf5uTD/b7gta2oDre3+5/Q5Mggf/CakxqJVoooCyHEmDXkYUur1WKz2YiLi+v1wr/4xS946qmnTuewI1Z/wpb8sSnGI/n8i/GqxeFm44EacgtsbDoYWEQ5NrSziHJGAgsnRkkRZSGEGGP6G7YGXFDkmmuuYcGCBWg0GmpqanqFrfb2dp555pkxF7aEEEKInsLMBq6ancRVs5OwuzqLKBfa+HB/NbWtTv7vs1L+77NSIiwGLp0RT3aGlSVSRFkIIcaVAYetCRMm8N5776GqKrNmzSI6OppZs2Yxa9Yszj33XL788ku5UF4IIcS4EmTUkd05aYbLo7D9SB15hTY+KKymvt3FP3aX84/d5QQbdSyf3llEeVocwVJEWQghxrTTHkZoMpnYtm0blZWV5Ofns2fPHr766isUReGxxx7jW9/61mC3dVjJMMKR4aabbqKpqYm33357uJvSp4suuojZs2fz3HPPDXdTzhr5/AtxYl5FZWdJA7kFNvIKbVQdV0R56ZRYsjOtrJgRT7jFMIwtFUIIMRBDNoywS3t7O3q97+lXXXXV6R5GCCGEGLN8xZGjOX9iND+/fCZfVjT7ankVVFFS38GH+6v5cH81eq2GRZOiycqwsjIjnrhQ+eFCCCHGgtMOW11BS4ih4HK5MBqNw90MIYQYNFqthtkpEcxOieAn2dM4WN3aGbxsHLC1svVwHVsP1/GzdwqYNyGSrM6ZDZMjpYiyEEKMVgOaHqm0tHRAB6+oqBjQ/mNNe3v7CReHw9Hvfe12e7/2Hai0tLRew91mz57Nww8/7H+s0Wh46aWXWLVqFUFBQaSnp/P3v//dv72kpASNRsNf//pXLrjgAsxmMxkZGWzatCnguPv27eOyyy4jJCSE+Ph41qxZQ11dnX/7RRddxN133829995LTEwMK1asOGnbH3nkEeLi4ggLC+OOO+7A5XL5tzmdTn7wgx8QFxeH2WxmyZIl7Ny507/91VdfJSIiIuB4b7/9dkB9nIcffpjZs2fz5z//mbS0NMLDw/nmN79Ja2urf5/29nb+8z//k5CQEBISEnj66adP2mYhhOii0WiYbg3jnkunknvPUjb++CLuXzWdWSkRqCrsLGnkl//ez5L/2cgVz29j3cYijtS2DXezhRBCDNCAwtb8+fO57bbb+Pzzz0+4T3NzM7/73e/IzMzkn//85xk3cDQLCQk54XLttdcG7BsXF3fCfVetWhWwb1paWp/7DZWf/exnXHvttezdu5f/+I//4IYbbmD//v0B+6xdu5b/+q//Ij8/nwsuuIArr7yS+vp6AKqqqli2bBmzZ89m165d5ObmUl1dzfXXXx9wjNdeew29Xs8nn3zCyy+/fML2fPTRR+zfv5+NGzfyxhtvsH79eh555BH/9vvuu4+33nqL1157jS+++ILJkyeTlZVFQ0PDgN73kSNHePvtt3nvvfd477332Lx5M08++WTAe964cSPr16/ngw8+YNOmTezevXtAryGEEADpMcHcuWwS79y1mO33X8zDV8xkYXoUWg18VdHMU3kHueTpzax4ZjNPf3CQwspmTvOSayGEEGfRgMYC7t+/n8cff5zs7GwMBgPz5s0jMTERs9lMY2Mj+/bto7CwkHnz5vHUU0/1CglidLruuuu49dZbAfjFL37Bhg0beP7553nxxRf9+9x9993+APnSSy+Rm5vLH/7wB+677z5eeukl5syZw+OPP+7f/49//CMpKSkcOnSIqVOnAjB58mR+9atfnbI9RqORP/7xj1gsFjIyMnj00UdZu3Ytv/jFL7Db7bz00ku8+uqr/s/f7373OzZs2MAf/vAH1q5d2+/3rSgKr776KqGhoQCsWbOGjz76iMcee4y2tjb+8Ic/8Kc//cnfC/faa6+RnJzc7+MLIURfEiOCuGlxOjctTqeuzcmGfb4iytuP1HG4po3DHxfx/MdFpER1FVFO4LyUCCmiLIQQI9CAwlZUVBS//vWv+eUvf8n777/P1q1bKSkpwW63ExMTw7e//W2ysrLIzMwcqvaOKm1tJx7yodMF1lmpqak54b5abWAHZElJyRm1a6AWLVrU6/GePXtOuI9er2fevHn+3q/du3ezcePGPnvfjhw54g9b8+bN61d7Zs2ahcXSfQ3DokWLaGtro6ysjObmZtxuN4sXL/ZvNxgMLFiwoFdv3KmkpaX5gxZAQkKC/7/TkSNHcLlcAe87KiqKadOmDeg1hBDiZGJCTNywIJUbFqTSbHfz8QFf8Np8qJayBju/21rM77YWExdqIivDyqpMKwvSo9BLEWUhhBgRBjzLxZo1a3j55Ze55ppruOaaa4aiTUOirKyMNWvWUFNTg16v52c/+xnXXXfdkL5mcHDwsO97MlqtttcwFLfb3a/n9ry+6VT7KIrCFVdcwf/8z//02qdnTbYzfV8ajcb/fo5vn6qq/nX9fd8GQ+A0zBqNBkVR/McTQoizKTzIwNXnJXP1ecl0uDxsPlhLbqGNj/fXUNPq5M+fHePPnx0jsrOI8qpzrCyeHINJL0WUhRBiuAz4p6/XX389oMfmjjvuoLGxMWCf/v7Bfjbp9Xqee+459u3bx4cffsiPfvSj05pUYiyJjY2lqqrK/7ilpYXi4uJe+3322We9Hk+fPv2E+3g8Hnbv3u3fZ86cORQWFpKWlsbkyZMDltMJWHv37g2YNOSzzz4jJCSE5ORkJk+ejNFoZNu2bf7tbrebXbt2MWPGDP/7bm1tDfjvf3xP3alMnjwZg8EQ8L4bGxs5dOjQgN+PEEIMlMWoZ9U5Cfzmm+ex62eX8v9ums835qUQaTHQ2OHm77vLufnVXcz9xYf84I183v+qinanZ7ibLYQQ486Aw9bxv+i/8cYbAWGruro6YOjVSJGQkMDs2bMB32QUUVFRA54wYay5+OKL+fOf/8zWrVspKCjgxhtv7DW8EeDvf/87f/zjHzl06BAPPfQQn3/+OXfffXfAPuvWrWP9+vUcOHCAu+66i8bGRm6++WYA7rrrLhoaGrjhhhv4/PPPOXr0KB988AE333wzXq93wO12uVzccsst7Nu3j5ycHB566CHuvvtutFotwcHBfPe732Xt2rXk5uayb98+brvtNjo6OrjlllsAWLhwIRaLhQcffJCioiJef/11Xn311QG1ISQkhFtuuYW1a9fy0UcfUVBQwE033dRryKcQQgw1k17H8ulx/M/Xz2XnTy/l9dsWcuOiCcSHmWhzenh3byXf+8sXzPnFBm7/0y7++UU5zR0j70dRIYQYi874L8O+hlP1nIa7v7Zs2cIVV1xBYmIiGo2Gt99+u9c+L774Iunp6ZjNZubOncvWrVtPp8ns2rULRVFISUk5reePFQ888ABLly7l8ssv57LLLmP16tVMmjSp136PPPIIf/3rXzn33HN57bXX+Mtf/sLMmTMD9nnyySf5n//5H2bNmsXWrVt55513iImJASAxMZFPPvkEr9frv6bvhz/8IeHh4acVTi655BKmTJnC0qVLuf7667niiisCpqt/8sknufbaa1mzZg1z5syhqKiIvLw8IiMjAd+1Vf/3f//H+++/zznnnMMbb7wR8Pz+euqpp1i6dClXXnkll156KUuWLGHu3LkDPo4QQgwWvU7LBZNieOSqTD69/xL++b0LuGPpRFKjLDg9Ch/sq+bev+1l7i838J9//JzXd5RS2+oc7mYLIcSYpVEHePGJVqvFZrMRFxcHQGhoKHv37mXixImAr2crMTFxwD0WOTk5fPLJJ8yZM4drr72W9evXs3r1av/2N998kzVr1vDiiy+yePFiXn75ZX7/+9+zb98+UlNTAZg7dy5OZ+9/ND744AMSExMBqK+v58ILL+T3v/89F1xwwQnb43Q6A47V0tJCSkoKzc3NhIWFBezrcDgoLi72B8GxRKPR9Ppv0VNJSQnp6enk5+f7ew7F+DKWP/9CjBWqqrK/qpXcQht5BTYOVnfXDNRoYH5aFNkZVrIyrSRFBA1jS4UQYnRoaWkhPDy8z2zQ04AnyADfdVtLly7lnHPOAfo3WcKprFq16qRTxT/zzDPccsst/inIn3vuOfLy8njppZd44oknAE5Z48jpdHL11VfzwAMPnDRoATzxxBMBtZuEEEKI0Uqj0TAzMYyZiWHcu2IqR2vb/MFrb3kznxc38HlxA4++t49ZyeFkZVrJzrAyMXboajgKIcR4MOCwtWTJEh566CFaW1sxGAx4PB4efPBBlixZwpw5c4iNjR30RrpcLnbv3s39998fsH7lypVs3769X8dQVZWbbrqJiy++mDVr1pxy/wceeIB7773X/7irZ0sIIYQY7SbGhvC9iybzvYsmU9FkJ6/ARm6hjZ0lDewtb2ZveTO/yj3ItPhQf/CakRA6KD+uCiHEeDLgsLVlyxYADh8+zO7du/niiy/YvXs3P/vZz2hqahqSL+K6ujq8Xi/x8fEB6+Pj47HZbP06xieffMKbb77Jueee678e7M9//rO/d+54JpMJk8l0Ru0eC041yjQtLU2mQRdCiFEsKSKIm5ekc/OSdGpbO4soF9rYXlTHwepWDla38r8fHWZCtMU/1HB2shRRFkKI/jitYYQAU6ZMYcqUKXzzm9/0rysuLmbXrl3k5+cPSuOOd7LaSaeyZMkSf42kgVi3bh3r1q07rVnzhBBCiNEkNtTEtxam8q2FqTR3uPmoRxHlY/UdvLzlKC9vOYo1zExWRjxZmVYWpEkRZSGEOJHTDlt9SU9PJz09fdCLBcfExKDT6Xr1YtXU1PTq7Rpsd911F3fddZf/IjghhBBiPAi3GLhmTjLXzEmm3elh86FacgtsfHygBluLg9c+PcZrnx4jKtjIihnxZGdauWBytBRRFkKIHgY1bA0Vo9HI3Llz2bBhA1dffbV//YYNG7jqqquGsWVCCCHE2Bds0nPZOQlcdk4CDreX7UfqyC2wsWFfNQ3tLt7cVcabu8oINem5eEYc2RlWlk2LxWIcFX9mCCHEkBkx34JtbW0UFRX5HxcXF7Nnzx6ioqJITU3l3nvvZc2aNcybN49FixbxyiuvUFpayp133jmMrRZCCCHGF7NBx8XT47l4ejwer8LnxQ3kFtrILbBR0+rknT2VvLOnErNBy7KpsWRnWrl4ejzhQYbhbroQQpx1A66zNVQ2bdrE8uXLe62/8cYbefXVVwFfUeNf/epXVFVVkZmZybPPPsvSpUuHtF09r9k6dOjQuKuzJcSpyOdfCAGgKCr5ZU3kFdrIKaiirMHu32bQabhgUgzZmVZWzIwnJkQmoBJCjG79rbM1YsLWSHeyEyp/bIrxTD7/QojjqarKvqoW8gps5BTYOFzT5t+m7SqinGklK8NKohRRFkKMQhK2BpmELSH6Jp9/IcSpFNW0kdc51PCriuaAbbNSIsjOsLIq00paTPAwtVAIIQamv2FrxFyzNVKNp6nf8/PzWbBgAeeffz5bt24d7uYIIYQYIybHhTA5bjJ3LZ9MeWMHeYXV5BZUsetYI3vLmthb1sT/5B5gujWU7Ewr2ZlWpsVLEWUhxOgnPVv9NB56ti688EKWLFnCunXraG5uln/kRL+Mlc+/EOLsq2l1+IooF9j49Eg9HqX7T5K0aAvZmQlkZ1qZlRwu/yYJIUYUGUY4yAYStlRVxe4enp6wIIPutP5Bev311/nrX//Kiy++SEpKCkVFRUyaNGkIWijGGglbQojB0NTh4sP9NeQW2NhyuBaXR/FvSwg3k5Xh6/GanxaFTivBSwgxvGQY4TCyu73M/HnesLz2vkezBlzXpL29nQcffJCcnBySk5MJDw9nz549EraEEEKcNREWI1+fm8zX5/qKKG886AteGw/UUNXs4NXtJby6vYToYCMrM+LJyrBywaQYjHrtcDddCCFOSMKW4LHHHiM7O5sZM2YAMHPmTPbs2cO11147zC0TQggxHgWb9Fx+biKXn5uIw+1l2+E6cgt9RZTr21288XkZb3xeRqhZz6UzfMFr2dRYgoy64W66EEIEkLB1CqczQUaQQce+R7OGsFUnf+2BOHr0KK+88goFBQX+dZmZmezZs2eQWyaEEEIMnNmg49KZ8Vw6Mx63V2HH0QZyC6vIK6ymttXJ+vwK1udXEGTQcdE0XxHl5dPjCDNLEWUhxPCTa7b6aaxOkHHVVVfx7rvvotN1hzRFUUhKSqKsrGwYWyZGi9H8+RdCjF6KovJFaSO5BTZyC22UN3YXUTbqtCyeHE12ppVLZ8QTLUWUhRCDTK7ZEqe0YcMGPvnkE/Lz89Hruz8KO3fu5Oabb6a+vp7o6OhhbKEQQgjRN61Ww7y0KOalRfHTr82gsLLFH7yKatrYeLCWjQdr0Wq+YmF6tL+IsjVcfhQSQpw9ErbGKY/Hww9/+EPWrl3L7NmzA7Z1pfM9e/ZwySWXDEPrhBBCiP7TaDRkJoWTmRTOj7OmUVTT6g9eBRUtfHq0nk+P1vPQu4Wcl+oropydaWVCtBRRFkIMLQlb49Tzzz9PfX09d999d69tKSkpWCwWCVtCCCFGpclxodx9cSh3XzyFsoYO8gpt5BbY2F3aSH5pE/mlTTyRc4AZCWH+4DU1PkRqeQkhBp1cs3UKPSfIOHTo0Ji7ZkuIMyWffyHEaFHT4iBvXzV5BTY+PVqPt0cR5YkxwWRlWsnOsHKuFFEWQpyCFDUeZGN1ggwhzpR8/oUQo1Fju4sP91eTV2hjy+G6gCLKieFmf/CaJ0WUhRB9kAkyhBBCCCFOIDLYyHXzUrhuXgptTg8bD9SQW+grolzZ7OD/fVLC//ukhJgQIytm+oYaLpoYLUWUhRADImFLCCGEEONaiEnPFbMSuWKWr4jy1sN15BbY+HB/NXVtLt74vJQ3Pi8lrLOIcnamlaVTYzEPsLalEGL8kbAlhBBCCNHJbNCxYmY8KzqLKH92tJ7cAht5hdXUtTn5Z34F/+wsorx8eizZmQksnxZLqBRRFkL0QcKWEEIIIUQfDDotF06J5cIpsTx6VWZ3EeUCGxVNdt7/ysb7X9kw6rQsmRJDdqaVFTPiiQw2DnfThRAjhIStU+g5G6EQQgghxiedVsP8tCjmp0Xx31+bQUFFC7mFVeQU2Dha287HB2r4+EANOq2GhelRrMq0sjLDSnyYTBwkxHgmsxH2k8xGKETf5PMvhBjvDld3F1EurGwJ2DYnNYJVmQlkZVhJjbYMUwuFEINNZiMUQgghhDgLpsSHMiU+lO9fMoXSel8R5ZyCKr4obfIvj72/n5kJYazK9M1sODlOiigLMR5I2BJCCCGEGCSp0RZuWzqR25ZOpLrFwQeFNnIKbOwobmBfVQv7qlp4esMhJsYG+4JXRgKZSWESvIQYo2QYYT/JMEIh+iaffyGEOLWGziLKuQU2th2uw+XtLqKcFBFEdmeP15zUSCmiLMQoIMMIxZBKS0ujpKRkuJshhBBCjApRwUaun5fC9fNSaHW42XiwltyCKjYeqKWiyc4fthXzh23FxISYyMrw1fI6f2I0Bp0UURZiNJOwJYQQQghxFoWaDVw5K5ErO4sobzlUS26BjQ37fbW8/rKjlL/sKCU8yMAlM+JYlZnAhVNipIiyEKOQhC0xIHfeeSefffYZlZWVzJ49myuvvJJHH310uJslhBBCjEpmg46VGb5p4l0eXxHlnAIbG/bZqGtz8c8vKvjnFxVYjDqWT4sjO9PK8ulxhJjkTzghRgO5ZusUetbZOnTo0Ji+Zis/P58FCxZw/vnns3Xr1pPuK8MIRZex8vkXQoiRxKuo7D7WSE5BFXkFNiqbHf5tRr2WCyf7iihfKkWUhRgW/b1mS8JWP42HCTIuvPBClixZwrp162hubj7pzEgStkSXsfL5F0KIkUpVVb6qaCanwEZugY3iunb/Np1Ww/kTo8jOTCBrZjxxUkRZiLNCJsgQA/L6668TGRnJXXfdxZNPPsnRo0eZNGnScDdLCCGEGPc0Gg3nJkdwbnIE92VN43BNGzlf+Yoo769q4ZOiej4pqufn7xQwJzWSVZlWsjKspERJEWUhhpuELUF7ezsPPvggOTk5JCcnEx4ezp49eyRsCSGEECOMRqNhanwoU+ND+eGlUzhW305ugS945Zc2sftYI7uPNfLLf+8nMymM7IyuIsqhw910IcYlCVtDQVXB3TE8r22wwAALIz722GNkZ2czY8YMAGbOnMmePXu49tprh6KFQgghhBgkE6KDuWPZJO5YNglbs4O8Qt9Qwx3F9RRUtFBQ0cKvPzjE5LgQf/DKSJQiykKcLRK2hoK7Ax5PHJ7XfrASjMH93v3o0aO88sorFBQU+NdlZmayZ8+ekz5v5syZp9tCIYQQQgwBa7iZGy9I48YL0qhvc3YXUS6qo6imjRdqinhhYxHJkUH+4DUnNRKtFFEWYshI2BrnfvSjH1FfX09ycrJ/naIoJCUlnfR577///lA3TQghhBCnKTrExDfmp/KN+am0ONxsPFBDboGNTQdrKW+08/ttxfx+WzGxoZ1FlDMSWDgxSoooCzHIJGwNBYPF18M0XK/dTxs2bOCTTz4hPz8fvb77o7Bz505uvvlm6uvriY6OHopWCiGEEOIsCTMbuGp2ElfNTsLu8rL5UC15hTY+3F9NbauT//uslP/7rJQIi4FLZ8STnWFliRRRFmJQSNgaChrNgIbyDQePx8MPf/hD1q5dy+zZswO2dU1fuWfPHi655JJhaJ0QQgghhkKQUUd2pm8IocujsP1IHXmFNj4orKa+3cU/dpfzj93lBBt1LJ/eWUR5WhzBUkRZiNMi/+eMU88//zz19fXcfffdvbalpKRgsVgkbAkhhBBjmFGv5aJpcVw0LY5frlbZWdJAboGNvEIbVc0O3vuyive+rMKo17J0SizZmVZWzIgn3GIY7qYLMWpIUeN+Gg9FjYU4HfL5F0KIsUVVVfaWN/umlC+ooqS+e4ZlvVbDoknRZGVYWZkRT1yofO+L8am/RY0lbPWThC0h+iaffyGEGLtUVeVgdWtn8LJxwNbq36bRwLwJkWR1zmyYHClFlMX40d+wJcMIT2HdunWsW7cOr9c73E0RQgghhDirNBoN061hTLeGcc+lUymuayev0EZOgY29ZU3sLGlkZ4mviPI5SeH+68EmxYYMd9OFGBGkZ6ufpGdLiL7J518IIcanyiY7H3QGr50lDSg9/qKcEhfiD14zE6SIshh7pGdLCCGEEEIMmcSIIG5anM5Ni9Opa3Py4b5qcgpsbD9Sx+GaNg5/XMTzHxeREtVVRDmB81IipIiyGFckbAkhhBBCiDMSE2LimwtS+eaCVJrtviLKOQVVbD5US1mDnd9tLeZ3W4uJCzWRlWFlVaaVBelR6KWIshjjJGwJIYQQQohBEx5kYPV5Saw+L4kOl4cth2rJKbDx8f4aalqd/PmzY/z5s2NEWgysmBlPdqaVxZNjMOmliLIYeyRsCSGEEEKIIWEx6snOTCA7MwGnx8v2I/XkfmVjw/5qGtpd/G1XOX/bVU6ISc/FnUWUl02NlSLKYsyQT7IQQgghhBhyJr2O5dPiWD4tjse8CjtLGsktqCKvsBpbi4N391by7t5KTHoty6b6iihfMl2KKIvRTcKWEEIIIYQ4q/Q6LYsmRbNoUjQPXZHB3vImcgt8MxuWNnTwwb5qPthXjV6r4YLJMWRnWFkxM57YUNNwN12IAZGp3/tJpn4Xom/y+RdCCDFYVFXlgK2VnAIbeQU2DlYHFlGenxZFdoaVrEwrSRFBw9hSMd7J1O9CCCGEEGJU0Wg0zEgIY0ZCGPeumMrR2jZyC33Ba295M58XN/B5cQOPvrePWcnhZGVayc6wMlGKKIsRSnq2+kl6toTom3z+hRBCnA0VTXbyCmzkFvqKKPf8C3ZafKg/eM1ICJUiymLISc+WGBPS0tIoKSkZ7mYIIYQQYpglRQRx85J0bl6STm2rkw37qskttLG9qI6D1a0crG7lfz86zIRoi3+o4exkKaIshpeELSGEEEIIMarEhpr41sJUvrUwleYONx8dqCa3wMbmQ7Ucq+/g5S1HeXnLUaxhZrIy4snKtLIgTYooi7NPwpYYke68804+++wzKisrmT17NldeeSWPPvrocDdLCCGEECNMuMXANXOSuWZOMh0uD5sO1pJbYOPjAzXYWhy89ukxXvv0GFHBRlbM8BVRvmBytBRRFmfFuLlmq7W1lYsvvhi3243X6+UHP/gBt912W7+fPx6u2crPz2fBggWcf/75bN26dUS8lgwjHPnGyudfCCHE2OL0ePmkqI7cAhsb9lXT2OH2bws16bl4RhzZGVaWTYvFYpT+BzEw/b1ma9yELa/Xi9PpxGKx0NHRQWZmJjt37iQ6Orpfzx8PYevCCy9kyZIlrFu3jubm5n5dXHrRRRdx0003cdNNNw3Ja0nYGvnGyudfCCHE2OXxKnxe3OCb2bDQRnWL07/NbOguonzx9HjCg6SIsji1/oatcTNwVafTYbFYAN8fh16vl3GSM/vl9ddfJzIykrvuuovW1laOHj06Jl5LCCGEEEKv03LB5BgevSqTT++/hLe+ewG3L51ISlQQDrdCXmE1P3pzL/N+uYEb//g5b3xeSl2b89QHFuIURkyf6ZYtW3jqqafYvXs3VVVVrF+/ntWrVwfs8+KLL/LUU09RVVVFRkYGzz33HBdeeGG/X6OpqYlly5Zx+PBhnnrqKWJiYgb5XfioqordYx+SY59KkD5owNOdtre38+CDD5KTk0NycjLh4eHs2bOHSZMmDXr7zuZrCSGEEEIcT6vVMHdCJHMnRPLAqunsq2rxTyl/qLqNzYdq2Xyolp+u/8pXRDnTSnamlYRwKaIsBm7EhK329nZmzZrFd77zHa699tpe2998803uueceXnzxRRYvXszLL7/MqlWr2LdvH6mpqQDMnTsXp7P3rxAffPABiYmJREREsHfvXqqrq7nmmmv4+te/Tnx8/KC/F7vHzsLXFw76cftjx7d2YDFYBvScxx57jOzsbGbMmAHAzJkz2bNnT5//Hc7U2XwtIYQQQoiT0Wg0ZCSGk5EYzr0rp3Gkto3cAhu5BTa+qmhmR3EDO4obeORf+5iVEsGqzlpeaTHBw910MUqMyGu2NBpNr56thQsXMmfOHF566SX/uhkzZrB69WqeeOKJAb/Gd7/7XS6++GKuu+66Prc7nc6A4NbS0kJKSkq/rtnqcHeMmrB19OhRFixYQEFBAVarFYDbb7+dqqoq/vWvf/Xa//HHH+fxxx/3P7bb7RgMBvT67tyek5PTZ4/jQF8L4LLLLuP999/v9/sRZ59csyWEEGIsKm/sIK+wmtyCKnYdawwoojzdGurv8ZoWL0WUx6MxVdTY5XKxe/du7r///oD1K1euZPv27f06RnV1NUFBQYSFhdHS0sKWLVv47ne/e8L9n3jiCR555JHTam+QPogd39pxWs89U0H6gXVx/+hHP6K+vp7k5GT/OkVRSEpK6nP/O++8k+uvv97/+Nvf/jbXXnst11xzjX/diZ470NcCJGgJIYQQYlgkR1q4ZUk6tyxJp6bV4SuiXGDj0yP1HLC1csDWynMfHiYt2kJ2ZgLZmVZmJYdL8BIBRkXYqqurw+v19hryFx8fj81m69cxysvLueWWW1BVFVVVufvuuzn33HNPuP8DDzzAvffe63/c1bPVHxqNZsBD+YbDhg0b+OSTT8jPzw/omdq5cyc333wz9fX1vWZrjIqKIioqyv84KCiIuLg4Jk+ePOivJYQQQggxEsSFmvn2wgl8e+EEmjpcfLS/hpwCG1sO11JS38FvNx/ht5uPkBBuJivD1+M1Py0KnVaC13g3KsJWl+N/KVBVtd+/HsydO5c9e/b0+7VMJhMmk2kgzRtVPB4PP/zhD1m7di2zZ88O2NbVFbpnzx4uueSSUfVaQgghhBBDKcJi5Nq5yVw7N5l2p6+Ick5BFRsP1FDV7ODV7SW8ur2E6GAjKzPiycqwcsGkGIz6cTMJuOhhVIStmJgYdDpdr16smpqaIZngoqd169axbt06vF7vkL7O2fb8889TX1/P3Xff3WtbSkoKFotl0ALQ2XwtIYQQQoizJdik52vnJvC1cxNwuH1FlHMKbHy4v5r6dhdvfF7GG5+XEWrWc+kMX/BaNjWWIKNuuJsuzpJRNUHG3LlzefHFF/3rZs6cyVVXXXVaE2QM1HgoaizE6ZDPvxBCCBHI3VlEOaegirzCampbuyddCzLouGiar4jy8ulxhJmliPJoNOomyGhra6OoqMj/uLi4mD179hAVFUVqair33nsva9asYd68eSxatIhXXnmF0tJS7rzzzmFstRBCCCGEEIEMOi2LJ8eweHIMj16ZSX5ZIzlf+Wp5lTfaySmwkVNgw6jTsnhyNNmZVi6dEU90yNi9hGW8GjFha9euXSxfvtz/uGtyihtvvJFXX32Vb3zjG9TX1/Poo49SVVVFZmYm77//PhMmTBjSdo3VYYRCCCGEEGLo+YooRzF3QhQ//doMCitbfLW8Cm0U1bSx8WAtGw/WotV8xcJ0X/DKyrBiDZfRImPBiBxGOBLJMEIh+iaffyGEEOL0FNW0+oNXQUVLwLbzUiPI7pzZcEK0FFEeafo7jFDCVj9J2BKib/L5F0IIIc5cWUMHeYU2cgts7C4NLKI8IyHMH7ymxodILa8RYNRdsyWEEEIIIcR4lRJl4dYLJ3LrhROpaXGQt6+avAIbnx6tZ39VC/urWnj2w0NMjAkmK9NKdoaVc6WI8ognYesU5JotIYQQQghxNsWFmVlz/gTWnO8rorxhXzV5hTa2HK7jaF07L206wkubjpAYbvYHr3lSRHlEkmGE/STDCIXom3z+hRBCiLOjzelh44EacgttbDxQQ4eruzMgJsTIiplWVmVaOX9itBRRHmIyjFAIIYQQQogxJMSk54pZiVwxKxGH28u2w91FlOvaXLzxeSlvfF5KWGcR5exMK0unxmI2SBHl4SJhSwghhBBCiFHGbNBx6cx4Lp0Zj9ursONodxHlujYn/8yv4J/5FQQZdCyfHkt2ZgLLp8USKkWUzyoJW6cg12wJIYQQQoiRzKDTsmRKDEumxPDoVZnklzaSU+Cb2bCiyc77X9l4/ytfEeUlU2LIzrSyYkY8kcHG4W76mCfXbPWTXLMF11xzDQ0NDWzatClg/U033cTvf/979Pru7H706FEKCwu54oorhrwNQ/Vaon/Gy+dfCCGEGG1UVaWwsoWcgipyCmwcrW33b9NpNSxMj2JVppWVGVbiw+Tf8IHo7zVbcuWc6Jf8/HxsNhsHDhwIWP/MM8+we/dufvrTn3Ls2DH/+pycnF77DlUbhuK1hBBCCCFGO41GQ2ZSOGuzpvPxf13Ehh8t5b9WTCUjMQyvorL9SD0/e6eQhY9/xDUvfsLvthyltL5juJs9pkjPVj+N956tyy+/nPvuu4/vfOc77Nixg5iYGP+243u2Nm/ezOrVq4mNjSU4OJjt27cTFBQ0JG0YqtcS/TcePv9CCCHEWFNa31lEudDG7mONAdtmJoSxKtNXRHlynBRR7ov0bIlBs2PHDpxOJ0uXLmXGjBns27fvpPsvW7aMzMxMPvroI/Lz8wcl/JyoDWf6WqqqcvvttxMVFYVGo2HPnj1n3FYhhBBCiJEuNdrCbUsn8tZ3L2DHg5fwi6syWDw5Gp1Ww76qFp7ecIgVz27hkmc281TeAb4qb0b6aAZOwtYprFu3jpkzZzJ//vzhbsqw+dnPfsajjz4K0K+wBVBeXk5KSspZacOZvFZubi6vvvoq7733HlVVVWRmZg5Ke7vcc889rF69elCPORAvvviiv8dp7ty5bN269aT7v/TSS5x77rmEhYURFhbGokWLyMnJOUutFUIIIcRwiA8zs2ZRGn+59Xx2/vRSfvX1c7lkehxGnZajte2s23iEK17YxpL/2cgv3tvHzpIGvIoEr/6Q2QhP4a677uKuu+7ydxX2h6qqqHb7ELesb5qgoEHt6t26dSvbtm3jhhtuAKCtrY1vfetbAft87WtfQ6vtzu3l5eUkJSUF7LN//3527NjBRRddRHx8/IB6oE7WhjN9rSNHjpCQkMAFF1zQ7/Ycz+VyYTT2PZvPzp07+drXvnbaxz4Tb775Jvfccw8vvvgiixcv5uWXX2bVqlXs27eP1NTUPp+TnJzMk08+yeTJkwF47bXXuOqqq8jPzycjI+NsNl8IIYQQwyAq2Mj181K4fl4KrQ43Gw/WkldgY+PBGiqa7PxhWzF/2FZMTIiJrAxfLa/zJ0Zj0EkfTl8kbA0B1W7n4Jy5w/La077YjcZiGdBzysrKuP/++3n//fcBWLVqFevWrSMyMpKf//znfPjhh/4w8sknn/DQQw8FPP+6664LeFxcXExiYmLAOkVReOGFF2hra+Puu+8etDaczmt1uemmm3jttdcA3wWkEyZMoKSkBKfTydq1a/nrX/9KS0sL8+bN49lnn/X3bl500UVkZmZiNBr505/+REZGBps3bw44ttvtJjg4GLfbzfbt2/npT3/KggUL2LFjxwnbM9ieeeYZbrnlFm699VYAnnvuOfLy8njppZd44okn+nzO8TM6PvbYY7z00kt89tlnEraEEEKIcSbUbODKWYlc2VlEecuhWnILbXy4z1fL6y87SvnLjlLCgwxcMiOOVZkJXDglRooo9yARdJwrKipi7ty5TJo0iU8//ZQPP/yQI0eOsHbtWj766CM0Gk1Ar8+UKVNOOYwwMzOTw4cPc8455/hnCTx27Bg33ngjEyZMoKMjcJabM2nDQF+rp9/85jc8+uijJCcnU1VVxc6dOwG47777eOutt3jttdf44osvmDx5MllZWTQ0NPif+9prr6HX6/nkk094+eWXex1bp9Oxbds2APbs2UNVVRV5eXknPW99efzxxwkJCTnp0tfQQJfLxe7du1m5cmXA+pUrV7J9+/Z+vbbX6+Wvf/0r7e3tLFq0aMBtF0IIIcTYYTboWJlh5ZnrZ7Prv1fwp5sX8K2FqcSEGGm2u/nnFxXc9qddzPnFBu76yxf8a28lbU7PcDd72MlshP00kNkIR9MwwksvvZTFixfzyCOP+Ne99dZbrF27lqNHjw5q22w2G1ar9ay04USvdbznnnuO5557jpKSEgDa29uJjIzk1Vdf9Q9VdLvdpKWlcc8997B27Vouuugimpubyc/PP+mx3377bW699Vbq6up6bVNVFY1Gw8MPP8zDDz/sf3y8hoaGgJDXl6SkpF5DJSsrK0lKSuKTTz4JCKqPP/44r732GgcPHjzh8b766isWLVqEw+EgJCSE119/ncsuu+yE+8tshEIIIcT45VVUdh9rJLfARl6hr4hyF6Ney4WTfUWULx1jRZT7OxuhDCMcAhqNZsBD+YbDsWPH+Oijj9i+fTtPP/20f73X6x3UyS269BV+hqoN/QlafTly5Ahut5vFixf71xkMBhYsWMD+/fv96+bNm3fKY+Xn5zNr1qw+t7300kvo9Xra29u5//77WbVqFcuWLeu1X1RUFFFRUafxTnyOD3AnCnU9TZs2jT179tDU1MRbb73FjTfeyObNm5k5c+Zpt0MIIYQQY5NOq2FBehQL0qP42eUz+KqimdwCG7kFNo7WtfPRgRo+OlCDTqth0cRosjKtZM2MJ26cFFGWsHUK69atY926dXi93uFuyqDbu3cvUVFRfV5HdLZqVY2ENvTU1dF7qpASHBx8ymPt2bPnhGHre9/7Hk899RT/+7//y8cffxwQ7np6/PHHefzxx0/6Ojk5OVx44YUB62JiYtDpdNhstoD1NTU1xMfHn/R4RqPRP0HGvHnz2LlzJ7/5zW/6HC4phBBCCNFFo9FwbnIE5yZHsDZrGodr2sgtsJFTYGN/VQvbiurYVlTHz98pYG5qJNmZVrIyrKREjfxOitMlYesUTmc2wtHCYDDQ2tpKQkJCv8LDWG1DT5MnT8ZoNLJt27aAYYS7du3innvuGdCxvvrqK66++uo+t/32t78lPDycH/zgB7z33nsoitIrMAHceeedXH/99Sd9neNnYwRfYJo7dy4bNmwIaMOGDRu46qqrBvQ+VFXF6XQO6DlCCCGEGN80Gg1T40OZGh/KDy6ZwrH6dvIKfcErv7SJXcca2XWskV/+ez+ZSWFkZ3QVUQ4d7qYPKglb49jChQsJCwtjzZo1/PznPyckJISioiJycnL4zW9+M27a0FNwcDDf/e53Wbt2LVFRUaSmpvKrX/2Kjo4ObrnllgEdS1EUvvzySyorKwkODg4I63fccUeva7b6cibDCO+9917WrFnDvHnzWLRoEa+88gqlpaXceeed/n1eeOEF1q9fz0cffQTAgw8+yKpVq0hJSaG1tZW//vWvbNq0idzc3NNqgxBCCCEEwIToYG5fOonbl07C1uzgg302cr6ysaO4noKKFgoqWvj1B4eYHBfiD14ZiWGDWtJoOEjYGseioqJ4//33+clPfsKyZctQVZXJkyezZs2acdWG4z355JMoisKaNWtobW1l3rx55OXlERkZOaDj/PKXv+QnP/kJzz77LPfee2/ANWldXxwPP/xwwOPB9I1vfIP6+noeffRRf8Hm999/nwkTJvj3qaur48iRI/7H1dXVrFmzhqqqKsLDwzn33HPJzc1lxYoVg94+IYQQQoxP1nAz/7kojf9clEZDu4sP91WTU1DFtqI6imraeKGmiBc2FpEcGeQPXnNSI9FqR1/wktkI+2kgsxEKMZ7I518IIYQQg6HF4WbjgRpyC2xsOliL3d09Z0JsqK+I8tfnpjA7JWL4GtlJZiMUQgghhBBCjBphZgNXzU7iqtlJ2F1ethyuJbfAxof7q6ltdfJ/n5WSGBE0IsJWf0nYEkIIIYQQQowoQUYdWRm+2QpdHoVPj9aTW1DFqsyE4W7agEjYEkIIIYQQQoxYRr2WZVNjWTY1dribMmDa4W7AWCKXv4nxSD73QgghhBB9k7B1CuvWrWPmzJnMnz//hPvodDoAXC7X2WqWECNGR0cH4KuZJoQQQgghuslshP10shlHVFWltLQUt9tNYmIiWq1kWDH2qapKR0cHNTU1REREkJAwusZQCyGEEEKcLpmN8CzSaDQkJCRQXFzMsWPHhrs5QpxVERERWK3W4W6GEEIIIcSII2FrkBiNRqZMmSJDCcW4YjAY/MNohRBCCCFEIAlbg0ir1UpRVyGEEEIIIQQgE2QIIYQQQgghxJCQsCWEEEIIIYQQQ0DClhBCCCGEEEIMAblmq5+6ZshvaWkZ5pYIIYQQQgghhlNXJjhVFS0JW/3U2toKQEpKyjC3RAghhBBCCDEStLa2Eh4efsLtUtS4nxRFobKyktDQUDQazbC2paWlhZSUFMrKyk5aRE2cHjm/Q0vO79CS8zu05PwOLTm/Q0vO79CS8zu0Rtr5VVWV1tZWEhMT0WpPfGWW9Gz1k1arJTk5ebibESAsLGxEfNjGKjm/Q0vO79CS8zu05PwOLTm/Q0vO79CS8zu0RtL5PVmPVheZIEMIIYQQQgghhoCELSGEEEIIIYQYAhK2RiGTycRDDz2EyWQa7qaMSXJ+h5ac36El53doyfkdWnJ+h5ac36El53dojdbzKxNkCCGEEEIIIcQQkJ4tIYQQQgghhBgCEraEEEIIIYQQYghI2BJCCCGEEEKIISBhSwghhBBCCCGGgIStEeDFF18kPT0ds9nM3Llz2bp160n337x5M3PnzsVsNjNx4kR++9vf9trnrbfeYubMmZhMJmbOnMn69euHqvkj3kDO7z//+U9WrFhBbGwsYWFhLFq0iLy8vIB9Xn31VTQaTa/F4XAM9VsZkQZyfjdt2tTnuTtw4EDAfvL57TaQ83vTTTf1eX4zMjL8+8jnt9uWLVu44oorSExMRKPR8Pbbb5/yOfL9238DPb/y/TswAz2/8v07MAM9v/L9OzBPPPEE8+fPJzQ0lLi4OFavXs3BgwdP+bzR+B0sYWuYvfnmm9xzzz389Kc/JT8/nwsvvJBVq1ZRWlra5/7FxcVcdtllXHjhheTn5/Pggw/ygx/8gLfeesu/z6effso3vvEN1qxZw969e1mzZg3XX389O3bsOFtva8QY6PndsmULK1as4P3332f37t0sX76cK664gvz8/ID9wsLCqKqqCljMZvPZeEsjykDPb5eDBw8GnLspU6b4t8nnt9tAz+9vfvObgPNaVlZGVFQU1113XcB+8vn1aW9vZ9asWbzwwgv92l++fwdmoOdXvn8HZqDnt4t8//bPQM+vfP8OzObNm7nrrrv47LPP2LBhAx6Ph5UrV9Le3n7C54za72BVDKsFCxaod955Z8C66dOnq/fff3+f+993333q9OnTA9bdcccd6vnnn+9/fP3116vZ2dkB+2RlZanf/OY3B6nVo8dAz29fZs6cqT7yyCP+x//v//0/NTw8fLCaOKoN9Pxu3LhRBdTGxsYTHlM+v93O9PO7fv16VaPRqCUlJf518vntG6CuX7/+pPvI9+/p68/57Yt8//ZPf86vfP+evtP5/Mr378DU1NSogLp58+YT7jNav4OlZ2sYuVwudu/ezcqVKwPWr1y5ku3bt/f5nE8//bTX/llZWezatQu3233SfU50zLHqdM7v8RRFobW1laioqID1bW1tTJgwgeTkZC6//PJev7yOB2dyfs877zwSEhK45JJL2LhxY8A2+fz6DMbn9w9/+AOXXnopEyZMCFgvn9/TI9+/Z5d8/w4N+f49O+T7d2Cam5sBev3/3tNo/Q6WsDWM6urq8Hq9xMfHB6yPj4/HZrP1+Rybzdbn/h6Ph7q6upPuc6JjjlWnc36P9/TTT9Pe3s7111/vXzd9+nReffVV3n33Xd544w3MZjOLFy/m8OHDg9r+ke50zm9CQgKvvPIKb731Fv/85z+ZNm0al1xyCVu2bPHvI59fnzP9/FZVVZGTk8Ott94asF4+v6dPvn/PLvn+HVzy/Xv2yPfvwKiqyr333suSJUvIzMw84X6j9TtYP2yv/P/bu7+Xpv44juNv3dlB88LE0sQww4t1VWTlZobehDfe+AckFl4JButul0O6MIiCQoJA62plPykIIqEUVBpenGFhEBVEoODVJOhG6N1FX4/fNX/sbJ6dbT4fsIsd3n722ZvP3vDa2IStrKws5b6qpl3bqf7f607XLGXZ9uLBgwcSjUblxYsXUldXZ18PhUISCoXs+x0dHdLa2iq3b9+WW7du7d7Gi4ST/gYCAQkEAvb99vZ2+fHjh1y/fl06OzuzWrPUZduL+/fvy/79+6W3tzflOuc3N8zf/GD+7j7mb/4wf50ZGhqShYUFmZmZ2bG2GGcwn2x56MCBA+Lz+dLS9srKSloqX3fo0KFN6w3DkNra2m1rtlqzVGXT33UTExMyMDAgjx49kvPnz29bW15eLmfOnNlz70zl0t//C4VCKb3j/P6VS39VVcbHx6Wvr09M09y2dq+e32wwf/OD+Zs/zN/dx/x15vLly/Ly5Ut59+6dHD58eNvaYp3BhC0PmaYpp06dksnJyZTrk5OTcvbs2U3/pr29Pa3+zZs3cvr0afH7/dvWbLVmqcqmvyJ/31G9ePGixGIx6enp2fFxVFUSiYQ0NDTkvOdikm1//2VZVkrvOL9/5dLf6elp+fLliwwMDOz4OHv1/GaD+es+5m9+MX93H/M3M6oqQ0ND8uzZM3n79q0cPXp0x78p2hmc39/jwL8ePnyofr9fx8bGdHFxUcPhsFZVVdm/XhOJRLSvr8+u//btm+7bt0+vXLmii4uLOjY2pn6/X588eWLXzM7Oqs/n05GREf306ZOOjIyoYRj6/v37vD8/rzntbywWU8MwdHR0VJeXl+1bMpm0a6LRqL5+/Vq/fv2qlmXppUuX1DAMjcfjeX9+XnPa35s3b+rz58/18+fP+vHjR41EIioi+vTpU7uG87vBaX/XXbhwQYPB4KZrcn43/Pz5Uy3LUsuyVET0xo0balmWfv/+XVWZv7ly2l/mrzNO+8v8dcZpf9cxfzMzODio1dXVOjU1lfJ6//Xrl11TKjOYsFUARkdH9ciRI2qapra2tqb87GV/f792dXWl1E9NTenJkyfVNE1tbm7WO3fupK35+PFjDQQC6vf79dixYynDdK9x0t+uri4VkbRbf3+/XRMOh7WpqUlN09SDBw9qd3e3zs3N5fEZFRYn/b127Zq2tLRoRUWF1tTU6Llz5/TVq1dpa3J+NzidD8lkUisrK/Xu3bubrsf53bD+U9hbvd6Zv7lx2l/mrzNO+8v8dSab+cD8zdxmvRURvXfvnl1TKjO4TPW/b5YBAAAAAHYN39kCAAAAABcQtgAAAADABYQtAAAAAHABYQsAAAAAXEDYAgAAAAAXELYAAAAAwAWELQAAAABwAWELAAAAAFxA2AIAAAAAFxC2AADIUDgclt7eXq+3AQAoEoQtAAAyND8/L21tbV5vAwBQJMpUVb3eBAAAhWxtbU2qqqpkbW3NvtbW1ibxeNzDXQEACp3h9QYAACh0Pp9PZmZmJBgMSiKRkPr6eqmoqPB6WwCAAkfYAgBgB+Xl5bK0tCS1tbVy4sQJr7cDACgSfGcLAIAMWJZF0AIAOELYAgAgA4lEgrAFAHCEsAUAQAY+fPggx48f93obAIAiQtgCACADv3//loWFBVlaWpLV1VWvtwMAKAKELQAAMnD16lWZmJiQxsZGGR4e9no7AIAiwP/ZAgAAAAAX8MkWAAAAALiAsAUAAAAALiBsAQAAAIALCFsAAAAA4ALCFgAAAAC4gLAFAAAAAC4gbAEAAACACwhbAAAAAOACwhYAAAAAuICwBQAAAAAuIGwBAAAAgAv+AGzs4XxBEnhsAAAAAElFTkSuQmCC", "text/plain": [ "Figure(PyObject