Self-organization during spherical torus formation by flux rope merging in the Mega Ampere Spherical Tokamak

Magnetic reconnection is a process for re-structuring magnetic fields and dissipating free magnetic energy which is of importance in a wide range of space and astrophysical plasmas, as well as magnetically-confined fusion plasmas in the laboratory [1–3]. Magnetic flux ropes are a basic building block of magnetic fields, and the merging of flux ropes through magnetic reconnection is a frequently-occurring process. Considerable understanding of flux-rope merging has been attained through purpose-built laboratory experiments e.g. [4–10].

Here, we report on the modelling of flux-rope merging in the Mega Ampere Spherical Tokamak (MAST), which is in a parameter regime more closely resembling astrophysical plasmas. MAST [11–13] is a tight-aspect ratio toroidal device, in which plasmas are created with typical major and minor radii R = 0.95 m, a = 0.60 m, plasma current I_p = 400–900 kA, toroidal field at the magnetic axis B_T = 0.40–0.58 T, and peak electron density and temperature n_e0 = 3 × 10¹⁹ m⁻³ and T_e0 = 1 keV. The MAST device has demonstrated the very promising potential of the spherical tokamak as a fusion power plant or Component Test Facility, especially in terms of confinement and high β operation, as well as providing insight into a range of physical processes relevant to conventional tokamaks. Several alternative plasma start-up methods can be used to reach these final parameter values, one of which is the merging-compression technique, whereby two flux ropes with parallel toroidal current move towards each other and then merge, forming a plasma torus with a single set of closed magnetic flux surfaces [2]. In addition to being an attractive scheme for plasma start-up and current drive without the use of a central solenoid, merging-compression in MAST also provides an opportunity to investigate the basic physics of flux-rope merging and magnetic reconnection in conditions of higher Lundquist number and lower plasma beta than those of other laboratory experiments, and in parameter regimes more closely resembling those of at least some astrophysical plasmas, such as the flaring solar corona. Moreover MAST is equipped with a comprehensive set of diagnostics, including Thomson scattering lasers (see right-hand frame of figure 1) that make it possible to obtain electron temperature and density profiles with sufficiently high time resolution to study the dynamics of reconnection during the merging process.

Zoom In Zoom Out Reset image size

The merging-compression formation scheme was pioneered on MAST's smaller predecessor START [14]. In this scenario, described more fully by [15–17], two toroidal flux ropes with parallel toroidal currents are produced around in-vessel poloidal field coils (the P3 coils shown in figure 1). Decreasing the current in the coils induces plasma current in the flux ropes and causes them to separate from the coils, mutually attracting due to their parallel currents; they eventually merge through reconnection of the poloidal field at the midplane (see figure 1). The process presently generates spherical tokamak plasmas with a current of up to 0.5 MA, which are rapidly heated (presumably, mainly by processes associated with the reconnection) to temperatures of up to 1 keV.

In section 2, we present an overview of recent numerical simulations of merging-compression in MAST using both resistive magnetohydrodynamic (MHD) and Hall-MHD equations, with the HiFi modelling framework [18, 19]. These are performed in 2D, both in Cartesian geometry and in tight-aspect ratio toroidal geometry; for more details, see [17]. Then, in section 3, we describe initial results from an analytical model explaining the self-organization of the flux ropes in terms of a helicity-conserving relaxation to a minimum energy state. Conclusions are presented, and the relevance for astrophysical plasmas is discussed, in section 4.

