<figure><img src="image_1.png" /></figure>

HAL
open science

# Holonomic equations and efficient random generation of binary trees

Pierre Lescanne

► To cite this version:

Pierre Lescanne. Holonomic equations and efficient random generation of binary trees. Discrete Mathematics and Theoretical Computer Science, 2023, vol. 25:2. hal-03674690v6

HAL Id: hal-03674690

https://hal.science/hal-03674690v6

Submitted on 22 Jan 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
# Holonomic equations and efficient random generation of binary trees

Pierre Lescanne

LIP, École Normale Supérieure de Lyon, France

received 13th Feb. 2023; accepted 7th July 2023.

Holonomic equations are recursive equations which allow computing efficiently numbers of combinatoric objects. Rémy showed that the holonomic equation associated with binary trees yields an efficient linear random generator of binary trees. I extend this paradigm to Motzkin trees and Schröder trees and show that despite slight differences my algorithm that generates random Schröder trees has linear expected complexity and my algorithm that generates Motzkin trees is in $O(n)$ expected complexity, only if we can implement a specific oracle with a $O(1)$ complexity. For Motzkin trees, I propose a solution which works well for realistic values (up to size ten millions) and yields an efficient algorithm.

**Keywords:** combinatorics, random generation, Motzkin number, Catalan number, binary tree, unary-binary tree

## 1 Introduction

In this paper, I am interested in efficient algorithms for randomly generating several sorts of binary trees. For this, I consider recurrences defining sequences of integer coefficients and, more precisely, I am interested in specific recurrences $F_n$ called “holonomic recurrence” where roughly speaking, “holonomic” means that $F_{n+s}$ is a combination, using polynomials in $n$, of the $F_i$'s, for $n \le i \le n+s$. More precisely, (see Flajolet and Sedgewick's book [15], Appendix B.4) the coefficients fulfill the following recurrence:

$$P_s(n)F_{n+s} + P_{s-1}(n)F_{n+s-1} + \dots + P_0(n)F_n = 0$$

for some $n \ge n_0$, where the $P_j(n)$ are polynomials in $n$. This kind of recurrence is called a *P-recurrence*. For instance, for Catalan numbers:

$$C_n - \sum_{k=0}^{n-1} C_k C_{n-k-1} = 0$$

is the classical recurrence that is used in general to defined them, but it is not a *P-recurrence*, whereas

$$(n+1)C_n - 2(2n-1)C_{n-1} = 0$$

is the *P-recurrence*, which will be considered later on. Notice that initial values should be added to this *P-recurrence*. This will be considered in the paper for each specific case.
<figure><img src="image_2.png" /><figcaption>Fig. 1: Rémy's right insertion of a leaf</figcaption></figure>

In this paper, I consider three families of binary trees (planar binary trees or Catalan trees, Motzkin trees aka unary-binary trees, Schröder trees) and their random generation. It turns out that holonomic recurrences play a key role in the design of efficient random generation algorithms.

The three examples: binary trees, Motzkin trees, Schröder trees are interesting because they have different holonomic equations, one (Catalan numbers) has one term on the right, one (Motzkin numbers) has a sum of two terms on the right and one (Schröder numbers) has a subtraction of two terms, on the right. These yield different random generation algorithms, as this will be explained further in this paper.

This paper is associated with a library of programs written in Haskell and in Python, available on GitHub. The reader who wants to read a better rendering of the programs of this paper is invited to get the version on my web page or the version available on GitHub.

## 2 Random binary trees

Rémy's algorithm [24] for generation of random binary trees is of linear complexity, i.e., $O(n)$. It is based on a constructive proof of the holonomic equation [25]:

$$ (n + 1)C_n = 2(2n - 1)C_{n-1} $$

Here "constructive" means that an explicit bijection between objects counted by the both sides of the equation is provided. In the case of binary trees, this holonomic equation is very peculiar since $C_n$ times a polynomial in $n$ is equal to $C_{n-1}$ times a polynomial in $n$. We will see that this is not the case for Motzkin numbers and Schröder numbers, but the paradigm can be extended. Rémy's algorithm is described by Knuth in [21] § 7.2.1.6 (pp. 18-19) and works on extended binary trees, or just binary trees in which we distinguish internal nodes and external nodes or leaves. The idea of the algorithm is that a random binary tree can be built by iteratively and randomly drawing an internal node or a leaf in a random binary tree and inserting, between it and its parent a new internal node and a new leaf either on the left or on the right (see Figure 1). An insertion is also possible at the root. In this case, the new inserted node becomes the root. The root can be seen as the child of a hypothetical node

A binary tree of size $n$ has $n-1$ internal nodes and $n$ leaves. We label binary trees with numbers between 0 and $2n-2$ such that internal nodes are labeled with odd numbers and leaves are labeled with even numbers. Inserting a node in a binary tree of size $n$ requires drawing randomly a number between 0 and $4n-3$. This process can be optimized by representing a binary tree as a list (a vector in Haskell), an idea sketched by Rémy and described by Knuth. In this vector, even values are for internal nodes and odd values are for leaves. The root is located at index 0. The left child of an internal node with label $2k+1$ is located at index $2k+1$ and its right child is located at index $2k+2$. Here is a vector representing a binary tree with 10 leaves and its drawing.
<table><tr><td>indices</td><td>0</td><td>1</td><td>2</td><td>3</td><td>4</td><td>5</td><td>6</td><td>7</td><td>8</td><td>9</td><td>10</td><td>11</td><td>12</td><td>13</td><td>14</td><td>15</td><td>16</td><td>17</td><td>18</td></tr><tr><td>values</td><td>1</td><td>13</td><td>0</td><td>2</td><td>5</td><td>9</td><td>7</td><td>8</td><td>4</td><td>11</td><td>17</td><td>12</td><td>10</td><td>15</td><td>3</td><td>16</td><td>14</td><td>18</td><td>6</td></tr></table>

<figure><img src="image_3.png" /><figcaption>This tree was built by inserting the node 17 together with the leaf 18 in the following vector.</figcaption></figure>

<table><tr><td>indices</td><td>0</td><td>1</td><td>2</td><td>3</td><td>4</td><td>5</td><td>6</td><td>7</td><td>8</td><td>9</td><td>10</td><td>11</td><td>12</td><td>13</td><td>14</td><td>15</td><td>16</td></tr><tr><td>values</td><td>1</td><td>13</td><td>0</td><td>2</td><td>5</td><td>9</td><td>7</td><td>8</td><td>4</td><td>11</td><td>6</td><td>12</td><td>10</td><td>15</td><td>3</td><td>16</td><td>14</td></tr></table>

which codes the tree

<figure><img src="image_4.png" /></figure>

