Next Article in Journal
Approach for 3D Cultural Relic Classification Based on a Low-Dimensional Descriptor and Unsupervised Learning
Next Article in Special Issue
Real Space Triplets in Quantum Condensed Matter: Numerical Experiments Using Path Integrals, Closures, and Hard Spheres
Previous Article in Journal
Effect of Different Local Vibration Frequencies on the Multiscale Regularity of Plantar Skin Blood Flow
Previous Article in Special Issue
Structural and Thermodynamic Peculiarities of Core-Shell Particles at Fluid Interfaces from Triangular Lattice Models
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Ultracold Bosons on a Regular Spherical Mesh

Dipartimento di Scienze Matematiche ed Informatiche, Scienze Fisiche e Scienze della Terra, Università degli Studi di Messina, Viale F. Stagno d’Alcontres 31, 98166 Messina, Italy
Submission received: 22 October 2020 / Revised: 8 November 2020 / Accepted: 11 November 2020 / Published: 13 November 2020
(This article belongs to the Special Issue Statistical Mechanics and Thermodynamics of Liquids and Crystals)

Abstract

:
Here, the zero-temperature phase behavior of bosonic particles living on the nodes of a regular spherical mesh (“Platonic mesh”) and interacting through an extended Bose-Hubbard Hamiltonian has been studied. Only the hard-core version of the model for two instances of Platonic mesh is considered here. Using the mean-field decoupling approximation, it is shown that the system may exist in various ground states, which can be regarded as analogs of gas, solid, supersolid, and superfluid. For one mesh, by comparing the theoretical results with the outcome of numerical diagonalization, I manage to uncover the signatures of diagonal and off-diagonal spatial orders in a finite quantum system.

1. Introduction

Gases of ultracold bosonic atoms loaded in an optical lattice provide the unique opportunity to study quantum many-body effects under controlled conditions [1,2]. To a very good approximation, the atoms can be described by a Bose–Hubbard (BH) Hamiltonian [3] whose parameters can be tuned by laser light [4,5]. By changing the configuration of the lasers, many lattice geometries can be explored, making optical lattices a powerful and versatile tool.
In a system at zero temperature ( T = 0 ), thermal fluctuations are frozen out and quantum fluctuations prevail. These microscopic fluctuations can induce phase transitions in the ground state of a many-body system, driven by a non-thermal control parameter, such as chemical potential, magnetic field, or chemical composition. As a concrete example, consider a dilute gas of bosons at temperatures low enough that a Bose-Einstein condensate is formed. The condensate is described by a wave function consisting of every particle in the same state spread over the entire volume of the system, and typically exhibits superfluidity. An interesting situation appears when the condensate is subject to a lattice potential in which particles can move from one lattice site to the next only by tunneling. If the lattice potential is increased smoothly, the system remains in the condensed phase as long as the repulsion between atoms is small compared to the tunnel coupling (assuming, by the way, that the range of the repulsion is much smaller than the lattice spacing). In this regime, where the tunneling term dominates the Hamiltonian, a delocalized wave function still minimizes the total energy of the many-body system. When the strength of the repulsion becomes large compared to the tunnel coupling, the total energy is made minimum when each lattice site is occupied by the same number of atoms; as a result, phase coherence is lost and the system becomes insulating. The addition of a longer-range repulsion will make the phase behavior richer, with the possibility of a non-superfluid density wave and a supersolid ground state where crystalline order coexists with superfluid behavior (see, e.g., Ref. [6]). Experimentally, a way to prepare ultracold gases of long-range interacting bosons is to use atoms (such as chromium [7] or dysprosium [8]) and molecules having a large magnetic or electric dipole moment.
The usual BH model predicts a T = 0 phase transition from a superfluid phase to a Mott insulator phase as the ratio of the hopping matrix element between adjacent sites (t, in absolute terms) to the on-site interaction (U) is reduced. The overall number density of particles is controlled via a chemical-potential parameter μ . As μ grows at fixed t and U, the lattice becomes increasingly filled with particles, but this can only occur outside the Mott regions since the insulator phase is incompressible [3]. The BH model has been studied in many lattice geometries and with several techniques (mean-field theory [9,10,11], perturbation theory [12,13,14,15,16], and quantum Monte Carlo simulation [17,18], to name but the most commonly employed). When a further repulsion V is introduced between nearest-neighbor atoms (“extended BH model”), new phases may arise, in primis a supersolid phase [19,20,21,22,23,24,25,26,27]. Another variant of the BH model is hard-core bosons, where site occupancy is restricted to zero or one, corresponding to taking the U limit [6,28,29,30,31,32].
I hereafter present the results of yet another investigation of the extended BH model, now choosing a finite graph as hosting space for bosons. Even though clearcut phase transitions (i.e., thermodynamic singularities) cannot occur in a few-particle system, a convenient choice of boundary conditions may alleviate the difference with an infinite system, making the study of a finite quantum system valuable anyway. A practical solution is to use spherical boundary conditions (SBCs), which have often been exploited in the past to discourage long-range triangular ordering at high density [33,34,35,36,37,38]. On the other hand, SBCs make it possible to observe novel forms of ordering, viz. into regular polyhedral structures, that are simply unknown to Euclidean space—see Refs. [39,40,41]. An added value of a spherical mesh of points is the possibility to vary the site coordination while keeping the overall geometry strictly two-dimensional. Bosons confined in spherical (bubble) traps have been produced experimentally [42,43] and are going to be studied soon under microgravity conditions [44,45]. In the present study, the extended BH model is considered on a finite mesh of points homogeneously distributed on the unit sphere, i.e., coincident with the vertices of a regular polyhedron inscribed in the sphere (we may call it a Platonic mesh). Despite consisting of a finite number of nodes, a regular spherical mesh shares an important feature in common with an ordinary lattice, namely all sites are equivalent; for this reason, the phase behavior of particles living on a Platonic mesh would not deviate much from that of an infinite system (this intuition will be checked in a particular case).
The plan of this paper is the following. After introducing the models in Section 2, the choice of the underlying mesh is discussed in detail, giving priority to those regular grids where a subset of nodes form a mesh that is also regular. Then, in Section 3, a mean-field (MF) analysis of the ground state is carried out for all the models considered; in one case, the indications of theory are validated against the results of exact diagonalization (Section 4). From this comparison, we find the artifacts of the MF approximation as applied to a finite system. Finally, some concluding remarks are given in Section 5.

2. Particle Models on a Spherical Mesh

2.1. Classical Models

