11institutetext:  Nonlinear Physics Centre and Centre for Ultra-high bandwidth
Devices for Optical Systems (CUDOS),
Australian National University, Canberra ACT 0200, Australia
Laboratory of Photonic Information Technology, School for Information and
Optoelectronic Science and Engineering,
South China Normal University, Guangzhou 510006, P.R. China

Wave propagation, transmission and absorption Optical bistability,
multistability, and switching, including local field effects Nonlinear
waveguides

# Nonlinear Mach-Zehnder-Fano interferometer

Yi Xu 1122 Andrey E. Miroshnichenko E-mail: 11<EMAIL_ADDRESS>1122

###### Abstract

We demonstrate that the interaction of loop and nonlinear Fano resonances
results in a formation of hybrid resonant states in Mach-Zehnder type
interferometers, providing with opportunities for an advanced phase
manipulation. The nonlinear response of such structures can be greatly
enhanced, leading to a low threshold 100$\%$ switching operation. We further
propose one of the possible realizations based on nonlinear photonic crystal
circuits, suitable for optimal all-optical switching.

###### pacs:

42.25.Bs

###### pacs:

42.65.Pc

###### pacs:

42.65.Wi

##  1 Introduction

Mach-Zehnder interferometer (MZI) is a key component in many branches of
physics because of its ability to manipulate a coherent signal [1]. By
coupling a resonator to the MZI can further increase the phase sensitivity of
the coherent manipulation [2, 3]. The enhanced all-optical switching [2] and
the bistability [4] have been demonstrated in a coupled ring-resonator Mach-
Zehnder interferometer, which provides the possibility for the effective and
coherent control by using a nonlinear resonator.

Figure 1:  (color online)(a) and (b) are the generic discrete models for the
system exhibiting Fano resonance. $A$ and $B$ are the incoming waves, $T$ and
$R$ are scattered amplitudes. $C_{i}$ and $D_{i}$ are the forward and backward
scattered waves. $V$ is the coupling strength between the nonlinear Fano
defect(located by $M$) and the site $A_{M}$ in the linear chain. $L$ and $N$
specify the number of sites. Without loss of generality, we choose L=1 which
provides the phase reference while $N$ represents the length of MZFI arm
including two junction sites.

Recently, we have introduced the concept of Mach-Zehnder-Fano interferometer
(MZFI) [5] providing with unique physical property that can not be found in a
macroscopic resonator enhanced MZI [2, 3, 4, 6]. The MZFI allows us to
manipulate the interaction of different types of resonances which leads to the
formation of a novel hybrid Fano-like resonant states [7]. Furthermore, the
counterpart of the ring-resonator-coupled Mach-Zehnder interferometer in the
microscopic scale, i.e. MZFI based on a photonics crystal (PhC) platform seems
to be more promising for future application owing to small volume compared
with the macroscopic resonators. Recent advantages in PhCs fabrication
technology [8], allow us to achieve ultra high-Q cavities facilitating low
threshold nonlinear bistability [9, 10, 11, 12, 13]. Indeed, the manifested
optical bistable state is the nonlinear Fano resonance [15]. As a result, the
systems supporting Fano resonances [14, 15, 16], associated with an asymmetric
scattering profile, attract a significant attention recently. It’s their
unique property that allows us to achieve an optimal high extinction ratio,
large modulation depth and lowest threshold nonlinear switching [17, 6].

The aim of this Letter is to introduce and demonstrate unique properties of
the nonlinear MZFI, which originate from the excitation of the nonlinear
hybrid Fano resonances. Such resonant states appear because of the interaction
between MZI loop’s resonance and nonlinear Fano resonances of the side-coupled
defect. It provides with an enhanced nonlinear response and optimal conditions
for low threshold dynamic bistability. As a particular realisation of our
model we provide with a PhC circuit example, which supports our results.

Figure 2:  (color online)(a) and (b) are the transmission of the linear MZFI
and the intensity of the defect, respectively. Here, $L=1$ and $N=7$. The
circles in (b) represent the normalized intensity in the nonlinear Fano defect
when the input power $I=0.01$ and $\chi^{(3)}=1$. (c) nonlinear response of
different resonances marked $\omega_{a-d}$ at Fig. 2 (a). Solid pink line
stands for the hybrid Fano resonance marked $\omega_{a}$ in Fig. 2 (a), dotted
red line presents the resonance of the MZI’s loop marked $\omega_{d}$ while
dashed green line and dashed-dotted blue line represents the resonance similar
with the one in resonator enhanced MZI [2, 3, 4, 6](marked $\omega_{c}$) and
conventional one(marked $\omega_{b}$), respectively. The shaded pink region
represents the area of the dynamic modulation instability initialled by a
pulse of $t_{0}=1.3\times 10^{4}$, $W=5\times 10^{3}$ and $I_{0}=0.25$. (d)
Double bistability at the frequency $\omega=0.43$ with detuning
$\Delta\omega=0.11$ respecting to $\omega_{a}$. Here, $t_{0}=1.3\times
10^{4}$, $W=5\times 10^{3}$ and $I_{0}=0.025$. The arrows indicate the upward
and downward nonlinear Fano bistable operation. The responses of
$\omega_{b-d}$ for the same detuning value are also presented.

##  2 Discrete nonlinear MZFI model

We start our analysis from a generic discrete model [15], which can describe
the dynamics of MZFI like structures. By using the modified Fano-Anderson
model, dynamical equations describing the nonlinear MZFI system shown in Fig.
1 (b) can be written as follows:

$\begin{array}[]{rcl}i\dot{\psi}_{n}=\sum_{k}\psi_{k}+\delta_{n,M}V\varphi_{d}\;,\\\
i\dot{\varphi}_{d}=E_{d}\varphi_{d}+\chi^{(3)}|\varphi_{d}|^{2}\varphi_{d}+V\psi_{M}\;,\end{array}$
(1)

where $\psi_{n}$ represents the linear chain with complex field amplitude,
$\varphi_{d}$ stands for the Fano defect, $M$ gives the location of the Fano
defect in the arm, $k$ is the total number of the nearest neighbour sites in
the chain($k=3$ in the Y-junction while $k=2$ at the others) and we neglect
the nonlocal interactions in this paper, $\chi^{(3)}$ is the cubic nonlinear
parameter, $V$ is the coupling strength between the chain and the Fano defect,
$E_{d}$ is the eigenfrequency of the Fano defect.

As it was demonstrated in Ref. [18], the side-coupled Fano defect acts as an
effective scattering potential, whose strength depends on the input frequency
$\omega$. In the presence of nonlinearity the effective scattering potential
of the Fano defect becomes input intensity dependent, which brings extra
flexibility to tune the Fano resonances of the MZFI.

To find static solutions of system (1) we employed S-matrix approach. The
S-matrices of the Y-junction can be obtained by applying scattering boundary
condition at three branches respectively [19]. A set of nonlinear equations
which describes the stationary nonlinear response of the system is as follow:

$\displaystyle\left(\begin{array}[]{c}R\\\ C_{1}\\\ C_{2}\end{array}\right)$
$\displaystyle={\mathbf{S}}^{Y}_{L}\left(\begin{array}[]{c}A\\\ D_{1}\\\
D_{2}\end{array}\right)\;,\;\left(\begin{array}[]{c}T\\\ D_{1}\\\
D_{3}\end{array}\right)={\mathbf{S}}^{Y}_{-N}\left(\begin{array}[]{c}B\\\
C_{1}\\\ C_{3}\end{array}\right),$ (14)
$\displaystyle\left(\begin{array}[]{c}D_{2}\\\ C_{3}\end{array}\right)$
$\displaystyle={\mathbf{S}}_{M}\left(\begin{array}[]{c}C_{2}\\\
D_{3}\end{array}\right),$ $\displaystyle X_{d}$
$\displaystyle=\frac{A_{M}V}{(\omega-E_{d}-\lambda|X_{d}|^{2})}$

where

${\mathbf{S}}^{Y}_{K}=\left(\begin{array}[]{ccc}e^{-2iKq}r_{Y}&t_{Y}&t_{Y}\\\
t_{Y}&e^{2iKq}r_{Y}&e^{2iKq}t_{Y}\\\
t_{Y}&e^{2iKq}t_{Y}&e^{2iKq}r_{Y}\end{array}\right),\\\
{\mathbf{S}}_{M}=\frac{1}{\varepsilon-2i\sin
q}\left(\begin{array}[]{cl}-\varepsilon e^{-2iMq}&-2i\sin q\\\ -2i\sin
q&-\varepsilon e^{2iMq}\end{array}\right),\\\
\varepsilon=\frac{V^{2}}{(\omega-E_{d}-\lambda|X_{d}|^{2})}$

As can be seen from the complex defect field $X_{d}$, the nonlinearity has an
effect on the phase and amplitude of the Fano defect and the corresponding
strengths in turn depend on the power of the Fano defect. At the same time,
the nonliear response of the Fano defect would give a feedback to the
scattering waves in the arms of MZFI (read as the nonlinear scattering
potential $\varepsilon$) of which forms a complex nonlinear resonance system.
Then, such nonlinear MZFI, whose nonlinear response can not be mapped from the
standard nonlinear Fano resonant system [18], is sophisticated even in the
linear region [7].

