Acoustic dispersion in a two-dimensional dipole system

Acoustic dispersion in a two-dimensional dipole system

Kenneth I. Golden

Department of Mathematics and Statistics, College of Engineering and Mathematical Sciences, University of Vermont, Burlington, Vermont 05401-1455, USA

Gabor J. Kalman

Department of Physics, Boston College, Chestnut Hill, Massachusetts 02467, USA Zoltan Donko and Peter Hartmann

Research Institute for Solid State Physics and Optics of the Hungarian Academy of Sciences, P.O. Box 49, H-1525 Budapest, Hungary 共Received 3 April 2008; published 8 July 2008兲

We calculate the full density response function and from it the long-wavelength acoustic dispersion for a two-dimensional system of strongly coupled point dipoles interacting through a 1/r³ potential at arbitrary degeneracy. Such a system has no random-phase-approximation共RPA兲limit and the calculation has to include correlations from the outset. We follow the quasilocalized charge共QLC兲approach, accompanied by molecular- dynamics 共MD兲 simulations. Similarly to what has been recently reported for the closely spaced classical electron-hole bilayer关G. J. Kalmanet al., Phys. Rev. Lett. 98, 236801共2007兲兴and in marked contrast to the RPA, we report a long-wavelength acoustic phase velocity that is wholly maintained by particle correlations and varies linearly with the dipole momentp. The oscillation frequency, calculated both in an extended QLC approximation and in the Singwi-Tosi-Land-Sjolander approximation关Phys. Rev. 176, 589共1968兲兴, is invari- ant in form over the entire classical to quantum domains all the way down to zero temperature. Based on our classical MD-generated pair distribution function data and on ground-state energy data generated by recent quantum Monte Carlo simulations on a bosonic dipole system关G. E. Astrakharchiket al., Phys. Rev. Lett. 98, 060405 共2007兲兴, there is a good agreement between the QLC approximation kinetic sound speeds and the standard thermodynamic sound speeds in both the classical and quantum domains.

DOI:10.1103/PhysRevB.78.045304 PACS number共s兲: 71.35.Ee, 52.27.Gr, 52.65.Yy, 05.30.Jp

I. INTRODUCTION

The formation of bound-electron-hole excitons in semi- conductors was predicted a long-time ago by Keldysh and co-workers and by Halperin and Rice.¹ Electron-hole bilay- ers共EHBs兲have created an especially promising medium for the formation of stable excitons.²In such systems the charges in the two layers have opposite polarities, and for sufficiently small layer separations, the positive and negative charges bind to each other in dipolelike excitonic formations. Recent quantum Monte Carlo共QMC兲 共Ref.3兲simulations and clas- sical Monte Carlo共MC兲 关Ref.4共a兲兴and molecular-dynamics 共MD兲 关Ref.4共b兲兴simulations have confirmed the emergence of the excitonic phase both in degenerate electron-hole³and in classical bipolar bilayers.⁴ The simulations have shown the existence of four phases in the strong-coupling regime:

Coulombic liquid and solid and dipole liquid and solid phases. The necessary existence of these phases was also pointed out in Ref.5. In electron-hole bilayers, the excitons may also form a Bose-Einstein condensate^3,6–8or possibly a supersolid.⁵

In a good approximation, the closely spaced EHB can be modeled as a two-dimensional共2D兲monolayer of interacting point dipoles, each of massm=me+mh. TheNpoint dipoles are free to move in the xy plane with dipolar moment ori- ented in the z direction; the interaction potential is accord- ingly given by ␾D共r兲=p²/r³, where p is the electric-dipole strength. The approximation that replaces the bound- electron-hole exciton by a point dipole has been considered by a number of investigators.^{9共a兲,9共b兲,10,11}

The coupling strength in the EHB is characterized at ar- bitrary degeneracy by ⌫˜=e²/a具Ekin典, where e is the elec- tronic charge anda= 1/

冑

_␲nis the average in-plane distance between particles. In the high-temperature classical domain, this becomes the customary coupling parameter ⌫=␤^e2/a, where ␤= 1/kBT, while at zero temperature, it becomes rs

=a/aB, wherea_B=ប/m_e,he². By the same token, for the di- pole system at arbitrary degeneracy,⌫˜

D=p²/a³具Ekin典can be taken as the coupling parameter, which becomes ⌫D

=␤p²/a³ in the high-temperature classical domain. At zero temperature, ⌫˜

D=rD=r₀/a is the appropriate measure of the coupling strength, where r₀=mp²/ប² is a characteristic length.^9共a兲 Here we focus on the strong-coupling regime⌫˜ Ⰷ1 that includes both the dipole liquid/solid phases. Since,D

in the symmetric 共me=mh=m/2兲 EHB, the Coulombic and dipole coupling parameters are related to each other by r_D

= 2共d²/a²兲rs共dis the spacing between layers兲, high coupling 共r_DⰇ1兲for point dipoles corresponds to the low-density re- gime in the closely spaced EHB, as dictated by the ordering a⬎dⰇaB.³

Various criteria have been put forward to determine whether the EHB can be considered as consisting of bound dipoles共excitons兲rather than individual electrons and holes.

De Paloet al.³defined the dipole phase on the basis of com- paring the energy of the assumed excitonic phase with that of a system of independent particles. In Ref. 4共a兲, the specific features of the correlation function showing a “correlation moat” surrounding the dipole were taken as the signature of the formation of permanent dipoles. Ranganathan and

Johnson^4共b兲 did not explicitly address the question of the Coulomb-dipole phase transition. In this paper, we will show that the 2D-point-dipole model as a representation of the system of dipoles of finite size is justified on the basis of the comparison of the values of the average potential in the point-dipole system with that in the EHB. All these criteria combined show that medium range d/a共sayd/a⬵0.6兲val- ues are sufficiently low to maintain the validity of the point- dipole model. Thus, the requirement for exceedingly highrs

values in the bilayer that seems to be required to balance closer layer separations共and which may be difficult to real- ize experimentally兲can be avoided.

