© ESO, 2011

1. Introduction

The mass of astrophysical objects – except maybe for stars – is generally difficult to determine with precision, mostly because of inappropriate tracers, relatively low spatial resolution, and a certain misunderstanding of the internal structure and physical processes involved. This is the case for giant disks orbiting supermassive black holes in active galactic nuclei (AGN). For some nearby objects, the cold gas rotating in the outermost regions (the subparsec scale typically from the center) is detected at radio wavelengths through water vapor emission (e.g. Miyoshi et al. 1995; Braatz et al. 2009). The inner regions, not accessible yet to current instruments, could host the bulk of the mass if the total surface density in the disk varies roughly with the cylindrical radius R as R^-2 or faster (Shakura & Sunyaev 1973; Collin-Souffrin & Dumont 1990). Estimating the disk mass is a complex task. It necessitates a global disk model capable of describing the dynamics of gas, its thermodynamics, its chemical complexity, as well as its interaction with radiation (lines and continuum). The disk mass is an important quantity in understanding the AGN phenomenon. Along with the accretion rate, turbulent viscosity, and black hole mass, it helps to put constraints on the activity of the AGN in terms of stability, lifetime, luminosity, and matter supplied from the host galaxy (Combes 2001; Collin & Zahn 2008).

The mass of disks can also be probed via the consequences of their gravity. All the material contained within a disk exerts gravitational forces on itself – the so-called “self-gravity” – which influences or even strongly governs (like in galaxies) internal orbital motions, sometimes up to instability (e.g. Mestel 1963; Binney & Tremaine 1987; Papaloizou & Lin 1995). Even in the presence of a massive central object, self-gravity may cause a slight deviation in Kepler’s law, which is interesting to analyze and to quantify. Obviously, non-Keplerian rotation can have other origins like pressure effects such as in slim/thick disks (Abramowicz et al. 1988) or magnetic fields (Heyvaerts & Priest 1989). Here, we focus on self-gravity, which is expected to play a role in geometrically thin disks (Shore & White 1982; Shlosman & Begelman 1987).

There are many articles that aims to establish the relation between the dynamics (through gravitational potentials) and mass density distribution, especially in the context of galactic dynamics (Binney & Tremaine 1987). Existing potential/density pairs do not however seem fully appropriate to gaseous disks surrounding a central object, probably because star/disk systems where the gas exhibit non-Keplerian motions are still marginal. As models and theories suggest (e.g. Shakura & Sunyaev 1973; Pringle 1981; Collin-Souffrin & Dumont 1990; Huré 1998), gaseous disks in AGN are large and are expected to exhibit a self-similar behavior over some radial range. Huré et al. (2008) determined an accurate formula for the gravitational potential in the midplane of a flat power-law disk with finite size and mass. It is valid for a wide range of the power index for surface density. In this article, we use this result to derive an algebraic relation between the orbital velocity of the gas, the disk parameters (surface density profile, mass, size), and the mass of the central object, assuming a pressure-less disk at centrifugal equilibrium, as commonly done. This relation furnishes a simple method for determining how the mass is shared between the disk and the central object. As expected, the “modified” rotation law is not, in this model, a power law of the radius as often considered in this context (Herrnstein et al. 2005). Although this study is valid for any kind of astrophysical star/disk system (like in circumstellar environments), we focus on AGN disks whose kinematics have been observed in VLBI through water-maser emission.

This article is organized as follows. In Sect. 2, we recall the model of a pressure-less disk at centrifugal equilibrium surrounding a central black hole. We introduce the formula for the gravitational potential in the disk midplane by Huré et al. (2008), and derive the general expression for the velocity of the orbiting gas as a function of the disk mass, black hole mass, and surface density profile through the power-law index. In Sect. 3, we show how this expression (or its zero-order version) can be used to estimate how the mass is shared between a central object and its surrounding disk by analyzing observational data in a “position-dynamical mass” diagram (instead of the classical position-velocity diagram). We first applied the method to IC 1481, thereby refining the disk parameters reported in Mamyoda et al. (2009). We discuss uncertainties in Sect. 4, and show how dispersion naturally goes against the method. Section 5 is devoted to applying the method to a sample of a few well known AGN hosting a masing outer disk. We conclude in Sect. 6.

