Development of an ESMF Based Flexible Coupling Application of ADCIRC and WAVEWATCH III for High Fidelity Coastal Inundation Studies

: To enable flexible model coupling in coastal inundation studies, a coupling framework based on the Earth System Modeling Framework ( ESMF ) and the National Unified Operational Prediction Capability ( NUOPC ) technologies under a common modeling framework called the NOAA Environmental Modeling System ( NEMS ) was developed. The framework is essentially a software wrapper around atmospheric, wave and storm surge models that enables its components communicate seamlessly, and efficiently to run in massively parallel environments. For the first time, we are introducing the flexible coupled application of the ADvanced CIRCulation model ( ADCIRC ) and unstructured fully implicit WAVEWATCH III including NUOPC compliant caps to read Hurricane Weather Research and Forecasting Model ( HWRF ) generated forcing fields. We validated the coupled application for a laboratory test and a full scale inundation case of the Hurricane Ike, 2008, on a high resolution mesh covering the whole US Atlantic coast. We showed that how nonlinear interaction between surface waves and total water level results in significant enhancements and progression of the inundation and wave action into land in and around the hurricane landfall region. We also presented that how the maximum wave setup and maximum surge regions may happen at the various times and locations depending on the storm track and geographical properties of the landfall area.


Introduction
A dynamically coupled system of storm surge, surface waves, inland hydrology and numerical weather prediction models is one of the most sophisticated tool to accurately simulate the total water level in a tropical land-falling hurricane inundation study. On top of that, based on the geographical location, other model components may need to be included. For instance, to set up an efficient coastal flooding prediction system for the Alaskan region, inclusion of a sea-ice model is essential [1].
In recent years, Earth System Models (ESM) were proven to be invaluable tools that enabled us to better understand and more accurately predict our environment. An ESM includes a coupled applications that consists of several model components to represent relevant physical processes. The model components are expected to interact with each other resembling what takes place in nature.
There are several Earth System Model software flavors that enable model components to communicate by importing and exporting data [2][3][4]. The Earth System Modeling Framework (ESMF) has been utilized to develop several earth system coupled applications worldwide (e.g., [5][6][7]). To increase ESMF interoperability, the National Unified Operational Prediction Capability (NUOPC) consortium developed a layer consisting of a set of generic components [8]. NUOPC layer is a software wrap around ESMF and was developed collaboratively by several research and operational centers. The primary objectives behind NUOPC design are to be reusable, extensible and a portable framework for ESM coupling.
The NOAA Environmental Modeling System (NEMS) is a coupled modeling infrastructure designed to address increasing needs for prediction of the earth environment at a range of time scales. NEMS includes several external model components that have a primary source code outside NOAA. Therefore, NOAA only needs to maintain and develop the coupling interfaces (so-called model caps) of a given modeling component. In turn, the NEMS ecosystem allows connecting various combinations of model components into a number of different coupled applications to address specific environmental phenomena at specific time scales.
For the first time, we developed a flexible and generic coupled application between ADvanced CIRCulation model (ADCIRC; [9]) and WAVEWATCH III (WW3; [10]) using NUOPC/ESMF technology. The coupled application provides an infrastructure to make future development and inclusion of various model components, such as river and inland flooding coupling, seamlessly possible. Most of the available coupling frameworks including the ADCIRC-SWAN coupled system [11] were developed by substantial changes to the source code of the both SWAN and ADCIRC models. Therefore, updating model components to newer versions or adding new models to the mix (e.g., ice model) would involve a tremendous amount of code development. The novelty of this work lies in the capability of flexible coupling both in terms of minimal changes to the original source code of the model components as well as in flexibility for using different meshes and domain decompositions for all the model components in the coupled application The current development of the NUOPC caps provides the possibility to perform dynamical coupling of ADCIRC and the unstructured version of WW3, as well as various atmospheric models ATM. The cap developed for ADCIRC is capable of importing atmospheric forcing and surface wave fields, and exporting water surface elevation and current velocity to the connected model components. Conversely, the cap developed for unstructured WW3 imports atmospheric forcing, water levels and current, and exports the wave radiation stresses required to force ADCIRC.
The first application of this new coupled system is the so-called Named Storm Event Model (NSEM), which is a high-fidelity model for simulating coastal inundation and total water level. It is being developed to meet the requirements of the NOAA's Consumer Option for an Alternative System to Allocate Losses (COASTAL) Act of 2012. This modeling system includes ADCIRC as the hydrodynamic component, WW3 as the wave model, the Hurricane Weather Research and Forecasting Model (HWRF) as the atmospheric component [12] and in the future the National Water Model (NWM) as the inland hydrological component [13], see Figure 1.
The structure of this paper is as follows. First, we describe the envisioned design of the NEMS ADCIRC-WW3 coupled application and the methodology. Then a detailed description of the ADCIRC cap implementation and available coupling options are given. This is followed by a similar description of the WW3 cap. Subsequently, we present the results of the coupled system. Finally we present verification of the coupled ADCIRC-WW3 application for the laboratory flume case of [14], as well as a full scale storm surge inundation event during Hurricane Ike, 2008 in the US Gulf of Mexico.

