# Quasi-crystalline order in vibrated granular matter

A. Plati  Université Paris-Saclay, CNRS, Laboratoire de Physique des Solides,
91405 Orsay, France  R. Maire  Université Paris-Saclay, CNRS, Laboratoire de
Physique des Solides, 91405 Orsay, France  E. Fayen  Université Paris-Saclay,
CNRS, Laboratoire de Physique des Solides, 91405 Orsay, France  F. Boulogne
Université Paris-Saclay, CNRS, Laboratoire de Physique des Solides, 91405
Orsay, France  F. Restagno  Université Paris-Saclay, CNRS, Laboratoire de
Physique des Solides, 91405 Orsay, France  F. Smallenburg  Université Paris-
Saclay, CNRS, Laboratoire de Physique des Solides, 91405 Orsay, France  G.
Foffi  Université Paris-Saclay, CNRS, Laboratoire de Physique des Solides,
91405 Orsay, France

###### Abstract

Quasi-crystals are aperiodic structures that present crystallographic
properties which are not compatible with that of a single unit cell. Their
revolutionary discovery in a metallic alloy, less than three decades ago, has
required a full reconsideration of what we defined as a crystal structure.
Surprisingly, quasi-crystalline structures have been discovered also at much
larger length scales in different microscopic systems for which the size of
elementary building blocks ranges between the nanometric to the micrometric
scale. Here, we report the first experimental observation of spontaneous
quasi-crystal self-assembly at the millimetric scale. This result is obtained
in a fully athermal system of macroscopic spherical grains vibrated on a
substrate. Starting from a liquid-like disordered phase, the grains begin to
locally arrange into three types of squared and triangular tiles that
eventually align, forming 8-fold symmetric quasi-crystal that has been
predicted in simulation but not yet observed experimentally in non-atomic
systems. These results are not only the proof of a novel route to
spontaneously assemble quasi-crystals but are of fundamental interest for the
connection between equilibrium and non-equilibrium statistical physics.

In 1982, Shechtman discovered the first alloy with a diffraction pattern for
which the Bragg peaks showed a symmetry that is forbidden by crystallography
in periodic solids [1]. This discovery was initially met with resistance: the
existence of structures in which atoms can be arranged in spatial structures
which lack long-range periodicity, while still preserving sufficient long-
range order to generate discrete Bragg peaks clashed with the elegant picture
of crystals as consisting of a repeating unit cell. Nonetheless, eventually
materials with this property, which were called quasi-crystals (QCs) [2],
changed the way in which scientists interpret the crystal state, by
disentangling the concept of order from the concept of periodicity, to the
point where the very definition of crystals had to be changed to include
aperiodic structures [3].

In the following years, quasi-crystalline structures have been observed in
several artificial alloys (see Refs. 4, 5 for a review) and were even
discovered in a natural occurring mineral, Icosahedrite [6], of probable
extraterrestrial origin [7]. More recently, quasi-crystals have been also
observed at much larger length scales in a wide range of soft matter systems
[8, 9, 10, 11, 12, 13, 14, 15, 16, Förster2013] revealing promising optical
properties for next-generation photonic devices [17, 18]. In soft quasi-
crystals [14], two fundamental questions arise: i) understanding up to which
length scales we can observe spontaneous quasi-crystalline order and ii)
identifying the key dynamical and interaction properties required for a soft-
matter system to form a quasi-crystal. For the former, the self-assembly of
QCs has been observed on length scales that are related to the nature of the
elementary building blocks that range from macromolecular structures [19] and
nanoparticles [10] to polymers aggregates made of micelles [12, 20], passing
by polymers [8, 13, 15, 16] and thin films [Förster2013]. To our knowledge,
the largest soft-matter quasi-crystals found to self-assemble consisted of
micrometer-size micelles [20, 21].

The second question is of a more fundamental nature and has been explored
extensively using numerical simulations. In particular, simple coarse-grained
models of interacting particles make it possible to simulate systems large
enough to display aperiodicity and to pin down its origin, from anisotropic
repulsive [22] and attractive [23] particles to simple isotropic potentials
[24, 25, 26] and hard spheres [27]. Interestingly, one of the simplest systems
leading to quasi-crystal formation in silico is a simple two-dimensional
binary mixture of non-additive hard disks undergoing equilibrium dynamics [28,
27]. This proved that, in the right geometrical conditions, quasi-crystalline
order can emerge from a purely entropy-driven self-assembly process. To our
knowledge, quasi-crystal self-assembly in dissipative non-equilibrium systems
where energy is constantly supplied from the environment and internally
dissipated (such as active or granular matter) is still unexplored.

A natural avenue to explore quasi-crystal formation beyond colloidal scale and
thermodynamic equilibrium is to use granular matter, which has proven itself
to be an ideal playground for the exploration of non-equilibrium phenomena
over the last three decades. Depending on how a granular system is confined
and the imposed external driving, it may exhibit either a fluid-like or a
solid-like behaviour [29, 30] and can undergo a variety of phase transitions
[31, 32, 33, 34, 35, 36]. However, spontaneous quasi-crystal formation in
systems driven out of thermodynamic equilibrium has not yet been observed in
either experiments or simulations.

In this paper, we report the experimental and numerical observation of quasi-
crystalline order in a binary mixture of millimeter-size spherical grains
vibrated on a substrate. Our findings demonstrate that quasi-crystals can be
formed also by macroscopic particles much beyond the scale at which thermal
agitation plays a role. Our system is indeed intrinsically out of equilibrium
due to dissipation arising from frictional forces and energy injection due to
external driving.

Figure 1:  Numerical results obtained with EDMD of the collisional model
defined in Eqs. 1. Here $\\{q,x_{S},\phi\\}=\\{0.5,0.68,0.85\\}$,
$N_{S}+N_{L}=5000$, $\alpha=0.95$, $\Delta=0.02$. a) Sketch of the granular
non-additive hard-core interactions, final configuration, reconstructed tiling
and structure factor. b) Bond orientation histrogram. c) Histogram of the area
fraction occupied by tiles of different types and orientations. The color of
the tiles is chosen according to their orientation but we also consider
defected regions which are not covered by the correct square-triangle tiling
(white areas) and tiles with ambiguous orientations or misaligned with the
dominant set of bond directions (yellow areas). Results are obtained with a
$5.3\times 10^{10}$ collision-long simulation.

In the following, we first discuss granular quasi-crystal formation in
numerical simulations of a coarse-grained collisional model. We then report
the main result of this paper namely the experimental self-assembly of a
granular quasi-crystal.

In both simulations and experiments, we consider a binary mixture of spherical
grains ($N_{S}$ with diameter $\sigma_{S}$ and $N_{L}$ with diameter
$\sigma_{L}$) lying on a horizontal substrate. In order to characterize the
geometrical properties of such a mixture, we use the following dimensionless
parameters: the size ratio $q=\sigma_{S}/\sigma_{L}$, the fraction of small
grains $x_{S}=N_{S}/(N_{S}+N_{L})$ and the area fraction
$\phi=(N_{S}\sigma_{S}^{2}+N_{L}\sigma_{L}^{2})\pi/4L^{2}$, where $L$ is the
side of the box. As shown in Fig. 1a, spheres on a substrate can be mapped
onto disks thanks to the introduction of a non-additivity parameter
$\delta=2\sqrt{q}(1+q)-1$ and letting them interact at a distance smaller than
$\sigma_{LS}=\frac{1}{2}(\sigma_{S}+\sigma_{L})(1+\delta)$, note that
$-1<\delta<0$. This is the main idea underlying the effective 2D model
considered in a previous numerical study focusing on elastic non-additive hard
disks following dynamics at thermodynamic equilibrium [27]. This previous
study revealed that, for sufficiently high area fractions and depending on the
specific $\\{q,x_{S}\\}$ combination, one can observe the self-assembly of
different crystals including a 12-fold and a 8-fold symmetric quasi-crystal.
These results are in agreement with the fact that for conservative dynamics,
geometrical constraints and hard-core interactions can be minimal ingredients
for thermodynamic stability of quasi-crystalline structures [37, 38]. In this
work, we focus our attention on the quasi-crystalline 8-fold symmetric phase
(QC8), which was observed to self-assemble significantly more rapidly than the
dodecagonal quasi-crystal [27].

The collisional model we use in our simulations represents the athermal/non-
equilibrium extension of the one considered in [27]. It has been shown to
embody dissipative and forcing mechanisms of experimental systems where
spherical grains are placed on vertically vibrating horizontal substrates [39,
40, 41]. In these systems, energy is injected along the vertical direction
through grain-substrate collisions and then transferred to the horizontal ones
through grain-grain collisions with an efficiency that depends on the impact
kinematics. In the model the dynamics is fully 2D: there is no vibrating plate
since its effect is coarse-grained out thanks to the introduction of
instantaneous grain-grain collisions which take into account both energy
injection and dissipation. In this model, a binary collision between grains of
mass $m_{i}$, $m_{j}$ obey the following rule for the velocity
($\mathbf{v}_{i},\mathbf{v}_{j}$) update:

$\begin{split}\mathbf{v}_{i}^{\prime}&=\mathbf{v}_{i}+\frac{m_{j}(1+\alpha)}{m_{i}+m_{j}}(\mathbf{v}_{ij}\cdot\boldsymbol{\hat{\sigma}}_{ij})\boldsymbol{\hat{\sigma}}_{ij}+\frac{2m_{j}}{m_{i}+m_{j}}\Delta\boldsymbol{\hat{\sigma}}_{ij}\\\
\mathbf{v}_{j}^{\prime}&=\mathbf{v}_{j}-\frac{m_{i}(1+\alpha)}{m_{i}+m_{j}}(\mathbf{v}_{ij}\cdot\boldsymbol{\hat{\sigma}}_{ij})\boldsymbol{\hat{\sigma}}_{ij}-\frac{2m_{i}}{m_{i}+m_{j}}\Delta\boldsymbol{\hat{\sigma}}_{ij},\end{split}$
(1)

where the primed letters refer to post-collisional variables, and
$\boldsymbol{\hat{\sigma}}_{ij}$ and $\mathbf{v}_{ij}$ are respectively the
unit vector joining particles $i$ and $j$ and the relative velocity between
them. The parameter $\alpha$ is the coefficient of restitution that embodies,
for $0\leq\alpha<1$, the dissipative nature of the collision. The last term in
the Eqs. 1, is responsible, via the parameter $\Delta$, for the velocity gain
arising from the non-planar collisions which are coarse-grained out in the
effective 2D description (see SI for a more detailed explanation). By
computing the energy change in a collision [39], it is possible to see that,
depending on the impact kinematics, one can have conditions in which the total
energy decreases or increases. In this simple granular model, the limit to the
equilibrium conservative case is recovered by setting $\alpha=1$ and
$\Delta=0$.

A granular system cannot attain thermodynamic equilibrium but it can reach a
non-equilibrium steady state thanks to a balance between injected and
dissipated energy. Nevertheless, one can often identify a conservative system
with the same geometrical properties and consider it as an equilibrium
counterpart. In this perspective, non-additive hard disks in the conservative
limit ($\alpha=1$ and $\Delta=0$) represent the equilibrium version of our
granular system. Of course, given the dissipative/athermal nature of the
dynamics, one cannot expect theoretical or numerical predictions for the
equilibrium counterpart to hold in the granular case. However, in some
specific conditions, vibrated granular materials have been shown to exhibit an
equilibrium-like phenomenology, such as in the case of tracer diffusion in
granular gases [42, 43] or hexagonal crystal formation in monodisperse
granular layers vibrated on a substrate [32, 33, 36]. The latter phenomenon is
particularly relevant for our study since many aspects of the equivalent
elastic hard-sphere crystallization have been observed also in the liquid-
solid granular phase transition. The main challenging aspect for the
observation of an equilibrium-like behaviour in our system is given by size
polydispersity because it usually triggers non-equilibrium effects in vibrated
granular materials. Size segregation [34, 44, 45, 46] and violation of energy
equipartition [47], are two examples of that. Thus, an important question
underlying the approach we propose is the following: can we tune the non-
equilibrium parameters i.e. related to forcing/dissipative mechanisms of the
granular system, such that it self-assembles into the quasi-crystalline
structures that have been observed in the equilibrium counterpart?

Results obtained through event-driven molecular dynamics (EDMD) simulations of
the model described by Eqs. 1 are reported in Fig. 1a. There we show the last
snapshot of a simulation and the associated scattering pattern computed from
the large grain positions. The final self-assembled structure exhibits no
periodic order but the scattering peaks reveal an underlying 8-fold symmetry.
This particular symmetry is forbidden by ordinary crystallographic order based
on repeated translations of a single unit cell. Indeed, in our system, the
final strcture can be decomposed into a combination of three different tiles
each one appearing with different orientations. In Fig. 1a, we also highlight
how small and large grains combine to form such tiles: we have i) small
squares made of four large grains surrounding one small grain, ii) isosceles
triangles made of three large grains surrounding one small grain and iii)
large squares made of four large grains surrounding a square of small grains.
The sides of the tiles coincide with bonds between large grain nearest
neighbours and we can identify long and short bonds. The former outlines large
square sides and triangle legs, the latter forms small square sides and
triangle bases. In addition to the scattering pattern, another piece of
evidence of the 8-fold symmetry is given by the histogram of bond orientations
over the entire system (Fig. 1b). Here we can clearly see that short and long
bonds are not uniformly distributed but are aligned along eight specific
directions. These dominant directions are spontaneously selected among a
continuum since PBC do not impose preferred orientations. Within a specific
8-fold set of bond directions, small and large squared tiles can only appear
with two specific orientations while triangular tiles can lay along eight
ones. One can then classify all the tiles in the system (see Methods) as shown
in Fig. 1c and reconstruct the overall tiling of the plane (bottom-left side
of Fig. 1a). The observed square-triangle tiling can cover aperiodically an
infinite plane with long-range 8-fold orientational order. What we observe in
our simulations is a finite portion of such an infinite quasi-crystal, the
fact that tiles with the same form but with different orientations cover a
similar area fraction is coherent with this picture [27]. Additional EDMD
simulations varying the geometrical parameters $\\{\phi,x_{S},q\\}$, the non-
equilibrium ones $\\{\alpha,\Delta\\}$ and considering different system sizes
confirmed the robustness of QC8 formation for this model (see SI).

We now turn our attention to the experimental realization. Our setup consists
of non-magnetic steel spherical grains confined in a quasi-2D squared
container (height $h$ and width $L\gg h$) which is vertically vibrated by an
electrodynamic shaker following a signal $z_{p}(t)$. The evolution of the
system is followed by a camera that allows to detect the horizontal position
of the grains ($xy$ coordinates) as a function of time. The relevant set of
explored geometrical parameters $\\{q,x_{S},\phi\\}$ for the experiments is
chosen close to those used in the simulation model. The realistic system
eventually attains a non-equilibrium steady state whose dynamical properties
sensitively depend on the driving force. To tune this, we performed
preliminary numerical simulations of the setup implemented through the
Discrete Element Method (DEM) [48, 49]. Such simulations implement granular
dynamics by means of accurate contact models which allow for studying in-
silico prototypes of realistic setups (see SI). From this analysis and
subsequent experimental tests (see SI), we found that to observe QC8 self-
assembly one generally needs sufficiently strong vibrations to allow for an
efficient vertical-to-horizontal energy transfer but not so strong that grains
pile up on each other (which alters the effective 2D packing fraction of the
system).

