Crossover Collision of Scroll Wave Filaments

Bernold Fiedler and Rolf M. Mantel

Received: December 20, 1999 Revised: December 19, 2000

Communicated by Alfred K. Louis

Abstract. Scroll waves are three-dimensional stacks of rotating spiral waves, with spiral tips aligned along filament curves. Such spatiotemporal patterns arise, for example, in reaction diffusion systems of excitable media type. We introduce and explore the crossover collision as the only generic possibility for scroll wave filaments to change their topological knot or linking structure. Our analysis is based on elementary singularity theory, Thom transversality, and abackwards uniqueness property of reaction diffusion systems. All phenomena are illustrated numerically by six mpeg movies downloadable at

http://www.mathematik.uni-bielefeld.de/documenta/vol-05/21.html

and, in the printed version, with six snapshots from each sequence.

1991 Mathematics Subject Classification: 35B05, 35B30, 35K40, 35K55, 35K57, 37C20 Keywords and Phrases: Parabolic systems, scroll wave patterns, scroll wave filaments, spirals, excitable media, crossover collision, singularity theory, Thom transversality, backwards uniqueness, video.

1 Introduction

Spatio-temporal scroll wave patterns have been observed both experimentally and in numerical simulations of excitable media in three space dimensions. See for example [36, 25, 20] and the references there. Typical experimental settings are the Belousov-Zhabotinsky reactions and its many variants. In two space dimensions, or in suitable planar sections through scroll wave patterns, rigidly rotating spiral wave patterns occur; see figure 1. For pioneering analysis motivated by propagation of electrical impulses in the heart muscle

Figure 1: Spiral wave patterns (model see section 6). Shown on the left is a rigidly rotating spiral wave with parameters as in section 6, on the right is a meandering spiral wave, with parameter a = 0.65 instead of a = 0.8. For color coding see section 7.

see [34], 1946. Meandering tip motions are also observed; see for example [35, 38, 5, 4] and the references there. There is some ambiguity in the definition of the tip of a spiral. It is an admissible definition in the sense of [13, sec.4], to associate tip positions (x_{1}, x_{2}) ∈ ℝ^{2} at time t ≥ 0 with the location of zeros of two components (u^{1}, u^{2}) of the solution describing the state of the system:

u = (u^{1}, u^{2})(t, x_{1}, x_{2}) = 0. (1.1)

In a typical excitable medium the values of (u^{1}, u^{2}) trace out a cycle as shown in figure 2, along x-circles around the spiral tip. In a singular perturbation setting, steep wave fronts are observed along these x-circles. Only near the spiral tip, these u-cycles shrink rapidly to the tip-value u = 0. This scenario, among other observations, motivated Winfree to attempt a phenomenological description in terms of states φ = u/|u| ∈ 𝕊^{1}, for (almost) all x ∈ ℝ^{2}, with remaining singularities of φ at the tip positions. In the present paper, we return to a reaction diffusion setting for u = u(t, x) ∈ ℝ^{2}, keeping in mind that the set u(t, x) = 0 is particularly visible, distinguished, and descriptively important – not as an “organizing center”which causes the global dynamics to follow its pace, but rather as a highly visible indicator of the global dynamics. In fact, defining tip positions by other nonzero levels (t, x) ≡ const., inside the cycle of figure 2, works just as well, and only reflects some of the ambiguity in the notion of “tip position”, as was mentioned above. With all our results below holding true, independently of such a shift of u-values, we proceed to work with u(t, x) = 0 as a definition of tip position. Scroll waves in three space dimensions x = (x_{1}, x_{2}, x_{3}) ∈ ℝ^{3} can be viewed

0.4 −0.4  0.4 −0.4  u^{1} u^{2} 0

Figure 2: A cycle of values (u^{1}, u^{2})(t, x_{0}) through a time-periodic wave front at a suitably fixed position x_{0} in an excitable medium (see section 6). Polar coordinates define a phase φ ∈ 𝕊_{1} along the dotted cycle.

as stacks of spiral waves with their tips aligned along a one-dimensional curve called the tip filament. As in the planar case, the tip filament may move around in ℝ^{3}, and the associated sectional spirals may continuously change their shapes and their mutual phase relations with time. Denoting by (u^{1}, u^{2}) two components of the solutions of the associated reaction diffusion systems, again, we can consider filaments φ^{t} as given by the zero sets

u = (u^{1}, u^{2})(t, x_{1}, x_{2}, x_{3}) = 0. (1.2)

We use two components here because the local dynamics of excitable media are essentially two-dimensional. More precisely, for each fixed time t > 0 the filaments φ^{t} describe the zeros x ∈ ℝ^{3} of the solution profile

x ↦ u(t, x). (1.3)

In other words, the filament φ^{t} is the zero level set of the solution profile u(t, ·) at time t. Suppose zero is a regular value of u(t, ·), that is, the x-Jacobian u_{x}(t, ·) possesses maximal rank 2 at any zero of u. Then the filaments φ^{t} consist of embedded curves in ℝ^{3}, by the implicit function theorem. Moreover the filaments depend as smoothly on t as smoothness of the solution u permits. Therefore, collision of filaments can occur only if the rank of u_{x}(t, ·) drops. To

Figure 3: A scroll wave and its filament. The band is tangential to the wave front at the filament.

a.) t < t_{0}  b.) t = t_{0}  c.) t > t_{0}

Figure 4: Crossover collision of oriented filaments at time t = t_{0}

analyze the simplest possible case, we assume

\begin{aligned}& u(t_{0}, x_{0}) & = & 0, \\ co-rank & u_{x}(t_{0}, x_{0}) & = & 1. \end{aligned} (1.4)

Let P denote a rank one projection along range u_{x}(t_{0}, x_{0}) onto any complement of that range. Let E = ker u_{x}(t_{0}, x_{0}) denote the two-dimensional null space of the 2 × 3 Jacobean matrix u_{x}. We assume the following non-degeneracy conditions for the time-derivative u_{t} and the Hessian u_{xx}, restricted to E:

P u_{t}(t_{0}, x_{0})  ≠ 0, and
P u_{xx}(t_{0}, x_{0})|_{E} is strictly indefinite.

(1.5)

A specific example u(t, x) satisfying assumptions (1.4), (1.5) at t = t_{0}, x_{0} = 0 is given by

\begin{aligned}u^{1}(t, x) & = & (t − t_{0}) + x_{1}^{2} − x_{2}^{2} \\ u^{2}(t, x) & = & x_{3}. \end{aligned} (1.6)

In figure 4 we observe the associated crossover collision of filaments in projection onto the null space E: at t = t_{0} two filaments collide, and then reconnect. Note that after collision the two filaments do not reconnect as before, re-establishing the previous filaments. Instead, they cross over, forming bridges between originally distinct filaments. Figure 4 describes the universal unfolding, by the time “parameter” t, of a standard transcritical bifurcation in x-space. In fact, suppose u(t, x) satisfies assumptions (1.4), (1.5). Then there exists a local diffeomorphism

\begin{aligned}τ & = & τ(t) \\ ξ & = & ξ(t, x) \end{aligned} (1.7)

mapping (t_{0}, x_{0}) to τ_{0} = t_{0}, ξ_{0} = 0, such that the original zero set transforms to that of example (1.6), rewritten in (τ, ξ)-coordinates. This follows from Lyapunov-Schmidt reduction and elementary singularity theory; see for example [15]. In an early survey, Tyson and Strogatz [31] hinted at topologically consistent changes of the connectivity of oriented tip filaments, as a theoretical possibility. The point of the present paper is to identify specific singularities, in the sense of singularity theory, which achieve such changes and which, in addition, are generic with respect to the initial conditions of general reaction diffusion systems. Genericity refers to topologically large sets. These sets contain countable intersections of open dense sets, and are dense. We caution our PDE readers here that we are not addressing issues like loss of regularity (smoothness) or development of singularities in a blow-up sense. Genericity is based on perturbations of only the initial conditions. We do not require any perturbations of the underlying partial differential equations themselves. We consider it a fundamental idea to study solutions u(t, x) of partial differential equations, qualitatively, by investigating the singularities of their level sets – possibly for all, or at least for generic initial conditions. Such an idea is already present in work by Schaeffer, [27], and more recently by Damon, [7], [8], [9] and the references there. In view of example (2.12) for linear scalar parabolic equations in one space dimension below, the first relevant example can even be attributed to Sturm [28], 1836. For present day relevance of Sturm’s observations, once motivated by Sturm-Liouville theory, see also [3], [12], [23]. The work by Schaeffer addresses level sets of strictly convex scalar hyperbolic conservation laws in one space dimension. His analysis is based on the variational formulation due to Lax: for almost every (t, x) the solution u(t, x) appears as the pointwise minimizer of a given function, which involves the initial conditions u_{0}(x) explicitly. The backwards uniqueness problem, a somewhat delicate technical point for our parabolic systems, is circumvented by the explicit Lax formula in his context.

