# 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]:

### 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]:
In [7]:
phi(T, z, t)(z,t).taylor((z,0), (t,0), 4)

Out[7]:
In [8]:
phi(phi(T, z, t), z, t)(z,t).taylor((z,0), (t,0), 4)

Out[8]:

## Motivation: Be Wary of Symbolics!¶

In [9]:
expr = 1/(1 + 1/x)
expr

Out[9]:
In [10]:
expr.series(x==0, 10)

Out[10]:
In [11]:
expr.simplify_rational()

Out[11]:
In [12]:
expr.simplify_rational().series(x==0, 6)

Out[12]:

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

In [13]:
expr = (2*exp(x) + exp(-x))
expr

Out[13]:
In [14]:
expr.factor()

Out[14]:
In [15]:
expr = (x + sqrt(x))/(x - sqrt(x))
expr

Out[15]:
In [16]:
expr.factor()

Out[16]:
In [17]:
(x + sqrt(x)).factor()

Out[17]:

## 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]:
A

Out[20]:
Asymptotic Ring <QQ^n * n^QQ * log(n)^QQ> over Symbolic Constants Subring
In [22]:
A.an_element()

Out[22]:
In [23]:
42 * (1/3)^n * n^(5/2) * log(n)^2  +  O((1/3)^n * n * log(n)^(1/2))

Out[23]:

### Straightforward Computations¶

In [24]:
log(n+1)

Out[24]:
In [25]:
1/(1 + n)

Out[25]:
In [26]:
(1 + 1/n)^n

Out[26]:

### 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]:
In [29]:
catalan_asy = binomial_2n_n / (n+1)
catalan_asy

Out[29]:

### 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]:

## 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]:
In [39]:
G_r(u, v).diff(v, 1)(v=1).simplify_rational()

Out[39]:

### 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]:

### 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]:
In [43]:
expected_tree_size = sum_of_remaining_nodes_asy / shifted_catalan_asy
expected_tree_size

Out[43]:
In [44]:
expected_tree_size.map_coefficients(lambda ex: ex.subs({1/sqrt(r^2 + 2*r + 1) : 1/(r+1)}).factor())

Out[44]: