The impact of ocean–atmosphere coupling and its possible seasonal dependence upon Weather Research and Forecasting (WRF) Model simulations of seven, wintertime cyclone events was investigated. Model simulations were identical aside from the degree of ocean model coupling (static SSTs, 1D mixed layer model, full-physics 3D ocean model). Both 1D and 3D ocean model coupling simulations show that SSTs following the passage of a nor’easter did tend to cool more strongly during the early season (October–December) and were more likely to warm late in the season (February–April). Model simulations produce SST differences of up to 1.14 K, but this change did not lead to significant changes in storm track (<100 km), maximum 10-m winds (<2 m s−1), or minimum sea level pressure (≤5 hPa). Simulated precipitation showed little sensitivity to model coupling, but all simulations did tend to overpredict precipitation extent (bias > 1) and have low-to-moderate threat scores (0.31–0.59). Analysis of the storm environment and the overall simulation failed to reveal any statistically significant differences in model error attributable to ocean–atmosphere coupling. Despite this result, ocean model coupling can reduce dynamical field error at a single level by up to 20%, and this was slightly greater (1%–2%) with 3D ocean model coupling as compared to 1D ocean model coupling. Thus, while 3D ocean model coupling tended to generally produce more realistic simulations, its impact would likely be more profound for longer-term simulations.
Computational power increases have driven the development of increasingly complex numerical weather prediction models, such as the Weather Research and Forecasting (WRF) Model (Skamarock et al. 2008). Despite the complexity and land surface coupling abilities of the WRF Model, it lacked any direct ocean–atmosphere coupling functionality until April 2009 (WRF version 3.2). This coupling is vital for tropical cyclone (Sutyrin and Khain 1984; Bender et al. 1993; Emanuel 1999; Warner et al. 2010; Olabarrieta et al. 2012) and coastal, midlatitude cyclone (Anthes et al. 1983; Kuo et al. 1991; Ren et al. 2004; Eckhardt and Stohl 2004; Knippertz and Wernli 2010) simulations because of those cyclones’ strong dependence upon ocean–atmosphere heat flux and momentum exchanges. Such dependencies motivate ocean–atmosphere coupled model development.
One early study (Bender et al. 1993) coupled the Geophysical Fluid Dynamics Laboratory tropical cyclone prediction model to a multilayer primitive equation ocean model to investigate wind stress–induced sea surface temperature (SST) cooling during Hurricane Gloria. Their simulations only varied cyclone propagation speed (slow, medium, fast) and produced maximum simulated SST coolings of 5.3°, 3.5°, and 1.8°C, respectively, and reduced ocean-to-atmosphere heat fluxes accordingly. These changes decreased simulated cyclone intensity and maximum winds by up to 7.3 hPa and 2.7 m s−1, respectively; storm track was not significantly altered.
More recently, Ren et al. (2004) coupled the Canadian Mesoscale Compressible Community atmospheric model to the full-physics, 3D Princeton Ocean Model to investigate ocean–atmosphere dynamics during the September 1998 extratropical transition of Hurricane Earl and an intense January 2000 cyclone dubbed the “Superbomb.” These systems were selected because of their comparable size and intensity, but the shallower mixed layer depth during Earl produced greater maximum SST reductions than for the Superbomb (5° versus 1°C). Despite this result, the slower propagation speed of the Superbomb led to greater intensity and maximum wind speed changes [4 hPa and 4 m s−1, respectively (1 hPa and 1 m s−1 more than Earl)]. Consistent with Bender et al. (1993), model coupling did not significantly alter storm track.
Using the Ren et al. (2004) coupled model, Yao et al. (2008) investigated 42 North Atlantic October midlatitude cyclones (locally dubbed “nor’easters”) off the U.S. East Coast. As compared to Ren et al. (2004), their simulations produced larger SST cooling (~6°C) and greater cyclone weakening (4–5 hPa), yet a similar decrease in 10-m winds (2–4 m s−1). Larger changes in Yao et al. (2008) were attributed to the thinner October ocean mixed layer depth and the wider case variety.
As of WRF version 3.6.1, three built-in SST options exist: static SST, SST update, or 1D ocean mixed layer coupling. For static SST, initial SST data are ingested from an external data source [e.g., 0.5° real-time global SSTs (RTGSST; Gemmill et al. 2007)] and are held constant. For SST update, SST values are updated at user-prescribed intervals using lower-boundary input files and are otherwise constant. Finally, the 1D ocean mixed layer model option couples WRF to the Pollard et al. (1972) 1D ocean mixed layer model. At initialization, standard ocean temperature profiles are affixed to SST values and a uniform, user-selected, mixed layer depth is prescribed over the entire model domain. During forward integration, atmospheric wind stress is applied to each simulated ocean column, but neither horizontal advection nor transport processes (e.g., Ekman transport) are simulated, which prevents accurate simulations of finescale SST gradients (e.g., the Gulf Stream) and the corresponding heat flux and momentum exchanges (Chelton et al. 2001; O’Neill et al. 2003; Tokinaga et al. 2005; LaCasse et al. 2008). Such limitations motivate the development of fully coupled modeling systems, such as the Coupled–Ocean–Atmosphere–Wave–Sediment Transport (COAWST) modeling system (Warner et al. 2010, hereafter W10).
The COAWST modeling system allows WRF to be directly coupled to the Regional Ocean Modeling System (ROMS; Shchepetkin and McWilliams 2005; Haidvogel et al. 2008), the Simulating Waves Nearshore (SWAN; Booij et al. 1999) wave model, and the Community Sediment Transport Modeling System (CSTMS; Warner et al. 2008) models via the Model Coupling Toolkit (MCT; Larson et al. 2004). As described by W10, COAWST is specifically designed “to better identify the significant processes affecting our coastlines and how those processes create coastal change.” Since its release, COAWST-related studies have focused primarily on two main themes: tropical cyclones (W10; Olabarrieta et al. 2012) and coastal process (rip currents, surf zone) simulations (Kumar et al. 2011; Olabarrieta et al. 2011). Applications of COAWST to baroclinically driven midlatitude cyclones (such as nor’easters) have yet to be addressed.
We define a nor’easter as a large (~2000 km), midlatitude cyclone occurring between October and April, bringing punishing winds, copious precipitation, and potential coastal flooding to the northeastern United States (Kocin and Uccellini 2004; Jacobs et al. 2005; Ashton et al. 2008). Total nor’easter-related damages can currently reach several billion U.S. dollars per event (NCDC 2015). With warmer temperatures and increased atmospheric water vapor content expected during the twenty-first century (Trenberth et al. 2007), the risk of damaging nor’easters is also expected to increase, primarily from sea level rise (Yin et al. 2009) and more intense precipitation (Lombardo et al. 2015). Since storm surge and precipitation crucially depend on air–sea interactions in the nor’easter environment, it is important to simulate these interactions accurately.
However, do we need a fully coupled model (such as COAWST) to capture these interactions, or is a less computationally demanding 1D ocean model sufficient? Furthermore, the ocean mixed layer is shallow in the fall, but deepens throughout the winter (Ren et al. 2004). Does this seasonality in the ocean’s thermal structure lead to a seasonal dependence in the impact of model coupling on nor’easter simulations? These are the questions this study aims to address.
a. COAWST overview
This study utilized COAWST revision 727 (May 2013 release), which offers four WRF-based configurations (one uncoupled, three coupled; see Fig. 1). Given our ocean–atmosphere focus, we used the WRF-only and WRF-ROMS configurations (Figs. 1a,b) and did not couple to either SWAN or CSTMS (Figs. 1c,d). The Spherical Coordinate Remapping Interpolation Package (SCRIP; Jones 1998) computed interpolation weights between the WRF (Fig. 2a) and ROMS (Fig. 2b) model grids needed to facilitate their coupling. COAWST did not allow for multidomain coupling when conducting this study, so ROMS domain 1 was only directly coupled to WRF domain 2 given their similar horizontal grid spacings (10 and 15 km, respectively) and spatial extent. Further WRF and ROMS model details are contained in the sections below.
b. WRF configuration
Version 3.4 of the Advanced Research WRF (hereafter W34) solves a set of fully compressible, nonhydrostatic, Eulerian equations in terrain-following coordinates (Skamarock et al. 2008). The three model domain grids have one-way feedback; 45-, 15-, and 5-km horizontal grid spacing, respectively; and 61 vertical levels, which afforded simulation of key precursor synoptic- and meso-α-scale phenomena (e.g., jet streaks, short- and long-wave troughs) on outer domains and smaller-scale phenomena (e.g., orographic forcing, condensational latent heating) on inner domains. Atmospheric lateral boundary conditions were derived from 1° × 1° Global Forecast System (GFS) model analysis (GMA), and the model top was set to 50 hPa. As in W10, RTGSST data were ingested for all three WRF domains, but for domain 2 these SST values were overwritten with ROMS-simulated SSTs within the ROMS coverage area (Fig. 2b). The three WRF domains had time steps of 90, 30, and 10 s, respectively.
shortwave radiation—new Goddard scheme (Chou and Suarez 1999),
land surface—Noah (Chen and Dudhia 2001),
cumulus parameterization—Kain–Fritsch (Kain 2004) (not applied to domain 3), and
microphysics—Goddard (Lang et al. 2007).
c. ROMS configuration
ROMS release 455 (Shchepetkin and McWilliams 2005; Haidvogel et al. 2008) numerically integrates 3D Reynolds-averaged Navier–Stokes equations using hydrostatic and Boussinesq approximations along vertically stretched, terrain-following (sigma) coordinates. The ROMS domain (Fig. 2b, shaded) has an average 10-km grid spacing (needed for finescale SST gradients), 16 sigma levels, and a 50-m sigma-coordinate surface–bottom-layer width. Time step size was fixed at 30 s to match WRF domain 2 and for consistency with W10.
Ocean initialization and open-boundary conditions required atmospheric, ocean, and tidal information. Atmospheric state variables are interpolated from WRF. Oceanic state variables (currents, temperature, salinity, sea surface height) originate from the Hybrid Coordinate Ocean Model (HyCOM; Bleck 2002). Tides were calculated from the Advanced Circulation (ADCIRC; Mukai et al. 2002) tidal dataset and included the five tidal constituents (M2, N2, S2, O1, and K1). All of these data sources were consistent with W10.
d. Cases for analysis and model run description
Research work focused on a small, but diverse sample of seven subjectively chosen nor’easter cases, which vary in both severity and time of year (Table 1). A modest sample size was needed to integrate computationally expensive, high-resolution COAWST simulations on available resources. To emphasize any potential seasonal dependencies, at least one case from every month (October–March) was selected, and the cases are sorted by month and day rather than by chronological order. Case severity was determined from the Northeast Snowfall Impact Scale (NESIS; Kocin and Uccellini 2004), which classifies nor’easters on a scale of 1 (notable) to 5 (extreme) based upon the population impacted, area affected, and snowfall severity. Early and late season storms (cases 1, 2, and 7) brought large impacts through heavy rain, strong winds, and coastal flooding, rather than snow, so they did not receive NESIS ratings.
Simulations are integrated for 180 h, starting roughly 72 h prior to the first precipitation impacts in the highly populated mid-Atlantic region. This starting point allowed sufficient model spinup time and for the establishment of both strong baroclinicity and latent heat flux off the U.S. East Coast. The high baroclinicity and heat fluxes are a by-product of the clash between the cold continental air and the expansive (>1000 km) northern edge of the Gulf Stream, which is an interaction crucial for accurate nor’easter simulations (Kuo et al. 1991; Mote et al. 1997; Yao et al. 2008). The first nor’easter-related precipitation is defined as the first 0.5-mm (~0.02 in.) precipitation reading from the New Jersey Weather and Climate Network (Robinson 2005). A smaller threshold is not used so that we avoid capturing isolated showers well ahead of the primary precipitation shield. A New Jersey–centric approach was chosen because of its high population density (461.6 km−2), significant contribution [$473 billion (U.S. dollars)] to the U.S. gross domestic product, and its central location in the mid-Atlantic (U.S. Census Bureau 2012).
To investigate the importance of ocean–atmosphere coupling, three model runs were completed for each case within the COAWST modeling system: a WRF-only run (CW34), a WRF run utilizing the Pollard 1D mixed layer model (CWPol), and a fully coupled WRF-ROMS simulation (CWR). All simulations were initialized identically, only varying by the ocean coupling option. Specifically, CW34 fixes SSTs, CWPol updates SST and mixed layer depth (200 m initially) every time step, and CWR exchanges ocean and atmosphere data between the component models every 10 min, as in W10.
e. Verification and analysis techniques
Atmospheric model validation compared WRF model output to both GMA and 4-km resolution stage IV precipitation data (Lin and Mitchell 2005). GMA data validated all model output (except precipitation) because of its extensive coverage, especially in data-sparse regions. Stage IV validated simulated precipitation because it is based on rain gauge and radar data, is gridded, and offers high resolution. Prior to any validations, all data were interpolated onto the coarsest grid spacing. Given our atmospheric focus, ROMS model data were not validated, but we may do so in future research.
Model output analysis is composed of several parts. First, skin temperatures (hereafter TSK) and their differences are analyzed; TSK is essentially SST at oceanic grid points. SST values are not compared because CWPol simulations do not update it and CWR simulations overwrite it with RTGSST data at lower-boundary update times. Storm track and intensity were derived from local minima in sea level pressure (SLP) via an objective, self-coded algorithm similar to that used at the Climate Prediction Center (Serreze 1995; Serreze et al. 1997). Using these tracks, the minimum SLP, storm maximum wind speed, and average TSK change were calculated and compared to GMA. Next, precipitation patterns and their distribution were evaluated against stage IV data and validated using bias and threat score calculations. Overall, the accuracy of the local-storm environment (i.e., within a 600-km-wide, storm-centered box that captures the storm and minimizes background contamination) and the large-scale environment (i.e., entire model domain and all times) was evaluated using the nonhydrostatic, moist, total energy norm (Kim and Jung 2009). Energy norm integrations were capped at ~100 hPa to limit large temperature errors near the model top:
In (1), Em is the moist total energy norm (J m2 kg−1); u′, υ′, and w′ are the zonal, meridional, and vertical wind perturbations (m s−1), respectively; p′ is the pressure perturbation (Pa); θ′ is the potential temperature perturbation (K); and q′ is the mixing ratio perturbation (kg kg−1). In addition, Nr, θr, ρr, Tr, and cs are the reference Brunt–Väisälä frequency (0.0124 s−1), reference potential temperature (270 K), reference air density (1.27 kg m−3), reference air temperature (270 K), and speed of sound (329.31 m s−1), respectively. Finally, cp is the specific heat at constant pressure (1005 J kg−1 K−1), g is the gravitational constant (9.8 m s−2), and ωq is a scaling factor (0.1).
Local storm environment energy norm calculations were split into two frameworks: model centered and GMA centered. The former centers the energy norm integration box at the model-indicated cyclone center location, whereas the latter integrates only at GMA-indicated nor’easter locations. For the large-scale environment analysis, the energy norm was calculated for each model domain and was complemented with root-mean-square error (RMSE) calculations of SLP, 850-hPa temperature, 500-hPa geopotential height, 300-hPa winds, and potential temperature on the dynamic tropopause (defined as the 2-PVU surface, where 1 PVU = 10−6 K m2 kg−1 s−1). The first four variables were referenced in Kocin and Uccellini (2004) for nor’easter analysis. Our decision to focus on the energy norm was influenced by Buizza et al. (2005), who made a compelling case for its usage at ECMWF for model validation given its total model volume integration, lack of single-layer sensitivity, and inclusion of temperature, wind, pressure, and moisture errors. For both RMSE and the energy norm, smaller values denote less error.
a. Simulated skin temperature comparison
Figure 3 displays TSK values from CW34 simulations at forecast hour 48 and the TSK differences for CWPol and CWR between forecast hour 48 and the end of the simulation (forecast hour 180) for cases 2, 4, and 7. These cases illustrate TSK variability during and following an early, mid-, and late-season nor’easter event and dramatically illustrate the greater TSK variability associated with 3D ocean–atmosphere coupling versus 1D coupling. Higher TSK differences in 3D ocean–atmosphere coupling are a product of shifts in the Gulf Stream and the loop current, wind-driven upwelling, and submerged warm water pockets as noted by Ren et al. (2004). To quantify how TSK varies in domain 2, we also show (Table 2) ocean TSK differences for all seven cases from GMA, CWPol, and CWR; differences for CW34 are exactly 0 K. Largely because the Gulf Stream’s location remains fixed (only temperature changes), average CWPol TSK differences in Table 2 vary by up to 1.18 K less than in CWR. Overall, the mean TSK differences vary up to 1.41 K (CWR, case 1). Table 2 and Fig. 3 also indicate a degree of seasonality in TSK differences. Early season cases (October–January) favor stronger overall cooling, whereas late-season cases (February–April) favor less cooling and even warming. These trends are not noted in previous nor’easter modeling studies because of their limited case or time focus and are associated with seasonal mixed layer depth [shallower mixed layer depth early in the season; Ren et al. (2004)] and solar heating variations. Notably, these trends are better handled by CWR than in CWPol because the former includes bathymetry rather than a set ocean depth (200 m), allows for advection between neighboring grid cells, and includes precipitation-based cooling.
b. Nor’easter track and property analysis
Figure 4 displays WRF- (colored lines) and GMA- (black lines) based storm tracks for each case. Figure 5 displays GMA-relative track errors (positive y axis = GMA cyclone propagation direction) every 6 h (small symbols, colors) and their mean (larger symbols, black) only when all simulations contained the nor’easter concurrently within domain 2. Average track error for all seven cases (km) is shown later (Table 4). Figure 5 shows a wide spread (up to 490 km) in average track error between cases, yet the track variability between experiments for each particular case remains fairly consistent (differences less than 100 km), which is consistent with previous studies (Ren et al. 2004; Yao et al. 2008; W10). Following from Fig. 5, Table 3 shows that no experiment exhibits a definitive (>60%) GMA-relative track bias over all 70 storm-track points (all cases). Even on a case-by-case basis, the track biases of the coupled runs only vary from CW34 at five (CWPol) and four (CWR) track points, and each is associated with SST variations along and near the Gulf Stream. Despite lacking definitive GMA-relative track biases, Fig. 5 does demonstrate that the accuracy and precision of storm tracks did vary between cases. For cases 1, 4, and 6, which had relatively low track error, each simulated nor’easter formed in a region of stronger differential cyclonic vorticity advection (CVA), whereas the other cases had higher track errors, and these nor’easters (especially case 7) underwent cyclogenesis in regions of less concentrated differential CVA.
In addition to average storm track error, Table 4 also displays various nor’easter properties including the lowest minimum SLP, highest sustained 10-m wind within 100 km of the cyclone center, and average TSK difference within 100 km of the cyclone center along the storm track for GMA and all three experiments. We compare each model to GMA because no in situ dataset exists for the data-sparse regions of the western Atlantic. With the exception of cases 2 and 5, simulated nor’easter intensity fell within 5 hPa of GMA. Notably for these two cases the nor’easter developed ahead of a 500-hPa cutoff where differential CVA was weak, as shown in Fig. 6, with 500-hPa geopotential height and relative vorticity fields from case 2. In GMA (not shown), weaker differential CVA values result in higher minimum SLP values. Relative to CW34, CWPol- and CWR-based SLP values are all within 5 hPa of each other with the exception of case 3. Further analysis of case 3 revealed a region of enhanced surface baroclinicity along the stronger SST gradient and a region of enhanced latent heat flux just west of the cyclone center and collocated with the SST gradient (Fig. 6). For CW34 or CWPol, both the baroclinicity and latent heat fluxes near the cyclone are weaker, which likely explains why CWR produced a stronger cyclone. Despite the enhanced SLP falls for CWR in case 3, for no case did the maximum 10-m wind between WRF simulations vary by more than 2 m s−1. Finally, similar to the SST differences averaged over domain 2, SST differences along each storm track also indicate a seasonal cycle, albeit generally more muted (lower magnitude) given the strong warm and cold SST pockets associated with the Gulf Stream.
c. Error analysis of the local-storm environment
Evaluation of the local-storm environment will focus exclusively on the energy norm (when storm center is at least 300 km from the domain 2 boundary) and will contain two parts: averaged energy norm total and its component parts at each model level (Fig. 7) and the vertical integrated energy norm at each relevant time (Fig. 8). Figure 7 shows the total energy norm (black) and its component parts (colors) for CW34 and then the percent change relative to CW34 represented by the difference between it and CWPol and CWR for cases 2–5. In Fig. 7, negative percentage values denote a decrease in the energy norm. Note that model error is dominated by the horizontal wind (U, V), temperature (T), and moisture (Q) components of the energy norm. These four cases were selected because they best illustrate situations where model coupling reduced the model error the most. These cases also cover four separate months. Among the missing cases, model coupling reduces error slightly in cases 1 and 7, but case 6 is slightly worse on average. In Fig. 7, model coupling is most beneficial to simulated water vapor (yellow line) and vertical velocity (brown line) values with decreases of up to 59% and 20% in the CW34 energy norm for a given level, respectively. Consistent with Ren et al. (2004) and Yao et al. (2008), ocean–atmosphere coupling leads to simulated differences throughout the entire troposphere.
Figure 8 shows the percentage change in the CW34 energy norm of the non-CW34 runs within both the GMA- (thin lines) and model-centered (thick lines) frameworks for all cases. In Fig. 8, a −5% value indicates a non-CW34 experiment has an energy norm 5% less than CW34, and data gaps represent times when the data integration box extended beyond the domain 2 boundary for one or more of the experiments. The CWPol and CWR energy norm differences range from 40% lower (CWPol, case 7) to 40% higher (CWPol, case 2) than the W34 energy norm and normally fall within 10% of the CW34 energy norm value. Taken as a whole, the average CWR energy norm is smaller than CW34 in all seven cases and CWPol in six of seven cases. Despite this success, on a time-by-time basis and combining all cases (64 total times), the bottom of Table 5 shows CWR only had model-relative energy norm values less than CW34 in 36 (56.3%) time periods, and less than CWPol in 40 (62.5%) time periods. In comparison, CWPol has energy norm values lower than CW34 in 32 (50.0%) time periods. Thus, CWR simulations do tend to give slightly better forecasts on average and in most situations, but two-tailed t tests did not reveal these energy norm differences as being statistically significant either overall or for any case (minimum p value 0.84; CWR–CW34 comparison). Perhaps errors in the GMA itself, which does not involve a fully coupled ocean, penalize the CWR enough in this metric to prevent statistical significance.
d. Stage IV precipitation analysis
WRF-simulated precipitation fields and their accuracy are strongly linked to storm track, forecast hour, and to a lesser degree ocean–atmosphere coupling. Figure 9 shows 72-h precipitation accumulations (forecast hours 48–120) from stage IV and CW34, CWPol, and CWR 72-h precipitation accumulations differenced from CW34, and the precipitation probability density function (PDF) and cumulative distribution function (CDF) of 72-h accumulated precipitation for stage IV and all models from cases 4 and 6. These two cases were selected for their low storm track error and relatively minimal track error spread. To supplement Fig. 9 and reference the other five cases, Table 6 shows the bias and threat scores for all seven cases assuming a 12.5-mm (0.5 in.) reference value. We exclude cases 2 and 3 from consideration because either the GMA or WRF model storm tracks lie outside the spatial coverage of the stage IV data, which results in extreme biases and low forecast skill.
Focusing on the remaining five cases, all WRF simulations tend to similarly overextend their precipitation extent from 1.17 (case 7, CW34) to 2.13 (case 4, CWPol) times that shown in the stage IV data. Threat score values show low (0.31; case 7, CWR) to moderate (0.59, case 5 all experiments) success in producing precipitation in the same locations as stage IV. Of the five cases, only in case 4 do both the bias and skill score improve slightly (relative to CW34) with model coupling. As seen in Fig. 9, this improvement is quite subtle and is visible as a small bump in the 10–15-mm bin on the CWR CDF curve and small reductions (10 mm) in precipitation in eastern Massachusetts; near Washington, D.C.; and east of the North Carolina coast. Figure 10 shows that these reductions are likely associated with changes within the lower troposphere given the similarity of the 500-hPa geopotential height and lower concentrations of differential CVA and other upper-level fields (i.e., 300-hPa winds and geopotential height). Instead, cooler SSTs along the mid-Atlantic coastline and a slight cooling (1–2 K) of the warm water tongue off North Carolina decreased the latent heat flux (50–100 W m−2) compared to nearby areas [see the asterisk (*) in Fig. 10], which increased local stability and, in turn, precipitation. The relatively minor improvements in QPF seen here are consistent with prior work that has shown QPF to be one of the most difficult forecast aspects for models to handle in a variety of situations (e.g., Fritsch and Carbone 2004; Wang and Clark 2010).
e. Error analysis over all WRF domains and times
Similar to Fig. 8, Fig. 11 displays the percentage change in the energy norm for the coupled experiments relative to CW34 (negative values denote reduced model error) for cases 4 and 6 across all three model domains. These two cases represent situations with high (case 4) and low (case 6) energy norm variability. Until cyclogenesis occurs, energy norm values tend to vary negligibly (<1%) because of a lack of strong surface wind stress and precipitation necessary to drive SST changes. Once generated, energy norm differences vary more rapidly (<5%) and sharply spike (up to 12%) when the nor’easter leaves the model domain because of varying fractions of the storm existing within the domain of interest until all simulations no longer contain the nor’easter within the model domain. As seen in Fig. 11, the percent change in the CW34 energy norm also tended to increase with high resolutions because the energy norm was summed over a smaller spatial extent so that small position errors in better-resolved updrafts and banded structures associated with the nor’easter contributed more significantly to the energy norm. Among the seven cases, 3D ocean–atmosphere coupling reduced model error in all model domains for at least four (domain 3) and up to six (all but case 6) out of seven cases. As compared to CWPol, CWR energy norm values were smaller in at least five (all but cases 5 and 6) of the seven cases. CWPol energy norm values also compared favorably (four or more cases with a lower energy norm) in all but domain 3, where CWPol energy norm values were often within 5% of CW34, as illustrated by Fig. 11, but more often than not were slightly higher than CW34.
To help quantify whether ocean–atmosphere coupling generally benefited the entire model simulation, Table 5 displays the total number and percentage of all 217 six-hour periods (31 six-hour periods per case, seven cases) where the shown inequality is valid for both the energy norm and five RMSE calculations. If only considering the energy norm, 3D ocean model coupling certainly produces more accurate nor’easter simulations than both CW34 and CWPol for domains 1 and 2 and arguably for domain 3 too. For the former, more than 62.7% of all 217 analysis time periods exhibit a lower energy norm from CWR simulations than either CW34 or CWPol. For the latter, this value drops to 51.6% of all 217 analysis times. This may be a function of model resolution and also a consequence of SST values not being passed from domain 2 to domain 3, where far coarser RTGSST data would be used. In comparison, CWR RMSEs show considerably weaker tendencies where regions above 500 hPa generally show a marginal to slight benefit (50.2%–60.8%) and a marginal to slightly negative (36.9%–48.4%) response for 850-hPa temperature and SLP. In comparison to CWR, 1D ocean model coupling produces energy norm and RMSE results that are comparatively less variable than CWR, but only two RMSE tests (domain 2, 300-hPa winds; domain 1, 850-hPa temperatures) show any distinct benefit (>55%) compared to the uncoupled simulation. As is, the RMSE results indicate that 1D ocean model coupling has greater benefits than 3D ocean model coupling in the lower atmosphere, yet in the mid- and upper atmosphere this situation reverses. However, overall model error (as indicated by the energy norm) indicates that 3D ocean model coupling produces better overall simulations than either 1D ocean model coupling or an uncoupled model. Despite the above results, two-tailed t tests did not reveal any significant variability among the energy norm results or any of the RMSE tests (minimum p value 0.73; CWR/CWPol domain 3 SLP).
An investigation into the benefits of ocean–atmosphere coupling and its potential seasonality was evaluated using seven sets of nor’easter simulations. Each set was composed of an uncoupled WRF (CW34), a coupled 1D ocean mixed layer model (CWPol), and a coupled 3D ocean model (CWR) experiment. All experiments were configured and initialized identically (except for the ocean coupling option) and ran for 180 h starting roughly 72 h prior to the first nor’easter-related precipitation measurements in the mid-Atlantic region. Model output was validated against stage IV and GMA data. Area-averaged (100-km radius), CWR-simulated skin temperature differences varied up to 1.41 K and were in all but one case (case 7) an order of magnitude larger than those from CWPol. CWR simulations (to a lesser extent CWPol) indicated a degree of seasonality in TSK differences where early season cases (October–January) favor stronger overall cooling, whereas late-season cases (February–April) favor less cooling and even warming. This seasonality is consistent with GMA, and also matches expectations based on the shallower mixed layer observed early in the cold season (Ren et al. 2004).
Although TSK differences could exceed 1 K, ocean model coupling did not produce significant intermodel changes to simulated storm track (<100 km), maximum 10-m wind speed (<0.79 m s−1), storm intensity (±5 hPa, exception CWR case 3, 5.25 hPa), or average, storm-track relative TSK (<1.27 K), or produce any distinctive track-bias tendencies. These results are consistent with those of Ren et al. (2004), Yao et al. (2008), and W10. Simulated precipitation fields showed little sensitivity to both 1D and 3D ocean model coupling, but bias and threat scores did show that WRF tended to overpredict the precipitation extent and have low (0.31) to moderate (0.59) model forecasting skill. Similar to previous modeling studies (i.e., Ridout et al. 2005; Dravitzki and McGregor 2011), most WRF simulations generated excess light precipitation and a dearth of heavy precipitation, with the notable exception of cases 2 and 5, where differential CVA aloft was weak.
Analysis of the local-storm environment and the overall simulation revealed the main sources of simulation error are associated with temperature and horizontal wind fields both within the planetary boundary layer and in the upper troposphere. With the exception of case 3, CWR simulations generated slightly lower energy norm values (20% reduction at a level) throughout the troposphere compared to corresponding CWPol simulations. Such reductions shed light upon how the error in CWR simulations decreased from CW34 in at least 34 out of 64 (56.3%) cases evaluated within both the model and GMA-relative frameworks. RMSE and energy norm analysis of the entire simulation demonstrated that ocean–atmosphere coupling does slightly to moderately (50.1%–71.9% of all analysis times) benefit from some degree of ocean model coupling. Benefits stemming from ocean–atmosphere coupling are not clear cut because neither CWR nor CWPol shows clear improvements in the RMSE calculations; however, energy norm results do clearly show CWR to have the best overall model simulation. Despite this claim, the overall impact to simulation outcome is modest given that simulation differences are shown to be nonsignificant (minimum p value 0.73; CWR/CWPol domain 3 SLP).
This study has investigated the merits of ocean–atmosphere coupling, compared 1D to 3D ocean modeling coupling, and determined whether any seasonal dependencies in SST differences exist. Although model coupling lacked any notable impact upon storm track, storm intensity, and storm maximum winds, and statistically significant RMSE and energy norm differences were not found, atmosphere–ocean coupling still does produce small, but visible differences in temperature, moisture, and precipitation, among other dynamical and physical fields. Despite these relatively meager differences, both 3D and to a limited extent 1D ocean model coupling did show a seasonal dependence where SSTs more strongly cooled and warmed depending upon the time of year. Because the SST differences from CWPol simulations tended to be an order of magnitude smaller than the CWR simulations, this does underscore the importance of advection and realistic bathymetry for realistic SST simulations. Even with nonsignificant error analysis results (both near the storm and throughout the domains), this study still demonstrates that ocean–atmosphere coupling has tangible benefits (albeit small) to nor’easter simulations by reducing model error. What remains less clear is the importance of 3D ocean modeling. This study has shown that CWR simulations did produce lower overall error, yet with a 3-day lead time, its simulations were not much improved upon the output from the 1D ocean model. Thus, it is likely that with longer model integration times the impacts from 3D ocean model coupling would become more apparent and crucial for model simulations.
5. Future work
This study investigated how ocean–atmosphere coupling impacts nor’easter simulations. Several expansions are possible. First, including the COAWST wave and sedimentation models (i.e., SWAN and CSTMS) would allow for further alterations in surface heat fluxes and wave-based variability in surface roughness values. These changes may lead to more realistic precipitation and cloud cover simulations. Second, these same simulations could be reexamined with the recently released (October 2014) COAWST, rev 828 or newer, which now includes multiple-nest coupling. Finally, this analysis could be expanded through comparing ROMS-simulated currents to Coastal Ocean Dynamics Applications Radar (CODAR; Paduan et al. 2004) and outgoing longwave radiation data to ex situ data from the Clouds and Earth’s Radiant Energy System (CERES; Wielicki et al. 1996).
This research was funded in part by the New Jersey Agricultural Experiment Station. We acknowledge both Dr. John Warner of the Woods Hole Oceanographic Institute and Dr. John Wilkin of the Institute of Marine and Coastal Sciences at Rutgers University for their invaluable and timely assistance. We also wish to thank the NASA Graduate Student Summer Program for the learning experience afforded by SDN’s internship at NASA Goddard Space Flight Center. We also acknowledge Dr. Wei-Kuo Tao of NASA Goddard Space Flight Center who generously supplied computational time for our research.
Current affiliation: NASA Postdoctoral Program Fellow, NASA Goddard Space Flight Center, Greenbelt, Maryland.