In the present project, the behavior of an arbitrary f(R) gravity in the presence of a thermal term was studied. The available code describes the phase space of the matter dominated era of the universe, however different cosmological eras can be easily studied by changing a single parameter. The study of the phase space is achieved through the creation of an autonomous system, hence the second time derivative of Hubble over Hubble to the power of 3 is assumed to be a constant value, different for each cosmological era. This particular code may help a beginner understand the process of phase space analysis. For any questions, feel free to send me an email.