Äîêóìåíò âçÿò èç êýøà ïîèñêîâîé ìàøèíû. Àäðåñ îðèãèíàëüíîãî äîêóìåíòà : http://chaos.phys.msu.ru/loskutov/PDF/BPHS631.pdf Äàòà èçìåíåíèÿ: Mon Jan 31 16:01:39 2011 Äàòà èíäåêñèðîâàíèÿ: Mon Oct 1 19:49:33 2012 Êîäèðîâêà:

ISSN 0006 3509, Biophysics, 2009, Vol. 54, No. 5, pp. 631­636. © Pleiades Publishing, Inc., 2009. Original Russian Text © E.P. Zemskov, A.Yu. Loskutov, 2009, published in Biofizika, 2009, Vol. 54, No. 5, pp. 908­915.

COMPLEX SYSTEMS BIOPHYSICS

Traveling Waves in a Piecewise Linear Reaction­Diffusion Model of Excitable Medium
E. P. Zemskova and A. Yu. Loskutovb
a

Computing Center, Russian Academy of Sciences, Moscow, 119333 Russia b Physical Faculty, Moscow State University, Moscow, 119992 Russia e mail: zemskov@ccas.ru, loskutov@chaos.phys.msu.ru
Received May 14, 2008

Abstract--One dimensional autowaves (traveling waves) in excitable medium described by a piecewise lin ear reaction­diffusion system have been investigated. Two main types of wave have been considered: a single pulse and a periodic sequence of pulses (wave trains). In a two component system, oscillations are due to the second component of the reaction­diffusion system, while in a one component system, they are caused by external periodic excitation (forcing). Using semianalytical solutions for the wave profile, the shape and velocity of autowaves have been found. It is shown that the dispersion relation for oscillating sequences of pulses has an anomalous character. Key words: excitable medium, autowaves, piecewise linear reaction­diffusion system DOI: 10.1134/S0006350909050145

INTRODUCTION A common problem in mathematical biophysics and in biophysics of ecological communities [1] is modeling of the ecosystem behavior. Therewith arises a task of describing the spatiotemporal dynamics of distributed objects consisting of many interacting ele ments. To solve this task, use is broadly made of the formalism of autowave processes [2], usually under stood as self sustaining processes (including stationary structures) retaining their period, velocity, and shape in an active nonlinear medium. An important role here is played by distributed systems admitting emer gence and interaction of various ordered spatial inho mogeneities and stable traveling waves. Active media are of three types: bistable, excitable, and oscillatory. The main type of structure character istic of bistable media is the switching wave (front). A single pulse is characteristic of an excitable medium; and a periodic pulse sequence (train), of an oscillatory medium [3]. All these structures are in essence travel ing waves and are mathematically described with cor responding nonlinear equations in partial derivatives of the reaction­diffusion type. A classical model of active medium based on the reaction­diffusion system includes two equations, with a cubic nonlinearity in the first one and a linear reaction term in the second one. Depending on the parameter values, all three regimes are realized in such a system. The equations of this system were proposed by FitzHugh [4] and Nagumo et al. [5] and were ini tially used as a model for excitation wave propagation in nerve tissue. All three types of nonlinear traveling

waves (fronts, pulses, and trains) have been intensely studied. It has turned out that along with waves of standard shape there can exist waves with profile oscil lations [6, 7]. In the latter case, the pulse train has an anomalous dispersion relation (wave velocity vs. period) [7, 8]. Autowave structures may form because of sponta neous emergence of traveling pulse sources as well as stationary dissipative structures under conditions when autowaves do not collide but move away from each other [9]. The developing structure, depending on conditions, can be in states characterized by the presence of one as well as several autowave processes. Thus in strongly excitable media, i.e., those with a low excitation threshold, there is a possibility of spontane ous appearance of a leading center and periodical birth of sources at the points of autowave decay [9]. Besides, there can appear stationary structures with a tempo rally pulsating spatial period, when the emerging sources of traveling pulses coexist with the stationary structure [10]. Such dynamics is observed in some real systems belonging to active media, in particular, car diac tissue, and is currently associated with arrhyth mias. In the present paper, with an example of an analyt ically solvable Rinzel­Keller model [11] qualitatively reproducing the FitzHugh­Nagumo dynamics, sev eral spatially oscillating solutions are found both for single pulses and for pulse trains. It is shown that an anomalous dispersion relation is admissible as well for a model presented with a single reaction­diffusion

