Finding apparent horizons in numerical relativity.

Abstract

This paper presents a detailed discussion of the ``Newton's method'' algorithm for finding apparent horizons in 3+1 numerical relativity. We describe a method for computing the Jacobian matrix of the finite differenced $H(h)$ function by symbolically differentiating the finite difference equations, giving the Jacobian elements directly in terms of the finite difference molecule coefficients used in computing $H(h)$. Assuming the finite differencing scheme commutes with linearization, we show how the Jacobian elements may be computed by first linearizing the continuum $H(h)$ equations, then finite differencing the linearized (continuum) equations. We find this symbolic differentiation method of computing the $H(h)$ Jacobian to be {\em much} more efficient than the usual numerical perturbation method, and also much easier to implement than is commonly thought. When solving the discrete $H(h) = 0$ equations, we find that Newton's method generally converges very rapidly. However, if the initial guess for the horizon position contains significant high-spatial-frequency error components, Newton's method has a small (poor) radius of convergence. This is {\em not} an artifact of insufficient resolution in the finite difference grid; rather, it appears to be caused by a strong nonlinearity in the continuum $H(h)$ function for high-spatial-frequency error components in $h$. Robust variants of Newton's method can boost the radius of convergence by O(1) factors, but the underlying nonlinearity remains, and appears to worsen rapidly with increasing initial-guess-error spatial frequency. Using 4th~order finite differencing, we find typical accuracies for computed horizon positions in the $10^{-5}$ range for $\Delta\theta = \frac{\pi/2}{50}$.