Figure 2:  Experimental results for
$\\{q,x_{S},\phi\\}=\\{0.5,0.675,0.848\\}$, $N_{S}+N_{L}=3840$,
$\sigma_{L}=2.381$ mm, $z_{p}(t)=A\sin(2\pi ft)$ where $A=22$ $\mu m$ and
$f=120$ Hz. Data from real-time measurements during QC8 formation (see also
supplemental video). a) Area fraction of tiles correctly aligned according to
the dominant orientational order, misaligned tiles and defects. Reconstructed
tiling (b), structure factor (c) and bond orientation histogram for early,
intermediate and late times. Yellow tiles are misaligned with the final long-
range order, green ones are in agreement with it. We note that the latter have
long bonds laying on the $xy$ axis. See also Movie 1.  Figure 3:  a)
Configuration, structure factor, reconstructed tiling obtained after $\sim
170$ hours of the experimental run already shown in Fig. 2. The specific
snapshot shown is the one with the highest QC8 area fraction (red arrow in
Fig. 2). We also show a sketch of the $xz$ projection of the setup and a
comparison between real space tiling and the reconstructed one in a zoomed
region on the sample. b) Bond orientation histogram. c) Tile area fraction
histogram with snapshot of single granular tiles forming the quasi-crystal.

For suitable choices of the driving parameters, our experimental system indeed
spontaneously forms the 8-fold quasi-crystal. The dynamics of this process is
shown in Fig. 2a (see also supplemental video) where we plot the occupied area
fraction as a function of time for three different groups of tiles: the ones
oriented according to the specific 8-fold set of directions which dominates at
the end of the experiment (green), misaligned tiles (yellow) and defects
(grey). These results are obtained from a seven-day-long experiment with
sinusoidal vibrations and quantify the two distinct processes contributing to
the QC8 self-assembly: single tile formation and tile orientational ordering.
We note that both aligned and misaligned tiles increase monotonously at very
short times. After that, their evolution becomes non-monotonous: aligned tiles
exhibit an increase interrupted by sudden drops while misaligned tiles exhibit
a decrease interrupted by sudden growths. Both of them eventually reach a
final plateau. Such behaviour reveals that, once tiles are formed and locally
ordered, they still need substantial rearrangements in order to form larger
quasi-crystalline domains and this makes the QC formation extremely slow.

From the evolution of reconstructed tilings shown in Fig. 2b, we observe that
QC8 with the dominant 8-fold symmetry is growing from the boundaries toward
the bulk. During this process, it is possible to observe the coexistence
between different misaligned QC domains ($t\sim t_{\text{max}}/2$). In Fig.
2c, we also show the evolution of the structure factor: the initial
configuration shows the typical rings of a liquid-like structure, the middle
one exhibits blurred peaks originating from the coexistence of multiply
oriented QC domains while the last one presents sharper peaks highlighting
long-range order according to a dominant set of 8 directions. Finally, from
the evolution of the bond orientation histogram (Fig. 2c), we can see that the
most favoured set of orientations is the one with large tile sides aligned
with the hard walls. Repetitions of the experiment revealed that this is a
reproducible feature. Indeed, an important difference with respect to the
numerical case discussed above is represented by the hard horizontal walls
which break the orientational symmetry. The enforced grain alignment at the
$xy$ boundaries favours two specific sets of 8-fold symmetric directions: one
with short bonds and one with long bonds aligned with walls. The dominance of
the latter can be explained considering that the forming quasi-crystalline
structure requires much more long bonds than short ones.

An interesting picture seems to emerge: tiles of the desired shapes form very
rapidly, probably as a consequence of their efficient local packing. However,
global alignment requires much more collective rearrangements that appear to
be rare events. It is important to point out that real-time measurements of
forming quasi-crystalline structures are extremely rare in experiments.
Generally, the techniques used for QC self-assembly (e.g. evaporation of
nanoparticle solutions) are not compatible with the observation of the
dynamics during the ordering process but only allows to analyze the final
structure. The experimental study of system configurations over time is then
another important novelty of macroscopic quasi-crystals since it can shed
light on unexplored dynamical properties of QC self-assembly.

In Fig. 3a, we show the spatial configuration, the related scattering pattern,
and the reconstructed tiling obtained from the experimental configuration with
the largest quasi-crystalline area fraction (see the red arrow in Fig. 2a). On
the right side of the same figure, we report the histograms of bond
orientation (b) and tile area fraction (c). We note that the key properties of
the quasi-crystalline structure observed in EDMD simulations are found also in
the experiment: we have the same 8-fold long-range orientational order,
confirmed by the scattering intensity, and the same square-triangle tiling. By
comparing Fig. 1 and 3, we realize that the main difference between EDMD
simulations and experiments is the presence of larger defects in the latter.
It is difficult to pin down the origin of this difference as the $z$-to-$xy$
energy transfer mechanism of the real experimental model depends, in a highly
non-trivial way, on the dissipative properties of the materials and the plate
roughness, while the energy injection in the EDMD simulations is fully
described by a simple collision rule that depends only on two parameters. We
argue that larger defects are also the reason why we observe lower large-
square area fraction in the experiment with respect to EDMD (compare Fig. 1c
and 3c). In fact, large squares are favoured by a large number of small grains
[27], but many of them are stuck in defects that appear mainly as clusters of
small grains.

In order to test robustness and reproducibility of the discussed results we
repeated the experiment for different initial random configurations, different
drivings (sawtooth vibrations) and slightly shifted state points. As reported
in the SI, the emergence of quasi-crystalline order proved to be a quite
robust feature of this granular system.

More than 30 years ago, the quasi-crystal self-assembly paradigm emerged in
purely atomistic systems but was later extended to the nanometric and
micrometric scales typical of soft matter. Here, we have shown that this idea
can be extended to much larger length scales in vibrofluidized granular
systems. Our study reports the first observation of quasi-crystalline order in
a physical system undergoing athermal dynamics in the visible scale where
real-time measurements of system configuration can be performed during the
self-assembly process. The passage from the micrometer to the millimeter scale
is not a trivial one, since in the latter the grains do not undergo thermal
agitation. This means that for granular systems, we cannot explain the
emergence of quasi-crystalline order based on entropy arguments alone. Despite
this fundamental difference, quasi-crystalline order emerges in the same
conditions as predicted at equilibrium, something that could not be trivially
expected.

## Acknoledgments

We thank Andrea Puglisi and Andrea Gnoli for their invaluable help in setting
up this project and for their comments on the manuscript. We also thank
Marianne Impéror-Clerc, Laura Filion, Anuradha Jagannathan and Franscesco
Sciortino for carefully reading and commenting on our paper and Stéphane
Cabaret for the design of the quasi-2D cell. This work has been done with the
support of Investissements d’Avenir of LabEx PALM (Grant No.
ANR-10-LABX-0039-PALM) and of the Agence Nationale de la Recherche (ANR),
grant ANR-18-CE09-0025.

## References

  * Shechtman and Blech [1985] D. Shechtman and I. A. Blech, Metallurgical Transactions A 16, 1005 (1985).
  * Levine and Steinhardt [1984] D. Levine and P. J. Steinhardt, Phys. Rev. Lett. 53, 2477 (1984).
  * IUC [1992] Acta Crystallographica Section A 48, 922 (1992).
  * DiVincenzo and Steinhardt [1999] D. DiVincenzo and P. Steinhardt,  _Quasicrystals: The State of the Art_, Series on directions in condensed matter physics (World Scientific, 1999).
  * Ranganathan and Chattopadhyay [1991] S. Ranganathan and K. Chattopadhyay, Annual Review of Materials Science 21, 437 (1991).
  * Bindi  _et al._ [2009] L. Bindi, P. J. Steinhardt, N. Yao, and P. J. Lu, Science 324, 1306 (2009).
  * Bindi  _et al._ [2012] L. Bindi, J. M. Eiler, Y. Guan, L. S. Hollister, G. MacPherson, P. J. Steinhardt, and N. Yao, Proc. Natl. Acad. Sci. U.S.A. 109, 1396 (2012).
  * Hayashida  _et al._ [2007] K. Hayashida, T. Dotera, A. Takano, and Y. Matsushita, Phys. Rev. Lett. 98, 195502 (2007).
  * Barkan  _et al._ [2011] K. Barkan, H. Diamant, and R. Lifshitz, Phys. Rev. B 83, 172201 (2011).
  * Talapin  _et al._ [2009] D. V. Talapin, E. V. Shevchenko, M. I. Bodnarchuk, X. Ye, J. Chen, and C. B. Murray, Nature 461, 964 (2009).
  * Takano  _et al._ [2005] A. Takano, W. Kawashima, A. Noro, Y. Isono, N. Tanaka, T. Dotera, and Y. Matsushita, J. Polym. Sci., Part B: Polym. Phys. 43, 2427 (2005).
  * Zeng  _et al._ [2004] X. Zeng, G. Ungar, Y. Liu, V. Percec, A. E. Dulcey, and J. K. Hobbs, Nature 428, 157 (2004).
  * Zhang and Bates [2012] J. Zhang and F. S. Bates, J. Am. Chem. Soc. 134, 7636 (2012).
  * Lifshitz and Diamant [2007] R. Lifshitz and H. Diamant, Philos. Mag. 87, 3021 (2007).
  * Lee  _et al._ [2010] S. Lee, M. J. Bluemle, and F. S. Bates, Science 330, 349 (2010).
  * Wasio  _et al._ [2014] N. A. Wasio, R. C. Quardokus, R. P. Forrest, C. S. Lent, S. A. Corcelli, J. A. Christie, K. W. Henderson, and S. A. Kandel, Nature 507, 86 (2014).
  * Jin  _et al._ [1999] C. Jin, B. Cheng, B. Man, Z. Li, D. Zhang, S. Ban, and B. Sun, Appl. Phys. Lett. 75, 1848 (1999).
  * Zoorob  _et al._ [2000] M. E. Zoorob, M. D. B. Charlton, G. J. Parker, J. J. Baumberg, and M. C. Netti, Nature 404, 740 (2000).
  * Liu  _et al._ [2022] Y. Liu, T. Liu, X.-Y. Yan, Q.-Y. Guo, H. Lei, Z. Huang, R. Zhang, Y. Wang, J. Wang, F. Liu, F.-G. Bian, E. Meijer, T. Aida, M. Huang, and S. Z. Cheng, Proc. Natl. Acad. Sci. U.S.A. 119, e2115304119 (2022).
  * Fischer  _et al._ [2011] S. Fischer, A. Exner, K. Zielske, J. Perlich, S. Deloudi, W. Steurer, P. Lindner, and S. Förster, Proc. Natl. Acad. Sci. U.S.A. 108, 1810 (2011).
  * Xiao  _et al._ [2012] C. Xiao, N. Fujita, K. Miyasaka, Y. Sakamoto, and O. Terasaki, Nature 487, 349 (2012).
  * Haji-Akbari  _et al._ [2009a] A. Haji-Akbari, M. Engel, A. S. Keys, X. Zheng, R. G. Petschek, P. Palffy-Muhoray, and S. C. Glotzer, Nature 462, 773 (2009a).
  * Noya  _et al._ [2021] E. G. Noya, C. K. Wong, P. Llombart, and J. P. K. Doye, Nature 596, 367 (2021).
  * Widom  _et al._ [1987] M. Widom, K. J. Strandburg, and R. H. Swendsen, Phys. Rev. Lett. 58, 706 (1987).
  * Dotera  _et al._ [2014] T. Dotera, T. Oshiro, and P. Ziherl, Nature 506, 208 (2014).
  * Malescio and Sciortino [2022] G. Malescio and F. Sciortino, J. Mol. Liq. 349, 118209 (2022).
  * Fayen  _et al._ [2023a] E. Fayen, M. Impéror-Clerc, L. Filion, G. Foffi, and F. Smallenburg, Soft Matter 19, 2654 (2023a).
  * Fayen  _et al._ [2020] E. Fayen, A. Jagannathan, G. Foffi, and F. Smallenburg, J. Chem. Phys. 152, 204901 (2020).
  * Jaeger  _et al._ [1996] H. M. Jaeger, S. R. Nagel, and R. P. Behringer, Rev. Mod. Phys. 68, 1259 (1996).
  * B. Andreotti and Pouliquen [2013] Y. F. B. Andreotti and O. Pouliquen,  _Granular Media. Between Fluid and Solid_ (Cambridge University Press, Cambridge, 2013).
  * Eshuis  _et al._ [2007] P. Eshuis, K. van der Weele, D. van der Meer, R. Bos, and D. Lohse, Phys. Fluids 19, 123301 (2007).
  * Olafsen and Urbach [1998] J. S. Olafsen and J. S. Urbach, Phys. Rev. Lett. 81, 4369 (1998).
  * Reis  _et al._ [2006] P. M. Reis, R. A. Ingale, and M. D. Shattuck, Phys. Rev. Lett. 96, 258001 (2006).
  * Aranson and Tsimring [2006] I. S. Aranson and L. S. Tsimring, Rev. Mod. Phys. 78, 641 (2006).
  * Panaitescu  _et al._ [2012] A. Panaitescu, K. A. Reddy, and A. Kudrolli, Phys. Rev. Lett. 108, 108001 (2012).
  * Komatsu and Tanaka [2015] Y. Komatsu and H. Tanaka, Phys. Rev. X 5, 031025 (2015).
  * Haji-Akbari  _et al._ [2009b] A. Haji-Akbari, M. Engel, A. S. Keys, X. Zheng, R. G. Petschek, P. Palffy-Muhoray, and S. C. Glotzer, Nature 462, 773 (2009b).
  * Fayen  _et al._ [2023b] E. Fayen, L. Filion, G. Foffi, and F. Smallenburg, “A hard-sphere quasicrystal stabilized by configurational entropy,” (2023b), arXiv:2306.03549 [cond-mat.soft] .
  * Brito  _et al._ [2013] R. Brito, D. Risso, and R. Soto, Phys. Rev. E 87, 022209 (2013).
  * Brito  _et al._ [2020] R. Brito, R. Soto, and V. Garzó, Phys. Rev. E 102, 052904 (2020).
  * Garzó  _et al._ [2021] V. Garzó, R. Brito, and R. Soto, Physics of Fluids 33, 023310 (2021).
  * D’Anna  _et al._ [2003] G. D’Anna, P. Mayor, A. Barrat, V. Loreto, and F. Nori, Nature 424, 909 (2003).
  * Sarracino  _et al._ [2010] A. Sarracino, D. Villamaina, G. Costantini, and A. Puglisi, J. Stat. Mech. 04013, P04013 (2010).
  * Aumaître  _et al._ [2001] S. Aumaître, C. A. Kruelle, and I. Rehberg, Phys. Rev. E 64, 041305 (2001).
  * Baldassarri  _et al._ [2015] A. Baldassarri, A. Puglisi, and A. Sarracino, C. R. Phys. 16, 291 (2015).
  * Kudrolli [2004] A. Kudrolli, Rep. Prog. Phys. 67, 209 (2004).
  * Feitosa and Menon [2002] K. Feitosa and N. Menon, Phys. Rev. Lett. 88, 198301 (2002).
  * Cundall and Strack [1979] P. A. Cundall and O. D. L. Strack, Géotechnique 29, 47 (1979).
  * Pöschel and Schwager [2005] T. Pöschel and T. Schwager,  _Computational Granular Dynamics_ (Springer, Berlin, 2005).


# One-loop corrections to light neutrino masses in gauged U(1) extensions of
the standard model

