# Phase transitions in the one-dimensional ionic Hubbard model

Myung-Hoon Chung<EMAIL_ADDRESS>College of Science and Technology,
Hongik University, Sejong 339-701, Korea

###### Abstract

We study quantum phase transitions by measuring the bond energy, the number
density, and the half-chain entanglement entropy in the one-dimensional ionic
Hubbard model. By performing the infinite density matrix renormalization group
with matrix product operator, we obtain ground states as the canonical form of
matrix product states. Depending on the chemical potential and the staggered
potential, the number density and the half-chain entanglement entropy shows
clear signatures of the Mott transition. Our results confirm the success of
the matrix product operator method for investigation of itinerant fermion
systems.

##  I Introduction

Intrinsic characteristics of quantum entanglements Osterloh  _et al._ (2002);
Osborne and Nielsen (2002) are essential elements in the study of quantum
phase transitions Sachdev (2011). In recent years, there has been a great deal
of interest in entanglement entropy Eisert  _et al._ (2010), which is a
quantum phase transition marker Cha and Chung (2018) in the interacting spin
Tagliacozzo  _et al._ (2008); Pollmann  _et al._ (2010); Pirvu  _et al._
(2012), boson Pino  _et al._ (2012), and fermion Gu  _et al._ (2004); Larsson
and Johannesson (2005); Iemini  _et al._ (2015); Cha (2018) systems.

Some cold atoms can obey the same statistical rules as electrons if atoms with
an odd number of neutrons are fermions Parsons  _et al._ (2015). Thus, some
cold atoms Cheuk  _et al._ (2015) in an optical lattice can mimic the behavior
of electrons in a real solid material. Since cold fermion systems are invented
in quasi one-dimensional optical lattices Schreiber  _et al._ (2015), the one-
dimensional fermion Hubbard model can be realized physically. The one-
dimensional fermion Hubbard model Essler  _et al._ (2005) was analytically
solved and became a prototype of a playground, where theoretical results could
be compared to physical reality.

In this paper, motivated by cold atoms, we study the one-dimensional Hubbard
model with an energy offset between even and odd lattice sites. The main
process in this study aims to see the role of different potential energies,
which is made by the second laser disturbance in the one-dimensional cold atom
systems. In fact, we impose the potential energy as

$V(x)=V_{\lambda}\sin(kx)+V_{2\lambda}\sin(\frac{k}{2}x),$ (1)

using the second laser of the doubly long wavelength. For example, we present
the staggered potential with typical values as shown in Fig. 1.

Figure 1: The shape of the potential energy made by two wavelength lasers
$\lambda$ and $2\lambda$. Explicitly, we plot $V(x)=5\sin(4\pi x)+\sin(2\pi
x)$. We notice the staggered potential $\Delta$.

The ionic Hubbard model Torbati  _et al._ (2014) stemmed from this staggered
potential. There have been many investigations into the ionic Hubbard model at
half filling, by using dynamical mean field theory Bag  _et al._ (2015); Lin
_et al._ (2015), determinant quantum Monte Carlo Bouadim  _et al._ (2007),
exact diagonalization and finite density matrix renormalization group Kampf
_et al._ (2003).

The purpose of this paper is to find a quantum phase diagram in the one-
dimensional ionic Hubbard model away from half filling. To do this, we propose
the ansatz of matrix product state (MPS) for the ground state, and perform the
infinite density matrix renormalization group (iDMRG) McCulloch (2008) method
with matrix product operator (MPO). By changing the chemical potential and the
staggered potential, we calculate the entanglement entropy, which shows
quantum phase transition. The occupation number density also indicates the
transition consistently. We determine the phase diagram in the parameter
space.

##  II Model and method

Introducing the staggered one-particle potential $\Delta$ into the Hubbard
model, we consider a simple quantum system of interacting fermions on a one-
dimensional lattice. The Hamiltonian for this model is written as

$\displaystyle H$ $\displaystyle=$ $\displaystyle-t\sum_{\langle
ij\rangle}(c^{\dagger}_{i\uparrow}c_{j\uparrow}+c^{\dagger}_{j\uparrow}c_{i\uparrow}+c^{\dagger}_{i\downarrow}c_{j\downarrow}+c^{\dagger}_{j\downarrow}c_{i\downarrow})$
(2)
$\displaystyle+U\sum_{i}(n_{i\uparrow}-\frac{1}{2})(n_{i\downarrow}-\frac{1}{2})-\mu\sum_{i}(n_{i\uparrow}+n_{i\downarrow})$
$\displaystyle+\frac{\Delta}{2}\sum_{i\in\text{even}}(n_{i\uparrow}+n_{i\downarrow})-\frac{\Delta}{2}\sum_{i\in\text{odd}}(n_{i\uparrow}+n_{i\downarrow}),$

where $\langle ij\rangle$ represents the nearest-neighbor hopping in a one-
dimensional lattice, and $n_{i\uparrow}$ and $n_{i\downarrow}$ are the spin-up
and the spin-down number operator, respectively. We let the hopping strength
$t$ be 1 and vary the strengths of the parameters: the on-site Coulomb
repulsion $U$, the chemical potential $\mu$, and the staggered potential
$\Delta$. The role of the chemical potential is to control the number of
fermions in the system. For the case of $\Delta=0$, the exact solution of the
original Hubbard model is known Yang  _et al._ (2000).

In order to handle $\mu$ and $\Delta$, we introduce two linear combinations
$\mu_{e}$ and $\mu_{o}$ for even and odd sites, denoting

$\mu_{e}\equiv\mu-\frac{\Delta}{2},~{}~{}~{}\mu_{o}\equiv\mu+\frac{\Delta}{2}.$
(3)

Because the Hamiltonian has the two-site translational symmetry as shown in
Fig. 1, we apply iDMRG with unit cell of two sites to determine the ground
state. Our strategy is to use MPO acting on MPS.

In relation to the basis, we set the physical index $\sigma_{i}$ in the MPS.
The state on the $i$-th site is represented by
$\sigma_{i}=(\alpha_{i},\beta_{i})$, such as $0=(0,0)$ for the empty,
$1=(0,1)$ for the spin-down, $2=(1,0)$ for the spin-up, and $3=(1,1)$ for both
of the spin-up and the spin-down. Thus, the physical index runs from $0$ to
$3$. The basis of the Fock space is written in terms of creation operators
$c^{\dagger}_{i\uparrow}$ and $c^{\dagger}_{i\downarrow}$ as follows:

$\displaystyle|\sigma_{0}\cdots\sigma_{L-1}\rangle$ $\displaystyle=$
$\displaystyle(c^{\dagger}_{0\uparrow})^{\alpha_{0}}(c^{\dagger}_{0\downarrow})^{\beta_{0}}(c^{\dagger}_{1\uparrow})^{\alpha_{1}}(c^{\dagger}_{1\downarrow})^{\beta_{1}}$
(4)
$\displaystyle\cdots(c^{\dagger}_{L-1\uparrow})^{\alpha_{L-1}}(c^{\dagger}_{L-1\downarrow})^{\beta_{L-1}}|0\rangle,$

where $\alpha_{i}(\beta_{i})=0$ or $1$ means that there is a spin-up (down)
fermion vacancy or occupancy at the $i$-th site, respectively. It is important
to maintain the order of fermions in the state of the Fock space to deal with
the negative sign caused by fermion exchange. We adopt the order of spin-up
first and then spin-down as described above.

In the process of iDMRG, we have to find the MPO for the Hamiltonian of Eq.
(2). We convince that the MPO of the $i$-th site is indeed given by the
following form:

$M_{i}=\left(\begin{array}[]{cccccc}1&c^{\dagger}_{i\uparrow}&c^{\dagger}_{i\downarrow}&c_{i\uparrow}&c_{i\downarrow}&b_{i}\\\
0&0&0&0&0&-tc_{i\uparrow}\\\ 0&0&0&0&0&-tc_{i\downarrow}\\\
0&0&0&0&0&tc^{\dagger}_{i\uparrow}\\\ 0&0&0&0&0&tc^{\dagger}_{i\downarrow}\\\
0&0&0&0&0&1\end{array}\right)$ (5)

where

