# Lattice Path Combinatorics and Interactions

Cyril Banderier, Christian Krattenthaler, Michael Wallner

► To cite this version:

Cyril Banderier, Christian Krattenthaler, Michael Wallner. Lattice Path Combinatorics and Interactions. Séminaire Lotharingien de Combinatoire., 87B, pp.1-302, 2024. hal-04799893

HAL Id: hal-04799893

https://hal.science/hal-04799893v1

Submitted on 23 Nov 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.

<figure><img src="image_1.png" /></figure>
Lattice Path Combinatorics and Interactions

Editors:

Cyril Banderier, Christian Krattenthaler, Michael Wallner

Séminaire Lotharingien de Combinatoire, volume 87B
This volume edited by Cyril Banderier, Christian Krattenthaler, and Michael Wallner, corresponds to the volume 87B of the journal Séminaire Lotharingien de Combinatoire, first published online in April 2024 at https://www.mat.univie.ac.at/~slc/.

©©⑈⑉: All the contributions to this volume are copyrighted by their authors under the CC-BY-SA license (International 4.0).
# PREFACE TO THE SPECIAL ISSUE FOR THE 2021 LATTICE PATH CONFERENCE

CYRIL BANDERIER, CHRISTIAN KRATTENTHALER, AND MICHAEL WALLNER

## 1. THE LATTICE PATH CONFERENCE

1.1. **History of the conference.** Following an interdisciplinary interest for lattice paths after the publication in 1979 of the two books *Lattice Path Combinatorics with Statistical Applications* by Tadepalli Venkata Narayana and *Lattice Path Counting and Applications* by Sri Gopal Mohanty, Mohanty initiated a series of conferences on this topic.

The first two events were held at McMaster University in 1984 and 1990, and then at the University of Delhi in 1994, University of Vienna in 1998, University of Athens in 2002, East Tennessee State University in 2007, University of Siena in 2010, Cal Poly Pomona University in 2015, Centre International de Rencontres Mathématiques in 2021. A new initiative started by dedicating the fourth conference in Wien to the memory of Germain Kreweras (1918–1998) and Tadepalli Venkata Narayana (1930–1987), both of whom made a significant contribution to the field. In the same spirit, the 2002 conference was dedicated to the memory of István Vincze (1912–1999). The 2015 conference was dedicated to Shreeram Shankar Abhyankar (1930–2013), Philippe Flajolet (1948–2011), and Lajos Takács (1924–2015).

For each conference there is the tradition to have a journal special issue where anyone could submit work related to lattice paths.<fnref n="1" /> We shall come back to the current special issue in Section 2. Before, let us say a few words about the event itself.

1.2. **The 2021 Lattice Path Conference.** The 9th Lattice Path Conference was held in hybrid form (given the circumstances of the global Covid pandemic) during 21–25 June 2021 at the Centre International de Rencontres Mathématiques (CIRM) in Luminy, France. This event attracted 227 registered participants with affiliations from the “four corners of the world” from 34 different countries of which 24 attended in person. Moreover, about 40 more people attended casually online without being registered (we made the event fully open), and some colleagues even organized some watch parties in their university, so the real number of participants in the event will remain unknown.

<figure><img src="image_2.png" /><figcaption>Sri Gopal Mohanty, father of
the Lattice Path Conference.</figcaption></figure>

<footnote n="1">¹See the website of the Lattice Path Conference for links to these special issues and further information.</footnote>
<figure><img src="image_3.png" /><figcaption>As visible via some details in the above photos, the 2021 Lattice Path Conference was a hybrid conference.</figcaption></figure>

The conference consisted of 27 invited talks and 22 accepted posters, leaving also time for scientific discussions. The presentations covered a wide range of topics, including

• algebraic combinatorics (Young tableaux, representation theory, symmetric polynomials, representation theory, ...),

• analytic combinatorics (kernel method, symbolic method, asymptotics, ...),

• enumerative combinatorics (walks in convex and non-convex cones, alternating sign matrices, generating functions, bijections, ...),

• probability theory (Brownian motion, Markov chains, random trees, ...),

• computer algebra (orthogonal polynomials, holonomy theory, ...),

• theoretical computer science (random generations, context-free grammars, ...),

• number theory (integer partitions, continued fractions, ...),

• theoretical physics (Liouville quantum gravity, square ice, ...).

The conference was a great success, with many stimulating talks and discussions. We would like to thank all of the participants, including the invited speakers and the attendees, for their contributions to the conference. We would also like to extend our gratitude to the staff at CIRM, who, like always, was warmly welcoming, even in this difficult period: CIRM reopened for our meeting after many months of closure due to the Covid containment and the interdiction of meetings. Finally, we would like to thank Kilian Raschel who was the main sponsor of this event via his ERC starting grant COMBINEPIC, which allowed us to help CIRM who financially suffered a lot of so many months of closure.

This event successfully brought together leading experts in the field and contributed to the development of new research directions. Several conjectures were even solved during the event! We hope that our special issue dedicated to this conference will serve as a valuable resource for researchers in combinatorics, probability theory, and statistical physics, and that they will inspire further research on lattice paths.

Let us now give the lists of organizers, speakers, and participants of this meeting.

### 1.3. Scientific and organizing committee.

• Cyril Banderier (CNRS, Université Sorbonne Paris Nord)

• Jehanne Dousse (CNRS, Université Lyon 1)

• Enrica Duchi (Université Paris Diderot)

• Christian Krattenthaler (Universität Wien)

• Greta Panova (University of Southern California)

• Kilian Raschel (CNRS, Université de Tours)

• Michael Wallner (TU Wien)
1.4. Invited talks.

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

George Andrews (Pennsylvania State University):
*Schmidt type partitions and partition analysis*

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

Andrei Asinowski (Alpen-Adria-Universität Klagenfurt):
*Vectorial kernel method and lattice paths with patterns*

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

Philippe Biane (CNRS, Université Paris-Est):
*Mating of discrete trees and walks in the quarter-plane*

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

Alin Bostan (INRIA, Saclay):
*How to prove or disprove the algebraicity of a generating function using a computer*

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

Mireille Bousquet-Mélou (CNRS, Université de Bordeaux):
*Invariants for walks avoiding a quadrant*

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

Timothy Budd (Radboud University):
*Winding of simple walks on the square lattice*

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

Philippe Di Francesco (University of Illinois at Urbana-Champaign and IPhT Saclay):
*Triangular ice: combinatorics and limit shapes*

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

Sergi Elizalde (Dartmouth College):
*Counting lattice paths by the number of crossings and major index*

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

Ilse Fischer (University of Vienna):
*The alternating sign matrices/descending plane partitions relation: n+3 pairs of equivalent statistics*
<figure><img src="image_13.png" /></figure>

Ira Gessel (Brandeis University):
*Redundant generating functions in lattice path enumeration*

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

Vadim Gorin (MIT):
*Addition of matrices at high temperature*

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

Tony Guttmann (Melbourne University):
*Extracting asymptotics from series coefficients*

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

Nina Holden (ETH Zürich):
*Random triangulations and bijective paths to Liouville quantum gravity*

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

Mourad E.H. Ismail (University of Central Florida):
*Orthogonal polynomials, moments, and continued fractions*

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

Satya Majumdar (CNRS, Université Paris Sud):
*Nonintersecting Brownian bridges in the flat-to-flat geometry*

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

Olya Mandelshtam (University of Waterloo):
*A Markov chain on tableaux that projects to a multispecies totally asymmetric zero range process*

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

Irène Marcovici (Université de Lorraine):
*Bjections between walks inside a triangular domain and Motzkin paths of bounded amplitude*

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

Stephen Melczer (University of Pennsylvania):
*Lattice walks and analytic combinatorics in several variables*
<figure><img src="image_22.png" /><figcaption>Robin Pemantle (University of Pennsylvania):
*Generating function technologies: applications to lattice paths*</figcaption></figure>

<figure><img src="image_23.png" /><figcaption>Bruno Salvy (INRIA / ENS Lyon):
*Computation of tight enclosures for Laplacian eigenvalues*</figcaption></figure>

<figure><img src="image_24.png" /><figcaption>Michael Singer (North Carolina State University):
*Differentially algebraic generating series for walks in the quarter plane*</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>Perla Sousi (University of Cambridge):
*The uniform spanning tree in 4 dimensions*</figcaption></figure>

<figure><img src="image_26.png" /><figcaption>Andrea Sportiello (CNRS, Université Paris Nord):
*Boltzmann sampling in linear time: irreducible context-free structures*</figcaption></figure>

<figure><img src="image_27.png" /><figcaption>Xavier Viennot (CNRS, Université de Bordeaux):
*Heaps and lattice paths*</figcaption></figure>

<figure><img src="image_28.png" /><figcaption>Karen Yeats (University of Waterloo):
*Łukasiewicz walks and generalized tandem walks*</figcaption></figure>

<figure><img src="image_29.png" /><figcaption>Doron Zeilberger (Rutgers University):
*Using symbolic dynamical programming in lattice paths combinatorics*</figcaption></figure>

<figure><img src="image_30.png" /><figcaption>Paul Zinn-Justin (Melbourne University):
*Generalized pipe dreams and lower-upper scheme*</figcaption></figure>

A recording of all the talks can be found at https://lipn.fr/~banderier/LPC/2021/.
1.5. Accepted posters.

* Ault Shaun, Charles Kicey: *From lattice paths to standard Young tableaux*

* Cyril Banderier, Marie-Louise Lackner, Michael Wallner: *Latticepathology and symmetric functions*

* Nicholas Beaton: *Walks obeying two-step rules on the square lattice*

* Swee Hong Chan, Igor Pak, Greta Panova: *Log-concavity in posets and random walks*

* Sergey Dovgal, Mohamed Lamine Lamali, Philippe Duchon: *A phase transition in non-deterministic walks with two or more variables*

* Andrew Elvey Price: *Enumeration of walks with small steps by winding angle*

* Rigoberto Flórez, José L. Ramírez, Fabio A. Velandia, Diego Villamiza: *Restricted Dyck paths*

* Xi Chen, Bishal Deb, Alexander Dyachenko, Tomack Gilmore, Alan Sokal: *Coefficientwise total positivity of some matrices defined by linear recurrences*

* Hans Höngesberg: *Weight-preserving bijections between integer partitions and a class of alternating sign trapezoids*

* Heba Ayeda, David Beecher, Alan Krinik, Jeremy J. Lin, David Perez, Thuy Vu Dieu Lu, Weizhong Wong: *Lattice paths with alternating probabilities*

* Josef Küstner, Michael Schlosser, Meesue Yoo: *Lattice paths and negatively indexed weight-dependent binomial coefficients*

* Florian Lehner, Christian Lindorfer, Wolfgang Woess: *The language of self-avoiding walks*

* Satya Majumdar, Francesco Mori, Gregory Schehr: *Distribution of the time between maximum and minimum of random walks*

* Stéphane Ouvry, Alexios Polychronakos, Shuang Wu: *Algebraic area counting for lattice closed random walks*

* Alan Krinik, Gerardo Rubino: *The exponential-dual matrix method: applications to Markov chain analysis*

* Andrei Asinowski, Benjamin Hackl, Sarah Selkirk: *Down-step statistics in generalized Dyck paths*

* Myrto Kallipoliti, Robin Sulzgruber, Elini Tzanaki: *Patterns in Shi tableaux and Dyck paths*

* Malvina Vamvakari: *On q-order statistics*

* Florian Aigner, Gabriel Frieden: *qRSt: A probabilistic Robinson-Schensted correspondence for Macdonald polynomials*

* Quang-Nhat Le, Sinai Robins, Christophe Vignat, Tanay Wakhare: *A continuous analogue of lattice path enumeration*

* Jisun Huh, Sun-Young Nam, Meesue Yoo: *LLT polynomials in a nutshell: on Schur expansion of LLT polynomials*

* Benjamin De Bruyne, Satya Majumdar, Gregory Schehr: *Generating discrete-time constrained random walks*.

The poster session was organized on the online platform Gather.town, where all the posters can still be perused.
1.6. **The Hotel Latticepathologia escape game.** The image below shows a part of the “Hotel Latticepathologia” (a virtual location designed by Cyril Banderier, Jehanne Dousse, and Michael Wallner), where the poster session was taking place.

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

