# Combining Synchrotron X-ray Diffraction, Mechanistic Modeling, and Machine
Learning for _In Situ_ Subsurface Temperature Quantification during Additive
Manufacturing

Rachel E. Lim  Tuhin Mukherjee  Chihpin Chuang  Thien Q. Phan  Tarasankar
DebRoy  Darren C. Pagan<EMAIL_ADDRESS>Pennsylvania State University,
University Park, PA 16802, USA  Argonne National Laboratory, Lemont, IL 60439,
USA  Lawrence Livermore National Laboratory, Livermore, CA 14850, USA

###### Abstract

Laser powder bed fusion additive manufacturing (AM) produces extreme gradients
of temperature in both space and time, which in turn influence microstructural
development in the material. Qualification and model validation of the process
itself and resulting material produced necessitates the ability to
characterize these temperature fields. However, well-established means to
directly probe material temperature below the surface of an alloy while it is
being processed are limited. To address this gap in characterization
capabilities, we present a novel means to extract subsurface temperature
distribution metrics, with uncertainty, from _in situ_ synchrotron X-ray
diffraction measurements to provide quantitative temperature evolution during
laser melting. Temperature distribution metrics are determined using Gaussian
Process Regression supervised machine learning surrogate models trained with a
combination of mechanistic modeling (heat transfer and fluid flow) and X-ray
diffraction simulation. Trained surrogate model uncertainties are found to
range from 5% to 15% depending on the metric and current temperature. The
surrogate models are applied to extracting temperature metrics from an Inconel
625 nickel superalloy wall specimen during AM. Maximum temperatures of the
solid phase in the diffraction volume during processing are found to reach the
solidus temperature as expected, with mean and minimum temperatures found to
be several hundred degrees less. The highest accuracy of extracted temperature
metrics is found to be near melting, while temperatures during cooling are
likely influenced by the development of residual stress.

###### keywords:

synchrotron X-ray diffraction, additive manufacturing , superalloys , machine
learning , heat transfer and fluid flow modeling , Gaussian Process Regression

††journal: arXiv

##  1 Introduction

Primary factors for controlling the microstructure, porosity, and residual
stress state during additive manufacture (AM) of engineering alloys include
heat input, resulting temperature, and temperature gradients through the
component. In response, significant efforts have been undertaken to
characterize temperature profiles during AM builds in order to guide the build
design process. These efforts include both thermal and optical imaging of the
specimen surfaces during a build to estimate both surface temperature and
melt-pool shape [1, 2, 3, 4, 5, 6, 7], along with predicting defect formation.
While valuable, these characterization efforts only provide information about
the temperature profile at the sample surface, precluding understanding
critical subsurface thermal profiles during initial melting and subsequent
reheating events as layers are added above a volume of material of interest.
Rather, different measurement modalities are needed to probe an alloy’s
evolving microstructure during repeated thermal cycling encountered during the
build process.

To address these challenges, new synchrotron X-ray imaging and diffraction
capabilities have been developed that can characterize structure at rapid time
scales ($\ll s$). These measurements probe subsurface thermomechanical state
(thermal and mechanical lattice strain) and microstructure evolution in
conditions mimicking a wide range of additive manufacturing processes [8, 9,
10, 11, 12, 13, 14]. During these experiments, average temperatures of the
crystalline phases within a diffraction volume are estimated from the shifts
in diffraction peak centroid positions due to thermal strains. With knowledge
of the coefficient of thermal expansion (CTE) of the material, and assuming
equilibrium CTEs are valid during rapid cooling, thermal strains are mapped
directly to temperatures. While valuable, accuracy of temperatures from this
peak centroid analysis can be compromised due to the inherent spatial
gradients of temperature, mechanical loading, and chemistry during the build
process. However, while a complication for data analysis, information
regarding the spatial gradients of temperature (along with mechanical strains
and chemistry variation) within a diffraction volume is encoded into each
diffraction peak. Unfortunately, the single projection of X-rays through the
diffraction volume during these experiments prevents the direct inversion or
reconstruction of the temperature field in a tomographic fashion.

Here we propose a novel path forward to extracting this temperature data in
which a mechanistic heat transfer and fluid flow model and X-ray simulations
provide a framework for interpreting and analyzing complex experimental
temperature distributions. The simulations are used to create a collection of
reference diffraction patterns (a ‘dictionary’) representing different thermal
states. A strategy is then adapted from [15] to use Gaussian Process
Regression (GPR) to learn a mapping between spatial temperature distribution
metrics within a diffraction volume (generated from mechanistic modeling) and
diffraction peak shapes. Once trained, the GPR surrogate models extract
evolving subsurface temperature distribution metrics through comparison to the
reference simulated diffraction patterns. Uncertainty for extracted
temperature metrics are derived from how similar the experimental diffraction
data is to the reference simulated diffraction patterns.

##  2 Material and Experiment Description

The specimens in both experiment and simulations used in this work were AM
Inconel 625 (IN625) nickel superalloy made using laser powder bed fusion
(LPBF). The experimental specimen was a thin wall 3 mm in height along the
build direction ($\bm{z^{S}}$) by 0.53 mm in thickness ($\bm{y^{S}}$) and 20
mm in length ($\bm{x^{S}}$). The experimental specimen was extracted using
electro discharge machining from a block (50 mm $\times$ 15 mm $\times$ 6 mm)
built in an EOS M290 machine using manufacturer recommended build parameters.
The build layer thickness was 40 $\mu$m, with a 110 $\mu$m hatch spacing. The
laser power and speed were 285 W and 960 mm/s respectively. After build, the
block was stress-relief heat treated at 800 ∘C for 2 hr.

The _in situ_ diffraction measurements during laser melting were performed at
beamline 1-ID at the Advanced Photon Source (APS). Fig. 1a shows a schematic
of the experimental geometry for the measurements including AM IN625 wall
specimen, heating laser, and incoming X-ray beam, and the orientations of the
sample (‘S’) and laboratory (‘L’) coordinate system. In the laboratory
coordinate system, the incoming X-ray beam travels in the $-\bm{z^{L}}$
direction while the heating laser was nominally aligned along $\bm{y^{L}}$.
During X-ray measurements, the specimen remained fixed in the laboratory
coordinate system as the laser traveled in the $\bm{x^{L}}$ / $\bm{x^{S}}$
direction. The angle between incoming and diffracted X-rays $2\theta$ is
labeled and is related to the spacing of diffracting sets of lattice planes.
The incoming X-ray beam was 61.332 keV and was focused vertically by a set of
Si sawtooth lenses to dimensions of 50 $\mu$m $\times$ 30 $\mu$m along
$\bm{x^{L}}$ and $\bm{y^{L}}$ respectively. X-ray diffraction images were
measured by a Pilatus3 X CdTe 2M detector sitting 752 mm downstream of the
sample. The detector has 1475$\times$1679 pixels and a pixel size of
172$\times$172 ${\mu m^{2}}$. Diffraction images were collected with an
exposure time of 1 ms and a frequency of 250 Hz throughout the experiment.

Figure 1: a) Schematic of the experimental setup with ‘S’ superscripts for the
sample reference frame associated with the heat transfer and fluid flow
modeling and ‘L’ superscripts for the laboratory reference frame associated
with the laser heating and with the X-ray diffraction simulations. b) Photo of
the experimental measurement set-up used to collect X-ray diffraction data for
developing the temperature extraction framework. Marked are (1) the sample
holder and sample, (2) the laser test system, (3) X-ray area detector, and (4)
incoming beam direction.

Laser melting was performed using an existing _in situ_ LPBF simulator at the
APS. A picture of the simulator in Sector 1-ID is shown in Fig. 1b and more
complete system details can be found in [16]. The LBPF simulator uses a
ytterbium fiber laser (IPG YLR-500-AC, USA) and an $intelliSCAN_{de}$® 30 for
laser motion. Prior to laser melting, the environment chamber was purged and
re-filled with high-purity argon gas. During testing, the specimen was placed
2.9 mm away from the laser focal plane to create a spot diameter of 100
$\mu$m. During X-ray diffraction measurements, the laser was rastered over the
wall specimen along $\bm{x^{L}}$ with a laser power of approximately 120 W and
speed of 0.05 m/s. We note that this linear power density (2,400 J/m) is
relatively high in comparison to standard LPBF parameters for IN625. These
parameters were chosen to ensure a relatively large melt-pool and extended
temperature gradient through the thickness of the specimen. Figure 2 shows the
diffracted intensity integrated azimuthally around the detector (along the
diffraction rings) vs time. In the figure, we can see the shifting of the
diffracted intensity to lower $2\theta$ at 125 ms. This shift to lower
$2\theta$ is due to an increase in lattice plane spacing as the diffraction
volume is rapidly heated due to the laser passing over the diffraction volume.

Figure 2: Evolution of azimuthally integrated diffracted intensity $I$ with
time $t$ during laser melting of the IN625 wall specimen. Each row corresponds
to the diffraction line profile (intensity vs 2$\theta$) for a given time step
with the color signifying the magnitude of diffracted intensity.

##  3 Methods

In this section, an overview of the various methods employed to build the
dictionary of reference diffraction patterns linked to underlying thermal
distributions is given. A description of the GPR model, training, and
hyperparameter selection used to learn the mapping between underlying
temperature distribution metrics and diffraction line profiles is also
provided. As previously described in the Introduction, as opposed to using
X-ray diffraction data to validate a simulation, we are using combinations of
heat transfer and fluid flow modeling, X-ray diffraction simulation, and
machine learning to interpret and extract information from the experimental
data. A schematic of the various components of the effort are shown in Fig. 3
and the components of the modeling and training effort are described in more
detail in the following subsections. As a short summary, synthetic X-ray data
are generated using thermal fields output from the mechanistic heat transfer
and fluid flow model. We adopt an approach outlined in Fig. 3 in which pairs
of synthetic X-ray line profiles (input) and underlying temperature metrics
(output) are used to train various GPR surrogate models. Then, the surrogate
models, trained using synthetic data (the Source Domain), are _transferred_
[17] to analyze experimental data (the Target Domain) to determine temperature
metrics through, essentially, comparisons with the previously generated
synthetic data.

Figure 3: Flow chart showing the steps included in this work to extract
temperature metrics from X-ray diffraction data. Within the Source Domain,
heat transfer and fluid flow modeling is used to inform X-ray diffraction
simulations which are then used to train GPR models. The trained GPR models
are transferred to the Target Domain and used to predict temperature metrics
from the X-ray experimental data.

###  3.1 Heat Transfer and Fluid Flow Modeling

A heat transfer and fluid flow model is utilized to calculate the 3D transient
temperature and velocity fields during the laser melting of IN625 specimens.
The model is discussed in detail in Mukherjee et al. [18, 19] and only the
important features are described here. The model was developed using an in-
house Fortran code and compiled with an Intel Compiler. The model calculates
the melt-pool size, temperature fields, and velocity fields during the LPBF
process taking the laser parameters, alloy, and environmental gas properties
as inputs. The model considers temperature-dependent thermophysical properties
for both powder and fully dense material. The thermophysical properties of
IN625 required in the model are calculated using the commercial package
JMatPro, while the uncertain parameters such as absorptivity and the power
distribution of the laser beam can be adjusted to match the experimental data
on thermal cycles and deposit geometry. The model iteratively solves the
equations of conservation of energy, mass, and momentum in a 3D computational
domain consisting of the substrate, power bed, deposited layers and hatches,
and the shielding gas. The equations are discretized in the computational
domain using a finite difference scheme and a traveling grid approach is used
to increase computational efficiency. The model provides accurate results on
melt-pool geometry and temperature fields by considering the effects of the
convective flow of molten metals.

Using this model, a series of single laser trace simulations were performed
around conditions similar to the previously performed _in situ_ synchrotron
experiments (§4.2). For all simulations, a unidirectional scan along
$\bm{x^{S}}$ of the laser beam is used, with the laser beam direction being
$-\bm{z^{S}}$. Positive $\bm{y^{S}}$ is perpendicular to the laser scanning
direction and represents the direction along the width of the wall specimen.
The laser power and velocity were varied around the nominal experimental
parameters (120 W laser power and 50 mm/s laser speed). In total, nine laser
melting simulations were performed according to the simulation test matrix in
Table 1. Each simulation captured the laser traveling over a 200 ms interval
with 100 µs time steps while temperature fields were output every 1 ms. The
heat transfer and fluid flow simulations were performed on the ROAR
supercomputer at Penn State using 40 cores (2.13 GHz) for each simulation,
with each thermal simulation taking one hr to complete.

P = 100 W | P = 100 W | P = 100 W
---|---|---
v = 0.04 m/s | v = 0.05 m/s | v = 0.06 m/s
P = 120 W | P = 120 W | P = 120 W
v = 0.04 m/s | v = 0.05 m/s | v = 0.06 m/s
P = 140 W | P = 140 W | P = 140 W
v = 0.04 m/s | v = 0.05 m/s | v = 0.06 m/s
Table 1: Test matrix showing the nine sets of power-velocity parameters which
were used in the thermal simulation.

###  3.2 Synthetic X-ray Data Generation

Each of the nine laser trace simulations was used to generate a series of
time-dependent area detector X-ray diffraction patterns using the framework
developed in Pagan et al. [20]. This simulator uses diffraction calculation
and projection algorithms contained in the Python-based HEXRD software package
[21]. In the diffraction simulation framework, X-ray diffraction is simulated
within a polychromatic (Laue) diffraction framework that employs a finite-
energy bandwidth to capture realistic beam conditions. This method provides
benefits over angular based diffraction solution methods in conditions where
the specimen is not rotating, such as the single laser trace experiments
described above. Here the average X-ray energy and bandwidth ($|\Delta E|/E$)
were chosen to match the experiment, 61.332 keV and of 5$\times 10^{-4}$
respectively.

In this framework, diffraction events are simulated from discretized volumes
in space, with the spatial positions of each volume being incorporated into
diffracted ray tracing calculations. For clarity, we will refer to these
discretized volumes as scattering volumes, while the total volume illuminated
by the X-ray beam is the diffraction volume. Inside each scattering volume,
individual grains (lattice orientations) are inserted from which diffraction
events are calculated. In these diffraction simulations, the grains have no
morphological features and their orientations are randomly generated. The
thermomechanical state of each scattering volume, and the accompanying changes
in lattice state, can vary spatially. Here, temperature fields produce
spatially varying thermal strain within the diffraction volume. Matching the
experiment, the simulated diffraction volume is 50 $\times$ 30 $\times$ 530
$\mu$m3, while each scattering volume was 20 $\times$ 20 $\times$ 20 $\mu$m3.
Each scattering volume contained two randomly oriented grains simulating
grains with approximately 25 $\mu$m equivalent diameter. Each grain contains
1∘ of lattice misorientation to provide some amount of diffraction peak
broadening.

To model the thermal strain’s effects on the measured diffraction data, the
lattice structure of an embedded grain is altered by stretching the reciprocal
lattice vectors, $\bm{g}$, of each grain isotropically (valid for a cubic
crystal):

$\bm{g}=(1-\varepsilon_{T})\underline{\bm{I}}\cdot\bm{g_{0}}$ (1)

where $\bm{g_{0}}$ is the unstrained reciprocal lattice vector in a crystal,
$\underline{\bm{I}}$ is the second order identity tensor, and the thermal
strain $\varepsilon_{T}$ is given by

$\varepsilon_{T}=\int_{T_{0}}^{T}\alpha(T)\,dT.$ (2)

The temperature-dependent coefficient of thermal expansion $\alpha(T)$ used in
this work was measured independently using bulk dilatometry measurements at
the Penn State Center for Innovative Sintered Products. As only crystalline
material will diffract, if the temperature in a scattering volume exceeds the
solidus temperature (1563 K [22]), no diffraction events are recorded. For the
high-energy transmission geometry used in the experiment, all diffracting
X-rays have nearly the same path length regardless of whether they diffract
from the upstream or downstream side of the specimen. As such, absorption is
not considered since it will only scale the total integrated intensity of the
diffraction peak and not change the peak shape.

The laser-trace simulations described in §3.1 are used as input for generating
synthetic diffraction images. As the thermal simulations use a traveling grid
formulation for calculation, thermal fields at each time step were mapped to a
regular grid of scattering volumes with 20 $\mu$m spacing in all three
directions. Each thermal simulation was used to generate three sets of
synthetic diffraction data capturing the evolution of different temperature
gradients by placing the X-ray beam at different portions of the sample. For
each thermal simulation time series, the center of the X-ray beam was placed
20, 40, and 60 $\mu$m below the top of the specimen and separate 2D X-ray
diffraction image sets were generated. In total, 27 sets of X-ray simulations
and 5400 diffraction images total were generated for surrogate model training.
Once the 2D diffraction patterns were simulated for the entire time series,
each pattern was integrated azimuthally around the diffraction rings to create
1D diffraction line profiles (intensity vs. 2$\theta$). The integrated
diffraction line profile data covers 2$\theta$ angles ranging from $5^{\circ}$
to $13^{\circ}$ and encompasses the first six sets of lattice planes (111),
(200), (220), (311), (222), (400). After integration, background noise
corresponding to scattering in the experimental station was added to the
synthetic data. A comparison of example experimental (blue) and synthetic
(dashed red) diffraction line profiles in the unheated conditions is given in
Fig. 4. Each diffraction simulation of a heating time series at a single beam
position took approximately 4 hr to complete, with diffraction calculations
from each scattering volume parallelized over 40 cores.

Figure 4: Comparison between experimental (blue) and synthetic (dashed red)
diffraction line profiles (intensity $I$ versus Bragg Angle $2\theta$) from
the AM IN625 wall specimens in the unheated condition.

###  3.3 Gaussian Process Regression Surrogate Model Description

Here we utilize the supervised machine learning technique, Gaussian Process
Regression [23] implemented via scikit-learn [24], to learn a mapping between
diffraction line profiles and various temperature metrics in the diffraction
volume. Figure 5 shows example synthetic diffraction patterns generated using
the mechanistic AM and X-ray diffraction modeling colored by a temperature
metric of interest in the diffraction volume (maximum temperature,
$T_{\mathrm{Max}}$). The goal of using GPR is to learn these mappings between
line profiles and underlying thermal distributions such that with a
diffraction line profile, a temperature metric of interest can be extracted.
The method was recently applied to developing a mapping between diffraction
line profiles and dislocation configurations within diffraction volumes [15].
GPR takes a Bayesian statistical approach to surrogate model prediction. The
GPR method creates a normal distribution of mapping functions, informed by
training data, with the mean of function distribution serving as a model
prediction. The variance of the function distribution can serve as a
confidence bound or to inform where more training data may be necessary. In
GPR, model output predictions (i.e., temperature metrics) are constructed from
linear combinations of transformations of the input data (i.e., diffraction
line profile data). The transformation and linear weights are fit according to
the input training data and a chosen covariance (kernel) function.

Figure 5: Simulated a) diffraction line profiles (intensity $I$ versus Bragg
Angle $2\theta$) and b) (111) diffraction peaks of IN625 colored by the
maximum temperature ($T_{\mathrm{Max}}$) in the diffraction volume which was
used to generate the pattern.

The most common kernel function, $k$, in GPR is the Gaussian kernel (also
often referred to as the squared exponential kernel, exponentiated quadratic,
or radial basis function kernel). The kernel dictates the weights used to make
predictions from training data points. The Gaussian kernel takes the form

$k(\bm{x_{a}},\bm{x_{b}})=\sigma^{2}\exp\left(-\frac{||\bm{x_{a}}-\bm{x_{b}}||^{2}}{2L^{2}}\right)$
(3)

where $\bm{x_{a}}$ and $\bm{x_{b}}$ are input data vectors and $\sigma^{2}$ is
the amplitude. For this work, input data vectors are intensity values in
diffraction line profiles. This kernel includes a length scale, $L$, which
controls the extent of influence of a data point, affecting the variance of
the function distribution. The Gaussian kernel is referenced because of its
representative behavior and familiarity with the shape of Gaussian functions.
As the distance between an input vector ($\bm{x_{a}}$) and a data point
($\bm{x_{b}}$) is decreased, the weight increases. Conversely, as the distance
increases the weights decay in an exponential fashion.

Here, we employ the related rational quadratic kernel which is equivalent to
the summation of many exponentiated quadratic kernels of different length
scales:

$k(\bm{x_{a}},\bm{x_{b}})=\left(1+\frac{||\bm{x_{a}}-\bm{x_{b}}||^{2}}{2\alpha
L^{2}}\right)^{-\alpha},\enspace\alpha>0$ (4)

where $\alpha$ is the relative weighting between large and small length
scales. Accordingly, increasing $\alpha$ reduces the amount of local variation
(slows the weighting decay rate), and when $\alpha\rightarrow\infty$, the
rational quadratic kernel converges to the exponentiated quadratic kernel. A
range of $L$ and $\alpha$ values were tested for model training, but the
closest fits to the training data (without overfitting) corresponded to $L=1$
and $\alpha=1$.

###  3.4 GPR Model Training and Temperature Metric Extraction

Prior to application of the GPR models to the experimental synchrotron data,
the accuracy of GPR predictions were evaluated using a set of reserved
simulations. GPR models were trained using 26 of the 27 synthetic diffraction
data sets (again nine laser parameter sets given in Table 1 with three beam
positions each) using a single processor. The reserved synthetic diffraction
time series data corresponds to conditions best matching the experiment: 120 W
laser power, 0.05 m/s laser speed, and placing the X-ray beam 20 $\mu$m below
the top of the sample. After GPR surrogate model training and testing of the
models against reserved simulated data, the trained GPR models were applied to
the X-ray diffraction data collected during the synchrotron experiment. In
this process, the diffraction data collected at each time step is treated as
an independent data point and used as input for the various GPR models. At the
end of the process, various temperature metric histories for the diffraction
volume probed during the experiment are generated.

##  4 Results

###  4.1 Surrogate Model Training

Figure 6 shows predictions of four temperature metrics within the diffraction
volume using trained GPR models (each metric has its own trained GPR model)
that use diffraction line profiles as input compared with the same temperature
metrics extracted directly from the reserved thermal simulations. The four
temperature metrics are mean $T_{\mathrm{Mean}}$ (Fig. 6a), maximum
$T_{\mathrm{Max}}$ (Fig. 6b), minimum $T_{\mathrm{Min}}$ (Fig. 6c), and median
$T_{\mathrm{Median}}$ (Fig. 6d) temperatures. As previously described, for
each GPR surrogate model prediction, the variance of the prediction associated
with the normal distribution of the mapping functions can also be calculated
and employed as an uncertainty. In Fig. 6, the uncertainty (square root of the
variance, standard deviation) of the GPR prediction is shown by the red error
bars. For each temperature metric, a linear regression line was fit to the GPR
predictions and is plotted with a black line. A blue dashed line corresponding
to perfect correlation between the GPR model predictions and the reserved
testing data is provided for comparison. The $R^{2}$ coefficient of
determination of the GPR predictions is also provided.

