1 Introduction

The Eastern Boundary Upwelling Systems (EBUSs) are among the major fishery grounds in the world, contributing more than 20% to global fish catches, while covering only 2% of the global ocean surface (Pauly and Christesen 1995). Coastal upwelling systems are generated by the action of the equatorward along-shore winds, which transport coastal surface water offshore, causing the upwelling of the cold and nutrient-rich deep waters by Ekman dynamics. These systems have a large ecological and economic importance (García-Reyes et al. 2015; Varela et al. 2015), hence knowing the response of coastal upwelling to changing climate is of uttermost importance. There are four EBUSs: Canary Current Upwelling System (CCUS), Benguela Current Upwelling System, Humboldt Current Upwelling System and California Current Upwelling System.

In the EBUSs, the large-scale pressure gradients between are the main drivers of upwelling favourable winds (Bakun 1990; García-Reyes et al. 2015). Bakun (1990) hypothesized that the strengthening of the ocean-land thermal gradient under the greenhouse warming would result in stronger alongshore winds, increasing the upwelling of the deeper water to the surface during the second half of the twentieth century. This hypothesis was supported by Wang et al. (2015), who found a robust relationship between the increase of the land-sea temperature differences and the upwelling intensity in the twenty-first century.

However, Ryckazewski et al. (2015) proposed an alternative hypothesis: changes in the magnitude, timing or location of upwelling winds could be associated with poleward migration and intensification of major atmospheric high-pressure cells in response to increased greenhouse gas concentrations. The upwelling systems may be more sensitive to this mechanism than to the increase of the ocean-land thermal gradient (García-Reyes et al. 2015).

Nevertheless, along-shore winds are not the sole driver of change in the EBUSs. Coastal surface warming may increase the water stratification, reducing the upward nutrient-rich transport to surface (Di Lorenzo et al. 2005; García-Reyes et al. 2015; Brady et al. 2019). For instance, Oyarzún and Brierley (2018) showed that the increase of the ocean stratification, as result of global warming, becomes the major mechanism of change during the twenty-first century in the Humboldt Current system. Therefore, the EBUS future behaviour is still uncertain and both stratification and wind changes may be complementary or competitive (Bonino et al. 2019).

In this study we focus on the CCUS, which is part of the North Atlantic subtropical gyre, extending from the northern tip of the Iberian Peninsula at 43° N to the south of Senegal at about 10° N (Fig. 1a). The coastal upwelling region between Cape Blanc and the Strait of Gibraltar is a year-round phenomenon that becomes stronger to the south of Cape Bojador (Cropper et al. 2014). In the western coast of the Iberian Peninsula (37º–43º N) there is a marked seasonality with along-shore upwelling favourable winds dominating during spring and summer months (Sousa et al. 2017a; Wooster et al. 1976). The Strait of Gibraltar divides both regions, representing a major discontinuity of the upwelling system due to an abrupt change in the coastline orientation in the Gulf of Cádiz (Kämpf and Chapman 2016).

Fig. 1
figure 1

a ROM coupling scheme, with sea surface temperature (SST) and 2 m air temperature (T2m). The masks used in SST section (red lines) and to calculate the T2mland-sea differences (blue lines) are overlaid. b MPIOM grid resolution and REMO domain for ROM (red boundary)

The resolution of the Global Climate Models (GCMs) is insufficient to reproduce mesoscale features or detailed coastal dynamics (Xiu et al. 2018). In fact, recent studies (García-Reyes et al. 2015; Wang et al. 2015; Sein et al. 2017; Gómez-Letona et al. 2017; Bindoff et al. 2019) stressed the need for much higher horizontal resolution for the explicit representation of mesoscale processes that are important for a correct simulation of the upwelling systems but are partly masked in the current GCMs. Atmosphere–ocean regional climate models (AORCMs) are able to provide a resolution that allows an explicit representation of these processes and represent a valuable tool for their study.

In this work we use a high-resolution atmosphere–ocean regional coupled model for the study of the CCUS. The main goals of the paper can be summarized as follows:

  • First, validating the representation of the inter-annual variability and the seasonal cycle of the CCUS.

  • Second, comparing the outputs of regionally coupled model with regional climate models (RCMs) and GCMs to highlight the advantages and disadvantages of the AORCMs in reproducing the CCUS dynamics.

2 Data and methods

ROM (REMO-OASIS-MPIOM; Mikolajewicz et al. 2005; Sein et al. 2015, 2020) is an atmosphere–ocean regional climate model in which the limited area atmospheric model REMO (REgional MOdel; Jacob 2001) is coupled to the global ocean model MPIOM (Max Plank Institute-Ocean Model; Marsland et al. 2003; Jungclaus et al. 2013). These components are coupled via the OASIS3.0 (Valcke 2013) coupler. ROM also simulates globally the lateral freshwater fluxes at the land surface through the Hydrological Discharge (HD, Hagemann and Dumenil-Gates 1998, 2001) model, which is run as part of REMO. Furthermore, the relevant carbon stocks of the atmosphere, the ocean and the sediments are included through the Hamburg Ocean Carbon Cycle (HAMOCC; Maier-Reimer et al. 2005), which is a MPIOM subsystem (Soares et al. 2018, 2019; Parras-Berrocal et al. 2020).

