# Global biasing using a Hardware-based artificial Zeeman term in Spinwave
Ising Machines

Victor H. González1, Artem Litvinenko1, Roman Khymyn1, and Johan Åkerman1
1Department of Physics, University of Gothenburg, Gothenburg, 41296, Sweden,
<EMAIL_ADDRESS>and<EMAIL_ADDRESS>
###### Abstract

A spinwave Ising machine (SWIM) is a newly proposed type of time-multiplexed
hardware solver for combinatorial optimization that employs feedback coupling
and phase sensitive amplification to map an Ising Hamiltonian into phase-
binarized propagating spin-wave RF pulses in an Yttrium-Iron-Garnet (YIG)
film. In this work, we increase the mathematical complexity of the SWIM by
adding a global Zeeman term to a 4-spin MAX-CUT Hamiltonian using a continuous
external electrical signal with the same frequency as the spin pulses and
phase locked with with one of the two possible states. We are able to induce
ferromagnetic ordering in both directions of the spin states despite
antiferromagnetic pairwise coupling. Embedding a planar antiferromagnetic spin
system in a magnetic field has been proven to increase the complexity of the
graph associated to its Hamiltonian and thus this straightforward
implementation helps explore higher degrees of complexity in this evolving
solver.

###### Index Terms:

combinatorial optimization problems, Ising machines, spinwaves, unconventional
computing, physical computing, spinwaves.

In the landscape of physical computation schemes, Ising machines (IM) have
attracted considerable attention and investment over the last decade, in both
academic and industrial research[1, 2, 3, 4, 5, 6], for their applicability to
combinatorial optimization problems, potential for scalability and progressive
increase in mathematical complexity. In this work, we contribute to the latter
as we implement a global bias to artificial spin states in a spinwave Ising
machine (SWIM).

A SWIM is a newly proposed [2] time-multiplexed hardware solver circuit for
problems in the NP (nonpolynomial time) complexity class that employs feedback
coupling and phase sensitive amplification to map an Ising Hamiltonian into
spinwave (SW) RF pulses propagating in an Yttrium-Iron-Garnet (YIG) film.
Spinwaves are suitable for Ising machines because of their GHz oscillation
frequencies, which permits the development of multiphysical systems using
cheap and efficient off-the-shelf microwave components for signal processing
[7, 8] and amplification which results in a small circuitry footprint and high
per-spin power efficiency. In this work, we use an external continuous
microwave signal to implement a global biasing to propagating spinwave RF
artificial spin states. Exploring the complexity limits of this circuit
implementation thus holds significant interest for technical and commercial
applications.

An IM operates as an mapping of an objective function into the Ising
Hamiltonian of a device composed of an array of $N$ binarized physical units
referred to as spins $s_{i}=\pm 1$:

$H=-\frac{1}{2}\sum_{i=1}^{N}\sum_{j=1}^{N}J_{ij}s_{i}s_{j}-\sum_{i=1}^{N}h_{i}s_{i}$
(1)

The objective function associated with the NP problem to solve is encoded into
the pairwise coupling $J_{ij}$ and external Zeeman bias $h_{i}$ such that
ground state of the system represents the solution to it. Combinatorial
problems of practical use, such as the traveling salesman and knapsack with
integer weights, have been shown to be encodable if one can add constraints to
the degrees of freedom of the system [9]. Adding a global bias (i.e. the
$h_{i}$ is identical for all $s_{i}$) to an antiferromagnetic planar graph has
been shown as a straightforward way to increase its mathematical complexity
[10]. Implementation of a global Zeeman term using software has also shown to
improve the stability and performance of time-multiplexed IMs with frustated
lattices [11]. Thus, exploring hardware-based alternative schemes can lend
versatility to small scale low-power devices.

Figure 1: Zeeman-biased spinwave Ising machine. PSA and LNA stand for phase-
sensitive and low noise amplifiers, respectively. The propagating RF pulses
have frequency of 3.13 GHz. The sign of the coupling is controlled by the
total phase accumulation in the coupling delay and the Zeeman field amplitude
and sign is controlled by the amplitude and phase of the injected signal
$\omega_{bias}$. Figure 2: Influence of Zeeman term at different signs of
magnetic field. (a) Measured potential landscape shift of the artificial spin
states. The sign of $h$ is given by the phase between $\omega_{ref}$ and
$\omega_{bias}$. The total phase sensitive amplification of the circuit
$A_{PSA}$ favors either 0 or $\pi$ depending on $h$ and thus changes the
magnetic ordering of the spins. (b) Power spectral densities (PSD) as a
functions of frequency for different $h$ signs. We observe that the modulation
harmonics can be used to identify the type of solution achieved, with both
ferromagnetic states having the same spectra. (c) Time traces of the RF pulses
colored with their respective instantaneous phase for different signs of $h$.
The non-zero signals have an amplitude of -15 dBm. $\omega_{bias}$ allows us
to change the ordering and direction of the artificial spin state.

The circuit of the SWIM with a hardware implemented Zeeman term is shown in
Fig. 1. The SWIM’s construction of the artificial spin state relies in phase-
sensitive amplification (PSA) of a reference signal $\omega_{ref}=$3.13 GHz.
The PSA binarizes the phase of the signal by amplifying only its in-phase
($\phi=0$) and out-of-phase ($\phi=\pi$) components. The signal is then
simultaneously injected into a YIG waveguide and a coaxial delay line using
coupler C2. The electrical signal excites spin waves within waveguide that
propagate at a much slower speed, allowing the coupling delay to re-inject a
shifted signal using coupler C1. Switch 1 (S1) then pulses the signal and a
low noise amplifier (LNA) compensates propagation losses as the cycle starts
again. The resulting time-multiplexed pulse train is our artificial spin
state, where each pulse is an artificial spin, with their electrical
amplitudes representing the norm of $s_{i}$ and their individual phases
representing its sign. The coupling delay’s length is such that every pulse
interferes (or couples) with the previous nearest neighbor. Phase shifter PS
ensures that the coupling term $J_{i,i+1}$ is antiferromagnetic and variable
attenuator VA controls its strength. The Zeeman term is implemented with
external signal $\omega_{bias}$ applied to the propagating pulses after PSA.

The role of $\omega_{bias}$ is to unbalance the potential landscape of the
artificial spin states and favor one phase over the other. Fig.2(a) shows the
changes in PSA amplitude ($A_{PSA}$) for different signs of h. The effective
sign of the Zeeman term is given by its relative phase with respect to the
reference signal, with negative $h$ being in-phase and positive out-of-phase.
The effective magnitude of $h$ is given by the signals’ amplitude, -15 dBm in
this case. The amplification imbalance allows us to change the phase
sensitivity of the circuit and globally bias the state of the spins. The
consequences of the bias are shown in Fig.2(b) and (c). In (b), we observe
that the modulation harmonics present in the spectra depend on the spin state,
with the biased signal containing a central carrier corresponding to
$\omega_{ref}$ that is absent for the unbiased solutions. $\omega_{bias}$,
whose frequency is the same as the reference signal’s, drives the oscillators
as they all acquire the same phase (i.e. align their spin direction). We
complement this picture in (c) using the time traces of the RF pulses colored
with their respective instantaneous phase and associated graph. It is clear
that despite antiferromagnetic coupling, we are able to induce a change in
ordering in both possible directions using $\omega_{bias}$ alone. Combining
time trace and spectrum analyses, we have different tools to study the
effectiveness of the biasing as well as develop differentiation and operation
protocols that allow us to program more complex problems. Although this is
very promising in terms of exploration of higher complexity schemes,
unexpected states also appear for intermediate amplitudes of $h$.

While it is expected that there will be a sudden breaking of up and down spin
symmetry at $h=2$, a gradual spin flipping occurs. Mixed spin states appear at
intermediate values of $\omega_{bias}$ amplitude ($\approx$-20 dBm), as seen
in Fig.3(a), resembling a chain with magnetic domains. Although the system is
indeed synchronized with $\omega_{bias}$, as shown in fig.3(b), these 3+1
states (three spins up and one down, or vice versa) are not a minimum energy
solution of eq.1 and suggest an unintended increase of the degrees of freedom
of the system.

In Fig.3(c) we show the phase diagram for a ring-shaped 4 spin system under an
external Zeeman field where one of the spin’s amplitude $|S_{1}|$ can be less
than one, i.e. one spin is shorter. The shortening of $S_{1}$ results in a
phase transition of the 4-spin system and appearance of the 3+1 states. In
fact, from the time traces shown in Fig.3, we can directly observe that the
spins have different electrical amplitudes and the odd one out has the
smallest one. Since electrical amplitude and spin amplitude are proportional
to each other, the emergence of the 3+1 states can give us information about
the operation of our circuit and its limitations.

Figure 3: Mixed artificial spin states for $\omega_{bias}$ with -20 dBm. (a)
Colored time traces of the RF spins. The appearance of 3+1 states is evidence
of additional degrees of freedom in the Hamiltonian of the system. (b) Phase
diagram of a ring-shaped 4 spin Ising machine with variable spin amplitude.
The 3+1 states are a consequence of phase-dependent amplification of the
spins. (c) PSD for both solutions with 3+1 states, we observe the same
synchronization as in fig.2(b), but the modulation peaks are characteristic to
these solutions.