A variety of collective modes can exist in the strongly coupled EHB system.^12,13 In the dipole approximation, the in-phase longitudinal mode can be identified with the density oscillation of the system of point dipoles. Prompted by this observation, there has been a recent flurry of activities di- rected at understanding the behavior of this collective mode both as a bilayer excitation and as a collective mode of the dipole system. A number of issues have emerged where re- sults obtained through different approaches are at variance with each other. First, there is the central question of the dependence of the collective-mode frequency on in-plane wave numberq. In Ref.5 an ␻共q→0兲⬀q^3/2dispersion has been proposed in the description of a Wigner 共super兲solid phase of dipoles. This behavior, in fact, can be regarded as the outcome of an application of the customary random- phase-approximation共RPA兲argument to a dipole system: the lack of validity of this approach is pointed out below. By contrast, all the investigators of Refs.9,10, and13reported an␻共q→0兲⬀q acoustic dispersion, albeit arrived at through different approaches. A second issue is the dependence of

␻共q→0兲on the layer separation distanced: the RPA analysis of the dipole system in Ref. 10 asserts that ␻共q→0兲⬀

冑

d, whereas the analysis in Ref. 13, which takes into account correlational effects from the outset, indicates that ␻共q

→0兲⬀d and so do the studies in Refs.9 of strongly corre- lated dipoles. Finally, there is the question of the precise numerical value of the phase velocity of the mode for which different values have been put forward in Refs.9,10, and13.

In this paper, we will rigorously show that the actual dis- persion for the density oscillations in the dipole liquid is acoustic, i.e., ␻共q→0兲⬃q, and that it can be attributed to dipole dynamics, i.e.,␻共q→0兲⬃d 共d being proportional to the dipole momentp兲; we will elucidate the correspondence between this result and a similar finding for the in-phase mode in the strongly correlated EHB.

One would expect that the RPA is an appropriate approach at least for the qualitative description of the collective modes both in the EHB and in the dipole systems. In fact, this expectation has not been borne out. The application of the RPA to an EHB leads to an acoustic q→0 behavior for the in-phase mode, but with an acoustic velocity s=␻0

冑

ad, where ␻02= 4e²/共ma³兲 is the nominal 2D plasma frequency.

This is not surprising, since the inappropriateness of the RPA for the analysis of the EHB has already been pointed out in Ref.13: the RPA treatment cannot reproduce the merging of the intrinsic dipole oscillation with the out-of-phase collec- tive mode. For the analysis of the dipole liquid, the inad- equacy of the RPA has a different origin and is more griev-

ous. An attempt to establish an RPA formalism fails entirely because the average Hartree field,

具␾D共r兲典H=n

冕

^{d2r␾D共r兲,}

of the dipole potential ␾D共r兲=p²/r³ diverges. Therefore, the Fourier transform of the dipole potential does not exist, im- plying that the 2D system of point dipoles interacting via this potential can have no RPA limit. As a consequence, the rou- tine argument that would generate the RPA collective-mode frequency via

␻²共q兲⬀q²␾D共q兲, ␾D共q兲= −q²␾COUL共q兲d²

关␾COUL共q兲= 2␲^e2/q兴 共Ref.5兲becomes invalid, ruling out the ensuing␻共q→0兲⬀q^3/2dispersion.

A correlational non-RPA study of the collective-mode spectrum of the EHBliquidin the classical domain was car- ried out in the recent work by Kalman et al.¹³ through a combined analytical/MD approach. The analytical portion of the study was based on the quasilocalized charge approxima- tion共QLCA兲,¹⁴which in thed→0 dipole limit led to a long- wavelength acoustic dispersion with phase velocity,

s共d→0兲=

冑

^2Km具␾D共r兲典=␻0d

冑

KI. 共1兲 This relationship was also corroborated by the accompanying MD simulations.¹³ With g₁₁共r¯兲, the inlayer pair distribution function, the average potential具␾D共r兲典is defined as

具␾D共r兲典=n

冕

^{d2rg11共r兲␾D共r兲= 2共p2/a3兲I,}

I=I共⌫,d兲=

冕

0

⬁

dr¯g₁₁共r¯兲1 r

¯², r

¯=r/a andK= 33/32 as calculated for the QLCA. For d/a

= 0.6, the integral I共⌫,d兲 is of the order of 0.8; its precise value and its dependence on ⌫ and d will be discussed in Sec. V. Expression共1兲is now in marked contrast to the RPA acoustic velocity s=␻0

冑

ad.

The analysis in Ref.13was extended to the solid phase of the EHB by applying the conventional harmonic approxima- tion for phonons and summing over lattice sites: the disper- sion of the longitudinal phonon in the d→0 dipole limit at T= 0 is given by a formula similar to Eq.共1兲, with a slightly different value of the integral; details will be given in Sec. III and Table I.

A similar calculation for the EHB lattice phonons was carried out by Kulakovskiiet al.¹²These authors also derived an acoustic dispersion but with a coefficient at variance with Eq. 共1兲:s共d→0兲=␲^1/4␻0d.

Turning now to the 2D-point-dipole system, Kachintsev and Ulloa¹⁰ were the first to analyze the collective excita- tions in a 2D fluid of bosonic dipoles modeled as point di- poles. They introduced a softened interaction potential

␾KU共r兲=␾D共r兲关1 − exp共−r²/d˜2兲兴, where˜d is of the order of the exciton radius; such a modification, which is Fourier transformable, makes it possible to realize an average Har-

tree field n兰d²r␾KU共r兲 and, consequently, an RPA limit re- sulting in an acoustic ␻共q→0兲⬀q dispersion with acoustic velocitys= 0.15␻0

冑

ad˜.

Convincing evidence for the acoustic behavior in the point-dipole system comes from recent QMC simulations of the degenerate bosonic dipole fluid in high^9共a兲 and rather weak to low^9共b兲 coupling regimes. These data provide indi- rect evidence for the acoustic dispersion through the appli- cation of the Feynman relation关see Eq.共45兲below兴with the input of computed static structure-function S共q兲 data. More will be said in Sec. V about quantitative comparisons of theory to the simulations in Ref. 9共a兲in the strong-coupling regime.

The QMC simulation of the degenerate bosonic dipole system was extended beyond the freezing point into the lat- tice phase.^9共a兲The results are not qualitatively different from those for the liquid phase and will be discussed in Sec. III.