To illustrate the main idea, let it initially be considered a system of classical point particles defined on the sites of a cubic mesh stretched over the surface of the unit sphere. Each of the eight nodes of the mesh has three nearest neighbors (NN, at chord distance r = 2 / 3 ) and three next-nearest neighbors (NNN, at chord distance 2 2 / 3 ). For simplicity, assume that the occupancy n i of site i can only be 0 or 1 ( i = 1 , , 8 ). Finally, choose the system Hamiltonian to be H [ n ] = V i , j n i n j with V > 0 (each NN pair in the sum is counted only once). The grand potential Ω ( T , μ ) of this system is
Ω = β 1 ln Ξ with Ξ = { n } e β ( H [ n ] μ i n i ) ,
where β is the inverse temperature. At T = 0 , the formula is simpler:
Ω = min { n } H [ n ] μ i n i .
Any of the points of absolute minimum in (2) represent the actual equilibrium state of the system for that μ . For particles living in a continuous space, minimization of H μ N is carried out among a selection of crystalline states that are thought to be relevant based on symmetry considerations (see, e.g., Ref. [46]). In the present case, where the total number of microstates is small ( 2 8 = 256 ), the T = 0 grand potential can be determined exactly for each μ by a scrutiny of all possible energies. While the mesh is empty for μ < 0 and completely filled with particles for μ > 3 V , in the interval from 0 to 3 V only half of the nodes are occupied, and these fall at the vertices of a regular tetrahedron. This twofold-degenerate state with checkerboard order can be viewed as the finite-size counterpart of a crystalline phase.
The rationale behind the choice of a cubic mesh is now clear: by introducing a repulsion between occupied NN sites, we promote the occurrence of a Platonic “crystal”, i.e., the regular tetrahedron (“CT model”, see Figure 1 left). There is only one other possibility to obtain a non-frustrated crystal-like state at T = 0 , which is using a regular dodecahedral mesh. We shall see that (i) by discouraging the occupancy of NN sites, a cubic “crystal” is stabilized for sufficiently small μ > 0 (“DC model”, Figure 1 center); (ii) if the repulsion is extended to embrace NNN sites too, then a tetrahedral “crystal” is stabilized in a range of positive μ values (“DT model”, Figure 1 right).
For a system of classical hard-core particles defined on a regular dodecahedral mesh, the total number of microstates is 2 20 = 1,048,576. By direct inspection, we see that the minimum- Ω “phase” for the DC model is the empty mesh for μ < 0 , a cube for 0 < μ < 3 V / 2 (a fivefold-degenerate state, since there are 5 ways to form a cube with the vertices of a regular dodecahedron), and the complement of a cube (“co-cube” in the following) for 3 V / 2 < μ < 3 V ; finally, for μ > 3 V the mesh is completely filled. For the DT model, the mesh is empty for μ < 0 , filled with particles located at the vertices of a regular tetrahedron for 0 < μ < 3 V / 2 (a tenfold degenerate state), and completely filled for μ > 9 V . For μ between 3 V / 2 and 9 V , a different ground state exists for each even number of occupied sites in the range 6 to 16, none of which corresponds to a simple geometric arrangement.
I have plotted in Figure 2 the evolution with temperature of a few properties of the CT and DC models. In addition to the grand potential Ω , the figure shows the total number of occupied sites (N), the total energy (E), and the order parameters for tetrahedral ( S t ) and cubic order ( S c ). The latter quantities are defined so as to discriminate the Platonic “phase” from the other T = 0 states. A proper order parameter should be insensitive to the orientation of the polyhedron, hence it can only depend on the relative angles between the vertices vector radii departing from the sphere center: for a regular tetrahedron all these angles are equal to α = arccos { 1 / 3 } (with sin α = 2 2 / 3 ), while they are π α , α , and π for a cube. With the idea to penalize configurations that do not match the wanted structure, my choice of the OPs is the following:
S t = 1 k t i < j cos θ i j + 1 3 2 , N = 4 0 , N 4
and
S c = 1 k c i < j sin θ i j sin θ i j 2 3 2 2 , N = 8 0 , N 8 ,
where the constants k t and k c are chosen, following the advice in Ref. [47], so that S vanishes for a random distribution of angles:
k t = 3 8 and k c = 36 7 59 π 128 2 1 .
Looking at Figure 2, it is clear that the only singularities are found at T = 0 (the system is finite) while all cusps and jumps are smoothened for T > 0 .
It is worth considering how the T = 0 phase diagram of the CT and DC models gets modified when the occupancy of a node is allowed to take any value. For simplicity, as relevant T = 0 states only the “cluster crystals” [48,49,50] originated from the previously identified ground states are considered. Call U > 0 the on-site energy and V the NN repulsion. The grand Hamiltonian now reads H [ n ] = U / 2 i n i ( n i 1 ) + V i , j n i n j μ i n i . For the CT model, the grand potential of a tetrahedral “cluster crystal” with n particles per site is Ω n ( 4 ) = 2 n ( n 1 ) U 4 n μ , while the grand potential of the “cluster crystal” with n particles per site is Ω n ( 8 ) = 4 n ( n 1 ) U + 12 n 2 V 8 n μ . For a given μ , the most stable “phase” corresponds to the minimum Ω . Up to μ = 0 the stable “phase” is still the empty mesh. As μ grows further, the site occupancy increases monotonically within each of the families Ω n ( 4 ) and Ω n ( 8 ) , but the exact sequence of stable “phases” depends on U / V . I show three cases in Figure 3. For U = 5 V (right panel), the behavior for small μ > 0 recalls that found in the hard-core limit; however, as μ grows, each site becomes increasingly populated, and a whole sequence of cluster states is found. The opposite occurs for U = 2 V (left panel), where only the cluster states with checkerboard order are stabilized for μ > 0 . A curious situation occurs for U = 3 V , where the competition between different cluster states becomes so stringent that “crystals” of the two families coexist at regular intervals along the μ axis (center panel).
As a matter of fact, hybrid ground states with features intermediate between those of the tetrahedral and cubic cluster states may also occur. For example, I have checked for U = 5 V that a microstate with two overlapping particles at the vertices of one tetrahedron and one particle at the vertices of the other tetrahedron is the most stable ground state in the μ interval from 8 V to 11 V .
Moving to the DC model, we have three different families of “cluster phases”: the cubic family with grand potential Ω n ( 8 ) = 4 n ( n 1 ) U 8 n μ ; the co-cubic family with grand potential Ω n ( 12 ) = 6 n ( n 1 ) U + 6 n 2 V 12 n μ ; and the family of those cluster states where every site is filled with the same number of particles, with grand potential Ω n ( 20 ) = 10 n ( n 1 ) U + 30 n 2 V 20 n μ . Without delving into details, the phase behavior as a function of U is similar to the CT model, with multiple coexistence between all three types of cluster crystals now occurring for U = 2 V .

2.2. Quantum Models

The previous analysis has served to set the stage for the forthcoming study of the extended BH model on a regular spherical mesh; only the quantum versions of the CT and DC models (denoted QCT and QDC, respectively) are considered below.
The system is defined by the following grand Hamiltonian:
H = t i , j a i a j + a j a i + U 2 i n i ( n i 1 ) + V i , j n i n j μ i n i .
In the present investigation, H describes a system of bosons on a spherical mesh of M sites, a i , a i are bosonic field operators ( i = 1 , , M ), and n i = a i a i is a number operator. Moreover, t 0 is the hopping amplitude between NN sites, U > 0 is the on-site repulsion, and V > 0 is the strength of the NN repulsion. In the hard-core limit U , the site occupancy is restricted to zero or one and the U term in (6) vanishes; it is only this limit that is treated hereafter.
In the BH model on a standard lattice, at T = 0 we observe an insulator–superfluid transition with increasing t for every positive μ . The addition of V may stabilize (depending on the lattice) a density wave at low t, as well as a supersolid phase at the boundary between the insulator and superfluid phases. These features are also present in the hard-core limit, as reported, e.g., in Refs. [30,32].

3. MF Investigation

MF theory is the method of choice when a new many-body problem is attacked; it has been frequently applied for continuous quantum systems as well, as an effective means to identify the ground states and quantum transitions between them (cf. Refs. [51,52,53]). Various versions of MF theory exist in discrete space; here the so-called decoupling approximation [11,29,32] is employed, which gives the advantage of a fully analytic treatment of the problem. Clearly, the accuracy of MF theory would be questionable when applied for a finite quantum system, even one lacking a boundary surface. This will make urgent an assessment of MF theory against exact results, which however is delayed until the next Section.
In the decoupling approximation, the two-site terms in the Hamiltonian (6) are linearized so that the eigenvalue problem becomes effectively one-site. This is accomplished by first writing the exact identity
a i a j = a i a j + a i a j a i a j + δ a i δ a j
with δ a i = a i a i (for i = 1 , , M ), and then neglecting the last term. Similarly,
n i n j n i n j + n i n j n i n j .
The averages a i ϕ i and n i ρ i (either ground-state averages or thermal averages) are to be determined self-consistently; ϕ i and ρ i (with 0 ρ i 1 ) represent the superfluid order parameter and the local density for site i, respectively. For hard-core bosons, we readily obtain H H MF with
H MF = t i F i a i + F i a i F i ϕ i + V 2 i 2 R i n i R i ρ i μ i n i .
In Equation (9), F i = j NN i ϕ j and R i = j NN i ρ j are sums over the nearest neighbors of the i-th site. For a bipartite mesh, sites are either of type A or B, hence the unknown parameters are four, namely ϕ A , ϕ B , ρ A , ρ B . The simplification is obvious: rather than working on a Fock space of dimensionality 2 M , for a bipartite mesh the basis states are just four, namely 0 , 0 , 0 , 1 , 1 , 0 , and 1 , 1 , corresponding to the possible occupancies of a pair of A and B sites.

