For each of the multiscale turbulence, reconnection, turbulent reconnection, and reconnection-driven turbulence projects funded by the Center, we are proceeding in three phases, partly in parallel. Very soon, we will begin experimental tests of the predictions from existing codes, at first in relatively simple contexts (reconnection and Alfvénic dynamics in LAPD, reconnection dynamics in VTF, ETG turbulence in DIII-D, NSTX and C-Mod) and later in more complex scenarios (NTM and sawtooth fluctuations, barrier formation and structure).

At the same time, a two-phase implementation of the multiscale algorithms will be undertaken. The first, easier stage will be the ``legacy'' stage. Here, wrappers around existing timesteppers will be used to perform projective integration, compute approximate slow manifolds, and transform the timesteppers to numerical continuation/bifurcation codes at the same level of description as the timestepper itself. This will validate the ability to transform (through a computational wrapper) GS2 and p3d into codes capable of performing new tasks (continuation/bifurcation, approximation of slow manifolds). We now describe an explicit plan of attack for exploring reconnection using the p3d code as the kinetic integrator with the coarse representation given by the set of MHD fields. For simplicity, we will start with a simple slab Harris current sheet in a 2-D system with a small initial magnetic perturbation to initiate magnetic reconnection. Since this system has been previously explored, its behavior is fairly well understood and it therefore will serve as an ideal test bed for the ``equation-free'' projective integration approach. Initially we will utilize a uniform grid although the BATS-R-US grid generation algorithm will be implemented as success with the simple system has been demonstrated. The MHD field variables will be used to generate a particle representation using p3d, which already has this initialization capability. The kinetic equations will be advanced through a specified number of time steps as discussed earlier and the results projected onto the coarse MHD variable set. Care must be taken to ensure that system has relaxed to the slow manifold. This can be checked by examining the convergence of several statistically independent kinetic representations of a given set of coarse data. We can then proceed with the projective integration along the ``slow manifold'' and repeat the process. The sensitivity of the results to the number of kinetic time steps will be tested and the desireability of ensemble averaging the kinetic representation will be explored. Success with this simple problem will demonstrate the viability of the ``equation-free'' approach and will already represent a major step forward for plasma physics.

In parallel with the equation-free approach we will also pursue a hybrid MHD/kinetic approach. The limitation of the MHD model is in the representation of the electric field, which typically comes from the one-fluid Ohm's law, and the heat fluxes of electrons and ions. The MHD model is a poor representation of these functions, especially at small spatial scales. Instead of advancing the MHD fields in the ``equation-free'' approach, we will simply evaluate the electric fields and heat fluxes from the kinetic model using p3d and advance the MHD equations in time using BATS-R-US. We will explore optimum time intervals for updating the electric field and heat fluxes. We note that although we are using BATS-R-US for this purpose, this hybrid approach, if it is successful, can readily be implemented with the standard MHD fusion codes, M3D and NIMROD.

Success with either model will allow us to move forward with computational advancements required to address the scientific foci of the Center, transport barrier dynamics, the sawtooth crash and neoclassical island growth. The multigrid algorithm in p3d will be altered to allow p3d to be run in the AMR environment of BARS-R-US. The FFT algorithm in GS2 will be similarly altered. The projective integration technique will then be updated so that traditional MHD time stepping can be run in parallel in spatial regions where a pure MHD rather than a kinetic representation can be utilized to advance the system in time.

As experience is gained from this stage, we will proceed to the second stage: the coarse, averaged computer-assisted analysis of the kinetic equations. Extensive numerical experimentation will be required here to establish the form of the coarse-graining of the kinetic equations that will lead to the most efficient convergence onto the ``slow manifold'' for efficient projective integration. Specifically, the MHD coarse grain variables may not be the optimum representation of the coarse data. Heat fluxes, for example, or temperature anisotropies might perhaps be an appropriate addition to the coarse data set. In extensive discussions with Dr. Kevrekidis we have identified a well-defined computational procedure for optimizing the coarse grain variables. Similar procedures have been successful in related applications, including the dynamics of bubbles in liquids. In the case of magnetic reconnection, for example, a quasi-steady X-line develops in the kinetic description (e.g., with p3d). This quasi-steady configuration will be perturbed and its relaxation measured. By carrying out this process with an ensemble of statistically independent perturbations, we can understand the spectrum of damped eigenvalues and eigenfunctions around the quasistationary X-line configuration. These eigenvalues control the rate of convergence of the fine scale representation to the ``slow manifold'' in projective integration. The procedure is to identify the most weakly damped eigenvalues and corresponding eigenfunctions so that these weakly damped eigenfunctions can be eliminated when creating the fine scale representations for projective integration. This procedure will allow us to create the optimimum coarse grain representation for the ``slow manifold'' projective integration. The expectation is that significant differences between the MHD coarse grain description will be discovered, producing an exciting advance for our field. The final selection of appropriate lifting and restriction operators will follow and lead to the optimization of the equation-free computational environment. This is a vital stage — once a good set of coarse observables is established, the remaining tasks are quite intensive, but in a sense clear and direct.

In the course of carrying out this research plan, we are thus extending the algorithmic advances (gyrokinetics, particle-in-cell) that have already made it possible to perform realistic simulations of nonlinear plasma dynamics, to incorporate the latest scientific advances in multiscale algorithm design. Equally importantly, we are pursuing detailed experimental tests of the simulation predictions along the way. These advances will contribute to a greater understanding of the key plasma physics problems that will be faced by a burning plasma experiment.