Economic and efficient energy production from enhanced geothermal systems in naturally fractured rocks depends on the stimulation of natural fractures by hydraulic fractures to generate complex connected fracture networks. This fracture network provides pathways for maintaining high fluid injection and production rates which lead to high heat extraction rates. However, little work has been conducted to describe the formation of these complex fracture networks formed by the interaction of multiple hydraulic fractures with natural fractures. This paper presents a dynamic (time dependent), fully coupled hydraulic fracturing model that is capable of capturing fracture propagation, fracture deflection, the interaction between propagating hydraulic fractures and pre-existing natural fractures and fluid flow inside the fractures. The proposed model is validated against an analytical model for a penny-shaped hydraulic fracture. This model is then applied to conduct a detailed parametric study to investigate the effect of pre-existing natural fracture networks on the connectivity and propagation of hydraulic fractures. The orientation, length distribution, width distribution and connectivity of the created fracture network is studied. Simulation results indicate that: (1) the hydraulic fracture network is oriented in the direction of the dominant natural fracture network when the length of the natural fractures is large, and this leads to a large connected fracture area; (2) If the hydraulic fracture propagates at a large angle to the natural fractures, it results in a smaller connected fracture area and a complicated hydraulic fracture network that is dominated by shear failure events; (3) Natural fracture density plays a dominant role in the propagation direction of hydraulic fractures and the resulting connected fracture area; (4) Increasing the injection rate increases the connected fracture area and width while creating more complex hydraulic fracture networks and leads to more fracture dilation and less shear opening of cracks. The detailed parametric study helps us better understand the creation of fracture networks and can serve to guide us in hydraulic fracture design and optimization in naturally fractured geothermal reservoirs.