text:
# A dilute mixture of atoms and molecules

## Abstract

We study a dilute gas with two species of Fermionic atoms of unequal concentrations, interacting via a short-range interaction with one deeply bound state. We study the properties of this system under the mean-field approximation. We obtain the effective interaction among the fermions and bosons, and discuss the collective modes of the system.

PACS numbers: 03.75.Fi, 67.60.-g

Recent experimental advances in trapping and cooling alkali atoms have already given us condensates of Bosons  and degenerate gases of Fermions . Much theoretical and experimental efforts have now turned to systems with trapped molecules, with interesting observations and predictions on the equilibrium and dynamical properties of this system.  These molecules can have their constituent atoms be either Bosons or Fermions. The latter case, on which the rest of this paper will concentrate, is also of particularly interest when it is considered together with superfluidity arising from Cooper pairing of Fermions. Indeed superfluid 3He and most known superconductors and have their superfluidity derived from Cooper pairs. It is also widely believed that Fermionic pairing in diatomic molecules and in Cooper pairs are in fact two limits of a continuum.  Much work has already been devoted to the cross-over between these two limits.  (and references therein) There is also a lot of interest in producing pairs (within any part of the above mentioned continuum) in trapped atomic Fermions.

We here would like to further our understanding of systems consisting of ”molecules”. We shall limit ourselves to ”molecules” made up of two distinguishable Fermions, interacting via a short range ($`<<a`$ below) interaction, in the s-wave state (Therefore by ’molecule’ we simply mean the Boson made up by a pair of Fermions). These two partners of the molecule can in general be atoms in different hyperfine spin projections or internal states, or different atomic species. For definiteness we shall restrict ourselves to the first case and for simplicity refer them hereafter as spin ”up” and ”down”, with concentrations $`n_{}`$ and $`n_{}`$. We shall study the properties of this system when the two constituents of the molecule are of different concentrations, say when $`n_{}>n_{}`$. Furthermore, we restrict ourselves to the case where the size of the bound state $`a`$ is small, such that $`n_{}a^3`$ and $`n_{}a^3`$ are both much smaller than $`1`$. Then at temperatures much below the dissociation temperature $`T_{\mathrm{disso}}`$ all down spins are paired into molecules but there are left over (or excess) spin-up Fermions. The Bosonic molecules and Fermionic atoms basically do not overlap and are well-defined entities. Thus one has a Bose-Fermi mixture where the Fermions are identical with one of the constituent particles of the Bosons, with the Boson and Fermion concentrations $`n_b=n_{}`$ and $`\nu =n_{}n_{}`$ respectively.

We shall obtain the properties of our Bose-Fermi mixture by applying many-body techniques on the (more basic) constituent spin-up and down Fermions (generalizing ref  —– we shall also be following the notations of the last reference quite closely, except that we shall use the more common sign conventions for the single particle propagators and self-energies). This system is attractive theoretically in that all interactions can be characterized by one interaction parameter, namely the s-wave scattering length $`a`$ ($`>0`$) between the spin up and down Fermions, or equivalently the binding energy $`ϵ_b`$ of the molecule ($`ϵ_b=\mathrm{}^2/ma^2`$, here $`m`$ is the mass of either of the Fermions). There are, however, several energy scales and hence temperature regimes. As mentioned we shall limit ourselves to where the temperature $`T`$ of the system is much lower than the dissociation temperature $`T_{\mathrm{disso}}`$ ($`ϵ_b`$ up to some entropic factors). The relevant scales are thus $`T`$, the Bose-Einstein condensation temperature $`T_{BEC}`$ for the Bosons and the Fermi temperature $`T_F`$ for the (left-over) Fermions ( $`T_{BEC}n_b^{2/3}\mathrm{}^2/m`$ and $`T_F\nu ^{2/3}\mathrm{}^2/m`$. Note they are much less than $`T_{\mathrm{disso}}`$). We shall study the some properties of this system at $`T=0`$ and compare them with those of a general dilute Bose-Fermi mixture. In particular, we shall be able to find the effective interaction, or equivalently the scattering length, between the Bosons and Fermions in terms of $`a`$. We further show that the Bogoliubov mode of the system is the same as that of a Bose-Fermi mixture.

In the case where the pairing between the spin-up and down particles are in the Cooper limit (weak attractive interaction so that $`a<0`$ and small), the system with $`n_{}>n_{}`$ has been studied before in the context of superconductivity under an applied Zeeman field $`h`$ . We shall compare our results with those of this better known system.

We now consider the system at $`T=0`$ under the mean-field approximation The necessary mean field equations can be obtained by generalizing those of ref  (which deals with $`n_{}=n_{}`$). We introduce three fields, the chemical potentials $`\mu _{}`$, $`\mu _{}`$ for the two spin species and a pairing field $`\mathrm{\Delta }`$. These fields obey the self-consistent equations