We study two cases, where the Fano defect is placed symmetrically and
asymmetrically. The detail parameters can be found in the caption of Fig. 2.
These two cases allow us to excite different sets of MZI loop’s modes, which
have nonzero overlap with the coupled site [7]. The scattering of the
conventional Fano resonance geometry [see Fig. 1(a)] is shown by dotted line
in Fig. 2(a) as a reference. In the symmetric case, the transmission in the
centre of the propagation band resembles a step-function. In such a case, the
eigenfrequency of the Fano defect is in the vicinity of MZI loop’s resonance,
where a hybrid Fano resonance is formed. At the same time, such resonance has
the highest intensity at the defect site among other resonances because of the
interference of two Fano-like resonances [7]. Intuitively, the high contrast
step function like transmission together with high resonant intensity of the
nonlinear defect are positive factors for further enhancment of the nonlinear
behaviour [20]. We, thus, investigate the nonlinear response of the system (1)
at these specific resonances.

The stationary nonlinear switching of four resonances marked
$\omega_{a-d}$(corresponding to transmission minima) in Fig. 2 (a) are shown
in Fig. 2 (c). It should be pointed out that the maximum of the Fano defect’s
intensity is in between the transmission dip and tip because of the sharp
asymmetric line shapes [15]. We emphasize that the resonance $\omega_{d}$ is
excited due to a symmetry breaking by side-coupled Fano defect [7] and it is
not the eigenfrequency of the defect. At the same time, the nonlinear response
at $\omega_{a}$ is greatly reduced by the step-function like linear
transmission compared to resonance $\omega_{b}$ and $\omega_{c}$. If we define
a figure of merit as $FoM=\Delta T/P_{th}$, where $\Delta
T=\max\limits_{P_{in}}[T(P_{in})-T(0)]$ refers to the maximum contrast of the
transmission and $P_{th}$ represents the necessary input power to pull the
transmission of the system up to $90\%$ of the $\Delta T$. The FoM of
resonance $\omega_{a}$, which describes both the enhanced transmission
contrast in the linear case and the reduction of the switching power, can be
enhanced more than 60 times compared to a given defect supporting conventional
Fano resonance $\omega_{b}$ [18]. The nonlinear responses at $\omega_{b}$ and
$\omega_{c}$ are similar because they are originated from the eigenfrequency
of the Fano defect.

To study the dynamical switching we employed high accuracy Crank-Nicolson
method [22] with suitable absorption boundary condition [23]. The dynamic
nonlinear response of the system at specified resonances is obtained by a long
pulse with $I=I_{0}\exp(-(t-t_{0})^{2}/W^{2})\sin(\omega t)$, where $W$ is the
pulse width, $\omega$ is set at the same frequencies as $\omega_{a-d}$.
Compared to the stationary one, the dynamic solution predicts similar
nonlinear response except for the case of $\omega_{a}$. One can see that
transmission cannot be well defined at frequency $\omega_{a}$ above a certain
threshold of the input power, indicated by the shaded region. This is caused
by the modulational instability of the Fano resonances [21], which also
indicates an enhanced nonlinear response. It is similar to the modulational
instability of waves scattered by a nonlinear centre [24]. Usually, it happens
in the finite interval of frequencies . The stability analysis of our system
with the nonlinear hybridization between Fano resonances gives similar results
to the conventional nonlinear Fano resonance with the same instability regions
[21]. Note here, that dynamical instability occurs independently from
bistability conditions.

Figure 3:  (color online) (a) Typical perfect bistability at the frequency
$\omega=-0.2$ with respect to $\omega_{a^{\prime}}$. Here, $t_{0}=1.3\times
10^{4}$, $W=5\times 10^{3}$ and $I_{0}=0.025$. (b) Intensity of the Fano
defect with detuning $\Delta\omega=0.11$ with respect to
$\omega_{a,a^{\prime}}$.

Figure 2 (d) demonstrates a butterfly like double bistability obtained by
dynamic pulse feeding, where the operation frequency is referred to
$\omega_{a}=0.43$. Such a hysteresis loop at $\omega_{a}$ shows interesting
properties that one loop is clockwise while the other is anti-clockwise. The
non-uniform shift of the normalized intensity in the nonlinear Fano defect
[see balls in Fig. 2 (b)] leads to the double bistable operation. As can be
seen from this figure, that such a detuning is not enough to initiate a
bistable response for other case of $\omega_{b}$, $\omega_{c}$, and
$\omega_{d}$. Thus, the nonlinear hybridization of various resonant states
offers a unique opportunity to manipulate the nonlinear Fano resonances.
According to the stationary normalized intensity of the defect [see Fig. 2
(b)], we can obtain further nonlinear enhancement by working with
$\omega_{a^{\prime}}$. Figure 3 (a) shows an enhanced bistability when the
operation frequency is set at $\omega=-0.21$. This kind of perfect bistable
state benefits from the sharp step function like linear transmission and the
further enhanced intensity of the nonlinear defect. The intensities of the
Fano defect [see Fig. 3(b)] which demonstrate the further enhancement of the
intensity in the Fano defect are in accord with the static result [see Fig.
2(b)]. The oscillation between two bistable state is the process that one
bistable state transfers to the other one and it is the properties of the
dynamical bistability [21].

Figure 4:  (color online)(a) Bistable operation obtain by pulse excitation in
nonlinear FDTD experiment; (b) electric field distribution ($|E|^{2}$) at the
upward dynamic switching point.

##  3 Photonic Crystal based nonlinear MZF interferometer

Using the modified Fano-Anderson model, we have proved that the nonlinear MZFI
can be considered as suitable candidate to realize the enhanced nonlinear
response and manipulate the dynamical bistability. We thus suggest a PhCs
platform as one of the possibilities to realize our idea, while the results
above can be applied to other varieties of nonlinear discrete system. The PhC
structure is shown in Fig. 4 (b) and It is formed by dielectric rods embedded
in air(the radius of dielectric rods is $r=0.19a$ where $a$ is the lattice
constant and refractive index $n=3.14$) except for the nonlinear Fano defect
(polymer rod $n=1.59$ and the third-order nonlinearity susceptibility of the
polymer is $\chi^{(3)}=1.14\times 10^{-12}cm^{2}/W$). We use the nonlinear
FDTD method to solve the Maxwell’s equations. The bistable response shown in
Fig. 4 (a) is obtained by the response of a Gaussian pulse with input
frequency $f=0.373$ $2\pi c/a$ and duration $30$ picoseconds. The profile is
similar to the theoretical result in Fig. 3. Figure 4 (b) shows the transient
electric field distribution ($|E|^{2}$) exactly at the downward switching
point (from on to off-state). Nearly perfect blocking of the input pulse
demonstrates the dynamical shutting down operation by a pulse and successful
suppression of the modulation instability. Both high intensities of the Fano
defect and the loop are the signature of nonlinear switching involving unique
hybrid Fano resonances in the MZFI [7]. The transmission contrast can be
further enhanced by careful engineering the interaction between the eigen-Fano
resonance and the MZI loop’s resonance which takes place exactly at the
propagation band center.

##  4 Conclusions

In conclusion, we have investigated the interaction of loop and side-coupled
Fano resonances in nonlinear MZFI type structures. By introducing nonlinear
Fano defects we are able to excite ”dark” states of pure MZI, which results in
the formation of very narrow hybrid resonant states. Both stationary and
dynamic studies convince the superiority of enhanced nonlinear response and
dynamical bistability in nonlinear MZFI. Our direct numerical simulations of
two-dimensional PhCs confirm the theoretical predictions of the dynamic
characteristics. We anticipate that such structures can be realized on
different platforms for ultra-high sensitive sensing operations. The idea of
nonlinear Fano resonances hybridization can also be generalized to the other
nanoscale structure. Particularly, in the new emerging field of plasmonics and
metamaterial where the electromagnetic wave is confined in subwavelength
scale, such nonlinear resonance interaction could be a basic physical idea for
designing the logic devices.

###### Acknowledgements.

The authors thank Prof. Yuri Kivshar and Dr. Anton Desyatnikov for useful
discussions. Y. Xu acknowledges the support from the China Scholarship Council
and the Nonlinear Research Centre at ANU for their hospitality. The work of A.
E. Miroshnichenko was supported by the Australian Research Council through
Future Fellowship program.

## References

  * [1] Ernst M. The Principles of Physical Optics Courier Dover Publications, Dover 2003
  * [2] Heebner J. E. Boyd R. W. Opt. Lett.241999847-849.
  * [3] Heebner J. E. et. al. Opt. Lett.292004769-771.
  * [4] Lu Y. et. al. Opt. Lett.3020053069-3071.
  * [5] Miroshnichenko A. E. Kivshar Y. S. Appl. Phys. Lett.952009121109.
  * [6] Mario L. Y., Darmawan S. Chin M. K. Opt. Express14200612770-12781.
  * [7] Xu Y. Miroshnichenko A. E. Phys. Rev. A842011033828.
  * [8] Akahane Y. et. al. Nature4252003944.
  * [9] Solja$\check{c}$i$\acute{c}$ M. et. al. Phys. Rev. E662002055601.
  * [10] Yanik M. F., Fan S. Soljacic M. Appl. Phys. Lett.1420032739-2741.
  * [11] Solja$\check{c}$i$\acute{c}$ M. et. al. Opt. Lett.282003637-639.
  * [12] Tanabe T. et. al. Opt. Lett.3020052575-2577.
  * [13] Yang X. et. al. Appl. Phys. Lett.912007051113.
  * [14] Fano U. Phys. Rev.12419611866.
  * [15] Miroshnichenko A. E., Flach S. Kivshar Y. S. Rev. Mod. Phys.8220102257.
  * [16] Luk’yanchuk B., Zheludev N. I., Maier S. A., Halas N. J., Nordlander P., Giessen H., Chong C. T. Nat. Mater.92010707-715.
  * [17] Fan S. Appl. Phys. Lett.802002908-910.
  * [18] Miroshnichenko A. E. et. al. Phys. Rev. E712005036626.
  * [19] Miroshnichenko A. E., Molina M. I. and Kivshar Y. S. Phys. Rev. E752007046602.
  * [20] Darmawan S., Landobasa Y. M. Chin M. K. Opt. Express152007437-448.
  * [21] Miroshnichenko A. E. et. al. Phys. Rev. A792009013809.
  * [22] Crank J. Nicolson P. Proc. Camb. Phil. Soc.43194750-67.
  * [23] Arnold A., Ehrhardt M. Sofronov I. Commun. Math. Sci.12003501.
  * [24] Malomed B. A. and Azbel M. Y. Phys. Rev. B47199310406.


