Fluid flow regimes affect the determination of hydraulic conductivity of fractured rocks, and the critical criteria for the onset of nonlinear fluid flow transitions in discrete fracture networks (DFNs) of rocks have yet to be established. First, the factors causing the fluid flow transition regime of fracture intersections and rough surface fractures are analyzed theoretically and numerically. This reveals that the fluid flow regime is governed by the fracture aperture, density of fracture intersections, surface roughness, and Reynolds number ( Re). Then, these identified parameters are redefined in DFN models, and their influence on the onset of nonlinear fluid flow is further investigated by performing computational fluid dynamic analysis. The results show that the fracture intersection and aperture play a more significant role in the linear-to-nonlinear fluid flow transition than the fracture aperture heterogeneity. With the increase in fracture aperture, unevenness of fracture surfaces and connectivity of DFNs, the onset of nonlinear fluid flow appeared at the lower flow velocity. With the Forchheimer equation, it is found that the critical hydraulic gradient Jc, defined as the hydraulic gradient at which inertial effects assume 10% of the total pressure loss, is highly correlated with the fracture aperture, fracture intersection and roughness of the surface. Finally, the mathematical expression of Jc and the Forchheimer coefficients are formulated based on the regression analysis of fluid dynamic computation results, which provides an approach to determine whether the cubic law should be applied as governing equations for the computation of fluid flow in DFNs.