# On positional representation of integer vectors

Edita Pelantová and  Tomáš Vávra Department of Mathematics, FNSPE Czech
Technical University in Prague
Trojanova 13, 120 00 Praha 2, Czech Republic  Department of Pure Mathematics,
University of Waterloo, Waterloo, Ontario, Canada N2L 3G1

###### Abstract.

We show that any $m\times m$ matrix $M$ with integer entries and $\det
M=\Delta\neq 0$ can be equipped by a finite digit set
$\mathcal{D}\subset\mathbb{Z}^{m}$ such that any integer $m$-dimensional
vector belongs to the set

${\rm Fin}_{\mathcal{D}}(M)=\Bigl{\\{}\sum_{k\in I}M^{k}{d}_{k}:\emptyset\neq
I\text{ finite subset of }\mathbb{Z}\text{ and }{d}_{k}\in\mathcal{D}\text{
for each }k\in
I\Bigr{\\}}\subset\bigcup\limits_{k\in\mathbb{N}}\frac{1}{\Delta^{k}}\mathbb{Z}^{m}\,.$

We also characterize the matrices $M$ for which the sets ${\rm
Fin}_{\mathcal{D}}(M)$ and
$\bigcup\limits_{k\in\mathbb{N}}\frac{1}{\Delta^{k}}\mathbb{Z}^{m}$ coincide.

Keywords: vector representation, number system, Jordan form

2000MSC: 11C20, 11A63, 15B36

##  1\. Introduction

The idea to represent $m$-dimensional vectors by a single string of digits can
be traced back to the work of A. Vince [18] and [19] who showed that for any
expansive matrix $M\in\mathbb{Z}^{m\times m}$ there exists a digit set
$\mathcal{D}\subset\mathbb{Z}^{m}$ such that any integer vector
$x\in\mathbb{Z}^{m}$ can be written in the form
$x=\sum_{k=0}^{n}M^{k}{d}_{k}$, where $d_{k}\in\mathcal{D}$. In other words,
the whole $\mathbb{Z}^{m}$ is representable in the matrix numeration system
$(M,\mathcal{D})$. On the other hand, if a matrix $M$ has an eigenvalue inside
the unit circle, no choice of the digit set $\mathcal{D}\subset\mathbb{Z}^{m}$
allows to represent all integer vectors as a combination of non-negative
powers of $M$ only.

Many works devoted to positional representations of elements of a commutative
finitely generated ring can be interpreted as a special case of the matrix
numeration systems. From this point of view, the history of the matrix
numeration systems had started several decades before the year 1993, the year
Vince’s results have been published. Positional number systems used to
represent the Gaussian integers or, more general, elements of a ring of
integers in a quadratic number field can be viewed as a matrix number system
given by a $2\times 2$ matrix. Number systems of this type were studied by
Penney [15], Kátai and Szabó [10], Kátai and B. Kovács, [8], [9], Gilbert [5].
This concept was extended to algebraic fields of higher order by B. Kovács
[11] and B. Kovács and Pethő [12]. A number system with an algebraic base
$\beta$ can be interpreted as a matrix number system with the base being the
companion matrix of the minimal polynomial of $\beta$. Hence the
characteristic polynomial of the matrix is irreducible over $\mathbb{Q}$. The
concept of the so-called canonical number system for irreducible polynomials
was further generalized to arbitrary polynomials from $\mathbb{Z}[x]$ by Pethő
[16]. The dynamic properties of canonical number systems can be well studied
in the formalism of the shift radix systems introduced in [1]. A very general
setting of the canonical number systems was recently considered in [2].

Our aim is to study the matrix numeration systems. The most important property
of a square matrix can be deduced from its Jordan form. Let us point out that
the Jordan form of the companion matrix associated to a polynomial has a
specific property: only one Jordan block corresponds to each eigenvalue.
Therefore the study of matrix numeration systems can display new phenomena and
is of its own interest. The matrix formalism for numeration systems – under
the name numeration systems in lattices – was systematically used by A. Kovács
in [13]. Of course, A. Kovács, just like Vince, considers integer matrices as
they map a lattice into itself. The Vince’s results on integer matrices were
recently generalized by J. Jankauskas and J. Thuswaldner to matrices
$M\in\mathbb{Q}^{m\times m}$ with rational entries and without eigenvalues in
modulus strictly smaller than 1, see [7].

In this paper, we generalize the Vince’s result in another direction: our
representation of an integer vector can use both positive and negative powers
of a matrix $M\in\mathbb{Z}^{m\times m}$ with $\det M=\Delta\neq 0$. Our aim
is to study the set

${\rm Fin}_{\mathcal{D}}(M)=\Bigl{\\{}\sum_{k\in I}M^{k}{d}_{k}:I\text{ finite
non-empty subset of }\mathbb{Z}\text{ and }{d}_{k}\in\mathcal{D}\text{ for
each }k\in
I\Bigr{\\}}\subset\bigcup\limits_{k\in\mathbb{N}}\frac{1}{\Delta^{k}}\mathbb{Z}^{m}\,.$

Notice that because the set of positions $I$ is finite, we can treat the
position not in $I$ as being occupied by the zero digit. Therefore we always
assume that zero is in the digit set $\mathcal{D}.$ We show that for any non-
singular matrix $M\in\mathbb{Z}^{m\times m}$ there exists a finite set
$\mathcal{D}\subset\mathbb{Z}^{m}$ such that $\mathbb{Z}^{m}\subset{\rm
Fin}_{\mathcal{D}}(M)$. For such $\mathcal{D}$, we also characterise the
matrices $M$ for which the sets ${\rm Fin}_{\mathcal{D}}(M)$ and
$\bigcup\limits_{k\in\mathbb{N}}\frac{1}{\Delta^{k}}\mathbb{Z}^{m}$ coincide.

##  2\. The main result

This section contains the main result and its proof that relies heavily on
Proposition 4. The proof of this proposition is, however, quite technical,
therefore we provide it in its own Section 3.

The method of proving the main theorem is distinct from the one used in [7].
Their setting allowed the problem to be reduced to a problem in canonical
number systems. This approach does not work in our case because we allow
eigenvalues that are smaller than one in absolute value.

######  Theorem 1.

Let $M\in\mathbb{Z}^{m\times m}$ be a non-singular matrix. Then there exists a
finite digit set $\mathcal{A}\subset\mathbb{Z}^{m}$ such that any vector
$z\in\mathbb{Z}^{m}$ can be written in the form

$z=\sum_{k\in I}M^{k}d_{k},\ \ \text{ where }\quad\emptyset\neq
I\subset\mathbb{Z},\ I\text{\ finite\ \ and }\ d_{k}\in\mathcal{A}\ \text{for
each }k\in I.$

The property $\mathbb{Z}^{m}\subseteq\mathrm{Fin}_{\mathcal{D}}(M)$ easily
implies a stronger property. Because $\mathrm{Fin}_{\mathcal{D}}(M)$ is closed
under multiplication by $M^{i}$ for $i\in\mathbb{Z},$ we easily obtain the
following corollary.

######  Corollary 2.

Let $\mathbb{Z}^{m}\subseteq\mathrm{Fin}_{\mathcal{D}}(M)$, then we have
$\bigcup_{n\in\mathbb{Z}}M^{n}\mathbb{Z}^{m}=\mathrm{Fin}_{\mathcal{D}}(M)$.

In particular, for each $M\in\mathbb{Z}^{m\times m}$ non-singular, there
exists a finite digit set $\mathcal{D}\subset\mathbb{Z}^{m}$, such that
$\mathrm{Fin}_{\mathcal{D}}(M)=\bigcup_{n\in\mathbb{Z}}M^{n}\mathbb{Z}^{m}$.

The full proof of Theorem 1 is rather technical. In order to simplify
readability, we put the technical part into Proposition 4 whose proof will be
the content of Section 3. Nevertheless, we can now present the main idea of
our work.

For $\mathcal{D}\subset\mathbb{Z}^{n}$ finite, define $T_{d}(x)=Mx-d$ and
$x_{k}=T_{d}(x_{k-1})$ with $x_{0}=x.$ We directly obtain

(1) $x=M^{-1}d_{1}+M^{-2}d_{2}+\dots+M^{-k}d_{k}+M^{-k}x_{k}.$

Assume that $S\subseteq\mathbb{Z}^{m}$ consists of vectors with the property
that there exist $d_{1},d_{2},\dots$ such that $\|x_{i}\|<C$ eventually (with
the constant being universal for the whole $S$). Then the digit set
$\mathcal{A}=\mathcal{D}+\\{d\in\mathbb{Z}^{n}:\|d\|<C\\}$ suffices for (1) to
be a finite representation of members of $S$ over $\mathcal{A}$. It then
suffices to show that for any $x\in\mathbb{Z}^{m}$, $M^{-n}x\in S$ for some
$n\in\mathbb{N}$. Indeed, we would have $M^{-n}x=\sum_{i=1}^{k}M^{-i}d_{i},$
i.e. $x=\sum_{i=1}^{k}M^{n-i}d_{i}.$

In the following, we will show that an appropriate choice of $S$ is all the
elements of $\mathbb{Z}^{m}$ whose projection into the expansive eigenspace is
bounded by a constant depending on $M$ only. Such a choice ensures that
$M^{-n}x\in S$ is achievable for all $x\in\mathbb{Z}^{m}$. It then remains to
be shown that with a suitable digit set, a sequence of iterations
$T_{d_{n}}T_{d_{n-1}}\cdots T_{d_{1}}(x)$ will have a bounded norm eventually
for any $x\in S$. The main obstacle are unimodular eigenvalues of $M$, in
particular those with different algebraic and geometric multiplicity.

Our proof uses the real Jordan form $J$ of the matrix $M$. If $M=P^{-1}JP$,

(2) ${\text{ the equality \ $z=\sum_{k\in I}M^{k}d_{k}$ \ \ gives \ \
$Pz=\sum_{k\in I}J^{k}Pd_{k}$}}.$

Instead of looking for a finite digit set from the lattice $\mathbb{Z}^{m}$
which is suitable for the matrix $M$ and representation of integer vectors, we
look for a digit set from the lattice $P\mathbb{Z}^{m}$, suitable for $J$ and
we represent lattice points of $P\mathbb{Z}^{m}$. The matrix $P$ transforming
$M$ to its Jordan form is not determined uniquely. We choose it carefully to
approximate in some sense the lattice $\mathbb{Z}^{m}$. For this purpose, we
introduce the perturbation set

$\mathcal{E}=\\{\varepsilon\in\mathbb{R}^{m}:\|\varepsilon\|_{\infty}<\tfrac{1}{3}\\}.$

In the article we work with two norms of $\mathbb{R}^{m}$:

$\|.\|_{\infty}$ denotes the norm defined by
$\|(x_{1},x_{2},\ldots,x_{m})\|_{\infty}=\max_{i}|x_{i}|$ .

$\|.\|_{2}$ denotes the Euclidean norm, i.e.
$\|(x_{1},x_{2},\ldots,x_{m})\|_{2}=\sqrt{\sum_{j=1}^{m}|x_{j}|^{2}}.$

######  Definition 3.

We say that a lattice $L\subset\mathbb{R}^{m}$ is close to the lattice
$\mathbb{Z}^{m}$, if $\mathbb{Z}^{m}\subset L+\mathcal{E}$.

For every $M\in\mathbb{R}^{m\times m}$, there exists a non-singular
$P\in\mathbb{R}^{m\times m}$ such that $PMP^{-1}M=J={\bigoplus_{k}J_{k}}$
where $J_{k}$ is a real Jordan block. Let us recall that the real Jordan block
to $\lambda\in\mathbb{C}$ is a matrix in the form

$\left(\begin{smallmatrix}R&I&&\\\\[-7.0pt] &R&\ddots&\\\
\vspace*{-1mm}&&\ddots&I\\\\[5.69054pt]
&&&R\end{smallmatrix}\right)\quad\qquad\text{with}\ R=\begin{cases}\
(\lambda)&\mbox{ if }\lambda\in\mathbb{R},\\\\[5.69054pt]
\left(\begin{smallmatrix}a&b\\\ -b&a\end{smallmatrix}\right)&\mbox{ if
}\lambda=a+ib\in\mathbb{C}\setminus\mathbb{R}\end{cases}$

Note that $I$ is a unit matrix of order 1 or 2, according to the order of $R$.

The vector space $\mathbb{R}^{m}$ can be decomposed into invariant subspaces
$\mathbb{R}^{m}=V_{e}\oplus V_{u}\oplus V_{c}$ of the matrix $J$ such that $J$
restricted to $V_{e}$ is expansive, $J$ restricted to $V_{c}$ is contractive
and $J$ restricted to $V_{u}$ is orthogonal. In other words, all eigenvalues
of $J$ restricted to $V_{e}$ are in modulus $>1$, all eigenvalues of $J$
restricted to $V_{u}$ are in modulus $=1$, and all eigenvalues of $J$
restricted to $V_{c}$ are in modulus $<1$. The subspaces $V_{e},V_{u}$ and
$V_{c}$ may be trivial. Note that $\dim V_{e}+\dim V_{u}>0$ as $M$ is an
integer non-singular matrix. Any $x\in\mathbb{R}^{m}$ can be uniquely written
as $x=x_{e}+x_{u}+x_{c}$ with $x_{e}\in V_{e}$, $x_{u}\in V_{u}$ and $x_{c}\in
V_{c}$. We will denote $x_{e}=\pi_{e}(x)$, $x_{u}=\pi_{u}(x)$, and
$x_{c}=\pi_{c}(x)$,

The proof of the main theorem will be a consequence of the following
proposition.

######  Proposition 4.

Let $J\in\mathbb{R}^{m\times m}$ be a non-singular matrix in the real Jordan
form and $\pi_{e}:\mathbb{R}^{m}\to V_{e}$ be the projection into the
expansive invariant subspace $V_{e}$ of $J$. Let $L\subset\mathbb{R}^{m}$ be a
lattice close to $\mathbb{Z}^{m}$. Then there exists a finite set
$\mathcal{D}\subset L$ and a constant $C$ such that

for any $x\in\mathbb{R}^{m}$ with $\|\pi_{e}(x)\|_{\infty}\leq 1$ and for any
sufficiently large $N\in\mathbb{N}$ we can write

(3) $x=J^{-1}d_{1}+J^{-2}d_{2}+\cdots+J^{-N}d_{N}+J^{-N}y,\ \text{ \ for some
}d_{1},d_{2},\ldots,d_{N}\in\mathcal{D}\,$ $\text{and}\ \ y\in\mathbb{R}^{m},\
\ \text{with}\ \ \|y\|_{\infty}\leq C\,.$