Structure of the Coupled Application
A typical NUOPC application includes a number of generic components that provide an interface to the underlying ESMF infrastructure for generating and operating a coupled application in a fairly straightforward and seamless manner. The generic components are defined as follows: • A Driver manages all the components to initialize, run, finalize and keep track of time for exchanging information among model components.

•
Connectors are used to execute field matching, grid remapping and data redistribution among model components. • A Model (cap) wraps each model component code (e.g., ADCIRC and WW3) to provide a generic interface and standard metadata suitable to be plugged into the Driver, and form a multi-model coupled application. • An optional Mediator wraps custom coupling code to calculate quantities which includes data from several model components or requires operations such as time averaging.
The system includes methods and utilities for time management, error handling, high performance inputs/outputs (I/O), grid remapping and field interpolation. Since NUOPC is a layer around ESMF library, function calls to both NUOPC and ESMF are possible and sometimes are necessary.
In this research, we developed a NUOPC application that includes a driver, three NUOPC enabled model components and four connectors. The components are not allowed to directly access each other's data. The only way the data moves in or out of a component is via instances of an ESMF_state class. The state is a container that wraps native data and also includes a metadata to let the other components know about name, coordinates and decomposition of the actual packed data.
The driver component accesses ADCIRC, WW3 and ATMesh model components via their SetServices methods. It reads basic information for how to initialize and run the model components from a configuration text file ( Figure 2). The configuration file contains information about the name of the model components, number of processes to be associated to each model component, the coupling time intervals and the order of data exchange among the components. The driver also initializes the number of connectors by providing the name of the sending and receiving model components. Therefore, for a dynamical two-way coupling between two model components, two connectors are required.
The connector component initializes at the run time by matching the list of available import and export fields advertised by the model components. The connector establishes the connection based on matched import and export fields. The connector also has access to the domain decomposition and computational domain discretization of the connected model components. It will generate a remapping and necessary weight matrices for interpolation of the fields among model components at the initialization phase. In other words, the connector receives exported data in the form of an ESMF_state in the native grid or mesh from the exported model component and passes it to the importing model components in their own native grid or mesh definition. This remapping facility allows the coupled system the freedom to use different meshes for, say, the circulation and wave modeling components, and/or different domain decompositions. This is useful in cases where the wave model component requires a different mesh optimization to resolve its distinct physics, or more computational cores for load balancing.
The ATMesh cap was developed as a placeholder interface for a full live atmospheric model, which was not included in our NSEM application due to scope limitations. This so-called data cap reads weather prediction outputs (from a NetCDF data file), initialize required NUOPC/ESMF objects and provide requested data and information to ADCIRC and WW3 caps via the NOUPC/ESMF backbone. The NWM hydrological model component and its associated connectors are being implemented (Figure 1).

ADCIRC Model
The ADvanced CIRCulation model (ADCIRC) is a finite element hydrodynamic community model originally developed by [9]. ADCIRC is undergoing continuous development by groups of scientists and engineers. Its natural finite element unstructured mesh capability, and several modules specifically addressing various aspects of the coastal flooding and tropical cyclone forcing, make it one of the best tools available for coastal inundation studies. ADCIRC operates in either two-dimensional depth-integrated depth-averaged (2D barotropic) and three-dimensional (baroclinic) modes. In the 2D mode, it solves equations for both water surface elevation and the depth-averaged velocity fields. For more details about ADCIRC governing equations, numerical methods and wave forcing implementation please see [9,11].
ADCIRC is written in modular FORTRAN and supports parallel execution on massive supercomputers using Message Passing Interface (MPI) architecture. The code structure is partitioned in three distinct initializing, running and finalizing phases ready for the ESMF coupling. The model initializes by a call to ADCIRC_Init which also receives an MPI communicator from the driver. The subroutine reads necessary input files for constructing the computational mesh including nodes location and connectivity. It also builds a local and global nodal map to reference which nodes reside on which MPI process, and to identify their global relationships. It reads input information to constrain the model such as bathymetry, meteorological forcing and freshwater inflow and open boundary conditions. As a part of the initialization, ADCIRC also checks and connects to all requested output files that will be used as containers to fill in the model results.
ADCIRC enters the run phase by a call to ADCIRC_Run subroutine, which also receives an argument for the number of time steps (NTIME_STP) for that specific run request. The start time and end time of the simulation is determined during the initialization phase. The model run takes place via a time loop in which, at every time step, a single call to the TIMESTEP subroutine occurs. All the computational steps for applying forcing and boundary conditions to produce the final results are being performed in this subroutine. The ADCIRC concludes its run by a call to ADCIRC_Final subroutine where some of the final post-processing and check for MPI finalizing are performed.

