Global 3D Simulations of Disc Accretion onto the classical T Tauri Star BP Tauri
MMon. Not. R. Astron. Soc. , 000–000 (0000) Printed 29 October 2018 (MN L A TEX style file v2.2)
Global 3D Simulations of Disc Accretion onto the classical TTauri Star BP Tauri
M. Long, (cid:63) , M. M. Romanova, † , A. K. Kulkarni and J.-F. Donati Center for Theoretical Astrophysics, Department of Physics, University of Illinois at Urbana-Champaign, Urbana, IL 61801-3080 Department of Astronomy, Cornell University, Ithaca, NY 14853-6801, USA Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, Cambridge, MA 02138, USA LATT-UMR 5572, CNRS& Univ. P. Sabatier, 14 Av. E. Belin, F-31400 Toulouse, France
29 October 2018
ABSTRACT
Recent spectropolarimetric observations of the classical T Tauri star BP Tau and anal-ysis of its surface magnetic field have shown that the magnetic field can be approx-imated as a superposition of slightly tilted dipole and octupole moments with re-spective strengths of the polar magnetic fields of 1.2 kG and 1.6 kG (Donati et al.2008, hereafter D08). We adopt the measured properties of BP Tau and model thedisc accretion onto the star by performing global three-dimensional magnetohydrody-namic simulations. We observed in simulations that the disc is disrupted by the dipolecomponent and matter flows towards the star in two funnel streams which form twoaccretion spots below the dipole magnetic poles. The octupolar component becomesdynamically important very close to the star and it redirects the matter flow to higherlatitudes and changes the distribution and shape of the accretion spots. The spots aremeridionally elongated and are located at higher latitudes, compared with the puredipole case, where crescent-shaped, latitudinally elongated spots form at lower lati-tudes. The position and shape of the spots are in good agreement with observations.The disk-magnetosphere interaction leads to the inflation of the field lines and tothe formation of magnetic towers above and below the disk. The magnetic field of BPTau is close to the potential inside the magnetospheric surface, where magnetic stressdominates over the matter stress. However, it strongly deviates from the potential atlarger distances from the star.A series of simulation runs were performed at different accretion rates. In one ofthem, the disk is truncated at r ≈ (6 − R (cid:63) which is close to the corotation radius, R cor ≈ . R (cid:63) . However, the accretion rate, . × − M (cid:12) yr − , is lower than thatobtained from most of the observations. In a sample model with a higher accretionrate . × − M (cid:12) yr − , the disk is truncated at r ≈ . R (cid:63) , but such a state can notbe a typical state for the slowly rotating BP Tau if it is in the rotational equilibrium.However, torque acting on the star is also small: it is about an order of magnitudelower than that which is required for the rotational equilibrium. We suggest that astar could lose most of its angular momentum at earlier stages of its evolution. Key words: accretion, accretion discs - magnetic fields - MHD - stars: magnetic fields.
Accretion-powered Classical T Tauri stars (CTTSs) are younglow-mass stars which often show signs of a strong magneticfield (e.g., Basri, Marcy & Valenti 1992; Johns-Krull, Valenti &Koresko 1999) which is expected to have a complex structure (cid:63)
E-mail:[email protected] † [email protected] (e.g., Johns-Krull 2007). The Zeeman-Doppler imaging tech-nique has proven very successful in obtaining surface mag-netic maps for many stars, and the external magnetic fieldsof the stars have been reconstructed from these maps underthe potential approximation (Donati & Collier Cameron 1997;Donati et al. 1999; Jardine et al. 2002, 2006; Gregory et al.2010). The magnetic field plays a crucial role in disc accretionby disrupting the inner regions of the disc and channeling the c (cid:13) a r X i v : . [ a s t r o - ph . S R ] S e p M. Long et al. matter onto the star, and hence it is important to know themagnetic field configurations in magnetized stars.D08 recently observed the CTTS BP Tau with the ES-PaDOnS and NARVAL spectropolarimeters and reconstructedthe surface magnetic field from the observations. They haveshown that the magnetic field of BP Tau can be approximatedby a combination of dipole and octupole components of 1.2kG and 1.6 kG, which are slightly (but differently) tilted aboutthe rotational axis.D07 analyzed the distribution of the accretion spots onthe stellar surface and found spots at high latitudes, whichcover about 8 per cent of the stellar surface (D08).Further, D08 extrapolated the surface magnetic field tolarger distances using the potential approximation, (i.e., as-suming that there are no currents outside the star and hencethe external matter does not influence the initial configura-tion of the field) and estimated the distance at which thedisk should be disrupted by the magnetosphere so that thematter flowing towards the star in funnel streams producesthe high-latitude spots. They concluded that this distanceshould be quite large, r (cid:38) R (cid:63) . However, this problem re-quires more complete analysis based on the MHD approach,where external currents can be taken into account, and thematter flow around the magnetosphere can be investigatedself-consistently, taking into account interaction of the exter-nal plasma with the magnetic field.In this paper, we investigate this problem using globalthree-dimensional MHD simulations. We solve the 3D MHDequations numerically in our simulation model to investigatethe structure of the external magnetic field, accretion flowsand location of accretion spots.In our previous work, we have performed global 3Dsimulations of accretion onto stars with misaligned dipolefields (Romanova et al. 2003, 2004a; Kulkarni & Romanova2005), and aligned or misaligned dipole plus quadrupolefields (Long, Romanova & Lovelace 2007, 2008). Recently,we were able to extend our method and to build a numer-ical model for stars with an octupolar component. The gen-eral properties of the model have been described in detail inLong, Romanova & Lamb (2010) and the model was appliedto another CTTS, V2129 Oph, with a strong octupole field inRomanova et al. (2010) which has been compared with ob-servations of V2129 (Donati et al. 2007, Donati et al. 2010).In this paper, we apply our 3D MHD model of stars withcomplex magnetic fields (Long, Romanova & Lamb 2010) tothe CTTS BP Tau. We investigate disc accretion onto the starby adopting the measured surface magnetic fields (D08) andother suggested properties of this star. The surface magneticfield is modeled as a superposition of a 1.2 kG dipole and 1.6kG octupole field, tilted by ◦ and ◦ with respect to therotational axis and located at opposite phases (the phase dif-ference is ◦ ). We also take into account other parametersof BP Tau: its mass M (cid:63) = 0 . M (cid:12) (Siess, Dufour & Forestini2000), and radius R (cid:63) = 1 . R (cid:12) (Gullbring et al. 1998). Itsage is about 1.5 Myr (D08), and its rotation period is . days(Vrba et al. 1986), which corresponds to a corotation radiusof R cor ≈ . R (cid:63) . The mass accretion rate derived from dif-ferent observations varies between ˙ M (cid:39) . × − M (cid:12) yr − (e.g. Gullbring et al. 1998) and × − M (cid:12) yr − (Schmittet al. 2005). We performed a series of simulation runs at dif-ferent mass accretion rates in order to investigate the caseswhere the disk stops at different distances from the star. We also calculated the torque on the star and compared it withstar’s age.To understand the role of the octupole field in channelingthe accreting matter, we compared our dipole plus octupolemodel of BP Tau with a similar model but with only the dipolecomponent.Thus, we focus on: (1) 3D MHD modeling of accretionflows around CTTS BP Tau modeled with dipole plus octupolemoments; (2) comparisons of accretion properties observedin simulations with observations of BP Tau, such as the shapeand distribution of hot spots, mass accretion rates and more;(3) deviation of the simulated magnetic field from the poten-tial field.Section § §
3. We end in § The global 3D MHD model originally developed by Koldoba etal. (2002) and used for modeling stars with dipole fields (Ro-manova et al. 2003, 2004a) was modified to include higherorder components (Long, Romanova & Lovelace 2007, 2008;Long, Romanova & Lamb 2010; Romanova et al. 2010) tosimulate disc accretion onto stars with complex fields. TheMHD equations are solved in a reference frame co-rotatingwith the star. A viscous term is incorporated into the MHDequations (only in the disc) to control the rate of matter flowthrough the disk. We use the α − prescription for viscosity with α = 0 . .(i) Initial conditions.
The simulation domain consists of acold, dense disc and a hot, low-density corona, which areinitially in rotational hydrodynamical equilibrium. The initialangular velocity in the disc is close to Keplerian. The angularvelocity in the corona at any given cylindrical radius is set tobe equal to that of the disk at that radius.(ii)
Boundary conditions.
At the inner boundary (the sur-face of the star), most of the variables A are set to have freeboundary conditions, ∂A/∂r = 0 . The initial magnetic fieldon the surface of the star is taken to be a superposition ofmisaligned dipole and octupole fields. As the simulation pro-ceeds, we assume that the normal component of the field re-mains unchanged, i.e., the magnetic field is frozen into thesurface of the star. At the outer boundary, free conditions aretaken for all variables. In addition, matter is not permitted toflow into the region from the outer boundary.(iii) Simulation region and grid.
We use the “cubed sphere”grid introduced by Koldoba et al. (2002) (see also Fig. 1 inLong, Romanova & Lamb 2010). The resolution of the gridis × × to simulate the accretion onto BP Tau in asimulation domain of R (cid:63) . We use a very high resolutionnear the star in order to resolve the complex structure of theoctupolar component of the field.(iv) Magnetic field configuration.
In our code we can modelthe magnetic field of the star by a superposition of three mul-tipole moments µ i ( i = 1 , , for dipole, quadrupole andoctupole respectively) which are tilted relative to the z − axis(which is aligned with the rotational axis Ω ) at different an-gles Θ i , and have different angles φ i between the xz planeand Ω − µ i planes. For simplicity, φ is set to be 0, whichmeans that the dipole moment µ is in xz plane. The general c (cid:13) , 000–000 imulations of CTTS BP Tau Table 1.
The reference values for CTTS BP Tau. The dimensional val-ues can be obtained by multiplying the dimensionless values fromsimulations by these reference values. B and the subsequent valuesbelow depend on (cid:101) µ .Reference Units (cid:101) µ = 1 (cid:101) µ = 2 (cid:101) µ = 3 M (cid:63) ( M (cid:12) ) 0 . – – R (cid:63) ( R (cid:12) ) 1 . – – B (cid:63) (G) – – R (cm) . × – – v (cm s − ) . × – – P (days) . – – B (G) . . . ρ (g cm − ) . × − . × − . × − ˙ M ( M (cid:12) yr − ) . × − . × − . × − F (erg cm − s − ) . × . × . × N (g cm s − ) . × . × . × magnetic field configuration is discussed in greater detail inLong, Romanova & Lamb (2010). The details of the magneticfield configuration in our BP Tau model are discussed in § Reference units.
The simulations are performed in di-mensionless variables (cid:101) A = A/A where A are reference val-ues. We choose the stellar mass M (cid:63) , radius R (cid:63) and the surfacedipole field strength B (cid:63) to build a set of reference values.The reference values are: length scale: R = R (cid:63) / . ; ve-locity: v = ( GM (cid:63) /R ) / ; time-scale: P = 2 πR /v . Thereference magnetic moments for dipole and octupole compo-nents are µ , = B R and µ , = B R respectively, where B is the reference magnetic field. Hence, the dimensionlessmagnetic moments are: (cid:101) µ = µ /µ , , (cid:101) µ = µ /µ , , wherethe dipole and octupole moments µ and µ of the star arefixed. We take one of the above relationships; for example,the one for the dipole; to obtain B = µ , R = 0 . B (cid:63) (cid:102) µ (cid:18) R (cid:63) R (cid:19) = 25 . (cid:18) B (cid:63) . (cid:19)(cid:18) (cid:102) µ (cid:19) G . (1)Hence, at fixed B (cid:63) , the reference magnetic field depends onthe dimensionless parameter (cid:101) µ . The reference density ρ andthe mass accretion rate ˙ M also depend on this parameter: ρ = B /v = 2 . × − (cid:18) B (cid:63) . (cid:19) (cid:18) (cid:102) µ (cid:19) gcm , (2) ˙ M = ρ v R = 1 . × − (cid:18) B (cid:63) . (cid:19) (cid:18) (cid:102) µ (cid:19) M (cid:12) yr . (3)The dimensional accretion rate then is ˙ M = (cid:102) ˙ M ˙ M , where (cid:102) ˙ M is the dimensionless accretion rate. One can see that thedimensional accretion rate depends on (cid:102) ˙ M which we find fromsimulations, and the dimensionless parameter (cid:101) µ , which weuse to vary the accretion rate. We change the dimension-less octupolar moment (cid:101) µ in same proportion so as to keepthe ratio µ /µ fixed. To find the ratio between the dimen-sionless moments, we use approximate formulae for alignedmoments: µ = 0 . B (cid:63) R (cid:63) and µ = 0 . B (cid:63) R (cid:63) (see Long,Romanova & Lamb 2010) and obtain for BP Tau the ratio (cid:101) µ / (cid:101) µ = B (cid:63) ˜ R (cid:63) / B (cid:63) ≈ . (where ˜ R (cid:63) = 0 . ).Other reference values are: angular momentum flux (atorque) ˙ N = ρ v R ; energy flux ˙ E = ρ v R ; temper-ature T = R p /ρ , where R is the gas constant; and the effective blackbody temperature T eff , = ( ρ v /σ ) / , where σ is the Stefan-Boltzmann constant. Tab. 1 shows the refer-ence values for CTTS V2129 Oph. In the subsequent sections,we show dimensionless values (cid:101) A for most of the variables anddrop the tildes( ∼ ). However, we keep them in (cid:101) µ , (cid:101) µ , and (cid:102) ˙ M because these are important parameters of the model.(vi) The magnetospheric radius.
The truncation radius, r t ,where the disc is truncated by the magnetosphere, could beestimated as (e.g., Elsner & Lamb 1977): r t = k ( GM (cid:63) ) − / ˙ M − / µ / , (4)where ˙ M is the accretion rate and µ is the dipole magneticmoment; k is a coefficient of order unity. For example, Long,Romanova & Lovelace (2005) obtained k ≈ . in numer-ical modeling of disc-accreting stars. For stars with known M (cid:63) , R (cid:63) , and dipole magnetic moment µ , such as BP Tau,the accretion rate determines where the disc stops and howthe matter flows onto the star. We performed a number of simulation runs at different accre-tion rates and observed that the disk is truncated at differentradii. We choose as our main case the one in which the diskstops sufficiently far away, but at the same time the accretionrate is not very low. Below, we describe this case in detail. Wealso show the results at higher accretion rates for comparisons(see § D08 decomposed the observed surface magnetic field of BPTau into spherical harmonics and found that the field ismainly poloidal with only of the total magnetic energyin the toroidal field. The poloidal component can be approxi-mated by dipole ( l = 1 ) and octupole ( l = 3 ) moments with and of the magnetic energy respectively. Other mul-tipoles (up to l < ) have only of the total magneticenergy. D08 concluded that the magnetic field of BP Tau isdominated by a 1.2 kG dipole and 1.6 kG octupole tilted by ◦ and ◦ . The meridional angle between the Ω − µ and Ω − µ planes is approximately ◦ . In our model, we onlyconsider the poloidal component.We convert the above parameters into dimensionless val-ues using our reference units and solve 3D MHD equations(see, e.g., Koldoba et al. 2002), in dimensionless form. Oneof the important parameters of the model is (cid:101) µ which is usedto vary the truncation radius in the dimensionless model (andthe accretion rate in the dimensional model, see eq. 1). Wefind the second dimensionless parameter (cid:101) µ from the rela-tionship: (cid:101) µ / (cid:101) µ = B (cid:63) ˜ R (cid:63) / B (cid:63) ≈ . (where ˜ R (cid:63) = 0 . ).This ratio is fixed for fixed values of the dipole and octupolecomponents, B (cid:63) = 1 . kG and B (cid:63) = 1 . kG. From a num-ber of simulation runs at different (cid:101) µ , we choose (cid:101) µ = 3 , (cid:101) µ = 0 . to obtain a large enough magnetosphere as sug-gested by D08 and investigate this case in detail. Other pa-rameters of the magnetic field configuration of BP Tau are Θ = 20 ◦ , Θ = 10 ◦ and φ = 180 ◦ .Fig. 1 shows the components of the simulated magneticfield ( B r , B φ and B θ ) in the polar projection down to thelatitude − ◦ . One can see that the distribution of the radial c (cid:13)000
The truncation radius, r t ,where the disc is truncated by the magnetosphere, could beestimated as (e.g., Elsner & Lamb 1977): r t = k ( GM (cid:63) ) − / ˙ M − / µ / , (4)where ˙ M is the accretion rate and µ is the dipole magneticmoment; k is a coefficient of order unity. For example, Long,Romanova & Lovelace (2005) obtained k ≈ . in numer-ical modeling of disc-accreting stars. For stars with known M (cid:63) , R (cid:63) , and dipole magnetic moment µ , such as BP Tau,the accretion rate determines where the disc stops and howthe matter flows onto the star. We performed a number of simulation runs at different accre-tion rates and observed that the disk is truncated at differentradii. We choose as our main case the one in which the diskstops sufficiently far away, but at the same time the accretionrate is not very low. Below, we describe this case in detail. Wealso show the results at higher accretion rates for comparisons(see § D08 decomposed the observed surface magnetic field of BPTau into spherical harmonics and found that the field ismainly poloidal with only of the total magnetic energyin the toroidal field. The poloidal component can be approxi-mated by dipole ( l = 1 ) and octupole ( l = 3 ) moments with and of the magnetic energy respectively. Other mul-tipoles (up to l < ) have only of the total magneticenergy. D08 concluded that the magnetic field of BP Tau isdominated by a 1.2 kG dipole and 1.6 kG octupole tilted by ◦ and ◦ . The meridional angle between the Ω − µ and Ω − µ planes is approximately ◦ . In our model, we onlyconsider the poloidal component.We convert the above parameters into dimensionless val-ues using our reference units and solve 3D MHD equations(see, e.g., Koldoba et al. 2002), in dimensionless form. Oneof the important parameters of the model is (cid:101) µ which is usedto vary the truncation radius in the dimensionless model (andthe accretion rate in the dimensional model, see eq. 1). Wefind the second dimensionless parameter (cid:101) µ from the rela-tionship: (cid:101) µ / (cid:101) µ = B (cid:63) ˜ R (cid:63) / B (cid:63) ≈ . (where ˜ R (cid:63) = 0 . ).This ratio is fixed for fixed values of the dipole and octupolecomponents, B (cid:63) = 1 . kG and B (cid:63) = 1 . kG. From a num-ber of simulation runs at different (cid:101) µ , we choose (cid:101) µ = 3 , (cid:101) µ = 0 . to obtain a large enough magnetosphere as sug-gested by D08 and investigate this case in detail. Other pa-rameters of the magnetic field configuration of BP Tau are Θ = 20 ◦ , Θ = 10 ◦ and φ = 180 ◦ .Fig. 1 shows the components of the simulated magneticfield ( B r , B φ and B θ ) in the polar projection down to thelatitude − ◦ . One can see that the distribution of the radial c (cid:13)000 , 000–000 M. Long et al.
Figure 1. : Polar projections of the magnetic field components atthe stellar surface in the dipole plus octupole model of BP Tau inspherical coordinates: radial magnetic field, B r ; azimuthal magneticfield, B φ ; and meridional magnetic field, B θ . The outer boundary,the bold circle and the two inner dashed circles represent the latitudeof − ◦ , the equator, and the latitudes of ◦ and ◦ respectively.The red and blue regions represent positive and negative polarities ofthe magnetic field. Figure 2.
The surface magnetic field in the dipole plus octupolemodel of BP Tau ( (cid:101) µ = 3 , (cid:101) µ = 0 . , Θ = 20 ◦ , Θ = 10 ◦ , φ = 180 ◦ ) as seen from the equatorial plane (left panel), the northpole (middle panel) and the south pole (right panel). The colors rep-resent different polarities and strengths of the magnetic field. component (left panel) is very similar to that obtained by D08in two observational epochs of Dec06 and Feb06 (see Fig. 14,left panels of D08). In both cases, there is a strong positivepole at colatitudes of − ◦ , a part of the negative octupolarbelt at colatitudes of ◦ − ◦ and a part of the positive oc-tupolar belt at colatitudes of ◦ − ◦ . The distribution ofthe meridional component, B θ , (right panel of Fig. 1) is qual-itatively similar to that of D08, though in the D08 plot theinner positive ring (red color in the plot) is weaker, while thenegative ring (blue) is stronger compared with our model.The azimuthal component of the field is weak and shows abutterfly pattern of the negative and positive polarities in boththe modeled and reconstructed fields.Fig. 2 shows three-dimensional views of the magneticfield distribution at the surface of the star. It can be seenthat there are two antipodal polar regions of opposite polar-ity where the field is strongest. They approximately coincidewith the octupolar high-latitude poles and their centers arelocated close to the octupolar axis µ . Next to these regions,there are negative (blue) and positive (red) octupolar belts.Their shapes are more complex than those of the belts in pureoctupole cases, where the belts are parallel to the magneticequator (see Long, Romanova & Lamb 2010). This is becausethe dipole component strongly distorts the “background” oc-tupolar field.Although the octupole field is strongest at the surface ofthe star, it decreases more rapidly than the dipole field with Figure 4.
3D view of matter flow and the magnetic field distributionin the main case at t = 10 . The left panel shows the distribution ofthe magnetic field lines. The middle panel shows one of the densitylevels, ρ = 0 . × − . The right panel shows the density distri-bution in the disc plane. The colors along the field lines representdifferent polarities and strength of the magnetic field. The thick cyan,white and orange lines represent the rotational axis and the dipoleand octupole moments respectively. distance from the star. To investigate the role of the dipoleand octupole components in channeling the accretion flow,we find the radius at which the dipole and octupole fieldsare equal. For this, we assume that both magnetic momentsare aligned with the rotational axis and take the magneticfield in the magnetic equatorial planes: B = µ /r and B = 3 µ / r (see Eqn. 1 in Long, Romanova & Lamb 2010).Noting that the strengths of the field at the magnetic poles are B (cid:63) = 2 µ /R (cid:63) and B (cid:63) = 4 µ /R (cid:63) , and equating B to B ,we find this radius: r eql = (cid:18) B (cid:63)B (cid:63) (cid:19) / R (cid:63) . (5)Substituting B (cid:63) = 1 . kG and B (cid:63) = 1 . kG, we obtain r eql ≈ R (cid:63) . Note that both dipole and octupole moments are tiltedabout the rotational axis, and hence the above formula givesonly an approximate value for this radius. We suggest thatthis radius should be located slightly above the surface of thestar to explain the dominance of the octupolar field seen inFig. 1 and Fig. 2.Fig. 3 shows the initial magnetic field distribution nearthe star in our model. One can see the octupolar field com-ponent in the vicinity of the star. It also modifies the dipolefield up to distances of r ≈ . R (cid:63) (above the surface of thestar). The dipole field dominates in the rest of the simulationregion. We should note that the dipole and octupole fields areequal at r eql ≈ R (cid:63) , but the octupole component disturbs thedipole field up to larger distances. Here, we show results of 3D MHD simulations of matter flowonto the dipole plus octupole model of BP Tau with param-eters corresponding to the main case ( (cid:101) µ = 3 , (cid:101) µ = 0 . )at time t = 10 when the system is in a quasi-stationary c (cid:13) , 000–000 imulations of CTTS BP Tau Figure 3.
The magnetic field of BP Tau modeled as a superposition of dipole and octupole components at t = 0 . The field shows an octupolarstructure only very close to the star, while the dipole field dominates in the rest of the simulation region. The color of the field lines representsthe polarity and strength of the field. The thick cyan, white and orange lines represent the rotational axis and the dipole and octupole momentsrespectively. state. Fig. 4 (middle panel) shows that the disk is truncatedby the dipole component of the field, and matter flows to-wards the star in two ordered funnel streams. The right panelshows the density distribution in the equatorial plane. Theleft panel shows that the magnetosphere is disturbed by thedisk-magnetosphere interaction.Fig. 5 shows slices of the density distribution and pro-jected field lines. Panel (b) shows the slice in the µ − Ω (or xz ) plane. One can see more clearly that the dipole compo-nent of the field is responsible for disk truncation and matteris channeled in the µ − Ω plane. Panel (c) shows that thedisk is stopped by the magnetosphere in the yz − plane. It isstopped at r ≈ (6 − R (cid:63) . Panel (d) shows that matter flow iscomplex in the equatorial plane, and that matter comes closerto the star in the yz − plane.The gray line in Fig. 5 shows the magnetospheric radius, r m , the distance where the matter and magnetic stresses areequal: β = ( p + ρv ) / ( B / π ) = 1 . The magnetic stress domi-nates at r < r m (Ghosh & Lamb 1978, 1979; Long, Romanova& Lamb 2010). Fig. 5 shows that r m ≈ R (cid:63) which is smallerthan the actual truncation radius r t ≈ (6 − R (cid:63) . We shouldnote that there are different criteria for the truncation radiusand the above criterion usually gives the smallest truncationradius (see Bessolaz et al. 2008 for details).The external magnetic field evolves and deviates fromthe initial configuration due to disk-magnetosphere interac-tion. The left panels of Fig. 6 show the initial field distribu-tion at t = 0 which corresponds to the potential field pro-duced by the star’s internal currents. As the magnetic fieldevolves, the currents produced by the motion of plasma out-side the star become important and they change the poten-tial magnetic field significantly. The middle and right panelsshow the evolved field from different directions at t = 10 . Itcan be seen that on a large scale the magnetic field lines wraparound the rotational axis of the star and form a magnetictower (e.g. Lynden-Bell 1996; Kato, Hayashi & Matsumoto2004; Romanova et al. 2004b). This is a natural result of the Figure 5.
Density distribution and selected field lines in differentslices. Panels (a) and (b) show xz slices at t = 0 and t = 10 . Panels(c) and (d) show yz and xy slices at t = 10 . The gray line shows thedistance at which the matter stress equals the magnetic stress. magnetic coupling between the disk and star. The foot-pointsof the field lines on the star rotate faster than the foot-pointsthreading the disk, which leads to differential rotation alongthe lines, and their stretching and inflation into the corona.The rotational energy is converted into the magnetic energyassociated with these field lines. The evolved field structureon the large scale significantly differs from the potential fieldshown in the left panels.In the vicinity of the star the situation is different. Themiddle top panel of Fig. 6 shows that the magnetic field dis-tribution at t = 10 is almost identical to the initial field distri-bution. We conclude that the potential approximation is validonly within the parts of the magnetosphere where the mag- c (cid:13) , 000–000 M. Long et al.
Figure 6.
Comparison of the initial (potential) magnetic field distribution at t = 0 (left panels) with the field distribution at t = 10 (middle andright panels). The top panels show the magnetic field in the vicinity of the star, while the bottom panels show the field distribution in the wholesimulation region. The left and middle columns show the side view, while the right column shows the axial view of the field The color on thestar’s surface shows different polarities and strengths of the magnetic field. Figure 7.
The simulated accretion spots viewed from different di-rections at t = 10 : from the equatorial plane (left-hand panel), thenorth pole (middle panel), and the south pole (right-hand panel).The color contours show the density distribution of the matter. Thesolid lines represent the magnetic moments of the dipole ( µ ) andoctupole ( µ ). netic stress dominates, that is, at r (cid:46) R (cid:63) in our case. Thefield distribution strongly departs from potential at larger dis-tances from the star. D08 analyzed accretion spots from the observed brightnessenhancements in chromospheric lines, such as Ca II IRT andHe I which presumably form in (or near) the shock front closeto the stellar surface. The random flaring component, as wellas the time-variable veiling component which reflects vari-ation of the intrinsic accretion rate were removed from themodeling. It was also suggested that matter flows in the vicin-ity of the strong magnetic field and hence Zeeman-splittingfeatures were used for analysis.The corresponding maps of the local surface brightnessare shown in Fig. 9 of D08 for two epochs of observations Figure 8.
The color background shows the distribution of the mag-netic field on the surface of the star. The blue contours show the en-ergy distribution in the accretion spot. The solid black lines show thedipole, octupole and rotational axes. (Feb06 and Dec06). It can be seen that in both plots there isone bright spot located at colatitudes of ◦ − ◦ and cen-tered at θ c ≈ ◦ , while a lower-brightness area spreads upto a colatitude of ◦ . In both cases the spots are elongatedin the meridional direction. The Feb06 spot also has an an-tipodal spot of weaker brightness and of similar shape. Theaccretion filling-factor is shown in Fig. 13 of D08. It showsa spot located at colatitudes of ◦ − ◦ and centered at θ c ≈ ◦ − ◦ .In our simulations the spots represent a slice of den-sity/energy taken across the funnel stream at the surface ofthe star. Hence, the spots show the distribution of these val-ues across the stream (Romanova et al. 2004a). The physicsof disk-magnetosphere interaction is expected to be morecomplex (e.g., Koldoba et al. 2008; Cranmer 2008, 2009;Brickhouse et al. 2010). However, this more complex physics c (cid:13) , 000–000 imulations of CTTS BP Tau Figure 9.
The simulated accretion spots in a polar projection down tocolatitude ◦ from the north pole. Left panel: density distribution;right panel: energy flux distribution. The equator is shown as a boldline. The dashed lines represent the latitude ◦ and ◦ respectively. strongly depends on the properties of funnel streams. Hencewe call these spots “accretion spots” and show the density dis-tribution in spots and also the energy flux: F = ρ v · ˆ r [( v − v (cid:63) ) / γp/ ( γ − ρ ] distribution (right panel) at t = 10 ,where v (cid:63) is the velocity of the star.Fig. 7 shows the density distribution at the surface of thestar at t = 10 . It can be seen that there are two antipodalspots which are centered in the µ − Ω plane slightly be-low the dipole magnetic pole. The spots are elongated in themeridional direction (unlike the spots in a pure dipole config-uration, see also § ◦ and spreadsbetween latitudes of ◦ and ◦ . These spots have a closeresemblance to the spots observed by D08. In both cases, thespots are stretched in the meridional direction and are locatedbetween latitudes of ◦ and ◦ . Analysis of our spots at dif-ferent moments of time has shown that the positions of thespots vary only slightly with time (within 0.1 in phase). How-ever, the spots observed by D08 (their Fig. 9) are located atdifferent phases. Feb06 observations show the main spot lo-cated at a higher phase compared with our spot, with a phasedifference of 0.15-0.2. This is quite good agreement, consid-ering the fact the D08 method of the dipole moment phasereconstruction has an error in phase of 0.1-0.2 or larger. Thephase difference is higher in the Dec06 observations wherethe simulated and observed spots are almost in anti-phase.To investigate the role of the dipole and octupole phasesin the spot’s position, we considered a few “extreme” cases,where the phase difference between the dipole and octupolemoments is φ = 0 ◦ , φ = 90 ◦ and φ = − ◦ (we have aphase difference of φ = 180 ◦ in the main case). We observedthat the accretion spot is similar in shape in all these cases,but is located at the phase corresponding to the phase of thedipole (in the µ − Ω plane). We conclude that the phase of Figure 11.
Fraction of the star’s surface covered by spots of density ρ and higher at t = 10 . The solid and dashed lines represent the dipoleplus octupole model and dipole model respectively. the accretion spots on the surface of the star gives a strongindication of the dipole component’s position. To understand the role of the octupolar field component inour dipole plus octupole model of BP Tau, we investigate asimilar model but with zero octupolar field, that is B (cid:63) = 1 . kG, Θ = 20 ◦ and B (cid:63) = 0 .Fig. 10 shows xz slices of the accretion flow in the dipoleplus octupole model (left panel) and pure dipole model (rightpanel). One can see that the flow is wider in the dipole plusoctupole model, because near the star the octupole field influ-ences the matter flow and redirects it in such a way that somematter flows towards the direction of the octupolar magneticpole located at anti-phase with the dipole pole. The streamsare narrower in the case of the pure dipole field.The bottom panels show that the accretion spots in thedipole plus octupole model are meridionally elongated, whilein the dipole model they have their typical crescent shape andare elongated in the azimuthal direction (as usually seen inthe pure dipole cases, e.g. Romanova et al. 2003, 2004a).We conclude that in BP Tau, the relatively weak octupolarmagnetic field which dominates only very close to the surfaceof the star, strongly influences the shape and position of theaccretion spots. This influence is not as dramatic as in anotherT Tauri star, V2129 Oph, which has a much weaker dipolecomponent, and where the octupolar field splits the funnelstream into polar and octupolar belt flows (Romanova et al.2010). D08 calculated the fraction of the star covered with visiblechromospheric spots on BP Tau and concluded that they arespread over up to about of the stellar surface, and coverabout of the surface, assuming one-fourth of each surfacepixel is subject to accretion (see Figs. 9 & 13 in D08). In thispaper we are interested only in the area of the spots spreadingthat is in the area covered by the funnel stream, and hencewe take the larger value (8%). Figs. 9 and 13 of D08 showthat the actual size of the spots depends on the brightness.Sometimes we take the brightest parts of spots (4-5)% andneglect the dimmer parts for comparison. In observations, the c (cid:13) , 000–000 M. Long et al.
Figure 10.
Comparison of matter flow and accretion spots in the dipole plus octuple model of BP Tau (left panels) with a pure dipole case (rightpanels) at t = 10 . The top panels show the density distribution (color background) and the sample poloidal field lines in the xz − plane. Thebottom panels show the density distribution in the accretion spots as seen along the y − axis (panels a,c) and along the rotational axis (panelsb,d). Figure 12.
The density distribution of spots with a cutoff at different density levels ρ , and the corresponding spot coverage f for the sum of thevisible and antipodal spots . spot’s coverage is calculated for one (visible) spot versus thewhole area of the star, while in the simulations we take intoaccount both spots, and hence the area in simulations is twiceas large compared with observations.We observed from the simulations that the size of theaccretion spot depends on the density (or energy) cutoff. Wecalculate the fraction of the star covered with the spots as f ( ρ ) = A ( ρ ) / πR (cid:63) , where A ( ρ ) is the area covered with (all)spots of density ρ or higher. Fig. 11 shows the distribution of f ( ρ ) versus ρ for the dipole plus octupole model and dipolemodel of BP Tau. It can be seen that for almost all chosendensity levels, the accretion spots occupy a larger surface areain the dipole plus octupole model than in the dipole model.This result is in good agreement with that obtained from Fig.10. We choose several cutoff densities and show the spotsin Fig. 12. One can see that in both simulations and obser-vations, spots cover different areas depending on the bright-ness/energy flux levels ranging from up to in simula-tions (for two spots). In both observations (see Fig. 9 of D08) and simulations (see our Fig. 10), the brightest parts of spotsare located at ◦ < θ c < ◦ .It was predicted that the area covered by spots in the caseof complex magnetic fields should be smaller than in the puredipole case (e.g., Mohanty & Shu 2008; Gregory et al. 2008,see also Calvet & Gullbring 1998). This paper and our previ-ous work (Long, Romanova & Lovelace 2008) show that thearea f could show a complex trend for mixed dipole and mul-tipole configurations, depending on how the multipole fieldredirects the accretion flow. The mass accretion rate of BP Tau obtained from observationsis not uniquely determined, and different authors give differ-ent results. For example, Gullbring et al. (1998) estimated themass accretion rate to be . × − M (cid:12) yr − using the lumi-nosity of L (cid:39) ( GM (cid:63) ˙ M/R (cid:63) )(1 − R (cid:63) /R in ) , where R in is theinner radius of the disk. Other estimates of the mass accre-tion rate include . × − M (cid:12) yr − (Valenti & Johns-Krull c (cid:13) , 000–000 imulations of CTTS BP Tau × − M (cid:12) yr − (Schmitt et al. 2005). Calvet et al.(2004) compared the mass accretion rates of intermediate-mass, (1 . − M (cid:12) , and low-mass, (0 . − . M (cid:12) , T Tauristars and concluded that the average mass accretion rate forintermediate-mass T Tauri stars is about . × − M (cid:12) yr − ,while for low-mass CTTSs (like BP Tau) it is about 5 timessmaller, that is, . × − M (cid:12) yr − .We obtained from simulations the dimensionless accre-tion rate of (cid:102) ˙ M ≈ . . Using the reference value for ˙ M fromTab. 1, ˙ M = 1 . × − M (cid:12) yr − (calculated for (cid:101) µ = 3 ),we obtain the dimensional accretion rate: ˙ M = ˙ M (cid:102) ˙ M ≈ . × − M (cid:12) yr − . This value is higher than the accretionrate suggested by Schmitt et al. (2005) but lower than theother estimates. Given that BP Tau has the low mass of ∼ . M (cid:12) (D08) (or . M (cid:12) Gullbring et al. 1998), the massaccretion rate from our simulation models may still be in thereasonable range according to the Calvet et al. (2004) analy-sis and discussion.Here, we should note that the “main model” consideredabove is one of the suggested models, where the disk is trun-cated at large distances. Below, we discuss additional models,where the disk is truncated at smaller distances and the ac-cretion rate is higher.
We performed two additional simulation runs at higher ac-cretion rates (at smaller values of parameter (cid:101) µ ) keepingthe ratio between the dipole and octupole moments fixed, (cid:101) µ / (cid:101) µ = 0 . . The parameters of these new models are thefollowing: (cid:101) µ = 2 , (cid:101) µ = 0 . and (cid:101) µ = 1 , (cid:101) µ = 0 . .Fig. 13 (middle and right panels) shows the matter flowand accretion spots for these two cases, while the left panelshows results for the main model ( (cid:101) µ = 3 , (cid:101) µ = 0 . ) forcomparison. One can see that in the new models, the disk istruncated at smaller distances from the star: r t ≈ (5 − R (cid:63) (for (cid:101) µ = 2 ) and r t ≈ . R (cid:63) (for (cid:101) µ = 1 ). The bottom panelsof Fig. 13 show that when the disk comes closer to the star,the spots are still at high latitudes, though they become longerin the meridional direction. This is probably because the oc-tupolar component has a stronger influence on the spots: thespots move closer to the octupolar magnetic pole. Note thatno octupolar ring spots appear (like in V2129 Oph, see R10).This is probably because the disk truncation radius is still faraway from the area where octupole has a strong influence onthe flow, which is closer to the star. It is interesting to notethat in the pure dipole case, the spots move towards lowerlatitudes when the disk comes closer to the star. Here, we seethe opposite: parts of the spots move to higher latitudes, be-cause the octupolar component becomes more significant atsmaller truncation radii.The dimensionless mass accretion rate which we ob-tain from our simulations is (cid:102) ˙ M ≈ . (for (cid:101) µ = 2 ) and (cid:102) ˙ M ≈ . (for (cid:101) µ = 1 ). We take reference values ˙ M for dif-ferent (cid:101) µ from Tab. 1 and obtain corresponding dimensionalaccretion rates ˙ M (see Tab. 2) which are . × − M (cid:12) yr − and . × − M (cid:12) yr − respectively. We see that the accretionrate for (cid:101) µ = 1 is high and close to many of the observedvalues. However, the disk comes too close to the star and thetruncation radius is much smaller than the corotation radius( R cor ≈ . R (cid:63) ), which would mean that the star is not in Figure 13.
Simulations at different accretion rates at t = 8 . (a)main case, (cid:101) µ = 3 , (cid:101) µ = 0 . ; (b) (cid:101) µ = 2 , (cid:101) µ = 0 . ; (c), (cid:101) µ = 1 , (cid:101) µ = 0 . . The top panels show the density distribution(color background), the poloidal field lines (red lines) and the β = 1 line (gray) in the xz − plane. The bottom panels show the energy fluxdistribution in the accretion spots in a polar projection. the rotational equilibrium state, which is unlikely. The modelwith (cid:101) µ = 2 is similar to the main model: the disk is truncatedfar away and the accretion rate is twice as high, but still it islower than many values given by observations. This model isslightly better than the “main” model.In the case of a pure dipole, we obtain from simulations (cid:102) ˙ M ≈ . and ˙ M ≈ . × − M (cid:12) yr − . 2 summarizes theresults obtained in our simulation models and some observa-tional properties of BP Tau for comparison. We calculated torque on the surface of the star for the abovethree models. The torque associated with the magnetic field ˜ N f is about 10-30 times larger than the torque associatedwith matter and thus dominates (see also Romanova et al.2002, 2004a; Long, Romanova & Lovelace 2007, 2008). Thedimensionless torque obtained from simulations, (cid:101) N , varies intime around some average value. Tab. 2 shows these valuesand values corresponding to these variations. We observed insimulations that the torque is negative in most of the modelsbecause the disk is truncated at the distances comparable withthe corotation radius and therefore the star spins-down. In amodel with (cid:101) µ = 1 , the disk is truncated at the distance of ∼ . R cor and it spins the star up. To estimate the dimensionaltorque, N f , we take the largest absolute value of (cid:101) N f fromTab. 2 and take into account the reference values of N fromTab. 1. Tab. 2 shows the dimensional torque.We estimate the time-scale of spinning-down of BP Tau. c (cid:13) , 000–000 M. Long et al.
Table 2.
Comparison of observational and modeled results shown in D08 and results from our simulations.Models r t (cid:102) ˙ M ˙ M ( M (cid:12) yr − ) (cid:102) N f N f (g cm s − ) τ (yr)D08 & other measurements > R (cid:63) – − − − –Model ( (cid:101) µ = 3 , (cid:101) µ = 0 . ) (6 − R (cid:63) .
13 1 . × − − . ± . − . × . × (spin-down)Model ( (cid:101) µ = 2 , (cid:101) µ = 0 . ) (5 . − . R (cid:63) .
075 1 . × − − . ± . − . × . × (spin-down)Model ( (cid:101) µ = 1 , (cid:101) µ = 0 . ) . R (cid:63) .
075 7 . × − +0 . ± . + . × . × (spin-up)Model ( (cid:101) µ = 3 , (cid:101) µ = 0 ) (6 − R (cid:63) .
12 1 . × − − . ± . − . × . × (spin-down) The star’s angular velocity is (at period P = 7 . days) Ω = 2 π/P ≈ . × − s − , its angular momentum is J = kM (cid:63) R (cid:63) Ω = 2 . × k gcm / s , where k < (wetake k = 0 . for estimations). The spin-down time-scale, τ = J/N f , was estimated for different models and is shown inTab. 2. One can see that the spin-up/down time-scale is aboutan order magnitude larger than the age of BP Tau ( . × years, see references in D08). Hence, the torque obtained insimulations is not sufficient to regulate the spin of the starat the considered epoch. It is possible that a star lost most ofits angular momentum at the earlier stages of its evolution.Note that similar low-torque situation has been observed inmodeling of another CTTS, V2129 Oph. We performed global 3D MHD simulations of accretion ontoa model star with parameters close to those of the classical TTauri star BP Tau, and with the magnetic field approximatedwith a superposition of slightly tilted dipole and octupole mo-ments with polar magnetic fields of 1.2kG and 1.6kG whichare in anti-phase (D08). In this star the dipole field dominatesalmost up to the surface of the star and determines the ma-jority of the observational properties.We performed a number of simulation runs for differenttruncation radii of the disk, and chose one of them where thedisk is truncated at r ≈ (6 − R (cid:63) which is sufficiently closeto the corotation radius R cor ≈ . R (cid:63) and investigated thiscase in greater detail. We observed that the dipole componentof the field truncates the disk, and matter flows in two funnelstreams towards the dipolar magnetic poles. However, nearthe star the flow is slightly redirected by the octupolar com-ponent towards higher latitudes, and this affects the shapeand position of the accretion spots: the spots are stretchedin the meridional direction and are centered at higher lati-tudes compared with spots in the pure dipole case, which arelatitudinally-elongated (crescent-shaped) and are centered atlower latitudes.The spots are located near the µ − Ω plane, where boththe dipole and octupole moments are situated, and are in anti-phase. Experiments with different relative phases between thedipole and octupole moments have shown that the spots arealways located in the meridional plane of the dipole moment.The spot’s position slightly varies in phase with the ac-cretion rate, but the variation is small, about ◦ at the most.Note that in the pure dipole case, the spot may rotate aboutthe magnetic pole (and hence may strongly change the phase)if the dipole is only slightly tilted about the rotational axis, Θ (cid:46) ◦ (Romanova et al. 2003, 2004a; Bachetti et al. 2010).However, in BP Tau, the tilt is larger, Θ ≈ ◦ and an oc- tupolar component is present and is significant near the star.Both factors lead to the restriction of such a motion.The meridional position and shape of spots observedin simulations is similar to those observed by Donati et al.(2008) (see their Fig. 9 for brightness distribution in spots re-constructed for the Dec06 or Feb06 epochs). However, theyare at different phases compared with our spots. This may beconnected with the relatively low accuracy of the phase recon-struction of the dipole component from the surface magneticfield, or some other reason.The accretion rates obtained in different models are inthe range of ˙ M ≈ (1 . − . × − M (cid:12) yr − . It is lower thanmost of the ˙ M values derived from observations which varyin the range of × − and . × − solar mass per yearand depend on the approaches used for the derivation. Thesmallest accretion rate obtained in simulations correspondsto the case where the disk is truncated at r t ∼ R cor , and astar slightly spins down, while the largest to the case where r t ∼ . R cor and a star spins up. We should note that atlarger accretion rates, say at . × − solar mass per year,the disk will come very close to the star, and the torque wouldbe much stronger. Such a state seems unlikely. Hence, if theaccretion rate is very high, then we should suggest that thedipole component should be larger than that derived by D08.The torque obtained in simulations is small and thetime-scale of the spinning-up/down is an order of magni-tude smaller compared with the age of BP Tau. This torqueis not sufficient to support a star in the rotational equilib-rium state, where a star spins up or spins down depend-ing on the accretion rate but has a zero torque on average(e.g., Ghosh & Lamb 1979, K¨onigl 1991, Cameron & Campbell1993, Long, Romanova & Lovelace 2005). We assume thatthis small torque matches a small accretion rate obtained insimulations because the torque generally correlates with theaccretion rate (e.g., Romanova et al. 2002).We suggest that a star may lose the majority of its an-gular momentum at earlier stages of its evolution due to the“propeller” effect (e.g., Romanova et al. 2005; Ustyugova etal. 2006), stellar winds (Matt & Pudritz 2005, 2008), or bysome other mechanism.Earlier, we performed global 3D MHD simulations of ac-cretion onto V2129 Oph (Romanova et al. 2010) which haveshown that in the case of a strong octupolar component, partsof the octupolar belt spots can be visible and can dominate atsufficiently high accretion rates. This is not the case in BP Tau,where the octupolar component dominates only in the closevicinity of the star.Disk-magnetosphere interaction leads to inflation of theexternal field lines and formation of a magnetic tower. Sim-ulations show that the potential approximation used in ex-trapolation of the magnetic field from the surface of the star c (cid:13) , 000–000 imulations of CTTS BP Tau to larger distances (e.g. Jardine et al. 2006, D08, Gregoryet al. 2010) is valid only inside the magnetospheric (Alfv´en)surface, where the magnetic stress dominates. At larger dis-tances, the magnetic field distribution strongly departs frompotential. ACKNOWLEDGMENTS
Resources supporting this work were provided by the NASAHigh-End Computing (HEC) Program through the NASA Ad-vanced Supercomputing (NAS) Division at Ames ResearchCenter and the NASA Center for Computational Sciences(NCCS) at Goddard Space Flight Center. The authors thankA.V. Koldoba and G.V. Ustyugova for the earlier develop-ment of the codes. The research was supported by NSF grantAST0709015. The research of MMR was supported by NASAgrant NNX08AH25G and NSF grant AST-0807129.
REFERENCES
Bachetti, M., Romanova, M. M., Kulkarni, A. K., Burderi, L.,di Salvo, T., 2010, MNRAS, 403, 1193Basri G., Marcy G.W., Valenti J.A., 1992, ApJ, 390, 622Bessolaz N., Zanni C., Ferreira J., Keppens R., Bouvier J.,2008, A&A, 478, 155Brickhouse N.S., Cranmer S.R., Dupree A.K., Luna G.J.M.,Wolk S., 2010, ApJ, 710, 1835Calvet N., Muzerolle J., Brice˜no C., Hern´ondez J., HartmannL., Saucedo J.S., Gordon K. D., 2004, AJ, 128, 1294Calvet N., Gullbring E., 1998, ApJ, 509, 802Cameron, A., Campbell, C, 1993, A&A, 274, 309Cranmer S. R., 2008, ApJ, 689, 316Cranmer S. R., 2009, ApJ, 706, 824Donati J.-F., Collier Cameron A. C., 1997, MNRAS, 291,1Donati J.-F., Collier Cameron A., Hussain G.A.J., Semel M.,1999, MNRAS, 302, 437Donati J.-F., Jardine M. M., Gregory S. G. et al., 2007, MN-RAS, 380, 1297Donati J.-F., Jardine M. M., Gregory S. G. et al., 2008, MN-RAS, 386, 1234Donati J.-F., Bouvier J., Walter F.M. et al., 2010, MNRAS, inpreparationElsner R. F., Lamb F. K., 1977, ApJ, 215, 897Ghosh P., Lamb F. K., 1978, ApJ, 223, L83Ghosh P., Lamb F. K., 1979, ApJ, 232, 259Gregory S. G., Matt S. P., Donati J.-F., Jardine M., 2008, MN-RAS, 389, 1839Gregory S. G., Jardine M., Gray, C.G., Donati J.-F., 2010, Rep.Prog. Phys., in pressGullbring E., Hartmann L., Briceno C., Calvet N., 1998, ApJ,492, 323Jardine M. M., Wood K., Collier Cameron A., Donati J.-F.,Mackay D.H., 2002, MNRAS, 336, 1364Jardine M. M., Collier Cameron A., Donati J.-F., GregoryS.G., Wood K., 2006, 367, 917Jardine M. M., Gregory S.G., Donati J. -F., 2008, MNRAS,386, 688Johns-Krull C., Valenti J.A., Koresko C., 1999, ApJ, 516, 900Johns-Krull C. M., 2007, ApJ, 664, 975Kato Y., Hayashi M.R., Matsumoto R., 2004, ApJ, 600, 338 Koldoba A. V., Romanova M. M., Ustyugova G. V., LovelaceR. V. E. 2002, ApJ, 576, L53Koldoba A. V., Ustyugova G. V., Romanova M. M., LovelaceR. V. E. 2008, MNRAS, 388, 357K¨onigl, A. 1991, ApJ, 370, L39Kulkarni A. K., Romanova M. M. 2005, ApJ, 633, 349Lamzin S. A., 1995, A&A, 295, L20Long M., Romanova M. M., Lovelace R. V. E., 2005, ApJ, 634,1214———. 2007, MNRAS, 374, 436———. 2008, MNRAS, 386, 1274Long M., Romanova M.M., Lamb F.K., 2010, MNRAS, inpress, (arxiv:0911:5455)Lynden-Bell D. 1996, MNRAS, 279, 389Matt, S., Pudritz, R.E. 2005, ApJ, 632, L135Matt, S., Pudritz, R.E. 2008, ApJ, 678, 1109Mohanty S., Shu F.H. 2008, ApJ, 687, 1323Romanova M. M., Kulkarni A. K., Lovelace R. V. E., 2008,ApJ, 673, L171Romanova M.M., Long M., Lamb F.K., Kulkarni A.K., DonatiJ.-F. 2010, MNRAS, in press (arXiv:0912.1681)Romanova M. M., Lovelace R. V. E., 2006, ApJ, 645, L73Romanova M. M., Ustyugova G. V., Koldoba A. V., LovelaceR. V. E., 2002, ApJ, 578, 420———. 2004a, ApJ, 610, 920———. 2004b, ApJ, 616, L151———. 2009, MNRAS, 399, 1802Romanova M. M., Ustyugova G. V., Koldoba A. V., Wick J. V.,Lovelace R. V. E., 2003, ApJ, 595, 1009Schmitt J. H. M. M., Robrade J., Ness J.-U., Favata F., StelzerB., 2005, A&A, 432, L35Siess L., Dufour E., Forestini M., 2000, A&A, 358, 593Ustyugova, G. V., Koldoba, A. V., Romanova, M. M., Lovelace,R. V. E., 2006, ApJ, 646, 304Vrba F. J., Rydgren A. E., Chugainov P. F., Shakovskaia N. I.,Zak D.S., 1986, ApJ, 306, 199Valenti J. A., Johns-Krull C.M., 2004, ApSS, 292, 619 c (cid:13)000