Periclase is the terminal component of the ferropericlase, and its chemical composition is MgO. It is well known that there exists a huge difference between the melting curves of MgO determined experimentally and theoretically. A feasible way to clarify the nature of the melting temperature is to investigate the possible new phase of MgO. Meanwhile, it is very important to study the new phase and the influence of temperature on structural stability of MgO in high-pressure condensed matter physics and geophysics. In the present work, we study in detail the phase stability and the possible existing structures of MgO, which include the structure predicted by particle swarm optimization algorithm through using the first-principles pseudopotential density functional method. We find that MgO crystallizes into a rocksalt structure in a pressure range from 0 to 580 GPa and that the CsCl-type structure is of a high-pressure phase at up to 800 GPa. Although an NiAs-type hexagonal phase perhaps explains the volume discontinuity at (170 ± 10) GPa along the MgO Hugoniot in a shock-compression experiment (Zhang L, Fei Y W 2008 Geophys. Res. Lett. 35L13302) and a wurtzite phase perhaps explains the huge difference between the melting curves of MgO determined experimentally and theoretically (Aguado A, Madden P A 2005 Phys. Rev. Lett. 94068501), neither of them is existent in the entire range of pressures studied, according to the thermodynamic stability calculations. The calculations of phonon spectra indicate that the B3, B4, B8 ₁, B8 ₂, and P3 m1 phases of MgO are dynamically stable at zero pressure. That is to say, all of the predicted structures are the metastable structures of MgO. In addition, the high-temperature structural stability of MgO is investigated by using very similar Lewis-Catlow and Stoneham-Sangster shell model potential based on the classical molecular dynamics (MD) simulations. In order to take into account the non-central force in crystal, the breathing shell model is also introduced in simulation. The thermodynamic melting curves are estimated on the basis of the thermal instability MD simulations and compared with the available experimental data and other theoretical results in the pressure range of 0－150 GPa.