text:
# Entropy of chains placed on the square lattice

## I Introduction

The term dimer was introduced in the thirties fr37  as an abbreviation for diatomic molecules in a model for their adsorption on crystal surfaces. Later dimer models were applied in the study of other physical systems such as ferroelectrics, and much is known about their thermodynamic properties nyb89 . One relevant question in these models is the entropy associated to placing dimers on a regular lattice. For the particular case of full covering of the square lattice by dimers, this question was answered exactly some time ago, using a technique based on pfaffians fis ; kast ; fisetemp . However, even the generalization of this problem for the case on partial covering of the square lattice is still an open question today, no exact result being known.

In this paper we address a generalization of the entropy of dimers problem, considering entropy related to covering the square lattice with chains with $`M`$ monomers each (we will call them $`M`$-mers), as a function of the fraction $`\rho `$ of sites of the lattice occupied by monomers. The chains will be considered flexible, so that there is no energy associated to bending them. Since the only energy in the model is the infinite excluded volume interaction, which forbids the presence of more than one monomer on the same lattice site, the problem is athermal. It may be a simple model for the adsorption of monodisperse flexible chains on the surface of a crystal. Besides the exact solution of the problem for $`M=2`$ and $`\rho =1`$ mentioned above, other cases were already considered in the literature. Rather precise transfer matrix calculations were performed in the polymer limit $`M\mathrm{}`$ for hamiltonian walks ($`\rho =1`$) dup . There are also mean-field approximations flor , Bethe - and Husimi lattice calculations jsmo , and series expansions in $`q^1`$, where $`q`$ is the number of first neighbors of each site in the lattice nem , and in those calculations approximate values for the entropy are obtained for both the full ($`\rho =1`$) and partial ($`\rho <1`$) coverage cases.

In this paper we obtain estimates for the entropy of flexible $`M`$-mers placed on the square lattice, using transfer matrix techniques. This is done calculating numerically exact values for the entropy of the problem on strips with finite widths $`m`$ and periodic boundary conditions and then using finite-size scaling to extrapolate the results to the two-dimensional limit $`m\mathrm{}`$. We separate the problem in cases where $`M`$ is finite or infinite (polymer limit). Also, the case of full coverage ($`\rho =1`$) may be treated separately from the general case. In the general case, it is convenient to address the problem in an ensemble which is grand-canonical with respect to the number of monomers placed on the lattice, whereas for full coverage it is easier to perform a microcanonical calculation.

The expressions we used to calculate the entropy are shown in section II. The model is discussed in more detail and the transfer matrices are described in section III. Our results for the entropies may be found in section IV, as well as the extrapolation procedure and their results. Section V presents final discussions and conclusions.

## II Determination of the entropy

For the case of full coverage, it is convenient to obtain the entropy directly from Boltzmann’s expression

$`s(\rho =1)=\underset{N\mathrm{}}{lim}{\displaystyle \frac{S}{Nk_B}}=\underset{N\mathrm{}}{lim}{\displaystyle \frac{1}{N}}\mathrm{ln}\mathrm{\Omega },`$ (1)
where $`\mathrm{\Omega }`$ is the number of ways to fill the lattice with $`N`$ sites completely with $`M`$-mers. In the polymer limit $`M\mathrm{}`$, we consider a single hamiltonian walk, that is, a self-avoiding walk (SAW) which visits all sites of the lattice.

In the general case where a fraction $`\rho `$ of lattice sites is occupied by monomers, we define the grand-canonical partition function

$`\mathrm{\Xi }(z)={\displaystyle \underset{p}{}}z^{pM}\mathrm{\Gamma }(M,N,p),`$ (2)
where $`z`$ is the activity of a monomer and $`\mathrm{\Gamma }(M,N,p)`$ is the number of ways to place $`p`$ chains with $`M`$ monomers each on the lattice with $`N`$ sites. For the polymer limit, again a single SAW is placed on the lattice and the partition function is defined as

$`\mathrm{\Xi }(z)={\displaystyle \underset{n}{}}z^n\mathrm{\Gamma }(n,N),`$ (3)
where $`\mathrm{\Gamma }(n,N)`$ is the number of ways to place a SAW with $`n`$ monomers on the $`N`$-site lattice. The density of monomers may now be written as

$`\rho (z)=z{\displaystyle \frac{d}{dz}}\varphi (z),`$ (4)
where the thermodynamic potential per lattice site is defined as

$`\varphi (z)=\underset{N\mathrm{}}{lim}{\displaystyle \frac{1}{N}}\mathrm{ln}\mathrm{\Xi }(z).`$ (5)
In the thermodynamic limit, a Legendre transformation allows us to rewrite the potential as

$`\varphi (z)\underset{\rho }{\mathrm{max}}\{\rho \mathrm{ln}z+s(\rho )\},`$ (6)
and thus the entropy may be written as

$`s(\rho )={\displaystyle _0^\rho }\mathrm{ln}z(\rho ^{})𝑑\rho ^{},`$ (7)
with $`s(0)=0`$.

## III Definition of the transfer matrix

We proceed defining a strip of width $`m`$ on the square lattice in the $`(x,y)`$ plane, so that $`1xm`$ and $`\mathrm{}y\mathrm{}`$, with periodic boundary conditions in both directions. A transfer matrix may be built for this problem, inspired on the prescription due to Derrida der  for infinite chains in strips. We thus consider the operation of including an additional step to the strip in the positive $`y`$ direction,adding $`m`$ new sites of the lattice. To properly take into account the statistical weight of the new step, we may define the state of the $`m`$ vertical bonds of the lattice which are incident to the new sites by specifying:

1. The number of monomers already present in the chain which passes through the vertical bond (it is equal to 0 if no chain is present). These numbers may be put into a vector $`|p`$, with $`m`$ components. It is necessary to keep track of this information so that we know when to end each chain.
2. The pairs of bonds which are connected to each other through a path lying entirely below the reference line (see figure 1). These pairs may be also be specified by a $`m`$-component vector $`|v`$, associating a different positive integer to each pair of connected bonds and 0 to the ones which are not connected to any other. This connectivity information prevents us from closing a ring at any level, since this configuration is not allowed in the model.

As an example of these definitions, we consider the case of pentamers ($`M=5`$) placed on a strip of width $`m=4`$. Among the possible configurations of a set of vertical bonds, the one depicted in figure 1 is described by the vectors $`|v=(0,1,1,0)`$ and $`|p=(1,2,2,0)`$. Elements of the line associated to this state of the transfer matrix $`𝒯`$ are obtained considering the possible continuations of the state $`\{|v,|p\}`$ one step upwards, as shown in two examples in figure 2. The resulting state in figure 2.a is described by the vectors $`|v=(0,0,0,0)`$ and $`|p=(2,0,4,0)`$, while the final state in figure 2.b corresponds to the vectors $`|v=(0,0,0,0)`$ and $`|p=(2,0,0,0)`$.

Each monomer placed on a site between the two sets of vertical bonds contributes with an activity $`z`$ to the partition function, so that the element of the transfer matrix which corresponds to the first configuration is equal to $`z^3`$, while the second configuration is associated to an element equal to $`z^4`$ in the transfer matrix. Only the second configuration contributes in the case of full occupancy.

For the polymer case $`M\mathrm{}`$, a single chain passes through the whole strip, so that it is enough to describe the connectivity at a particular set of $`m`$ vertical bonds by indicating the bond which is connected to the initial monomer of the chain (in $`y\mathrm{}`$) and the pairs of bonds connected to each other, exactly as was done in the original work of Derrida der . Thus, a single vector $`|v`$ is enough to describe the state in this limit.

Once the transfer matrix $`𝒯`$ is obtained, the entropy of the model on the strip in the thermodynamic limit is related to the largest eigenvalue of the matrix. For the case of full occupancy ($`\rho =1`$), the number of configurations is given by

$`\mathrm{\Omega }=Tr(𝒯_{}^{}{}_{}{}^{l}),`$ (8)
where $`N=ml`$ is the number of sites and the elements of the matrix $`𝒯^{}`$ are defined by the limit

$`𝒯_{}^{}{}_{i,j}{}^{}=\underset{z\mathrm{}}{lim}{\displaystyle \frac{𝒯_{i,j}}{z^m}}.`$ (9)
The entropy is then related to the largest eigenvalue $`\lambda ^{^{}}`$ of this matrix, so that

$`s(\rho =1)={\displaystyle \frac{1}{m}}\mathrm{ln}\lambda ^{^{}}.`$ (10)
For the general case, where a fraction $`\rho `$ of lattice sites are occupied by monomers, the grand-canonical partition function is related to the transfer matrix through

$`\mathrm{\Xi }(z)=Tr(𝒯^l),`$ (11)
and thus the density $`\rho (z)`$ will be

$`\rho (z)=\underset{N\mathrm{}}{lim}{\displaystyle \frac{z}{N}}{\displaystyle \frac{d}{dz}}ln\mathrm{\Xi }(z)={\displaystyle \frac{z}{m}}{\displaystyle \frac{d}{dz}}ln\lambda ,`$ (12)
where $`\lambda `$ is the largest eigenvalue of the transfer matrix $`𝒯`$. This relation may be inverted to obtain the entropy as a function of the density using equation 7.

## IV Numerical results

The size of the transfer matrix increases very fast with both the molecular weight $`M`$ and the width $`m`$ of the strip grow. This sets an upper limit to the widths we were able to consider for each chain of a given molecular weight. Figure 3 shows this effect, and in the inset one may appreciate that the growth of the transfer matrix is roughly exponential.

Furthermore, as was already observed in similar calculations for polymers hite , the values of the entropy for each class of $`M`$-mers are split into subsets with different finite-size scaling behavior in each subset according to the width of the strips, so that extrapolations must be done within each subset. These splitting seems to be related to frustration effects in the limit of the fully occupied lattice and the subsets are indicated in table 1, where all the widths we considered are given.

The data of each subset were extrapolated to the two-dimensional limit $`m\mathrm{}`$ using the Shanks transformation barb , since we expect finite-size corrections to be exponential. Since at least three values for the entropy are needed in a subset to obtain an estimate for the two-dimensional entropy and its confidence interval, not all subsets may be extrapolated, so that, for example, no estimate could be found for hexamers ($`M=6`$), where we calculated entropies for widths up to $`m=7`$. The final estimate was chosen to be the highest possible extrapolant, and the error associated to it was obtained from the previous generation of extrapolants, through