###### Proof of Theorem 1.

We use the statement and notation of Proposition 4. Let $M=P^{-1}JP$, where
$J$ is the real Jordan canonical form of $M$. The matrix $P$ is not given
uniquely. For example, any matrix $\alpha P$ with a non-zero $\alpha$
transforms $M$ to its Jordan form, as well. Therefore we can assume without
loss of generality that the lattice $L=P\mathbb{Z}^{m}$ is close to the
lattice $\mathbb{Z}^{m}$ in the sense of Definition 3. As $P^{-1}JP=M$, we get
$JP\mathbb{Z}^{m}=PM\mathbb{Z}^{m}\subset P\mathbb{Z}^{m}$, i.e., $J$ maps the
lattice $L$ into $L$. For such $J$ and $L$, we find by Proposition 4, the
digit set $\mathcal{D}\subset L$.

Due to (2), we have to show that there exists a finite digit set
$\mathcal{A}\subset L$ such that

(4) $z\in L\qquad\Longrightarrow\quad\ \ z=\sum_{k\in I}J^{k}d_{k}$, where
$I\subset\mathbb{Z}$ is finite and $d_{k}\in\mathcal{A}$ for each $k\in I$.

Define

$\mathcal{A}:=\mathcal{D}+\mathcal{B},\quad\text{ where }\
\mathcal{B}:=\\{d\in L:\|d\|_{\infty}\leq C\\}.$

Obviously, $\mathcal{B}$ is finite and $0\in\mathcal{B}$. Therefore
$\mathcal{A}$ is finite and $\mathcal{D}\subset\mathcal{A}\subset L$. Let us
show that our choice of the digit set $\mathcal{A}$ has the property stated in
(4).

Let $z\in L$. Since $V_{e}$ is a contractive subspace of the matrix $J^{-1}$,
there exists $j\in\mathbb{N}$ such that $\|\pi_{e}(J^{-j}z)\|_{\infty}<1$.
Applying (3) to $x=J^{-j}z$ with $N\geq j$, we find
$d_{1},\ldots,d_{n}\in\mathcal{D}\subset\mathcal{A}$ and $y$ with
$\|y\|_{\infty}\leq C$ such that

(5) $J^{-j}z=J^{-1}d_{1}+J^{-2}d_{2}+\cdots+J^{-N}d_{n}+J^{-N}y.$

Multiplying (5) by $J^{N}$, we deduce
$y+d_{N}=J^{N-j}z-\sum_{k=1}^{N}J^{N-k}d_{k}$. As $J$ maps the lattice $L$
into $L$, we have $y\in L$ and $y\in\mathcal{B}$. Obviously,
$y+d_{n}\in\mathcal{A}$. Altogether,
$z=\sum_{k=1}^{j-1}J_{-k}d_{k}+M^{j-N}(d_{N}+y)$, and all the coefficients
$d_{1}d_{2},\ldots,d_{N-1}$ and $d_{N}+y$ belong to $\mathcal{A}$.

∎

##  3\. Proof of Proposition 4

The proof of Proposition 4 will be done in following way. To the matrix $J$ in
the real Jordan form, we find a finite digit set $\mathcal{D}\subset L$, such
that the iterations of the transformations $Jx-d$ for some $d\in\mathcal{D}$,
have small norm eventually.

First we find a digit set in $\mathbb{Z}^{m}$ and then we replace each integer
digit by a close element from the lattice $L$. When working with the matrix
$J$, the integer lattice has an important advantage. It can be decomposed into
the direct sum $\mathbb{Z}^{m}=\oplus_{j}\mathbb{Z}^{m_{j}}$ where each
lattice $\mathbb{Z}^{m_{j}}$ is contained in an invariant subspace of a real
Jordan block. Therefore, one can treat each Jordan block separately.

###  3.1. A real Jordan block to an eigenvalue $\lambda$ on the unit circle.

Let us assume first that $\lambda$ is not real. The real Jordan block to such
a $\lambda$ has an even size, say $2\ell$, and the form

(6) $J=\left(\begin{smallmatrix}R&I&&\\\\[-7.0pt] &R&\ddots&\\\
\vspace*{-1mm}&&\ddots&I\\\\[5.69054pt]
&&&R\end{smallmatrix}\right)\in\mathbb{R}^{2\ell\times
2\ell}\quad\qquad\text{with }\
R=\left(\begin{smallmatrix}\cos\varphi&-\sin\varphi\\\
\sin\varphi&\cos\varphi\end{smallmatrix}\right)\ \ \ \text{and }\
I=\left(\begin{smallmatrix}1&0\\\ 0&1\end{smallmatrix}\right).$

######  Definition 5.

Let $K_{1},K_{2},\ldots,K_{\ell}$ be positive constants and $J$ the matrix
given in (6). We define the index ${\rm
ind}:\mathbb{R}^{2\ell}\mapsto\\{0,1,\cdots,\ell\\}$ in the following way. For
$x=(x_{1},x_{2},\ldots,x_{2\ell-1},x_{2\ell})^{T}$, we put

  * •

${\rm ind}(x)=0$, if $\|(x_{2i-1},x_{2i})\|_{2}<K_{i}$ for all
$i=1,2,\ldots,\ell$;

  * •

otherwise, ${\rm ind}(x)=$ the maximal index $j$ such that
$\|(x_{2j-1},x_{2j})\|_{2}\geq K_{j}$.

First we show two simple claims.

######  Claim 6.

Given $q\in\mathbb{N},$ denote
$\mathcal{B}=\\{(0,0)^{T},(3q,0)^{T},(0,3q)^{T},(-3q,0)^{T},(0,-3q)^{T}\\}\subset\mathbb{R}^{2}$.
Then for each $z\in\mathbb{R}^{2}$ there exists $b\in\mathcal{B}$ such that

$\|z-b\|_{2}\leq 6q\qquad\text{or}\qquad\|z-b\|_{2}\leq\|z\|_{2}-q.\ $

###### Proof.

Because of the symmetry of the digit set $\mathcal{B}$, it is enough to
consider $z=(z_{1},z_{2})^{T}\in\mathbb{R}^{2}$ with $z_{1}\geq z_{2}\geq 1$.

If $z_{1}\geq 3q$, put $b=(3q,0)^{T}$. Then
$\|z-b\|_{2}=\sqrt{(z_{1}-3q)^{2}+z_{2}^{2}}$. The inequality
$\|z-b\|_{2}\leq\|z\|_{2}-q$ we want to show is equivalent to
$4q+\sqrt{z_{1}^{2}+z_{2}^{2}}<3z_{1}$. The last inequality can be easily
checked, since $4q<(3-\sqrt{2})z_{1}$ and
$\sqrt{z_{1}^{2}+z_{2}^{2}}\leq\sqrt{2}z_{1}$.

If $z_{1}<3q$, put $b=(0,0)^{T}$. Then
$\|z-b\|_{2}=\sqrt{z_{1}^{2}+z_{2}^{2}}\leq\sqrt{(3q)^{2}+(3q)^{2}}\leq 6q$. ∎

######  Claim 7.

Let a constant $c_{1}\geq 0$ and
$R=\left(\begin{smallmatrix}\cos\varphi&-\sin\varphi\\\
\sin\varphi&\cos\varphi\end{smallmatrix}\right)\in\mathbb{R}^{2\times 2}$ be
given. Then there exist a constant $c_{2}>0$ and
$\mathcal{B}\subset\mathbb{Z}^{2}$ with $\\#\mathcal{B}=5$ such that for any
$u,v\in\mathbb{R}^{2}$, $\|v\|_{2}\leq c_{1},$ there exists $b\in\mathcal{B}$
such that

$\|Ru+v-b\|_{2}<c_{2}-\tfrac{1}{2}\qquad\text{or}\qquad\|Ru+v-b\|_{2}<\|u\|_{2}-1.$

###### Proof.

Find $q\in\mathbb{N},q\geq c_{1}+1$. For this $q$, Claim 6 gives us
$\mathcal{B}$. Denote $c_{2}=6q+\tfrac{1}{2}$. Any rotation $R$ preserves the
euclidean norm, i.e., $\|Ru\|_{2}=\|u\|_{2}$. According to Claim 6, for each
vector $z=Ru+v$ we find $b\in\mathcal{B}$ such that

$\|Ru+v-b\|_{2}<6q=c_{2}-\tfrac{1}{2}\quad\text{or}\quad\|Ru+v-b\|_{2}<\|Ru+v\|_{2}-q\leq\|u\|_{2}+\|v\|_{2}-q\leq\|u\|_{2}+\underbrace{c_{1}-q}_{\leq-1}\,.$

∎

######  Lemma 8.

Let $J$ be the matrix given in (6). There exist a finite digit set
$\mathcal{D}\subset\mathbb{Z}^{2\ell}$ and constants $K_{1},\dots,K_{\ell}$,
such that the function ${\rm ind}$ defined by $K_{1},\dots,K_{\ell}$ has the
following property.

For each $x\in\mathbb{R}^{2\ell}$ there exists a digit $d\in\mathcal{D}$ such
that for all $\varepsilon\in\mathcal{E}$, the vector $y=Jx-d+\varepsilon$
satisfies

  * •

${\rm ind}(y)\leq{\rm ind}(x)$;

  * •

If ${\rm ind}(y)={\rm ind}(x)=j$, then
$\|(y_{2j-1},y_{2j})\|_{2}\leq\|(x_{2j-1},x_{2j})\|_{2}-\tfrac{1}{2}$.

###### Proof.

Consider the matrix $J\in\mathbb{R}^{2\ell\times 2\ell}$ given in (6). Using
$\ell$ times Claim 7 , we will construct $\ell$ constants
$K_{1},\ldots,K_{\ell}$ to determine the function ${\rm ind}$:

by $K_{\ell}$ and $\mathcal{B}_{\ell}$ we denote the constant $c_{2}$ and the
digit set found by Claim 7 for $c_{1}=0$;

by $K_{\ell-1}$ and $\mathcal{B}_{\ell-1}$ we denote the constant $c_{2}$ and
the digit set found by Claim 7 for $c_{1}=K_{\ell}$;

by $K_{\ell-2}$ and $\mathcal{B}_{\ell-2}$ we denote the constant $c_{2}$ and
the digit set found by Claim 7 for $c_{1}=K_{\ell-1}$;

etc.

The digit set $\mathcal{D}$ is defined by

(7) $d\in{\mathcal{D}}\qquad\text{ if and only if }\qquad
d=\left(\begin{smallmatrix}b^{(1)}\\\ \vspace*{-2mm}b^{(2)}\\\
\vspace*{-1mm}\vdots\\\\[5.69054pt]
b^{(\ell)}\end{smallmatrix}\right)\in\mathbb{Z}^{2\ell},\ \ \text{where}\ \
b^{(k)}\in\mathcal{B}_{k},\text{ for }k=1,2,\ldots,\ell.$

We show that the function ${\rm ind}$ defined by the constants
$K_{1},\ldots,K_{\ell}$ and the digit set $\mathcal{D}$ have the property
declared in the statement of the lemma.

Let $x\in\mathbb{R}^{2\ell}$. Then

$\left(\\!\\!\\!\begin{array}[]{l}(Jx)_{2i-1}\\\
(Jx)_{2i}\end{array}\\!\\!\\!\right)=R\left(\\!\\!\\!\begin{array}[]{l}x_{2i-1}\\\
x_{2i}\end{array}\\!\\!\\!\right)+\left(\\!\\!\\!\begin{array}[]{l}x_{2i+1}\\\
x_{2i+2}\end{array}\\!\\!\\!\right),\ \ \ \text{ for $i=1,2,\ldots,\ell$,
where we put}\ \ \left(\\!\\!\\!\begin{array}[]{l}x_{2\ell+1}\\\
x_{2\ell+2}\end{array}\\!\\!\\!\right)=\left(\\!\\!\\!\begin{array}[]{l}0\\\
0\end{array}\\!\\!\\!\right).$

To this $x$, we define $d\in\mathcal{D}$ by determining all its components
$b^{(i)}$.

  * •

If $i<{\rm ind}(x)$, then we choose $b^{(i)}\in\mathcal{B}_{i}$ randomly.

  * •

If $i\geq{\rm ind}(x)$, then $\|(x_{2i+1},x_{2i+2})\|_{2}<K_{i+1}$ and
$b^{(i)}\in\mathcal{B}_{i}$ is determined by Claim 7 for
$u=(x_{2i-1},x_{2i})^{T}$ and $v=(x_{2i+1},x_{2i+2})^{T}$. Recall that we used
Claim 7 to find $K_{i}:=c_{2}$ for $c_{1}:=K_{i+1}$. Therefore,

(8) $\left\|\left(\\!\\!\\!\begin{array}[]{l}(Jx)_{2i-1}\\\
(Jx)_{2i}\end{array}\\!\\!\\!\right)-b^{(i)}\right\|_{2}<K_{i}-\tfrac{1}{2}\quad\text{or}\quad\left\|\left(\\!\\!\\!\begin{array}[]{l}(Jx)_{2i-1}\\\
(Jx)_{2i}\end{array}\\!\\!\\!\right)-b^{(i)}\right\|_{2}<\left\|\left(\\!\\!\\!\begin{array}[]{l}x_{2i-1}\\\
x_{2i}\end{array}\\!\\!\\!\right)\right\|_{2}-1.$

For the digit $d$ we have described above and an $\varepsilon\in\mathcal{E}$,
we focus on $y=Jx-d+\varepsilon$. Let us realize that the inequality
$\|\varepsilon\|_{\infty}<\tfrac{1}{3}$ means
$|\varepsilon_{j}|\leq\tfrac{1}{3}$ for each coordinate $j$ and thus
$\|(\varepsilon_{2i-1},\varepsilon_{2i})\|_{2}\leq\tfrac{1}{2}$.

For each $i>{\rm ind}(x)$, we have $\|(x_{2i-1},x_{2i})\|_{2}<K_{i}$, and thus
both inequalities in (8) imply

$\left\|\left(\\!\\!\\!\begin{array}[]{l}(Jx)_{2i-1}\\\
(Jx)_{2i}\end{array}\\!\\!\\!\right)-b^{(i)}+\left(\\!\\!\\!\begin{array}[]{l}\varepsilon_{2i-1}\\\
\varepsilon_{2i}\end{array}\\!\\!\\!\right)\right\|_{2}<K_{i}.$

Consequently, ${\rm ind}(y)\leq{\rm ind}(x)$.

