Capturing the cascade: a transseries approach to delayed bifurcations

Transseries expansions build upon ordinary power series methods by including additional basis elements such as exponentials and logarithms. Alternative summation methods can then be used to"resum"series to obtain more efficient approximations, and have been successfully widely applied in the study of continuous linear and nonlinear, single and multidimensional problems. In particular, a method known as transasymptotic resummation can be used to describe continuous behaviour occurring on multiple scales without the need for asymptotic matching. Here we apply transasymptotic resummation to discrete systems and show that it may be used to naturally and efficiently describe discrete delayed bifurcations, or"canards", in singularly-perturbed variants of the logistic map which contain delayed period-doubling bifurcations. We use transasymptotic resummation to approximate the solutions, and describe the behaviour of the solution across the bifurcations. This approach has two significant advantages: it may be applied in systematic fashion even across multiple bifurcations, and the exponential multipliers encode information about the bifurcations that are used to explain effects seen in the solution behaviour.


Introduction
Transseries are a natural extension to classical asymptotic power series which are used to study systems in which the solution behaviour depends on multiple distinct exponential scales. A transseries represents the solution to a system as the sum of multiple power series, each multiplied by a different exponential prefactor [28]. The value in this approach is that a transseries developed in one region of parameter space can typically be extended into regimes in which the solution depends on different scales, simply be applying different summation methods to the transseries itself, without the need to rebalance the equation terms and apply matched asymptotic expansions to connect the regimes.
The transasymptotic method introduced in [13,14,15] consists of constructing a transseries in terms of some small parameter ε. The transseries terms are then reordered, and higher-order exponential terms at each order of the small parameter are summed. This change of summation order, or "resummation", captures the behaviour of the system in regions where different exponential terms dominate the solution. Transasymptotic methods have been been used to determine the location of moveable poles in Painlevé equations directly from asymptotic solutions [5,15] (see also [16,17]). In this work, we show that these transseries and summation approaches can be used in a systematic fashion to identify complicated bifurcations in discrete systems that typically require careful application of multiple scales [40] or renormalisation methods [7].
In this study, we demonstrate that transseries resummation can be used systematically to accurately capture the behaviour of the solution to a discrete equation containing a periodic-doubling cascade with delayed bifurcations. Delayed bifurcations may occur in dynamical systems where an underlying parameter is itself slowly varying and the solution initially clings to a metastable branch of the solution before eventually jumping to the stable branch. They have been studied widely in systems of ordinary differential equations (see, for example, [55]). These "slow-fast" systems have behaviour occurring on two (or more) distinct timescales, with the solution trajectory remaining near to an unstable solution for a significant distance after stability has been lost; solutions containing this behaviour have been termed "canards".
There has been a significant volume of work studying the asymptotic behaviour of canards in continuous settings, see for example, using composite asymptotic expansions in [9,32,27], and Borel summation methods in [24]. Borel summation methods are closely connected to transseries resummation methods (see [1]), and have been used to study discrete problems, as in [47]. This motivates the idea that transseries resummation methods could be a useful technique for studying delayed bifurcation behaviour; in this study, we will focus on delayed bifurcations appearing in discrete systems, and in particular, singularly perturbed variants of the logistic map.
We will show that period doubling bifurcations depend on the interaction between different exponential factors, and it is therefore advantageous to represent them explicitly using transseries. By expanding in the asymptotic limit, we may determine terms in the algebraic power series to determine the initially stable non-periodic solution. The next step will be to reorder the transseries terms and perform a transasymptotic resummation, which will produce an accurate description of the doubling phenomena. This approach has the additional advantage that it allows us to determine further subdominant exponential scales in the transseries explicitly which dictate subsequent doubling bifurcations present in the solution.
By incorporating a multiple scales ansatz into the transseries expression, we will shown that transseries resummation -which was developed to describe continuous behaviour -can be used to calculate discrete variation without any further analysis to the transseries method.
We study here two variants of the ubiquitous (and generic) standard logistic map y(n + 1) = λy(n) [1 − y(n)] , 0 < y(0) < 1, where λ is a dimensionless bifurcation parameter 0 < λ ≤ 4. This system contains a period-doubling route to chaos, found by allowing the parameter λ to vary. In the range 1 < λ < 3, this system tends to a stable equilibrium without periodic effects. In the range 3 < λ < 1 + √ 6, the system tends to a 2-periodic stable equilibrium. Increasing λ beyond 1 + √ 6 produces systems that tend to stable equilibria with higher periodicity.
The earliest study of the delayed bifurcations in the slowly-varying logistic map is [7], who applied renormalisation methods to derive asymptotic scaling laws for the delays between period doubling, and performed analysis and numerical experiments to determine the location of the bifurcation points. In addition to establishing specific results about the slowly-varying logistic map, this study established that delayed bifurcations can play an essential role in the behaviour of discrete systems. Similar methods were used to study delayed bifurcations in more general unimodal maps [20], as well as discrete maps with noise [6,21,22].
Further studies of this system appeared in subsequent years. In [29,30,32], the existence of canard solutions was rigorously proven in general classes of discrete maps that include the slowly-varying logistic map. Further discussions of canard solutions to both discrete continuous and discrete dynamical systems are given in [31,33].
In more recent years, this system was studied using matched asymptotic expansions and multiple scales methods [40]. The purpose of this previous work was to show that the method of multiple scales could be used to combine a "fast" discrete timescale with a slow time variable that could be treated as continuous, while still capturing the essentially discrete-scale behaviour present in the problem. By carefully balancing terms, the authors were able to identify the bifurcation points and produce accurate asymptotic approximations to the solution behaviour on both sides of the delayed bifurcation. In the works described here, the slowlyvarying logistic equation has provided a useful testing ground for treatments of discrete systems, due to the complicated behaviour that it produces.
It has been shown in [42] that transseries approaches may be used to improve upon asymptotic results obtained using matched asymptotic expansions. In that study, transseries resummation methods were used to obtain a uniform approximation to a continuous problem that had been previously solved using multiple scales methods. The transseries approach was able to naturally incorporate higher exponential terms, and thereby improve on the accuracy of the results, even for values of the perturbation parameter that were not extremely small. Motivated by this result, we will show that transseries resummation can be used to improve on existing multiple scales results in discrete systems.
In section 2 we will first study the "static" logistic map, corresponding to λ = 3 + ε for a fixed choice of (small) real parameter ε. This will be used to introduce and to demonstrate the technical details of the transseries process The second variant, studied in section 3, is the "dynamical" logistic equation, which has a slowly increasing bifurcation parameter λ = 3 + εn. In this case, the bifurcation parameter grows, with different solutions becoming stable and unstable as n increases. In Figure 1, Apart from the pedagogical aspects of demonstrating the applicability of transseries to such problems to recover and extend existing results, in this work we demonstrate four main enhancements.
Firstly, we show that transseries resummation can capture the behaviour of solutions to both a static and slowly-varying logistic map in a systematic and efficient fashion. The multiple scales method used by [40] was able to produce asymptotically valid approximations to the solution, but the process requires a careful The period-doubling cascade is apparent; the transition between non-periodic and 2-periodic behaviour is visible, as is the transition between 2-periodic and 4-periodic behaviour. As the solution continues, it eventually becomes chaotic. The 2-periodic behaviour in the solution begins to contribute immediately, but is not immediately visibly apparent due to the delay in the bifurcation behaviour. Similarly, from the analysis in Section 3, we will determine that solution begins to display 4-periodic behaviour at n ≈ 3455, shown as a red line, but it is also not immediately visibly apparent.
expansion and asymptotic matching each time a bifurcation occurs. The transseries resummation approach used here is systematic, without any need for matching, and can be applied in largely identical fashion to capture each successive bifurcation. Secondly, transseries resummation allows us to capture behaviour when the bifurcation parameter is not necessarily small. In [40], the static logistic map was studied in the limit that the bifurcation parameter was close to three, leading to 2-periodic behaviour in the static map. In this study, the transseries resummation approach allows for the study of larger values of the bifurcation parameter, describing 4-and 8-periodic behaviour in the static map.
Thirdly, we demonstrate that transseries can capture and control the onset of period-doubling behaviour through in terms of exponential weights in the transseries coupled with their resummation. By studying these exponential terms, we are able to determine precisely when higher-periodicity behaviour appears in the solution, and when the bifurcation starts to grow. We will show that transseries resummation can approximate the behaviour of 2-and 4-periodic solutions, and that calculating the transseries exponential terms can explain the onset of higher 4-or even 8-periodic behaviour as the bifurcation parameter grows. In principle the method could be continued in the same fashion to determine this behaviour as the dynamic map sweeps through subsequent bifurcations.
Fourthly, we are able to use transseries summation to significantly improve the approximation accuracy of the solutions over and above that afforded by matched asymptotics in several parameter regimes, with minimal, if any, additional effort. The increase in accuracy is a consequence of the inclusion of multiple exponential scales in the solution approximation, and is most apparent in parameter regimes in which different exponential scales all contribute to the solution behaviour.

