Robust and reliable hydraulic fracturing models that appropriately account for random initiation of fractures, strongly nonlinear coupling among deformation, fracturing and fluid flow in fracture apertures and leakage into porous rock matrix, would be a key step toward developing a better understanding of physics associated with hydraulic fracturing process. In this paper, we present a physics-based hydraulic fracturing simulator based on coupling a quasi-static discrete element model (DEM) for deformation and fracturing with conjugate lattice network flow model for fluid flow in both fractures and porous matrix. The coupled DEM-network flow model reproduces a variety of realistic growth patterns of hydraulic fractures. The effects of in situ stress, fluid viscosity, heterogeneity of rock mechanical properties and injection rate on the fracture patterns will be presented and discussed. In particular, simulation results of multistage horizontal wellbore with multiple perforations clearly demonstrate that elastic interactions among multiple propagating fractures, strong coupling between fluid pressure fluctuations within fractures and fracturing, and lower length scale heterogeneities, collectively lead to complicating fracturing patterns.