2. The basic model

2.1. A pressure-less disk at centrifugal equilibrium

We consider a gaseous disk with inner edge a_in and outer edge a_out ≫ a_in, orbiting a central black hole with mass M_BH. This disk is assumed to be axially symmetrical, flat (i.e. no vertical thickness), pressure-less, and steady. At centrifugal equilibrium, the rotation velocity v of material at cylindrical distance R in the midplane of the disk, in the reference frame of the black hole, is given by the standard relation: (1)where ψ_d is the gravitational potential of the disk. This latter function critically depends on the surface density profile Σ(R) through the Poisson integral. It is generally not easily accessed by analytical means, even in the actual one-dimensional case.

There is a broad literature devoted to determining potential-density pairs (ψ_d,Σ) for axially symmetric systems (e.g. Binney & Tremaine 1987; Evans & Collett 1993). Here, we consider the class of flat, power-law distributions where the surface density varies according to (2)where R = a_outϖ, Δ = a_in/a_out is the axis ratio, Σ_out the surface density at the outer edge, and s is a constant. Such a profile seems well-suited for large, gaseous disks in AGN, at least in the framework of geometrically thin disk models that predict s ≈  − 1 typically (Shakura & Sunyaev 1973; Pringle 1981; Collin-Souffrin & Dumont 1990; Huré 1998). Potential-density pairs for flat power-law disks, including Mestel’s solution (s =  −1), are summarized in Evans & Read (1998). These correspond to infinite disks (i.e. a_in = 0 and a_out → ∞) whose mass is infinite as soon as s >  −2. Conway (2000) has produced formal solutions corresponding to finite disks and no inner edge (i.e. a_in = 0), but for even positive indexes (i.e. s = 0,2,...). Unfortunately, accounting for edges increases the mathematical difficulties. Huré et al. (2008) have recently produced a reliable approximation for ψ_d associated with Eq. (2), namely (see their Eq. (53)) (3)where¹ψ_out = 2πGΣ_outa_out and (4)where C is the Catalan’s constant. This approximation is accurate within a few percent, provided the disk is large enough (Δ ≪ 1) and  − 3 ≲ s ≲ 0. These conditions are probably met in most astrophysical disks, especially in AGN disks (Shakura & Sunyaev 1973; Pringle 1981; Collin-Souffrin & Dumont 1990; Huré 1998), and others (e.g. Dubrulle 1992). Actually, we have Δ ~ 3R_S/a_out ≈ 10^-5 for a 10⁸   M_⊙ AGN black hole accreting a parsec size disk (R_S being the Schwarzschild radius). Besides, Eq. (3) does not suffer from the edge singularities expected when considering sharp edges.

2.2. The dynamical mass

As Eq. (1) shows, a “good” variable to measure masses in this model is (e.g. Yamauchi et al. 2004): (5)and the disk makes its own contribution through ψ_d. In the case of a spherical distribution (where R → r), this dynamical mass would represent the enclosed mass at a given radius (Herrnstein et al. 2005). Here, things are quite different since matter is gathered in a plane. There is no obvious use of the Gauss theorem, although the monopole approximation gives the right order of magnitude. Inserting in Eq. (5) the velocity v given by Eq. (1), and using Eq. (3) for the disk potential ψ_d, we finally get (6)where M_d is the disk mass (at a_out), and (7)is a second-order polynomial in s. It is displayed in Fig. 1. In the range of interest, P₂(s) is always positive and has a limited range of variation, since 0.14 ≲ P₂(s) ≲ 1.11. We also have P₂(−2) = P₂(−1) = 1 and a mean value of  ~0.8.

Fig. 1 P2(s) in the range of validity of Eq. (3).

