Saddle Slow Manifold: Definition and Computation
The results in this chapter have appeared in . As was explained in chapter 1, slow manifolds in slow-fast systems are perturbations of the associated invariant manifolds of the fast subsystem. The slow manifolds together with invariant manifolds of equilibria and periodic orbits organise the local and global dynamics. For instance, it is well known that the interaction of attracting and repelling slow manifolds can lead to canard explosions [8, 69, 96]. More recently, it has been estab-lished that saddle slow manifolds (SSMs) and their stable and unstable manifolds can play important roles in the dynamics of a system: the number of spikes in a bursting periodic orbit is organised by the intersection of the stable and unstable manifolds of an SSM . In , the e↵ect of changing a parameter on the number of spikes in the response of a system is investigated when a short-time stimulus is applied with fixed amplitude; a transition between solutions with di↵erent numbers of spikes occurs in an exponentially small parameter interval, and SSMs and their stable manifolds are an integral part of the mechanism for spike adding. Hence, the ability to compute accurate approximations to slow manifolds, including SSMs and their stable and unstable manifolds, is of significant interest. Slow manifolds experience extremely strong attraction or repulsion because of the fast dynamics normal to the manifolds. Hence, their numerical approximation is a challenge, and shooting methods are often unhelpful, because small errors in the initial conditions grow exponentially quickly. The computation of SSMs is perhaps even more challeng-ing because these manifolds have both repelling and attracting properties. There are well-established numerical methods for computing attracting and repelling slow mani-folds [22, 23, 25, 45] but methods for the approximation of SSMs are scarce. The first method for the computation of an SSM and its associated (un)stable manifolds was pre-sented in 2009 by Guckenheimer and Kuehn . There, the SSM is approximated using a collocation method; the corresponding stable (unstable) manifold is then computed as the union of trajectory segments integrated backward (forward) in time starting a small distance from the computed SSM in the direction of the stable (unstable) eigenvectors of the corresponding branch of equilibria of the fast subsystem. This idea is the same as in  and defines the (un)stable manifold of the SSM in terms of the corresponding manifold in the singular limit » = 0. Kristiansen  introduces an iterative method for computing slow manifolds and particularly an SSM by enforcing the invariance condition of a slow manifold as an algebraic equation and solving it from an initial guess. The asso-ciated (un)stable manifold of an SSM is computed by another iterative method through a projection onto the SSM. Basically, the computations are split into two nonsti ↵ parts: one on the SSM and the other as the connection to and from the SSM. The method has some advantages; for example, one only needs to know the vector field and its Jacobian, but its convergence is guaranteed only for small enough « . In this thesis, we compute an approximation of an SSM using pseudo-arclength continuation in Auto [27, 29] which is similar to the collocation method in . Our method for computing the associated (un)stable manifold of the SSM is based on a two-point boundary value problem (2PBVP) approach in the full system and does not rely on information from the fast subsystem. This is quite di↵erent from the methods presented in [24, 46]. Our method is fast and very accurate because of the set-up in Auto. In this chapter, we assume that the critical manifold is folded twice, resulting in a middle branch that is of saddle type. Corresponding to this saddle branch, there will be family of one-dimensional SSMs each with a corresponding family of two-dimensional stable and unstable manifolds . A trajectory started on such a stable manifold in a small neighbourhood of the SSM converges very quickly towards the SSM, follows it for time interval of O(1) on the slow time scale, and then diverges from the SSM along an unstable manifold, again very quickly . We approximate the stable manifold of an SSM as a one-parameter family of trajectory segments. The method is implemented in the software package Auto [27, 29] using pseudo-arclength continuation and a 2PBVP set-up. This chapter is organised as follows. In section 4.1, we define an approximation of an SSM and its stable and unstable manifolds. The implementation of the 2PBVP set-up for the computation of the stable manifold of an SSM in Auto is explained in section 4.2. We validate the computed stable manifold of an SSM of a polynomial example model  against the definition of these manifolds in section 4.3. Here, we also discuss the accuracy of the computed manifold in terms of its convergence properties as » ! 0.
Saddle slow manifolds and their (un)stable man-ifolds
If the critical manifold (1.5) is normally hyperbolic, Fenichel theory [40, 41] guarantees that C perturbs to a family of locally invariant manifolds with compatible stability properties, each of which is Cr and lies in an O(« )-neighbourhood of C for » suﬃciently small. Typically, C has folds with respect to the slow variable z, which means that there exist values of z for which the fast subsystem exhibits a saddle-node bifurcation. At such points, C is not normally hyperbolic, but we can divide C into several isolated branches, so that each of these branches is normally hyperbolic and gives rise to a corresponding family of slow manifolds. We are interested in branches of C that are of saddle type. We define a compact, connected submanifold S of C extending from z = zin to z = zout, where we choose zin < zout and assume that > 0 in (1.1) for all points p 2 S. Note that each point p 2 S is a hyperbolic saddle equilibrium of the fast subsystem with fixed z 2 [zin, zout], and p has exactly one negative and one positive eigenvalue. Hence, in the fast subsystem, p has a one-dimensional stable manifold, denoted W s(p), consisting of two trajectories that converge to p in forward time. Similarly, p has a one-dimensional unstable manifold, Figure 4.1(b) shows a sketch of S, together with (local) manifolds W s(S) and W u(S). The Stable Manifold Theorem and the smoothness of system (1.2) guarantee that W s(S) and W u(S) are also Cr-smooth; see [40, 58, 60]. Associated with S, provided » is small enough, Fenichel theory guarantees the ex-istence of a family of SSMs S« that are each locally invariant. Local invariance means that a solution started from a point in S« with z-coordinate z 2 (zin, zout) stays in S« until z = zout. While the theory does not guarantee uniqueness of S« , all manifolds in the family are exponentially close to one another [41, 70]. Fenichel theory [40, 41, 60, 70] also implies the existence of locally invariant Cr-smooth stable and unstable manifolds, denoted W s(S« ) and W u(S« ), associated with W s(S) and W u(S), respectively. Locally near S« , these manifolds lie in an O(« )-neighbourhood of their unperturbed counterparts. The (un)stable manifolds W s(S« ) and W u(S« ) are also not unique and exist as families of manifolds that lie exponentially close to one another. As mentioned in the introduction, each chosen pair W s(S« ) and W u(S« ) intersect in an SSM S« . Furthermore, there are trajectories that enter a small neighbourhood of each S« close to W s(S« ), and follow S« for a certain length of time, after which they leave close to W u(S« ). We approximate W s(S« ) by selecting a one-parameter family from those trajectories that follow S« up to z = zout; the unstable manifold W u(S« ) can be approximated in the same way by reversing time and considering z = zin. Figure 4.2 is a sketch of the local stable and unstable manifolds of S« . These surfaces are perturbations of the stable and unstable manifolds of S in Figure 4.1(b). For ease of visualization, we show just one sheet of each of W s(S« ) and W u(S« ). As indicated, solutions on W s(S« ) approach S« very fast at di↵erent values of the slow variable z along S« . The same thing happens for trajectories on W u(S« ) in reversed time.
Selecting a saddle slow manifold
We first provide a suitable definition of S« as a trajectory segment along S. Here, we are inspired by the theory of nonautonomous systems. We used [33, 103] as references, but see also the first papers on the subject [79, 92, 93, 95]. Let Bδ (z) denote a two-dimensional closed disk in the plane z = z with radius δ and centre (x, z) 2 S. Here, δ is small, but it must be at least of order « . We define which is a tubular compact set around S; see Figure 4.3(a). For » small enough, the intersection between Bδ (S) and the family of SSMs, which lies O(« ) from S, is not empty. Moreover, we can choose δ such that there is a set of trajectories, including S« , that enter Bδ (S) at zin, and leave Bδ (S) at zout. This is illustrated in Figure 4.3(a), where two such trajectories are sketched inside Bδ (S). We approximate S« \ Bδ (S) by the specific trajectory from this set that spends the longest time in Bδ (S); we denote this approximation by S« x. It is possible that there exists more than one trajectory with this property; we simply choose one of them. Note that, by definition, S« x can be parameterised by z 2 [zin, zout].
The (un)stable manifold of S« x
We now proceed with defining the approximation W s(S« x) of W s(S« ). Note that we only approximate one sheet of W s(S« ); the approximation of the other sheet is similar. As shown in Figure 4.3(b), we define a tubular neighbourhood of S« x similar to the way we defined Bδ (S). Specifically, we define where B (z« ) is now a disk of radius centred at a point (x« , z« ) = (x« (z« ), z« ) 2 S« x and is small. We formulate a definition of W s(S« x) as follows. The distance between a hyperbolic trajectory and the orbits on its stable (unstable) manifold goes to 0 as t ! 1 (t !1 ); see for example . For a finite-time trajectory, the stable (unstable) manifold is nonuniquely defined as the family of orbits with a strictly decreasing (increasing) distance to the trajectory in a small neighbourhood . Hence, the family of stable manifolds of S« x corresponds to the set of trajectories that enter B (S« x) (not necessarily from zin) and, as long as they are in the -neighbourhood of S« x, come closer to S« x. In other words, a trajectory φt(z« ) that enters B (S« x) at z« 2 [zin, zout) lies on a member of the family of stable manifolds if dz(‘t(z« ), S« x), the Euclidean distance between the intersections with the plane z = constant of ‘t(z« ) and S« x, decreases in forward time (increasing z) until φt(z« ) reaches the disk z = zout. We refer to such a trajectory as a converging trajectory. We emphasise that this definition also embraces all trajectories that enter B (S« x) from zin. Since we assume that the vector field is Cr-smooth, the distance function is also Cr-smooth, which is important for our definition of W s(S« x). We approximate a representative W s(S« x) of the stable manifold family of S« x as follows. We consider all the converging trajectories that enter B (S« x) at some fixed z-value, say z = z« , with zin < z« < zout. Of all these trajectories, φt(z« ) is chosen such that dz(φt(z« ), S« x) at z = zout is minimal, that is, at z = zout, the converging trajectory φt(z« ) lies closest to S« x. As for the definition of S« x, it is possible that there is more than one trajectory with the minimum distance; here, we also select only one of them. The union of all of these trajectories over di↵erent values of z« gives the manifold W s(S« x). The Cr-smoothness of the vector field enables us to choose the minimum-distance trajectory for di↵erent z« -values in a continuous manner. The extension of the selected trajectories Remark 4.1.1 We remark that we can define both sheets of W s(S« x) and S« x together as a single surface if we also include trajectories on W s(S« x) that enter B (S« x) from the disk z = zin. For example, this is possible by including all trajectories ‘t(zin, r), r 2 [0, ), that enter B (S« x) at z = zin with distance r to S« x, and have minimal dz(‘t(zin, r), S« x) at z = zout. As we will see in the next section, this extension is convenient from a computational perspective, because it means that W s(S« x) can be computed in its entirety as a single computation run.
When approximating W s(S« x), we would like to avoid the dependence on the choice of zin and zout and base the computations on a maximal possible interval [zin, zout]. In what follows, we begin with a set-up that follows the definition as given in section 4.1, but our final selection of trajectory segments that lie on W s(S« x) is determined by a maximum in total integration time, rather than a minimum in distance with respect to zout; this approach is similar to the numerical methods presented in chapter 3 for the detection of spike onset in the bursting patterns of neurons [82, 84]. Specifically, this means that we do not need to fix zout. Furthermore, zin is e↵ectively used to define a particular region of interest. We compute an approximation to W s(S« x) with the pseudo-arclength continuation package Auto [27, 29]. To this end, we set up a 2PBVP in section 4.2.1 that achieves the required approximation in three steps. In section 4.3, we show that our numerical approach leads to a good approximation of W s(S« x) as defined in section 4.1.
1.1 Background on geometric singular perturbation theory
1.2 Numerical methods for slow-fast systems
1.3 Examples of spike adding: an overview
1.4 Outline of the thesis
2 Background on spike adding
3 Spike-adding mechanisms in a polynomial model
3.2 The model and the transient response
3.3 Bifurcation analysis
3.4 The organisation of curves of spike onset
3.5 Spike-adding mechanism in the presence of extra equilibria
3.6 E↵ect of the injected current and its duration
3.7 Spike adding in another bursting pattern
3.8 Discussion and remarks
4 Saddle Slow Manifold:Definition and Computation
4.1 Saddle slow manifolds and their (un)stable manifolds
4.2 The algorithm
4.3 Accuracy of the method
4.4 Discussion and remarks
5 Numerical examples
5.1 A polynomial model
5.2 A thalamic neuron model
5.3 Chay–Cook model
5.4 Bursting periodic orbits in the Morris–Lecar model
5.5 Discussion and remarks
6 Conclusions and outlook
GET THE COMPLETE PROJECT
A thesis submitted in fulfilment of the requirements for the degree of Doctor of Philosophy in Mathematics, the University of Auckland, 2017.