Determination of broken KAM surfaces for particle orbits in toroidal confinement systems

Only through resonance can orbits be significantly modified in an integrable conservative dynamical system. Without resonance the trajectories in phase space occupy Kolmogorov–Arnold–Moser [1] (KAM) surfaces that prohibit diffusion without particle collisons, and in the presence of collisions only relatively slow neoclassical diffusion results. Isolated resonances modify particle distributions by phase mixing in islands and overlapping resonances allow stochastic transport. Often individual perturbations of an integrable system are so small that the system is very far from Chirikov overlap, and thus the destruction of the last KAM surface is not a good paradigm for the examination of diffusive transport. Even if nearby resonances overlap the random phase approximation is not well suited to describe the transport in many cases. It can be necessary to understand the domains of broken KAM surfaces for each perturbation separarately, and to find the resulting transport due to the action of all modes together. In this work we compare three methods of determining the stability of orbits in a Hamiltonian system, the method of frequency analysis [2–4], that of the fast Lyapunov indicator [5] (FLI) and that of phase vector rotation [6, 7].

Although the methods apply to any conservative Hamiltonian system, we consider here the dynamics of charged particles in an axisymetric toroidal magnetic confinement system, perturbed by time dependent magnetohydrodynamic (MHD) modes of a given spectrum. As well as providing an interesting Hamiltonian system of two variables, it is of important practical interest in the pursuit of controlled thermonuclear fusion. In the absense of perturbations the system is integrable, the orbits forming a phase space of invariant tori with two degrees of freedom, corresponding to poloidal and toroidal motion in a topological torus. Typically in experiments several small amplitude modes are present, each with a different frequency and toroidal and poloidal mode numbers [8]. The modes are driven unstable by the free energy of a high energy particle distribution such as fusion produced alpha particles and thus are expected in any magnetic fusion experiment. Although Poincaré sections can be made for each mode individually, such plots are often well below Chirikov threshold [9], individual modes not producing diffusion because of the Nekhoroshev theorem [10] and only the synergistic action of the collection of modes leads to particle transport [6, 7].

Using the guiding center drift approximation a particle orbit in an axisymmetric toroidal system is completely described by the values of the toroidal canonical momentum ${{P}_{\zeta}}$, the energy E and the magnetic moment μ. Particle spatial coordinates are given by ${{\psi}_{p}},\theta ,\zeta$, respectively the poloidal flux coordinate, a label of the unperturbed topologically toroidal flux surfaces, and the poloidal and toroidal angles. The magnetic field is given by

$$\begin{eqnarray}&&\vec{B}=g\nabla \zeta +I\nabla \theta +\delta \nabla {{\psi}_{p}}.\end{eqnarray} \tag{ 1 }$$

where $g({{\psi}_{p}})$, $I({{\psi}_{p}})$ and $\delta ({{\psi}_{p}},\theta )$ are functions determining the form of the equilibrium. The toroidal field strength is given by g, the poloidal field by I, and δ is related to the nonorthogonality of the coordinate system. The guiding center Hamiltonian is

$$\begin{eqnarray}&&H=\rho _{\parallel}^{2}{{B}^{2}}/2+\mu B+\Phi,\end{eqnarray} \tag{ 2 }$$

where ${{\rho}_{\parallel}}={{v}_{\parallel}}/B$ is the normalized parallel velocity, ${{v}_{\parallel}}$ is the particle velocity parallel to the magnetic field, μ is the magnetic moment, and $\Phi$ the electric potential. The field magnitude B and the potential may be functions of ${{\psi}_{p}}$, θ and also ζ if axisymmetry is broken. We consider perturbations with the frequency well below the cyclotron frequency and ignore particle collisions, so μ is conserved and may be considered a constant parameter. The construction of equilibrium fields involves the solution of the Grad-Shafranov equation and is well known. The equilibrium used for most of this paper is a simple circular large aspect ratio equilibrium with the field line helicity $q=\text{d}\theta /\text{d}\zeta$ quadratic in the minor radius and with q on axis of 0.8 and q  =  4 at the plasma edge. The major radius is 100 cm, the minor radius 25 cm and the magnetic field on axis is 20 kG. The details of the equilibrium and the high energy particle distribution determine the mode spectrum. In the following we will use a simple ad hoc mode spectrum, sufficient to display the properties of each method.