From Eq. (6), we expect two extreme behaviors of the function μ(ϖ) depending on the disk properties. If s ≤  −2 (the case of “centrally peaked” distributions), we have (8)In the absence of radial gradient of μ, it is not possible to separate the disk and the central black hole, regardless of the disk mass. A Keplerian rotation curve result from a massive disk (without central object) or to a massive central object (with a light disk, as often considered). If s >  −2 (the bulk of the disk mass stands in the outermost regions), we have (9)which is essentially an increasing function of the radius (positive gradient). Disk rotation is therefore always super-Keplerian. The larger the disk mass, the larger the deviation from Kepler’s law.

We conclude that, if the disk mass is significant with respect to the central mass, v²(R) is the sum of two power laws of the radius, and this sum is not a power law (see below). In other words, fitting the global rotation curve of a system containing a massive disk and a black hole with a single power law (e.g. v ∝ R^γ), cannot give any quantitative information about the mass distribution in the framework of Newtonian gravity. This approach is often considered when a massive disk is suspected (see e.g. Herrnstein et al. 2005; Kondratko et al. 2008; McCallum et al. 2009). See, however, Appendix A for a short discussion of the velocity index γ.

For s =  −1, Eq. (6) reads (10)which is to be compared to the case of a Mestel disk (Mestel 1963) (11)where is the cumulative disk mass at the actual radius (linear with the radius ϖ). As already pointed out elsewhere (Binney & Tremaine 1987), this equation “happens to give the same answer” as what is deduced from the Gauss theorem for a spherically symmetric distribution. In contrast, Eq. (10) corresponds to a Mestel disk truncated on both sides. It includes edge effects and total disk mass, and explicitly contains the central point mass (which is not part of Mestel’s disk model). In the limit Δ → 0, we then recover Mestel’s solution.

3. AGN disk/black hole mass determination method. The case of IC 1481

3.1. Position-dynamical mass diagram

We immediately see from above that, if the rotation curve of the disk is partly known in the form of N observational points  { (R_i,v_i) } , then some constraints can be set on the disk mass, black hole mass, and surface density profile by fitting the data  { (R_i,μ_i) }  through Eqs. (8) or (9). Obviously, this procedure does not guarantee that the triplet (M_BH,M_d,s) is physically meaningful given the simplicity of our model and assumptions. Uncertainties in data also fragilize the inversion. Thus, there are three different possibilities.

• A:

  Data points  { (R_i,μ_i) } _N show no noticeable variation around a constant value, only a certain dispersion. The systems thus appears in Keplerian rotation. We deduce that either there is a light disk surrounding a massive black hole or the disk is rather massive but the gas is distributed such that s ≤  −2. The diagram only gives the quantity M_BH + M_d, which can be identified with the so-called “binding mass” M_b (or enclosed mass). There is no way to separate the black hole and the disk in this analysis.

• B:

  Data points show a significant variation, still with a certain dispersion (see Sect. 4). The gas rotates faster than Keplerian. If μ increases faster than ϖ, then s >  −1, otherwise s <  −1. In either case, fitting the data points through Eq. (9) can yield a triplet (M_BH, M_d, s).

• C:

  Data cannot be fitted by Eq. (9), or inferred parameters are non physical. In this case, our model is inappropriate. Various reasons can be invoked (see Sect. 6).

3.2. Zeroth-order: disk mass and central mass

The zeroth-order treatment of the non linear formula is interesting and instructive because it gives the orders of magnitude. Actually, if we consider that astrophysical disks are characterized by s ≈  − 1, we can expand Eq. (9) around s =  −1. We find (see also Eq. (10) with Δ → 0) (12)We conclude that, if observational data plotted on a position-dynamical mass diagram are almost linearly distributed in a position-dynamical mass diagram, then the slope is the disk mass M_d and the intercept is the black hole mass. In the following, we discuss both approaches in the context of AGN disks whose rotation, for some of them, is known from maser emission.

Fig. 2 Position-velocity data of maser spots in the disk of IC 1481 (top), and corresponding position-dynamical mass diagram (bottom). A non linear fit of the data (green curve) gives the disk mass, the black hole mass and the s-parameter. The linear fit (black line) is also shown. Data are from Mamyoda et al. (2009).

3.3. Scale-free formula. Scaling to AGN disks

