A Formation Mechanism of Collapsar Black Hole -- early evolution phase
aa r X i v : . [ a s t r o - ph ] J un A Formation Mechanism of Collapsar Black Hole – early evolutionphase –
Yuichiro
Sekiguchi and Masaru
Shibata
Graduate School of Arts and Sciences, University of Tokyo, Tokyo 153-8902, Japan
The latest studies of massive star evolution indicate that an initially rapidly rotatingstar with sufficiently low metallicity can produce a rapidly rotating, massive stellar core thatcould be a progenitor of long-soft gamma-ray bursts (LGRBs). Motivated by these studies,we follow the collapse of a rapidly rotating massive stellar core to a ’collapsar’ black hole (BH)surrounded by a massive, hot accretion disk performing fully general relativistic simulations.We focus on the general relativistic dynamics of the collapse, and the relevant microphysics istreated in a qualitative manner. The simulations are performed until the system consisting ofthe BH and the disk has relaxed to a quasi-stationary state. A novel mechanism found in thisstudy is that strong shock waves are formed at the inner part of the disk after the formationof the BH. These shock waves propagate mainly along the rotation axis, heating the diskand sweeping materials around the rotational axis, and thereby forming a low density region.The temperature of the disk is high enough for copious neutrino emission. All these featuresindicate that the direct formation of a rapidly rotating BH is a promising source of LGRBseven in the absence of strong magnetic fields. §
1. Introduction
The observed association of afterglows of long gamma-ray bursts (LGRBs) withstar forming regions in galaxies, such as the connections between GRB980425and the Type Ib/c supernova SN1998bw, between GRB030329 and the Type Icsupernova SN2003dh, and between GRB060218 and the Type Ic supernovaSN2006aj, (see also Ref. 11) for other suggested associations between LGRBsand core-collapse supernovae), suggested that at least some LGRBs are associatedwith supernovae and the core collapse of massive stars. Supported by these evidences,the collapsar scenario, in which the central engine of LGRBs is composed of a rotatingblack hole (BH) surrounded by a massive, hot accretion disk, is currentlythe most favored model for LGRBs. It requires the progenitors to be rotating rapidlyenough that an accretion disk can be formed around the BH. Also, the relativisticjets that are to produce the LGRB have to reach the stellar surface. This constrainsthe size of the progenitor; the progenitor should not have an extended hydrogenenvelope.
Thus, the progenitors of LGRBs are now believed to be rapidlyrotating massive Wolf-Rayet (WR) stars.However, WR stars are known to be accompanied by strong stellar winds drivenby radiation pressure which lead to a rapid spin-down of the cores of the WR stars.
Some authors have proposed binary-evolution models for producing rapidlyrotating progenitors of LGRBs. On the other hand, Yoon and Langer and Woosleyand Heger have recently discovered that even a single star can fulfill therequirements of the collapsar scenario if it is initially rapidly rotating ( &
50% of theKeplerian velocity at the equatorial surface) and of low metallicity (
Z/Z ⊙ . . typeset using PTP
TEX.cls h Ver.0.9 i Y. Sekiguchi and M. Shibata
The mass loss is suppressed by the low metallicity.
The rapid rotation results in ashort mixing timescale, which leads to a chemically homogeneous state throughoutthe hydrogen burning phase.
In this case, a single star can become a rapidlyrotating WR star without losing the hydrogen envelope through the stellar wind,avoiding the red giant phase that otherwise would cause a significant decrease of thecore angular momentum due to magnetic torques.
In addition to the above single-star evolutionary scenario to the LGRB pro-genitor, recent growing empirical evidence indicates that LGRBs may prefer a lowmetallicity environment.
At least some LGRBs are likely to be formedfrom a rapidly rotating WR star which is born in a low metallicity environment andexperiences chemically homogeneous evolution.Motivated by these latest studies, we perform axisymmetric simulations of rapidlyrotating stellar core collapse to a BH in full general relativity. We explore the out-come of the collapse of very massive, rapidly rotating WR stars in the context of thecollapsar scenario. As a first step toward more realistic simulations, microphysicalprocesses are treated only in a qualitative manner, focusing on the general relativis-tic hydrodynamics of the collapse. Throughout this paper, c and G denote the speedof light and the gravitational constant, respectively. §
2. Setting
To model the core of WR stars just before the collapse, we adopt the polytropicequation of state (EOS) P = Kρ Γ , with Γ = 4 /
3, taking account of degenerateelectron and radiation pressures. Here K is a constant. Rigid rotation is adopted asthe rotation profile for simplicity. The maximum angular velocity (which is equal tothe Kepler angular velocity of the equatorial surface) is assumed, since the progenitor(a WR star) of LGRBs is likely to be rapidly rotating. The central density andtemperature of a pre-collapse iron core of massive metal-poor stars are ρ c & g/cm and T c . K (according to Umeda and Nomoto). We set the centraldensity of the initial conditions to ρ c ≈ × g/cm . Note that the chemicalpotential of the electrons at such a density is larger than the temperature T ≈ K, and hence the electrons are degenerate even at such a high temperature.According to the latest models of stellar evolution, rapid rotation resultsin a chemically homogeneous state and leads to a large CO core which could producean iron core with mass ≫ M ⊙ . The CO core mass M CO can be approximately 10,17, 25 M ⊙ for models of initial mass 20, 30, 40 M ⊙ , according to the result of Ref. 20).Taking this fact into account, we choose models with a core mass of M ≈ . . M ⊙ . This is achieved by setting K ≈ . × cm /s /g / . [Note that themass is approximately given by 4 . K/G ) / for the Γ = 4 / M = 4 . M ⊙ in the following. In this model, the ratio of the rotational kinetic energy T rot to thegravitational potential energy W is T rot / | W | ≈ . q ≡ cJ/GM is ≈ .
98, where M and J denote the mass and angular ormation Mechanism of Collapsar Black Hole Table I. The parameter set ( Γ , Γ , ρ nuc , Γ th ), the maximum gravitational mass M max , and the time t AH at which the first apparent horizon is formed.Model Γ Γ ρ nuc (g/cm ) Γ th M max ( M ⊙ ) t AH (ms)A 1.315 2.6 2.0 × × × momentum of the core, respectively.For the evolution, we adopt a hybrid EOS (see, e.g., Ref. 29)) in which thepressure consists of the sum of cold and thermal parts, as P = P cold + P th , where P cold = (cid:26) K ρ Γ , ρ ≤ ρ nuc ,K ρ Γ , ρ ≥ ρ nuc , (2.1)( ρ nuc denoting the nuclear density) and P th = ( Γ th − ρ ( ε − ε cold ). Here, ρ is the rest-mass density, ε and ε cold are the total and cold parts of the specific internal energy, Γ , Γ , and Γ th are constants, K = 5 × in cgs units, and K = K ρ Γ − Γ nuc ,respectively. By choosing a value of Γ slightly smaller than 4 /
3, the effect of thedepletion of the degenerate electron pressure due to photo-dissociation and electroncapture is qualitatively taken into account. The specific internal energy is set to thesame value as in the case of the Γ = 4 / ε = 3 Kρ / .The values of Γ are chosen as 1.315, 1.32, and 1.325 in this paper. For thecase of ordinary supernova simulations (i.e. the collapse of a stellar core with amass smaller than those in the present work), the value of Γ may be in the range1 . ≤ Γ ≤ . For the collapse of cores with larger masses, it has been pointedout that the value of Γ tends to be larger than in the lower mass case, due to theirlarger entropies (see, e.g., Ref. 30)). For this reason, we believe that the valuesof Γ adopted in this paper may be appropriate. The values of Γ and ρ nuc arechosen so that the maximum allowed gravitational mass of the cold neutron staris M max ≈ . M ⊙ , which is approximately equal to the highest pulsar mass evermeasured, i.e., M = 2 . ± . M ⊙ . We set Γ th = Γ for simplicity. With suchsmall values, the low efficiency of the shock heating due to the energy loss throughphoto-dissociation and neutrino emission is qualitatively taking in to account . Theadopted values for the parameter set ( Γ , Γ , ρ nuc , Γ th ) are listed in Table I.The temperature T is estimated from ε by solving the equation ε = ε gas + ε rad + ε ν + ε cold , (2.2)where ε gas , ε rad , and ε ν denote the specific internal energy of the gas, the radiation,and the neutrinos and are functions of ρ and T . We use the same forms for thesequantities as in Ref. 32), setting ε gas = 3 ρkT m p X nuc , (2.3) ε rad = 114 aT , (2.4) Y. Sekiguchi and M. Shibata ε ν = 218 aT Θ ( τ − τ crit ) , (2.5) ε cold = hc π (cid:18) π ρY e m p (cid:19) / , (2.6)where k is the Boltzmann constant, h Planck’s constant, m p the nucleon mass, X nuc the mass fraction of free nucleons, a the radiation constant, and Y e the electronfraction. Further, Θ denotes the Heaviside step function; i.e., for an optical depth τ larger than the critical optical depth τ crit ≈ /
3, we assume that the neutrinos areopaque with respect to nucleons and electrons and take into account their energydensities. (The definition of τ is given in § The so-called simple excision technique was implemented to follow the evolutionof the BH spacetime in a stable manner. (For recent calculations, see, e.g., Refs. 35)and 36).) The regridding technique was adopted during the collapse. The grid sizeand grid spacing in the final regridding were (3300 , ω, z ) and ∆x ≈ M/ ω = p x + y is the radius in cylindrical coordinates. §
3. Results
In the early stage, during which the central density of collapsing core is smallerthan the nuclear density, the collapse proceeds in an approximately homologous man-ner. Then, as the central density approaches the nuclear density, the collapse aroundthe equatorial plane is decelerated, due to the centrifugal force, while the collapsealong the rotational axis is relatively accelerated to form a flattened structure. Sincethe mass of the stellar core is much larger than the maximum allowed mass for thegiven EOS, a BH is formed directly. The formation of the BH is ascertained by find-ing the apparent horizon. After the apparent horizon is formed (for the formationtime, see Table I), a centrifugally supported, thin accretion disk is formed aroundthe BH (see the top left panel of Fig. 1): The disk rotates at approximately the Ke-plerian velocity. (see the bottom panel of Fig. 2). The evolution of the irreduciblemass of the BH, M irr ≡ p A/ π , and the disk mass M disk are shown in Fig. 3. Here A is the area of the apparent horizon.After the disk formation, shocks are formed at the inner part of the disk,converting the kinetic energy of the infall into thermal energy (See Ref. 38) fordiscussion of a similar phenomenon.) The gravitational potential energy releasedis E ∼ GM BH M disk /R ISCO ≈ × ergs, where M BH , M disk ∼ . M ⊙ , and R ISCO ≈ Gc − M BH are the black hole mass, the disk mass, and the radius of theinnermost stable circular orbit around the BH, respectively. As the thermal energyis stored, the disk height H increases. For H < R
ISCO , H in the vicinity of the BHis approximately determined by the balance relation P disk − P ram H ∼ GM BH ρ s HR , (3.1)where the left-hand and right-hand sides represent the pressure gradient and thegravity of the BH. The quantities ρ s , P disk , and P ram are the disk density near the ormation Mechanism of Collapsar Black Hole Fig. 1. Snapshots of density contours and velocity fields of the central region of the core in the x - z plane at t ≈ . x - z plane at 290.14 ms.All figures are for model B. surface, the pressure inside the disk, and the ram pressure of the infalling matter,respectively. Equation (3.1) provides the order of magnitude of the pressure as( P disk − P ram ) ∼ GM BH ρ s H R ∼ (cid:18) ρ s g / cm (cid:19) (cid:18) HR ISCO (cid:19) dyn / cm . (3.2) Y. Sekiguchi and M. Shibata
Fig. 2. Various quantities characterizing the disk at 290.14 ms for model B. The top panel plotsthe surface density Σ ≡ Σ/ g/cm (solid curve) and the optical depth τ ≈ κ ν Σ (dottedcurve). The middle panel plots the temperature evaluated in the equatorial plane. The bottompanel plots the angular velocity Ω (solid line) and the Keplerian angular velocity (dotted curve).The vertical dashed lines at R ≈
19 km indicate the expected location of the marginally stablecircular orbit around the BH.
The ram pressure can be expressed as P ram ∼ ρ f v ∼ (cid:18) ρ f g / cm (cid:19) dyn / cm (3.3)where ρ f and v f ∼ (2 GM BH /R ISCO ) / ∼ . c are the density and velocity of theinfalling matter.The density ( ρ disk ) and the temperature ( T disk ) inside the disk eventually increaseto ∼ g/cm and ∼ K (and hence, P disk ∼ dyn/cm ), while the rampressure ( P ram ) decreases to . . P disk ( ≪ P disk ), since the density of the infallingmatter ( ρ f ) decreases to . g/cm . Consequently, H increases to ∼ R ISCO for ρ s ∼ g/cm . Then, the approximate balance relation becomes( P disk − P ram ) ∼ GM BH ρ s H . (3.4)Since the binding due to the gravitational force decreases as H increases, the diskexpands preponderantly in the z -direction, forming strong shock waves (see the topright and middle left panels of Fig. 1). The propagation speed of the shock wavesis ≈ . c , i.e., mildly relativistic. As the shock waves propagate, the mass accretionrate to the BH significantly decreases (see Fig. 3). After the shock propagation, alow density funnel region is formed around the rotational axis (see the middle right ormation Mechanism of Collapsar Black Hole Fig. 3. Evolutions of the irreducible mass M irr and the disk mass M disk . The long dashed, solid,and dashed lines represent the results for models A, B, and C, respectively. The horizontal lineat M ≈ . M ⊙ indicates the total mass. and bottom left panels of Fig. 1). In addition, the funnel region is surrounded by adense, hot wall (see the bottom right panel of Fig. 1). The temperature in the funnelincreases to approximately the same level as that in the wall through shock heating,and hence such a low density funnel can be supported by the thermal pressure.Figure 2 displays the surface density, the optical depth, the temperature, andthe angular velocity of the disk at t = 290 .
14 ms. Due to the shock heating, thetemperature increases to T ≈ × K, which is about a factor of 4–5 larger thanthat in the unshocked regions (see the middle panel of Fig. 2). The surface densityaround the rotational axis is Σ ≈ × g/cm , while it is Σ ∼ g/cm in thedisk. This reflects the formation of the funnel structure formation (see the top panelof Fig. 2). Assuming that neutrino opacity is κ ν = 7 × − ( T / K) (see, e.g.,Di Matteo et al. 2002), we compute the optical depth τ = κ ν Σ . Then we find thatthe disk is optically thick for a radius of R ≈ Y. Sekiguchi and M. Shibata
In our models, a highly relativistic fireball will be produced by the neutrino pair( ν ¯ ν ) annihilation around the rotational axis. Here, let us approximately estimatethe energy deposition rate ˙ E ν ¯ ν through ν ¯ ν annihilation. In the diffusion limit, theneutrino flux is given by F ν ≈ N ν σT κ ν Σ , (3.5)where N ν = 3 is the number of neutrino species and σ is the Stefan-Boltzmannconstant. The neutrino luminosity from an optically thick disk is then L ν ≈ πR F ν ∼ × T Σ − (cid:18) R disk
70 km (cid:19) ergs / s (3.6)where T ≡ T / K and Σ ≡ Σ/ g/cm . According to the results ofSetiawan et al., the expected energy deposition rate through the ν ¯ ν annihilationwould then be ˙ E ν ¯ ν ∼ ergs/s. The low density funnel region above the BH willbe a favorable place for an efficient production of ee + pairs through ν ¯ ν annihilation.The formation of a dense, hot wall surrounding the funnel will play an impor-tant role in collimating relativistic outflows. In the absence of this wall, energeticoutflows fail to be collimated, as the outflows travel through the stellar mantle. Theshock waves formed at the birth of the BH are likely to sweep the matter alongthe rotational axis, reaching the stellar surface before the main relativistic outflows,which will be produced in the accretion phase (e.g., Ref. 39)), are driven.For larger values of | Γ − / | , the time at the onset of the shock wave propagationis delayed, because the rate of increase of the disk pressure as well as the rate ofdecrease of the ram pressure are smaller, reflecting the less homologous nature of thecollapse. This also results in the fact that it takes a longer time for the BH to relaxto a stationary state (see Fig. 3). However this ’time delay’ is at most ∼
10 ms for Γ & . α ( GM ˙ m/R ) ≫ L ν ( ∼ ergs / s), where α . m is the massaccretion rate. This gives ˙ m ≫ L ν R/αGM ≈ α − M ⊙ s − . According to our results,˙ m ≫ M ⊙ s − at least for ( t − t AH ) .
10 ms. Thus, the neutrino cooling doesnot play an important role in the thermal-energy-store phase, unless the conversionefficiency is extremely low, i.e., α ≪ .
1. To summarize, our results are universal, atleast qualitatively, for any of the EOSs and microphysical processes.We also performed simulations for 0 . . q . .
2, as well as for weakly differen-tially rotating cases and found essentially the same results for all these parametervalues. For larger values of q ≈ .
2, more collimated shock waves are formed. Forsmaller values of q ≤ .
8, a massive disk, which is essential for shock formation, isnot formed. For sufficiently large values of q > .
2, on the other hand, a BH is notpromptly formed, due to the effect of the centrifugal force. In this case, a BH maybe formed after a sufficient amount of angular momentum is transported. ormation Mechanism of Collapsar Black Hole M irr ≈ . M ⊙ is born at first for any of the EOSs. This value is approximatelydetermined by the maximum allowed mass of the rigidly rotating configuration forthe adopted EOSs. The figure also reveals that approximately 95% of the totalmass is eventually swallowed by the BH for any of the adopted EOSs. Since thetotal angular momentum is conserved, the angular momentum of the BH can beindirectly estimated from that of the disk, which is ∼
20% of the total angularmomentum. Thus, the final spin parameter of the BH is estimated as ≈ . §
4. Summary
Motivated by recently developed single-star evolutionary models of LGRB pro-genitors, we performed fully general relativistic simulations for rapidly rotating stel-lar core collapse to a BH. We found that a BH is directly formed as a result of thecollapse of a sufficiently massive progenitor. Soon after the BH formation, a disk ofdensity ∼ g / cm is formed around the BH. The subsequent infall of matter fromoutside produces shocks at the surface of the disk, and thermal energy is generated,heating the disk. The thermal pressure eventually reaches ∼ dyn / cm , whichis much larger than the ram pressure of the infalling matter, and hence it is used tosustain the vertical gravitational force. Because of the increasing thermal energy, thedisk and shock front gradually expand, and when the scale height is comparable tothe disk radius, the shock front starts to expand with mildly relativistic speed. Thestrong shock waves sweep matter around the rotational axis and heat up the diskmatter to ∼ K. Then, a low-density ( ρ . g / cm ) funnel region surroundedby a hot, dense wall is formed. Formation of such a structure in advance of thesubsequent quasi-stationary accretion from the hot disk and the resulting relativisticoutflows will be quite favorable for producing the fireball and LGRBs.In the present simulation, we did not incorporate magnetic fields. It shouldbe noted that even in the absence of magnetic stress, funnel structure, which willbe essential for producing collimated jets, is formed. However, magnetic fields canfurther promote the formation of such a funnel if the field strength is larger than ∼ G (i.e., the magnetic pressure is larger than ∼ dyn / cm , which iscomparable to the ram pressure of the infalling matter). Study of magnetic fieldeffects is left as a future project.Although our treatments of the EOS and microphysical processes in the presentwork are approximate, we believe that the qualitative features of the BH formationand subsequent shock formation found here will be universal. To confirm this processmore rigorously, more detailed quantitative studies with a realistic EOS, includingrelevant microphysics, is needed. Such a study is in progress, and the results will bereported in the near future.0 Y. Sekiguchi and M. Shibata
Acknowledgements
We thank R. Takahashi and K. Maeda for valuable comments and discussions.Numerical computations were performed on the FACOM VPP5000 at the data anal-ysis center of NAOJ and on the NEC SX-8 at YITP. This work is supported by JSPSResearch Grant for Young Scientists (No. 1611308) and by Monbukagakusho Grant(No. 17030004).
References
1) J. Gorosabel et al., Astron. Astrophys. (2003), 127.2) L. Christensen J. Hjorth, and J. Gorosabel, Astron. Astrophys. (2004), 913.3) P. Jakobsson et al., Mon. Not. R. Astron. Soc. (2005), 245.4) T. J. Galama et al., Nature (1998), 670.5) S. R. Kulkarni et al., Nature (1998), 663.6) J. Hjorth et al., Nature (2003), 847.7) K. Z. Stanek et al., Astrophys. J. (2003), L17.8) K. S. Kawabata et al., Astrophys. J. (2003), L19.9) M. Modjaz et al., Astrophys. J. (2006), L21.10) J. Sollerman et al., Astron. Astrophys. (2006), 503.11) A. Zeh, S. Klose and D. H. Hartmann, Astrophys. J. (2004), 952; in
Proc. 22nd TexasSymposium on Relativistic Astrophysics , ed. P. Chen et al. (Stanford, SLAC, 2005).12) S. E. Woosley, Astrophys. J. (1993), 273.13) B. Paczy´nski, Astrophys. J. (1998), L45.14) A. MacFadyen and S. E. Woosley, Astrophys. J. (1999), 262.15) W. Zhang and S. E. Woosley, Astrophys. J. (2004), 365.16) N. Langer, Astron. Astrophys. (1998), 551.17) P. Podsiadlowski et al., Astrophys. J. (2004), L17.18) C. L. Fryer and A. Heger, Astrophys. J. (2005), 302.19) S.-C. Yoon and N. Langer, Astron. Astrophys. (2005), 643; in
ASP Conf. Ser., StellarEvolution at Low Metallicity: Mass-Loss, Explosions, Cosmology , ed. H. Lamers et al. (SanFrancisco: ASP, 2006); astro-ph/0511222.20) S.-C. Yoon, N. Langer and C. Norman, Astron. Astrophys. (2006), 199.21) S. E. Woosley and A. Heger, Astrophys. J. (2006), 914; in
AIP Conf. Proc., GammaRay Bursts in the Swift Era , ed. S. S. Holt, N. Gehrels, and J. Nousek (San Francisco:ASP, 2006); astro-ph/0604131.22) J. S. Vink and A. de Koter, Astron. Astrophys. (2005), 587.23) A. Maeder, Astron. Astrophys. (1987), 159.24) J. P. U. Fynbo et al., Astron. Astrophys. (2003), L63.25) R. L. C. Starling et al., Astron. Astrophys. (2005), L21.26) J. Gorosabel et al., Astron. Astrophys. (2005), 711.27) K. Z. Stanek et al., Acta. Astron. (2006), 333.28) S. L. Shapiro and S. A. Teukolsky, Black Holes, White Dwarfs, and Neutron Stars (Wi1eyinterscience, New York, 1983).29) T. Zwerger and E. M¨uller, Astron. Astrophys. (1997), 209.30) K. Nakazato, K. Sumiyoshi, and S. Yamada, Astrophys. J. (2006), 519.31) D. J. Nice et al., Astrophys. J. (2005), 1242.32) T. Di Matteo, R. Perna and R. Narayan, Astrophys. J. (2002), 706.33) M. Shibata, Phys. Rev. D (2003), 024033.34) M. Alcubierre and B. Br¨ugmann, Phys. Rev. D (2001), 104006.35) M. Shibata and K. Taniguchi, Phys. Rev. D (2006), 064027.36) M. Shibata, M. D. Duez, Y. T. Liu, S. L. Shapiro, and B. C. Stephens, Phys. Rev. Lett. (2006), 031102.37) M. Shibata and S. L. Shapiro, Astrophys. J. (2002), L39.38) W. H. Lee and E. Ramirez-Ruiz, Astrophys. J. (2006), 961.39) W. H. Lee, E. Ramirez-Ruiz and D. Page, Astrophys. J. (2005), 421.40) D. Proga et al., Astrophys. J. (2003), L5. ormation Mechanism of Collapsar Black Hole
41) T. Takiwaki et al., Astrophys. J. (2004), 1086.42) S. Setiawan, M. Ruffert and H.-Th. Janka, Mon. Not. R. Astron. Soc. (2004), 753;Astron. Astrophys.458