## 2 Milestones

Numerical solving is a thousand-year-old art, which developed into modern numerical analysis several decades ago with the advent of modern computers and supercomputers. For a compelling account of the early history of numerical analysis and computing we refer the reader to Goldstine [359, 360].It is impossible to summarize all the important work on the subject in this review, but we find it instructive to list a chronogram of several relevant milestones taking us to 2014, in the context of GR. The following is a list – necessarily incomplete and necessarily biased – of works which, in our opinion, have been instrumental to shape the evolution of the field. A more complete set of references can be found in the rest of this review.

⋅ 1910 – The analysis of finite difference methods for PDEs is initiated with Richardson [648].

⋅ 1915 – Einstein develops GR [293, 294*].

⋅ 1916 – Schwarzschild derives the first solution of Einstein’s equations, describing the gravitational field generated by a point mass. Most of the subtleties and implications of this solution will only be understood many years later [687*].

⋅ 1917 – de Sitter derives a solution of Einstein’s equations describing a universe with constant, positive curvature . His solution would later be generalized to the case [255].

⋅ 1921, 1926 – In order to unify electromagnetism with GR, Kaluza and Klein propose a model in which the spacetime has five dimensions, one of which is compactified on a circle [463*, 476*].

⋅ 1928 – Courant, Friedrichs and Lewy use finite differences to establish existence and uniqueness results for elliptic boundary-value and eigenvalue problems, and for the initial-value problem for hyperbolic and parabolic PDEs [228].

⋅ 1931 – Chandrasekhar derives an upper limit for white dwarf masses, above which electron degeneracy pressure cannot sustain the star [193*]. The Chandrasekhar limit was subsequently extended to NSs by Oppenheimer and Volkoff [590*].

⋅ 1939 – Oppenheimer and Snyder present the first dynamical collapse solution within GR [589].

⋅ 1944 – Lichnerowicz [515*] proposes the conformal decomposition of the Hamiltonian constraint laying the foundation for the solution of the initial data problem.

⋅ 1947 – Modern numerical analysis is considered by many to have begun with the influential work of John von Neumann and Herman Goldstine [763], which studies rounding error and includes a discussion of what one today calls scientific computing.

⋅ 1952 – Choquet-Bruhat [327*] shows that the Cauchy problem obtained from the spacetime decomposition of the Einstein equations has locally a unique solution.

⋅ 1957 – Regge and Wheeler [641*] analyze a special class of gravitational perturbations of the Schwarzschild geometry. This effectively marks the birth of BH perturbation theory, even before the birth of the BH concept itself.

⋅ 1958 – Finkelstein understands that the surface of the Schwarzschild geometry is not a singularity but a horizon [320]. The so-called “golden age of GR” begins: in a few years there would be enormous progress in the understanding of GR and of its solutions.

⋅ 1961 – Brans and Dicke propose an alternative theory of gravitation, in which the metric tensor is non-minimally coupled with a scalar field [128*].

⋅ 1962 – Newman and Penrose [575*] develop a formalism to study gravitational radiation using spin coefficients.

⋅ 1962 – Bondi, Sachs and coworkers develop the characteristic formulation of the Einstein equations [118*, 667*].

⋅ 1962 – Arnowitt, Deser and Misner [47*] develop the canonical 3 + 1 formulation of the Einstein equations.

⋅ 1963 – Kerr [466*] discovers the mathematical solution of Einstein’s field equations describing rotating BHs. In the same year, Schmidt identifies the first quasar (quasi-stellar radio source) [681]. Quasars are now believed to be supermassive BHs, described by the Kerr solution.

⋅ 1963 – Tangherlini finds the higher-dimensional generalization of the Schwarzschild solution [740*].

⋅ 1964 – Chandrasekhar and Fock develop the post-Newtonian theory [194, 325].

⋅ 1964 – First documented attempt to solve Einstein’s equations numerically by Hahn & Lindquist [385]. Followed up by Smarr & Eppley about one decade later [710, 311].

⋅ 1964 – Seymour Cray designs the CDC 6600, generally considered the first supercomputer. Speeds have increased by over one billion times since.

⋅ 1964 – Using suborbital rockets carrying Geiger counters new sources of cosmic X-rays are discovered. One of these X-ray sources, Cygnus X-1, confirmed in 1971 with the UHURU orbiting X-ray observatory, is soon accepted as the first plausible stellar-mass BH candidate (see, e.g., [110]). The UHURU orbiting X-ray observatory makes the first surveys of the X-ray sky discovering over 300 X-ray “stars”.

⋅ 1965 – Penrose and Hawking prove that collapse of ordinary matter leads, under generic conditions, to spacetime singularities (the so-called “singularity theorems”) [608*, 401]. A few years later, Penrose conjectures that these singularities, where quantum gravitational effects become important, are generically contained within BHs – The cosmic censorship conjecture [610*, 767*].

⋅ 1965 – Weber builds the first GW detector, a resonant alluminium cylinder [771, 772].

⋅ 1966 – May and White perform a full nonlinear numerical collapse simulation for some realistic equations of state [543].

⋅ 1967 – Wheeler [661*, 778] coins the term black hole (see the April 2009 issue of Physics Today, and Ref. [779] for a fascinating, first-person historical account).

⋅ 1967, 1971 – Israel, Carter and Hawking prove that any stationary, vacuum BH is described by the Kerr solution [453, 188*, 403, 406*]. This result motivates Wheeler’s statement that “a BH has no hair” [661].

⋅ 1968 – Veneziano proposes his dual resonance model, which will later be understood to be equivalent to an oscillating string [759]. This date is considered the dawn of SMT.

⋅ 1969 – Penrose shows that the existence of an ergoregion allows to extract energy and angular momentum from a Kerr BH [610]. The wave analogue of the Penrose process is subsequently shown to occur by Zeldovich, who proves that dissipative rotating bodies (such as Kerr BHs, for which the dissipation is provided by the horizon) amplify incident waves in a process now called superradiance [827*, 828*].

⋅ 1970 – Zerilli [829, 830*] extends the Regge–Wheeler analysis to general perturbations of a Schwarzschild BH. He shows that the problem can be reduced to the study of a pair of Schrödinger-like equations, and applies the formalism to the problem of gravitational radiation emitted by infalling test particles.

⋅ 1970 – Vishveshwara [762] studies numerically the scattering of GWs by BHs: at late times the waveform consists of damped sinusoids (now called ringdown waves, or quasi-normal modes).

⋅ 1971 – Davis et al. [250*] carry out the first quantitative calculation of gravitational radiation emission within BH perturbation theory, considering a particle falling radially into a Schwarzschild BH. Quasi-normal mode (QNM) ringing is excited when the particle crosses the maximum of the potential barrier of the Zerilli equation, located close to the unstable circular orbit for photons.

⋅ 1973 – Bardeen, Carter and Hawking derive the four laws of BH mechanics [74].

⋅ 1973 – Teukolsky [743*] decouples and separates the equations for perturbations in the Kerr geometry using the Newman–Penrose formalism [575*].

⋅ 1973 – York [808, 809] introduces a split of the extrinsic curvature leading to the Lichnerowicz–York conformal decomposition, which underlies most of the initial data calculations in NR.

⋅ 1973 – Thorne provides a criterium for BH formation, the hoop conjecture [750*]; it predicts collapse to BHs in a variety of situations including very high-energy particle collisions, which were to become important in TeV-scale gravity scenarios.

⋅ 1974 – Hulse and Taylor find the first pulsar, i.e., a radiating neutron star (NS), in a binary star system [447]. The continued study of this system over time has produced the first solid observational evidence, albeit indirect, for GWs. This, in turn, has further motivated the study of dynamical compact binaries and thus the development of NR and resulted in the 1993 Nobel Prize for Hulse and Taylor.

⋅ 1975 – Using quantum field theory in curved space, Hawking finds that BHs have a thermal emission [405]. This result is one of the most important links between GR and quantum mechanics.

⋅ 1977 – NR is born with coordinated efforts to evolve BH spacetimes [708, 287, 711].

⋅ 1978 – Cunningham, Price and Moncrief [229, 230, 231] study radiation from relativistic stars collapsing to BHs using perturbative methods. QNM ringing is excited.

⋅ 1979 – York [810*] reformulates the canonical decomposition by ADM, casting the Einstein equations in a form now commonly (and somewhat misleadingly) referred to as the ADM equations.

⋅ 1980 – Bowen & York develop the conformal imaging approach resulting in analytic solutions to the momentum constraints under the assumption of maximal slicing as well as conformal and asymptotic flatness [121*].

⋅ 1983 – Chandrasekhar’s monograph [195*] summarizes the state of the art in BH perturbation theory, elucidating connections between different formalisms.

⋅ 1985 – Stark and Piran [724] extract GWs from a simulation of rotating collapse to a BH in NR.

⋅ 1985 – Leaver [504, 505, 506*] provides the most accurate method to date to compute BH QNMs using continued fraction representations of the relevant wavefunctions.

⋅ 1986 – McClintock and Remillard [547] show that the X-ray nova A0620-00 contains a compact object of mass almost certainly larger than , paving the way for the identification of many more stellar-mass BH candidates.

⋅ 1986 – Myers and Perry construct higher-dimensional rotating, topologically spherical, BH solutions [565*].

⋅ 1987 – ’t Hooft [736*] argues that the scattering process of two point-like particles above the fundamental Planck scale is well described and calculable using classical gravity. This idea is behind the application of GR for modeling trans-Planckian particle collisions.

⋅ 1989 – Echeverria [290] estimates the accuracy with which one can estimate the mass and angular momentum of a BH from QNM observations. The formalism is substantially refined in Refs. [97, 95*].

⋅ 1992 – The LIGO detector project is funded by the National Science Foundation. It reaches design sensitivity in 2005 [6]. A few years later, in 2009, the Virgo detector also reaches its design sensitivity [10].

⋅ 1992 – Bona and Massó show that harmonic slicing has a singularity-avoidance property, setting the stage for the development of the “1+log” slicing [115].

⋅ 1992 – D’Eath and Payne [256*, 257*, 258*, 259*] develop a perturbative method to compute the gravitational radiation emitted in the head-on collision of two BHs at the speed of light. Their second order result will be in good agreement with later numerical simulations of high-energy collisions.

⋅ 1993 – Christodoulou and Klainerman show that Minkowski spacetime is nonlinearly stable [219].

⋅ 1993 – Anninos et al. [37*] first succeed in simulating the head-on collision of two BHs, and observe QNM ringing of the final BH.

⋅ 1993 – Gregory and Laflamme show that black strings, one of the simplest higher-dimensional solutions with horizons, are unstable against axisymmetric perturbations [367*]. The instability is similar to the Rayleigh–Plateau instability seen in fluids [167*, 162]; the end-state was unclear.

⋅ 1993 – Choptuik finds evidence of universality and scaling in gravitational collapse of a massless scalar field. “Small” initial data scatter, while “large” initial data collapse to BHs [212*]; first use of mesh refinement in NR.

⋅ 1994 – The “Binary Black Hole Grand Challenge Project”, the first large collaboration with the aim of solving a specific NR problem (modeling a binary BH coalescence), is launched [542, 213].

⋅ 1995, 1998 – Through a conformal decomposition, split of the extrinsic curvature and use of additional variables, Baumgarte, Shapiro, Shibata and Nakamura [695*, 78*] recast the ADM equations as the so-called BSSN system, partly building on earlier work by Nakamura, Oohara and Kojima [569*].

⋅ 1996 – Brügmann [140*] uses mesh refinement for simulations of BH spacetimes in 3 + 1 dimensions.

⋅ 1997 – Cactus 1.0 is released in April 1997. Cactus [154] is a freely available environment for collaboratively developing parallel, scalable, high-performance multidimensional component-based simulations. Many NR codes are based on this framework. Recently, Cactus also became available in the form of the Einstein Toolkit [521, 300*].

⋅ 1997 – Brandt & Brügmann [126*] present puncture initial data as a generalization of Brill–Lindquist data to the case of generic Bowen–York extrinsic curvature.

⋅ 1997 – Maldacena [536*] formulates the AdS/CFT duality conjecture. Shortly afterward, the papers by Gubser, Klebanov, Polyakov [372*] and Witten [798*] establish a concrete quantitative recipe for the duality. The AdS/CFT era begins. In the same year, the correspondence is generalized to non-conformal theories in a variety of approaches (see [15] for a review). The terms “gauge/string duality”, “gauge/gravity duality” and “holography” appear (the latter had been previously introduced in the context of quantum gravity [737*, 734*]), referring to these generalized settings.

⋅ 1998 – The hierarchy problem in physics – the huge discrepancy between the electroweak and the Planck scale – is addressed in the so–called braneworld scenarios, in which we live on a four-dimensional subspace of a higher-dimensional spacetime, and the Planck scale can be lowered to the TeV [46*, 40*, 638*, 639*].

⋅ 1998 – First stable simulations of a single BH spacetime in fully dimensional NR within a “characteristic formulation” [508, 362], and two years later within a Cauchy formulation [23*].

⋅ 1998 – The possibility of BH formation in braneworld scenarios is first discussed [45, 69*]. Later work suggests BH formation could occur at the LHC [279*, 353*] or in ultra-high energy cosmic ray collisions [315, 33, 304].

⋅ 1999 – Friedrich & Nagy [335] present the first well-posed formulation of the initial-boundary-value problem (IBVP) for the Einstein equations.

⋅ 2000 – Brandt et al. [127] simulate the first grazing collisions of BHs using a revised version of the Grand Challenge Alliance code [227].

⋅ 2000 – Shibata and Uryū [698] perform the first general relativistic simulation of the merger of two NSs. More recent simulations [62], using a technique developed by Baiotti and Rezzolla that circumvents singularity excision [64], confirm that ringdown is excited when the merger leads to BH formation. In 2006, Shibata and Uryū perform NR simulations of BH-NS binaries [699].

⋅ 2001 – Emparan and Reall provide the first example of a stationary asymptotically flat vacuum solution with an event horizon of non-spherical topology – the “black ring” [307*].

⋅ 2001 – Horowitz and Maeda suggest that black strings do not fragment and that the end-state of the Gregory–Laflamme instability may be an inhomogeneous string [440], driving the development of the field. Non-uniform strings are constructed perturbatively by Gubser [371*] and numerically by Wiseman, who, however, shows that these cannot be the end-state of the Gregory–Laflamme instability [789*].

⋅ 2003 – In a series of papers [479*, 452*, 480*], Kodama and Ishibashi extend the Regge–Wheeler–Zerilli formalism to higher dimensions.

⋅ 2003 – Schnetter et al. [684*] present the publically available Carpet mesh refinement package, which has constantly been updated since and is being used by many NR groups.

⋅ 2005 – Pretorius [629*] achieves the first long-term stable numerical evolution of a BH binary. Soon afterwards, other groups independently succeed in evolving merging BH binaries using different techniques [159*, 65*]. The waveforms indicate that ringdown contributes a substantial amount to the radiated energy.

⋅ 2007 – First results from NR simulations show that spinning BH binaries can coalesce to produce BHs with very large recoil velocities [363, 161].

⋅ 2007 – Boyle et al. [122*] achieve unprecedented accuracy and number of orbits in simulating a BH binary through inspiral and merger with a spectral code that later becomes known as “SpEC” and uses multi-domain decomposition [618*] and a dual coordinate frame [678*].

⋅ 2008 – The first simulations of high-energy collisions of two BHs are performed [719*]. These were later generalized to include spin and finite impact parameter collisions, yielding zoom-whirl behavior and the largest known luminosities [697*, 720*, 717*, 716*].

⋅ 2008 – First NR simulations in AdS for studying the isotropization of a strongly coupled supersymmetric Yang–Mills plasma through the gauge/gravity duality [205*].

⋅ 2009 – Dias et al. show that rapidly spinning Myers–Perry BHs present zero-modes, signalling linear instability against axially symmetric perturbations [272*], as previously argued by Emparan and Myers [305*]. Linearly unstable modes were subsequently explored in Refs. [271*, 270*].

⋅ 2009 – Shibata and Yoshino evolve Myers–Perry BHs nonlinearly and show that a non-axisymmetric instability is present [701*].

⋅ 2009 – Collisions of boson stars show that at large enough energies a BH forms, in agreement with the hoop conjecture [216*]. Subsequent investigations extend these results to fluid stars [288*, 647*].

⋅ 2010 – Building on previous work [215], Lehner and Pretorius study the nonlinear development of the Gregory–Laflamme instability in five dimensions, which shows hints of pinch-off and cosmic censorship violation [511*].

⋅ 2010, 2011 – First nonlinear simulations describing collisions of higher-dimensional BHs, by Zilhão et al., Witek et al. and Okawa et al. [841*, 797*, 587*].

⋅ 2011 – Bizoń and Rostworowski extend Choptuik’s collapse simulations to asymptotically AdS spacetimes [108*], finding evidence that generic initial data collapse to BHs, thereby conjecturing a nonlinear instability of AdS.

⋅ 2013 – Collisions of spinning BHs provide evidence that multipolar structure of colliding objects is not important at very large energies [716*].