631


632

ZEMSKOV, LOSKUTOV

equation where oscillations are enforced by a periodi cal external stimulus. UNPERTURBED REACTION­DIFFUSION SYSTEM The two component system of the reaction­diffu sion type, with which we start the consideration of oscillatory solutions, can be composed of the follow ing equations: u u = f(u) ­ v + 2 , t x
2

When ±, the pulse tends to a constant value. The general expressions for pulse (3) are supple mented with linking conditions for the functions un(), vn(), n = 1, 2, 3 and their derivatives dun()/d, dvn()/d in two points = 0 and = * where two 0 neighboring parts of the solution are joined: u1 ( 0 ) = u2 ( 0 ) , du 1 ( 0 ) du 2 ( 0 ) , = d d v1 ( 0 ) = v2 ( 0 ) , d v1 ( 0 ) d v2 ( 0 ) , = d d * * u2 ( 0 ) = u3 ( 0 ) , du 2 ( * ) du 3 ( 0 = d d v2 ( * ) = v3 ( * ) 0 0 *) 0 , , (4)

(1) 2 v v = (u ­ v) + 2 . t x The first one describes the dynamics of the activa tor; the second one, of the inhibitor. The reaction term f(u) in the first equation is a piecewise linear function f(u) = ­u ­ 1 + 2(u ­ u0) approximating the cubic nonlinearity f(u) = u ­ u3 + u0. Here (u ­ u0) is a Heaviside step function, , , and u are model param eters (constants). We are interested only in such solutions of (1) that propagate in space without change in shape or velocity. These solutions are called autowaves. Introducing the autowave variable = x ­ ct, where c is wave velocity, we can write (1) as: du d
2 2 2

d v2 ( * ) d v3 ( * ) 0 0 = . d d

Since we know the u() values at the linking points, hence we can obtain two more additional equations u1(0) = u3( * ) = u0. 0 Owing to the translational invariance of the model equations, the position of one linking points can be taken arbitrarily, so for simplicity we take 0 = 0. Then the linking equations are solved and the integration constants An and pulse velocity c are determined there from. Figure 1 displays the solitary pulses for activator u = u() and inhibitor v = v(), and the (u, v) phase diagram. Calculations show that two pulses can coex ist, fast and slow. This reproduces the situation with non oscillating pulses as obtained elsewhere [11], where it was shown that the fast pulse is a stable solu tion, and the slow one is unstable. Since the pulse velocities take positive values, the waves in Fig. 1a,b propagate from left to right, i.e. vividly pronounced oscillations in the wave profile are in the tail while the front is visually smooth. In the (u, v) phase dia gram such profile oscillations appear as a spiral. Trains of Pulses

+c

(2) dv + c + (u ­ v) = 0, 2 d d where u = u() and v = v(). Spatially oscillating solutions of such a system in general form have been obtained by us earlier [6]. Using these results, let us build solutions for single pulses and trains thereof. For convenience, the main relationships from [6] are given in the Appendix. dv Single Traveling Pulses Let us first build single (solitary) pulses. Such solu tions consist of three parts: peak and two (rising and fading) "tails." Using the expressions for each part (Appendix), we can present the resulting pulse as: * u1 ( ) = u1 ( ) + u1 ,
+ ­ u 2 ( ) = u 21 ( ) + u 22 ( ) + u * , 2 ­ u3 ( ) = u3 ( ) + u* , 1 + v1 ( ) = v1 ( ) + v* , 1 + ­ v 2 ( ) = v 21 ( ) + v 22 ( ) + v * , 2 ­ v3 ( ) = v3 ( ) + v* . 1 +

du + f(u) ­ v = 0, d

(3)

