Algebraic Walk Theory

P.-L Giscard

► To cite this version:

P.-L Giscard. Algebraic Walk Theory. Combinatorics [math.CO]. Université du Littoral - Côte d'Opale, 2023. tel-04323927v2

HAL Id: tel-04323927

https://hal.science/tel-04323927v2

Submitted on 8 Jan 2024

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
Document de synthèse en vue de
**L'Habilitation à Diriger des Recherches**

Discipline:

**Mathématiques**

Ecole doctorale Sciences, Technologie, Santé (ED585)

Présentée par

**Pierre-Louis GISCARD**

**Algebraic Walk Theory**

Soutenance: 15 Décembre 2023

Devant le jury composé de :

<table><tr><td>M. Cyril BANDERIER</td><td>CNRS - Université Paris-Nord, Rapporteur</td></tr><tr><td>M. Christian BONHOMME</td><td>Sorbonne Université, Rapporteur</td></tr><tr><td>M. Shalom ELIAHOU</td><td>Université du Littoral Cte d'Opale, Examinateur</td></tr><tr><td>M. Loïc FOISSY</td><td>Université du Littoral Cte d'Opale, Examinateur</td></tr><tr><td>M. Guoniu HAN</td><td>CNRS - Université de Strasbourg, Examinateur</td></tr><tr><td>M. Dominique MANCHON</td><td>CNRS - Université Clermont-Auvergne, Examinateur</td></tr><tr><td>M. Paul ROCHET</td><td>École Nationale de l'Aviation Civile, Examinateur</td></tr><tr><td>Mme. María RONCO</td><td>Universidad de Talca (Chile), Rapporteure</td></tr></table>
<figure><img src="image_2.png" /></figure>
# Foreword

## Acknowledgements

I am heavily indebted to my family, in particular to my wife, for the profound support I received from her throughout the years that we shared.

I am grateful to the people of the Laboratoire de Mathématiques Pures et Appliquées Joseph Liouville: firstly, they gave me food and shelter in the form of a permanent research position, ending years of wandering and anxiety. Secondly, they have proven to be marvellous colleagues, great Magic players and mathematical luminaries, making the LMPA into a remarkably enjoyable work environment on top of a place of constant mathematical learning.

Beyond the confines of the LMPA, I would like to warmly thank Christian Bonhomme for contacting me out of the blue back in 2017, starting an ongoing scientific adventure and to me, a friendship. Christian is one of the nicest human beings I have met. Just as remarkable is Stefano Pozza: it takes some incredible daring, confidence and expertise to jump into ★-Lanczos, invent the ★-product and orient your career from there, all after only one hour's worth of discussion in a Scottish visitor's office. I am glad to have the chance to be among Stefano's collaborators, witnessing in awe his regular achievements.

María Ronco deserves a special mention for her patience, humility and brilliance: *patience* because I first contacted her in 2016 about a research project which we still have not finished since every step along the way I need to learn something anew; *humility* because I was a postdoctoral researcher and she is a researcher of international stature who not only answered my unclear questions and even came to visit but she also consistently took the time and effort to come down to my level so the project would progress; *brilliance* because she regularly discovers fundamental objects and laws of algebra yet seems oblivious to the size of her achievements. María is uniquely fluent in mathematics and languages alike and a constant traveler: I hope to be part of her voyage in the years to come.

## On real authorship

Almost all of the results presented in this *Habilitation à Diriger des Recherches* are the fruits of extended collaborations. These results should by no means be attributed to myself alone and the present thesis should certainly not be taken to be a claim of this. The researchers who, by their work, contributed to ‘walk theory’ and are therefore effectively authors of this thesis are (in alphabetical order):

Franka Bause, Tolga Birkandan, Christian Bonhomme, Zheng Choo, Loïc Foissy, Mohammadali Foroozandezeh, Jean Fromentin, Dieter Jaksch, Théo Karaboghossian, Nils Kriege, Kevin Lui, Cécile Mammez, Stefano Pozza, Paul Rochet, María Ronco, Manon Ryckebusch, Aditya Tamar, Simon James Thwaite, Richard C. Wilson.

Finally, there are those discreet people who shared their profound knowledge with me over discussions, leading to advances in the present research (and diminishing my vast ignorance a bit), notably Abderrahman Bouhamidi (always enthusiastic about distributions), Shalom Eliahou (who consistently delivers the most captivating seminars I have attended), Jean Fromentin (the living embodiment of mathematics *and* computer science within a single person), Loïc Foissy (who apparently thinks it normal to solve an arduous problem and write a detailed account of the solution plus context and further results in LaTeX while also heading the maths department, *all in the space of a few days*) and, of course, Christian Miebach (whom I secretly believe to be the best mathematician around).

## Mode of presentation

Mathematicians probably prefer proofs to theorems. Proofs are often—though not always—either trivial or lengthy and technical. For this reason and given my goal of presenting an overarching view of a theory of walks, I have decided in all but a few cases to either briefly sketch or completely omit proofs while concentrating on results. The proofs are all available and published (hence accessible to the inquisitive reader) but the tentative vision of a coherent theory of walks offered by a mass of results concentrated in a single report can only be found here.
<figure><img src="image_3.png" /></figure>
# Contents

## Introduction

<table><tr><td>1</td><td>The two branches of walk theory</td><td>9</td></tr><tr><td>1.1</td><td>NAP-co-preLie bialgebra on walks</td><td>9</td></tr><tr><td>1.2</td><td>Hike monoids</td><td>10</td></tr></table>

## Part 1: The fundamental bialgebraic structure on walks

<table><tr><td>2</td><td>NAP-co-preLie bialgebra on graph walks</td><td>12</td></tr><tr><td>2.1</td><td>Co-preLie co-algebra of walks</td><td>13</td></tr><tr><td>2.2</td><td>Non-associative permutative algebra of walks</td><td>15</td></tr><tr><td>2.3</td><td>NAP-co-preLie bialgebra</td><td>16</td></tr><tr><td>2.4</td><td>Path-Sum theorem</td><td>17</td></tr></table>

<table><tr><td>3</td><td>Implications and applications to differential calculus</td><td>19</td></tr><tr><td>3.1</td><td>The ⋆-product</td><td>20</td></tr><tr><td>3.2</td><td>Path-sum representation of the solution of systems of non-autonomous ODEs</td><td>26</td></tr><tr><td>3.3</td><td>Tridiagonalization of differential systems: ⋆-Lanczos algorithm</td><td>29</td></tr></table>

<table><tr><td>4</td><td>Further applications of the path-sum theorem</td><td>32</td></tr><tr><td>4.1</td><td>Statistical inference</td><td>32</td></tr><tr><td>4.2</td><td>Quantum dynamics</td><td>33</td></tr><tr><td>4.3</td><td>Representation of Heun functions</td><td>38</td></tr></table>

## Part 2: Extending number theory to walks on graphs

<table><tr><td>5</td><td>Walk theory: combinatorics on hike monoids</td><td>42</td></tr><tr><td>5.1</td><td>On the necessity of a theory of walks</td><td>42</td></tr><tr><td>5.2</td><td>Hike monoids and number theory</td><td>44</td></tr><tr><td>5.3</td><td>Algebraic relations between hikes series</td><td>48</td></tr></table>

<table><tr><td>6</td><td>Number-theoretic sieves on hikes</td><td>52</td></tr><tr><td>6.1</td><td>Asymptotic enumeration of walks by their last erased loop</td><td>52</td></tr><tr><td>6.2</td><td>Proof of the infinite sieve theorem</td><td>55</td></tr><tr><td>6.3</td><td>Relation with the self-avoiding polygons counting problem</td><td>60</td></tr></table>

<table><tr><td>7</td><td>Exact enumeration formulas for simple cycles and simple paths</td><td>62</td></tr><tr><td>7.1</td><td>Counting from weakly connected subgraphs</td><td>62</td></tr><tr><td>7.2</td><td>Algorithmic analysis of the counting formulae</td><td>64</td></tr><tr><td>7.3</td><td>Counting via convolutions over induced subgraphs</td><td>67</td></tr><tr><td>7.4</td><td>An Hopf algebra structure for the hikes</td><td>69</td></tr></table>

<table><tr><td>8</td><td>On the tenuous relation between walks and their graphs</td><td>72</td></tr><tr><td>8.1</td><td>Realizable cycle structures in digraphs</td><td>72</td></tr><tr><td>8.2</td><td>Simple graphs do not relate well with their simple cycles</td><td>76</td></tr><tr><td>8.3</td><td>Invariants of hike monoids</td><td>77</td></tr></table>
<table><tr><td>9</td><td>Further applications of hike monoids</td><td>80</td></tr><tr><td>9.1</td><td>Cycle counting in social analysis and computer science</td><td>80</td></tr><tr><td>9.2</td><td>Sieves on hikes as a cycle centrality measure</td><td>82</td></tr></table>

Conclusion

<table><tr><td>10</td><td>And from there...</td><td>86</td></tr><tr><td>10.1</td><td>Exotic umbral calculi for path-sum solutions to non-linear ODEs and PDEs</td><td>86</td></tr><tr><td>10.2</td><td>Sieves on hikes</td><td>88</td></tr></table>

<table><tr><td>List of Symbols</td><td>90</td></tr></table>

<table><tr><td>Bibliography</td><td>92</td></tr></table>
# Publications

This is a list of published material I authored or co-authored since completion of my DPhil in August 2014.

## Articles in peer-reviewed journals

1. P.-L. Giscard and M. Foroozandeh. Exact solutions for the time-evolution of quantum spin systems under arbitrary waveforms using algebraic graph theory, Computer Physics Communications (2022), p. 108561.

2. T. Birkandan, P.-L. Giscard, and A. Tamar. Computations of general Heun functions from their integral series representations. Proceedings of the 2021 Days on Diffraction (DD), 2021, pp. 1218.

3. J. Fromentin, P.-L. Giscard, and T. Karaboghossian. Realizable cycle structures in digraphs, European Journal of Combinatorics, Vol. 113, 103748 (2023).

4. P.-L. Giscard. Counting walks by their last erased self-avoiding polygons using sieves, Discrete Mathematics 344.4 (2021), p. 112305.

5. P.-L. Giscard and S. Pozza. Tridiagonalization of systems of coupled linear differential equations with variable coefficients by a Lanczos-like method, Linear Algebra and its Applications 624 (2021), pp. 153173.

6. P.-L. Giscard and A. Tamar. Elementary Integral Series for Heun Functions. With an Application to Black-Hole Perturbation Theory. Editor's pick at Journal of Mathematical Physics, Vol.63, Issue 6: 063501 (2022).

7. P.-L. Giscard. On the solutions of linear Volterra equations of the second kind with sum kernels. In: Journal of Integral Equations and Applications 32 (4 2020), pp. 429445.

8. P.-L. Giscard and C. Bonhomme. Dynamics of quantum systems driven by time-varying Hamiltonians: Solution for the Bloch-Siegert Hamiltonian and applications to NMR. In: Phys. Rev. Research 2 (2 04/2020), p. 023081.

9. P.-L. Giscard and S. Pozza. Lanczos-like algorithm for the time-ordered exponential: The $\star$-inverse problem. In: Applied Mathematics (2020), pp. 121.

10. P.-L. Giscard, N. Kriege, and R. C. Wilson. A general purpose algorithm for counting simple cycles and simple paths of any length, Algorithmica 81 (7 2019), pp. 27162737.

11. P.-L. Giscard, S. Pozza, Lanczos-like algorithm for the time-ordered exponential: the $\*$-inverse problem, Applied Mathematics, 1-21 (2020).

12. N. Kriege, P.-L. Giscard, F. Bause, R. C. Wilson, Computing Optimal Assignments in Linear Time for Graph Matching, 19th IEEE International Conference on Data Mining (01/2019).

13. P.-L. Giscard and P. Rochet. Enumerating Simple Paths from Connected Induced Subgraphs, Graphs and Combinatorics 34(6): 1197-1202, (2018).

14. P.-L. Giscard, P. Rochet, and R. C. Wilson. A Hopf algebra for counting cycles, Discrete Mathematics 341.5 (2018), pp. 14391448.

15. P.-L. Giscard and R. C. Wilson. A centrality measure for cycles and subgraphs II, Applied Network Science 3.1 (06/2018), p. 9.

16. P.-L. Giscard and R. C. Wilson. Cycle-Centrality in Economic and Biological Networks, in: *Complex Networks & Their Applications VI*. Ed. by C. Cherifi et al. Cham: Springer International Publishing, 2018, pp. 1428.
17. P.-L. Giscard and P. Rochet. Algebraic Combinatorics on Trace Monoids: Extending Number Theory to Walks on Graphs, SIAM Journal on Discrete Mathematics 31.2 (2017), pp. 14281453.

18. P.-L. Giscard, P. Rochet, and R. C. Wilson. Evaluating balance on social networks from their simple cycles, Journal of Complex Networks 5.5 (2017), pp. 750775.

19. P.-L. Giscard et al. Exact Inference on Gaussian Graphical Models of Arbitrary Topology using Path-Sums, Journal of Machine Learning Research 17.71 (2016), pp. 119.

20. N. M. Kriege, P.-L. Giscard, and R. C. Wilson. On Valid Optimal Assignment Kernels and Applications to Graph Classification, Proceedings of the 30th International Conference on Neural Information Processing Systems. NIPS16. Barcelona, Spain: Curran Associates Inc., 2016, pp. 16231631.

21. P.-L. Giscard, K. Lui, S. J. Thwaite, and D. Jaksch. An exact formulation of the time-ordered exponential using path-sums, Journal of Mathematical Physics 56.5 (2015), p. 053503.

## Presentations in peer-reviewed conferences

1. P.-L. Giscard and C. Bonhomme. "Undamped and Damped Bloch Equations: Exact Solutions by Path-Sum". 61st Experimental Nuclear Magnetic Resonance Conference, Baltimore, Maryland. 2020.

2. C. Bonhomme and P.-L. Giscard. Theory of Solid state NMR: from Dyson, Magnus, Feynman to Path-Sum, Alpine Conference on Magnetic Resonance in Solids. 2019.

3. P.-L. Giscard and C. Bonhomme. A new approach of time-ordered exponential in NMR: the path-sum. In: 60th Experimental Nuclear Magnetic Resonance Conference. 2019.

## Preprints

1. L. Foissy, P.-L. Giscard, C. Mammez. A co-preLie structure from chronological loop erasure in graph walks, arXiv preprint 2306.17605 (2023)

2. P.-L. Giscard and R. C. Wilson. The All-Paths and Cycles Graph Kernel, preprint arXiv 1708.01410 (2017).
# Chapter 1

# The two branches of walk theory

All the results presented here realize a coherent albeit incomplete body of work on a *theory of walks*, its necessity, applications and *implications*. Though they may seem disparate at first sight, they fall into two distinct but closely related branches, which differ only in the role played by the initial vertex of closed walks.

## 1.1 NAP-co-preLie bialgebra on walks

The first branch concerns a bialgebraic structure on graph walks and its consequences throughout mathematics. The structure arises from a simple procedure, now known as *Lawler's loop erasing* [123], first conceived in the context of percolation theory to randomly generate simple paths–walks where all vertices are distinct– from a sample of random walks [115, 121, 119, 116]. The procedure consists of a chronological removal of cycles (called loops in Lawler's original work) as one walk along on the graph: consider for instance the complete graph $K_4$ on 4 vertices and label these vertices with integers 1 through 4. Walking along the path $1 \to 2 \to 1 \to 3 \to 4 \to 3 \to 1 \to 3$ on the graph and removing cycles whenever they appear, we are left with the simple path $1 \to 3$ after having successively ‘erased’ the cycles $1 \to 2 \to 1$, then $3 \to 4 \to 3$ and finally $1 \to 3 \to 1$. Note how $1 \to 3 \to 1$ does not appear contiguously in the original walk. Once terminated, Lawler’s loop-erasing has eliminated a set of cycles, all of whose internal vertices are distinct, leaving a possibly trivial walk-skeleton behind. If the initial walk was itself a cycle, this skeleton is the empty walk on the initial vertex (also called length-0 walk) and otherwise it is a self-avoiding path. Remark that because the loop-removal occurs in a chronological fashion, Lawler’s process is strongly non-Markovian: complete knowledge of all the past steps of a walk is required to decide the current and future erased sections at any point of the walk.

In Chapter 2 we show that Lawler’s procedure, which formally sends a walk to a family of closed walks and at most one open path, correspond to a pre-Lie co-product. This dual of this co-structure is a non-associative permutative algebra based on a product called nesting, which consists of inserting walks into another one according to specific rules. Nesting induces a remarkable property which led to its discovery even before its relation to Lawler’s process was uncovered: any walk admits a unique factorisation into nesting products of simple cycles and simple paths. These, it turns out, are exactly the primitive elements of the co-preLie structure and precisely those loops and walk skeleton produced by Lawler’s loop-erasing.

The unique factorisation property sustains extremely powerful combinatorial results, allowing in particular for a sophisticated resummation of all walks between any two vertices of a graph into a single branched continued fraction of finite depth and breadth over the simple cycles and simple paths of the graph, known as path-sum. Path-sums, which were first devised from purely combinatorial arguments during my PhD thesis, have since found an astounding number of applications, in quantum mechanics, black hole theory, chemistry and in various fields of mathematics including probability theory, linear algebra, differential calculus and numerical analysis. Somme of these are presented in Chapter 4. Surprisingly, this line of ideas has proven particularly fertile in its *implications* even more than through its applications. Indeed, in the case of differential equations, the very fact that a path-sum formulation of the solution of a problem *ought to* exist implied the existence of novel mathematical results. As an example, we shall cite in particular the $\star$-product on distributions, the tridiagonalization of systems of non-autonomous Ordinary Differential Equations (ODEs) and more; subjects which, at first sight, might seem very far from walks indeed. This area of research, the results of which are in Chapter 3, is now especially active notably in the pursuit of novel formal solutions to systems of non-autonomous, fractional, non-linear ODEs as well as Partial Differential Equations (PDEs).
1.2 Hike monoids

The second branch of research results from quotienting out the initial vertex (i.e. the root) of the walks in the bialgebraic structure mentioned above. This is a rather dramatic step as in Lawler's process, which is inherently chronological, the initial vertex is the start of history and the point from which all erasing decisions are taken. Dually, for nesting, the initial vertex of a closed walk determines the few places where a walk may be inserted into another. To properly describe the objects obtained from walks and sets of walks on removing the information about the starting point necessitates a return to a lower level presentation: words on the formal alphabet directed edges $w_{i \to j} \in E$ of the graph $G = (V, E)$. Changing the starting point, we see that the backtrack $w_{1 \to 2}w_{2 \to 1}$ must represent the same object as $w_{2 \to 1}w_{1 \to 2}$. Since exchanging $w_{1 \to 2}$ and $w_{2 \to 1}$ passes from one presentation to the other, this implies in general that two edges with different starting points must be allowed to commute in the word representations of the same object. Meanwhile, walks $1 \to 2 \to 1 \to 3 \to 1$ and $1 \to 3 \to 1 \to 2 \to 1$ remain distinct. So the words $w_{1 \to 2}w_{2 \to 3}w_{3 \to 1}$ and $w_{1 \to 3}w_{3 \to 2}w_{2 \to 1}$ represent distinct objects. This implies that edges with identical starting points–here $w_{1 \to 2}$ and $w_{1 \to 3}$–must not be allowed to commute.

These two rules define the Cartier-Foata monoids. We shall see that a purely cycle-based vision also exists for these monoids and closed walks are now really heaps of cycles on the graph, the combinatorial properties of which are described by semi-commutative monoids called *hike monoids*. Algebraically and combinatorially these monoids obey a semi-commutative extension of number theory presented in Chapter 5. This approach is particularly well suited to enumeration purposes, in Chapter 7 we present novel formulas for counting the simple cycles and simple paths as well as the best general purpose algorithm for this task. This is because all the theoretical tools developed for number theory continue to exist in walk theory, where they share the same relations with one-another. We may therefore exploit these well-understood tools in the novel context of walks. In this framework, the semi-commutative extension of the Prime Number Theorem is a famous open problem, namely that of counting asymptotically the number of self-avoiding polygons of length $\ell$ as $\ell \to \infty$ on infinite regular lattices. Pursuing this idea, in Chapter 6 we use number theoretic sieves to obtain closed form finite, explicit results for enumerating walks per their last erased loops on infinite lattices.

The need for a dedicated mathematical field of study devoted to walks independently from graph theory is exemplified by results from this branch. In Chapter 8 we show that the hike monoidal structure of walks is essentially independent from the topology of the graphs that sustains them. In particular, deciding whether an arbitrary Cartier-Foata monoid is a hike monoid, i.e. if some arrangement of cycles on a graph exists at all, turns out to be extremely difficult.

Applications of the results of this branch beyond pure mathematics are presented in Chapter 9. Here we use the number-theoretic sieves to rate the importance of cycles in real-world networks. As these cycles represent feedback processes in the exchanges between the entities forming the complex system modeled by the network, ranking their importance might lead to uncovering properties of the system itself. For example we show how the sieves are sensitive enough to predict which plant proteins will be targeted by pathogens and which proteins participate in biological complexes.
# Part 1: The fundamental bialgebraic structure on walks
# Chapter 2

## NAP-co-preLie bialgebra on graph walks

We begin by showing that a co-preLie structure on graph walks arises naturally from Lawler's loop erasing. In addition, slightly relaxing the chronological constraints by allowing simultaneous erasures under some compatibility conditions leads to Hopf algebra structures on the tensor and symmetric algebras of graph walks.

In §2.2 we show that Lawler's process is also naturally associated with a non-associative permutative product, known as nesting [76], which satisfies the Livernet compatibility condition [125] with the co-preLie co-product defined here. This provides the very first concrete example of the NAP - co-preLie operad in a 'living' context. This construction appears to be of paramount importance given the pervasive use of graph-walks in mathematics and mathematical-physics. In particular this leads to a useful bridge between formal sums over infinite families of walks and branched continued fractions, the path-sum theorem §2.4.

### Notations and definitions

We consider a non-empty digraph $G$ as defined throughout this work. The set of all walks on a graph $G$ is denoted $W(G)$.

Consider a walk $\omega = w_0 \dots w_\ell$. A subwalk of a walk $\omega = w_0 \dots w_\ell$ is any walk $w_k \dots w_{k'}$ where $0 \le k \le k' \le \ell$. If $k \ne k'$ and $w_k = w_{k'}$, we designate by $\omega^{k,k'} := w_k w_{k+1} \dots w_{k'}$ the closed subwalk of $\omega$ with root $w_k$. In a complementary way, we define the remainder section $\omega_{k,k'} := w_0 \dots w_k w_{k'+1} \dots w_\ell$ to be what remains of $\omega$ after removal of the section $\omega^{k,k'}$. Note, for convenience we denote $\omega_{k,k'}^{l,l'}$ for $(\omega_{k,k'})^{l,l'}$, the section $w_l \dots w_{l'}$ erased from the remainder $\omega_{k,k'} = w_0 \dots w_k w_{k'+1} \dots w_\ell$. This means in particular that in $\omega_{k,k'}^{l,l'}$, integers $k, k', l$ and $l'$ all refer to indices from $\omega$.

A rooted walk in which all vertices are distinct is said to be a *simple path* or *self-avoiding walk*. The set of all such walks on a digraph $G$ is denoted $SAW(G)$. The set of all simple cycles on $G$ is $SAP(G)$. For $G$ any (directed multi)graph, to ease the notation, we also denote by $W(G)$ the $K$-vector space spanned by all walks, $K$ being a field of characteristic 0. For a walk $\omega \in W(G)$, we designate by $V(\omega)$ the support of $\omega$, that is the set of *distinct* vertices visited by $\omega$; and by $E(\omega)$ the *multiset* of directed edges visited by $\omega$.

As stated above, Lawler's loop-erasing procedure consists in erasing all cycles from a walk $\omega$ in the *chronological* order in which they appear. Formally, it is a selection-quotient process which transforms a walk into its self-avoiding skeleton. To construct the algebraic structures associated with Lawler's procedure we must not only consider its end product but also what it produces during its intermediary stages and what it removes from the walk, in its original context:

**Definition 2.0.1 (Loop-erased sections).** Let $G$ be a digraph and consider $\omega = w_0 \dots w_\ell \in W(G)$. The set $\text{LES}(\omega)$ of loop-erased sections is the set of all *closed subwalks* of $\omega$ *erased* by Lawler's procedure.

**Example 2.0.1.** On the complete graph $K_5$ on 5 vertices (including self-loops), consider the walk $\omega = 1234522$. The simple cycles erased by Lawler's procedure are $\omega^{1,3} = 232$, $\omega^{3,6} = 2452$ and $\omega^{6,7} = 22$ and the set of erased closed subwalks of $\omega$ is therefore,

$$ \text{LES}(1234522) = \{\omega^{1,3}, \omega^{3,6}, \omega^{1,6}, \omega^{6,7}, \omega^{3,7}, \omega^{1,7}\} = \{232, 2452, 232452, 22, 24522, 2324522\}. $$

**Definition 2.0.2 (Loop-erased walks).** Let $G$ be a digraph and $\omega = w_0 \dots w_\ell \in W(G)$ of length $\ell$. For $0 \le k \le \ell$, we designate $\text{LEW}_k(\omega)$, called loop-erased walk $\omega$ at step $k$, to be what is left of $\omega$ after its first $k$ steps while performing Lawler's procedure.

By the definitions of $\text{LES}(\omega)$ and $\text{LEW}(\omega)$ we obtain what was remarked above, namely that loop-erased sections may not straddle over one-another, a consequence of their step-by-step erasure in chronological order:
**Lemma 2.0.1.** Let $G$ be a digraph and $\omega = w_0 \dots w_\ell \in W(G)$. Then $\omega^{k,k'} \in \text{LES}(\omega)$ if and only if there does not exist a pair of integers $0 \le l < k < l' < k' \le \ell$ with $w_k = w_{k'} \ne w_l = w_{l'}$ and $\omega^{l,l'} \in \text{LES}(\omega)$

We remark that the notion of loop-erased walks allows for an alternative but equivalent definition of that of loop-erased section:

**Remark** (A recursive procedure for constructing $\text{LES}(\omega)$). Let $G$ be a digraph and consider $\omega = w_0 \dots w_\ell \in W(G)$. The set $\text{LES}(\omega)$ of loop-erased sections of $\omega = w_0 \dots w_\ell$ is constructed recursively as follows. Initialize with $\text{LES}(\omega) = \emptyset$. Then for $k \in \{1, \dots, \ell-1\}$, if $w_{k+1} \in V(\text{LEW}_k(\omega))$ denote $k'$, the greatest integer such that $0 \le k' \le k$ and $w_{k'} = w_{k+1}$. If $k'$ exists, then:

1. add the closed walk $\omega^{k',k+1} = w_{k'} \dots w_{k+1}$ to $\text{LES}(\omega)$;

2. if there exists $\omega^{k'',k'} \in \text{LES}(\omega)$, add the closed walk $\omega^{k'',k+1} = w_{k''} \dots w_{k+1}$ to $\text{LES}(\omega)$ as well.

While equivalent to Definition 2.0.1, the above formulation is more formal in flavor and recursive in nature, thus better suited to algorithm designs and easier to wield in proofs.

## 2.1 Co-preLie co-algebra of walks

From a walk $\omega$, Lawler's process, once terminated, produces a set of erased simple cycles and one self-avoiding skeleton (possibly trivial). It is therefore natural to seek a co-product which to the walk $\omega$ would associate a sum over erased sections $\omega^{k,k'}$ and associated remainders $\omega_{k,k'}$, so that $\omega$ could be obtained back from these through grafting of the former onto the latter. The 'grafting' product appropriate to that end, known as *nesting* is presented in §2.2 is permutative non-associative reflecting Lawler's process' chronological constraints. It is difficult to maintain any form of compatibility with nesting via such an indiscriminate procedure as cutting out all loop-erased sections however, as not all pairs $(\omega^{k,k'}, \omega_{k,k'})$ can be consistently grafted back to form the original walk; and when grafting is possible, it may be so in more than one way. These problems arise from certain ladders and all corollas, respectively.

Consider first an instance of the former, $\omega = 1233231$, which is a ladder in the sense that the self-loop 33 is attached 'on top of' cycle 232, itself attached to the 'base' triangle 1231. Here $\omega^{2,3} = 33$ is a valid loop-erased section of $\omega$, yet can be grafted back onto $\omega_{2,3}$ in two distinct ways: one producing $\omega$ and the other yielding the walk $\omega' = 1232331$. Remark how in $\omega'$, the self-loop 33 occurs one level below its original location in $\omega$ since it is now attached directly to the 'base' triangle 1231. Algebraically such instances correspond to cases where the nesting product fails to be associative. Second, for the issue with corollas, i.e. bouquets of closed walks with the same root, consider e.g. $\omega = 12131$. Here both $121, 131 \in \text{LES}(\omega)$; yet cutting e.g. $\omega^{0,2} = 121$ and grafting it back onto $\omega_{0,2} = 131$ either gives back the walk $\omega = 12131$ or the completely different one $\omega' = 13121$. Algebraically, these instances translate into cases where the nesting product fails to be commutative.

### Admissible cuts

To resolve the difficulties mentioned above, which become extensive when taken together in arbitrary long walks, we must refine the set of loop-erased sections that can be cut out of the original walk by the co-product. Here, as earlier, the major hurdle is due to the chronological constraints inherent to Lawler's process. Because of this, special attention must be paid to erased sections that appear within longer erased sections, the latter providing the temporal context of the former:

**Definition 2.1.1** (Temporal context of an erased section). Let $G$ be a digraph, $\omega \in W(G)$ and $\omega^{k,k'} \in \text{LES}(\omega)$. We denote $\text{LES}(\omega)_{<_{k,k'}} \subset \text{LES}(\omega)$ the subset of loop-erased sections $\omega^{l,l'}$ which *strictly* include $\omega^{k,k'}$ as left subwalk, i.e. $l \le k < k' < l'$. Because we require $k' < l'$ strictly, $\text{LES}(\omega)_{<_{k,k'}}$ may be empty. Otherwise, we denote $\omega^{\min}_{k,k'}$ the smallest element of $\text{LES}(\omega)_{<_{k,k'}}$ for inclusion.

By construction, if $\omega^{\min}_{k,k'}$ exists, it is the tightest erased section which comprises $\omega^{k,k'}$ entirely. It provides the relevant temporal context for $\omega^{k,k'}$ since anything outside of $\omega^{\min}_{k,k'}$ creates no further chronological constraints on $\omega^{k,k'}$ beyond those on $\omega^{\min}_{k,k'}$. We can now control the loop-erased sections that a co-product may extract by admitting only those cuts which are corollas within their relevant temporal context and only if those cuts are contiguous subwalks including the last petals of the corolla:

**Definition 2.1.2** (Admissible cuts). Let $G$ be a digraph and $\omega = w_1 \dots w_\ell \in W(G)$. A non-empty loop-erased section $\omega^{k,k'} := w_k w_{k+1} \dots w_{k'} \in \text{LES}(\omega)$ is an *admissible cut* of $\omega$ when $\omega^{k,k'} \ne \omega$ and either $\omega^{l,l'} := \omega^{\min}_{k,k'}$ does not exist or $w_k$ does not appear in $w_{k'+1} \dots w_l$. The set of admissible cuts of $\omega$ is denoted $\text{AdC}(\omega)$.
The condition that for $\omega^{k,k'} \in \text{LES}(\omega)$, $\omega^{l,l'} := \omega_{k,k'}^{\min}$ either does not exist or $w_k$ does not appear in $w_{k'+1} \cdots w_{l'}$ implies that admissible cuts can only be made right to left in the walk, that is from the latest to the earliest, in reverse chronological order. The notion of admissible cut is well defined because the property of being admissible does not depend on the order in which admissible cuts are considered and removed from the original walk. In particular, if a loop-erased section is an admissible cut of an admissible cut of a walk or of its remainder, then it is an admissible of that walk and vice-versa. This is significant because it indicates that, in spite of the strong chronological constraints created by Lawler's process, cutting out admissible cuts does not alter the other cuts relevant temporal context and thence, their admissibility:

**Proposition 2.1.1.** Let $G$ be a digraph and $\omega \in W(G)$.

**Case 1.** If $k < k' < l < l'$ or $l < l' < k < k'$ then,

$$ \omega^{k,k'} \in \text{AdC}(\omega) \text{ and } \omega^{l,l'} \in \text{AdC}(\omega_{k,k'}) \iff \omega^{l,l'} \in \text{AdC}(\omega) \text{ and } \omega^{k,k'} \in \text{AdC}(\omega_{l,l'}). $$

**Case 2.** If $k < l < l' \le k'$ then,

$$ \omega^{k,k'} \in \text{AdC}(\omega) \text{ and } \omega^{l,l'} \in \text{AdC}(\omega^{k,k'}) \iff \omega^{l,l'} \in \text{AdC}(\omega) \text{ and } \omega_{l,l'}^{k,k'} \in \text{AdC}(\omega_{l,l'}). $$

This Proposition is the central result on which all other proofs pertaining to the co-structure are built. The rather long proof of the Proposition, given in [57], works entirely by contradiction, treating every possible configuration of loops in a walk on a case-by-case basis. It is close in spirit to the original proof of the path-sum theorem [76], which relies on a case-by-case analysis of the configuration of loops in walks in conjunction with combinatorial arguments.

## Co-product definition

With the notion of admissible cut, we may now formally define the co-product associated to Lawler's process, by mapping a walk to a sum over all its admissible cuts tensored with their remainders:

**Definition 2.1.3 (Co-product).** Let $G$ be a digraph. The co-product associated to Lawler's process is the linear map $\Delta_{\text{CP}}$ defined by,

$$ \Delta_{\text{CP}} : \left\{ \begin{array}{rcl} \mathcal{W}(G) & \to & \mathcal{W}(G) \otimes \mathcal{W}(G) \\ \omega & \mapsto & \Delta_{\text{CP}}(\omega) = \sum_{\omega^c \in \text{AdC}(\omega)} \omega_c \otimes \omega^c. \end{array} \right. $$

An essential property of this co-product is that a walk is primitive for it if and only if it is a simple path or a simple cycle. This is the ultimate algebraic explanation for the role these objects play in the path-sum theorem.

**Proposition 2.1.2.** Let $G$ be a digraph and $\omega \in W(G)$. Then, $\Delta_{\text{CP}}(\omega) = 0 \iff \omega \in \text{SAW}(G) \cup \text{SAP}(G)$.

## Co-preLie property

Having established the definition of the co-product associated to the Lawler process and identified its primitive walks, we now turn to the co-algebraic structure it gives to the walk vector space $W(G)$. Recall that:

**Definition 2.1.4.** A co-preLie co-algebra is a couple $(\mathcal{V}, \Delta)$ where $\mathcal{V}$ is a vector space and $\Delta : \mathcal{V} \to \mathcal{V} \otimes \mathcal{V}$ is a linear map such that for any $v \in \mathcal{V}$ the following relation is satisfied

$$ (\Delta \otimes \text{Id} - \text{Id} \otimes \Delta) \circ \Delta(v) = (\text{Id} \otimes \tau) \circ (\Delta \otimes \text{Id} - \text{Id} \otimes \Delta) \circ \Delta(v) $$

where $\text{Id}$ is the identity map and $\tau$ is the twisting linear map, $\tau : \mathcal{V} \otimes \mathcal{V} \to \mathcal{V} \otimes \mathcal{V}$, $\tau(u \otimes v) = v \otimes u$.

**Theorem 2.1.1.** The vector space $W(G)$, equipped with the coproduct $\Delta_{\text{CP}}$, is a co-preLie (but not co-unital) co-algebra.

There exists two proofs of this result. The first is a direct verification of the co-preLie property based on the properties of admissible cuts. The second proof obtains the theorem as a corollary of the Hopf structure on the tensor algebra generated by $W(G)$ via a brace coalgebra construction. Note that the Hopf structure on the tensor algebra is a corollary of the properties of admissible cuts, in particular Proposition 2.1.1.
<figure><img src="image_4.png" /><figcaption>Figure 2.1: An example of nesting: the walk $w = \alpha_1\alpha_2\alpha_3\alpha_4\alpha_5\alpha_6\alpha_7\alpha_8\alpha_9\alpha_7\alpha_4$ is obtained upon inserting the triangle $c_1 = \alpha_7\alpha_8\alpha_9\alpha_7$ into the square $c_2 = \alpha_4\alpha_5\alpha_6\alpha_7\alpha_4$ and then into the simple path $p = \alpha_1\alpha_2\alpha_3\alpha_4$, that is $\omega = p \odot (c_2 \odot c_1)$.</figcaption></figure>

## Chronological structure of walks

Lemma 2.0.1 and Proposition 2.1.1 strongly suggest that any walk on any graph is chronologically equivalent to a tree where the root node is the self-avoiding skeleton of the walk and each non-root node stands for a simple cycle, see Theorem 2.1.2 below. In that tree, Lawler’s procedure erases nodes from the leaves down to the root and operates on the branches from left to right (or more precisely along the direction given to time). That is, time totally orders the walk’s tree structure.

**Definition 2.1.5** (Time-ordering of the admissible cuts). Let $G$ be a digraph and $\omega \in W(G)$. Assuming that $\mathrm{AdC}(\omega) \neq \emptyset$ we define the relation $\leq_{\otimes}$ on $\mathrm{AdC}(\omega)$ as follows:

$$ \omega^{k,k'}, \omega^{l,l'} \in \mathrm{AdC}(\omega) : \omega^{k,k'} \leq_{\otimes} \omega^{l,l'} \iff l \leq k < k' \leq l' \text{ or } k < k' < l < l'. $$

That is, $\omega^{k,k'} \leq_{\otimes} \omega^{l,l'}$ if and only if either $\omega^{k,k'}$ is erased prior to $\omega^{l,l'}$ or, if both are erased simultaneously, $\omega^{k,k'}$ began after $\omega^{l,l'}$.

With this definition, the set of admissible cuts of a walk is totally ordered chronologically. As a consequence, all walks are temporal trees in the following sense. For a walk $\omega$ one can always construct an unique tree $t(\omega)$ whose nodes and leaves represent simple cycles eliminated by Lawler’s process from $\omega$. Two nodes are linked in the tree if and only if the child node is a cycle that is found inside of the parent node in the walk. The theorem makes this intuitive presentation rigorous:

**Theorem 2.1.2.** Let $G$ be a digraph and $\omega \in W(G)$. Then $\omega$ has the temporal-structure of a tree $t(\omega)$ whose nodes are totally ordered by $\leq_{\otimes}$ according to a reverse depth-first order.

Although the tree $t(\omega)$ to which a walk is mapped depends on the walk $\omega$, an universal tree can be constructed for all walks of a given digraph $G$. Considering only the trees $t(\kappa)$ obtained from walks with no repeated sections produces a finite number of structurally distinct trees from all walks on $G$. These trees can be ordered partially by inclusion and the resulting poset always admits an unique maximum. This maximum tree is of paramount importance to $G$: it is one of the few invariants of its hike monoid (see Chapter 5 and §8.3) and dictates the shape of the branched continued fraction in the path-sum theorem §2.4.

## 2.2 Non-associative permutative algebra of walks

### Nesting product

Nesting was first introduced in [76] during my PhD thesis. The motivation for doing so was solely combinatorial: the definition aimed at ensuring the existence and uniqueness of the factorisation of walks into products of “cycle-irreducible” objects. The nesting product crafted in consequence *turns out* to be closely related to Lawler’s process in that the factors of a walk per nesting are exactly the admissible loop-erased sections of that walk. This is remarkable since nesting was devised independently from Lawler’s ideas, implying that existence and uniqueness of walk factorisations is necessarily tied with the chronological ordering of the edges of a walk. In this section we give an updated, streamlined definition of nesting. After detailing its relation with the co-preLie coproduct introduced earlier, we discuss results on the NAP-co-preLie bialgebraic structure on walks.

Let $\omega = w_1w_2\cdots w_n$ be a walk. By the loop-erasing procedure, we obtain that $w$ can be seen as a self-avoiding skeleton on top of which are attached closed-walks (loops) $\ell_1, \dots, \ell_n$, i.e.

$$ \omega = \omega_1 \ell_1 \omega_2 \ell_2 \cdots \ell_n \omega_{n+1} \qquad (2.2.1) $$

where $\omega_j$ are (possibly trivial) self-avoiding vertex sequences. If $\omega$ is open we denote $\mathrm{rad}(\omega) := \omega_1\omega_2\cdots\omega_{n+1}$ its self-avoiding skeleton. If instead $\omega$ is closed, we define $\mathrm{rad}(\omega) := \mathrm{rad}(\omega')$ to be the self-avoiding skeleton of the open walk $\omega' := w_1w_2\cdots w_{n-1}$ obtained from $\omega$ by removing its last vertex.
<figure><img src="image_5.png" /><figcaption>Figure 2.2: Illustration of the inherent non-associativity of walk construction from Lawler’s erased loops reflected in nesting’s lack of associativity. Left: walk $\omega = (x \odot y) \odot z$. Right: walk $\omega = x \odot (y \odot z)$.</figcaption></figure>

The decomposition expounded in Eq. 2.2.1 holds recursively: we can decompose all $\ell_i$ closed walks by seeing them as comprising closed walks on their radicals rad($\ell_i$). We can continue in this fashion until all closed walks invoked in the decomposition are simple cycles. The set of the simple cycles so obtained is unique, coincides with the loops erased in Lawler’s procedure. Furthermore there is a unique way to build $w$ back from these loops [76].

**Definition 2.2.1 (Nesting product).** Consider $(a, b) \in W(G)^2$ with $b$ a cycle with root $\beta$. Then the couple $(a, b)$ is nestable if and only if $\beta$ occurs in rad($a$), no closed walk is already attached to $\beta$ in rad($a$) and no vertex that is visited by rad($a$) before $\beta$ is also visited by $b$.

For a nestable couple $(a, b)$, the nesting product $a \odot b$ is defined as the walk whose vertex sequence is formed by replacing the last appearance of $\beta$ in $a$ by the entire vertex sequence of $b$.

Nesting is neither commutative nor associative, however for $x, y, z \in W(G)$, if $(x, y)$ and $(x, z)$ are both nestable couples then $(a \odot b) \odot c = (a \odot c) \odot b$ since rad$(a \odot b) = \text{rad}(a)$. As a consequence:

**Proposition 2.2.1.** The nesting product is a Non-Associative Permutative (NAP) product and $(W(G), \odot)$ is a NAP-algebra

A compatibility relation between nesting and the co-preLie coproduct $\Delta_{CP}$ exists as a corollary of the properties of admissible cuts. Indeed, take $a, b \in W(G)$, if $(a, b)$ is nestable then $\{b\} \cup \text{AdC}(b) \subset \text{AdC}(a \odot b)$ and $\text{AdC}(a \odot b) = \{a, b\} \cup \text{AdC}(a) \cup \text{AdC}(b)$ since by Proposition 2.1.1 no admissible cut of $a \odot b$ may straddle over $a$ and $b$. Thus:

**Corollary 2.2.1.** Let $G$ be a digraph, $a, b \in W(G)$, then $(\Delta_{CP}, \odot)$ share the Livernet compatibility condition,

$$ \Delta_{CP} \circ \odot = \mathbf{Id}^{\otimes 2} + (\odot \otimes \mathbf{Id}) \circ (\mathbf{Id} \otimes \tau)(\Delta \otimes \mathbf{Id}) + (\odot \otimes \mathbf{Id}) \circ (\mathbf{Id} \otimes \Delta) $$

that is, using Sweedler’s notation,

$$ \Delta_{CP}(a \odot b) = a \otimes b + a \odot b_{(1)} \otimes b_{(2)} + a_{(1)} \odot b \otimes a_{(2)} $$

and $(W(G), \odot, \Delta_{CP})$ forms a NAP-co-preLie bialgebra.

## 2.3 NAP-co-preLie bialgebra

Corollary 2.2.1 states that walks obey a NAP-co-preLie bialgebra structure. They are therefore merely an example of realization of this bialgebra, albeit a natural and eminently important one. A preLie-coNAP structure, dual of that of interest here was first studied by Livernet [125] yet much general algebraic and operadic results on the structure remain to be established. The example offered by walks is currently the main stimulus for such studies. Since the problem of identifying the simple cycles and simple paths on a graph is paramount to many applications, see e.g. Chapters 5–9 and given that these are the primitive elements of the co-algebra of walks, it is natural to seek the projection map onto these elements within the bialgebraic structure. The formula for this projection was discovered by Ronco, then proven independently by Ronco and Foissy.

Let $(A, \bullet, \Delta)$ be a connected NAP-co-preLie bialgebra, that is $(A, \mu)$ is a NAP algebra, $(A, \Delta)$ a co-preLie co-algebra and $\bullet, \Delta$ share the Livernet compatibility condition. For any $n \in \mathbb{N}$ define the maps $\Delta_n : A \to A^{\otimes(n+1)}$ recursively by setting $\Delta_0 = \mathbf{Id}_A$ and, for $n \ge 1$,

$$ \Delta_n := \sum_{i=1}^{n} \left( \mathbf{Id}_A^{\otimes(i-1)} \otimes \Delta \otimes \mathbf{Id}_A^{\otimes(n-i)} \right) \circ \Delta_{n-1}. $$

Dually define the maps $\mu_n : A^{\otimes(n+1)} \to A$ by

$$ \mu_n(a_0 \otimes a_1 \otimes \dots \otimes a_n) := (\dots ((a_0 \bullet a_1) \bullet a_2) \dots) \bullet a_n. $$

**Definition 2.3.1.** Let $E : A \to A$ be the linear map $E := \sum_{n=0}^{\infty} \frac{(-1)^n}{n!} \mu_n \circ \Delta_n$.
**Theorem 2.3.1** (Foissy, Ronco, unpublished). Let $(A, \bullet, \Delta)$ be connected NAP-co-preLie bialgebra. Let $\text{Prim}(A)$ be the set of primitive elements of $A$. Then the linear map $E$ is the Eulerian idempotent of $A$, in particular it satisfies that $a \in \text{Prim}(A) \iff E(a) = a$ and for any $a, b \in A$, $E(a \bullet b) = 0$.

This theorem proven by direct combinatorial arguments notably implies the existence of explicit formulas for counting simple cycles on graph. Foissy offers an alternative form for $E$ in terms of trees which, by Theorem 2.1.2, is also immediately interpretable in terms of walks. The Eulerian idempotent enters a formula for the identity map on $A$, which reconstructs all the elements of $A$ out of its primitive elements. This formula can be partially resummed leading to the path-sum theorem. This line of proof for the path-sum theorem, purely algebraic, is the subject of ongoing work. Without it the only proof of the path-sum theorem is a combinatorial one based on the uniqueness of the factorisation of walks into nesting products of erased loops.

## 2.4 Path-Sum theorem

The main interest of the existence and uniqueness of the factorisation of walks into nestings of primitive elements of its co-algebraic structure is that it allows resummations of families of terms in series of walks. These turn the series into an expression that involves only the primitive elements of the bialgebraic structure, that is the simple cycles and simple paths of $G$. This is analogous to how the fundamental theorem of arithmetic leads to the existence of Euler products for the Riemann zeta function and other totally multiplicative functions on the integers. This observation is made rigorous in Part 2 of this work where the relation with number theory is fully elucidated.

Let $\alpha$ and $\omega$ be two vertices on $G$, let $W(G)_{\alpha\omega}$ be the set of all walks on $G$ from $\alpha$ to $\omega$, $\Pi(G)_{\alpha\omega}$ be the set of simple paths from $\alpha$ to $\omega$ and $\Gamma(G)_{\alpha}$ the set of simple cycles from $\alpha$ to itself. The formal series of all walks from $\alpha$ to $\omega$ on $G$ is the formal sum over $W(G)_{\alpha\omega}$,

$$ \Sigma_{G; \alpha\omega} := \sum_{w \in W(G)_{\alpha\omega}} w. $$

In other words, the coefficient of $w$ in $\Sigma_{G; \alpha\omega}$, denoted is 1 if $w \in W(G)_{\alpha\omega}$ and 0 otherwise.

By using the fact that every open walk can be factorised into a simple path and a collection of nested cycles, we rewrite $\Sigma_{G; \alpha\omega}$ as a series over simple paths by modifying each path in the series to include all collections of cycles that can be nested off the vertices it visits. To preserve the vertex-edge notation of walks, we implement this modification by replacing each vertex $\alpha$ in a simple path by a 'dressed vertex' $(\alpha)'_G$ defined to represent the characteristic series of all cycles that can be nested off $\alpha$ on $G$:

$$ (\alpha)'_G := \sum_{c \in W(G)_{\alpha\alpha}} c = \Sigma_{G; \alpha\alpha}. \qquad (2.4.1) $$

We rewrite this characteristic series as a series over simple cycles $\gamma \in \Gamma(G)_\alpha$ by replacing each vertex $\mu$ visited by a simple cycle $\gamma$ by a dressed vertex representing the characteristic series of all the cycles that can be nested off $\mu$ on the appropriate subgraph of $G$. Applying this procedure recursively yields a representation of the formal series $\Sigma_{G; \alpha\omega}$ as a branched continued fraction of finite depth and breadth which only involves simple paths and simple cycles of $G$, the path-sum.

To present the path-sum theorem we need the vertex-edge notation for walks, where e.g. walk $w = \alpha\beta\omega$ is written $w = (\alpha)(\alpha\beta)(\beta)(\beta\omega)(\omega)$. Then:

**Theorem 2.4.1** (Path-sum theorem). The formal series over all walks from $\alpha$ to $\omega$ on $G$ has the following expression,

$$ \Sigma_{G; \alpha\omega} = \sum_{\Pi(G)_{\alpha\omega}} (\alpha)'_G (\alpha\nu_1) (\nu_1)'_{G\setminus\{\alpha\}} \cdots (\nu_{\ell(p)-1}\omega) (\omega)'_{G\setminus\{\alpha, \nu_1, \dots, \nu_{\ell(p)-1}\}}, \qquad (2.4.2a) $$

where $p = (\alpha\nu_1 \cdots \nu_{\ell(p)-1}\omega)$ is a simple path of length $\ell(p)$ from $\alpha \equiv \nu_0$ to $\omega \equiv \nu_{\ell(p)}$, and $(\alpha)'_G$ denotes the dressed vertex $\alpha$ on $G$, defined as the formal series of all cycles off $\alpha$ on $G$ and given explicitly by

$$ (\alpha)'_G = \left[ (\alpha) - \sum_{\gamma \in \Gamma(G)_\alpha} (\alpha) (\alpha\mu_1) (\mu_1)'_{G\setminus\{\alpha\}} (\mu_1\mu_2) \cdots (\mu_{\ell(\gamma)-1})'_{G\setminus\{\alpha, \mu_1, \dots, \mu_{\ell(\gamma)-1}\}} (\mu_{\ell(\gamma)-1}\alpha)(\alpha) \right]^{-1}, \qquad (2.4.2b) $$

with $\gamma = (\alpha\mu_1 \cdots \mu_{\ell(\gamma)-1}\alpha)$ a simple cycle of length $\ell(\gamma)$ off $\alpha$.

The formal series $\Sigma_{G; \alpha\omega}$ is expressed recursively in terms of formal series on subgraphs of $G$. We term these formal series the dressed vertices, and denote them by e.g. $(\mu_j)'_{G\setminus\{\alpha, \mu_1, \dots, \mu_{j-1}\}}$. These subseries are
in turn obtained through Eq. (2.4.2b), but on the subgraphs of G (e.g. $G\{\alpha, \mu_1, \dots, \mu_{j-1}\}$ in the case of $(\mu_j)'_{G\{\alpha, \mu_1, \dots, \mu_{j-1}\}}$). The recursion stops when vertex $\mu_j$ has no neighbour on this subgraph. In this case the dressed vertex is given by

$$ (\mu_j)'_{G\{\alpha, \mu_1, \dots, \mu_{j-1}\}} = \sum_{n \ge 0} (\mu_j \mu_j)^n = \begin{cases} [(\mu_j) - (\mu_j \mu_j)]^{-1} & \text{if the loop } (\mu_j \mu_j) \text{ exists,} \\ (\mu_j) & \text{otherwise,} \end{cases} \qquad (2.4.3) $$

where $(\mu_j)$ is the trivial walk off the vertex $\mu_j$. Note that this trivial walk is a local unit in the $(W(G), \odot)$ algebra. The recursive nature of Eq. (2.4.2b) implies that the result of the path-sum theorem for $\Sigma_{G; \alpha\omega}$ yields a formal continued fraction involving only simple cycles and simple paths. On finite digraphs, the depth of this continued fraction is necessarily finite.

The path-sum theorem has proven to have an extremely wide range of applications, e.g. for expressing analytic matrix functions [181, 77, 63]. At first I thought this was owing to the ubiquity of graph walks throughout mathematics, physics, chemistry and, via network analysis, into computer science and biology. The problem with this view is that it contrasts with the mathematically extremely limited scope of the theorem. Indeed, recall that the fundamental result of algebraic graph theory is that powers of the ordinary adjacency matrix A of a graph G counts the walks on that graph [19], for $n \in \mathbb{N}$, $(A^n)_{i,j}$ is the number of walks of length n from vertex i to vertex j on G. This result generalises naturally to labelled adjacency matrix $W_{i,j} := w_{i,j}A_{i,j}$ where $w_{i,j}$ are formal labels (see also p. 44). Then it becomes readily apparent that, formally,

$$ \Sigma_{G; \alpha\omega} = \mathrm{Id} + W + W^2 + \dots = (\mathrm{Id} - W)^{-1}. $$

In other terms, all the path-sum theorem does is express a resolvent. Alarmingly, resolvents are only to be found in linear problems since one verifies easily that if $r = (1 - a)^{-1}$ then r solves the linear equation $a.r = r - 1$. Having paid attention to that, it becomes apparent that the theorem's sole reason for being useful is that resolvents appear far more often than could reasonably be expected. This in turn, I now see as a manifestation of the existence of exotic umbral calculi on finite and infinite dimensional vector spaces, see the discussion in §10.1.
# Chapter 3

## Implications and applications to differential calculus

The most promising field of applications of the path-sum theorem is that of differential equations. While re-investigating more closely the original [77] path-sum formulation of the ordinary matrix-exponential, it appeared that a path-sum expression for a more general object, the *time-ordered exponential* was possible. This object arises in coupled systems of linear ordinary differential equations with non-constant coefficients (also known as non-autonomous systems of ODEs),

$$ \frac{d}{dt} U(t) = A(t)U(t), \quad U(0) = \mathbf{I}, $$

where $A$ is a matrix of smooth functions of time and $U$, known as the evolution operator, is the solution that is sought after. Under the assumption that $A$ commutes with itself at all times, i.e. $A(t')A(t) - A(t)A(t') = \mathbf{0}$ for all $t', t \in I$, then $U$ is the ordinary matrix exponential of $A$, $U(t) = \exp\left(\int_0^t A(\tau) d\tau\right)$. In general, however, the evolution operator has no known explicit form in terms of $A$ and is rather expressed as

$$ U(t) = \mathfrak{T} \exp \left( \int_{0}^{t} A(\tau) d\tau \right). $$

Here $\mathfrak{T}$ represents the time-ordering (or path-ordering) operator defined as

$$ \mathfrak{T}[A(t')A(t)] = \begin{cases} A(t')A(t), & \text{if } t' \ge t, \\ A(t)A(t'), & \text{otherwise.} \end{cases} $$

The action of this operator transforms the series representation of the ordinary matrix-exponential function into a matrix-valued Peano-Baker series (also known as Dyson series [47]) involving an infinite sum of nested integrals. This series is frankly unwieldy, hinders theoretical analyses of the solution, is computationally costly to evaluate and purely perturbative from a physical point of view.<fnref n="1" /> Furthermore the time-ordered exponential cannot be evaluated via direct diagonalisation of $A$. Only two analytical methods have been devised to evaluate the time-ordered exponential:<fnref n="2" />

* If the matrix $A$ is periodic, i.e. for all $t \in I$ there exists $T > 0$ with $A(t + T) = A(t)$, then 19th century Floquet theory dictates that the evolution operator takes on the form $U(t) = P(t) \exp(Ft)$. In this expression $P(t)$ a periodic time-dependent matrix and $F$ a constant matrix, both of which are determined in the Fourier domain from the solution to an *infinite-dimensional* but autonomous linear system of ODEs. This system can only be solved approximately by truncating it to a finite size either analytically [23], or otherwise via numerical procedures [86, 82]. Floquet formalism was first used by Shirley [185], who successfully applied it to study the quantum dynamics of a linearly polarised excitation in magnetic resonance. Shirley obtained low orders analytical approximations for the now celebrated Bloch-Siegert effect [24].

* The second analytical approach, discovered in 1954, is that of Magnus series [129]. It posits that the evolution operator can be expressed as an ordinary exponential $U = \exp\Omega(t)$ and gives an expression for $\Omega(t)$ as

<footnote n="1">¹That is, the terms of the series involve powers of a parameter that should be small to guarantee rapid enough convergence for an acceptable approximation to be reached by truncation of the series in applications.</footnote>

<footnote n="2">²A plethora of numerical methods have been devised: (i) Fer and Magnus-Floquet hybrids [187, 133], (ii) Zassenhaus and Suzuki-Trotter propagator approximations [30, 45, 140]; (iii) Flow-equation approaches [197]; (iv) Average Hamiltonian Theory [199, 89]; and many more case-specific methods. The expansions presented above all suffer from various drawbacks including: the divergence of the representation at long times; the perturbative nature of the numerical or theoretical approach; the non-avoidable propagation of errors at long time; the failure to find exact solutions even in 2 × 2 cases.</footnote>
an increasingly intricate infinite series of nested commutators of A with itself at different times. This series is so complicated to write explicitly that it is rarely calculated beyond its first two terms. Most importantly, Magnus expansion suffers from severe and incurable divergences as already mentioned by Magnus [129] and Fel'dman [53], now well studied [34, 136, 99, 176]. Ultimately, these result from the inherently faulty assumption that $\Omega(t)$ even exists, see Example 3.2.1 page 27 for an explicit case where U presents a negative real eigenvalue. While not restricted to periodic matrices it is well known that Magnus expansion suffers from a further two limitations for such matrices: the stroboscopic detection of events, and the impossibility to take into account more than one characteristic period. In spite of all this, Magnus series are very much in use nowadays in quantum physics [23], notably because they guarantee the unitarity of the approximated solutions even when the series diverges (!).

The observation, first made in 2015 [67] that the time-ordered exponential is amenable to a path-sum formulation is thus striking on two accounts:

1) It constitutes the only exact expression of the time-ordered exponential that involves a finite number of (★-)operations;

2) It both necessitates and implies the existence of a novel product, the ★-product, which induces a Fréchet-Lie group structure on distributions [175].

The existence of formal expressions for the solution of differential problems in ★-algebras opens the floodgates to classical linear algebra methods for solving problems in the differential context. This is exemplified by the ★-Lanczos algorithm reported here in §3.3 and consequent theorem regarding the tridiagonalization of systems of non-autonomous ODEs. These results are among the unexpected implications of the path-sum theorem, yet by no mean corollaries of it.

**Remark** (Fiber bundles and space-time graphs). Before we delve in the detailed treatment of non-autonomous ODEs, let us sketch the reason why there should even exist a path-sum formulation of time-ordered exponentials.

Consider a graph G whose edges are endowed with smooth time-dependent weight functions of $C^\infty[I]$. This may be encoded by a time-dependent adjacency matrix $A(t)$, $t \in I$. At a given time $t$, some of these weights may be 0 so the graph $G(t)$ so-obtained may be topologically different from another $G(t')$ for $t' \neq t$. This view is rather inconvenient because one cannot decide whether a walk passing from a vertex $i$ to a vertex $j$ from time $t$ to time $t'$ exists at all, unless one knows the continuous family of graphs $\{G(\tau)\}_{\tau \in [t,t']}$ It is thus more convenient to encode this family into a single larger, semi-continuous graph-like object $\mathfrak{G}$ so that time slicing $\mathfrak{G}$ at time $t$ recovers $G(t)$ while space slicing $\mathfrak{G}$ at edge $w_{i \to j}$ yields the time-dependent weight function for this edge. Now we verify that if the edge set of $\mathfrak{G}$ is just the trivial fiber-bundle $E \times C^\infty[I]^{|E|}$ then $A(t)$ commutes with itself at all times and the sum of all space-time walks on $\mathfrak{G}$ is the ordinary matrix exponential of $A(t)$.<fnref n="3" /> If instead the edge set of $\mathfrak{G}$ is a nontrivial fiber bundle, then $A$ fails to commute with itself at different times and the sum of all space-time walks on $\mathfrak{G}$ is now a time-ordered exponential. As a sum of walks it must admit a path-sum representation as per the path-sum theorem. In this view the ★-product defined below appears naturally from the continuous extension of the matrix-product, continuously generating the walking process along the time dimension, just as the ordinary matrix product generates discrete walk steps along the space dimension.

## 3.1 The ★-product

To present the path-sum formulation of the time-ordered exponential, we must first introduce the algebraic structure in which this formulation makes sense and the product on which it is based. While we introduced the ★-product in [67, 69], we here follow the more recent presentation by Ryckebusch2023 [175].

### Definition

Let $I$ be a compact subspace of $\mathbb{R}$. We write $C^\infty(I^2)$ for the set of smooth functions on an open neighborhood of $I^2$ and $C^\infty(I)$ for the set of smooth functions on an open neighborhood of $I$. Let $\mathcal{D}$ be the set of distributions of the form

$$d(t', t) = \tilde{d}(t', t)\Theta(t' - t) + \sum_{i=0}^{+\infty} \tilde{d}_i(t', t)\delta^{(i)}(t' - t),$$

where $\tilde{d}, \tilde{d}_i \in C^\infty(I^2)$ are complex valued functions defined on $I^2$, $\Theta(\cdot)$ stands for the Heaviside theta function with the convention $\Theta(0) = 1$ and $\delta^{(i)}(t' - t)$ is the $i$th Dirac delta derivative in the sense of distributions evaluated in $t' - t$. Then we have:

<footnote n="3">³YES: the sum of all walks on a purely discrete graph is a resolvent $\sum_n A^n = (\mathbf{Id}-A)^{-1}$ but the sum of walks on a semi-continuous graph is a matrix exponential! This is seemingly little known.</footnote>
**Theorem 3.1.1** (Ryckebusch 2023 [175]). The $\star_I$-product, defined by

$$ (f \star_I g)(x, y) := \int_I f(x, \tau)g(\tau, y)d\tau, \qquad (3.1.1) $$

is well defined on the weak closure $\overline{C}^\infty(I^2)$ of $C^\infty(I^2)$, that is for any $f,g \in \overline{C}^\infty(I^2)$, $f \star_I g$ exists and is in $\overline{C}^\infty(I^2)$. Furthermore the $\star_I$-product is associative over $\overline{C}^\infty(I^2)$, $D \subset \overline{C}^\infty(I^2)$ and for $t', t \in I^2$, $1_* := \delta(t' - t)$ is the unit of the $\star_I$ product. The $\star_\mathbb{R}$-product is a well defined associative product on $D$, that is for any $d, e \in D$,

$$ (d \star_{\mathbb{R}} e)(t', t) := \int_{-\infty}^{\infty} d(t', \tau)e(\tau, t)d\tau, $$

is a distribution of $D$.

From now-on the $\star_\mathbb{R}$ product will be generically designated as ‘the $\star$-product’. The above theorem is non-trivial and follows from a lengthy and careful sequential approach to distributions in the spirit of Mikusiski’s pioneering work [143]. From now on, we omit the $t' - t$ arguments of the Heaviside $\Theta$ functions and Dirac deltas $\delta^{(i)}$ to alleviate the equations.

**Remark.** The $\star$-product extends directly to time-dependent matrices $A, B \in D^{n \times n}$ in the natural way,

$$ (A \star B)_{i,j}(t', t) = \sum_{k=1}^{n} (B_{i,k} \star A_{k,j})(t', t) = \int_{-\infty}^{\infty} (A(t', \tau)B(\tau, t))_{i,j} d\tau. $$

## Properties

The $\star$-product is non-commutative, induces the convolution $f \star g \equiv f \star g$ whenever both $f$ and $g$ depend only on the difference between their variables; but also the Volterra compositions of first and second kind; the pointwise product on $C^\infty(I^2)$; the Schwartz bracket on $\overline{C}^\infty(I)$ and can be seen as ‘continuum version of the matrix product.

Basic $\star$-multiplicative identities on $D$ include $\Theta \star \delta' = \delta$, which indicates that $\delta^{\star-1} = \Theta$ and equivalently $\Theta^{\star-1} = \delta' = \delta^{(1)}$ since $\delta \equiv 1_*$ acts as the unit of the $\star$-product. In addition since, for $n \in \mathbb{N} \setminus \{0\}$, $\Theta^{\star n} = (t' - t)^{n-1}/(n-1)!$ $\Theta$ then $\Theta^{\star n} \propto \Theta$, implying that all negative $\star$-powers of $\delta'$ are included in $D$. Conveniently, differentiation coincide with $\star$-powers of $\delta'$, $\delta^{(j)} = (\delta')^{\star j}$. In view of our earlier remarks this in fact holds for $j \in \mathbb{Z}$, understanding that $\delta^{(-|j|)} := (\delta')^{*-|j|} = \Theta^{*-|j|}$. This means that iterated integration also coincides with $\star$-powers. We may therefore summarily write, for all $i, j \in \mathbb{Z}$,

$$ \delta^{(i)} \star \delta^{(j)} = (\delta')^{\star i} \star (\delta')^{\star j} = \delta^{\star i+j} = \delta^{(i+j)}. $$

These results follow easily from the reduction of the $\star$-product to convolutions $\star$ but there are more general $\star$-multiplicative identities which cannot be obtained this way. Let $\tilde{f} \in C^\infty(I^2)$, and denote $\tilde{f}^{(k,l)}(\tau, \rho)$ the $k$th $x$-derivative and $l$th $y$-derivative of $\tilde{f}$ evaluated at $x = \tau$, $y = \rho$ with the conventions that $k = 0$ or $l = 0$ means no derivative is taken and $k = -1$ or $l = -1$ denotes integration. By associativity of the $\star$-product, $(\delta^{(k)} \star \tilde{f}) \star \delta^{(l)} = \delta^{(k)} \star (\tilde{f} \star \delta^{(l)}) = \tilde{f}^{(k,l)}$ is well defined. We showed in [69] that, for any $i, j \ge -1$,

$$ \delta^{(j)} \star (\tilde{f}(t', t)\delta^{(i)}) = \tilde{f}^{(j,0)}(t', t)\delta^{(i)} + \sum_{k=1}^{j} \tilde{f}^{(j-k,0)}(t, t)\delta^{(i+k)}, \qquad (3.1.2a) $$

$$ (\tilde{f}(t', t)\delta^{(i)}) \star \delta^{(j)} = (-1)^j \tilde{f}^{(0,j)}(t', t)\delta^{(i)} + \sum_{k=1}^{j} (-1)^{j+k} \tilde{f}^{(0,j-k)}(t', t')\delta^{(i+k)}. \qquad (3.1.2b) $$

Observe that the smooth function’s partial derivatives are evaluated in $(t, t)$ and $(t', t')$ in all but the first term.

**Remark.** The representation of $\star$-products between elements of $D$ can be deceptive. Consider for example $\tilde{f} \in C^\infty(I^2)$ such that all its derivative of order up to $k$ are null whenever $t' = t$, i.e. $\exists k \in \mathbb{N}: \forall j_1, j_2 \in \mathbb{N}, j_1 + j_2 \le k$, $\tilde{f}^{(j_1, j_2)}(x, x) = 0$. Then $\tilde{f}(x, y)\delta^{(k)} \equiv 0_D$ is null both as a linear functional $C^\infty(I^2) \to \mathbb{C}$ whose action is defined by the Schwartz bracket and as endomorphism of $\overline{C}^\infty(I^2)$, something which is not readily apparent from the notation alone. To further illustrate the notational difficulties consider calculating $(\tilde{f}(x, y)\delta') \star (\tilde{g}(x, y)\delta')$. Relying on the $\star$-action of the leftmost $\delta'$, Eq. (3.1.2a) yields

$$ (\tilde{f}(x, y)\delta') \star (\tilde{g}(x, y)\delta') = \tilde{f}^{(0,1)}(x, x)\tilde{g}(x, y)\delta' + \tilde{f}(x, x)\tilde{g}^{(1,0)}(x, y)\delta' + \tilde{f}(x, x)\tilde{g}(x, y)\delta^{(2)}. \qquad (3.1.3) $$

But we could equally well calculate this relying on the $\star$-action of the rightmost $\delta'$. Then Eq. (3.1.2b) gives

$$ (\tilde{f}(x, y)\delta') \star (\tilde{g}(x, y)\delta') = -\tilde{f}^{(0,1)}(x, y)\tilde{g}(x, y)\delta' - \tilde{f}(x, y)\tilde{g}^{(1,0)}(x, y)\delta' + \tilde{f}(x, y)\tilde{g}(x, y)\delta^{(2)}. \qquad (3.1.4) $$
In spite of appearances the two results are equal: their action as linear functionals $C^\infty(I^2) \to \mathbb{C}$ are the same on any test function. This is because by [180, eqs. II,1; 57, p. 35] for any smooth function $\tilde{h}$ and $k \in \mathbb{N}$,

$$ \tilde{h}(x)\delta^{(k)} = (-1)^k (\tilde{h}(y)\delta)^{(0,k)}, \quad \tilde{h}(y)\delta^{(k)} = (\tilde{h}(x)\delta)^{(k,0)}, $$

so defining e.g. $\tilde{h}(\tau) := \tilde{f}^{(0,1)}(x, \tau)\tilde{g}(\tau, y)$ and so on, one turns Eq. (3.1.4) into Eq. (3.1.3).

## On ★-inverses

The calculation of ★-inverses of functions $f(t', t)$ carries the gist of the difficulty inherent in obtaining explicit expressions for time-ordered exponentials. Contrary to ★-resolvents of elements of $\mathrm{Sm}_\Theta$ (see p. 24) which can be always represented by the Peano-Baker series; given an arbitrary smooth function $\tilde{f}(t', t)$ and barring any further assumption, the ★-inverse of $f(t', t) = \tilde{f}(t', t)\Theta(t' - t)$ cannot in general be given explicitly. In this section, we show that the ★-inverse $f^{★-1}$ exists and is accessible explicitly provided that $\tilde{f}(t', t)$ is a separable function that is smooth in $t, t'$ and not identically null. A function $\tilde{f}(t', t)$ is separable if and only if there exist ordinary functions $\tilde{a}_i$ and $\tilde{b}_i$ with

$$ \tilde{f}(t', t) = \sum_{i=1}^{k} \tilde{a}_i(t') \tilde{b}_i(t). $$

This result is of central importance because it implies two highly non-trivial results: i) that the ★-invertible elements of $D$ are dense in $D$; and ii) that the ★-Lanczos algorithm for the tridiagonalization of differential systems suffers from no further breakdowns than the classical non-Hermitian Lanczos algorithm.

We begin our results on ★-inverses with that of functions of a single time variable times a Heaviside function. These inverses are found by direct calculation of the solution of degree-one ordinary differential equations:

**Proposition 3.1.1.** Let $a(t', t) := \tilde{a}(t')\Theta(t' - t)$ and $b(t', t) := \tilde{b}(t)\Theta(t' - t)$ so that $\tilde{a}$ and $\tilde{b}$ are differentiable, and not identically null over $I$. Then

$$ a^{★-1}(t', t) = \frac{\partial}{\partial t'} \frac{\delta(t' - t)}{\tilde{a}(t')}, \quad b^{★-1}(t', t) = -\frac{\partial}{\partial t} \frac{\delta(t' - t)}{\tilde{b}(t)}. $$

The method employed in the proof of the proposition above generalises to polynomials in at least one time variable, here taken to be $t'$. An analogous result can be given for functions that are polynomials in $t$.

**Proposition 3.1.2.** Let $p(t', t) = \tilde{p}(t', t)\Theta(t' - t)$ be so that $\tilde{p}(t', t)$ is a polynomial of degree $k \ge 1$ in $t'$ and is smooth in $t$. If $p(t, t)$ is not identically null over $I$, then

$$ p(t', t)^{★-1} = x(t', t) \star \delta^{(k+1)}(t' - t), $$

where $x(t', t) = \tilde{x}(t', t)\Theta(t' - t)$ is the solution of the linear homogeneous ordinary differential equation in $t \sum_{j=0}^{k} (-1)^j \tilde{p}^{(k-j,0)}(t, t) \tilde{x}^{(0,j)}(t', t) = 0$, with the boundary conditions $\tilde{x}^{(0,k-1)}(t', t') = (-1)^{k-1} \tilde{p}(t', t')^{-1}$, $\tilde{x}^{(0,k-2)}(t', t') = 0, \dots, \tilde{x}(t', t') = 0$.

**Remark.** The condition $\tilde{p}(t, t) \ne 0$ for $p^{★-1}(t', t)$ to exist shortens the proof and facilitates the presentation of the results but can be relaxed. What is truly necessary is that $p(t', t)$ itself must not be identically zero on $I^2$.

A technique similar to the one used in the proof of Proposition 3.1.2 can be applied to a more general class of functions. For instance, whenever differentiating leads to an expression like $\delta^{(k)} \star f(t', t) = \tilde{h}(t)f(t', t) + g(t', t)$, the expression can be rewritten as

$$ (\delta^{(k)} - \tilde{h}(t)\delta) \star f(t', t) = g(t', t). $$

Then we can go on with a further combination of differentiations until there is no Heaviside function left on the right-hand side of the above equality. In particular, such a technique can be used when dealing with commonly encountered exponential or trigonometric functions.

The strategy used in the proof of Proposition 3.1.1 can be extended to give ★-inverses in the much more general case of functions which are separable and piecewise smooth in both time variables over the interval $I$. The proof is quite technical so it is omitted, see [68].

**Theorem 3.1.2.** Consider a function $f(t', t) := \tilde{f}(t', t)\Theta(t' - t)$ with $\tilde{f}(t', t) \in C^\infty[I^2]$ and so that $\tilde{f}(t, t)$ is not identically null. Assume that there exists a distribution $L(t', t) := \sum_{j=0}^{k+1} \tilde{g}_j(t')\delta^{(j)} \in D$, with $\tilde{g}_j(t') \in C^\infty[I]$ and $\tilde{g}_{k+1} \ne 0$, such that $L \star \tilde{f} = 0$. Then, if $k > 0$, the ★-inverse of $f$ is

$$ f^{★-1} = \tilde{r}_{-1}(t', t)\Theta + \sum_{m=0}^{k} \tilde{r}_{m}(t')\delta^{(m)}, $$
with the separable smooth function $\tilde{r}_{-1}(t', t) := \sum_{j=0}^{k+1} (-1)^j \tilde{y}_j^{(0,j)}(t', t)$, and the smooth functions $\tilde{r}_{m \ge 0}(t') := \sum_{j=m+1}^{k+1} (-1)^j \tilde{y}_j^{(0,j-1-m)}(t', t')$, where $\tilde{y}_j(t', t) := \tilde{x}(t', t)\tilde{g}_j(t)$ and $\tilde{x}(t', t)$ is the solution of the linear homogeneous ordinary differential equation in $t$, $\sum_{m=0}^{k} \tilde{h}_m(t)\tilde{x}^{(0,m)}(t', t) = 0$, with boundary conditions $\tilde{x}^{(0,k-1)}(t', t') = -\tilde{h}_k(t')^{-1}$, $\tilde{x}^{(0,k-2)}(t', t') = 0$, ..., $\tilde{x}(t', t') = 0$. In these expressions, $\tilde{h}_m(t)$ are smooth functions given by

$$ \tilde{h}_m(t) := \sum_{j=m+1}^{k+1} \sum_{\ell=m}^{j-1} \binom{\ell}{m} (-1)^{\ell} \tilde{f}^{(j-\ell-1,0)}(t,t) \tilde{g}_j^{(\ell-m)}(t). $$

If instead $k=0$, the $\star$-inverse of $f$ is trivially given by $f^{\star-1}(t', t) = (\tilde{g}_1(t')\tilde{f}(t', t'))^{-1}L(t', t)$.

Inverting the role of $t'$ and $t$, a completely similar theorem is proven by changing all left $\star$-multiplications by $\delta^{(j)}$ with right multiplications and vice-versa. In this situation, $\tilde{x}$ satisfies a linear homogeneous ordinary differential equation in $t'$, and the boundary conditions involve the variable $t$.

**Remark.** The assumption $f(t,t) = 0$ is not necessary. We can reformulate the theorem statement so that the condition is $f$ not identically zero on $I$.

The most stringent condition imposed by Theorem 3.1.2 is the existence of the differential operator $L$ with coefficients that depend only on $t'$. This condition can be made more transparent upon relating it to the class of separable functions. Let $\tilde{y}_1(t')$, ..., $\tilde{y}_{k+1}(t')$ be smooth functions of $t'$, and $\tilde{a}_1(t)$, ..., $\tilde{a}_{k+1}(t)$ be functions of $t$. If $\tilde{y}_1(t')$, ..., $\tilde{y}_{k+1}(t')$ are linearly independent, equivalently, the related Wronskian $W(\tilde{y}_1, \dots, \tilde{y}_{k+1})$ is not identically null, i.e.

$$ W(\tilde{y}_1, \dots, \tilde{y}_{k+1}) := \begin{vmatrix} \tilde{y}_1 & \tilde{y}_2 & \dots & \tilde{y}_{k+1} \\ \tilde{y}'_1 & \tilde{y}'_2 & \dots & \tilde{y}'_{k+1} \\ \vdots & \vdots & \ddots & \vdots \\ \tilde{y}^{(k)}_1 & \tilde{y}^{(k)}_2 & \dots & \tilde{y}^{(k)}_{k+1} \end{vmatrix} \neq 0, $$

then there exist $L$ as in Theorem 3.1.2 so that $L \star \tilde{f} = 0$ for every separable function

$$ \tilde{f}(t', t) = \tilde{a}_1(t)\tilde{y}_1(t') + \dots + \tilde{a}_{k+1}(t)\tilde{y}_{k+1}(t'). \quad (3.1.5) $$

Indeed, the conditions $L \star \tilde{y}_j = 0$, for $j = 1, \dots, k+1$, give the system

$$ \begin{bmatrix} \tilde{y}_1 & \tilde{y}'_1 & \dots & \tilde{y}_1^{(k)} \\ \tilde{y}_2 & \tilde{y}'_2 & \dots & \tilde{y}_2^{(k)} \\ \vdots & \vdots & \ddots & \vdots \\ \tilde{y}_{k+1} & \tilde{y}'_{k+1} & \dots & \tilde{y}_{k+1}^{(k)} \end{bmatrix} \begin{bmatrix} \tilde{g}_0 \\ \tilde{g}_1 \\ \vdots \\ \tilde{g}_k \end{bmatrix} = -\tilde{g}_{k+1} \begin{bmatrix} \tilde{y}_1^{(k+1)} \\ \tilde{y}_2^{(k+1)} \\ \vdots \\ \tilde{y}_{k+1}^{(k+1)} \end{bmatrix}, $$

whose solutions exist since the Wronskian is not identically null. In particular, at least one of the solutions has smooth coefficients. Theorem 3.1.2 thus yields the following corollary for separable functions:

**Corollary 3.1.1.** Let $f(t', t) := \tilde{f}(t', t)\Theta(t' - t)$ with $\tilde{f}(t', t) \in C^\infty[I^2]$ so that $\tilde{f}(t, t)$ is not identically null. Then $f^{\star-1}$ exists and is given as in Theorem 3.1.2.

Thanks to the previous results, we can also show that the set comprising all the separable distributions in $D$ is closed under $\star$-inversion. Here a distribution $d = \sum_{i=-1}^{N} \tilde{d}_i(x, y)\delta^{(i)} \in D$ is said to be separable if and only if all of its smooth coefficients $\tilde{d}^{(i)}$ are separable. If $N \in \mathbb{N} \cup \{+\infty\}$ is finite, we say that $d$ is separable of order $N$.

**Corollary 3.1.2.** Let $d \in D$ a separable distribution of order $N$. Then the $\star$-inverse of $d$ exists and can be written as

$$ d^{\star-1}(t', t) = \Theta(t' - t)^{\star(N+1)} \star h^{\star-1}(t', t), $$

where $h(t', t)$ is a separable, hence $\star$-invertible function. Furthermore, $d^{\star-1}$ is a separable distribution.

Before we turn to the consequences of these results, let us present one example of $\star$-inverses since these are much less commonly encountered than the $\star$-resolvents.

**Example 3.1.1.** Let us determine the $\star$-inverse of $f(t', t) = (t'^2 - t/t')\Theta(t' - t)$. Since $\tilde{f}(t', t)$ is separable, smooth in both variables, and $\tilde{f}(t, t)$ is not identically null, Theorem 3.1.2 applies immediately. Setting $L(t', t) := \tilde{g}_0(t')\delta+\tilde{g}_1(t')\delta'+\tilde{g}_2(t')\delta''$ with $\tilde{g}_0(t') := 1$, $\tilde{g}_1(t') := 0$ and $\tilde{g}_2(t') := -t'^2/2$, we have $k = 1$ and $L(t', t)\star\tilde{f}(t', t) = 0$. This leads to

$$ \tilde{h}_0(t) := 3t/2, \quad \tilde{h}_1(t) := (t^4 + t^2)/2, $$
which are the only non-identically null functions $\tilde{h}_m$. The function $\tilde{x}$ is thus the solution of

$$3\tilde{x}(t', t) + t(1 + t^2)\tilde{x}^{(0,1)}(t', t) = 0, \quad \tilde{x}(t', t') \frac{t'^2}{2}(t'^2 + 1) = -1.$$

We find

$$\tilde{x}(t', t) = - \frac{2t'(t^2 + 1)^{3/2}}{(t'^2 + 1)^{5/2}t^3}.$$

We verify that

$$x \star L \star f = \tilde{h}_1(t)\tilde{x}(t', t')\delta = \frac{t^2(t^2 + 1)}{t'^2(t'^2 + 1)}\delta = \delta,$$

indicating that indeed

$$f^{\star-1} = x \star L = \left( \frac{2t'(t^2 + 1)^{3/2}}{(t'^2 + 1)^{5/2}t^3} \Theta \right) \star \left( \delta - \frac{t'^2}{2}\delta'' \right).$$

In conjunction with the Stone-Weierstrass theorem, study of the topology of $\mathcal{D}$ as well asconvergence arguments in $\mathcal{D}$, the result on $\star$-invertibility of separable distributions leads to:

**Theorem 3.1.3** (Ryckebusch2023 [175]). *The set Inv($\mathcal{D}$) of $\star$-invertible distributions of $\mathcal{D}$ is a dense subset of $\mathcal{D}$ that is made of automorphisms of $C^\infty(I^2)$. Furthermore, $(\text{Inv}(\mathcal{D}), \star)$ is a Fréchet Lie-group.*

## Properties of $\star$-resolvents

As a case of special importance in applications is the subset $\text{Sm}_\Theta \subset \mathcal{D}$ comprising those distributions of the form $f(t', t) = \tilde{f}(t', t)\Theta(t' - t)$, $\tilde{f} \in C^\infty(I^2)$. For $f_1, f_2 \in \text{Sm}_\Theta$, the $\star$-product between $f_1, f_2$ simplifies to the Volterra composition of the second kind between $\tilde{f}_1$ and $\tilde{f}_2$,

$$(f_2 \star f_1)(t', t) = \int_{-\infty}^{\infty} \tilde{f}_2(t', \tau) \tilde{f}_1(\tau, t) \Theta(t' - \tau) \Theta(\tau - t) d\tau = \Theta(t' - t) \int_{t}^{t'} \tilde{f}_2(t', \tau) \tilde{f}_1(\tau, t) d\tau,$$

showing that $(\text{Sm}_\Theta, \star)$ is a monoid. The primary role played by $\text{Sm}_\Theta$ comes from $\star$-resolvents of elements of $\text{Sm}_\Theta$, which solve Volterra integral equations of the second kind.

Let $f \in \text{Sm}_\Theta$, the $\star$-resolvent of $f$ is the distribution $R_f := (1_\star - f)^{\star-1} \in \text{Sm}_\Theta \cup \{1_\star\}$. The $\star$-resolvent of an element of $\text{Sm}_\Theta$ is guaranteed to exist [84], as shown by the $\star$-Neumann series (also known as Peano-Baker series and Picard iteration)

$$(1_\star - f)^{\star-1} = \sum_{n=0}^{\infty} f^{\star n} = 1_\star + f + f \star f + f \star f \star f + \dots$$

Indeed, since $\tilde{f} \in C^\infty(I^2)$, $\|\tilde{f}|_{I^2}\|_\infty$ exists on any compact $I^2$ of $\mathbb{R}^2$. Then $\|f^{\star n}\| \le \|\tilde{f}\|_\infty(t' - t)^n/n!$ and the $\star$-Neumann series converges. Let now $\tilde{R}_f \Theta := R_f - 1_\star$. Since $f \star R_f = R_f - 1_\star$, we find that $\tilde{R}_f$ solves the inhomogenous linear Volterra integral equation of the second kind with kernel $\tilde{f}$ [124, 88],

$$\tilde{f}(t', t) + \int_{t}^{t'} \tilde{f}(t', \tau) \tilde{R}_f(\tau, t) d\tau = \tilde{R}_f(t', t).$$

In general, such equations are usually solved indirectly through transformations mapping the equation into a system of non-autonomous ODEs [36, 84, 162]. The solutions of such systems are in fact themselves non-obvious since they involve time-ordered exponentials. To make matter worse this process is often circular in applications, as the ODEs are usually the equations meant to be solved! Restricting the problem to separable kernels, i.e. assuming $\tilde{f}(t', t) = \sum_{i}^{K} \tilde{a}_i(t')\tilde{b}_i(t)$, is unfortunately not sufficient. Indeed, the only available general result is for $K = 1$,

**Proposition 3.1.3.** Let $\tilde{f}(t', t) = \tilde{a}(t')\tilde{b}(t) \in C^\infty(I^2)$. Then $\tilde{R}_f = \tilde{a}(t')\tilde{b}(t)e^{\int_{t}^{t'} a(\tau)b(\tau)d\tau}$.

To the best of our knowledge, whenever $K > 1$ there is no known closed-form expression for $\star$-resolvents with separable kernels, save in a few very special cases [162]. Examples where $\tilde{R}_f$ is a higher transcendental function are known [203]. All that can be said in general is that if $\tilde{f}$ is separable then $\tilde{R}_f$ is separable as well [161]. Consequently, the best possible analytical approach is to provide a series representation of the solution with 'good' convergence and truncation properties. Taylor polynomial expansion have been proposed [182] but the $\star$-Neumann series gives a simpler starting point in this quest, as it is amenable to series acceleration [60]:
**Theorem 3.1.4** (Accelerated ★-Neumann series). Let $\tilde{f} := \sum_{i=1}^{N} \tilde{f}_i$ with $\tilde{f}_i \in C^{\infty}(I^2)$ and let $\tilde{R}_i$ be the smooth part of the ★-resolvent of $\tilde{f}_i$. Then for any permutation $\sigma \in S_N$,

$$ \tilde{R}_f = \sum_{n=0}^{\infty} T^{\star n} \star_{i=1}^{N} \tilde{R}_{\sigma(i)} - 1_{\star} = \star_{i=1}^{N} \tilde{R}_{\sigma(i)} + T \star \star_{i=1}^{N} \tilde{R}_{\sigma(i)} + T^{\star 2} \star \star_{i=1}^{N} \tilde{R}_{\sigma(i)} \dots $$

with

$$ T := 1_{\star} - \star_{i=1}^{d} \tilde{R}_{\sigma(i)} \star (1_{\star} - \tilde{f}) $$

This acceleration is particularly effective when one of the $\tilde{f}_i$ dominates over the others or when $T$ is very small for other reasons, see §4.2.

## On numerical computations of ★-products

We turn to the question of evaluating ★-products and ★-resolvents numerically. Let $I = [a, b]$ be an interval of interest for computations and let $\{t_i \in I\}_{0 \le i \le N-1}$ be the discrete values at which the numerical evaluations are to be made. For simplicity, suppose that the step $\Delta t := |t_{i+1} - t_i|$ is the same for all $0 \le i \le N-2$. This assumption is not necessary but alleviates the notation. Now for $f \in \text{Sm}_{\Theta}, \tilde{f} \in C^{\infty}(I^2)$, define a matrix $F$ with entries

$$ F_{i,j} := f(t_i, t_j) = \tilde{f}(t_i, t_j) \Theta(t_i - t_j). $$

By construction $F$ is lower triangular owing to the Heaviside step function. Now let $f, l \in \text{Sm}_{\Theta}$ with $F$ and $L$ the corresponding matrices. Once $I$ is discretised the ★-product—which is a Volterra composition of second kind on elements of $\text{Sm}_{\Theta}$—turns into an ordinary matrix product,

$$ (f \star l)(t_i, t_j) = \int_{t_j}^{t_i} \tilde{f}(t_i, \tau) \tilde{l}(\tau, t_j) d\tau \ \Theta(t_i - t_j) \rightarrow \sum_{t_j \le t_k \le t_i} \tilde{f}(t_i, t_k) \tilde{l}(t_k, t_j) \Delta t = (F \cdot L)_{i,j} \Delta t, $$

since $\lim_{\Delta t \to 0} (F \cdot L)_{i,j} \Delta t = \int_{t_j}^{t_i} \tilde{f}(t_i, \tau) \tilde{l}(\tau, t_j) d\tau$. This line of reasoning extends to ★-resolvents,

$$ (1_{\star} - f)^{-1}(t_i, t_j) = \lim_{\Delta t \to 0} \frac{1}{\Delta t} (\mathrm{Id} - \Delta t \mathrm{F})^{-1}_{i,j}, $$

which in practice would mean taking $\Delta t \ll 1$ and computing an ordinary matrix resolvent. We improve on the above by noting that these results correspond to using the rectangular rule of integration. Using the trapezoidal rule instead leads to more accurate results. To follow this rule, the usual matrix product $F \cdot L \Delta t$ representing $f \star l$ must be replaced by

$$ \frac{\Delta t}{2} (F - dF) \cdot L + \frac{\Delta t}{2} F \cdot (L - dL), $$

where $dF$ is the diagonal of $F$. Now the ★-resolvent of $f$ becomes

$$ (1_{\star} - f)^{-1}(t_i, t_j) \approx R_f := \frac{1}{\Delta t} \left( \mathrm{Id} - \Delta t \mathrm{F} + \frac{\Delta t}{2} \mathrm{dF} \right)^{-1}_{i,j}. $$

with an error scaling quadratically in $O(\Delta t^2)$. In applications the time-ordered exponential $U$ is rarely necessary at all $t', t \in I$, rather only $U(t, 0)$ is required. This means that only the first column of $R_f$ is useful. This observation is profitably exploited numerically, as it avoids the need to even compute the matrix inverse in $R_f$, rather asking for the solution of the triangular system $(\mathrm{Id} - \Delta t \mathrm{F} + \frac{\Delta t}{2} \mathrm{dF}) \cdot \vec{x} = \vec{v}$, where $v = (1, 0, 0, \dots)^{\mathrm{T}}$, which is faster and requires less memory. From there further matrix multiplications can be implemented vectorially on $\vec{x}$.

Other, better, quadrature rules than these trapezoidal one can be implemented to improve the error scaling with $\Delta t$. Simpson's quadrature for example leads to an error scaling in $O(\Delta t^3)$. Taken together these results completely bypass the ★-Neumann series, which is not desirable numerically, rather evaluating directly what the series converges to. The entire complexity inherent in the calculation of ★-resolvent is reduced to multiplying and inverting triangular, well-conditioned matrices. For the latter assertion, consider the diagonal elements of $\mathrm{Id} - \Delta t \mathrm{F} + (\Delta t/2) \mathrm{dF}$: these can be made as far from 0 as desired by tuning $\Delta t$.

While these ideas are intuitive and stem from the deeper view of the ★-product as a continuum matrix product, they are not the best way to proceed numerically. Instead it is better to express all functions to be manipulated in the basis of Legendre polynomials. In this basis, inherently discrete, a matrix representation for elements of $\text{Sm}_{\Theta}$ also arises while the error scaling, stability and speed of the computations are profoundly improved. This now represents an extensive line of research conducted mainly by S. Pozza's group with ground-breaking performances. We refer the reader to [165, 164, 166, 167, 39, 25] and references therein.
## 3.2 Path-sum representation of the solution of systems of non-autonomous ODEs

In its most general form Theorem 2.4.1 stems from a fundamental algebraic property of the set of all walks on any graph. Consequently the representation of the formal series of all walks with fixed endpoints as a finite branched continued fraction is universal and continues to hold *mutatis mutandis* for series of weighted walks as well. It is therefore sufficient to show that the solution of a problem takes the form of a resolvent to imply the existence of path-sum formulation for it. This bypasses the space-time view sketched earlier but is equivalent to it [67].

In the case of non-autonomous systems of coupled linear ordinary differential equations this is a consequence of the ⋆-product formalism. Indeed consider the system

$$ \frac{d}{dt} U(t) = \tilde{A}(t) U(t), \quad U(0) = \mathrm{Id}, $$

where $\tilde{A}$ is a time-dependent matrix whose entries are smooth functions of time over some interval $I$ of interest; and $U$, called the evolution operator, is the solution of the system. Now remark that the system above is equivalent to the two-variable system

$$ \frac{d}{dt'} \{ U(t') \} \Theta(t' - t) = A(t', t) \cdot U(t') \Theta(t' - t), \quad U(t, t) = \mathrm{Id}. $$

where we introduced $A(t', t) := \tilde{A}(t')\Theta(t' - t)$. Now define the Green's function $G := \delta' \star U\Theta = \frac{d}{dt'} \{U(t')\} \Theta(t' - t) + \mathrm{Id}_\star$ with $\mathrm{Id}_\star := \mathrm{Id} \mathrm{Id}_\star$. Observe also that $A \cdot U\Theta = A \star G$. Thus the Green's function obeys the equation

$$ G - \mathrm{Id}_\star = A \star G \Rightarrow G = (\mathrm{Id}_\star - A)^{-1}, \qquad (3.2.1) $$

while $U = \Theta \star G$. This shows that the solution of non-autonomous systems of coupled linear ODEs stems from the ⋆-resolvent of $A = \tilde{A}\Theta$. This in turn entails the existence of a path-sum formulation for $G$ and $U$.

### Path-sum formulation of the time-ordered exponential

We consider $\tilde{A}(t) \in C^\infty[I]^{n \times n}$ a time-dependent matrix whose entries are smooth on an open neighborhood containing the compact $I \subset \mathbb{R}$. As earlier we denote $A := \tilde{A}\Theta$.

Let $G = (V, E)$ be the graph with weighted adjacency matrix $A$, i.e. if there exists $t \in I$ with $A_{\beta\alpha}(t) \neq 0$, then $(\alpha, \beta) \in E$. We designate $G \setminus \{\alpha_1, \dots, \alpha_n\}$ the induced subgraph of $G$ with vertex set $V \setminus \{\alpha_1, \dots, \alpha_n\}$. Recall that $\Pi(G)_{\alpha\omega}$ and $\Gamma(G)_\alpha$ are the set of simple paths from $\alpha$ to $\omega$ on $G$ and the set of simple cycles from $\alpha$ to itself on $G$, respectively. If $G$ has finitely many vertices and edges, these two sets are finite.

**Theorem 3.2.1.** *The time-ordered exponential U of A(t) is given by*

$$ U_{\omega\alpha} = \Theta \star \sum_{p \in \Pi(G)_{\alpha\omega}} G_{G \setminus \{\alpha, \nu_1, \dots, \nu_{\ell-1}\}}; \omega \star A_{\omega\nu_{\ell-1}} \star \dots \star A_{\nu_2\nu_1} \star G_{G \setminus \{\alpha\}}; \nu_1 \star A_{\nu_1\alpha} \star G_{G;\alpha}, \qquad (3.2.2a) $$

where $p = (\alpha, \nu_1, \dots, \nu_{\ell-1}, \omega)$ is a simple path of length $\ell$ from $\alpha$ to $\omega$ on $G$ and $G_{G;\alpha}$ is a Green's function given by

$$ G_{G;\alpha} = \left[ \mathrm{Id}_{\star} - \sum_{\gamma \in \Gamma(G)_{\alpha}} G_{G \setminus \{\alpha, \mu_1, \dots, \mu_{\ell'-2}\}; \mu_{\ell'-1}} \star \dots \star G_{G \setminus \{\alpha\}}; \mu_1 \star H_{\mu_1\alpha} \right]^{\star-1}, \qquad (3.2.2b) $$

with $\gamma = (\alpha, \mu_1, \dots, \mu_{\ell'-1}, \alpha)$ is a simple cycle of length $\ell'$ from $\alpha$ to itself on $G$.

Note, in these expressions the time arguments $(t', t)$ have been omitted for the sake of clarity.

**Remark.** As expected of a branched continued fraction, the Green's function $G_{G;\alpha}$ is obtained recursively through Eq. (3.2.2b). Indeed, $G_{G;\alpha}$ is expressed in terms of Green's functions such as $G_{G \setminus \{\alpha, \mu_1, \dots, \mu_{\ell-2}\}; \mu_{j-1}}$, which is in turn defined through Eq. (3.2.2b) but on the smaller subgraph $G \setminus \{\alpha, \dots, \mu_{j-1}\}$ of $G$. The recursion stops when vertex $\mu_j$ has no neighbour on this subgraph, in which case

$$ G_{G \setminus \{\alpha, \mu_1, \dots, \mu_{\ell-2}\}; \mu_{j-1}} = \begin{cases} [\mathrm{Id}_\star - A(t)_{\mu_j \mu_j}]^{\star-1} & \text{if } A(t)_{\mu_j \mu_j} \neq 0 \text{ for some } t \in I, \\ \mathrm{Id}_\star & \text{otherwise,} \end{cases} $$

and $\mathrm{Id}_\star$ is the ⋆-algebra matrix identity of appropriate dimension.
**Example 3.2.1.** In this example, we consider the following system of differential equations on the interval $I = [0, T]$, $T > 0$, for $t \in I$

$$ \tilde{A}(t)v(t) = \dot{v}(t), \quad \tilde{A}(t) = \begin{pmatrix} 0 & t & 0 \\ 0 & 0 & 1 \\ 1 & 0 & 0 \end{pmatrix}, \quad v(0) = v_0. $$

Matrix $\tilde{A}(t)$ does not commute with itself at different times: $[\tilde{A}(t'), \tilde{A}(t)] \neq 0$ for $t' \neq t$. Furthermore the Magnus series representation of the solution of this system is guaranteed to be convergent only for $t < \sqrt{2\pi}$ [23] and indeed the path-sum solution reveals the appearance of a negative real eigenvalue in the time-ordered exponential of $\tilde{A}$ at $t \simeq 1.15\sqrt{2\pi}$, implying the divergence of the Magnus series from that point on.

To present the path-sum formulation of the solution we first note that the graph corresponding to $\tilde{A}$ is an oriented triangle, which we designate $T$:

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

Let us detail the calculation of entry $U_{11}$ of the time-ordered exponential of $\tilde{A}$. There is a single simple cycle from vertex 1 to itself on $T$, namely $1 \to 3 \to 2 \to 1$. Furthermore $T$ is acyclic as soon as a vertex is removed from it. Thus, following Theorem 3.2.1, the path-sum formulation for $U_{11}$ is

$$ U_{11} = \int_{0}^{t} G_{T;11}(\tau) d\tau, \quad G_{T;11}(\tau) = \left( 1_{\star} - \underbrace{\begin{pmatrix} \text{Edge } 1 \leftarrow 2 & \text{Edge } 2 \leftarrow 3 & \text{Edge } 3 \leftarrow 1 \end{pmatrix}}_{\text{Triangle}} \right)^{\star-1} (\tau, 0). $$

here we omitted writing the $\Theta$ altogether to facilitate the identification of the edge weights but $t$ above should be understood as $t'\Theta(t'-t)$, 1 and as $1\Theta(t'-t)$. The $\star$-resolvent of $f(t', t) := (t\star 1\star 1)(t', t) = (t'-t)^2t'/2\Theta(t'-t)$ is obtained from the $\star$-Neumann series $\sum_n f^{*n}$, i.e. the Peano-Baker series for $\tilde{f}$. This gives

$$ G_{T;11}(t) = \delta(t) + \frac{1}{2} t^3 {}_0F_2 \left( ; \frac{5}{4}, \frac{3}{2}; \frac{t^4}{64} \right) \Rightarrow U_{11}(t) = \int_0^t G_{T;11}(\tau) d\tau = {}_0F_2 \left( ; \frac{1}{4}, \frac{1}{2}; \frac{t^4}{64} \right). $$

Similarly we get

$$ U_{22}(t) = {}_0F_2 \left( ; \frac{1}{2}, \frac{3}{4}; \frac{t^4}{64} \right), \quad U_{33}(t) = {}_0F_2 \left( ; \frac{1}{4}, \frac{3}{4}; \frac{t^4}{64} \right). $$

Note how these entries differ from $U_{11}(t)$, despite arising from circular permutations of the same triangle (cycles 1321, 2132 and 3213). This is because $t$ and 1 do not $\star$-commute, that is $t \star 1 \neq 1 \star t$.

Since there is a single simple path between any two vertices of $T$, path-sum formulation for the off-diagonal elements is, for example,

$$ \begin{aligned} U_{12}(t) &= \left( \Theta \star \overbrace{t}^{\text{Edge } 2 \leftarrow 1} \star G_{T;22} \right) (t) = \int_0^t d\tau \int_0^\tau d\tau' \tau G_{T;22}(\tau'), && \text{Path} \\ U_{13}(t) &= \left( \Theta \star \overbrace{t}^{\text{Edge } 1 \leftarrow 2} \star \overbrace{1}^{\text{Edge } 2 \leftarrow 3} \star G_{T;22} \right) (t) = \int_0^t d\tau \int_0^\tau d\tau_1 \int_{\tau_1}^\tau d\tau_2 \tau G_{T;33}(\tau_1), && \text{Path} \end{aligned} $$

where again $t \equiv t'\Theta(t' - t)$ etc. Similarly,

$$ U_{21}(t) = \int_{0}^{t} d\tau (1 \star 1 \star G_{T;11})(\tau), \quad \text{Path} $$

$$ U_{23}(t) = \int_{0}^{t} d\tau (1 \star G_{T;33})(\tau), \quad \text{Path} $$

$$ U_{31}(t) = \int_{0}^{t} d\tau (1 \star G_{T;11})(\tau), \quad \text{Path} $$

$$ U_{32}(t) = \int_{0}^{t} d\tau (1 \star t \star G_{T;22})(\tau), \quad \text{Path} $$
Calculating these (straightforward) integrals, the time-ordered exponential of $\tilde{A}$ is found to be

$$ U(t) = \begin{pmatrix} Q_{\frac{1}{4}, \frac{1}{2}}(t) & \frac{1}{2} t^2 Q_{\frac{3}{4}, \frac{3}{2}}(t) & \frac{1}{2} t^3 Q_{\frac{3}{4}, \frac{5}{4}}(t) - \frac{1}{6} t^3 Q_{\frac{1}{4}, \frac{7}{4}}(t) \\ t^2 Q_{\frac{1}{2}, \frac{5}{4}}(t) - \frac{1}{2} t^2 Q_{\frac{1}{4}, \frac{3}{2}}(t) & Q_{\frac{1}{2}, \frac{3}{4}}(t) & t Q_{\frac{3}{4}, \frac{5}{4}}(t) \\ t Q_{\frac{1}{2}, \frac{5}{4}}(t) & \frac{1}{2} t^3 Q_{\frac{3}{4}, \frac{3}{2}}(t) - \frac{1}{3} t^3 Q_{\frac{1}{2}, \frac{7}{4}}(t) & Q_{\frac{1}{4}, \frac{3}{4}}(t) \end{pmatrix}, $$

where $Q_{a,b}(t) := {}_0F_2(a,b;t^4/64)$. While the solution is closed form here, we emphasise that this very rarely the case as time-ordered exponentials typically produce hitherto undescribed transcendent functions. In the case of 2 × 2 systems, Kovacic's algorithm allows for a rigorous a priori decision on the nature of the solution (transcendent versus non-transcendent) relying on Picard-Vessiot theory [110]. The question is non-trivial, profound difficulties must still be overcome before it can effectively be extended to larger systems.

**Example 3.2.2** (Time-dependent Dyson equation). Let $\tilde{A}(t) \in C^\infty[I]^{n \times n}$, $A := \tilde{A}(t')\Theta(t' - t)$ and consider a 2 × 2 partition of A into 4 blocks. Let $A_{\alpha\alpha}$ and $A_{\omega\omega}$ be the two diagonal blocks and $G \equiv K_2$ be the corresponding graph on vertices $\alpha$ and $\omega$. This means that $A_{\alpha\alpha}$ and $A_{\omega\omega}$ are the time-dependent weights of the self-loops on vertices $\alpha$ and $\omega$, respectively. Similarly, $A_{\omega\alpha}$ and $A_{\alpha\omega}$ designate the off-diagonal blocks and are the weights of the edges $(\alpha\omega)$ and $(\omega\alpha)$, respectively. Then Theorem 3.2.1 states that the Green's function at $\alpha$ is given by

$$ G_{K_2; \alpha} = [\mathrm{Id}_* - A_{\alpha\alpha} - A_{\alpha\omega} \star G_{K_2 \setminus \{\alpha\}; \omega} \star A_{\omega\alpha}]^{*-1}, \quad (3.2.5) $$

and $G_{K_2 \setminus \{\alpha\}; \omega} = (\mathrm{Id}_* - A_{\omega\omega})^{*-1}$. The above result is the celebrated Dyson equation. To see this, let $G := G_{K_2; \alpha}$, $G_0 := (\mathrm{Id}_* - A_{\alpha\alpha})^{*-1}$ and $\Sigma := A_{\alpha\omega} \star G_{K_2 \setminus \{\alpha\}; \omega} \star A_{\omega\alpha}$. Then Eq. (3.2.5) indicates that the Green's function $G$ is $G = (G_0^{*-1} - \Sigma)^{-1}$, or equivalently,

$$ G = G_0 + G \star \Sigma \star G_0. $$

This is the time-dependent Dyson equation [6], which arises naturally from resummations of the Dyson series for the time-ordered exponential in the context of quantum many-body physics. This equation appears for example when considering a physical entity (for example a particle, or an ensemble of sites in a solid) in contact with a larger system. In this situation, the Hamiltonian driving the system + entity is naturally partitioned into four submatrices (blocks): $A_{\alpha\alpha}$, which drives the isolated entity, $A_{\omega\omega}$ which drives the rest of the system without the entity and $A_{\omega\alpha}(t)$ and $A_{\alpha\omega}(t)$, which represent the interactions system-entity. In this context, $\Sigma$ is known as the self-energy.

Hence we see that the Dyson equation stems from Theorem 3.2.1 on the complete graph on two vertices. In general, Theorem 3.2.1 can be seen as extending the Dyson equation to an arbitrary number of systems/entities in contact with each other. It therefore provide an explicit non-perturbative formula for the self-energy that involves finitely many terms for finite systems.

Decay properties of the time-ordered exponential

In the last 20 years, a number of significant results have established exponentially decaying bounds for the magnitude of the entries of holomorphic functions of sparse matrices [16, 17, 18]. These results have given rise to a flurry of applications in linear algebra and physics as they underlie efficient approximation techniques, see e.g. [16, 40, 183]. The techniques used to prove these results do not extend to the time-ordered exponential function however. Relying instead on $\star$-resolvent formalism, we establish super-exponential decaying bounds for the magnitude of the entries of the time-ordered exponential of sparse matrices.

We have seen that the time-ordered exponential U of any-time-dependent matrix A is $U = \Theta \star (\mathrm{Id}_* - A)^{-1}$. Furthermore, the Peano-Baker series for $\tilde{A}$ is unconditionally convergent so that $(\mathrm{Id}_* - A)^{-1} = \sum_{n \ge 0} A^{*n}$. Using the Taylor series remainder theorem then leads to the following bound for the entries of the time-ordered exponential of a sparse matrix:

**Proposition 3.2.1.** Let $I \subset \mathbb{R}$, $\tilde{A}(t) \in C^\infty[I]^{n \times n}$ and $A := \tilde{A}\Theta$. Let $G$ be the graph encoding the sparsity structure of A. Let $d := d(\alpha, \omega)$ be the length of the shortest path from vertex $\alpha$ to vertex $\omega$ on G. Define $h := \sup_{t \in I} \max_{\alpha, \beta} |\tilde{A}_{\alpha\beta}(t)|$ and let $|W_{G;\alpha\omega}|$ be the number of walks of length $k$ from $\alpha$ to $\omega$ on G. Then

$$ |U_{\omega\alpha}(t', t)| \le \sum_{k \ge d} \frac{h^k}{k!} (t' - t)^k |W_{G;\omega\omega}|, \quad (3.2.6a) $$

with equality when $\tilde{A}$ is a constant multiple of the ordinary adjacency matrix of G. Let $\Delta$ be the maximum degree of any vertex of G. If $\Delta$ is finite, we have

$$ |U_{\omega\alpha}(t', t)| \le e^{\Delta h(t' - t)} \frac{(\Delta h(t' - t))^d}{d!}. \quad (3.2.6b) $$
The bound of Eq. (3.2.6b) demonstrates the super-exponential decay of the ordered-exponential function of any time-dependent sparse matrix, contingent on the assumption that $\Delta$ is finite$^4$. Furthermore, the result of the proposition is non-trivial only when the maximum distance $D := \max_{\alpha, \omega} d(\alpha, \omega)$ between any two vertices $\alpha$ and $\omega$ on $G$ is infinite. Otherwise, a super-exponentially decaying bound can always be found for any matrix, by choosing a large enough multiplying constant. We emphasise the validity of these results in the case of infinite matrices of finite norm. Indeed, convergence of the Peano-Baker series necessitates only $\sup_{t \in I} \|\tilde{A}(t)\|$ to be finite.

## 3.3 Tridiagonalization of differential systems: ⋆-Lanczos algorithm

The ⋆-resolvent formulation of the solution to systems of non-autonomous linear ODEs, Eq. (3.2.1), shows that the ⋆-algebra formed by matrices of $Sm_{\Theta}^{N \times N}$ is the relevant algebraic context for dealing with such differential systems. This algebra allows for the import of linear algebraic methods in the differential context. In this section, we report on the first such method: an extension to differential systems of the celebrated Lanczos algorithm. Given a matrix $M \in \mathbb{C}^{N \times N}$ and two vectors $\mathbf{w}$ and $\mathbf{v} \in \mathbb{C}^N$ as inputs and barring breakdowns, this algorithm outputs a tridiagonal matrix $T_N \in \mathbb{C}^{N \times N}$ such that for all $n \in \mathbb{N}$, $(T_N^n)_{1,1} = \mathbf{w}^H A^n \mathbf{v}$. Such moment matching property, if it could be extended from the ordinary matrix product to the ⋆-product, would allow for the computation of time-ordered exponentials from tridiagonal matrices.

Indeed, let $A \in Sm_{\Theta}^{N \times N}$. We have seen that the time-ordered exponential $U(t', t)$ of $A$ is $U = \Theta \star G_A$, where the Green's function $G_A$ is the ⋆-resolvent of $A$, that is

$$ G_A := (\mathrm{Id}_\star - A)^{\star-1} = \sum_{n \ge 0} A^{\star n}, \qquad (3.3.1a) $$

the series on the right-hand side converging when $\tilde{A}$ elements are bounded over $I^2$. Note that here we wrote $\mathrm{Id}_\star$ for $\mathrm{Id} \times \mathrm{Id}_\star$. In particular then the time-ordered exponential $U$ of $A$ obeys

$$ \mathbf{w}^H \cdot \mathbf{U} \cdot \mathbf{v} = \Theta \star \sum_{n \ge 0} \mathbf{w}^H \cdot \mathbf{A}^{\star n} \cdot \mathbf{v}. \qquad (3.3.1b) $$

This confirms that if a tridiagonal $T_N$ could be found with a ⋆-moment matching property, then $\mathbf{w}^H \cdot \mathbf{U} \cdot \mathbf{v}$ could be evaluated from it. As with the original Lanczos method, this presents the major advantage that $T_N$ is extremely sparse while $A$ can be full.

Now we can present our results on the matter [70]: *baring breakdowns*—which we will characterise below—the ⋆-Lanczos algorithm presented in Table 3.1 produces a sequence of tridiagonal matrices $T_n, 1 \le n \le N$, of the form

$$ T_n := \begin{bmatrix} \alpha_0 & \delta & & & \\ & \alpha_1 & \delta & & \\ & & \ddots & \ddots & \\ & & & \delta & \\ & & & & \alpha_{n-1} & \delta \\ & & & & & \alpha_{n-1} \end{bmatrix}, \qquad (3.3.2) $$

and such that the matching ⋆-moment property is achieved:

**Theorem 3.3.1.** Let $A, \mathbf{w}, \mathbf{v}$ and $T_n$ be as described above, then

$$ \mathbf{w}^H (\mathbf{A}^{\star j}) \mathbf{v} = \mathbf{e}_1^H (T_n^{\star j}) \mathbf{e}_1, \quad \text{for } j = 0, \dots, 2n-1. \qquad (3.3.3) $$

Combining this with Eqs. (3.3.1) we have, for $n = N$, the exact expression

$$ \mathbf{w}^H \mathbf{U} \mathbf{v} = \Theta \star (\mathrm{Id}_\star - T_n)^{\star-1}_{1,1}, $$

while for $n < N$, the right-hand side yields an approximation to the time-ordered exponential. The path-sum theorem then gives,

$$ (\mathrm{Id}_\star - T_n)^{\star-1}_{1,1} = (1_\star - \alpha_0 - (1_\star - \alpha_1 - (1_\star - (\dots(1_\star - \alpha_{n-1})^{\star-1} \dots))^{\star-1} \star \beta_2)^{\star-1} \star \beta_1)^{\star-1}. $$

The ⋆-resolvents in the expression above are accessible either numerically or analytically from the ⋆-Neumann series.

<footnote n="4">⁴In [67] we give an example where $\Delta$ infinite leads to 'only' exponential decay.</footnote>
Input: A complex time-dependent matrix A, and complex vectors w, v such that $w^H v = 1$.

Output: Coefficients $\alpha_0, \dots, \alpha_{n-1}$ and $\beta_1, \dots, \beta_{n-1}$ defining the matrix $T_n$ of Eq. (3.3.2) which satisfies Eq. (3.3.3).

Initialize: $\mathbf{v}_{-1} = \mathbf{w}_{-1} = 0, \mathbf{v}_0 = \mathbf{v} 1_\star, \mathbf{w}_0^H = \mathbf{w}^H 1_\star$.

$$ \alpha_0 = \mathbf{w}^H A \mathbf{v}, $$

$$ \mathbf{w}_1^H = \mathbf{w}^H A - \alpha_0 \mathbf{w}^H, $$

$$ \hat{\mathbf{v}}_1 = A \mathbf{v} - \mathbf{v} \alpha_0, $$

$$ \beta_1 = \mathbf{w}^H A^{*2} \mathbf{v} - \alpha_0^{*2}, $$

If $\beta_1$ is not $\star$-invertible, then stop, otherwise,

$$ \mathbf{v}_1 = \hat{\mathbf{v}}_1 \star \beta_1^{*-1}, $$

For $n = 2, \dots$

$$ \alpha_{n-1} = \mathbf{w}_{n-1}^H \star A \star \mathbf{v}_{n-1}, $$

$$ \mathbf{w}_n^H = \mathbf{w}_{n-1}^H \star A - \alpha_{n-1} \star \mathbf{w}_{n-1}^H - \beta_{n-1} \star \mathbf{w}_{n-2}^H, $$

$$ \hat{\mathbf{v}}_n = A \star \mathbf{v}_{n-1} - \mathbf{v}_{n-1} \star \alpha_{n-1} - \mathbf{v}_{n-2}, $$

$$ \beta_n = \mathbf{w}_n^H \star A \star \mathbf{v}_{n-1}, $$

If $\beta_n$ is not $\star$-invertible, then stop, otherwise,

$$ \mathbf{v}_n = \hat{\mathbf{v}}_n \star \beta_n^{*-1}, $$

end.

<caption>Table 3.1: The $\star$-Lanczos algorithm.</caption>

**Remark.** The described tridiagonalization of the system of ODEs effected by Theorem 3.3.1 can also be seen as a $\star$-factorization of the matrix A. Consider the matrices $W_N = [\mathbf{w}_0, \dots, \mathbf{w}_{N-1}]$ and $V_N = [\mathbf{v}_0, \dots, \mathbf{v}_{N-1}]$ composed of the vectors computed by the $\star$-Lanczos algorithm. Then

$$ T_N = W_N^H \star A \star V_N, \quad W_N^H \star V_N = \mathrm{Id}_\star, \qquad (3.3.4) $$

$$ \text{and } G_A = V_N \star (\mathrm{Id}_\star - T_N)^{-1} \star W_N^H $$

A crucial assumption underlying these results is that the algorithm suffers no breakdown. This is related to the nature of the $\alpha_j$ and $\beta_j$ distributions appearing in the $T_n$ matrices and which are produced by the $\star$-Lanczos procedure through recurrence relations. These necessitate the $\star$-inversion of the $\beta_j$ without which the algorithm breaks down. If $\beta_j$ is not identically null, the existence of $\beta_j^{*-1}$ is established by Corollary 3.1.1 assuming *ad minima* that $\alpha_j, \beta_j \in \mathrm{Sm}_\theta$. Proving that this assumption is always verified turns out to be surprisingly difficult. In the course of proof, we also established that the remaining breakdowns of the $\star$-Lanczos-arising when a $\beta_j$ is identically null-coincide exactly with the breakdowns of the non-Hermitian Lanczos algorithm.

Before we state the main theorems on the tridiagonalization of systems of coupled linear differential equations with non-constant coefficients, we begin by exhibiting a relation between breakdowns in the $\star$-Lanczos procedure and breakdowns in the ordinary non-Hermitian Lanczos procedure. This shows that the feasibility of tridiagonalization does not depend on the nature of the entries of the original matrix nor on the kind of product between these entries. Rather breakdowns in tridiagonalization must be topological in origin, i.e. they depend on the structure and the edge weights of the graph whose adjacency matrix is A. We recall that, given as inputs a time-independent matrix B and time-independent vectors p, q, the ordinary non-Hermitian Lanczos algorithm aims to compute the (time-independent) vectors $\mathbf{p}_0, \dots, \mathbf{p}_n$ and $\mathbf{q}_0, \dots, \mathbf{q}_n$, respectively biorthonormal bases of the Krylov subspaces

$$ \mathrm{span}\{\mathbf{p}, \mathbf{B}\mathbf{p}, \dots, \mathbf{B}^n\mathbf{p}\}, \quad \mathrm{span}\{\mathbf{q}, \mathbf{B}^H\mathbf{q}, \dots, (\mathbf{B}^H)^n\mathbf{q}\}. $$

There is a serious breakdown at the $n$th iteration of the ordinary non-Hermitian Lanczos algorithm when $\mathbf{q}_n^H \mathbf{p}_n = 0$ and $\mathbf{p}_n, \mathbf{q}_n \neq 0$; see [155, Remark 3.2]. Assuming no breakdown arises before the $n$th iteration, the non-Hermitian Lanczos algorithm also gives the $n \times n$ tridiagonal matrix

$$ S_n := [\mathbf{q}_0, \dots, \mathbf{q}_{n-1}]^H \mathbf{B} [\mathbf{p}_0, \dots, \mathbf{p}_{n-1}]. \qquad (3.3.5) $$

The first key result relating breakdowns of the $\star$-Lanczos algorithm and of its classical non-Hermitian counterpart is the following Lemma:
# Chapter 4

## Further applications of the path-sum theorem

As explained in §2.4, the first field of application of the path-sum theorem is in linear algebra, in particular for the representation of matrix functions. The most immediate generalizations stem from mathematical problems that are naturally mappable to graphs. This is for example the case of Gaussian Markov random fields (GMRF), discussed here briefly. We then summarise example of applications of path-sum in conjunction with ⋆-calculus in quantum mechanics and in the representation of transcendent functions.

### 4.1 Statistical inference

A GMRF is a random vector that follows a multivariate normal distribution and satisfies conditional independence assumptions. If $X_1, X_2, X_3$ are random variables with a joint probability density function (or joint probability mass function in a discrete case), we say that $X_1$ is conditionally independent of $X_2$ given $X_3$, denoted $X_1 \perp X_2|X_3$, if $f(x_1, x_2|x_3) = f(x_1|x_3)f(x_2|x_3)$. Here we use $f$ as a generic symbol for the probability density function of the random variables corresponding to its arguments. GMRFs have a simple interpretation and find their applications, for example, in image analysis, spatial statistics, structural time series analysis and analysis of longitudinal and survival data [174].

Consider now a random vector $\mathbf{X} = (X_1, X_2, \dots, X_n) \sim N(\mu, \Sigma)$ following a multivariate normal distribution with mean $\mu$ and covariance matrix $\Sigma$. The probability density function of $\mathbf{X}$ is given as

$$f(x) = \frac{1}{\sqrt{(2\pi)^n \det(\Sigma)}} \exp \left[ -\frac{1}{2} (x - \mu)^T \Sigma^{-1} (x - \mu) \right] = g(x) \exp \left[ -\frac{1}{2} x^T J x + h^T x - k(\mu, \Sigma) \right]. \quad (4.1.1)$$

Here $\Sigma$ is a symmetric and positive definite matrix, $J = \Sigma^{-1}$ is called information matrix, $h = J\mu$ is called the potential vector, $g(x) = (2\pi)^{-\frac{n}{2}}$ and $k(\mu, \Sigma) = \frac{1}{2}\mu^T J\mu - \frac{1}{2}\ln(\det(J))$. We may now construct a graphical model for statistical inference on GMRF as follows.

Let $G = (V, E)$ be a graph and denote $X_{ij}$ the set of variables with $X_i$ and $X_j$ removed from $\mathbf{X}$. If $J = (J_{ij})_{i,j \in V}$ is positive definite, then for $i, j \in V$, where $i \neq j$, we have [114, Proposition 5.2] $X_i \perp X_j | X_{ij} \Leftrightarrow J_{ij} = 0$.

**Definition 4.1.1** ([174]). A random vector $\mathbf{X}$ is called a GMRF with respect to a graph $G = (V, E)$ with information matrix $J$ and potential vector $h$ if and only if its density has the form of (4.1.1) and $J_{ij} \neq 0 \Leftrightarrow (i, j) \in E$, for all $i$ and $j$.

It is known that $X$ satisfies the Markov property on $G$ [174, Theorem 2.4]<fnref n="1" />. By definition there is a one-to-one correspondence between the structure of $J$ and the structure of $G$. Most information matrices $J$ for GMRFs are sparse with only $O(n)$ non-zero entries in $J$ a key property facilitating the simulation of GMRFs through their information matrix.

Estimations of $\mathbf{X}$ given noisy observations $\mathbf{Y}$ are also readily accessible in this framework. Assume that $\mathbf{Y} = C\mathbf{X} + \varepsilon$, where $C$ is an $n \times n$ real matrix and $\varepsilon \sim N(0, M)$. Then the conditional distribution of $\mathbf{X}|\mathbf{Y}$ is given by $f(x|y) = \frac{f(y|x)f(x)}{f(y)} \propto \exp \left[ -\frac{1}{2}x^T \tilde{J}x + \tilde{h}^T x \right]$, where $\tilde{J} = J + C^T M^{-1}C$ and $\tilde{h} = h + C^T M^{-1}y$. Thus given noisy observations, one only needs to update the information matrix and the potential vector to construct a graphical model for $f(x|y)$. For simplicity, we shall still use $J$ and $h$ to denote the parameters after absorption of observations.

<footnote n="1">¹For a GMRF, the pairwise Markov property, the local Markov property and the global Markov property are equivalent. This is proven by using [114, Proposition 3.8], in conjunction with the Hammersley-Clifford Theorem [114, Theorem 3.9].</footnote>
As apparent from Eq. 4.1.1 one needs both the covariance matrix $\Sigma = J^{-1}$ and the mean vector $\mu = J^{-1}h$ to obtain the marginal distributions of $\mathbf{X}$ or $\mathbf{X}|\mathbf{Y}$. Since knowing $J^{-1}$ is sufficient to recover $\mu$, we focus our efforts on calculating $J^{-1}$. Direct inversion of $J$ has complexity $O(n^3)$ and does not exploit the sparsity of $J$.<fnref n="2" /> In a simple situation where the graph $G$ of a GMRF is a tree, belief propagation (BP) efficiently calculates the correct marginals [132, 158]. For graphs with cycles (also called loopy graphs), the method of Loopy Belief Propagation (LBP) can be used to efficiently approximate the marginals yet, estimates of the covariance matrices it provides are generally incorrect [200].

As earlier we denote $G\{\alpha, \beta \dots \}$ the subgraph of $G$ obtained by deleting from $G$ the vertices $\{\alpha, \beta \dots \} \subset V$ and the edges incident to them. The path-sum theorem immediately gives

**Theorem 4.1.1.** Let $J$ be a finite GMRF information matrix. Let $\Pi(G)_{\alpha\omega}$ and $\Gamma(G)_{\alpha}$ be the sets of simple paths from $\alpha$ to $\omega$ on $G$ and the set of simple cycles from $\alpha$ to itself on $G$, respectively. Then each entry of the covariance matrix $\Sigma = J^{-1}$ admits the path-sum representation

$$ J_{\omega\alpha}^{-1} = \sum_{p \in \Pi(G; \alpha\omega)} (-1)^{\ell(p)} \prod_{j=1}^{\ell(p)+1} \left\{ (J_{G\{\alpha, \nu_2, \dots, \nu_{j-1}\}})_{\nu_j \nu_j}^{-1} J_{\nu_{j+1} \nu_j} \right\} J_{\alpha\alpha}^{-1}, \\ J_{\alpha\alpha}^{-1} = \left( \sum_{\gamma \in \Gamma(G; \alpha\alpha)} (-1)^{\ell(\gamma)+1} J_{\mu_1 \mu_{\ell(\gamma)}} \prod_{j=2}^{\ell(\gamma)} \left\{ (J_{G\{\alpha, \mu_2, \dots, \mu_{j-1}\}})_{\mu_j \mu_j}^{-1} J_{\mu_j \mu_{j-1}} \right\} \right)^{-1} $$

where the products are right-to-left (i.e. $\prod_{i=1}^m a_i = a_m \cdots a_1$), $p = (\nu_1, \nu_2, \dots, \nu_{\ell(p)+1})$ is a simple path of length $\ell(p)$ with $\alpha \equiv \nu_1$ and $\omega \equiv \nu_{\ell(p)+1}$ for convenience; and $\gamma = (\mu_1, \mu_2, \dots, \mu_{\ell(\gamma)}, \mu_1)$ is a simple cycle of length $\ell(\gamma)$ from $\alpha \equiv \mu_1$ to itself.

**Remark.** The above path-sum formulation is valid for all matrix-partitions of $J$. For a GMRF $\mathbf{X}$, a partition of $J$ is equivalent to a partition of the set of random variables $\mathbf{X}$ into $B > 1$ disjoint subsets $\mathbf{X}_1, \dots, \mathbf{X}_B$. Let $\mathbf{X}' = (\mathbf{X}_1, \dots, \mathbf{X}_B)$ be a GMRF with respect to a new graph $G' = (V', E')$ with information matrix $J'$. Note that each $\mathbf{X}_i$ is now a random vector instead of a random variable. Following Definition 4.1.1 we define $\mathbf{X}_i \perp \mathbf{X}_j | \mathbf{X}'_{\setminus ij} \iff J'_{ij} = 0 \iff (j, i) \notin E'$, with $\mathbf{X}'_{\setminus ij}$ the set of variables with $\mathbf{X}_i$ and $\mathbf{X}_j$ removed from $\mathbf{X}'$. Note that $\mathbf{X}_i \perp \mathbf{X}_j | \mathbf{X}'_{\setminus ij}$ implies the global Markov property (GMP), which, for GMRFs, is equivalent to the pairwise Markov property (PMP).

**Remark.** The authors of [132] provided a walk-based derivation for Gaussian belief propagation on trees. The equations they obtain follow immediately from an application of the above theorem under the assumption that $G$ is a tree [63]. Another approach, presented in [102], proved that some entries of $J^{-1}$ are sums over the simple paths of $G$ with determinantal weights. This follows immediately from the above theorem on using the adjugate formula for the quantities $(J_{G\{\alpha, \nu_2, \dots, \nu_{j-1}\}})_{\nu_j \nu_j}^{-1}$ when using the standard partition of $J$ into its entries.

## 4.2 Quantum dynamics

A non-relativistic isolated quantum system evolves in time according to Schrdinger's equation, which we may synthesise as $(\hbar = 1)$

$$ \dot{U} = -iH.U $$

Here $H$ is the possibly time-dependent Hamiltonian operator, which describes the system's energy. If the physical system's degrees of freedom are purely discrete then $H$ is a (possibly infinite) matrix. Then Schrdinger's equation is a system of non-autonomous linear ODEs. We may consequently use the $\star$-product formalism in conjunction with the path-sum theorem to solve this equation. These purely mathematical observations also make sense physically, as we now briefly discuss.<fnref n="3" />

### Quantum evolution, walks on graphs and scale invariance

Quantum systems with discrete degrees of freedom such as spin systems, obey a discrete analog to Feynman's path integrals. To illustrate this from a physical standpoint, define one history of a quantum system as a temporal succession of orthogonal quantum states $h: |s_1\rangle \mapsto |s_2\rangle \mapsto |s_3\rangle \cdots$, each transition $|s_i\rangle \mapsto |s_{i+1}\rangle$ happening at a specified time $t_i$. Overall the history $h$ acquires a complex weight which is the product of the

<footnote n="2">²Algorithms that compute some entries of a symmetric $n \times n$ sparse matrix with complexity less than $O(n^3)$ do exist. The path-sum representation achieves this as well, computing the covariance of a pair of variables with complexity $O(n)$ whenever $G$ is a tree [59].</footnote>

<footnote n="3">³See also Chapter 4 of my PhD thesis for extended discussion and a physics-based proof of the walk-sum lemma, a weaker form of the path-sum theorem [59].</footnote>
weights of all the transitions in the history. The weight of an individual transition $|s_i\rangle \mapsto |s_{i+1}\rangle$ is dictated by the Hamiltonian as $\langle s_{i+1}|H(t_i)|s_i\rangle$.

A natural representation of such discrete histories is as walks on a graph. Define $G_t$ the graph where each vertex $v_i$ corresponds to one member $|s_i\rangle$ of an orthonormal basis for the entire state-space and give the directed edge $v_i \mapsto v_j$ the time-dependent weight $\langle s_j|H(t)|s_i\rangle$. In this picture, a system history as defined earlier is a walk on $G_t$ and $H(t)$ is the adjacency matrix of $G_t$. Because the Hamiltonian is time-dependent, the graph itself is dynamical.

Now just as for Feynman's path-integrals, the exact evolution of the system is obtained from the superposition of all its possible histories. Equivalently, every element $\langle s_j|U(t)|s_i\rangle$ of the evolution operator $U(t)$ is given by the sum over all walks from $v_i$ to $v_j$ on $G_t$, including all possible jumping times for each transition between vertices. While individual walks are the discrete counterpart of Feynman diagrams, their algebraic structure is better understood (Chapter 2). Since all walks are uniquely representable as nestings of simple cycles and paths, so walks series have a formulation involving only these: the path-sum. Then any $\langle s_j|U(t)|s_i\rangle$ is given as branched continued fractions comprising only the weights of the simple cycles and paths of the graph. Because the graph $G_t$ is finite, there are finitely many such cycles and paths and the fraction is finite in both depth and breadth. It is thus unconditionally convergent. As presented in Chapter 3, the same principles apply regardless of whether the Hamiltonian depends on time or not.

Combinatorially, path-sums stem from formal resummations of families of walks. This principle does not depend on what those walks represent. In particular, it remains unchanged by the nature of the evolving system. To exploit this observation, consider a more general type of system histories made of temporal successions of orthogonal vector spaces $\tilde{h}: V_1 \mapsto V_2 \mapsto V_3 \dots$. Physically such histories can describe an evolving subsystem, such as a group of protons in a large molecule. Mathematically they correspond to walks on a coarse-grained representation of the quantum state space, a subgraph $\tilde{G}_t$ of $G_t$. To see this, take a complete family of orthogonal spaces, i.e. $\bigoplus_{i=1} V_i = V$, where $V$ is the entire quantum state space. To each $V_i$ associate a vertex $v_i$ and give the edge $v_i \mapsto v_j$ the time-dependent weight $P_{V_j} H(t) P_{V_i}$. Here $P_{V_k}$ is the projector onto $V_k$. Observe then that these edge weights are generally non-Abelian. Yet, because path-sums fundamentally retain the order and time of the transitions in histories when performing resummations of walks, this setup poses no further difficulty. It follows that the submatrix $P_{V_j} U(t', t) P_{V_i}$ of the evolution operator is again given as a matrix-valued branched continued fraction of finite depth and breadth. While the shape of this fraction depends on the particular choice of vector spaces, its existence and convergence properties do not. If the vector spaces are chosen so that the shape of the fraction itself is unchanged, and such a choice is always possible, then the path-sum formulation is truly invariant under scale changes in the quantum state space.

An immediate consequence of scale-invariance is that there is always a path-sum calculation rigorously relating the global evolution of a system to that of any ensemble of its subsystems, such as clusters of spins in a large molecule. In this scheme, we can evolve each subsystem *separately* from one-another using any preferred method; only to then combine these isolated evolutions *exactly* via a path-sum to generate the true system evolution.

## Solutions for small non-autonomous systems

We give the explicit path-sum solutions for all $2 \times 2$ and $3 \times 3$ problems, solutions to an $8 \times 8$ and $14 \times 14$ problem are presented in [64] and [62], respectively.

### General $2 \times 2$ case

Consider the time-dependent matrix

$$ \tilde{H}(t) = \begin{pmatrix} \tilde{H}_{11}(t) & \tilde{H}_{12}(t) \\ \tilde{H}_{21}(t) & \tilde{H}_{22}(t) \end{pmatrix}, $$

where all $\tilde{H}_{ij} \in C^\infty[I]$ for some $I \subset \mathbb{R}$ compact. In the standard partition of $H$ into its entries, the corresponding graph $G$ is the complete graph $K_2$. Thus, the path-sum theorem gives the time-ordered matrix exponential $U$ of $H := \tilde{H}\Theta$, solution of $\dot{U} = \tilde{H}.U$, as

$$ \begin{aligned} U_{ii} &= \Theta \star G_i, \quad G_i := (1_* - H_{ii} - H_{ij} \star (1_* - H_{jj})^{*-1} \star H_{ji})^{*-1}, \\ U_{ji} &= \Theta \star (1_* - H_{jj})^{*-1} \star H_{ji} \star G_i, \end{aligned} $$

where $i, j = 1, 2$ and $i \neq j$. Furthermore $(1_* - H_{jj})^{*-1} = 1_* + H_{jj} \exp(\Theta \star H_{jj})$, $j = 1, 2$ by virtue of Proposition 3.1.3, p. 24. This notably implies $H_{ij} \star (1_* - H_{jj})^{*-1} = H_{ij} \exp(\Theta \star H_{jj})$ and $\Theta \star (1_* - H_{jj})^{*-1} = \exp(\Theta \star H_{jj})$. Since the graph comprises only two vertices, the path-sum theorem reduces to a use of the Schur complement in the $\star$-algebra $(C^\infty[I]^{2\times2}, \star)$, as may be recognised above.
General 3 × 3 case

We consider now a time-dependent matrix $\tilde{H} \in C^\infty[I]^{3\times3}$. The corresponding graph is the $K_3$ graph with self-loops. The path-sum representation of the time-ordered exponential of $H := \tilde{H}\Theta$ is $U = \Theta \star G$ with $G$ a matrix whose entries $G_{ij}(t', t)$, $i, j, k \in \{1, 2, 3\}$, $i \neq j \neq k$, are explicitly given by

$$
\begin{align*}
G_{ii} &= (1_{\star} - H_{ii} - (H_{ij} + H_{ik} \star G_{K_3 \setminus \{i,j\};k} \star H_{kj}) \star G_{K_3 \setminus \{i\};j} \star H_{ji} \\
&\qquad - (H_{ik} + H_{ij} \star G_{K_3 \setminus \{i,k\};j} \star H_{jk}) \star G_{K_3 \setminus \{i\};k} \star H_{ki})^{*-1}, \\
G_{ij} &= G_{K_3 \setminus \{j\};i} \star H_{ij} \star G_{jj} + G_{K_3 \setminus \{j,k\};i} \star H_{ik} \star G_{K_3 \setminus \{j\};k} \star H_{kj} \star G_{jj}.
\end{align*}
$$

In these expressions,

$$
\begin{align*}
G_{K_3 \setminus \{j\};i} &:= (1_{\star} - H_{ik} \star H_{ki} - H_{ik} \star G_{K_3 \setminus \{i,j\};k} \star H_{ki})^{*-1}, \\
G_{K_3 \setminus \{i,j\};k} &= 1_{\star} + H_{kk} \exp(\Theta \star H_{kk}).
\end{align*}
$$

Bloch-Siegert dynamics

The Bloch-Siegert Hamiltonian, here denoted $\tilde{H}_{BS}(t)$, is possibly the simplest model to exhibit non-trivial
physical effects (such as the Bloch-Siegert shifts [185]) due to time-dependencies in the driving radio-frequency
fields. The detailed study of these effects is of paramount importance in the broad field of quantum computing,
as they have a deleterious impact on qubit driving and stored quantum information [207]. The Hamiltonian
reads, in the $\{|\uparrow\rangle, |\downarrow\rangle\}$ basis,

$$
\tilde{H}_{BS}(t) = \begin{pmatrix} \omega_0/2 & 2\beta \cos(\omega t) \\ 2\beta \cos(\omega t) & -\omega_0/2 \end{pmatrix}. \qquad (4.2.1)
$$

This is usually given in the interaction picture $\tilde{H}_{BS}(t) = 2\beta \cos(\omega t) \cos(\omega_0 t)\sigma_x - 2\beta \cos(\omega t) \sin(\omega_0 t)\sigma_y$. In these expressions, the coupling parameter $\beta$ is the amplitude of the radio-frequency field. Beyond the rotating wave approximation—which omits the field’s counter-rotating terms and is limited to near resonant $\omega \sim \omega_0$ ultra-weak couplings $\beta/\omega \ll 1$ and in spite of the much theoretical efforts [7, 8, 204, 208], a non-perturbative truly analytical solution at all orders over the entire coupling range, on and off resonance, is ultimately lacking.

Path-sum solution

Since in the rotating frame, $(\mathbb{H}_{BS})_{11} = (\mathbb{H}_{BS})_{22} = 0$, the graph representing $\mathbb{H}_{BS}$ in its standard partition is $K_2$ but with no self-loops. Then the general solution given above yields the time-ordered exponential of $-i\mathbb{H}_{BS}$

$$
\begin{align*}
U(t)_{11} &= \int_0^t G_\uparrow(\tau, 0) d\tau, \quad U(t)_{22} = \int_0^t G_\downarrow(\tau, 0) d\tau, \\
U(t)_{21} &= -2i\beta \int_0^t \int_0^{\tau_1} \cos(\omega\tau_1) e^{i\omega_0\tau_1} G_\uparrow(\tau_0, 0) d\tau_0 d\tau_1, \quad U(t)_{12} = -2i\beta \int_0^t \int_0^{\tau_1} \cos(\omega\tau_1) e^{-i\omega_0\tau_1} G_\downarrow(\tau_0, 0) d\tau_0 d\tau_1,
\end{align*}
$$

while $G_\uparrow(t', t) = (1_\star - K_\uparrow\Theta)^{-1}$, $G_\downarrow = (1_\star - K_\downarrow\Theta)^{-1}$ with

$$
\begin{align*}
K_{\uparrow}(t', t) &= \frac{4\beta^2}{\omega^2 - \omega_0^2} \cos(\omega t') (k_{\uparrow}(t) e^{-i\omega_0(t'-t)} - k_{\uparrow}(t')), \\
K_{\downarrow}(t', t) &= \frac{i\beta^2}{\omega^2 - \omega_0^2} (1 + e^{-2i\omega t'}) (k_{\downarrow}(t') - k_{\downarrow}(t) e^{i(\omega+\omega_0)(t'-t)}),
\end{align*}
$$

where $k_\uparrow(t) = i\omega_0 \cos(\omega t) + \omega \sin(\omega t)$ and $k_\downarrow(t) = e^{2i\omega t}(\omega + \omega_0) - (\omega - \omega_0)$. In spite of the apparent divergences in the resonant case $\omega_0 \to \omega$, the kernels $K_\uparrow$ and $K_\downarrow$ are actually well defined in this limit. The quantity $G_\uparrow$ as obtained from $K_\uparrow$ has no closed form, rather it is a hitherto unknown higher special function. It is nonetheless analytically available thanks to unconditionally convergent ★-Neumann expansion $G_\uparrow = 1_\star + \sum_{k>0} K_\uparrow^{\star k}$.<fnref n="4" /> The ★-Neumann series is well suited to analytical computations, yielding (here displaying only the first two orders on resonance $\omega_0 = \omega$),

$$
G_{\uparrow}(t, 0) = \delta(t) - \frac{\beta^2}{\omega} e^{-i\omega t} \cos(\omega t) (-ie^{2i\omega t} + 2\omega t + i) \\
+ \frac{\beta^4}{24\omega^3} e^{-3i\omega t} \cos(\omega t) (3ie^{6i\omega t} + 6e^{4i\omega t} (-2i\omega^2 t^2 + 2\omega t + i) + e^{2i\omega t} (8\omega^3 t^3 + 12i\omega^2 t^2 + 12\omega t - 15i) - 12\omega t + 6i) + \dots
$$

<footnote n="4">⁴This observation holds for all $N \times N$ time-dependent Hamiltonians treated by path-sum.</footnote>
Of particular interest for qubit-driving experiments is the evolution of the transition probability $P_{↑→↓}(t) := |U_{21}(t)|^2$ from state $|↑⟩$ to state $|↓⟩$ [4, 126, 206]. This quantity is usually found perturbatively in the weak coupling regime using Floquet theory [185] as Magnus series again suffer from divergences [136]. It is now easily accessible–$U_{21}$ being given above. We find that $P_{↑→↓}(t)$ takes on the form of a Fourier-like series

$$P_{\uparrow \to \downarrow}(t) = \sum_{k=0}^{\infty} \sin(2k\omega t)S_{2k}(\beta, t) + \cos(2k\omega t)C_{2k}(\beta, t), \quad (4.2.2)$$

with $S_{2k}$ and $C_{2k}$ functions of $\beta$ and $t$, a representation of which is analytically available [62]. This form of $P_{↑→↓}(t)$ is due to the path-sum integral of $U_{21}$, which resembles a Fourier transform. We emphasise that this is not a general feature of path-sum nor of $2 \times 2$ Hamiltonians, but solely of the present Hamiltonian with linearly polarised driving.

## Anomalous resonances and accelerated ⋆-Neumann series

One of the striking effects related to Bloch-Siegert dynamics is an anomalous resonance known as Coherent Destruction of Tunneling (CDT) [82, 126]. Concretely this effect manifests itself when $4\beta/\omega$ falls on a root of the 0th Bessel function $J_0$. At this point there is a strong suppression of the transition between states $|\psi_{\pm}\rangle := \frac{1}{\sqrt{2}}(|\uparrow\rangle \pm |\downarrow\rangle)$ in spite of a very strong coupling $\beta/\omega_0 \gg 1$ between the $|\uparrow\rangle$ and $|\downarrow\rangle$ states in the Hamiltonian. In this situation we found it highly advantageous to see $-i\mathcal{H}_{BS}$ as comprising two parts,

$$K_1(t) = -2i\beta \begin{pmatrix} 0 & \cos(\omega t) \\ \cos(\omega t) & 0 \end{pmatrix}, \quad K_2 = -i\omega_0 \begin{pmatrix} \frac{1}{2} & 0 \\ 0 & -\frac{1}{2} \end{pmatrix}.$$

In the trivial partition where $\mathcal{H}_{BS}$ is a single block, the ⋆-resolvent formulation indicates that $U = \Theta \star (\mathrm{Id}_* + i\mathcal{H}_{BS}\Theta)^{-1}$, which is also the path-sum formulation for the graph with a single vertex and a self-loop. We calculate the ⋆-resolvent of $\mathcal{H}$ directly with an accelerated ⋆-Neumann series (see §3.1), suitably truncated for approximations. We expect this to work well because in strong coupling regime $\beta/\omega_0 \gg 1$, $K_1$ dominates $K_2$ and the ⋆-resolvents of $K_{1,2}$ are both known exactly. For reasons unknown we found this scheme to work beyond well in that it gave exact expressions at order 0 as all subsequent orders of the accelerated series vanished when on the CDT resonance. We emphasise that this is non-trivial and not understood. Anyhow at the 0th order of the accelerated series, we get the ⋆-resolvent as $G^{(acc,0)} = (\mathrm{Id}_* - K_1\Theta)^{-1} \star (\mathrm{Id}_* - K_2\Theta)^{-1}$ so that

$$U^{(acc,0)}(t) = \begin{pmatrix} \cos\left(\frac{2\beta}{\omega}\sin(\omega t)\right) + e^{-\frac{1}{2}i\omega_0 t} - 1 & -i\sin\left(\frac{2\beta}{\omega}\sin(\omega t)\right) \\ -i\sin\left(\frac{2\beta}{\omega}\sin(\omega t)\right) & \cos\left(\frac{2\beta}{\omega}\sin(\omega t)\right) + e^{\frac{1}{2}i\omega_0 t} - 1 \end{pmatrix} + \int_0^t \begin{pmatrix} i\omega_0 e^{-\frac{1}{2}i\omega_0 \tau} \sin^2\left(\frac{2\beta}{\omega}(\sin(\omega\tau) - \sin(\omega t))\right) & -\frac{1}{2}\omega_0 e^{\frac{1}{2}i\omega_0 \tau} \sin\left(\frac{4\beta}{\omega}(\sin(\omega\tau) - \sin(\omega t))\right) \\ \frac{1}{2}\omega_0 e^{-\frac{1}{2}i\omega_0 \tau} \sin\left(\frac{4\beta}{\omega}(\sin(\omega\tau) - \sin(\omega t))\right) & -i\omega_0 e^{\frac{1}{2}i\omega_0 \tau} \sin^2\left(\frac{2\beta}{\omega}(\sin(\omega\tau) - \sin(\omega t))\right) \end{pmatrix} d\tau.$$

The integrals in $U^{(acc,0)}(t)$ have no closed form as revealed by standard expansions over Bessel functions, e.g.

$$\sin(\alpha + z \sin(\phi)) = \sin(\alpha) \left( J_0(z) + 2 \sum_{n=1}^{\infty} J_{2n}(z) \cos(2n\phi) \right) + 2 \cos(\alpha) \sum_{n=0}^{\infty} J_{2m+1}(z) \sin((2n+1)\phi).$$

All physical quantities of interest are now accessible, for example, the return probability to the $|↑⟩$ state is found to be

$$P_{\uparrow \to \uparrow}^{(acc,0)}(t) = \left| \cos\left(\frac{2\beta}{\omega}\sin(\omega t)\right) + e^{-\frac{1}{2}i\omega_0} - 1 + \int_0^t i\omega_0 e^{-\frac{1}{2}i\tau\omega_0} \sin^2\left(\frac{\beta}{\omega}(\sin(\omega\tau) - \sin(\omega t))\right) d\tau \right|^2, \quad (4.2.3)$$

This formula becomes exact when either $\omega_0 \to 0$ or $\beta \to 0$, as expected from the acceleration procedure. In general, it provides excellent approximations when $\beta/\omega_0$ is large, see Fig. (4.1 a, b, c). The time-average of the return probability is therefore found to be

$$\langle P_{\uparrow \to \uparrow}^{(acc,0)}(t) \rangle_t = \frac{1}{2} \left( 1 + J_0\left(\frac{4\beta}{\omega}\right) \right), \quad (4.2.4)$$

which is exactly $1/2$ on CDT resonances where $J_0(4\beta/\omega) = 0$, consistent with the understanding of CDT offered by Floquet analysis. To be more precise let us study CDT directly by considering the states $|\psi_{\pm}\rangle = \frac{1}{\sqrt{2}}(|\uparrow\rangle \pm |\downarrow\rangle)$.
<figure><img src="image_7.png" /><figcaption>Figure 4.1: Coherent destruction of tunneling: Top line, return probability $P_{↑→↑}(t)$ in the ultra-strong coupling regime $\beta/\omega_0 = 30$ for (a) $\omega = 4\omega_0$; (b) $\omega = 20\omega_0$; and (c) $\omega = 100\omega_0$. Shown here are $P_{↑→↑}^{(acc,0)}(t)$ as given by Eq. (4.2.3) (solid blue line), the numerical solution (dashed black line), and its predicted time-average Eq. (4.2.4) (solid red straight line, indistinguishable from the numerically computed time-average). Middle line: transition probability $P_{\psi_{-} \to \psi_{+}}(t)$ for a system starting in the $|\psi_{-}\rangle$ state at $t = 0$ with: (d) $4\beta/\omega = 2.404...$, first zero of $J_0(4\beta/\omega)$; (e) $4\beta/\omega = 11.79...$, fourth zero of $J_0(4\beta/\omega)$; and (f) $4\beta/\omega = 27.49...$, ninth zero of $J_0(4\beta/\omega)$. Note the changes of scales. Shown here are the formula of Eq. (4.2.5) (solid blue line) and the numerical solution (dashed black line), these two being completely indistinguishable. Bottom line: far off-resonance $\omega = 100\omega_0$ expectation value of $\sigma_x$ for a system starting in the $|\uparrow\rangle$ state at $t = 0$ with: (g) $4\beta/\omega = 2.404...$, first zero of $J_0(4\beta/\omega)$; (h) $4\beta/\omega = 11.79...$, fourth zero of $J_0(4\beta/\omega)$; and (i) $4\beta/\omega = 27.49...$, ninth zero of $J_0(4\beta/\omega)$. Note the changes of scales in $\langle\sigma_x\rangle$. Shown here are the formula of Eq. (4.2.6) (solid blue line) and the numerical solution (dashed black line), these two being indistinguishable.</figcaption></figure>

The probability of transition between these states, denoted $P_{\psi_{-} \to \psi_{+}}(t)$, is found from $U^{(acc,0)}(t)$ in the situation where $\omega_0 \ll (\beta/\omega)^{1/2}$, as

$$P_{\psi_{-} \to \psi_{+}}^{(acc,0)}(t) = \frac{\omega_0^2}{4} \left( \int_0^t \sin \left( \frac{4\beta}{\omega} (\sin(\omega t) - \sin(\omega \tau)) \right) d\tau \right)^2 + \frac{\omega_0^2}{4} \left( \int_0^t \cos \left( \frac{4\beta}{\omega} (\sin(\omega t) - \sin(\omega \tau)) \right) d\tau \right)^2. \quad (4.2.5)$$

This expression *flawlessly* reproduces the numerical solution in its finest details, details which had hitherto not been captured with such accuracy [126]. Minimizing the time-average of this formula confirms that the CDT condition is exactly $J_0(4\beta/\omega) = 0$, i.e. this is not changed by the non-perturbative corrections. Mathematically, the reason for this is simple: the $J_0$ function is quadratically dominant over the other terms of the Bessel-series expansion of Eq. (4.2.5) because it stems from the sole term of that expansion which does not depend on $\tau$ in both integrals.

While these results are as expected from the standard theory of CDT, it is not so for all physical quantities. Consider for example, the expectation value of the observable $\sigma_x$ for a system initially prepared in the $|\uparrow\rangle$ state. As observed by [191], $\langle\sigma_x\rangle$ presents anomalous fluctuations on CDT resonances, a fact that was interpreted as a hallmark of and resulting from a crossing Floquet states. This interpretation is in fact not correct. Indeed,
at order 0 of the accelerated expansion of the path-sum solution we get,

$$ \langle \sigma_x \rangle^{(acc,0)} = \omega_0 \int_0^t \cos\left(\frac{1}{2}\omega_0\tau\right) \sin\left(\frac{2\beta}{\omega}\sin(\omega\tau)\right) d\tau \\ + 2\omega_0 \sin\left(\frac{1}{4}\omega_0 t\right) \int_0^t \sin\left(\frac{1}{4}\omega_0(t-2\tau)\right) \sin\left(\frac{2\beta}{\omega}(\sin(\omega t) - \sin(\omega\tau))\right) d\tau. $$

In the regime $\omega_0 \ll \omega$, both $\cos(\omega_0 t/2)$ and $\sin(\omega_0 t/4)$ are essentially equal to their initial $t = 0$ values, leading to the simpler

$$ \langle \sigma_x \rangle^{(acc,0)} = \omega_0 \int_0^t \sin \left( \frac{4\beta}{\omega} \sin(\omega\tau) \right) d\tau. \qquad (4.2.6) $$

This simple expression fits once again absolutely flawlessly with the numerically computed expectation $\langle \sigma_x \rangle$, see Fig. (4.1 d, e, f). Now evaluating the integral remaining in Eq. (4.2.6) via Bessel functions shows that the time average of $\langle \sigma_x \rangle$ is

$$ \langle\langle \sigma_x \rangle\rangle^{(acc,0)}_t = \frac{2\omega_0}{\omega} \sum_{n=0}^{\infty} J_{2n+1} \left( \frac{4\beta}{\omega} \right) \frac{1}{2n+1}, $$

whose extrema are reached whenever

$$ 1 - \frac{\pi}{2} H_1 \left( \frac{4\beta}{\omega} \right) = 0, \qquad (4.2.7) $$

with $H_1(.)$ the first Struve function. Remarkably, the difference $\Delta_n$ between the location of the nth zero of $J_0(.)$ and of the nth zero of Eq. (4.2.7) tends asymptotically to 0 as $\Delta_n \sim 1/(2\pi n)$ for $n \gg 1$. This asymptotics develops quite quickly: while $\Delta_1 \simeq 0.4$, already $\Delta_2 \simeq 0.03$. The fact that the anomalous fluctuations in the expectation value of $\sigma_x$ peak at the zeroes of Eq. (4.2.7) rather than on CDT resonances is confirmed by the numerical simulations. This analysis indicates that while $\langle \sigma_x \rangle$ does indeed seem to fluctuate the most on CDT resonances, it is in fact not true and the phenomenon driving these fluctuations is subtly different from that behind CDT.

These results demonstrate the power of various expansions of the path-sum solution, enabling very precise and hitherto unequaled analytical analysis of subtle phenomena, e.g. $P_{\psi_{-} \to \psi_{+}}(t)$ is on the order of $10^{-5}$ on CDT resonances and reproduced exactly by the formulae provided. This is not because of special features of the Bloch-Siegert Hamiltonian. Rather, the path-sum approach and acceleration of ★-Neumann series are both generally valid as shown on p. 34. In particular the same treatment is valid for dissipative non-Hermitian operators.

## 4.3 Representation of Heun functions

Heun differential equations are the most general second order Fuchsian equations with four regular singularities. Series representations for the Heun functions have finite radiuses of convergence causing difficulties notably for black hole perturbation theory [135], where conditions at the black hole horizon ($z = 1$) cannot easily be inferred from Earth-bound observations ($z = +\infty$). These issues were noted in the recent review [97] on Heun's functions and are well-understood mathematically to be inherently tied with the use of series representations. This has been clearly identified as a major obstacle when extracting physical meaning in areas of black holes astrophysics [97], yet remains unaddressed in the mathematical literature. The quest for an unconditionally-convergent representation has led physicists to seek integral transformations relating Heun functions with other functions, all of which are transcendent. In spite of these efforts, [97] observes that: "No example has been given of a solution of Heun's equation expressed in the form of a definite integral or contour integral involving only functions which are, in some sense, simpler. [...] This statement does not exclude the possibility of having an infinite series of integrals with 'simpler' integrands".

Thanks to the path-sum theorem and ★-algebra, we provide such elementary integral representations of the solutions of all equations of the Heun class: general, confluent, bi-confluent, doubly-confluent and triconfluent. For brevity's sake we here only state the results in the general Heun function case. Our results lead to the first representation of the solution to the Teukolsky radial equation<fnref n="5" /> governing the metric perturbations of rotating black holes that is convergent everywhere from the black hole horizon up to spatial infinity.

<footnote n="5">⁵The Teukolsky Equation [189] is a gauge invariant equation that governs the curvature perturbations of the Kerr black hole [147].</footnote>
## Heun's equation

In its canonical form, the General Heun Equation (GHE) is given as follows:

$$ \frac{d^2y(z)}{dz^2} + \left[ \frac{\gamma}{z} + \frac{\delta}{z-1} + \frac{\epsilon}{z-t} \right] \frac{dy(z)}{dz} + \frac{\alpha\beta z - q}{z(z-1)(z-t)} y(z) = 0. \quad (4.3.1) $$

In this equation, $q \in \mathbb{C}$ while the other parameters satisfy the Fuch's condition: $1 + \alpha + \beta = \gamma + \delta + \epsilon$. The GHE has four singular points at $z = 0, 1, t, \infty$. Concerning its solutions, Maier-completing a task initiated by Heun [96] himself-has shown that there are 192 solutions that can be generated using the symmetries of $D_4$ [131]. For specific parameter values the GHE reduces to other well-known equations of importance: e.g. setting $\epsilon = 0, \gamma = \delta = 1/2$ yields the Mathieu equation [42].

In order to produce an integral representation for Heun's functions, remark that the GHE (and, in fact, all other Heun equations) take the form

$$ y''(z) - B_1(z)y'(z) - B_2(z)y(z) = 0, \quad (4.3.2) $$

We may thus focus on obtaining the integral representation of the solution of the above equation in terms of integrals involving $B_1$ and $B_2$. To this end, remark that we may transform the ODE Eq. (4.3.2) into a matrix equation:

**Proposition 4.3.1.** Let $y(z)$ be a solution of Eq. (4.3.2) with initial conditions $y(z_0) = y_0$ and $y'(z_0) = y'_0$. Let

$$ M(z) = \begin{pmatrix} 1 & 1 \\ B_1(z) + B_2(z) - 1 & B_1(z) - 1 \end{pmatrix}, \quad (4.3.3) $$

and let $U(z, z_0) := Pe^{\int_{z_0}^z M(\zeta)d\zeta}$ be the path-ordered exponential of $M$, that is the solution of $\dot{U}(z, z_0) = M(z)U(z, z_0)$. Then

$$ y(z) = y_0 U_{11}(z, z_0) + (y'_0 - y_0) U_{12}(z, z_0). $$

Hence the problem of determining any Heun function is equivalent with that of determining the path-ordered exponential of a (rather simple) matrix. Using the results of Chapter 3 and the path-sum theorem, one gets as a consequence:

**Corollary 4.3.1.** Let $H_G(z)$ be solution of the General Heun Equation,

$$ \frac{d^2 H_G(z)}{dz^2} + \left[ \frac{\gamma}{z} + \frac{\delta}{z-1} + \frac{\epsilon}{z-t} \right] \frac{dH_G(z)}{dz} + \frac{\alpha\beta z - q}{z(z-1)(z-t)} H_G(z) = 0, $$

with initial conditions $H_G(z_0) = H_0$ and $\dot{H}_G(z_0) = H'_0$, assuming that $z_0 \in \mathbb{R}$ is not a singular point of $H_G$. Denote $I$ the largest real interval that contains $z_0$ and does not contain any singular point of $H_G$. Then, for any $z \in I$,

$$ H_G(z) = H_0 + H_0 \int_{z_0}^{z} G_1(\zeta, z_0) d\zeta + (H'_0 - H_0) \left( e^{z-z_0} - 1 + \int_{z_0}^{z} (e^{z-\zeta} - 1) G_2(\zeta, z_0) d\zeta \right), $$

where $G_i = (1_i - \tilde{K}_i \Theta)^{i-1} = \sum_{n=1}^{\infty} \tilde{K}_i^{*n} \Theta$ and

$$ \tilde{K}_1(z, z_0) := 1 + e^{-z} \int_{z_0}^{z} \frac{\zeta \gamma (\zeta_1 - 1)^\delta (t - \zeta_1)^\epsilon}{z^\gamma (z-1)^\delta (t-z)^\epsilon} e^{\zeta_1} \left( \frac{q - \alpha\beta\zeta_1}{(\zeta_1 - 1)\zeta_1(\zeta_1 - t)} - \frac{\epsilon}{t - \zeta_1} - \frac{\gamma}{\zeta_1} - \frac{\delta}{\zeta_1 - 1} - 1 \right) d\zeta_1, $$

$$ \tilde{K}_2(z, z_0) := \left( \frac{q - \alpha\beta z}{(z-1)z(z-t)} - \frac{\epsilon}{t-z} - \frac{\gamma}{z} - \frac{\delta}{z-1} - 1 \right) e^{z-z_0} - \frac{q - \alpha\beta z}{(z-1)z(z-t)}. $$

In particular both $\tilde{K}_1$ and $\tilde{K}_2$ are bounded over $I$, implying convergence of the $\star$-Neumann series for $G_1$ and $G_2$ for all $z \in I$.

**Example 4.3.1** (Elementary integral series converging to a general Heun function). In order to illustrate concretely the above corollary, consider the following General Heun equation (here with arbitrary parameters),

$$ \frac{d^2 H_G(z)}{dz^2} + \left[ \frac{2}{z} + \frac{7}{z-1} + \frac{(-1)}{z-4} \right] \frac{dH_G(z)}{dz} + \frac{(3/2)z - 1}{z(z-1)(z-4)} H_G(z) = 0, \quad (4.3.4) $$

with initial conditions $H_G(6) = H'_G(6) = 1$. Here, the largest real interval containing 6 and none of the singular points 0, 1 and $t = 4$ is $I = ]4, +\infty[$. Thus Corollary 4.3.1 indicates that for any $z \in ]4, +\infty[$,

$$ H_G(z) = (\Theta \star G_1)(z, z_0) = 1 + \int_{z_0}^{z} G_1(\zeta, z_0) d\zeta $$
with $G_1 = (1_* - \tilde{K}_1\Theta)^{-1}$ and the kernel $\tilde{K}_1$ given by

$$ \tilde{K}_1(z, z_0) = 1 - e^{-z} \frac{(z-4)}{z^2(z-1)^7} \int_{z_0}^{z} e^{\zeta_1} \frac{\zeta_1 (\zeta_1 - 1)^6}{2 (\zeta_1 - 4)^2} (2\zeta_1^3 + 10\zeta_1^2 - 67\zeta_1 + 14) d\zeta_1. $$

All the other types of Heun functions are treated similarly [75]. Of particular importance is the solution of the Confluent Heun equation which is that modelling the event horizon of a Kerr black-hole. Since the singular points of this equation are located at $z = 0, 1, +\infty$, given any initial conditions for $H(z_0)$ and $\dot{H}(z_0)$ at $z_0 \in ]1, +\infty[$, the $\star$-Neumann series representation for the confluent Heun function resulting from the path-sum theorem is guaranteed to converge on the entire domain $]1, +\infty[$, that is from the black hole horizon up to spatial infinity. This crucial property stands in stark contrast with the hypergeometric and Coulomb series, which converge close to 1 and to $+\infty$, respectively. Algorithmic implementations of the novel representations of the Heun function reported are discussed in [20].
# Part 2: Extending number theory to walks on graphs
# Chapter 5

## Walk theory: combinatorics on hike monoids

### 5.1 On the necessity of a theory of walks

The precise nature of the relation between graphs and their walks has, to the best of our knowledge, not yet been thoroughly scrutinised. It seems to be often assumed that walks are slave objects to graphs, in that once the graph is specified its walks can be studied and, in principle, be perfectly known. As a corollary it is expected that properties of the graph leave indelible imprints on its walks, from which the former can thus be inferred. These simple arguments mask the subtle nature of the relation between walks and graphs. For example, consider the following bidirected, vertex-transitive, hence regular, bipartite graph

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

and let $W_{G:○→○'}$ be the set of all walks from any vertex $○$ to any vertex $○'$ on $G$. Now consider the following two digraphs

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

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

Remark that $G'$ is directed and neither vertex-transitive nor regular. Yet, the set of all walks from any gray vertex to any gray vertex on $G'$ is in bijection with the corresponding set $W_{G:○→○'}$ on $G$. The correspondence between both walk sets is simple: the white transient vertices of $G'$ do not add any new cycle to the digraph as compared to $G$ and so only cause the lengths of all of these walks to be multiplied by 2. The structure of every individual walk, i.e. the way it is composed of cycles and of at most one open simple path, is clearly preserved. Rigorously, this makes the bijection into an isomorphism between monoids on walks and walk-like objects, called hikes, something we discuss in details later. Similarly, sets of walks between gray vertices on $G''$ are in bijection with the corresponding sets on $G$. In this case however the bijection does not act plainly on the length of individual walks but it continues to be a monoid isomorphism preserving the internal structure of walks. Remarkably here, $G''$ is not even bipartite. That is, from the point of view of graph theory the transformations from $G$ to $G'$ and from $G$ to $G''$ are nontrivial and subtle. They correspond to the loss of respectively three and four fundamental graph properties namely bidirectedness, vertex-transivity, regularity and bipartiteness. Yet walks sets are essentially unchanged by these alterations.

At the very least, these observations suggest a rather loose relationship between graph properties and walk properties. It raises the question of whether there are more graph transformations-excluding graph automorphisms-that preserve walk sets, and if so are they all somewhat trivial, as above with the addition of transient vertices? This question indirectly sprung up in the fields of applied network analysis and machine learning. For example it was noted that vertex centralities, mathematical quantities designed to grasp the relative importance of nodes in a graph, sometimes fail to do so meaningfully [169, 153]. They can for instance
predict that an outlying vertex is just as 'central' as another one which appears to be at the heart of a graph:

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

On this digraph with adjacency matrix $A$, the gray and black vertices are deemed to be equally central by the subgraph centrality measure $(e^A)_\bullet = (e^A)_\bullet$, a widely used index of vertex importance in real-world networks [50]. Similarly, measures aimed at grasping the degree of “similarity” between networks are prone to failures. In the field of machine learning these measures, called graph-kernels, are used primarily for tasks of automatic graph classification. In particular it is known that walk-based graph-kernels suffer from the existence of differing graph structures with similar walk counts or arrangement of walks, hindering the automatic distinction of these graph structures, see e.g. [170, 130, 113] and references therein. To illustrate this observation, consider the following four directed graphs

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

While digraphs $G_1$ and $G_2$ might seem most similar with one another among all four, sets of walks between any pairs of gray vertices on graphs $G_3$ and $G_4$ are in bijection with sets of walks on $G_1$. As in the earlier example, this bijection preserves the structure of individual walks. At the opposite, there is no such correspondence between the set of all walks on graph $G_2$ and those on $G_1$.

When questioned, these failures have at times been put down to the underlying non-rigorous notions of “importance of a node” and of “similarity between graphs” as we *intuitively* understand them. It is argued that such an understanding, imprecise and fraught with preconceptions, is difficult to express rigorously in mathematical terms thus leading to seemingly absurd results. But the failures could also run deeper and be yet more manifestations of the misappreciated mathematical relation between graphs and walks. Indeed almost all of the measures proposed so far for quantifying such notions as network centrality and graph similarity are algebraic quantities which translate into statements about closed walks. If indeed the cause of the problem is that graphs and their walks give relatively little or at least indirect information about one another, then this issue should be made mathematically rigorous. Considering in particular the case of closed walks on strongly connected digraphs, is perfect knowledge of the former sufficient to determine fundamental properties of the later, possibly up to some triviality?

We answer the preceding question in the negative: not only are there many transformations between strongly connected digraphs that profoundly alter their properties yet leave their sets of closed walks invariant, but these transformations are diverse and far from trivial. Some produce bijective mappings between walk sets that not only preserve their monoidal structure as in the above examples, but also the length of all closed walks. Algebraic quantities routinely used to characterise graphs are left invariant under such mappings. A complete classification of all walk-preserving graph transformations seems to be particularly arduous. The problem finds a wider context in the assertion that a digraph can be reconstructed based solely on structural information about its simple cycles, something which, we find, fails to hold in many a strange way. The complementary question, namely deciding whether a digraph exhibiting certain structural relations between simple cycles exists at all turns out to be unexpectedly difficult. We show at least that this question is decidable and give a meaning for and examples of algebraically closed sets of walks and walk-like objects that cannot, by themselves, be drawn on digraphs. In one more unexpected result we find that these ‘undrawable’ walks do exist on larger digraphs where they are accompanied by a host of algebraically unrelated walks.

In order to talk about closed walks and digraphs in a clearly separate manner so as to untangle their thorny relationship, it is necessary to have a mathematical language for describing sets of closed walks and walk-like objects independently from the digraphs that sustain them.
## 5.2 Hike monoids and number theory

Several schools of thoughts have emerged from the literature concerned such a task, aiming at studying walks as algebraic objects. Among the numerous structures proposed over the years are those based on walk concatenation such as the path-algebra [29] and later, nesting [76] or the cycle space [43]. A promising approach using the theory of partially commutative monoids (also called trace monoids), consists in viewing the arcs (i.e. the directed edges) of a graph as letters forming an alphabet and walks as words on this alphabet. A crucial idea in this approach, proposed by [32], is to define a specific commutation rule on the alphabet: two arcs commute if and only if they initiate from different vertices. This construction yields a semi-commutative monoid which allows for a great flexibility in the walk structure while preserving the ability to distinguish between different walks composed of the same arcs. A remarkable consequence of this construction is the existence of a stable subset of traces, formed by collections of cycles: the hikes. More specifically, hikes constitute a simplified trace monoid that carries most of the information pertaining to the graph structure and, in the case of undirected graphs, all the information. In this trace monoid, the simple cycles form the alphabet while the independence relation is characterised by vertex-disjointness.

Of fundamental importance for the trace monoid of hikes is the hitherto underappreciated prime-property satisfied by the simple cycles. Recall that an element of a monoid is prime if and only if, whenever it is factor of the product of two elements, then it is a factor of at least one of the two. The importance of the prime property lies in that because of it, the partially ordered set formed by the hikes ordered by divisibility is host to a plethora of algebraic relations in direct extension to number theory. This includes identities involving many more objects beyond the well-studied Möbius function [32, 173], such as the von Mangoldt and Liouville functions. In this respect hikes are natural objects to consider, as most of their algebraic properties follow from analytical transformations of the weighted adjacency matrix. The study of the algebraic structures associated with hikes is the main subject of the present work. These structures provide an extended semi-commutative framework to number theory from which both well-known and novel relations in general combinatorics are derived as particular consequences.

### Notations for graphs and walks

While we begin by recalling standard definitions for graphs, we introduce somewhat less common concepts for walks, of which we advise the reader to take special notice.

A graph $G = (V, E)$ is a finite set of vertices $V$ and a finite set $E$ of distinct paired vertices, called edges or arcs, denoted $\{i, j\}$, $i, j \in V$. A digraph $G = (V, E)$ is a finite set of vertices $V$ and a finite set $E \subseteq V^2$ of directed edges (or arcs), denoted $(i, j)$ for the arc from $i$ to $j$. A directed multigraph (or multidigraph) is defined the same way as a digraph, except that $E$ is a multiset. An edge of $E$ is then denoted $(i, j)_k$, the integer $k$ specifying which edge from $i$ to $j$ we consider. We denote by $A_G$ or simply by $A$ the adjacency matrix of $G$ defined as $A_{ij} := n$ with $n \ge 0$ the number of directed edges $(i, j)_k \in E$ from vertex $i \in V$ to vertex $j \in V$.

The labeled adjacency matrix $W = (w_{ij})_{i,j=1,...,N}$ of the graph $G$ is built by attributing a formal variable $w_{ij}$ to every pair $(v_i, v_j) \in V^2$ and setting $w_{ij} = 0$ whenever there is no arc from $v_i$ to $v_j$ in $G$. In the sequel, we identify each arc of a digraph $G$ to the corresponding non-zero variable $w_{ij}$.

We say that $G' = (V', E')$ is a subgraph of a multidigraph $G$, denoted by $G' \subseteq G$, if $V' \subseteq V$ and $E' \subseteq E$. We say that $G'$ is an induced subgraph of a digraph $G$ if furthermore we have $E' = V'^2 \cap E$, that is $G'$ is obtained from $G$ by deleting some of its vertices and only the edges adjacent to them.

A rooted walk, or rooted path, of length $\ell$ from vertex $i$ to vertex $j$ on a multi directed graph $G$ is a contiguous sequence of $\ell$ arcs starting from $i$ and ending in $j$, e.g. $w = (i, i_1)_{k_1} (i_1, i_2)_{k_2} \cdots (i_{\ell-1}, j)_{k_\ell}$ (a sequence of arcs is said to be contiguous if each arc but the first one starts where the previous ended). The rooted walk $w$ is open if $i \neq j$ and closed otherwise, in which case it is also called rooted cycle. A rooted cycle $(i_0, i_1)_{k_1} (i_1, i_2)_{k_2} \cdots (i_{\ell-1}, i_0)_{k_\ell}$ of non-zero length for which all vertices $i_t$ are distinct is said to be simple or self-avoiding. A self-loop $(i, i)_k$ is considered a rooted simple cycle of length one. On digraphs we may also represent walks unambiguously as ordered sequences of vertices $w = i, i_1, \cdots, i_{\ell-1}, j$.

Discarding the piece of information regarding the position of the root turns rooted simple cycles into simple cycles. An induced cycle is a simple cycle for which no pair of visited vertices is linked by an arc that does not belong to the cycle.

A strongly connected component of a digraph $G$ is a maximal digraph $G' \subseteq G$ such that for every pair of vertices $v_1, v_2$ in $G'$, there is a rooted path in $G'$ from $v_1$ to $v_2$. A digraph $G$ is then said to be strongly connected if it is its sole strongly connected component.
## Partially commutative structure on the arcs

Let $\Sigma$ be an alphabet then its Kleene star $\Sigma^*$ designates the set of all finite words on the letters of $\Sigma$. Now let $\mathcal{I} \subseteq \Sigma^2$ be a set of pairs of letters. This set defines a rule, called independence relation, which affirms that pairs of letters in $\mathcal{I}$ are independent and can be commuted when they occur next to each other in a word. Thus $\mathcal{I}$ induces an equivalence relation $\sim_{\mathcal{I}}$ between words $w_1, w_2 \in \Sigma^*$ with $w_1 \sim_{\mathcal{I}} w_2$ if and only if it is possible to pass from $w_1$ to $w_2$ by commuting adjacent pairs of independent letters. Then the trace monoid $\mathcal{T} = \Sigma^*/\sim_{\mathcal{I}}$, is the free partially commutative monoid formed by the $\sim_{\mathcal{I}}$ equivalence classes on $\Sigma^*$. These classes are generically called traces. Different elements of an equivalence class, i.e. trace $t \in \mathcal{T}$ will be referred to as representations of a trace.

At the heart of the trace monoid $\mathcal{T}$ is the partially commutative structure induced by $\mathcal{I}$. This structure is best represented as a graph, called the dependency graph $H$ of $\mathcal{T}$. It is the graph for which vertices represent all letters of the alphabet $\Sigma$ and two vertices $v_i$ and $v_j$ are joined by an undirected edge if and only if letters $i, j \in \Sigma$ are not allowed to commute per $\mathcal{I}$. Formally, $H$ is the complement graph of $(\Sigma, \mathcal{I})$.

Trace monoids where introduced by Pierre Cartier and Dominique Foata in their quest for a purely combinatorial proof of MacMahon's master theorem [32, 127]. They considered more specifically the trace monoid $\mathcal{T}_{CF}$ on the alphabet of labeled directed edges $\Sigma = E$ of a digraph $G = (V, E)$ with independence relation

$$ \mathcal{I}_{CF} := \{ \{(i, j), (k, l)\} : i \neq k \}, $$

where $(i, j) \in E$ and $(k, l) \in E$ are arcs of $G$. This rule implies that two directed edges with different starting points are allowed to commute. For example, words $(1, 2)(2, 3)(1, 4)$ and $(2, 3)(1, 2)(1, 4)$ belong to the same Cartier-Foata trace, while $(2, 3)(1, 4)(1, 2)$ belongs to a distinct trace due to the forbidden commutation of $(1, 4)$ with $(1, 2)$ required to pass from the former trace to the latter. The trace monoid $\mathcal{T}_{CF}$ is sometimes called the Cartier-Foata monoid. In this monoid, a walk (that is a sequence of contiguous arcs) may have non-contiguous representations. For instance, the walk $w_{12}w_{23}$ from $v_1$ to $v_3$ can be rewritten as $w_{23}w_{12}$ since $w_{23}$ and $w_{12}$ start from different vertices. In fact, an open walk always has a unique contiguous representation, as any allowed permutations of arcs would break the contiguity. Surprisingly, the uniqueness of the contiguous representation no longer holds for closed walks. This consequence is an important feature of the partially commutative structure on the arcs: two closed walks starting from different vertices define the same object if they can be obtained from one another by permuting arcs with different starting points.

To illustrate this statement, consider the example pictured in Figure 5.1. The only closed walk starting from $v_1$ that covers every arc exactly once is $c_1 := w_{12}w_{23}w_{34}w_{45}w_{53}w_{31}$. Since the only non-commuting arcs are $w_{31}$ and $w_{34}$, the cycle can be rewritten as starting from $v_3$ by $c_1 = w_{34}w_{45}w_{53}w_{31}w_{12}w_{23}$. On the other hand, there are two closed walks starting from $v_3$ covering every arc once, one is $c_1$ and the other is $c_2 := w_{31}w_{12}w_{23}w_{34}w_{45}w_{53}$. One cannot go from $c_1$ to $c_2$ without permuting $w_{31}$ and $w_{34}$, thus $c_1$ and $c_2$ are different elements in $\mathcal{T}_{CF}$. Here, $w_{12}w_{23}w_{31}$ and $w_{31}w_{12}w_{23}$ are equivalent representations of the same cycle since the permutations of the arcs to go from one to the other are allowed. More generally, the starting vertex of a simple cycle never influences its value.

<figure><img src="image_13.png" /><figcaption>Figure 5.1: The closed walks $c_1 = w_{34}w_{45}w_{53}w_{31}w_{12}w_{23}$ and $c_2 = w_{31}w_{12}w_{23}w_{34}w_{45}w_{53}$ are different although composed of the same arcs. Both are achievable starting from $v_3$ but only $c_1$ is achievable from $v_1$.</figcaption></figure>

## Multiplication and factorization of hikes

A closed walk can be characterised as a contiguous sequence of arcs comprising the same number of ingoing and outgoing arcs for each vertex. A hike is obtained upon relaxing the contiguity condition:

**Definition 5.2.1.** A hike is a trace $h = w_{i_1j_1} \cdots w_{i_{\ell}j_{\ell}} \in \mathcal{M}$ for which the numbers of outgoing and ingoing arcs for each vertex are equal. Formally, the indices $i_k, j_k, k = 1, \dots, \ell$ of the arcs of $h$ satisfy,

$$ \forall i = 1, \dots, N, \quad \#\{k : i_k = i\} = \#\{k : j_k = i\}, \qquad (5.2.1) $$

where $\#$ stands for the cardinality.

**Remark.** Hikes have been studied in various forms in the literature. They originated under the name "circuits" in the seminal paper [32], written in french. We here use a different term to avoid confusion, as circuit may refer to other common objects in graph theory. Hikes can also be defined as the particular heaps of pieces formed from the simple cycles of the graph, see [196].
We denote by $\mathcal{H}$ the set of hikes, which is a subset of $\mathcal{T}_{\text{CF}}$. By convention, the trivial walk 1 viewed as the empty sequence is considered to be a hike. We emphasise that, since hikes are elements of $\mathcal{T}_{\text{CF}}$, they obey the partially commutative structure on the arcs: two hikes $h$ and $h'$ are equal if, and only if, $h'$ can be obtained from $h$ by permuting arcs in $h$ with different starting point. In particular, while every closed walk is a hike, a hike is a closed walk only if it has a contiguous representation. Moreover, we say that a hike is self-avoiding if all its arcs are different and commute.

The multiplication of two hikes $h, h'$, simply defined as the concatenation, yields a hike and shall be denoted by $h.h'$ or simply $hh'$ in the sequel. We define hike division as the reverse operation: $d \in \mathcal{H}$ divides $h \in \mathcal{H}$, which we write $d|h$, if there exists $h' \in \mathcal{H}$ such that $h = d.h'$. We shall use the standard division notation

$$h = d.h' \iff h' = \frac{h}{d}.$$

Here the choice of left-division, rather than right-division, is only a matter of convention. Remark that because the multiplication of hikes is not commutative, $d|h$ does not necessarily imply that $h/d$ divides $h$.

A main concern of this chapter is to treat the decomposition of a hike into simple cycles as a *prime decomposition*, seeing the simple cycles as prime factors. Rigorously, an element $p$ of a monoid is prime if and only if, whenever $p$ is a factor of $a.b$, then $p$ is a factor of $a$ or $b$. Thus, in this context, the primes are indeed given by the simple cycles. We emphasise that, because of the lack of commutativity, the prime factors of $h$, i.e. the elements of the prime decomposition, are different from its prime divisors.

An important consequence of the commutation rule on the arcs is that the prime decomposition $h = c_1 \cdots c_k$ is unique up to permutations of consecutive vertex-disjoint simple cycles. Indeed, switching two different consecutive cycles in the prime decomposition $h = c_1 \cdots c_k$ violates the commutation rule as soon as $V(c_i) \cap V(c_{i+1}) \neq \emptyset$. This property highlights that $\mathcal{H}$ forms a sub-monoid of $\mathcal{T}_{\text{CF}}$, whose alphabet is the set of prime hikes (the simple cycles) $\Sigma_{\mathcal{H}} := \{c_1, \dots, c_k\}$ and with independence relation defined by

$$I_{\mathcal{H}} = \{(c_i, c_j) : V(c_i) \cap V(c_j) = \emptyset\}. \qquad (5.2.2)$$

A geometric interpretation of the prime decomposition can be found in Viennot's theory of heaps of pieces [196]. In this case, the simple cycles are pieces "piled up" in such a way that two simple cycles can only be put to the same level if they share no vertex in common. In fact, heaps of pieces provide a geometric construction for the Cartier-Foata clique decomposition of a trace.

**Definition 5.2.2.** The independence graph of $\mathcal{H}$ is the undirected graph whose vertices are the simple cycles of $G$ and with an edge between two simple cycles $c, c' \in \Sigma_{\mathcal{H}}$ if they share no vertex in common. The hike dependency graph of a hike monoid is the complement of its independence graph.

<figure><img src="image_14.png" /><figcaption>Figure 5.2: An example of a digraph $G$ (left), the simple cycles on $G$ (middle) and the corresponding hike dependency graph $H$ (right).</figcaption></figure>

A clique of the independence graph of $\mathcal{H}$ can be identified with a product of pairwise non-intersecting simple cycles, that is, a self-avoiding hike. The Cartier-Foata clique decomposition of a hike $h$ can then be built as follows.

i) The maximal self-avoiding divisor of a hike $h$ (i.e. the product $s(h)$ of its prime divisors), is the first clique in the Cartier-Foata decomposition of $h$.

ii) If $h$ is self-avoiding, then $s(h) = h$ and $h$ is its own Cartier-Foata decomposition. All simple cycles composing $h$ are vertex disjoint so that they can be put to the same level, forming a heap of height 1.

iii) Otherwise, consider a collection of self-avoiding hikes $s_k$, initiated by $s_1 = s(h)$, and setting

$$s_{k+1} = s\left(\frac{h}{s_1 \cdots s_k}\right)$$

until all simple cycles of $h$ are made part of a clique $s_k$. Each clique of the Cartier-Foata decomposition defines a layer in the heap of pieces.
We refer to [111] for a detailed explanation of the representation of heaps of pieces and its equivalence with Cartier-Foata monoids.

In the sequel, $\ell(h)$ represents the length of a closed hike $h$ while the number of elements in its prime-decomposition (counted with multiplicity) is denoted by $\Omega(h)$. By convention, the trivial hike 1 is not prime and thus $\Omega(1) = 0$.

## Hikes incidence algebra

The hikes ordered by division form a locally finite partially ordered set, or poset, which we denote $P_G$. The reduced incidence algebra on this poset is the set $\mathcal{F}$ of real-valued functions on $\mathcal{H}$ endowed with the Dirichlet convolution

$$f * g(h) := \sum_{d|h} f(d)g\left(\frac{h}{d}\right), h \in \mathcal{H}.$$

Here, the sum is taken over all left-divisors $d$ of $h$, including $h$ itself and the trivial hike 1. One verifies easily that the Dirichlet convolution is associative and distributive over addition. However, it is not commutative since $d$ can divide $h$ without it being the case for $h/d$.

**Proposition 5.2.1.** The reduced incidence algebra $(\mathcal{F},*)$ is isomorphic to the algebra of formal series

$$\sum_{h \in \mathcal{H}} f(h)h, f \in \mathcal{F}$$

endowed with hike multiplication.

Important functions of the reduced incidence algebra include the identity $\delta(.)$ equal to one for $h = 1$ and zero otherwise, the constant function $1(h) = 1$, $\forall h \in \mathcal{H}$ or the Möbius function, the inverse of 1 through the Dirichlet convolution. We refer to [173] for a more comprehensive study. It is one of the main results of the present work that many more number-theoretic functions beyond 1 and $\mu$ have generalised analogs in the reduced incidence algebra $(\mathcal{F},*)$ and that these analogs satisfy the same relations as their number-theoretic counterparts, see §5.3.

The next theorem gives the expression of the Möbius function on $\mathcal{H}$. This result is discussed in Remark 3.6 in [32].

**Proposition 5.2.2.** The Möbius function on $\mathcal{H}$ is given by

$$\mu(h) := \begin{cases} 1 & \text{if } h = 1 \\ (-1)^{\Omega(h)} & \text{if } h \text{ is self-avoiding} \\ 0 & \text{otherwise.} \end{cases} \qquad (5.2.3)$$

Proposition 5.2.2 confirms the characterization of $\mathcal{H}$ as the trace monoid generated by the alphabet of simple cycles $\Sigma_{\mathcal{H}} = \{c_1, \dots, c_k\}$ with independence relation defined in Equation (5.2.2) (see Eq. (56), Chapter 2.5 in [177]). The formal series associated to the Möbius function for $\mathcal{H}$ then appears in the identity

$$\det(\mathrm{Id} - W) = \sum_{h \in \mathcal{H}} \mu(h)h, \qquad (5.2.4)$$

where we recall that $W$ denotes the labelled adjacency matrix of $G$. A proof of this identity can be found in Theorem 1 of [163] on noting that for self-avoiding hikes, the concatenation of arcs coincides with the ordinary multiplication. Proposition 5.2.2 thus provides a determinant formula for the Möbius function of $\mathcal{H}$ and the formal series associated to 1 (i.e. the analogue of the zeta function) is obtained via the formal inversion

$$\det(\mathrm{Id} - W)^{-1} = \frac{1}{\sum_{h \in \mathcal{H}} \mu(h)h} = \sum_{h \in \mathcal{H}} h.$$

**Remark (Coprimality).** The Möbius function is multiplicative on vertex-disjoint hikes,

$$V(h) \cap V(h') = \emptyset \implies \mu(hh') = \mu(h)\mu(h'). \qquad (5.2.5)$$

This identity is reminiscent of the multiplicative property of the number-theoretic Möbius function $\mu_N$ for which $\mu_N(nm) = \mu_N(n)\mu_N(m)$ whenever $n$ and $m$ are coprime integers. The fact that (5.2.5) only holds for vertex-disjoint hikes suggests a more general notion of coprimality on $\mathcal{H}$: two hikes are coprime if they share no vertex in common. In particular, coprime hikes have different prime factors, but contrary to natural integers, this condition is in general not sufficient. The two notions of coprimality coincide on a class of graphs where $\mu_N$ is recovered from $\mu$, see §5.3.
## 5.3 Algebraic relations between hikes series

In this section we show that a plethora of number theoretic relations find natural extensions on the trace monoid of hikes. These provide powerful algebraic tools in a novel graph theoretic context and yield further insights into well established results. For example, we find in §5.3 that MacMahon's master theorem and the Dirichlet inverse of totally multiplicative functions over the integers both originate from the same general result about series of hikes. Throughout this section, $G$ designates a digraph and $P_G$ is the poset of hikes on $G$ ordered by divisibility.

**Definition 5.3.1.** The formal series associated to a function $f \in \mathcal{F}$ is defined as

$$Sf(s) := \sum_{h \in \mathcal{H}} e^{-s\ell(h)} f(h)h$$

for $s$ a complex variable. The particular case $f = 1$ is the zeta function $\zeta(s) := S1(s)$.

Recall that because of the lack of commutativity between hikes, Dirichlet convolution typically acts non-commutatively on functions on hikes $g * f \neq f * g$ and thus hike-series also multiply non-commutatively, i.e. $Sf.Sg = S(f * g) \neq Sg.Sf = S(g * f)$. For convenience, we write $Sf/Sg$ for the right multiplication with the inverse: $Sf/Sg = Sf.(Sg)^{-1}$. We begin with two simple relations counting the left divisors and left prime divisors of a hike.

**Proposition 5.3.1.** Let $\tau(h)$ be the number of left divisors of $h \in \mathcal{H}$. Then $S\tau(s) = \zeta^2(s)$. Let $1_p$ be the indicator function on primes and $\omega(h)$ the number of prime divisors of $h$. Then $S\omega(s) = S1_p(s) \cdot \zeta(s)$.

The results follow immediately from combinatorial arguments on the reduced incidence algebra of $P_G$. The functions $\Omega$ and $\omega$ coincide over self-avoiding hikes, unlike their number theoretic counterparts for which $\Omega(n) = \omega(n)$ if and only if $n$ is a square-free integer. This is due to the stronger characterization of the co-primality for hikes, which requires that the prime factors be not only different but also vertex-disjoint. Indeed, if two different prime factors of a hike $h$ do intersect, then at least one of them is not a left-divisor, resulting in $\Omega(h)$ being greater than $\omega(h)$. The parallel with number theory is however accurate if all the simple cycles commute, in which case a prime factor is always a divisor.

### Walk von Mangoldt function

We begin with a hike version of the number theoretic von Mangoldt function.

**Definition 5.3.2.** The walk von Mangoldt function $\Lambda: \mathcal{H} \to \mathbb{N}$ is defined as the number of contiguous representations of a hike that is, $\Lambda(h)$ is the number of possible contiguous rearrangements of the arcs in $h$, obtained without permuting two arcs with the same starting point.

By convention, we set $\Lambda(1) = 0$. Remark that by this definition, $\Lambda(h) = 0$ whenever $h$ is not a walk. Since different contiguous representations of the same walk start from different vertices, the series associated to $\Lambda$ is obtained by

$$S\Lambda(s) = \sum_{h \in \mathcal{H}} e^{-s\ell(h)} \Lambda(h)h = \mathrm{Tr}(e^{-sW} + e^{-2sW} + \dots) = \mathrm{Tr}((\mathrm{Id} - e^{-sW})^{-1}) - N. \quad (5.3.1)$$

The heaps of pieces point of view provides a remarkable characterization of closed walks as heaps of cycles with a unique high-most element (such heaps are called pyramids). Using our terminology, this translates into a (non-trivial) closed walk being a hike with a unique prime right-divisor, whose length is precisely the von Mangoldt function. To go even further, Viennot remarked that each piece in a heap can be associated to the pyramid formed by the pieces below it, including itself, thus revealing the bijection existing between the pieces composing a heap and its sub-heaps containing only one maximal element (this fact is discussed in the proof of Proposition 6 in [196]). In the context of hikes, this signifies that the prime factors a hike $h$ can be put in bijection with the non-trivial walks dividing it. In particular, summing the von Mangoldt function over all divisors of a hike $h$ reduces to summing the length of its prime factors:

$$\forall h \in \mathcal{H}, \sum_{d|h} \Lambda(d) = \ell(h) \iff \Lambda * 1 = \ell \iff \Lambda = \ell * \mu.$$

This powerful observation yields the following result as an immediate consequence.
**Proposition 5.3.2.** The von Mangoldt function is linked to the zeta function by the following relation

$$ S\Lambda(s) = \sum_{h \in \mathcal{H}} e^{-s\ell(h)} \Lambda(h)h = -\frac{\zeta'(s)}{\zeta(s)} \qquad (5.3.2) $$

$$ \text{where } \zeta'(s) := d\zeta(s)/ds = -\sum_{h \in \mathcal{H}} e^{-s\ell(h)}\ell(h)h. $$

Equation (5.3.2), which is none other than the series version of the identity $\Lambda = \ell * \mu$, coincides with the number theoretic formula of the von Mangoldt series obtained as the logarithmic derivative of the zeta function. The indefinite integral with respect to $s$ yields the logarithmic identity,

$$ \log \zeta(s) := \int \frac{\zeta'(s)}{\zeta(s)} ds = \sum_{h \in \mathcal{H}} e^{-s\ell(h)} \frac{\Lambda(h)}{\ell(h)} h, \qquad (5.3.3) $$

also reminiscent of its number theoretic counterpart, with the length of a hike playing the role of the logarithm of an integer. Further analogies with number theory are discussed in Section 5.3.

**Example 5.3.1.** To illustrate the relation $\Lambda = \ell * \mu$, consider the following graph on 4 vertices:

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

Let $p_1$ be the backtrack and $p_2$ the triangle and let us calculate $\Lambda(p_1p_2)$ and $\Lambda(p_2p_1)$ from $\ell * \mu$. Since the left divisors of $p_1p_2$ are 1, $p_1$ and $p_1p_2$, we have

$$ \Lambda(p_1p_2) = \ell(1)\mu(p_1p_2) + \ell(p_1)\mu(p_2) + \ell(p_1p_2)\mu(1) = 0 \times 0 + 2 \times (-1) + 5 \times 1 = 3. $$

We proceed similarly for $\Lambda(p_2p_1)$:

$$ \Lambda(p_2p_1) = \ell(1)\mu(p_2p_1) + \ell(p_2)\mu(p_1) + \ell(p_2p_1)\mu(1) = 0 \times 0 + 3 \times (-1) + 5 \times 1 = 2. $$

Let us now compare these results with a direct calculation of $\Lambda$, by way of counting all the contiguous sequences in the equivalence classes $p_1p_2$ and $p_2p_1$. We find

$$ \begin{aligned} w_{21}w_{12}w_{23}w_{34}w_{42} &\approx p_1p_2, & w_{12}w_{23}w_{34}w_{42}w_{21} &\approx p_2p_1, \\ w_{42}w_{21}w_{12}w_{23}w_{34} &\approx p_1p_2, & w_{23}w_{34}w_{42}w_{21}w_{12} &\approx p_2p_1, \\ w_{34}w_{42}w_{21}w_{12}w_{23} &\approx p_1p_2, \end{aligned} $$

This confirms that $\Lambda(p_1p_2) = 3$ and $\Lambda(p_2p_1) = 2$, as expected.

## Totally additive functions on hikes

The von Mangoldt identity is a particular case of a more general result concerning totally additive functions over hikes. A function $f: \mathcal{H} \to \mathbb{R}$ is said to be totally additive if

$$ \forall h, h' \in \mathcal{H}, f(h.h') = f(h) + f(h'). $$

If the totally additive function $f$ also respects the divisibility order, i.e. $h \le h' \Rightarrow f(h) \le f(h')$, then we say that $f$ is a rank function.

**Proposition 5.3.3.** Let $f$ be a totally additive function over hikes, then

$$ f * \mu(h) = \begin{cases} f(c) & \text{if } h \text{ is a non-trivial walk, with } c \text{ its unique prime right-divisor,} \\ 0 & \text{otherwise.} \end{cases} $$

The proof of this result is similar to that effected in number-theory, relying on properties of Mbius inversion and the totally additive property, see [71]. Now the von Mangoldt identity can be obtained as an application of Proposition 5.3.3 to the length function $\ell: \mathcal{H} \to \mathbb{N}$, which is obviously totally additive. Another notable application concerns a relation between the prime factors counting function $\Omega$ and the indicator function over walks.
**Corollary 5.3.1.** For all $h \in \mathcal{H}$, $\Omega * \mu(h) = 1_w(h)$ where $1_w$ is the indicator function over non-trivial walks.

The proof is omitted as it is a direct application of Proposition 5.3.3 to the totally additive function $\Omega$. Alternatively, it follows from the fact that the number of prime factors of a hike $h$ (counted with multiplicity) equals the numbers of non-trivial walks dividing it, i.e. $\Omega = 1_w * 1$.

Observe that a similar relation holds in number theory, namely that for $n \in \mathbb{N}$,

$$ \Omega_{\mathbb{N}} * \mu_{\mathbb{N}}(n) = \sum_{d|n} \Omega_{\mathbb{N}}(d) \mu_{\mathbb{N}}\left(\frac{n}{d}\right) = \begin{cases} 1 & \text{if } n = p^k \text{ with } p \text{ prime and } k \in \mathbb{N}, \\ 0 & \text{otherwise.} \end{cases} $$

Here, $\Omega_{\mathbb{N}}(n)$ equals the number of prime factors (counted with multiplicity) of $n \in \mathbb{N}$ and $\mu_{\mathbb{N}}$ is the number theoretic Möbius function. This parallel indicates that powers of primes are in fact the “walks” of number theory, being the only numbers with a unique prime divisor, ergo a unique prime right-divisor.

## Totally multiplicative functions on hikes

A consequence of the Möbius inversion between $\Lambda$ and $\ell$, $\Lambda = \ell * \mu$, concerns totally multiplicative functions on hikes $f \in \mathcal{F}$. We say that $f$ is totally multiplicative if

$$ \forall h, h' \in \mathcal{H}, f(h.h') = f(h)f(h'). $$

**Lemma 5.3.1.** Let $f$ be a totally multiplicative function. The inverse of $f$ through the Dirichlet convolution is given by

$$ f^{-1} = \mu f : h \mapsto \mu(h)f(h), \quad h \in \mathcal{H}. \qquad (5.3.4) $$

The proof is obtained by direct calculation. Lemma 5.3.1 is a semi-commutative extension of the number-theoretic result on the inverse of totally multiplicative functions, $f^{-1}(n) = \mu_{\mathbb{N}}(n)f(n)$, $n \ge 0$, with $\mu_{\mathbb{N}}$ the number-theoretic Möbius function, see below for a rigorous statement on the relation between both results. In a similar vein, the above Lemma leads to another result with a strong number-theoretic flavor:

**Corollary 5.3.2.** Let $f$ be a totally multiplicative function on hikes, $F(s) = Sf(s)$ and $F'(s) = dF(s)/ds$. Then,

$$ \frac{F'(s)}{F(s)} = - \sum_{h \in \mathcal{H}} e^{-s\ell(h)} \Lambda(h) f(h) h. $$

An important extension to MacMahon's master theorem stems from the formal series version of Lemma 5.3.1:

$$ Sf(s) = \sum_{h \in \mathcal{H}} e^{-s\ell(h)} f(h)h = \frac{1}{\sum_{h \in \mathcal{H}} e^{-s\ell(h)} \mu(h) f(h)h}. $$

To see this, consider first a weighted version of the graph $G$ where all arcs pointing to a vertex $i$ are given a formal weight $t_i$. The weighted adjacency matrix of this weighted graph is TW, with $T$ the diagonal matrix where $T_{ii} = t_i$. Now observe that a totally multiplicative function on hikes is completely determined by its value on the primes (since $f(hh') = f(h)f(h')$ regardless of the commutativity of $h$ and $h'$). We may therefore consider the totally multiplicative function which associates any prime $p$ with its weight,

$$ f(p) = \mathrm{weight}(p) = t_{i_2} \cdots t_{i_{\ell(p)}} t_{i_1}. \qquad (5.3.5) $$

where $\{i_1, \cdots, i_{\ell(p)}\}$ is the set of vertices visited by $p$. Then, Lemma 5.3.1 yields

$$ Sf(0) = \sum_{h \in \mathcal{H}} f(h)h = \frac{1}{\sum_{h \in \mathcal{H}} \mu(h)f(h)h} = \frac{1}{\det(\mathbf{Id} - \mathrm{TW})}. \qquad (5.3.6) $$

This is the non-commutative generalization of MacMahon's theorem discovered by Cartier and Foata [32]. MacMahon's original result [127] is then recovered upon replacing $W$ by the adjacency matrix $A$, thus attributing the value 1 to every hike.

In general, totally multiplicative functions on hikes do not have to take on the extremely restricted form of Eq. (5.3.5). In these cases Lemma 5.3.1 goes beyond even the non-commutative generalization of MacMahon's theorem. A striking example of this observation is given by the Liouville function. In number theory, this function is defined as $\lambda(n) = (-1)^{\Omega(n)}$, where we recall $\Omega(n)$ is the number of prime factors of the positive integer $n$, counted with multiplicity. We define the walk Liouville function similarly.

**Definition 5.3.3.** The walk Liouville function $\lambda(h) : \mathcal{H} \to \{-1, 1\}$ is defined by $\lambda(h) := (-1)^{\Omega(h)}$, where $\Omega(h)$ is the number of prime factors of $h$, counted with multiplicity.
The series $S\lambda(s) := \sum_{h \in \mathcal{H}} e^{-s\ell(h)}\lambda(h)h$ associated to the walk Liouville function has a remarkably simple expression showing that calculating it is #P-complete on arbitrary graphs.

**Proposition 5.3.4.** *The formal series of the walk Liouville function $\lambda$ satisfies*

$$S\lambda(s) = \sum_{h \in \mathcal{H}} e^{-s\ell(h)}(-1)^{\Omega(h)}h = \frac{1}{\operatorname{perm}(\mathbf{Id} + e^{-s}\mathbf{W})},$$

where $\operatorname{perm}$ designates the permanent.

The permanent of $\mathbf{Id} + e^{-s}\mathbf{W}$ is the series associated to the indicator function on self-avoiding hikes, that is, the absolute value of $\mu$:

$$\operatorname{perm}(\mathbf{Id} + e^{-s}\mathbf{W}) = \sum_{h \in \mathcal{H}} e^{-s\ell(h)}|\mu(h)|h = S|\mu|(s).$$

Hence, the walk Liouville function $\lambda$ is the inverse of $|\mu|$ through the Dirichlet convolution, similarly as its number theoretic counterpart.

## Relation to number theory

The unique factorization of hikes into products of hikes satisfying the prime property is reminiscent of the fundamental theorem of arithmetic. The difference between these two results stems from the non-commutativity of the product operation between hikes. Unsurprisingly then, on a graph where all prime cycles commute, the prime factorization of hikes identifies with that of the integers and the poset $P_G$ becomes isomorphic to the poset of integers ordered by divisibility, which we denote $P_N$.

**Theorem 5.3.1.** Let $G$ be an infinite directed graph formed by a countable union of vertex-disjoint oriented cycles. Then, the hike poset $P_G$ is isomorphic to the poset $P_N$ of natural integers ordered by divisibility. In particular, the reduced incidence algebra of $P_G$, $(\mathcal{F},*)$, is isomorphic to the algebra of Dirichlet series equipped with ordinary multiplication.

The theorem follows from the existence of a bijection between the simple cycles of $G$ and a family of prime numbers. Then all results obtained earlier on hike series and their relations yield valid number theoretic results when applied to the infinite digraph $G$ composed of a countable union of vertex-disjoint simple cycles. For example, on such digraph $G$:

* Two simple cycles are vertex-disjoint if, and only if, they are different. Thus, the Möbius function on $\mathcal{H}$ coincides with the number theoretic Möbius function. More generally, the notion of co-primality which extends over hikes to vertex-disjointness is here equivalent to the usual definition of co-primality over $\mathbb{N}$.

* All prime factors of a hike are also divisors. It follows that $\tau(h) := 1*1(h)$ and $\omega(h) := 1_p*1(h)$ (where $1_p$ is the prime indicator function) give respectively the number of divisors and the number of prime factors of $h$, similarly as for the arithmetic versions of these functions.

* Closed walks take the form $h = p^k$, for $p$ a prime hike (i.e. a simple cycle) and $k \in \mathbb{N}$. Thus, Definition 5.3.2 yields the von Mangoldt function

$$\Lambda(h) = \begin{cases} \ell(p), & \text{if } h = p^k, \ p \text{ prime} \\ 0, & \text{otherwise.} \end{cases}$$

This recovers the number-theoretic von Mangoldt function, provided we identify the length of a hike with the logarithm of an integer.

* More generally, Proposition 5.3.3 shows that, for any totally additive function $f$, $f * \mu$ has its support on non-trivial walks. On a digraph $G$ where all primes commute, the walks are the powers of primes thus recovering the number theoretic version of the result: for all $f: \mathcal{H} \to \mathbb{R}$ totally additive,

$$f * \mu(h) = \sum_{d|h} f(d)\mu\left(\frac{h}{d}\right) = \begin{cases} f(p) & \text{if } h = p^k, \ p \text{ prime} \\ 0 & \text{otherwise.} \end{cases}$$
# Chapter 6

## Number-theoretic sieves on hikes

The “widely open problem of counting Self-Avoiding Walks (SAWs) and Self-Avoiding Polygons (SAPs) on lattices” (quoting Flajolet & Sedgewick [55]) was first conceived from the study of polymer chemistry in 1947 [154]. Mathematically speaking, SAWs are walks that do not self-intersect, which mimic well actual polymer molecules. As the formal study of SAWs and their closed counterpart the SAPs started in earnest, it was quickly realised that self-avoiding objects arise in a wide range of physical and mathematical problems; e.g. as phase boundaries [186] and in percolation clusters [201, 120, 15] or, as B. Mandelbrot observed [134], from the outer frontier of Brownian motion [117]. Because in such models SAWs and SAPs are invariably realised through a random process, the problem of studying and counting them has so far only been attacked with tools from statistical physics and probability theory. Works along these directions have yielded deep insights into renormalisation and conformal mappings such as the relation between self-avoiding curves and the Schramm-Loewner Evolution (SLE) [118, 14], and the result of H. Duminil-Copin and S. Smirnov who proved the value of the connective constant of the honeycomb lattice [46]. One of the most active subfield of this research concerns loop-erased random walks (LERWs) [121, 119, 116], of which walk theory is the deterministic pendant.

We have seen in the preceding chapter that walks obey a semi-commutative extension of number theory where hikes extend the natural integers, self-avoiding hikes extend the square-free integers, walks extend prime powers, i.e. integers of the form $p^k$ with $p$ prime and $k \in \mathbb{N}$ and finally simple paths and simple cycles extend the primes themselves. This last observation implies that in the framework of walk theory, the extension of the Prime Number Theorem (PNT) describes the asymptotic growth of the number of SAPs on regular lattices as their length goes to infinity. It therefore appears possible to work out an alternative deterministic approach to the 70-year old problem of SAP counting by using extending to a semi-commutative setting the number-theoretic tools originally developed for proving the PNT.

In this spirit, we use combinatorial sieves to prove an exact, explicit and compact formula for the asymptotic fraction of all closed walks on any finite or infinite vertex-transitive graph whose last erased loop is any chosen SAP. In stark contrast with approaches based on probability theory, we proceed via purely deterministic arguments relying on Viennot's theory of heaps of pieces seen as a semi-commutative extension of number theory. Our approach sheds light on the origin of the discrepancies between exponents stemming from loop-erased walk and self-avoiding polygon models, and suggests a natural route to bridge the gap between both. Our results are illustrated by calculations on the infinite square lattice.

### 6.1 Asymptotic enumeration of walks by their last erased loop

The main result, which we prove in the remainder of this work is as follows:

**Infinite sieve theorem.** Let $G$ be an infinite vertex transitive graph of bounded degree $\lambda$. Let $p$ be a self-avoiding polyong on $G$. Let $\{G_N^{Tor}\}$ be a sequence of vertex-transitive graphs on $N$ vertices converging to $G$ as $N \to \infty$. Then the fraction of all hikes, i.e. heaps of cycles, which are closed walks whose last erased loop is $p$ is given asymptotically for $N \gg 1$ by

$$ \frac{\alpha^N F_p}{N \lambda^{\ell(p)}} $$

In this expression $\alpha := \lim_{z \to 1/\lambda^{-}} \exp(\int \frac{1}{z}(R(z) - 1)dz)$, $\alpha \in ]0, 1[$, is well defined and $F_p \lambda^{-\ell(p)}$ designates the fraction of all closed walks defined up to translation on $G$ which are multiples of $p$. This fraction is explicitly given by

$$ \frac{F_p}{\lambda^{\ell(p)}} = \frac{1}{\lambda^{\ell(p)+1}} \deg^T \operatorname{adj}(\mathrm{Id} + C_G|_p \cdot \mathbf{B}_p) \cdot \mathbf{1}, \quad (6.1.1) $$

where $\operatorname{adj}(\cdot)$ designates the adjugate operator, $\mathbf{B}_p$ is the adjacency matrix of the graph $G_p$ induced by $p$ and its immediate neighbours on $G$, $\mathbf{1}$ designates the vector full of 1 and $\deg = \operatorname{diag}(\mathbf{B}_p^2)$ is the vector of vertex-degrees