Fully unconstrained cell optimizations of uranium dioxide (UO2) bulk model systems using the Gaussian plane waves method are presented. Different effective values are employed for various initial uranium 5f orbital occupation patterns in the framework of a spin polarized DFT+U approach based on the generalized gradient approximation proposed by Perdew, Burke, and Ernzerhof. The detected low-lying states are presented and their properties are analyzed. A new lowest energy state has been detected with a mixed 5f orbital occupation pattern different from the proposed ones so far which shows also favorable properties. The presented results are relevant for future investigations using DFT+U or hybrid functional methods for the study of actinide systems like UO2, especially defective ones.