This was done by drawing a node (internal node or leaf, here the node with label 6, right child of the node with label 9) and a direction (here right) and by inserting above this node a new internal node (labeled 17) and, below the new inserted internal node, a new leaf of the left (labeled 18). This double action (inserting the internal node and attaching the leaf) is done by choosing a number in the interval [0..33] (in general, in the interval [0..(4n - 3)]). Assume that in this case the random generator returns 21. 21 contains two pieces of information : its parity (a boolean) and floor of its half. Half of 21 is 10, which tells that the new node 17 must be inserted above the 11th node (in the vector) namely 6. Since 21 is odd, the rest of the tree (here reduced to the leaf 6) is inserted on the right (otherwise it would be inserted on the left). A new leaf 18 is inserted on the left (otherwise it would be inserted on the right).
Consider the same tree and suppose that the random value is 8. Half of 8 is 4. Hence the new internal node labeled by 17 is inserted above the node labeled by 5

<figure><img src="image_5.png" /></figure>

and, since 8 is even, the rest of the tree is inserted on the left and a new leaf (labeled 18) is inserted on the right.

<figure><img src="image_6.png" /></figure>

The algorithm (Figure 2) works as follows. If $n = 0$, Rémy's algorithm returns the vector starting at 0 and filled with anything, since the whole algorithm works on the same vector with the same size. In general, say that, for $n - 1$, Rémy's algorithm returns a vector $v$ (vector is the concept used in Haskell for arrays that can be changed in place). In our Haskell implementation the function yields an object of type `Gen (Vector Int)` which returns vector and carries a hidden random number generator. One accesses to the generator by `get` and stores the new generator by `put`. One draws a random integer $x$ between 0 and $4n - 3$. Let $k$ be half of $x$. In the vector $v$ one replaces the $k^{th}$ position with $2n - 1$ and one appends two elements, namely the $k^{th}$ item of $v$ followed by $2n$ if $x$ is even and $2n$ followed by the $k^{th}$ item of $v$ if $x$ is odd.