Probing into the origin of this amplitude mismatch, we propose that the
emergence of the 3+1 states depends on the non-linearity of the saturation of
the LNA. Since we are injecting additional power to the system with the Zeeman
signal, the LNA saturates and thus gain compression occurs beyond its linear
range. Even if they are very close in the power transfer curve[12], its non-
linearity results in spins of different signs being amplified differently with
the minority spin shortening. An amplifier with a higher linear regime would
mitigate this state degeneration. If compression gain is unavoidable, a
stronger coupling between spins and a smoother saturation curve for the LNA
would suppress 3+1 states. Digital feedback using a field programmable gate
array (FPGA) has been employed successfully in previous time-multiplexed Ising
machines to improve amplitude stability [13] and can also be used instead to
the delay line to modify pairwise or all-to-all coupling. These findings can
be implemented in future circuit designs to improve the quality and complexity
of the solutions with larger amounts of spins.

Figure 4: Five spin states at two different $\omega_{bias}$ amplitudes. (a)
Time traces at -15 dBm and -20 dBm. We can observe that although the bias
manages to stabilize the phase, the solution is not phase binarized. (b) Power
spectral densities (PDS) of the solutions at the same amplitudes. We see that
each solution is synchronized to $\omega_{bias}$ has a characteristic
spectrum.

It is worthwhile to mention that spectral analysis can help us understand the
synchronization dynamics of the system as well as use for soution
differentiation. As we see in fig.2(b) and fig.3(b), the peaks at
$\omega_{ref}$ show that $\omega_{bias}$ drives the oscillators. Additionally,
both biased solutions have their own characteristic spectrum. Thus, including
this information in the digital feedback can allow us to design
differentiation metrics and stopping conditions for relaxation and annealing
protocols in bigger and more complex systems.

Finally, we tried to recreate a stable globally biased 3+2 solution with five
spins as it would appear in a spin ring. Evidently, a ring with an odd number
of spins will not have a stable unbiased solution as the phase difference on
each circulation period will never be zero. It was clear that for a large
enough bias, we do see all spins parallel to each other, as shown for -15 dBm
in fig.4(a). The phase transition mentioned before and shown in fig.3 could be
an indication that a 3+2 state could be viable and could allow us to construct
a magnetic system analog with two clearly defined domains. An amplitude of -23
dBm (fig.4(a)), is unable to produce such solution because the spins do not
synchronize with the external field (we can observe that the highest peak in
fig.4(b) is not at $\omega_{ref}$). Instead, the resulting state’s phase is
not binarized and produces spins do not comply with the definition of eq.1. We
believe that frustration is responsible for this phase slip and stable
solutions are achievable with digital coupling and individual bias, both
implementable with the aforementioned FPGA.

We have shown the broad features of a globally biased artificial spin state
space composed of RF pulses in a YIG waveguide. Despite antiferromagnetic
coupling between nearest neighbors in a 4 spin ring, we are able to induce
ferromagnetic ordering by injecting an external signal of the same frequency
to emulate the role of the Zeeman term in the Ising Hamiltonian. Intermediate
values of this Zeeman signal introduce degeneracy in the amplification of the
pulses and, consequently, 3+1 spin states. These effects be mitigated by
alternative amplification schemes whose implementation would guide future work
in enabling all-to-all spin coupling for tackling non-trivial optimization
tasks. The present work improves upon the emerging technology of commercially
feasible IM hardware accelerators. The SWIM concept has a high potential for
further scaling in terms of spin capacity, physical size and low-power low-
footprint circuits for applied combinatorial optimization.

## References

  * [1] N. Mohseni, P. L. McMahon, and T. Byrnes, “Ising machines as hardware solvers of combinatorial optimization problems,” _Nature Reviews Physics_ , vol. 4, no. 6, pp. 363–379, 2022.
  * [2] A. Litvinenko, R. Khymyn, V. H. González, A. A. Awad, V. Tyberkevych, A. Slavin, and J. Åkerman, “A spinwave ising machine,” _arXiv preprint arXiv:2209.04291_ , 2022.
  * [3] D. I. Albertsson, M. Zahedinejad, A. Houshang, R. Khymyn, J. Åkerman, and A. Rusu, “Ultrafast ising machines using spin torque nano-oscillators,” _Applied Physics Letters_ , vol. 118, no. 11, p. 112404, 2021.
  * [4] A. Houshang, M. Zahedinejad, S. Muralidhar, R. Khymyn, M. Rajabali, H. Fulara, A. A. Awad, J. Åkerman, J. Checiński, and M. Dvornik, “Phase-binarized spin hall nano-oscillator arrays: Towards spin hall ising machines,” _Physical Review Applied_ , vol. 17, no. 1, p. 014003, 2022.
  * [5] T. Honjo, T. Sonobe, K. Inaba, T. Inagaki, T. Ikuta, Y. Yamada, T. Kazama, K. Enbutsu, T. Umeki, R. Kasahara, K. ichi Kawarabayashi, and H. Takesue, “100,000-spin coherent ising machine,” _Science Advances_ , vol. 7, no. 40, p. eabh0952, 2021.
  * [6] K. Tatsumura, M. Yamasaki, and H. Goto, “Scaling out ising machines using a multi-chip architecture for simulated bifurcation,” _Nature Electronics_ , vol. 4, no. 3, pp. 208–217, 2021.
  * [7] A. Litvinenko, R. Khymyn, V. Tyberkevych, V. Tikhonov, A. Slavin, and S. Nikitov, “Tunable magnetoacoustic oscillator with low phase noise,” _Physical Review Applied_ , vol. 15, no. 3, p. 034057, 2021.
  * [8] A. Litvinenko, S. Grishin, Y. P. Sharaevskii, V. Tikhonov, and S. Nikitov, “A chaotic magnetoacoustic oscillator with delay and bistability,” _Technical Physics Letters_ , vol. 44, no. 3, pp. 263–266, 2018.
  * [9] A. Lucas, “Ising formulations of many np problems,” _Frontiers in physics_ , vol. 2, p. 5, 2014.
  * [10] F. Barahona, “On the computational complexity of ising spin glass models,” _Journal of Physics A: Mathematical and General_ , vol. 15, no. 10, p. 3241, oct 1982. [Online]. Available: https://dx.doi.org/10.1088/0305-4470/15/10/028
  * [11] Y. Inui, M. D. S. H. Gunathilaka, S. Kako, T. Aonishi, and Y. Yamamoto, “Control of amplitude homogeneity in coherent ising machines with artificial zeeman terms,” _Communications Physics_ , vol. 5, no. 1, p. 154, Jun 2022\. [Online]. Available: https://doi.org/10.1038/s42005-022-00927-x
  * [12] MiniCircuits. Coaxial low noise amplifier zx60-83ln12+. [Online]. Available: https://www.minicircuits.com/pdfs/ZX60-83LN12+.pdf
  * [13] H. Takesue, K. Inaba, T. Inagaki, T. Ikuta, Y. Yamada, T. Honjo, T. Kazama, K. Enbutsu, T. Umeki, and R. Kasahara, “Simulating ising spins in external magnetic fields with a network of degenerate optical parametric oscillators,” _Phys. Rev. Appl._ , vol. 13, p. 054059, May 2020. [Online]. Available: https://link.aps.org/doi/10.1103/PhysRevApplied.13.054059


# Primitive Pythagorean triples and neutrino mixing

Yuta Hyodo<EMAIL_ADDRESS>Graduate School of Science and
Technology, Tokai University, 4-1-1 Kitakaname, Hiratsuka, Kanagawa 259-1292,
Japan  Micro/Nano Technology Center, Tokai University, 4-1-1 Kitakaname,
Hiratsuka, Kanagawa 259-1292, Japan  Teruyuki Kitabayashi<EMAIL_ADDRESS>Department of Physics, Tokai University,
4-1-1 Kitakaname, Hiratsuka, Kanagawa 259-1292, Japan

###### Abstract

The primitive Pythagorean triples are the three natural numbers $(a,b,c)$ that
satisfy $c^{2}=a^{2}+b^{2}$ in a right triangle. We constructed a neutrino
mixing models related to primitive Pythagorean triples that satisfy the
observed values within the $3\sigma$ region for the reactor, solar, and
atmospheric neutrino mixing angles, as well as the Dirac CP phase.

##  I Introduction

The primitive Pythagorean triplesWilliams2017 ; MartensarXiv2112 ; Kocik2007 ;
SchmelzerarXiv2021 ; YekutieliarXiv2021 ; Cha2018 ; DeriyarXiv2023 ;
PricearXiv2008 ; AlperinarXiv2000  are the three natural numbers $(a,b,c)$
that satisfy $c^{2}=a^{2}+b^{2}$ in the right triangle (Fig. 1). The most
famous primitive Pythagorean triple is the combination $(3,4,5)$. The internal
angles of a triangle with the combination $(3,4,5)$ are
$(36.87^{\circ},53.13^{\circ},90.00^{\circ})$.

Figure 1: Right triangle.

The two internal angles $(36.87^{\circ},53.13^{\circ})$ closely approximate
the upper limit of the $3\sigma$ region for the so-called solar neutrino
mixing angle, $35.74^{\circ}$, and the upper limit of the $3\sigma$ region for
the so-called atmospheric neutrino mixing angle,
$51.0^{\circ}~{}\left(51.5^{\circ}\right)$, in the case of normal mass
ordering (NO) [inverted mass ordering (IO)]. In this paper, we constructed a
neutrino mixing model for the solar and atmospheric mixing angles with
primitive Pythagorean triples.

The rest of the paper is organized as follows. In Section. II, we review the
Pontecorvo-Maki-Nakagawa-Sakata (PMNS) mixing matrix and observed values to
compare the mixing model with observed values. In Section. III describes the
construction of a neutrino mixing model related to primitive Pythagorean
triples. In Section. IV presents the numerical calculation results of the
neutrino mixing model constructed in Section. III. In Section. V describes the
investigation of if the neutrino mass matrix from the proposed neutrino mixing
model is invariant under $Z_{2}$ symmetry. Section. VI presents the summary.