$b_{i}=U(n_{i\uparrow}-\frac{1}{2})(n_{i\downarrow}-\frac{1}{2})-\bar{\mu}(n_{i\uparrow}+n_{i\downarrow}),$
(6) $\bar{\mu}=\left\\{\begin{array}[]{cc}\mu_{e}&\text{for }i=\text{even,}\\\
\mu_{o}&\text{for }i=\text{odd.}\\\ \end{array}\right.$ (7)

Here, we omit to present the boundary operators. We note that the MPO for even
sites is slightly different from the MPO for odd sites, where $\mu_{e}$ must
be replaced with $\mu_{o}$. Because $c$ and $c^{\dagger}$ are fermion
operators, the commutation relations between $M_{i}$ and
$c^{\dagger}_{j\sigma}$ ($i\neq j,~{}\sigma=\uparrow,\downarrow$) are given by
$M_{i}c_{j\sigma}^{\dagger}=c_{j\sigma}^{\dagger}M_{i}^{*}$, where

$M_{i}^{*}=\left(\begin{array}[]{cccccc}1&-c^{\dagger}_{i\uparrow}&-c^{\dagger}_{i\downarrow}&-c_{i\uparrow}&-c_{i\downarrow}&b_{i}\\\
0&0&0&0&0&tc_{i\uparrow}\\\ 0&0&0&0&0&tc_{i\downarrow}\\\
0&0&0&0&0&-tc^{\dagger}_{i\uparrow}\\\
0&0&0&0&0&-tc^{\dagger}_{i\downarrow}\\\ 0&0&0&0&0&1\end{array}\right)$ (8)

When the Hamiltonian acts on the basis, the sign caused by fermion exchange
must be taken care of such as

$\displaystyle M_{0}\cdots
M_{L-1}(c^{\dagger}_{0\uparrow})^{\alpha_{0}}(c^{\dagger}_{0\downarrow})^{\beta_{0}}(c^{\dagger}_{1\uparrow})^{\alpha_{1}}(c^{\dagger}_{1\downarrow})^{\beta_{1}}$
$\displaystyle\cdots(c^{\dagger}_{L-1\uparrow})^{\alpha_{L-1}}(c^{\dagger}_{L-1\downarrow})^{\beta_{L-1}}|0\rangle~{}~{}~{}~{}~{}~{}~{}$
$\displaystyle~{}~{}~{}=M_{0}(c^{\dagger}_{0\uparrow})^{\alpha_{0}}(c^{\dagger}_{0\downarrow})^{\beta_{0}}~{}\bar{M}_{1}(c^{\dagger}_{1\uparrow})^{\alpha_{1}}(c^{\dagger}_{1\downarrow})^{\beta_{1}}\cdots$
$\displaystyle~{}\bar{M}_{L-1}(c^{\dagger}_{L-1\uparrow})^{\alpha_{L-1}}(c^{\dagger}_{L-1\downarrow})^{\beta_{L-1}}|0\rangle,~{}~{}~{}~{}~{}~{}~{}~{}$
(9)

where the total number of fermions in the left side of the $i$-th site,
$N_{i}\equiv\sum_{k=0}^{i-1}(\alpha_{k}+\beta_{k})$, determines the $6\times
6$ matrix:

$\bar{M}_{i}=\left\\{\begin{array}[]{cc}M_{i}&\text{for }N_{i}=\text{even,}\\\
M_{i}^{*}&\text{for }N_{i}=\text{odd.}\\\ \end{array}\right.$ (10)

In order to manage the sequence of Eq. (9), we note that the even and odd
structure of the fermion numbers $(E~{}~{}O)$ is kept recursively by
$e\equiv\\{|0\rangle,c^{\dagger}_{\uparrow}c^{\dagger}_{\downarrow}|0\rangle\\}$
and
$o\equiv\\{c^{\dagger}_{\uparrow}|0\rangle,c^{\dagger}_{\downarrow}|0\rangle\\}$
such as

$\displaystyle(E~{}~{}O)\left(\begin{array}[]{cc}e&o\\\ o&e\\\
\end{array}\right)\left(\begin{array}[]{cc}e&o\\\ o&e\\\
\end{array}\right)\cdots\left(\begin{array}[]{cc}e&o\\\ o&e\\\
\end{array}\right)$ (17)
$\displaystyle=(E~{}~{}O)\left(\begin{array}[]{cc}e&o\\\ o&e\\\
\end{array}\right)\cdots\left(\begin{array}[]{cc}e&o\\\ o&e\\\
\end{array}\right).$ (22)

Using the relation of Eq. (22), we double the MPO for the $i$-th site,
$\bar{M}_{i}$, which is modified into $12\times 12$ matrix, $D_{i}$, such as

$D_{i}e=\left(\begin{array}[]{cc}M_{i}e&0\\\ 0&M_{i}^{*}e\\\
\end{array}\right),~{}~{}D_{i}o=\left(\begin{array}[]{cc}0&M_{i}o\\\
M_{i}^{*}o&0\\\ \end{array}\right),$ (23)

where $M_{i}^{*}$ is used in the second row because there are odd number of
fermions in the left side of the $i$-th site. By including the $4$-dimensional
physical index, we present the effective MPO as the graphical diagram with the
four-leg tensor $D$ having $4\times 12\times 12\times 4$ elements, of which
the number of non-zero elements is small:

$\begin{array}[]{c}\\\ \cdots\\\ \\\ \end{array}\begin{array}[]{c}\\\ -\\\ \\\
\end{array}\begin{array}[]{c}|\\\ D\\\ |\\\ \end{array}\begin{array}[]{c}\\\
-\\\ \\\ \end{array}\begin{array}[]{c}|\\\ D\\\ |\\\
\end{array}\begin{array}[]{c}\\\ -\\\ \\\ \end{array}\begin{array}[]{c}|\\\
D\\\ |\\\ \end{array}\begin{array}[]{c}\\\ -\\\ \\\
\end{array}\begin{array}[]{c}\\\ \cdots\\\ \\\ \end{array}$ (24)

It is remarkable that the so-called fermionic matrix product states Bultinck
_et al._ (2017) are introduced in a similar fashion of doubling the MPO. We
have already tested this MPO approach Chung  _et al._ (2019), and the
numerical results by doubling the MPO were compared with those obtained by the
Jordan-Wigner transformation in good agreement.

By adopting the two-site iDMRG algorithm, we optimize two tensors in our MPS,
$A_{ab}^{\sigma}$ and $B_{ab}^{\sigma}$. The physical index $\sigma$ takes a
value of $0$ to $3$, and the indices $a$ (left) and $b$ (right) run from $0$
to $\chi-1$, where $\chi$ is the internal bond dimension. The Schmidt
coefficients are denoted by $\lambda^{AB}_{b}$ between $A_{ab}^{\sigma}$ and
$B_{bc}^{\rho}$, and by $\lambda^{BA}_{b}$ between $B_{ab}^{\sigma}$ and
$A_{bc}^{\rho}$. A state in the MPS space is represented by graphic notation
such as

$\begin{array}[]{c}\\\ \cdots\\\ \end{array}\begin{array}[]{c}\\\ -\\\
\end{array}\begin{array}[]{c}|\\\ A\\\ \end{array}\begin{array}[]{c}\\\ -\\\
\end{array}\begin{array}[]{c}\\\ \lambda^{AB}\\\
\end{array}\begin{array}[]{c}\\\ -\\\ \end{array}\begin{array}[]{c}|\\\ B\\\
\end{array}\begin{array}[]{c}\\\ -\\\ \end{array}\begin{array}[]{c}\\\
\lambda^{BA}\\\ \end{array}\begin{array}[]{c}\\\ -\\\
\end{array}\begin{array}[]{c}|\\\ A\\\ \end{array}\begin{array}[]{c}\\\ -\\\
\end{array}\begin{array}[]{c}\\\ \cdots\\\ \end{array}$ (25)

By evaluating $4\times 12\times 12\times 4$ elements of $D$, we perform iDMRG
by acting the MPO of Eq. (24) on the MPS of Eq. (25). Regardless of the values
of $t$, $U$, $\mu_{e}$, and $\mu_{o}$ used in our calculations, we have
observed a smooth convergence. Our numerical iDMRG results show that the
ground-state solutions are divided into two classes: the MPS are either of the
form $\cdots ABABAB\cdots$, which we will identify as a Mott phase with unit
cell of two sites, near half filling, or uniform $\cdots AAAAAA\cdots$, which
we will identify as metallic, further away from half filling.

##  III Numerical Results

By setting $t=1$ and $U=4$, for a given $\mu_{e}$ and $\mu_{o}$, we calculate
the ground state by the two-site iDMRG. With the ground state, we extract
three physical quantities corresponding to even and odd sites: the bond energy
$\langle h_{e}\rangle$ and $\langle h_{o}\rangle$, the number density $\langle
n_{e}\rangle$ and $\langle n_{o}\rangle$, and the half-chain entanglement
entropy $S_{e}$ and $S_{o}$. Here, the bond Hamiltonian
$h_{e}~{}(i=\text{even})$ is written as

$\displaystyle h_{e}$ $\displaystyle=$
$\displaystyle-t(c^{\dagger}_{i\uparrow}c_{j\uparrow}+c^{\dagger}_{j\uparrow}c_{i\uparrow}+c^{\dagger}_{i\downarrow}c_{j\downarrow}+c^{\dagger}_{j\downarrow}c_{i\downarrow})$
(26)
$\displaystyle+\frac{U}{2}(n_{i\uparrow}-\frac{1}{2})(n_{i\downarrow}-\frac{1}{2})-\frac{\mu_{e}}{2}(n_{i\uparrow}+n_{i\downarrow})$
$\displaystyle+\frac{U}{2}(n_{j\uparrow}-\frac{1}{2})(n_{j\downarrow}-\frac{1}{2})-\frac{\mu_{o}}{2}(n_{j\uparrow}+n_{j\downarrow}).$

Also, the occupation number density operator $n_{e}~{}(i=\text{even})$ is
written as

$n_{e}=n_{i\uparrow}+n_{i\downarrow}=c^{\dagger}_{i\uparrow}c_{i\uparrow}+c^{\dagger}_{i\downarrow}c_{i\downarrow}.$
(27)

Using the Schmidt coefficients $\lambda_{a}$ between $A$ and $B$ in the MPS,
we determine $S_{e}$. Keeping the normalization with
$\sum^{\chi-1}_{a=0}\lambda^{2}_{a}=1$, we calculate the half-chain
entanglement entropy $S_{e}$, which is given by

$S_{e}=-\sum^{\chi-1}_{a=0}\lambda^{2}_{a}\log_{2}\lambda^{2}_{a}.$ (28)

Obviously, we choose $i=\text{odd}$ for $h_{o}$ and $n_{o}$. We use the
Schmidt coefficients $\lambda^{BA}$ between $B$ and $A$ in the MPS for
$S_{o}$.

Figure 2: Ground-state bond energy $(\langle h_{e}\rangle+\langle
h_{o}\rangle)/2$ versus $\mu_{o}$ showing a level-crossing between Mott-
insulating and metallic solutions, at $\mu_{e}=0.65$ and $\chi=100$. We set
$t=1$ and $U=4$. Inset: the number density of $\langle n_{e}\rangle$, $\langle
n_{o}\rangle$ and $(\langle n_{e}\rangle+\langle n_{o}\rangle)/2$ versus
$\mu_{o}$ showing the Mott-insulating states, at $\mu_{e}=0$ and $\chi=100$.
Although the difference of $(\langle n_{e}\rangle-\langle n_{o}\rangle)$
increases as the staggered potential $\Delta=(\mu_{o}-\mu_{e})$ becomes
bigger, the average number density $(\langle n_{e}\rangle+\langle
n_{o}\rangle)/2$ remains fixed as $1$ in the Mott-insulating states.

In order to verify the correctness of the calculation, we use an important
analytic result Essler  _et al._ (2005) from the critical point of $\mu_{+}$
at which the Mott insulator - metallic phase transition takes place in the
restricted parameter space of $\mu_{e}=\mu_{o}$ or $\Delta=0$. In fact, the
transition point $\mu_{+}$ is given by

$\mu_{+}=-2+\frac{U}{2}+2\int_{0}^{\infty}\frac{d\omega}{\omega}\frac{J_{1}(\omega)\exp(-\omega
U/4)}{\cosh(\omega U/4)}.$ (29)

For $U=4$, we obtain the exact critical point,
$\mu_{+}=\mu_{e}=\mu_{o}=0.643364$. We present the average bond energy
$(\langle h_{e}\rangle+\langle h_{o}\rangle)/2$ as shown in Fig. 2, where the
level-crossing takes place at $\mu_{e}=\mu_{o}\sim 0.65$.

In Fig. 2, to determine the level-crossing of the average bond energy, we
compute the ground state starting from $\mu_{o}=0$. On the other hand, we can
also calculate the ground state of the uniform solution by decreasing
$\mu_{o}$ from large values. In each iDMRG calculation, the tensors of the
initial environment are given by the previous solution of the different
$\mu_{o}$.

In order to determine the Mott transition in the parameter space of $\mu_{e}$
and $\mu_{o}$, we compute the average number density $(\langle
n_{e}\rangle+\langle n_{o}\rangle)/2$ and the average entanglement entropy
$(S_{e}+S_{o})/2$. Fixing $\mu_{o}$, we present $(\langle n_{e}\rangle+\langle
n_{o}\rangle)/2$ and $(S_{e}+S_{o})/2$ by varying $\mu_{e}$, as shown in Fig.
3 and 4. We note that the abrupt behaviour corresponds to the Mott transition.

Figure 3: The number density $(\langle n_{e}\rangle+\langle n_{o}\rangle)/2$
versus $\mu_{e}$ showing the abrupt change near $\mu_{e}=0.35$, at
$\mu_{o}=1$, by setting $t=1$, $U=4$ and $\chi=100$. We note that the phase
transition takes place near $\mu_{e}=0.35$. Inset: By setting $\mu_{o}=-1$,
the number density versus $\mu_{e}$ showing the phase transition near
$\mu_{e}=-0.35$.  Figure 4: The half-chain entanglement entropy
$(S_{e}+S_{o})/2$ versus $\mu_{e}$ showing the jump near $\mu_{e}=0.35$, at
$\mu_{o}=1$, by setting $t=1$, $U=4$ and $\chi=100$. We note that the phase
transition takes place near $\mu_{e}=0.35$. Inset: By setting $\mu_{o}=-1$,
the half-chain entanglement entropy versus $\mu_{e}$ showing the phase
transition near $\mu_{e}=-0.35$.

As the main numerical result, we present the phase diagram in Fig. 5. One
notices the robustness of the Mott state. In fact, we find that a small
difference of energy offsets, which is the staggered potential
$\Delta=\mu_{o}-\mu_{e}$, can not destroy the Mott phase.

Figure 5: Numerical phase diagram of the one-dimensional ionic Hubbard model.
We have set $t=1$ and $U=4$. In the region of the Mott-insulating phase, the
number density is given by $1$, and the entanglement entropy is almost a
constant. The red line is a guide to the eyes. The region of the Mott-
insulating phase shrinks for smaller $U$.

##  IV Conclusion

In summary, we have used iDMRG with MPO to obtain the ground state in the one-
dimensional ionic Hubbard model. We calculate the occupation number density
and the half-chain entanglement entropy to determine the Mott transition. We
draw the phase diagram in the two-dimensional parameter space of the chemical
potential and the staggered potential.

It is of interest to extend our method to manage the number of fermions in
two-dimensional systems. A typical topic of interest may be the two-
dimensional fermion Hubbard model, for which our first task is to build a
local tensor product operator like MPO. Although we do not encounter the
notorious sign problem here, we should overcome the sign problem in the two-
dimensional Hubbard model with the tensor product operator. As we double MPO
as explained above, we can double the tensor product operator in order to
support the number of fermions. We anticipate progress in the two-dimensional
case.

###### Acknowledgements.

This work was partially supported by the Basic Science Research Program
through the National Research Foundation of Korea (NRF) funded by the Ministry
of Education, Science and Technology (Grant No. NRF-2017R1D1A1A0201845). The
author would like to thank M. C. Cha for helpful discussions.

## References

  * Osterloh  _et al._ (2002) A. Osterloh, L. Amico, G. Falci, and R. Fazio, Nature 416, 608 (2002).
  * Osborne and Nielsen (2002) T. J. Osborne and M. A. Nielsen, Phys. Rev. A 66, 032110 (2002).
  * Sachdev (2011) S. Sachdev,  _Quantum Phase Transitions_ (Cambridge University Press, 2011).
  * Eisert  _et al._ (2010) J. Eisert, M. Cramer, and M. B. Plenio, Rev. Mod. Phys. 82, 277 (2010).
  * Cha and Chung (2018) M. C. Cha and M. H. Chung, Physica B 536, 701 (2018).
  * Tagliacozzo  _et al._ (2008) L. Tagliacozzo, T. R. de Oliveira, S. Iblisdir, and J. I. Latorre, Phys. Rev. B 78, 024410 (2008).
  * Pollmann  _et al._ (2010) F. Pollmann, A. M. Turner, E. Berg, and M. Oshikawa, Phys. Rev. B 81, 064439 (2010).
  * Pirvu  _et al._ (2012) B. Pirvu, G. Vidal, F. Verstraete, and L. Tagliacozzo, Phys. Rev. B 86, 075117 (2012).
  * Pino  _et al._ (2012) M. Pino, J. Prior, A. M. Somoza, D. Jaksch, and S. R. Clark, Phys. Rev. A 86, 023631 (2012).
  * Gu  _et al._ (2004) S.-J. Gu, S.-S. Deng, Y.-Q. Li, and H.-Q. Lin, Phys. Rev. Lett. 93, 086402 (2004).
  * Larsson and Johannesson (2005) D. Larsson and H. Johannesson, Phys. Rev. Lett. 95, 196406 (2005).
  * Iemini  _et al._ (2015) F. Iemini, T. O. Maciel, and R. O. Vianna, Phys. Rev. B 92, 075423 (2015).
  * Cha (2018) M. C. Cha, Phys. Rev. B 98, 235161 (2018).
  * Parsons  _et al._ (2015) M. F. Parsons, F. Huber, A. Mazurenko, C. S. Chiu, W. Setiawan, K. Wooley-Brown, and M. G. S. Blatt, Phys. Rev. Lett. 114, 213002 (2015).
  * Cheuk  _et al._ (2015) L. W. Cheuk, M. A. Nichols, M. Okan, T. Gersdorf, V. V. Ramasesh, W. S. Bakr, T. Lompe, and M. W. Zwierlein, Phys. Rev. Lett. 114, 193001 (2015).
  * Schreiber  _et al._ (2015) M. Schreiber, S. S. Hodgman, P. Bordia, H. P. Lüschen, M. H. Fischer, R. Vosk, E. Altman, U. Schneider, and I. Bloch, Science 21, 842 (2015).
  * Essler  _et al._ (2005) F. H. L. Essler, H. Frahm, F. Göhmann, A. Klümper, and V. E. Korepin,  _The one-dimensional Hubbard model_ (Cambridge University Press, 2005).
  * Torbati  _et al._ (2014) M. H. Torbati, N. A. Drescher, and G. S. Uhrig, Phys. Rev. B 89, 245126 (2014).
  * Bag  _et al._ (2015) S. Bag, A. Garg, and H. R. Krishnamurthy, Phys. Rev. B 91, 235108 (2015).
  * Lin  _et al._ (2015) H. F. Lin, H. D. Liu, H. S. Tao, and W. M. Liu, Scientific Reports 5, 9810 (2015).
  * Bouadim  _et al._ (2007) K. Bouadim, N. Paris, F. Hebert, G. G. Batrouni, and R. T. Scalettar, Phys. Rev. B 76, 085112 (2007).
  * Kampf  _et al._ (2003) A. P. Kampf, M. Sekania, G. I. Japaridze, and P. Brune, J of Phys.: Condensed Matter 15, 5895 (2003).
  * McCulloch (2008) I. P. McCulloch, e-print arXiv:0804.2509 (2008).
  * Yang  _et al._ (2000) C. Yang, A. N. Kocharian, and Y. L. Chiang, J. Phys.: Condens. Matter 12, 7433 (2000).
  * Bultinck  _et al._ (2017) N. Bultinck, D. J. Williamson, J. Haegeman, and F. Verstraete, Phys. Rev. B 95, 075108 (2017).
  * Chung  _et al._ (2019) M. H. Chung, E. Orignac, D. Poilblanc, and S. Capponi, e-print arXiv:1912.10203 (2019).


CTPU-PTC-21-14

$\geq 4\mu$ signal from a vector-like lepton decaying to a muon-philic
$Z^{\prime}$ boson at the LHC

Junichiro Kawamuraa,b<EMAIL_ADDRESS>and Stuart Rabyc<EMAIL_ADDRESS>

a Center for Theoretical Physics of the Universe, Institute for Basic Science
(IBS), Daejeon 34051, Korea

bDepartment of Physics, Keio University, Yokohama 223-8522, Japan

c Department of Physics, Ohio State University, Columbus, Ohio, 43210, USA

We propose a novel possibility to detect a very distinctive signal with more
than four muons originating from pair-produced vector-like leptons decaying to
a muon-philic $Z^{\prime}$ boson. These new particles are good candidates to
explain the anomalies in the muon anomalous magnetic moment and the $b\to
s\ell\ell$ processes. The doublet (singlet) vector-like leptons lighter than
1.3 (1.0) TeV are excluded by the latest data at the LHC if $\mathrm{BR}(E\to
Z^{\prime}\mu)=1$. We also show that the excess in the signal region with more
than five leptons can be explained by this scenario if the vector-like lepton
is a weak singlet, with mass about 400 GeV and $\mathrm{BR}(E\to
Z^{\prime}\mu)=0.25$. The future prospects at the HL-LHC are discussed.

##  1 Introduction

The Large Hadron Collider (LHC) explores new physics beyond the Standard Model
(SM) at TeV-scale. The SM has been established as the theory just above the
electroweak (EW) scale, particularly by the discovery of the 125 GeV Higgs
boson at the LHC [1, 2]. Although most of the experiments are consistent with
the predictions of the SM, there are 2-3 $\sigma$ discrepancies in the
measurements of rare semi-leptonic $B$ meson decays [3, 4, 5, 6, 7, 8, 9, 10,
11, 12, 13, 14, 15, 16, 17, 18, 19], $b\to s\ell\ell$, and the $4.2\sigma$
discrepancy in the anomalous magnetic moment of muon, $\Delta a_{\mu}$ [20,
21, 22]. An interesting coincidence here is that both anomalies are found in
physics related to muons, and hence these could be explained by the same
origin. One way to establish these discrepancies as evidence of new physics is
to increase the significance by reducing the uncertainties in the experimental
measurements and in the predictions of the SM. Another way is by directly
discovering new particles at the LHC, which we pursue in this paper.

It was shown in Refs. [23, 24] that both anomalies in $b\to s\ell\ell$ and
$\Delta a_{\mu}$ are addressed by introducing vector-like (VL) fermions and a
${Z^{\prime}}$ boson associated with an additional gauge symmetry
$U(1)^{\prime}$ 333 See Refs. [25, 26, 27, 28, 29, 30] for models with VL
fermions and $U(1)^{\prime}$ for the anomalies. The VL lepton explanation for
$\Delta a_{\mu}$ is studied in e.g. Refs. [31, 32, 33, 34, 35, 36] . The
former anomaly is explained by ${Z^{\prime}}$ exchange at the tree-level 444
The $b\to s\ell\ell$ anomaly can be explained by loop corrections involving VL
families [37, 38, 39, 40, 41, 42, 43, 44, 45]. , while the latter is explained
by loop corrections involving the VL leptons and the ${Z^{\prime}}$ boson. In
this paper, we point out the possibility that pair productions of VL leptons
provide very distinctive signals with more than four muons. We shall discuss
the current limits from the recent ATLAS data [46] and future prospects at the
HL-LHC in a simplified model with a VL lepton and ${Z^{\prime}}$ boson. We
also discuss the possible explanation for the excess in the more than
$5$-lepton signal found in Ref. [46].

The rest of this paper is organized as follows. The simplified model is
defined and then the relation to the anomalies are discussed in Section 2. In
Section 3, we discuss limits from the high-multiplicity lepton signal at the
LHC. Section 4 is devoted to summary. The model proposed in Refs. [23, 24] are
reviewed in Appendix A as a UV completion of the simplified model.

##  2 Simplified model

We shall consider the simplified model with a VL lepton $E$, which is weak
singlet-like, $E_{1}$, or doublet-like, $L=(E_{2},N)$, where $N$ is the
$SU(2)_{L}$ partner of $E_{2}$. The ${Z^{\prime}}$ boson couplings to the
leptons are given by

$\displaystyle\mathcal{L}_{Z^{\prime}}=$ $\displaystyle\
Z^{\prime}_{\mu}\begin{pmatrix}\overline{\mu}&\overline{E}\end{pmatrix}\gamma^{\mu}\left[\begin{pmatrix}g^{L}_{\mu\mu}&g^{L}_{\mu
E}\\\ g^{L}_{\mu E}&g^{L}_{EE}\\\
\end{pmatrix}P_{L}+\begin{pmatrix}g^{R}_{\mu\mu}&g^{R}_{\mu E}\\\ g^{R}_{\mu
E}&g^{R}_{EE}\\\ \end{pmatrix}P_{R}\right]\begin{pmatrix}{\mu}\\\
{E}\end{pmatrix}$ (2.1) $\displaystyle\
+Z^{\prime}_{\mu}\begin{pmatrix}\overline{\nu}&\overline{N}\end{pmatrix}\gamma^{\mu}\left[\begin{pmatrix}g^{L}_{\nu\nu}&g^{L}_{\nu
N}\\\ g^{L}_{\nu N}&g^{L}_{NN}\\\ \end{pmatrix}P_{L}+\begin{pmatrix}0&0\\\
0&g^{R}_{NN}\\\ \end{pmatrix}P_{R}\right]\begin{pmatrix}{\nu}\\\
{N}\end{pmatrix},$

where $E=E_{1}$ or $E_{2}$ and the interactions with $N$ in the second line
are absent in the case of weak singlet VL lepton. We assume that the off-
diagonal couplings of the SM bosons to the SM and VL leptons are negligible,
such that the dominant decay modes of the VL leptons are the decays to a
${Z^{\prime}}$ boson and SM lepton. In fact, this is achieved in the model
proposed in Refs. [23, 24].

The loop corrections involving the VL leptons and the ${Z^{\prime}}$ boson
contribute to the anomalous magnetic moment of the muon. It is known that the
chiral-flip effect should be sizable to explain the current discrepancy of
$\mathcal{O}\left({10^{-9}}\right)$ with the new particles above the EW scale.
In models with VL leptons, the chiral-flip effects may come from the non-zero
VEV of the SM Higgs doublet. Hence the size of the loop correction is
estimated as

$\displaystyle\Delta a_{\mu}\sim$ $\displaystyle\ -\frac{m_{\mu}\kappa
v_{H}}{8\pi^{2}m_{Z^{\prime}}^{2}}g^{L}_{\mu E}g^{R}_{\mu E}C_{\Delta
a_{\mu}}$ (2.2) $\displaystyle\sim$ $\displaystyle\ 2.9\times
10^{-9}\times\left(\frac{500~{}\mathrm{GeV}}{m_{Z^{\prime}}}\right)^{2}\left(\frac{\kappa}{0.5}\right)\left(\frac{\sqrt{g^{L}_{\mu
E}g^{R}_{\mu E}}}{0.25}\right)^{2}\left(\frac{C_{\Delta
a_{\mu}}}{0.1}\right),$

where $\kappa$ is the Yukawa couping constant for
$\tilde{H}\overline{L}_{R}E_{L}$. $C_{\Delta a_{\mu}}$ is the factor from loop
functions which is typically of $\mathcal{O}\left({0.1}\right)$, see Appendix
A for the explicit form in the example model.

The ${Z^{\prime}}$ boson couplings to muons, $g_{\mu\mu}^{L}$ and
$g_{\mu\mu}^{R}$, directly relate to the Wilson coefficients for the $b\to
s\ell\ell$ decay. The effective Hamiltonian is given by [47, 48]

$\displaystyle\mathcal{H}_{\mathrm{eff}}=-\frac{4G_{F}}{\sqrt{2}}\frac{\alpha_{e}}{4\pi}V_{tb}V_{ts}^{*}\left(C_{9}\mathcal{O}_{9}+C_{10}\mathcal{O}_{10}\right),$
(2.3)

where

$\displaystyle\mathcal{O}_{9}:=\bigl{[}\overline{s}\gamma^{\mu}P_{L}b\bigr{]}\bigl{[}\overline{\mu}\gamma_{\mu}\mu\bigr{]},\quad\mathcal{O}_{10}:=\bigl{[}\overline{s}\gamma^{\mu}P_{L}b\bigr{]}\bigl{[}\overline{\mu}\gamma_{\mu}\gamma_{5}\mu\bigr{]}.$
(2.4)

The coefficients in this model are given by

$\displaystyle C_{9}=$ $\displaystyle\
-\frac{\sqrt{2}}{4G_{F}}\frac{4\pi}{\alpha_{e}}\frac{1}{V_{tb}V_{ts}^{*}}\frac{g^{L}_{sb}}{2m^{2}_{Z^{\prime}}}\left(g^{R}_{\mu\mu}+g^{L}_{\mu\mu}\right),$
(2.5) $\displaystyle C_{10}=$ $\displaystyle\
-\frac{\sqrt{2}}{4G_{F}}\frac{4\pi}{\alpha_{e}}\frac{1}{V_{tb}V_{ts}^{*}}\frac{g^{L}_{sb}}{2m^{2}_{Z^{\prime}}}\left(g^{R}_{\mu\mu}-g^{L}_{\mu\mu}\right),$
(2.6)

where $g^{L}_{sb}$ is the ${Z^{\prime}}$ couplings to $\overline{s}b$ in the
left-current. The value of $C_{9}$ is estimated as

$\displaystyle\left|{C_{9}}\right|\sim
0.87\times\left(\frac{m_{Z^{\prime}}}{500~{}\mathrm{GeV}}\right)^{2}\left(\frac{g^{L}_{sb}}{0.0007}\right)\left(\frac{g^{L}_{\mu\mu}+g^{R}_{\mu\mu}}{0.5}\right).$
(2.7)

Note that the ${Z^{\prime}}$ couplings to quarks are tiny to explain the $b\to
s\ell\ell$ anomaly, while those to muons are large so that $\Delta a_{\mu}$ is
explained when $g_{\mu\mu}^{L,R}\sim g_{\mu E}^{L,R}$, which is true in the
sample model. This feature ensures that the ${Z^{\prime}}$ mass of
$\mathcal{O}\left({100~{}\mathrm{GeV}}\right)$ is not excluded by the di-
lepton resonance search at the LHC [49] 555 See Refs. [50, 51] for general
discussions for ${Z^{\prime}}$ boson responsible for $b\to s\ell\ell$ .

In the ${Z^{\prime}}$ boson explanation, the ratio of the coefficients are
given by

$\displaystyle\frac{C_{10}}{C_{9}}=\frac{g^{R}_{\mu\mu}-g^{L}_{\mu\mu}}{g^{R}_{\mu\mu}+g^{L}_{\mu\mu}}.$
(2.8)

The recent analyses [52, 53] including the measurement of $R_{K}$ based on the
full run-2 data at the LHCb [19] 666 See Refs. [54, 55, 56, 57, 58, 59, 60,
61, 62, 63, 64, 65, 66, 67, 68, 69, 70] for the analyses before the Moriond
2021.  favor $C_{9}$-only, $C_{10}$-only and $C_{9}=-C_{10}$ scenarios among
the one dimensional analyses, which correspond to
$g_{\mu\mu}^{R}=g_{\mu\mu}^{L}$, $g_{\mu\mu}^{R}=-g_{\mu\mu}^{L}$ and
$g^{R}_{\mu\mu}=0$, respectively. Among these three cases, the explanation by
$C_{9}=-C_{10}$ is not preferred to explain the $\Delta a_{\mu}$ anomaly
because $\Delta a_{\mu}\propto g_{\mu\mu}^{L}g_{\mu\mu}^{R}$, see Eq. (2.2).
From this observation, we shall consider the case with
$\left|{g^{L}_{\mu\mu}}\right|=\left|{g^{R}_{\mu\mu}}\right|$ which predicts

$\displaystyle\mathrm{BR}\left({{Z^{\prime}}}\to{\mu\mu}\right):=\frac{\Gamma\left({{Z^{\prime}}}\to{\mu\mu}\right)}{\Gamma\left({{Z^{\prime}}}\to{\nu\nu}\right)+\Gamma\left({{Z^{\prime}}}\to{\mu\mu}\right)}\simeq\frac{\left|{g^{L}_{\mu\mu}}\right|^{2}+\left|{g^{R}_{\mu\mu}}\right|^{2}}{2\left|{g^{L}_{\mu\mu}}\right|^{2}+\left|{g^{R}_{\mu\mu}}\right|^{2}}=\frac{2}{3}.$
(2.9)

Here, we assume $g_{\mu\mu}^{L}=g_{\nu\nu}^{L}$ as expected from the
$SU(2)_{L}$ symmetry and the ${Z^{\prime}}$ boson decay to quarks are
negligible as expected from Eq. (2.7). Studies for
$\left|{g^{L}_{\mu\mu}}\right|\neq\left|{g^{R}_{\mu\mu}}\right|$, as preferred
by the two dimensional analyses on $(C_{9},C_{10})$ plane, are interesting but
are beyond the scope of this paper.

##  3 LHC signals

Figure 1:  Processes produces more than four muons.

In this paper, we study signals from pair produced VL leptons decaying to the
second generation leptons and the ${Z^{\prime}}$ boson. This can be realized
when the VL leptons are heavier than the ${Z^{\prime}}$ boson. If the
${Z^{\prime}}$ boson is heavier, the VL lepton may decay to a SM boson and a
lepton. The limits for VL leptons in such a case are studied in Refs. [71, 72,
73, 74, 75, 76, 77]. It is also possible that the VL lepton decays to a new
boson, such as the physical mode of the $U(1)^{\prime}$ breaking scalar. Thus
we treat the branching fraction of $E/N\to{Z^{\prime}}\mu/\nu$ as a free-
parameter. We further assume
$\mathrm{BR}\left({N}\to{{Z^{\prime}}\nu}\right)=\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)$
for simplicity. Figure 1 shows the relevant processes which can generate
signals for more than four muons 777 We used TikZ-FeynHand to draw these
figures [78, 79]. Only the left process is relevant for the singlet-like case.

We recast the limits obtained in Refs. [46] and [80]. The former searches for
signals with more than four leptons, and the latter searches for signals with
exactly two leptons with large missing transverse energy,
$E_{T}^{\mathrm{miss}}$. We have generated events using
MadGraph5$\\_$2$\\_$8$\\_$2 [81] based on a UFO [82] model file generated with
FeynRules$\\_$2$\\_$3$\\_$43 [83, 81]. The events are showered with PYTHIA8
[84] and then run through the fast detector simulator Delphes3.4.2 [84]. We
used the default ATLAS card for the detector simulation, but the threshold on
$p_{T}$ for the muon efficiency formula is changed to $5$ GeV from $10$ GeV
since muons with $p_{T}>5~{}\mathrm{GeV}$ are counted as signal muons in Ref.
[46].

Table 1:  The number of events observed (data), fitted SM backgrounds (SM) and $95\%$ C.L. upper bound on the number of signal events $(S^{95})$ in the signal regions [46, 80].  | $\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$ | $\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ | SR5L | SR2L
---|---|---|---|---
data | 11 | 1 | 21 | 37
SM | $11.5^{+2.9}_{-2.2}$ | $3.5^{+2.0}_{-2.2}$ | $12.4\pm 2.3$ | $37.3\pm 3.0$
$S^{95}$ | $9.79$ | $3.87$ | $17.88$ | $14.3$

We recast the experimental limits on the signal regions without $Z$ boson,
b-jet and hadronic $\tau$ defined in Ref. [46].888We recast the analysis such
that light leptons, ($e,\mu$), in Ref. [46], for our simple model these are
just muons. These are named $\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$,
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ and SR5L. The requirements
for the events, in addition to the b-jet veto and hadronic $\tau$-veto, common
in the signal regions are as follows. To meet the trigger thresholds, $p_{T}$
of the leading muon, ordered by $p_{T}$, must be larger than
$27~{}\mathrm{GeV}$, or $p_{T}$’s of the leading and next-to leading muons are
required to be larger than $(15,15)~{}\mathrm{GeV}$ or
$(23,9)~{}\mathrm{GeV}$. If an opposite-sign (OS) muon pair whose invariant
mass $m_{\mathrm{OS}}$ is less than $4~{}\mathrm{GeV}$ or
$8.4<m_{\mathrm{OS}}<10.4~{}\mathrm{GeV}$, both leptons are discarded. If two
muons are found in $\Delta R<0.6$ and one of them has
$p_{T}<30~{}\mathrm{GeV}$, both leptons are discarded. The first (second) $Z$
candidate is found from a pair of OS muons whose $m_{\mathrm{OS}}$ is the
(second) closest to the $Z$ boson mass $m_{Z}=91.2~{}\mathrm{GeV}$. A pair is
identified as a $Z$ boson if $m_{\mathrm{OS}}\in[81.2,101.2]$ $\mathrm{GeV}$.
Further, the event is considered to have a $Z$ boson if any system of
$\mu^{+}\mu^{-}\mu^{\pm}$ or $\mu^{+}\mu^{-}\mu^{+}\mu^{-}$ has invariant mass
in $[81.1,101.2]$ GeV. In the signal regions
$\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$ and
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$, there must be more than four
muons after the selections above, and the event must not have any combinations
of muons which is identified as a $Z$ boson. Further, the effective mass of
the event $m_{\mathrm{eff}}$, defined as the scalar sum of
$E_{T}^{\mathrm{miss}}$, $p_{T}$ of signal leptons and $p_{T}$ of the jets
with $p_{T}>40$ GeV, is required to be larger than 600 (1250)$~{}\mathrm{GeV}$
in the $\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$
($\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$). In the SR5L, the
requirement is simply the lepton number to be larger than five, and no further
selection applied.

We also study the limits from the SUSY slepton search [80] which requires
exactly two leptons and large $E_{T}^{\mathrm{miss}}$. The most relevant
signal region for our scenario is with same flavor (SF) two leptons without
any jet. There must be exactly two OSSF leptons, both with
$p_{T}>25~{}\mathrm{GeV}$. Events are rejected if there are more muons with
$p_{T}>10$ GeV and $\left|{\eta}\right|<2.7$ or the two leading leptons are
not opposite sign. The missing energy $E_{T}^{\mathrm{miss}}$ and invariant
mass of two leptons $m_{\mathrm{OS}}$ must be larger than $110~{}\mathrm{GeV}$
and $121.2~{}\mathrm{GeV}$, respectively. The stransverse mass
$m_{\mathrm{T2}}$ [85, 86] is required to be larger than $160$ GeV 999 We used
the code provided by Ref. [87] to calculate the stransverse mass. . We name
this signal region as SR2L.

Figure 2:  Pair-production cross section of VL leptons at $\sqrt{s}=13$ and 14
TeV. The productions of VL neutrinos are included in the doublet case.

The number of observed events, fitted SM backgrounds and $95\%$ C.L. upper
bounds on the signal events in the signal regions are shown in Table 1. We see
that there is an excess over the SM background in SR5L for which the local
significance is $1.9\sigma$. Figure 2 shows the production cross sections of
VL leptons at $\sqrt{s}=13$ TeV and $14$ TeV calculated by MadGraph5. We
calculated the probability of how many events pass the cuts in each signal
region from the VL lepton pair production at $\sqrt{s}=13$ TeV by generating
25000 (50000) events at each point on the $(m_{E},m_{{Z^{\prime}}})$ plane for
the singlet-like (doublet-like) VL leptons.

###  3.1 Current limits

Figure 3:  Limits on $\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)$ for the
doublet-like VL lepton.  Figure 4:  Limits on
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)$ for the singlet-like VL
lepton.

