Rapid generation of fully relativistic extreme-mass-ratio-inspiral waveform templates for LISA data analysis

Our paper on rapidly computing fully relativistic EMRI waveforms is out today at arXiv:2008.06071. This is the first time that waveforms with full harmonic content can be computed on timescales useful for LISA data analysis. With our new model you can compute a year-long EMRI waveform in 10s of seconds on a CPU and in < 1s on a GPU. This is the same timing as the kludge models but without making any weak-field approximations. All the code is available in the BHPToolkit. It’s written in Python, with C++ and CUDA backends, which makes it very easy to use — see the screenshot below.


At the moment the model computes adiabatic inspirals into a Schwarzschild black hole, but the techniques used and the modular code is extensible to post-adiabatic inspirals into a Kerr black hole (watch this space). The key observation with post-adiabatic EMRI waveforms (i.e., those that include self-force corrections) is that although the waveform phase needs to be known to post-adiabatic order, the instantaneous waveform amplitudes only need to be known at leading order. Self-force corrections also add a short orbital timescale to the equations of motion. This can greatly slowdown the calculation of the inspiral but we know how to over this with either near-identity transformations or a two-timescale framework.


The main bottleneck to rapidly computing EMRI waveforms is thus interpolating the mode amplitudes and summing over thousands of modes at each time step (and each year-long waveform will have millions of time steps at LISAs sampling rate). For the interpolation we used order-reduction methods and a neural network. Neural networks are not great at high precision interpolation but we don’t need that as the error in the amplitudes does not accumulate over the inspiral. Both the neural network for the waveform amplitudes and the splines we use to sample the waveform can be evaluated extremely efficiently on a GPU. This results in year-long waveforms that can be computed in 100s of milliseconds.


The waveforms are also very accurate, with a mismatch of less than 5e-4 (relative to slow to compute high precision waveforms) across the strong-field Schwarzschild inspiral parameter space (e <= 0.7)


Finally, let me emphasize again the modular framework of the code. This should make it easy to add post-adiabatic corrections as well as environmental effects or beyond-GR/standard model physics that effect the inspiral phasing.

The (virtual) Black Hole Perturbation Toolkit Workshop

Before the lockdown we applied for and had been awarded GWVerse COST funding to host a workshop for the Black Hole Perturbation Toolkit in Prague in the Czech Republic. Unfortunately we had to put that face-to-face meeting on hold – hoping to host the workshop in March 2021 instead – but this opened a new opportunity: a virtual workshop. For the original workshop we had expected 25-30 participants but by the time the virtual workshop came around we had almost 200 participants registered! We had a range of great talks from contributors to the Toolkit and many of them also prepared awesome exercises, notebooks and extra material help people get to grips with various pieces of the Toolkit. It was amazing to see ~110 people in the Zoom room and ~20-25 people watching live on YouTube on the first day. All the talks were recorded on YouTube and can be watched again on the BHPToolkit YouTube channel and the talk material (slides, notebooks) can be found on the workshop website.

At the time of writing the videos have already amassed over 700 views. The first two days of the workshop were devoted to introducing the tools to new users with an emphasis on giving worked examples and exercises showing how the code works and how to use the code in your own projects. We also had short contributed talks from 8-9 people who spoke for 5 mins each on a topic related to, or already using, the Toolkit. It was awesome to meet new people using our code to make interesting calculations. The third day was the developer day and we set up a Slack channel to facilitate interactions between developers and potential new developers. This was a great success and already new code has been added to the Toolkit and I can see a lot more coming by looking at the spike of activity in the development branches.

The virtual workshop format worked far better than I could have ever imagined. A huge thank-you to our “local” organizing committee (Vojtech Witzany and Georgios Lukes-Gerakopoulos) for all their work behind the scenes; a huge thank you to all the speakers who gave consistently excellent talks that they had clearly put a lot of effort into preparing; and finally a huge thank you for all the attendees who turns up on mass and engaged throughout. It’s clear that virtual workshops give an opportunity to reach a large audience: we had 84 attendees from Europe, 37 from North America, 16 from Asia, and 5 from South America (rough stats based on an analysis of the email addresses people used to sign up – usually their host institution addresses but I had to remove all the .coms etc so not a perfect analysis). The map below shows the countries, based on my rough analysis, that attendees came from.



Given the huge success of the virtual workshop I think that even once the lockdowns are all over we will give serious consideration to running a virtual workshop again.

(virtual) Talk at the Gravity Exploration Institute at Cardiff University

I gave my second virtual talk of the lockdown today at the Gravity Exploration Institute of Cardiff University. Screenshot 2020-05-19 19.01.42


