| |
ISSN 1320-0682 | |||
| Volume 02 | April 1995 | |||
Mark L. Spano,
William L. Ditto
School of Physics
John Lindner
The College of Wooster
U.S. Naval Surface Warfare Center,
Silver Spring, Maryland 20903
Email:mark@chaos.nswc.navy.mil
Georgia Institute of Technology
Atlanta, Georgia 30332, USA
Email:wditto@acl.gatech.edu
Wooster, Ohio 44691, USA
Although the theory [1] and the implementation [2] of the control of chaos have been demonstrated on many experimental systems since 1990, we still lack a general understanding of the control of spatiotemporal chaos. Many papers have been written [3] on the characterisation of spatiotemporal chaos, but only recently has interest turned to its control [4].
In this paper, we present our efforts to simulate arrays of coupled oscillators for two cases: a spherical array of Duffing oscillators and a planar array of heart cells.
In an attempt to understand this most complicated of problems, we have developed a program to simulate a one- or two-dimensional array of coupled "oscillators". The code runs on advanced personal computers under the Microsoft Windows/Windows NT operating systems at a quite reasonable integration rate.
The choices of host computer and of operating system are nontrivial and set the basis for the entire simulation. The governing criteria are speed and cost. Since the speed of modern PC processors such as the Intel Pentium rivals that of many traditional workstations and since PC's may be obtained for a fraction of the cost of such a workstation, we chose to work on a 60 MHz Pentium PC with 64 megabytes of memory and a large hard disk. (The cost of such a machine is currently between $6,000 and $8,000.) For the operating system, we decided on Microsoft Windows NT.
The major concerns here were speed, memory and portability. Windows NT is a fully 32-bit operating system and provides a flat address space, unlike other PC operating systems such as DOS. The 32-bit feature essentially doubles the speed of calculation as compared to DOS. NT also provides full access to the 64 Mb of memory on the machine, a feature which is vital when using the huge arrays necessary for such modelling. In addition, NT provides virtual memory, using the hard disk as an alternative to physical memory, if needed. If a bit of care is taken during programming to avoid using certain esoteric system calls, the resultant program may also be run under Windows 3.1 with little degradation in speed (after installation of a library called Win32s). Thus, the program achieves a high degree of portability. The results shown in this paper on the van Capelle-Durrer model were obtained on this computer/operating system platform. (Another operating system that provides equivalent functionality is the NextStep operating system. The program has been ported to that operating system and the results presented here on the Duffing model on the sphere were obtained on that platform.)
Another advantage of the Windows NT operating system is the fact that, like many modern operating systems, it has a graphical user interface. Thus, it is possible to have available a number of view windows at any given time. We have taken full advantage of this fact to provide plots of position and velocity versus time, of position and velocity delay co-ordinate embeddings (the variable versus the same variable at some previous time) and of phase space (velocity versus position) at each site of the array. In addition, each of these is available for both the flow data and, where appropriate, for section data (data strobed at the driving frequency). The phase space section is the plot of interest for implementing chaos control at a given site according to such prescriptions as the algorithm of Ott, Grebogi and Yorke [1]. Of course, no completely general algorithm exists for implementing global control of the spatiotemporal chaos.
This program generally treats an 
array using a fourth order
Runge-Kutta method. Although faster methods of numerical integration
exist, the fourth order Runge-Kutta method is the fastest and most
stable
non-variable step size technique, commonly available. We have found that,
since the frequency content of chaotic systems generally spans a broad
range and may vary as a function of time, variable step size routines may
lead to considerable inaccuracies in the integration. The usual choice
for integration time step is one fiftieth of the period of the driving
force. We find that, with these choices, we can integrate an entire
array of oscillators at about a time step per second rate on our
existing hardware.
Figure 1:
Development of spatiotemporal complexity on the surface of a sphere. The
coupling 
increases from left to right: 
(left), 
(center),
and 
(right).
The first coupled oscillator system of interest is an array of generalised Duffing oscillators:
Here, x
may be taken as the position of the oscillator, 
is the velocity,
is
the strength of the damping, 
is the coefficient of the linear term,
is the coefficient of the nonlinear
spatial term, and A
and B are the coefficients of the driving terms. is
the strength of the coupling, which is diffusive in nature. The
development of order in this system has been studied for different array
sizes 
and shapes (spheroidal, toroidal or planar), different cell
tiling symmetries (rectangular and hexagonal), various boundary conditions
(free and periodic) and over different ranges of values for the
parameters.
An example of Duffing oscillators tiled with hexagonal symmetry on a sphere is shown in Figure 1. Note how the scale of the spatial patterns changes as the coupling is increased. This shows the development of regions in which the oscillators are strongly correlated with each other. For toroidal arrays, we find that the regions can be characterised by a certain length scale (correlation length) which is related to the group velocity of waves travelling through the array. The correlation lengths of the spatial patterns scale with the coupling as a power law with an exponent near 1/2.
This can be understood by realising that, in the continuum limit, the system is described by coupled partial differential equations representing a dissipative, dispersive, nonlinear wave medium in which the group velocity can be identified as the square root of the coupling. Thus, the size of the correlated regions would appear to be proportional to the group velocity. This is reasonable since the group velocity is the speed with which information and signals travel through the medium. A region is correlated only if its parts can communicate with each other. Thus, the development of regions of a given characteristic size (as opposed to random and highly variable sizes) is a natural result of a system in which the coupling is uniform.
Figure 2:
Spiral wave sequence in the van Capelle-Durrer model. In the 1st
picture at the upper left we show the array shortly after we deliver a
stimulus to a single cell. This models the way a "normal" beat might be
initiated in the heart. In the 2nd and 3rd frames, the wave of excitation
from this stimulus propagates through the medium and, if no other event
were to intervene, would reach the edge of the medium and die out.
However, in the 4th frame we introduce a second stimulus, such as might be
caused by an electrical shock to or pathology in the heart. This second
stimulus cannot propagate to the right, top or bottom because the passage
of the earlier wave has lowered the excitability of those cells. However,
the cells directly to the left of the second stimulus have had time to
recover and, thus, the new wave propagates to the left, as shown in the 5th
and 6th frame. Note that in the 6th frame, the new wavefront has started
to curl around. This is the beginning of a self-sustaining spiral wave.
Frames 7 and 8 detail the evolution of the spiral and the subsequent
formation of a new leftward traveling wave, shown in Frame 9 (which is
essentially the same as Frame 5). Shown here is a 
array of
cells using free boundary conditions.
The second system of interest is based on the van Capelle-Durrer model
[5]
of the human heart. Unlike the Duffing oscillator, this is not a driven
system. The steady state is quiescent and, therefore, the system evolves
only if a stimulus is delivered (such as that given by the pacemaker node
of the heart). However, under suitable conditions, a self-sustaining wave
may be set up in the array. The presence of such a wave
(fibrillation)
interferes with the normal beat-rest-beat-rest- 
rhythm of the heart.
The model is:
where the quantities 
are the electrical potential of each cell and the 
are the excitability of each cell. (Once a cell "beats", it requires a
resting period before it may "beat" again; this waiting period is related
to the time a cell needs to recover its "excitability".) 
and
are
phenomenologically derived functions that describe currents flowing into
and out of the cell, 
is similarly derived to model the observed
excitability and C and T
are constants. The coupling to neighbouring cells
is expressed by an external current 
due to the nearest neighbour cells:
where 
and 
are the parallel and transverse resistances of the heart
tissue, respectively. A typical value of the ratio 
is 2. 
provides
for the introduction of an external electrical stimulus into the system at
any site [6].
We have chosen free boundary conditions for this model. This mimics the fact that the atrium of the heart is surrounded by a non-conducting ring of tissue.
The normal response of this system to a stimulus (such as that from the heart's own pacemaker node) is to propagate a wave of excitation out from the stimulus site. However, if a second stimulus is introduced before the first wave has progressed to the edge of the tissue and died out, the conditions resulting from the passage of the first wavefront may induce the second wavefront to produce self-sustaining spiral waves, as shown in Figure 2. Spiral waves may also be generated by a single wave of excitation encountering an obstacle, such as tissue that has died during a heart attack. The wave splits to move around the dead tissue and eventually meets itself on the opposite side of the obstruction, sometimes leading to the formation of spirals. These results are very similar to spiral waves observed experimentally in in vivo tissue [6, 7].
Spatial order can appear in a wide class of coupled oscillator models. The type and strength of the coupling between oscillators is often the dominant factor governing the appearance of spatial order. It has been shown that an understanding of these systems can be obtained in finite time using only a conventional PC to integrate the model equations. Future work on these systems will concentrate on understanding how to control the spatiotemporal chaos. Some work has already been started but, as yet, a general theory of control in these situations is lacking.
Spatiotemporal Chaos: Visual Simulations of Duffing Oscillators and Heart Tissue
This document was generated using the LaTeX2HTML translator Version 96.1 (Feb 5, 1996) Copyright © 1993, 1994, 1995, 1996, Nikos Drakos, Computer Based Learning Unit, University of Leeds.
The command line arguments were:
latex2html spano.tex.
The translation was initiated by Pam Milliken on Tue Nov 12 11:27:04 EST 1996