As a whole, there is very good agreement between the reserved testing data and
the GPR model predictions. No aphysical predictions are found across the
metrics such as predictions of temperature significantly below room
temperature or above the solidus temperature. It is noted that there is
generally more training and testing data in the lower temperature regions due
to the laser passing rapidly over the specimen, leading to generally increased
uncertainty at higher temperatures (see larger red error bars). This feature
is most notable for the maximum temperature predictions near the solidus
temperature which has the largest uncertainties (on the order of 200 K or
15%). Conversely, the mean temperature predictions have the smallest
uncertainty (on the order of 20 K to 40 K or 6% to 12%) which is not
surprising. The mean temperature is strongly correlated to the point of
highest intensity on the diffraction peak. This can be contrasted to the
minimum and maximum temperatures which will generally correspond to small
volumes of material contributing to the tails of the diffraction peaks. Again,
we emphasize here that these temperature metric predictions take into account
contributions to peak broadening from the spatial distribution of temperature
within the diffraction volume. In other words, the hottest and coldest regions
may not necessarily correspond to the most extreme tail positions on the
diffraction peak depending on the spatial location of the diffraction event.

Figure 6: Comparison of prediction of temperature metrics from trained GPR
surrogate models using simulated diffraction data input reserved from model
training. The metrics are the a) mean $T_{\mathrm{Mean}}$, b) maximum
$T_{\mathrm{Max}}$, c) minimum $T_{\mathrm{Min}}$, and d) median
$T_{\mathrm{Median}}$ of the temperature distribution present in the
diffraction volume. The red error bars correspond to the square root of the
variance (standard deviation) of the GPR predictions. A linear regression line
has been fit to the testing data and GPR predictions and is shown with a black
line. The dashed blue line corresponds to perfect correlation between reserved
testing data and GPR model predictions.

###  4.2 Application of GPR Surrogate Models to Experimental Data

Figure 7 shows the evolution of mean $T_{\mathrm{Mean}}$ (Fig. 7a), maximum
$T_{\mathrm{Max}}$ (Fig. 7b), minimum $T_{\mathrm{Min}}$ (Fig. 7c), and median
$T_{\mathrm{Median}}$ (Fig. 7d) temperatures within the diffraction volume
versus time. Again, the red error bars associated with each temperature metric
measurement correspond to the uncertainty in the extracted quantity as given
by the square root of variance of the GPR prediction.

The point where the moving laser passes over the diffraction volume is the
clear peak in all four metrics. At its highest point, the prediction for
$T_{\mathrm{Max}}$ is close to the solidus temperature for IN625 as expected,
since melted material will not contribute to the diffraction peaks. Similar to
the cross-validation with the simulated data, $T_{\mathrm{Max}}$ has the
highest uncertainty (largest error bars). We also observe across all four
metrics, that the temperature remains relatively high well after the laser has
passed over the diffraction volume (on the order of 200 K above room
temperature. As will be discussed, this may be due to mechanical stress
developing within the specimen upon rapid cooling, as a positive mean stress
in the volume probed will ‘appear’ as an elevated temperature. As previously
described, the uncertainty in the temperature predictions is related to how
close input diffraction line profiles are to data used for GPR model training.
With the transfer learning approach, if the simulations in the Source Domain
used to generate the training data are missing physics, such as the
development of stress due to thermal gradients, the accuracy in the Target
Domain will decrease. Taking this into account, temperature metric values
closer to the solidus temperature are likely the most accurate since at this
point thermal expansion is largest and the mechanical stresses are lowest.
This will be further discussed in §5.2.

Figure 7: The evolving a) mean $T_{\mathrm{Mean}}$, b) maximum
$T_{\mathrm{Max}}$, c) minimum $T_{\mathrm{Min}}$, and d) median
$T_{\mathrm{Median}}$ of the distribution of temperature within the
experimental X-ray diffraction volume with respect to time $t$ extracted using
the trained GPR surrogate models. The red error bars correspond to the square
root of the variance (standard deviation) of the GPR surrogate model
predictions.

##  5 Discussion

Here, we have described and demonstrated a novel approach for quantifying
temperature distributions within alloys during extreme heating processes that
utilizes _in situ_ synchrotron X-ray diffraction, mechanistic modeling and
X-ray simulation, and supervised machine learning. The approach was applied to
quantifying temperature metrics in the bulk of an IN625 specimen during high-
speed laser melting mimicking AM. The development of approaches such as that
presented in this paper are important for quantifying and controlling
temperature distributions during the AM build process. In turn, this
information can be used to accelerate process certification and optimization
of build routines to control microstructure and minimize defects. The
presented method is a significant advance from current X-ray diffraction based
approaches in the literature that are only capable of estimating the average
temperature in an illuminated volume [11, 12, 13]. Usually analytical
functions (e.g., Gaussian or Lorentzian) are fit to the diffraction line
profiles and temperature is calculated from shifts of diffraction peaks.
Often, to improve temperature extraction from peak fittings, experiments are
conducted such that melt pools are much larger and scanning conditions are
unrealistic. During this process numerous assumptions are made regarding the
temperature distribution present and X-ray interaction with the sample
increase the uncertainty of the estimation, most notably that all X-rays are
emitted from a point source and the shape of the function used to fit the
diffraction peak. In addition to decreasing accuracy with these assumptions,
all information about the temperature distribution in the illuminated volume
is lost.

Our effort addresses these shortcomings by directly accounting for realistic
spatial thermal gradients. With these gradients accounted for, temperature
metrics determined from experimental data will have increased accuracy, and we
can access information about thermal gradients, a major driver of
microstructure formation. Now we will examine the temperature distribution
that was probed during the _in situ_ measurements. While the method presented
is a major advance forward for bulk temperature quantification during AM, we
will also discuss means to further increase accuracy and to extend the method
to other quantities, such as melt-pool volume.

###  5.1 Temperature Distribution Evolution

A primary benefit of the approach presented is the ability to explore the
evolving _distributions_ of temperature present within the diffraction volume.
As an example, we can analyze the distribution of temperature within the
diffraction volume during the _in situ_ synchrotron experiment. Figure 8 shows
the evolution of the mean, maximum, minimum, and median temperature metrics
extracted using the various GPR surrogate models together. With regards to the
temperature distribution, of most interest is the range (difference of maximum
and minimum) of temperatures in the crystalline phase throughout the
diffraction volume. From the figure, we can see that there is nearly a 600 K
difference between the maximum and minimum temperatures in the solid phase as
the laser passes over the diffraction volume. We can also see the temperature
difference remains throughout cooling and is still over 100 K at the end of
measurement. As the mean and median temperatures are closer to the minimum
temperature than the maximum, we can infer that the bulk of the diffraction
volume remains cooler through thickness, matching intuition regarding
localized melting. Using these data we can also begin to establish lower
bounds for the temperature gradient present. As the dimension of the
diffraction volume corresponding to the specimen thickness (530 $\mu$m) is
significantly larger than those defined by the incoming beam (50 $\mu$m and 30
$\mu$m), we can assume that spread of temperature is primarily along the
thickness of the specimen. With this being the case, and the hottest portion
of the specimen being in the center of the specimen, the lower bound of the
temperature gradient can be estimated to be approximately 2250 K/mm (600K /
0.265 mm). With any liquid phase present being significantly hotter, the
gradient will be larger, but a lower bound is of value for process design.

Figure 8: Comparison of the evolution of temperature $T$ metrics with time $t$
extracted from the experimental diffraction data using the GPR surrogate
models.

###  5.2 Effects of Mechanical Elastic Strain and Stress

In §4.2, we observed that in Fig. 7, all temperature metrics appear to be
converging towards predicted temperatures greater than room temperature at the
end of data collection. As mentioned, this may be due to accumulation of
mechanical (residual) stress within the specimen. Due to large thermal
gradients in the material during laser melting, mechanical elastic strain and
accompanying stress distributions in the material are formed to maintain
deformation compatibility. These stresses can become large enough to drive
plastic flow and set-up residual elastic strain and stress fields that remain
in the material even upon cooling [25, 26, 27]. In a cubic alloy such as
IN625, tensile stresses and strains distort the crystal lattice, and
subsequently diffraction patterns, in the same fashion as increased
temperature does. In the case of a lattice distortion that would be equivalent
to an approximately 100 K increase in temperature (lattice strain of
$1.5\times 10^{-3}$), the same lattice strain would correspond to
approximately 300 MPa of tensile stress applied (assuming a Young’s modulus of
approximately 200 GPa). As this stress magnitude is well below the yield
stress of Inconel 625, it is reasonable to assume that temperatures determined
from GPR in the low temperature regime are artificially high due to the
presence of mechanical strains. Residual strain measured in as-built IN625
artifacts have been reported to approach 3$\times 10^{-3}$ [26].

However, it should be noted that elastic strains during most mechanical
deformation modes are primarily deviatoric while thermal strains in cubic
materials are solely volumetric. With regards to the diffraction data, most
stress states will distort diffraction rings around the area detector, while
heating will cause uniform contraction or expansion of the diffraction rings.
A recent experimental effort [28] has demonstrated that with a large panel
area detector and capturing full diffraction rings, the effects of mechanical
stresses could largely be decoupled from the effects of heating. Keeping this
in mind, the method presented here could be extended to train GPR surrogate
models with diffraction patterns including thermomechanical effects rather
than just thermal effects. The GPR surrogate models would also need to be
provided with independent diffraction line profiles from different azimuthal
angles around the detector, but this should be readily possible. If
successful, low temperature measurement accuracy will be significantly
enhanced.

###  5.3 Melt-Pool Volume Estimation

For the extraction of temperature distribution metrics, the GPR surrogate
models that have been presented are primarily learning connections between
diffraction peak shape and position with the temperature distributions present
in the illuminated diffraction volume. As has been mentioned, when a volume of
material melts, that volume will no longer contribute intensity to the
measured diffraction peaks. The total integrated intensity will therefore
reflect the volume of unmelted alloy and as well as the converse, the relative
volume of the melt-pool. To test this idea, a final GPR surrogate model was
trained in the same fashion as described above, but trained to connect
diffraction line profiles to the volume fraction of the melt-pool present.
This surrogate model was then provided the experimental X-ray data similarly
to §4.2. Figure 9 shows the evolution of the melt-pool volume fraction during
the experiment. We can see that the melt-pool volume fraction has a maximum of
0.3 when the laser passes over the diffraction volume. With the relatively
large laser beam size (100 $\mu$m) and high power density laser conditions, a
melt pool volume extending through approximately 30% of the diffraction volume
(150 $\mu$m) appears reasonable. While high-speed X-ray radiography has been
used extensively in the literature to image melt-pool size, diffraction and
the approach presented may provide a complimentary means to characterize sub-
surface melt-pool dynamics in alloys in addition to samples with relatively
little density difference between liquid and solid phases and their minimal
radiographic imaging contrast.

Figure 9: Estimated experimental melt-pool volume fraction (VF) versus time
$t$ determined from the experimental diffraction data using a trained GPR
surrogate model.

##  6 Summary and Conclusions

For the first time, bulk maximum and minimum temperatures of the solid phase
in an engineering alloy (along with temperature distribution information) are
extracted from _in situ_ AM measurements. This was made possible by using
high-fidelity mechanistic and supervised machine learning modeling to
determine quantities from the experimental data, as opposed to taking a
traditional approach of using experimental data to calibrate the models. The
approach consists of training Gaussian Process Regression surrogate models
using a combination of heat transfer and fluid flow modeling and X-ray
diffraction modeling. Each surrogate model links diffraction line profiles to
metrics describing the temperature distribution present within the diffraction
volumes including maximum, minimum, mean, and median temperature. The smallest
uncertainties determined from the GPR models are approximately 5% for the
minimum, median, and mean temperatures. In contrast, the largest for the
maximum temperature (approximately 15%), near the solidus temperature. The
trained surrogate models were successfully applied to extracting these metrics
from _in situ_ high-energy synchrotron diffraction data collected during
additive manufacturing (laser melting) of a thin wall of Inconel 625. Despite
the larger uncertainties from the GPR surrogate model output, temperature
metric accuracy is believed to be greatest in the high temperature regime and
decreases in the low temperature regime when applied to experimental data, as
residual stress forms and distorts the material present. The discussion
described future efforts to account for residual stress formation to increase
model accuracy and extend the approach to extract other information about the
volume probed including the volume of the melt-pool in the diffraction volume.