The talk was on three techniques for rapidly computing accurate EMRI waveforms. By combining these techniques we will be able to compute year-long EMRI waveforms in under a second. The three techniques are:

  • Non-linear black hole perturbation theory. This is needed in order to compute the post-adiabatic corrections to the waveform phase that will leave a residual error in the phase that scales with the mass ratio (so it the error in the waveform phase should be tiny for an EMRI). Our first paper on this topic presenting a result for a binary was recently published in Physical Review Letters. In this talk I also presented some exciting new (preliminary) results where we have computed the second-order (in the mass-ratio) flux for the first time. There are still many checks to go (and for some reason the l=2,m=2 mode is not working yet) but it feels great to be reaching a major goal with this long term project.
  • Near-identity (averaging) transformations. Once the oscillatory pieces of the self-force are included in the inspiral models the phase space trajectory and the orbital phases oscillate (with the number of oscillations scaling with the mass ratio). Numerically computing the phase space trajectory and the phases thus becomes very slow as the integrator has to resolve tens to hundreds of thousands of small oscillations. This can be overcome using so-called near-identity (averaging) transformations, or NITs. These average out the short timescale physics to produce equations of motion that capture the correct long-term behaviour of the system but that do not oscillate on the orbital timescale. The resulting NIT’ed equations of motion can be solved in milliseconds rather than minutes. Maarten van de Meent and I published a paper on this in Classical and Quantum Gravity, which looked at the general transformation and applied to to inspirals in Schwarzschild spacetimes. Since then one of my PhD students, Philip Lynch, has explicitly derived the NIT for generic inspiral into a Kerr black hole and has a numerical model implemented for equatorial orbits in Kerr spacetime.
  • Fast EMRI waveforms. To compute an EMRI waveform you need to sum over many thousands of harmonics at each time step. A straightforward implementation of this approach is quite slow, resulting in year-long waveforms that take 10s-100s of second to produce. In collaboration with Alvin Chua, Michael L. Katz, and Scott Hughes, I presented our new method for computing year-long EMRI waveforms in ~200ms. This involves a computation of a reduced-order model (ROM), neural networks and acceleration using graphics processor units (GPUs).

For those who are interested here’s a link to the slides from my talk.

The show must go on

Recently I was invited to give a seminar as part of the Maynooth Theoretical Physics Seminar Series and with the lockdown it went online, so earlier today I gave my first virtual seminar. On the whole it worked pretty well. My talk was titled ‘Gravitational Waves from Galactic Cores’ and reviewed modelling of small mass ratio binaries. Towards the end I highlighted our recent result on second-order calculations (now published in Physical Review Letters). I also included material highlighting our work from last spring on gravitational waves from objects orbiting Sgr A*, which unfortunately I didn’t get time to talk about during the seminar.

Screenshot 2020-04-17 18.14.04

For those who are interested here’s a link to the slides from the talk.

Talk at Caltech

I recently visited Caltech to work with Alvin Chua (and Michael Katz, who was also visiting) on a new fast method to compute full relativistic EMRI waveforms that can incorporate all known self-force corrections to the phase. We made great progress during the week and are now writing up that work. I’ll post on that again when the paper is out.

Whilst at Caltech I gave a talk on our recent second-order in the mass-ratio result. That talk was recorded in case you are interested:

One week, two papers

It’s been a busy end to the year and this week I had two papers out on the arXiv. I briefly review them below.

The location of the last stable orbit in Kerr spacetime

The first was with Leo Stein. One of the key differences between Newtonian and Einstein gravity is the presence of unstable and plunging orbits around compact objects. The surface in parameter space which separates stable and unstable orbits is called the ‘separatrix’. The location of the separatrix is non-trivial and in this work we present a new approach to finding it based upon finding roots of a high-order polynomial. This allows for analytic study of the separatrix. This allows us to easily derive many limits in previously calculated in the literature (such as the ISCO, ISSO, IBSO, etc). The polynomial is of order 12 in ‘p’ (the semi-latus rectum) so there is are closed form solutions in general but we provide a robust numerical scheme for finding the roots.


Screenshot 2019-12-18 11.13.02

The above figure shows the an extended view in (p,e,x) space of the solutions of the separatrix polynomial. The physical separatrix is one of the sheets (restricted to 0 <= e <=1).


Dissipation in extreme-mass ratio binaries with a spinning secondary

In this work we examine the gravitational-wave flux from a compact binary where the secondary is spinning. We derive the balance law in this case, which relates the flux to the local change in the secondary’s orbital constants of motion (energy, angular momentum). We then perform an explicit calculation for the case of a spin-aligned binary where the primary is a non-rotating (Schwarzschild) black hole. We make this calculation using the Teukolsky formalism reconstructing in the radiation gauge for the local calculation (we make a post-Newtonian and numerical calculation using this method), and a numerical Lorenz-gauge calculation. We find perfect agreement between the approaches and show explicitly that our balance law holds.

Screenshot 2019-12-20 08.23.53Screenshot 2019-12-20 08.24.09

The above two plots show the comparison of the flux to infinity between the PN and numerical results, and the Lorenz-gauge metric perturbation for the (2,2)-mode for a spinning body. All the PN and numerical results can be found in the Black Hole Perturbation Toolkit.

Second-order self-force calculation of gravitational binding energy

It is with great pleasure that we present today the first second-order (in the mass-ratio) self-force result. From the monopole piece of the second-order metric perturbation, and specializing to circular orbits in Schwarzschild spacetime, we compute the binding energy of the binary. We have compared the results of our calculation with the result from the first law of binary mechanics (FLBM) and find great agreement — see the figure below.

