# Extending the Broad Histogram Method for Continuous SystemsOral presentation on CCP 1998, 2-5 September 1998, Granada, Spain

## Abstract

We propose a way of extending the Broad Histogram Monte Carlo method (BHMC) to systems with continuous degrees of freedom, and we apply these ideas to investigate the three-dimensional XY-model. Our method gives results in excellent agreement with Metropolis and Histogram Monte Carlo simulations and calculates for the whole temperature range $`1.2<T<4.7`$ using only $`2`$ times more computer effort than the Histogram method for the range $`2.1<T<2.2`$. Our way of treatment is general, it can also be applied to other systems with continuous degrees of freedom.

Usually one wants to calculate average values $`<Q>_T`$ of a certain quantity $`Q`$ for a system in equilibrium at temperature $`T`$ according to, for instance, the canonical ensemble

$$<Q>_T=\frac{\underset{E}{}<Q>_Eg(E)\mathrm{exp}(E/k_BT)}{_Eg(E)\mathrm{exp}(E/k_BT)},$$
(1)   
where $`g(E)`$ is the number of states with energy $`E`$, $`<Q>_E`$ denotes the micro-canonical average of $`Q`$ at energy $`E`$, and $`k_B`$ is the Boltzmann constant (in the rest $`k_B=1`$ is taken). It is clear that given $`g(E)`$ and $`<Q>_E`$ one can calculate $`<Q>_T`$ at any desired temperature.

The Broad Histogram Monte Carlo Method (BHMC) developed by de Oliveira et. al.  calculates $`g(E)`$ and $`<Q>_E`$ directly. It chooses a microreversible protocol of allowed virtual movements in the space of states of the system, and counts $`N_{\mathrm{u}p}(X)`$ ($`N_{\mathrm{d}n}(X)`$) as the number of allowed movements that increases (decreases) the energy of the configuration $`X`$ by a fixed amount $`\mathrm{\Delta }E_{\mathrm{f}ix}`$. These movements are virtual, in the sense that they are never performed. They are introduced only to estimate $`g(E)`$ and should not be mixed with the dynamics employed to take the samples. It has been shown  that the total number of possible ways to go up, summed over all states with energy $`E`$, equals the total number of ways to go down, summed over all states with energy $`E+\mathrm{\Delta }E_{\mathrm{f}ix}`$, i.e.

$$g(E)<N_{\mathrm{u}p}(E)>=g(E+\mathrm{\Delta }E_{\mathrm{f}ix})<N_{\mathrm{d}n}(E+\mathrm{\Delta }E_{\mathrm{f}ix})>,$$
(2)   
where $`<N_{\mathrm{u}p}(E)>`$ ($`<N_{\mathrm{d}n}(E)>`$) is the micro-canonical average of $`N_{\mathrm{u}p}(X)`$ ($`N_{\mathrm{d}n}(X)`$) at energy $`E`$. By taking logarithms of both sides of Eq. (2) and by subsequent dividing by $`\mathrm{\Delta }E_{\mathrm{f}ix}`$ one can obtain expressions for $`\mathrm{ln}g(E+\mathrm{\Delta }E_{\mathrm{f}ix})\mathrm{ln}g(E)`$ and $`\beta (E)d\mathrm{ln}g(E)/dE`$ that can be either added up or integrated numerically to obtain $`\mathrm{ln}g(E)`$. The micro-canonical averages involved in the method can be calculated either analytically or by Monte Carlo simulations. In this case four histograms are required, i.e. $`N_{up}(E)`$, $`N_{dn}(E)`$, $`Q(E)`$ and the number of visits $`V(E)`$. This procedure gives accurate results with small computer efforts for many discrete systems .

To extend the BHMC method to continuous systems we propose to choose a protocol of virtual random movements in the space of states of the system, such that for each allowed movement the probability to perform it equals the probability to revert it, i.e.

$`P\left(X_{\mathrm{o}ld}X_{\mathrm{n}ew}\right)=P\left(X_{\mathrm{n}ew}X_{\mathrm{o}ld}\right).`$ (3)   
Next, let us define the probability density function (p.d.f.) $`f^X(\mathrm{\Delta }E)`$ i.e. the probability to obtain an energy change between $`\mathrm{\Delta }E`$ and $`\mathrm{\Delta }E+d\mathrm{\Delta }E`$ starting from the configuration $`X`$. We propose