After this lengthy preamble, we can state in precise terms the purpose of this paper: it is to approach the question of the small-qdispersion of the strongly coupled excitonic fluid by studying directly the collective-mode dispersion of the strongly coupled 2D dipolar liquid. We will show that the dispersion is acoustic, and we will study its characteristic sound velocity and its relationship to the corresponding quantity in the EHB. This analysis is the central objective of the present work.

We propose to calculate the collective-mode behavior by invoking two well-tested and rather different approaches:共i兲 the QLCA dynamical equation-of-motion/collective coordi- nate approach¹⁴ and 共ii兲 the Singwi-Tosi-Land-Sjolander 共STLS兲 kinetic equation approach.^15–17 There is no need to explicitly specify the degree of degeneracy in either formal- ism 共Secs. II and VI兲. However, the more involved STLS analysis共Sec. VI兲is carried out first in the high-temperature classical domain and then in the quantum domain at arbitrary temperatures. We contend that in the strong-coupling regime the QLCA is superior to the STLS approach; this is borne out by the comparison of the model theoretical results both with standard thermodynamic results and with those generated from computer simulations. Nevertheless, we follow this strategy in order to illustrate that in the domain of interest, quantum effects have no bearing on the architecture of the collective-mode dispersion. There remains one open ques- tion, namely, whether the formation of a condensate would affect the mode dispersion. That this indeed may be the case is known from the Bogoliubov analysis of the excitation spectrum of weakly interacting bosons. Here, however, the Bose-Einstein condensate fraction can be considered to be negligibly small since strong dipole-dipole interactions tend to destroy coherence. This observation is borne out by the QMC simulation in Ref.9共a兲.

The essential point to be noted in all the theoretical cal- culations is that, as discussed above, the dipole potential does not admit an RPA-like approximation. Hence, the den- sity response function for the 1/r³interaction can be calcu- lated only through the introduction of correlations in the for- malism from the outset. This is why it is crucial to rely on a nonperturbative calculational method, such as the proposed QLCA and STLS theoretical approaches.

Our theoretical analysis is accompanied by a molecular- dynamics study of the strongly coupled classical dipole liq- uid. While the full scope of this work will be reported elsewhere,¹⁸here we will cite the conclusions that pertain to the low-q behavior of the density oscillation mode.

As to the plan of the rest of the paper, in Sec. II, we reformulate the QLCA, which was originally created for the classical charged liquid,¹⁴ into an approximation method suitable for the description of the 2D system of strongly interacting point dipoles. This will be done through a two- stage procedure: 共i兲 first develop a classical QLCA theory along the lines in Ref.14to be followed by共ii兲its extension into the quantum domain along the lines in Ref. 19. In Sec.

III, we use the harmonic approximation to calculate the long- wavelength acoustic behavior of the longitudinal phonon in the 2D dipole crystal; we also compare the results with those of a similar calculation for the phonons in the EHB crystal.¹³ In Sec. IV, we calculate thermodynamic sound speeds of the 2D dipole liquid both in the classical and in the zero- temperature quantum domains. These may serve as standards for comparison in Sec. V where we establish linkages be- tween the QLCA sound velocity in Sec. II and the corre- sponding data generated from our classical MD and the quantum MC simulations in Ref.9共a兲; we also provide com- parisons with the results of our earlier work¹³ on the sound velocity in the EHB. In Sec. VI, in order to see how quantum effects may or may not alter the semiclassical results, we adapt the classical and quantum STLS kinetic equation de- scriptions of the three-dimensional共3D兲electron gas^15,16共us- ing the more tractable quantum kinetic equation formalism in Ref.20兲to the strongly coupled 2D-point-dipole liquid. Con- clusions are drawn in Sec. VII.

II. QLCA DESCRIPTION

We turn now to the formulation of a QLCA scheme for the model 2D monolayer of Nstrongly interacting point di- poles. LetAbe the large but bounded area of the monolayer and n=N/A the average density. In the two-stage develop- ment of the extended QLCA, we begin with the derivation in the classical domain.

The quasilocalized charge 共QLC兲 method has already been established and successfully applied to strongly coupled charged particle systems. Here we follow the paradigm of the original derivation, focusing on the differences that distin- guish the point-dipole system from a system of point charged particles. Similarly to what has been established for charged particle systems, the observation that serves as the basis of the QLC theory is that the dominating feature of the physical state of a classical dipolar liquid with coupling parameter

⌫DⰇ1 is the quasilocalization of the point dipoles. The en- suing model closely resembles a disordered solid where the dipoles occupy randomly located sites and undergo small- amplitude oscillations about them. However, the site posi- tions also change and a continuous rearrangement of the un- derlying quasiequilibrium configuration takes place. Inherent in the model is the assumption that the two-time scales are well separated and that it is sufficient to consider the time average 共converted into ensemble average兲 of the drifting quasiequilibrium configuration.

In the first stage, we wish to calculate the linear response to a weak perturbing external potential energy⌽Dext. Follow- ing the procedure in Ref.14, let Xi共t兲=xi+␰i共t兲 be the mo- mentary position of the ith point dipole, x_i its quasiequilib- rium site position, and ␰i共t兲 the perturbed amplitude of its small excursion; ␰i共␻兲 is its Fourier transform. In the equa- tions that follow,iandjsubscripts enumerate particles and␮ and␯ are vector indices; Einstein summation convention of the repeated vector indices is understood. The microscopic equation of motion for theith dipole is

−m␻²␰i,␮共␻兲+

兺

j

K_ij,␮␯␰j,␯共␻兲= − ⳵

⳵^xi,␮⌽D

ext共x_i,␮,␻兲, 共2兲

K_ij,␮␯=共1 −␦ij兲 ⳵²␾ij

⳵^x_i,␮⳵^x_j,␯⁻␦ij

兺

_ᐉ ^{共1 −␦iᐉ兲}_⳵xi,␮^⳵2␾_⳵^iᐉ_xᐉ,␯^,

共3兲 where␾ij=p²/兩x_i−x_j兩³is the 2D projection for fixed dipoles of the tensorial point-dipole potential,

