Vol.0 (20xx) No.0, 000–000

11institutetext: School of Physics and Astronomy, Sun Yat-sen University,
Zhuhai 519082, China; tanbxuan<EMAIL_ADDRESS>
22institutetext: CSST Science Center for the Guangdong-Hong Kong-Macau Greater
Bay Area, Sun Yat-Sen University, Zhuhai 519082, China
Received 20xx month day; accepted 20xx month day

# Revealing Double White Dwarf Mergers with Multi-messenger Signals

He-Wen Yang 1122 Pak-Hin Thomas Tam 1122 Lili Yang 1122

###### Abstract

A significant number of double white dwarfs (DWDs) are believed to merge
within the Hubble time due to the gravitational wave (GW) emission during
their inspiraling phase. The outcome of a DWD system is either a type Ia
Supernova as the double-degenerate model, or a massive, long-lasting merger
remnant. Expected multi-messenger signals of these events will help us to
distinguish detailed merging physical processes. In this work, we aim to
provide a generic scenario of DWD merging, investigate the emission of all
major messengers, with a focus on GWs and neutrinos. Our goal is to provide
some guidance for current and future (collaborative) efforts of multi-
messenger observations. Throughout the merging evolution of a DWD system,
different messengers (GW, neutrino and electromagnetic wave) will dominate at
different times. In this work, we show that DWD merger events located at the
distance of $1\ \mathrm{kpc}$ can indeed produce detectable signals of GWs and
neutrinos. The GW frequency are in $0.3-0.6\ \mathrm{Hz}$ band around $10$
days before tidal disruption begin. We estimate that in optimistic situations,
the neutrino number detected by upcoming detectors such as JUNO and Hyper-
Kamiokande can reach O(1) for a DWD merging event at $\sim$1 kpc.

###### keywords:

white dwarfs — binaries: close — neutrinos — gravitational waves

##  1 Introduction

Many evidences show that most stars in the Universe are not single isolated
stars, but in binary systems (Sana et al. 2012; Duchêne & Kraus 2013; El-Badry
et al. 2021). Stellar evolution models suggest that most stars will become
white dwarfs (WDs). Therefore, the most common compact object mergers form out
of double white dwarf (DWD) binaries.

In the past few years, some studies of DWD mergers show that every close, DWD
binary may eventually merge due to energy loss and angular momentum loss
through gravitational wave (GW) radiation (Shen 2015). DWD mergers can lead to
various kinds of observable astronomical phenomena including some most
luminous transient events in the Universe. On one hand, DWDs are treated as
one of the progenitors of type Ia Supernovae (SNe Ia) for a long time, known
as the double-degenerate (DD) scenario (Yungelson & Kuranov 2016; Livio &
Mazzali 2018). Other suggested manifestations of DWD mergers are transient
electromagnetic signals such as fast blue optical transients (Lyutikov &
Toonen 2019) and less-luminous nova-type explosions (Roy et al. 2022). On the
other hand, if the merger core fails to ignite, the merger product is not
clear, suggestions include a rapidly rotating WD with a strong magnetic field
(Ji et al. 2013; Rueda et al. 2019) which latter collapses to a neutron star
(NS) (Ruiter et al. 2019; Liu & Wang 2020), or a heavy, long-lasting WD
(Schwab 2021; Wu et al. 2022). Therefore, a full understanding of the outcomes
of DWD mergers not only provides us with the merger rates of the most common
compact objects in the Universe, it will also help us to understand the origin
of a broad range of transient phenomena and the formation of a class of
specific WDs and NSs.

DWD systems are important gravitational wave sources (Nelemans et al. 2001;
Huang et al. 2020). During the inspiral phase, DWD systems are stable, low
frequency ($1-10\ \mathrm{mHz}$) GW sources, and are expected to be the first
detected GW sources for the space-based GW observatories like LISA (Robson et
al. 2019), Taiji and TianQin (Gong et al. 2021). In the Galaxy, the high
number distribution of low-mass (CO-He) DWDs makes them the dominant GW
sources at frequency $\geq 5\ \mathrm{mHz}$, while massive (CO-CO or ONeMg)
DWDs are less numerous, and the former evolve to merge more quickly than the
latter after common-envelope phases or stable Roche lobe overflow phase (Yu &
Jeffery 2010). However, if they successfully evolve to very close binary
(i.e., they tend to merge soon), the GW frequency of the system can reach
$0.1\ \mathrm{Hz}$ or even $1\ \mathrm{Hz}$ (Maselli et al. 2020; Zou et al.
2020).

Neutrinos are believed to be generated from astrophysical sources such as
stellar cores, novae, core-collapse supernovae and compact binary mergers at
$\mathrm{MeV}$ energies or from the cosmic rays accelerators at high and
ultra-high energies (Halzen & Hooper 2002; Katz & Spiering 2012; Vitagliano et
al. 2020). However, the merging physics of DWDs is still unclear, and few
works involving numerical simulations focus on neutrino emission. During the
merging stage, the interaction region of DWD can create density and
temperature conditions which are high enough for efficient neutrino production
such that neutrinos can dominate the carried-away energy (Aznar-Siguán et al.
2014). In recent years, some single-degenerate models indicate that “normal”
SNe Ia bursts can produce significant neutrino flux that is four orders of
magnitude smaller than that of core-collapse SNe, regardless of the exact
explosion mechanism (Kunugise & Iwamoto 2007; Odrzywolek & Plewa 2011;
Seitenzahl et al. 2015; Wright et al. 2017). However, the neutrino flux
associated with double-degenerate model has rarely been studied by
simulations. It has previously been suggested that neutrino signals could help
to distinguish these two progenitor models of SNe Ia (Raj 2020). It is
commonly thought that SNe Ia come from the detonation of WD carbon cores when
the temperature and density reach the ignition condition. Besides, the
neutrino production is efficient when the temperature reaches above $10^{9}\
\mathrm{K}$, and carbon burning can reach such temperature (Itoh et al. 1996),
and in such case neutrino production will dominate the energy release.

In this work, we investigate the multi-messenger detection prospects of DWDs,
focusing on gravitational waves and neutrinos. We mainly focus on the CO DWD
systems which are detached binaries. In Section 2, we provide a simple
description of DWD merging processes and a schematic overview of possible
multi-messenger signals as well as their emission times for a merging event.
In Section 3, we assume a toy model for gravitational wave and neutrino
production when DWDs merge, introducing some physical quantities to describe
the two types of signals. In Section 4, we calculate the amplitude spectral
density for the final inspiraling evolution and the neutrino number events
from the merging process, and discuss their detectability using current and
planned detectors. Finally, we present the discussions and conclusion of our
results in Section 5.

##  2 The mechanism of double white dwarf mergers

###  2.1 Physical processes

In recent years, the merging dynamics of DWDs have been investigated by
smoothed particle hydrodynamics (SPH) simulations, which reveal the details of
the merging processes (Lorén-Aguilar et al. 2009; Pakmor et al. 2010, 2012;
Raskin et al. 2012; Ji et al. 2013; Sato et al. 2015, 2016). In general, it is
thought that the less massive WD could be tidally disrupted by the more
massive WD that remains. A small amount of mass of the disrupted one, which
carries most of the angular momentum, are fastly accreted onto the surface of
the undisrupted one, forming a thick Keplerian disk.

Depending on the initial masses and mass transfer rate of the DWD, it may
result in the immediate ignition of the carbon core of the more massive WD
that explodes (a violent merger), or that the product fails to be detonated
after they merge. In the latter case, a magnetized corona above the disk and a
strongly magnetized bi-conical jets perpendicular to the disk can form because
of the development of the magneto-rotational instability within the disk (Ji
et al. 2013). And a viscous outflow is driven at the interface of the corona
and the jet. It is also argued that the merging may go through some viscous
processes with a timescale longer than the dynamical timescale, producing SN
Ia explosion. Otherwise, if there is no explosion, the merge product is a more
massive (ONeMg) WD. Even when the remnant mass exceeds the Chandrasekhar
limit, the remnant will undergo a long-term evolution before it finally
collapses to form a NS (Schwab 2021).

###  2.2 Multi-messenger Signals

In Figure 1, we depict anticipated emission signals from DWD merger events. A
plausible timeline of the signals associated with various processes is
outlined as follows:

  * •

During the inspiral phase, the emission is dominated by GW. It is believed
that GW emission reaches a maximum when tidal disruption occurs and thereafter
is largely suppressed.

  * •

When tidal disruption starts, the dynamical timescale is thought to be $\sim
1\ \mathrm{s}$ as estimated by the free-fall motion of the accreting mass
across the tidal disruption radius of the secondary.

  * •

Subsequently, a neutrino burst with a time duration of $\mathcal{O}(1\
\mathrm{s})$ follows immediately the dynamical merging process.

  * •

If SNe Ia explosion occurs, most materials of the DWD is thrown out, giving
rise to a supernova light curves that may differ from (and indeed are likely
broader than) those in Chandrasekhar-mass explosion models (Fryer et al. 2010;
Moll et al. 2014). The luminosity reaches the maximum at $\sim 30$ days after
they merge.

  * •

If there is no explosion, then in the viscous outflow process, cosmic rays are
being accelerated via magnetic dissipation in the outflow region, subsequently
leading to high-energy neutrino emission $\sim 10^{4}\ \mathrm{s}$ after the
merge (Xiao et al. 2016). It is also argued that gamma rays might be produced
too, but can hardly be detected because of the surrounding optically thick
environment. About $\sim 7$ days, after the merge optical photons are emitted
due to the adiabatic cooling of ejecta. In addition, fallback accretion and
spindown of the newly-formed central WD will produce X-rays about $150-200$
days post-merging (Rueda et al. 2019).

Figure 1: Plausible signals timeline before and after DWD merging in a multi-
messenger perspective

###  2.3 Merger Rate

The merger rate of DWDs is the highest among all types of merger events in the
Universe. In the Milky Way, the merger rate of DWD is expected to be
$(1-80)\times 10^{-13}\ \mathrm{yr}^{-1}M_{\odot}^{-1}\ (2\sigma)$ or $(7\pm
2)\times 10^{-13}\ \mathrm{yr}^{-1}M_{\odot}^{-1}\ (1\sigma)$ (Maoz &
Hallakoun 2017). For comparison, by using the Milky Way stellar mass
$6.4\times 10^{10}\ M_{\odot}$ and the extrapolating factor of Milky Way-type
galaxies $1.16\times 10^{-2}\ \mathrm{Mpc}^{-3}$ (Kalogera et al. 2001) for
DWDs, the merger rate of double neutron stars (DNSs) and double black holes
(DBHs) are shown in Table 1, which are smaller than DWDs by $2-4$ orders of
magnitude. Additionally, previous studies also show that DWD merger rate are
consistent and on the same order of SNe Ia rates deduced from observations
(Toonen et al. 2012; Maoz et al. 2014, 2018).

Table 1: DWD, DNS and DBH merger rate in the local Universe and in the Milky
Way

System | $\mathcal{R}\ [\mathrm{Gpc}^{-3}\ \mathrm{yr}^{-1}]$ |  MW Event Rate $[(100\ \mathrm{yr})^{-1}]$ | References
---|---|---|---
DWD | $(0.74-5.94)\times 10^{6}\ (2\sigma)$ | $0.64-51$ | (Maoz & Hallakoun 2017)
DWD | $(5.2\pm 1.5)\times 10^{5}\ (1\sigma)$ | $4.5\pm 1.3$
DNS | $320^{+490}_{-240}$ | $2.8^{+4.2}_{-2.1}\times 10^{-3}$ |  (Abbott et al. 2021)
DBH | $23.9^{+14.3}_{-8.6}$ | $2.1^{+1.2}_{-0.7}\times 10^{-4}$ |  (Abbott et al. 2021)

##  3 A multi-messenger perspective of double white dwarf mergers

In this section, we focus on the anticipated GW and neutrino emission from a
DWD merging event.

###  3.1 Gravitational Waves

In this work, we assume that there is no mass-transfer nor interaction between
the two WDs until tidal disruption happens, because we assume that both
components of the binary are of the type CO-WD or even more massive (ONeMg)
WD. It is therefore reasonable to assume a two point-mass approximation during
the inspiral phase.

Considering a DWD system with component masses $m_{1}$ and $m_{2}$ and an
orbital separation of $a$, the orbital frequency is derived from its Keplerian
motion as:

$\displaystyle f_{\mathrm{orb}}=\sqrt{\frac{G(m_{1}+m_{2})}{2\pi^{2}a^{3}}},$
(1)

where $G$ is the gravitational constant. In the case of a circular orbit, the
GW is radiated only in the second harmonics, and hence the GW frequency $f$ is
twice that of the orbital frequency $f=2f_{\mathrm{orb}}$. During the
inspiraling phase, the two independent polarization states (or waveforms) $+$
and $\times$ of GW radiation from DWD in the quadrupole approximation (Landau
& Lifshitz 1962) is

$\displaystyle h_{+}$
$\displaystyle=\frac{\left(G\mathcal{M}_{c}\right)^{5/3}}{c^{4}d}\left(\pi
f\right)^{2/3}2(1+\cos^{2}{i})\cos(2\pi ft),$ (2) $\displaystyle h_{\times}$
$\displaystyle=\frac{\left(G\mathcal{M}_{c}\right)^{5/3}}{c^{4}d}\left(\pi
f\right)^{2/3}4\cos{i}\sin(2\pi ft).$ (3)

where $c$ is the speed of
light,$\mathcal{M}_{c}=(m_{1}m_{2})^{3/5}/(m_{1}+m_{2})^{1/5}$ is the chirp
mass and $d$ is the distance to the source. The GW emission luminosity is:

$\displaystyle
L_{\mathrm{GW}}=\dot{E}=\frac{32}{5}\frac{G^{4}}{c^{5}}\frac{m_{1}^{2}m_{2}^{2}(m_{1}+m_{2})}{a^{5}}.$
(4)

Alternatively, it can be expressed as a function of chirp mass and GW
frequency

$\displaystyle L_{\mathrm{GW}}=\frac{32}{5}\frac{G^{7/3}}{c^{5}}\left(\pi
f\mathcal{M}_{c}\right)^{10/3}.$ (5)

Due to GW emission loss, the GW frequency evolves as:

$\displaystyle\dot{f}=\frac{96}{5\pi}\frac{(G\mathcal{M}_{c})^{5/3}}{c^{5}}(\pi
f)^{11/3}.$ (6)

Introducing the strain amplitude $h_{0}$ (or the instantaneous root-mean-
square amplitude),

$\displaystyle h_{0}$
$\displaystyle=\sqrt{\left<h_{+}^{2}\right>+\left<h_{\times}^{2}\right>}=\sqrt{\frac{32}{5}}\frac{\left(G\mathcal{M}_{c}\right)^{5/3}}{c^{4}d}\left(\pi
f\right)^{2/3}.$ (7)

Averaging over a full orbital period and all inclination angle, it is related
to the GW luminosity by (Postnov & Yungelson 2014)

$\displaystyle h_{0}^{2}=\frac{1}{(\pi
d)^{2}}\frac{G}{c^{3}}\frac{L_{\mathrm{GW}}}{f^{2}}.$ (8)

The characteristic strain amplitude $h_{c}$ for inspiraling binaries is given
by (Finn & Thorne 2000; Moore et al. 2014)

$\displaystyle
h_{c}^{2}=\left(\frac{2f^{2}}{\dot{f}}\right)h_{0}^{2}=\frac{2}{3}\frac{(G\mathcal{M}_{c})^{5/3}}{c^{3}d^{2}}f\left(\pi
f\right)^{-4/3}.$ (9)

The amplitude spectral density (ASD) is defined as
$\sqrt{S_{h}}=h_{c}f^{-1/2}$, where $S_{h}$ is the power spectral density of
sources (Moore et al. 2014). Similarly, for GW detectors, $\sqrt{S_{N}(f)}$ is
the effective strain spectral density of the detector noise used in
sensitivity curves (Robson et al. 2019; Huang et al. 2020).

Inspiraling of compact binaries can be well-modeled in theoretical analysis
while the merging of binaries (especially those bursting sources) is poorly
modeled (Moore et al. 2014). As mentioned above, we consider the binary
evolution during the inspiral phase, during which GW frequency is slowly
changing $\Delta f=\dot{f}\Delta t$, and the GW signal is weak and remains
almost unchanged for many years, then the corresponding signal is accumulated
at detectors over a long time $\Delta t$, increasing the signal-to-noise ratio
(SNR). The ASD can be described as $\sqrt{S_{h}}=\sqrt{2}h_{0}\sqrt{\Delta
t}$. However, for latter phases of DWD evolution (i.e., right before the tidal
disruption), the frequency and strain is enhanced within a very short time
(e.g., several orbital periods), they can be treated as burst sources. These
systems produce large amplitude signals that could be much higher than the
detector noise and specific waveform models are not necessary in signal
detection. In practice, the signal amplitude described for burst sources is
the root-sum-square amplitude

$\displaystyle
h_{\mathrm{rss}}^{2}=\int\mathrm{d}t\left(\left|h_{+}(t)\right|^{2}+\left|h_{\times}(t)\right|^{2}\right).$
(10)

It can be approximated to $h_{\mathrm{rss}}\simeq|\tilde{h}(f)|\sqrt{\Delta
f}$ by assuming that the GW mode is linearly polarized and $|\tilde{h}(f)|$ is
almost constant within a frequency band $\Delta f$, where $\tilde{h}(f)$ is
the Fourier transform of $h(t)$. However, the quantity $h_{\mathrm{rss}}$ is
on the same order as the characteristic strain $h_{c}$ as defined before
(Moore et al. 2014). In this work, the main phase we focused on is the last
stages of inspiraling before tidal disruption, and we therefore use $h_{c}$ to
represent GW emission.

###  3.2 Neutrinos

####  3.2.1 neutrinos production

We now consider neutrino emission from the merger event. When the density and
temperature are high enough, thermal neutrino processes begin and produce
$\mathcal{O}(1)\ \mathrm{MeV}$ neutrinos, at the same time neutrino cooling
becomes efficient. In general, five thermal neutrino processes are taken into
account: electron–positron annihilation, plasmon decay, photoemission,
neutrino Bremsstrahlung and recombination, which produce neutrinos with all
kinds of flavors (Itoh et al. 1996). Once the temperature exceeds $10^{9}\
\mathrm{K}$, electron–positron annihilation

$\displaystyle e^{+}+e^{-}$ $\displaystyle\rightarrow\nu+\bar{\nu},$ (11)

starts to dominate the thermal neutrino production (Itoh et al. 1996).

When the merging process is violent, the temperature can easily reach $\geq
5\times 10^{9}\ \mathrm{K}$ such that weak processes including electron and
positron capture becomes more efficient than thermal neutrino production
channels. These weak processes mainly produce electron neutrinos and anti-
electron neutrinos via,

$\displaystyle p+e^{-}$ $\displaystyle\rightarrow n+\nu_{e},$ (12)
$\displaystyle n+e^{+}$ $\displaystyle\rightarrow p+\bar{\nu}_{e},$ (13)
$\displaystyle(Z,A)+e^{-}$ $\displaystyle\rightarrow(Z-1,A)+\nu_{e},$ (14)
$\displaystyle(Z-1,A)+e^{+}$ $\displaystyle\rightarrow(Z,A)+\bar{\nu}_{e}.$
(15)

During the DWDs merging process, the material is proton-rich so that the
electron capture happens mainly between free protons/nuclei and electrons,
resulting in a dominant $\nu_{e}$ emission.

We assume that neutrinos follow the Fermi-Dirac distribution with zero
chemical potential, and the neutrino number emitted per energy per time per
area could be written as

$\displaystyle\Phi_{\nu}(E)=\frac{L_{\nu}}{4\pi
d^{2}}\frac{120}{7\pi^{2}}\frac{E^{2}}{T^{4}}\frac{1}{1+e^{E/k_{B}T}},$ (16)

where $L_{\nu}$ is the neutrino luminosity, $T$ the neutrino temperature, $E$
the neutrino energy and $d$ the distance to the earth. Here we take the
average neutrino energy to be $\left<E\right>\sim 3.15k_{B}T$, a result from
the Fermi-Dirac distribution.

Neutrino luminosity and average energy strongly depend on the merging
processes. Dynamical interactions of DWD failing to ignite can result in
optimistic values of neutrino luminosity with the total energy $\sim 10^{48}\
\mathrm{erg}$ (Aznar-Siguán et al. 2014). As mentioned before, SNe Ia can
produce a large amount of neutrinos (Kalogera et al. 2001; Wright et al.
2017), so it should be reasonable to assume that SNe Ia from the double-
degenerate channel can also generate a neutrino burst with the luminosity
$\sim 10^{49}\ \mathrm{erg\ s^{-1}}$ or even $\sim 10^{50}\ \mathrm{erg\
s^{-1}}$. Here we adopt three characteristic temperature values in our study:

  * •

$\geq 10^{9}\ \mathrm{K}$: this is above the carbon burning temperature,
electron-positron pair annihilation starts to dominate the thermal production
process;

  * •

around $4\times 10^{9}\ \mathrm{K}$: nuclear statistical equilibrium (NSE)
starts, and the electron capture process become efficient;

  * •

around $10^{10}\ \mathrm{K}$: it is the typical core neutrino temperature in a
SN Ia explosion, resulting in the neutrino average energy $\left<E\right>\sim
3\ \mathrm{MeV}$.

During the neutrino burst, the weak or thermal production processes happen.
For weak processes, the timescale of NSE can be approximated as
$\tau_{\mathrm{NSE}}\sim 1\ \mathrm{s}$ for SNe Ia explosion (Kunugise &
Iwamoto 2007), which is the same order as the explosion timescale. Assuming
that DWDs merging could also reach the NSE condition (i.e., cases 1, 2, 3(a),
and 3(b)), it would be reasonable to assume that it also lasts for about 1
second. While in the thermal production, which is related to density and
temperature during the merge processes, we assume that it lasts $\sim 1\
\mathrm{s}$, to be consistent with the case of explosion. The neutrino
luminosity $L_{\nu}$, emission duration $\Delta t$ and temperature $T$ are
assumed in this work as shown in Table 2. In Table 2, we argue that Cases 1
and 2 likely represent situations which will lead to SNe Ia, and neutrinos
carry away $\sim(1-10)\%$ of the total released energy of $10^{51}\
\mathrm{erg}$. For Cases 3(a), 3(b) and 3(c), the luminosity is lower compared
with the SNe Ia case, and we consider cases with temperature values (i.e.,
$4\times 10^{9}\ \mathrm{K}$ or even $10^{10}\ \mathrm{K}$) which are somewhat
optimistic.

Table 2: The neutrino luminosity, duration time and temperature we adopt in
this work

Cases | $L_{\nu_{e}}\ [\mathrm{erg\ s^{-1}}]$ | $\Delta t\ [\mathrm{s}]$ | $T\ [\mathrm{K}]$ | Reference
---|---|---|---|---
Case 1 | $10^{50}$ | $1$ | $\left<E\right>=3\ \mathrm{MeV}$ |  (Kunugise & Iwamoto 2007)
Case 2 | $10^{49}$ | $2$ | $4\times 10^{9}$ |  (Wright et al. 2017)
Case 3(a) | $10^{48}$ | $1$ | $10^{10}$ |
Case 3(b) | $10^{48}$ | $1$ | $4\times 10^{9}$ |  (Aznar-Siguán et al. 2014)
Case 3(c) | $10^{48}$ | $1$ | $10^{9}$ |

Based on the assumed parameter values as in Table 2, the respective neutrino
flux values are calculated with Eq. (16). From Figure 2, for higher
temperature (averaged energy), the distribution of neutrino (flux) has broader
energy range.

Figure 2: Neutrino fluxes of five cases presented in Table 2, and the distance
from the source is $d=1\ \mathrm{kpc}$.

Here, the material is proton-rich (electron fraction $Y_{e}>0.5$) during DWD
merge, it implies that electron capture dominates neutrino production process,
and $\nu_{e}$ will be dominant during the emission. 111While in DNS merge
$\bar{\nu}_{e}$ dominates due to positron capture in neutron-rich matters.
Nevertheless, it is a neutronization process which losses lepton number by
neutrino release (nucleosynthesis that produce heavier elements also occurs).
On the other hand, if the thermal burning process dominates, neutrinos of all
flavors will be produced with similar proportion. However, if weak processes
become significant then the production rate of thermal processes become sub-
dominant at the same physical conditions, thus the thermal contribution to the
neutrino flux can be safely ignored. We note also that no neutrino trapping is
expected as occurred in core-collapse supernovae, because the matter density
of DWD merger is not high enough222Core-collapse supernovae are triggered by
massive star collapse. Before the explosion, the center core of star is so
dense ($\geq 10^{11}\ \mathrm{g\ cm^{-3}}$) that neutrinos cannot escape and
are trapped in the inner regions for some time. Only after the explosion
starts, the core density decreases and the trapping conditions do not exist
any more, leading to a neutrino burst..