Periodic sequences of pulses may be observed in the same excitable medium where single pulses are seen. Unlike the latter, the number of trains can differ depending on parameter values. Thus upon changing the wave period, one can observe one or several trains with different velocities. [Here we imply waves of com plex shape with profile oscillations.] The diagrams for the wave velocity vs. period are called dispersion rela tions. For waves of a common profile, the dispersion relation curves are monotonic. For waves with profile oscillations the dispersion relations are anomalous [7]. Now let us build pulse trains. Such solutions for one wave period in the piecewise linear model consist of two parts, left and right (in the phase diagram). Using the corresponding expressions from the Appen dix, the resulting wave can be written as:
BIOPHYSICS Vol. 54 No. 5 2009


TRAVELING WAVES IN A MODEL EXCITABLE MEDIUM
+ ­ u 1 ( ) = u 11 ( ) + u 12 ( ) + u * , 1 + ­ u 2 ( ) = u 21 ( ) + u 22 ( ) + u * , 2 ­ + v 1 ( ) = v 11 ( ) + v 12 ( ) + v * , 1 + ­ v 2 ( ) = v 21 ( ) + v 22 ( ) + v * . 2

633

(a)

u 15 10 5

(5)

It must be taken into account that the trajectory for a pulse train on the (u, v) phase plane is a smooth closed line, so here we again have two linking points. The conditions in the first point are the same as for single pulses. However, in the second point we have to link two parts of solution with different values: the 0 curve "sets out" from point u = u0 at = 0 , passes again through this point at = 0 and "winds up" in the same point at = * . Thus, the linking conditions 0 are described as follows: u1 ( 0 ) = u2 ( 0 ) , du 1 ( 0 ) du 2 ( 0 ) = , d d v1 ( 0 ) = v2 ( 0 ) , d v1 ( 0 ) d v2 ( 0 ) = , d d
0 u2 ( * ) = u3 ( 0 ) , 0 0 du 2 ( * ) du 1 ( 0 ) 0 = , d d 0 v2 ( * ) = v1 ( 0 ) , 0

-80

-60

-40

-20 0 -5 -10

20

40

(b)

v 1.5 1.0 0.5

(6)

-80

-60

-40

-20 -0.5

20

40

d v2 ( * ) d v1 ( ) 0 = . d d

0 0

-1.0 -1.5 (c) v 1.5 1.0

0 Hence it follows that L = * ­ 0 gives the wave 0 period. Comparing the number of linking equations with the number of unknown constants, we can con clude that period L is now a new additional parameter of the solution. Figures 2 and 3 give several examples of trains of pulses u(), v() and the (u, v) diagrams for L = 57. At this value there exist four waves, three of which are shown here (the fourth is very close in profile and velocity to the intermediate one in the Figs. 2b, 3b). All the waves travel from left to right. The fast wave has the simplest profile, while the slow waves have a more complex shape because of oscillations. The latter in the phase diagram deform the closed curve represent ing the train. The dispersion relation for a similar two component model [7] exhibits a complex oscillatory behavior, pointing to the possible multiplicity of waves at certain period values.

0.5 -10 -5 0 -0.5 5 10 u

-1.5
Fig. 1. Pulse profiles for (a) activator u() and (b) inhibitor v(), and (c) the (u, v) phase diagram. The thick line is the fast and the thin line is the slow pulse. Waves plotted at u0 = 2.5, = 0.0107, = 0.01. Isoclines f(u, v) = 0 and g(u, v) = 0 in (c) are shown as dashed lines.

PERTURBED REACTION­DIFFUSION SYSTEM: PULSE TRAINS The wave profile oscillations described above were caused by the presence of the second component in the model equations. However, spatial oscillations in solution can also be observed in a system of one equa tion if it is subjected to an external perturbation of spe cial kind. Such an external stimulus resembles a trav eling wave and is expressed as follows (see [12]):
BIOPHYSICS Vol. 54 No. 5 2009

I ( x, t ) = h cos ( t + k x ) + I 0 , (7) where h, k, and I0 are some constants, so that the model equation acquires the form: u = f ( u ) + I ( x, t ) + u . 2 t x
2

(8)


634 u, v (a) 8 6 4 2 10 0 -2 -4 (b) 10 5 10 0 -5 -10 40 30 20 10 0 -10 -20 10 20 30 40 50 (c) 20 30 40 50 20 30 40 50