␾ij,␳␴=p_i,␳p_j,␴

r_ij³

冋

^{␦␳␴− 3rij,r␳r}ij^ij,␴

2

册

^{, rij=兩xi−xj兩.}

Similarly, ⌽D

ext共xi,t兲 here represents a 共fictitious兲 external scalar potential 共energy兲 with an arbitrary space and time dependence derivatively coupled to the particle displacement

␰_i,␮共t兲. Equation共3兲shows the characteristic separation of the potential energy共1/2兲兺i,jK_ij,␮␯␰i,␮␰j,␯into diagonal共␦ij兲and off-diagonal 关共1 −␦ij兲兴 contributions: the former originates from the displacement of a dipole in a fixed environment of the other dipoles, while the latter originates from the fluctu- ating environment.

We next introduce collective coordinates␰k via the Fou- rier representation,

␰_i,␮共␻兲= 1

冑

mN

兺

_k ^{␰k,␮共␻兲exp共ik·xi兲. 共4兲}

Substituting Eq.共4兲into Eq.共2兲and following the procedure in Ref. 14, one ultimately obtains the ensemble-averaged equation of motion in terms of the dynamical tensorC_␮␯共q兲,

关␻²␦␮␯−C_␮␯共q兲兴␰q,␯共␻兲= inq_␮

冑

mN⌽D

ext共q,␻兲, 共5兲

C_␮␯共q兲= 1 mN

兺

i,j

具K_ij,␮␯exp关−iq·共xi−xj兲兴典

=3np²

m

冕

^d2rr15g共r兲关exp共iq·r兲− 1兴

冋

^{␦␮␯− 5r␮rr2␯}

册

^,

共6兲 whereg共r兲is the equilibrium pair distribution function. Pro- jecting out the longitudinal共LL: with respect to q兲element of the dynamical tensor, we derive an equation for the aver- age density responsen共q,␻兲by using the relationn共q,␻兲=

−共iqN/

冑

mN兲␰q,L共␻兲,

关␻^2−C共q兲兴n共q,␻兲=nq²

m ⌽Dext共q,␻兲, 共7兲

C共q兲 ⬅C_LL共q兲=3␲^np2 m

冕

0

⬁

dr1

r⁴g共r兲关3 − 3J₀共qr兲+ 5J₂共qr兲兴, 共8兲 whereJ₀共x兲andJ₂共x兲are Bessel functions of the first kind. It is useful to introduce the notation

⌿共q兲= m nq²C共q兲.

Equation共7兲and the constitutive relation

n共q,␻兲=␹共q,␻兲⌽Dext共q,␻兲 共9兲 then give the QLCA density response function,

␹共q,␻兲= nq²/m␻²

1 −⌿共q兲nq²/m␻^{2, 共10兲}

⌿共q兲=3␲^p2 q²

冕

0

⬁

dr1

r⁴g共r兲关3 − 3J0共qr兲+ 5J₂共qr兲兴. 共11兲 Note that the usual procedure of splitting⌿共q兲 into RPA and correlational parts by replacingg共r兲by 1 +h共r兲does not work because both of the separated terms would be repre- sented by divergent integrals.

We observe that the derivation of Eq.共10兲is predicated on the reasonable assumption that thermal motions are negli- gible in the high coupling regime. Nevertheless, the effects of random motion of the particles can be incorporated in the formalism²¹by replacing thenq²/m␻²factors in Eq.共10兲by the Vlasov density response function␹0

V共q,␻兲 共in the classi- cal domain兲 or by the Lindhard density response function

␹0

L共q,␻兲 共in the quantum domain兲.

In the quantum domain where the fluctuations are more important, the second-stage reformulation of the QLCA par- allels the procedure in Ref.19. One may accordingly assume that the effect of random quantum fluctuations is well repre- sented by replacing the nq²/m␻² factors by the Lindhard function,

␹0

L共q,␻兲= 1

A

兺

_k

冋

^{បn␻k−q−共ប/2−2/m兲knk+q/·2q}

册

^{, 共12兲}

resulting in the density response function,

␹共q,␻兲= ␹0 L共q,␻兲 1 −⌿共q兲␹0

L共q,␻兲, 共13兲 which replaces Eq. 共10兲 for arbitrary degeneracy; nk is the momentum distribution function for particles with energy spectrum ␧k=ប²k²/共2m兲, andN=兺knkwithn=N/A the av- erage 2D density. At strong coupling and in theq→0 limit,

␹0

L共q→0 ,␻兲⬇nq²/m␻², so that Eqs.共10兲and共13兲coincide and both the extended quantum QLC and classical QLC ap- proximations lead to

␹共q→0,␻兲 nq²/m␻²

1 −⌿共q→0兲nq²/m␻^{2, 共14兲}

⌿共q→0兲=33 8 ␲^p2

冕

0

⬁

dr1

r²g共r兲=33

16␲^a2具␾D共r兲典. 共15兲 Defining the dipole oscillation frequency 共equivalent of the 2D plasma frequency兲

␻D2= 2␲^p2n ma³ ,

the acoustic mode oscillation frequency then follows from setting the denominator of Eq.共14兲equal to zero,

␻²共q→0兲=33 8

␲^np2 m q²

冕

0

⬁

dr1

r²g共r兲= 33 16J共⌫˜

D兲␻D2a²q². 共16兲 We thus obtain the phase velocity as

s=␻Da

冑

2KJ共⌫˜

D兲, 共17兲

withK= 33/32 and J共⌫˜

D兲=

冕

0

⬁

dr¯1 r

¯²g共r¯兲, 共18兲 where¯r=r/a. We note thatJ共⌫˜

D兲is identical toI共⌫˜,d兲in Eq.

共1兲, except for the difference in the correlation functions un- der the integral:I共⌫˜,d兲is defined through the inlayer corre- lation function of the bilayer and has a weak dependence on d, while inJ共⌫˜

D兲the correlation function is that of the point dipoles. Tabulated values of the J共⌫˜

D兲 andI共⌫,d兲 integrals are displayed in TablesIIandIV, respectively.