2.1 ROM configuration

The atmospheric component of ROM, REMO, is the only model component that is run in a regional configuration. The REMO domain (Fig. 1b) extends to the North Atlantic, the eastern tropical Pacific and the Mediterranean Sea regions, including the whole CCUS region, and has a constant horizontal resolution of 25 km with a rotated grid (Parras-Berrocal et al. 2020).

MPIOM is discretized on a curvilinear grid with variable spatial resolution ranging from 5 km at the West African coast to 100 km in the Southern oceans. In the area of study, the MPIOM resolution is not coarser than 10 km. The MPIOM horizontal resolution in CCUS is sufficiently high to allow for the study of local scale processes in the region while maintaining a global domain with an acceptable computing cost. A similar MPIOM configuration was successfully applied in a process study of the Mediterranean Outflow Waters (Izquierdo and Mikolajewicz 2019). The model has 40 vertical levels with increasing level thickness towards the ocean bottom (Sein et al. 2015; Parras-Berrocal et al. 2020). The spin-up of MPIOM was done according to the procedure described in Sein et al. (2015). MPIOM is started with climatological temperature and salinity data (Levitus et al. 1998). Subsequently, it is integrated six times through the 1958–2002 period forced by ERA-40 and one time by ERA-Interim reanalysis (1979–2012) and with 60 min of the coupling frequency.

In this work, we assess the ROM performance in the CCUS using a simulation forced by ERA-Interim (Dee et al. 2011). ERA-Interim provides lateral boundary conditions to REMO and surface forcing to MPIOM outside the coupling region. Both models are hydrostatic and solve the Navier–Stokes equations using the Boussinesq approximation. In this run HAMOCC was disabled.

2.2 Assessment strategy and data sets

The model outputs were assessed in terms of seasonal cycle and inter-annual variability, focusing on the larger scale processes, the latitudinal variability and the mesoscale dynamics that characterizes the CCUS. For the validation, we use daily and monthly data of the most representative atmosphere and ocean variables: near-surface air temperature (T2m), sea surface temperature (SST), wind stress and ocean temperature. This validation is carried out against several observational and reanalysis data sets (Table 1). Additionally, we use the Max Plank Institute—Earth System Model (MPI-ESM), seven RCMs from AFRICA-CORDEX44 (AFR44; Table 2), three RCMs from AFRICA-CORDEX22 (AFR22; Table 3) and four GCMs from CMIP5 (Table 4) to assess ROM performance. The validation strategy addressed:

  • Larger scale: The CCUS climate at the larger scale was evaluated by means of OISSTv2 (SST) and ERA5 (wind stress) data. The NOAA 1/4° daily OISST (Reynolds et al. 2007) is an analysis constructed by combining observations from different platforms (satellites, ships, buoys, and Argo floats) on a regular global grid. ERA5 (Copernicus Climate Change Services 2017) is a product of the European Centre for Medium Range Weather Forecast (ECMWF) and covers the Earth on a 30 km grid and resolves the atmosphere on 137 levels from the surface up to a height of 80 km. We also used the European Space Agency (ESA) SST Climate Change Initiative (CCI) and C3S global Sea Surface Temperature Reprocessed product (Good et al. 2019), which provides daily averaged SST with a resolution of 5 km (Merchant et al. 2019).

  • Latitudinal variability: Between the four major EBUSs, the CCUS stands out as the most spatially and seasonally variable one in terms of primary production (Sylla et al. 2019). Based on the seasonality, there are different arguments to divide the CCUS in sub-regions but without a broad consensus (Aristegui et al. 2009; Pardo et al. 2011; Cropper et al. 2014; Gómez-Letona et al. 2017; Fischer et al. 2019). Here we will focus on the seasonal and latitudinal variability of T2m, wind and SST over the coastal band, allowing an integrated view of the coastal upwelling. The local wind was assessed using an upwelling index (UI) based on the offshore wind-driven Ekman transport (Q). The UI proposed by Bakun (1973) was calculated following the approach used by Gomez-Gesteira et al. (2006), Cropper et al. (2014) and Sousa et al. (2017a, b):

    $$Q_{x} = \frac{{\tau_{y} }}{f\rho }$$
    $$Q_{y} = - \frac{{\tau_{x} }}{f\rho }$$
    $$UI = - \sin \left( {\theta - \frac{\pi }{2}} \right)Q_{x} + \cos \left( {\theta - \frac{\pi }{2}} \right)Q_{y}$$

    where Qx, Qy and \({\tau }_{x}\), \({\tau }_{y}\) are the zonal and meridional components of the horizontal Ekman transport and the wind stress vector, respectively; \(\rho\) is the reference sea water density (1025 kg m−3); \(f\) is the Coriolis parameter and \(\theta\) is the angle between the coastline and the equator.

  • Thermal vertical structure: The cross-shelf thermal vertical structure is relevant for the characterization of the coastal upwelling intensity. The vertical structure of the ocean temperature was assessed through SODA (Carton et al. 2018) and GLORYS12V1 (Drèvillon et al. 2018) gridded products and WOD18 (Boyer et al. 2019) in situ data. SODA3.4.1 is the latest version of SODA with 1/4º × 1/4º resolution in 50 vertical levels. The GLORYS12V1 product is the Copernicus Marine Environment Monitoring Service (CMEMS) global ocean eddy-resolving reanalysis with observations assimilated from along track altimeter data, satellite sea surface temperature, sea ice concentration and in situ temperature and salinity vertical profiles. The WOD18 includes in situ measurements of temperature, salinity, dissolved oxygen and nutrients from 1773 to the present. For the evaluation we used a cross-shelf transect along 31º N (Cape Ghir), a region where the spatial resolution could play a role due to the bathymetry features and coastline.

  • Mesoscale events: Coastal upwelling filaments are mesoscale structures with a great importance in the transport of organic carbon and nutrient-rich upwelled waters several hundred kilometers offshore and have typical lifetime of a few weeks (Menna et al. 2016), fueling the biological activity of the oligotrophic open waters (Lovecchio et al. 2018). We use two high resolution datasets able to detect the filaments: MODIS Aqua (Hall and Riggs 2007) and GLORYS12V1 (Drèvillon et al. 2018). The Moderate Resolution Imaging Spectroradiometer (MODIS) is a key instrument aboard Aqua satellite, which views the entire Earth’s surface every 2 days, acquiring data in 36 spectral bands (Hall and Riggs 2007).