Here, we describe results of simulations using both single-fluid MHD and Hall-MHD. The Hall terms in Ohm's law are included here because, for typical MAST parameters (see below), the ion skin-depth is comparable with the device size (normalized d_i = 0.145) and is much larger than the width of the Sweet–Parker current sheet predicted by resistive MHD. A hyper-resistive term [1] representing anomalous electron viscosity due to effects of micro-turbulence and 3D instabilities is also included, which sets a dissipation scale for whistler and kinetic Alfvén waves and provides a parallel electric field to allow reconnection in the 2D simulations. Other terms, such as electron inertia and off-diagonal pressure tensor terms, are likely to be less significant, but will be included in future modelling. We solve the nonlinear compressible Hall-MHD equations in normalized form [1]. The normalization is with respect to typical MAST start-up plasma values: density n₀ = 5 × 10¹⁸ m⁻³, length-scale L₀ = 1 m, magnetic field B₀ = 0.5 T, and (equal) ion and electron temperatures T₀ = 10 eV. Velocities are normalized by the Alfvén speed v₀ = B₀ (μ₀n₀m_i)^−1/2. The full set of equations is

$$\begin{eqnarray} &&\fl \partial_t n + {\boldsymbol \nabla} \cdot (n{\bit v}_{\rm i}) = 0, \end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray} &&\fl \partial_t (n {\bit v}_{\rm i}) + {\boldsymbol \nabla} \cdot (n {\bit v}_{\rm i} {\bit v}_{\rm i} + p\mathbb{I} + {\boldsymbol \pi}_{\rm i}) = {\bit j} \times {\bit B}, \end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray} &&\fl {\bit E} = -\partial_t {\bit A} = -{\bit v}_{\rm e}\times {\bit B} - \frac{d_{\rm i}}{n}{\boldsymbol \nabla}p_{\rm e} + \eta {\bit j} - \eta_{\rm H} \nabla^2 {\bit j}, \end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray} &&\fl \partial_t {\bit B} = - {\boldsymbol \nabla} \times {\bit E} \end{eqnarray} \tag{ 4 }$$

$$\begin{eqnarray} &&\fl (\gamma - 1)^{-1}\left[\partial_t p + {\bit v}_{\rm i} \cdot {\boldsymbol \nabla}p + \gamma p{\boldsymbol \nabla}\cdot {\bit v}_{\rm i}\right]\nonumber\\ &&= \eta j^2 + \eta_{\rm H} ({\boldsymbol \nabla}{\bit j})^2 -{\boldsymbol \pi}_{\rm i}:{\boldsymbol \nabla}{\bit v}_{\rm i} - {\boldsymbol \nabla} \cdot {\bit q}, \end{eqnarray} \tag{ 5 }$$

where n is the density, j = ∇ × B the current density, v_i the ion velocity (with electron velocity v_e = v_i − d_i j/n), B the magnetic field, p = p_i + p_e the total (sum of the ion and electron) thermal pressure, E the electric field and A the magnetic vector potential. The ion stress tensor is ${\boldsymbol \pi}_{\rm i} = -\mu_{\rm i}({\boldsymbol \nabla}{\bit v}_{\rm i} + {\boldsymbol \nabla}{\bit v}_{\rm i}^{\rm T})$, and the heat-flux vector q has anisotropic form ${\bit q} = -\kappa^\parallel_{\rm e} {\boldsymbol \nabla}_\parallel T - \kappa^\perp_{\rm i} {\boldsymbol \nabla}_\perp T$ where ${\boldsymbol \nabla}_\parallel = {\hat{\bit b}}({\hat{\bit b}}\cdot {\boldsymbol \nabla})$. The coefficients are: normalized ion skin-depth, $d_{\rm i}=c(n_0 e^2/\epsilon_0m_{\rm i})^{-1/2}L_0^{-1}$, resistivity in terms of the parallel Spitzer value η = (μ₀v₀L₀)⁻¹η_Sp,∥, hyper-resistivity η_H, viscosity $\mu_{\rm i} =(\mu_0 v_0 L_0)^{-1}\mu_{\rm i}^\parallel$, and the normalized parallel electron, $\kappa_{\rm e}^\parallel$, and perpendicular ion, $\kappa_{\rm i}^\perp$, heat conductivities; here, we set the normalized values to be $\kappa_{\rm e}^\parallel=10^{-1}$ and $\kappa_{\rm i}^\perp =10^{-7}$. The final term in equation (3) is the hyper-resistive diffusion. On the basis of Braginskii values based on pre-merging MAST temperatures, η = 10⁻⁵ and μ_i = 10⁻³; but the values of η, μ_i and η_H are also sometimes varied in order to study the scaling effects of collisions on the merging. With these standard values, the Lundquist number S = 10⁵ and the magnetic Prandtl number P_m = 100.

