Dynamics of Non-polar Solutions to the Discrete Painlevé I Equation ¹¹1This work was partially funded by NSF under grant DMS-1615921 to N. Ercolani.

This article is concerned with the following discrete-time non-autonomous dynamical system on the plane $\mathbb{R}^{2}$:

where $r$ and $N$ are parameters. Equation (1) defines a sequence of birational maps, $\phi_{n}$ on the $(x,y)$-plane, for $n\geq 1$, which can be extended to $n\leq 0$. For $n\neq 0$ these maps are undefined along the $y$-axis; however, this does not affect our principal considerations, as we will explain. In many areas of mathematics (1) is well known as the discrete Painlevé I equation (or just dP1). There is a vast literature on the relevance of Painlevé equations, both continuous and discrete, to various subfields of mathematics and we refer the reader to [CM20], [KNY17], and [VA18] for recent reviews of many of these connections.

Painlevé equations in general were originally characterized by the special, restricted behavior of the singularities that their solutions may have. Indeed, in the autonomous limit of (1), which we may specify by setting $N=\alpha^{-1}\,n$ for some non-zero real constant $\alpha$, the system is analytically completely integrable, by which we mean that its solutions lie on closed curves and may be explicitly written in terms of meromorphic functions. For instance, in the case of autonomous dP1, solutions are expressed in terms of rational functions on an elliptic curve or one of its degenerations [BE20]. For general continuous Painlevé equations, this translates to solutions being expressible in terms of functions whose global analytic continuations are restricted only to have poles as singularities. This criterion is known as the Painlevé property. For general discrete Painlevé equations, this property can be re-expressed dynamically in terms of singularity confinement: solutions may become arbitrarily large in a finite number of steps before returning to a bounded region of the plane [GRP91, LG03]. That number of steps is fixed for almost all solutions and excursions to infinity can only take place along a fixed set of directions in the plane.

System (1) is not known to be analytically completely integrable and the existence of a closed form analytic expression of its solutions remains an open problem. Nevertheless, dP1 continues to possess many if not all of the well-known properties typically associated with integrable systems, such as singularity confinement of its generic solutions [GRP91, OTGR99, GR04], zero algebraic entropy [BV99, OTGR99], as well as Lax pair [FIK91, PNGR92] and Hirota birational [KNY17] formulations. We will refer to solutions that leave the first quadrant and exhibit singularity confinement as polar solutions. There are, however, certain initial conditions leading to solutions that increase without bound for all time. We will refer to these as non-polar solutions, similar to the pole-free terminology originally used in [Jos97]. We do not claim that this dichotomy is exhaustive although numerical explorations (a few examples of which are provided in this manuscript) suggest this is the case. The primary focus of this paper is on non-polar solutions.

An important example of a non-polar solution, which originally motivated our interest in this work, stems from the connection of dP1 to approximation theory, more explicitly to analyzing the structure of orthogonal polynomials with exponential weights. Recall that a family of orthonormal polynomials $\{p_{m}\}$ associated to an exponential weight $w(\lambda)=\exp(-V(\lambda))$ is a complete basis of polynomials for the weighted $L^{2}$ space such that

It is known that such polynomials are algebraically specified as solutions to a three-term recurrence relation with real coefficients. When $V$ is even, this recurrence takes the form

and the recurrence coefficients $\{b_{n}\}$ themselves are solutions to a nonlinear difference equation (see for instance [VA18]). In the particular case where $V(\lambda)=N\left(\frac{1}{2}\lambda^{2}+\frac{r}{4}\lambda^{4}\right)$ [Fre76, Mag99, VA18], this difference equation reads