## Acknowledgements

REL and DCP were supported by the National Institute of Standards and
Technology under award 70NANB20H208. The work of TQP was performed under the
auspices of the U.S. Department of Energy by Lawrence Livermore National
Laboratory under Contract DE-AC52-07NA27344. Computations for this research
were performed on the Pennsylvania State University’s Institute for
Computational and Data Sciences’ Roar supercomputer. This research used
resources of the Advanced Photon Source, a U.S. Department of Energy (DOE)
Office of Science user facility operated for the DOE Office of Science by
Argonne National Laboratory under Contract No. DE-AC02-06CH11357. Mr. Ian
Wietecha-Reiman is thanked for help performing the dilatometry measurements.

## References

  * [1] S. Moylan, E. Whitenton, B. Lane, J. Slotwinski, Infrared thermography for laser-based powder bed fusion additive manufacturing processes, in: AIP Conference Proceedings, Vol. 1581, American Institute of Physics, 2014, pp. 1191–1196.
  * [2] J. C. Fox, B. M. Lane, H. Yeung, Measurement of process dynamics through coaxially aligned high speed near-infrared imaging in laser powder bed fusion additive manufacturing, in: Thermosense: thermal infrared applications XXXIX, Vol. 10214, SPIE, 2017, pp. 34–50.
  * [3] B. A. Fisher, B. Lane, H. Yeung, J. Beuth, Toward determining melt pool quality metrics via coaxial monitoring in laser powder bed fusion, Manufacturing letters 15 (2018) 119–121.
  * [4] M. Montazeri, A. R. Nassar, C. B. Stutzman, P. Rao, Heterogeneous sensor-based condition monitoring in directed energy deposition, Additive Manufacturing 30 (2019) 100916.
  * [5] A. J. Dunbar, A. R. Nassar, Assessment of optical emission analysis for in-process monitoring of powder bed fusion additive manufacturing, Virtual and Physical Prototyping 13 (1) (2018) 14–19.
  * [6] J.-B. Forien, N. P. Calta, P. J. DePond, G. M. Guss, T. T. Roehling, M. J. Matthews, Detecting keyhole pore defects and monitoring process signatures during laser powder bed fusion: A correlation between in situ pyrometry and ex situ x-ray radiography, Additive Manufacturing 35 (2020) 101336.
  * [7] A. Ashby, G. Guss, R. K. Ganeriwala, A. A. Martin, P. J. DePond, D. J. Deane, M. J. Matthews, C. L. Druzgalski, Thermal history and high-speed optical imaging of overhang structures during laser powder bed fusion: A computational and experimental analysis, Additive Manufacturing 53 (2022) 102669\.
  * [8] C. Kenel, D. Grolimund, J. L. Fife, V. A. Samson, S. Van Petegem, H. Van Swygenhoven, C. Leinenbach, Combined in situ synchrotron micro X-ray diffraction and high-speed imaging on rapidly heated and solidified Ti-48Al under additive manufacturing conditions, Scripta Materialia 114 (2016) 117–120.  doi:10.1016/j.scriptamat.2015.12.009.
  * [9] N. P. Calta, J. Wang, A. M. Kiss, A. A. Martin, P. J. Depond, G. M. Guss, V. Thampy, A. Y. Fong, J. N. Weker, K. H. Stone, et al., An instrument for in situ time-resolved x-ray imaging and diffraction of laser powder bed fusion additive manufacturing processes, Review of Scientific Instruments 89 (5) (2018) 055101.
  * [10] R. Cunningham, C. Zhao, N. Parab, C. Kantzos, J. Pauza, K. Fezzaa, T. Sun, A. D. Rollett, Keyhole threshold and morphology in laser melting revealed by ultrahigh-speed x-ray imaging, Science 363 (6429) (2019) 849–852.
  * [11] S. Hocine, H. Van Swygenhoven, S. Van Petegem, C. S. T. Chang, T. Maimaitiyili, G. Tinti, D. Ferreira Sanchez, D. Grolimund, N. Casati, Operando X-ray diffraction during laser 3D printing, Materials Today 34 (2020) 30–40.  doi:10.1016/j.mattod.2019.10.001.
  * [12] S. A. Oh, R. E. Lim, J. W. Aroh, A. C. Chuang, B. J. Gould, B. Amin-Ahmadi, J. V. Bernier, T. Sun, P. C. Pistorius, R. M. Suter, et al., High speed synchrotron x-ray diffraction experiments resolve microstructure and phase transformation in laser processed ti-6al-4v, Materials Research Letters 9 (10) (2021) 429–436.
  * [13] S. A. Oh, R. E. Lim, J. W. Aroh, A. C. Chuang, B. J. Gould, J. V. Bernier, N. Parab, T. Sun, R. M. Suter, A. D. Rollett, Microscale observation via high-speed x-ray diffraction of alloy 718 during in situ laser melting, JOM 73 (1) (2021) 212–222.
  * [14] V. Thampy, A. Y. Fong, N. P. Calta, J. Wang, A. A. Martin, P. J. Depond, A. M. Kiss, G. Guss, Q. Xing, R. T. Ott, et al., Subsurface cooling rates and microstructural response during laser based metal additive manufacturing, Scientific reports 10 (1) (2020) 1–9.
  * [15] D. Bamney, A. Tallman, L. Capolungo, D. E. Spearot, Virtual diffraction analysis of dislocations and dislocation networks in discrete dislocation dynamics simulations, Computational Materials Science 174 (2020) 109473.
  * [16] C. Zhao, K. Fezzaa, R. W. Cunningham, H. Wen, F. De Carlo, L. Chen, A. D. Rollett, T. Sun, Real-time monitoring of laser powder bed fusion process using high-speed x-ray imaging and diffraction, Scientific reports 7 (1) (2017) 1–11.
  * [17] K. Weiss, T. M. Khoshgoftaar, D. Wang, A survey of transfer learning, Journal of Big data 3 (1) (2016) 1–40.
  * [18] T. Mukherjee, H. L. Wei, A. De, T. DebRoy, Heat and fluid flow in additive manufacturing—Part I: Modeling of powder bed fusion, Computational Materials Science 150 (April) (2018) 304–313.  doi:10.1016/j.commatsci.2018.04.022.
  * [19] T. Mukherjee, H. L. Wei, A. De, T. DebRoy, Heat and fluid flow in additive manufacturing – Part II: Powder bed fusion of stainless steel, and titanium, nickel and aluminum base alloys, Computational Materials Science 150 (April) (2018) 369–380.  doi:10.1016/j.commatsci.2018.04.027.
  * [20] D. C. Pagan, K. K. Jones, J. V. Bernier, T. Q. Phan, A Finite Energy Bandwidth-Based Diffraction Simulation Framework for Thermal Processing Applications, Jom 72 (12) (2020) 4539–4550.  doi:10.1007/s11837-020-04443-7.
  * [21] J. V. Bernier, N. R. Barton, U. Lienert, M. P. Miller, Far-field high-energy diffraction microscopy: a tool for intergranular orientation and strain analysis, Journal of Strain Analysis for Engineering Design 46 (7) (2011) 527–547.  doi:10.1177/0309324711405761.
  * [22] R. Pawel, R. Williams, Survey of physical property data for several alloys, Tech. rep., Oak Ridge National Lab. (1985).
  * [23] C. K. Williams, C. E. Rasmussen, Gaussian processes for machine learning, Vol. 2, MIT press Cambridge, MA, 2006.
  * [24] F. Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, et al., Scikit-learn: Machine learning in python, Journal of machine learning research 12 (Oct) (2011) 2825–2830.
  * [25] Z. Wang, E. Denlinger, P. Michaleris, A. D. Stoica, D. Ma, A. M. Beese, Residual stress mapping in Inconel 625 fabricated through additive manufacturing: Method for neutron diffraction measurements to validate thermomechanical model predictions, Materials and Design 113 (2017) 169–177.  doi:10.1016/j.matdes.2016.10.003.
  * [26] T. Q. Phan, M. Strantza, M. R. Hill, T. H. Gnaupel-Herold, J. Heigel, C. R. D’Elia, A. T. DeWald, B. Clausen, D. C. Pagan, J. P. Ko, et al., Elastic residual strain and stress measurements and corresponding part deflections of 3d additive manufacturing builds of in625 am-bench artifacts using neutron diffraction, synchrotron x-ray diffraction, and contour method, Integrating Materials and Manufacturing Innovation 8 (3) (2019) 318–334.
  * [27] J. L. Bartlett, X. Li, An overview of residual stresses in metal powder bed fusion, Additive Manufacturing 27 (January) (2019) 131–149.  doi:10.1016/j.addma.2019.02.020.
  * [28] F. Schmeiser, E. Krohmer, N. Schell, E. Uhlmann, W. Reimers, Internal stress evolution and subsurface phase transformation in titanium parts manufactured by laser powder bed fusion—an in situ x-ray diffraction study, Advanced Engineering Materials 23 (11) (2021) 2001502.


# Revisiting the cosmic string origin of GW190521

Josu C. Aurrekoetxea<EMAIL_ADDRESS>Astrophysics,
University of Oxford, Oxford OX1 3RH, United Kingdom  Charlie Hoy
<EMAIL_ADDRESS>University of Portsmouth, Portsmouth, PO1 3FX, United
Kingdom  Mark Hannam  Gravity Exploration Institute, Cardiff University,
Cardiff, United Kingdom

###### Abstract

For the first time we analyse gravitational-wave strain data using waveforms
constructed from strong gravity simulations of cosmic string loops collapsing
to Schwarzschild black holes; a previously unconsidered source. Since the
expected signal is dominated by a black-hole ringdown, it can mimic the
observed gravitational waves from high-mass binary black hole mergers. To
illustrate this, we consider GW190521, a short duration gravitational-wave
event observed in the third LIGO–Virgo–KAGRA observing run. We show that
describing this event as a collapsing cosmic string loop is favoured over
previous cosmic string analyses by an approximate log Bayes factor of $22$.
The binary black hole hypothesis is still preferred, mostly because the cosmic
string remnant is non-spinning. It remains an open question whether a spinning
remnant could form from loops with angular momentum, but if possible, it would
likely bring into contention the binary black hole preference. Finally, we
suggest that searches for ringdown-only waveforms would be a viable approach
for identifying collapsing cosmic string events. This work opens up an
important new direction for the cosmic-string and gravitational-wave
communities.

Introduction.— The observation of gravitational-waves (GWs) [1, 2, 3, 4, 5]
has paved the way to search for new physics. Cosmic strings [6, 7, 8, 9, 10,
11] are a well-motivated example that naturally arise when the rapid cooling
of the universe triggers a phase transition [12, 13, 14]. Cosmic strings may
manifest themselves through several channels, such as imprints via lensing on
the Cosmic Microwave Background (CMB) [15, 16], and a stochastic background of
gravitational waves (SGWB) [17, 18, 19, 20, 19], which is the total integrated
power of incoherent GWs from all individual emissions that are too weak to be
detected. The LIGO–Virgo–KAGRA (LVK) collaboration currently searches for the
SGWB [21], and places constraints on the dimensionless string tension
$G\mu/c^{2}$; a key property that sets their gravitational coupling strength
and the energy scale of the phase transition, providing a unique link to the
early universe. Localized coherent events of cosmic strings can also be
searched for, if they are energetic enough to be directly detected.