$$1=\frac{4\pi \mathrm{}^2a}{m}_\stackrel{}{k}\left[\frac{1}{2ϵ_k}\frac{1}{2E_k}\left(1f(E_kh)f(E_k+h)\right)\right]$$
(1)
$$n_\sigma =_\stackrel{}{k}n_\sigma (\stackrel{}{k}),\sigma =\mathrm{or}$$
(2)
with

$$n_{}(\stackrel{}{k})=u_k^2f(E_kh)+v_k^2(1f(E_k+h))$$
(3)
$$n_{}(\stackrel{}{k})=u_k^2f(E_k+h)+v_k^2(1f(E_kh))$$
(4)
where $`f`$ is the Fermi function. We have defined for convenience the chemical potential average and difference $`\mu (\mu _{}+\mu _{})/2`$, $`h(\mu _{}\mu _{})/2`$ ($`>0`$), the energies $`ϵ_k\mathrm{}^2k^2/2m`$, $`E_k[(ϵ_k\mu )^2+|\mathrm{\Delta }|^2]^{1/2}`$, and the ”coherence factors” $`u_k^2=\frac{1}{2}(1+\frac{ϵ_k\mu }{E_k})`$, $`v_k^2=\frac{1}{2}(1\frac{ϵ_k\mu }{E_k})`$, $`u_kv_k=\frac{\mathrm{\Delta }}{2E_k}`$ are exactly as in BCS theory for superconductivity. Eq (1) is the self-consistent equation for $`\mathrm{\Delta }`$ whereas eqs (2) are the constraints for the Fermion concentrations, with $`n_\sigma (\stackrel{}{k})`$ the occupation number at wavevector $`\stackrel{}{k}`$. Notice that if we view the system as a Bose-Fermi mixture, the chemical potential for the Bosons is $`\mu _b=2\mu `$ and that for the (left-over) Fermions is $`\mu _{}`$. In the dilute limit $`0<x_{}n_{}a^3<<1`$ and $`0<x_{}n_{}a^3<<1`$ which we are studying (in contrast to the Cooper limit where $`a<0`$ ), it is convenient to study the solution to eqs (1)-(2) as an expansion in $`x_{}`$ and $`x_{}`$. At $`T=0`$, $`f(E_k+h)=0`$ and $`f(E_kh)`$ is non-vanishing only for $`k<k_c`$ where $`E_{kc}=h`$. From eq (3) and (4) we get immediately

$$\nu =k_c^3/6\pi ^2.$$
(5)
$`k_c`$ plays the role of the Fermi wavevector for the excess spin-up Fermions (which form a ‘Fermi sea’, see also below)

We shall see (as in the corresponding case of ref ) that $`\mu <0`$ and $`|\mathrm{\Delta }|<<|\mu |`$. Using this inequality, eq(1) and eq(4), we obtain two equations determinating $`\mathrm{\Delta }`$ and $`\mu `$. In the dilute limit we get

$$\mu =\frac{\mathrm{}^2}{2ma^2}+\frac{2\pi (n_b+2\nu )\mathrm{}^2a}{m}$$
(6)
$$|\mathrm{\Delta }|^2=8\pi n_bϵ_b^2a^3[1\mathrm{}]$$
(7)
(and hence $`|\mathrm{\Delta }|^2/\mu ^2n_ba^3<<1`$ as promised). Eq (6) can be interpreted easily by recalling $`\mu _b=2\mu `$. The first term arises from the binding energy of the molecule. The two parts of the second term represent the mean-field corrections to the chemical potential due to the interaction of the molecules among themselves and with the Fermions respectively. Equating these two expressions with ($`1/2`$ times) $`n_bg_{bb}`$ and $`\nu g_{bf}`$, we obtain the interaction constants

$`g_{bb}`$ $`=`$ $`{\displaystyle \frac{4\pi \mathrm{}^2a}{m}}`$ (8)
$`g_{bf}`$ $`=`$ $`{\displaystyle \frac{8\pi \mathrm{}^2a}{m}}`$ (9)
or equivalently the scattering lengths $`a_{bb}=2a`$ and $`a_{bf}=8a/3`$.

We now discuss some properties associated with the Fermions. The occupation numbers of the constituent Fermions $`n_\sigma (\stackrel{}{k})`$ are sketched in Fig. 1 for $`T=0`$. $`n_{}(k)=1`$ whereas $`n_{}(k)=0`$ for $`k<k_c`$. $`n_{}(k)`$ and $`n_{}(k)`$ are equal but small for $`k>k_c`$. These results can be understood by considering the wavefunction of the sytem .

It is also interesting to study the spectral function $`A_{}(\stackrel{}{k},\omega )`$ ($`=\frac{1}{\pi }\mathrm{Im}G_{11}^R(\stackrel{}{k},\omega )`$ given below) corresponding to the density of states for addition/removal of a spin-up Fermion at momentum $`\stackrel{}{k}`$ and energy $`\omega `$. $`A_{}(\stackrel{}{k},\omega )`$ has the same formal expression as in the case of Cooper pairing (c.f. eq (12) below): $`A(\stackrel{}{k},\omega )=u_k^2\delta (\omega (E_kh))+v_k^2\delta (\omega +(E_k+h))`$. There are thus two delta-function peaks, one at positive and one at negative energies. Since $`|\mathrm{\Delta }|<<\mu `$, $`v_k^2<<1`$ and thus the weight at the negative frequency is very small. To a very good approximation, $`A_{}\delta (\omega (ϵ_k\mu _{}))`$, that of a free Fermi gas. We shall compare these results with those of case (B) below.

