  The Hubbard model arises naturally when electron-electron interactions are
added to the tight-binding descriptions of many condensed matter systems. For
instance, the two-dimensional Hubbard model on the honeycomb lattice is central
to the ab initio description of the electronic structure of carbon
nanomaterials, such as graphene. Such low-dimensional Hubbard models are
advantageously studied with Markov chain Monte Carlo methods, such as Hybrid
Monte Carlo (HMC). HMC is the standard algorithm of the lattice gauge theory
community, as it is well suited to theories of dynamical fermions. As HMC
performs continuous, global updates of the lattice degrees of freedom, it
provides superior scaling with system size relative to local updating methods.
A potential drawback of HMC is its susceptibility to ergodicity problems due to
so-called exceptional configurations, for which the fermion operator cannot be
inverted. Recently, ergodicity problems were found in some formulations of HMC
simulations of the Hubbard model. Here, we address this issue directly and
clarify under what conditions ergodicity is maintained or violated in HMC
simulations of the Hubbard model. We study different lattice formulations of
the fermion operator and provide explicit, representative calculations for
small systems, often comparing to exact results. We show that a fermion
operator can be found which is both computationally convenient and free of
ergodicity problems.
