The Kelvin-Helmholtz instability in weakly ionised plasmas: Ambipolar dominated and Hall dominated flows
aa r X i v : . [ a s t r o - ph . GA ] J u l Mon. Not. R. Astron. Soc. , 1–12 (—-) Printed 13 November 2018 (MN L A TEX style file v2.2)
The Kelvin-Helmholtz instability in weakly ionisedplasmas : Ambipolar dominated and Hall dominated flows.
A.C. Jones , and T. P. Downes , , ⋆ School of Mathematical Sciences, Dublin City University, Glasnevin, Dublin 9, Ireland School of Cosmic Physics, Dublin Institute for Advanced Studies, 31 Fitzwilliam Place, Dublin 2, Ireland National Centre for Plasma Science and Technology, Dublin City University, Glasnevin, Dublin 9, Ireland
Accepted –. Received –; in original form –
ABSTRACT
The Kelvin-Helmholtz instability is well known to be capable of converting well-ordered flows into more disordered, even turbulent, flows. As such it could representa path by which the energy in, for example, bowshocks from stellar jets could be con-verted into turbulent energy thereby driving molecular cloud turbulence. We presentthe results of a suite of fully multifluid magnetohydrodynamic simulations of this in-stability using the HYDRA code. We investigate the behaviour of the instability in aHall dominated and an ambipolar diffusion dominated plasma as might be expectedin certain regions of accretion disks and molecular clouds respectively.We find that, while the linear growth rates of the instability are unaffected bymultifluid effects, the non-linear behaviour is remarkably different with ambipolardiffusion removing large quantities of magnetic energy while the Hall effect, if strongenough, introduces a dynamo effect which leads to continuing strong growth of themagnetic field well into the non-linear regime and a lack of true saturation of theinstability.
Key words: mhd – instabilities – ISM:clouds – ISM:kinematics and dynamics
The Kelvin-Helmholtz (KH) instability is an important in-stability in almost any system involving fluids: it can occuranywhere that has a velocity shear. In astrophysical plas-mas the KH instability can provide the means of producingturbulence in a medium or the mixing of material betweentwo boundary layers.The KH instability has been studied in a variety ofastrophysical systems, from solar winds (Amerstorfer et al.2007; Bettarini et al. 2006; Hasegawa et al. 2004) and pul-sar winds (Bucciantini & Del Zanna 2006) to thermal flares(Venter & Meintjes 2006). Due to its ability to drive mix-ing and turbulence, the KH instability has been con-sidered relevant in protoplanetary disks (Johansen et al.2006; G´omez & Ostriker 2005), accretion disks and magne-tospheres (Li & Narayan 2004), and other jets and outflows(Baty & Keppens 2006). Generally speaking, the assump-tions of ideal magnetohydrodynamics (MHD) have beenused in order to simplify the system of equations to besolved. These assumptions are, however, not always valid.Weakly ionised plasmas, for example, contain a large frac-tion of neutral particles as well as a number of charged parti- ⋆ E-mail: [email protected] (TPD) cle fluids with differing physical characteristics. Interactionsbetween the various species can introduce non-ideal effects.Ambipolar dissipation and the Hall effect are two non-idealeffects that can greatly influence the development of theKH instability in a system by altering the dynamics of theplasma and the evolution of the magnetic field. Astrophys-ical examples of such weakly ionised systems include densemolecular clouds (e.g. Ciolek & Roberge 2002) and accre-tion disks around young stellar objects (e.g. Wardle 1999).In these systems, the relevant length scales are such thatnon-ideal effects can play an important role (Wardle 2004a;Downes & O’Sullivan 2009, 2011).Many authors have investigated the role of the KHinstability in both magnetised and unmagnetised astro-physical flows (e.g. Frank et al. 1996; Malagoli et al. 1996;Hardee et al. 1997; Downes & Ray 1998; Keppens et al.1999). Most of these studies have investigated the KH in-stability in the context of either hydrodynamics or idealMHD. We know that non-ideal effects are important inmolecular clouds at length scales below about 0.2 pc (e.g.Oishi & Mac Low 2006; Downes & O’Sullivan 2009) andhence it is of interest to explore the KH instability in thecontext of either non-ideal MHD or, preferably, fully multi-fluid MHD. In more recent years the emphasis of KH stud-ies has been on including non-ideal effects. Keppens et al. c (cid:13) —- RAS A. C. Jones and T. P. Downes (1999) studied both the linear growth and subsequent non-linear saturation of the KH instability using resistive MHDnumerical simulations. The inclusion of diffusion allowed formagnetic reconnection and non-ideal effects were observedthrough tearing instabilities and the formation of magneticislands.Palotti et al. (2008) also carried out a series of simu-lations using resistive MHD. They found that, following itsinitial growth, the KH instability decays at a rate that de-creases with decreasing plasma resistivity, at least within therange of of resistivities accessible to their simulations. Theyalso found that magnetisation increases the efficiency of mo-mentum transport, and that the transport increases withdecreasing resistivity. Birk & Wiechen (2002) examined thecase of a partially ionised dusty plasma, using a multifluidapproach in which collisions could be included or ignored.They found that collisions between the neutral fluid anddust particles could lead to the stabilisation of KH modesof particular wavelengths. The unstable modes led to a sig-nificant local amplification of the magnetic field strengththrough the formation of vortices and current sheets. In thenonlinear regime they observed the magnetic flux being re-distributed by magnetic reconnection. It was suggested thatthis could be applicable to dense molecular clouds and haveimportant implications for the magnetic flux loss problem(Umebayashi & Nakano 1990).A comprehensive study was carried out by Wiechen(2006) which demonstrated the effect of dealing with theplasma using a multifluid scheme. This study focused onthe effect of varying the properties of the dust grains. Theresults of the simulations led to the conclusions that moremassive dust grains have a stabilising effect on the systemwhile higher charged numbers have a destabilizing effect. Itwas found that there is no significant dependence on thecharge polarity of the dust.In a linear study of stellar outflows Watson et al. (2004)described how the charged and neutral fluids are affecteddifferently by the presence of a magnetic field. This studyis carried out using parameters chosen to reflect those ofmolecular clouds, and so is particularly relevant to our ownstudy. The principal result of this paper is that for muchof the relevant parameter space, neutrals and ions are suf-ficiently decoupled that the neutrals are unstable while theions are held in place by the magnetic field. Since the mag-netic field is frozen to the ionised plasma, it is not tangledby the turbulence in the boundary layer. The authors pre-dict that with well-resolved observations, there should be adetectably narrower line profile in ionised species tracing thestellar outflow compared with neutral species, since ionisedspecies are not participating in the turbulent interface withthe ambient interstellar medium. The paper also includes astudy of the growth rate of the instability. It is found thatat short length scales, the growth rate is well approximatedby the growth rate of the hydrodynamic system. At largerscales and for super-Alfv´enic flows, the fastest growing modeis equal to that of the ideal MHD case.Shadmehri & Downes (2008) carried out an analyticalstudy of the Kelvin-Helmholtz instability in dusty and par-tially ionised outflows. They investigated primarily the ef-fect of the presence of dust particles by varying their mass,charge and charge polarity. It was found that as the charge ofthe grain increased, the growth timescales also increased, im- plying a stabilising effect on the system. The stability of thesystem was also examined for dependence on the mass of thedust particles. It was found that for stronger magnetic fields,this did not affect the stability of the system. However, forweaker magnetic fields, the larger dust particles had a sta-bilising effect on the growing modes. This was in agreementwith previous laboratory experiments (Luo et al. 2001) andnumerical simulations (Wiechen 2006). Finally, as the mag-netic field strength increased, the growth timescale of theunstable modes at a particular perturbation wavelength de-creased. By examining the combinations of the wavelength ofthe perturbation used, and the resultant growth timescalesof the instability, Shadmehri & Downes (2008) concludedthat the Kelvin-Helmholtz instability is a possible candidatefor causing the formation of some of the physical structuresobserved in molecular outflows from young stars.In this paper we perform numerical simulations of thecomplete evolution of the Kelvin-Helmholtz instability in aweakly ionised, multifluid plasma including both its lineardevelopment, saturation and its subsequent behaviour. Weinclude the physics of the Hall effect, ambipolar diffusionand parallel resistivity in these simulations and we analysethe role each of the former two effects on the developmentof the instability.The aim of this work is to investigate the influence ofmultifluid effects on the growth and saturation of the KHinstability. In order to develop a full understanding of theroles of the various non-ideal effects, in particular the Halleffect and ambipolar diffusion, we run simulations with pa-rameters chosen to simulate very high, medium and verylow magnetic Reynolds number systems and with parame-ters chosen to ensure ambipolar-dominated flows and Hall-dominated flows. We focus on gaining an understanding ofthe general characteristics of the KH instability in weaklyionised flows. A detailed study of this instability in the spe-cific context of molecular clouds, and which is of interestfrom the point of view of turbulence generation by stellaroutflows, is the subject of a future work.In section 2 we outline the numerical and physical modelemployed, in section 3 we discuss how we analyse the resultsof our simulations while in section 4 we detail the results inboth the linear and non-linear regimes, separating out theeffects of ambipolar diffusion and the Hall effect in order tomore fully understand the influence of each.
The simulations described in this work are performed us-ing the HYDRA code (O’Sullivan & Downes 2006, 2007)for multifluid magnetohydrodynamics in the weakly ionisedregime. We further assume the flow is isothermal. The as-sumption of weak ionisation allows us to ignore the inertiaof the charged species and allows us to derive a (relatively)straightforward generalised Ohm’s law. The resulting sys-tem of equations, given below, incorporates finite parallel,Hall and Pederson conductivity and is valid in, for example,molecular clouds. In such regions the viscous lengthscalesare much smaller than those over which nonideal effects areimportant. This leads to high Prandtl numbers and plasmaflows in these regions can be considered to be effectivelyinviscid. In this work we examine low Mach number flows c (cid:13) —- RAS, MNRAS , 1–12 he Kelvin-Helmholtz instability in weakly ionised plasmas which, taken in concert with the isothermal assumption,means that features in the flow such as shocks are unlikelyto create regions of high ionisation. The HYDRA code solves the following equations for a sys-tem of N fluids. The simulations described in this paperconsist of three fluids, indexed by i = 0 for the neutral fluidand i = 1 and i = 2 for the electron and ion fluids respec-tively. The equations to be solved are ∂ρ i ∂t + ∇ · ( ρ i q i ) = 0 (0 i N − , (1) ∂ρ q ∂t + ∇ · ( ρ q q + a ρ I ) = J × B , (2) ∂ B ∂t + ∇ · ( q B − Bq ) = −∇ × E ′ , (3) α i ρ i ( E + q i × B )+ ρ i ρ K i ( q − q i ) = 0 (1 i N − , (4) ∇ · B = 0 , (5) ∇ × B = J , (6) N − X i =1 α i ρ i = 0 , (7)where ρ i , q i , B , and J are the mass densities, velocities,magnetic field and current density, respectively. a denotesthe sound speed, and α i and K i are the charge-to-mass ra-tios and the collision coefficients between the charged speciesand the neutral fluid, respectively.These equations lead to an expression for the electricfield in the frame of the fluid, E ′ , given by the generalisedOhm’s Law E ′ = E O + E H + E A , (8)where the components of the field are given by E O = ( J · a O ) a O , (9) E H = J × a H , (10) E A = − ( J × a H ) × a H , (11)using the definitions a O ≡ f O B , a H ≡ f H B , a A ≡ f A B ,where f O ≡ √ r O /B , f H ≡ r H /B and f A ≡ √ r A /B . Theresistivities given here are the Ohmic, Hall and ambipolarresistivities, respectively, and are defined by r O ≡ σ O , (12) r H ≡ σ H σ + σ , (13) r A ≡ σ A σ + σ , (14)where the conductivities are given by σ O = 1 B N − X i =1 α i ρ i β i , (15) σ H = 1 B N − X i =1 α i ρ i β i , (16) σ A = 1 B N − X i =1 α i ρ i β i β i , (17)where the Hall parameter β i for a charged species is givenby β i = α i BK i ρ . (18)To solve these equations numerically we use three dif-ferent operators:(i) solve equations (1), (2), (3), including the restrictionof equation (5) and for i = 0, using a standard second or-der, finite volume shock-capturing scheme. Note that for thisoperator the resistivity terms in equation (3) are not incor-porated. Equation (5) is incorporated using the method ofDedner (Dedner et al. 2002).(ii) Incorporate the resistive effects in equation (3) us-ing super-time-stepping to accelerate the ambipolar diffu-sion term and the Hall Diffusion Scheme to deal with theHall term.(iii) Solve equations (4) for the charged species velocitiesand use these to update equation (1) (with i = 1 , . . . , N − The simulations are carried out on a 2.5 D slab grid in the xy -plane. The grid consists of 6400 × × x , y , and z directions respectively. This resolution was chosenon the basis that it reproduces the initial linear growth ofthe ideal MHD system in Keppens et al. (1999). Resolutionstudies were performed to confirm the resolution as beingappropriate (see Sect. 4.2.1 and 4.3.1).The initial set-up used was that of two plasmas flowinganti-parallel side-by-side on a grid of size x = [0 , L ] and y = [0 , L ]. The plasma velocities are given by + V and − V in the y -direction, with a tangential shear layer of width 2 a at the interface at x = 16 L . This velocity profile is describedby v = V (cid:16) x − La (cid:17) ˆ y . (19)The width of the shear layer is chosen to be aL = 0 .
05, orapproximately 20 grid zones. The magnetic field is initiallyset to be uniform and aligned with the plasma flow.The initial background for all three fluids in the systemis now an exact equilibrium. The initial neutral velocity field, V is then augmented with a perturbation given by δv x = δV sin( − k y y ) exp (cid:18) − ( x − L ) σ (cid:19) . (20)where δV is set to 10 − V . The wavelength of the per-turbation is set equal to the characteristic length scale, λ = πk y = L , so that a single wavelength fits exactly intothe computational domain. This maximises the possibilityof resolving structures that are small relative to the initialperturbed wavelength, (Frank et al. 1996). The perturba-tion attenuation scale is chosen so that it is larger than the c (cid:13) —- RAS, MNRAS , 1–12 A. C. Jones and T. P. Downes
Fluid Density K i α i . × − × − × × − . × . × Table 1.
The density, collision coefficient ( K i ) and charge-to-mass ratios ( α i ) for each of the charged fluids in simulation full-low-hr (see table 2). These parameters are modified to vary theresistivities as necessary for the other simulations in this work.See text. shear layer, but small enough so that the instability can beassumed to interact only minimally with the x -boundaries(see Palotti et al. 2008), and is set using σL = 0 . k y ischosen to be 2 π in order to maximise the growth rate of theinstability (Keppens et al. 1999). This normalises the lengthscale of the simulation so that L = 1. The timescale is thennormalised by the sound speed, so that c s = LT = 1. Themass scale is chosen such that the initial mass density is setto unity, ρ = 1. In the isothermal case the adiabatic index γ = 1, and this gives us a sound speed c s ≡ p γp ρ = 1,and the initial pressure is therefore also equal to unity, p = 1. A transonic flow is chosen with sonic Mach number M s = V c s = 1, so that the plasma has velocity ± V = ± . M A = V v A = 10 is cho-sen. This sets the Alfv´en velocity, v A = q B ρ = 0 .
1, andthe initial magnetic field strength to B = 0 . y boundaries. Since we wish to study not only the initialgrowth phase of the instability, but also its subsequent non-linear behaviour we must ensure that waves interacting withthe high and low x boundaries do not reflect back into thedomain to influence the dynamics. Several test simulationsfor various parameters have shown that a large width of 32is necessary to ensure this. We use gradient zero boundaryconditions at the high and low x and z boundaries.Finally we must choose the parameters describing theproperties of our charged fluids in our multifluid system.Our basic parameter set is contained in table 1. See Jones(2011) for more details.Table 2 contains the nomenclature we will use for therest of this paper when referring to the simulations. Eachsimulation is denoted by xxxx-yyyy-zz where the first setof characters denote the dominant resistivity (Hall or am-bipolar or “full” if both resistivities have the same mag-nitudes), the second set denote the level of the resistivity(low, medium or high) and the final two characters denotethe resolution (low, medium or high resolution denoted bylr, mr and hr respectively). The high resolution simulations(6400 × Simulation Resolution r H r A mhd-zero-hr 400 ×
200 0 0hd-zero-hr 400 ×
200 0 0ambi-high-lr 1600 ×
50 3 . × − . × − ambi-high-mr 3200 ×
100 3 . × − . × − ambi-high-hr 6400 ×
200 3 . × − . × − full-low-hr 6400 ×
200 3 . × − . × − ambi-med-hr 6400 ×
200 3 . × − . × − hall-high-lr 1600 ×
50 3 . × − . × − hall-high-mr 3200 ×
100 3 . × − . × − hall-high-hr 6400 ×
200 3 . × − . × − hall-med-hr 6400 ×
200 3 . × − . × − Table 2.
The resolutions and (initial) resistivities of each of thesimulations presented in this work. The Hall resistivity is denotedby r H while the ambipolar resistivity is denoted by r A . Note thatthe Pederson conductivity is always very large and hence its as-sociated resistivity is always several orders of magnitude smallerthan both r H and r A . In order to study the growth of the instability, the evolutionof a number of parameters can be measured with time. Inparticular, we measure the transverse kinetic energy E k ,x ≡ Z Z ρv x dx dy and the magnetic energy E b ≡ Z Z (cid:8)(cid:2) B x + B y + B z (cid:3) − B (cid:9) dx dy in the system where B is the magnitude of the magneticfield at t = 0. As the entire plasma flow is initially in the y -direction, with only a very small perturbation in the x -direction, any growth of E k ,x is due to the growth of theinstability. It is possible to determine the growth rate of theinstability directly from the growth of the transverse kineticenergy (Keppens et al. 1999): the transverse kinetic can beexpressed as E k ,x = ( ρ + δρ )( | v x | + δv x ) (21)= ( ρ + δρ )( δv x ) (22) ≈ ρ δv x (23) ∝ exp[2 i ( k · r − ωt )] (24)presuming the initial perturbation is proportional to exp[ i ( k · r − ωt )]. Hence the kinetic energy grows at a rate of 2 ω , where ω is the growth rate of the instability. We start by describing the validation of the set-up used bycomparing an ideal MHD simulation run using HYDRA withpreviously published literature. We then go on to discuss thebehaviour of the KH instability in ambipolar-dominated andHall-dominated flows respectively. c (cid:13) —- RAS, MNRAS , 1–12 he Kelvin-Helmholtz instability in weakly ionised plasmas s -25-20-15-10-50 l og ( / ρ v x ) Figure 1.
Plot of the log of the transverse kinetic energy withtime for mhd-zero-hr. A linear function (dashed line) with slope2.63 has been over-plotted. s / B Figure 2.
Plot of the evolution of magnetic energy with time forsimulation mhd-zero-hr.
In order first to validate our set-up we examine our mhd-zero-hr simulation with HYDRA and determine the growthrate using the kinetic energy of motions in the x directionas described in Sect. 3. Figure 1 contains a plot of the log ofthe transverse kinetic energy as a function of time. At earlytimes this growth is clearly exponential and can be fittedwith a line of slope 2.63 implying a growth rate, normalisedby the width of the shear layer and the initial relative veloc-ity, for the dominant mode of the KH instability of 0.1315.We can compare this with the value of the growth rate calcu-lated analytically by Miura & Pritchett (1982) (their Figure4) at this wavenumber of 0.13. While comparisons betweenlinear studies of incompressible flows, and numerical stud-ies of compressible flows are bound to differ to some extent,these results are seen to agree exceptionally well.We wish to examine not only the linear regime but alsothe non-linear regime. We compare our results for the growthof magnetic energy with those of Malagoli et al. (1996) (theupper panel of their Figure 5). Figure 2 contains a plot of themagnetic energy, calculated as R R (cid:0) B x + B y + B z (cid:1) dx dy , asa function of time. The maximum magnetic energy reachedin our simulations matches that of Malagoli et al. (1996) towithin 10%. Our simulation reaches saturation at a latertime but the exact time of saturation depends on the initialamplitude of the perturbation and so this is not a concern.We are therefore confident of the behaviour of HYDRAin simulating the KH instability. We now move on to investi-gating the influence of multifluid MHD effects on the growth,saturation and non-linear behaviour of this instability. s / ρ v x Figure 3.
Plot of the evolution of the transverse kinetic energyagainst time in the case of high ambipolar resistivity for ambi-high-lr, ambi-high-mr and ambi-high-hr.
We begin our study of the multifluid KH instability bychoosing our fluid parameters to ensure our ambipolar resis-tivity is dynamically significant while minimising the Hallresistivity. This allows us to isolate the influence of the am-bipolar resistivity on the instability. In order to increase theambipolar resistivity we change the value of the collisioncoefficient for species 2 so that K , is decreased by 3 or-ders of magnitude from that given in table 1 for simulationambi-med-hr and by 4 orders of magnitude for ambi-high-hr.These alterations of K , give values of r A of 3 . × − and3 . × − respectively, and (ambipolar) magnetic Reynold’snumbers, Re m , of 2 . × and 2 . × respectively.These simulations are examined in comparison to the full-low-hr simulation. With a formal magnetic Reynolds num-ber 2 . × , the diffusion in this set-up is predominantlynumerical, and as such, it is effectively an ideal MHD simu-lation. In non-ideal MHD we must ensure that the length scalesover which the diffusion of the magnetic field (or the whistlerwaves in the case of Hall dominated flows) must be resolvedin order to properly track the dynamics of the system. Tothis end we perform a resolution study using simulationsambi-high-lr, ambi-high-mr and ambi-high-hr (see Table 2).Figures 3 and 4 contain plots of the evolution of E k ,x and E b for each of the simulations in our resolution study. Itcan be seen that the linear growth in ambi-high-lr is signif-icantly lower than the two other simulations. However, thelinear behaviour is almost identical for ambi-high-mr andambi-high-hr. The subsequent non-linear behaviour is simi-lar with only relatively small variations after t ∼ ×
200 is sufficient to capture the initial growth andsaturation of the instability. Subsequently, the dynamics iscaptured at least qualitatively.
Generally speaking, the evolution of the KH instability inideal MHD leads to a wind-up of both the plasma and themagnetic field at the interface between the two fluids, re-sulting in the “Kelvin’s cat’s eye” vortex. Multifluid effects c (cid:13) —- RAS, MNRAS , 1–12 A. C. Jones and T. P. Downes s / B - / B Figure 4.
As figure 3 but for the perturbed magnetic energy.
15 16 17 x0.00.20.40.60.81.0 y y Figure 5.
Plots of the magnitude and vector field of the magneticfield for full-low-hr, at time t = 8 t s (upper panel) and for ambi-high-hr (lower panel). The magnetic field is wound up in a similarfashion to the velocity in the upper panel, but not the lower panel. alter the nature of this significantly. Figure 5 contains plotsof the magnetic field at t = 8 t s which is the time at whichthe instability saturates for both full-low-hr (effectively idealMHD) and ambi-high-hr. The difference in the morphologyis clear.Given these striking differences at saturation it is in-teresting to investigate whether the linear growth rate isinfluenced by the addition of ambipolar diffusion. Figure 6contains plots of the evolution of E k ,x with time for thevarious simulations. The linear growth rate remains almostunchanged with the addition of ambipolar diffusion. Onthe other hand, figure 7 contains plots of E b as a func-tion of time. It is clear that the perturbed magnetic energyis strongly influenced, even well before saturation, by thepresence of ambipolar diffusion. We will discuss this in moredetail in section 4.2.3. s -25-20-15-10-50 l og ( / ρ v x ) Figure 6.
Plots of the log of the transverse kinetic energyfor varying ambipolar resistivity. The thicker lines represent themodels with higher ambipolar resistivity. The dashed line repre-sents hd-zero-hr. There is very little difference between the lineargrowths for the three magnetised cases. s / B - / B Figure 7.
The perturbed magnetic energy in the system is plot-ted against time in the three cases of varying ambipolar resistivity.The thicker lines represent the models with higher ambipolar re-sistivity. Higher ambipolar resistivity causes the magnetic field toexperience less amplification, through diffusion and decoupling ofthe bulk fluid from the magnetic field.
It is clear from figure 7 that the magnetic energy is stronglyinfluenced by the presence of ambipolar diffusion. This isnot too surprising as ambipolar diffusion, being a genuinelydiffusive process (unlike the Hall effect) allows the magneticfield to diffuse relative to the bulk flow.As the collision rate between the ion and neutral fluidsis decreased, the ion fluid decouples from the bulk fluid, andthus the magnetic field becomes decoupled from the bulkflow: the frozen-in approximation of ideal MHD is broken.As a result, the magnetic field is able to diffuse through thebulk fluid rather than being tied to it. In figure 7, diffusioncan be identified as the cause of the decrease of the am-plification of the magnetic energy with time for increasingambipolar resistivity. A cursory examination of the topologyof the magnetic field for the low resistivity case (full-low-hrin figure 5) demonstrates clearly that there are regions in theflow which will be susceptible to diffusion: regions in whichthe magnetic field lines have been compressed and amplified.It is clear, again from figure 7, that there is a signifi-cant decrease in the growth of magnetic energy as a resultof diffusion for even a moderate amount of ambipolar re-sistivity. The diffusion has an influence on the dynamics ofthe neutral fluid: when the magnetic energy has not been c (cid:13) —- RAS, MNRAS , 1–12 he Kelvin-Helmholtz instability in weakly ionised plasmas s -10 -9 -9 -9 / ( ρ v x ) i on Figure 8.
Plots of the transverse kinetic energy of the ion fluidas a function of time for each of full-low-hr, ambi-med-hr, ambi-high-hr and hd-zero-hr. The thicker lines represent the modelswith higher ambipolar resistivity. The dashed line represents thehydrodynamic case. as strongly amplified, the field can no longer exert the sameeffect on the neutral and ion fluid. Examination of figure 6reveals that the peak reached by E k ,x increases with increas-ing ambipolar resistivity. The value reached tends towardthe hydrodynamic limit for two reasons: increasing resistiv-ity implies increasing diffusion and hence a decreasing fieldstrength, and it also implies less coupling between the mag-netic field and the neutral fluid with increasing resistivity.As expected, the weaker magnetic field strength allows thevortex to become more rolled up (Faganello et al. 2009).It can be shown that for the ambi-med-hr simulation,the increase in the wind-up of the bulk fluid is due solelyto the first source: the magnetic diffusion. The slight de-coupling of the ion and neutral fluid is sufficiently high toallow magnetic diffusion, while still being sufficiently low toforce the charged fluids to behave in a manner similar tothe bulk fluid. This can be seen in a plot of the transversekinetic energy of the ion fluid, as shown in figure 8. Withmoderate amounts of ambipolar diffusion, the transverse ki-netic energy of the ion fluid (and electron fluid) reaches ahigher maximum, meaning that the fluid experiences morewind-up.However, as the ambipolar resistivity is increased fur-ther, the ion fluid becomes more decoupled from the neutralfluid, and is instead influenced primarily by the dynamics ofthe magnetic field. This significant decoupling occurs onlyfor magnetic Reynolds number lower than Re m ≈ We now turn our attention to the likely influence of theHall effect on the KH instability in multifluid MHD flows.In order to attain the resistivities we want (see table 2)we reduce the charge-to-mass ratio of the ion fluid, causingits Hall parameter to become smaller, causing the ion fluiddynamics to more closely resemble the bulk (neutral) fluiddynamics. The electrons are, however, still well-tied to thefield lines and so a relative drift emerges between the ion s / ρ v x Figure 9.
Plot of the transverse kinetic energy against time forsimulations hall-high-lr, hall-high-mr and hall-high-hr. and electron fluids leading to a current perpendicular to themagnetic field and hence the Hall effect.
As in the ambipolar resistivity study, a resolution study iscarried out to ensure that the small-scale non-ideal dynam-ics are captured. For this purpose, simulations are again runat three different resolutions (see Table 2). These simula-tions are run with the highest level of Hall resistivity toensure that the smallest-scale dispersive effects are in placewhen examining whether they are sufficiently resolved. Theinclusion of the Hall term in the dispersion relation in prin-ciple allows for the introduction of waves with a signal speedwhich tends towards infinity as their wavelength tends to-wards zero. While the Hall term is handled in the equa-tions by the HYDRA code using the explicit Hall DiffusionScheme (HDS) (O’Sullivan & Downes 2006, 2007), the codenaturally does not resolve these waves of vanishing wave-length.As has been seen, the introduction of non-ideal effectsdoes not tend to greatly influence the linear growth rateof the instability. As a result, the growth rate does not pro-vide a good means of measuring convergence with increasingresolution. The nonlinear evolution of the transverse kineticenergy, E k ,x is, however, strongly influenced by the non-idealeffects. We do not examine the evolution of the perturbedmagnetic energy as, in the case of high Hall resistivity, it nolonger demonstrates a simple growth to an initial maximum.The evolution of the transverse kinetic energy for each sim-ulation is plotted in figure 9. It can clearly be seen that thesimulations have started to converge at higher resolutions.While there is a notable difference between the two simula-tions of lower resolution, the gap closes significantly in thecomparison between the two simulations of higher resolu-tion. In particular, the initial maxima of E k ,x are almostidentical in simulations hall-high-mr and hall-high-hr.We are, therefore, confident that the dynamics in thehall-hr are well resolved and that our conclusions as to thephysical processes occurring are well-founded. Figures 10 and 11 contain plots of E k ,x and E b as a functionof time for simulations hall-med-hr and full-low-hr. In thelinear regime, the influence of the Hall effect on the growth c (cid:13) —- RAS, MNRAS , 1–12 A. C. Jones and T. P. Downes s -25-20-15-10-50 l og ( / ρ v x ) Figure 10.
The log of the transverse kinetic energy in the sys-tem is plotted against time for the full-low-hr (thin line) andhall-med-hr (thick line). The dashed line represents hd-zero-hr.There is very little difference between the linear growths for themagnetised systems. s / B - / B Figure 11.
The perturbed magnetic energy in the system forboth full-low-hr (thin line) and hall-med-hr (thick line) are plot-ted against time. No significant difference can be seen betweenthe maxima reached in each case. rate of the instability, and on the magnetic and kinetic en-ergy in the system at saturation can be seen to be negligible.Hence we can conclude that neither ambipolar diffusion northe Hall effect influence the energetics of the system in thelinear regime.We expect the Hall effect to re-orient the magnetic fieldout of the xy -plane in which it resides at t = 0. Figure 12contains plots of E b , the perturbed magnetic energy in the xy -plane and the perturbed magnetic energy in the z direc-tion as functions of time. It is clear that there is some growthof the magnetic field in the z direction. Interestingly, at satu-ration the magnetic energy in the xy -plane is noticeably lessthan E b and yet the overall value of E b is the same as thatderived from the full-low-hr (i.e. quasi-ideal MHD) simula-tion. An interesting point to note about this is that, whereasit has been shown (Jones et al. 1997) that the strength ofthe field in the direction of the initial flow is what is impor-tant in determining the effect of the magnetic field on theKH instability, here we can see that the strength of the fieldperpendicular to this plane appears to play a role also. Following the initial linear growth of the instability, the sys-tem experiences a period of transferring energy back andforth between the magnetic field and velocity field, as in s / B - / B Figure 12.
The total perturbed magnetic energy (solid line) inthe system is plotted against time for hall-med-hr. Also plottedfor this simulation are the evolutions of the perturbed magneticenergies in the xy -plane only (dashed line) and in the z -direction(dot-dashed line). While the maximum reached by the total mag-netic energy is the same as in full-low-hr, the magnetic field hasexperienced some re-orientation into the z -direction. s / ρ v x and / B - / B Figure 13.
Plot of the total perturbed magnetic energy (dot-dashed lines) and transverse kinetic energy (solid lines) with timefor the full-low-hr (thin lines) and hall-med-hr (thick lines) sim-ulations. full-low-hr (see figure 13). Interestingly, the amplitude ofthe oscillations - i.e. the amount of energy being transferredbetween motion and the magnetic field - is larger in hall-med-hr than in full-low-hr. This may be due to re-orientationof the magnetic field out of the plane of the instability andhence a reduction in (numerical) reconnection. It is worth re-calling here that the Hall effect, although it appears similarto a diffusion term in the induction equation, is a dispersiveeffect which conserves magnetic energy.Perhaps the most important difference between full-low-hr and hall-med-hr is that while the kinetic energy in the y -direction gradually tends towards a constant value in thelow resistivity case, it continues to steadily decay in hall-med-hr. The conclusion is that in full-low-hr, the instabilityhas completed its growth and is returning to a quasi-steadystate. In hall-med-hr however, the instability is continuingto consume the parallel kinetic energy available to it andthe instability undergoes a further stage of development asa result of the inclusion of the Hall effect.Since the magnetic field gains a component in the z direction in hall-med-hr due to the Hall effect it is interest-ing to examine the behaviour of the kinetic energy in the z direction also. Figure 14 contains plots of the transversekinetic energy and the kinetic energy in the z direction (i.e. E k ,z ≡ ρv z ). Clearly the transverse kinetic energy grows c (cid:13) —- RAS, MNRAS000
Plot of the total perturbed magnetic energy (dot-dashed lines) and transverse kinetic energy (solid lines) with timefor the full-low-hr (thin lines) and hall-med-hr (thick lines) sim-ulations. full-low-hr (see figure 13). Interestingly, the amplitude ofthe oscillations - i.e. the amount of energy being transferredbetween motion and the magnetic field - is larger in hall-med-hr than in full-low-hr. This may be due to re-orientationof the magnetic field out of the plane of the instability andhence a reduction in (numerical) reconnection. It is worth re-calling here that the Hall effect, although it appears similarto a diffusion term in the induction equation, is a dispersiveeffect which conserves magnetic energy.Perhaps the most important difference between full-low-hr and hall-med-hr is that while the kinetic energy in the y -direction gradually tends towards a constant value in thelow resistivity case, it continues to steadily decay in hall-med-hr. The conclusion is that in full-low-hr, the instabilityhas completed its growth and is returning to a quasi-steadystate. In hall-med-hr however, the instability is continuingto consume the parallel kinetic energy available to it andthe instability undergoes a further stage of development asa result of the inclusion of the Hall effect.Since the magnetic field gains a component in the z direction in hall-med-hr due to the Hall effect it is interest-ing to examine the behaviour of the kinetic energy in the z direction also. Figure 14 contains plots of the transversekinetic energy and the kinetic energy in the z direction (i.e. E k ,z ≡ ρv z ). Clearly the transverse kinetic energy grows c (cid:13) —- RAS, MNRAS000 , 1–12 he Kelvin-Helmholtz instability in weakly ionised plasmas s / ρ v z and / ρ v x Figure 14.
Plot of the transverse kinetic energy (solid line) andkinetic energy in the z -direction (dashed line) with time for hall-med-hr. It can be seen that at later times, the growth of energyin the z -direction has become significant. rapidly during the linear development of the instability but E k ,z also grows and, eventually, even becomes larger than E k ,x . Figure 15 contains plots of the total energy (kineticand magnetic) in full-low-hr, hall-med-hr and hall-high-hr.Somewhat surprisingly we find that the hall-med-hr simula-tion loses energy somewhat faster than full-low-hr. This is,on the face of it, a little puzzling since the Hall effect doesnot itself dissipate magnetic energy. If we consider the hall-high-hr simulation we find that the total energy is roughlyconstant - it behaves roughly the same as hd-zero-hr. Infact, the magnetic energy in the hall-high-hr simulation be-haves qualitatively differently to the hall-med-hr. Figure 16contains plots of the perturbed magnetic energy as a func-tion of time for hall-high-hr, and full-low-hr. It is clear thatsomething dramatic is happening in the hall-high-hr simula-tion. To gain insight into this, in figure 17 we plot the totalperturbed magnetic energy and the magnetic energy in the xy -plane and in the z -direction for simulation hall-high-hr.The growth of the magnetic field in the z direction is signifi-cant, as we expect from a system with the Hall effect. In fact,the growth of the magnetic field in the xy -plane is faster inthe hall-med-hr case. This isn’t too surprising as, in order toincrease the Hall effect, the coupling between the electronsand the neutrals is much weaker than that between the ionsand neutrals. Hence, while the neutrals and ions generatethe usual cat’s eye vortex, the electrons (and the magneticfield) do not. Figure 19 demonstrates these morphologicaldifferences between the various fluids.Furthermore, in the case of full-low-hr, the magneticfield in the xy -plane grows to such an extent that it op-poses further wind-up of the fluid in the vortex. In the hall-high-hr case this does not happen as magnetic energy is re-distributed to the z component of the field which does notoppose this wind-up. Hence the vortex is not destroyed in thenon-linear regime in the hall-high-hr case. This is the maindifference between the hall-med-hr and hall-high-hr simula-tions. In the non-linear regime, then, the ion fluid remainsspinning in the KH vortex while the electron fluid remainstied to the magnetic field. This maintains a velocity differ-ence between the electrons and ions which, in turn, causesthe magnetic energy to be further re-distributed to the z component. In this way the Hall effect, if strong enough, in-troduces strong dynamo action into the KH instability. Thisdynamo behaviour, which is not possible in a 2.5D, ideal s -0.10-0.050.000.05 C hange i n t o t a l ene r g y s -0.10-0.050.000.05 C hange i n t o t a l ene r g y Figure 15.
The change in total energy of the system is plottedagainst time for hall-med-hr (thick line, upper panel) and hall-high-hr cases (thick line, lower panel). In each case the change intotal energy is plotted for full-low-hr (thin line) for comparisonand also for hd-zero-hr (dashed line). s / B - / B Figure 16.
The total perturbed magnetic energy in the systemis plotted against time for full-low-hr (thin line) and hall-high-hr(thick line) cases. The magnetic energy is seen to experience verydifferent growth in hall-high-hr.
MHD system is made possible by the Hall term introducinga handedness into the flow (e.g. Wardle 1999).If we follow the dynamics further we find that the KHinstability, in the presence of high Hall resistivity does notsaturate to a quasi-steady state as it does in, for example,the full-low-hr case. As the z components of the current andmagnetic field continue to grow, the Hall effect now actson the non-parallel currents and magnetic fields that havearisen between the z -directions and the xy -plane. This hasthe result of re-orienting some of the magnetic field, andelectron fluid flow, back onto the xy -plane. During this pro-cess the electron fluid obtains a velocity away from the KHvortex, which results in a broader volume of plasma beingdisturbed. This feeds the continuous growth of the magneticenergy in the xy -plane, and thus causes continuous growth c (cid:13) —- RAS, MNRAS , 1–12 A. C. Jones and T. P. Downes s / B - / B Figure 17.
The perturbed magnetic energy in the system is plot-ted against time for hall-high-hr (solid line). Also plotted for com-parison are the growth of the magnetic energies in the xy -plane(dashed line) and in the z -direction (dot-dashed line). It can beseen that there is significant growth in the z -direction.
15 16 17 x0.00.20.40.60.81.0 y Figure 18.
Plot of the magnitude and vector field of the magneticfield in the xy -plane in hall-high-hr at time 8 t s . It can be seenthat the magnetic field does not undergo as much wind-up as isseen in the bulk velocity field (upper panel of figure 19). of the electron transverse kinetic energy, as can be seen infigure 20. In this particular simulation, the size of the simu-lation grid is set so that it provides a sufficiently large regionof ordered flow that can be transformed into disordered flowduring the time of the simulation. The instability will nat-urally saturate when it has exhausted the area of orderedflow available to it.To summarise, through their strong decoupling from themagnetic field, the dynamics of the bulk fluid and ion fluiddemonstrate behaviour very similar to that of hd-zero-hr inwhich the KH vortex remains intact. The high Hall casecan in fact be thought of as two separate systems occur-ring simultaneously; the bulk fluid demonstrating hydrody-namic behaviour and the continuously widening volume ofperturbed fluid perpetually feeding the growth of the mag-netic field through the Hall effect. These two systems are in-trinsically entwined through the requirement of charge neu-trality, by which the electron fluid causes a widened areaof perturbed ion fluid, and thus the bulk fluid. Both sys-tems are relatively energy efficient and, following the initialgrowth of the KH instability, the overall system experienceslittle further loss of energy. The supply of energy to feedthe magnetic field is limited only to the physical size of thecomputational domain over which the simulation is run.This dynamo action occurs only under certain condi-tions. As the initial Hall resistivity is increased from moder-
15 16 17 x0.00.20.40.60.81.0 y y y Figure 19.
Plot of the magnitude and vector field of the bulkvelocity (upper panel), ion velocity (middle panel) and electronfield (lower panel) for hall-high-hr. The dynamics of the electronfluid can be seen to be very different to the other two fluids, andmore similar to the magnetic field (figure 18). s -18 -18 -18 -18 -18 / ( ρ v x ) e l e c t r on Figure 20.
The transverse kinetic energy of the electron fluidis plotted against time for full-low-hr (thin line) and hall-high-hr(thick line). Also plotted in the electron transverse kinetic energyfor hd-zero-hr (dashed line). c (cid:13) —- RAS, MNRAS000
The transverse kinetic energy of the electron fluidis plotted against time for full-low-hr (thin line) and hall-high-hr(thick line). Also plotted in the electron transverse kinetic energyfor hd-zero-hr (dashed line). c (cid:13) —- RAS, MNRAS000 , 1–12 he Kelvin-Helmholtz instability in weakly ionised plasmas ate to high, the electron fluid, and thus the magnetic field,becomes increasingly decoupled from the neutrals. As a re-sult, the neutral fluid tends toward hydrodynamic behavior.This dynamo action is observed only when the Hall resis-tivity is increased to a sufficiently high value that the KHvortex formed by the neutral fluid is no longer constrainedby the magnetic field. Unlike the pure hydrodynamic case,the presence of charged fluids undergoing different dynamicsleads to the dynamo action observed. If the Hall resistivity isnot sufficiently large, the magnetic field eventually leads tothe destruction of the KH vortex in the neutral fluid throughreconnection as in the ideal MHD case.Even a moderate amount of Hall resistivity results ina wider volume of fluid being disturbed by the instability.This agrees with previous studies of the Hall effect on theKH instability (e.g. Huba 1994). The re-orientation of themagnetic field lines within the KH instability has also beenobserved in studies of MRI in accretion disks (e.g. Wardle1999). Kunz (2008) investigated a simple model of accretiondisks without rotation and demonstrated that the combinedactions of the shear instability and the Hall effect leads toincreased stretching of the magnetic field lines. This wasshown to result in continued growth of the instability, whichcorresponds well to the dynamo action observed in our sim-ulations with high Hall resistivity. It is important to note,though, that studies by both Kunz (2008) and Wardle (1999)are linear studies and, as such, do not extend to the nonlin-ear regime of the instability. A more recent numerical studyby Nykyri & Otto (2004) demonstrates the twisting of mag-netic field lines, but due to the inclusion of magnetic re-connection, doesn’t produce the magnetic dynamo observedhere. We have presented the results of a suite of fully multifluidMHD simulations of the KH instability in weakly ionised flu-ids such as, for example, molecular cloud material. Throughvarying the collision coefficients between the various chargedspecies and the neutrals we were able to investigate systemsin which ambipolar diffusion dominates the multifluid effectsand ones in which the Hall effect dominates. We validatedour KH simulations through comparison of an ideal MHDsimulation with previously published results and performedresolution studies for each of these cases to ensure that ourconclusions are not unduly effected by our numerical reso-lution.Our findings can be summarised as follows: • The multifluid effects do not significantly influence thelinear growth rate of the instability. • Ambipolar diffusion dramatically reduces the energy as-sociated with the perturbed magnetic field. This happensthrough diffusion for moderate ambipolar resistivity, butthrough both diffusion and decoupling of the magnetic fieldfrom the bulk flow for higher resistivity. • The Hall effect, as expected, rotates the magnetic fieldout of the initial xy -plane. • In contrast to both ambipolar dominated and idealMHD flows, the Hall effect causes the system to fail to settleto a quasi-steady state after saturation of the instability. • For moderate Hall resistivity the perturbed magneticfield contains higher energy than in the ideal MHD case,presumably due to a field topology which impedes (numeri-cal) reconnection. • For high Hall resistivity strong dynamo action is seenas energy associated with the magnetic field grows withoutany apparent signs of saturation.
ACKNOWLEDGEMENTS
The research of A.C.J. has been part supported by the Cos-moGrid project funded under the Programme for Researchin Third Level Institutions (PRTLI) administered by theIrish Higher Education Authority under the National Devel-opment Plan and with partial support from the EuropeanRegional Development Fund.The authors wish to acknowledge the SFI/HEA IrishCentre for High-End Computing (ICHEC) for the provisionof computational facilities and support.
REFERENCES
Amerstorfer U.V., Erkaev N.V., Langmayr D., BiernatH.K., 2007, P&SS, 55, 1811Baty H., Keppens R., 2006, A&A, 447, 9Bettarini L., Landi S., Rappazzo F.A., Velli M., Opher M.,2006, A&A, 452, 321Birk G.T., Wiechen H., 2002, PhPl, 9, 964Bucciantini N., Del Zanna L., 2006, A&A, 454, 393Chandrasekhar, S., 1961, Hydrodynamic and hydromag-netic stability. Clarendon, OxfordCiolek G.E., Roberge W.G., 2002, ApJ, 567, 947Dedner A., Kemm F., Kr¨oner D., Munz C.-D., SchnitzerT., Wesenberg M., 2002, JCP, 175, 645Downes T.P., Ray T.P., 1998, A&A, 331, 1130Downes T.P., O’Sullivan S., 2009, ApJ, 701, 1258Downes T.P., O’Sullivan S., 2011, ApJ, In press.Dyson J.E., 1984, Ap & SS, 106, 181Faganello M., Califano F., Pegoraro F., 2009, NJPh, 11,3008Frank A., Jones T.W., Ryu D., Gaalaas J.B., 1996, ApJ,460, 777G´omez G.C., Ostriker E.C., 2005, ApJ, 630, 1093Hardee P., Stone J., Rosen A., 1997, IAUS, 182, 132Hasegawa H., Fujimoto M., Phan T.-D., R`eme H., BaloghA., Dunlop M.W., Hashimoto C., TanDokoro R., 2004,Nature, 430, 755Huba J.D., 1994, PhRvL, 72, 2033Jones, A.C., PhD thesis, Dublin City UniversityJones T.W., Gaalaas J.B., Ryu D., Frank A., 1997, ApJ,482, 230Johansen A., Henning T., Klahr H., 2006, ApJ, 643, 1219Keppens R., T´oth G., Westermann R.H.J., Goedbloed J.P.,1999, JPlPh, 61, 1Kunz M.W., 2008, MNRAS, 385, 1494Li L.-X., Narayan R., 2004, ApJ, 601, 414Luo Q.Z., D’Angelo N., Merlino R.L., 2001, PhPl, 8, 31Malagoli A., Bodo G., Rosner R., 1996, ApJ, 456, 708Miura A., Pritchett P.L., 1982, JGR, 87, 7431Norman C., Silk J., 1980, ApJ, 238, 158 c (cid:13) —- RAS, MNRAS , 1–12 A. C. Jones and T. P. Downes
Nykyri, K., Otto, A., 2004, Annales Geophysicae, 22, 935Oishi, J.S., Mac Low, M.-M. 2006, ApJ, 638, 281O’Sullivan S., Downes T. P., 2006, MNRAS, 366, 1329O’Sullivan S., Downes T. P., 2006, MNRAS, 376, 1648Palotti M.L., Heitsch F., Zweibel E.G., Huang Y.-M., 2008,ApJ, 678, 234Shadmehri M., Downes T.P., 2008, MNRAS, 387, 1318Umebayashi T., Nakano T., 1990, MNRAS, 243, 103Venter L.A., Meintjes P.J., 2006, MNRAS, 366, 557Wardle M., 1999, MNRAS, 307, 849Wardle M., 2004, Ap&SS, 292, 317Watson C., Zweibel E.G., Heitsch F., Churchwell E., 2004,ApJ, 608, 274Wiechen, H.M., 2006, PhPl, 13, 2104 c (cid:13) —- RAS, MNRAS000