We next study the dynamic properties of our system (at $`T=0`$) by evaluating the two particle vertex function (matrix) $`𝚪(Q,\mathrm{\Omega })`$ for repeated scattering between a spin-up and spin-down particle (hole) at small wavevector $`Q`$ and frequency $`\mathrm{\Omega }`$. We shall see that, as in the case of $`n_{}=n_{}`$,  $`𝚪`$ is related to the propagator for the Bosons. $`𝚪`$, in Matsubara frequencies (we shall drop the superscript $`M`$ for simplicity. The desired $`𝚪`$ in real frequency can be obtained via $`i\mathrm{\Omega }_\nu \mathrm{\Omega }`$ with $`\mathrm{\Omega }_\nu >0`$ as usual) obeys the Bethe-Salpeter equation (generalizing that of )

$$\mathrm{\Gamma }_{\alpha 1,\alpha 2}^1(Q,\mathrm{\Omega }_\nu )=\frac{m}{4\pi \mathrm{}^2a}\delta _{\alpha 1,\alpha 2}+M_{\alpha 1,\alpha 2}(Q,\mathrm{\Omega }_\nu )$$
(10)
where

$$M_{\alpha 1,\alpha 2}(Q,\mathrm{\Omega }_\nu )_\stackrel{}{k}\left[T\underset{ϵ_n}{}G_{\alpha 1,\alpha 2}(\stackrel{}{Q}\stackrel{}{k},\mathrm{\Omega }_\nu ϵ_n)G_{\alpha 1,\alpha 2}(\stackrel{}{k},ϵ_n)\frac{m}{\mathrm{}^2k^2}\delta _{\alpha 1,\alpha 2}\right]$$
(11)
and $`G_{\sigma \alpha 1,\alpha 2}`$ are the matrix Green’s function in particle ($`\alpha =1`$) and hole ($`\alpha =2`$) space for spin-$`\sigma `$ defined in a similar manner as in Gorkov’s theory of superconductivity , and $`\mathrm{\Omega }_\nu `$ and $`ϵ_n`$ are Bosonic and Fermionic Matsubara frequencies respectively. In principle there is a contribution to the self-energy for $`𝐆`$ proportional to $`𝚪`$ and therefore $`𝐆`$ and $`𝚪`$ have to be found self-consistently (see  ). However, at $`T=0`$ and under our diluteness condition we can simply use the mean-field solution for $`G`$, i.e.,

$$G_{11}(\stackrel{}{k},ϵ_n)=\frac{u_k^2}{iϵ_n(E_kh)}+\frac{v_k^2}{iϵ_n+(E_k+h)}$$
(12)
and

$$G_{12}(\stackrel{}{k},ϵ_n)=u_kv_k\left[\frac{1}{iϵ_n(E_kh)}\frac{1}{iϵ_n+(E_k+h)}\right]$$
(13)
where $`u_k`$, $`v_k`$ are the coherence factors introduced before. $`G_{}`$’s are given by similar expressions except $`hh`$ and $`\mathrm{\Delta }\mathrm{\Delta }`$. $`𝚪^1(Q=0,\mathrm{\Omega }_\nu =0)`$ can be calculated easily and is of order $`\mathrm{\Delta }^2`$ due to eq (1). However, evaluation of $`𝚪^1`$ at finite $`Q`$ and $`\mathrm{\Omega }`$ needs special care. In $`M_{11}`$ there arises the contribution