####  3.2.2 neutrino oscillation effects

From the production region to surface, neutrino flavors will be regulated by
the Mikheyev-Smirnov-Wolfenstein (MSW) effect (Wolfenstein 1978; Mikheyev &
Smirnov 1985; Blennow & Smirnov 2013). These regulations are different
depending on the mass order hierarchy: the normal mass hierarchy (NH) or the
inverted mass hierarchy (IH) (Dighe & Smirnov 2000; Fogli et al. 2005). If we
denote the initial flux as $\Phi_{\nu_{\alpha}}^{0}$ (
$\Phi_{\bar{\nu}_{\alpha}}^{0}$) where $\alpha=e,\mu,\tau$, then with MSW
resonance effect, the flux at the surface is $\Phi_{\nu_{\alpha}}$
($\Phi_{\bar{\nu}_{\alpha}}$). When neutrinos propagate from the source
surface to the Earth, neutrino oscillation in vacuum has to be taken into
account. The neutrino fluxes observed $\phi_{\nu_{\alpha}}$
($\phi_{\bar{\nu}_{\alpha}}$) are different from $\Phi_{\nu_{\alpha}}$ (
$\Phi_{\bar{\nu}_{\alpha}}$). In general, the transition probabilities between
$\phi_{\nu_{\alpha}}$ ( $\phi_{\bar{\nu}_{\alpha}}$) and $\Phi_{\nu_{\alpha}}$
( $\Phi_{\bar{\nu}_{\alpha}}$) is:

$\displaystyle\phi_{\nu_{\beta}}(E)=\sum_{\alpha}\Phi_{\nu_{\alpha}}(E)P_{\nu_{\alpha}\rightarrow\nu_{\beta}}\left({L}/{E}\right),$
(17)
$\displaystyle\phi_{\bar{\nu}_{\beta}}(E)=\sum_{\alpha}\Phi_{\bar{\nu}_{\alpha}}(E)P_{\bar{\nu}_{\alpha}\rightarrow\bar{\nu}_{\beta}}\left({L}/{E}\right),$
(18)

i.e., it depends on the energy of the neutrinos and their travel distance
(evolution time). However, the astrophysical sources are so far way from the
Earth, so the probabilities should be averaged along the way and can be
described as

$\displaystyle P_{\alpha\beta}$
$\displaystyle=\left<P_{\nu_{\alpha}\rightarrow\nu_{\beta}}\right>=\sum_{i}^{3}\left|U_{\alpha
i}\right|^{2}\left|U_{\beta i}\right|^{2},$ (19)
$\displaystyle\bar{P}_{\alpha\beta}$
$\displaystyle=\left<P_{\bar{\nu}_{\alpha}\rightarrow\bar{\nu}_{\beta}}\right>=P_{\alpha\beta},$
(20)

where $U_{\alpha i}$ is the element of Pontecorvo-Maki-Nakagawa-Sakata (PMNS)
mixing matrix (Pontecorvo 1957; Maki et al. 1962)

$\displaystyle
U=\left(\begin{array}[]{ccc}c_{12}c_{13}&s_{12}c_{13}&s_{13}e^{-i\delta}\\\
-s_{12}c_{23}-c_{12}s_{23}s_{13}e^{i\delta}&c_{12}c_{23}-s_{12}s_{23}s_{13}e^{i\delta}&s_{23}c_{13}\\\
s_{12}s_{23}-c_{12}c_{23}s_{13}e^{i\delta}&-c_{12}s_{23}-s_{12}c_{23}s_{13}e^{i\delta}&c_{23}c_{13}\end{array}\right),$
(24)

where $s_{ij}=\sin\theta_{ij}$, $c_{ij}=\cos\theta_{ij}$, in which
$\theta_{ij}$ is the lepton flavor mixing angle and $\delta$ is the CP
violation phase. The above parameters are taken from Zyla et al. (2020).

For the situation of electron neutrino dominant,
$\Phi_{\nu_{e}}^{0}:\Phi_{\nu_{x}}^{0}:\Phi_{\bar{\nu}_{e}}^{0}:\Phi_{\bar{\nu}_{x}}^{0}=1:0:0:0$
, it leads to
$\phi_{\nu_{e}}:\phi_{\nu_{x}}:\phi_{\bar{\nu}_{e}}:\phi_{\bar{\nu}_{x}}=0.23:0.77:0:0$
for NH case and $0.32:0.68:0:0$ for IH case. Here, the subscript $e$ and $x$
denote the electron flavor and other ($\mu$ and $\tau$) flavors, and the bar
denote the antineutrinos with different types. If three flavors of neutrinos
are equally emitted, i.e.
$\Phi_{\nu_{e}}^{0}:\Phi_{\nu_{x}}^{0}:\Phi_{\bar{\nu}_{e}}^{0}:\Phi_{\bar{\nu}_{x}}^{0}=1/6:1/3:1/6:1/3$,
it has trivial effect
$\phi_{\nu_{e}}:\phi_{\nu_{x}}:\phi_{\bar{\nu}_{e}}:\phi_{\bar{\nu}_{x}}=1/6:1/3:1/6:1/3$.
In this case, the neutrino mass order hierarchy cannot be distinguished.

##  4 Detectability

###  4.1 Gravitational waves

The power of GW emission reaches the maximum when tidal disruption occurs and
is largely suppressed when the secondary mass starts to be disrupted. We
assume that the merging process occurs when the orbital separation $a_{m}$ is
equal to the tidal disruption radius:

$\displaystyle a_{m}=R_{\mathrm{td}}=qR_{2}$ (25)

where $q$, $R_{2}$ is the mass-ratio ($m_{2}/m_{1}$) and the radius of the
secondary WD, respectively. We assume the WD radius to be:

$\displaystyle R\approx\left(\frac{3m}{4\pi\rho}\right)^{\frac{1}{3}},$ (26)

where $\rho$ is the mean density of WD. The final cut-off GW frequency is then
determined by Eq. (1).

In this work, we consider four hypothetical DWD systems: $(1.0+0.8)M_{\odot}$,
$(1.0+1.0)M_{\odot}$, $(1.1+0.9)M_{\odot}$ and $(1.2+1.0)M_{\odot}$ which are
usually assumed in simulations and might mostly lead to SNe Ia, which also
have a high probability to produce neutrinos with significant luminosity. For
various masses, the mean densities are taken to be $10^{6}\ \mathrm{g\
cm^{-3}}$ (for $0.8\ M_{\odot}$), $2\times 10^{6}\ \mathrm{g\ cm^{-3}}$ (for
$0.9\ M_{\odot}$) and $3\times 10^{6}\ \mathrm{g\ cm^{-3}}$ (for $1.0\
M_{\odot}$), which we refer to Figure 1 of Zou et al. (2020).

Here, we use the Python package
LEGWORK333https://legwork.readthedocs.io/en/latest/ (Wagg et al. 2021) to
calculate the strain amplitude $h_{0}$ and characteristic strain amplitude
$h_{c}$ of binaries systems and the evolution of their orbital frequency
$f_{\mathrm{orb}}$. LEGWORK is also used to calculate SNRs of GW emitted from
inspiraling binary systems as detected by certain GW detectors. However, in
this work, we mainly focus on the latter stage of DWD merging process. We also
take the initial GW frequency to be $5\ \mathrm{mHz}$ for four systems, and
follow the evolution of GW frequency based on the two-point approximation as
described in Sec. 3.1.

Figure 3 show the luminosities of GW emission $L_{\mathrm{GW}}$ for these
systems from $10$ days before tidal disruption, $t_{\mathrm{td}}$, to the time
of disruption. It can be seen that for more massive systems, the luminosities
are higher as expected. With the same total mass, the equal mass-ratio (e.g.,
$(1.0+1.0)M_{\odot}$) system has a smaller luminosity than the unequal mass-
ratio (e.g., $(1.1+0.9)M_{\odot}$) one. Integrating luminosities from
$t_{\mathrm{td}}-10\ \mathrm{day}$ to $t_{\mathrm{td}}$, it follows that the
energy release by GW emission is on the order of $10^{48}-10^{49}\
\mathrm{erg}$.

Figure 3: The GW luminosity in final evolution phase (i.e., beginning $10\
\mathrm{days}$ before the tidal disruption). The line colors denote systems
with different component masses: $(1.0+0.8)M_{\odot}$, $(1.0+1.0)M_{\odot}$,
$(1.1+0.9)M_{\odot}$ and $(1.2+1.0)M_{\odot}$.

We assume the distance from source to the Earth to be $1\ \mathrm{kpc}$, and
calculate $h_{0}$ and $h_{c}$ which could represent the instantaneous moment
and _accumulated_ strength of GW emission. We show the evolution of $h_{0}$
and GW frequency of these systems since $10\ \mathrm{days}$ before tidal
disruption in Figure 4. It can be seen that $h_{0}$ and GW frequency are not
evolving with the same rate of change, i.e., the latter is changing more
rapidly. It implies that the strain amplitude varies across the frequency band
in the final stage of merging. Besides, the $h_{0}$ and $f$ of DWDs are
_slightly_ increasing whereas in the cases of DNSs and DBHs they change
suddenly. This contrast can be understood because the components of binaries
are more compact and have smaller radii in the latter cases. Thus, the cut-off
frequencies and the evolution rates of late stage are larger. On the other
hand, the merge of DWD seems gentle comparing with DNSs and DBHs.
Nevertheless, the GW signals will be detected by planned space-based
interferometers assuming a distance of $1\ \mathrm{kpc}$.

Figure 4: The final evolution of strain amplitude $h_{0}$ and GW frequency $f$
beginning $10$ days before the tidal disruption. The line colors represent the
same systems as in Fig. 3.

Concerning the detectability of the GW signal, it is a common practise to
compare the ASD $\sqrt{S_{h}}$ with the sensitivity curves of detectors where
$S_{h}$ represents the energy density. Because we are only concerned about the
final stage of inspiraling, we can safely ignore the galactic confusion noise
which is related to detector’s observing time. From Figure 5, if the sources
is at the distance of $1\ \mathrm{kpc}$, it is clear that the four systems we
assume have a significant signal strength that is easily detected by LISA
(Robson et al. 2019) and TianQin (Huang et al. 2020) in $0.3-0.6\ \mathrm{Hz}$
frequency band within final $10$ days before tidal disruption. We do not
calculate the related SNRs because they are so large and the noise is
unimportant due to the short time duration in our consideration.

Figure 5: The GW amplitude spectral density (ASD) of the four systems as
presented in Fig. 3. The sensitivity curves of LISA and TianQin are also
showed in this figure. GW frequency in these cases are in $0.3-0.6\
\mathrm{Hz}$ band which can be well detected by space-based GW detectors. Here
we assume the source is $1\ \mathrm{kpc}$ away from the Earth.

###  4.2 Neutrinos

