# PLUE: Language Understanding Evaluation Benchmark
for Privacy Policies in English

Jianfeng Chi1,2 Wasi Uddin Ahmad3 Yuan Tian3 Kai-Wei Chang3
1Meta AI, 2University of Virginia, 3University of California, Los Angeles
jianfengchi<EMAIL_ADDRESS>Work done while at
UCLA.

###### Abstract

Privacy policies provide individuals with information about their rights and
how their personal information is handled. Natural language understanding
(NLU) technologies can support individuals and practitioners to understand
better privacy practices described in lengthy and complex documents. However,
existing efforts that use NLU technologies are limited by processing the
language in a way exclusive to a single task focusing on certain privacy
practices. To this end, we introduce the Privacy Policy Language Understanding
Evaluation (PLUE) benchmark, a multi-task benchmark for evaluating the privacy
policy language understanding across various tasks. We also collect a large
corpus of privacy policies to enable privacy policy domain-specific language
model pre-training. We demonstrate that domain-specific pre-training offers
performance improvements across all tasks. We release the benchmark to
encourage future research in this domain.111 https://github.com/JFChi/PLUE

##  1 Introduction

Privacy policies are legal documents that outline how an organization
collects, uses, shares, and protects the personal information of individuals.
Without a clear understanding of privacy policies, individuals may not know
how their personal information is being used or who it is being shared with,
and the privacy violation might cause potential harm to them. However, privacy
policies are lengthy and complex, prohibiting users from reading them in
detail (Commission et al., 2012; Gluck et al., 2016; Marotta-Wurgler, 2015).

Various natural language understanding (NLU) technologies have recently been
developed to understand privacy policies, such as text classification (Wilson
et al., 2016a; Harkous et al., 2018), question answering (Ravichander et al.,
2019; Ahmad et al., 2020; Parvez et al., 2022), semantic parsing (Ahmad et
al., 2021), and named-entity recognition (Bui et al., 2021). These tasks focus
on understanding specific privacy practices at different syntax or semantics
levels and require significant effort for data annotations (e.g., domain
experts). It is hard to develop generic pre-trained language models (e.g.,
BERT (Devlin et al., 2019)) with task-specific fine-tuning using limited
annotated data. Besides, the unique characteristics of privacy policies, such
as reasoning over ambiguity and vagueness, modality, and document structure
(Ravichander et al., 2021), make it challenging to directly apply generic pre-
trained language models to the privacy policy domain.

To address these problems and encourage research to develop NLU technologies
in the privacy policy domain, we introduce the Privacy Policy Language
Understanding Evaluation (PLUE) benchmark, a multi-task benchmark for
evaluating the privacy policy language understanding across six tasks,
including text classification, question answering, semantic parsing, and
named-entity recognition. PLUE also includes a pre-training privacy policy
corpus that we crawl from the websites to enable privacy policy domain-
specific language model pre-training. We use this corpus to pre-train BERT
(Devlin et al., 2019), Electra (Clark et al., 2020), and SpanBERT (Joshi et
al., 2020) and fine-tune them for different tasks. We demonstrate that domain-
specific pre-training provides performance improvements across all tasks. We
release the benchmark, including the pre-training codes and corpus and
privacy-policy-related tasks and datasets, to assist natural language
processing (NLP) researchers and practitioners in future exploration of
privacy-policy-related NLP.

##  2 Policy Language Understanding Evaluation (PLUE) Benchmark

Dataset | Task | Sub-domain | |Policy| |  $|$Train$|$ |  $|$Dev$|$ |  $|$Test$|$ | Metric
---|---|---|---|---|---|---|---
OPP-115 | Classification | Websites | 115 | 2,771 | 395 | 625 | F1
APP-350 | Classification | Mobile Apps | 350 | 10,150 | 2,817 | 2,540 | F1
PrivacyQA | QA | Mobile Apps | 35 | 1,350 | – | 400 | P / R / F1
PolicyQA | QA | Mobile Apps | 115 | 17,056 | 3,809 | 4,152 | F1 / EM
PolicyIE |  Intent Classification Slot Filling |  Websites Mobile Apps | 31 | 4,209 | – | 1,041 | F1 / EM
PI-Extract | NER | Websites | 30 | 3,034 | – | 1,028 | F1
Table 1: Statistics of the PLUE datasets and tasks.

PLUE is centered on six English privacy policy language understanding tasks.
The datasets and tasks are selected based on the following principles: (1)
usefulness: the selected tasks can help practitioners in the privacy policy
domain quickly understand particular privacy practices without reading the
whole privacy policy; (2) task diversity: the selected tasks focus on
different semantic levels, e.g., words (phrases), sentences, and paragraphs;
(3) task difficulty: the selected tasks should be adequately challenging for
more room for improvement; (4) training efficiency: all tasks can be trainable
on a single moderate GPU (e.g., GeForce GTX 1080 Ti) for no more than ten
hours; (5) accessibility: all datasets are publicly available.

###  2.1 Datasets and Tasks

PLUE includes six tasks in four categories. Table 1 presents an overview of
the datasets and tasks within PLUE. We describe each task as follows.

#### OPP-115

Wilson et al. (2016a) presented 115 Online Privacy Policies (OPP-115). The
dataset comprises website privacy policies with annotations that specify
privacy practices. Each text segment in privacy policies is annotated with one
or more privacy practices from ten categories (e.g., First Party
Collection/Use, Data Retention, User Access, Edit and Deletion). We train a
multi-label classifier to predict the privacy practices given a sentence from
a policy document.

#### APP-350

Zimmeck et al. (2019) presented APP-350, a collection of mobile application
privacy policies annotating data-type-related privacy practices and could also
be used for text classification. Like OPP-115, each text segment in a privacy
policy is annotated with zero or more privacy practices. Still, APP-350
focuses on what types of data mobile applications collect or share from users.
For example, the privacy practices in APP-350 include Demographic_Age and
Location_GPS (i.e., describes the collection/sharing of the user’s age and GPS
location data). In total, there are 30 data-type-related classes in APP-350,
and we assign one more class, No_Mention, to those text segments that do not
contain any data-type-related classes.

#### PrivacyQA

Ravichander et al. (2019) proposed a question-answering dataset, PrivacyQA,
comprised of 35 mobile application privacy policies. Given a question from a
mobile application user and a sentence from a privacy policy, the task is to
predict whether the sentence is relevant to the question. PrivacyQA includes
unanswerable and subjective questions and formulates the QA task as a binary
classification task.

#### PolicyQA

Ahmad et al. (2020) proposed a reading comprehension (Rajpurkar et al., 2016)
style dataset, PolicyQA. The dataset is derived from OPP-115 annotations that
include a set of fine-grained attributes and evidence text spans that support
the annotations. Considering the annotated spans as the answer spans, PolicyQA
generates diverse questions relating to the corresponding privacy practices
and attributes. The task is to predict the answer text span given the
corresponding text segment and question.

#### PolicyIE

Ahmad et al. (2021) proposed a semantic parsing dataset composed of two tasks:
intent classification and slot filling. Given a sentence in a privacy policy,
the task is to predict the sentence’s intent (i.e., privacy practice) and
identify the semantic concepts associated with the privacy practice. For
example, the sentence “We do not sell your personal information to third
parties” has the intent Data Sharing/Disclosure and slots-value pairs (1) Data
Sharer: First Party Entity–We, (2) Polarity: Negation–not, (3) Action–sell,
(4) Data Shared: General Data–personal information, and (5) Data Receiver:
Third Party Entity–third parties. Based on the role of the slots in privacy
practices, PolicyIE groups them into type-I and type-II slots. In total, there
are four intent labels and 14 type-I and four type-II slot labels. We
individually train a text classifier and sequence taggers to perform intent
classification and slot filling, respectively.

#### PI-Extract

Bui et al. (2021) presented PI-Extract, a named-entity recognition (NER)
dataset. It aims to identify what types of user data are (not) collected or
shared mentioned in the privacy policies. It contains 4 types of named
entities: COLLECT, NOT_COLLECT, SHARE and NOT_SHARE. As the name suggests,
COLLECT and NOT_COLLECT correspond to the types of data that are collected and
not collected mentioned in the privacy policies.

###  2.2 Pre-training Corpus Collection

The existing pre-trained language models (PLMs) mostly use data from
BooksCorpus (Zhu et al., 2015) and English Wikipedia. Language models pre-
trained on text from those sources might not perform well on the downstream
privacy policy language understanding tasks, as privacy policies are composed
of text written by domain experts (e.g., lawyers). Gururangan et al. (2020)
suggested that adapting to the domain’s unlabeled data (domain-adaptive pre-
training) improves the performance of domain-specific tasks. Therefore, we
collect a large privacy policy corpus for language model pre-training. In
order to achieve broad coverage across privacy practices written in privacy
policies (William, 2020; Ahmad et al., 2021), we collect the privacy policies
from two sources: mobile application privacy policies and website privacy
policies.

We use MAPS, the mobile application privacy policy corpus presented by Zimmeck
et al. (2019). MAPS consists of the URLs of 441K mobile application privacy
policies, which were collected from April to May 2018 from the Google Play
store. We remove the duplicated URLs, crawl the privacy policy documents in
HTML/PDF format, convert them to raw text format, and filter out the documents
with noise (e.g., empty documents resulting from obsolete URLs). Finally, we
ended up with 64K privacy policy documents. For website privacy policies, we
use the Princeton-Leuven Longitudinal Corpus of Privacy Policies (Amos et al.,
2021).222The corpus is publicly available at https://github.com/citp/privacy-
policy-historical. The Princeton-Leuven Longitudinal Corpus of Privacy
Policies contains 130K website privacy policies spanning over two decades. We
use the documents with the latest date and convert them (from markdown format)
into text format. Combining these two corpora, we obtain our pre-training
corpus with 332M words.

Models | OPP-115 | APP-350 | PrivacyQA | PolicyQA | PI-Extract
---|---|---|---|---|---
2-6 | F1 | F1 | P / R / F1 | F1 / EM | F1
Human | - | - | 68.8 / 69.0 / 68.9 | - | -
BERT | 75.3 | 59.6 | 44.6 / 35.9 / 36.3 | 55.1 / 27.7 | 62.3
SpanBERT | 62.8 | 32.8 | 24.8 / 24.8 / 24.8 | 55.2 / 27.8 | 61.3
Electra | 74.0 | 49.3 | 42.7 / 36.0 / 36.1 | 57.5 / 29.9 | 68.1
PP-BERT | 78.0 | 62.8 | 44.8 / 36.9 / 37.7 |  58.3 / 30.0 | 68.5
PP-SpanBERT | 78.1 | 61.9 | 43.4 / 36.4 / 36.8 | 55.8 / 27.5 | 65.7
PP-Electra | 73.1 | 57.1 |  48.3 / 38.8 / 39.3 | 58.0 / 30.0 | 66.2
LEGAL-BERT | 76.0 | 57.4 | 45.6 / 37.6 / 38.2 | 55.1 / 27.7 | 66.3

Table 2:  Performance comparison of pre-trained models on text classification,
question answering, and named entity recognition tasks. We fine-tune all the
models three times with different seeds and report average performances. Human
performances are reported from the respective works.

