# Discrete coherent states for $n$ qubits

Carlos Muñoz  Departamento de Física, Universidad de Guadalajara, 44420
Guadalajara, Jalisco, Mexico  Andrei B. Klimov  Departamento de Física,
Universidad de Guadalajara, 44420 Guadalajara, Jalisco, Mexico
<EMAIL_ADDRESS>Luis L. Sánchez-Soto  Departamento de Óptica, Facultad
de Física, Universidad Complutense, 28040 Madrid, Spain
<EMAIL_ADDRESS>Gunnar Björk  School of Information and Communication
Technology, Royal Institute of Technology (KTH), Electrum 229, SE-164 40
Kista, Sweden
<EMAIL_ADDRESS>
(Day Month Year; Day Month Year)

###### Abstract

Discrete coherent states for a system of $n$ qubits are introduced in terms of
eigenstates of the finite Fourier transform. The properties of these states
are pictured in phase space by resorting to the discrete Wigner function.

###### keywords:

Discrete phase space; coherent states; Wigner function.

##  1 Introduction

Discrete quantum systems were studied originally by Weyl[1] and Schwinger,[2]
and later by many authors.[3, 4] However, many concepts that appear sharp for
continuous systems become fuzzy when one tries to apply them to discrete ones.
The reason is that the in the continuum we have only one harmonic oscillator,
while for finite systems, there are a number of candidates for that role, each
one with its own virtues and drawbacks.[5]

Coherent states constitute an archetypical example of the situation: for the
standard harmonic oscillator they are well understood, and sensible
generalizations have been devised to deal with systems with more general
dynamical groups.[6] However, in spite of interesting advances,[7, 8] the
discrete counterparts are still under heated discussion.

Number states are eigenstates of the Fourier transform, but they are not the
only ones: coherent states, with their associated Gaussian wave functions,
also are. In our opinion, this subtle, yet obvious, observation has not been
taken in due consideration in this field. From this perspective, a decisive
step was given by Mehta,[9] who obtained the eigenstates of the finite Fourier
transform. The purpose of this paper is to further explore this path, showing
how we can define physically feasible coherent states for $n$ qubits that are
also Fourier eigenstates.

Since the natural arena of these discrete coherent states is the phase space,
we use the discrete Wigner function to picture the corresponding results. The
problem of generalizing the Wigner function to finite systems has a long
story. Using the notions presented in the comprehensive review Bjork2008, we
construct a Wigner function for these coherent states and discuss some of
their properties.

##  2 Coherent states for $n$ qubits

We consider $n$ identical qubits (i.e., $n$ noninteracting spin 1/2 systems).
We recall that the Dicke states, belonging to the symmetric subspace of the
representation of SU(2)⊗n, are given by

$|n,k\rangle=\sqrt{\frac{k!(n-k)!}{n!}}\sum_{k}P_{k}(|1_{1},1_{2},\ldots,1_{k},0_{k+1},\ldots,0_{n}\rangle)\,,$
(1)

where $\\{P_{k}\\}$ denotes the complete set of all the possible permutations
of the qubits. These states can be expressed in terms of the elements of the
Galois field $\mathrm{GF}(2^{n})$ using the standard decomposition in the
self-dual basis (a short review of the concepts of finite fields needed in
this paper is presented in the appendix)

$|1_{1},1_{2},\ldots,1_{k},0_{k+1},\ldots,0_{n}\rangle\mapsto|1\sigma_{1}+1\sigma_{2}+\ldots+1\sigma_{k}+0\sigma_{k+1}+\ldots+0\sigma_{n}\rangle\,.$
(2)

Now, consider a SU(2) coherent state[6]

$|\xi\rangle=\frac{1}{(1+|\xi|^{2})^{n/2}}\sum_{k=0}^{n}\sqrt{\frac{n!}{k!(n-k)!}}\;\xi^{k}|n,k\rangle\,,$
(3)

where the complex number $\xi$ is related with the angular coordinates
$(\vartheta,\varphi)$ on the Bloch sphere by

$\xi=\cot(\vartheta/2)\,e^{-i\varphi}\,.$ (4)

Using the previous correspondence, $|\xi\rangle$ can be recast as

$|\xi\rangle=\frac{1}{(1+|\xi|^{2})^{n/2}}\sum_{\gamma\in\mathrm{GF}(2^{n})}\xi^{h(\gamma)}\,|\gamma\rangle\,,$
(5)

where the function $h(\gamma)$, when applied to the field element
$\gamma=\sum_{k=1}^{n}\gamma_{k}\sigma_{k}$ indicates the number of nonzero
coefficients $\gamma_{k}$.

In this case, the Fourier operator is

$F=\frac{1}{2^{n/2}}\sum_{\mu,\nu\in\mathrm{GF}(2^{n})}\chi(\mu\nu)\,|\mu\rangle\langle\nu|\,,$
(6)

$\chi$ being an additive character defined in (36). As it is well known,
$F^{2}=\leavevmode\hbox{\small 1\normalsize\kern-3.30002pt1}$, so if we impose
that the states $|\xi\rangle$ are also eigenstates of $F$ we are lead to

$F|\xi\rangle=\pm|\xi\rangle\,.$ (7)

This immediately implies (all the spins are pointing in the same direction)
that there are two SU(2) coherent states (with $\xi_{\pm}=\pm\sqrt{2}-1$) that
simultaneously are eigenstates of the Fourier operator: they are precisely our
candidates to be coherent states for $n$ qubits. In particular,
$|\xi_{\pm}\rangle$ satisfy the following condition

$\frac{1}{2^{n/2}}\sum_{\mu,\gamma\in\mathrm{GF}(2^{n})}\xi^{h(\gamma)}\chi(\mu\gamma)|\mu\rangle=\pm\sum_{\mu\in\mathrm{GF}(2^{n})}\xi^{h(\mu)}|\mu\rangle\,,$
(8)

or, equivalently,

$\frac{1}{2^{n/2}}\sum_{\gamma\in\mathrm{GF}(2^{n})}\xi^{h(\gamma)}\chi(\mu\gamma)=\pm\xi^{h(\mu)}\,,$
(9)

and the minus sign may appear only for odd number of qubits.

Equation (5) is the abstract form of the SU(2) coherent state. It factorizes
in a product of single-qubit states when represented in the self-dual basis,
i.e.,

$|\xi\rangle=\frac{1}{(1+|\xi|^{2})^{n/2}}\sum_{c_{1},\ldots,c_{n}\in\mathbb{Z}_{2}}\xi^{h\left(\sum_{k=1}^{n}c_{k}\sigma_{k}\right)}|c_{1}\rangle\ldots|c_{n}\rangle=\prod_{j=1}^{n}\frac{(|0\rangle+\xi|1\rangle)_{j}}{(1+|\xi|^{2})^{1/2}}\,,$
(10)

$c_{k}$ being the expansion coefficients of $\gamma$ in that basis. The
operator transforming from the arbitrary basis
$\\{\varepsilon_{1},\varepsilon_{2},\ldots,\varepsilon_{n}\\}$ into a
factorized form is always a permutation given by

$P=\sum_{\mu\in\mathrm{GF}(2^{n})}(|\mu_{1}^{\prime}\rangle\ldots|\mu_{n}^{\prime}\rangle)\,(\langle\mu_{1}|\ldots\langle\mu_{n}|)\,,$
(11)

where

