Table Of Content

Astronomy & Astrophysics manuscript no. bhlflow February 2, 2008 (DOI: will be inserted by hand later) High Mach-number Bondi–Hoyle–Lyttleton flow around a small accretor 1 Richard G. Edgar Stockholms Observatorium, AlbanovaUniversitetscentrum,SE-106 91, Stockholm, Sweden 5 Received / Accepted 0 0 2 Abstract. In this paper, we discuss a two-dimensional numerical study of isothermal high Mach number Bondi– Hoyle–Lyttleton flow around a small accretor. The flow is found to beunstable at high Mach numbers, with the n instability appearing even for a larger accretor. The instability appears to be the unstable radial mode of the a accretion column predicted by earlier analytic work. J 5 Key words.Accretion, accretion disks – Hydrodynamics– Instabilities 1 v 0 1. Introduction studied, due to the largecomputationalexpenseinvolved. 7 Edgar & Clarke (2004) observed a flow instability while 0 The Bondi–Hoyle–Lyttleton (BHL) accretion problem studing the effect of radiative feedback on Bondi–Hoyle– 1 concerns the accretion rate of a point mass moving at Lyttleton flow, and this paper looks more closely at this 0 constantvelocitythroughauniformgascloud.Thisprob- 5 instability.Thecomputationalmethodisdiscussedinsec- lem was first considered by Hoyle & Lyttleton (1939). 0 tion2andtheresultsinsection3.Finally,theconclusions ConsideringNewtoniantrajectoriesinapressurelessfluid, / are presented in section 4. h Hoyle & Lyttleton concluded that material with an im- p pact parameter (ζ) less than - o 2. Method 2GM r ζ <ζ = (1) st HL v∞2 The Zeus-2d code of Stone & Norman (1992) was used a to perform the calculations presented here. The author : would accrete onto the point mass, M (v is the velocity v ∞ parallelised the code using OpenMP, and found that it of the gas relative to the point mass at infinity). If the i X scaled well to up to eight processors. Extra routines to density of the gas is ρ (at infinity) then the accretion ∞ monitor accretion rates and drag were added to the main r rate onto the point mass will be a code, and calculations were performed in spherical polar 4πG2M2ρ (r,θ) co-ordinates. M˙ =πζ2 v ρ = ∞ (2) HL HL ∞ ∞ v3 A linear grid of 100 cells was used in the θ direction, ∞ while the radialgrid was logarithmically spaced. The size For a full discussion of the problem, see Edgar (2004). of the innermost radial cells was chosen to make the in- Numericalstudies haveshown(see Edgar(2004),andref- ner portions of the grid as ‘square’ as possible (to avoid erences therein) that, although the accretion rate pre- possible truncation errors). For the inner boundary con- dicted by equation 2 is fairly accurate, the flow geome- dition, three alternatives were tried: Zeus’ own ‘outflow’ try is somewhat more complicated, due to the effects of boundary, a constant (small) density, and the density set pressure. toasmallfractionoftheinnermost‘active’cell.Theouter Bondi–Hoyle–Lyttleton accretion is of particular in- boundary lay at 10ζ (as defined by equation 1) with a HL terest in star formation studies, since protostars in clus- uniforminflowimposedontheupstreamside,andoutflow ters are thought to undergo competitive accretion, at conditions on the downstream side. rates similar to that given by equation 2 (Bonnell et al. The gridwas initialised to a uniform flow,and the ac- 2001). However, under such conditions, v is low, and ∞ cretor introduced over a number of timesteps at the start so ζ is very large compared to the radius of the ac- HL of the run. The simulations were run for as long as possi- cretor (R ). This region of parameter space is poorly acc ble, to ensure that any ‘switch-on’ effects died away. The Send offprint requests to: Richard Edgar, e-mail: goal was a crossing time of the whole grid, but this was [email protected] not possible for the runs with the small accretor. 2 R. G. Edgar: High M BHL Flow around a small accretor Fig.1.AccretionratesasafunctionoftimeforlowMach number flow around a large accretor Fig.2.Detailoffigure1toshowdifferencesduetobound- Inthesimulations,dimensionlessunits wereused.The ary conditions incoming material had v = 1, and the sound speed was ∞ set according to the required Mach number (M). The point mass was set to M = 1 (as was the gravitational constant, G). Using equation 1 we see that ζ = 2 for HL these choices. A value of ρ =0.01 was chosen (arbitrar- ∞ ily,sincetheflowwasisothermalandnotself-gravitating), giving M˙ =0.126 in these units. HL 3. Results Three sets of simulations were performed. For eachof the sets,onerunforeachoftheinnerboundaryconditionswas made. First, a low Mach number (M = 3) flow around a large(R =0.005ζ )accretorwasusedtotestthecode. acc HL The second set used the same size of accretor, but with the Mach number increased to M=20. For the final set, the accretor size was reduced to 0.0005ζHL. Fig.3. Gravitational drag force in units of M˙HLv∞ for We will now discuss the resultsof eachof these sets in low Mach number flow around a large accretor turn. thissetofruns,wecanconcludethatthecodeisbehaving 3.1. Low Mach Number, Large Accretor well. The accretion rate as a function of time is presented in figure1.Afteraninitialturn-onperiod,theaccretionrate 3.2. High Mach Number, Large Accretor settles to a fairly steady value, slightly larger than that predicted by equation 2. Changing the inner boundary The accretionrates for the three runs of high Mach num- condition does not affect the results significantly, as is ber flow arounda large accretorare presented in figure 4. shown in figure 2. The different boundary conditions are The behaviour is very different to that shown in figure 1, givingslightlydifferentcurves,butareveryclosetogether. with a large oscillating accretion rate. We can see that The gravitational drag force on the accretor for this changingtheboundaryconditiondoeschangetheflowhis- setofrunsisplottedinfigure3.Thisisplottedinunitsof tory, but the change does not appear to be significant - M˙ v whichis the naturalscalingfor the problem.The the instability is similar in all cases. HL ∞ drag force rises smoothly to around 18M˙HLv∞. Figure 5 shows the corresponding gravitational drag Inspection of the density field shows a wake with a for this set of runs. As might be expected from the accre- fairlybroadopeningangle,butnobowshock-asistypical tion rate curves, it is quite noisy. Again, there are slight for isothermalBHL flow.None of these results are partic- differencesduetotheinnerboundarycondition,butthese ularly surprising, when compared to the work of Ruffert do not seem to be significant compared to the general in- (1996)(model‘GS’isclosesttothatpresentedhere).From stability of the flow. R. G. Edgar: High M BHL Flow around a small accretor 3 Fig.4.AccretionratesasafunctionoftimeforhighMach Fig.6. Sample radial velocity in the wake for high Mach number flow around a large accretor number flow around a large accretor trinsic to the flow along the line. Settling this dispute is beyond the scope of the present work. Looking at a single accretion rate history curve from figure 4, around two to three oscillations occur per simu- lationtime unit.Thisisinaccordancewiththeprediction of Cowie, who pointed out that the instability must grow faster than the characteristicaccretiontime (ζ /v =2 HL ∞ in this paper), or it will be swallowed by the accretor be- fore it can grow. Examining the velocity structure of the wake (figure 6 for a sample curve), oscillations in vr are observed superimposed on the freefall curve, as predicted byCowie.Theoscillationwavelengthdoesnotvarymuch, as required by the WKB approximation. The stagnation Fig.5.Gravitationaldrag(inunitsofM˙ v )asafunc- point is close to r =3.5,whichis unsurprising,giventhat HL ∞ the accretion rates are slightly higher than M˙ (which tion of time for high Mach number flow around a large HL would place the stagnation point at r = 2 – cf Edgar accretor (2004)). Unfortunately, computation of the power spectrum of Thewake’sopeningangleismuchlessthanforthelow the accretionrates (figure 7) is rather unenlightening due Mach number case, as one would expect. The oscillating to noise. The spectrum peaks at a frequency of around 2 accretionrateappearstobethe‘overstability’ofradialos- to3(itistoonoisytobemoreaccurate),andhasapower cillationsdescribedbyCowie(1977)andSoker(1990).The lawfalloffathigherfrequencies.Theremayalsobe aside analysiswas performed for a pressure-free(ballistic) flow, peak at a frequency of around 0.5, but this is only cov- andthehighMachnumberisothermalnatureofthesecal- ered by a few points, and hence may just be noise. The culationsisacloseapproximationofthiscondition.Cowie very intense bursts of accretion seen in figure 4 seem to extendedtheanalysisofBondi & Hoyle(1944)(seeEdgar beassociatedwithlarger‘blobs’ofmatterbeingaccreted. (2004) for a slightly more detailedexposition) to time de- However, these are relatively rare, with only a couple ap- pendentflow.Perturbingthesteadystatesolutionandap- pearing over the entire duration of the simulation (and plyingtheWKBapproximation(thatis,thewavelengthof one of those was probably a ‘switch-on’ transient). the perturbation is assumed to be short compared to the scale on which it varies), Cowie found that the accretion 3.3. High Mach Number, Small Accretor column should be unstable. Soker showed that viscosity couldsuppressthisinstability(whichhasnothappenedin Finally, we come to the simulations of high Mach num- thepresentwork).However,theydisagreedonthephysical ber flow around a small (0.0005ζ ) accretor. These sim- HL natureoftheinstability:Cowieattributedittotheincom- ulations were painfully expensive computationally. Each ingmaterial.Sokerassertedthattheequationshowingthe simulation in the first two sets required only about six possibilty for instability (the dispersion relation from the days of wall time (on two processors). These final three WKB approximation) did not contain terms relating to each occupied two processors for eight months, but still the inflowing material, and hence the instability was in- did not achieve the same elapsed simulation time. As the 4 R. G. Edgar: High M BHL Flow around a small accretor Fig.7. Sample accretion rate power spectrum for high Fig.9.Gravitationaldrag(inunitsofM˙ v )asafunc- HL ∞ Mach number flow around a large accretor tion of time for high Mach number flow around a small accretor Table 1. Mean Accretion Rates in units of M˙ for high HL Machnumberflowfordifferentinnerboundaryconditions Accretor Size Zeus Const. ρ Const. ρ Frac. Large 1.47 1.50 1.48 Small 1.46 1.46 1.47 ble 1. Shrinking the accretor does not seem to have made anysignificantdifferencetotheresults(beyondincreasing the requiredcomputing resourcesby well overan orderof magnitude). 4. Conclusions Fig.8.AccretionratesasafunctionoftimeforhighMach number flow around a small accretor Based on the results presented above, we see that high MachnumberBondi–Hoyle–Lyttletonflowaroundasmall accretor is unstable. Although the average accretion gridisextendedtowardsthe pointmass,thegravitational rate is generally similar to that originally predicted by accelerationgetsstronger.Therearemorecells,theveloc- Hoyle & Lyttleton (1939), shortincreases to rates several ities are higher, and the cell size is smaller. These three timesthisarepossible.Theinstabilityseemstoberelated effects combine to cause the consumption of prodigious to the high Mach number of the flow - at lower Mach computing resources. numbers,pressuremusthavea stabilisingeffect. So far as The accretion rate onto the point mass as a function these simulations can determine, the instability is the ra- of time is plotted in figure 8. The behaviour is generally dialoscillationoftheaccretioncolumndescribedbyCowie similar to figure 4, with large fluctuations once the flow (1977) and Soker (1990). becomes established. There are differences in detail be- Observationally, the instability may be visible for tween the three boundary conditions, but the qualitative deeply embedded young protostars.For a 1M⊙ protostar behaviourisobviouslythe same.Figure9showsthe grav- orbiting in a protocluster, the time unit in the present itational drag for this set of runs. Note that the gravita- paper would correspond to approximately one year. The tional force reaches larger values - as one might expect fluctuations in the accretion rate would correspond to a given the smaller size of the point mass (cf the Coulomb brightening of a factor of a few over a period of days to logarithm which appears in the calculation of dynamical weeks,anddimmingoverasimilartimescale.Thisisquite friction.Seee.g.Binney & Tremaine(1987)).Iftheflowis different from the FU Orionis phenomenon which is at- dense enough to make M˙HL interestingly large, then this tributedtoaccretiondiscactivity(cfBell et al.(2000)and dragforcewillalterv significantly.However,thegeneral referencestherein),wheretheoutburstsriseonatimescale ∞ patternissimilartofigure5,showingfluctuationscompa- of months and last for years. rable to those of the accretion rate. In the author’s opinion, there is little to be gained The mean accretion rates for the two sets of simula- by the further study of the pure Bondi–Hoyle–Lyttleton tions with high Mach number flow are presented in ta- problem. Shrinking the accretor size revealed no new be- R. G. Edgar: High M BHL Flow around a small accretor 5 haviour, at a tremendous cost in computer time. Instead, theBondi–Hoyle–Lyttletonmodelshouldberetainedasa useful ‘benchmark’ and work focused on improving simulations of more realistic scenarios. References Bell, K. R., Cassen, P. M., Wasson, J. T., & Woolum, D. S. 2000, Protostars and Planets IV, 897 Binney, J. & Tremaine, S. 1987, Galactic dynamics (Princeton, NJ, Princeton University Press, 1987, 747 p.) Bondi, H. & Hoyle, F. 1944, MNRAS, 104, 273 Bonnell,I.A.,Bate,M.R.,Clarke,C.J.,&Pringle,J.E. 2001, MNRAS, 323, 785 Cowie, L. L. 1977, MNRAS, 180, 491 Edgar, R. 2004,New Astronomy Review, 48, 843 Edgar, R. & Clarke, C. 2004, MNRAS, 349, 678 Hoyle,F. & Lyttleton,R. A. 1939,Proc.Cam.Phil. Soc., 35, 405 Ruffert, M. 1996, A&A, 311, 817 Soker, N. 1990, ApJ, 358, 545 Stone, J. M. & Norman, M. L. 1992, ApJS, 80, 753 Acknowledgements. The author is supported by EU-RTN HPRN-CT-2002-00308,“PLANETS”andthecalculationspre- sentedherewereperformedonthemachinesofHPC2N,Ume˚a and NSC, Linköping. Comments from an anonymous referee helped improve thediscussion of theinstability