# Ultra-coherent nanomechanical resonators based on density phononic crystal
engineering

Dennis Høj<EMAIL_ADDRESS>Center for Macroscopic Quantum States (bigQ),
Department of Physics, Technical University of Denmark, Fysikvej, 2800 Kgs.
Lyngby, Denmark  Ulrich Busk Hoff<EMAIL_ADDRESS>Center for
Macroscopic Quantum States (bigQ), Department of Physics, Technical University
of Denmark, Fysikvej, 2800 Kgs. Lyngby, Denmark  Ulrik Lund Andersen
<EMAIL_ADDRESS>Center for Macroscopic Quantum States (bigQ),
Department of Physics, Technical University of Denmark, Fysikvej, 2800 Kgs.
Lyngby, Denmark

###### Abstract

Micro- and nanomechanical systems with exceptionally low dissipation rates are
enabling the next-generation technologies of ultra-sensitive detectors and
quantum information systems. New techniques and methods for lowering the
dissipation rate have in recent years been discovered and allowed for the
engineering of mechanical oscillators with phononic modes that are extremely
well isolated from the environment and thus possessing quality factors close
to and beyond 1 billion. A powerful strategy for isolating and controlling a
single phononic mode is based on phononic crystal engineering. Here we propose
a new method for phononic crystal engineering of nanomechanical oscillators
that is based on a periodic variation of the material density. To circumvent
the introduction of additional bending losses resulting from the variation of
material density, the added mass constitutes an array of nanopillars in which
the losses will be diluted. Using this novel technique for phononic crystal
engineering, we design and fabricate corrugated mechanical oscillators with
quality factors approaching one billion in a room temperature environment. The
flexibility space of these new phononic crystals is large and further
advancement can be attained through optimized phononic crystal patterning and
strain engineering via topology optimization. This will allow for the
engineering of mechanical membranes with quality factors approaching 10
billion. Such extremely low mechanical dissipation rates will enable the
development of radically new technologies such as quantum-limited atomic force
microscopy at room-temperature, ultra-sensitive detectors of dark matter,
spontaneous waveform collapses, gravity, and high-efficiency quantum
information transducers.

Phononic crystal engineering is an extremely powerful method for controlling
the phononic properties of material systems in much the same way as photonic
crystals are able to control light [1, 2]. The method is being used to
engineer material systems for controlling heat and sound both at the macro-
and the micro-scale, for example with applications in the design of
loudspeakers and microphones to control the direction of sound [3], in
buildings and vehicles to isolate sound, and in communication technologies to
enable signal processing at 5G frequencies. In particular, in recent years,
phononic crystals have been used to control and isolate specific phononic
modes of nano- and micromechanical resonators, protecting them from the
surrounding thermal heat bath [4, 5, 6, 7, 8, 9]. They greatly outperform
other resonator designs [10, 11, 12, 13] and have given rise to a completely
new generation of ultra-coherent mechanical oscillators with quality factors
reaching one billion at room temperature.

The effect of enhancing the quality factor of mechanical resonators using
phononic crystals is two-fold. First, as a result of the periodic phononic
structure, an acoustic bandgap is formed in which phononic modes are trapped
and thus unable to radiate into the surrounding environment that otherwise
would cause dissipation [4, 5, 14]. Second, the phononic crystal structure can
be also engineered to reduce the phononic mode curvature of the resonator at
the clamping points to the substrate [15, 16], which will lead to a decrease
in the intrinsic dissipation rate for highly stressed systems – a technique
that is known as soft-clamping [7, 8, 9, 17, 18, 19]. Phononic crystal
membranes and strings exploiting these effects in harmony have been
fabricated, and the membranes have been used in a series of optomechanical
experiments [20, 21, 22].

A phononic crystal in a mechanical membrane is traditionally engineered by
mapping a periodic pattern of holes into a pre-stressed silicon membrane. This
produces a periodic variation in the stress of the material across the
membrane, which in turn creates a phononic bandgap with isolated modes that
are associated with the vibrations at a central ”defect” of the pattern. As
mentioned above, this leads to acoustic isolation as well as soft-clamping of
a particular mode, and thus very large quality factors. However, the stress of
the material at the location of the high-quality mode is low which means that
the quality-enhancing effect of dissipation dilution[23, 10, 24] is not
optimally used, and therefore eventually bounds the quality factor of the
system. The technique of strain engineering can in some geometries be applied
to co-locate a high stress region with the mechanical motion [8] but in this
particular phononic membrane, it seems not possible.

In this Letter, we propose a new type of phononic crystal engineering of a
nano- and micromechanical resonator that offers a larger flexibility in the
design of periodic patterns and in which the stress can be engineered. Our
strategy for phononic crystal nano-engineering is based on the creation of a
periodic variation in the mass density across the membrane. More specifically,
we propose and demonstrate the mapping of a periodic hexagonal pattern with a
central defect onto a membrane that consists of an array of nanopillars,
resulting in a membrane with periodic pattern of cell domains each containing
an array of nanopillars as illustrated in Fig.1a-c (note that the central
defect is not shown in the figure). Such a density phononic membrane enables
the formation of a frequency bandgap with modes that are highly shielded from
the substrate (acoustic filtering) and it reduces the curvature of the mode
amplitude towards the clamping points, thereby lowering the intrinsic
dissipation losses (soft-clamping). Moreover, the use of a fine mesh of
nanopillars for controlling the mass density ensures low intrinsic bending
dissipation at the pillars. We therefore expect the density phononic membranes
to exhibit ultra-low dissipation rates.

Intuitively, the simplest strategy for phononic crystal engineering via
material density variation is to periodically add mass by depositing large
pillars representing the cells of the phononic crystal pattern. Such a
structure will in principle create a bandgap (similar to a holey structure)
but the structure will be too stiff and it will be strongly affected by
boundary bending losses at each of the domain walls. However, this stiffness
and dissipation effect can be largely reduced by patterning a fine mesh of
nanopillars within each of the cells. With such a strategy, the density
contrast required for phononic engineering is kept while the bending loss is
significantly reduced. This reduction is caused by the dilution of the loss
into many nanopillars each of which produces almost negligible loss due to
their small diameter in comparison to the wavelength of the phononic mode. We
simulate the potential loss enhancement (see methods) by considering a
hexagonal pattern of pillars on a membrane with the effective density

$\rho_{\textrm{eff}}=\rho\left(1+\frac{\pi\sqrt{3}h_{\textrm{pil}}d_{\textrm{pil}}^{2}}{6ta_{\textrm{pil}}^{2}}\right)$
(1)

where $\rho$ is the density of the material (assuming the same material for
membrane and pillars), while $a_{\textrm{pil}}$ and $h_{\textrm{pil}}$ are the
diameter and height of the pillars, $d_{\textrm{pil}}$ is the distance between
them (see Fig. 1b) and $t$ is the thickness of the membrane. The results of
the simulations are presented in Fig. 1d where we plot the loss enhancement
$\alpha_{\gamma}$ of a membrane patterned with a hexagonal array of pillars
against the periodicity, $a_{\textrm{pil}}$, or pillar size,
$d_{\textrm{pil}}$ (assuming $a_{\textrm{pil}}=2d_{\textrm{pil}}$). It is
clear that the intrinsic loss of the membrane is dramatically reduced by the
introduction of a nanopillar array, and for very small pillars, the loss
contribution is effectively negligible.

Since the phononic crystal pattern is created by material density variation,
there exists a large flexibility in defining its form and shape. Two examples
are a binary and a sinusoidal pattern where the material density is either
changed abruptly or smoothly across the membrane following a certain
periodicity as shown in Fig. 2c. Analytical expressions for these two patterns
can be found in the supplementary material, and in the following we
interrogate the binary pattern while the results of the sinusoidal pattern are
discussed in the supplementary material. We assume the cell periodicity to be
given by $a_{\textrm{ph}}$, the diameter of each cell to be $d_{\textrm{ph}}$
and the ratio $\alpha_{w}=d_{\textrm{ph}}/a_{\textrm{ph}}$. Furthermore we
define the density contrast parameter, $g=\rho_{\textrm{eff}}/\rho$.

The optimization of the density phononic crystal is a trade-off between the
size of the bandgap, which should be as large as possible to encapsulate the
mode of interest, and the dissipation of an isolated mode, which should be as
low as possible. In Fig. 2a-b we present the results of simulating the
relative bandgap width and the quality factor of the isolated mode assuming
that intrinsic losses are dominating. It is clear that the bandgap becomes
wider as the density contrast, $g$, increases as expected, and that an optimum
is obtained for a specific relative distribution width, $\alpha_{w}$, of
around 0.4 (corresponding to $a_{\textrm{ph}}=2.5d_{\textrm{ph}}$). However,
contrarily, the quality factor for the isolated mode is decreasing as the
density contrast is increased. This is caused by the fact that a larger
material density will lead to a smaller phase velocity and therefore a shorter
wavelength of the propagating phononic mode. A shorter wavelength relative to
the fixed size of the nanopillars will naturally lead to larger bending losses
at the nanopillars, and thus a reduced quality factor. In the following, we
choose a contrast of $g=5$ which is a good compromise between having bandgap
of reasonable size and a large quality factor. In addition, based on these
simulations, in the following we choose $\alpha_{w}=0.4$.

To trap an isolated mode in the membrane, a defect is introduced simply by
removing one of the cell domains. Such a disturbance in the phononic crystal
allows some modes to occur in the bandgap which will be spatially well
confined at the defect. However, the resonance frequency of the high quality
mode is not centered in the bandgap but at the lower end of it which might
lead to some mode leakage. To avoid this effect, we increase the frequency of
the mode by lowering the material density of the phononic cells located just
next to the defect. By choosing the contrast of these inner cells to $g=3.9$
or the relative size to $\alpha_{w}=0.32$, the frequency will be centered in
the bandgap. In Fig. 2, we show the contour of the density contrast (Fig. 2d),
the displacement amplitude of the high quality mode (Fig. 2e) and the
dispersion relation that clearly exhibits the bandgap (Fig. 2f).
Interestingly, we note that it generates a full bandgap in contrast to stress-
induced phononic membranes which produces a pseudo-bandgap with in-plane modes
[7]. Finally, in Fig. 2c we present the results of the simulated quality
factor (based on intrinsic losses) for different number of layers of the
phononic crystal cells as well as three different values for the periodicity
of the phononic pattern ($a_{\textrm{ph}}={199,221,253}\ \mu m$). In these
simulations we assume no added losses from the pillars ($\alpha_{\gamma}=1$) .
The effect of soft-clamping is clear. As more layers are added, the mode field
curvatures at the boundary is gradually reduced which in turn leads to lower
intrinsic dissipation rates, thereby increasing the quality factor. We see
that at least 7 layers of phononic cells are required for saturating the
quality factor to a value above $10^{9}$.