##  II PMNS Matrix

We review the PMNS mixing matrixPontecorvo1957 ; Pontecorvo1958 ; Maki1962PTP
; PDG  and observed values.

PMNS matrix is written as

$\displaystyle U$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{ccc}U_{e1}&U_{e2}&U_{e3}\\\ U_{\mu
1}&U_{\mu 2}&U_{\mu 3}\\\ U_{\tau 1}&U_{\tau 2}&U_{\tau 3}\end{array}\right)$
(4) $\displaystyle=$
$\displaystyle\left(\begin{array}[]{ccc}c_{12}c_{13}&s_{12}c_{13}&s_{13}e^{-i\delta}\\\
-s_{12}c_{23}-c_{12}s_{23}s_{13}e^{i\delta}&c_{12}c_{23}-s_{12}s_{23}s_{13}e^{i\delta}&s_{23}c_{13}\\\
s_{12}s_{23}-c_{12}c_{23}s_{13}e^{i\delta}&-c_{12}s_{23}-s_{12}c_{23}s_{13}e^{i\delta}&c_{23}c_{13}\\\
\end{array}\right),$ (8)

where, $c_{ij}=\cos\theta_{ij}$, $s_{ij}=\sin\theta_{ij}$ ($i,j$=1,2,3),
$\theta_{ij}$ is a mixing angle, and $\delta$ is the Dirac CP phase. The sines
and cosines of the three neutrino mixing angles of the PMNS matrix, $U$, are
given by

$\displaystyle s_{12}^{2}$ $\displaystyle=$
$\displaystyle\frac{|U_{e2}|^{2}}{1-|U_{e3}|^{2}},~{}~{}s_{23}^{2}=\frac{|U_{\mu
3}|^{2}}{1-|U_{e3}|^{2}},~{}~{}s_{13}^{2}=|U_{e3}|^{2},$ $\displaystyle
c_{12}^{2}$ $\displaystyle=$
$\displaystyle\frac{|U_{e1}|^{2}}{1-|U_{e3}|^{2}},~{}~{}c_{23}^{2}=\frac{|U_{\tau
3}|^{2}}{1-|U_{e3}|^{2}}.$ (9)

The Jarlskog rephasing invariantJarlskog ,

$\displaystyle J=\rm{Im}(U_{e1}U_{e2}^{\ast}U_{\mu 1}^{\ast}U_{\mu
2})=s_{12}s_{23}s_{13}c_{12}c_{23}c_{13}^{2}\sin{\delta}$ (10)

can be used to calculate $\delta$.

A global analysis of current data reveals the following the best-fit values of
the neutrino mixing angles in the case of NO, $m_{1}<m_{2}<m_{3}$, where
$m_{1},m_{2},$ and $m_{3}$ are the neutrino mass eigenvalues, NuFit :

$\displaystyle s_{12}^{2}$ $\displaystyle=$ $\displaystyle
0.303^{+0.012}_{-0.012}\quad(0.270\sim 0.341),$ $\displaystyle s_{23}^{2}$
$\displaystyle=$ $\displaystyle 0.451^{+0.019}_{-0.016}\quad(0.408\sim
0.603),$ $\displaystyle s_{13}^{2}$ $\displaystyle=$ $\displaystyle
0.02225^{+0.00056}_{-0.00059}\quad(0.02052\sim 0.02398),$
$\displaystyle\delta/^{\circ}$ $\displaystyle=$ $\displaystyle
232^{+36}_{-26}\quad(144\sim 350),$ (11)

where the $\pm$ denotes the $1\sigma$ region and the parentheses denote the
$3\sigma$ region. For the IO, $m_{3}<m_{1}<m_{2}$, we have

$\displaystyle s_{12}^{2}$ $\displaystyle=$ $\displaystyle
0.303^{+0.012}_{-0.012}\quad(0.270\sim 0.341),$ $\displaystyle s_{23}^{2}$
$\displaystyle=$ $\displaystyle 0.569^{+0.016}_{-0.021}\quad(0.412\sim
0.613),$ $\displaystyle s_{13}^{2}$ $\displaystyle=$ $\displaystyle
0.02223^{+0.00058}_{-0.00058}\quad(0.02048\sim 0.02416),$
$\displaystyle\delta/^{\circ}$ $\displaystyle=$ $\displaystyle
276^{+22}_{-29}\quad(194\sim 344).$ (12)

##  III Primitive primitive Pythagorean triples and neutrino mixing model

###  III.1 Primitive primitive Pythagorean triples

Primitive Pythagorean triples are the combinations $(a,b,c)$ of three natural
numbers that satisfy $a^{2}+b^{2}=c^{2}$. The relationship among $a,b,c$ is
described by

$\displaystyle a<b<c.$ (13)

$a,b,c$ can be expressed in terms of two natural numbers, $m$ and $n$, using
either

$\displaystyle\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$
(14)

or

$\displaystyle\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$
(15)

MartensarXiv2112 ; Kocik2007 ; SchmelzerarXiv2021 . In the case of Eq. (14)
(Eq. (15)), $a$ is an odd (even) number, and $b$ is an even (odd) number.

The relationship between $m$ and $n$ is as followsMartensarXiv2112 ; Kocik2007
; SchmelzerarXiv2021 :

  * •

$m>n$

  * •

$m$ and $n$ are natural numbers.

  * •

The parities of $m$ and $n$. (one is an odd number; and the other is an even
number)

  * •

$m$ and $n$ are coprime.

###  III.2 Relationship between the solar mixing angle, atmospheric mixing
angle, and primitive Pythagorean triples

The mixing matrix predicting a reactor mixing angle ($\theta_{13}$) of $0$ can
be written as

$\displaystyle U=R_{23}R_{12}.$ (16)

$R_{23}$ $(R_{12})$ is a rotation matrix for the 2-3 (1-2) planePDG ; it is
denoted by

$\displaystyle R_{23}=\left(\begin{array}[]{ccc}1&0&0\\\
0&\cos{\theta_{23}}&\sin{\theta_{23}}\\\
0&-\sin{\theta_{23}}&\cos{\theta_{23}}\end{array}\right),$ (20) $\displaystyle
R_{12}=\left(\begin{array}[]{ccc}\cos{\theta_{12}}&\sin{\theta_{12}}&0\\\
-\sin{\theta_{12}}&\cos{\theta_{12}}&0\\\ 0&0&1\end{array}\right).$ (24)

The relationship between primitive Pythagorean triples and neutrino mixing
angles, $\theta_{12}$ and $\theta_{23}$, is as follows:

$\displaystyle s_{12}=\frac{a}{c},~{}~{}s_{23}=\frac{b}{c}.$ (25)

With the relationship provided in Eq. (25), the mixing matrix predicting a
reactor mixing angle of $0$ can be written as

$\displaystyle
U_{\rm{PM}}=\left(\begin{array}[]{ccc}\frac{b}{c}&\frac{a}{c}&0\\\
-\frac{a^{2}}{c^{2}}&\frac{ab}{c^{2}}&\frac{b}{c}\\\
-\frac{ab}{c^{2}}&-\frac{b^{2}}{c^{2}}&\frac{a}{c}\end{array}\right).$ (29)

We call this mixing matrix primitive Pythagorean triple neutrino mixing (PM).

The mixing matrix predicts the vanishing reactor mixing angle. Next, we
performed modifications to improve the prediction of the reactor mixing angle.

###  III.3 Improved reproducibility of reactor mixing angle

To improve the reproducibility of the reactor mixing angle, we employed the
method of changing from tri-bimaximal (TBM) mixing to trimaximal (TM)
mixingXing2007PLB ; Harrison2006PRD . Using this method, the first column
(second column) of the $U_{\rm{PM}}$ remains unchanged, whereas the second and
third columns (first and third) are adjusted.

####  III.3.1 Modifications of the second and third columns

As shown in Ref. Xing2007PLB , the first column of the $U_{\rm{PM}}$ remains
unchanged, whereas the second and third columns are adjusted as follows:

$\displaystyle\left(U_{\rm{PM}}\right)_{\rm{C1}}$ $\displaystyle=$
$\displaystyle U_{\rm{PM}}U_{23}$ (33) $\displaystyle=$
$\displaystyle\left(\begin{array}[]{ccc}\frac{b}{c}&\frac{a}{c}\cos{\theta}&\frac{a}{c}\sin{\theta}e^{-i\phi}\\\
-\frac{a^{2}}{c^{2}}&\frac{b}{c^{2}}\left(a\cos{\theta}-c\sin{\theta}e^{i\phi}\right)&\frac{b}{c^{2}}\left(c\cos{\theta}+a\sin{\theta}e^{-i\phi}\right)\\\
\frac{ab}{c^{2}}&-\frac{1}{c^{2}}\left(b^{2}\cos{\theta}+ac\sin{\theta}e^{i\phi}\right)&\frac{1}{c^{2}}\left(ac\cos{\theta}-b^{2}\sin{\theta}e^{-i\phi}\right)\end{array}\right),$

where

$\displaystyle U_{23}=\left(\begin{array}[]{ccc}1&0&0\\\
0&\cos{\theta}&\sin{\theta}e^{-i\phi}\\\
0&-\sin{\theta}e^{i\phi}&\cos{\theta}\end{array}\right),$ (37)

$\theta$ denotes a rotation angle, and $\phi$ denotes a phase parameter.

From Eq. (9), the neutrino mixing angles of
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$ are as follows:

$\displaystyle s_{12}^{2}$ $\displaystyle=$
$\displaystyle\frac{a^{2}\cos^{2}{\theta}}{c^{2}-a^{2}\sin^{2}{\theta}},$
$\displaystyle s_{23}^{2}$ $\displaystyle=$
$\displaystyle\frac{b^{2}c^{2}\cos^{2}{\theta}+ab^{2}\left(a\sin^{2}{\theta}+c\sin{2\theta}\cos{\phi}\right)}{c^{4}-a^{2}c^{2}\sin^{2}{\theta}},$
$\displaystyle s_{13}^{2}$ $\displaystyle=$
$\displaystyle\frac{a^{2}}{c^{2}}\sin^{2}{\theta}.$ (38)

In term of $s_{13}^{2}$, detailing $s_{12}^{2}$ and $s_{23}^{2}$ leads to

$\displaystyle s_{12}^{2}$ $\displaystyle=$
$\displaystyle\frac{a^{2}}{c^{2}}\times\frac{\cos^{2}{\theta}}{1-s_{13}^{2}},$
$\displaystyle s_{23}^{2}$ $\displaystyle=$
$\displaystyle\frac{b^{2}}{c^{3}}\times\frac{a\sin{2\theta}\cos{\phi}+c\left(\cos^{2}{2\theta}+s_{13}^{2}\right)}{1-s_{13}^{2}}.$
(39)

From Eq. (10), the Dirac CP phase, $\delta$, is written as follows:

$\displaystyle\tan{\delta}=\frac{a(c^{2}+b^{2}+a^{2}\cos{2\theta})\sin{\phi}}{a\\{a^{2}+(c^{2}+b^{2})\cos{2\theta}\\}\cos{\phi}+c(a-b)(a+b)\sin{2\theta}},$
$\displaystyle(a=m^{2}-n^{2},b=2mn,c=m^{2}+n^{2}),$ (40)

and

$\displaystyle\tan{\delta}=\frac{a(2c^{2}-a^{2}+a^{2}\cos{2\theta})\sin{\phi}}{-a\\{a^{2}+(2c^{2}-a^{2})\cos{2\theta}\\}\cos{\phi}+c(b-a)(b+a)\sin{2\theta}},$
$\displaystyle(a=2mn,b=m^{2}-n^{2},c=m^{2}+n^{2}).$ (41)

Figure 2: Dependence of $m$ on $n$ of $\left(U_{\rm{PM}}\right)_{\rm{C1}}$ in
the case of $\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$.
$\theta$, and $\phi$ varied within the ranges of $0$ to $\frac{\pi}{2}$ and
$0$ to $2\pi$, respectively.

Figure 3: Dependence of $s_{12}^{2}$ and $s_{23}^{2}$ on $s_{13}^{2}$ of
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$ in the case of
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$. $\theta$,
$\phi$ and $m$ and $n$ varied within the ranges of $0\sim\frac{\pi}{2}$, the
range of $\phi$ is varied within $0\sim 2\pi$ and the range of $m$ and $n$ is
varied within $0$ to $\frac{\pi}{2}$ and $0$ to $2\pi$, and $0\sim 1000$,
respectively.

Figure 4: Dependence of $s_{12}^{2}$ and $s_{13}^{2}$ on $\theta$ in the cases
of $m=11,n=6$, $m=502,n=265$ and $m=1000,n=551$.

Figure 5: Dependence of $s_{23}^{2}$ on $\theta$ in the cases of $m=11,n=6$,
$m=502,n=265$ and $m=1000,n=551$. We selected
$0$,$\frac{\pi}{2}$,$\pi$,$\frac{5\pi}{4}$ for $\phi$.

Figure 6: Dependence of $s_{23}^{2}$ on $\phi$ in the cases of $m=11,n=6$,
$m=502,n=265$ and $m=1000,n=551$. We selected $\theta$ such that $s_{13}^{2}$
is approximately $0.022$.

Figure 7: Dependence of Dirac CP phase $\delta$ on $\theta$ in the cases of
$m=11,n=6$, $m=502,n=265$ and $m=1000,n=551$. We selected
$0$,$\frac{\pi}{2}$,$\pi$,$\frac{5\pi}{4}$ for $\phi$.

Figure 8: Dependence of Dirac CP phase $\delta$ on $\phi$ in the cases of
$m=11,n=6$, $m=502,n=265$ and $m=1000,n=551$. We selected $\theta$ such that
$s_{13}^{2}$ is approximately $0.022$.

####  III.3.2 Modifications first and third columns

As shown in Ref. Harrison2006PRD , the second column of the $U_{\rm{PM}}$
remains unchanged, whereas the first and third columns are adjusted as
follows:

$\displaystyle\left(U_{\rm{PM}}\right)_{\rm{C2}}$ $\displaystyle=$
$\displaystyle U_{\rm{PM}}U_{13}$ (45) $\displaystyle=$
$\displaystyle\left(\begin{array}[]{ccc}\frac{b}{c}\cos{\theta}&\frac{a}{c}&\frac{b}{c}\sin{\theta}e^{-i\phi}\\\
-\frac{1}{c^{2}}\left(a^{2}\cos{\theta}+bc\sin{\theta}e^{i\phi}\right)&\frac{ab}{c^{2}}&\frac{1}{c^{2}}\left(bc\cos{\theta}-a^{2}\sin{\theta}e^{-i\phi}\right)\\\
\frac{a}{c^{2}}\left(b\cos{\theta}-c\sin{\theta}e^{i\phi}\right)&-\frac{b^{2}}{c^{2}}&\frac{a}{c^{2}}\left(c\cos{\theta}+b\sin{\theta}e^{-i\phi}\right)\end{array}\right),$

where

$\displaystyle
U_{13}=\left(\begin{array}[]{ccc}\cos{\theta}&0&\sin{\theta}e^{-i\phi}\\\
0&1&0\\\ -\sin{\theta}e^{i\phi}&0&\cos{\theta}\end{array}\right),$ (49)

where $\theta$ denotes a rotation angle, and $\phi$ denotes a phase parameter.

From Eq. (9), the neutrino mixing angles of
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$ are as follows:

$\displaystyle s_{12}^{2}$ $\displaystyle=$
$\displaystyle\frac{a^{2}}{c^{2}-b^{2}\sin^{2}{\theta}},$ $\displaystyle
s_{23}^{2}$ $\displaystyle=$
$\displaystyle\frac{b^{2}c^{2}\cos^{2}{\theta}-a^{2}bc\sin{2\theta}\cos{\phi}+a^{4}\sin^{2}{\theta}}{c^{4}-b^{2}c^{2}\sin^{2}{\theta}},$
$\displaystyle s_{13}^{2}$ $\displaystyle=$
$\displaystyle\frac{b^{2}}{c^{2}}\sin^{2}{\theta}.$ (50)

In term of $s_{13}^{2}$, detailing $s_{12}^{2}$ and $s_{23}^{2}$ leads to

$\displaystyle s_{12}^{2}$ $\displaystyle=$
$\displaystyle\frac{a^{2}}{c^{2}}\times\frac{1}{1-s_{13}^{2}},$ $\displaystyle
s_{23}^{2}$ $\displaystyle=$
$\displaystyle\frac{b^{2}c^{2}\cos^{2}{\theta}-a^{2}bc\sin{2\theta}\cos{\phi}+a^{4}\sin^{2}{\theta}}{c^{4}(1-s_{13}^{2})}.$
(51)

From Eq. (10), $\delta$ is written as follows:

$\displaystyle\tan{\delta}=\frac{b(2c^{2}-b^{2}+b^{2}\cos{2\theta})\sin{\phi}}{-b\\{b^{2}+(2c^{2}-b^{2})\cos{2\theta}\\}\cos{\phi}+c(a-b)(a+b)\sin{2\theta}},$
$\displaystyle(a=m^{2}-n^{2},b=2mn,c=m^{2}+n^{2}),$ (52)

and

$\displaystyle\tan{\delta}=\frac{b(c^{2}+a^{2}+b^{2}\cos{2\theta})\sin{\phi}}{b\\{b^{2}+(c^{2}+a^{2})\cos{2\theta}\\}\cos{\phi}+c(b-a)(b+a)\sin{2\theta}},$
$\displaystyle(a=2mn,b=m^{2}-n^{2},c=m^{2}+n^{2}).$ (53)

Eqs. (40) and (53), as well as Eqs. (41) and (52), are identical. We will
algebraically demonstrate this topic in a future study.

##  IV Numerical calculations

First, we performed numerical calculations for $m$ and $n$, $\theta$ and
$\phi$ within the ranges of $1\sim 1000$, $0$ to $\frac{\pi}{2}$, and $0$ to
$2\pi$, searching for combinations that satisfy the observed values.
Subsequently, we selected three combinations of $m$ and $n$ and conducted
numerical calculations to investigate the dependence of neutrino parameters on
$\theta$ and $\phi$. We investigated the relationship between $m$ and $n$
using the Pearson product-moment correlation coefficient (PPMCC)Pearson1895
as follows:

$\displaystyle
r_{mn}=\frac{{\displaystyle\sum_{k=1}^{l}(m_{k}-\overline{m})(n_{k}-\overline{n})}}{\sqrt{{\displaystyle\sum_{k=1}^{l}(m_{k}-\overline{m})^{2}}}\sqrt{{\displaystyle\sum_{k=1}^{l}(n_{k}-\overline{n})^{2}}}}.$
(54)

PPMCC measures the linear correlation between two sets of data within the
range of

$\displaystyle-1\leq r_{mn}\leq 1.$ (55)