Sho Iwamoto<EMAIL_ADDRESS>Institute for Theoretical Physics, ELTE
Eötvös Loránd University, Pázmány Péter sétány 1/A, 1117 Budapest, Hungary
Timo J. Kärkkäinen<EMAIL_ADDRESS>Institute for Theoretical
Physics, ELTE Eötvös Loránd University, Pázmány Péter sétány 1/A, 1117
Budapest, Hungary  Zoltán Péli<EMAIL_ADDRESS>ELKH-DE Particle
Physics Research Group, 4010 Debrecen, PO Box 105, Hungary  Zoltán Trócsányi
<EMAIL_ADDRESS>Institute for Theoretical Physics, ELTE Eötvös Loránd
University, Pázmány Péter sétány 1/A, 1117 Budapest, Hungary
ELKH-DE Particle Physics Research Group, 4010 Debrecen, PO Box 105, Hungary

###### Abstract

We consider gauged U(1) extensions of the standard model of particle physics
with three right-handed sterile neutrinos and a singlet scalar. The neutrinos
obtain mass via the type I seesaw mechanism. We compute the one loop
corrections to the elements of the tree level mass matrix of the light
neutrinos and show explicitly the cancellation of the gauge dependent terms.
We present a general formula for the gauge independent, finite one-loop
corrections for arbitrary number new U(1) groups, new complex scalars and
sterile neutrinos. We estimate the size of the corrections relative to the
tree level mass matrix in a particular extension, the super-weak model.

gauged U(1) extension, neutrino mass, seesaw mechanism, gauge independence

##  I Introduction

The standard model (SM) of particle interactions is one of the most successful
physics models with unprecedented precision for predicting physical
quantities, for instance for the anomalous magnetic moment of the electron.
However, it does not contain right handed neutrinos as they are sterile under
the SM gauge group. This in turn leads to the prediction that neutrinos are
massless which is in conflict of the now well established experimental result,
that at least two neutrinos are massive [1, 2], and therefore, signals that
the SM requires an extension to explain the origin of the neutrino masses.
There are lots of models attempting to explain neutrino masses. Among those
perhaps the most economical one that requires the least extension of the SM,
is the type I seesaw mechanism where neutrinos acquire masses after
spontaneous symmetry breaking (SSB) of one or more scalar fields [3, 4, 5, 6,
7, 8, 9, 10].

Recently there has been a lot of interest in gauged U(1)-extended models in
particle physics phenomenology, motivated by the observed difference between
the measured and SM predicted values of the anomalous magnetic moment of the
muon [11] and also anomalies in short-baseline neutrino oscillations [12].
Gauged $B-L$, $B-3L_{\tau}$, $L_{e}-L_{\mu}$ and $L_{\mu}-L_{\tau}$ have been
considered [13, 14, 15, 16, 17, 18, 19], as well as a general gauged U(1) not
related to flavour [20]. In these models, both seesaw and radiative one-loop
neutrino mass generation mechanisms have been considered.

As the effects of new physics are typically much smaller than those of the SM
interactions, computations in theories beyond the SM are often considered only
at tree level. Yet, the loop corrections may be sizable and can affect
significantly the validity region in the parameter space of the model. For
instance, the lightness of active neutrinos requires that the loop corrections
to the mass matrix of those particles must also be small in order to have a
phenomenologically viable model. Computations of such one-loop corrections
have been carried out previously in Refs. [21, 22] for the canonical seesaw
case, and in the context of multi-Higgs doublet models [23, 24, 25, 26]. In
the cases of gauged U(1) models we are not aware of a computation of the one-
loop corrections to active neutrino mass matrix.

In this article we consider gauged U(1) extensions of the SM and derive a
general formula for the one-loop corrections of the mass matrix of the active
neutrinos. The mass matrix of the active neutrinos emerges after SSB due to
the type I seesaw mechanism. Our goal is to derive the one-loop corrections to
that mass matrix and estimate their sizes relative to the tree level for a
particular example, called the super-weak force [27]. The super-weak model
contains three additional right-handed sterile (under the SM interactions)
neutrinos and one complex scalar field in addition to the fields of the SM.
The loop corrections involve all the gauge and scalar bosons which couple to
neutrinos.

In order to obtain the one-loop corrections to the elements of the light
neutrino mass matrix, we perform our computations in the $R_{\xi}$ gauge and
show explicitly the intricate cancellation of the gauge fixing parameters from
the corrections. In addition, we shall also demonstrate the cancellation of
the $\epsilon$ poles when the loop integrals are regulated by dimensional
regularization in $d=4-2\epsilon$ dimensions. These cancellations are highly
non-trivial, and therefore provide strong checks on the correctness of the
computations.

The paper is composed as follows. We introduce the model to the extent needed
for the present work in Sect. II. We define and compute the one loop
correction to mass matrix of the active neutrinos in Sect. III. In Sect. IV we
provide numerical estimates of the one-loop corrections and show that those
are very small. Finally we summarize our findings in Sec. V. We collect
auxiliary formulas in the appendices and also provide an auxiliary zip file
containing the SARAH model, parameter and particle files.

##  II Particle model, mixings and interactions

We consider an extension of the standard model by a U(1)z gauge group with
particle content and charge assignment defined in Ref. [27]. The super-weak
model is an economical extension of the standard model that provides a
framework to explain the origin of (i) neutrino masses and oscillations [28],
(ii) dark matter [29], (iii) cosmic inflation and stabilization of the
electroweak vacuum [30], (iv) matter-antimatter asymmetry of the universe. The
complete model including Feynman rules in the unitary gauge was presented
fully in Ref. [27]. As we are to compute one-loop corrections to neutrino
masses, we recall the details relevant to such computations, with Feynman
rules in the $R_{\xi}$ gauge. We generated those Feynman rules with SARAH[31,
32, 33, 34] but here we present simpler forms for the rules needed in our
computations to make those more comprehensive. We also recall some of the
conventions that are different in SARAH and the original definition of the
model. We stick to the SARAH conventions throughout this work.111We present
the model files in a separate file SuperWeak.zip.

###  II.1 Mixing of neutral gauge bosons

The particle content of the standard model is extended by 3 right-handed
neutrinos $\nu_{\mathrm{R}}$, a new scalar $\chi$, and the U(1)z gauge boson
$B^{\prime}$. As the field strength tensors of the U(1) gauge groups are gauge
invariant, kinetic mixing is allowed between the gauge fields belonging to the
hypercharge U(1)y and the new U(1)z gauge symmetries, whose strength is
measured by $\epsilon$ in

$\begin{split}\mathcal{L}&\supset-\frac{1}{4}F^{\mu\nu}F_{\mu\nu}-\frac{1}{4}F^{\prime\mu\nu}F^{\prime}_{\mu\nu}-\frac{\epsilon}{2}F^{\mu\nu}F^{\prime}_{\mu\nu}\,,\\\
\mathcal{D}^{\mathrm{U}(1)}_{\mu}&=-\mathrm{i}(yg_{y}B_{\mu}+zg_{z}B^{\prime}_{\mu})\end{split}$
(II.1)

where $B^{\mu}$ is the U(1)y gauge field. However, equivalently, we can choose
the basis—the convention in SARAH—in which the gauge-field strengths do not
mix, while the couplings are given by a $2\times 2$ coupling matrix in the
covariant derivative

$D^{\mathrm{U}(1)}_{\mu}=-\mathrm{i}\begin{pmatrix}y&z\end{pmatrix}\begin{pmatrix}\hat{g}_{yy}&\hat{g}_{yz}\\\
\hat{g}_{zy}&\hat{g}_{zz}\end{pmatrix}\begin{pmatrix}\hat{B}_{\mu}\\\
\hat{B}^{\prime}_{\mu}\end{pmatrix}$ (II.2)

where $y$ and $z$ are the U(1)y and U(1)z charges. We can parametrize the
coupling matrix as

$\hat{\textbf{g}}=\begin{pmatrix}\hat{g}_{yy}&\hat{g}_{yz}\\\
\hat{g}_{zy}&\hat{g}_{zz}\end{pmatrix}=\begin{pmatrix}g_{y}&-\eta
g^{\prime}_{z}\\\
0&g^{\prime}_{z}\end{pmatrix}\begin{pmatrix}\cos\epsilon^{\prime}&\sin\epsilon^{\prime}\\\
-\sin\epsilon^{\prime}&\cos\epsilon^{\prime}\end{pmatrix}\,.$ (II.3)

The coupling mixing matrix containing $\eta$ is equivalent to the kinetic
mixing in the Lagrangian (II.1) and the parameters of the two representations
are related by $g^{\prime}_{z}=g_{z}/\sqrt{1-\epsilon^{2}}$ and $\eta=\epsilon
g_{y}/g_{z}$. In this paper, it will be convenient to use the kinetic mixing
representation defined by (II.1).

The rotation with angle $\epsilon^{\prime}$ is unphysical as it can be
absorbed into the mixing of the neutral gauge fields $B^{\mu}$,
$B^{\prime\mu}$ and $W^{3\,\mu}$ to the mass eigenstates $A^{\mu}$, $Z^{\mu}$
and $Z^{\prime\mu}$, which then can be described by a rotation matrix

$\left(\begin{array}[]{c}\hat{B}^{\mu}\\\ W^{3\,\mu}\\\
\hat{B}^{\prime\mu}\end{array}\right)=\left(\begin{array}[]{ccc}\cos\theta_{\rm
W}&-\cos\theta_{Z}\sin\theta_{\rm W}&-\sin\theta_{Z}\sin\theta_{\rm W}\\\
\sin\theta_{\rm W}&\cos\theta_{Z}\cos\theta_{\rm W}&\cos\theta_{\rm
W}\sin\theta_{Z}\\\
0&-\sin\theta_{Z}&\cos\theta_{Z}\end{array}\right)\left(\begin{array}[]{c}A^{\mu}\\\
Z^{\mu}\\\ Z^{\prime\mu}\end{array}\right).$ (II.4)

This matrix depends on two mixing angles: $\theta_{\rm W}$ is the weak mixing
(or Weinberg) angle and $\theta_{Z}$ is the $Z-Z^{\prime}$ mixing angle
222Note the opposite sign convention for $\theta_{Z}$ in this work and in Ref.
[27] where this mixing angle was denoted as $\theta_{T}$, so
$\theta_{Z}=-\theta_{T}$.. In terms of the coupling parameters

$\kappa=\cos\theta_{\rm
W}(\gamma^{\prime}_{y}-2\gamma_{z}^{\prime})\text{~{}~{}and~{}~{}}\tau=2\cos\theta_{\rm
W}\gamma_{z}^{\prime}\tan\beta\,,$ (II.5)

introduced in Ref. [27], this new mixing angle is given implicitly by
$\tan(2\theta_{Z})=2\kappa/(1-\kappa^{2}-\tau^{2})$. In Eq. (II.5)
$\tan\beta=w/v$ is the ratio of the vacuum expectation values (VEVs) of the
scalar fields (see below) and
$\gamma^{\prime}_{y}=(\epsilon/\sqrt{1-\epsilon^{2}})(g_{y}/g_{\rm L})$,
$\gamma^{\prime}_{z}=g^{\prime}_{z}/g_{\rm L}$, i.e. the couplings are
normalized by the SU(2)L coupling.

We can express the elements of the $Z-Z^{\prime}$ mixing matrix explicitly,

$\begin{split}\sin\theta_{Z}&=\operatorname{sgn}\left(\kappa\right)\left[\frac{1}{2}\left(1-\frac{1-\kappa^{2}-\tau^{2}}{\sqrt{(1+\kappa^{2}+\tau^{2})^{2}-4\tau^{2}}}\right)\right]^{1/2},\\\
\cos\theta_{Z}&=\left[\frac{1}{2}\left(1+\frac{1-\kappa^{2}-\tau^{2}}{\sqrt{(1+\kappa^{2}+\tau^{2})^{2}-4\tau^{2}}}\right)\right]^{1/2}\,,\end{split}$
(II.6)

which also appear in the neutral currents
$\Gamma^{\mu}_{V\bar{f}f}=-\mathrm{i}e\gamma^{\mu}(C_{V\bar{f}f}^{R}P_{R}+C_{V\bar{f}f}^{L}P_{L})$
where $e$ is the electromagnetic coupling and $P_{R/L}\equiv
P_{\pm}=\frac{1}{2}(1\pm\gamma^{5})$ are the usual chiral projections. In
particular, for neutrinos

$\begin{split}eC_{Z\nu\nu}^{L}&=\frac{g_{\rm L}}{2\cos\theta_{\rm
W}}\Big{[}\cos\theta_{Z}-(\gamma_{y}^{\prime}-\gamma_{z}^{\prime})\sin\theta_{Z}\cos\theta_{\rm
W}\Big{]}\,,\quad eC_{Z\nu\nu}^{R}=-\frac{g_{\rm
L}}{2}\gamma_{z}^{\prime}\sin\theta_{Z}\,,\\\
eC_{Z^{\prime}\nu\nu}^{L}&=\frac{g_{\rm L}}{2\cos\theta_{\rm
W}}\Big{[}\sin\theta_{Z}+(\gamma_{y}^{\prime}-\gamma_{z}^{\prime})\cos\theta_{Z}\cos\theta_{\rm
W}\Big{]}\,,\quad eC_{Z^{\prime}\nu\nu}^{R}=\frac{g_{\rm
L}}{2}\gamma_{z}^{\prime}\cos\theta_{Z}\,,\end{split}$ (II.7)

i.e. $C_{Z^{\prime}\nu\nu}^{L/R}$ can be obtained from $C_{Z\nu\nu}^{L/R}$ by
the replacement

$(Z\to
Z^{\prime})\Rightarrow(\cos\theta_{Z},\sin\theta_{Z})\to(\sin\theta_{Z},-\cos\theta_{Z})\,.$
(II.8)

###  II.2 Mixings of scalar and Goldstone bosons

In addition to the usual $SU(2)_{\mathrm{L}}$-doublet Brout-Englert-Higgs
(BEH) field

$\phi=\left(\\!\\!\begin{array}[]{c}\phi^{+}\\\
\phi^{0}\end{array}\\!\\!\right)=\frac{1}{\sqrt{2}}\left(\\!\\!\begin{array}[]{c}\phi_{1}+\mathrm{i}\phi_{2}\\\
\phi_{3}+\mathrm{i}\phi_{4}\end{array}\\!\\!\right)\,,$ (II.9)

there is another complex scalar $\chi$ in the model, with charges specified in
[27]. The Lagrangian of the scalar fields contains the potential energy

$V(\phi,\chi)=V_{0}-\mu_{\phi}^{2}|\phi|^{2}-\mu_{\chi}^{2}|\chi|^{2}+\left(|\phi|^{2},|\chi|^{2}\right)\left(\\!\\!\begin{array}[]{cc}\lambda_{\phi}&\frac{\lambda}{2}\\\
\frac{\lambda}{2}&\lambda_{\chi}\end{array}\\!\\!\right)\left(\\!\\!\begin{array}[]{c}|\phi|^{2}\\\
|\chi|^{2}\end{array}\\!\\!\right)\subset-\mathcal{L}$ (II.10)

where $|\phi|^{2}=|\phi^{+}|^{2}+|\phi^{0}|^{2}$. In the $R_{\xi}$ gauge we
parametrize the scalar fields after spontaneous symmetry breaking as

