# The curvature-induced gauge potential and the geometric momentum for a
particle on a hypersphere

Z. Li  School for Theoretical Physics, School of Physics and Electronics,
Hunan University, Changsha 410082, China  L. Q. Lai  School for Theoretical
Physics, School of Physics and Electronics, Hunan University, Changsha 410082,
China  Y. Zhong<EMAIL_ADDRESS>School for Theoretical Physics, School of
Physics and Electronics, Hunan University, Changsha 410082, China  Q. H. Liu
<EMAIL_ADDRESS>School for Theoretical Physics, School of Physics and
Electronics, Hunan University, Changsha 410082, China

###### Abstract

A particle that is constrained to freely move on a hyperspherical surface in
an $N\left(\geq 2\right)$ dimensional flat space experiences a curvature-
induced gauge potential, whose form was given long ago (J. Math. Phys.
34(1993)2827). We demonstrate that the momentum for the particle on the
hypersphere is the geometric one including the gauge potential and its
components obey the commutation relations $\left[p_{i},p_{j}\right]=-i\hbar
J_{ij}/r^{2}$, in which $\hbar$ is the Planck’s constant, and $p_{i}$
($i,j=1,2,3,...N$) denotes the $i-$th component of the geometric momentum, and
$J_{ij}$ specifies the $ij-$th component of the generalized angular momentum
containing both the orbital part and the coupling of the generators of
continuous rotational symmetry group $SO(N-1)$ and curvature, and $r$ denotes
the radius of the $N-1$ dimensional hypersphere.

curvature, gauge potential, hypersphere, geometric momentum

##  I Introduction

In quantum mechanics, a constrained dynamical system is usually associated
with a gauge structure, and it is quite well-understood in, for instance,
gravitational field RMP1957 ; Hehl , condensed matter physics 2010PR ; Iorio ,
quantum fields weinberg  and particle physics Qiugg . We are recently
interested in the constrained motion, i.e., a particle remains and freely
moves on a hypersurface liu07 ; liu11 ; liu13 ; liu132 ; liu133 ; liu14 ;
liu15 ; liu16 ; liu18 ; liu19 ; liu20 ; liu21 . There are also a lot of papers
paying attention to the curvature-induced gauge structure for the system, and
the form of the gauge potential is well-known japan1993 ; japan1993-2 ;
japan1995 ; italy1995 ; japan1997 ; japan1997-2 ; japan1997-3 ; japan1998 ;
usa2003 ; pra14 ; soiqm ; japan2007 ; yang . The form of the gauge potential
is clearly a coupling of the curvature and all generators of a rotational
symmetry group, which in quantum mechanics is not necessarily understood to
represent the spin though it is the important situation. Even so, the gauge
potential is far beyond fully understood. For instance, once a classical
bracket (c.f. Eq. (13)) involves both the momentum and the orbital angular
momentum, the corresponding quantum commutator (c.f. Eq. (19)) based on the
quantization of the classical bracket contains both the momentum operator and
the angular momentum operator, and it appears transparent. However, the
quantum commutator is highly non-trivial because the momentum explicitly
contains the extrinsic curvature, and the previous studies unintentionally
missed this geometric nature of the momentum till 2011 when the geometric
momentum came into sight liu11 . The main aim of the present study is to show
that the relation (19) holds true with reasonable inclusion of the gauge
potential into the geometric momentum.

This paper is organized in the following. In section II, how Ohnuki and
Kitakado obtained the gauge potential is outlined and commented. In section
III, the Dirac formalism of quantization for the particle on the sphere is
invoked, where the Dirac brackets between the momentum components and the
orbital angular momentum components play central role. To show that the
fundamental quantum conditions in the Dirac scheme of quantization are
completely compatible with the gauge potential, we must utilize the proper
form of the momentum. In section IV and V, we explicitly prove two sets of
fundamental quantum conditions in the Dirac scheme, respectively. The final
section VI is a brief conclusion and discussion.

##  II Ohnuki and Kitakado gauge potential on hyperspherical surface

By a hyperspherical surface $\mathit{S}^{N-1}$ ($N\geqslant 2$), we mean a
spherical surface in $N$ dimensional flat space $E^{N}$. The surface equation
is given by,

$x_{i}x_{i}=r^{2},(i=1,2,3,...,N),\text{ and }n_{i}=\frac{x_{i}}{r},$ (1)

where $r$ denotes the radius of the sphere, and $x_{i}$ is the $i$-th
coordinate and $n_{i}$ is the $i$-th component of the normal vector.
Throughout the paper, the Einstein’s summation convention is adopted, which
implies the indices repeated twice in a term are summed over the range of the
index unless specified, and the Roman indices $i,j,k,l,m,n$ run from $1$ to
$N$ to specify the coordinates of the surface in $E^{N}$, while the Greek ones
$\alpha,\beta,\lambda,\mu$ run from $1$ to $N-1$ which can not be confused
with the $N-1$ local coordinates of the surface.

Ohnuki and Kitakado started from a _fundamental algebra_ on the flat $E^{N}$
space japan1993 , the corresponding Euclidean group $E(N)$ whose generators
$\left\\{x_{i},J_{jk}\right\\}$ satisfy following relations,

$\displaystyle[x_{i},x_{j}]$ $\displaystyle=$ $\displaystyle
0,\left[x_{k},J_{ij}\right]=i\hbar\left(x_{i}\delta_{kj}-x_{j}\delta_{ki}\right),$
(2) $\displaystyle\left[J_{jk},J_{lm}\right]$ $\displaystyle=$ $\displaystyle
i\hbar\left(\delta_{jl}J_{km}-\delta_{jm}J_{kl}+\delta_{km}J_{jl}-\delta_{kl}J_{jm}\right).$
(3)

The $SO(N)$ (3) is a subgroup of $E(N)$, whose generators are the generalized
angular momentum,

$J_{ij}=L_{ij}+f_{ij},$ (4)

where $L_{ij}$ $=x_{i}(-i\hbar\partial_{j})-x_{j}(-i\hbar\partial_{i})$ is
$ij-$th component of the usual orbital angular momentum, and $f_{ij}$ is
related to $S_{\alpha\beta}$ that obey commutation relations (3) as well and
constitute the irreducible representation of $SO(N-1)$ group japan1993 ,

$\displaystyle f_{\alpha\beta}$ $\displaystyle=$ $\displaystyle-
f_{\beta\alpha}=S_{\alpha\beta},$ (5) $\displaystyle f_{\alpha N}$
$\displaystyle=$ $\displaystyle-
f_{N\alpha}=-\frac{r}{r+x_{N}}S_{\alpha\beta}n_{\beta}.$ (6)

Ohnuki and Kitakado obtained the form of $f_{ij}$ by means of consideration of
the transportation of a state from any point to another on $\mathit{S}^{N-1}$
by successive infinitesimal unitary transformation japan1993 . The gauge
structure is rooted in a coupling between the generators $S_{\alpha\beta}$ of
rotational symmetry group $SO(N-1)$ and curvature via relation (6). To be
explicit, the general form of the gauge potential $\mathbf{A}$ is in the
following japan1993 ,

