Various in-situ datasets are presented alongside model output, with the aim of assessing the model performance at multiple time-scales, starting from seasonal climatologies.

The seasonal variability of the water column structure across the continental shelf is taken from quasi monthly Conductivity-Temperature-Depth (CTD) casts along the St Helena Bay Monitoring Line (SHBML) from March 2000 to August 2017 (Oceans and Coastal Research, 2017). The transect is comprised of twelve stations orientated approximately perpendicular to the isobaths off St Helena Bay ( Figure 1 ). A detailed description of the data is provided in Lamont et al. (2015). Measurements from all available cruises were used to calculate monthly climatologies of temperature and salinity at each station on a 1 m vertical grid.

The inter-annual and event-scale upwelling/relaxation dynamics in the nearshore region of St Helena Bay is taken from daily averaged temperature observations from a fixed mooring in Eland’s Bay near the surface (∼5 m water depth) and at the seabed (∼20 m water depth) over the three year period from November 2008 to November 2011 (‘EB’ in Figure 1 ). These data are described in more detail in Pitcher et al. (2014).

The super-diurnal variability of the system is informed by concurrent deployments of Wirewalker wave-powered profilers (Rainville and Pinkel, 2001; Pinkel et al., 2011) and bottom-mounted Acoustic Doppler Current Profilers (ADCP) providing horizontal velocity components and temperature at a temporal resolution of 10 min and vertical resolutions of 1 m and 0.25 m, respectively (‘WW’ in Figure 1 ). The full dataset is described in detail in Lucas et al. (2014), although only data from the month of March 2011 are revisited in this study for model validation purposes. This period was specifically selected so that the model output could be compared with the reduced physics experiments of Fearon et al. (2020; 2022). The data were filtered in time to provide a two hour running mean at 30 min intervals, sufficient for analysing the super-diurnal processes of interest.