Canonical momenta are

$$\begin{eqnarray}&&{{P}_{\zeta}}=g{{\rho}_{\parallel}}-{{\psi}_{p}},\qquad {{P}_{\theta}}=\psi +{{\rho}_{\parallel}}I,\end{eqnarray} \tag{ 3 }$$

where ψ is the toroidal flux, with $\text{d}\psi /\text{d}{{\psi}_{p}}=q({{\psi}_{p}})$, the field line helicity. The equations of motion in Hamiltonian form are

$$\begin{eqnarray}&&\begin{array}{*{20}{l}}\overset{\cdot}{\mathop{\theta}}=\frac{\partial H}{\partial {{P}_{\theta}}},\,\overset{\cdot}{{{P_{\theta}}}}=-\frac{\partial H}{\partial \theta}\\ \overset{\cdot}{\mathop{\zeta}}=\frac{\partial H}{\partial {{P}_{\zeta}}},\,\overset{\cdot}{{{P_{\zeta}}}}=-\frac{\partial H}{\partial \zeta}.\end{array}\end{eqnarray} \tag{ 4 }$$

Hamiltonian equations for advancing particle positions in time, also in the presence of time dependent flute-like perturbations of the form $\delta \vec{B}=\nabla \times \alpha \vec{B}$ with $\vec{B}$ the equilibrium field and $\alpha =\sum\nolimits_{m,n}{{\alpha}_{m,n}}({{\psi}_{p}})\sin \left(n\zeta -m\theta -{{\omega}_{n}}t\right)$ can easily be derived [11]. Since the equilibrium is axisymmetric, the toroidal mode number n is well defined, so the mode frequency ${{\omega}_{n}}$ depends only on n. Each mode normally possesses several poloidal harmonics given by m. The eigenfunctions of the perturbation ${{\alpha}_{mn}}({{\psi}_{p}})$ have the dimensions of length, and are zero at the axis and the plasma edge. We will cite amplitudes giving the maximum value in terms of the major radius R. The guiding center equations including MHD perturbations are realized using a fourth order Runge–Kutta algorithm in the code ORBIT [12]. The units are conveniently defined by the on-axis cyclotron frequency ${{\omega}_{0}}$ (time) and the major radius R (distance). Another characteristic time is the toroidal transit time for a particle of a particular energy and $\mu =0$ at the magnetic axis.

The variables E and ${{P}_{\zeta}}$ are constant in an unperturbed system, but ${{P}_{\zeta}}$ and ${{P}_{\theta}}$ are not action variables so $\text{d}\theta /\text{d}t$ and $\text{d}\zeta /\text{d}t$ are not constant in time. However, in an axisymmetric system, without perturbations, the system is integrable, and orbits close in the poloidal plane in a time ${{T}_{\theta}}$, which depends on the values of E and ${{P}_{\zeta}}$ and thus because of the toroidal precession, frequencies of the system are given by ${{\omega}_{\theta}}=2\pi /{{T}_{\theta}}$, and ${{\omega}_{\zeta}}=\Delta \zeta /{{T}_{\theta}}$, with $\Delta \zeta$ the toroidal motion of the orbit in time ${{T}_{\theta}}$. One can also use the the helicity of a particle orbit in the mode frame, $h=\left(\Delta \zeta -{{\omega}_{n}}t/n\right)/\Delta \theta$, where $\Delta \zeta$ and $\Delta \theta$ are the toroidal and poloidal angles traversed in time t, and the time is taken to be large compared to ${{T}_{\theta}}$, useful because it must be equal to a low order rational for a resonance to occur. An orbit closing upon itself after a finite number of poloidal and toroidal transits is a necessary but not sufficient condition for resonance to occur. In figure 1 is shown a mapping of points in the ${{P}_{\zeta}}$, E plane describing confined orbits into the plane of ${{\omega}_{\theta}}$, ${{\omega}_{\zeta}}$, demonstrating the topological equivalance of these variables to action angle variables within this range. The energy is in units of keV. The frequencies are normalized to the toroidal transit frequency of an on-axis particle with energy 12 KeV and $\mu =0$. In the equilibrium used they have a value of 10 kHz. The introduction of symmetry breaking perturbations destroys the one to one character of the map from the ${{P}_{\zeta}}$, E plane into the plane of ${{\omega}_{\theta}}$, ${{\omega}_{\zeta}}$, normally producing folds along singular surfaces.