ZEMSKOV, LOSKUTOV (a) v 1.5 1.0 0.5 -30 -20 -10 10 -0.5 -1.0 v 1.5 1.0 0.5 -30 -20 -10 10 -0.5 -1.0 -1.5 (c) v 3 2 1 -30 -20 -10 10 20 30 u -1 -2 -3
Fig. 3. Phase diagrams for trains shown in Fig. 2 (L = 57): (a) slow, (b) intermediate, and (c) fast waves. Dashed are isoclines f(u, v) = 0 and g(u, v) = 0.

20

30 u

-1.5

(b)

20

30 u

Fig. 2. Profiles for trains of one period (L = 57): (a) slow, (b) intermediate, and (c) fast waves. Thick and thin lines are for activator u and inhibitor v, respectively. Waves plot ted at u0 = 2, = = 0.01.

This effect is called autowave modulation of excit able medium. Here we will use a special case of function (7) when the velocities of modulation and of the traveling wave (the solution sought for) coincide, i.e., in the traveling wave equation d u + c du + f ( u ) + I ( ) = 0 2 d d
2

u() =

n=1



2

An e

n

+ u ( ) + u *,

A n, u * = const , (10)

where the eigenvalues of n are specified as 1, 2 = ­c/2 ± c /4 + 1 . Since the simplest case of an external stimulus of type (7) is described by function I() = h cos(k), the partial solution can be presented as u () = R cos(k) + Q sin(k), where R, Q = const. Now, to determine R and Q we must substitute this partial solution into (9) and collect all terms containing cos(k) and sin(k). As the result we obtain: R=h k +1 (k + 1) + c k
2 2 22 2 2

(9)

the function of external stimulus I depends only on the traveling wave coordinate and has the form of (7). Choosing as before the piecewise linear function f(u) = ­u ­ 1 + 2(u ­ u0) as the reaction term, we write down the general solution u() as

,

Q = ­h

. 2 22 ( k + 1 ) + c k (11)
2

ck

BIOPHYSICS

Vol. 54

No. 5

2009


TRAVELING WAVES IN A MODEL EXCITABLE MEDIUM u 0.2 0.1 0 -0.1 u 0.5 0 -0.5 -1.0 u 1.5 1.0
c>0

635

(a)

c 3

2 1 2 3 4 5 6

1

0 (b) -1 1
c>1

2

4

6

8

12

14 L

2

3

45

0
6

7

Fig. 5. Dependence of velocity on the period of the pulse train (dispersion relation).



CONCLUSIONS At present there are two distinct methods for math ematical description of active media: axiomatic and dynamic. The goal of the axiomatic theory initially consisted in describing qualitatively the processes of wave propagation in biological excitable media, such as pulse propagation in nerve and muscle tissues. Application of the axiomatic theory does not require detailed knowledge of the kinetics of real processes, which allows considering a broad class of problems in general form. However, with this approach it is impos sible to observe finer effects or to achieve quantitative fit to experimental data. In the dynamic approach, an active medium can be described in much more detail using the correspond ing equations in partial derivatives. The dimensional ity of the medium in such a system is determined by the number of independent spatial variables. Such reaction­diffusion systems prove quite complicated for analytical investigation. Thus no exact solutions in general form have yet been found for spiral waves. Therefore, approximate methods are applied. The best known among such asymptotic methods is the so called kinematic approach designed specially for excitable media. Kinematics remains a powerful tool for description and analysis of autowave structure in excitable media. Its development involves account of the effects of refractoriness and recovery. In the framework of the kinematic theory one can describe many processes inherent in active media, including some types of arrhythmia. Since arrhythmias are caused by certain disturbances in the cardiac mus cle, modeling of such systems is of great practical importance and can help solve the question of the feasi bility of controlling their behavior by external stimuli.

(c)

0.5 0 -0.5 -1.0 2 4 6 8 10

-0.5 < c < 0

c < -0.5

Fig. 4. Profiles for trains of period (a) L = 1.8, L = 4, L = 5; (b) L = 6.2, and (c) L = 10.