Static Logistic Equation
First we consider the static logistic equation, given by We will write the solution as a continuous transseries in terms of ε > 0. We first produce an asymptotic expansion in the limit 0 < ε 1 but as we shall see later, the transseries approach will be used to extend this result to produce an accurate approximation for ε = O(1). We will then show that this continuous transseries is capable of capturing discrete period-doubling effects seen in this system, and approximating higher periodicity behaviour for values of ε that lead to 2-, 4-and even 8-periodic solutions.
In [40], the authors studied the asymptotic behaviour of this system for small ε using multiple scales methods. This showed the manner in which the behaviour approached the 2-periodic stable manifold associated with λ > 3. Using transseries methods, we can extend this approach to consider systems in which ε is not asymptotically small, and demonstrate the manner in which the solution approaches the stable solution for higher periodicities.
In order to first determine the non-periodic and 2-periodic solutions, we ignore the initial condition and solve (2) with the condition that y(n + 2) = y(n). This gives three unique solutions. One solution is non-periodic, and is given by The remaining two solutions are 2-periodic, and are given by For ε > 0, the non-periodic solution is unstable. For 0 < ε < √ 6 − 2, the 2-periodic solution is stable. If ε exceeds √ 6 − 2, the 2-periodic solution is unstable, and the stable solution to the system becomes 4-periodic, and can be identified by solving y(n + 4) = y(n), however this solution cannot be expressed in closed form. Continuing to increase ε leads to the periodicity of the stable solution increasing until chaotic behaviour is eventually obtained.

Transseries ansatz
We begin by applying a transseries ansatz, including a continuous variable x and the small parameter ε. We set x = εn, assuming for now ε 1, and R(x) = y(n) obeys We formulate an ansatz for the solution behaviour in terms of both ε and a transseries parameter σ0. The preliminary ansatz is given as the standard (non-logarithmic) transseries: where βm will be chosen such that Rm,0 takes a nonzero value. Applying this expression to (5) allows for the system to be matched in powers of σ0. At leading order, we find This expression gives the non-periodic manifold, which is stable for 1 < λ < 3. For values of λ greater than three, corresponding to positive ε, we expect to see additional periodic effects emerge from the transseries expression. Continuing to the next order in σ0, we find that By expanding R1 as a power series in ε, (6) gives Expanding as a Taylor series in ε gives A(x + ε) = A(x) + εA (x) + · · · . Matching equal powers of ε on both sides in (9) therefore requires exp(−A (x)) = −1, or The arbitrary constant in A(x) may be absorbed into the transseries parameter σ0, and is therefore set to be zero for convenience. The solution will only be evaluated for n ∈ Z, which corresponds to x/ε ∈ Z. The choice of p has no effect at these values, and we therefore set p = 0 without loss of generality. Using the expressions for R0 in (7) and A in (10) in (5) and tracking orders of σ0, we obtain a recurrence relation for Rm: for m ≥ 1, where we take the convention that the summation term is zero when m = 1. It is straightforward to show by direct substitution that the solution to this recurrence takes the form where the functions Rm(ε) can be computed directly. R1(ε) may be chosen to be an arbitrary constant, which again can be absorbed into σ0. For algebraic convenience, we can select R1(ε) = ε. The subsequent terms may hence be obtained using the recursion (11), which gives the first few terms as . (13) By analysing the form of the recurrence solution, it is straightforward to determine βm in general, giving βm = (m + 1)/2 for m even, and βm = (m + 2)/2 for m even. It can also be seen by direct calculation that