Zoom In Zoom Out Reset image size

In section II we discuss kinetic Poincaré plots as a means of examining the destruction of KAM surfaces and give two cases to be examined with the three methods. Section III considers the method of frequency analysis, section IV, the fast Lyapunov indicator, first introduced in the analysis of the stability of the solar system, and section V that of phase vector rotation. Conclusions are given in section VI.

We are interested in the case of the interaction of particles of arbitrary pitch with low amplitude modes of nonzero frequency. It is fairly easy to assess the effect of a particular mode of a single n, ω, but possibly with many poloidal harmonics m, on a particle distribution by examining a Poincaré plot for a particular choice of either co-moving and trapped or of counter-moving particles, which we refer to as a kinetic Poincaré plot to distinguish it from a plot of the magnetic field. Points are plotted in the poloidal cross section whenever $n\zeta -{{\omega}_{n}}t=2\pi k$ with k integer.

The toroidal motion then gives successive Poincaré points in the poloidal cross section ${{\psi}_{p}}$, θ, or better, since ${{P}_{\zeta}}$ and E are constant in the absence of perturbations, the ${{P}_{\zeta}},\theta$ plane or the $E,\theta$ plane. ${{P}_{\zeta}}$, E, and ${{\psi}_{p}}$ are simply related for a mode of a single frequency and toroidal mode number. Individual modes produce islands in the phase space of the particle orbits, which through phase mixing produce local flattening of the particle distribution. In addition, overlap of these islands, the Chirikov criterion, leads to stochastic transport of particles [9, 13]. Such a plot shows the canonical division of orbits into those following good KAM surfaces, isolated islands bounded by separatrices, and stochastic domains.

To obtain a kinetic Poincaré plot the orbits must be initiated with a fixed value of μ and ${{\omega}_{n}}{{P}_{\zeta}}-nE=c$. Here $E-{{\omega}_{n}}{{P}_{\zeta}}/n$ is simply the particle energy in the frame of the rotating mode and is conserved in the presence of this mode, an orbit moving in the ${{P}_{\zeta}}$, E plane only along a line with this slope. A plot of orbits with fixed μ and energy E does not give a coherent plot; it contains intersecting surfaces, since it is really an overlaying of plots with different values of c. Thus a plot can be obtained only for a perturbation consisting of modes with a single value of $\omega /n$, typically meaning only a single frequency and n value, although it may contain many poloidal harmonics m.

Orbits in a general toroidal equilibrium are classified according to whether they are poloidally trapped or passing, whether they circle the magnetic axis, and the direction in which they circle toroidally, dividing the space of ${{P}_{\zeta}},E,\mu$ into distinct domains [11]. In the following we will restrict ourselves for simplicity to a single domain, that of co-moving passing orbits. We take a simple example, consisting of a perturbation with a single harmonic, with m  =  6, n  =  5, a frequency of 10 kHz, and a simple global radial structure for α. We will examine first KAM stability along a single line in the ${{P}_{\zeta}},E,\mu$ volume, with $\mu {{B}_{0}}=5$ keV, ${{\omega}_{n}}{{P}_{\zeta}}-nE=c$ and the energy chosen to be 12 keV at the plasma edge, ${{\psi}_{p}}={{\psi}_{w}}$. This simple example will serve as a test case for the three methods of determining the breaking of KAM surfaces in the space of particle orbits.