Our next step is to experimentally validate the expected quality of our
proposed density phononic crystal. Towards this end, we fabricate membrane
structures for different pillar periodicity and for different numbers of
phononic cell layers. We fabricate the membranes from two thin layers of high-
stress (about 1GPa) silicon nitride films grown on a silicon wafer using low-
pressure chemical vapor deposition (CVD) and plasma-enhanced CVD,
respectively. Using UV lithography and reactive ion etching the nanopillar
pattern was transferred to the front-side while a window for releasing the
resonator was transferred to the back-side. More details on the fabrication
process can be found in the method section.

We characterize the fabricated membranes by determining their frequency power
spectrum as well as their quality factors using an interferometric ringdown
measurement setup operating in ultrahigh vacuum to avoid air damping (see
methods). For the ringdown measurements, we excite the membrane by the
exertion of radiation pressure force from a light beam modulated at the
frequency of the mode under investigation, turn off the modulation and
subsequently record the decay of the membrane’s motion using the same laser
mode. An example of the measured power spectrum is illustrated in Fig. 3c and
shows the clear existence of a bandgap around the confined mode at 1.4 MHz.
The result of the ringdown measurement of this particular mode is presented in
Fig. 3d and yields a quality factor of 5.65$\times$108.

In Fig. 3e we present the results for the measured quality factors of the
bandgap mode with different number of layers of phononic cells. As the number
of cell layers increases, the quality factor improves since the increasing
number of layers will gradually strengthen the effects of soft-clamping and
acoustic shielding against phonon tunneling loss. When the acoustic shielding
is weak (that is, when only a few cell layers are present), the coupling to
the surrounding substrate mode becomes stochastic resulting in a large spread
of the quality factors. For a large number of cell layers ($N_{ph}\geq 9$),
the acoustic shielding effect becomes more effective and the spread is low. In
Fig. 3e, we compare the measured data with theory that accounts for intrinsic
losses but neglecting the phonon tunneling losses, and we clearly observe
reasonable agreements to the best measured values.

Results for the quality factor of the bandgap mode for different nanopillar
periodicity, $a_{pil}$, are presented in Fig. 3f. As we kept the phononic
crystal pattern constant, a change in the periodicity effectively leads to a
change in the nanopillar diameter ($d_{pil}=a_{pil}/2$) and density. As
expected, the quality factor increases as we lower the periodicity (and thus
the size of the nanopillars) as this will dilute the effect of bending losses
at the pillars. We expect the quality factors to be limited by the intrinsic
bending losses since we have chosen a large number of cell layers, $N_{ph}=8$,
diminishing the phonon tunneling losses. Indeed, by using a model that
accounts only for intrinsic losses we obtain good agreement with the measured
data as shown in Fig. 3f.

The smallest pillars that we were able to fabricate with the present
fabrication method was around one micron which yields the maximum quality
factor of $9.1\times 10^{8}$, and a Qf-product of $1.3\times 10^{15}$
(corresponding to about 216 quantum coherent oscillations) which we believe is
the largest reported Qf-product at room temperature for membranes and similar
to the best measured perimeter mode of a polygon-shaped resonator [19].
Moreover, according to our model, the quality factor will increase to beyond
one billion by further decreasing the periodicity. Furthermore, due to the
structural functionality of our system, it is possible to enhance the strain
of the membrane and thereby further decreasing the intrinsic losses through
dissipation dilution. One particular strategy for increasing the stress is to
use the design depicted in Fig. 4 where the stress of the central pad (on
which the density phononic crystal could be located) is reaching approximately
3 GPa. In addition to improving the quality factor further by strain
engineering, it is also possible to improve it by choosing optimised
variations of the material density. E.g. in the supplementary material we have
considered a sinusoidal variation of the material density by which the quality
factor is slightly improved, and we expect even further improvements by
optimizing the density variation using for example the mathematical framework
of topology optimization [25].

In conclusion, we proposed and demonstrated a new method for phononic crystal
engineering of micro- and nanomechanical membranes that realizes acoustic
shielding and soft-clamping, and that leaves a lot of flexibility for
parameter optimization and strain-engineering. The fabricated devices exhibits
quality factors as high a $9.1\times 10^{8}$ at room temperature, and with
further optimization of pillar size (see Fig. 2c), pattern design (see Fig. 8b
in supplementary materials) and strain-engineering (see Fig. 4b), we
anticipate quality factors approaching 10 billion at room temperature. The
method is material independent and we expect it can be combined with
crystalline materials for even greater quality factors [26, 9]. With such
quality factors, it is possible to bring the membrane into the quantum ground
state at room temperature via feedback cooling (possibly without the use of a
cavity) and to eventually prepare non-classical mechanical states such as
squeezed and entangled states with applications in quantum repeaters [27] and
sensors. Even more intriguingly, the high-quality of our membranes might
eventually facilitate the preparation of quantum non-Gaussian mechanical
states through the interaction with two-level systems such as defect centers
in diamond [28] or superconducting qubits [29], and the realization of quantum
transduction to coherently connect two disparate systems with close to unity
efficiency [30, 31, 32]. Finally, in the current work, we have focused on the
improvement of the quality factor, and the Qf-product, which is of high
relevance for many of the applications mentioned above. However, in a future
work, it would be interesting to optimize the Q/f-fraction for applications in
atomic force microscopy or optimize the $Q\times m$-product (where m is the
effective mass) as required for testing alternative models of quantum
mechanics [33, 34, 35].

Figure 1: Concept illustration: (a) The desired phononic crystal pattern of
high and low material density regions. To approximate this behaviour the
distribution is used as a mask on-top an even distribution of nano-pillars
shown in (b). The end result is an engineered distribution of nano-pillars
based on any desired material distribution as illustrated in (c). The zoomed-
out inset of (b) highlights the relevant pillar parameters. (d) The increased
losses relative to a hypothetical membrane (without pillars) of increased
material density while keeping the total mass constant. It assumes a 20 nm
thick stoichiometric silicon nitride membrane vibrating at 1.1 MHz. The
results are practically identical to the other studied frequencies. (e)
Illustration of evanescent bending losses at the pillar base, estimated using
3D simulations. Figure 2: Phononic analysis: (a) and (b) show the relative
bandgap, $\Delta f/f_{c}$, and quality factor, $Q$, respectively, versus the
phononic crystal contrast, $g$, and the relative domain size, $\alpha_{w}$
(fixing the periodicity to $a_{ph}=253\ \textrm{\textmu{}m}$ membrane. (c)
Predicted quality factors versus number of isolating phononic crystal cell
layers. The insert illustrates a membrane with $N=4$. In (b) and (c), the Q
factor was numerically estimated without accounting for the additional pillar
loss (assuming very small pillars) and without accounting for acoustic
radiation loss. (d) and (e) illustrate the density distribution and the
modeshape of the phononic membrane, respectively. (f) The band structure of
the phononic crystal, where the shaded area highlights the bandgap. Figure 3:
Experimental results: (a) Microscope image of the fabricated density phononic
membrane. (b) SEM image of the nanopillars of a specific membrane. (c)
Spectrum of a density phononic membrane where the band-gap region is
highlighted by a shaded grey. (d) An example of a ringdown measurement of a
specific membrane exhibiting a Q factor of 5.65$\times$108. (e) Measured
quality factors versus number of phononic cells $N$ for $a_{\textrm{pil}}=2\
\textrm{\textmu m}$. (f) Measured quality factor versus pillar periodicity
$a_{\textrm{pil}}$ for $N=8$. The solid lines correspond to fitted simulation
results based on the best samples. Figure 4: 2D strain engineering: (a)
Theoretical design of a stoichiometric silicon nitride membrane consisting of
two different thicknesses (20nm and 200nm) which enables stress enhancement of
two-dimensional membranes. (b) Illustration of the resulting von Mises tensile
stress distribution.

##  I Methods

Simulations: The simulations of the membranes and nanopillars were performed
in COMSOL Multiphysics. The nanopillars were simulated in full 3D for a single
pillar using Floquet periodic boundary conditions. The wavenumber was defines
as $k=2\pi f\sqrt{g\rho/\sigma}$ (where $\sigma$ is the tensile stress and $f$
is the frequency) and the damping was extracted in post-processing by
integrating the strains over the entire volume. The loss factor and other
material parameters were assumed constant and isotropic. Based on these
results, the membranes were modelled in 2D assuming a uniform tensile stress
distribution and by adjusting the material density using
$\rho\to\rho(x,y)=g(x,y)\rho$. The losses in the density-increased domains are
scaled based on the chosen pillar periodicity $a_{\textrm{pil}}$ and pillar
simulation results. The quality factor $Q=2\pi W/\Delta W$ is then estimated
in post-processing using the kinetic energy

$W=2\pi^{2}f^{2}t\iint\rho(x,y)u^{2}(x,y)\ \mathrm{d}x\mathrm{d}y$ (2)

and dissipation via bending

$\Delta W=\frac{\pi Et^{3}}{12(1-\nu^{2})}\iint\phi(x,y)\\\
\cdot\left(\frac{\partial^{2}u(x,y)}{\partial
x^{2}}+\frac{\partial^{2}u(x,y)}{\partial y^{2})}\right)^{2}\
\mathrm{d}x\mathrm{d}y$ (3)

where $t,E,\nu$ are the thickness, Young’s modulus and Poisson’s ratio of the
resonator material. $u(x,y)$ is the mode shape and $\phi(x,y)$ is the
effective loss angle of the membrane which is affected by the choice of
material, its thickness, and the presence of pillars.

Fabrication: Double-sided polished Silicon wafers (of 500 µm thickness) were
initially deposited with 20 nm stoichiometric silicon nitride via low-pressure
chemical vapor deposition, followed by a layer of 1000 nm silicon nitride
using plasma-enhanced chemical vapor deposited (PECVD). The pillar pattern was
transferred onto a spincoated layer of photo-resist via UV-lithography
followed by reactive ion ethcing almost through the PECVD layer. Another round
of lithography and reactive ion etch was performed on the backside to define
the membrane window. Finally, a short etch in buffered hydrofluoric acid
followed by an etch in potassium hydroxide ensured the removal of the residual
PECVD layer and release of the membrane. The membranes were then cleaned in
hydrochloric acid followed by sulfiric acid mixed with ammoniumpersulfate and
later dried in ethanol fumes.