ADCIRC Coupling Interface (Cap)
The ADCIRC NUOPC cap performs the coupling in all the three phases: initialize, run and finalize. In the development of the NUOPC cap for ADCIRC, extreme care and attention were paid to minimize changes to the original ADCIRC code. At the initialization of the NUOPC application, a global MPI communicator is created by ESMF infrastructure and a dedicated set of processes passed to ADCIRC via a MPI communicator based on the number of processes requested for ADCIRC in the configuration file. At the initialization, ADCIRC cap also gets connected to available import and export field matches accepted by the communicators.
After information exchange among the model components, the ModelAdvance subroutine of the ADCIRC cap calls the ADCIRC_Run subroutine to perform the next run interval. Table 1 shows the list of the exported and imported fields currently accepted by the ADCIRC cap. The naming conventions of these variables are defined in the NUOPC field dictionary to allow interoperability with other NUOPC components. We modified and tested ADCIRC preprocessing and main model code to accommodate various coupling arrangements.
ADCIRC uses drag law to calculate wind stress.
where U 10 is the wind speed at 10m elevation, C d is the drag coefficient and ρ is the water density.
The NWS input parameters in fort.15 input file are described in Table 2 [15]. NWS is the parameter controlling whether wind velocity or stress, wave radiation stress and atmospheric pressure are used to force ADCIRC.

WAVEWATCH III
WAVEWATCH III (WW3) [10] is a third-generation spectral wave model that solves the wave action balance equation that accounts for the growth, propagation, non-linear interaction and dissipation of wind waves in the ocean by: where N(k, θ) is the wave action density spectrum, related to the wave energy density spectrum F(t, x, y, k, θ) where N(t, x, y, k, θ) = F(k, θ)/σ. Here ∇ h is horizontal derivative, x and y are horizontal coordinates, k is the wavenumber andk its propagation speed due to depth-or current-induced Doppler shifting, θ is the wave direction andθ its propagation speed due to depth-or current refraction, σ is the wave frequency, c g is the group velocity vector and U is the depth-averaged current velocity vector. On the right-hand side, S represents the sum of source terms, including wave growth, nonlinear interaction and dissipation. WW3 was originally developed on a regular grid for global operational wave forecasting, with 2-way nesting for regional applications. More recently, it has been extended to curvilinear grids for Arctic applications [16], as well as unstructured meshes for high-resolution coastal application [17]. Most recently, its traditional "card deck" MPI parallel implementation [18] has been supplemented with a more conventional domain-decomposition approach using ParaMETIS [19] equipped with an optional implicit equation solver. Along with these improvements in numeric, source terms suitable for nearshore application have been added, including depth-induced breaking, reflection, three-wave nonlinear interaction and wave-ice interaction, amongst others. The WW3 code is written in modular FORTRAN, similar to ADCIRC, and is broken up into a collection of sub-programs which are run in sequence to carry out a simulation. The most important of these are ww3_grid, which compiles the computational mesh, and physics and numerics settings into mod_def.ww3, a binary resource file, ww3_prep, which preprocesses all forcing files, ww3_multi, the multi-grid core wave model and ww3_ounf and ww3_ounp, which are NetCDF post processing routines. To comply with the ESMF protocol, the core wave model ww3_multi has been broken up into w3init, w3wave and w3final to perform the main steps of model initialization, model advancing and model finalization, respectively [20]. During the initialization step with w3init, the configuration of the computational mesh file, including node indices, geographical location and the mesh connectivity is read from the binary resource file mod_def.ww3. Using the domain decomposition from the PDLIB library, a local and global nodal map is built that references which nodes reside on which MPI process, and how they are related globally. During the model advance step w3wave, forcing fields such as water depth, wind velocity and currents as well as boundary conditions are updated, followed by the solution of the wave action equation equation by means of fractional stepping. Output is written to a set of binary output files for later post-processing (using ww3_ounf and ww3_ounp). Model finalization is completed by calling w3final.