The resulting Poincaré maps are shown in figure 2, spanning almost the full domain of confined orbits along this line in the ${{P}_{\zeta}}$, E, μ volume, determined by the value of ${{P}_{\zeta}}$. Principal island chains are seen at the resonances given by $\sin \left({{m}^{\prime}}\theta -n\zeta \right)$ with, starting from the right, ${{m}^{\prime}}/n=6/5,7/5,8/5,9/5,10/5,$ and 11/5. The range of q covered does not include q  =  1, otherwise there would also be visible a ${{m}^{\prime}}/n=5/5$ island on the far right. In between these island chains there are also higher order Fibonacci sequence chains, obtained from modes m/n and ${{m}^{\prime}}/{{n}^{\prime}}$ through $\left(m+{{m}^{\prime}}\right)/\left(n+{{n}^{\prime}}\right)$ not always visible without using higher resolution Poincaré plots. Clearly seen is only the 19/10 resonance between 10/5 and 9/5. The quantity ${{m}^{\prime}}$ is used to distinguish the number of islands in the chain from the poloidal mode number m of the perturbation, m  =  6.

Zoom In Zoom Out Reset image size

This plot illustrates the difficulty in analytically predicting resonance location and strength in this system. Conventional wisdom has it that a perturbation with poloidal and toroidal mode numbers m/n will produce resonances due to the m  =  1 character of the orbital shift due to drift motion with islands at ${{m}^{\prime}}/n=\left(m\pm 1\right)/n$. We see that the actual situation is much more complex, with ${{m}^{\prime}}$ ranging from 6 to 11 even with the perturbation far from Chirikov overlap, and with these resonance islands all practically of the same order in size. Resonance locations can be found by numerical integral methods, [14–16], but while these methods give the location, they do not provide the resonance width and thus the extent of the effect of the mode on transport.

The second case, shown in figure 3 has the same perturbation but a larger amplitude, by a factor of 8, producing some chaotic domains with remnant islands imbedded in them. The perturbations used are smooth and global radially, not localized, so they overlap everywhere. Although some good KAM surfaces still exist, transport in such a system is very complicated. In addition, as has been well demonstrated in simple systems [17], transport even well above the breaking of the last KAM surface involves many long time correlations, and is not well described by a random phase approximation. The transport is characterized by Levy flights, and can be of the form ${{r}^{2}}\sim {{t}^{\beta}}$ with $\beta <1$ subdiffusive or $\beta >1$ superdiffusive [18]. However in the present case the domains of broken KAM surfaces are too small to allow significant long time averages. Probably the only useful measure is the time necessary for the flattening of a distribution in the domain, given by the island rotation frequency.

Zoom In Zoom Out Reset image size

These examples are sufficient for the examination of methods for the determination of broken KAM surfaces, containing resonant islands of various size, regions of unperturbed surfaces, regions with highly perturbed surfaces with no islands, and stochastic domains. A Poincaré plot is useful for a single line in the ${{P}_{\zeta}}$, E, μ volume, but not useful to characterize points in the volume according to whether they belong to good KAM surfaces or are broken, as it depends on visual inspection and gives information along one line at a time.

Both the method of frequency analysis and FLI given in the references had the luxury of knowing the phases of the O-points in the Poincaré plots of the system, and therefore being able to choose initial conditions for orbits so that they intersected the O-points at the resonant values of the action variable. This provides for a very clear demarcation of a resonance and of the Arnold web because it is impossible to miss the center of the resonance. This is not the case in general, where one cannot know these phases, as is clear from the plot of the resonances in figures 2 and 3. Higher order Fibonacci resonances will not in general appear at the same phases as the principal ones. We compensate for this lack of knowledge where possible by initiating orbits with different phase values, to maximize the probability of intersecting the resonance.