Optical characterization: The phononic spectra as well as the Q factors are
measured using optical interferometry at 1550 nm. The membrane under study is
placed inside a vacuum chamber with pressure $<$10-7mBar and excited with a
modulated laser beam reflecting off its central part. Once excited, the laser
modulation is turned off, and the phase shift of the reflected beam is
measured as a function of time with a phase-locked homodyne detector. The
signal is used to record the phononic mode spectrum as well as the amplitude
decay of the phase modulation with specific examples shown in Fig. 3c and Fig.
3d.

## References

  * Maldovan [2013] M. Maldovan, Sound and heat revolutions in phononics, Nature 503, 209 (2013).
  * Vasileiadis et al [2021] T. Vasileiadis et al, Progress and perspectives on phononic crystals, J. Appl. Phys. 129, 160901 (2021).
  * Song  _et al._ [2020] A. Song, J. Li, C. Shen, T. Chen, and S. A. Cummer, Switchable directional sound emission with improved field confinement based on topological insulators, Appl. Phys. Lett. 117, 043503 (2020).
  * Yu et al [2014] P.-L. Yu et al, A phononic bandgap shield for high-Q membrane microresonators, Applied Physics Letters 104, 023510 (2014).
  * Tsaturyan et al [2014] Y. Tsaturyan et al, Demonstration of suppressed phonon tunneling losses in phononic bandgap shielded membrane resonators for high-Q optomechanics, Optics Express 22, 6810 (2014).
  * Ghadimi et al [2017] A. Ghadimi et al, Radiation and Internal Loss Engineering of High-Stress Silicon Nitride Nanobeams, Nano Letters 17, 3501 (2017).
  * Tsaturyan  _et al._ [2017] Y. Tsaturyan, A. Barg, E. S. Polzik, and A. Schliesser, Ultracoherent nanomechanical resonators via soft clamping and dissipation dilution, Nature Nanotechnology 12, 776 (2017), arXiv:1608.00937 .
  * Ghadimi  _et al._ [2017] A. H. Ghadimi, D. J. Wilson, and T. J. Kippenberg, Radiation and Internal Loss Engineering of High-Stress Silicon Nitride Nanobeams, Nano Letters 17, 3501 (2017).
  * Beccari  _et al._ [2021] A. Beccari, D. A. Visani, S. A. Fedorov, M. J. Bereyhi, V. Boureau, N. J. Engelsen, and T. J. Kippenberg, Strained crystalline nanomechanical resonators with ultralow dissipation, Nature Physics , 436 (2021), arXiv:2107.02124v1 .
  * Verbridge  _et al._ [2006] S. S. Verbridge, J. M. Parpia, R. B. Reichenbach, L. M. Bellan, and H. G. Craighead, High quality factor resonance at room temperature with nanostrings under high tensile stress, Journal of Applied Physics 99, 10.1063/1.2204829 (2006).
  * Reinhardt  _et al._ [2016] C. Reinhardt, T. Müller, A. Bourassa, and J. C. Sankey, Ultralow-noise SiN trampoline resonators for sensing and optomechanics, Physical Review X 6, 1 (2016), arXiv:1511.01769 .
  * Norte et al [2016] R. Norte et al, Mechancial resonators for quantum optomechanical experiments at room temperature, Phys. Rev. Lett. 116, 147202 (2016).
  * Pratt  _et al._ [2021] J. R. Pratt, A. R. Agrawal, C. A. Condos, C. M. Pluchar, S. Schlamminger, and D. J. Wilson, Nanoscale torsional dissipation dilution for quantum experiments and precision measurement,  (2021), arXiv:2112.08350 .
  * de Jong  _et al._ [2022] M. H. J. de Jong, M. A. ten Wolde, A. Cupertino, P. G. Steeneken, and R. A. Norte, Mechanical dissipation by substrate-mode coupling in SiN resonators,  , 1 (2022), arXiv:2206.00990 .
  * Schmid  _et al._ [2011] S. Schmid, K. D. Jensen, K. H. Nielsen, and A. Boisen, Damping mechanisms in high-Q micro and nanomechanical string resonators, Physical Review B - Condensed Matter and Materials Physics 84, 1 (2011).
  * Villanueva and Schmid [2014] L. G. Villanueva and S. Schmid, Evidence of surface loss as ubiquitous limiting damping mechanism in SiN micro- and nanomechanical resonators, Physical Review Letters 113, 1 (2014), arXiv:arXiv:1405.6115v1 .
  * Bereyhi et al [2022a] M. Bereyhi et al, Hierarchical tensile structures with ultralow mechanical dissipation, Nature Comm. 13, 3097 (2022a).
  * Shin et al [2022] D. Shin et al, Spiderweb nanomechanical resonators via bayesian optimization: Inspired by Nature and guided by machine learning, Adv. Mater. 34, 2106248 (2022).
  * Bereyhi et al [2022b] M. Bereyhi et al, Perimeter Modes of Nanomechanical Resonators Exhibit Quality Factors Exceeding $10^{9}$ at Room Temperature, Physical Review X 12, 021036 (2022b).
  * Rossi et al [2018] M. Rossi et al, Measurement-based quantum control of mechanical motion, Nature 563, 7729 (2018).
  * Mason et al [2019] D. Mason et al, Continuous force and displacement measurement below the standard quantum limit, Nature Physics 15, 745 (2019).
  * Chen et al [2020] J. Chen et al, Entanglement of propagating optical modes via a mechanical interface, Nature communications 11, 943 (2020).
  * Gonzalez and P.R. [1994] G. Gonzalez and S. P.R., Brownian motion of a mass suspended by an anelastic wire, The Journal of the Acoustical Society of America 96, 207 (1994).
  * Federov et al [2019] S. A. Federov et al, Generalized dissipation dilution in strained mechanical resonators, Physical Review B 9, 054107 (2019).
  * Høj  _et al._ [2021] D. Høj, F. Wang, W. Gao, U. B. Hoff, O. Sigmund, and U. L. Andersen, Ultra-coherent nanomechanical resonators based on inverse design, Nature Communications 2021 12:1 12, 1 (2021), arXiv:2103.15601 .
  * Romero  _et al._ [2020] E. Romero, V. M. Valenzuela, A. R. Kermany, L. Sementilli, F. Iacopi, and W. P. Bowen, Engineering the Dissipation of Crystalline Micromechanical Resonators, Physical Review Applied 13, 1 (2020), arXiv:1911.03835 .
  * Fiaschi et al [2021] N. Fiaschi et al, Optomechanical quantum teleportation, Nat. Photon. 15, 817 (2021).
  * Lee et al [2017] D. Lee et al, Topical review:spins and mechanics in diamond, J. Opt. 19, 033001 (2017).
  * Clerk et al [2020] A. Clerk et al, Hybrid quantum systems with circuit quantum electrodynamics, Nat. Phys. 16, 257 (2020).
  * Barzanjeh et al [2022] S. Barzanjeh et al, Optomechanics for quantum technologies, Nature Physics 18, 15 (2022).
  * Lauk et al [2020] N. Lauk et al, Perspectives on quantum transduction, Quantum Sci. Technol. 5, 020501 (2020).
  * Delaney et al [2022] R. Delaney et al, Superconducting-qubit readout via low-backaction electro-optic transduction, Nature 606, 489 (2022).
  * Manley et al [2021] J. Manley et al, Searching for vector dark matter with an optomechanical accelerometer, Physical Review Letters 126, 061301 (2021).
  * Nimmrichter et al [2014] S. Nimmrichter et al, Optomechnical sensing of spontaneous wavefunction collapse, Physical Review Letters 113, 020405 (2014).
  * Westphal et al [2021] T. Westphal et al, Measurement of gravitational coupling between millimeter-sized masses, Nature 591, 225 (2021).


# Charged Higgs contribution to rare top-quark decay $t\to b\bar{b}bW^{+}$

S. Slabospitskii  NRC “Kurchatov Institute” - IHEP, Protvino, Moscow Region,
Russia<EMAIL_ADDRESS>
###### Abstract

The calculation of the rare $t$-quark decay $t\to b\,\bar{b}\,b\,W^{+}$ within
the Standard Model as well as the charged Higgs contribution to this decay is
presented. The role of possible background processes is discussed. It is shown
that this decay provides an additional possibility for charged Higgs searches.
PACS: 12.38.-t, 14.54Ha

###### keywords:

top-quark , rare decay , charged Higgs

##  1 Introduction

The physics of the top-quark is one the main areas of study in current high
energy physics. Indeed almost all characteristics of the processes with top-
quark can be calculated with high theoretical accuracy. The decay $t\to bW$ is
by far the dominant one with the Standard Model (SM). All other decay channels
have very small branching fractions and are predicted to be smaller by several
orders of magnitude in the SM [1, 2].

On the other hand, rare decays of the $t$-quark are very sensitive to the
manifestation of New Physics beyond the SM. A very well known example is top-
quark decays due to anomalous flavor changing neutral currents [1].

In this article we study a rare top-quark decay:

$\displaystyle t\,\to\,b\,W^{+}\,\bar{b}\,b$

as well as the charged Higgs boson contribution to this decay. It is shown
that this decay has well identified final states and has the relatively large
branching fraction, $\mathsf{Br}(t\to
b\,\bar{b}\,b\,W^{+})\sim\mathcal{O}(10^{-3})$. We demonstrate that this four-
body decay channel provides an additional way to search for charged
$H^{\pm}$-boson.

Charged Higgs bosons appear in the scalar sector of several SM extensions, and
is the object of various beyond the Standard Model searches at the LHC. In
this article we explore a generic two-Higgs-doublet model (2HDM), which is one
of the simplest SM extensions featuring a charged scalar [3, 4]. Within this
class of models, two isospin doublets are introduced to break the $SU(2)\times
U(1)$ symmetry, leading to the existence of five physical Higgs bosons, two of
which are charged particles ($H^{\pm}$). The latest constraints the allowed
$H^{\pm}$ mass range as a function of $\tan\beta$ can be found in [5, 6].

