This paper presents a new method for handling the phase change process in numerical simultations of freezing and thawing soils. Moisture and heat transfer in soils subjected to both freezing and thawing are discussed. Simulated freezing results of temperature and total water content (water + ice) are compared with experimental data reported by Jame (1977). Simulated and experimental results were similar. The effects of different time-dependent temperature boundary conditions were evaluated and discussed. The method was used both with abrupt and smooth temperature boundary conditions and both resulted in stable numerical solutions. Finally, results from a simulation of a freezing and thawing cycle are presented and discussed qualitatively.