Table 1 Observational and reanalysis data products selected to assess ROM performance in the CCUS
Table 2 AFR44 simulations used to compare ROM with RCMs
Table 3 AFR22 simulations used to compare ROM with RCMs
Table 4 CMIP5 simulations used to compare ROM with GCMs

2.3 Model intercomparison

We use the Max Plank Institute for Meteorology—Earth System Model (MPI-ESM; Giorgetta et al. 2013) in order to assess the impact of the high resolution of ROM, particularly regarding mesoscale dynamics. MPI-ESM has the same ocean component (MPIOM) than ROM. Here we use the two versions of MPI-ESM contributing to CMIP5. The first, Low Resolution (-LR) uses for the ocean a bipolar grid with 1.5º resolution. The Medium Resolution (-MR) has a tripolar grid with higher ocean horizontal resolution of 0.4º (Giorgetta et al. 2013).

Additionally, we compare ROM output to the ensemble mean of several Regional Climate Models (RCMs) from the AFRICA-CORDEX, with boundary conditions forced by ERA-Interim. We should note that in AFRICA-CORDEX the available simulations are atmosphere-only and the horizontal resolution available were 0.22º and 0.44º. Following Sousa et al. (2017b), we chose seven AFR44 (Table 2) and three AFR22 (Table 3) in the study.

Finally, also following Sousa et al. (2017a), ROM is compared with an ensemble of four Global Climate Models (GCMs) from the CMIP5 project (Table 4).

3 Results

3.1 Larger scale

In this section we assess the larger scale climate of the CCUS through the SST, the wind stress and the T2m from the Iberian Peninsula to Cape Blanc. Besides the coastal area, we also analyze an extensive offshore area that contains the CCUS (see the area of study in Fig. 2).

Fig. 2
figure 2

SST biases (ºC) in DJF (upper row) and in JJA (lower row) between ROM (a, e), MPI-ESM-LR (b, f), MPI-ESM-MR (c, g) and CMIP5 (d, h) with OISST

The mean SST field in the CCUS for the 1982–2012 period presents a year-round meridional gradient and a colder patch by the NW African and Iberian coasts because of the upwelled waters (Fig. 1a). For the evaluation of SST, we compare the output of ROM with OISST, MPI-ESM-LR, MPI-ESM-MR and CMIP5 for the 1982–2012 period (Fig. 2). In DJF, ROM provides a good agreement with OISST, displaying biases smaller than 1.0 °C over most of domain except for the northwest region, where the differences reach 2.5 °C (Fig. 2a). MPI-ESM-LR and MPI-ESM-MR present biases up to 2.0 °C with larger differences at local regions such as Cape Bojador and Cape Blanc (Fig. 2b, c). CMIP5 ensemble mean shows an overall cold bias with a maximum in the Gulf of Cádiz (Fig. 2d).