To see thatg共r→0兲tends to zero sufficiently fast to guar- antee the convergence of the integral in Eqs. 共15兲,共16兲, and 共18兲, we observe that this is indeed the case in the high- temperature classical domain where one would expect that g共r→0兲⬀exp共−␤^p2/r³兲. In fact, this has been verified by our MD simulation. To make the case for convergence in the zero-temperature quantum domain, we observe that when two point dipoles are in close proximity to each other, the pair wave function ␺共r兲 and consequently the pair distribu- tion functiong共r兲⬀兩␺共r兲兩² are determined by the solution to the two-particle Schrödinger equation in ther→0 limit. Par- alleling Kimball’s electron-gas calculation,²² one readily finds that

␺共r→0兲=K₀共2

冑

r₀/r兲 ⬇

冑

_␲

2

冋

^rr0

册

^{1/4exp共− 2}

^冑

^r0/r兲,

共19兲 with the characteristic length r₀=mp²/ប² introduced above andK₀共2

冑

r₀/r兲is the modified Bessel function of the second kind. This small-r behavior has also been reported in Ref.

9共a兲. Consequently,

g共r→0兲⬀␲

4

冑

r^r₀ exp共− 4

冑

r₀/r兲, 共20兲 again guaranteeing the convergence of the integral in Eqs.

共15兲,共16兲, and共18兲.

III. DIPOLE SOLID

The philosophy of the QLCA scheme not being substan- tially different from that of the harmonic approximation for lattice phonons, the results of Sec. II can be converted,mu- tatis mutandis, into a description of the q→0 dispersion of lattice phonons.¹³ Based on Eq. 共5兲, the lattice dispersion relation can be written as

储␻²␦␮␯−C_␮␯共q兲储= 0, 共21兲

C_␮␯共q兲=3np² m

兺

_i

⫽0

关exp共iq·ri兲− 1兴1

r_i⁵

冋

^{␦␮␯− 5ri␮rri2i␯}

册

^.

共22兲 Since toO共q²兲the triangular lattice exhibits isotropic behav- ior, one can focus on the longitudinal dispersion and obtain toO共q²兲,

␻²共q→0兲=C共q→0兲=33

32M␻D2a²q² 共23兲 or

s_solid=␻Da

冑

³³32M, 共24兲

where

M=

兺

i

1

¯r_i³ 共25兲

is the lattice sum over the triangular lattice with¯r_i⬅ri/a. In effect,M/2 replaces the integralJ共⌫˜

D兲in Eq.共16兲. The value of M has been calculated by a number of workers^23–26 with slightly different results; the most recent semianalytic calcu- lation is due to Rozenbaum,²⁶according to which

M 2 =1

211.341

冋 ^冑

²␲³

册

^{3/2= 0.821. 共26兲}

Our own lattice sum computation for the 2D dipole crystal involving 1.9⫻10⁹particles provides

M

2 = 0.7985. 共27兲

From Eq.共24兲, the corresponding sound speed is then

s_solid= 1.283␻Da. 共28兲

This can be compared with the quantum MC formula for the potential energy of the dipole crystal quoted from Ref.9共a兲 as

E_triang= 4.446共nr₀²兲^3/2E₀= p²

2a³M, 共29兲 where E₀ is defined in Eq. 共31兲 below. From Eq. 共29兲 we calculate

M

2 = 4.446

␲^{3/2 = 0.7984, 共30兲} so that the sound speed共28兲again results.

As to the linkage with the classical zero-temperature EHB crystal, our results关Eqs.共27兲and共28兲兴can be compared with our corresponding EHB lattice sumsL共d兲, which, of course, depend on the layer separation,

L共d兲=

兺

_i

冋

^¯r1i−

^冑

^¯ri21+¯d2

册

^{. 共31兲}

TheL共d兲values, together with their associated sound veloci- ties, are tabulated in Table I. The tabulatedL共d兲 and sound speed values are quite close to the关Eqs.共27兲and共28兲兴M/2 and sound speed values, even ford/a= 0.6 as expected.

IV. DIPOLE LIQUID-PHASE THERMODYNAMICS We turn next to the straightforward derivation of the ther- modynamic sound speed in the 2D-point-dipole liquid. We consider first the classical domain. Starting from the correla- tion energy per particle of the dipole system,

E_corr=n

2

冕

^{d2r␾d共r兲g共r兲= pa32J共⌫D兲 共32兲}

or more succinctly

␤E_corr=⌫DJ共⌫D兲. 共33兲 The J共⌫D兲 integral is defined in Eq. 共18兲 above. The total thermodynamic pressure Pis calculated to be

␤P n = 1 +3

2⌫DJ共⌫D兲, 共34兲 and the dipole sound speed formula

␤ms²共⌫D兲=␤⳵^P

⳵^{n = 1 +} 15

4 ⌫DJ共⌫D兲+9 4⌫D

2J

⬘

共⌫D兲 共35兲 readily follows from Eq.共34兲. To facilitate comparison with the QLCA sound speed关Eq.共17兲兴, Eq.共35兲can also be recast ina␻Dunits,

s=␻Da

冑

_2⌫¹_D^{+ 15}8 J共⌫D兲+9

8⌫DJ

⬘

共⌫D兲. 共36兲 Isothermal sound speed values generated from Eq. 共36兲are tabulated in Table IIover a wide range of liquid-phase cou- pling strengths. We note that the QLCA values are a few percentabovethe corresponding thermodynamic values.

We turn next to the derivation of the thermodynamic sound speed in the 2D bosonic dipole liquid at zero tempera- ture. The starting point for the calculation is the ground-state energy fitting formula given by the quantum MC simulation work^9共a兲 in the strong-coupling regime of the degenerate di- pole liquid,

E=关a₁共nr₀²兲^3/2+a₂共nr₀²兲^5/4+a₃共nr₀²兲^1/2兴E₀

共4ⱕnr₀²ⱕ256兲, 共37兲

where E₀=ប²/mr₀² is the dipole equivalent of the Rydberg energy, r₀=mp²/ប²,

a₁= 4.536, a₂= 4.38, and

a₃= 1.2.

The thermodynamic pressure and sound speed are then cal- culated to be

TABLE I. EHB lattice sum and sound velocity共ina␻Dunits兲as functions of layer separationd.