$\begin{split}\phi=\frac{1}{\sqrt{2}}\binom{-\mathrm{i}\sqrt{2}\sigma^{+}}{v+h^{\prime}+\mathrm{i}\sigma_{\phi}}\,,\quad\chi=\frac{1}{\sqrt{2}}(w+s^{\prime}+\mathrm{i}\sigma_{\chi})\end{split}$
(II.11)

where $v$ and $w$ denotes the vacuum expectation values (VEVs) of the fields,
whose values are

$v=\sqrt{2}\sqrt{\frac{2\lambda_{\chi}\mu_{\phi}^{2}-\lambda\mu_{\chi}^{2}}{4\lambda_{\phi}\lambda_{\chi}-\lambda^{2}}}\,,\qquad
w=\sqrt{2}\sqrt{\frac{2\lambda_{\phi}\mu_{\chi}^{2}-\lambda\mu_{\phi}^{2}}{4\lambda_{\phi}\lambda_{\chi}-\lambda^{2}}}\,.$
(II.12)

Using the VEVs, we can express the quadratic couplings as

$\mu_{\phi}^{2}=\lambda_{\phi}v^{2}+\frac{\lambda}{2}w^{2}\,,\qquad\mu_{\chi}^{2}=\lambda_{\chi}w^{2}+\frac{\lambda}{2}v^{2}\,.$
(II.13)

The fields $h^{\prime}$ and $s^{\prime}$ are two real scalars and
$\sigma_{\phi}$ and $\sigma_{\chi}$ are the corresponding Goldstone bosons
that are weak eigenstates. We shall denote the mass eigenstates with $h$, $s$
and $\sigma_{Z}$, $\sigma_{Z^{\prime}}$. These different eigenstates are
related by the rotations

$\binom{h}{s}=\textbf{Z}_{S}\binom{h^{\prime}}{s^{\prime}}\equiv\begin{pmatrix}\cos\theta_{S}&-\sin\theta_{S}\\\
\sin\theta_{S}&~{}~{}\cos\theta_{S}\end{pmatrix}\binom{h^{\prime}}{s^{\prime}}$
(II.14)

and

$\binom{\sigma_{Z}}{\sigma_{Z^{\prime}}}=\textbf{Z}_{\rm
G}\binom{\sigma_{\phi}}{\sigma_{\chi}}\equiv\begin{pmatrix}\cos\theta_{\rm
G}&-\sin\theta_{\rm G}\\\ \sin\theta_{\rm G}&~{}~{}\cos\theta_{\rm
G}\end{pmatrix}\binom{\sigma_{\phi}}{\sigma_{\chi}}$ (II.15)

where $\theta_{S}$ and $\theta_{\rm G}$ are the scalar and Goldstone mixing
angles that can be determined by the diagonalization of the mass matrix of the
real scalars and that of the neutral Goldstone bosons.

The scalar mixing angle $\theta_{S}$ is related to the potential parameters by
[27]

$\tan(2\theta_{S})=-\frac{\lambda
vw}{\lambda_{\phi}v^{2}-\lambda_{\chi}w^{2}}\,.$ (II.16)

The condition $\theta_{S}\in(-\frac{\pi}{4},\frac{\pi}{4})$ implies that the
scalar mass eigenstates are not labeled by mass hierarchy.

The mass matrix of the Goldstone bosons is given in principle by the sum of
gauge-independent and gauge-dependent terms. However, the gauge-independent
terms vanish by Eq. (II.13):

$\left(\begin{array}[]{cc}\frac{1}{2}\lambda
w^{2}+\lambda_{\phi}v^{2}-\mu_{\phi}^{2}&0\\\ 0&\frac{1}{2}\lambda
v^{2}+\lambda_{\chi}w^{2}-\mu_{\chi}^{2}\end{array}\right)=\textbf{0}\,,$
(II.17)

so the mass matrix contains only gauge-dependent terms,

$\textbf{m}^{2}_{A}=\xi_{Z}\textbf{m}^{2}_{A_{Z}}+\xi_{Z^{\prime}}\textbf{m}^{2}_{A_{Z^{\prime}}}$
(II.18)

where $\xi_{Z}$ and $\xi_{{Z^{\prime}}}$ are the gauge parameters. The mass
matrix is symmetric, so we can write it formally as

$\textbf{m}^{2}_{A_{x}}=\left(\begin{array}[]{cc}m^{2}_{A_{x},11}&m^{2}_{A_{x},12}\\\
m^{2}_{A_{x},12}&m^{2}_{A_{x},22}\end{array}\right),$ (II.19)

for both $x=Z$ and ${Z^{\prime}}$. Explicitly,

$\begin{split}m^{2}_{A_{Z},11}&=v^{2}e^{2}\Big{(}C^{L}_{Z\nu\nu}-C^{R}_{Z\nu\nu}\Big{)}^{2}=\left(\frac{M_{W}}{\cos\theta_{\rm
W}}\right)^{2}(\cos\theta_{Z}-\kappa\sin\theta_{Z})^{2}\,,\\\
m^{2}_{A_{Z},12}&=2vwe^{2}\Big{(}C^{L}_{Z\nu\nu}-C^{R}_{Z\nu\nu}\Big{)}C^{R}_{Z\nu\nu}=\left(\frac{M_{W}}{\cos\theta_{\rm
W}}\right)^{2}(\cos\theta_{Z}-\kappa\sin\theta_{Z})(-\tau\sin\theta_{Z})\,,\\\
m^{2}_{A_{Z},22}&=w^{2}e^{2}\Big{(}2C^{R}_{Z\nu\nu}\Big{)}^{2}=\left(\frac{M_{W}}{\cos\theta_{\rm
W}}\right)^{2}(-\tau\sin\theta_{Z})^{2}\end{split}$ (II.20)

where $M_{W}=\frac{vg_{\rm L}}{2}$ is the mass of the W bosons, and the
elements of $\textbf{m}^{2}_{A_{Z^{\prime}}}$ can be obtained by the
replacement $Z\to Z^{\prime}$ in the chiral couplings, which implies the
replacement (II.8) in the second forms of the matrix elements. The latter are
the most convenient ones for the diagonalization of the mass matrix. Using Eq.
(II.6), one can check that the matrix

$\textbf{Z}_{\rm G}\textbf{m}^{2}_{A}\textbf{Z}_{\rm
G}^{T}=\textbf{m}^{2}_{\text{diag},A}$

is indeed diagonal provided we have for the Goldstone mixing angle

$\cos\theta_{\rm
G}=\frac{\cos\theta_{Z}-\kappa\sin\theta_{Z}}{\sqrt{(\cos\theta_{Z}-\kappa\sin\theta_{Z})^{2}+(\tau\sin\theta_{Z})^{2}}}$
(II.21)

and

$\sin\theta_{\rm
G}=\frac{\tau\sin\theta_{Z}}{\sqrt{(\cos\theta_{Z}-\kappa\sin\theta_{Z})^{2}+(\tau\sin\theta_{Z})^{2}}}\,.$
(II.22)

###  II.3 Masses of neutral gauge bosons

As expected, the elements of the diagonal matrix
$\textbf{m}^{2}_{\text{diag},A}$ coincide with the squares of the masses of
the neutral gauge bosons [27],

$M_{Z}^{2}=\left(\frac{M_{W}}{\cos\theta_{\rm
W}}\right)^{2}\Big{[}(\cos\theta_{Z}-\kappa\sin\theta_{Z})^{2}+(\tau\sin\theta_{Z})^{2}\Big{]}$
(II.23)

and

$M_{Z^{\prime}}^{2}=\left(\frac{M_{W}}{\cos\theta_{\rm
W}}\right)^{2}\Big{[}(\sin\theta_{Z}+\kappa\cos\theta_{Z})^{2}+(\tau\cos\theta_{Z})^{2}\Big{]}\,,$
(II.24)

which can also be expressed conveniently with the chiral couplings and
Goldstone mixing angle. First we note that using Eq. (II.23), we find the
simple relation

$\sin\theta_{\rm G}=\tau\frac{\sin\theta_{Z}}{\cos\theta_{\rm
W}}\frac{M_{W}}{M_{Z}}$ (II.25)

between the Goldstone and neutral boson mixing angles, and also

$\cos\theta_{\rm G}=\tau\frac{\cos\theta_{Z}}{\cos\theta_{\rm
W}}\frac{M_{W}}{M_{Z^{\prime}}}\,.$ (II.26)

Next, we can substitute the relations found in Eq. (II.20) into Eqs. (II.23)
and (II.24) together with the definition of the right handed couplings defined
in Eq. (II.7), resulting in

$M_{Z}^{2}=v^{2}e^{2}\Big{(}C^{L}_{Z\nu\nu}-C^{R}_{Z\nu\nu}\Big{)}^{2}+w^{2}g_{z}^{\prime\,2}\sin^{2}\theta_{Z}$
(II.27)

and also using Eq. (II.8),

$M_{Z^{\prime}}^{2}=v^{2}e^{2}\Big{(}C^{L}_{Z^{\prime}\nu\nu}-C^{R}_{Z^{\prime}\nu\nu}\Big{)}^{2}+w^{2}g_{z}^{\prime\,2}\cos^{2}\theta_{Z}\,.$
(II.28)

From Eq. (II.25) and (II.26) we can express

$wg_{z}^{\prime}\sin\theta_{Z}=M_{Z}\sin\theta_{\rm
G}\textrm{~{}~{}and~{}~{}}wg_{z}^{\prime}\cos\theta_{Z}=M_{Z^{\prime}}\cos\theta_{\rm
G}\,,$ (II.29)

which after substitution and simple rearrangement leads to

$M_{Z}^{2}=\frac{v^{2}e^{2}}{\cos^{2}\theta_{\rm
G}}\Big{(}C^{L}_{Z\nu\nu}-C^{R}_{Z\nu\nu}\Big{)}^{2}\,,\quad
M_{Z^{\prime}}^{2}=\frac{v^{2}e^{2}}{\sin^{2}\theta_{\rm
G}}\Big{(}C^{L}_{Z^{\prime}\nu\nu}-C^{R}_{Z^{\prime}\nu\nu}\Big{)}^{2}\,.$
(II.30)

###  II.4 Mass terms and mixing of neutrinos

The masses of the neutrinos are generated by the leptonic Yukawa terms in the
Lagrangian [27],

$-\mathcal{L}_{Y}^{\ell}=\frac{1}{2}\overline{\nu_{R}^{c}}\>\textbf{Y}_{N}\>\nu_{R}\chi+\overline{L_{L}}\>\phi^{c}\>\textbf{Y}_{\nu}\>\nu_{R}+\text{h.c.}$
(II.31)

where $\overline{L_{L}}$ is the Dirac adjoint of the left handed lepton
dublet, $\textbf{Y}_{N}$ and $\textbf{Y}_{\nu}$ are $3\times 3$ matrices, the
superscript $c$ denotes charge conjugation,
$\nu^{c}=-\mathrm{i}\gamma_{2}\nu^{*}$. After SSB this Lagrangian becomes

$-\mathcal{L}_{Y}^{\ell}=\frac{w+s^{\prime}+\mathrm{i}\sigma_{\chi}}{2\sqrt{2}}\overline{\nu_{R}^{c}}\>\textbf{Y}_{N}\>\nu_{R}+\frac{v+h^{\prime}-\mathrm{i}\sigma_{\phi}}{\sqrt{2}}\overline{\nu_{L}}\>\textbf{Y}_{\nu}\>\nu_{R}+\text{h.c.}$
(II.32)

and the terms proportional to the VEVs provide the mass matrices

$\textbf{M}_{N}=\frac{w}{\sqrt{2}}\textbf{Y}_{N}\,,\quad\textbf{M}_{D}=\frac{v}{\sqrt{2}}\textbf{Y}_{\nu}$
(II.33)

where the Majorana mass matrix $\textbf{M}_{N}$ is real and symmetric, while
the Dirac mass matrix $\textbf{M}_{D}$ is complex and Hermitian.

In flavour basis the $6\times 6$ mass matrix for the neutrinos that can be
written in terms of $3\times 3$ blocks as

$\textbf{M}^{\prime}=\begin{pmatrix}\textbf{0}_{3}&\textbf{M}_{D}^{T}\\\
\textbf{M}_{D}&\textbf{M}_{N}\end{pmatrix}.$ (II.34)

The weak (flavour) eigenstates
$(\nu_{e},\,\nu_{\mu},\,\nu_{\tau},\,\nu_{R,1},\,\nu_{R,2},\,\nu_{R,3})$ can
be transformed into the basis of $\nu_{i}$ ($i=1-6$) mass eigenstates with a
$6\times 6$ unitary matrix 333This matrix coincides with the unitary matrix
$U^{V,\dagger}$ used by SARAH. U where the mass matrix is diagonal,

$\textbf{U}^{T}\textbf{M}^{\prime}\textbf{U}=\textbf{M}=\operatorname{diag}(m_{1},m_{2},m_{3},m_{4},m_{5},m_{6})\,.$
(II.35)

It is helpful to decompose the matrix U into two $3\times 6$ blocks
$\textbf{U}_{L}$ and $\textbf{U}_{R}^{*}$,

$\textbf{U}=\binom{\textbf{U}_{L}}{\textbf{U}_{R}^{*}}\,,$ (II.36)

so $\textbf{U}^{T}=(\textbf{U}_{L}^{T},\textbf{U}_{R}^{\dagger})$ where both
blocks are $6\times 3$ matrices. It may be worth to emphasize that in spite of
what might be implied by the notation, the matrices $\textbf{U}_{L}$ and
$\textbf{U}_{R}^{*}$ are only semi-unitary. Useful relations of these matrices
are collected in Appendix A.

###  II.5 Gauge boson – neutrino interactions

As the neutral currents are written in terms of flavour eigenstates, the
interactions between the neutral gauge bosons and the propagating mass
eigenstate neutrinos include also the neutrino mixing matrices:

${\mbox{\boldmath$\Gamma$}}^{\mu}_{V\nu_{i}\nu_{j}}=-\mathrm{i}e\gamma^{\mu}\Big{(}{\mbox{\boldmath$\Gamma$}}^{L}_{V\nu\nu}P_{L}+{\mbox{\boldmath$\Gamma$}}^{R}_{V\nu\nu}P_{R}\Big{)}_{ij}$
(II.37)

where

${\mbox{\boldmath$\Gamma$}}^{L}_{V\nu\nu}=C^{L}_{V\nu\nu}\textbf{U}_{L}^{\dagger}\textbf{U}_{L}-C^{R}_{V\nu\nu}\textbf{U}_{R}^{T}\textbf{U}_{R}^{*}$
(II.38)

and

${\mbox{\boldmath$\Gamma$}}^{R}_{V\nu\nu}=-C^{L}_{V\nu\nu}\textbf{U}_{L}^{T}\textbf{U}^{*}_{L}+C^{R}_{V\nu\nu}\textbf{U}_{R}^{\dagger}\textbf{U}_{R}=-\Big{(}{\mbox{\boldmath$\Gamma$}}^{L}_{V\nu\nu}\Big{)}^{*}$
(II.39)

for both $V=Z$ and $V=Z^{\prime}$.

###  II.6 Scalar boson – neutrino and Goldstone boson – neutrino interactions

The terms containing the scalar and Goldstone bosons in Eq. (II.32) provide
interactions between those and the neutrinos. These interactions have the same
structure with small differences. For the propagating scalar states $S_{k}$ or
$\sigma_{k}$ ($k=1$ denoting $h$ or the Goldstone boson belonging to $Z$ and
$k=2$ referring to $s$ or the Goldstone boson belonging to $Z^{\prime}$) such
interactions can be decomposed into left and right chiral terms

