ResearchinAstron.Astrophys. 20xxVol.9No.XX,000–000 R esearchin http://www.raa-journal.org http://www.iop.org/journals/raa A stronomyand A strophysics Detailed Analysis of Fan-Shaped jets in Three Dimensional Numerical Simulation 1 1 0 R.L.Jiang1,2,4,K.Shibata2,H.Isobe3 andC.Fang1,4 2 1 DepartmentofAstronomy,NanjingUniversity,Nanjing210093,China;[email protected] n 2 KwasanandHidaObservatories,KyotoUniversity,Yamashina,Kyoto607-8471,Japan a J 3 UnitofsynergeticStudiesforSpace,KyotoUniversity,Yamashina,Kyoto607-8471,Japan 4 4 KeyLaboratoryofModernAstronomyandAstrophysics(NanjingUniversity),Ministryof 2 Education,ChinaReceived[year][month][day];accepted[year][month][day] ] R S Abstract We performedthreedimensionalresistivemagnetohydrodynamicsimulations . to study the magneticreconnectionusingan initially shearingmagneticfield configura- h tion (force free field with a currentsheet in the middle of the computationalbox). It is p - shown that there are two types of reconnection jets: the ordinary reconnection jets and o fan-shapedjets,whichareformedalongtheguidemagneticfield.Thefan-shapedjetsare r muchdifferentfromtheordinaryreconnectionjetswhichareejectedbymagnetictension t s force. There are two driving forces for accelerating the fan-shaped jets. The one is the a [ Lorentzforcewhichdominatesthemotionoffluidelementsatfirstandthenthegaspres- sure gradientforce accelerates the fluid elements in the later stage. The dependenceon 1 magnetic reconnectionangle and resistivity value has also been studied. The formation v and evolution of these jets provide a new understanding of dynamic magnetohydrody- 8 namicjets. 9 5 4 Keywords: magnetohydrodynamics:MHD—methods:numerical . 1 0 1 1 INTRODUCTION 1 : Magneticreconnectionplaysaveryimportantroleinsolarflares,coronalmassejections,andothersolar v activities.Theproblemofflareenergyreleasehasbeenpuzzlingpeopleformanyyearsbeforethedevel- i X opmentofmagneticreconnectionmechanism.Thefirstpioneersinthisfieldofmagneticreconnection areSweet(1958)andParker(1957),whodevelopedatheorywellknownasSweet-Parkermechanism. r a However,it can notexplainthe solar flares. The biggestproblemin their theoryis thatthe time scale of energy release is much longer than the realistic one of flares. Several years later, Petschek (1964) improvedthis work by consideringa pair of slow mode shocksoutside a verysmall diffusionregion, whichgreatlyincreasesthereconnectionrate.Thefastreconnectionmodelissuccessful,thoughsome peopledonotbelievethatthelargeenergyreleaseprocesshassuchasmalldiffusionregion.Thefamous CSHKP model(Carmichael1964; Sturrock 1966; Hirayama 1974; Kopp&Pneuman 1976) based on magneticreconnectionforflareshasbeendevelopedsince1960s. Besides observations, magnetohydrodynamic(MHD) simulations play a key role in studying the magnetic reconnection (Ugai&Tsuda 1977; Forbes&Priest 1983; Forbes 1990; Forbes&Malherbe 1991;Yokoyama&Shibata1994;Magaraetal.1996;Ugai1996;Chenetal.1999).Withthedevelop- ment of computerand observationtechnology,the CSHKP flare modelhas been extended. More and ⋆E-mail:[email protected] 2 R.L.Jiang,K.Shibata,H.Isobe&C.Fang moresimulationsshowthattheplasmoidejectionmaygreatlyincreasethemagneticreconnectionrate (Shibataetal. 1995; Shibata 1996, 1997; Magaraetal. 1997; Nishidaetal. 2009). The plasmoidscan forminananti-parallelmagneticfield.Beforeejectionmagneticenergyisrestoredinthediffusionre- gion,andsimultaneouslytheplasmoidsmergewitheachother.Finally,therestoredmagneticenergyis releasedinaveryshorttimescaleaftertheejectionofthebigplasmoid. Recently,thethreedimensional(3D)simulationisoneofthehottestpointsinastrophysics.Inthis case,thereconnectingprocessismuchmorecomplicatedthanthatintwodimensional(2D)case.The3D simulationscanprovideusverycomplexandrealisticstructures(Yokoyama&Shibata1995;Isobeetal. 2005;Ugai&Zheng2005;Isobeetal.2006;Shimizuetal.2009).Itisalsoagoodtooltoexplainthe mechanismofmanysmallscalesolaractivities(Cirtainetal.2007;Shibataetal.2007;Katsukawaetal. 2007)whichhavebeenobservedbyHinode(Kosugietal.2007).Someoftheseobservedjetlikefeatures are believedto be a result of magnetic reconnectionin the so-called interlocking-comblike magnetic configurationin sunsportpenumbra(i.e.,penumbraljets). Themotivationofthispaperisto studythe effect of the guide field using 3D numericalsimulations. In our new simulation, we found that some partsofthereconnectionjetswhichareejectedfromthediffusionregioncanmovealongthemagnetic guidefieldlinesintheinterlocking-comblikemagneticconfiguration.Moreover,thereconnectionjets whichmovealongtheguidefieldlinesarealmostperpendiculartotheambientmagneticfield.Inorder todistinguishthejetsmovingindifferentdirections,wecallthejetsmovingalongtheambientmagnetic fieldasordinaryreconnectionjetsandthejetsmovingalongtheguidefieldasfan-shapedjets(wewill seetheshapeofthejetsissimilartoafaninthenextsection).Besidesoursimulations,twoveryrecent paperbyPontinetal.(2005)andUgai(2010)alsodiscussedtheeffectofguidefield.Theirsimulations showedthattheguidefieldcandistortthepropagationof3Dplasmoidsandtheplasmaisalmostejected alongtheguidefield.Itseemsthattheguidefieldplaysaveryimportantroleforformingthesejetsin the3Dmagneticreconnection.Inourfirstpaper(Jiangetal.2011)webrieflyreportedthemainresults of our 3D simulation, while in this paper the detailed analysis and discussion of fan-shaped jets are presented. This paper is organizedas follows. A detailed description of the magnetohydrodynamics(MHD) equationsandtheinitialconditionaregiveninthenextsections.InSection3weshowthe3D,2Dand 1Dstructuresofourtypicalsimulationresults.Moreover,theanalysisoftheLagrangianfluidelements, reconnectionrateandslowmodeshockarealsogiven.Section4showsdependenceondifferentrecon- nectionangles(definedintheinitialconditon)anddifferentvaluesofresistivityinthediffusionregion. Finally,discussionandsummaryaregiveninthelastSection. 2 NUMERICALMETHOD Wesolvethethree-dimensional,nonlinear,resistive,compressibleMHDequationsnumerically.Therole ofgravity,thermalconduction,radiationorsomeothereffectsarenotsoimportantinoursimulations, becausewhatwearereallyinterestedinisthedynamicsofthemagneticreconnection.Theadditional considerationmayaffectthefinalresultsbutnotessential.ThesimplifiedMHDequationsaregivenin cgsunitsasfollows: ∂ρ +v·∇ρ=−ρ∇·v, (1) ∂t ∂v ∇p 1 1 2 +v·∇v=− − ∇B + (B·∇)B, (2) ∂t ρ 8πρ 4πρ ∂p 2 +v·∇p=−γp∇·v+(γ−1)ηj , (3) ∂t ∂B =−c∇×E, (4) ∂t 1 E=ηj− v×B, (5) c Fan-Shapedjets 3 c j= ∇×B, (6) 4π κ B p= ρT, (7) m where the eight independent variables are the velocity (υ ,υ ,υ ), density (ρ), magnetic field x y z (B ,B ,B ), andpressure(p).η isthe resistivity,jthecurrentdensity,T the temperature,andE the x y z electricfield.Thespecificheatratioofmonoatomicperfectgasisγ =5/3.Thenormalizedunitsofthe variablesare listed in Table 1. Hereafter,all the valuesused in ouranalysis anddiscussion arein non dimensionalformforsimplicity. Table1 NormalizationUnits. Variable Quantity Unit x,y,z Length L0 T Temperature T0 ρ Density ρ0 p Pressure p0 =ρ0T0κB/m t Time t0 =L0/v0 η Resistivity η0 =(1/c2)L0v0 v Velocity v0 =(p0/ρ0)1/2 B Magneticfield B0 =p10/2 E Electricfield E0 =(1/c)v0B0 j Currentdensity j0 =(c)B0/L0 F Forceperunitmass F0 =p0/L0 vA Alfve´nvelocity vA0 =B0/(ρ0)1/2 τA Alfve´ntimescale τA0=L0/vA EA Alfve´nelectricfield EA0=(1/c)vAB0 Notes: Inthistable,m,κ andcmeanmolecularmass,Boltzmannconstantandlightspeed,respectively. B Fig.1 Initialcondition,andtheboxsizeis−4 < x < 4,−16 < y < 16,−16 < z < 16in nondimensionalunit. ThecomputationalboxisshowninFigure1.Thedomainofthisboxis−4<x<4,−16<y <16 and−16<z <16inallofourcases(Notethattheimagesofsimulationresultsonlyshowthecentral 4 R.L.Jiang,K.Shibata,H.Isobe&C.Fang partofthe computationalboxand the scale forx, y andz coordinatesis differentbecause ofthe thin currentsheetinxdirection).Themostimportantconditionistheconfigurationoftheforcefreemagnetic fieldwhichisshownintherightpanelofFigure1.Themathematicalformofthisfieldis B =0, (8) x −B for x<−△h ini x B = B sin{[θx−(π−θ)△h ]/2△h } for |x|<△h , (9) y ini x x x B sin[(2θ−π)/2] for x>△h ini x 0 for x<−△h x B = B cos{[θx−(π−θ)△h ]/2△h } for |x|<△h , (10) z ini x x x B cos[(2θ−π)/2] for x>△h ini x whereB isaconstant,whichmeansthemagneticfieldhasaninitialuniformstrength.Thedirection ini ofthisfieldisshowninFigure1.The△h isthehalfwidthofthecurrentsheet.θ isthereconnection x angleshowninFigure1,andithasbeentakentobeoneofthefourvaluesinoursimulationcases,i.e., π,3π/4,π/2andπ/3,aslistedinTable2. Table2 ComputationalCases. Case β η Reynoldsnumber(v L /η ) Reconnectionangle(θ) Totalgridpoints ini A x ini 1 1.0 0.025 226 π 300×240×240 2 0.8 0.050 126 π 240×480×480 3 0.8 0.050 126 3π/4 240×480×480 4 0.8 0.050 126 π/2 240×480×480 5 0.8 0.050 126 π/3 240×480×480 6 1.0 0.025 226 π 240×480×480 7 1.0 0.050 113 π 240×480×480 8 1.0 0.075 75 π 240×480×480 Notes: Inthistable,Lxisthehalflengthofthecomputationalboxsizeinxdirection. In order to get fast reconnection process, we put a localized diffusion region at the center of the currentsheet(Ugai1986;Yokoyama&Shibata1994).Theadoptedfunctionalformfortheanomalous resistivityisasfollows: xπ yπ zπ η =η cos[ ]cos[ ]cos[ ] for |x|<△h , |y|<△h , |z|<△h , (11) ini x y z 2△h 2△h 2△h x y z where△h = 0.4(itisthe samevalueasthe widthofthecurrentsheet),△h = 1.6,△h = 1.6is x y z thehalfwidthofthediffusionregioninx,y andz direction,respectively.η isthemaximumofthe ini resistivityinthediffusionregionandistakenasaconstant. In the numerical procedures, we solve the MHD equations using CIP-MOCCT scheme (Yabe&Aoki1991;Yabeetal.1991;Kudohetal.1999)withanartificialviscosity(Stone&Norman 1992) and libraries of CANS (Coordinated Astronomical Numerical Softwares). The calculation do- main is discretizedinto 240×480×480nonuniformgridsin mostof cases. The minimumgridsize is △x = 0.00667,△y = 0.0667and △z = 0.0667.Since there is nogravityin our simulations,the freeboundaryconditionissimplyappliedinallofsixboundaries(namely,anequivalentextrapolation was applied for all quantities). The non-divergencecondition of the magnetic field is guaranteed by theconstrainttransport(CT)algorithmofEvans&Hawley(1988)inCIP-MOCCTscheme.Forinitial condition,thedensityandpressureareuniform(ρ = 1,p = 1)inthecomputationalboxandtheother parametersareshowninTable2. Fan-Shapedjets 5 3 TYPICALCASE Our typical case is the one in which the reconnection angle is π with a plasma beta β = 0.8 and η = 0.050 (case 2 in Table 2). Figure 2 shows the 3D visualization of gas pressure distribution at ini different times. The bright colors and high opacity stand for the high values. The volume rendering in the 3D visualization does not show the low gas pressure volumes which exist everywhere in the computationalbox. Magnetic field lines are drawn from Lagrangianfluid elements located outside of thediffusionregion,andthesame fieldlinesareshownin theupperpanelsofFigure2. We cansee a fan-shapedstructure in the currentsheet. For the lowerpanels, no magnetic field lines are plotted but we adda cross section near z = 0 and a velocity field plane at x = 0. The velocityarrowsshow the fan-shapedstructurealmostalongtheguidefieldlines(thereisnolegendforthevelocityfield,thescale forthevelocityisuselessinthese3Dimagesbecauseoftheprojectioneffect).Thecrosssectionlocated nearz =0.0showstheordinarymagneticreconnectionjetsinthisplanewhichisverysimilartothe2D or 2.5D (two dimension with three components) simulations (Chenetal. 2001; Yokoyama&Shibata 2001; Jiangetal. 2010). We will discuss this in detail in the next paragraph. Figure 3 shows the 3D resultsatdifferentviewpointsatthetime= 11.5.Fromthisfigurewecanseethefan-shapedstructure jetsmoreclearly. The2DcrosssectionsatdifferentpositionsanddifferenttimesareshowninFigures4-6.Figure4 depictsthegaspressure,magneticpressureandtotalpressuredistributionsinx−yplaneatz =0.0(we canalsobrieflyunderstandthepositionbythesmallbosabovethecolorbarinthisfigure).Theguide fieldisperpendiculartothisplane,whichmeanswhatwepresentinthisplaneissimilartosimulations in 2D or 2.5D case as we mentioned above. As reconnection goes on, the central diffusion region is heated by Joule dissipation, and a pair of ordinary reconnection jets are ejected by magnetic tension force.Thecentrallowmagneticpressureregionisduetothemagneticreconnectionprocess,whichcan becomparedtothehighgaspressureregion.Forasteadystate,thetwopressuresarealmostbalanced witheachotherinthexdirectionoftheinflowregionshownbythetotalpressurein thelowerpanels of Figure 4. The X-shaped structures (we can find them in the panels showing the gas and magnetic pressures)almostdisappear. Figure5showsthegaspressure,magneticpressureandtotalpressuredistributionsiny−zplaneat x=0.0.Inthisplane,theordinaryreconnectionjetsarealmostperpendiculartothemagneticfieldlines (inotherwords,movingalongthey and−y directions)andthejetsalmostmovingalongthemagnetic fieldlines(alongthezand−zdirections)areshownonlyinacross-sectionofthefan-shapedjets.The velocityof the fan-shapedjets is smaller than thatof the ordinaryreconnectionjets whose velocityis approximatelytheAlfve´nspeed(about1.58).Highmagneticpressureregionshavealsobeenejectedas shown in the magnetic pressuredistribution imagesin Figures4 and 5 (high magneticpressure parts, yellowcolor).ThetemperatureanddensitydistributionsinthisplaneareshowninFigure6.Similarto thegasandmagneticpressuredistributions,thedistributionsofthedensityandtemperatureareroughly oppositewitheachother. Figures7 and 8 depictpressure distributionsin y −z plane scanningat x = −0.1, x = 0.0 and x=0.1atthetime=8.5and11.5,respectively.Becausethedomainsofthesefiguresareinthecurrent sheet,wecanseethatthemagneticfieldisshearingbetweendifferentcolumns.Theinterestingthingis thatthejetsalsohavedifferentdirectionsfordifferentxpositions.Inotherwords,thedirectionofthe jetscontinuouslychangesfollowingthedirectionofthemagnetidfieldwhichisshearinginthecurrent sheetforformingfan-shapedjetsas shownin the3D visualizationin Figures2 and3. Thecartoonin Figure9showsbothordinaryreconnectionjetsandthefan-shapedjets.Onecanclearlyunderstandthe shapeandthestuctureofthejetsbythiscartoon. Theinteresting1Dvelocitydistributionsalongx(y = 0.0,z = 0.0),y (x = 0.0,z = 0.0)andz (x = 0.0,y = 0.0)directionsattime= 8.5and11.5areshowninFigure10.Notethatthevelocities in x, y and z directions stand for inflow speed, ordinary reconnection outflow speed and fan-shaped jet speed, respectively.The inflow speed (x direction)becomesa little smaller at the later stage (time =11.5),whichindicatesthatthemagneticreconnectionrate(asdiscussedlater)reachesthemaximum aroundthetime=8.5.Atthelatterstagetheordinaryreconnectionoutflowspeedisapproximatelythe 6 R.L.Jiang,K.Shibata,H.Isobe&C.Fang Fig.2 3Dvisualizationofthegaspressuredistributionsatdifferenttimes(typicalcase).The brightcolorsandhighopacitystandforthehighgaspressureandthelowvalueistransparent. Solidtubesaremagneticfieldlinesintheupperpanelsandarrowsarevelocityinthelower panels. Since the current sheet is very thin, a differentscale in x direction is used in these panels. ambientAlfve´nspeed(1.58).Moreover,themaximumvelocityofthefan-shapedjetsalongz direction isabout0.6.Thus,thedrivingforcesof thefan-shapedjetsandtheordinaryreconnectionoutfloware different. Inordertoinvestigatethedrivingforceofthefan-shapedjets,weexaminethetrajectoriesandthe forcesof Lagrangianfluidelements(testparticles).Thisis also a methodto showhowplasma moves and what kind of force drives plasma to move (Shibata&Uchida 1986; Katoetal. 2002). Figure 11 showsthetrajectoryofatypicalfluidelementinthefan-shapedjet.Thisfluidelementislocatedatthe edgeofthediffusionregionattheinitialtime andthemagneticfieldisperpendiculartothetrajectory (i.e.,perpendiculartothevelocitydirection).Afterthat,thiselementisejectedasajetatthetime=8.5 andeventuallymovesalongtheguidefieldline.Thecosinevalueoftheanglebetweenvelocity(v)and magneticfield(B)vectorsisshownintheleftpanelofFigure12.Attheinitialtime,thevelocityofthe particle is perpendicularto the magnetic field and almost parallel to that at the time = 14. The right Fan-Shapedjets 7 Fig.3 3Dvisualizationofthegaspressuredistributionsatdifferentviewpoints(typicalcase). Thebrightcolorsandhighopacitystandforthehighgaspressureandthelowvalueistrans- parent.Solidtubesaremagneticfieldlines. panelof Figure 12 shows the forcesin z directionon this element. It shows that there are two stages foracceleration.Inthefirststage(time=8−10),themagnetictensionforcedrivesthiselementsince the gas pressure gradient force balances the magnetic pressure gradient force in the inflow region as we mentioned before. After that it is dominated by the gas pressure gradient in the later stage (time = 10−11.5). The trajectory in the x− z plane is shown in Figure 13. Note that the y direction is differentatdifferenttimestofollowthetrajectory.Atthetime= 8.5,theelementisdrivenasainflow andatthetime=11.5asafan-shapedoutflow. Becauseofthecomplexityinmeasuringthereconnectionratein3Dsimulations,wegivetwokinds ofmeasurementforthereconnectionrateusingtheinflowspeedandelectricfieldasshowninFigure14. Theinflowspeedisthemaximumoneinxdirection(atthepositionsy =0andz =0)andtheelectric fieldistakenatthepointneartheedgeofthediffusionregion(E andE arealmostzero,sowecan x y notsee them in this figure).By comparingthese two, one can have a overallview abouthow fast the magneticfield dissipatesin oursimulations.Theupperpanelsin Figure14show thatthese two kinds ofreconnectionrate havea similar profile.Bothof themreachthemaximumalmostatthe same time 8 R.L.Jiang,K.Shibata,H.Isobe&C.Fang Fig.4 2Ddistributions(typicalcase)ofthegaspressure,magneticpressureandtotalpressure inx−yplaneatz =0.0.Thearrowsintheupperpanelsindicatethemagneticfieldandothers thevelocityfield. Fan-Shapedjets 9 Fig.5 2Ddistributions(typicalcase)ofgaspressure,magneticpressureandtotalpressurein y−zplaneatx=0.0.Thearrowsintheupperpanelsindicatethemagneticfieldandothers thevelocityfield. 10 R.L.Jiang,K.Shibata,H.Isobe&C.Fang Fig.6 2Ddistributions(typicalcase)ofdensityandtemperatureiny−z planeatx = 0.0. Thearrowsintheupperpanelsindicatethemagneticfieldandothersthevelocityfield. (around8.5)asweanalyzedin1Ddistributionsandthatisalsothereasonwhywechoosethetime=8.5 andalatertime= 11.5toshowthe3D,2Dand1Ddistributions.ThelowerleftpanelinFigure14is themassandenergyconservationswithconsideringthefluxcominginoroutofthecomputationalbox. EvenifCIP-MOCCTschemeisnon-conservativeforvariables,ournumericalerrorisstillsmallenough to get a realistic result. W , W and W in the lower right panel are the total thermal energy, total t k m kineticenergyandtotalmagneticenergy,respectively.Themagneticenergyreleaseratebecomesfaster andfasterafterthereconnectionratereachingthemaximumbecauseoftheformationandextensionof theslowmodeshock(Nittaetal.2001).Actually,mostofmagneticenergyisreleasedbytheslowmode shock(Petschek1964). TheupperpanelsofFigure15depictthegaspressureandcurrentdensitydistributionsinzdirection, whichisonlythecentralpartofFigure4.Otherpanelsshowonedimensionalplotalongthewhiteline indicated in the upper left panel at the time = 14.0. We can see a pair of slow mode shocks in the current density distribution and the 1D panels. The value jumps between two sides of the shock are shown by the shadow partsin the 1D distributions.Maybe the shocksare notsharp enough,which is due to the lack of grid points (about 65 grid points along the white line). However, we can compare