Selfgravity at the scale of the polar cell
J.M. Huré^{1,2}  A. Pierens^{3}  F. Hersant^{1,2}
1  Université de Bordeaux, Observatoire Aquitain des Sciences de l'Univers, 2 rue de l'Observatoire, BP 89, 33271 Floirac cedex, France
2 
CNRS/INSU/UMR 5804/LAB, BP 89, 33271 Floirac Cedex, France
3 
LALIMCCE/USTL, 1 impasse de l'Observatoire, 59000 Lille, France
Received 6 February 2009 / Accepted 5 March 2009
Abstract
We present the exact calculus of the gravitational potential and acceleration along the symmetry axis of a plane, homogeneous, polar cell as a function of mean radius ,
radial extension ,
and opening angle
.
Accurate approximations are derived in the limit of high numerical resolution at the geometrical mean
of the inner and outer radii (a keyposition in current FFTbased Poisson solvers). Our results are the full extension of the approximate formula given in the textbook of Binney & Tremaine to all resolutions. We also clarify definitely the question about the existence (or not) of selfforces in polar cells. We find that there is always a selfforce at radius
except if the shape factor
,
asymptotically. Such cells are therefore well suited to build a polar mesh for high resolution simulations of selfgravitating media in two dimensions. A byproduct of this study is a newly discovered indefinite integral involving complete elliptic integral of the first kind over modulus.
Key words: accretion, accretion disks  gravitation  methods: analytical  methods: numerical
1 Introduction
The structure and dynamical evolution of astrophysical discs is governed mainly by gravity. In this context, the numerical computation of the gravitational potential and forces of discs is of fundamental importance and the improvement in accuracy, resolution and computing time remains an interesting challenge (e.g., Li et al. 2008; Londrillo 2004; Matsumoto & Hanawa 2003; Huré 2005; Mueller & Steinmetz 1995; Jusélius & Sundholm 2007). In many (if not all) hydrodynamical simulations of flat, selfgravitating discs (e.g., Baruteau & Masset 2008; Huber & Pfenniger 2001; Zhang et al. 2008), the potential is derived everywhere in the disc by means of Fast Fourier Transforms (e.g., Binney & Tremaine 1987). This is an efficient technique, which can be extended to deduce directly accelerations (Baruteau & Masset 2008). A weak point however is in the determination of the selfgravitating component (i.e., the effect of a cell on itself), which involves an improper integral. An approximation of this selfpotential was proposed in Binney & Tremaine (1987), but for a cell in which the surface density varies with the polar radius R as R^{3/2}. In the same conditions, Baruteau & Masset (2008) argued that the radial acceleration is expected to vanish at the center of each computational cell. This assertion is certainly correct asymptotically for homogeneous cells as their size becomes smaller and smaller. We stress that because Newton's law goes like R^{2}, the selfgravitating component generally provides a significant contribution to the total potential or force, and must therefore be treated as properly as possible, especially at high resolution. Any error, even relatively small, may introduce artifacts or biases in models.
In this paper, we report the exact expressions for both the potential and radial acceleration due to a homogeneous polar cell inside the cell itself and study the highresolution limit. We recall in Sect. 2 the integral expression for the potential along the axis of a polar cell which can be decomposed into two parts: a potential created by an homogeneous disc and the other a potential due to horseshoelike surface. Sections 3 and 4 are devoted to the determinatons of these two contibutions. The results for the polar cell (potential and acceleration) are obtained in Sect. 5. We discuss the case of high resolutions in Sect. 6 and show that there is always a selfacceleration except if the cell has a special shape. These results are important since they contribute to the reliability and accuracy of current simulations of selfgravitating media on a polar mesh.
2 The polar mesh
We consider a planar, homogeneous, polar cell with inner polar radius a_{1}, outer radius
,
lower azimuth ,
and upper azimuth
as shown in Fig. 1. At any position
in the plane of this cell, the gravitational potential is given by Newton's generalized formula:
(1) 
where is the surface density. This cell has a single axis of symmetry defined by . If we introduce the quantity:
as well as the new variable such that , the potential along its axis of symmetry is given by:
(3) 
where . The integral over can be decomposed into two integrals with 0 as lower bound. We then have:
where
(5) 
is the incomplete elliptic integral of the first kind, m is the modulus, is the amplitude, and is the complete elliptic integral of the first kind. The interpretation of this relation is presented in Fig. 2 and follows from the superposition principle: the first part of the integral of Eq. (4) is the potential due to a flat disc, namely:
and the second part is that of a flat, horseshoelike surface:
Figure 1: A polar cell with inner radius a_{1}, outer radius a_{2}, radial extension , and opening angle . 

Open with DEXTER 
Figure 2: A polar cell is the superposition of a disc and a horseshoe with negative surface density. 