Equation (9) is totally scale-free if we divide μ by the black hole mass. Actually, if q = M_d/M_BH denotes the mass ratio in the system, Eq. (9) writes (13)so that the model applies to any kind of star/disk system. In the context of AGN, black holes are supermassive (several million solar masses), the central accretion disk typically has the parsec size, and rotational velocities are, at such distances, hundreds of km s^-1. With this scaling, the formula for the dynamical mass becomes² where v_i is, following our model, the disk rotational velocity measured relative to the systemic/receding velocity of the system, and R_i relative to the rotation axis of the disk. In the following, masses are expressed in units of 10⁷   M_⊙, and denoted . Then, Eq. (9) becomes (17)and its linear version is (18)

3.4. An example: the case of IC 1481

The top panel in Fig. 2 shows the rotational velocity of maser spots observed in the outer disk of the active nucleus in galaxy IC 1481 versus the distance from the center (Mamyoda et al. 2009). The radius of the outermost maser spot is set to a_out (although the gas disk can extend farther away³). The bottom panel displays the same data once converted into a (ϖ,μ)-diagram, as well as the result of the fit by Eq. (17) (non-linear) and by Eq. (18) (linear). For each point, the square of the radial deviation between the data and the linear fit is shown (see below). The parameters of the two fits are gathered in Table 1. These are in good agreement, especially because the solution of the non linear fit gives an s-parameter close to  − 1. In such a case, the curvature is difficult to detect by eye. As announced in Mamyoda et al. (2009)⁴, the disk mass is higher than the black hole mass by a factor of 2 − 3.

4. Uncertainties and data dispersion

Position-velocity data deduced from observations generally suffer from dispersion and uncertainties, which can have different origins: physical (e.g. variability, non uniform dynamics, geometry and deprojection) and instrumental (i.e. lack of resolution). In particular, locating the position of emitters precisely is a critical point (e.g. Uscanga et al. 2007). In order to check whether the method is “robust” enough to infer some reliable information about masses, we generated a sample of N points , obeying exactly Eq. (17) for a given reference triplet . To these synthetic data, we added uncertainties on the radius (without presuming their origin) with three different levels of “noise”. We considered a Gaussian noise, with various standard deviations σ. The perturbed radii ϖ_i were then all rescaled so that finally ϖ_N = 1. Unfortunately, this rescaling procedure introduce a slight bias by tending to overestimate both the black hole mass mass and the disk mass.

Fig. 3 Synthetic data set made with N = 20 points generated from Eq. (17) with (top) and with (bottom). Uncertainties corresponding to a Gaussian noise with a standard deviation σ = 0.03 have been added to all the radii. Data were rescaled in radius such that the outer point has a normalized radius ϖN = 1.

Here, we report a typical experiment obtained for

• a massive disk configuration with corresponding to mass ratio q = 4;

• a massive black hole configuration with corresponding to ;

with s^ref =  −1.5 and N = 20 points in both cases. Initially, data were randomly spread over the range  [0.5,1 ] , and we considered 3 levels of noise: σ =  {0.03,0.05,0.1} . The synthetic -diagram obtained for σ = 0.03, which is presumed to mimic dispersed observationnal data, is shown in Fig. 3 with the exact rotation law given by Eq. (17). We tried to determine the best triplet in three ways:

• by the linear approach;

• by the non linear approach, but forcing s = s^ref;

• by the non linear approach.

The second case is only illustrative, as it can not be used in practice since s is not known a priori. For each configuration, the relative distance between the reference triplet and the one deduced by fitting the synthetic noisy data is measured by the parameter δ, with (19)As a result, δ → 0 when the method works, whereas δ ≳ 1 when it fails (i.e. when the output and the input triplets differ significantly⁵). Practically, the linear approach is performed from a least-square procedure. The results are gathered in Table 2. As expected, the lower the dispersion, the better the method. For the massive disk configuration and σ = 0.03 (low dispersion), we recover the disk mass within a few percent by the non linear fit, while the error in the black hole mass is about 30%. By the linear approach, we get the disk mass within 30%, while there is a factor 2 − 3 for the black hole mass. This is not all that surprising since the curvature of the function is noticeable in this case, so the error on the intercept can be large. As the dispersion increases, the linear approach becomes more and more reliable for getting the disk mass (within a few percent), but the black hole mass is still poorly determined (a factor 2 typically). Regarding the massive black hole configuration and σ = 0.03 (low dispersion), the black hole mass and the disk mass are deduced correctly within a few percent typically, especially using the non linear method. As the dispersion increases, the uncertainty in both two quantities rises, for the same reason as mentioned hereabove; however, the black hole mass is determined with about 15%.

