GWFrames compatibility layer

GWFrames is an old python package that was too hard to upgrade to python 3, and therefore became impossible to maintain. Many of its capabilities have been superseded by the scri and sxs packages, but not the post-Newtonian capabilities. This package now supersedes the latter. For convenience, the PostNewtonian.GWFrames submodule provides a simple function to mimic the arguments used with the original GWFrames package to obtain a PN waveform (with a couple extra arguments). We won't bother to provide a return type that can fully mimic the object returned by the original GWFrames package, though it should be similar to the ones used by scri and sxs.

This interface can be used directly from Julia, but is more likely to be wanted by python users. To borrow from the example in the standard python interface documentation:

# Any python imports you need go here
import numpy as np
import matplotlib.pyplot as plt

# Start the Julia session
from sxs.julia import GWFrames

# Declare the essential parameters
delta = 0.0
chi1 = [0.1, 0.2, 0.3]
chi2 = [-0.2, 0.1, -0.3]
Omega_i = 0.01

# Call into Julia to run some function
w = GWFrames.PNWaveform("TaylorT1", delta, chi1, chi2, Omega_i)

# Plot the magnitudes of all the modes as functions of time
plt.semilogy(w.t, np.abs(w.data))

This should produce precisely the same plot as the python example, but using the GWFrames interface. Note that the parameters given here are required, but optional parameters are the same as in the original GWFrames interface — which are listed in full in the docstring below.

When called directly from Julia, the result will just be a NamedTuple with fields described in the docstring below. When called from python as shown above, the resulting w will be an sxs.WaveformModes object, with all the usual methods and properties, as well as several fields relevant to post-Newtonian waveforms:

  • M1 (array): The primary mass as a function of time.
  • M2 (array): The secondary mass as a function of time.
  • chi1 (array): The primary spin as a function of time.
  • chi2 (array): The secondary spin as a function of time.
  • frame (array): The quaternionic frame as a function of time.
  • v (array): The orbital velocity as a function of time.
  • orbital_phase (array): The orbital phase as a function of time.
PostNewtonian.GWFrames.PNWaveformFunction
PNWaveform(Approximant, delta, chi1_i, chi2_i, Omega_orb_i; kwargs...)

Compute a PN waveform, with the same call signature as GWFrames.PNWaveform

This is essentially a compatibility layer for the corresponding function in the original GWFrames Python package, with several additional optional arguments: inertial, dt, quiet, ell_min, ell_max, Lambda1, and Lambda2 (see below). Also, this function accepts optional arguments either as positional arguments (which the original GWFrames requires) or as keyword arguments.

Warning

We do not expect the result of this function to be identical to the result from the GWFrames Python package. In particular, this package uses more general expressions for the tidal-heating terms, fixes an error in the 2PN quadratic-spin terms for the waveform modes, uses more accurate (and efficient) ODE integration, and uses a more accurate method to compute the number of steps per orbit (by default).

Also note that there are differences in the order of operations for computing intermediate variables. Cancellation and roundoff error caused by these differences can have surprisingly large effects on the orbital evolution in many cases — particularly for precessing systems. Results from the two packages have been painstakingly analyzed, leading to the conclusion that all differences are caused by such errors or the differences in formulations mentioned above.

The Julia interface is more detailed, flexible, and efficient than the simple GWFrames interface that this function emulates. In particular, orbital_evolution takes essentially all the same arguments that DifferentialEquations.solve takes, and returns a solution that provides dense output and more details about the ODE solution itself. For example, one reason this function is more efficient than GWFrames is that we can use dense output to solve with fewer timesteps, while accurately and efficiently interpolating to the requested timesteps. While orbital_evolution solves for the dynamics, coorbital_waveform or inertial_waveform provides the actual waveform; both are returned by this function.

Required arguments

  • Approximant: Currently, supported arguments are "TaylorT1", "TaylorT4", and "TaylorT5".
  • delta: Fractional mass difference $(M₁-M₂)/(M₁+M₂)$
  • chi1_i: Normalized spin vector $S⃗₁/M₁²$
  • chi2_i: Normalized spin vector $S⃗₂/M₂²$
  • Omega_orb_i: Orbital angular frequency at initial instant

Optional arguments

As mentioned above, the following may be given either as positional arguments in this order (though any number of them may be omitted from the end), or as keyword arguments.

  • Omega_orb_0=Omega_orb_i: Orbital angular frequency at first instant found in data. If this is less than Omega_orb_i, the system is integrated backwards in time from the latter value to this value.
  • R_frame_i=Rotor(1): Initial orientation of the frame.
  • MinStepsPerOrbit=32: Number of time steps in the output data per orbit. Because the waveform modes go as high as $m=8$, this number should be at least 16 to avoid Nyquist aliasing in those modes. Note that this value may be overridden by dt (see below).
  • PNWaveformModeOrder=4.0: Maximum PN order of terms in the waveform formulas.
  • PNOrbitalEvolutionOrder=4.0: Maximum PN order of terms in the orbital-evolution formulas.
  • inertial=false: If true, transform waveform to the inertial frame; otherwise, the waveform will be in the co-orbital frame.
  • dt=0: Uniform time step size of the output. If this is not a strictly positive number, MinStepsPerOrbit will be used instead.
  • quiet=true: If false, show informational messages about the reasons for terminating the ODE integration. In either case, warnings will still be issued if terminating for bad or suspicious reasons. See the documentation of orbital_evolution for an example of how to filter warnings also.
  • ell_min=2: The lowest ℓ value in the output waveform.
  • ell_max=8: The highest ℓ value in the output waveform.
  • Lambda1=0: Tidal-coupling parameter of object 1.
  • Lambda2=0: Tidal-coupling parameter of object 2.

Returned values

This function returns a NamedTuple with the following keys:

  • t: The vector of time steps at which the data are evaluated. The time $t=0$ corresponds to the initial values that are arguments to this function.
  • data: Matrix of complex values of the mode weights. The shape is length(t) x 77. The first dimension enumerates the values at each instant of time. The second dimension enumerates the modes, starting with $(2,-2)$, then $(2,-1)$, up to $(2,2)$, followed by $(3,-3)$, and so on up to $(8,8)$. This is the opposite ordering as results from GWFrames, but the same as the ordering used by the sxs and scri packages. However, also note that certain conversions between Julia and Python may transpose matrices, because Julia is Fortran-ordered by default, whereas numpy is C-ordered. It is best to check the shape manually to be sure which dimension is which.
  • frame: Matrix of shape length(t) x 4 representing the frame-orientation quaternion as a function of time t.
  • M1, M2: Vectors of the respective masses as functions of time t. Note that only at the time corresponding to Omega_orb_i will the total mass be precisely 1. Generally, tidal heating will lead to time-dependent masses.
  • chi1, chi2: Matrices of shape length(t) x 3 representing the spins as functions of time t.
  • v: PN velocity parameter as a function of time t.
  • Phi: Orbital phase as a function of time t.

Because this is a NamedTuple, the fields can be accessed much like the fields of a WaveformModes object in the scri or sxs Python packages — as in w.t and w.data, where w is the object returned by this function.

source