${\mbox{\boldmath$\Gamma$}}_{S_{k}/\sigma_{k}\,\nu_{i}\nu_{j}}=\Big{(}{\mbox{\boldmath$\Gamma$}}_{S_{k}/\sigma_{k}\,\nu\nu}^{L}P_{L}+{\mbox{\boldmath$\Gamma$}}_{S_{k}/\sigma_{k}\,\nu\nu}^{R}P_{R}\Big{)}_{ij}$
(II.40)

where the matrices ${\mbox{\boldmath$\Gamma$}}^{L/R}$ contain both the mixing
matrix of the neutrinos and the mixing matrix of the scalar or Goldstone
bosons. The left chiral coefficients are

$\displaystyle{\mbox{\boldmath$\Gamma$}}^{L}_{S_{k}\nu\nu}$
$\displaystyle=-\mathrm{i}\left[\Big{(}\textbf{M}\textbf{U}_{L}^{\dagger}\textbf{U}_{L}+\textbf{U}_{L}^{T}\textbf{U}_{L}^{*}\textbf{M}\Big{)}\frac{(\textbf{Z}_{S})_{k1}}{v}+\textbf{U}_{R}^{\dagger}\textbf{M}_{N}\textbf{U}_{R}^{*}\frac{(\textbf{Z}_{S})_{k2}}{w}\right]\,,$
(II.41)

and

$\displaystyle{\mbox{\boldmath$\Gamma$}}^{L}_{\sigma_{k}\nu\nu}$
$\displaystyle=-\left[\Big{(}\textbf{M}\textbf{U}_{L}^{\dagger}\textbf{U}_{L}+\textbf{U}_{L}^{T}\textbf{U}_{L}^{*}\textbf{M}\Big{)}\frac{(\textbf{Z}_{\rm
G})_{k1}}{v}+\textbf{U}_{R}^{\dagger}\textbf{M}_{N}\textbf{U}_{R}^{*}\frac{(\textbf{Z}_{\rm
G})_{k2}}{w}\right]$ (II.42)

and the right chiral ones are related by complex conjugation,
${\mbox{\boldmath$\Gamma$}}^{R}_{S_{k}/\sigma_{k}\,\nu\nu}=-\Big{(}{\mbox{\boldmath$\Gamma$}}^{L}_{S_{k}/\sigma_{k}\,\nu\nu}\Big{)}^{*}$.

##  III Neutrino mass matrix at one-loop order

We are interested in the one-loop correction $\delta\textbf{M}_{L}$ to the
tree-level mass matrix of the light neutrinos. In perturbation theory we deal
with propagating states which are mass eigenstates. Hence, we can compute loop
corrections to self energies of mass eigenstates of neutrinos. The neutrino
mass matrix at one-loop order is then obtained from Eq. (II.35), with diagonal
mass matrix substituted at one loop, $\textbf{M}+\delta\textbf{M}$ where

$\delta\textbf{M}=\operatorname{diag}(\delta m_{1},\delta m_{2},\delta
m_{3},\delta m_{4},\delta m_{5},\delta m_{6})\,.$ (III.1)

Hence, the correction is obtained by

$\delta\textbf{M}^{\prime}=\begin{pmatrix}\delta\textbf{M}_{L}&\delta\textbf{M}_{D}^{T}\\\
\delta\textbf{M}_{D}&\delta\textbf{M}_{N}\end{pmatrix}=\textbf{U}^{*}\delta\textbf{M}\textbf{U}^{\dagger}\,.$
(III.2)

Using Eq. (II.36), we can compute the $3\times 3$ blocks as

$\delta\textbf{M}_{L}=\textbf{U}_{L}^{*}\delta\textbf{M}\textbf{U}_{L}^{\dagger},\quad\delta\textbf{M}_{D}=\textbf{U}_{R}\delta\textbf{M}\textbf{U}_{L}^{\dagger},\quad\delta\textbf{M}_{N}=\textbf{U}_{R}\delta\textbf{M}\textbf{U}_{R}^{T}\,.$
(III.3)

In the following subsections we prove that the one-loop correction to the mass
matrix of the active neutrinos have the form

$\delta\textbf{M}_{L}=\frac{1}{16\pi^{2}}\sum_{k=1,2}\left[3(\textbf{Z}_{\rm
G})_{k1}^{2}\frac{M^{2}_{V_{k}}}{v^{2}}\textbf{F}(M^{2}_{V_{k}})+(\textbf{Z}_{\rm
S})_{k1}^{2}\frac{M^{2}_{S_{k}}}{v^{2}}\textbf{F}(M^{2}_{S_{k}})\right]$
(III.4)

where we introduced the finite matrix valued function

$\textbf{F}_{ij}(M^{2})=\sum_{a=1}^{6}(\textbf{U}_{L}^{*})_{ia}(\textbf{U}_{L}^{\dagger})_{aj}\frac{m_{a}^{3}}{M^{2}}\frac{\ln\frac{m_{a}^{2}}{M^{2}}}{\frac{m_{a}^{2}}{M^{2}}-1}$
(III.5)

of dimension mass and with summation running over all neutrinos.

###  III.1 Self-energy decomposition

The neutrino self energy is a $6\times 6$ matrix that can be decomposed as

$\mathrm{i}{\mbox{\boldmath$\Sigma$}}(p)=\textbf{A}_{L}(p^{2})\not{p}P_{L}+\textbf{A}_{R}(p^{2})\not{p}P_{R}+\textbf{B}_{L}(p^{2})P_{L}+\textbf{B}_{R}(p^{2})P_{R}\,.$
(III.6)

Using this decomposition, $\delta\textbf{M}_{L}$ is given by [23]

$\delta\textbf{M}_{L}=\textbf{U}_{L}^{*}\textbf{B}_{L}(0)\textbf{U}_{L}^{\dagger}\,.$
(III.7)

The matrix $\textbf{B}_{L}(0)$ receives contributions involving a neutrino and
either a neutral vector boson Z, Z’, or a scalar boson $\sigma_{Z}$,
$\sigma_{Z^{\prime}}$ (Goldstone boson), $h$, $s$ (Higgs-like scalar) in the
loop. The relevant Feynman graphs that give contributions to the neutrino self
energies at one-loop order are shown in Fig. 1. There are also tadpole
contributions to $\textbf{B}_{L}(0)$. Those are proportional to the scalar-
neutrino coupling ${\mbox{\boldmath$\Gamma$}}^{L}_{S_{k}\nu_{i}\nu_{j}}$ given
in Eq. (II.40), which vanishes when sandwiched between $\textbf{U}_{L}^{*}$
and $\textbf{U}_{L}^{\dagger}$, see Eq. (A.5). The charged vector boson
together with a charged lepton in the loop (bottom right diagram in Fig. 1)
contributes only to $\textbf{A}_{L/R}$. Thus we compute the first three graphs
explicitly. For a given boson $x$ in the loop, the matrix $\textbf{B}_{L}(0)$
depends on the mass $M_{x}$ and also the tree-level masses of the neutrinos
$\\{m_{a}\\}$, $\textbf{B}_{L}(0)=\textbf{B}_{L}^{x}(M_{x},\\{m_{a}\\})$.

Figure 1:  One-loop graphs contributing to the neutrino self energy. Top left:
Goldstone boson contribution. Top right: scalar contribution. Bottom left:
neutral gauge boson contribution. Bottom right: charged gauge boson
contribution. Note that $W$ boson loop does not contribute to the matrix BL.

###  III.2 Contributions with neutral gauge bosons in the loop

The contribution of the neutral gauge boson $V$ is

$\Big{(}B_{L}^{V}(M_{V},\\{m_{a}\\};\xi_{V})\Big{)}_{ij}P_{L}=\mathrm{i}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\sum_{a=1}^{6}\Gamma^{\mu}_{V\nu_{i}\nu_{a}}\frac{\not{p}-\not{\ell}+m_{a}}{(p-\ell)^{2}-m_{a}^{2}}\Gamma^{\nu}_{V\nu_{a}\nu_{j}}P_{\mu\nu}(\ell,M_{V}^{2};\xi_{V})$
(III.8)

where $\xi_{V}$ is the gauge parameter and

$P_{\mu\nu}(\ell,M_{V}^{2};\xi_{V})=\frac{g_{\mu\nu}}{\ell^{2}-M_{V}^{2}}-(1-\xi_{V})\frac{\ell_{\mu}\ell_{\nu}}{(\ell^{2}-M_{V}^{2})(\ell^{2}-\xi_{V}M_{V}^{2})}\,.$
(III.9)

Introducing the $6\times 6$ matrix

$\textbf{m}_{\ell}^{(n)}=\operatorname{diag}\left(\frac{m_{1}^{n}}{\ell^{2}-m_{1}^{2}},\dots,\frac{m_{6}^{n}}{\ell^{2}-m_{6}^{2}}\right),$
(III.10)

and using the result of Appendix B, we obtain the following expression for a
neutral vector boson in the loop:

$\delta\textbf{M}_{L}^{V}=\mathrm{i}e^{2}\Big{(}C^{L}_{V\nu\nu}-C^{R}_{V\nu\nu}\Big{)}^{2}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\textbf{U}_{L}^{*}\bigg{[}\frac{d\>\textbf{m}_{\ell}^{(1)}}{\ell^{2}-M_{V}^{2}}+\frac{\textbf{m}_{\ell}^{(3)}}{M_{V}^{2}}\left(\frac{1}{\ell^{2}-\xi_{V}M_{V}^{2}}-\frac{1}{\ell^{2}-M_{V}^{2}}\right)\bigg{]}\textbf{U}_{L}^{\dagger}.$
(III.11)

###  III.3 Contributions with neutral Goldstone bosons in the loop

The contribution of the neutral Goldstone boson $\sigma_{V}$ ($V=1$ means the
Goldstone boson belonging to the $Z$ field and $V=2$ refers to the
$Z^{\prime}$ field) is

$\Big{(}B_{L}^{\sigma_{V}}(m_{\sigma_{V}},\\{m_{a}\\};\xi_{V})\Big{)}_{ij}P_{L}=-\mathrm{i}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\sum_{a=1}^{6}{\mbox{\boldmath$\Gamma$}}_{\sigma_{V}\nu_{i}\nu_{a}}\frac{m_{a}}{\ell^{2}-m_{a}^{2}}{\mbox{\boldmath$\Gamma$}}_{\sigma_{V}\nu_{a}\nu_{j}}\frac{1}{\ell^{2}-\xi_{V}M^{2}_{V}}\,.$
(III.12)

Using the matrix notation, we can write

$\textbf{U}_{L}^{*}\textbf{B}_{L}^{\sigma_{V}}(0)\textbf{U}_{L}^{\dagger}P_{L}=-\mathrm{i}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\textbf{U}_{L}^{*}{\mbox{\boldmath$\Gamma$}}_{\sigma_{V}\nu\nu}\textbf{m}_{\ell}^{(1)}{\mbox{\boldmath$\Gamma$}}_{\sigma_{V}\nu\nu}\textbf{U}_{L}^{\dagger}\frac{1}{\ell^{2}-\xi_{V}M^{2}_{V}}\,.$
(III.13)

Substituting the vertex functions of Eq. (II.40) and employing the matrix
relations in Eqs. (A.2) and (A.5), we obtain the correction to the mass matrix
as

$\delta\textbf{M}_{L}^{\sigma_{V}}=-\mathrm{i}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\textbf{U}_{L}^{*}\textbf{M}\textbf{m}_{\ell}^{(1)}\textbf{M}\textbf{U}_{L}^{\dagger}\left(\frac{(\textbf{Z}_{\rm
G})_{V1}}{v}\right)^{2}\frac{1}{\ell^{2}-\xi_{V}M^{2}_{V}}\,.$ (III.14)

We now substitute
$\textbf{M}\textbf{m}_{\ell}^{(1)}\textbf{M}=\textbf{m}_{\ell}^{(3)}$ and
using Eq. (II.30), we obtain

$\delta\textbf{M}_{L}^{\sigma_{V}}=-\mathrm{i}e^{2}\Big{(}C^{L}_{V\nu\nu}-C^{R}_{V\nu\nu}\Big{)}^{2}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\textbf{U}_{L}^{*}\frac{\textbf{m}_{\ell}^{(3)}}{M_{V}^{2}}\textbf{U}_{L}^{\dagger}\frac{1}{\ell^{2}-\xi_{V}M^{2}_{V}}\,.$
(III.15)

###  III.4 Contributions with scalar bosons in the loop

The scalar – neutrino vertex is very similar to the Goldstone boson neutrino
vertex, so the contribution with a scalar boson $S_{k}$ in the loop can be
written immediately in analogy with Eq. (III.14):

$\begin{split}\delta\textbf{M}_{L}^{S_{k}}&=\mathrm{i}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\textbf{U}_{L}^{*}\textbf{M}\textbf{m}_{\ell}^{(1)}\textbf{M}\textbf{U}_{L}^{\dagger}\left(\frac{(\textbf{Z}_{S})_{k1}}{v}\right)^{2}\frac{1}{\ell^{2}-M^{2}_{S_{k}}}\\\
&=\mathrm{i}\left(\frac{(\textbf{Z}_{S})_{k1}}{v}\right)^{2}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\textbf{U}_{L}^{*}\textbf{m}_{\ell}^{(3)}\textbf{U}_{L}^{\dagger}\frac{1}{\ell^{2}-M^{2}_{S_{k}}}\,.\end{split}$
(III.16)

###  III.5 The complete one-loop mass correction

Combining Eqs. (III.11), (III.15) and (III.16), we find that that the gauge-
dependent pieces of the vector boson contribution cancel exactly with the
Goldstone boson contribution, and obtain

$\begin{split}\delta\textbf{M}_{L}&=\sum\limits_{V=Z,Z^{\prime}}\left(\delta\textbf{M}_{L}^{V}+\delta\textbf{M}_{L}^{\sigma_{V}}\right)+\sum\limits_{k=1,2}\delta\textbf{M}_{L}^{S_{k}}\\\
&=\sum\limits_{V=Z,Z^{\prime}}\mathrm{i}e^{2}\Big{(}C^{L}_{V\nu\nu}-C^{R}_{V\nu\nu}\Big{)}^{2}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\textbf{U}_{L}^{*}\bigg{[}\frac{d\>\textbf{m}_{\ell}^{(1)}}{\ell^{2}-M_{V}^{2}}-\frac{\textbf{m}_{\ell}^{(3)}}{M_{V}^{2}}\frac{1}{\ell^{2}-M_{V}^{2}}\bigg{]}\textbf{U}_{L}^{\dagger}\\\
&\quad+\sum\limits_{k=1,2}\mathrm{i}\left(\frac{(\textbf{Z}_{S})_{k1}}{v}\right)^{2}\int\\!\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\textbf{U}_{L}^{*}\frac{\textbf{m}_{\ell}^{(3)}}{\ell^{2}-M^{2}_{S_{k}}}\textbf{U}_{L}^{\dagger}\,.\end{split}$
(III.17)

Introducing the integral

$I_{0}(M^{2},m_{a}^{2};\mu^{2},\epsilon)=\mu^{2\epsilon}\int\frac{\mathrm{d}^{d}\ell}{(2\pi)^{d}}\frac{1}{(\ell^{2}-M^{2})(\ell^{2}-m_{a}^{2})}\,,$
(III.18)

the matrix $\textbf{I}^{(n)}$ with elements