In JJA, the cooling of coastal surface waters evidences the intensification of the upwelling system. This mechanism seems to be magnified by ROM, which despite of having small differences with OISST throughout the domain presents a cold bias of 2.0 °C in the coastal strip between 24º N and 34º N (Fig. 2e), corresponding to the permanent upwelling zone with summer intensification (Arístegui et al. 2009). MPI-ESMs show a patchy bias distribution along the coast (− 2.0 to + 2.0 °C) presenting a cold bias (3.0 °C) in the northwest corner of the domain (Fig. 2f, g). CMIP5 biases present a similar patchy pattern but with lower values along the NW African coast (Fig. 2h), where the horizontal resolution seems to play an important role.

The interannual variability and the seasonal cycle of the spatially averaged SST were evaluated using OISST, ERA5 and ESA data sets (Fig. 3). ERA5 and ESA were added to the analysis to account for uncertainty on SST observation-derived data. The time series of yearly spatially averaged ROM SST for the period 1980–2012 shows warm and cold biases ranging from + 0.59 to − 0.28 °C compared to OISST (Fig. 3a), being mostly within the observation data sets spread range.

Fig. 3
figure 3

a Times series of yearly mean (a) and seasonal cycle (b) (1982–2012) SST (ºC) averaged over the domain

Figure 3b shows the seasonal cycle of ROM SST compared to OISST, ERA5 and ESA. Minimum (maximum) temperatures are reached in February (September), being the amplitude of the seasonal cycle 5.4 °C for ESA and 5.5 °C for ERA5, OISST and ROM.

Wind stress is the main driver of the CCUS (Fischer et al. 2019); erroneous intensities or directions of wind stress can have a strong impact on the seasonal upwelling system. To assess the ability of ROM reproducing the CCUS we analyze the seasonal wind stress intensity and direction compared to ERA5, MPI-ESM-LR, MPI-ESM-MR and CMIP5 during 1980–2012 period (Fig. 4). In DJF, when the anticyclonic circulation is weaker, the mean ERA5 wind stress intensity increases to north and south from 34º N where a zonal band with minimum values is located. Moreover, low wind stress values extend along the coasts of the Iberian Peninsula and NW Africa, increasing gradually to the south of the Canary islands (Fig. 4a). ROM overestimates the ERA5 wind stress due to an overestimation of the Azores high during the winter months (Fig. 4b). In both MPI-ESM runs (Fig. 4c, d) this overestimation is even larger. On the contrary, the CMIP5 ensemble mean underestimates the strength but shows a spatial pattern closer to the wind stress depicted by ERA5 (Fig. 4e). Despite the larger scale differences, ROM adequately represents the spatial pattern induced by coast in the wind stress because of its high horizontal resolution.

Fig. 4
figure 4

a Wind stress intensity and direction for the winter (upper row) and summer (lower row) months. It is shown ERA5 (a, f), ROM (b, g), MPI-ESM-LR (c, h), MPI-ESM-MR (d, i) and CMIP5 (e, j)

In JJA, the Azores high strengthens and migrates to the NW, where ERA5 shows the lower intensities of wind stress (Fig. 4f). In ERA5, the wind stress increases to the south showing the largest values in the coastal strip from Cape Ghir to the south, and downwind from Canary Islands. By the Iberian and NW African coastlines the wind direction is along the coast, i.e. upwelling favourable. As in winter, MPI-ESMs strongly overestimate the wind stress field strength (Fig. 4h, i), while CMIP5 underestimates it, but both properly reproduce the wind stress field directions (Fig. 4j). Improving on the GCMs performance, ROM reproduces the JJA wind stress field remarkably well, including smaller scale features such as the local maximum off Cape Ghir, the Madeira and Canary Islands shadowing effect, the wind intensification in the passage between Canary Islands and Africa or the coastal weakening effect (Fig. 4g).

3.2 Latitudinal variability over the coastal band

The CCUS is a dynamically complex system with a large spatial and seasonal variability. Here, we will assess SST, wind and T2m variability focusing on the coastal band as the upwelling front is located within a few degrees offshore from the shelf break (Pelegrí and Bennazzouz 2015), supplementing the analysis of the previous section.

The coastal cold SST is a key indicator to determinate the intensity of the upwelling system. The seasonal cycle of ROM was compared with the seasonal cycle of OISST over the coastal band (the mask is shown in Fig. 1a with red lines) and with the GCMs (MPI-ESM-LR, MPI-ESM-MR and CMIP5 ensemble). ROM performs better than GCMs with biases smaller than 0.5ºC in winter. It is also evident the effect of the coarse resolution in GCMs, especially in MPI-ESM-LR and CMIP5 (Fig. 5a–d).

Fig. 5
figure 5