d/a L共d兲 s_solid共a␻D兲

0.1 0.7974 1.282

0.2 0.7944 1.280

0.3 0.7895 1.276

0.4 0.7828 1.271

0.5 0.7744 1.264

0.6 0.7646 1.256

TABLE II. 2D-point-dipole liquid: QLCA共s_QLCA兲, MD共s_MD兲, and thermodynamic共s_COMP兲sound speeds as functions of the classical coupling parameter⌫D. Columns 3–5 is in units ofa␻D; columns 6–8 are in units of 1/冑␤m.

⌫D J共⌫D兲 s_QLCA共a␻D兲 s_MD共a␻D兲 s_comp共a␻D兲 s_MD共1/冑_␤m兲 s_comp共1/冑_␤m兲 s_QLCA共1/冑_␤m兲

20 0.8865 1.352 1.339 1.272 8.47 8.042 8.55

40 0.8515 1.325 1.268 1.251 11.34 11.19 11.85

60 0.8367 1.314 1.281 1.243 14.03 13.62 14.39

80 0.8298 1.308 1.284 1.239 16.24 15.67 16.55

100 0.8245 1.304 1.285 1.235 18.17 17.47 18.44

P=n²⳵^E

⳵^{n =n共nr0}

2兲

冋

^{32a1共nr02兲1/2+54a2共nr02兲1/4}

+1

2a₃共nr₀²兲^−1/2

册

^{E0, 共38兲}

ms²=⳵^P

⳵^{n =nr02}

冋

^{154 a1共nr02兲1/2+4516a2共nr02兲1/4}

+3

4a₃共nr₀²兲^−1/2

册

^{E0. 共39兲}

In terms of the convenienta␻Dunits, Eq.共35兲becomes s²= a²␻D2

2␲^3/2

冋

^{154 a1+4516a2共nr02兲−1/4+ 34nra3}02

册

共4ⱕnr₀²ⱕ256兲. 共40兲

Sound speeds generated from Eqs. 共39兲 and 共40兲 are tabu- lated in TableIIIover a wide range of liquid-phase coupling strengths. In contrast to the classical case, the QLCA values are a few percent below the corresponding thermodynamic values.

V. SIMULATIONS

In the strong-coupling regime, we expect that the collective-mode behavior is well emulated by a classical model. In order to further study the collective-mode behavior and to assess the validity of the QLCA, we have performed a classical MD simulation of the 2D dipole liquid. Details of this method and of the result for the full dispersion of the entire mode spectrum will be described elsewhere.¹⁸Here we quote the relevant q→0 results for the longitudinal collec- tive mode. The values of the integral J共⌫D兲 as calculated from MD simulated pair distribution functions and the QCLA as well as MD sound velocities are given in TableII.

Using the thermodynamic sound velocity关Eq.共36兲兴as a ref- erence, the discrepancy between QLCA and thermodynamic sound velocities ranges from 5.55% at⌫D= 100 to 6.34% at

⌫D= 20. As expected, the MD sound velocity is somewhat closer to the thermodynamic sound velocity.

As to correspondence with the classical EHB liquid, we see from TablesIIandIVthat the EHB and 2D-point dipole respective I共⌫,d兲 andJ共⌫D兲 values are quite close. In gen- eral, the EHB sound speeds are lower than their point-dipole

counterparts. For example, at⌫= 121 andd/a= 0.6, the EHB sound speed based in Eq.共1兲is⬃1.2% lower than that com- puted for the 2D-point-dipole liquid at the corresponding

⌫D= 43.6 coupling strength, but this is within the uncertain- ties of the determination of the sound speed from the MD data. We can also compare the point-dipole QLCA sound velocities 共TableII, column 3兲 with the “exact”关i.e., calcu- lated without the assumption of linear dependence on d in

␻共q,d兲兴EHB QLCA phase velocities 共TableIV, column 5兲.

Again the EHB sound speeds are lower, with a somewhat larger difference, ranging from 7.6% at ⌫D= 20 to 6.2% at

⌫D= 100. Thus, we may conclude that the 2D-point-dipole model reasonably well emulates the in-phase mode of the EHB.

At the present time, only the classical MD simulations generate a direct description of the dynamics of the collec- tive modes. The recent quantum MC simulations in Ref.9共a兲 of a strongly coupled bosonic dipole system at zero tempera- ture, however, provide some indirect insight into the collective-mode structure atT= 0. The comparison to the pre- viously obtained QLCA and MD results can be afforded on three different levels: first, the average potentialJ共⌫˜

D兲in Eq.

共18兲can be replaced by itsT= 0 equivalent; second, the fitted ground-state energy equation of state formula, as given in Ref.9共a兲, can be employed to find the thermodynamic sound speed; and third, in order to obtain the collective-mode dis- persion from static structure-function data, the Feynman con- struction, as used in Ref.9共a兲, can be invoked. The first step TABLE III. 2D-point-dipole liquid: QLCA共s_QLCA兲, quantum MC共s_QMC兲, and thermodynamic共s_COMP兲

sound speeds as a function of the zero-temperature coupling parameterr_D=r₀/a. Columns 3–5 are in units of a␻D. Columns 6 and 7 are in units of

冑

E₀/m, whereE₀=ប²/mr₀².

nr₀² r_D s_QMC共a␻D兲 s_QLCA共a␻D兲 s_comp共a␻D兲 s_QLCA共

冑

E₀/m兲 s_comp共

冑

E₀/m兲

32 10.0 ⬃1.76 1.3 1.14 58.2 63.3

64 14.2 ⬃1.70 1.3 1.39 97.9 105

128 20.0 ⬃1.69 1.3 1.36 164.7 172.7

256 28.4 ⬃1.80 1.3 1.34 277.7 286.2

TABLE IV. EHB: IntegralI and QLCA sound speeds共in a␻D

units兲as a function of⌫ford/a= 0.6. The QLCA sound speeds in column 4 are calculated from Eq.共1兲 with the input of column 3.

The more exact sound speed values in column 5 are extracted from the EHB in-phase oscillation frequency共2兲in Ref.13which is valid for arbitraryqanddvalues.