Computing the terms in the resummed transseries
The alternating behaviour of (14) combined with the general form of Rm(x, ε) described in (12) suggests that the ansatz may be conveniently re-written to incorporate these elements explicitly. In particular, the fact that we have two different sub-series in (14) depending on the parity of m suggests that we should split the series into odd and even powers of m. We therefore write the ansatz (6) as where we have switched the order of summation, noting that the exponential terms are all of the same order in ε. Note that in the static case both sums in (15) are convergent. We define a new series parameter τ0, as well as odd and even power series in this parameter, such that The transseries expression is now given by (the x, σ0 dependence is encoded in τ0): We may apply this "resummed" transseries to the logistic equation (5) and equate powers of ε to obtain expressions for Ω o,k and Ω e,k . This process is somewhat technical, and is shown explicitly in Appendix A. Equating terms of order ε and ε 3/2 respectively gives Solving the ordinary differential equation gives Ωo,0(τ0) = ±τ0(C + 9τ 2 0 ) −1/2 . The sign and constant may be chosen arbitrarily, as this choice again can be absorbed into σ0, which is yet to be determined. In order to maintain consistency with (13), we select the positive sign and C = 1, giving Ωe,0(τ0) = − 3 2 Continuing this process and equating terms of order ε 2 and ε 5/2 respectively gives Solving these equations gives where yet again the arbitrary constant of integration can be chosen arbitrarily with the choice being absorbed into σ0, and was therefore selected in order to maintain consistency with (13). In principle, this process can be continued to higher transasymptotic order indefinitely by matching at higher orders of ε.
We can compare these expressions to the multiple scales expansion obtained in [40]. A straightforward comparison shows that the expression Ωo,1 gives the first correction to the composite expansion, denoted P (s), from [40]. The expression for P (s) contains an exponential multiplier e s , which corresponds to the exponential scaling in τ0. This confirms that the resummed transseries identifies the region in which the power series expression breaks down, and allows for the expansion to be continued past this region even without the use of matched asymptotic expansions.

Initial value problem
In order to complete the approximation, we must determine the value of σ0 using the initial condition in (5), which requires that R(x = 0, ε; σ0) = 2/3. Noting that τ0(x = 0, ε) = √ εσ0, this corresponds to solving By expanding σ0 as a power series in ε such that Matching powers of ε in (22) now gives This process may be continued as Ω o,k and Ω e,k are computed for higher values of k, giving the increasingly complete expression for the 2-periodic behaviour.

Transseries ansatz
To obtain the 4-periodic solution requires an adaptation of the previous process. We now take a four-periodic perturbation about the 2-periodic solution obtained in (22). This allows us to form a transseries that can be used to capture solutions which tend to a four-periodic stable manifold. In [40], this would have required solving a challenging multiple scales problem, as the asymptotic solution obtained therein is only valid for small ε. Using the transseries approach, we obtain a significant more general result. We write the solution as a perturbation around the non-period behaviour and the 2-periodic behaviour captured by the transseries expression (17), in terms of the variable τ0, here written asR(τ0, ε).
We will then show that this perturbation starts contributing for values of ε large enough. We can see by direct substitution into (5) that In order to identify the correct scaling for S(x), we note the form of the 2-periodic manifold, given in (4). We represent the 2-periodic behaviour in terms of the continuous variable x by writing (−1) n as α = sign(cos(πx/ε)), such that where the asymptotic order of this expression can be obtained by rewriting (22) in powers of τ −1 0 and equating terms. We may now follow similar methods to the 2-periodic case, and formulate an ansatz for the solution in terms of ε and a new transseries parameter, denoted σ1. Motivated by (12), we choose the ansatz noting that we now allow an ε dependence on the exponential scale. The exponential scaling B(x, ε) may then be determined by considering the large-x behaviour. This is convenient, as we have the form for R(x) in this limit, given in (27), and we know from the asymptotic order of this expression that the solution approaches this limit exponentially as x becomes large. Using (27) in (26) and matching powers of ε in an identical fashion to (9) gives This may be solved to give where where the exponential weight B(x) is given in (30). If Re[f (ε)] > 0, the 4-periodic exponential contribution is exponentially small in ε, while if Re[f (ε)] < 0, the contribution is large, and must be incorporated into any approximation in order to accurately describe the system behaviour. In parameter regime , this exponential contribution is not present in the transseries, and is therefore denoted as a dashed curve. In regime , the 4-periodic exponential terms appear, but are exponentially small. In regime , the exponential contributions become large, and 4-periodic behaviour becomes apparent in the solution. The 4-periodicity of the solution arises due to Im[f (ε)]. This represents a multiplicative factor in B(x) of −i, corresponding to 4-periodic behaviour in the exponential term. and g(x, ε) is a bounded function that vanishes for x = εn for n ∈ Z, given by As g(nε, ε) = 0 for n ∈ Z, this expression could be ignored and the result will still give the correct value of B(x, ε), and hence the correct exponential scaling, on x = εn. This therefore suggests that we can capture the 4-periodic solution by defining a new variable τ1, such that In subsequent analysis, it will be useful to have a convenient expression for the value of τ1 at x + ε and x + 2ε. Through direct substitution, we find that