which, setting $x_{n}=b_{n}^{2}$, is precisely equivalent to the dP1 system (1). Freud [Fre76] went on to pin down that $x_{n}\propto n^{1/2}$, using properties of orthogonal polynomials. We call the corresponding orbit $\{x_{n}=b_{n}^{2}\}$, the Freud orbit. Subsequently, Máté, Nevai and Zaslavsky [MNZ85] built on Freud’s work to show the existence of a full asymptotic expansion of the orthogonal polynomial solution in powers of $n^{1/2}$, although they did not describe the coefficients in that expansion. Later Ercolani, McLaughlin, and Pierce [EMP08] developed an alternative approach to the asymptotic expansion of the recursion coefficients using Riemann-Hilbert analysis. This work provided a geometric characterization of the recursion coefficients, in terms of a graphical enumeration problem related to diagrammatic expansions in mathematical physics. Since $x_{n}=b_{n}^{2}>0$, Freud’s special orbit remains in the first quadrant for all time. This feature motivated Lew and Quarles [LQ83] to seek all solutions of

that possess this property. For $x_{n}\neq 0$, (4) is equivalent to (1) without the term in $1/r$ on the right-hand side of the equation for $x_{n+1}$. In [LQ83], Lew & Quarles established, by an elegant contraction mapping argument, that there is a one parameter family of such non-polar solutions. According to [ANSVA15], addressing the existence and uniqueness of the positive solution of (4) with $y_{1}=0$ was a problem posed by Nevai, and independently solved by him in [N83]. In Appendix B, we extend the Lew-Quarles construction to describe a broad class of non-polar orbits of dP1, which we refer to as the Lew-Quarles orbits.

In this paper, we adopt a dynamical systems perspective for (1). Doing so enables us to simultaneously describe and compare polar and non-polar solutions, thereby illuminating novel features of the fuller class of non-polar solutions of dP1. To this end, we use the Painlevé property to complete the phase space at infinity with an asymptotic change of variables based on the Riemann-Hilbert scalings used in [EMP08]. A further transformation allows us to define the asymptotic change of coordinates $(x,y,n)\to(s,f,u)$,

which reveals the existence of an invariant plane at infinity ($u=0$). The dynamics in this invariant plane, which organizes the various asymptotic behaviors of the full system, is shown in Figure 1 and corresponds to the reduced discrete system

The picture we will develop is that the Freud orbit, when extended to negative discrete times, is a singular heteroclinic connection from the origin (which we call $P_{-\infty}$) of the $(s,f,u)$ space to the fixed point $P_{\infty}$ with coordinates $s=f=2$ and $u=0$. In addition, the Lew-Quarles orbits form a family of non-singular heteroclinic connections that converge to the Freud orbit, both as $n\to-\infty$ and $n\to\infty$. The asymptotic points $P_{-\infty}$ and $P_{\infty}$ lie in the plane $u=0$ and are marked as black dots in Figure 1. Given the phase plane structure suggested by this figure, heteroclinic connections from $P_{-\infty}$ to $P_{\infty}$ are trajectories that “escape” into the third dimension before asymptotically returning to the invariant plane. Our goal is to characterize the structure of these connections.

The rest of this article is organized as follows. Section 2 defines the Freud orbit, explains how it is initialized, and presents a high precision numerical simulation of forward and backward iterates from the selected initial condition, illustrating convergence to $P_{-\infty}$ as $n\to-\infty$, and $P_{\infty}$ as $n\to\infty$. Section 3 introduces the change of variables that transforms system (1) into an autonomous system in $(s,f,u)$ coordinates, and reviews the basic properties of the resulting discrete dynamical system. Section 4 presents a detailed study, both analytical and numerical, of the Lew-Quarles orbits, and characterizes these solutions as heteroclinic connections between $P_{-\infty}$ and $P_{\infty}$. In particular, we provide strong evidence that these trajectories live on the center and center-stable manifolds of these two points respectively, and converge exponentially to the Freud orbit as $n\to\infty$. We also build on these results to propose a unique characterization of the Freud orbit as a singular limit of the Lew-Quarles orbits. In Section 5, we prove that orbits that converge to $P_{-\infty}$ and $P_{\infty}$ along specific invariant curves track sequences of points formed by the period-2 points (as $n\to-\infty)$ and fixed points (as $n\to\infty$) of the autonomous dP1 system, in which $n$ appears as a parameter. This system is further described in Appendix C. In addition, we obtain expansions of these orbits in powers of $|n|^{-1/2}$, valid as $n\to-\infty$ and $n\to\infty$. Applying these results to the Freud orbit, as legitimized by our numerical observations, provides a simple means to obtain its asymptotic expansions to arbitrary order in powers of $|n|^{-1/2}$ as $n\to-\infty$ or $n\to\infty$. Finally, Section 6 summarizes our findings and identifies future directions that build on the dynamical perspective introduced in this work. Numerical methods are presented in Appendix D, and asymptotic expansions for the invariant curves near $P_{-\infty}$ and $P_{\infty}$ are given in Appendix E.