Positive and negative linear correlations are observed as the value approaches
$1$ and $-1$, respectively. If the value is close to $0$, no correlation is
observed.
$(m_{1},n_{1}),(m_{2},n_{2}),\cdots,(m_{k},n_{k}),\cdots,(m_{l},n_{l})$
represents $l$ sets of data comprising $m$ and $n$, and it becomes
$(m_{1},n_{1}),(m_{2},n_{2}),\cdots,(m_{k},n_{k}),\cdots,(m_{9716},n_{9716})=(11,6),(15,8),\cdots,(1000,551)$
(Fig. 2). Furthermore, $\overline{m}$ and $\overline{n}$ represent the
average, as described in

$\displaystyle\overline{m}=\frac{{\sum_{k=1}^{l}m_{k}}}{l},$ (56)

and

$\displaystyle\overline{n}=\frac{{\sum_{k=1}^{l}n_{k}}}{l}.$ (57)

In the case of Fig. 2, $l$ is $9716$, and the averages of $m$ and $n$ are
$\overline{m}=665$ and $\overline{n}=350$.

If PPMC $r_{mn}$ is close to 1 (-1), we employed the linear regression model,

$\displaystyle n=wm+\beta,$ (58)

where $w$ and $\beta$ denote the slope and intercept, respectively.

Figure 9: Similar to Fig. 2, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$.

Figure 10: Similar to Fig. 3, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$.

Figure 11: Similar to Fig. 4, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. We selected
$(m,n)=(10,3),(500,163),(1000,331)$.

Figure 12: Similar to Fig. 5, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. We selected
$(m,n)=(10,3),(500,163),(1000,331)$.

Figure 13: Similar to Fig. 6, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. We selected
$(m,n)=(10,3),(500,163),(1000,331)$.

Figure 14: Similar to Fig. 7, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. We selected
$(m,n)=(10,3),(500,163),(1000,331)$.

Figure 15: Similar to Fig. 8, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. We selected
$(m,n)=(10,3),(500,163),(1000,331)$.

###  IV.1 $\left(U_{\rm{PM}}\right)_{\rm{C1}}$

####  IV.1.1 The case of the
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$

Fig. 2 illustrates the dependence of the $m$ on $n$ for
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$ with the use of $\left(a,b,c\right)=$

$\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$. We varied $\theta$ and $\phi$
within the ranges of $0$ to $\frac{\pi}{2}$, and $0$ to $2\pi$, respectively.

Fig. 2 shows $m$ and $n$ satisfying the $3\sigma$ region for mixing angles,
$s_{ij}^{2}$, and the Dirac CP phase, $\delta$. The left and right panels
shows the NO and IO cases, respectively. The key observations from Fig. 2 are
as follows:

  * •

In the NO and IO cases, the PPMCs were both $0.997$. Considering that both
PPMCs were close to $1$, a very strong positive linear correlation was
observed between $m$ and $n$.

  * •

The relationship between $m$ and $n$ exhibits strong linearity. Thus, using
Eq. (58), the slope, $w$, was $0.527$ and $0.526$ in the NO and IO cases,
respectively. If $m$ increased by $10$, $n$ increased by $5$. The intercept,
$\beta$, was $-0.0869$ and $-0.104$ in the NO and IO cases, respectively.
Slope $w$ and intercept $\beta$ did not significantly vary between the NO and
IO cases.

Fig. 3 illustrates the relationship between $s_{12}^{2}$ and $s_{13}^{2}$, in
addition to the relationship between $s_{23}^{2}$ and $s_{13}^{2}$, for
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$ using
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$, satisfying the
$3\sigma$ regions of mixing angles $s_{ij}^{2}$ and the Dirac CP phase,
$\delta$. We varied $\theta$, $\phi$, and $m$ and $n$ within the ranges of $0$
to $\frac{\pi}{2}$, $0$ to $2\pi$, and $0\sim 1000$, respectively.

The upper (lower) panel shows the relationship between $s_{12}^{2}$ and
$s_{13}^{2}$ ($s_{23}^{2}$ and $s_{13}^{2}$). The upper and lower left panels
show the NO case, and the upper and lower right panels show the IO case. The
red line in the figure represents the best-fit values of $s_{12}^{2}$,
$s_{13}^{2}$, and $s_{23}^{2}$, and the black line shows the $3\sigma$ region
for $s_{12}^{2}$, $s_{13}^{2}$, and $s_{23}^{2}$.

  * •

In the NO and IO cases, the predicted values of $s_{12}^{2}$, $s_{13}^{2}$,
and $s_{23}^{2}$ from the model fell within the $3\sigma$ region.

  * •

If we appropriately select $m$, $n$, and $\theta$, the predicted values of
$s_{12}^{2}$,$s_{13}^{2}$, and $s_{23}^{2}$ from the model can simultaneously
satisfy the best-fit values.

A benchmark point,

$\displaystyle(m,n,\theta,\phi)=(11,6,0.2892~{}\rm{rad},3.364~{}\rm{rad})$
(59)

yields

$\displaystyle(s_{12}^{2},s_{13}^{2},s_{23}^{2},\delta)=(0.2759,0.02384,0.4735,191.6^{\circ}).$
(60)

Next, using $(m,n)=(11,6),(502,265),(1000,551)$, we investigated the
relationship between neutrino mixing angles, $\delta$, rotation angle
$\theta$, and phase parameter $\phi$.

Fig. 4 illustrates the dependence of the $s_{12}^{2}$ and $s_{13}^{2}$ on
$\theta$ in the cases of $m=11,n=6$, $m=502,n=265$ and $m=1000,n=551$. The red
line represents the best-fit values of $s_{12}^{2}$ and $s_{13}^{2}$, and the
black line in the upper panel indicates the $3\sigma$ region for $s_{12}^{2}$.
The enlarged view in the lower panel illustrates the dependence of
$s_{13}^{2}$ on $\theta$ when expanded into the 3$\sigma$ region of
$s_{13}^{2}$. The key observations from Fig. 4 are as follows:

  * •

$\theta$ exists to satisfy the 3$\sigma$ region of $s_{13}^{2}$ and
$s_{12}^{2}$.

  * •

The range of $\theta$ satisfying the 3$\sigma$ region of $s_{13}^{2}$ was
limited to $0.255~{}\rm{rad}\leq\theta\leq 0.296~{}\rm{rad}$.

  * •

If the range of $\theta$ was $0.255~{}\rm{rad}\leq\theta\leq
0.296~{}\rm{rad}$, $s_{12}^{2}$ was within the 3$\sigma$ region of
$s_{12}^{2}$.

Fig. 5 illustrates the dependence of $s_{23}^{2}$ on $\theta$ in the cases of
$m=11,n=6$, $m=502,n=265$ and $m=1000,n=551$. We selected
$0$,$\frac{\pi}{2}$,$\pi$,$\frac{5\pi}{4}$ for $\phi$. The red line represents
the best-fit values of $s_{23}^{2}$, and the black line shows the $3\sigma$
region for $s_{23}^{2}$. The left and right panels show the NO and IO cases,
respectively.

Fig. 6 illustrates the dependence of $s_{23}^{2}$ on $\phi$ in the cases of
$m=11,n=6$, $m=502,n=265$ and $m=1000,n=551$. We selected $\theta$ such that
$s_{13}^{2}$ was approximately $0.022$. The left and right panels show the NO
and IO cases, respectively. The key observations from Figs.5 and 6 are as
follows:

  * •

$\theta$ and $\phi$ exist to satisfy the 3$\sigma$ region of $s_{23}^{2}$.

  * •

This combination of $m$ and $n$ does not satisfy the best-fit value for NO;
however, it satisfies the best-fit value for IO.

  * •

The other $m$ and $n$ combinations satisfy the best-fit value of NO. For
example, a benchmark point,
$(m,n,\theta,\phi)=(44,23,0.2743~{}\rm{rad},3.188~{}\rm{rad})$, yields
$(s_{12}^{2},s_{13}^{2},s_{23}^{2},\delta)=$

$(0.3093,0.02390,0.4512,182.5^{\circ})$.

  * •

As shown in Fig. 6, the range of $\phi$ satisfying the 3$\sigma$ region of
$s_{23}^{2}$ was limited to $1.8~{}\rm{rad}\leq\phi\leq 4.5~{}\rm{rad}$.

Fig. 7 illustrates the dependence of $\delta$ on $\theta$ in the cases of
$m=11,n=6$, $m=502,n=265$ and $m=1000,n=551$. We selected
$0$,$\frac{\pi}{2}$,$\pi$,$\frac{5\pi}{4}$ for $\phi$. The red line represents
the best-fit values of $\delta$, and the black line shows the $3\sigma$ region
for $\delta$. The left and right panels show the NO and IO cases,
respectively.

Fig. 8 illustrates the dependence of $\delta$ on $\phi$ in the cases of
$m=11,n=6$, $m=502,n=265$ and $m=1000,n=551$. We selected $\theta$ such that
$s_{13}^{2}$ was approximately $0.022$. The left and right panels show the NO
and IO cases, respectively. The key observations from Figs.7 and 8 are as
follows:

  * •

No difference was observed in the $\delta$ for the different combinations of
$m$ and $n$.

  * •

As shown in Fig. 8, the range of $\phi$ satisfying the 3$\sigma$ region of
$\delta$ was limited to $2.4~{}\rm{rad}\leq\phi\leq 6.2~{}\rm{rad}$.

####  IV.1.2 The case of the
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$

Fig. 9 is the same as Fig. 2, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. The key observations from Fig. 9 are as
follows:

  * •

In the NO and IO cases, the PPMCs were both $0.994$. The PPMCs were close to
$1$; thus, a very strong positive linear correlation was observed between $m$
and $n$.

  * •