To date, all observed GW signals are consistent with binary black hole (BBH)
and/or neutron star mergers [1, 2, 3, 4, 5], with signal-to-noise ratios that
are mainly gained during the inspiral phase (around $20-50\mathrm{Hz}$). The
absence or modification of this stage on the GW signal may serve as smoking
gun for new physics. A candidate with these properties is GW190521 [22, 23],
which featured as a short transient in the LIGO–Virgo–KAGRA (LVK) detection
pipelines due to its low frequency nature. The origin of this event has been
extensively discussed in the literature as a massive black hole binary merger
[24, 25, 26, 27, 28, 29, 30], eccentric encounter [31, 32, 33, 34, 35, 36, 37,
38] or new physics [39, 40, 41, 42, 43]. Another studied (and disfavoured)
hypothesis of interest to this letter is a cosmic string cusp [23], a
gravitational-wave burst released when a fragment of a string doubles back on
itself and moves at the speed of light. This waveform inherently lacks of a
quasi-normal mode ringdown phase, as no black hole is formed during the
process [44].

In this letter we revisit the cosmic string hypothesis of GW190521 using novel
GW waveforms of cosmic string loops collapsing to black holes [45]. Guided by
fully general relativistic field-theory simulations of Abelian Higgs strings,
we construct an example waveform for a loop with dimensionless string tension
$G\mu/c^{2}=10^{-7}$ that collapses to form a black hole, radiating $3.4\%$ of
its initial mass in gravitational waves, see Fig. 1. At a fiducial distance
$d_{L}=4\,\mathrm{Gpc}$, this corresponds to a cosmic string loop with radius
$R_{0}=2.85\,\mathrm{AU}$ observed with an inclination angle
$\iota=34^{\circ}$. We show that this collapsing cosmic string yields a
signal-to-noise ratio $\rho\approx 12$, greatly improving on the cusp result
($\rho\approx 10$). This maps to an approximate log Bayes’ factor of
$\log\mathcal{B}^{\mathrm{loop}}_{\mathrm{cusp}}=22$ in favour of the
collapsing cosmic string over existing cusp results. Although the collapsing
cosmic string is disfavoured relative to the BBH hypothesis ($\rho\approx
14$), a full Bayesian analysis that encompasses the full configuration space
could potentially show comparable agreement. The purpose of our analysis is to
show that while cusps are easily distinguishable from black-hole mergers,
collapsing loops can be excellent mimickers. This makes it difficult to
conclusively identify a high-mass binary black hole merger if collapsing
cosmic string loops occur in nature.

Figure 1: Gravitational waveform $h^{+}_{\ell m}$ of a collapsing cosmic
string loop with radius $R_{0}=2.85\,\mathrm{AU}$ and dimensionless string
tension $G\mu/c^{2}=10^{-7}$. This corresponds to an initial mass of
$M_{0}\approx 181{M_{\odot}}$ that collapses to form a remnant black hole of
mass $M\approx 175{M_{\odot}}$ in the source frame, radiating $\approx
6{M_{\odot}}$ in gravitational waves. The vertical dashed line depicts the
transition between the infall and black hole formation stages, where there is
a sudden change in the frequency content of the signal.

Collapsing string loop waveforms.— After cosmic strings form and the universe
evolves, long strings self-intersect producing a network with abundant closed
loops of a vast range of sizes [46, 47, 48, 49, 50, 51, 52, 53]. These loops
can oscillate and contribute to the SGWB [54, 55, 56, 57, 58, 59] and/or may
eventually circularize and collapse. The linear mass density $\mu$ and radius
$R_{0}$ of circular string loops set their rest mass via $M_{0}=2\pi
R_{0}\mu$. The Schwarzschild radius of such configurations can then be
expressed as

$R_{\mathrm{schw}}\equiv\frac{2GM_{0}}{c^{2}}=4\pi
R_{0}\frac{G\mu}{c^{2}}~{}.$ (1)

If during the loop’s lifetime its mass is enclosed within its Schwarzschild
radius, it will form a black hole [60, 61, 62, 63, 64]. Using the hoop
conjecture, this condition becomes rather generic when the loops are
circular111This is less likely for non self-intersecting loops [60], although
details about their circularization timescales remain unknown., as the initial
radius needed to form a black hole is very small

$R_{0}\gtrsim
10^{-26}\,\mathrm{m}\,\left(\frac{G\mu/c^{2}}{10^{-7}}\right)^{-3/2}$ (2)

with masses

$M_{0}\approx
100\,M_{\odot}\left(\frac{R_{0}}{1\,\mathrm{AU}}\right)\left(\frac{G\mu/c^{2}}{10^{-7}}\right)\,.$
(3)

Here $1\,\mathrm{AU}\approx 1.49\times 10^{8}\,\mathrm{km}$ denotes
astronomical units, so these solar-system-sized loops are small when compared
to the size of the cosmological network.

The gravitational waveform of a circular cosmic string loop can be described
in three stages [45]: infall, black hole formation and ringdown:

$h(t)=\begin{cases}h_{\mathrm{[infall]}}\qquad t<t_{\mathrm{bh}}\,,\\\
h_{\mathrm{[bh]}}\qquad t_{\mathrm{bh}}<t<t_{\mathrm{qnm}}\,,\\\
h_{\mathrm{[qnm]}}\qquad t_{\mathrm{qnm}}<t\,.\end{cases}$ (4)

Analogous to BBHs [65, 66], string loop waveforms can be constructed using
semi-analytical and numerical techniques.

Infall: During the early stages of the loop collapse, when its radius is much
larger than its Schwarzschild radius, local backreaction effects can be
neglected and the dynamics can be treated within the weak-field limit,

$h_{ij}^{[\mathrm{infall}]}(t)=\frac{4G}{rc^{4}}\int_{-\infty}^{\infty}d^{3}x~{}T_{ij}\left(t_{\mathrm{ret}}+\mathbf{x}\cdot\mathbf{n},\mathbf{x}\right),$
(5)

where $t_{\mathrm{ret}}\equiv t-r/c$, and $r$ and $\mathbf{n}$ are the
distance and direction of the observer [67]. In the above expression, $T_{ij}$
is the energy momentum tensor, which for a (infinitesimally thin) circular
cosmic string loop oscillating on the $z=0$ plane is given by

$T^{\alpha\beta}(t,\mathbf{x})=\mu\,\delta\left[\sqrt{x^{2}+y^{2}}-R(t)\right]\delta[z]\,U^{\alpha}U^{\beta},$
(6)

where $U^{\alpha}=\gamma\,(c,\,V(t)\sin(\phi),\,V(t)\cos(\phi),\,0)$ is the
four-velocity and $\gamma$ the Lorentz contraction factor. The radius and
velocity during the infall are given by [68, 62]

$\displaystyle R(t)$
$\displaystyle=R_{0}\cos\left(\frac{c\,t}{R_{0}}\right),\quad
V(t)=c\sin\left(\frac{c\,t}{R_{0}}\right)\,.$ (7)

Even if the periodicity of these solutions seems to suggest the existence of
ever-oscillating circular loops (via the trigonometric functions), this
description breaks down when the loop radius is comparable to the string
width. In this work we will focus on the last oscillation of a planar,
circular loop, where the endpoint is the formation of a black hole.
Geometrically, non-circular loops can still form black holes as long as the
loops are planar. We expect non-planar and non-self-intersecting loops will
eventually circularize via radiation of scalars, vectors, and gravitational
waves, but this timescale is not known. The maximum radius of a circular loop
$R_{0}$ sets the oscillation frequency and thus the frequency content of the
infall signal, which is approximately

$f_{\mathrm{[infall]}}\approx
10^{-3}\,\mathrm{Hz}\left(\frac{1\,\mathrm{AU}}{R_{0}}\right)\,.$ (8)

From Eqn. (5), it can be shown that,

$\displaystyle rh^{+}_{[\mathrm{infall}]}(t,\theta,\phi)$
$\displaystyle=R_{0}\frac{G\mu}{c^{2}}\mathcal{I}(t,\theta,\phi)$ (9)
$\displaystyle rh^{\times}_{[\mathrm{infall}]}(t,\theta,\phi)$
$\displaystyle=0$ (10)

where $h^{\times}_{[\mathrm{infall}]}=0$ due to the axial symmetry of the
collapse and the details of the numerical integral
$\mathcal{I}(t,\theta,\phi)$ can be found in the appendices of Ref. [45]. A
template bank can be constructed using the scaling relations

$\displaystyle h^{+}_{[\mathrm{infall}]}$
$\displaystyle\propto\frac{R_{0}}{r}\frac{G\mu}{c^{2}}\propto\frac{1}{r}\frac{GM_{0}}{c^{2}}\,,$
(11) $\displaystyle\Delta t_{[\mathrm{infall}]}$
$\displaystyle\propto\frac{R_{0}}{c}\,.$ (12)

Black hole formation: At early times, the weak-field GW description of the
collapsing loop during its infall stage is valid. However, this linearized
treatment breaks down when the loop’s size is comparable to its Schwarzschild
radius, and gravitational effects become significant. We define this time
$t_{\mathrm{bh}}$, when $R(t_{\mathrm{bh}})\approx R_{\mathrm{schw}}$ and can
be estimated via Eqns. (1) and (7), to be

$t_{\mathrm{bh}}\approx\frac{R_{0}}{c}\arccos\left(4\pi\frac{G\mu}{c^{2}}\right)~{}.$
(13)

This quantity sets the timescale of the infall, after which the waveform
$h^{+}_{\mathrm{[infall]}}$ is no longer valid and needs to be matched to the
black hole formation stage $h^{+}_{\mathrm{[bh]}}$. This is inherently a
strong-gravity phenomenon and thus relies on general relativistic simulations.
These become prohibitive when there exist several physical scales that need to
be resolved, as is the case of extreme mass ratio inspirals of black holes. In
the cosmic string scenario, we would need to resolve the string width
$\delta\approx 10^{-35}\mathrm{km}$ (for $G\mu/c^{2}=10^{-7}$), and the
Schwarzschild radius of the black hole $R_{\mathrm{schw}}\approx
500\mathrm{km}$ (for GW190521). In addition, the Lorentz contraction of the
loop along the collapsing direction at the latest stages of the infall is
$\gamma\propto 1/G\mu\approx 10^{11}$. This results in a separation of scales
$\delta/R_{\mathrm{schw}}\ll 10^{-37}$, beyond the capabilities of adaptive
mesh refinement techniques [69]. However, results from numerical relativity
simulations of collapsing cosmic string loops with $G\mu\approx 10^{-3}$ and
$\delta/R_{\mathrm{schw}}\approx 10^{-2}$ exhibit a featureless intermediate
stage between the infall and the ringdown, and we will thus model
$h^{+}_{\mathrm{[bh]}}$ as a simple time-interpolating function between both
phases. We choose this interpolation such that the radiated energy in
gravitational waves agrees with the balance between the initial loop and final
mass of the remnant black hole.

Ringdown: Slightly after the black hole has formed $t_{\mathrm{qnm}}\gtrapprox
t_{\mathrm{bh}}+15GM/c^{3}$, the remnant enters the ringdown stage,
$h^{+}_{\mathrm{[qnm]}}$. This is modelled as a linear sum of damped sinusoids