The works of Laskar [2–4] discuss various averaging means of obtaining frequencies associated with orbits. Since we are interested in resonances, the orbit helicity h is a quantity which is heuristically useful. Orbits are initiated along a single line in the ${{P}_{\zeta}},E,\mu$ volume, with $\mu {{B}_{0}}=5$ keV, ${{\omega}_{n}}{{P}_{\zeta}}-nE=c$ and the energy chosen to be 12 keV at the plasma edge, ${{\psi}_{p}}={{\psi}_{w}}$, ie the same set of orbits used for the kinetic Poincaré plots shown above. If an orbit is followed for a sufficient length of time the bounce motion of a particle trapped in an island is negligible, and the observed helicity of the orbit is simply the helicity of the resonance. The interior of an island thus produces a range of parameter values in which h is constant. At the edges of this domain there is a jump in the value of h and hence an easily observed singularity in the second derivative. In our case μ is fixed and the initial value of ${{P}_{\zeta}}$ characterizes the orbit. We typically used a total of one thousand orbits, spanning the range of confined orbits in the system used in figure 2, and a time of one thousand transits.

The results of this method depend on the phase between the orbit trajectory and the mode phase. In figure 4 is shown how the results vary for a single resonance, with plots made with different phases between orbit and mode. If the trajectory intersects the mode at an elliptic point the range of constant h is as broad as the island width, but if it intersects the island at a hyperbolic point it has zero width. For a given phase the plot of helicity is approximately symmetric about the rational surface. The phase dependent island width shown by this method complicates the analysis if the mode phases are not known or if they vary from resonance to resonance, giving widths of regions with constant h that do not represent the maximum island width for all resonances.

Zoom In Zoom Out Reset image size

III.A. Case I

In figures 5 and 6 are shown the results for the same orbits shown in figure 2. Following Laskar, we varied the number of orbits used, and hence the spacing in ${{P}_{\zeta}}$ as well as the the time step in the integration process to obtain good results. Flat domains of h are seen for the major island chains, and the numerical second derivative of the helicity [$\delta \delta h=h(k+1)+h(k-1)-2h(k)$] clearly indicate the edges of major island chains as well as three of the first order Fibonacci chains. The red lines in figure 5 are at helicities of $n/{{m}^{\prime}}=5/6,5/7,5/8,5/9,5/10,5/11$, in excellent agreement with the flat domains of the function h. In the domains of good unperturbed KAM surfaces $h\left({{P}_{\zeta}}\right)$ is seen to be very smooth. The only drawback of this method is that the flat regions of h do not give the full island width unless the orbit happens to intersect the island at an elliptic point. We have managed to choose a phase so that the islands are well indicated, so the widths are apparent. Note that the same phase must be used for all orbits, using random values of phase for each point would make h practically a random function of ${{P}_{\zeta}}$.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

III.B. Case II

