{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "### Example 2.29 (Simple Lattice Walks in $\\mathbb{N}^2$)\n", "Determining connection coefficients for lattice walks on the steps $\\{(\\pm1,0),(0,\\pm1)\\}$ restricted to $\\mathbb{N}^2$. \n", "*Requirements: [ore_algebra package](https://github.com/mkauers/ore_algebra/), internet access (optional)*" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# Define the rings of differential and shift operators\n", "from ore_algebra import *\n", "Ind. = PolynomialRing(QQ); Shift. = OreAlgebra(Ind)\n", "Pols. = PolynomialRing(QQ); Diff. = OreAlgebra(Pols)" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(16*z^4 - z^2)*Dz^3 + (128*z^3 + 8*z^2 - 6*z)*Dz^2 + (224*z^2 + 28*z - 6)*Dz + 64*z + 12" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# The generating function of these walks satisfies the differential equation\n", "diff = z^2*(4*z-1)*(4*z+1)*Dz^3+2*z*(4*z+1)*(16*z-3)*Dz^2+(2*(112*z^2+14*z-3))*Dz+64*z+12\n", "diff" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(-n^3 - 9*n^2 - 26*n - 24)*Sn^2 + (8*n^2 + 36*n + 40)*Sn + 16*n^3 + 80*n^2 + 128*n + 64" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# The counting sequence satisfies the recurrence defined by the following operator\n", "rec = diff.to_S(Shift)\n", "rec" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[1, 2, 6, 18, 60, 200, 700, 2450, 8820, 31752]" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# We verify the first terms match what we expect and check with OEIS\n", "LST = rec.to_list([1,2],10)\n", "LST" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0: A005566: Number of walks of length n on square lattice, starting at origin, staying in first quadrant." ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Search the inital sequence terms in the OEIS *this requires an internet connection -- it is not required in the following*\n", "oeis(LST)" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ], "text/latex": [ "\\begin{math}\n", "\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left[(4)^{n}n^{-1}\\Bigl(1 - \\frac{3}{2}n^{-1} + O(n^{-2})\\Bigr), (-4)^{n}n^{-3}\\Bigl(1 - \\frac{9}{2}n^{-1} + O(n^{-2})\\Bigr)\\right]\n", "\\end{math}" ], "text/plain": [ "[4^n*n^(-1)*(1 - 3/2*n^(-1) + O(n^(-2))),\n", " (-4)^n*n^(-3)*(1 - 9/2*n^(-1) + O(n^(-2)))]" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Asymptotically, the counting sequence is a C-linear combination of the following\n", "show(rec.generalized_series_solutions(n=2))" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[z^(-2) - 4*z^(-1)*log(z) - 8*log(z) - 16*z*log(z) - 8*z, z^(-1), 1 + 2*z]" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# There is a basis of solutions for the C-vector space defined by the differential operator, with the following expansions at the origin\n", "diff.local_basis_expansions(0, order = 4)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "# Note the final basis element is the generating function, so in this basis the generating function has coordinates\n", "ini = [0,0,1]" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "[log(z - 1/4) - 6*(z - 1/4)*log(z - 1/4) + 31*(z - 1/4)^2*log(z - 1/4) - 150*(z - 1/4)^3*log(z - 1/4) - 2/3*(z - 1/4)^3,\n", " 1 - 14*(z - 1/4)^2 + 108*(z - 1/4)^3,\n", " (z - 1/4) - 15/2*(z - 1/4)^2 + 43*(z - 1/4)^3]" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# There is a basis of solutions for the C-vector space defined by the differential operator, with the following expansions at z=1/4\n", "diff.local_basis_expansions(1/4, order = 4)" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ], "text/latex": [ "\\begin{math}\n", "\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left(\\begin{array}{rrr}\n", "\\verb|[-5.1|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|0.0138]|\\phantom{\\verb!x!}\\verb|+|\\phantom{\\verb!x!}\\verb|[+/-|\\phantom{\\verb!x!}\\verb|5.45e-3]*I| & \\verb|[+/-|\\phantom{\\verb!x!}\\verb|2.25e-4]|\\phantom{\\verb!x!}\\verb|+|\\phantom{\\verb!x!}\\verb|[+/-|\\phantom{\\verb!x!}\\verb|1.98e-4]*I| & \\verb|[-1.3|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|0.0268]|\\phantom{\\verb!x!}\\verb|+|\\phantom{\\verb!x!}\\verb|[+/-|\\phantom{\\verb!x!}\\verb|2.81e-5]*I| \\\\\n", "\\verb|[4.1e+1|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|0.306]|\\phantom{\\verb!x!}\\verb|+|\\phantom{\\verb!x!}\\verb|[1.6e+1|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|0.0109]*I| & \\verb|[4.0|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|2.68e-4]|\\phantom{\\verb!x!}\\verb|+|\\phantom{\\verb!x!}\\verb|[+/-|\\phantom{\\verb!x!}\\verb|3.34e-4]*I| & \\verb|[-2.4|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|9.34e-3]|\\phantom{\\verb!x!}\\verb|+|\\phantom{\\verb!x!}\\verb|[4.0|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|3.9e-5]*I| \\\\\n", "\\verb|[-2.5e+2|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|0.0646]|\\phantom{\\verb!x!}\\verb|+|\\phantom{\\verb!x!}\\verb|[-9.6e+1|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|0.0489]*I| & \\verb|[-1.6e+1|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|1.12e-3]|\\phantom{\\verb!x!}\\verb|+|\\phantom{\\verb!x!}\\verb|[+/-|\\phantom{\\verb!x!}\\verb|1.41e-3]*I| & \\verb|[1.3e+1|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|0.110]|\\phantom{\\verb!x!}\\verb|+|\\phantom{\\verb!x!}\\verb|[-2.4e+1|\\phantom{\\verb!x!}\\verb|+/-|\\phantom{\\verb!x!}\\verb|1.6e-4]*I|\n", "\\end{array}\\right)\n", "\\end{math}" ], "text/plain": [ "[ [-5.1 +/- 0.0138] + [+/- 5.45e-3]*I [+/- 2.25e-4] + [+/- 1.98e-4]*I [-1.3 +/- 0.0268] + [+/- 2.81e-5]*I]\n", "[ [4.1e+1 +/- 0.306] + [1.6e+1 +/- 0.0109]*I [4.0 +/- 2.68e-4] + [+/- 3.34e-4]*I [-2.4 +/- 9.34e-3] + [4.0 +/- 3.9e-5]*I]\n", "[[-2.5e+2 +/- 0.0646] + [-9.6e+1 +/- 0.0489]*I [-1.6e+1 +/- 1.12e-3] + [+/- 1.41e-3]*I [1.3e+1 +/- 0.110] + [-2.4e+1 +/- 1.6e-4]*I]" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# We can compute the change of basis matrix between the basis with series at the origin and the basis with series at 1/4\n", "# We take a low numeric accuracy here to fit on screen\n", "show(diff.numerical_transition_matrix([0,1/4],1e-1))" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "([-1.2732395447351627+/-2.54e-17]+[+/-2.65e-22]*I)*log(z - 1/4) + ([7.6394372684109761+/-5.21e-17]+[+/-1.59e-21]*I)*(z - 1/4)*log(z - 1/4) + ([-2.3906971441245563+/-1.99e-17]+[4.0000000000000000+/-2.78e-17]*I)*1 + ([12.890661954217663+/-2.55e-16]+[-24.000000000000000+/-1.12e-16]*I)*(z - 1/4)" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Thus the generating function has a singular expansion at z=1/4 which begins\n", "bas = diff.numerical_transition_matrix([0,1/4]) * vector(ini)\n", "loc = diff.local_basis_expansions(1/4,2)\n", "add([a*b for [a,b] in zip(bas,loc)])" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "([1.2732395447351627 +/- 2.54e-17])*4^n/n" ] }, "execution_count": 12, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# The dominant singular behaviour is given by the term with log(z-1/4)\n", "# Using the coefficient extraction [z^n]log(1-4*z) = -4^n/n dominant asymptotic behaviour is\n", "var('n')\n", "-bas[0].real() * 4^n/n" ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1.27323954473516" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# We will determine the leading constant to be 4/pi\n", "(4/pi).n()" ] }, { "cell_type": "code", "execution_count": 0, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "SageMath 9.2", "language": "sage", "name": "sagemath" }, "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.8.5" } }, "nbformat": 4, "nbformat_minor": 4 }