The relationship between $m$ and $n$ exhibits strong linearity. Thus, using
Eq. (58), slope $w$ was $0.311$ in the NO and IO cases, respectively. If $m$
increased by $10$, $n$ increased by $3$. The intercept $\beta$ was $-0.150$
and $-0.196$ in the NO and IO cases, respectively. Slope $w$ and intercept
$\beta$ did not significantly vary between the NO and IO cases.

  * •

The slope of Fig. 9 is smaller than that of Fig. 2.

Fig. 10 is the same as Fig. 3, except, for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. The key observations from Fig. 10 are as
follows:

  * •

In the NO and IO cases, the predicted values of $s_{12}^{2}$, $s_{13}^{2}$,
and $s_{23}^{2}$ from the model fell within the $3\sigma$ region.

  * •

If we appropriately selected $m$, $n$ and $\theta$, the predicted values of
$s_{12}^{2}$, $s_{13}^{2}$, and $s_{23}^{2}$ from the model can simultaneously
satisfy the best-fit values.

A benchmark point,

$\displaystyle(m,n,\theta,\phi)=(10,3,0.2785~{}\rm{rad},3.952~{}\rm{rad}),$
(61)

yields

$\displaystyle(s_{12}^{2},s_{13}^{2},s_{23}^{2},\delta)=(0.2867,0.02290,0.5327,318.1^{\circ}).$
(62)

Next, using $(m,n)=(10,3),(500,163),(1000,331)$, we investigated the
relationship between the neutrino mixing angles, $\delta$, rotation angle
$\theta$ and phase parameter $\phi$.

Fig. 11 is the same as Fig. 4, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. The key observations from Fig. 11 are as
follows:

  * •

$\theta$ exists to satisfy the 3$\sigma$ region of $s_{13}^{2}$ and
$s_{12}^{2}$.

  * •

The range of $\theta$ satisfying the 3$\sigma$ region of $s_{13}^{2}$ was
limited to $0.241~{}\rm{rad}\leq\theta\leq 0.285~{}\rm{rad}$.

  * •

If the range of $\theta$ is $0.241~{}\rm{rad}\leq\theta\leq 0.285~{}\rm{rad}$,
$s_{12}^{2}$ is within the 3$\sigma$ region of $s_{12}^{2}$.

Fig. 12 is the same as Fig. 5, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. Fig. 13 is the same as Fig. 6, except
for $\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. The key observations from Figs. 12 and
13 are as follows:

  * •

$\theta$ and $\phi$ exist to satisfy the 3$\sigma$ region of $s_{23}^{2}$.

  * •

There is a combination of $m$ and $n$ that satisfies the best-fit values of NO
and IO.

  * •

As shown in Fig. 13, the range of $\phi$ satisfying the 3$\sigma$ region of
$s_{23}^{2}$ was limited to $1.7~{}\rm{rad}\leq\phi\leq 4.6~{}\rm{rad}$.

Fig. 14 is the same as Fig. 7, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. Fig. 15 is the same as Fig. 8, except
for $\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. The key observations from Fig. 14 and 15
are as follows:

  * •

No difference was observed in the $\delta$ for different combinations of $m$
and $n$.

  * •

As shown in Fig. 15, the range of $\phi$ satisfying the 3$\sigma$ region of
$\delta$ was limited to $0~{}\rm{rad}\leq\phi\leq 0.5~{}\rm{rad}$ and
$3.3~{}\rm{rad}\leq\phi\leq 2\pi~{}\rm{rad}$.

###  IV.2 $\left(U_{\rm{PM}}\right)_{\rm{C2}}$

####  IV.2.1 The case of the
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$

Figure 16: Similar to Fig. 2, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$.

Figure 17: Similar to Fig. 3, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$.

Figure 18: Similar to Fig. 4, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. We selected
$(m,n)=(19,10),(500,261),(1000,553)$.

Figure 19: Similar to Fig. 5, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. We selected
$(m,n)=(19,10),(500,261),(1000,553)$ and
$\phi=0,\frac{\pi}{2},\pi,\frac{7\pi}{4}$.

Figure 20: Similar to Fig. 6, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$.We selected
$(m,n)=(19,10),(500,261),(1000,553)$.

Figure 21: Similar to Fig. 7, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. We selected
$(m,n)=(19,10),(500,261),(1000,553)$ and
$\phi=0,\frac{\pi}{2},\pi,\frac{7\pi}{4}$.

Figure 22: Similar to Fig. 8, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. We selected
$(m,n)=(19,10),(500,261),(1000,553)$.

Fig. 16 is the same as Fig. 2, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. The key observations from Fig. 16 are as
follows:

  * •

In the NO and IO cases, the PPMCs were both $0.999$. Both PPMCs were close to
$1$; thus, a very strong positive linear correlation was observed between $m$
and $n$.

  * •

The relationship between $m$ and $n$ exhibits strong linearity. Thus, using
Eq. (58), $w$ was $0.527$ in both cases. If $m$ increased by $10$, $n$
increased by $5$. $\beta$ was $-0.106$ and $-0.0698$ in the NO and IO cases.
Slope $w$ and intercept $\beta$ did not significantly vary between the NO and
IO cases.

Fig. 17 is the same as Fig. 3, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. The key observations from Fig. 17 are as
follows:

  * •

In the NO and IO cases, the predicted values of $s_{12}^{2}$, $s_{13}^{2}$,
and $s_{23}^{2}$ from the model fell within the $3\sigma$ region.

  * •

The predicted values from the model exceeded or were equal to $0.3$ and $0.55$
for $s_{12}^{2}$ and $s_{23}^{2}$, respectively, in both cases. Thus, it was
impossible to simultaneously satisfy the best-fit values for $s_{12}^{2}$,
$s_{13}^{2}$ and $s_{23}^{2}$ in the case of NO.

A benchmark point,

$\displaystyle(m,n,\theta,\phi)=(19,10,0.1751~{}\rm{rad},0.4214~{}\rm{rad})$
(63)

yields

$\displaystyle(s_{12}^{2},s_{13}^{2},s_{23}^{2},\delta)=(0.3273,0.02062,0.5914,157.1^{\circ}).$
(64)

Next, using $(m,n)=(19,10),(500,261),(1000,553)$, we investigated the
relationship between the neutrino mixing angles, $\delta$, rotation angle
$\theta$, and phase parameter $\phi$.

Fig. 18 is the same as Fig. 4, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. The key observations from Fig. 18 are as
follows:

  * •

$\theta$ exists to satisfy the 3$\sigma$ region of $s_{13}^{2}$ and
$s_{12}^{2}$.

  * •

The range of $\theta$ satisfying the 3$\sigma$ region of $s_{13}^{2}$ was
limited to $0.173~{}\rm{rad}\leq\theta\leq 0.187~{}\rm{rad}$.

  * •

If the range of $\theta$ was $0.173~{}\rm{rad}\leq\theta\leq
0.187~{}\rm{rad}$, $s_{12}^{2}$ was within the 3$\sigma$ region of
$s_{12}^{2}$. However, it exceeded the best-fit value.

Fig. 19 is the same as Fig. 5, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. Fig. 20 is the same as Fig. 6, except
for $\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. The key observations from Figs.19 and 20
are as follows:

  * •

$\theta$ and $\phi$ exist to satisfy the 3$\sigma$ region of $s_{23}^{2}$.

  * •

This combination of $m$ and $n$ satisfies the best-fit value for IO.

  * •

This combination of $m$ and $n$ does not satisfy the best-fit value for NO.

  * •

If $s_{13}^{2}$ aligns with the best-fit value, $s_{23}^{2}$ was close to the
upper limit of the 3$\sigma$ region.

  * •

As shown in Fig. 20, the range of $\phi$ satisfying the 3$\sigma$ region of
$s_{23}^{2}$ was limited to $0~{}\rm{rad}\leq\phi\leq 1~{}\rm{rad}$ and
$5.3~{}\rm{rad}\leq\phi\leq 2\pi~{}\rm{rad}$.

Fig. 21 is the same as Fig. 7, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. Fig. 22 is the same as Fig. 8, except
for $\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. The key observations from Figs.21 and 22
are as follows:

  * •

No difference was observed in the $\delta$ for different combinations of $m$
and $n$.

  * •

As shown in Fig. 22, the range of $\phi$ satisfying the 3$\sigma$ region of
the $\delta$ was limited to $0~{}\rm{rad}\leq\phi\leq 0.6~{}\rm{rad}$ and
$3.3~{}\rm{rad}\leq\phi\leq 2\pi~{}\rm{rad}$.

####  IV.2.2 The case of the
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$

Figure 23: Similar to Fig. 2, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$.

Figure 24: Similar to Fig. 3, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$.

Figure 25: Similar to Fig. 4, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$.We selected
$(m,n)=(16,5),(500,153),(1000,317)$.

Figure 26: Similar to Fig. 5, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. We selected
$(m,n)=(16,5),(500,153),(1000,317)$ and
$\phi=0,\frac{\pi}{2},\pi,\frac{7\pi}{4}$.

Figure 27: Similar to Fig. 6, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. We selected
$(m,n)=(16,5),(500,153),(1000,317)$.

Figure 28: Similar to Fig. 7, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. We selected
$(m,n)=(16,5),(500,153),(1000,317)$ and
$\phi=0,\frac{\pi}{2},\pi,\frac{7\pi}{4}$.

Figure 29: Similar to Fig. 8, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. We selected
$(m,n)=(16,5),(500,153),(1000,317)$.

Fig. 23 is the same as Fig. 2, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. The key observations from Fig. 23 are as
follows:

  * •

