Alternative possibility of GW190521: Gravitational waves from high-mass black hole-disk systems

We evolve high-mass disks of mass $15-50M_\odot$ orbiting a $50M_\odot$ spinning black hole in the framework of numerical relativity. Such high-mass systems could be an outcome during the collapse of rapidly-rotating very-massive stars. The massive disks are dynamically unstable to the so-called one-armed spiral-shape deformation with the maximum fractional density-perturbation of $\delta \rho/\rho \gtrsim 0.1$, and hence, high-amplitude gravitational waves are emitted. The waveforms are characterized by an initial high-amplitude burst with the frequency of $\sim 40-50$ Hz and the maximum amplitude of $(1-10)\times 10^{-22}$ at the hypothetical distance of 100 Mpc and by a subsequent low-amplitude quasi-periodic oscillation. We illustrate that the waveforms in our models with a wide range of the disk mass resemble that of GW190521. We also point out that gravitational waves from rapidly-rotating very-massive stars can be the source for 3rd-generation gravitational-wave detectors for exploring the formation process of rapidly-rotating high-mass black holes of mass $\sim 50-100M_\odot$ in an early universe.


I. INTRODUCTION
In their third observational run, Advanced LIGO and Advanced Virgo announced the discovery of a high-mass binary black hole, GW190521 [1]. The estimated mass of two black holes are m 1 = 85 +21 −14 M and m 2 = 66 +17 −18 M with the estimated distance to the source of 5.3 +2. 4 −2.6 Gpc. Stellar evolution calculations predict that black holes of mass gap ∼ 50-120M should not be formed from the massive stars because of the presence of pulsational pair instability and pair-instability explosion for the massivestar evolution with a range of the initial stellar mass [2][3][4]. For this reason, the announcement of GW190521 by the LIGO-Virgo collaboration has been attracting much attention and several possibilities for forming the black holes in the mass gap have been suggested (see, e.g., Refs. [5][6][7][8][9]). Subsequently Nitz and Capano [10] reanalyzed the data with different priors and waveform models from those in Ref. [1]. They found that more probable component mass might be m 1 = 168 +15 −61 M and m 2 = 16 +33 −3 M , and thus, the mass of them might not be in the mass gap range of 50-120M . A very-massive black hole of mass 150M could be naturally interpreted as a remnant of the collapse of a very-massive star (see below). Another interesting point in their work is that their maximum-likelihood waveform indicates that after the burst waveform of a few cycles with the highest amplitude, gravitational waves might be emitted: They show that with their maximum-likelihood waveform in which the post-burst waves are present, the signal-to-noise ratio is higher than in the absence of the post-burst waves. For binary black hole mergers, the gravitational waveforms after a peak amplitude is reached should be character-ized by a ringdown waveform and the amplitude should approach zero quickly. The possible presence of gravitational waves after the burst-waves emission suggests that the source of GW190521 might not be the merger of high-mass binary black holes.
Besides the frequency of gravitational waves, the waveforms such as that of GW190521, which are composed initially of burst waves of a few wave cycles and subsequently of a small-amplitude oscillation, have been often reported in the results of numerical simulations for the systems dynamically unstable to non-axisymmetric deformation. For example, during the collapse of very rapidly rotating stars, the non-axisymmetric deformation can occur and the resulting gravitational waveform can be qualitatively similar to that of GW190521 (besides the frequency: see, e.g., Fig. 15(b) of Ref. [11]). Thus, it would not be non-sense to explore other astrophysical possibilities for interpreting the observational results of GW190521, supposing that post-burst gravitational waves indicated in Ref. [10] may indeed be real signals.
Based on this motivation, we perform a numericalrelativity simulation for the system of a massive black hole and a massive disk. We choose the black hole of initial mass M BH,0 = 50M and the disk of rest mass M disk = 15-50M in this paper. Here the blackhole mass of 50M is chosen in order to reproduce the typical frequency of gravitational waves for GW190521 which is ∼ 50 Hz [1] (note that the frequency is approximately proportional to M −1 BH,0 ). Such a high-mass black hole-disk system could be formed during the collapse of very-massive metal-poor stars of initial mass larger than ∼ 190-260M [12][13][14][15][16][17][18]. Here the critical initial mass depends on the angular momentum of the progenitor stars [15,17]. For these stars, the mass of the carbon-oxygen core formed at the final stage of the stellar evolution is M CO 125M and the collapse is triggered by the pair instability. The collapse of such heavy cores leads to the direct formation of a high-mass black hole, avoiding the pair-instability explosion, while for 70M M CO 125M , a pair-instability supernova is the final outcome, leaving no compact objects at the center [15,17].
Here, it is important to emphasize that the mass of the black hole during its formation may not be always equal to the entire mass of the carbon-oxygen core and it could be smaller than M CO in the presence of rapid rotation. One reason for this is that a fraction of the mass can form disks (and ejecta from the disks). Indeed, our previous work in numerical relativity shows that this is the case: Reference [19] illustrates that after the collapse of a moderately rapidly-rotating very-massive star of the initial mass 320M and M CO ≈ 150M , which is obtained by a stellar evolution calculation [15], the outcome is a black hole of mass ∼ 130M surrounded by a disk of mass ∼ 20M . Our subsequent work further illustrates that in the presence of the neutrino cooling effect which was absent in Ref. [19], the initial mass of the black hole formed during the early stage of the collapse is likely to be substantially smaller than M CO [20] because the collapse of the central region is significantly accelerated in a runaway manner toward the earlier formation of a black hole. The temporal mass of the black hole and disk formed during the collapse could also depend strongly on the angular velocity profile of the pre-collapse carbonoxygen core, which is not well understood by the stellar evolution calculations. Thus, in this paper, we suppose the case in which a high angular momentum is present in the carbon-oxygen core just prior to the collapse, and a system of a black hole and a disk, for both of which the mass is several tens of solar mass, could be formed during the collapse of rapidly-rotating very-massive stellar cores.
By numerically evolving high-mass black hole-disk systems as a plausible outcome formed during the collapse of rapidly-rotating very-massive stars, we show that the massive disks in such systems are dynamically unstable to the one-armed spiral-shape deformation as often found in previous works for low-mass disks [21][22][23][24][25] (but see also Ref. [26] for the possible importance of the magnetic field effects for low-mass disks). Then, we demonstrate that gravitational waves emitted from the nonlinearly perturbed disks are indeed composed of a highamplitude initial burst and a subsequent low-amplitude quasi-periodic oscillation. For M BH,0 = 50M , the frequency of the initial burst gravitational waves is ∼ 40-50 Hz for a wide range of the disk mass (for higher disk mass, the frequency is lower), and we find that the waveforms for several models are similar to the waveforms shown in Refs. [1,10], depending only weakly on the disk mass.
We note that during an early formation stage of a black hole and a massive disk, the configuration of the system is likely to be approximately axisymmetric until the onset of the one-armed spiral-shape deformation. The amplitude of gravitational waves emitted in the (approximately) axisymmetric stage is not as high as that emitted from the non-axisymmetric system even at the moment of the black-hole formation (see, e.g., Ref. [20]). The non-axisymmetric instability in the disk is likely to set in when the growth timescale of the dynamical instability becomes as short as the formation timescale of the disk, which would be approximately equal to the local free-fall timescale. The amplitude of gravitational waves emitted by the non-axisymmetric deformation of the massive disk should be by ∼ 1 orders of magnitude higher than that emitted during the axisymmetric stage (e.g., Refs. [11,27]). Thus, for low signal-tonoise events, only gravitational waves emitted in the nonaxisymmetric stage would be observable. This motivates us to perform simulations focusing only on the black holedisk system for exploring gravitational waves from the collapse of rapidly-rotating very-massive stars as a first step.
This paper is organized as follows. In Sec. II, we summarize the set-up of the present numerical-relativity simulations. In Sec. III, the unstable evolution of the high-mass disk surrounding the high-mass black hole and resulting gravitational waveforms are presented. Section IV is devoted to discussions. Throughout this paper, G and c denote the gravitational constant and the speed of light, respectively.