As mentioned above, Freud’s asymptotic analysis of the recurrence coefficients for orthogonal polynomials with quartic weight (2) determines a particular solution of (1) with the property that in forward time the orbit always lies in the first quadrant of the phase plane and, further, that $x_{n}\propto\sqrt{n}$. Using the orthonormality of the $p_{n}$ (see Appendix A), one finds that the initial conditions for this orbit are

where the $\mu_{i}$ are the moments:

We will refer to iterates of $(x_{1}^{F},y_{1}^{F})$ under (1) as the Freud orbit and use the superscript $F$ to distinguish this specific set of points. Another form of initial conditions is given by $(x_{2}^{F},x_{1}^{F})$ (See Appendix A).

Figure 2 shows a high-precision numerical simulation of the Freud orbit in the $(x,y,n)$ (top left) and $(s,f,u)$ (top right) coordinates, for $n\in[-224,225]$ and parameter values $r=N=1$. The color scale goes from dark blue (large, negative values of $n$) to dark red (large, positive values of $n$), with lighter colors corresponding to negative and positive values of $n$ that are smaller in magnitude. The bottom row shows projections of the orbit on the $(x,y)$ plane (left), and on the $(s,f)$ plane (right). In the right column, the points, $P_{-\infty}=(0,0,0)$ and $P_{\infty}=(2,2,0)$, as well as their projections on the $(s,f)$ plane, are marked with black dots. In the bottom left panel, the point $(x_{1}^{F},y_{1}^{F})$ is shown in bright green on the $x$-axis with $x_{1}^{F}\simeq 0.47$. The image of $(x_{1}^{F},y_{1}^{F},n=1)$ in the $(s,f,u)$ space (top right panel) is the point of approximate coordinates $(3.1,4.6,-2.1)$ in the bottom right corner of the plot. As detailed in Section 4, we will define the pre-image of $(x_{1}^{F},y_{1}^{F})$ as the point on the $y$-axis of coordinates $x_{0}^{F}=0$, $y_{0}^{F}\simeq-1.5$. In the $(s,f,u)$ space, iterates blow up when $n=0$ but are well defined for $n\neq 0.$ The forward part ($n>0$) of Freud’s orbit is in green, yellow, and red, with $x\simeq y\simeq n^{1/2}$. Its backward iterates (negative values of $n$) alternate between the second and fourth quadrants (points shown in green as well as light and dark blue). In the asymptotic coordinate system (right column), the Freud orbit moves away from $P_{-\infty}$ while alternating between each side of the plane $u=0$, and converges to $P_{\infty}$ as $n\to\infty$. The thin solid curves correspond to Equations (LABEL:eq:CMP1) and (38) (see Section 4), which capture the dynamics near $P_{\infty}$ and $P_{-\infty}$ respectively.

As is clearly illustrated in Figure 2, the $(s,f,u)$ coordinate system provides a natural framework to analyze the properties of the Freud orbit. We now explain the origins of this change of coordinates.

The non-autonomous dP1 mapping (1) may be transformed into an autonomous 3-dimensional system by introducing the variable $\alpha_{n}=n/N$. The choice of denominator here is motivated by a fundamental re-scaling used in the Riemann-Hilbert analysis of [EMP08], which amounts to considering a limit in which $n$ and $N$ go to infinity together at a fixed rate. We first set