3.1. QCT Model

For the QCT model, the mesh is bipartite and formed by two tetrahedral sub-meshes with four points each. A point of grid A has three neighbors, all belonging to grid B, and vice versa. Hence,
F A = 3 ϕ B , F B = 3 ϕ A , R A = 3 ρ B , and R B = 3 ρ A .
Using a subscript A (B) for operators relative to a single A (B) site, the MF Hamiltonian reads:
H MF = E 0 12 t ϕ B a A + ϕ B a A + ϕ A a B + ϕ A a B + 4 ( 3 V ρ B μ ) n A + 4 ( 3 V ρ A μ ) n B
with
E 0 = 12 t ( ϕ A ϕ B + ϕ A ϕ B ) 12 V ρ A ρ B .
On the 4-vector basis B = 0 , 0 , 0 , 1 , 1 , 0 , 1 , 1 , the Hamiltonian (11) is represented by a 4 × 4 Hermitian matrix:
H MF = E 0 12 t ϕ A 12 t ϕ B 0 12 t ϕ A E 0 + 4 ( 3 V ρ A μ ) 0 12 t ϕ B 12 t ϕ B 0 E 0 + 4 ( 3 V ρ B μ ) 12 t ϕ A 0 12 t ϕ B 12 t ϕ A E 0 + 12 V ( ρ A + ρ B ) 8 μ .
Non-superfluid phases have ϕ A = ϕ B = 0 . In this case, the matrix becomes diagonal, meaning that the basis vectors are all (grand-)energy eigenstates. From the self-consistency equations ρ A = n A and ρ B = n B , it readily follows that ρ A = ρ B = 0 for 0 , 0 (empty mesh), ρ A = 0 and ρ B = 1 for 0 , 1 (tetrahedral-B “crystal”), ρ A = 1 and ρ B = 0 for 1 , 0 (tetrahedral-A “crystal”), and ρ A = ρ B = 1 for 1 , 1 (fully occupied mesh). The eigenvalue gives the grand potential Ω , which is zero for the empty mesh, 4 μ for the twofold-degenerate tetrahedral “crystal”, and 12 V 8 μ for the filled mesh. These grand potentials are exactly the same as in the CT model.
Superfluid and supersolid “phases” have non-zero, possibly distinct, complex values of ϕ A and ϕ B . For sure, ϕ A and ϕ B have equal phases since only the magnitude of the order parameter can be spatially modulated; hence, the arbitrary phase can be taken as zero. Without loss of generality, we may assume ϕ A and ϕ B to be positive quantities. We shall see in the next section that the exact eigenstates of H for t > 0 do not distinguish between A and B; hence, our search can be restricted to homogeneous (superfluid) solutions: ϕ A = ϕ B = ϕ and ρ A = ρ B = ρ . We are thus led to the following characteristic equation (with E 0 = 24 t ϕ 2 12 V ρ 2 ):
E 0 λ 12 t ϕ 12 t ϕ 0 12 t ϕ E 0 + 4 ( 3 V ρ μ ) λ 0 12 t ϕ 12 t ϕ 0 E 0 + 4 ( 3 V ρ μ ) λ 12 t ϕ 0 12 t ϕ 12 t ϕ E 0 + 8 ( 3 V ρ μ ) λ = 0 .
Doing the simplifications, we arrive at
( λ b ) 2 ( λ a ) ( λ 2 b + a ) 4 u 2 = 0
with a = E 0 , b = E 0 + 4 ( 3 V ρ μ ) , and u = 12 t ϕ . The minimum root of (15) is
E = b ( a b ) 2 + 4 u 2 = 4 6 t ϕ 2 3 V ρ 2 + 3 V ρ μ ( 3 V ρ μ ) 2 + 36 t 2 ϕ 2 .
A real eigenvector of E is
ψ E = E 2 b + a 2 u 0 , 0 + 0 , 1 + 1 , 0 + E a 2 u 1 , 1 ,
or, in explicit terms,
ψ E = 3 V ρ μ + ( 3 V ρ μ ) 2 + 36 t 2 ϕ 2 6 t ϕ , 1 , 1 , ( 3 V ρ μ ) + ( 3 V ρ μ ) 2 + 36 t 2 ϕ 2 6 t ϕ .
Now imposing the conditions
ρ A , B = ψ E | n A , B | ψ E ψ E | ψ E and ϕ A , B = ψ E | a A , B | ψ E ψ E | ψ E ,
we arrive at the two coupled equations
1 2 ρ = 3 V ρ μ ( 3 V ρ μ ) 2 + 36 t 2 ϕ 2 and 1 = 3 t ( 3 V ρ μ ) 2 + 36 t 2 ϕ 2 ,
which are easily solved to give
ρ = μ + 3 t 3 V + 6 t and ϕ = ( μ + 3 t ) ( 3 V + 3 t μ ) 3 V + 6 t .
The above ϕ solution only exists provided that 3 t < μ < 3 V + 3 t , which is a necessary condition for the existence of the superfluid. Plugging these ρ and ϕ in Equation (16), we finally obtain the superfluid grand potential:
Ω SF = 4 ( μ + 3 t ) 2 3 V + 6 t .
This outcome can also be obtained from the general relation (see Equations (6)–(8))
H MF = 2 t i , j ϕ i ϕ j + V i , j ρ i ρ j μ i ρ i .
For the superfluid, we have ϕ i = ϕ and ρ i = ρ , with ϕ and ρ given by Equation (21); hence,
H MF = 4 6 t ϕ 2 + 3 V ρ 2 2 μ ρ = 4 ( μ + 3 t ) 2 3 V + 6 t .
By comparing the grand potentials of all the “phases”, we arrive at the ground-state diagram in Figure 4. Along the straight lines μ = 3 t and μ = 3 V + 3 t , separating the superfluid from the insulator “phases” at small and large μ , the μ derivative of the grand potential (= N ) is continuous. Instead, along the line
μ = 3 V 2 ± 3 2 V 2 4 t 2 ,
separating the “crystalline” lobe from the superfluid, the average number of occupied sites shows a jump discontinuity. Only at the lobe vertex ( V / 2 , 3 V / 2 ) the tetrahedral-superfluid transition is continuous. We will see in Section 4 to what extent the indications of MF theory are accurate, considering that the cubic mesh consists of 8 sites only.

3.2. QDC Model