Figures 3 and 4 show upper bounds on
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)$ in the signal regions, where
the signal cross section is proportional to the branching fraction squared. In
the gray region, $m_{Z^{\prime}}>m_{E}$ and hence the decay
$E\to{Z^{\prime}}\mu$ is kinematically forbidden. The white region is not
excluded by the current data even if
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)=1$.

For the doublet-like VL lepton, we see that
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ gives the strongest bound if
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)\gtrsim 0.2$, because of fewer
backgrounds satisfying the tighter $m_{\mathrm{eff}}$ cut. The current limit
is about 1350 GeV if $\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)=1$. If
the branching fraction is smaller, then
$\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$ gives the stronger bound,
since the cut by $m_{\mathrm{eff}}>1250~{}\mathrm{GeV}$ of the
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ is too tight for
$m_{E_{2}}\lesssim 600$ GeV. The limit from SR5L is weaker because of the
excess and that from SR2L is also weaker due to the larger backgrounds. Note
that the limits do not change much as the mass difference between the
${Z^{\prime}}$ boson and VL lepton decreases since the signal muons can
originate from ${Z^{\prime}}$ decays.

For the singlet-like VL lepton, the strongest bound of about 1000 GeV is again
from the $\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ if
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)=1$. The
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ is not sensitive to the cases
when $\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)\lesssim 0.3$. The
difference comes from smaller production cross section of the singlet-like
case. The $\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$ gives the strongest
constraint for smaller branching fractions. The limits from SR5L and SR2L are
not shown, since the limits are much weaker than those from
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ and
$\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$ for the same reason as the
doublet-like case. In particular, SR2L gives no bounds for the singlet-like
case.