$A_{i}=\frac{1}{r}f_{ij}n_{j}=\frac{1}{r+x_{N}}\left\\{\begin{array}[]{cc}S_{i\beta}n_{\beta},&i\neq
N\\\ 0,&i=N\end{array}\right..$ (7)

In the following, we point out three elementary properties of the gauge field
$\mathbf{A}$. 1. Commutators $\left[A_{\alpha},A_{\beta}\right]$ between the
two different components of the gauge potential $\mathbf{A}$ indicate that the
gauge field is non-Abelian,

$\displaystyle\left[A_{\alpha},A_{\beta}\right]$ $\displaystyle=$
$\displaystyle\frac{n_{\alpha}n_{\beta}}{\left(r+x_{N}\right)^{2}}\left[S_{\alpha\mu},S_{\beta\nu}\right]$
(8a) $\displaystyle=$
$\displaystyle\frac{i\hbar}{\left(r+x_{N}\right)^{2}}\left(\left(1-\frac{x_{N}^{2}}{r^{2}}\right)S_{\alpha\beta}+\left(n_{\alpha}S_{\beta\mu}-n_{\beta}S_{\alpha\mu}\right)n_{\mu}\right)$
(8b) $\displaystyle=$ $\displaystyle
i\hbar\left(\frac{\left(r-x_{N}\right)}{r^{2}\left(r+x_{N}\right)}S_{\alpha\beta}+\frac{1}{r+x_{N}}\left(n_{\alpha}A_{\beta}-n_{\beta}A_{\alpha}\right)\right)$
(8c) $\displaystyle=$ $\displaystyle
i\hbar\left(\frac{\left(r-x_{N}\right)}{r^{2}\left(r+x_{N}\right)}S_{\alpha\beta}+\frac{1}{r+x_{N}}\left(\mathbf{n}\wedge\mathbf{A}\right)_{\alpha\beta}\right)\neq
0.$ (8d) The field strength and the gauge transformation of the gauge field
$\mathbf{A}$ are discussed in Ref. japan1997-3 .

2\. It is easily to check that the $N$-th component of the gauge potential
$\mathbf{A}$ is zero, i.e.,
$A_{N}=n_{\alpha}S_{\alpha\beta}n_{\beta}/(r+x_{N})=0$, which amounts to
choose a convenient gauge condition. More importantly, the gauge potential
$\mathbf{A}$ is defined on the tangential plane on $\mathit{S}^{N-1}$ for we
have,

$\mathbf{A\cdot n}=0\text{, or, }A_{\alpha}n_{\alpha}=0.$ (9)

3\. A very important case is: Once $S_{\alpha\beta}$ is understood as spin
matrix, the gauge field $\mathbf{A}$ can be interpreted as spin-curvature
coupling, where the curvature part $n_{\beta}/\left(r+x_{N}\right)$ in (7) is
related to spin-connections, which is in detail discussed in Ref. japan1997-2
. Once the radius $r$ is infinitely large, we have
$n_{\beta}/\left(r+x_{N}\right)\rightarrow 0$. Thus, the spin is not directly
the usual intrinsic angular momentum of the particle defined in flat space.

##  III Dirac formalism of quantization on the sphere $\mathit{S}^{N-1}$ and
an $SO(N,1)$ algebra

Let us now consider following equation of the hyperspherical surface equation,

$f\left(x\right)=\frac{1}{2r}\left(x_{i}^{2}-r^{2}\right)=0.$ (10)

From it we know $\mathbf{n}=\nabla f$. Dirac formalism for constrained motion
gives following Dirac brackets liu13 ,

$\displaystyle\left[x_{i},x_{j}\right]_{D}$ $\displaystyle=$ $\displaystyle
0,$ (11) $\displaystyle\left[x_{i},p_{j}\right]_{D}$ $\displaystyle=$
$\displaystyle\left(\delta_{ij}-n_{i}n_{j}\right),$ (12)
$\displaystyle\left[p_{i},p_{j}\right]_{D}$ $\displaystyle=$
$\displaystyle-\frac{L_{ij}}{r^{2}}=-\frac{x_{i}p_{j}-x_{j}p_{i}}{r^{2}},$
(13)

where $L_{ij}=x_{i}p_{j}-x_{j}p_{i}$ is the usual $ij$-component of the
orbital angular momentum satisfying,

$\displaystyle\left[x_{k},L_{ij}\right]_{D}$ $\displaystyle=$
$\displaystyle\left(x_{i}\delta_{kj}-x_{j}\delta_{ki}\right),$ (14)
$\displaystyle\left[p_{k},L_{ij}\right]_{D}$ $\displaystyle=$
$\displaystyle\left(p_{i}\delta_{kj}-p_{j}\delta_{ki}\right),$ (15)
$\displaystyle\left[L_{jk},L_{lm}\right]_{D}$ $\displaystyle=$
$\displaystyle\left(\delta_{jl}L_{km}-\delta_{jm}L_{kl}+\delta_{km}L_{jl}-\delta_{kl}L_{jm}\right).$
(16)

When transition to quantum mechanics, we utilize the full Dirac formalism of
quantization procedure to define a quantum commutator for any pair of
operators $f$ and $g$ is given by $\left[f,g\right]=i\hbar$
$O\left\\{\left[f,g\right]_{D}\right\\}$. Thus, the proper form of quantum
conditions from (11) to (16) must be, together with (2) and (3),

$\displaystyle\left[x_{i},x_{j}\right]$ $\displaystyle=$ $\displaystyle 0,$
(17) $\displaystyle\left[x_{i},p_{j}\right]$ $\displaystyle=$ $\displaystyle
i\hbar\left(\delta_{ij}-n_{i}n_{j}\right),$ (18)
$\displaystyle\left[p_{i},p_{j}\right]$ $\displaystyle=$
$\displaystyle-i\hbar\frac{J_{ij}}{r^{2}}\neq-\frac{i\hbar}{r^{2}}\left(x_{i}p_{j}-x_{j}p_{i}\right),$
(19) $\displaystyle\left[p_{k},J_{ij}\right]$ $\displaystyle=$ $\displaystyle
i\hbar\left(p_{i}\delta_{kj}-p_{j}\delta_{ki}\right).$ (20)

Here, we replace the symbol $L_{ij}$ by $J_{ij}$, which denotes the
generalized angular momentum, then the quantum momentum $p_{i}$ in our problem
is not the usual one. It is clear that operators ($p_{i},J_{ij}$) form the
$SO(N,1)$ algebra from relations (3), (19) and (20).

We will show that with $\mathbf{p}=\mathbf{\Pi-A}$, relations (19) and (20)
are valid, in which the $\mathbf{A}$-independent part $\mathbf{\Pi}$ is a
well-defined geometric momentum liu07 ; liu11 ; liu13 ; liu132 ; liu133 ;
liu14 ; liu15 ; liu16 ; liu18 ; liu19 ; liu20 ; liu21 . In general, the
geometric momentum takes the form, $\mathbf{\Pi=}\Pi_{i}\mathbf{e}_{i}=$
$-i\hbar({\nabla_{\Sigma}}+{M{\mathbf{n/}}}2)$ with
$\nabla_{S}\equiv\nabla-\mathbf{n}\left(\mathbf{n\cdot}\nabla\right)=\nabla-\mathbf{n}\partial_{n}$
denoting the gradient operator on $\mathit{S}^{N-1}$, and the mean curvature
${M\equiv-\nabla_{\Sigma}}\cdot\mathbf{n}$ is defined by the sum of the all
principal curvatures. Notice that the mean curvature ${M}$ is an extrinsic
curvature, this form of momentum $\mathbf{p}=\mathbf{\Pi-A}$ is fundamentally
different from the canonical ones in curvilinear coordinates for $\mathbf{p}$
depends on the geometric invariants whereas the the canonical ones do not.
However, for a particle on the hyperspherical surface, the geometric momentum
assumes a simpler form which can be obtained via a simple quantization of
relation $p_{i}=n_{j}L_{ji}/r$ which in quantum mechanics turns out to be,

$\Pi_{i}=\frac{1}{2r}\left(n_{j}L_{ji}+L_{ji}n_{j}\right),\text{ i.e.,
}\mathbf{\Pi}=\mathbf{-}i\hbar\left(\nabla_{S}-(N-1)\frac{\mathbf{n}}{2r}\right),$
(21)

where the fact $\left[n_{j},L_{ji}\right]\neq 0$ is noted. Now we have
$J_{ij}=L_{ij}+f_{ij}$ (4). The similar quantization of the classical relation
$p_{i}=n_{j}L_{ji}/r$ gives explicitly how the gauge potential enters in the
geometric momentum,

$p_{i}=\frac{1}{2r}\left(n_{j}J_{ji}+J_{ji}n_{j}\right)=\Pi_{i}-A_{i},$ (22)

where $A_{i}=f_{ij}n_{j}/r$ comes from its definition (7). What is more, the
momentum still lies on the tangential plane on the surface for we have,

$\mathbf{p\cdot n+n\cdot p}{\mathbf{=0,}}$ (23)

where $A_{i}n_{i}=0$ directly arises from its definition (7) as well, and
$\mathbf{\Pi\cdot n+n\cdot\Pi=0}$ is easily verifiable with definition of
$\mathbf{\Pi}$ (21).

There is a curious fact. The form $\left[p_{i},p_{j}\right]_{D}=-L_{ij}/r^{2}$
(13) is seldom used, and it is usually written in the form
$\left[p_{i},p_{j}\right]_{D}=-\left(x_{i}p_{j}-x_{j}p_{i}\right)/r^{2}$ whose
quantization assumes the form
$\left[p_{i},p_{j}\right]=-i\hbar\left(x_{i}p_{j}-x_{j}p_{i}\right)/r^{2}$,
see for instance Weinberg’s lectures on quantum mechanics Weinberg1  and also
Ohnuki and Kitakado paper japan1997-2 , in which only the orbital angular
momentum presents. However, if starting from
$\left[p_{i},p_{j}\right]_{D}=-L_{ij}/r^{2}$ (13) which is explicitly
expressed in terms of orbital angular momentum, the quantum angular momentum
allows the spin-related angular momentum to enter because the angular momentum
in quantum mechanics is determined by the algebra. In other words, the quantum
conditions assume the form $\left[p_{i},p_{j}\right]=-i\hbar J_{ij}/r^{2}$
(19) in which $J_{ij}$ is known, then the proper form the $p_{i}$ is sought.

In section IV, we will show $\left[p_{i},p_{j}\right]=-i\hbar J_{ij}/r^{2}$
(19) with $p_{i}$ defined in (22); and in section V, we will prove
$\left[p_{k},J_{ij}\right]=i\hbar\left(p_{i}\delta_{kj}-p_{j}\delta_{ki}\right)$
(20).

##  IV A proof of $\left[p_{i},p_{j}\right]=-i\hbar J_{ij}/r^{2}$

We start from $\mathbf{p}=\mathbf{\Pi-A}$, from which the commutators
$\left[p_{i},p_{l}\right]$ are,

$\left[p_{i},p_{l}\right]=\left[\Pi_{i}-A_{i},\Pi_{l}-A_{l}\right].$ (24)

We split $N$ components $p_{i}$ of the momentum $\mathbf{p}$ into two
categories of $p_{\alpha}$ $\left(\alpha=1,2,\cdots,N-1\right)$ and $p_{N}$,

$\displaystyle p_{\alpha}$ $\displaystyle=$
$\displaystyle\Pi_{\alpha}+\frac{1}{r^{2}}x_{i}f_{i\alpha}=\Pi_{\alpha}-A_{\alpha},$
(25) $\displaystyle p_{N}$ $\displaystyle=$
$\displaystyle\Pi_{N}-A_{N}=\Pi_{N},$ (26)

where we used $A_{N}=0$. The calculations of $\left[p_{i},p_{l}\right]$ will
be done separately, and we first study $\left[p_{\alpha},p_{\beta}\right]$ and
secondly deal with $\left[p_{\alpha},p_{N}\right]$.

In order to compute $\left[p_{\alpha},p_{\beta}\right]$, we split the
commutator $\left[p_{\alpha},p_{\beta}\right]$ into four parts,

$\left[p_{\alpha},p_{\beta}\right]=\left[\Pi_{\alpha}-A_{\alpha},\Pi_{\beta}-A_{\beta}\right]=\left[\Pi_{\alpha},\Pi_{\beta}\right]-\left[\Pi_{\alpha},A_{\beta}\right]-\left[A_{\alpha},\Pi_{\beta}\right]+\left[A_{\alpha},A_{\beta}\right].$
(27)

In the right-hand-side of above expression, the first part is known liu13  and
the result is given by,

$\left[\Pi_{i},\Pi_{j}\right]=-\frac{i\hbar}{r^{2}}L_{ij}.$ (28)

The second and the third part is essentially the same for
$\left[\Pi_{\alpha},A_{\beta}\right]=-\left[A_{\beta},\Pi_{\alpha}\right]$,
and the last one is $\left[A_{\alpha},A_{\beta}\right]$ which is also
understood via Eq. (8c). To carry out $\left[\Pi_{\alpha},A_{\beta}\right]$,
we start from $\left[\mathbf{\Pi},A_{\beta}\right]$ and then project onto the
$\alpha$-direction. To note a relation
$\left[\mathbf{\Pi},A_{\beta}\right]=-i\hbar\left[\nabla-\mathbf{n}\left(\mathbf{n\cdot}\nabla\right),A_{\beta}\right]$,
we need thus to compute the commutator $\left[\nabla,A_{\beta}\right]$,

$\displaystyle\left[\nabla,A_{\beta}\right]$ $\displaystyle=$
$\displaystyle\nabla A_{\beta}$ (29) $\displaystyle=$
$\displaystyle\nabla\left(\frac{n_{\mu}}{r+x_{N}}\right)S_{\beta\mu}$
$\displaystyle=$
$\displaystyle\left(-\left(\frac{n_{\mu}\left(n_{N}+2\right)}{\left(r+x_{N}\right)^{2}}\right)\mathbf{n-}\left(\frac{n_{\mu}}{\left(r+x_{N}\right)^{2}}\right)\mathbf{e}_{N}+\left(\frac{1}{r\left(r+x_{N}\right)}\right)\mathbf{e}_{\mu}\right)S_{\beta\mu}.$

Then $\mathbf{n\cdot}\nabla A_{\beta}$ is with use of $\mathbf{n\cdot
e}_{N}=n_{N}$ and $\mathbf{n\cdot e}_{\mu}=n_{\mu}$,

$\mathbf{n\cdot}\nabla
A_{\beta}=\left(-\frac{2n_{\mu}\left(n_{N}+1\right)}{\left(r+x_{N}\right)^{2}}\mathbf{+}\frac{n_{\mu}}{r\left(r+x_{N}\right)}\right)S_{\beta\mu}.$
(30)

Thus, we obtain,

$\displaystyle\left[\nabla-\mathbf{n}\left(\mathbf{n\cdot}\nabla\right),A_{\beta}\right]$
$\displaystyle=$ $\displaystyle\nabla
A_{\beta}-\mathbf{n}\left(\mathbf{n\cdot}\nabla\right)A_{\beta}$ (31)
$\displaystyle=$
$\displaystyle\left(\left(\frac{n_{\mu}n_{N}}{\left(r+x_{N}\right)^{2}}\mathbf{-}\frac{n_{\mu}}{r\left(r+x_{N}\right)}\right)\mathbf{n}-\left(\frac{n_{\mu}}{\left(r+x_{N}\right)^{2}}\right)\mathbf{e}_{N}+\left(\frac{1}{r\left(r+x_{N}\right)}\right)\mathbf{e}_{\mu}\right)S_{\beta\mu}$
$\displaystyle=$
$\displaystyle\left(\mathbf{-}\frac{n_{\mu}}{\left(r+x_{N}\right)^{2}}\mathbf{n}-\frac{n_{\mu}}{\left(r+x_{N}\right)^{2}}\mathbf{e}_{N}+\frac{1}{r\left(r+x_{N}\right)}\mathbf{e}_{\mu}\right)S_{\beta\mu}.$

Projecting it onto the $\alpha$-direction, we get
$\left[\Pi_{\alpha},A_{\beta}\right]$, with
$\mathbf{e}_{\alpha}\cdot\mathbf{n=}$ $n_{\alpha}$,
$\mathbf{e}_{\alpha}\cdot\mathbf{e}_{N}=0$, and
$\mathbf{e}_{\alpha}\cdot\mathbf{e}_{\mu}=\delta_{\alpha\mu}$,

$\displaystyle\left[\Pi_{\alpha},A_{\beta}\right]$ $\displaystyle=$
$\displaystyle-i\hbar\mathbf{e}_{\alpha}\cdot\left[\nabla-\mathbf{n}\left(\mathbf{n\cdot}\nabla\right),A_{\beta}\right]$
(32) $\displaystyle=$
$\displaystyle-i\hbar\left(\mathbf{-}\frac{n_{\mu}}{\left(r+x_{N}\right)^{2}}n_{\alpha}+\frac{1}{r\left(r+x_{N}\right)}\delta_{\alpha\mu}\right)S_{\beta\mu}$
$\displaystyle=$
$\displaystyle-i\hbar\left(-\frac{n_{\alpha}}{r+x_{N}}A_{\beta}+\frac{S_{\beta\alpha}}{r\left(r+x_{N}\right)}\right).$

So, two terms in right-hand-side of (27)
$\left[\Pi_{\alpha},A_{\beta}\right]+\left[A_{\alpha},\Pi_{\beta}\right]$ is,

$\displaystyle\left[\Pi_{\alpha},A_{\beta}\right]+\left[A_{\alpha},\Pi_{\beta}\right]$
$\displaystyle=$
$\displaystyle\left[\Pi_{\alpha},A_{\beta}\right]-\left[\Pi_{\beta},A_{\alpha}\right]$
(33) $\displaystyle=$
$\displaystyle-i\hbar\left(-\frac{n_{\alpha}}{r+x_{N}}A_{\beta}+\frac{S_{\beta\alpha}}{r\left(r+x_{N}\right)}-\left(-\frac{n_{\beta}}{r+x_{N}}A_{\alpha}+\frac{S_{\alpha\beta}}{r\left(r+x_{N}\right)}\right)\right)$
$\displaystyle=$
$\displaystyle-i\hbar\left(-\frac{1}{r+x_{N}}\left(n_{\alpha}A_{\beta}-n_{\beta}A_{\alpha}\right)+\frac{2S_{\beta\alpha}}{r\left(r+x_{N}\right)}\right).$

Combining results (8c) and (33), we find a simple result for last three
commutators in (27),

$\displaystyle-\left[\Pi_{\alpha},A_{\beta}\right]-\left[A_{\alpha},\Pi_{\beta}\right]+\left[A_{\alpha},A_{\beta}\right]$
$\displaystyle=$ $\displaystyle
i\hbar\left(\frac{2S_{\beta\alpha}}{r\left(r+x_{N}\right)}-\frac{1}{r+x_{N}}\left(n_{\alpha}A_{\beta}-n_{\beta}A_{\alpha}\right)\right)$
(34)
$\displaystyle+i\hbar\left(\frac{\left(r-x_{N}\right)}{r^{2}\left(r+x_{N}\right)}S_{\alpha\beta}+\frac{1}{r+x_{N}}\left(n_{\alpha}A_{\beta}-n_{\beta}A_{\alpha}\right)\right)$
$\displaystyle=$ $\displaystyle
i\hbar\left(\frac{\left(r-x_{N}\right)}{r^{2}\left(r+x_{N}\right)}-\frac{2}{r\left(r+x_{N}\right)}\right)S_{\alpha\beta}$
$\displaystyle=$ $\displaystyle
i\hbar\left(-\frac{1}{r^{2}}+\left(\frac{1}{r^{2}}-\frac{1}{r\left(r+x_{N}\right)}-\frac{x_{N}}{r^{2}\left(r+x_{N}\right)}\right)\right)S_{\alpha\beta}$
$\displaystyle=$ $\displaystyle-i\hbar\frac{S_{\alpha\beta}}{r^{2}}.$

Thus we have from Eqs. (28) and (34),

$\left[p_{\alpha},p_{\beta}\right]=-\frac{i\hbar}{r^{2}}\left(L_{\alpha\beta}+S_{\alpha\beta}\right)=-\frac{i\hbar}{r^{2}}J_{\alpha\beta}.$
(35)

Similarly, we can prove,

$\left[p_{\alpha},p_{N}\right]=-\frac{i\hbar}{r^{2}}J_{\alpha N}.$ (36)

Combination of two equations (35) and (36) above gives the result (19). Q.E.D.

##  V A proof of
$\left[p_{k},J_{ij}\right]=i\hbar\left(p_{i}\delta_{kj}-p_{j}\delta_{ki}\right)$

The proof of fundamental quantum conditions
$\left[p_{k},J_{ij}\right]=i\hbar\left(p_{i}\delta_{kj}-p_{j}\delta_{ki}\right)$
(20) is also straightforward. To do it, we split $N$ components $p_{i}$ of the
momentum $\mathbf{p}$ into two categories of $p_{\alpha}$
$\left(\alpha=1,2,\cdots,N-1\right)$ and $p_{N}$. In consequence, the
commutators $\left[p_{k},J_{ij}\right]$ are divided into two categories as
$\left[p_{\lambda},J_{ij}\right]$ and $\left[p_{N},J_{ij}\right]$. To process
$\left[p_{\lambda},J_{ij}\right]$, we compute $[p_{\lambda},J_{\alpha\beta}]$
and $\left[p_{\lambda},J_{iN}\right]$ one by one. The detailed steps of
calculation of commutators $[p_{\lambda},J_{\alpha\beta}]$ is given in the
following.

We split commutators $[p_{\lambda},J_{\alpha\beta}]$ into four parts,

$[p_{\lambda},J_{\alpha\beta}]=[\Pi_{\lambda}-A_{\lambda},L_{\alpha\beta}+S_{\alpha\beta}]=[\Pi_{\lambda},L_{\alpha\beta}]+[\Pi_{\lambda},S_{\alpha\beta}]-[A_{\lambda},L_{\alpha\beta}]-[A_{\lambda},S_{\alpha\beta}].$
(37)

The first part of the commutators in the right-hand-side was already treated
liu13  and the results are given by,

$[\Pi_{\lambda},L_{\alpha\beta}]=i\hbar(\Pi_{\alpha}\delta_{\beta\lambda}-\Pi_{\beta}\delta_{\alpha\lambda}).$
(38)

The second part vanishes because of,

$[\Pi_{\lambda},S_{\alpha\beta}]=0.$ (39)

The third part $[A_{\lambda},L_{\alpha\beta}]$ is with $L_{\alpha\beta}\equiv
x_{\alpha}\Pi_{\beta}-x_{\beta}\Pi_{\alpha},$

$\displaystyle[A_{\lambda},L_{\alpha\beta}]$ $\displaystyle=$
$\displaystyle[\frac{S_{\lambda\tau}n_{\tau}}{r+x_{N}},x_{\alpha}\Pi_{\beta}-x_{\beta}\Pi_{\alpha}]$
(40a) $\displaystyle=$ $\displaystyle
x_{\alpha}[\frac{S_{\lambda\tau}n_{\tau}}{r+x_{N}},\Pi_{\beta}]-x_{\beta}[\frac{S_{\lambda\tau}n_{\tau}}{r+x_{N}},\Pi_{\alpha}]$
(40b) $\displaystyle=$ $\displaystyle
i\hbar(\frac{S_{\lambda\beta}n_{\alpha}}{r+x_{N}}-\frac{S_{\lambda\alpha}n_{\beta}}{r+x_{N}}),$
(40c) where we used results (32). The fourth part
$[A_{\lambda},S_{\alpha\beta}]$ in (37) can be easily carried out with use of
the algebraic relations (3),
$[A_{\lambda},S_{\alpha\beta}]=\frac{n_{\tau}}{r+x_{N}}[S_{\lambda\tau},S_{\alpha\beta}]=-i\hbar(\delta_{\alpha\lambda}A_{\beta}-\delta_{\beta\lambda}A_{\alpha}+\frac{n_{\alpha}S_{\lambda\beta}}{r+x_{N}}-\frac{n_{\beta}S_{\lambda\alpha}}{r+x_{N}}).$
(41)

Combination of results for the third and fourth part (40c) and (41) yields a
simple result,

$\displaystyle-[A_{\lambda},L_{\alpha\beta}]-[A_{\lambda},S_{\alpha\beta}]$
$\displaystyle=$
$\displaystyle-i\hbar(\frac{S_{\lambda\beta}n_{\alpha}}{r+x_{N}}-\frac{S_{\lambda\alpha}n_{\beta}}{r+x_{N}})+i\hbar(\delta_{\alpha\lambda}A_{\beta}-\delta_{\beta\lambda}A_{\alpha}+\frac{n_{\alpha}S_{\lambda\beta}}{r+x_{N}}-\frac{n_{\beta}S_{\lambda\alpha}}{r+x_{N}})$
(42a) $\displaystyle=$ $\displaystyle
i\hbar(\delta_{\alpha\lambda}A_{\beta}-\delta_{\beta\lambda}A_{\alpha}).$
(42b) Thus, the commutators $[p_{\lambda},J_{\alpha\beta}]$ give,
$\displaystyle[p_{\lambda},J_{\alpha\beta}]$ $\displaystyle=$ $\displaystyle
i\hbar(\Pi_{\alpha}\delta_{\beta\lambda}-\Pi_{\beta}\delta_{\alpha\lambda})+i\hbar(\delta_{\alpha\lambda}A_{\beta}-\delta_{\beta\lambda}A_{\alpha})$
(43a) $\displaystyle=$ $\displaystyle
i\hbar(p_{\alpha}\delta_{\beta\lambda}-p_{\beta}\delta_{\alpha\lambda}).$
(43b)

Similarly, we can prove following two relations,

$[p_{\lambda},J_{\alpha N}]=-i\hbar p_{N}\delta_{\alpha\lambda},$ (44)

and

$\left[p_{N},J_{ij}\right]=i\hbar\left(p_{i}\delta_{Nj}-p_{j}\delta_{Ni}\right).$
(45)

In final, from results (43b)-(45), we see that the proof of the fundamental
quantum conditions (20) is complete. Q.E.D.

##  VI Conclusions and discussions

For a particle that is constrained on an $\left(N-1\right)$-dimensional
hypersphere, the classical bracket is
$\left[p_{i},p_{j}\right]_{D}=-L_{ij}/r^{2}$. In contrast to the usual
understanding of the corresponding quantum condition remains
$\left[p_{i},p_{j}\right]=-i\hbar L_{ij}/r^{2}$ in which $L_{ij}$ represents
orbital angular momentum, we argue that it must be replaced by
$\left[p_{i},p_{j}\right]=-i\hbar J_{ij}/r^{2}$ in which $J_{ij}$ includes two
parts, and one of them is the usual orbital angular momentum $L_{ij}$ and
another comes from the gauge potential $f_{ij}$. In other words, we adapt the
definition of momentum such that $\left[p_{i},p_{j}\right]=-i\hbar
J_{ij}/r^{2}$ holds true. What is more, we see that another set of the
fundamental quantum conditions
$\left[p_{k},J_{ij}\right]=i\hbar\left(p_{i}\delta_{kj}-p_{j}\delta_{ki}\right)$
comes true as well. The momentum is the geometric one which came into our
sight recently.

The present study gives us an insight into the noncommutativity between
different components of momentum. Our general conjecture is: where there is
such a noncommutativity, there is a generalized angular momentum which is
possibly associated with both the gauge field angular momentum and the orbital
one. The most familiar system is a charged particle in $E^{3}$ under influence
of uniform magnetic field $\mathbf{B}=B_{0}\mathbf{e}_{z}$, where we have
$[p_{x},p_{y}]=i\hbar qB_{0}$ in which $p_{x}$ and $p_{y}$ are $x$\- and
$y$-component of the kinetic momentum $\mathbf{p}=-i\hbar\nabla-q\mathbf{A}$,
and $q$ is the charge. The magnetic field is a gauge field, but where is the
angular momentum? A usually overlooked fact is that the $z$-component of the
field angular momentum is prl , $L_{z}^{field}=qB_{0}\rho^{2}/2$ with
$\rho=\sqrt{x^{2}+y^{2}}$, we immediately find that $[p_{x},p_{y}]=2i\hbar
L_{z}^{field}/\rho^{2}$. Thus, the conjecture is true for the well-known
system. With this understanding at hand, we can boldly guess that the angular
momentum part $f_{ij}$ in (4) to be closely related to the angular momentum of
the gauge field, which is still under investigations.

###### Acknowledgements.

This work is financially supported by National Natural Science Foundation of
China under Grant No. 11675051.

## References

  * (1) D. R. Brill, and J. A. Wheeler, Rev. Mod. Phys. 29(1957)465-479.
  * (2) M. Blagojevic, F. W. Hehl, Gauge Theories of Gravitation, Imperial College Press, London, 2013, pp. 365–380.
  * (3) M. A. H. Vozmediano, M. I. Katsnelson, F. Guinea, Phys. Rep. 496(2010)109-148.
  * (4) A. Iorio and G. Lambiase, Phys. Rev. D 90(2014)025006.
  * (5) S. Weinberg, The quantum theory of fields, Vol. 1, Cambridge University Press, Cambridge, 1995, pp. 369–371.
  * (6) C. Quigg, Gauge of theories of strong, weak and electromagnetic interactions, 2nd Ed. Princeton University Press, 2014, pp. 57–70.
  * (7) Q. H. Liu, C. L. Tong, and M. M. Lai, J. Phys. A 40(2007)4161-4168.
  * (8) Q. H. Liu, L. H. Tang, D. M. Xun, Phys. Rev. A 84(2011)042101.
  * (9) Q. H. Liu, J. Math. Phys. 54(2013)122113.
  * (10) Q. H. Liu, J. Phys. Soc. Jap. 82(2013)104002.
  * (11) D. M. Xun, Q. H. Liu, X. M. Zhu, Ann. Phys. (NY), 338(2013)123-133.
  * (12) D. M. Xun, Q. H. Liu, Ann. Phys. (NY), 341(2014)132-141.
  * (13) Z. S. Zhang, S. F. Xiao, D. M. Xun, and Q. H. Liu, Commun. Theor. Phys. 63(2015)19-24.
  * (14) L. D. Hu, D. K. Lian, Q. H. Liu, Eur. Phys. J. C, 76(2016)655.
  * (15) D. K. Lian, L. D. Hu, Q. H. Liu, Ann. Phys.(Berlin), 530(2018)1700415.
  * (16) Q. H. Liu, Z. Li, X. Y. Zhou, Z. Q. Yang, W. K. Du, Eur. Phys. J. C, 79(2019)712.
  * (17) Z. Q. Yang, X. Y. Zhou, Z. Li, W. K. Du and Q. H. Liu, Phys. Lett. A, 384(2020)126604.
  * (18) Z. Li, X. Yang and Q. H. Liu, Commun. Theor. Phys. 73(2021)025104.
  * (19) Y. Ohnuki and S. Kitakado, J. Math. Phys. 34(1993)2827.
  * (20) K. Fujii and N. Ogawa, Prog. Theor. Phys. 89(1993)575.
  * (21) D. McMullan and I. Tsutsui, Ann. Phys. (NY), 237(1995)269.
  * (22) P. Maraner, Ann. Phys.(NY), 246(1996)325-346.
  * (23) K. Fujii, N. Ogawa, S. Uchiyama, N. M. Chepilko, Int. J. Mod. Phys. A, 12(1997)5235-5277.
  * (24) N. Ogawa, Mod. Phys. Lett. A, 12(1997)1583-1588.
  * (25) M. Hirayama, H. -M. Zhang, T. Hamada, Prog. Theoret. Phys. 97(1997)679–689.
  * (26) H. Ikemori, S. Kitakdo, H. Nakatani, H. Otsu, and T. Sato, Mod. Phys. Lett. A. 13(1998)15-21.
  * (27) P. C. Schuster, R. L. Jaffe, Ann. Phys. (NY), 307(2003)132-143.
  * (28) J. Stockhofe and P. Schmelcher, Phys. Rev. A. 89(2014)033630.
  * (29) C. Cohen-Tannoudji, B. Diu and F. Laloe, Quantum mechanics, Vol. Two, John Wiley & Sons, NY, 1977, pp. 315–328.
  * (30) N. Hatano, R. Shirasaki and H. Nakamura, Phys. Rev. A, 75(2007)032107.
  * (31) C. N. Yang, Phys. Today, 33(1980)42-49.
  * (32) S. Weinberg, Lectures on Quantum Mechanics, 2nd ed. , Cambridge University Press, Cambridge, 2015, pp. 335–340.
  * (33) C. R. Greenshields, R. L. Stamps, S. Franke-Arnold, and S. M. Barnett, Phys. Rev. Lett. 113(2014)240404.


#  Multipolar ordering in the three-orbital Hubbard model

Naoya Chikano  Department of Physics, Saitama University, Saitama 338-8570,
Japan  Shintaro Hoshino  Department of Physics, Saitama University, Saitama
338-8570, Japan  Hiroshi Shinaoka  Department of Physics, Saitama University,
Saitama 338-8570, Japan

###### Abstract

The ground-state phase diagrams of the three-orbital $t_{\mathrm{2g}}$ Hubbard
model are studied using a Hartree-Fock approximation. First, a complete set of
multipolar order parameters for $t_{\mathrm{2g}}$ models defined in terms of
the effective total angular momentum $j_{\mathrm{eff}}$ are theoretically
derived. These order parameters can classify off-diagonal orders between
$j_{\mathrm{eff}}$ = 1/2 and $j_{\mathrm{eff}}$ = 3/2 manifolds. Second,
through extensive Hartree-Fock calculations, the ground-state phase diagrams
in the space of (1) the onsite Coulomb repulsion $U$, (2) the spin-orbit
coupling (SOC) $\lambda$, and (3) the number of electrons $n$ are mapped out.
A variety of nontrivial quantum phases with $j_{\mathrm{eff}}$–diagonal and
$j_{\mathrm{eff}}$–off-diagonal multipole orders are found. Finally, future
studies using more numerically expensive methods, such as dynamical mean-field
theory are discussed.

##  I Introduction

Novel phenomena arising from competing spin-orbit coupling (SOC) and
electronic correlation are matters of considerable interest in condensed
matter physics. In particular, $5d$ transition metal oxide provides a suitable
playground to study such phenomena because these interactions are comparable.

$5d$ transition metal oxides exhibit a variety of exotic quantum phases. For
example, the layered iridate Sr2IrO4 with a $t_{\mathrm{2g}}^{5}$ ($n$=5)
configuration exhibits a nontrivial spin-orbital-entangled Mott insulator[1,
2, 3]. The pyrochlore oxide Cd2Re2O7, with a $t_{\mathrm{2g}}^{2}$
configuration, exhibits various multipolar ordered states [4],
superconductivity at ambient pressure, and under high pressures [5, 6]. For
$n=3$, Cd2Os2O7 shows a finite-temperature metal-insulator transition into a
non-collinear magnetically ordered state [7, 8, 9]. The emergence of novel
magnetic states has been reported in various materials such as Eu2Ir2O7 [10,
11], Ba2YIrO6 [12], and Na2IrO3 [13, 14].

On the theoretical side, the $t_{\mathrm{2g}}$ three-orbital Hubbard model
with SOC is the simplest model for clarifying the novel phenomena, arising
from a competing SOC and electron interaction. Recently, Sato et al. [15, 16,
17] analyzed a model with infinite spatial dimensions for $n=4$ and $5$ using
the dynamical mean-field theory (DMFT), which can accurately describe local
strong correlation effects [18]. They showed that this model hosts many
intriguing quantum phases, such as a spin-orbital-entangled Mott insulating
phase at $n=5$ [15] and a nonmagnetic excitonic insulator [16] with quadrupole
ordering at $n=4$ [17].

However, further exploration of the phase diagram at arbitrary fillings is
computationally expensive. In particular, the hybridization-expansion
continuous-time quantum Monte Carlo method [19, 20] , which was used to solve
the impurity model in previous studies, suffers from a severe sign problem
when the SOC is strong [21, 22]. Although efforts to alleviate the sign
problem in CT-QMC continue [21], it is still difficult to compute the low-$T$
properties of a multi-orbital Hubbard model under SOC, especially in cases
which are slightly away from filling [21]. Thus, it is important to map the
phase diagram using a simpler and computationally feasible method.

Another theoretical issue is the complete classification of the (local)
multipolar order parameters. Multipole representation is suitable for the
description of order parameters in spin-orbit entangled systems as developed
in $f$-electron systems [23, 24, 25, 26]. Although a complete clarification
should involve 36 (= $6\times 6$) distinct order parameters for
$t_{\mathrm{2g}}$ systems, only a subset of order parameters have been used to
classify quantum phases in previous studies [17]. Moreover, these order
parameters do not correctly detect the off-diagonal orders of the so-called
$j_{\mathrm{eff}}=1/2$ and $j_{\mathrm{eff}}=3/2$ manifolds.

In this study, a complete set of multipolar order parameters designated for
classifying spin-orbital entangled states, similar to Ref. [26], are derived.
In particular, the diagonal and off-diagonal matrix elements with respect to
$j_{\mathrm{eff}}$ are focused upon. These order parameters can distinguish
diagonal and off-diagonal orders correctly. Based on this result, the ground-
state phase diagram of the $t_{\mathrm{2g}}$ Hubbard model using the
unrestricted Hartree-Fock approximation is systematically explored and all
possible (particle-conserving) local symmetry-breaking patterns are
considered. The phase diagram in the space of three parameters: the onsite
intra-orbital Coulomb repulsion $U$, the strength of the SOC $\lambda$, and
the number of electrons $n$, are mapped out, thus revealing the presence of
various spin-orbital entangled states and quantum multi-critical points.
Furthermore, the intensity maps of multipolar order parameters as a function
of $U$ and $\lambda$, which reveal several interesting features about the non-
trivial quantum phases, are studied.

The remainder of this paper is organized as follows. In Secs. II and III, the
model and numerical method are introduced. In Sec. IV, the multipole expansion
described in a complete basis set is defined. The relation between the order
parameters introduced in this study and the conventional ones used in previous
studies are also discussed. The phase diagram computed for $n=4$ is
represented in Sec. V. In Sec. VI, cases of general electronic filling and
intensity maps of the order parameters are discussed. In Sec. VII, this paper
is summarized.

##  II Model

Here, a three-orbital $t_{2{\rm g}}$ Hubbard model is considered, whose
Hamiltonian is given by

$\displaystyle\mathcal{H}$
$\displaystyle=\mathcal{H}_{0}+\mathcal{H}_{\mathrm{int}}+\mathcal{H}_{\mathrm{SOC}},$
(1)

where $\mathcal{H}_{0}$ is the non-interacting part,
$\mathcal{H}_{\mathrm{int}}$ is the on-site interaction part, and
$\mathcal{H}_{\mathrm{SOC}}$ is the spin-orbit coupling (SOC). The non-
interacting Hamiltonian $\mathcal{H}_{0}$ corresponds to an orbital-diagonal
semi-circular DOS with full bandwidth $W=4t$, that is,
$D(\omega)=\sqrt{4t^{2}-\omega^{2}}/\pi$. $t=1$ is taken as the unit of
energy.

For the interaction term, the standard, fully rotationally invariant Slater-
Kanamori interactions are taken:

$\displaystyle\mathcal{H}_{{\rm int}}$
$\displaystyle=U\sum_{i\alpha}n_{i\alpha\uparrow}n_{i\alpha\downarrow}+\sum_{i,\alpha>\beta,\sigma}U^{\prime}n_{i\alpha\sigma}n_{i\beta\sigma^{\prime}}$
$\displaystyle\hskip
11.38109pt+\sum_{i,\alpha>\beta,\sigma}(U^{\prime}-J_{{\rm
H}})n_{i\alpha\sigma}n_{i\beta\sigma}$ $\displaystyle\hskip 11.38109pt+J_{{\rm
H}}\sum_{i,\alpha\neq\beta}(c_{i\alpha\downarrow}^{\dagger}c_{i\beta\uparrow}^{\dagger}c_{i\alpha\uparrow}c_{i\beta\downarrow}+c_{i\alpha\uparrow}^{\dagger}c_{i\alpha\downarrow}^{\dagger}c_{i\beta\uparrow}c_{i\beta\downarrow})$
$\displaystyle=\frac{1}{2}\sum_{\alpha\beta\alpha^{\prime}\beta^{\prime}\sigma\sigma^{\prime}}U_{\alpha\beta\alpha^{\prime}\beta^{\prime}}c_{i\alpha\sigma}^{\dagger}c_{i\beta\sigma^{\prime}}^{\dagger}c_{i\beta^{\prime}\sigma^{\prime}}c_{i\alpha^{\prime}\sigma},$
(2)

where $i$ is the site index, $U$ is the intra-orbital Coulomb repulsion,
$U^{\prime}(=U-2J_{{\rm H}})$ is the inter-orbital Coulomb repulsion, and
$J_{{\rm H}}$ is Hund’s coupling. Here, $\alpha(\beta)$ and
$\sigma(\sigma^{\prime})$ are the orbital and spin indices, respectively. The
site index $i$ in the following section is omitted for simplicity. The Coulomb
tensor is defined as
$U_{\alpha\alpha\alpha\alpha}=U,U_{\alpha\beta\alpha\beta}=U-2J_{{\rm
H}},U_{\alpha\beta\beta\alpha}=U_{\alpha\alpha\beta\beta}=J_{{\rm
H}}(\alpha\neq\beta)$. Throughout the present study, $J_{{\rm H}}=0.15U$ is
taken which is motivated by a first-principles estimate for a typical $5d$
compound Na2 IrO3 ($U=2.72$ eV, $J_{\mathrm{H}}=0.23$ eV)[27].

The SOC term is written as

$\displaystyle\mathcal{H}_{{\rm
SOC}}=\lambda\sum_{\alpha,\beta,\sigma,\sigma^{\prime}}\expectationvalue{\alpha|\hat{\bm{l}}|\beta}\cdot\expectationvalue{\sigma|\hat{\bm{s}}|\sigma^{\prime}}c_{\alpha\sigma}^{\dagger}c_{\beta\sigma^{\prime}}$
(3)

where $\lambda$ is the SOC strength, and $\hat{\bm{l}}$ and $\hat{\bm{s}}$ are
the angular and spin momenta, respectively. By choosing the index on the order
of $xy,yz,zx$, the matrix elements of $\hat{\bm{l}}$ in the $t_{2{\rm g}}$
system and $\hat{\bm{s}}$ are given by

$\displaystyle l_{x}=\begin{pmatrix}0&0&-i\\\ 0&0&0\\\ i&0&0\\\
\end{pmatrix},$ (4) $\displaystyle l_{y}=\begin{pmatrix}0&i&0\\\ -i&0&0\\\
0&0&0\\\ \end{pmatrix},$ (5) $\displaystyle l_{z}=\begin{pmatrix}0&0&0\\\
0&0&i\\\ 0&-i&0\\\ \end{pmatrix},$ (6) $\displaystyle
s_{x}=\frac{1}{2}\begin{pmatrix}0&1\\\ 1&0\\\ \end{pmatrix},$ (7)
$\displaystyle s_{y}=\frac{1}{2}\begin{pmatrix}0&-i\\\ i&0\\\ \end{pmatrix},$
(8) $\displaystyle s_{z}=\frac{1}{2}\begin{pmatrix}1&0\\\ 0&-1\\\
\end{pmatrix}.$ (9)

Under the SOC, the effective total angular momentum $\hat{\bm{J}}_{{\rm
eff}}=-\hat{\bm{L}}+\hat{\bm{S}}$ is a good quantum number for the local
Hamiltonian. Here, $\hat{\bm{L}}$ ($\hat{\bm{S}}$) is the angular (spin)
momentum in the spin-orbital space and is expressed as

$\displaystyle\hat{\bm{L}}=\sum_{\alpha,\beta,\sigma}\expectationvalue{\alpha|\hat{\bm{l}}|\beta}c_{\alpha\sigma}^{\dagger}c_{\beta\sigma},$
(10)
$\displaystyle\hat{\bm{S}}=\sum_{\alpha,\sigma,\sigma^{\prime}}\expectationvalue{\sigma|\hat{\bm{s}}|\sigma^{\prime}}c_{\alpha\sigma}^{\dagger}c_{\alpha\sigma^{\prime}}.$
(11)

The single-particle eigenstates of $\hat{\bm{J}}_{\rm eff}^{2}$ and
$\hat{J}_{\rm eff}^{z}$ forms a complete local $j_{\mathrm{eff}}$ basis set as
$|j_{\mathrm{eff}},j_{\mathrm{eff}}^{z}\rangle=(|\tfrac{3}{2},\tfrac{3}{2}\rangle,|\tfrac{3}{2},\tfrac{1}{2}\rangle,|\tfrac{3}{2},-\tfrac{1}{2}\rangle,|\tfrac{3}{2},-\tfrac{3}{2}\rangle,|\tfrac{1}{2},\tfrac{1}{2}\rangle,|\tfrac{1}{2},-\tfrac{1}{2}\rangle)$.
The basis transformation matrix from the $t_{2{\rm g}}$ basis (by choosing the
index on the order of
$xy\uparrow,xy\downarrow,yz\uparrow,yz\downarrow,zx\uparrow,zx\downarrow$) to
the $j_{\mathrm{eff}}$ basis is given by

$\displaystyle
V=\begin{pmatrix}\frac{1}{\sqrt{3}}&0&0&\frac{2}{\sqrt{6}}&0&0\\\
0&-\frac{1}{\sqrt{3}}&0&0&\frac{2}{\sqrt{6}}&0\\\
0&\frac{1}{\sqrt{3}}&\frac{1}{\sqrt{2}}&0&\frac{1}{\sqrt{6}}&0\\\
\frac{1}{\sqrt{3}}&0&0&-\frac{1}{\sqrt{6}}&0&\frac{1}{\sqrt{2}}\\\
0&-\frac{i}{\sqrt{3}}&\frac{i}{\sqrt{2}}&0&-\frac{i}{\sqrt{6}}&0\\\
\frac{i}{\sqrt{3}}&0&0&-\frac{i}{\sqrt{6}}&0&-\frac{i}{\sqrt{2}}\\\
\end{pmatrix},$ (12)

where each column denotes the expansion coefficients of the corresponding
$j_{\mathrm{eff}}$ basis function on the $t_{\mathrm{2g}}$ basis.

##  III Method

In this study, a model (1) using the unrestricted Hartree-Fock approximation
at zero temperature ($T$) is solved by considering all possible local mean
fields. To investigate effects of local electronic correlations, only the
uniform solutions are considered. In this approximation, the interaction term
is decoupled as

$\displaystyle\mathcal{H}_{\mathrm{int}}=\frac{1}{2}\sum_{pqrs}U_{pqrs}c^{\dagger}_{p}c_{q}^{\dagger}c_{s}c_{r}$
$\displaystyle\hskip
22.76219pt\simeq\sum_{pqrs}(U_{pqrs}-U_{pqsr})D_{pr}c^{\dagger}_{q}c_{s}$
$\displaystyle\hskip
22.76219pt-\frac{1}{2}\sum_{pqrs}(U_{pqrs}-U_{pqsr})D_{pr}D_{qs},$ (13)

where, composite indices ($p$, $q$, $r$, $s$) for spin orbitals, and a
$6\times 6$ density matrix
$D_{pr}\equiv\expectationvalue{c^{\dagger}_{p}c_{r}}$ were introduced. Note
that the site index $i$ in Eq. (13) is dropped. Self-consistent calculations
using different initial density matrices are performed to obtain the lowest-
energy states.

As will be seen later, the ground state of this model may be degenerated
because the spin and orbital moments can be completely decoupled without SOC.
Thus a very small SOC $\lambda=10^{-5}$ is introduced to lift the degeneracy.
Furthermore, the spin moment is always aligned along the $z$ axis.

##  IV Complete set of order parameters for $t_{\mathrm{2g}}$ systems

Multipoles | label |  rank |  time-reversal |  $j$-parity
---|---|---|---|---
$j_{{\rm eff}}$-diagonal |  $N^{\rm 3/2,even}$ |  $a$ | 0 (monopole) |  $+$ (electric) | +
$M_{\mu}^{\rm 3/2,odd}$ |  $a$ | 1 (dipole) |  $-$ (magnetic) | +
$Q_{\lambda}^{\rm 3/2,even}$ |  $a$ | 2 (quadrupole) |  $+$ | +
$T_{\xi}^{\rm 3/2,odd}$ |  $a$ | 3 (octupole) |  $-$ | +
$N^{\rm 1/2,even}$ |  $b$ | 0 |  $+$ | +
$M_{\mu}^{\rm 1/2,odd}$ |  $b$ | 1 |  $-$ | +
$j_{{\rm eff}}$-offdiagonal |  $M_{\mu}^{\rm offd,odd}$ |  $c$ | 1 |  $-$ | $-$
$Q_{\lambda}^{\rm offd,even}$ |  $c$ | 2 |  $+$ | $-$
$M_{\mu}^{\rm offd,even}$ |  $d$ | 1 |  $+$ | $-$
$Q_{\lambda}^{\rm offd,odd}$ |  $d$ | 2 |  $-$ | $-$
Table 1: List of the complete multipole operators in three-orbital model with
spin-orbit coupling. The suffices are $\mu=x,y,z$,
$\lambda=xy,yz,zx,3z^{2}-r^{2},x^{2}-y^{2}$,
$\xi=xyz,x(5x^{2}-3r^{2}),y(5y^{2}-3r^{2}),z(5z^{2}-3r^{2}),x(y^{2}-z^{2}),y(z^{2}-x^{2}),z(x^{2}-y^{2})$
with $r^{2}=x^{2}+y^{2}+z^{2}$.

The symmetry-breaking information is fully encoded in the single-particle
density matrix. Multipole is a useful tool for describing the entanglement of
spin and orbit degrees of freedom. The multipole moment is described by a
polynomial form of the effective total angular momentum $\hat{\bm{J}}_{{\rm
eff}}=-\hat{\bm{L}}+\hat{\bm{S}}$. Such multipoles are introduced by means of
crystallographic point groups.

The multipole representation used so far for the spin-orbit-coupled three-
orbital Hubbard model [16, 17] are reviewed first. In a $t_{\mathrm{2g}}$
orbital system, a conventional multipole is constructed as a polynomial of
$\hat{\bm{J}}_{\rm eff}$, and is described by 16 order parameters. The density
matrix $D$ is expanded as follows:

$\displaystyle D=\sum^{16}_{\xi=1}C_{\xi}O_{\xi},$ (14)

where $C$ is the weight of $O$ and $\xi$ is the index of the basis set. The
breakdown is one order parameter(OP) for the electric monopole $N$, three OPs
for the magnetic dipole $M$, five OPs for the electric quadrupole $Q$, and
seven OPs for the magnetic octupole $T$. The higher-rank tensors are trivially
zero. Operator $O$ is imposed the orthonormality

$\displaystyle{\rm
Tr}\left[O_{\xi}O_{\eta}^{\dagger}\right]=\delta_{\xi\eta}.$ (15)

The weight of the order parameter $C$ can be computed as

$\displaystyle C_{\xi}={\rm Tr}\left[DO_{\xi}^{\dagger}\right].$ (16)

However, this definition does not produce a complete basis set because it
requires $6\times 6=36$ order parameters.

To construct a complete set, the density matrix is decomposed into 36 one-
particle tensor operators $O$ using the local projection operator. Note that
multipoles for this complete set are not constructed by a polynomial of
$\hat{\bm{J}}_{{\rm eff}}$, but $\bm{K}=\hat{\bm{L}}+\hat{\bm{S}}$ (see
Appendix A). which makes it possible to account for a component perpendicular
to $\hat{\bm{J}}_{\rm eff}$. These tensor operators are categorized into four
blocks: diagonal components $\ket{j_{{\rm
eff}}=3/2}\otimes\ket{j_{\mathrm{eff}}=3/2}$, $\ket{1/2}\otimes\ket{1/2}$,
off-diagonal components $\ket{3/2}\otimes\ket{1/2}$, and
$\ket{1/2}\otimes\ket{3/2}$. They are labelled as $a,b,c$, and $d$ in the
definition order. A schematic illustration of this process is presented in
Fig. 1. Components $c$ and $d$ show that there is hybridization between the
$j_{\mathrm{eff}}=3/2$ and $1/2$ orbitals. Note that the OPs in blocks $c$ and
$d$ are given by a linear combination of the off-diagonal components to make
it Hermitian. A time-reversal even and odd is imposed on the $c$ and $d$
components, respectively.

Figure 1: Schematic illustration of multipolar order parameters. The density
matrix is decomposed into four blocks in terms of hybridization between the
$j_{\mathrm{eff}}$ orbitals. Block $a$ has one monopole, three dipoles, five
quadrupoles, and seven octupoles. Block $b$ contains one monopole and three
dipoles. Blocks $c$ and $d$ have three dipoles and five quadrupoles. For more
details, see the text.

The complete classification of multipoles are derived in terms of the quantum
numbers of rank, time reversal, and $j$-parity. The results are summarized in
Table 1 (see Appendix A for details of the classification). The advantage of
the present classification is that one can identify diagonal and off-diagonal
orders in terms of $j_{\mathrm{eff}}$ manifolds separately using the four
labels ($a$, $b$, $c$, $d$). The density matrix $D$ is decomposed into two
monopoles $N$, 12 dipoles $M$, 15 quadrupoles $Q$, and seven Octupoles $T$.
Here, the $j$-parity distinguishes the diagonal and off-diagonal components
and is similar to $sp$-hybridized systems where the parity operator
distinguishes the angular momentum $\ell=0$ or $\ell=1$. In previous studies
[16, 17], they found quantum phases with off-diagonal orders, which was called
“excitonic phases”. They identified the existence of these phases by the off-
diagonal elements in the hybridization functions on the $j_{\mathrm{eff}}$
basis. The present scheme allows to directly detect and classify such phases
even further.

Figure 2: $U$–$\lambda$ phase diagram for $n=4$. PM, PI, MM, MI, and QCP
represent the paramagnetic metal phase, paramagnetic insulator phase, magnetic
metal phase, magnetic insulator phase, and quantum critical point,
respectively. The numbers (I) and (II) of the MM phase denote the absence and
presence of the spontaneous breaking of the orbital degeneracy, respectively.
The dashed line denotes a metal-insulator transition. Figure 3: Conventional
OPs computed as a function of $U$ at $n=4$. The background colors correspond
to those of the phases shown in Figs. 2).  Figure 4: Angular momentum
$\hat{\Braket{\bm{L}}}_{x},\hat{\Braket{\bm{L}}}_{y}$ and
$\hat{\Braket{\bm{L}}}_{z}$ computed as a function of $U$ at $n=4$ and
$\lambda=0$. The background colors correspond to those of the phases shown in
Figs. 2.

