  The canonical one-band Hubbard model is studied using a computational method
that mixes the Monte Carlo procedure with the mean field approximation. This
technique allows us to incorporate thermal fluctuations and the development of
short-range magnetic order above ordering temperatures, contrary to the crude
finite-temperature Hartree-Fock approximation, which incorrectly predicts a
N\'eel temperature $T_N$ that grows linearly with the Hubbard $U/t$. The
effective model studied here contains quantum and classical degrees of freedom.
It thus belongs to the "spin fermion" model family widely employed in other
contexts. Using exact diagonalization, supplemented by the traveling cluster
approximation, for the fermionic sector, and classical Monte Carlo for the
classical fields, the Hubbard $U/t$ vs. temperature $T/t$ phase diagram is
studied employing large three and two dimensional clusters. We demonstrate that
the method is capable of capturing the formation of local moments in the normal
state without long-range order, the non-monotonicity of $T_N$ with increasing
$U/t$, the development of gaps and pseudogaps in the density of states, and the
two-peak structure in the specific heat. Extensive comparisons with determinant
quantum Monte Carlo results suggest that the present approach is qualitatively,
and often quantitatively, accurate, particularly at intermediate and high
temperatures. Finally, we study the Hubbard model including plaquette diagonal
hopping (i.e. the $t-t^\prime$ Hubbard model) in two dimensions and show that
our approach allows us to study low temperature properties where determinant
quantum Monte Carlo fails due to the fermion sign problem. Future applications
of this method include multi-orbital Hubbard models such as those needed for
iron-based superconductors.