The singling out of the coordinates $\theta_{1}=n/(Nx)^{2}$ and $\theta_{2}=y/x$ stems from an alternative approach to extending the phase space at infinity developed in [KNY17], based on a resolution of singularities (see also [Tip20]). The above change of variables leads to

where $\gamma=r/N.$ This system has two fixed points,

corresponding to positive ($P_{\infty};\ \theta_{1}>0$) and negative ($P_{-\infty};\ \theta_{1}<0$) values of $n$. The associated values of $Z$ are $1$ and $-1$, respectively. The eigenvalues and eigenvectors of the linearization of (9) about these fixed points are

We now introduce a change of coordinates centered on the fixed point $P_{-\infty}$ and consistent with the basis of eigenvectors of the linearization of (9) near $P_{-\infty}.$ Specifically, setting

transforms (9) into the following discrete dynamical system

with which we will now work.

System (11) has exactly two fixed points, which in the $(s,f,u)$ coordinates are given by

and therefore lie in the invariant plane $u=0$. Looking for periodic orbits of higher order, we find by direct calculation that there is only one genuine period-2 orbit, given by

and a line of period-3 orbits such that $u_{n}=0$ and $s_{n}+f_{n}=1$, with $s_{n}\neq 0$ and $f_{n}\neq 0$. They are of the form

where $s_{0}\in\mathbb{R}\setminus\{0,1\}.$ Moreover, these are the only period-3 orbits in the invariant plane $u=0$, other than the two fixed points. The points $s_{n}=0$ and $f_{n}=0$ are special, in the sense that they are part of a singular period-3 orbit of the form

Our numerical explorations, in which we set $\gamma=1$, indicate that polar orbits (see examples in Appendix D) track the phase plane structure shown in Figure 1 as soon as $|x|$ is larger than a few units (thus corresponding to $|u|<0.5$). We also observe almost periodic orbits associated with large values of $x$ and $y$ ($|x|,\ |y|>1000$), which in $(s,f,u)$ coordinates correspond to solutions close to (12) on the line $s+f=1$; in that case, $u$ remains small ($|u|<10^{-3}$) and a slight drift is observed across the line $s+f=1$. Polar orbits for which there exists a value of $n$ such that $x_{n}={\mathcal{O}}(\epsilon)$ and $y_{n}={\mathcal{O}}(1)$ display generic singularity confinement: iterates of $x_{n}$ become large before returning to values of order $y_{n}$. Specifically,

In $(s,f,u)$ coordinates, this translates to

as $\epsilon\to 0$. In other words, during singularity confinement, iterates of dP1 visit neighborhoods of the points $(1,0,0)$ and $(0,1,0)$ of the singular period-3 orbit (13).

The linearization about the line of period-3 orbits (12) is given by

for $s_{n}=s_{0}+\mu_{n}$, $f_{n}=1-s_{0}+\nu_{n}$, $u_{n}=O(\epsilon)$, $\mu_{n}=O(\epsilon)$, and $\nu_{n}=O(\epsilon).$ The matrix

has eigenvalues $\lambda_{A}=0$ with eigenvector $\zeta_{0}=(0,-1,1)^{T}$ and $\lambda_{A}=1$ with eigenvector $\zeta_{1}=(1,-1,0)^{T}$. In addition, since $\lambda_{A}=1$ has algebraic multiplicity 2 but geometric multiplicity 1, we define the generalized eigenvector $\zeta_{2}$ such that

With these conventions, any initial condition of the form $\zeta=a\,\zeta_{0}+b\,\zeta_{1}+c\,\zeta_{2}$ evolves according to

so that a perturbation transverse to the line of period-3 orbits, of the form $\zeta=a\,\zeta_{0}+c\,\zeta_{2}$ will linearly drift along this line according to $s_{3p}=s_{0}+3\,p\,c-c/3$, $f_{3p}=1-s_{0}-3\,p\,c$, $u_{3p}=0$, and eventually reach the vicinity of the singular orbit (13).

Finally, we note that the change of coordinates from $(s,f,u)$ to $(x,y,n)$ is given by

with $\alpha=n/N$.