##  V Results for $n=4$

In this section, results for $n=4$ are discussed. Figure 2 shows a
$U$–$\lambda$ phase diagram computed for $n=4$, which involves four different
phases: a paramagnetic metallic (PM) phase, a magnetic metallic (MM) phase, a
paramagnetic insulating (PI) phase, and a magnetic insulating (MI) phase. The
labels (I) and (II) denote the absence and presence of non-zero higher-rank
multipoles $Q$ and $T$, respectively, as shown in Fig. 3. The detailed
analysis of multipole order parameters will be given in the next subsections.
In the PM and PI phases, no spontaneous symmetry breaking takes place. At a
larger $U$, the time-reversal symmetry is broken in the MM and MI phases. In
the phase diagram, there is a dome of MM-II at moderate values of $U$, while
the MI is always stable in the strong $U$ limit. In these magnetic phases, the
electric quadrupole is also active.

Note that the PI phase corresponds to the relativistic band insulator at
$n=4$, where the $j_{\mathrm{eff}}=3/2$ and $j_{\mathrm{eff}}=1/2$ bands are
completely separated in energy by the strong SOC. The critical value
$\lambda_{\mathrm{c}}$ of the metal-insulator transition between the PM and PI
is $\lambda_{\mathrm{c}}=8/3$ at $U=0$. As seen in Fig. 2,
$\lambda_{\mathrm{c}}$ decreases as $U$ increases until it reaches the
magnetic transition line. On the other hand, the boundary between PI and MI-II
is determined by the competing $\lambda$ and $J_{\mathrm{H}}$, and the slope
of the boundary is roughly proportional to $J_{\mathrm{H}}$ (not shown).