$\displaystyle\textbf{I}^{(n)}_{ij}(M^{2})=\mathrm{i}\sum\limits_{a=1}^{6}(\textbf{U}_{L}^{*})_{ia}m_{a}^{n}(\textbf{U}_{L}^{\dagger})_{aj}I_{0}(M^{2},m_{a}^{2};\mu^{2},\epsilon)\,,$
(III.19)

and using the relations (II.30) allows us to recast Eq. (III.17) into a neatly
condensed form

$\displaystyle\delta\textbf{M}_{L}$
$\displaystyle=\sum\limits_{k=1,2}\left\\{\left(\frac{(\textbf{Z}_{\rm
G})_{k1}}{v}\right)^{2}\left(dM_{V_{k}}^{2}\textbf{I}^{(1)}(M_{V_{k}}^{2})-\textbf{I}^{(3)}(M_{V_{k}}^{2})\right)+\left(\frac{(\textbf{Z}_{S})_{k1}}{v}\right)^{2}\textbf{I}^{(3)}(M_{S_{k}}^{2})\right\\}$
(III.20)

with $V_{1}=Z$ and $V_{2}=Z^{\prime}$. In Eq. (III.18) $2\epsilon=d-4$ and
$\mu$ is the regularization scale.

###  III.6 Finiteness and scale independence of $\delta\textbf{M}_{L}$

We show here that the one loop mass correction $\delta\textbf{M}_{L}$ is
finite and independent of the scale $\mu$. Evaluating the integral (III.18)
yields

$\displaystyle I_{0}(M^{2},m^{2};\mu^{2},\epsilon)$
$\displaystyle=I^{(\mathrm{s})}_{0}(\epsilon)+I^{(\mathrm{f})}_{0}\left(\frac{m^{2}}{M^{2}},\frac{\mu^{2}}{M^{2}}\right)+\mathrm{O}(\epsilon)$
(III.21)

where ‘s’ stands for the singular and ‘f’ for the finite functions

$\displaystyle
I^{(\mathrm{s})}_{0}(\epsilon)=\frac{\mathrm{i}}{16\pi^{2}}\left(\frac{1}{\epsilon}-\gamma_{\mathrm{E}}+\ln
4\pi+1\right)$ $\displaystyle,\quad
I^{(\mathrm{f})}_{0}(x,x_{\mu})=\frac{\mathrm{i}}{16\pi^{2}}\left(\frac{x\ln
x}{1-x}+\ln x_{\mu}\right).$ (III.22)

with $\gamma_{\mathrm{E}}\simeq 0.5722$ being the Euler-Mascheroni constant.
It is also convenient to split the matrix (III.19) in a similar fashion

$\textbf{I}^{(n)}(M^{2})=\textbf{I}^{(\mathrm{s},n)}+\textbf{I}^{(\mathrm{f},n)}(M^{2})$
(III.23)

such that

$\textbf{I}^{(\mathrm{s},n)}=\mathrm{i}\textbf{U}_{L}^{*}\textbf{M}^{n}\textbf{U}_{L}^{\dagger}I^{(\mathrm{s},n)}_{0}\,,\quad\Big{(}\textbf{I}^{(\mathrm{f},n)}(M^{2})\Big{)}_{ij}=\mathrm{i}\sum_{a=1}^{6}(\textbf{U}_{L}^{*})_{ia}m_{a}^{n}(\textbf{U}_{L}^{\dagger})_{aj}I^{(\mathrm{f},n)}_{0}\left(\frac{m_{a}^{2}}{M^{2}},\frac{\mu^{2}}{M^{2}}\right).$
(III.24)

Then the one-loop correction to the mass matrix of the light neutrinos can
also be decomposed as

$\delta\textbf{M}_{L}=\delta\textbf{M}^{(\mathrm{s})}_{L}+\delta\textbf{M}^{(\mathrm{f})}_{L}+\mathrm{O}(\epsilon)$
(III.25)

where

$\displaystyle\delta\textbf{M}^{(\mathrm{s})}_{L}$
$\displaystyle=\sum_{k=1,2}\left[dM_{V_{k}}^{2}\left(\frac{(\textbf{Z}_{\rm
G})_{k1}}{v}\right)^{2}\textbf{I}^{(\mathrm{s},1)}+\frac{(\textbf{Z}_{S})^{2}_{k1}-(\textbf{Z}_{\rm
G})^{2}_{k1}}{v^{2}}\textbf{I}^{(\mathrm{s},3)}\right]$ (III.26)

and

$\displaystyle\delta\textbf{M}^{(\mathrm{f})}_{L}$
$\displaystyle=\sum\limits_{k=1,2}\left[\left(\frac{(\textbf{Z}_{\rm
G})_{k1}}{v}\right)^{2}\left(dM_{V_{k}}^{2}\textbf{I}^{(\mathrm{f},1)}(M_{V_{k}}^{2})-\textbf{I}^{(\mathrm{f},3)}(M_{V_{k}}^{2})\right)+\left(\frac{(\textbf{Z}_{S})_{k1}}{v}\right)^{2}\textbf{I}^{(\mathrm{f},3)}(M_{S_{k}}^{2})\right].$
(III.27)

In order to prove that $\delta\textbf{M}_{L}$ is finite, one has to show that
$\delta\textbf{M}^{(\mathrm{s})}_{L}$ is free from $\epsilon$ poles. We prove
that it in fact vanishes because the matrix $\textbf{I}^{(\mathrm{s},1)}$ is
zero matrix due to the identity (A.5), while the coefficient in the second
term cancels because the matrices ZS and ZG are orthogonal, so

$\sum_{k=1}^{2}(\textbf{Z}_{S})_{k1}^{2}-\sum_{k=1}^{2}(\textbf{Z}_{G})_{k1}^{2}=0\,.$
(III.28)

Hence the mass independent terms, including the divergent pieces of the light
neutrino one-loop mass correction cancel, and we can set $\epsilon=0$, which
yields $\delta\textbf{M}_{L}=\delta\textbf{M}^{(\mathrm{f})}_{L}$.
Furthermore, the terms depending on the regularization scale in
$\delta\textbf{M}^{(\mathrm{f})}_{L}$ cancel in an identical way as the second
term does in $\delta\textbf{M}^{(\mathrm{s})}_{L}$ (using Eq. (III.28)).

The remaining finite terms give the final, regularization-scale independent
and finite one-loop correction to the light neutrinos as given in Eq. (III.4).
It is the linear combination of the matrix valued function F given in Eq.
(III.5) with different arguments and coefficients corresponding to the
different one-loop contributions. The function F gives the mass correction
corresponding to a one loop diagram, before coupling suppression; see Fig. 4
in Sect. IV for details where we shall also give a numerical estimate for its
eigenvalues $\delta m^{\nu}_{i,0}$. It is well defined for any non-negative
$x$ because

$\lim\limits_{x\rightarrow 0}\frac{x\ln
x}{1-x}=0\text{~{}~{}and~{}~{}}\lim\limits_{x\rightarrow 1}\frac{x\ln
x}{1-x}=-1\,.$ (III.29)

###  III.7 Generalization to arbitrary number of neutral bosons and neutrinos

Our predictions for the one-loop correction to the light neutrino mass matrix
can easily be generalized to any number of $n_{V}$ massive neutral gauge
bosons, $n_{S}$ neutral real scalars coupling to $n_{a}$ active and $n_{s}$
sterile neutrinos. Clearly, the matrix form of gauge-dependent parts in Eq.
(III.11) and Eq. (III.15) is unchanged, and they cancel in the same way.

The correction without gauge parameters $\xi_{V}$ in Eq. (III.4) is
straightforwardly generalized to a case where the sums go over an arbitrary
positive integer $n_{V}$ and $n_{S}$.

The neutrino mass and mixing matrices with arbitrary $n_{a}$ and $n_{s}$ are
written identically in the block form, differing only on the block shape: UL
is a $n_{a}\times(n_{a}+n_{s})$ matrix and UR is a $n_{s}\times(n_{a}+n_{s})$
matrix. The finite correction derived in Eq. (III.4) is then immediately
generalized to

$\begin{split}\delta\textbf{M}_{L}=\frac{1}{16\pi^{2}}\left[3\sum_{k=1}^{n_{V}}(\textbf{Z}_{\rm
G})_{k1}^{2}\frac{M^{2}_{V_{k}}}{v^{2}}\textbf{F}(M_{V_{k}}^{2})+\sum_{k=1}^{n_{S}}(\textbf{Z}_{\rm
S})_{k1}^{2}\frac{M^{2}_{S_{k}}}{v^{2}}\textbf{F}(M_{S_{k}}^{2})\right]\end{split}$
(III.30)

where the upper limit in the summation in the matrix F is $n_{a}+n_{s}$. The
factor 3 in front of the first term in the bracket of Eq. (III.30) stems from
the three polarization states of the propagating massive neutral gauge bosons.
The corresponding factor is of course one in the case of the scalars. This
formula is also independent of the new U(1) charge assignments.

##  IV Numerical estimate of the corrections

We are now ready to estimate the order of magnitude of the corrections. We
assume large mixing in the scalar sector: $\theta_{S}=\mathrm{O}(1)$. The
$Z^{\prime}$ mass and mixing angle $\theta_{G}$ are fixed by the gauge
couplings $g_{y}^{\prime}=\gamma_{y}^{\prime}g_{L}$ and $g_{z}^{\prime}$ and
ratio of VEVs, $\tan\beta\equiv w/v$. We plot their magnitudes in Fig. 2,
scanning the parameters $g_{y}^{\prime},g_{z}^{\prime}\in[10^{-6},1]$ and for
$w=100,750$ GeV. Note that larger $\tan\beta$ corresponds to a larger
Goldstone angle. Smaller $\tan\beta$ distorts the $M_{Z^{\prime}}$ contours so
that the same $Z^{\prime}$ mass can be achieved with larger gauge couplings
$g_{y}^{\prime}$ and $g_{z}^{\prime}$ compared to large $\tan\beta$. In
addition, we set $M_{s}/v=\mathrm{O}(1)$, that is, only the mass of the
$Z^{\prime}$ boson is free, and may be far from electroweak scale. The
relevant gauge couplings can then be estimated as from Fig. 2 after
identifying the region in $(g_{y}^{\prime},g_{z}^{\prime})$ plane
corresponding to $M_{Z^{\prime}}\in[20,200]$ MeV, which is the relevant mass
region for the super-weak model to reproduce the dark matter relic density,
allowed by experimental constraints [29].

Figure 2:  Absolute values of $\sin\theta_{G}$ (top) and mass of $Z^{\prime}$
boson (bottom) in logarithmic $(g_{y}^{\prime},g_{z}^{\prime})$ plane. For
$\theta_{G}$ the contour labels $n$ correspond to value $10^{n}$. Left plots:
$w=100$ GeV; right plots: $w=750$ GeV.

Then we identify the order-of-magnitude estimate for $|\sin\theta_{G}|$ by
comparing the regions relevant to the mass range of $M_{Z^{\prime}}$. For
$w=100$ GeV, we have $|\sin\theta_{G}|<10^{-6}$, which we take as a
conservative upper limit. Then the prefactors in gauge boson contributions to
$\delta\textbf{M}_{L}$ are

$e^{2}(C^{L}_{Z\nu\nu}-C^{R}_{Z\nu\nu})^{2}=\cos^{2}\theta_{G}\frac{M_{Z}^{2}}{v^{2}}\sim\mathrm{O}(10^{-1})$
(IV.1)

and

$e^{2}(C^{L}_{Z^{\prime}\nu\nu}-C^{R}_{Z^{\prime}\nu\nu})^{2}=\sin^{2}\theta_{G}\frac{M_{Z^{\prime}}^{2}}{v^{2}}\sim\mathrm{O}(10^{-19})\times\left(\frac{M_{Z^{\prime}}}{100\,\text{MeV}}\right)^{2}.$
(IV.2)

Then the numerical estimate for the total correction in Eq. (III.4) can be
written as

$\displaystyle(\delta\textbf{M}_{L})_{ij}<\mathrm{O}(10^{-7})\,\text{eV}+\mathrm{O}(10^{-21})\times\left(\frac{M_{Z^{\prime}}}{100\text{
MeV}}\right)^{2}\textbf{F}_{ij}(M_{Z^{\prime}}^{2})\,.$ (IV.3)

The elements of the matrix F are plotted as a function of the mass of the
boson of the loop $m_{\rm loop}$ in Fig. 3 and the eigenvalues of the matrix
corresponding to corrections to active neutrino species in Fig. 4. The
eigenvalues of F themselves exceed the active neutrino tree-level masses, as
the latter are at most about 10 eV for the MeV scale $Z^{\prime}$ boson.
However, the coupling suppressions in Eq. (IV.3) are sufficient to tame the
relative correction to the tree level mass below the per cent level. Assuming
the active neutrino masses to be $\mathrm{O}(10^{-3})$ eV, a rough estimate
for the relative correction to active neutrino masses is of
$\mathrm{O}(10^{-4})$.

We may maximize the effect of $Z^{\prime}$ loop by allowing $Z^{\prime}$ mass
to be free and setting large $|\sin\theta_{G}|=\mathrm{O}(10^{-1})$, which is
obtained when $g_{y}^{\prime}$ and $g_{z}^{\prime}$ are O($10^{-1}$). This
corresponds to $M_{Z^{\prime}}=\mathrm{O}(M_{Z})$, which is of course,
excluded. Yet, even in this case, the correction from $Z$ and $Z^{\prime}$
loops are small, have the same order of magnitude, $\mathrm{O}(10^{-7})$ eV.
Thus, the individual contributions from BSM loops cannot be significantly
larger than the SM contributions.

Figure 3:  Matrix elements Fij as a function of the mass $m_{\text{loop}}$ of
the boson in the loop are confined to the blue band, assuming normal neutrino
mass hierarchy. We have highlighted with vertical bands the relevant mass
regions where the masses of the bosons in the loop lie. The scalar $s$ is
required to have mass between 144 and 558 GeV requiring stability of the
vacuum [38]. Left plot: $m_{1}^{\text{tree}}=0.01$ eV,
$m_{4}^{\text{tree}}=30$ keV, $m_{5}^{\text{tree}}\approx
m_{6}^{\text{tree}}=2.5$ GeV; right plot: $m_{1}^{\text{tree}}=0.001$ eV,
$m_{4}^{\text{tree}}=7.1$ keV, $m_{5}^{\text{tree}}\approx
m_{6}^{\text{tree}}=3.0$ GeV. Figure 4:  Eigenvalues of the matrix F as a
function of the mass of the boson in the loop $m_{\text{loop}}$, assuming
$m_{1}^{\text{tree}}=0.01$ eV, $m_{4}^{\text{tree}}=30$ keV,
$m_{5}^{\text{tree}}\approx m_{6}^{\text{tree}}=2.5$ GeV, and normal neutrino
mass hierarchy.

##  V Conclusions

In this paper, we have computed the one-loop corrections to the mass matrix of
the active neutrinos in a gauged U(1) extension of the standard model of
particle interactions. The field content of the model consists of a new
complex scalar field and three right-handed neutrinos—sterile under the
standard model interactions—in addition to the fields in the standard model.
The neutrino masses are generated by Dirac and Majorana type Yukawa terms,
which after spontaneous symmetry breaking of both scalar fields give rise to
neutrino masses in the way of the type I see-saw mass generation. We used
$R_{\xi}$ gauge and have shown that the one-loop corrections are (i)
independent of the gauge fixing parameters, (ii) finite and (iii) independent
of the regularization scale. We also demonstrated how the formula for the one-
loop mass corrections can be generalized to the case of arbitrary number of
new U(1) groups, complex scalars and right-handed neutrinos.