Neutrinos from DWD merger events are in the $\mathrm{MeV}$ energy scale and
are observable via inverse beta decay (IBD) and electron-neutrino elastic
scattering (ES) channels. The IBD channel is typically used to detect
$\bar{\nu}_{e}$, while the ES channel is used to detect $\nu_{e}$ and other
types of neutrinos. The threshold energy of the IBD channel is $\sim 1.8\
\mathrm{MeV}$, while that of the ES channel heavily depends on the material
used in the detectors. For water Cherenkov detectors such as Super-Kamiokande
(Super-K) (Abe et al. 2016) and Hyper-Kamiokande (Hyper-K) (Abe et al. 2018),
the threshold energy of recoil electrons is $\sim 0.77\ \mathrm{MeV}$,
corresponding to a minimum neutrino energy of $\sim 0.97\ \mathrm{MeV}$. The
fiducial volume of the inner detector of Super-K and Hyper-K is $22.5\
\mathrm{kt}$ and $220\ \mathrm{kt}$, respectively. The corresponding target
particle number is
$N_{\mathrm{tar}}\approx(M_{\mathrm{tar}}/m_{p})\times(2/18)=6.7\times
10^{30}(M_{\mathrm{tar}}/1\ \mathrm{kt})$, where $M_{\mathrm{tar}}$ is the
effective mass of water and $m_{p}$ is the proton mass. Therefore, the target
electron number of Super-K and Hyper-K are, respectively, $1.51\times 10^{32}$
and $1.47\times 10^{33}$. For JUNO – a liquid organic scintillation, the
target number is $1.29\times 10^{33}$ (An et al. 2016) with a rather low
detection energy threshold of $0.2\ \mathrm{MeV}$ (Fang et al. 2020) that
corresponds to the neutrino threshold energy of $0.35\ \mathrm{MeV}$. The
cross-section of IBD and ES interaction channels are taken from the code
package provided by Kate et al. (2017), and the formulae of IBD and ES are
taken from Strumia & Vissani (2003); Marciano & Parsa (2003). In this work,
$\nu_{e}$ dominate the emission, so we mainly focus on the ES channel.

The neutrino event number can be calculated with

$\displaystyle N_{\nu}=N_{\mathrm{tar}}\Delta
t\int^{E_{\mathrm{max}}}_{E_{\mathrm{th}}}\phi(E)\sigma(E)\eta(E)\mathrm{d}E,$
(27)

where $N_{\mathrm{tar}}$ is the target particle number of protons, $\Delta t$
denotes the event duration ( observed time window), $\sigma(E)$ is the
interaction-channel cross section, $\phi(E)$ is neutrino fluxes at earth,
$\eta(E)$ is the efficiency of neutrino detector and is usually taken as
$\eta\simeq 0.9$.

We could calculate the neutrino numbers and show in Table 3.

Table 3: The electron neutrino $\nu_{e}$ (other neutrino $\nu_{x}$) events
number calculated in detectors sensitive to $\mathcal{O}(1\ \mathrm{MeV})$. NH
and IH denotes the normal and inverted mass order hierarchy for neutrino
respectively. In this table, the values are presented in scientific notation
and values in bracket are for ${\nu}_{x}$

Detector | hierarchy | Case 1 | Case 2 | Case 3(a) | Case 3(b) | Case 3(c)
---|---|---|---|---|---|---
JUNO | NH | 1.2E-00 | 2.0E-01 | 1.2E-02 | 9.8E-03 | 1.6E-03
(7.2E-01) | (1.3E-01) | (7.2E-03) | (6.7E-03) | (1.2E-03)
IH | 1.6E-00 | 2.7E-01 | 1.6E-02 | 1.4E-02 | 2.2E-03
(6.4E-01) | (1.2E-01) | (6.4E-03) | (6.0E-03) | (1.1E-03)
Super-K | NH | 1.4E-01 | 1.8E-02 | 1.3E-03 | 8.9E-04 | 5.0E-06
(8.4E-02) | (1.2E-02) | (8.3E-04) | (6.0E-04) | (3.5E-06)
IH | 1.9E-01 | 2.5E-02 | 1.9E-03 | 1.2E-03 | 7.0E-06
(7.4E-02) | (1.1E-02) | (7.3E-04) | (5.3E-04) | (3.1E-06)
Hyper-K | NH | 1.3E-00 | 1.7E-01 | 1.3E-02 | 8.7E-03 | 4.9E-05
(8.2E-01) | (1.2E-01) | (8.1E-03) | (5.8E-03) | (3.4E-05)
IH | 1.8E-00 | 2.4E-01 | 1.8E-02 | 1.2E-02 | 6.8E-05
(7.2E-01) | (1.0E-01) | (7.2E-03) | (5.1E-03) | (3.0E-05)

From Table 3, one can see that in the most optimistic situation (Case 1), the
neutrino number $N_{\nu_{e}}$ ($N_{\nu_{x}}$) can reach about $2$ (at the
order of $\mathcal{O}(1)$). $N_{\nu}$ is different between NH and IH, i.e.,
the ratio of NH and IH ($N_{\mathrm{NH}}/N_{\mathrm{IH}}$) is $0.72$ for
$\nu_{e}$ and $1.13$ for $\nu_{x}$. The neutrino event number of Case 3(c) is
very small because the average energy of this case are around the threshold of
water Cherenkov detectors. It is clear that the higher temperature of
neutrinos, the more events we observe in detectors. Furthermore, if the
evolution of luminosity is specified, one could expect the event numbers
evolves with time, which reflect the physics of merge processes, particularly
the explosion mechanism if it leads to SN Ia. However, these processes and
mechanism are still poorly studied, and we simply take the average of
luminosity with $\Delta t\sim 1\ \mathrm{s}$ time window. Nevertheless, it
means that the neutrino production processes are important ingredients in
understanding the merging physics.

##  5 Discussion and Conclusions

In this work, we show that DWD merger events are able to produce detectable GW
and neutrino signals if they are located at distance $1\ \mathrm{kpc}$. From
GW point of view, the GW frequency DWD systems are within $0.3-0.6\
\mathrm{Hz}$ band in $10$ days before the tidal disruption in the four example
DWD systems that we have studied. Unlike DNSs and DBHs, which are the main
potential sources of ground-based GW detectors such as the Advanced LIGO and
Virgo, DWD systems are more easily to be detected by space-based GW detectors
such as LISA and TianQin at the last stage of inspiraling. Besides, the
evolution of strains and GW frequency are moderate while for DNSs and DBHs
they are more dramatic. This can be understood since WD are not so dense and
the GW signals reach the maximum at tidal disruption radius and thereafter the
two bodies merge. The amplitude spectral density of these systems are large
comparing to the sensitive curves of LISA and TianQin for the final short time
interval, implying that DWDs merging could be well detected by space-based GW
detectors.

From neutrino point of view, given that the average energy from these events
is $\mathcal{O}(1)\ \mathrm{MeV}$, the upper limit of neutrino event numbers
is on the order of $\mathcal{O}(1)$ for current and upcoming detectors such as
Super-K, Hyper-K and JUNO. We found that the merger remnant cases (i.e., case
3(a), 3(b) and 3(c)) are harder to detect because of the lower neutrino
luminosity. In $\mathcal{O}(1)\ \mathrm{MeV}$ energy, the background mainly
come from solar neutrinos but it can be easily removed given the short
duration time. We also show that the detected neutrino flux could help to
distinguish the neutrino mass order hierarchy when taking the neutrino
oscillation into account. Therefore, we need detectors with larger volume
(i.e., more target number) and lower energy threshold, which will also extend
the detection horizon distance. More theoretical studies about merging physics
and explosion mechanisms of DWDs are also necessary to provide the neutrino
luminosity, spectra and oscillation mode.

In the era of multi-messenger astronomy, the DWD merging events are important
sources and can be studied using various multi-messenger emissions. Apart from
GW and neutrino signals, the subsequent electromagnetic signals such as
optical-infrared, X-rays and bolometric light curves are also interesting and
important observable signatures (Rueda et al. 2019; Moll et al. 2014; Aznar-
Siguán et al. 2014). Even if the electromagnetic signal could be similar
between different explosion mechanisms (e.g., DD and SD SNe Ia) and/or merging
processes, the neutrino detection can help to distinguish theoretical
scenarios once the neutrino production models of DWDs merging are specified.
The detection of these sources can join the next-generation supernova early
warning system (SNEWS 2.0) (Al Kharusi et al. 2021), and will guide the
observation of DWD merging events in the future.

###### Acknowledgements.

This work is supported by the National Natural Science Foundation of China
(NSFC) grants 11633007, 12005313 and U1731136, Guangdong Major Project of
Basic and Applied Basic Research (Grant No. 2019B030302001), Key Laboratory of
TianQin Project (Sun Yat-sen University) of the Ministry of Education, and the
China Manned Space Project (No. CMS-CSST-2021-B09).