Models | Intent | Slot Filling
---|---|---
Classification | Type-I Slots | Type-II Slots
2-6 | F1 | F1 | EM | F1 | EM
Human | 96.5 | 84.3 | 56.6 | 62.3 | 55.6
BERT | 73.7 | 55.2 | 19.7 | 34.7 | 29.8
SpanBERT | 71.9 | 44.0 | 10.8 | 29.7 | 17.5
Electra | 73.7 | 56.4 | 22.8 | 36.5 | 30.7
PP-BERT | 76.9 | 56.7 | 22.8 | 38.7 | 32.5
PP-SpanBERT | 75.0 | 54.1 | 19.8 | 33.6 | 26.7
PP-Electra | 77.1 | 58.2 | 24.1 | 37.8 | 32.9
LEGAL-BERT | 72.6 | 53.8 | 19.5 | 36.1 | 29.7

Table 3:  Performance comparison of pre-trained models on intent
classification and slot filling tasks (PolicyIE). We fine-tune all the models
three times with different seeds and report average performances. Human
performances are reported from the respective works.

##  3 Models

We benchmark a few pre-trained language models as baselines to facilitate
future works.

#### BERT

Devlin et al. (2019) proposed transformer-based language model pre-trained on
BooksCorpus and English Wikipedia using masked language modeling (MLM) and
next sentence prediction.

#### SpanBERT

Joshi et al. (2020) shares the same architecture and pre-training corpus as
BERT but differs in the pre-training objectives. It extends BERT by masking
contiguous spans instead of single tokens and training the span boundary
representations to predict the masked spans.

#### Electra

Clark et al. (2020) pre-trains a generator and a discriminator on the same
corpus as BERT, where the generator takes a masked text as input and is
trained using the MLM objective. The discriminator takes the predictions from
the generator and detects which tokens are replaced by the generator. After
pre-training, the generator is discarded, and the discriminator is used as the
language model for the downstream tasks.

#### LEGAL-BERT

Chalkidis et al. (2020) pre-trained BERT using 12 GB of the English text (over
351K documents) from several legal fields (e.g., contracts, legislation, court
cases) scraped from publicly available resources. Since privacy policies serve
as special legal documents to protect company and consumers’ privacy rights,
we study LEGAL-BERT’s effectiveness on the PLUE tasks.

#### Domain-specific Training

We continue to train BERT, SpanBERT, and Electra on the pre-training corpus
described in Section 2.2. We refer to them as PP-BERT, PP-SpanBERT, and PP-
Electra, respectively. Table 5 in the Appendix presents the hyperparameters we
use to train the PLMs.

#### Task-specific Fine-tuning

Given the pre-trained language models, we fine-tune the models for each task.
Table 5 in the Appendix lists the hyperparameters for all the downstream
tasks.

##  4 Experiment Results

Tables 3 and 3 present the results for all the experiment models for PLUE
tasks. We observe that the language models (PP-BERT, PP-SpanBERT, PP-Electra)
adapted to the privacy policy domain via continual pre-training outperform the
general language models consistently in all the tasks. In particular, PP-BERT
performs the best for APP-350, PolicyQA, and PI-Extract, while PP-Electra
performs the best for PolicyIE and PrivacyQA.333Ravichander et al. (2019)
reported 39.8% F1 score for BERT model; however, we were able to achieve
36.3%. On the other hand, PP-SpanBERT achieves the best performance for
OPP-115. In contrast, LEGAL-BERT performs comparably or shows moderate
improvements over BERT, SpanBERT, and Electra. This indicates that pre-
training on the general legal corpus does not necessarily help privacy policy
language understanding.

It is interesting to see continual pre-training of the language models using
the privacy policy domain data benefits them differently. For example, the
performance difference between SpanBERT and PP-SpanBERT are most significant
in the text classification tasks (i.e., OPP-115 and APP-350). Lastly, even
though domain-specific pre-training helps boost the performance for all tasks,
the performance of some tasks and datasets (e.g., APP-350, PrivacyQA, slot
filling in PolicyIE) remains low, which indicates much potential for further
work on NLP for the privacy policy domain.

##  5 Related Work

#### Privacy Policy Benchmarks

Usable Privacy Project Sadeh et al. (2013) is the most significant effort to
date, resulting in a large pool of works Wilson et al. (2016a, b); Sathyendra
et al. (2016); Mysore Sathyendra et al. (2017); Bhatia and Breaux (2015);
Bhatia et al. (2016); Hosseini et al. (2016); Zimmeck et al. (2019) to
facilitate the automation of privacy policy analysis. Natural language
processing (NLP) techniques including text alignment (Liu et al., 2014;
Ramanath et al., 2014), text classification (Wilson et al., 2016a; Harkous et
al., 2018; Zimmeck et al., 2019), question answering (QA) (Shvartzshanider et
al., 2018; Harkous et al., 2018; Ravichander et al., 2019; Ahmad et al.,
2020), named entity recognition Bui et al. (2021), and semantic parsing Ahmad
et al. (2021) have been explored. In this work, we aggregate some of these
efforts under one umbrella and benchmark pre-trained language models.

#### Pre-trained Language Models

In the last few years, natural language processing (NLP) research has
witnessed a radical change with the advent of pre-trained language models like
ELMo (Peters et al., 2018) and BERT (Devlin et al., 2019). These pre-trained
language models achieved state-of-the-art results in many language
understanding benchmarks. In particular, PLMs found to be very effective in
low-resource settings. Consequently, PLMs have been developed for NLP in a
wide range of domains, giving birth to SciBERT Beltagy et al. (2019) for
scientific NLP, BioBERT Lee et al. (2020), Med-BERT Rasmy et al. (2021), and
ClinicalBERT Alsentzer et al. (2019) for medical NLP, LEGAL-BERT Chalkidis et
al. (2020) for legal NLP, CyBERT Ranade et al. (2021) and CySecBERT Bayer et
al. (2022) for Cybersecurity domain. This work particularly investigates the
adaptation of PLMs to facilitate NLP research in the privacy policy domain.

##  6 Conclusion

Reliable aggregation of datasets benefits NLP research for low-resource
domains. Moreover, benchmarking foundation models on these datasets could
facilitate future research. In this work, we present PLUE, a benchmark to
facilitate developing, evaluating, and comparing new security and privacy
policy models. Our contributions will help privacy and NLP communities
benchmark language understanding tasks more easily under a unified setup.

## Limitations

The pre-training privacy policy corpus and the downstream task datasets are
unlikely to contain toxic or biased content. Therefore, they should not
magnify toxicity or bias in the pre-trained and fine-tuned models, although
the models may exhibit such behavior due to their original pre-training. The
pre-training and benchmark datasets are formed based on privacy policies
crawled in the past; as a result, they could be outdated by now. This work
focuses on the English language only, and the findings may not apply to other
languages.

## Ethics Statement

#### License

The OPP115 and APP350 datasets are made available for research, teaching, and
scholarship purposes only, with further parameters in the spirit of a Creative
Commons Attribution-NonCommercial License (CC BY-NC). The PolicyQA and PI-
Extract datasets are derived from OPP115 datasets. The PrivacyQA and PolicyIE
dataset are released under an MIT license. The pre-training corpus, MAPS
Policies Dataset, is released under CC BY-NC. We strictly adhere to these
licenses and will release the PLUE benchmark resources under CC BY-NC-SA 4.0.

#### Carbon Footprint

We avoided using large models for continual training on the privacy policy
domain, reducing their environmental impacts. The PP-BERT, PP-SpanBERT, and
PP-Electra models were trained for 100k steps on Tesla V100 GPUs that took 1-2
days. Therefore, the training would emit only 9kg of carbon into the
environment.444Calculated using https://mlco2.github.io/impact, based on a
total of 100 hours of training on Tesla V100 and Amazon Web Services as the
provider. All fine-tuning experiments were very lightweight due to the small
size of the datasets, resulting in approximately 12kg of carbon emission.