The algorithm builds a uniformly random *decorated binary tree*, i.e., a binary tree with its leaves numbered 0, 2, ..., $2n$. We notice that the construction of a tree with such labels is unique, the labels of the internal nodes are a consequence of the construction, hence are deduced from the labels of the leaves. If
1 rbt :: Int -> Int -> Gen (Vector Int)
2 rbt seed 0 = do put (mkStdGen seed)
3 return initialVector // [(0,0)]
4 rbt seed n =
5 do v <- rbt seed (n-1)
6 generator <- get
7 let (rand, newGenerator) = randomR (0::Double, 1) generator
8 put newGenerator
9 let x = floor (rand * fromIntegral (4*n-3))
10 -- x is a random value between 0 and 4n-3 --
11 k = x 'div' 2
12 case even x of
13 True -> return(v // [(k,2*n-1),(2*n-1,v!k),(2*n,2*n)])
14 False -> return(v // [(k,2*n-1),(2*n-1,2*n),(2*n,v!k)])

<caption>Fig. 2: Haskell program for Rémy’s algorithm</caption>

we ignore the leaves, we get a uniform distribution for the undecorated binary trees (i.e., with no labels on the leaves).

In the program, `rands` is a vector of random floating numbers between 0 and 1.

## 3 Motzkin trees

*Motzkin trees* are also called *unary-binary trees*. This paper proposes an algorithm for random generation of Motzkin trees. The algorithm takes the same paradigm as this of Rémy’s linear algorithm for random generation of *binary trees* [24]. Assume *n* is the size of the trees. My algorithm for random generation of Motzkin trees is based on a bijective proof due to Dulucq and Penaud [13] of the inductive equality:

$$ (n + 2)M_n = (2n + 1)M_{n-1} + 3(n - 1)M_{n-2} $$

where the $M_n$'s are the Motzkin numbers. At some point of the algorithm, an “oracle” choices which subprogram to call, based on $M_{n-1}$ and $M_{n-2}$. Since $M_{n-1}$ and $M_{n-2}$ are big numbers, this induces potentially a not $O(1)$ computation. A preprocessing allows a constant time computation for the oracle.

## 4 Motzkin numbers and Motzkin trees

The $n^{th}$ Motzkin number $M_n$ is the number of different ways of drawing non-intersecting chords between $n$ points on a circle (not necessarily touching every point by a chord). Motzkin numbers count also well parenthesized expressions with a constant c, called *Motzkin words*. They are words of length $n$ in the language generated by the grammar $M$.

$$ M = \varepsilon | c M | (M) M $$
<figure><img src="image_7.png" /><figcaption>Fig. 3: The 9 non-intersecting chords between 4 points on a circle</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Fig. 4: The 9 parenthesized words with constant c.</figcaption></figure>

The bijection between sets of non intersecting chords and well parenthesized words with constant c is
as follows: first one numbers nodes on the circle counterclockwise, as follows:

<figure><img src="image_9.png" /></figure>

A position at the beginning of a chord on the circle corresponds to an opening parenthesis. A position at the end of a chord on the circle corresponds to a closing parenthesis. A position which is neither of those corresponds to the constant c.

Motzkin numbers count also routes in the upper quadrant from (0, 0) to (0, 4) with move *up*, *down* and *straight*.

The bijection is as follows: an opening parenthesis corresponds to an *up*, a closing parenthesis corresponds to a *down* and the constant *c* corresponds to a *straight*.

Motzkin number $M_n$ counts also the number in unary-binary planar trees with $n$ edges, that are tree structures with nodes of arity one or two and with $n$ edges. Let us call this number $n$ of edges the *size* of the Motzkin tree. Notice that the number of nodes of a Motzkin tee of size $n$ is $n + 1$, i.e., the size plus one. Figure 6 gives the trees for $n = 4$. The bijection $f$ from well parenthesized expressions with constant
<figure><img src="image_10.png" /><figcaption>Fig. 5: The 9 routes in the upper right quadrant from (0, 0) to (0, 4) with move *up*, *down* and *straight*</figcaption></figure>

c to Motzkin trees is as follows. Its inverse $f^{-1}$ is also given.

$$
\begin{align*}
f[\varepsilon] &= \bullet & f^{-1}[\bullet] &= \varepsilon \\
f[\mathbf{c} w] &= f[w] & f^{-1}\left[\begin{smallmatrix} \uparrow \\ t \end{smallmatrix}\right] &= \mathbf{c} f^{-1}[t] \\
f[(w_1) w_2] &= f[w_1] \bigwedge f[w_2] & f^{-1}\left[\begin{smallmatrix} \uparrow \\ t_1 \\ \downarrow \\ t_2 \end{smallmatrix}\right] &= (f^{-1}[t_1]) f^{-1}[t_2]
\end{align*}
$$

Motzkin numbers fulfill the equation:

$$
M_{n+1} = M_n + \sum_{i=0}^{n-1} M_i M_{n-1-i}
$$

5 Dulucq-Penaud bijection proof

As I said, Motzkin numbers fulfill the holonomic equation [15]:

$$
(n + 2)M_n = (2n + 1)M_{n-1} + 3(n - 1)M_{n-2}.
$$

Together with the equalities $M_0 = 1$ and $M_1 = 1$, Motzkin numbers can be computed and form the sequence A001006 in the *online encyclopedia of integer sequences* [19]. In this section, I present Dulucq and Penaud’s proof of this equality [13]. This proof relies on the exhibition of a bijection between the objects counted by the left-hand side and those counted by the right-hand side. The first idea is to consider specific binary trees called *slanting binary trees* and divide those trees into 7 subclasses.

5.1 Slanting binary trees

Following Dulucq and Penaud, I represent *Motzkin trees* as specific *binary trees* in which leaves □ are added. In such binary trees, only the three first configurations below are allowed and the fourth and rightmost one is not.

<figure><img src="image_11.png" /></figure>
<figure><img src="image_12.png" /><figcaption>Fig. 6: 9 Motzkin trees with 4 edges and their slanting trees</figcaption></figure>
The first configuration corresponds to a binary node, the second configuration corresponds to a unary node and the third configuration corresponds to an end node in the classical presentation (for instance in Figure 6). I call such trees slanting trees or binary slanting trees, from the french arbres binaires penchés.

Figure 6 shows the 9 Motzkin trees with 4 edges and their corresponding slanting trees. Let us label each node of a slanting tree with a number between 1 and $2m+1$, where $m$ is the number of internal nodes of the slanting tree. Let us call such a labeled tree a *labeled slanting tree*. Now consider labeled slanting trees with one marked leaf. Let us call it a *leaf-marked slanting tree*. Below there is a labeled slanting tree of size 4 and a leaf-marked labeled slanting tree, where the mark is on the leaf labeled 8.

<figure><img src="image_13.png" /></figure>

This corresponds to the vector $[3, 0, 2, 5, 4, 7, 1, 6, 8]$.

How nodes and leaves are labeled by numbers will be explained below and is essentially like binary trees. Like binary trees, just notice that *internal nodes* have odd labels and *leaves* have even labels. From now on, let us forget the labels, but let us mark one of the leaves. In such trees with a marked leaf, we can distinguish 7 general patterns of subtrees containing the marked leaf (Figure 7 first column). The marked leaf is denoted by a star in a square, namely $\star$. In the first group of 4 slanting trees, there are the patterns where the marked leaf is a right child and in the second group of 3 slanting trees, there are the patterns where the marked leaf is a left child, hence, due to the constraints on slanting trees, the other child (a right child) is a leaf as well.

Let us call *node-marked* a slanting tree in which one internal node is marked. Let us call *marked* tree, a slanting tree in which either a leaf or an internal node is marked.

## How many leaves in a Motzkin tree?

The slanting tree associated with a Motzkin tree of size $n$ (number of its edges) has $n+2$ leaves. This can be shown by induction.

**Basic case:** If the Motzkin tree is $\bullet$, its size is 0 and its associated slanting tree $\Box$ has 2 leaves.

**Adding a unary node:** Assume we add a unary node above a Motzkin tree $t$ of size $n$, this yields a Motzkin tree $t'$ of size $n+1$. The slanting tree associated with $t'$ has $n+2$ leaves (the number of the leaf of the slanting tree $u$ associated with $t$) plus a new one added, then all together $n+3$.

$$t \rightarrow \overset{\oplus}{t}$$

$$u \rightarrow u \Box$$

**Adding a binary node:** Assume we add a binary node above two Motzkin trees $t_1$ and $t_2$ of size $n_1$ and
<table><caption>Fig. 7: The 7 patterns of leaf-marked slanting trees</caption><thead><tr><th></th><th>leaf-marked slanting trees</th><th>marked slanting trees</th><th>Choice and node-marked slanting trees</th></tr></thead><tbody><tr><td>1.</td><td></td><td></td><td></td></tr><tr><td>2.</td><td></td><td></td><td>LR,</td></tr><tr><td>3.</td><td></td><td></td><td></td></tr><tr><td>4.</td><td></td><td></td><td>RR,</td></tr><tr><td>5.</td><td></td><td></td><td>RL,</td></tr><tr><td>6.</td><td></td><td></td><td>LL,</td></tr><tr><td>7.</td><td></td><td></td><td></td></tr></tbody></table>

$n_2$, this yields a Motzkin tree $t'$ of size $n_1 + n_2 + 2$. The slanting trees associated with $t_1 \stackrel{\cdot}{\rightleftharpoons} t_2$ have $n_1 + n_2 + 4$ leaves.

## 5.2 A taxonomy of slanting trees

Recall Rémy's algorithm which consists in inserting, at a marked position (internal node or leaf), a leaf in a marked tree (see Figure 1). After insertion, the formerly marked node becomes unmarked and the inserted leaf becomes marked. Here, since we are interested in Motzkin trees, we insert a leaf on the right, above the marked node in the marked slanting tree and like in *Rémy's algorithm*, a leaf insertion on a marked slanting tree is performed, but unlike *Rémy's algorithm* the insertion is performed only on the right and a leaf-marked slanting tree is produced. This corresponds to what is done to *pattern1*, *pattern3*
and *pattern4* in the middle column of Figure 7. The leaf is inserted at a marked position in the marked tree of the middle column producing the leaf-marked tree of the left column. But as we will see for the other patterns, there are other ways to increase a slanting tree when it does not correspond to one of these three patterns.

# 6 The bijection

Beside the fact one works on slanting trees with constraints, what makes also the random generation of Motzkin trees trickier than Rémy's algorithm is the structure of the holonomic equation:

$$ (n + 2)M_n = (2n + 1)M_{n-1} + 3(n - 1)M_{n-2} $$

when compared to the equation:

$$ (n + 1)C_n = 2(2n - 1)C_{n-1} $$

First, if we use a construction based on that equation, a Motzkin tree of size $n$ can be built from a Motzkin tree of size $n - 1$ or from a Motzkin tree of size $n - 2$. Thus there are at least two cases to consider. Actually 7 cases as we will see, since the construction of a Motzkin tree of size $n - 1$ splits in 4 cases and the construction of a Motzkin tree of size $n - 2$ splits in 3 cases. Notice that $M_n$ counts both the number of Motzkin trees of size $n$ and the slanting trees with $n + 2$ leaves.

## Interpreting the holonomic equation

We conclude that $(n+2)M_n$ counts the number of *leaf-marked* slanting trees of size $n$, that $(2n+1)M_{n-1}$ counts the number of *marked* slanting trees of size $n - 1$ and that $(n-1)M_{n-2}$ counts the number of *node-marked* slanting trees of size $n - 2$. Therefore looking at the equation, we see that we should be able to build a leaf-marked slanting tree of size $n$ from either a marked slanting tree of size $n - 1$ or from a pair made of an item which can take one of three values and of a *node-marked* slanting tree of size $n - 2$. Let us see how Dulucq and Penaud propose to proceed.

## A taxonomy of leaf-marked slanting trees

Leaf-marked slanting trees can be sorted according to the position of their mark. This is done in the first column of Figure 7. This column has two parts.

### The upper part

In the upper part, we have four patterns in which the marked leaf is a right child. Let us call them *pattern1*, *pattern2*, *pattern3* and *pattern4* (Figure 7). Three of them *pattern1*, *pattern3* and *pattern4* are obtained by Rémy's right insertion of a leaf in a marked slanting tree. The other *pattern2* is not. Indeed if the marked leaf is removed, the tree that is obtained has a leaf on the left and a node on the right, which is forbidden. This pattern will be obtained another way.

### The lower part

In the lower part, there are three patterns which correspond to the case where the marked leaf is a left child, hence the sibling of a leaf (Figure 9); *pattern7* cannot be obtained by a right insertion of a leaf (Figure 7). I annotate it with a 😊. But Dulucq and Penaud noticed that pattern ☐ among marked slanting trees is not taken into consideration. Thus they propose to associate this pattern ☐ with *pattern7*, as shown in Figure 7.
<figure><img src="image_14.png" /><figcaption>**Fig. 8:** Contribution of marked slanting tree of size $n-1$. 😊 for *pattern7* marks a specific case explained in Section *The lower part*</figcaption></figure>

## The bijection by cases

### Case $(2n+1)M_{n-1}$

The previously explained contribution to leaf-marked slanting trees of size $n$ from marked slanting tree of size $n-1$ is summarized in Figure 8. All the patterns of marked slanting trees are taken into account.

### Case $3(n-1)M_{n-2}$

Let us now look at the three remaining patterns: *pattern2*, *pattern5* and *pattern6*; forming the lines of Figure 9. Those three patterns have the same model, namely an internal node with two children: one child is an internal node and the other child is an internal node whose children are two leaves, one of which is marked, the other is not. Depending on the position of the marked leaf, we distinguish three cases.

* **LL** corresponds to the case where the marked node is on the left of the top node and on the left of its parent node.

* **LR** corresponds to the case where the marked node is on the left of the top node and on the right of its parent node.

* **RL** corresponds to the case where the marked node is on the right of the top node and on the left of its parent node.

One notices that there is no case **RR**, because this would correspond to *pattern4* considered in the previous section. As a matter of fact, the three cases **LL**, **LR** and **RL** correspond to the multiplicative factor 3 in the holonomic equation.
<figure><img src="image_15.png" /><figcaption>Fig. 9: Contribution of node-marked slanting tree of size $n-2$</figcaption></figure>

## Forgetting the marks on leaves

As Rémy noticed for binary trees, since we generate the leaf-marked slanting trees of size $n$ uniformly, we also get a uniform distribution of slanting trees of size $n$. Thus we can forget the marks, which we do in the concrete algorithm.

# 7 A concrete algorithm for random generation of Motzkin trees

The Haskell program of Figure 10 presents the algorithm for random generation of Motzkin trees. In what follows, I make no distinction between the algorithm and the program and I consider the Haskell program as an executable specification. The main function is called `rMt` and returns an object of type `Gen (Vector Int)` like `rbt` in Figure 2.

Assume that there is a function `motzkin` that returns the $n^{th}$ Motzkin number. Like for Rémy's algorithm, one represents a labeled slanting tree by a vector. In this vector, odd labels are for internal and even labels are for leaves. Notice that the algorithm preserves two properties:

1. The vector codes a slanting tree.

2. The vector of a Motzkin tree of size $n$ has a length $2n + 3$.

In order to choose which case to consider, namely $(2n+1)M_{n-1}$ (case1) or $3(n-1)M_{n-2}$ (case2), the algorithm `rMt` requires a random value between 0 and $(n+2)M_n$ which we call $r$. If $r$ is less than or equal to $(2n+1)M_{n-1}$, we are in case1, else we are in case2. Said otherwise, given a random value $c$ between 0 and 1 if $c \le \frac{(2n+1)M_{n-1}}{(n+2)M_n}$ we choose case1, if not we choose case2.

More abstractly, let us forget Motzkin numbers and consider an *oracle* which, given a random number $r$ between 0 and 1 and an $n$ chooses between case1 and case2. If the oracle runs in constant time and returns a boolean according to the distribution given by the above inequality then the algorithm has a linear time complexity and returns a random Motzkin tree distributed according to the above inequality.
rMt :: Int -> Int -> Gen (Vector Int)
rMt seed 0 = do put (mkStdGen seed)
             return $ initialVector // [(0,1),(1,0),(2,2)]
rMt seed 1 = do put (mkStdGen seed)
             return $ initialVector // [(0,1),(1,3),(2,0),(3,2),(4,4)]

rMt seed n =
  do generator <- get

let (rand, newGenerator) = randomR (0::Double, 1) generator
  put newGenerator
  case oracle n rand of
    True -> case1 seed n
    False -> case2 seed n

case1 seed n =

do generator <- get

let (rand, newGenerator) = randomR (0::Double, 1) generator
  k = floor (rand * (fromIntegral (2*n)))
v <- rMt seed (n-1)
put newGenerator
case odd k || odd (v!k) || odd (v!(k-1)) of
  True -> return $ v // [(k,2*n+1),(2*n+1,v!k),(2*n+2,2*n+2)]
  False -> return $ v // [(k-1,2*n+1),(2*n+1,v!(k-1)),(2*n+2,2*n+2)]

case2 seed n =

do generator <- get

let (rand, newGenerator) = randomR (0::Double, 1) generator
  r = floor (rand * (fromIntegral (3*n-6)))
  k = r 'div' 3
c = r 'rem' 3
v <- rMt seed (n-2)
put newGenerator
case c < 2 of
  True -> return $ v // [(2*k+1,2*n-1),(2*k+2,2*n+1),
                          (2*n-1,2*n),(2*n,2*n+2),
                          (2*n+1,v!(2*k+1)),(2*n+2,v!(2*k+2))]
  False -> return $ v // [(2*k+1,2*n-1),(2*k+2,2*n+1),
                            (2*n-1,v!(2*k+1)),(2*n,v!(2*k+2)),
                            (2*n+1,2*n),(2*n+2,2*n+2)]