Searches for $H^{\pm}$have been performed at LEP [7], at the Fermilab Tevatron
[8, 9]. The ATLAS and CMS Collaborations have covered several $H^{\pm}$ decay
channels, such as $\tau\nu_{\tau},\;tb,\;cs,\;cb$ (see [10, 11, 12, 13, 14,
15, 16] and references therein). Note, that $H^{\pm}$-boson (like a SM Higgs)
has couplings proportional to fermion masses.

The large value of the $t\,H^{\pm}\,b$ coupling leads to the fact that
processes involving top quarks are considered for the charged Higgs production
processes. The experiments explore three scenarios for the search for a
charged Higgs boson. For “light” $H^{\pm}$-boson ($M(H^{\pm})<m_{t}$) top-
quark decay channel into charged Higgs is used. For the case when $M(H^{\pm})$
is greater than $t$-quark mass, it is considered the charged Higgs production
with subsequent decay into $t\bar{b}$ state. As a result the charged Higgs
with masses

$\displaystyle M(H^{\pm})<160\;\;\hbox{GeV}\quad\hbox{and}\quad
M(H^{\pm})>180\;\;\hbox{GeV}$ (1)

with a wide range of $\tan\beta$ are excluded [10, 11, 12, 13, 14, 15, 16].

The third scenario explores the production of the Higgs boson in association
with a top quark (see [11]) with subsequent $H^{\pm}\to\tau^{\pm}\nu_{\tau}$
decay. As a result, the charged Higgs with mass range $(90-2000)$ GeV and
$\tan\beta\raisebox{-3.0pt}{$\>\stackrel{{\scriptstyle>}}{{\scriptstyle\sim}}\>$}1$
is excluded.

Note, that in all experiments, a search is made for the direct production of
charged Higgs with subsequent decays in observed fermions. In the present
articles we propose an additional “indirect” way for $H^{\pm}$ searches by
study of a rare four-body top-quark decay channel ($t\to bW^{+}b\bar{b}$). It
is shown, that $t$-quark decay ($t\to bW^{+}b\bar{b}$) provides a reasonable
possibility to detect charged Higgs contribution.

Bearing in mind the experimental constraint (1) in the present article we
investigate the role of the charged $H^{\pm}$-boson with mass range and
$\tan\beta$ as follows:

$\displaystyle M(H^{\pm})=(160\,\div
180)\;\;\hbox{GeV}\quad\hbox{and}\quad\tan\beta\raisebox{-3.0pt}{$\>\stackrel{{\scriptstyle<}}{{\scriptstyle\sim}}\>$}1$
(2)

Throughout of this article we follow [17] for the notations, the SM vertices
and SM parameters. For numerical calculations of the decay widths and
production processes the C++ version of the TopReX package [18] is used. All
calculations were done the following quark masses:

$\displaystyle\left.\begin{array}[]{lrlr}m_{d}=&0.33\;\;\hbox{GeV},&m_{u}=&0.33\;\;\hbox{GeV},\\\
m_{s}=&0.5\;\;\hbox{GeV},&m_{c}=&1.5\;\;\hbox{GeV},\\\
m_{b}=&4.8\;\;\hbox{GeV},&m_{t}=&172.5\;\;\hbox{GeV},\end{array}\right.$ (6)

We explore the $b$-jet tagging with $b$-quark identification efficiency of
70%, with probability to misidentify $c$ quark, and light-flavor quark and
gluon jets as $b$-jets of approximately 10% and 1.%, respectively (see, e.g
[19]):

$\displaystyle\epsilon_{b}=0.7,\epsilon_{c}=0.1,\epsilon_{q,\;g}=0.01$ (7)

The article is organized as follows. Charged Higgs interaction model is
described in the Section 2. The constraints on $\tan\beta$ evaluated with
taking into account for $t$-quark total decay width, as well as the branching
fractions of $t$-quark decays, is presented in the Section 3. The width of the
top-quark decay $t\to b\,\bar{b}\,b\,W$ is calculated in the Section 4. A
brief comment on $t$-quark decay to a multilepton final state is given in the
Section 5. The different background sources are considered in the Sections 6
and 7. The kinematic of the decay $t\,\to\,b\,\bar{b}\,b\,W^{+}$ is discussed
in the Section 8. The last Section 9 summarize the obtained results.

##  2 Charged Higgs interaction Lagrangian

The interaction Lagrangian describing the $H^{\pm}$-boson-fermions
interactions in the MSSM is (see [1], [4], [3]):

$\displaystyle{\cal L}$ $\displaystyle=$
$\displaystyle\frac{g}{\sqrt{2}M_{W}}H^{+}\\{V_{ud}\overline{u}(m_{u}\cot\beta
P_{L}+m_{d}\tan\beta P_{R})d+\,\overline{\nu}(\tan\beta
m_{\ell}P_{R})\ell\\},$ (8)

where $P_{L/R}=\frac{1}{2}\left(1\mp\gamma^{5}\right)$, symbols $u$ and $d$
stand for ”up” ($u,c,t$) and ”down” ($d,s,b$) quarks, $\nu$ and $\ell$ for
neutrino and charged leptons, $V_{ud}$ is the Cabibbo-Kobayashi-Maskawa matrix
element.

At tree level the corresponding decay widths are equal to [1]:

$\displaystyle\Gamma\left(H^{+}\rightarrow l\nu\right)$ $\displaystyle=$
$\displaystyle\frac{g^{2}M_{H}}{32\pi M^{2}_{W}}m^{2}_{l}\tan^{2}\beta$ (9)
$\displaystyle\Gamma\left(H^{+}\rightarrow q\bar{q}\right)$ $\displaystyle=$
$\displaystyle\frac{3g^{2}}{32\pi
M^{2}_{W}M_{H}}\left|V_{q\bar{q}}\right|^{2}\lambda^{1/2}\left(1,\frac{m^{2}_{q}}{M^{2}_{H}},\frac{m^{2}_{\bar{q}}}{M^{2}_{H}}\right)\times$
(10)
$\displaystyle\left[\left(M^{2}_{H}-m^{2}_{\bar{q}}-m^{2}_{q}\right)\left(m^{2}_{q}\cot^{2}\beta+m^{2}_{\bar{q}}\tan^{2}\beta\right)-4m^{2}_{q}m^{2}_{\bar{q}}\right]$
$\displaystyle\Gamma\left(t\rightarrow bH^{+}\right)$ $\displaystyle=$
$\displaystyle\frac{g^{2}}{64\pi
M^{2}_{W}m_{t}}\left|V_{tb}\right|^{2}\lambda^{1/2}\left(1,\frac{m^{2}_{b}}{m^{2}_{t}},\frac{M^{2}_{H}}{m^{2}_{t}}\right)\times$
(11)
$\displaystyle\left[\left(m^{2}_{t}+m^{2}_{b}-M^{2}_{H}\right)\left(m^{2}_{t}\cot^{2}\beta+m^{2}_{b}\tan^{2}\beta\right)+4m^{2}_{t}m^{2}_{b}\right]$

where

$\lambda\left(a,b,c\right)=a^{2}+b^{2}+c^{2}-2\left(ab+ac+bc\right)$

Note, that for small $\tan\beta$ the additional decay channel ($H^{+}\to
t^{*}\bar{b}\to W^{+}b\bar{b}$) can get a noticeable contribution to the total
charged Higgs decay width [20] (see diagram in Fig. 1).

Figure 1: Diagram describing $H^{+}\to t^{*}\bar{b}\to W^{+}b\bar{b}$ decay
due to virtual $t$-quark contribution.

In the Fig. 2 we present the behavior of the branching fractions for three
$H^{\pm}$ decay modes ($c\bar{s},\,\tau\nu_{tau},\;b\bar{b}W$) for four values
of $m_{H^{+}}$ versus $\tan\beta$. As it seen for small $m_{H^{+}}$ values
($m_{H^{+}}\sim 100$ GeV) this three-body decay ($H^{+}\to b\bar{b}W^{+}$) is
almost negligible. Consequently, for these $m_{H^{+}}$-values one expects
$H^{+}\to\tau^{+}\nu$ (favoured for large $\tan\beta$) and/or $H^{+}\to
c\bar{s}$ (favoured for small $\tan\beta$) to be the dominant decays.

On the other hand, if
$\tan\beta\raisebox{-3.0pt}{$\>\stackrel{{\scriptstyle<}}{{\scriptstyle\sim}}\>$}2$
and
$m_{H^{+}}\raisebox{-3.0pt}{$\>\stackrel{{\scriptstyle>}}{{\scriptstyle\sim}}\>$}150$
GeV, the large mass (or coupling) of the $t$-quark causes
$\mathsf{Br}(H^{+}\to b\bar{b}W^{+})$ to exceed $\mathsf{Br}(H^{+}\to
c\bar{s})$. Moreover, for $\tan\beta<1$ this three-body decay mode become
dominant [1, 20] (see Fig. 2).

Figure 2: Branching fractions for three $H^{\pm}$ decay modes for four values
of $m_{H^{+}}$ vs. $\tan\beta$.

##  3 Constraints on charged Higgs parameters

As explained in the introduction, we consider the charged Higgs boson with
mass from 160 GeV to 180 GeV and small values of $\tan\beta$ (this area has
not yet been studied in detail in the experiments).

In this section we present very rough estimates for the limits of the allowed
region of $\tan\beta$ values. The charged Higgs (if it exists) should make a
contribution into the total $t$-quark decay width as well as to the branching
fractions. These values are taken from “Review of Particle Properties” [2]).

$\displaystyle\Gamma_{t}$ $\displaystyle=$ $\displaystyle
1.42^{+0.19}_{-0.15}\;\;\;\hbox{GeV}$ (12) $\displaystyle\mathsf{Br}(t\to
bq\bar{q})$ $\displaystyle=$ $\displaystyle(66.5\pm 1.4)\%$ (13)
$\displaystyle\mathsf{Br}(t\to b\,\tau\nu_{\tau})$ $\displaystyle=$
$\displaystyle(11.1\pm 0.9)\%$ (14)

The first constraints could be obtained from requirement that the contribution
to the total top-quark decay width could not exceed
$\Delta\Gamma_{t}^{H}\approx 0.2$ GeV. Additional constrains could be deduced
be using of the top-quark branching fractions (13) and (14).

Using the errors of these branching fractions we require that charged Higgs
contribution to branching fractions of these decays ($t\to bq\bar{q}$ and
$t\to b\,\tau\nu_{\tau}$) should be less then the errors from (13-14).
Therefore, we explore the constraints as follows:

$\displaystyle(1)$ $\displaystyle\Delta\Gamma_{t}^{H}<0.2\;\;\hbox{GeV}$ (15)
$\displaystyle(2)$
$\displaystyle\left.\begin{array}[]{lcl}\Delta\mathsf{Br}^{(H)}(t\to
bq\bar{q})=\displaystyle\frac{\Gamma^{(H)}(t\to
bq\bar{q})}{\Gamma_{LO}+\Gamma^{(H)}(t\to bH^{+})}&<&0.014\\\
\Delta\mathsf{Br}^{(H)}(t\to
b\tau\nu_{tau})=\displaystyle\frac{\Gamma^{(H)}(t\to
bq\bar{q})}{\Gamma_{LO}+\Gamma^{(H)}(t\to
bH^{+})}&<&0.009\end{array}\right\\}$ (18)

where $\Gamma_{LO}$ is $t$-quark decay width evaluated at the leading order:
$\Gamma_{LO}=\Gamma(t\to b\,W^{+})=1.47$ GeV [1]

The resulted allowed parameter ranges are presented in the Table 1 and in Fig.
3. It follows form this table and figure that the more “narrow” limits are
resulted from the branching fraction constraints.

Table 1:  The allowed charged Higgs parameters range. The symbols (1) and (2) corresponds to the constraints (15) and (18). The mass of $H^{\pm}$-boson is in GeV.  $M(H^{\pm})$ | $\tan\beta$ | $M(H^{\pm})$ | $\tan\beta$
---|---|---|---
| (1) | (2) |  | (1) | (2)
$120$ | $1.50\div 24.$ | - | $160$ | $0.37\div 93$ | $0.55\div 30.$
$125$ | $1.4\div 27.$ | $5.9\div 6.1$ | $165$ | $0.20\div 115$ | $0.26\div 47$
$130$ | $1.27\div 29.$ | $5.\div 8.$ | $170$ | $0.05\div 135$ | $0.089\div 74$
$140$ | $0.98\div 37$ | $3.3\div 10.7$ | $175$ | $0.045\div 150$ | $0.079\div 84$
$150$ | $0.70\div 52$ | $2.0\div 15.8$ | $180$ | $0.032\div 165$ | $0.071\div 92$
Figure 3: The allowed charged Higgs parameters range. The wide area correspond
to the constraint (15), while the small one is resulted by using of the
constraints (18).

Certainly, these results should be considered as very rough estimates. To get
more correct results one needs to take into account the NLO calculations and
the running masses of the quarks.

Moreover, we do not compare these limits (see Table 1) with results from [5,
6], because, in fact, we use a rather narrow ”working” area of $\tan\beta$
parameter (see Section 4 and Fig. 5). We use these range for $M(H^{\pm})$ (2)
and $\tan\beta$ in further calculations. In addition, for illustration, we
present the total decay widths of the charged Higgs boson
$\Gamma_{tot}(H^{\pm})$ calculated within these limits (see Table 2).

Table 2:  The total charged Higgs boson decay width $\Gamma_{tot}(H^{\pm})$ $\tan\beta$ | $160$ | $165$ | $170$ | $175$ | $180$
---|---|---|---|---|---
$0.1$ | - | - | $0.62$ | $1.18$ | $6.96$
$0.2$ | - | - | $0.154$ | $0.30$ | $1.75$
$0.35$ | - | $0.031$ | $0.051$ | $0.096$ | $0.051$
$0.5$ | $0.010$ | $0.0154$ | $0.025$ | $0.048$ | $0.028$
$1$ | $0.0028$ | $0.004$ | $0.0065$ | $0.0132$ | $0.022$
$10$ | $0.046$ | $0.047$ | $0.049$ | $0.051$ | $0.058$
$50$ | - | $1.18$ | $1.22$ | $1.26$ | $1.51$

##  4 Top-quark decay $t\to b\,W^{+}\,b\,\bar{b}$

In this article we consider charged Higgs contribution to rare four-body top-
quark decay:

$\displaystyle t\to bW^{+}\;b\bar{b}$ (19)

Within SM this decay process is described by 28 Feynman diagrams (see Fig. 4).
The picture $(a)$ corresponds to three diagrams with virtual quarks $q_{U}$
exchange ($q_{U}=t,\;u$ or $q_{U}=c$).

Figure 4: The diagrams($(a),(b),(c)$ and $(d)$) describing $t\to
bW^{+}\;b\bar{b}$ decay within SM. The diagram $(h)$ describes charged Higgs
contribution to this decay.

In the same fashion the pictures $(b)$ and $(c)$ correspond to eight diagrams
with virtual bosons $B$ exchange ($B=g$, $B=\gamma$, $B=Z$ or $B=H$). The
picture $(d)$ corresponds to three diagrams with virtual bosons $V$ exchange
($V=\gamma$, $V=Z$ or $V=H$). The picture $(h)$ corresponds to three diagrams
with charged Higgs and with virtual quarks $q_{U}$ exchange
($q_{U}=t,\,u,\,c$). Note, that the light virtual quarks ($u$ and $c$) as well
virtual photon, $Z$ and $H$ bosons give very small contributions to the $t\to
bW^{+}\;b\bar{b}$ decay width.

The calculated decay width and branching fraction are as follows:

$\displaystyle\left.\begin{array}[]{lcl}\Gamma(t\to bW^{+}b\bar{b})&=&(9.30\pm
0.03)\times 10^{-4}\;\;\;\hbox{GeV}\\\ \mathsf{Br}(t\to
bW^{+}b\bar{b})&=&(6.29\pm 0.02)\times 10^{-4}\end{array}\right.$ (22)

In the Table 3 we present the results for several charged Higgs mass values.
As it follows from calculations the partial decay width for three variants of
the parameters ($\huge\\{M(H^{+})=160$ GeV, $\tan\beta<0.35\huge\\}$,
$\huge\\{M(H^{+})=165$ GeV, $\tan\beta<0.2\huge\\}$, and
$\huge\\{M(H^{+})=170$ GeV, $\tan\beta<0.05\huge\\}$) is exceed value 0.2 GeV.
Consequently, this parameter space range should also be excluded from our
further discussion.

Table 3:  The partial decay widths of the top-quark decay $t\to b\bar{b}bW$ with account of $H^{\pm}$ contribution. All width values (in GeV) are multiplied by $10^{4}$. The “SM” means that charged Higgs contribution is less the 5% and the corresponding width equals $\Gamma_{SM}\approx 9.3\times 10^{-4}$ GeV. For $tan\beta>50$ and $M_{H}\geq 170$ GeV the decay widths are equal to $\Gamma_{SM}$ $\tan\beta$ | $160$ | $165$ | $170$ | $175$ | $180$ | $185$
---|---|---|---|---|---|---
$0.1$ | - | - | $282\pm 2.4$ | $122\pm 1.2$ | $75.1\pm 0.6$ | $50.6\pm 0.4$
$0.2$ | - | - | $28.6\pm 0.0.2$ | $17.7\pm 0.23$ | $14.5\pm 0.20$ | $13.1\pm 0.20$
$0.25$ | - | $1631\pm 211$ | $17.8\pm 0.13$ | $13.1\pm 0.1$ | $11.7\pm 0.09$ | $11.0\pm 0.09$
$0.5$ | $1110\pm 340$ | $388\pm 47$ | $10.3\pm 0.19$ | SM | SM | SM
$1$ | $254\pm 84$ | $114\pm 37$ | SM | SM | SM | SM
$30$ | $9.9\pm 0.09$ | SM | SM | SM | SM | SM
$50$ | - | $10.3\pm 0.11$ | SM | SM | SM | SM

The Fig. 5 presents the same plot as in Fig. 3 (the allowed charged Higgs
parameters range). The hatched area corresponds to cases when the charged
Higgs contribution to the $t\to b\,\bar{b}\,b\,W^{+}$ decay width exceeds 5%.

Figure 5: The same plot as in Fig. 3 (the allowed charged Higgs parameters
range). The hatched area corresponds to cases when the charged Higgs
contribution to the $t\to b\,\bar{b}\,b\,W^{+}$ decay width exceeds 5%.

There is one important comment concerning the final $W$-boson decay into
observed particles. For $W$-boson decaying into lepton pair ($W\to\ell\nu$) a
sharp peak near the $W$-boson should be expected, and the corresponding decay
width is

$\displaystyle\Gamma(t\to b\bar{b}b\,\ell\nu)=\Gamma(t\to
b\bar{b}b\,W^{+})\times\mathsf{Br}(W^{+}\to\ell^{+}\nu)$

However, this is not true for the $W$-boson decay into quarks. For example,
there are additional diagrams, where the $c\,\bar{s}$ pair is produced not
from $W$-boson (see Fig. 6).

Figure 6: The representative diagrams describing the contribution to $t\to
b\bar{b}b\,c\bar{s}$ decay, where the $c\,\bar{s}$ pair is not from the
$W$-boson.

Here we present few details concerning top-quark decay channel: $t\to
b\bar{b}b\,c\bar{s}$. The Fig. 7 presents the distribution over the invariant
mass of the $c\bar{s}$ pair in the $t\to b\bar{b}b\,c\bar{s}$ decay.

Figure 7: The distribution over the invariant mass of the $c\bar{s}$ pair in
the $t\to b\bar{b}b\,c\bar{s}$ decay. The right sharp peak corresponds to
virtual $W$-boson decay into $c\bar{s}$ pair, while the left wide region
corresponds to the case when the $s\,\bar{s}$ pair is not from the $W$-boson.

As it follows from calculations, the contribution to the decay width from this
“non-resonant” (i.e. $c\bar{s}$ is not from the $W*$-boson) region
approximately coincides with the contribution from $W\to q\bar{q}^{\prime}$
region. As a result, the decay with of the 5-body top-quark decay channel is
approximately twice the $t\to b\bar{b}b\,W$ decay width:

$\displaystyle\left\\{\begin{array}[]{lcl}\Gamma(t\to
b\bar{b}b\,\ell^{+}\nu)&=&(1.0\pm 0.2)\times 10^{-4}\;\;\hbox{GeV}\\\
\Gamma(t\to b\bar{b}b\,c\bar{s})&=&(6.8\pm 0.7)\times
10^{-4}\;\;\hbox{GeV}\end{array}\right.$ (25) $\displaystyle\Rightarrow$
$\displaystyle\Gamma(t\to
b\bar{b}b\,f\bar{f}^{\prime})=3\,\Gamma(\ldots\,\ell^{+}\nu)+2\,\Gamma(\ldots\,c\bar{s})=(16.6\pm
1.0)\times 10^{-4}\;\;\hbox{GeV}$ (26)

Therefore, the decay width value (26) is about twice as large than
$\Gamma(t\to b\bar{b}b\,W^{+})=9.3\times 10^{-4}$ GeV from (22).

This circumstance should be taken into account when searching for such a rare
decay of the top quark. Therefore, we think that the most efficient way is the
search for this decay in lepton mode $t\,\to\,b\,\bar{b}\,b\,\ell^{+}\nu$.
Therefore, in what follows the decay of the $W$-boson into leptons will be
assumed in this article.

##  5 Top-quark decay $t\,\to\,\ell^{+}_{1}\ell^{-}_{1}\ell_{2}^{+}\nu_{2}$

It was established many years ago, that at tree level, the decay $t\,b\,Z\,W$
has some peculiar features, since the process occurs near the kinematic
threshold ($m_{t}\sim M_{Z}+M_{W}+m_{b}$) [21, 22, 23, 24, 25, 26, 27]. A
brief discussion of this decay channel can be found also in [1]. Note, that
charged Higgs (at large $\tan\beta$ values) should also make a significant
contribution to the top-quark decay into four leptons:

$\displaystyle
t\,\to\,b\,Z^{*}(\ell^{+}_{1}\ell^{-}_{1})\;W^{+*}(\ell_{2}^{+}\nu_{2})\,\;\;+\;\;b\,Z^{*}(\ell^{+}_{1}\ell^{-}_{1})\,\widetilde{H^{+*}}(\ell_{2}^{+}\nu_{2})$
(27)

where $\widetilde{H^{+*}}$ corresponds to all charged Higgs contributions,
like $t\to b\,H^{+\,*}(\to\ell^{+}_{1}\,\ell^{-}_{1}\,\ell_{2}^{+}\,\nu_{2})$,
$t\to bZ^{*}\;H^{+*}(\to\ell_{2}^{+}\nu_{2})$, $\ldots$

Here, for the sake of completeness, we recalculate the width of this decay.
The results of calculations for the SM case and taking into account the
$H^{\pm}$-boson contribution, are given in the Table 4.

Table 4:  The partial widths for $t\to b\,\ell^{+}_{1}\,\ell^{-}_{1}\,\ell_{2}^{+}\,\nu_{2}$ decay. All width values (in GeV) are multiplied by $10^{8}$.  $M(H^{\pm}),\;\tan\beta$ | $\mu^{+}\mu^{-}e^{+}\nu_{e}$ | $\mu^{+}\mu^{-}\tau^{+}\nu_{\tau}$ | $\tau^{+}\tau^{-}\mu^{+}\nu_{\mu}$ | $\tau^{+}\tau^{-}\tau^{+}\nu_{\tau}$
---|---|---|---|---
SM | $6.1\pm 0.3$ | $6.1\pm 0.3$ | $6.1\pm 0.3$ | $5.7\pm 0.24$
$160,\;\;30$ | $6.13\pm 0.4$ | $11.8\pm 0.9$ | $28.\pm 4.$ | $34.7\pm 3.$
$170,\;\;75$ | $6.23\pm 0.5$ | $7.6\pm 0.5$ | $14.9\pm 0.7$ | $17.4\pm 0.6$

As can be seen from this table, for large values of $\tan\beta$ the
$H^{\pm}$-boson decay into $\tau$-leptons gives a noticeable contribution to
this partial width. However, the branching fraction is too small to for
experimental searches and will not be considered further.

##  6 Backgrounds 1. Four-body top-decays

In this section the main sources of the background processes to the
investigated decay $t\to bW^{+}b\bar{b}$ are discussed. We consider two types
of background processes: top-quark decay into other final states and
$t\bar{t}$ production in hadronic collisions accompanied by others partons
(for example $pp\,\to\,t\bar{t}\,b\bar{b}$).

We start with the other four-body top-quark decay channels:

$\displaystyle t\to bW^{+}\,gg\;\;\hbox{and}\;\;t\to
bW^{+}\,q\bar{q}^{\prime}$ (28)

There are about 30 decay channels with $(W\,q_{1}\bar{q}_{2}q_{3})$ final
state and three decay channels with gluons: $t\to qggW^{+},\;q=d,s,b$. The
most part of these decays have very small decay widths. Therefore, we present
the results for such decays that have widths at least 10% of the $t\to
b\,W\,b\,\bar{b}$ decay width. The values of the calculated widths are given
in the Table 5.

Table 5:  The partial widths for $t\to b\,W\,\bar{q}_{1}q_{2}(gg)$ decay channels. The widths are in GeV. The values in third column are evaluated with requirement that invariant mass of any two particle in the final state should be greater then 20 GeV. The fourth column (“ratio”) presents the ratios of the corresponding decay widths to the $\Gamma(bWb\bar{b})$ width, taking into account the $b$-tagging efficiency from (7).  $t\to b\,W^{+}\,\ldots$ | no cuts |  $M(i\,j)>20$ GeV | ratio
---|---|---|---
$b\,\bar{b}$ | $(9.30\pm 0.03)\times 10^{-4}$ | $(1.34\pm 0.03)\times 10^{-4}$ | 1
$g\,g$ |  | $(1.8\pm 0.04)\times 10^{-3}$ | 0.003
$c\,\bar{c}$ | $(4.2\pm 0.3)\times 10^{-3}$ | $(1.3\pm 0.3)\times 10^{-4}$ | 0.02
$s\,\bar{s}$ | $(1.1\pm 0.15)\times 10^{-2}$ | $(1.3\pm 0.3)\times 10^{-4}$ | 0.0002
$d\bar{d}\;(u\,\bar{u})$ | $(1.6\pm 0.4)\times 10^{-2}$ | $(1.3\pm 0.3)\times 10^{-4}$ | 0.0002

As it seen from third column of this table the decay widths of these channels
are comparable (or even greater) then the width of the main process
$\Gamma(t\to b\bar{b}b\,W^{+})=9.3\times 10^{-4}$ GeV from (22). However, the
$b$-tagging application provides a substantial suppression of the background
decay channels (see fourth column in the Table 5).

##  7 Backgrounds 2. The production processes

In this section we consider the backgrounds from $t\bar{t}$ hadronic
production accompanied by two partons:

$\displaystyle p\,p\,\to\,t\,\bar{t}\,j_{1}j_{2},\quad
j_{1}j_{2}=b\,\bar{b},\;\;g\,g,\;\;g\,q,\;\;q\,\bar{q}^{\prime}$ (29)

As it shown latter the main background comes from $t\bar{t}b\bar{b}$
production in $pp$-collisions. For $pp$-collisions at $\sqrt{s}=13$ TeV the
cross-section values were calculated [28, 29, 30, 31, 32, 33, 34] and were
measured by the ATLAS and CMS Collaborations [10, 35, 13].

Following these results (for $t\bar{t}b\bar{b}$) for further calculations we
use the values for cross section production of $t\bar{t}$ reaction
($\sigma(t\bar{t})=803\pm 2(\hbox{stat})\pm 25(\hbox{syst})\pm
20(\hbox{lumi})$ pb (see [36] and reference therein) and for
$t\bar{t}b\bar{b}$ process

$\displaystyle\left.\begin{array}[]{lcl}\sigma(t\bar{t})&=&800\;\;\;\hbox{pb}\\\
\sigma(t\bar{t}b\bar{b})&=&4\;\;\;\hbox{pb}\end{array}\right\\}$ (32)

The evaluated cross-section production values for reactions (29) at
$\sqrt{s}=13$ TeV are given in Table 6 below. These values are evaluated with
the following kinematic cuts:

$\displaystyle p_{\top}(j)\geq 20\;\;\hbox{GeV},\quad M(j_{1}j_{2})\geq
20\;\;\hbox{GeV}$ (33)

As it seen the requirement to have two $b$-tagged jets ($jj\to$ ($b$-tagging)
$\to j_{b}j_{b}$) provides the substantial suppression of the all background
processes except the $pp\to t\bar{t}b\bar{b}$ production.

Table 6:  The cross section values for the processes $pp\to t\bar{t}j_{1}j_{2}$ (where $j$ stands for quarks or gluon) with kinematic cuts $p_{\top}(j)\geq 20$ GeV and $M(j_{1}j_{2})\geq 20$ GeV. The third column presents cross section values times $b$-tagging efficiencies: $\epsilon_{b}=0.7,\epsilon_{c}=0.1,\epsilon_{q,\;g}=0.01$ from (7).  reaction |  $\sigma$ (pb) |  $\sigma\times\epsilon_{1}\epsilon_{2}$ (pb) | $\sigma/\sigma(t\bar{t}\,b\bar{b})$
---|---|---|---
$pp\to t\,\bar{t}\,b\,\bar{b}$ | $0.76\pm 0.08$ | $0.372\pm 0.040$ | $1$
$pp\to t\,\bar{t}\,g\,g$ | $57.3\pm 12.8$ | $(5.73\pm 1.28)\times 10^{-3}$ | $1.5\times 10^{-2}$
$pp\to t\bar{t}u\bar{u}$ | $0.61\pm 0.08$ | $(6.1\pm 0.8)\times 10^{-5}$ | $1.6\times 10^{-3}$
$pp\to t\bar{t}uu$ | $0.64\pm 0.28$ | $(6.4\pm 2.8)\times 10^{-5}$ | $1.7\times 10^{-4}$
$pp\to t\bar{t}d\bar{d}$ | $0.53\pm 0.08$ | $(5.3\pm 0.8)\times 10^{-5}$ | $1.4\times 10^{-3}$
$pp\to t\bar{t}dd$ | $0.26\pm 0.2$ | $(2.6\pm 2.)\times 10^{-5}$ | $7.0\times 10^{-5}$
$pp\to t\,\bar{t}\,c\,\bar{c}$ | $0.56\pm 0.11$ | $(5.6\pm 1.0)\times 10^{-3}$ | $1.5\times 10^{-2}$
$pp\to t\bar{t}\,b\,g$ | $0.37\pm 0.14$ | $(3.7\pm 0.14)\times 10^{-3}$ | $1.0\times 10^{-2}$
$pp\to t\bar{t}\,\bar{b}\,g$ | $0.44\pm 0.12$ | $(4.4\pm 1.2)\times 10^{-3}$ | $1.2\times 10^{-2}$

Thus, the main background process comes from the $t\bar{t}b\bar{b}$ production
in $pp$-collisions. This background can be suppressed by applying kinematic
cuts to the final partons. Indeed, the processes of $t\bar{t}$ production with
the subsequent four-body top-quark decay and $t\bar{t}b\bar{b}$ production
have the same final state:

$\displaystyle
p\,p\,\to\,t\,(\to\,b\bar{b}bW^{+})\;\;\bar{t}(\to\,\bar{b}W^{-})\;\;\;\;\;\;\to
b\bar{b}bW^{+}\bar{b}W^{-}$ (34) $\displaystyle
p\,p\,\to\,t\,(\to\,bW+)\;\;\;\;\;\bar{t}\,(\to\,\bar{b}W^{-})\;b\bar{b}\;\to
b\bar{b}bW^{+}\bar{b}W^{-}$ (35)

In what follows the $b$-quark and $\bar{b}$-quark will be considered as a
$b$-jet. Therefore, the reactions have to following (identical) final states:
$4j_{b}W_{1}W_{2}$. The pair $j_{b}$ and $W$, with a mass closest to the mass
of the top quark is treated as a top-quark decaying into $b$ and $W$. The
system system consisting of the remaining particles (three $b$-jets and $W$),
is considered as “top”-quark $T_{3b}$. Inside this $T_{3b}$-system we select a
pair of $j_{b}$ and $W$ with maximum mass. This pair is treated as a “virtual
top-quark” ($t^{*}$).

Note, that the distributions over the invariant masses of two and three
$b$-quarks and the virtual top mass are very different for the reactions (34)
and (35). The final state $3j_{b}W$ in the reaction (34) has a sharp peak (due
to top-quark decay). These distributions are presented in the Figs. 8.

Figure 8: $d\,N/M_{inv}$ distributions. The mass of two (three) $b$-quarks is
shown in top left (top right) plot. The distributions over $3j_{b}\,W$
(virtual top $t^{*}$) is shown in the bottom left (bottom right) plot.

Therefore, to suppress the background from $t\bar{t}b\bar{b}$ production (35)
we apply the kinematic cuts as follows:

$\displaystyle\left.\begin{array}[]{lcclcc}M(j_{1}\,j_{2})&\geq&20\;\;\hbox{GeV},&p_{\top}(j)&\geq&20\;\;\hbox{GeV}\\\
M(j_{b}\,j_{b})&\leq&60\;\;\hbox{GeV},&M(3j_{b})&\leq&90\;\;\hbox{GeV}\\\
M(3j_{b}\,W)&\leq&200\;\;\hbox{GeV},&M(t^{*})&=&100\div
150\;\;\hbox{GeV}\end{array}\right\\}$ (39)