###  3.2 Explanation for the excess in SR5L

Figure 5:  Maximum values of the number of signal events in SR5L consistent
with the limits of $\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$ and
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ for the doublet-like
(singlet-like) VL lepton in the left (right) panel. Background colors are the
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)$.

Figure 6:  Histograms of $E_{T}^{\mathrm{miss}}$ (left) and $m_{\mathrm{eff}}$
(right). In the left panel,
$(m_{Z^{\prime}},m_{E_{1}})=(390,400)~{}\mathrm{GeV}$ and
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)=0.25$. In the right panel,
the values are normalized such that the total of the bins is unity.

Figure 5 shows the upper bound on the number of signal events in the SR5L
allowed by the limits from the other signal regions. The background colors
represent maximum values of $\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)$.
Since the limits from the four lepton signal regions are severe, the branching
fraction should be so small that the limits from
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ and
$\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$ are relaxed by the fewer
events passing the $m_{\mathrm{eff}}$ cut. The excess is explained on the
solid green line when the SM background is at the central value shown in
Table. 1. The yellow band corresponds to the uncertainty of the background
estimation. The singlet-like VL lepton can more easily explain the excess. The
limits from $\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ and
$\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$ are much stronger for the
doublet-like case, since the production cross section is larger and there are
fewer muon signals originating from the VL neutrino production.

The left panel of Fig. 6 shows the $E_{T}^{\mathrm{miss}}$ distribution after
the selection of SR5L at the benchmark point with
$m_{Z^{\prime}}=390~{}\mathrm{GeV}$, $m_{E_{1}}=400~{}\mathrm{GeV}$ and
$\mathrm{BR}\left({E_{1}}\to{{Z^{\prime}}\mu}\right)=0.25$. The SM
contributions are represented by yellow bars and the black dots (bars) show
the data (its error bar), which are read from Fig. 8 of Ref. [46]. We see that
the $E_{T}^{\mathrm{miss}}$ distribution is well described by our scenario.
The benchmark point in our example model which realizes these masses and
branching fractions and is also consistent with the muon anomalies at the same
time is shown Appendix A.3. The right panel of Fig. 6 shows the
$m_{\mathrm{eff}}$ distribution after the selection of the number of muons to
be larger than four. The values of the bars are normalized such that the sum
of all the bins is unity. The green bars are for the same benchmark point as
the left panel, and the red hatched bars are for another benchmark point with
$m_{Z^{\prime}}=400~{}\mathrm{GeV}$ and $m_{E_{1}}=1000~{}\mathrm{GeV}$. We
see that the peak of the distribution in $m_{\mathrm{eff}}$ is about
$2m_{E_{1}}$, and hence the strong constraint from
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$, which requires
$m_{\mathrm{eff}}>1250$ GeV, is avoided and the tightest bound of
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)\lesssim 0.25$ is from
$\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$.

We emphasize that this model can only explain the SR5L excess by muons. Thus
the SR5L excess can not be explained, in this scenario, if it includes signals
with electrons. The limits from the data would be significantly tightened if
lepton flavors are specified in the signal regions. Thus the information of
lepton flavor is crucial to test this model with a muon-philic ${Z^{\prime}}$
and VL leptons.