$`M_{11}^a(Q,\mathrm{\Omega })`$ $``$ $`{\displaystyle _\stackrel{}{k}}u_{\stackrel{}{Q}\stackrel{}{k}}^2v_\stackrel{}{k}^2T{\displaystyle \underset{ϵ_n}{}}\left[{\displaystyle \frac{1}{iϵ_n+(E_{\stackrel{}{Q}\stackrel{}{k}}h)\mathrm{\Omega }}}\times {\displaystyle \frac{1}{iϵ_n+(E_\stackrel{}{k}h)}}\right]`$ (14)
$`=`$ $`{\displaystyle _\stackrel{}{k}}u_{\stackrel{}{Q}\stackrel{}{k}}^2v_\stackrel{}{k}^2{\displaystyle \frac{[1sgn((E_{\stackrel{}{Q}\stackrel{}{k}}h)\mathrm{\Omega })sgn(E_\stackrel{}{k}h)]/2}{|E_{\stackrel{}{Q}\stackrel{}{k}}E_\stackrel{}{k}\mathrm{\Omega }|}}.`$ (15)
The integrand is singular as $`Q`$ and $`\mathrm{\Omega }`$ approaches $`0`$. To evalute $`M_{11}^a`$, we first note that under our diluteness approximation, we can replace $`E_k`$ by $`ϵ_k\mu `$ and thus $`E_kh=E_kE_{kc}ϵ_kϵ_{kc}`$ etc. The integrand is therefore finite only if $`ϵ_{\stackrel{}{Q}\stackrel{}{k}_c}ϵ_{k_c}\mathrm{\Omega }`$ and $`ϵ_kϵ_{k_c}`$ are of opposite signs, and thus for $`Q<<k_c`$, $`\mathrm{\Omega }<<k_c^2/m`$, $`\stackrel{}{k}`$ is restricted to be near the ”Fermi surface”, i.e., $`|\stackrel{}{k}|k_c`$. We can then use $`u_{\stackrel{}{Q}\stackrel{}{k}}^2v_\stackrel{}{k}^2u_{kc}^2v_{kc}^2\left(\frac{|\mathrm{\Delta }|}{ϵ_b}\right)^2`$, and replace the $`\stackrel{}{k}`$ integral by integration over $`\widehat{k}`$ and $`\xi _kϵ_kϵ_{k_c}`$. We finally obtain $`M_{11}^a(Q,\mathrm{\Omega })\frac{mk_c}{2\pi ^2}(\frac{|\mathrm{\Delta }|}{ϵ_b})^2`$, and hence

$$\mathrm{\Gamma }_{11}^1(Q,\mathrm{\Omega })=\frac{1}{8\pi ϵ_b^2a^3}\left\{\mathrm{\Omega }+\frac{Q^2}{4m}+\frac{|\mathrm{\Delta }|^2}{2ϵ_b}8\pi |\mathrm{\Delta }|^2a^3\left(\frac{mk_c}{2\pi ^2}\right)\right\}.$$
(16)
Similarly

$$\mathrm{\Gamma }_{12}^1(Q,\mathrm{\Omega })=\frac{1}{8\pi ϵ_b^2a^3}\left\{\frac{\mathrm{\Delta }^2}{2ϵ_b}8\pi \mathrm{\Delta }^2a^3\left(\frac{mk_c}{2\pi ^2}\right)\right\}.$$
(17)
As in the case of equal number of spin-up and down particles, $`𝚪`$ is thus found to be proportional to the connected part of the Bosonic propagator matrix $`𝐆_b`$. $`𝐆_b^1`$, under the random phase approximation, reads (c.f. )

$$G_{b11}^1(Q,\mathrm{\Omega })=\mathrm{\Omega }+\frac{Q^2}{2m_b}+g_{bb}n_bg_{bf}^2n_b\chi _f(Q,\mathrm{\Omega })$$
(18)
$$G_{b12}^1(Q,\mathrm{\Omega })=+g_{bb}n_bg_{bf}^2n_b\chi _f(Q,\mathrm{\Omega })$$
(19)
where $`\chi _f`$ is the density response of the Fermions and we have taken the gauge where the phase of the condensate is zero. To further check our claim that $`𝚪𝐆_b`$ above, we note that $`m_b=2m`$, $`g_{bb}n_b=\frac{|\mathrm{\Delta }|^2}{2ϵ_b}`$, $`g_{bf}^2n_b=8\pi |\mathrm{\Delta }|^2a^3`$ (using eqs (7), (8) and (9)), $`\chi _f(Q,\mathrm{\Omega })=\left(\frac{mk_c}{2\pi ^2}\right)`$. In this last expression, we have used the value of $`\chi _f(Q,\mathrm{\Omega })`$ for $`\mathrm{\Omega }/Q0`$. That this is the correct limit can be seen as follows. The bare sound velocity for the Bosons is given by $`c_b^0=(n_bg_{bb}/m_b)^{1/2}`$ whereas the bare Fermi velocity is $`v_f^0=k_c/m`$. Thus $`(c_b^0/v_f^0)^2(n_ba/k_c^2)<<1`$ under our diluteness condition. For the mode under consideration, $`(\mathrm{\Omega }/v_fQ)(c_b^0/v_f^0)<<1`$.

In conclusion, we have studied, via many-body technique, some properties of an atom-molecule mixture where the Fermionic atoms are identical with one of the constituent atoms of the diatomic molecules. If $`T<<T_{BEC}`$, our system behaves like a general Bose-Fermi mixture. In this case we have also obtained the scattering lengths $`a_{bb}`$ and $`a_{bf}`$ for our mixture in terms of $`a`$.

I would like to thank Darwin Chang for pointing out and making available to me ref . This research was supported by the National Science Council of Taiwan under grant number 89-2112-M-001-061.


format: text

text:
# Scaling laws for the 2d 8-state Potts model with Fixed Boundary Conditions

## I Introduction