## References

  * Ahmad et al. (2021) Wasi Ahmad, Jianfeng Chi, Tu Le, Thomas Norton, Yuan Tian, and Kai-Wei Chang. 2021\.  Intent classification and slot filling for privacy policies.  In _Proceedings of the 59th Annual Meeting of the Association for Computational Linguistics and the 11th International Joint Conference on Natural Language Processing (Volume 1: Long Papers)_ , pages 4402–4417, Online. Association for Computational Linguistics.
  * Ahmad et al. (2020) Wasi Ahmad, Jianfeng Chi, Yuan Tian, and Kai-Wei Chang. 2020.  PolicyQA: A reading comprehension dataset for privacy policies.  In _Findings of the Association for Computational Linguistics: EMNLP 2020_ , pages 743–749, Online. Association for Computational Linguistics.
  * Alsentzer et al. (2019) Emily Alsentzer, John Murphy, William Boag, Wei-Hung Weng, Di Jindi, Tristan Naumann, and Matthew McDermott. 2019.  Publicly available clinical BERT embeddings.  In _Proceedings of the 2nd Clinical Natural Language Processing Workshop_ , pages 72–78, Minneapolis, Minnesota, USA. Association for Computational Linguistics.
  * Amos et al. (2021) Ryan Amos, Gunes Acar, Eli Lucherini, Mihir Kshirsagar, Arvind Narayanan, and Jonathan Mayer. 2021.  Privacy policies over time: Curation and analysis of a million-document dataset.  In _Proceedings of the Web Conference 2021_ , WWW ’21, page 2165–2176, New York, NY, USA. Association for Computing Machinery.
  * Bayer et al. (2022) Markus Bayer, Philipp Kuehn, Ramin Shanehsaz, and Christian Reuter. 2022.  Cysecbert: A domain-adapted language model for the cybersecurity domain.  _arXiv preprint arXiv:2212.02974_.
  * Beltagy et al. (2019) Iz Beltagy, Kyle Lo, and Arman Cohan. 2019.  SciBERT: A pretrained language model for scientific text.  In _Proceedings of the 2019 Conference on Empirical Methods in Natural Language Processing and the 9th International Joint Conference on Natural Language Processing (EMNLP-IJCNLP)_ , pages 3615–3620, Hong Kong, China. Association for Computational Linguistics.
  * Bhatia and Breaux (2015) Jaspreet Bhatia and Travis D Breaux. 2015.  Towards an information type lexicon for privacy policies.  In _2015 IEEE eighth international workshop on requirements engineering and law (RELAW)_ , pages 19–24. IEEE.
  * Bhatia et al. (2016) Jaspreet Bhatia, Morgan C Evans, Sudarshan Wadkar, and Travis D Breaux. 2016.  Automated extraction of regulated information types using hyponymy relations.  In _2016 IEEE 24th International Requirements Engineering Conference Workshops (REW)_ , pages 19–25. IEEE.
  * Bui et al. (2021) Duc Bui, Kang G. Shin, Jong-Min Choi, and Junbum Shin. 2021.  Automated extraction and presentation of data practices in privacy policies.  _Proceedings on Privacy Enhancing Technologies_ , 2021(2):88–110.
  * Chalkidis et al. (2020) Ilias Chalkidis, Manos Fergadiotis, Prodromos Malakasiotis, Nikolaos Aletras, and Ion Androutsopoulos. 2020.  LEGAL-BERT: The muppets straight out of law school.  In _Findings of the Association for Computational Linguistics: EMNLP 2020_ , pages 2898–2904, Online. Association for Computational Linguistics.
  * Clark et al. (2020) Kevin Clark, Minh-Thang Luong, Quoc V Le, and Christopher D Manning. 2020.  Electra: Pre-training text encoders as discriminators rather than generators.  In _International Conference on Learning Representations_.
  * Commission et al. (2012) Federal Trade Commission et al. 2012.  Protecting consumer privacy in an era of rapid change.  _FTC report_.
  * Devlin et al. (2019) Jacob Devlin, Ming-Wei Chang, Kenton Lee, and Kristina Toutanova. 2019.  BERT: Pre-training of deep bidirectional transformers for language understanding.  In _Proceedings of the 2019 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, Volume 1 (Long and Short Papers)_ , pages 4171–4186, Minneapolis, Minnesota. Association for Computational Linguistics.
  * Gluck et al. (2016) Joshua Gluck, Florian Schaub, Amy Friedman, Hana Habib, Norman Sadeh, Lorrie Faith Cranor, and Yuvraj Agarwal. 2016.  How short is too short? implications of length and framing on the effectiveness of privacy notices.  In _Twelfth Symposium on Usable Privacy and Security ( $\\{$SOUPS$\\}$ 2016)_, pages 321–340.
  * Gururangan et al. (2020) Suchin Gururangan, Ana Marasović, Swabha Swayamdipta, Kyle Lo, Iz Beltagy, Doug Downey, and Noah A. Smith. 2020.  Don’t stop pretraining: Adapt language models to domains and tasks.  In _Proceedings of the 58th Annual Meeting of the Association for Computational Linguistics_ , pages 8342–8360, Online. Association for Computational Linguistics.
  * Harkous et al. (2018) Hamza Harkous, Kassem Fawaz, Rémi Lebret, Florian Schaub, Kang G Shin, and Karl Aberer. 2018.  Polisis: Automated analysis and presentation of privacy policies using deep learning.  In _27th $\\{$USENIX$\\}$ Security Symposium ($\\{$USENIX$\\}$ Security 18)_, pages 531–548.
  * Hosseini et al. (2016) Mitra Bokaei Hosseini, Sudarshan Wadkar, Travis D Breaux, and Jianwei Niu. 2016\.  Lexical similarity of information type hypernyms, meronyms and synonyms in privacy policies.  In _2016 AAAI Fall Symposium Series_.
  * Joshi et al. (2020) Mandar Joshi, Danqi Chen, Yinhan Liu, Daniel S. Weld, Luke Zettlemoyer, and Omer Levy. 2020.  SpanBERT: Improving pre-training by representing and predicting spans.  _Transactions of the Association for Computational Linguistics_ , 8:64–77.
  * Lee et al. (2020) Jinhyuk Lee, Wonjin Yoon, Sungdong Kim, Donghyeon Kim, Sunkyu Kim, Chan Ho So, and Jaewoo Kang. 2020.  Biobert: a pre-trained biomedical language representation model for biomedical text mining.  _Bioinformatics_ , 36(4):1234–1240.
  * Liu et al. (2014) Fei Liu, Rohan Ramanath, Norman Sadeh, and Noah A. Smith. 2014.  A step towards usable privacy policy: Automatic alignment of privacy statements.  In _Proceedings of COLING 2014, the 25th International Conference on Computational Linguistics: Technical Papers_ , pages 884–894, Dublin, Ireland. Dublin City University and Association for Computational Linguistics.
  * Marotta-Wurgler (2015) Florencia Marotta-Wurgler. 2015.  Does “notice and choice” disclosure regulation work? an empirical study of privacy policies,”.  In _Michigan Law: Law and Economics Workshop_.
  * Mysore Sathyendra et al. (2017) Kanthashree Mysore Sathyendra, Shomir Wilson, Florian Schaub, Sebastian Zimmeck, and Norman Sadeh. 2017.  Identifying the provision of choices in privacy policy text.  In _Proceedings of the 2017 Conference on Empirical Methods in Natural Language Processing_ , pages 2774–2779, Copenhagen, Denmark. Association for Computational Linguistics.
  * Parvez et al. (2022) Md Rizwan Parvez, Jianfeng Chi, Wasi Uddin Ahmad, Yuan Tian, and Kai-Wei Chang. 2022\.  Retrieval enhanced data augmentation for question answering on privacy policies.  _arXiv preprint arXiv:2204.08952_.
  * Peters et al. (2018) Matthew E. Peters, Mark Neumann, Mohit Iyyer, Matt Gardner, Christopher Clark, Kenton Lee, and Luke Zettlemoyer. 2018.  Deep contextualized word representations.  In _Proceedings of the 2018 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, Volume 1 (Long Papers)_ , pages 2227–2237, New Orleans, Louisiana. Association for Computational Linguistics.
  * Rajpurkar et al. (2016) Pranav Rajpurkar, Jian Zhang, Konstantin Lopyrev, and Percy Liang. 2016.  SQuAD: 100,000+ questions for machine comprehension of text.  In _Proceedings of the 2016 Conference on Empirical Methods in Natural Language Processing_ , pages 2383–2392, Austin, Texas. Association for Computational Linguistics.
  * Ramanath et al. (2014) Rohan Ramanath, Fei Liu, Norman Sadeh, and Noah A. Smith. 2014.  Unsupervised alignment of privacy policies using hidden Markov models.  In _Proceedings of the 52nd Annual Meeting of the Association for Computational Linguistics (Volume 2: Short Papers)_ , pages 605–610, Baltimore, Maryland. Association for Computational Linguistics.
  * Ranade et al. (2021) Priyanka Ranade, Aritran Piplai, Anupam Joshi, and Tim Finin. 2021.  Cybert: Contextualized embeddings for the cybersecurity domain.  In _2021 IEEE International Conference on Big Data (Big Data)_ , pages 3334–3342. IEEE.
  * Rasmy et al. (2021) Laila Rasmy, Yang Xiang, Ziqian Xie, Cui Tao, and Degui Zhi. 2021.  Med-bert: pretrained contextualized embeddings on large-scale structured electronic health records for disease prediction.  _NPJ digital medicine_ , 4(1):1–13.
  * Ravichander et al. (2021) Abhilasha Ravichander, Alan W Black, Thomas Norton, Shomir Wilson, and Norman Sadeh. 2021.  Breaking down walls of text: How can NLP benefit consumer privacy? In _Proceedings of the 59th Annual Meeting of the Association for Computational Linguistics and the 11th International Joint Conference on Natural Language Processing (Volume 1: Long Papers)_ , pages 4125–4140, Online. Association for Computational Linguistics.
  * Ravichander et al. (2019) Abhilasha Ravichander, Alan W Black, Shomir Wilson, Thomas Norton, and Norman Sadeh. 2019.  Question answering for privacy policies: Combining computational and legal perspectives.  In _Proceedings of the 2019 Conference on Empirical Methods in Natural Language Processing and the 9th International Joint Conference on Natural Language Processing (EMNLP-IJCNLP)_ , pages 4947–4958, Hong Kong, China. Association for Computational Linguistics.
  * Sadeh et al. (2013) Norman Sadeh, Alessandro Acquisti, Travis D Breaux, Lorrie Faith Cranor, Aleecia M McDonald, Joel R Reidenberg, Noah A Smith, Fei Liu, N Cameron Russell, Florian Schaub, et al. 2013.  The usable privacy policy project.  _Technical report, Technical Report, CMU-ISR-13-119_.
  * Sathyendra et al. (2016) Kanthashree Mysore Sathyendra, Florian Schaub, Shomir Wilson, and Norman Sadeh. 2016\.  Automatic extraction of opt-out choices from privacy policies.  In _2016 AAAI Fall Symposium Series_.
  * Shvartzshanider et al. (2018) Yan Shvartzshanider, Ananth Balashankar, Thomas Wies, and Lakshminarayanan Subramanian. 2018.  RECIPE: Applying open domain question answering to privacy policies.  In _Proceedings of the Workshop on Machine Reading for Question Answering_ , pages 71–77, Melbourne, Australia. Association for Computational Linguistics.
  * William (2020) Henry William. 2020.  Do web apps and mobile apps need separate privacy policies?
  * Wilson et al. (2016a) Shomir Wilson, Florian Schaub, Aswarth Abhilash Dara, Frederick Liu, Sushain Cherivirala, Pedro Giovanni Leon, Mads Schaarup Andersen, Sebastian Zimmeck, Kanthashree Mysore Sathyendra, N. Cameron Russell, Thomas B. Norton, Eduard Hovy, Joel Reidenberg, and Norman Sadeh. 2016a.  The creation and analysis of a website privacy policy corpus.  In _Proceedings of the 54th Annual Meeting of the Association for Computational Linguistics (Volume 1: Long Papers)_ , pages 1330–1340, Berlin, Germany. Association for Computational Linguistics.
  * Wilson et al. (2016b) Shomir Wilson, Florian Schaub, Rohan Ramanath, Norman Sadeh, Fei Liu, Noah A Smith, and Frederick Liu. 2016b.  Crowdsourcing annotations for websites’ privacy policies: Can it really work?  In _Proceedings of the 25th International Conference on World Wide Web_ , pages 133–143. International World Wide Web Conferences Steering Committee.
  * Zhu et al. (2015) Yukun Zhu, Ryan Kiros, Rich Zemel, Ruslan Salakhutdinov, Raquel Urtasun, Antonio Torralba, and Sanja Fidler. 2015.  Aligning books and movies: Towards story-like visual explanations by watching movies and reading books.  In _Proceedings of the IEEE international conference on computer vision_ , pages 19–27.
  * Zimmeck et al. (2019) Sebastian Zimmeck, Peter Story, Daniel Smullen, Abhilasha Ravichander, Ziqi Wang, Joel Reidenberg, N Cameron Russell, and Norman Sadeh. 2019.  Maps: Scaling privacy compliance analysis to a million apps.  _Proceedings on Privacy Enhancing Technologies_ , 2019(3):66–86.

| PP-BERT | PP-SpanBERT | PP-Electra
---|---|---|---
Learning Rate | 1e-4 | 1e-4 | 1e-4
Train Steps | 100,000 | 100,000 | 100,000
batch Size | 256 | 256 | 256
Learning Rate Schedule | linear | polynomial_decay | linear
# warm-up steps | 1000 | 1000 | 1000
Optimizer | AdamW | AdamW | AdamW

Table 4: Hyperparameters for pre-training language models.

|  Text Classification |  Question Answering |  Semantic Parsing |  NER
---|---|---|---|---
Dropout | 0.1 | 0.1 | 0.1 | 0.1
Weight decay | 0.0 | 0.0 | 0.0 | 0.0
Optimizer | AdamW | AdamW | AdamW | AdamW
Batch Size | 32 | 32 | 32 | 32
Learning rate | 3e-5 | 2e-5 | 3e-5 | 5e-5
Learning Rate Schedule | Linear | Linear | Linear | Linear
Warm-up Ratio | 0.05 | 0.0 | 0.05 | 0.05
# epoch | 20 | 3 | 20 | 20

Table 5: Hyperparameters for fine-tuning pre-trained language models on
different PLUE tasks.


# Where binary neutron stars merge: predictions from IllustrisTNG

Jonah C. Rose Department of Astronomy, University of Florida, Gainesville, FL
32611, USA  Paul Torrey Department of Astronomy, University of Florida,
Gainesville, FL 32611, USA  K.H. Lee Department of Physics, University of
Florida, Gainesville, FL 32611, USA  I. Bartos Department of Physics,
University of Florida, Gainesville, FL 32611, USA

(Received June 1, 2019; Revised January 10, 2019)

###### Abstract