The excess with electrons, might be explained by VL leptons decaying to $Z$ or
$W$ boson, where the SM bosons decay leptonically. If
$m_{Z^{\prime}}>m_{E_{1}}$ which is the opposite case to our scenario, the VL
lepton will decay to a SM boson, including the Higgs boson. In addition,
electrons may come from the decays of the heavier VL leptons, such as
$E_{2}\to E_{1}Z$. We note that the roughly degenerate mass of the VL leptons
are favored to explain the sizable $\Delta a_{\mu}$ in the model [24]. These
possibilities are interesting, but beyond the scope of this paper.

###  3.3 Future prospects

Figure 7:  Future prospects of the upper bound on
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)$ for the doublet-like VL
lepton at the HL-LHC.  Figure 8:  Future prospects of the upper bound on
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)$ for the singlet-like VL
lepton at the HL-LHC.

We shall discuss the discovery and exclusion potential at the HL-LHC with
$3~{}\mathrm{ab}^{-1}$ data. We consider the two signal regions,
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ and SR5L and propose two more
new signal regions, SRZp and $\texttt{SR5L}^{\prime}$.

We rescale the backgrounds to $\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$
and SR5L by simply multiplying the ratio of integrated luminosity, $3000/139$.
For the signal events, we use the same efficiency times acceptance factor as
those used in the analyses for the current limits. These are then multiplied
by the integrated luminosity and the production cross section at
$\sqrt{s}=14~{}\mathrm{TeV}$ shown in Fig. 2 to calculate the number of signal
events.

We define the two new signal regions, named SRZp and $\texttt{SR5L}^{\prime}$.
In the SRZp, at least four muons are required and the $Z$-veto is applied.
Then, two ${Z^{\prime}}$ candidates are chosen from any OS pair of muons,
where the reference ${Z^{\prime}}$ mass, $m_{{Z^{\prime}}}^{\mathrm{ref}}$, is
set at $500~{}\mathrm{GeV}$, in the same manner as the $Z$ candidates. The
(next-to) leading ${Z^{\prime}}$ candidate, $m_{\mathrm{OS}}$ is the (second)
closest to $m_{Z^{\prime}}^{\mathrm{ref}}$, must satisfy
$\left|{m_{\mathrm{OS}}-m_{Z^{\prime}}^{\mathrm{ref}}}\right|<100$ ($250$)
GeV. We take relatively large range for the selection, because we do not know
the ${Z^{\prime}}$ mass. The sensitivity can be improved by requiring strict
range for $\left|{m_{\mathrm{OS}}-m_{Z^{\prime}}^{\mathrm{ref}}}\right|$ and
scan over $m_{Z^{\prime}}^{\mathrm{ref}}$ in the analysis, as in the
${Z^{\prime}}$ searches [49, 88]. In the $\texttt{SR5L}^{\prime}$, more than 5
muons are required. Then, the $Z$-veto for any OS pair of muons and
$m_{\mathrm{eff}}>1000~{}\mathrm{GeV}$ cuts are applied. In these two signal
regions, we assume that there are 10 SM background events per
$3~{}\mathrm{ab}^{-1}$ data. This may be a conservative assumption, since the
cut is very tight and not so many background events will survive, c.f. the
rescaled backgrounds in $\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ and
SR5L are $76.6$ and $272$, respectively.

We quantify the future discovery and exclusion limits by the $p$-values
proposed in Ref. [89],

$\displaystyle p_{\mathrm{disc}}=\frac{\gamma(s+b,b)}{\Gamma(s+b)},\quad
p_{\mathrm{excl}}=\frac{\Gamma(b+1,s+b)}{\Gamma(b+1)},$ (3.1)

where $s$ and $b$ are the number of signals and backgrounds. $\Gamma(z)$,
$\gamma(a,z)$ and $\Gamma(a,z)$ are the ordinary, lower incomplete and upper
complete Gamma functions. The discovery (exclusion) limit corresponds to
$p_{\mathrm{disc}}<2.867\times 10^{-7}$ ($p_{\mathrm{excl}}<0.05$) where the
significance is $>5$ ($>1.645$). Here, we do not consider uncertainties in the
signals and backgrounds for simplicity.

The future prospects at the HL-LHC for the doublet-like and the singlet-like
VL leptons are shown in Fig. 7 and 8, respectively. The background colors are
the exclusion limits ($p_{\mathrm{excl}}<0.05$) for the branching fraction
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)$. The white lines show the
discovery potential for a given branching fraction as labeled on the lines
($p_{\mathrm{disc}}<2.867\times 10^{-7}$). Assuming
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)=1$, the doublet-like
(singlet-like) VL lepton will be discovered up to $m_{E}\lesssim 1.5$ ($1.15$)
TeV by $\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$. The limits from SR5L
are weaker. The SRZp may cover a wider parameter range than that of
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ at $m_{Z^{\prime}}\sim
500~{}\mathrm{GeV}$ even if we set the relatively large range for
$\left|{m_{\mathrm{OS}}-m_{Z^{\prime}}^{\mathrm{ref}}}\right|$. The
$\texttt{SR5L}^{\prime}$ may also cover a wider parameter range independent of
the ${Z^{\prime}}$ mass, based on our assumption of the background. It is
interesting that the entire parameter range with
$\mathrm{BR}\left({E}\to{{Z^{\prime}}\mu}\right)\sim 0.25$, which can explain
the excess in SR5L, can be discovered in SR5L, SRZp and
$\texttt{SR5L}^{\prime}$.

##  4 Discussions

In this paper, we study the signal with more than four muons originating from
the pair-production of VL leptons decaying to a ${Z^{\prime}}$ boson which
couples to muons and/or muon neutrinos. These particles may provide a way to
resolve the tensions in the $b\to s\ell\ell$ decays and $\Delta a_{\mu}$. The
current limits can reach about $1~{}\mathrm{TeV}$ when the VL lepton decays to
the ${Z^{\prime}}$ boson exclusively owing to the very low backgrounds. We
showed that the excess in the signal region with five leptons or more may be
explained in our model if the excess is given only by muons. A benchmark point
in our example model is given in Appendix A.3 which simultaneously explains
the muon anomalies as well as the excess in SR5L. If the excess is not only
muons, then the cascade decay of the heavier VL lepton might be a nice
candidate to explain the excess with electrons and muons in this kind of
model. The information of lepton flavor is crucial to test these new physics
models which explain the muon anomalies.

Note added. While finalizing this manuscript, the new experimental data of the
muon anomalous magnetic moment was announced from the FNAL. The discrepancy
from the SM prediction reaches to $4.2~{}\sigma$ [22]. Models with VL leptons
and a ${Z^{\prime}}$ boson coupling to muons nicely explains the discrepancy
of $\Delta a_{\mu}$, as well as the anomalies in $b\to s\ell\ell$. These new
particles could be confirmed by the LHC as discussed in this paper.

## Acknowledgment

The work of J.K. is supported in part by the Institute for Basic Science
(IBS-R018-D1), and the Grant-in-Aid for Scientific Research from the Ministry
of Education, Science, Sports and Culture (MEXT), Japan No. 18K13534. The work
of S.R. is supported in part by the Department of Energy (DOE) under Award No.
DE-SC0011726.

##  Appendix A Review of vector-like $U(1)^{\prime}$ model

Table 2:  Matter contents. Electric charge of fermion $f$ is $Q_{f}=T_{f}^{3}+Y_{f}/2$.  | $\ell_{L}$ | $\overline{\mu}_{R}$ | $H$ | $L_{L}$ | $\overline{E}_{R}$ | $\overline{L}_{R}$ | $E_{L}$ | $Z^{\prime}$ | $\Phi$
---|---|---|---|---|---|---|---|---|---
$SU(2)_{L}$ | $\mathbf{2}$ | $\mathbf{1}$ | $\mathbf{2}$ | $\mathbf{2}$ | $\mathbf{1}$ | $\mathbf{2}$ | $\mathbf{1}$ | $\mathbf{1}$ | $\mathbf{1}$
$U(1)_{Y}$ | $-1$ | $2$ | $-1$ | $-1$ | $2$ | $1$ | $-2$ | $0$ | $0$
$U(1)^{\prime}$ | $0$ | $0$ | $0$ | $-1$ | $1$ | $1$ | $-1$ | $0$ | $-1$

In this Appendix, we review the model proposed in Refs. [23, 24] as an example
of a UV completion of the simplified model. The matter contents of our model
is given by Table 2. The $SU(2)_{L}$ doublets are defined as

$\displaystyle\ell_{L}=(\nu_{L},\mu_{L}),\quad H=(H_{0},H_{-}),\quad
L_{L}=(N_{L}^{\prime},E_{L}^{\prime}),\quad\overline{L}_{R}=(-\overline{E}_{R}^{\prime},\overline{N}^{\prime}_{R}).$
(A.1)

We only consider muons, and assume that the couplings with the other leptons
are negligible for simplicity. The masses of VL states and Yukawa interactions
are given by

$\displaystyle\mathcal{L}\supset$ $\displaystyle\
-m_{L}\overline{L}_{R}L_{L}-m_{E}\overline{E}_{R}E_{L}$ (A.2) $\displaystyle\
+y_{\mu}\overline{\mu}_{R}\ell_{L}H+{\kappa^{\prime}}\overline{E}_{R}L_{L}H-\kappa\overline{L}_{R}\tilde{H}E_{L}+\lambda_{L}\Phi\overline{L}_{R}\ell_{L}-\lambda_{E}\Phi\overline{\mu}_{R}E_{L}+h.c.,$

where $\tilde{H}:=i\sigma_{2}H^{*}=(H_{-}^{*},-H_{0}^{*})$. The $SU(2)_{L}$
indices are contracted via $i\sigma_{2}$. After the symmetry breaking by
$v_{H}:=\langle{H_{0}}\rangle$ and $v_{\Phi}:=\langle{\Phi}\rangle$, the mass
matrix for the leptons are given by

$\displaystyle\overline{\mathbf{e}}_{R}\mathcal{M}_{e}\mathbf{e}_{L}:=$
$\displaystyle\
\begin{pmatrix}\overline{\mu}_{R}&\overline{E}_{R}&\overline{E}_{R}^{\prime}\end{pmatrix}\begin{pmatrix}y_{\mu}v_{H}&0&\lambda_{E}v_{\Phi}\\\
0&{\kappa^{\prime}}v_{H}&m_{E}\\\ \lambda_{L}v_{\Phi}&m_{L}&\kappa
v_{H}\end{pmatrix}\begin{pmatrix}\mu_{L}\\\ E^{\prime}_{L}\\\
E_{L}\end{pmatrix},$ (A.3)
$\displaystyle\overline{\mathbf{n}}_{R}\mathcal{M}_{n}\mathbf{n}_{L}:=$
$\displaystyle\
\overline{N}_{R}^{\prime}\begin{pmatrix}\lambda_{L}v_{\Phi}&m_{L}\end{pmatrix}\begin{pmatrix}\nu_{L}\\\
N_{L}^{\prime}\end{pmatrix}.$ (A.4)

The mass basis is defined as

$\displaystyle\hat{\mathbf{e}}_{L}:=U_{L}^{\dagger}\mathbf{e}_{L},\quad\hat{\mathbf{e}}_{R}:=U_{R}^{\dagger}\mathbf{e}_{R},\quad\hat{\mathbf{n}}_{L}:=V_{L}^{\dagger}\mathbf{n}_{L},\quad\hat{\mathbf{n}}_{R}:=\mathbf{n}_{R},$
(A.5)

where unitary matrices diagonalize the mass matrices as

$\displaystyle
U_{R}^{\dagger}\mathcal{M}_{e}U_{L}=\mathrm{diag}\left(m_{\mu},m_{E_{2}},m_{E_{1}}\right),\quad\mathcal{M}_{n}V_{L}=\begin{pmatrix}0&m_{N}\end{pmatrix},$
(A.6)

where $E_{1}$ ($E_{2}$) is the singlet-like (doublet-like) VL lepton 101010 We
restrict cases which $m_{E_{1}}\ll m_{E_{2}}$ or $m_{E_{1}}\gg m_{E_{2}}$, so
we can always identify the VL lepton is singlet-like or doublet-like. . The
non-zero mass of the SM neutrino will be explained by introducing the right-
handed counterparts, but these are irrelevant for the present discussion.

We define the Dirac fermions as

$\displaystyle\mathbf{e}:=\left(\mu,E_{2},E_{1}\right),\quad\mathbf{n}:=\left(\nu,N\right),$
(A.7)

where

$\displaystyle\left[\mathbf{e}\right]_{i}:=\left(\left[\hat{\mathbf{e}}_{L}\right]_{i},\left[\hat{\mathbf{e}}_{R}\right]_{i}\right),\quad\nu:=$
$\displaystyle\ \left(\left[\hat{\mathbf{n}}_{L}\right]_{1},0\right),\quad
N:=\left(\left[\hat{\mathbf{n}}_{L}\right]_{2},N^{\prime}_{R}\right),$ (A.8)

with $i=1,2,3$.

###  A.1 Interactions

The gauge interactions with the ${Z^{\prime}}$ boson in the mass basis are
defined as

$\displaystyle\mathcal{L}_{V}=$ $\displaystyle\
{Z^{\prime}}_{\mu}\sum_{f=\mathbf{e},\mathbf{n}}\overline{f}\gamma^{\mu}\left({g}^{Z^{\prime}}_{f_{L}}P_{L}+{g}^{Z^{\prime}}_{f_{R}}P_{R}\right)f,$
(A.9)