$`ϵ=\underset{l^{}\mathrm{}}{lim}2|s_{l^{}1}s_{l^{}+1}|.`$ (13)
The extrapolated values of the entropies for $`\rho =1`$ and their uncertainties are displayed in table 2, together with values obtained with other techniques and best values found in the literature.

Our results may be compared with other values in the literature. One may notice that the mean-field estimates are systematically smaller than the values obtained here, but no such general trend is apparent for the Bethe- and Husimi lattice results. The estimate for dimers agrees with the exact value obtained by Fisher, Temperley and Kasteleyn fis ; kast ; fisetemp , and the entropy for hamiltonian walks $`(M\mathrm{})`$ is consistent with both transfer matrix calculations dup  and the result of series expansions up to third order in $`q^1`$ lise , which is $`s_{\mathrm{}}0.38629`$.

Another relevant question is the value of the molecular weight which maximizes the entropy at full occupancy. Mean-field and Bethe lattice results show maximum entropy at $`M=8`$ for a lattice with coordination number $`q=4`$ flor ; jsmo , while series up to second order in $`q^1`$ on the square lattice result in a maximum entropy at $`M=7`$. Our results suggest that this maximum actually occurs at $`M=4`$ on the square lattice, if we suppose that that only one maximum exists in the curve $`s(\rho =1)\times M`$ and also disregard the value obtained for $`M=7`$ due to the large uncertainty associated to it.

For partial occupancy of the lattice, the results are similar to the ones shown in figure 5 for the entropy of dimers as a function of the fraction of occupied lattice sites $`\rho `$. For all the cases we considered the entropy displays a single maximum. The density at which this maximum occurs increases with $`M`$, getting closer to $`\rho 0.79`$, the value found in the polymer limit. The densities and maximum entropies are listed in table 3, and it may be noticed that the largest value for the maximum entropy occurs for tetramers, as was also found for $`\rho =1`$.

In the polymer limit $`M\mathrm{}`$, the model of a polymer placed on a strip in the grand-canonical ensemble displays a first order phase transition at a critical activity $`z_c`$, with the coexistence of a non-polymerized phase ($`\rho =0`$) and a polymerized phase ($`\rho =\rho _c>0`$) pw83 . As the width $`m`$ of the strip is increased, the discontinuity in the density at the transition becomes smaller and in the two-dimensional limit $`m\mathrm{}`$ a continuous transition is found at $`z_c0.3790522`$ gut . The entropy of a polymer on a strip of finite width is therefore not defined for $`\rho <\rho _c`$, as may be seen in figure 6, where the extrapolated value of the entropy as a function of the density for polymers is depicted. As larger widths are considered, the step in the entropy decreases and in the two-dimensional limit we have $`s(\rho =0)=0`$.

## V Conclusion

In this paper we estimate the entropy of chains with $`M`$ monomers each placed on the square lattice as a function of the fraction $`\rho `$ of lattice sites occupied by monomers. The estimates were obtained by extrapolating numerically exact values for the entropy on strips of finite widths $`m`$, calculated using a transfer matrix approach, to the two-dimensional limit $`m\mathrm{}`$.

A particular case where results may easily be obtained analytically is the one dimensional problem ($`m=1`$). In this case, the dimension of the transfer matrix $`𝒯`$ is equal to $`M`$, and we have

$$𝒯_{i,j}=\delta _{i,1}\delta _{j,1}+z(\delta _{i+1,j}+\delta _{i,M}\delta _{j,1}),$$
(14)
where $`1i,jM`$. One may easily obtain the secular equation of this matrix, which is

$$\lambda ^M\lambda ^{M1}z^M=0.$$
(15)
The density may then be found as a function of the largest eigenvalue $`\lambda `$ through using the secular equation above and also expression 12. One obtains

$$\lambda =\frac{1\alpha \rho }{1\rho },$$
(16)
where $`\alpha =(M1)/M`$. Then the entropy may be found by performing the integration in equation 7 changing the integration variable from $`\rho `$ to $`\lambda `$. The result is

$`s`$ $`=`$ $`(1\alpha \rho )\mathrm{ln}(1\alpha \rho )(1\rho )\mathrm{ln}(1\rho )+`$ (17)
$``$ $`\rho (1\alpha )\mathrm{ln}[\rho (1\alpha )].`$
This result is equal to the expression which is obtained if the coordination number of the Bethe lattice result (expression 22 in jsmo ) is taken equal to 2. The entropy in the one-dimensional case vanishes, as expected, for $`\rho =1`$, and the maximum is located at a value of the density which is equal to 1/2 for monomers ($`\alpha =0`$) increasing monotonically with $`\alpha `$ and approaching 1 in the polymer limit $`\alpha 1`$ (figure 7), where the entropy vanishes for all values of $`\rho `$. The value of the maximum entropy is a decreasing function of $`\alpha `$.

For a given value of $`M`$, the density which maximizes the entropy is obtained through he equation

$`(1\alpha \rho )^\alpha (1\rho )[\rho (1\alpha )]^{(1\alpha )}=1`$ (18)
On the square lattice, our calculations indicate that the absolute maximum of the entropy $`s_M(\rho )`$ occurs for tetramers ($`M=4`$). In the one-dimensional case, the maximum entropy of $`s_M(\rho )`$ is a monotonically decreasing function of $`M`$, the absolute maximum $`s_1(1/2)=\mathrm{ln}(2)`$ being obtained for monomers and approaching 0 as $`M\mathrm{}`$.

Finally, the problem of $`M`$-mers confined inside strips of finite width with closed boundary conditions, is an interesting extension of earlier work done in the polymer limit s99 .

###### Acknowledgements.

This research was partially supported by the Brazilian agencies CAPES, FAPERJ and CNPq, whose assistance is gratefully acknowledged.


format: text

text:
# Geometric phase shift in quantum computation using superconducting nanocircuits: nonadiabatic effects

Quantum computation is now attracting increasing interest both theoretically and experimentally. So far, a number of systems have been proposed as potentially viable quantum computer models, including trapped ions, cavity quantum electrodynamics, nuclear magnetic resonce(NMR), etc . In particular, a kind of solid state qubits using controllable low-capacitance Josephson juctions has been paid considerable attention. A two-qubit gate in many experimental implementations is the controlled phase shift, which may be achieved using either a conditional dynamic or geometric phase. A remarkable feature of the latter lies in that it depends only on the geometry of the path executed, and therefore provides a possibility to perform quantum gate operations by an intrinsically fault-tolerant way.

Recently, several basic ideas of adiabatic geometric quantum computation by using NMR, superconducting nanocircuits or trapped ions were proposed. However, since some of the quantum gates are quite sensitive to perturbations of the phase factor of the computational basis states, control of the phase factor becomes an important issue for both hardware and software. Moreover, the adiabatic evolution appears to be quite special, and thus the nonadiabatic correction on the phase shift may need to be considered in some realistic systems as it may play a significant role in a whole process . In this paper, we focus on the nonadiabatic geometric phase in superconducting nanocircuits. We indicate that the adiabatic Berry’s phase, as well as the single qubit gate controlled by this phase, may hardly be implemented in the present experimental setup. On the other hand, since the 2-qubit operations are about $`10^2`$ times slower than the 1-bit operations, the conditional adiabatic phase is extremely difficult to be achieved. A serious disadvantange of the adiabatic conditional phase shift is that the adiabatic condition requires that the evolution time must be much longer than the typical operation time $`\tau _0`$ ($`=\mathrm{}/E_J`$ with $`E_J`$ as the Josephson energy), which leads to an intrinsical time limitation on the operation of quantum gate. Therefore, a generalization to nonadiabatic cases is important in controlling the quantum gates. We find that the nonadiabatic geometric phase shift can also be used to achieve the phase shift in quantum gates.

We first consider a single qubit using Josephson juctions described in Ref. ( see the j-th qubit in Fig.1). The qubit consists of a superconducting electron box formed by an asymmetric SQUID with the Josephson coupling $`E_1`$ and $`E_2`$, pierced by a magnetic flux $`\mathrm{\Phi }`$ and subject to an applied gate voltage $`V_x=2en_x^e/C_x`$ (here we omit the subscript $`j`$, and $`2en_x^e`$ is the offset charge). In the charging regime (where $`E_{1,2}`$ are much smaller than the charging energy $`E_{ch}`$) and at low temperatures, the system behaves as an artificial spin-$`1/2`$ particle in a magnetic field, and the effective Hamiltonian reads

$$\widehat{H}=\frac{1}{2}𝐁\stackrel{}{\sigma },$$
(1)
where $`\sigma _{x,y,z}`$ are Pauli matrices, and the fictitious field

$$𝐁=\{E_Jcos\alpha ,E_Jsin\alpha ,E_{ch}(12n_x^e)\}$$
(2)
with $`E_J=\sqrt{(E_1E_2)^2+4E_1E_2cos^2(\pi \mathrm{\Phi }/\mathrm{\Phi }_0)}`$, $`tan\alpha =(E_1E_2)tan(\pi \mathrm{\Phi }/\mathrm{\Phi }_0)/(E_1+E_2)`$, and $`\mathrm{\Phi }_0=h/2e`$. In this qubit Hamiltonian, charging energy is equivalent to the $`B_z`$ field whereas the Josephson term determines the fields in the $`xy`$ plane. By changing $`V_x`$ and $`\mathrm{\Phi }`$ the qubit Hamiltonian describes a curve in the parameter space $`\{𝐁\}`$. Therefore by adiabatically changing $`\widehat{H}`$ around a circuit in $`\{𝐁\}`$, the eigenstates will accumulate a Berry’s phase $`\gamma _B=\mathrm{\Omega }/2`$, where the signs $`\pm `$ depend on whether the system is in the eigenstate aligned with or against the field. The solid angle $`\mathrm{\Omega }`$, which represents the magnetic field trajectory subtends at $`𝐁=0`$, is derived as

$$\mathrm{\Omega }=_0^\tau \frac{B_x_tB_yB_y_tB_x}{|𝐁|(B_z+|𝐁|)}𝑑t,$$
(3)
under the condition $`𝐁(\tau )=𝐁(0)`$.