Exponential Weights
The form of B(x, ε) provides insight into the behaviour of the solution as ε increases outside of the range of validity of the original small-ε transseries. The behaviour of this term is shown in Figure 2, which identifies three distinct ranges of ε which must be considered separately. From the form of (33), we see that B(x, ε) is the exponential controlling factor for S, and therefore determines how this series will contribute as x grows. If Re[B] > 0, corresponding to Re[f (ε)] > 0, the exponential contribution will decay as x grows, while if Re[B] < 0, corresponding to Re[f (ε)] < 0, the exponential part will grow and become the most significant contribution for large x.
This change in sign occurs at ε = −2+ √ 6. At this point, Re[f (ε)] becomes negative, and the exponentials in (28) therefore grow as x becomes large, rather than decaying. This means that in Region the S term is no longer a small decaying perturbation aroundR, but rather plays a significant role in the limiting behaviour as x → ∞. If S is ignored, this behaviour is not captured in the transseries, and the resultant expression for R is an inaccurate description of the solution behaviour.
We note that Im[f (ε)] = −3π/2 for ε > −2 + √ 5. This has the effect of making τ m 1 4-periodic in m ∈ Z + , due to having a factor of −i, rather than the 2-periodic behaviour associated with a factor of −1. Hence, this exponential behaviour represented by B(x, ε) corresponds to 4-periodic effects in the solution.
In order to include this behaviour in the transseries expression, we cannot simply expand the solution about ε = 0. We must instead expand S about some point ε0 such that the 4-periodic behaviour is present in the expansion. This requirement suggests that ε0 = −2 + √ 6 is a sensible choice, as 4-periodic effects are apparent in the solution at this value.
Finally, we must consider the region in which the series terms obtained by expanding about ε0 are valid. If we examine the behaviour of B(x, ε) for ε > ε0, we see that the real part of f (ε) becomes infinite as ε → −2 + √ 5. This corresponds to the exponentials disappearing, as every exponential term tends to zero. The series expansion about ε0 is not valid for ε smaller than this value. Consequently, region contains exponentially small 4-periodic behaviour, while no such behaviour exists in region . In Figure 2 we have represented f (ε) in region as a dashed curve, to indicate that it does not have any effect on the transseries.
Consequently, simply by studying B(x, ε), we are able to describe the onset of 4-periodicity in the solution. In region , there are no 4-periodic effects present. In region , there are 4-periodic effects caused by the appearance of new exponential terms, but they are exponentially subdominant compared to the 2-periodic behaviour. In region , these effects grow to become the most significant effect in the solution behaviour. Note that the switching of the 4-period exponentials is independent of the initial data, the latter only determines how quickly they grow to dominate the solution. For higher values of ε, there must be values for which higher-periodicity behaviour appears. We discuss the onset of 8-periodic behaviour in Section 2.4.
Finally, we note that the change in exponential contribution has a parallel with a Borel transform approach to asymptotic expansions. Borel transforms encode the different exponential weights of an asymptotic series as singularities in a complex domain known as the "Borel plane". There, a change in the number of exponential contributions corresponds to singularities moving across a branch cut onto a different Riemann sheet of the Borel plane, giving rise to behaviour known as the "Stokes phenomenon" [10] as the number of exponential contributions in an asymptotic series abruptly changes. A similar, but not identical, behaviour occurs in this system at ε = −2 + √ 5, where Re[f (ε)] becomes infinite and Im[f (ε)] changes instantaneously, corresponding to a branch point in the f -plane. For more information on Borel transform methods, and their links to transseries and transasymptotic summations, see [1,13,15,43,48,49].

Computing the terms in the resummed transseries
Writing an appropriate form for the 4-periodic ansatz is slightly more involved than in the 2-periodic case, given in (15). Recall from Section 2.2.2 that the significant change in the behaviour of the exponential contribution occurs for values of ε greater than ε0 = −2 + √ 6. We therefore define a new series variable √ 6η = ε − ε0, where the √ 6 term is included for subsequent algebraic convenience. In analogous fashion to (17), we again divide the ansatz up into separate power series. In the 2-periodic case, it was clear from the form of the previously calculated terms that splitting the odd and even powers of τ0 would capture the discrete variation effectively. From the analysis in Section 2.2.2, we determine that the power series for the 4-periodic solution should instead be split into four parts, such that Consequently, we now write each split power series as functions Θ j,k , j = 1, 2, 3, 4, giving Noting that each series consists only of powers τ m 1 with the same m mod 4, and comparing this with the expression for τ1 in (33) indicates that the functions Θ j,k for j = 1, . . . , 4 must have the symmetries At this stage, it might be expected that we should express the governing equation (26) in terms of η, and perform an expansion in this variable; however, a comparison of the terms in (34) and (35) suggests that iterating the map once leads to a simplification. Writing the x dependence explicitly, the equation becomes: This expression does not contain any S(x + ε, ε) terms, and instead only contains the double iteration term, S(x + 2ε, ε). This is convenient, as the expression for τ1(x + 2ε, ε) is substantially simpler than τ1(x + ε, ε), as it does not contain g(x, ε). This simplifies significantly the subsequent analysis.
Expressing the left-hand term in (40) in terms of τ1 and η gives Rewriting (40) in terms of τ1 and η therefore gives S(−(1 + 12η + 6η 2 )τ1, η) = S(τ1, η) 1−α 2 + 12η + 6η 2 + (1 + √ 6(1 + η))S(τ1, η) Analogously to the analysis of the 2-periodic case in Appendix A, the next step is to expand this expression as a power series in η, and apply the series expression for S(τ1, η) given in (36). Matching powers of η j/2 for j = 1, . . . 4 produces a system of four equations -two of these equations are algebraic, and two are nonlinear ordinary differential equations in τ1. We omit the details of this step here, as it requires only algebraic manipulations, and the intermediate mathematical expressions are quite lengthy. These four equations may be simplified using the symmetry relations in (38)- (39), resulting in the following system of equations where By substituting the power series (36) into the governing equation (42), it can be seen at leading order as η → 0 and τ1 → 0 that S1,0 = −S3,0, providing one initial condition for the system (43)- (46). The second initial condition may be chosen arbitrarily, as this choice may be absorbed into the expression for σ1, in the same manner as the constant C in (18). For algebraic convenience, and without loss of generality, we select S1,0 = 1. These conditions are sufficient to uniquely solve (43)- (46). The solution to this system is given by In principle, we can match the expansion of (42) at higher powers of η in order to obtain Θ j,k for j = 1, . . . 4 with k > 0. For the purposes of this example, however, the first four terms of the series will produce a useful approximation for the solution behaviour. The final step is to determine the behaviour of the transseries parameter σ1. This is slightly more complicated than in the 2-periodic problem, as we must incorporate the behaviour ofR(x, ε) into the calculations. We include the details of this process in Appendix B, where we show that We have now determined enough transseries terms to accurately approximate the solution behaviour in the 4-periodic regime.