This virtual place gave the online participants the opportunity to socialize and to partic-
ipate in an escape game! In its first part, you have to identify a mysterious mathematician
in each room in the Hotel Latticepathologia (https://tinyurl.com/s6njbsku), using clues
scattered throughout. Will you succeed?

In its second part, you have to solve a collection of puzzles designed in collaboration
with Vivien Ripoll; see https://lipn.fr/~cb/LPC/2021/Puzzles/. It includes a musical
concert by 9 conference participants. Don't miss it!

Below is a sample puzzle from the escape game. Can you solve it?

Confusing Dream

I had the strangest dream last night: a Zoom meeting with some of our mathematical heroes. Luckily I was able to take a screenshot before waking up! They look utterly **confused** by the situation, though I note they all have one extraordinary feature.

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

They were discussing the discovery of a new lattice path. I managed to copy it there:

https://www.youtube.com/watch?v=0000000000
1.7. **Registered participants.** Most of the 227 registered participants were online. A few of them had the opportunity to enjoy CIRM's beautiful surroundings.

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

List of registered participants: David Adame-Carrillo, Mohammed Ageel, Florian Aigner, Marie Albenque, Seamus Albion, Ian Alevy, Irha Ali, George Andrews, Omer Angel, Margaret Archibald, David Ash, Andrei Asinowski, Shaun Ault, Jean-Christophe Aval, Heba Ayeda, Arvind Ayyer, Beáta Bényi, Cyril Banderier, Josaphat Baolahy, Elena Barcucci, Jean-Luc Baril, Erik Bates, Nicholas Beaton, David Beecher, Chiheb Ben Bechir, Sudip Bera, Olivier Bernardi, Antonio Bernini, David Bevan, Philippe Biane, Arthur Blanc-Renaudie, Aubrey Blecher, Alin Bostan, Mireille Bousquet-Mélou, Cédric Boutillier, Jérémie Bouttier, Timothy Budd, Théophile Buffière, Ariane Carrance, Giulio Cerbai, Swee Hong Chan, Linxiao Chen, Shaoshi Chen, Frédéric Chyzak, Lapo Cioni, Alice Contat, Michael Coopman, Sylvie Corteel, Logan Crew, Cesar Cuenca, Nicolas Curien, Stéphane Dartois, Benjamin De Bruyne, Bishal Deb, Nachum Dershowitz, Hiranya Kishore Dey, Philippe Di Francesco, Lucia Di Vizio, Ruiwen Dong, Robert Donley, Jehanne Dousse, Sergey Dovgal, Thomas Dreyfus, Michael Drmota, Enrica Duchi, Dennis Eichhorn, Sergi Elizalde, Andrew Elvey Price, Sen-Peng Eu, Wenjie Fang, Valentin Féray, Luca Ferrari, Ilse Fischer, Rigoberto Flórez, Luis Fredes, Gabriel Frieden, Éric Fusy, Ira Gessel, Sudhir Ghorpade, Juan Gil, Tomack Gilmore, Vadim Gorin, Adam Gregory, Tony Guttmann, Hans Höngesberg, Benjamin Hackl, Aliakbar Haghighi, Eva-Maria Hainzl, Charlotte Hardouin, Kilian Hermann, Clemens Heuberger, Pawel Hitczenko, Hung Hoang, Nina Holden, Sam Hopkins, Yueyun Hu, Justin Hua, Mourad E. H. Ismail, Svante Janson, Helen Jenne, Frédéric Jouhet, Josef Küstner, Wonwoo Kang, Manuel Kauers, Rinat Kedem, Ghizlane Kettani, Mikhail Khristoforov, Charles Kicey, Donghyun Kim, Sergey Kirgizov, Victor Kleptsyn, Arnold Knopfmacher, Isaac Konan, Irina Kourkova, Christian Krattenthaler, Alan Krinik, Nishu Kumari, Raunak Kundagrami, Florian Lehner, Helder Lima, Zhicong Lin, Christian Lindorfer, Martin Loebl, Baptiste Louf, Torsten Mütze, Satya Majumdar, Pritam Majumder, Olya Mandelshtam, Jean-Francois Marckert, Irène Marcovici, Barbara Margolius, Hana Melánová, Stephen Melczer, Laurent Menard, Sri Gopal Mohanty, Derrick Mohlala, Francesco Mori, Lukas Nabergall, Philippe Nadeau, Victor Nador, Mehdi Naima, Hiroshi Naruse, Andreas Nessmann, David Nguyen, Hadrien Notarantonio, Soichi Okada, Stéphane Ouvry, J. E. Paguyo, Nimisha Pahuja, Greta Panova, Jay Pantone, Eveliina Peltola, Robin Pemantle, Karol Penson, Leonid Petrov, Renzo Pinzani, Thomas Prellberg, Helmut Prodinger, Sanjay Ramassamy, José Ramírez, Kilian Raschel, Vivien Ripoll, Tom Roby, Martin Rubey, Gerardo Rubino, Tiadora Ruza, Nasser Saad, Bruce Sagan, Manjil Pratim Saikia, Bruno Salvy, Yoshio Sano, Gilles Schaeffer, Gregory Schehr, Michael Schlosser, Jeanne Scott, Blair Seidler, Sarah Selkirk, Timo Seppalainen, Michael F. Singer, Alexandros Singh, Erik Slivken, Rebecca Smith, Alan Sokal, U-Keun Song, Perla Sousi, Andrea Sportiello, Richard Stanley, Dennis Stanton, Benedikt Stufler, Adrian Tanasa, Benjamin Terlat, Vasu Tewari, Paul Thévenin, Mikhail Tikhonov, Jordan Tirrell, Jessica Tomasko, Joonas Turunen, Eleni Tzanaki, Malvina Vamvakari, Roger Van Peski, Zoé Varin, Ekaterina Vassilieva, Fabio Velandia, Xavier Viennot, Christophe Vignat, Diego Villamizar, Michael Voit, Trung Vu, David Wahiche, Tanay Wakhare, Michael Wallner, Harriet Walsh, Guoliang Wang, Sebastian Wild, Mark Wilson, Peter Winkler, Wolfgang Woess, Elaine Wong, Zaidan Wu, Karen Yeats, Meesue Yoo, Sergey Yurkevich, Doron Zeilberger, Noam Zeilberger, Jiang Zeng, Yan Zhuang, Paul Zinn-Justin.
## 2. OUR SPECIAL ISSUE IN THE “SÉMINAIRE LOTHARINGIEN DE COMBINATOIRE”

The “Séminaire Lotharingien de Combinatoire” is an international biannual seminar, cofounded in 1980 by Dominique Foata (Strasbourg), Adalbert Kerber (Aachen and Bayreuth), and Volker Strehl (Erlangen). The name of the seminar comes from the fact that these cities were almost covered by Lotharingia, a part of the Carolingian Empire. In 1994, an eponymous journal was launched and quickly gathered a wider international audience, welcoming also articles independent of any participation to the actual seminar. On some occasions, the journal has had special issues dedicated to Festschriften or conference proceedings. All volumes are freely accessible at https://www.mat.univie.ac.at/~slc/.

We are pleased to add another volume to this collection of the journal of the Séminaire Lotharingien de Combinatoire. Our volume, dedicated to the themes of the 2021 Lattice Path Conference, contains the following contributions:

1. Helmut Prodinger: *A walk in my lattice path garden*

2. Anthony J. Guttmann and Václav Kotěšovec: *A numerical study of L-convex polyominoes and 201-avoiding ascent sequences*

3. David W. Ash: *Introducing DASEP: the doubly asymmetric simple exclusion process*

4. Stéphane Ouvry and Alexios P. Polychronakos: *Signed area enumeration for lattice paths*

5. Aubrey Blecher and Arnold Knopfmacher: *Left-to-right maxima in Dyck paths*

6. Rigoberto Flórez, Toufik Mansour, José L. Ramírez, Fabio A. Velandia, and Diego Villamizar: *Restricted Dyck paths on valleys sequence*

7. Sergi Elizalde: *Counting lattice paths by crossings and major index II: tracking descents via two-rowed arrays*

8. Jang Soo Kim and Dennis Stanton: *Three families of q-Lommel polynomials*

9. Malvina Vamvakari: *On q-order statistics*

10. Thomas Dreyfus: *Differential algebraic generating series of weighted walks in the quarter plane*

11. Rodolphe Garbit and Kilian Raschel: *The generating function of the survival probabilities in a cone is not rational*

12. Rafik Aguech, Asma Althagafi, and Cyril Banderier: *Height of walks with resets, the Moran model, and the discrete Gumbel distribution*

We would like to thank all the authors for their patience during the interactions with the editors, which (hopefully!) resulted in a pleasantly polished volume... Last but not least, we would like to also thank all the referees for their excellent work. We hope that the reader will enjoy this volume and will eventually be motivated to contribute to the next Lattice Path Conference!
3. A PANORAMA ON LATTICE PATHS

Let us end this preface with two beautiful photos taken during a random walk in the proximity of the conference location (warm thanks to Andreas Nessmann and Sergey Dovgal for sharing them with us).

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

A path... what else could end our preface?
# Prelude

<figure><img src="image_35.png" /><figcaption>A famous lattice path related to Bach...</figcaption></figure>
<figure><img src="image_36.png" /></figure>

Jackson Park, Windsor, Ontario
Wikimedia Commons

# A WALK IN MY LATTICE PATH GARDEN

HELMUT PRODINGER¹

¹STELLENBOSCH UNIVERSITY AND NITHECS (NATIONAL INSTITUTE FOR THEORETICAL AND COMPUTATIONAL SCIENCES), SOUTH AFRICA; https://www.math.tugraz.at/~prodinger

**Abstract.** Various lattice path models are reviewed. The enumeration is done using generating functions. A few bijective considerations are woven in as well. The kernel method is often used. Computer algebra was an essential tool. Some results are new, some have appeared before, but all are interesting.

The lattice path models considered are Hoppy walks and several models involving skew Dyck paths, Schröder paths, hex-trees, decorated ordered trees, multi-edge trees, etc., related to the sequence A002212 in the On-line Encyclopedia of Integer Sequences (created by N. Sloane). Weighted unary-binary trees also occur and we there improve on our old paper on Horton-Strahler numbers [P. Flajolet and H. Prodinger, 1986], by using a different substitution. Some material on Motzkin numbers and paths is also discussed. Some new results on ‘Deutsch paths’ in a strip are included as well. During the Covid period, I spent much time with this beautiful concept that I dare to call Deutsch paths, since Emeric Deutsch stands at the beginning with a problem that he posted in the American Mathematical Monthly some 20 years ago. Peaks and valleys, studied by Rainer Kemp 40 years ago under the names MAX-turns and MIN-turns, are revisited with a more modern approach, streamlining the analysis, relying on the ‘subcritical case’ (named so by Philippe Flajolet), the adding a new slice technique and once again the kernel method.

**Keywords:** Skew Dyck paths, decorated Dyck paths, generating functions, Motzkin paths, kernel method.
CONTENTS

<table><tr><td>1. Introduction</td><td>22</td></tr><tr><td>2. Hoppy walks</td><td>23</td></tr><tr><td>3. Combinatorics of the OEIS sequence A002212</td><td>29</td></tr><tr><td>4. Binary trees and Horton-Strahler numbers</td><td>30</td></tr><tr><td>5. Marked ordered trees</td><td>36</td></tr><tr><td>6. A bijection between multi-edge trees and 3-coloured Motzkin paths</td><td>39</td></tr><tr><td>7. The combinatorics of skew Dyck paths</td><td>42</td></tr><tr><td>8. More about Motzkin paths</td><td>53</td></tr><tr><td>9. Oscillations in Dyck paths revisited</td><td>59</td></tr><tr><td>10. Deutsch-paths in a strip</td><td>63</td></tr><tr><td>11. Conclusion</td><td>67</td></tr><tr><td>References</td><td>68</td></tr></table>

1. INTRODUCTION

Around 20 years ago I published a collection of examples about applications of the kernel method in the present journal. The success of this enterprise was unexpected and came as a very pleasant surprise. My current plan is to present again a collection of subjects, loosely related as they all have a lattice path flavour (trees are also allowed in my private book when lattice paths are mentioned). The subjects cover my last 2 or 3 years of research; some results were only posted on arXiv, and some are completely new.

As in the predecessor paper, the kernel method plays a role again, but also analytic techniques like singularity analysis and Mellin transform, as well as bijective results.

Let me emphasize that this is *not* a survey about lattice path enumeration, but a *personal survey*, that is, a very *personal* account of some of *my* interests and recent activities. Of course, I hope that some people will like it, and that some readers will even contact me to further investigate the models presented in this article. It is *not* a long novel, I organized the material rather as a collection of short stories, roughly in the order as I worked on them.

While some methods might be intimidating for the uninitiated, I tried to provide an accessible introduction, also by explaining and simplifying older proofs and applying old and powerful tools to many beautiful and attractive up-to-date problems. Even with traditional bricks, traditional timber, traditional paint one can build a beautiful house, a beautiful home, a beautiful garden!

People who want to properly learn the subject (via different approaches all keeping some intimate links with enumerative combinatorics) can go to Christian Krattenthaler's survey [29], or the older books by Mohanty [31] and Narayana [32] or to many deep and sophisticated articles by Mireille Bousquet-Mélou and Philippe Flajolet. I apologize to all those that I did not mention although I should have.
The sections are arranged in roughly the form they were conceived. They have their own little introductions so that a casual reader can look at various parts at his/her leisure. Like a walk in a real garden, you can just cross it or you can explore, in butterfly style, various flowers and other beauties.

I spent probably a year of work on that project, and so it can be hoped that the interested enthusiast can learn something from it. In the first place, I would like to mention the large variety of combinatorial objects, perhaps not known to everybody. Then I mention the "right" substitution of auxiliary variables that one finds in various places of the paper. Without them, certain computations/considerations would be (almost) impossible, in particular since computer algebra systems are not good to deal with roots. There are a few asymptotic calculations woven in, both related to my old friends, the height of trees and the register function (Horton-Strahler numbers). As mentioned, I do not distinguish much (cum grano salis) between lattice paths and trees. There are also a few attractive bijections to be found. I spent myriads of hours with Maple to guess the quantities of interest, especially in the context of Deutsch paths in a strip and bounded marked trees. Many explicit formulæ were found with the kernel method, often in combination with the adding-a-new-slice procedure.

My motivation for this work, apart from numerous improvements related to the existing literature, is a personal credo about how to proceed in various instances of this personal survey. I like generating functions and explicit expressions and I am not so thrilled by bounds and estimates. I know that I open a can of worms with such a statement, but nobody is forced to proceed like me. But those who do can take home something beautiful. It is clear to me that many researchers will not agree with 'beautiful' and replace it with something else. Fortunately, different views are allowed, and I mention freely that Philippe Flajolet's work, especially his younger period, had a lasting effect on the way I consider things. I should also mention Emeric Deutsch here, who is responsible for skew Dyck paths, marked trees, Deutsch paths and more.

So I send this *swan song* into the world. May it encourage younger people to pick up some stones in the garden in the hope that they are actually rough diamonds.

## 2. HOPPY WALKS

This section (which is extending our unpublished preprint arXiv:2009.13474) can be seen as a warm-up, introducing all the typical techniques, before a much longer section.

Deng and Mansour [10] introduce a rabbit named Hoppy and let him move according to certain rules. While the story about Hoppy is charming and entertaining, we do not need this here and move straight ahead to the enumeration issues. Eventually, the enumeration problem is one about *k*-Dyck paths ($k \ge 1$). The up-steps are $(1, k)$ and the down-steps are $(1, -1)$. The model that has $(1, 1)$ as up-step and the down-step are $(1, -k)$ will also be called *k*-Dyck paths.

The question is about the length of the last sequence of down-steps (shown in red in Figure 1). Or, phrased differently, how many *k*-Dyck paths end on level *j* with an up-step. Note that such paths have length $n = m + km - j$. The recent paper [49] consider these paths, although without the restriction that the last step must be an up-step.
<figure><img src="image_37.png" /></figure>

<figure><img src="image_38.png" /><figcaption>FIGURE 1. The rabbit Hoppy thinking in the lattice path garden... on the number of final down-steps in paths with steps (1, −1) and (1, k).</figcaption></figure>

Source: https://en.wikipedia.org/wiki/File:Oryctolagus_cuniculus_Tasmania_2.jpg

The original description of Deng and Mansour is a reflection of Figure 1, with up-steps of size 1 and down-steps of size $-k$, but we prefer it as given here, since we are going to use the adding-a-new-slice method; see [18, 38]. A slice is here a run of down-steps, followed by an up-step. So, for each path, one begin with an up-step, and then $m-1$ new slices are added. We keep track of the level after each slice, using a variable $u$. The variable $z$ is used to count the number of up-steps.

Deng and Mansour work out a formula which comprises $O(m)$ terms. For our walks, we obtain a more compact sum of only $O(j)$ terms (recall that $j$ is the level of the last point).

We start with the following substitution which encodes that one adds a new slice

$$u^j \rightarrow z \sum_{0 \le h \le j} u^{h+k} = \frac{zu^k}{1-u}(1-u^{j+1}).$$

Now let $F_m(z, u)$ be the generating function of paths having $m$ runs of down-steps. The substitution leads to

$$F_{m+1}(z, u) = \frac{zu^k}{1-u}F_m(z, 1) - \frac{zu^{k+1}}{1-u}F_m(z, u), \quad F_0(z, u) = zu^k.$$

Let $F = \sum_{m \ge 0} F_m$, so that we do not care about the number $m$ anymore; then

$$F(z, u) = zu^k + \frac{zu^k}{1-u}F(z, 1) - \frac{zu^{k+1}}{1-u}F(z, u),$$

or, more conveniently,

$$F(z, u)(1 - u + zu^{k+1}) = zu^k(1 - u + F(z, 1)).$$

The equation $1 - u + zu^{k+1} = 0$ (the so-called kernel equation) is famous when enumerating $(k+1)$-ary trees (or $k$-Dyck paths). Its relevant combinatorial solution (also the only one being analytic at the origin) is

$$\bar{u} = \sum_{\ell \ge 0} \frac{1}{1 + \ell(k+1)} \binom{1 + \ell(k+1)}{\ell} z^{\ell}.$$
Now, since $u := \bar{u}$ cancels the kernel and the left-hand side, $(u - \bar{u})$ must be a factor of the right-hand side (which is a polynomial in $u$); this gives

$$zu^k(1 - u + F(z, 1)) = -zu^k(u - \bar{u}).$$

Cancelling the kernel equation is thus a method which brings additional equations, allowing us to identify $F(z, 1)$, and then $F(z, u)$ which is given by

$$F(z, u) = zu^k \frac{\bar{u} - u}{1 - u + zu^{k+1}}.$$

The first factor has even a combinatorial interpretation, as a description of the first step of the path. It is also clear from this that the level reached is at least $k$ after each slice. We do not care about the factor $zu^k$ anymore, as it produces only a simple shift. The main interest is now how to get to the coefficients of

$$\frac{\bar{u} - u}{1 - u + zu^{k+1}}$$

in an efficient way. First we deal with the denominators $(j \ge k + 1)$

$$S_j := [u^j] \frac{1}{1 - u + zu^{k+1}} = \sum_{0 \le m \le j/k} (-1)^m \binom{j - km}{m} z^m.$$

One way to see this formula is to prove by induction that the sums $S_j$ satisfy the recursion $S_j - S_{j-1} + zS_{j-k-1} = 0$ and initial conditions $S_0 = \cdots = S_k = 1$. In [49] such expressions also appear as determinants. Summarizing,

$$\frac{1}{1 - u + zu^{k+1}} = \sum_{m \ge 0} (-1)^m z^m \sum_{j \ge km} \binom{j - km}{m} u^j.$$

Now we read off coefficients:

$$[u^j] \frac{\bar{u}}{1 - u + zu^{k+1}} = \sum_{0 \le m \le j/k} (-1)^m \binom{j - km}{m} z^m \sum_{\ell \ge 0} \frac{1}{1 + \ell(k + 1)} \binom{1 + \ell(k + 1)}{\ell} z^\ell$$

and further

$$[z^n][u^j] \frac{\bar{u}}{1 - u + zu^{k+1}} = \sum_{0 \le m \le j/k} \frac{(-1)^m}{1 + (n-m)(k+1)} \binom{j-km}{m} \binom{1 + (n-m)(k+1)}{n-m}.$$

The final answer to the Deng–Mansour enumeration (without the shift) is

$$\left( \sum_{0 \le m \le j/k} \frac{(-1)^m}{1 + (n-m)(k+1)} \binom{j-km}{m} \binom{1 + (n-m)(k+1)}{n-m} \right) - (-1)^n \binom{j-1-kn}{n} \quad (2.1)$$

If one wants to take care of the factor $zu^k$ as well, one needs to do the replacements $n \to n + 1$ and $j \to j + k$ in the formula just derived. That enumerates then the $k$-Dyck paths ending at level $j$ after $n$ up-steps, where the last step is an up-step.

The main contribution of this section is this equation (2.1); let us now discuss three variations about these walks.
**An application.** In [1] the authors considered the total number of down-steps of the last down-run in all k-Dyck paths. For $k = 2, 3, 4$, this corresponds to the sequences A334680, A334682, A334719 in the OEIS<fnref n="1" />, respectively. So, if the path ends on level j, the contribution to the total is j.

All we have to do here is to differentiate

$$F(z, u) = zu^k \frac{\bar{u} - u}{1 - u + zu^{k+1}}$$

with respect to $u$, and then replace $u$ by 1. The result is

$$\frac{\bar{u}}{z} - \bar{u} - \frac{1}{z},$$

and the coefficient of $z^m$ therein is

$$\frac{1}{1 + (m + 1)(k + 1)} \binom{1 + (m + 1)(k + 1)}{m + 1} - \frac{1}{1 + m(k + 1)} \binom{1 + m(k + 1)}{m}.$$

The bivariate generating function does this enumeration cleanly and quickly.

**Hoppy's early adventures.** Now we investigate what Hoppy does after his first up-step; he might follow with 0, 1, ..., k down-steps. Eventually, we want to sum all these steps (red in the picture).

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

A new slice is now an up-step, followed by a sequence of down-steps. The substitution of interest is:

$$u^i \rightarrow z \sum_{0 \le h \le i+k} u^h = \frac{z}{1-u} - \frac{zu^{i+k+1}}{1-u}.$$

Furthermore

$$F_{h+1}(z, u) = \frac{z}{1-u} F_h(z, 1) - \frac{zu^{k+1}}{1-u} F_h(z, u),$$

and $F_0 = u^h$, the starting level. We have

$$H(z, u) = \sum_{h \ge 0} F_h(z, u) = u^h + \frac{z}{1-u} H(z, 1) - \frac{zu^{k+1}}{1-u} H(z, u)$$

or

$$H(z, u)(1 - u + zu^{k+1}) = u^h(1 - u) + zH(z, 1).$$

<footnote n="1">¹The On-line Encyclopedia of Integer Sequences (OEIS) is a database available at https://oeis.org/.</footnote>
Plugging in $\bar{u}$ into the right-hand side gives 0, thus one has

$$zH(z, 1) = -\bar{u}^h(1 - \bar{u}),$$

which itself implies

$$H(z, u) = \frac{u^h(1-u) - \bar{u}^h(1-\bar{u})}{1 - u + zu^{k+1}}.$$

But we only need $H(z, 0)$, since we return to the x-axis at the end:

$$H(z, 0) = [h = 0] + \bar{u}^{h+1} - \bar{u}^h.$$

(The Iverson notation $\llbracket P \rrbracket$ which is 1 when $P$ is true and 0 otherwise is the notation of choice for combinatorialists, [21].) The total contribution of red steps is then

$$k + \sum_{h=0}^{k} (k-h)(\bar{u}^{h+1} - \bar{u}^h) = \sum_{h=1}^{k} \bar{u}^h;$$

the coefficient of $z^m$ in this is the total contribution. Since $\bar{u} = 1 + z\bar{u}^{k+1}$, there is the further simplification

$$-1 + \frac{1}{z} + \frac{1}{1-\bar{u}} = \sum_{m \ge 1} \frac{k}{m+1} \binom{(k+1)m}{m} z^m.$$

Indeed, for $m \ge 1$, we have

$$\begin{align*}
[z^m]\left(-1 + \frac{1}{z} + \frac{1}{1-\bar{u}}\right) &= -[z^m]\frac{1}{z\bar{u}^{k+1}} \\
&= -[z^{m+1}]\sum_{\ell \ge 0} \frac{-(k+1)}{(k+1)\ell - (k+1)} \binom{(k+1)\ell - (k+1)}{\ell} z^{\ell} \\
&= [z^{m+1}]\sum_{\ell \ge 0} \frac{(k+1)}{(k+1)(\ell-1)} \binom{(k+1)(\ell-1)}{\ell} z^{\ell} \\
&= \frac{(k+1)}{(k+1)m} \binom{(k+1)m}{m+1} = \frac{k}{m+1} \binom{(k+1)m}{m}.
\end{align*}$$

We did not expect such a simple answer $\frac{k}{m+1} \binom{(k+1)m}{m}$ to this question about Hoppy's early adventures! This analysis of Hoppy's early adventures covers sequences A007226, A007228, A124724 of [52], with references to [1].

**Hoppy walks into negative territory.** Hoppy is now adventurous and allows himself to go to level $-1$ as well, but not deeper. The setup with generating functions is the same, but the $u$-variable counts the level relative to the $-1$ level, so this has to be corrected later.

Hoppy, after some initial frustration discovers that he can now start with an up-step or a down-step. First, let us start Hoppy with an up-step:

$$F(z, u) = zu^{k+1} + \frac{zu^k}{1-u}F(z, 1) - \frac{zu^{k+1}}{1-u}F(z, u),$$

which we conveniently rewrite as

$$F(z, u)(1 - u + zu^{k+1}) = zu^{k+1}(1 - u) + zu^k F(z, 1).$$
Since the left-hand side cancels for $u = \bar{u}$, we get that $u = \bar{u}$ is also cancelling the right-hand side (which is a polynomial in $u$), and this implies that

$$ \bar{u}(1 - \bar{u}) + F(z, 1) = 0. $$

This finally gives

$$ F(z, u) = \frac{zu^k}{1 - u + zu^{k+1}} \left( u(1 - u) - \bar{u}(1 - \bar{u}) \right). $$

But Hoppy can also start with a downstep. So we have to add the result of the previous computation, and get finally

$$ G(z, u) = \frac{zu^k}{1 - u + zu^{k+1}} \left( u(1 - u) - \bar{u}(1 - \bar{u}) \right) + \frac{zu^k}{1 - u + zu^{k+1}} \left( \bar{u} - u \right) $$

or better

$$ G(z, u) = \frac{zu^k}{1 - u + zu^{k+1}} (\bar{u}^2 - u^2). $$

Now we need

$$ \frac{\partial}{\partial u} G(z, 1) - G(z, 1). $$

This subtraction is necessary, since the contribution of $u^j$ is not $j$ as before but only $j-1$. The result is

$$ \frac{\bar{u}^2}{z} - 2\bar{u}^2 - \frac{1}{z}. $$

Hoppy knows that $\bar{u}^d$ has beautiful coefficients:

$$ \bar{u}^d = \sum_{\ell \ge 0} \binom{d-1+(k+1)\ell}{\ell} \frac{d}{k\ell+d} $$

and he inserts $k=2$ which gives A030983:

$$ 3z + 16z^2 + 83z^3 + 442z^4 + 2420z^5 + \dots, $$

$k=3$ which gives A334608:

$$ 5z + 34z^2 + 236z^3 + 1714z^4 + 12922z^5 + \dots, $$

$k=4$ which gives A334610:

$$ 7z + 58z^2 + 505z^3 + 4650z^4 + 44677z^5 + \dots. $$

In general, we have

$$ \frac{\bar{u}^2}{z} - 2\bar{u}^2 - \frac{1}{z} = \sum_{\ell \ge 0} \left[ \binom{1+(k+1)(\ell+1)}{\ell+1} \frac{2}{k(\ell+1)+2} - 2 \binom{1+(k+1)\ell}{\ell} \frac{2}{k\ell+2} \right] z^{\ell}. $$

Happy Hoppy decides to stop this line of computations here.
### 3. COMBINATORICS OF THE OEIS SEQUENCE A002212

The following (sub)sections give some (mostly new) results about the sequence

1, 1, 3, 10, 36, 137, 543, 2219, 9285, 39587, 171369, 751236, 3328218, 14878455, ...,

which is A002212 in the OEIS. In the following four sections we consider different combinatorial structures enumerated by this sequence.

* Hex-trees are identified as weighted unary-binary trees, with weight one (see the article by Hana Kim and Richard Stanley [26]). Apart from left and right branches, as in binary trees, there are also unary branches, and they can come in different colours, here in just one colour. Unary-binary trees played a role in the present authors scientific development, as documented in [17], a paper written with the late and great Philippe Flajolet, about the register function (Horton-Strahler numbers) of unary-binary trees. In Section 4, we offer an improvement, using a “better” substitution than in [17]. The results can now be made fully explicit. As a by-product, this provides a definition and analysis of the Horton-Strahler numbers of hex-trees.

* Then we move to skew Dyck paths, as considered by Emeric Deutsch, Emanuele Munarini, and Simone Rinaldi in [12]. They are like Dyck paths, but allow for an extra step (−1, −1), provided that the path does not intersect itself. An equivalent model, defined and described using a bijection, is from [12]: marked ordered trees; see Section 5. They are like ordered trees, with an additional feature, namely each rightmost edge (except for one that leads to a leaf) can be coloured with two colours. Since we find this class of trees to be interesting, we analyze two parameters of them: number of leaves and height. While the number of leaves for ordered trees is about $n/2$, it is only $n/10$ in the new model. For the height, the leading term $\sqrt{\pi n}$ drops to $\frac{2}{\sqrt{5}}\sqrt{\pi n}$. Of course, many more parameters of this new class of trees could be investigated, which we encourage to do.

* The next two classes of structures are multi-edge trees; see Section 6. Our interest in them was already triggered in an earlier publication, together with Clemens Heuberger and Stephan Wagner [24]. They may be seen as ordered trees, but with weighted edges. The weights are integers $\ge 1$, and a weight $a$ may be interpreted as $a$ parallel edges. The other class are 3-Motzkin paths. They are like Motzkin paths (Dyck paths plus horizontal steps); however, the horizontal steps come in three different colours. A bijection is described. Since 3-Motzkin paths and multi-edge trees are very much alike (using a variation of the classical rotation correspondence), all the structures that are discussed in this paper can be linked via bijections. Since these trees are not so common in the combinatorics community, details and examples are presented for the readers' benefit.

* Skew Dyck paths are finally discussed in more detail in Section 7.
Amazingly, the recursion for the generating functions $R_p(z)$ can be solved explicitly!
The substitution

$$z = \frac{u}{(1+u)^2}$$

that de Bruijn, Knuth, and Rice [9] also used, produces the nice expression

$$R_p(z) = \frac{1-u^2}{u} \frac{u^{2p}}{1-u^{2p+1}}.$$

Of course, once this is known, it can be proved by induction, using the recursive formula.
For the readers benefit, this will be sketched now. We start with the auxiliary formula

$$\sum_{0 \le j < p} \frac{u^{2j}}{1 - u^{2j+1}} = \frac{u}{1-u} - \frac{u^{2p}}{1 - u^{2p}},$$

which we will prove now by induction: For $p = 0$, the formula $0 = \frac{u}{1-u} - \frac{u}{1-u}$ is correct, and then

$$\begin{aligned} \sum_{0 \le j < p+1} \frac{u^{2j}}{1 - u^{2j+1}} &= \frac{u}{1-u} - \frac{u^{2p}}{1 - u^{2p}} + \frac{u^{2p}}{1 - u^{2p+1}} \\ &= \frac{u}{1-u} - \frac{u^{2p}(1+u^{2p})}{1 - u^{2p+1}} + \frac{u^{2p}}{1 - u^{2p+1}} = \frac{u}{1-u} - \frac{u^{2p+1}}{1 - u^{2p+1}}. \end{aligned}$$

Now the formula for $R_p(z)$ can also be proved by induction. First, $R_0(z) = \frac{1-u^2}{u} \frac{u}{1-u^2} = 1$, as it should, and

$$\begin{aligned} & z R_{p-1}^2(z) + 2z R_p(z) \sum_{j<p} R_j(z) \\ &= \frac{u}{(1+u)^2} \frac{(1-u^2)^2}{u^2} \frac{u^{2p}}{(1-u^{2p})^2} + \frac{2u}{(1+u)^2} R_p(z) \sum_{j<p} \frac{1-u^2}{u} \frac{u^{2j}}{1-u^{2j+1}} \\ &= \frac{u^{2p-1}(1-u)^2}{(1-u^{2p})^2} + \frac{2(1-u)}{(1+u)} R_p(z) \sum_{j<p} \frac{u^{2j}}{1-u^{2j+1}}. \end{aligned}$$

Solving

$$R_p(z) = \frac{u^{2p-1}(1-u)^2}{(1-u^{2p})^2} + \frac{2(1-u)}{(1+u)} R_p(z) \left[ \frac{u}{1-u} - \frac{u^{2p}}{1-u^{2p}} \right]$$

leads to

$$R_p(z) \frac{1-u}{1+u} \left[ 1 + 2 \frac{u^{2p}}{1-u^{2p}} \right] = \frac{u^{2p-1}(1-u)^2}{(1-u^{2p})^2},$$

or, simplified

$$R_p(z) = \frac{u^{2p-1}(1-u^2)}{(1-u^{2p})(1+u^{2p})} = \frac{1-u^2}{u} \frac{u^{2p}}{1-u^{2p+1}},$$

which is the formula that we needed to prove. Alternatively, this formula can also be proved by converting the sum into a telescoping sum, by extending the numerator and denominator by $(1 - u^{2j})$, and using a partial fraction decomposition. □
**Weighted unary-binary trees and Horton–Strahler numbers.** The family of unary-binary trees $\mathcal{M}$ might be defined by the symbolic equation

$$ \mathcal{M} = \Box + \mathcal{O} + \mathcal{O} + \mathcal{M} \setminus \{\Box\} + \mathcal{M} \setminus \mathcal{M} $$

The equation for the generating function is

$$ M(z) = 1 + z(M(z) - 1) + zM(z)^2 $$

with the solution

$$ M(z) = \frac{1 - z - \sqrt{1 - 6z + 5z^2}}{2z} = 1 + z + 3z^2 + 10z^3 + 36z^4 + \dots; $$

the coefficients form again sequence A002212 in [52] and enumerate Schröder paths, among many other things. We will come to equivalent structures a bit later.

In the instance of unary-binary trees, we can also work with a substitution. Set $z = \frac{u}{1+3u+u^2}$, then $M(z) = 1 + u$. Unary-binary trees and the register function were investigated in [17], but the present favourable substitution was not used. Therefore, in this previous paper, asymptotic results were available but no explicit formulæ.

This works also with a weighted version, where we allow unary edges with a different colours. Then

$$ \mathcal{N} = \Box + a \cdot \mathcal{O} + \mathcal{O} + \mathcal{N} \setminus \{\Box\} + \mathcal{N} \setminus \mathcal{N} $$

and with the substitution $z = \frac{u}{1+(a+2)u+u^2}$, the generating function is beautifully expressed as $N(z) = 1 + u$. For $a = 0$, this covers also binary trees.

We will consider the Horton–Strahler numbers of unary-binary trees in the sequel. The definition is naturally extended by

$$ \mathrm{reg}\left(\mathcal{O} \middle| t\right) = \mathrm{reg}(t). $$

Now we can move again to $R_p(z)$, the generating function of (generalized) unary-binary trees with Horton–Strahler number $p$. The recursion (for $p \ge 1$) is

$$ \mathcal{R}_p = \mathcal{R}_{p-1} \mathcal{R}_{p-1} + \mathcal{R}_p \sum_{j<p} \mathcal{R}_j + \sum_{j<p} \mathcal{R}_j \mathcal{R}_p + a \cdot \mathcal{R}_p $$

In terms of generating functions, these equations read as

$$ R_p(z) = zR_{p-1}^2(z) + 2zR_p(z) \sum_{j<p} R_j(z) + azR_p(z), \quad p \ge 1; \quad R_0(z) = 1. $$

Amazingly, with the substitution $z = \frac{u}{1+(a+2)u+u^2}$, formally we get the same solution as in the binary case:

$$ R_p(z) = \frac{1 - u^2}{u} \frac{u^{2p}}{1 - u^{2p+1}}. $$
Then we introduce $S_p(z) = R_p(z)+R_{p+1}(z)+R_{p+2}(z)+\cdots$, the generating function of trees with Horton-Strahler number $\ge p$. Using the summation formula proved earlier, we get

$$S_p(z) = \frac{1-u^2}{u} \frac{u^{2p}}{1-u^{2p}} = \frac{1-u^2}{u} \sum_{k \ge 1} u^{k2p}.$$

**Asymptotics.** We start by deriving asymptotics for the number of (generalized) unary-binary trees with $n$ (internal) nodes. This is a standard application of singularity analysis of generating functions, as described in [16] and [20].

We start from the generating function

$$N(z) = \frac{1 - az - \sqrt{1 - 2(a + 2)z + a(a + 4)z^2}}{2z}$$

and determine the singularity closest to the origin, which is the value making the square root disappear: $z = \frac{1}{a+4}$. After that, the local expansion of $N(z)$ around this singularity is determined:

$$N(z) \sim 2 - \sqrt{a+4}\sqrt{1-(a+4)z}.$$

The translation lemmas given in [16] and [20] provide the asymptotics:

$$[z^n]N(z) \sim [z^n]\left(2 - \sqrt{a+4}\sqrt{1-(a+4)z}\right) \\ = -\sqrt{a+4}(a+4)^n \frac{n^{-3/2}}{\Gamma(-\frac{1}{2})} = (a+4)^{n+1/2} \frac{1}{2\sqrt{\pi n^{3/2}}}.$$

Just note that $a=0$ is the well-known formula for binary trees with $n$ nodes.

Now we move to the generating function for the average number of registers. Apart from normalization it is

$$\sum_{p \ge 1} p R_p(z) = \sum_{p \ge 1} S_p(z) = \frac{1-u^2}{u} \sum_{p \ge 1} \sum_{k \ge 1} u^{k2p} = \frac{1-u^2}{u} \sum_{n \ge 1} v_2(n)u^n,$$

where $v_2(n)$ is the highest exponent $k$ such $2^k$ divides $n$.

This has to be studied around $u=1$, which, upon setting $u = e^{-t}$, means around $t=0$. Eventually, and that is the only thing that is different here, this is to be retranslated into a singular expansion of $z$ around its singularity, which depends on the parameter $a$.

For the reader's convenience, we also repeat the steps that were known before. The first factor is elementary:

$$\frac{1-u^2}{u} \sim 2t + \frac{1}{3}t^3 + \cdots.$$

For

$$\sum_{p \ge 1} \sum_{k \ge 1} e^{-k2pt},$$

one applies the Mellin transform, with the result

$$\frac{\Gamma(s)\zeta(s)}{2^s - 1}.$$

Applying the inversion formula, one finds

$$\sum_{p \ge 1} \sum_{k \ge 1} e^{-k2pt} = \frac{1}{2\pi i} \int_{2-i\infty}^{2+i\infty} t^{-s} \frac{\Gamma(s)\zeta(s)}{2^s - 1} ds.$$
Shifting the line of integration to the left, the residues at the poles $s = 1, s = 0, s = \chi_k = \frac{2k\pi i}{\log 2}, k \neq 0$ provide enough terms for our asymptotic expansion.

$$ \frac{1}{t} + \frac{\gamma}{2\log 2} - \frac{1}{4} - \frac{\log \pi}{2\log 2} + \frac{\log t}{2\log 2} + \frac{1}{\log 2} \sum_{k \neq 0} \Gamma(\chi_k)\zeta(\chi_k)t^{-\chi_k}. $$

Combined with the elementary factor, this leads to

$$ 2 + \left( \frac{\gamma}{\log 2} - \frac{1}{2} - \frac{\log \pi}{\log 2} + \frac{\log t}{\log 2} \right) t + \frac{2t}{\log 2} \sum_{k \neq 0} \Gamma(\chi_k)\zeta(\chi_k)t^{-\chi_k} + O(t^2 \log t). $$

Now we want to translate into the original z-world. Since $z = \frac{u}{1+(a+2)u+u^2}$, $u = 1$ translates into the singularity $z = \frac{1}{a+4}$. Further,

$$ t \sim \sqrt{a+4} \cdot \sqrt{1-z(a+4)}, $$

let us abbreviate $A = a + 4$, and we now want to get the asymptotic behaviour of the coefficients in the power series expansion of

$$ \frac{\sqrt{A} \cdot \sqrt{1 - zA \log(1 - zA)}}{2 \log 2} + \left( \frac{\gamma}{\log 2} - \frac{1}{2} - \frac{\log \pi}{\log 2} + \frac{\log A}{2 \log 2} \right) \sqrt{A} \cdot \sqrt{1 - zA} + \frac{2}{\log 2} \sum_{k \neq 0} \Gamma(\chi_k) \zeta(\chi_k) A^{\frac{1-\chi_k}{2}} (1-zA)^{\frac{1-\chi_k}{2}}. $$

By singularity analysis (see [16, 20]), one has

$$ [z^n](1-z)^\alpha \sim \frac{n^{-\alpha-1}}{\Gamma(-\alpha)} $$

and

$$ [z^n] \log(1-z) \sqrt{1-z} \sim \frac{n^{-3/2} \log n}{2\sqrt{\pi}} + \frac{n^{-3/2}}{2\sqrt{\pi}} (-2 + \gamma + 2 \log 2). $$

We start with the most complicated term:

$$ \frac{[z^n] \frac{\sqrt{A} \cdot \sqrt{1-zA \log(1-zA)}}{2 \log 2}}{[z^n] N(z)} \sim \frac{\sqrt{A} A^n \left( \frac{n^{-3/2} \log n}{2\sqrt{\pi}} + \frac{n^{-3/2}}{2\sqrt{\pi}}(-2+\gamma+2\log 2) \right)}{2 \log 2 A^{n+1/2} \frac{1}{2\sqrt{\pi n^{3/2}}}} \\ = \log_4 n + 1 + \frac{\gamma}{2 \log 2} - \frac{1}{\log 2}. $$

The next term we consider is

$$ \left( \frac{\gamma}{\log 2} - \frac{1}{2} - \frac{\log \pi}{\log 2} + \frac{\log A}{2 \log 2} \right) \sqrt{A} \frac{[z^n] \sqrt{1-zA}}{[z^n] N(z)} \\ \sim \left( \frac{\gamma}{\log 2} - \frac{1}{2} - \frac{\log \pi}{\log 2} + \frac{\log A}{2 \log 2} \right) \sqrt{A} \frac{[z^n] \sqrt{1-zA}}{-\sqrt{A}[z^n] \sqrt{1-zA}} \\ = -\frac{\gamma}{\log 2} + \frac{1}{2} + \frac{\log \pi}{\log 2} - \frac{\log A}{2 \log 2}. $$
The last term we consider is

$$ \frac{2}{\log 2} \Gamma(\chi_k) \zeta(\chi_k) A^{\frac{1-\chi_k}{2}} \frac{[z^n](1-zA)^{\frac{1-\chi_k}{2}}}{-\sqrt{A}[z^n]\sqrt{1-zA}} \sim -\frac{4\sqrt{\pi}}{\log 2} \frac{\Gamma(\chi_k)\zeta(\chi_k)}{\Gamma(\frac{\chi_k-1}{2})} A^{\frac{1-\chi_k}{2}} n^{\chi_k/2}. $$

Eventually we have evaluated the average value of the Horton–Strahler numbers:

**Theorem 4.1.** The average Horton–Strahler number of weighted unary-binary trees with $n$ nodes is given by the asymptotic formula

$$ \log_4 n - \frac{\gamma}{2 \log 2} - \frac{1}{\log 2} + \frac{3}{2} + \frac{\log \pi}{\log 2} - \frac{\log A}{2 \log 2} - \frac{4\sqrt{\pi A}}{\log 2} \sum_{k \neq 0} \frac{\Gamma(\chi_k)\zeta(\chi_k)}{\Gamma(\frac{\chi_k-1}{2})} A^{-\frac{\chi_k}{2}} n^{\chi_k/2} \\ = \log_4 n - \frac{\gamma}{2 \log 2} - \frac{1}{\log 2} + \frac{3}{2} + \frac{\log \pi}{\log 2} - \frac{\log A}{2 \log 2} + \psi(\log_4 n), $$

with a tiny periodic function $\psi(x)$ of period 1.

These oscillations are usually bounded by $10^{-5}$, say. See [19] for some explicit error bounds.

## 5. MARKED ORDERED TREES

In [12] we find the following variation of ordered trees. Each rightmost edge might be marked or not, if it does not lead to an endnode (leaf). We depict a marked edge by the red colour and draw all of them of size 4 (4 nodes):

<figure><img src="image_40.png" /><figcaption>FIGURE 2. All 10 marked ordered trees with 4 nodes.</figcaption></figure>

Accordingly, the marked ordered trees satisfy the following symbolic equation (where $\mathcal{A} \cdots \mathcal{A}$ refers to $\ge 0$ copies of $\mathcal{A}$):

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

In terms of generating functions, this gives the functional equation

$$ A(z) = z + \frac{z}{1-A(z)}z + \frac{z}{1-A(z)}2(A(z)-z), $$

whose solution is

$$ A(z) = \frac{1 - z - \sqrt{1 - 6z + 5z^2}}{2} = z + z^2 + z^3 + 3z^3 + 10z^4 + 36z^5 + \dots $$
In fact, as proved in [12], these trees are in bijection with an interesting family of lattice paths, called skew Dyck paths. The bijection performs a walk around the contour of the tree (that is, a depth-first search traversal) and translates it into a skew Dyck path as follows

• black or red edges on the way down become a (+1, +1) step,

• black edges on the way up become a (+1, -1) step,

• red edges on the way up become a (-1, -1) step.

Thus, the 10 trees of Figure 2 translate as follows into skew Dyck paths of length 6:

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

We will analyze several parameters of skew Dyck paths in Section 7. But, as the present author believes that trees are more manageable (than these paths) when it comes to enumeration issues, let us now investigate these marked trees in more detail.

**Parameters of marked ordered trees.** There are many parameters, usually considered in the context of ordered trees, that can be considered for marked ordered trees. Of course, we cannot be encyclopedic about such parameters. We just consider a few parameters and leave further analysis to the future.

*The number of leaves.* To get this, it is most natural to use an additional variable $u$ when translating the symbolic equation, so that $z^nu^k$ refers to trees with $n$ nodes and $k$ leaves. One obtains

$$F(z, u) = zu + \frac{z}{1 - F(z, u)} (zu + 2(F(z, u) - zu)),$$

with the solution

$$F(z, u) = -z + \frac{zu}{2} + \frac{1}{2} - \frac{1}{2}\sqrt{4z^2 - 4z + z^2u^2 - 2zu + 1} \\ = zu + z^2u + (2u + u^2)z^3 + (4u + 5u^2 + u^3)z^4 + \cdots$$

The factor $4u + 5u^2 + u^3$ corresponds to distribution of leaves in the 10 trees of Figure 2.

Of interest is also the average number of leaves, when all marked ordered trees of size $n$ are considered to be equally likely. For that, we differentiate $F(z, u)$ with respect to $u$, and set $u := 1$, with the result

$$\frac{z}{2} + \frac{z - z^2}{2\sqrt{1 - 6z + 5z^2}} = \frac{z}{1 - v}, \quad \text{with the parametrization} \quad z = \frac{v}{1 + 3v + v^2}. \quad (5.1)$$

Since $F(z, 1) = z(1 + v)$, it follows that the average is asymptotic to

$$\frac{[z^{n+1}]_{\frac{z}{1-v}}}{[z^{n+1}]z(1+v)} = \frac{[z^n]_{\frac{1}{1-v}}}{[z^n](1+v)} = \frac{[z^n]_{\frac{1}{\sqrt{5}}} \frac{1}{\sqrt{1-5z}}}{5^{n+\frac{1}{2}} \frac{1}{2\sqrt{\pi}} n^{3/2}} = \frac{\frac{n^{-1/2}}{\Gamma(\frac{1}{2})}}{5^{n+\frac{1}{2}} \frac{1}{2\sqrt{\pi}} n^{3/2}} = \frac{n}{10}. \quad (5.2)$$

Note that the corresponding number for ordered trees (unmarked) is $\frac{n}{2}$, so we have significantly less leaves here.
*The height.* As in the seminal paper [9], we define the height in terms of the longest chain of nodes from the root to a leaf. Further, let $p_h = p_h(z)$ be the generating function of marked ordered trees of height at least $h$. From the symbolic equation, one has

$$p_{h+1} = z + \frac{z^2}{1-p_h} + \frac{2z(p_h-z)}{1-p_h} = -z + \frac{2z-z^2}{1-p_h} \quad (\text{for } h \ge 1) \text{ and } p_1 = z.$$

By some creative guessing, separating numerator and denominator, we find the solution (where we use the auxiliary algebraic function $v$, implicitly defined in (5.1)):

$$p_h = z(1+v) \frac{(1+2v)^{h-1} - v^h(v+2)^{h-1}}{(1+2v)^{h-1} - v^{h+1}(v+2)^{h-1}}.$$

This formula is in fact proved by induction (we start with $p_1 = z(1+v)^{\frac{1-v}{1-v^2}} = z$ and, then, the induction step is best checked using a computer).

The limit of $p_h$ for $h \to \infty$ is $z(1+v)$, the generating function of all marked ordered trees, as expected. Taking differences, we get the generating functions of trees of height at least $h$:

$$p_{\infty} - p_h = z(1-v^2) \frac{(v+2)^{h-1}v^h}{(1+2v)^{h-1} - v^{h+1}(v+2)^{h-1}}.$$

From this, the average height can be worked out, as in the model paper [24]. We sketch the essential steps. For the average height, one needs

$$\sum_{h \ge 0} z(1-v^2) \frac{(v+2)^{h-1}v^h}{(1+2v)^{h-1} - v^{h+1}(v+2)^{h-1}}$$

and its behaviour around $v = 1$, viz.

$$2z(1-v) \sum_{h \ge 0} \frac{3^{h-1}v^h}{3^{h-1} - v^{h+1}3^{h-1}} \sim 2z(1-v) \sum_{h \ge 1} \frac{v^h}{1-v^h}.$$

The behaviour of the series can be taken straight from [24]. We find there

$$\sum_{h \ge 1} \frac{v^h}{1 - v^h} = - \frac{\log(1 - v)}{1 - v}$$

and

$$\sum_{h \ge 0} z(1 - v^2) \frac{(v + 2)^{h-1} v^h}{(1 + 2v)^{h-1} - v^{h+1} (v + 2)^{h-1}} \sim -2z \log(1 - v).$$

Thus, the coefficient of $z^{n+1}$ is asymptotic to $-2[z^n] \log(1-v)$. Since $1-v \sim \sqrt{5\sqrt{1-5z}}$, $-2z \log(1-v) \sim -2z \log \sqrt{1-5z} = -z \log(1-5z)$,

and the coefficient of $z^{n+1}$ in it is asymptotic to $\frac{5^n}{n}$. This has to be divided (as derived inside Formula (5.2)) by

$$5^{n+\frac{1}{2}} \frac{1}{2\sqrt{\pi n^{3/2}}}$$

with the result

$$2 \frac{5^n}{n} \frac{1}{5^{n+\frac{1}{2}}} \sqrt{\pi n^{3/2}} = \frac{2}{\sqrt{5}} \sqrt{\pi n}.$$

Note that the constant in front of $\sqrt{\pi n}$ for ordered trees is $\frac{2}{\sqrt{4}} = 1$, so the average height for marked ordered trees is indeed a bit smaller thanks to the extra markings.
## 6. A BIJECTION BETWEEN MULTI-EDGE TREES AND 3-COLOURED MOTZKIN PATHS

Multi-edge trees are like ordered (planar, plane, ...) trees, but instead of edges there are multiple edges. When drawing such a tree, instead of drawing, say 5 parallel edges, we just draw one edge and put the number 5 on it as a label. These trees were studied in [14, 24]. We also considered this model in our unpublished preprint arXiv:2105.03350; the bijection presented hereafter is new.

The generating function $F(z)$ (where one counts edges) satisfies

$$F(z) = \sum_{k \ge 0} \left( \frac{z}{1-z} F(z) \right)^k = \frac{1}{1 - \frac{z}{1-z} F(z)},$$

whence

$$F(z) = \frac{1 - z - \sqrt{1 - 6z + 5z^2}}{2z} = 1 + z + 3z^2 + 10z^3 + 36z^4 + 137z^5 + 543z^6 + \dots$$

The coefficients form sequence A002212 in the OEIS [52].

A Motzkin path consists of up-steps, down-steps, and horizontal steps; see sequence A091965 in [52] and the references given there. As Dyck paths, they start at the origin and end, after $n$ steps again at the $x$-axis, but are not allowed to go below the $x$-axis. A 3-coloured Motzkin path is built as a Motzkin path, but there are 3 different types of horizontal steps, which we call red, green, blue. The generating function $M(z)$ satisfies

$$M(z) = 1 + 3zM(z) + z^2M(z)^2 = \frac{1 - 3z - \sqrt{1 - 6z + 5z^2}}{2z^2}, \quad \text{or} \quad F(z) = 1 + zM(z).$$

So multi-edge trees with $N$ edges (counting the multiplicities) correspond to 3-coloured Motzkin paths of length $N-1$.

The purpose of this section is to describe a bijection. It transforms trees into paths, but all steps are reversible.

**The details.** As a first step, the multiplicities will be ignored, and the tree then has only $n$ edges. The standard translation of such tree into the world of Dyck paths, which is in every book on combinatorics, leads to a Dyck path of length $2n$. Then the Dyck path will be transformed bijectively to a 2-coloured Motzkin path of length $n-1$ (the colours used are red and green). This transformation plays a prominent role in [13], but is most likely much older. I believe that people like Viennot know this for 40 years. I would be glad to get a proper historic account from the gentle readers.

The last step is then to use the third colour (blue) to deal with the multiplicities.

The first up-step and the last down-step of the Dyck path will be deleted. Then, the remaining $2n-2$ steps are coded pairwise into a 2-Motzkin path of length $n-1$:

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

The last step is to deal with the multiplicities. If an edge is labelled with the number $a$, we will insert $a-1$ horizontal blue steps in the following way. Since there are currently $n-1$ symbols in the path, we have $n$ possible positions to enter something (in the beginning, in the end, between symbols). We go through the tree in pre-order, and enter the multiplicities one by one using the blue horizontal steps.
To illustrate this procedure, we give in Table 1 the list of 10 multi-edge trees with 3 edges, and the corresponding 3-Motzkin paths of length 2, with intermediate steps completely worked out.

<table><caption>TABLE 1. First row is a multi-edge tree with 3 edges, second row is the standard Dyck path (multiplicities ignored), third row is cutting off first and last step, and then translated pairs of steps, fourth row is inserting blue horizontal edges, according to multiplicities.</caption><thead><tr><th>Multi-edge tree</th><th>Dyck path</th><th>2-Motzkin path</th><th>blue edges added</th></tr></thead><tbody><tr><td>1<br>1<br>1</td><td></td><td></td><td></td></tr><tr><td>2<br>1</td><td></td><td></td><td></td></tr><tr><td>1<br>2</td><td></td><td></td><td></td></tr><tr><td>3</td><td></td><td></td><td></td></tr><tr><td>1<br>1</td><td></td><td></td><td></td></tr><tr><td>2<br>1</td><td></td><td></td><td></td></tr><tr><td>1<br>2</td><td></td><td></td><td></td></tr><tr><td>1<br>1</td><td></td><td></td><td></td></tr><tr><td>1<br>1</td><td></td><td></td><td></td></tr></tbody></table>

**Connecting unary-binary trees with multi-edge trees.** This is not too difficult: We start from multi-edge trees, and ignore the multiplicities at the moment. Then we apply the classical rotation correspondence (also called: natural correspondence). Then we add vertical edges, if the multiplicity is higher than 1. To be precise, if there is a node, and an edge with multiplicity $a$ leads to it from the top, we insert $a - 1$ extra unary nodes in a chain on the top, and connect them with unary branches.
<table><caption>In Table 2 below, we illustrate this procedure on 10 objects.</caption><thead><tr><th>Multi-edge tree</th><th>Binary tree (rotation)</th><th>vertical edges added</th></tr></thead><tbody><tr><td><pre>  1
  |
  1
  |
  1</pre></td><td><pre>  ●
  |
  ●
  |
  ●</pre></td><td><pre>  ●
  |
  ●
  |
  ●</pre></td></tr><tr><td><pre>  2
  |
  1</pre></td><td><pre>  ●
  |
  ●
  |
  ●</pre></td><td><pre>  ●
  |
  ●
  |
  ●</pre></td></tr><tr><td><pre>  1
  |
  2</pre></td><td><pre>  ●
  |
  ●
  |
  ●</pre></td><td><pre>  ●
  |
  ●
  |
  ●</pre></td></tr><tr><td><pre>  3</pre></td><td><pre>  ●</pre></td><td><pre>  ●</pre></td></tr><tr><td><pre>  1
  |  1
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td></tr><tr><td><pre>  2
  |  1
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td></tr><tr><td><pre>  1
  |  2
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td></tr><tr><td><pre>  1
  |  1
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td></tr><tr><td><pre>  1
  |  1
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td><td><pre>  ●
  |  ●
  ●   ●
  |   |
  ●   ●</pre></td></tr></tbody></table>

<caption>TABLE 2. First row is a multi-edge tree with 3 edges, second row the corresponding binary tree, according to the classical rotation correspondence, ignoring the unary branches. Third row is inserting extra horizontal edges when the multiplicities are higher than 1.</caption>
7. THE COMBINATORICS OF SKEW DYCK PATHS

Let us come back to skew Dyck paths, which we introduced in Section 5 as objects in bijection with marked ordered trees. As we saw, skew Dyck paths are a variation of Dyck paths, where additionally to steps (1, 1) and (1, -1) a south-west step (-1, -1) is also allowed, provided that the path does not intersect itself. Also, like for Dyck paths, it must never go below the x-axis and end eventually (after 2n steps) on the x-axis. These paths were considered in [6, 12, 26, 44]. We extend here on our paper [46], giving here additional results on the prefixes of these paths. The enumerating sequence is 1, 1, 3, 10, 36, 137, 543, 2219, 9285, 39587, 171369, . . . , which is A002212 in the OEIS [52].

Let us now give a more thorough analysis of skew Dyck paths, using generating functions and the kernel method. Here is the list of the 10 skew Dyck paths consisting of 6 steps:

<figure><img src="image_44.png" /><figcaption>FIGURE 3. All 10 skew Dyck paths of length 6.</figcaption></figure>

We prefer to work with the equivalent model (resembling more traditional Dyck paths) where we replace each step (-1, -1) by (1, -1) but label it red (see Figure 4, and compare with Figure 3):

<figure><img src="image_45.png" /><figcaption>FIGURE 4. The 10 paths of length 6 redrawn, with red south-east edges instead of south-west edges.</figcaption></figure>

The rules to generate such decorated Dyck paths are: each edge (1, -1) may be black or red, but ∧ and ∨ are forbidden.

Our interest is in particular in *partial* decorated Dyck paths, ending at level *j*, for fixed *j* ≥ 0; the instance *j* = 0 is the classical case. The analysis of partial skew Dyck paths was recently started by Baril et al. in [6] (using the notion ‘prefix of a skew Dyck path’) using Riordan arrays instead of our kernel method. The latter gives us *bivariate* generating functions, from which it is easier to draw conclusions. Two variables, *z* and *u*, are used, where *z* marks the length of the path and *j* marks the end-level. We briefly mention that one can, using a third variable *w*, also count the number of red edges.

Again, once all generating functions are explicitly known, many corollaries can be derived in a standard fashion. We only do this in a few instances. But we would like to emphasize that the substitution $z = \frac{v}{1+3v+v^2}$, which was used in [24, 44] allows to write *explicit enumerations*, using the notion of a (weighted) trinomial coefficient:

$$ \binom{n; 1, 3, 1}{k} := [t^k](1 + 3t + t^2)^n. $$
**Generating functions and the kernel method.** We catch the essence of a decorated Dyck path using a state-diagram:

<figure><img src="image_46.png" /><figcaption>FIGURE 5. Three layers of states according to the type of steps leading to them (up, down-black, down-red).</figcaption></figure>

It has three types of states, with $j$ ranging from 0 to infinity; in the drawing, only $j = 0..8$ is shown. The first layer of states refers to an up-step leading to a state, the second layer refers to a black down-step leading to a state and the third layer refers to a red down-step leading to a state. We will work out generating functions describing all paths leading to a particular state. We will use the notations $f_j, g_j, h_j$ for the three respective layers, from top to bottom. Note that the syntactic rules of forbidden patterns $\wedge$ and $\vee$ can be clearly seen from the picture. The functions depend on the variable $z$ (marking the number of steps), but mostly we just write $f_j$ instead of $f_j(z)$, etc.

The following recursions can be read off immediately from the diagram:

$$f_0 = 1, \quad f_{i+1} = zf_i + zg_i, \quad i \ge 0,$$

$$g_i = zf_{i+1} + zg_{i+1} + zh_{i+1}, \quad i \ge 0,$$

$$h_i = zg_{i+1} + zh_{i+1}, \quad i \ge 0.$$

And now it is time to introduce the promised bivariate generating functions:

$$F(z, u) = \sum_{i \ge 0} f_i(z)u^i, \quad G(z, u) = \sum_{i \ge 0} g_i(z)u^i, \quad H(z, u) = \sum_{i \ge 0} h_i(z)u^i.$$

Again, often we just write $F(u)$ instead of $F(z, u)$ and treat $z$ as a 'silent' variable. Summing the recursions leads to

$$\sum_{i \ge 0} u^i f_{i+1} = \sum_{i \ge 0} u^i z f_i + \sum_{i \ge 0} u^i z g_i,$$

$$\sum_{i \ge 0} u^i g_i = \sum_{i \ge 0} u^i z f_{i+1} + \sum_{i \ge 0} u^i z g_{i+1} + \sum_{i \ge 0} u^i z h_{i+1},$$

$$\sum_{i \ge 0} u^i h_i = \sum_{i \ge 0} u^i z h_{i+1} + \sum_{i \ge 0} u^i z g_{i+1}.$$

This can be rewritten as

$$\frac{1}{u}(F(u) - 1) = zF(u) + zG(u),$$

$$G(u) = \frac{z}{u}(F(u) - 1) + \frac{z}{u}(G(u) - G(0)) + \frac{z}{u}(H(u) - H(0)),$$

$$H(u) = \frac{z}{u}(G(u) - G(0)) + \frac{z}{u}(H(u) - H(0)).$$

Such systems of equations having more unknowns than equations can be solved with the kernel method (see [37] for a gentle example-driven introduction to this method).
We begin by rewriting our system as

$$F(u) = \frac{z^2 u G(0) + z^2 u H(0) + z^2 u - u - z^3 + 2z}{-z^3 - u + 2z + zu^2 - z^2 u}, \qquad (7.1)$$

$$G(u) = \frac{z(H(0) - uzH(0) + z^2 + G(0) - zuG(0) - zu)}{-z^3 - u + 2z + zu^2 - z^2 u}, \qquad (7.2)$$

$$H(u) = \frac{z(-uzH(0) - z^2 - zuG(0) + G(0) - z^2H(0) + H(0) - z^2G(0))}{-z^3 - u + 2z + zu^2 - z^2 u}, \qquad (7.3)$$

The denominator is the same for each equation and it factors as $z(u - r_1)(u - r_2)$, with

$$r_1 = \frac{1 + z^2 + \sqrt{1 - 6z^2 + 5z^4}}{2z}, \quad r_2 = \frac{1 + z^2 - \sqrt{1 - 6z^2 + 5z^4}}{2z}.$$

Consider Equation (7.1), since $F(u)$ is a power series in $z$, the factor $u - r_2$ in the denominator is “bad”$^2$, thus this factor must also be a factor of the numerator (seen as a polynomial of degree 1 in $u$). This implies

$$G(0) = -H(0) - 1 + \frac{1}{z^2} + \frac{z - 2/z}{r_2}.$$

Applying the same principle to either (7.2) or (7.3), we get after simplification

$$H(0) = \frac{1 - 4z^2 + z^4 + (z^2 - 1)\sqrt{1 - 6z^2 + 5z^4}}{2 - z^2}.$$

Thus, with $W = \sqrt{1 - 6z^2 + 5z^4} = \sqrt{(1 - z^2)(1 - 5z^2)}$, one has

$$F(u) = \frac{-1 - z^2 - W}{2z(u - r_1)} = \frac{1 + z^2 + W}{2zr_1(1 - u/r_1)},$$

$$G(u) = \frac{-1 + z^2 + W}{2z(u - r_1)} = \frac{1 - z^2 - W}{2zr_1(1 - u/r_1)},$$

$$H(u) = \frac{-1 + 3z^2 + W}{2z(u - r_1)} = \frac{1 - 3z^2 - W}{2zr_1(1 - u/r_1)}.$$

The total generating function (summing the 3 cases that lead to the same level) is

$$S(u) = F(u) + G(u) + H(u) = \frac{3 - 3z^2 - W}{2zr_1(1 - u/r_1)}.$$

The coefficient of $u^jz^n$ in $S(u)$ counts the partial paths of length $n$, ending at level $j$. We will write $s_j = [u^j]S(u)$. At this stage, we are only interested in

$$s_j = f_j + g_j + h_j = [u^j] \frac{3 - 3z^2 - W}{2zr_1(1 - u/r_1)} = \frac{3 - 3z^2 - W}{2zr_1^{j+1}},$$

which is the generating function of all (partial) paths ending at level $j$. Parity considerations give us that only coefficients $[z^n]s_j$ are non-zero if $n \equiv j \bmod 2$. To make this more transparent, we set

$$P(z) = zr_1 = \frac{1 + z^2 + \sqrt{1 - 6z^2 + 5z^4}}{2}.$$

<footnote>$^2$In an identity having the form $A(z, u) = B(z, u)/(u - r(z))$, we say that the factor $(u - r(z))$ is bad if $A(u, z)$ is a power series in $z$ while $1/(u - r(z))$ is not.</footnote>
We thus get

$$s_j = f_j + g_j + h_j = z^j \frac{3 - 3z^2 - W}{2P^{j+1}}.$$

Now we read off coefficients. We do this using residues and contour integration. The path of integration, in both variables $x$ resp. $v$ is a small circle or an equivalent contour. The following computation is abbreviated.

$$[z^{2m+j}]s_j = \frac{1}{2\pi i} \oint \frac{dx}{x^{m+1}} \frac{(1+v)(1+2v)}{v^{j+1}(v+2)^{j+1}} (1+3v+v^2)^j \\ = [v^{m+j+1}] \frac{(1+v)^2(1+2v)(1-v)}{(v+2)^{j+1}} (1+3v+v^2)^{m-1+j}.$$

Note that

$$(1 + v)^2(1 + 2v)(1 - v) = -9 + 27(v + 2) - 29(v + 2)^2 + 13(v + 2)^3 - 2(v + 2)^4;$$

consequently

$$[v^k] \frac{(1+v)^2(1+2v)(1-v)}{(v+2)^{j+1}} = -9 \frac{1}{2^{j+1+k}} \binom{-j-1}{k} + 27 \frac{1}{2^{j+k}} \binom{-j}{k} - 29 \frac{1}{2^{j-1+k}} \binom{-j+1}{k} \\ + 13 \frac{1}{2^{j-2+k}} \binom{-j+2}{k} - 2 \frac{1}{2^{j-3+k}} \binom{-j+3}{k} =: \lambda_{j;k}.$$

With this abbreviation $\lambda_{j;k}$ we find

$$[v^{m+j+1}] \frac{(1+v)^2(1+2v)(1-v)}{(v+2)^{j+1}} (1+3v+v^2)^{m-1+j} = \sum_{k=0}^{m+j+1} \lambda_{j;k} \binom{m-1+j; 1, 3, 1}{m+j+1-k}.$$

This is not extremely pretty but it is explicit and as good as it gets. Here are the first few generating functions:

$$\bullet \quad s_0 = 1 + z^2 + 3z^4 + 10z^6 + 36z^8 + 137z^{10} + 543z^{12} + \dots,$$

$$\bullet \quad s_1 = z + 2z^3 + 6z^5 + 21z^7 + 79z^9 + 311z^{11} + 1265z^{13} + \dots.$$

We could also give such lists for the functions $f_j$, $g_j$, $h_j$, if desired. We summarize the essential findings of this section in the following theorem.

**Theorem 7.1.** *The generating function of decorated (partial) Dyck paths, consisting of n steps, ending on level j, is given by*

$$S(z, u) = \frac{3 - 3z^2 - \sqrt{1 - 6z^2 + 5z^4}}{2zr_1(1 - u/r_1)},$$

with

$$r_1 = \frac{1 + z^2 + \sqrt{1 - 6z^2 + 5z^4}}{2z}.$$

Furthermore

$$[u^j]S(z, u) = \frac{3 - 3z^2 - \sqrt{1 - 6z^2 + 5z^4}}{2zr_1^{j+1}}.$$
**Open ended paths.** If we do not specify the end of the paths, in other words we sum over all $j \ge 0$, then at the level of generating functions this is very easy, since we only have to set $u := 1$. We find

$$S(1) = - \frac{(z+1)(z^2+3z-2) + (z+2)\sqrt{1-6z^2+5z^4}}{2z(z^2+2z-1)} \\ = 1 + z + 2z^2 + 3z^3 + 7z^4 + 11z^5 + 26z^6 + 43z^7 + 102z^8 + 175z^9 + 416z^{10} + \dots$$

**Counting red edges.** We can use an extra variable, $w$, to count additionally the red edges that occur in a path. We use the same letters for generating functions. Eventually, the coefficient $[z^n u^j w^k]S$ is the number of (partial) paths consisting of $n$ steps, leading to level $j$, and having passed $k$ red edges. The endpoint of the original skew path has then coordinates $(n - 2k, j)$. The computations are very similar, and we only sketch the key steps.

$$f_0 = 1, \quad f_{i+1} = zf_i + zg_i, \quad i \ge 0,$$

$$g_i = zf_{i+1} + zg_{i+1} + zh_{i+1}, \quad i \ge 0,$$

$$h_i = wzg_{i+1} + wzh_{i+1}, \quad i \ge 0;$$

$$\frac{1}{u}(F(u) - 1) = zF(u) + zG(u),$$

$$G(u) = \frac{z}{u}(F(u) - 1) + \frac{z}{u}(G(u) - G(0)) + \frac{z}{u}(H(u) - H(0)),$$

$$H(u) = \frac{wz}{u}(G(u) - G(0)) + \frac{wz}{u}(H(u) - H(0));$$

$$F(u) = \frac{z^2 u G(0) + z^2 u H(0) + z^2 u - u - wz^3 + z + wz}{-wz^3 - u + z + wz + zu^2 - wz^2 u},$$

$$G(u) = \frac{z(H(0) - uzH(0) + wz^2 + G(0) - zuG(0) - zu)}{-wz^3 - u + z + wz + zu^2 - wz^2 u},$$

$$H(u) = \frac{wz(-uzH(0) - z^2 - zuG(0) + G(0) - z^2H(0) + H(0) - z^2G(0))}{-wz^3 - u + z + wz + zu^2 - wz^2 u}.$$

The denominator factors as $z(u - r_1)(u - r_2)$, with

$$r_1 = \frac{1 + wz^2 + \sqrt{1 - (4 + 2w)z^2 + (4w + w^2)z^4}}{2z},$$

$$r_2 = \frac{1 + wz^2 - \sqrt{1 - (4 + 2w)z^2 + (4w + w^2)z^4}}{2z}.$$

Note the factorization $1 - (4 + 2w)z^2 + (4w + w^2)z^4 = (1 - z^2w)(1 - (4 + w)z^2)$. Since the factor $u - r_2$ in the denominator is “bad,” it must also cancel in the numerators. From this we eventually find, with the abbreviation $W = \sqrt{1 - (4 + 2w)z^2 + (4w + w^2)z^4}$:

$$F(u) = \frac{-1 - wz^2 - W}{2z(u - r_1)},$$

$$G(u) = \frac{-1 + wz^2 + W}{2z(u - r_1)},$$

$$H(u) = \frac{-1 + (2 + w)z^2 + W}{2z(u - r_1)}.$$
The total generating function is

$$S(u) = F(u) + G(u) + H(u) = \frac{-2 - w + z^2(w + w^2) + wW}{2z(u - r_1)}.$$

The special case $u = 0$ (return to the $x$-axis) is to be noted:

$$S(0) = \frac{-2 - w + z^2(w + w^2) + wW}{-2zr_1} = \frac{1 - wz^2 - W}{2z^2}.$$

Since there are only even powers of $z$ in this function, we replace $x = z^2$ and get

$$S(0) = \frac{1 - wx - \sqrt{1 - (4 + 2w)x + (4w + w^2)x^2}}{2x} \\ = 1 + x + (w + 2)x^2 + (w^2 + 4w + 5)x^3 + (w^3 + 6w^2 + 15w + 14)x^4 + \cdots$$

Compare the factor $(w^2 + 4w + 5)$ with the earlier drawing of the 10 paths.

There is again a substitution that allows for better results:

$$z = \frac{v}{1 + (2 + w)v + v^2}, \quad \text{then} \quad S(0) = 1 + v.$$

Reading off coefficients can now be done using modified trinomial coefficients:

$$\binom{n; 1, 2 + w, 1}{k} = [t^k](1 + (2 + w)t + t^2)^n.$$

Again, we use contour integration to extract coefficients:

$$\begin{aligned} [x^n](1 + v) &= \frac{1}{2\pi i} \oint \frac{dx}{x^{n+1}}(1 + v) \\ &= \frac{1}{2\pi i} \oint \frac{dx}{v^{n+1}} \frac{1 - v^2}{(1 + (2 + w)v + v^2)^2} (1 + (2 + w)v + v^2)^{n+1}(1 + v) \\ &= [v^n](1 - v)(1 + v)^2(1 + (2 + w)v + v^2)^{n-1} \\ &= \binom{n-1; 1, 2 + w, 1}{n} + \binom{n-1; 1, 2 + w, 1}{n-1} \\ &\quad - \binom{n-1; 1, 2 + w, 1}{n-2} - \binom{n-1; 1, 2 + w, 1}{n-3}. \end{aligned}$$

Now we want to count the average number of red edges. For that, we differentiate $S(0)$ with respect to $w$, and set $w := 1$. This leads to

$$\frac{-1 + 6x - 5x^2 + (1 + 3x)\sqrt{1 - 6x + 5x^2}}{2(1 - x)(1 - 5x)}.$$

A simple application of singularity analysis leads to $\frac{\frac{1}{2\sqrt{5}}[x^n]\frac{1}{\sqrt{1-5x}}}{-\sqrt{5}[x^n]\sqrt{1-5x}} \sim \frac{n}{5}$.

So, a random path consisting of $2n$ steps has about $n/5$ red steps, on average. For readers who are not familiar with singularity analysis of generating functions [16, 20], we just mention that one determines the local expansion around the dominating singularity, which is at $z = \frac{1}{5}$ in our instance. In the denominator, we just have the total number of skew Dyck paths, according to the sequence A002212 in the OEIS [52]. In the example of Figure 2, the exact average is 6/10, which curiously is exactly the same as 3/5.
We finish the discussion by considering fixed powers of $w$ in $S(0)$, counting skew Dyck paths consisting of zero, one, two, three, ... red edges. We find

$$[w^0]S(0) = \frac{1 - \sqrt{1 - 4x}}{2x}, \quad [w^1]S(0) = \frac{1 - 2x - \sqrt{1 - 4x}}{2\sqrt{1 - 4x}},$$

$$[w^2]S(0) = \frac{x^3}{(1 - 4x)^{3/2}}, \quad [w^3]S(0) = \frac{x^4(1 - 2x)}{(1 - 4x)^{5/2}}, \quad \&c.$$

The generating function $[w^0]S(0)$ is of course the generating function of Catalan numbers, since no red edges just means: ordinary Dyck paths. We can also conclude that the asymptotic behaviour is of the form $n^{k-3/2}4^n$, where the polynomial contribution gets higher, but the exponential growth stays the same: $4^n$. This is compared to the scenario of an arbitrary number of red edges, when we get an exponential growth of the form $5^n$.

**Dual skew Dyck paths.** The mirrored version of skew Dyck paths with two types of up-steps, $(1, 1)$ and $(-1, 1)$ are also cited among the objects in A002212 in the OEIS [52]. We call them dual skew paths and drop the ‘dual’ when it isn’t necessary. When the paths come back to the $x$-axis, no new enumeration is necessary, but this is no longer true for paths ending at level $j$.

Here is a list of the 10 skew paths consisting of 6 steps:

<figure><img src="image_47.png" /><figcaption>FIGURE 6. All 10 dual skew Dyck paths of length 6 (consisting of 6 steps).</figcaption></figure>

We prefer to work with the equivalent model (resembling more traditional Dyck paths) where we replace each step $(-1, -1)$ by $(1, -1)$ but label it blue. Here is the list of the 10 paths again (Figure 2):

<figure><img src="image_48.png" /><figcaption>FIGURE 7. All 10 dual skew Dyck paths of length 6 (consisting of 6 steps).</figcaption></figure>

The rules to generate such decorated Dyck paths are: Each edge $(1, -1)$ may be black or blue, but $\vee$ and $\wedge$ are forbidden.

Our interest is in particular in *partial* decorated Dyck paths, ending at level $j$, for fixed $j \ge 0$; the instance $j = 0$ is the classical case.

As before, two variables, $z$ and $u$, are used, where $z$ marks the length of the path and $j$ marks the end-level. We briefly mention that one can, using a third variable $w$, also count the number of blue edges. The substitution $x = \frac{v}{1 + 3v + v^2}$ is again the key to the success.
**Generating functions and the kernel method.** We catch the essence of a decorated (dual skew) Dyck path using a state-diagram:

<figure><img src="image_49.png" /><figcaption>FIGURE 8. Three layers of states according to the type of steps leading to them (down, up-black, up-blue).</figcaption></figure>

It has three types of states, with $j$ ranging from 0 to infinity; in the drawing, only $j = 0..8$ is shown. The first layer of states refers to an up-step leading to a state, the second layer refers to a black down-step leading to a state and the third layer refers to a blue down-step leading to a state. We will work out generating functions describing all paths leading to a particular state. We will use the notations $c_j, a_j, b_j$ for the three respective layers, from top to bottom. Note that the syntactic rules of forbidden patterns $\wedge$ and $\vee$ can be clearly seen from the picture. The functions depend on the variable $z$ (marking the number of steps), but mostly we just write $a_j$ instead of $a_j(z)$, etc.

The following recursions can be read off immediately from the diagram:

$$ 
\begin{aligned}
a_0 &= 1, & a_{i+1} &= za_i + zb_i + zc_i, & i \ge 0, \\
b_i &= za_{i+1} + zb_{i+1}, & i \ge 0, \\
c_{i+1} &= za_i + zc_i, & i \ge 0.
\end{aligned}
$$

And now it is time to introduce the bivariate generating functions:

$$ A(z, u) = \sum_{i \ge 0} a_i(z) u^i, \quad B(z, u) = \sum_{i \ge 0} b_i(z) u^i, \quad C(z, u) = \sum_{i \ge 0} c_i(z) u^i. $$

Summing the recursions leads to

$$ 
\begin{aligned}
\sum_{i \ge 0} u^i a_i &= 1 + u \sum_{i \ge 0} u^i (za_i + zb_i + zc_i) = 1 + uzA(u) + uzB(u) + uzC(u), \\
\sum_{i \ge 0} u^i b_i &= \sum_{i \ge 0} u^i (za_{i+1} + zb_{i+1}) = \frac{z}{u} \sum_{i \ge 1} u^i a_i + \frac{z}{u} \sum_{i \ge 1} u^i b_i, \\
\sum_{i \ge 1} u^i c_i &= uz \sum_{i \ge 0} u^i a_i + uz \sum_{i \ge 0} u^i c_i.
\end{aligned}
$$

This can be rewritten as

$$ 
\begin{aligned}
A(u) &= 1 + uzA(u) + uzB(u) + uzC(u), \\
B(u) &= \frac{z}{u}(A(u) - a_0) + \frac{z}{u}(B(u) - b_0), \\
C(u) &= c_0 + uzA(u) + uzC(u).
\end{aligned}
$$

Note that $a_0 = 1, c_0 = 0$. Simplification leads to

$$ C(u) = \frac{uzA(u)}{1 - uz} \quad \text{and} \quad B(u) = \frac{z(A(u) - 1 - B(0))}{u - z}. $$
This leaves us with just one equation

$$A(u) = \frac{(z - u + uz^2 + uz^2B(0))(uz - 1)}{u^2z^3 + uz^2 - 2u^2z - z + u}.$$

This is again a typical application of the kernel method: One writes

$$u^2z^3 + uz^2 - 2u^2z - z + u = z(z^2 - 2)(u - s_1)(u - s_2).$$

The denominator thus factors as $2z(z^2 - 2)(u - s_1)(u - s_2)$, with

$$s_1 = \frac{1 + z^2 + \sqrt{1 - 6z^2 + 5z^4}}{2z(2 - z^2)}, \quad s_2 = \frac{1 + z^2 - \sqrt{1 - 6z^2 + 5z^4}}{2z(2 - z^2)}.$$

Note that $s_1s_2 = \frac{1}{2-z^2}$. Since the factor $u - s_2$ in the denominator is “bad,” it must also cancel in the numerators. We get $B(0) = \frac{zs_2}{1-2zs_2}$ and, again with the abbreviation $W = \sqrt{1 - 6z^2 + 5z^4}$,

$$A(u) = \frac{(1 - uz)(1 + z^2 + W)}{2z(z^2 - 2)(u - s_1)}, \quad B(u) = \frac{1 - 2z^2 - W}{z(2 - z^2)(u - s_1)}, \quad C(u) = \frac{1 + z^2 + W}{2(z^2 - 2)} \frac{u}{u - s_1},$$

and for the function of main interest

$$G(u) = A(u) + B(u) + C(u) = \frac{3z^2 - 3 + W}{2z(2 - z^2)(u - s_1)}.$$

Since one has

$$\frac{1}{s_1} = \frac{1 + z^2 - \sqrt{1 - 6z^2 + 5z^4}}{2z} = zS \quad \text{and} \quad \frac{1}{s_2} = \frac{1 + z^2 + \sqrt{1 - 6z^2 + 5z^4}}{2z},$$

we then get

$$[u^j]G(u) = [u^j]\frac{3z^2 - 3 + W}{2z(z^2 - 2)s_1(1 - u/s_1)} = \frac{3z^2 - 3 + W}{2z(z^2 - 2)s_1^{j+1}} = \frac{3z^2 - 3 + W}{2(z^2 - 2)}z^j S^{j+1}.$$

So $[u^j]G(u)$ contains only powers of the form $z^{j+2N}$. Now we continue

$$\begin{align*}
[z^{j+2N} u^j] G(u) &= [z^{2N}] \frac{3z^2 - 3 + W}{2(z^2 - 2)} S^{j+1} \\
&= [x^N] \frac{3x - 3 + \sqrt{1 - 6x + 5x^2}}{2(x - 2)} \left( \frac{1 + x - \sqrt{1 - 6x + 5x^2}}{2x} \right)^{j+1} \\
&= [x^N] (v + 1)(v + 2)^j
\end{align*}$$

which is the generating function of all (partial) paths ending at level j.

Now we read off coefficients as before:

$$\begin{align*}
[z^{j+2N} u^j] G(u) &= [x^N](v+1)(v+2)^j \\
&= \frac{1}{2\pi i} \oint \frac{dx}{x^{N+1}} (v+1)(v+2)^j \\
&= \frac{1}{2\pi i} \oint \frac{dv}{v^{N+1}} (1+3v+v^2)^{N+1} \frac{(1-v^2)}{(1+3v+v^2)^2} (v+1)(v+2)^j \\
&= [v^N](1+3v+v^2)^{N-1}(1-v)(1+v)^2(v+2)^j.
\end{align*}$$
Note that

$$ (1 - v)(1 + v)^2 = 3 - 7(v + 2) + 5(v + 2)^2 - (v + 2)^3; $$

consequently

$$ [z^{j+2N} u^j] G(u) = [v^N](1 + 3v + v^2)^{N-1} [3 - 7(v + 2) + 5(v + 2)^2 - (v + 2)^3](v + 2)^j. $$

We abbreviate:

$$ \begin{aligned} \mu_{j;k} &= [v^k] [3(v+2)^j - 7(v+2)^{j+1} + 5(v+2)^{j+2} - (v+2)^{j+3}] \\ &= 3 \binom{j}{k} 2^{j-k} - 7 \binom{j+1}{k} 2^{j+1-k} + 5 \binom{j+2}{k} 2^{j+2-k} - \binom{j+3}{k} 2^{j+3-k}. \end{aligned} $$

With this notation we get

$$ [z^{j+2N} u^j] G(u) = \sum_{0 \le k \le N-1} \mu_{j;k} \binom{N-1; 1, 3, 1}{N-k}. $$

Here are the first few generating functions:

$$ \bullet \quad G_0 = 1 + z^2 + 3z^4 + 10z^6 + 36z^8 + 137z^{10} + 543z^{12} + 2219z^{14} + \dots, $$

$$ \bullet \quad G_1 = 2z + 3z^3 + 10z^5 + 36z^7 + 137z^9 + 543z^{11} + 2219z^{13} + 9285z^{15} + \dots, $$

$$ \bullet \quad G_2 = 4z^2 + 8z^4 + 29z^6 + 111z^8 + 442z^{10} + 1813z^{12} + 7609z^{14} + 32521z^{16} + \dots, $$

$$ \bullet \quad G_3 = 8z^3 + 20z^5 + 78z^7 + 315z^9 + 1306z^{11} + 5527z^{13} + 23779z^{15} + 103699z^{17} + \dots. $$

We could also give such lists for the functions $a_j$, $b_j$, $c_j$, if desired. We summarize the essential findings of the rest of this section:

**Theorem 7.2.** The generating function of decorated (partial) dual skew Dyck paths, consisting of n steps, ending on level j, is given by

$$ G(z, u) = \frac{3z^2 - 3 + \sqrt{1 - 6z^2 + 5z^4}}{2z(2 - z^2)(u - s_1)}, $$

with

$$ s_1 = \frac{2z}{1 + z^2 - \sqrt{1 - 6z^2 + 5z^4}}. $$

Furthermore

$$ [u^j] G(z, u) = \frac{3z^2 - 3 + \sqrt{1 - 6z^2 + 5z^4}}{2(z^2 - 2)} z^j S^{j+1}, $$

with

$$ S = \frac{1 + z^2 - \sqrt{1 - 6z^2 + 5z^4}}{2z^2}. $$

**Open ended paths.** If we do not specify the end of the paths, in other words we sum over all $j \ge 0$, then at the level of generating functions this is very easy, since we only have to set $u := 1$. We find

$$ \begin{aligned} G(1) &= \frac{(1+z)(1-3z)}{2z(z^2+2z-1) - \sqrt{1-6z^2+5z^4}} \\ &= 1 + 2z + 5z^2 + 11z^3 + 27z^4 + 62z^5 + 151z^6 + 354z^7 + 859z^8 + 2036z^9 + \dots. \end{aligned} $$
**Counting blue edges.** We use an extra variable $w$ as before to count additionally the blue edges that occur in a path. Eventually, the coefficient $[z^n u^j w^k]S$ is the number of (partial) paths consisting of $n$ steps, leading to level $j$, and having passed $k$ blue edges. The endpoint of the original skew path has then coordinates $(n-2k, j)$. The computations are very similar, and we only sketch the key steps.

$$ 
\begin{aligned}
a_0 &= 1, \quad a_{i+1} = za_i + zb_i + zc_i, \quad i \ge 0, \\
b_i &= za_{i+1} + zb_{i+1}, \quad i \ge 0, \\
c_{i+1} &= wza_i + wzc_i, \quad i \ge 0.
\end{aligned}
$$

This leads to

$$ 
\begin{aligned}
A(u) &= 1 + uzA(u) + uzB(u) + uzC(u), \\
B(u) &= \frac{z}{u}(A(u) - a_0) + \frac{z}{u}(B(u) - b_0), \\
C(u) &= c_0 + wuzA(u) + wuzC(u).
\end{aligned}
$$

Solving,

$$ S(u) = A(u) + B(u) + C(u) = \frac{u - wuz^2 - zA(0) - zB(0) + uwz^2A(0) + uwz^2B(0)}{u^2z^3w + u - wu^2z - u^2z - z + wuz^2}. $$

The denominator factors as $-z(1+w-z^2w)(u-s_1)(u-s_2)$, with

$$ 
\begin{aligned}
s_1 &= \frac{1 + z^2w + \sqrt{1 - 2z^2w + z^4w^2 - 4z^2 + 4z^4w}}{2z(1 + w - z^2w)}, \\
s_2 &= \frac{1 + z^2w - \sqrt{1 - 2z^2w + z^4w^2 - 4z^2 + 4z^4w}}{2z(1 + w - z^2w)}.
\end{aligned}
$$

Note the factorization $1 - (4+2w)z^2 + (4w+w^2)z^4 = (1-z^2w)(1-(4+w)z^2)$. Since the factor $u-r_2$ in the denominator is “bad,” it must also cancel in the numerators. From this we eventually find, with the abbreviation $W = \sqrt{1-(4+2w)z^2+(4w+w^2)z^4}$

$$ G(0) = \frac{1 - z^2w - W}{2z^2}, $$

and further

$$ G(u) = \frac{w - z^2w^2 - wW + 2 - 2z^2w}{2z(-w - 1 + z^2w)(u - s_1)}. $$

The special case $u=0$ (return to the $x$-axis) is to be noted:

$$ G(0) = 1 + z^2 + (w + 2) z^4 + (w^2 + 4w + 5) z^6 + (w + 2) (w^2 + 4w + 7) z^8 + \dots $$

Compare the factor $(w^2+4w+5)$ with the earlier drawing of the 10 paths. There is again a substitution that allows for better results:

$$ z = \frac{v}{1 + (2 + w)v + v^2}, \quad \text{then} \quad G(0) = 1 + v. $$

Since $S(u) = G(u)$ with $S(u)$ from the first part of the paper, as it means the same objects, read from left to right resp. from right to left, no new analysis is required.
## 8. MORE ABOUT MOTZKIN PATHS

Several other interesting problems related to Motzkin paths were considered in the literature. To wet the readers appetite, we briefly mention below the Retakh paths, the amplitude of paths, and then analyze in more detail skew Motzkin paths.

**Retakh's Motzkin paths.** V. Retakh [15] introduced a restricted class of Dyck paths: Peaks are only allowed on level 1 and on even-numbered levels. For the analysis of this class using generating functions, including also the average height and the number of leaves; see [43].

**The amplitude of Motzkin paths.** Another interesting parameter is the amplitude of Motzkin paths that was recently analyzed in [45]. Here, we want to give a few introductory remarks. The average height of a random Motzkin path of length $n$ was considered in an early paper [35], it is asymptotically given by $\sqrt{\frac{\pi n}{3}}$.

In the recent paper [8] an interesting new concept was introduced: the *amplitude*. It is basically twice the height, but with a twist. If there exists a horizontal step on level $h$, which is the height, the amplitude is $2h + 1$, otherwise it is $2h$. To clarify the concept, we created a list of all 9 Motzkin paths of length 4 with height and amplitude given.

<table><thead><tr><th>Motzkin path</th><th>horizontal on maximal level</th><th>height</th><th>amplitude</th></tr></thead><tbody><tr><td></td><td>Yes</td><td>0</td><td>1</td></tr><tr><td></td><td>No</td><td>1</td><td>2</td></tr><tr><td></td><td>No</td><td>1</td><td>2</td></tr><tr><td></td><td>Yes</td><td>1</td><td>3</td></tr><tr><td></td><td>Yes</td><td>1</td><td>3</td></tr><tr><td></td><td>Yes</td><td>1</td><td>3</td></tr><tr><td></td><td>No</td><td>1</td><td>2</td></tr><tr><td></td><td>No</td><td>1</td><td>2</td></tr><tr><td></td><td>No</td><td>2</td><td>4</td></tr></tbody></table>

The goal of an extended analysis is to investigate this new parameter, using trinomial coefficients $\binom{n,3}{k} = [t^k](1 + t + t^2)^n$ (notation following Comtet's book [7]). The intuitive result that the average amplitude is about twice the average height, can be confirmed.

**Skew Motzkin paths.** This section was written to provide a more complete analysis related to [30]. The methods are of course by now familiar to the readers.

As seen before, Motzkin paths allow additional flat (horizontal) steps of unit length. A skew path allows ‘left’ step $(-1, -1)$ as well, but the path is not allowed to intersect itself. We prefer ‘red’ steps $(1, -1)$; see our analysis in [46]. For Motzkin paths, some analysis was provided in [30]. Here, we provide further analysis that allows to consider partial paths as well, so we do not need to land at the $x$-axis. It uses the kernel method [37].
Apart from being not below the x-axis, the restrictions are that a left (red) step cannot follow or precede an up-step. The situation is best described by a graph (state-diagram); see Figure 9.

<figure><img src="image_50.png" /><figcaption>FIGURE 9. Four layers of states according to the type of steps leading to them. Traditional up-steps and down-steps are black, level-steps are blue, and left steps are red.</figcaption></figure>

In further sections, the asymptotic equivalent for the number of skew Motzkin paths of given size is derived, as well as the height, meaning that the generating function of paths with a bounded height (bounded by H) is given, as well as the average height, which is approximately $\text{const} \cdot \sqrt{n}$, which is typical for families of paths.

*Generating functions for skew Motzkin paths.* We translate the state diagram accordingly; $f_j$, $g_j$, $h_j$, $k_j$ are generating functions in the variable $z$ (marking the length of the path), ending at level $j$. The four families are related to the four layers of states.

$$ 
\begin{aligned}
f_{j+1} &= zf_j + zg_j + zh_j, \quad f_0 = 1, \\
g_j &= zf_{j+1} + zg_{j+1} + zh_{j+1} + zk_{j+1}, \\
h_j &= zf_j + zg_j + zh_j + zk_j, \\
k_j &= zg_{j+1} + zh_{j+1} + zk_{j+1}.
\end{aligned}
$$

Now we introduce bivariate generating functions, namely

$$ F(u) := \sum_{j \ge 0} f_j u^j, \quad G(u) := \sum_{j \ge 0} g_j u^j, \quad H(u) := \sum_{j \ge 0} h_j u^j, \quad K(u) := \sum_{j \ge 0} k_j u^j. $$

The recursions then take this form:

$$ 
\begin{aligned}
F(u) &= 1 + zuF(u) + zuG(u) + zuH(u), \\
uG(u) &= zF(u) + zG(u) + zH(u) + zK(u) - z - zg_0 - zh_0 - zk_0, \\
H(u) &= zF(u) + zG(u) + zH(u) + zK(u), \\
uK(u) &= zG(u) + zH(u) + K(u) - zg_0 + zh_0 + zk_0.
\end{aligned}
$$
Solving the system,

$$
\begin{align*}
F(u) &= \frac{\mathcal{F}}{2z - u + zu - z^2u + zu^2 - z^3 - z^3u}, \\
G(u) &= \frac{\mathcal{G}}{2z - u + zu - z^2u + zu^2 - z^3 - z^3u}, \\
H(u) &= \frac{\mathcal{H}}{2z - u + zu - z^2u + zu^2 - z^3 - z^3u}, \\
K(u) &= \frac{\mathcal{K}}{2z - u + zu - z^2u + zu^2 - z^3 - z^3u},
\end{align*}
$$

with

$$
\mathcal{F} = -z^3 + 2z - u + zu + z^2u + z^2ug_0 + z^2uh_0 + z^2uk_0 + z^3ug_0 + z^3uh_0 + z^3uk_0,
$$

$$
\begin{align*}
\mathcal{G} = -z^2h_0 + z^4 + z^4k_0 - z^2ug_0 - z^2u - z^2k_0 - z^2 - z^2uh_0 - z^2uk_0 + z^4h_0 + z^3 - z^2g_0 + zh_0 \\
\qquad + zg_0 + zk_0 + z^4g_0,
\end{align*}
$$

$$
\mathcal{H} = -z^4 + 2z^2g_0 + 2z^2h_0 + 2z^2k_0 + 2z^2 - zu - z^4g_0 - z^4h_0 - z^4k_0 - z^3ug_0 - z^3uh_0 - z^3uk_0,
$$

$$
\mathcal{K} = zg_0 - z^3 - z^2g_0 - z^2h_0 + zk_0 - z^2k_0 - z^3g_0 - z^3h_0 - z^3k_0.
$$

One cannot immediately insert $u = 0$ to identify the constants, but one can use the kernel method. For that, one factorizes the denominator:

$$
2z - u + zu - z^2u + zu^2 - z^3 - z^3u = z(u - u_1)(u - u_2)
$$

with

$$
u_1 = \frac{1 - z + z^2 + z^3 + (1 + z)W}{2z}, \quad u_2 = \frac{1 - z + z^2 + z^3 - (1 + z)W}{2z}
$$

and

$$
W = \sqrt{(1-z)(1-3z-z^2-z^3)} = \sqrt{1-4z+2z^2+z^4}.
$$

Since $u_2 \sim 2z$ for small $z$, $u - u_2$ is a 'bad' factor and must cancel from both, numerator and denominator. This yields

$$
F(u) = \frac{-1 + z + z^2 + z^2g_0 + z^2h_0 + z^2k_0 + z^3g_0 + z^3h_0 + z^3k_0}{z(u - u_1)},
$$

$$
G(u) = \frac{-z^2 - z^2 g_0 - z^2 h_0 - z^2 k_0}{z(u - u_1)},
$$

$$
H(u) = \frac{-z - z^3 g_0 - z^3 h_0 - z^3 k_0}{z(u - u_1)},
$$

$$
K(u) = \frac{-z^2 g_0 - z^2 h_0 - z^2 k_0}{z(u - u_1)}.
$$

Now we can plug in $u = 0$ and identify the constants:

$$
g_0 = \frac{-z^5 + z^3 W - z^2 - z W + 3z - 1 + W}{2z^2 (-2 + z^2)},
$$

$$
h_0 = - \frac{-z^2 + 2z - 1 + W}{2z},
$$

$$
k_0 = - \frac{-z^4 + z^3 + z^2 W + z W - 3z + 1 - W}{2z^2 (-2 + z^2)}.
$$
Adding these quantities yields

$$1 + g_0 + h_0 + k_0 = -\frac{-z^2 + 2z - 1 + W}{2z^2},$$

which is the generating function of the number of skew Motzkin paths (returning to the x-axis); the series expansion is

$$1+z+2z^2+5z^3+13z^4+35z^5+97z^6+275z^7+794z^8+2327z^9+6905z^{10}+20705z^{11}+\cdots,$$

as already given in [30], the coefficients are the sequence A082582 in the OEIS [52].

We further get

$$F(u) = \frac{-1 + z - z^2 - z^3 + u_2 z}{z(u - u_1)},$$

$$G(u) = \frac{(z - u_2)}{(u - u_1)(1 + z)},$$

$$H(u) = \frac{-1 - z + 2z^2 + z^3 - zu_2}{(u - u_1)(1 + z)},$$

$$K(u) = \frac{z^2 + 2z - u_2}{(u - u_1)(1 + z)}.$$

Altogether,

$$F(z) + G(z) + H(z) + K(z) = \frac{-1 - z + 2z^2 + z^3 - zu_2}{z(u - u_1)(1 + z)}$$

and

$$[u^j](F(z) + G(z) + H(z) + K(z)) = \frac{1 + z - 2z^2 - z^3 + zu_2}{z(1 + z)u_1^{j+1}},$$

which is the generating function of partial skew Motzkin paths, landing on level j. Here are the examples for j = 1, 2, 3, 4 (leading terms only):

$$z + 2z^2 + 5z^3 + 13z^4 + 36z^5 + 102z^6 + 295z^7 + 866z^8 + 2574z^9 + 7730z^{10} + 23419z^{11},$$

$$z^2 + 3z^3 + 9z^4 + 26z^5 + 77z^6 + 230z^7 + 694z^8 + 2110z^9 + 6459z^{10} + 19890z^{11} + 61577z^{12},$$

$$z^3 + 4z^4 + 14z^5 + 45z^6 + 143z^7 + 451z^8 + 1421z^9 + 4478z^{10} + 14129z^{11} + 44654z^{12},$$

$$z^4 + 5z^5 + 20z^6 + 71z^7 + 242z^8 + 806z^9 + 2653z^{10} + 8670z^{11} + 28213z^{12}.$$

One can also substitute $u = 1$, which means that all partial skew Motzkin paths are counted with respect to length, regardless on which level they end:

$$\frac{2 - 3z - 7z^2 - z^3 + z^4 - (2 + z)(1 + z)W}{2z(1 + z)(2z^2 + 3z - 1)}.$$

The series expansion is

$$1 + 2z + 5z^2 + 14z^3 + 40z^4 + 117z^5 + 348z^6 + 1049z^7 + 3196z^8 + 9823z^9 + 30413z^{10} + \cdots$$
Counting flat and left (red) steps. Using two extra variables $t$ and $w$, we can count the number of flat resp. left steps in a skew Motzkin path. The recursions are self-explanatory.

$$
\begin{align*}
f_{j+1} &= zf_j + zg_j + zh_j, \quad f_0 = 1, \\
g_j &= zf_{j+1} + zg_{j+1} + zh_{j+1} + zk_{j+1}, \\
h_j &= ztf_j + ztg_j + zth_j + ztk_j, \\
k_j &= zwg_{j+1} + zwh_{j+1} + zwk_{j+1}.
\end{align*}
$$

Again, here is the system for the multi-variate generating functions;

$$
\begin{align*}
F(u) &= 1 + zuF(u) + zuG(u) + zuH(u), \\
uG(u) &= zF(u) + zG(u) + zH(u) + zK(u) - z - zg_0 - zh_0 - zk_0, \\
H(u) &= ztF(u) + ztG(u) + ztH(u) + ztK(u), \\
uK(u) &= zwG(u) + zwH(u) + zwK(u) - zwg_0 + zwh_0 + zwk_0.
\end{align*}
$$

And following a similar procedure as before we get

$$
\begin{align*}
1 + g_0 + h_0 + k_0 &= \frac{-zw + u_2}{z(1 + wt)} \\
&= 1 + tz + (t^2 + 1)z^2 + (tw + 3t + t^3)z^3 + (2 + 6t^2 + w + 3wt^2 + t^4)z^4 + \dots;
\end{align*}
$$

the quantity $u_2$ is now

$$
u_2 = \frac{1 - tz + wz^2 + twz^3 - \sqrt{(1 - z^2 w)(1 - 2tz + (t^2 - 4 - w)z^2 - 2twz^3 - wt^2 z^4)}}{2z}.
$$

Quantities like $F(u)$, $G(u)$, $H(u)$, $K(u)$ can also be computed easily, following the approach
from the previous section.

Asymptotics for the number of skew Motzkin paths. We must analyze the generating
function

$$
\mathcal{S}M = \frac{(1-z)^2 - \sqrt{(1-z)(1-3z-z^2-z^3)}}{2z^2}
$$

which is of the sqrt-type [16, 20] around the singularity $\rho$ closest to the origin, which we
call $\rho$. It is a solution of $1 - 3z - z^2 - z^3 = 0$ and can be expressed as

$$
\rho = \frac{\sqrt[3]{26 + 6\sqrt{33}}}{3} - \frac{8}{3\sqrt[3]{26 + 6\sqrt{33}}} - \frac{1}{3} \approx 0.295597.
$$

Expanding the generating function around $z = \rho$, we get

$$
\mathcal{S}M \sim \frac{(1-\rho)^2 + 2\sqrt{(1-\rho-\rho^3)(z-\rho)}}{2\rho^2}.
$$

Singularity analysis of generating function [16, 20] gives the estimate

$$
[z^n] \mathcal{S} \mathcal{M} \sim \frac{\sqrt{1 - \rho - \rho^3}}{2\sqrt{\pi \rho^3}} \rho^{-n} n^{-3/2}. \quad (8.1)
$$

The error at $n = 100$ is about 3%. This is to be expected by this type of approximation.
Skew Motzkin paths of bounded height. Now we introduce a parameter $h$ and do not allow the path to reach any level higher than $h$. We can still work with the system

$$
\begin{align*}
f_{j+1} &= zf_j + zg_j + zh_j, \quad 0 \le j \le h-1, \quad f_0 = 1, \\
g_j &= zf_{j+1} + zg_{j+1} + zh_{j+1} + zk_{j+1}, \quad 0 \le j < h, \\
h_j &= zf_j + zg_j + zh_j + zk_j, \quad 0 \le j \le h, \\
k_j &= zg_{j+1} + zh_{j+1} + zk_{j+1}, \quad 0 \le j < h.
\end{align*}
$$

This is now a finite linear system, and we are only interested in paths that return to the $x$-axis. For a given $h$, we write $s[h] = f_0 + g_0 + h_0 + k_0$ for the generating function of path of height $\le h$. It can be proved that both the numerator and the denominator of $s[h]$ satisfy the recursion

$$X_{h+2} + (-1 + z - z^2 - z^3)X_{h+1} + (2z^2 - z^4)X_h = 0.$$

Thus, adjusting this to the initial conditions, we get

$$s[h] = \frac{A_o(1+z^3+z^2-z+\omega)^h + B_o(1+z^3+z^2-z-\omega)^h}{A_u(1+z^3+z^2-z+\omega)^h + B_u(1+z^3+z^2-z-\omega)^h}$$

with

$$
\begin{align*}
\omega &= \sqrt{z^6 + 2z^5 + 3z^4 - 5z^2 - 2z + 1}, \\
A_o &= (z^3 + z^2 + 3z - 1)(z + 1) + (z - 1)\omega, \\
B_o &= (z^3 + z^2 + 3z - 1)(z + 1) - (z - 1)\omega, \\
A_u &= (1 - z^2)(z^3 + z^2 + 3z - 1) + \frac{z^3 - z^2 + 3z - 1}{1 - z}\omega, \\
B_u &= (1 - z^2)(z^3 + z^2 + 3z - 1) - \frac{z^3 - z^2 + 3z - 1}{1 - z}\omega.
\end{align*}
$$

When $h$ goes to infinity, the second terms go away, and we are left with

$$s[\infty] = \frac{A_o}{A_u} = \frac{(1-z)^2 - \sqrt{(1-z)(1-3z-z^2-z^3)}}{2z^2} = \mathcal{S}\mathcal{M},$$

as expected. Now we consider $s[>h]$, the generating function of skew Motzkin paths of height $> h$. Taking differences, we find

$$
\begin{align*}
s[>h] = s[\infty] - s[h] &= \frac{A_o B_u - A_u B_o}{A_u} \frac{(1+z^3+z^2-z-\omega)^h}{A_u(1+z^3+z^2-z+\omega)^h + B_u(1+z^3+z^2-z-\omega)^h} \\
&\sim \frac{A_o B_u - A_u B_o}{A_u^2} \frac{\left(\frac{1+z^3+z^2-z-\omega}{1+z^3+z^2-z+\omega}\right)^h}{1 - \left(\frac{1+z^3+z^2-z-\omega}{1+z^3+z^2-z+\omega}\right)^h}.
\end{align*}
$$

A computer computation leads to (always in the neighbourhood of $z = \rho$)

$$\frac{A_o B_u - A_u B_o}{A_u^2} \sim 18.854986275200314363\sqrt{\rho - z}.$$
Now we approximate and write for convenience:

$$ \frac{1+z^3+z^2-z-\omega}{1+z^3+z^2-z+\omega} \sim 1-5.2213516788791457598\sqrt{\rho-z} \\ \sim \exp(-5.2213516788791457598\sqrt{\rho-z}) = e^{-t}. $$

For the average height, we need apart from the leading factor,

$$ \sum_{h \ge 0} s[>h] \sim \sum_{h \ge 0} \frac{e^{-th}}{1-e^{-th}}. $$

Since we only compute the leading term of the asymptotics of the average height, we might start the sum at $h \ge 1$, and expand the geometric series:

$$ \sum_{h \ge 1} s[>h] \sim \sum_{h,k \ge 1} e^{-thk} = \sum_{k \ge 1} d(k)e^{-kt} \sim -\frac{\log t}{t}, $$

with $d(k)$ being the number of divisors of $k$. This type of analysis, although having been done often before, has been described in much detail in [24]. Together with the factor in front, we are at

$$ -18.854986275200314363 \frac{\log \sqrt{\rho-z}}{5.2213516788791457598\sqrt{\rho-z}} \\ \sim -1.8055656307800996608 \log(1 - z/\rho). $$

Singularity analysis [16] gives the following estimate for the coefficient of $z^n$:

$$ 1.8055656307800996608 \frac{\rho^{-n}}{n}. $$

For the average height we need to normalize, that is, we divide by the total number of skew Motzkin numbers of size $n$ given by (8.1):

$$ \frac{1.8055656307800996608 \frac{\rho^{-n}}{n}}{5.1256244361431546460 \frac{1}{2\sqrt{\pi}} \rho^{-n} n^{-3/2}} = 0.70452513767814089508\sqrt{\pi n}. $$

## 9. OSCILLATIONS IN DYCK PATHS REVISITED

This section in honour of Rainer Kemp was written for this personal survey.

Rainer Kemp's paper [25] was unfortunately largely overlooked. An extension was published quickly [28], and then it fell into oblivion. We want to come back to this gem, with modern methods, in particular, the kernel method and singularity analysis. Kemp was interested in peaks and valleys of Dyck paths, which he called MAX-turns and MIN-turns, probably motivated by computer science applications. The peaks/valleys are enumerated from left to right, and the height of the j-th one is analyzed. In the corresponding ordered (plane) tree, the peaks correspond to the leaves.

Very precise information is available for leaves of binary trees [22, 27, 33, 34] but the situation is a bit different for Dyck paths since the number of peaks/valleys isn't directly related to the length of the Dyck path. (Narayana numbers enumerate them.) Kemp's results in a nutshell are: The average height of the m-th peak/valley is $\sim 4\sqrt{2m/\pi}$ (it is asymptotically independent of the length $n$ of the path), and the difference between the height of the peak and the next valley is about 2, with more terms being available in principle.
**A trivariate generating function for heights of valleys.** The goal is to derive an expression for $\Phi(u) = \Phi(u; z, w)$, where $z$ is used for the length of the path, $w$ for the enumeration of the *valleys* ($w^m$ corresponds to the $m$-th valley), and $u$ is used to record the height of the $m$-th (and last) valley of a partial Dyck path (the path does not need to return to the $x$-axis). We could think about it continued in any possible fashion, as in the following figure. We will figure out the generating function of partial Dyck paths with $m$ valleys, and the generating function of the 'rest', which (if it is not empty) is a partial Dyck path starting with an up-step and ends on the $x$-axis, where the number of valleys is immaterial. The enumeration of the rest is easy, when one thinks about it from right to left, since then it is just a Dyck path ending on a prescribed level $j$ with a down-step. This can be obtained from the first part by setting $w = 1$, i.e., not counting the valleys.

<figure><img src="image_51.png" /><figcaption>FIGURE 10. The third valley at level $j$.</figcaption></figure>

Our goal is, as often, to use the adding-a-new-slice technique, namely adding another 'mountain' (a maximal sequence of up-steps, followed by a maximal sequence of down-steps), or going from the $m$-th valley to the $(m+1)$-st valley. We investigate what can happen to $u^j$:

$$ \sum_{l \ge 1} \sum_{i=1}^{j+l} z^l u^{j+l} z^i u^{-i}. $$

Working this out, the following substitution is essential for our problem:

$$ u^j \rightarrow \frac{z^2 u^k}{(u-z)(1-zu^k)} u^j - \frac{z^{k+2}}{(u-z)(1-z^{k+1})} z^j. $$

Working this into a generating function of the type

$$ \Phi(u) = \sum_{m \ge 0} w^m \varphi_m(u), $$

where the variable $w$ keeps track of the number of mountains, we find from the substitution

$$ \Phi(u) = 1 + \frac{wz^2u}{(u-z)(1-zu)}\Phi(u) - \frac{wz^3}{(u-z)(1-z^2)}\Phi(z), $$

where 1 stands for the empty path having no mountains. Rearranging,

$$ \Phi(u) \frac{z(u-s_1)(u-s_2)}{(u-z)(zu-1)} = 1 - \frac{wz^3}{(u-z)(1-z^2)} \Phi(z), $$

and

$$ \Phi(u) = \frac{(zu-1)}{z(u-s_1)(u-s_2)} \left[ u - z - \frac{wz^3}{(1-z^2)} \Phi(z) \right]. $$
Here,

$$s_2 = \frac{z^2 + 1 - wz^2 - \sqrt{z^4 - 2z^2 - 2z^4w + 1 - 2wz^2 + w^2z^4}}{2z} \quad \text{and} \quad s_1 = \frac{1}{s_2}.$$

In the spirit of the kernel method, the factor $u - s_2$ is 'bad' and must cancel out. That leads first to

$$\Phi(z) = \frac{(1 - z^2)(s_2 - z)}{wz^3}$$

and further to

$$\begin{aligned} \Phi(u) &= \frac{(zu - 1)}{z(u - s_1)} = \frac{s_2(1 - zu)}{z(1 - us_2)} \\ &= 1 + wz^2 + wuz^3 + (w^2 + w + wu^2)z^4 + (2w^2u + wu + wu^3)z^5 + \dots \end{aligned}$$

From this it is easy to read off coefficients in general:

$$[u^j]\Phi(u) = [u^j]\frac{s_2(1 - zu)}{z(1 - us_2)} = \frac{1}{z}s_2^{j+1} - s_2^j.$$

Note that setting $w = 1$ ignores the number of mountains, and the generating function would then be enumerating partial Dyck paths ending on level $j$ with a down-step. The answer could then be derived by combinatorial means as well.

For Kemp's problem, we need

$$S = \sum_{j \ge 0} j \left( \frac{1}{z} s_2^{j+1} - s_2^j \right) \cdot \left( \frac{1}{z} s_2^{j+1} - s_2^j \right) \Big|_{w=1}.$$

Recall that the two parts of the Dyck path, according to our decomposition, are glued together, which just means multiplication of generating functions. The factor $j$ comes in because of the average value of the height of the valley, the first factor is what we just worked out, and the third factor is the rest, which, when read from right to left, is just what we discussed, since the number of valleys or mountains in the rest is irrelevant. Thanks to computer algebra (not available when Kemp worked on the oscillations), we get

$$S = 4 \frac{(-3z + W_1z - W_1 + 1)(-W_2 + wzW_2 + 1 + z^2w^2 - wz^2 - 2wz - z)}{z(-3z - W_1z + 1 - W_1 - wz + wzW_1 - W_2 + W_2W_1)^2}$$

with

$$W_1 = \sqrt{1 - 4z} \quad \text{and} \quad W_2 = \sqrt{z^2 - 2z - 2z^2w + 1 - 2wz + w^2z^2}.$$

Note carefully that $z^2$ was replaced by $z$, since Dyck paths (returning to the x-axis) must have an even number of steps. Their enumeration is classical:

$$D(z) = \frac{1 - \sqrt{1 - 4z}}{2z} \sim 2 - 2\sqrt{1 - 4z},$$

for $z$ close to the (dominant) singularity $z = \frac{1}{4}$. We are in the regime of the subcritical case; see [20, Section IX-3].
The function $S$ has a similar local expansion:

$$S \sim C_1(w) - C_2(w)\sqrt{1 - 4z},$$

and the function $\frac{C_2(w)}{2}$ is the resulting generating function. Working out the details,

$$S \sim \frac{w + \sqrt{(1 - w)(9 - w)} - 3}{-1 + w} - \sqrt{1 - 4z} \left( \frac{w^2 + 2w - 3 + (1 + w)\sqrt{(1 - w)(9 - w)}}{(1 - w)^2} \right) + \dots$$

Eventually we are led to

$$\text{Valley}(w) := \frac{w^2 + 2w - 3 + (1 + w)\sqrt{(1 - w)(9 - w)}}{2(1 - w)^2}.$$

To say it again, the coefficient of $w^m$ in this is the average value of the $m$-th valley in a 'very long' Dyck path. To say more about it, we can use singularity analysis again and expand (this time around $w = 1$, which is dominant):

$$\text{Valley}(w) \sim \frac{2\sqrt{2}}{(1-w)^{3/2}} - \frac{2}{1-w} - \frac{7}{8}\frac{\sqrt{2}}{\sqrt{1-w}}.$$

The traditional translation theorems [16, 20] lead to the average value of the height of the $m$-th valley:

$$4\sqrt{2}\sqrt{\frac{m}{\pi}} - 2 + \frac{5\sqrt{2}}{8\sqrt{\pi m}} + \dots$$

**From valleys to peaks.** We do not need too many new computations, as we can modify the previous results. If one adds an arbitrary non-empty number of up-steps after the $m$-th valley, one has reached the $(m+1)$-st peak! This is basically a substitution!

<figure><img src="image_52.png" /><figcaption>FIGURE 11. The third peak at level $j$.</figcaption></figure>

Start from

$$\Phi(u) = \frac{s_2(1 - zu)}{z(1 - us_2)}$$

and attach a sequence of up-steps: $u^j \to \frac{zu}{1-zu}u^j$. A factor $w$ is also important, since the $m$-th valley corresponds to the $(m+1)$-st peak. Now

$$\frac{zuw}{1-zu} \frac{s_2(1-zu)}{z(1-us_2)} = \frac{us_2w}{1-us_2} = w \sum_{j \ge 1} u^j s_2^j.$$
The computation

$$w \sum_{j \ge 0} j s_2^j \cdot s_2^j \Big|_{w=1}$$

was basically done before, and the local expansion leads to

$$\frac{2w}{1-w} - \frac{2w\sqrt{(1-w)(9-w)}}{(1-w)^2}\sqrt{1-4z},$$

and the generating function of the average values of the m-th peak is

$$\text{Peak}(w) = \frac{w\sqrt{(1-w)(9-w)}}{(1-w)^2}.$$

A local expansion of this results in

$$\text{Peak}(w) \sim \frac{2\sqrt{2}}{(1-w)^{3/2}} - \frac{15}{8} \frac{\sqrt{2}}{\sqrt{1-w}}.$$

Taking differences:

$$\text{Peak}(w) - \text{Valley}(w) \sim \frac{2}{1-w} - \frac{\sqrt{2}}{\sqrt{1-w}},$$

and translating into asymptotics:

$$2 - \frac{\sqrt{2}}{\sqrt{\pi m}}.$$

The formula $2+O(m^{-1/2})$ was already known to Kemp [25]. As Kemp stated in [25], which was confirmed in [28], the generating functions Peak(w) and Valley(w) can be expressed by Legendre polynomials at special values. This is a bit artificial and not too useful in itself.

## 10. DEUTSCH-PATHS IN A STRIP

Emeric Deutsch [11] had the idea to consider a variation of ordinary Dyck paths, by augmenting the usual up-steps and down-steps by one unit each, by down-steps of size 3, 5, 7, ...; the set of down-steps is $\{(1, -1), (1, -3), (1, -5), ...\}$. This leads to ternary equations, as can be seen for instance in [41].

The present author started to investigate a related but simpler model of down-steps 1, 2, 3, 4, ... and investigated it (named Deutsch paths in honour of Emeric Deutsch) in a series of papers [39, 40, 42]. This simpler model can also be seen in the context of Łukasiewicz paths, except that horizontal steps are not allowed; see [5]. Another relevant paper that deals with infinite step sets is [3].

This section is a further member of this series (and extends our unpublished preprint arXiv:2108.12797): The condition that (as with Dyck paths) the paths cannot enter negative territory is relaxed by introducing a negative boundary $-t$. Here are two recent publications about such a negative boundary: [51] and [49].

Instead of allowing negative altitudes, we think about the whole system shifted up by $t$ units, and start at the point $(0, t)$ instead. This is much better for the generating functions that we are going to investigate. Eventually, the results can be re-interpreted as results about enumerations with respect to a negative boundary.
The setting with flexible initial level $t$ and final level $j$ allows us to consider the Deutsch paths also from right to left (they are not symmetric!), without any new computations.

The next sections achieves this, using the celebrated kernel method. An additional upper bound is introduced, so that the Deutsch paths live now in a strip. The way to attack this is linear algebra. Once everything has been computed, one can relax the conditions and let lower/upper boundary go to $\mp\infty$.

**Generating functions and the kernel method.** As discussed, we consider Deutsch paths starting at $(0, t)$ and ending at $(n, j)$, for $n, t, j \ge 0$. First we consider univariate generating functions $f_j(z)$, where $z^n$ stays for $n$ steps done, and $j$ is the final destination. The recursion is immediate:

$$f_j(z) = \llbracket t = j \rrbracket + z f_{j-1}(z) + z \sum_{k>j} f_k(z),$$

where $f_{-1}(z) = 0$. Next, we consider $F(z, u) := \sum_{j \ge 0} f_j(z) u^j$, and get

$$
\begin{align*}
F(z, u) &= u^t + zuF(z, u) + z \sum_{j \ge 0} u^j \sum_{k>j} f_k(z) = u^t + zuF(z, u) + z \sum_{k>0} f_k(z) \sum_{0 \le j < k} u^j \\
&= u^t + zuF(z, u) + z \sum_{k \ge 0} f_k(z) \frac{1-u^k}{1-u} = u^t + zuF(z, u) + \frac{z}{1-u} [F(z, 1) - F(z, u)] \\
&= \frac{u^t(1-u) + zF(z, 1)}{z - zu + zu^2 + 1 - u}.
\end{align*}
$$

Since the critical value is around $u = 1$, we write the denominator as

$$z(u-1)^2 + (u-1)(z-1) + z = z(u-1-r_1)(u-1-r_2),$$

with

$$r_1 = \frac{1 - z + \sqrt{1 - 2z - 3z^2}}{2z}, \quad r_2 = \frac{1 - z - \sqrt{1 - 2z - 3z^2}}{2z}.$$

The factor $(u-1-r_2)$ is bad, so the numerator must vanish for $[u^t(1-u)+zF(z, 1)]|_{u=1+r_2}$, therefore $zF(z, 1) = (1+r_2)^t r_2$. Furthermore

$$F(z, u) = \frac{\frac{u^t(1-u)+zF(z, 1)}{u-r_2}}{z(u-r_1)}.$$

The expressions become prettier using the substitution $z = \frac{v}{1+v+v^2}$; then $r_1 = \frac{1}{v}$, $r_2 = v$. It can be proved by induction (or computer algebra) that

$$\frac{u^t(1-u) + v(1+v)^t}{u-1-v} = -v \sum_{k=0}^{t-1} (1+v)^{t-1-k} - u^t.$$

Furthermore

$$\frac{1}{z(u-1-r_1)} = - \frac{1}{z(1+r_1)(1-\frac{u}{1+r_1})},$$

and so

$$f_j(z) = [u^j] F(z, u) = [u^j] \left[ v \sum_{k=0}^{t-1} (1+v)^{t-1-k} u^k + u^t \right] \sum_{\ell \ge 0} \frac{u^\ell}{z(1+r_1)^{\ell+1}}.$$
Les us mention two interesting special cases: the case $t = 0$ (which was also studied before [39])

$$f_j = \frac{(1 + v + v^2)v^j}{(1 + v)^{j+1}},$$

and the case $j = 0$ for general $t$, as it may be interpreted as Deutsch paths read from right to left, starting at level 0 and ending at level $t \ge 1$ (for $t = 0$, the previous formula applies). It gives

$$f_0(z) = [u^0] \left[ v \sum_{k=0}^{t-1} (1+v)^{t-1-k} u^k + u^t \right] \sum_{\ell \ge 0} \frac{u^\ell}{z(1+r_1)^{\ell+1}} \\ = v(1+v)^{t-1} \frac{1}{z(1+r_1)} = v(1+v+v^2)(1+v)^{t-2}.$$

The next section will present a simplification of the expression for $f_j(z)$, which could be obtained directly by distinguishing cases and summing some geometric series.

**Refined analysis: lower and upper boundary.** Now we consider Deutsch paths bounded from below by zero and bounded from above by $m-1$; they start at level $t$ and end at level $j$ after $n$ steps. For that, we use generating functions $\varphi_j(z)$ (the quantity $t$ is a silent parameter here). The recursions that are straight-forwarded are best organized in a matrix, as the following example shows.

$$\left( \begin{array}{ccccccccc} 1 & -z & -z & -z & -z & -z & -z & -z \\ -z & 1 & -z & -z & -z & -z & -z & -z \\ 0 & -z & 1 & -z & -z & -z & -z & -z \\ 0 & 0 & -z & 1 & -z & -z & -z & -z \\ 0 & 0 & 0 & -z & 1 & -z & -z & -z \\ 0 & 0 & 0 & 0 & -z & 1 & -z & -z \\ 0 & 0 & 0 & 0 & 0 & -z & 1 & -z \\ 0 & 0 & 0 & 0 & 0 & 0 & -z & 1 \end{array} \right) \left( \begin{array}{c} \varphi_0 \\ \varphi_1 \\ \varphi_2 \\ \varphi_3 \\ \varphi_4 \\ \varphi_5 \\ \varphi_6 \\ \varphi_7 \end{array} \right) = \left( \begin{array}{c} 0 \\ 0 \\ 0 \\ 1 \\ 0 \\ 0 \\ 0 \\ 0 \end{array} \right) \left\{ t \right\}$$

The goal is now to solve this system. For that the substitution $z = \frac{v}{1+v+v^2}$ is used throughout. The method is to use Cramer's rule, which means that the right-hand side has to replace various columns of the matrix, and determinants have to be computed. At the end, one has to divide by the determinant of the system.

Let $D_m$ be the determinant of the matrix with $m$ rows and columns. The recursion

$$(1 + v + v^2)^2 m_{n+2} - (1 + v + v^2)(1 + v)^2 D_{m+1} + v(1 + v)^2 D_m = 0$$

appeared already in [39] and is not difficult to derive and to solve:

$$D_m = \frac{(1 + v)^{m-1}}{(1 + v + v^2)^m} \frac{1 - v^{m+2}}{1 - v}.$$
To solve the system with Cramer's rule, we must compute a determinant of shape

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

where the various rows are replaced by the right-hand side. While it is not impossible to solve this recursion by hand, it is very easy to make mistakes, so it is best to employ a computer. Let $D(m; t, j)$ the determinant according to the drawing.

It is not unexpected that the results are different for $j < t$ resp. $j \ge t$. Here is what we found:

$$D(m; t, j) = \frac{(1 + v)^{t-j-3+m}(1 - v^{j+1})v(1 - v^{m-t})}{(1 - v)^2(1 + v + v^2)^{m-1}}, \quad \text{for } j < t,$$

$$D(m; t, j) = \frac{v^{j-t}(1 - v^{t+2})(1 - v^{1-j+m})}{(1 - v)^2(1 + v + v^2)^{m-1}(1 + v)^{j-t+3-m}}, \quad \text{for } j \ge t.$$

To solve the system, we have to divide by the determinant $D_m$, with the result

$$\varphi_j = \frac{D(m; t, j)}{D_m} = \frac{(1 + v)^{t-j-2}(1 - v^{j+1})v(1 - v^{m-t})(1 + v + v^2)}{(1 - v)(1 - v^{m+2})}, \quad \text{for } j < t,$$

$$\varphi_j = \frac{D(m; t, j)}{D_m} = \frac{v^{j-t}(1 - v^{t+2})(1 - v^{1-j+m})(1 + v + v^2)}{(1 - v)(1 + v)^{j-t+2}(1 - v^{m+2})}, \quad \text{for } j \ge t.$$

We found all this using computer algebra. Some critical minds may argue that this is only experimental. One way of rectifying this would be to show that indeed the functions $\varphi_j$ solve the system, which consists of summing various geometric series; again, a computer could be helpful for such an enterprise.

Of interest are also the limits for $m \to \infty$, i.e., no upper boundary:

$$\varphi_j = \lim_{m \to \infty} \frac{D(m; t, j)}{D_m} = \frac{(1 + v)^{t-j-2}(1 - v^{j+1})v(1 + v + v^2)}{(1 - v)}, \quad \text{for } j < t,$$

$$\varphi_j = \frac{v^{j-t}(1 - v^{t+2})(1 + v + v^2)}{(1 - v)(1 + v)^{j-t+2}}, \quad \text{for } j \ge t.$$

The special case $t = 0$ appeared already in the previous section:

$$\varphi_j = \frac{v^j(1 + v + v^2)}{(1 + v)^{j+1}}.$$

Likewise, for $t \ge 1$,

$$\varphi_0 = v(1 + v + v^2)(1 + v)^{t-2}.$$

In particular, the formulae show that the expression from the previous section can be simplified in general, which could have been seen directly, of course.
**Theorem 10.1.** The generating function of Deutsch path with lower boundary 0 and upper boundary $m-1$, starting at $(0, t)$ and ending at $(n, j)$ is given by

$$ \frac{(1+v)^{t-j-2}(1-v^{j+1})v(1-v^{m-t})(1+v+v^2)}{(1-v)(1-v^{m+2})} \quad \text{for } j < t, $$ $$ \frac{v^{j-t}(1-v^{t+2})(1-v^{1-j+m})(1+v+v^2)}{(1-v)(1+v)^{j-t+2}(1-v^{m+2})} \quad \text{for } j \ge t, $$

with the substitution $z = \frac{v}{1+v+v^2}$.

By shifting everything down, we can interpret the results as Deutsch walks between boundaries $-t$ and $m-1-t$, starting at the origin $(0,0)$ and ending at $(n, j-t)$.

**Theorem 10.2.** The generating function of Deutsch path with lower boundary $-t$ and upper boundary $h$, starting at $(0,0)$ and ending at $(n,i)$ with $-t \le i \le h$ is given by

$$ \frac{(1+v)^{i-2}(1-v^{i+t+1})v(1-v^{h+1})(1+v+v^2)}{(1-v)(1-v^{h+t+3})} \quad \text{for } i < 0, $$ $$ \frac{v^i(1-v^{t+2})(1-v^{2-i+h})(1+v+v^2)}{(1-v)(1+v)^{i+2}(1-v^{h+t+3})} \quad \text{for } i \ge 0. $$

It is possible to consider the limits $t \to \infty$ and/or $h \to \infty$ resulting in simplified formulæ.

**Theorem 10.3.** The generating function of Deutsch path with lower boundary $-t$ and upper boundary $\infty$, starting at $(0,0)$ and ending at $(n,i)$ with $-t \le i$ is given by

$$ \frac{(1+v)^{i-2}(1-v^{i+t+1})v(1+v+v^2)}{(1-v)} \quad \text{for } i < 0, $$ $$ \frac{v^i(1-v^{t+2})(1+v+v^2)}{(1-v)(1+v)^{i+2}} \quad \text{for } i \ge 0. $$

**Theorem 10.4.** The generating function of Deutsch path with lower boundary $-\infty$ and upper boundary $h$, starting at $(0,0)$ and ending at $(n,i)$ with $\le i \le h$ is given by

$$ \frac{(1+v)^{i-2}v(1-v^{h+1})(1+v+v^2)}{(1-v)} \quad \text{for } i < 0, $$ $$ \frac{v^i(1-v^{2-i+h})(1+v+v^2)}{(1-v)(1+v)^{i+2}} \quad \text{for } i \ge 0. $$

**Theorem 10.5.** The generating function of unbounded Deutsch path starting at $(0,0)$ and ending at $(n,i)$ is given by

$$ \frac{(1+v)^{i-2}v(1+v+v^2)}{(1-v)} \quad \text{for } i < 0, $$ $$ \frac{v^i(1+v+v^2)}{(1-v)(1+v)^{i+2}} \quad \text{for } i \ge 0. $$

## 11. CONCLUSION

After this personal survey was completed (it will never be complete!) a few more papers about skew Dyck path papers were written; see, e.g., Baril, Kirgizov, Maréchal, Vajnovszki [4] and [47].
REFERENCES

[1] Andrei Asinowski, Benjamin Hackl, and Sarah J. Selkirk. Down-step statistics in generalized Dyck paths. Discrete Math. Theor. Comput. Sci., 24(1):Paper No. 17, 40, 2022.

[2] Alexandru T. Balaban, Roy O. Davies, Frank Harary, Anthony Hill, and Roy Westwick. Cubic identity graphs and planar graphs derived from trees. J. Austral. Math. Soc., 11:207-215, 1970.

[3] Cyril Banderier. Limit laws for basic parameters of lattice paths with unbounded jumps. In Mathematics and Computer Science II (Versailles, 2002), pages 33-47. Birkhäuser, 2002.

[4] Jean-Luc Baril, Sergey Kirgizov, Rémi Maréchal, and Vincent Vajnovszki. Enumeration of Dyck paths with air pockets. J. Integer Seq., 26:Article 23.3, 2023.

[5] Jean-Luc Baril and Helmut Prodinger. Enumeration of partial Łukasiewicz paths. Enumer. Comb. Appl., 3(1):Paper No. S2R2, 13, 2023.

[6] Jean-Luc Baril, José L. Ramírez, and Lina M. Simbaqueba. Counting prefixes of skew Dyck paths. J. Integer Seq., 24(8):Art. 21.8.2, 15, 2021.

[7] Louis Comtet. Advanced Combinatorics. D. Reidel Publishing Co., enlarged edition, 1974.

[8] Julien Courtiel, Andrew Elvey Price, and Irène Marcovici. Bijections between walks inside a triangular domain and Motzkin paths of bounded amplitude. Electron. J. Combin., 28(2):Paper No. 2.6, 2021.

[9] Nicolaas G. de Bruijn, Donald E. Knuth, and Stephen O. Rice. The average height of planted plane trees. In R. C. Read, editor, Graph Theory and Computing, pages 15-22. Academic Press, 1972.

[10] Eva Y. P. Deng and Toufik Mansour. Three Hoppy path problems and ternary paths. Discrete Appl. Math., 156(5):770-779, 2008.

[11] Emeric Deutsch. Problem 10751. Amer. Math. Monthly, 106(7):686, 1999. Solution in AMM vol. 108(9):872, 2001.

[12] Emeric Deutsch, Emanuele Munarini, and Simone Rinaldi. Skew Dyck paths. J. Statist. Plann. Inference, 140(8):2191-2203, 2010.

[13] Emeric Deutsch and Louis W. Shapiro. A bijection between ordered trees and 2-Motzkin paths and its many consequences. Discrete Math., 256(3):655-670, 2002. LaCIM 2000 Conference on Combinatorics, Computer Science and Applications (Montreal, QC).

[14] Maciej Dziemiańczuk. Counting lattice paths with four types of steps. Graphs Combin., 30(6):1427-1452, 2014.

[15] Shalosh B. Ekhad and Doron Zeilberger. Automatic counting of restricted Dyck paths via (numeric and symbolic) dynamic programming. arXiv, 2020.

[16] Philippe Flajolet and Andrew Odlyzko. Singularity analysis of generating functions. SIAM J. Discrete Math., 3(2):216-240, 1990.

[17] Philippe Flajolet and Helmut Prodinger. Register allocation for unary-binary trees. SIAM J. Comput., 15(3):629-640, 1986.

[18] Philippe Flajolet and Helmut Prodinger. Level number sequences for trees. Discrete Math., 65(2):149-156, 1987.

[19] Philippe Flajolet, Jean-Claude Raoult, and Jean Vuillemin. The number of registers required for evaluating arithmetic expressions. Theoret. Comput. Sci., 9(1):99-125, 1979.

[20] Philippe Flajolet and Robert Sedgewick. Analytic Combinatorics. Cambridge University Press, 2009.

[21] Ronald L. Graham, Donald E. Knuth, and Oren Patashnik. Concrete Mathematics. Addison-Wesley, second edition, 1994.

[22] Walter Gutjahr and Georg C. Pflug. The asymptotic distribution of leaf heights in binary trees. Graphs Combin., 8(3):243-251, 1992.

[23] Frank Harary and Ronald C. Read. The enumeration of tree-like polyhexes. Proc. Edinburgh Math. Soc. (2), 17:1-13, 1970.

[24] Clemens Heuberger, Helmut Prodinger, and Stephan Wagner. The height of multiple edge plane trees. Aequationes Math., 90(3):625-645, 2016.

[25] Rainer Kemp. On the average oscillation of a stack. Combinatorica, 2(2):157-176, 1982.

[26] Hana Kim and Richard P. Stanley. A refined enumeration of hex trees and related polynomials. European J. Combin., 54:207-219, 2016.
[27] Peter Kirschenhofer. On the height of leaves in binary trees. *J. Combin. Inform. System Sci.*, 8(1):44–60, 1983.

[28] Peter Kirschenhofer and Helmut Prodinger. On the average hyperoscillations of planted plane trees. *Combinatorica*, 2(2):177–186, 1982.

[29] Christian Krattenthaler. Lattice path enumeration. In *Handbook of Enumerative Combinatorics*, Discrete Math. Appl., pages 589–678. CRC Press, 2015.

[30] Qing Lin Lu. Skew Motzkin paths. *Acta Math. Sin. (Engl. Ser.)*, 33(5):657–667, 2017.

[31] Sri Gopal Mohanty. *Lattice Path Counting and Applications*. Academic Press [Harcourt Brace Jovanovich, Publishers], 1979.

[32] Tadepalli V. Narayana. *Lattice Path Combinatorics with Statistical Applications*. Mathematical Expositions, No. 23. University of Toronto Press, 1979.

[33] Alois Panholzer and Helmut Prodinger. Descendants and ascendants in binary trees. *Discrete Math. Theor. Comput. Sci.*, 1(1):247–266, 1997.

[34] Alois Panholzer and Helmut Prodinger. Moments of level numbers of leaves in binary trees. *J. Statist. Plann. Inference*, 101(1-2):267–279, 2002. Special issue dedicated to the conference Lattice Path Combinatorics and Applications (Vienna, 1998).

[35] Helmut Prodinger. The average height of a stack where three operations are allowed and some related problems. *J. Combin. Inform. System Sci.*, 5(4):287–304, 1980.

[36] Helmut Prodinger. Some analytic techniques for the investigation of the asymptotic behaviour of tree parameters. *EATCS Bulletin*, pages 180–199, 1992.

[37] Helmut Prodinger. The kernel method: a collection of examples. *Sém. Lothar. Combin.*, 50:Art. B50f, 19, 2003/04.

[38] Helmut Prodinger. Analytic methods. In *Handbook of Enumerative Combinatorics*, Discrete Math. Appl., pages 173–252. CRC Press, 2015.

[39] Helmut Prodinger. Deutsch paths and their enumeration. *Open J. Discrete Appl. Math.*, 4(1):12–18, 2021.

[40] Helmut Prodinger. Enumeration of Deutsch paths by the adding-a-new-slice method and applications. *Matimyás Mat.*, 44(1):15–22, 2021.

[41] Helmut Prodinger. Generating functions for a lattice path model introduced by Deutsch. *Spec. Matrices*, 9:217–225, 2021.

[42] Helmut Prodinger. Nondecreasing Deutsch paths. *Fibonacci Quart.*, 59(3):232–236, 2021.

[43] Helmut Prodinger. Retakh's Motzkin paths and some combinatorial comments. *Enumer. Comb. Appl.*, 1(1):Paper No. S2R4, 5, 2021.

[44] Helmut Prodinger. Weighted unary-binary trees, hex-trees, marked ordered trees, and related structures. *arXiv*, 2021.

[45] Helmut Prodinger. The amplitude of Motzkin paths. *Notes Number Theory Discrete Math.*, 28:41–47, 2022.

[46] Helmut Prodinger. Counting ternary trees according to the number of middle edges and factorizing into (3/2)-ary trees. *Asian-Eur. J. Math.*, 15(11):Paper No. 2250187, 5, 2022.

[47] Helmut Prodinger. Partial Dyck paths with air pockets. *Integers*, 22:Paper No. A94, 8, 2022.

[48] Helmut Prodinger. Philippe Flajolet's early work in combinatorics. *Enumer. Comb. Appl.*, 2(1):Paper No. S1H1, 11, 2022.

[49] Helmut Prodinger. On *k*-Dyck paths with a negative boundary. *J. Comb. Math. Comb. Comput.*, 2023. To appear.

[50] Helmut Prodinger. Introduction to Philippe Flajolet's work on the register function and related topics, Chapter introduction, *Collected Works of Philippe Flajolet*, to appear.

[51] Sarah J. Selkirk. On a generalisation of *k*-Dyck paths. Master's thesis, Stellenbosch University, 2019.

[52] Neil J. A. Sloane and The OEIS Foundation Inc. *The On-Line Encyclopedia of Integer Sequences*, 2023.

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

<figure><img src="image_55.png" /></figure>
# A NUMERICAL STUDY OF L-CONVEX POLYOMINOES AND 201-AVOIDING ASCENT SEQUENCES

ANTHONY J. GUTTMANN¹ ID AND VÁCLAV KOTĚŠOVEC² ID

¹SCHOOL OF MATHEMATICS AND STATISTICS, THE UNIVERSITY OF MELBOURNE, VIC. 3010, AUSTRALIA; https://blogs.unimelb.edu.au/tony-guttmann

²PRAGUE, CZECH REPUBLIC; http://www.kotesovec.cz/math.htm

**Abstract.** For L-convex polyominoes we give the conjectured asymptotics of the generating function coefficients, obtained by analysis of the coefficients derived from the functional equation given by Castiglione et al.

For 201-avoiding ascent sequences, we conjecture the solution, obtained from the first twenty-three coefficients of the generating function. This solution is D-finite, indeed algebraic. The conjectured solution then correctly generates all subsequent coefficients. We also obtain the asymptotics, both from direct analysis of the coefficients, and from the conjectured solution.

As well as presenting these new results, our purpose is to illustrate the methods used, so that they may be more widely applied.

**Keywords:** L-convex polyominoes, ascent sequences.

## 1. INTRODUCTION

In [4], Castiglione et al. gave a functional equation for the number of L-convex polyominoes. These are defined as polyominoes with the property that any two cells may be joined by an L-shaped path, that is to say, a path with at most one right-angle bend. An example is shown in Figure 1. It can be seen that such polygons can be described as a stack polyomino placed atop an upside-down stack polyomino. A stack polyomino is just a row-convex bargraph polyomino. The perimeter generating function of L-convex polyominoes has a simple, rational expression,

$$P(x) = \frac{(1-x)^2}{2(1-x)^2 - 1} = 1 + 2x + 7x^2 + 24x^3 + \dots,$$

and is the sequence A003480 in the On-line Encyclopaedia of Integer Sequences (OEIS), [9]. Accordingly, one has

$$[x^n]P(x) = \frac{(2 + \sqrt{2})^{n+1} - (2 - \sqrt{2})^{n+1}}{4\sqrt{2}} \sim \frac{1 + \sqrt{2}}{4}(2 + \sqrt{2})^n.$$

The area generating function is given by Castiglione et al. [4]

$$A(q) = 1 + \sum_{k \ge 0} \frac{q^{k+1} f_k(q)}{(1-q)^2 (1-q^2)^2 \dots (1-q^k)^2 (1-q^{k+1})} = 1 + q + 2q^2 + 6q^3 + 15q^4 + \dots, \qquad (1.1)$$
where

$$f_k(q) = 2f_{k-1}(q) - (1 - q^k)^2 f_{k-2},$$

with initial conditions $f_0(q) = 1$, and $f_1(q) = 1 + 2q - q^2$. We used this expression to generate 2000 terms of the sequence, and these are given in the OEIS as sequence A126764. Analysis of this sequence allowed us to derive the conjectured asymptotics as

$$[q^n]A(q) \sim \frac{13\sqrt{2}}{768 \cdot n^{3/2}} \exp(\pi \sqrt{13n/6}).$$

In the next section we will describe how this estimate was obtained.

<figure><img src="image_56.png" /><figcaption>FIGURE 1. An L-convex polyomino.</figcaption></figure>

The second problem we are considering is that of 201-avoiding ascent sequences, defined below. Given a sequence of non-negative integers, $n_1n_2n_3...n_k$, the number of *ascents* in this sequence is

$$asc(n_1n_2n_3...n_k) = |\{1 \le j < i : n_j \le n_{j+1}\}|.$$

The given sequence is an *ascent sequence* of length $k$ if it satisfies $n_1 = 0$ and

$$n_i \in [0, 1 + asc(n_1n_2n_3...n_{k-1})] \text{ for all } 2 \le i \le k.$$

For example, $(0, 1, 0, 2, 3, 1, 0, 2)$ is an ascent sequence, but $(0, 1, 2, 1, 4, 3)$ is not, as $4 > asc(0121) + 1 = 3$.

Ascent sequences came to prominence when Bousquet-Mélou et al. [2] related them to $(2+2)$-free posets. They have subsequently been linked to other combinatorial structures. See [10] for a number of examples. Later, Duncan and Steingrímsson [6] studied *pattern-avoiding ascent sequences*.
A pattern is simply a word on non-negative integers (repetitions allowed). Given an ascent sequence $(n_1n_2n_3...n_k)$, a pattern $p$ is a subsequence $n_{i_1}n_{i_2}...n_{i_j}$, where $j$ is just the length of $p$, and where the letters appear in the same relative order of size as those in $p$. For example, the ascent sequence $(0, 1, 0, 2, 3, 1)$ has three occurrences of the sequence 001, namely 002, 003 and 001. If an ascent sequence does not contain a given pattern, it is said to be *pattern avoiding*.

The connection between pattern-avoiding ascent sequences and other combinatorial objects, such as set partitions, is the subject of [6], while the connection between pattern-avoiding ascent sequences and a number of stack sorting problems is explored in [5].

Considering patterns of length three, the number of ascent sequences of length $n$ avoiding the patterns 001, 010, 011, and 012 is $2^{n-1}$ (the sequence A000079 in the OEIS). For the pattern 102 the number is $(3^n+1)/2$ (OEIS A007051), while for 101 and 021 the number is just given by the $n^{th}$ Catalan number (OEIS A000108).

For the pattern 201, the first twenty-eight terms of the generating function are given in the OEIS as sequence A202062, and it is this sequence that we have used in our investigation. First, we found, experimentally, that the coefficients given in the OEIS satisfied a recurrence relation, given in Section 3 below. This recurrence can be converted to a second-order inhomogeneous ODE, or, as we prefer, a third-order homogeneous ODE. The smallest root of the polynomial multiplying the third derivative in the ODE is $x = 0.1370633395...$ and is the radius of convergence of the generating function, and of course the reciprocal of the growth constant $\mu = 7.295896946...$.

This ODE, readily converted into differential operator form, can be factored into the direct sum of two differential operators, one of first order and one of second order. The solution of the first order ODE is a rational function while the solution of the second turns out to satisfy a cubic algebraic equation. This can be solved by one's favourite computer algebra package (we give the solution below), and expanding this solution, and adding it to the expansion of the solution of the first-order ODE, gives the required expansion.

This analysis required only the first 24 terms given in OEIS, so the correct prediction of the next four terms gives us confidence that this is indeed the exact solution. Expanding this solution and analysing the coefficients, as described in Section 3, leads to us conjecturing the following asymptotic behaviour for these coefficients:

$$u(n) \sim C \frac{\mu^n}{n^{9/2}},$$

where

$$\mu = \frac{14}{3} \cos \left( \frac{\arccos\left(\frac{13}{14}\right)}{3} \right) + \frac{8}{3}$$

and

$$C = \frac{35}{16} \left( \frac{4107}{\pi} - \frac{84}{\pi} \sqrt{9289} \cos \left( \frac{\pi}{3} + \frac{1}{3} \arccos \left[ \frac{255709\sqrt{9289}}{24653006} \right] \right) \right)^{1/2}.$$

In the next two sections we give the derivation of the results given above.
## 2. L-CONVEX POLYOMINOES

As mentioned above, a typical L-convex polyomino can be considered as a stack polyomino placed atop an upside-down stack polyomino. Stack polyominoes counted by area have generating function

$$S(q) = \sum s_n q^n = \sum_{n \ge 1} \frac{q^n}{(q)_{n-1}(q)_n},$$

where $(q)_n := \prod_{k=1}^n (1 - q^k)$, and, as first shown by Auluck [1], one has

$$s_n \sim \frac{\exp(2\pi\sqrt{n/3})}{8 \cdot 3^{3/4} \cdot n^{5/4}}.$$

Thus putting two such objects together, one would expect a similar expression for the asymptotic form of the coefficients of the generating function (1.1), that is to say, an expression of the form

$$l_n \sim \frac{\exp(a\pi n^\beta)}{cn^\delta}, \qquad (2.1)$$

where we write $L(x) = \sum l_n x^n$ for the ordinary generating function of L-convex polyominoes. We expect both exponents $\beta$ and $\delta$ to be simple rationals, as for stack polyominoes, and the constants $a$ and $c$ to be products of integers and small fractional powers.

The analysis of series with asymptotics of this type is described in detail in [8] and we will not repeat that discussion here, but simply apply the methods described there.

First, we consider the ratios of successive coefficients, $r_n = l_n/l_{n-1}$. For a power-law singularity, one expects the sequence of ratios to approach the growth constant linearly when plotted against $1/n$. In our case the growth constant is 1. That is to say, there is no exponential growth. From the asymptotic behaviour (2.1), which is called of *stretched exponential* type, it follows that the ratio of coefficients behaves as

$$r_n = \frac{l_n}{l_{n-1}} = 1 + \frac{a\beta\pi}{n^{1-\beta}} + O\left(\frac{1}{n}\right), \qquad (2.2)$$

so we expect the ratios to approach a limit of 1 linearly when plotted against $1/n^{1-\beta}$, and to display curvature when plotted against $1/n$. We show the ratios plotted against $1/n$ and $1/\sqrt{n}$ in Figure 2.

<figure><img src="image_57.png" /><figcaption>FIGURE 2. L-convex ratios $r_n$ plotted against $1/n$ (left) and against $1/\sqrt{n}$ (right).</figcaption></figure>
These plots are behaving as expected, with the plot against $1/n$ displaying considerable curvature, while the plot against $1/\sqrt{n}$ is visually linear. This is strong evidence that $\beta = 1/2$, just as is the case for stack polyominoes.

In fact we can easily refine this estimate. From Equation (2.2), one sees that

$$r_n - 1 = a\beta\pi \cdot n^{\beta-1} + O\left(\frac{1}{n}\right).$$

Accordingly, a plot of $\log(r_n - 1)$ versus $\log n$ should be linear, with gradient $\beta - 1$. We would expect an estimate of $\beta$ close to that which linearised the ratio plot. In Figure 3 we show the log-log plot, and in Figure 4 we show the local gradient plotted against $1/\sqrt{n}$. The linearity of the first plot is obvious, while the second is convincingly going to a limit of $-0.5$ as $n \to \infty$.

<figure><img src="image_58.png" /><figcaption>FIGURE 3. Log-log plot of $r_n - 1$ against $n$.</figcaption></figure>

<figure><img src="image_59.png" /><figcaption>FIGURE 4. Gradient of log-log plot.</figcaption></figure>

Having convincingly established that $\beta = 1/2$, just as for stack polyominoes, it remains to determine the other parameters. There are several ways one might proceed, but here is one that works quite well. From the conjectured asymptotic form, we write

$$\lambda_n := \frac{\log(l_n)}{\pi\sqrt{n}} \sim a - \frac{\delta \log n}{\pi\sqrt{n}} - \frac{\log c}{\pi\sqrt{n}},$$

so one can readily fit successive triple of coefficients $\lambda_{k-1}, \lambda_k, \lambda_{k+1}$, to the linear equation

$$\lambda_n = e_1 + e_2 \frac{\log n}{\pi\sqrt{n}} + e_3 \frac{1}{\pi\sqrt{n}},$$

with $k$ increasing until one runs out of known coefficients. Then $e_1$ should give an estimator of $a$, $e_2$ should give an estimator of $-\delta$ and $e_3$ should give an estimator of $-\log(c)$. The result of doing this is shown for $e_1$ and $e_2$ in Figures 5 and 6 respectively.
<figure><img src="image_60.png" /><figcaption>FIGURE 5. Plot of $e_1$ against $1/\sqrt{n}$.</figcaption></figure>

<figure><img src="image_61.png" /><figcaption>FIGURE 6. Plot of $e_2$ against $1/\sqrt{n}$.</figcaption></figure>

We estimate the limits as $n \to \infty$ of $e_1$ as approximately 1.472, and $e_2$ as -1.5. From the asymptotic expression for $s_n$, we expect $a$ to likely involve a square root. So we look at $e_1^2 = 2.16678$, which we conjecture to be 13/6. The exponent $\delta$ is expected to be a simple rational, and 3/2 is indeed a simple rational! We don't show the plot for $e_3$, as it does not give a precise enough estimate to conjecture the value of $\log(c)$ with any precision.

So at this stage we can reasonably conjecture that

$$l_n \sim \frac{\exp(\pi\sqrt{13n/6})}{c \cdot n^{3/2}}. \qquad (2.3)$$

We reached this stage based on only 100 terms in the expansion. In order to both gain more confidence in the conjectured form, and to calculate the constant, we needed more terms, and eventually generated 2000 terms from Equation (1.1).

With hindsight, an arguably more elegant way to analyse this series is to consider only the coefficients $l_{n^2}$. Denote $\ell_n := l_{n^2}$. The conjectured form (2.3) then becomes $\ell_n \sim \frac{\exp(n\pi\sqrt{13/6})}{c \cdot n^3}$. We have 44 coefficients of the series $\ell_n$ available, and these grow in the usual power-law manner, that is, $\ell_n \sim D \cdot \mu^n \cdot n^g$.

We now analyse this sequence assuming its asymptotic form to be

$$\ell_n \sim \frac{\exp(n\pi\sqrt{a})}{c \cdot n^b}, \qquad (2.4)$$

with $a$, $b$, and $c$ to be determined. Then we form the ratios,

$$r_n^{(sq)} = \ell_n / \ell_{n-1} = \mu(1 - b/n + o(1/n)),$$

where $\mu = \exp(\pi\sqrt{a})$ and where the superscript (sq) is a mnemonic recalling that we here consider sequences derived from the subsequence $\ell_n = l_{n^2}$ of square indices. Plotting the ratios $r_n^{(sq)}$ against $1/n$ should give a linear plot with gradient $-b\mu$ and ordinate $\mu$. For a pure power-law the term $o(1/n)$ is $O(1/n^2)$, and the estimate of $\mu$ can thus be refined by plotting the linear intercepts $\ell_n^{(sq)} = n \cdot r_n - (n-1) \cdot r_{n-1}$ against $1/n^2$. The results of doing this are shown in Figures 7 and 8 for the ratios and linear intercepts respectively.
<figure><img src="image_62.png" /><figcaption>FIGURE 7. Plot of ratios $r_n^{(sq)}$ against $1/n$.</figcaption></figure>

<figure><img src="image_63.png" /><figcaption>FIGURE 8. Plot of linear intercepts $\ell_n^{(sq)}$ against $1/n^2$.</figcaption></figure>

It can be seen that the linear intercepts have a faster convergence. We can go further and eliminate the $O(1/n^2)$ term by forming the sequence $t_n = (n^2 \cdot \ell_n^{(sq)} - (n-1)^2 \cdot \ell_{n-1}^{(sq)})/(2n-1)$, and these are shown in Figure 9. From this we estimate that the intercept of the plot with the ordinate is about 101.931. This is the growth constant $\mu = \exp(\pi\sqrt{a})$, from which we find $a \approx 2.16666$, which strongly suggests that $a = 13/6$ exactly.

To estimate the exponent $b$ in the asymptotic form (2.4), we introduce

$$g_n = (r_n^{(sq)} / \mu - 1) \cdot n,$$

noting that $\lim_n g_n = g = -b$. Then, using the estimate of $\mu$ just given, we obtain the plot shown in Figure 10. This is rather convincingly approaching $g = -3$.

<figure><img src="image_64.png" /><figcaption>FIGURE 9. Plot of sequence $t_n$ against $1/n^3$.</figcaption></figure>

<figure><img src="image_65.png" /><figcaption>FIGURE 10. Plot of exponent estimates $g_n$ against $1/n$.</figcaption></figure>
We can do better by calculating the linear intercepts $g2_n := n \cdot g_n - (n-1) \cdot g_{n-1}$. A plot of $g2_n$ against $1/n^2$ is shown in Figure 11. The result $g = -3$ is totally convincing.

<figure><img src="image_66.png" /><figcaption>FIGURE 11. Plot of sequence $g2_n$ against $1/n^3$.</figcaption></figure>

In order to calculate the constant $c$ in the asymptotic form (2.4), we introduce the sequence

$$c_n := \frac{\exp(\pi\sqrt{13n/6})}{l_n \cdot n^{3/2}},$$

and extrapolate the sequence $c_n$ using any of a variety of standard methods.

For this extrapolation, we used the Bulirsch-Stoer method (see [12, Chapter 3.5] or [3] for more details), applied to the coefficient sequence {$l_n$}, with parameter 1/2, and 44 terms in the sequence (corresponding to $44^2 = 1936$ terms in the original series). This gave the estimate $c \approx 0.023938510821419$. This unknown number is likely to involve a square root, cube root or fourth root of a small integer, just as did $s_n$.

We investigated this by dividing by various powers of small integers, and tried to identify the result. Fortuitously, dividing the approximate value by $\sqrt{2}$ gave a result that the Maple command *identify* reported as 13/768. This implies $c = 13\sqrt{2}/768 = 0.023938510821419577...$, agreeing to all quoted digits with the approximate value. The occurrence of 13 in this fraction, as well as in the exponent square-root, is a reassuring feature, as is the factorisation of 768 as $3 \cdot 2^8$.

Thus we conclude with the confident conjecture that the asymptotic form of the coefficients of L-convex polyominoes is

$$l_n \sim \frac{3 \cdot 2^8 \cdot \exp(\pi\sqrt{13n/6})}{13\sqrt{2} \cdot n^{3/2}}.$$
### 3. 201-AVOIDING ASCENT SEQUENCE

From the coefficients $u(n)$ for $n = 0, \dots, 27$ (this is the sequence A202062 in the OEIS), we used the `gfun` package of Maple [11] and immediately found that the coefficients satisfy the recurrence relation

$$ (2n^2 + n)u(n) + (6n^2 + 45n + 60)u(n + 1) + (-34n^2 - 263n - 480)u(n + 2) \\ + (44n^2 + 421n + 984)u(n + 3) + (-20n^2 - 235n - 684)u(n + 4) \\ + (2n^2 + 31n + 120)u(n + 5) = 0, \\ \text{with } u(0) = 1, u(1) = 1, u(2) = 2, u(3) = 5, u(4) = 15. $$

This recurrence can be converted to a second-order inhomogeneous ODE, or, as we prefer, a third-order homogeneous ODE, using the `gfun` command `diffeqtohomdiffeq`, giving

$$ P_3(x)f'''(x) + P_2(x)f''(x) + P_1(x)f'(x) + P_0(x)f(x) = 0, $$

where

$$ \begin{aligned} P_3(x) &= -2x^2(x^3 + 5x^2 - 8x + 1)(4x^4 - 30x^3 + 48x^2 - 36x + 15)(x - 1)^2, \\ P_2(x) &= -3x(x - 1)(12x^8 - 30x^7 - 652x^6 + 2734x^5 - 4767x^4 + 4758x^3 - 2843x^2 + 870x - 85), \\ P_1(x) &= -24x^9 + 30x^8 + 2754x^7 - 13278x^6 + 28884x^5 - 38106x^4 + 32436x^3 - 16620x^2 + 4350x - 420, \\ P_0(x) &= 30(3x - 2)(3x^5 - 10x^4 + 19x^3 - 28x^2 + 24x - 7). \end{aligned} $$

The smallest root of the cubic factor in $P_3(x)$ is $x = 0.1370633395...$ and is the radius of convergence of the solution. Accordingly, the growth constant $\mu$ thus satisfies

$$ \mu = \frac{1}{x} = \frac{14}{3} \cos \left( \frac{\arccos \left( \frac{13}{14} \right)}{3} \right) + \frac{8}{3} = 7.295896946... $$

This ODE can then be studied using the `Maple` package `DEtools`. We first convert the ODE to differential operator form through the command `de2diffop`, then factor this into the direct sum of two differential operators by the command `DFactorLCLM`. One of these operators is first order and one is second order.

The solution of the first order ODE is immediately given by the `dsolve` command, and is the rational function

$$ y_1(x) = \frac{x^4 + 26x^3 - 45x^2 + 18x + 1}{12(x - 1)x^3}. $$

To solve the second-order ODE, we obtain a series solution, the first term of which is $O(x^{-3})$. We multiply the solution by $x^3$ to obtain a regular power series, then use the `gfun` command `seriestoalgeq` to discover the cubic equation,

$$ \begin{aligned} & 4(x - 1)^3 y_2(x)^3 \\ & - 3(x - 1)(x^2 - x + 1)(x^6 - 235x^5 + 1430x^4 - 1695x^3 + 270x^2 + 229x + 1)y_2(x) \\ & + x^{12} + 510x^{11} - 14631x^{10} + 80090x^9 - 218058x^8 + 316290x^7 - 253239x^6 \\ & + 131562x^5 - 70998x^4 + 37950x^3 - 8955x^2 - 522x + 1 = 0. \end{aligned} \quad (3.1) $$
This can be solved by Maple's, solve command, giving three solutions. Inspection of their expansion reveals the appropriate one, and simplifying this gives the following rather cumbersome solution: Let

$$P_1 = x^{12} + 510x^{11} - 14631x^{10} + 80090x^9 - 218058x^8 + 316290x^7 - 253239x^6 + 131562x^5 - 70998x^4 + 37950x^3 - 8955x^2 - 522x + 1 - 24\sqrt{3x(x-1)(x^3 + 5x^2 - 8x + 1)^7},$$

$$P_2 = (x^2 - x + 1)(x - 1)^4(x^6 - 235x^5 + 1430x^4 - 1695x^3 + 270x^2 + 229x + 1),$$

$$P_3 = (3^{5/6}i + 3^{1/3}), \text{ and } P_4 = (-3^{5/6}i + 3^{1/3}). \text{ Then, one has}$$

$$y_2(x) = \frac{-3^{2/3} (P_4 (-P_1 \cdot (x-1)^6)^{2/3} + P_2 \cdot P_3)}{12 (-P_1 \cdot (x-1)^6)^{1/3} (x-1)^3}.$$

The solution to the original ODE is then

$$y(x) = \frac{y_2(x)}{12x^3} - y_1(x) = 1 + x + 2x^2 + 5x^3 + 15x^4 + \dots$$

This analysis required only the first 24 terms given in the OEIS, so the correct prediction of the next four terms gives us confidence that this is indeed the exact solution.

We next obtained the first 5000 terms in only a few minutes of computer time by expanding this solution. We used these terms to calculate the amplitude. That is to say, we now know that the coefficients behave asymptotically as $u(n) \sim C\mu^n n^{-9/2}$. Equivalently, the generating function behaves as

$$U(x) = \sum u(n)x^n = A(1 - \mu \cdot x)^{7/2},$$

where $C = A/\Gamma(-7/2) = 105A/(16\sqrt{\pi})$. We estimate $C$ by assuming a pure power law, so that

$$\frac{u(n) \cdot n^{9/2}}{\mu^n} = C(1 + \sum_{k \ge 1} a_k/n^k).$$

We calculated the first twenty coefficients of this expansion, which allowed us to estimate $C = 13.4299960869...$ with 74-digit accuracy (as checked later). Unless one is very fortunate (for example, when the Maple command *identify* determines an expression for this constant, which it doesn't in our case), to identify this constant requires some experience-based guesswork.

Such constants in favourable cases are a product of rational numbers and square roots of small integers, sometimes with integer or half-integer powers of $\pi$. These powers of $\pi$ usually arise from the conversion factor in going from the generating function amplitude $A$ to the coefficient amplitude $C$. That is to say, we might expect the amplitude $A$ to be simpler than $C$. And, to eliminate square-roots, we will try and identify $A^2$ rather than $A$.

We do this by seeking the minimal polynomial with root $A^2$, using the command *MinimalPolynomial* in either *Maple* or *Mathematica*. In fact, one only requires 20 digit accuracy in the estimate of $A^2$ to establish the minimal polynomial, $A^6 - 1369A^4 + 17839A^2 + 1$, which can be solved to give

$$C = \frac{35}{16} \left( \frac{4107}{\pi} - \frac{84}{\pi} \sqrt{9289} \cos \left( \frac{\pi}{3} + \frac{1}{3} \arccos \left[ \frac{255709\sqrt{9289}}{24653006} \right] \right) \right)^{1/2}.$$
This derivation includes a degree of hindsight. In fact we searched for the minimal polynomial for the amplitude $C$, by including various powers of $\pi$, and then choose the polynomial of minimal degree. This required a much greater degree of precision in our estimate of $C$ to ensure we found the correct minimal polynomial.

It has been pointed out to us by Jean-Marie Maillard that the amplitude $A$ can be obtained directly from the solution of the cubic equation (3.1), by extracting the coefficient of $(1 - \mu \cdot x)^{7/2}$, as explained in [7, Chapter VII.7.1]. This gives the minimal polynomial that we obtained by numerical experimentation. This alternative way to derive asymptotic expansions is a more elegant method, as it is automatic, but it only works for algebraic functions. There are thus many sequences for which it is not applicable, as in the case of $L$-convex polyominoes (for which the generating function is not algebraic, as its radius of convergence is not algebraic), while our numerical approach can still yield conjecturally exact results.

## 4. CONCLUSION

We have shown how experimental mathematics can be used to conjecture exact asymptotics, in the case of $L$-convex polyominoes, and to conjecture an exact solution, in the case of 201-avoiding ascent sequences. We hope that the results will be of interest, and that the methods will be more widely applied, as there are many outstanding combinatorial problems that lend themselves to such an approach.

We recognise that these results are conjectural. We leave proofs to those more capable, and in the hope that the maxim of the late lamented J. M. Hammersley to the effect that "it is much easier to prove something when you know that it is true" will aid that endeavour.

**Acknowledgements.** AJG wishes to acknowledge helpful discussions with Jean-Marie Maillard, and with Paolo Massazza on the topic of $L$-convex polyominoes, and to thank the ARC Centre of Excellence for Mathematical and Statistical Frontiers (ACEMS) for support.

## REFERENCES

[1] Faqir Chand Auluck. On some new types of partitions associated with generalized Ferrers graphs. *Proc. Cambridge Philos. Soc.*, 47:679–686, 1951.

[2] Mireille Bousquet-Mélou, Anders Claesson, Mark Dukes, and Sergey Kitaev. (2+2)-free posets, ascent sequences and pattern avoiding permutations. *J. Comb. Theory, Ser. A*, 117(7):884–909, 2010.

[3] Roland Bulirsch and Josef Stoer. Fehlerabschätzungen und Extrapolation mit rationalen Funktionen bei Verfahren vom Richardson-Typus. *Numer. Math.*, 6:413–427, 1964.

[4] Giusi Castiglione, Andrea Frosini, Emanuele Munarini, Antonio Restivo, and Simone Rinaldi. Combinatorial aspects of $L$-convex polyominoes. *Eur. J. Comb.*, 28(6):1724–1741, 2007.

[5] Giulio Cerbai, Anders Claesson, and Luca Ferrari. Stack sorting with restricted stacks. *J. Comb. Theory, Ser. A*, 173:19 pp., 2020.

[6] Paul Duncan and Einar Steingrímsson. Pattern avoidance in ascent sequences. *Electron. J. Comb.*, 18(1):17 pp., 2011.

[7] Philippe Flajolet and Robert Sedgewick. *Analytic Combinatorics*. Cambridge University Press, 2009.

[8] Anthony J. Guttmann. Analysis of series expansions for non-algebraic singularities. *J. Phys. A*, 48(4):045209, 33, 2015.

[9] OEIS Foundation Inc. The On-Line Encyclopedia of Integer Sequences. https://oeis.org, 2022.

[10] Sergey Kitaev. *Patterns in Permutations and Words*. Springer, 2011.

[11] 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.

[12] Josef Stoer and Roland Bulirsch. *Introduction to Numerical Analysis*. Springer, 2002.
# INTRODUCING DASEP:
THE DOUBLY ASYMMETRIC SIMPLE EXCLUSION PROCESS

DAVID W. ASH¹

¹OPTIN INC., SAN FRANCISCO BAY AREA, CA, USA; https://orcid.org/0000-0001-8970-5784

**Abstract.** Research in combinatorics has often explored the asymmetric simple exclusion process (ASEP). The ASEP, inspired by examples from statistical mechanics, involves particles of various species moving around a lattice. With the traditional ASEP particles of a given species can move but do not change species. In this paper a new combinatorial formalism, the DASEP (doubly asymmetric simple exclusion process), is explored. The DASEP is inspired by biological processes where, unlike the ASEP, the particles can change from one species to another. The combinatorics of the DASEP on a one dimensional lattice are explored, including the associated generating function. The stationary probabilities of the DASEP are explored, and results are proven relating these stationary probabilities to those of the simpler ASEP.

**Keywords:** ASEP, DASEP, lattice, algebraic combinatorics, steady state probabilities, species, lattice paths.

## 1. INTRODUCTION

The ASEP (asymmetric simple exclusion process) is a structure that has frequently been referred to in the combinatorics literature. In its simplest form, the ASEP consists of a one dimensional infinite lattice, with each point on the lattice being populated with either a particle or a hole. At random intervals, each particle attempts to move either to the left or the right with different but fixed probabilities (hence the term ‘asymmetric’). The ASEP can be thought of as a form of Markov process as noted in [4] by Corteel et al. Multiline queues [5] were introduced by Ferrari et al. as a combinatorial approach to the analysis of the ASEP. Originally the ASEP particles were thought of as all belonging to a single species. More recent work by Cantini et al. [3] generalized the concept to multiple species and uncovered a link with Macdonald polynomials. Although we focus on the homogeneous ASEP (transition probabilities do not depend on position in the lattice), several researchers (Lam et al. [7], Ayyer et al. [1], Cantini [2], Mandelshtam [9], and Kim et al. [6]) have explored the inhomogeneous ASEP in which transition probabilities do depend on lattice position.

## 2. DEFINITIONS

Following [4], a partition can be defined as follows:

**Definition 2.1.** A partition $\lambda$ is a nonincreasing sequence of $n$ nonnegative integers $\lambda = (\lambda_1 \ge \lambda_2 \ge \dots \ge \lambda_n \ge 0)$.
We will start by working through a simple example of the ASEP before introducing the new concept of the DASEP. We will ordinarily write a partition as defined above as an n-tuple: $\lambda = (\lambda_1, \lambda_2, ..., \lambda_n)$.

**Definition 2.2.** We write $S_n(\lambda)$ to mean the set of all permutations of $\lambda$.

**Example 2.3.** So, for $\lambda = (2, 2, 1)$, $S_3(\lambda) = \{(2, 2, 1), (2, 1, 2), (1, 2, 2)\}$.

The multispecies asymmetric simple exclusion process ASEP($\lambda$) is then defined to be a Markov process on $S_n(\lambda)$ with certain specific transition probabilities:

**Definition 2.4.** For all partitions $\lambda$ as defined in Definition 2.1, ASEP($\lambda$) is a Markov process on $S_n(\lambda)$. We let $t$ be a constant with $0 \le t \le 1$. The transition probability $P_{\mu,\nu}$ between two permutations $\mu \in S_n(\lambda)$ and $\nu \in S_n(\lambda)$ is given by:

* If $\mu = (\mu_1, ..., \mu_k, i, j, \mu_{k+2}, ..., \mu_n)$ and $\nu = (\mu_1, ..., \mu_k, j, i, \mu_{k+2}, ..., \mu_n)$, with $i \neq j$, then $P_{\mu,\nu} = \frac{t}{n}$ if $i > j$ and $P_{\mu,\nu} = \frac{1}{n}$ if $j > i$.

* If $\mu = (i, \mu_2, \mu_3, ..., \mu_{n-1}, j)$ and $\nu = (j, \mu_2, \mu_3, ..., \mu_{n-1}, i)$ with $i \neq j$, then $P_{\mu,\nu} = \frac{t}{n}$ if $j > i$ and $P_{\mu,\nu} = \frac{1}{n}$ if $i > j$.

* If neither of the above conditions apply but $\nu \neq \mu$ then $P_{\mu,\nu} = 0$. If $\nu = \mu$ then $P_{\mu,\mu} = 1 - \sum_{\nu \neq \mu} P_{\mu,\nu}$.

It is possible to compute steady state probabilities for ASEP($\lambda$). For the purposes of the example that we will develop as we introduce DASEP, we are primarily interested in ASEP($\lambda$) for $\lambda = (2, 2, 0)$, $\lambda = (2, 1, 0)$, and $\lambda = (1, 1, 0)$, so we will focus mostly on these three processes as we work through the computation of the steady state probabilities. Continuing to follow [4] as we develop this example, to compute these probabilities we need to define the concept of a multiline queue.

**Definition 2.5.** A ball system $B$ is an $L \times n$ matrix each element of which is either 0 or 1. Moreover for all $i$ the number of 1's in row $i+1$ is less than or equal to the number of 1's in row $i$.

**Definition 2.6.** Given a ball system $B$ a multiline queue $Q$ is obtained by augmenting $B$ with a labeling and matching system. Each cell in $B$ will be labelled with a number from 0 to $L$ inclusive, and each cell with a 1 element in row $i+1$, for $i \ge 1$, will be matched to a cell with a 1 element in row $i$. Such a matching must be obtained through an application of the following algorithm:

* Step 1: Find the highest numbered row with unlabelled 1 elements. Label each of those elements with the number of the row. If this is row 1, or there are no remaining unlabelled 1 elements in the matrix, exit.

* Step 2: Find the row with labelled but unmatched elements. If this is row 1, go back to step 1. If it is row $i+1$, for $i \ge 1$, first match each labelled but unmatched element that can be matched to an unlabelled element directly below it to that element. This is considered a *trivial* match. Then proceed from right to left (highest to lowest numbered columns) matching each remaining labelled but unmatched element to an unlabelled element in the row below–these are the *nontrivial* matches. Give all newly matched elements in row $i$ the same label as the element it has just been matched to. Repeat step 2.
A multiline queue is often visualized as a ball system with an element with a 1 value being shown as a ball and a 0 value by the absence of a ball. Matches between elements (balls) are drawn by lines between the matched balls. The following shows a multiline queue associated with ASEP($\lambda$) where $\lambda = (2, 2, 0)$. Note that the line matching the ball at upper right to the one at the lower middle wraps around to the right.

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

The labels in the bottom row determine the partition of the associated ASEP. The above multiline queue has $\lambda = (2, 2, 0)$ since the bottom row includes two 2's and a 0-by convention an element without a ball is assumed to be labeled with a 0. Likewise the following would be a multiline queue with $\lambda = (2, 1, 0)$:

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

Each multiline queue is also associated with a permutation $\alpha \in S_n(\lambda)$ corresponding to the labels of its bottom row in unsorted order. For example, for the above multiline queue, $\lambda = (2, 1, 0)$ but $\alpha = (0, 1, 2)$. We will write $\lambda(Q) = \lambda$ and $\alpha(Q) = \alpha$.

### 3. STEADY STATE PROBABILITIES WITH EXAMPLE

To determine steady state probabilities-and continue with the example started in the introduction-we next assign to each nontrivial matching $p$ in $Q$ two values $f(p)$ and $s(p)$. $f(p)$ is the number of choices that were available for the match when the match was made. $s(p)$ is the number of legal matches that were skipped, if we imagine ourselves considering possible matches from left to right and wrapping around the end if needed, before the actual choice was made. We can then define a weight on $p$ as $wt(p) = \frac{(1-t)t^{s(p)}}{1-tf(p)}$. Here we are proceeding from [4] but with the simplifying assumption that $q = 1$, since in the sequel we will rely on a theorem that requires $q = 1$. Next we can define a weight on the entire multiline queue $wt(Q) = \prod_{p \in Q} wt(p)$ where the product is taken over all nontrivial matches $p$ in $Q$. A theorem due to Martin [10] then gives the required steady state probabilities:

$$ \mathrm{Pr}(\alpha) = \frac{\sum_{\alpha(Q)=\alpha} \mathrm{wt}(Q)}{\sum_{\lambda(Q)=\lambda} \mathrm{wt}(Q)}. $$

Before moving on to the DASEP, we need to evaluate the steady state probabilities for the examples that we will ultimately use to develop the DASEP. For the above multiline queue, there is exactly one nontrivial pair $p$. When this pair is matched, there are two available options so $f(p) = 2$. As we picked the second available option, $s(p) = 1$. So $wt(Q) = \frac{(1-t)t}{1-t^2}$. As noted above, $\alpha = (0, 1, 2)$ and the only other multiline queue with $\alpha = (0, 1, 2)$ is as follows:

<figure><img src="image_69.png" /></figure>
Here there is no nontrivial matching pair, so $wt(Q) = 1$. Hence:

$$ \sum_{\alpha(Q)=(0,1,2)} \mathrm{wt}(Q) = 1 + \frac{(1-t)t}{1-t^2} = \frac{1+2t}{1+t}. $$

For reasons of symmetry:

$$ \sum_{\alpha(Q)=(0,1,2)} \mathrm{wt}(Q) = \sum_{\alpha(Q)=(1,2,0)} \mathrm{wt}(Q) = \sum_{\alpha(Q)=(2,0,1)} \mathrm{wt}(Q) = \frac{1+2t}{1+t}. $$

Next we look at $\alpha = (0, 1, 2)$, for which there are also two multiline queues. The first of these is as follows:

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

Here there is one nontrivial matching pair $p$. When this pair is matched, there are two available options so $f(p) = 2$. As we picked the first available option, $s(p) = 0$. So $wt(Q) = \frac{1-t}{1-t^2}$. The other multiline queue with $\alpha = (2, 1, 0)$ is as follows:

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

Again there is no nontrivial matching pair, so $wt(Q) = 1$. Hence:

$$ \sum_{\alpha(Q)=(2,1,0)} \mathrm{wt}(Q) = 1 + \frac{1-t}{1-t^2} = \frac{2+t}{1+t}. $$

For reasons of symmetry, one has

$$ \sum_{\alpha(Q)=(2,1,0)} \mathrm{wt}(Q) = \sum_{\alpha(Q)=(1,0,2)} \mathrm{wt}(Q) = \sum_{\alpha(Q)=(0,2,1)} \mathrm{wt}(Q) = \frac{2+t}{1+t}. $$

So we get

$$ \sum_{\lambda(Q)=(2,1,0)} \mathrm{wt}(Q) = 3(\frac{1+2t}{1+t}) + 3(\frac{2+t}{1+t}) = 9. $$

We are now ready to give the steady state probabilities

$$ \mathrm{Pr}(0, 1, 2) = \mathrm{Pr}(1, 2, 0) = \mathrm{Pr}(2, 0, 1) = \frac{1+2t}{9(1+t)} $$

and

$$ \mathrm{Pr}(2, 1, 0) = \mathrm{Pr}(1, 0, 2) = \mathrm{Pr}(0, 2, 1) = \frac{2+t}{9(1+t)}. $$

Trivial computations also give

$$ \mathrm{Pr}(0, 1, 1) = \mathrm{Pr}(1, 1, 0) = \mathrm{Pr}(1, 0, 1) = \frac{1}{3} $$

and

$$ \mathrm{Pr}(0, 2, 2) = \mathrm{Pr}(2, 2, 0) = \mathrm{Pr}(2, 0, 2) = \frac{1}{3}. $$

This concludes our computation for the steady state probabilities of this model; in the next section we introduce the DASEP model.
<figure><img src="image_72.png" /><figcaption>FIGURE 1. An example of DASEP(n, p, q): DASEP(3,2,2).</figcaption></figure>

Each of the 12 small triangles represents a state of the DASEP in the circular lattice with $n = 3$ sites with $q = 2$ balls (i.e., the nonzero labels), each nonzero label is $\le p = 2$ (i.e., one has $p = 2$ species). Each state corresponds to a permutation of the partition $(2, 2, 0)$, $(2, 1, 0)$, or $(1, 1, 0)$. The transitions are explained in Definition 4.1 hereafter.
## 4. DOUBLY ASYMMETRIC SIMPLE EXCLUSION PROCESS

We are now ready to introduce the DASEP (doubly asymmetric simple exclusion process). While the ASEP is inspired by statistical mechanics where particles do not change species, the DASEP, by contrast, is inspired by biological processes where particles can change species, which we denote by DASEP(n, p, q) where n is the number of positions on the lattice, p is the number of types of species, and q is the number of particles.

**Definition 4.1.** For all positive integers $n, p$, and $q$ with $n > q$, DASEP$(n, p, q)$ is a Markov process on the set $\bigcup_{\lambda_1 \le p, \lambda'_1 = q} S_n(\lambda)$, where one uses the notation of Definition 2.2, and where $\lambda'_1 = q$ refers to the dual partition [8] of $\lambda$, namely $\lambda'$, and uses the fact that $\lambda'_1$ gives the number of nonzero terms in the original partition $\lambda$. The transition probability $P_{\mu,\nu}$ on two permutations $\mu$ and $\nu$ is as follows:

* If $\mu = (\mu_1, \dots, \mu_k, i, j, \mu_{k+2}, \dots, \mu_n)$ and $\nu = (\mu_1, \mu_2, \dots, \mu_k, j, i, \mu_{k+2}, \dots, \mu_n)$ with $i \neq j$, then $P_{\mu,\nu} = \frac{t}{3n}$ if $i > j$ and $P_{\mu,\nu} = \frac{t}{3n}$ if $j > i$.

* If $\mu = (i, \mu_2, \dots, \mu_{n-1}, j)$ and $\nu = (j, \mu_2, \mu_3, \dots, \mu_{n-1}, i)$ with $i \neq j$, then $P_{\mu,\nu} = \frac{t}{3n}$ if $j > i$ and $P_{\mu,\nu} = \frac{t}{3n}$ if $i > j$.

* If $\mu = (\mu_1, \dots, \mu_k, i, \mu_{k+2}, \dots, \mu_n)$ and $\nu = (\mu_1, \dots, \mu_k, i+1, \mu_{k+2}, \dots, \mu_n)$ with $i \ge 1$, then $P_{\mu,\nu} = \frac{u}{3n}$.

* If $\mu = (\mu_1, \dots, \mu_k, i+1, \mu_{k+2}, \dots, \mu_n)$ and $\nu = (\mu_1, \dots, \mu_k, i, \mu_{k+2}, \dots, \mu_n)$ with $i \ge 1$, then $P_{\mu,\nu} = \frac{1}{3n}$.

* If none of the above conditions apply but $\nu \neq \mu$ then $P_{\mu,\nu} = 0$. If $\nu = \mu$ then $P_{\mu,\mu} = 1 - \sum_{\nu \neq \mu} P_{\mu,\nu}$.

Figure 1 shows the simple example of the DASEP that we are working through. All possible transitions within a single ASEP (the first and second bullet points in the definition above) are shown with blue arrows on this diagram. To keep the diagram relatively clean in appearance, only selected transitions between different ASEPs (the third and fourth bullet points) are shown (with red arrows). Other ASEPs such as ASEP(1, 0, 0) or ASEP(2, 0, 0) are not shown since these are not part of DASEP(3, 2, 2). This is because, per Definition 4.1, for DASEP(3, 2, 2) we always have $\lambda'_1 = 2$, whereas for ASEP(1, 0, 0) and ASEP(2, 0, 0), we would have $\lambda'_1 = 1$.

Similar to with the ASEP, with the DASEP we wish to compute steady state probabilities for permutations $\alpha$ which we will call $Pd(\alpha)$. We will focus on continuing to develop the example we have been working on which turns out to be DASEP(3, 2, 2). Here $n = 3$ means that the particles move on the circular lattice with 3 sites, $p = 2$ means that each particle is allowed to take on the value 0, 1, or 2, and $q = 2$ means that each permutation $\alpha$ has exactly 2 nonzero values. We therefore find ourselves interested in the following 12 steady state probabilities:

$Pd(0, 1, 1), Pd(0, 1, 2), Pd(0, 2, 1), Pd(0, 2, 2), Pd(1, 0, 1), Pd(1, 0, 2), Pd(1, 1, 0), Pd(1, 2, 0), Pd(2, 0, 1), Pd(2, 0, 2), Pd(2, 1, 0), Pd(2, 2, 0).$
Note here that particles in the DASEP are allowed to switch back and forth between species 1 and 2, but not back and forth from 0 to anything else. That is because a value of 0 is understood to not so much be a species but the absence of a species. Due to symmetries we can now focus on solving for the following four probabilities:

$$w = \mathrm{Pd}(0, 1, 1), \quad x = \mathrm{Pd}(0, 1, 2), \quad y = \mathrm{Pd}(0, 2, 1), \quad z = \mathrm{Pd}(0, 2, 2).$$

From the above transition probabilities, this reduces to solving the system

$$\begin{cases} 2uw = x + y \\ (2+t)x + x + ux = (1+2t)y + z + uw \\ (1+2t)y + y + uy = (2+t)x + uw + z \\ 2z = u(x + y) \end{cases}$$

which in turn implies the relation

$$(5 + 2t + u)x = (3 + 4t + u)y.$$

We can then ask ourselves the question of when the proportions of steady state probabilities for the DASEP are the same as for the previous ASEP. Noting that $\mathrm{Pr}(0, 1, 2) = \frac{1+2t}{9(1+t)}$ and $\mathrm{Pr}(2, 1, 0) = \frac{2+t}{9(1+t)}$ such equality will happen if $(5+2t+u)(1+2t) = (3+4t+u)(2+t)$, or $5+2t+u+10t+4t^2+2tu = 6+8t+2u+3t+4t^2+tu$, or $t(1+u) = 1+u$. So this will happen if and only if $t=1$. We have therefore proven the following proposition.

**Proposition 4.2.** If $D = \mathrm{DASEP}(3, 2, 2)$ is parameterized as described above by $t$ and $u$, then the following two statements are equivalent:

* $t=1$.

* For all partitions $\lambda$ with $S_n(\lambda) \subseteq D$ and all permutations $\mu, \nu \in S_n(\lambda)$, the following equality holds: $\frac{\mathrm{Pr}(\mu)}{\mathrm{Pr}(\nu)} = \frac{\mathrm{Pd}(\mu)}{\mathrm{Pd}(\nu)}$. That is, the ratio between steady state probabilities does not change in moving from the ASEP to the DASEP.

In fact, we conjecture the following more general statement.

**Conjecture 4.3.** If $D = \mathrm{DASEP}(n, p, q)$ is parameterized as described above by $t$ and $u$, then the following two statements are equivalent:

* $t=1$.

* For all partitions $\lambda$ with $S_n(\lambda) \subseteq D$ and all permutations $\mu, \nu \in S_n(\lambda)$, the following equality holds: $\frac{\mathrm{Pr}(\mu)}{\mathrm{Pr}(\nu)} = \frac{\mathrm{Pd}(\mu)}{\mathrm{Pd}(\nu)}$. That is, the ratio between steady state probabilities does not change in moving from the ASEP to the DASEP.

*Partial proof.* We will prove this only in the $\implies$ direction. If $t=1$ we can replace $\lambda$ with a similar partition but with species of the same type being replaced by similar distinct species. For example, if $\lambda = (3, 3, 3, 2, 1, 0, ...)$ we would map this to $\hat{\lambda} = (3_1, 3_2, 3_3, 2, 1, 0, ...)$ and allow adjacent species originally of the same type to be exchanged with the same transition probability. This will create a completely symmetric situation, so all steady state probabilities are equal. As an equal number of $\hat{\lambda}$'s can be derived from each $\lambda$ this means all original steady state probabilities are equal as well, so $\frac{\mathrm{Pr}(\mu)}{\mathrm{Pr}(\nu)} = \frac{\mathrm{Pd}(\mu)}{\mathrm{Pd}(\nu)} = 1$. This completes the proof in the $\implies$ direction. $\square$
Let us motivate Conjecture 4.3 by showing that it holds on one example. Following are the nine values (of the nine steady state probabilities) we must solve for to prove this conjecture for DASEP(3, 3, 2):

$$
\begin{align*}
a_1 &= \mathrm{Pd}(0, 1, 1), & a_2 &= \mathrm{Pd}(0, 2, 2), & a_3 &= \mathrm{Pd}(0, 3, 3), \\
b_1 &= \mathrm{Pd}(0, 2, 3), & b_2 &= \mathrm{Pd}(0, 1, 3), & b_3 &= \mathrm{Pd}(0, 1, 2), \\
c_1 &= \mathrm{Pd}(0, 3, 2), & c_2 &= \mathrm{Pd}(0, 3, 1), & c_3 &= \mathrm{Pd}(0, 2, 1).
\end{align*}
$$

These values can be obtained by solving the following set of nine equations:

$$
\left\{
\begin{array}{l}
2ua_1 = b_3 + c_3 \\
(2+t)b_3 + ub_3 + ub_3 + b_3 = (1+2t)c_3 + b_2 + a_2 + ua_1 \\
(1+2t)c_3 + uc_3 + uc_3 + c_3 = (2+t)b_3 + c_2 + a_2 + ua_1 \\
a_2 + a_2 + ua_2 + ua_2 = b_1 + c_1 + ub_3 + uc_3 \\
(2+t)b_2 + ub_2 + b_2 = (1+2t)c_2 + b_1 + ub_3 \\
(1+2t)c_2 + uc_2 + c_2 = (2+t)b_2 + c_1 + uc_3 \\
2a_3 = ub_1 + uc_1 \\
(2+t)b_1 + ub_1 + b_1 + b_1 = (1+2t)c_1 + a_3 + ub_2 + ua_2 \\
(1+2t)c_1 + uc_1 + c_1 + c_1 = (2+t)b_1 + a_3 + ua_2 + uc_2.
\end{array}
\right.
$$

Without working through all the details, this can be solved to give

$$
\begin{align*}
& (4u^3 + 36u^2t + 90ut^2 + 72t^3 + 32u^2 + 206ut + 270t^2 + 108u + 322t + 120)c_3 \\
& = (4u^3 + 24u^2t + 54ut^2 + 36t^3 + 44u^2 + 190ut + 198t^2 + 160u + 350t + 200)b_3.
\end{align*}
$$

As previously discussed, Pr(0, 1, 2) = $\frac{1+2t}{9(1+t)}$ and Pr(2, 1, 0) = $\frac{2+t}{9(1+t)}$, so for $b_3 = \mathrm{Pd}(0, 1, 2)$ and $c_3 = \mathrm{Pd}(2, 1, 0)$ to be in the same ratio we would require $b_3 = k(1+2t)$ and $c_3 = k(2+t)$ for some $k$. It follows, after also dividing through by 2, that

$$
\begin{align*}
& (2u^3 + 18u^2t + 45ut^2 + 36t^3 + 16u^2 + 103ut + 135t^2 + 54u + 161t + 60)(t + 2) \\
& = (2u^3 + 12u^2t + 27ut^2 + 18t^3 + 22u^2 + 95ut + 99t^2 + 80u + 175t + 100)(2t + 1).
\end{align*}
$$

This can be expanded to

$$
\begin{align*}
& 2u^3t + 18u^2t^2 + 45ut^3 + 36t^4 + 4u^3 + 52u^2t + 193ut^2 \\
& \qquad + 207t^3 + 32u^2 + 260ut + 431t^2 + 108u + 382t + 120 \\
&= 4u^3t + 24u^2t^2 + 54ut^3 + 36t^4 + 2u^3 + 56u^2t + \\
& \qquad 217ut^2 + 216t^3 + 22u^2 + 255ut + 449t^2 + 80u + 375t + 100.
\end{align*}
$$

This can be reduced to

$$
2u^3t + 6u^2t^2 + 9ut^3 - 2u^3 + 4u^2t + 24ut^2 + 9t^3 - 10u^2 - 5ut + 18t^2 - 28u - 7t - 20 = 0.
$$

This can be factored as

$$
(t - 1)(2u^3 + 6u^2t + 9ut^2 + 10u^2 + 33ut + 9t^2 + 28u + 27t + 20) = 0.
$$

Since $u \ge 0$ and $t \ge 0$, it follows that $t = 1$. This completes the proof in the $\iff$ direction for the DASEP(3, 3, 2) case.
5. PROOF OF THE CONJECTURE FOR DASEP(3, p, 2)

It would be an endless game to prove the conjecture “case by case”, with more and more cumbersome computations, so let us now prove it for an infinite family of models. More precisely, we now prove Conjecture 4.3 for DASEP(3, p, 2) (our previous examples covered the cases p = 2 and p = 3). To solve this case we essentially need to solve for each of p² prior probabilities $p_{i,j} = \text{Pd}(0, i, j)$ for $1 \le i, j \le p$. The steady state probabilities can be obtained by solving a set of p² linear equations each of which essentially demands equilibrium for each of the possible states of the process. The generic form of such an equation, for $i < j$, is given by

$$ (4 + t + 2u)p_{i,j} = (1 + 2t)p_{j,i} + p_{i+1,j} + p_{i,j+1} + up_{i-1,j} + up_{i,j-1}. \quad (5.1) $$

For $i > j$ the equation is

$$ (3 + 2t + 2u)p_{i,j} = (2 + t)p_{j,i} + p_{i+1,j} + p_{i,j+1} + up_{i-1,j} + up_{i,j-1} $$

For $i = j$ the equation simplifies to

$$ (2 + 2u)p_{i,i} = p_{i+1,i} + p_{i,i+1} + up_{i-1,i} + up_{i,i-1}. $$

The equation may be similarly simplified for other edge cases such as $i = 1 < j$, $i < j = p$, $i = 1 < j = p$, $i > j = 1$, $i = p > j$, $i = p > j = 1$, $i = j = 1$, and $i = j = p$. For the sake of brevity we do not list all such cases in detail.

From the first above equation we can define a polynomial $A_{i,j}$ by gathering all terms on the left:

$$ A_{i,j} := (4 + t + 2u)p_{i,j} - (1 + 2t)p_{j,i} - p_{i+1,j} - p_{i,j+1} - up_{i-1,j} - up_{i,j-1}. $$

We can similarly define $A_{i,j}$ under the conditions stated for the various edge cases. We next define a $p^2 \times p^2$ matrix $B$ as follows:

$$ B_{p(i_1-1)+j_1, p(i_2-1)+j_2} = [p_{i_1,j_1}] A_{i_2,j_2}. $$

The next step is to prove that the rank of $B$ is $p^2 - 1$. To see this, we first observe that the sum of all rows of $B$ is identically zero, meaning that the rank cannot be $p^2$. For the rank to then be $p^2 - 1$, we would then need to show that no nontrivial linear combination of a proper subset of the rows can be zero. If we let row $i, j$ be $R_{i,j}$ and for some coefficients $c_{i,j}$ we have $\sum_{i,j} c_{i,j} R_{i,j} = 0$, then we need to show that if any $c_{i,j} = 0$, then all $c_{i,j} = 0$. The only rows with a $t$ term in column $i, j$ will be $R_{i,j}$ and $R_{j,i}$. Hence if $c_{i,j} = 0$, it follows that $c_{j,i} = 0$.

We next show that if $c_{i,j} = 0$ it follows that $c_{i-1,j-1} = 0$. We can do this by first showing that $c_{i-1,j}$ and $c_{i,j-1}$ must be negations of one another. The only rows with a $u$ term in column $i, j$ will be $R_{i,j}$, $R_{i-1,j}$, and $R_{i,j-1}$, with the latter two having the same coefficient. Hence the following two statements are equivalent: $c_{i,j} = 0$ and $c_{i-1,j} + c_{i,j-1} = 0$. We can similarly show that $c_{i,j} = 0$ and $c_{i+1,j} + c_{i,j+1} = 0$ are equivalent. So from $c_{i,j} = 0$ we can derive $c_{i-1,j-1} = 0$. By repeated application of the same argument we will get $c_{k,1} = 0$ or $c_{1,k} = 0$ for some $k$.
Likewise, using the equations for the edge cases $i = 1 < j$ and $i > j = 1$, the only rows with a $u$ term in column 1, $k$ will be $R_{1,k}$ and $R_{1,k-1}$ and the only rows with a $u$ term in column $k$, $1$ will be $R_{k,1}$ and $R_{k-1,1}$. So from $c_{k,1} = 0$ we can derive $c_{k-1,1} = 0$ and from $c_{1,k} = 0$ we can derive $c_{1,k-1} = 0$. By repeated application of this we will get to $c_{1,1} = 0$. By reversing the above arguments it follows that $c_{i,j} = 0$ for any $i, j$ and we have proven:

**Lemma 5.1.** The rank of the matrix $B$ as defined above is $p^2 - 1$.

We next prove a result about the values of the $p_{i,j}$.

**Proposition 5.2.** One has

$$p_{i,j} + p_{j,i} = \frac{2u^{i+j-2}}{\left(\sum_{k=0}^{n-1} u^k\right)^2} \quad \text{and} \quad p_{i,i} = \frac{u^{2i-2}}{\left(\sum_{k=0}^{n-1} u^k\right)^2}.$$

*Proof.* This can be proven by eliminating the variable $t$ from the set of linear equations above. For example, if we add the equations for $i < j$ and $j < i$ we get the following:

$$\begin{aligned} & (4 + t + 2u)p_{i,j} + (3 + 2t + 2u)p_{j,i} \\ & = (2 + t)p_{i,j} + (1 + 2t)p_{j,i} + p_{i+1,j} + p_{j,i+1} + p_{i,j+1} + p_{j+1,i} \\ & \quad + up_{i-1,j} + up_{j,i-1} + up_{i,j-1} + up_{j-1,i}. \end{aligned}$$

If we let $q_{i,j} = p_{i,j} + p_{j,i}$ the above can be simplified to

$$(2 + 2u)q_{i,j} = q_{i+1,j} + q_{i,j+1} + uq_{i-1,j} + uq_{i,j-1}.$$

If we substitute in the values for $q_{i,j}$ from the theorem we are attempting to prove to the above equation, we see that it does satisfy the above equation. Therefore the values of $q_{i,j}$ given in the theorem represent one possible feasible solution to the set of equations. Moreover, via Lemma 5.1 about the rank of $B$, the solution must be unique. This completes the proof. $\square$

To continue with the proof of Conjecture 4.3 in the $\Leftarrow$ direction, we note that from $\frac{\mathrm{Pr}(\mu)}{\mathrm{Pr}(\nu)} = \frac{\mathrm{Pd}(\mu)}{\mathrm{Pd}(\nu)}$ it follows that $\frac{\mathrm{Pr}(0,2,1)}{\mathrm{Pr}(0,1,2)} = \frac{\mathrm{Pd}(0,2,1)}{\mathrm{Pd}(0,1,2)}$ or $\frac{2+t}{1+2t} = \frac{p_{2,1}}{p_{1,2}}$. This expands as $(2+t)p_{1,2} = (1+2t)p_{2,1}$. From the above theorem we know that

$$p_{1,2} + p_{2,1} = \frac{2u}{\left(\sum_{k=0}^{n-1} u^k\right)^2}.$$

We can then solve for $p_{1,2}$ giving

$$p_{1,2} = \frac{2(1 + 2t)u}{3(1 + t) \left(\sum_{k=0}^{n-1} u^k\right)^2}.$$

From the equation (5.1) for $i = 1 < j$ we get

$$(3 + t + 2u)p_{1,2} = (1 + 2t)p_{2,1} + p_{2,2} + p_{1,3} + up_{1,1}.$$

Substitute in to get

$$\frac{2(3 + t + 2u)(1 + 2t)u}{3(1 + t) \left(\sum_{k=0}^{n-1} u^k\right)^2} = \frac{2(2 + t)(1 + 2t)u}{3(1 + t) \left(\sum_{k=0}^{n-1} u^k\right)^2} + \frac{3(1 + t)(1 + u)u}{3(1 + t) \left(\sum_{k=0}^{n-1} u^k\right)^2} + p_{1,3}.$$
This simplifies to

$$p_{1,3} = \frac{(5ut + u + t - 1)u}{3(1 + t) \left(\sum_{k=0}^{n-1} u^k\right)^2}.$$

A similar argument to that used to produce the above equation for $p_{1,2}$ will give us

$$p_{1,3} = \frac{2(1 + 2t)u^2}{3(1 + t) \left(\sum_{k=0}^{n-1} u^k\right)^2}.$$

Equating the last two equations and solving gives us $t = 1$. This completes the proof of **Theorem 5.3.** Conjecture 4.3 holds for $D = \text{DASEP}(3, p, 2)$.

## 6. FUTURE WORK

Three main potential directions for future work are indicated. One is that further results should be obtained with a view to eventually proving Conjecture 4.3. We proved it for DASEP(3, p, 2) and the suggestion would be to prove it for DASEP(n, 2, 2) and DASEP(n, 2, q) before eventually proceeding to DASEP(n, p, q). Similarly considering the case where 0 represents a ball with species 0 rather than the absence of a species is a variant that should be explored. The other, and more ambitious, possible goal for future research would be to come up with a complete combinatorial characterization of the steady state probabilities for the DASEP. For the ASEP, this has been done in [4] and [10] leading to a deep relationship being discovered between the ASEP and Macdonald polynomials.

**Acknowledgments.** Thanks to the editors of this special issue as well as to the anonymous referees for their many helpful suggestions in the review and publication of this paper.

**Funding.** Thanks to Real Time Agents Inc., Pleasant Hill, CA, USA, for partially funding this work.

## REFERENCES

[1] Arvind Ayyer and Svante Linusson. An inhomogeneous multispecies TASEP on a ring. *Advances in Applied Mathematics*, 57:21–43, 2014.

[2] Luigi Cantini. Inhomogeneous multispecies TASEP on a ring with spectral parameters. *arXiv*, 2016.

[3] Luigi Cantini, Jan de Gier, and Michael Wheeler. Matrix product formula for Macdonald polynomials. *Journal of Physics A: Mathematical and Theoretical*, 48:38–64, 2015.

[4] Sylvie Corteel, Olya Mandelshtam, and Lauren Williams. From multiline queues to Macdonald polynomials via the exclusion process. *Proceedings of the 31st Conference on Formal Power Series and Algebraic Combinatorics*, 82B:1–12, 2019.

[5] Pablo Ferrari and James Martin. Stationary distributions of multi-type totally asymmetric exclusion processes. *Annals of Probability*, 35:807–832, 2007.

[6] Donghyun Kim and Lauren Williams. Schubert polynomials and the inhomogeneous TASEP on a ring. *Séminaire Lotharingien de Combinatoire*, 20:1–12, 2021.

[7] Thomas Lam and Lauren Williams. A Markov chain on the symmetric group that is Schubert positive? *Experimental Mathematics*, 21:189–192, 2012.

[8] Ian G. Macdonald. *Symmetric Functions and Hall Polynomials*. Oxford University Press, 1995.

[9] Olya Mandelshtam. Toric tableaux and the inhomogeneous two-species TASEP on a ring. *Advances in Applied Mathematics*, 113:101958, 2020.

[10] James Martin. Stationary distributions of the multi-type ASEP. *Electronic Journal of Probability*, 43:1–41, 2020.
# SIGNED AREA ENUMERATION FOR LATTICE WALKS

STÉPHANE OUVRY<fnref n="1" /> ID AND ALEXIOS P. POLYCHRONAKOS<fnref n="2" /> ID

**Abstract.** We give a summary of recent progress on the signed area enumeration of closed walks on planar lattices. Several connections are made with quantum mechanics and statistical mechanics. Explicit combinatorial formulae are proposed which rely on sums labelled by the multicompositions of the length of the walks.

**Keywords:** lattice walks, signed area enumeration, Hofstadter model, exclusion statistics.

## 1. INTRODUCTION

The seminal problem of the signed area enumeration of walks on planar lattices of various kinds has been around for a long time. It is well known that this purely combinatorial problem can be equivalently reformulated in the realm of Hofstadter-like quantum mechanics models (note that, in physics, the “signed area” is often called the “algebraic area”). Recently, in [16], this problem has been given a boost in the form of an explicit enumeration formula which in turn could be reinterpreted [5, 14, 15] in terms of statistical mechanics models with exclusion statistics, again a purely quantum concept. It is a striking fact that an enumeration quest regarding classical random walks should be in the end so intimately connected to quantum physics.

In this note we give a summary of this recent progress starting with the original signed area enumeration problem for closed walks on a square lattice and then enlarging the perspective to other kinds of lattices and walks via the statistical mechanics reinterpretation. The first question we address is: Among the $\binom{N}{N/2}^2$ closed $N$-step walks that one can draw on a square lattice starting from and returning to a given point (note that $N$ is then necessarily even), how many of them enclose a given signed area $A$?

The signed area enclosed by a directed walk is weighted by its *winding number*: If the walk moves around a region in a counterclockwise direction, its area counts as positive, otherwise it counts as negative; if the walk winds around more than once, the area is counted with multiplicity. These regions inside the walk are called winding sectors.

<footnote n="1">¹LPTMS, CNRS, UNIVERSITÉ PARIS-SACLAY, FRANCE;
http://lptms.u-psud.fr/fr/blog/stphane-ouvry</footnote>

<footnote n="2">²PHYSICS DEPARTMENT, CITY COLLEGE OF NEW YORK, USA;
GRADUATE CENTER, CITY UNIVERSITY OF NEW YORK, USA;
https://www.ccny.cuny.edu/profiles/alexios-polychronakos</footnote>
<figure><img src="image_73.png" /><figcaption>FIGURE 1. A closed walk of length $N = 36$ starting from and returning to the same bullet red point with its various winding sectors $m = +2, +1, 0, -1, -1$ (containing respectively 2, 14, 1, 1, 1 unit lattice cells). Note the double arrow on the horizontal link (above the +2 sector) which indicates that the walk has moved twice on this link, here in the same left direction.</figcaption></figure>

In Figure 1, the 0-winding sector inside the walk arises from a superposition of a +1 and a -1 winding. Summing the areas of each sector, with the corresponding multiplicative weight, gives the signed area $A = (-1) \times 2 + (+0) \times 1 + (+1) \times 14 + (+2) \times 2 = 16$.

More formally, if $\gamma : [0, 1] \to \mathbb{R}^2$ is a closed path that begins and ends at the origin, the signed area of this path is

$$A = \int_{\mathbb{R}^2} \eta(\gamma, \mathbf{x}) d\mathbf{x} = \sum_{m=-\infty}^{\infty} m S_m,$$

where $\eta(\gamma, \mathbf{x})$ is the winding number of $\gamma$ around the point $\mathbf{x} \in \mathbb{R}^2$, and where $S_m$ denotes the classical area of the $m$-winding sectors inside the path (i.e. the number of unit lattice cells it encloses with winding number $m$, where $m$ can be positive or negative).

Winding sectors for continuous Brownian curves as well as for discrete lattice walks have been the subject of study for a long time. In this respect, we note in the last few years some advances in [2] where<fnref n="1" /> an explicit formula for the expected area $\langle S_m \rangle$ of the $m$-winding sectors inside square-lattice walks is proposed, to the exception of the 0-winding sector, for the simple reason that the latter is difficult to distinguish from the outside (i.e. 0-winding again) sector, which is of infinite size. Taking the continuous limit allows us to recover the results previously obtained in [4] for Brownian curves. One notes that for Brownian curves the expected area $\langle S_0 \rangle$ of the 0-winding sectors is also known by other means thanks to the SLE machinery [7]. However, it remains an open problem for discrete lattice walks.

<footnote n="1">¹Note that Timothy Budd gave a talk on his article [2] at the conference Lattice Paths, Combinatorics and Interactions, at CIRM, in 2021. A video of his talk is available on the website of this conference.</footnote>
Counting the number of closed walks of length $N$ on the square lattice enclosing a signed area $A$ can be achieved in a most straightforward way by introducing two lattice hopping operators (symbols) $u$ and $v$ respectively in the right and up directions, as well as operators $u^{-1}$ and $v^{-1}$ corresponding to hops in the left and down direction. A directed walk on the square lattice starting at the origin is then represented by the ordered product of the hopping operators corresponding to its individual steps. By convention we order the operators from right to left as we trace the steps of the walk: $vu$ corresponds to a up step $u$, followed by a right step $v$. Clearly the set of all walks of length $N$ on the lattice is reproduced by the $4^N$ terms in the expansion of

$$ (u + u^{-1} + v + v^{-1})^N \qquad (1.1) $$

into monomials of products of symbols, each with $N$ factors. The operator $u+u^{-1}+v+v^{-1}$ can be considered as the generator of walks.

We are interested in closed walks, and in counting their multiplicity according to their signed area. To this end, we endow the above operators with the relations

$$ uu^{-1} = u^{-1}u = vv^{-1} = v^{-1}v = 1, $$

to which we add a non-commutativity relation (which expresses the fact that the elementary walk circling one lattice cell in the counterclockwise direction has signed area 1):

$$ v^{-1} u^{-1} v u = Q, $$

where $Q$ is a central element (that is, $Q$ commutes with all operators). This entails

$$ vu = Quv, \quad vu^{-1} = Q^{-1}u^{-1}v, \quad v^{-1}u = Q^{-1}uv^{-1}, \quad v^{-1}u^{-1} = Qu^{-1}v^{-1}, \qquad (1.2) $$

which allows us to reduce all terms in (1.1) into monomials of the form $u^m v^n$, $(m, n) \in \mathbb{Z}^2$ being the lattice coordinates of the end of the walk, with coefficients which are powers of $Q$. In particular, closed walks correspond to the monomial $u^0 v^0$ in (1.1).

The non-commutativity relation $vu = Quv$ has the effect of flipping a right-up two-step segment into an up-right one, producing a factor of $Q$, and similarly for each relation in (1.2). In each case, the coefficient is $Q$ to the power of the signed area of the unit lattice cells left behind by the exchange. Repeated application of these relations reduces each walk to hook-shape walk $u^m v^n$ with an overall coefficient $Q^A$, with $A$ the signed area of the original walk prolonged into a closed walk by joining its end to the origin with a vertical and a horizontal straight walk. In particular, closed walks correspond to monomials $Q^A u^0 v^0$, with $A$ the signed area of the walk. This area is maximal if the walk forms a square of width $N/4$, and then the signed area is $\pm N^2/16$. Therefore, using the notation $[u^0 v^0]$ for the extraction of the constant term in a Laurent polynomial in $u$ and $v$, the distribution of the signed area of closed walks of length $N$ is given by

$$ [u^0 v^0] (u + u^{-1} + v + v^{-1})^N = \sum_{A=-\lfloor N^2/16 \rfloor}^{\lfloor N^2/16 \rfloor} C_N(A) Q^A, \qquad (1.3) $$

where $C_N(A)$ counts the closed walks of length $N$ enclosing a signed area $A$. For example, one easily checks that $[u^0 v^0](u + u^{-1} + v + v^{-1})^4 = 28 + 4Q + 4Q^{-1}$, indicating that among the $\binom{4}{2}^2 = 36$ closed walks making 4 steps $C_4(0) = 28$ enclose a signed area $A = 0$ and $C_4(1) = C_4(-1) = 4$ enclose a signed area $A = \pm 1$.
## 2. THE HOFSTADTER MODEL

In any irreducible representation of the operator relation $vu = Quv$, the central element $Q$ will be represented by a number. Restricting to unitary representations, for which $u^\dagger = u^{-1}$, $v^\dagger = v^{-1}$, $Q$ will necessarily be a complex number of norm unity, i.e., a phase. This provides a mapping between the $u, v$ representation for walks and quantum mechanics, interpreting $u$ and $v$ as unitary operators acting on a quantum Hilbert space, and the operator $u + u^{-1} + v + v^{-1}$ as the Hamiltonian of a quantum system.

In fact, such a quantum system exists and corresponds to a well-known model in physics. Interpreting $u$ and $v$ as operators that generate hops of a quantum particle by one link on the square lattice, the non-commutativity relation $vu = Quv$ indicates that translations of the particle in the horizontal and vertical directions do not commute. This can be interpreted as that the particle is charged and coupled to a homogeneous magnetic field perpendicular to the lattice. The magnetic flux associated to this (constant) vector field is $\Phi a^2$ for any surface of area $a^2$.

Let us now take $Q = e^{i2\pi\Phi/\Phi_o}$, where $\Phi$ is the magnetic flux through any unit lattice cell (i.e. for the surface of the square of width $a = 1$) and where $\Phi_o = h/c$ is the flux quantum ($h$ is the Planck constant and $c$ the particle's charge). The Hermitian operator

$$H = u + u^{-1} + v + v^{-1} \qquad (2.1)$$

then becomes a Hamiltonian modelling a quantum particle hopping on a square lattice and coupled to a perpendicular magnetic field. This model is known as the Hofstadter model [8].

To make the physics connection completely explicit, we note that in quantum mechanics the hopping operators $u$ and $v$ are written as

$$u = e^{i(p_x - cA_x)/\hbar} \quad \text{and} \quad v = e^{i(p_y - cA_y)/\hbar},$$

where $A_x = -\Phi y$ and $A_y = 0$ are the two components of the vector potential of the magnetic field in the Landau gauge and $p_x = -i\hbar\partial_x$ and $p_y = -i\hbar\partial_y$ those of the momentum operator (where we use the standard notation $\hbar = h/(2\pi)$). Thus, assuming (1), the relation $vu = Quv$ follows from the Baker-Campbell-Hausdorff formula, using the Heisenberg commutators $[x, p_x] = [y, p_y] = i\hbar$ (this identity is often called the “canonical commutation relation”).

One now introduces the quantum state $\Psi_{m,n}$ representing the probability amplitude of the particle being at lattice site $(m, n)$, on which hopping operators act as

$$u\Psi_{m,n} = e^{icn\Phi/\hbar}\Psi_{m+1,n}, \quad v\Psi_{m,n} = \Psi_{m,n+1}.$$

As $c\Phi/\hbar = 2\pi\Phi/(h/c) = 2\pi\Phi/\Phi_o$, the factor appearing in the action of $u$ on $\psi_{m,n}$ is $Q^n$. Using translation invariance in the horizontal direction we can further choose $\Psi_{m,n}$ to be an eigenstate of $p_x$, that is, $\Psi_{m,n} = e^{imk_x}\Phi_n$. The action of $u, v$ on $\Phi_n$ becomes

$$u\Phi_n = e^{ik_x}Q^n\Phi_n, \quad v\Phi_n = \Phi_{n+1}.$$

For the Hofstadter model, the Schrödinger equation $H\Psi = E\Psi$ (which determines the eigenvalue $E$ of the spectrum) can be rewritten as

$$(u + u^{-1} + v + v^{-1})\Psi_{m,n} = E\Psi_{m,n} \Rightarrow \Phi_{n+1} + \Phi_{n-1} + (Q^n e^{ik_x} + Q^{-n}e^{-ik_x})\Phi_n = E\Phi_n.$$
Going a step further, a simplification arises when the flux is rational (i.e. when one has $Q = e^{i2\pi p/q}$ with $p, q$ two coprime integers): It induces a $q$-periodicity of the Schrödinger equation in the vertical direction. Then, as Bloch's theorem states that solutions to the Schrödinger equation in a periodic potential take the form of a plane wave modulated by a periodic function $\tilde{\Phi}_n$, we can write

$$ \Phi_n = e^{ink_y} \tilde{\Phi}_n, \quad \tilde{\Phi}_{n+q} = \tilde{\Phi}_n. \qquad (2.2) $$

Indeed, for this rational flux $Q^q = 1$, thus $u^q, v^q$ become Casimir (a physicist's term for central elements), and the choice of Bloch states (2.2) can be interpreted mathematically as choosing an irreducible representation of the $u, v$ algebra. Acting on such states, $u^q$ and $v^q$ become $u^q = e^{iqk_x}$ and $v^q = e^{iqk_y}$. One ends up with $u$ and $v$, acting on $\tilde{\Phi}_n$, becoming the $q \times q$ matrices

$$ u = e^{ik_x} \begin{pmatrix} Q & 0 & 0 & \cdots & 0 & 0 \\ 0 & Q^2 & 0 & \cdots & 0 & 0 \\ 0 & 0 & Q^3 & \cdots & 0 & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & 0 & \cdots & Q^{q-1} & 0 \\ 0 & 0 & 0 & \cdots & 0 & 1 \end{pmatrix} \quad \text{and} \quad v = e^{ik_y} \begin{pmatrix} 0 & 1 & 0 & \cdots & 0 & 0 \\ 0 & 0 & 1 & \cdots & 0 & 0 \\ 0 & 0 & 0 & \cdots & 0 & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & 0 & \cdots & 0 & 1 \\ 1 & 0 & 0 & \cdots & 0 & 0 \end{pmatrix} $$

involving two real quantities $k_x$ and $k_y$. Finding the energy spectrum (which depends on $k_x$ and $k_y$) reduces to computing the eigenvalues $E_1, \dots, E_q$ of the $q \times q$ Hamiltonian matrix $H_q := u + u^{-1} + v + v^{-1}$, i.e.

$$ H_q = \begin{pmatrix} Qe^{ik_x} + Q^{-1}e^{-ik_x} & e^{ik_y} & 0 & \cdots & 0 & e^{-ik_y} \\ e^{-ik_y} & Q^2e^{ik_x} + Q^{-2}e^{-ik_x} & e^{ik_y} & \cdots & 0 & 0 \\ 0 & e^{-ik_y} & () & \cdots & 0 & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & 0 & \cdots & () & e^{ik_y} \\ e^{ik_y} & 0 & 0 & \cdots & e^{-ik_y} & Q^q e^{ik_x} + Q^{-q}e^{-ik_x} \end{pmatrix}. $$

All the machinery of quantum mechanics is now at our disposal. Selecting as in (1.3) the $u^0v^0$ monomial of $(u+u^{-1}+v+v^{-1})^N$ translates in the quantum world to computing the trace of $H_q^N$. The quantum trace is defined as

$$ \mathbf{Tr} H_q^N := \frac{1}{q} \int_{-\pi}^{\pi} \int_{-\pi}^{\pi} \frac{dk_x}{2\pi} \frac{dk_y}{2\pi} \mathrm{tr} H_q^N = \frac{1}{q} \int_{-\pi}^{\pi} \int_{-\pi}^{\pi} \frac{dk_x}{2\pi} \frac{dk_y}{2\pi} \sum_{i=1}^{q} E_i^N, \qquad (2.3) $$

that is, one sums over the $q$ eigenvalues $E_i$ of $H_q$ (yielding the standard matrix trace $\mathrm{tr} H_q^N$) and integrates over $k_x$ and $k_y$ while enforcing a continuous normalization in $k_x, k_y$ and one rescales by a factor $1/q$ (thus, if one considers for example the $q \times q$ identity matrix $I_q$, one has $\mathrm{tr} I_q = q$, while $\mathbf{Tr} I_q = 1$). Under this definition of the trace, $\mathbf{Tr} u^m v^n = \delta_{m,0}\delta_{n,0}$ (integration over $k_x, k_y$ eliminates the traces of terms involving $u^{qm}$ and $v^{qn}$). We thus get a first noteworthy result (also obtained via another approach by Bellisard et al. in [1]):

**Theorem 2.1.** Assuming $q > \frac{N^2}{8}$, the signed area enumeration of closed paths is given by

$$ \sum_A C_N(A) Q^A = \mathbf{Tr} H_q^N. \qquad (2.4) $$