$\mu=\sum_{i=1}^{n}\mu_{i}\,\varepsilon_{i}=\sum_{i=1}^{n}\mu_{i}^{\prime}\,\sigma_{i}\,,\qquad\mu\in\mathrm{GF}(2^{n}),\quad\mu_{i},\mu_{i}^{\prime}\in\mathbb{Z}_{2}\,.$
(12)

Let us examine the simple yet illustrative example of a two-qubit coherent
state. In its abstract form it reads as

$|\xi\rangle=\frac{1}{1+\xi^{2}}(|0\rangle+\xi|\sigma\rangle+\xi|\sigma^{2}\rangle+\xi^{2}|\sigma^{3}\rangle).$
(13)

In the self-dual basis ($\sigma,\sigma^{2}$) we have the representation

$|0\rangle=|00\rangle=\left(\begin{array}[]{c}0\\\ 0\\\ 0\\\
1\end{array}\right),\,|\sigma\rangle=|10\rangle=\left(\begin{array}[]{c}0\\\
0\\\ 1\\\
0\end{array}\right),\,|\sigma^{2}\rangle=|01\rangle=\left(\begin{array}[]{c}0\\\
1\\\ 0\\\
0\end{array}\right),\,|\sigma^{3}\rangle=|11\rangle=\left(\begin{array}[]{c}1\\\
0\\\ 0\\\ 0\end{array}\right),$ (14)

in such a way that

$|\xi\rangle=\frac{1}{1+\xi^{2}}\left(\begin{array}[]{c}\xi^{2}\\\ \xi\\\
\xi\\\
1\end{array}\right)=\frac{1}{\sqrt{1+\xi^{2}}}\left(\begin{array}[]{c}\xi\\\
1\end{array}\right)\otimes\frac{1}{\sqrt{1+\xi^{2}}}\left(\begin{array}[]{c}\xi\\\
1\end{array}\right)\,.$ (15)

In a non self-dual basis ($\sigma,\sigma^{3}$) we have

$|0\rangle=|00\rangle=\left(\begin{array}[]{c}0\\\ 0\\\ 0\\\
1\end{array}\right),\,|\sigma\rangle=|10\rangle=\left(\begin{array}[]{c}0\\\
0\\\ 1\\\
0\end{array}\right),\,|\sigma^{3}\rangle=|01\rangle=\left(\begin{array}[]{c}0\\\
1\\\ 0\\\
0\end{array}\right),\,|\sigma^{2}\rangle=|11\rangle=\left(\begin{array}[]{c}1\\\
0\\\ 0\\\ 0\end{array}\right)\,,$ (16)

and

$|\xi\rangle=\frac{1}{1+\xi^{2}}\left(\begin{array}[]{c}\xi\\\ \xi^{2}\\\
\xi\\\ 1\end{array}\right)\,,$ (17)

which cannot be factorized. The transition operator for this case is

$P=\left(\begin{array}[]{cccc}0&1&0&0\\\ 1&0&0&0\\\ 0&0&1&0\\\
0&0&0&1\end{array}\right)\,,$ (18)

and it is nothing but a CNOT gate performing the operation

$|00\rangle+|01\rangle\rightarrow|00\rangle+|11\rangle\,.$ (19)

##  3 Discrete Wigner function

To gain further insights into the coherent states $|\xi_{\pm}\rangle$ we
proceed to picture them in phase space. To this end, we first note that, while
in the continuous case it is possible to translate a state by an infinite
distance, this is clearly not possible if the space is finite. To “prevent” a
state from “escaping” the finite phase space it is natural and convenient to
use the field $\mathrm{GF}(2^{n})$ in the representation of the states.

In consequence, we denote by $|\alpha\rangle$, with
$\alpha\in\mathrm{GF}(2^{n})$, an orthonormal basis in the Hilbert space of
the system. Operationally, the elements of the basis can be labeled by powers
of a primitive element, and the basis reads

$\\{|0\rangle,\,|\sigma\rangle,\ldots,\,|\sigma^{2^{n}-1}=1\rangle\\}\,.$ (20)

These vectors are eigenvectors of the operators $Z_{\beta}$ belonging to the
generalized Pauli group, whose generators are now defined as

$Z_{\beta}=\sum_{\alpha\in\mathrm{GF}(2^{n})}\chi(\alpha\beta)\,|\alpha\rangle\langle\alpha|\,,\qquad\qquad
X_{\beta}=\sum_{\alpha\in\mathrm{GF}(2^{n})}|\alpha+\beta\rangle\langle\alpha|\,,\\\
$ (21)

so that

$Z_{\alpha}X_{\beta}=\chi(\alpha\beta)\,X_{\beta}Z_{\alpha}\,.$ (22)

The operators (21) can be factorized into tensor products of powers of single-
particle Pauli operators $\sigma_{z}$ and $\sigma_{x}$, whose expression in
the standard basis of the two-dimensional Hilbert space is

$\hat{\sigma}_{z}=|1\rangle\langle 1|-|0\rangle\langle
0|\,,\qquad\qquad\hat{\sigma}_{x}=|0\rangle\langle 1|+|1\rangle\langle 0|\,.$
(23)

This factorization can be carried out by mapping each element of
$\mathrm{GF}(2^{n})$ onto an ordered set of natural numbers. As we have
already seen, a convenient choice for this is the self-dual basis, since the
finite Fourier transform factorizes then into a product of single-particle
Fourier operators, which leads to

$Z_{\alpha}=\hat{\sigma}_{z}^{a_{1}}\otimes\ldots\otimes\hat{\sigma}_{z}^{a_{n}}\,,\qquad\qquad
X_{\beta}=\hat{\sigma}_{x}^{b_{1}}\otimes\ldots\otimes\hat{\sigma}_{x}^{b_{n}}\,,$
(24)

where $(a_{1},\ldots,a_{n})$ and $(b_{1},\ldots,b_{n})$ are the expansion
coefficients of $\alpha$ and $\beta$, respectively, in the self-dual basis.

It was shown that the operators

$D(\alpha,\beta)=\phi(\alpha,\beta)\,Z_{\alpha}X_{\beta}\,,$ (25)

where $\phi(\alpha,\beta)$ is a phase, form an operational basis in the
discrete phase space.[10] The unitarity condition imposes the condition
$\phi^{2}(\alpha,\beta)=\chi(-\alpha\beta)$. These displacement operators (or
phase-point operators in the notation of Wootters[11]) allows us to introduce
a Hermitian kernel

$\Delta(\alpha,\beta)=\frac{1}{2^{n}}\sum_{\mu,\nu\in\mathrm{GF}(2^{n})}\chi(\alpha\nu-\beta\mu)\,D(\mu,\nu)\,,$
(26)

in terms of which we can define a well-behaved Wigner function as

$W_{\varrho}(\alpha,\beta)=\mathop{\mathrm{Tr}}\nolimits[\varrho\,\Delta(\alpha,\beta)]\,,$
(27)

where $\varrho$ is the density matrix of the system.

Figure 1: Wigner function for a coherent state of a system of three identical
qubits.

After some calculations, the Wigner function for our coherent states turns out
to be

$W_{|\xi_{+}\rangle}(\alpha,\beta)=\frac{1}{2^{n}}\frac{1}{(1+|\xi|^{2})^{n}}\sum_{\mu,\nu,\gamma\in\mathrm{GF}(2^{n})}\xi^{h(\gamma)}\overline{\xi}^{h(\gamma+\nu)}\chi(\alpha\nu+\beta\mu+\mu\nu+\mu\gamma)\phi(\mu,\nu)\,.$
(28)