The rate and location of Binary Neutron Star (BNS) mergers are determined by a
combination of the star formation history and the Delay Time Distribution
(DTD) function. In this paper, we couple the star formation rate histories
(SFRHs) from the IllustrisTNG model to a series of varied assumptions for the
BNS DTD to make predictions for the BNS merger host galaxy mass function.
These predictions offer two outcomes: (i) in the near term: influence BNS
merger event follow-up strategy by scrutinizing where most BNS merger events
are expected to occur and (ii) in the long term: constrain the DTD for BNS
merger events once the host galaxy mass function is observationally well
determined. From our fiducial model analysis, we predict that 50% of BNS
mergers will occur in host galaxies with stellar mass between
$10^{10}-10^{11}$ $M_{\odot}$, 68% between $4\times 10^{9}-3\times 10^{11}$
$M_{\odot}$, and 95% between $4\times 10^{8}-2\times 10^{12}$ $M_{\odot}$. We
find that the details of the DTD employed does not have a strong effect on the
peak of the host mass function. However, varying the DTD provides enough
spread that the true DTD can be determined from enough electromagnetic
observations of BNS mergers. Knowing the true DTD can help us determine the
prevalence of BNS systems formed through highly eccentric and short separation
fast-merging channels and can constrain the dominant source of r-process
material.

neutron star merger, gravitational waves, methods: numerical, stars: neutron,
binaries: close

††journal: AJ

##  1 Introduction

Since the first discovery of gravitational waves by LIGO (Aasi et al., 2015),
a growing number of compact object mergers have been detected. To date, two
detections have been confirmed as binary neutron star (BNS) mergers (Abbott et
al., 2017a, 2020). Of these two, the BNS event GW170817 was detected across
the electromagnetic spectrum (Abbott et al., 2017b), beginning the age of
multi-messenger astronomy.

Within the next few years we expect tens of new BNS mergers events to be
observed by LIGO, Virgo (Acernese et al., 2014) and KAGRA (Akutsu et al.,
2019), broadening our understanding of BNS systems and their host galaxies
Abbott et al. (2018). Developing a clearer understanding of the link between
BNS merger events and their host galaxies is useful for multiple reasons in
both the short and long term.

In the short term, having clear predictions for the BNS host galaxy mass
function could inform follow-up strategies for future BNS merger events
detected by LIGO/Virgo/KAGRA. Locating the electromagnetic counterpart of GW
events is difficult owing to narrowly peaked observability windows and large
localization areas (Smartt et al., 2017; Metzger & Berger, 2012). GW
localizations can extend tens-to-hundreds of square degrees, making them
impractical to completely cover in a reasonable time after the initial event
with most telescopes (Gehrels et al., 2016; Bartos et al., 2013, 2015, 2018,
2019a). Long-term radio emission could allow sufficient time for follow-up
observations, but this will only be possible for nearby events within dense
circum-merger media (Bartos et al., 2019b; Lee et al., 2020; Grandorf et al.,
2020).

GW follow-up strategies have taken two approaches to search the localization
area more efficiently: covering the entire area or targeting galaxies (e.g.
Bartos et al., 2014; Arcavi et al., 2017; Chan et al., 2018; Antier et al.,
2020). Covering the entire localization area increases the likelihood of
imaging the correct host galaxy, but risks missing the transient owing to the
limited exposure times. Soares-Santos et al. (2017) used this method to
successfully locate the kilonova after GW170817 by covering 80.7% of the
probability weighted localization area. In contrast, targeted follow-ups use
galaxy catalogs to preferentially search galaxies based on select criteria
(e.g. galaxy blue luminosity), potentially reducing the required number of
pointings by a factor of 10 to 100 and increasing exposure times (Gehrels et
al., 2016; Ducoin et al., 2020). This strategy was used in the first
successful detection of the optical counterpart of GW170817 (Coulter et al.,
2017). However, it is more likely that targeted strategies will miss the event
if the BNS merger takes place in a less-massive galaxy. Efficient follow-up
strategies are important for maximising the chance of identifying the
electromagnetic counterpart with limited observations. One way to achieve this
is to build a clearer understanding of the expected host galaxy mass function
for BNS mergers.

Longer term, the link between BNS merger events and their host galaxies can
help determine the dominant formation channel of $r$-process material by
constraining the true BNS Delay Time Distribution (DTD) (Marchant et al.,
2016; Barrett et al., 2018; Mapelli et al., 2019; Santoliquido et al., 2020;
McCarthy et al., 2020). Specifically, while core-collapse supernovae (SNe) and
BNS mergers have been proposed as $r$-process formation channels, the dominant
$r$-process production channel must be able to recreate the observed
decreasing trend in Eu/Fe vs Fe/H (Matteucci et al., 2015). BNS mergers must
produce $r$-process elements in less than 100 Myr to dominate $r$-process
element production (Hotokezaka et al., 2018; Côté et al., 2017), and in less
than 1 Myr – with a steep cutoff slope – to be the source of all $r$-process
material (Matteucci et al., 2015). While it is possible for BNS systems to
merge in this time (Safarzadeh et al., 2019b), they require highly eccentric
orbits from high-velocity kicks or low initial separation from case BB mass
transfer, both of which may not occur in BNS formation (Tauris et al., 2017).
These models also predict a shallower DTD based on the current understanding
of BNS formation channels (Giacobbo & Mapelli, 2019; Simonetti et al., 2019;
Safarzadeh & Berger, 2019).

Current stellar populations synthesis models suggest BNS DTDs are best
modelled by power law distributions with an exponent between -1 and -1.5
(Simonetti et al., 2019) and the minimum time from creation of the binary
system to merger ($t_{\mathrm{min}}$) between 1 Myr to 1 Gyr (e.g. Simonetti
et al., 2019; Safarzadeh et al., 2019b). These models assert that the main
formation channel that forms BNS systems begin with two OB stars that are
close enough to undergo mass transfer (Tauris et al., 2017; Giacobbo &
Mapelli, 2018; Safarzadeh et al., 2019b). The rest of the systems are born
through so-called fast-merging channels where a binary system forms with
either a high eccentricity through large natal kicks, or with a small initial
separation through unstable case BB mass transfer (Tauris et al., 2017). If
the BNS DTD could be observationally constrained, it would not only shed light
on the physical formation channels (fast-merging vs. OB star mass transfer),
but could also help constrain progenitor metallicity, common-envelope
efficiency, natal kicks, mass ratio, and initial binary separation through
comparisons with population synthesis codes (Giacobbo & Mapelli, 2018;
Belczynski et al., 2018; Dominik et al., 2012). Constraining the BNS DTD may
be one of the best and most practical ways to constrain the physical origin
and implications of BNS merger events.

Taken over a whole galaxy, the rate of BNS merger events can be determined by
convolving the DTD with the star formation rate history (SFRH). Metallicity is
also accounted for in some models, but has been found to play a minor role in
influencing the DTD (Mapelli et al., 2019; Giacobbo & Mapelli, 2019; Côté et
al., 2017; Giacobbo & Mapelli, 2018). Previous studies have measured the BNS
merger rate using SFRHs derived from EAGLE and Illustris cosmological
simulations (Artale et al., 2019; Mapelli et al., 2018, 2019), the FIRE zoom-
in simulation (Lamberts et al., 2018), dark matter only simulations (Marassi
et al., 2019; Cao et al., 2018), or from semi-analytical models (Adhikari et
al., 2020; Toffano et al., 2019) with population synthesis codes or an assumed
DTD. Each of these models provide a different SFRH for the galaxies in that
simulation which provide different distributions for BNS mergers given the
same DTD.

In this paper, we use the IllustrisTNG simulations to make predictions for the
BNS merger host galaxy mass function. This extends the work that which has
been presented in (Mapelli et al., 2018) and (Artale et al., 2019) by focusing
on the uncertainty/variation introduced by variations in the assumed DTD.
Moreover, we consider here how in the future an observed BNS host galaxy mass
function could be used to constrain the real/underlying DTD. To do this, we
take as input the IllustrisTNG galactic SFRHs – which are known to match a
wide range of observed galaxy properties and galaxy scaling relations – and
employ varied assumptions about the BNS DTD. Our chosen SFRHs and DTDs allow
us to demonstrate the galactic masses at which we expect most BNS mergers to
occur, as well as to identify the level of variation that would be induced
based on changes to the BNS DTD.

The structure of this paper is as follows. In Section §2 we outline our
methods including a brief description of the IllustrisTNG simulations, our
adopted DTDs, and our methods for calculating the galaxy-by-galaxy BNS merger
rate. In §3 we present our main results including predictions for the BNS
merger host galaxy mass function and the sensitivity of this prediction to the
assumed DTD. In §4 we discuss the implications of our results. Finally, in §5
we discuss our results and conclude.

##  2 Methods

In this paper, we make predictions for the BNS merger host galaxy mass
function by adopting SFRHs from cosmological simulations and DTDs from basic
stellar population synthesis models.

###  2.1 Delay Time Distributions

Generally, the current BNS merger rate for any galaxy is given by convolving
its SFRH with the appropriate DTD. The BNS merger rate for any collection of
material (e.g. galaxy, volume, etc.) is given by

$r(t)=\int_{0}^{t}\psi(\tau)\Gamma(t-\tau,\,Z)d\tau$ (1)

where $\psi$ is the star formation rate, $\Gamma(t-\tau,\,Z)$ is the DTD, and
the integration is performed from the Big Bang ($t=0$) to the time of
observation, $t$ (e.g. Maoz et al., 2012). The metallicity dependence, $Z$, in
$\Gamma$ is only present in some DTDs, otherwise the DTD takes the form
$\Gamma(t-\tau)$. In the case of cosmological galaxy formation simulations,
this can be reduced to a sum over contributions from all relevant stellar
populations

$r(t)=\sum_{j}M_{j}\Gamma(t-t_{j},\,Z_{j})$ (2)

where the sum is performed over all stars (or stellar populations) in the
region of interest (e.g. within a specific galaxy), $M_{j}$ is the mass of
each stellar particle or stellar population, $t_{j}$ is birth time of that
stellar population such that $t-t_{j}$ is the age of the stellar population,
and $Z_{j}$ is the metallicity of the stellar population. For any cosmological
galaxy formation simulation, the BNS merger rate is easily evaluated once a
BNS DTD is specified.

Figure 1: (left) The two fiducial DTDs shown as the BNS merger rate vs time.
The BPASS DTD is split into two lines to show the range covered by the DTD as
the metallicity of the host star changes. (right) The normalized BNS merger
rate for the two fiducial DTDs as a function of the host galaxy stellar mass.
The solid line shows the merger rate given the power law DTD with an exponent
of s=-1 and $t_{\rm cut}$=0.01Gyr. The dashed line shows the merger rate given
the BPASS DTD. Both merger rates have been normalized individually such that
the total merger rate across the simulation for the given DTD is unity. The
shaded bands show the mass range which contain 50, 68, and 95 percent of the
mergers around the peak merger rate for the fiducial power law DTD. The arrow
points to the host galaxy mass of the only BNS merger with a detected
electromagnetic counterpart so far.

###  2.2 Power Law DTD

In this paper, we adopt two fiducial DTDs: (i) a simple paramaterized power
law and (ii) the metalicity dependent DTD from BPASS (Eldridge & Stanway,
2016). The power law DTD is given by

$\Gamma(t-t_{j})=\begin{cases}0&t-t_{j}\leq t_{\mathrm{cut}}\\\
\Gamma_{0}(t-t_{j})^{s}&t-t_{j}>t_{\mathrm{cut}}\end{cases}$ (3)

where $\Gamma_{0}$ is a normalization coefficient, $s$ is the power law index,
and $t_{\mathrm{cut}}$ is the minimum time/age before the first BNS merger
event occurs. Figure 1 shows our fiducial power law DTD ($s=-1$;
$t_{\mathrm{cut}}=10^{7}\,\mathrm{yrs}$). In addition to our fiducial power
law DTD, we also consider DTDs that have varied power law exponents ranging
from $s=-2$ to $s=2$ and cutoff times ranging from
$t_{\mathrm{cut}}=0.001\,\mathrm{Gyr}$ to $t_{\mathrm{cut}}=10\,\mathrm{Gyr}$.