The decoupling approximation for the QDC model works similarly as for the QCT model. Again, the starting point is Equation (9) and the mesh is partitioned into a cube (A, 8 sites) and a co-cube (B, 12 sites). We have:
F A = 3 ϕ B , F B = 2 ϕ A + ϕ B , R A = 3 ρ B , and R B = 2 ρ A + ρ B .
The MF Hamiltonian then reads:
H MF = E 0 24 t ( ϕ B a A + ϕ B a A ) 12 t ( 2 ϕ A + ϕ B ) a B + ( 2 ϕ A + ϕ B ) a B + 4 ( 6 V ρ B 2 μ ) n A + 4 ( 6 V ρ A + 3 V ρ B 3 μ ) n B
with
E 0 = 24 t ( ϕ A ϕ B + ϕ A ϕ B ) + 12 t | ϕ B | 2 24 V ρ A ρ B 6 V ρ B 2 .
On the B basis, the Hamiltonian (27) is represented by the matrix:
H MF = E 0 12 t ( 2 ϕ A + ϕ B ) 24 t ϕ B 0 12 t ( 2 ϕ A + ϕ B ) E 0 + 4 ( 6 V ρ A + 3 V ρ B 3 μ ) 0 24 t ϕ B 24 t ϕ B 0 E 0 + 4 ( 6 V ρ B 2 μ ) 12 t ( 2 ϕ A + ϕ B ) 0 24 t ϕ B 12 t ( 2 ϕ A + ϕ B ) E 0 + 4 ( 6 V ρ A + 9 V ρ B 5 μ ) .
For phases with ϕ A = ϕ B = 0 , the matrix is diagonal and, like in the QCT model, the basis vectors are all (grand-)energy eigenstates. From the self-consistency conditions ρ A = n A and ρ B = n B it follows that ρ A = ρ B = 0 for 0 , 0 (empty mesh), ρ A = 0 and ρ B = 1 for 0 , 1 (co-cubic “crystal”), ρ A = 1 and ρ B = 0 for 1 , 0 (cubic “crystal”), and ρ A = ρ B = 1 for 1 , 1 (filled mesh). The grand potential Ω is zero for the vacuum, 8 μ for the cubic “crystal”, 6 V 12 μ for the co-cubic “crystal”, and 30 V 20 μ for the filled mesh. These values are the same as for the DC model, hence the same sequence of “phases” as a function of μ is observed in the QDC model for t = 0 .
Moving to phases with ϕ A , ϕ B 0 , I first consider the possibility of a superfluid ( ϕ A = ϕ B = ϕ > 0 and ρ A = ρ B = ρ ). In this case, the characteristic equation takes the form:
a λ u v 0 u b λ 0 v v 0 c λ u 0 v u d λ = 0
with E 0 = 60 t ϕ 2 30 V ρ 2 , a = E 0 , b = E 0 + 12 ( 3 V ρ μ ) , c = E 0 + 8 ( 3 V ρ μ ) , d = E 0 + 20 ( 3 V ρ μ ) , u = 36 t ϕ , v = 24 t ϕ . Observing that a + d = b + c s , the minimum eigenvalue of H MF turns out to be (see Appendix A):
E = 1 2 s ( a b ) 2 + 4 u 2 ( a + b s ) 2 + 4 v 2 = 10 6 t ϕ 2 3 V ρ 2 + 3 V ρ μ ( 3 V ρ μ ) 2 + 36 t 2 ϕ 2 .
Notice that this energy is exactly 5 / 2 times that of the QCT model (see Equation (16)).
Using p = ( b c ) / 2 = ( d a ) / 10 = 2 ( 3 V ρ μ ) and q = u / 3 = v / 2 = 12 t ϕ , the coordinates ( x , y , z , w ) of an eigenvector ψ E of E satisfy the following linear system:
5 p x + 3 q y + 2 q z = 5 p 2 + q 2 x 3 q x + p y + 2 q w = 5 p 2 + q 2 y 2 q x p z + 3 q w = 5 p 2 + q 2 z 2 q y + 3 q z + 5 p w = 5 p 2 + q 2 w ,
implying that ψ E is also an eigenvector of the simpler matrix
5 p 3 q 2 q 0 3 q p 0 2 q 2 q 0 p 3 q 0 2 q 3 q 5 p
with eigenvalue 5 p 2 + q 2 . For t > 0 , one such vector is
ψ E = p + p 2 + q 2 q , 1 , 1 , p p 2 + q 2 q ,
or
ψ E = 3 V ρ μ + ( 3 V ρ μ ) 2 + 36 t 2 ϕ 2 6 t ϕ , 1 , 1 , ( 3 V ρ μ ) + ( 3 V ρ μ ) 2 + 36 t 2 ϕ 2 6 t ϕ .
This is identical to the state (18) describing the superfluid in the QCT model. The reason of this equivalence is that, in both cubic and dodecahedral meshes, every site has three neighbors; hence, in the superfluid, the quantities F A , F B , R A , R B are the same. No surprise, then, if also the expressions of ρ and ϕ turn out to be equal for the two models (the consistency equations are identical). As already commented, the energy of the QDC superfluid is instead a factor 5 / 2 larger (in absolute terms) than in the QCT model:
Ω SF = 10 ( μ + 3 t ) 2 3 V + 6 t .
The main novelty with respect to the QCT model is the existence of a stable supersolid “phase” in the QDC model. To seek for MF solutions having the character of a supersolid, I have first derived the exact equations for the MF parameters ρ A , ρ B , ϕ A , ϕ B , without a priori assuming them to be equal in pairs (this is done in Appendix B). Then, I have solved these equations numerically so as to find the region where the supersolid is more stable than the other phases. This task is made simpler by noting that, thanks to a symmetry property of the equations, from a supersolid solution with μ > 3 V / 2 it is possible to obtain another solution with μ < 3 V / 2 . Indeed, it easily follows from Equation (A15) that, if ( ρ A , ρ B , ϕ A , ϕ B ) is a solution for a certain μ , then ( 1 ρ A , 1 ρ B , ϕ A , ϕ B ) is a solution for 3 V μ . Moreover, for μ = 3 V / 2 the two solutions share the same grand potential. It turns out that the ( t , μ ) region where the supersolid is stable is symmetric about the μ = 3 V / 2 axis, lying across the boundary between the solid “phases” and the superfluid (see Figure 5). To all evidence, the boundary between the supersolid and the superfluid lies at t = 1 / 3 . In fact, the supersolid consists of two distinct “phases”, SS1 below μ = 3 V / 2 (where ρ B < 1 / 2 < ρ A ) and SS2 above μ = 3 V / 2 (where ρ A < 1 / 2 < ρ B ), with ϕ A < ϕ B < 1 / 2 in both. In SS1, cubic sites are more occupied, on average, than co-cubic sites, whereas the opposite occurs in SS2. The two supersolids coexist for μ = 3 V / 2 , for all t in the range 0.25 to 1 / 3 . Hard-core bosons on the triangular lattice have a similar phase diagram [32], but in that case the supersolid extends down to t = 0 and the solid lobes do not overlap each other.
The full T = 0 phase diagram of the QDC model according to the decoupling approximation is reported in Figure 5. Along the straight lines μ = 3 t and μ = 3 V + 3 t , the μ derivative of the grand potential is continuous. Instead, along the curves
μ = 3 t + 6 V ± 6 V 2 V t 6 t 2 5 and μ = 3 t + 9 V ± 6 V 2 V t 6 t 2 5 ,
which respectively separate the cubic and co-cubic regions from the superfluid region, the total number of particles jumps discontinuously. Only at the lobe vertices ( V / 3 , V ) and ( V / 3 , 2 V ) the two “crystal”—superfluid transitions are continuous. The line μ = 3 V / 2 separating the two “crystalline” regions ends at the point ( 3 V / 10 , 3 V / 2 ) where the two boundary curves (37) cross each other. However, this triple point is only metastable since, in the whole region bounded by the blue circuit of Figure 5, the stable phase is actually supersolid.

4. Assessment of MF Theory