#  $B$ Anomalies and Dark Matter in an $L_{\mu}-L_{\tau}$ Model with General
Kinetic Mixing

Douglas Tuckler<EMAIL_ADDRESS>Department of Physics, Carleton
University, Ottawa, ON K1S 5B6, Canada

###### Abstract

We revisit the $L_{\mu}-L_{\tau}$ extension of the Standard Model that can
simultaneously address anomalies in semileptonic $B$ meson decays and the
nature of dark matter (DM). In the region favored by the $B$ anomalies, this
scenario is excluded by a combination of low-energy flavor constraints and
stringent DM direct detection constraints if the kinetic mixing between
$Z^{\prime}$ gauge boson and the photon vanishes at high energy scales, since
this leads to a sizable coupling between DM and the SM in low-momentum
scattering processes. However, this is no longer the case if the kinetic
mixing vanishes at low energy scales instead. In this scenario, the low-
momentum DM scattering rate can be suppressed and constraints from direct
detection experiments can be substantially relaxed. As a result, we find a re-
opening of the region of parameter space where the $B$ anomalies and DM can be
simultaneously explained.

##  I Introduction

The $Z^{\prime}$ gauge boson of a broken $U(1)_{L_{\mu}-L_{\tau}}$ symmetry is
a simple extension of the Standard Model (SM) He:1991qd ; He:1990pn  that can
address the tension in the anomalous magnetic moment of the muon $(g-2)_{\mu}$
Baek:2001kca , can serve as mediator for dark matter (DM) interactions with
the SM Altmannshofer:2016jzy ; Arcadi:2018tly ; Kamada:2018zxi ;
Foldenauer:2018zrz ; Borah:2021jzu ; Holst:2021lzm , and can address anomalies
observed in the decays of $B$ mesons Altmannshofer:2014cfa ; Crivellin:2015mga
; Altmannshofer:2016jzy ; Crivellin:2016ejn ; Ko:2017yrd ; Arcadi:2018tly . In
fact, a virtue of the $U(1)_{L_{\mu}-L_{\tau}}$ extension of the SM is that in
can address more than one of these questions at the same time.

Previously, it was shown that the $Z^{\prime}$ can address both the $B$
anomalies and DM Altmannshofer:2016jzy . By introducing a flavor-violating
coupling of the $Z^{\prime}$ to left-handed $b$ and $s$ quarks the
$Z^{\prime}$ is able to modify the $b\to s\mu\mu$ transition as preferred by
global fits to rare $B$ decays Altmannshofer:2014cfa ; Altmannshofer:2015mqa ;
Altmannshofer:2021qrr , and the kinetic mixing between the $Z^{\prime}$ and
the SM photon allows for $L_{\mu}-L_{\tau}$ charged DM to interact with SM
particles. Strong constraints from DM direct detection experiments allowed for
only a narrow window where the DM annihilates nearly on resonance in the early
Universe Altmannshofer:2016jzy . As we will see, the latest DM direct
detection results completely exclude this scenario.

However, the interpretation of these direct detection constraints depends
strongly on the energy dependence of the kinetic mixing. Typically, it is
assumed that the kinetic mixing vanishes at energy scales much higher than the
$\tau$ lepton mass, which leads to a sizable coupling at low energies between
DM and nuclei of direct detection experiments. This assumption is by no means
general and other boundary conditions are equally viable.

Recently, it was shown that DM direct detection constraints can be
substantially weakened if the kinetic mixing is suppressed at energy scales
lower than the muon mass Hapitas:2021ilr . It was shown in Hapitas:2021ilr
that this scenario can open up new parameter space for a simultaneous
explanation of DM and the $(g-2)_{\mu}$ anomaly.

In this paper, we revisit the explanation of the $B$ anomalies and DM in the
$U(1)_{L_{\mu}-L_{\tau}}$ model by considering this alternative treatment of
the kinetic mixing between the $Z^{\prime}$ and the SM photon. We will show
that stringent constrains from DM direct detection experiments can be relaxed
if the kinetic mixing is suppressed at low momentum transfer. We find a
reopening of viable parameter space that addresses the $B$ anomalies and DM
simultaneously, consistent with all other existing constraints from neutrino
trident production, collider searches, and DM indirect detection.

This paper is organized as follows. In Sec. II we introduce the
$U(1)_{L_{\mu}-L_{\tau}}$ model and discuss the momentum dependence of the
kinetic mixing between $Z^{\prime}$ and the photon. In Sec. III we discuss the
connection between the model under consideration and the $B$ anomalies, as
well as flavor constraints that arise from the flavor-violating quark
coupling. Sec. IV focuses on the DM phenomenology, updating the model of
Altmannshofer:2016jzy  with the latest DM direct detection constraints and
showing the effect of the momentum dependence of the kinetic mixing on the
these constraints. In Sec. V we compute the constraints from shifts in the $Z$
boson mass induced by the kinetic mixing with the $Z^{\prime}$ at high scales.
We conclude in Sec. VI.

##  II The $L_{\mu}-L_{\tau}$ Model

Following Altmannshofer:2016jzy , we consider the gauged
$U(1)_{L_{\mu}-L_{\tau}}$ extension of the SM, whose symmetry breaking gives
rise to a $Z^{\prime}$ gauge boson with mass $m_{Z^{\prime}}$.111We remain
agnostic about the details of this symmetry breaking and assume that the
physics involved is sufficiently decoupled, having no effect on the
phenomenology discussed in this paper. In addition, we introduce a Dirac
fermion DM candidate $\chi$ with mass $m_{\chi}$ and charged under the
$U(1)_{L_{\mu}-L_{\tau}}$ symmetry. We also introduce a flavor-violating
interaction of the $Z^{\prime}$ to left-handed $b$ and $s$ quarks. The
Lagrangian describing this model is

$\begin{split}\mathcal{L}&{}_{L_{\mu}-L_{\tau}}\supset\
-\frac{1}{4}Z^{\prime}_{\alpha\beta}Z^{\prime\alpha\beta}+\frac{1}{2}m^{2}_{Z^{\prime}}Z^{\prime}_{\alpha}Z^{\prime\alpha}+\frac{\varepsilon_{0}}{2}Z^{\prime}_{\alpha\beta}F^{\alpha\beta}\\\
&+q_{\ell}g^{\prime}\big{(}\bar{\mu}\gamma^{\alpha}\mu-\bar{\tau}\gamma^{\alpha}\tau+\bar{\nu}_{\mu}\gamma^{\alpha}P_{L}\nu_{\mu}-\bar{\nu}_{\tau}\gamma^{\alpha}P_{L}\nu_{\tau}\big{)}Z^{\prime}_{\alpha}\
\\\ &+i\bar{\chi}\gamma^{\alpha}\partial_{\alpha}\chi-
m_{\chi}\bar{\chi}\chi+q_{\chi}g^{\prime}\bar{\chi}\gamma^{\alpha}\chi
Z^{\prime}_{\alpha}\ \\\
&+g_{bs}\bar{s}\gamma^{\alpha}P_{L}bZ^{\prime}_{\alpha}\ ,\end{split}$ (1)

where $g^{\prime}$ is the $U(1)_{L_{\mu}-L_{\tau}}$ coupling constant,
$q_{\ell}$ and $q_{\chi}$ are the charge of the leptons and DM, respectively,
under $U(1)_{L_{\mu}-L_{\tau}}$, and $g_{bs}$ is the flavor-violating coupling
of $Z^{\prime}-b-s$ interaction. Note that $q_{\ell,\chi}$ are free parameters
but we will simply assume $q_{\chi}=q_{\ell}=1$.

In Eq. (1) we have introduced a bare kinetic mixing parameter
$\varepsilon_{0}$ between $Z^{\prime}$ and the photon field strengths. This is
a free parameter of the model and depends on unknown ultraviolet (UV) physics.
However, we will see that it can be fixed by considering the asymptotic
boundary conditions of the total kinetic mixing in, for example, DM-nucleon
scattering.

###  II.1 Momentum Dependence in Kinetic Mixing

Figure 1: Momentum dependent kinetic mixing $\varepsilon_{\text{tot}}(Q)$ as a
function of momentum transfer $Q$. The solid (dashed) curve corresponds to the
boundary condition $\varepsilon_{\text{IR}}=0$ of Eq. (4)
($\varepsilon_{\text{UV}}=0$ of Eq. (3)). The vertical, dashed purple lines
indicate where the momentum transfer is equal to $m_{\mu}$ and $m_{\tau}$.

In this model, the $Z^{\prime}$ gauge boson does not couple directly to
electrons and quarks of the SM, and the scattering of DM (and neutrinos) with
nucleons of a detector is induced by the kinetic mixing between the
$Z^{\prime}$ and the photon. The total kinetic mixing has a tree-level
contribution $\varepsilon_{0}$ and a one-loop contribution from virtual $\mu$
and $\tau$ leptons. The total kinetic mixing with explicit momentum transfer
dependence is given by

$\varepsilon_{\rm
tot}(Q)=\varepsilon_{0}-\frac{eg^{\prime}}{2\pi^{2}}\int_{0}^{1}dx\,x(1-x)\log\left[\frac{m_{\tau}^{2}+x(1-x)Q^{2}}{m_{\mu}^{2}+x(1-x)Q^{2}}\right]\
,$ (2)