<caption>Fig. 10: Haskell program for random generation of Motzkin trees</caption>
* **case1:** one draws at random a leaf or an internal node in a slanting tree of size $n-1$. This means choosing at random an index $k$ in the vector $v$. We get *pattern7* if three conditions are fulfilled.

1. **The marked item, should be a right child.** This means that $k$ is even, since the left child of a node of index $2p+1$ is $2p+1$ and the right child of this node is $2p+2$.

2. **The marked item is a leaf.** This means that $v[k]$ is even, since leaves have even labels. Notice that Haskell uses the notation $v!k$ for our mathematical notation $v[k]$.

3. **The sibling item of the marked item is a leaf** (a left child by the way). This means that $v[k-1]$ is even.

In this case ($k$ is even, $v[k]$ is even and $v[k-1]$ is even) one inserts a node and a leaf as shown by Figure 8, which corresponds in the code to:

```haskell
v // [(k-1, 2*n+1), (2*n+1, v!(k-1)), (2*n+2, 2*n+2)]
```

In Haskell, the operator `//` updates vectors at once, it is called a **bulk update**. $(2*n+1, v!(k-1))$ means that the left child is a new node, at index $2*n + 1$, which points to the former value of $v!(k-1)$. The right child is a new leaf.

The other cases (*pattern1*, *pattern3*, *pattern4*) correspond to cases when one of $k$, $v[k]$ or $v[k-1]$ is odd. The update is then