If ${\rm ind}(y)={\rm ind}(x)=:j$, then necessarily
$\|(y_{2j-1},y_{2j})\|_{2}=\left\|\left(\\!\\!\\!\begin{array}[]{l}(Jx)_{2j-1}\\\
(Jx)_{2j}\end{array}\\!\\!\\!\right)-b^{(j)}+\left(\\!\\!\\!\begin{array}[]{l}\varepsilon_{2j-1}\\\
\varepsilon_{2j}\end{array}\\!\\!\\!\right)\right\|_{2}>K_{j}$ and the left
most inequality in (8) cannot hold true. The validity of the right most
inequality implies
$\|(y_{2i-1},y_{2i})\|_{2}\leq\|(x_{2i-1},x_{2i})\|_{2}-\tfrac{1}{2}$.

∎

In case of real $\lambda$, that is, $\lambda=\pm 1,$ the approach is
analogous, only easier. The digit set
$\\{(a_{1},\dots,a_{m})^{T}:a_{i}\in\\{-1,0,1\\}\\}$ now suffices for the
analogy Lemma 8.

######  Corollary 9.

Let $J=J_{m}(\lambda)$ with $|\lambda|=1$. Then there exists
$\mathcal{D}\subset\mathbb{Z}^{m}$ finite with the following property. For
each $x\in\mathbb{R}^{m}$ there exists a sequence $(d_{i+1})_{i\in\mathbb{N}}$
of digits from $\mathcal{D}$ such that for any sequence
$(\varepsilon_{i+1})_{i\in\mathbb{N}}$ from the perturbation set
$\mathcal{E}$, the sequence defined recursively

$x^{(0)}=x$ and $x^{(n+1)}:=Jx^{(n)}-d_{n+1}+\varepsilon_{n+1}$ for
$n\in\mathbb{N}$

satisfies ${\rm ind}(x^{(N)})=0$ for all sufficiently large $N\in\mathbb{N}$.
In particular, there exists a constant $C$ such that $\|x^{(N)}\|_{\infty}<C$
for all sufficiently large $N\in\mathbb{N}$.

###  3.2. A real Jordan block to an eigenvalue $\lambda$ strictly inside the
unit circle

######  Lemma 10.

Let $J\in\mathbb{R}^{m\times m}$ be the real Jordan block to $\lambda$ of
modulus $<1$ and $\mathcal{E}$ be the perturbation set. Then there exist a
norm $\|.\|_{c}$ of $\mathbb{R}^{m}$ and a constant $\gamma$ such that

for any $x\in\mathbb{R}^{m}$ and any $\varepsilon\in\mathcal{E}$, the vector
$y=Jx+\varepsilon$ satisfies:

  * •

If $\|x\|_{c}\leq\gamma$, then $\|y\|_{c}\leq\gamma$;

  * •

If $\|x\|_{c}\geq\gamma$, then $\|y\|_{c}\leq\|x\|_{c}-\tfrac{1}{2}$.

###### Proof.

Let us choose $\beta$ such that $1>\beta>|\lambda|$. By Theorem 3 from [6],
there exists a norm $\|x\|_{c}$ of $\mathbb{R}^{m}$ such that
$\|Jz\|_{c}\leq\beta\|z\|_{c}$ for all $z\in\mathbb{R}^{m}$. Then

$\|Jx+\varepsilon\|_{c}\leq\|Jx\|_{c}+\|\varepsilon\|_{c}\leq\beta\|x\|_{c}+E$
where $E=\max\\{\|\varepsilon\|_{c}:\varepsilon\in\mathcal{E}\\}$.

Now, it is enough to check that

$\|x\|_{c}\leq\frac{1}{1-\beta}\bigl{(}\tfrac{1}{2}+E\bigr{)}$, implies
$\beta\|x\|_{c}+E\leq\frac{1}{1-\beta}\bigl{(}\tfrac{1}{2}+E\bigr{)}$

and

$\|x\|_{c}\geq\frac{1}{1-\beta}\bigl{(}\tfrac{1}{2}+E\bigr{)}$, implies
$\beta\|x\|_{c}+E\leq\|x\|_{c}-\tfrac{1}{2}$.

Therefore, we put $\gamma=\frac{1}{1-\beta}\bigl{(}\tfrac{1}{2}+E\bigr{)}$.

∎

######  Corollary 11.

Let $J\in\mathbb{R}^{m\times m}$ be the real Jordan block to $\lambda$ of
modulus $<1$. Then there exists a constant $\gamma$ and a norm $\|x\|_{c}$ of
$\mathbb{R}^{m}$ such that for each $x\in\mathbb{R}^{m}$ and for any sequence
$(\varepsilon_{i+1})_{i\in\mathbb{N}}$ from $\mathcal{E}$, the sequence
defined recursively

$x^{(0)}=x\ \ and\ \ \ x^{(n+1)}:=Jx^{(n)}+\varepsilon_{n+1}\quad\text{ for
}\quad n\in\mathbb{N}$

satisfies $\|x^{(N)}\|_{c}<\gamma$ for all sufficiently large
$N\in\mathbb{N}$.

###  3.3. A real Jordan block to an eigenvalue $\lambda$ strictly outside the
unit circle

######  Lemma 12.

Let $J\in\mathbb{R}^{m\times m}$ be the real Jordan block to $\lambda$ of
modulus $>1$ and $\mathcal{E}$ be the perturbation set. Then there exists a
digit set $\mathcal{D}\subset\mathbb{Z}^{m}$ with the property:

for each $x\in\mathbb{R}^{m}$ with $\|x\|_{\infty}\leq 1$ there exists a digit
$d\in\mathcal{D}$ such that for all $\varepsilon\in\mathcal{E}$, the vector
$y=Jx-d+\varepsilon$ satisfies $\|y\|_{\infty}\leq 1$.

###### Proof.

Since $\|Jx\|_{\infty}\leq(2|\lambda|+1)\|x\|_{\infty}$ for each
$x\in\mathbb{R}^{m}$, we put
$\mathcal{D}=\\{d\in\mathbb{Z}^{m}:\|d\|_{\infty}\leq 2|\lambda|+2+E\\}$,
where $E=\max\\{\|\varepsilon\|_{\infty}:\varepsilon\in\mathcal{E}\\}$.

∎

######  Corollary 13.

Let $J\in\mathbb{R}^{m\times m}$ be the real Jordan block to $\lambda$ of
modulus $>1$, and $\mathcal{D}$ and $\mathcal{E}$ be as in Lemma 12. Then for
each $x\in\mathbb{R}^{m}$ there exists a sequence $(d_{i+1})_{i\in\mathbb{N}}$
of digits from $\mathcal{D}$ such that for any sequence
$(\varepsilon_{i+1})_{i\in\mathbb{N}}$ from $\mathcal{E}$ the sequence defined
recursively

$x^{(0)}=x$ and $x^{(n+1)}:=Jx^{(n)}-d_{n+1}+\varepsilon_{n+1}$ for
$n\in\mathbb{N}$,

satisfies $\|x^{(N)}\|_{\infty}<1$ for all sufficiently large
$N\in\mathbb{N}$.

###  3.4. Completion of proof of Proposition 4

Now we combine properties of all types of the real Jordan blocks.

###### Proof of Proposition 4.

Let $J_{1},J_{2},\ldots,J_{s}$ be the real Jordan blocks of the size
$m_{1},\ldots,m_{s}$ respectively, such that $J=J_{1}\oplus
J_{2}\oplus\cdots\oplus J_{s}$.

To each block $J_{k}\in\mathbb{R}^{m_{k}}$ corresponding to $\lambda$ with
$|\lambda|\geq 1$ we find a constant $C_{k}$ and a finite digit set
$\mathcal{D}_{k}\subset\mathbb{Z}^{m_{k}}$ with the properties described in
Corollaries 13 and 9.

To each block $J_{k}$ corresponding to $\lambda$ with $|\lambda|<1$, we assign
the constant $\gamma_{k}$ found in Corollary 11 for the norm $\|.\|_{c}$. As
two norms on a finite dimensional space are equivalent, there exists a
constant $C_{k}$ such that $\|x\|_{c}<\gamma_{k}$ implies
$\|x\|_{\infty}<C_{k}$. We put
$\mathcal{D}_{k}=\\{0\\}\subset\mathbb{R}^{m_{k}}$.

We use the digit sets $\mathcal{D}_{k}$ for construction of a new digit set
$\widetilde{\mathcal{D}}\subset\mathbb{Z}^{m}$

(9) $\tilde{d}\in\widetilde{\mathcal{D}}\qquad\text{ if and only if
}\qquad\tilde{d}=\left(\begin{smallmatrix}d^{(1)}\\\ \vspace*{-2mm}d^{(2)}\\\
\vspace*{-1mm}\vdots\\\\[5.69054pt]
d^{(s)}\end{smallmatrix}\right)\in\mathbb{Z}^{m},\ \ \text{where}\ \
d^{(k)}\in\mathcal{D}_{k},\text{ for }k=1,2,\ldots,s.$

Put $C=\max_{k}{C_{k}}$. In this notation, we get the following.

Claim A: for each $x\in\mathbb{R}^{m}$ with $\|\pi_{e}(x)\|_{\infty}\leq 1$
there exists a sequence $(\tilde{d}_{i})_{i\in\mathbb{N}}$ of digits from
$\widetilde{\mathcal{D}}$ such that for any sequence
$(\varepsilon_{i})_{i\in\mathbb{N}}$ from $\mathcal{E}$ the sequence defined
recursively

(10) $x^{(0)}=x$ and $x^{(n+1)}:=Jx^{(n)}-\tilde{d}_{n+1}+\varepsilon_{n+1}$
for $n\in\mathbb{N}$,

satisfies $\|x^{(N)}\|_{\infty}<C$ for all sufficiently large
$N\in\mathbb{N}$.

Since $L$ is close to $\mathbb{Z}^{m}$, we can assign to each
$\tilde{d}\in\widetilde{\mathcal{D}}$ a lattice point $d\in L$ such that
$d\in\tilde{d}+\mathcal{E}$. We obtain a new digit set $\mathcal{D}\subset L$
and the size of $\mathcal{D}$ does not exceed the size of
$\widetilde{\mathcal{D}}$. We rewrite Equation (10) for the specific choice of
the sequence $(\varepsilon_{i})_{i\in\mathbb{N}}$, namely for
$\varepsilon_{i}:={\tilde{d}}_{i}-d_{i}\in\mathcal{E}$. We get a modification
of the previous claim

Claim B: for each $x\in\mathbb{R}^{m}$ with $\|\pi_{e}(x)\|_{\infty}\leq 1$
there exists a sequence $({d}_{i})_{i\in\mathbb{N}}$ of digits from
${\mathcal{D}}$ such that the sequence defined recursively

$x^{(0)}=x$ and $x^{(n+1)}:=Jx^{(n)}-{d_{n+1}}$ for $n\in\mathbb{N}$,

satisfies $\|x^{(N)}\|_{\infty}<C$ for all sufficiently large
$N\in\mathbb{N}$.

Multiplying the equality $x^{(n+1)}=Jx^{(n)}-{d_{n+1}}$ by $J^{-n-1}$ and
summing up for $n=0,1,\ldots,N-1$ we get

$\sum_{n=0}^{N-1}J^{-n-1}x^{(n+1)}=\sum_{n=0}^{N-1}J^{-n}x^{(n)}-\sum_{n=0}^{N-1}J^{-n-1}d_{n+1}\quad\text{and
thus}\quad\sum_{n=1}^{N}J^{-n}x^{(n)}=\sum_{n=0}^{N-1}J^{-n}x^{(n)}-\sum_{n=1}^{N}J^{-n}d_{n}\,.$

It implies $J^{-N}x^{(N)}=x^{(0)}-\sum_{n=1}^{N}J^{-n}d_{n}$ with
$\|x^{(N)}\|_{\infty}<C$, as required in the statement of the proposition.

∎

##  4\. ${\rm Fin}_{\mathcal{D}}(M)$ versus
$\bigcup\limits_{k\in\mathbb{N}}\frac{1}{\Delta^{k}}\mathbb{Z}^{m}$

In general, we have that ${\rm
Fin}_{\mathcal{D}}(M)\subseteq\bigcup\limits_{k\in\mathbb{N}}\frac{1}{\Delta^{k}}\mathbb{Z}^{m}$.
In this section we describe when the inclusion turns into an equality. The
inverse to a matrix $M\in\mathbb{Z}^{m\times m}$ with the determinant
$\Delta\neq 0$ belongs to $\frac{1}{\Delta}\mathbb{Z}^{m\times m}$. Therefore,
the definition of the set ${\rm Fin}_{\mathcal{D}}(M)$ directly implies the
following properties.

######  Lemma 14.

Let $M\in\mathbb{Z}^{m\times m}$ with $\det M=\Delta\neq 0$ and
$\mathcal{D}\subset\mathbb{Z}^{m}$, $\mathcal{D}$ finite. Then

  1. (1)

$M^{k}\Bigl{(}{\rm Fin}_{\mathcal{D}}(M)\Bigr{)}={\rm Fin}_{\mathcal{D}}(M)$
for each $k\in\mathbb{Z}$ .

  2. (2)

${\rm
Fin}_{\mathcal{D}}(M)\subseteq\bigcup\limits_{k\in\mathbb{N}}\frac{1}{\Delta^{k}}\mathbb{Z}^{m}$
.

  3. (3)

If $\mathbb{Z}^{m}\subset{\rm Fin}_{\mathcal{D}}(M)$, then ${\rm
Fin}_{\mathcal{D}}(M)$ is closed under addition and subtraction.

Any classical $b$-ary numeration system, with base $b\in\mathbb{N},b\geq 2$
and the canonical digit set $\mathcal{D}=\\{0,1,\ldots,b-1\\}$ can be
considered as a matrix system with $M=b\in\mathbb{Z}^{1\times 1}$ and the
determinant $\Delta=b$. The canonical digit set allows to represent only non-
negative numbers and obviously ${\rm
Fin}_{\mathcal{D}}(M)=\bigcup\limits_{k\in\mathbb{N}}\frac{1}{b^{k}}\mathbb{N}$.

In the next auxiliary claim we denote by $e_{i}\in\mathbb{Z}^{m}$ the vector
$(\delta_{i1},\delta_{i2},\ldots,\delta_{im})^{T}$, where $\delta_{ij}$ is the
Kronecker symbol, i.e. $\delta_{ij}=0$, if $i\neq j$ and $\delta_{ii}=1$.

######  Claim 15.