The introduction of computer simulation methods has been a breakthrough in the study of phase transitions in lattice models. The rapid increase on the computers power has enabled to analyze with great accuracy the scaling laws, governed by the critical exponents, and even the corrections to these scaling laws. In such analysis, finite size effects must be taken into account carefully Binder . One of these effects, the disturbance from the boundary, has usually been dismissed by the adoption of periodic boundary conditions (PBC). Nevertheless, in some situations the adoption of periodic lattices may not be adequate, either for practical or theoretical reasons. This is the case of the free boundary conditions used in the analysis of free surfaces, or the so called boundary fields, used in the analysis of wetting phenomena. In the present paper we will focus on a particular election of the boundary conditions, the so called fixed boundary conditions (FBC), which have been recently applied to spin models BEJ ; BEJ2  and gauge models BF .

Second order phase transitions exhibit universality. For this reason, all the details of the system near the phase transition point become irrelevant for the critical exponents. By contrast, first order phase transitions are not universal and, hence, all details of a simulation must be considered carefully. This includes, in particular, the choice of boundary conditions. What are the appropriate set of scaling laws for a first order phase transition with FBC? This is the question we address in this paper. Starting from the theoretical analysis of Borgs and Kotecky BK ; BK2  and the diploma Thesis of Medved M  on the dynamics of first order transitions, we present the finite size scaling laws applicable to the case of FBC and we check them performing a numerical simulation of the 2-d 8-state spin Potts model Potts  with FBC.

The paper is divided as follows. In section II a brief summary of some recent simulations where FBC have been adopted serves as a motivation for a detailed analysis of the finite size scaling laws which are also presented and discussed. Section III is devoted to a discussion of our numerical simulation and the results that we have obtained. In section IV we analyze our results in the light of the scaling laws presented in Sec. II. Finally, in section V we give some concluding remarks.

## II Fixed Boundary Conditions

### II.1 The motivation for FBC

Recently, the so called gonihedric spin modelshave been proposed as a laboratory to study discrete versions of string theories SW . All these simulations have been performed imposing standard periodic boundary conditions on a three-dimensional lattice. Nevertheless, in some cases three intersecting inner planes of spins were fixed to break the large energy degeneracy of the hamiltonianBEJ ; BEJ2 . Due to the periodicity of the boundaries, this is equivalent to fix the spins belonging to the six planes of the 2-d boundary of the 3-d cube formed by the spins. Since for a certain range of the coupling parameter, in particular for $`\kappa =0`$, the transition is clearly of first order, the analysis of the finite size effects should have been done using the scaling laws presented in this paper. We expect that the application of the FBC scaling laws may overcome some anomalies recently observed in the analysis of this transition BV .

Another situation where the knowledge of the FBC scaling laws seems to be crucial is the issue of the triviality of lattice QED. Indeed, it has been claimed that the formation of artificial monopole structures, which close over the boundaries, in a simulation of the 4-d U(1) gauge model may be responsible for turning the phase transition of this model from second to first order LJ . To avoid this problem, originated probably by an incorrect choice of the boundaries, it was suggested to perform the Monte Carlo simulations on a lattice with a topology of an sphere. Along these lines, Baig and FortBF  proposed the adoption of FBC to simulate an spherical topology. Effectively, to fix all the variables belonging to the 3-d border to unity is the higher dimensional equivalent of converting a 2-d plane square lattice to the 2-d surface of a sphere by collapsing the lines of the border to a single point. Nevertheless, to discriminate between a first or a second order nature for a transition, an accurate analysis of data produced is necessary, and, in particular, this will be only possible if one knows for certain the applicable scaling laws.

### II.2 The scaling laws

Although speaking properly no critical exponents can be defined for first order phase transitions, it is usual to define a set of characteristic exponents, together with a set of scaling laws borrowed from those of the second order phase transitions. The pioneering work of Privman and Fisher Privman , Binder and Landau Binder2  and Challa et al. Challa , provided a phenomenological understanding of the scaling for first order transitions. A more rigorous theoretical justification for these first order scaling laws was presented by Borgs and Kotecky BK ; BK1 . The formulation of its applicability to finite size scaling expressions in terms of the lattice size was the work of Borgs, Kotecky and Miracle-Solé BKS  and, independently, of Janke Janke . But in all these developments the existence of periodic boundary conditions was assumed. Recently, though, Borgs and Kotecky BK2  have extended their analysis to include surface effects in addition to the standard volume effects which govern first order transitions. Following this work, Mendev M  has deduced the scaling laws for the spin Potts model in the presence of surface effects, in particular adopting boundary conditions other than the periodic ones.

Following the general analysis of Mendev M , finite size scaling laws in terms of the lattice size for the case of fixed boundary conditions can easily be deduced. They are summarized in Table 1, together with the standard laws for periodic conditions. In the rest of this paper we will check these modified scaling laws with the results of our numerical simulation.

It should be noticed that the suggestion that in the case of free boundary conditions every transition is shifted by a $`1/L`$ correction term caused by surface effects is quite old. Binder Binder3 , for instance, reports on a series of experimental results Experiments  supporting this conclusion.

## III Numerical simulation