We have provided a numerical estimate of the size of the mass corrections in
the context of the super-weak model, in which the new neutral gauge boson
$Z^{\prime}$ is much lighter than the $Z$ boson of the standard model. We have
found that in the mass range of $M_{Z^{\prime}}\in[20,200]$ MeV, motivated by
a possible explanation of the relic density of dark matter in the Universe,
the relative mass corrections to the tree-level mass matrix elements do not
exceed the per mill level. Hence the model is stable against higher-order
corrections in the neutrino sector, which motivates further studies to explore
the viable parameter space of the model regarding the mixings between the
active and sterile neutrinos [28].

### Acknowledgments

We are grateful to Josu Hernández-García for discussions on this project. This
work was supported by grant K 125105 of the National Research, Development and
Innovation Fund in Hungary.

##  Appendix A Some properties of neutrino mass and mixing matrices

In this appendix, we derive some useful relations among the neutrino mass and
mixing matrices. The matrix U that diagonalizes the neutrino mass matrix is
unitary, hence

$\textbf{UU}^{\dagger}=\binom{\textbf{U}_{L}}{\textbf{U}_{R}^{*}}(\textbf{U}_{L}^{\dagger},\textbf{U}_{R}^{T})=\begin{pmatrix}\textbf{U}_{L}\textbf{U}_{L}^{\dagger}&\textbf{U}_{L}\textbf{U}_{R}^{T}\\\
\textbf{U}_{R}^{*}\textbf{U}_{L}^{\dagger}&\textbf{U}_{R}^{*}\textbf{U}_{R}^{T}\end{pmatrix}=\begin{pmatrix}{\mbox{\boldmath$1$}}_{3}&\textbf{0}_{3}\\\
\textbf{0}_{3}&{\mbox{\boldmath$1$}}_{3}\end{pmatrix},$ (A.1)

from which we obtain the following important relations:

$\textbf{U}_{L}\textbf{U}_{L}^{\dagger}={\mbox{\boldmath$1$}}_{3}\,,\quad\textbf{U}_{R}\textbf{U}_{R}^{\dagger}={\mbox{\boldmath$1$}}_{3}$
(A.2)

and

$\textbf{U}_{L}\textbf{U}_{R}^{T}=\textbf{U}_{R}^{*}\textbf{U}_{L}^{\dagger}=\textbf{0}_{3}$
(A.3)

where ${\mbox{\boldmath$1$}}_{n}$ denotes the $n\times n$ unit matrix. The
second unitarity conditions gives

$\textbf{U}^{\dagger}\textbf{U}=(\textbf{U}_{L}^{\dagger},\textbf{U}_{R}^{T})\binom{\textbf{U}_{L}}{\textbf{U}_{R}^{*}}=\textbf{U}_{L}^{\dagger}\textbf{U}_{L}+\textbf{U}_{R}^{T}\textbf{U}_{R}^{*}={\mbox{\boldmath$1$}}_{6}\,.$
(A.4)

Using Eq. (II.35) we derive

$\displaystyle\textbf{U}_{L}^{*}\textbf{M}\textbf{U}_{L}^{\dagger}$
$\displaystyle=\textbf{U}_{L}^{*}(\textbf{U}_{L}^{T},\textbf{U}_{R}^{\dagger})\begin{pmatrix}\textbf{0}_{3}&\textbf{M}_{D}^{T}\\\
\textbf{M}_{D}&\textbf{M}_{N}\end{pmatrix}\binom{\textbf{U}_{L}}{\textbf{U}_{R}^{*}}\textbf{U}_{L}^{\dagger}$
$\displaystyle=\textbf{U}_{L}^{*}\textbf{U}_{L}^{T}\textbf{M}_{D}\textbf{U}_{R}^{*}\textbf{U}_{L}^{\dagger}+\textbf{U}_{L}^{*}\textbf{U}_{R}^{\dagger}\textbf{M}_{D}\textbf{U}_{L}\textbf{U}_{L}^{\dagger}+\textbf{U}_{L}^{T}\textbf{U}_{R}^{\dagger}\textbf{M}_{N}\textbf{U}_{R}^{*}\textbf{U}_{L}^{\dagger}\,,$

and then with relations in Eq. (A.3) we obtain

$\textbf{U}_{L}^{*}\textbf{M}\textbf{U}_{L}^{\dagger}={\mbox{\boldmath$0$}}_{3}\,.$
(A.5)

Analogous calculations yield

$\textbf{U}_{R}\textbf{M}\textbf{U}_{L}^{\dagger}=\textbf{M}_{D}\,.$ (A.6)

Multiplying Eq. (A.6) with $\textbf{U}_{R}^{\dagger}$ from the left and using
Eq. (A.4), we find

$\textbf{U}_{R}^{\dagger}\textbf{M}_{D}=\textbf{U}_{R}^{\dagger}\textbf{U}_{R}\textbf{M}\textbf{U}_{L}^{\dagger}=({\mbox{\boldmath$1$}}_{6}-\textbf{U}_{L}^{T}\textbf{U}_{L}^{*})\textbf{M}\textbf{U}_{L}^{\dagger}=\textbf{M}\textbf{U}_{L}^{\dagger}-\textbf{U}_{L}^{T}\textbf{U}_{L}^{*}\textbf{M}\textbf{U}_{L}^{\dagger}$
(A.7)

where the last term vanishes by Eq. (A.5), so

$\textbf{U}_{R}^{\dagger}\textbf{M}_{D}=\textbf{M}\textbf{U}_{L}^{\dagger}\,.$
(A.8)

Finally,

$\displaystyle\textbf{U}_{R}\textbf{M}\textbf{U}_{R}^{T}$
$\displaystyle=\textbf{U}_{R}(\textbf{U}_{L}^{T},\textbf{U}_{R}^{\dagger})\begin{pmatrix}\textbf{0}&\textbf{M}_{D}^{T}\\\
\textbf{M}_{D}&\textbf{M}_{N}\end{pmatrix}\binom{\textbf{U}_{L}}{\textbf{U}_{R}^{*}}\textbf{U}_{R}^{T}$
$\displaystyle=\textbf{U}_{R}\Big{(}\textbf{U}_{L}^{T}\textbf{M}_{D}^{T}\textbf{U}_{R}^{*}+\textbf{U}_{R}^{\dagger}\textbf{M}_{D}\textbf{U}_{L}+\textbf{U}_{R}^{\dagger}\textbf{M}_{N}\textbf{U}_{R}^{*}\Big{)}\textbf{U}_{R}^{T}$

Expanding the factors into the parenthesis, the first two terms give vanishing
contribution by Eq. (A.4), while utilizing Eq. (A.2), the last one is simply
$\textbf{M}_{N}$, so

$\textbf{U}_{R}\textbf{M}\textbf{U}_{R}^{T}=\textbf{M}_{N}\,.$ (A.9)

Now Eq. (A.4) allows us to derive

$\textbf{U}_{R}^{\dagger}\textbf{M}_{N}=\textbf{U}_{R}^{\dagger}\textbf{U}_{R}\textbf{M}\textbf{U}_{R}^{T}=({\mbox{\boldmath$1$}}_{6}-\textbf{U}_{L}^{T}\textbf{U}_{L}^{*})\textbf{M}\textbf{U}_{R}^{T}=\textbf{M}\textbf{U}_{R}^{T}-\textbf{U}_{L}^{T}\textbf{U}_{L}^{*}\textbf{M}\textbf{U}_{R}^{T}\,.$
(A.10)

where the second term on the right does not vanish this time.

##  Appendix B Evaluation of the vector boson exchange diagram

The vector boson exchange diagrams, shown in the bottom row of Fig. 1
contribute gauge-dependent terms to the neutrino self energy. In order to show
that the gauge-dependent terms cancel once contributions from all particles
are considered, it is useful to eliminate the loop momentum from the one loop
integral corresponding to the vector boson exchange diagram. A decomposition
to achieve this was used in Ref. [39], and we shall derive it here as well. In
Ref. [40], Eq (4.4) contains the self energy:

$\mathrm{i}\Sigma_{V}(p)=-\int\frac{d^{d}\ell}{(2\pi)^{d}}{\mbox{\boldmath$\Gamma$}}^{\mu\dagger}\textbf{P}(p-\ell){\mbox{\boldmath$\Gamma$}}^{\nu}P_{\mu\nu}(\ell,M_{V}^{2};\xi_{V}),$
(B.1)

where the $6\times 6$ matrices are defined as follows. The matrix P is the
fermion propagator, diagonal in the mass eigenstates,

$\textbf{P}(p-\ell)=[(\cancel{p}-\cancel{\ell}){\mbox{\boldmath$1$}}-\textbf{M}]^{-1}\,,\text{~{}~{}while~{}~{}}{\mbox{\boldmath$\Gamma$}}^{\mu}=-\mathrm{i}e\gamma^{\mu}\textbf{A}\,,\text{~{}~{}with~{}~{}}\textbf{A}={\mbox{\boldmath$\Gamma$}}^{L}P_{L}+{\mbox{\boldmath$\Gamma$}}^{R}P_{R}.$
(B.2)

In the following we shall write P for $\textbf{P}(p-\ell)$. The matrix A is
self-adjoint, $\textbf{A}^{\dagger}=\textbf{A}$, and so is
${\mbox{\boldmath$\Gamma$}}^{L/R}$. We also introduce the abbreviation

$\tilde{\textbf{A}}={\mbox{\boldmath$\Gamma$}}^{R}P_{L}+{\mbox{\boldmath$\Gamma$}}^{L}P_{R}\,,$
(B.3)

which will simplify our calculations. In order to compute the loop integral
easily containing the neutral vector boson propagator in the neutrino self-
energy loop, in this Appendix we perform tensor reduction of the matrix
product

$\cancel{\ell}\textbf{A}^{\dagger}\textbf{P}\cancel{\ell}\textbf{A}$ (B.4)

such that the numerator factor be at most linear in the loop momentum $\ell$.

When the fermion momentum $p$ appears as $\cancel{p}$ at the extreme left or
right of the expression, it satisfies the Dirac equation
$\cancel{p}{\mbox{\boldmath$1$}}_{6}=\textbf{M}$ (both Dirac and Majorana
fermions do so), thus we can replace formally $\cancel{p}$ with M. Let us
first write the identity

$\displaystyle\cancel{\ell}\textbf{A}^{\dagger}$
$\displaystyle=(\cancel{\ell}-\cancel{p}+\textbf{M})\textbf{A}^{\dagger}=\textbf{MA}^{\dagger}-(\cancel{p}-\cancel{\ell})\textbf{A}^{\dagger}.$
(B.5)

The chiral coupling matrix A anticommutes with the Dirac matrices
$\gamma_{\mu}$, hence

$\displaystyle\cancel{\ell}\textbf{A}^{\dagger}=\textbf{M}\textbf{A}^{\dagger}-\tilde{\textbf{A}}^{\dagger}\Big{(}(\cancel{p}-\cancel{\ell}){\mbox{\boldmath$1$}}-\textbf{M}+\textbf{M}\Big{)}$
$\displaystyle=\textbf{M}\textbf{A}^{\dagger}-\tilde{\textbf{A}}^{\dagger}\textbf{P}^{-1}-\tilde{\textbf{A}}^{\dagger}\textbf{M}\,,$
(B.6)

and similarly,

$\displaystyle\textbf{P}\cancel{\ell}\textbf{A}=-\textbf{A}-\textbf{P}\textbf{MA}+\textbf{P}\tilde{\textbf{A}}\textbf{M}\,.$
(B.7)

Multiplying Eqs. (B.6) and (B.7), we obtain the expression (B.4), and its
expansion yields

$\begin{split}\cancel{\ell}\textbf{A}^{\dagger}\textbf{P}\cancel{\ell}\textbf{A}=&-\textbf{MA}^{\dagger}\textbf{A}-\textbf{MA}^{\dagger}\textbf{P}\textbf{MA}+\textbf{MA}^{\dagger}\textbf{P}\textbf{AM}+\tilde{\textbf{A}}^{\dagger}\textbf{P}^{-1}\textbf{A}\\\
&+\tilde{\textbf{A}}^{\dagger}\textbf{MA}-\tilde{\textbf{A}}^{\dagger}\tilde{\textbf{A}}\textbf{M}+\tilde{\textbf{A}}^{\dagger}\textbf{MA}+\tilde{\textbf{A}}^{\dagger}\textbf{M}\textbf{P}\textbf{MA}-\tilde{\textbf{A}}^{\dagger}\textbf{M}\textbf{P}\tilde{\textbf{A}}\textbf{M}\,.\end{split}$
(B.8)

Using that $\cancel{p}\textbf{A}=\tilde{\textbf{A}}\cancel{p}$, the fourth
term can rearranged as

$\begin{split}\tilde{\textbf{A}}^{\dagger}\textbf{P}^{-1}\textbf{A}&=\tilde{\textbf{A}}^{\dagger}\Big{(}(\cancel{p}-\cancel{\ell})\textbf{1}-\textbf{M}\Big{)}\textbf{A}=\frac{1}{2}\tilde{\textbf{A}}^{\dagger}\cancel{p}\textbf{A}+\frac{1}{2}\tilde{\textbf{A}}^{\dagger}\cancel{p}\textbf{A}-\tilde{\textbf{A}}^{\dagger}\cancel{\ell}\textbf{A}-\tilde{\textbf{A}}^{\dagger}\textbf{MA}\\\
&=\frac{1}{2}\cancel{p}\textbf{A}^{\dagger}\textbf{A}+\frac{1}{2}\tilde{\textbf{A}}^{\dagger}\tilde{\textbf{A}}\cancel{p}-\tilde{\textbf{A}}^{\dagger}\cancel{\ell}\textbf{A}-\tilde{\textbf{A}}^{\dagger}\textbf{MA}\end{split}$
(B.9)

The $\cancel{p}$ is on extreme left and right, hence can be replaced with M,
giving

$\tilde{\textbf{A}}^{\dagger}\textbf{P}^{-1}\textbf{A}=\frac{1}{2}\textbf{M}\textbf{A}^{\dagger}\textbf{A}+\frac{1}{2}\tilde{\textbf{A}}^{\dagger}\tilde{\textbf{A}}\textbf{M}-\tilde{\textbf{A}}^{\dagger}\cancel{\ell}\textbf{A}-\tilde{\textbf{A}}^{\dagger}\textbf{MA}\,.$
(B.10)

Substituting Eq. (B.10) into Eq. (B.8), we obtain

$\cancel{\ell}\textbf{A}^{\dagger}\textbf{P}\cancel{\ell}\textbf{A}=\textbf{M}_{0}+\textbf{M}_{1}+\textbf{M}_{2}$
(B.11)

where we introduced the abbreviations

$\begin{split}\textbf{M}_{0}&=-\tilde{\textbf{A}}^{\dagger}\cancel{\ell}\textbf{A}\,,\\\
\textbf{M}_{1}&=-\frac{1}{2}\textbf{M}\textbf{A}^{\dagger}\textbf{A}-\frac{1}{2}\tilde{\textbf{A}}^{\dagger}\tilde{\textbf{A}}\textbf{M}+\tilde{\textbf{A}}^{\dagger}\textbf{MA}\end{split}$
(B.12)

and