```haskell
v // [(k, 2*n+1), (2*n+1, v!k), (2*n+2, 2*n+2)]
```

which corresponds to the first lines of Figure 8.

* **case2:** In this case we consider a random node-marked slanting trees of size $n-2$ and a random values among **LR**, **RL**, **LL**. For that we draw a number $r$ between 0 and $3n-6$, from which $r \div 3$ gives a random number between 0 and $n-2$ (a random node) and $r$ mode 3 yields a random number among 0, 1 and 2. We notice that **LR** and **LL** correspond to the same transformation, while **RL** corresponds to another transformation. In each case one adds four nodes, with labels $2n-1, 2n, 2n+1$ and $2n+2$. Thus,

```haskell
v // [(2*k+1, 2*n-1), (2*k+2, 2*n+1), (2*n-1, 2*n), (2*n, 2*n+2),
     (2*n+1, v!(2*k+1)), (2*n+2, v!(2*k+2))]
```

is the transformation for **LR** and **LL** and

```haskell
v // [(2*k+1, 2*n-1), (2*k+2, 2*n+1), (2*n-1, v!(2*k+1)), (2*n, v!(2*k+2)),
     (2*n+1, 2*n), (2*n+2, 2*n+2)]
```

corresponds to **RL**. We let the reader check that the code of Figure 10 corresponds to the pictures of Figure 9.
1 oracleMotzkin :: Int -> Double -> Bool

2 oracleMotzkin n rand =

3 let r = floor (rand * fromIntegral (fromIntegral (n+2) * (motzkin n)))

4 in r <= (fromIntegral (2*n+1)) * (motzkin (n-1))

<caption>Fig. 11: The plain oracle</caption>

1 oracleRatioM :: Int -> Double -> Bool

2 oracleRatioM n rand = case n <= 100004 of

3 True -> ratioM!n <= rand

4 False -> 0.66666666667 <= rand

<caption>Fig. 12: The approximation oracle</caption>

## 8 Linearity and oracle

The linearity of the algorithm depends on an oracle which should decide an inequality in constant time. In the implementation is the plain implementation of Figure 11, the complexity of the algorithm $rMt$ is not asymptotically linear, due to computations on big numbers.

However the inequality $rand < \frac{(2n+1)Mn-1}{(n+2)Mn}$ of the oracle can be approximated by a precomputed table which I call $ratioM$ (of actual size 10004 in my case), and which is used in the oracle in Figure 12. For sizes larger than this bound, I take just the value 0.66666666667, which is a good approximation, since the fraction goes down to the limit $\frac{2}{3}$. The actual implementation runs up to size $9 \times 10^6$ and never requires (in the benchmarks) more than 10 digits of precision (see next section). But despite the asymptotic linearity is not guaranteed, the algorithm is linear in practice.

For a better efficiency, the program $rMt$ can run without recursive calls. In this case, the stack of the calls is first computed and the construction of the vector is performed by popping the stack. I wrote a Python program, which generates random Motzkin trees of size 9 millions in less that 50 seconds on a laptop. In all the benchmarks, I checked that the difference between the drawn *rand* and the fraction is always larger than $10^{-10}$. Hence, digits that are not 6 are never checked, which shows that the oracle called *oracleRatioM* corresponds to the oracle called *oracleMotzkin* in this size interval $[0.. 9 \times 10^6]$ and runs in constant time.

## 9 Schröder trees

In this section, I define Schröder trees and a proof by bijection of the holonomic equation for defining numbers that count Schröder trees, aka *Schröder numbers* or *Schröder-Hipparchus numbers*. This proof is the translation for Knuth's definition of Schröder trees of this due to Foata and Zeilberger [17]. From this proof I derive a quasi linear algorithm for random generation of Schröder trees.

### 9.1 Definition of Schröder trees

We take the definition of Knuth [21] § 7.2.1.6 (pp. 41): "A Schröder tree is a binary tree in which every nonnull right link is colored either white or black". We represent black links with `\_` and white links
<figure><img src="image_16.png" /><figcaption>Fig. 13: The 11 Schröder trees with 4 leaves.</figcaption></figure>

with $\mathcal{R}$. But unlike Knuth, we follow Foata and Zeilberger [17] and we say that a Schröder tree has size $n$ if it has $n$ leaves, hence $n-1$ nodes. For instance, Schröder trees with 4 leaves are given in Figure 13.

Schröder trees are counted by numbers $S_n$, which form the sequence A001003 in the online encyclo-
pedia of integer sequences [19] and which fulfill the holonomic equations:

$$ (n + 1)S_{n+1} - 3(2n - 1)S_n + (n - 2)S_{n-1} = 0 $$

that we will use on the form

$$ (n + 1)S_{n+1} = 3(2n - 1)S_n - (n - 2)S_{n-1} $$

Notice that Foata and Zeilberger [17] use for the constructive proof:

$$ 3(2n - 1)S_n = (n + 1)S_{n+1} + (n - 2)S_{n-1} \quad (1) $$

Schröder trees count also ordinary trees with $n$ leaves and no node of degree one. One “shrinks” the white edges, merging the nodes that they connects (see [21] Exercise 66). For trees with 4 leaves, the reader can see the bijection in the lines of Figure 13 and in the odd lines of Figure 14. Trees in correspondence are on the same row. Schröder numbers count also partitions of polygons [14]. The bijection between ordinary trees with 4 leaves and no node of degree one and partitions of pentagons (polygons with 5 edges) is illustrated by Figure 14. First one associates with each tree a parenthesized expression in $\alpha, \beta, \gamma$ and $\delta$, in this order. Then each parenthesized expression is associated with a unique partition of the pentagon. For this purpose, a basic edge $\equiv$ is chosen and other edges are named clockwise by letters $\alpha, \beta, \gamma, \delta$. If one proceeds starting at edge $\alpha$, each diagonal closes a polygon (with less edges) in which the other edges
<figure><img src="image_17.png" /><figcaption>Fig. 14: 11 Trees having no node of degree one and 4 leaves, associated with 11 partitions of pentagons</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>Fig. 15: Insertion of a leaf in Foata and Zeilberger bijection</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>Fig. 16: Three impossible left insertions of leaves</figcaption></figure>

