We present different methods to increase the performance of Hybrid Monte Carlo simulations of the Hubbard model in two-dimensions. Our simulations concentrate on a hexagonal lattice, though can be easily generalized to other lattices. It is found that best results can be achieved using a flexible GMRES solver for matrix