The number of limit cycles for regularized piecewise polynomial systems is unbounded

In this paper, we extend the slow divergence-integral from slow-fast systems, due to De Maesschalck, Dumortier and Roussarie, to smooth systems that limit onto piecewise smooth ones as $\epsilon\rightarrow 0$. In slow-fast systems, the slow divergence-integral is an integral of the divergence along a canard cycle with respect to the slow time and it has proven very useful in obtaining good lower and upper bounds of limit cycles in planar polynomial systems. In this paper, our slow divergence-integral is based upon integration along a generalized canard cycle for a piecewise smooth two-fold bifurcation (of type visible-invisible called $VI_3$). We use this framework to show that the number of limit cycles in regularized piecewise smooth polynomial systems is unbounded.


Introduction
In this paper, we consider smooth systems of the forṁ z = Z(z, φ(h(z) −1 )), (1.1) for z ∈ R n , 0 < 1 and where h : R n → R is regular, φ is a regularization function: and where Z is affine in its second component: These systems have recently received a great deal of attention [8,9,36,35,39,40,41,42,43]. The motivation is three-fold. Firstly, in the limit → 0 the system (1.1) becomes piecewise smooth (PWS)ż = Z + (z) for h(z) > 0, Z − (z) for h(z) < 0, (1.4) with Σ := {z : h(z) = 0} being a discontinuity/switching manifold, see Fig. 1. Figure 1. A planar piecewise smooth system (1.4), having Σ as a switching manifold. Here regular orbits of Z + and Z − reach Σ in finite time. The point Σ T is a tangency point of Z + (a visible fold if the tangency is quadratic). In the present case, where Z − is transverse to Σ it divides Σ into sliding (to the left of Σ T ) where Z ± are in opposition relative to Σ and crossing (to the right of Σ T ) where Z ± point in the same direction relative to Σ. The situation is different if Z − also has a fold at Σ T (called a two-fold). Then there can be sliding (stable and unstable) on each side of Σ T , see Fig. 2.
For 0 < 1, the system (1.1) is therefore a regularized PWS system [18,28]. The reason for restricting to (1.3) is that in this case, one can show [55] that the singular limit system is a Filippov system [18,28]. Lately, there has been a growing interest in understanding how PWS phenomena (folds, grazing, boundary equilibria,... [45]) unfold in the smooth version [35,36,40,41,42]. For this purpose methods from Geometric Singular Perturbation Theory (GSPT) and blowup have been refined to deal with resolving the special singular limit of (1.1) [40,41]. Finally, the interest in systems of the form (1.1) is also motivated by applications. For example, in biology switches [18,43] are frequently modeled by functions (1.2) and friction is also inherently piecewise smooth [3].
Mathematically, piecewise smooth system has also received a great deal of attention over the past decades. Starting from the groundbreaking work of Filippov [28] and Utkin [57], there was an effort to extend Peixito's program of structural stability to PWS systems [10,55]. Subsequently, there has been a focus on characterizing and interpreting the lack of uniqueness of solutions in PWS systems [34].
Parallel to this effort, there has been an attempt to bound the number limit cycles in PWS systems in the plane where n = 2. In contrast to the smooth linear setting, limit cycles can exist for piecewise linear systems and J. Llibre and co-workers have obtained upper bounds for a number of cases [26,49,52]. Of course, the interest in bounding the number of limit cycles, comes from Hilbert's 16th problem [48] which seeks to bound the number of limit cycles of polynomial systems:ẋ = P N (x, y), y = Q N (x, y), (1.5) with P N and Q N of fixed degree N . Hilbert's 16th problem remains unsolved to this day. Whereas general progress has been made on N = 2 [1,4,5,20,24,25,31,53] and on Smale's version of the problem where (1.5) is restricted to classical Liénard type: P N (x, y) = y −p N (x), Q N = −x, see [11,46,50,54], there has been an emphasis on obtaining good lower bounds on the number of limit cycles (see [12,29,30,51] and references therein). Following the work of De Maesschalck, Dumortier and Roussarie, see [14,15,19,22,23], a key tool in this effort has been the slow divergence-integral from slow-fast systems and canard theory; in particular, the roots of the slow divergence-integral provide candidates for limit cycles. For example, using this tool good lower bounds on the number of limit cycles in Liénard equations can be found (see [16,17,21,33,59]).
1.1. Main result. In this paper, we work at the interface of these research fields. In particular, we consider (1.1) with n = 2, put z = (x, y) and restrict attention to the case h(z) = y so that the switching manifold is Σ = {(x, y) : y = 0} and then ask the following question: Does there exist polynomial vector-fields Z ± such that the number of limit cycles of Z is unbounded?
We prove that this is in fact true, even for quadratic vector-fields. More precisely we prove the following. Theorem 1.1. There exists a quadratic vector-field Z + (·, λ) and a linear vector-field Z − (·, λ), depending smoothly on a parameter λ ∈ R, such that the following holds true in a compact domain U : For every k ∈ N there exist: (a) k > 0, (b) a regularization function φ k : R → R, and (c) a continuous function λ k c : [0, k [→ R such that the regularized vector-field: ), has at least k limit cycles contained in U for all ∈]0, k [.
We give examples of Z + and Z − later on, see (4.13)(4.14). We emphasize that the unboundedness of limit cycles stems from the regularization and not from the vector-fields Z ± . We use smooth regularization functions in order to find an unbounded number of limit cycles. It is known that boundedness of limit cycles is closely related to the notion of o-minimality in function spaces (see e.g. [37]). Our smoothings are taken from a family that does not have this o-minimality property. From this viewpoint it is not surprising that we find that the number of limit cycles is unbounded.
At the same time, Theorem 1.1 also illustrates a certain degree of deficiency with smoothing piecewise smooth systems (since the result may depend upon how we regularize). On the other hand, there are other complementary results, see [9,35,36,40], that show that smoothing play little role (at least on a macroscopic-level, i.e. at O(1)) for different types of PWS singularities and bifurcations. In [40] for example, it was shown that the regularization of the visible-invisible fold in R 3 , with Σ being two-dimensional, is independent of the smoothing function. In fact, for the system in Theorem 1.1 it is also only in an exponentially small parameter regime that a different number of limit cycles can be realized for different regularization functions.
To prove Theorem 1.1, we will follow the approach of [14] and use a slow divergenceintegral. But seeing that our system is nonsmooth (as opposed to slow-fast) in the singular limit → 0 we will first have to develop this framework within the setting of (1.1). For slowfast systems, the slow divergence-integral is defined along a canard trajectory, i.e. along a folded critical manifold with an equilibrium at the fold in such way that the reduced problem goes from the attracting sheet to the repelling one with nonzero speed. In the setting of (1.1), our slow divergence-integral will be based upon the PWS two-fold bifurcation [7,42], which is reminiscent of the standard canard [22,44]. In particular, Z ± in Theorem 1.1 will be chosen so that the PWS system has a two-fold bifurcation (of type visible-invisible called V I 3 [45]). Proposition 3.2 then describes the structure of the difference map near the associated canard-like limit periodic sets (see Section 3.2). Proposition 3.2 is not only relevant and important for proving Theorem 1.1, but also for studying bifurcations of limit cycles inside such visible-invisible two-folds (see Remark 3.4). This proposition is therefore also one of our main results, but we delay the detailed statement to later sections after having introduced the two-fold bifurcation model (see Section 2).
Our approach for constructing an unbounded number of limit cycles, does not work for the piecewise linear case. It would be interesting to study the linear case more carefully in future work.

1.2.
Overview. The paper is organized as follows: In Section 2, we define a planar PWS two-fold and revisit some results from [7,42] on canards of (1.1) for 0 < 1. Next in Section 3 we define the slow divergence-integral and prove that simple roots of this function lead to hyperbolic limit cycles (Theorem 3.1). In the proof of Theorem 3.1, we describe the difference map in terms of the slow divergence integral in Proposition 3.2. For the proof of this statement, we also use Appendix A and Appendix B. In Section 4 we then prove Theorem 1.1, using Theorem 3.1, see also Theorem 4.3, and finally in Section 5 we illustrate our approach with numerical examples.

The two-fold bifurcation
We consider (1.1) with h(z) = y:ż for z = (x, y) ∈ R 2 . In comparison with (1.1) we have also included λ ∼ λ 0 ∈ R as an additional unfolding parameter. Notice also that we write −2 in (2.1) rather than just −1 , since this will be convenient later on (see Section 3.1). The basic assumption is that the right hand side Z is smooth in each entry (in this paper, by "smooth" we mean differentiable of class C ∞ ). In particular we suppose that it is affine in the second component, i.e., where Z ± = (X ± , Y ± ) are smooth in (z, λ). The function φ : R → R is a smooth sigmoidal function satisfying the following assumptions: (A1) The function φ has the following asymptotics when s → ±∞: (A2) The function φ is strictly monotone, i.e., φ (s) > 0 for all s ∈ R. By assumption (A1), the system (2.1) is piecewise smooth (PWS) in the limit → 0: the set Σ defined by (x, 0) being the discontinuity set/switching manifold, for each λ ∼ λ 0 . In fact, from assumption (A3) we have that (2.1) is a regular perturbation of Z + or Z − outside any fixed neighborhood of y = 0. In particular: smoothly and uniformly with respect to → 0.
By plugging the expression for 3), we may write X sl as The vector-field (2.5) is known as the Filippov sliding vector-field [28] and PWS systems with this vector-field prescribed on Σ sl are called Filippov systems.
2.1. Folds. Clearly, Σ = Σ sl (λ) ∪ Σ T (λ) ∪ Σ cr (λ) for each λ. We further classify the points in Σ T as follows (see also [18]): (4) A point q ∈ Σ T (λ) is a fold point from "above" if the orbit of Z + (·, λ) through q has a quadratic tangency with Σ at q. In terms of Lie-derivatives Z ± (h)(·, λ) := ∇h · Z ± (·, λ), with h(x, y) = y, the last condition becomes: We define a fold point from "below" in terms of Z − in a similar way. (5) A fold point q ∈ Σ T (λ) from "above" is said to be visible, if the orbit of Z + (·, λ) through q is contained within y > 0 in neighborhood of q. It is said to be invisible otherwise. In terms of Lie-derivatives, we clearly have Z 2 + (h)(q, λ) > 0 iff q satisfying Z + (q, λ) = 0, Z + (h)(q, λ) = 0 is visible. Fold points from below are classified in a similar way. In particular, Fold points that are only PWS singularities on one side of Σ are persistent by the implicit function theorem, in the following sense: If Σ T (λ 0 ) consists of a fold point q(λ 0 ) (from above or below), then Σ T (λ) also consists of a fold point q(λ) (from above or below, respectively) for any λ ∼ λ 0 . In fact, q(λ) then also depends smoothly on λ ∼ λ 0 .
2.2. Two-folds. Now, we finally arrive at the concept of two-folds in PWS systems, which will play the role of a canard point in our analysis of (2.1).
(6) A two-fold q ∈ Σ T (λ) is a point with quadratic tangencies from above and from below. In terms of Lie-derivatives we have: with these equations understood to hold for both ±. (7) A two-fold is said to be visible-visible, visible-invisible, invisible-invisible according to the "visibility" of the fold from above and below, respectively, see item (5) above. The three distinct cases are illustrated in Fig. 2. The further details depend on the direction of the flow. In fact, according to [45] there are 7 cases, two visible-visible cases (called V V 1,2 ), three visible-invisible cases (called V I 1−3 ) and two invisible-invisible cases II 1,2 ). We refer to [45] as well as [7,42] for further details here. They will not be needed in the present manuscript.
In contrast to a fold, a two-fold is a co-dimension one (PWS) bifurcation [7]. Consequently, Here we denote by () the partial derivative with respect to x, a convention we will continue to adopt in the following. invisible-invisible (c). We have deliberately not put arrows on the orbits of Z − (red) and Z + (green), because Σ sl and Σ cr depend on this direction. Notice Σ (orange) is the x-axis in all figures. Following [45] there are 7 cases, two visible-visible cases (called V V 1,2 ), three visible-invisible cases (called V I 1−3 ) and two invisible-invisible cases II 1,2 ). The case V I 3 , which will be our main focus, is illustrated separately in Fig. 3.
In the present paper, we will focus on the visible-invisible two-fold. In this case, [7, Lemma 2.8] shows that if q is a visible-invisible two-fold for λ = λ 0 , then locally Consequently, X sl (x, λ 0 ) is in this case locally defined for all points on Σ(λ 0 ) except q (see Theorem 2.2). Notice in particular from the form (2.5) that X sl (x, λ 0 ) has a "0/0" at the two-fold. However, by (2.6) and (2.8) we also have that at (q, λ 0 ), and consequently from (2.5) we see that X sl (x, λ 0 ) can be extended locally to all of Σ by L'Hospital in this case. We collect the findings in the following proposition (fixing q = 0 for simplicity).
Proposition 2.3. Consider a PWS system (2.2) in a sufficiently small neighborhood of the origin. Suppose furthermore that (2.10) Then the following holds about system (2.2) for λ = λ 0 : Henceforth, we suppose that (2.10) holds and that X sl (x, λ 0 ) > 0 for all x ∈ Σ, so that the flow of X sl takes points from stable sliding to unstable sliding. These conditions -which following Proposition 2.3 item (iv) and (2.5) imply that det Z > 0 (2.11) -correspond to the specific visible-invisible two-fold called V I 3 in [45]. See an illustration of this case in Fig. 3. We collect these assumptions and (2.7) into the following hypothesis.
(A4) Suppose that (2.7) and (2.10) both hold and that there are µ − < 0 and µ + > 0 such that the PWS system (2.2) for λ = λ 0 has stable sliding for all x ∈ [µ − , 0[ and unstable sliding for is the x-value of the first intersection with the x-axis of the forward flow of (x, 0) following Z − for λ = λ 0 .
2.3. Canards of (2.1). In [7,42], it was independently shown that under the assumption (A4), the two invariant manifolds for x < 0 and x > 0 (which are slow manifolds within the scaling regime defined by y = 2 y 2 , recall the proof of Theorem 2.2) intersect along some λ ∼ λ 0 for all 0 < 1. Such orbits are also called canards. The reference [42] used the blowup method, which will also form the basis of our analysis.
for all x locally so that the stable sliding region x < 0 is connected to the unstable sliding region by the flow of X sl (extended through x = 0). ξ and ξ ± are used in relation to the slow divergence-integral in Section 3. .

The slow divergence-integral and canard limit cycles
Consider (2.1) satisfing (A1)-(A4). For λ = λ 0 the singular limit (Filippov) system is shown in Fig. 3. The situation is clearly reminiscent of the classical canard situation. In particular, at the level λ = λ 0 , we denote by Γ x for x ∈]0, µ + ], the limit periodic set consisting of the segment [ξ(x), x] ⊂ Σ and the regular orbit of Z − connecting (x, 0) and (ξ(x), 0). We call Γ x a canard cycle. We then define the associated slow divergence-integral along the segment [ξ(x), x]: The slow divergence-integral is the integral of the divergence of the vector field (2.4), for = 0, computed along the critical manifold w.r.t. the slow time τ defined by dτ = dx X sl (x,λ 0 ) . It follows from (A4) that I in (3.1) is well-defined. The following result plays a crucial role in proving Theorem 1.1.
The periodic orbit O i is isolated, hyperbolic and Hausdorff close to the canard cycle Γ x i , for each i = 1, . . . k.
A result similar to Theorem 3.1 for smooth planar slow-fast systems can be found in [17,19].
Notice that the statement of Theorem 3.1 deals only with limit cycles of size O(1) in the (x, y)-phase space. Once the positive simple zeros of the slow divergence integral I are detected, the related canard cycles Γ x 1 , . . . (and hence the limit cycles O 1 , . . . ) are of size O(1). Thus, the limit cycles born from the origin (x, y) = (0, 0) are not covered by Theorem 3.1.
We divide the proof of Theorem 3.1 into three parts. In the first part we consider the extended fast-time system (z , ) = ( 2 Z, 0) and then gain smoothness by applying a cylindrical blow-up (see Section 3.1). Using the cylindrical blow-up we replace the discontinuity line Σ of the PWS system (2.2) with a half-cylinder and we show that near the canard trajectories on the top of the cylinder, we are in the framework of [14]. In [14] a very general smooth planar slow-fast model has been studied containing a normally attracting branch of singularities, a normally repelling branch of singularities and a turning point between them (an additional critical curve passing through the turning point is possible). One usually uses the results of [14] for specific slow-fast families by checking the assumptions in [14] (see for example [25,32,47]). We do the same here. In the second part (Section 3.2) we find the structure of the difference map of (2.1) near Γ x using [14] and Proposition B.1 (Appendix B) near the hyperbolic edge of the cylinder. In the third part (Section 3.3) we establish a one-to-one correspondence between simple zeros of the slow divergence-integral (3.1) and simple zeros of the difference map by choosing a suitable control function λ = λ c ( ), following [2,14].
3.1. Cylindrical blow-up. First we introduce the following scaling: where λ ∼ 0 is called a regular breaking parameter. We study the system Z given in (2.1) in nonsmooth limit ε → 0 in the classical way, see e.g. [40]. We consider the extended fast-time system (z , ) = ( 2 Z, 0) and apply the cylindrical blow-up with r ≥ 0, (ȳ,¯ ) ∈ S 1 and¯ ≥ 0. Let F denote the vector field on (x, r, (ȳ,¯ )), i.e. the pullback of ( 2 Z, 0) under (3.2). We then perform desingularization by division of the right hand side by¯ 2 . In other words, it is F :=¯ −2 F that we shall study. To study the dynamics of F in a neighborhood of the cylinder, we use different charts. Based upon Section 3.1.1 and Section 3.1.2, we illustrate the transformation and the properties of F in Fig. 4(a).
3.1.1. Dynamics in the scaling chart¯ = 1. We consider the chart-specific coordinate y 2 defined by y = 2 y 2 , with (x, y 2 ) kept in a large compact subset of R 2 , > 0 and ∼ 0. Inserting this into the extended system (z , ) = ( 2 Z, 0) produces the following equations: The critical set of (3.4) is given by the union of two critical manifolds:  and the curve C given by Notice that at the point an intersection of H and C appears. All the singularities on H are nilpotent except for p 0 which is linearly zero. In the rest of this section we show that the slow-fast system (3.3) satisfies Assumptions T0-T6 in [14] along the critical curve C. Then we can use [14,Theorem 4] and prove Proposition 3.2 in Section 3.2. Theorem 4 says that the leading term of the integral of divergence of the vector field (3.3), computed along canard orbits near C between ξ(x) < 0 and x > 0, is I(x) 2 with I defined in (3.1) (see for example the exponent of the exponential term in (3.10)). This term remains dominant in the expression for the difference map of (2.1) near Γ x , see (3.8).
The singularities on C are normally attracting when x < 0, normally repelling when x > 0 and the slow dynamics X sl (x, λ 0 ) -given in (2.5) -is regular, pointing from the attracting part to the repelling part of C. Thus, if we denote the vector field in (3.3) by F S and if M x is any local C n center manifold of F 1 S := F S + 0 ∂ ∂ at normally hyperbolic singularity x ∈ C, then 1 2 F 1 S | Mx is a local flow box containing C and pointing from the left to the right. (The exponent 2 in the term 2 is often called the order of degeneracy.) This implies that Assumptions T0-T2 of [14] are satisfied. It remains to show that (3.3) satisfies Assumptions T3-T6 of [14] in an ( , λ)-uniform neighborhood of the turning point p 0 . In order to do that, we have to blow up the degenerate line H = {x = 0}, inside the slow-fast system (3.3), to a half-cylinder (see Fig. 4(b)). For the sake of readability of Section 3, we prove that Assumptions T3-T6 are satisfied in Appendix A.
Proposition 3.2. The transition maps ∆ ± have the following form: tend to zero as → 0, uniformly in (y, λ), and where Proof. We treat the forward transition map ∆ − (the backward transition map ∆ + can be studied in similar fashion). We split up the forward transition map ∆ − between S 0 and S 3 in three parts (see sketch of sections in Fig. 4(a)): (a) We define a section S 1 ⊂ {r 1 = r 10 } parametrized by x ∈ J ⊂ [µ − , 0[, J being a segment, and 1 ∈ [0, 0 r 10 ] where r 10 > 0 is a small constant and (x, r 1 , 1 ) are the coordinates of (3.5). The segment J is chosen large enough such that the transition map x = ∆ 01 (y, , λ) between S 0 and S 1 is well defined. Notice that 1 = r 10 . Since Z − has no singularities between S 0 and S 1 and the passage between S 0 and S 1 is located outside a fixed neighborhood of y = 0, it is clear that ∆ 01 is smooth in (y, , λ) (see also Lemma 2.1).
where f − , κ 1 and κ 2 are smooth, including at = 0, κ 1 = O( ) andĪ(x) < 0 is the slow divergence-integral of the form (3.9) computed betweenx and 0. We have the negative sign in front of the exponential term due to the chosen parametrization of S 2 and S 3 . Combining (a), (b) and (c), we obtain (3.8). We use that g 12 (∆ 01 (y, 0, 0), 0) = ξ − (y). Since Assumption T6 of [14] is satisfied (see Section A), the function f − − f + has the property given in Proposition 3.2, where f + is obtained in a similar way by studying the backward transition map ∆ + .

Proof of Theorem 1.1
To prove Theorem 1.1 we now use Theorem 3.1. We consider Z ± (·, λ) and suppose that (A4) holds with λ 0 = 0. Moreover, we will suppose that Z − is invariant under the symmetry (x, t) → (−x, −t) for λ = 0: Based upon the following simple result, this leads to a significant simplification of the calculations that follow. Proof. From the symmetry, we have that if (ξ − (y), 0) is the first intersection with Σ by the forward flow of (0, y) then (ξ + (y), 0) with ξ + (y) = −ξ − (y) is the first intersection with Σ by the backward flow.
In the following, while we continue to use () to denote the partial derivative with respect to x evaluated at (x, y, λ) = (0, 0, 0), we will also use () to indicate the second order partial derivative with respect to x also evaluated at (x, y, λ) = (0, 0, 0).
Proof. For simplicity write 0, 0). (g(u)). Notice that g and h both have "0/0" at x = 0, but each has a smooth extension to x = 0 due to the assumption of the two-fold by L'Hospital, recall also Proposition 2.3. In particular, Moreover, using assumption (A6).
In the same way, h(0) = 0 and We compute the partial derivatives of i(x) of even degree using the Faá di Bruno rule: where B m,n are the Bell polynomials. Here we have used that h(0) = 0. Each g (l) (0) can be written in terms of φ (y 2c ), . . . , φ (l) (y 2c ) (as well as the partial derivatives of Y ± ). This follows from the rule of inverse differentiation. To show the explicit expression for the coefficient of φ (2k) (y 2c ), including the expression for C 2k , we consider the term in (4.3) with n = m = 2k−1: using that B n,n (x) = x n . By the Leibniz integral rule, the result -including the stated properties of J 2k−1 -then follows.
Suppose that β > χ > 0, ξ = 0. Then it is a simple calculation to show that assumptions (A4), (A5) and (A6) all hold true for Z ± with Z + being quadratic of the form whereas Z − is linear of the form . (4.14) Upon invoking Theorem 3.1, we then conclude that for each k ∈ N Theorem 4.3 gives the existence of k limit cycles of z = Z(z, φ k (y −2 ), λ k c ( )) for all 0 < 1. In turn, this then completes the proof of Theorem 1.1.
The roots are very sensitive with respect to δ; increasing δ only slightly in each of our cases k = 4, 6, and 8 lead to fewer roots. For example for k = 8 we only find 5 roots for δ 9.449 × 10 −5 . In any case, δ > 0 has to be taken quite small to realize the desired number of roots. In turn, this implies that I is extremely small on the relevant domains; for k = 8 for example, we find (see Fig. 6(c)) that I(x) ∼ 10 −82 − 10 −85 ! We therefore expect -in line with [13] -that the desired number of limit cycles for 0 < 1 can also only be realized for extremely small values of > 0 and that these are therefore extremely difficult (if not impossible) to detect in numerical computations.
Appendix A. Blowing up the degenerate line H In this section we show that the system (3.3) from Section 3.1.1 satisfies Assumptions T3-T6 of [14] near the intersection p 0 of H with C. We write where Z ± = ( X ± , Y ± ), like in Section 2.2, and λ ∼ 0 is introduced in Section 3.1. We blow up H to a cylinder through the following blow-up transformation where (x,˜ ) ∈ S 1 and˜ ≥ 0. Again we will work with different charts. Let's first consider the end points of the normally attracting part and the normally repelling part of C on the edge of the cylinder {ρ =˜ = 0}. A.1. Dynamics in the phase directional chartsx = ±1. In the phase directional chart x = 1 we have (x, ) = (ρ, ρ˜ ). In these coordinates system (3.3) becomes, after division by ρ > 0,ẏ When ρ =˜ = 0, system (A.2) has a semi-hyperbolic singularity The eigenvalues of the linearization are the first eigenvalue being positive by Proposition 2.3 item (iv). Two-dimensional center manifolds of (A.2) at this singularity are transverse to the unstable manifold given by the y 2 -axis. Thus, the end point of the repelling part of C is normally hyperbolic (Assumption T3). Moreover, each center manifold of (A.2) at the singularity is the graph of Using the (ρ,˜ )-component of (A.2) we easily find the center behavior: Since X sl (0, λ 0 ) > 0, this system has, after division by˜ 2 , an isolated hyperbolic saddle (ρ,˜ ) = (0, 0) (Assumption T4). Notice that the exponent in˜ 2 is equal to the order of degeneracy mentioned in Section 3.1.1. Notice also that the center manifold, restricted to ρ = 0, is unique because (A.2) is of saddle type inside ρ = 0.
Remark A.1. In the framework of [14] a turning point is usually replaced with a sphere S 2 and Assumptions T3-T4 have to be satisfied at the end points of normally hyperbolic branches of the critical curve on the equator of the sphere. In our slow-fast setting (3.3) it is more convenient to use a cylindrical blow-up. This is not a problem because locally near the end points, located on the edge of the cylinder, one can use the normal linearization theorem of [6], like in [14].
A.2. Dynamics in the family chart˜ = 1. In the scaling chart we obtain x = x 2 . The system (3.3) changes intȯ upon desingularization (dividing the right hand side by ) and (subsequently) setting = 0. We used (A.1). For λ = 0, (A.3) has an invariant line γ defined by y 2 = y 2c with the dynamicṡ x 2 = X sl (0, λ 0 ) on it. The line γ is a heteroclinic connection on the cylinder connecting the end point of the attracting part of C to the end point of the repelling part of C (Assumption T5). See also Fig. 4(b).
To show that the invariant line γ breaks in a regular way as we vary λ ∼ 0 (Assumption T6), we follow [42, section 6.2] and extend the system by augmenting˙ λ = 0. In this formulation the center manifolds fromx = ±1 -that extend the attracting and repelling parts of C onto H -become two-dimensional and γ belongs to the intersection of these within ρ = 0 (where the manifolds are overflowing and unique). Write {(A.3),˙ λ = 0} in terms of dy 2 dx 2 , d λ dx 2 and consider the variational equations around the solution y 2 = y 2c , λ = 0 (corresponding to γ): du dx 2 = Ax 2 u + Bv, dv dx 2 = 0, Notice that A > 0, see Proposition 2.3 item (iv), and that B = 0 by (A4), see specifically (2.7). It is then straightforward to show, using the asymptotics of the error-function erf, see also [42,Lemma 6.2], that there are two linearly independent solutions of (A.4): with exponential growth for x 2 → ∞ (resp. −∞) and algebraic growth for x 2 → −∞ (resp. ∞). By [56, Proposition 4.2] the extended center manifolds therefore intersect transversally along γ, which completes the verification of Asssumption T6. 23 Appendix B. Transition maps near the hyperbolic edges In this section we study the transition map near the line of singularities {r 1 = 1 = 0} oḟ x = r 2 1 X(x, r 1 , 1 ), r 1 = −r 1 Y (x, r 1 , 1 ), where X and Y are smooth functions. We assume that Y (x, 0, 0) > 0 for each x ∈ J, J being a compact set. Notice that (x, 0, 0) is a set of equilibria and the linearization has ∓Y (x, 0, 0) as two nonzero eigenvalues.
Proposition B.1. Fix n ∈ N then there are constants r 10 > 0 and 10 > 0 small enough such that x → π x Q 1 (x, r 10 , 1 ) is C n uniformly and continuously in 1 . In particular, π x Q 1 (x, r 10 , 1 ) = g 0 (x) + O( 1 log −1 1 ), 1 → 0, with g 0 smooth and this expression can be differentiated with respect to x without changing the order of the remainder.
If X and Y in (B.1) depend smoothly on a parameter α, then π x Q 1 will also depend C n -smoothly on this parameter. This also follows from studying (B.2). We simply study the variational equations obtained by differentiating with respect to α and apply a similar estimation.