have been already named by an expression. Therefore one can associate a partition with each expression as this is presented in Figure 14. One can process the same way for hexagons, heptagons etc.

Didier Arquès and Alain Giorgetti have shown that planar rooted hypergraphs with vertices only of the outer face are also counted by Schröder numbers [3].

## 9.2 Foata and Zeilberger bijection

In identity (1):

* $(2n-1)S_n$ counts marked Schröder trees of size $n$, where nodes and leaves can be marked, therefore $3(2n-1)S_n$, counts pairs made of one of three values and of a marked Schröder tree,

* $(n+1)S_{n+1}$ counts leaf-marked Schröder trees of size $n+1$,

* $(n-2)S_{n-1}$ counts node-marked Schröder trees of size $n-1$.

In order to prove identity (1), Foata and Zeilberger build a one-to-one function which associates, with a pair of a number 0, 1 or 2 (coded by Foata and Zeilberger as $L_1$, $L_2$ and $R_1$) and a marked Schröder tree of size $n$, either a leaf-marked Schröder tree of size $n+1$ or a node-marked Schröder tree of size $n-1$.

Foata and Zeilberger bijection [17] that proves identity (1) is based on constructions similar to Rémy's insertion (Figure 1), with the difference that there are two left insertions, one that corresponds to a black link and one that corresponds to a white link (Figure 15). Due to the constraints on white links in Schröder trees, there are three patterns (Figure 16), with impossible insertion, when one tries to insert a left leaf to become a left sibling of a leaf connected with a white link. Two patterns ($a$ and $c$ in Figure 16) can be "recovered". The third pattern ($b$ in Figure 16) is shrunk toward a Schröder tree of size $n-1$. This will correspond to a fail in the algorithm.

Similarly there are two unreachable cases by left or right insertion. They are patterns which are the results of a right or left insertion on a leaf which are a right child by a white link.
<figure><img src="image_20.png" /><figcaption>Fig. 17: Two unreachable</figcaption></figure>

Let us call $L_1, L_2$ and $R_1$ the three labels. Therefore Figure 18 gives the correspondence between pairs of a label from $\{L_1, L_2, R_2\}$ and of a marked Schröder tree, with either a leaf-marked tree of size $n+1$ or a node marked tree of size $n-1$.

**Cases $L_1$ and $R_1$.** One inserts a leaf with a black link. This is exactly like Rémy's insertion.

**Case $L_2$.** This case deals with a left insertion with a white right link.

**First line:** The mark is on a node. The left insertion with a white right link is possible.

**Second line:** The mark is on a left leaf and the right link is black. The left insertion yields a forbidden pattern (a right white link toward a leaf). But by twisting the tree and swapping the colors of the right links, one reaches the first unreachable pattern of Figure 17.

**Third line:** The mark is on a left leaf and the right link is white. Then the left child is not a leaf. The left insertion is forbidden as well. Therefore by removing the left leaf one gets a node marked tree of size $n-1$.

**Fourth line:** The mark is on a right leaf. The left insertion yields a forbidden pattern, but by swapping the colors of the right link one gets the second unreachable pattern of Figure 17.

# 10 A concrete algorithm for random generation of Schröder trees

The program of Figure 19 presents the algorithm for random generation of Schröder trees. Like in previous cases one uses a vector (an array) of size $2n+1$. But in addition to the indices for the next nodes, one adds a boolean. This boolean says that the right link that starts from the node corresponding to this index is white. Hence each cell of the array contains a pair $(k, b)$ where $k$ is an index and $b$ is a boolean. According to the constraints induced on the Schröder trees, the pair corresponding to a boolean *True* has, as a first component, an odd number, since this first component corresponds to a link to a node. Since the node is a right child, it is located at an even index. If these constraints are not fulfilled, the second component must be a *False*. Therefore in the program when she adds a pair $(k, True)$ at a position $m$, she has to check that $k$ is odd and $m$ is even. In another hand, if she writes $(m, (k, True))$ for such an operation, she has to check that it is of the form $(2p, (2q+1, True))$. This constraints is an invariant of the program. The case when Foata and Zeilberger produce a Schröder tree of size $n-2$ corresponds in my algorithm to a "failure", that is a "retry": the subprogram **body** is called with a new random generator.

There are six cases. Assume one draws a number $x$ between 0 and $6n-4$ and let us call $k$ the number $x \div 3$. The value $x \bmod 3$ discriminates among $L_1, L_2$ and $R_2$: 0 for $L_1$, 1 for $L_2$ and 2 for $R_1$. The cases $L_1$ and $R_1$ are easy. Assume that at the $k^{th}$ index the array contains $(h, b)$.

$L_1$ corresponds in the code to
<figure><img src="image_21.png" /><figcaption>Fig. 18: Foata and Zeilberger isomorphism in pictures</figcaption></figure>
type Gen = State StdGen

rst :: Int -> Int -> Gen (Vector (Int, Bool))

rst seed 0 = do put (mkStdGen seed)

return ( initV // [(0,(0, False ))])

rst seed n =

let g = rst seed (n-1)

body generator =

do v <- g

let (rand, newGenerator) = randomR (0::Double, 1) generator

x = floor (rand * fromIntegral (6*n-4))

-- x is a random value between 0 and 6n - 4

k = x 'div' 3

put newGenerator

case x 'mod' 3 of

0 {\-L1\-} -> return (v // [(k,(2*n-1, False )),
                 (2*n-1,(2*n, False )),
                 (2*n,v!k)]) -- L1

1 {\-L2\-} -> case odd(fst (v!k)) of

True {\- not a leaf\-} -> return (v// [(k,(2*n-1, False )),
                 (2*n-1,(2*n, False )),
                 (2*n,(fst (v!k), True ))])
False -> case odd k of

False {\- the other link is black\-} ->
return (v //[(k,v!(k+1)),

(k+1,(2*n-1, True )),
(2*n-1,v!k),(2*n,(2*n, False ))])

{\-Failure\-} True {\- the other link is white \-} -> body newGenerator

False {\- the leaf is on the right \-} ->
return (v //[(k,(2*n-1, True )),

(2*n-1,v!k),
(2*n,(2*n, False ))]

2 {\-R1\-} -> return (v // [(k,(2*n-1, False )),

(2*n-1,v!k),
(2*n,(2*n, False ))])