## References

  * Abbott et al. (2021) Abbott, R., Abbott, T. D., Abraham, S., et al. 2021, The Astrophysical Journal Letters, 913, L7
  * Abe et al. (2016) Abe, K., Haga, Y., Hayato, Y., et al. 2016, Astroparticle Physics, 81, 39
  * Abe et al. (2018) Abe, K., Abe, K., Aihara, H., et al. 2018, arXiv e-prints, arXiv:1805.04163
  * Al Kharusi et al. (2021) Al Kharusi, S., BenZvi, S. Y., Bobowski, J. S., et al. 2021, New Journal of Physics, 23, 031201
  * An et al. (2016) An, F., An, G., An, Q., et al. 2016, Journal of Physics G: Nuclear and Particle Physics, 43, 030401
  * Aznar-Siguán et al. (2014) Aznar-Siguán, G., García-Berro, E., Magnien, M., & Lorén-Aguilar, P. 2014, MNRAS, 443, 2372
  * Blennow & Smirnov (2013) Blennow, M., & Smirnov, A. Y. 2013, Advances in High Energy Physics, 972485
  * Dighe & Smirnov (2000) Dighe, A. S., & Smirnov, A. Y. 2000, Physical Review D, 62, 033007
  * Duchêne & Kraus (2013) Duchêne, G., & Kraus, A. 2013, Annual Review of Astronomy and Astrophysics, 51, 269
  * El-Badry et al. (2021) El-Badry, K., Rix, H.-W., & Heintz, T. M. 2021, Monthly Notices of the Royal Astronomical Society, 506, 2269
  * Fang et al. (2020) Fang, X., Zhang, Y., Gong, G., et al. 2020, Journal of Instrumentation, 15, P03020
  * Finn & Thorne (2000) Finn, L. S., & Thorne, K. S. 2000, Phys. Rev. D, 62, 124021
  * Fogli et al. (2005) Fogli, G. L., Lisi, E., Mirizzi, A., & Montanino, D. 2005, J. Cosmology Astropart. Phys, 2005, 002
  * Fryer et al. (2010) Fryer, C. L., Ruiter, A. J., Belczynski, K., et al. 2010, ApJ, 725, 296
  * Gong et al. (2021) Gong, Y., Luo, J., & Wang, B. 2021, Nature Astronomy, 5, 881
  * Halzen & Hooper (2002) Halzen, F., & Hooper, D. 2002, Reports on Progress in Physics, 65, 1025
  * Huang et al. (2020) Huang, S.-J., Hu, Y.-M., Korol, V., et al. 2020, Phys. Rev. D, 102, 063021
  * Itoh et al. (1996) Itoh, N., Hayashi, H., Nishikawa, A., & Kohyama, Y. 1996, ApJS, 102, 411
  * Ji et al. (2013) Ji, S., Fisher, R. T., García-Berro, E., et al. 2013, The Astrophysical Journal, 773, 136
  * Kalogera et al. (2001) Kalogera, V., Narayan, R., Spergel, D. N., & Taylor, J. H. 2001, ApJ, 556, 340
  * Kate et al. (2017) Kate, S., Joshua, B., Justin, V., Felix, M., & Erin, C. 2017, SuperNova Observatories with GLoBES, https://github.com/SNOwGLoBES/snowglobes, (accessed April 26, 2022)
  * Katz & Spiering (2012) Katz, U. F., & Spiering, C. 2012, Progress in Particle and Nuclear Physics, 67, 651
  * Kunugise & Iwamoto (2007) Kunugise, T., & Iwamoto, K. 2007, PASJ, 59, L57
  * Landau & Lifshitz (1962) Landau, L. D., & Lifshitz, E. M. 1962, The classical theory of fields; 2nd ed., Course of theoretical physics (London: Pergamon), trans. from the Russian
  * Liu & Wang (2020) Liu, D., & Wang, B. 2020, Monthly Notices of the Royal Astronomical Society, 494, 3422
  * Livio & Mazzali (2018) Livio, M., & Mazzali, P. 2018, Phys. Rep., 736, 1
  * Lorén-Aguilar et al. (2009) Lorén-Aguilar, P., Isern, J., & García-Berro, E. 2009, A&A, 500, 1193
  * Lyutikov & Toonen (2019) Lyutikov, M., & Toonen, S. 2019, MNRAS, 487, 5618
  * Maki et al. (1962) Maki, Z., Nakagawa, M., & Sakata, S. 1962, Progress of Theoretical Physics, 28, 870
  * Maoz & Hallakoun (2017) Maoz, D., & Hallakoun, N. 2017, Monthly Notices of the Royal Astronomical Society, 467, 1414
  * Maoz et al. (2018) Maoz, D., Hallakoun, N., & Badenes, C. 2018, Monthly Notices of the Royal Astronomical Society, 476, 2584
  * Maoz et al. (2014) Maoz, D., Mannucci, F., & Nelemans, G. 2014, Annual Review of Astronomy and Astrophysics, 52, 107
  * Marciano & Parsa (2003) Marciano, W. J., & Parsa, Z. 2003, Journal of Physics G Nuclear Physics, 29, 2629
  * Maselli et al. (2020) Maselli, A., Marassi, S., & Branchesi, M. 2020, A&A, 635, A120
  * Mikheyev & Smirnov (1985) Mikheyev, S. P., & Smirnov, A. Y. 1985, Yadernaya Fizika, 42, 1441
  * Moll et al. (2014) Moll, R., Raskin, C., Kasen, D., & Woosley, S. E. 2014, ApJ, 785, 105
  * Moore et al. (2014) Moore, C. J., Cole, R. H., & Berry, C. P. L. 2014, Classical and Quantum Gravity, 32, 015014
  * Nelemans et al. (2001) Nelemans, G., Yungelson, L. R., & Portegies Zwart, S. F. 2001, A&A, 375, 890
  * Odrzywolek & Plewa (2011) Odrzywolek, A., & Plewa, T. 2011, A&A, 529, A156
  * Pakmor et al. (2010) Pakmor, R., Kromer, M., Röpke, F. K., et al. 2010, Nature, 463, 61
  * Pakmor et al. (2012) Pakmor, R., Kromer, M., Taubenberger, S., et al. 2012, ApJ, 747, L10
  * Pontecorvo (1957) Pontecorvo, B. 1957, Zh. Eksp. Teor. Fiz., 34, 247
  * Postnov & Yungelson (2014) Postnov, K. A., & Yungelson, L. R. 2014, Living Reviews in Relativity, 17
  * Raj (2020) Raj, N. 2020, Phys. Rev. Lett., 124, 141802
  * Raskin et al. (2012) Raskin, C., Scannapieco, E., Fryer, C., Rockefeller, G., & Timmes, F. X. 2012, ApJ, 746, 62
  * Robson et al. (2019) Robson, T., Cornish, N. J., & Liu, C. 2019, Classical and Quantum Gravity, 36, 105011
  * Roy et al. (2022) Roy, N., Tiwari, V., Bobrick, A., et al. 2022, arXiv e-prints, arXiv:2204.09683
  * Rueda et al. (2019) Rueda, J., Ruffini, R., Wang, Y., et al. 2019, Journal of Cosmology and Astroparticle Physics, 2019, 044
  * Ruiter et al. (2019) Ruiter, A. J., Ferrario, L., Belczynski, K., et al. 2019, Monthly Notices of the Royal Astronomical Society, 484, 698
  * Sana et al. (2012) Sana, H., de Mink, S. E., de Koter, A., et al. 2012, Science, 337, 444
  * Sato et al. (2015) Sato, Y., Nakasato, N., Tanikawa, A., et al. 2015, The Astrophysical Journal, 807, 105
  * Sato et al. (2016) Sato, Y., Nakasato, N., Tanikawa, A., et al. 2016, ApJ, 821, 67
  * Schwab (2021) Schwab, J. 2021, ApJ, 906, 53
  * Schwab (2021) Schwab, J. 2021, The Astrophysical Journal, 906, 53
  * Seitenzahl et al. (2015) Seitenzahl, I. R., Herzog, M., Ruiter, A. J., et al. 2015, Phys. Rev. D, 92, 124013
  * Shen (2015) Shen, K. J. 2015, The Astrophysical Journal, 805, L6
  * Strumia & Vissani (2003) Strumia, A., & Vissani, F. 2003, Physics Letters B, 564, 42
  * Toonen et al. (2012) Toonen, S., Nelemans, G., & Portegies Zwart, S. 2012, A&A, 546, A70
  * Vitagliano et al. (2020) Vitagliano, E., Tamborra, I., & Raffelt, G. 2020, Reviews of Modern Physics, 92, 045006
  * Wagg et al. (2021) Wagg, T., Breivik, K., & de Mink, S. E. 2021, arXiv e-prints, arXiv:2111.08717
  * Wolfenstein (1978) Wolfenstein, L. 1978, Phys. Rev. D, 17, 2369
  * Wright et al. (2017) Wright, W. P., Kneller, J. P., Ohlmann, S. T., et al. 2017, Phys. Rev. D, 95, 043006
  * Wu et al. (2022) Wu, C., Xiong, H., & Wang, X. 2022, MNRAS, 512, 2972
  * Xiao et al. (2016) Xiao, D., Mészáros, P., Murase, K., & Dai, Z.-G. 2016, ApJ, 832, 20
  * Yu & Jeffery (2010) Yu, S., & Jeffery, C. S. 2010, A&A, 521, A85
  * Yungelson & Kuranov (2016) Yungelson, L. R., & Kuranov, A. G. 2016, Monthly Notices of the Royal Astronomical Society, 464, 1607
  * Zou et al. (2020) Zou, Z.-C., Zhou, X.-L., & Huang, Y.-F. 2020, Research in Astronomy and Astrophysics, 20, 137
  * Zyla et al. (2020) Zyla, P., et al. 2020, PTEP, 2020, 083C01


###### Proof.

That the sectional nonassociativity is nonnegative means that
$L_{\circ}(x\circ x)\geq L_{\circ}(x)^{2}\geq 0$ for all $x\in\mathbb{A}$.
This shows
$\mathbb{Idem}(\mathbb{A},\circ)\subset\mathbb{Q}(\mathbb{A},\circ,h)\subset\mathbb{C}(\mathbb{A},\circ,g)$.
In a Euclidean metrized commutative algebra every multiplication endomorphism
$L_{\circ}(x)$ is self-adjoint, so diagonalizable with real eigenvalues. If
$x,y\in\mathbb{A}$ and $t\in[0,1]$ then
$L_{\circ}(tx+(1-t)y)=tL_{\circ}(x)+(1-t)L_{\circ}(y)\geq 0$ because the set
of nonnegative definite self-adjoint endomorphisms of a Euclidean vector space
is a closed convex cone. ∎

With the hypotheses of Lemma 8.7, are the extreme rays of
$\square(\mathbb{A},\circ)$ generated by elements of
$\mathbb{Idem}(\mathbb{A},\circ)$? This is true for the cones of squares in
the simple Euclidean Jordan algebras over $\mathbb{R}$ or $\mathbb{C}$; these
are the cones of nonnegative definite Hermitian matrices and their extremal
rays are generated by the rank one idempotents. For the Griess algebras of OZ
VOAs (Example 19) it seems the question has not been studied.