where $Q^{2}=\sqrt{(-q^{2})}>0$, with $q^{\mu}$ being the momentum transfer in
the scattering process. The bare kinetic mixing parameter $\varepsilon_{0}$
can be fixed by considering boundary conditions of the total kinetic mixing.
Here we consider two boundary conditions:

$\displaystyle\varepsilon_{\rm tot}(Q)$ $\displaystyle\to
0~{}\text{as}~{}Q\to\infty~{}~{}\Rightarrow~{}~{}\varepsilon_{0}=0$ (3)
$\displaystyle\varepsilon_{\rm tot}(Q)$ $\displaystyle\to 0~{}\text{as}~{}Q\to
0~{}~{}~{}\Rightarrow~{}~{}\varepsilon_{0}=\frac{eg^{\prime}}{12\pi^{2}}\log\frac{m^{2}_{\tau}}{m^{2}_{\mu}}\
,$ (4)

which fixes the bare kinetic mixing parameter $\varepsilon_{0}$ to the values
after the double right arrows in Eqs. (3) and (4). For brevity, we will denote
the boundary condition in Eq. (3) by $\varepsilon_{\text{UV}}=0$ (i.e
$\varepsilon_{\rm tot}$ vanishes at high or UV scales), and the boundary
condition in Eq. (4) by $\varepsilon_{\text{IR}}=0$ (i.e $\varepsilon_{\rm
tot}$ vanishes at low or infrared (IR) scales).

The boundary condition in Eq. (3) is one that is typically assumed in the
literature, but it is by no means a generic choice. Recently, it was shown
that varying the boundary condition of the total kinetic mixing has non-
trivial effects the phenomenology of neutrino and DM scattering in the
$L_{\mu}-L_{\tau}$ model Hapitas:2021ilr . In particular, the boundary
condition in Eq. (4) results in a suppression of the total kinetic mixing at
low momentum transfer, as depicted in Fig. (1).

An important consequence of this behavior is that constraints from processes
that occur at low momentum transfer, such as those from coherent elastic
neutrino-nucleus scattering (CE$\nu$NS) and DM direct detection experiments,
can be substantially weaker and regions of parameter space that would
otherwise be severely constrained become viable. We will explore this for DM
direct detection in Sec. IV.

##  III $B$ Anomalies and Constraints

In this section we discuss the how the $Z^{\prime}$ of the
$U(1)_{L_{\mu}-L_{\tau}}$ model can address the $B$ anomalies. We also briefly
discuss constraints on the parameter space from meson mixing, neutrino trident
production, and $Z\to 4\mu$ searches.

###  III.1 $B$ Anomalies

Over the past decade, $B$ factories such as Belle, LHCb, ATLAS, and CMS have
observed anomalies in rare semileptonic decays of $B$ mesons. Discrepencies
between experimental results and SM predictions have been seen in differential
branching ratios, angular observables, and lepton flavor universality (LFU)
ratios of various $b\to s$ transitions such as $B^{0}\to K^{\ast
0}\ell^{+}\ell^{-}$, $B^{+}\to K^{+}\ell^{+}\ell^{-}$, and
$B_{s}\to\phi\ell^{+}\ell^{-}$ Belle:2009zue ; CDF:2011buy ; BaBar:2012mrf ;
LHCb:2014cxe ; LHCb:2015tgy ; LHCb:2015wdu ; CMS:2015bcy ; LHCb:2016ykl ;
Belle:2016fev ; LHCb:2017avl ; CMS:2017rzx ; ATLAS:2018gqc ; LHCb:2019hip ;
Belle:2019oag ; BELLE:2019xld ; LHCb:2020lmf ; LHCb:2020gog ; LHCb:2021trn ;
LHCb:2021zwz ; LHCb:2021lvy .

Global fits to all rare $B$ decay data find a strong preference for new
physics in the form a four-fermion contact interaction described by the
effective Hamiltonian

$\mathcal{H}^{\text{NP}}_{\text{eff}}=-\frac{4G_{F}}{\sqrt{2}}V_{tb}V^{\ast}_{ts}\frac{e^{2}}{16\pi^{2}}C^{\text{NP}}_{9}(\bar{s}\gamma_{\alpha}P_{L}b)(\bar{\mu}\gamma^{\alpha}\mu)$
(5)

with a best fit value for the Wilson coefficient of $C_{9}^{\text{NP}}=-0.73$
Altmannshofer:2021qrr . This kind of interaction can be accommodated in the
$U(1)_{L_{\mu}-L_{\tau}}$ model since the $Z^{\prime}$ has the required vector
coupling to muons, but does not couple to electrons to leading order.
Integrating out the $Z^{\prime}$ in Eq. (1) leads to the contact interaction
in Eq. (5) and the Wilson coefficient is determined by $g^{\prime}$,
$m_{Z^{\prime}}$, and $g_{bs}$ as

$C^{\text{NP}}_{9}=-\frac{g^{\prime}g_{bs}}{2m^{2}_{Z^{\prime}}}\frac{\sqrt{2}}{4G_{F}}\frac{16\pi^{2}}{e^{2}}\frac{1}{V_{tb}V^{\ast}_{ts}}\
.$ (6)

We will fix the value of $g_{bs}$ such that this always holds i.e. we have

$g_{bs}=-\frac{4G_{F}}{\sqrt{2}}V_{tb}V^{\ast}_{ts}\frac{e^{2}}{16\pi^{2}}C^{\text{NP}}_{9}\frac{2m^{2}_{Z^{\prime}}}{g^{\prime}}\simeq\frac{10^{-9}}{\text{GeV}^{2}}\times\frac{m^{2}_{Z^{\prime}}}{g^{\prime}}\
,$ (7)

so that the phenomenology of this model is entirely determined by $g^{\prime}$
and $m_{Z^{\prime}}$. To successfully address the $B$ anomalies
$m_{Z^{\prime}}$ cannot be too small. If $m_{Z^{\prime}}$ is smaller than the
$B$ meson mass, it could be produced on-shell in $B$ decays and would appear
as a resonance in the dimuon invariant mass distribution. The non-observation
of such resonance by LHCb LHCb:2015nkv  restricts the $Z^{\prime}$ mass to be
above the $B$ meson mass, and we will conservatively require
$m_{Z^{\prime}}\gtrsim 10$ GeV to explain the $B$ anomalies. In Fig. (4) this
bound is denoted by the vertical dashed, black line at $m_{Z^{\prime}}=10$
GeV. Note, the choice of $m_{Z^{\prime}}\gtrsim 10$ GeV also allows us to use
the EFT framework described in Eq. (5).

###  III.2 Flavor Constraints

In addition to the $b\to s\mu\mu$ transition discussed above, the flavor
violating coupling in the last line of Eq. (1) will contribute at tree-level
to low energy flavor observables. In particular, it will lead to a tree-level
contribution to the $B_{s}$ mass difference $\Delta M_{s}$. The mass
difference is parameterized by the mixing amplitude $M_{12}$ as

$\Delta M_{s}=\Delta
M^{\text{SM}}_{s}\Bigg{|}1+\frac{M^{Z^{\prime}}_{12}}{M^{\text{SM}}_{12}}\Bigg{|}\
,$ (8)

where $\Delta M^{\text{SM}}_{s}$ is the SM prediction and
$M^{Z^{\prime}}_{12}$ ($M_{12}^{\text{SM}}$) is the $Z^{\prime}$ (SM)
contribution to the mixing amplitude. The modification due the $Z^{\prime}$
contribution to $M_{12}$ is given by

$\frac{M_{12}^{Z^{\prime}}}{M_{12}^{\text{SM}}}=\frac{e^{4}}{2\pi^{2}}\frac{|C_{9}^{\text{NP}}|^{2}}{M_{W}^{2}S_{0}(x_{t})}\frac{m^{2}_{Z^{\prime}}}{(q_{\ell}g^{\prime})^{2}}$
(9)

where $S_{0}(x_{t})$ us an Inami-Lim function Inami:1980fz  whose value is
$\simeq 2.3$, with $x_{t}=m_{t}/m_{W}$ being the ratio of the top quark mass
to the $W$ boson mass. Bounds on generic new physics contributions to $\Delta
M_{s}$ allow for $\Big{|}M_{12}^{Z^{\prime}}/M_{12}^{\text{SM}}\Big{|}\lesssim
0.12$ Charles:2020dfl , and will produce an upper bound on $m_{Z^{\prime}}$ as
a function of $g^{\prime}$ once $C_{9}^{\text{NP}}$ is fixed. This bound is
shown in Fig. (4) by the gray shaded region in the bottom-right corner of the
plots labeled “$B_{s}$ Meson Mixing”.

###  III.3 Collider Constraints

A powerful probe of the $Z^{\prime}$ parameter space is neutrino trident
production – the production of muon pairs when a muon neutrino scatters with
the Coulomb field of a target nucleus. This most precise measurement of this
process is from the CCFR CCFR:1991lpl  collaboration which finds good
agreement with the SM Altmannshofer:2014cfa

$\sigma_{\text{CCFR}}/\sigma_{\text{SM}}=0.82\pm 0.28\,.$ (10)

This result can set strong constraints on the parameter space of the
$U(1)_{L_{\mu}-L_{\tau}}$ model. Following the procedure in
Altmannshofer:2014pba , the region of parameter space excluded by neutrino
trident production is depicted in Fig. (4) by the gray shaded region above the
line labeled CCFR.