Let us again reconsider the QCT model at T = 0 . Due to the relatively small dimensionality of its Hilbert space (=256), the model can be solved numerically, determining (among others) the exact ground state g and its eigenvalue (i.e., the grand potential) as a function of t and μ . This is obtained by representing the grand Hamiltonian H on the Fock basis { x 1 , x 2 , , x 8 } (with x i = 0 or 1 for all i) and diagonalizing the ensuing matrix. An extensive mapping of a few quantum averages enables us to clarify the nature of the “phases” present.
As usual, the sites of the mesh are classified according to what tetrahedral sub-mesh, A or B, they belong to. Then, I compute the average occupancies of A and B sites (corresponding to the MF parameters ρ A and ρ B , respectively), the average values of a A and a B (corresponding to the MF parameters ϕ A and ϕ B , respectively), and the superfluid density ρ SF (see, e.g., Ref. [13,22]). In the present case, the latter quantity reads:
ρ SF 1 8 g a ˜ 0 a ˜ 0 g = 1 64 i , j = 1 8 g a i a j g = 1 8 j = 1 8 g a 1 a j g ,
where a ˜ 0 = ( 1 / 8 ) i = 1 8 a i is the zero-momentum field operator. Notice that, in a large lattice of M sites, a ˜ 0 a ˜ 0 = N 0 is the average number of condensate particles, hence ρ SF = N 0 / M is indeed the condensate density.
As far as the occupancies are concerned, exact diagonalization shows that they are always equal for an A and a B site, with the only exception of t = 0 and 0 < μ < 3 V where the occupancies are as in MF theory (i.e., either n A = 1 and n B = 0 or vice versa). In particular, for t > 0 the equivalence n A = n B holds in the whole ( t , μ ) region pertaining, according to MF theory, to the tetrahedral “phase”. Hence, no spontaneous symmetry breaking does really occur in the QCT model, except for t = 0 . Indeed, for t > 0 the Fourier coefficients of g relative to any pair of Fock states equal by A-B inversion are the same. I show in Figure 6 the average site occupancy in the whole space of parameters. We see that the ( t , μ ) plane is divided in zones where the site occupancy, which overall grows with μ , takes the same constant value. The possible values are k / 8 , for k = 0 , 1 , , 8 . This outcome is not entirely unexpected considering that the Hamiltonian commutes with the N operator, implying that the ground state (which is non-degenerate for t > 0 ) should also be an eigenstate of N. Like in MF theory, the occupancy is zero below μ = 3 t and 1 above μ = 3 V + 3 t ; in a whole region around μ = 3 V / 2 the occupancy is 0.5 for all t, with no abrupt transition from “crystal” to superfluid values.
Another difference with MF theory concerns the ground-state averages of a A and a B : these are identically zero for all t and μ values, which may seem in stark contrast to the behavior of the MF parameter ϕ . In fact, it is not these averages that should be monitored in an exact treatment but rather the ρ SF quantity defined at Equation (38). In the top panel of Figure 7, I have reported the ρ SF values computed along a few iso-t lines; each jump discontinuity of ρ SF occurs in coincidence with one of the site occupancy. For comparison, in the bottom panel of Figure 7, I show the MF values of ϕ 2 along the same lines. We see a clear correlation between the two behaviors, which demonstrates the existence of signatures of superfluidity in a finite quantum system. However, intriguing differences also exist: (i) first observe that ρ SF = 0.125 in the filled mesh, which is an oddity for a perfectly insulating state. In fact, 0.125 is 1 / M for M = 8 , suggesting that this is an artifact of the finite system size. (ii) Another finite-size effect is the non-zero values of ρ SF in the purported tetrahedral region. However, ρ SF is here sufficiently small so that its peaks in the μ gaps between the tetrahedral and insulator regions remain well visible—in MF theory these gaps fall in the superfluid region.
Finally, in Figure 8, I plot the QCT grand potential as a function of μ for a number of t values in the range 0 to 0.5. In the MF curves, cusp singularities are associated with the crossing of first-order transition lines. We see that MF theory systematically overestimates exact values, the more so the larger t is.

5. Conclusions

I have worked out the zero-temperature phase diagram of two systems of hard-core bosons defined on the nodes of a regular spherical mesh. The interaction is of Bose–Hubbard type, with a further repulsion between neighboring particles. Choosing a suitable mesh, bosons are pushed to form a Platonic “crystal” in a range of chemical potentials. In the QCT model, the mesh is cubic and a tetrahedral “crystal” is formed; in the QDC model, the mesh is dodecahedral and a cubic “crystal” is formed instead.
Using a mean-field approximation, I have obtained fully analytic results for the thermodynamic properties of the two models. In addition to a number of insulating phases, both systems also exhibit a superfluid ground state. In the QDC model, triple coexistence between two solids and the superfluid is superseded by the occurrence of a more stable supersolid phase. Clearly, while the predictions of mean-field theory are generally accurate for an infinite Bose–Hubbard system, deviations will unavoidably be observed in a small system, where no true singularity can occur. However, discrepancies are probably less strong for bosons on a regular spherical mesh, which has equivalent sites and is devoid of natural boundaries.
To check this expectation for the QCT model, I have diagonalized its Hamiltonian exactly, finding the ground state and a number of ground-state averages. Overall, the exact T = 0 behavior of the finite system does not depart much from theory, though obviously phase-transition lines are now reduced to simple crossovers; however, clear marks of superfluidity have been detected in the same region of parameters where MF theory predicts it to occur. In conclusion, I hope that the present work may stimulate research aimed at the realization of a new experimental platform for ultracold bosons that is akin to a regular spherical mesh.

Funding

This research received no external funding.

Conflicts of Interest

The author declares no conflict of interest.

Appendix A. On the Solutions to an Algebraic Equation

In this Appendix, I show how to determine the eigenvalues of a matrix like (29), namely the four real solutions λ to the characteristic Equation (30) with a + d = b + c s .
Upon doing the algebra, the characteristic equation turns out to be:
( λ a ) ( λ b ) ( λ s + a ) ( λ s + b ) u 2 ( λ a ) ( λ b ) + ( λ s + a ) ( λ s + b ) v 2 ( λ a ) ( λ s + b ) + ( λ s + a ) ( λ b ) + ( u 2 v 2 ) 2 = 0 .
Let us now reshuffle the lhs of (A1) term by term:
( i ) ( λ a ) ( λ b ) ( λ s + a ) ( λ s + b ) = 1 16 ( 2 λ s ) 2 ( 2 a s ) 2 ( 2 λ s ) 2 ( 2 b s ) 2 = 1 16 ( 2 λ s ) 4 ( 2 a s ) 2 + ( 2 b s ) 2 ( 2 λ s ) 2 + ( 2 a s ) 2 ( 2 b s ) 2 = 1 16 ( 2 λ s ) 4 2 ( a b ) 2 + ( a + b s ) 2 ( 2 λ s ) 2 + ( a b ) 2 ( a + b s ) 2 2 ;
( i i ) u 2 ( λ a ) ( λ b ) + ( λ s + a ) ( λ s + b ) = 1 2 u 2 ( 2 λ s ) 2 + ( a + b s ) 2 ( a b ) 2 ;
( i i i ) v 2 ( λ a ) ( λ s + b ) + ( λ s + a ) ( λ b ) = 1 2 v 2 ( 2 λ s ) 2 + ( a b ) 2 ( a + b s ) 2 .
Plugging Equations (A2)–(A4) in (A1), the latter equation becomes:
( 2 λ s ) 4 2 ( a b ) 2 + ( a + b s ) 2 + 4 ( u 2 + v 2 ) ( 2 λ s ) 2 + ( a b ) 2 ( a + b s ) 2 2 + 8 ( u 2 v 2 ) ( a b ) 2 ( a + b s ) 2 + 16 ( u 2 v 2 ) 2 = 0
or equivalently
( 2 λ s ) 4 2 ( a b ) 2 + ( a + b s ) 2 + 4 ( u 2 + v 2 ) ( 2 λ s ) 2 + ( a b ) 2 ( a + b s ) 2 + 4 ( u 2 v 2 ) 2 = 0 .
Upon defining the two quantities X = ( a b ) 2 + 4 u 2 and Y = ( a + b s ) 2 + 4 v 2 , the lhs of Equation (A6) can be factorized as follows:
( 2 λ s ) 4 2 ( X + Y ) ( 2 λ s ) 2 + ( X Y ) 2 = ( 2 λ s ) 2 ( X + Y ) 2 4 X Y = ( 2 λ s ) 2 X Y + 2 X Y ( 2 λ s ) 2 X Y 2 X Y = ( 2 λ s ) 2 X Y 2 ( 2 λ s ) 2 X + Y 2 .
Hence, two independent second-order equations are obtained from (A6):
( 2 λ s ) 2 = X Y 2 and ( 2 λ s ) 2 = X + Y 2 ,
whose solutions are the searched eigenvalues:
λ 1 = 1 2 s X + Y ; λ 2 = 1 2 s + X Y ; λ 3 = 1 2 s X Y ; λ 4 = 1 2 s + X + Y .
The minimum eigenvalue is solution no. 3 above.

