Molecular dynamics (MD) simulations of an MFI-type zeolite have been performed at various temperatures. Crystal structures and vibrational spectra obtained from the simulations were in good agreement with the experimental data. The temperature induced phase transition between orthorhombic and monoclinic phases, which was already known from experimental measurements, was also reproduced. The potential model used in this study was thus sufficiently applicable to molecular simulations of this type of zeolites. Our MD simulations suggest Pn symmetry for the monoclinic phase instead of experimental P2(1)/n. A small change in the thermal expansion coefficient was found on the phase transition. In the higher temperature range, another remarkable change was found : the thermal expansion coefficient became negative. This indicates the existence of a higher temperature orthorhombic phase.