This paper develops a numerical model to study fracture network evolution during the nitrogen fracturing process in shale reservoirs. This model considers the differences of incompressible and compressible fluids, shear and tensile failure modes, shale heterogeneity, and the strength and permeability of both shale matrix and bedding planes through the coupling of mechanical-seepage-damage during fracturing fluid injection. The results show that nitrogen fracturing has a lower breakdown pressure and larger seepage zone than hydraulic fracturing under the same injection pressure. Tensile failure was identified as the major reason for the initiation and propagation of fractures. Ignoring the effect of bedding planes, the fracture initiation pressure, breakdown pressure, and fracturing effectiveness reached their maxima when the stress ratio is 1. Under the same strength ratio, the propagation path of the fractures was controlled by the stronger effect that was casused by the bedding angle and stress ratio. With increasing the strength ratio, the fracture number and shearing of the bedding plane increased significantly and the failure pattern changed from tensile-only mode to tensile-shear mode. These analyses indicated that the fracture network of bedding shale was typically induced by the combined impacts of stress ratio, bedding angle and strength ratio.Energies 2018, 11, 2503 2 of 22 more effective if N 2 is applied for the hydraulic fracturing [22,23]. Different fracturing fluids lead to the different pore pressure distributions during fluid flow due to the various properties such as viscosity and compressibility. According to the principle of effective stress, the difference in pore pressure distribution induces the difference in the stress distribution in rocks, which plays an important role in fracturing behavior during the fracturing process [24][25][26][27]. Lower breakdown pressure, larger seepage zone and more complex fracture network were observed during nitrogen fracturing [24][25][26][27]. However, these studies were performed on simple rock samples and the influence of other factors during nitrogen fracturing has not been well evaluated. It is essential to investigate fracture network evolution during the nitrogen fracturing process.Several numerical models have been developed to investigate fracture behaviors during the hydraulic fracturing process. These models include the boundary element method [28], finite element method [29][30][31], and discrete element method [32,33]. Wang et al. [29] applied a new mathematical model to investigate fracture network evolution of multi-stranded fractures with pre-existing natural fractures. Shalev and Lyakhovsky [31] studied the damage evolution and seismicity phenomenon in hydraulic fracturing using an improved damage model based on the finite element method. However, there are still only a few studies on a numerical model for nitrogen fracturing and the comparison between water fracturing and nitrogen fracturing. On the other hand, shale reservoirs occur at di...