The simulations use a grid with N_R = 180, N_z = 360 finite elements in the R, z directions, with a polynomial of order 4 in each element, giving an effective resolution of 720 × 1440. The spacing is non-uniform, with a finer grid near the current sheet location around z = 0. Convergence studies have been undertaken.

The initial conditions are taken as two localized distributions of toroidal current representing the plasma when the current rings have detached from the in-vessel coils (which are not modelled). The toroidal field within the flux ropes is calculated so as to give local force-balance, but an attractive force between the parallel currents remains and triggers the merging process.

Figure 2 shows a sequence of fields for a resistive MHD simulation (d_i = 0, η_H = 0; η = 10⁻⁵, μ = 10⁻³), in Cartesian geometry. It can be seen that as the flux ropes approach, a current sheet initially forms, which becomes the site of reconnection. The pile-up of the poloidal field (B_R), associated with a layer of reverse toroidal current (see second panel of figure 2), opposes the merging and causes the reconnection to stall. This process continues, leading to oscillations in the reconnection rate or 'sloshing' [20, 21] which can be seen in some traces on figure 3. The 'sloshing' is known to occur for sufficiently small values of η, and can be interpreted as being due to the increase in induced current between approaching magnetic islands, on a timescale more rapid than resistive decay, leading to a repulsive force [22]; thus, the frequency scales with the Alfven time. Nevertheless, full merging of the flux ropes is eventually achieved, as shown in the final panel of figure 2. For the resistive MHD case, analysis of the scalings with μ and η [17] shows that the reconnection rate is consistent with the visco-resistive Sweet–Parker reconnection [23, 24].

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The time dependence of the reconnection rate for Hall-MHD is shown in figure 3, along with a comparable single-fluid MHD case in which the hyper-resistivity has been included, in order to distinguish between the effects of finite ion skin-depth and hyper-resistive diffusion. Note that the reconnection rate for the single-fluid case, with hyper-resistivity—performed largely as a numerical test of the effects of hyper-resistivity—is very close to that for resistive MHD, since the breaking of the frozen-in condition is largely provided by the resistivity η at this low value of η_H = 10⁻¹⁰). Evidently, inclusion of the Hall terms in Ohm's law can significantly increase the reconnection rate. As the hyper-resistivity decreases from 10⁻⁶ to 10⁻⁷, both peak and average reconnection rates decrease; as η_H approaches 10⁻¹⁰, the peak reconnection rate becomes larger (see below). Note that the contribution to the reconnection electric field by the hyper-resistivity dominates over the resistive contribution for the Hall-MHD simulations.

The structure of the reconnection region also changes significantly as the Hall terms are included: the symmetry is broken, leading to a tilted current sheet and tilted ion outflow jets [17]. The effects of varying collisionality are explored by changing the hyper-resistivity η_H [17], showing a transition between various reconnection regimes. For higher collisionality, η_H = 10⁻⁶, a broad current sheet forms which does not fragment. Reducing the collisionality leads to break-up of the sheet due to tearing-type instabilities. At low values of η_H, the thin current sheet becomes unstable to secondary island formation, and a series of small islands (regions of closed magnetic flux) are generated, separated by additional X-points; see figure 4. At the lowest values simulated, η_H = 10⁻¹⁰, the current sheet spreads asymmetrically across the separatrix between open and closed flux, and the outflow region opens, leading to fast reconnection.

Zoom In Zoom Out Reset image size