⌫ ⌫D I s_QLCA^a s_QLCA^b

52 18.72 0.8506 1.325 1.251

60 21.60 0.8449 1.320 1.247

105 37.80 0.8309 1.304 1.236

121 43.56 0.8285 1.307 1.234

212 76.32 0.8186 1.299 1.226

243 87.48 0.8140 1.296 1.224

280 100.8 0.8126 1.295 1.223

aEquation共1兲.

bReference13.

is made possible by recalling that 共n/2兲兰d²r␾D共r兲g共r兲

=共p²/a³兲J共rD兲is the dipole-dipole interaction energy per par- ticle, E_int共rD兲, where rD=r₀/a is the effective coupling pa- rameter共defined above兲for the zero-temperature 2D dipolar fluid. From the fitting formula in Ref. 9共a兲 for the ground- state energy, we identify the interaction energy as the leading 共nr₀²兲^3/2term in the series. Introducing the dipole equivalent of the Rydberg energy^9共a兲

E₀=p² r₀³ = ប⁶

p⁴m³= ប²

mr₀², 共41兲 one finds the relation

E_int=a₁共nr₀²兲^3/2E₀= p²

a³J共rD兲. 共42兲 Taking a₁= 4.536 from the fitting formula in Ref.9共a兲, one can then calculate

J共rD兲= a₁

␲^{3/2= 0.8146 共43兲}

as the共apparently independent ofr_D兲value of the integral, in reasonable agreement with its classical equivalent for ⌫D

ⱖ60. The QLCA sound velocity

s= 1.30␻Da 共44兲

then results from Eqs.共17兲and共42兲. This value can be com- pared with the other data in TableII; note in particular that it is quite close to thes= 1.314␻Da value at⌫D= 60.

Addressing the second comparison, the simulation in Ref.

9共a兲 reports that the 2D dipolar fluid crystallizes at nr₀²

= 290共corresponding to rD= 30.18兲. We therefore selectnr₀²

= 32, 54, 128, and 256共r_D= 10, 14.2, 20, and 28.4兲as repre- sentative of the strongly coupled fluid phase for our compari- son. These values, when plugged into Eq. 共40兲, result in the thermodynamic sound speeds tabulated in Table III. The thermodynamic sound speeds are higher than the QLCA ki- netic sound speed= 1.3a␻D. This discrepancy decreases from 7.8% to 3% as the coupling parameter increases from rD

= 10 to 28.4.

Addressing the third comparison, from Ref. 9共a兲we cal- culate the quantum MC sound speed based on Feynman’s relation,

ប␻共q→0兲= ប²q²

2mS共q→0兲, 共45兲

It appears that the phase velocity extracted from Eq. 共45兲 compares less satisfactorily with the thermodynamic sound speed关Eq.共40兲兴and with the corresponding QLCA value. In this connection, we refer to Fig. 3 in Ref. 9共a兲showing the dispersion curves generated from Eq. 共45兲 with input of QMC S共q兲 data 共note that the vertical axis in that figure is incorrectly labeled; it should read mEk/nប²兲. Choosing the wave-number value q= 0.5

冑

n in the acoustic domain and 共roughly兲 reading the mEk/nប² values off the nr₀²= 32, 64, 128, and 256 fluid-phase curves of Fig. 3 in Ref.9共a兲result in the quantum MC sound speeds tabulated in TableIII. We find that the thermodynamic sound speed is 22%–34% lower than the QMC sound speed. While it is true that the Feynman

excitation spectrum关Eq.共45兲兴constitutes an upper bound to the actual collective-mode dispersion, it is also the case that Eq. 共45兲 should reasonably well describe the dispersion in the acoustic regime, especially for zero-temperature bosons.

This would imply that the sizable discrepancies could be due to possible inaccuracies in the input S共q→0兲 data in Ref.

9共a兲. In any case, the resolution of this issue is not in the purview of the present work.

We can also compare the classical and quantum QLCA sound speeds as given in Tables II andIII, and we observe the marked closeness of the two.

VI. STLS DESCRIPTION

We have already stated the philosophy that leads us to pursue the STLS calculation, in addition to the QLC analysis already carried out, even though the comparison with various simulation results convincingly demonstrates the reliability of the QLCA scheme. What makes the STLS method attrac- tive is that it has both a classical and a quantum formulation and that the latter is derived from first principles, without recourse to the heuristic arguments exploited in the quantum generalization of the QLCA. Thus, what we are interested in here is less the actual value of the sound velocity, as pre- dicted by the STLS scheme, but rather seeing whether the calculation corroborates the conclusion we have arrived at through the QLCA analysis, namely, that in the q→0 limit there is no difference between the classical and quantum ar- chitectures of the point-dipole system’s longitudinal collec- tive mode.

First, we adapt the STLS kinetic equation approximation scheme¹⁵ to the calculation of the density response function and long-wavelength dispersion for the 2D dipolar fluid in the high-temperature classical domain. The starting point for the calculation is the Fourier-transformed linearized kinetic equation for the perturbed one-particle distribution function f^共1兲共v,r,␻兲, which, in the presence of a weak external dipole potential energy⌽D

ext共q,␻兲, is given by 关␻^−q·v兴f^共1兲共v,q,␻兲+ 1

mq·⳵^f0共1兲共v兲

⳵^{v ⌽Dext共q,}␻兲

= i m

⳵

⳵^v·

冕

^{d2rexp共−iq·r兲}

冕

^d2r

^⬘

⫻

冕

^d2v

^⬘

^{f共2兲共v,r;v}

^⬘

^,r

^⬘

^;␻兲ⵜ兩r−^p2r

⬘

兩³, 共46兲 wheref^共2兲共v,r;v

⬘

^,r

⬘

^;␻兲is the perturbed two-particle veloc- ity distribution function,f₀^共1兲共v兲=n共␤m/2␲兲exp共−␤^mv2/2兲is the Maxwellian distribution normalized to the average 2D density n, andn共q,␻兲=兰d²vf^共1兲共q,␻兲is the average density response. Introducing the equilibrium pair distribution func- tion g共兩r−r

⬘

兩兲, we then make use of the STLS closure hy- pothesis,