After applying of all cuts from (39) the total efficiencies
($\varepsilon_{cuts}$) of event detection for the $T_{3b}(b\bar{b}bW)$ system
are as follows:

$\displaystyle\left.\begin{array}[]{lcl}pp\to
t\bar{t}&:&\varepsilon_{cuts}=0.086\\\ pp\to
t\bar{t}b\bar{b}&:&\varepsilon_{cuts}=0.00034\end{array}\right\\}$ (42)

As a result the application of all cuts (39) provides a rather well
suppression of the background process:

$\displaystyle
R=\frac{\sigma(t\bar{t}b\bar{b})\times\epsilon_{cuts}}{\sigma(t\bar{t})\times\epsilon_{cuts}\times\mathsf{Br}(t\to
3bW)}=0.031$ (43)

where $\sigma(t\bar{t}b\bar{b})$ and $\sigma(t\bar{t})$ are taken from (32).
Therefore, in what follows we consider the signal process (34) only.

Then, for estimation the expected number of events we use the following
options:
$\bullet$ the LHC Run-3 integrated luminosity equals $L_{tot}=300$ fb-1;
$\bullet$ $\mathsf{Br}(t\to b\bar{b}bW)=6.29\times 10^{-4}$;
$\bullet$ the kinematic cuts from (39);
$\bullet$ $W^{+}W^{-}$ decay into $\ell\,\nu$ and $q\bar{q}^{\prime}$ :
$W^{+}W^{-}\to e(\mu)\,\nu\;\;q\bar{q^{\prime}}$.

As a result, at LHC Run-3 option the expected number of events for
$t\bar{t}$-pair production with subsequent decays (with the following
$W^{+}W^{-}$ decays) are as follows:

$\displaystyle t\,\bar{t}\to
b\bar{b}bW^{+}\,\bar{b}W^{-}\quad\to\;\;2b\,2\bar{b}\,\ell\nu\,q\,\bar{q}^{\prime}$
$\displaystyle N(2b\,2\bar{b}\,\ell\nu\,q\,\bar{q}^{\prime})\approx 2000$ (44)

This number of events is large enough that the search for such rare decay
looks like a promising goal for the experimental searches.

##  8 Kinematics

In this section we study the obvious question: how one get a manifestation of
the charged Higgs boson contribution to this rare $t\to bW^{+}b\bar{b}$ decay
of the top-quark ?.

In this section , the symbol ”SM” corresponds to the case when top-quark
decays within the SM, while ”H” denotes the inclusion (exchange) of charged
Higgs boson:

$\displaystyle t_{SM}\to b\bar{b}b\,W^{+}$ $\displaystyle:$ (SM only) (45)
$\displaystyle t_{H}\to b\bar{b}b\,W^{+}$ $\displaystyle:$
$\displaystyle\quad\hbox{(SM}\;+\;H^{\pm})$ (46)

As it seen from the Table 3 for a small $\tan\beta$ values $(\tan\beta\leq
0.2)$ the decay width due to $H^{\pm}$ contribution exceeds the SM value by at
least one (or even two) order of magnitude.

Therefore, one could expect a significant increase in the number of events in
this case. Additional manifestation of the $H^{\pm}$ contribution can be seen
in different kinematic distributions. We have found one variable, whose
distribution demonstrates clear difference in decays (45) and (46). This is
the cosine of the angle between the 3-momenta of the final $b$-quarks in the
$t$-quark rest frame, $\cos\theta_{bb}$. These distributions are shown in the
Figs. 9 and 10.

Figure 9: The distribution over $\cos\theta_{bb}$ for two charged Higgs masses
and for several $\tan\beta$ values. The black curve corresponds to SM $t\to
b\,\bar{b}\,b\,W^{+}$ decay.

Figure 10: The same distributions as in Fig. 9, but for $M(H^{\pm})=175$ and
$180$ GeV.

As can be seen from these plots, for small $\tan\beta$ values ($\tan\beta\leq
0.2$) these distributions for $t_{H}$ decay (46) have a shape different from
the SM decay.

The next plots (see Fig. 11) present the same distributions, but evaluated for
top-quark production in $pp$-collisions at $\sqrt{s}=13$ TeV after application
all cuts from (39). It is seen an obvious difference in distributions over
$\cos\theta_{bb}$.

Figure 11: The distribution over $\cos\theta_{bb}$ for for-body top-quark
decay. Here the $t$-quark is produced in $pp$-collisions at $\sqrt{s}=13$ TeV.
The black curve corresponds to SM decay $t\to b\,\bar{b}\,b\,W^{+}$.

Certainly, this conclusion is nothing more than a demonstration of different
kinematic in (45) and (46) decays. More realistic simulation could destroy
this pictures and conclusions. Nevertheless, we believe that advanced methods
of the experimental data analysis will make it possible to search for this
rare $t$-quark four-body decay and (in principle) to detect manifestation of a
charged Higgs boson (if it exists) to this $t\to bW^{+}b\bar{b}$ decay.

##  9 Conclusion

In this article the rare top decay $t\to b\,\bar{b}\,b\,W^{+}$ is considered.
The relatively large partial decay width, $\Gamma(t\to
b\,\bar{b}\,b\,W^{+})=(9.30\pm 0.03)\times 10^{-4}$ GeV (and $\mathsf{Br}(t\to
b\,\bar{b}\,b\,W^{+})=6.29\times 10^{-4}$), makes possible the searches for
this rare $t$-quark decay at LHC.

The role of the charged Higgs boson (with the mass $M(H^{\pm})=(160\div 180)$
GeV and
$\tan\beta\raisebox{-3.0pt}{$\>\stackrel{{\scriptstyle<}}{{\scriptstyle\sim}}\>$}1$)
in top-quark decays is considered. Taking in account the current value of the
$t$-quark total decay width as well as the branching fractions of $t$-quark
decays, the allowable region of $\tan\beta$ versus $M(H^{\pm})$ is evaluated.

In this article we propose an additional method for $H^{\pm}$-boson searches.
This scenario explores only the increase of the partial width of the $t\to
b\,\bar{b}\,b\,W^{+}$ decay and does not require the appearance of the peak in
the invariant mass distributions around the mass of the charged Higgs boson.

It is shown that the $H^{\pm}$-boson contribution to $\Gamma(t\to
b\,\bar{b}\,b\,W^{+})$ with
$\tan\beta\raisebox{-3.0pt}{$\>\stackrel{{\scriptstyle<}}{{\scriptstyle\sim}}\>$}1$
could increase this width by up to two orders of magnitude (depending on
$\tan\beta$ values).

The role of the possible background processes (for $t$-quark decay as well as
for production) is investigated. The exploration of the $b$-tagging
($\varepsilon_{b}\approx 70\%$) provides sufficiently good suppression of all
other decay channels. The main background to $t\bar{t}$ production with
subsequent $t\to b\,\bar{b}\,b\,W^{+}$ decay comes from four heavy quarks
production ($pp\to t\bar{t}b\bar{b})$. This background also could be
suppressed by exploring of $b$-tagging and proper kinematic cuts.

As a result one may expect to get about 2000 events at LHC energies with total
luminosity of ${\mathscr{L}}_{tot}=300$ fb-1 for the SM $t$-quark decay $t\to
b\,\bar{b}\,b\,W^{+}(\to e/\mu\,\nu)$. Such a number of expected events makes
it realistic challenge for this study is in the LHC experiments.

## Acknowledgments

In conclusion the author is grateful to P.S. Mandrik, V.F. Obraztsov, R.N.
Rogalyov, and A.M. Zaitsev for multiple and fruitful discussions.