where the coupling matrices are given by

$\displaystyle g^{Z^{\prime}}_{\mathbf{e}_{L}}=$ $\displaystyle\
{g^{\prime}}U_{L}^{\dagger}Q^{\prime}_{e}U_{L},\quad
g^{Z^{\prime}}_{\mathbf{e}_{R}}={g^{\prime}}U_{R}^{\dagger}Q^{\prime}_{e}U_{R},\quad
g^{Z^{\prime}}_{\mathbf{n}_{L}}={g^{\prime}}V_{L}^{\dagger}Q^{\prime}_{n}V_{L},\quad
g^{Z^{\prime}}_{\mathbf{n}_{R}}={g^{\prime}}Q^{\prime}_{n}.$ (A.10)

$P_{L}$ ($P_{R}$) are the chiral projections onto the left- (right-)handed
fermions. ${g^{\prime}}$ is the gauge coupling constant for $U(1)^{\prime}$.

We expand the neutral scalar fields as

$\displaystyle
H_{0}=v_{H}+\frac{1}{\sqrt{2}}\left(h+ia_{h}\right),\quad\Phi=v_{\Phi}+\frac{1}{\sqrt{2}}\left(\chi+ia_{\chi}\right),$
(A.11)

where $h$ and $\chi$ are the physical real scalar fields, while the pseudo-
scalar components $a_{h}$ and $a_{\chi}$ are absorbed by the $Z$ and
${Z^{\prime}}$ bosons, respectively. The Yukawa interactions are given by

$\displaystyle-\mathcal{L}_{Y}=\frac{1}{\sqrt{2}}\sum_{S=h,\chi}\sum_{f=\mathbf{e},\mathbf{n}}S\overline{f}\;Y^{S}_{f}P_{L}f+h.c.,$
(A.12)

where

$\displaystyle Y^{h}_{\mathbf{e}}=U_{R}^{\dagger}\begin{pmatrix}y_{\mu}&0&0\\\
0&{\kappa^{\prime}}&0\\\ 0&0&\kappa\end{pmatrix}U_{L},\quad
Y^{\chi}_{\mathbf{e}}=U_{R}^{\dagger}\begin{pmatrix}0&0&\lambda_{E}\\\
0&0&0\\\ \lambda_{L}&0&0\end{pmatrix}U_{L},\quad Y^{h}_{\mathbf{n}}=0_{2\times
2},\quad Y^{\chi}_{\mathbf{n}}=\begin{pmatrix}0&0\\\
\lambda_{L}&0\end{pmatrix}V_{L},$ (A.13)

Let us define the approximate masses of the VL leptons as

$\displaystyle M_{L}:=\sqrt{m_{L}^{2}+\lambda_{L}^{2}v_{\Phi}^{2}},\quad
M_{E}:=\sqrt{m_{E}^{2}+\lambda_{E}^{2}v_{\Phi}^{2}}.$ (A.14)

Assuming $\kappa v_{H}\ll\left|{M_{L}-M_{E}}\right|$, the diagonalization
matrices for the charged lepton mass matrix are given by

$\displaystyle U_{L}\sim\begin{pmatrix}c_{L}&s_{L}&-\delta_{L}s_{L}\\\
-s_{L}&c_{L}&-\delta_{L}c_{L}\\\
0&\delta_{L}&1\end{pmatrix}+\mathcal{O}\left({\delta_{L}^{2}}\right),\quad
U_{R}\sim\begin{pmatrix}c_{R}&s_{R}\delta_{R}&s_{R}\\\
-s_{R}&c_{R}\delta_{R}&c_{R}\\\
0&1&-\delta_{R}\end{pmatrix}+\mathcal{O}\left({\delta_{R}^{2}}\right),$ (A.15)

where

$\displaystyle c_{L}:=$ $\displaystyle\ \frac{m_{L}}{M_{L}},\quad
s_{L}:=\frac{\lambda_{L}v_{\Phi}}{M_{L}},\quad\delta_{L}:=\frac{\kappa
v_{H}M_{L}}{M_{L}^{2}-M_{E}^{2}},$ (A.16) $\displaystyle c_{R}:=$
$\displaystyle\ \frac{m_{E}}{M_{E}},\quad
s_{R}:=\frac{\lambda_{E}v_{\Phi}}{M_{E}},\quad\delta_{R}:=\frac{\kappa
v_{H}M_{E}}{M_{L}^{2}-M_{E}^{2}}.$ (A.17)

The diagonalized mass matrix is given by

$\displaystyle
U_{R}^{\dagger}\mathcal{M}_{e}U_{L}=\begin{pmatrix}\left(y_{\mu}c_{L}c_{R}+{\kappa^{\prime}}s_{L}s_{R}\right)v_{H}&\mathcal{O}\left({m_{\mu}}\right)&0\\\
0&M_{L}+\mathcal{O}\left({\kappa v_{H}\delta_{L,R}}\right)&0\\\
\mathcal{O}\left({m_{\mu}}\right)&\mathcal{O}\left({\kappa
v_{H}\delta^{2}_{L,R}}\right)&M_{E}+\mathcal{O}\left({\kappa
v_{H}\delta_{L,R}}\right)\end{pmatrix}+\mathcal{O}\left({m_{\mu}\delta_{L,R}}\right),$
(A.18)

where we assume $y_{\mu}v_{H},{\kappa^{\prime}}v_{H}\lesssim m_{\mu}$ to
explain the muon mass without fine-tuning.

The ${Z^{\prime}}$ couplings in the mass basis are approximately given by

$\displaystyle g^{Z^{\prime}}_{\mathbf{e}_{L}}=$ $\displaystyle\
-{{g^{\prime}}}\begin{pmatrix}s_{L}^{2}&-c_{L}s_{L}&c_{L}s_{L}\delta_{L}\\\
-c_{L}s_{L}&c_{L}^{2}&s_{L}^{2}\delta_{L}\\\
c_{L}s_{L}\delta_{L}&s_{L}^{2}\delta_{L}&1\end{pmatrix}+\mathcal{O}\left({\delta_{L}^{2}}\right),$
(A.19) $\displaystyle\quad g^{Z^{\prime}}_{\mathbf{e}_{R}}=$ $\displaystyle\
-{{g^{\prime}}}\begin{pmatrix}s_{R}^{2}&-c_{R}s_{R}\delta_{R}&-c_{R}s_{R}\\\
-c_{R}s_{R}\delta_{R}&1&-s_{R}^{2}\delta_{R}\\\
-c_{R}s_{R}&-s_{R}^{2}\delta_{R}&c_{R}^{2}\end{pmatrix}+\mathcal{O}\left({\delta_{R}^{2}}\right).$
(A.20)

Hence, the effective couplings defined in Eq. (2.1) are given by

$\displaystyle\begin{pmatrix}g_{\mu\mu}^{L}&g_{\mu E}^{L}\\\ g_{\mu
E}^{L}&g_{EE}^{L}\end{pmatrix}\sim-{g^{\prime}}\begin{pmatrix}s_{L}^{2}&-s_{L}c_{L}\\\
-s_{L}c_{L}&c_{L}^{2}\end{pmatrix},\quad\begin{pmatrix}g_{\mu\mu}^{R}&g_{\mu
E}^{R}\\\ g_{\mu
E}^{R}&g_{EE}^{R}\end{pmatrix}\sim-{g^{\prime}}\begin{pmatrix}s_{R}^{2}&-c_{R}s_{R}\delta_{R}\\\
-c_{R}s_{R}\delta_{R}&1\end{pmatrix},$ (A.21)

in the doublet-like case, and these are given by

$\displaystyle\begin{pmatrix}g_{\mu\mu}^{L}&g_{\mu E}^{L}\\\ g_{\mu
E}^{L}&g_{EE}^{L}\end{pmatrix}\sim-{g^{\prime}}\begin{pmatrix}s_{L}^{2}&c_{L}s_{L}\delta_{L}\\\
c_{L}s_{L}\delta_{L}&1\end{pmatrix},\quad\begin{pmatrix}g_{\mu\mu}^{R}&g_{\mu
E}^{R}\\\ g_{\mu
E}^{R}&g_{EE}^{R}\end{pmatrix}\sim-{g^{\prime}}\begin{pmatrix}s_{R}^{2}&-c_{R}s_{R}\\\
-c_{R}s_{R}&c_{R}^{2}\end{pmatrix},$ (A.22)

in the singlet-like case.

The Yukawa couplings with $\chi$ are given by

$\displaystyle
Y^{\chi}_{\mathbf{e}}\sim\begin{pmatrix}0&\lambda_{E}c_{R}\delta_{L}&\lambda_{E}c_{R}\\\
\lambda_{L}c_{L}&\lambda_{L}s_{L}&\lambda_{E}s_{R}\delta_{R}-\lambda_{L}s_{L}\delta_{L}\\\
-\lambda_{L}c_{L}\delta_{R}&\lambda_{E}s_{R}\delta_{L}-\lambda_{L}s_{L}\delta_{R}&\lambda_{E}s_{R}\end{pmatrix},\quad
Y^{\chi}_{\mathbf{n}}\sim\lambda_{L}\begin{pmatrix}0&0\\\
c_{L}&s_{L}\end{pmatrix},$ (A.23)

where the coupling of $\chi$ with $\mu\mu$ is as small as $m_{\mu}/m_{E}$. The
couplings to the SM bosons are the SM-like up to
$\mathcal{O}\left({m_{\mu}/m_{E}}\right)$.

###  A.2 Muon anomalies

The ${Z^{\prime}}$ and $\chi$ boson contribution to $\Delta a_{\mu}$ is given
by [90, 33]

$\displaystyle\Delta a_{\mu}\sim\frac{m_{\mu}\kappa
v_{H}}{64\pi^{2}v_{\Phi}^{2}}s_{2L}s_{2R}C_{LR},$ (A.24)

with

$\displaystyle
C_{LR}:=\sqrt{x_{L}x_{E}}\frac{G_{Z}(x_{L})-G_{Z}(x_{E})}{x_{L}-x_{E}}+\frac{1}{2}\sqrt{y_{L}y_{R}}\frac{{y_{L}}G_{{S}}(y_{L})-{y_{R}}G_{{S}}(y_{R})}{y_{L}-y_{R}},$
(A.25)

where $x_{L}:=M_{L}^{2}/m_{Z^{\prime}}^{2}$,
$x_{E}:=M_{E}^{2}/m_{Z^{\prime}}^{2}$, $y_{L}:=M_{L}^{2}/m_{\chi}^{2}$ and
$y_{E}:=M_{E}^{2}/m_{\chi}^{2}$. Here,
$m_{Z^{\prime}}^{2}=2{g^{\prime}}^{2}v_{\Phi}^{2}$ is used. The loop functions
are given by

$\displaystyle G_{Z}(x):=$ $\displaystyle\
\frac{x^{3}+3x-6x\ln{(x)}-4}{2(1-x)^{3}},\quad
G_{S}(y):=\frac{y^{2}-4y+2\ln{(y)}+3}{(1-y)^{3}}.$ (A.26)

The contribution from the scalar $\chi$ is included since it is sizable unless
$m_{\chi}$ is very heavy which requires very large quartic couplings.

For the $b\to s\mu\mu$ anomaly, the Wilson coefficients are given by

$\displaystyle C_{9}\sim$ $\displaystyle\
-\frac{\sqrt{2}}{4G_{F}}\frac{4\pi}{\alpha_{e}}\frac{1}{V_{tb}V^{*}_{ts}}\frac{1}{4v_{\Phi}^{2}}(s_{R}^{2}+s_{L}^{2})\epsilon_{Q_{e}}\epsilon_{Q_{3}},$
(A.27) $\displaystyle C_{10}\sim$ $\displaystyle\
-\frac{\sqrt{2}}{4G_{F}}\frac{4\pi}{\alpha_{e}}\frac{1}{V_{tb}V^{*}_{ts}}\frac{1}{4v_{\Phi}^{2}}(s_{R}^{2}-s_{L}^{2})\epsilon_{Q_{e}}\epsilon_{Q_{3}},$
(A.28)

where the ${Z^{\prime}}$ boson couplings to the SM doublet quarks are
parametrized as

$\displaystyle\left[g^{{Z^{\prime}}}_{d_{L}}\right]_{ij}\sim\left[g^{{Z^{\prime}}}_{u_{L}}\right]_{ij}\sim-{g^{\prime}}\epsilon_{Q_{i}}\epsilon_{Q_{j}}.$
(A.29)

$\epsilon_{Q_{i}}$ is the similar quantity as
$s_{L}:=\lambda_{L}v_{\Phi}/m_{L}$, but we now consider the couplings with the
second and third generation quarks and these are typically small in contrast
to that for muon.

From Eq. (A.24),

$\displaystyle\Delta a_{\mu}\sim 2.9\times
10^{-9}\times\left(\frac{1.0~{}\mathrm{TeV}}{v_{\Phi}}\right)^{2}\left(\frac{\kappa}{1.0}\right)\left(\frac{s_{2L}s_{2R}}{1.0}\right)\left(\frac{C_{LR}}{0.1}\right).$
(A.30)

For the $b\to s\mu\mu$ anomaly,