Simulations have also been performed in tight-aspect-ratio toroidal geometry, with a confining vertical field in addition to the dominant toroidal field (with 1/r dependence), more closely representing MAST. The flux-rope merging proceeds as before, leading to the formation of a single 'spherical tokamak' state with closed flux surfaces and a realistic, monotonically-increasing profile of the safety factor q, with q > 1 everywhere. There are some differences in the dynamics from the Cartesian case, especially due to the inboard-outboard asymmetry; in particular, the islands which form in the current sheet are more easily ejected. For further details, see [17, 25].

The simulations also predict the plasma density profiles. It has been shown that, in the case of Hall-MHD simulations in toroidal geometry, a double-peaked radial density structure in the midplane is predicted, which agrees with experimental measurements using Thomson scattering [15, 17]. This is discussed in more detail in [17], where it is shown that explaining the observed double-peaked density structure requires both two-fluid and toroidal effects.

The concept of helicity-conserving relaxation to a minimum energy state was pioneered for reverse field pinch plasmas by Taylor [26] and has subsequently been applied with great success in a range of a laboratory plasmas (see review [27]) as well as astrophysical plasmas such as the solar corona (e.g. [28, 29]). The basic idea is that a plasma subject to some initial disturbance or instability relaxes towards the state of lowest magnetic energy, with the constraint that the total magnetic helicity is conserved. The latter constraint is appropriate in a highly-conducting plasma in which dissipation is mostly due to reconnection within thin current layers; the basic idea is that reconnection enhances energy dissipation, but not helicity dissipation (whilst both slowly dissipate due to global resistive diffusion) [27, 30]. In this case, the relaxed state is a linear force-free field,

$$\begin{equation} \nabla\times{\bit B} = \mu{\bit B}, \end{equation} \tag{ 6 }$$

where μ = μ₀B · j/B² is spatially constant.

We propose, following an earlier suggestion for START (M Bevir and C G Gimblett, private communication), that flux-rope merging in MAST can be modelled as a Taylor relaxation. Initial results of such a model are presented here; for simplicity, a Cartesian geometry is assumed, with a rectangular cross-section infinite-aspect ratio tokamak ('toroidal' invariant direction z). The model could be easily adapted to a tight-aspect-ratio system in cylindrical geometry, following [31]; but, as demonstrated by the numerical simulations in section 2, this should not have a major effect on the results. Furthermore, the aim of the relaxation model is to draw out some essential underlying physics, and this would not be greatly enhanced by a more accurate choice of geometry.

3.1. Relaxation model

The initial state, shown in the left-hand panel of figure 5, consists of two adjacent, twisted flux ropes, with square outer cross-sections given by x = 0, a (for both flux ropes) and y = 0, a for one flux-rope and y = a, 2a for the other; hence, the flux ropes touch along the surface y = a. We assume each initial flux rope is force-free. The field then undergoes helicity-conserving relaxation to a single flux rope contained within the same overall volume (0 ⩽ x ⩽ a; 0 ⩽ y ⩽ 2a). The relaxed state is described by (6). For ease of calculation, we assume that the initial flux ropes are individually linear force-free fields, also described by (6); there is no particular reason to assume such a μ profile, but as we have no detailed knowledge of the current distribution, this is not unreasonable (note therefore that the initial current profile here is somewhat different from that in section 2). However, the initial state is not a minimum energy state, because it contains a current sheet at the boundary between the two flux ropes (at y = a), due to the discontinuity in the poloidal field B_x which reverses direction. Hence the model initial state represents the fields at the point in time at which the flux ropes have been brought together by the attractive force, but have not yet commenced reconnection (somewhat later than the initial state of the simulations). This is a force-free equilibrium, since the magnetic pressure balances (by symmetry) across the interface y = 0; but is not a relaxed (constant-μ) state because μ has a negative value, with delta-function dependence, in the current sheet.

Zoom In Zoom Out Reset image size

The fields for the general case of a rectangular boundary a × b are found using a flux function ψ, where B = (∂ψ/∂y, −∂ψ/∂x, μψ), so (6) leads to

$$\begin{equation} \nabla^2\psi+\mu^2\psi=0. \end{equation} \tag{ 7 }$$