Appendix B. Self-Consistency Conditions in the QDC Model

I here derive the general equations obeyed by MF parameters in the QDC model. These are obtained by first determining the leading eigenvector of the matrix representing the MF grand Hamiltonian on B . Then, the self-consistency conditions are written in the same way as Equation (19).
The MF Hamiltonian is Equation (29) with real ϕ A and ϕ B . Therefore, the eigenvalue equation is as in (30), with a = E 0 = 48 t ϕ A ϕ B + 12 t ϕ B 2 24 V ρ A ρ B 6 V ρ B 2 , b = E 0 + 12 ( 2 V ρ A + V ρ B μ ) , c = E 0 + 8 ( 3 V ρ B μ ) , and d = E 0 + 4 ( 6 V ρ A + 9 V ρ B 5 μ ) ; moreover, u = 12 t ( 2 ϕ A + ϕ B ) and v = 24 t ϕ B . The minimum eigenvalue is like in the first line of Equation (31), but its explicit form in terms of MF parameters is now
E = E 0 + 2 ( 6 V ρ A + 9 V ρ B 5 μ ) 6 ( 2 V ρ A + V ρ B μ ) 2 + 4 t 2 ( 2 ϕ A + ϕ B ) 2 4 ( 3 V ρ B μ ) 2 + 36 t 2 ϕ B 2 .
Using the latter formula, the linear system for the coordinates ( x , y , z , w ) of an eigenvector ψ E of E are:
( p 1 + p 2 ) x + q 1 y + q 2 z = p 1 2 + q 1 2 + p 2 2 + q 2 2 x q 1 x + ( p 1 p 2 ) y + q 2 w = p 1 2 + q 1 2 + p 2 2 + q 2 2 y q 2 x ( p 1 p 2 ) z + q 1 w = p 1 2 + q 1 2 + p 2 2 + q 2 2 z q 2 y + q 1 z + ( p 1 + p 2 ) w = p 1 2 + q 1 2 + p 2 2 + q 2 2 w
with
p 1 = 3 ( 2 V ρ A + V ρ B μ ) ; q 1 = 6 t ( 2 ϕ A + ϕ B ) ; p 2 = 2 ( 3 V ρ B μ ) ; q 2 = 12 t ϕ B .
Therefore, ψ E is also an eigenvector of a matrix simpler than the original one, with eigenvalue p 1 2 + q 1 2 p 2 2 + q 2 2 . When q 1 , q 2 0 , one such vector is:
ψ E = p 1 + p 1 2 + q 1 2 q 1 p 2 + p 2 2 + q 2 2 q 2 , p 2 + p 2 2 + q 2 2 q 2 , p 1 + p 1 2 + q 1 2 q 1 , 1 .
Now we are ready to impose self-consistency, which eventually leads to:
ρ A = 1 2 p 2 2 p 2 2 + q 2 2 ; ρ B = 1 2 p 1 2 p 1 2 + q 1 2 ; ϕ A = q 2 2 p 2 2 + q 2 2 ; ϕ B = q 1 2 p 1 2 + q 1 2 ,
or explicitly:
1 2 ρ A = 3 V ρ B μ ( 3 V ρ B μ ) 2 + 36 t 2 ϕ B 2 ; 1 2 ρ B = 2 V ρ A + V ρ B μ ( 2 V ρ A + V ρ B μ ) 2 + 4 t 2 ( 2 ϕ A + ϕ B ) 2 ; ϕ A = 3 t ϕ B ( 3 V ρ B μ ) 2 + 36 t 2 ϕ B 2 ; ϕ B = t ( 2 ϕ A + ϕ B ) ( 2 V ρ A + V ρ B μ ) 2 + 4 t 2 ( 2 ϕ A + ϕ B ) 2 .
Probably, the simpler method to solve these four coupled non-linear equations in four unknowns is to minimize a suitable non-negative function constructed in such a way as to vanish when the (A15) are all fulfilled. This is easy to do numerically with a computer. By this method, I have drawn the supersolid boundaries in Figure 5.