$h^{+}_{\mathrm{[qnm]}}(t)=\frac{1}{r}\sum_{\ell m}A_{\ell
m}\exp\left[{{i\mkern 1.0mu}\omega_{\ell m}\frac{c^{3}}{GM}t}\right]$ (14)

where $\omega_{\ell m}=\mathrm{Re}(\omega_{\ell m})+{i\mkern
1.0mu}\mathrm{Im}(\omega_{\ell m})$ are a set of complex frequencies known as
quasi-normal modes. The real and imaginary parts set the oscillation frequency
and damping rates of the sinusoids and can be obtained from perturbative
calculations [70, 71], whilst the amplitudes $A_{\ell m}$ depend on the
strong-field regime details of the event itself. In our case, the remnant
black hole is non-spinning and predominantly radiates in the $(\ell,m)=(2,0)$
multipole.

The radiation of energy-flux generally results in a permanent displacement
$\Delta h=h(\infty)-h(-\infty)$ so-called gravitational-wave memory [72, 73,
74, 75, 76]. The magnitude can be estimated using the Christodolou formula
$\Delta h\sim 4GE_{\mathrm{rad}}/c^{4}r$, where $E_{\mathrm{rad}}$ is the
total energy radiated. From general relativistic simulations [45] and
geometrical arguments by Hawking [77], collapsing circular cosmic string loops
are expected to radiate between $1\%-29\%$ of their initial mass in
gravitational waves. The expected memory222The significant memory obtained in
[45] was due to the large anisotropic emission of matter, which is not the
case for our waveforms here as $\delta\ll R_{\mathrm{Schw}}$. for sources that
feature in the LVK frequency band is then given by $\Delta h\lessapprox
10^{-27}\,\left(1\,\mathrm{Gpc}/r\right)$, so we will ignore it when
constructing the waveforms.

Figure 2: The time-domain gravitational wave strain data for the LIGO-
Livingston detector [78, 79] around the time of GW190521 [22]. The top, middle
and bottom panels plot the best-fit detector-frame waveforms assuming a cosmic
string cusp, loop collapse and BBH merger respectively. The cosmic string cusp
waveform is the best-fitting template from the LVK’s cosmic string matched
filter search pipeline [80, 81] and the BBH waveform is the reconstructed
maximum likelihood waveform from the LVK analysis [22, 82]. The strain data
and waveforms are whitened by the noise power spectral density provided as
part of the LVK data release [82]. We shift the data and waveforms by the
coalescence time in the LIGO-Livingston detector $t_{\mathrm{c}}$ as reported
in the LVK BBH analysis.

Data analysis.— Due to limited numerical relativity simulations for collapsing
cosmic string loops and the difficulty of studying low string tensions, we do
not construct a parameterized model, and hence are unable to perform a
template based matched filter search or a Bayesian analysis to determine the
parameters that maximises the signal-to-noise ratio (SNR). Instead, we
calculated the matched-filter SNR for an example cosmic string loop with
$G\mu/c^{2}=10^{-7}$ and $R_{0}=2.85\,\mathrm{AU}$ that radiates $3.4\%$ of
its initial mass in gravitational waves. Although these may not correspond to
the parameters that best fit GW190521, we demonstrate that our example
waveform is still favoured when compared to previously studied cosmic string
cusps.

The construction of waveforms is generally done in the source frame, whilst
the analysis is in the detector frame. The transition between source-frame and
detector-frame waveforms relies on additional parameters describing the sky
position of the source, as well as its inclination $\iota$ and redshift $z$.
When these sources are at cosmological distances, we infer a larger remnant
mass in the detector frame $(1+z)M$ owing the redshifted waveforms. Similarly,
the amplitude of the observed GWs is suppressed by a factor of
$\sin^{2}\iota/d_{L}$ where $d_{L}$ is the luminosity distance to the
source333Note that this is different from the expression used for BBHs, since
the amplitude suppression depends on the multipole considered. BBHs
predominantly emit in the ($\ell,m)=(2,2)$ while circular collapsing loops in
the ($\ell,m)=(2,0)$ multipole.. As a result of the short timescales over
which a significant fraction of the initial mass is radiated in gravitational
waves, the amplitude of the waveform is large, see Fig. 1. Suppressing the
amplitude to typical GW strain levels restricts the allowed values of
$\sin^{2}\iota/d_{L}$, and due to the degeneracy between $\iota$ and $d_{L}$,
this can take any value between $0$ and $1/d_{L}$. We are therefore free to
choose a distance and inclination angle. For this analysis, we use a fiducial
distance $d_{L}=4\,\mathrm{Gpc}$, chosen to match the LVK BBH analysis [23],
and inclination angle $\iota=34^{\circ}$.

In Fig. 2 we plot the best-fit cusp [80, 81] and loop waveforms with the
strain data in the LIGO–Livingston GW detector. For comparison, we also add
the the BBH waveform, which is the reconstructed maximum likelihood waveform
from the LVK analysis [22, 82]. We whiten the data and waveforms using the
publicly available power spectral densities (PSDs) provided as part of the LVK
data release [82]. Owing to the short duration of GW190521 (approximately
$0.1\,\mathrm{s}$), we additionally bandpass the strain data between
$(10-256\,\mathrm{Hz})$ to suppress high frequencies. At early stages of the
loop waveform, we see that the infall is restrained in comparison to Fig. 1.
This is due to its low frequency nature, which is below the sensitivity of the
LVK detectors [78, 83, 84]. However, it retains the late-time ringdown, which
is similar to the BBH waveform444The best-fit BBH waveform favours a spinning
remnant, which has a lower damping rate of the ringdown. This allows it to
capture more quasi-normal mode cycles than the non-spinning remnant black hole
from the loop collapse.. We see that the loop waveform is dominated by the
black hole formation and ringdown stages, which solely depend on the mass of
the loop (and thus black hole). This introduces a degeneracy between the
radius and string tension via Eqn. (1), which can only be broken in the
presence of the infall signal as this timescale depends on $R_{0}$. This means
that comparable waveforms with different string tensions and radii can be
generated, as long as the mass of the loops remain unchanged.

To assess the probability that GW190521 was formed from a collapsing cosmic
string, we calculate the matched-filter SNR for our example string loop
waveform described above. This can be considered as analogous to a matched-
filter search [85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97] with only
one template. For comparison, we repeat the calculation for the best-fit BBH
and cusp models. We calculate the SNR $\rho$ by evaluating the noise-weighted
inner product [98]. We use the publicly available PyCBC [90, 93] software to
perform the analysis. Often signal-consistency tests are additionally
performed to distinguish between genuine astrophysical signals, and non-
Gaussian noise artefacts [99]. We did not perform any signal-consistency tests
in this work as we are ultimately interested in comparing the collapsing
cosmic string, BBH and cusp models. In Fig. 3 we plot the SNR timeseries for
our example collapsing cosmic string loop and best-fit BBH and cusp models. We
see that the loop waveform greatly out performs the cusp model with an
increase in SNR of $20\%$: $10$ to $12$. Assuming that the likelihood follows
$\exp(-\rho^{2}/2)$, we can approximate the $\log$ Bayes’ factor by taking the
ratio of likelihoods. We find that a $\log$ Bayes factor in favour of our loop
waveform over cusps555When calculating odd ratios priors in favour of one
hypothesis over the other must be specified. We do not take this into account
in our analysis but note that cusps generically form, whereas the
circularization timescale for generic loops is unknown. of $22$. This implies
that the loop hypothesis is strongly preferred over the cusp model, and is
therefore the most accurate representation of GW190521 assuming a cosmic
string origin to date. We see that the BBH hypothesis obtains a larger SNR,
meaning that it is still the preferred progenitor of GW190521. Without further
details about the low-frequency nature of GW190521 and a more complete
template bank of loop waveforms, we are unable of drawing further conclusions.

Figure 3: Network signal-to-noise ratios (assuming a network of the LIGO-
Hanford, LIGO Livingston and Virgo detectors) for the best-fit BBH, loop and
cusp waveforms, $\rho=\\{14,\,12,\,10\\}$ respectively. As with Fig. 2, we
shift the timeseries by the coalescence time in the LIGO-Livingston detector
$t_{\mathrm{c}}$ as reported in the LVK BBH analysis.

Conclusion.— In this letter we have revisited the origin of GW190521 using,
for the first time, waveforms of circular cosmic string loops collapsing to
black holes; a previously unconsidered source. We have shown that GW190521 can
be described as collapsing string loop with a signal-to-noise ratio of
$\rho\approx 12$. This greatly improves upon previous cosmic string cusp
analyses, where $\rho\approx 10$ is obtained. By approximating the $\log$
Bayes’ factor, we show that the collapsing cosmic string hypothesis is
favoured over existing cusp results by
$\log\mathcal{B}^{\mathrm{loop}}_{\mathrm{cusp}}=22$.

Our cosmic string loop waveforms are constructed from circular, planar loops
that have no angular momentum, and thus collapse to form Schwarzschild black
holes. Indeed, this is a very special case – in general loops will be non-
planar [100], or carry initial spin in its internal field configuration, or
possess traveling kinks. A loop with initial angular momentum might collapse
to form a black hole with significant spin. This would reduce the damping rate
of the ringdown phase and improve the SNR, so it is entirely conceivable that
such signals would be indistinguishable from binary black holes in an
observation like GW190521. However, the details about the angular momentum
loss rate during the loop collapse process are not known, so we leave this
study for future work.

The timescales shown in Fig. 1 highlight the multi-frequency nature of these
events. The data favours the absence of an oscillatory infall signal,
suggesting that it must be below the low-frequency sensitivity of the LVK
detectors, $f_{\mathrm{[infall]}}\lessapprox 10\,\mathrm{Hz}$. This implies
that $R_{0}\gtrapprox 10^{-4}\,\mathrm{AU}$ via Eqn. (8), and assuming the
progenitor was a cosmic string loop of $M_{0}\approx 181M_{\odot}$, we can put
a weak constraint on the dimensionless string tension of this event to be
$G\mu/c^{2}\lessapprox 10^{-2}$. This is not competitive with those obtained
from CMB [101] or SGWB [21] searches. Note however that these are sensitive to
details of the loop distribution model and how the string network evolves.
This is at present a point of debate [21] as it depends on numerical
simulations of large cosmological networks. An alternate channel of detection
via coherent signals of individual loop collapse can shed observational light
on this distribution and guide theoretical modeling of the cosmic string
network.

Based on our results, we propose collapsing string loops as a potential
hypothesis for future analyses. The degeneracy between high mass binary black
hole and string collapse signals will make a definitive detection of this
source difficult, although it might be possible to use higher modes to more
conclusively distinguish them. The absence of additional power prior to black-
hole formation may be a smoking gun in LVK observations, but it will be
necessary to understand how easily string collapse can be distinguished from
dynamical capture mergers that also lack significant pre-merger power. We
suggest using ringdown-only waveform templates as a simple and viable approach
to search for these events. Future GW detectors that target lower frequencies,
such as LISA [102, 103, 104], will be able to probe the early-phase of GW
events that feature at the edge of the LVK sensitivity, helping to distinguish
between high mass BBHs or signals arising from new physics.

_Acknowledgements._ — JCA would like to thank Thomas Helfer, Eugene Lim,
Samaya Nissanke and Andrew Williamson for earlier collaboration in cosmic
string simulations and discussions related to using higher modes to
distinguish exotica. We also thank Jade Powell, Juan Calderon Bustillo, and
Salvatore Vitale for comments on this manuscript. We are grateful for helpful
discussions with Jose Juan Blanco-Pillado, Katy Clough, Pedro Ferreira, Ian
Harry, Max Isi, Laura Nuttall, and Ken Olum. JCA acknowledges funding from the
Beecroft Trust and The Queen’s College via an extraordinary Junior Research
Fellowship (eJRF). CH thanks the UKRI Future Leaders Fellowship for support
through the grant MR/T01881X/1. MH thanks the Science and Technology
Facilities Council (STFC) grant ST/V00154X/1.