Additional constraints on the parameter space of this model can be obtained
from measurements of $Z\to 4\mu$ decays when the $Z^{\prime}$ is resonantly
produced. A search for this process was performed by the CMS collaboration
CMS:2018yxg  and sets strong constraints on the $g^{\prime}$ coupling for 5
GeV $\gtrsim m_{Z^{\prime}}\gtrsim$ 70 GeV. These constraints are shown in
Fig. (4) by the gray curve labelled CMS. We see that CMS sets stronger
constraints than CCFR in this region of parameter space.

Taking into account existing constraints from CCFR Altmannshofer:2014pba , CMS
CMS:2018yxg , and $B_{s}$ meson mixing, and also requiring that
$m_{Z^{\prime}}>$ 10 GeV, the region favored by the $B$ anomalies is the
central white region in Fig. (4) enclosed by the gray shaded regions above
$m_{Z^{\prime}}\simeq 50$ GeV.222Note, that there are additional constraints
from $e^{+}e^{-}\to 4\mu$ searches from BaBar TheBABAR:2016rlg  but these
apply outside of region favored by the $B$ anomalies and we omit them from
Fig. (4).

##  IV Dark Matter Phenomenology

In this section we discuss the relic density, direct detection, and indirect
detection of DM in the model, motivated by the region of parameter where the
$Z^{\prime}$ can explain the $B$ anomalies. In particular, we pay special
attention to the role of the momentum dependence of the kinetic mixing in DM-
nucleon scattering and the effect this has on DM direct detection constraints.

###  IV.1 The thermal dark matter target

Figure 2: Regions of the $g^{\prime}$ vs $m_{Z^{\prime}}$ parameter space
where the DM obtains the observed relic abundance. The red, coral, blue, and
green curves are for $m_{\chi}$ = 10, 25, 50, and 100 GeV, respectively.

We consider a DM candidate $\chi$ whose relic abundance is obtained via freeze
out. The relic abundance of DM today is determined by DM annihilation into
$Z^{\prime}$ bosons, or annihilation into SM leptons via off-shell
$Z^{\prime}$ exchange. The corresponding annihilation cross sections are

$\displaystyle(\sigma v)_{\chi\bar{\chi}\to Z^{\prime}Z^{\prime}}$
$\displaystyle=\frac{g^{\prime 4}}{16\pi
m_{\chi}^{2}}\left(1-\frac{m_{Z^{\prime}}^{2}}{m_{\chi}^{2}}\right)^{3/2}\left(1-\frac{m_{Z^{\prime}}^{2}}{2m_{\chi}^{2}}\right)^{-2}\
,$ (11a) $\displaystyle(\sigma
v)_{\chi\bar{\chi}\to\nu_{\ell}\bar{\nu}_{\ell}}$
$\displaystyle=\frac{g^{\prime
4}m_{\chi}^{2}}{2\pi[(4m_{\chi}^{2}-m_{Z^{\prime}}^{2})^{2}+m_{Z^{\prime}}^{2}\Gamma_{Z^{\prime}}^{2}]}\
,$ (11b) $\displaystyle(\sigma v)_{\chi\bar{\chi}\to\ell^{+}\ell^{-}}$
$\displaystyle=\frac{g^{\prime
4}m_{\chi}^{2}}{\pi[(4m_{\chi}^{2}-m_{Z^{\prime}}^{2})^{2}+m_{Z^{\prime}}^{2}\Gamma_{Z^{\prime}}^{2}]}$
$\displaystyle~{}~{}\times\left(1+\frac{m_{\ell_{a}}^{2}}{2m_{\chi}^{2}}\right)\sqrt{1-\frac{m_{\ell_{a}}^{2}}{m_{\chi}^{2}}}\
,$ (11c)

where the flavor index $\ell=\mu$ or $\tau$. The $Z^{\prime}$ decay width is
Kelly:2020pcy

$\Gamma_{Z^{\prime}}=\frac{g^{\prime
2}m_{Z^{\prime}}}{12\pi}\left[1+\sum_{\alpha=\mu}^{\tau}(1+2r_{\alpha})(1-4r_{\alpha})^{1/2}\Theta(1-4r_{\alpha})\right]\
,$ (12)

where $r_{\alpha}=m_{\alpha}^{2}/m_{Z^{\prime}}^{2}$ and $\Theta$ is the
Heaviside theta function. Given the annihilation cross sections, the DM relic
abundance can be found by numerically solving the Boltzmann equation.

In Fig. 2 we show the region of parameter space where DM obtains the observed
relic abundance of $\Omega h^{2}\simeq 0.12$ Planck:2018vyg . The red, coral,
blue and green curves are for $m_{\chi}$ = 10, 25, 50, and 100 GeV,
respectively. The behavior of these curves can be understood by analyzing the
annihilation cross sections in different regimes. For $m_{\chi}\gg
m_{Z^{\prime}}$ the annihilation cross section is proportional to $g^{\prime
4}/m^{2}_{\chi}$ and $g^{\prime}$ is a constant for a fixed DM mass, as shown
in the left part of the curves of Fig. (2). Then, we see the effect of
annihilation near a resonance at $m_{\chi}\simeq 2m_{Z^{\prime}}$, where the
cross section is enhanced and a much smaller value of $g^{\prime}$ can produce
the correct relic abundance. Finally, when $m_{\chi}\ll m_{Z^{\prime}}$ the
annihilation cross section scales as $(g^{\prime}/m_{Z^{\prime}})^{4}$ and
larger couplings are needed for larger $Z^{\prime}$ masses.

###  IV.2 Direct Detection

Since DM does not couple directly to quarks, its scattering with nuclei of DM
direct detection experiments occurs through the kinetic mixing between
$Z^{\prime}$ and the SM photon, and depends on the choice of boundary
condition for the total kinetic mixing. Here we consider DM direct detection
constrains for the two different choices $\varepsilon_{\text{UV}}=0$ and
$\varepsilon_{\text{IR}}=0$ discussed in Sec. II.

The DM-proton scattering cross section is given by

$\sigma_{\chi p}\approx\frac{g^{\prime 2}\varepsilon_{\rm
tot}^{2}e^{2}\mu_{\chi p}^{2}}{\pi m_{Z^{\prime}}^{4}}\ ,$ (13)

where $\mu_{\chi p}=m_{\chi}m_{p}/(m_{\chi}+m_{p})$ is the reduced mass of DM-
proton system, and we’ve neglected the momentum transfer in the $Z^{\prime}$
propagator. To account for the momentum transfer dependence, we consider the
nucleus level DM scattering rate derived using the standard halo model

$\displaystyle R$ $\displaystyle=\frac{N_{A}M_{A}\rho_{\odot}g^{\prime
2}Z^{2}e^{2}}{2\pi m_{\chi}}$ $\displaystyle\times\int_{E_{\rm
th}}^{E_{R}^{\rm max}}dE_{R}\frac{\varepsilon^{2}_{\rm tot}(Q)F_{\rm
Helm}^{2}(Q)}{(Q^{2}+m_{Z^{\prime}}^{2})^{2}}\int_{v_{\rm min}(E_{R})}^{v_{\rm
esc}}d^{3}v\frac{f(\vec{v})}{v}\ ,$ (14)

where $N_{A}$ is the number of target nucleus, $M_{A}$ and $Ze$ are the mass
and electric charge of the nucleus, $F_{\rm Helm}$ is the Helm form factor
Helm:1956zz ; Jungman:1995df ; Duda:2006uk ; Hoferichter:2020osn , and
$Q\simeq\sqrt{2M_{A}E_{R}}$ for non-relativistic scattering. The local DM mass
density is $\rho_{\odot}=0.3\,{\rm GeV/cm^{3}}$ and the velocity distribution
is Maxwellian
$f(\vec{v})=C\exp(-|\vec{v}+\vec{v}_{\odot}|^{2}/v_{0}^{2})\,\Theta(v_{\rm
esc}-|\vec{v}+\vec{v}_{\odot}|)$, where $\vec{v}$ is the DM velocity in the
rest frame of the solar system, $v=|\vec{v}|$, $v_{\odot}=220\,{\rm km/s}$,
$v_{0}=235\,{\rm km/s}$, $v_{\rm esc}=550\,{\rm km/s}$ Lin:2019uvt , and $C$
is a normalization factor such that $\int d^{3}vf(\vec{v})=1$.

The limits of integration over velocity and recoil energy in Eq. (IV.2) are

$v_{\rm min}(E_{R})=\sqrt{\frac{M_{A}E_{R}}{2\mu_{\chi A}^{2}}}\ ,\quad\quad
E_{R}^{\rm max}=\frac{2\mu_{\chi A}^{2}(v_{\rm esc}+v_{\odot})^{2}}{M_{A}}\ ,$
(15)

where $\mu_{\chi A}=m_{\chi}M_{A}/(m_{\chi}+M_{A})$ is the reduced mass of DM-
nucleus system, and $E_{\rm th}$ is the energy threshold of the detector under
consideration.

To determine the direct detection bounds on the coupling $g^{\prime}$ we begin
by first considering the $\varepsilon_{\text{UV}}=0$. Because we are
interested in the region of parameter space where $m_{Z^{\prime}}\gtrsim 1$
GeV, we can neglect the momentum transfer in propagator of the scattering
cross section and directly translate the experimental upper limit on the DM-
nucleon scattering cross section into limits on $g^{\prime}$ as a function of
$m_{Z^{\prime}}$, for a for a fixed value of $m_{\chi}$, using Eq. (13).

Taking into account the strongest constraints from the LUX-ZEPELIN (LZ)
experiment LUX-ZEPLIN:2022qhg , the upper limit for
$\varepsilon_{\text{UV}}=0$ is depicted in Fig. (4) by the dashed blue lines.
We can see that direct detection constraints completely exclude the region of
parameter space where the DM gets the relic abundance for all four DM masses
that we consider, and that the model originally discussed in
Altmannshofer:2016jzy  can no longer simultaneously address the $B$ anomalies
and DM.

