Electron and boson clusters in confined geometries: Symmetry breaking in quantum dots and harmonic traps
See allHide authors and affiliations

Edited by A. Welford Castleman, Jr., Pennsylvania State University, University Park, PA, and approved January 26, 2006 (received for review October 20, 2005)
Abstract
We discuss the formation of crystalline electron clusters in semiconductor quantum dots and of crystalline patterns of neutral bosons in harmonic traps. In a first example, we use calculations for two electrons in an elliptic quantum dot to show that the electrons can localize and form a molecular dimer. The calculated singlet–triplet splitting (J) as a function of the magnetic field (B) agrees with cotunneling measurements with its behavior reflecting the effective dissociation of the dimer for large B. Knowledge of the dot shape and of J(B) allows determination of the degree of entanglement. In a second example, we study strongly repelling neutral bosons in twodimensional harmonic traps. Going beyond the Gross–Pitaevskii (GP) meanfield approximation, we show that bosons can localize and form polygonalringlike crystalline patterns. The total energy of the crystalline phase saturates in contrast to the GP solution, and its spatial extent becomes smaller than that of the GP condensate.
Explorations of the sizedependent evolution of the properties of materials are at the frontier of modern condensed matter and materials research. Indeed, investigations of clusters containing a finite, well defined number of elementary building units (atoms, molecules, electrons, or other elementary constituents) allow investigations of the transition from the atomic or molecular regime to the finite nanoaggregate domain, and ultimately of the convergence with increasing size to the condensed phase, extended system category. Moreover, investigations of clusters provide opportunities for discovery of novel properties and phenomena that are intrinsic properties of finite systems, distinguishing them from bulk materials (1).
Commonly, studies of materials clusters involve atoms or molecules interacting through electrostatic or electromagnetic potentials, with the heavier nuclei being the “structural skeleton” and the much lighter electrons serving as the “glue” that binds the atoms together. In this article we focus on novel, somewhat exotic, types of clusters. In particular, we discuss clusters of electrons in manmade (artificial) quantum dots (QDs) created through lithographic and gatevoltage techniques at semiconductor interfaces, and clusters of neutral atoms in traps under conditions that may relate to formation of Bose–Einstein condensates (BECs). We illustrate that these cluster systems reveal interesting emergent physical behavior arising from spontaneous breaking of spatial symmetries; symmetry breaking is defined as a circumstance where a lower energy solution of the Schrödinger equation is found that is characterized by a lower symmetry than that of the Hamiltonian of the system. Such symmetry breaking is exhibited through the formation of clusters of localized electrons (often called Wigner molecules, or WMs) in 2D QDs [see Fig. 1 and in particular our discussion of twoelectron WMs below (2)]. Symmetry breaking is also manifested in the transition (3), induced by increasing the interatomic repulsive interaction strength, of the BEC state of neutral atoms confined by a parabolic 2D trap to a crystalline cluster state.
TwoElectron WMs
Electron localization leading to formation of molecularlike structures (the aforementioned WMs) within a single circular 2D QD at zero magnetic field (B) has been theoretically predicted to occur (4–11), as the strength of the e–e repulsive interaction relative to the zeropoint energy increases, as expressed through an increasing value of the Wigner parameter R
_{W}, defined as R
_{W} = Z
^{2}
e
^{2}/(ℏω_{0}
l
_{0}), with l
_{0} =
Because of the finite size, WMs are expected to show new behavior that differs from the classical Wigner crystal familiar from solid state physics. The limit of a classical Wigner crystal is expected to be reached for a higher number of electrons N and very large R _{W} (4). In the following we use the term “Wigner molecule” even in the case of only two localized electrons. For this case the WM exhibits close analogies to an H_{2} natural molecule, as described below.
Here, we focus on a twoelectron (2e) WM, in light of the current experimental effort (19, 20) aiming at implementation of a spinbased (21) solidstate quantum logic gate that employs two coupled oneelectron QDs (double dot). We present an exact diagonalization (EXD) and an approximate (generalized Heitler–London, GHL) microscopic treatment for two electrons in a single elliptic QD specified by the parameters of a recently investigated experimental device (22). While formation of WMs in circular QDs requires weak confinement (that is, small ω_{0} in the expression for R _{W} given above), and thus large dots of lower densities (so that the interelectron repulsion dominates), we show that formation of such WMs is markedly enhanced in highly deformed (e.g., elliptic) dots due to their lower symmetry. The calculations provide a good description of the measured J(B) curve [the singlet–triplet (ST) splitting] when screening (23, 24) due to the metal gates and leads is included (in addition to the dielectric constant of the semiconductor, GaAs). In particular, our results reproduce the salient experimental findings pertaining to the vanishing of J(B) for a finite value of B ∼ 1.3 T [associated with a change in sign of J(B) indicating a ST transition], as well as the flattening of the J(B) curve after the ST crossing. These properties, and in particular the latter one, are related directly to the formation of an electron molecular dimer and its effective dissociation for large magnetic fields. The effective dissociation of the electron dimer is most naturally described through the GHL approximation, and it is fully supported by the more accurate, but physically less transparent, EXD.
Of special interest for quantum computing (21) is the degree of entanglement exhibited by the twoelectron molecule in its singlet state. Entanglement is a purely quantum mechanical phenomenon in which the quantum state of two or more objects cannot be described independently of each other, even when the individual objects are spatially separated. The highest degree of entanglement occurs at full separation, as discussed in the celebrated Einstein, Podolsky, and Rosen paper (25). Electrons confined in a QD are not necessarily spatially separated from each other, and consequently their degree of entanglement may be lower than the maximal one, as shown by us below.
Here, in relation to the microscopic calculations, we investigate two different measures of entanglement. The first, known as the concurrence (𝒞) for two indistinguishable fermions (26, 27), has been used in the analysis of the experiment in ref. 22 [this measure is related to the operational cycle of a twospinqubit quantum logic gate (26, 27)]. The second measure, referred to as the von Neumann entropy (𝒮) for indistinguishable particles, was developed in ref. 28 and used in ref. 29. We show that the present wavefunctionbased methods, in conjunction with the knowledge of the dot shape and the J(B) curve, enable theoretical determination of the degree of entanglement, in particular for the elliptic QD of ref. 22. The increase in the degree of entanglement (for both measures) with stronger magnetic fields correlates with the dissociation of the 2e molecule. This supports the experimental assertion (22) that cotunneling spectroscopy can probe properties of the electronic wave function of the QD, and not merely its lowenergy spectrum. Our methodology can be straightforwardly applied to other cases of strongly interacting devices, e.g., double dots with strong interdot tunneling.
Clusters of Neutral Bosons in Harmonic Traps
BECs in harmonic traps (30, 31) are normally associated with weakly interacting neutral atoms, and their physics is described adequately by the Gross–Pitaevskii (GP) meanfield theory (32). Lately, however, experimental advances in controlling the interaction strength (33–36) permit the production of novel bosonic states in the regime of strong interparticle repulsions. Theoretical efforts motivated by this capability include studies of the Bose–Hubbard model (37, 38), and investigations about the “fermionization” limit of a 1D gas of trapped impenetrable bosons (39–41), often referred to as the Tonks–Girardeau (TG) regime (39, 42). Here, we address the problem of strongly repelling (impenetrable) bosons in higher dimensions. In particular, we discuss 2D interacting bosons in a circular harmonic trap, with the extension to 3D systems being straightforward. To this end, we use computational methods that go beyond the GP method.
We explore the transition from a BEC (diffuse cloud) state to a crystalline phase, in which the trapped localized bosons form crystalline patterns. At the meanfield level, these crystallites are static and are portrayed directly in the singleparticle densities. After restoration of rotational symmetry, the singleparticle densities are circularly symmetric, and thus the crystalline symmetry becomes “hidden”; however, it can be revealed in the conditional probability distribution (CPD) (anisotropic pair correlation), P(r, r _{0}), which expresses the probability of finding a particle at r given, that the “observer” (i.e., reference point) is riding on another particle at r _{0} (10, 16).
Methods
TwoElectron QD: Microscopic Treatment.
The Hamiltonian for two 2D interacting electrons is where the last term is the Coulomb repulsion, κ is the dielectric constant, and r _{12} = r _{1} − r _{2}. H(r) is the singleparticle Hamiltonian for an electron in an external perpendicular magnetic field B and an appropriate confinement potential. When positiondependent screening is included, the last term in Eq. 1 is modified by a function of r _{12} (see below). For an elliptic QD, the singleparticle Hamiltonian is written as where T = (p − e A/c)^{2}/2m*, with A = 0.5(−By, Bx, 0) being the vector potential in the symmetric gauge. m* is the effective mass, and p is the linear momentum of the electron. The second term is the external confining potential; the last term is the Zeeman interaction with g* being the effective g factor, μ_{B} is the Bohr magneton, and s is the spin of an individual electron.
The GHL method for solving the Hamiltonian (1) consists of two steps. In the first step, we solve selfconsistently the ensuing unrestricted Hartree–Fock (UHF) equations allowing for lifting of the doubleoccupancy requirement (imposing this requirement gives the restricted HF method, RHF). For the S _{z} = 0 solution, this step produces two singleelectron orbitals u _{L,R}(r) that are localized left (L) and right (R) of the center of the QD [unlike the RHF method that gives a single doublyoccupied elliptic (and symmetric about the origin) orbital]. At this step, the manybody wave function is a single Slater determinant Ψ_{UHF} (1 ↑, 2 ↓) ≡  u _{L}(1 ↑)u _{R}(2 ↓)〉 made out of the two occupied UHF spinorbitals u_{L}(1 ↑) ≡ u_{L}(r _{1})α (1) and u _{R}(2 ↓) ≡ u _{R}(r _{2})β(2), where α(β) denotes the up (down) [↑ (↓)] spin. This UHF determinant is an eigenfunction of the projection S _{z} of the total spin S = s _{1} + s _{2} but not of S ^{2} (or the parity spacereflection operator).
In the second step, we restore the broken parity and totalspin symmetries by applying to the UHF determinant the projection operator (43, 44) P ^{s} ^{,t} = 1 ∓ ϖ_{12}, where the operator ϖ_{12} interchanges the spins of the two electrons; the upper (minus) sign corresponds to the singlet. The final result is a GHL twoelectron wave function Ψ_{GHL} ^{s,t}(r _{1}, r _{2}) for the groundstate singlet (index s) and firstexcited triplet (index t), which uses the UHF localized orbitals, where χ ^{s} ^{,t} = (α(1)β(2) ∓ α(2)β(1)) is the spin function for the 2e singlet and triplet states. The general formalism of the 2D UHF equations and of the subsequent restoration of broken spin symmetries can be found in refs. 11 and 43–45.
The use of optimized UHF orbitals in the GHL is suitable for treating single elongated QDs. The GHL is equally applicable to double QDs with arbitrary interdottunneling coupling (43, 44). In contrast, the Heitler–London treatment (46) (known also as valence bond), where nonoptimized “atomic” orbitals of two isolated QDs are used, is appropriate only for the case of a double dot with small interdottunneling coupling (21).
The orbitals u _{L,R}(r) are expanded in a real Cartesian harmonicoscillator basis, i.e., where the index j ≡ (m, n) and ϕ _{j} (r) = X_{m} (x)Y_{n} (y), with X_{m} (Y_{n} ) being the eigenfunctions of the 1D oscillator in the x(y) direction with frequency ω _{x} (ω _{y} ). The parity operator 𝒫 yields 𝒫X_{m} (x) = (−1) ^{m}X_{m} (x), and similarly for Y_{n} (y). The expansion coefficients C _{j} ^{L,R} are real for B = 0 and complex for finite B. In the calculations, we use K = 79, yielding convergent results.
In the exactdiagonalization method, the manybody wave function is written as a linear superposition over the basis of noninteracting twoelectron determinants, i.e., where ψ(1; i) = ϕ _{i} (1 ↑) if 1 ≤ i ≤ K and ψ(1; i) = ϕ _{i} _{−K}(1 ↓) if K + 1 ≤ i ≤ 2K [and similarly for ψ(2; j)]. The total energies E _{j} ^{L,R} and the coefficients Ω_{ij} ^{s,t} are obtained through a “brute force” diagonalization of the matrix eigenvalue equation corresponding to the Hamiltonian in Eq. 1 . The EXD wave function does not immediately reveal any particular form, although our calculations below show that it can be approximated by a GHL wave function in the case of the elliptic dot under consideration.
TwoElectron QD: Measures of Entanglement.
To calculate the concurrence 𝒞 (26, 27), one needs a decomposition of the GHL wave function into a linear superposition of orthogonal Slater determinants. Thus, one needs to expand the nonorthogonal u ^{L,R}(r) orbitals as a superposition of two other orthogonal ones. To this effect, we write u ^{L,R}(r) ∝ Φ^{+}(r) ± ξΦ^{−}(r), where Φ^{+}(r) and Φ^{−}(r) are the parity symmetric and antisymmetric (along the x axis) components in their expansion given by Eq. 4 . Subsequently, with the use of Eq. 3 , the GHL singlet can be rearranged as follows: where the socalled interaction parameter (27), η = ξ^{2}, is the coefficient in front of the second determinant. Knowing η allows a direct evaluation of the concurrence of the singlet state, since 𝒞 ^{s} = 2η/(1 + η^{2}) (27). Note that Φ^{+}(r) and Φ^{−}(r) are properly normalized. It is straightforward to show that η = (1 − S _{LR})/(1 + S _{LR}), where S _{LR} (with S _{LR} ≤ 1) is the overlap of the original u ^{L,R}(r) orbitals.
For the GHL triplet, one obtains an expression independent of the interaction parameter η, i.e., which is a maximally (𝒞 ^{t} = 1) entangled state. Note that underlying the analysis of the experiments in ref. 22 is a conjecture that wave functions of the form given in Eqs. 6 and 7 describe the two electrons in the elliptic QD.
To compute the von Neumann entropy, one needs to bring both the EXD and the GHL wave functions into a diagonal form [the socalled “canonical form” (28, 47)], i.e., with the Φ(i)s being appropriate spin orbitals resulting from a unitary transformation of the basis spin orbitals ψ(j)s (see Eq. 5 ); only terms with z_{k} ≠ 0 contribute. The upper bound M can be smaller (but not larger) than K (the dimension of the singleparticle basis); M is referred to as the Slater rank. One obtains the coefficients of the canonical expansion from the fact that the z_{k} ^{2} are eigenvalues of the hermitian matrix Ω^{†}Ω [Ω (see Eq. 5 ) is antisymmetric]. The von Neumann entropy is given by with the normalization Σ_{k=1} ^{M} z_{k} ^{2} = 1. Note that the GHL wave functions in Eqs. 6 and 7 are already in canonical form, which shows that they always have a Slater rank of M = 2. One finds 𝒮_{GHL} ^{s} = log_{2}(1 + η^{2}) − η^{2} log_{2}(η^{2})/(1 + η^{2}), and 𝒮_{GHL} ^{t} = 1 for all B. For large B, the overlap between the two electrons of the dissociated dimer vanishes, and thus η → 1 and 𝒮_{GHL} ^{s} → 1.
Neutral Repelling Bosons in Harmonic Traps: GP Equation.
Meanfield symmetry breaking for bosonic systems has been discussed earlier in the context of twocomponent condensates, where each species is associated with a different space orbital (48, 49). We consider here one species of bosons but allow each particle to occupy a different space orbital φ _{i} (r _{i} ). The permanent Φ _{N} 〉 = Perm[φ_{1}(r _{1}), …, φ _{N} (r _{N} )] serves as the manybody wave function of the unrestricted Bose–Hartree–Fock (UBHF) approximation. This wave function reduces to the GP form with the restriction that all bosons occupy the same orbital φ_{0}(r), i.e., Φ_{N} ^{GP}〉 = Π_{i=1} ^{N}φ_{0}(r _{i} ), and φ_{0}(r) is determined selfconsistently at the restricted Bose–Hartree–Fock (RBHF) level via the equation (50) Here, V(r _{1}, r _{2}) is the twobody repulsive interaction, which is taken to be a contact potential, V _{δ} = gδ(r _{1} − r _{2}), for neutral bosons. The singleparticle Hamiltonian is given by H _{0}(r) = − ℏ^{2}▿^{2}/(2m) + mω_{0} ^{2} r ^{2}/2, where ω_{0} characterizes the harmonic confinement.
Neutral Repelling Bosons in Harmonic Traps: Symmetry Breaking.
We simplify the solution of the UBHF problem by considering explicit analytic expressions for the space orbitals φ _{i} (r _{i} ). In particular, since the bosons must avoid occupying the same position in space to minimize their mutual repulsion, we take all of the orbitals to be of the form of displaced Gaussians, namely, φ _{i} (r _{i} ) = π^{−1/2}σ^{−1} exp[−(r _{i} − a _{i} )^{2}/(2σ^{2})]. The positions a _{i} describe the vertices of concentric regular polygons, with both the width σ and the radius a = a _{i}  of the regular polygons determined variationally through minimization of the total energy E _{UBHF} = 〈Φ _{N} HΦ _{N} 〉/〈Φ _{N} Φ _{N} 〉, where H = Σ_{i=1} ^{N} H _{0}(r _{i} ) + Σ_{i<j} ^{N} V(r _{i} , r _{j} ) is the manybody Hamiltonian.
With the above choice of localized orbitals, the unrestricted permanent Φ _{N} 〉 breaks the continuous rotational symmetry. However, the resulting energy gain becomes substantial for stronger repulsion. Controlling this energy gain (the strength of correlations) is the ratio R _{δ} between the strength of the repulsive potential and the zeropoint kinetic energy. Specifically, for a 2D trap, one has R _{δ} = gm/(2π ℏ^{2}) for a contact potential.
Neutral Repelling Bosons in Harmonic Traps: Restoration of Broken Symmetry.
Although the optimized UBHF permanent Φ _{N} 〉 performs exceptionally well regarding the total energies of the trapped bosons, in particular in comparison with the restricted wave functions (e.g., the GP anzatz), it is still incomplete. Indeed, due to its localized orbitals, Φ _{N} 〉 does not preserve the circular (rotational) symmetry of the 2D manybody Hamiltonian H. Instead, it exhibits a lower pointgroup symmetry, i.e., a C _{2} symmetry for N = 2 and a C _{5} one for the (1, 5) structure of N = 6 (see below). As a result, Φ _{N} 〉 does not have a good total angular momentum. This paradox is resolved through a postHartree–Fock step of restoration of broken symmetries via projection techniques (45, 51), yielding a new wave function Ψ_{N,L} ^{PRJ}〉 with a definite angular momentum L, that is where Φ _{N} (γ)〉 is the original UBHF permanent having each localized orbital rotated by an azimuthal angle γ, with L being the total angular momentum. The projection yields wave functions for a whole rotational band. Note that the projected wave function Ψ_{N,L} ^{PRJ}〉 in Eq. 11 may be regarded as a superposition of the rotated permanents Φ _{N} (γ)〉, thus corresponding to a “continuousconfigurationinteraction” solution.
Here, we are interested in the projected groundstate (L = 0) energy, which is given by
Results and Discussion
TwoElectron QD.
To model the experimental elliptic QD device, we take, following ref. 22, ℏω _{x} = 1.2 meV and ℏω _{y} = 3.3 meV. The effective mass of the electron is taken as m* = 0.067 m_{e} (GaAs). Since the experiment did not resolve the lifting of the triplet degeneracy caused by the Zeeman term, we take g* = 0. Using the twostep method, we calculate the GHL ST splitting J _{GHL}(B) = E _{GHL} ^{s}(B) − E _{GHL} ^{t}(B) as a function of the magnetic field in the range 0 ≤ B ≤ 2.5 T. Screening of the e–e interaction due to the metal gates and leads must be considered to reproduce the experimental J(B) curve. This screening can be modeled, to first approximation, by a positionindependent adjustment of the dielectric constant κ (52). Indeed, with κ = 22.0 (instead of the GaAs dielectric constant; i.e., κ = 12.9), good agreement with the experimental data is obtained (see Fig. 2). In particular, we note the ST crossing for B ≈ 1.3 T and the flattening of the J(B) curve beyond this crossing.
We have also explored, particularly in the context of the EXD treatment, a positiondependent screening using the functional form, (e ^{2}/κr _{12})[1 − (1 + 4d ^{2}/r _{12} ^{2})^{−1/2}], proposed in ref. 24, with d as a fitting parameter. The J _{EXD}(B) result for d = 18.0 nm is depicted in Fig. 2 (dotted line) and is in very good agreement with the experimental measurement.
The singletstate electron densities from the GHL and the EXD treatments at B = 0 and B = 2.5 T are displayed in Fig. 3. These densities illustrate the dissociation of the electron dimer with increasing magnetic field. The asymptotic convergence (beyond the ST point) of the energies of the singlet and triplet states, i.e., [J(B) → 0 as B → ∞], is a reflection of the dissociation of the 2e molecule, since the groundstate energy of two fully spatially separated electrons (zero overlap) does not depend on the total spin.
In contrast, the singletstate RHF electron densities fail to exhibit formation of an electron dimer for all values of B. This underlies the failure of the RHF method to describe the behavior of the experimental J(B) curve. In particular, J _{RHF}(B = 0) has the wrong sign, while J _{RHF}(B) diverges for high B as is the case for the RHF treatment of double dots (see ref. 43).
For the GHL singlet, using the overlaps of the left and right orbitals, we find that starting with η = 0.46 (𝒞 ^{s} = 0.76) at B = 0, the interaction parameter (singletstate concurrence) increases monotonically to η = 0.65 (𝒞 ^{s} = 0.92) at B = 2.5 T. At the intermediate value corresponding to the ST transition (B = 1.3 T), we find η = 0.54 (𝒞 ^{s} = 0.83).
Our B = 0 theoretical results for η and 𝒞 ^{s} are in remarkable agreement with the experimental estimates (22) of η = 0.5 ± 0.1 and 𝒞 ^{s} = 0.8, which were based solely on conductance measurements below the ST transition (i.e., near B = 0). We note that, for the RHF, 𝒞 _{RHF} ^{s} = 0, since a single determinant is unentangled for both the two measures considered here.
Since the EXD singlet has obviously a Slater rank M > 2, the definition of concurrence is not applicable to it. The von Neumann entropy for the EXD singlet (𝒮_{EXD} ^{s}) is displayed in Fig. 4, along with that (𝒮_{GHL} ^{s}) of the GHL singlet. 𝒮_{EXD} ^{s} and 𝒮_{GHL} ^{s} are rather close to each other for the entire B range, and it is remarkable that both remain close to unity for large B, although the maximum allowed mathematical value is log_{2}(K) [as aforementioned, we use K = 79; i.e., log_{2}(79) = 6.3]; this maximal value applies for both the EXD and GHL approaches. The saturation of the entropy for large B to a value close to unity reflects the dominant (and roughly equal at large B) weight of two configurations in the canonical expansion (see Eq. 8 ) of the EXD wave function, which are related to the two terms (M = 2) in the canonical expansion of the GHL singlet (Eq. 6 ). This is illustrated by the histograms of the z _{k} ^{s}^{2} coefficients for B = 1.3 T at the top of Fig. 4. These observations support the GHL approximation, which is computationally less demanding than the exact diagonalization, and can be used easily for larger N.
Neutral Repelling Bosons in Harmonic Traps.
In Fig. 5, we display as a function of the parameters R _{δ} the total energies for N = 6 bosons calculated at several levels of approximation. In both cases the lowest UBHF energies correspond to a (1,5) crystalline configuration, namely one boson is at the center and the rest form a regular pentagon of radius a. Observe that the GP total energies are slightly lower than the E _{RBHF} ^{G} ones; however, both exhibit an unphysical behavior since they diverge as R _{δ} → ∞. This behavior contrasts sharply with that of the UHF energies, E _{UBHF} and PRJ (see below), which saturate as R _{δ} → ∞; in fact, a value close to saturation is achieved already for R _{δ} ∼ 10. We have checked that for all cases with N = 2 − 7, the total energies exhibit a similar behavior. For a repulsive contact potential, the saturation of the UBHF energies is associated with the ability of the trapped bosons (independent of N) to minimize their mutual repulsion by occupying different positions in space, and this is one of our central results. For N = 2, the two bosons localize at a distance 2a apart to form an antipodal dimer. For N ≤ 5, the preferred UBHF crystalline arrangement is a single ring with no boson at the center [usually denoted as (0, N)]. N = 6 is the first case having one boson at the center [designated as (1, N − 1)], and the (0, 6) arrangement is a higher energy isomer. The structural parameters (e.g., the width of the Gaussian orbitals and the radii of the polygonal ring, calculated via the UBHF method) show a saturation behavior similar to that illustrated above for the energy of the system (3). In contrast, the width of the condensate cloud (i.e., the GP solution) diverges with increasing repulsion strength (R _{δ}).
The saturation found here for 2D trapped bosons interacting through strong repelling contact potentials is an illustration of the “fermionization” analogies that appear in strongly correlated systems in all three dimensionalities. Indeed such energy saturation has been shown for the TG 1D gas (42, 39) and has also been discussed for certain 3D systems [i.e., three trapped bosons (53) and an infinite boson gas (54)]. Saturation of the energy and the length of the trapped atom cloud (and thus of the interparticle distance) has been measured recently for the 1D TG gas (see in particular figures 3 and 4 in ref. 36, and compare with the similar trends predicted here for the 2D case in Fig. 5).
For N = 6 2D bosons, Fig. 5 shows that the E _{0} ^{PRJ} energies share with the UBHF ones the saturation property for the case of a contactpotential repulsion. However, the projection brings further lowering of the total energies compared with the UBHF ones. [The projected ground state is always lower in energy than the original brokensymmetry one (55).] Thus, for strong interactions (large values of R _{δ}) the restorationofbrokensymmetry step yields an excellent approximation of both the exact manybody wave function and the exact total energy.
The transformations of the singleparticle densities (displayed in Fig. 6 for N = 6 neutral bosons interacting via a contact potential and R _{δ} = 25) obtained from application of the successive approximations provide an illustration of the twostep method of symmetry breaking with subsequent symmetry restoration. Indeed, the GP singleparticle density (Fig. 6 a) is circularly symmetric, but the UBHF one (Fig. 6 b) explicitly exhibits a (1,5) crystalline configuration. After symmetry restoration (Fig. 6 c) the circular symmetry is reestablished, but the singleparticle density is radially modulated unlike the GP density. In addition, the crystalline structure in the projected wave function is now hidden; however, it can be revealed through the use of the conditional probability distribution (10, 16) (see Fig. 6 d), which resembles the (crystalline) UBHF singleparticle density, but with one of the humps on the outer ring missing (where the observer is located). In particular, P(r _{0}, r _{0}) ≈ 0 and the boson associated with the observer is surrounded by a “hole” similar to the exchangecorrelation hole in electronic systems. This is another manifestation of the “fermionization” of the strongly repelling 2D bosons. However, here as in the 1D TG case (42, 39), the vanishing of P(r _{0}, r _{0}) results from the impenetrability of the bosons. For the GP condensate, the conditional probability distribution is independent of r _{0}, i.e., P _{GP}(r, r _{0}) ∝ φ_{0}(r)^{2}, reflecting the absence of any space correlations.
It is of importance to observe that the radius of the BEC (GP case; Fig. 6 a) is significantly larger than the actual radius of the stronglyinteracting crystalline phase (projected wave function; Fig. 6 c). This is because the extent of the crystalline phase saturates, whereas that of the GP condensate grows with no bounds as R _{δ} → ∞. Such dissimilarity in size (between the condensate and the strongly interacting phase) has been also predicted (40) for the trapped 1D Tonks–Girardeau gas and indeed observed experimentally (36). In addition, the 2D singleparticle momentum distributions for neutral bosons have a onehump shape with a maximum at the origin (a behavior exhibited also by the trapped 1D TG gas). The width of these momentum distributions versus R _{δ} increases and saturates to a finite value, whereas that of the GP solution vanishes as R _{δ} → ∞.
Summary
In this article we explored symmetrybreaking transitions predicted to occur in confined fermionic and bosonic systems when the strength of the interparticle repulsive interactions exceeds an energy scale that characterizes the degree of confinement. For two electrons in an elliptic QD, we predicted formation and effective dissociation (with increasing magnetic field) of an electron dimer, which is reflected in the behavior of the computed ST splitting, J(B), that agrees well (Fig. 2) with measurements (22).
Furthermore, we showed that, from a knowledge of the dot shape and of J(B), theoretical analysis along the lines introduced here allows probing of the correlated groundstate wave function and determination of its degree of entanglement. This presents an alternative to the experimental study where determination of the concurrence used conductance data (22). Such information is of interest to the implementation of spinbased solidstate quantum logic gates.
For the case of 2D trapped bosonic clusters, we found with increasing repulsive twobody interaction localization of the bosons in the trap, resulting in formation of crystalline patterns made of polygonal rings; although we have focused here on repulsive contact interactions, similar results were obtained also for a Coulomb repulsion (3).
These results provide the impetus for experimental efforts to access the regime of strongly repelling neutral bosons in two dimensions. To this end we anticipate that extensions of methodologies developed for the recent realization of the Tonks–Girardeau regime in one dimension [using a finite small number of trapped ^{87}Rb and optical lattices, with a demonstrated wide variation of R _{δ} from 5 to 200 (35) and from 1 to 5 (36)] will prove most promising. Control of the interaction strength via the use of the Feshbach resonance may also be considered (33).
Acknowledgments
This work was supported by U.S. Department of Energy Grant FG0586ER45234 and National Science Foundation Grant DMR0205328.
Footnotes
 ^{†}To whom correspondence may be addressed. Email: constantine.yannouleas{at}physics.gatech.edu or uzi.landman{at}physics.gatech.edu

Author contributions: C.Y. and U.L. designed research, performed research, and wrote the paper.

Conflict of interest statement: No conflicts declared.

This paper was submitted directly (Track II) to the PNAS office.
 Abbreviations:
 BEC,
 Bose–Einstein condensate;
 EXD,
 exact diagonalization;
 GHL,
 generalized Heitler–London;
 GP,
 Gross–Pitaevskii;
 QD,
 quantum dot;
 ST,
 singlet–triplet;
 UHF,
 unrestricted Hartree–Fock;
 WM,
 Wigner molecule.
Abbreviations:
 © 2006 by The National Academy of Sciences of the USA
References

↵
 Landman U.

↵
 Yannouleas C. ,
 Landman U.
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Elzerman J. M. ,
 Tarucha S. ,
 Kouwenhoven L. P.
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Xe L. ,
 Bester G. ,
 Zunger A.
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵

↵
 Kinoshita G. T. ,
 Wenger T. ,
 Weiss D. S.
 ↵
 ↵

↵
 Girardeau M. D. ,
 Wright E. M.
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
Citation Manager Formats
Article Classifications
 Physical Sciences
 Physics