###  2.3 BPASS DTD

In addition to using a simple power law DTD, we adopt our second fiducial DTD
from BPASS (Eldridge & Stanway, 2016) as shown in Figure 1. BPASS is a stellar
population synthesis code which follows the evolution of a large suite of
varying binary stars (Eldridge et al., 2008; Eldridge & Stanway, 2016). The
critical feature of BPASS that is important for this paper is that the
simulated stellar population matches the observed population of binaries in
abundance along with supernovae progenitors and rates (Eldridge et al., 2008,
2013, 2015; Eldridge & Stanway, 2016). The remnants of supernova can have a
mass in the full range between .1 and 300 $M_{\odot}$, allowing for more
realistic evolution of these systems. The simulation also encompass a wide
range of stellar metalicities which has been shown to potentially affect
(albeit weakly) the DTD of BNS mergers (e.g. Mapelli et al., 2018). The final
time for the BNSs to merge in BPASS is then the sum of the progenitor stars’
evolution and the in-spiral time once the BNS system has formed (Eldridge &
Stanway, 2016). The tabulated BPASS BNS merger rates are a function of stellar
population age and metallicity, $\Gamma(t-t_{j},Z_{j})$, which can be employed
in conjunction with Equation 2 to determine the total BNS merger rate.

We note that the overall normalization for all of our DTDs ($\Gamma_{0}$, in
the case of the powerlaw DTD) can be specified to match the expected global
rate of BNS mergers. However, its exact value is not important to the present
work as we are only interested in the relative/normalized distribution of BNS
merger events as a function of galaxy mass. We therefore normalize the total
BNS merger rate across the entire simulation box to unity for each DTD
individually. To achieve this, we divide the rate for an individual galaxy,
$j$, by the rate of the entire box for the given DTD, $\Gamma$. The normalized
form of Equation 2 becomes

$R_{i}(t)=\frac{r_{i}(t)}{\sum_{k}M_{k}\Gamma(t-t_{k},\,Z_{k})}$ (4)

where $r_{i}$ represents the BNS merger rate for an individual galaxy given by
equation 2, and the denominator sums over the rates of all galaxies ($k$) in
the simulation box.

We linearly interpolate across the ages and metallicities presented in
Eldridge & Stanway (2016). We do not extrapolate outside of the provided
metallicity values; all stars with $Z\leq 0.0001$ follow the DTD for
$Z=0.0001$ and all stars with $Z\geq 0.014$ follow the DTD for $Z=0.014$. The
BNS merger rate is then calculated using Equation 4. The BPASS DTD includes
information on natal kicks from the initial supernovae. For each supernova,
the kick velocity and direction are determined from Hobbs et al. (2005). For
more information, see Eldridge et al. (2011). There are no natal kicks
included in the calculation of the power law delay time distributions,
including the fiducial power law model. Thus, the power law DTDs are fully
specified with two parameters controlling (i) the time of BNS mergers onset
and (ii) the subsequent BNS merger rate evolution.

###  2.4 The IllustrisTNG Simulation Suite

In order to calculate the BNS merger rates, we adopt SFRHs from IllustrisTNG
simulation (Pillepich et al., 2018; Nelson et al., 2018a; Marinacci et al.,
2018; Springel et al., 2018; Naiman et al., 2018). IllustrisTNG is a suite of
cosmological hydrodnamical simulations which includes a comprehensive galaxy
formation model (Pillepich et al., 2018; Weinberger et al., 2017) and builds
upon the original Illustris model (Vogelsberger et al., 2013; Torrey et al.,
2014). The critical feature of the IllustrisTNG simulations for this paper is
that the simulations have been shown to broadly reproduce the cosmic star
formation rate density and the redshift-dependent galaxy stellar mass function
(Pillepich et al., 2018; Springel et al., 2018; Nelson et al., 2018b). These
tests give confidence that both the global and galaxy-by-galaxy star formation
histories produced by IllustrisTNG reasonably match those of the Universe.

While it is possible that the constraints for any given galaxy differ from a
“real” galaxy, when averaged over a large enough sample the trends are well
matched. Specifically, we employ the TNG-100-1 simulation which includes a
$100\;\rm{Mpc}$ cubed volume with hundreds-of-thousands of galaxies with
varied SFRHs and self-consistently evolved metallicity distributions. All
stellar particles are used in equation 4 to calculate the normalization, but
results are only presented for galaxies down to a stellar mass of
$10^{7}M_{\odot}$. We use the full information from the simulated galaxy
populations including the age and metallicity distribution to evaluate the BNS
merger rate.

##  3 Results

Figure 1 shows the DTDs (left) and host galaxy mass function (right) for our
two fiducial setups. The DTDs follow the same general form with mergers being
most likely shortly after $t_{\mathrm{cut}}$ then dropping with increasing
time. The main difference comes from the metallicity dependence in the BPASS
DTD. Owing to the similarities in the DTDs, the resulting BNS merger host
galaxy mass functions are remarkably similar. In particular, we find that the
peak of BNS mergers occurs in host galaxies with stellar masses around
$M_{*}=5\times 10^{10}$ $M_{\odot}$ and that half of the mergers occur in
galaxies with masses between $10^{10}M_{\odot}<M_{*}<10^{11}$ $M_{\odot}$.
Despite the significant added complexity of the BPASS models, the predicted
BNS host galaxy mass function is not significantly different from the power
law DTD. Additionally, the mass of the host galaxy from GW170817 (Blanchard et
al., 2017) is indicated with a downward facing arrow in the right panel of
Figure 1. While it is a sample size of one and should not be over-interpreted,
we note that GW170817’s host galaxy mass is in the peak region of expected BNS
host galaxy masses for both of our fiducial DTD models.

Figure 2: The individually normalized BNS merger rates as a function of
stellar mass for varied power law exponents (left) and varied
$t_{\mathrm{cut}}$ values (right). There is significant variation in the
predicted host galaxy mass functions when the DTD is perturbed from the
fiducial values.

Figure 2 shows the host galaxy mass functions for the power law DTDs with
varied power law exponents and $t_{\mathrm{cut}}$. For completeness, we
explore a large range of values for the power law exponents and cutoff times
$t_{\mathrm{cut}}$ which go beyond what is believed to be physically correct
(Côté et al., 2017). These values are included to demonstrate the variation in
host galaxy mass functions which would result from varied DTDs. Even with this
large spread of DTDs, we find that the results shown in Figure 1 are broadly
stable. In particular, despite the very significant variation in the DTDs, all
cases show a peak BNS merger rate that occurs in galaxies with a stellar mass
in the range $10^{10}$-$10^{11}$ $M_{\odot}$.

The DTD for which $s=0$ is of particular interest in Figure 2 because it
tracks the total stellar mass found in each mass bin – nearly independent of
star formation history.111There is a dependence on the amount of stellar mass
that formed in the past $t_{\mathrm{cut}}=10\mathrm{Myrs}$, but this is
expected to be only a $\sim 0.1\%$ correction. Owing to the shape of the
simulated (and observed) galaxy stellar mass functions, the peak of the
stellar mass distribution is in galaxies with stellar masses between $10^{10}$
and $10^{11}$ $M_{\odot}$. Therefore, the majority of BNS mergers for this DTD
are also found in that mass range. Importantly, because the predicted host
galaxy mass function for a DTD with power law index of $s=0$ is nearly
independent of galactic formation history, this result is not very sensitive
to the detailed SFRHs predicted by IllustrisTNG, but only the shape of the
galaxy stellar mass function. Insofar as other simulations or models reproduce
the same galaxy stellar mass functions, their predicted host galaxy mass
function for $s=0$ would be nearly identical.

Changing the power law exponent to values away from $s=0$ introduces a direct
dependence on the assumed star formation history by placing emphasis either on
the older or younger stellar populations. Specifically, power law exponents
higher than $s=0$ lead to systematic changes in which the host galaxy mass
function is biased toward galaxies with older stellar masses. This naturally
results in a shift of the peak in the host galaxy mass function toward older,
more massive systems. Conversely, changing the power law exponent to values
lower than $s=0$ (which is the more physical case) biases the host galaxy mass
function toward systems with younger stellar populations. Thus, as the power
law exponent is decreased, there is an expectation that an increasing number
of BNS mergers occur in low mass galaxies with current or recent ongoing star
formation. For a fixed DTD, the detailed shape that we predict for the BNS
merger host galaxy mass function is dependent on the IllustrisTNG galaxy
stellar mass function and SFRHs, and therefore should be checked against other
models (e.g. EAGLE, SIMBA, etc.). However, owing to observational constraints
provided by the cosmic star formation rate density and redshift dependent
galaxy stellar mass functions, we do not expect these results to substantially
change.

Despite the stability in the the peak of the host mass function across
different DTDs, there is still significant spread in the resulting host mass
functions at other masses. For example, when examining the BNS merger rate in
galaxies with a host mass near $10^{9}$, the BNS merger rates differ by 1-dex
between the $s=2$ and $s=-2$ DTDs with a fixed $t_{\mathrm{cut}}$. There is an
even greater spread in the highest mass systems where the merger rate differs
by 2-dex between the $s=2$ and $s=-2$ DTDs. A similar range in host mass
functions occurs across the different $t_{\mathrm{cut}}$ values at a fixed
value of $s$. The predicted variability in host galaxy mass functions suggests
that as GW BNS detections with EM follow-up observations mount, a careful
comparisons of observed and predicted host galaxy mass functions could be used
to constrain the true DTD for BNS mergers. These results agree with Safarzadeh
& Berger (2019) and running KS-tests on our power law host mass functions also
results in $\mathcal{O}(1000)$ observations being required to determine a true
DTD.

When optimizing BNS merger event follow-up, a question arises of which
galaxies should be targeted first. While Figures 1 and 2 indicate that most
BNS mergers will occur in roughly Milky Way mass galaxies, Figure 3 shows the
host mass function normalized by the number of galaxies in each mass bin which
indicates the predicted number of BNS mergers per galaxy. Here, the host mass
function no longer peaks in the $10^{10}$-$10^{11}$ $M_{\odot}$ range, but
instead peaks at larger masses (in the $10^{12}$-$10^{13}$ $M_{\odot}$ range).
This indicates that while our analysis predicts that most BNS merges will
occur in $\sim$ Milky Way mass galaxies, the highest likelihood of finding a
BNS merger based on a single observation still favors more massive systems,
simply because they have more mass. This conclusion is somewhat dependent on
the detailed assigned prescription for the BNS merger DTD. In particular,
while the fiducial values (see the magenta line in the left panel of Figure 3)
clearly peaks at the highest mass bin resolved in the IllustrisTNG volume, the
steeper exponent cases ($s=-1.5$ and $s=-2$) are much flatter above
$M_{*}=10^{10.5}M_{\odot}$.

A closer examination of how the BNS merger rate correlates with different
observables indicates a dependence on the DTD. This analysis was conducted by
comparing, through the Pearson correlation coefficient, the BNS merger rate
for each galaxy to one of three observable properties: its star formation rate
(SFR), its blue luminosity, and its stellar mass. For very steep and negative
DTDs, $s=-2$, we find that SFR is best correlated with merger rate (R=0.978).
For less steep and negative DTDs, $s=-1$, we find that blue luminosity is best
correlated (R=0.993). For flat or increasing DTDs, $s=0+$ we find that stellar
mass is best correlated (R=1.0, 0.994, 0.985 for s=0,1,2 respectively).