Open with DEXTER 
and
where is the complementary modulus, +m' is for (i.e., matter is located right to the position R) and m' is for (i.e., matter is located to the left). From Eqs. (2) and (9), we deduce that the selfgravitating potential in the cell along its symmetry axis is:
We note that this expression is rarely seen in the literature.
3 The homogeneous disc
3.1 Classical derivation
To calculate
,
one classically changes^{} the modulus m of
by setting
where .
The new modulus x is either
or
depending on R.
So, the potential due to an homogeneous disc is given by:
where region (I) is for , region (II) is for , and region (III) is for (see Fig. 2). Since the indefinite integrals in Eqs. (12) are known^{}, a closeform expression for can finally be deduced. It is:
3.2 A new indefinite integral?
From a numerical point of view, it is preferable to use u and v as modulus of the complete elliptic integrals rather than m (see Sect. 6). However, an interesting issue is raised if we perform a change of modulus in Eqs. (14), to restore m as the modulus of
and .
We find:
If we now compare these expressions with the first part of Eq. (10), we conclude that:
for any modulus k and complementary modulus k'. This indefinite integral is probably new (Jeffrey 2009).
4 The homogeneous horseshoe
To calculate
,
we have attempted to proceed as for the disc, by applying a change of modulus. We failed, mainly because, as mentioned, the new
amplitude
of F resulting from Eq. (11) now becomes a function of the modulus. This introduces a severe mathematical difficulty that we were unable to circumvent. The answer is however obtained by considering
Eq. (16) with incomplete elliptic integrals. Using the partial derivative^{} of F and E with respect to their modulus (keeping the amplitude constant), we find that:
where
and . In fact, can be expressed exactly in terms of basic functions. Actually, after some algebra, we have:
It follows from Eqs. (10), (18), (19), and (20) that the potential along the axis of the homogeneous horseshoe is given by:
As for Eq. (14), the expression for region (II) can be simplified.
Figure 3: The potential along the symmetry axis of a homogeneous polar cell with geometrical parameters a_{1}=0.95, and ( and ). The minimum lies at . Also shown are the potential of the associated disc and horseshoe. 

Open with DEXTER 
5 Potential and acceleration in the polar cell
The potential inside a polar mesh and along its symmetry axis is found from Eqs. (15b) and (21b). It is:
where
with .
This expression is exact. Figure 3 displays
,
and their difference
versus R for a typical cell. By deriving this formula with respect to R, we obtain the relation:
which is the opposite of the gravitational acceleration g_{R} along the cell axis. If we set and where
is the potential due to the cell at the origin R=0, and
then Eq. (23) takes the form:
This is the generalization of the Ordinary Differential Equation (ODE) derived by Huré & Hersant (2007) to a polar cell with opening angle (which becomes a disc for ). This ODE can be easily solved numerically using standard schemes since is known both at the center (see Eq. (24)) and at infinity where . The full generalization to polar cells where the surface density is a powerlaw of the radius, i.e., , seems straightforward.
6 Approximations at high numerical resolutions and the ``magic'' shape factor
At high resolution, polar cells tends to become segments, squares, or arcs depending on the shape factor
of the cell defined by:
(27) 
where is the center of physical cells. When both (i.e. ) and , we can expand , , , and over m' (e.g., Gradshteyn & Ryzhik 1965; Carlson & Gustafson 1985) to obtain an approximation for and g_{R} as a function of R and the cell geometrical parameters . A radius plays a key role in current FFTbased Poisson solvers (see references in Sect. 1): this is the geometrical mean of the inner and outer radii, namely . This radius is the center of computational cells when a radial logarithmic mapping is applied to the entire physical mesh. For , the two modulus m_{1} and m_{2} are equal (and so are the complementary), and we have to first order. After some algebra, we find:
and
Figure 4 compares exact values of and given by Eqs. (22) and (23) with approximate values obtained from Eqs. (28) and (29) respectively, for three factors . We note that, for , relative errors increase. This is due to the ``low dynamic'' of the variable m as it approaches unity (and ). In this case, it is advisable to use u and v as the moduli of elliptic integrals instead of m to restore the accuracy in Eq. (22).
Figure 4: Top panels: exact potential and acceleration in a polar cell at radius for three shape factors ( circles) as well as their approximations ( lines) by Eqs. (28) and (29). Bottom panels: decimal logarithm of the relative error. Also displayed are the results obtained with the ``magic'' value . 

