Reduced drift-kinetics with thermal velocity distribution across magnetic field MykolaGordovskyy and PhilippaBrowning 6 ∗ 1 JodrellBank CentreforAstrophysics,UniversityofManchester, ManchesterM139PL, UK 0 2 n Abstract. The goal of this study is to develop an approximate self-consistent description of particle a J motion in strongly magnetised solar corona. We derive a set of reduced drift-kinetic equations based on 1 theassumption thatthegyro-velocity distribution isMaxwellian. Theequations aretested using simple1D 3 models. ] R 1 Introduction S . h p Generally, solving even the 2D drift-kinetic problem would require at least a 4-dimensional phase space - (2D2V). This is practically difficult: taking into account a number of processors normally available (up to o r 100) and a ‘reasonable’ wall-time (up to 3-5 days), it would be possible to consider a domain with up to st ∼107-108gridpoints. Realistically,thephase-spacecanhave3dimensionswithareasonableresolution(e.g. a ∼ 2D1V), or 4 dimensions (e.g. 2D2V)with very scarce resolution along at least one of dimensions. Hence, [ fullkinetic(ordrift-kinetic) treatmentwouldbetoonumerically expensiveincaseof2Dor3Dgeometry. 1 Variousanalyticalstudiesandnumericalsimulations(e.g.Gordovskyy etal.,2010;Gordovskyy&Browning, v 1 2011) show that particle acceleration by quasi-stationary electric fields in the solar corona affects the par- 4 allel component of particle velocity v = ~v ~b (where~b = B~/B is the magnetic field direction), while the 3 || · gyro-velocities v normally remain nearly thermal. This is natural, taking into account that in the corona 0 g 0 magnetic field curvature is small compared to particle Larmor radii and collisional times are longer than . 2 acceleration times, i.e. there arenostrong scattering mechanisms, andaccelerated particles areexpected to 0 remaincollimated alongmagneticfield. 6 Kinetic description of particles with small pitch-angles (i.e. v v ) has been discussed in number 1 g : ofpapers inthe laboratory plasma context (e.g. Pfefferleetal.,2015||).≫Generally, anapproximation of zero v i Larmor velocities, at least, for energetic particles, would substantially simplify the kinetic equations, re- X moving a dimension from the phase space. However, in context of large-scale particle kinetics in the solar r a corona, Larmor radii, even small, may play an important role in some cases, for instance, during particle mirroring from the strongly converging magnetic field at the bottom of the corona. Therefore, it might be morerealistic toassume thatparticle gyro-velocities are, generally, non-zero, buttheir distributions remain Maxwellian. Assuming the gyro-velocities always have Maxwellian distribution, one could reduce the phase space by looking for a distribution function F(~r,v ;t) and the perpendicular temperature τ(~r,v ;t) instead of the || || distribution function f(~r,~v;t). Here, the perpendicular temperature defines the width of the Maxwellian distribution of gyro-velocities for a given particle specii with given parallel velocity at a given location. Hence,usingthisformalismtheproblemcanbereducedtocalculationofF(x,y,v ;t)andτ(x,y,v ;t)in2D || || case, andtocalculation of F(x,y,z,v ;t)and τ(x,y,z,v ;t)in3Dcase. Belowthis formalism willbecalled || || ‘reduced kinetics’. e-mail:mykola.gordovskyy[AT]manchester.ac.uk ∗ 1 We derive a set of reduced drift kinetic equations by averaging the Larmor gyration velocity at each location of the phase space (~r,v ). Theaveraging is done assuming that the distribution in respect of gyro- || velocity always remains Maxwellian. Theresulting equations are required toconserve the particle number andenergy. 2 Reduced equation derivation 2.1 Full kineticequation Considerakineticequation inthefollowingform: ∂f ∂f dv ∂f dv2 ∂f Lˆf = + ~u+v ~b + || + g = 0, (1) ∂t || ∂~r dt ∂v dt ∂v2 (cid:16) (cid:17) g || where the distribution function is f = f(~r,v ,v2;t), ~u is the guiding centre drift velocity, ~b = B~/B is the g || magnetic field direction vector. This equation is adopted from astandard form of drift-kinetic equation for particles withnon-zero magneticmoments(Kulsrud,1983,seee.g.). Thedriftvelocityconsistsofthefollowingterms: ~u=~u +~u , (2) ∗ B ∇ where~u =~u +~u isthesumofExBandcurvaturedrifts ∗ E C E~ B~ ~u = × , E B2 mv2 ~uC = ||[~b (~b ~)~b], qB × ·∇ whichdon’tdependonthegyrovelocity, and mB~ ~B ~u = ×∇ v2. ∇B q 2B2 g Parallel velocity canbeaffected bytheparallel electric fieldandbymagnetic fieldgradient along mag- neticfieldlines: dv q ~B ~b || = ~ ~b v2∇ · , (3) dt mE· − g 2B where ~iselectric field,and Band~baretheabsolute valueanddirection (~b = B~/B)ofthemagneticfield. E Variation of gyro-velocity can be derived from the magnetic moment conservation v2/B = const. Dif- g ferentiating thisinrespect oftime,andsubstituting dB = ∂B +(v ~b+~u) ~Bgives dt ∂t || ·∇ dv2 ~B ~b ~B ~u g = ∇ · v2v + ∇ · v2. (4) dt B g || B g The first term in the RHS corresponds to the magnetic mirroring effect, so that 2v dv + dv2g = 2v dv + dt|| dt dt|| || || ~B~bv2v = 0. (We ignore terms containing ∂E and ∂B, assuming that field variation timescale is much ∇B· g ∂t ∂t || longerthan1/ω ofconsidered particlespecies.) g 2 2.2 Integrated distribution function and averagegyro-velocity In terms of the parallel velocity v and squared gyro-velocity v2 the thermal distrbution with total specific g + + || + + energy (~r)= ∞ ∞(v2+v2)Fdv2dv andtotalparticle number (~r)= ∞ ∞Fdv2dv is Eth R R0 || g g || Nth R R0 g || −∞ v2+v2 −∞ where f = th andthe”equivalefntht(tve||m,vp2ge)ra=tufrteh”0eτxp=−2 ||τtht.h g, (5) th0 √Nπτ3/2 th 3 Eth th N Weassumethatthedistribution function canbewritteninthefollowingform: v2 v2 f(~r,v ,v2;t)= P(~r,v ;t) 0 exp g = P(~r,v ;t) (v2,τ(~r,v ;t)), (6) || g || τ −τ(~r,v ;t) || S g || wherevg isgyro-velocity andv0 issomecharacteristic co||nstant velocity. Letusintroduce anewdistribution function integratedinrespectofv2 g ∞ v2 ∞ v2 F(~r,v ;t) = f(~r,v ,v2;t)dv2 = P 0 exp g dv2. (7) || Z || g g τ Z −τ(~r,v ;t) g 0 0  ||  Theoriginalkineticequation1cannotbeexactlyintegratedinrespectofv2ingeneralcasebecausesome g coefficients depend on the gyro-velocity. The idea is to substitute the gyro-velocity by the ”perpendicular temperature”, whichisthe average gyro-velocity. Using thedistribution function form 6,itiseasytoshow that ∞v2f(~r,v ,v2;t)dv2 g g g || v2 = R0 = τ(~r,v ;t). (8) h gi ∞f(~r,v ,v2;t)dv2 || g g || R0 2.3 Integrated kinetic equation Hereweintegrate thekinetic equation, eacheffectisconsidered separately. Firstly,severaltermsintheequation1don’tdependonv2 andtheirintegration istrivial: g ∞ ∂f ∞ ∂f ∞ q ∂f dv2 + (~u +v ~b) dv2 + ~ ~b dv2 = Z ∂t! g Z ∗ || ∂~r g Z mE· ∂v g 1 0 0 0 || ∂F ∂F q ∂F +(~u +v ~b) + ~ ~b . (9) ∗ ∂t ! || ∂~r mE· ∂v 1 || Next,weusetheaveragesquared gyro-velocity fortheperpendicular Bdrift: ∇ ∞ ∂f ∞ ∂f dv2 + ~u dv2 Z ∂t! g Z ∇B∂~r g ≈ 2 0 0 ∞ ∂f ∞ ∂f ∂F ∂F dv2 +U~ dv2 = +U~ , (10) Z ∂t! g ∇BZ ∂~r g ∂t ! ∇B ∂~r 2 2 0 0 3 where mB~ ~B U~ = ×∇ τ. (11) ∇B q 2B2 Finally,weuseaveragesquaredvelocity todescribe thev variation duetothemagneticmirroring: || ∞ ∂f ∞ ~B ~b ∂f dv2 + v2∇ · dv2 Z ∂t! g Z − g 2B ∂v g ≈ 0 3 0   || ∞ ∂f ∞ ~B ~b ∂f ∂F 1 ∂F dv2 + τ∇ · dv2 = τ , (12) Z ∂t! g Z − 2B ∂v g ∂t ! − 2G ∂v 0 3 0   || 3 || where ~B ~b = ∇ · . (13) G B Integrating the last term containing the ∂f in respect of v2 is, obviously, zero. Therefore, one can write an ∂v2 g g approximate drift-kinetic equation basedonthev -averaging: g ∂F ∂F q ∂F 1 ∂F Lˆ F = +(~u +U~ +v ~b) + ~ ~b τ = 0. (14) R ∗ B ∂t ∇ || ∂~r mE· ∂v − 2G ∂v || || 2.4 Energy equation Nowweneedanequationgoverningtheevolutionofthe”perpendicular temperature”τ. Sinceitisassumed thatthemagneticmomentv2/Bisconserved, similartoapproximations above,weassumethatthe”average g moment”isconserved aswell ∞v2g f(~r,v ,v2;t)dv2 2B g g || τ = R0 = . M 2B ∞f(~r,v ,v2;t)dv2 g g || R0 Mathematically, thiscanbewrittenas 1 τ Lˆ F = Lˆ F = 0. (15) R R M 2 (cid:18)B (cid:19) Expanding theaboveequationyields ∂F τ ∂F τ q ∂F τ 1 ∂F τ +(~u +U~ +v ~b) + ~ ~b τ + ∗ B ∂t B ∇ || ∂~r B mE· ∂v B − 2G ∂v B || || ∂τF ∂τF Fτ q ∂τ F 1 ∂τ F +(~u +U~ +v ~b) (~u +U~ +v ~b)~B + ~ ~b τ = 0. ∂t B ∗ ∇B || ∂~r B − ∗ ∇B || ∇ B2 mE· ∂v B − 2G ∂v B || || Substracting τLˆ F andmultiplying by B yieldsthefollowing: B R F ∂τ ∂τ τ q ∂τ 1 ∂τ +(~u +U~ +v ~b) (~u +U~ +v ~b)~B + ~ ~b τ = 0. (16) ∗ B ∗ B ∂t ∇ || ∂~r − ∇ || ∇ B mE· ∂v − 2G ∂v || || 4 2.5 Fieldequations Electricandmagneticfieldevolution canbedescribed usingtheMaxwellequations ∂B~ = ~ E~ (17) ∂t −∇× ∂E~ 1 1 = ~ B~ ~j (18) ∂t ǫ µ ∇× − ǫ 0 0 0 ρ ~ E~ = (19) ∇· ǫ 0 ~ B~ = 0, (20) ∇· where ρ = Σsqs F(~r,v ;t)dv  and ZV || || ~j= Σsqs (~u∗+U~ B+v ~b)F(~r,v ;t)dv . ZV ∇ || || || 3 Numerical tests Hereweinvestigatethemagneticmirroring,whichwouldproducethehighestsystematicerrorinreducedki- neticapproximation. Theresultsofsimplenumericaltestscomparingthefullkineticsolutionswithreduced kinetic solutions in stationary magnetic field are shown below. The magnetic field distribution is shown in Figure 1. The first test (Figure 2) shows the evolution of the distribution functions for particles moving through the weak magnetic mirrors (i.e. most particles are in the loss cone), the second test (Figure3) is for particles oscillating between two strong mirrors (i.e. most particles are outside the loss cone). All test modelshereareone-dimensional, withparticlesmovingalongtheaxesofcylindricallysymmetricmagnetic configurations (i.e., B~ ~B= 0). ×∇ Itcanbeseenthatthereducedkineticsolutionissimilartothe’fullkinetic’solution,althoughtheformer ismorecompact. Thisisnotsurprising,becauseoftheintroducedaveragingoverthegyro-velocity: asingle valueofdv /dtresults inalowerdispersion inv and,hence,in x. || || References GordovskyyM.,BrowningP.K.&Vekstein, G.E.,2010,Astrophys. J GordovskyyM.&BrowningP.K.,2011, Astrophys.J KulsrudR.,1983,HandbookofPlasmaPhysics,eds.M.N.RosenbluthandR.Z.Sagdeev,v.1: BasicPlasma PhysicsI. Pfefferle,D.,Graves,J.P.,&Cooper, W.A.,2015,PlasmaPhys.Contr.Fus.,57,54017 5 Figure 1: Magnetic field distribution in the 1D test simulations. The value of magnetic field at the bound- aries,max,isequal1inthemodelwithweakconvergence and32inthemodelwithstrongconvergence. 6 Figure 2: Particle number versus position (left panels) and velocity (right panels). The system is periodic in x. The magnetic mirrors are located at the left and right boundaries. Magnetic convergence ratio is 2, magnetic ’cork’ thickness is 16. The initial particle velocity is 4, the initial velocity dispersion is 2, and the initial gyro-temperature is 4 (corresponding to the velocity dispersion of 2). Panels from top tobottom correspond totime0,32,64,128. 7 Figure3: SameasinFigure1,butthemagneticconvergence is32. 8