$`N_{\mathrm{u}p}(X)f^X(\mathrm{\Delta }E_{fix})`$ ; $`N_{\mathrm{d}n}(X)f^X(\mathrm{\Delta }E_{fix}).`$ (4)   
Eq. (2) is still valid, and the method proceeds as before .

To test our ideas we chose the XY-model  with the Hamiltonian

$$=J\underset{<\mathrm{𝑖𝑗}>}{}\stackrel{}{\sigma }_i\stackrel{}{\sigma }_j=J\underset{<\mathrm{𝑖𝑗}>}{}\mathrm{cos}(\theta _i\theta _j),$$
(5)   
where the summations $`<ij>`$ are taken over all pairs of nearest-neighbor sites, $`\theta _i`$ denotes the planar angle of the spin at site $`i`$ and $`J`$ is the maximal energy per bond (in the rest $`J=1`$). As protocol of virtual random movements, we take at random one site $`i`$ with angle $`\theta _{i}^{}{}_{\mathrm{o}ld}{}^{}`$ and choose a new angle $`\theta _{i}^{}{}_{\mathrm{n}ew}{}^{}`$ for it with uniform probability on $`[\pi ,+\pi ]`$. We obtain 

$`f^X\left(\mathrm{\Delta }E\right)={\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}f^i\left(\mathrm{\Delta }E\right)`$ ; $`f^i\left(\mathrm{\Delta }E\right)=\{\begin{array}{cc}\frac{1}{\pi }\frac{1}{\sqrt{A_i^2(\epsilon _i+\mathrm{\Delta }E)^2}}& \left|\mathrm{\Delta }E+\epsilon _i\right|<A_i\\ 0& otherwise\end{array},`$ (8)   
where $`N`$ denotes the number of sites $`i`$, and

$`\epsilon _i{\displaystyle \underset{j}{}}\mathrm{cos}(\theta _{i}^{}{}_{\mathrm{o}ld}{}^{}\theta _j)`$ ; $`A_i\sqrt{\left(_j\mathrm{cos}\theta _j\right)^2+\left(_j\mathrm{sin}\theta _j\right)^2}`$ (9)   
The summation over $`j`$ is always performed on all nearest neighbors of site $`i`$. Finally $`N_{\mathrm{u}p}(X)`$ ($`N_{\mathrm{d}n}(X)`$) equals the average of $`f^i(\mathrm{\Delta }E_{\mathrm{f}ix})`$ ($`f^i(\mathrm{\Delta }E_{\mathrm{f}ix})`$) on all sites.

To estimate the micro-canonical averages we binned the whole energy range of positive temperatures (i.e. energies per bond between $`1.0`$ and $`0.0`$) in windows of size $`\mathrm{\Delta }E_{fix}`$ and performed a micro-canonical simulation on each window as follows . Starting from a configuration inside the window one spin was rotated to a new angle chosen at random with uniform probability over the interval $`[\pi ,\pi ]`$. If the new system energy falls inside the window, the change is performed, otherwise it is rejected. It is clear that this sampling process maintains a detailed balance condition with equal probabilities for the old an the new configurations and, therefore, it samples all the states inside the window with the same probability. An entire lattice sweep was performed by repeating this procedure on all sites.

Fig. 1 compares the specific heat we obtained for the SC $`10\times 10\times 10`$ XY-model with periodic boundary conditions by using the BHMC method (circles), the Histogram Monte Carlo method (HMC)  (squares) and Metropolis  simulations (diamonds).

For the BHMC method $`\mathrm{\Delta }E_{fix}=6.0`$ ($`\mathrm{\Delta }E_{fix}=0.002`$ in energy-per-bond units) was chosen.

The HMC method took $`188`$ sec. per run, a Metropolis for one temperature point $`28.5`$ sec. and the BHMC method $`383`$ sec. per run on a Digital Alpha workstation, i.e. only $`2.0`$ times more computer effort than the HMC method and less than $`14`$ Metropolis simulations. It can be observed that all methods coincide on the temperature range $`2.1<T<2.2`$, i.e. the classical range of validity for this HMC simulation , but the BHMC method is still precise on the whole temperature range $`1.2<T<4.7`$.

In conclusion, our proposal to apply the BHMC method to continuous systems gives excellent agreement with Metropolis and HMC simulations. The strategy proposed is completely general, it can also be applied to other systems with continuous degrees of freedom.