In figures 7 and 8 are shown the results for the same orbits shown in figure 3, with $\alpha =4\times {{10}^{-4}}R$. In this plot h is flat at most of the remaining imbedded islands and also at higher order Fibonacci islands in between them. Island 8/5 is missing although present in the Poincaré plot, only an extended chaotic domain is seen. The numerical second derivative of the helicity [$\delta \delta h=h(k+1)+h(k-1)-2h(k)$] shows the remaining large islands as well as Fibonacci islands in between each of the principal resonances, with large resonances at 10/13 and 10/19. There are three clear domains of chaos visible, near ${{P}_{\zeta}}/{{\psi}_{w}}=.05$, from 0.1 to 0.2, and from 0.24 to 0.27, surrounding the remaining isolated islands.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The fast Lyapunov indicator (FLI) is another means to determine the destabilization of orbits [5]. Consider two very nearby orbits and examine the distance d between them in the ${{P}_{\zeta}},\theta$ plane. The normal Lyapunov is the limit of ln(d)/t for large time t. The FLI is the value $\lambda =ln(d)$ at fixed time t, which keeps track of the topological differences between resonant regular motion and KAM tori. Using the usual Lyapunov number, for $t\to \infty$ both domains give zero, but at finite time the linear dependence of d on time for points at good KAM surfaces distinguishes these points from island interiors, where d remains relatively constant. To account for the fact that the phase location of the islands is unknown, several different orbit pairs were used for each value of ${{P}_{\zeta}}$ with different phase with respect to the mode and the minimum value λ taken for each ${{P}_{\zeta}}$. Thus the value of λ is given with the inital conditions of the orbit closest in phase to being at an island elliptic point. The time t was chosen to optimize the resolution of the calculation. The pair separation was chosen to be $\text{d}{{\psi}_{p}}=2\times {{10}^{-4}}{{\psi}_{w}}$, limiting the detection of islands to those larger than this.

IV.A. Case I

Shown in figure 9 is the fast Lyapunov indicator for the perturbation shown in figure 2. Five orbit pairs were initiated at each value of ${{P}_{\zeta}}$, and the minimum value of the separation betwen the points taken for λ. We find the distance between the points d to be a more sensative indicator than the logarithm of d. The islands are seen as regions where the indicator is near zero. Domains with good KAM surfaces show linear separation of the orbit pairs in time and thus nonzero λ. This method roughly indicates the width of the islands, but not exactly because there is also some spread of orbit pairs located inside an island due to relative rotation rates about the island elliptic point.

Zoom In Zoom Out Reset image size

IV.B. Case II

Shown in figure 10 is the fast Lyapunov indicator for the perturbation shown in figure 3. This indicator gives only a very qualitative indication of the stochastic domains and remnant islands. One can roughly make out regions where the remnant islands exist, but with very limited accuracy.

Zoom In Zoom Out Reset image size

Another method for numerically determining the existence of or the destruction of good KAM surfaces can be obtained using the method of phase vector rotation [6, 7]. Consider following two orbits located nearby one another. Examine a Poincaré section in ${{P}_{\zeta}},\theta$ and define the angle χ to give the orientation of the vector joining them in this plane. If good KAM surfaces exist χ can change by at most an angle of π, due to their relative velocity in the angular coordinate.

However two orbits within an island rotate around one another with χ increasing with the rotation about the island, also referred to as the bounce frequency of a particle trapped in the wave, which increases with the size of the island. The rate of change of χ is a function of distance from the island O-point, dropping to zero at the separatrix. Pairs of orbits are illustrated in figure 11, showing vectors between the nearby points in the ${{P}_{\zeta}},\theta$ plane on good KAM surfaces and in a resonance.

Zoom In Zoom Out Reset image size

Thus we determine the nonexistence of good KAM surfaces by examining nearby pairs of orbits, looking for phase vector rotation χ exceeding π. Because local perturbation of orbits can cause χ to excede π by some small amount, whereas χ continues to rotate in an unbounded manner in an island, we record the time at which $|\chi |\,=2\pi$, thus giving the rotation frequency ${{\omega}_{\chi}}$ as a function of the initial position of the orbit. The pair separation was chosen to be $\text{d}{{\psi}_{p}}=2\times {{10}^{-4}}{{\psi}_{w}}$, limiting the detection of islands to those larger than this, and the time used was 500 toroidal transit times, making the minimum observable frequency to be ${{\omega}_{\chi}}=2\times {{10}^{-3}}$, normalized to the toroidal transit frequency. The results depend on the phase relation of the orbit with respect to the mode, but using this method one can simply launch a few pairs of orbits at a given value of ${{P}_{\zeta}}$ with different phases and select the maximum value of ${{\omega}_{\chi}}$, corresponding to the pair which is closest in phase to the island elliptic point.

V.A. Case I