SST (ºC) biases in the closest grid-points to the coast (red line Fig. 1a; 100 km from Coast to Offshore), between ROM (a, e), MPI-ESM-LR (b, f), MPI-ESM-MR (c, g) and CMIP5 (d, h) with OISST, in DJF (upper row) and JJA (lower row)

In summer boreal months, when the Canary upwelling intensifies, ROM shows very small warm biases north of the Strait of Gibraltar and cold biases in the region of the NW African coast (Fig. 5e). Remarkably at 25º N ROM cold bias extends uniformly offshore. MPI-ESM-LR shows the largest warm biases (up to 3.0 °C) in Cape Ghir and the Gulf of Cádiz. MPI-ESM-MR and CMIP5 present cold biases in the Gulf of Cádiz and Cape Bojador (2.0 °C) and warm biases in Cape Ghir and Cape Blanc (Fig. 5g, h).

From these results it could be concluded that ROM magnifies in summer the upwelling in the zone of weak permanent upwelling located between 25° N and 33° N (Cropper et al. 2014; Gómez-Letona et al. 2017). To clarify the role of the high resolution in this issue we compared the ROM biases with respect to OISST, ERA5 and ESA in the grid points that are closer to the coast (Fig. 6). In DJF ROM biases are similar for the three datasets, with cold biases smaller than 1.0 °C along the African coast. From Cape Ghir to Cape Blanc, the biases with ESA are notably smaller (Fig. 6a). In JJA, the differences in the biases increase, being close to 0.5 °C in Cape Ghir and Cape Bojador for ESA and reaching the 4.0 °C for ERA5 at those locations (Fig. 6b). In general, ROM differences with ESA are smaller than with OISST and ERA5 along the African coast. The reason for these discrepancies can be related to the OISST and ERA5 resolution, which does not allow for a clear representation of the SST coastal pattern of the upwelling.

Fig. 6
figure 6

ROM SST (ºC) biases in the closest grid-point to coast with OISST, ESA and ERA5, in DJF (a) and JJA (b)

The coastal wind stress was assessed through UI (Fig. 7) in order to quantify the upwelling intensity from Ekman transport. Positive (negative) UI corresponds to upwelling (downwelling) conditions. In DJF (Fig. 7a), ERA5 presents nearly neutral conditions along the coasts of the Iberian Peninsula. From the Gulf of Cádiz to Cape Ghir ERA5-based UI is positive increasing towards the south, where from Cape Bojador to Cape Blanc the upwelling becomes intense. This latitudinal variability is well represented by ROM and AFRICA-CORDEX (AFR44 and AFR22), highlighting higher values in ROM for Cape Ghir and in southern regions and weaker in AFRICA-CORDEX for the region between the Gulf of Cádiz to Cape Blanc. CMIP5 also represents reasonably well the latitudinal variability, although it is smoother and with a general tendency to underestimate the upwelling index.

Fig. 7
figure 7

UI (m3 s−1 km−1) averaged over the closest grid-points to the coast (c), for DJF (a) and JJA (b). ROM is assessed with ERA5 and compared with MPI-ESM-LR, MPI-ESM-MR, CMIP5, AFR44 and AFR22

In JJA, ERA5-based UI shows upwelling conditions along the Iberian coast and a clear intensification of the upwelling by the NW African coast, appearing a local maximum between Cape Beddouza and Cape Ghir (30.5º N–32.5º N). ROM and AFRICA-CORDEX present the same latitudinal pattern as ERA5, while CMIP5 reproduces the latitudinal variability, but with lower values of the UI and MPI-ESMs show different local maxima. In general, AFRICA-CORDEX shows a slightly better performance than ROM and CMIP5, with MPI-ESMs having the worst coefficients of determination (Table 5).

Table 5 Latitudinal coefficient of determination for UI averaged over the closest grid-points to coast in DJF and JJA, comparing ERA5 with ROM, MPI-ESM-LR, MPI-ESM-MR, CMIP5, AFR44 and AFR22

The T2m land-sea difference lies in the ground of Bakun’s (1990) hypothesis of change in upwelling systems under global warming. The land-sea gradient simulated by ROM is validated against ERA5 and compared with MPI-ESM-LR, MPI-ESM-MR, CMIP5 and AFRICA-CORDEX simulations. For each latitude, we calculate the difference between the 2 m air temperature zonally averaged over 100 km inshore (bounded by the red line, land) and zonally averaged over 100 km offshore (bounded by the blue line, ocean) as shown in the Fig. 8c.

Fig. 8
figure 8

Along-shore land (c; blue)—sea (c; red) temperature difference (ºC) for DJF (a) and JJA (b). ROM is assessed with ERA5 and compared with MPI-ESM-LR, MPI-ESM-MR, CMIP5, AFR44 and AFR22

