The creep and consolidation behaviors of clays subjected to thermal cycles are of fundamental importance in the application of energy geostructures. This study aims to numerically investigate the physical mechanisms for the temperature-triggered volume change of saturated clays. A recently developed thermodynamic framework is used to derive the thermomechanical constitutive model for clays. Based on the model, a fully coupled thermo-hydro-mechanical (THM) finite element (FE) code is developed. Comparison with experimental observations shows that the proposed FE code can well reproduce the irreversible thermal contraction of normally consolidated and lightly overconsolidated clays, as well as the thermal expansion of heavily overconsolidated clays under drained heating. Simulations reveal that excess pore pressure may accumulate in clay samples under triaxial drained conditions due to low permeability and high heating rate, resulting in thermally induced primary consolidation. Results show that four major mechanisms contribute to the thermal volume change of clays: (i) the principle of thermal expansion, (ii) the decrease of effective stress due to the accumulation of excess pore pressure, (iii) the thermal creep, and (iv) the thermally induced primary consolidation. The former two mechanisms mainly contribute to the thermal expansion of heavily overconsolidated clays, whereas the latter two contribute to the noticeable thermal contraction of normally consolidated and lightly overconsolidated clays. Consideration of the four physical mechanisms is important for the settlement prediction of energy geostructures, especially in soft soils.