To test the scaling laws of Table 1, we have performed a numerical simulation of the 2-d 8-state spin Potts model defined by the partition function

$$Z_{potts}=\underset{\{\sigma _i\}}{}e^{\beta E},$$
(1)
where the energy is

$$E=\underset{ij}{}\delta _{\sigma _i\sigma _j},(\sigma _i=1,\mathrm{},8),$$
(2)
with $`\beta =J/kT`$ in natural units. It is well known that this model exhibits a first order phase transition Wu  and for this reason it has been chosen as a test model in several previous studies.

Fixed boundary conditions have been implemented along the lines stated by Baig and Fort BF . In a 2-d grid with points labeled by $`(n_x,n_y)`$, all spins corresponding to the lattice points $`(1,n_y)`$ and $`(n_x,1)`$ for $`n_x,n_y=1,\mathrm{},L`$ have been fixed during all the simulation at its initial values $`\sigma =1`$. With this precaution, the structure of the program, that implements PBC, assures the persistence of the frozen boundary.

We have performed the lattice updating applying a well tested head bath algorithm. During the simulation we recorded time series files for the energy $`E`$ and the magnetization $`M`$ defined as

$$M=\frac{q\text{max}\{n_i\}L^d}{q1},$$
(3)
where $`q=8`$ and $`n_i`$ is the number of spins in a given orientation.

Table 2 summarizes the details of the simulations that have performed from $`L=70`$ up to $`L=350`$. The number of production Monte Carlo sweeps varies from $`n_{\mathrm{prod}}=6000000`$ for $`L=70`$, to $`n_{\mathrm{prod}}=32700000`$ for $`L=350`$. Since we took measurements only every $`n_{\mathrm{flip}}=8`$ sweeps, the number of total measurements per run is $`n_{\mathrm{meas}}=n_{\mathrm{prod}}/n_{\mathrm{flip}}`$. We left at least $`20n_{\mathrm{flip}}\tau _\mathrm{e}`$ thermalization sweeps before taking measurementsSokal ; Sokal2 ; Janke2 . To estimate the autocorrelation time of energy measurements $`\tau _\mathrm{e}`$, we have applied two different methods. First, we use the fact that $`\tau _\mathrm{e}`$ enters the error estimate $`ϵ_{\mathrm{JK}}=\sqrt{2\tau _\mathrm{e}/n_{\mathrm{meas}}}ϵ_{\mathrm{naive}}`$ for the mean energy $`<E>`$ of $`n_{\mathrm{meas}}`$ correlated energy measurements of variance

$$ϵ_{\mathrm{naive}}^2=\underset{j=1}{\overset{n_{\mathrm{meas}}}{}}(<E>E_j)^2/(n_{\mathrm{meas}}1).$$

The “true” error estimate $`ϵ_{\mathrm{JK}}`$ is obtained splitting the energy time-series into 50 bins, which were in their turn jackknived jack ; jack2  to decrease the bias in the analysis. The second way of obtaining $`\tau _\mathrm{e}`$ is by a direct computation of the integrated autocorrelation time

$$\tau _e^{\mathrm{int}}=\frac{1}{2}+\underset{j=1}{\overset{k_{\mathrm{max}}1}{}}(1j/k_{\mathrm{max}})\frac{\frac{1}{k_{\mathrm{max}}(j1)}_{i=1}^{k_{\mathrm{max}}(j1)}(E_i<E>)(E_{i+j}<E>)}{<E^2><E>^2},$$

where $`k_{\mathrm{max}}`$ is a suitable cut-off Sokal  around $`6\tau _e^{\mathrm{int}}<k_{\mathrm{max}}<10\tau _e^{\mathrm{int}}`$. The corresponding error in $`\tau _e^{\mathrm{int}}`$ is derived from the a priori formula $`\sqrt{2(2k_{\mathrm{max}}+1)/n_{\mathrm{meas}}}\tau _e^{\mathrm{int}}`$.

In Fig. 1 we present the energy time-series for the $`L=300`$ and $`\beta _{MC}=1.34146`$ simulation run. The expected characteristic behaviour for a first order phase transition can be clearly seen. The system remains in one of the two coexisting phases for a long period of time. The energy histogram for the full series is also presented in this figure. The similar height of the two peaks confirms that the simulation was performed very near the pseudo-critical inverse temperature.