## References

  * [1] A. A. Albert, _On the nuclei of a simple Jordan algebra_ , Proc. Nat. Acad. Sci. U.S.A. 50 (1963), 446–447.
  * [2] B. N. Allison, _A class of nonassociative algebras with involution containing the class of Jordan algebras_ , Math. Ann. 237 (1978), no. 2, 133–156.
  * [3] T. Ashihara and M. Miyamoto, _Deformation of central charges, vertex operator algebras whose Griess algebras are Jordan algebras_ , J. Algebra 321 (2009), no. 6, 1593–1599.
  * [4] K. M. R. Audenaert, _Variance bounds, with an application to norm bounds for commutators_ , Linear Algebra Appl. 432 (2010), no. 5, 1126–1143.
  * [5] A. M. Bloch and A. Iserles, _Commutators of skew-symmetric matrices_ , Internat. J. Bifur. Chaos Appl. Sci. Engrg. 15 (2005), no. 3, 793–801.
  * [6] A. Böttcher and D. Wenzel, _How big can the commutator of two matrices be and how big is it typically?_ , Linear Algebra Appl. 403 (2005), 216–228.
  * [7] by same author, _The Frobenius norm and the commutator_ , Linear Algebra Appl. 429 (2008), no. 8-9, 1864–1885.
  * [8] M. Bremner and V. Dotsenko, _Associator dependent algebras and Koszul duality_ , Ann. Mat. Pura Appl. (2022), Advance online publication.
  * [9] R. B. Brown and A. Gray, _Vector cross products_ , Comment. Math. Helv. 42 (1967), 222–236.
  * [10] C.-M. Cheng, S.-W. Vong, and D. Wenzel, _Commutators with maximal Frobenius norm_ , Linear Algebra Appl. 432 (2010), no. 1, 292–306.
  * [11] S. S. Chern, M. do Carmo, and S. Kobayashi, _Minimal submanifolds of a sphere with second fundamental form of constant length_ , Functional Analysis and Related Fields (F. E. Browder, ed.), Springer, New York, 1970, pp. 59–75.
  * [12] J. H. Conway, _A simple construction for the Fischer-Griess monster group_ , Invent. Math. 79 (1985), no. 3, 513–540.
  * [13] B. Eckmann, _Stetige Lösungen linearer Gleichungssysteme_ , Comment. Math. Helv. 15 (1943), 318–339.
  * [14] A. Elduque, _Composition algebras_ , arXiv:1810.09979.
  * [15] A. Elduque and H. C. Myung, _On flexible composition algebras_ , Comm. Algebra 21 (1993), no. 7, 2481–2505.
  * [16] by same author, _Composition algebras satisfying the flexible law_ , Comm. Algebra 32 (2004), no. 5, 1997–2013.
  * [17] J. Faraut and A. Korányi, _Analysis on symmetric cones_ , Oxford Mathematical Monographs, The Clarendon Press Oxford University Press, New York, 1994.
  * [18] D. J. F. Fox, _Commutative algebras with nondegenerate invariant trace form and trace-free multiplication endomorphisms_ , arXiv:2004.12343.
  * [19] by same author, _Geometric structures modeled on affine hypersurfaces and generalizations of the Einstein-Weyl and affine sphere equations_ , Extended abstracts Fall 2013—geometrical analysis, type theory, homotopy theory and univalent foundations, Trends Math. Res. Perspect. CRM Barc., vol. 3, Birkhäuser/Springer, Cham, 2015, pp. 15–19.
  * [20] by same author, _Harmonic cubic homogeneous polynomials such that the norm-squared of the Hessian is a multiple of the Euclidean quadratic form_ , Anal. Math. Phys. 11 (2021), no. 1, 43.
  * [21] E. Frenkel and D. Ben-Zvi, _Vertex algebras and algebraic curves_ , Mathematical Surveys and Monographs, vol. 88, American Mathematical Society, Providence, RI, 2001.
  * [22] I. Frenkel, J. Lepowsky, and A. Meurman, _Vertex operator algebras and the Monster_ , Pure and Applied Mathematics, vol. 134, Academic Press Inc., Boston, MA, 1988.
  * [23] H. Freudenthal, _Oktaven, Ausnahmegruppen und Oktavengeometrie_ , Geom. Dedicata 19 (1985), no. 1, 7–63.
  * [24] J. Ge, F. Li, and Y. Zhou, _Some generalizations of the DDVV and BW inequalities_ , Trans. Amer. Math. Soc. 374 (2021), no. 8, 5331–5348.
  * [25] R. W. Gebert, _Introduction to vertex algebras, Borcherds algebras and the Monster Lie algebra_ , Internat. J. Modern Phys. A 8 (1993), no. 31, 5441–5503.
  * [26] M. Goze and E. Remm, _Lie-admissible algebras and operads_ , J. Algebra 273 (2004), no. 1, 129–152.
  * [27] by same author, _Poisson algebras in terms of non-associative algebras_ , J. Algebra 320 (2008), no. 1, 294–317.
  * [28] A. A. Ivanov, _The Monster group and Majorana involutions_ , Cambridge Tracts in Mathematics, vol. 176, Cambridge University Press, Cambridge, 2009.
  * [29] A.A. Ivanov, D.V. Pasechnik, Á. Seress, and S. Shpectorov, _Majorana representations of the symmetric group of degree 4_ , Journal of Algebra 324 (2010), no. 9, 2432 – 2463.
  * [30] N. Jacobson, _Structure and representations of Jordan algebras_ , American Mathematical Society Colloquium Publications, Vol. XXXIX, American Mathematical Society, Providence, R.I., 1968.
  * [31] V. Kac, _Vertex algebras for beginners_ , 2 ed., University Lecture Series, vol. 10, American Mathematical Society, Providence, RI, 1998.
  * [32] M.-A. Knus, A. Merkurjev, M. Rost, and J.-P. Tignol, _The book of involutions_ , American Mathematical Society Colloquium Publications, vol. 44, American Mathematical Society, Providence, RI, 1998.
  * [33] M. Koecher, _The Minnesota notes on Jordan algebras and their applications_ , Lecture Notes in Mathematics, vol. 1710, Springer-Verlag, Berlin, 1999.
  * [34] F. Kosier, _On a class of nonflexible algebras_ , Trans. Amer. Math. Soc. 102 (1962), 299–318.
  * [35] C. H. Lam, _On VOA associated with special Jordan algebras_ , Comm. Algebra 27 (1999), no. 4, 1665–1681.
  * [36] T. Y. Lam, _Introduction to quadratic forms over fields_ , Graduate Studies in Mathematics, vol. 67, American Mathematical Society, Providence, RI, 2005.
  * [37] H. T. Laquer, _Invariant affine connections on Lie groups_ , Trans. Amer. Math. Soc. 331 (1992), no. 2, 541–551.
  * [38] Z. Lu, _Normal scalar curvature conjecture and its applications_ , J. Funct. Anal. 261 (2011), no. 5, 1284–1308.
  * [39] Z. Lu and D. Wenzel, _The normal Ricci curvature inequality_ , Recent advances in the geometry of submanifolds—dedicated to the memory of Franki Dillen (1963–2013), Contemp. Math., vol. 674, Amer. Math. Soc., Providence, RI, 2016, pp. 99–110.
  * [40] M. Markl and E. Remm, _Algebras with one operation including Poisson and other Lie-admissible algebras_ , J. Algebra 299 (2006), no. 1, 171–189.
  * [41] A. Matsuo and K. Nagatomo, _Axioms for a vertex algebra and the locality of quantum fields_ , MSJ Memoirs, vol. 4, Mathematical Society of Japan, Tokyo, 1999.
  * [42] W. Meyer and W. Neutsch, _Associative subalgebras of the Griess algebra_ , J. Algebra 158 (1993), no. 1, 1–17.
  * [43] M. Miyamoto, _Griess algebras and conformal vectors in vertex operator algebras_ , J. Algebra 179 (1996), no. 2, 523–548.
  * [44] H. C. Myung, _Some classes of flexible Lie-admissible algebras_ , Trans. Amer. Math. Soc. 167 (1972), 79–88.
  * [45] by same author, _Mal ′cev-admissible algebras_, Progress in Mathematics, vol. 64, Birkhäuser Boston, Inc., Boston, MA, 1986.
  * [46] S. Norton, _The Monster algebra: some new formulae_ , Moonshine, the Monster, and related topics (South Hadley, MA, 1994), Contemp. Math., vol. 193, Amer. Math. Soc., Providence, RI, 1996, pp. 297–306.
  * [47] S. Okubo, _Pseudo-quaternion and pseudo-octonion algebras_ , Hadronic J. 1 (1978), no. 4, 1250–1278.
  * [48] by same author, _Introduction to octonion and other non-associative algebras in physics_ , Montroll Memorial Lecture Series in Mathematical Physics, vol. 2, Cambridge University Press, Cambridge, 1995.
  * [49] S. Okubo and H. C. Myung, _Some new classes of division algebras_ , J. Algebra 67 (1980), no. 2, 479–490.
  * [50] by same author, _Adjoint operators in Lie algebras and the classification of simple flexible Lie-admissible algebras_ , Trans. Amer. Math. Soc. 264 (1981), no. 2, 459–472.
  * [51] S. Okubo and J. M. Osborn, _Algebras with nondegenerate associative symmetric bilinear forms permitting composition_ , Comm. Algebra 9 (1981), no. 12, 1233–1261.
  * [52] G. Ovando, _Two-step nilpotent Lie algebras with ad-invariant metrics and a special kind of skew-symmetric maps_ , J. Algebra Appl. 6 (2007), no. 6, 897–917.
  * [53] H. P. Petersson, _Eine Identität fünften Grades, der gewisse Isotope von Kompositions-Algebren genügen_ , Math. Z. 109 (1969), 217–238.
  * [54] E. Remm and K. Morand, _Weakly associative algebras, Poisson algebras and quantization_ , arXiv:2005.04273.
  * [55] H. Röhrl and M. B. Wischnewsky, _Subalgebras that are cyclic as submodules_ , Manuscripta Math. 19 (1976), no. 2, 195–209.
  * [56] A. A. Sagle, _Malcev algebras_ , Trans. Amer. Math. Soc. 101 (1961), 426–458.
  * [57] R. D. Schafer, _On the algebras formed by the Cayley-Dickson process_ , Amer. J. Math. 76 (1954), 435–446.
  * [58] by same author, _An introduction to nonassociative algebras_ , Pure and Applied Mathematics, Vol. 22, Academic Press, New York, 1966.
  * [59] D. B. Shapiro, _Compositions of quadratic forms_ , De Gruyter Expositions in Mathematics, vol. 33, Walter de Gruyter & Co., Berlin, 2000.
  * [60] H. Shima, _The geometry of Hessian structures_ , World Scientific Publishing Co. Pte. Ltd., Hackensack, NJ, 2007.
  * [61] H. Suzuki, _Automorphism groups of multilinear maps_ , Osaka J. Math. 20 (1983), no. 3, 659–673.
  * [62] È. B. Vinberg, _The theory of homogeneous convex cones_ , Trudy Moskov. Mat. Obšč. 12 (1963), 303–358.
  * [63] by same author, _Invariant norms in compact Lie algebras_ , Funkcional. Anal. i Priložen 2 (1968), no. 2, 89–90.
  * [64] M. Whybrow, _An infinite family of axial algebras_ , J. Algebra 577 (2021), 1–31.
  * [65] Y.-D. Wu and X.-Q. Liu, _A short note on the Frobenius norm of the commutator_ , Math. Notes 87 (2010), no. 6, 903–907.
  * [66] X. Zhan, _Matrix theory_ , Graduate Studies in Mathematics, vol. 147, American Mathematical Society, Providence, RI, 2013.
  * [67] H. Zhao, _On vertex operator algebras associated to Jordan algebras of Hermitian type_ , Comm. Algebra 48 (2020), no. 1, 331–350.
  * [68] A. Zohrabi and P. Zusmanovich, _On Hermitian and skew-Hermitian matrix algebras over octonions_ , J. Nonlinear Math. Phys. 28 (2021), no. 1, 108–122.


# On-sky measurements of atmospheric dispersion:
II. Atmospheric models characterization

B. Wehbe,1,2 A. Cabral,3,4 L. Sbordone5, and G.Ávila
1Instituto de Astrofísica e Ciências do Espaço, Universidade do Porto, CAUP,
Rua das Estrelas, 4150-762 Porto, Portugal
2Departamento de Física e Astronomia, Faculdade de Ciências, Universidade do
Porto, Rua Campo Alegre, 4169-007 Porto, Portugal
3Instituto de Astrofísica e Ciências do Espaço, Universidade de Lisboa, Campus
do Lumiar, Estrada do Paço do Lumiar 22, Edif. D, PT1649-038 Lisboa, Portugal
4Departamento de Física, Faculdade de Ciências, Universidade de Lisboa, Campo
Grande 1749-016 Lisboa Portugal
5European Southern Observatory, Alonso de Córdova 3107, Vitacura, Región
Metropolitana, Chile E-mail<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

Differential atmospheric dispersion is a wavelength-dependent effect
introduced by Earth’s atmosphere that affects astronomical observations
performed using ground-based telescopes. It is important, when observing at a
zenithal angle different from zero, to use an Atmospheric Dispersion Corrector
(ADC) to compensate this atmospheric dispersion. The design of an ADC is based
on atmospheric models that, to the best of our knowledge, were never tested
against on-sky measurements. We present an extensive models analysis in the
wavelength range of 315-665 nm. The method we used was previously described in
the paper I of this series. It is based on the use of cross-dispersion
spectrographs to determine the position of the centroid of the spatial profile
at each wavelength of each spectral order. The accuracy of the method is 18
mas. At this level, we are able to compare and characterize the different
atmospheric dispersion models of interest. For better future ADC designs, we
recommend to avoid the Zemax model, and in particular in the blue range of the
spectra, when expecting residuals at the level of few tens of milli-
arcseconds.

###### keywords:

atmospheric effects - instrumentation: spectrographs - methods: data analysis

††pubyear: 2020††pagerange: On-sky measurements of atmospheric dispersion:
II. Atmospheric models characterization–On-sky measurements of atmospheric
dispersion:  II. Atmospheric models characterization

##  1 Introduction

Differential atmospheric dispersion is a wavelength-dependent effect
introduced by Earth’s atmosphere that affects astronomical observations
performed using ground-based telescopes. It is important to take this
atmospheric dispersion into consideration, when observing at a zenithal angle
different from zero. The common way to compensate this dispersion, is to
implement an Atmospheric Dispersion Corrector (ADC) to compensate the on-sky
dispersion. Most of ADCs are counter-rotating prisms used to correct the
effects of atmospheric dispersion (Cabral et al., 2012). Each prism is a set
of two prisms at least, made from different glasses. The choice of glasses is
based on their refractive index, and their ability to reproduce the same
amount of dispersion that needs to be corrected. The design of an ADC is based
on atmospheric models used to estimate the on-sky dispersion to be corrected.
During this phase, several pairs of glasses are chosen and tested. The final
selection of glass is usually based on the minimum amount of residuals (on-sky
dispersion - ADC dispersion) they can provide. This can be considered a good
function of merit if the level of accuracy of the atmospheric dispersion
models is lower than the level of residuals desired. High-precision
astronomical instruments, like state-of-the-art spectrographs, are now
equipped with ADCs that are supposed to deliver residuals at the level of few
tens of milli-arcseconds (mas). Having residuals below 100 milli-arcseconds
(mas), for objects that have a full width at half maximum (FWHM) of 1000 mas
is quite challenging (Wehbe et al., 2020a). We notice that the difference
between the most used atmospheric models, Filippenko (1982), hereafter
Filippenko, and Hohenkerk & Sinclair (1985), hereafter Zemax, can reach $\sim$
50 mas, the limit of residuals needed (see Figure 1). Designing ADCs to return
residuals at the same level of the difference between the models will
introduce errors that will affect flux losses and the final radial velocity
(RV) computed as shown in Wehbe et al. (2020a). Therefore, it is important to
characterize different atmospheric models and quantify their accuracy in terms
of reproducing the on-sky dispersion down to the level of few mas.
Several attempts to theoretically compare different atmospheric models can be
found in the literature (e.g., Spanò, 2014). To the best of our knowledge,
these models have never been tested on-sky. In the first paper of this series
(Wehbe et al., 2020b), we presented a method to measure on-sky the atmospheric
dispersion. In Section 2, we describe the models we are interested in. We also
present a qualitative analysis of the different models. The method is
summarized in Section 3, together with the different observations used in this
work. The data analysis is presented in Section 4, and the conclusions in
Section 5.

Figure 1: Differential atmospheric dispersion at Paranal for different models
between 315 nm and 665 nm. The average between the models was used as
reference.

##  2 Models

The effects of atmospheric dispersion are well known and have been extensively
discussed by several authors (e.g., Filippenko, 1982; Skemer et al., 2009;
Bechter et al., 2018; Blackman et al., 2019; Wehbe et al., 2020a). Several
attempts were done to measure the refractivity of dry air, $n$ , using an air
refractometer (for example Bönsch & Potulski, 1998, hereafter B&P) , and
derive a formula to compute the refractivity as a function of wavelength.
Using the air refractivity formula, we can derive the amount of atmospheric
dispersion using the following:

$\displaystyle\Delta R(\lambda)$ $\displaystyle=R(\lambda)-R(\lambda_{\rm
ref})$ $\displaystyle\Delta R(\lambda)$ $\displaystyle\approx
206265\left[n(\lambda)-n(\lambda_{\rm ref})\right]\times\tan Z,$ (1)