In the NO and IO cases, the PPMCs were both $0.999$. The PPMCs were close to
$1$; thus, a very strong positive linear correlation was observed between $m$
and $n$.

  * •

The relationship between $m$ and $n$ exhibits strong linearity. Thus, using
Eq. (58), $w$ was $0.311$ and $0.310$ in the NO and IO cases, respectively. If
$m$ increased by $10$, $n$ increased by $3$. $\beta$ was $-0.0995$ and
$-0.0370$ in the NO and IO cases, respectively. $w$ and $\beta$ did not
significantly vary between the NO and IO cases.

  * •

The slop of Fig. 23 is smaller than that of Fig. 16.

Fig. 24 is the same as Fig. 3, except for
$\left(a,b,c\right)=\left(m^{2}-n^{2},2mn,m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. The key observations from Fig. 24 are as
follows:

  * •

In the NO and IO cases, the predicted values of $s_{12}^{2}$, $s_{13}^{2}$,
and $s_{23}^{2}$ from the model fell within the $3\sigma$ region.

  * •

The predicted values from the model exceeded or were equal to $0.3$ for
$s_{12}^{2}$ and $0.55$ for $s_{23}^{2}$ in both case. Thus, it was impossible
to simultaneously satisfy the best-fit values for $s_{12}^{2}$, $s_{13}^{2}$,
and $s_{23}^{2}$ in the case of NO.

A benchmark point,

$\displaystyle(m,n,\theta,\phi)=(16,5,0.1881~{}\rm{rad},5.529~{}\rm{rad})$
(65)

yields

$\displaystyle(s_{12}^{2},s_{13}^{2},s_{23}^{2},\delta)=(0.3321,0.02363,0.5986,319.2^{\circ}).$
(66)

Next, using $(m,n)=(16,5),(500,153),(1000,317)$, we investigated the
relationship between the neutrino mixing angles, $\delta$, rotation angle
$\theta$, and phase parameter $\phi$.

Fig. 25 is the same as Fig. 11, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. The key observations from Fig. 25 are as
follows:

  * •

$\theta$ exists to satisfy the 3$\sigma$ region of $s_{13}^{2}$ and
$s_{12}^{2}$.

  * •

The range of $\theta$ satisfying the 3$\sigma$ region of $s_{13}^{2}$ was
limited to $0.170~{}\rm{rad}\leq\theta\leq 0.190~{}\rm{rad}$.

  * •

If the range of $\theta$ was $0.170~{}\rm{rad}\leq\theta\leq
0.190~{}\rm{rad}$, $s_{12}^{2}$ was within the 3$\sigma$ region of
$s_{12}^{2}$. However, it exceeded the best-fit value.

Fig. 26 is the same as Fig. 5, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. Fig. 27 is the same as Fig. 6, except
for $\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$. The key observations from Fig. 26 and
Fig. 27 are as follows:

  * •

$\theta$ and $\phi$ exist to satisfy the 3$\sigma$ region of $s_{23}^{2}$.

  * •

This combination of $m$ and $n$ satisfies the best-fit value for IO.

  * •

This combination of $m$ and $n$ does not satisfy the best-fit value for NO.

  * •

If $s_{13}^{2}$ aligns with the best-fit value, $s_{23}^{2}$ was close to the
upper limit of the 3$\sigma$ region.

  * •

As shown in Fig. 27, the range of $\phi$ satisfying the 3$\sigma$ region of
$s_{23}^{2}$ was limited to $0~{}\rm{rad}\leq\phi\leq 1.0~{}\rm{rad}$ and
$5.3~{}\rm{rad}\leq\phi\leq 2\pi~{}\rm{rad}$.

Fig. 28 is the same as Fig. 7, except for
$\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. Fig. 29 is the same as Fig. 8, except
for $\left(a,b,c\right)=\left(2mn,m^{2}-n^{2},m^{2}+n^{2}\right)$ and
$\left(U_{\rm{PM}}\right)_{\rm{C2}}$. The key observations from Figs.28 and 29
are as follows:

  * •

No difference was observed in the $\delta$ for different combinations of $m$
and $n$.

  * •

As shown in Fig. 29, the range of $\phi$ satisfying the 3$\sigma$ region of
$\delta$ was limited to $2.6~{}\rm{rad}\leq\phi\leq 6.0~{}\rm{rad}$.

##  V $Z_{2}$ Symmetry

The topic of symmetry is crucial in neutrino physics. As an example, a texture
of the neutrino mass matrix based on the $A_{n}$ or $S_{n}$ symmetry has been
proposed Ref. Altarelli2010RMP . Discussions regarding the $Z_{2}$ symmetry
are conducted in the context of TBM or TM mixing models LamPRD2006 ;
LamPLB2007 ; LamPRL2008 . Ref. DicusPRD2011  discusses the $Z_{2}$ symmetry of
the mixing matrix proposed in Ref. DicusPRD2011 . Our PM exhibits $Z_{2}$
symmetry. Herein, we discuss the $Z_{2}$ symmetry of PM.

Assuming that the mass matrix of charged leptons is diagonal and real, the
neutrino mass matrix, $M_{\nu}$, is as follows:

$\displaystyle
M_{\nu}=U^{\ast}\textrm{diag}(m_{1},m_{2},m_{3})U^{{\dagger}}=\left(\begin{array}[]{ccc}A&B&C\\\
B&D&E\\\ C&E&F\end{array}\right),$ (70)

where $m_{1}$, $m_{2}$, and $m_{3}$ are the neutrino mass eigenvalues.

If $M_{\nu}$ satisfies the transformation,

$\displaystyle G^{T}M_{\nu}G$ $\displaystyle=$ $\displaystyle M_{\nu},$
$\displaystyle G^{2}$ $\displaystyle=$ $\displaystyle\textrm{diag}(1,1,1).$
(71)

$M_{\nu}$ is invariant under $Z_{2}$ symmetryLamPRD2006 . In this context, $G$
is expressed as follows:

$\displaystyle G$ $\displaystyle=$ $\displaystyle
g_{1}v_{1}v_{1}^{{\dagger}}+g_{2}v_{2}v_{2}^{{\dagger}}+g_{3}v_{3}v_{3}^{{\dagger}},$
(72)

where $v_{1}$, $v_{2}$, and $v_{3}$ are column vectors of the mixing matrix,
and one of the eigenvalues among $g_{1}$, $g_{2}$, and $g_{3}$ was equal to
$-1~{}(1)$. The remaining two eigenvalues were equal to $1~{}(-1)$LamPRD2006 ;
LamPLB2007 ; LamPRL2008 .

The column vectors of the PM are

$\displaystyle v_{1}$ $\displaystyle=$
$\displaystyle\left(\frac{b}{c},-\frac{a^{2}}{c^{2}},\frac{ab}{c^{2}}\right)^{T},$
$\displaystyle v_{2}$ $\displaystyle=$
$\displaystyle\left(\frac{a}{c},-\frac{ab}{c^{2}},-\frac{b^{2}}{c^{2}}\right)^{T},$
$\displaystyle v_{3}$ $\displaystyle=$
$\displaystyle\left(0,\frac{b}{c},\frac{a}{c}\right)^{T}.$ (73)

Based on Eq. (72), $G_{i}~{}(i=1,2,3)$ can be expressed as

$\displaystyle
G_{1}=-v_{1}\left(v_{1}\right)^{{\dagger}}+v_{2}\left(v_{2}\right)^{{\dagger}}+v_{3}\left(v_{3}\right)^{{\dagger}}$
$\displaystyle=\left(\begin{array}[]{ccc}\frac{(a-b)(a+b)}{c^{2}}&\frac{2a^{2}b}{c^{3}}&-\frac{2ab^{2}}{c^{3}}\\\
\frac{2a^{2}b}{c^{3}}&\frac{-a^{4}+a^{2}b^{2}+b^{2}c^{2}}{c^{4}}&\frac{ab(a^{2}-b^{2}+c^{2})}{c^{4}}\\\
-\frac{2ab^{2}}{c^{3}}&\frac{ab(a^{2}-b^{2}+c^{2})}{c^{4}}&\frac{b^{4}+a^{2}(c-b)(c+b)}{c^{4}}\end{array}\right),$
(77) $\displaystyle
G_{2}=v_{1}\left(v_{1}\right)^{{\dagger}}-v_{2}\left(v_{2}\right)^{{\dagger}}+v_{3}\left(v_{3}\right)^{{\dagger}}$
$\displaystyle=\left(\begin{array}[]{ccc}\frac{(b+a)(b-a)}{c^{2}}&-\frac{2a^{2}b}{c^{3}}&\frac{2ab^{2}}{c^{3}}\\\
-\frac{2a^{2}b}{c^{3}}&\frac{a^{4}-a^{2}b^{2}+b^{2}c^{2}}{c^{4}}&\frac{ab(-a^{2}+b^{2}+c^{2})}{c^{4}}\\\
\frac{2ab^{2}}{c^{3}}&\frac{ab(-a^{2}+b^{2}+c^{2})}{c^{4}}&\frac{-b^{4}+a^{2}(b^{2}+c^{2})}{c^{4}}\end{array}\right),$
(81)

and,

$\displaystyle
G_{3}=v_{1}\left(v_{1}\right)^{{\dagger}}+v_{2}\left(v_{2}\right)^{{\dagger}}-v_{3}\left(v_{3}\right)^{{\dagger}}$
$\displaystyle=\left(\begin{array}[]{ccc}\frac{a^{2}+b^{2}}{c^{2}}&0&0\\\
0&\frac{a^{4}+a^{2}b^{2}-b^{2}c^{2}}{c^{4}}&-\frac{ab(a^{2}+b^{2}+c^{2})}{c^{4}}\\\
0&-\frac{ab(a^{2}+b^{2}+c^{2})}{c^{4}}&\frac{b^{4}+a^{2}(b-c)(b+c)}{c^{4}}\end{array}\right),$
(85)

where $G_{1}^{2}=G_{2}^{2}=G_{3}^{2}=\textrm{diag}(1,1,1)$.

The neutrino mass matrix obtained from $\left(U_{\rm{PM}}\right)_{\rm{C1}}$ is
written as follows:

$\displaystyle
M_{\left(U_{\rm{PM}}\right)_{\rm{C1}}}=\left(U_{\rm{PM}}\right)_{\rm{C1}}^{\ast}\textrm{diag}(m_{1},m_{2},m_{3})\left(U_{\rm{PM}}\right)_{\rm{C1}}^{{\dagger}}$
$\displaystyle=\left(\begin{array}[]{ccc}A&B&C\\\
B&D&\frac{a^{4}B-b^{4}B-a^{3}bC-a^{2}b\left\\{bB+c\left(A-D\right)\right\\}}{ab^{2}c}\\\
C&\frac{a^{4}B-b^{4}B-a^{3}bC-a^{2}b\left\\{bB+c\left(A-D\right)\right\\}}{ab^{2}c}&\frac{a^{5}B+ab^{3}\left(cA-
bB\right)-a^{4}bC-b^{5}C-a^{3}b\left\\{2bB+c\left(A-D\right)\right\\}}{ab^{3}c}\end{array}\right).$
(89)

The neutrino mass matrix obtained from $\left(U_{\rm{PM}}\right)_{\rm{C2}}$ is
written as follows:

$\displaystyle M_{\left(U_{\rm{PM}}\right)_{\rm{C2}}}$ $\displaystyle=$
$\displaystyle\left(U_{\rm{PM}}\right)_{\rm{C2}}^{\ast}\textrm{diag}(m_{1},m_{2},m_{3})\left(U_{\rm{PM}}\right)_{\rm{C2}}^{{\dagger}}$
$\displaystyle=\left(\begin{array}[]{ccc}A&B&C\\\
B&D&\frac{a^{3}B+b^{3}C+abc\left(D-A\right)}{b^{2}c}\\\
C&\frac{a^{3}B+b^{3}C+abc\left(D-A\right)}{b^{2}c}&\frac{a^{5}B+ab^{3}\left(cA+bB\right)+a^{4}bC+2a^{2}b^{3}C-b^{5}C+a^{3}bc\left(D-A\right)}{ab^{3}c}\end{array}\right).$
(93)

Considering that these neutrino mass matrices,
$M_{\left(U_{\rm{PM}}\right)_{\rm{C1}}}$ and
$M_{\left(U_{\rm{PM}}\right)_{\rm{C2}}}$, were satisfied through
transformation

$\displaystyle
G_{1}^{T}M_{\left(U_{\rm{PM}}\right)_{\rm{C1}}}G_{1}=M_{\left(U_{\rm{PM}}\right)_{\rm{C1}}},$
(94)

and

$\displaystyle
G_{2}^{T}M_{\left(U_{\rm{PM}}\right)_{\rm{C2}}}G_{2}=M_{\left(U_{\rm{PM}}\right)_{\rm{C2}}},$
(95)

they were invariant under $Z_{2}$ symmetry.

##  VI Summary

Here, we focused on primitive Pythagorean triples Williams2017 ;
MartensarXiv2112 ; Kocik2007 ; SchmelzerarXiv2021 ; YekutieliarXiv2021 ;
Cha2018 ; DeriyarXiv2023 ; PricearXiv2008 ; AlperinarXiv2000  to construct a
neutrino mixing model. Primitive Pythagorean triples refer to sets of three
natural numbers $(a,b,c)$ satisfying $c^{2}=a^{2}+b^{2}$. For instance, in the
triangle, $(3,4,5)$, the three internal angles are
$(36.87^{\circ},53.13^{\circ},90.00^{\circ})$. Among these,
$(36.87^{\circ},53.13^{\circ})$ closely approximate the upper limit of
$3\sigma$ for the solar neutrino mixing angle, $35.74^{\circ}$, and the upper
limit of the $3\sigma$ atmospheric mixing angle,
$51.0^{\circ}~{}\left(51.5^{\circ}\right)$, in the case of NO (IO).

First, the relationship between primitive Pythagorean triples and neutrino
mixing angles, $\theta_{12}$ and $\theta_{23}$, is as follows:

$\displaystyle s_{12}=\frac{a}{c},~{}~{}s_{23}=\frac{b}{c},$ (96)

where the relationship among $a$, $b$, and $c$ is described by $a<b<c$. From
this relationship, we constructed the neutrino mixing matrix $U_{\rm{PM}}$.
However, the mixing matrix predicts the vanishing reactor mixing angle. To
improve the reproducibility of the reactor mixing angle, we modified
$\left(U_{\rm{PM}}\right)_{\rm{C1}}=U_{\rm{PM}}U_{23}$
($\left(U_{\rm{PM}}\right)_{\rm{C2}}=U_{\rm{PM}}U_{13}$) where $U_{23}$ and
$U_{13}$ are

$\displaystyle U_{23}=\left(\begin{array}[]{ccc}1&0&0\\\
0&\cos{\theta}&\sin{\theta}e^{-i\phi}\\\
0&-\sin{\theta}e^{i\phi}&\cos{\theta}\end{array}\right)$ (100)

and

$\displaystyle
U_{13}=\left(\begin{array}[]{ccc}\cos{\theta}&0&\sin{\theta}e^{-i\phi}\\\
0&1&0\\\ -\sin{\theta}e^{i\phi}&0&cos{\theta}\end{array}\right).$ (104)

Consequently, we obtained a new neutrino mixing model related to primitive
Pythagorean triples that satisfies the observed values.

We showed that the neutrino mass matrix,
$M_{\left(U_{\rm{PM}}\right)_{\rm{C1}}}$
$\left(M_{\left(U_{\rm{PM}}\right)_{\rm{C2}}}\right)$, obtained from
$\left(U_{\rm{PM}}\right)_{\rm{C1}}$
$\left(\left(U_{\rm{PM}}\right)_{\rm{C2}}\right)$ was invariant under $Z_{2}$
symmetry.

Finally, in this paper, we modified $U_{\rm{PM}}$ into $U_{\rm{PM}}U_{23}$
($U_{\rm{PM}}U_{13}$). We would like to consider the charged lepton mixing,
$U_{e}$, as in Ref. PetcovNPB2015  and the texture of the neutrino mass matrix
in the future. In the present numerical calculation results, the simplest
primitive Pythagorean triple, $(3,4,5)$, did not match the observed values.
The simplest primitive Pythagorean triple $(3,4,5)$ shared relations with
other primitive Pythagorean triples, as observed in Ref. Cha2018 ;
DeriyarXiv2023 ; PricearXiv2008 ; AlperinarXiv2000 . Such relationships are
referred to as the trees of primitive Pythagorean triples. We would like to
investigate the relationship between the trees of primitive Pythagorean
triples and the neutrino mixing in the future.

## References

  * (1) K. R. Williams, TRIPLES: Applictions of Pythagorean Triples (Inspiration Books, 2017)
  * (2) G. J. Martens, arXiv:2112.09553 [math.GM].
  * (3) J. Kocik,Advances in Applied Clifford Algebras. 17, 71-93 (2007).
  * (4) A. Schmelzer, and S. Chetty, arXiv:2102.03451 [math.NT].
  * (5) A. Yekutieli, arXiv:2101.12166 [math.NT].
  * (6) B. Cha, E. Nguuyen, and B. Tauber, J. of Number Theor. 185, 218-256 (2018).
  * (7) P. Deriy, arXiv:2310.15174 [math.HO].
  * (8) H. Lee. Price, arXiv:0809.4324 [math.HO].
  * (9) R. C. Alperin, arXiv:0010281 [math.HO].
  * (10) B. Pontecorvo, Sov. Phys. JETP. 6, 429 (1957).
  * (11) B. Pontecorvo, Sov. Phys. JETP. 7, 172 (1958).
  * (12) Z. Maki, M. Nakagawa and S. Sakata, Prog. Theor. Phys. 28, 870 (1962).
  * (13) M. Tanabashi et al. (Particle Data Group), Phys. Rev. D 98, 030001 (2018).
  * (14) C. Jarlskog, Phys. Rev. Lett. 55, 1039 (1985).
  * (15) I. Esteban, M. C. Gonzalez-Garcia, M. Maltoni, T. Schwetz, and A. Zhou, J. High Energy Phys. 09, 178 (2020). See also, NuFIT 5.2 (2022), www.nu-fit.org.
  * (16) Z. Z. Xing and S. Zhou, Phys. Lett. B 653, 278 (2007).
  * (17) J. D. Bjorken, P. F. Harrison, and W. G. Scott, Phys. Rev. D 74, 073012 (2006).
  * (18) K. Pearson, Proc. R. Soc. Lond. 58, 240 (1895).
  * (19) G. Altarelli and F. Feruglio, Rev. Mod. Phys. 82, 2701 (2010).
  * (20) C. S. Lam, Phys. Rev. D 74, 113004 (2006).
  * (21) C. S. Lam, Phys. Lett. B 656, 193 (2007).
  * (22) C. S. Lam, Phys. Rev. Lett. 101, 121602 (2008).
  * (23) D. A. Dicus, S. -F. Ge, and W. W. Repko, Phys. Rev. D 83, 093007 (2011).
  * (24) S. T. Petcov, Nucl. Phys. B 892, 400-428 (2015).