We now restrict our attention to the family of Lew-Quarles orbits, whose proof of existence is summarized in Appendix B. The goal of this section is to show that these orbits, which are initiated in the first quadrant, can be extended to heteroclinic connections between $P_{-\infty}$ and $P_{\infty}$, and that the Freud orbit is a singular limit of such solutions. For reference, Figure 3 shows the Lew-Quarles (in brown) and Freud (in color) orbits in the $(x,y)$ plane.

This section provides the conceptual framework for understanding a striking observation: invariant curves near the fixed points of (11) track periodic points of the associated autonomous limits of dP1. This conceptual framework is grounded in a novel process that realizes a simple and elegant mechanism for generating classical asymptotic expansions of the Freud orbit.

In Section 4, we made the conjecture that near $P_{-\infty}$ and $P_{\infty}$, the Lew-Quarles orbits lie on invariant curves whose expansions are given by (38) and (LABEL:eq:CMP1) respectively. We now use this assumption to determine the behavior of $u_{n}$ as a function of $n$ along the Freud orbit near $P_{-\infty}$ and $P_{\infty}$. Near $P_{-\infty}$, we look for a sequence $\{u_{n}\}$ that converges to $0$ as $n\to-\infty$ and is consistent with the last equation of (14). We thus require that

Since per (18) $s_{-\infty}(u_{n})+f_{-\infty}(u_{n})={\mathcal{O}}(u_{n}^{4})$, we see that at dominant order, the iterates $\{u_{n}\}$ solve $\gamma\,n\,u_{n}^{2}-u_{n}+1=0$, which is the equation defining the period-2 solutions of the autonomous dP1 system in $(s,f,u)$ coordinates. As further described in Appendix C, this system is obtained from dP1 by setting $\alpha=n/N$ and then assuming that $\alpha$ is a constant parameter. We call the resulting autonomous map $\alpha$-dP1. Figure 9 of Appendix C.2 illustrates the numerical convergence of the Freud orbit to the sequence of period-two points of $\alpha$-dP1 defined in Equation (35), as $n\to-\infty$. Equation (19) also provides an expansion of $u_{n}$ as a function of $n$,

which is obtained by writing $u_{-\infty,n}$ as a Laurent expansion in powers of $\sqrt{-n}$, substituting into the right-hand equation of (19), and solving term by term. At this point, this expansion is formal because we do not have a proof of the existence of the smooth functions $s_{-\infty}$ and $f_{-\infty}$ that appear in expansions (18), and therefore do not have enough control to bound the remainder in (20). The approach however, is very general. The definition of $\alpha$ in terms of $s$, $f$, and $u$, together with Equation (38), immediately leads to two results: that solutions of dP1 on the invariant curve associated with (38) track the sequence of period-two points of the autonomous dP1 system as $n\to-\infty$, and the expansion (20).

Near $P_{\infty}$, we proceed in a similar fashion, although in that case the existence of the center manifold $\mathcal{C}$ makes the resulting expansions asymptotic. Equation (19) is replaced by

Using (16), we see that at leading order, the iterates $\{u_{n}\}$ solve $\gamma\,n\,u_{n}^{2}+u_{n}-3=0$, which is the equation defining the fixed point solutions of the autonomous dP1 system in $(s,f,u)$ coordinates (see Appendix C). Solving for $u_{n}$ as a function of $n$, we find

We note that control of the big $\mathcal{O}$ term follows from the center manifold theorem that guarantees the existence of all Taylor coefficients of $s_{\infty}$ and $f_{\infty}$, and then a straightforward application of Taylor’s remainder theorem. Plots of $u_{n}$ as a function of $n$ for the numerically estimated Freud orbit, and of the expansions $u_{-\infty,n}^{\pm}$ and $u_{\infty,n}^{-}$ given above are provided in Figure 7; they show excellent agreement, even for values of $n$ near 0. This figure dramatically reinforces the fact that the Freud orbit is singled out among the Lew-Quarles orbits by its singularity at $n=0$.