Figure 5: $U$ dependence of the computed DOS projected onto $l_{z}$=0, $\pm 1$
orbitals computed at $n=4$ and $\lambda=10^{-5}$. The Fermi energy was
$\omega=0$.  Figure 6: $\lambda$ dependence of the computed partial DOS for
$n=4$ and $U=6$. The DOS is projected onto (a) the $j_{\mathrm{eff}}$=1/2 and
3/2 orbitals and (b) $l_{z}$=0, $\pm 1$ orbitals. The Fermi energy is located
at $\omega=0$.

###  V.1 Results for $\lambda=0$

First, the results of conventional OPs computed for $\lambda=0$ are discussed.
In this study, the result of $\lambda=10^{-5}$ is treated as $\lambda=0$
because the spin and orbital moments are completely decoupled at $\lambda=0$
in the MI phase. As seen in Fig. 3, all the conventional OPs are zero up to
$U\simeq 2.6$. Since the magnetic moment is aligned along the $z$ axis, only
$M_{z},Q_{u}(3z^{2}-r^{2})$ and $T^{\alpha}_{z}(z(5z^{2}-3r^{2}))$ are non-
zero. A transition point separates the PM phase and the MM-I phase around
$U_{{\rm c}1}\simeq 2.6$. This transition is characterized by the spontaneous
breaking of the time-reversal symmetry through the emergence of $M_{z}\neq 0$.
This is followed by a subsequent first-order transition at $U_{{\rm c}2}\simeq
5.86$, which is characterized by the emergence of $Q_{u}$ and
$T_{z}^{\alpha}$. This corresponds to a spontaneous breaking of the orbital
degeneracy. As shown in Fig. 4, the angular momentum $L_{z}$ becomes finite
for $U\geq U_{{\rm c}2}$. With a further increase in $U$, the metal-insulator
transition occurs, and $L_{z}$ reaches one at $U_{{\rm c}3}\simeq 7.2$. The
symmetry of the system does not change across this transition point.