${\rm
Fin}_{\mathcal{D}}(M)=\bigcup\limits_{k\in\mathbb{N}}\frac{1}{\Delta^{k}}\mathbb{Z}^{m}$
if and only if for each $i\in\\{1,2,\ldots,m\\}$ there exists
$\ell_{i}\in\mathbb{N}$ such that $\frac{1}{\Delta}{e}_{i}\in
M^{-\ell_{i}}(\mathbb{Z}^{m})\,.$

###### Proof.

$(\Leftarrow)$ Due to Item (2) of Lemma 14, it is enough to show that
$\frac{1}{\Delta^{k}}\mathbb{Z}^{m}\subset{\rm Fin}_{\mathcal{D}}(M)$ for each
$k\in\mathbb{N}$. We show it by induction on $k\in\mathbb{N}$.

The assumption $\frac{1}{\Delta}{e}_{i}\in M^{-\ell_{i}}(\mathbb{Z}^{m})$
gives
$\tfrac{1}{\Delta}\mathbb{Z}^{m}\subset\sum_{i=1}^{m}M^{-\ell_{i}}(\mathbb{Z}^{m})$.
Multiplying this inclusion by $\tfrac{1}{\Delta^{k}}$ and using the induction
hypothesis we obtain

$\tfrac{1}{\Delta^{k+1}}\mathbb{Z}^{m}\subset\sum_{i=1}^{m}M^{-\ell_{i}}\bigl{(}\tfrac{1}{\Delta^{k}}\mathbb{Z}^{m}\bigr{)}\subset\sum_{i=1}^{m}M^{-\ell_{i}}\bigl{(}{\rm
Fin}_{\mathcal{D}}(M)\bigr{)}={\rm Fin}_{\mathcal{D}}(M)\,.$

The last equality follows from Items (1) and (3) of Lemma 14.

$(\Rightarrow)$ This implication is obvious. ∎

######  Theorem 16.

Let $M\in\mathbb{Z}^{m\times m}$, $\det M=\Delta\neq 0$. Let
$\mathcal{D}\subset\mathbb{Z}^{m}$ be a finite digit set such that
$\mathbb{Z}^{m}\subset{\rm Fin}_{\mathcal{D}}(M)$. Then

${\rm
Fin}_{\mathcal{D}}(M)=\bigcup\limits_{k\in\mathbb{N}}\frac{1}{\Delta^{k}}\mathbb{Z}^{m}\qquad\Longleftrightarrow\qquad\exists\ell\in\mathbb{N}\text{
such that }M^{\ell}=\Theta\mod\Delta\,.$

In particular, if $\det M=\pm 1$, then $\mathbb{Z}^{m}={\rm
Fin}_{\mathcal{D}}(M)$ .

###### Proof.

By Claim 15 we have to study when $\frac{1}{\Delta}e_{i}\in
M^{-\ell_{i}}(\mathbb{Z}^{m})$ for some $\ell_{i}\in\mathbb{N}$. It can be
rewritten into $M^{\ell_{i}}e_{i}\in\Delta\mathbb{Z}^{m}$, or equivalently,
$M^{\ell_{i}}e_{i}=0\mod\Delta$. We denote the $i^{th}$ column of a matrix $A$
by $A_{\bullet i}$. Obviously, $M^{\ell_{i}}e_{i}=(M^{\ell_{i}})_{\bullet i}$.

Let us realize that if $(M^{\ell_{i}})_{\bullet i}=0\mod\Delta$, then
$(M^{\ell_{i}})_{\bullet i}=0\mod\Delta$ for all
$\ell\in\mathbb{N},\ell>\ell_{i}$. Indeed, $(M^{\ell_{i}})_{\bullet
i}=M^{\ell-\ell_{i}}(M^{\ell_{i}})_{\bullet i}=M^{\ell-\ell_{i}}0\mod\Delta$.
Therefore, we can look for a common exponent $\ell$ for all columns of the
matrix $M$. ∎

######  Example 17.

Let $M_{1}=\left(\begin{array}[]{rr}2&1\\\ 2&2\end{array}\right)$ and
$M_{2}=\left(\begin{array}[]{rr}1&2\\\ 2&2\end{array}\right)$. Then $\det
M_{1}=2$ and $\det M_{2}=-2$.

$M_{1}=\left(\begin{array}[]{rr}0&1\\\ 0&0\end{array}\right)\mod\Delta,\qquad
M_{1}^{2}=\Theta\mod\Delta\qquad\text{and thus \ }{\rm
Fin}_{\mathcal{D}}(M)=\bigcup\limits_{k\in\mathbb{N}}\frac{1}{2^{k}}\mathbb{Z}^{2}$

for a suitable $\mathcal{D}\subset\mathbb{Z}^{2}$.

On the other hand,

$M_{2}=\left(\begin{array}[]{rr}1&0\\\
0&0\end{array}\right)\mod\Delta,\qquad\text{and thus \ }\ \
M_{2}^{\ell}\neq\Theta\mod\Delta\ \text{ for all $\ell\in\mathbb{N}$\ }\,.$

In particular, $\Bigr{(}\\!\\!\\!\begin{array}[]{r}\tfrac{1}{2}\\\
0\end{array}\\!\\!\\!\Bigl{)}\notin{\rm Fin}_{\mathcal{D}}(M_{2})$ for any
choice $\mathcal{D}\subset\mathbb{Z}^{2}$.

##  5\. Comments

The main idea of this paper comes from [17] concerning number systems with
algebraic base without any restrictions on its Galois conjugates. We do not
give an explicit bound on the size of the digit set, even though our proof is
a constructive one. For matrices without eigenvalues on the unit circle, the
method used in [4],[3] can be used to obtain explicit digit sets, as well as
arithmetic algorithms.

To find a digit set $\mathcal{D}\subset\mathbb{Z}^{m}$ of the minimal
cardinality may be a hard problem. It is known that if $M$ is expanding, then
we need a digit set of size at least $|\det M|$ to be able to represent
integer vectors in the form $\sum_{k=0}^{n}M^{k}d_{k}$. In our case it seems
that an appropriate lower bound would be around the product of absolute values
of the eigenvalues outside the unit circle.

In any case, the digit set needs to have at least one digit beside zero. Two
digits might be already enough for some matrices of dimension two or more, as
is illustrated by the following examples.

The first example is a case of a matrix having only one eigenvalue outside the
unit circle. It can be seen from the proof of Proposition 4 that this dominant
eigenvalue will determine the size of the digit set. We will, however, utilize
the connection between matrix numeration systems and those with real base.

######  Example 18.

Let $T=\left(\\!\\!\\!\begin{array}[]{rrr}0&0&-1\\\ 1&0&1\\\
0&1&-1\end{array}\\!\\!\\!\right)$. As $\det T=-1$, we have ${\rm
Fin}_{\mathcal{D}}(T)=\mathbb{Z}^{3}$ for a suitable digit set guaranteed by
Theorem 1. The characteristic polynomial of $T$ is $x^{3}+x^{2}-x+1$. $T$ has
one negative eigenvalue $\beta=-1,839\dots$ and a pair of complex conjugate
eigenvalues in modulus smaller than 1. As shown in [14], any element of the
ring $\mathbb{Z}[\beta]=\\{a+b\beta+c\beta^{2}:a,b,c\in\mathbb{Z}\\}$ can be
written as

(11) $a+b\beta+c\beta^{2}=\sum_{i=k}^{n}a_{i}\beta^{i},\quad\text{where
}k,n\in\mathbb{Z},k\leq n\text{ and \ }a_{i}\in\\{0,1\\}\,.$

We use the isomorphism $\psi:\mathbb{Z}[\beta]\to\mathbb{Z}^{3}$ of two
additive groups given by
$\psi(a+b\beta+c\beta^{2})=\left(\\!\\!\\!\begin{array}[]{r}a\\\ b\\\
c\end{array}\\!\\!\\!\right)$. Since
$\beta(a+b\beta+c\beta^{2}))=-c+(a+c)\beta+(b-c)\beta^{2}$, we have

$\psi(\beta(a+b\beta+c\beta^{2}))=\left(\\!\\!\\!\begin{array}[]{c}-c\\\
a+c\\\
b-c\end{array}\\!\\!\\!\right)=\left(\\!\\!\\!\begin{array}[]{rrr}0&0&-1\\\
1&0&1\\\
0&1&-1\end{array}\\!\\!\\!\right)\left(\\!\\!\\!\begin{array}[]{r}a\\\ b\\\
c\end{array}\\!\\!\\!\right)=T\left(\\!\\!\\!\begin{array}[]{r}a\\\ b\\\
c\end{array}\\!\\!\\!\right).$

Applying this rule to Equation (11) we deduce that any element of
$\mathbb{Z}^{3}$ can be expressed as

$\left(\\!\\!\\!\begin{array}[]{r}a\\\ b\\\
c\end{array}\\!\\!\\!\right)=\psi(a+b\beta+c\beta^{2})=\sum_{i\in
I}\psi(\beta^{i}a_{i})=\sum_{i=k}^{n}T^{i}\psi(a_{i})=\sum_{i=k}^{n}T^{i}\left(\\!\\!\\!\begin{array}[]{c}a_{i}\\\
0\\\ 0\end{array}\\!\\!\\!\right)\,.$

In other words, ${\rm Fin}_{\mathcal{A}}(T)=\mathbb{Z}^{3}$ already for the
small digit set $\mathcal{A}=\\{(0,0,0)^{T},(1,0,0)^{T}\\}$.

The second example is a rotation matrix in $\mathbb{R}^{2}.$ Notice that in
this case the link to real/complex base numeration systems is missing, for
there is no eigenvalue of modulus $>1.$ Moreover, the actual minimal digit set
is significantly smaller that the one ensured by the the proof of Proposition
4.

######  Example 19.

Consider $M=\begin{pmatrix}0&1\\\ -1&0\end{pmatrix}$ with the digit set
$\mathcal{D}=\\{(0,0)^{T},(1,0)^{T}\\}.$

We have the identities

$M^{4k}=\begin{pmatrix}1&0\\\ 0&1\end{pmatrix},\quad
M^{4k+1}=\begin{pmatrix}0&1\\\ -1&0\end{pmatrix},\quad
M^{4k+2}=\begin{pmatrix}-1&0\\\ 0&-1\end{pmatrix},\quad
M^{4k+3}=\begin{pmatrix}0&-1\\\ 1&0\end{pmatrix}.$

Hence

$M^{4k}\begin{pmatrix}1\\\ 0\end{pmatrix}=\begin{pmatrix}1\\\
0\end{pmatrix},\quad M^{4k+1}\begin{pmatrix}1\\\
0\end{pmatrix}=\begin{pmatrix}0\\\ -1\end{pmatrix},\quad
M^{4k+2}\begin{pmatrix}1\\\ 0\end{pmatrix}=\begin{pmatrix}-1\\\
0\end{pmatrix},\quad M^{4k+3}\begin{pmatrix}1\\\
0\end{pmatrix}\begin{pmatrix}0\\\ 1\end{pmatrix}.$

It is now easy to see how the representations are created. For instance, if
$a,b>0,$ then

$\begin{pmatrix}a\\\ b\end{pmatrix}=\sum_{k=1}^{a}M^{4k}\begin{pmatrix}1\\\
0\end{pmatrix}+\sum_{k=1}^{b}M^{4k+3}\begin{pmatrix}1\\\ 0\end{pmatrix}.$

## Acknowledgements

E.P. acknowledges support of the project CZ.02.1.01/0.0/0.0/16_019/0000778.
T.V. was supported by Charles University, projects PRIMUS/20/SCI/002 and
UNCE/SCI/022.

## References

  * [1] S. Akiyama, T. Borbély, H. Brunotte, A. Pethő, and J. M. Thuswaldner, Generalized radix representations and dynamical systems. I, Acta Math. Hungar., 108 (2005), pp. 207–238.
  * [2] J.-H. Evertse, K. Győry, A. Pethő, J. Thuswaldner, Number systems over general orders, Acta Mathematica Hungarica 159(1) (2019), 187–-205.
  * [3] Ch. Frougny, P. Heller, E. Pelantová, and M. Svobodová, k-block parallel addition versus 1-block parallel addition in non-standard numeration systems, Theor. Comp. Sci., 543 (2014), 52-67.
  * [4] Ch. Frougny, E. Pelantová, and M. Svobodová. Parallel addition in non-standard numeration systems, Theor. Comp. Sci., 412 (2011), 5714-5727.
  * [5] W. J. Gilbert, Radix representations of quadratic Fields, J. Math. Anal. Appl., 83 (1981), 264–274.
  * [6] E. Isaacson and H. B. Keller, Analysis of numerical methods, John Wiley& Sons 1966.
  * [7] J. Jankauskas and J. Thuswaldner, Characterization of rational matrices that admit finite digit representations, Linear Algebra Appl., 557 (2018), 350–358.
  * [8] I. Kátai and B. Kovács, Kanonische Zahlensysteme in der Theorie der quadratischen algebraischen Zahlen, Acta Sci. Math. (Szeged), 42 (1980), 99–107.
  * [9] I. Kátai and B. Kovács, Canonical number systems in imaginary quadratic fields, Acta Math. Acad. Sci. Hungar., 37 (1981), 159–164.
  * [10] I. Kátai and J. Szabó, Canonical number systems for complex integers, Acta Sci. Math. (Szeged) 37:3-4(1975), 255–260.
  * [11] B. Kovács. Canonical number systems in algebraic number fields. Acta Math. Acad. Sci. Hungar., 37(4):405–407, 1981.
  * [12] B. Kovács and A. Pethő. Canonical systems in the ring of integers. Publ. Math. Debrecen, 30(1-2):39–45, 1983.
  * [13] A. Kovács. Number expansions in lattices. Math. Comput. Modelling, 38(7-9):909–-915, 2003. Hungarian applied mathematics and computer applications.
  * [14] Z. Krčmáriková, W. Steiner and T. Vávra, Finite beta-expansions with negative bases, Acta Math. Hungarica, vol. 152, 485–-504, 2017.
  * [15] W. Penney, A binary system for complex numbers, J. Assoc. Comput. Mach. 12 (1965) 247–248.
  * [16] A. Pethő, On a polynomial transformation and its application to the construction of a public keycryptosystem, in: A. Pethő, M. Pohst, H.G. Zimmer, H.C. Williams (Eds.), Computational NumberTheory, Proc., Walter de Gruyter, Berlin, 1991, pp. 31–-44.
  * [17] T. Vávra. Periodic representations in Salem bases, to appear in Israel Journal of Math., https://arxiv.org/abs/1812.08228.
  * [18] A. Vince. Radix representation and rep-tiling. In Proceedings of the Twenty-fourth South- eastern International Conference on Combinatorics, Graph Theory, and Computing (Boca Raton, FL, 1993), volume 98, pages 199–212, 1993.
  * [19] A. Vince. Replicating tessellations. SIAM J. Discrete Math., 6(3):501–521, 1993.


