AsymptoticRing SageMath Demonstration¶

Power Series Expansion¶

$T(z,t)$: Bivariate Generating function for combinatorial class of rooted plane trees ($z$ counts inner nodes, $t$ counts leaves)

In [3]:

def T(z, t):
    return 1/2 * (1 - (z-t) - sqrt(1 - 2*(z+t) + (z-t)^2))

In [4]:

T(z,t).taylor((z,0), (t,0), 4)

Out[4]:

$\newcommand{\Bold}[1]{\mathbf{#1}}t^{3} z + 3 \, t^{2} z^{2} + t z^{3} + t^{2} z + t z^{2} + t z + t$

Expansion Operator $\Phi$¶

In [5]:

def phi(f, z, t):
    return lambda z, t: (1-t) * f(z/(1-t)^2, z*t/(1-t)^2)

In [6]:

phi(T, z, t)(z,t).simplify_rational()

Out[6]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} \, {\left(t - 1\right)} \sqrt{\frac{t^{2} - 2 \, {\left(t + 1\right)} z + z^{2} - 2 \, t + 1}{t^{2} - 2 \, t + 1}} - \frac{1}{2} \, t - \frac{1}{2} \, z + \frac{1}{2}$

In [7]:

phi(T, z, t)(z,t).taylor((z,0), (t,0), 4)

Out[7]:

$\newcommand{\Bold}[1]{\mathbf{#1}}t^{3} z + 3 \, t^{2} z^{2} + t z^{3} + t^{2} z + t z^{2} + t z$

In [8]:

phi(phi(T, z, t), z, t)(z,t).taylor((z,0), (t,0), 4)

Out[8]:

$\newcommand{\Bold}[1]{\mathbf{#1}}3 \, t^{2} z^{2} + t z^{3} + t z^{2}$

Motivation: Be Wary of Symbolics!¶

In [9]:

expr = 1/(1 + 1/x)
expr

Out[9]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{\frac{1}{x} + 1}$

In [10]:

expr.series(x==0, 10)

Out[10]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{2} x + \mathcal{O}\left(x^{10}\right)$

In [11]:

expr.simplify_rational()

Out[11]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{x}{x + 1}$

In [12]:

expr.simplify_rational().series(x==0, 6)

Out[12]:

$\newcommand{\Bold}[1]{\mathbf{#1}}1 x + {(-1)} x^{2} + 1 x^{3} + {(-1)} x^{4} + 1 x^{5} + \mathcal{O}\left(x^{6}\right)$

factor() is currently broken (and will be fixed soon, hopefully).

In [13]:

expr = (2*exp(x) + exp(-x))
expr

Out[13]:

$\newcommand{\Bold}[1]{\mathbf{#1}}e^{\left(-x\right)} + 2 \, e^{x}$

In [14]:

expr.factor()

Out[14]:

$\newcommand{\Bold}[1]{\mathbf{#1}}3 \, e^{x}$

In [15]:

expr = (x + sqrt(x))/(x - sqrt(x))
expr

Out[15]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{x + \sqrt{x}}{x - \sqrt{x}}$

In [16]:

expr.factor()

Out[16]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{x + 1}{x - 1}$

In [17]:

(x + sqrt(x)).factor()

Out[17]:

$\newcommand{\Bold}[1]{\mathbf{#1}}x + 1$

AsymptoticRing and Friends¶

In [18]:

SCR = SR.subring(no_variables=True) # complex numbers + symbolic constants
A = AsymptoticRing(growth_group='QQ^n * n^QQ * log(n)^QQ', 
                   coefficient_ring=SCR, 
                   default_prec=5)
n = A.gen()

In [20]:

Out[20]:

Asymptotic Ring <QQ^n * n^QQ * log(n)^QQ> over Symbolic Constants Subring

In [22]:

A.an_element()

Out[22]:

$\newcommand{\Bold}[1]{\mathbf{#1}}O\!\left(\left(\frac{1}{2}\right)^{n} n^{\frac{1}{2}} \log\left(n\right)^{\frac{1}{2}}\right)$

In [23]:

42 * (1/3)^n * n^(5/2) * log(n)^2  +  O((1/3)^n * n * log(n)^(1/2))

Out[23]:

$\newcommand{\Bold}[1]{\mathbf{#1}}42 \left(\frac{1}{3}\right)^{n} n^{\frac{5}{2}} \log\left(n\right)^{2} + O\!\left(\left(\frac{1}{3}\right)^{n} n \log\left(n\right)^{\frac{1}{2}}\right)$

Straightforward Computations¶

In [24]:

log(n+1)

Out[24]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\log\left(n\right) + n^{-1} - \frac{1}{2} n^{-2} + \frac{1}{3} n^{-3} - \frac{1}{4} n^{-4} + \frac{1}{5} n^{-5} + O\!\left(n^{-6}\right)$

In [25]:

1/(1 + n)

Out[25]:

$\newcommand{\Bold}[1]{\mathbf{#1}}n^{-1} - n^{-2} + n^{-3} - n^{-4} + n^{-5} + O\!\left(n^{-6}\right)$

In [26]:

(1 + 1/n)^n

Out[26]:

$\newcommand{\Bold}[1]{\mathbf{#1}}e - \frac{1}{2} \, e n^{-1} + \frac{11}{24} \, e n^{-2} - \frac{7}{16} \, e n^{-3} + \frac{2447}{5760} \, e n^{-4} + O\!\left(n^{-5}\right)$

Example: Catalan Numbers¶

In [27]:

try:
    factorial(n)
except TypeError as e:
    print(e)

no common canonical parent for objects with parents: 'Asymptotic Ring <QQ^n * n^QQ * log(n)^QQ> over Symbolic Constants Subring' and 'Asymptotic Ring <(e^(n*log(n)))^QQ * (e^n)^QQ * n^QQ * log(n)^QQ> over Symbolic Constants Subring'

In [28]:

binomial_2n_n = asymptotic_expansions.Binomial_kn_over_n(n, k=2, precision=5)
binomial_2n_n

Out[28]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{\sqrt{\pi}} 4^{n} n^{-\frac{1}{2}} - \frac{1}{8 \, \sqrt{\pi}} 4^{n} n^{-\frac{3}{2}} + \frac{1}{128 \, \sqrt{\pi}} 4^{n} n^{-\frac{5}{2}} + \frac{5}{1024 \, \sqrt{\pi}} 4^{n} n^{-\frac{7}{2}} - \frac{21}{32768 \, \sqrt{\pi}} 4^{n} n^{-\frac{9}{2}} - \frac{399}{262144 \, \sqrt{\pi}} 4^{n} n^{-\frac{11}{2}} + \frac{869}{4194304 \, \sqrt{\pi}} 4^{n} n^{-\frac{13}{2}} + O\!\left(4^{n} n^{-\frac{15}{2}}\right)$

In [29]:

catalan_asy = binomial_2n_n / (n+1)
catalan_asy

Out[29]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{\sqrt{\pi}} 4^{n} n^{-\frac{3}{2}} - \frac{9}{8 \, \sqrt{\pi}} 4^{n} n^{-\frac{5}{2}} + \frac{145}{128 \, \sqrt{\pi}} 4^{n} n^{-\frac{7}{2}} - \frac{1155}{1024 \, \sqrt{\pi}} 4^{n} n^{-\frac{9}{2}} + \frac{36939}{32768 \, \sqrt{\pi}} 4^{n} n^{-\frac{11}{2}} + O\!\left(4^{n} n^{-\frac{13}{2}}\right)$

Extract coefficients from GF¶

In [30]:

def catalan(z):
    return (1 - sqrt(1 - 4*z)) / (2*z)

In [31]:

catalan_asy = A.coefficients_of_generating_function(catalan, 
                                                    singularities=(1/4,),
                                                    precision=5)
catalan_asy

Out[31]:

Preparation for later: There are $C_{n-1}$ plane trees with $n$ nodes.

In [32]:

shifted_catalan_asy = catalan_asy.subs(n=n-1)
shifted_catalan_asy

Out[32]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{1}{4 \, \sqrt{\pi}} 4^{n} n^{-\frac{3}{2}} + \frac{3}{32 \, \sqrt{\pi}} 4^{n} n^{-\frac{5}{2}} + \frac{25}{512 \, \sqrt{\pi}} 4^{n} n^{-\frac{7}{2}} + \frac{105}{4096 \, \sqrt{\pi}} 4^{n} n^{-\frac{9}{2}} + \frac{1659}{131072 \, \sqrt{\pi}} 4^{n} n^{-\frac{11}{2}} + O\!\left(4^{n} n^{-6}\right)$

Leaf-Cutting: Expected Tree Size after $r$ Rounds¶

In [33]:

var('u v')
var('r', domain='integer')
assume(r > 0)

In [34]:

A.<n> = AsymptoticRing('QQ^n * n^QQ', SR)
A_func.<Z> = AsymptoticRing('Z^QQ * log(Z)^QQ', SR)

Some Preparations¶

In [35]:

# solve z = u/(1+u)^2, get u as a function of z
def upsilon(z):
    return (1 - sqrt(1 - 4*z))/(2*z) - 1

In [36]:

# results from Proposition
def G_r(u, v):
    return (1-u^(r+2)) / ((1-u^(r+1)) * (1+u)) * \
            T(u * (1-u^(r+1))^2 / (1-u^(r+2))^2 * v, 
              u^(r+1) * (1-u)^2 / (1-u^(r+2))^2 * v)

Sum of remaining tree sizes: $[z^n] \partial_v G_r(u,v)|_{v=1}$.

Function Preparation¶

Expand $\partial_v G_r(u,v)|_{v=1}$ at $z = 1/4$ in terms of $Z$ where $Z = (1 - 4z)^{-1}$.

In [37]:

z_asy = 1/4 * (1-Z^(-1))
u_asy = upsilon(z_asy) + O(Z^(-3))

In [38]:

# Expansion of u(z) in terms of Z
u_asy

Out[38]:

$\newcommand{\Bold}[1]{\mathbf{#1}}1 - 2 Z^{-\frac{1}{2}} + 2 Z^{-1} - 2 Z^{-\frac{3}{2}} + 2 Z^{-2} - 2 Z^{-\frac{5}{2}} + O\!\left(Z^{-3}\right)$

In [39]:

G_r(u, v).diff(v, 1)(v=1).simplify_rational()

Out[39]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\frac{{\left(u^{3} - u^{2}\right)} u^{2 \, r} + {\left(u^{3} - 2 \, u^{2} + u\right)} u^{r} - u^{2} + {\left({\left(u^{3} + u\right)} u^{r} - u - u^{2 \, r + 3}\right)} \sqrt{\frac{{\left(u^{4} - 2 \, u^{3} + u^{2}\right)} u^{2 \, r} - 2 \, {\left(u^{3} - 2 \, u^{2} + u\right)} u^{r} + u^{2} - 2 \, u + 1}{u^{2 \, r + 4} - 2 \, u^{r + 2} + 1}} + u}{2 \, {\left({\left(u^{4} + u^{3}\right)} u^{2 \, r} - {\left(u^{3} + 2 \, u^{2} + u\right)} u^{r} + u + 1\right)} \sqrt{\frac{{\left(u^{4} - 2 \, u^{3} + u^{2}\right)} u^{2 \, r} - 2 \, {\left(u^{3} - 2 \, u^{2} + u\right)} u^{r} + u^{2} - 2 \, u + 1}{u^{2 \, r + 4} - 2 \, u^{r + 2} + 1}}}$

Function Expansion¶

In [40]:

%%time
sum_of_remaining_nodes_func = fast_callable(G_r(u, v).diff(v, 1)(v=1)\
                                                     .simplify_rational(),
                                            vars=(u, r))(u_asy, r)

CPU times: user 8.2 s, sys: 151 ms, total: 8.35 s
Wall time: 7.76 s

In [41]:

sum_of_remaining_nodes_func

Out[41]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\left(\frac{{\left(r + 2\right)} {\left(r + 1\right)}}{4 \, {\left(r^{2} + 3 \, r + 2\right)} \sqrt{r^{2} + 2 \, r + 1}}\right) Z^{\frac{1}{2}} - \frac{{\left(4 \, r^{2} + 12 \, r + \frac{{\left({\left(2 \, r + 3\right)}^{2} - 2 \, r^{2} - 8 \, r - 9\right)} \sqrt{r^{2} + 2 \, r + 1}}{r + 2} + 8\right)} {\left(r + 2\right)}}{8 \, {\left(r^{2} + 3 \, r + 2\right)} \sqrt{r^{2} + 2 \, r + 1}} + \frac{{\left(r^{3} + 5 \, r^{2} + 8 \, r + 4\right)} {\left(r + 2\right)} {\left(r + 1\right)}}{2 \, {\left(r^{2} + 3 \, r + 2\right)}^{2} \sqrt{r^{2} + 2 \, r + 1}} + O\!\left(Z^{-\frac{1}{2}}\right)$

Coefficient Extraction¶

In [42]:

sum_of_remaining_nodes_asy = sum_of_remaining_nodes_func\
                             ._singularity_analysis_(n, 
                                                     zeta=1/4, 
                                                     precision=2)
sum_of_remaining_nodes_asy

Out[42]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\left(\frac{{\left(r + 2\right)} {\left(r + 1\right)}}{4 \, \sqrt{\pi} {\left(r^{2} + 3 \, r + 2\right)} \sqrt{r^{2} + 2 \, r + 1}}\right) 4^{n} n^{-\frac{1}{2}} + O\!\left(4^{n} n^{-\frac{3}{2}}\right)$

In [43]:

expected_tree_size = sum_of_remaining_nodes_asy / shifted_catalan_asy
expected_tree_size

Out[43]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\left(\frac{{\left(r + 2\right)} {\left(r + 1\right)}}{{\left(r^{2} + 3 \, r + 2\right)} \sqrt{r^{2} + 2 \, r + 1}}\right) n + O\!\left(1\right)$

In [44]:

expected_tree_size.map_coefficients(lambda ex: ex.subs({1/sqrt(r^2 + 2*r + 1) : 1/(r+1)}).factor())

Out[44]:

$\newcommand{\Bold}[1]{\mathbf{#1}}\left(\frac{1}{r + 1}\right) n + O\!\left(1\right)$