References

  1. Bloch, I.; Dalibard, J.; Zwerger, W. Many-body physics with ultracold gases. Rev. Mod. Phys. 2008, 80, 885–964. [Google Scholar] [CrossRef] [Green Version]
  2. Amico, L.; Boshier, M.; Birkl, G.; Minguzzi, A.; Miniatura, C.; Kwek, L.-C.; Aghamalyan, D.; Ahufinger, V.; Andrei, N.; Arnold, A.S.; et al. Roadmap on Atomtronics. arXiv 2020, arXiv:2008.04439. [Google Scholar]
  3. Fisher, M.P.A.; Weichman, P.B.; Grinstein, G.; Fisher, D.S. Boson localization and the superfluid-insulator transition. Phys. Rev. B 1989, 40, 546–570. [Google Scholar] [CrossRef] [Green Version]
  4. Jaksch, D.; Bruder, C.; Cirac, J.I.; Gardiner, C.W.; Zoller, P. Cold Bosonic Atoms in Optical Lattices. Phys. Rev. Lett. 1998, 81, 3108–3111. [Google Scholar] [CrossRef]
  5. Greiner, M.; Mandel, O.; Esslinger, T.; Hänsch, T.W.; Bloch, I. Quantum phase transition from a superfluid to a Mott insulator in a gas of ultracold atoms. Nature 2002, 415, 39–44. [Google Scholar] [CrossRef]
  6. van Otterlo, A.; Wagenblast, K.-H.; Baltin, R.; Fazio, R.; Schön, G. Quantum phase transitions of interacting bosons and the supersolid phase. Phys. Rev. B 1995, 52, 16176–16186. [Google Scholar] [CrossRef] [Green Version]
  7. Griesmaier, A.; Werner, J.; Hensler, S.; Stuhler, J.; Pfau, T. Bose-Einstein Condensation of Chromium. Phys. Rev. Lett. 2005, 94, 160401. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  8. Lu, M.; Burdick, N.Q.; Youn, S.H.; Lev, B.L. Strongly Dipolar Bose-Einstein Condensate of Dysprosium. Phys. Rev. Lett. 2011, 107, 190401. [Google Scholar] [CrossRef] [PubMed]
  9. Rokhsar, D.S.; Kotliar, B.G. Gutzwiller projection for bosons. Phys. Rev. B 1991, 44, 10328–10332. [Google Scholar] [CrossRef] [PubMed]
  10. Krauth, W.; Caffarel, M.; Bouchaud, J.-P. Gutzwiller wave function for a model of strongly interacting bosons. Phys. Rev. B 1992, 45, 3137–3140. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  11. Sheshadri, K.; Krishnamurthy, H.R.; Pandit, R.; Ramakhrishnan, T.V. Superfluid and Insulating Phases in an Interacting-Boson Model: Mean-Field Theory and the RPA. Europhys. Lett. 1993, 22, 257–263. [Google Scholar] [CrossRef]
  12. Freericks, J.K.; Monien, H. Phase diagram of the Bose-Hubbard Model. Europhys. Lett. 1994, 26, 545–550. [Google Scholar] [CrossRef] [Green Version]
  13. van Oosten, D.; van der Straten, P.; Stoof, H.T.C. Quantum phases in an optical lattice. Phys. Rev. A 2001, 63, 053601. [Google Scholar] [CrossRef] [Green Version]
  14. Schroll, C.; Marquardt, F.; Bruder, C. Perturbative corrections to the Gutzwiller mean-field solution of the Mott-Hubbard model. Phys. Rev. A 2004, 70, 053609. [Google Scholar] [CrossRef] [Green Version]
  15. dos Santos, F.E.A.; Pelster, A. Quantum phase diagram of bosons in optical lattices. Phys. Rev. A 2009, 79, 013614. [Google Scholar] [CrossRef] [Green Version]
  16. Kübler, M.; Sant’Ana, F.T.; dos Santos, F.E.A.; Pelster, A. Improving mean-field theory for bosons in optical lattices via degenerate perturbation theory. Phys. Rev. A 2019, 99, 063603. [Google Scholar] [CrossRef] [Green Version]
  17. Batrouni, G.G.; Rousseau, V.; Scalettar, R.T.; Rigol, M.; Muramatsu, A.; Denteneer, P.J.H.; Troyer, M. Mott Domains of Bosons Confined on Optical Lattices. Phys. Rev. Lett. 2002, 89, 117203. [Google Scholar] [CrossRef] [Green Version]
  18. Capogrosso-Sansone, B.; Prokof’ev, N.V.; Svistunov, B.V. Phase diagram and thermodynamics of the three-dimensional Bose-Hubbard model. Phys. Rev. B 2007, 75, 134302. [Google Scholar] [CrossRef] [Green Version]
  19. Góral, K.; Santos, L.; Lewenstein, M. Quantum Phases of Dipolar Bosons in Optical Lattices. Phys. Rev. Lett. 2002, 88, 170406. [Google Scholar] [CrossRef] [Green Version]
  20. Kovrizhin, D.L.; Pai, G.V.; Sinha, S. Density wave and supersolid phases of correlated bosons in an optical lattice. Europhys. Lett. 2005, 72, 162–168. [Google Scholar] [CrossRef] [Green Version]
  21. Sengupta, P.; Pryadko, L.P.; Alet, F.; Troyer, M.; Schmid, G. Supersolids versus Phase Separation in Two-Dimensional Lattice Bosons. Phys. Rev. Lett. 2005, 94, 207202. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  22. Yamamoto, K.; Todo, S.; Miyashita, S. Successive phase transitions at finite temperatures toward the supersolid state in a three-dimensional extended Bose-Hubbard model. Phys. Rev. B 2009, 79, 094503. [Google Scholar] [CrossRef] [Green Version]
  23. Pollet, L.; Picon, J.D.; Büchler, H.P.; Troyer, M. Supersolid Phase with Cold Polar Molecules on a Triangular Lattice. Phys. Rev. Lett. 2010, 104, 125302. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  24. Ng, K.-K. Thermal phase transitions of supersolids in the extended Bose-Hubbard model. Phys. Rev. B 2010, 82, 184505. [Google Scholar] [CrossRef] [Green Version]
  25. Iskin, M. Route to supersolidity for the extended Bose-Hubbard model. Phys. Rev. A 2011, 83, 051606(R). [Google Scholar] [CrossRef] [Green Version]
  26. Kimura, T. Gutzwiller study of phase diagrams of extended Bose-Hubbard models. J. Phys. Conf. Ser. 2012, 400, 012032. [Google Scholar] [CrossRef] [Green Version]
  27. Ohgoe, T.; Suzuki, T.; Kawashima, N. Commensurate Supersolid of Three-Dimensional Lattice Bosons. Phys. Rev. Lett. 2012, 108, 185302. [Google Scholar] [CrossRef]
  28. Wessel, S.; Troyer, M. Supersolid Hard-Core Bosons on the Triangular Lattice. Phys. Rev. Lett. 2005, 95, 127205. [Google Scholar] [CrossRef] [Green Version]
  29. Kurdestany, J.M.; Pai, R.V.; Pandit, R. The Inhomogeneous Extended Bose-Hubbard Model: A Mean-Field Theory. Ann. Phys. 2012, 524, 234–244. [Google Scholar] [CrossRef] [Green Version]
  30. Zhang, X.-F.; Dillenschneider, R.; Yu, Y.; Eggert, S. Supersolid phase transitions for hard-core bosons on a triangular lattice. Phys. Rev. B 2011, 84, 174515. [Google Scholar] [CrossRef] [Green Version]
  31. Yamamoto, D.; Masaki, A.; Danshita, I. Quantum phases of hardcore bosons with long-range interactions on a square lattice. Phys. Rev. B 2012, 86, 054516. [Google Scholar] [CrossRef] [Green Version]
  32. Gheeraert, N.; Chester, S.; May, M.; Eggert, S.; Pelster, A. Mean-Field Theory for Extended Bose-Hubbard Model with Hard-Core Bosons. In Selforganization in Complex Systems: The Past, Present, and Future of Synergetics; Pelster, A., Wunner, G., Eds.; Springer: Zurich, Switzerland, 2016; pp. 289–296. [Google Scholar]
  33. Post, A.J.; Glandt, E.D. Statistical thermodynamics of particles adsorbed onto a spherical surface. I. Canonical ensemble. J. Chem. Phys. 1986, 85, 7349–7358. [Google Scholar] [CrossRef]
  34. Prestipino Giarritta, S.; Ferrario, M.; Giaquinta, P.V. Statistical geometry of hard particles on a sphere. Physica A 1992, 187, 456–474. [Google Scholar] [CrossRef]
  35. Prestipino Giarritta, S.; Ferrario, M.; Giaquinta, P.V. Statistical geometry of hard particles on a sphere: Analysis of defects at high density. Physica A 1993, 201, 649–665. [Google Scholar] [CrossRef]
  36. Prestipino, S.; Speranza, C.; Giaquinta, P.V. Density anomaly in a fluid of softly repulsive particles embedded in a spherical surface. Soft Matter 2012, 8, 11708–11713. [Google Scholar] [CrossRef]
  37. Vest, J.-P.; Tarjus, G.; Viot, P. Glassy dynamics of dense particle assemblies on a spherical substrate. J. Chem. Phys. 2018, 148, 164501. [Google Scholar] [CrossRef] [Green Version]
  38. Guerra, R.E.; Kelleher, C.P.; Hollingsworth, A.D.; Chaikin, P.M. Freezing on a sphere. Nature 2018, 554, 346–350. [Google Scholar] [CrossRef]
  39. Franzini, S.; Reatto, L.; Pini, D. Formation of cluster crystals in an ultra-soft potential model on a spherical surface. Soft Matter 2018, 14, 8724–8739. [Google Scholar] [CrossRef] [Green Version]
  40. Prestipino, S.; Giaquinta, P.V. Ground state of weakly repulsive soft-core bosons on a sphere. Phys. Rev. A 2019, 99, 063619. [Google Scholar] [CrossRef] [Green Version]
  41. Prestipino, S.; Sergi, A.; Bruno, E.; Giaquinta, P.V. A variational mean-field study of clusterization in a zero-temperature system of soft-core bosons. EPJ Web Conf. 2020, 230, 00008. [Google Scholar] [CrossRef]
  42. Zobay, O.; Garraway, B.M. Atom trapping and two-dimensional Bose-Einstein condensates in field-induced adiabatic potentials. Phys. Rev. A 2004, 69, 023605. [Google Scholar] [CrossRef] [Green Version]
  43. Garraway, B.M.; Perrin, H. Recent developments in trapping and manipulation of atoms with adiabatic potentials. J. Phys. B Mol. Opt. Phys. 2016, 49, 172001. [Google Scholar] [CrossRef]
  44. Elliott, E.R.; Krutzik, M.C.; Williams, J.R.; Thompson, R.J.; Aveline, D.C. NASA’s Cold Atom Lab (CAL): System development and ground test status. NPJ Microgravity 2018, 4, 16. [Google Scholar] [CrossRef]
  45. Lundblad, N.; Carollo, R.A.; Lannert, C.; Gold, M.J.; Jiang, X.; Paseltiner, D.; Sergay, N.; Aveline, D.C. Shell potentials for microgravity Bose-Einstein condensates. NPJ Microgravity 2019, 5, 30. [Google Scholar] [CrossRef] [Green Version]
  46. Prestipino, S.; Saija, F. Phase diagram of Gaussian-core nematics. J. Chem. Phys. 2007, 126, 194902. [Google Scholar] [CrossRef] [Green Version]
  47. Errington, J.R.; Debenedetti, P.G. Relationship between structural order and the anomalies of liquid water. Nature 2001, 409, 318–321. [Google Scholar] [CrossRef]
  48. Likos, C.N.; Lang, A.; Watzlawek, M.; Löwen, H. Criterion for determining clustering versus reentrant melting behavior for bounded interaction potentials. Phys. Rev. E 2001, 63, 031206. [Google Scholar] [CrossRef] [Green Version]
  49. Mladek, B.M.; Charbonneau, P.; Likos, C.N.; Frenkel, D.; Kahl, G. Multiple occupancy crystals formed by purely repulsive soft particles. J. Phys. Condens. Matter 2008, 20, 494245. [Google Scholar] [CrossRef]
  50. Prestipino, S.; Saija, F. Hexatic phase and cluster crystals of two-dimensional GEM4 spheres. J. Chem. Phys. 2014, 141, 184502. [Google Scholar] [CrossRef] [Green Version]
  51. Kunimi, M.; Kato, Y. Mean-field and stability analyses of two-dimensional flowing soft-core bosons modeling a supersolid. Phys. Rev. B 2012, 86, 060510(R). [Google Scholar] [CrossRef] [Green Version]
  52. Macrì, T.; Maucher, F.; Cinti, F.; Pohl, T. Elementary excitations of ultracold soft-core bosons across the superfluid-supersolid phase transition. Phys. Rev. A 2013, 87, 061602(R). [Google Scholar] [CrossRef] [Green Version]
  53. Prestipino, S.; Sergi, A.; Bruno, E. Freezing of soft-core bosons at zero temperature: A variational theory. Phys. Rev. B 2018, 98, 104104. [Google Scholar] [CrossRef] [Green Version]