A plot of this function for the case of three qubits is shown in Fig. 1. We
also note that the marginal distributions take a very simple form:

$\sum_{\alpha}W_{|\xi_{+}\rangle}(\alpha,\beta)=\frac{|\xi^{h(\beta)}|^{2}}{(1+|\xi|^{2})^{n}},\qquad\qquad\sum_{\beta}W_{|\xi_{+}\rangle}(\alpha,\beta)=\frac{|\xi^{h(\alpha)}|^{2}}{(1+|\xi|^{2})^{n}}\,.$
(29)

To conclude we wish to mention that it is also possible to introduce the
notion of squeezing for these states[12]. In the basis of the eigenstates of
$Z_{\alpha}$, such an operator has the following form

$S_{\lambda}=\sum_{\kappa\in\mathrm{GF}(2^{n})}|\kappa\rangle\langle\lambda\kappa|\,.$
(30)

The following relations hold

$S_{\lambda}^{\dagger}Z_{\alpha}S_{\lambda}=Z_{\alpha\lambda^{-1}}\,,\qquad\qquad
S_{\lambda}^{\dagger}X_{\alpha}S_{\lambda}=X_{\alpha\lambda}\,,$ (31)

so that

$\langle\xi_{\pm}|S_{\lambda}^{\dagger}X_{\alpha}S_{\lambda}|\xi_{\pm}\rangle=\langle\xi_{\pm}|S_{\lambda}^{\dagger}Z_{\alpha\lambda^{2}}S_{\lambda}|\xi_{\pm}\rangle\,.$
(32)

##  4 Conclusions

In summary, we have formulated a new sensible approach to deal with coherent
states for a system of $n$ qubits. The associated discrete Wigner function has
also been worked out. Some related problems, as the behavior under time
evolution or the extension to systems of qudits, will be addressed elsewhere.

##  Appendix A Galois fields

We briefly recall the minimum background of finite fields needed to proceed
through this paper. The reader interested in more mathematical details is
referred, e.g., to the excellent monograph by Lidl and Niederreiter.[13]

A commutative ring is a set $R$ equipped with two binary operations, called
addition and multiplication, such that it is an Abelian group with respect the
addition, and the multiplication is associative. Perhaps, the motivating
example is the ring of integers $\mathbb{Z}$ with the standard sum and
multiplication. On the other hand, the simplest example of a finite ring is
the set $\mathbb{Z}_{n}$ of integers modulo $n$, which has exactly $n$
elements.

A field $F$ is a commutative ring with division, that is, such that 0 does not
equal 1 and all elements of $F$ except 0 have a multiplicative inverse (note
that 0 and 1 here stand for the identity elements for the addition and
multiplication, respectively, which may differ from the familiar real numbers
0 and 1). Elements of a field form Abelian groups with respect to addition and
multiplication (in this latter case, the zero element is excluded).

The characteristic of a finite field is the smallest integer $p$ such that

$p\,1=\underbrace{1+1+\ldots+1}_{\mbox{\scriptsize$p$ times}}=0$ (33)

and it is always a prime number. Any finite field contains a prime subfield
$\mathbb{Z}_{p}$ and has $d=p^{n}$ elements, where $n$ is a natural number.
Moreover, the finite field containing $p^{n}$ elements is unique and is called
the Galois field $\mathrm{GF}(p^{n})$.

Let us denote as $\mathbb{Z}_{p}[x]$ the ring of polynomials with coefficients
in $\mathbb{Z}_{p}$. Let $P(x)$ be an irreducible polynomial of degree $n$
(i.e., one that cannot be factorized over $\mathbb{Z}_{p}$). Then, the
quotient space $\mathbb{Z}_{p}[X]/P(x)$ provides an adequate representation of
$\mathrm{GF}(p^{n})$. Its elements can be written as polynomials that are
defined modulo the irreducible polynomial $P(x)$. The multiplicative group of
$\mathrm{GF}(p^{n})$ is cyclic and its generator is called a primitive element
of the field.

As a simple example of a nonprime field, we consider the polynomial
$x^{2}+x+1=0$, which is irreducible in $\mathbb{Z}_{2}$. If $\sigma$ is a root
of this polynomial, the elements
$\\{0,1,\sigma,\sigma^{2}=\sigma+1=\sigma^{-1}\\}$ form the finite field
$\mathrm{GF}(2^{2})$ and $\sigma$ is a primitive element.

A basic map is the trace

$\mathop{\mathrm{tr}}\nolimits(\alpha)=\alpha+\alpha^{2}+\ldots+\alpha^{p^{n-1}}\,,$
(34)

which satisfies

$\mathop{\mathrm{tr}}\nolimits(\alpha+\beta)=\mathop{\mathrm{tr}}\nolimits(\alpha)+\mathop{\mathrm{tr}}\nolimits(\beta),$
(35)

and leaves the prime field invariant. In terms of it we define the additive
characters as

$\chi(\alpha)=\exp\left[\frac{2\pi
i}{p}\mathop{\mathrm{tr}}\nolimits(\alpha)\right],$ (36)

and posses two important properties:

$\chi(\alpha+\beta)=\chi(\alpha)\chi(\beta),\qquad\qquad\sum_{\alpha\in\mathrm{GF}(p^{n})}\chi(\alpha\beta)=p^{n}\delta_{0,\beta}.$
(37)

Any finite field $\mathrm{GF}(p^{n})$ can be also considered as an
$n$-dimensional linear vector space. Given a basis $\\{\theta_{k}\\}$,
($k=1,\ldots,n$) in this vector space, any field element can be represented as

$\alpha=\sum_{k=1}^{n}a_{k}\,\theta_{k},$ (38)

with $a_{k}\in\mathbb{Z}_{p}$. In this way, we map each element of
$\mathrm{GF}(p^{n})$ onto an ordered set of natural numbers
$\alpha\Leftrightarrow(a_{1},\ldots,a_{n})$.

Two bases $\\{\theta_{1},\ldots,\theta_{n}\\}$ and
$\\{\theta_{1}^{\prime},\ldots,\theta_{n}^{\prime}\\}$ are dual when

$\mathop{\mathrm{tr}}\nolimits(\theta_{k}\theta_{l}^{\prime})=\delta_{k,l}.$
(39)

A basis that is dual to itself is called self-dual.

There are several natural bases in $\mathrm{GF}(p^{n})$. One is the polynomial
basis, defined as

$\\{1,\sigma,\sigma^{2},\ldots,\sigma^{n-1}\\},$ (40)

where $\sigma$ is a primitive element. An alternative is the normal basis,
constituted of

$\\{\sigma,\sigma^{p},\ldots,\sigma^{p^{n-1}}\\}.$ (41)

The choice of the appropriate basis depends on the specific problem at hand.
For example, in $\mathrm{GF}(2^{2})$ the elements $\\{\sigma,\sigma^{2}\\}$
are both roots of the irreducible polynomial. The polynomial basis is
$\\{1,\sigma\\}$ and its dual is $\\{\sigma^{2},1\\}$, while the normal basis
$\\{\sigma,\sigma^{2}\\}$ is self-dual.