Part of this work used the DiRAC@Durham facility managed by the Institute for
Computational Cosmology on behalf of the STFC DiRAC HPC Facility
(www.dirac.ac.uk) under DiRAC RAC15 Grant ACTP316. The equipment was funded by
BEIS capital funding via STFC capital grants ST/P002293/1, ST/R002371/1 and
ST/S002502/1, Durham University and STFC operations grant ST/R000832/1. DiRAC
is part of the National e-Infrastructure.

This research has made use of data or software obtained from the Gravitational
Wave Open Science Center (https://www.gw-openscience.org), a service of LIGO
Laboratory, the LIGO Scientific Collaboration, the Virgo Collaboration, and
KAGRA. LIGO Laboratory and Advanced LIGO are funded by the United States
National Science Foundation (NSF) as well as the Science and Technology
Facilities Council (STFC) of the United Kingdom, the Max-Planck-Society (MPS),
and the State of Niedersachsen/Germany for support of the construction of
Advanced LIGO and construction and operation of the GEO600 detector.
Additional support for Advanced LIGO was provided by the Australian Research
Council. Virgo is funded, through the European Gravitational Observatory
(EGO), by the French Centre National de Recherche Scientifique (CNRS), the
Italian Istituto Nazionale di Fisica Nucleare (INFN) and the Dutch Nikhef,
with contributions by institutions from Belgium, Germany, Greece, Hungary,
Ireland, Japan, Monaco, Poland, Portugal, Spain. KAGRA is supported by
Ministry of Education, Culture, Sports, Science and Technology (MEXT), Japan
Society for the Promotion of Science (JSPS) in Japan; National Research
Foundation (NRF) and Ministry of Science and ICT (MSIT) in Korea; Academia
Sinica (AS) and National Science and Technology Council (NSTC) in Taiwan. This
material is based upon work supported by NSF’s LIGO Laboratory which is a
major facility fully funded by the National Science Foundation.

_Software_ – We used Matplotlib [105] NumPy [106], pycbc [93] and SciPy [107]
for analysis and plotting

## References

  * Abbott  _et al._ [2021a] R. Abbott _et al._ (LIGO Scientific, VIRGO, KAGRA), (2021a), arXiv:2111.03606 [gr-qc] .
  * Nitz  _et al._ [2023] A. H. Nitz, S. Kumar, Y.-F. Wang, S. Kastha, S. Wu, M. Schäfer, R. Dhurkunde, and C. D. Capano, Astrophys. J. 946, 59 (2023), arXiv:2112.06878 [astro-ph.HE] .
  * Venumadhav  _et al._ [2020] T. Venumadhav, B. Zackay, J. Roulet, L. Dai, and M. Zaldarriaga, Phys. Rev. D 101, 083030 (2020), arXiv:1904.07214 [astro-ph.HE] .
  * Zackay  _et al._ [2019] B. Zackay, T. Venumadhav, L. Dai, J. Roulet, and M. Zaldarriaga, Phys. Rev. D 100, 023007 (2019), arXiv:1902.10331 [astro-ph.HE] .
  * Zackay  _et al._ [2021] B. Zackay, L. Dai, T. Venumadhav, J. Roulet, and M. Zaldarriaga, Phys. Rev. D 104, 063030 (2021), arXiv:1910.09528 [astro-ph.HE] .
  * Kibble [1976] T. W. B. Kibble, J. Phys. A 9, 1387 (1976).
  * Vilenkin [1981a] A. Vilenkin, Phys. Rev. D 23, 852 (1981a).
  * Vilenkin [1981b] A. Vilenkin, Phys. Rev. D 24, 2082 (1981b).
  * Vilenkin [1985] A. Vilenkin, Phys. Rept. 121, 263 (1985).
  * Turok and Brandenberger [1986] N. Turok and R. H. Brandenberger, Phys. Rev. D 33, 2175 (1986).
  * Hindmarsh and Kibble [1995] M. B. Hindmarsh and T. W. B. Kibble, Rept. Prog. Phys. 58, 477 (1995), arXiv:hep-ph/9411342 .
  * Vachaspati and Vilenkin [1984] T. Vachaspati and A. Vilenkin, Phys. Rev. D 30, 2036 (1984).
  * Jeannerot  _et al._ [2003] R. Jeannerot, J. Rocher, and M. Sakellariadou, Phys. Rev. D 68, 103514 (2003), arXiv:hep-ph/0308134 .
  * Copeland and Kibble [2010] E. J. Copeland and T. W. B. Kibble, Proc. Roy. Soc. Lond. A 466, 623 (2010), arXiv:0911.1345 [hep-th] .
  * Vilenkin [1984] A. Vilenkin, Astrophys. J. Lett. 282, L51 (1984).
  * Ade  _et al._ [2014] P. A. R. Ade _et al._ (Planck), Astron. Astrophys. 571, A25 (2014), arXiv:1303.5085 [astro-ph.CO] .
  * Vilenkin [1981c] A. Vilenkin, Phys. Lett. B 107, 47 (1981c).
  * Vachaspati and Vilenkin [1985] T. Vachaspati and A. Vilenkin, Phys. Rev. D 31, 3052 (1985).
  * Allen and Shellard [1992] B. Allen and E. P. S. Shellard, Phys. Rev. D 45, 1898 (1992).
  * Hogan and Rees [1984] C. J. Hogan and M. J. Rees, Nature 311, 109 (1984).
  * Abbott  _et al._ [2021b] R. Abbott _et al._ (LIGO Scientific, Virgo, KAGRA), Phys. Rev. Lett. 126, 241102 (2021b), arXiv:2101.12248 [gr-qc] .
  * Abbott  _et al._ [2020a] R. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. Lett. 125, 101102 (2020a), arXiv:2009.01075 [gr-qc] .
  * Abbott  _et al._ [2020b] R. Abbott _et al._ (LIGO Scientific, Virgo), Astrophys. J. Lett. 900, L13 (2020b), arXiv:2009.01190 [astro-ph.HE] .
  * Fishbach and Holz [2020] M. Fishbach and D. E. Holz, Astrophys. J. Lett. 904, L26 (2020), arXiv:2009.05472 [astro-ph.HE] .
  * Nitz and Capano [2021] A. H. Nitz and C. D. Capano, Astrophys. J. Lett. 907, L9 (2021), arXiv:2010.12558 [astro-ph.HE] .
  * Estellés  _et al._ [2022] H. Estellés _et al._ , Astrophys. J. 924, 79 (2022), arXiv:2105.06360 [gr-qc] .
  * Capano  _et al._ [2022] C. D. Capano, J. Abedi, S. Kastha, A. H. Nitz, J. Westerweck, Y.-F. Wang, M. Cabero, A. B. Nielsen, and B. Krishnan, (2022), arXiv:2209.00640 [gr-qc] .
  * Xu and Hamilton [2023] Y. Xu and E. Hamilton, Phys. Rev. D 107, 103049 (2023), arXiv:2211.09561 [gr-qc] .
  * Sberna  _et al._ [2022] L. Sberna _et al._ , Phys. Rev. D 106, 064056 (2022), arXiv:2205.08550 [gr-qc] .
  * Siegel  _et al._ [2023] H. Siegel, M. Isi, and W. M. Farr, Phys. Rev. D 108, 064008 (2023), arXiv:2307.11975 [gr-qc] .
  * Romero-Shaw  _et al._ [2020] I. M. Romero-Shaw, P. D. Lasky, E. Thrane, and J. C. Bustillo, Astrophys. J. Lett. 903, L5 (2020), arXiv:2009.04771 [astro-ph.HE] .
  * Gayathri  _et al._ [2022] V. Gayathri, J. Healy, J. Lange, B. O’Brien, M. Szczepanczyk, I. Bartos, M. Campanelli, S. Klimenko, C. O. Lousto, and R. O’Shaughnessy, Nature Astron. 6, 344 (2022), arXiv:2009.05461 [astro-ph.HE] .
  * Tagawa  _et al._ [2021] H. Tagawa, B. Kocsis, Z. Haiman, I. Bartos, K. Omukai, and J. Samsing, Astrophys. J. Lett. 907, L20 (2021), arXiv:2010.10526 [astro-ph.HE] .
  * Gamba  _et al._ [2023] R. Gamba, M. Breschi, G. Carullo, S. Albanesi, P. Rettegno, S. Bernuzzi, and A. Nagar, Nature Astron. 7, 11 (2023), arXiv:2106.05575 [gr-qc] .
  * Iglesias  _et al._ [2022] H. L. Iglesias _et al._ , (2022), arXiv:2208.01766 [gr-qc] .
  * Guo  _et al._ [2022] W. Guo, D. Williams, I. S. Heng, H. Gabbard, Y.-B. Bae, G. Kang, and Z.-H. Zhu, Mon. Not. Roy. Astron. Soc. 516, 3847 (2022), arXiv:2203.06969 [gr-qc] .
  * Ramos-Buades  _et al._ [2023] A. Ramos-Buades, A. Buonanno, and J. Gair, (2023), arXiv:2309.15528 [gr-qc] .
  * Morton  _et al._ [2023] S. Morton, S. Rinaldi, A. Torres-Orjuela, A. Derdzinski, M. P. Vaccaro, and W. Del Pozzo, (2023), arXiv:2310.16025 [gr-qc] .
  * Calderón Bustillo  _et al._ [2021] J. Calderón Bustillo, N. Sanchis-Gual, A. Torres-Forné, J. A. Font, A. Vajpeyi, R. Smith, C. Herdeiro, E. Radu, and S. H. W. Leong, Phys. Rev. Lett. 126, 081101 (2021), arXiv:2009.05376 [gr-qc] .
  * Sakstein  _et al._ [2020] J. Sakstein, D. Croon, S. D. McDermott, M. C. Straight, and E. J. Baxter, Phys. Rev. Lett. 125, 261105 (2020), arXiv:2009.01213 [gr-qc] .
  * Clesse and Garcia-Bellido [2022] S. Clesse and J. Garcia-Bellido, Phys. Dark Univ. 38, 101111 (2022), arXiv:2007.06481 [astro-ph.CO] .
  * Abedi  _et al._ [2023] J. Abedi, L. F. Longo Micchi, and N. Afshordi, Phys. Rev. D 108, 044047 (2023), arXiv:2201.00047 [gr-qc] .
  * Calderon Bustillo  _et al._ [2022] J. Calderon Bustillo, N. Sanchis-Gual, S. H. W. Leong, K. Chandra, A. Torres-Forne, J. A. Font, C. Herdeiro, E. Radu, I. C. F. Wong, and T. G. F. Li, (2022), arXiv:2206.02551 [gr-qc] .
  * Damour and Vilenkin [2000] T. Damour and A. Vilenkin, Phys. Rev. Lett. 85, 3761 (2000), arXiv:gr-qc/0004075 .
  * Aurrekoetxea  _et al._ [2020] J. C. Aurrekoetxea, T. Helfer, and E. A. Lim, Class. Quant. Grav. 37, 204001 (2020), arXiv:2002.05177 [gr-qc] .
  * Martins and Shellard [1996] C. J. A. P. Martins and E. P. S. Shellard, Phys. Rev. D 54, 2535 (1996), arXiv:hep-ph/9602271 .
  * Ringeval  _et al._ [2007] C. Ringeval, M. Sakellariadou, and F. Bouchet, JCAP 02, 023 (2007), arXiv:astro-ph/0511646 .
  * Martins and Shellard [2006] C. J. A. P. Martins and E. P. S. Shellard, Phys. Rev. D 73, 043515 (2006), arXiv:astro-ph/0511792 .
  * Lorenz  _et al._ [2010] L. Lorenz, C. Ringeval, and M. Sakellariadou, JCAP 10, 003 (2010), arXiv:1006.0931 [astro-ph.CO] .
  * Blanco-Pillado  _et al._ [2011] J. J. Blanco-Pillado, K. D. Olum, and B. Shlaer, Phys. Rev. D 83, 083514 (2011), arXiv:1101.5173 [astro-ph.CO] .
  * Blanco-Pillado  _et al._ [2014] J. J. Blanco-Pillado, K. D. Olum, and B. Shlaer, Phys. Rev. D 89, 023512 (2014), arXiv:1309.6637 [astro-ph.CO] .
  * Allen and Shellard [1990] B. Allen and E. P. S. Shellard, Phys. Rev. Lett. 64, 119 (1990).
  * Martins and Shellard [2002] C. J. A. P. Martins and E. P. S. Shellard, Phys. Rev. D 65, 043514 (2002), arXiv:hep-ph/0003298 .
  * Blanco-Pillado and Olum [2017] J. J. Blanco-Pillado and K. D. Olum, Phys. Rev. D 96, 104046 (2017), arXiv:1709.02693 [astro-ph.CO] .
  * Afzal  _et al._ [2023] A. Afzal _et al._ (NANOGrav), Astrophys. J. Lett. 951, L11 (2023), arXiv:2306.16219 [astro-ph.HE] .
  * Ellis and Lewicki [2021] J. Ellis and M. Lewicki, Phys. Rev. Lett. 126, 041304 (2021), arXiv:2009.06555 [astro-ph.CO] .
  * Blasi  _et al._ [2021] S. Blasi, V. Brdar, and K. Schmitz, Phys. Rev. Lett. 126, 041305 (2021), arXiv:2009.06607 [astro-ph.CO] .
  * Blanco-Pillado  _et al._ [2021] J. J. Blanco-Pillado, K. D. Olum, and J. M. Wachter, Phys. Rev. D 103, 103512 (2021), arXiv:2102.08194 [astro-ph.CO] .
  * Ringeval and Suyama [2017] C. Ringeval and T. Suyama, JCAP 12, 027 (2017), arXiv:1709.03845 [astro-ph.CO] .
  * Hawking [1989] S. W. Hawking, Phys. Lett. B 231, 237 (1989).
  * Polnarev and Zembowicz [1991] A. Polnarev and R. Zembowicz, Phys. Rev. D 43, 1106 (1991).
  * Helfer  _et al._ [2019] T. Helfer, J. C. Aurrekoetxea, and E. A. Lim, Phys. Rev. D 99, 104028 (2019), arXiv:1808.06678 [gr-qc] .
  * Bramberger  _et al._ [2015] S. F. Bramberger, R. H. Brandenberger, P. Jreidini, and J. Quintin, JCAP 06, 007 (2015), arXiv:1503.02317 [astro-ph.CO] .
  * James-Turner  _et al._ [2020] C. James-Turner, D. P. B. Weil, A. M. Green, and E. J. Copeland, Phys. Rev. D 101, 123526 (2020), arXiv:1911.12658 [astro-ph.CO] .
  * Khan  _et al._ [2016] S. Khan, S. Husa, M. Hannam, F. Ohme, M. Pürrer, X. Jiménez Forteza, and A. Bohé, Phys. Rev. D 93, 044007 (2016), arXiv:1508.07253 [gr-qc] .
  * Bohé  _et al._ [2017] A. Bohé _et al._ , Phys. Rev. D 95, 044028 (2017), arXiv:1611.03703 [gr-qc] .
  * Maggiore and Press [2008] M. Maggiore and O. U. Press,  _Gravitational Waves: Volume 1: Theory and Experiments_, Gravitational Waves (OUP Oxford, 2008).
  * Nambu [1969] Y. Nambu, (1969).
  * Radia  _et al._ [2022] M. Radia, U. Sperhake, A. Drew, K. Clough, P. Figueras, E. A. Lim, J. L. Ripley, J. C. Aurrekoetxea, T. França, and T. Helfer, Class. Quant. Grav. 39, 135006 (2022), arXiv:2112.10567 [gr-qc] .
  * Kokkotas and Schmidt [1999] K. D. Kokkotas and B. G. Schmidt, Living Rev. Rel. 2, 2 (1999), arXiv:gr-qc/9909058 .
  * Berti  _et al._ [2009] E. Berti, V. Cardoso, and A. O. Starinets, Class. Quant. Grav. 26, 163001 (2009), arXiv:0905.2975 [gr-qc] .
  * Zel’dovich and Polnarev [1974] Y. B. Zel’dovich and A. G. Polnarev, Sov. Astron. 18, 17 (1974).
  * Braginsky and Thorne [1987] V. B. Braginsky and K. S. Thorne, Nature 327, 123 (1987).
  * Christodoulou [1991] D. Christodoulou, Phys. Rev. Lett. 67, 1486 (1991).
  * Thorne [1992] K. S. Thorne, Phys. Rev. D 45, 520 (1992).
  * Jenkins and Sakellariadou [2021] A. C. Jenkins and M. Sakellariadou, Class. Quant. Grav. 38, 165004 (2021), arXiv:2102.12487 [gr-qc] .
  * Hawking [1990] S. W. Hawking, Phys. Lett. B 246, 36 (1990).
  * Aasi  _et al._ [2015] J. Aasi _et al._ (LIGO Scientific), Class. Quant. Grav. 32, 074001 (2015), arXiv:1411.4547 [gr-qc] .
  * Abbott  _et al._ [2023] R. Abbott _et al._ (LIGO Scientific, VIRGO, KAGRA), (2023), arXiv:2302.03676 [gr-qc] .
  * Abbott  _et al._ [2018] B. P. Abbott _et al._ (LIGO Scientific, Virgo), Phys. Rev. D 97, 102002 (2018), arXiv:1712.01168 [gr-qc] .
  * Collaboration and Collaboration [2022a] L. S. Collaboration and V. Collaboration, “GW190521 parameter estimation samples and figure data,”  (2022a).
  * Collaboration and Collaboration [2022b] L. S. Collaboration and V. Collaboration, “GWTC-2.1: Deep Extended Catalog of Compact Binary Coalescences Observed by LIGO and Virgo During the First Half of the Third Observing Run - Parameter Estimation Data Release,”  (2022b).
  * Acernese  _et al._ [2014] F. Acernese, M. Agathos, K. Agatsuma, D. Aisa, N. Allemandou, A. Allocca, J. Amarni, P. Astone, G. Balestri, G. Ballardin, _et al._ , Classical and Quantum Gravity 32, 024001 (2014).
  * Somiya [2012] K. Somiya (KAGRA), Class. Quant. Grav. 29, 124007 (2012), arXiv:1111.7185 [gr-qc] .
  * Cannon  _et al._ [2012] K. Cannon, R. Cariou, A. Chapman, _et al._ , Astrophys. J. 748, 136 (2012), arXiv:1107.2665 [astro-ph.IM] .
  * Privitera  _et al._ [2014] S. Privitera, S. R. P. Mohapatra, P. Ajith, _et al._ , Physical Review D 89, 024003 (2014), arXiv:1310.5633 [gr-qc] .
  * Messick  _et al._ [2017] C. Messick, K. Blackburn, P. Brady, _et al._ , Physical Review D 95, 042001 (2017), arXiv:1604.04324 [astro-ph.IM] .
  * Hanna  _et al._ [2020] C. Hanna, S. Caudill, C. Messick, A. Reza, S. Sachdev, L. Tsukada, K. Cannon, K. Blackburn, J. D. E. Creighton, H. Fong, and et al., Physical Review D 101, 022003 (2020), arXiv:1901.02227 [gr-qc] .
  * Sachdev  _et al._ [2019] S. Sachdev, S. Caudill, H. Fong, R. K. L. Lo, C. Messick, _et al._ , arXiv:1901.08580 (2019).
  * Usman  _et al._ [2016] S. A. Usman _et al._ , Class. Quant. Grav. 33, 215004 (2016), arXiv:1508.02357 [gr-qc] .
  * Nitz  _et al._ [2017] A. H. Nitz, T. Dent, T. Dal Canton, S. Fairhurst, and D. A. Brown, Astrophys. J. 849, 118 (2017), arXiv:1705.01513 [gr-qc] .
  * Nitz  _et al._ [2018] A. H. Nitz, T. Dal Canton, D. Davis, and S. Reyes, Physical Review D 98, 024050 (2018), arXiv:1805.11174 [gr-qc] .
  * Nitz  _et al._ [2019] A. Nitz, I. Harry, D. Brown, C. M. Biwer, J. Willis, T. D. Canton, C. Capano, L. Pekowsky, T. Dent, A. R. Williamson, _et al._ , “gwastro/pycbc: Pycbc release v1.15.2,” Zenodo (2019).
  * Adams  _et al._ [2016] T. Adams, D. Buskulic, V. Germain, G. M. Guidi, F. Marion, M. Montani, B. Mours, F. Piergiovanni, and G. Wang, CQGra 33, 175012 (2016), arXiv:1512.02864 [gr-qc] .
  * Hooper  _et al._ [2012] S. Hooper, S. K. Chung, J. Luan, D. Blair, Y. Chen, and L. Wen, Physical Review D 86, 024012 (2012), arXiv:1108.3186 [gr-qc] .
  * Liu  _et al._ [2012] Y. Liu, Z. Du, S. K. Chung, S. Hooper, D. Blair, and L. Wen, CQGra 29, 235018 (2012).
  * Guo  _et al._ [2018] X. Guo, Q. Chu, S. K. Chung, Z. Du, L. Wen, and Y. Gu, CoPhC 231, 62 (2018).
  * Brown [2004] D. A. Brown,  _Searching for gravitational radiation from binary black hole MACHOs in the galactic halo_ , Other thesis (2004), arXiv:0705.1514 [gr-qc] .
  * Davies  _et al._ [2020] G. S. Davies, T. Dent, M. Tápai, I. Harry, C. McIsaac, and A. H. Nitz, Phys. Rev. D 102, 022004 (2020), arXiv:2002.08291 [astro-ph.HE] .
  * Matsunami  _et al._ [2019] D. Matsunami, L. Pogosian, A. Saurabh, and T. Vachaspati, Phys. Rev. Lett. 122, 201301 (2019), arXiv:1903.05102 [hep-ph] .
  * Lizarraga  _et al._ [2016] J. Lizarraga, J. Urrestilla, D. Daverio, M. Hindmarsh, and M. Kunz, JCAP 10, 042 (2016), arXiv:1609.03386 [astro-ph.CO] .
  * Auclair  _et al._ [2020] P. Auclair _et al._ , JCAP 04, 034 (2020), arXiv:1909.00819 [astro-ph.CO] .
  * Auclair  _et al._ [2023] P. Auclair, S. Babak, H. Quelquejay Leclere, and D. A. Steer, Phys. Rev. D 108, 043519 (2023), arXiv:2305.11653 [gr-qc] .
  * Afshordi  _et al._ [2023] N. Afshordi _et al._ (LISA Consortium Waveform Working Group), (2023), arXiv:2311.01300 [gr-qc] .
  * Hunter [2007] J. D. Hunter, CSE 9, 90 (2007).
  * Harris  _et al._ [2020] C. R. Harris _et al._ , Nature 585, 357 (2020).
  * Virtanen  _et al._ [2020] P. Virtanen _et al._ , Nature Methods 17, 261 (2020).