$\displaystyle
C_{9}\sim-0.62\times\left(\frac{1.0~{}\mathrm{TeV}}{v_{\Phi}}\right)^{2}\left(\frac{s_{L}^{2}+s_{R}^{2}}{1}\right)\left(\frac{\epsilon_{Q_{2}}\epsilon_{Q_{3}}}{-0.002}\right).$
(A.31)

Assuming $s_{L}=s_{R}=1/\sqrt{2}$, i.e. $\lambda_{L}v_{\Phi}=m_{L}$ and
$\lambda_{E}v_{\Phi}=m_{E}$, the quark mixing angles are given by

$\displaystyle\epsilon_{Q_{2}}\epsilon_{Q_{3}}\sim-0.003\times\left(\frac{C_{9}}{-0.82}\right)\left(\frac{2.51\times
10^{-9}}{\Delta
a_{\mu}}\right)\left(\frac{\kappa}{1.0}\right)\left(\frac{C_{LR}}{0.1}\right),$
(A.32)

when the both anomalies are explained. With such small couplings with quarks,
${Z^{\prime}}$ boson is sufficiently suppressed to be consistent with the
constraints from the resonant di-lepton signal search at the LHC, unless
$\epsilon_{Q_{2}}\sim 1$ or $\epsilon_{Q_{3}}\sim 1$ to have large production
cross section from $s\overline{s}$ or $b\overline{b}$, respectively.

###  A.3 Benchmark

We show a benchmark scenario which explains the anomalies in $\Delta a_{\mu}$
and $b\to s\ell\ell$ and the excess in SR5L simultaneously. As discussed in
the main text, singlet-like VL lepton is more suitable to explain the excess
in SR5L. We take 111111 With these values, $v_{\Phi}=1103~{}\mathrm{GeV}$
which is sufficiently large to evade the bound from the neutrino trident
process [91, 92, 93, 94, 95, 96].

$\displaystyle m_{Z^{\prime}}=$ $\displaystyle\ 390~{}\mathrm{GeV},\quad
M_{L}=1.1~{}\mathrm{TeV},\quad M_{E}=404~{}\mathrm{GeV},\quad
m_{\chi}=365~{}\mathrm{GeV},$ (A.33) $\displaystyle s_{L}=$ $\displaystyle\
s_{R}=1/\sqrt{2},\quad
y_{\mu}v_{H}=2m_{\mu},\quad{\kappa^{\prime}}=0,\quad\kappa=-0.821,\quad{g^{\prime}}=0.25.$

The VL lepton masses are $400$ and $1111~{}\mathrm{GeV}$. The correction to
the anomalous magnetic moment of the muon is $\Delta a_{\mu}=2.51\times
10^{-9}$. $s_{L}=s_{R}$ realizes the $C_{9}$-only scenario, and
$C_{9}\sim-0.81$ is explained if
$\epsilon_{Q_{2}}\epsilon_{Q_{3}}\sim-0.0032$.

The partial decay widths of the singlet VL lepton $E_{1}$ are approximately
given by

$\displaystyle\Gamma\left(E_{1}\to{Z^{\prime}}\mu\right)\sim$ $\displaystyle\
\frac{M_{E}^{3}}{64\pi v_{\Phi}^{2}}c_{R}^{2}s_{R}^{2}(1-z)^{2}(1+2z),$ (A.34)
$\displaystyle\Gamma\left(E_{1}\to\chi\mu\right)\sim$ $\displaystyle\
\frac{M_{E}^{3}}{64\pi v_{\Phi}^{2}}c_{R}^{2}s_{R}^{2}(1-x)^{2},$ (A.35)

where $z:=m_{Z^{\prime}}^{2}/M_{E}^{2}$ and $x:=m_{\chi}^{2}/M_{E}^{2}$. Hence
the branching fraction of $E_{1}$, assuming no other decay modes, is
approximately given by

$\displaystyle\mathrm{BR}\left({E_{1}}\to{{Z^{\prime}}\mu}\right)\sim\frac{(1-z)^{2}(1+2z)}{(1-z)^{2}(1+2z)+(1-x)^{2}}.$
(A.36)

At the benchmark point,
$\mathrm{BR}\left({E_{1}}\to{{Z^{\prime}}\mu}\right)\simeq 0.25$. The $\chi$
boson predominantly decays to VL fermions as far as these are kinematically
allowed. If this is not the case, it should decay to a pair of SM leptons or
quarks. For the lepton coupling, as seen from Eq. (A.23), the coupling to 2
muons are strongly suppressed by the muon mass. Hence, the dominant decay mode
of $\chi$ may be to a pair of top quarks due to the weaker suppression if
$m_{\chi}>2m_{t}$ which is true at the benchmark point. In this case, the
processes with $\chi$ decays will not contribute to the
$\texttt{SR0}^{\mathrm{tight}}_{\mathrm{bveto}}$ and
$\texttt{SR0}^{\mathrm{loose}}_{\mathrm{bveto}}$ due to the b-jet veto, and
thus the results in the main text will not be changed. If $m_{\chi}<2m_{t}$,
$\chi$ decays to a pair of bottom quarks, where the relevant Yukawa is
estimated as $\sim\epsilon_{3}m_{b}/M_{Q}\sim 10^{-4}$ for $\epsilon_{3}\sim
0.1$ and $M_{Q}\sim 4~{}\mathrm{TeV}$. This would be comparable to the decay
to a pair of muons which the relevant Yukawa coupling is estimated as $\sim
m_{\mu}/M_{L}\sim 10^{-4}$ for $M_{L}\sim 1~{}\mathrm{TeV}$. Thus, there will
be additional contributions with six muons on top of the decays from
${Z^{\prime}}$ boson.