## References

  * [1] H. Weyl, The Theory of Groups and Quantum Mechanics (Dover, New York, 1950).
  * [2] J. Schwinger, Unitary operator basis, Proc. Natl. Acad. Sci. USA 46 (1960) 570–576.
  * [3] A. Vourdas, Quantum systems with finite Hilbert space, Rep. Prog. Phys. 67 (2004) 267–320.
  * [4] A. Vourdas, Quantum systems with finite Hilbert space: Galois fields in quantum mechanics, J. Phys. A 40 (2007) R285–R331.
  * [5] M. Ruzzi, Jacobi $\vartheta$ functions and discrete Fourier transforms, J. Math. Phys. 47 (2006) 063507.
  * [6] A. Perelomov, Generalized Coherent States and their Applications (Springer, Berlin, 1986).
  * [7] D. Galetti and M. A. Marchiolli, Discrete coherent states and probability distributions in finite-dimensional spaces, Ann. Phys. (N Y) 249 (1996) 454–480.
  * [8] M. Ruzzi, M. Marchiolli and D. Galetti, Extended Cahill-Glauber formalism for finite-dimensional spaces: I. fundamentals, J. Phys. A 38 (2005) 6239–6245.
  * [9] M. L. Mehta, Eigenvalues and eigenvectors of the finite Fourier transform, J. Math. Phys. 28 (1987) 781–785.
  * [10] G. Björk, A. B. Klimov and L. L. Sánchez-Soto, The discrete Wigner function, Prog. Opt. 51 (2008) 470–514
  * [11] W. K. Wootters, A Wigner-function formulation of finite-state quantum mechanics, Ann. Phys. (N Y) 176 (1987) 1–21.
  * [12] M. A. Marchiolli, M. Ruzzi and D. Galetti, Discrete squeezed states for finite-dimensional spaces, Phys. Rev. A 76 (2007) 032102.
  * [13] R. Lidl, H. Niederreiter, Introduction to finite fields and their applications (Cambridge University Press, Cambridge, 1986).


# Scattering of matter-wave soliton from a narrow barrier: transmission
coefficient and induced collapse

F. Lorenzi  Dipartimento di Fisica e Astronomia "Galileo Galilei", Università
di Padova, Via Marzolo 8, 35131 Padova, Italy  Istituto Nazionale di Fisica
Nucleare (INFN), Sezione di Padova, via Marzolo 8, 35131 Padova, Italy
<EMAIL_ADDRESS>L. Salasnich  Dipartimento di Fisica e
Astronomia "Galileo Galilei", Università di Padova, Via Marzolo 8, 35131
Padova, Italy  Istituto Nazionale di Fisica Nucleare (INFN), Sezione di
Padova, via Marzolo 8, 35131 Padova, Italy  Padua Quantum Technology Research
Center, Università di Padova, Via Gradenigo 6/A, 35131 Padova, Italy  Istituto
Nazionale di Ottica (INO) del Consiglio Nazionale delle Ricerche (CNR), via
Nello Carrara 1, 50019 Sesto Fiorentino, Italy

###### Abstract

We report systematic numerical simulations of the collision of a bright
matter-wave soliton made of Bose-condensed alkali-metal atoms through a narrow
potential barrier by using the three-dimensional Gross-Pitaevskii equation. In
this way, we determine how the transmission coefficient depends on the soliton
impact velocity and the barrier height. Quite remarkably, we also obtain the
regions of parameters where there is the collapse of the bright soliton
induced by the collision. We compare these three-dimensional results with the
ones obtained by three different one-dimensional nonlinear Schrödinger
equations. We find that a specifically modified nonpolynomial Schrödinger
equation is able to capture all the main features of the three-dimensional
findings. In particular, this simplified but very effective one-dimensional
model takes into account the transverse width dynamics of the soliton with an
ordinary differential equation coupled to the partial differential equation of
the axial wavefunction of the Bose-Einstein condensate.

Localized soliton-like structures, known as bright matter-wave solitons, can
be generated in Bose-Einstein Condensates (BEC) with attractive interatomic
interactions. Since the first experimental realization of such structure about
two decades ago [2], dynamics of matter-wave bright solitons in an attractive
BEC have been intensely studied, both at the quantum level and using the
Gross-Pitaevskii equation (GPE). Current experimental capabilities offer an
unprecedented opportunity to test many-body theories in an ultracold Bose gas,
and matter-wave solitons are an excellent target for the predictions [3, 4].
Moreover, several technological applications are based on the possibility of
generating and manipulating such kind of coherent structures: some of them are
interferometry [5] even beyond the quantum limit, and quantum-enhanced
metrology [6, 7]. The remarkable analogy of models based on the 3D-GPE with
the equations of motion occurring in optics with Kerr media allowed to study
common aspects on the same ground, such as the Hong-Ou-Mandel experiment [8].

In a typical setup of an atomic interferometer using solitons, a matter-wave
soliton is prepared in a quasi-1D trap, i.e. a confining potential made of a
strong radial component and a weak or absent axial component. This setup
experimentally allows the creation of cigar-shaped condensates in the case of
repulsive interparticle interaction and in the noninteracting case. By using
attractive interparticle interactions that are obtainable for example by using
Feshbach resonances, one can generate a matter-wave soliton. In the latter
case, solitons loaded into quasi-1D traps have typical axial widths that are
comparable to the radial potential characteristic length . To achieve
interference, the soliton is set into motion by phase imprinting, and it
collides with a narrow potential barrier set by a narrow laser beam, acting
analogously to an optical beam splitter, designed to be able to split the
number of atoms into two even solitonic packets. The resulting two solitons
are then recombined in a later stage, through the same barrier. After the
first splitting, split solitons may achieve a differential phase shift, thus
allowing the observation of interference in the recombined packet.

From the theoretical point of view, the study of quantum matter-wave solitons
was carried out mostly in 1D, where the many-body wave function is well known
to have an exact solution by Bethe ansatz [9]. This procedure relies on having
very strong radial confinement, and it is not sensible to the 3D dynamics,
thus losing details about the transverse degrees of freedom that are
especially important near the point of GPE collapse in which the reduction of
the size of the condensate brings it to a regime, in which other interaction
effects start to be non-negligible, like three-body interactions causing
depletion from the trap [10].

The stability of solitons in quasi-1D harmonic traps is highly nontrivial, as
the cubic 1D-GPE has a critical dimension $2$. Instability is in the form of a
collapse, also known in the mathematical literature as nonlinear blow-up [11].
In one dimension, the collapse is prevented by the Vlasov-Petrishev-Talanov
theorem [11]. GPE collapse due to an arbitrary attractive interaction
potential can be triggered in this context by loading into the trap a suitably
high number of particles [12]. Moreover, the radial anisotropy of the trap can
play a role in the critical number of particles for collapse [13, 14]. It is
fundamental to remark that, even with a purely 1D model, by adding a barrier-
like external potential to the 1D-GPE the problem becomes non-integrable, and
requires approximate methods to be tackled.

The 3D-GPE dynamics represents a useful tool not only as an approximation of
the full quantum dynamics but also as providing signatures of soliton
entanglement across the barrier [15]. In fact, the discontinuity in the
reflection coefficient indicates the possibility of creating “Schrödinger cat"
states, exploring quantum entanglement phenomena [15]. The 3D-GPE model
predicts a peculiar behavior of the transmission coefficient with the barrier:
at low values of the velocity and the barrier height it is a discontinuous
function of the parameters. This aspect has been investigated [15] and is
frequently referred to as the particle behavior of the impinging soliton.