Damon’s work is motivated by Gaussian blurring and by applications of the linear heat equation to image processing, but applies to a large class of differential operators. Unfortunately, the partial differential equations are viewed as purely local constraints on the k-jet of “solutions”. Neither initial nor boundary conditions are imposed on these “solutions”. Genericity is understood purely in the space of smooth such “solutions”. The important nonlocal PDE issue of genericity in terms of initial conditions, as addressed in our present paper, has not been resolved by Damon’s approach. In contrast to these abstract results, strongly in the spirit of pure singularity theory, our motivation is the global qualitative dynamics of reaction diffusion systems. In particular, we do require our solutions u = u(t, x) to not only satisfy the underlying partial differential equations near (t_{0}, x_{0}) but also the respective initial and boundary conditions. For a technically detailed statement see our main result, theorem 2.1 below. As a consequence, the crossover of filaments just described is the one and only non-destructive collision of filaments possible – for a generic set of initial conditions. See theorem 2.2. The remaining sections are organized as follows. Preparing for the proof of theorem 2.1, we provide an abstract jet perturbation lemma in section 3 which is based on backwards uniqueness results for linear, non-autonomous parabolic systems. In section 4, we prove theorem 2.1 using Thom’s jet transversality theorem. Moreover we present a generalization to the vector case u ∈ ℝ^{m}, m ≥ 2, in corollary 4.2. Theorem 2.2 is proved in section 5. Section 6 summarizes a fast numerical method, due to [11, 22], for time integration of a specific excitable medium with steep fronts in three space dimensions. In section 7 we adapt this method to compute filaments and their associated local isochrone phase bands. We conclude with numerical examples illustrating crossover collisions in autonomous and periodically forced reaction diffusion systems, including the unlinking of linked twisted scroll rings and the unknotting of a trefoil torus knot filament; see section 8. Acknowledgment. Both authors are grateful to the Institute of Mathematics and its Applications (IMA), Minneapolis, Minnesota. The main part of this work was completed there during a PostDoc stay of the second author and several visits of the first author as senior visiting scientist during the special year ”Emerging Applications of Dynamical Systems”, 1997/98. We are indebted to Jim Damon for helpful discussions, and to the referee for additional references. We thank Martin Rumpf and Peter Serocka for help with visualization. Support by the Deutsche Forschungsgemeinschaft is also gratefully acknowledged.

2 Main Results

For a technical setting we consider semilinear parabolic systems

u_{t}^{i} = div_{x}(d^{i}(t, x)∇_{x}u^{i}) + f^{i}(t, x, u, ∇_{x}u) (2.1)

throughout the present paper. Here u = (u^{1}, …, u^{m}) ∈ ℝ^{m}, x = (x_{1}, …, x_{N}) ∈ Ω ⊂ ℝ^{N}. The data d^{i}, f^{i} are smooth with uniformly positive definite diffusion matrices d^{i}. The bounded open domain Ω is assumed to have smooth boundary. Inhomogeneous mixed linear boundary conditions

α_{i}(x)u^{i}(t, x) + β_{i}(x)∂_{ν}u^{i}(t, x) = γ(x) (2.2)

with smooth data and α_{i}, β_{i} ≥ 0, α_{i}^{2} + β_{i}^{2} ≡ 1 are imposed. Periodic boundary conditions are also admissible, as well as uniformly parabolic semilinear equations on compact manifolds with smooth boundaries, if any. The solutions

u = u(t, x; u_{0}) (2.3)

of (2.1), (2.2) with initial condition

u(0, x; u_{0}) := u_{0}(x) (2.4)

define a local semi-evolution system in the phase space X of profiles u_{0}(·) in any of the Sobolev spaces W^{k′,p}(Ω), k′ > N/p, which satisfy the boundary conditions (2.2); see [16] for a reference. By the smoothing property of the parabolic system, solutions are in fact smooth in their maximal open intervals of existence t ∈ (0, t_{+}(u_{0})) and depend smoothly on u_{0} ∈ X, both when viewed pointwise and when viewed as x-profiles u(t, ·; u_{0}) ∈ X. To address the issue of singularities u(t_{0}, x_{0}) = 0, in the sense of singularity theory, we consider the jet space J_{x}^{k} of Taylor-polynomials in x = (x_{1}, …, x_{N}) ∈ ℝ^{N} of degree at most k, with real coefficients and vector values u ∈ ℝ^{m}. Defining the k-jet j_{x}^{k}u with respect to x at (t_{0}, x_{0}) as

(j_{x}^{k}u)(t_{0}, x_{0}) := (u, ∂_{x}u, …, ∂_{x}^{k}u)(t_{0}, x_{0}), (2.5)

Taylor expansion at x_{0} allows us to interpret j_{x}^{k}u(t_{0}, x_{0}) as an element of our linear jet space J_{x}^{k} satisfying

u(t_{0}, x_{0}) = 0. (2.6)

Here and below, we assume that k′ > k + N/p so that the evaluation

u ↦ j_{x}^{k}u(t_{0}, x_{0}) (2.7)

becomes a bounded linear map from X to J_{x}^{k}, by Sobolev embedding. On the level of k-jets, a notion of equivalence is induced by the action of local C^{k}-diffeomorphisms x ↦ Φ(x), u ↦ Ψ(u) fixing the origins of x ∈ ℝ^{N}, u ∈ ℝ^{m}, respectively. Indeed, for any polynomial p(x) ∈ J_{x}^{k} with p(0) = 0, we may consider the transformed polynomial

j_{x}^{k}(Ψ ∘ p ∘ Φ) ∈ J_{x}^{k}. (2.8)

We call the jet (2.8) contact equivalent to j_{x}^{k}p = p; see for example [15].

By a variety S ⊂ ℝ^{ℓ} we here mean a finite disjoint union

S = ⋃_{j=0}^{j_{0}} S_{j} (2.9)

of embedded submanifolds S_{j} ⊂ ℝ^{ℓ} with strictly decreasing dimensions such that S_{j_{1}} ∪ … ∪ S_{j_{0}} is closed for any j_{1}. We call codim_{ℝ^{ℓ}} S_{0} the codimension of the variety S in ℝ^{ℓ}. Similarly, by a singularity (in the sense of singularity theory) we mean a variety S ⊂ J_{x}^{k} in the sense of (2.9), which satisfies u = 0 and is invariant under any of the contact equivalences (2.8). Let codim_{J_{x}^{k}} S denote the codimension of S, viewed as a subvariety of J_{x}^{k}. Shifting codimension by N = dim x for convenience we call

(2.10) codim S := (codim_{J_{x}^{k}} S) − N

the codimension of the singularity S. For example, a typical map (t_{0}, x_{0}) ↦ j_{x}^{k}u(t_{0}, x_{0}) with x_{0} ∈ ℝ^{N}, u ∈ ℝ^{m} will miss singularities of codimension 2 or higher. In contrast, the map can be expected to hit singularities S of codimension 1 at isolated points t = t_{0}, and for some x_{0} ∈ ℝ^{N}. Having shifted codimension by N in (2.10) therefore conveniently allows us to observe that typical profiles of functions u(t, ·) miss singularities of codimension 2 entirely, and encounter such singularities of codimension 1, anywhere in x ∈ ℝ^{N}, only at discrete times t. We aim to show that this simple arithmetic also works for PDE solutions u(t, x) under generic initial conditions. Since the geometrically simple issue of codimension is overloaded with – sometimes conflicting – definitions in singularity theory, we add some examples which illustrate our terminology. First consider the simplest case

S = {u = 0} ⊂ J_{x}^{k}. (2.11)

where u(t, ·) : ℝ^{N} → ℝ^{m}. Then codim S = m − N. For systems of m = 2 equations in N = 0 space dimensions, that is, for ordinary differential equations in the plane, typical trajectories fail to pass through the origin in finite time: codim S = 2. For N = 1, we can expect the solution curve profile u(t, ·) to pass through the origin at certain discrete times t_{0} and positions x_{0}, because codim S = 1. For N = 2 we have codim S = 0. We therefore expect isolated zeros to move continuously with time: see our intuitive description of planar spiral waves in section 1 and figure 1. Since codim S = −1 for N = 3, we expect zeros of u(t_{0}, ·) to occur along one-dimensional filaments, even for fixed t_{0}. This is the case of scroll wave filaments φ^{t_{0}} in excitable media. Next we consider a scalar one-dimensional equation, m = N = 1. Multiple zeros are characterized by

S = {u = 0, u_{x} = 0}, (2.12)

x_{0} x′_{0} x  t t_{0}  t′_{0}

Figure 5: Saddle-node singularities of codimension 1.

a set to which we ascribe codimension 1. Indeed, we can typically expect a pair of zeros to coalesce and disappear as in (t_{0}, x_{0}) of figure 5. The opposite case, a pair creation of zeros as in (t′_{0}, x′_{0}), does not occur for scalar nonlinearities f satisfying f(t, x, 0, 0) = 0. This observation, going back essentially to Sturm [28], conveys considerable global consequences for the associated semiflows; see for example [12] and the references there. Passing to planar 2-systems, m = N = 2, the same saddle-node bifurcations of figure 5 could for example correspond to annihilation and creation of a pair of tips of counter-rotating spirals, respectively. We conclude our series of motivating examples with the singularity (1.4) of filament collision in systems satisfying N = m + 1:

S = {u = 0, co-rank u_{x} ≥ 1}. (2.13)

Note that codim S = 1. For the stratum S_{0} of S with lowest codimension we can assume that the quadratic form P u_{xx}|_{E} is indeed nondegenerate, in the notation of (1.5). Under the additional transversality assumption P u_{t} ≠ 0, the strictly indefinite case was discussed in section 1. It leads to crossover collisions, which are our main applied motivation here. The strictly definite case, positive or negative, leads to creation/annihilation of small circular filaments. For a numerical realization of the associated scroll ring annihilation we refer to the simulation in figure 8. After our intermezzo on singularities we now address genericity. We say that a property of solutions u(t, x; u_{0}) of our semilinear parabolic system (2.1) – (2.4) holds for generic initial conditions u_{0} ∈ X if it holds for a generic subset of initial conditions. Here subsets are generic (or residual) if they contain a countable intersection of open dense subsets of X. Recall that generic subsets and countable intersections of generic subsets are dense in complete metric spaces X, by Baire’s theorem; see [10, ch. 12]. With these preparations we can now state our main result concerning solutions u(t, x) of our parabolic system (2.1) – (2.4) with generic initial conditions u_{0} ∈

t t′_{0} t_{0} x′_{0} x_{0} x

Figure 6: Annihilation (left) and creation (right) of closed filaments

X ⊂ W^{k′,p} ↪ C^{k}. As before 0 ≤ t < t_{+}(u_{0}) denotes the maximal interval of existence. Finally, we recall that a map ρ : V → J between Banach spaces is transverse to a variety S = S_{0} ∪ … ∪ S_{j_{0}}, in symbols:

ρ ∩ ⊤ S, (2.14)

if ρ(v) ∈ S_{j} implies

T_{ρ(v)}S_{j} + range Dρ(v) = J; (2.15)

see for example [1, 19].