WAVEWATCH III NUOPC cap
The WW3 NUOPC cap carries out the coupling in the three phases (initialization, advancing and finalizing) described above. A regular grid version of this cap was developed by [20] for global and regional-scale NUOPC applications. In the present work, this cap was extended to support unstructured meshes and domain decomposition for high-resolution coastal modeling. During the initialization of the NUOPC application, import and export meshes are defined based on the PDLIB decomposition, a global MPI communicator is created by the ESMF infrastructure, and a dedicated set of processes are passes to WW3 via a MPI communicator based on the number of processes requested for WW3 in the configuration file ( Figure 2). During this initialization step, the WW3 cap is also connected to available import and export field matches accepted by the communicators.
After the information exchange among the model components, the ModelAdvance subroutine of the WW3 cap calls the w3wave subroutine to perform the next run interval. Table 1 shows the list of the exported and imported fields accepted by the WW3 cap for the current application. It is noted that a larger set of import and export variables is supported by the WW3 cap in general, including surface roughness variables, Stokes drift and bed roughness, used in other NUOPC applications featuring wind waves. For more details, see [20].

Wave-Induced Stresses
Breaking waves transfer their momentum to ocean currents. From spectral wave models such as WW3, the radiation stress vectors can be evaluated from the computed wave energy density as follows [10]: where ρ w is the water density, g is the gravitational acceleration, the directional wave energy density spectrum F(t, x, y, k, θ) = σN(t, x, y, k, θ), d is the water depth and n is the ratio of the wave group c g to wave phase speed c p for a given depth and frequency, given by: n = 1 2 + kd sinh(2kd) . In order to account for the impact of the waves on the mean circulation, the spatial gradient of wave radiation stress per unit area τ s,waves is calculated as: and incorporated as additional depth-integrated stresses alongside wind stresses and bottom stresses into ADCIRC's Generalized Wave Continuity Equation and vertically-integrated momentum equations, following [11].

Model Setup
We utilized the existing NOAA's Hurricane Surge On-demand Forecast System (HSOFS) unstructured triangular mesh as the base of the computational domain for the coupled model setup. The HSOFS mesh covers the entire Gulf of Mexico and extends into the Atlantic Ocean to the approximate longitude of 65°W, allowing for appropriate generation of storm surge from atmospheric effects over a large region. The HSOFS mesh has 1.8 M nodes and covers the shallow coastal regions up to a topographic height of 10 m above local mean sea level with the finest mesh resolution of approximately 250 m (See Section 7.2).

Results
The successful implementation of the ESMF based coupled application and the importance of dynamical coupling of surge and surface waves on the spatial extent of the inundation and active wave action area were investigated for Hurricane Ike, 2008.
Hurricane Ike was a powerful tropical cyclone that swept through portions of the Greater Antilles and Northern America in September 2008, wreaking havoc on infrastructure and agriculture, particularly in Cuba and Texas. As the ninth tropical storm, fifth hurricane and third major hurricane of the 2008 Atlantic hurricane season, Ike developed from a tropical wave west of Cape Verde on 1 September and strengthened to a peak intensity as a Category 4 hurricane over the open waters of the central Atlantic on 4 September as it tracked westward. Several fluctuations in strength occurred before Ike made landfall on eastern Cuba on 8 September. The hurricane weakened prior to continuing into the Gulf of Mexico, but increased its intensity by the time of its final landfall on Galveston, Texas on 13 September.
The wave-surge coupled application (hereafter "Fully coupled") and stand alone models (hereafter "Stand alone") were forced with an identical HWRF meteorological forcing (See Section 7.2.1). As a reference, the ADCIRC model results forced with tidal boundary condition is also presented (hereafter "Only tide").
The maximum surge level during the simulation is calculated by subtracting maximum water level of the "tide only run" from the maximum water level of the "fully coupled run". The results reveal that the most severe inundation during Ike, 2008 with more than 6 [m] above the maximum tide level was took place on the east side of the hurricane track in the region between the Galveston Bay, TX and the Sabine Lake, TX.
The maximum wave contribution to total water level is calculated by subtracting water level from Stand Alone run (ADCIRC forced only by wind and atmospheric surface pressure) from Fully coupled run results (Figure 3b). It is shown that contribution of the wave induced stress on inundation extends the maximum inundation limit landwards where the combined hurricane wind stress and pressure deficit are close to their maximum influence. This is significant because it shows how wave induced momentum released from breaking waves increases the total water level which in turn causes the wave active breaking region to advance further into the land and release wave action on the structures further landwards.
Wave height significantly enhanced due to the dynamical coupling of the surge and wave components. Figure 4a,b represents the maximum wave height difference during the whole storm between Fully coupled and Stand alone (which includes no tidal or surge/inundation effects) cases. Wave height increases more than 2 m along the track as well as in the east of the Galveston Bay at the coastal and landfall region. The comparison at the 6 quick deployed wave gauges shows significant contribution of the surge on the eastern side of the hurricane track in the nearshore region (See Section 7.2.3).
It should be noted that wave setup contribution seen in Figure 3b at the Mississippi and Atchafalaya rivers delta region occurred hours before the actual landfall therefore the maximum surge and wave setup in this region did not happen at the same time. However, it also points to the possibility of experiencing large swells and rip-current events hours before the hurricane makes its actual landfall.
To further analyze this mechanism, we plot changes of the total water level (TWL) and wave height for the Fully coupled and Stand alone cases at a transect shown in Figure 5 (The location of the transect is shown in Figure 6). We also shown the topobathy values along the transect in a positive upward vertical coordinate system where Z = 0 m located at the local mean sea level (black line in 5b). We plotted the only available High Water Mark observation in the 1 km radius from the transect with red squares (location shown in Figure 5). The TWL and wave height line plots start from the land (kilometer 0) in which ground level is almost 6 m above the local mean sea level (latitude ∼29.80 N) and continues towards the ocean to 8 m water depth (kilometer 30) in the ocean side (latitude ∼29.55 N). The Fully coupled solution for the total water level continuously shows enhancement over the Stand alone results (Figure 5a). Around the shoreline where the Z ∼0 m and farther off-shore, the Fully coupled model shows greater total water level in comparison with the Stand alone solution. This directly relates to wave height and its dissipation shown in Figure 5b.
We see that Stand alone model (forced without water surface elevation from surge component) shows wave height of almost zero right after waves cross the shoreline (landwards of the Z = 0 m; Kilometer ∼24). On the other hand, the wave height from the coupled model shows significant wave height of ∼2.2 m in the same region (landwards of the Z = 0 m; Kilometer ∼24). This pattern is visible for the wave height evolution even more landwards. This mechanism in which greater total water level (∼0.5 m) leads to the potential for more active local wave generation and propagation. Therefore greater wave dissipation and release of wave action trigger further enhancement of total water level in the inundated. We also looked at a time series of water level at "P" shown as a blue triangle in the Figure 6. This point is located very close to the maximum spatial extent of the inundated region ( Figure 7) which is helpful to examine time variation of the total water level for Fully coupled and Stand alone cases. The ground level is also plotted by a black dashed-line which shows that P is not wet before the storm as it is located ∼3 m above the local mean sea level. During the land fall both Fully coupled and Stand alone cases produced total water level that inundated the area. Fully coupled shows inundation of ∼1 m above the ground level which is ∼0.5 m above the water level resulted by the Stand alone case at the time of the peak of the storm.      Figure 6). Black dahsed-line is the ground level. All model configurations and results are pre-decisional and for official use only.

Model System Validation
We verified the coupled ADCIRC-WW3 application in a step-by-step manner. In the first step, we verified all the ESMF intermediate exchange ESMF_state fields before sending and after receiving by the other model component. Then we performed a basic verification using a small setup to make sure the coupled ADCIRC-WW3 application runs smoothly. Finally, we switched to the full scale inundation test case for Hurricane Ike, 2008.

Laboratory Case (Boers 1996)
To validate the ADCIRC-WW3 coupled system, we first compared the coupled system results against the laboratory flume experiment of [14]. This experiment was carried out at Delft University of Technology in Spring 1993 to investigate the interaction between wind waves and the mean circulation via wave dissipation and radiation stress transfer. The geometry of the flume is shown in Figure 8c. The bathymetry represents an immovable (concrete) profile of a typical barred beach along the western Dutch coastline. The [14] experiment features three test cases, 1A, 1B and 1C, in which 1A and 1B represent violently-breaking, locally-generated wind waves, and 1C represents a mildly-breaking swell. Table 3 shows the wave height and peak period parameters of the imposed wave conditions at the wave maker, while Figure 9 shows their corresponsing energy spectra. Figure 8a,b shows the domain decompositions of the WW3 and ADCIRC model components, respectively. In both cases the decomposition was generated using METIS [21], but recall that in NUOPC/ESMF they do not necessarily need to match-for WW3 the decomposition was done for 24 cores, while for ADCIRC it was done for 47 cores. The wave component of the coupled ADCIRC-WW3 was forced with the observed spectra at the upstream wave maker, and the water level in the flume is initially set to rest. The WW3 component is configured with a time step of 1 s. The relevant wave physics processes are depth-induced breaking and three-wave nonlinear interaction, using the source term formulations of, respectively, Battjes and Janssen (1978), with γ BJ = 0.80 and α = 1, and [22]. The ADCIRC model component was run at a time step of 1 s, and forced from rest by only the coupled wave radiation stresses. The NEMS coupling time step for the two model components is 1 s. Figure 10 shows that the modeled significant wave heights and the wave setup produced by coupled ADCIRC-WW3 application are in good agreement in terms of significant wave heights and the maximum water surface for all three cases. For the more energetic wave case 1B, wave heights are somewhat underestimated offshore of the outer bar. This overestimated dissipation appears to result in the overestimation of the wave setup over this region leading up to the outer bar. Conversely, in the milder-breaking swell case 1C, we observe an overestimation of wave heights inshore of the outer bar, presumably due to insufficient depth-induced breaking. This results in a slight underestimation of nearshore water levels. Nevertheless, these results show good skill of our coupled ADCIRC-WW3 application in capturing the wave height and wave setup evolution correctly.    Table 3.
All model configurations and results are pre-decisional and for official use only.  Table 3. All model configurations and results are pre-decisional and for official use only.