This study makes use of the Group for High Resolution Sea Surface Temperature (GHRSST) dataset for the Eastern Atlantic Region from the Spinning Enhanced Visible and InfraRed Imager (SEVIRI, https://podaac.jpl.nasa.gov/dataset/SEVIRI_SST-OSISAF-L3C-v1.0) (OSISAF, 2015). The geostationary orbit of the satellite allows data to be made available at hourly intervals from June 2004 to present, on a 0.05° regular grid. The high temporal resolution of the data is used to assess day vs night-time SST’s over the area of interest. Only data with a quality flag of 3 and greater are used in the presented analysis, as recommended in the scientific validation report for the product (SAF, 2018).

The ocean model employed in this study is the V1.0 official release of the Coastal and Regional Ocean COmmunity model (CROCO, http://www.croco-ocean.org/), an ocean modelling system built upon ROMS AGRIF (Shchepetkin and McWilliams, 2005). CROCO is a free-surface, terrain-following coordinate oceanic model which solves the Navier-Stokes primitive equations by following the Boussinesq and hydrostatic approximations. The model solves equations governing the conservation of horizontal momentum, hydrostatic balance, incompressibility and the conservation of tracers (temperature and salinity). A curvilinear Arakawa C-grid is used for the discretisation of the horizontal plane, while the vertical grid is discretised using a terrain-following (σ) coordinate reference system.

The choice of the vertical turbulent viscosity scheme is particularly important within the context of this study, as we attempt to highlight diapycnal mixing driven by super-diurnal wind variability. Initial simulations indicated that CROCO’s default K-profile parameterization (LMD-KPP) scheme resulted in an over-estimation of vertical mixing within the nearshore regions of St Helena Bay when compared with the observations. The results presented in this study therefore adopt the k-ϵ turbulent closure scheme within the Generic Length-Scale (GLS) formulation (Umlauf and Burchard, 2003; 2005), which provided a better fit to the observations. Horizontal eddy viscosity and diffusivity are typically ignored in regional configurations using ROMS (e.g. Marchesiello et al., 2003; Veitch et al., 2010). Instead, horizontal dissipation at the grid scale is determined by dissipation associated with a split and rotated third-order upstream biased horizontal advection scheme (Marchesiello et al., 2009). Bottom friction is parametrised using a quadratic drag law, where the bottom roughness length parameter is taken as 0.001 m. A nonlinear equation of state is used for the computation of density (Jackett and Mcdougall, 1995). The model configuration, as described below, effectively constitutes a dynamical downscaling of a global reanalysis product to ∼1 km resolution over St Helena Bay, so that the high-frequency bay-scale processes can be better assessed.

Each of the three types of wind forcing ( Table 1 ) are used to integrate the model over the period 01 November 2005 to 31 December 2012 (7 years), spanning the duration over which all forcing data are available. Three hourly ‘snapshots’ of temperature generated over the entire simulation period are used to assess day- vs night-time sea surface temperature differences in the model. We use model output times of 14:00 and 02:00 (local time, UTC+2) for this purpose, and limit the analysis to the nominally upwelling months of October to March, which also coincides with enhanced land-sea breeze forcing in the region of interest (Fearon et al., 2020).

Daily averaged output from the 7 year simulations are used to compute climatologies of temperature and horizontal currents, allowing for an assessment of the impact of the land-sea breeze on the mean state of the vertical water column structure and circulation over the area of interest. Again, this analysis is confined to the nominally upwelling months of October to March. Upwelling and relaxation wind regimes are considered separately, as different sub-diurnal processes are at play in each. To distinguish between these regimes we consider the time-series of the 3 day average alongshore wind stress for a location in the centre of St Helena Bay, as shown in Figure 2 . We take active upwelling as times where the 75 ^th percentile 3 day average alongshore wind stress is exceeded (i.e. the black time-series above the red line), while relaxation periods are identified as times where the 3 day average alongshore wind stress falls below the 25 ^th percentile value (i.e. the black time-series below the blue line). Composite means are then presented based on these wind regimes.

Given that one of the main implications of the land-sea breeze in EBUS is its effect on vertical mixing, we present an experiment specifically designed to examine the event-scale mixing processes in the model over the weekly time-scale associated with upwelling events. To do this, we compare hourly output from two simulations, one forced by τ d a i l y and one forced by τ d a i l y + a c , each initialised on 25 January 2007 from the simulation forced with τ d a i l y (this ensures that the influence of the land-sea breeze is absent from the initial condition in both simulations). The initialisation date was subjectively chosen to be at the end of a relaxation event during which a shallow stratified surface layer developed within St Helena Bay, as shown in Figures 3A, B . The vertical isotherms at the offshore extent of this shallow surface layer (60 to 80 km offshore) are associated with the upwelling cell off Cape Columbine, and associated equatorward flowing Columbine Jet, which represents a dynamic offshore boundary of the bay. The model is integrated over a 7 day period which is characterised by a local amplification of land-sea breeze, as indicated by the amplitude of the diurnal anticlockwise rotary component of the winds ( τ a c 0 ) ( Figure 3C ). It is also a period of upwelling-favourable winds, as indicated by the mean wind stress ( Figure 3D ). Note that we use τ a c 0 to refer to the amplitude of the rotating wind stress component, while τ a c denotes the rotating wind stress vector. Figure 3E provides the time-series of wind stress components for τ d a i l y and τ d a i l y + a c for a location in the centre of St Helena Bay.

The analysis of diapycnal mixing in the model is aided by initialising the model with a subsurface passive tracer ( C , in arbitrary tracer units per volume, ATU m^-3). C is initialised to a value of 1 for temperatures of 12°C and less, 0 for temperatures of 13°C and above, while linear interpolation is applied for temperatures between 12°C and 13°C. This approach is also adopted for specifying the passive tracer at the lateral boundaries of the model. The cumulative diapycnal mixing in the model is quantified by integrating the passive tracer from the depth of the 12°C isotherm to the surface:

It is informative to isolate the effect of the land-sea breeze forcing on C s , which is done by computing the difference in C s for simulations which both include and exclude the diurnal anticlockwise rotary component of the wind stress ( τ a c ):

Δ C s therefore represents a diagnostic aimed at highlighting the increase in cumulative diapycnal mixing due to the inclusion of land-sea breeze forcing.

We further use the vertical displacement of the 12°C isotherm to diagnose internal wave generation and propagation in the model due to the inclusion of land-sea breeze forcing. To this end, the depth of the 12C isotherm is presented as the displacement from the daily running average, allowing us to isolate the super-diurnal variability induced by the applied land-sea breeze forcing.

The purpose of the model validation is to assess the skill of the model at representing the variability of the system at multiple time-scales, including seasonal fluctuations ( Figure 4 ), upwelling/relaxation events ( Figure 5 ), and higher frequency process associated with the land-sea breeze ( Figure 6 ). Preference is given to in-situ data over satellite observations for the model validation owing to the nearshore and high-frequency nature of the processes of interest.

Figure 4 presents a comparison of modelled and observed seasonal temperature climatologies along the St Helena Bay Monitoring Line (‘SHBML’ in Figure 1 ). The vertical water column characteristics for the offshore region (Station 7 and offshore) is shown to be distinct from those in the nearshore waters of St Helena Bay (Stations 1 to 4). During austral summer (December to February), Stations 3 and 4 (in the middle of St Helena Bay) are characterised by a highly stratified shallow surface layer overlaying cool subsurface waters. In contrast, winter months at these stations are characterised by a nearly fully mixed water column. Nearshore temperatures are on average cooler than those offshore, particularly during summer, in line with the seasonal variability in upwelling-favourable winds. The climatological model bias within St Helena Bay (Stations 1 to 4) is below 1°C for both summer and winter seasons, indicating that the seasonal variability in the mean state of water column over the area of interest is well represented in the model.

Figure 5 presents a 3 year time-series of measured and modelled daily averaged temperature at the Elands Bay fixed mooring (‘EB’ in Figure 1 ), while summary statistics for the nominally winter (April-September) and summer (October-March) periods are provided in the Supplementary Table S1 . The data indicate a well-mixed water column during winter with relatively little variability, while summer reflects the variability associated with upwelling/relaxation events. Extended periods of upwelling result in a drop in temperature in both surface and subsurface waters, while relaxation of equatorward winds drive stratification events through the warming of surface waters which are advected inshore and poleward (Fawcett et al., 2007; Fawcett et al., 2008; Pitcher et al., 2014). Relaxation events can lead to the subduction of warm water to the bottom, as reflected in the time-series at 20 m depth. The model is shown to under-predict surface warming for some relaxation events (top panel of Figure 5 ), suggesting that the surface heat flux and/or the inshore and poleward advection of heat is not adequately captured for these events. This results in a cool bias of ~1.2C in the surface waters during summer. Surface model temperatures achieve correlation coefficients of 0.8 and 0.92 for the summer and winter periods, respectively. Subsurface bias is low throughout the year (<0.6C) with correlation coefficients of 0.88 and 0.89 for the summer and winter periods, respectively. Overall, these results indicate that the model shows acceptable skill at capturing the salient features of the upwelling/relaxation dynamics in the nearshore waters of St Helena Bay.

For reference, Supplementary Table S2 provides summary statistics at the Elands Bay mooring, but for the simulation forced with daily averaged wind stress ( τ d a i l y ). The statistics confirm that the model performance is significantly reduced when the super-diurnal winds are excluded in the model forcing. For example, the correlation coefficient for surface temperature reduces from 0.8 to 0.68 during summer, and from 0.92 to 0.83 during winter, highlighting the importance of super-diurnal winds in representing upwelling-relaxation dynamics in the model.

We now further analyse the extent to which the observed super-diurnal variability, embedded within the sub-diurnal upwelling dynamics, is captured by the model forced with τ r e a l . Figure 6 presents the temporal evolution of temperature and horizontal velocity components through the water column for a 7 day period at the Wirewalker mooring in 60 m water depth (‘WW’ in Figure 1 ), both in the observations and in hourly output from the model. For reference, this is the same 7 day period considered in the reduced physics experiments of Fearon et al. (2020; 2022). A time-series of wind stress components over this time period at the Wirewalker location can be found in Figure 2 of Fearon et al. (2022). The considered period begins with a relaxation event characterised by a highly stratified two layer system. The onset of an upwelling event is shown to be accompanied by strong variability at the diurnal-inertial frequency in both the temperature and horizontal velocity profiles, with a 180° phase shift between the surface and subsurface currents. The model performance over the full month of March 2011 is summarised in the Supplementary Figure S1 . Modelled temperatures are shown to be under-estimated near the surface, although here temperature correlation is highest (∼0.8 in the surface layer of the model), reflecting general agreement in the super-diurnal variability in surface temperature. Correlation coefficients of 0.4-0.6 and 0.6-0.7 are attained for near-bottom and near-surface current velocities, respectively. Overall, the comparison indicates agreement between the modelled and observed super-diurnal variability, and provides some confidence that the model can be used to gain further insight into the spatial and temporal variability of the processes of interest.

The shaded time in Figure 6 denotes a 24 hour period characterised by the onset of strong diurnal variability in surface temperature, both in the observations and in the model. The evolution of the spatial variability in sea surface temperature (SST) and surface currents over St Helena Bay for this 24 hour period is further assessed in Figure 7 . It is worth noting that the shown period displays some typical features of the region, as documented in the literature. For instance, upwelling along the coastline south of St Helena Bay sets up a shelf edge front and the associated Goodhope Jet, a strong equatorward current in geostrophic balance with the cross-shore pressure gradient (Strub et al., 1998; Veitch et al., 2018). North of Cape Columbine, the Goodhope Jet tends to bifurcate due to shelf widening, resulting in offshore and alongshore components (Shannon and Nelson, 1996; Veitch et al., 2018), the latter corresponding to the Columbine Jet which effectively defines a dynamic boundary between the nearshore waters of St Helena Bay and those offshore (Lamont et al., 2015). Also evident in the model is upwelling along a narrow belt on the eastern periphery of the bay and the formation of an inner shelf upwelling front (Jury, 1985; Taunton-Clark, 1985). Mean surface flow on the inner shelf during active upwelling is typically equatorward, in agreement with observations (Fawcett et al., 2008).

Figure 7 highlights the super-diurnal variability in surface currents and temperature embedded within these features, particularly within St Helena Bay. Here, the oscillatory nature of the surface currents is clearly evident, whereby the cross-shore component in the surface current is offshore at night, turning to onshore during the day. This phasing is tied to that of the land-sea breeze, indicating they are driven primarily by the diurnal-inertial resonance phenomenon. The surface currents drive the outcropping of the inner shelf upwelling front within St Helena Bay during the night (e.g. refer to the 08:30 time-step), however six hours later the newly upwelled waters are once again subducted by onshore surface currents. These processes contribute to temperature variability of greater than 5°C within the inner shelf regions of St Helena Bay over the considered 24 hour period.

Given this stark difference between day- and night-time surface temperature for this particular 24 hour period, it is prudent to assess the persistence of this effect. This is achieved by computing day- and night-time climatologies of SST, considering only the nominally upwelling months of October to March. Figure 8 presents this analysis for the GHRSST-SEVIRI satellite data over the 7 year duration of the model integration, while that of the realistically forced model is presented in Figure 9 . The satellite data confirm that indeed the observed day-time SST’s are on average warmer than those at night, and that the effect increases towards the coast, where the difference exceeds 1°C ( Figure 8C ). A similar pattern is seen in the model, although the effect is more exaggerated, particularly in the nearshore region of St Helena Bay where the mean day-night difference in SST exceeds 2°C over much of the bay ( Figure 9C ). SST’s are also notably higher in the satellite data than in the model ( Figures 8A, B vs Figures 9A, B ). It should however be noted that flagging techniques used to ‘de-cloud’ thermal infrared satellite data are known to contribute to large warm biases in upwelling systems by attributing cold upwelled water as missing data (Dufois et al., 2012; Smit et al., 2013; Meneghesso et al., 2020; Carr et al., 2021). Our high frequency model output is therefore seen as a preferred reference for assessing diurnal variability in SST’s over the area of interest.

Although a strong day-night signal is evident in Figure 9C , it isn’t clear how much of this signal is due to the movement of the front by the diurnal-inertial oscillations versus surface heat fluxes driven by solar irradiance (the land-sea breeze is itself a product of diurnal heating/cooling over the land and ocean). To distil these two effects we make use of the 7 year simulation forced with τ^daily ( Table 1 ), which maintains hourly heat fluxes but removes the diurnal-inertial oscillations. The mean day-night difference in SST for this simulation is shown in Figure 9D , while the effect of the diurnal-inertial oscillations is estimated by subtracting this from the overall signal ( Figure 9E ). The influence of the land-sea breeze on the diurnal fluctuations in SST is shown to be particularly enhanced in the nearshore region of St Helena Bay, where the mean difference exceeds 1°C. Diurnal warming due to solar irradiance is also enhanced in the southern extent of the bay, leading to a pronounced combined effect.

One of the main implications of the high amplitude diurnal-inertial oscillations is their impact on vertical mixing. We now present the results from an experiment specifically designed to examine these processes over an upwelling event accompanied by an enhanced land-sea breeze (refer to Section 2.5 for the experiment setup). Figure 10 shows the temporal evolution of the cross-shore component of velocity and the passive tracer at a location in the centre of St Helena Bay (‘P’ in Figure 3 ). The comparison of the simulations forced with both τ d a i l y + a c (left panels) and τ d a i l y (right panels) clearly demonstrates how the inclusion of τ a c leads to shear-driven deepening of the thermocline primarily in response to the diurnal-inertial resonance phenomenon and the associated injection of subsurface waters into the surface layer, all of which are embedded in the sub-diurnal upwelling dynamics.

Figure 10 also indicates the depth of the 12°C isotherm, which reveals a complex pattern of thermocline displacements when land-sea breeze forcing is included in the model. These thermocline displacements can be interpreted from the Hovmöller diagram shown in Figure 11A . The model reveals how the super-diurnal onshore and offshore advection of the surface layer generates large thermocline displacements at both the landward and offshore extents of the bay. Offshore surface flow is associated with upward thermocline displacements at the land boundary and downward thermocline displacements at the offshore extent of the bay, and vice versa for onshore surface flow. These forced diurnal vertical motions trigger offshore and onshore propagating internal waves which are evanescent, given that the latitude is poleward of 30°S (Zhang et al., 2010). In addition to the expected thermocline pumping and internal wave generation at the land boundary (Millot and Crepon, 1981; Fearon et al., 2022), the results point to an additional offshore source of internal wave generation associated with the Cape Columbine upwelling plume, where thermocline pumping is in fact larger than at the land boundary. This is attributed largely to the spatial heterogeneity in the depth of the Ekman layer which changes abruptly at the upwelling plume (see Figure 3B ) and introduces large spatial variability in the forced diurnal-inertial surface currents. Enhanced vorticity associated with the alongshore jet would further serve to introduce spatial variability in the diurnal-inertial motions in this region (Weller, 1982; Klein et al., 2004; Elipot et al., 2010). The spatial variability in the forced diurnal-inertial oscillations drive convergences and divergences in the surface flow which are responsible for the thermocline pumping and internal waves seen in the model.

Figure 11B illustrates the development of enhanced diapycnal mixing across St Helena Bay due to the inclusion of super-diurnal winds (refer to Equation 2 for how this diagnostic is computed). We note that the relative contribution of the locally forced oscillations and the evanescent internal waves on the cumulative diapycnal mixing can be difficult to separate, although we would expect mixing to be largely driven by the locally forced oscillations, while the internal waves are expected to play a secondary role (Fearon et al., 2022). The region of maximum enhanced diapycnal mixing is shown to be associated with the shallow stratified surface layer within St Helena Bay. Within this region, the cumulative diapycnal mixing varies considerably at time-scales of less than a day, in response to the onshore and offshore advection of the surface layer, and the large thermocline displacements. Regions of the largest enhancement of cumulative diapycnal mixing are typically associated with large downward thermocline displacements, given that C s represents the integration of C from the 12°C isotherm to the surface. This super-diurnal variability is embedded within the sub-diurnal offshore advection of the region of maximum enhanced mixing associated with the surface layer. The spatial variability in enhanced diapycnal mixing is further explored in Figure 12 , which presents the increase in C s due to the inclusion of land-sea breeze forcing, averaged over the fifth day of the simulation (29-30 January 2007). The analysis clearly identifies a region of particularly enhanced diapycnal mixing relative to the surrounding areas, associated with the shallow surface layer which had formed within St Helena Bay ( Figure 3 ).

We now consider how the presented super-diurnal processes may affect the mean state of the vertical water column structure and circulation over the area of interest. This is achieved by comparing 7 year climatologies for simulations forced with τ d a i l y and τ d a i l y + a c , thereby isolating the effect of the diurnal anticlockwise rotary component of the wind stress ( τ a c ) ( Table 1 ). Climatological comparisons are provided in Figures 13 , 14 for upwelling and relaxation wind regimes, respectively (refer to Section 2.4 for how upwelling and relaxation wind regimes are defined for this analysis).

The inclusion of τ a c during upwelling conditions is shown to drive a notable deepening of the thermocline within St Helena Bay ( Figure 13B ), consistent with the effects of enhanced vertical mixing driven by the diurnal-inertial oscillations. For example, at an offshore distance of 30 km, the mean 12°C isotherm depth is ~14 m for the simulation forced with τ d a i l y , increasing to a depth of ~21 m for the simulation forced with τ d a i l y + a c . The net warming of the entire water column in the inner shelf region (right panel of Figure 13B ) however points to effects which can’t be explained purely by the vertical exchange of surface and subsurface waters. Fearon et al. (2022) used a 2D model to show how a deepened thermocline leads to a reduction in the offshore advection of the surface layer, thereby maintaining the upwelling front closer to the coast and driving a net warming of nearshore waters. This is invoked as an explanation for the nearshore warming seen in the 3D model during upwelling conditions.

Figure 13C reveals mean alongshore currents which are largely geostrophically driven, considering the cross-shore structure of the isotherms ( Figure 13B ). The strong equatorward current offshore of St Helena Bay (mean surface current >0.4 m s^-1 at an offshore distance of ~55 km) corresponds to the Columbine Jet, in geostrophic balance with the cross-shore pressure gradient set up by the Columbine upwelling plume. The model also indicates an equatorward mean flow associated with the inner shelf upwelling front (<~15 km offshore), in agreement with observations (Fawcett et al., 2008), while a poleward mean flow is associated with the offshore extent of the shallow surface layer inside St Helena Bay (refer to depths of ~10-30 m and offshore distances of ~20-35 km). Mean near-bottom currents within St Helena Bay are shown to be poleward, which may in part be linked to barotropic continental shelf waves (Holden, 1987), cyclonic barotropic flow in the lee of Cape Columbine (Penven et al., 2000), and negative wind stress curl through the Sverdrup relation (Veitch et al., 2010). The inclusion of land-sea breeze forcing ( τ a c ) is shown to amplify the geostrophically driven alongshore currents associated with either end of the shallow surface layer within St Helena Bay. For example, the mean equatorward currents associated with the inner shelf upwelling front are roughly doubled in the area of maximum impact, increasing from 0.09 m s^-1 to 0.17 m s^-1 (right panel of Figure 13C ). The mean poleward flow at the offshore end of the surface layer increases from 0.02 m s^-1 to 0.05 m s^-1 in the area of maximum impact. These impacts are consistent with a locally deepened thermocline, which serves to steepen cross-shore pressure gradients and therefore alongshore geostrophic flow.

During periods of wind relaxation, the inclusion of land-sea breeze forcing ( τ a c ) is shown to drive cooler SST’s over most of the model domain (right panel of Figure 14A ). This effect is particularly evident within the nearshore regions of St Helena Bay, where the mean difference between the two simulations exceeds 2°C. Inspection of the vertical water column structure ( Figure 14B ) reveals that the lower SST’s driven by τ a c are compensated by warmer subsurface waters due to a deeper thermocline, consistent with the effects of enhanced diapycnal mixing by the diurnal-inertial oscillations. It is interesting to note that the net warming of the nearshore water column in response to τ a c during upwelling periods ( Figure 13B ) is not at all present during periods of relaxation. This further supports the suggested mechanism for the net warming to be linked to reduced offshore advection of the surface layer, which would be largely absent during relaxation of equatorward winds.

Periods of wind relaxation are associated with a weakened Columbine Jet, and a poleward mean flow within St Helena Bay, again consistent with inner shelf observations (Fawcett et al., 2008). During these periods, the only notable impact of τ a c on the mean flow is again on the inner shelf, where a geostrophically driven poleward flow is associated with the subduction of the thermocline. In this case, the deepened thermocline due to the inclusion of τ a c serves to reduce the cross-shore pressure gradients and weakens this poleward flow.