{ "cells": [ { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "## Exercises: An Invitation to Analytic Combinatorics in Several Variables\n", "\n", "Created by Stephen Melczer\n", "\n", "This notebook complements the exercises found on the [ALEA ACSV course webpage](https://melczer.ca/ALEA22/)\n", "\n", "A quick Sage tutorial can be found [here](https://melczer.ca/files/SageIntro.html) (try an interactive version in your browser [here](https://mybinder.org/v2/git/https%3A%2F%2Fgit.uwaterloo.ca%2Fsmelczer%2Fintro-to-sage/HEAD?filepath=A%20Brief%20Introduction%20to%20Sage.ipynb)).\n", "\n", "See [https://melczer.ca/textbook](https://melczer.ca/textbook) for further Sage notebooks solving problems in analytic combinatorics in several varibles. In particular, [this notebook](https://melczer.ca/files/TextbookCode/Chapter5/Example5-SmoothASM.ipynb) (also available as a [static HTML page](https://melczer.ca/files/TextbookCode/Chapter5/Example5-SmoothASM.html)) gives an algorithm to compute asymptotic terms. Don't use it to solve these exercises, but use it to further check your answers! (Or to compute asymptotics for other problems!)" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "# Helper function to compute the Hessian matrix M from the function H(vars)\n", "# in the direction R at the point CP, specified by a list of substitutions.\n", "# Copied from melczer.ca/textbook/\n", "def getHes(H,R,vars,CP):\n", " dd = len(vars)\n", " V = zero_vector(SR,dd)\n", " U = matrix(SR,dd)\n", " M = matrix(SR,dd-1)\n", "\n", " for j in range(dd):\n", " V[j] = R[j]/R[-1]\n", " for i in range(dd):\n", " U[i,j] = vars[i]*vars[j]*diff(H,vars[i],vars[j])/vars[-1]/diff(H,vars[-1])\n", " for i in range(dd-1):\n", " for j in range(dd-1):\n", " M[i,j] = V[i]*V[j] + U[i,j] - V[j]*U[i,-1] - V[i]*U[j,-1] + V[i]*V[j]*U[-1,-1]\n", " if i == j: M[i,j] = M[i,j] + V[i]\n", " return M.subs(CP)\n", "\n", "\n", "# Helper function to compute leading asymptotics of the R-diagonal of G(vars)/H(vars)\n", "# determined by the Main Asymptotic Theorem of Smooth ACSV at the point CP, specified\n", "# by a list of substitutions. We take det(M) as an input that can be computed by the\n", "# above function.\n", "var('n')\n", "def leadingASM(G,H,detM,R,vars,CP):\n", " dd = len(R)\n", " lcoeff = -G/vars[-1]/H.diff(vars[-1])\n", " exp = 1/mul([vars[k]^R[k] for k in range(dd)])^n\n", "\n", " ASM = exp * (2*pi*n*R[-1])^((1-dd)/2) / sqrt(detM) * lcoeff\n", " return ASM.subs(CP)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "These functions can be used to compute the matrix M appearing in asymptotics, as well as the leading asymptotic term in an asymptotic expansion. \n", "\n", "Here is an example of their use to find asymptotics for the main diagonal of $1/(1-x-y-z)$." ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left[\\left[x = \\left(\\frac{1}{3}\\right), y = \\left(\\frac{1}{3}\\right), z = \\left(\\frac{1}{3}\\right)\\right]\\right]\\]" ], "text/latex": [ "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left[\\left[x = \\left(\\frac{1}{3}\\right), y = \\left(\\frac{1}{3}\\right), z = \\left(\\frac{1}{3}\\right)\\right]\\right]$$" ], "text/plain": [ "[[x == (1/3), y == (1/3), z == (1/3)]]" ] }, "execution_count": 2, "metadata": { }, "output_type": "execute_result" } ], "source": [ "# Introduce variables x, y, and z\n", "var('x y z')\n", "H = 1 - x - y - z\n", "\n", "# In the main diagonal direction this has a critical point at (1/3,1/3,1/3)\n", "CPs = solve([H,x*diff(H,x) - y*diff(H,y), x*diff(H,x) - z*diff(H,z)],[x,y,z])\n", "show(CPs)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\verb|M|\\verb| |\\verb|=| \\left(\\begin{array}{rr}\n", "2 & 1 \\\\\n", "1 & 2\n", "\\end{array}\\right)\\]" ], "text/latex": [ "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\verb|M|\\verb| |\\verb|=| \\left(\\begin{array}{rr}\n", "2 & 1 \\\\\n", "1 & 2\n", "\\end{array}\\right)$$" ], "text/plain": [ "'M = ' [2 1]\n", "[1 2]" ] }, "execution_count": 3, "metadata": { }, "output_type": "execute_result" } ], "source": [ "# Let CP be the critical point, defined as a list [x == 1/3, y == 1/3, z == 1/3]\n", "CP = CPs[0]\n", "\n", "# Get the matrix M \n", "M = getHes(H,[1,1,1],[x,y,z],CP)\n", "show(\"M = \", M)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The dominant asymptotic behaviour of the main diagonal is\n" ] }, { "data": { "text/html": [ "\\[\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{\\sqrt{3}}{2 \\, \\pi \\left(\\frac{1}{27}\\right)^{n} n}\\]" ], "text/latex": [ "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\frac{\\sqrt{3}}{2 \\, \\pi \\left(\\frac{1}{27}\\right)^{n} n}$$" ], "text/plain": [ "1/2*sqrt(3)/(pi*(1/27)^n*n)" ] }, "execution_count": 4, "metadata": { }, "output_type": "execute_result" } ], "source": [ "# Get and print leading asymptotics\n", "ASM = leadingASM(1,H,M.determinant(),[1,1,1],[x,y,z],CP)\n", "print(\"The dominant asymptotic behaviour of the main diagonal is\")\n", "show(ASM)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "We can check our asymptotic approximation by computing series terms." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "0.978030899006393" ] }, "execution_count": 5, "metadata": { }, "output_type": "execute_result" } ], "source": [ "# First, define the ring of formal power series (more efficient for computations)\n", "# We use capital letters to denote the formal power series variables\n", "P. = QQ[['X,Y,Z']] \n", "\n", "# Computes the series expansion up to precision 3*N\n", "N = 10\n", "ser = 1/(1-X-Y-Z + O(X^(3*N+1)))\n", "\n", "# Check ratio of asymptotic formula to actual coefficients -- this should go to 1!\n", "(ser.coefficients()[X^N*Y^N*Z^N]/ASM.subs(n=N)).n()" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "## Question 1: Delannoy Numbers\n", "The *Delannoy number* $d_{a,b}$ is the number of paths from the origin $(0,0)$ to the point $(a,b)$ using only the steps $\\textsf{N}=(0,1)$, $\\textsf{E} = (1,0)$, and $\\textsf{NE}=(1,1)$.\n", "\n", "**(a)** Prove the recurrence \n", "$$ d_{a,b} = \n", "\\begin{cases} \n", "1 &: \\text{ if $a=0$ or $b=0$} \\\\\n", "d_{a-1,b} + d_{a,b-1} + d_{a-1,b-1} &:\\text{ otherwise}\n", "\\end{cases}\n", "$$\n", "Conclude that \n", "$$ D(x,y) = \\sum_{a,b\\geq0}d_{a,b}x^ay^b = \\frac{1}{1-x-y-xy}. $$\n", "\n", "**(b)** Use the Main Theorem of Smooth ACSV to find asymptotics of $d_{n,n}$ as the $(1,1)$-diagonal of $D(x,y)$. What are the critical points in the $(1,1)$ direction? Which are minimal?\n", "\n", "**(c)** Use the Main Theorem of Smooth ACSV to find asymptotics of the $(r,s)$-diagonal of $D(x,y)$ for any $r,s>0$." ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Function to numerically compute terms in the expansion\n", "\n", "Check your computed asymptotics against this function!" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The coefficient [x^10y^20]D(x,y) = 4354393801\n" ] } ], "source": [ "# First, define the ring of formal power series (more efficient for computations)\n", "P. = QQ[['X,Y']] \n", "\n", "# Code to compute the coefficient of x^(N*R) * y^(N*S) where R, S, and N are positive integers\n", "R, S, N = 1, 2, 10\n", "N2 = (R+S)*N\n", "ser = 1/(1-X-Y-X*Y + O(X^(N2+1)))\n", "coef = ser.coefficients()[X^(R*N)*Y^(S*N)]\n", "print(\"The coefficient [x^{}y^{}]D(x,y) = {}\".format(N*R,N*S,coef))" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "## Question 2: Apéry Asymptotics\n", "\n", "Recall from lecture that a key step in Apéry's proof of the irrationality of $\\zeta(3)$ is determining the exponential growth of the sequence that can be encoded as the main diagonal\n", "of \n", "$$ F(x,y,z,t) = \\frac{1}{1 - t(1+x)(1+y)(1+z)(1+y+z+yz+xyz)}. $$\n", "Use the Main Theorem of Smooth ACSV to find dominant asymptotics of this sequence." ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Function to numerically compute terms in the expansion\n", "\n", "Check your computed asymptotics against this function!" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The coefficient [(xyzt)^(30)]F(x,y,z,t) = 11320115195385966907843180411829810312080825\n" ] } ], "source": [ "# Numerically compute coefficient of (xyzt)^n (can take a long time for large N)\n", "P. = QQ['X,Y,Z'] \n", "N = 30\n", "ser = ((1+X)*(1+Y)*(1+Z)*(1+Y+Z+Y*Z+X*Y*Z))^N\n", "coef = ser[X^N*Y^N*Z^N]\n", "print(\"The coefficient [(xyzt)^({})]F(x,y,z,t) = {}\".format(N,coef))" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "## Question 3: Pathological Directions\n", "**(a)** Find asymptotics of the $(r,s)$-diagonal of $F(x,y) = \\frac{1}{1-x-xy}$\n", "for any $0 = QQ[['X,Y']] \n", "\n", "# Code to compute the coefficient of x^(N*R) * y^(N*S) where R, S, and N are positive integers\n", "R, S, N = 2, 1, 10\n", "N2 = (R+S)*N\n", "ser = 1/(1-X-X*Y + O(X^(N2+1)))\n", "coef = ser.coefficients()[X^(R*N)*Y^(S*N)]\n", "print(\"The coefficient [x^{}y^{}]F(x,y) = {}\".format(N*R,N*S,coef))" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "## Question 4: A Composition Limit Theorem\n", "\n", "An *integer composition* of size $n\\in\\mathbb{N}$ is an ordered tuple of positive integers (of any length) that sum to $n$. A composition can be viewed as an integer partition where the order of the summands matters. Let $c_{k,n}$ denote the number of compositions of size $n$ that contain $k$ ones. \n", "\n", "**(a)** If you know the symbolic method, species theory, or similar enumerative constructions, prove that \n", "$$ C(u,x) = \\sum_{n,k\\geq0}c_{k,n}u^kx^n = \\frac{1-x}{1-2x-(u-1)x(1-x)}. $$\n", "\n", "**(b)** Prove that the distribution for the number of ones in a composition of size $n$ satisfies a local central limit theorem as $n\\rightarrow\\infty$. More precisely, find a constant $t>0$ and normal density $\\nu_n(s)$ such that \n", "$$ \\sup_{s \\in \\mathbb{Z}} |t^nc_{s,n} - \\nu_n(s)| \\rightarrow 0 $$\n", "as $n\\rightarrow\\infty$." ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Function to plot series coefficients\n", "\n", "Check your computed distribution against this function!" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The following plot shows the coefficients of [x^(200)]C(u,x)\n" ] }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAksAAAGFCAYAAADtt7dbAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAva0lEQVR4nO3de3xcZZ348c/TW8otoUKBQlMot5ZAsVwiIJcWEbmJy8W7gFcUdVcQcfdXVldRMOuKFRdmraK7gqyCiKIoKIh1ASm0lEJbLi1IS1taSkvbpGKbtun5/fEkzaXJNDPJzJnMfN6v13l1njPnzPlmcmby7fM853tCkiRIkiSpe4PSDkCSJKmUmSxJkiRlYbIkSZKUhcmSJElSFiZLkiRJWZgsSZIkZWGyJEmSlIXJkiRJUhYDIlkKUXUIIaQdiyRJqixDUj5+r8qHNzY2UlNTQ2NjY6HjkSRJ5SnvDpcB0bMkSZKUlpJOljKZDHV1ddTX16cdiiRJqlAh5Rvp9urgTU1N24bhqqurCx2TJEkqP3kPw6U9Z0nqtZ/+FJ5+Gt7xDjjttLSjkSRVipIehpPafPvb8KEPwX/8B5x+Otx/f9oRSZIqRZ+SpRDClBBCEkK4YQfbTQohzA4hbAwhvBRCuKwvx1XlWL0abrkFfvzj9nVJAvfem1pIkqQKk/cwXAihHvgkMHcH240F7gVuBi4CTgT+K4SwKuX5Uipxr78O9fWwePH2zzU3w69+Be98JwwdWvTQJEkVJK+epRDCrsD/ApcCa3ew+WXAkiRJrkiS5LkkSX4I/DdwVT7HVuX4/e87J0pDh8bk6dBDYdo0uOCCmCxt3ZpaiJKkCpDvMFwG+F2SJH/sxbYnAF1nmPwBOHbz5s3d7tDc3ExTU1OnRZVnn306tw84AG6/HRYubF93//3w/PNFDUuSVGFyTpZCCO8Hjgam9HKXfYCVXdatBIasXr262x0aGhqoqanZttTW1uYapsrAaafBl78MI0bA+PHws59BTU3nYbchQ2D33VMLUZJUAXJKlkIItcB3gYuSJNmYw65dJyeF1tfrduMpU6bQ2Ni4bVm6dGkuYaqMfO1rsGYNPPccHHMM7LFHnOz9pjfFJOkHP4B99007SklSOcupKGUI4TzgV0BLh9WDicnQVqAqSZKWLvs8BMxJkuTyDuvOB36+adOmIUN7MTvXopSSJKmPilaU8kFgQpd1/wM8D3yza6LUagZwbpd17wCeGDp06PE5Hl+SJKmochqGS5JkfZIk8zsuwBvA662PCSE0hBBu7bDbNGD/EMLUEMJhIYSPAR8Hru+vH0Ll5emn4dhjYf/94TvfSTsaSVKlK0QF71HAmLZGkiSLgLOBycBTwJeBzyVJclcBjq0y8J73wOzZsGQJXHklzJiRffvVq6Gluz5NSZL6QZ+TpSRJJidJckWH9keSJJncZZv/S5Lk6CRJqpIkGZskybS+Hlfla8mSzu2e5vevXw8nnggjR0JtLczNWh5VkqT8eG84lZwPf7j98ejRcOqp3W93003w6KPx8YoV8IUvFD42SVLlyft2J1KhTJsGb3tbvN3JBRfEnqPubNiQvS1JUn/IqXRAAWQ9eCaTIZPJ0NLSwsKFCy0doE6WLYvDcEuWwPDh8V5xZ56ZdlSSpBKVd+mAkk6W2lhnST1pbIxXz40dG+ctSZLUg6LVWZJKSk0NnHJK2lFIksqZE7wlSZKyMFmSJEnKwmRJkiQpC5MlSZKkLEyWVDKWLYuX/y9YkHYkkiS1K+lkKZPJUFdXR319fdqhqMDmzYMjjohFKI88Eu67L+2IJEmKrLOkkvC5z8GNN7a3zzwz94Rp5UpYtQrGj4chFsWQJHWWd52lku5ZUuXYfffO7REjctv/jjtiUcoJE+KtUjZu7LfQJEkVzmRJJeGLX4TTToMQYOJE+OY3c9v/qqtg8+b4+OGH4c47+z1ESVKFcrBCJWG33eCPf4StW2FQHil8133yeQ1JkrrjnxSVlHyTnO9+N95MF+Dtb4f3vKf/YpIkVTZ7llQWzjsPli+HtWvhgAPsWZIk9R+TJZWNESNynxguSdKOlPT/v62zJEmS0madJUmSVAmssyRJklQIJkuSJElZmCxJkiRlYbIkSZKUhcmSJElSFiZLkiRJWZgsKXWLF8PSpWlHIUlS90o6WbIoZfn7p3+CsWNhzBj46lfTjkaSpO1ZlFKpef55OOywzuteew1Gjuzb627cCNOnx1ufHH98315LklQ2LEqpgae7PL2vufvGjTBpEpx9NpxwAlx9dd9eT5IkkyWl5rDD4B//sb39la/AXnv17TUffBBmzmxvf+tb0NLSt9eUJFW2nJKlEMKnQwhzQwhNrcuMEMJZWbafHEJIulnG9z10lYMbb4S//hWWLOmfOUs1NZ3bu+4Kgwf3/XUlSZUr156lZcD/A45tXf4E/DqEcPgO9hsHjOqwvJDjcVXGDjwQamv757VOOgmuuiomSDU1cNtt/fO6kqTK1ecJ3iGENcAXkyT5UTfPTQamAyOSJFnXze5O8FZBbN4MQ4ZAyHs6nySpzOT9F2FI3kcMYTDwHmAXYMYONp8TQhgOPAtcmyTJ9GwbNzc309zcvK3d1NSUb5iqUEOHph2BJKlc5DzBO4QwIYTwN6AZmAacnyTJsz1svgL4JHAhcAGwAHgwhHBKtmM0NDRQU1OzbantrzEaSZKkHOU8DBdCGAaMAXYnJkGfACZlSZi67n8PkCRJ8i56GIbrrmeptrbWYThJkpSv4g3DJUmyCXixtflECKEeuBz4VC9f4jHgomwbVFVVUVVVlWtokiRJ/a4/6iwFIJfM5iji8JwkSVLJy6lnKYTwDeA+YCmwG/B+YDJwZuvzDcB+SZJc0tq+AlgMPAMMI/YoXdi6SJIklbxch+H2Bn5CrJXUCMwFzkyS5IHW50cR5zO1GQZcD+wHbCAmTeckSXJvX4KWJEkqFm+kK0mSKoE30pUkSSoEkyVJkqQsSjpZymQy1NXVUV9fn3YokiSpQjlnSZIkVQLnLEmSJBWCyZJSsWEDLF0KW7emHYkkSdmZLKnoHn8cRo+GMWPghBOgqalwx3r9dZg5ExobC3cMSVJ5M1lS0V11FaxZEx/PnAnTphXmOE8+CYccAscdB+PHw8KFhTmOJKm8mSyp6LZs6dzevLkwx2logLVr4+NXX4WpUwtzHElSeTNZUtFdcw3svHN8PG4cfPKThTnO0KGd28OGFeY4kqTyVtKlAzKZDJlMhpaWFhYuXGjpgDKycmWc4F1X15449bcXX4S3vx1efjkOwz34IOy7b2GOJUkqeXmXDijpZKmNdZaUr82b4bXXYO+9YUiut42WJJWTvJMl/3yorA0dCvvtl3YUkqSBzDlLkiRJWZgsSZIkZWGyJEmSlIXJkiRJUhYmS5IkSVmUdLKUyWSoq6ujvr4+7VAkSVKFss6SJEmqBHnXWSrpniVJkqS0mSxJkiRlYbIkSZKUhcmSJElSFiZLkiRJWZgsSZIkZWGyJEmSlEVJJ0sWpZQkSWmzKKUkSaoEFqXUwNDUBC0txT/u/Plwzjlw2mkwfXrxjy9JGrhMllQUW7bAhRdCTQ2MHAkPPVTcY59xBtx7L/zpT3DuubBiRfGOL0ka2HJKlkIInw4hzA0hNLUuM0IIZ+1gn0khhNkhhI0hhJdCCJf1LWQNRHfcAb/8ZXy8di1cVsSzYM0aWL68vf3GG/DSS8U7viRpYMu1Z2kZ8P+AY1uXPwG/DiEc3t3GIYSxwL3Aw8BRwDeA/wwhXJh3xBqQ3ngje7uQRo6EY49tb48eDRMmFO/4kqSBrc8TvEMIa4AvJknyo26e+ybwriRJDuuwbhrw5iRJTsAJ3hVj3To48UR49lkYNAhuvhk+9rHiHv+GG2DjRvjMZ2DMmOIdW5JUEvKe4D0k7yOGMBh4D7ALMKOHzU4A7u+y7g/Ax0MIQ3tK1Jqbm2lubt7WbmpqyjdMlYjdd4eZM+Oy774wblzxj//Vrxb3mJKk8pDzBO8QwoQQwt+AZmAacH6SJM/2sPk+wMou61YSk7Q9ezpGQ0MDNTU125ba2tpcw1QJ2mUXOPXU4idKkiT1RT5Xwy0AJgLHA98Dbgkh1GXZvmv3Uehh/TZTpkyhsbFx27J06dI8wpQkSeq7nIfhkiTZBLzY2nwihFAPXA58qpvNXyX2LnW0F7AFeL2nY1RVVVFVVZVraJIkSf2uP+osBaCnzGYGcHqXde8AnkiSZHM/HFuSJKmgcq2z9I0QwskhhANa5y5dB0wG/rf1+YYQwq0ddpkG7B9CmBpCOCyE8DHg48D1/RS/JElSQeU6DLc38BNgFNAIzAXOTJLkgdbnRwHbLspOkmRRCOFs4DvAZ4HlwOeSJLmrr4FLkiQVgzfSlSRJlcAb6UqSJBVCSSdLmUyGuro66uvr0w5FkiRVKIfhJElSJXAYTpIkqRBMliRJkrIwWZIkScrCZEmSJCkLkyVJkqQsTJYkSZKyMFmSJEnKoqSTJYtSSpKktFmUUhXnuuvg5z+HAw+EadNg773TjkiSVAR5F6Uc0p9RSKXurrvgS1+Kj+fOhZYW+M1v0o1JklTaSnoYTupvCxdmb0uS1JXJkgrqgQdg9GioroZvfjPtaODss2H48Pb2hRemF4skaWBwzpIKJklgjz1g7dr2dbNnw9FHpxcTwJw5cejtoIPgoovSjUWSVDTOWVLp2bQJ1q3rvG7VqlRC6eSoo+IiSVJvOAyngqmqgksvbW9PnAgnn5xaOJIk5aWkh+EymQyZTIaWlhYWLlzoMNwAdd990NQE73wn7LJL2tFIkipU3sNwJZ0stXHOkiRJ6qO8kyWH4SRJkrIwWZIkScrCZEmSJCkLkyVJkqQsTJYkSZKyMFmSJEnKwmRJkiQpi5JOljKZDHV1ddTX16cdiiRJqlAWpZQkSZXAopSSJEmFYLIkSZKURU7JUghhSghhVghhfQjhtRDC3SGEcTvYZ3IIIelmGd+30CVJkgov156lSUAGOB44HRgC3B9C6M295McBozosL+R4bEmSpKIbksvGSZKc2bEdQvgo8BpwDPDQDnZ/LUmSdTlFJ0mSlLKckqVu1LT+u6YX284JIQwHngWuTZJkek8bNjc309zcvK3d1NTUpyAlSZLylfcE7xBCAKYCjyRJMj/LpiuATwIXAhcAC4AHQwin9LRDQ0MDNTU125ba2tp8w5SyeuUVmDkTNm5MOxJJUqnKu85SCCEDnAOclCTJshz3vQdIkiQ5t7vnu+tZqq2ttc6S+tUvfwkf+ABs2gQTJsBDD8Huu6cdlSSpQIpbZymEcCPwLuDUXBOlVo8Bh/T0ZFVVFdXV1Z0Wqb996UsxUQKYNw9+8pN045EklaZcSweEEMJNxOG0tyVJsijP4x5FHJ6TUjNsWPa2JEmQe89SBrgI+CCwPoSwT+uyU9sGIYSGEMKtHdpXhBDOCyEcEkI4PITQQJy/dFN//ABSvr77XahpvUTh1FPhwx9ONx5JUmnK9Wq4T7f+++cu6z8K/Lj18ShgTIfnhgHXA/sBG4BngHOSJLk3x2NL/WrSJFixAtauhVGjIOQ9mi1JKmfeSFf9btMmuPxyePhhOP54uPFG2GmnHe8nSVIB5f1f4r7WWZK209AA06bFx888AyNHxnWSJA1EJX0j3UwmQ11dHfX19WmHohy88EL2tiRJA4nDcOp3v/oVXHghtJ1a//u/8MEPphuTJKniOQyn0nH++fDAA/CXv8Bxx8EZZ6QdkSRJ+bNnSZIkVYLiVvCWJEmqFCZLkiRJWZgsSZIkZWGyJEmSlIXJkiRJUhYlnSxZlFKSJKXN0gGSJKkSWDpAkiSpEEyWJEmSsjBZkiRJysJkSZIkKQuTJUmSpCxMliRJkrIo6WTJOksqlscegwsvhIsvhsWL045GklRKrLOkirdiBYwbB+vXx/Yhh8CCBRDyrsghSSpB1lmS8vX88+2JEsALL8C6damFI0kqMSZLqnhHHAF77NHePvJIGDEivXgkSaVlSNoBSGkbORIeegi++13YeWeYMiXtiCRJpcQ5S5IkqRI4Z0mSJKkQTJYkSZKyKOlkyTpLA8u8eXDppXD55bByZdrRSJLUP5yzpH7x2mtw2GGwZk1sH3kkPP10ujFJktSBc5aUrmeeaU+UAObOhbVr04tHkqT+YrKkfjF+POy2W3v70ENh991TC0eSpH6TU7IUQpgSQpgVQlgfQngthHB3CGFcL/abFEKYHULYGEJ4KYRwWf4hqxSNGgUPPADvfjdcckl87O1CJEnlIKc5SyGE3wO3A7OIBS2vAyYAdUmSvNHDPmOB+cDNwPeBE4H/Aj6QJMkvenNc5yxJkqQ+yvu/8DlV8E6S5MxORw3ho8BrwDHAQz3sdhmwJEmSK1rbz4UQjgWuyi1USZKk4uvrnKWa1n/XZNnmBOD+Luv+ABy7efPmPh5ekiSpsPK+N1wIIQBTgUeSJJmfZdN9gK5Vd1YCQ1avXs2oUaO226G5uZnm5uZt7aampnzDlCRJ6pO+9CzdBBwJfKAX23adGBUAQg8zgBsaGqipqdm21NbW9iFMSZKk/OVVlDKEcCNwHnBKkiSLdrDtQ8CcJEku77DufODnmzZtGjJ06NDt9umuZ6m2ttYJ3pIkKV/FmeDdOvR2I3A+MHlHiVKrGcC5Xda9A3hi6NChx3e3Q1VVFVVVVbmEJkmSVBC5DsNlgIuADwLrQwj7tC47tW0QQmgIIdzaYZ9pwP4hhKkhhMNCCB8DPg5c39fgJUmSCi3XZOnTxCvg/gys6LC8r8M2o4AxbY3W3qezgcnAU8CXgc8lSXJXnjFLkiQVjTfSlSRJlcAb6UqSJBWCyZLUwfr1cO65MGIEnHUWrFuXdkSSpLSVdLKUyWSoq6ujvr4+7VBUIb7+dfjtb2OS9Pvfw5e/nHZEkqS0lXSy9NnPfpZnn32WWbNmpR2KKsSrr2ZvS5IqT0knS1KxffSjMGxYfDxkCHz84+nGI0lKX973hpPK0amnwuzZMHMmHHMMvPnNaUckSUqbyZLUxRFHxEWSJHAYTpIkKSuTJUmSpCxMliRJkrIwWZIkScqipJMli1IODH/9K/z5z/C3v6UdiSRJ/c8b6apPfvYzuOQS2LIFxo2Dv/wF9tgj7agkSdqON9JVOq65JiZKAAsWwG23pRuPJEn9zWRJfbLzzp3bu+ySThySJBWKyZL6JJNpH3Y766w4JCdJUjmxgrf65IQTYOVKWL8edt897WgkSep/9iypzwYPNlGSJJUvkyVJkqQsSjpZss6SJElKm3WWJElSJbDOkiRJUiGYLEmSJGVh6QCpB1u2xIrk69fD+98PI0emHZEkKQ0mS1IPPvhBuPPO+PiGG+DJJ6GmJtWQJEkpcBhO6kZzc3uiBPDSSzBjRnrxSJLSY7IkdaOqCkaNam8PGgS1tenFI0lKT0knS9ZZUpp+8xs45hg4+GC4+WY4/PC0I5IkpcE6S5IkqRJYZ0mSJKkQTJYkSZKyyDlZCiGcEkK4J4SwPISQhBDO28H2k1u367qMzztqSZKkIsmnztIuwNPA/wB35bDfOKCpQ3tVHseWJEkqqpyTpSRJ7gPuAwghp7lSryVJsi7X40mSJKWpmHOW5oQQVoQQHgwhnFrE40qSJOWtGLc7WQF8EpgNVAEXAw+GECb3VLagubmZ5ubmbe2mpqZut5MkSSq0gidLSZIsABZ0WDUjhFALXNXTPg0NDVxzzTWFDk19sHAh3HsvHHQQnHtu2tFIklQ4fSpKGUJIgPOTJLk7x/3+FbgoSZJur4jrrmeptrbWopQl4rnn4LjjYP362L7uOrj66nRjkiRpBwZcUcqjiMNz3aqqqqK6urrTotJx993tiRLAbbelFookSQWX8zBcCGFX4OAOq8aGECYCa5IkWRJCaAD2S5LkktbtrwAWA88Aw4CLgAtbFyd6D0D779+5PWZMOnFIklQM+cxZOhaY3qE9tfXfW4CPAKOAjn8+hwHXA/sBG4hJ0zlJktybx7FVAj7wAXjqKbjjjjhn6eab045IkqTC8Ua6kiSpEgy4OUuSJEkDgsmSJElSFiWdLGUyGerq6qivr087FFW4P/4RPvIR+PKX4e9/TzsaSVIxOWdJ2oEnn4x1pbZsie33vQ9uvz3dmCRJOXPOklQoM2a0J0oADz2UXiySpOIzWZJ24JhjYFCHT8pb3pJeLJKk4ivGjXSlAe344+EXv4Bbb4XaWvj619OOSJJUTM5ZkiRJlcA5S5IkSYVgsiRJkpRFSSdL1lmSJElpc86SJEmqBM5ZkiRJKgSTJUmSpCxMliRJkrIwWZIkScrCZEmSJCkLkyX12vz5cM45cNppMH162tFIklQc3htOvbJlC5xxBixfHtuPPw4vvACjRqUblyRJhVbSPUsWpSwdr7/enigBvPEGvPRSevFIklQsFqVUryQJvOUt8MQTsT16dByWq6lJNy5Jknop76KUDsOpV0KABx6AG26AjRvhs5+tzERp/Xq46ip4/nk47zz4/OfTjkiSVGj2LEk5+MhH4JZb2ts//zm85z2phSNJ6j1vdyIVw9NPd27Pm5dOHJKk4jFZknJw1lntjwcNgtNPTy8WSVJxOGdJysG118KYMXHO0rveBSefnHZEkqRCc86SJEmqBOU5Z8k6S5IkKW32LEmSpEpQnj1LkiRJaTNZkiRJyiLnZCmEcEoI4Z4QwvIQQhJCOK8X+0wKIcwOIWwMIbwUQrgsr2glSZKKLJ+epV2Ap4F/7M3GIYSxwL3Aw8BRwDeA/wwhXJjHsSVJkooq5zpLSZLcB9wHEEKv5kpdBixJkuSK1vZzIYRjgatyPbYkSVKxFWPO0gnA/V3W/QE4dvPmzUU4vCRJUv6KUcF7H2Bll3UrgSGrV69m1KhR2+3Q3NxMc3PztnZTU1NBA5QkSepJsa6G61pPKUDPw3gNDQ3U1NRsW2prawsdnyRJUreKkSy9Suxd6mgvYMsee+zR7Q5TpkyhsbFx27J06dJCx6gsvv99+MQn4Cc/STsSSZKKrxjJ0gyg673Z3wE8MXTo0G53qKqqorq6utOidHznO3DZZfCjH8Ell8Ctt6YdUen44Q/hTW+CvfaCO+9MOxpJUqHkU2dp1xDCxBDCxNZVY1vbY1qfbwghdPyTOg3YP4QwNYRwWAjhY8DHgev7GrwK709/6tyePj2dOErNkiUxiVy7FlatgosvhsbGtKOSJBVCPj1LxwJzWheAqa2Pv9baHgWMads4SZJFwNnAZOAp4MvA55IkuSuviFVUxx6bvV2p1qyBlpb2dnMzeB2CJJUnb6SrrLZsgWuvhVmzYNIk+OIXoXfltcrbli1w+unw5z/H9nnnwS9/6XsjSSUs729okyUpT83NcM89MGQInHsuDB6cdkSSpCxMliRJkrLIO1kqVp2lvGQyGerq6qivr087FEmSVKHsWZIkSZWgPHuWJEmS0mayJEmSlIXJkiRJUhYmS5IkSVmYLEmSJGVhsiRJkpRFSSdL1lnSQLFyJdx3HyxalHYkkqT+Zp0lqY+eew5OPhlefx2GD4+3QHn729OOSpLUhXWW1L+SBJ59FpYtSzuS0jdtWkyUADZuhG9/O914JEn9y2RJ22lpgfPOg8MPh/33h//6r7QjKm1dOzvt/JSk8uIwnLbzhz/AmWe2t4cPhzfegEGm1t1qbIR3vhMeeQTGjYtzl8aOTTsqSVIXeQ/DDenPKFQeuiZFIcRF3aupgYcfhg0bYKed0o5GktTf7CvQdk47Dd773vh4yBC46SaTpd4wUZKk8uQwnHq0eDHsuivsuWfakUiS1GcOw6n/HXBA2hFIkpS+kh6GsyilJElKm8NwkiSpEliUUpIkqRBMliRJkrIwWZL60erV8OlPw7vfDQ8+mHY0kqT+4NVwUj+64IJYoBLgt7+Fp56C8eNTDUmS1Ef2LEn96LHH2h83N8dkSZI0sJksSf3olFPaH++0Exx7bHqxSJL6R0kPw2UyGTKZDC0tLWmHUhFaWuC//xuWL4f3vc/ho3z84hdw7bVx7tKll8LBB6cdkSSpr6yzpG0+9Sn4wQ/i4+pqmDMHDjww3ZgkSeon1llS3/3yl+2Pm5q8mkuSJDBZUgeHHtq5fcgh6cQhSVIpyStZCiF8JoSwKISwMYQwO4RwcpZtJ4cQkm4WZ8SUmJ/9DM4+GyZOhJtugsmT045IkqT05TzBO4TwPuAG4DPAX4BPAfeFEOqSJFmSZddxQFOH9qpcj63CGjMGfve7tKOQJKm05NOzdCXwoyRJfpgkyXNJklwBLAU+vYP9XkuS5NUOi5e4qext2gTpXkMhSeqrnJKlEMIw4Bjg/i5P3Q+8dQe7zwkhrAghPBhCODXbhs3NzTQ1NXVapIEkSeLVhcOHw557OllekgayXHuW9gQGAyu7rF8J7NPDPiuATwIXAhcAC4AHQwin9LA9DQ0N1NTUbFtqa2tzDFNK1733xjIMSQJr1sBHP5p2RJKkfOVblLLrwELoZl3cMEkWEBOkNjNCCLXAVT29+JQpU7jyyiu3tZuamkyYNKCsX9+5beeoJA1cufYsrQZa2L4XaS+2723K5jGgxwvTq6qqqK6u7rRIA8k558CRR7a3r746vVgkSX2TU89SkiSbQgizgdOBX3V46nTg1zm81FHE4TnLB6gs7bYbPPooPPIIjBwJRx+ddkSSpHzlMww3FfhJCOEJYAZxPtIYYBpACKEB2C9Jkkta21cAi4FngGHARcT5SxcCWSd6qziefBLuvx/q6uBd70o7mvKxyy5wxhlpRyFJ6quck6UkSe4IIewB/BswCpgPnJ0kycutm4wiJk9thgHXA/sBG4hJ0zlJktzbl8DVPx5/HE45JV7iDjB1Knz+8+nGJElSKfFGuhVuyhT4939vbx93HDz2WHrxSJJUIN5IV/k58MDsbfWPdetg3jzYsCHtSCRJucq3dIDKxMc/DgsWwD33xDlLN96YdkTlZ+ZMOPNMWLsWDj4YHnoIRo1KOypJUm+V9DBcJpMhk8nQ0tLCwoULHYbTgHTWWfD737e3/+VfOg99SpKKojyH4T772c/y7LPPMmvWrLRDkfIWQva2JKm0lXSyJJWDa6+N94cDGD8errgi1XAkSTlyzpJUYEcfDYsXw/LlsP/+MGxY2hFJknJhsiQVwS67wCE93uBHklTKHIarUA88EO9ddsQR8Nvfph2NJEmly2SpAjU2wvnnx7o/zzwD730vrFqVdlSV4e674YADYMwY+PnP045GktQbJksVaNUqeOON9vaGDbByZXrxVIq1a+EDH4CXX4alS+Hii+HVV9OOSpK0IyWdLGUyGerq6qivr087lLIydiy89a3t7WOOgXHj0ounUqxdCxs3trc3bYLVq9OLR5LUOyVdlLKN94brf3//O9xyC2zdCpdcArvtlnZE5W/rVnjHO+DBB2P7pJNg+nQY4mUWklQMeVe5M1mSiqi5Ge68MyZO730vDB+edkSSVDFMlqSBaN06qK6GQSU9IC5JZaE8b3ei/vfoo/D978Pzz6cdSWX7299g0iQYMSIWqpw/P+2IJEk9MVmqILfdFufJXHZZrCr9xBNpR1S5Mhl46KH4eNkyuPLKdOORJPXMZKmC/OhH0DbqumED/PSn6cZTyTqWbuiuLUkqHSZLFWTffbO3VTyXXgqjR8fHVVVw9dXpxiNJ6pkTvCvIypXwwQ/Gyt1nnBF7mrypa3rWrYM5c+DAA+OcpeXL4ZxzTGIlqUDK82q4TCZDJpOhpaWFhQsXmiypLF19NTQ0xMf77AOzZ5swSVIBlGey1Maepb654w743vdgr73gO9+B/fZLOyJ1tPfe8Npr7e2bb4ZPfCK9eCSpTOWdLFk7uMzNnh2H3rZuje1ly2L5AJWOffftnCzZqyRJpcUJ3mXumWfaEyWI85VUWm67DY46KvYwve99MHUqvPvdsHhx2pFJksCepbJ30knxvm/r18f2mWemG4+2d/jh8OSTsHAhTJgQb7ALsXCoxSolKX0mS2Vq/Xr493+Pd7WfNi0Ox+29N1x+edqRqSfPPdeeKEHsFdyyxRvtSlLanOBdps44A+6/Pz7eZZc4/DZ2bLoxKbsVK+CII2DNmtjeY4/47+TJcOutsPPOqYUmSeXAe8Ops+nT2x+/8QbMnJleLOqdUaPgkUfg85+H+np4/fW43HVXe2kBSVLxlXSylMlkqKuro76+Pu1QBoxbboGLLup8RdXQoXEujErfYYfFCd5dr4i75RZ4+9s7J8GSpOJwGK6M/OxnsUxAm4kT4y01PvMZOOus1MJSHn79a7jwQmhp6bx+p53gT3+Kk8J32y2d2CRpgLLOUiVbsgT+8z/hwQc7r99zT7jnnnRiUt/8wz/AY4/Bn/8MX/xi+/oNG+CEE6C6OiZUkyenFaEkVQ6TpQFs9ep4ufn73w9Ll27//HHHFT8m9Z9jj4VjjokV2J94ovNzTU3xPnI77xx7E7/1rXjV3KCSHliXpIHJYbgBZv16mDEjJkeXXx4nb3d1xhlw4okwZYqXnZeDdevibWoefrjnOUuDBsHw4XDNNTBiREyyJk4sZpSSVPK8N1w52Lw5TsZesiT+4Vu3LtZHGjUK/ud/Yk/S/Pnx+Z6MGgWLFkFVVdHCVpGsWweTJsHcuRACZPvoDh4cb8rb0hKrgm/aFOevnXxyrAw+aRI0NsZketw4WLkyDtsOGRJrOw0dWqyfSpKKpnyTpe99D775zSZefrmGU05pZNGiasaPh+XLYePGePXQ3Lnx5rBJEte/+c3w7LMx4dh331gJ+dBDYdUq+Nvf4uTYuXNjkcahQ2MvzRFHwAsvxD8W++8fCwIeeGAc7li3Ll5NNndu/F97dTW89BLU1cU/PFu3wsEHx0Rm//3jvJJVq2Icc+fGibgjRsCLL8Y/TEuWxNgPOgjmzIkxNjXFP1gHHBCTnRBiLJs37/gPI8TXHTcOrrsu/iwqT83N8Xy+6y74+tfzf51hw9oLYI4cGc/XPfeM59v69fHcff55eNOb4vpFi+JnYMmSeL6PH99+vm/eHO9t9+Y3x8/NrrvG13zxRbr9rO67bzynly+Pr/n88zG533dfWLAADjkk/sdg/fp4Ls+dG28CXVUVj3/EEfG1Bw2Kn5dnnomfpaamWKPqyCO3/6wefnj8GbZuja8/fz6MGRPfz7bY582Ln9U994zfBYcdBq+8ErcZP37775kjj4yFRHv6njniCHj66R3H3vY9s3ZtjOPpp/see9v7nkvsbe97T7FPmBBfu6fYe3rfFy+OSXtb7LW18dzrGvvIkXFaQcfYDzssxtF2dWjH7/eqqvgz9Sb2Hb3vHWOvqYG//rVz7IceGuPsGvv8+fF873rOdD3f22Jve9+7xv766+3n+9NPx/di+PC+nTM7et+7xt7TZ7W72EePjo/7Gnu+58z8+bF+4MiRMfZTTok3IN9ppx1+9Q28ZCmEEBobG7d291xzczPNzc3Mmwdnnw2wHqgDlgLl27OUq0GD4hfoxInwm9/ED7kqx/Tp8cvl9tutoyWpsv3zP8O//mv2bWpqamqA9UkeiU+ayVI10JjKwSVJUiWqSZKkKdedSrpnad26eOPXV15ZAbwFeBbYr8vrbD9EtaNhq+zP1wOzctwnnzjqGTx4Vrd1dDZsiI8PPDAOER51VOx63LwZLr00DjkMHx5ff0fq6+uZNWv7n6e/FeM4TU1N1NbWsnTp0qLMXRvI793WrbFsxPr1sTv717+GmTPrGTNmFosWxW7tBQvitm03Wh46NJ5j0HmYrk3H8zn7+d77z1D+n6smoJbuepsL8VkNYVY/xp7bcfr2mt09X7zjdDwP+vf333FN53OhN9MWenOc7p9v/5n6//e/489Q32Lv/XEG2t/VQYPglluauPji7H8f+tKzlNq1Ur0JdvToOMH5Bz+Af/s3uPPO3Vi2rJrjjovjpc3N8Ja3xFtEHHRQ/KOweHGcxPrYY/Gy6kMOgVmz4Oij4eWX41jp8cfDX/4SX3+nneI4+Yknxsuzv/a1wdx1VzWPPx7HU1eujGOyJ5wAjz4ax8H32COOvR5/fBxv3bo1vv6jj8Z5TE1NcYz3pJPilWsjRsQ5BnPmxMv5X3wRvvCFwUyfXs3MmTEp2m23OMZ74onw6qsx9v337/v7PHjw4KIkFsU6DkB1dXVZ/UyFOs7FF7c//qd/grq6wcyf336c116LE8GHDIlzQA46KH5GXn89nvuPPBLn9e28c5wj8Ja3wFNPxc/ZxInxMzR+fJwovnw5vPWt8Xy/7rrB3H57NU8+GW/bsmBB/A/AccfFK/oOOih+2S1aFM/3xx+Pxzj00Dic2PZZbWyMn7u2z+rOO8fXeutbY2xXXgl//GM18+ZVc+SR8edZvbrzZ3XPPeNckxNO2P6zOn58nOPyyisxjhkzYPfdYyI5Z078eRcsgKuvHsxvf1u97Xtm69b4PZNv7LNnx/d8woR4zAkT4jynK64YzO9+V82jj8a5GCNH9j32tve9Y+yf/GT8eWbO3HHs++0X54b0FHvb+75qVXy+Y+wf+tBg7r67mnnz+h77gQfGc6bj+77TTnGe5kMPxSt/v//9ajZtqu6X2OfNi9/vXWOfOnUwU6dW88or8fv90Ufj9Ie27/eusf/lL+335Ox4vrfFPnNm/M/wkiWd3/dvfWswP/xhNc8/v+PYTzwx7jNyZDzn586NrzNvXrxY45hjen7fGxoG89OfVjNnTvysLlyYf+xtn8v99otzodpif/JJuOaawfzqV9X9GntNTfwb+eST8Tjjx8dYs/19yKdHqePOaS69snTp0gRIli5d2ttd8nbTTTcV/BgeJz+NjY0JkDQ2Nhb8WElSXu9duR2nmOdCOb1v5XgczwWPkyS9Pg/yzldK/mo4gGXLlm0bfhk9enShY1KJqpQSEtoxzwW18VwQ9Po8yPtquAFR77eqtWhQlcWDKlpVVRVf+cpXPA/kuaBtPBcEhT8PBkTPkv9zkCRJfTTw6izlokOZgZqkLxO0JEmScjRQkqUA7Eael/xJkiTla0AkS5IkSWkZEBO8JUmS0mKypJITQvhqCCHpsrza4fnQus3yEMKGEMKfQwiHpxmz+i6EcEoI4Z7W32sSQjivy/M7/L2HEKpCCDeGEFaHEN4IIfwmhGC9kQGmF+fCj7v5jnisyzaeCwNcCGFKCGFWCGF9COG1EMLdIYRxXbYpyveCyZJK1TPAqA7LhA7P/TNwJfCPxDr6rwIPhBB2K3aQ6le7AE8Tf6/d6c3v/QbgfOD9wEnArsBvQwiDCxSzCmNH5wLA7+n8HXF2l+dvwHNhoJsEZIDjgdOJdx25P4SwS4dtivO90JeKli4uhViArwJP9fBcAFYA/9JhXRWwDvhU2rG79Ns5kADn5fJ7B2qATcD7OmyzL9ACnJH2z+TSP+dC67ofA3dn2cdzoQwXYGTr+XBKa7to3wv2LKlUHdLarboohHB7COHA1vVjgX2A+9s2TJKkGfg/4K0pxKni6M3v/RhgaJdtlgPz8dwoR5Nbh2YWhhBuDiHs1eE5z4XyVNP675rWf4v2vWCypFL0OHAJcAZwKfHD8GgIYY/WxwAru+yzssNzKj+9+b3vA2xKkmRtlm1UHu4DPgS8DfgCcfjlTyGEtvLNngtlprWE0FTgkSRJ5reuLtr3wpDcwpUKL0mS+zo054UQZgB/BT4MtE3i7FrzInSzTuUnn9+750aZSZLkjg7N+SGEJ4CXgXOAX2bZ1XNh4LoJOJI456irgn8v2LOkkpckyRvAPOAQ4uQ92P5/BHux/f8uVD5683t/FRgWQhiRZRuVoSRJVhCTpUNaV3kulJEQwo3Au4BTkyRZ1uGpon0vmCyp5LV2rR9GnMi3iHjyn97h+WHEqyYeTSVAFUNvfu+zgc1dthkFHIHnRllrHaKvJX5HgOdCWWgtC3ATcAHwtiRJFnXZpGjfCw7DqeSEEK4H7gGWELP/LwHVwC1JkiQhhBuAq0MILwAvAFcDfwd+mk7E6g8hhF2BgzusGhtCmAisSZJkyY5+70mSNIYQfgR8O4TwOnES6PXEXsk/Fu0HUZ9lOxdal68CdxGTowOAbwCrgV+B50IZyQAfBP4BWB9CaOtBakySZENv/h7027mQ9qWALi5dF+B2YDnxcs9XiF+KdR2eD8QvyxXARuKVD0ekHbdLn3/vk4lzCLouP+7t7x0YDtwIvN76hXkPUJv2z+bSf+cCsBPwB+C11u+Il1vX13Z5Dc+FAb70cA4kwEc6bFOU7wXvDSdJkpSFc5YkSZKyMFmSJEnKwmRJkiQpC5MlSZKkLEyWJEmSsjBZkiRJysJkSZIkKQuTJUmSpCxMliRJkrIwWZIkScrCZEmSJCkLkyVJkqQs/j8CDq/Lp2iscgAAAABJRU5ErkJggg==", "text/plain": [ "Graphics object consisting of 200 graphics primitives" ] }, "execution_count": 9, "metadata": { }, "output_type": "execute_result" } ], "source": [ "# Plot series terms versus computed density\n", "K. = QQ['U']\n", "P. = K[['X']]\n", "\n", "# Set the value of n to test\n", "N = 200\n", "mser = (1 - X)/(1 - 2*X - (U-1)*X*(1-X) + O(X^(N+1)))\n", "uvals = mser[N]\n", "\n", "plt = point([])\n", "for k in range(N):\n", " plt += point([k,uvals[k]])\n", "\n", "print(\"The following plot shows the coefficients of [x^({})]C(u,x)\".format(N))\n", "plt" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ ] } ], "metadata": { "kernelspec": { "display_name": "SageMath 9.5", "language": "sagemath", "metadata": { "cocalc": { "description": "Open-source mathematical software system", "priority": 10, "url": "https://www.sagemath.org/" } }, "name": "sage-9.5", "resource_dir": "/ext/jupyter/kernels/sage-9.5" }, "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.5" } }, "nbformat": 4, "nbformat_minor": 4 }