Error comparison
As a consequence of the preceding analysis, we are able to derive an approximation for the solution to the logistic equation in the 2-periodic and 4-periodic parameter regimes, which we denote as R2,app(x) and R4,app(x) respectively. Combining (16), (17), (19), (21), and (24), we find that in the 2-periodic parameter regime where τ0 and σ0 are approximated as A comparison of the exact solution against the approximation is shown for ε = 0.05 in Figure 3(a). The exact solution is shown as red circles, while the approximation is shown as blue dots. The two curves are visually indistinguishable. The approximation error is shown in Figure 3 (53), against the exact solution for ε = 0.51, or √ 6η ≈ 0.0605. The approximation errors, given by the difference between the exact solution R(x) and the approximations are shown in (c) and (d). The 2-periodic approximation has maximum error in the region just before reaching the 2-periodic steady solution. The 4-periodic approximation has maximum error in the initial region; this is to be expected, as the initial condition was obtained directly from the 2-periodic solution, and is not expected to be highly accurate in the 4-periodic regime.
that the error has a peak at the end of the transition region, just before the solution settles into the stable 2-periodic behaviour.
In the 4-periodic parameter regime ε > −2 + √ 6, the approximated transseries is given combining the expressions in (33), (48)- (50), and the previous approximation (51), to give R(x) ≈ R4,app(x) = R2,app(x) + √ η(Θ1,0(τ1)+Θ3,0(τ1)) + η(Θ2,0(τ1) + Θ4,0(τ1)), where τ1 and σ1 are approximated as Note that we do not include the term containing g(x, ε) in B(x, ε) from (30). This term disappears for integer values of n, and therefore can be omitted at this stage without altering the approximation. In Figure 4(a), we show the approximation error for a range of values of ε, where the error is measured as the maximum difference between the exact solution and the transseries approximation, shown as a blue curve. This error measure was chosen to allow for direct comparison with equivalent results from [40], which are shown as a red curve. The transseries approximation is more accurate than the multiple scales approximation in this parameter regime, and the error decays faster in the limit that ε → 0. The reason for this behaviour is that the transseries approach allowed for higher-order exponential corrections to be easily computed and retained. The maximum approximation error occurs at the end of the transition region between non-periodic and 2-periodic behaviour, where the exponential contributions contribute significantly to the solution behaviour. Computing these exponential corrections using multiple scales methods would be an algebraically significantly more demanding task, requiring matched asymptotic expansions to be applied at higher orders of the expansion.

8-periodic solution
We may continue this process to understand the emergence of the next period doubling bifurcation. While we will not include a full explicit, algebraic analysis here, we will show that the exponential factor can be used to identify the appearance of 8-periodic stable solutions as ε is increased further.  The plot in (a) shows the resummed transseries approximation error in blue, corresponding to the maximum difference between the approximated and exact value. This measure of the error was chosen to be consistent with the error measure provided in [40]; this error is shown as a red curve. Due to the ease with which the transseries method captures higher-order exponential behaviour, which plays an important role in the transition region between non-periodic and 2-periodic behaviour, it outperforms the multiple scales approximation. In (b), we show the error of the four solution branches as n → ∞, approximated by taking |R − R 4,app | on each of the four branches for a value of n sufficiently large that the error is not visibly changing. For each of the four branches, the error decreases as η → 0, as would be expected.
The method from Section 2.2 can be applied again in order to obtain approximations for solutions with even higher periodicity. We can now write the next term in the transseries such that R(x, ε) = R(x, ε) + S(x, ε) + T (x, ε). The quantity T (x, ε) is defined in terms of a new transseries parameter σ2 to be The transseries termsR + S capture the 4-periodic solution behaviour, and therefore must tend to the 4-periodic solution in the limit that τ0 and τ1 become large. We denote this solution as R4(ε). Hence, we apply the expression R(x, ε) = R4(ε) + T (x, ε) to the governing equation (5) and find an expression for the exponential weights, in similar fashion to the process for obtaining (10) or (29). The exponential weights may again be written in the form F (x, ε) = f (ε)x + εg(x, ε), where g disappears on n ∈ Z. The behaviour of f (ε) is illustrated in Figure 5. A very similar set of inferences may be drawn from this image as for Figure 2. In region , the 8-periodic behaviour does not contribute to the solution, as discussed for the 4-periodic case in Section 2.2.2. This 8-periodic contribution appears in the transseries as ε moves into region . In this range of ε, there are 8-periodic contributions to the solution, but they are smaller than the 4-periodic solution contribution, as the exponential term is relatively small compared to those in S(x, ε), decaying exponentially as x → ∞. Finally, in region , the 8-periodic solution grows exponentially, and the behaviour of T (x, ε) dominates the solution behaviour.
It is therefore clear that we can explain the onset of these higher periodicity solutions by explicitly studying the exponential weights of the transseries solution; while the algebraic complexity of the process increases after each doubling, the steps for identifying this behaviour remain essentially the same. The resummed transseries therefore provides a systematic approach to studying bifurcations even for larger values of the bifurcation parameter, where classical asymptotic methods typically fail.