The boundary condition is ψ = ψ_b = constant on the boundary, where ψ_b, which is effectively a normalization constant for the fields, must be non-zero. With the choice B_z = μψ made above, with the allowable arbitrary constant set to zero, this is equivalent to requiring non-zero toroidal field on the boundary. (For ψ_b = 0, the only non-trivial solutions to (7) are spheromak-like solutions for discrete eigenvalues of μ, which are not relevant to spherical tokamaks). Following [32] and [31], the solution is expressed as

$$\begin{equation} \psi=\psi_{\rm b}\bigg[1 + \sum\limits_{m,n} {a_{nm}\sin\left(\frac{m\pi x}{a}\right)\sin\left(\frac{n\pi y}{b}\right)}\bigg] \end{equation} \tag{ 8 }$$

which automatically satisfies the required boundary condition. Substitution into (7) gives

$$\begin{eqnarray} &&\fl \sum\limits_{m,n}\left(-\frac{m^2\pi^2}{a^2}-\frac{n^2\pi^2}{b^2}+\mu^2\right) a_{nm}\sin\left(\frac{m\pi x}{a}\right)\nonumber\\ &&\times \sin\left(\frac{n\pi y}{b}\right)=-\mu^2; \end{eqnarray} \tag{ 9 }$$

using the orthogonality of the eigenfunctions, sin (mπ x/a) sin (nπy/b), within the summation then yields

$$\begin{eqnarray} &&\fl a_{mn}=\frac{16\mu^2}{mn\pi^2[\pi^2(m^2/a^2+n^2/b^2)-\mu^2]}, \quad {\rm for}~m,n~{\rm odd};\nonumber\\ \end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray} &&\fl a_{mn}=0 , \tqs {\rm for}~m~or~n~{\rm even;} \end{eqnarray} \tag{ 11 }$$

see also [33]. Results are presented in dimensionless form with lengths scaled with respect to a. We also normalize so that the toroidal flux Φ_t ≡ ∫_AB_z dx dy = 1 (thereby ensuring conservation of flux during the relaxation process described below), so that ψ_b is a function of μ. It should be noted that, with this normalization, ψ_b, and all related quantities such as energy and helicity, diverge as μ approaches the first eigenvalue of (7) [32, 34]; however, since we are interested here in tokamak-like magnetic configurations with strong toroidal fields, this occurs outside the relevant parameter range.

The fields from (8), (10) are used both to represent the individual flux ropes in the initial state (with b = a) and the single relaxed state (with b = 2a). The final value of μ, μ_f, is determined from the initial value by the constraint that helicity, K, and total toroidal flux are conserved (where the dimensionless total toroidal flux initially is 2, because of the two flux ropes). A root-finding process is used to find μ_f so as to conserve $K/\Phi_{\rm t}^2$. The helicity is given by

$$\begin{equation} K=\int_{V} {\bit A} \cdot {\bit B} \mathrm{d}V, \end{equation} \tag{ 12 }$$

where we fix A_z = 0 on the boundary so that the gauge-correction term arising from the flux through the torus [35] vanishes. It can be shown that this is related to the total magnetic energy, W, by

$$\begin{equation} K=2\mu_{0}W/\mu-\psi_{\rm b}\Phi_{\rm t}, \end{equation} \tag{ 13 }$$

where the energy (per unit length) is obtained by integration of the fields, as

$$\begin{eqnarray} &&\fl W(\mu;a,b) = \frac{ab\psi_{\rm b}^2}{2\mu_0}\nonumber\\ &&\fl \times\left[\mu^2+ \sum\limits_{m,n {\rm ~odd}}\left[\frac{a_{mn}^2}{4}\left(\frac{m^2\pi^2}{a^2}+\frac{n^2\pi^2}{b^2}+\mu^2\right) +\frac{8a_{mn}\mu^2}{\pi^2mn}\right]\right],\nonumber\\ \end{eqnarray} \tag{ 14 }$$