in do generator <- get
body generator

<caption>Fig. 19: Haskell program for random generation of Schröder trees.</caption>
```
(v // [(k, (2*n-1, False)), (2*n-1, (2*n, False)), (2*n, v!k)]
```

The links are black hence the booleans are *False*.

**L₂ and h is odd.** This means that the mark is on a node. This corresponds to the code:

```
v// [(k, (2*n-1, False)), (2*n-1, (2*n, False)), (2*n, (fst(v!k), True))]
```

Clearly, 2*n*, (fst(v!k), True) fulfills the constraints (2p, (2q+1, True)) since fst(v!k), which we called h, is odd.

**L₂ and h is even and k is odd and the second component of v! (k+1) is False:** This means that the mark is on a leaf (h is even) which is a left child (k is odd) and the link that goes to the sibling leaf is black. This corresponds to the code:

```
v// [(k, v! (k+1)), (k+1, (2*n-1, True)), (2*n-1, v!k), (2*n, (2*n, False))]
```

The reader may check that the code corresponds to the picture of Figure 18, line 2. (k+1, (2*n-1, True)) fulfills the constraints (2p, (2q + 1, True)) since k + 1 is even and 2n - 1 is odd.

**L₂ and h is even and k is odd and the other link is white and the second component of v! (k+1) is True:** This corresponds to a *failure*, and the program loops with a new random generator.

**L₂ and h is even and k is even:** This means that the mark is a leaf (h is even) which is a right child. This corresponds to the code:

```
v// [(k, (2*n-1, True)), (2*n-1, v!k), (2*n, (2*n, False))]
```

(k, (2*n-1, True)) fulfills the constraints (2p, (2q + 1, True)) since k is even and 2n - 1 is odd.

$R_2$ corresponds in the code to

```
v // [(k, (2*n-1, False)), (2*n-1, v!k), (2*n, (2*n, False))]
```

The links are black hence the booleans are *False*.

One may notice that in the array, only values at odd indices require to carry a boolean, and this boolean is required only when the first value of the pair is odd. This suggests a better data structure which may save little space, but I did not implement it.
<table><caption>Fig. 20: Benchmarks for Schröder trees</caption><thead><tr><td>size</td><td>time</td><td>ratio</td></tr></thead><tbody><tr><td>1 000</td><td>0.012s</td><td>0.024</td></tr><tr><td>5 000</td><td>0.031s</td><td>0.0288</td></tr><tr><td>10 000</td><td>0.064s</td><td>0.025</td></tr><tr><td>50 000</td><td>0.200s</td><td>0.0269</td></tr><tr><td>100 000</td><td>0.290s</td><td>0.02707</td></tr><tr><td>500 000</td><td>1.295s</td><td>0.027762</td></tr><tr><td>1 000 000</td><td>3.065s</td><td>0.027883</td></tr><tr><td>5 000 000</td><td>15.183s</td><td>0.0276378</td></tr><tr><td>10 000 000</td><td>30.738s</td><td>0.0275827</td></tr></tbody></table>

## 10.1 Ratio of failures and average linear complexity

The loop in `rst` can be sketched by the following program P:

```python
1 while True:
2   r = random.random(0,2)
3   if r == 0:
4     break
5   elif r == 1:
6     continue
7   else:
8     break
```

despite in `rst`, the loop happens less often. In P, one breaks with probability $2/3$, in the first run of the loop, with probability $2/3^2$ in the second run of the loop, with probability $2/3^3$ in the third run of the loop, etc. Therefore if the cost of executing the statements of the loop is $c$, the total average cost of the loop is:

$$c \sum_{i=1}^{\infty} \frac{2}{3^i} = c.$$

## 10.2 Benchmarks

Benchmarks are done on a Python program, with a `while` loop, which mimics the Haskell program. See Figure 20. A third column yields the ratio of number of fails over the size of the produced tree. Computations are run on a laptop with 4 Gb of memory.

# 11 Related Works

## On Motzkin trees

Our work was inspired by Jean-Luc Rémy's algorithm [24]. Laurent Alonso [1] proposed an algorithm for generating uniformly Motzkin trees. His method consists in generating the number $k$ of binary nodes with the correct probability law; he uses then standard techniques to generate a unary-binary tree with $k$
binary nodes among $n$ nodes. The number of trees with $k$ binary nodes is over-approximated by values that follow a binomial distribution: choosing $k$ is therefore done using random generations for a binomial law and rejections. Therefore his algorithm is *linear on the average*, with possible but extremely rare long sequences of rejection. Dominique Gouyou-Beauchamps and Cyril Nicaud [18] propose a random generation for color Motzkin trees which is linear on the average and Srečko Brlek et al. [9] propose an extension of Alonso's algorithm to generalized versions of Motzkin trees.

Axel Bacher, Olivier Bodini and Alice Jacquot [4, 20, 5] propose an algorithm with similar ideas. Especially their Figure 2 shares similarity with our Figure 9. There “operations” $G_3$, $G_4$ and $G_5$ are connected with our cases **RL**, **LR** and **LL** respectively. ⊠ corresponds to ⌢ and ⊗ corresponds to ⌣. Like Alonso’s, the algorithm they propose has a linear expected complexity, due to failures similar to those of our Schröder tree generation algorithm. Let us also mention generic Boltzmann’s samplers with exact-size which apply among others to Motzkin trees [7, 22] and generic algorithms [11, 23] with linear expected complexity.

Denise and Zimmermann [10] discuss what can be done on floating-point arithmetic when generating random structures. The authors focus on decomposable labeled structures [16] and address the problem of choice (which I call oracle), with a specific section on Motzkin trees.

## On Schröder trees

After studying the random generation of Motzkin trees and reading Foata and Zeilberger paper, I started the implementation of a random generation of Schröder trees, which turns out to be of expected linear complexity. Actually Laurent Alonso, René Schott and Jean-Luc Rémy [2] proposed another linear algorithm for random generation of Schröder trees, on the same principle as Alonso’s algorithm for the generation of Motzkin trees. Like this quasi-linear algorithm, it proceeds in two steps: first, it chooses randomly the number $k$ of leaves with an adequate probability and a rejection technique, second, it generates a random Schröder tree with $k$ leaves. In comparison, my algorithm is direct. I deal only with Schröder trees, not with Schröder trees with $k$ leaves. I am very closed to Rémy’s algorithm and to my algorithm for random generation of Motzkin trees.

# 12 Conclusion