# Spin-orbit-coupling-induced phase separation in trapped Bose gases

Zhiqian Gui  Department of Physics, Shanghai University, Shanghai 200444,
China  Zhenming Zhang  CAS Key Laboratory of Quantum Information, University
of Science and Technology of China, Hefei 230026, China  Jin Su  Department of
Basic Medicine, Changzhi Medical College, Changzhi 046000, China  Hao Lyu
Quantum Systems Unit, Okinawa Institute of Science and Technology Graduate
University, Onna, Okinawa 904-0495, Japan  Yongping Zhang
<EMAIL_ADDRESS>Department of Physics, Shanghai University, Shanghai
200444, China

###### Abstract

In a trapped spin-$1/2$ Bose-Einstein condensate with miscible interactions, a
two-dimensional spin-orbit coupling can introduce an unconventional spatial
separation between the two components. We reveal the physical mechanism of
such a spin-orbit-coupling-induced phase separation. Detailed features of the
phase separation are identified in a trapped Bose-Einstein condensate. We
further analyze differences of phase separation in Rashba and anisotropic
spin-orbit-coupled Bose gases. An adiabatic splitting dynamics is proposed as
an application of the phase separation.

##  I Introduction

Phase separation is a generic phenomenon from classical physics to quantum
physics, for example, the oil-water separation and spin Hall effect [1]. Two-
component atomic Bose-Einstein condensates (BECs) provide a tunable platform
for the investigations of phase separation [2, 3, 4, 5, 6, 7, 8]. The two
components can be realized by using different atomic species or same species
with different electric hyperfine states. Such a system features intra- and
inter-component interactions. When the inter-component interactions dominate
over the intra-component interactions, two components prefer to be phase-
separated in order to minimize the inter-component interactions [9, 10]. The
interactions for phase separation are called immiscible. The immiscibility of
two-component BECs are completely tunable in experiments. Phase separation
effect induces rich physics in quantum gases, such as the formation of vector
solitons and vortex-soliton structures, coherent spin dynamics, and pattern
formations [11, 12, 13, 14, 15, 16, 17, 18, 19].

In a two-component BEC, an artificial spin-orbit coupling can be synthesized
between different hyperfine states via Raman lasers [20, 21, 22, 23]. Such a
Raman-induced spin-orbit coupling is one dimensional. Rashba spin-orbit
coupling, which is two dimensional, has also been experimentally realized in
BECs [24, 25]. The implementation of spin-orbit coupling in BECs gives rise to
exotic quantum phases and rich superfluid properties, which opens an avenue
for simulating topological matters and exploring superfluid dynamics [26, 27,
28, 29, 30, 31, 25, 32]. In a Raman-type spin-orbit-coupled BEC, a stripe
phase [33] can exist in miscible interactions [26, 28]. In contrast, for a
Rashba spin-orbit-coupled BEC, the stripe phase may appear in the immiscible
regime [34]. Very interestingly, Refs. [35, 36] have numerically found a
spatially phase-separated ground state in a Rashba-coupled and harmonically
trapped BEC with miscible interactions. Such a ground state in the miscible
regime is unexpected for a usual two-component BEC without spin-orbit
coupling. Reference [35] identifies that the exotic phase separation satisfies
a combined symmetry of parity and a spin flip. The existence of this state is
attributed by Refs. [36, 37] to a spin-dependent force. The force is intrinsic
in the presence of Rashba spin-orbit coupling and drives the two components
moving in opposite directions. The force concept provides an intuitive picture
for the unexpected phase separation. However, its weakness is obvious. The
force is proportional to the square of Rashba spin-orbit coupling strength.
Therefore, a large strength is expected to generate a larger spatial
separation. In contrast, numerical results show that the separation decreases
with an increasing strength [36, 37]. So far, the physical origin of the
unconventional phase separation in the miscible regime is yet to be addressed.
It calls for an unambiguous interpretation, since the phase separation has
already found a broad application in other excited states. In Refs. [38, 39],
a spin-orbit-coupled bright soliton is found to be spatially separated in
center-of-mass between the two components. Dynamics of the separation in
bright solitons is analyzed by varying spin-orbit coupling strength [40]. A
spin-orbit-coupled single-vortex state, in which each component carries a
singly quantized vortex, shows spatial separation between two components, and
the separation is inversely proportional to spin-orbit coupling strength [41].
Recently, dynamics of the separation is triggered by a sudden quench of spin-
orbit coupling strength in a trapped BEC [42]. All mentioned separations occur
in the miscible regime, causing a counterintuitive expectation.

In this paper, we provide the physical mechanism for the unconventionally
spin-orbit-coupling-induced phase separation. Eigenstates of a two-dimensional
spin-orbit coupling have a momentum-dependent relative phase
$\varphi(\vec{k})$ between the two components. Closely around a fixed momentum
$\vec{k}_{0}$, the relative phase may present a linear dependence
$\varphi(\vec{k})\propto(\vec{k}-\vec{k}_{0})\cdot\vec{r}_{0}$ with a constant
$\vec{r}_{0}$. The linear dependence is a momentum kick to move two components
relatively. The superposition of these eigenstates distributing around
$\vec{k}_{0}$ constitutes a spatially separated wave packet. The separation,
whose amplitude can be calculated, is a completely single-particle effect of
spin-orbit coupling. A weakly trapped BEC with two-dimensional spin-orbit
coupling is a perfect platform to simulate the phase separation. The miscible
interactions force atoms to condense at a certain spin-orbit-coupled momentum
state with a momentum-dependent relative phase. Meanwhile, weak traps broaden
the condensed momentum so that the condensation occupies momentum states in a
narrow regime, which give rise to the linear dependence of the relative phase.
We numerically identify detailed features of spin-orbit-coupling-induced phase
separation in a trapped BEC with miscible interactions by analyzing its ground
states. The phase separation matches with the single-particle prediction when
spin-orbit coupling strength dominates. We also compare the separation
differences between Rashba and an anisotropic spin-orbit coupling. Finally, as
an application of the phase separation, we propose an adiabatic splitting
dynamics.

The paper is organized as follows. In Sec. II, the physical mechanism of spin-
orbit-coupling-induced phase separation is unveiled. The separation amplitudes
are predicted. From the mechanism, we know that the separation is a single-
particle effect. In Sec. III, we identify separated features of ground states
in a trapped BEC with Rashba spin-orbit coupling by the imaginary-time
evolution method and the variational method. In Sec. IV, we reveal the effect
of the anisotropy of spin-orbit coupling on the phase separation. In Sec. V,
we propose an adiabatic dynamics to dynamically split two components basing on
the phase separation. For the completeness of our discussion, immiscible-
interaction-induced phase separation is shown in Sec. VI. The conclusion
follows in Sec. VII.

##  II Spin-orbit-coupling-induced phase separation

Rashba spin-orbit-coupling-induced phase separation is a completely single-
particle effect. We reveal the physical origin of such phase separation. The
Rashba spin-orbit-coupled Hamiltonian is

$H_{\text{SOC}}=\frac{p_{x}^{2}+p^{2}_{y}}{2}+\lambda(p_{x}\sigma_{y}-p_{y}\sigma_{x}),$
(1)

where $p_{x}$ and $p_{y}$ are the momenta along the $x$ and $y$ directions
respectively, $\lambda$ is the spin-orbit coupling strength, and
$\sigma_{x,y}$ are spin-1/2 Pauli matrices. The eigenenergy of the Hamiltonian
has two bands. The lower band is

$E=\frac{k_{x}^{2}+k_{y}^{2}}{2}-\lambda\sqrt{k_{x}^{2}+k_{y}^{2}},$ (2)

with associated eigenstates being

$\Phi=\frac{1}{\sqrt{2}}e^{ik_{x}x+ik_{y}y}\begin{pmatrix}e^{i\frac{\varphi}{2}}\\\
e^{-i\frac{\varphi}{2}}\end{pmatrix}.$ (3)

Since the Hamiltonian possesses continuously translational symmetry, the
eigenstates are plane waves with $k_{x,y}$ being the quasimomenta along the
$x$ and $y$ directions, respectively. The outstanding feature is that Rashba
spin-orbit coupling generates a relative phase $\varphi$ between the two
components, which satisfies

$\tan(\varphi)=\frac{k_{x}}{k_{y}}.$ (4)

It is noted that $(k_{x},k_{y})=(0,0)$ is a singularity, closely around which
the relative phase cannot be defined. Therefore, the eigenstate in Eq. (3)
works beyond the regime around the singularity.

We construct a wave packet by superposing these eigenstates,

$\Psi=\int_{-\infty}^{\infty}dk_{x}dk_{y}\mathcal{G}({\mathbf{k}}-\mathbf{\bar{k}})\Phi,$
(5)

with the superposition coefficient $\mathcal{G}$ being a momentum-dependent
localized function centering around $\mathbf{\bar{k}}$. For a straightforward
illustration, we take a Gaussian distribution as an example,

$\mathcal{G}({\mathbf{k}}-\mathbf{\bar{k}})=\frac{1}{2\pi\sqrt{\Delta_{x}\Delta_{y}}}e^{-\frac{(k_{x}-\bar{k}_{x})^{2}}{2\Delta_{x}}-\frac{(k_{y}-\bar{k}_{y})^{2}}{2\Delta_{y}}}.$
(6)

The Gaussian distributed superposition coefficient is centered at
$\mathbf{\bar{k}}=(\bar{k}_{x},\bar{k}_{y})$ with the packet widths
$\sqrt{\Delta_{x,y}}$ along $x$ and $y$ directions. If the widths are narrow,
the superposition mainly happens around $\mathbf{\bar{k}}$. Therefore, we
analyze the eigenstates around $\mathbf{\bar{k}}$, and the relative phase
becomes

$\varphi({\mathbf{k}})\approx\varphi(\mathbf{\bar{k}})+\left.(k_{x}-\bar{k}_{x})\frac{\partial\varphi}{\partial
k_{x}}\right|_{\mathbf{\bar{k}}}+\left.(k_{y}-\bar{k}_{y})\frac{\partial\varphi}{\partial
k_{y}}\right|_{\mathbf{\bar{k}}},$ (7)

which is linearly dependent on the momenta $k_{x,y}$. This is true since
expanding any continuous function around a certain parameter point leads to
dominant linear-dependence. Such linear dependence in Eq. (7) induced a
momentum kick, generating the relative motion between the two components.
After substituting the Gaussian distribution in Eq. (6) and $\varphi$ in Eq.
(7) into Eq. (5) and performing integration, we get the wave packet,

$\displaystyle\Psi=$
$\displaystyle\frac{1}{\sqrt{2}}e^{i\bar{k}_{x}x+i\bar{k}_{y}y}$
$\displaystyle\times\begin{pmatrix}e^{-\frac{\Delta_{x}}{2}\left[x+\frac{1}{2}\frac{\partial\varphi(\mathbf{\bar{k}})}{\partial
k_{x}}\right]^{2}-\frac{\Delta_{y}}{2}\left[y+\frac{1}{2}\frac{\partial\varphi(\mathbf{\bar{k}})}{\partial
k_{y}}\right]^{2}+i\frac{\varphi(\mathbf{\bar{k}})}{2}}\\\
e^{-\frac{\Delta_{x}}{2}\left[x-\frac{1}{2}\frac{\partial\varphi(\mathbf{\bar{k}})}{\partial
k_{x}}\right]^{2}-\frac{\Delta_{y}}{2}\left[y-\frac{1}{2}\frac{\partial\varphi(\mathbf{\bar{k}})}{\partial
k_{y}}\right]^{2}-i\frac{\varphi(\mathbf{\bar{k}})}{2}}\end{pmatrix}.$ (8)

The outstanding feature of the resultant wave packet is that the two
components have a relative position displacement. The displacements along the
$x$ and $y$ directions are

$\left.\frac{\partial\varphi(\mathbf{k})}{\partial
k_{x}}\right|_{\mathbf{\bar{k}}}=\frac{\bar{k}_{y}}{\bar{k}_{x}^{2}+\bar{k}_{y}^{2}},\
\ \left.\frac{\partial\varphi(\mathbf{k})}{\partial
k_{y}}\right|_{\mathbf{\bar{k}}}=-\frac{\bar{k}_{x}}{\bar{k}_{x}^{2}+\bar{k}_{y}^{2}}.$
(9)

The nonzero displacements give rise to a phase separation between two
components. From the construction of the phase-separated wave packets, we can
see that the origin of the phase separation is the existence of the momentum-
dependent relative phase in eigenstates and the occupation of these
eigenstates confined in a narrow momentum regime.

Rashba spin-orbit coupled BEC is an ideal platform to generate such a phase-
separated state. The lower band in Eq. (2) has infinite energy minima which
locate at the quasimomenta satisfying $k_{x}^{2}+k_{y}^{2}=\lambda^{2}$;
therefore, $k_{x}=\lambda\cos(\theta)$ and $k_{y}=\lambda\sin(\theta)$ with
$\theta$ being an angle. The interacting atoms spontaneously choose one of
energy minima to condense and form a BEC [33, 35, 36]. This means that
$\theta$ is spontaneously chosen to be a value $\bar{\theta}$. In real atomic
BEC experiments, traps are inevitable. A weak harmonic trap naturally broadens
the BEC momentum giving rise to a Gaussian distribution centered at
$(\bar{k}_{x},\bar{k}_{y})=\lambda(\cos(\bar{\theta}),\sin(\bar{\theta}))$.
Furthermore, the broadening is narrow so that Eq. (7) is satisfied.
Consequently, the Rashba-coupled BEC presents as a phase separated state in
Eq. (II) with $\partial\varphi(\mathbf{k})/\partial
k_{x}|_{\mathbf{\bar{k}}}=\sin(\bar{\theta})/\lambda$ and
$\partial\varphi(\mathbf{k})|/\partial
k_{y}|_{\mathbf{\bar{k}}}=-\cos(\bar{\theta})/\lambda$. The position
displacement is inversely proportional to the spin-orbit coupling strength
$\lambda$, which clearly indicates Rashba spin-orbit-coupling-induced phase
separation. When the strength goes to zero ($\lambda\approx 0$), the momentum
$(\bar{k}_{x},\bar{k}_{y})\approx(0,0)$ becomes a singularity so that the
eigenstate in Eq. (3) is not physical. Without the spin-orbit coupling, the
BEC becomes the conventional one, and there is no phase separation between two
components. If the strength is enhanced gradually from zero, the position
displacements should continuously increase from zero to catch up with the
predicted value $(\sin(\bar{\theta})/\lambda,-\cos(\bar{\theta})/\lambda)$.
When the strength $\lambda$ is large enough, the position displacements
decrease towards zero again since they are inversely proportional to
$\lambda$. In this case, the plane-wave phase $(\bar{k}_{x}x+\bar{k}_{y}y)$
dominates, while the relative phase in the eigenstates [Eq. (3)] is
independent of $\lambda$, i.e.,
$\tan(\varphi)=\bar{k}_{x}/\bar{k}_{y}=\cot(\bar{\theta})$. Consequently, the
effect of the relative phase is obliterated by the plane-wave phase, and phase
separation disappears.