f^共2兲共v,r;v

⬘

^,r

⬘

^;␻兲=关f0共1兲共v兲f^共1兲共v

⬘

^,r

⬘

^,␻兲

+f^共1兲共v,r,␻兲f₀^共1兲共v

⬘

^{兲兴g共兩r−r}

⬘

^兩兲, 共47兲 which, when substituted into Eq. 共46兲, gives

关␻^−q·v兴f^共1兲共v,q,␻兲+ 1

mq·⳵^f0共1兲共v兲

⳵^{v ⌽Dext共q,}␻兲

= −3ip² m

⳵^f0共1兲共v兲

⳵^{v ·n共q,}␻兲

冕

^{d2RRR5g共R兲exp共−iq·R兲,}

共48兲 whereR=r−r

⬘

. Solving forf^共1兲共v,q,␻兲and taking the den- sity moment, one readily obtains

n共q,␻兲= ␹0 V共q,␻兲 1 −⌳共q兲␹0

V共q,␻兲⌽D

ext共q,␻兲, 共49兲

⌳共q兲=3p²

q²

冕

^{d2rqr·5rg共r兲exp共−iq·r兲}

=6␲^p2 q

冕

0

⬁

dr1

r³g共r兲J1共qr兲. 共50兲 Note that

⌳共q→0兲= 3␲^p2

冕

0

⬁

dr1

r²g共r兲=3

2␲^a2具␾D共r兲典. 共51兲 One can observe the natural emergence of the Vlasov func- tion,

␹0

V共q,␻兲= − 1

m

冕

^d2vq·␻^{⳵f−0共1兲q共v·兲/v⳵v 共52兲} in the formalism. Comparing Eq. 共49兲 and the constitutive relation共9兲, we obtain

␹共q,␻兲= ␹0 V共q,␻兲 1 −⌳共q兲␹0

V共q,␻兲. 共53兲 When compared with the QLCA, the replacement of ⌿共q兲 with⌳共q兲 is the hallmark of the STLS approach. Then fol- lowing the same pattern of reasoning, one finds, similarly to Eq. 共16兲,

␻²共q→0兲= 3␲^np2q2 m

冕

0

⬁

dr1

r²g共r兲=3

2J共⌫D兲␻D 2a²q²

共54兲 or

s=␻Da

冑

³2J共⌫d兲. 共55兲 We therefore recover the acoustic phase velocity 关Eq. 共17兲兴 with the QLCA K= 33/32 value therein replaced by the STLSK= 3/4 value.

We turn now to the STLS description of the collective mode in the quantum domain. The analysis is facilitated by adapting Niklasson’s quantum kinetic equation for the 3D electron fluid²⁰ to the 2D dipolar bosonic fluid at arbitrary temperature. Referring to the definitions of the distribution functions provided in Ref.20, the starting point for our cal- culation is the linearized kinetic equation for the perturbed one-particle Wigner distribution function 共WDF兲, f_k^共1兲共q,␻兲,

which, in the presence of the weak external dipole potential energy ⌽Dext共q,␻兲, is given by

关ប␻⁻共ប²/m兲k·q兴f_k^共1兲共q,␻兲− 1

A关n_k−q/2−n_k+q/2兴⌽Dext共q,␻兲

= 1

A

冕

^{d2r␾D共r兲}_q

^兺

_⬘,k⬘^exp共i−q

^⬘

^{·r兲关fk−q共2兲 ⬘/2,k⬘共q−q}

^⬘

^,q

^⬘

^;␻兲

−f_k+q^共2兲 _⬘/2,k⬘,共q−q

⬘

^,q

⬘

^;␻兲兴, 共56兲 where f_k−q^共2兲 _⬘/2,k⬘共q−q

⬘

^,q

⬘

^;␻兲 is the perturbed two-particle WDF,nkis the momentum distribution function for particles with energy spectrum␧k=ប²k²/共2m兲,n=共1/A兲兺kn_k=N/Ais the average 2D density, andn共q,␻兲=兺kf_k^共1兲共q,␻兲is the per- turbed density response to ⌽D

ext共q,␻兲. Introducing the pair distribution function g共r兲 with Fourier transform g共q兲, we now invoke the linearized STLS hypothesis,¹⁷

f_k^共2兲_{⫾q⬘/2,k⬘}共q−q

⬘

^,q

⬘

^;␻兲=1

A关nk⫾q⬘^/2f_k^共1兲_⬘共q,␻兲g共兩q−q

⬘

^兩兲 +n_k⬘f_k^共1兲_⫾q⬘/2共q,␻兲g共q

⬘

兲兴, 共57兲 which, when substituted into Eq. 共56兲, gives

f_k^共1兲共q,␻兲= 1

A

冋

^បn␻^k−q−共ប^/2−2/m兲k^nk+q/2·q

册

^{⌽Dext共q,␻兲}

+n共q,␻兲

冕

^{d2r␾D共r兲A12}

^兺

_q⬘

冋

^{បn␻k−q−⬘共ប/2−2/mnk+q兲k⬘·/2q}

册

⫻g共兩q−q

⬘

兩兲exp共−iq

⬘

^·r兲. 共58兲 Upon taking the density moment of Eq.共58兲and comparing the result with constitutive relation 共9兲, one readily obtains the quantum STLS density response function for the 2D di- polar liquid,

␹共q,␻兲=

␹0 L共q,␻兲

1 −

冕

^{d2r␾D共r兲A1}

^兺

_q⬘ ^g共兩q−q

^⬘

^{兩兲␹0L共q,q}

^⬘

^{,␻兲exp共−iq}

^⬘

^·r兲.

共59兲 Note the natural emergence of the inhomogeneous Lindhard function,^16,17

␹0

L共q,q

⬘

^,␻兲=1 A

兺

k

冋

^{បn␻k−q−⬘共ប/2−2/m兲knk+q⬘·/2q}

册

^{, 共60兲}

where␹0

L共q,␻兲=␹0

L共q,q,␻兲. At long wavelengths and in the strong-coupling regime,

␹0

L共q→0,q

⬘

^,␻兲=n共q·q

⬘

^兲/共m␻²兲. 共61兲 Consequently,