1 Introduction

TwoPunctures computes conformally flat initial data for two non-spinning or spinning, non-moving or moving black holes using the moving puncture approach [2]. The central task is to solve the Hamiltonian constraint for the conformal factor \(\psi \), which is decomposed as \begin {equation} \psi = 1 + \frac {m_+}{2 r_+} + \frac {m_-}{2 r_-} + u, \end {equation} where \(m_\pm \) are the bare masses of the two punctures located at \((\pm b, 0, 0)\) (or offset by center_offset), \(r_\pm \) are the coordinate distances from each puncture, and \(u\) is a smooth correction function that satisfies a nonlinear elliptic equation derived from the full Hamiltonian constraint. The Bowen-York extrinsic curvature [1] is used to incorporate linear momenta and spins.

The nonlinear elliptic equation for \(u\) is solved with a spectral method described in Ansorg, Brügmann & Tichy (2004) [3]. Newton iteration is used to handle the nonlinear terms. Once \(u\) is known, all ADMBase metric and extrinsic-curvature components are filled on the 3D Cactus grid by interpolation from the spectral solution.

2 Physical System

2.1 Conformal Decomposition

The spatial metric is assumed to be conformally flat: \begin {equation} \gamma _{ij} = \psi ^4 f_{ij}, \end {equation} where \(f_{ij}\) is the flat metric. The extrinsic curvature is decomposed into a trace-free part \(A_{ij}\) and trace \(K\): \begin {equation} K_{ij} = \psi ^{-2}\,A_{ij} + \tfrac {1}{3}\,\gamma _{ij}\,K. \end {equation} For the time-symmetric slice (\(A_{ij}=0\) and \(K=0\)) the Hamiltonian constraint reduces to a flat Laplace equation for \(\psi \). For the general Bowen-York case (\(K=0\) is still enforced on the slice but \(A_{ij}\neq 0\)), the Hamiltonian constraint becomes a nonlinear equation for \(u\).

2.2 Bowen-York Extrinsic Curvature

The Bowen-York extrinsic curvature for a single puncture with linear momentum \(P^i\) and spin \(S^i\) is \begin {equation} A_{ij}^{\rm BY} = \frac {3}{2r^2}\left [ P_i n_j + P_j n_i - (f_{ij} - n_i n_j) P^k n_k + \frac {2}{r}\left (\epsilon _{kli} S^k n^l n_j + \epsilon _{klj} S^k n^l n_i\right ) \right ], \end {equation} where \(n^i\) is the unit outward normal to a sphere of radius \(r\) centered on the puncture. The total extrinsic curvature is the sum of the contributions from both punctures.

2.3 Masses

The thorn distinguishes between bare masses \(m_\pm \) (the parameters \(\texttt {par\_m\_plus}\) and \(\texttt {par\_m\_minus}\)) and ADM masses \(M_\pm \) measured at the respective asymptotic infinities. When give_bare_mass = yes the user sets bare masses directly. When give_bare_mass = no the thorn iterates to find the bare masses that reproduce the requested target ADM masses target_M_plus and target_M_minus to within the tolerance adm_tol.

3 Numerical Implementation

3.1 Spectral Grid

The equation for \(u\) is solved on a spectral grid parameterized by coordinates \((A, B, \phi )\), where \(A \in [0,1]\) is a compactified radial coordinate covering both punctures and spatial infinity, \(B \in [-1,1]\) is an angular coordinate, and \(\phi \in [0, 2\pi )\) is the azimuthal angle. The number of spectral coefficients in each direction is set by npoints_A, npoints_B, and npoints_phi, respectively.

3.2 Newton Solver

Newton’s method is applied to the nonlinear equation for \(u\) with tolerance Newton_tol and at most Newton_maxit iterations. The linear system at each Newton step is solved by a banded direct solver. OpenMP parallelism is used within the Newton solver to accelerate the construction and application of the Jacobian.

3.3 Grid Interpolation

Once the spectral coefficients are known, the solution is transferred to the Cactus 3D grid in one of two ways controlled by grid_setup_method:

• Taylor expansion (default): the solution is Taylor-expanded about the nearest spectral collocation point. This is fast but may be slightly less accurate far from collocation points.

• evaluation: each grid point is evaluated using all spectral coefficients. More accurate but significantly slower for large grids.

The interpolation loop is parallelized with OpenMP.

3.4 Singularity Regularization

Two small parameters help avoid numerical issues near the punctures:

• TP_epsilon: replaces \(r_\pm \) with \((r_\pm ^4 + \epsilon ^4)^{1/4}\), smoothing the singularity.

• TP_Tiny: a floor applied to \(r_\pm \) after the above regularization.