Dynamical Logistic Equation
In the previous section, we studied the classical logistic equation, and showed that the higher periodicity solutions may be obtained directly using a transseries approach. In this section, we consider a more complicated variant of this problem, known as the slowly-varying logistic equation. with ε > 0. The bifurcation parameter is given by λ = λ0 + εn, and in this case, it changes slowly over time. In previous studies [31,33,29], this has been shown as an example of a "canard" solution, in which the behaviour appears to remain near the unstable solution for an extended period of time, before rapidly jumping to approach the stable solution with higher periodicity. As n increases, this parameter will pass through values across which the solution stability is known to change. When λ0 + εn = 3, the 1-periodic equilibrium becomes unstable, and the 2-periodic equilibrium becomes stable. As n increases further, eventually λ exceeds 1 + √ 6, and the 2-periodic equilibrium becomes unstable, with the 4-periodic equilibrium becoming stable. This process continues until the bifurcation parameter becomes sufficiently large that the solution becomes chaotic. For the problem studied here, we will set λ0 = 3 and y(0) = 2/3.
In [40], it was shown that a discrete multiple scales approach can be used to describe this behaviour asymptotically. This approach required balancing several different timescales, and using asymptotic matching to connect the solutions in each different asymptotic region.
In this section, we will show that this process can be described using a transseries approach, with the resulting expansion to be valid even as the solution behaviour changes dramatically, and increases in periodicity. We will now show that transseries provide a systematic and generally more accurate approach than the multiple scales procedure of [40] in describing the solution behaviour as it transitions from an unstable to a stable manifold; this will demonstrate that transseries expansions can be used to effectively capture canard behaviour in discrete systems.
We will show the first stability transition in detail. We will subsequently provide an outline of how this method can be extended to describe the second transition, together with some results; however, the algebraic manipulations for this process are quite involved, and the precise details will be omitted.

Transseries ansatz
The difference from (5) above is that in the prefactor of the r.h.s. the perturbative parameter ε is replaced now by x. We again begin by applying a multiple scales ansatz, and expanding as a transseries in a continuous variable x. Setting x = εn and R(x) = y(n) gives We again formulate an ansatz for the solution in terms of ε and a transseries parameter σ0. The ansatz is identical to that given in (6), but has been included below for convenience: where βm will again be chosen such that Rm,0 takes nonzero value. It is straightforward to compute the first few terms of the algebraic portion of the series expression, corresponding to m = 0 in (59), which gives a power series expression for the non-periodic manifold. The recursion relation is given obtained by expanding R(x + ε) using a power series in ε, and matching powers of ε in the resultant expression. This process gives The first few iterations of this recurrence relation give (62) This process may be continued indefinitely in order to continue calculating terms in the power series for the non-periodic manifold. This process will not, however, capture the transition to the 2-periodic manifold. In order to obtain an approximation for this behaviour, we are required to consider terms in the ansatz (59) with m = 0. Continuing to the next order in σ0, we find that As before, the argument of the exponential may be determined by expanding R1 as a power series in ε, as well as expanding A(x + ε) = A(x) + εA (x) + · · · . At leading order in ε, this gives the differential equation Hence, we obtain where we follow the same reasoning as the analysis used to determine (10), and absorb the constant into the series parameter. We may again set p = 0; this choice will have no effect on the behaviour of the solution for integer values of n.
It is now possible to apply the second part of the ansatz in (59) and to match powers of ε in this expression. By direct substitution, we find that βm = m gives the result that Rm,0 is nonzero. By subsequently matching terms which are O(ε) in the small ε limit, it is possible to generate an equation for R1,0 and a recurrence relation for Rm,0 for m ≥ 2. We find that The initial condition in (67) may be chosen arbitrarily, as this choice may be absorbed into the transseries parameter. Choosing R1,0(0) = 1 gives .
The recurrence relation for subsequent terms is given by Rn,0(x)Rm−n,0(x), m ≥ 2.
Continuing to match higher powers of ε allows for the direct computation of terms further terms such as R m,k , obtained by matching terms which are O(ε k ) in the small ε limit. The direct computation of further terms is not required for the present analysis.

Computing terms in the resummed transseries
Motivated by the analysis of the static system, and in particular the form of (15), we switch the order of summation in the transseries (59), writing it as A main difference from the static one is that the expansion in powers of eps is asymptotic, while the sum over the exponentials (m ≥ 0) is convergent. Thus (70) is a formal expansion. As for the static system, we define a new series parameter τ0, and new quantities Ω k (τ0) such that It will be helpful later to note that The transseries expression in (70) is now given by We can now apply this expression to (58) and match orders of ε. At leading order, we find that where (72) was used to obtain the leading-order on the left-hand side. At this stage, we could mechanically obtain the function Ω0 as a Taylor series in τ0, which is convergent, with some finite radius of convergence. It happens, however, that there exists a particularly convenient variable transformation that converts the right-hand side from a dilation to a translation. If we define a new variable y such that y = −x log(τ0)/A (x), the expression in (74) becomes Ω0(y + x) = (3 + x)Ω0(y)(1 − Ω0(y)).
This expression has the same form as the static logistic map equation, given in (5), with x in place of ε. Furthermore, since x = εn, it is valid to apply the asymptotic solution derived for this expression in Section 2.1. As we are interested in capturing the first transition, across which the solution switches from having no periodic component to having a 2-periodic component, we can directly apply the transseries expression for the 2-periodic solution given in (17). In order to take into account the form of (75), we must replace ε and x with x and y respectively in (17). We must also replace the τ0 in this expression with a new transseries parameter τ 0, in which ε and x are again replaced with x and y respectively. This gives where σ0 is a new transseries parameter that remains to be determined. Making the appropriate substitutions in (17) now gives the form of Ω0(y) as where Ω o,k and Ω e,k are defined in (16), and Ω o,k and Ω e,k for k = 0 and k = 1 are given explicitly in (19) and (21) respectively. In a typical problem of this form, σ0 would be determined using the fact that Ω0 = 2/3 at y = 0; however, this is enforced by the transformation y = −x log(τ0)/A (x), which forces x to be zero if y = 0. Consequently, the initial condition cannot be used to determine σ0. This is to be expected, as the the initial condition will instead be used to determine the original transseries parameter σ0.
Instead, we expand (77) as a Taylor series about x = 0 using the form of Ωo,0 and Ωe,0 given in (19). This gives where the omitted terms are O(x 3 τ 2 0 ). We may now match powers of τ0 with (70) to determine that xσ0 = R1,0, which was explicitly calculated in (68). We therefore find that .   (82) and exact solution of (57) for ε = 0.001. The difference between these is shown in (b). The points labeled , and will be referenced below in Figure 7. We see that the error reaches a maximum at the start of the 2-periodic regime. It then decreases, although will eventually increase as n grows, due to the increasing influence of the 4-periodic region which was not computed. Note that the small error at the point labelled corresponds to the value where the approximation crosses the exact solution. This occurs at some point in the 2-periodic region for any choice of ε, and therefore does not signify a special parameter choice. It is an artifact of the error calculation.