Figure 5 shows the partial density of states (DOS) projected onto the angular
momentum $l_{z}$=$0,\pm 1$ orbitals for typical values of $U$ in the PM, MM-I,
and MM-II phases. Note that the mean fields are diagonal in the $l_{z}$ basis
for these parameters. The six basis vectors are labelled by the spin
$(\uparrow,\downarrow)$ and angular momentum $l_{z}$=$0,\pm 1$. A magnetic
transition happens at $U_{c\rm{1}}\sim 2.6$. At $U\simeq 3$, in the MM-I
phase, a Lifshitz transition occurs, and then the spin becomes fully polarized
as $\hat{\Braket{\bm{S}}}_{z}=1$. As seen in Fig. 5, the three up-spin
orbitals are completely filled, leaving one electron in the down-spin
orbitals. The partially filled three down-spin orbitals are still degenerate,
and thus $\hat{\Braket{\bm{L}}}_{z}=0$. After the spontaneous breaking of the
orbital degeneracy ($U\geq U_{{\rm c}2}$), the partial $l_{z}$=+1$\uparrow$
DOS is separated from the partial $l_{z}$=$0$$\uparrow$ and partial
$l_{z}$=$-1$$\uparrow$ DOS. At $U=U_{c3}$, a finite gap finally opens at the
Fermi energy, and the system becomes a magnetic insulator with
$\hat{\Braket{\bm{S}}}_{z}=1$ and $\hat{\Braket{\bm{L}}}_{z}=1$. Note that the
spin and orbital moments are coupled in an antiparallel manner owing to the
small SOC. In the fully spin-polarized situation, the self-consistent equation
is simplified. See Appendix B for further details.