Various dimensional reduction schemes for the 3D-GPE have been proposed, [3,
16, 17, 18]. In this work, we compare three schemes of dimensional reduction
with full 3D simulations. 1D effective equations have a great computational
advantage in the description of the dynamics and are routinely used in studies
of atomic interferometers. The simplest one is the 1D-GPE, obtained by
imposing a fixed transverse wavefunction as the eigenstate of the transverse
harmonic potential. An improved model is called nonpolynomial Schrödinger
equation (NPSE) [18], which is based on assuming the transverse width of the
trial as a variational parameter and obtaining the equation of motion as
Euler-Lagrange equations. Furthermore, in the original NPSE formulation,
derivatives of the transverse width parameter present in the Lagrangian are
approximated to zero, thus obtaining an algebraic equation of motion for the
transverse width parameter. So we also consider the non-approximated version
of the NPSE, which we call NPSE+ for brevity. Previous work [19] highlighted
the behavior of the transmission coefficient and the barrier-induced collapse
with the 3D-GPE and the NPSE, studying regions in the barrier strength vs.
number of particles plane. We investigate the differences in the description
of the collision process, focusing on the transmission coefficient and the
onset of barrier-induced collapse.

### Gross-Pitaevskii equations

The model is based on the Hartree approximation for bosons [9], using which it
is possible to derive the following Lagrangian, called Gross-Pitaevskii
Lagrangian, for the field $\psi({\bf r},t)$, representing the wavefunction of
the Hartree product state with all the particles in the same single-particle
quantum state,

$\mathcal{L}=\int d^{3}\mathbf{r}\
\psi^{*}\left[i\hbar\frac{\partial}{\partial
t}+\frac{\hbar^{2}}{2m}\nabla^{2}-U-\frac{g}{2}(N-1)|\psi|^{2}\right]\psi,$
(1)

where $U$ is the external potential, $N$ is the number of particles, and $g$
is the contact potential, which can be linked to the s-wave scattering length
$a_{s}$ with the expression

$g=\frac{4\pi\hbar^{2}a_{s}}{m}.$ (2)

The associated Euler-Lagrange (EL) equation is the 3D-GPE:

$i\hbar\dfrac{\partial}{\partial
t}\psi=\left[-\dfrac{\hbar^{2}}{2m}\nabla^{2}+U+g(N-1)|\psi|^{2}\right]\psi.$
(3)

Standard dimensional reduction of the 3D-GPE in a tight transverse harmonic
potential relies on the assumption that the transverse degree of freedom of
the wavefunction is frozen to the ground state of the harmonic potential, as
we will briefly review now. Let the external potential be written as

$U(x,y,z)=\frac{1}{2}m\omega_{\perp}^{2}(y^{2}+z^{2})+V(x),$ (4)

where $\omega_{\perp}$ is the isotropic strength of the potential, and $V$ is
a weaker axial potential. The role of anisotropy on the transverse potential
was studied in [14, 13]. The harmonic potential strength $\omega_{\perp}$
naturally sets a characteristic length scale
$l_{\perp}=\sqrt{\hbar/(m\omega_{\perp})}$. Let us assume that the
wavefunction is composed of a constant Gaussian transverse part $\phi$, which
is the ground state of the transverse harmonic potential, and a time-varying
axial component $f$ as

$\psi(\mathbf{r},t)=f(x,t)\phi(y,z),$ (5)

where

$\phi(y,z)=\dfrac{1}{\sqrt{\pi}l_{\perp}}\exp\left[-\frac{y^{2}+z^{2}}{2l_{\perp}^{2}}\right].$
(6)

This is physically justified when the interaction energy is much smaller than
the energy difference between the first excited state and the ground state of
the transverse potential. Inserting this ansatz into Eq. (3), and integrating
along the transverse coordinates, one obtains the corresponding wave equation,
called the 1D-GPE:

$i\hbar\dfrac{\partial}{\partial
t}f=\left[-\dfrac{\hbar^{2}}{2m}\dfrac{\partial^{2}}{\partial
x^{2}}+V(z)+\hbar\omega_{\perp}+g_{1D}|f|^{2}\right]f$ (7)

where we defined $g_{1D}=g(N-1)/(2\pi l_{\perp})$. The expression of $g_{1D}$
can be used to define a normalized nonlinear parameter
$\gamma=(N-1)|a_{s}|/l_{\perp}$. Using Eq. (7), for $g_{1D}<0$ the ground
states of the axial problem are constituted by stable solitons, and no
collapse is expected for any interaction strength. Instead, for the 3D-GPE
case, there exists a critical nonlinear parameter for the existence of stable
solitons. The value, above which static wavefunction collapse is expected, is
about $\gamma_{c}\approx 0.67$.

(a) Soliton solutions compared. (b) Zoom applied to the top part of the
soliton solutions.

Figure 1: Comparison of the ground state of the 3D-GPE, the 1D-GPE, the NPSE,
and the NPSE with derivatives. $x$ is in units of
$l_{\bot}=\sqrt{\hbar/(m\omega_{\bot})}$, the characteristic length of
transverse harmonic confinement of frequency $\omega_{\bot}$. The nonlinear
parameter is set to $\gamma=0.65$. The three-dimensional Gross-Pitaevskii
equation is the red solid line, the nonpolynomial Schödinger equation without
the corrections is the green dotted line, the one with the corrections is the
green dashed line. The one-dimensional Gross-Pitaevskii equation is in solid
grey.

### Variational ansatz and NPSE

As shown in [18], a better approximation is to consider the separation of the
total wavefunction in a transverse Gaussian component with non-constant
transverse width $\sigma(x,t)$ and to find the equation of motion using a
variational principle. The resulting equation is the NPSE. The function $\phi$
in the ansatz Eq. (5) is substituted by a more general

$\phi(y,z,\sigma(x,t))=\dfrac{1}{\sqrt{\pi}\sigma(x,t)}\exp\left[-\frac{y^{2}+z^{2}}{2\sigma(x,t)^{2}}\right]$
(8)

where $\sigma$ is a function to be determined as a variational parameter. In
[18], the calculations were done assuming that derivatives of $\sigma$ are
negligible. By keeping these derivatives terms, it is possible to write the
following effective 1D Lagrangian (detailed calculations are given in the
Methods section)

$\mathcal{L}=\int dx\ f^{*}\bigg{[}i\hbar\frac{\partial}{\partial
t}+\frac{\hbar^{2}}{2m}\frac{\partial^{2}}{\partial
x^{2}}-V-\frac{\hbar^{2}}{2m}\frac{1}{\sigma^{2}}\left(1+\left(\frac{\partial}{\partial
x}\sigma\right)^{2}\right)-\frac{m\omega_{\perp}}{2}\sigma^{2}-\frac{2\hbar^{2}a_{s}(N-1)}{m\sigma^{2}}|f|^{2}\bigg{]}f.$
(9)

The corresponding EL equations for $f$ and $\sigma$ provide the solution to
the variational problem and are obtained as shown in [20]:

$i\hbar\frac{\partial}{\partial
t}f=\bigg{[}-\frac{\hbar^{2}}{2m}\frac{\partial^{2}}{\partial
x^{2}}+V+\frac{\hbar^{2}}{2m}\frac{1}{\sigma^{2}}\left(1+\left(\frac{\partial}{\partial
x}\sigma\right)^{2}\right)+\frac{m\omega_{\perp}}{2}\sigma^{2}+\frac{2\hbar^{2}a_{s}(N-1)}{m\sigma^{2}}|f|^{2}\bigg{]},$
(10)
$\sigma^{4}-l_{\perp}^{4}\left[1+2a_{s}|f|^{2}\right]+l_{\perp}^{4}\left[\sigma\frac{\partial^{2}}{\partial
x^{2}}\sigma-\left(\frac{\partial}{\partial x}\sigma\right)^{2}\right]=0.$
(11)

We will refer to the above coupled equations as NPSE+. By neglecting the
derivative of $\sigma$, one can obtain the same effective 1D Lagrangian as in
Ref. [18], whose EL equations correspond to

$\begin{split}i\hbar\frac{\partial}{\partial
t}f=\bigg{[}-\frac{\hbar^{2}}{2m}\frac{\partial^{2}}{\partial
x^{2}}+V+\frac{\hbar^{2}}{2m}\frac{1}{\sigma^{2}}+\frac{m\omega_{\perp}}{2}\sigma^{2}+\frac{2\hbar^{2}a_{s}(N-1)}{m\sigma^{2}}|f|^{2}\bigg{]}f,\end{split}$
(12) $\sigma^{2}=l_{\perp}^{2}\sqrt{1+2a_{s}(N-1)|f|^{2}}.$ (13)

Eq. (12) is called NPSE. We remark that the NPSE+, as opposed to the NPSE,
respects the variational principle, so the corresponding ground state energy
is bound to be greater or equal than the true ground state energy of the
3D-GPE.

We will use the 3D-GPE as a reference equation, and compare the predictions of
the 1D-GPE, the NPSE, and the NPSE+. The axial densities of the ground state
solutions are shown in Fig. 1, where we have set the nonlinear parameter to a
very high value $\gamma=0.65$, near the 3D-GPE static collapse value
$\gamma_{c}$. We notice that, as the 1D-GPE fails to represent accurately the
axial wavefunction, NPSE+ is more accurate with respect to the NPSE, but has
the additional cost of the solution of the transverse width differential
equation coupled to the axial wavefunction partial differential equation.

### Generalization of bound on splitting energy

The soliton splitting event can be verified only for specific ranges of the
transmission coefficient [21], depending on the initial soliton velocity.
These values can be computed by analyzing energy conservation in the splitting
event. The interplay of kinetic energy and internal energy of the solitons
during the scattering event has been discussed in [22, 23, 21], by using Lieb-
Liniger[9] energies $E_{G}$ pertaining to the soliton internal degrees of
freedom in the total Hamiltonian. Imposing energy conservation, the kinetic
energy of the initial soliton must satisfy

$E_{k}>E_{G}\left(N-n\right)+E_{G}\left(n\right)-E_{G}(N),$ (14)

where $N$ is the number of atoms in the initial soliton, and $n$ is the one in
the transmitted soliton. This is the condition that must be satisfied for
having a splitting event of transmission coefficient $T=n/N$. In our case, the
internal energy of the soliton can be computed either numerically or
analytically. The chemical potential of the stationary solution for the NPSE
has been obtained [18] as satisfying the implicit relation

$(1-\mu)^{3/2}-\dfrac{3}{2}(1-\mu)^{1/2}\dfrac{3}{2\sqrt{2}}\gamma=0$ (15)

where $\gamma=(N-1)|a_{s}|/l_{\perp}$, and selecting only the stable branch of
the solutions, i.e. the one satisfying the Vakhitov-Kolokolov criterion
$\frac{\partial}{\partial n}\mu<0$. For the other equations used in this work,
it is possible to obtain numerically the value of $\mu(n)$ from stationary
state solutions. Using the values of the chemical potential, we are able to
write the ground state energy of the nonlinear wave equation corresponding to
an $N$-particle soliton as

$E_{G}(n)=\int_{0}^{n}dn^{\prime}\mu(n^{\prime}),$ (16)

it is possible to obtain different ranges of transmission coefficients that
are accessible for a given value of the initial kinetic energy.

##  1 Results

Figure 2: Chemical potential dependence on the nonlinear parameter.
$\gamma=(N-1)|a_{s}|/l_{\perp}$. $\mu$ is in units of $\hbar\omega_{\bot}$.The
three-dimensional Gross-Pitaevskii equation is the red solid line, the
nonpolynomial Schödinger equation without the corrections is the green dotted
line, and the one with the corrections is the green dashed line. The one-
dimensional Gross-Pitaevskii equation is in solid grey.

Figure 3: Comparison of the transverse width parameter. Both $\sigma$ and $x$
are in units of $l_{\bot}=\sqrt{\hbar/(m\omega_{\bot})}$. The three-
dimensional Gross-Pitaevskii equation is the red solid line, the nonpolynomial
Schödinger equation without the corrections is the green dotted line, and the
one with the corrections is the green dashed line. The one-dimensional Gross-
Pitaevskii equation is in solid grey.

### Soliton solutions

We start to review some properties of the solitonic ground state of the
equations, comparing them. We study the highly nonlinear regime, in which
$\gamma=0.65$. Soliton solutions in this case are stable for the 3D-GPE and
the NPSE and NPSE+ for $\gamma<\gamma_{c}$ [17, 11]. Soliton solutions are
shown for all the equations in Fig. 1. The simulations show a better agreement
of the NPSE+ equation with respect to the 3D-GPE, whereas both the 1D-GPE and
the NPSE are inaccurate. The accuracy of the NPSE+ ground state is reflected
also in the transverse width, as shown in Fig. 3. The computation of the NPSE+
transverse width is obtained by iteratively solving Eq. (10) and then Eq.
(11). In the solution of the ordinary differential equation Eq. (11) we apply
Dirichlet boundary conditions corresponding to the vanishing of the axial
wavefunction at infinity. The computation of the 3D-GPE transverse width is
done by a least square fit on the radial distribution of the wavefunction,
namely defined as

$\sigma^{2}(x)=\int dydz\,(y^{2}+z^{2})|\psi(x,y,z)|^{2}.$ (17)

The results are shown in Fig. 3. The chemical potential of the solitons is
shown in Fig. 3. We remark that being the NPSE only an approximation of the
true variational solution, its chemical potential is allowed to be less than
the 3D-GPE chemical potential, thus becoming less than the bound set by the
variational principle. As expected, the difference becomes more pronounced at
high nonlinearities, implying a more localized wavefunction, where the terms
proportional to the derivatives of $\sigma$ in Eq. (10) and (11) become more
relevant.

### Scattering from a narrow barrier

We are interested in computing the transmission coefficient for various
velocities and barriers. In our simulations, the velocity ranges in
$v\in[0.1,1.0]$, and the barrier in $b\in[0.0,1.0]$. In particular, by setting
a sufficiently high $\gamma$, for example near to $\gamma_{c}$, namely
$\gamma=0.65$, we analyze the onset of barrier-induced collapse, happening for
high soliton velocity and high barrier energy. Our results show that the
vanishing of the transverse width predicted by the NPSE, and suggesting a
barrier-induced collapse, is a very weak indicator of an actual collapse, as
can be seen by comparing Fig.4 (a) and (c) panels. Instead, the NPSE+
collapsing region is not due to a vanishing of the, but to a sudden
concentration of the axial density in smaller and smaller regions, like in the
3D-GPE case. In fact, we have set the numerical threshold of the collapse to
the detection of a single probability per site greater than $0.3$. The abrupt
change in the local maximum density we observe between stable solutions and
collapsing ones justifies the validity of this criterion. By comparing Fig.4
(a) and (b) panels, it is clear that NPSE+ represents a better effective
equation in the description of the barrier-induced collapse. We also reported
results obtained by the familiar 1D-GPE in Fig.4 panel (d). As expected, this
show no collapse region.