It is instructive to compare the energy histograms corresponding to the adoption of fixed or periodic boundary conditions. To this end we have performed two different Monte Carlo runs, close to the respective pseudo-critical inverse temperatures, which are $`\beta _{MC}^{PBC}=1.342027`$ and $`\beta _{MC}^{FBC}=1.3378`$ for a lattice size $`L=100`$. These simulations has been done using $`8000000`$ production sweeps, with $`n_{\mathrm{flip}}=8`$, discarding the initial $`250000`$ ($`150000`$) sweeps in the case of PBC (FBC) for the thermalization of the system. Both histograms can be seen in Fig. 2. They show the characteristic two-peaks structure. Nevertheless, the latent heat, i.e., the separation between the maximum of the two peaks, is clearly smaller for Fixed Boundary Conditions. This qualitative observation suggest that a simple analysis of the energy histograms of a true first order phase transition simulated with Fixed Boundary Conditions might be misleading. Effectively, if the lattice size is not large enough the energy histogram could show (apparently) a single peak and, in consequence, one can get the erroneous conclusion that the model exhibit a second order phase transition. Nevertheless, even with FBC the evolution of the energy histograms when the size of the system increases shown in Fig. 2 ($`L=100`$), and in Fig.1 ($`L=300`$), exhibit the expected behaviour of a first order transition. This observation may be relevant in the interpretation of the analysis of Baig and Fort BF , where a disappearance of a two peaks structure was observed when FBC were imposed to the system N1 .

In addition to the qualitative analysis of the histograms, we have computed the specific heat, magnetic susceptibility and the Binder kurtosis parameter at nearby values of $`\beta _{MC}`$ by means of standard reweighting techniques Swendsen . They are defined as

$`C(\beta )`$ $`=`$ $`{\displaystyle \frac{\beta ^2}{L^2}}(E^2E^2),`$ (4)
$`\chi (\beta )`$ $`=`$ $`{\displaystyle \frac{\beta ^2}{L^2}}(M^2M^2),`$ (5)
$`B(\beta )`$ $`=`$ $`1{\displaystyle \frac{E^4}{3E^2^2}}.`$ (6)
In table 3 we show the extrema of the magnitudes above defined, together with their pseudo-critical inverse temperatures. The error bars of these quantities have been estimated splitting the time-series data into 50 bins, which were jackknived to decrease the bias in the analysis of reweighted data.

## IV Scaling laws analysis

Once we have the results from the numerical simulation on finite lattices, we can proceed to analyze the data imposing the scaling laws of Table 1.

### IV.1 Analysis of the pseudo-critical inverse temperature

In Table 4 we present the results of fitting the pseudo-critical betas of $`C_{\mathrm{max}}`$, $`\chi _{\mathrm{max}}`$ and $`B_{\mathrm{min}}`$ to the ansatz $`\beta _c+a_1/L+a_2/L^2`$ suggested by the finite-size scaling laws presented in Table 1. Notice that we have performed two set of fits, one for the full range $`84L350`$, and a second including only results from the lattice sizes $`100L350`$. Notice that the fits are extremely good even for the initial range $`84L350`$, but they improve slightly if $`L=84`$ is discarded. Remember that reasonable fits should have a goodness-of-fit NR , $`Q`$, above 0.05. Fig. 3 depicts the fit for $`\beta _\chi ^{\mathrm{max}}(L)`$ in the range $`84L350`$. The exact critical inverse temperature for the 2d 8-state Potts model is $`\beta _c(exact)=\mathrm{ln}(1+\sqrt{(}8))=1.342454\mathrm{}`$. Our results of Table 4 are in perfect agreement with this value.

We have also fitted our data to the ansatz $`\beta _c^{peaks}=\beta (\mathrm{})+\theta _1/L^2+\theta _2/L^4`$ corresponding to the PBC finite-size scaling law. Even though the goodness-of-fit, $`Q`$, obtained does not allow to discard the fits, the infinite volume $`\beta _c(\mathrm{})`$ resulting from them do not coincide with the exactly known value, showing that this ansatz is unsuitable. I.e. for the $`\beta _{\mathrm{max}}^C(L)`$ in the range $`84L350`$, the fit produces $`Q=0.10`$ and $`\beta _c(\mathrm{})=1.342063(11)`$ and for $`\beta _{\mathrm{max}}^\chi (L)`$ in the range $`100L350`$, the results are $`Q=0.21`$ and $`\beta _c(\mathrm{})=1.342079(13)`$.

### IV.2 Analysis of $`C_{\mathrm{max}}`$, $`\chi _{\mathrm{max}}`$ and $`B_{\mathrm{min}}`$

The results of the fits to the specific heat and susceptibility maxima, $`C_{\mathrm{max}}`$ and $`\chi _{\mathrm{max}}`$, together with the kurtosis minimum are summarized in Table 5. As before, we show the fits for two ranges of lattice sizes. Notice that the linear correction coefficients, $`c_1`$ and $`e_1`$, are two orders of magnitude larger than the coefficients, $`c_2`$ and $`e_2`$, of the dominant contribution $`L^2`$. This makes necessary to adjust the data to the ansatz $`C_{\mathrm{max}}(L)=c_0+c_1L+c_2L^2`$, and allows to estimate the corrections to the leading term.

In simulations with PBC, the correction to the leading term is of the order $`\gamma _1/L^2`$. If we fit our specific heat data in the range $`L=126350`$ to the ansatz $`C_{\mathrm{max}}(L)=\gamma _0+\gamma _1/L^2+\gamma _2L^2`$, the goodness-of-fit is $`Q=0.0003`$ with an absurdly high value for $`\gamma _1`$. On the other hand, if we do not allow for a correction term and fit the data to $`C_{\mathrm{max}}(L)=\gamma _0+\gamma _2L^2`$, the goodness-of-fit turns out to be 0.