Screenshot 2019-08-21 09.34.25

In comparing against the FLBM we find a slight discrepancy in the very strong field. This discrepancy is not unexpected as the FLBM only applies to a conservative spacetime (no inspiral) whereas we include the effects of dissipation in our calculation.

This work builds upon almost 7 years of effort since the foundational second-order papers. It’s been a long road with 10+ papers already out and I can think of at least another 5 in preparation.

This first concrete calculation really marks the end of the beginning. As we compute higher modes we will get the binding angular momentum, the radiated flux, gauge invariants, local force, and then the first complete post-adiabatic waveform. Watch this space!


Traveling by train

I think my favourite way to travel is by train so when the opportunity arises, if possible, I will try to take the train there, back or both. A good example of this is that I was fortunate enough to be invited to lecture at a summer school in Beijing in 2016. Whilst I flew out there the only logical choice was to take the train back (or at least as far as St. Petersburg).

One thing I’ve been thinking of doing for a while now is making map of all the rail journey’s I’ve been on. This has resulted in my Rail Journey Map which shows all the intercity rail journeys I have been on totalling ~24,000 km of travel.

Compiling this map was rather challenging and it took me some time to research how to do it. Below is the method I used. It might not be the best, but it worked for me.

Creating the rail journey map

To create the map made use of a piece of open-source software called OpenRailRouting. The instructions are pretty clear on GitHub page about how to install and run it but a key piece of information is missing, at least for those not used to using OpenStreetMap (OSM) data.

To get the OSM data I had to download the 44GB planet PBF file and then use the Osmium tool to extract all the railway information. The command for this is:

osmium tags-filter -o planet-rail.osm.pbf planet.osm.pbf nw/railway

This results in a much more manageable ~350MB file. With this in hand you can run OpenRailRouting as described on the software’s webpage. This will setup a GraphHopper server locally which you can reach at localhost:8989. With this GUI you can select start and end points on a rail journey and the software will provide the route in between. You can also set intermediate points in case the software returns another route you didn’t take. Once you’ve settled on a route, you can export the GPX file by clicking the little ‘gpx’ button.

Finally, to plot all the routes together I used Mathematica as this made it very easy to change the map projection.

One useful thing to note is that the default configuration for OpenRailRouting only routes on standard gauge (1435mm) tracks. Many countries use this rail gauge but some, such as Ireland, Russia and Mongolia, use different gauges. You can modify the config file to add additional gauges easily.

Gravitational waves from bodies orbiting the Galactic Center black hole and their detectability by LISA

Today we put out a paper investigating bodies orbiting Sgr A* as a possible source of gravitational waves for LISA. It just so happens that the 4.1 million solar masses of Sgr A* places the gravitational wave (GW) frequency for bodies on a circular orbit in the strong-field right in the sweet spot for LISA. Our paper then addresses two questions: i) what sort of objects can get close enough to Sgr A* without being tidally disrupted and ii) for objects that are not tidally disrupted what would the resulting GWs and signal-to-noise (SNR) look like in LISA.

For the first question, clear candidates (from EMRI research) are stellar mass black holes, neutron stars and white dwarfs all of which can cross the innermost stable circular orbit (ISCO) without being tidal disrupted. We also investigate the Roche limit for planets (rocky and gaseous), low mass stars, red dwarfs and brown dwarfs. These results are nicely summarised in Fig. 1 of the paper:

LISA sensitivity curve and various gravitational wave frequencies from circular orbits around Sgr A*

This shows that detecting Jupiter-like or rocky planets is unlikely but low-mass stars or brown dwarves seems possible. We then go on to show that such object are all detectable in one year of LISA data with a signal-to-noise ratio above 10. We do not in this work attempt any estimate of event rates for these objects to be in the LISA band but hopefully it helps to motivate those who know how to calculate these things to look into this!

Lots more details can be found in our paper at: arXiv:1903.02049 where we also consider sources from another nearby massive black hole in M32 (see appendix C).

This work was in collaboration with Éric Gourgoulhon, Frederic Vincent and Alexandre Le Tiec.


LISA Waveform Working Group meeting

The first stand-alone meeting of the LISA Waveform working group (WavWG) will take place at the AEI Potsdam from 13 – 15 May 2019.

The WavWG is dedicated to the development of waveform templates for the LISA mission and connects the broader scientific community to the LISA Consortium. As such, the WavWG coordinates efforts to model a broad spectrum of gravitational wave sources, establishing links with LISA’s astrophysics, cosmology, fundamental physics and simulation working groups, and bridge between communities employing different methodologies to prepare waveforms.

This meeting is the first stand-alone meeting of the WavWG and aims at getting together leading international experts and young scientists to identify pressing tasks concerning LISA waveform modelling in its broadest sense and to foster new international collaborations. The meeting will have few talks per day and provide ample time for discussions.

While the WavWG is part of the LISA Consortium, you do not need to be in the Consortium to attend.

You can find additional information here: https://workshops.aei.mpg.de/lisawavwg/

The registration will end on *30 April 2019*.

We have limited funds available to support PhD students. Should you require financial support, please follow the instructions on our website and apply until *30 March 2019*.