## References

  * [1] M. Beneke, et al., Top quark physics, in: Workshop on Standard Model Physics (and more) at the LHC (First Plenary Meeting), 2000, pp. 419–529 (3 2000).  arXiv:hep-ph/0003033.
  * [2] P. A. Zyla, et al., Review of Particle Physics, PTEP 2020 (8) (2020) 083C01, https://pdg.lbl.gov/ (2020).  doi:10.1093/ptep/ptaa104.
  * [3] V. D. Barger, R. J. N. Phillips, COLLIDER PHYSICS, Vol. 71, Addison-Wesley, 1987 (1987).
  * [4] G. C. Branco, P. M. Ferreira, L. Lavoura, M. N. Rebelo, M. Sher, J. P. Silva, Theory and phenomenology of two-Higgs-doublet models, Phys. Rept. 516 (2012) 1–102 (2012).  arXiv:1106.0034, doi:10.1016/j.physrep.2012.02.002.
  * [5] A. G. Akeroyd, et al., Prospects for charged Higgs searches at the LHC, Eur. Phys. J. C 77 (5) (2017) 276 (2017).  arXiv:1607.01320, doi:10.1140/epjc/s10052-017-4829-2.
  * [6] A. Arbey, F. Mahmoudi, O. Stal, T. Stefaniak, Status of the Charged Higgs Boson in Two Higgs Doublet Models, Eur. Phys. J. C 78 (3) (2018) 182 (2018).  arXiv:1706.07414, doi:10.1140/epjc/s10052-018-5651-1.
  * [7] G. Abbiendi, et al., Search for Charged Higgs bosons: Combined Results Using LEP Data, Eur. Phys. J. C 73 (2013) 2463 (2013).  arXiv:1301.6065, doi:10.1140/epjc/s10052-013-2463-1.
  * [8] T. Aaltonen, et al., Search for Higgs bosons predicted in two-Higgs-doublet models via decays to tau lepton pairs in 1.96-TeV p anti-p collisions, Phys. Rev. Lett. 103 (2009) 201801 (2009).  arXiv:0906.1014, doi:10.1103/PhysRevLett.103.201801.
  * [9] V. M. Abazov, et al., Search for Higgs bosons of the minimal supersymmetric standard model in $p\bar{p}$ collisions at $\sqrt{s}=1.96$ TeV, Phys. Lett. B 710 (2012) 569–577 (2012).  arXiv:1112.5431, doi:10.1016/j.physletb.2012.03.021.
  * [10] M. Aaboud, et al., Measurements of inclusive and differential fiducial cross-sections of $t\overline{t}$ production with additional heavy-flavour jets in proton-proton collisions at $\sqrt{s}$ = 13 TeV with the ATLAS detector, JHEP 04 (2019) 046 (2019).  arXiv:1811.12113, doi:10.1007/JHEP04(2019)046.
  * [11] M. Aaboud, et al., Search for charged Higgs bosons decaying via $H^{\pm}\to\tau^{\pm}\nu_{\tau}$ in the $\tau$+jets and $\tau$+lepton final states with 36 fb-1 of $pp$ collision data recorded at $\sqrt{s}=13$ TeV with the ATLAS experiment, JHEP 09 (2018) 139 (2018).  arXiv:1807.07915, doi:10.1007/JHEP09(2018)139.
  * [12] G. Aad, et al., Search for charged Higgs bosons decaying into a top quark and a bottom quark at $\sqrt{\mathrm{s}}$ = 13 TeV with the ATLAS detector, JHEP 06 (2021) 145 (2021).  arXiv:2102.10076, doi:10.1007/JHEP06(2021)145.
  * [13] A. M. Sirunyan, et al., Measurement of the $\mathrm{t\bar{t}}\mathrm{b\bar{b}}$ production cross section in the all-jet final state in pp collisions at $\sqrt{s}=$ 13 TeV, Phys. Lett. B 803 (2020) 135285 (2020).  arXiv:1909.05306, doi:10.1016/j.physletb.2020.135285.
  * [14] A. M. Sirunyan, et al., Search for a light charged Higgs boson in the H± $\to$ cs channel in proton-proton collisions at $\sqrt{s}=$ 13 TeV, Phys. Rev. D 102 (7) (2020) 072001 (2020).  arXiv:2005.08900, doi:10.1103/PhysRevD.102.072001.
  * [15] A. M. Sirunyan, et al., Search for a charged Higgs boson decaying into top and bottom quarks in events with electrons or muons in proton-proton collisions at $\sqrt{\mathrm{s}}$ = 13 TeV, JHEP 01 (2020) 096 (2020).  arXiv:1908.09206, doi:10.1007/JHEP01(2020)096.
  * [16] A. M. Sirunyan, et al., Search for charged Higgs bosons decaying into a top and a bottom quark in the all-jet final state of pp collisions at $\sqrt{s}$ = 13 TeV, JHEP 07 (2020) 126 (2020).  arXiv:2001.07763, doi:10.1007/JHEP07(2020)126.
  * [17] V. I. Borodulin, R. N. Rogalyov, S. R. Slabospitskii, CORE 3.1 (COmpendium of RElations, Version 3.1) (2 2017).  arXiv:1702.08246.
  * [18] S. R. Slabospitsky, L. Sonnenschein, TopReX generator (version 3.25): Short manual, Comput. Phys. Commun. 148 (2002) 87–102 (2002).  arXiv:hep-ph/0201292, doi:10.1016/S0010-4655(02)00471-X.
  * [19] A. M. Sirunyan, et al., Identification of heavy-flavour jets with the CMS detector in pp collisions at 13 TeV, JINST 13 (05) (2018) P05011 (2018).  arXiv:1712.07158, doi:10.1088/1748-0221/13/05/P05011.
  * [20] E. Ma, D. P. Roy, J. Wudka, Enhanced three-body decay of the charged Higgs boson, Phys. Rev. Lett. 80 (1998) 1162–1165 (1998).  arXiv:hep-ph/9710447, doi:10.1103/PhysRevLett.80.1162.
  * [21] R. Decker, M. Nowakowski, A. Pilaftsis, Dominant three-body decays of a heavy Higgs and top quark, Z. Phys. C 57 (1993) 339–348 (1993).  arXiv:hep-ph/9301283, doi:10.1007/BF01565067.
  * [22] G. Mahlon, S. J. Parke, Finite width effects in top quark decays, Phys. Lett. B 347 (1995) 394–398 (1995).  arXiv:hep-ph/9412250, doi:10.1016/0370-2693(95)00083-W.
  * [23] E. E. Jenkins, The Rare top decays $t\to bW^{+}Z$ and $t\to cW^{+}W^{-}$, Phys. Rev. D 56 (1997) 458–466 (1997).  arXiv:hep-ph/9612211, doi:10.1103/PhysRevD.56.458.
  * [24] G. Altarelli, L. Conti, V. Lubicz, The t —$>$ WZ b decay in the standard model: A Critical reanalysis, Phys. Lett. B 502 (2001) 125–132 (2001).  arXiv:hep-ph/0010090, doi:10.1016/S0370-2693(00)01333-2.
  * [25] A. Papaefstathiou, G. Tetlalmatzi-Xolocotzi, Rare top quark decays at a 100 TeV proton–proton collider: $t\rightarrow bWZ$ and $t\rightarrow hc$, Eur. Phys. J. C 78 (3) (2018) 214 (2018).  arXiv:1712.06332, doi:10.1140/epjc/s10052-018-5701-8.
  * [26] P. Onyisi, A. Webb, Impact of rare decays $t\to\ell^{\prime}\nu b\ell\ell$ and $t\to qq^{\prime}b\ell\ell$ on searches for top-associated physics, JHEP 02 (2018) 156 (2018).  arXiv:1704.07343, doi:10.1007/JHEP02(2018)156.
  * [27] N. Quintero, J. L. Diaz-Cruz, G. Lopez Castro, Lepton pair emission in the top quark decay $t\to bW^{+}\ell^{-}\ell^{+}$, Phys. Rev. D 89 (9) (2014) 093014 (2014).  arXiv:1403.3044, doi:10.1103/PhysRevD.89.093014.
  * [28] A. Bredenstein, A. Denner, S. Dittmaier, S. Pozzorini, NLO QCD corrections to t anti-t b anti-b production at the LHC: 1. Quark-antiquark annihilation, JHEP 08 (2008) 108 (2008).  arXiv:0807.1248, doi:10.1088/1126-6708/2008/08/108.
  * [29] F. Cascioli, P. Maierhöfer, N. Moretti, S. Pozzorini, F. Siegert, NLO matching for $t\bar{t}b\bar{b}$ production with massive $b$-quarks, Phys. Lett. B 734 (2014) 210–214 (2014).  arXiv:1309.5912, doi:10.1016/j.physletb.2014.05.040.
  * [30] M. V. Garzelli, A. Kardos, Z. Trócsányi, Hadroproduction of $t\bar{t}b\bar{b}$ final states at LHC: predictions at NLO accuracy matched with Parton Shower, JHEP 03 (2015) 083 (2015).  arXiv:1408.0266, doi:10.1007/JHEP03(2015)083.
  * [31] G. Bevilacqua, M. Czakon, C. G. Papadopoulos, M. Worek, Hadronic top-quark pair production in association with two jets at Next-to-Leading Order QCD, Phys. Rev. D 84 (2011) 114017 (2011).  arXiv:1108.2851, doi:10.1103/PhysRevD.84.114017.
  * [32] G. Bevilacqua, M. Worek, On the ratio of $t\overline{t}b\overline{b}$ and $t\overline{t}jj$ cross sections at the CERN Large Hadron Collider, JHEP 07 (2014) 135 (2014).  arXiv:1403.2046, doi:10.1007/JHEP07(2014)135.
  * [33] G. Bevilacqua, M. V. Garzelli, A. Kardos, $t\bar{t}b\bar{b}$ hadroproduction with massive bottom quarks with PowHel (9 2017).  arXiv:1709.06915.
  * [34] T. Ježo, J. M. Lindert, N. Moretti, S. Pozzorini, New NLOPS predictions for $t\bar{t}+b$ -jet production at the LHC, Eur. Phys. J. C 78 (6) (2018) 502 (2018).  arXiv:1802.00426, doi:10.1140/epjc/s10052-018-5956-0.
  * [35] A. M. Sirunyan, et al., Measurement of the cross section for $\text{t}\bar{\text{t}}$ production with additional jets and b jets in pp collisions at $\sqrt{s}=$ 13 TeV, JHEP 07 (2020) 125 (2020).  arXiv:2003.06467, doi:10.1007/JHEP07(2020)125.
  * [36] S. Grancagnolo, Top quark pair-production cross section measurements at LHC, PoS DIS2019 (2019) 152 (2019).  doi:10.22323/1.352.0152.