Generating Motzkin trees and Schröder trees has many potential applications [12, 6, 2]. My algorithm for generation of Motzkin trees has a simple code and is linear and my algorithm for Schröder trees is direct, which means it deals only with Schröder trees. Among the possible extensions of my method which could be explored, there is the generation of extended versions of Motzkin structures like Motzkin trees with colored leaves [18] or Motzkin paths with $k$ long straights [9]. $k = 1$ corresponds to Motzkin paths and $k = 2$ to Schröder paths. On another hand, Bracucci et al. [6] study a family of sets of permutations: $M_1$, $M_2$, ..., $M_\infty$, in which $M_1$ is for Motzkin permutations (that are Motzkin trees up to a bijection) and $M_\infty$ is for Catalan permutations (that are binary trees up to a bijection) an interpolation of our method seems doable.

Deriving linear algorithms for the generation of random objects applies to other structures. Indeed the production of a holonomic function can be mechanized by the software GFUN [26]. However, for my purpose, this method has limitations, since the size of the equation can blow up as shown by Figure 3 of [8] for instance.
# Acknowledgments

I thank Laurent Alonso, Maciej Bendkowski, Alain Giorgetti and Jean-Luc Rémy for interesting discussions and suggestions.

# References

[1] Laurent Alonso. Uniform generation of a Motzkin word. *Theor. Comput. Sci.*, 134(2):529–536, 1994.

[2] Laurent Alonso, Jean-Luc Rémy, and René Schott. Uniform generation of a Schröder tree. *Inf. Process. Lett.*, 64(6):305–308, 1997. doi:10.1016/S0020-0190(97)00174-9.

[3] Didier Arquès and Alain Giorgetti. Une bijection géométrique entre une famille d'hypercartes et une famille de polygones énumérées par la série de Schröder. *Discret. Math.*, 217(1-3):17–32, 2000. doi:10.1016/S0012-365X(99)00253-8.

[4] Axel Bacher, Olivier Bodini, and Alice Jacquot. Efficient random sampling of binary and unary-binary trees via holonomic equations. *CoRR*, abs/1401.1140, 2014.

[5] Axel Bacher, Olivier Bodini, and Alice Jacquot. Efficient random sampling of binary and unary-binary trees via holonomic equations. *Theor. Comput. Sci.*, 695:42–53, 2017. doi:10.1016/j.tcs.2017.07.009.

[6] Elena Barcucci, Alberto Del Lungo, Elisa Pergola, and Renzo Pinzani. From Motzkin to Catalan permutations. *Discret. Math.*, 217(1-3):33–49, 2000. doi:10.1016/S0012-365X(99)00254-X.

[7] Maciej Bendkowski, Olivier Bodini, and Sergey Dovgal. Polynomial tuning of multiparametric combinatorial samplers. In Markus E. Nebel and Stephan G. Wagner, editors, *Proceedings of the Fifteenth Workshop on Analytic Algorithmics and Combinatorics, ANALCO 2018, New Orleans, LA, USA, January 8-9, 2018*, pages 92–106. SIAM, 2018. doi:10.1137/1.9781611975062.9.

[8] Maciej Bendkowski and Pierre Lescanne. On the enumeration of closures and environments with an application to random generation. *Log. Methods Comput. Sci.*, 15(4), 2019. doi:10.23638/LMCS-15(4:3)2019.

[9] Srecko Brlek, Elisa Pergola, and Olivier Roques. Non uniform random generation of generalized Motzkin paths. *Acta Informatica*, 42(8-9):603–616, 2006. doi:10.1007/s00236-006-0008-x.

[10] Alain Denise and Paul Zimmermann. Uniform random generation of decomposable structures using floating-point arithmetic. *Theor. Comput. Sci.*, 218(2):233–248, 1999. doi:10.1016/S0304-3975(98)00323-5.

[11] Luc Devroye. Simulating size-constrained Galton-Watson trees. *SIAM J. Comput.*, 41(1):1–11, 2012. doi:10.1137/090766632.
[12] R. Donaghey and L. W. Shapiro. Motzkin numbers. *Journal of Combinatorial Theory, Series A*, 23(3):291–301, 1977.

[13] Serge Dulucq and Jean-Guy Penaud. Interprétation bijective d'une récurrence des nombres de Motzkin. *Discret. Math.*, 256(3):671–676, 2002. doi:10.1016/S0012-365X(02)00342-4.

[14] I. M. H. Etherington. Some problems of non-associative combinations. *The Edinburgh math. notes*, 32:i-vi, 1940. doi:https://doi.org/10.1017/S0950184300002639.

[15] Philippe Flajolet and Robert Sedgewick. *Analytic Combinatorics*. Cambridge University Press, 2008.

[16] Philippe Flajolet, Paul Zimmermann, and Bernard Van Cutsem. A calculus for the random generation of labelled combinatorial structures. *Theor. Comput. Sci.*, 132(2):1–35, 1994. doi:10.1016/0304-3975(94)90226-7.

[17] Dominique Foata and Doron Zeilberger. A classic proof of a recurrence for a very classical sequence. *J. Comb. Theory, Ser. A*, 80(2):380–384, 1997. doi:10.1006/jcta.1997.2814.

[18] Dominique Gouyou-Beauchamps and Cyril Nicaud. Random generation using binomial approximations. In *21st International Meeting on Probabilistic, Combinatorial, and Asymptotic Methods in the Analysis of Algorithms (AofA'10)*, pages 359–372, Vienna, Austria, 2010.

[19] OEIS Foundation Inc. The on-line encyclopedia of integer sequences. Published electronically at https://oeis.org/, 2022.

[20] Alice Jacquot. *Construction par greffes, combinatoire analytique et génération aléatoire*. PhD thesis, U. Paris-Nord – Sorbonne Paris Cité, April 2014.

[21] Donald E. Knuth. *The Art of Computer Programming, Volume 4, Fascicle 4: Generating All Trees, History of Combinatorial Generation*. Addison-Wesley Publishing Company, 2006.

[22] Pierre Lescanne. Boltzmann samplers for random generation of lambda terms. *CoRR*, abs/1404.3875, 2014. URL: http://arxiv.org/abs/1404.3875, arXiv:1404.3875.

[23] Konstantinos Panagiotou, Leon Ramzews, and Benedikt Stufler. Exact-size sampling of enriched trees in linear time, 2021. URL: https://arxiv.org/abs/2110.11472, doi:10.48550/ARXIV.2110.11472.

[24] Jean-Luc Rémy. Un procédé itératif de dénombrement d'arbres binaires et son application à leur génération aléatoire. *RAIRO Theor. Informatics Appl.*, 19(2):179–195, 1985. doi:10.1051/ita/1985190201791.

[25] Olinde Rodrigues. Sur le nombre de moyens d'effectuer un produit de n facteurs. *Journal de mathématiques pures et appliquées*, page 549, 1938.

[26] Bruno Salvy and Paul Zimmermann. GFUN: a Maple package for the manipulation of generating and holonomic functions in one variable. *ACM Trans. Math. Softw.*, 20(2):163–177, 1994. doi:10.1145/178365.178368.