(a) 3D-GPE (b) NPSE+ (c) NPSE (d) 1D-GPE

Figure 4: Comparison of the transmission coefficient and the collapse regions
for the various equations. $v$ is in units of $\sqrt{\hbar\omega_{\perp}/m}$,
$b$ is in units of $\hbar\omega_{\perp}$.

Moreover, we compare the transmission coefficient versus the barrier strength
for two values of the impinging velocity, namely $v=0.5$ and $v=1$ in Fig.5.
Quite remarkably, the NPSE+ is able to describe in a very accurate way the
transmission coefficient in the non-collapsing region.

(a) $v=0.5$. (b) $v=1$.

Figure 5: Transmission coefficient versus barrier strength with fixed
velocity. $v$ is in units of $\sqrt{\hbar\omega_{\perp}/m}$, $b$ is in units
of $\hbar\omega_{\perp}$. The three-dimensional Gross-Pitaevskii equation is
the red solid line, the nonpolynomial Schödinger equation without the
corrections is the green dotted line, and the one with the corrections is the
green dashed line.

##  2 Discussion

In this article, we have presented a numerical study of the collision of a
bright matter-wave soliton with a narrow potential barrier using the three-
dimensional Gross-Pitaevskii and three dimensionally reduced versions of it.
We investigated how the choice of dimensional reduction impacts the
description of some features of the process, namely the transmission
coefficient and the onset of barrier-induced collapse, also using the familiar
one-dimensional Gross-Pitaevskii. We first reviewed the ground state
properties given by all the schemes, highlighting the role of the variational
transverse width. Then we compared the scattering properties: our results show
that by using the NPSE in a regime of high barrier energy and high velocity it
fails to describe the 3D dynamics due to the the vanishing of the transverse
width of the solution. In such cases the collapse phenomena in the 3D
solutions are absent, and the NPSE is not capturing the correct dynamics.

Our main result is that by adopting a slight modification of the transverse
width using the true variational solution with the NPSE ansatz, the efficacy
of the method is restored, showing good agreement with the 3D-GPE. We believe
the present work to be a valuable contribution to the field of matter-wave
soliton interferometry and quantum measurement, as results can be used in
predicting the dynamics of experimentally accessible scenarios. For example,
the NPSE+ can be used for modeling interferometric experiments in highly
nonlinear regimes. In the setting of a quasi 1D harmonic trap corresponding to
a transverse frequency $\omega_{\perp}/2\pi=254\text{Hz}$, loaded with about
$28\times 10^{3}$ atoms of 7Li, analogously to a past experiment [24], the
nonlinear regime we have studied is achieved for an s-wave scattering length
of $a_{s}\approx-5.52\times 10^{-11}\text{m}$. By using the 1D-GPE or the NPSE
as one-dimensional effective models, it is not possible to accurately describe
the 3D dynamics.

##  3 Methods

### Numerical method

The time-marching scheme we use for all the simulations is the SSFM, adopting
Strang splitting of the nonlinear and linear part of the evolution operator.
The SSFM is well-known to be accurate to the second order in time and to every
order in space, thus being highly efficient in the spatial discretization [25,
26]. The drawback of the method - or the feature, depending on the point of
view - is to natively implement PBC. The implementation of absorption
boundaries is still possible but not straightforward [27, 28]. We assume the
field to be localized away from the boundary in order to neglect this problem.
In our setup, we use a unit of energy of $\hbar\omega_{\perp}$, a unit of time
of $\omega_{\perp}^{-1}$ and a unit of length of $l_{\perp}$, constituting the
natural units for the (isotropic) harmonic confinement. In these units, we
consider for the 1D simulations a total length of $L=40$, with a grid of
$N=512$ points. In the 3D simulations, we use a grid of
$(N_{x},N_{y},N_{z})=(512,40,40)$ points, with total lengths of
$(L_{x},L_{y},L_{z})=(40,10,10)$. The time step in both setups is chosen to be
$h_{t}=0.01$. These parameters have been proven to give a total truncation
error in the $L_{\infty}$ norm of the order of $10^{-4}$, and allow for a
reasonable computation time of all the calculations.

The ground state solutions are computed using an imaginary-time propagation
method. We point out that some modifications of this method are available
under the name of normalized gradient-flow methods [29].

###  3.1 Derivation of the NPSE+

Following [17], we write the 3D Lagrangian

$\mathcal{L}=\int dx\int dy\,dz\
f^{*}\phi^{*}\left[i\hbar\frac{\partial}{\partial
t}+\right.\left.\frac{\hbar^{2}}{2m}\nabla^{2}-U-\frac{1}{2}g(N-1)|f\phi|^{2}\right]f\phi.$
(18)

We are interested in integrating along the transverse coordinates without
neglecting the terms proportional to $\frac{\partial}{\partial z}\sigma$ and
$\frac{\partial^{2}}{\partial z^{2}}\sigma$. The novel terms arise from
$i\hbar\frac{\partial}{\partial t}(f\phi)$ and
$\frac{\hbar}{2m}\nabla^{2}(f\phi)$. By separating the derivatives, we have

$\mathcal{L}=\int dz\int
dxdyf^{*}\phi^{*}\biggl{[}i\hbar\phi\frac{\partial}{\partial t}f+i\hbar
f\phi\left(\dfrac{y^{2}+z^{2}}{\sigma^{3}}-\dfrac{1}{\sigma}\right)\frac{\partial}{\partial
t}\sigma+\frac{\hbar^{2}}{2m}\left(f\nabla_{\perp}^{2}\phi+f\frac{\partial^{2}}{\partial
z^{2}}\phi+\phi\frac{\partial^{2}}{\partial
z^{2}}f\right)-Uf\phi-\frac{1}{2}g(N-1)|f\phi|^{2}f\phi\biggr{]}.$ (19)

Integrating the term proportional to $\frac{\partial}{\partial t}\sigma$ gives
$0$, as one may realize by looking at the symmetry of its prefactor. However,
the term proportional to $\frac{\partial^{2}}{\partial z^{2}}\phi$ gives a
non-null contribution to the 1D Lagrangian. The integration gives

$\mathcal{L}=\int dz\,f^{*}\bigg{[}i\hbar\frac{\partial}{\partial
t}+\frac{\hbar^{2}}{2m}\frac{\partial}{\partial
z}^{2}-V-\frac{\hbar^{2}}{2m\sigma^{2}}\left(1+\left(\frac{\partial}{\partial
z}\sigma\right)^{2}\right)-\frac{m\omega_{\perp}^{2}}{2}\sigma^{2}-\frac{1}{2}\frac{g(N-1)}{2\pi\sigma^{2}}|f|^{2}\bigg{]}f.$
(20)

Let us now consider the corresponding Euler-Lagrange equations. We will have a
simple corrective term in the equation for $f$, and $\sigma$