Figure 3:  The individually normalized BNS merger rates as a function of
stellar mass per galaxy for varied power law exponents (left) and varied
$t_{\mathrm{cut}}$ values (right). While the host galaxy mass function (Figure
2) predicts most BNS mergers will happen in roughly Milky Way mass galaxies
when averaged over the whole galaxy population, this host galaxy specific-mass
function (this figure) indicates that the rate of BNS merger rate is higher is
higher in higher mass galaxies, when compared on an individual basis.

##  4 Discussion

The ability to connect LIGO-detected BNS merger events to their host galaxy
opens new scientific opportunities. Specifically, while transient event
detection and host galaxy association is well-established in astronomy,
traditional methods for kilonova detection yield little direct information
about the progenitor system. In contrast, wave form fitting of LIGO detected
compact object merger events provides detailed information about the
progenitor system including the masses of the merging objects. This new
information links mergers of specific object types to host galaxies with a
limited level of ambiguity or uncertainty that was not previously possible. As
we have discussed in this paper, this opens up the possibility of developing a
more intimate link between galactic star formation rate histories (SFRHs), BNS
delay time distributions, and the observed host galaxy stellar mass function.

In this paper, we have leveraged the galactic SFRHs from the IllustrisTNG
cosmological galaxy formation model. At some level, these SFRHs are likely not
a perfect reflection of real galactic SFRHs. However, the model is able to
broadly match the cosmic star formation rate history as well as the redshift
dependent galaxy stellar mass function. This gives us a reasonable level of
confidence that these simulated SFRHs provide good approximations to those of
real galaxies. Moreover, since fairly different physical models (e.g. those of
Illustris model, the EAGLE model, and semi-analytical models) yield similar
BNS merger rates, we believe the more holistic analysis obtained from the
TNG-100 simulation can further our understanding of the host galaxy mass
function and its dependence on the DTD.

An important feature of the IllustrisTNG simulation is that it self-
consistently tracks gas- and stellar-phase metallicities. The stellar
metallicities, in turn, impact the results of the metallicity-dependent DTDs,
such as the BPASS DTD. It has been shown that IllustrisTNG’s stellar
metallicity vs stellar mass relation generally agrees with observations but is
too flat leading to higher metallicities at lower masses (Nelson et al.,
2018b). The largest discrepancy is $\sim 0.5$ dex near $10^{10.5}\;M_{\odot}$.
To understand how this uncertainty affects our results, we use two host galaxy
mass functions with the BPASS DTD. Each host galaxy mass function is made by
setting all stellar metallicities to either $Z=0.0001$ or $Z=0.014$. Given the
large difference in metallicity, $\sim 2.5$ dex, between these host galaxy
mass functions, we expect any differences to be larger than those introduced
from uncertainties in the IllustrisTNG stellar metallicities. We find very
little variation between the host galaxy mass functions between the lowest and
highest metallicities across all galaxy masses. Given the large spread in
metallicities used in creating these DTDs, it is unlikely that our results
would be changed significantly if we instead employed stellar metallicities
from a different galaxy formation model. Given that our results do not change
when accounting for the more complicated BPASS DTD, additional credibility can
be given to results derived from power law DTDs.

Some efforts have already begun to explore the new connection between merger
events and their host galaxy using population synthesis codes and various star
formation histories. These studies investigate how the host galaxy mass
function is affected by the merger progenitors. The unique contribution of
this paper is to focus on systematic variations of the employed DTD coupled to
cosmologically motivated SFRHs. A similar set of DTD variations was employed
in Safarzadeh & Berger (2019), albeit with analytically simplified SFRHs. They
concluded that the host galaxy mass function peaks at high masses and that
$\mathcal{O}(10^{3})$ observations are required to constrain the a true power
law distribution. This paper agrees with their conclusions for the host galaxy
mass range that they cover, $10^{9}-10^{11.25}M_{\odot}$. This study extends
Safarzadeh & Berger (2019) by pairing a broad set of assumed DTDs to SFRHs
naturally derived in a cosmological environment and examining how the assumed
DTD affects the host galaxy mass function over the large range of host masses
allowed by IllustrisTNG.

Similar results to those presented in this paper have also been discussed in
Artale et al. (2019); Safarzadeh et al. (2019a); McCarthy et al. (2020).
Artale et al. (2019) uses the EAGLE simulation to create a stellar mass vs
specific BNS merger rate plot similar to Figure 3. They find that stellar mass
is an excellent tracer for specific merger rate. This result is consistent
with our result up to $\sim 10^{10.5}M_{\odot}$. However, at higher masses we
find a dependence on the DTD, causing faster merging times to not depend on
stellar mass.

Safarzadeh et al. (2019a); Adhikari et al. (2020); McCarthy et al. (2020) also
present host galaxy mass functions using different SFRH models. Safarzadeh &
Berger (2019) uses an analytic model with the set us DTDs used in this paper
to understand how the host galaxy mass function is affected by the DTD. Our
results are generally consistent with theirs, but our extended range of host
masses allow us to see that most BNS mergers do not happen in galaxies with
highest mass, but in galaxies with masses between $10^{10}$ and
$10^{11}M_{\odot}$. McCarthy et al. (2020) also uses an analytic model but
paired with SDSS observations to explore the host galaxy mass function along
with other host observables. For the host galaxy mass function, our results
are consistent with theirs. However, our larger range of DTDs presented in
Figure 2 reveal the large spread between different assumed DTDs and the stable
peak near $10^{10.5}M_{\odot}$. Adhikari et al. (2020) also find that other
host observables paired with stellar mass are necessary to obtain a better
understanding where BNS merge. Overall, we find that future explorations of
this topic will need to consider a wide range of DTDs and the full range of
the observable in question.

While the work we present here continues our understanding of what we can
learn from observations of BNS host galaxies, further investigations are
necessary to fully understand how BNS form and evolve. One example of such an
investigation is to expand the set of DTDs examined using the methods in this
paper. The set of DTDs we examine are broad, covering those most commonly
referenced (e.g. Safarzadeh & Berger, 2019; Eldridge & Stanway, 2016), but we
do not exhaustively search the full range of DTDs proposed (e.g. Simonetti et
al., 2019; Dominik et al., 2012). Also, our convolution of IllustrisTNG’s SFRH
with our DTDs does not include any form of natal kicks. If these kicks are
strong enough to dislodge the binary from smaller galaxies, it is possible
their addition would weight the host mass functions toward higher mass
galaxies. With a greater range of DTDs examined and a more detailed
convolution, we will gain a clearer picture of where BNS mergers are located,
which delay times can be distinguished using the host galaxy mass function,
and the most likely places they will be observed. Another way to incorporate a
more complete set of DTDs would be to use a varied set of population synthesis
models which cover a wide range of binary separations, kick velocities,
initial mass functions, etc.

Including other star formation histories could also provide a more detailed
look at the spread in possible host galaxy mass functions. While IllustrisTNG
is broadly consistent with the cosmic star formation rate density and redshfit
dependent galaxy stellar mass functions (Pillepich et al., 2018), its accuracy
should not be over interpreted and different simulations will surely produce
somewhat varied star formation histories that could impact our results.
However, we can say that up to their mass cutoff, our results align with
Artale et al. (2019) who found no significant difference when comparing
results from Illustris and Eagle. The lack of variation in Artale et al.
(2019) most likely indicates that – while there is some variation – the SFRHs
in Illustris and EAGLE are sufficiently similar to not significantly impact
the results. Thus, by adopting the SFRHs from galaxy formation simulations and
assumptions about the functional form of the DTD, predictions can be made
about the BNS host galaxy mass function. Additionally, similarities between
the different simulations suggest that uncertainties in the poorly constrained
DTDs are likely larger than the uncertainties introduced from the SFRHs. In
particular, the detailed shape of the BNS host galaxy mass function will be
sensitive to assumptions about the DTD.

##  5 Conclusion

We presented predictions for the host galaxy mass function and host galaxy
specific-mass function for BNS mergers. Our predictions were generated by
convolving a set of power law and BPASS DTDs with the star formation histories
from the IllustrisTNG cosmological simulation. Our main conclusions are as
follows:

  1. 1.

We find almost no difference between the host galaxy mass functions produced
by our fiducial power law (slope of $s=-1$, minimum time of $t_{\rm
cut}=10^{7}$ yrs) and the BPASS DTDs (Figure 1).

  2. 2.

The peak of the host galaxy mass function occurs around the Milky Way mass
scale, with roughly $\sim 50\%$ of BNS mergers happening in the
$10^{10}M_{\odot}<M_{*}<10^{11}M_{\odot}$ mass range for our fiducial DTDs
(Figure 1). This mass bin includes NGC 4993, the host galaxy of GW170817.

  3. 3.

While the detailed shape of the host galaxy mass function is sensitive to
details of the adopted DTD, the peak does not change significantly when
varying over a broad range of DTDs (Figure 2). The peak of the host galaxy
specific-mass function is similarly insensitive to changes in the adopted DTD
(Figure 3).

  4. 4.

The peak of the host galaxy specific-mass function is located in the highest
mass bin for the fiducial power law DTD and BPASS model. Thus, while we expect
most BNS mergers to happen in somewhat lower mass systems for our fiducial
DTDs, high mass galaxies are more likely to host a BNS merger on a per-galaxy
basis (Figures 1 and 3).

  5. 5.

Host galaxy mass functions constructed from different DTDs vary up to one dex
at low masses and up to two dex at high masses. This provides an opportunity
through which an observationally reconstructed host galaxy mass function can
be used to constrain the true BNS DTD (Figure 2).

  6. 6.

The observable galactic property (or properties) that is expected to provide
the best correlation with the BNS merger rate depends on the true DTD.

In the short term, the results found here in both the host mass and specific-
mass functions paint an interesting picture on how astronomers should
structure electromagnetic follow-ups for BNS events. The peak of the host
galaxy specific-mass function laying in the highest mass bin suggests that the
optimal way to quickly find the resulting kilonova from a BNS merger would be
to search the highest mass galaxies first. This agrees with the current method
most follow-up strategies use in locating BNS mergers (e.g. Gehrels et al.,
2016; Arcavi et al., 2017; Singer et al., 2016). However, the peak of the host
mass function laying in the mass range $M_{*}=10^{10}-10^{11}\;M_{\odot}$
suggests that this method will miss, or take a longer to locate, most of the
BNS mergers. Determining the true DTD would allow for more efficient
electromagnetic follow-up by determining which observable: SFR, blue
luminosity, or stellar mass, best correlates with BNS merger rate.

In the long term, LIGO/Virgo/KAGRA will create a host mass function which can
be used to determine the true BNS DTD. With this DTD, the minimum delay time,
$t_{\rm cut}$, can constrain the proportion of BNS systems which form through
highly eccentric and low separation fast-merging channels. Understanding this
proportion will place constraints on natal kick velocity and common envelop
efficiency. The minimum delay time can also determine whether BNS mergers are
the dominant source of r-process elements. The overall shape of the true DTD
allows various physical parameters of BNS systems to be constrained, such as
the progenitor’s metallicity, masses, mass ratio, common envelope efficiency,
natal kicks, and initial binary separation through comparisons with resulting
DTDs from population synthesis codes.

## Acknowledgements

The authors thank Steve Eikenberry for his useful ideas and comments. JCR
acknowledges support from the University of Florida Graduate School’s Graduate
Research Fellowship. PT acknowledges support from NSF grant AST-1909933, NASA
ATP Grant 19-ATP19-0031. IB acknowledges support from the Alfred P. Sloan
Foundation and the University of Florida.