II. SET-UP OF NUMERICAL SIMULATIONS
The set-up of the present numerical simulations is essentially the same as that in our previous papers on axisymmetric simulations for the massive black hole-disk systems [28,29], but in this paper, we perform threedimensional hydrodynamics simulations. We numerically solve Einstein's equation and the hydrodynamics equations. Einstein's equation is solved using the original version of the Baumgarte-Shapiro-Shibata-Nakamura formalism [30] together with the puncture formulation [31], Z4c constraint propagation prescription [32], and 6thorder Kreiss-Oliger dissipation. For hydrodynamics, we do not consider the neutrino emission and weak interaction throughout this paper, although our implementation can take into account the neutrino effects in the same manner as in Refs. [33,34]. In the presence of the neutrino emission, the disk shrinks by the neutrino cooling and its compactness is increased. This effect makes the dependence of the growth timescale of the disk instability and resulting gravitational waveforms on the disk morphology obscure. Thus we decided to switch off the neutrino effects in this work. In this prescription, only the advection part for the equations of the electron fraction, Y e , is taken into account.
The quantities of black holes (mass and spin) are determined from their area and circumferential radii of ap-TABLE I. Initial conditions for the numerical simulation. Described are the model name, initial black-hole mass (MBH,0), rest mass of the disk (M disk ), dimensionless spin of the black hole (χ), the coordinate radii at the inner and outer edges of the disk (rin and rout), entropy per baryon (s/k with k the Boltzmann's constant), the maximum density of the disk (ρmax), and the orbital period at the maximum density (P orb ). The units are M for the mass, GMBH,0/c 2 ≈ 73.84(MBH,0/50M ) km for rin and rout, 10 10 g/cm 3 for ρmax, and GMBH,0/c 3 for P orb , respectively.
Model MBH,0 M disk χ rin rout s/k ρmax P orb parent horizons [35], assuming that these quantities are written as functions of the mass and spin in the same formulation as in the vacuum Kerr black hole. Gravitational waves are derived through the extraction of the outgoing component of the complex Weyl scalar [36]. As in our previous work [28,29], we employ a tabulated equation of state based on the DD2 equation of state [37] for a relatively high-density part and the Timmes (Helmholtz) equation of state for the low-density part [38]. In this tabulated equation of state, thermodynamics quantities such as ε, P , and h are written as functions of ρ, Y e , and T where ε, P , h(= c 2 +ε+P/ρ), ρ, and T are the specific internal energy, pressure, specific enthalpy, rest-mass density, and matter temperature, respectively. We choose the lowest rest-mass density to be 0.1 g/cm 3 in the table.
Axisymmetric equilibrium states for black hole-disk systems are prepared as the initial conditions [28,29] using the method of Ref. [39]. As in Refs. [28,29], we determine the angular velocity profile from the relation where j = c −2 hu ϕ is the specific angular momentum. Ω is the angular velocity defined by u ϕ /u t with u µ the fluid four-velocity. n is a constant that determines the profile of the angular velocity, for which we choose n = 1/7 following our previous papers [28,29] to obtain a profile of Ω close to the Keplerian one that does not result in the geometrically thick disks. For obtaining the initial conditions, we assume a plausible relation between ρ and Y e for the high-density matter in the same form of ρ(Y e ) as in Ref. [28]. For this model, the value of Y e is 0.07 in a high density region with ρ ≥ 10 11 g/cm 3 , 0.07 − (43/400) log 10 [ρ/(10 11 g cm −3 )] for 10 7 g/cm 3 ≤ ρ ≤ 10 11 g/cm 3 , and settles to 1/2 for ρ ≤ 10 7 g/cm 3 , for which the effect of the electron degeneracy to the neutron richness is weak. In addition, we assume that the specific entropy, s, is initially constant, in order to obtain the first integral of the Euler equation easily. For the high-mass disks, we choose it s/k = 11-17 in this paper where k is the Boltzmann's constant: By adjusting this value, the compactness for given mass is controlled (for the higher entropy, the disk becomes less compact for a given value of the disk mass).
We note that our setting for the disk profile is idealized. For more realistic work, it is obviously necessary to perform a simulation started from rapidly-rotating verymassive progenitor stellar cores. However, the collapse simulation with a wide variety of the parameters (e.g., angular momentum profile) is computationally very expensive. The purpose of the present work is to understand the nature of the non-axisymmetric instability of the massive disks orbiting a black hole and resulting gravitational waves emitted. For this, we believe that the present setting is acceptable.
The initial black-hole mass is fixed to be M BH,0 = 50M while a wide range of the disk mass, M disk ≈ 15-50M , is employed. The initial dimensionless spin of the black hole is set to be approximately 0.8 supposing that the very-massive stellar cores, which form a massive disk during the collapse, should be rapidly rotating. We set the inner coordinate radius of the disk to be close to the radius of the innermost stable circular orbit around the rapidly spinning black hole as r in ≈ 2GM BH,0 /c 2 , and the outer edge is varied in the range of r out /(GM BH,0 c −2 ) = 21-53 (see Table I). These are the plausible sizes of the dense region of the disk formed during the collapse of rotating very-massive stellar cores [19]. With these setups, the maximum density of the disks is of the order of 10 10 g/cm 3 (see Table I). Thus, the matter in the high-density region is mildly neutronrich.
Numerical simulations are performed with a fixed mesh-refinement implementation used for the neutronstar mergers [40][41][42]. We prepare 10 refinement levels with each region composed of 241×241×121 uniform grid points for x-y-z (the reflection symmetry with respect to the equatorial plane is imposed). For the finest level, the grid spacing is 0.018GM BH,0 /c 2 , and the outer boundary along each axis is located at L ≈ 1106GM BH,0 /c 2 for all the simulations. We checked that with such setting the black hole can be accurately evolved within the fractional change of its mass and dimensional spin of 0.1% before the onset of the mass accretion resulting from the disk instability.
FIG. 1. Snapshots of the rest-mass density profile of the disk in the equatorial plane for model M14. P orb denotes the initial orbital period of the disk at the density maximum (see Table I). The color bar shows log ρ with the unit of ρ being g/cm 3 . See also http://www2.yukawa.kyoto-u.ac.jp/˜masaru.shibata/BHdiskM14.mp4 for the corresponding animation.
We do not add any perturbation initially, except for numerical noises. By a small random perturbation automatically introduced in numerical evolution, the dynamically unstable deformation of the massive disks is enhanced for all the simulations.

A. Evolution of dynamically unstable disks
Irrespective of the initial conditions prepared in this paper, the disk is dynamically unstable to the one-armed spiral-shape deformation. Figure 1 plots the snapshots of the rest-mass density profile of the disk in the equatorial plane for model M14 (see Table I for the initial condition of this model). It is found that a small numerical noise introduced during the simulation triggers the subsequent exponential growth of the one-armed spiral-shape deformation mode. We note that in the present context (i.e., in fully general relativistic computation), the enhanced density perturbation in the disk enforces the black hole to slightly shift to the higher-density side of the disk, because the center of mass of the system should be preserved. This enhances the gravitational attraction between the black hole and the high-density part of the disk, and accelerates the growth of the one-armed density perturbation. Note that this effect cannot be taken into account in the simulations with fixed background spacetime, and thus, a fully general relativistic simulation is needed to explore the one-armed spiral-shape instability for high-mass disks.
After the significant growth of the spiral-shape deformation, the angular momentum transport inside the disk, which is caused by the gravitational torque from the non-axisymmetric density perturbation as well as by the black-hole tidal field, is enhanced. As a result, a significant mass accretion onto the black hole proceeds. During the initial growth of the one-armed spiral-shape deformation until its saturation, 20-30% of the total rest mass of the disk falls into the black hole. Subsequently, a weakly perturbed disk is formed, and from such a disk, quasi-steady mass inflow to the black hole continues (e.g., Refs. [24,25]). For all the cases, the black hole spins up by the mass accretion and the dimensionless spin always becomes higher than 0.85 after a substantial fraction of the mass falls into the black hole. In particular for the massive and compact disk models such as M12, H14, E15, and E16, the dimensionless spin exceeds 0.92.
Because of the presence of the non-axisymmetric structure in the disk, the gravitational torque is exerted continuously to the matter in the outer part of the disk, which is expanded significantly. Even in the present sim- . The time is shown in units of P orb , which denotes the initial orbital period of the disks at the maximum rest-mass density (see Table I).
ulations of relatively short duration, we find that 1-2% of the total mass of the disk is ejected from the system. As a result of the mass ejection (as well as of the mass accretion into the black hole), the maximum density of the disk decreases (compare the first and last panels of Fig. 1). With a longer-term simulation taking into account the viscous effect, we could find that a more fraction of the disk matter will be ejected from the system. However, exploring such a long-term mass ejection process is beyond the scope of this paper. Figure 2 displays the evolution of the mode amplitude of the density perturbation for the selected models. To define the mode amplitude, we first extract each nonaxisymmetric mode of the density perturbation by calculating with m = 0-4 and ϕ = tan −1 (y/x). Here, the volume integral is performed only for the outside of the apparent horizon. We then define the normalized mode amplitude by In the following, we refer to δ m as the mode amplitude.
For all the simulations, δ m initially increases by the numerical noise. However, after δ 1 reaches ∼ 10 −6 , the one-armed spiral instability, which is developed from the initial seed perturbation that results from the numerical noise, makes the dominant perturbation. Thus, in the following, we focus only on the growth rate for δ 1 ≥ 10 −5 .  II. The growth timescale of the dynamical instability in units of P orb , the maximum value of δ1, and the signal-to-noise ratio (SNR) for gravitational waves in the detection at the hypothetical distance to the source of 100 Mpc with the most optimistic direction assuming the detection by a single advanced LIGO detector with the designed sensitivity. The last row describes the type of the waveform (see the text for the classification).  Figure 2(b)-(d) display the m = 1 mode amplitude, δ 1 , as a function of time for several models. These show that δ 1 increases exponentially with time for the stage of δ 1 10 −5 , for which we extract the growth timescale, τ dyn , by fitting δ 1 in the form of ∝ exp(t/τ dyn ). It is found that the exponential growth timescale for the models chosen in this paper is 0.5-2.5 times of the initial orbital period, P orb , of the matter at the highest-density region of the disk (cf. Table II). This implies that the instability found here is indeed the dynamical instability. Figure 2(b) compares the m = 1 mode amplitude as a function of time for models M12-M15 for which the compactness of the disk defined by GM BH /(c 2 r out ) is different for a fixed value of the disk mass. In Fig. 2(c) we also compare the m = 1 mode amplitudes between E15 and E16, between K13 and K14, and between L11 and L12 to see the dependence of the growth timescale of the instability on the disk compactness for a fixed value of the disk mass. It is always found that the growth timescale of the dynamical instability (in units of P orb ) is shorter for the more compact disk models. We note that if we measure the growth timescale in units of GM BH,0 /c 3 , this trend is further remarkable. Thus, we can conclude that more compact disks are more subject to the dynamical instability. On the other hand, this result also suggests that a sufficiently less compact disk may be dynamically stable. Figure 2(b) also shows that the maximum value of δ 1 is larger for more compact disks (see also Table II). Figure 2(d) compares the m = 1 mode amplitude for models L12, J13, K13, M13, H14, and E15 for which c 2 r out /(GM BH ) is ∼ 30. This shows that the growth timescale of the dynamical instability is shorter for more massive disks with a given disk compactness. In particular for M disk /M BH,0 0.6 the growth timescale can be shorter than P orb . On the other hand, this figure suggests that for a disk of mass smaller than a critical value the disk may be dynamically stable [24]. This figure also illustrates that the maximum value of δ 1 is larger for more massive disks (see also Table II).
Irrespective of the models, after δ 1 reaches the maximum of δ 1 0.1, i.e., the density perturbation becomes non-linear, the growth of the density perturbation saturates. Then, the degree of the deformation is approximately preserved with the density perturbation of δ 1 = 0.01-0.1. This result is qualitatively the same as, e.g., in the previous results [24,25]. Because the gravitational torque resulting from the one-armed density perturbation transports the angular momentum outward and thus the disk compactness decreases, the value of δ 1 will become smaller in the long-term evolution of t 10P orb . In Table II, we summarize the exponential growth timescale, τ dyn , and the maximum value of δ 1 , δ 1,max . Broadly speaking, for more compact or more massive disks, τ dyn /P orb is smaller and δ 1,max is larger. As shown in Sec. III B, the dependence of these quantities on the compactness and mass of the disks is reflected in the gravitational waveforms. We will find that three different types of gravitational waveforms are produced by varying r out and the disk mass.
Before closing this subsection, we note the following point. In the presence of magnetic fields, the magnetorotational instability should occur [43] and subsequently, a turbulence could be excited in the disk [26,44]. Then the turbulent viscosity is enhanced, and the angular momentum in the disk would be redistributed. If the viscous angular momentum transport efficiently works in the disk, the growth of the one-armed spiral instability studied in this paper should be suppressed [26]. The timescale of the viscous angular momentum transport is where R and ν denote the cylindrical radius in the disk and the shear viscous coefficient, respectively. For the geometrically thin disks, ν can be rewritten as ν = α ν c s H where α ν is the dimensionless alpha parameter [45], c s is the sound velocity, and H is the scale height of the disk. Numerical simulations for accretion disks have shown α ν = O(10 −2 ) (e.g., Ref. [44]). Then the ratio of τ vis to P orb is where we used the approximate relation of the force balance with respect to the vertical direction of the disk as HΩ ≈ c s . Because R/H 2 for the disk employed in this work, it is found that τ vis /P orb > 10 as long as α ν is not unusually large, i.e., α ν 0.05. Thus, for the massive and compact disks of τ dyn 2.5P orb which we study in this paper, the timescale of the viscous angular momentum transport is much longer than the growth timescale of the one-armed spiral-shape instability. This is in particular the case for the models that result in the type I waveforms which we are most interested in (see Sec III B). Therefore, the magnetic/viscous effects could be safely neglected for the models employed in this paper. However, as shown in Ref. [26], the magnetic/viscous effects are likely to play an important role for low-mass or or less compact disks for which τ dyn 10P orb .

B. Gravitational waves
Figures 3 and 4 display gravitational waves observed along the rotation axis (i.e., perpendicular to the orbital plane) for all the models employed in this paper. Figure 3 is for relatively less massive disks with 15M ≤ M disk ≤ 25M and Fig. 4 For many of the models studied in this paper (the models except for L12, J13, K14, and M15), the waveforms are characterized by initial burst waves of a high amplitude and subsequent quasi-periodic waves of a low amplitude. We refer to this waveform type as type I. The initial burst waves are most efficiently emitted when the value of δ 1 becomes the maximum. The number of the wave cycle of the burst waves is 1-3 for this type, and is smaller for more compact or more massive disks for a given value of c 2 r out /(GM BH,0 ). The maximum amplitude of the burst gravitational waves is 10 −22 -10 −21 at the hypothetical distance to the source of D = 100 Mpc and is higher for more compact disks  with smaller values of c 2 r out /(GM BH ) or for more massive disks. On the other hand, the maximum amplitude of the quasi-periodic waves does not depend on the value of c 2 r out /(GM BH ) as strongly as that of the initial burst waves and is always 10 −22 for D = 100 Mpc. For the waveforms of models L12, J13, K14 and M15, the amplitude of the initial burst gravitational waves is only by a factor of 2-3 higher than the subsequent quasiperiodic waves in contrast to the models of the type I waveform. We refer to this waveform type as type II. For the type II waveform, the maximum amplitude is 2 × 10 −22 for D = 100 Mpc irrespective of the disk mass. In terms of δ 1,max , the type I waveform is obtained for δ 1,max 0.15 while the type II waveform is for δ 1,max 0.15 (see Table II).
We can consider that there are two types in the type I waveform. For the models of the very compact or massive disks, the burst waves are characterized by only one wave of a very high amplitude (see the waveforms for models M12, H14, E15, and E16). We refer to this waveform as type Is. The short duration of the highest amplitude phase is due to the fact that the peak values of δ 1 for these cases are relatively large and hence by the gravitational torque exerted from the non-axisymmetric structure, a large amount of the disk matter is ejected from the highdensity region of the disk settling to a less compact disk in a short timescale. For other type I models (L11, J12, K13, M13, M14, H15, H16, E17), the gravitational torque at the moment that δ 1 is highest is not large enough to expel a large fraction of matter from the high-density region of the disk. As a result, the burst waves are emitted for a couple of wave cycles. We refer to the waveform type as type Iw for this case. In terms of δ 1,max , the type Iw waveform is obtained only for δ 1,max 0.25. We note that the waveforms for models M13 and H15 may be classified as type Is. Thus, we consider that these models are located near the boarder for distinguishing types Iw and Is (note that the waveform should change continuously with the change of the compactness and mass of the disk, and hence, it would not be possible to definitely classify the waveform into a particular type).
The numerical results presented here illustrate that the compactness and mass of the disk (which are determined by the angular momentum profile of the progenitor star prior to the collapse) are the key parameters for determining the type of the waveforms. Specifically, for the models employed in this paper, the type II waveforms are obtained broadly for M disk /M < −10 + c 2 r out /(GM BH,0 ), and otherwise, the type I waveform is the result (see Fig. 5). Also, broadly for M disk /M > 3 + c 2 r out /(GM BH,0 ), the waveform becomes the type Is. We note that this would be quantitatively valid only for c 2 r out /(GM BH,0 ) ≥ 20 and for the initial angular momentum profile employed in this paper. However, we infer that qualitatively similar relations would be satisfied irrespective of the disk configuration for the massive disks.
For the type I waveform, the amplitude of the burstwave part is appreciably larger than the quasi-periodic waves subsequently emitted. Thus, if gravitational waves of this type are detected by gravitational-wave detectors with a low signal-to-noise ratio, it would not be very easy to confirm the detection of the quasi-periodic wave part, and essentially, the attention would be paid only to the burst-waves part. For the case of GW190521, the detection with a sufficiently high signal-to-noise ratio appears to be done only for the first ∼ 2-3 wave cycles that have the highest amplitude [1,10]. This waveform is similar to the type Iw waveform. This indicates a possibility that the source of GW190521 may not be a merger of highmass binary black holes but an outcome of the collapse of a rapidly-rotating very-massive stellar core leading to a black hole and a massive disk with moderate compactness. We note that for a variety of the mass and the compactness of the disk, the type Iw waveform resembles GW190521; that is, as long as a relation between the mass and the compactness of the disk is satisfied (see Fig. 5), no artificial fine-tuning is necessary to reproduce a waveform similar to GW190521.
By contrast, for the type II waveform, i.e., for the case of less compact and less massive disks, the waveforms do not resemble that of GW190521: For these models, the amplitude of the initial burst waves is not much higher than the subsequent quasi-periodic waves. Also, for the type Is waveform, i.e., for high-mass and compact disk models such as M12, H14, E15, and E16, the wave cycle of the high-amplitude burst waves is ∼ 1, which is not as many as that of GW190521. Figure 6 displays the dimensionless effective amplitude of gravitational waves for the selected models: The left panel compares the spectra for given values of the disk mass with different compactness and the right panel compares the spectra for approximately the same value of c 2 r out /(GM BH,0 )(∼ 30) with different disk mass. Here, Here the peak frequency is slightly lower, ∼ 40 Hz, for the models with more massive disks such as H14-H16 and E15-E17. Our interpretation for this is that the total mass inside the region for which the disk has the highest density, M tot , is appreciably larger than M BH,0 , and thus, the angular velocity of the one-armed spiral mode is slower (in inversely proportional to M tot ) for a give compactness of the disk.
The peak amplitude of h eff is higher for more compact disks (i.e., for smaller values of c 2 r out /(GM BH,0 )) with a given value of the disk mass and for more massive disks, as expected from Figs. 3 and 4. The noteworthy features found from Fig. 6 are as follows: Irrespective of the waveform type, (i) the effective amplitude increases with the frequency for f ≤ f peak and (ii) the effective amplitude steeply drops for the frequency beyond f peak . Although the feature (ii) is similar to those of binary black-hole mergers, the feature (i) is totally different from them. Thus, in the detection of gravitational waves with a high signal-to-noise ratio, it would be easy to distinguish the waveforms of the massive black hole-disk systems from those of binary black hole mergers. However, for the detection with a low signal-to-noise ratio (which is in particular the case for the low-frequency region of 30 Hz with the gravitational-wave detectors operated in the thirdobservational run [46]), such difference is unlikely to be confirmed. . This feature is indeed similar to that of GW190521 presented in Ref. [1] under the condition that the signal-to-noise ratio is not very high. We note that as shown in Fig. 6, for models H16 and E17 for which the disk mass is comparable to the black-hole mass, the prominent spot appears at a relatively low frequency with f ∼ 40 Hz while for other models, it is ∼ 50 Hz. On the other hand, for model M15 of the type II waveform, not one prominent spot but an extend bright region appears. This feature is universally found for models L12, J13, and K14, i.e., for the lowmass or less compact models, and is different from that of GW190521.
The signal-to-noise ratio of gravitational waves with respect to the designed sensitivity of advanced LIGO at the hypothetical distance of 100 Mpc is listed in Table II. Here, the signal-to-noise ratio is calculated for the detection by a single detector. Note that the numerical simulations were stopped at 0.3-0.5 s after the peak gravitational-wave amplitude is reached. Since longlasting quasi-periodic waves with the duration longer than 0.5 s is likely to be present, the signal-to-noise ratio shown here should be considered as the minimum value.
Taking into account that the sensitivity of advanced LIGO in the third observational run was by a factor of ∼ 2 lower than the designed sensitivity [46], the horizon distance to the source for gravitational waves from the massive black hole-disk systems would be ∼ 100-300 Mpc for advanced LIGO at the third observational run. Here, for the models with more compact and more massive disks, the horizon distance is larger. Thus, if the source of GW190521 was a massive black hole-disk system that results from a collapse of a rapidly-rotating very-massive star, it should occur in a fairly nearby galaxy. Then, we could have a question whether a metal-poor star, which can form a very massive star, could exist in the rela- tively late-time universe. However, currently, we do not have rich observational information on this, although this question may be answered by the future observation [47].
It is also worthy to note that in the future detectors currently proposed, such as Einstein Telescope [48] and Cosmic Explorer [49], the signal-to-noise ratio for the events studied in this paper will be more than 10 times as high as the detectors currently in operation. For such future detectors, the horizon distance will be deeper than 1 Gpc and the formation rate of the very-massive stars is likely to be higher than in the late-time universe. Considering that the typical frequency of gravitational waves is approximately proportional to the inverse of the blackhole mass and that the future detectors will have a high sensitivity down to 20 Hz [48,49], the detectability will be also high for black holes of higher mass of ∼ 100M . Thus, in the future, gravitational waves not only from binary black hole mergers but also from high-mass black hole-disk systems could be an interesting target for exploring the formation process and formation history of massive black holes of mass ∼ 50-100M .

IV. DISCUSSION
Just prior to the formation of the systems composed of a massive black hole and a massive disk which we considered in this paper, a stellar core collapse of very-massive star, triggered by the pair instability, should occur. Since the black hole is likely to be formed immediately after the collapse (e.g., Refs. [19,20]), electromagnetic counterparts associated with the usual core-collapse supernovae are unlikely to accompany with this. However, after the formation of a black hole, the material could be ejected from the massive disk surrounding the black hole for the long-term evolution, and such material injects the outgoing momentum and energy into the extended envelop of radius ∼ 10 14 cm surrounding the carbon-oxygen core [19]. For a sufficiently high energy injection, the envelop will explode as in the core-collapse supernovae, and the ejecta may be the source for an electromagnetic counterpart in the optical-infrared bands with a high luminosity. The brightness of the supernova-like explosion is likely to depend strongly on the kinetic energy, E kin , of the mass outflow from the disk [19]. Our previous work [29] shows that 10-20% of the disk material could be ejected by a viscous process with the typical velocity of 0.05c from massive disks surrounding central black holes. Assuming these values, we have where M eje and v eje are the rest mass and velocity of the ejecta. Thus for the ejecta mass of 1-10M , the injection energy is broadly considered to be 10 51 -10 52 erg. We note that the ejecta energy may be increased by subsequent nucleosynthesis.
Using the Arnett's model [50], Ref. [19] estimated that for E kin = 10 51 -10 52 erg and the envelop mass of ∼ 150M , the absolute luminosity of an explosion would be ∼ 10 42 -10 43 erg/s with the duration of O(1 yr). Thus, for a hypothetical distance of 100 Mpc, the apparent magnitude is ∼ 17-19 mag. Such a source is observable or excludable by the current optical-infrared telescopes for the transient electromagnetic counterpart search as indicated in Ref. [51]. Unfortunately, the electromagnetic counterparts of GW190521 were not seriously searched by the optical-infrared telescopes, because GW190521 was announced as a candidate for binary black holes.
The collapse of rapidly-rotating very-massive stars of initial mass larger than ∼ 200M would be rare and hence the typical distance to this collapse event would be large 10 Mpc. Thus, the possible electromagnetic counterparts of this kind of the event are not likely to be very luminous, and it might not be easy to detect it in the absence of the information of the sky localization. As we studied in this paper, gravitational waveforms from the collapse of rapidly-rotating very-massive stars could be similar to that of the merging binary black holes with a value of high chirp mass 60M . Thus an alert for such an event could be issued as a candidate for a binary black hole merger as in the case of GW190521 in the future. Although a careful follow-up observation by electromagnetic telescopes is usually absent for the candidates announced as binary black holes, in the future, it will be interesting to perform follow-up observations for the candidates of very high-mass binary black holes. The problem is that quite unfortunately the LIGO-Virgo collaboration never announce the chirp mass (i.e., typical frequency) of the candidate events of gravitational waves. Thus it will not be possible to distinguish the interesting events among a huge amount of the "binary black hole" candidates unless the policy of the LIGO-Virgo collaboration is changed in the future.
In this paper, we focus only on the excitation of an unstable mode in massive disks surrounding a high-mass black hole and associated gravitational waves emitted. Another interesting aspect of this system is that a large amount of matter could be ejected from the massive disk during the enhancement of the one-armed spiral density wave and probably through the subsequent long-term viscous evolution. For a reliable estimation of E kin mentioned above, it is important to perform a long-term simulation for exploring the mass ejection process [28,29] and to quantitatively explore the value of E kin . In the ejecta from the disk, nucleosynthesis could also occur and the radio-active energy could be released for the additional source of the electromagnetic counterparts. If an amount of Ni with mass larger than ∼ 1M is synthesized in the ejecta, the event may be similar to the luminous supernovae (e.g., Ref. [52]). The detailed exploration of the mass ejection and nucleosynthesis is an interesting topic for the subsequent work. It is also interesting to explore how much mass can be ejected from the system. If a significant fraction of the matter forms a massive disk and is subsequently ejected from the system, the final mass of the formed black hole can be smaller than the progenitor mass, and may be smaller than 120M , i.e., may come into the mass-gap range.
There is another possible phenomenon for emitting gravitational waves similar to GW190521. In this paper, we suppose that a black hole is formed during the collapse. If the centrifugal force in the central region of the collapsing very-massive stars is strong enough, not a black hole but a very massive spheroidal or toroidal object can be formed (see, e.g., Refs. [11,22] for other contexts). If such an outcome is very compact and dynamically unstable to one-armed spiral mode or bar-mode deformation, gravitational waves with a high amplitude, which would have the same order of magnitude as that from massive black hole-disk systems, can be emitted. Assuming the formation of a compact object, the estimated frequency of gravitational waves is where M obj is the mass of the spheroidal/toroidal object and R obj is the typical radius for the highest-density part of such an object at the onset of the dynamical instability, respectively. Thus, the formation of a compact spheroid/toroid from a rapidly rotating very-massive stellar core can be another scenario for GW190521.
Finally, we note that for forming a rapidly-rotating very-massive star, a metal-poor star would be necessary. Such a star is unlikely to be formed frequently in the late-time universe, and thus, it would be rare that the collapse of the rapidly-rotating very-massive star is observed at the distance of several hundred Mpc. However, for given uncertainty on the formation rate of such stars in metal-poor galaxies [47], it is important to impose the constraint using the observational results. From this perspective, the gravitational-wave observation together with the electromagnetic-counterpart search can provide us a valuable opportunity.
To summarize, in this paper, we have suggested that burst gravitational waves of frequency ∼ 50 Hz, which are observed with a small signal-to-noise ratio and are composed only of a few cycles of high-amplitude waves, can be interpreted not only by those of binary black holes but by those by other scenarios. We suspect that the source of GW190521 might not be a merger of binary black holes but a stellar collapse of a very massive star leading temporarily to a black hole of mass ∼ 50M and a massive disk of several tens of solar mass that is dynamically unstable to the one-armed spiral-shape deformation. The point of this scenario is that the progenitor should be heavy enough to form a high-mass black hole of mass ∼ 50M and a high-mass disk, which can emit gravitational waves of low frequency ∼ 50 Hz, although we do not need the fine-tuning of the disk mass. One should keep in mind that it is not easy to definitely conclude the source of this-type of burst events only from a gravitational-wave observation with a low signal-to-noise ratio. In the future, electromagnetic counterpart searches will play a key role for pining down the most probable scenario and an alert from gravitational-wave observation community suitable for such searches is obviously necessary.
In the present paper, we prepare an initial condition of a massive black hole and a massive disk supposing that it would be the outcome formed during the collapse of a rapidly-rotating very-massive stellar core. Obviously, this setting is idealized. For more realistic work, it is necessary to perform a simulation started from a rapidlyrotating very-massive progenitor star. We plan to perform such simulations in the subsequent work.