Fig. 4 Position-dynamical mass diagram for each AGN listed in Table 3 (left) and deviation in radius with the linear fitting (right). Observational data are circles, and the red dashed line is the averaged value of the considered. The linear fit through Eq. (18) is the bold line (see Table 4).

Since the errors on the radius propagate to the variable μ, we made the same fit in the (R,v)-diagram and found very similar results. In brief, we see that when the “global slope” is steep (massive disk configuration), the disk mass is determined with precision. Conversely, when the “global slope” is almost zero (massive black hole configuration), the black hole mass (value at the intercept) is well determined. When dispersion is too large, the correlation becomes too weak for interpretation, and the method fails to give the triplet, but it is important to note that the linear approach always gives the value of the most massive component within a factor less than about 15%, and the less massive by a factor 2 − 3.

5. Application to a few AGN masing disks

In the case where position-velocity diagrams are available in the literature, we analyzed a few AGN masing disks to try to put constraints on masses. This is not a new problem regarding the mass of black hole. The mass of disks is, however, much a subject of debate, since there is almost no systematic or generic study like the one presented here, except isolated attempts (see e.g. Huré 2002; Lodato & Bertin 2003; Herrnstein et al. 2005; Kondratko et al. 2008; McCallum et al. 2009). Table 3 lists the systems considered here. Position-velocity data points were obtained by digitalizing graphs when published, without special treatment (in particular, we fully trust in the analysis by the authors to furnish deprojected position-velocity data in the reference frame of the black hole). The table also contains masses previously proposed for the black hole and disk, mostly through models of spherical distributions. In each case listed in Table 3, we have tried to fit the data by the linear formula and by the non linear formula as well. Position-dynamical mass diagrams, together with the linear fit and averaged value, are displayed in Fig. 4. The results are given in Table 4. In particular, from the linear regression, we estimated the degree of dispersion by computing the standard deviation σ (Col. 8). From the simple analysis performed in Sect. 4, we immediately expect no precise values for NGC 4945 and Circinus. For NGC 4258 and NGC 1068, the dispersion is much weaker with σ ≈ 0.2−0.3. For these two objects, the linear approach should give the mass of the most massive component quite correctly (i.e. see Table 2). Finally, for IC 1481, UGC 3789 and NGC 3393, dispersion is less than  ~10%, so both mass components should be obtained with “accuracy”.

It was not possible to derive a triplet (M_BH, M_d, s) for all these systems through the non linear approach (negative mass or non convergence of the fitting procedure), even for UGC 3789 and NGC 3393 whose data, like for IC 1481, are spread over a large radial domain and little dispersed. In two cases (NGC 4945 and Circinus), the method is inappropriate: no correlation really exists in the position-dynamical mass and inferred parameters are non physical (in particular, M_d < 0 for NGC 4945). We are aware that some disks have been predicted to be geometrically thick, so the model of flat disk used here should naturally fail. In contrast, results obtained through the linear approach are satisfactory for UGC 3786, NGC 3393, and NGC 1068 where disk mass is similar to, or higher than the mass of the central black hole. For these objects, the massive disk should be prone to instabilities. For NGC 4258, the slope is very weak and a mass ratio of about 0.05 is expected.

6. Concluding remarks

In this paper, we have reported a simple method to estimate the mass of a disk surrounding a central black hole and assuming a flat, pressure-less disk at centrifugal balance. This model obviously applies to other kinds of systems containing a disk and (possibly) a central star. Instead of the conventional position-velocity diagram, the position-dynamical mass diagram makes the measurement of possible deviations to Kepler’s law easier. Using recent calculations for the gravitational potential of “truncated self-similar” disks, we have shown that i) the rotation law is generally not a power law of the radius, and ii) such deviations are, at zero order at least, directly related to the cumulative disk mass.