where ψ_b(μ, a, b) is given by the normalization Φ_t = 1 (and μ₀ = 1 in dimensionless units).

3.2. Relaxation results

A typical relaxed state is shown in the right-hand panel of figure 5. Note that this represents the post-merging state, and is indeed very similar to the final states of the numerical simulations as shown in figure 2.

The energy released during relaxation is calculated as ΔW = W_i−W_f where the initial energy W_i = 2W(μ_i;a, a) and the final energy is W_f = 4W(μ_f;a, 2a) (the factors 2 and 4 are to account for, respectively, two flux ropes with Φ_t = 1 initially, and a single flux rope with Φ_t = 2 finally. In figure 6, the final value of μ is plotted against the initial value, as well as the non-dimensionalized energy release ΔW.

Zoom In Zoom Out Reset image size

Since the toroidal plasma current is I_t = μΦ_t, figure 6 may also be interpreted, re-scaling the axes by the appropriate μ values, as showing the dependence of final plasma current on initial current (in one flux rope). However, the initial state also includes a current sheet, and in order to calculate the total initial plasma current, we must add the (negative) net current in this sheet, which may be simply obtained as $I_{\rm t,cs}=-(2/\mu_{0})\int\limits{_0^a}B_{x}(x,0)\,\rmd x$.

It can be shown analytically that for small μ_i, $\Delta W \sim \mu_{\rm i}^2$, μ_f ∼ μ_i. Hence the final plasma current depends linearly on the initial current; whilst the energy release depends quadratically on μ_i. The dependence of μ_f on μ_i shown in figure 6 appears as a straight line, but the scaling is in fact stronger than linear; for larger μ_i, approaching the first eigenvalue, the linearity clearly breaks down.

As shown in section 2, the magnetic energy released by reconnection may both heat the plasma and drive plasma flows—but the latter may also be dissipated to provide further heating. We obtain an estimate (strictly an upper bound) on the temperature increase—here simply within single-fluid MHD—by assuming that the released magnetic energy is converted fully into thermal energy. For illustration with typical MAST parameters, we take a = 0.4 m and B_t = 0.5 T ⇒ Φ_t = 0.08 Wb. Taking $\mu_{\rm i}^{*} = 1.25$ gives $I_{\rm t,cs}^{*}=-0.625$ (here, an asterisk denotes dimensionless quantities), leading to a total initial plasma current of 300 kA, compatible with experiment and the modelling described above (but with somewhat larger currents in the plasma rings). The energy release is ΔW^* = 0.0022, which, equating to the thermal energy (3/2)(2n)kT (2a²) gives a temperature rise of about 150 eV, which is compatible with what is experimentally observed—of course, the relaxation model cannot predict any spatial distribution of temperature. Note that the temperature increase depends (to a very good approximation) quadratically on μ_i and hence on the initial plasma current and on the pre-merger poloidal field, in agreement with experiment [15].

3.3. Comparison of relaxation model with numerical simulations

As noted above, the relaxation model captures the self-organization of two flux ropes into a single one, and the predicted relaxed state is in very good qualitative agreement with the final state of the simulations (comparing final panels of figures 2 and 5). As a further test, we calculate the time evolution of magnetic helicity (see equation (12) as well as the magnetic, kinetic and internal energies, in the resistive MHD simulation simialr to that described in section 2. This is shown in figure 7. The helicity is indeed subject to global resistive diffusion and decays, as expected; however, the helicity evolution is seen to be insensitive to the dynamics of the reconnection and shows no sign of the 'sloshing' oscillations which are evident in the energy traces. The energy traces show that the reconnection causes an increase in kinetic energy (associated with reconnection outflows), while the magnetic energy falls. The semi-oscillatory nature of the reconnection is clearly seen in all the energy traces, but most clearly in the kinetic energy. Note also that, in agreement with the relaxation picture, the kinetic energy eventually falls to a very low value, so that most of the dissipated magnetic energy is converted into thermal energy.

Zoom In Zoom Out Reset image size

The drop in normalized helicity, about 1.2 × 10⁻³, is somewhat lower than the corresponding drop in magnetic energy, 6.0 × 10⁻³. This is an indication that the magnetic energy release due to reconnection is not strongly dominating over the dissipation due to global Ohmic diffusion. In simulations with multiple reconnection sites, where dissipation due to reconnection is dominant, the disparity between energy and helicity dissipation is much greater [30]. Nevertheless, the relaxation picture still captures the essential features of the final state and the energy conversion. We would also expect the resistive dissipation of helicity to be far lower in more highly-conducting plasmas such as in the solar corona.

We have described two approaches to modelling the merging of flux ropes in the MAST device. In both cases, the self-organization through magnetic reconnection of two initial flux ropes into a single flux rope is predicted. Resistive MHD and two-fluid simulations using the HiFI framework, in both 2D Cartesian and toroidal geometries, provide detailed information about the dynamics of the reconnection process and its dependence on the collisionality and other parameters. A simple analytical model based on Taylor relaxation theory captures some aspects of the self-organization, predicting the final field state and the plasma heating for given initial conditions.

The numerical simulations presented here are all in 2D. The relaxation model could in principle account for a non-axisymmetric relaxed state, but this would only be attained for far higher values of μ than are relevant here [26, 27]. The twisted flux ropes have strong toroidal fields (typically 5 times the poloidal field), and hence are very unlikely to be subject to kink instabilities. The final state of the simulations also has q > 1 everywhere, indicating stability. (Note that the stability will depend also on the aspect ratio, which is larger for the initial individual flux ropes than for the final state.) Nevertheless, there may be some interest in considering 3D effects in future, particularly on the post-merger evolution, where experimental results sometimes indicate the development of non-axisymmetric structure.

One motivation for studying reconnection in MAST is the relevance to astrophysical and space plasmas. In the solar corona, reconnection occurs with dramatic consequences in solar flares, as well as being likely to contribute to the overall heating of the coronal plasma, and occurring frequently during emergence of flux ropes from the interior. A comparison between MAST and coronal parameters is shown in table 1. Whilst the coronal value of the Lundquist number S, at least on global scales, is very much larger than the MAST value, MAST has a higher S than any other currently operating laboratory reconnection experiment [15]. Coronal reconnection is widely modelled using single-fluid MHD. However, the predicted length-scales of reconnecting current sheets in solar flares are comparable with the ion skin-depth, provided that the dissipation-scale is not anomalously enhanced: thus, Hall terms (at least) should be included and the parameter regime is quite similar to that described here.

Work is currently in progress to extend the two-fluid simulations to include separate electron and ion temperature distributions, which may be compared with Thomson scattering measurements and ion temperature measurements, where available. It is also planned to explore the effects of changing the initial flux-rope currents. This modelling may also be adapted to study other devices with merging flux ropes.

The relaxation-based model will be extended to toroidal geometry in future. It will also be interesting to investigate further the comparison between the simulation results and the relaxation model, for example, by checking the scalings of heating with plasma current predicted by the latter. It will be particularly interesting to do this both for the resistive MHD and the Hall cases. In the latter case, alternatives to the standard Taylor relaxation approach should also be considered, with additional invariants [36, 37]. The relaxation model may also be easily applied to astrophysical phenomena; in this context, it would be interesting to consider merging and self-organization of multiple flux ropes e.g. [34, 38], with implications for solar coronal heating.

This work was funded by STFC, the US DoE Experimental Plasma Research program, the RCUK Energy Programme under grant EP/I501045, and by the European Communities under the Contract of Association between EURATOM and CCFE. The views and opinions expressed herein do not necessarily reflect those of the European Commission. Simulations were run at the Solar-Terrestrial Environment Laboratory (STEL) cluster at Nagoya University,the US National Energy Research Scientific Computing Centre and the St-Andrews MHD cluster funded by STFC. A S would like to thank Kanya Kusano for access to the STEL computer.