According to the above mechanism of the phase separation, if there is no weak
trap to broaden the condensed momentum, the spin-orbit-coupled BEC can not
present the position displacement. This is why a spatially homogeneous BEC
with spin-orbit coupling does not show phase separation as studied in most
literature. Nevertheless, in order to broaden the condensed momentum without
traps, we may consider spatially localized excitation states, such as bright
solitons and vortices. These spatially self-trapped states naturally broadens
the condensed momentum. Therefore, the resultant phase separation between two
components in Rashba spin-orbit-coupled bright solitons and quantum vortices,
which have been numerically revealed in Refs. [38, 39, 41], can be understood
by a generalization of our mechanism. Interestingly, the position displacement
of quantum vortex is inversely proportional to the spin-orbit coupling
strength as uncovered numerically in [41], can be explained unambiguously.

We emphasize that the spin-orbit-coupling-induced phase separation only works
for a two-dimensional spin-orbit coupling. For an one-dimensional spin-orbit
coupling, i.e., the Raman-induced one, the single-particle Hamiltonian is
$H^{\prime}=p_{x}^{2}/2+\lambda p_{x}\sigma_{z}+\Omega\sigma_{x}$ with
$\Omega$ being the Rabi frequency due to Raman lasers [20, 23]. The lower
energy band of this system is
$E=k_{x}^{2}/2-\sqrt{\lambda^{2}k_{x}^{2}+\Omega^{2}}$ with eigenstates being
$\Phi=e^{ik_{x}x}(-\sin(\Theta),\cos(\Theta))^{T}$. Here,
$\tan(\Theta)=\Omega/(\lambda k_{x})$, and $T$ is the transpose operator. It
is noted that there is no momentum-dependent relative phase in the
eigenstates. Therefore, the Raman-induced spin-orbit coupling, in principle,
can not generate the phase separation.

In above, we have revealed the physical mechanism of Rashba spin-orbit-
coupling-induced phase separation. We demonstrate that a weakly trapped spin-
orbit-coupled BEC satisfies requirements of the mechanism. Quantum phase in
trapped spin-orbit-coupled BECs may be phase separated states. The spin-orbit-
coupling-induced phase separation is a single-particle effect. The role of
nonlinearity in the BEC is to spontaneously choose one energy minimum for
condensation. In the following, we study ground states of a trapped spin-
orbit-coupled BEC and identify features of spin-orbit-coupling-induced phase
separation.

##  III Rashba spin-orbit-coupling-induced phase separation in trapped BECs

We consider a quasi-two-dimensional spin-1/2 BEC with Rashba spin-orbit
coupling. The trap frequency $\omega_{z}$ along the $z$ direction is assumed
to be very large so that the dynamics is completely frozen into the ground
state of the $z$-directional harmonic trap. Such the strong trap can be
implemented by an optical lattice in the $z$ direction in experiments. After
integrating the atomic state along $z$ direction, we are left with a quasi-
two-dimensional system. Rashba spin-orbit coupling can be artificially
implemented by an optical Raman lattice [24], generating the Hamiltonian
$H_{\text{SOC}}$ shown in Eq. (1). The spin-orbit-coupled BEC is described by
the following Gross-Pitaevskii (GP) equation,

$\displaystyle i\frac{\partial\Psi}{\partial
t}=\left(H_{\text{SOC}}+V+H_{\text{int}}\right)\Psi.$ (10)

with $\Psi=\left(\Psi_{1},\Psi_{2}\right)^{T}$ being the two-component wave
function. The harmonic trap in the $x-y$ plane is
$V=\frac{1}{2}\omega^{2}(x^{2}+y^{2})$ with $\omega$ the dimensionless trap
frequency. $H_{\text{int}}$ denotes nonlinear interactions,

$\displaystyle
H_{\text{int}}=\left(\begin{matrix}g\left|\Psi_{1}\right|^{2}+g_{12}\left|\Psi_{2}\right|^{2}&0\\\
0&g_{12}\left|\Psi_{1}\right|^{2}+g\left|\Psi_{2}\right|^{2}\end{matrix}\right),$
(11)

The GP equation is dimensionless, and the units of length, time, momentum and
energy are chosen as $l_{z}=\sqrt{\hbar/(m\omega_{z})}$, $1/\omega_{z}$,
$\hbar/l_{z}$ and $\hbar\omega_{z}$, respectively. With the units, the inter-
and intra-component interaction coefficients become $g=Na\sqrt{8\pi}/l_{z}$
and $g_{12}=Na_{12}\sqrt{8\pi}/l_{z}$. Here, $N$ is the atom number, and $a$
and $a_{12}$ are corresponding $s$-wave scattering lengths, respectively. The
wave functions satisfy the normalization condition, $\int
dxdy(|\Psi_{1}|^{2}+|\Psi_{2}|^{2})=1$. In numerical calculations,
experimentally accessible parameters are used. The typical trap frequency is
$\omega_{z}=2\pi\times 200$ Hz, leading to the unites of length and time
$l_{z}=0.76\mu m$ and $1/\omega_{z}=0.8ms$ respectively. $a\sim 100~{}a_{0}$
with $a_{0}$ being the Bohr radius, and $N\sim 300$, lead to $g\sim 10$. The
spin-orbit coupling strength can be changed by tuning parameters of Raman
lasers in experiments [24].

When $g>g_{12}$, the interactions are miscible. We first study ground states
of the system in this regime by performing the imaginary-time evolution of the
GP equation. The evolution is numerically implemented by the split-step
Fourier method. The window of two-dimensional space is chosen as
$(x,y)\in[-6\pi,6\pi]$ and is discretized into a 256 × 256 grid.

A typical result is shown in Fig. 1. As expected from the prediction in the
previous section, the ground state is phase-separated. The two components are
spatially separated along the $x$ direction, as shown by Figs. 1(a) and 1(b).
The ground state spontaneously chooses $\bar{\theta}=-\pi/2$ so that the atoms
condense at $(\bar{k}_{x},\bar{k}_{y})=(0,-\lambda)$, which can be seen from
the momentum-space density distributions in Figs. 1(c) and 1(d). In this case,
according to Eq. (9), the position displacement occurs along the $x$
direction, and the first component shifts by $1/(2\lambda)$ on the right side
and the second component shifts oppositely by $1/(2\lambda)$ on the left side.

Figure 1: Ground state of a trapped Rashba spin-orbit-coupled BEC with
miscible interactions. (a) and (b) Density distributions $|\Psi_{1}|^{2}$ and
$|\Psi_{2}|^{2}$ in the coordinate space. (c) and (d) Density distributions
$|\Psi_{1}|^{2}$ and $|\Psi_{2}|^{2}$ in the momentum space. The parameters
are $\omega/\omega_{z}=0.1$, $l_{z}\lambda/\hbar=0.2$, $g=12$ and $g_{12}=8$.
Figure 2:  Rashba spin-orbit-coupling-induced phase separation in a trapped
BEC. The parameters are $\omega/\omega_{z}=0.1$, $g=12$ and $g_{12}=8$. (a)
The center of mass for the two components along the $x$ direction as a
function of the spin-orbit coupling strength $\lambda$. The solid lines are
from the variational method, and the circles are obtained by the imaginary-
time evolution of the GP equation. The dashed lines are $\pm 1/(2\lambda)$
predicted from the single-particle model. The red (blue) color represents the
first (second) component. (b) The condensed momentum $\bar{k}^{\prime}_{y}$ as
a function of $\lambda$. The red solid line and blue circles are obtained by
the imaginary-time evolution and the variational method, respectively. The
variational parameters $\Delta^{\prime}_{x,y}$ are shown in (c) and (d).

In the presence of interactions, it is impossible to construct analytical wave
function of ground state from the procedure demonstrated in the previous
section. Nevertheless, the single-particle wave functions in Eq. (II) and
phase-separated results shown in Fig. 1 stimulate us to use a trial wave
function to study the phase separation by the variational method [37]. The
trial wave function is assumed to be

$\displaystyle\Psi(x,y)=\frac{(\Delta^{\prime}_{x}\Delta^{\prime}_{y})^{\frac{1}{4}}}{\sqrt{2\pi}}e^{i\bar{k}^{\prime}_{y}y}\left(\begin{array}[]{cc}e^{-\frac{\Delta^{\prime}_{x}}{2}(x-\delta_{x})^{2}-\frac{\Delta^{\prime}_{y}}{2}y^{2}}\\\
e^{-\frac{\Delta^{\prime}_{x}}{2}(x+\delta_{x})^{2}-\frac{\Delta^{\prime}_{y}}{2}y^{2}}\end{array}\right).$
(14)

Here, we have assumed that the atoms spontaneously condenses at
$(0,\bar{k}^{\prime}_{y})$ in momentum space and therefore the phase
separation only happens along the $x$ direction with the relative position
displacement $2\delta_{x}$. $1/\sqrt{\Delta^{\prime}_{x,y}}$ characterize the
widths of the wave packet along the $x,y$ directions. The unknown parameters
$\bar{k}^{\prime}_{y},\delta_{x},\Delta^{\prime}_{x,y}$ are to be determined
by minimizing the energy functional,

$\displaystyle\mathcal{E}$ $\displaystyle=\int
dxdy\Psi^{*}(H_{\text{SOC}}+V)\Psi$ $\displaystyle\phantom{={}}+\int
dxdy\left[\frac{g}{2}(|\Psi_{1}|^{4}+|\Psi_{2}|^{4})+g_{12}|\Psi_{1}|^{2}|\Psi_{2}|^{2}\right],$
(15)

Substituting the trial wave function into the energy functional $\mathcal{E}$
leads to

$\displaystyle\mathcal{E}$
$\displaystyle=\frac{\bar{k}_{y}^{{}^{\prime}2}}{2}+\frac{\Delta^{\prime}_{x}+\Delta^{\prime}_{y}}{4}\left(1+\frac{\omega^{2}}{\Delta^{\prime}_{x}\Delta^{\prime}_{y}}\right)+\frac{1}{2}\omega^{2}\delta_{x}^{2}$
$\displaystyle+\lambda(\Delta^{\prime}_{x}\delta_{x}-\bar{k}^{\prime}_{y})e^{-\Delta^{\prime}_{x}\delta_{x}^{2}}+\frac{\sqrt{\Delta^{\prime}_{x}\Delta^{\prime}_{y}}}{8\pi}\left(g+g_{12}e^{-2\Delta^{\prime}_{x}\delta_{x}^{2}}\right).$
(16)

By minimizing the energy functional with respect to the unknown parameters,
$\partial\mathcal{E}/\partial X=0$
$(X=\bar{k}^{\prime}_{y},\delta_{x},\Delta^{\prime}_{x,y})$, we obtain all
information of the trial wave function. The phase separation can be
characterized by the center of mass of each component,

$\displaystyle\bar{\bm{r}}_{1,2}=\int\bm{r}|\Psi_{1,2}(\bm{r})|^{2}d\bm{r},$
(17)

with $\bm{r}=(x,y)$. In Fig. 2(a), the solid lines show
$\bar{x}_{1,2}=\pm\delta_{x}$ calculated from the variational method, while
the results obtained by the imaginary-time evolution of the GP equation are
demonstrated by the circles. We find that the results from the two calculation
methods agree very well. Without spin-orbit coupling ($\lambda=0$), the
conventional BEC has $\bar{x}_{1,2}=0$ and condensates at
$\bar{k}^{\prime}_{y}=0$, as shown in Fig. 2(b). With the growth of $\lambda$,
$\bar{k}^{\prime}_{y}$ always increases linearly [see Fig. 2(b)]. The
displacement $\bar{x}$ first increases drastically to a maximum value and then
declines to the predicted $\pm 1/(2\lambda)$ obtained by the single-particle
model [see the dashed lines in Fig. 2(a)]. The dependence of the displacement
on $\lambda$ exactly follows the expectation in the previous section. In the
dramatic increase regime for $\bar{x}$, the variational parameters
$\Delta^{\prime}_{x,y}$ also change dramatically [see Figs. 2(c) and 2(d)].

Figure 3:  The center of mass of two components for a non-dominant
$l_{z}\lambda/\hbar=0.5$. The solid lines are obtained by the variational
method and the circles are from the imaginary-time evolution of the GP
equation. The red and blue colors represent the first and second components,
respectively. (a) The center-of-mass as a function of the inter-component
interaction coefficient $g_{12}$. $\omega/\omega_{z}=0.1$ and $g=12$. (b) The
center-of-mass as a function of the trap frequency $\omega$. $g=12$ and
$g_{12}=8$.

Rashba spin-orbit coupling introduces an intrinsic force,

$\displaystyle\bm{F}=\frac{d\bm{p}}{dt}$
$\displaystyle=-\big{[}[\bm{r},H_{\text{SOC}}],H_{\text{SOC}}\big{]}$
$\displaystyle=2\lambda^{2}(\bm{p}\times\bm{e}_{z})\sigma_{z},$ (18)

with $\bm{e}_{z}$ being the unit vector along the $z$ direction and $\bm{p}$
the atom momentum. The force originates from spin-orbit-coupling-induced
anomalous velocity [43, 44, 45]. Considering the ground states shown in Fig.
2, the force operator in momentum space becomes
$F_{x}=2\lambda^{2}\bar{k}^{\prime}_{y}\sigma_{z}$ and $F_{y}=0$. The two
components feel opposite force $F_{x}$ along the $x$ direction. Ground states
must compensate the intrinsic force to reach equilibrium. It can be
implemented by displacing two component opposite to the force. Since
$\bar{k}^{\prime}_{y}<0$ in the case shown in Fig. 2, the first component is
displaced towards to the right side and the second to the left side. The force
concept has been used in Refs. [36, 37] to explain the phase separation. Since
the force is proportional to $\lambda^{2}$, it seems that a large displacement
would be induced for a large $\lambda$. However, as shown in Fig. 2(a), the
dependence of the displacement on $\lambda$ does not follow the force. We can
see that the intrinsic force cannot explain the phase separation in the large
$\lambda$ regime.