$i\hbar\frac{\partial}{\partial
t}f=\left[-\dfrac{\hbar^{2}}{2m}\frac{\partial}{\partial
z}^{2}+V+\dfrac{\hbar^{2}}{2m}\dfrac{1}{\sigma^{2}}+\dfrac{m\omega_{\perp}^{2}}{2}\sigma^{2}+\dfrac{N-1}{2\pi\sigma^{2}}g|f|^{2}\right]f$
(21)

$-m\omega_{\perp}^{2}\sigma+\left[\dfrac{\hbar^{2}}{m}+\dfrac{N-1}{2\pi}g|f|^{2}\right]\sigma^{-3}+\frac{\hbar^{2}}{m}\sigma^{-3}\left(\sigma\frac{\partial}{\partial
z}^{2}\sigma-\left(\frac{\partial}{\partial z}\sigma\right)^{2}\right)=0.$
(22)

We obtain a differential equation for $\sigma$ that can be turned into an
algebraic equation as in the standard derivation of the NPSE [18] by
neglecting the derivatives of $\sigma$.

## Acknowledgements

F.L. and L.S. acknowledge a National Grant of the Italian Ministry of
University and Research for the PRIN 2002 project "Quantum Atomic Mixtures:
Droplets, Topological Structures, and Vortices". L.S. is partially supported
by the BIRD grant “Ultracold atoms in curved geometries” of the University of
Padova, by the “Iniziativa Specifica Quantum” of INFN, by the European Quantum
Flagship project PASQuanS 2, and by the European Union-NextGenerationEU within
the National Center for HPC, Big Data and Quantum Computing (Project No.
CN00000013, CN1 Spoke 10: “Quantum Computing”).

## References

  * [1]
  * [2] Khaykovich, L. _et al._ Formation of a matter-wave bright soliton.  _Science_ 296, 1290–1293 (2002).
  * [3] Weiss, C. & Castin, Y.  Creation and detection of a mesoscopic gas in a nonlocal quantum superposition.  _Phys. Rev. Lett._ 102, 010403 (2009).
  * [4] Streltsov, A. I., Alon, O. E. & Cederbaum, L. S.  Scattering of an attractive bose-einstein condensate from a barrier: Formation of quantum superposition states.  _Phys. Rev. A_ 80, 043616 (2009).
  * [5] Helm, J. L., Cornish, S. L. & Gardiner, S. A.  Sagnac interferometry using bright matter-wave solitons.  _Phys. Rev. Lett._ 114, 134101 (2015).
  * [6] Dunningham, J. A. & Burnett, K.  Sub-shot-noise-limited measurements with bose-einstein condensates.  _Phys. Rev. A_ 70, 033601 (2004).
  * [7] Dunningham, J. A.  Using quantum theory to improve measurement precision.  _Contemp. Phys._ 47, 257–267 (2006).
  * [8] Sun, Z.-Y., Kevrekidis, P. G. & Krüger, P.  Mean-field analog of the hong-ou-mandel experiment with bright solitons.  _Phys. Rev. A_ 90, 063612 (2014).
  * [9] Lieb, E. H. & Liniger, W.  Exact analysis of an interacting bose gas. i. the general solution and the ground state.  _Phys. Rev._ 130, 1605–1616 (1963).
  * [10] Cornish, S. L., Thompson, S. T. & Wieman, C. E.  Formation of bright matter-wave solitons during the collapse of attractive bose-einstein condensates.  _Phys. Rev. Lett._ 96, 170401 (2006).
  * [11] Sulem, C. & Sulem, P.-L.  _The Nonlinear Schrödinger Equation: Self-Focusing and Wave Collapse_.
  * [12] Bergé, L.  Wave collapse in physics: principles and applications to light and plasma waves.  _Phys. Rep._ 303, 259–370 (1998).
  * [13] Mazzarella, G. & Salasnich, L.  Collapse of triaxial bright solitons in atomic bose–einstein condensates.  _Phys. Lett. A_ 373, 4434–4437 (2009).
  * [14] Gammal, A., Tomio, L. & Frederico, T.  Critical numbers of attractive bose-einstein condensed atoms in asymmetric traps.  _Phys. Rev. A_ 66, 043619 (2002).
  * [15] Gertjerenken, B. _et al._ Generating mesoscopic bell states via collisions of distinguishable quantum bright solitons.  _Phys. Rev. Lett._ 111, 100406 (2013).
  * [16] Khaykovich, L. & Malomed, B. A.  Deviation from one dimensionality in stationary properties and collisional dynamics of matter-wave solitons.  _Phys. Rev. A_ 74, 023607 (2006).
  * [17] Salasnich, L., Parola, A. & Reatto, L.  Condensate bright solitons under transverse confinement.  _Phys. Rev. A_ 66, 043603 (2002).
  * [18] Salasnich, L., Parola, A. & Reatto, L.  Effective wave equations for the dynamics of cigar-shaped and disk-shaped bose condensates.  _Phys. Rev. A_ 65, 043614 (2002).
  * [19] Cuevas, J., Kevrekidis, P. G., Malomed, B. A., Dyke, P. & Hulet, R. G.  Interactions of solitons with a gaussian barrier: splitting and recombination in quasi-one-dimensional and three-dimensional settings.  _New J. Phys._ 15, 063006 (2013).
  * [20] Salasnich, L., Malomed, B. A. & Toigo, F.  Matter-wave vortices in cigar-shaped and toroidal waveguides.  _Phys. Rev. A_ 76, 063614 (2007).
  * [21] Gertjerenken, B., Billam, T. P., Khaykovich, L. & Weiss, C.  Scattering bright solitons: Quantum versus mean-field behavior.  _Phys. Rev. A_ 86, 033608 (2012).
  * [22] Helm, J. L., Rooney, S. J., Weiss, C. & Gardiner, S. A.  Splitting bright matter-wave solitons on narrow potential barriers: Quantum to classical transition and applications to interferometry.  _Phys. Rev. A_ 89, 033610 (2014).
  * [23] Wang, C.-H., Hong, T.-M., Lee, R.-K. & Wang, D.-W.  Particle-wave duality in quantum tunneling of a bright soliton.  _Opt. Expr._ 20, 22675–22682 (2012).
  * [24] Nguyen, J. H., Dyke, P., Luo, D., Malomed, B. A. & Hulet, R. G.  Collisions of matter-wave solitons.  _Nature Physics_ 10, 918–922 (2014).
  * [25] Taha, T. R. & Ablowitz, M. J.  Analytical and numerical aspects of certain nonlinear evolution equations. i. analytical.  _J. Comput. Phys._ 55, 192–202 (1984).
  * [26] Taha, T. R. & Ablowitz, M. I.  Analytical and numerical aspects of certain nonlinear evolution equations. ii. numerical, nonlinear schrödinger equation.  _J. Comput. Phys._ 55, 203–230 (1984).
  * [27] Antoine, X., Bao, W. & Besse, C.  Computational methods for the dynamics of the nonlinear schrödinger/gross–pitaevskii equations.  _Comput. Phys. Comm._ 184, 2621–2633 (2013).
  * [28] Antoine, X., Geuzaine, C. & Tang, Q.  Perfectly matched layer for computing the dynamics of nonlinear schrödinger equations by pseudospectral methods. application to rotating bose-einstein condensates.  _Commun. Nonlinear. Sci. Numer. Simulat._ 90, 105406 (2020).
  * [29] Bao, W. & Du, Q.  Computing the ground state solution of bose-einstein condensates by a normalized gradient flow.  _SIAM J. Sci. Comput._ 25, 1674–1697 (2004).