Figure 3: Differential scattering rate of DM scattering with Xenon for the
boundary condition $\varepsilon_{\text{IR}}=0$ (sold curves) and
$\varepsilon_{\text{UV}}=0$ (dashed curves), and fixed $m_{Z^{\prime}}=$50
GeV. The red, coral, blue, and green curves are for $m_{\chi}$ = 10, 25, 50,
and 100 GeV, respectively. Note, the dip at $E_{R}\simeq 100$ keV is due to
the behavior of the Helm form factor $F_{\text{Helm}}$ in the scattering rate
and is independent of the model under consideration.

Figure 4: Constraints of the $g^{\prime}$ vs $m_{Z^{\prime}}$ parameter space
for $m_{\chi}=$ 10 GeV (top left), 25 GeV (top right), 50 GeV (bottom left),
and 100 GeV (bottom right). The red, coral, dark blue, and green curves are
the same as those of Fig. (2), where the DM obtains the correct relic
abundance. The gray shaded regions enclosed by solid gray lines are existing
constraints from neutrino trident production at CCFR Altmannshofer:2014pba , a
CMS search for $Z\to 4\mu$ CMS:2018yxg , and from $B_{s}$ meson mixing. The
dashed gray line is the upper bound from a shift in the $Z$ boson mass for the
$\varepsilon_{\text{IR}}=0$ case. The purple shaded regions are constraints
from energy injection into the CMB. The solid blue and dashed blue curves are
DM direct detection constraints from LZ for $\varepsilon_{\text{IR}}=0$ and
$\varepsilon_{\text{UV}}=0$, respectively.

However, it can not be overemphasized that this conclusion holds only for the
case where the kinetic mixing vanishes at high scales i.e.
$\varepsilon_{\text{UV}}=0$. This boundary condition of the kinetic mixing, as
stated in Sec. II, is not generic and the DM direct detection constraints can
change significantly if the boundary condition is changed.

To this end, we now consider the case $\varepsilon_{\text{IR}}=0$ where the
kinetic mixing instead vanishes as low scales. To determine the upper bound on
$g^{\prime}$ we need to rescale the bounds on the $\varepsilon_{\text{UV}}$ =
0 case by

$\Bigg{(}\frac{R_{\varepsilon_{\text{UV}}=0}}{R_{\varepsilon_{\text{IR}}=0}}\Bigg{)}^{1/4}$
(16)

where $R$ is the scattering rate given by Eq. (IV.2) evaluated for a
particular boundary condition. Note, this ratio is independent of
$m_{Z^{\prime}}$ once the DM mass is fixed. Furthermore, because the
scattering rate for $\varepsilon_{\text{IR}}=0$ is suppressed at low momentum
transfer compared to $\varepsilon_{\text{UV}}=0$, the ratio in Eq. (16) is
always greater than 1. This is illustrated in Fig. (3) where we show the
differential scattering rate $dR/dE_{R}$ for $\varepsilon_{\text{IR}}=0$
(solid curves) and $\varepsilon_{\text{UV}}=0$ (dashed curves). We fix
$m_{Z^{\prime}}=50$ GeV, and the different colored curves correspond to
different DM masses. The suppression in the scattering rate can be seen
clearly and implies that DM direct detection constrains will be weaker for
$\varepsilon_{\text{IR}}=0$ .

The results of this rescaling on DM direct detection constraints are depicted
in Fig. (4) by the solid blue lines. We see that the DM direct detection
constraints can be about an order of magnitude weaker for
$\varepsilon_{\text{IR}}=0$ compared to $\varepsilon_{\text{UV}}=0$ case.
Parameter space where DM obtains the correct relic abundance becomes viable
for $m_{\chi}=$ 25, 50, and 100 GeV, and does not have to occur near the
resonance $m_{\chi}\simeq 2m_{Z^{\prime}}$. For $m_{\chi}=10$ GeV we see that
the relic abundance curve is excluded by a combination of existing constraints
from CCFR Altmannshofer:2014pba  and CMS CMS:2018yxg , but we keep it as an
illustration of how the direct detection bounds vary when moving from
$\varepsilon_{\text{UV}}=0$ to $\varepsilon_{\text{IR}}=0$.

To summarize, Fig. (4) shows that, by changing the boundary condition of the
total kinetic mixing, we are able to open up parameter space where the
$Z^{\prime}$ gauge boson can address the $B$ anomalies and, at the same time,
serves as a mediator for DM that can obtain the correct relic abundance via
thermal freeze out. Furthermore, the annihilation of DM in the early universe
does not have to occur at the resonance, allowing a much wider parameter space
for thermal DM.

###  IV.3 Indirect Detection

The parameter space for $U(1)_{L_{\mu}-L_{\tau}}$ charged DM can also be
probed by indirect detection. In particular, for $m_{\chi}>m_{\mu,\tau}$ DM
will annihilate to visible particles with a cross section given by Eq. (11c),
and this can inject energy into the cosmic microwave background (CMB) and
change the ionization history of the universe. Precise measurements of the CMB
power spectrum can set strong constraints on the DM annihilation cross section
Slatyer:2015jla ; Planck:2018vyg ; Kawasaki:2021etm .

These constraints are shown in Fig. (4) by the purple shaded regions. We see
that for $m_{\chi}=10$ GeV the CMB is probing almost all of the relic
abundance curve, while for $m_{\chi}=$ 25, 50, and 100 GeV energy injection
into the CMB is constraining only the region where DM annihilates on resonance
$m_{\chi}\simeq 2m_{Z^{\prime}}$.

However, it is important to note that DM annihilation near a resonance is
sensitive to the temperature at which the annihilation occurs. The
annihilation cross section at the temperature of freeze-out can be different
than when the CMB is formed Griest:1990kh ; Gondolo:1990dk ; Edsjo:1997bg ;
Ibe:2008ye ; Guo:2009aj . Therefore, the CMB constraints near the resonance
should only be taken as an order of magnitude estimate. A more precise
treatment is beyond the scope of this work, but will ultimately not have an
effect on the CMB constrains outside of the resonance.

##  V Electroweak Precision Constraints

In the previous section we considered the kinetic mixing between the
$Z^{\prime}$ and SM photon since we are interested in DM direct detection,
which is a low momentum process. Generally, however, the $Z^{\prime}$ will mix
with the $Z$ boson at high scales and can shift the values of $m_{Z}$ and the
$Z$ coupling to fermions away from their SM values. The shifts in the
properties of the $Z$ boson can have an effect on predictions of electroweak
(EW) precision observables, such as SM gauge boson masses and widths $m_{Z,W}$
and $\Gamma_{Z,W}$, the $Z$-peak hadronic cross section
$\sigma^{0}_{\text{had}}$, $Z$ boson partial width ratios $R^{0}_{\ell,c,b}$,
and forward-backward asymmetries $A^{f,0}_{FB}$.

Measurements of these observables from LEP and the LHC are in good agreement
with SM predictions Baak:2014ora ; Haller:2018nnx , and can be used to
constrain the parameter space of this model. A full treatment of the effect of
$Z-Z^{\prime}$ mixing on EW precision observables is beyond the scope of this
paper, but approximate constraints can be obtained by considering the effect
on the predicted value of $Z$ boson mass.

At high energies, where the $Z$ boson mass is precisely measured, the total
kinetic mixing is approximately given by

$\varepsilon_{\text{tot}}(Q)\simeq\varepsilon_{0}+\frac{eg^{\prime}}{2\pi^{2}}\frac{m^{2}_{\tau}}{Q^{2}}\
.$ (17)

We can ignore the second term since $Q^{2}\gg m^{2}_{\tau}$ at LEP and LHC
energies. For the $\varepsilon_{\text{UV}}=0$ case, we have
$\varepsilon_{0}=0$ and there is no constraint from EW precision observables.
On the other hand, $\varepsilon_{0}$ is nonzero for the
$\varepsilon_{\text{IR}}=0$ case, and we can place constraints on
$g^{\prime}$.

To find the physical masses of the $Z$ and $Z^{\prime}$ gauge bosons, we need
to diagonalize the kinetic terms in Eq. (1) and the resulting neutral gauge
boson mass-squared matrix. The physical masses are given by Curtin:2014cca

$\displaystyle m^{2}_{Z,Z^{\prime}}$
$\displaystyle=\frac{m^{2}_{Z,0}}{2}\Bigg{(}1+\delta^{2}+\eta^{2}\sin^{2}\theta_{W}$
$\displaystyle\pm\text{Sign}(1-\delta^{2})\sqrt{(1+\delta^{2}+\eta^{2}\sin^{2}\theta_{W})^{2}-4\delta^{2}}\Bigg{)},$
(18)

where $m^{2}_{Z,0}=(g^{2}_{1}+g^{2}_{2})v^{2}/4$ is the $Z$ boson mass before
mixing, with $g_{1},g_{2}$ and $v\simeq 246$ GeV being the $U(1)_{Y}$ gauge
coupling, $SU(2)_{L}$ gauge coupling, and the SM Higgs vacuum expectation
value, respectively. The parameter $\delta$ is given by the relation
$m^{2}_{Z^{\prime},0}\equiv m^{2}_{Z,0}\times\delta^{2}$ where
$m^{2}_{Z^{\prime},0}$ is the $Z^{\prime}$ mass before mixing, and
$\eta\equiv\varepsilon_{0}/\sqrt{1-\varepsilon_{0}^{2}}$. Finally,
$\theta_{W}$ is the Weinberg angle.