Figure 2(a) shows that the separation follows the single-particle prediction
$\pm 1/(2\lambda)$ when $\lambda$ dominates. When $\lambda$ is weak, the
displacement also depends on other parameters, such as nonlinear coefficients
and the harmonic trap. In Fig. 3(a), we plot the displacement $\bar{x}$ as a
function of the inter-component interaction coefficient $g_{12}$ for a non-
dominant $\lambda$. The displacement slightly rises with an increasing
$g_{12}$, and it reaches the maximum when $g_{12}=g$. If $g_{12}>g$, the
interactions become immiscible, leading to ground states different from the
trial wave function in Eq. (14). The dependence of the displacement on the
trap frequency is shown in Fig. 3(b). We find that the displacement decreases
as the trap frequency increases. This is because the displacement requires
more kinetic energy in a tight trap. It is notice that there is a slight
mismatching between the results from the variational method (the solid line)
and the imaginary-time evolution (the circles) in Fig. 3. The origin of such
the mismatching is that the Gaussian profile in the trial wave function in Eq.
(14) cannot exactly describe the imaginary-time-evolution-generated wave
function as shown in Fig. 1.

Figure 4:  The anisotropic spin-orbit-coupling-induced phase separation in a
trapped BEC. The parameters are $\omega/\omega_{z}=0.1$, $g=12$, and
$g_{12}=8$. (a1) The lower band of $H^{\prime}_{\text{SOC}}$ in Eq. (19) with
$l_{z}\lambda_{1}/\hbar=0.3$ and $l_{z}\lambda_{2}/\hbar=0.6$. (a2) and (a4)
show corresponding ground-state density distributions of the first component
$|\Psi_{1}|^{2}$ in coordinate and momentum spaces, respectively. (a3) and
(a5) are for the second component $|\Psi_{2}|^{2}$. (b1)-(b5) are same as
(a1)-(a5) but with $l_{z}\lambda_{1}/\hbar=0.6$ and
$l_{z}\lambda_{2}/\hbar=0.3$. In (a2), (a3), (b2), and (b3), white stars
represent the center of wave packets predicted by the single-particle model.

##  IV The anisotropic spin-orbit-coupling-induced phase separation in trapped
BECs

Rashba-spin-orbit-coupling-induced phase separation has been analyzed in the
previous section. In the two-dimensional spin-orbit-coupled BEC experiment
[24], the spin-orbit coupling strengths are tunable, which leads to an
anisotropic coupling. It has been revealed that the anisotropic spin-orbit
coupling has a great impact on ground states of a spatially homogeneous BEC
[46]. In this section, we study anisotropic-spin-orbit-coupling-induced phase
separation. The single-particle Hamiltonian of the anisotropic spin-orbit
coupling is

$H^{\prime}_{\text{SOC}}=\frac{p_{x}^{2}+p^{2}_{y}}{2}+\lambda_{1}p_{x}\sigma_{y}-\lambda_{2}p_{y}\sigma_{x},$
(19)

with the anisotropic strengths $\lambda_{1}\neq\lambda_{2}$. The lower band of
$H^{\prime}_{\text{SOC}}$ is

$\displaystyle
E=\frac{k_{x}^{2}+k_{y}^{2}}{2}-\sqrt{\lambda_{1}^{2}k_{x}^{2}+\lambda_{2}^{2}k_{y}^{2}}.$
(20)

with the associated eigenstates being the same as Eq. (3) but having the
different relative phase which can be written as

$\tan(\varphi)=\frac{\lambda_{1}k_{x}}{\lambda_{2}k_{y}}.$ (21)

According to the mechanism of the spin-orbit-coupling-induced phase
separation, the anisotropic coupling can generate position displacements
related to the derivatives of the relative phase. The displacements along the
$x$ and $y$ directions are

$\displaystyle\left.\frac{\partial\varphi(\mathbf{k})}{\partial
k_{x}}\right|_{\mathbf{\bar{k}}}=\frac{\lambda_{1}\lambda_{2}\bar{k}_{y}}{\lambda_{1}^{2}\bar{k}_{x}^{2}+\lambda_{2}^{2}\bar{k}_{y}^{2}},$
$\displaystyle\left.\frac{\partial\varphi(\mathbf{k})}{\partial
k_{y}}\right|_{\mathbf{\bar{k}}}=-\frac{\lambda_{1}\lambda_{2}\bar{k}_{x}}{\lambda_{1}^{2}\bar{k}_{x}^{2}+\lambda_{2}^{2}\bar{k}_{y}^{2}}.$
(22)

Here, $\mathbf{\bar{k}}=(\bar{k}_{x},\bar{k}_{y})$ is the momentum at which
the atoms condense. The lowest energy minima of the lower band depend on the
anisotropy. When $\lambda_{1}<\lambda_{2}$, the two minima locate at
$(\bar{k}_{x},\bar{k}_{y})=(0,\pm\lambda_{2})$ [see Fig. 4(a1)]. They locate
at $(\bar{k}_{x},\bar{k}_{y})=(\pm\lambda_{1},0)$ when
$\lambda_{1}>\lambda_{2}$ [see Fig. 4(b1)]. With the miscible interactions,
the BEC spontaneously chooses one of these two minima to condense. The ground
state that spontaneously condenses at
$(\bar{k}_{x},\bar{k}_{y})=(0,-\lambda_{2})$ for $\lambda_{1}<\lambda_{2}$ is
demonstrated in Figs. 4(a2)-(a5). We obtain ground states by the imaginary-
time evolution of the GP equation with $H^{\prime}_{\text{SOC}}$. From the
single-particle prediction in Eq. (22), the phase separation of this ground
state happens only along the $x$ direction, and the center-of-mass of the
first component is $\lambda_{1}/(2\lambda_{2}^{2})$ and that of the second
component is $-\lambda_{1}/(2\lambda_{2}^{2})$ [see the white stars in Figs.
4(a2) and 4(a3)]. Density distributions shown in Figs. 4(a2) and 4(a3) clearly
indicate the phase separation following the predictions. The ground state that
spontaneously condenses at $(\bar{k}_{x},\bar{k}_{y})=(-\lambda_{1},0)$ for
$\lambda_{1}>\lambda_{2}$ is demonstrated in Figs. 4(b2)-(b5). The single-
particle mechanism in Eq. (22) predicts that for this ground state the
separation happens along the $y$ direction and the center-of-mass are
$\mp\lambda_{2}/(2\lambda_{1}^{2})$ for two components [see the white stars in
Figs. 4(b2) and 4(b3)]. The results from the imaginary-time evolution shown in
Figs. 4(b2) and 4(b3) match with the single-particle predictions.

These analyses have shown that the center of mass of each component strongly
depends on the ratio of the spin-orbit coupling strengths. To reveal the
dependence of phase separation on $\lambda_{2}/\lambda_{1}$, we calculate
ground states with a fixed $\lambda_{1}$ and a changeable $\lambda_{2}$ by
using the imaginary-time evolution. The results are summarized in Fig. 5,
where the circles (crosses) represent the center of mass for the first
(second) component. For $\lambda_{2}<\lambda_{1}=1$, the phase separation
occurs along the $y$ direction and $|\bar{y}|$ increases with the increase of
$\lambda_{2}$ [see red circles and crosses in Fig. 5], while $\bar{x}$ is zero
[see blue circles and crosses in Fig. 5]. When $\lambda_{2}=0$, the spin-orbit
coupling becomes one-dimensional, there is no phase separation due to the
absence of the relative phase. The results change for
$\lambda_{2}>\lambda_{1}=1$ and the phase separation along the $x$ direction
is observed. In this case, the separation decreases with $\lambda_{2}$
increasing. For a very large $\lambda_{2}$, the separation disappears since
the spin-orbit coupling effectively turns to be one-dimensional. The results
in Fig. 5 demonstrate that the maximum separation happens for
$\lambda_{1}=\lambda_{2}$ which is Rashba spin-orbit coupling. This is also
expected from the single-particle prediction in Eq. (22).

Figure 5: Anisotropic-spin-orbit-coupling-induced phase separation as a
function of $\lambda_{2}$ with a fixed $l_{z}\lambda_{1}/\hbar=1$. Circles are
for the first component and crosses are for the second component. The blue
(red) color represents separation along the $x$ ($y$) direction. Other
parameters are $\omega/\omega_{z}=0.1$, $g=12$ and $g_{12}=8$.  Figure 6:
Adiabatic splitting dynamics of a trapped BEC with Rashba spin-orbit coupling
by slowly switching off the linear coupling. The parameters are
$\omega/\omega_{z}=0.1$, $g=12$, $g_{12}=8$, $\Omega_{0}/\omega_{z}=3$, and
$\omega_{z}\tau_{q}=150$. The red (blue) dots represent the center-of-mass of
the first (second) component. Insets (a,b) [(c,d)] are density distributions
of the first and second components at $t=0$ [$t=\tau_{q}$], respectively.

##  V Adiabatic splitting dynamics

Figure 7: Immiscible-interaction-induced phase-separated ground states in a
trapped BEC with Rashba spin-orbit coupling. The parameters are $g=4$,
$g_{12}=8$ and $\omega/\omega_{z}=0.1$. (a1)-(a4) When
$l_{z}\lambda/\hbar=0.5$ the ground state is a half-quantum vortex state. (a1)
[a(3)] and (a2) [(a4)] are the coordinate (momentum) space density
distributions of the first and second components respectively. (b1)-(b4) When
$l_{z}\lambda/\hbar=1.5$ the ground state is a stripe state. (b1) [b(3)] and
(b2) [(b4)] are the coordinate [momentum] space density distributions of the
first and second components respectively.

We have shown that ground states of a trapped BEC with two-dimensional spin-
orbit coupling and miscible interactions are phase-separated. As an important
application, we study adiabatic dynamics of the phase separation. As pointed
out by previous works, a linear coupling between two component favors
miscibility regardless of interactions [47, 48]. Therefore, a miscible-to-
immiscible transition may occur by decreasing the coupling. The adiabatic
dynamics is stimulated by slowly switching off the linear coupling.
Theoretically, the process is described by the time-dependent GP equation,

$\displaystyle i\frac{\partial\Psi}{\partial
t}=\left[H_{\text{SOC}}+\Omega(t)\sigma_{x}+V+H_{\text{int}}\right]\Psi.$ (23)

Here, $\Omega(t)\sigma_{x}$ represents the linear coupling between the two
components, and can be experimentally achieved by using a radio-frequency
coupling [6]. The time-dependent Rabi frequency is

$\Omega(t)=\Omega_{0}(1-t/\tau_{q}),$ (24)

with $\Omega_{0}$ being the initial value of the linear coupling and
$\tau_{q}$ is the quench duration. At $t=0$, the presence of $\Omega_{0}$
greatly suppresses the ground-state phase separation. We obtain ground state
by the imaginary-time evolution of Eq. (23) with $\Omega(t)=\Omega_{0}$. A
typical ground state is shown in insets (a) and (b) of Fig. 6, and the
separation between two components is not obvious. Using this ground state as
initial state, we evolve the time-dependent GP equation. The center-of-mass
$\bar{x}$ for two components is recorded during the time evolution in Fig. 6.
By decreasing the linear coupling adiabatically, the separation between two
component gradually increases. When it is completely switched off, i.e.,
$t=\tau_{q}$, the separation is maximized [see corresponding density
distributions in insets (c) and (d)]. The two components can realize a
dynamically spatial splitting, which move along opposite directions. Such
adiabatic splitting dynamics are reminiscent of a kind of “atomic spin Hall
effect” [49].

##  VI Immiscible interactions induced phase separation

In all above, the interactions are miscible ($g>g_{12}$), which support atoms
to condense at a particular momentum state. On the other hand, immiscible
interactions ($g<g_{12}$) prefer a spatial separation between two components
in order to minimize the inter-component interactions proportional to
$g_{12}$. In the presence of spin-orbit coupling, the immiscible-interaction-
induced phase separation presents interesting features [36, 35, 50]. In Fig.
7, we show two different kinds of immiscible-interaction-induced phase
separated ground states with different values of spin-orbit coupling strength
$\lambda$. For $\lambda=0.5$, the ground state obtained by the imaginary-time
evolution is a half-quantum vortex state which has been first revealed in
Refs. [35, 50]. The first component distribution has a Gaussian shape [see
Figs. 7(a1) and 7(a3)], and the second component is a vortex with a winding
number $w=1$ [see Figs. 7(a2) and 7(a4)]. The first component is filled in the
density dip of the second one, forming a spatial separation along the radial
direction. For $\lambda=1.5$, the ground state becomes a stripe state which
has been first revealed in Refs. [33, 34]. The ground state condenses
simultaneously at two different momenta [see Figs. 7(b3) and 7(b4)]. Such
momentum occupation generates spatially periodic modulations in density
distributions [Figs. 7(b1) and 7(b2)]. Meanwhile, stripes of the two
components are spatially separated.

We emphasize that phase separations induced by spin-orbit coupling and
immiscible interaction have different physical origins. The spin-orbit-
coupling-induced phase separation only works for a two-dimensional spin-orbit
coupling. However, phase separations have been also studied for a BEC with a
one-dimensional spin-orbit coupling, the mechanism of which is different. In
the pioneered spin-orbit-coupled experiment, experimentalists observed a
spatial separation between two dressed states with a Raman-induced spin-orbit
coupling [20]. The spin-orbit coupling generates two energy minima, whose
occupations can be considered as two dressed states. In the dressed state
space, atomic interactions turn to be immiscible between two dressed states in
the presence of the Rabi frequency. The phase separation happens in the
dressed state space due to immiscibility. In addition, Ref. [47] reveals the
existence of phase separation in a spin-1 BEC with the Raman-induced spin-
orbit coupling. The single-particle Hamiltonian of the system is
$H=(p_{x}+\lambda^{\prime}F_{z})^{2}/2+\Omega^{\prime}F_{x}+\epsilon
F_{z}^{2}$. Here, $F_{x,y,z}$ are the spin-1 Pauli matrices,
$\lambda^{\prime}$ is the spin-orbit coupling strength, $\Omega^{\prime}$ is
the Rabi frequency, and $\epsilon$ is the quadratic Zeeman shift. The spinor
interactions include density-density part with the coefficient $c_{0}$ and
spin-spin part with the coefficient $c_{2}$. In particular, a very negative
quadratic Zeeman shift $\epsilon=-\lambda^{2}/2$ was considered. With such a
large negative $\epsilon$, the occupation in the second component can be
eliminated. The spinor only occupies the first and third components.
Interestingly, the spinor interactions between the first and third components
are immiscible for a negative spin-spin interaction ($c_{2}<0$). Different
phase-separated states between the first and third component are due to
immiscible interactions [47] .