Figure 1. The two regular meshes considered in this work (the circumscribed sphere is not shown). The dots (some are hidden) are the sites/nodes of the mesh; their number is 8 for the cubic mesh and 20 for the regular dodecahedral mesh. (Left) Cubic mesh with a regular tetrahedron inside (CT model). (Center) Regular dodecahedral mesh with a cube inside (DC model). (Right) Regular dodecahedral mesh with a regular tetrahedron inside (DT model).
Figure 1. The two regular meshes considered in this work (the circumscribed sphere is not shown). The dots (some are hidden) are the sites/nodes of the mesh; their number is 8 for the cubic mesh and 20 for the regular dodecahedral mesh. (Left) Cubic mesh with a regular tetrahedron inside (CT model). (Center) Regular dodecahedral mesh with a cube inside (DC model). (Right) Regular dodecahedral mesh with a regular tetrahedron inside (DT model).
Entropy 22 01289 g001
Figure 2. Thermal averages for the CT model (a) and the DC model (b), plotted as a function of μ for four values of T in units of V / k B , where k B is Boltzmann’s constant (see legend in the top left panel). (Top left): Grand potential. (Top right): Total number of occupied sites. (Bottom left): Total energy. (Bottom right): Order parameter (see text).
Figure 2. Thermal averages for the CT model (a) and the DC model (b), plotted as a function of μ for four values of T in units of V / k B , where k B is Boltzmann’s constant (see legend in the top left panel). (Top left): Grand potential. (Top right): Total number of occupied sites. (Bottom left): Total energy. (Bottom right): Order parameter (see text).
Entropy 22 01289 g002
Figure 3. CT model with multiple occupancy allowed. The grand potential Ω (in V units) is plotted as a function of μ for two families of “cluster crystals”, Ω n ( 4 ) (blue) and Ω n ( 8 ) (red), and for three values of U / V . Red and blue numbers are values of n.
Figure 3. CT model with multiple occupancy allowed. The grand potential Ω (in V units) is plotted as a function of μ for two families of “cluster crystals”, Ω n ( 4 ) (blue) and Ω n ( 8 ) (red), and for three values of U / V . Red and blue numbers are values of n.
Entropy 22 01289 g003
Figure 4. Phase diagram of the QCT model according to the decoupling approximation (see text). The blue dashed lines are second-order transition lines, while the red full line is a first-order line (only at the lobe vertex the transition is continuous).
Figure 4. Phase diagram of the QCT model according to the decoupling approximation (see text). The blue dashed lines are second-order transition lines, while the red full line is a first-order line (only at the lobe vertex the transition is continuous).
Entropy 22 01289 g004
Figure 5. Phase diagram of the QDC model according to the decoupling approximation (see text). The blue dashed lines are second-order transition lines, while red full lines are first-order lines (only at the lobe vertices the transition is continuous). Inside the blue circuit, the stable phase is supersolid (SS1 below μ = 3 V / 2 ; SS2 above μ = 3 V / 2 ). The full lines through the data points (open black squares) are a guide for the eye.
Figure 5. Phase diagram of the QDC model according to the decoupling approximation (see text). The blue dashed lines are second-order transition lines, while red full lines are first-order lines (only at the lobe vertices the transition is continuous). Inside the blue circuit, the stable phase is supersolid (SS1 below μ = 3 V / 2 ; SS2 above μ = 3 V / 2 ). The full lines through the data points (open black squares) are a guide for the eye.
Entropy 22 01289 g005
Figure 6. Exact average site occupancy in the QCT model (full symmetry occurs between the A and B sub-meshes). The plotted numbers are the occupancies in the whole regions delimited by the colored lines. These values reflect a perfect particle–hole symmetry around μ = 3 V / 2 . The MF boundary between the tetrahedral and superfluid “phases” is also shown for comparison (dashed line).
Figure 6. Exact average site occupancy in the QCT model (full symmetry occurs between the A and B sub-meshes). The plotted numbers are the occupancies in the whole regions delimited by the colored lines. These values reflect a perfect particle–hole symmetry around μ = 3 V / 2 . The MF boundary between the tetrahedral and superfluid “phases” is also shown for comparison (dashed line).
Entropy 22 01289 g006
Figure 7. Superfluid density in the QCT model: comparison between exact diagonalization (top panel) and decoupling approximation (bottom panel). The superfluid density is plotted as a function of μ for a number of t values in units of V, increasing from bottom to top (see legends).
Figure 7. Superfluid density in the QCT model: comparison between exact diagonalization (top panel) and decoupling approximation (bottom panel). The superfluid density is plotted as a function of μ for a number of t values in units of V, increasing from bottom to top (see legends).
Entropy 22 01289 g007
Figure 8. Grand potential of the QCT model as a function of μ , for a number of t values (from top to bottom, t = 0.05 , 0.1 , 0.15 , 0.2 , 0.3 , 0.4 , 0.5 ). MF data (thick full lines) are compared with exact values (full dots). To improve figure readability, the dots are joined by thin straight-line segments and the data are displaced vertically by 1 with respect to one another, starting from t = 0.1 .
Figure 8. Grand potential of the QCT model as a function of μ , for a number of t values (from top to bottom, t = 0.05 , 0.1 , 0.15 , 0.2 , 0.3 , 0.4 , 0.5 ). MF data (thick full lines) are compared with exact values (full dots). To improve figure readability, the dots are joined by thin straight-line segments and the data are displaced vertically by 1 with respect to one another, starting from t = 0.1 .
Entropy 22 01289 g008
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Prestipino, S. Ultracold Bosons on a Regular Spherical Mesh. Entropy 2020, 22, 1289. https://0-doi-org.brum.beds.ac.uk/10.3390/e22111289

AMA Style

Prestipino S. Ultracold Bosons on a Regular Spherical Mesh. Entropy. 2020; 22(11):1289. https://0-doi-org.brum.beds.ac.uk/10.3390/e22111289

Chicago/Turabian Style

Prestipino, Santi. 2020. "Ultracold Bosons on a Regular Spherical Mesh" Entropy 22, no. 11: 1289. https://0-doi-org.brum.beds.ac.uk/10.3390/e22111289

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop