Compound flood modelling framework for surface-subsurface water interactions

. Compound floods are an active area of research where the complex interaction between pluvial, fluvial, coastal or groundwater flooding are analyzed. A number of studies have simulated the compound flooding impacts of precipitation, river 15 discharge and storm surge variables with different numerical models and linking techniques. However, groundwater flooding is often neglected in flood risk assessments due to its sporadic frequency - as most regions have water tables sufficiently low that do not exacerbate flooding conditions -, isolated impacts and considerably less severity in respect to other types of flooding. This paper presents a physics-based, loosely-coupled modelling framework using FLO-2D and MODFLOW-2005 that is capable of simulating surface-subsurface water interactions. FLO-2D, responsible for the surface hydrology and 20 infiltration processes, transfers the infiltration volume as recharge to MODFLOW-2005 until the soil absorption capacity is exceeded, while MODFLOW-2005 returns exchange flow to the surface when the groundwater heads are higher than the surface depth. Three events characterized by short-duration intense precipitation, average tide levels and unusually high water table levels are used to assess the relevance of groundwater flooding in the Arch Creek Basin, a locality in North Miami particularly prone to flooding conditions. Due to limitations in water level observations, the model was calibrated based on 25 properties that have experienced repetitive flooding losses and validated using image-based volunteer geographic information (VGI). Results suggest that groundwater-induced flooding is localized, and high groundwater heads influence pluvial flooding, as the shallow water table undermines the soil infiltration capacity. Understanding groundwater flood risk is of particular interest to low-elevation coastal karst environments as the sudden emergence of the water table at ground surface can result in social disruption, adverse effects to essential services and damage of infrastructure. Further research should assess the


Introduction
Flood inundation modelling is of critical importance for better planning, forecasting and decision-making practices (Teng et al., 2017). Scientific and technological innovations in numerical algorithms have continuously improved the performance of 35 physically-based hydrologic, ocean circulation and hydraulic modelling packages to simulate faster and more accurate flood physical processes over the computational domain at various scales and resolutions (Devia et al., 2015). However, most flood inundation models are designed to simulate specific flood hazards (i.e., pluvial, fluvial, coastal, groundwater) independently and are unable to assess complex flood dynamics per se due to code limitations and burdensome compatibility. To address these numerical constraints, some models have the ability to operate as linked units or groups by using coupling schemes (i.e. 40 one-way, loosely, tightly, fully) to build compound models capable of simulating multiple flood drivers (Santiago-Collazo et al., 2019).
Compound floods (CF) are high-impact low-probability events characterized by a non-linearity behavior resulted from the complex interactions of interrelated flood drivers triggered at the same spatial and temporal scales (Field et al., 2012;45 Seneviratne et al., 2012;van Westen and Greiving, 2017;Zscheischler et al., 2018). Research on CF has received increasing attention in recent years due to their adverse impacts at the global scale. Deterministic and probabilistic approaches are preferred frameworks to analyze CF events. Stochastic models through copula-based probability analysis and extreme value theory examine the interrelationship between flood drivers, while physically-based numerical simulations provide a tangible depiction of the flood dynamics for current and future climate projections. Several compound flooding studies have used 50 physically-based hydrodynamic models as the reference model to simulate the combined effects of rainfall-runoff and storm surge (Christian et al., 2015;Gori et al., 2020;Ikeuchi et al., 2017;Karamouz et al., 2015;Kumbier et al., 2018;Olbert et al., 2017). Failure to consider the compound interactions of flood drivers can result in significant uncertainties in the magnitude, timing, and estimation of flood risk (Wahl et al., 2015). Therefore, the transition from traditional univariate approaches to a multivariate perspective is necessary to improve flood hazard understanding and predictions (Bates et al., 2021). 55 The significance of groundwater flooding is rarely disputed as it is only relevant to geographical regions sitting on top of permeable rock that are prone to groundwater emergence (i.e., Miami, Yucatán Peninsula, United Kingdom). Groundwater floods are events limited to prolonged rainfall in low-elevation karst watersheds characterized by unconfined aquifers that experience sudden increases of already high-water table levels above normal conditions (Finch et al. 2004). Although there 60 has been a substantial increase in groundwater flooding literature since the 2000s as well as advances in understanding surface water/groundwater interactions (Brunner et al., 2017;Sophocleous, 2002), relevant knowledge gaps and lack of understanding of this phenomenon persist from the complex relationship between topography and hydrogeology (Bradford, 2002;Hughes et al., 2011;Ó Dochartaigh et al., 2019). The water table response time to hydrological events is controlled by the soil, vegetation and aquifer properties, which influence the infiltration capacity, recharge rate and response time (Nalesso, 2009). Similarly, 3 the groundwater dynamics are influenced by spatial-temporal variations of single or compound flood drivers (i.e. precipitation events, high river levels, above-average tides and sea level rise conditions) over long or repetitive periods of time (Ascott et al., 2017). Thus, the water table response to hydrological mechanisms (García-Gil et al., 2015), system fluctuations and residence time (MacDonald et al., 2014) determine the severity of groundwater flooding.

70
While probabilistic and empirical approaches have contributed to the development of regional groundwater flood maps (Cobby et al., 2009;Jacobs, 2007), physically-based models are scarce. Abboud et al. (2018) found that the June 2013 compound flood disaster in the Elbow River (Canada) was induced by steady precipitation and increased river flow discharges from upstream basins resulting in basement flooding due to the rise of the water table. The combined effects of fluvial and groundwater flooding were not considered in that study since the MODFLOW river package focused exclusively on groundwater flow. 75 Similarly, Yu et al. (2019) applied the coupled surface-subsurface model PIHM to produce a comprehensive groundwater flood risk and damage assessment over the Koiliaris River (Greece). Yang & Tsai (2020) investigated the impacts of water  (Fish and Stewart, 1991), estimation of aquifer parameters to calculate groundwater flow (Cunningham et al., 2004), and statistical analysis of hydrological measurements Melesse, 2012, 2011;Prinos and Dixon, 2016). Similarly, Hughes and 85 White (2016) investigated the effect of pump practices and sea level rise on surface water routing and groundwater interactions in Miami-Dade County (MDC) using MODFLOW. Currently this is the main reference model for MDC regional research and planning purposes in hydrologic, ecologic, and environmental fields. Regarding the study area, Sukop et al. (2018) developed a MODFLOW model that analyzed the current and future response of the water table to rainfall events in a portion of the Arch Creek Basin. The study highlighted precipitation as the main trigger for groundwater-induced flooding, with tidal fluctuations 90 and sea level rise increasing the shallow water table. Researching the flood risk potential from surface-subsurface water interactions in MDC where the water table is near to the ground surface is critical as it could reveal hidden risks from the compound impact of major storms and coastal forcing variables for present and future scenarios.
The main purpose of this study is to present a loosely-coupled modeling framework capable of simulating surface and 95 subsurface water interactions to advance flood vulnerability assessments in regions prone to groundwater-induced flooding and complex compound flooding phenomena. To better understand the effects of the water table in a low elevation coastal zone, a methodology is developed to couple the 2D hydrodynamic software FLO-2D and the groundwater model MODFLOW-2005. The Arch Creek Basin in North Miami was selected as an ideal test site due to its unique hydrogeomorphology, low-lying topography, and high vulnerability to flood events. For the purpose of this analysis, three events characterized by short-100 lived heavy precipitation, regular tide levels and unusually high-water tables were selected to demonstrate the importance of simulating surface-subsurface water interactions in urbanized karst coasts, as high groundwater heads may exacerbate flooding conditions. In the context of this paper, compound flooding is defined as the interaction of overland flow and groundwater emergence, while surge levels are normal and have a minimal influence in the inundation beyond the coast. Finally, the coupled model results were calibrated based on official database from FEMA and validated using volunteered geographic information 105 (VGI) flood observations from the study area. The paper is organized as follows: a complete description of the study area is introduced (Section 2), followed by data collection and the methodology presented in Sections 3 and 4. Model calibration and results illustrate the main findings (Section 5); the discussion compares the results with similar work in the region (Section 6); and the conclusion section includes the advantages, limitations, and future research (Section 7). 110 2 Study Area

Site description
The Arch Creek Basin is located in the northeastern part of MDC, along the coast of Biscayne Bay in the city of North Miami, Florida. Prior to anthropogenic interventions, the Arch Creek River served as an important flow corridor that connected the Everglades to Biscayne Bay, controlling the flood pulse dynamics in the tropical wetland system (Fig. 1). 115 The gradual modifications in land use and the construction of the Biscayne Canal in the 1920s marked the transition of the natural environment to agricultural lands. Variations in the soil moisture conditions and infiltration levels due to changes in the streamflow and drainage patterns in the area caused unsustainable farming practices that lead to a shift to residential development (Fig. 2). The urbanization process along Biscayne Bay required considerable cut and fill earthworks to create 120 ideal urban development conditions (Miami-Dade, 2016).
The Arch Creek Basin (16.95 km 2 ) is a low-lying coastal zone predominantly urbanized (90.1%) and economically diverse.
The population is distributed within five jurisdictions, primarily concentrated in North Miami and North Miami Beach (Table  5 from year to year (1000 -2000 mm/yr), due to tropical storms and extreme hydrometeorological events which highly influence rainfall amounts. A reported increasing trend in rainfall of 2.1 mm/yr from 1906 to 2016, mainly attributable to an increase in wet season rainfall (Abiy et al., 2019), underscores that MDC is under a continued threat from flooding.

Hydrogeology and groundwater
The Arch Creek Basin sits atop one of the most permeable aquifers in the world, known as the Biscayne Aquifer. The Biscayne 135 Aquifer stores 34 billion m 3 of water and spans an area of 10,000 km 2 (Price et al., 2020) tapering from near the center of peninsula Florida towards the eastern coastline where its maximum thickness is about 38 meters (Parker and Cooke, 1944) and hydraulic conductivities exceeding 3,000 m/day (Fish and Stewart, 1991).
The stratigraphy of Biscayne aquifer consists entirely of unconfined permeable limestones of the Fort Thompson and Miami 140 Limestone Formations and contains numerous solution conduits, resulting in rapid infiltration and recharge to the aquifer (Cunningham and Florea, 2009;Hoffmeister et al., 1967;Parker and Cooke, 1944). Recharge via precipitation occurs primarily in the Everglades and groundwater flows eastward towards the shore where it discharges to Biscayne Bay (Cunningham and Florea, 2009).

Flood risk and vulnerability 145
Floods resulting from extreme weather and climate events represent a major threat to low-lying neighborhoods and housing infrastructure in the Arch Creek Basin. Historically, frontal systems and summer cloudbursts are responsible for most of the significant pluvial flooding events in the study area compared to strong tropical systems, with Hurricane Irene (1999), Katrina areas. The capacity of these communities to respond to hydrometeorological phenomena is limited or non-existent, resulting in repetitive negative impacts on livelihoods and residential property, expanding the socio-economic gap and inequality of MDC communities (Keenan et al., 2018).
Frameworks to integrate flood risk mitigation and climate change adaptation strategies are a main component in Miami Dade

Data description
This section presents the data sets required to build the 2D surface-subsurface flood modelling study, including the topographic 165 input, and hydrologic monitoring stations that provide rainfall, tide and well gauge records, as well as verified flood observations.

Topography
The Light Detection and Ranging (LiDAR) digital elevation model (DEM) is a 2-meter spatial resolution produced by Miami-Dade County, Florida. The LiDAR scanner corresponds to the actual bare-earth surface, removing tops of vegetation, 170 buildings, and vehicles, and the project coordinate system is UTM zone 17N Horizontal Datum WGS84. In terms of elevation, the North American Vertical Datum of 1988 (NAVD 88) was assigned as the reference geodetic vertical datum for this study, substituting the original measurements based on the National Geodetic Vertical Datum of 1929 (NGVD 29).

Hydrologic input
Hydrologic modeling included hydrologic conditions of the time periods 1-4 October 2000, 6-8 June 2013, and 23-26 May 175 2020. Boundary and initial hydrologic inputs such as precipitation, tide and ocean-side water levels, and groundwater heads over the specified time periods were obtained from the following sources.

Rainfall
The NEXRAD Radar Rainfall Application is a scientific web map interface developed by the South Florida Water Management District (SFWMD) on which rainfall data is reported based on spatial coverage configurations in the form of the entire district, 180 counties, Arch Hydro Enhanced Database (AHED) watersheds, or Rain Grid. The NEXRAD Rain Grid Layer is a 2 km grid resolution that provides an accurate representation of precipitation every 15 minutes. Rainfall Grid cell 10044042 was selected to characterize the Arch Creek Basin's rainfall conditions.

Tides and ocean-side water levels
DBHYDRO is the official SFWMD repository for climate, hydrologic, and environmental databases 185 (https://www.sfwmd.gov/science-data/dbhydro). Ocean-side water levels were obtained from stations S28_H and S28_T, located in the Biscayne Canal Number C-8 on the Arch Creek southern boundary edge.
The NOAA Tides & Currents website (https://tidesandcurrents.noaa.gov/) provides local water levels, tides, current predictions, and other oceanographic and meteorological conditions. The closest coastal sensor to the Arch Creek Basin is located at the Virginia Key, Biscayne Bay Station (ID #8723214).

Groundwater heads
The Unites States Geological Survey (USGS) National Water Information System (https://waterdata.usgs.gov/nwis/gw), in cooperation with the SFWMD, records daily summary data of maximum groundwater levels in the south Florida region. The groundwater level data was obtained from well G-852 adjacent to the outer western boundary of the study area (Fig. 3). Daily 195 field water level measurements have been recorded since 1973, and 15-minute intervals since October 2007.

Repetitive flood claims
FEMA's severe repetitive loss properties program is designed to provide grants and financial assistance to residential properties that have experienced frequent flood losses over the years (FEMA, 2021). Currently, seventy-five properties have requested 200 financial assistance for property acquisition or to recoup with some of their investments due to flood damages in the Arch Creek Basin (Miami-Dade, 2017). The database stores detailed information on the date of loss, building type, flood zone designation, type of insurance and claim payments between 1995 to 2015, providing a clear footprint of flooding risk hotspots and flood prone communities. This dataset will be used to calibrate the flood inundation maps.

Hydraulic Model: FLO-2D
FLO-2D is a physically-based volume conservation model that combines hydrology and hydraulics to simulate the propagation of water dynamics in urban, riverine, and coastal environments for flood hazard mapping, floodplain delineation, flood vulnerability assessments and mitigation planning (O'Brien et al., 1993). The flood routing model applies the dynamic wave 210 approximation to the momentum equation to calculate the average flow velocity across the square grid system one direction at a time in eight potential flow directions over the floodplain. Hydrological processes are represented as rainfall data over the computational domain or as input hydrographs that can be specified in the channel, floodplain, or along the coasts. Various attributes (elevations, roughness coefficient), components (channel, infiltration, storm drain) and features (streets, hydraulic structures) can be incorporated into the FLO-2D model to produce more refined simulations (O'Brien, 2011). Details are

MODFLOW-2005
MODFLOW-2005 is a fully distributed model developed by the USGS that simulates groundwater flow in aquifer layers (confined or unconfined) using a block-centered finite-difference approach (Harbaugh, 2005). The spatial discretization of the aquifer(s) into grid elements computes the horizontal and vertical flow stresses of the hydrogeological system (water heads, 220 recharge, zetas) at the center of the cell. Similarly, the model offers several solvers for matrix equations, as well as subsidence, observations, surface-water routing, and transport packages. Technical documentation on the model description and groundwater flow equations is presented in Harbaugh (2005).

Coupling surface-groundwater models 225
The main factors determining the coupling process between FLO-2D and MODFLOW-2005 include the algorithms' mathematical solver compatibility to calculate and transfer the exchanged volumes in opposite directions within a fully integrated framework and share consistent spatial and temporal scales.
In terms of the spatial scale, a perfect match between FLO-2D and MODFLOW-2005 surface elevation layers is necessary for 230 the surface and subsurface water interactions to happen. This agreement is subject to identical geographical position, reference system, size resolution, and topographic cell elevations (Fig. 4) pressure and the assumption that the piezometric head is similar to the datum elevation in unconfined aquifers (Nalesso, 2009). 305 The soil saturation percentage is determined based on the surface flow and water table levels. The infiltration calculation continues as long as the water table levels are lower than the terrain elevation. Conversely, the water exchange can also occur in the opposite direction due to a sudden rise in the water table. If the groundwater heads calculated in MODFLOW-2005 are higher than the surface depth in FLO-2D, the depth of water from groundwater will be added to the surface depth. The infiltration calculation is switched off at each node as long as the saturation condition persists, meaning that infiltration will 310 not be calculated until the soil absorption capacity is reestablished.

Model configuration and set-up
The FLO-2D hydraulic model requires a grid of square cells to represent the topography of the floodplain domain. The structured grid size of the computational domain defines the hydraulic model resolution. The LIDAR DTM was used as source 315 floodplain topographic information, and an interpolation algorithm was implemented to produce a resampled DTM floodplain model to be used as input elevation of the hydraulic model. The nearest neighbor interpolation method was selected to resample data from the high-resolution 2 m LiDAR to a 20m resolution (42,621 cells).
In addition to the topographic features, a detailed representation of the built environment is relevant for urban flood modeling 320 in order to simulate the flow wave propagation dynamics realistically. All buildings in the domain (7827 features where the grid element surface area is considered impervious and is removed from potential water interactions.
Rainfall and tides were considered for the hydrologic forcing, setting the precipitation over the grid system and tide levels in 325 the easternmost cells to represent the Biscayne Bay's coastal conditions. Both time series are structured on a one-hour basis and are presented in the following section. The inclusion of the storm drain system, French drains, surface water control structures and pump stations in the modelling framework is beyond the scope of this study.
The infiltration method selected for the case study was the Green & Ampt. Global soil parameters correspond to the urbanized 330 and permeable surfaces characteristics. Considering that MDC is characterized by the water table response to rainfall events, conservative infiltration estimates for the impermeable surfaces were selected to account for the influence of the French drains in the system. For simplicity, the Manning roughness coefficient was assumed as 0.40 for green land cover areas and 0.04 for the impervious urbanized environment, canal bed, and Biscayne coast.

335
Bathymetric measures were available for the Little Arch Creek River. A 1D hydraulic model with natural cross-sections was imported into FLO-2D extending from NE 143 rd Street to structure G-58 located downstream of the Enchanted Forest Elaine Gordon Park. Official bathymetry from the Biscayne shore, Keystone Island, and Sans Souci canals was not available for this study due to jurisdiction restrictions. To compensate for the missing geometry, aerial imagery Google Earth was used to measure the canal's width, while a 10-meter bottom elevation was used as constant depth based on the Miami Florida 340 Intracoastal Topography database from the Oleta River. The boundary area applicable to the Arch Creek Basin was extracted from the regional model using the ModelMuse graphical 350 user interface (Winston, 2009), and the grid spacing across the model was regenerated to a 20 meters resolution. The spatial discretization of the model on the horizontal axis consists of 265 columns and 285 rows for a total of 75,525 cells. The Biscayne aquifer is simplified to be a one-layer 35 meters thickness compared to the three layer units of the regional hydrogeological system. Taking the upper aquifer parameters as reference, the hydraulic conductivities (Kx≈1,890 meter/day), specific storage  (Fig. 7).

Flood events 370
Three flood events characterized by similar high intensity rainfall, tide levels, and unusually high-water table levels with different response times were selected to compare the surface-subsurface model results (Fig. 8). Tropical Storm Leslie (2-4 October 2000) was responsible for one of the most severe events of North Miami in recent history in terms of flooding and property damages, with an accumulated rainfall of 454 mm over 65 hours and an estimated return period of 50 years (Franklin et al., 2001). Similarly, Tropical Storm Andrea (6-8 June 2013) was a short-lived storm that formed in the Gulf of Mexico 375 which produced very heavy precipitation across Broward and MDC (Beven II, 2013), and a total rainfall of 317 mm in the Arch Creek Basin. The 25 May 2020 event is categorized as a 25-year storm with a total daily rainfall depth of 263 mm, producing localized rainfall in the North Biscayne Bay watershed, specifically in the Arch Creek Basin, due to antecedent rainfall conditions since mid-April 2020.

Calibrated coupled surface-subsurface model
Simulating surface-subsurface water physical processes through physics-based flood modelling frameworks is relevant and meaningful to better assess the severity of groundwater-induced flooding in low elevation coastal environments characterized 385 by porous permeable soil. Fig. 9 illustrates the simulated maximum inundation depths corresponding to the magnitudes of Tropical Storm Leslie, Tropical Storm Andrea, and the 25 May 2020 storm. Tide levels per se do not pose significant threats to infrastructure as the coastal waters remain within the channels. Fig. 10 illustrates the emergence of the groundwater heads to the surface as a result of the increase in the water table. The simulation proves reasonable in terms of maximum flood depth and extent due to the similarities in the hydrologic conditions, being Tropical Storm Leslie the most severe of all three storms. 390 FEMA's records on properties subject to frequent flooding were used as a calibration approach to verify a match between the model results with flood observations. Although the available records do not specify the observed inundation depths, an agreement between the property locations and maximum water levels may offer sufficient evidence that the model provides reasonable results (Fig. 11). The calibrated results and display of the water table timeseries in selected locations for Tropical Storm Leslie are shown in Fig. 12-13. 395

Identification of flooding hotspots
The groundwater flood maps for Tropical Storm Leslie (37.17%), Tropical Storm Andrea (13.87%) and the May 2020 event (20.82%) are showed in Fig. 10. The simulation demonstrates that slight variations in the water table depth (Fig. 8) can exacerbate groundwater emergence extent, resulting in ≈ 10 cm across the Arch Creek Basin. Interestingly heavy precipitations 400 scenarios with very high water tables over extended periods of time (May 2020 event) are more likely to trigger groundwater induced flooding compared to very high precipitation with high water table levels (Tropical Storm Andrea). Fig. 11 presents reasonable results between the reported claims and localized flooding, indicating that the housing infrastructure in these neighborhoods are likely to experience additional flood losses at some point in the future. The simulated storm events illustrate that most of the properties experienced moderate to high flood depths (> 0.5 meters) in predefined locations. Although rainfall-405 runoff is the primary source of flooding in the urbanized Arch Creek Basin, abnormally high groundwater levels triggered groundwater-induced flooding near historic waterways and zones below the County's land elevation flood criteria, with flood depths ≈ 1 meter (Fig. 12a -12b). The groundwater plots illustrate the effect of tidal and groundwater boundary conditions on the behavior of the simulated water table, in turn demonstrating the importance of both variables in the modeling set-up and influence in subsurface dynamics, as a cyclic high-low pattern characterizes the tide fluctuations of the Biscayne Bay (Fig.  410 12b -12e) compared to the defined water heads behavior from well G-852 in the western boundary of the domain (Fig. 12a,   12f). In terms of residential damage, Tropical Storm Leslie and Tropical Storm Andrea may be considered the costliest events in the Arch Creek Basin as both account for 60% of the reported claims (25 and 17 respectively) ( Table 2).
Sources of uncertainty in the coupled numerical model could be reduced by increasing the model's resolution and incorporating 415 storm-water infrastructure features (i.e., French drains). For example, the increase of the water table levels could challenge the ability of the storm drain system to convey water towards the Bay, resulting in prolonged flooding conditions, or anti-flood pump stations may alleviate the impacts of flooding by draining water from the streets and swales back to the ocean.
Nevertheless, the repetitive loss records only reflect a small percentage of the damaged infrastructure and cannot be generalized at the Basin scale as the property owners may not meet the criteria to file the claim. Therefore, the presented modelling results 420 fall more on the conservative side and might overestimate the real flooding conditions.

Validation using crowdsourced data from Tropical Storm Andrea
A limited number of real-time crowdsourced flooding observations in the Arch Creek Basin were available for Tropical Storm Andrea (Fig. 13). The visual comparison indicates a spatial agreement between the maximum flood depth of the coupled 425 simulation and the interpreted depth of the crowdsourced data (Table 3)

Discussion
The results of this investigation determined that areas in the Arch Creek Basin below 1.0 meter elevation are potentially vulnerable to groundwater-induced flooding (Fig 10, 12a, 12b). Similar results were obtained by Sukop et al. (2018) who found 440 that precipitation as the main trigger for rainfall-induced and groundwater-induced flooding in elevations below 0.9 meters and 1.5 meters respectively, with tidal fluctuations and sea level rise increasing the shallow water table, contributing to the reduction of the storm drain capacity. The present study also determined that antecedent rainfall events were important in the height of the water table at the start of the rainfall events investigated. A simple groundwater model was approximated to be 2D in the horizontal axis and 1D in the vertical axis. Considering that most of the water table interactions occurred in the upper aquifer layer of the regional model (≈ 7 meters) and the short simulation time of the selected events (64 and 84 hours), we presume that differences in the modelling set up are not significant compared to the regional model and can be considered adequate for the purpose of this study. Additional work may be necessary for the coupled model to be fully operational as the groundwater model should represent the heterogeneous aquifer 450 system to assess the sensitivity of the water table dynamics.
Seasonal water table fluctuations are expected throughout the year, presenting a higher level frequency during the winter and spring seasons due to climate variability and hydrological forcing (Gurdak et al., 2009;Taylor and Alley, 2001). Nevertheless, as we observed with Tropical Storm Leslie and Tropical Storm Andrea, the potential rise of groundwater levels to the surface 455 during dry season cannot be ruled out since the hydraulically non-restrictive nature of the carbonate strata in MDC allows for rapid infiltration and high recharge rates during heavy precipitation events. The hydrologic forcing input and modeling results suggest that the joint occurrence of a high-intensity short-duration precipitation (> 50 mm peak, 250 mm total) with already high groundwater levels (> 1 meter) result in a CF event. Further research on linking multivariate statistical analysis with coupled hydrodynamic modeling frameworks may prove beneficial to identify thresholds that trigger CF conditions (Couasnon 460 et al., 2018;Jane et al., 2020;Moftakhari et al., 2019;Saksena et al., 2019;Sebastian et al., 2017;Serafin et al., 2019).
Although this investigation determined that rainfall and tide levels alone did not produce significant flooding, the modeling efforts did not include storm surge flooding that are often accompany by large hurricanes (Zhang et al., 2013). Nonetheless induced storm surge flooding conditions and sea level rise projections are beyond the scope of this study, future work on 465 assessing the impact of high tide and storm surge induced flooding are fundamental to assess CF events and future flood risk scenarios (Obeysekera et al., 2019).

Conclusions 470
Surface-subsurface water interactions are increasing in coastal cities due to multiple factors related to climate change. The Arch Creek Basin in North Miami, which served as a vital flow corridor that connected the Everglades to the Biscayne Bay, is an appropriate location to study the influence of high water tables in flood conditions. Results corroborate that groundwaterinduced flooding is localized; thus, becoming an underlying condition that must be considered in low elevation coastal karst environments where the water table dynamics are subject to swift fluctuations caused by rainfall events. 475 A knowledge gap regarding a consolidated groundwater modelling framework was identified and addressed by proposing a loosely-coupled flood model that integrates surface hydrology and groundwater. The ability to produce more comprehensive flood hazard mapping from coupled surface and subsurface water interactions is scientifically relevant to professionals in hydroinformatics since it improves the replicability of flood dynamics, setting the path to improve the understanding, 480 prediction, and response time of groundwater levels as a potential trigger to compound flooding phenomena that can exacerbate floodwater depth and areal extent. This work opens new horizons on the development of CF models from a holistic perspective.
The quality and accuracy of flood hazard mapping in urban areas are strictly related to the model spatial resolution considering that the vertical datum and built-up environment influence flow propagation dynamics. A 20-meters grid resolution was 485 selected to balance the computational demands with a certain level of precision without compromising the quality of the simulation. However, the investigation of higher and coarser resolutions in surface-subsurface modelling studies might yield insights into the estimation of inundated areas and time performance at different scales.
Considering Miami's hydrogeomorphology is one of the most complex globally, the compounding effects of flood drivers may 490 respond differently in diverse geographic settings. Therefore, further research should consider the proposed modeling framework to assess the CF risk in different geographical regions prone to multiple flood drivers, specifically in areas that have access to post-event flooding maps in the form of remote sensing products or VGI data for calibration and validation purposes.

495
The contributions of this research are substantial and go beyond the numerical simulation scope, as it supports numerous fields and real applications including flood management, urban planning and design, flood mapping and zoning, disaster risk reduction, flood insurance policies and policy making. The ability to simulate rising groundwater levels may be of great interest to Miami-Dade authorities on the impact of flooded septic systems and pollutants from a water quality, ecological and public health perspective. Ultimately, this research is a small piece of multidisciplinary work that analyzes the ripple effects of 500 flooding in a wide range of fields (such as socio-economic costs, urban and ecological degradation, and health) and can set the basis for prevention, protection, accommodation, and even retreat/relocation policies.

Author contributions.
FP and JO jointly conceptualized the research experiment, from the design of the procedure to the presentation of results. FP 505 gathered and processed the case study data, developed the coupled model framework, calibrated, and validated the simulation results, wrote initial version of manuscript, and produced all figures and tables. NGR provided the technical expertise to achieve the coupling between FLO-2D and MODFLOW-2005. FN, JO, and AM provided guidance and supervised the work of FP. RP and FC shared ideas to improve the results and discussion sections. FN, AM, JO, RP, FC, and TC contributed to the paper revisions. 510 Competing interests. The authors declare that they have no conflict of interest.
Acknowledgements. We gratefully acknowledge Marcia Steelman from MDC for the kind support throughout this research, including the provision of detailed background of the study area, documentation, historic imagery, shapefiles, and 515 crowdsourced data. We thank Angela Montoya from MDC for her helpful assistance on understanding MDC's regional groundwater model using MODFLOW, and Ruben Arteaga from SFWMD for sharing flood protection and planning drainage reports. We would also like to thank our colleagues Michael C. Sukop