Figure 7: Complete OPs computed as a function of $U$ at $n=4$. The
superscripts $a,b$, and $c$ denote $j_{\mathrm{eff}}=3/2$ components,
$j_{\mathrm{eff}}=1/2$ components, and $j_{\mathrm{eff}}=3/2$–$1/2$ entangled
components, respectively. The background colors correspond to those of the
phases shown in Figs. 2).

###  V.2 Results for $\lambda>0$

Once $\lambda$ is turned on, $Q_{u}$ and $T_{z}^{\alpha}$ become coupled to
$M_{z}$. Thus, the MM-II phase is characterized by
$M_{z},Q_{u},T_{z}^{\alpha}\neq 0$. As shown in the inset of Fig. 2, the
first-order transition terminates at a quantum critical endpoint at
$\lambda\sim 0.003$ and $U\sim 5.68$. On the other hand, the critical value
$U$ of the magnetic transitions increases as $\lambda$ increases, and the SOC
suppresses magnetization. Figure 6 shows the $\lambda$ dependence of the
partial DOS projected onto the $j_{\mathrm{eff}}$ basis and the $l_{z}$ basis
at $U=6$. In the weak SOC limit (MM-II), the partial $l_{z}$ DOS is not
distorted, so the $l_{z}$ scheme is a better representation. As $\lambda$ is
increased, the $j_{\mathrm{eff}}=1/2$ contributions become dominant for the
unoccupied DOS, whereas those of $j_{\mathrm{eff}}=3/2$ become dominant for
the occupied DOS. In this regime, the $j_{\mathrm{eff}}$ representation is
better. Figure 7 shows the complete OPs computed as a function of $U$ at
$n=4$. $M_{z}^{a},M_{z}^{b},M_{z}^{c},Q_{u}^{a},Q_{u}^{c}$ and $T^{\alpha
a}_{z}$ are non-zero. The difference between conventional (Fig. 3), and the
complete OPs are clearly seen in the magnetic insulator phase. With increasing
$U$, one of the old OPs, $M_{z}$, decreases toward zero and vanishes
completely for $U\geq U_{c3}$. In terms of the complete OPs, however, the off-
diagonal component of the magnetic dipole $M^{c}_{z}$ remains the most
dominant.

##  VI Results for general filling

In this section, the filling $n$ dependence of the ground-state phase diagrams
are discussed.

###  VI.1 $n$–$U$ phase diagram at $\lambda=0$

Figure 8(a) shows the $n$–$U$ ground-state phase diagram computed for
$\lambda=0$. The PM exits only for a small $U$ at $\lambda=0$. With an
increment in $U$, the ground state turns into the magnetic metallic phase
(MM-I), where only the magnetic dipole $M_{z}$ is finite. With a further
increase in $U$, the MM-I phase turns into the MM-II phase through a first-
order transition. In the MM-II phase, the symmetry of the orbital is broken,
and the angular momentum $L$ and higher-order multipole order parameters are
finite as $L_{z}\neq 0$, $M_{z}\neq 0$, $Q_{u}\neq 0$, $T_{z}^{\alpha}\neq 0$
(see Figs. 3, and 4). At $n=3$, symmetry breaking in the orbital sector does
not occur because $\hat{\Braket{\bm{L}}}_{z}=0$ and
$\hat{\Braket{\bm{S}}}_{z}=3/2$, as illustrated in Fig. 9.

Furthermore, a new distinct phase (MM-III) emerges only for $1<n<2$ and
$4<n<5$. The transition between MM-II and MM-III is second-order. The three MM
phases can be distinguished by different degeneracies of the DOS projected
onto the $l_{z}$ orbitals for the major spin. Figure 10 plots typical data of
DOS at $U=7.5$. In the MM-I phase, all three orbitals are degenerate for each
spin. In the MM-II phase, only two of the three orbitals remain degenerate. In
contrast, in the MM-III phase, all three orbitals become non-degenerate. The
MM-I, MM-II, and MM-III phases meet at $U\sim 5.5$, $n=1.5$, and $4.5$.
Moreover, the first-order transitions between MM-II and MM-III become
continuous only at these critical points.

The insulating states can be stable only if $n$ is an integer. Except for
$n=3$, the metal-insulator transition occurs at $U_{{\rm c}3}\sim 7.2$. At
$n=3$, the metal-insulator transition occurs at a smaller $U\sim 3$ because a
high-spin configuration is stable. Figure 9 illustrates the electron
configurations projected onto the $l_{z}$ orbitals in the insulting phases at
$n=1,2,\cdots,5$. At $n=3$, 4, and 5, the spin-up orbitals are fully occupied
for a large $U$. Thus, the system can be effectively regarded as a spinless
system of spin-down orbitals. This gives rise to an interesting emergent
symmetry between $n=4$ and $n=5$: A particle-hole transformation for the spin-
down orbitals connects these two states. Consequently, the MM-II phase is
symmetric with respect to $n=4.5$, as shown in Fig. 8(a). This symmetry
originates from the fact that the spin-up orbitals are fully occupied, which
is due to the ignorance of quantum spin fluctuations in the zero-$T$ Hartree-
Fock calculations. See Appendix B for a more detailed discussion. In addition,
the phase is symmetric with respect to $n=3$ owing to the particle-hole
symmetry for the spinful system. Furthermore, the quantum critical end point
is located at $U_{{\rm c}2}\sim 6.8$ and $\lambda\sim 0.003$ for $n=1,2,4$,
and $5$.

Figure 8: The $n$–$U$ phase diagram for (a) $\lambda=0$, (b) $\lambda=0.003$,
(c) $\lambda=0.5$ and (d) $\lambda=1.0$. PM and MM represent a paramagnetic
metal phase and magnetic metal phase, respectively. The numbers I and II(III)
indicate the absence and appearance of the higher-rank multipoles $Q$ or $T$,
respectively. The solid and broken lines indicate the first-order and second-
order transitions, respectively. The blue lines on $n=1,2,3,4$, and $5$
indicate the insulating phases. The broken circles denote multiple critical
points where the first-order transition lines merge (see text).  Figure 9:
Schematic illustration of electronic configurations projected onto $l_{z}$=0,
$\pm 1$ orbitals at $\lambda=0$ and $U=7.5$. Note that the mean fields and
density matrices are diagonal in this basis.

###  VI.2 $n$–$U$ phase diagrams at $\lambda>0$

Figures 8(b), 8(c), and 8(d) show the $n$–$U$ ground-state phase diagrams for
$\lambda=0.003,0.5,1.0$, respectively. As seen in Fig. 8(b), even for an
infinitesimal value of $\lambda$, the MM-I phase changes into the MM-II phase
because $Q$ and $T$ emerge owing to the coupling between $M$ and $Q$/$T$
through $\lambda$ (refer to the discussion in Sec. V and Fig. 2).

For finite $\lambda$, the MM-III phase survives up to $\lambda\sim 0.0001$
(not shown). With further increasing $\lambda$, the MM-I($b$) phase and the
MM-II($a$) phase emerge. The MM-I($b$) phase is located around $n=5$. In this
phase, the $b$ component $(M_{z}^{b})$ is finite, and the other components are
zero. Comparing the phase diagrams at $\lambda=0$ and $0.003$ in Fig. 8, one
can see that the MM-II phase around $n=5$ is replaced by MM-I($b$). With an
increment in $\lambda$, the MM-I($b$) phase broadens around $n=5$. In this
phase, the two $j_{\mathrm{eff}}$ = 3/2 orbitals are completely filled, and
the $j_{\mathrm{eff}}$ = 1/2 orbital remains partially filled. As a
consequence, the hybridization between $j_{\mathrm{eff}}$ = 3/2 and 1/2
vanishes.

In contrast, the MM-II($a$) phase was located for $n\leq 1$. In this phase,
the $a$ components $(M_{z}^{a},Q_{u}^{a}$ and $T^{\alpha a}_{z})$ are finite,
and the other components are zero. The MM-II($a$) phase appears when the
energy gap between $j_{\mathrm{eff}}$ = 3/2 and 1/2 is large. The
$j_{\mathrm{eff}}$ = 1/2 orbital is empty and its contribution is negligible.

Figure 10: $n$ dependence of the computed DOS projected onto $l_{z}$=0, $\pm
1$ orbitals for three different phases ($U=7.5$, $\lambda=0$). The Fermi
energy is located at $\omega=0$.

###  VI.3 $U$–$\lambda$ phase diagrams

Figure 11 shows the ground-state $U-\lambda$ phase diagrams for $n=1,2,3$ and
$5$ (The phase diagram for $n=4$ was discussed in Sec. V).

At $n=5$, the phase boundaries between PM, MM-I$(b)$, and MI-I$(b)$ are
vertical for $\lambda\gtrsim 0.7$. In this regime, the $j_{\mathrm{eff}}$ =
3/2 orbitals are completely filled and the $j_{\mathrm{eff}}$ = 1/2 orbital is
half-filled owing to the gap between these two manifolds induced by $\lambda$.
Thus, these transitions can be regarded as phase transitions in the effective
single-orbital model of $j_{\mathrm{eff}}$ = 1/2. This explains why the
critical values of $U$ do not depend on $\lambda$.

Figure 11: Phase diagrams computed for $n=1,2,3,5$. The data are the same as
those shown in Fig. 8). See the main text and captions of Fig. 2, to define
the phases (labels).

For $n=3$, the boundary between MM-II and MI-II runs linearly to the upper
right. This may be due to the competition between $\lambda$ and local
interactions ($U$, $J_{\mathrm{H}}$).

The phase diagrams for $n=1$ and $n=2$ are similar in nature. This is because
the electrons partially occupy $j_{\rm eff}=3/2$, and the fillings are less
than half in both cases. The partially filled orbitals always have magnetic
moments regardless of the value of $\lambda$, which is magnetically ordered by
interaction effects. In contrast, at $n=4$, a sufficiently large $\lambda$
suppresses the magnetic order because active moments are absent.

Next, the intensity maps of the OPs computed for $n=1,2,3,4$, and $5$ are
shown in Fig. 12. The absolute values of the computed OPs are plotted. The
phase transition lines are marked with white lines. $j_{\mathrm{eff}}$ = 1/2
diagonal components denoted by $a$ are dominant for $n=5$. For $n=3$ and $4$,
$a$ and $c$ components ($j_{\mathrm{eff}}$ = 3/2–1/2 entangled) are enhanced
depending on the balance between $U$ and $\lambda$. As shown in Fig. 8, a
first-order transition with spontaneous symmetry breaking is observed, except
for $n=3$. In addition, as shown in the case of $n=4$ (Sec. V), higher-order
multipoles $Q$ and $T$ are enhanced when crossing the quantum critical line or
its crossover line in the direction in which $U$ grows.

Figure 12: Intensity maps for $n=1,2,3,4,5$. The solid and broken lines
indicate the magnetic phase transitions and metal-insulator transitions,
respectively. See Table 1 for the definition of the order parameters. The
superscripts $a,b$, and $c$ denote $j_{\mathrm{eff}}=3/2$ components,
$j_{\mathrm{eff}}=1/2$ components, and $j_{\mathrm{eff}}=3/2$–$1/2$ entangled
components, respectively.

The characteristic features of the multipole order parameters at each filling
are discussed now. At $n=1$, symmetry breaking of the orbital enhances
$Q_{u}^{a}$. Because the electron is mainly stored in the $j_{\mathrm{eff}}$ =
3/2 orbitals regardless of $U$ or $\lambda$, the $a$ component is dominant.
When $U$ is sufficiently large (i.e., $U\geq 7$), all the $b$ and $c$
components become inactive and only the $a$ components remain active. At the
same time, the metal-insulator transition occurs, and the MI-II$(a)$ phase
appears. A valley of OPs, such as $T_{z}^{\alpha a}$ at $U=5$, indicates a
sign change similar to $Q_{u}$, as shown in Fig. 3).

At $n=2$, $Q_{u}^{c}$ and $T^{\alpha a}_{z}$ are notable. The former is
enhanced by the symmetry breaking in the orbital sector. On the other hand,
the latter is amplified by the SOC because $j_{\mathrm{eff}}$ = 3/2 and 1/2
are split, and the $a$ component is activated. $T^{\alpha a}_{z}$ has a peak
at $n=2$, probably because it belongs to the same irreducible representation
$\Gamma_{4u}$ as $M_{z}$.