• TP_Extend_Radius: inside this radius around each puncture the singular \(1/r\) term is replaced by a smooth polynomial extension, useful when the puncture lies inside the computational domain.

4 Using This Thorn

4.1 Obtaining This Thorn

TwoPunctures is part of the EinsteinInitialData arrangement, which can be checked out via the Einstein Toolkit release mechanism or directly from the repository. The thorn requires ADMBase, StaticConformal, and the Cactus flesh.

4.2 Basic Usage

Add EinsteinInitialData/TwoPunctures to your thorn list and set the following parameters in your parameter file:

ADMBase::initial_data  = "twopunctures"
ADMBase::initial_lapse = "twopunctures-antisymmetric"

TwoPunctures::par_b       = 3.0   # half separation (M)
TwoPunctures::par_m_plus  = 0.5
TwoPunctures::par_m_minus = 0.5
TwoPunctures::par_P_plus[1]  =  0.084
TwoPunctures::par_P_minus[1] = -0.084
TwoPunctures::npoints_A   = 30
TwoPunctures::npoints_B   = 30
TwoPunctures::npoints_phi = 16

The punctures are placed at \((+b, 0, 0)\) and \((-b, 0, 0)\) unless center_offset shifts the center of mass. To separate the holes along the \(z\)-axis instead, set swap_xz = yes.

4.3 Lapse Options

The initial lapse is selected via ADMBase::initial_lapse:

• twopunctures-antisymmetric: \(\alpha = (1 - \psi _{\rm sing}/\psi ) / (1 + \psi _{\rm sing}/\psi )\), ranges in \([-1, +1]\).

• twopunctures-averaged: average of the antisymmetric lapse and unity, ranges in \([0, +1]\).

• psiˆn: \(\alpha = \psi ^n\) where \(n\) is set by initial_lapse_psi_exponent (typically \(-2\)).

• W: \(\alpha = \psi ^{-2}\), suitable as a slow-start lapse for the W-formulation [5].

• brownsville: \(\alpha = 2/(1 + \psi ^n)\) [4].

4.4 Conformal Storage

When ADMBase::metric_type = "static conformal", the thorn populates the StaticConformal grid functions (psi, psix, etc.) in addition to the physical metric, with the level of detail controlled by StaticConformal::conformal_storage.

4.5 ADM Mass Iteration

Setting give_bare_mass = no enables an outer iteration that finds bare masses \(m_\pm \) such that the resulting ADM masses match target_M_plus and target_M_minus to within adm_tol. The iteration uses the analytic inversion of the ADM mass relation.

4.6 Interaction With Other Thorns

TwoPunctures writes to the ADMBase grid functions gxx, gyy, gzz, gxy, gxz, gyz, kxx, kyy, kzz, kxy, kxz, kyz, and alp. It exposes the scalar grid functions mp, mm, mp_adm, mm_adm, E, J1, J2, and J3 for use by other thorns (e.g., Multipole, QuasiLocalMeasures).

Matter source terms can be supplied via the optional function interfaces Set_Rho_ADM and Set_Initial_Guess_for_u, enabling neutron star or BH-NS initial data when use_sources = yes.

5 History

The spectral algorithm implemented in this thorn was developed by Marcus Ansorg and is described in [3]. The thorn was originally written for the Cactus framework by Marcus Ansorg and Erik Schnetter. It has since been maintained and extended by the Einstein Toolkit community.

References

6 Parameters

initial_lapse_psi_exponent	Scope: restricted	REAL
Description: Exponent n for psin initial lapse profile
Range		Default: -2.0
(*:*)	Should be negative

7 Interfaces

General

Implements:

twopunctures

Inherits:

admbase

staticconformal

grid

Grid Variables

7.0.1 PRIVATE GROUPS

7.0.2 PUBLIC GROUPS

Adds header:

TwoPunctures.h

TP_utilities.h

TP_params.h

TP_Guts.h

8 Schedule

This section lists all the variables which are assigned storage by thorn EinsteinInitialData/TwoPunctures. Storage can either last for the duration of the run (Always means that if this thorn is activated storage will be assigned, Conditional means that if this thorn is activated storage will be assigned for the duration of the run if some condition is met), or can be turned on for the duration of a schedule function.

Storage

Scheduled Functions

CCTK_PARAMCHECK (conditional)

  twopunctures_paramcheck

  check parameters and thorn needs

ADMBase_InitialData (conditional)

  twopunctures_group

  twopunctures initial data group

CCTK_INITIAL (conditional)

  twopunctures_group

  twopunctures initial data group

TwoPunctures_Group (conditional)

  twopunctures

  create puncture black hole initial data

TwoPunctures_Group (conditional)

  twopunctures_metadata

  output twopunctures metadata