ERA5 T2m land-sea differences present mostly negative values in DJF, changing sign in the southern part from 21º N to 23º N (Cape Blanc). AFICA-CORDEX reproduces accurately the latitudinal variability, while ROM and CMIP5 fail in simulating the sign change in the southernmost region. The MPI-ESMs do not reproduce properly the latitudinal variability and show too low values of the T2m land-sea difference (around 0.5 °C).

In JJA the insolation increases and the land becomes warmer than the sea, so T2m land-sea differences are mostly positive (Fig. 8b). ERA5 presents a high latitudinal variability, with two regions showing the largest positive differences (Cape Ghir and from Cape Bojador to Cape Blanc). ROM and AFRICA-CORDEX again reproduce well the T2m land-sea differences. CMIP5 does not reproduce properly the T2m land-sea differences latitudinal pattern, showing a spurious maximum in the Gulf of Cádiz, and the MPI-ESMs showed very low values compared to those of ERA5. AFRICA-CORDEX and ROM show a better performance, especially in summer when CMIP5 and MPI-ESMs present a rather low coefficient of determination (Table 6).

Table 6 Latitudinal coefficient of determination for T2m land-sea differences averaged over the closest grid-points to coast in DJF and JJA, comparing ERA5 with ROM, MPI-ESM-LR, MPI-ESM-MR, CMIP5, AFR44 and AFR22

3.3 Thermal vertical structure

The cross-shelf thermal vertical structure is commonly used to characterize coastal upwelling. The thermal vertical structure of the upper 200 m simulated by ROM was compared to SODA and GLORYS reanalysis and to in-situ data from WOD18. We choose a cross-shelf transect at Cape Ghir (31.5º N). All observations used from WOD18 were taken between 1980 and 2012, and the isotherms were plotted as a qualitative reference (Fig. 9a, b).

Fig. 9
figure 9

Temperature (ºC) transect for JJA (left) and DJF (right) in Cape Ghir (1980–2012) for WOD18 (a, b), GLORYS (c, d), ROM (e, f) and SODA (g, h)

At Cape Ghir WOD18 presents a clear thermal stratification, more evident in summer, with isotherms sloping up to the coast and a temperature range between 14 °C and 22 °C. In summer the upwelling front outcrops in the surface near 10.3º W (see 17 °C isotherm), while it is relaxed in winter with the 17 ºC isotherm outcropping more to the east by the coast (Fig. 9a, b). The cross-shelf structure and seasonality is properly represented by GLORYS and SODA reanalysis, but ROM is more accurate in the representation of the vertical gradients and the isotherm sloping by the coast (Fig. 9c–h).

3.4 Upwelling filaments

The MPIOM high horizontal resolution in the central region of the CCUS allows marginally the representation of the mesoscale variability because the climatological first baroclinic Rossby radius of deformation is around 30 km in this area (Chelton et al. 1998). Upwelling filaments are elongated mesoscale structures of upwelled water extending offshore in the upper surface layer (Brink 1983). Alvarez-Salgado et al. (2007) and Lovecchio et al. (2018) highlighted the relevant contribution of upwelling filaments to offshore organic carbon in the CCUS. The most prominent filaments in the CCUS are located at Cape Ghir, Cape Juby, Cape Bojador and Cape Blanc. Cape Ghir filament, one of the largest and existing nearly all year round (Hagen et al. 1996), exports large amounts of organic material into the open ocean (García-Muñoz et al. 2005; Pelegrí et al. 2005, 2006). To assess ROM skills in reproducing these relevant mesoscale features we use two filament events in August 2006 and August 2009, the latter studied in detail by Sangrá (2015) and Sangrá et al. (2015). We compare the observed MODIS Aqua SST with ROM and GLORYS output averaged for 21st to 28th August 2006 and 13th to 20th August 2009 (Fig. 10).

Fig. 10
figure 10

Averaged SST (ºC) in 21–28 August of 2006 (upper row) and 13–20 August of 2009 (lower row) for MODIS Aqua (a, d), ROM (b, e) and GLORYS (c, f)

As described by Sangrá et al. (2015), ROM reproduces the filament cold core, with SST below 19 °C and a broader cool embracing region with temperatures between 19 °C and 21 °C. ROM reproduces properly both filament events, quite accurately the offshore extension of the embracing region and a overextended cold core, while GLORYS filament is excessively limited to the coast.

4 Discussion