At $n=3$, either the high-spin state or the low-spin state is stable,
depending on $U$ and $\lambda$. For $U\gg\lambda$, the high-spin state is
stable. In this state, the angular momentum $L$ vanishes, and thus $Q$ and $T$
are zero. Owing to this circumstance, the magnetic insulating state without
$Q$ and $T$ (MI-I) is realized at $n=3$. For $U\ll\lambda$, a low-spin state
with three electrons in the $j_{\mathrm{eff}}=3/2$ orbitals becomes stable. In
this state, both, the spin and orbital degrees of freedom remain; thus,
higher-order multipoles can be activated.

At $n=4$, as already discussed in Sec. V, the $c$ component of the magnetic
dipole $M_{z}^{c}$ is dominant. In the large $\lambda$ region, the
paramagnetic insulating (PI) state, which is induced by the strong SOC,
appears for $\lambda\geq 1.3$.

At $n=5$ and large $\lambda$, the two $j_{\mathrm{eff}}$ = 3/2 orbitals are
almost filled, leaving the $j_{\mathrm{eff}}$ = 1/2 orbital half-filled.
Consequently, all the $a$ components are inactive, while the $b$ component
$(M_{z}^{b})$ is enhanced. At small $U$ and $\lambda$, the $a$ and $c$
components are active because of the hybridization between $j_{\mathrm{eff}}$
= 1/2 and 3/2. However, with further increments of $U$ or $\lambda$, the $a$
and $c$ components become inactive and the MM-I($b$) and MI-I($b$) phases
emerge.

##  VII SUMMARY

In conclusion, the $t_{2{\rm g}}$ Hubbard model using the unrestricted
Hartree-Fock approximation with all possible (particle-conserving) local
symmetry breaking patterns has been studied.

The main results of this work are twofold. First, a complete set of multipole
order parameters are constructed that can identify entanglements between
$j_{\mathrm{eff}}$= 1/2 and 3/2 manifolds. Second, through extensive Hartree-
Fock calculations, the ground-state phase diagrams in the parameter space of
the onsite Coulomb repulsion $U$, the strength of the spin-orbit coupling
(SOC) $\lambda$, and filling $n$ are systematically investigated. It has been
determined that this model hosts many non-trivial quantum phases with
multipole ordering as well as peculiar phase structures, such as multi-
critical points. Furthermore, the intensity maps of the multipolar order
parameters for the phase diagrams have been computed.

The results of the present study clearly show that the simple
$t_{\mathrm{2g}}$ model with a semi-circular density of states can host a
variety of phenomena. An interesting future study may be the analyses of the
model using a more elaborate method, such as the dynamical mean-field (DMFT)
theory, although DMFT calculations suffer from a sign problem. Another future
direction may be comparisons with $4d$ and $5d$ materials by considering
realistic band structures.

###### Acknowledgements.

N.C. and H.S. were supported by JSPS KAKENHI Grants No. 18H01158 and No.
21H01003, and S.H. by No. 19H01842 and No. 21K03459.

## References

  * Shimura  _et al._ [1995] T. Shimura, Y. Inaguma, T. Nakamura, M. Itoh, and Y. Morii, Phys. Rev. B 52, 9143 (1995).
  * Cao  _et al._ [1998] G. Cao, J. Bolivar, S. McCall, J. E. Crow, and R. P. Guertin, Phys. Rev. B 57, R11039 (1998).
  * Kim  _et al._ [2008] B. J. Kim, H. Jin, S. J. Moon, J.-Y. Kim, B.-G. Park, C. S. Leem, J. Yu, T. W. Noh, C. Kim, S.-J. Oh, J.-H. Park, V. Durairaj, G. Cao, and E. Rotenberg, Phys. Rev. Lett. 101, 076402 (2008).
  * Harter  _et al._ [2017] J. W. Harter, Z. Y. Zhao, J.-Q. Yan, D. G. Mandrus, and D. Hsieh, Science 356, 295 (2017).
  * Hanawa  _et al._ [2001] M. Hanawa, Y. Muraoka, T. Tayama, T. Sakakibara, J. Yamaura, and Z. Hiroi, Phys. Rev. Lett. 87, 187001 (2001).
  * C. Kobayashi  _et al._ [2011] T. C. Kobayashi, Y. Irie, J.-i. Yamaura, Z. Hiroi, and K. Murata, Journal of the Physical Society of Japan 80, 023715 (2011).
  * Sleight  _et al._ [1974] A. Sleight, J. Gillson, J. Weiher, and W. Bindloss, Solid State Communications 14, 357 (1974).
  * Mandrus  _et al._ [2001] D. Mandrus, J. R. Thompson, R. Gaal, L. Forro, J. C. Bryan, B. C. Chakoumakos, L. M. Woods, B. C. Sales, R. S. Fishman, and V. Keppens, Phys. Rev. B 63, 195104 (2001).
  * Yamaura  _et al._ [2012] J. Yamaura, K. Ohgushi, H. Ohsumi, T. Hasegawa, I. Yamauchi, K. Sugimoto, S. Takeshita, A. Tokuda, M. Takata, M. Udagawa, M. Takigawa, H. Harima, T. Arima, and Z. Hiroi, Phys. Rev. Lett. 108, 247205 (2012).
  * Wang  _et al._ [2017] Y. Wang, H. Weng, L. Fu, and X. Dai, Phys. Rev. Lett. 119, 187203 (2017).
  * Chun  _et al._ [2018] S. H. Chun, B. Yuan, D. Casa, J. Kim, C.-Y. Kim, Z. Tian, Y. Qiu, S. Nakatsuji, and Y.-J. Kim, Phys. Rev. Lett. 120, 177203 (2018).
  * Fuchs  _et al._ [2018] S. Fuchs, T. Dey, G. Aslan-Cansever, A. Maljuk, S. Wurmehl, B. Büchner, and V. Kataev, Phys. Rev. Lett. 120, 237204 (2018).
  * Singh and Gegenwart [2010] Y. Singh and P. Gegenwart, Phys. Rev. B 82, 064412 (2010).
  * Liu  _et al._ [2011] X. Liu, T. Berlijn, W.-G. Yin, W. Ku, A. Tsvelik, Y.-J. Kim, H. Gretarsson, Y. Singh, P. Gegenwart, and J. P. Hill, Phys. Rev. B 83, 220403 (2011).
  * Sato  _et al._ [2015] T. Sato, T. Shirakawa, and S. Yunoki, Phys. Rev. B 91, 125122 (2015).
  * Kim  _et al._ [2017] A. J. Kim, H. O. Jeschke, P. Werner, and R. Valentí, Phys. Rev. Lett. 118, 086401 (2017).
  * Sato  _et al._ [2019] T. Sato, T. Shirakawa, and S. Yunoki, Phys. Rev. B 99, 075117 (2019).
  * Anisimov  _et al._ [1997] V. I. Anisimov, A. I. Poteryaev, M. A. Korotin, A. O. Anokhin, and G. Kotliar, Journal of Physics: Condensed Matter 9, 7359–7367 (1997).
  * Werner and Millis [2006] P. Werner and A. J. Millis, Phys. Rev. B 74, 155107 (2006).
  * Werner  _et al._ [2006] P. Werner, A. Comanac, L. de’ Medici, M. Troyer, and A. J. Millis, Phys. Rev. Lett. 97, 076405 (2006).
  * Kim  _et al._ [2020] A. J. Kim, P. Werner, and R. Valentí, Phys. Rev. B 101, 045108 (2020).
  * Shinaoka  _et al._ [2019] H. Shinaoka, Y. Motome, T. Miyake, S. Ishibashi, and P. Werner, Journal of Physics: Condensed Matter 31, 323001 (2019).
  * J. Ohkawa  _et al._ [1983] F. J. Ohkawa, H. Fukuyama, and K. Yosida, Journal of the Physical Society of Japan 52, 1701 (1983).
  * Shiina  _et al._ [1997] R. Shiina, H. Shiba, and P. Thalmeier, Journal of the Physical Society of Japan 66, 1741 (1997).
  * Kusunose [2008] H. Kusunose, Journal of the Physical Society of Japan 77, 064710 (2008).
  * Kusunose  _et al._ [2020] H. Kusunose, R. Oiwa, and S. Hayami, Journal of the Physical Society of Japan 89, 104704 (2020).
  * Yamaji  _et al._ [2014] Y. Yamaji, Y. Nomura, M. Kurita, R. Arita, and M. Imada, Phys. Rev. Lett. 113, 107201 (2014).
  * J. Ohkawa [1983] F. J. Ohkawa, Journal of the Physical Society of Japan 52, 3897 (1983).
  * Kuramoto  _et al._ [2009] Y. Kuramoto, H. Kusunose, and A. Kiss, Journal of the Physical Society of Japan 78, 072001 (2009).
  * Hayami  _et al._ [2018] S. Hayami, M. Yatsushiro, Y. Yanagi, and H. Kusunose, Phys. Rev. B 98, 165110 (2018).
  * Tamura  _et al._ [2020] S. Tamura, S. Iimura, and S. Hoshino, Phys. Rev. B 102, 024505 (2020).

##  Appendix A Complete basis set with spin-orbit coupling

For systems with spin-orbit coupling, the energy level is split into lower
$j_{\rm eff}=3/2$ and higher $j_{\rm eff}=1/2$ states. At the filling $n=4$,
the lower orbital is fully occupied, so that no degrees of freedom are left
within the $j_{\rm eff}=3/2$ multiplet. Then, we need to consider the excited
$j_{\rm eff}=1/2$ state, which can, in principle, mix through thermal and/or
interaction effects. Such degrees of freedom are called “excitonic” in the
sense that the transition matrix involves the process from the ground state to
the energetically excited states, and are naturally described through the
$j_{{\rm eff}}$-off-diagonal matrix elements. Here the terminology “$j_{{\rm
eff}}$-off-diagonal” is introduced, in a manner similar to the pairing
amplitude in superconductivity which is called the off-diagonal long-range
order. The purpose of this section is to construct a proper complete set of
operators to describe the order parameters for multi-orbital systems with
spin-orbit coupling.

Let us begin with the total angular momentum

$\displaystyle\hat{\bm{J}}_{{\rm eff}}$
$\displaystyle=-\hat{\bm{L}}+\hat{\bm{S}}$ (17)

which is a good quantum number for the local Hamiltonian. The eigenstates form
a complete local basis set $|j_{{\rm eff}},j_{{\rm
eff},z}\rangle=(|\tfrac{3}{2},\tfrac{3}{2}\rangle,|\tfrac{3}{2},\tfrac{1}{2}\rangle,|\tfrac{3}{2},-\tfrac{1}{2}\rangle,|\tfrac{3}{2},-\tfrac{3}{2}\rangle,|\tfrac{1}{2},\tfrac{1}{2}\rangle,|\tfrac{1}{2},-\tfrac{1}{2}\rangle)$.
For multi-orbital systems, the operators of these states are classified by
utilizing the concept of multipole expansion. The concept of multipoles had
been introduced originally for a description of the local degrees of freedom
of $f$ electrons in terms of the total angular momentum $J_{{\rm eff}}$ [28,
24, 25, 29, 30].

High-rank multipoles are needed that can be introduced as a polynomial form of
the $\hat{\bm{J}}_{{\rm eff}}$ operator [17]. However, the $j_{\rm eff}$-off-
diagonal components, which are the essential quantities for the order
parameters at $n=4$, cannot be described because the amplitude of the total
angular momentum $\hat{\bm{J}}^{2}_{{\rm eff}}$ is a conserved quantity. The
corresponding set of multipole operators is incomplete. Hence, one needs to
consider an operator that includes the transition between different $j_{\rm
eff}$ s. Namely, another angular momentum is defined:

$\displaystyle\bm{K}=\alpha\hat{\bm{L}}+\beta\hat{\bm{S}}$ (18)

where $\alpha$ and $\beta$ ($\in\mathbb{R}$) are constants, which in general
includes the “perpendicular” component with respect to $\hat{\bm{J}}_{{\rm
eff}}$. This operator also commutes with the local Hamiltonian, but does not
commute with $\hat{\bm{J}}_{{\rm eff}}$. Multipole operators can be
constructed based on the polynomial expressions of $\bm{K}$ as:

$\displaystyle N$ $\displaystyle=1$ (19) $\displaystyle\bm{M}$
$\displaystyle=\bm{K}$ (20) $\displaystyle Q_{xy}$
$\displaystyle=\overline{K_{x}K_{y}}$ (21) $\displaystyle Q_{yz}$
$\displaystyle=\overline{K_{y}K_{z}}$ (22) $\displaystyle Q_{zx}$
$\displaystyle=\overline{K_{z}K_{x}}$ (23) $\displaystyle Q_{3z^{2}-r^{2}}$
$\displaystyle=3K_{z}^{2}-\bm{K}^{2}\ \ (=Q_{u})$ (24) $\displaystyle
Q_{x^{2}-y^{2}}$ $\displaystyle=\overline{K_{x}^{2}-K_{y}^{2}}\ \ (=Q_{v})$
(25) $\displaystyle T_{xyz}$ $\displaystyle=\overline{K_{x}K_{y}K_{z}}$ (26)
$\displaystyle T_{x(5x^{2}-3r^{2})}$
$\displaystyle=\overline{K_{x}(5K_{x}^{2}-3\bm{K}^{2})}\ \ (=T_{x}^{\alpha})$
(27) $\displaystyle T_{y(5y^{2}-3r^{2})}$
$\displaystyle=K_{y}(5K_{y}^{2}-3\bm{K}^{2})\ \ (=T_{y}^{\alpha})$ (28)
$\displaystyle T_{z(5z^{2}-3r^{2})}$
$\displaystyle=\overline{K_{z}(5K_{z}^{2}-3\bm{K}^{2})}\ \ (=T_{z}^{\alpha})$
(29) $\displaystyle T_{x(y^{2}-z^{2})}$
$\displaystyle=\overline{K_{x}(K_{y}^{2}-K_{z}^{2})}\ \ (=T_{x}^{\beta})$ (30)
$\displaystyle T_{y(z^{2}-x^{2})}$
$\displaystyle=\overline{K_{y}(K_{z}^{2}-K_{x}^{2})}\ \ (=T_{y}^{\beta})$ (31)
$\displaystyle T_{z(x^{2}-y^{2})}$
$\displaystyle=\overline{K_{z}(K_{x}^{2}-K_{y}^{2})}\ \ (=T_{z}^{\beta})$ (32)

where the overline symmetrizes the expression as
$\overline{ABC}=(ABC+ACB+BAC+BCA+CAB+CBA)/3!$, for example. These operators
are referred to as monopole ($N$), dipole ($M_{\mu}$), quadrupole
($Q_{\lambda}$), and octupole ($T_{\xi}$) in accordance with the number of
multiplied angular momenta [26, 31]. The further high-rank tensors are zero.
$6\times 6=36$ multipoles are expected, but there are only 16 in the above,
which is not enough. The complete matrix basis can be constructed using the
local projection operators $P_{3/2}$ and $P_{1/2}$, each of which singles out
the $j_{\mathrm{eff}}=3/2$ and $j_{\mathrm{eff}}=1/2$ components, respectively
($P_{3/2}+P_{1/2}=1$). First, the operators in the $j_{\mathrm{eff}}=3/2$
subspace are introduced as

$\displaystyle N^{\rm 3/2,even}(\equiv N^{a})$ $\displaystyle\propto P_{3/2}$
(33) $\displaystyle M^{\rm 3/2,odd}_{\mu}(\equiv M^{a}_{\mu})$
$\displaystyle\propto P_{3/2}M_{\mu}P_{3/2}$ (34) $\displaystyle Q^{\rm
3/2,even}_{\lambda}(\equiv Q^{a}_{\lambda})$ $\displaystyle\propto
P_{3/2}Q_{\lambda}P_{3/2}$ (35) $\displaystyle T^{\rm 3/2,odd}_{\xi}(\equiv
T^{a}_{\xi})$ $\displaystyle\propto P_{3/2}T_{\xi}P_{3/2}=T_{\xi}$ (36)

where $\mu,\lambda,\xi$ represent polynomials, and the ones in
$j_{\mathrm{eff}}=1/2$ subspace are defined by

$\displaystyle N^{\rm 1/2,even}(\equiv M^{b}_{\mu})$ $\displaystyle\propto
P_{1/2}$ (37) $\displaystyle M^{\rm 1/2,odd}_{\mu}(\equiv M^{b}_{\mu})$
$\displaystyle\propto P_{1/2}M_{\mu}P_{1/2}$ (38)

The “even” and “odd” symbols represent the sign from the time-reversal
operation $\mathscr{T}=\exp(-i\pi J_{y})\mathscr{K}$ with the complex
conjugation $\mathscr{K}$. For examples, one can confirm

$\displaystyle\mathscr{T}M^{\rm 3/2,odd}_{\mu}\mathscr{T}^{-1}$
$\displaystyle=-M^{\rm 3/2,odd}_{\mu}$ (39) $\displaystyle\mathscr{T}Q^{\rm
3/2,even}_{\lambda}\mathscr{T}^{-1}$ $\displaystyle=+Q^{\rm
3/2,even}_{\lambda}$ (40)

and so on. The quadrupoles and octupoles do not exist for the
$j_{\mathrm{eff}}=1/2$ subspace:
$P_{1/2}Q_{\lambda}P_{1/2}=P_{1/2}T_{\xi}P_{1/2}=0$. We used a short-hand
notation such as $N^{a}$ and $M^{b}$ in the main text (see Tab. I).

The $j_{{\rm eff}}$-off-diagonal (the symbol “offd” assigned) operators are
also obtained, which are classified by the time-reversal operation. The time-
reversal odd $j_{\rm eff}$-off-diagonal operators are

$\displaystyle M^{\rm offd,odd}_{\mu}(\equiv M^{c}_{\mu})$
$\displaystyle\propto P_{3/2}M_{\mu}P_{1/2}+P_{1/2}M_{\mu}P_{3/2}.$ (41)
$\displaystyle Q^{\rm offd,odd}_{\lambda}(\equiv Q^{d}_{\lambda})$
$\displaystyle\propto-i(P_{3/2}Q_{\lambda}P_{1/2}-P_{1/2}Q_{\lambda}P_{3/2}).$
(42)

which are magnetic dipoles and magnetic quadrupoles. In addition, there are
time-reversal even ones

$\displaystyle M^{\rm offd,even}_{\mu}(\equiv M^{d}_{\mu})$
$\displaystyle\propto-i(P_{3/2}M_{\mu}P_{1/2}-P_{1/2}M_{\mu}P_{3/2}),$ (43)
$\displaystyle Q^{\rm offd,even}_{\lambda}(\equiv Q^{c}_{\lambda})$
$\displaystyle\propto P_{3/2}Q_{\lambda}P_{1/2}+P_{1/2}Q_{\lambda}P_{3/2}.$
(44)

which are electric dipoles and electric quadrupoles. In particular, for the
electric dipoles, there is another symbolic expression

$\displaystyle\bm{M}^{\rm offd,even}\propto\bm{L}\times\bm{S}$ (45)

which is identified as time-reversal even.

Thus, the complete 36 basis matrices have been constructed. Once the matrix
representation of these operators is normalized by the trace of the squared
matrices, these matrices do not depend on the constants $\alpha$ and $\beta$.
The simple choice is then $\alpha=\beta=1$ in the present case. The above
matrices are also orthogonal as

$\displaystyle{\rm
Tr\,}O_{\xi}O^{\dagger}_{\xi^{\prime}}=\delta_{\xi\xi^{\prime}}$ (46)

where the indices $\xi$ take the 36 types of multipoles. With these setups,
the magnitudes of the expectation values of each multipole operator can be
compared and the primary order parameter that is largest in magnitude can be
found.

The classification based on the point group is also possible. The ranks
0,1,2,3 respectively correspond to the monopoles ($N^{\rm 3/2,even},N^{\rm
1/2,even}$), dipoles ($M^{\rm 3/2,odd},M^{\rm 1/2,odd},M^{\rm offd,odd},M^{\rm
offd,even}$), quadrupoles ($Q^{\rm 3/2,even},Q^{\rm offd,odd},Q^{\rm
offd,even}$), and octupoles ($N^{\rm 3/2,odd}$), respectively, which are
defined above. Each rank is further decomposed based on the cubic harmonics.
The obtained results are summarized in Table 1 in the main text. $j$-parity
for the complete classification of the multipole moments is introduced. This
is similar to the case of $sp$ hybridized systems, where the spatial parity is
$(-1)^{\ell}$, giving $+1$ for $s$-electrons ($\ell=0$) and $-1$ for
$p$-electrons ($\ell=1$). In the same manner, a similar parity transformation
is considered to distinguish $j_{\rm eff}=3/2$ from $j_{\rm eff}=1/2$, where
the total angular momentum differs by 1. Namely, the transformation matrix is
defined using the projection operators as

$\displaystyle\mathscr{P}$ $\displaystyle=P_{3/2}-P_{1/2}$ (47)

which transforms the wave function as

$\displaystyle\mathscr{P}|j_{\rm eff}j_{{\rm eff},z}\rangle$
$\displaystyle=(-1)^{j_{\rm eff}+1/2}|j_{\rm eff}j_{{\rm eff},z}\rangle$ (48)

and the transformed operators are, for example,

$\displaystyle\mathscr{P}M_{\mu}^{3/2,\rm odd}\mathscr{P}^{-1}$
$\displaystyle=+M_{\mu}^{3/2,\rm odd}$ (49)
$\displaystyle\mathscr{P}M_{\mu}^{\rm offd,odd}\mathscr{P}^{-1}$
$\displaystyle=-M_{\mu}^{\rm offd,odd}$ (50)

The $j$-parity gives $-1$ for the $j_{\mathrm{eff}}$-off-diagonal components.
Thus, every multipole is uniquely classified in terms of rank, time reversal,
and $j$parity. Note that all the multipoles in this study are even under real
spatial parity transformation because the $d$ electrons are considered. It is
also noted that the concept of $j$-parity is introduced for the classification
of the multipole moments, and the interacting Hamiltonian is not invariant
under this transformation. This is why $j_{\rm eff}$-off-diagonal multipoles
mix in general, as shown in the numerical results.

The multipole expansion is regarded as the choice of a set of basis matrices
and is not unique. Our multipole basis based on the total angular momentum is
different from the previously proposed ones [10, 26] in that $j_{\rm
eff}$-diagonal and $j_{\rm eff}$-off-diagonal components are classified. If we
restrict ourselves to the $j_{\mathrm{eff}}=3/2$ diagonal subspace, the
matrices are the same as those used in Ref. [31].

##  Appendix B Symmetry in orbital space

Here a comment on the symmetry in orbital space is made. In the case without
spin-orbit coupling, the Hamiltonian with the hopping term and Slater-Kanamori
interaction has SO(3) symmetry. Namely, the Hamiltonian is invariant under the
transformation

$\displaystyle
c_{i\alpha\sigma}\longrightarrow\sum_{\beta}V_{\alpha\beta}c_{i\beta\sigma}$
(51)

where $V_{\alpha\beta}\in\mathbb{R}$ is an orthogonal $3\times 3$ matrix.

If the full-spin polarized situation is considered, as in phase MM(II) in Fig.
8, then, only the $\sigma=\uparrow$ components need to be looked at. Then, the
symmetry in the orbital space is elevated to SU(3), where the transformation
matrix $V$ is unitary with complex variables. This emergent symmetry is broken
once the down spin components are mixed. In addition to the above spin-
polarized situation, when the spin-orbit coupling is turned on, we have an
additional term $\propto L_{z}S_{z}$, which acts as the magnetic field in the
orbital space.

It is shown that a full-spin polarization significantly simplifies the self-
consistent equations. The Hamiltonian is

$\displaystyle\mathscr{H}$
$\displaystyle=\sum_{\bm{k}\alpha}(\varepsilon_{\bm{k}}-\mu)c^{\dagger}_{\bm{k}\alpha\uparrow}c_{\bm{k}\alpha\uparrow}+\frac{\lambda}{2}\sum_{i\alpha\beta}\ell^{z}_{\alpha\beta}c^{\dagger}_{i\alpha\uparrow}c_{i\beta\uparrow}$
$\displaystyle+\frac{U-J}{2}\sum_{\alpha\beta}c^{\dagger}_{i\alpha\uparrow}c_{i\alpha\uparrow}c^{\dagger}_{i\beta\uparrow}c_{i\beta\uparrow}$
(52)

where $\alpha,\beta=xy,yz,zx$. It is notable that the Hamiltonian has a
particle-hole symmetry within the $\uparrow$-spin sector. The $\alpha$-basis
is changed to the $\ell_{z}(=m)$basis by a unitary transformation ($m=0,\pm
1$). Then the Hamiltonian becomes

$\displaystyle\mathscr{H}$
$\displaystyle=\sum_{\bm{k}m}(\varepsilon_{\bm{k}}-\mu)c^{\dagger}_{\bm{k}m\uparrow}c_{\bm{k}m\uparrow}+\frac{\lambda}{2}\sum_{im}mc^{\dagger}_{im\uparrow}c_{im\uparrow}$
$\displaystyle+\frac{U-J}{2}\sum_{mm^{\prime}}c^{\dagger}_{im\uparrow}c_{im\uparrow}c^{\dagger}_{im^{\prime}\uparrow}c_{im^{\prime}\uparrow}$
(53)

The density-type mean-field is assumed and the self-consistent equations
obtained, as follows :

$\displaystyle\Delta_{m}$ $\displaystyle=(U-J)\sum_{m^{\prime}\neq
m}F(\mu-\Delta_{m^{\prime}}-\frac{\lambda}{2}m^{\prime})$ (54) $\displaystyle
n$ $\displaystyle=\sum_{m}F(\mu-\Delta_{m}-\frac{\lambda}{2}m)$ (55)

at zero temperature, where

$\displaystyle F(\varepsilon)$
$\displaystyle=\int_{-\infty}^{\varepsilon}D(\varepsilon^{\prime})d\varepsilon^{\prime}$
(56)

is the integrated density of states. The above argument relies only on the
spin-polarized situation and is applicable to cases with general filling,
interaction, and small spin-orbit coupling. It is confirmed that the solution
of these simplified equations is consistent with the ones discussed in the
main text inside the magnetic phases at small $\lambda$.