To determine the constraints on the $U(1)_{L_{\mu}-L_{\tau}}$ parameter space,
we fix $m_{Z,0}$ to the PDG fit value of $m_{Z}=91.1876\pm 0.0021$
Workman:2022ynf , that is we assume that the SM perfectly agrees with the
experimental measurement. Allowing for a 2$\sigma$ shift in the measured value
of the $Z$ boson mass, we obtain the upper limits shown in Fig. (4) by the
dashed gray curves labeled “EWPO”. We see that constraints from the shift in
$m_{Z}$ are generally weaker than constraints from CCFR and CMS, except in the
region where $m_{Z^{\prime}}\simeq m_{Z}$.

##  VI Conclusions

In this work, we revisited the simultaneous explanation of the $B$ anomalies
and DM in the $U(1)_{L_{\mu}-L_{\tau}}$ extension of the SM presented in
Altmannshofer:2016jzy  with a more general treatment of the kinetic mixing
between the $Z^{\prime}$ and the photon. In general, the kinetic mixing
depends on the momentum transfer of a particular process and different
boundary conditions can lead to drastically different behavior of the total
kinetic mixing.

Typically, it is assumed that the total kinetic mixing vanishes at high energy
scales and we have seen in Fig. (4) that direct detection constraints,
together with constraints from $B_{s}$ meson mixing, completely exclude the DM
relic abundance curves, and a simultaneous explanation of $B$ anomalies and DM
in this scenario is not possible. However, we have seen that there is a
freedom to choose the boundary condition of the total kinetic mixing, and the
picture changes significantly if we instead require that the kinetic mixing
vanishes at low scales.

The scattering rate of low momentum transfer processes such as DM-nucleon
scattering can be substantially reduced, and we saw that strong direct
detection constraints on the parameter space of this model can be lifted, as
shown in Fig. (4). We found that the parameter space that is favored by both
the $B$ anomalies and DM becomes viable when considering this alternative
boundary condition for the kinetic mixing. These results emphasize the utility
of treating the momentum dependence of the kinetic mixing in a way that is
different than what is typically studied.

###### Acknowledgements.

We thank Tim Hapitas for collaboration in the early stages of this work. We
are grateful to Wolfgang Altmannshofer, Stefania Gori, and Yue Zhang for
helpful discussions and valuable feedback on the manuscript. The work of DT is
supported by the Arthur B. McDonald Canadian Astroparticle Physics Research
Institute.