Open with DEXTER 
The potential always has its minimum inside the cell at a certain radius reagardless of the opening angle . This radius is found by equating the righthand side of Eq. (23) to zero. It depends on the parameters , , and . We expect for large opening angles (or low azimuthal resolutions) and for small ones (high azimuthal resolutions). When , as the polar cell becomes a rectangle. Since , there is in general a selfforce at at high azimuthal resolution. However, the gravitational acceleration exceptionally cancels out if the shape factor has the remarkable value . This is the single root of Eq. (29). For , we have . Thus, if the computational mesh contains cells all with shape factor , then there is no need to compute the gravitational force at , since it is zero in the high resolution limit.
More generally, since the two moduli m_{1} and m_{2} appearing in
Eq. (23) are functions of the ratio a_{2}/a_{1}, it is always possible to find a shape factor for which the acceleration vanishes at a fixed value of R/a_{1}. It can therefore be advantageous to build the polar mesh accordingly. For instance, let N_{R} and
be the numbers of cells in the radial and azimuthal directions, respectively, and
be the radius of the disc inner edge. Once the shape factor
has been numerically determined such that g_{R}=0 (this is
in the high resolution limit), the radius
of the outer edge is found from the relation:
(30) 
With such a mesh, there is no gravitational acceleration at the center of computational cells.
7 Concluding remarks
We have reported exact expressions for the selfgravitating potential and acceleration in homogeneous polar cells whatever their shape. The high resolution limit has been analyzed and reliable approximations have been derived. We confirm the validity of Binney & Tremaine (1987)'s prescription for the singular term K(0,0) to first order. We have shown that there is always a selfacceleration at the center of computational cells. However, this acceleration can cancel out by an appropriate choice of cell shape. This work^{} contributes to the treatment of selfgravity in hydrodynamical simulations, and especially for discs.
Acknowledgements
It is a pleasure to thank P. Barge, C. Baruteau and C. Surville for valuable discussions, as well as A. Jeffrey, coeditor of ``Tables of integrals, series and products'' by Gradshteyn & Ryzhik (1965).
References
 Baruteau, C., & Masset, F. 2008, ApJ, 678, 483 [NASA ADS] [CrossRef]
 Binney, J., & Tremaine, S. 1987, Galactic dynamics (Princeton, NJ, Princeton University Press) (In the text)
 Carlson, B., & Gustafson, J. L. 1985, SIAM J. Math. Anal., 16, 1072 [CrossRef]
 Gradshteyn, I. S., & Ryzhik, I. M. 1965, Table of integrals, series and products, ed. Yu. V. Geronimus, & M. Yu. Tseytlin 4th Ed. (New York: Academic Press)
 Huber, D., & Pfenniger, D. 2001, A&A, 374, 465 [NASA ADS] [CrossRef] [EDP Sciences]
 Huré, J.M. 2005, A&A, 434, 1 [NASA ADS] [CrossRef] [EDP Sciences]
 Huré, J.M., & Hersant, F. 2007, A&A, 467, 907 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
 Jeffrey, A. 2009, private communication (In the text)
 Jusélius, J., & Sundholm, D. 2007, J. Chem. Phys., 126, 094101 [NASA ADS] [CrossRef]
 Li, S., Buoni, M. J., & Li, H. 2008, ArXiv eprints
 Londrillo, P. 2004, Mem. Soc. Astron. Ital. Suppl., 4, 69
 Matsumoto, T., & Hanawa, T. 2003, ApJ, 583, 296 [NASA ADS] [CrossRef]
 Mueller, E., & Steinmetz, M. 1995, Comput. Phys. Commun., 89, 45 [NASA ADS] [CrossRef]
 Zhang, H., Yuan, C., Lin, D. N. C., & Yen, D. C. C. 2008, ApJ, 676, 639 [NASA ADS] [CrossRef]
Footnotes
 ... changes^{}
 The transformation is (e.g. Gradshteyn & Ryzhik 1965):
 ... known^{}
 In particular, we have for any
(Gradshteyn & Ryzhik 1965):
where is the complete elliptic integral of the second kind and is the complementary modulus.  ... derivative^{}
 In particular, we have (e.g. Gradshteyn & Ryzhik 1965):
where is the incomplete elliptic integral of the second kind.  ... work^{}
 A Fortran 90 package called PolarCELL is available at the following address: http://www.obs.ubordeaux1.fr/radio/JMHure/intro2polarcell.html
All Figures
Figure 1: A polar cell with inner radius a_{1}, outer radius a_{2}, radial extension , and opening angle . 

Open with DEXTER  
In the text 
Figure 2: A polar cell is the superposition of a disc and a horseshoe with negative surface density. 

Open with DEXTER  
In the text 
Figure 3: The potential along the symmetry axis of a homogeneous polar cell with geometrical parameters a_{1}=0.95, and ( and ). The minimum lies at . Also shown are the potential of the associated disc and horseshoe. 

Open with DEXTER  
In the text 
Figure 4: Top panels: exact potential and acceleration in a polar cell at radius for three shape factors ( circles) as well as their approximations ( lines) by Eqs. (28) and (29). Bottom panels: decimal logarithm of the relative error. Also displayed are the results obtained with the ``magic'' value . 

Open with DEXTER  
In the text 
Copyright ESO 2009