where $R(\lambda)$ is the refraction angle, $\lambda_{\rm ref}$ is the
reference wavelength, and $Z$ is the zenithal angle of observation. Since
observatories are usually located at an altitude higher than sea level, the
atmospheric parameters, such as temperature, pressure and relative humidity
(RH), are also taken into consideration while deriving the air refractivity
formula. Different attempts to measure the refractivity using in-lab
measurements, give different results. This difference is minimized when taking
into account the variation of pressure, temperature and RH. Another important
difference is the limits of the wavelength range where their formulaf are
trustable.
In this section, we will briefly describe the different models of interest and
compare them in a qualitative way.

###  2.1 Edlen

A dispersion formula for standard air (15∘ and 760 mbar) was derived by Edlen
(1953) based on measurements done by Barrell & Sears (1939). In 1966, an
improved version of the dispersion formula was provided by Edlén (1966),
hereafter Edlen, that took into account the carbon dioxide and water vapor.
The main dispersion formula is given by equation 2, followed by the effects of
carbon dioxide, temperature and pressure, and the water vapor, respectively.

$(n-1)_{\rm s}\times
10^{8}=8342.13+\frac{2406030}{130-(1/\lambda)^{2}}+\frac{15997}{38.9-(1/\lambda)^{2}}$
(2)

$(n-1)_{\rm x}=[1+0.54(x-0.0003)](n-1)_{\rm s}$ (3)

$(n-1)_{\rm tp}=\frac{p(n-1)_{\rm
s}}{720.775}\times\frac{1+p(0.817-0.0133t)\times 10^{-6}}{1+0.003661t}$ (4)

$n_{\rm tpf}=n_{\rm tp}-f(5.722-0.0457(1/\lambda)^{2})\times 10^{-8},$ (5)

where $x$ is the percentage of carbon dioxide, $f$ is the water vapor
pressure, $t$ is the temperature, and $p$ is the atmospheric pressure.
Combining equations 2, 3, 4, and 5 with equation 1, result with the
atmospheric dispersion using Edlen’s formula.

###  2.2 Filippenko

Filippenko (1982) introduced a new set of formula to compute the air
refractivity of standard air (15∘ and 760 mbar). He started using the 1953
Edlen’s formula (Edlen, 1953) and updated the dependence on pressure,
temperature, and water vapor (Equations 6 to 8).

$n_{\rm f}=f\times\frac{0.0624-0.00068(1/\lambda)^{2}}{1+0.003661t}$ (6)

$(n-1)_{\rm s}\times
10^{6}=64.328+\frac{29498.1}{146-(1/\lambda)^{2}}+\frac{255.4}{41-(1/\lambda)^{2}}-n_{\rm
f}$ (7)

$(n-1)_{\rm tp}=(n-1)_{\rm s}\times\frac{p[1+(1.049-0.0157t)p\times
10^{-6}]}{720.883(1+0.003661t)}$ (8)

###  2.3 Bönsch & Potulski

Bönsch & Potulski (1998) updated the Edlen’s formula to adapt the change of
systems where the standard air is now considered at 20∘ and 750 mbar. They
also built an air refractometer and were able to measure the refractivity.
They refined their formula using the best fit of their results. In a similar
way, we will list the main equations of B&P to compute the atmospheric
dispersion (Equations 9 to 12).

$(n-1)_{\rm s}\times
10^{8}=8091.37+\frac{2333983}{130-(1/\lambda)^{2}}+\frac{15518}{38.9-(1/\lambda)^{2}}$
(9)

$(n-1)_{\rm x}=[1+0.5327(x-0.0004)](n-1)_{\rm s}$ (10)

$(n-1)_{\rm
tp}=\frac{p(n-1)_{s}}{93214.6}\times\frac{1+p(0.5953-0.009876t)\times
10^{-8}}{1+0.003661t}$ (11)

$n_{\rm tpf}=n_{\rm tp}-f(3.802-0.0384(1/\lambda)^{2})\times 10^{-10},$ (12)

###  2.4 Zemax

The Zemax model of atmospheric dispersion is based on the work of Hohenkerk &
Sinclair (1985). The main refractivity formula is based on the one published
by Barrell & Sears (1939). The dependence on temperature, pressure, and RH, is
computed in a complex way to take into account all the variations in the
troposphere and stratosphere layers of the atmosphere. Therefore, we will only
mention the main refractivity formula here (Equation 13). For more details,
the reader may refer to Barrell & Sears (1939), and Hohenkerk & Sinclair
(1985).

$(n-1)_{\rm s}\times
10^{8}=27263.24+\frac{154.4}{\lambda^{2}}+\frac{1.29}{\lambda^{4}}$ (13)

###  2.5 Zemax updated

Spanò (2014) introduced an improved version of the Zemax atmospheric
dispersion model. The Zemax model underestimates the atmospherisc dispersion
in the UV and the in the blue part of the visible spectrum by 0.07 %.
Therefore, the author introduced a thin prism before the Zemax model to better
describe the air refraction with a “gradient 5" surface. The prism apex angle
is scaled as a function of the observed zenithal angle Z (see Equation 14).
The reader may refer to Spanò (2014) for more details.

$\tan\delta=\frac{1}{G}\left[A\tan Z+B(\tan Z)^{3}\right],$ (14)

where G is a scale factor used to define the prism index N, $\rm
A~{}=~{}1.00472$, and B = -0.00103.

###  2.6 Qualitative analysis

First, it is important to qualitatively address the behavior of each model.
This will help us determine the sensitivity of each model in terms of
atmospheric parameters inputs: temperature, pressure, and relative humidity.
It is clear from Figure 1, that all of the models behave in a similar way
except for Zemax. In order to test the sensitivity of the models, we compare
the results of a set of two models by varying the inputs as follow: i)
temperature range from -5 to 25∘C; ii) pressure range from 700 to 800 mbar;
and iii) relative humidity from 0 to 70 %. We also assume the following
default inputs: T = 13∘C, P = 740 mbar, RH = 20%, and Z = 60∘. These ranges
are the same ranges where the VLT telescopes operate. For demonstration
purposes, we will show the comparison between the B & P model as well as the
Zemax model with respect to the Filippenko model (the other models return
similar behavior as the ones presented in this subsection).

####  2.6.1 B & P vs Filippenko

The left column of Figure 2, represents the results of the comparison between
the B & P and Filippenko models regarding the temperature, pressure, and
relative humidity variation. From these figures, we notice that the general
difference between these models, due to the main refractivity formula
(Equations 7, 9), is negligible ($\sim$ 2.5 mas). We also notice that they
have similar sensitivities in terms of temperature and pressure, while there
is a variation of 1 mas (in the bluest part of the spectrum) due to the RH
variation.

####  2.6.2 Zemax vs Filippenko

The right column of Figure 2, represents the results of the comparison between
the Zemax and Filippenko models regarding the temperature, pressure, and
relative humidity. It is clear that the sensitivities of the models are
similar, and the general difference is due to the main refractivity formula
(Equations 7, 13).

Figure 2: Left column: difference (in dispersion from the models) between the
B & P and Filippenko assuming the default inputs parameters: P = 740 mbar, RH
= 20%, and Z = 60∘. The temperature varies from -5 to 25 ∘C as shown in the
colorbar. The wavelength reference is assumed to be 490 nm; right: similar to
left but comparing Zemax with Filippenko. We note that the scale on the y-axis
between for B&P and Zemax, is different.

##  3 Method & Observations

###  3.1 Method

The method we developed to measure the atmospheric dispersion on-sky, is
described in details in Wehbe et al. (2020b). We limit ourselves here to
mention the main steps only. We used the cross-dispersed slit spectrograph
UVES (Dekker et al., 2000), with the slit directed along the dispersion
direction and we removed the ADC from the optical path. This setup, allowed us
to directly measure the atmospheric dispersion as a function of wavelength.
The method is based on the determination of the position of the centroid of
the spatial profile at each wavelength of each spectral order. Using the flat
field frames as a reference point, we overlap the science frames and by
applying cuts perpendicular to the spectral orders, we extract the centroids
using a Gaussian and a Sigmoid function (for more details see Wehbe et al.,
2020b). In Figure 2 of Wehbe et al. (2020b), we represent the overlapping of
the two frames as well as the centroids extraction. By repeating the same
procedure to all the spectral orders, through all the cuts, we get a plot that
represents the variation of atmospheric dispersion in pixel as a function of
wavelength. By using one of the two methods explained in Wehbe et al. (2020b),
we convert the atmospheric dispersion from pixels to mas in the sky. Using the
method, we also were able to detect a non-negligible instrumental dispersion
that should be corrected from the final results. The method shows an accuracy
(see Table 1), at the 1 $\sigma$ confidence level, of $\pm$ 18 mas in the blue
and red ranges of UVES (315 to 665 nm).

Table 1: Uncertainties budget of the method, showing the main sources of errors in on-sky measurements (adopted from Table 3 of Wehbe et al., 2020b). Source of uncertainty | Contribution
---|---
Pixel scale | 0.87 mas
Instrumental dispersion | 8.30 mas
Atmospheric dispersion | 16.18 mas
Accuracy of the method: 18 mas

###  3.2 Observations

In order to perform an analysis over the optical range of the spectrum, we
obtained observation from UVES, using a dichroic to cover the wavelength range
of 303 nm to 665 nm. Our data was acquired based on observations made with
UVES (program ID 4103.L-0942(A); PI: B. Wehbe). The spectrograph is equipped
with two arms: i) blue in the 303-384 nm range, ii) red in the 487-567 nm and
590-665 nm (the red arm detector is a mosaic of two CCDs); the wavelength
ranges are following our setup using a dichroic. In all observations, the slit
of the spectrograph was oriented along the atmospheric dispersion direction,
the ADC was off (outside of the optical path), and using dichroic 1 along with
the two cross-dispersers (CD) 1 & 3\. In order to achieve the desired accuracy
of $\sim$ 1 % (on the dispersion measurements), 19 exposures of 20 seconds
each were done for each target. We observed 3 targets between zenithal angles
of 50∘ and 60∘, and one target between 14∘ and 19∘. Table 2 summarizes the
atmospheric conditions at the time of observation of each target. The targets
chosen are O-type stars (HD17490 and HD150574), and B-type stars (HD143449 and
HD165320). We choose O and B-type stars in order to achieve a higher signal-
to-noise (S/N) ratios in particular in the blue range of the spectrum where
the atmospheric dispersion is more severe. We also state that another
important factor when choosing the targets was the zenithal angle at which it
can be observed. The higher the angle, the bigger the atmospheric dispersion
is. Two of these targets (HD117490 and HD143449) were observed with increasing
zenithal angles between the first and last exposures. The other two targets
(HD150574 and HD165320) were observed with decreasing zenithal angles between
the first and last exposures. This is in particular importance because the
atmospheric dispersion is in opposite directions between increasing and
decreasing zenithal angles; an information important to take into
consideration while correcting from the instrumental dispersion as it is
considered as a fixed dispersion in quantity and direction.

Table 2: Atmospheric conditions at the time of observation of each target. Target | Zenithal angle | Temperature | Pressure | RH
---|---|---|---|---
| ∘ |  ∘C | mbar | %
HD 117490 | 57.4 - 59.9 | 11.7 - 12.3 | 742.6 | 4.5 - 6
HD 150574 | 55.9 - 59.4 | 13.4 | 745 | 13.5 - 17
HD 143449 | 50 - 52.3 | 9.5 - 9.7 | 741.8 | 18.5 - 23
HD 165320 | 13.8 - 18.5 | 13.6 -13.8 | 744.7 | 14.5 - 16.5

##  4 Data Analysis

As mentioned before, UVES is equipped with two arms: a blue one and a red one.
Unfortunately, these two ranges do not overlap in our setup. This means that
we can’t evaluate the full range from 303 nm to 665 nm in one exposure.
Therefore, we will divide our analysis in two subsections focusing on each arm
separately. In order to reach the accuracy needed, we required 19 exposures of
each target. It is important, when analyzing the results, to remove any
outliers from the measurements. We first subtract the mean value from the
residuals of each exposure. This will allow us to center the residuals around
zero. This in fact is similar to using the guiding camera to center the target
in the fiber. Then, we compute the half-peak-to-valley (HPTV) of the residuals
of each exposure, for all the models of interest. Using the HPTV of all the
exposures, we can identify the outliers using the Z-score method (Equation
15). Any exposure that return a Z-score higher than the threshold (Equation
16), is considered as an outlier. Note that we repeat this procedure for each
model. In the following part of the paper, we only take into account the
exposures that return a Z-score below the threshold: all the exposures of
HD117490 and HD165320, 11 exposures of HD150574, and 11 exposures of HD143449.