##  VII Conclusion

In summary, we have revealed the physical mechanism of spin-orbit-coupling-
induced phase separation. The mechanism, which is very different from the
conventional immiscible-interaction-induced separation, is a complete single-
particle effect of spin-orbit coupling. We have analyzed separation features
in a trapped BEC with Rashba spin-orbit coupling and miscible interactions and
studied the effects of the anisotropy of spin-orbit coupling on the
separation. All features can be explained by the single-particle mechanism. As
an interesting application of the phase separation, we propose an adiabatic
dynamics that can dynamically split two components spatially.

## Acknowledgments

This work was supported by National Natural Science Foundation of China with
Grants No.12374247 and 11974235. H.L. acknowledges support from Okinawa
Institute of Science and Technology Graduate University.

## References

  * Sinova  _et al._ [2015] J. Sinova, S. O. Valenzuela, J. Wunderlich, C. H. Back, and T. Jungwirth, Spin Hall effects, Rev. Mod. Phys. 87, 1213 (2015).
  * Shin  _et al._ [2006] Y. Shin, M. W. Zwierlein, C. H. Schunck, A. Schirotzek, and W. Ketterle, Observation of Phase Separation in a Strongly Interacting Imbalanced Fermi Gas, Phys. Rev. Lett. 97, 030401 (2006).
  * Papp  _et al._ [2008] S. B. Papp, J. M. Pino, and C. E. Wieman, Tunable Miscibility in a Dual-Species Bose-Einstein Condensate, Phys. Rev. Lett. 101, 040402 (2008).
  * Thalhammer  _et al._ [2008] G. Thalhammer, G. Barontini, L. De Sarlo, J. Catani, F. Minardi, and M. Inguscio, Double Species Bose-Einstein Condensate with Tunable Interspecies Interactions, Phys. Rev. Lett. 100, 210402 (2008).
  * Tojo  _et al._ [2010] S. Tojo, Y. Taguchi, Y. Masuyama, T. Hayashi, H. Saito, and T. Hirano, Controlling phase separation of binary Bose-Einstein condensates via mixed-spin-channel Feshbach resonance, Phys. Rev. A 82, 033609 (2010).
  * Nicklas  _et al._ [2011] E. Nicklas, H. Strobel, T. Zibold, C. Gross, B. A. Malomed, P. G. Kevrekidis, and M. K. Oberthaler, Rabi Flopping Induces Spatial Demixing Dynamics, Phys. Rev. Lett. 107, 193001 (2011).
  * Jiménez-García  _et al._ [2019] K. Jiménez-García, A. Invernizzi, B. Evrard, C. Frapolli, J. Dalibard, and F. Gerbier, Spontaneous formation and relaxation of spin domains in antiferromagnetic spin-1 condensates, Nat. Commun. 10, 1422 (2019).
  * He  _et al._ [2020] L. He, P. Gao, and Z.-Q. Yu, Normal-Superfluid Phase Separation in Spin-Half Bosons at Finite Temperature, Phys. Rev. Lett. 125, 055301 (2020).
  * Trippenbach  _et al._ [2000] M. Trippenbach, K. Góral, K. Rzazewski, B. Malomed, and Y. B. Band, Structure of binary Bose-Einstein condensates, J. Phys. B 33, 4017 (2000).
  * Petrov [2015] D. S. Petrov, Quantum Mechanical Stabilization of a Collapsing Bose-Bose Mixture, Phys. Rev. Lett. 115, 155302 (2015).
  * Shrestha  _et al._ [2009] U. Shrestha, J. Javanainen, and J. Ruostekoski, Pulsating and Persistent Vector Solitons in a Bose-Einstein Condensate in a Lattice upon Phase Separation Instability, Phys. Rev. Lett. 103, 190401 (2009).
  * Law  _et al._ [2010] K. J. H. Law, P. G. Kevrekidis, and L. S. Tuckerman, Stable Vortex–Bright–Soliton structures in Two-Component Bose-Einstein Condensates, Phys. Rev. Lett. 105, 160405 (2010).
  * Lee  _et al._ [2016] K. L. Lee, N. B. Jørgensen, I.-K. Liu, L. Wacker, J. J. Arlt, and N. P. Proukakis, Phase separation and dynamics of two-component Bose-Einstein condensates, Phys. Rev. A 94, 013602 (2016).
  * Roy  _et al._ [2021] A. Roy, M. Ota, A. Recati, and F. Dalfovo, Finite-temperature spin dynamics of a two-dimensional Bose-Bose atomic mixture, Phys. Rev. Res. 3, 013161 (2021).
  * Eto  _et al._ [2016] Y. Eto, M. Takahashi, M. Kunimi, H. Saito, and T. Hirano, Nonequilibrium dynamics induced by miscible–immiscible transition in binary Bose–Einstein condensates, New J. Phys. 18, 073029 (2016).
  * Bernier  _et al._ [2014] N. R. Bernier, E. G. Dalla Torre, and E. Demler, Unstable Avoided Crossing in Coupled Spinor Condensates, Phys. Rev. Lett. 113, 065303 (2014).
  * Shirley  _et al._ [2014] W. E. Shirley, B. M. Anderson, C. W. Clark, and R. M. Wilson, Half-Quantum Vortex Molecules in a Binary Dipolar Bose Gas, Phys. Rev. Lett. 113, 165301 (2014).
  * Mistakidis  _et al._ [2018] S. I. Mistakidis, G. C. Katsimiga, P. G. Kevrekidis, and P. Schmelcher, Correlation effects in the quench-induced phase separation dynamics of a two species ultracold quantum gas, New J. Phys. 20, 043052 (2018).
  * Lous  _et al._ [2018] R. S. Lous, I. Fritsche, M. Jag, F. Lehmann, E. Kirilov, B. Huang, and R. Grimm, Probing the Interface of a Phase-Separated State in a Repulsive Bose-Fermi Mixture, Phys. Rev. Lett. 120, 243403 (2018).
  * Lin  _et al._ [2011] Y.-J. Lin, K. Jiménez-García, and I. B. Spielman, Spin–orbit-coupled Bose–Einstein condensates, Nature (London) 471, 83 (2011).
  * Goldman  _et al._ [2014] N. Goldman, G. Juzeliūnas, P. Öhberg, and I. B. Spielman, Light-induced gauge fields for ultracold atoms, Rep. Prog. Phys. 77, 126401 (2014).
  * Zhai [2015] H. Zhai, Degenerate quantum gases with spin–orbit coupling: a review, Rep. Prog. Phys. 78, 026001 (2015).
  * Zhang  _et al._ [2016] Y. Zhang, M. E. Mossman, T. Busch, P. Engels, and C. Zhang, Properties of spin–orbit-coupled Bose–Einstein condensates, Front. Phys. 11, 118103 (2016).
  * Wu  _et al._ [2016] Z. Wu, L. Zhang, W. Sun, X.-T. Xu, B.-Z. Wang, S.-C. Ji, Y. Deng, S. Chen, X.-J. Liu, and J.-W. Pan, Realization of two-dimensional spin-orbit coupling for Bose-Einstein condensates, Science 354, 83 (2016).
  * Valdés-Curiel  _et al._ [2021] A. Valdés-Curiel, D. Trypogeorgos, Q.-Y. Liang, R. P. Anderson, and I. B. Spielman, Topological features without a lattice in Rashba spin-orbit coupled atoms, Nat. Commun. 12, 593 (2021).
  * Li  _et al._ [2012] Y. Li, L. P. Pitaevskii, and S. Stringari, Quantum Tricriticality and Phase Transitions in Spin-Orbit Coupled Bose-Einstein Condensates, Phys. Rev. Lett. 108, 225301 (2012).
  * Manchon  _et al._ [2015] A. Manchon, H. C. Koo, J. Nitta, S. M. Frolov, and R. A. Duine, New perspectives for Rashba spin–orbit coupling, Nat. Mater. 14, 871 (2015).
  * Li  _et al._ [2017] J.-R. Li, J. Lee, W. Huang, S. Burchesky, B. Shteynas, F. Ç. Top, A. O. Jamison, and W. Ketterle, A stripe phase with supersolid properties in spin–orbit-coupled Bose–Einstein condensates, Nature (London) 543, 91 (2017).
  * Khamehchi  _et al._ [2017] M. A. Khamehchi, K. Hossain, M. E. Mossman, Y. Zhang, T. Busch, M. M. Forbes, and P. Engels, Negative-Mass Hydrodynamics in a Spin-Orbit–Coupled Bose-Einstein Condensate, Phys. Rev. Lett. 118, 155301 (2017).
  * Kartashov and Zezyulin [2019] Y. V. Kartashov and D. A. Zezyulin, Stable Multiring and Rotating Solitons in Two-Dimensional Spin-Orbit-Coupled Bose-Einstein Condensates with a Radially Periodic Potential, Phys. Rev. Lett. 122, 123201 (2019).
  * Hasan  _et al._ [2022] M. Hasan, C. S. Madasu, K. D. Rathod, C. C. Kwong, C. Miniatura, F. Chevy, and D. Wilkowski, Wave Packet Dynamics in Synthetic Non–Abelian Gauge Fields, Phys. Rev. Lett. 129, 130402 (2022).
  * Frölian  _et al._ [2022] A. Frölian, C. S. Chisholm, E. Neri, C. R. Cabrera, R. Ramos, A. Celi, and L. Tarruell, Realizing a 1D topological gauge theory in an optically dressed BEC, Nature (London) 608, 293 (2022).
  * Wang  _et al._ [2010] C. Wang, C. Gao, C.-M. Jian, and H. Zhai, Spin-Orbit Coupled Spinor Bose-Einstein Condensates, Phys. Rev. Lett. 105, 160403 (2010).
  * Ho and Zhang [2011] T.-L. Ho and S. Zhang, Bose-Einstein Condensates with Spin-Orbit interaction, Phys. Rev. Lett. 107, 150403 (2011).
  * Hu  _et al._ [2012] H. Hu, B. Ramachandhran, H. Pu, and X.-J. Liu, Spin-Orbit Coupled Weakly Interacting Bose-Einstein Condensates in Harmonic Traps, Phys. Rev. Lett. 108, 010402 (2012).
  * Zhang  _et al._ [2012] Y. Zhang, L. Mao, and C. Zhang, Mean-Field Dynamics of Spin-Orbit Coupled Bose-Einstein Condensates, Phys. Rev. Lett. 108, 035302 (2012).
  * Song  _et al._ [2013] S.-W. Song, Y.-C. Zhang, L. Wen, and H. Wang, Spin–orbit coupling induced displacement and hidden spin textures in spin-1 Bose–Einstein condensates, J. Phys. B 46, 145304 (2013).
  * Sakaguchi  _et al._ [2014] H. Sakaguchi, B. Li, and B. A. Malomed, Creation of two-dimensional composite solitons in spin-orbit-coupled self-attractive Bose-Einstein condensates in free space, Phys. Rev. E 89, 032920 (2014).
  * Xu  _et al._ [2015] Y. Xu, Y. Zhang, and C. Zhang, Bright solitons in a two-dimensional spin-orbit-coupled dipolar Bose-Einstein condensate, Phys. Rev. A 92, 013633 (2015).
  * Wang  _et al._ [2020] Y.-J. Wang, L. Wen, G.-P. Chen, S.-G. Zhang, and X.-F. Zhang, Formation, stability, and dynamics of vector bright solitons in a trapless Bose–Einstein condensate with spin–orbit coupling, New J. Phys. 22, 033006 (2020).
  * Kato  _et al._ [2017] M. Kato, X.-F. Zhang, and H. Saito, Vortex pairs in a spin-orbit-coupled Bose-Einstein condensate, Phys. Rev. A 95, 043605 (2017).
  * Ravisankar  _et al._ [2021] R. Ravisankar, H. Fabrelli, A. Gammal, P. Muruganandam, and P. K. Mishra, Effect of Rashba spin-orbit and Rabi couplings on the excitation spectrum of binary Bose-Einstein condensates, Phys. Rev. A 104, 053315 (2021).
  * Mardonov  _et al._ [2015a] S. Mardonov, M. Modugno, and E. Y. Sherman, Dynamics of spin-orbit coupled bose-einstein condensates in a random potential, Phys. Rev. Lett. 115, 180402 (2015a).
  * Mardonov  _et al._ [2015b] S. Mardonov, E. Y. Sherman, J. G. Muga, H.-W. Wang, Y. Ban, and X. Chen, Collapse of spin-orbit-coupled bose-einstein condensates, Phys. Rev. A 91, 043604 (2015b).
  * Mardonov  _et al._ [2018] S. Mardonov, V. V. Konotop, B. A. Malomed, M. Modugno, and E. Y. Sherman, Spin-orbit-coupled soliton in a random potential, Phys. Rev. A 98, 023604 (2018).
  * Ozawa and Baym [2012] T. Ozawa and G. Baym, Ground-state phases of ultracold bosons with Rashba-Dresselhaus spin-orbit coupling, Phys. Rev. A 85, 013612 (2012).
  * Gautam and Adhikari [2014] S. Gautam and S. K. Adhikari, Phase separation in a spin-orbit-coupled Bose-Einstein condensate, Phys. Rev. A 90, 043619 (2014).
  * Merhasin  _et al._ [2005] I. M. Merhasin, B. A. Malomed, and R. Driben, Transition to miscibility in a binary Bose–Einstein condensate induced by linear coupling, J. Phys. B 38, 877 (2005).
  * Beeler  _et al._ [2013] M. C. Beeler, R. A. Williams, K. Jiménez-García, L. J. LeBlanc, A. R. Perry, and I. B. Spielman, The spin Hall effect in a quantum gas, Nature (London) 498, 201 (2013).
  * Ramachandhran  _et al._ [2012] B. Ramachandhran, B. Opanchuk, X.-J. Liu, H. Pu, P. D. Drummond, and H. Hu, Half-quantum vortex state in a spin-orbit-coupled Bose-Einstein condensate, Phys. Rev. A 85, 023606 (2012).