Full Scale Inundation Case
Here we present the results of the full scale inundation study for Hurricane Ike, 2008. We utilized the unstructured triangular mesh used by NOAA's operational Hurricane Surge On-demand Forecast System (HSOFS). We applied 8 tidal constituents (M2, K1, O1, P1, Q1, N2, S2, K2) at the model open boundaries (Figure 11). We also included detailed information regarding ADCIRC, WW3 and coupled application in Table 4.
For the Ike test case, we associated 720 processes to each of the ADCIRC and WW3 models and 1 process to atmospheric data component for reading the forcing data. The domain decomposition algorithm performed within both models as they initialize and start their run sequence independently. Due to the flexibility that NUOPC/ESMF coupling framework provides, there are no restrictions on the wave and surge models subdomain number or decompositions (see Figure 12). Based on the current domain decomposition configuration, it takes ∼4 h of computation for a day of real time simulation. Please note that we are working on the load balancing aspect of the simulation. For instance, here we use a similar number of processes for wave and surge models and this is not necessarily the best approach when it comes to optimization of the spectral wave model and storm surge model coupled application.

Atmospheric Forcing
The atmospheric forcing for this study is provided by NOAA's Hurricane Weather Research and Forecasting (HWRF) modeling system, coupled to the MPIPOM ocean model and empowered by a movable multi-level nesting technology [23,28]. The model grid is triple-nested using telescopic, two-way interactive horizontal grid resolutions from synoptic with 0.18 • resolution as the outer domain (spanning about 75 • × 75 • ), to moving storm nest with 0.06 • resolution (10 • × 10 • ) and core of about 6 • × 6 • with 0.02 • resolution. These nests follow the hurricane best track, ensuring the highest resolution around the eye of a hurricane. In this study, we have interpolated the hourly HWRF model outputs from multiple cycles initiated with analysis data and 9 forecast time steps. Every 6 h, reanalysis data from the next cycle are smoothly ramped into the wind and pressure fields. The atmospheric forcing has been validated against National Data Buoy Center (NDBC) and satellite altimeter data. We extracted wind velocity at 10 m height and surface pressure from the original GRIB2 output files and saved them in NetCDF format. The ATMesh NUOPC data cap reads the meteorological forcing from NetCDF file and provides it to ADCIRC and WW3 caps at every coupling time step. The HWRF model was forced with initial and boundary conditions provided by NOAA's Global Forecast System (GFS) with 0.5 • spatial grid resolution. Figure 13 gives an impression of the quality of the HWRF atmospheric forcing fields used to force the coupled ADCIRC-WW3 model, by comparing 10 m wind speed and direction against NDBC buoy observations in the Gulf of Mexico (see Figure 4 for locations). We can see that the model and observation are in agreement, in particular for wind directions. However, a tendency to overestimate the wind speeds at the storm peak is found at the mid-Gulf NDBC buoys 42001 and 42002, and the shelf buoy 42019. The high bias in the wind speed is particularly evident at landfall, as seen at NDBC 42035, located on the shelf just offshore of Galveston. It would be expected that these overestimated winds would lead to a degree of overestimation of the locally-generated significant wave heights. It is also worth mentioning that the NOAA buoy 42035 broke free during the storm (See https: //www.ndbc.noaa.gov/hurricanes/2008/ike/). Therefore, the observation time series might not be at the same location of extraction coordinate from the models.

ADCIRC Model Component Validation
We validated the total water level results against two sets of water level observations. The tidal gauge time-series were measured by the NOAA's Center for Operational Oceanographic Products and Services (CO-OPS) and High Water Marks (HWM) were measured and provided by the United States Geographical Survey (USGS).
Comparison between the total water surface elevation from the coupled ADCIRC-WW3 application and tidal gauges time series for four stations close to the actual hurricane track (location of the tidal gauges in Figure 6 are presented in Figure 14a-d). The modeled water surface elevation is in agreement with the observation in terms of the level and timing of the peak of the storm inundation. The forerunner surge generated a freely propagating continental shelf wave with greater than 1 [m] elevation in Ike, 2008 case that traveled coherently along the coast to Southern Texas, and was 300 km in advance of the storm track at the time of landfall. We diagnose the forerunner surge as being generated by Ekman setup on the wide and shallow shelf. The longer forerunner time scale additionally served to increase water levels significantly in narrow-entranced coastal bays. We attribute the majority of the mismatch prior to the main peak of the storm to the fact that the current 2D ADCIRC implementation is not resolving the Ike forerunner. It should be noted that our concentration here is to ascertain general performance of the coupled model. In terms of the Ike forerunner, we reference our reader to [29].
High Water Marks are an important source of observations for validation and enhancement of the storm surge and flood inundation studies. After significant flooding due to a land-falling hurricane, a rapid High Water Mark (HWM) data collection by USGS takes place to document the event and to help improving future disaster preparedness activities. Comparison of the total water level from the coupled ADCIRC-WW3 application and HWMs observation reveal that both are in general agreement in particular around the hurricane track in the land-fall region ( Figure 15).
We also compared scatter plots and a number of statistical metrics of the total water levels and HWM data ( Figure 16). Both Fully coupled and Stand alone results show underestimation of the model in comparison to observations. However, upon inclusion of the wave forcing, an improvement in the error statistics is found, such that a reduction in the relative bias (RB) from −0.594 m to −0.392 m, and a reduction root mean square error (RMSE) from 0.899 m to 0.832 m were presented (See Appendix A).

WW3 Component Verification
The coupled ADCIRC-WW3 coupled application performance with respect to the wave modeling aspects was firstly assessed in the offshore and over the shelf of the Gulf of Mexico using NDBC monitoring buoys, and secondly at the location of landfall by using rapid-deployment pressure sensors, placed by [29].
Considering the NDBC monitoring buoys, four groups can be distinguished. At mid-Gulf (NDBC 42001 and 42002), we find some overestimation of significant wave heights during the passing of the hurricane in response to the somewhat overestimated wind speeds seen above (Figure 17a). Similarly, the peak wave periods are overestimated, reflecting excessive wind-sea growth (Figure 17b). However, wave direction is reproduced well (Figure 17c). Interestingly, in the far-field, towards the eastern half of the Gulf (NDBC 42036, 42039 42040), significant wave heights tended to be underestimated during the passing of the storm, while the peak period and direction were captured well. Moving onto the shelf in the region of landfall (NDBC 42019 and 42020) the agreement of the model results with the observations improved in general, although there is still overestimation in the significant wave height and peak period. Finally, at NDBC 42035, just offshore of the landfall location at Galveston, all wave model parameters are in good agreement with the observations at the storm peak.
The significant wave height is underestimated at this shallow water location just ahead of the storm peak which could be related to the forerunner not being resolved by the surge model component. Regarding the differences between the Fully coupled and Stand alone versions of the wave model, these can only be seen, as expected, at the two shallower stations NDBC 42020 and 42035. As expected, the increased water levels seen above increase the modeled significant wave height at these stations ( Figure 17). Moving more towards shoreline, model results at the pressure sensors deployed in the surf zone by [29] show a significant sensitivity to the inclusion of variable water levels from the coupled ADCIRC model ( Figure 18). Referring to Figure 4b for the station locations, the most significant of these are the stations ANDKNDY-X, Y and Z located under the eastern half of the landfalling hurricane. For these three stations, we can see a large influence of the added surge level on the significant wave height at the storm peak, in all cases improving the agreement with observations. We note that, as discussed above, the stations U to Z also might relate to underestimation of the water surface elevation during the Ike forerunner. Stations S and R, located to the south of the landfall location received mostly offshore winds and a water level set-down. The significant wave heights at these stations are reproduced well.  Figure 19 provides a more in-depth view of the results at the [29] nearshore stations by considering their spectrograms. This figure shows a number of important features of the nearshore wave field transformation under both the Fully coupled and Stand alone models. Up to 13 September, the simulated variance density at the northerly stations ANDKNDY-X, Y and Z is underestimated by both the stand-alone ( Figure 19b) and coupled (Figure 19c) models in connection with non-resolved Ike forerunner effects.
Due to this water level underestimation, in turn, the water depths are underestimated, and hence bed friction and depth-induced breaking are overestimated in the wave model component. After 13 September, the variance density abruptly increases, as the main storm peak arrives with landfall. At this time, the coupled model (panel c) shows greater levels of variance density than the Stand-alone model (panel b), which agrees better with the observations. A final important feature of the nearshore spectra is the frequency upshift of the spectral peak, due to the nonlinear three-wave interactions. This strong upshift at stations V, W, X, Y, Z is not seen as intensely in the model results, indicating an underestimation of the magnitude of this process by the LTA nonlinear interaction source term. Interestingly, the Stand-alone model captured this upshift process somewhat better than the coupled model, since the depth underestimation in the former fortuitously enhances the computed nonlinear interaction (see stations X, Y, Z).

Summary and Conclusions
We developed a flexible coupling application for coastal inundation studies in the framework of NEMS using NUOPC/ESMF infrastructure. This application includes NUOPC model interfaces (or caps) for ADCIRC and WW3, and also a data cap to read and provide atmospheric forcing fields from HWRF model results. The application is examined using a standard laboratory flume case for set-up and wave dissipation, and validated for Hurricane Ike, 2008 on NOAA's HSOFS mesh, a 1.8 M node triangular mesh with a nominal resolution of ∼500 m. The model skills and improvement due to wave effects on the final inundation were examined and discussed using time series from tide gauges and High Water Mark observations. We conclude that the coupling improves the performance of the ADCIRC and WW3 model components. The [14] case shows that the coupling between these two models behaves correctly under laboratory conditions.
In the Ike field case, the simulated water levels of ADCIRC showed generally better agreement with observations upon inclusion of the coupled wave effects from WW3. Accounting for the forerunner effect is expected to improve the results further. Conversely, the coastal wave field simulated by WW3 improves significantly upon coupling with ADCIRC, in particular in the surf zone and in inundated regions.
In this research our main goal is to introduce the flexible coupling framework which has been developed based on NUOPC/ESMF to dynamically couple ADCIRC and WW3 for the first time. For Both ADCIRC and WW3, we employed standard physics packages and settings. More work is in progress to find optimized settings for the coupled system by simulating more recent storms (See: [30][31][32]).
As outlined in the Introduction, the first application of the NUOPC model combination (or App) described here is the Named Storm Event Model (NSEM), a high-resolution hindcast model which is being developed to meet the requirements of the COASTAL Act (2012). However, we anticipate that these flexible and generic NUOPC coupling caps for ADCIRC and WW3 will enable future development and seamless inclusion of various additional model components for forecasting applications. For example, the current cap development provides the possibility to perform three-way dynamical coupling of ADCIRC, WW3 and atmospheric prediction models. In the hindcast application described here, we limited the atmospheric model, run offline, to a one-way forcing to ADCIRC-WW3 via a data cap. However, the ADCIRC and WW3 caps are capable of dynamically importing atmospheric model forcing, and exporting current velocity, water surface elevation (accounting for inundated zones which alter the surface roughness) and enhanced sea surface roughness due to the presence of waves (Charnock parameter) back to the connected atmospheric model component.
Further future coastal applications of the presented flexible NUOPC modeling frameworks are incorporating processes such as river and inland flooding and sea ice. As part of the Named Storm Event Model, work is currently underway to incorporate river flow into this coupled system by including a NUOPC cap for WRFHydro. This extension will allow the modeling of complex coastal flooding events such as Hurricane Harvey (2017), which featured extreme rainfall and inland flooding alongside the surge and wave forcing from the ocean. Sea ice has significant and complex impacts on coastal surge [33] as well as dissipative and dispersive effects on the coastal wave field [16]. The flexible coastal modeling system presented here can be extended to include these coupled processes via an existing NUOPC cap for the sea ice model CICE [20]. Funding: This work has been supported by The Consumer Option for an Alternative System to Allocate Losses (COASTAL) Act Program and Water Initiative project within the National Oceanic and Atmospheric Administration (NOAA). The WaveWatch III development has been supported by the US. Army Corps of Engineers (USACE).

Acknowledgments:
The authors acknowledge Andrew Kennedy for providing the in-situ data. The results presented here are preliminary and for presenting the incremental developments of coupled application which is under development to serve above mentioned projects. The authors also would like to acknowledge Chris Massey for his constructive communications.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A. Metrics for the Evaluation of Data-Model Agreement
In order to assess model performance for water levels, root mean square error (RMSE), BIAS, relative BIAS (RB), Correlation (Cor), Index of Agreement (IA) and peak error (Peak) were used.
The RMSE is given by where M i is the modeled data, O i is the measured data and N is the total number of data. BIAS shows the systematic deviation from the observations and is given by Relative BIAS (RB) shows relative systematic deviation from the observations and is given by Peak error is calculated by The Index of Agreement (IA) is formulated as where brackets, · , denote time averaging. I A = 1 shows perfect agreement and I A = 0 means complete disagreement. The Pearson correlation (Cor) coefficient is calculated by It has a value between +1 and −1, where 1 is total positive linear correlation, 0 is no linear correlation and −1 is total negative linear correlation.