$Z=\frac{x-\mu}{\sigma},$ (15)

where x is the HPTV of each exposure, $\rm\mu$ is the PV average of all the
exposures, and $\rm\sigma$ is the standard deviation.

$t=\frac{n-1}{\sqrt{n}},$ (16)

where n is the number of exposures of each target (19).

###  4.1 The blue arm

Similar to paper I, in our analysis we use the blue range starting 315 nm,
since the points below are affected by low atmospheric transmission and low
S/N ratios. In Figure 3, we show the case of one exposure of the target
HD117490. Panel a of Figure 3 shows the measured atmospheric dispersion as
well as the computed one using the Filippenko’s model. In the other panels (b
to f) of Figure 3, we show the residuals (i.e., data - model) for all the
models of interest. The different models show similar behavior in terms of
residuals, except for the Zemax model which seems to be underestimating the
atmospheric dispersion in the bluest part of the range. To investigate more,
we should look at the results of the 19 exposures of each target. After
removing the outlier exposures, we plot the residuals of all the models. The
Edlen model is excluded from these plots as it return similar results to the B
& P model. Figure 4 represents the results. It is clear from Figure 4, that
the models B & P, Filippenko, and Zemax updated, are returning similar
residuals for the 4 targets. It is also clear that they are able to return
residuals between -10 and 10 mas, for all the zenithal angles tested (the 4
targets). Also from Figure 4, we can see that the Zemax model is
underestimating the atmospheric dispersion, especially in the bluest part of
the spectrum. This is particularly clear in the highest zenithal angles cases
(HD117490 and HD150574), where the dispersion is more severe. The updated
version of the Zemax model introduced by Spanò (2014), is taking these
underestimation into consideration and returning better results in terms of
residuals. The Zemax model should be avoided when dealing with the blue range,
especially when expecting residuals below the 20 mas level. At lower zenithal
angles, hence lower atmospheric dispersion (the case of HD165320 in Figure 4),
the 4 models return similar results between -10 and 10 mas.

Figure 3: Panel a: measured atmospheric dispersion for one exposure of
HD117490, overplotted by the expected dispersion from the Filippenko’s model;
panels b to f: atmospheric dispersion residuals for all the modesl of
interest. Figure 4: Residuals of all the exposures of all the targets for all
the models of interest in the blue arm (315 - 385 nm). Each row represents a
target, hence a zenithal angle. Each column represents a different atmospheric
model. The two dashed lines represent the accuracy of the method: $\pm$ 18
mas.

###  4.2 The red arm

As it was mentioned before, the red arm detector is a mosaic of two CCDs
resulting in two separate ranges: 487-567 nm, and $\rm 590-665~{}nm$. In this
subsection, we will show the results of the red arm in a similar way as
presented in subsection 4.1. In Figures 5, and 6, we show the case of one
exposure of the target HD117490. Panel a of both figures shows the measured
atmospheric dispersion as well as the computed one using the Filippenko’s
model, for the corresponding range. In the other panels (b to f), we show the
residuals for all the models of interest. The different models show similar
behavior in terms of residuals. In a similar manner to subsection 4.1, we plot
the residuals of all the exposures of all the targets. These residuals are
shown in Figures 8 and 9 for the wavelength ranges of $\rm 487-567~{}nm$ and
590-665 nm, respectively.

Figure 5: Similar to Figure 3, for the wavelength range 487-567 nm. Figure 6:
Similar to Figure 3, for the wavelength range 590-665 nm.

####  4.2.1 The red range of 487-567 nm

The residuals of all the models of all the targets are within the uncertainty
level as seen in Figure 8. In a first analysis, we can consider them
negligible. Nevertheless, looking at Figure 8 and not taking into
consideration the uncertainty level, a slope in the residuals is clear. As
this effect is withing our uncertainty level, it does not affect our results.
However, we tried to investigate the possible cause.
This effect is model independent which rule out the idea of an error in the
models. It is also target independent. The four targets were observed at
different days, different zenithal angles, but with the same optical setup of
the spectrograph. This indicates that it is an instrumental effect, and not
target related, that we were not able to detect using the airmass 1 data.
One possible cause to explain this slope could be the derotator, affecting our
calibration between the science and flat frames used. To investigate this
effect, we requested a set of flat frames111This test was executed as an
instrumental test and as such it’s not available from the ESO archive. The
data is available upon request from the author. while rotating the derotator
of UVES between 0∘ and 360∘, with steps of 30∘. The angles of 180∘ and 210∘
were excluded from the set due to technical issues. By measuring the
atmospheric dispersion using different flat frames with different derotator
positions, we were able to quantify the effect of the derotator on the
instrumental dispersion. We use the frame of 0∘ as a reference point. Figure 7
shows the effect of the derotator position on the measured atmospheric
dispersion.

Figure 7: Effect of derotator angle on the measured atmospheric dispersion.
The frame at derotator position 0∘ is used as a reference frame for the
comparison. Top: for the blue range of 315-384 nm, middle: for the red range
of 487-567 nm, bottom: for the red range of 590-665 nm.

The effect of the derotator is at the level of $\pm$ 2 mas in the blue range,
and $\pm$ 1 mas in the red range (as seen in Figure 7) which rule out the
responsibility of the derotator. Another possible cause could be a small miss-
alignment variation in the optics. It is important to state here that this
effect is so small that it will not affect the science observations of UVES.
However, since we are trying to measure residuals at the level of few mas, we
are able to detect it. The set of the flat frames was taken in two nights (27
and 28 of December 2020) where there was an earthquake near Paranal where UVES
is installed. By comparing the zero flat frame before and after the
earthquake, we can quantify a linear shift in the dispersion at the level of
10 mas in the blue range and 3 mas in the red range. This shift, even though
at a different level to what we are seeing in Figure 8, shows a similar
behavior. This might explain the residuals linear behavior. In fact, the
airmass 1 set of data used to quantify the instrumental dispersion was taken
in 2015, while our observations were done in 2019. Any small misalignment
variation on the red mosaic chip, or even an earthquake, might introduce such
an effect. Another possible cause could be related to thermal variations from
observation to observation that might create some offsets on the data.
Unfortunately, we were not able to find a single cause that can confirm the
presence of this slope. We suspect that it could be due to a sum of several
small effects, that don’t affect the UVES science observations, but clear in
our analysis.
Even though there is an extra dispersion that it is unclear for us from where,
we are confident from our results as this slope effect is within our
uncertainty level. It has no impact on our conclusions. The residuals of all
the targets are between the -10 and 10 mas, similar to the blue range.

####  4.2.2 The red range of 590-665 nm

The residuals of the range 590-665 nm are between -5 and 5 mas for all the
targets and all the models, as seen in Figure 9. The residuals are smaller in
this range from the fact that the atmospheric dispersion is smaller and the
difference between the models at this level is almost negligible. All the
models are able to return residuals this low, independent of the zenithal
angles. The slope detected in HD143449 in Figure 9, is due to errors in our
measurements. They are at the level few mas and can be considered negligible.
The periodicity seen in the case of HD165320 is due to the order separation.
It is particularly seen in this range since the level of atmospheric
dispersion, hence the residuals, is really small compared to the other red
range and to the blue range.

Figure 8: Similar to Figure 4, but for the wavelength range of 487-567 nm.
Figure 9: Similar to Figure 4, but for the wavelength range of 590-665 nm.

##  5 Conclusions

We present an analysis of different atmospheric dispersion models used during
the design phase of an ADC. Using the method developed by Wehbe et al.
(2020b), we were able to compare the models of interest to on-sky
measurements. To the best of our knowledge, this is the first time that
atmospheric models are being compared to on-sky measurements, in particular in
the blue part of the spectra, where the atmospheric dispersion is more severe.
Even though the method we used has an accuracy of 18 mas, we were able to
compare the different models. We can notice that in the blue part of the
range, the Zemax model is returning higher residuals than the other models.
While in the red ranges, the results are similar. From this paper, we find the
following:

  1. 1.

The validation of the method is successfully extended from paper I to four
different targets covering different ranges of atmospheric parameters and in
particular different zenithal angles.

  2. 2.

It is clear that the Zemax updated model is better than the Zemax model, in
terms of estimating the on-sky dispersion. We conclude this from the fact that
the residuals are lower, in particular in the blue range of the spectra where
it is well known by now that the Zemax model is underestimating the
atmospheric dispersion.

  3. 3.

In the blue range of the spectra, and when expecting residuals at the level of
few tens of mas, we recommend to avoid using the Zemax model.

  4. 4.

In the red range of the spectra, all the models we tested return similar
results. Any of the models can be safely used when expecting residuals at the
level of few tens of mas.

  5. 5.

The effect seen in Figure 8, is to our estimation, due to an extra
instrumental dispersion not detected in the data of airmass 1. Even that the
amplitude is below our method accuracy, the effect is clear. Ruling out the
effect of the derotator position, strengthen the accuracy of our method. We
suspect that it could be to a miss-alignment in the red mosaic chip. When
dealing with the range of 487-567 nm, we reached the method’s limits ($\pm$ 18
mas) . Even by using UVES, one of the state-of-the-art spectrographs, it is
not easy to quantify the atmospheric dispersion below the accuracy of our
method. This only shows that measuring on-sky the atmospheric dispersion is a
challenging task.

## Acknowledgements

The first author is supported by a Fundação para a Ciência e Tecnologia (FCT)
fellowship (PD/BD/135225/2017), under the FCT PD Program PhD::SPACE
(PD/00040/2012). This work was supported by FCT/MCTES through national funds
and by FEDER - Fundo Europeu de Desenvolvimento Regional through COMPETE2020 -
Programa Operacional Competitividade e Internacionalização by these grants:
UID/FIS/04434/2019; PTDC/FIS-AST/32113/2017 & POCI-01-0145-FEDER-032113.

## Data availability

The data underlying this article are based on observations made with ESO
Telescopes at the Paranal Observatory under program ID 4103.L-0942(A) (with
the UVES spectrograph at the ESO VLT UT2 telescope) and available in ESO
Science Archive Facility at http://archive.eso.org/, and can be accessed with
the program ID. The flat frames used in the derotator test are available upon
request from the author.

## References

  * Barrell & Sears (1939) Barrell H., Sears J., 1939, Philosophical Transactions of the Royal Society of London Series A, 238, 1
  * Bechter et al. (2018) Bechter A. J., Bechter E. B., Crepp J. R., King D., Crass J., 2018, in Proc. SPIE. p. 107026T (arXiv:1812.02704), doi:10.1117/12.2313658
  * Blackman et al. (2019) Blackman R. T., Ong J. M. J., Fischer D. A., 2019, AJ, 158, 40
  * Bönsch & Potulski (1998) Bönsch G., Potulski E., 1998, Metrologia, 35, 133
  * Cabral et al. (2012) Cabral A., et al., 2012, in Ground-based and Airborne Telescopes IV. p. 84444F, doi:10.1117/12.926093
  * Dekker et al. (2000) Dekker H., D’Odorico S., Kaufer A., Delabre B., Kotzlowski H., 2000, Design, construction, and performance of UVES, the echelle spectrograph for the UT2 Kueyen Telescope at the ESO Paranal Observatory. pp 534–545, doi:10.1117/12.395512
  * Edlen (1953) Edlen B., 1953, Journal of the Optical Society of America (1917-1983), 43, 339
  * Edlén (1966) Edlén B., 1966, Metrologia, 2, 71
  * Filippenko (1982) Filippenko A. V., 1982, PASP, 94, 715
  * Hohenkerk & Sinclair (1985) Hohenkerk C. Y., Sinclair A. T., 1985, NAO Tech. Note No. 63
  * Skemer et al. (2009) Skemer A. J., et al., 2009, PASP, 121, 897
  * Spanò (2014) Spanò P., 2014, in Advances in Optical and Mechanical Technologies for Telescopes and Instrumentation. p. 915157, doi:10.1117/12.2057072
  * Wehbe et al. (2020a) Wehbe B., Cabral A., Martins J. H. C., Figueira P., Santos N. C., Ávila G., 2020a, MNRAS, 491, 3515
  * Wehbe et al. (2020b) Wehbe B., Cabral A., Ávila G., 2020b, MNRAS, 499, 183