We present numerical simulations of elastic wave propagation, scattering and attenuation in two-dimensional fractured media. Natural fracture systems following a power law length scaling are modeled by the discrete fracture network approach for the geometrical representation of fracture distributions and the displacement discontinuity method for the mechanical computation of fracture-wave interactions. The model is validated against analytical solutions for wave reflection, transmission and scattering by single fractures, after which we apply it to solve the spatiotemporal wavefield evolution in various synthetic fracture networks. We find that the dimensionless angular frequency ῶ plays a crucial role in governing wave transport. When ῶ is smaller than the critical frequency ῶ_c (≈ 5), waves are in the extended mode, either propagating (for small ῶ) or diffusing by multiple scattering (for intermediate ῶ); as ῶ exceeds ῶ_c, the wave energy becomes trapped, entering either the Anderson localization regime (kl^* ≈ 1) in well-connected fracture systems or the weak localization regime (kl^* > 1) in poorly-connected fracture systems, where k is the incident wavenumber and l^* is the mean free path length. Consequently, the inverse quality factor Q^-1 scales with ῶ obeying a two-branch power law dependence, showing significant frequency dependence when ῶ < ῶ_c and almost frequency independence when ῶ > ῶ_c. In addition, when ῶ < ῶ_c, the wavefield exhibits a weak dependence on fracture network geometry, whereas when ῶ > ῶ_c, the fracture network connectivity has an important impact on the wave behavior such that strong attenuation occurs in well-connected fracture systems.