Theorem 2.1 For some fixed k ≥ 1, consider a finite collection of singularities S^{i} ⊂ J_{x}^{k}, each of codimension at least 1. Then the following holds true for solutions u(t, x) of (2.1) – (2.4) with generic initial conditions u_{0} ∈ X. Singularities S^{i} with

codim S^{i} ≥ 2 (2.16)

are not encountered at any (t_{0}, x_{0}) ∈ (0, t_{+}(u_{0})) × Ω. In other words, j_{x}^{k}u(t_{0}, x_{0}) ∈ S^{i} for some 0 < t_{0} < t_{+}(u_{0}), x_{0} ∈ Ω implies codim S^{i} = 1. The map

\begin{matrix}(0, t_{+}(u_{0})) × Ω & → & J_{x}^{k} \\ (t_{0}, x_{0}) & ↦ & j_{x}^{k}u(t_{0}, x_{0}) \end{matrix} (2.17)

is in fact transverse to each of the varieties S^{i}. In particular, the points (t_{0}^{n}, x_{0}^{n}) where the solution u(t, x) encounters singularities S^{i} of codimension 1 are isolated in the domain [0, t_{+}(u_{0}))×Ω of existence. Although there can be countably many singular points (t_{0}^{n}, x_{0}^{n}) accumulating to the boundary t_{+}(u_{0}) or ∂Ω, the values t_{0}^{n} are pairwise distinct.

Theorem 2.2 For some fixed k ≥ 1, consider solutions u(t, x) of (2.1) – (2.4) with N = 3, m = 2, that is with x ∈ Ω ⊂ ℝ^{3} and u(t, x) ∈ ℝ^{2}. Then for generic initial conditions u_{0} ∈ X the following holds true. Except for at most countably many times t = t_{0}^{n} ∈ (0, t_{+}(u_{0})), the filaments

{x ∈ Ω | u(t, x) = 0} (2.18)

are curves embedded in Ω, possibly accumulating at the boundary. At each exceptional value t = t_{0}^{n}, exactly one of the following occurs at a unique location x_{0}^{n} ∈ Ω:

(i) a creation of a closed filament, or

(ii) an annihilation of a closed filament, or

(iii) a crossover collision of filaments.

For cases (i),(ii) see figures 6, 8; for case (iii) see figures 4, 9–13, and (1.4) – (1.6).

3 Jet Perturbation

In this section we prove a perturbation result, lemma 3.1, which is crucial to our proof of theorem 2.1. We work in the technical setting of semilinear parabolic systems (2.1) – (2.4) with associated evolution

u = u(t, x; u_{0}) (3.1)

on the phase space X of W^{k′,p}(Ω)-profiles u(t, ·,; u_{0}) satisfying Robin boundary conditions (2.2). Let k′ − \frac{N}{p} > k ≥ 1, to ensure the Sobolev embedding X ↪ C^{k}(Ω). Let

D := {(t, x, u_{0}) | x ∈ Ω, u_{0} ∈ X, 0 < t < t_{+}(u_{0})} (3.2)

denote the interior of the domain of definition.

Lemma 3.1 The map

\begin{matrix}j_{x}^{k}u: D & → & J_{x}^{k} \\ (t, x, u_{0}) & ↦ & j_{x}^{k}u(t, x; u_{0}) \end{matrix} (3.3)

is a C^{κ} map, for any κ. For any (t, x, u_{0}) ∈ D, the derivative

D_{u_{0}}j_{x}^{k}u(t, x; u_{0}): X → J_{x}^{k} (3.4)

is surjective.

Proof: The regularity claim follows from smoothness of the data d^{i}, f^{i}, α_{i}, β_{i} and the smoothing action of parabolic systems; see for example [16, 26, 29, 14, 21]. To prove surjectivity of the linearization (3.4) with respect to the initial condition, we essentially follow [16]. First observe that for any fixed x_{0} ∈ Ω the linear evaluation map

\begin{matrix}j_{x}^{k}: X & → & J_{x}^{k} \\ v & ↦ & j_{x}^{k}v(x_{0}) \end{matrix} (3.5)

is bounded, because X ↪ C^{k}(Ω), and trivially surjective. Moreover, the jet space J_{x}^{k} is finite-dimensional. It is therefore sufficient to show that the linearization

\begin{matrix}D_{u_{0}}u(t, ·; u_{0}) : X & → & X \\ v_{0} & ↦ & v(t·) \end{matrix} (3.6)

possesses dense range, for all u_{0} ∈ X, 0 < t_{0} < t_{+}(u_{0}). Here v(t, ·) satisfies the linearized parabolic system

v_{t}^{i} = div_{x}(d_{i}(t, x)∇_{x}v^{i}) + f_{p}^{i} · ∇_{x}v + f_{u}^{i} · v (3.7)

with boundary conditions (2.2) for v and initial condition v(0, ·) = v_{0}. The partial derivatives f_{p}^{i}, f_{u}^{i} of the nonlinearity f = f(t, x, u, p) are to be evaluated along (t, x, u(t, x), ∇_{x}u(t, x)). To show the density of range D_{u_{0}}u(t, ·; u_{0}) in X, we now proceed indirectly. Suppose

clos_{X}D_{u_{0}}u(t_{0}, ·; u_{0})X ≠ X. (3.8)

Then X contains a nonzero element w(t_{0}, ·) in the L^{2}-orthogonal complement of D_{u_{0}}u(t, ·; u_{0})X in X. Consider the associated solution w(t, ·) ∈ X of the formal adjoint equation

w_{t}^{i} = −div_{x}(d_{i}(t, x)^{T} ∇_{x}w^{i}) + ∑_{j} div_{x}(w^{j}f_{p_{i}}^{j}) − (f_{u}^{T} w)_{i} (3.9)

for 0 ≤ t ≤ t_{0}, still with boundary conditions (2.2) but with “initial” condition w(t_{0}, ·) at t = t_{0}. We again use the notation f_{p_{i}}^{j} for the partial derivative of f^{j} with respect to ∇u^{i}, here. Direct calculation shows that scalar products ⟨·, ·⟩ between solutions v(t, ·) of the linearization (3.7) and solutions w(t, ·) of its formal adjoint (3.9) in L^{2}(Ω) are time-independent. Therefore, by construction of w(t_{0}, ·)

⟨v(t, ·), w(t, ·)⟩_{L^{2}(Ω)} = ⟨v(t_{0}, ·), w(t_{0}, ·)⟩ = 0, (3.10)

for all 0 ≤ t ≤ t_{0}. Evaluating at t = 0, v(0, ·) = v_{0} ∈ X, we conclude

⟨v_{0}, w(0, ·)⟩_{L^{2}(Ω)} = 0 (3.11) for all v_{0} ∈ X, and hence

w(0, ·) = 0. (3.12)

In other words, the backwards parabolic system (3.9) possesses a solution w(t, ·) which starts nonzero at t = t_{0} > 0 but ends up zero at t = 0. This is a contradiction to the so-called backwards uniqueness property of parabolic equations. See for example [14], [16] and the references there. By contradiction, we have therefore proved that

clos_{X}D_{u_{0}}u(t_{0}, ·; u_{0})X = X, (3.13)

contrary to our indirect assumption (3.8). This completes the indirect proof of the perturbation lemma. ▷◁

4 Proof of Theorem 2.1

Our proof of theorem 2.1 is based on Thom’s transversality theorem [30, 1]. For convenience we first recall a modest adaptation of the transversality theorem, fixing notation. We use the concept of transversality of a map ρ to a variety S as explained in (2.9), (2.14), (2.15). The proof is based on Sard’s theorem and is not reproduced here.

Theorem 4.1 [Thom transversality] Let X be a Banach space, D ⊆ ℝ^{ℓ} × X open and

\begin{matrix}ρ : D & → & ℝ^{ℓ′} \\ (y, u_{0}) & ↦ & ρ(y, u_{0}) \end{matrix} (4.1)

a C^{κ}-map. Let S ⊂ ℝ^{ℓ′} be a variety and assume

ρ ∩ ⊤ S, (4.2)

(4.3) κ > max{0, ℓ − codim_{ℝ^{ℓ′}} S}.

Then the set

(4.4) X_{S} := {u_{0} ∈ X | ρ(·, u_{0}) S, where defined}

is generic in X (that is: contains a countable intersection of open dense sets).

The point of the theorem is, of course, that in X_{S} transversality to S is achieved, for fixed u_{0}, by varying only y in ρ(y, u_{0}). For example, u_{0} ∈ X_{S} and codim_{ℝ^{ℓ′}} S > ℓ imply

ρ(y, u_{0}) ∉ S (4.5) whenever y is such that (y, u_{0}) ∈ D. This follows immediately from condition (2.15) on transversality. In other words, for generic u_{0} the image of ρ(·, u_{0}) misses varieties of sufficiently high codimension. We now use theorem 4.1 to prove our main result, theorem 2.1. We consider the jet evaluation map

ρ(t, x, u_{0}) := j_{x}^{k}u(t, x; u_{0}) (4.6)

of the evolution u(t, ·; u_{0}) associated to our parabolic system; see (2.1) – (2.5). We choose D to be the (open) domain of definition

D = {(t, x, u_{0}) | 0 < t < t_{+}(u_{0}), x ∈ Ω, u_{0} ∈ X} (4.7)

of the evolution; clearly y = (t, x) ∈ ℝ^{N+1} so that ℓ = N + 1. For the variety S we choose, successively, any of the finitely many singularities S^{i} ⊂ J_{x}^{k} of theorem (2.1). Their codimensions as subvarieties of J_{x}^{k} ≅ ℝ^{ℓ′} are

(4.8) codim_{J_{x}^{k}} S^{i} = N + codim S^{i};

see (2.10). Note that assumptions (4.2) and (4.3) both hold, independently of the choice of k for the varieties S^{i} ⊆ J_{x}^{k}, by lemma 3.1. Claim (2.17) about transversality of (t_{0}, x_{0}) ↦ u(t_{0}, x_{0}; u_{0}) to any singularity S^{i} is now just the statement of theorem 4.1. Next, we prove that singularities S^{i} with codim S^{i} ≥ 2 are missed altogether, for generic initial conditions u_{0} ∈ X, as was claimed in (2.16). We evaluate (4.8) to yield