Initial Condition
We have now explicitly calculated all of the required quantities for the transseries approximation except for σ0, which must be determined from the initial condition at x = 0. At x = 0, it follows that τ0 = σ0ε. Consequently, the initial condition is given by R(τ0 = σ0ε, x = 0) = 2/3, which we apply to the first expression in (59). We then express σ0 as a power series in ε, where the series terms are functions of R m,k for various values of m and k, giving Using the second expression from (59) and matching powers of ε allows us to compute σ0,j. We have obtained enough R m,k terms to solve for σ0,0, giving Computing subsequent series terms for σ0 requires values of R m,k that are not presented in this study, as even this first order approximation is sufficiently accurate as we now show.
The most useful feature of this approximation is that it is valid before, during, and after the transition region in the slowly-varying logistic equation. We illustrate an example comparison in Figure 6(a), corresponding to ε = 0.001. The approximation is shown as blue dots, and overlaid on top of the exact solution, shown as red circles. The two solutions are visually almost indistinguishable.   [40] as belonging to the inner region, transition region, and outer region, shown as points , and in Figure 6. In each case, the error is shown as a blue curve. This curve becomes smaller as ε decreases. The point at which the error dips is an artifact of the observation that the approximation crosses the exact solution at the identified point for this choice of ε, and does not represent any significant phenomenon within the transseries approximation. The cause of this behaviour is explained in more detail within the description of Figure 6. We have chosen a similar range of small parameter to the analysis in [40], shown in red. The transseries outperforms the multiple scales method in both the transition region and the outer region, in which the exponential terms play an important role in describing the solution behaviour. These terms are more naturally captured using transseries methods, leading to an improved approximation.
The approximation error for this example is shown in Figure 6(b), calculated by y(n) − Rapp(εn). The error reaches a peak following the transition region, at the beginning of the stable 2-periodic behaviour. The error does grow in this region as n becomes large, and continues to do so until the transition to 4-periodic behaviour occurs. This behaviour is not depicted in Figure 6(b).
In order to obtain a more complete picture of the accuracy of the transseries approximation, we determined the approximation error at three selected values of n. These values were tested in [40] relative to other methods, to obtain representative computations of the approximation error in important parts of the solution domain. The first point is n = 1/ √ ε . This point is found in the early non-periodic region before the transition from non-periodic behaviour to 2-periodic behaviour occurs. It is labelled in the example solution from Figure 6.
For comparison, we need to identify the remaining representative points used in [40], which required the computation of an intermediate quantity K, satisfying This quantity was derived in [40] although it has been adjusted to take into account the slightly different form for the slowly-varying logistic equation considered here. The second point falls within the transition region between the non-periodic unstable manifold and the 2-periodic stable manifold, and is given by This point is labelled in the example solution from Figure 6. Finally, we also determine the error at a point in the region where the solution has completed its transition to 2-periodic stable behaviour. This point is given by n = K + 15K −1 / √ ε , and is labelled in the example solution from Figure 6. The error for each of these three points was studied in [40] allowing for direct comparison between the transseries approximation and the multiple scales approximation errors.
The error for each of the three representative points over a range of ε values may be seen in Figure 7, shown in blue. The error for the approximation from the multiple scales approximation in [40] is shown in red for each point. In each region, both approximations are relatively accurate. In the non-periodic region, the multiple scales approximation outperforms the transseries approximation, while in the transition and 2-periodic region, the transseries approximation is substantially more accurate. This outcome is sensible; the transseries approximation tracks the contribution of exponentials in the solution, and accurately incorporates them into the solution behaviour. In the non-periodic region, the solution is best represented by an algebraic power series in ε. The multiple scales approach involves calculating this power series to several terms, while our transseries approximation relies only on the leading-order . For z < √ 5−2, the weight B(z) contains an imaginary term −πi, which corresponds to 2-periodic behaviour. After z exceeds √ 5 − 2, the slope of the imaginary term changes to −3πi/2, which leads to the appearance of 4-periodic behaviour. This behaviour is not immediately apparent, as the contribution is exponentially small if Re[z] > 0, corresponding to z < z 0 , where z 0 ≈ 0.9951. For z > z 0 , the 4-periodic terms become significant in the solution behaviour. We note that, due to the bifurcation delay, this behaviour is not immediately visibly apparent in the solution; however, a careful analysis of the corresponding transseries terms will identify the transition between 2-periodic and 4-periodic behaviour.
behaviour of this series. In the transition and 2-periodic region, however, these exponential contributions become more significant, and this corresponds to the transseries approximation becoming more accurate than the multiple scales approximation. While the multiple scales approximation is able to capture some of the exponential behaviour, the transseries approximation is able to incorporate several exponential corrections in a straightforward fashion, producing greater accuracy in the solution regions where these corrections play an important role. Furthermore, increasing the accuracy of the transseries approximation in the non-periodic region can be done systematically by including higher corrections in ε.
Finally, we note that there are points in Figure 7(b)-(c) where the error appears to drop to zero. This corresponds to a coincidental crossing between the approximation and actual solution occurring at this value of n. The crossing may be seen in the example solution from 6 at n ≈ 250. Any solution with a reasonable amount of accuracy will have some value of n where this crossing occurs; this does not provide any added insight into the accuracy of the approximation.