The work of Medved M  shows that the coefficient of $`L^2`$ in the finite size scaling of $`C_{\mathrm{max}}`$ is related to the latent heat, $`\mathrm{\Lambda }_{\mathrm{FBC}}`$, via $`c_2=(\mathrm{\Lambda }_{\mathrm{FBC}}\beta _c/2)^2`$. In fact, it is the same relationship that holds for periodic boundary conditions Challa ; BKS ; Janke . If we use our estimation $`c_2=0.0389(29)`$ from Table 5 and $`\beta _c=\mathrm{ln}(1+\sqrt{(}8))`$, we obtain for the latent heat

$$\mathrm{\Lambda }_{\mathrm{FBC}}=0.294(11).$$
(7)
Another way of estimating the latent heat is from the direct calculation, right at the transition, of the internal energies per site of the ordered and disordered phases, $`e_{\mathrm{ord}}=E_{\mathrm{ord}}/V`$, and $`e_{\mathrm{dis}}=E_{\mathrm{dis}}/V`$. Of course, the latent heat is just $`\mathrm{\Lambda }=e_{\mathrm{dis}}e_{\mathrm{ord}}`$. Lee and Kosterlitz proposed Kosterlitz  to reweight a given energy histogram until both peaks have equal height. The locations of the two maxima in the histogram can be taken as finite size estimates, $`e_\mathrm{o}(L)`$ and $`e_\mathrm{d}(L)`$, for the infinite-volume limits at $`\beta _c`$ of $`e_{\mathrm{ord}}`$ and $`e_{\mathrm{dis}}`$. The scaling of $`e_\mathrm{o}(L)`$ and $`e_\mathrm{d}(L)`$ for fixed boundary conditions M  as well as periodic boundary conditions Kosterlitz  is $`e_\mathrm{o}(L)=e_{\mathrm{ord}}+O(1/L)`$ and $`e_\mathrm{d}(L)=e_{\mathrm{dis}}+O(1/L)`$.

We smoothed NR  our energy histograms to reduce the noise and searched for $`e_\mathrm{o}(L)`$ and $`e_\mathrm{d}(L)`$. Table 6 shows the estimations that we found. Fitting them to the ansatz $`e_\mathrm{o}(L)=e_{\mathrm{ord}}+k_1/L`$ and $`e_\mathrm{d}(L)=e_{\mathrm{dis}}+k_2/L`$, we obtained $`e_{\mathrm{ord}}=1.6032(48)`$ and $`e_{\mathrm{dis}}=1.3114(92)`$, with goodness-of-fit NR  $`Q=1`$ and $`Q=0.9`$ respectively. Consequently another estimation for the latent heat is

$$\mathrm{\Lambda }=0.292(10).$$
(8)
The agreement with our previous estimation could not be better: it is quite comforting.

R.J. Baxter Baxter1 ; Baxter2  derived an analytical expression for the latent heat of the $`q`$-state Potts model assuming periodic boundary conditions. Numerical evaluations of his expression are tabulated in Wu Wu  and Janke Janke . For $`q=8`$, the latent heat for the Potts model with periodic boundary conditions is $`\mathrm{\Lambda }_{\mathrm{PBC}}=0.486358\mathrm{}`$ Obviously our estimations of the latent heat do not coincide with this value, but it should not be so surprising in view of Fig. 2, where it can be seen that, for $`L=100`$, the distance between peaks for P.B.C. is so different from the distance between peaks for F.B.C. Although such differences could tend towards the same value with $`L\mathrm{}`$, our analysis indicates that in fact they do not.

Notice that, unlike the latent heat, the analytically known infinite-volume critical inverse temperature, $`\beta _c=\mathrm{ln}(1+\sqrt{(}q))`$, for the $`q`$-state Potts model is derived Potts ; Wu2 ; Hintermann ; Baxter2  using the self-duality property of the model, which is independent of boundary conditions when $`L\mathrm{}`$. Let us recall that our estimations of $`\beta _c`$ are consistent with $`\beta _c=1.342454\mathrm{}`$

## V Conclusions

The first order phase transition finite-size scaling laws for Fixed Boundary Condition lattices of Borgs-Kotecky-Medved have been presented, tested and shown to be the only ones that hold for the 2d 8-state Potts model.

It is clear from our analysis that Monte Carlo simulations for FBC are necessarily going to be much more time consuming than those for PBC, since for PBC the system sets into the finite-size scaling region as $`\beta _c(L)=\beta _c(\mathrm{})+\theta _1/L^d`$, while for FBC it does it at the slower pace of $`\beta _c(L)=\beta _c(\mathrm{})+a_1/L`$. Besides, we have found that the latent heat is affected by the boundaries.

It is a pleasure to thank W. Janke and A. Salas for very stimulating discussions. Financial support from CICYT contracts AEN98-0431, AEN99-0766, is acknowledged.