(4.9) codim_{J_{x}^{k}} S^{i} = N + codim S^{i} ≥ N + 2 > N + 1 = ℓ

In view of example (4.5), this proves our claim (2.16): generically, only singularities S^{i} with codim S^{i} = 1 are encountered. Now we prove that the positions (t_{0}^{n}, x_{0}^{n}), where singularities S^{i} with codim S^{i} = 1 are encountered, are generically isolated in [0, t_{+}(u_{0})) × Ω. Indeed assuming j_{x}^{k}u_{0} ∉ S^{i}, we have t_{0}^{n} > 0 without loss of generality. Since the lowerdimensional strata S_{j}^{i}, j ≥ 1 of the singularity S^{i} are of (singularity) codimension ≥ 2, they are missed by solutions entirely, for generic initial conditions u_{0}. Therefore

j_{x}^{k}u(t_{0}^{n}, x_{0}^{n}; u_{0}) ∈ S_{0}^{i} (4.10)

only hit the maximal strata, staying away from the closed union of lowerdimensional strata, uniformly in compact subsets of [0, t_{+}(u_{0})) × Ω. Because the S_{0}^{i} are finitely many embedded submanifolds of codimension N +1 in J_{x}^{k} and because the crossings (4.10) are transverse, the corresponding crossing points (t_{0}^{n}, x_{0}^{n}) are also isolated in [0, t_{+}(u_{0})) × Ω, as claimed. It remains to show that the values t_{0}^{n} are mutually distinct for generic initial conditions u_{0} ∈ X. To this end we consider the augmented map

\begin{matrix}ρ̃: D̃ → J_{x}^{k} × J_{x}^{k} \\ (t, x_{1}, x_{2}, u_{0}) → (j_{x}^{k}u(t, x_{1}; u_{0}), j_{x}^{k}u(t, x_{2}; u_{0})) \end{matrix} (4.11)

on the open domain

(4.12)  D̃ := {(t, x_{1}, x_{2}, u_{0}) | 0 < t < t_{+}(u), x_{1}, x_{2} ∈ Ω, x_{1} ≠ x_{2}, u_{0} ∈ X}.

To apply Thom’s transversality theorem 4.1, we only need to check the transversality assumption (4.2). In fact we show

ρ̃ ∩ ⊤ {0} ∈ J_{x}^{k} × J_{x}^{k}. (4.13)

This follows, analogously to lemma 3.1, from x_{1} ≠ x_{2} and the fact that the linearization D_{u_{0}}u(t_{0}, ·; u_{0}) possesses dense range in X; see (3.6) – (3.13). We can therefore apply theorem 4.1 to ρ̃ with respect to the varieties

S̃ := S^{i_{1}} × S^{i_{2}}. (4.14)

In J_{x}^{k} × J_{x}^{k}, these varieties have codimension

(4.15) codim_{J_{x}^{k}×J_{x}^{k}} S̃ = 2N + codim S^{i_{1}} + codim S^{i_{2}} = 2N + 2

Since this number exceeds

dim(t, x_{1}, x_{2}) = 2N + 1, (4.16)

the variety S̃ is missed by ρ̃(·, ·, ·; u_{0}), for generic u_{0} ∈ X. See example (4.5) again. Therefore the times t_{0}^{n} where singularities S^{i} can occur are pairwise distinct for generic initial conditions, completing the proof of theorem 2.1. ▷◁

Reviewing the proof of theorem 2.1, which hinges crucially on the transversality statement (3.4) of our jet perturbation lemma 3.1, we state an easy generalization which is important from an applied viewpoint. Suppose that only m′ ≤ m profiles (or m′ linear combinations) out of the m profiles u = (u^{1}, ..., u^{m})(t, x) are observable:

û := P̂u, (4.17)

for some linear rank m′ projection of ℝ^{m}. Then û(t, x; u_{0}) may encounter certain singularities Ŝ^{i} in the space Ĵ_{x}^{k} of k-jets with values in range P̂.

Corollary 4.2 Under the assumptions of theorem 2.1 and in the above setting, theorem 2.1 remains valid, verbatim, for singularities Ŝ^{i} ⊂ Ĵ_{x}^{k} of the k-jets j_{x}^{k}û(t, x) of the observables û := P̂u. We emphasize that codimensions of Ŝ^{i} are then to be computed in Ĵ_{x}^{k}.

Proof: Acting on the dependent variables (u^{1}, …, u^{m}), only, the projection P̂ lifts to a projection P̂_{k} from J_{x}^{k} onto Ĵ_{x}^{k} such that

j_{x}^{k}P̂u(t, x; u_{0}) = P̂_{k}j_{x}^{k}u(t, x; u_{0}) (4.18)

Therefore the surjectivity property (3.4) of lemma 3.1 remains valid for

D_{u_{0}}j_{x}^{k}û(t, x; u_{0}): X → Ĵ_{k}. (4.19)

Repeating the proof of theorem 4.1, now on the level of û, Ĵ_{x}^{k}, Ŝ^{i}, proves the corollary. ▷◁

5 Proof of Theorem 2.2

To prove theorem 2.2 we invoke theorem 2.1 for x ∈ Ω ⊂ ℝ^{3}, u(t, x) ∈ ℝ^{2}, and appropriate singularities S^{i} ⊂ J_{x}^{k} of singularity codimension 1, in the sense of (2.10). We first consider the case that 0 is a regular value of u(t, ·) on Ω, that is

rank u_{x}(t_{0}, x_{0}) = 2 (5.1)

is maximal, whenever u(t_{0}, x_{0}) = 0, 0 < t_{0} < t_{+}(u_{0}), x_{0} ∈ Ω. Then the filament

{x ∈ Ω | u(t_{0}, x) = 0} (5.2)

is an embedded curve in Ω, as claimed in (2.18). Next consider the case

rank u_{x}(t_{0}, x_{0}) ≤ 1. (5.3)

Let S ⊂ J_{x}^{k=2} be the set of those 2-jets (u, u_{x}, u_{xx}) ∈ J_{x}^{k=2} satisfying u = 0 and rank u_{x} = 1. Clearly S is a singularity in the sense of (2.9), (2.10) and

codim S = 1 (5.4)

as was discussed in example (2.13). We recall that the maximal stratum S_{0} of S, determining the codimension, is given by the conditions

rank u_{x} = 1,
P u_{xx}|_{E} nondegenerate.

(5.5)

Here E := ker u_{x} denotes the kernel and P denotes a projection in ℝ^{2} onto a complement of the range of the Jacobian u_{x}. In view of example (2.13) and section 1, nondegeneracy of P u_{xx}|_{E} gives rise to the three cases (i) - (iii) of corollary 2.2, via theorem 2.1, if only we show that

P u_{t}(t_{0}, x_{0}) ≠ 0 (5.6)

whenever j_{x}^{2}u(t_{0}, x_{0}) ∈ S. By theorem 2.1, we have

j_{x}^{2}u(·, ·) ∩ ⊤ S (5.7) in J_{x}^{2}, at (t_{0}, x_{0}). Evaluating only transversality in the first component u = 0 of j_{x}^{2}u = (u, u_{x}, u_{xx}) ∈ J_{x}^{2}, we see that

rank (u_{t}, u_{x}) = 2 (5.8)

at (t_{0}, x_{0}). Since P u_{x} = 0 by definition of P, this implies

P u_{t}(t_{0}, x_{0}) ≠ 0 (5.9)

and the proof of corollary 2.2 is complete. ▷◁

6 Numerical Model and Methods

For our numerical simulations, we use two-variable N = 2 reaction-diffusion equations

\begin{matrix}∂_{t}ũ^{1} = △ũ^{1} + f(ũ^{1}, ũ^{2}) \\ ∂_{t}ũ^{2} = D△ũ^{2} + g(ũ^{1}, ũ^{2}) \end{matrix} (6.1)

on a square or cube Ω with Neumann boundary conditions. The functions f(ũ^{1}, ũ^{2}) and g(ũ^{1}, ũ^{2}) express the local reaction kinetics of the two variables ũ^{1} and ũ^{2}. The diffusion coefficient for the ũ^{1} variable has been scaled to unity, and D is the ratio of diffusion coefficients. For the reaction kinetics we use

\begin{matrix}f(ũ^{1}, ũ^{2}) = ϵ^{−1}ũ^{1}(1 − ũ^{1})(ũ^{1} − u_{th}(ũ^{2})) \\ g(ũ^{1}, ũ^{2}) = ũ^{1} − ũ^{2}, \end{matrix} (6.2)