Figure 4 shows the wave profiles u = u() for differ ent wave periods. At small periods (L = 1, ..., 6) there is only one train for each period (Fig. 4a), whereas at large periods (L > 6) there are several trains (Fig. 4b,c) differing both in profile and in velocity. At that, the dis persion relation exhibits anomalous behavior (Fig. 5) as required by the oscillatory profile of the waves. The complex shape of the dispersion curve is caused by the interaction of pulses within the train via their oscillatory tails. This gives rise to regions of insta bility [13] in the dispersion curve, its rupture and for mation of gaps [14]. Besides, so called isolas can form for the train [13]. In the 2D case the described wave dynamics leads to a complex process of spiral wave for mation [8, 15].

BIOPHYSICS

Vol. 54

No. 5

2009


636

ZEMSKOV, LOSKUTOV

ACKNOWLEDGMENTS We are grateful to H. Engel and G. Bordyugov for extensive discussions and useful comments. The work was supported by the Russian Founda tion for Basic Research (07 01 00295). APPENDIX Presented are the main relations of the type of trav eling waves with profile oscillations [6]. The general form of u() and v() solutions for set (2) is a sum of exponentials: u() = v() =

where r= z= ( c b + d + c b ) /2 , ( c b + d ­ c b ) /2 , c b = c /4 + b
± ^± ± ^± and A = aA± ± dB , B = aB± ± dA , d = 2 2 2 2 2

(A5)

­a .

2

n=1


4

4

REFERENCES
1. A. D. Bazykin, Mathematical Biophysics of Interacting Populations (Nauka, Moscow, 1985) [in Russian]. 2. V. A. Vasil'ev, Yu. M. Romanovskii, and V. G. Yakhno, Autowave Processes (Nauka, Moscow, 1987) [in Rus sian]. 3. A. Yu. Loskutov and A. S. Mikhailov, Basic Theory of Complex Sysems (RKhD, Moscow, 2008) [in Russian]. 4. R. FitzHugh, Biophys. J. 1, 445 (1961). 5. J. Nagumo, S. Arimoto, and Y. Yoshizawa, Proc. IRE 50, 2061 (1962). 6. E. P. Zemskov, V. S. Zykov, K. Kassner, and S. C. MÝller, Nonlinearity 13, 2063 (2000). 7. G. Bordiougov and H. Engel, Phys. Rev. Lett. 90, 148302 (2003). 8. N. Manz, C. T. Hamik, and O. Steinbock, Phys. Rev. Lett. 92, 248301 (2004). 9. I. M. Tsyganov, M. A. Tsyganov, A. B. Medvinsky, and G. R. Ivanitsky, Dokl. RAN 346, 825 (1996). 10. I. M. Tsyganov, R. R. Aliev, and G. R. Ivanitsky, Dokl. RAN 352, 699 (1997). 11. J. Rinzel and J. B. Keller, Biophys. J. 13, 1313 (1973). 12. S. Zykov, V. S. Zykov, and V. Davydov, Europhys. Lett. 73, 335 (2006). 13. G. RÆder, G. Bordyugov, H. Engel, and M. Falcke, Phys. Rev. E 75, 036202 (2007). 14. M. Falcke, M. Or Guil, and M. BÄr, Phys. Rev. Lett. 84, 4753 (2000). 15. O. Steinbock, Phys. Rev. Lett. 88, 228302 (2002).

An e

n

+ u *, + v *, (A1)

n=1



Bn e

n

u *, v * = const , where An, Bn, u* and v* are constants defined in each of the regions u < u0 and u > u0. Constants Bn are expressed through An. The n eigenvalues are of the form:
1, 2

= ­ c /2 + c /4 + b ± a ­ ,
2 2

2

2

(A2)

3, 4 = ­ c /2 ­ c /4 + b ± a ­ , with a = ( ­ )/2 and b = ( + )/2. Parameters and are positive constants, but ( ­ )2/4 may be smaller ­ + ± than . This is true when im < < im , where im = + 2 ± 2 + 1 . Then n has an imaginary part so that solutions u ( ) = u ( ) + u ( ) + u *, v ( ) = v ( ) + v ( ) + v *, contain cosines and sines: u () = e
± ± ( ­ c /2 ± r ) ( ­ c /2 ± r ) + ­ + ­

(A3)

v () = e

[ A cos ( z ) + B sin ( z ) ] , ^± ^± [ A cos ( z ) + B sin ( z ) ] ,

±

±

(A4)

BIOPHYSICS

Vol. 54

No. 5

2009