1 ｜ INTRODUCTION
Divertor heat load mitigation is one of the most crucial issues for realizing nuclear fusion reactors. For this purpose, development of stable divertor detachment operation scenario is very important. In helical devices, it has been observed that the edge magnetic field structure plays a key role in the detachment control. In LHD, it has been found that the application of resonant magnetic perturbation (RMP) field provides stable sustainment of divertor detachment  (1.2) . In W7-AS, the size of the edge magnetic island is found to be a key parameter to stabilize detachment [3] . The recent experiments of W7-X with island divertor configuration have also realized a stable operation of complete detachment [4] . The mechanisms of the detachment stabilization have been discussed, such as core-edge coupling through neutral penetration [3] , or competition between impurity confinement time and radiation penetration time  [5] , etc. However, the analyses have been limited to rather simplified models so far, while in order to understand the stability of the detachment, time-dependent analysis with appropriate geometrical effects is essential. In these contexts, analyses with time-dependent transport simulation with edge magnetic island structure have been prompted, in particular, geometrical effects of the magnetic field structure on the impurity radiation distribution and its consequences for the detachment stabilization have been investigated. The rest of this paper is structured as follows. In the next section, the transport model, geometry, and computation scheme are described. In section 3, the results of the computation are presented, where it is shown that dynamics of radiation layer formation is strongly affected by the magnetic island. The paper is summarized in section 4.
2 ｜ TRANSPORT MODEL
Figure 1 (a) shows the edge magnetic field structure of LHD with RMP application, which induces edge magnetic island as shown in the figure. Due to the inherent stochasticity at the edge region, the separatrix of the island is no longer clear, but the structure is visible with the bundle of the long flux tubes indicated with the dashed lines. In order to obtain basic physics insight into the geometrical effects of the magnetic island on the impurity radiation dynamics, as a first step, we intend to utilize the simple magnetic field structure in the present analysis. For this purpose, a two dimensional (2D) computation grid is constructed as shown in Fig.1   (b), where the island structure is modeled in a cylindrical geometry. An RMP field of m/n=1/1 with  = 10 −3 is applied, which induces an edge magnetic island of ～ 11cm width according to R B the island width formula, 4 where R 0 = 3.90 (m), m = 1, r = q = 5.1 (1/m), respectively ( l = 1/q). The island structure, which rotates in poloidal direction (θ) as moving in toroidal direction (φ), is mapped back to φ = 0 plane to create the 2D grid. The grid of the magnetic island is colored red, where the island is resolved with 20 radial grid surfaces. The "normal" flux surfaces outside of the island, which are colored green, are discretized with 23 grid surfaces in total, 12 and 11 surfaces at radially inner and outer than the island, respectively. The poloidal direction is discretized by 200 grid points. The magnetic field line lengths (pitch) between the poloidal grid points changes according to the iota profile. Especially, the change is significant near the X and O-points due to the singularity. These changes are stored in the grid points, and used for parallel transport. In the present model, the divertor legs are degenerated to a boundary condition at the outer-most surface, where Bohm boundary condition is applied as described later in this section. This assumption is introduced accounting for the experimental observations of the low density at the divertor legs, ～ 10 18 m − 3 (6.7) , which is smaller by an order of magnitude than the density at the stochastic layer, ≥ 1 × 10 19 m −3 . Thus impurity radiation is usually larger in the stochastic layer than in the divertor legs as analyzed in ref. [5,8] . The effects of the stochastic layer are not included in the present simplified geometry, but rather we focus on the effects of island structure as mentioned above. The effects of the stochasticity will be addressed in future works with such as enhanced cross-field transport coefficients etc. The following time-dependent equations for particle, parallel momentum, and energy transport, are solved in the code:
electron and ion densities, particle and heat diffusivities perpendicular to the field lines. η , = m,n,D, (kg/m/s), is a perpendic- ular viscosity. K ε , i = κ ερ io T ε j  (1 /m/s) is parallel heat conductivities for electron and ion, with κ ερ io = 2.1 × 10 22 , 6.5 × 10 20 , respectively for T ej in eV. S q , S m, S cj are the source terms for particle, parallel momentum, and energy transport, which are defined as,
(6)
(8)
Here k ion and k C X are ionization and charge exchange (CX) rate coefficients ( m ³ /s ). K represents an equilibration coefficient for energy exchange between electrons and ions, K = 3mznc with τ ci = 2.03 × 10 10 T 1.5 /n c (sec, T c in eV). The term, n c n imp L(T c ), is cooling of electron due to impurity radiation, where n imp is an impurity density summed over all charge states, n imp = 1 and L(T,) ( m 3W ) is a cooling rate of the impurity. In the present computations carbon is used as an impurity, of which the divertor plates is made and is dominant radiator in LHD, and non-coronal model with n,τ = 10 17 ( sm −3 ) is used for L(T,). [9] The neutrals are introduced as a "cold" atoms at the outer most surface of computation domain, n m ( a ), and decays in radial direction according to,
T α = 3 eV, assuming Flanck-Condon atoms dissociated from molecules. The collision between the atoms "hot" CX neutrals, n Cx , transport of which is modeled with fluid approximation as follows, and ions produces
T . In the present analyses, quasi-neutrality, n ς ≈ n is assumed. As we checked the ratio n imp /n ς a posteriori in the with V, = results of computations, it is a few percent at maximum. So that the present approximation is considered to be valid within the error in this range. For impurity, by summing up for all charge states, the particle transport equations is reduced to,
The impurity density is initially set zero everywhere at a lowest plasma density computation, where the impurity radiation is reasonably small. Then the neutral impurity flux, Γ ( θ ) ( α ), is launched at the outer-most surface proportional to the divertor TABLE 1  Input parameters for present computations.
particle flux, Γ ρ , with multiplication factor, sputtering coefficient, C sput , such that Γ ( ω ) = C sput Γ ρ . The impurity flux then decays radially as,
with T (0) = 2 eV is assumed. The obtained where k imp(0) is an ionization rate coefficient of neutral impurity, and n imp distribution is used as a initial condition for the next computations with increased plasma density ( n α t ( a )). The parallel impurity velocity, V impli , is obtained by assuming force balance between the friction force and the ion thermal force (ion temperature gradient force) [10]
(ln Λ = 17), the impurity slowing down time. It is noted that the above equation becomes charge independent due to the cancellation of Z imp between C i and τ s . D imp = D 1 is assumed for the present computations. Carbon, sputtered from the divertor target plates, is considered as the dominant impurity species The boundary conditions are defined as follows. At the inner-most surface, power flux to the computation domain, PSOL- is given, which is uniformly distributed in a poloidal direction, and is equally divided to electrons and ions. The neutral flux that penetrates into the confinement region through the inner-most surface is given to plasma particle as a net flux through the boundary. As mentioned above, the divertor legs are degenerated to a boundary condition at the outer-most surface, where Bohm boundary condition is applied for particle and energy transport with multiplication factor β, that is, Γ ρ = β nc, and , T + T , sound speed, and γ ej being sheath heat transmission coefficients for electron and Γ ç = β(γ ç Τ ç + γ ç Τ ç ) nc ç , with c ç = Br ≈ 10 −3 is introduced to convert the particle and energy flux at the divertor plates to the ion, respectively. Here β = r · b = perpendicular flux at the outer-most surface. The input parameters used in the present computations are summarized in table 1  . The transport equations are integrated by applying an implicit finite volume scheme on the 2D grid. The time step of the computations are set to 2 × 10 −7 ～ 2 × 10 −8 (sec). The atom density at the outer-most surface, n α ( α ), is gradually increased to reproduce a prescribed level of the electron density at the core boundary, and thus to induce detachment transition. This method mimics gas puff scan in experiments.
3 ｜ RESULTS OF COMPUTATIONS
The computations are started with low plasma density where the impurity radiation is reasonably small, and density is gradually increased by increasing n α ( a ) stepwise. Usually, the computation reaches steady state within physical time of 0.1 - 0.2 sec after the stepwise increase of n ,,,( a ). When n ,,,( a ) was increased to 2.9× 10 17m −3 , the radiated power from the impurity starts to exceed the power to the divertor plates. At the convergence, the power balance with respect to P SOL is as follows. The power to the divertor plates: 25%, the power of impurity radiation: 54%, energy loss via hydrogen (excitation and ionization): 11%, energy loss via CX: 10%. When n n ( a ) was increased to 3.0 × 10 17m −3 , rapid increase of radiated power and concomitant decrease of the divertor power/particle load are encountered after 47 millisecond in physical time. The time traces of radiated power and divertor particle flux integrated over the computation domain are plotted in Fig. 2   (a). The onset of detachment transition, t = 0 sec, is defined at the timing of decrease of the divertor particle flux as shown in Fig. 2 (a). It is seen that the radiated power fraction starts increasing significantly after the onset, even slightly exceeding 100%, during which significant cooling of plasma takes place in an unsteady state. The temporal evolution of T c and n c at various locations are also plotted in Figs. 2   (b-e). T , at the downstream (outermost radial surface) of X-point region (Fig. 2   (b), "down") starts to decrease at the same time of the divertor particle flux reduction. On the other hand, T c at at the downstream of O-point region (Fig. 2   (c), "down") decreases with delay. This indicates that the detachment is initiated around the X-point. n , at the downstream of X-point region (Fig. 2   (d), "down") shows significant increase right after the detachment onset. This is due to the thermal condensation instability caused by the magnetic field structure of the island as explained below. Figure 3   shows 2D distribution of plasma parameters during the detachment transition. From Fig.3  (a1-3), it is observed that the radiated power by impurity, P rad , is enhanced around/at the X-point. This is found to be due to the geometrical effect of the island structure on the transport, as explained in the following. Along the island separatrix, which connects radially inner and outer regions of O-point, parallel temperature gradient is developed through parallel transport, as indicated by arrows in Fig.3 (b1), where the arrow directions indicate the direction of energy flux. The energy flux that reaches the outer edge of the O-point diffuses out perpendicularly and, at the same time, transported toward the side of X-point through parallel transport, as indicated by the outer-most arrows in Fig.3  (b1), and thus develops parallel temperature gradient there. Through the parallel momentum balance, these temperature gradients develop parallel flow, V ill , in the same direction of the arrows in Fig.3 (b1). The impurity flow, V impll , was found to have similar flow pattern as V ill . The resulting flow field of V impli is plotted in Fig.3  (d) as a poloidal projection of the parallel flow. The flow at the outer-most region is clearly visible, which is directed toward the side of X-point. This outer-most flow is most effective to drive the impurity toward the side of X-point, since the impurity comes from the outside of plasma and has the largest density at the periphery. Due to the increased impurity at the X-point side, P og is increased. This geometrical effect is observed already at the attached phase, and is enhanced as approaching the detachment onset since temperature gradient develops further due to the lowering of the temperature. After the detachment onset, the radiation region moves radially inward as shown in Fig.3  (a2), and the temperature is fur- ther cooled down, Fig.3  (b2). The impurity flow towards the X-point side and localization of impurity there are enforced as seen in Figs.3  (d2), and (e2). The density peaking at the radiation region is also observed through the momentum conservation, Fig.3  (c2). These processes facilitate each other, leading to the radiation increase at the X-point region, exhibiting thermal con- densation instability. At t = 12 ms, the peaked radiation reaches the X-point. Indication of the X-point radiation is also observed in experiments with AXUV and bolometer measurements (1,11). The process described above may explain the experimental observation. The computations were further continued and the radiation layer penetrates further into the island, where the temperature of the island decreases down to ～ 1 eV. In the experiments, on the other hand, the temperature at the island is kept at around 10 to 20 eV and the radiation region seems to stay around the X-point. This discrepancy between the modeling and the experiments will be further studied by investigating various aspects of the model, such as dynamical change of impurity sputtering coefficient 6 IRE 3  2D distributions of (α) impurity radiation P nas , (b) electron temperure T nas electron density nas (d) represent director the site of dash not transition. The arrow
4 | SUMMARY
ACKNOWLEDGMENTS