with u_{th}(ũ^{2}) = (ũ^{2} + b)/a. This choice differs from traditional FitzHugh-Nagumo equations, but facilitates fast computer simulations [11]. In nonautonomous simulations, we periodically force the excitability threshold b = b(t) = b_{0} + A cos(ωt). We keep most model parameters fixed at a = 0.8, b_{0} = 0.01, ϵ = 0.02, and D = 0.5. Without forcing, the medium is strongly excitable, see figure 1. See figure 2 for the dynamics of a wave train. In two space dimensions, the equations generate rigidly rotating spirals with small cores. These spirals are far from the meander instability, and appropriate initial conditions quickly converge to rotating waves. We map the coordinates (ũ^{1}, ũ^{2}) into the (u^{1}, u^{2})-coordinates of theorem 2.1 by setting u^{1} = ũ^{1} − 0.5 and u^{2} = ũ^{2} − (a/2 − b_{0}). We have remarked in the introduction, already, that our results are not effected by such a shift of level sets. In the autonomous cases we choose a forcing amplitude A = 0, of course. For collision of spirals in two dimensions, we choose A = 0.01, ω = 3.21. For collision of scroll wave filaments in three dimensions, we choose A = 0.01, ω = 3.92. The challenging aspect of computing wave fronts in excitable media is the resolution of both spatial and temporal details of the wave fronts while the interesting global phenomena occur on a much slower time scale. Since both spatial and temporal resolutions have to be high, the main computational speedup is achieved by minimizing the number of operations necessary per time step and space point. Simulations with cellular automata encounter problems due to grid isotropies [17, 32, 33]. The existence of persistent spatial wave fronts impedes algorithms with variable time steps. Due to linearity of the spatial operator, methods with fixed, small time steps are feasible. Moreover, ũ^{1} and ũ^{2} can be updated in place away from the wave front. We use a third-order semi-implicit stepping routine to time step f, combined with explicit Euler time stepping for g and the Laplacian term. In the evaluation of f and in the diffusion of ũ^{1}, we take into account that ũ^{1} ≈ 0 in a large part of the domain, and that f(0, ũ^{2}) = 0. This allows a cheap update of approximately half of the grid elements and, even with a straightforward finite-difference method, enables simulation on a workstation. The extra effort of an adaptive grid with frequent re-meshing has been avoided. In three space dimensions N = 3, we use a 19-point stencil with good numerical properties (isotropic error, mild time-step constraint) for approximating the Laplacian operator. In two dimensions N = 2, we use the analogous 9-point stencil. Neumann boundary conditions are imposed on all boundaries. For specific simulation runs in this paper, we take 125^{3} grid points. The domain Ω is chosen sufficiently large, in terms of diffusion length, to exhibit scroll wave collision phenomena. The time step △t is chosen close to maximal: let h denote grid size, σ = 3/8 the stability limit of the Laplacian stencil, and choose △t := 0.784σh^{2}. This results in the following numerical parameters: domain Ω = −[15, 15]^{3}, grid spacing h = 30/124 ≈ 1/4, time step △t = 0.0172086, giving △t/ϵ = 0.86043. For high-accuracy studies of the collision of scroll waves, we use a higher resolution of Ω = [−10, 10]^{3}, h = 20/124 ≈ 1/6, △t = 0.00764828, giving △t/ϵ = 0.3882414. Note that △t/ϵ < 1 in both cases, which means that the temporal dynamics are well resolved. Further numerical details for the three-dimensional simulations are given in [11].

7 Filament Visualization

