ABSTRACT: This study presents a coupled hydro-mechanical model for unsaturated fractured rock masses. We solve the Richard's equation to compute unsaturated groundwater flow dynamics, incorporating the van Genuchten retention model for determining relative permeability and capillary pressure. We model the mechanical behavior of rock masses based on linear poroelasticity, where the hydraulic and mechanical processes are coupled for both partially and fully saturated conditions. The model has been examined based on a series of benchmark tests, covering fully saturated poroelasticity, static and dynamic unsaturated flow, and solid rock mass deformation. Our model shows a good agreement with analytical solutions, experimental observations, and other independent numerical simulations. We further use this model to simulate a synthetic fractured rock with three distinct regimes identified: fracture-dominated flow at higher saturation, a transitional flow regime, and matrix-dominated flow at lower saturation. The model also demonstrates that effective saturation influences mechanical properties, with Young's modulus increasing and Poisson's ratio decreasing as the saturation reduces.