## References

  * [1] ATLAS Collaboration, G. Aad et al., Observation of a new particle in the search for the Standard Model Higgs boson with the ATLAS detector at the LHC, Phys. Lett. B 716 (2012) 1–29, [arXiv:1207.7214].
  * [2] CMS Collaboration, S. Chatrchyan et al., Observation of a New Boson at a Mass of 125 GeV with the CMS Experiment at the LHC, Phys. Lett. B 716 (2012) 30–61, [arXiv:1207.7235].
  * [3] LHCb Collaboration, R. Aaij et al., Test of lepton universality using $B^{+}\rightarrow K^{+}\ell^{+}\ell^{-}$ decays, Phys. Rev. Lett. 113 (2014) 151601, [arXiv:1406.6482].
  * [4] 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].
  * [5] 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].
  * [6] Belle Collaboration, A. Abdesselam et al., Test of lepton flavor universality in ${B\to K^{\ast}\ell^{+}\ell^{-}}$ decays at Belle, arXiv:1904.02440.
  * [7] LHCb Collaboration, R. Aaij et al., Differential branching fraction and angular analysis of the decay $B_{s}^{0}\to\phi\mu^{+}\mu^{-}$, JHEP 07 (2013) 084, [arXiv:1305.2168].
  * [8] BaBar Collaboration, J. P. Lees et al., Measurement of the $B\to X_{s}l^{+}l^{-}$ branching fraction and search for direct CP violation from a sum of exclusive final states, Phys. Rev. Lett. 112 (2014) 211802, [arXiv:1312.5364].
  * [9] 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].
  * [10] 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].
  * [11] LHCb Collaboration, R. Aaij et al., Measurement of Form-Factor-Independent Observables in the Decay $B^{0}\to K^{*0}\mu^{+}\mu^{-}$, Phys. Rev. Lett. 111 (2013) 191801, [arXiv:1308.1707].
  * [12] 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. B753 (2016) 424–448, [arXiv:1507.08126].
  * [13] LHCb Collaboration, R. Aaij et al., Angular analysis of the $B^{0}\to K^{*0}\mu^{+}\mu^{-}$ decay using 3 fb-1 of integrated luminosity, JHEP 02 (2016) 104, [arXiv:1512.04442].
  * [14] Belle Collaboration, A. Abdesselam et al., Angular analysis of $B^{0}\to K^{\ast}(892)^{0}\ell^{+}\ell^{-}$, in Proceedings, LHCSki 2016 - A First Discussion of 13 TeV Results: Obergurgl, Austria, April 10-15, 2016, 2016\.  arXiv:1604.04042.
  * [15] 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].
  * [16] ATLAS Collaboration, T. A. collaboration, Angular analysis of $B^{0}_{d}\to K^{*}\mu^{+}\mu^{-}$ decays in $pp$ collisions at $\sqrt{s}=8$ TeV with the ATLAS detector, .
  * [17] CMS Collaboration, C. Collaboration, Measurement of the $P_{1}$ and $P_{5}^{\prime}$ angular parameters of the decay $\mathrm{B}^{0}\to\mathrm{K}^{*0}\mu^{+}\mu^{-}$ in proton-proton collisions at $\sqrt{s}=8~{}\mathrm{TeV}$, .
  * [18] LHCb Collaboration, R. Aaij et al., Angular analysis of the $B^{+}\rightarrow K^{\ast+}\mu^{+}\mu^{-}$ decay, arXiv:2012.13241.
  * [19] LHCb Collaboration, R. Aaij et al., Test of lepton universality in beauty-quark decays, arXiv:2103.11769.
  * [20] Muon g-2 Collaboration, G. W. Bennett et al., Final Report of the Muon E821 Anomalous Magnetic Moment Measurement at BNL, Phys. Rev. D73 (2006) 072003, [hep-ex/0602035].
  * [21] T. Aoyama et al., The anomalous magnetic moment of the muon in the Standard Model, Phys. Rept. 887 (2020) 1–166, [arXiv:2006.04822].
  * [22] Muon g-2 Collaboration, B. Abi et al., Measurement of the Positive Muon Anomalous Magnetic Moment to 0.46 ppm, Phys. Rev. Lett. 126 (2021) 141801, [arXiv:2104.03281].
  * [23] J. Kawamura, S. Raby, and A. Trautner, Complete vectorlike fourth family and new U(1)’ for muon anomalies, Phys. Rev. D 100 (2019), no. 5 055030, [arXiv:1906.11297].
  * [24] J. Kawamura, S. Raby, and A. Trautner, Complete vectorlike fourth family with U(1)’ : A global analysis, Phys. Rev. D 101 (2020), no. 3 035026, [arXiv:1911.11075].
  * [25] B. Allanach, F. S. Queiroz, A. Strumia, and S. Sun, $Z^{\prime}$ models for the LHCb and $g-2$ muon anomalies, Phys. Rev. D93 (2016), no. 5 055045, [arXiv:1511.07447]. [Erratum: Phys. Rev.D95,no.11,119902(2017)].
  * [26] W. Altmannshofer, M. Carena, and A. Crivellin, $L_{\mu}-L_{\tau}$ theory of Higgs flavor violation and $(g-2)_{\mu}$, Phys. Rev. D94 (2016), no. 9 095026, [arXiv:1604.08221].
  * [27] E. Megias, M. Quiros, and L. Salas, $g_{\mu}-2$ from Vector-Like Leptons in Warped Space, JHEP 05 (2017) 016, [arXiv:1701.05072].
  * [28] S. Raby and A. Trautner, Vectorlike chiral fourth family to explain muon anomalies, Phys. Rev. D97 (2018), no. 9 095006, [arXiv:1712.09360].
  * [29] S. F. King, Flavourful $Z^{\prime}$ models for ${R}_{K^{\left(\ast\right)}}$, JHEP 08 (2017) 019, [arXiv:1706.06100].
  * [30] L. Darmé, K. Kowalska, L. Roszkowski, and E. M. Sessolo, Flavor anomalies and dark matter in SUSY with an extra U(1), JHEP 10 (2018) 052, [arXiv:1806.06036].
  * [31] A. Czarnecki and W. J. Marciano, The Muon anomalous magnetic moment: A Harbinger for ’new physics’, Phys. Rev. D64 (2001) 013014, [hep-ph/0102122].
  * [32] K. Kannike, M. Raidal, D. M. Straub, and A. Strumia, Anthropic solution to the magnetic muon anomaly: the charged see-saw, JHEP 02 (2012) 106, [arXiv:1111.2551]. [Erratum: JHEP10,136(2012)].
  * [33] R. Dermisek and A. Raval, Explanation of the Muon g-2 Anomaly with Vectorlike Leptons and its Implications for Higgs Decays, Phys. Rev. D88 (2013) 013017, [arXiv:1305.3522].
  * [34] Z. Poh and S. Raby, Vectorlike leptons: Muon g-2 anomaly, lepton flavor violation, Higgs boson decays, and lepton nonuniversality, Phys. Rev. D96 (2017), no. 1 015032, [arXiv:1705.07007].
  * [35] J. Kawamura, S. Okawa, and Y. Omura, Current status and muon $g-2$ explanation of lepton portal dark matter, JHEP 08 (2020) 042, [arXiv:2002.12534].
  * [36] Y. Bai and J. Berger, Muon g-2 in Lepton Portal Dark Matter, arXiv:2104.03301.
  * [37] B. Gripaios, M. Nardecchia, and S. A. Renner, Linear flavour violation and anomalies in B physics, JHEP 06 (2016) 083, [arXiv:1509.05020].
  * [38] P. Arnan, L. Hofer, F. Mescia, and A. Crivellin, Loop effects of heavy new scalars and fermions in $b\to s\mu^{+}\mu^{-}$, JHEP 04 (2017) 043, [arXiv:1608.07832].
  * [39] B. Grinstein, S. Pokorski, and G. G. Ross, Lepton non-universality in $B$ decays and fermion mass structure, JHEP 12 (2018) 079, [arXiv:1809.01766].
  * [40] P. Arnan, A. Crivellin, M. Fedele, and F. Mescia, Generic Loop Effects of New Scalars and Fermions in $b\to s\ell^{+}\ell^{-}$ and a Vector-like $4^{\rm th}$ Generation, arXiv:1904.05890.
  * [41] C.-W. Chiang and H. Okada, A simple model for explaining muon-related anomalies and dark matter, arXiv:1711.07365.
  * [42] J. M. Cline and J. M. Cornell, $R({K^{(*)}})$ from dark matter exchange, Phys. Lett. B782 (2018) 232–237, [arXiv:1711.10770].
  * [43] J. Kawamura, S. Okawa, and Y. Omura, Interplay between the b$\to s\ell\ell$ anomalies and dark matter physics, Phys. Rev. D96 (2017), no. 7 075041, [arXiv:1706.04344].
  * [44] D. G. Cerdeño, A. Cheek, P. Martín-Ramiro, and J. M. Moreno, B anomalies and dark matter: a complex connection, Eur. Phys. J. C 79 (2019), no. 6 517, [arXiv:1902.01789].
  * [45] G. Arcadi, L. Calibbi, M. Fedele, and F. Mescia, Muon $g-2$ and $B$-anomalies from Dark Matter, arXiv:2104.03228.
  * [46] ATLAS Collaboration, G. Aad et al., Search for supersymmetry in events with four or more charged leptons in $139\,\mbox{fb${}^{-1}$}$ of $\sqrt{s}=13$ TeV $pp$ collisions with the ATLAS detector, arXiv:2103.11684.
  * [47] A. J. Buras and M. Munz, Effective Hamiltonian for $B\to X(s)e^{+}e^{-}$ beyond leading logarithms in the NDR and HV schemes, Phys. Rev. D52 (1995) 186–195, [hep-ph/9501281].
  * [48] C. Bobeth, M. Misiak, and J. Urban, Photonic penguins at two loops and $m_{t}$ dependence of $BR[B\to X_{s}l^{+}l^{-}]$, Nucl. Phys. B574 (2000) 291–330, [hep-ph/9910220].
  * [49] ATLAS Collaboration, G. Aad et al., Search for high-mass dilepton resonances using 139 fb-1 of $pp$ collision data collected at $\sqrt{s}=$13 TeV with the ATLAS detector, Phys. Lett. B 796 (2019) 68–87, [arXiv:1903.06248].
  * [50] M. Kohda, T. Modak, and A. Soffer, Identifying a $Z^{\prime}$ behind $b\to s\ell\ell$ anomalies at the LHC, Phys. Rev. D97 (2018), no. 11 115019, [arXiv:1803.07492].
  * [51] B. C. Allanach, J. M. Butterworth, and T. Corbett, Collider Constraints on $Z^{\prime}$ Models for Neutral Current $B-$Anomalies, arXiv:1904.10954.
  * [52] L.-S. Geng, B. Grinstein, S. Jäger, S.-Y. Li, J. Martin Camalich, and R.-X. Shi, Implications of new evidence for lepton-universality violation in $b\to s\ell^{+}\ell^{-}$ decays, arXiv:2103.12738.
  * [53] W. Altmannshofer and P. Stangl, New Physics in Rare B Decays after Moriond 2021, arXiv:2103.13370.
  * [54] W. Altmannshofer, C. Niehoff, P. Stangl, and D. M. Straub, Status of the $B\rightarrow K^{*}\mu^{+}\mu^{-}$ anomaly after Moriond 2017, Eur. Phys. J. C77 (2017), no. 6 377, [arXiv:1703.09189].
  * [55] W. Altmannshofer, P. Stangl, and D. M. Straub, Interpreting Hints for Lepton Flavor Universality Violation, Phys. Rev. D96 (2017), no. 5 055008, [arXiv:1704.05435].
  * [56] A. K. Alok, B. Bhattacharya, A. Datta, D. Kumar, J. Kumar, and D. London, New Physics in $b\to s\mu^{+}\mu^{-}$ after the Measurement of $R_{K^{*}}$, Phys. Rev. D96 (2017), no. 9 095009, [arXiv:1704.07397].
  * [57] B. Capdevila, A. Crivellin, S. Descotes-Genon, J. Matias, and J. Virto, Patterns of New Physics in $b\to s\ell^{+}\ell^{-}$ transitions in the light of recent data, JHEP 01 (2018) 093, [arXiv:1704.05340].
  * [58] M. Ciuchini, A. M. Coutinho, M. Fedele, E. Franco, A. Paul, L. Silvestrini, and M. Valli, On Flavourful Easter eggs for New Physics hunger and Lepton Flavour Universality violation, Eur. Phys. J. C77 (2017), no. 10 688, [arXiv:1704.05447].
  * [59] G. D’Amico, M. Nardecchia, P. Panci, F. Sannino, A. Strumia, R. Torre, and A. Urbano, Flavour anomalies after the $R_{K^{*}}$ measurement, JHEP 09 (2017) 010, [arXiv:1704.05438].
  * [60] L.-S. Geng, B. Grinstein, S. Jäger, J. Martin Camalich, X.-L. Ren, and R.-X. Shi, Towards the discovery of new physics with lepton-universality ratios of $b\to s\ell\ell$ decays, Phys. Rev. D96 (2017), no. 9 093006, [arXiv:1704.05446].
  * [61] D. Ghosh, Explaining the $R_{K}$ and $R_{K^{*}}$ anomalies, Eur. Phys. J. C77 (2017), no. 10 694, [arXiv:1704.06240].
  * [62] A. Arbey, T. Hurth, F. Mahmoudi, and S. Neshatpour, Hadronic and New Physics Contributions to $b\to s$ Transitions, Phys. Rev. D98 (2018), no. 9 095027, [arXiv:1806.02791].
  * [63] J. Aebischer, W. Altmannshofer, D. Guadagnoli, M. Reboud, P. Stangl, and D. M. Straub, $B$-decay discrepancies after Moriond 2019, arXiv:1903.10434.
  * [64] M. Algueró, B. Capdevila, A. Crivellin, S. Descotes-Genon, P. Masjuan, J. Matias, M. Novoa Brunet, and J. Virto, Emerging patterns of New Physics with and without Lepton Flavour Universal contributions, Eur. Phys. J. C 79 (2019), no. 8 714, [arXiv:1903.09578]. [Addendum: Eur.Phys.J.C 80, 511 (2020)].
  * [65] A. K. Alok, A. Dighe, S. Gangal, and D. Kumar, Continuing search for new physics in $b\to s\mu\mu$ decays: two operators at a time, arXiv:1903.09617.
  * [66] M. Ciuchini, A. M. Coutinho, M. Fedele, E. Franco, A. Paul, L. Silvestrini, and M. Valli, New Physics in $b\to s\ell^{+}\ell^{-}$ confronts new data on Lepton Universality, arXiv:1903.09632.
  * [67] A. Datta, J. Kumar, and D. London, The $B$ Anomalies and New Physics in $b\to se^{+}e^{-}$, arXiv:1903.10086.
  * [68] K. Kowalska, D. Kumar, and E. M. Sessolo, Implications for New Physics in $b\to s\mu\mu$ transitions after recent measurements by Belle and LHCb, arXiv:1903.10932.
  * [69] A. Arbey, T. Hurth, F. Mahmoudi, D. Martinez Santos, and S. Neshatpour, Update on the b$\to$s anomalies, arXiv:1904.08399.
  * [70] D. Kumar, K. Kowalska, and E. M. Sessolo, Global Bayesian Analysis of new physics in $b\to s\mu\mu$ transitions after Moriond-2019, 2019.  arXiv:1906.08596.
  * [71] ATLAS Collaboration, G. Aad et al., Search for heavy lepton resonances decaying to a $Z$ boson and a lepton in $pp$ collisions at $\sqrt{s}=8$ TeV with the ATLAS detector, JHEP 09 (2015) 108, [arXiv:1506.01291].
  * [72] R. Dermisek, J. P. Hall, E. Lunghi, and S. Shin, Limits on Vectorlike Leptons from Searches for Anomalous Production of Multi-Lepton Events, JHEP 12 (2014) 013, [arXiv:1408.3123].
  * [73] CMS Collaboration, A. M. Sirunyan et al., Search for vector-like leptons in multilepton final states in proton-proton collisions at $\sqrt{s}$ = 13 TeV, arXiv:1905.10853.
  * [74] N. Kumar and S. P. Martin, Vectorlike Leptons at the Large Hadron Collider, Phys. Rev. D92 (2015), no. 11 115018, [arXiv:1510.03456].
  * [75] A. Falkowski, D. M. Straub, and A. Vicente, Vector-like leptons: Higgs decays and collider phenomenology, JHEP 05 (2014) 092, [arXiv:1312.5329].
  * [76] S. A. R. Ellis, R. M. Godbole, S. Gopalakrishna, and J. D. Wells, Survey of vector-like fermion extensions of the Standard Model and their phenomenological implications, JHEP 09 (2014) 130, [arXiv:1404.4398].
  * [77] P. N. Bhattiprolu and S. P. Martin, Prospects for vectorlike leptons at future proton-proton colliders, Phys. Rev. D 100 (2019), no. 1 015033, [arXiv:1905.00498].
  * [78] J. Ellis, TikZ-Feynman: Feynman diagrams with TikZ, Comput. Phys. Commun. 210 (2017) 103–123, [arXiv:1601.05437].
  * [79] M. Dohse, TikZ-FeynHand: Basic User Guide, arXiv:1802.00689.
  * [80] ATLAS Collaboration, G. Aad et al., Search for electroweak production of charginos and sleptons decaying into final states with two leptons and missing transverse momentum in $\sqrt{s}=13$ TeV $pp$ collisions using the ATLAS detector, Eur. Phys. J. C 80 (2020), no. 2 123, [arXiv:1908.08215].
  * [81] J. Alwall, R. Frederix, S. Frixione, V. Hirschi, F. Maltoni, O. Mattelaer, H. S. Shao, T. Stelzer, P. Torrielli, and M. Zaro, The automated computation of tree-level and next-to-leading order differential cross sections, and their matching to parton shower simulations, JHEP 07 (2014) 079, [arXiv:1405.0301].
  * [82] C. Degrande, C. Duhr, B. Fuks, D. Grellscheid, O. Mattelaer, and T. Reiter, UFO - The Universal FeynRules Output, Comput. Phys. Commun. 183 (2012) 1201–1214, [arXiv:1108.2040].
  * [83] N. D. Christensen and C. Duhr, FeynRules - Feynman rules made easy, Comput. Phys. Commun. 180 (2009) 1614–1641, [arXiv:0806.4194].
  * [84] DELPHES 3 Collaboration, J. de Favereau, C. Delaere, P. Demin, A. Giammanco, V. Lemaître, A. Mertens, and M. Selvaggi, DELPHES 3, A modular framework for fast simulation of a generic collider experiment, JHEP 02 (2014) 057, [arXiv:1307.6346].
  * [85] C. G. Lester and D. J. Summers, Measuring masses of semiinvisibly decaying particles pair produced at hadron colliders, Phys. Lett. B 463 (1999) 99–103, [hep-ph/9906349].
  * [86] A. Barr, C. Lester, and P. Stephens, m(T2): The Truth behind the glamour, J. Phys. G 29 (2003) 2343–2363, [hep-ph/0304226].
  * [87] C. G. Lester and B. Nachman, Bisection-based asymmetric MT2 computation: a higher precision calculator than existing symmetric methods, JHEP 03 (2015) 100, [arXiv:1411.4312].
  * [88] CMS Collaboration, A. M. Sirunyan et al., Search for resonant and nonresonant new phenomena in high-mass dilepton final states at $\sqrt{s}=$ 13 TeV, arXiv:2103.02708.
  * [89] P. N. Bhattiprolu, S. P. Martin, and J. D. Wells, Criteria for projected discovery and exclusion sensitivities of counting experiments, Eur. Phys. J. C 81 (2021), no. 2 123, [arXiv:2009.07249].
  * [90] F. Jegerlehner and A. Nyffeler, The Muon g-2, Phys. Rept. 477 (2009) 1–110, [arXiv:0902.3360].
  * [91] W. Altmannshofer, S. Gori, M. Pospelov, and I. Yavin, Quark flavor transitions in $L_{\mu}-L_{\tau}$ models, Phys. Rev. D89 (2014) 095033, [arXiv:1403.1269].
  * [92] 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].
  * [93] G. Magill and R. Plestid, Neutrino Trident Production at the Intensity Frontier, Phys. Rev. D95 (2017), no. 7 073004, [arXiv:1612.05642].
  * [94] S.-F. Ge, M. Lindner, and W. Rodejohann, Atmospheric Trident Production for Probing New Physics, Phys. Lett. B772 (2017) 164–168, [arXiv:1702.02617].
  * [95] P. Ballett, M. Hostert, S. Pascoli, Y. F. Perez-Gonzalez, Z. Tabrizi, and R. Zukanovich Funchal, Neutrino Trident Scattering at Near Detectors, JHEP 01 (2019) 119, [arXiv:1807.10973].
  * [96] W. Altmannshofer, S. Gori, J. Martín-Albo, A. Sousa, and M. Wallbank, Neutrino Tridents at DUNE, arXiv:1902.06765.