Abstract <p>A generalized multiscale finite element method (GMsFEM) is proposed for the numerical simulation of coupled thermo-poroelastic systems with phase transitions in heterogeneous porous media. The model accounts for interactions between heat transfer, fluid flow, and solid deformation, incorporating ice-water phase transitions via a regularized phase function. To resolve fine-scale heterogeneities while reducing computational cost, we employ an offline model reduction framework with localized multiscale basis functions constructed in selected regions. The governing system includes energy balance with latent heat effects, mass conservation with porosity evolution, and momentum equations with thermal stress coupling. Numerical experiments demonstrate that the proposed approach provides accurate approximations of temperature, pressure, and displacement fields compared to fine grid solutions. The method offers a flexible and efficient framework for simulating multiscale thermo-poroelastic problems relevant to geotechnical applications in permafrost-affected areas.</p>