## References

  * (1) X.-G. He, G. C. Joshi, H. Lew, and R. R. Volkas, “Simplest Z-prime model,” Phys. Rev. D 44 (1991) 2118–2132.
  * (2) X. G. He, G. C. Joshi, H. Lew, and R. R. Volkas, “NEW Z-prime PHENOMENOLOGY,” Phys. Rev. D 43 (1991) 22–24.
  * (3) S. Baek, N. G. Deshpande, X. G. He, and P. Ko, “Muon anomalous g-2 and gauged L(muon) - L(tau) models,” Phys. Rev. D 64 (2001) 055006, arXiv:hep-ph/0104141.
  * (4) W. Altmannshofer, S. Gori, S. Profumo, and F. S. Queiroz, “Explaining dark matter and B decay anomalies with an $L_{\mu}-L_{\tau}$ model,” JHEP 12 (2016) 106, arXiv:1609.04026 [hep-ph].
  * (5) G. Arcadi, T. Hugle, and F. S. Queiroz, “The Dark $L_{\mu}-L_{\tau}$ Rises via Kinetic Mixing,” Phys. Lett. B 784 (2018) 151–158, arXiv:1803.05723 [hep-ph].
  * (6) A. Kamada, K. Kaneta, K. Yanagi, and H.-B. Yu, “Self-interacting dark matter and muon $g-2$ in a gauged U$(1)_{L_{\mu}-L_{\tau}}$ model,” JHEP 06 (2018) 117, arXiv:1805.00651 [hep-ph].
  * (7) P. Foldenauer, “Light dark matter in a gauged $U(1)_{L_{\mu}-L_{\tau}}$ model,” Phys. Rev. D 99 (2019) no. 3, 035007, arXiv:1808.03647 [hep-ph].
  * (8) D. Borah, M. Dutta, S. Mahapatra, and N. Sahu, “Muon $(g-2)$ and XENON1T Excess with Boosted Dark Matter in $L_{\mu}-L_{\tau}$ Model,” arXiv:2104.05656 [hep-ph].
  * (9) I. Holst, D. Hooper, and G. Krnjaic, “The Simplest and Most Predictive Model of Muon $g-2$ and Thermal Dark Matter,” arXiv:2107.09067 [hep-ph].
  * (10) W. Altmannshofer, S. Gori, M. Pospelov, and I. Yavin, “Quark flavor transitions in $L_{\mu}-L_{\tau}$ models,” Phys. Rev. D 89 (2014) 095033, arXiv:1403.1269 [hep-ph].
  * (11) A. Crivellin, G. D’Ambrosio, and J. Heeck, “Explaining $h\to\mu^{\pm}\tau^{\mp}$, $B\to K^{*}\mu^{+}\mu^{-}$ and $B\to K\mu^{+}\mu^{-}/B\to Ke^{+}e^{-}$ in a two-Higgs-doublet model with gauged $L_{\mu}-L_{\tau}$,” Phys. Rev. Lett. 114 (2015) 151801, arXiv:1501.00993 [hep-ph].
  * (12) A. Crivellin, J. Fuentes-Martin, A. Greljo, and G. Isidori, “Lepton Flavor Non-Universality in B decays from Dynamical Yukawas,” Phys. Lett. B 766 (2017) 77–85, arXiv:1611.02703 [hep-ph].
  * (13) P. Ko, T. Nomura, and H. Okada, “Explaining $B\to K^{(*)}\ell^{+}\ell^{-}$ anomaly by radiatively induced coupling in $U(1)_{\mu-\tau}$ gauge symmetry,” Phys. Rev. D 95 (2017) no. 11, 111701, arXiv:1702.02699 [hep-ph].
  * (14) W. Altmannshofer and I. Yavin, “Predictions for lepton flavor universality violation in rare B decays in models with gauged $L_{\mu}-L_{\tau}$,” Phys. Rev. D 92 (2015) no. 7, 075022, arXiv:1508.07009 [hep-ph].
  * (15) W. Altmannshofer and P. Stangl, “New Physics in Rare B Decays after Moriond 2021,” arXiv:2103.13370 [hep-ph].
  * (16) T. Hapitas, D. Tuckler, and Y. Zhang, “General kinetic mixing in gauged U(1)L$\mu$-L$\tau$ model for muon g-2 and dark matter,” Phys. Rev. D 105 (2022) no. 1, 016014, arXiv:2108.12440 [hep-ph].
  * (17) Belle Collaboration, J. T. Wei et al., “Measurement of the Differential Branching Fraction and Forward-Backward Asymmetry for $B\to K^{(*)}\ell^{+}\ell^{-}$,” Phys. Rev. Lett. 103 (2009) 171801, arXiv:0904.0770 [hep-ex].
  * (18) CDF Collaboration, T. Aaltonen et al., “Observation of the Baryonic Flavor-Changing Neutral Current Decay $\Lambda_{b}\to\Lambda\mu^{+}\mu^{-}$,” Phys. Rev. Lett. 107 (2011) 201802, arXiv:1107.3753 [hep-ex].
  * (19) BaBar Collaboration, J. P. Lees et al., “Measurement of Branching Fractions and Rate Asymmetries in the Rare Decays $B\to K^{(*)}l^{+}l^{-}$,” Phys. Rev. D 86 (2012) 032012, arXiv:1204.3933 [hep-ex].
  * (20) LHCb Collaboration, R. Aaij et al., “Differential branching fractions and isospin asymmetries of $B\to K^{(*)}\mu^{+}\mu^{-}$ decays,” JHEP 06 (2014) 133, arXiv:1403.8044 [hep-ex].
  * (21) LHCb Collaboration, R. Aaij et al., “Differential branching fraction and angular analysis of $\Lambda^{0}_{b}\rightarrow\Lambda\mu^{+}\mu^{-}$ decays,” JHEP 06 (2015) 115, arXiv:1503.07138 [hep-ex]. [Erratum: JHEP 09, 145 (2018)].
  * (22) LHCb Collaboration, R. Aaij et al., “Angular analysis and differential branching fraction of the decay $B^{0}_{s}\to\phi\mu^{+}\mu^{-}$,” JHEP 09 (2015) 179, arXiv:1506.08777 [hep-ex].
  * (23) CMS Collaboration, V. Khachatryan et al., “Angular analysis of the decay $B^{0}\to K^{*0}\mu^{+}\mu^{-}$ from pp collisions at $\sqrt{s}=8$ TeV,” Phys. Lett. B 753 (2016) 424–448, arXiv:1507.08126 [hep-ex].
  * (24) LHCb Collaboration, R. Aaij et al., “Measurements of the S-wave fraction in $B^{0}\rightarrow K^{+}\pi^{-}\mu^{+}\mu^{-}$ decays and the $B^{0}\rightarrow K^{\ast}(892)^{0}\mu^{+}\mu^{-}$ differential branching fraction,” JHEP 11 (2016) 047, arXiv:1606.04731 [hep-ex]. [Erratum: JHEP 04, 142 (2017)].
  * (25) Belle Collaboration, S. Wehle et al., “Lepton-Flavor-Dependent Angular Analysis of $B\to K^{\ast}\ell^{+}\ell^{-}$,” Phys. Rev. Lett. 118 (2017) no. 11, 111801, arXiv:1612.05014 [hep-ex].
  * (26) LHCb Collaboration, R. Aaij et al., “Test of lepton universality with $B^{0}\rightarrow K^{*0}\ell^{+}\ell^{-}$ decays,” JHEP 08 (2017) 055, arXiv:1705.05802 [hep-ex].
  * (27) CMS Collaboration, A. M. Sirunyan et al., “Measurement of angular parameters from the decay $\mathrm{B}^{0}\to\mathrm{K}^{*0}\mu^{+}\mu^{-}$ in proton-proton collisions at $\sqrt{s}=$ 8 TeV,” Phys. Lett. B 781 (2018) 517–541, arXiv:1710.02846 [hep-ex].
  * (28) ATLAS Collaboration, M. Aaboud et al., “Angular analysis of $B^{0}_{d}\rightarrow K^{*}\mu^{+}\mu^{-}$ decays in $pp$ collisions at $\sqrt{s}=8$ TeV with the ATLAS detector,” JHEP 10 (2018) 047, arXiv:1805.04000 [hep-ex].
  * (29) LHCb Collaboration, R. Aaij et al., “Search for lepton-universality violation in $B^{+}\to K^{+}\ell^{+}\ell^{-}$ decays,” Phys. Rev. Lett. 122 (2019) no. 19, 191801, arXiv:1903.09252 [hep-ex].
  * (30) Belle Collaboration, A. Abdesselam et al., “Test of Lepton-Flavor Universality in ${B\to K^{\ast}\ell^{+}\ell^{-}}$ Decays at Belle,” Phys. Rev. Lett. 126 (2021) no. 16, 161801, arXiv:1904.02440 [hep-ex].
  * (31) BELLE Collaboration, S. Choudhury et al., “Test of lepton flavor universality and search for lepton flavor violation in $B\rightarrow K\ell\ell$ decays,” JHEP 03 (2021) 105, arXiv:1908.01848 [hep-ex].
  * (32) LHCb Collaboration, R. Aaij et al., “Measurement of $CP$-Averaged Observables in the $B^{0}\rightarrow K^{*0}\mu^{+}\mu^{-}$ Decay,” Phys. Rev. Lett. 125 (2020) no. 1, 011802, arXiv:2003.04831 [hep-ex].
  * (33) LHCb Collaboration, R. Aaij et al., “Angular Analysis of the $B^{+}\rightarrow K^{\ast+}\mu^{+}\mu^{-}$ Decay,” Phys. Rev. Lett. 126 (2021) no. 16, 161802, arXiv:2012.13241 [hep-ex].
  * (34) LHCb Collaboration, R. Aaij et al., “Test of lepton universality in beauty-quark decays,” Nature Phys. 18 (2022) no. 3, 277–282, arXiv:2103.11769 [hep-ex].
  * (35) LHCb Collaboration, R. Aaij et al., “Branching Fraction Measurements of the Rare $B^{0}_{s}\rightarrow\phi\mu^{+}\mu^{-}$ and $B^{0}_{s}\rightarrow f_{2}^{\prime}(1525)\mu^{+}\mu^{-}$\- Decays,” Phys. Rev. Lett. 127 (2021) no. 15, 151801, arXiv:2105.14007 [hep-ex].
  * (36) LHCb Collaboration, R. Aaij et al., “Tests of lepton universality using $B^{0}\to K^{0}_{S}\ell^{+}\ell^{-}$ and $B^{+}\to K^{*+}\ell^{+}\ell^{-}$ decays,” Phys. Rev. Lett. 128 (2022) no. 19, 191802, arXiv:2110.09501 [hep-ex].
  * (37) LHCb Collaboration, R. Aaij et al., “Search for hidden-sector bosons in $B^{0}\\!\to K^{*0}\mu^{+}\mu^{-}$ decays,” Phys. Rev. Lett. 115 (2015) no. 16, 161802, arXiv:1508.04094 [hep-ex].
  * (38) T. Inami and C. S. Lim, “Effects of Superheavy Quarks and Leptons in Low-Energy Weak Processes k(L) —$>$ mu anti-mu, K+ —$>$ pi+ Neutrino anti-neutrino and K0 $<$—$>$ anti-K0,” Prog. Theor. Phys. 65 (1981) 297. [Erratum: Prog.Theor.Phys. 65, 1772 (1981)].
  * (39) J. Charles, S. Descotes-Genon, Z. Ligeti, S. Monteil, M. Papucci, K. Trabelsi, and L. Vale Silva, “New physics in $B$ meson mixing: future sensitivity and limitations,” Phys. Rev. D 102 (2020) no. 5, 056023, arXiv:2006.04824 [hep-ph].
  * (40) CCFR Collaboration, S. R. Mishra et al., “Neutrino tridents and W Z interference,” Phys. Rev. Lett. 66 (1991) 3117–3120.
  * (41) W. Altmannshofer, S. Gori, M. Pospelov, and I. Yavin, “Neutrino Trident Production: A Powerful Probe of New Physics with Neutrino Beams,” Phys. Rev. Lett. 113 (2014) 091801, arXiv:1406.2332 [hep-ph].
  * (42) CMS Collaboration, A. M. Sirunyan et al., “Search for an $L_{\mu}-L_{\tau}$ gauge boson using Z$\to 4\mu$ events in proton-proton collisions at $\sqrt{s}=$ 13 TeV,” Phys. Lett. B 792 (2019) 345–368, arXiv:1808.03684 [hep-ex].
  * (43) BaBar Collaboration, J. P. Lees et al., “Search for a muonic dark force at BABAR,” Phys. Rev. D 94 (2016) no. 1, 011102, arXiv:1606.03501 [hep-ex].
  * (44) K. J. Kelly, M. Sen, W. Tangarife, and Y. Zhang, “Origin of sterile neutrino dark matter via secret neutrino interactions with vector bosons,” Phys. Rev. D 101 (2020) no. 11, 115031, arXiv:2005.03681 [hep-ph].
  * (45) Planck Collaboration, N. Aghanim et al., “Planck 2018 results. VI. Cosmological parameters,” Astron. Astrophys. 641 (2020) A6, arXiv:1807.06209 [astro-ph.CO]. [Erratum: Astron.Astrophys. 652, C4 (2021)].
  * (46) R. H. Helm, “Inelastic and Elastic Scattering of 187-Mev Electrons from Selected Even-Even Nuclei,” Phys. Rev. 104 (1956) 1466–1475.
  * (47) G. Jungman, M. Kamionkowski, and K. Griest, “Supersymmetric dark matter,” Phys. Rept. 267 (1996) 195–373, arXiv:hep-ph/9506380.
  * (48) G. Duda, A. Kemper, and P. Gondolo, “Model Independent Form Factors for Spin Independent Neutralino-Nucleon Scattering from Elastic Electron Scattering Data,” JCAP 04 (2007) 012, arXiv:hep-ph/0608035.
  * (49) M. Hoferichter, J. Menéndez, and A. Schwenk, “Coherent elastic neutrino-nucleus scattering: EFT analysis and nuclear responses,” Phys. Rev. D 102 (2020) no. 7, 074018, arXiv:2007.08529 [hep-ph].
  * (50) T. Lin, “Dark matter models and direct detection,” PoS 333 (2019) 009, arXiv:1904.07915 [hep-ph].
  * (51) LUX-ZEPLIN Collaboration, J. Aalbers et al., “First Dark Matter Search Results from the LUX-ZEPLIN (LZ) Experiment,” arXiv:2207.03764 [hep-ex].
  * (52) T. R. Slatyer, “Indirect dark matter signatures in the cosmic dark ages. I. Generalizing the bound on s-wave dark matter annihilation from Planck results,” Phys. Rev. D 93 (2016) no. 2, 023527, arXiv:1506.03811 [hep-ph].
  * (53) M. Kawasaki, H. Nakatsuka, K. Nakayama, and T. Sekiguchi, “Revisiting CMB constraints on dark matter annihilation,” JCAP 12 (2021) no. 12, 015, arXiv:2105.08334 [astro-ph.CO].
  * (54) K. Griest and D. Seckel, “Three exceptions in the calculation of relic abundances,” Phys. Rev. D 43 (1991) 3191–3203.
  * (55) P. Gondolo and G. Gelmini, “Cosmic abundances of stable particles: Improved analysis,” Nucl. Phys. B 360 (1991) 145–179.
  * (56) J. Edsjo and P. Gondolo, “Neutralino relic density including coannihilations,” Phys. Rev. D 56 (1997) 1879–1894, arXiv:hep-ph/9704361.
  * (57) M. Ibe, H. Murayama, and T. T. Yanagida, “Breit-Wigner Enhancement of Dark Matter Annihilation,” Phys. Rev. D 79 (2009) 095009, arXiv:0812.0072 [hep-ph].
  * (58) W.-L. Guo and Y.-L. Wu, “Enhancement of Dark Matter Annihilation via Breit-Wigner Resonance,” Phys. Rev. D 79 (2009) 055012, arXiv:0901.1450 [hep-ph].
  * (59) Gfitter Group Collaboration, M. Baak, J. Cúth, J. Haller, A. Hoecker, R. Kogler, K. Mönig, M. Schott, and J. Stelzer, “The global electroweak fit at NNLO and prospects for the LHC and ILC,” Eur. Phys. J. C 74 (2014) 3046, arXiv:1407.3792 [hep-ph].
  * (60) J. Haller, A. Hoecker, R. Kogler, K. Mönig, T. Peiffer, and J. Stelzer, “Update of the global electroweak fit and constraints on two-Higgs-doublet models,” Eur. Phys. J. C 78 (2018) no. 8, 675, arXiv:1803.01853 [hep-ph].
  * (61) D. Curtin, R. Essig, S. Gori, and J. Shelton, “Illuminating Dark Photons with High-Energy Colliders,” JHEP 02 (2015) 157, arXiv:1412.0018 [hep-ph].
  * (62) Particle Data Group Collaboration, R. L. Workman and Others, “Review of Particle Physics,” PTEP 2022 (2022) 083C01.