## References

  * Aasi et al. (2015) Aasi, J., et al. 2015, Class. Quantum Grav., 32, 074001, doi: 10.1088/0264-9381/32/7/074001
  * Abbott et al. (2017a) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2017a, Phys. Rev. Lett., 119, 161101, doi: 10.1103/PhysRevLett.119.161101
  * Abbott et al. (2017b) —. 2017b, ApJ, 848, L12, doi: 10.3847/2041-8213/aa91c9
  * Abbott et al. (2018) Abbott, B. P., et al. 2018, Living Rev. Relativ., 21, 3, doi: 10.1007/s41114-018-0012-9
  * Abbott et al. (2020) Abbott, B. P., Abbott, R., Abbott, T. D., et al. 2020, ApJ, 892, L3, doi: 10.3847/2041-8213/ab75f5
  * Acernese et al. (2014) Acernese, F., Agathos, M., Agatsuma, K., et al. 2014, Classical and Quantum Gravity, 32, 024001, doi: 10.1088/0264-9381/32/2/024001
  * Adhikari et al. (2020) Adhikari, S., Fishbach, M., Holz, D. E., Wechsler, R. H., & Fang, Z. 2020, arXiv e-prints, arXiv:2001.01025.  https://arxiv.org/abs/2001.01025
  * Akutsu et al. (2019) Akutsu, T., Ando, M., Arai, K., et al. 2019, Nature Astronomy, 3, 35–40, doi: 10.1038/s41550-018-0658-y
  * Antier et al. (2020) Antier, S., Agayeva, S., Aivazyan, V., et al. 2020, MNRAS, 492, 3904, doi: 10.1093/mnras/stz3142
  * Arcavi et al. (2017) Arcavi, I., McCully, C., Hosseinzadeh, G., et al. 2017, ApJ, 848, L33, doi: 10.3847/2041-8213/aa910f
  * Artale et al. (2019) Artale, M. C., Mapelli, M., Giacobbo, N., et al. 2019, MNRAS, 487, 1675, doi: 10.1093/mnras/stz1382
  * Barrett et al. (2018) Barrett, J. W., Gaebel, S. M., Neijssel, C. J., et al. 2018, MNRAS, 477, 4685, doi: 10.1093/mnras/sty908
  * Bartos et al. (2013) Bartos, I., Brady, P., & Márka, S. 2013, Class. Quantum Grav., 30, 123001, doi: 10.1088/0264-9381/30/12/123001
  * Bartos et al. (2019a) Bartos, I., Corley, K. R., Gupte, N., et al. 2019a, MNRAS, 490, 3476, doi: 10.1093/mnras/stz2848
  * Bartos et al. (2015) Bartos, I., Crotts, A. P. S., & Márka, S. 2015, ApJ, 801, L1, doi: 10.1088/2041-8205/801/1/L1
  * Bartos et al. (2019b) Bartos, I., Lee, K. H., Corsi, A., Márka, Z., & Márka, S. 2019b, MNRAS, 485, 4150, doi: 10.1093/mnras/stz719
  * Bartos et al. (2014) Bartos, I., Veres, P., Nieto, D., et al. 2014, MNRAS, 443, 738, doi: 10.1093/mnras/stu1205
  * Bartos et al. (2018) Bartos, I., Di Girolamo, T., Gair, J. R., et al. 2018, MNRAS, 477, 639, doi: 10.1093/mnras/sty602
  * Belczynski et al. (2018) Belczynski, K., Bulik, T., Olejak, A., et al. 2018, arXiv e-prints, arXiv:1812.10065.  https://arxiv.org/abs/1812.10065
  * Blanchard et al. (2017) Blanchard, P. K., Berger, E., Fong, W., et al. 2017, ApJ, 848, L22, doi: 10.3847/2041-8213/aa9055
  * Cao et al. (2018) Cao, L., Lu, Y., & Zhao, Y. 2018, MNRAS, 474, 4997, doi: 10.1093/mnras/stx3087
  * Chan et al. (2018) Chan, M. L., Messenger, C., Heng, I. S., & Hendry, M. 2018, Phys. Rev. D, 97, 123014, doi: 10.1103/PhysRevD.97.123014
  * Côté et al. (2017) Côté, B., Belczynski, K., Fryer, C. L., et al. 2017, ApJ, 836, 230, doi: 10.3847/1538-4357/aa5c8d
  * Coulter et al. (2017) Coulter, D. A., Foley, R. J., Kilpatrick, C. D., et al. 2017, Science, 358, 1556, doi: 10.1126/science.aap9811
  * Dominik et al. (2012) Dominik, M., Belczynski, K., Fryer, C., et al. 2012, ApJ, 759, 52, doi: 10.1088/0004-637X/759/1/52
  * Ducoin et al. (2020) Ducoin, J. G., Corre, D., Leroy, N., & Le Floch, E. 2020, MNRAS, 492, 4768, doi: 10.1093/mnras/staa114
  * Eldridge et al. (2015) Eldridge, J. J., Fraser, M., Maund, J. R., & Smartt, S. J. 2015, MNRAS, 446, 2689, doi: 10.1093/mnras/stu2197
  * Eldridge et al. (2013) Eldridge, J. J., Fraser, M., Smartt, S. J., Maund, J. R., & Crockett, R. M. 2013, MNRAS, 436, 774, doi: 10.1093/mnras/stt1612
  * Eldridge et al. (2008) Eldridge, J. J., Izzard, R. G., & Tout, C. A. 2008, MNRAS, 384, 1109, doi: 10.1111/j.1365-2966.2007.12738.x
  * Eldridge et al. (2011) Eldridge, J. J., Langer, N., & Tout, C. A. 2011, MNRAS, 414, 3501, doi: 10.1111/j.1365-2966.2011.18650.x
  * Eldridge & Stanway (2016) Eldridge, J. J., & Stanway, E. R. 2016, MNRAS, 462, 3302, doi: 10.1093/mnras/stw1772
  * Gehrels et al. (2016) Gehrels, N., Cannizzo, J. K., Kanner, J., et al. 2016, ApJ, 820, 136, doi: 10.3847/0004-637X/820/2/136
  * Giacobbo & Mapelli (2018) Giacobbo, N., & Mapelli, M. 2018, MNRAS, 480, 2011, doi: 10.1093/mnras/sty1999
  * Giacobbo & Mapelli (2019) —. 2019, MNRAS, 482, 2234, doi: 10.1093/mnras/sty2848
  * Grandorf et al. (2020) Grandorf, C., McCarty, J., Rajkumar, P., et al. 2020, arXiv e-prints, arXiv:2008.05330.  https://arxiv.org/abs/2008.05330
  * Hobbs et al. (2005) Hobbs, G., Lorimer, D. R., Lyne, A. G., & Kramer, M. 2005, MNRAS, 360, 974, doi: 10.1111/j.1365-2966.2005.09087.x
  * Hotokezaka et al. (2018) Hotokezaka, K., Beniamini, P., & Piran, T. 2018, International Journal of Modern Physics D, 27, 1842005, doi: 10.1142/S0218271818420051
  * Lamberts et al. (2018) Lamberts, A., Garrison-Kimmel, S., Hopkins, P. F., et al. 2018, MNRAS, 480, 2704, doi: 10.1093/mnras/sty2035
  * Lee et al. (2020) Lee, K. H., Bartos, I., Privon, G. C., Rose, J. C., & Torrey, P. 2020, arXiv e-prints, arXiv:2007.00563.  https://arxiv.org/abs/2007.00563
  * Maoz et al. (2012) Maoz, D., Mannucci, F., & Brandt, T. D. 2012, MNRAS, 426, 3282, doi: 10.1111/j.1365-2966.2012.21871.x
  * Mapelli et al. (2019) Mapelli, M., Giacobbo, N., Santoliquido, F., & Artale, M. C. 2019, MNRAS, 487, 2, doi: 10.1093/mnras/stz1150
  * Mapelli et al. (2018) Mapelli, M., Giacobbo, N., Toffano, M., et al. 2018, MNRAS, 481, 5324, doi: 10.1093/mnras/sty2663
  * Marassi et al. (2019) Marassi, S., Graziani, L., Ginolfi, M., et al. 2019, MNRAS, 484, 3219, doi: 10.1093/mnras/stz170
  * Marchant et al. (2016) Marchant, P., Langer, N., Podsiadlowski, P., Tauris, T. M., & Moriya, T. J. 2016, A&A, 588, A50, doi: 10.1051/0004-6361/201628133
  * Marinacci et al. (2018) Marinacci, F., Vogelsberger, M., Pakmor, R., et al. 2018, MNRAS, 480, 5113, doi: 10.1093/mnras/sty2206
  * Matteucci et al. (2015) Matteucci, F., Romano, D., Arcones, A., Korobkin, O., & Rosswog, S. 2015, MNRAS, 447, 326, doi: 10.1093/mnras/stu2463
  * McCarthy et al. (2020) McCarthy, K. S., Zheng, Z., & Ramirez-Ruiz, E. 2020, MNRAS, 499, 5220, doi: 10.1093/mnras/staa3206
  * Metzger & Berger (2012) Metzger, B. D., & Berger, E. 2012, ApJ, 746, 48, doi: 10.1088/0004-637X/746/1/48
  * Naiman et al. (2018) Naiman, J. P., Pillepich, A., Springel, V., et al. 2018, MNRAS, 477, 1206, doi: 10.1093/mnras/sty618
  * Nelson et al. (2018a) Nelson, D., Pillepich, A., Springel, V., et al. 2018a, MNRAS, 475, 624, doi: 10.1093/mnras/stx3040
  * Nelson et al. (2018b) —. 2018b, MNRAS, 475, 624, doi: 10.1093/mnras/stx3040
  * Pillepich et al. (2018) Pillepich, A., Nelson, D., Hernquist, L., et al. 2018, MNRAS, 475, 648, doi: 10.1093/mnras/stx3112
  * Pillepich et al. (2018) Pillepich, A., Springel, V., Nelson, D., et al. 2018, MNRAS, 473, 4077
  * Safarzadeh & Berger (2019) Safarzadeh, M., & Berger, E. 2019, ApJ, 878, L12, doi: 10.3847/2041-8213/ab24df
  * Safarzadeh et al. (2019a) Safarzadeh, M., Berger, E., Leja, J., & Speagle, J. S. 2019a, ApJ, 878, L14, doi: 10.3847/2041-8213/ab24e3
  * Safarzadeh et al. (2019b) Safarzadeh, M., Ramirez-Ruiz, E., Andrews, J. J., et al. 2019b, ApJ, 872, 105, doi: 10.3847/1538-4357/aafe0e
  * Santoliquido et al. (2020) Santoliquido, F., Mapelli, M., Bouffanais, Y., et al. 2020, arXiv e-prints, arXiv:2004.09533.  https://arxiv.org/abs/2004.09533
  * Simonetti et al. (2019) Simonetti, P., Matteucci, F., Greggio, L., & Cescutti, G. 2019, MNRAS, 486, 2896, doi: 10.1093/mnras/stz991
  * Singer et al. (2016) Singer, L. P., Chen, H.-Y., Holz, D. E., et al. 2016, ApJ, 829, L15, doi: 10.3847/2041-8205/829/1/L15
  * Smartt et al. (2017) Smartt, S. J., Chen, T. W., Jerkstrand, A., et al. 2017, Nature, 551, 75, doi: 10.1038/nature24303
  * Soares-Santos et al. (2017) Soares-Santos, M., Holz, D. E., Annis, J., et al. 2017, ApJ, 848, L16, doi: 10.3847/2041-8213/aa9059
  * Springel et al. (2018) Springel, V., Pakmor, R., Pillepich, A., et al. 2018, MNRAS, 475, 676, doi: 10.1093/mnras/stx3304
  * Tauris et al. (2017) Tauris, T. M., Kramer, M., Freire, P. C. C., et al. 2017, ApJ, 846, 170, doi: 10.3847/1538-4357/aa7e89
  * Toffano et al. (2019) Toffano, M., Mapelli, M., Giacobbo, N., Artale, M. C., & Ghirlanda, G. 2019, MNRAS, 489, 4622, doi: 10.1093/mnras/stz2415
  * Torrey et al. (2014) Torrey, P., Vogelsberger, M., Genel, S., et al. 2014, MNRAS, 438, 1985, doi: 10.1093/mnras/stt2295
  * Vogelsberger et al. (2013) Vogelsberger, M., Genel, S., Sijacki, D., et al. 2013, MNRAS, 436, 3031, doi: 10.1093/mnras/stt1789
  * Weinberger et al. (2017) Weinberger, R., Springel, V., Hernquist, L., et al. 2017, MNRAS, 465, 3291, doi: 10.1093/mnras/stw2944


