ABSTRACT Aquifer storage and recovery (ASR) is used world-wide to supplement available water supplies by storing surplus water in aquifers and recovering it during periods of drought and increased demand. The use of ASR as an option for increasing available municipal irrigation and fresh water supplies is threatened as a result of the mobilization of arsenic in some aquifers during ASR. Arsenic is liberated from arsenic-bearing sulfide minerals as a result of the mixing of oxidizing injected water with reducing insitu groundwater. Fracture networks can have significant influence on the migration and distribution of arsenic in the Upper Floridan Aquifer (UFA) during ASR operations through effects on fluid flow, chemical reactions, and transport characteristics. To characterize fracture flow and associated mass transport, numerical three-dimensional models constructed with MODFLOW and FracMan are used to represent fractures in equivalent continuum, discontinuum, and stochastic discontinuum or discrete fracture network (DFN) dual porosity or hybrid models. The geochemical reaction (PHREEQC-2) and transport (MT3DMS) models are coupled to the three dimensional numerical flow model (MODFLOW 2000) as PHT3D- 2003, and utilized to simulate the flow, transport, and inorganic reactions associated with the injection of oxidized water into the UFA of Southwest Florida during ASR cycles. The discrete fractures, implicitly simulated in MODFLOW as high flow zones, are model layers of varying thicknesses with uniform hydraulic conductivity and storage parameters, and as stochastically-generated horizontal and vertical fractures with varying physical attributes including orientation, aperture widths, fracture intensity, and fracture distributions, distributed within a lower conductivity matrix. Discrete fracture networks are simulated with FracMan and the results imported into MODFLOW. Although each fracture zone layer is assigned a unique stochastic distribution of hydraulic conductivity, each model represents a single realization. The FracMan output of stochastic distributions of hydraulic conductance and storage parameters is "upscaled" for use in MODFLOW. The vertical migration of solute due to variations in the density of injectate and groundwater does not appear to be a significant characteristic of the modeled flow system. The modeling results support the hypothesis that arsenopyrite, which is stable under reducing conditions, liberates arsenic during recharge cycles as a result of oxidation. The results also indicate that fracture flow significantly controls the distribution of all solutes affected by the ASR flow system due to the significantly higher transmissivity of the fractures compared to the matrix. The simulated distribution of arsenic in the matrix is significantly less than in the fractures as a result of the limited penetration of oxidized recharge waters into the inter-fracture matrix. Under the simulated aquifer and geochemical conditions, arsenic travels farther from the injection well via fractures than is observed in monitor wells, suggesting that the partially-penetrating monitoring well network does not intercept many of the fractures. The modeled increases in concentrations of arsenic in the ASR wells during the recovery cycles are also consistent with observations. Explicit representation of fracture zones in numerical transport models provides an increased understanding of the flow system and the potential occurrence and distribution of arsenic in groundwater.