When Bakun (1990) formulated his hypothesis highlighting that under global warming the increase in ocean-land temperature contrast may oppose the general tendency for ocean-basin scale circulation to slow down, the enormous uncertainty about the evolution of coastal ocean upwellings under global climate change conditions was evidenced. Additionally, the importance of the mesoscale in the decadal changes of the CCUS is widely acknowledged (e. g. Relvas et al. 2009), as it is the role of the mesoscale in transporting nutrients and phytoplankton (e.g. Gruber et al. 2011). Since then many attempts have been done to get rid of such uncertainty and to deepen the knowledge about the behaviour of such fundamental marine ecosystems under global warming conditions, and different observational and model products have been used to resolve the different scales relevant for these questions. The study of the upwelling systems in the CMIP5 simulations (e. g. Oyarzun and Brierley et al. 2018) allowed identifying the climate change drivers and the main limitations of the GCMs to reproduce the mesoscale processes in the EBUSs (Garcia-Reyes et al. 2015; Bindoff et al. 2019). These local and mesoscale features may have greater impact on EBUSs than large-scale wind patterns (Renault et al. 2016; Xiu et al. 2018). The difficulties found in the representation of the CCUS by global models can be significantly resolved with the help of high resolution AORCMs (e. g. Xiu et al. 2018). In this work, the use of ROM in the CCUS region allows us to represent the fine-scale atmosphere–ocean feedback, reproducing the spatial and temporal structure needed in the regional climates (Li et al. 2012; Sein et al. 2015). Moreover, ROM configuration allows a high resolution in the region of interest maintaining a global ocean domain (Soares et al. 2018; Parras-Berrocal et al. 2020).

The provided climate assessment shows that ROM reproduces the basin scale climate and the seasonal signal affecting the CCUS (Sect. 3.1) better than the GCMs (MPI-ESMs and CMIP5) and is in good agreement with OISST and ERA5. ROM presents also good skills in reproducing the interannual variability of the basin scale SST when compared to ERA5, OISST and ESA data sets (Fig. 3). However, in winter ROM has a warm bias in the northwest corner of the domain (Fig. 2a). We have also identified an overestimation of the Azores high in winter, originating a wrong wind stress field representation along the West Iberian coast (Fig. 4b). Parras-Berrocal et al. (2020) attributed those differences to the role played by the deficiencies of ROM to simulate the ocean circulation in the North Atlantic, with a cold bias centred east of Flemish Cap. Other studies stand out the difficulties of the ocean models to simulate the Gulf Stream in winter (Eden and Greatbatch 2003; Bryan et al. 2007) and Cabos et al. (2020) detected those deficiencies not only in the ocean models, but also in the Climate Forecasting System Reanalysis. However, according to our results the impact of those winter basin-scale inaccuracies on the CCUS seems to be quite limited, as they impact mostly in the West Iberian coast during the no upwelling season.

ROM also presents a cold bias compared to OISST along the NW African coast in summer, when the upwelling intensifies (Figs. 2 and 5). Mason et al. (2011) and Santana-Falcón et al. (2020) found a similar cold bias in their ocean regional models over the CCUS, attributing it to the warm bias in SST forcing data set identified by Dufois et al. (2012) in a number of modelling studies in the world EBUS. According to Dufois et al. (2012) the cause of the SST warm coastal bias in monthly Pathfinder data during summer was a flagging method based on an OISST reference test, which is explained by strong coastal SST gradients in these regions, which cannot be satisfactorily represented by the large scale OISST product. This drawback of OISST can partly explain the strong cold bias in ROM with respect to this observational product. Indeed, when we compared ROM SST with a higher resolution dataset (ESA) the cold bias was notably reduced (Fig. 6), demonstrating the importance of the high resolution in mesoscale regions as the CCUS.

This high ocean resolution is key for the representation of the CCUS mesoscale, which is determined by the accuracy in the representation of the coastal wind stress field and the land–ocean temperature contrast and characterized by the cross shelf thermal structure and the generation of events such as upwelling filaments.

ROM properly distinguishes the three main CCUS sub-regions (e.g. Aristegui et al. 2009, Cropper et al. 2014; Gómez-Letona et al. 2017) from the UI (Fig. 9): from 21º N to 26º N (permanent upwelling zone), from 26º N to 33º N (weak permanent upwelling zone) and from 33º N to 43º N (seasonal upwelling zone), although in part of the West Iberian coast ROM produces a nearly year-round upwelling (Fig. 7) as a consequence of the misrepresentation of the wind field in that area as previously discussed. Despite this fact, ROM performs better than the MPI-ESMs and the CMIP5 ensemble, and slightly worse than the AFRICA-CORDEX ensemble in terms of the UI (which is calculated from the wind stress). It must be taken into account that, the lowest resolution of AFRICA-CORDEX configuration (AFR44) reproduces better the latitudinal variability than higher resolution (AFR22). This fact may be related to the smaller ensemble in AFR22 (only 3 members vs 7 members in AFR44) and the different regional models in both ensembles. Moreover, although AFR44 models have lower horizontal resolution than ROM atmosphere (0.44º vs 0.22º) their west and north boundaries, forced by ERA-Interim, are pretty close to the CCUS and that the SST is also prescribed from ERA-Interim. The UI values along the Iberian Peninsula simulated by ROM are similar to other studies, with a localised downwelling by the Galician coasts (− 200 m3 s−1 km−1) in winter and weak upwelling (500 m3 s−1 km−1) in the summer months (Pardo et al. 2011; Sousa et al. 2017a). The upwelling in the African coast also is in the range obtained by different studies (Pardo et al. 2011; Sousa et al. 2017b), with year-round permanent upwelling (from 600 to 1000 m3 s−1 km−1) in winter, increasing to 1500 m3 s−1 km−1 in the African coastline.