# RelativisticDynamics.jl: Relativistic Spin-Orbital Dynamics in Julia

Tom Kimpson

††margin:  DOI: 10.21105/joss.04992 Software • Review • Repository • Archive
Submitted: 15 November 2022
Published: 10 October 2023 License
Authors of papers retain copyright and release the work under a Creative
Commons Attribution 4.0 International License (CC BY 4.0).

## Summary

Relativistic binaries composed of a millisecond pulsar (MSP) orbiting a much
more massive ($\gtrsim 10^{3}M_{\odot}$), spinning black hole (BH) are
exceptional probes for investigating key questions of fundamental physics and
astrophysics. Such systems are natural sources of gravitational waves (GWs) in
the mHz regime, expected to be detectable by the next generation of space-
based GW detectors such as LISA (Thorpe et al., 2019). The associated radio
emission from the companion pulsar raises the possibility of an
electromagnetic (EM) counterpart, enabling high precision multimessenger
measurements to be made. The description of the orbital dynamics of these
systems, and the influence on the resultant observed EM and GW signals, is
non-trivial. A proper treatment of the spin-orbital dynamics can be derived
from the conservation of the energy-momentum tensor

${T^{\mu\nu}}_{;\nu}=0$ (1)

which when expanded into a set of infinite multipole moments leads to a
description of the momentum vector $p^{\mu}$ and the spin tensor $s^{\mu\nu}$

$\frac{Dp^{\mu}}{d\lambda}=-\frac{1}{2}{R^{\mu}}_{\nu\alpha\beta}u^{\nu}s^{\alpha\beta}$
(2) $\frac{Ds^{\mu\nu}}{d\lambda}=p^{\mu}u^{\nu}-p^{\nu}u^{\mu}$ (3)

for affine parameter $\lambda$, 4-velocity $u^{\nu}$ and Riemann curvature
tensor ${R^{\mu}}_{\nu\alpha\beta}$. The system is closed by providing a spin
supplementary condition, equivalent to specifying the observer-dependent
centre of mass. For this work we take the Tulczyjew-Dixon condition (Dixon,
1964; Tulczyjew, 1959)

$s^{\mu\nu}p_{\nu}=0$ (4)

Together, equations 2 \- 4 form the Mathisson-Papetrou-Dixon (MPD) equations
(Dixon, 1964; Mathisson, 1937; Papapetrou, 1951), and describe the spin-
orbital evolution in a fully consistent way that is applicable to strong field
regimes.

## Statement of need

RelativisticDynamics.jl is an open-source Julia package for relativistic spin-
orbital dynamics in the gravitational strong field for a Kerr spacetime.
Existing codes for modelling the dynamics of spinning objects like pulsars in
the strong-field regime are generally lacking, since such systems occupy an
intermediate regime that is generally overlooked. At the “low” end of this
regime there are post-Newtonian or geodesic descriptions (e.g. Damour &
Taylor, 1992) which neglect the influence of the pulsar spin on the underlying
spacetime metric (“spin-curvature” coupling). At the “high” end there is the
full Numerical Relativity (NR) solutions (e.g. Andrade et al., 2021) which are
primarily applicable to two BHs with a mass ratio $\mathcal{O}(1)$, and are
computationally intractable for these MSP systems which are observed over a
large number of orbital cycles.

RelativisticDynamics.jl aims to bridge this gap by providing a modern, fast
code for accurate numerical evolution of spinning relativistic systems, via
the MPD formalism. Julia is a modern language that solves the “two language
problem”, enabling fast dynamic typing and JIT compilation in conjunction with
petaflop performance, comparable with numerical languages that are better
known in the astrophysics community such as C or Fortran. As a modern
language, it also provides a dedicated package manager and a large catalogue
of _composable_ packages for scientific computing. This enables
RelativisticDynamics.jl to easily leverage and interface with other scientific
computing packages. The author and collaborators have used the general methods
and mathematics described in this package for multiple research projects (e.g.
Kimpson et al., 2019, 2020a, 2020b; Li et al., 2019) with a particular focus
on the radio signals from spinning pulsar systems. This package represents an
attempt to create a documented, well-tested, open source resource for public
use in this area, that can also be used as a computational playground for
exploring techniques that could be applicable to more advanced numerical
models. The package has been formulated in terms of ODE integration, rather
than using e.g. action-angle variables (Witzany, 2022), to allow for extension
to general spacetime metrics and straightforward computation of quantities
relevant for pulsar observations e.g. spin axis orientation.

In addition to providing a fast, modern package for strong field spin
dynamics, RelativisticDynamics.jl has two additional important features from
the perspective of modern relativistic astrophysics. Firstly, it is fully type
flexible, being able to support arbitrary number formats. By making use of
Julia’s type-flexibility the model is written in such a way so as to be able
to support hardware accelerated, low precision arithmetic and alternative
rounding methods such as stochastic rounding. This enables rapid prototyping
and exploration of reduced precision numerical techniques in astrophysics, an
approach common in other numerical fields such as weather and climate
modelling (e.g. Váňa et al., 2017). Secondly, RelativisticDynamcis.jl is
written to be fully differentiable via automatic differentiation. This enables
the package to be used for differentiable physics applications in
astrophysics, for example gravitational waveform modelling and parameter
estimation or training neural networks based on the model. Automatic
differentiation also provides a potential avenue for extension of the package
to general (i.e. non-Kerr) spacetimes, whereby a user can specify the metric
and the associated Christoffel symbols and Riemann tensors - which are simply
linear combinations of the metric derivatives - are calculated automatically.

Future potential extensions of this code include taking the dynamics beyond
second order in the multipole expansion, and the inclusion of alternative spin
conditions and spacetime metrics. The inclusion of a diagnostics tool for
extracting gravitational waveforms in the time domain via a numerical kludge
method would also be a worthwhile addition. Moreover, we have considered only
bound dynamical systems - the ability to also explore hyberbolic systems would
also be an interesting development.

Figure 1: Example orbital trajectories for a ms-pulsar with eccentricity
$e=0.1$ (left panels), $e=0.8$ (right panels), orbiting a massive BH with
extremal spin, $a=0.998$. The orbital motion is presented in the $x-y$ plane
(top panels) and $x-z$ plane (bottom panels). The pulsar is initialised in the
orbital plane with zero inclination. In the absence of spin-curvature coupling
the particle would remain in the plane ($z=0$). Note the $z$-motion is on the
scale of km, not gravitational radii.

## Acknowledgements

This work exploring the spin-evolution of relativistic systems via the the MPD
equations was originally motivated through interesting discussions with Kinwah
Wu. The port to a modern, precision-flexible model in Julia was heavily
inspired by Milan Klöwer. Our thanks to both.

## References

reAndrade, T., Salo, L. A., Aurrekoetxea, J. C., Bamber, J., Clough, K.,
Croft, R., Jong, E. de, Drew, A., Duran, A., Ferreira, P. G., Figueras, P.,
Finkel, H., França, T., Ge, B.-X., Gu, C., Helfer, T., Jäykkä, J., Joana, C.,
Kunesch, M., … Wong, K. (2021). GRChombo: An adaptable numerical relativity
code for fundamental physics. _Journal of Open Source Software_ , _6_(68),
3703. https://doi.org/10.21105/joss.03703

preDamour, T., & Taylor, J. H. (1992). Strong-field tests of relativistic
gravity and binary pulsars. _Phys. Rev. D_ , _45_ , 1840–1868.
https://doi.org/10.1103/PhysRevD.45.1840

preDixon, W. G. (1964). A covariant multipole formalism for extended test
bodies in general relativity. _Il Nuovo Cimento_ , _34_(2), 317–339.
https://doi.org/10.1007/BF02734579

preKimpson, T., Wu, K., & Zane, S. (2019). Pulsar timing in extreme mass ratio
binaries: a general relativistic approach. _486_(1), 360–377.
https://doi.org/10.1093/mnras/stz845

preKimpson, T., Wu, K., & Zane, S. (2020a). Orbital spin dynamics of a
millisecond pulsar around a massive BH with a general mass quadrupole.
_Monthly Notices of the RAS_ , _497_(4), 5421–5431.
https://doi.org/10.1093/mnras/staa2103

preKimpson, T., Wu, K., & Zane, S. (2020b). Radio timing in a millisecond
pulsar - extreme/intermediate mass ratio binary system. _Astronomy and
Astrophysics_ , _644_ , A167. https://doi.org/10.1051/0004-6361/202038561

preLi, K. J., Wu, K., & Singh, D. (2019). Spin dynamics of a millisecond
pulsar orbiting closely around a massive black hole. _Monthly Notices of the
Royal Astronomical Society_ , _485_(1), 1053–1066.
https://doi.org/10.1093/mnras/stz389

preMathisson, A. (1937). _Acta Phys. Pol._ , _6_ , 163.
https://doi.org/10.1007/s10714-010-0939-y

prePapapetrou, A. (1951). Spinning test-particles in general relativity. i.
_Proceedings of the Royal Society of London A: Mathematical, Physical and
Engineering Sciences_ , _209_(1097), 248–258.
https://doi.org/10.1098/rspa.1951.0200

preThorpe, J. I., Ziemer, J., Thorpe, I., Livas, J., Conklin, J. W., Caldwell,
R., Berti, E., McWilliams, S. T., Stebbins, R., Shoemaker, D., Ferrara, E. C.,
Larson, S. L., Shoemaker, D., Key, J. S., Vallisneri, M., Eracleous, M.,
Schnittman, J., Kamai, B., Camp, J., … Wass, P. (2019). The Laser
Interferometer Space Antenna: Unveiling the Millihertz Gravitational Wave Sky.
_Bulletin of the American Astronomical Society_ , _51_ , 77.
https://arxiv.org/abs/1907.06482

preTulczyjew, W. M. (1959). Motion of multipole particles in general
relativity theory binaries. _Acta Phys. Pol._ , _18_ , 393.

preVáňa, F., Düben, P., Lang, S., Palmer, T., Leutbecher, M., Salmond, D., &
Carver, G. (2017). Single precision in weather forecasting models: An
evaluation with the IFS. _Monthly Weather Review_ , _145_(2), 495–502.
https://doi.org/10.1175/MWR-D-16-0228.1

preWitzany, V. (2022). Action-angle coordinates for black-hole geodesics I:
Spherically symmetric and Schwarzschild. _arXiv e-Prints_ , arXiv:2203.11952.
https://doi.org/10.48550/arXiv.2203.11952

p