$\begin{split}\textbf{M}_{2}&=-\textbf{MA}^{\dagger}\textbf{P}\textbf{MA}+\textbf{MA}^{\dagger}\textbf{P}\tilde{\textbf{A}}\textbf{M}+\tilde{\textbf{A}}^{\dagger}\textbf{M}\textbf{P}\textbf{MA}-\tilde{\textbf{A}}^{\dagger}\textbf{M}\textbf{P}\tilde{\textbf{A}}\textbf{M}\\\
&=(\textbf{MA}^{\dagger}-\tilde{\textbf{A}}^{\dagger}\textbf{M})\textbf{P}(\tilde{\textbf{A}}\textbf{M}-\textbf{MA})\,,\end{split}$
(B.13)

which correspond to constant, linear and quadratic terms in the neutrino mass
matrix M. We now discuss the contribution from each term in Eq. (B.11)
separately.

The first, constant term gives vanishing contribution to the loop integral as
it is odd in the loop momentum. The other two terms can be decomposed into
left and right chiral pieces:

$\textbf{M}_{i}=\textbf{M}_{i}^{L}P_{L}+\textbf{M}_{i}^{R}P_{R}\,.$ (B.14)

Our goal is to compute the one-loop correction (III.7) to the tree-level mass
matrix of the light neutrinos. In order to obtain it, one sandwiches the left
handed pieces $\textbf{M}_{i}^{L}$ between the matrices $\textbf{U}_{L}^{*}$
and $\textbf{U}_{L}^{\dagger}$. Using the properties of the neutrino mixing
matrices of Appendix A, we immediately see that

$\textbf{U}_{L}^{*}\textbf{M}_{1}^{L}\textbf{U}_{L}^{\dagger}=0\,,$ (B.15)

while lengthy computations yield

$\begin{split}\textbf{U}_{L}^{*}\textbf{M}_{2}^{L}\textbf{U}_{L}^{\dagger}=&-(C^{L}_{V\nu\nu}-C^{R}_{V\nu\nu})^{2}\textbf{U}_{L}^{*}\textbf{M}\textbf{P}\textbf{MU}_{L}^{\dagger}\\\
&+\text{terms that do not contribute to {\bf B}}_{L}(0)\,.\end{split}$ (B.16)

Here we outline the steps needed to reach Eq. (B.16).

Firstly, in order to find the left-chiral part $\textbf{M}_{2}^{L}$, we
substitute A and $\tilde{\textbf{A}}$ into Eq. (B.13). We write the
denominator of the fermion propagator as

$(\textbf{P})_{ij}=\delta_{ij}\frac{\not{p}-\not{\ell}+m_{i}}{(p-\ell)^{2}-m_{i}^{2}}$
(B.17)

and use the following relations for the Dirac projectors:

$P_{L/R}(\not{q}+m)P_{L/R}=mP_{L/R}\,,\quad
P_{L/R}(\not{q}+m)P_{R/L}=\not{q}P_{R/L}\,,$ (B.18)

valid for any momentum $q$ and mass $m$. Hence

$\begin{split}P_{L/R}(\textbf{P})_{ij}P_{L/R}&=\left((\textbf{P})_{ij}-\delta_{ij}\frac{\not{p}-\not{\ell}}{(p-\ell)^{2}-m_{i}^{2}}\right)P_{L/R}\,,\\\
P_{L/R}(\textbf{P})_{ij}P_{R/L}&=\delta_{ij}\frac{\not{p}-\not{\ell}}{(p-\ell)^{2}-m_{i}^{2}}P_{R/L}\,,\end{split}$
(B.19)

and therefore, we obtain

$\textbf{M}_{2}^{L}=-\textbf{M}{\mbox{\boldmath$\Gamma$}}^{L\dagger}\textbf{P}\textbf{M}{\mbox{\boldmath$\Gamma$}}^{L}+\textbf{M}{\mbox{\boldmath$\Gamma$}}^{L\dagger}\textbf{P}{\mbox{\boldmath$\Gamma$}}^{R}\textbf{M}+{\mbox{\boldmath$\Gamma$}}^{R\dagger}\textbf{M}\textbf{P}\textbf{M}{\mbox{\boldmath$\Gamma$}}^{L}-{\mbox{\boldmath$\Gamma$}}^{R\dagger}\textbf{M}\textbf{P}{\mbox{\boldmath$\Gamma$}}^{R}\textbf{M}+\textbf{D}$
(B.20)

where the last term is proportional to $(\not{p}-\not{\ell})$:

$\textbf{D}=[\textbf{M}({\mbox{\boldmath$\Gamma$}}^{L\dagger}-{\mbox{\boldmath$\Gamma$}}^{R\dagger})-({\mbox{\boldmath$\Gamma$}}^{L\dagger}-{\mbox{\boldmath$\Gamma$}}^{R\dagger})\textbf{M})](\not{p}-\not{\ell})[(p-l)^{2}{\mbox{\boldmath$1$}}_{6}-\textbf{M}^{2}]^{-1}({\mbox{\boldmath$\Gamma$}}^{R}\textbf{M}-\textbf{M}{\mbox{\boldmath$\Gamma$}}^{L})\,.$

Then using the matrix relations derived in Appendix A, we can compute the
following identities:

$\displaystyle\textbf{U}_{L}^{*}\textbf{M}{\mbox{\boldmath$\Gamma$}}^{L\dagger}$
$\displaystyle=-C_{V\nu\nu}^{R}\textbf{U}_{L}^{*}\textbf{M}\,,\quad\textbf{U}_{L}^{*}{\mbox{\boldmath$\Gamma$}}^{R\,\dagger}=-C^{L}_{V\nu\nu}\textbf{U}_{L}^{*}\,,$
(B.21) $\displaystyle{\mbox{\boldmath$\Gamma$}}^{R}\textbf{MU}_{L}^{\dagger}$
$\displaystyle=C^{R}_{V\nu\nu}\textbf{MU}_{L}^{\dagger}\,,\qquad\;{\mbox{\boldmath$\Gamma$}}^{L}\textbf{U}_{L}^{\dagger}=C^{L}_{V\nu\nu}\textbf{U}_{L}^{\dagger}\,.$
(B.22)

Finally sandwiching Eq. (B.20) gives us

$\begin{split}\textbf{U}_{L}^{*}\textbf{M}_{2}^{L}\textbf{U}_{L}^{\dagger}&=-(C^{L}_{V\nu\nu}-C^{R}_{V\nu\nu})^{2}\textbf{U}_{L}^{*}\textbf{M}\textbf{P}\textbf{MU}_{L}^{\dagger}+\textbf{U}_{L}^{*}\textbf{D}\textbf{U}_{L}^{\dagger}\,.\end{split}$
(B.23)

As mentioned, the last term is proportional to $(\not{p}-\not{\ell})$, but
only the term with $\not{\ell}$ contributes to $\textbf{B}_{L}(p=0)$. That
piece, being an odd function of $\ell$, vanishes upon integration, which
completes the proof of Eq. (B.16).

The charged vector bosons $W^{\pm}$ also contribute to the neutrino self-
energy. The corresponding Feynman rules are

${\mbox{\boldmath$\Gamma$}}^{\mu}_{W^{-}\bar{\ell}\nu}=-\mathrm{i}e\gamma^{\mu}{\mbox{\boldmath$\Gamma$}}^{L}_{W^{-}\bar{\ell}\nu}P_{L},\quad{\mbox{\boldmath$\Gamma$}}^{\mu}_{W^{+}\nu\ell}=-\mathrm{i}e\gamma^{\mu}{\mbox{\boldmath$\Gamma$}}^{L}_{W^{+}\nu\ell}P_{L}.$
(B.24)

where

${\mbox{\boldmath$\Gamma$}}^{L}_{W^{-}\bar{\ell}\nu}=C_{W\ell\nu}(\textbf{U}_{L}^{\ell}\textbf{U}_{L})_{ij}$
(B.25)

with $\textbf{U}_{L}^{\ell}$ being the charged lepton mixing matrix and
${\mbox{\boldmath$\Gamma$}}^{L}_{W^{-}\bar{\ell}\nu}=({\mbox{\boldmath$\Gamma$}}^{L}_{W^{+}\nu\ell})^{\dagger}$.
The charged vector boson contribution to (III.7) is proportional to
$\textbf{U}_{L}^{*}\textbf{M}\textbf{U}_{L}^{\dagger}$, which vanishes
identically as shown in Appendix A.

## References

  * [1] Y. Fukuda et al.  Measurements of the solar neutrino flux from super-kamiokande’s first 300 days.  Phys. Rev. Lett., 81:1158–1162, Aug 1998.
  * [2] Q. R. Ahmad et al.  Measurement of the rate of ${\nu}_{e}+d\rightarrow p+p+e^{-}$ interactions produced by ${}^{8}b$ solar neutrinos at the sudbury neutrino observatory.  Phys. Rev. Lett., 87:071301, Jul 2001.
  * [3] H. Fritzsch, M. Gell-Mann, and P. Minkowski.  Vectorlike weak currents and new elementary fermions.  Physics Letters B, 59(3):256 – 260, 1975.
  * [4] Peter Minkowski.  $\mu\to e\gamma$ at a Rate of One Out of $10^{9}$ Muon Decays?  Phys. Lett., 67B:421–428, 1977.
  * [5] Murray Gell-Mann, Pierre Ramond, and Richard Slansky.  Complex Spinors and Unified Theories.  Conf. Proc., C790927:315–321, 1979.
  * [6] Tsutomu Yanagida.  Horizontal symmetry and masses of neutrinos.  Conf. Proc., C7902131:95–99, 1979.
  * [7] Rabindra N. Mohapatra and Goran Senjanović.  Neutrino mass and spontaneous parity nonconservation.  Phys. Rev. Lett., 44:912–915, Apr 1980.
  * [8] J. Schechter and J. W. F. Valle.  Neutrino Masses in SU(2) $\times$ U(1) Theories.  Phys. Rev., D22:2227, 1980.
  * [9] M. Magg and C. Wetterich.  Neutrino Mass Problem and Gauge Hierarchy.  Phys. Lett., 94B:61–64, 1980.
  * [10] S. L. Glashow.  The Future of Elementary Particle Physics.  NATO Sci. Ser. B, 61:687, 1980.
  * [11] B. Abi et al.  Measurement of the Positive Muon Anomalous Magnetic Moment to 0.46 ppm.  Phys. Rev. Lett., 126(14):141801, 2021.
  * [12] A. A. Aguilar-Arevalo et al.  Observation of a Significant Excess of Electron-Like Events in the MiniBooNE Short-Baseline Neutrino Experiment.  2018\.
  * [13] Ann E. Nelson and Jonathan Walsh.  Short Baseline Neutrino Oscillations and a New Light Gauge Boson.  Phys. Rev. D, 77:033001, 2008.
  * [14] Julian Heeck and Werner Rodejohann.  Gauged $L_{\mu}-L_{\tau}$ and different Muon Neutrino and Anti-Neutrino Oscillations: MINOS and beyond.  J. Phys. G, 38:085005, 2011.
  * [15] Ernest Ma.  Gauged B - 3L(tau) and radiative neutrino masses.  Phys. Lett. B, 433:74–81, 1998.
  * [16] Ernest Ma, D. P. Roy, and Sourov Roy.  Gauged L(mu) - L(tau) with large muon anomalous magnetic moment and the bimaximal mixing of neutrinos.  Phys. Lett. B, 525:101–106, 2002.
  * [17] Kento Asai.  Predictions for the neutrino parameters in the minimal model extended by linear combination of U(1)${}_{L_{e}-L_{\mu}}$, U(1)${}_{L_{\mu}-L_{\tau}}$ and U(1)B-L gauge symmetries.  Eur. Phys. J. C, 80(2):76, 2020.
  * [18] Disha Bhatia, Sabyasachi Chakraborty, and Amol Dighe.  Neutrino mixing and $R_{K}$ anomaly in U(1)X models: a bottom-up approach.  JHEP, 03:117, 2017.
  * [19] Rathin Adhikari, Jens Erler, and Ernest Ma.  Seesaw Neutrino Mass and New U(1) Gauge Symmetry.  Phys. Lett. B, 672:136–140, 2009.
  * [20] Debasish Borah, Lopamudra Mukherjee, and Soumitra Nandi.  Low scale U(1)X gauge symmetry as an origin of dark matter, neutrino mass and flavour anomalies.  JHEP, 12:052, 2020.
  * [21] D. Aristizabal Sierra and Carlos E. Yaguna.  On the importance of the 1-loop finite corrections to seesaw neutrino masses.  JHEP, 08:013, 2011.
  * [22] J. Lopez-Pavon, S. Pascoli, and Chan-fai Wong.  Can heavy neutrinos dominate neutrinoless double beta decay?  Phys. Rev. D, 87(9):093007, 2013.
  * [23] Walter Grimus and Luis Lavoura.  One-loop corrections to the seesaw mechanism in the multi-Higgs-doublet standard model.  Phys. Lett. B, 546:86–95, 2002.
  * [24] W. Grimus and M. Löschner.  Renormalization of the multi-Higgs-doublet Standard Model and one-loop lepton mass corrections.  JHEP, 11:087, 2018.
  * [25] Vytautas Dūdėnas and Thomas Gajdosik.  Gauge dependence of tadpole and mass renormalization for a seesaw extended 2HDM.  Phys. Rev. D, 98(3):035034, 2018.
  * [26] Ansgar Denner, Laura Jenniches, Jean-Nicolas Lang, and Christian Sturm.  Gauge-independent $\overline{MS}$ renormalization in the 2HDM.  JHEP, 09:115, 2016.
  * [27] Zoltán Trócsányi.  Super-weak force and neutrino masses.  Symmetry, 12(1):107, 2020.
  * [28] Timo J. Kärkkäinen and Zoltán Trócsányi.  prepared for submission.
  * [29] Sho Iwamoto, Károly Seller, and Zoltán Trócsányi.  Sterile neutrino dark matter in a U(1) extension of the standard model.  4 2021.
  * [30] Zoltán Péli, István Nándori, and Zoltán Trócsányi.  Particle physics model of curvaton inflation in a stable universe.  Phys. Rev. D, 101(6):063533, 2020.
  * [31] F. Staub.  Sarah, 2012.
  * [32] Florian Staub.  From superpotential to model files for feynarts and calchep/comphep.  Computer Physics Communications, 181(6):1077–1086, 2010.
  * [33] Florian Staub.  Automatic calculation of supersymmetric renormalization group equations and loop corrections.  Computer Physics Communications, 182(3):808–833, 2011.
  * [34] Florian Staub.  Sarah 4: A tool for (not only susy) model builders.  Computer Physics Communications, 185(6):1773–1790, 2014.
  * [35] We present the model files in a separate file SuperWeak.zip.
  * [36] Note the opposite sign convention for $\theta_{Z}$ in this work and in Ref. [27] where this mixing angle was denoted as $\theta_{T}$, so $\theta_{Z}=-\theta_{T}$.
  * [37] This matrix coincides with the unitary matrix $U^{V,\dagger}$ used by SARAH.
  * [38] Zoltán Péli and Zoltán Trócsányi.  Stability of the vacuum as constraint on $U$(1) extensions of the standard model.  2 2019.
  * [39] Maximilian Löschner.  Renormalization and one-loop corrections of lepton masses.  PhD thesis, Vienna U., 2018.
  * [40] Steven Weinberg.  Perturbative Calculations of Symmetry Breaking.  Phys. Rev. D, 7:2887–2910, 1973.