Given their uniform exponential convergence to the Freud orbit exhibited above, it is natural to further conjecture that all the Lew-Quarles orbits defined in Theorem 1 have this same asymptotic expansion, and so differ only in terms that are beyond all orders with respect to the algebraic gauge $n^{-1/2}$.

By introducing the asymptotic change of coordinates (5), we transformed the dP1 mapping (1) into the 3-dimensional discrete autonomous dynamical system (11), whose two fixed points, $P_{-\infty}$ and $P_{\infty}$, correspond to solutions $\{x_{n}\}$ of dP1 that grow without bounds as $|n|\to\infty$. This transformation, together with a combination of analytical and numerical investigations, allowed us to characterize known non-polar orbits of dP1 as heteroclinic connections between these two fixed points. In particular, we described the Freud orbit as a singular limit of the Lew-Quarles orbits. By understanding how these solutions leave $P_{-\infty}$ to converge to $P_{\infty}$ as $n$ increases, we discovered that they track sequences of points constructed from period-2 (near $P_{-\infty}$) and period-1 (near $P_{\infty}$) points of the autonomous counterpart of dP1. Moreover, our results are consistent with and in many aspects support the conjecture that the Lew-Quarles orbits are on center-stable manifolds of $P_{\infty}$ and, as $n\to\infty$, exponentially converge to the Freud orbit, which itself lives on a center manifold of $P_{\infty}$. The presence of invariant curves that contain the Lew-Quarles (as $n\to\infty$) and Freud orbits provides a method to find explicit expansions of these solutions in powers of $|n|^{-1/2}$ as $|n|\to\infty$. These expansions are asymptotic near $P_{\infty}$ and formal near $P_{-\infty}$.

Compared to their continuous counterparts, many aspects of discrete dynamical systems remain relatively unexplored, especially in the non-autonomous case. In this environment, the study of concrete examples takes on an added value. Due to its analytically completely integrable discrete limit, dP1 is of particular interest. The resulting rich inherent structure therefore provides an auspicious environment for the exploration of explicit features of non-autonomous discrete dynamics that cannot be directly attacked in a general setting. This paper does that in the context of what we have termed non-polar orbits and relates its findings to previously known aspects of Painlevé systems. This approach also leads to interesting applications regarding the behavior of specific orbits as $|n|\to\infty$.

The autonomous dP1 system has orbits that are degenerations of its quasi-periodic elliptic solutions and correspond to trigonometric solutions along a separatrix. This separatrix, and the solutions along it, were examined in detail from a dynamical systems perspective in [BE20] and used to better understand their relation to the combinatorial problem of geodesic distance in the enumeration of planar graphs. That paper also examined a related system, $x_{n}(x_{n+1}+x_{n-1}+1/r)=n/Nr$, which has relevance to the enumeration of labelled trees and super-Brownian excursions. Both of these are instances of the application of dP1 type systems to other areas of mathematics. The present manuscript reveals a subtle connection between the separatrices of autonomous dP1 discussed in [BE20] and Freud’s orbit. Specifically, we discover a novel and deep connection between the mapping (1) and its autonomous counterpart, as the sequence $\{x_{n}\}$ in Freud’s recurrence follows the fixed points of a collection of autonomous mappings. In Appendix C, we provide a geometric description of the important role played by these separatrices in making it possible for the Freud orbit to grow without bounds as a solution of dP1. In addition, the principal message of this article is that Freud’s orbit coincides with a center manifold of a fixed point at infinity. This property of having, effectively, a limiting fixed point at infinity is highly atypical. We do not believe that such a connection between autonomous fixed points and an asymptotic non-autonomous fixed point has been noticed before in the literature.

Another novel discovery is the realization that recursive constructions of invariant curves provide an elegant mechanism for generating asymptotic expansions of non-polar solutions of dP1 as $|n|\to\infty$ (see (20) and (21)). For $n>0$, these have relevance for another combinatorial problem related to random tilings of Riemann surfaces (or geometric foams and quantum gravity in the physics literature). Together with the results of [BE20], we therefore now have three examples of dynamical systems structures (specifically two instances of a stable manifold of a hyperbolic fixed point in [BE20], and a center manifold here) associated with solutions to combinatorial problems. Such connections between different areas of mathematics are highly intriguing. In a subsequent paper we will develop some of these connections using our recursive construction of invariant curves and their associated asymptotic expansions. Methods developed in [Tip20] will enable us to demonstrate the key feature of uniform validity in $r$ and $N$, within appropriate ranges, for these expansions.

