A two-dimensional discrete fracture network model was established based on the characteristics of the fractured rock mass of a water conveyance tunnel. The constitutive law for two-phase flow of Bingham fluids in a single saturated fracture was extended to fracture networks to simulate the diffusion process of grouts in fractured rock masses. The obtained grouting zones with different grouting times were input into a three-dimensional discrete fracture network model taking into account the coupled hydro-mechanical effect, and the distribution characteristics of stress, deformation and hydraulic pressure in the rock mass were numerically calculated. The results show that the geometry of the grouting zone is anisotropic due to the distribution characteristics of the primary fracture set. The lining can effectively control the deformation of the rock mass, and the maximum displacement is reduced by around 50% compared with the case without lining. The formation of grouting zone can effectively reduce the compressive stress acting on the lining and relieve the stress concentration caused by the external water pressure. These results can provide technical support for the reinforcement design of underground engineering such as tunnels and caverns in water-rich rock masses.