After discretization in the cube domain Ω, and time integration, the solution data u(t, x) ∈ ℝ^{2} are given as values u(t_{i}, x_{i}) at time steps t_{i} and at positions x_{i} on a Cartesian lattice. In our two-dimensional examples, figure 1 and example 8.2, we show the vector field (ũ^{1}, ũ^{2}) = (u^{1} + 0.5, u^{2} + (a/2 − b_{0}), choosing for each point a color vector in RGB space of (u^{1}, 0.73 ∗ (u^{2})^{2}, 1.56 ∗ u^{2}). We also mark the (past) trace of the tip path in white, to keep track of the movements of the spiral tip. In figure 3 and example 8.3, we depict the wave front in x ∈ Ω as the surface u^{1} = 0. To determine the filament location, alias the level set

φ^{t} := {x ∈ Ω | u^{1}(t, x) = u^{2}(t, x) = 0}, (7.1)

we use a simplicial algorithm in the spirit of [2, ch. 12].

As in section 6, let Q ⊆ Ω be any of the small discretization cubes. We triangulate its faces by bisecting diagonals, denoting the resulting closed triangles by τ. The corners of τ are vertices of Q. We orient τ according to the induced orientation of ∂Q by its outward normal ν and the right hand rule applied to (τ, ν). By linear interpolation, u(t, τ) ⊂ ℝ^{2} is also an oriented triangle. The filament φ^{t} passes through τ, on the discretized level, if and only if 0 ∈ u(t, τ). Inverting the linear approximation u on τ defines an approximation φ_{ι}^{t} ∈ τ to φ^{t} ∩ τ. We orient φ^{t} to leave Q through τ, if the orientation of the triangle u(t, τ) is positive (”door out”). In the opposite case of negative orientation we say that φ^{t} enters Q through τ (”door in”). By elementary degree theory, the numbers of in-doors and of out-doors coincide for any small discretization cube Q. Matching in-doors φ_{ι}^{t} and out-doors φ_{ι′}^{t} in pairs defines a piecewise linear, oriented approximation to the filament φ^{t}. For orientations before and after crossover-collision see figure 4. Note that here and below, we freely discard certain degenerate, non-generic situations from our discussion which complicate the presentation and tend to confuse the simple issue. In fact, due to homotopy invariance of Brouwer degree, this piecewise linear (PL) method is robust with respect to perturbations of degeneracies like filaments touching a face of the cube Q or repeatedly threading through the same triangle τ. To indicate the phase near the filament φ_{t}, we compute a tangential approximation to the accompanying somewhat arbitrary isochrone

χ^{t} := {x ∈ Ω | u^{1}(t, x) ≥ 0 = u^{2}(t, x)} (7.2)

as follows. The values (u^{1}, u^{2})(t, x) = (α, 0) with α > 0 define a local half line in the face triangle x ∈ τ through the filament point φ_{ι}^{t} ∈ τ. Together with a filament point φ_{ι−1}^{t} in another cube face, this half line also defines a half space which approximates the isochrone χ^{t}, locally . We choose a point φ̃_{ι}^{t} in this half space, a fixed distance from φ_{ι}^{t} and such that the line from φ_{ι}^{t} to φ̃_{ι}^{t} is orthogonal to the filament line from φ_{ι−1}^{t} to φ_{ι}^{t}. The sequence of triangles (φ_{ι−1}^{t}, φ̃_{ι−1}^{t}, φ̃_{ι}^{t}), (φ_{ι−1}^{t}, φ̃_{ι}^{t}, φ_{ι}^{t}) then define a triangulated isochrone band approximating χ^{t} near the filament φ^{t}. In practical computations shown in the next section, we distinguish an absolute front and back of the isochrone band by color, independently of camera angle and position. This difference reflects the absolute orientation of filaments, introduced above, which induces an absolute orientation and an absolute normal for the accompanying isochrone χ^{t}. The absolute normal of the isochrone χ^{t} also points into the propagation direction of the isochrone, by our choice of orientation.

8 Examples

In this section we present four simulations of three-dimensional filament dynamics, both in autonomous and in periodically forced cases. All examples are based on equations (6.1) with the set of nonlinearities and parameters specified there. We use a cube Ω = [−15, 15]^{3} as a spatial domain, together with Neumann boundary conditions. Only in example 8.4, we use a smaller cube Ω = [−10, 10]^{3}. Reflecting the solutions through the boundaries we obtain an extension to the larger cube 2Ω with periodic boundary conditions. Viewing this system on the flat 3-torus T^{3}, equivalently, eliminates all boundary conditions and avoids the issue of ∂Ω not being smooth. In the paper version, each of the spatio-temporal singularities at (t_{0}, x_{0}) is illustrated by a series of still shots: t ⪆ 0, t ⪅ t_{0}, t = t_{0}, t ⪆ t_{0} and t = tend for the respective run. In the Internet version, each sequence is replaced by a downloadable animation in MPEG-1 format; see

              http://www.math.fu-berlin.de/~Dynamik/

For possible later, updated and revised versions, please contact the authors. Discretization was performed by 125^{3} cubes and a time step of △t = 0.0172086 (△t = 0.00764828 in example 8.4); see section 6. Autonomous cases refer to the forcing amplitude A = 0, whereas A = 0.01 switches on non-autonomous additive forcing.

8.1 Initial Conditions

Prescribing approximate initial conditions for colliding scroll waves in three space dimensions is a somewhat delicate issue. We describe the construction in 8.1.1, 8.1.2 below. We discuss our four examples in sections 8.3-8.6.

8.1.1 Two-dimensional spirals

According to our numerical simulations, planar spiral waves are very robust objects. In fact, sufficiently separated nondegenerate zeroes of the planar “vector field” (u_{0}^{1}, u_{0}^{2})(x_{1}, x_{2}) of initial conditions typically seemed to converge into collections of single-armed spiral waves. Their tips were located nearby the prescribed zeroes of u_{0}. To prepare for our construction of scroll waves below, we nevertheless construct u_{0} as a composition of two maps,

(8.1)  u_{0} = σ ∘ γ
(8.2) γ: ℝ^{2} ⊇ Ω → ℂ
(x_{1}, x_{2}) ↦ z
(8.3) σ: ℂ → ℝ^{2}
z ↦ (u_{0}^{1}, u_{0}^{2})

Here γ prescribes the geometric location of the spiral tip and wave fronts. The scaling map σ is chosen piecewise linear. It adjusts for the appropriate range of u-values to trace out a wave front cycle in our excitable medium, see fig. 2. Specifically, we choose

σ(z) = (u^{1}, u^{2}) = (Re z, Im z/4) (8.4) near the origin. Further away, we cut off by constants as follows:

u^{1} := \begin{cases}−0.5, & when Re(z) < −0.5 \\ Re(z), & when Re(z) ∈ [−0.5, 0.5] \\ 0.5, & when Re(z) > 0.5 \end{cases}

u^{2} := \begin{cases}−0.4, & when Im(z) < −1.6 \\ 0.25Im(z), & when Im(z) ∈ [−1.6, 1.6] \\ 0.4, & when Im(z) > 1.6 \end{cases}

(8.5)

In the following, we will sometimes further decompose σ = σ_{2} ∘ σ_{1} where

σ_{1}(z) = (Re(z), Im(z)/4) (8.6)

is linear and the clamping σ_{2} : ℝ^{2} → ℝ^{2} is the cut-off

(8.7) (u^{1}, u^{2}) ↦ (sign(u^{1}) min{|u^{1}|, 0.5}, sign(u^{2}) min{|u^{2}|, 0.4}).

For example, this choice of σ, combined with the simplest geometry map γ(x_{1}, x_{2}) = x_{1} + ix_{2}, results in a spiral wave rotating clockwise around the origin, with wave front at x_{1} = 0, x_{2} < 0, initially, and wave back at x_{1} = 0, x_{2} > 0. A possible initial condition for a spiral — antispiral pair as in example 8.2 below would be

\begin{matrix}γ: [−15, 15]^{2} & → & ℂ \\ (x_{1}, x_{2}) & ↦ & |x_{1}| − 6 + ix_{2}. \end{matrix}

This reflection symmetric initial condition creates a pair of spirals rotating around (±6, 0). The spiral at (6, 0) rotates clockwise and the symmetric spiral around (−6, 0) rotates anti-clockwise.

8.1.2 Three-dimensional scrolls

It is useful to visualize a three-dimensional scroll wave as a stack foliated by twodimensional slices which contain planar spirals. Initial conditions u_{0} = σ ∘γ for scroll waves then contain the following ingredients: a mapping γ : ℝ^{3} → ℂ that stacks the spirals into the desired three-dimensional geometry, and a scaling σ : ℂ → ℝ^{2}. For planar γ: ℝ^{2} → ℂ as in (8.2), the scaling σ of (8.3)– (8.7) generates a spiral whose tip is at the origin in ℝ^{2}. For γ: ℝ^{3} → ℂ, the preimage in ℝ^{3} of the origin under the stacking map γ will therefore comprise the filament of the three-dimensional scroll wave. For example, it is easy to find a stacking map γ that gives rise to a single straight scroll wave with vertical filament: γ(x_{1}, x_{2}, x_{3}) := x_{1} + ix_{2}. As soon as filaments are required to form rings, linked rings or knots, however, the design of stacking maps γ with the appropriate zero set becomes more difficult.

For the generation of more complicated stacking maps γ, we largely follow the method pioneered by Winfree et al [37, 18, 39]. This approach uses a standard method of embedding an algebraic knot in 3-space [6]. For convenience of our readers, we briefly recall the construction here. We construct stacking maps γ: ℝ^{3} → ℂ with prescribed, possibly linked or knotted zero set as a composition

γ = p ∘ s. (8.8)

Here the embedding s: ℝ^{3} → ℂ^{2} will be related to the map

s̃: ℝ^{3} → 𝕊_{ε}^{3} ⊂ ℝ^{4} = ℂ^{2} (8.9)

denoting the inverse of the standard stereographic projection from the standard 3- sphere 𝕊_{ε}^{3} of radius ε to ℝ^{3}; see (8.11) below. The map

p: ℂ^{2} → ℂ (8.10)

is a complex polynomial p = p(z_{1}, z_{2}) in two complex variables z_{1}, z_{2}. The zero set of p describes a real, two-dimensional variety V in ℂ^{2}. Consider the intersection φ̃ of V with the small 3-sphere 𝕊_{ε}^{3}, that is φ̃ := V ∩ 𝕊_{ε}^{3}. Typically, φ := s^{−1}(φ̃) ⊂ ℝ^{3}, the zero set of γ, will be a one-dimensional curve or a collection of curves: the desired filament of our scroll wave. In the simplest case φ may be a circle embedded into the 3-sphere 𝕊_{ε}^{3}. If however zero is a critical point of the polynomial p, then the filament φ need not be a topological circle. And even if φ̃ happens to be a topological circle, it may be embedded as a knot in 𝕊_{ε}^{3}. The inverse stereographic map s̃ is given explicitly by

s̃(x_{1}, x_{2}, x_{3}) = \frac{1}{R^{2} + ε^{2}} (\begin{matrix}2ε^{2}x_{1} \\ 2ε^{2}x_{2} \\ 2ε^{2}x_{3} \\ (R^{2} − ε^{2})ε \end{matrix}) ≅ \frac{2ε^{2}}{R^{2} + ε^{2}} (\begin{matrix}x_{1} + ix_{2} \\ x_{3} + i\frac{(R^{2}−ε^{2})}{(2ε)} \end{matrix}) (8.11)

where R^{2} ≡ x_{1}^{2} + x_{2}^{2} + x_{3}^{2}. Note that points inside 𝕊_{ε}^{2} ⊂ ℝ^{3} are mapped to the lower hemisphere, points outside 𝕊_{ε}^{2} to the upper hemisphere of 𝕊_{ε}^{3}. In our construction (8.8) of the stacking map γ, we now replace the inverse stereographic map s̃ by the embedding

s(x_{1}, x_{2}, x_{3}) ≅ c (\begin{matrix}x_{1} \\ x_{2} \\ x_{3} \\ (cR^{2} − \frac{1}{4c}) \end{matrix}) ≅ (\begin{matrix}cx_{1} + icx_{2} \\ cx_{3} + i(c^{2}R^{2} − \frac{1}{4}) \end{matrix}) (8.12)

with a suitable scaling factor c. Clearly x → ∞ in ℝ^{3} implies s(x) → ∞ in ℂ^{2}. In the examples 8.5, 8.6 of a pair of linked rings and of a torus knot below, the filaments φ = γ^{−1}(0) ∩ Ω, φ̃ = s(φ) = p^{−1}(0) ∩ s(Ω) do not intersect the compact boundaries of the cube ∂Ω, s(∂Ω), respectively. Therefore, the embedded paraboloid s(ℝ^{3}) can in fact be modified outside s(Ω) without changing the filaments in Ω. We modify s such that clos s(ℝ^{3}) closes up to a diffeomorphically embedded 3-shere S diffeotopic to 𝕊^{3} in ℂ^{2} \ {0}, by a family s_{ϑ} of embeddings 0 ≤ ϑ ≤ 1. Moreover, we will choose p = p(z_{1}, z_{2}) such that z_{1} = z_{2} = 0 is the only critical point of p in ℂ^{2}. If the embedding s_{ϑ}(ℝ^{3}) remains transverse to p^{−1}(0) in ℂ^{2} \ {0} throughout the diffeotopy, then the variety p^{−1}(0) is an embedded real surface in ℂ^{2}, outside z = 0. The filament φ̃ = s(φ) = p^{−1}(0) ∩ s(Ω) is diffeotopic to some components of p^{−1}(0) ∩ 𝕊_{ε}^{3}, which in turn are described classically in algebraic geometry. The same remarks apply, slightly more generally, if we replace s by a composition

s ∘ ℓ (8.13)

where ℓ denotes a nondegenerate affine transformation in ℝ^{3}. In summary, we generate our initial conditions by applying the following composition of mappings:

u_{0} = σ ∘ γ = (σ_{2} ∘ σ_{1}) ∘ (p ∘ s). (8.14)

Here the scaling σ is given by (8.5)–(8.7). The modified stereographic projection s is given by (8.12) with ℓ = id, except in example 8.5, and with appropriate scaling constant c. The polynomial p is chosen according to the desired topology of the filament. The initial conditions thus created do not necessarily respect the boundary conditions; however any intersection of a filament with the boundary is transverse. Anyways, such intersections only occur in example 8.4. Neumann boundary conditions can be enforced artificially, by standard implementation, without introducing additional filaments.

8.2 Two-dimensional spiral pair annihilation

As a preparation to visualizing the three-dimensional behavior, we begin with the collision of a pair of counter-rotating planar spirals. We use a domain Ω = [−15, 15]^{2} and discretize with 125^{2} grid points, resulting in the same spatial and temporal resolution as with our three-dimensional experiments. In the movie and pictures, we show the subdomain [−15, 15] × [−11.25, 11.25] to get the 3:4 size ratio typical for video. For initial conditions, we take the fully developed rigidly rotating spiral of figure 1 with origin at (−6, 0), for the half-plane x_{1} ≤ 0, and reflect at the vertical x_{2}-axis. Near-resonant periodic forcing with an amplitude A = 0.01 and ω = 3.21 causes the spirals to drift towards each other until they collide. The forcing makes the spiral tips drift on an almost straight, epicyclic trajectory, until they reach interaction distance at time t = 19.2. The paths of the tips show that the forcing is strong enough to move the spirals by approximately twice their tip radius per rotation (which is small in comparison to their wave 718 Bernold Fiedler and Rolf M. Mantel

t = 0  t = 19.60

t = 31.77  t = 35.08

t = 37.67  t = 38.46

Figure 7: Interaction and collision of a pair of spiral waves in the plane.

MPEG-Movie [26.4MB,gzipped]

http://www.mathematik.uni-bielefeld.de/documenta/vol-05/21.mpeg/twodim.mpg.gz

length). During the interaction time of the spiral tips, the u^{2} gradients are much shallower than at other times. This can be seen by the fact that the bright red part of the wave front is further away from the tip location. The spirals then wander along the vertical axis, the excited center getting smaller with every revolution. Finally the center is too small to sustain excitation (t_{0} = 39.825) and disappears; the spirals annihilate. The purely local interaction between the spiral tips shortly before collision from time t^{∗} = 19.2 up to the extinction at t_{0} = 39.825, x_{0} = (0, 1.4) is clearly visible from the tip paths. In view of theorem 2.1, this annihilation illustrates the left saddle-node singularity of fig. 5 for dim u = dim x = 2.

8.3 Scroll ring annihilation

Our first three-dimensional example shows the disappearance of a closed circular filament as described, from an abstract singularity theory point of view, in theorem 2.2,(ii), and as illustrated in figure 6. The example is autonomous, A = 0. Viewed in a vertical planar slice through the center, the dynamics is reminiscent of the two-dimensional spiral pair annihilation 8.2. Instead of periodic forcing, this time, the curvature of the three-dimensional filament seems to be responsible for the filament contraction and annihilation [24]. The simplest initial conditions to create a scroll ring would be via the polynomial p(z_{1}, z_{2}) = z_{2}, resulting in the vertical axis s̃(Rez_{2}) = x_{3} = 0 being a symmetry axis both for u_{0} and for Ω ⊂ ℝ^{3}. In order for the initial conditions to be less symmetric with respect to the boundaries of the domain Ω, we apply the translation ℓx = x − x_{∗} with x_{∗} = (−1.5, 3, 0), and we choose a polynomial p that also depends on z_{1}. Our initial conditions are prescribed by (8.14), using

\begin{aligned}p & = & z_{2} + 0.1 iz_{1}, \\ c & = & 8/21. \end{aligned} (8.15)

Under discretization, scroll ring annihilation occurs at

t_{0} = 9.10; x_{0} = (−1.5, 3.5, −0.5). (8.16)

For illustration/animation see figs. 8.

8.4 Crossover collision of scroll waves

We now return to the motivating phenomenon of this paper, outlined in the introduction; see (1.6) and figure 4. For finer spatial resolution, we choose a smaller domain, Ω = [−10, 10]^{3}, with discretization into 125^{3} cubes. Due to the finer space discretization of 20/124 instead of 30/124, we choose a smaller time step of △t = 0.00764828. The example is non-autonomous, with forcing amplitude A = 0.01 and frequency ω = 3.92. Circumventing the polynomial construction γ = p ∘ s, we take

γ(x̃/c) = ((x_{3} + π/6) + i sin(x_{1}))(sin(x_{2}) − i(x_{3} − π/6)), (8.17)

720 Bernold Fiedler and Rolf M. Mantel

t = 4.16  t = 4.58

t = 15.32  t = 15.67

t = 15.74  t = 16.38

Figure 8: Scroll ring annihilation. By t = 4.2, a spiral-like cross-section has formed. The scroll ring emits ball shaped target waves twice per revolution, starting at approximately t = 4.58. After scroll ring annihilation at t_{0} = 23.45, the surface u^{1} = 0 largely follows a concentric target wave pattern rather than a scroll ring pattern. The remaining target waves move outwards, and the medium becomes quiescent.

MPEG-Movie [10.7MB,gzipped]

http://www.mathematik.uni-bielefeld.de/documenta/vol-05/21.mpeg/ring.mpg.gz

Documenta Mathematica 5 (2000) 695–731

Crossover Collision of Scroll Wave Filaments 721

t = 6.13  t = 17.875

t = 31.60  t = 35.654

t = 36.446  t = 39.918

Figure 9: Collision of scroll waves: Two scroll wave filaments drift towards each other. After t = 17, they start interacting visibly. Around t = 34, the filaments have found a common tangent plane and start lining up for collision. The crossover collision occurs at t_{0} = 35.83, x_{0} = (−3.25, 3.25, 0). After collision, the filaments connect adjacent faces of the cube rather than opposite faces.

MPEG-Movie [10.6MB,gzipped]

http://www.mathematik.uni-bielefeld.de/documenta/vol-05/21.mpeg/scrolls.mpg.gz

which has zeros in [−π/2, π/2]^{3} at (0, x_{2}, −π/6) and at (x_{1}, 0, π/6). Taking u_{0} = σ ∘ γ, we then start with explicit initial conditions

u_{0}^{1}(x) = sin(ax_{1})(ax_{3} − π/6) + sin(ax_{2})(ax_{3} + π/6),

u_{0}^{2}(x) = 0.25 ∗ (sin(ax_{1}) sin(ax_{2}) − (ax_{3} + π/6)(ax_{3} − π/6)). (8.18)

The spatial scaling factor a is chosen as π/20. This example was selected because (8.17) has zeroes in [−π/2, π/2]^{3} at (0, x_{2}, −π/6) and at (x_{1}, 0, π/6). Then (u_{0}^{1}, u_{0}^{2}) has zeroes at (0, x̃_{2}, −10/3) and at (x̃_{1}, 0, 10/3). Therefore, at t = 0, filaments are at right angles to each other. Near resonant forcing with amplitude A = 0.01 and frequency ω = 3.92 is chosen, together with an appropriate initial phase, such that the filaments drift towards each other and eventually interact. Under discretization, crossover collision occurs at

t_{0} = 35.83; x_{0} = (−3.25, 3.25, 0). (8.19)

For illustration/animation see figures 9 and 10. Naively, there would be at least two options for non-destructive collision of the two scroll wave filaments. In figures 4, 9 and 10, the two primary filaments are seen to touch, forming a crossing with four emanating semi-branches. Keeping their orientation, the semi-branches could either simply re-connect, as before the collision. Alternatively, they could separate and connect with that semibranch of matching orientation which they were not attached to previously. The first scenario of a crossing collision may be more intuitive at first: the two incoming semi-branches simply reconnect to their previous outgoing partners without exchanging their pairing. Such a crossing clearly would not change the global connectivity of the filaments. Viewed in projection onto the tangent plane E at collision time t_{0}, however, the filament branches would then have to remain crossing immediately before and after collision time t_{0}, in contradiction to both theorem 2.2 and numerical observation in figures 9 and 10. Note that the filaments, albeit initially straight lines, have to bend out of their way considerably in order to accommodate a generic crossover collision in the tangent plane E. Indeed, initial conditions, periodic forcing, and boundary conditions are all chosen invariant under a rotation by 180^{o} around the axis A which diagonally connects the mid-edge points (−10, 10, 0) and (10, −10, 0) of the domain Ω. This rotation invariance is preserved by the solution u(t, .). Because rotation initially maps one filament into the other, the collision point x_{0} must occur on the axis A – and it does, see (8.19). Similarly, the tangent plane E must be orthogonal to A, forming angles of 45^{o} with the straight line initial conditions. We found it fascinating to watch the numerical filaments obey all these predictions. We caution the reader here that theorems 2.1 and 2.2, as they stand, do not directly apply within restricted classes of symmetric initial conditions. In full generality, the necessary modifications require a restriction to, and analysis of, invariant singularities and their codimensions in spaces of symmetry invariant

t = 35.764  t = 35.918

tangent plane

side view

top view

Figure 10: Details of the crossover collision: breaking and reconnecting scroll wave filaments, consistently with theorem 2.2. The two incoming semi-branches exchange their pairing with the two outgoing semi-branches at t = t_{0}, x = x_{0}. Each incoming semi-branch crosses over to its opposite outgoing semi-branch. The projected branches, when viewed locally in the tangent plane E = ker u_{x} to the collision configuration at t = t_{0}, x = x_{0}, neither cross before nor after collision.

MPEG-Movie [10.7MB,gzipped]

http://www.mathematik.uni-bielefeld.de/documenta/vol-05/21.mpeg/scrollzoom.mpg.gz

Documenta Mathematica 5 (2000) 695–731

k-jets, again based on transversality, lemma 3.1. The present example and its codimension, however, comply with our simple rotation symmetry. In the coordinates (1.6) of crossover collision this can be seen from invariance under the 180^{o} rotation (x_{1}, x_{2}) ↦ (−x_{1}, −x_{2}) around the x_{3} axis.

8.5 Collision of linked twisted scroll rings

In the previous non-autonomous example we have seen how crossover collisions change the local connectivity of tip filaments. We now present an autonomous example, with forcing amplitude A = 0, where two linked filaments merge into a single filament. After collision the resulting filament is neither knotted nor self-linked but is isotopic to a circle. We start with initial conditions u_{0} prescribed by (8.14), with the polynomial p = z_{1}^{2} − z_{2}^{2} and stereographic scaling factor c = 8/21 in (8.12).

\begin{aligned}p & = & z_{1}^{2} − z_{2}^{2}, \\ c & = & 8/21. \end{aligned} (8.20)

Under discretization, crossover collision occurs at

t_{0} = 4.90, x_{0} = (0, 0, −2.14) (8.21)

For illustration/animation see fig. 11. We comment on the changes of the global topological characteristics of twist and linking which occur at the crossover collision in this example. See figure 12 for a caricature of the essential features. To determine the twist of a non self-intersecting closed oriented filament φ^{t}, we first orient the tip filament φ^{t} as described in section 7. Then we count the integer winding number of the accompanying isochrone band χ^{t} around φ^{t}, according to the right hand rule. The integer twist can be positive, negative, or zero. Next suppose the filament φ^{t} spans an embedded disk, as all filaments in figures 11, 12 do. The orientation of φ^{t} induces an orientation of the disk which, again by the right hand rule, we can represent by a field of vectors ν normal to the disk. To any other oriented filament crossing the disk transversely, we associate a crossing sign +1, if the crossing is in the direction of ν, and −1 otherwise. Following [40], the sum of crossing signs on the disk adds up to the twist of the boundary filament φ^{t}. Applied to the schematic representation of figure 11 in figure 12, we conclude that the two filaments φ_{ℓ}^{t}, φ_{r}^{t} for t < t_{0} each have twist −1. After collision the single remaining filament is untwisted. Our example therefore indicates that one can hope, at best, for a conservation of the parity of the total twist.

8.6 Unknotting the trefoil knot by crossover collision

In the previous example two linked but unknotted filaments merged into a single filament. Also, the initial conditions were far from a long-term solution

t = 0.00  t = 1.665

t = 4.278  4.881

t = 4.910  t = 10.598

Figure 11: Crossover collision of two linked twisted filaments at t_{0} = 4.90, x_{0} = (0, 0, −2.14) into a single untwisted filament.

MPEG-Movie [3.3MB,gzipped]

http://www.mathematik.uni-bielefeld.de/documenta/vol-05/21.mpeg/rings.mpg.gz

φ_{ℓ}^{t}  φ_{r}^{t}

Figure 12: A caricature of the crossover collision of two linked, simply twisted filaments φ_{ℓ}^{t}, φ_{r}^{t} at t = t_{0}. Before collision, each scroll ring possesses a twist of −1. After collision, the resulting scroll ring is untwisted, globally.

of the equations. In contrast, we now take a trefoil knot as an initial condition that already exhibits fully developed scroll waves. We then rescale space, which is equivalent to a change of diffusion constants. This brings the filaments into sufficiently close contact for interaction. The initial conditions for this autonomous example, A = 0, are the numerical end state of a coarser simulation on a domain Ω_{1} = [−25, 25]^{3}, also running on a numerical grid of 125^{3} grid points. The initial condition for the coarser simulation (starting at time t = −10) is created using the polynomial p = z_{1}^{2}−z_{2}^{3} with stereographic scaling factor c = 1/5 in (8.12):

z_{1} = 1/5(x_{1} + ix_{2});
z_{2} = 1/5x_{3} + i((x_{1}^{2} + x_{2}^{2} + x_{3}^{2})/5^{2} − 1/4);
u_{0}^{1}(x) = Re(z_{1}^{2} − z_{2}^{3})  clamped by (8.7);
u_{0}^{2}(x) = 0.25Im(z_{1}^{2} − z_{2}^{3})  clamped by (8.7).

(8.22)

At time t = 0, we stop the simulation, keeping the same numerical data at grid points but rescaling the domain to Ω = [−15, 15]^{3}. This is the initial condition at t = 0. Under discretization, crossover collision from a trefoil knot to two linked rings is observed at

t_{0} = 8.94, x_{0} = (0, 0, −9.28) (8.23)

For illustration/animation see figure 13. Again we provide a caricature in figure 14.

8.7 Discussion of examples

We conclude our series of examples with some remarks. Concerning example 8.3 of scroll ring annihilation we observe that only untwisted scroll rings can be directly annihilated. This follows from the normal form of the corresponding singularity with positive definite quadratic form ⟨P u_{t}, P u_{xx}|_{E}⟩ at (t_{0}, x_{0}); see section 1. More globally, it also follows from the observation that the shrinking disk spanned by a circular filament near annihilation is not traversed by other filaments. Indeed a filament shrinking around another, large filament would

t = 0.00  t = 4.905

t = 7.985  t = 8.914

t = 8.966  t = 12.597

Figure 13: Decomposing the trefoil knot into two linked twisted unknotted filaments by crossover collision at t_{0} = 8.94, x_{0} = (0, 0, −9.28). As explained in example 8.3, we see in figures 13, 14 how the trefoil knot with twist ±3 decomposes into two unknotted, but mutually linked twisted filaments, each of twist −1.

MPEG-Movie [8.9MB,gzipped]

http://www.mathematik.uni-bielefeld.de/documenta/vol-05/21.mpeg/knot.mpg.gz

Figure 14: A caricature of the unknotting of the trefoil knot, showing the orientations of all filaments.

require a three-dimensional kernel and hence a singularity of codimension at least six. We have not presented an example for the process opposing annihilation: the creation of a circular filament by a negative definite quadratic form ⟨P u_{t}, P u_{xx}|_{E}⟩ at (t_{0}, x_{0}). Since the definiteness required for P u_{xx}|_{E} does not predetermine the direction of △u, we could construct initial conditions u_{0}(x) corresponding to scroll ring creation at t_{0} = 0, x_{0} ∈ Ω. Although we expect scroll ring creation to be feasible also for large positive times t_{0}, we did not observe this phenomenon in our simulations so far. Our results provide specific examples of the “internal” collision type, which [31] have described as topologically viable; furthermore, we show that crossover collision is the only generic way for scroll waves to change their topological linking type. From a modeling point of view, experimental systems may require substantially more than just two dependent variables u^{1}, u^{2} for an adequate description by parabolic reaction diffusion systems. We repeat that theorem 4.2 predicts the described two-variable phenomena to occur in any projection setting, where only two combinations of the relevant quantities u^{1}, ..., u^{m} are observable, for example by color shading. We emphasize that this observation neither requires, nor corresponds to, a dynamic reduction of the full underlying reaction diffusion system by inertial manifolds or related techniques of dimension reduction. Aiming at the ubiquitous wealth of phenomena of pattern formation and pattern transformation, our paper has detected and addressed just a few elementary dynamic effects peculiar to systems of two equations in three space dimensions. Clearly, the theoretical framework supports significantly more complicated spatio-temporal effects than were presented here. Applicability hopefully also will reach far beyond the specific motivating context of Belousov-Zhabotinsky patterns or excitable media.

References

[1] R. Abraham and J. Robbin. Transversal Mappings and Flows. Benjamin, New York, 1967.

[2] E. L. Allgower and K. Georg. Numerical Continuation Methods. An Introduction. Springer-Verlag, Berlin, 1990.

[3] S. Angenent. The zero set of a solution of a parabolic equation Crelle J. reine angew. Math., 390:79–96, 1988.

[4] D. Barkley. Spiral meandering. In Chemical Waves and Patterns, R. Kapral and K. Showalter (eds.), 163–190. Kluwer, 1995.

[5] D. Barkley, M. Kness, and L. S. Tuckerman. Spiral-wave dynamics in a simple model of excitable media - the transition from simple to compound rotation. Phys. Rev. A, 42:2489–2492, 1990.

[6] K. Brauner. Zur Geometrie der Funktionen zweier komplexer Veränderlicher iii,iv. Abh. Math. Sem. Hamburg, 6:8–54, 1928.

[7] J. Damon. Local Morse theory for solutions to the heat equation and Gaussian blurring. J. Diff. Eqns., 115:368–401, 1995.

[8] J. Damon. Generic properties of solutions to partial differential equations. Arch. Rat. Mech. Analysis, 140:353–403, 1997.

[9] J. Damon. Singularities with scale threshold and discrete functions exhibiting generic properties. Mat. Contemp., 12:45–65, 1997.

[10] J. Dieudonné. Elements d’Analyse 2. Gauthier-Villars, Paris, 1969.

[11] M. Dowle, R. M. Mantel, and D. Barkley. Fast simulations of waves in three-dimensional excitable media. Int. J. Bifur. Chaos, 7:2529–2546, 1997.

[12] B. Fiedler and C. Rocha. Orbit equivalence of global attractors of semilinear parabolic differential equations. Trans. Amer. Math. Soc., 352:257-284, 2000.

[13] B. Fiedler and D. Turaev. Normal forms, resonances. and meandering tip motions near relative equilibria of Euclidean group actions. Arch. Rat. Mech. Analysis, 145:129–159, 1998.

[14] A. Friedman. Partial Differential Equations. Holt, Reinhart & Winston, New York, 1969.

[15] M. Golubitsky and D. G. Schaeffer. Singularities and Groups in Bifurcation Theory I. Appl. Math. Sc. 51, Springer-Verlag, 1985.

[16] D. Henry. Geometric Theory of Semilinear Parabolic Equations. Springer-Verlag, Berlin, 1981.

[17] C. Henze and J. J. Tyson. Cellular-automaton model of 3-dimensional excitable media. J. Chem. Soc. Farad. Transac., 92:2883–2895, 1996.

[18] C. Henze and A.T. Winfree. A stable knotted singularity in an excitable medium. Int. J. Bif. Chaos, 1:891–922, 1991.

[19] M. W. Hirsch. Differential Topology. Springer-Verlag, New York, 1976.

[20] W. Jahnke, C. Henze, and A. T. Winfree. Chemical vortex dynamics in 3-dimensional excitable media. Nature, 336:662–665, 1988.

[21] O. A. Ladyzhenskaja, V. A. Solonnikov, and N. N. Ural´ ceva. Linear and Quasilinear Equations of Parabolic Type. Moscow, 1967. English: AMS, Providence, 1968.

[22] R. M. Mantel and D. Barkley. Periodic forcing of scroll-wave patterns in three-dimensional exciatable media. Preprint, University of Minnesota, IMA #1568, 1998.

[23] H. Matano. Nonincrease of the lap-number of a solution for a onedimensional semi-linear parabolic equation. J. Fac. Sci. Univ. Tokyo Sec. IA, 29:401–441, 1982.

[24] A. V. Panfilov, A. N. Rudenko, and V. I. Krinsky. Turbulent rings in 3-dimensional active media with diffusion by 2 components. Biofizika, 31:850–854, 1986.

[25] A. V. Panfilov and A. T. Winfree. Dynamical simulations of twisted scroll rings in 3-dimensional excitable media. Physica D, 17:323–330, 1985.

[26] A. Pazy. Semigroups of Linear Operators and Applications to Partial Differential Equations. Springer-Verlag, New York, 1983.

[27] D.G. Schaeffer. A regularity theorem for conservation laws. Adv. Math., 11:368–386, 1973.

[28] C. Sturm. Sur une classe d’equations a difference partielles. J. Math. Pure Appl., 1:373–444, 1836.

[29] H. Tanabe. Equations of Evolution. Pitman, London, 1979.

[30] R. Thom. Un lemma sur les applications differentiables. Bol. Soc. Mat. Mex., 2:59–71, 1956.

[31] J. J. Tyson and S. H. Strogatz. The differential geometry of scroll waves. Int. J. Bif. Chaos, 1:723–744, 1991.

[32] J. R. Weimar, J. J. Tyson, and L. T. Watson. 3rd generation cellular automaton for modeling excitable media. Physica D, 55:328–339, 1992.

[33] J. R. Weimar, J. J. Tyson, and L. T. Watson. Diffusion and wavepropagation in cellular automaton models of excitable media. Physica D, 55:309–327, 1992.

[34] N. Wiener, A. Rosenblueth. The mathematical formulation of the problem of conduction of impulses in a network of connected excitable elements, specifically in cardiac muscle. Arch. Inst. Cardiol. Mexico, 16:205-265, 1946.

[35] A. T. Winfree. Spiral waves of chemical activity. Science, 175:634–636, 1972.

[36] A. T. Winfree. Scroll-shaped waves of chemical activity in three dimensions. Science, 181:937–939, 1973.

[37] A. T. Winfree. When Time Breaks Down. Princeton Univ. Press, Princeton, NJ, 1987.

[38] A. T. Winfree. Varieties of spiral wave behavior: An experimentalist’s approach to the theory of excitable media. Chaos, 1:303–334, 1991.

[39] A. T. Winfree. Persistent tangles of vortex rings in excitable media. Physica D, 84:126–147, 1995.

[40] A. T. Winfree, E. M. Winfree, and M. Seifert. Organizing centers in a cellular excitable medium. Physica D, 17:109–115, 1985.

Bernold Fiedler Institut für Mathematik I Freie Universität Berlin Arnimallee 2-6 D-14195 Berlin, Germany fiedler@math.fu-berlin.de Rolf M. Mantel Philipps-Universität Marburg Fachbereich Mathematik und Informatik Hans-Meerwein Str. 35032 Marburg mantel@mathematik.uni-marburg.de 732

Documenta Mathematica 5 (2000)