It is clear that the disk model is very simple and can be improved in several ways. At the same time, it is difficult to make the method robust and universal, since position-velocity data generally show a certain dispersion that can have several origins like thickness effects, instabilities (warps), etc. Through a simple analysis of uncertainties, we have shown that the linear approach gives quite correctly the most massive component (typically with a few tens of percent) if data dispersion is large (standard deviation larger than 0.1). For weak dispersion, both the disk mass and the black hole mass should be accessible through non linear data fitting, as well as the surface density power index. Moreover, the systemic velocity plays a major role, and the asymmetry often observed between the redshifted part and the blueshifted part of the rotation curve must be accounted for. This method must therefore be seen as a first step in the analysis of masses in star/disk systems based on gravity.

Acknowledgments

We thank S. Collin and F. Herpin for useful comments. The referees are acknowledged for their comments and suggestions, in particular about the problem of uncertainties.

References

Appendix A: Note on the power index of the rotation curve

Let γ be the local index of the rotation velocity, namely (A.1)where M₀ is a reference mass. In the present model, it is not a constant, but a function of the radius. We easily find from Eq. (5) (A.2)with for a Keplerian motion associated with a central point mass M₀. From Eq. (9), we also deduce (A.3)and then γ can be written in the form: (A.4)where (A.5)is the deviation to the Keplerian index. For s =  −2, we have δγ_kep. = 0: the rotation curve resembles a Keplerian curve due to a point mass (see above). This result was already known. The mean deviation is (A.6)where q = M_d/M_BH is the disk-to-black hole mass ratio. The above expression can be integrated exactly for some values of the s-parameter. For instance, with s =  −1, we find . It means that a deviation to Kepler’s law as low as 0.02 on the velocity index (i.e. γ =  −0.48) could imply a mass ratio q = 10%.

Fig. A.1 Relation between the mean index of the super-Keplerian rotation law and the disk mass to black hole mass ratio for some values of the s-parameter. To a given mean index  ⟨ γ ⟩ , it corresponds a minimum mass ratio q (red). To a given mass ratio, it corresponds a maximal mean index (blue).

This value is, in magnitude, consistent (within a factor  ~ 2) with the monopole approximation that predicts a mass ratio  ~ 0.4. For systems containing a disk significantly less massive than the central object (i.e. q ≲ 1) and Δ ≪ 1, we find (A.7)Figure A.1 displays the mean velocity index  ⟨ γ ⟩  versus q for a few values of the s-parameter. This plot enables either to predict

the maximum departure to Kepler’s law for a given mass or to bound the disk mass given a mean velocity index. Practically, if deviations remain of small amplitude for low-mass disks, these may not be exploitable (data dispersion, thickness effects).

All Tables

All Figures

	Fig. 1 P2(s) in the range of validity of Eq. (3).
In the text

	Fig. 2 Position-velocity data of maser spots in the disk of IC 1481 (top), and corresponding position-dynamical mass diagram (bottom). A non linear fit of the data (green curve) gives the disk mass, the black hole mass and the s-parameter. The linear fit (black line) is also shown. Data are from Mamyoda et al. (2009).
In the text

	Fig. 3 Synthetic data set made with N = 20 points generated from Eq. (17) with (top) and with (bottom). Uncertainties corresponding to a Gaussian noise with a standard deviation σ = 0.03 have been added to all the radii. Data were rescaled in radius such that the outer point has a normalized radius ϖN = 1.
In the text

	Fig. 4 Position-dynamical mass diagram for each AGN listed in Table 3 (left) and deviation in radius with the linear fitting (right). Observational data are circles, and the red dashed line is the averaged value of the considered. The linear fit through Eq. (18) is the bold line (see Table 4).
In the text

	Fig. A.1 Relation between the mean index of the super-Keplerian rotation law and the disk mass to black hole mass ratio for some values of the s-parameter. To a given mean index  ⟨ γ ⟩ , it corresponds a minimum mass ratio q (red). To a given mass ratio, it corresponds a maximal mean index (blue).
In the text