However, the adiabatic evolution is quite special, and thus the generalization to nonadiabatic noncyclic cases is of significance. We now recall how to calculate the Pancharatnam phase. For a spin-$`1/2`$ particle subject to an arbitrary magnetic field, each spin state $`|\psi =[e^{i\phi /2}cos(\theta /2),e^{i\phi /2}sin(\theta /2)]^T`$ may be mapped into a unit vector $`𝐧=(sin\theta cos\phi ,sin\theta sin\phi ,cos\theta )`$, with $`𝐧`$ a unit sphere $`S^2`$, via the relation $`𝐧=\psi |\stackrel{}{\sigma |}\psi `$, where $`T`$ represents the transposition of matrix. By changing the magnetic field, the evolution of spin state is a curve on $`S^2`$ from an initial state $`(\theta _i,\phi _i)`$ to a final state $`(\theta _f,\phi _f)`$, and the Pancharatnam phase accumulated in this evolution was found to be

$$\gamma =\frac{1}{2}_C(1cos\theta )𝑑\phi +arctan\frac{sin(\phi _f\phi _i)}{cot\frac{\theta _f}{2}cot\frac{\theta _i}{2}+cos(\phi _f\phi _i)},$$
(4)
where $`C`$ is along the actual evolution curve on $`S^2`$, and is determined by the equation: $`_t𝐧(t)=𝐁(t)\times 𝐧(t)/\mathrm{}`$. This $`\gamma `$ phase recovers the Aharonov-Anandan (AA) phase (Berry’s phase) in a cyclic (adiabatic) evolution.

At this stage, we propose how to detect the nonadiabatic or adiabatic geometric phase in the charge qubit system. The system is prepared in the ground state of the Hamiltonian at $`n_x^e=0`$ and $`\mathrm{\Phi }=0`$, and then changes to the fictitious field $`𝐁(\mathrm{\Phi }(t),n_x^e(t))`$, which is a periodic function of time $`t`$ with the period $`\tau `$. We consider the process where a pair of orthogonal states $`|\psi _\pm `$ evolve cyclically(but not necessary adiabatically). This process can be realized in the present system. Noting that the adiabatic approximation is merely a sufficient but not necessary condition for the above cyclic evolution, we here focus on a nonadiabatic generalization. In this evolution, the initial state is given by

$$|\psi _i=a_+|\psi _+(\theta _i,\phi _i)+a_{}|\psi _{}(\theta _i,\phi _i),$$

where

$`|`$ $`\psi _+(\theta ,\phi )=[e^{i\phi /2}cos(\theta /2),e^{i\phi /2}sin(\theta /2)]^T,`$
$`|`$ $`\psi _{}(\theta ,\phi )=[e^{i\phi /2}sin(\theta /2),e^{i\phi /2}cos(\theta /2)]^T,`$
$`a`$ $`{}_{+}{}^{}=cos[(\eta \theta _i)/2]cos\phi _i/2icos[(\eta +\theta _i)/2]sin\phi _i/2,`$
$`a`$ $`{}_{}{}^{}=sin[(\eta \theta _i)/2]cos\phi _i/2+isin[(\eta +\theta _i)/2]sin\phi _i/2,`$
with $`tan\eta =E_J(\mathrm{\Phi }=0)/E_{ch}`$, $`tan\theta _i=[E_J(t)/B_z(t)]|_{t=0}`$, and $`tan\phi _i=[B_y(t)/B_x(t)]|_{t=0}`$. A phase difference between $`|\psi _{+/}`$ can be introduced by changing $`\widehat{H}`$. The phases acquired in this way will have both geometrical and dynamical components. But the dynamical phase accumulated in the whole procedure may be removed, thus only the geometric phase remains. By taking into account the cyclic condition $`𝐧(0)=𝐧(\tau )`$ for $`|\psi _\pm `$, the final state in this case is given by

$$|\psi _f=a_+e^{i\gamma }|\psi _+(\theta _i,\phi _i)+a_{}e^{i\gamma }|\psi _{}(\theta _i,\phi _i),$$
(5)
where $`\gamma `$ can be calculated from Eq.(4). The contribution from the second term of Eq.(4) vanishes simply because $`𝐧(0)=𝐧(\tau )`$. Thus the geometric phase considered here is the cyclic AA phase. The probability of measuring a charge $`2e`$ ($`n=1`$) in the box at the end of this procedure is derived as

$$P_1=|a_+sin\frac{\theta _i}{2}+a_{}cos\frac{\theta _i}{2}e^{2i\gamma }|^2.$$
(6)
This probability can be simplified to

$$P_1=[1cos(\eta \theta _i)cos\theta _i+sin(\eta \theta _i)sin\theta _icos2\gamma ]/2.$$
(7)
when $`\mathrm{\Phi }(0)=0`$. Note that Eq.(7) recovers $`sin^2(\gamma )`$ in Ref. even in a nonadiabatic but cyclic evolution. Thus the nonadiabatic phase may be determined by the probability of the charge state in the box at the end of this process. It is worth pointing out that the parameters $`\eta `$ and $`\theta _i`$ in Eq.(6)(or(7)) are fully determined by the experimentally controllable parameters $`\mathrm{\Phi }`$ and $`n_x^e`$, as in the adiabatic Berry’s phase case .

It is remarkable that the probability obtained in Eq.(6)(or(7)) may be directly detected by the dc current through the probe junction $`C_b^j`$ under a finite bias voltage $`V_b^j`$. Assume that we have achieved one SQUID qubit as well as the detector circuit, as shown in Fig.1. By changing $`V_{jx}`$ and $`\mathrm{\Phi }_j`$ in time $`[0,\tau ]`$, the system oscillates between $`|0`$ and $`|1`$, and the final state would be determined by the geometric phase. The measurable dc current through the probe junction formulates by the processes: $`|1`$ emits two electrons to the probe, while $`|0`$ does nothing. Consequently, the probability described by Eq.(6)(or (7)) as well as the geometric phase may be detected by the dc current.

The single qubit gate may be realized by this geometric phase. For example, it is straightforward to check that the unitary evolution operator defined by $`|\psi _f=U_1^{sq}|\psi _i`$, is given by

$$U_1^{sq}(\gamma )=\left(\begin{array}{cc}cos\gamma \hfill & isin\gamma \hfill \\ isin\gamma \hfill & cos\gamma \hfill \end{array}\right),$$
(8)
when $`\theta _i=0`$ and $`\phi _i=0`$. Clearly, the operation depends on the geometric phase $`\gamma `$; $`\gamma =\pi /2`$ and $`\gamma =\pi /4`$ procduce a spin flip (NOT-operation) and an equal-weight superposition of spin states, respectively. On the other hand, the phase-flip gate $`U_2^{sq}=exp(2i\gamma |11|)`$ (up to an irrelevant over phase ) is derived by $`\theta _i=0`$ and $`\phi _i=0`$. The noncommutable $`U_1^{sq}`$ and $`U_2^{sq}`$ gates are the two well-known universal gates for single-qubit operation. The Berry’s phase may be used to achieve intrinsical fault-tolerant quantum computation since it depends only on the evolution path in the parameter space. The nonadiabatic cyclic phase is also rather universal in a sense that it is the same for a infinite number of possible ways of motion along the curves in the projective Hilbert space. Consequently, the nonadiabatic phase may also be used as a tool for some fault-tolerant quantum computation.

We now illustrate how to achieve the cyclic state for quantum gates in two processes. The parameters $`(\mathrm{\Phi }(t),n_x^e(t))`$ in process I change as