The T2m land-sea difference and its latitudinal variability is relevant for CCUS representation as Bakun hypothesis is based on the global warming induced changes in the T2m land-sea differences. Several authors (Sydeman et al. 2014; Wang et al. 2015) supported Bakun’s proposition, but remains unclear the relationship between coastal wind intensification and global warming (García-Reyes et al. 2015). Our results have shown a good performance of ROM representing the T2m land-sea contrast variability in CCUS as compared to ERA5. The latitudinal variability of the land sea temperature contrast is strongly reduced in the GCMs, indicating that global models have difficulties to represent the land-sea thermal differences in coastal regions (Ning et al. 2014; Roxy et al. 2015; Ward et al. 2020). ROM does not reproduce the sign change of T2m land-sea difference in Cape Blanc, while it is properly captured by AFRICA-CORDEX. In JJA, ROM overestimates the T2m land-sea differences in Cape Ghir as a consequence of a strong warming in the Atlas mountains.

Therefore, the ROM improvement over the GCMs and good performance as compared to AFRICA-CORDEX shows the impact of the higher resolution and coupling for CCUS climate modelling. Coming to details, we have identified the Cape Ghir as an intense upwelling location both in DJF as JJA, and as a place where upwelling filaments are a nearly permanent mesoscale feature. ROM properly reproduces the summer increase in ocean thermal stratification in the upper 100 m at Cape Ghir and the increase in the slope of the isotherms in the immediate proximity to the coast (Fig. 9), comparing well with WOD18 hydrographic data and with GLORYS reanalysis. The improvement over SODA is remarkable, especially over the shelf where the high resolution plays a key role.

The mesoscale events play an important factor over the CCUS, where only organic carbon export by filaments in the subtropical northeast Atlantic contributes to the 63% of the annual primary production associated with the coastal upwelling (Santana-Falcón et al. 2016, 2020). ROM reproduces the two selected events in Cape Ghir, exporting the water masses from coast to 150 km offshore (Fig. 10). Quite often limited area models are used for reproducing these events with the difficulties of dealing with the open boundary conditions. This is even more important when considering the biogeochemical relevance of the upwelling filament dynamics. Troupin et al. (2012) and Santana-Falcón et al. (2020) approached the study of upwelling filaments using nested domains and a one-way offline coupling. However, when the focus is on the basin-scale, approaches accounting for the influence of the mesoscale on the larger scale are needed. Lovecchio et al. (2017, 2018) used such an approach by means of a telescopic curvilinear grid with strong refinement in the NW African coast, being able to highlight the key role of mesoscale processes in the offshore transport of organic carbon and concluding that a great part of this flux out of the upwelling regions is not accounted for in coarse global models. Here, we present a model approach capable of reproducing the large-scale climate signal accounting for the relevant upwelling mesoscale dynamics, which is of uttermost importance to assess the future evolution of CCUS and its socio-economic consequences under climate change scenarios.

5 Conclusions

In the present work, we assess the ability of the atmosphere–ocean regionally coupled model ROM to represent the large-scale climate and mesoscale processes involved in the CCUS dynamics, and compare it to two MPI-ESMs configurations and ensembles of CMIP5 GCMs and AFRICA-CORDEX RCMs. Our findings can be summarized as follows:

  • ROM shows a better performance than the MPI-ESMs and CMIP5 models in representing the larger-scale wind stress and SST fields, although a relatively overestimation of the Azores high in winter affects slightly the CCUS over the West Iberian coast. Besides, ROM reproduces adequately the seasonal and interannual variability of the ERA5, ESA and OISST.

  • High resolution is key to reproduce the latitudinal variability of the CCUS, as ROM represents the observed coastal SST with higher accuracy than the GCMs. Moreover, ROM successfully reproduces the coastal UI as well as the T2m land-sea differences, highlighting the impact of the higher resolution against the GCMs with a performance comparable to AFRICA-CORDEX.

  • The mesoscale processes in the CCUS are well simulated by ROM, which is able to transfer the coastal upwelling waters to the open ocean in a realistic way. Thus, it successfully represents two coastal upwelling filaments off Cape Ghir, which are not accounted for in most global models.

In conclusion, ROM is a powerful atmosphere–ocean model system able to reproduce with accuracy the CCUS, performing better than the GCMs. The improvement is related to a much higher horizontal resolution, which allows a better simulation of the dominant mesoscale coastal dynamics. The results here give ground to the future use of ROM to gain a deeper insight into the CCUS by the end of twenty-first century.