Results are shown in figure 12 for a perturbation amplitude of $\alpha =5\times {{10}^{-5}}R$. This method clearly indicates the presence of phase space islands. In addition to the major resonances seen in figure 2 the first order Fibonacci resonances, 13/10, 15/10, 17/10, and 19/10 are all visible. The magnitude of the phase vector rotation ${{\omega}_{\chi}}$ gives the time scale for phase mixing of orbits trapped inside the island. Note that for larger islands one even obtains the internal profile of the rotation frequency. To achieve this several orbit pairs are necessary, with different phases with respect to the mode, as was used for the FLI. Note that the phase vector rotation and the FLI are obtained with the same calculation. 2000 orbits, 500 transits.

Zoom In Zoom Out Reset image size

V.B. Case II

Results are shown in figure 13 for the larger perturbation. The major remnant islands are seen, and it is also possible to see overlap of these chains by noting that there is no region in which ${{\omega}_{\chi}}=0$ in between some chains, indicating that the island is surrounded by a stochastic sea, with no intervening KAM surfaces. Note that the frequencies for phase mixing are three times as rapid as in the lower amplitude case shown in figure 12, scaling as the square root of the perturbation amplitude. It is interesting that the stochastic domains do not produce more rapid phase vector rotation than the coherent rotation present in the remnant island centers. Taking the maximum values of ${{\omega}_{\chi}}$ over the mode-orbit phase relations always clearly shows the centers. The edge of the remnant islands can even be inferred by the minimum values of ${{\omega}_{\chi}}$ surounding it, where overlap with another island chain and chaos intervenes. Again for these calculations 1000 orbits were followed for 500 transits.

Zoom In Zoom Out Reset image size

The Arnold web consists of resonances in the plane of the action variables. It is an open set, dense but with small measure if the perturbations are small. In many cases the resonant lines form complicated web like patterns in this plane. Arnold diffusion is diffusive motion along a resonance line this plane.

Examples of Arnold web construction for the methods of frequency determination and FLI are given in the cited references, so we will not reproduce these methods here. The analysis of the single line in the ${{P}_{\zeta}},E$ plane given in the preceding sections is sufficient to judge the capabilities and shortcomings of the methods without examining the whole plane of the Arnold web. In the system under study the Arnold web is simpler than in many Hamiltonian systems, consisting almost entirely of nonintersecting lines of resonances, occasionally broken in intervals. However there are often several perturbations present of different frequencies, which thus cannot be observed in a single Poincaré plot, and the action of these modes together produces significant chaotic transport. In the present work the resonances are always observed to be nearly vertical in the plane of ${{P}_{\zeta}}$, E, almost parallel to the bounding surfaces, given on the left by contact with the plasma boundary, and on the right by the magnetic axis. No resonance lines with nearly constant energy E are observed. Thus in the particular system of particle orbits in a toroidal confinement device, the Arnold web does not lead to the possibility of significant transport. It could provide some movement in the energy variable, but nothing leading to paticle loss or profile flattening.

Resonance is the only means by which modes can modify particle distributions, through Landau mixing and stochastic motion due to the overlap of resonances. Thus the determination of the location and extent of resonances due to perturbations is important for understanding induced transport. Three methods of determining the destruction of KAM surfaces were compared, that of frequency analysis, that of the fast Lyapunov incicator, and that of phase vector rotation. The method of phase vector rotation appears to be a more efficient and reliable means of finding resonances. It contains more information, giving the extent of island structures, even when immersed in a stochastic sea, and also gives a measure of the mixing time within the islands. The mixing time, given in figures 12 and 13 by $2\pi /{{\omega}_{\chi}}$, directly gives the time scale for profile flattening within the domain of broken KAM surfaces, and thus provides transport information even when other quantitative measures such as diffusion rates are lacking.

The author acknowledges useful discussions with J Laskar. This work was partially supported by the U.S. Department of Energy Grant DE-AC02-09CH11466.