$$(I)\{\begin{array}{cc}\{\frac{4\mathrm{\Phi }_mt}{\tau },\frac{1}{2}\},\hfill & t[0,\frac{\tau }{4})\hfill \\ \{\mathrm{\Phi }_m,\frac{1}{2}+4(n_{xm}^e\frac{1}{2})(\frac{t}{\tau }\frac{1}{4})\},\hfill & t[\frac{\tau }{4},\frac{1\tau }{2})\hfill \\ \{4\mathrm{\Phi }_mt/\tau +3\mathrm{\Phi }_m,n_{xm}^e\},\hfill & t[\frac{1\tau }{2},\frac{3\tau }{4})\hfill \\ \{0,n_{xm}+4(\frac{1}{2}n_{xm}^e)(\frac{t}{\tau }\frac{3}{4})\}.\hfill & t[\frac{3\tau }{4},\tau )\hfill \end{array}$$
(9)
The path in the parameter space $`\{𝐁\}`$ swept out in this case is exactly the same as that proposed in Ref., Since the evolution in this process is cyclic only under the adiabatic condition, we need to answer a key question: whether the adiabatic approximation is valid for the given parameters? As for process II , the parameters $`(\mathrm{\Phi }(t),n_x^e(t))`$ change as

$$(II)\{\begin{array}{c}\mathrm{\Phi }(t)=\frac{\mathrm{\Phi }_0}{\pi }atan[\frac{E1+E2}{E1E2}tan(\omega t)],\\ n_x^e(t)=\frac{1}{2}(1\frac{E_Jctg\chi _0+\mathrm{}\omega }{E_{ch}}).\end{array}$$
(10)
The fictitious field described by Eq.(10) guarantees that the angle $`\chi _0=arctan[E_J/(B_z(t)\mathrm{}\omega )]`$ (and $`n_z`$) is time-independent. It is found that the state described by the vectors $`𝐧(\chi _0,\omega t)`$ in this process evolves cyclically with period $`\tau =2\pi /\omega `$ , and the AA phase for one cycle is given by $`\gamma =\pi (1cos\chi _0)`$, which may be used to achieve the mentioned single-qubit gates geometrically. For the present system, the dynamic phase can be removed by simply choosing $`\omega =4(E_1+E_2)E_k[4E_1E_2/(E_1E_2)^2]/\pi sin(2\chi _0)`$ with $`E_k(x)`$ the complete elliptic integral of the first kind.

The nonadiabatic effect should be important if $`\tau `$ is not short. We first consider the evolutions described by Eq.(9). Figure 2 shows $`n_z(t)`$ and $`\widehat{B}_z(t)=B_z(t)/|𝐁(t)|`$ versus time, with the parameters being the same as those in Ref.. The deviation of $`𝐧(t)`$ from $`\widehat{𝐁}(t)(=𝐁(t)/|𝐁(t)|)`$ indicates clearly whether or not the adiabatic approximation is valid because $`𝐧(t)`$ almost follows the trajectory of the magnetic field $`\widehat{𝐁}(t)`$ under this approximation. It is seen from Fig.2 that the adiabatic approximation is satisfied in the first case when $`\tau >500\tau _0`$, where $`\tau _0=\mathrm{}/(E_1+E_2)84ps`$. The adiabatic condition for process II is in the same order of magnitude (see Fig.3). It is worth pointing out that the coherence time achieved in a single SQUID is merely about $`3040\tau _0`$, which is not long enough for the adiabatic evolution, implying that the adiabatic condition is not satisfied in the above two processes for realistic systems. But, fortunately, the nonadiabatic phase can be measured and used in achieving geometric quantum gates.

Conditional geometric phase accumulated in one sub-system evolution depends on the quantum state of another sub-system, which may be realized by coupling capacitively two asymmetric SQUIDS (see any neighboring pair of qubits in Fig.1.). If the coupling capacitance $`C_{ij}`$ is smaller than the others, the Hamiltonian reads

$$\widehat{H}=\underset{i=1}{\overset{N}{}}\widehat{H}_i+\underset{i=1}{\overset{N1}{}}(\widehat{H}_{i,i+1}+H.C.),$$
(11)
where $`H_i`$ refer to the uncoupled qubits defined in Eq.(1) and $`\widehat{H}_{i,i+1}=E_{i,i+1}(n_i^en_{x,i}^e)(n_{i+1}^en_{x,i+1}^e)`$ with $`E_{i,i+1}=E_{ch}C_{i,i+1}/C`$. The gate voltage and magnetic flux can be independently fixed for all qubits. We address firstly a two-qubit operation, e.g., $`i`$ and $`j`$ qubits are two neighbour qubits with the i-th as the control qubit and the j-th as the target qubit. The fictitious field on the target qubit is $`[E_J(\mathrm{\Phi }_j)cos\alpha _j,E_J(\mathrm{\Phi }_j)sin\alpha _j,B_z^l]`$ with $`B_z^l=E_{ch}(12n_{x,j}^e)+E_{i,j}(n_{x,i}^el)`$, where $`l`$ represent the control qubit state $`0`$ or $`1`$. Obviously, the geometric phase $`\gamma _j`$ for $`j`$-th qubit in decoupled case is different from $`\gamma _j^l`$ even the changeings of $`(\mathrm{\Phi }_j,n_{x,j}^e)`$ are the same, where $`\gamma _j^l`$ is the geometric phase of the target qubit when the charge state of the control qubit is $`l`$. $`\gamma _j^l`$ may be directly derived from Eq.(4). It is worth to pointing out that the state described by the vector $`𝐧(\chi ^l,\omega t)`$ with $`\chi _0^l=arctan[E_J/(B_z^l\mathrm{}\omega )]`$ is still a cyclic evolution, and may be used to achieve the two-qubit operation. In terms of the basis $`\{|00,|01,|10,|11\}`$, the unitary operator to describe the two-qubit gate is given by

$$U_{(\gamma _j^0,\gamma _j^1)}=diag(e^{i\gamma _j^0},e^{i\gamma _j^0},e^{i\gamma _j^1},e^{i\gamma _j^1}).$$
(12)
The combination with single-bit operations allow us to perform the XOR gate. The unitary operation for XOR gate can be obtained by $`U_{XOR}=\left[IU_1^{sq}(\pi /4)\right]U_{(0,\frac{3\pi }{2})}\left[IU_1^{sq}(\pi /4)\right]^{}`$ with $`I`$ as a $`2\times 2`$ unit matrix. This XOR gate together with single qubit gates constitutes a universality: they are sufficient for all manipulations required for quantum computation . Therefore, all the elements of quantum computation may be achievable by (nonadiabatic) geometric phase.

We now compute the geometric phases required by the spin flip operation ($`\gamma =3\pi /2`$) and NOT operation ($`\gamma =3\pi /4`$) accumulated in the second process. The comparison of the nonadiabatic geometric phase $`\gamma `$ with $`\gamma _a`$ is shown in Fig.3, where the $`\gamma _a`$ is the phase calculated under the adiabatic approximation. It is seen that the $`\gamma _a`$ deviates evidently from the $`\gamma `$ for $`\tau <150\tau _0`$. Thus the operation time required by the adiabatic condition in both processes I and II is in the same order of magnitude. The dynamic phases can be removed when $`\tau 3.57\tau _0`$ for $`\gamma =3\pi /2`$ and $`\tau 4.12\tau _0`$ for $`\gamma =3\pi /4`$, respecitively. Therefore, by accurately controlling the parameters $`\mathrm{\Phi }`$ and $`n_x^e`$, we may control the state in the projective Hilbert space. It is striking that the present operation time for nonadiabatic geometric gates is much shorter than that with the adiabatic scheme. Note that the coherence time achieved in a single SQUID by current technology is about $`3040\tau _0`$, implying that tens of geometric NOT operation may be achieved experimentally. Therefore, the generalization of the adiabatic phase to the nonadiabatic case is of significance since the coherence time achieved in charge qubit in Josephson Junctions is short. Moreover, the large number qubits required for useful computation may be devised by a network similar to Fig.1.

In conclusion, we study how to detect the nonadiabatic phase in superconducting nanocircuits, and the possibility to use the nonadiabatic phase as a tool to achieve the quantum computation.

We wish to acknowledge valuable discussions with Dr. L. M. Duan and Dr. X. B. Wang. This work was supported by the RGC grant of Hong Kong under Grant Nos. HKU7118/00P and HKU7114/02P, the Ministry of Science and Technology of China under Grant No. G1999064602, and the URC fund of HKU. S. L. Z. was supported in part by the SRF for ROCS, SEM, the NSF of Guangdong under Grant No. 021088, and the NNSF of China under Grant No. 10204008.


text:
# Magnetic Field Effects on the Transport Properties of One-sided Rough Wires

## Abstract

We present a detailed numerical analysis of the effect of a magnetic field on the transport properties of a ‘small-$`N`$’ one-sided surface disordered wire. When time reversal symmetry is broken due to a magnetic field $`B`$, we find a strong increase with $`B`$ not only of the localization length $`\xi `$ but also of the mean free path $`\mathrm{}`$ caused by boundary states. Despite this, the universal relationship between $`\mathrm{}`$ and $`\xi `$ does hold. We also analyze the conductance distribution at the metal-insulator crossover, finding a very good agreement with Random Matrix Theory with two fluctuating channels within the Circular Orthogonal(Unitary) Ensemble in absence(presence) of $`B`$.

The physics of disordered media has attracted the attention of scientists from different areas since long agoReviews , and it is an up-to-date topic due to its relevance in wave transport, be it in the form of quantum (electrons) or classical (electromagnetic or acustic) waves. Transport in disordered media, in particular in the mesoscopic regime, presents a rich panoply of phenomena, such as *Anderson localization*Anderson , *Weak localization*Weak , and *Universal Conductance Fluctuations*UCF , just to mention a few.

When dealing with disordered systems, it is of fundamental importance to know about the full statistical properties of the relevant transport quantities, i.e. the conductance $`G`$, the eigenchannels $`\tau `$ and the transmission $`T_i,T_{ij}`$ and reflection $`R_i,R_{ij}`$ coefficients. Much is known about statistics of both the transmissionLangen ; PRL1  and the reflectionPRL2  coefficients, since the full distribution function has been obtained from a Random Matrix Theory (RMT) approachLangen ; PRL2  and numerical techniquesPRL1 ; PRL2  as well as experimentally in the microwave regimeStoytchev .

Several approaches have been considered to analyze the effect of disorder on the transport properties, ranging from analytical techniques using RMTDMPK ; RMT ; Pichard ; muttalib  to numerical onesNumerical ; gatoprl . Most of them consider that the disorder is inside the material where the wave travels through (bulk disorder). However, the recent developments in nanotechnology open the way to new sources of scattering. In particular, as system sizes are shrunk down to the nanometer scale, the surface-to-volume ratio becomes larger, and surface effects (for example disorder) become very important. This is more relevant in the case of few-modes quantum wires, where the effect of surface disorder (roughness) can dominate their transport properties. In addtion, due to improvements in fabrication techniques bulk disorder is almost absent. Many statistical properties of these surface disordered wires have been shown to be in very good agreement with those predicted by RMTAPL ; PRL1 ; PRL2 ; surface ; gatoprl .

The consequences of time reversal symmetry breaking in disordered media, such as the doubling of the localization length as a consequence of the phase mismatch of time reversed paths, are conceptually known since long ago (see for instance Reviews ; RMT ). Nevertheless, a fundamental question was still not answered: how does the transition take place? This issue has recently been addressed, both theoreticallykolesnikov  and numerically, by means of an Anderson model (bulk disorder)schomerus , finding that the expected doubling of the localization length takes place in a rather smooth way, while the mean free path remains constant.

Another conceptual issue concerns the shape of the full conductance distribution $`P(G)`$ in the diffusion-localization transition. It is known that $`P(G)`$ evolves from Gaussian, in the diffusive regime, to log-normal in the localization regime. It has been recently predictedmuttalib  that, in the crossover regime, $`P(G)`$ exhibits a highly non trivial shape; the subsequent numerical analysis revealed a direct relation between the number of eigenchannels that actually describe the distribution and its universal characteristicsgatoprl .

In this paper we analyze the consequences that the application of a magnetic field has on the transport properties of a one-sided surface disordered wire (see Fig. 1). We show that breaking time reversal symmetry leads to drastic effects on the electron propagation through rough wires in the small-number-of-channels regime. The dependence of the mean free path $`\mathrm{}`$ and localization length $`\xi `$ on the magnetic field strength is studied by means of an analysis of $`1/G`$ and $`\mathrm{ln}G`$ as a function of both the length and the magnetic field strength. This analysis reveals a peculiar behavior of this two relevant length scales: $`\xi `$ increases, but much more than the expected factor of twokolesnikov ; schomerus , while at the same time $`\mathrm{}`$, which is expected to remain constant (for the low fields considered), also increases. We find that the latter phenomenon is related to the precursors of the edge states (we are far away from the quantum Hall regime), which propagate almost freely due to the specific type of disorder. We show that the unexpected large increase of $`\xi `$ does not contradict the predictions of RMT when the anomalous behavior of $`\mathrm{}`$ is taken into account. The validity of RMT is confirmed by the analysis of the numerical conductance distribution $`P(G)`$ in the vicinity of the diffusion-localization transition. This analysis shows that the universal features of the distribution are not affected by the lifting of time reversal invariance.

The system which we are dealing with (depicted schematically in Fig. 1) is made of two semi-infinite two-dimensional clean leads of width $`W`$ and a surface-disordered region of length $`L`$ in between. The corrugated region is composed of slices of length $`\mathrm{}`$ laterally bounded by hard walls, whose random width is uniformly distributed in the interval $`[W\delta ,W+\delta ]`$ around the mean value $`W`$. We will take $`W/\lambda =2.6`$ ($`\lambda `$ being the electron Fermi wavelength), which allows five propagating channels in the clean part $`N=5`$, and a disorder strength $`\delta `$ so that $`W/\delta =7`$note0 . We consider *One-sided Surface Disorder* (OSD), i.e., the roughness is located only on one surface. In order to study how the presence or absence of time reversal symmetry affects the transport, a magnetic field $`B`$ perpendicular to the plane of the wire is considered.

The conductance (in units of $`2e^2/h`$) of the total system is simply given by $`G=\mathrm{trace}\{tt^{}\}`$. The transmission matrix $`t`$ is calculated for each realization of the disorder by a generalized scattering matrix techniqueAPL ; WeissJosanMetodo  where the magnetic field has been introduced in the formalism as in the method of Ref.gvrboese . This consists in choosing a gauge, in which only the component of the vector potential along the transverse direction is non zero; the system is then divided in discretization slices (separated by thin lines in Fig. 1), and a constant vector potential is assigned to each slice. To minimize errors due to the discretization, the slices have to be chosen so that they contain less than a flux quantum, furthermore convergence can be checked by increasing the number of them. Then, in order to obtain the statistical features introduced by the disorder, we perform configuration averages $`\mathrm{}`$. To calculate the distributions we average over 10000 different configurations of the disorder.

In the upper panel of Fig 2 we plot $`\mathrm{ln}G`$ versus $`L/W`$ for four different values $`B`$. The magnetic field strength is given by the quantity $`\mathrm{\Phi }_\xi =\xi _0WB`$ that measures the amount of magnetic flux trapped in a portion of wire with a length equal to the localization length at zero field ($`\xi _0`$). Throughout this work all magnetic fluxes are given in units of the flux quantum $`\mathrm{\Phi }_0`$.

The localization length $`\xi _B`$ for the different values of $`B`$ is obtained from the linear part of the plot since $`\mathrm{ln}GL/\xi _B`$. The value of $`\xi _B`$ evolves as expected, the higher $`B`$ the larger $`\xi _B`$, but instead of saturating to the expected value $`\xi _B=2\xi _0`$RMT ; Pichard ; kolesnikov ; schomerus  when time reversal invariance is broken, for OSD this does not occur giving rise to values $`\xi _B2\xi _0`$.

The origin of this lies in the fact that, whereas usually $`\mathrm{}`$ is assumed to be constant for the small values of $`B`$ that are considered (see, for example, the results for a two-dimensional Anderson model with on-site disorder of Ref. schomerus ), in wires with OSD one has to take the magnetic field dependence of $`\mathrm{}`$ into account. This effect can be seen in lower panel of Fig. 2, where the average resistance $`1/G`$ is plotted as a function of $`L/W`$. By virtue of $`1/G1/N+L/(N\mathrm{})`$ (as long as $`\mathrm{}<L<\xi `$) the best linear fit in the diffusive regime gives $`\mathrm{}_B`$.

The general relationship between mean free path and localization length predicted by RMT from general symmetry considerations can be written as

$$\xi _\text{B}=\beta N\mathrm{}_\text{B},$$
(1)
with $`\beta `$ being 1 in the presence of time reversal symmetry, and 2 when it is broken (by the magnetic field). In Eq. (1) the magnetic field dependence of the mean free path has to be accounted for.

Eq. (1), predicts a transition from $`\beta =1`$ to $`\beta =2`$, as soon as time reversal symmetry is removed, i.e. when the magnetic flux comprised in a section of wire of length of the order of $`\xi _0`$ is about one flux quantum. As shown in Fig. 3, both $`\xi _B`$ and $`\mathrm{}_B`$ increase smoothly when $`B`$ increases and in the same qualitative way. In fact, the tremendous increase of $`\mathrm{}_B`$, allows the possibility to have huge values of $`\xi _B`$ without violating Eq. (1) in any case (see inset in Fig. 3). The huge increase of these two quantities gives rise to the actual possibility of changing the length of the diffusive regime with a magnetic field. The reason for that lays in the presence of one clean surface. In fact, even for really small $`B`$, the developing edge states may propagate along the clean surface suffering almost no scattering. We have checked that this is indeed the case by doing the same analysis for a wire with *two-sided surface disorder*, with uncorrelated disorder of equal strength on both surfaces (not shown here). In the absence of the clean surface, the non fully developed edge states still suffer from scattering, and no appreciable variation of the mean free path is found. It is important to stress once more that the magnetic fields that give rise to the aforementioned effect are small (cyclotron radius much bigger than wire width). Thus, we are *far away* from the quantum Hall regime, where it is natural that edge states propagate freely without suffering any scattering. We note in passing that the magnetoresistance of these one-sided rough wires may offer interesting technical possibilities.

It is worth to mention that our numerical simulations show that the lifting of time reversal invariance is a gradual process, with $`\beta =2`$ reached for around 100 flux quanta per localization length times system width, as it is shown in the inset of Fig. 3. This result is in perfect agreement with the values presented in Ref. schomerus  even though several characteristics of the transport depend on the kind of disorder. We also see no evidence of the ’two-scale localization’ predicted in Ref. kolesnikov  and not found in Ref. schomerus .

Now we go beyond the first moments of the distribution, and devote our attention to the entire distribution function $`P(G)`$ at the transition point. It has been recently reportedgatoprl  that $`P(G)`$ contains some universal characteristics that depend only on the value of $`G`$, and in particular, that the transition point corresponds to a transition in the number of active eigenchannels. We would like to address these issues, when the action of the magnetic field eliminates time reversal invariance, by analyzing $`P(G)`$ in the vicinity of the diffusion-localization transition. Let us focus on the characteristic features of $`P(G)`$ just before the localization onset, more explicitly when $`G=1`$. At this value of the mean conductance the distribution has a well define cusp at $`G=1`$. In Fig. 4 we show that $`P(G)`$ for $`G=1`$ shrinks with increasing $`B`$ and at the same time the cusp becomes higher and more pronounced. This confirms that general characteristics of $`P(G)`$ at the transition point actually does not depend on the time reversal invariance conditions.

It has been shown that the RMT pictureBaran94  within the maximum entropy Ansatz (this leads to the Circular Orthogonal Ensemble (COE) when time reversal symmetry is present to the Circular Unitary Ensemble (CUE) when it is broken) gives an almost quantitative description of the distributions when time reversal symmetry is preservedgatoprl . Following the ideas of Ref. gatoprl  we obtain $`P(G)`$ for two active eigenchannels (Eq. 2):

$`P(G,G)`$ $`{\displaystyle \frac{2}{3}}(1{\displaystyle \frac{G}{2}})^3\left({\displaystyle \frac{G}{2}}\right)^{2\alpha }\times `$ (2)
$`{}_{2}{}^{}F_{1}^{}({\displaystyle \frac{3}{2}},\alpha ,{\displaystyle \frac{5}{2}},\left({\displaystyle \frac{2}{G}}1\right)^2)\mathrm{\Theta }(G1)`$
$`+`$ $`{\displaystyle \frac{\sqrt{\pi }}{2}}\left({\displaystyle \frac{G}{2}}\right)^{2\alpha +3}{\displaystyle \frac{\mathrm{\Gamma }(1+\alpha )}{\mathrm{\Gamma }(\frac{5}{2}+\alpha )}}\mathrm{\Theta }(1G).`$
where

$$\alpha =\frac{4(G1)}{2G},$$
(3)
depends only on the value of the mean conductancenote1 .

This distribution is shown in in the insets of Fig. 4, together with Eq. 3 from Ref. gatoprl  for COE. The qualitative agreement with the numerical data is striking.

In conclusion, we have presented extensive numerical simulations for nanowires with one-sided surface disorder in the presence of magnetic fields. The transition between the orthogonal and unitary case has been demonstrated to happen smoothly. Although our model of disorder exhibits a large increase of the mean free path leading to giant wave delocalization, the relationship between mean free path and localization length predicted by RMT is fulfilled. The behavior of both mean free path and localization length suggest the possibility of manipulating the diffusion and localization onsets by applying a magnetic field. The conductance distribution around the metal-insulator transition has confirmed its universal nature, being very well described by RMT with the maximum entropy Ansatz and two fluctuating channels.

We are gratefully indebted to F. Evers, A.D. Mirlin, J.J. Sáenz and A. Rosch for valuable discussions. AGM acknowledges support from the Emmy-Noether program of the DFG under Grant No. Bu 1107/2-1. This work was partially supported by the DFG Center for Functional Nanostructures at the University of Karlsruhe.


text:
# GROUND STATE OF THE DOUBLE-WELL CONDENSATES INTERACTING WITH TRAP OSCILLATIONS

## I Introduction

The realization of two-component condensates Myatt et al. (1997) as well as condensates in a double-well potential Andrews et al. (1997) has attracted considerable interest, both experimental and theoretical (see, e.g., Leggett (2001) and references therein) to the phenomena caused by phase coherence of two condensate modes. In the work Andrews et al. (1997) spatial quantum coherence was observed by means of interference pattern in two overlapping condensates. This interference pattern was conformed in the work Rohrl et al. (1997) with the numerical simulation of the Gross-Pitaevskii equation. In the work Hall et al. (1998) coherent oscillations of the relative populations have been observed in driven two-component condensates with different internal states. As is well known, a bright manifestation of phase coherence is Josephson effect. In numerous studies devoted to the Josephson effect in systems of two condensates in different internal states Hall et al. (1998), Williams et al. (1999) or in a double-well potential Zapata et al. (1998)-Milburn et al. (1997) coherent Josephson oscillations are considered for various dynamical regimes caused by the competition between tunneling and intracondensate interaction (nonlinearity). In the work Zapata et al. (1998) the Josephson coupling energy is calculated for small-amplitude oscillations in a double-well potential. There are estimated damping effects due to the normal currents at finite temperature. In the work Milburn et al. (1997) it is shown that for a relatively weak interaction, the particle number oscillations between the condensates are complete. They are suppressed when the total number of atoms in the condensates exceeds a critical value and the behavior of the system is governed by the nonlinearity. Nonlinear Josephson-type oscillations in the relative oscillations of driven two-component condensates are studied in the work Williams et al. (1999). Decoherence effects and quantum corrections to mean-field solutions have been considered in Lewenstein and You (1996), Imamoglu et al. (1997). In Meir and Zweger (2001) the damping effects (even at zero temperature) of the Josephson current are derived by means of an functional integral approach. The detailed treatment of the nonlinear classical dynamic of the condensates in a double-well potential has been performed in the works Raghavan et al. (1999),Smerzi and Fantoni (1997). In the work Pitaevskii and Stringari (2001) the quantum and thermal fluctuations of the phase are studied for the condensates in the double well potential.

It should be emphasized that experimental observation of the Josephson effect is difficult because the small energy splitting associated with Josephson coupling means that thermal and quantum fluctuations will destroy the phase coherence between two condensates even at the lowest achievable temperatures Meir and Zweger (2001), Pitaevskii and Stringari (2001). While the energy splitting can be increased, for instance, by lowering of the barrier height, it then becomes comparable with that of motion states of the condensates.

However, the problem of the interaction between the Josephson degrees of freedom and states of motion (oscillations) of the trapped condensate has yet to be analyzed. The present paper focuses on a mechanism of formation of the bound state of the Josephson degrees of freedom and trap oscillations due to their interaction. The mechanism proposed in what follows may be important for detection of the excited Josephson states.

The results achieved in this paper may be listed as follows.

1. To formulate the problem of the interaction between the Josephson and oscillation degrees of freedom in an adequate manner, the quantizied spectrum of the particle number generated by Josephson coupling is derived. The states of this spectrum represent the quantum analogue of the nonlinear coherent Josephson oscillations considered in Raghavan et al. (1999). In what follows, the states of the quantizied spectrum are named as Josephson states. The spectrum is highly nonequidistant and has the logarithmic singularity in the density of states at the $`2E_J`$ energy, where $`E_J`$ is the Josephson coupling energy. As is shown in what follows, any Josephson state can be realized by means of a given initial disbalance of the particle number in two condensates.

2. There is considered the interaction between the trap oscillations and excited Josephson states corresponding to the large enough initial disbalance of the particle number. It is shown that this interaction is responsible for the formation of the bound state of $`\overline{n}_m1`$ oscillation quanta with the Josephson state corresponding to the definite initial disbalance. In the Thomas-Fermi approximation at $`\mu \omega _0`$, where $`\mu `$ is the chemical potential of the condensates, $`\omega _0`$ is the characteristic frequency of a trap, $`h=1`$, the bound state arises in the region of the dense oscillation spectrum. In this region the level separations are small in comparison with the harmonic oscillation frequency $`\omega _0`$. The equilibrium values of the oscillation quanta and initial disbalance of the particle number are coupled self-consistently and can realize an absolute minimum of the thermodynamic energy at the large enough interaction. The thermodynamic average $`\overline{n}_m0`$ generates the equilibrium distortion of the condensate shape. It allows us to detect the Josephson states by observing the change in the condensate shape.

## II The quantum spectrum of the particle number

A Josephson coupling is realized for the condensates in an symmetrical/asymmetrical double-well potential formed by two different traps with a barrier between them Andrews et al. (1997). The barrier is created by the laser light, and its height is directly proportional to the laser power and, thus, can be varied with ease. The proposed mechanism is also suitable for the condensates in different internal states in the same trap. Experimentally, this may be a superposition of two $`\mathrm{Rb}^{87}`$ condensates in the $`|F=1,m_F=1,|F=2,m_F=1`$ states Hall et al. (1998), Williams et al. (1999). In the case of a weak Josephson coupling, the basis states are the self-consistent ground states in two condensates separately. The wave function of the condensate with Josephson coupling takes the form of a superposition of these states, namely, $`\mathrm{\Psi }(𝐫,t)=\psi _1(𝐫)a_1(t)+\psi _2(𝐫)a_2(t)`$, where $`\psi _i(𝐫)`$ are the normalized solutions of the Gross-Pitaevskii equation, $`a_i(t)=N_i^{1/2}(t)e^{i\theta _i(t)}`$, $`i=1,2`$, $`N_{1,2}`$, $`\theta _i(t)`$ are the particle numbers and phases of each condensate.

As is well known Leggett (2001), Zapata et al. (1998), Milburn et al. (1997), Smerzi and Fantoni (1997) the Hamiltonian of two condensates with Josephson coupling has the form

$$H_JE_0=E_C(\mathrm{\Delta }N)^22E_J\mathrm{cos}\varphi ,$$
(1)
here $`\mathrm{\Delta }N=N_1N_2`$, $`E_C=\mu /N`$, $`\mu \mu _1=\mu _2`$ are the chemical potentials, $`N`$ are the total particle numbers, $`\varphi =\theta _1\theta _2`$ is the relative phase of the condensates. The quantities $`E_C`$, $`E_J`$ depend on the total particle number $`N`$. $`(\mathrm{\Delta }N)`$, $`\varphi `$ are canonically conjugate variables. In (1) the energy origin is the mean-field summary energy of the condensates, namely, $`E_0\mu N`$.

The quantization of the Hamiltonian (1) produces the spectrum of the particle number in the Josephson potential $`E_J\mathrm{cos}\varphi `$. As is shown in what follows, any Josephson state can be realized by means of a given initial disbalance of the particle number. For this reason, it is interesting to obtain the complete spectrum generated by the Hamiltonian of (1), and to represent it as a function of the initial disbalance.

The Schrodinger equation for the Hamiltonian (1) is derived by means of the quantization rule: $`(\mathrm{\Delta }N)i/\varphi `$. As a result, we obtain Mathieu’s equation

$$\left[\lambda \frac{d^2}{d\varphi ^2}2E_J\mathrm{cos}\varphi _r\right]\mathrm{\Psi }=\epsilon \mathrm{\Psi },\epsilon H_JE_0.$$
(2)
For $`\epsilon >2E_J`$ this equation has a continuous spectrum. The states of this spectrum correspond to the classical states with unlimited phase change: $`\mathrm{}<\varphi <+\mathrm{}`$ that are named as the self-trapping states in the works Raghavan et al. (1999); Smerzi and Fantoni (1997). In the region $`2E_J<\epsilon <2E_J`$ Eq.(2) has the discrete spectrum. It corresponds to the region of finite motion of Hamiltonian (1), where the relative phase changes within $`\mathrm{arccos}(\epsilon /2E_J)<\varphi <\mathrm{arccos}(\epsilon /2E_J)`$ for every $`\epsilon `$. In the Josephson regime at $`E_JE_C`$ Leggett (2001), the number of levels in a well is large and the discrete spectrum is specified by the Bohr-Sommerfeld formula:

$$\begin{array}{cc}& \nu (\epsilon _\nu )=\frac{d\varphi _r}{\pi }\mathrm{\Delta }N(\varphi _r;\epsilon _\nu )\hfill \\ \hfill =& \frac{d\varphi _r}{\pi }\left[\frac{1}{E_C}\left(\epsilon _\nu +2E_J\mathrm{cos}\varphi \right)\right]^{1/2}\hfill \\ & =\nu _c\left[E(\kappa )(1\kappa ^2)K(\kappa )\right];\hfill \end{array}$$
(3)
$$\nu _c=\frac{8}{\pi }\left(\frac{E_J}{E_C}\right)^{1/2};\kappa ^2=\frac{\epsilon _\nu +2E_J}{4E_J};\nu _c1;$$

$`K(\kappa )`$ and $`E(\kappa )`$ are complete elliptic integrals of the first and second kinds. The density of states $`\rho _d(\epsilon _\nu )`$ results from Eq.(5) and is equal to

$$\rho _d(\epsilon _\nu )=\frac{1}{2\pi }\frac{d\nu }{d\epsilon _\nu }=\frac{K(\kappa )}{\pi ^2\omega _m}.$$
(4)
The level separations in (3) are given by the following expression

$$\omega (\epsilon _\nu )=\frac{d\epsilon }{d\nu }=\frac{\pi \omega _m}{2K(\kappa )},\omega _m=2(E_CE_J)^{1/2}.$$
(5)
At $`(E_J/E_C)^{1/2}1`$ we have the relation $`\omega _m2E_J`$. The quantity $`\omega _m`$ determines the maximum splitting of the levels in the Josephson well. In what follows, the states with $`\nu <\nu _c`$ are denoted as ”libration” states.

In the region of $`\nu \nu _c`$ the $`\nu (\epsilon )`$ dependence and density of states are determined by the following expressions

$$\nu (\epsilon )=\frac{4}{\pi ^2}\left(\frac{E_J}{E_C}\right)^{1/2}\kappa E(\kappa ^1);\rho _c(\epsilon )=\frac{1}{2\pi ^2}\frac{\kappa ^1K(\kappa ^1)}{\omega _m}.$$
(6)
The states (6) with $`\nu >\nu _c`$ are named as the self-trapping states.

Eqs.(3), (6) imply that

$$\begin{array}{cc}& \epsilon (\nu )2E_J+\omega _m\nu ;\mathrm{\hspace{0.33em}\hspace{0.33em}\hspace{0.33em}1}\nu \nu _c,\hfill \\ & \epsilon (\nu )4\pi ^2E_C\nu ^2;\nu \nu _c.\hfill \end{array}$$
(7)
At the same time, it is easy to show that $`d^2\epsilon (\nu )/d\nu ^2<0`$ at $`\nu <\nu _c`$ and $`d^2\epsilon (\nu )/d\nu ^2>0`$ at $`\nu >\nu _c`$. At $`\nu =\nu _c`$ the curve $`\epsilon (\nu )`$ has a flex point.

Since the energy is conserved, the state with a given $`\nu `$ value can be realized by means of the definition of the initial values of $`(\mathrm{\Delta }N)_0`$ and $`\varphi (0)`$. Namely, $`\epsilon (\nu )=E_C(\mathrm{\Delta }N)_0^22E_J\mathrm{cos}\varphi (0)`$. Supposing that $`\varphi (0)=0`$, we obtain the following relation between $`\nu `$ and $`(\mathrm{\Delta }N)_0`$

$$\epsilon (\nu )=2E_J+E_C(\mathrm{\Delta }N)_0^2$$
(8)
Using Eqs.(7), (8), we arrive at the expressions

$$\begin{array}{cc}& \nu =\frac{1}{2\pi }|(\mathrm{\Delta }N)_0|;\text{for the self-trapping (sf) states,}\hfill \\ & \nu =\left(\frac{E_C}{\omega _m}\right)(\mathrm{\Delta }N)_0^2;\text{ for the “libration” (l) states.}\hfill \end{array}$$
(9)
Combining (4), (6), we come to the depedence

$$\rho _{d,c}(\epsilon )\omega _m^1\mathrm{ln}|1\epsilon /2E_J|^1;\epsilon 2E_J^\pm $$
(10)
Thus, there appears a new logarithmic singularity at the boundary separating the libration and self-trapping spectra.

## III The interaction of the Josephson states and condensate oscillations

In this and next sections it is shown that the spectrum of the system can change drastically due to the interaction between the excited Josephson states (7)-(9) with large enough $`\nu `$ values and oscillations of the condensate.

The interaction can be realized by the following mechanisms. First, the interaction is secured if we allow for the dependence of $`E_C`$ in Eq.(8) on the atom displacements. The latter are generated by the condensate oscillation. Second, the interaction can be realized by applying two-photon traveling-wave laser pulse with the Rabi frequency $`\mathrm{\Omega }`$. The pulse both creates the condensates with the different particle numbers and induces the interaction of atom’s displacements with the excited Josephson states corresponding to the particle number disbalance created by the pulse. General description proposed in what follows is independent on the specific mechanism producing the interaction.

Let us consider the classical states of motion of the condensate. These states may be described in terms of the complex amplitudes $`a^{},a=n^{1/2}e^{\pm i\phi _1}`$, where $`n=a|\widehat{a}^+\widehat{a}|a=|a|^2`$ is the average number of quanta in the coherent state $`|a`$. The variables $`n,\phi _1`$ are canonical. By classical state of motion we mean that it’s number of quanta is very large, $`n1`$. It is convenient to specify the relation between the $`a,a^{}`$ amplitudes and $`\widehat{a},\widehat{a}^+`$ operators by the following way: $`a=N^{1/2}\widehat{a}`$. At this, the commutator of $`a,a^{}`$ is equal to zero with macroscopic accuracy: $`[a,a^{}]=1/N0`$. The Hamiltonian of the motion states can be written in the form: $`Nϵ(n)`$.

For the quasiclassical Josephson states with $`\nu 1`$, the $`c_\nu ,c_\nu ^{}`$ amplitudes may be written in the form $`c_\nu =\nu ^{1/2}e^{i\phi _2}`$. However, it is convenient to rewrite $`\epsilon (\nu )`$, $`c_\nu `$ in terms of the variable $`x=|(\mathrm{\Delta }N)_0|/N^{1/2}1`$. Combining this inequality with the requirement that $`x=|(\mathrm{\Delta }N)_0|/N1`$ we arrive at the conditions for the $`x`$ values:

$$1xN^{1/2}.$$
(11)
Using Eqs.(7),(8), we find that $`\nu =\nu (x)`$, $`\epsilon (\nu )=N(E_J/N+E_Cx^2)`$. In general case, the $`\nu (x)`$ dependence is implicit. It is determined by Eqs.(3), (6), (8). However, in the particular cases of the ”libration” ($`\epsilon (\nu )E_J`$) and self-trapping ($`\epsilon (\nu )E_J`$) states, the relations between $`\nu `$ and $`(\mathrm{\Delta }N)_0`$ can be represented in a simple form, as is seen from (9). By means of Eqs.(8), (9), we come to the following expressions

$$\begin{array}{cc}\hfill H_0N& \epsilon _0(n,x)=N[\frac{E_J}{N}+ϵ(n)+E_Cx^2]\hfill \\ & =N[\frac{E_J}{N}+ϵ(n)+\epsilon _J(x)],\hfill \end{array}$$
(12)
$$\begin{array}{cc}& c_\nu ^{(sf)}=N^{1/4}x^{1/2}e^{i\phi _2}N^{1/4}c_x^{(sf)};\hfill \\ & c_\nu ^{(l)}=N\left(\frac{E_C}{\omega _m}\right)^{1/2}xe^{i\phi _2}N\left(\frac{E_C}{\omega _m}\right)^{1/2}c_x^{(l)}.\hfill \end{array}$$
(13)
For any mechanism producing the interaction between two subsystems, it can be written in the form of a multiple Fourier series in $`\phi _1`$, $`\phi _2`$

$$\begin{array}{cc}& H_{int}=N\underset{k_1k_2}{}[g_{k_1k_2}(N)a^{k_1}c_\nu ^{k_2}+c.c.]\hfill \\ & =N\underset{k_1,k_2}{}[g_{k_1k_2}^{(sf,l)}(N)n^{k_1/2}x^{\alpha k_2}e^{i(k_1\phi _1k_2\phi _2)}+c.c.].\hfill \end{array}$$
(14)
Here $`\alpha _{sf}=1/2,\alpha _l=1`$, $`k_1,k_2`$ are integer. For the sake of simplicity, we disregard the phase-independent interaction. Using Eq.(13), we obtain

$$g_{k_1k_2}^{(sf)}(N)=gN^{1+k_2/4};g_{k_1k_2}^{(l)}(N)=gN^{1+k_2}\left(\frac{E_C}{\omega _m}\right)^{k_2/2}.$$
(15)
The constant $`g`$ is specified by the concrete mechanism producing the interaction. Let us assume that the term with the phase $`\varphi _k^r=k_{1r}\phi _1k_{2r}\phi _2`$, which varies anomalously slowly with time, can be set off in sum (14). It is possible to make under two conditions

$$\begin{array}{cc}\hfill k_{1r}\left(\frac{dϵ(n)}{dn}\right)=k_{2r}& \left(\frac{d\epsilon _J(x)}{dx}\right),ϵ_n^{}\left(\frac{dϵ(n)}{dn}\right),\hfill \\ & x_m=\frac{k_{1r}\epsilon _n^{}}{2k_{2r}E_C}.\hfill \end{array}$$
(16)
$$\left(\frac{d^2H_0}{dx^2}\right)_{x=x_m}(\mathrm{\Delta }x)_{max}\left(\frac{H_{int}}{x}\right)_{x=x_m}.$$
(17)
Eq.(17) is written with taking into account that $`H_0`$, $`H_{int}`$ are functions of a single dynamic variable, e.g., x. The $`ϵ_n^{}`$ quantity defines the level separations of an oscillation spectrum.

As is shown in what follows, the condition (16) is equivalent to that of the minimum of the function $`\epsilon _0(n,x)`$ over x. In turn, when the minimum exists, it can provide the principal contribution into the thermodynamic functions.

In addition, the condition (16) implies that the phase $`\varphi _k^r`$ is an approximate integral of motion in the absence of changing $`x`$ near $`x_m`$: $`d\varphi _k^r/dtk_{1r}H_0/nk_{2r}H_0/x0`$.

The inequality (17) implies that the width of the near-minimum region is large at the characteristic interaction variation scale. From Eqs.(16), (17) one can obtain that time changing the $`\varphi _k^r`$ phase is proportional to $`(d^2H_0/dx^2)_m\mathrm{\Delta }x`$, where $`\mathrm{\Delta }x`$ is the change in $`x`$ near the $`x_m`$ value. The maximum value $`(\mathrm{\Delta }x)_{max}`$ specifies the width of the near-minimum region in that $`d\varphi _k^r/dt\mathrm{\Delta }x`$. The estimation for $`(\mathrm{\Delta }x)_{max}`$ is given in what follows. Thus, the leading term in sum (14) has the form

$$H_{int}^{(r)}=Ng_k^{(sf,l)}(N;n,x)\mathrm{cos}\varphi _k,\varphi _k=k_{1r}\phi _1k_{2r}\phi _2;$$
(18)
All remaining terms in this sum are rapidly oscillating perturbations and will be disregarded in this work. Here and below the index $`k`$ in $`g_k`$ and $`\varphi _k`$ denotes a set of $`k_{1r},k_{2r}`$. One can easy to show that aside from the energy $`H=N[\epsilon _0(n,x)+g_k(N;n,x)\mathrm{cos}\varphi _k)]`$ the system in question has an additional integral of motion: $`n_0=n/k_{1r}+x/k_{2r}`$, $`dn_0/dt=0`$. Owing to this, the condition (16) is equivalent to that of the minimum of $`\epsilon _0(n_0,x)`$ over x at a given $`n_0`$ value, as is mentioned above. Using Eqs.(16)-(18), it is straightforward to write the Hamiltonian $`H_m=H_0+H_{int}^{(r)}`$ near the minimum to the first nonvanishing order in $`\mathrm{\Delta }x`$:

$$H_m=N\left[\epsilon _0(n_m;x_m)+\left(\frac{d^2\epsilon _0}{dx^2}\right)_m(\mathrm{\Delta }x)^2g_{km}\mathrm{cos}\varphi _k\right],$$
(19)
$`n_m=n_0k_{1r}x_m/k_{2r}`$, $`(d^2\epsilon _0/dx^2)_m=2E_c`$, $`g_{km}=g_k^{(sf,l)}(N;n_m,x_m)`$. Terms with the derivatives of $`H_{int}`$ are absent in Eq.(19) due to the condition (17).

Using that $`E_C\omega _0(a/a_0)^{2/5}N^{3/5}`$ in the Thomas-Fermi approximation Dalfovo et al. (1999) (here $`a`$, $`a_0`$ are the scattering and oscillator lengths, respectively), we can represent the range of $`1x_mN^{1/2}`$ in the form

$$\frac{1}{N}\left(\frac{Na}{a_0}\right)^{2/5}\frac{k_{1r}}{k_{2r}}\left(\frac{ϵ_n^{^{}}}{\omega _0}\right)\frac{1}{N^{1/2}}\left(\frac{Na}{a_0}\right)^{2/5}.$$
(20)
As is known Dalfovo et al. (1999), the relation $`(Na/a_0)1`$ occurs in the Thomas-Fermi approximation. However, $`N^{1/2}(Na/a_0)^{2/5}1`$. Owing to this relation, the condition (20) (or, what is the same, the condition (16)) specifies the region of the dense oscillation spectrum, where $`ϵ_n^{^{}}\omega _0`$. Both here and in what follows, we suppose that $`k_{1r}=k_{2r}=1`$ for the sake of simplicity.

From Eq.(19) the $`(\mathrm{\Delta }x)_{max}`$ can be estimate to be $`(\mathrm{\Delta }x)_{max}(g_{km}/E_C)^{1/2}`$. Hence, the condition (17) takes the form:

$$\frac{|x_mn_m|}{x_mn_m}\left(\frac{E_C}{g_{km}(N)}\right)^{1/2},$$
(21)
In what follows (see Eq.(34)), it will be shown that the relation (21) is fulfilled with macroscopic accuracy.

## IV The ground state

At a fixed $`n_0`$ value the principal contribution to the partition function comes from the neighbourhood of the minimum at $`x=x_m`$. The expression for $`Z(n_0;x_m;T)`$ is equal to

$$\begin{array}{cc}& Z(n_0;x_m;T)=const\underset{\mathrm{}}{\overset{\mathrm{}}{}}d(\mathrm{\Delta }x)\underset{\pi }{\overset{\pi }{}}𝑑\varphi _ke^{\beta H_m(n_0,(\mathrm{\Delta }x),\varphi _k)}\hfill \\ & =\frac{const}{(\beta NE_C)^{1/2}}\mathrm{exp}\left[\beta N\epsilon _0(N;n_m;x_m)+\mathrm{ln}I_0(\beta Ng_{km})\right],\hfill \end{array}$$
(22)
where $`\beta =1/T`$, $`T`$ is temperature, $`I_0(x)`$ is the modified Bessel function. Eq.(22) implies that the expression for the free energy of the system has the form

$$F=N\epsilon _0(n_m;x_m)+\frac{1}{2}T\mathrm{ln}(\beta NE_c)T\mathrm{ln}I_0(\beta Ng_{km}),$$
(23)
Using Eq.(23), we come to the following equation for the $`\overline{n}_m`$ value realizing the minimum of the free energy

$$\left(\frac{d\epsilon _{0m}}{dn_m}\right)_{n_m=\overline{n}_m}=\left(\frac{dg_{km}}{dn_m}\right)_{n_m=\overline{n}_m}\frac{I_1(\beta Ng_{km})}{I_0(\beta Ng_{km})},$$
(24)
here $`I_1(x)=I_0^{}(x)`$. In addition to $`\overline{n}_m`$, the thermodynamic average of $`\mathrm{cos}\varphi _k`$ may be determined from Eqs.(22) or (23). This average is equal to

$$\mathrm{cos}\varphi _k_T=\frac{\mathrm{ln}Z}{(\beta Ng_{km})}=\frac{F}{T(\beta Ng_{km})}=\frac{I_1(\beta Ng_{km})}{I_0(\beta Ng_{km})},$$
(25)
The order parameters $`\overline{n}_m`$, $`\mathrm{cos}\varphi _k_T`$ describe new coherent state. There is the bound state of the $`\overline{n}_m`$ oscillation quanta and Josephson state generated by the initial disbalance of the particle number that corresponds to the $`x_m`$ value. In addition, this state has the equilibrium phase coherence factor $`\mathrm{cos}\varphi _k_T`$. The $`\overline{n}_m0`$ value provides the equilibrium distortion of the condensate shape. The aboveobtained equations imply that the shape distortion is self-consistently coupled to the $`x_m`$ value defining the equilibrium initial disbalance of the particle number.

At $`T=0`$ the $`\overline{n}_m`$ value realizes the minimum of the thermodynamic energy

$$E_m=N[\epsilon _0(N;n_m,x_m)g_{km}(N;n_m,x_m)].$$
(26)
To determine $`\overline{n}_m(T=0)`$, it is suitable to use the following consideration. As is well known, the level separations $`ϵ_n^{}`$ change slowly in dependence on $`n`$ within the range of dense (quasiclassical) spectrum. For this reason, we can suppose that $`ϵ_n^{^{}}const\omega _b\omega _0`$. On this assumption, the $`x_m`$ value does not depend on $`n_m`$ and the $`\overline{n}_m`$ quantity is equal to

$$\overline{n}_m^{1/2}=\frac{g^{(sf,l)}(N)x_m^\alpha }{2\omega _b},\alpha _{sf,l}=\frac{1}{2},1,$$
(27)
(It is worth noting that $`(^2E_m/n_m^2)>0`$). Both here and in what follows, the denotation $`g_{k=1}^{(sf,l)}(N)g^{(sf,l)}(N)`$ is used. Taking Eq.(27) into account, one readily gets

$$\overline{g}_{1m}^{(sf)}=\frac{g^{(sf)2}(N)x_m}{2\omega _b},\overline{g}_{1m}^{(l)}=\frac{g^{(l)2}(N)x_m^2}{2\omega _b},$$

$$\begin{array}{cc}& \frac{E_m^{(sf)}}{N}=\frac{E_J}{N}\frac{g^{(sf)2}(N)}{8E_C}\left(1\frac{8E_c^2x_m^2}{g^{(sf)2}(N)}\right),\hfill \\ & \frac{E_m^{(l)}}{N}=\frac{E_J}{N}\frac{g^{(l)2}(N)x_m}{16E_C}\left(1\frac{16E_c^2x_m}{g^{(l)2}(N)}\right).\hfill \end{array}$$
(28)
The expressions for energies imply, first, that we obtain the minimum in the region of the dense enough oscillation spectrum, which satisfies the condition (20). The minimum corresponds to the formation of the bound state for the $`\overline{n}_m,x_m`$ values. Second, as is seen from Eq.(28), the absolute minimum of $`E_m`$ can be realized within the ranges of

$$1x_m<\frac{g^{(sf)}(N)}{8E_C};\mathrm{\hspace{0.33em}\hspace{0.33em}\hspace{0.33em}1}x_m<\frac{g^{(l)2}(N)}{16E_C^2}.$$
(29)
These conditions are met when the interaction matrix elaments $`g^{(sf,l)}(N)`$ are large enough. Let us estimate the condensate parameters that are required for existence of the absolute minimum. In the Thomas-Fermi approximation, the inequality $`g^{(sf)}(N)E_C`$ occurs provided the total particle number is not very large, namely,

$$N\frac{g}{\omega _0}\left(\frac{a_0}{a}\right)^{8/3}.$$
(30)
In turn, the relation $`g^{(l)2}(N)E_C^2`$ holds true within the range of

$$N^{0.1}\frac{g^2}{\mathrm{\Omega }^{1/2}\omega _0^{3/2}}\left(\frac{a_0}{a}\right)^{0.6}1,$$
(31)
Here we use that $`E_J=\mathrm{\Omega }N`$. The condition (31) is fulfilled for all admissible parameters, if $`g^2/(\mathrm{\Omega }^{1/2}\omega _0^{3/2})1`$. It should be emphasized that the right-side inequalities in Eq.(29) are much stronger than the condition $`x_mN^{1/2}`$.

The transition to the state with $`\overline{n}_m0`$, $`\mathrm{cos}\varphi _k_T0`$ is second-order. By imposing that $`\overline{n}_m0`$ at the transition temperature, one gets from Eq.(24).

$$T_c^{(sf,l)}=\frac{g^{(sf,l)2}(N)x_m^\alpha }{\omega _b}N.$$
(32)
The dependences of the transition temperatures (32) on the total particle number are given by the expressions

$$T_c^{(sf)}N^{0.2},T_c^{(l)}N^{0.8}.$$
(33)
Thus, the transition temperature $`T_c^{(sf)}`$ has the macroscopic smallness in comparison with the $`T_c^{(l)}`$ temperature. Along with the conditions (30), (31), this fact implies that the libration Josephson state forms the bound state with the condensate oscillation rather than the self-trapping state.

## V Concluding remarks

We have found that the interaction between the Josephson and oscillation states results in new coherent ground state of the double-well trapped condensate. There is the bound state of the anharmonic condensate oscillation and excited Josephson state. The latter is specified by the definite initial disbalance of the particle number. The $`x_m`$ value defining the disbalance is self-consistently coupled to the equilibrium number of the oscillation quanta entering the bound state. Along with the $`\overline{n}_m`$ order parameter, this state has the equilibrium phase coherence factor $`\mathrm{cos}\varphi _k_T`$, where $`\varphi _k`$ is the relative phase of the Josephson and oscillation modes.

It should be emphasized that the bound state arises in the neighbourhood of the minimum of $`\epsilon _0(n_0;x)`$ over $`x`$, where the conditions (16), (17) are met. The condition (16) specifies the type of the interaction between the Josephson and oscillation degrees of freedom. In addition, it imposes the definite restrictions on the spectrum of the oscillation states which may effectively interact with the Josephson degrees of freedom. Namely, the density of oscillation states should be large enough in order to satisfy the conditions (20). For instance, in an asymmetrical double-well potential the energy $`ϵ(n)`$ of it’s classical oscillation states has three branches. Two branches have energies of $`ϵ_{1,2}(n)V_b`$, here $`V_b`$ is the barrier height. There are the maximums at $`ϵ_{1,2}(n_{max})=V_b`$. The third branch has both the energy $`ϵ_3(n)V_b`$ and the minimum at $`ϵ_3(n_{min})=V_b`$. Thus, the regions of the dense spectrum exist in the neighbourhood of the extremums of the functions $`ϵ_i(n)`$ ($`i=1,2,3`$). In the other words, in the neighbourhood of barrier top.

As is shown in section III, the condition (17) can be represent in the form (21). Substituting the aboveobtained expressions for $`\overline{n}_m`$, $`\overline{g}_{1m}`$ into (21), we find that it takes the form

$$|x_mn_m|^{(sf,l)}\left(\frac{g_{1m}^{(sf,l)}(N)}{E_C}\right)^{1/2}.$$
(34)
This inequality is satisfied with macroscopic accuracy.

It is worth noting that the conditions (30), (31) of the absolute minimum of the thermodynamic energy are experimentally controlled by means of either the interaction matrix element $`g`$ or the particle number $`N`$.

The existence of the bound state generates the equilibrium distortion of the condensate shape specified by the $`\overline{n}_m`$ value. This mechanism can provide the experimental detection of the excited Josephson states. The latter can be observed by changing the condensate shape. In addition, the phase transition to the new ground state occurs at $`T=T_c`$, where transition temperatures $`T_c^{(sf,l)}`$ are defined by Eq.(32).

###### Acknowledgements.

I am grateful to Yu.M.Kagan and L.A.Maksimov for helpful discussions. This work was supported by the Russian Foundation for Basic Research.