4-periodic solution
From Figure 1, we see that as n increases, eventually the bifurcation parameter becomes sufficiently large that the solution becomes 4-periodic. As discussed in Section 2.2.2, this higher periodicity must be encoded in the transseries solution as the weights of new exponential scales. We will not perform a full explicit calculation here, but we will demonstrate that these exponential weights do, in fact, predict the emergence of stable behaviour with higher periodicity.
We include a new term with transseries parameter σ1, and the analysis suggests that it is natural to define a new scaled variable z = 2nε. This new contribution to the transseries, denoted S(z, ε; σ1), is given by By adding S as a perturbation to the 2-period solution approximated by (82) and balancing terms in (58) in a similar fashion to Section 3.1.1, we obtain an equation for B(z) that gives where the constant of integration is picked to set B(0) = 0 for convenience, though this choice may be absorbed into the parameter σ1. The behaviour of B(z) is depicted in Figure 8. There are two significant conclusions that may be drawn from this figure. In Figure 8(b), we see that Noting the format of (85), we see that this exponent changes from 2-periodic behaviour to 4-periodic behaviour when crossing the value z = √ 5 − 2. This change in exponent gives rise to 4 different branches in the solution, and therefore explains the onset of 4-periodic behaviour in the solution to the dynamical logistic equation.
The second important observation is that this 4-periodicity is not immediately apparent in the solution, due to the behaviour of Re[z]. In Figure 8(a), we see that there is a value of z, denoted z0 and located at z0 ≈ 0.9951, at which the real part of B(z) changes sign from positive to negative. From the form of (85), we see that this corresponds to the 4-periodic transseries contribution being exponentially small for z < z0, before growing to have a significant impact on the solution behaviour for z > z0.
This value of z0 corresponds to 4-periodic behaviour becoming apparent at n ≈ 0.4975/ε. For example, in Figure 1, we would expect that 4-periodic solution to become significant at n ≈ 3455, which is consistent with the appearance of the second transition region in this image.
A more detailed transseries analysis would permit us to calculate a series approximation for the 4periodic behaviour; however, as we expected from the transseries approach, a straightforward analysis of the exponential weights in the transseries is sufficient to explain the onset of the higher periodicity, and identify the location in z (and hence, in n) where this transition to dominant 4-periodic behaviour takes place.
Finally, we note that the points where the periodicity changes correspond to values of n where the real part of the exponential weights changes sign, or z0 in Figure 8. In asymptotic analysis, this corresponds to the crossing of a curve known as an anti-Stokes curve. This suggests that the Stokes phenomenon plays a role in this system behaviour. In fact, the solution does contain Stokes curves, which are responsible for appearance of exponential factors in the solution; however, finding these Stokes curves requires continuing the solution in the negative-n direction, and was therefore not presented here. Nonetheless, the study the Stokes phenomenon in the dynamic logistic map is an interesting and rich subject which is beyond the scope of the present work.

Discussion
We have obtained transasymptotic approximations for the solutions to both the standard and slowly varying logistic equation. In each case, we were able not only to reproduce the results calculated in [40] using multiple scales asymptotic methods, but to go significantly further.
As, a priori, transseries methods allow for the straightforward calculation of higher-order exponentials, the transseries approximation was able to represent the solution more accurately than the multiple scales method both during and after the delayed bifurcation, as seen in Figure 7; during and after the birfurcation, the initially subdominant exponentials contribute significantly to the solution, so it should be expected that the transseries approximation would be particularly accurate compared to other methods in these regions.
Furthermore, the transseries approach can still provide a useful approximation when the parameter ε is not particularly small, as the solution can simply be rescaled to determine the next asymptotic weight.
We considered the dynamic logistic equation with ε > 0, producing a cascade of delayed bifurcations. If ε < 0, causing the bifurcation parameter to decrease rather than increase, bifurcations appear earlier than the solution stability would suggest, rather than later [7]. A transseries approach could be used in almost identical fashion to the present study in order to approximate these accelerated bifurcations; however, that analysis is beyond the scope of this study.
There are several significant and general advantages to the transseries resummation approach. The first is that the method we have described can be applied in systematic fashion to a wide range of problems, including both discrete and continuous systems. Whilst such advantages have already been seen elsewhere, in the context of the logistic equation it has been instructive to compare this to the multiple scales approach from [40], which required the careful comparison of asymptotic terms up to several orders. In order to capture the fast discrete scale, as well as both the inner and outer continuous scales near the bifurcation, asymptotic matching was performed through three scales. The transseries approach here was able to reproduce this behaviour by resumming the series in order to ensure that the transasymptotic approximation contained behaviour encoded in the subdominant exponentials; this behaviour contained all of the information found using asymptotic matching methods. Furthermore, improving approximation accuracy by computing more terms of a multiple scales expansion requires comparing the asymptotic behaviour of more terms and checking to determine when the relative dominance of terms changes, while obtaining a more accurate transseries expression simply requires the systematic calculation of further series terms in the transseries. While these calculations can prove challenging, the steps required to obtain the subdominant exponentials, and the associated solution behaviour, follow the same consistent process at each stage which it is applied.
Computing the subdominant exponentials is not valuable simply in that it produces a more accurate approximation. In fact, a second major advantage of the transseries method is that the exponential weights have a significant effect on the system behaviour, and computing just these weights can tell us the form of the solution as parameters in the problem vary. In our analysis of the standard logistic map, we showed that 2-, 4-, and 8-periodic behaviour can be determined simply by carefully studying the subdominant weights. This explains the appearance of higher periodicities in the solution, and suggests that if this process is continued, it can be used to study further bifurcations in the period doubling process.
In our subsequent analysis of the period doubling cascade found in the slowly varying logistic equation, we were able to predict the onset of 2-periodic and 4-periodic behaviour in the solution, simply by studying the relative size of the exponential weights associated with the 2-and 4-periodic contributions. It would be particularly interesting to continue to investigate how the full period doubling cascade is encoded in the exponential weights of this system, and whether this can provide (at least theoretically) further insight into the period doubling route to chaos.
Finally the full analysis of the movements of exponential contributions between Riemann sheets, seen in the dynamic logistic map, also merits further investigation. Examples of such phenomena have been observed recently in novel features of aeroacoustic flows [53,54]. Initial explorations appear to suggest this is commonly found in other physical and mathematical contexts.