This paper presents a stress-induced variable aperture model to characterise the effect of polyaxial stress conditions on the fluid flow in three-dimensional (3D) persistent fracture networks. Geomechanical modelling of the fractured rock is achieved by the finite-discrete element method (FEMDEM), which can capture deformability of matrix blocks, heterogeneity of stress fields as well as sliding and opening of pre-existing fractures. Propagation of new cracks is not required for this study of persistent fracture systems. The deformed fracture network topologies include details of dilation, opening and closing of fracture apertures, from which the local variations in hydraulic apertures are derived. Stress-controlled distribution of fracture apertures is modelled with both fracture-scale and network-scale effects considered. Under a geomechanical condition with low differential stress ratio, fracture porosity is dominated by the fracture-scale roughness. However, with the increase of stress ratio, some favourably oriented fractures are reactivated for shearing, and matrix blocks are promoted to rotate and generate large openings along their boundaries, which tend to be the key contributors to the aperture field in such persistent systems. The flow behaviour is then considered for these stressed but static solid skeletons and is investigated using a finite element solution to the Laplace problem of single-phase fluid flow. The equivalent permeability tensor of each cube-shaped rock mass is computed based on a series of flow simulations under a macroscopic pressure differential applied at opposite model boundaries with no-flow conditions on the remaining boundaries. Components of the permeability tensor are found to vary more than three orders of magnitude with respect to the change of stress ratio. Large aperture channels formed under a critical stress state accommodate significant localisation features in the flow structure of the network. The results of this study have important implications for upscaling permeability to grid block properties for reservoir flow simulation.