As mentioned in the introduction, Máté, Nevai, and Zaslavsky in 1985 [MNZ85] established the existence of an asymptotic expansion for the Freud orbit in powers of $n^{1/2}$, whose leading behavior can be shown to be [Tip20]

We also note that in [Jos97] Joshi carried out a formal Painlevé dominant balance analysis, seeking asymptotic fixed points of dPI in terms of an asymptotic constant of motion. From the perspective of this paper one can see that this formal approach correctly captured the leading order behavior at $P_{\infty}$ and $P_{-\infty}$. It would be interesting to determine if and where this approach might match the “level set” structure evident in Figure 1, and possibly even approximate our exact trajectories in the vicinity of $P_{\pm\infty}$. The mechanisms just referred to in the previous paragraph reproduce and extend this type of asymptotic analysis. Indeed, in these connections, our results may say more about the asymptotics of orthogonal polynomials.

Our subsequent work will also make central use of the realization that the leading order term in the large $n$ asymptotic expansion of [EMP08] in fact coincides with our expression for the fixed points of autonomous dP1. This, together with the analysis developed in Section 5, enables us to devise an efficient scheme for explicitly counting topologically distinct quadrangulations, involving a fixed number of tiles, of compact Riemann surfaces. Other, higher dimensional, Painlevé systems have orbits of Freud type corresponding to degree $2\nu$ potentials that in turn yield information about the enumeration of more general polygonal tilings. For potentials of odd dominant degree traditional methods of orthogonal polynomial theory break down and one must consider generalizations such as non-Hermitian orthogonal polynomials [EP12], which presents obstacles to asymptotic analysis. The particular case of a cubic potential, which relates to topological triangulations of surfaces, is of special interest and involves a non-autonomous planar Painlevé dynamical system. The dynamical systems approach developed in this paper may help to overcome some of these obstacles.

Finally, we note that our work has focused on the dP1 regime with $r>0$. However, there is also interest in the singular regime where $r\to-1/12$, the so-called “Boutroux regime” that has been explored in [JL15] and [DK06]. This is the so-called double scaling limit of random matrix theory, which provides a bridge between the discrete and continuous versions of Painlevé I. A problem of general physical importance is to study and relate the behavior of a correlation function for a statistical mechanical system at small values of a scaling parameter to its behavior at large values. This is called a connection problem. For a class of explictly solvable statistical mechanical problems this is related to connection problems for continuous Painlevé equations. For Painlevé I and II a systematic study of this problem was carried out in [JK92]. A natural extension of our work is to explore the potentially interesting connection problem between the limits $r\to\infty$ and $r\to-1/12$ from a dynamical systems perspective.

As stated at the outset of these Conclusions, the work presented here combines theory and asymptotics with high-precision numerical simulations to arrive at a detailed and compelling picture for the structure of non-polar dynamics in the dP1 system. The main work needed to convert this picture into fully rigorous statements centers on establishing our stated conjectures about center-stable and center manifolds laid out in Section 4. This challenge is of independent interest. Its resolution will have significance for discrete dynamical systems theory broadly speaking. Within that scope particular open questions of importance include

1. 1.

   The rigorous verification of the numerically evident fact that all of the Lew-Quarles orbits converge exponentially to the Freud orbit at a uniform rate as $n\to\infty$.

2. 2.

   A more general conceptual understanding of the unique singularity formation in the Freud orbit as compared to the other Lew-Quarles orbits. In other (continuum) examples of center manifolds (see e.g. [M07] Section 5.6) one observes that singularities form in finite time along those center manifolds whose asymptotic expansions do not have any corrections beyond all orders. Could this be the case for Freud’s orbit?