Original Research ARTICLE
Initiation of Rotors by Fast Propagation Regions in Excitable Media: A Theoretical Study
- 1Max-Planck-Institute for Dynamics and Self-Organization, Göttingen, Germany
- 2School of Physics and Information Technology, Shaanxi Normal University, Xi'an, China
- 3German Center for Cardiovascular Research, Göttingen, Germany
- 4Institute for Nonlinear Dynamics, University of Göttingen, Göttingen, Germany
- 5Cornell University, Ithaca, NY, United States
We study the effect of geometry of a fast propagation region (FPR) in an excitable medium on the rotor initiation using a generic two-dimensional reaction-diffusion model. We find that, while the flat boundary of a rectangularly shaped FPR may block the propagation of the excitation wave, a large local curvature at the rounded corners of the FPR would prevent the blockage and thus initiate a rotor. Our simulations demonstrate that the prerequisites for the rotor initiation are the degree of the heterogeneity, its shape and size. These results may explain the incidence of arrhythmias by local heterogeneities induced, for example, by a cardiac tissue remodeling.
Rotors, also known as spiral waves, are observed in many systems, including the Belousov-Zhabotinsky chemical reactions [1–4], autocatalytic reactions of carbon monoxide on a platinum surface , aggregations of Dictyostelium discoideum amoebae , Xenopus oocytes , disinhibited mammalian neocortex , chicken retina , and especially cardiac tissue [10, 11]. Rotors, resulting in reentries in heart tissue, are known to cause cardiac arrhythmias and even sudden death [12–14]. To understand the mechanism of the rotor initiation and to eliminate the consequential malignant arrhythmias, the effects of the electrophysiological heterogeneity are thought to be one of the major causes and have attracted much attention [15–20].
Destabilization of wave fronts and the subsequent initiation of reentrant excitation can result from both intrinsic and dynamical heterogeneity. For example, a possible result of a multiple pacing of cardiac tissue is a dynamically induced heterogeneities of repolarization leading to a destabilization of a propagating wave and initiation of a self-sustained activity [21–25]. The heterogeneities of the electrical coupling and automaticity also might lead to the appearance of fragmented ectopic waves . Furthermore, the boundary layer between the well-coupled and uncoupled cardiac tissues would create a rich set of phenomena associated with self-organized spiral waves and ectopic waves . Transient rotors could be also initiated in complicated and rare situations [28, 29]. An abrupt transition of the coupling gradient would block the wave propagation, but nearby parts with a smooth transition would not and therefore cause a reentry. The wave blockage was also found in a model of human ventricular tissue due to an abrupt transition of the anisotropic coupling .
As exemplary mentioned above, there are many situations, which can lead to the initiation of a self-sustained excitation wave. One novel scenario, which was found recently in a generic model for the excitable system, is that a region with the fast propagation of an excitation wave might cause a unidirectional block of the wave propagation . This unidirectional block is realized based on a phenomenon termed source-sink mismatch in the cardiology literature . The further study demonstrated that rotors could be nucleated in the presence of a localized fast propagation region (FPR) after application of one stimulus only . It was shown also that various geometrical factors play an important role in rotor initiation .
Here we show that in the two-dimensional medium the flat boundary of a rectangularly shaped FPR may block the propagation of the outside excitation wave. However, a large local curvature at the rounded corner of the FPR would prevent the blockage and thus let the outside excitation penetrate into the FPR causing a rotor initiation. We demonstrate that the rotor initiation critically depends on the size of the FPR and the degree of the heterogeneity. If the FPR size is below a certain threshold, the initiated rotor would vanish eventually when it approaches the medium's boundary. The critical size of the FPR depends on the degree of the heterogeneity.
2. Model and Method
Although some aspects of the complex electrical activity in the cardiac tissue need to be studied using the reaction-diffusion equations with the detailed ionic channel model, many general spatiotemporal features of cardiac dynamics can be reproduced by a relatively simple but universal two-component system as follows
where u and v are the activator and inhibitor, respectively. The local kinetics of u and v are specified by the nonlinear functions F(u, v) and G(u, v). Let us consider a widely-used computationally-efficient model proposed by Barkley . In this generic model, the two nonlinear functions read as
To simulate a relatively quick recovery of the excitability after a pulse generation, the original Barkley model is slightly modified by introducing an additional parameter kϵ > 1.
Below a FPR is considered to be a rounded rectangular region of the length L with two rounded corners of the radius R, as illustrated in Figure 1A. Inside this region, the values of D and/or A are larger than outside. It is inserted into a square shaped medium of size 450 × 450 in space unit, and the no-flux boundary conditions at the boundaries are implied. The parameters a = 1, b = 0.44, ϵ = 0.00011 and kϵ = 10 are fixed in our simulations. We use the explicit finite difference method in the Cartesian coordinates. The staircase approximation is used at the rounded corner. The spatial step dx = 0.3, and the time step dt = 0.01, when R ≥ 15. The finer spatial and time steps are used for smaller R. For instance, dx = 0.2 and dt = 4.44 × 10−3, when R = 10, and dx = 0.1 and dt = 1.11 × 10−3, when R = 5.
Figure 1. Rotor initiated due to a rectangularly shaped fast propagation region (FPR). (A) The FPR marked by a black dashed line has length L = 300 and two rounded corners with the radius of R = 15. Within the FPR D = 1, A = 2 and outside the FPR D = 1 and A = 1. A plane wave propagates from the left toward the FPR. (B) The wave is blocked at the flat boundary of the FPR but penetrates into it at the rounded corner. (C–E) A rotor is initiated and stably rotates. The white dot and line are the rotors phase change point (PCP) and its trajectory, respectively.
The variables u and v in Equations (1)-(4) are vary within the range 0 < u < 1 and 0 < v < a−2b. The spatiotemporal dynamics of u and v is represented in Figures 1–4, 6, 7 by color-coded distribution of the excitation phase ϕ, where −π < ϕ < π. The phase is defined as ϕ = α + 3π/4, in which an angle α is determined by the direction of the vector with components (u − 1/2) and (v − a/2 + b)/(a − 2b) on the (u, v) phase plane. According to this definition, ϕ = 0 corresponds to the resting state of the medium (green areas in Figures 1–4, 6, 7), narrow yellow (dark blue) regions represent the wave front (wave back), and red areas correspond to a wave plateau, whereas blue ones represent the refractory regions.
3.1. Rotor Initiated from a Rectangularly Shaped FPR
To illustrate the phenomenon of the rotor initiation from a rounded rectangular heterogeneous region, we set L = 300, R = 15, as shown in Figure 1. Inside this region D = 1 and A = 2, while outside D = 1 and A = 1. Due to an increased value of A, this rounded rectangular region could be considered as a FPR since the propagation velocity inside it is larger than outside. For such parameter choice, the flat boundary of the FPR would unidirectionally block a plane wave propagating through the medium outside FPR. However, the local curvature at the rounded corner of the FPR 1/R ≈ 0.067 is large enough to prevent the blockage and let the excitation penetrate into the FPR. Then, a phase change point (PCP) emerges, and a self-sustained rotor is initiated. The process is similar to the scenario described in Zykov et al. [33, 34].
However, if the local curvature at the rounded corner of the FPR is below some critical value, as illustrated in Figure 2, there would be no penetration into the FPR at the corner. Thus, the FPR would act for propagating waves as an obstacle. The transient rotor starts to circulate around the FPR and vanishes eventually when it approaches the medium boundary. The critical value of the corner curvature depends on D and A within the FPR.
Figure 2. No wave penetration occurs at the rounded corner of a FPR if R is above some critical value. (A) No penetration is observed at the corner of radius R = 30, as an example. (B–D) The PCP (white dot) would travel along the FPR boundary and vanishes eventually, as shown by its trajectory (white line). Other parameters are the same as in Figure 1.
Another scenario appears when the values of D and/or A within the FPR are below some critical values. As illustrated in Figure 3, in this case, the flat boundary of the FPR would not block the propagating wave. The plane wave would become curved, propagates through the FPR, and vanishes eventually when it reaches the medium's boundary. No self-sustained rotors are initiated.
Figure 3. No blockage occurs at the flat boundary of a FPR if D and/or A inside it are below some critical values. (A) No PCP would emerge for D = 0.9 and A = 1.8, as an example. (B) A slightly deformed plane wave propagates through the FPR and disappears at the right boundary of the medium. Other parameters are the same as in Figure 1.
3.2. Critical Length of a Rectangularly Shaped FPR
For a given D and A within the FPR, its length L is also a critical parameter to initiate a self-sustained rotor. If L is shorter than a certain critical length, as illustrated in Figure 4, the PCP would approach the medium boundary and vanishes eventually. After the PCP has disappeared, a curved wave propagates through the medium and vanishes at its boundary. Thus, no rotors are initiated.
Figure 4. No stably rotating rotor appears if L is shorter than a critical length. (A) Although a PCP would be initiated like in the case shown in Figures 1A–C, the PCP (white dot) would be too close to the medium boundary and vanishes eventually, as shown by its trajectory (white line). (B) After the PCP vanishing the curved wave would propagate to the right boundary of the medium, and no stably rotating rotor would exist. The FPR length L = 200 and other parameters are the same as in Figure 1.
We investigate in detail the critical length of the FPR needed to initiate a self-sustained rotor and determine the boundaries of the non-penetration and non-blockage regions for various A and D, as shown in Figure 5. Here the radius of the rounded corner is fixed as R = 15. As shown in Figure 5A, the non-penetration occurs when D and/or A are above some critical values, and the non-blockage occurs when D and/or A are below other critical values. Between these two boundaries, the initiation of a self-sustained rotor is possible with the color-coded critical length Lc of the rounded rectangular FPR. The larger D and/or A are, the shorter Lc would be.
Figure 5. (A) The phase diagram shows the regions of non-penetration, non-blockage, and rotor initiation due to the FPR with critical length Lc for different values of D and A but fixed R = 15. (B) The detailed dependence of Lc on A for fixed D = 1 and different R. (C) The detailed dependence of Lc on D for fixed A = 2 and different R.
This dependence of Lc on D and A is further confirmed in Figures 5B,C, where we show how Lc shrinks as D or A increases for different corner radius R. Figures 5B,C also demonstrate how Lc changes as R increases.
We also investigate the impact of the width of the rectangular FPR. The simulation results show that the FPR width has no significant effect on the rotor initiation if it is larger than 2R.
To analyze the conditions for the rotor initiation by a rounded rectangular FPR, we simplify the two-component reaction-diffusion equations by taking ϵ = 0 and setting v = 0. In this limiting case, the initial Equations (1)–(4) can be reduced to
This equation describes a bistable extended system, where the resting state u = 0, the excited state u = 1, and the unstable steady state u = b/a exist. The value β = b/a is the excitation threshold. The bistable equation has been widely used to analyze the propagation of the wave front when ϵ ≪ 1 . It is also useful to establish fundamental mechanism behind the blockage and penetration at the FPR boundary.
4.1. Analysis of the Non-blockage and Non-penetration Boundary
To analyze the conditions for the blockage of the initial plane wave at the flat boundary of the FPR, we could further simplify Equation (5) to consider a stationary wave profile for a one-dimensional bistable system as follows
where A(x) = 1, D(x) = 1 for x ≤ 0 and A(x) = A, D(x) = D for x > 0. The boundary conditions and the continuity conditions at x = 0 read as
Multiplying Equation (6) by du/dx, integrating over x from −∞ to 0 and from 0 to ∞, using Equations (7) and (8), we obtain the following equation
which determines the value of u(0) at the point of the parameter jump as a function of the product DA. Note that the front could be stopped only if u(0) < β. Thus, the Equation (9) for u(0) = β gives the critical value of the product DA, above which the propagation blockage could be observed. Therefore, the non-blockage boundary in the phase diagram has an analytical form as follows
This expression represents a modification of another one mentioned already in Zykov et al. [33, 34]. It looks more simple because of the normalized values of the parameters A = 1 and D = 1 in the part of the medium outside an FPR. This normalization performed by rescaling of time and space variables in Equation (5) is made without loss of generality. It is worth to note also that this expression generalizes a similar analysis for the case of a steep rise of the parameter D under constant value of A performed earlier in Pauwelussen  and Mornev .
It is important to stress that the obtained analytical expression (10) gives a very precise estimate of the non-blockage boundary obtained by numerical computations illustrated by Figure 5A. The deviations do not exceed one percent.
To analyze the conditions to prevent wave penetration at the rounded corner of the FPR, we look in detail into the process of the (non-)penetration, as illustrated in Figure 6A. The initial plane wave would become curved at the rounded corner of the FPR. It is analogous to a circular wave penetrating into a circular FPR with a radius r, as shown in Figure 6B.
Figure 6. Analogy between the non-penetration at the corner of a rectangularly shaped FPR and the non-penetration of a circular wave into a circular shaped FPR. (A) The non-penetration at the rounded corner of the FPR (zoomed from Figure 2A). (B) A circular wave blocked by the circular shaped FPR of the radius r = 60. (C) The non-penetration boundaries for the rectangularly shaped FPR for R = 15 (black dots) and for the circular wave and the circular shaped FPR with the radii rmin = 27.6 (red line) and rmax = 29.7 (blue line). (D) The dependence of 1/r = 1/(2R) is proofed to be valid for a large range of R.
To verify the analogy, we compare the non-penetration boundaries for the rounded rectangular FPR with the circular FPR in a A − D diagram. As shown in Figure 6C, the non-penetration boundary for the rectangularly shaped FPR with the corner radius of R is located between two curves corresponding to non-penetration boundaries for the circular FPR with the radius rmin and rmax, as shown in Figure 6C. Note, that in a vicinity of the rounded corner the boundary curvature jumps from zero to 1/R. Hence, it is natural to assume that the non-penetration boundaries for the circular FPR with the curvature about an averaged curvature of the rounded corner, 1/r ≈ 1/(2R), will approximate the non-penetration boundary for the rounded rectangular FPR. This approximation is working well for 10≤R≤50, as demonstrated in Figure 6D.
Therefore, we can use the results of the simulations for the circular FPR and the radius relation r ≈ 2R to approximate the non-penetration boundary for the rounded rectangular FPR. Since the circular FPR in the polar coordinates (ρ, θ) has a rotational symmetry, this allows us to transform the two-dimensional Equations (1) and (2) to the one-dimensional ones as follows
This one-dimensional equations considerably simplifies the analysis. The corresponding computations have been performed by use of the explicit finite difference method with the spatial step dρ = 0.3 and the time step dt = 0.005. In order to simulate a circular wave approaching a circular shaped FPR, a part of the medium with ρ > ρext is assumed to be in the excited state at t = 0, as illustrated in Figure 6B.
4.2. Analysis of the Critical Length Lc of a Rounded Rectangular FPR
To understand the mechanism behind the dependence of Lc on the characteristics of the FPR, i.e., D, A and R, we separate Lc into three parts, as shown in Figure 7. The first part is the distance from the rounded corner of the FPR, where the penetration of the excitation starts, to the position where the PCP emerges for the first time. This part, named Lexc, should be determined by the characteristics of the FPR since it describes the propagation of the excitation inside the FPR. The second part is the distance from the initial position of the PCP to the highest position in its trajectory. This part, named Lr, describes the range of the PCP trajectory along the FPR but located outside the FPR. This trajectory part should be practically independent of the FPR characteristics. The third part is the distance from the highest PCP position in its trajectory to the top medium boundary. This part, named Lmin, should be above some minimum distance toward the top medium's boundary. Otherwise, the PCP would be too close to the boundary and vanishes eventually. The value of Lmin should only depend on the given characteristics of the medium outside the FPR, and thus is fixed in our simulations. Therefore, the value of Lc is the sum of Lmin and Lr, which are fixed, and Lexc, which is determined by D and A inside the FPR, and R at the FPR corner.
Figure 7. The composition of the critical length Lc of a rectangularly shaped FPR to initiate a rotor. The white dot and line are the initial location and the following trajectory of the PCP, respectively. Three components of Lc are: Lmin which describes the minimum distance of a persistent PCP from the top medium boundary, Lr which describes the vertical range of the PCP trajectory in the medium outside the FPR, and Lexc which describes the propagation of the penetrated excitation within the FPR. The FPR length L = Lc = 229 and other parameters are the same as in Figure 1.
Hence, Lexc is the only part which is varied and depends on D, A, and R. Its value may read as
where tp is the time when the penetration of the excitation at the rounded corner of the FPR starts, tr is the time when the PCP initially emerges, and cv is the propagation velocity of the excitation along the flat border of the FPR. As shown in Figure 8, tr remains constant for different sets of D, A, and R, while tp varies. The velocity cv changes with time and also depends on D, A, and R.
Figure 8. The propagation speed cv of the penetrated excitation along the vertical flat boundary of a FPR over time. tp is the time when the penetration at the FPR corner starts. tr is the time when the PCP initially emerges. Larger D, A, or R lead to a delayed tp, but nearly the same tr. The subfigure shows the change of the ratio between cv and cp over time, where cp is the plane wave velocity for the medium's parameters established inside the FPR.
Based on these results, three conclusions can be made. First, larger D, A, or R would delay tp. Second, tr is nearly the same in all cases since it is determined by the time when the wave back of the plane wave reaches the left flat boundary of the FPR. Therefore, it is determined by the fixed characteristics of the medium outside the FPR. Third, for the most part of the trajectory, cv is larger than the plane wave velocity cp in a homogeneous medium where the parameters D and A are the same as inside the FPR, as shown in the subfigure of Figure 8. Obviously, cv is accelerated since the value of the activator u > 0 in the vicinity of the left flat boundary of the FPR due to a diffusive influence of the blocked plane wave. Such increase of the propagation velocity is a general effect in bistable models of one-dimensional excitable media if ahead of the wave front the activator value exceeds the resting state . In the context of a flame propagation, which also can be described by Equation (6), this phenomenon is named as preheating effect .
4.2.1. Mechanism of Delay of tp
The delay of tp occurs at the rounded corner of a FPR. As shown above, the penetration of the plane wave into the rounded corner with the radius of R is analogous to the penetration of a circular wave into a circular FPR with the radius r ≈ 2R. Thus, considering the analogy and just focusing on the propagation of the excitation wave front, we could use the bistable version of Equations (11) and (12) to investigate the variation of u at the rounded FPR corner with time. The bistable equation for the circular FPR expressed in the polar coordinates (ρ, θ) reads as
Using the finite difference method with space step Δρ, at the circular FPR boundary r = 2R, Equation (14) could be expanded as
When the circular wave reaches the circular FPR boundary at r, we have u|r−Δρ < u|r < u|r+Δρ, and u|r is still smaller than the excitation threshold β. Thus, u|r(u|r − 1)(u|r − β) > 0. Then we can divide the right hand side of Equation (15) into two terms. The term which would cause an increase of u|r with time is named as the source term. It reads as
The other term which would cause a decrease of u|r with time is named the sink term. It reads as
From the above expressions of two terms, we find that larger D would enhance the source term (Equation 16) but enhances the sink term (Equation 17) even more. Larger A would not affect the source term but enhance the sink term. Larger r, i.e., 2R in the analogy, would reduce the source term, but enhance the sink term.
Therefore, the conclusion is that the larger D, A and R are, the stronger the sink term would be, and the later u|r reaches the excitation threshold. This is the cause of the delay of the start time of the excitation penetration near the corner of the rounded rectangular FPR, i.e., tp in Equation (13). As shown in Figure 9, the numerical simulation results prove our explanation of the delay effect by plotting the value of u at the corner of the rounded rectangular FPR with time.
Figure 9. Temporal dynamics of the value of u at the junction point between the flat boundary and the corner of a rectangularly shaped FPR. An increase of D, A, or R inside the FPR delays the time tp when the excitation threshold is reached.
4.2.2. Influence of “Preheating” on cv
Inside the rounded rectangular FPR, the accelerating effect on the propagation velocity cv occurs near its vertical flat boundary. When the initial plane wave is blocked at the flat boundary, although it does not penetrate inside the FPR, it yet increases the value of u in a vicinity of the FPR boundary. This is quite similar to a preheating effect in the flame propagation when the fuel temperature ahead of the flame front is increased [39, 40]. This “preheated” medium would accelerate the propagation velocity cv of the excitation wave front along the vertical flat boundary of the FPR.
The mechanism of the accelerated wave front could be analytically understood from Equation (5) for the bistable distributed system. If a preheated part of the FPR near its flat boundary is assumed as a nearly one-dimensional medium, we can establish a comoving frame as z = x + ct, where c is the propagation velocity of the wave front. Thus, Equation (5) would be simplified as
The preheating effect increases the value of u to some preheated state up, and makes the excitation start from up > 0 instead of the resting state u = 0. Based on the theory described in Keener and Sneyd , the propagation velocity of the excitation wave front could be expressed as
The analytical expressions of Equation (19) would be obtained for two limiting cases. The first one is the unpreheated case at which up is equal to the resting state. That gives
The second is the fully preheated case at which up is equal to the excitation threshold β. This gives
The above two analytical expressions apparently demonstrate that c(β) > c(0), since β > 0.
We also investigate the preheated propagation velocity in the numerical simulations of Equation (5). As shown in Figure 10, the numerical results elucidate the acceleration of the propagation velocity c(up) as a function of the preheated state up. The analytical results from Equations (20) and (21) perfectly describe both limiting cases following from these numerical data.
Figure 10. The propagation speed cp corresponding to the preheated state u = up. The values of D = 1 and A = 2 are taken as an example. The solid line represents the numerical simulation results in a one-dimensional medium described by Equation (5). The value of u ahead of the wave front is set to be up. The open circle and square are the analytic results from Equations (20) and (21) at the limiting cases of up = 0 and up = β, i.e., the resting state and the excitation threshold, respectively.
5. Conclusions and Applications
Our results demonstrate that a self-sustained rotor could be initiated from the spatial heterogeneity, i.e., a rectangularly shaped FPR. We use a generic model to parameterize the heterogeneity with three parameters D, A, and R. In the D − A diagram at a given R, the region of the rotor initiation is located between the non-blockage and non-penetration regions. The two boundaries of the rotor initiation region could be estimated by the analytical equation for the bistable distributed system and the simulations in a one-dimensional medium for a circular FPR, respectively. We also show that to initiate the self-sustained rotor the length of the rounded rectangular FPR should be larger than the critical Lc. The critical value Lc depends on the parameters D and A, within the FPR, as well as on the radius R of a rounded corner.
Our findings in the generic model might be applicable to describe the electrophysiological dynamics of cardiac tissue. Indeed, the distribution of transmembrane potential V in a two-dimensional tissue could be described by the reaction-diffusion equation as follows 
where χ is the surface-to-volume ratio of the cardiac cells, Cm is the membrane capacitance, σ is the tensor of electric conductivity, and Iion is the sum of ion channel currents. Intensity of each separate current is determined by corresponding component of the vector . Equations (22) and (23) can be generalized into a two-component reaction-diffusion system as follows
where the effective diffusion coefficient tensor D = σ/χCm and the description of the ion currents is reduced to a scalar value h. In an isotropic tissue, we can simplify the tensors σ and thus D to be scalars. Thus, the reduced system which describes electrophysiological properties of the cardiac tissue looks similar to the reaction-diffusion model we use.
Nowadays many detailed models of human atria incorporate both structural and electrophysiological heterogeneities leading to differences in conduction velocity between the neighboring regions [42–44]. It is also well known that atrial fibrosis in the aging heart can result in spatial variations in the electrical conductivity of a part of the cardiac muscle . If some regions within this part remain unchanged, they can resemble fast propagation regions introduced in our model. Note, that a similar nonhomogeneity in the electrical conductivity can appear, for instance, when fresh stem cells aggregates implanted in strongly remodeled cardiac tissue form gap junctions with adult cardiac myocytes . Moreover, some cardiac diseases cause ion channel remodeling . This remodeling can be represented as a variation of the term Iion in Equation (22). This is to some extent equivalent to a variation of the parameter A in our model. If this remodeling occurs non-uniformly in space, one can expect the creation of some spots with a negligible variation of this parameter in comparison to its strong decrease in the surrounding regions. Thus, nonhomogeneous remodeling can results in a creation of fast propagation regions considered in this study.
Of course, the model used above is aimed to reproduce only most generic features of electrical activity in myocardial tissue. Investigation of specific dynamical features can be done by application of more detailed models widely used in the literature [48–50]. It is important to note that our recent results based on the Fenton-Karma model  indicate that all scenarios of rotor initiation obtained with the Barkley model are perfectly reproducible . An obvious reason for this is that the restitution of action potential duration in detail reproduced in the Fenton-Karma model plays only a restricted role in the described scenarios, where spiral waves are generated after application of a single excitation stimulus. Of course, the following dynamics of the initiated rotors is strongly influenced by many other factors and specific features of cardiac tissue, which are not reproduced in the framework of the generic model used in the study.
Therefore, computer simulations of a real tissue in the framework of much more detailed models and most importantly experimental investigations definitely can help to verify the role of the observed scenario for generation of cardiac arrhythmias.
All authors conceived of the presented idea. XG, AK, VZ performed the computations. All authors discussed the results and contributed to the final manuscript.
This work was supported by the National Natural Science Foundation of China under Grants No. 11447026 and No. 11705114, the Max Planck Society, the German Center for Cardiovascular Research (DZHK), and Grant DFG-SFB937 Collective Behavior of Soft and Biological Matter.
Conflict of Interest Statement
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
5. Jakubith S, Rotermund HH, Engel W, Vonoertzen A, Ertl G. Spatiotemporal concentration patterns in a surface-reaction - propagating and standing waves, rotating spirals, and turbulence. Phys Rev Lett. (1990) 65:3013–6. doi: 10.1103/PhysRevLett.65.3013
9. Yu Y, Santos LM, Mattiace LA, Costa ML, Ferreira LC, Benabou K, et al. Reentrant spiral waves of spreading depression cause macular degeneration in hypoglycemic chicken retina. Proc Natl Acad Sci USA (2012) 109:2585–9. doi: 10.1073/pnas.1121111109
22. Cherry EM, Fenton FH. Suppression on alternans and conduction blocks despite steep APD restitution: electrotonic, memory, and conduction velocity restitution effects. Am J Physiol Heart Circ Physiol. (2004) 286:H2332–41. doi: 10.1152/ajpheart.00747.2003
23. Gelzer ARM, Koller ML, Otani NF, Fox JJ, Enyeart MW, Hooker GJ, et al. Dynamic mechanism for initiation of ventricular fibrillation in vivo. Circulation (2008) 118:1123–9. doi: 10.1161/CIRCULATIONAHA.107.738013
24. Bueno-Orovio A, Hanson BM, Gill JS, Taggart P, Rodriguez B. In vivo human left-to-right ventricular differences in rate adaptation transiently increase pro-arrhythmic risk following rate acceleration. PLoS ONE (2012) 7:e52234. doi: 10.1371/journal.pone.0052234
25. Bueno-Orovio A, Cherry EM, Evans SJ, Fenton FH. Basis for the unduction of tissue-level phase-2 reentry as a repolarization disorder in the Brugada syndrome. BioMed Res Int. (2015) 2015:197586. doi: 10.1155/2015/197586
30. Kudryashova NN, Kazbanov IV, Panfilov AV, Agladze KI. Conditions for waveblock due to anisotropy in a model of human ventricular tissue. PLoS ONE (2015) 10:e0141832–849. doi: 10.1371/journal.pone.0141832
38. Mornev OA. Elements of optics of autowaves. In: Krinsky VI, editor. Self-Organization: Autowaves and Structures Far from Equilibrium Springer Series in Synergetics, Vol. 28, Berlin; Heidelberg: Springer (1984). p. 111–8.
42. Seemann G, Höper C, Sachse FB, Dössel O, Holden AV, Zhang H. Heterogeneous three-dimensional anatomical and electrophysiological model of human atria. Philos Trans R Soc A (2006) 364:1465–81. doi: 10.1098/rsta.2006.1781
43. Xie Y, Garfinkel A, Camelliti P, Kohl P, Weiss JN. Effects of fibroblast-myocyte coupling on cardiac conduction and vulnerability to reentry: a computational study. Heart Rhythm (2009) 6:1641–9. doi: 10.1016/j.hrthm.2009.08.003
44. Sanchez C, Bueno-Orovio A, Pueyo E, Rodriguez B. Atrial fibrillation dynamics and ionic block effects in six heterogeneous human 3D virtual atria with distinct repolarization dynamics. Front Bioeng Biotechnol. (2017) 5:29. doi: 10.3389/fbioe.2017.00029
47. Nattel S. Effects of heart disease on cardiac ion current density versus current amplitude: important conceptual subtleties in the language of arrhythmogenic ion channel remodeling. Circul Res. (2008) 102:1298–300. doi: 10.1161/CIRCRESAHA.108.178087
Keywords: excitable media, rotor initiation, arrhythmia, source-sink mismatch, fast propagation region
Citation: Gao X, Krekhov A, Zykov V and Bodenschatz E (2018) Initiation of Rotors by Fast Propagation Regions in Excitable Media: A Theoretical Study. Front. Phys. 6:8. doi: 10.3389/fphy.2018.00008
Received: 01 November 2017; Accepted: 26 January 2018;
Published: 13 February 2018.
Edited by:Simonetta Filippi, Università Campus Bio-Medico, Italy
Reviewed by:Alfonso Bueno-Orovio, University of Oxford, United Kingdom
Winfried Mayr, Medizinische Universität Wien, Austria
Copyright © 2018 Gao, Krekhov, Zykov and Bodenschatz. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Vladimir Zykov, firstname.lastname@example.org