Modelling the impact of flow-driven turbine power plants on great wind-driven ocean currents and the assessment of their energy potential

1Institut des Géosciences de l’Environnement, Centre National de la Recherche Scientifique, Université Grenoble Alpes, CNRS/UGA/IRD/G-INP, Grenoble, France. 2Shirshov Institute of Oceanology, Russian Academy of Sciences, Moscow, Russia. 3School of Environmental Sciences, University of Liverpool, Liverpool, UK. 4Laboratoire des Ecoulements Geophysiques et Industriels, Grenoble, France. 5Ocean Next, Grenoble, France. ✉e-mail: bernard.barnier@cnrs.fr Mitigating human impacts on the climate and the environment has become a pressing challenge that makes ocean sources of renewable energy particularly attractive. Great wind-driven ocean currents flowing along the western boundaries of ocean basins (for example the Gulf-Stream or Kuroshio) are among the largest renewable energy resources on the planet1,2. These currents are relatively persistent in intensity and direction, which suggests that turbine power plants (TPPs) immersed into these currents could have a higher capacity of utilization than those exploiting tidal or wind energies2,3. With the advances in the development of power-generating technologies for ocean turbines2–7, the proper identification of the most favourable locations for TPPs and the accurate assessment of their energy potential are becoming pressing issues8–10. A variety of approaches are used to evaluate the potential exploitation of great ocean currents, particularly the Gulf-Stream and Kuroshio, as a source of electrical power. Computations of the kinetic energy flux through cross-sections where main current characteristics are observed have yielded crude estimates of a total power reserve of ~20 GW for the Florida Current11,12 and of ~30 GW for the Kuroshio Current13. Attempts to improve power resource estimates by accounting for the time variability of ocean currents by employing realistic model simulations4,13–16 faced limitations associated with possible underestimations of the modelled currents’ magnitudes16. This limitation was partly overcome in a few cases thanks to re-scaling of the modelled currents using local measurements (where they existed)13,15. Attempts to characterize the power extractable without creating drastic alterations of flow energetics (rather than estimating the global resource) have been proposed. These studies use simplified models9,17,18 (which provide much-idealized representations of western boundary currents (WBCs)) or regional realistic models19 and account for the impact of TPPs by increasing dissipation locally in the WBC. Varying the strength of dissipation, they search for the ‘acceptable dissipation limit’ that the current system can support, and interpret it as the maximum extractable power. San (2016)18 for example applied such an additional friction over 100 km along an idealized WBC comparable with the Gulf-Stream and estimated that about 10 GW could be extracted with only a slight altering of the mean flow. Using a realistic model of the North Atlantic, Haas et al. (2018)19 found an important sensitivity of the main flow path and the extracted power to various scenarios of turbine distribution in the Florida Current. Another aspect of the problem is how to precisely identify those sites that are favourable for energy extraction. Current velocities obtained from satellite or surface drifter observations are used for this purpose8,10,13,14 because of their global coverage. However, they cannot provide an accurate estimation of the overall power of the current because they do not sample the vertical profile. The power density, which provides an estimate of the power available per unit area of an energy production device, is used instead. These approaches, which developed relevant multi-parameter criteria for site selection8, are somewhat adequate for identifying regions and providing preliminary estimates of energy resources in the case of a small number of devices that induce negligible changes on the existing flow, but are not pertinent for large TPPs that are able to cause large flow changes. In the present study we go beyond the above limitations by using the time-dependent currents produced by a realistic, state-of-theart, global eddy-resolving ocean circulation model. Our strategy is threefold. First, following the approach used by, for example, Chang et al. (2015)8, we use a multi-parameter criterion to search for the most favourable sites. Further, we implement virtual TPPs, all identical, at the selected sites located in the WBCs of the world ocean. Because the constraint of the TPP on the current is the same at every location, the sensitivity of the various sites can be compared. The formulation of the effect of the turbines is comparable with that followed by Haas et al. (2014)19 or used in the works analysing Modelling the impact of flow-driven turbine power plants on great wind-driven ocean currents and the assessment of their energy potential

M itigating human impacts on the climate and the environment has become a pressing challenge that makes ocean sources of renewable energy particularly attractive. Great wind-driven ocean currents flowing along the western boundaries of ocean basins (for example the Gulf-Stream or Kuroshio) are among the largest renewable energy resources on the planet 1,2 . These currents are relatively persistent in intensity and direction, which suggests that turbine power plants (TPPs) immersed into these currents could have a higher capacity of utilization than those exploiting tidal or wind energies 2,3 . With the advances in the development of power-generating technologies for ocean turbines 2-7 , the proper identification of the most favourable locations for TPPs and the accurate assessment of their energy potential are becoming pressing issues [8][9][10] .
A variety of approaches are used to evaluate the potential exploitation of great ocean currents, particularly the Gulf-Stream and Kuroshio, as a source of electrical power. Computations of the kinetic energy flux through cross-sections where main current characteristics are observed have yielded crude estimates of a total power reserve of ~20 GW for the Florida Current 11,12 and of ~30 GW for the Kuroshio Current 13 . Attempts to improve power resource estimates by accounting for the time variability of ocean currents by employing realistic model simulations 4,13-16 faced limitations associated with possible underestimations of the modelled currents' magnitudes 16 . This limitation was partly overcome in a few cases thanks to re-scaling of the modelled currents using local measurements (where they existed) 13,15 .
Attempts to characterize the power extractable without creating drastic alterations of flow energetics (rather than estimating the global resource) have been proposed. These studies use simplified models 9,17,18 (which provide much-idealized representations of western boundary currents (WBCs)) or regional realistic models 19 and account for the impact of TPPs by increasing dissipation locally in the WBC. Varying the strength of dissipation, they search for the 'acceptable dissipation limit' that the current system can support, and interpret it as the maximum extractable power. San (2016) 18 for example applied such an additional friction over 100 km along an idealized WBC comparable with the Gulf-Stream and estimated that about 10 GW could be extracted with only a slight altering of the mean flow. Using a realistic model of the North Atlantic, Haas et al. (2018) 19 found an important sensitivity of the main flow path and the extracted power to various scenarios of turbine distribution in the Florida Current.
Another aspect of the problem is how to precisely identify those sites that are favourable for energy extraction. Current velocities obtained from satellite or surface drifter observations are used for this purpose 8,10,13,14 because of their global coverage. However, they cannot provide an accurate estimation of the overall power of the current because they do not sample the vertical profile. The power density, which provides an estimate of the power available per unit area of an energy production device, is used instead. These approaches, which developed relevant multi-parameter criteria for site selection 8 , are somewhat adequate for identifying regions and providing preliminary estimates of energy resources in the case of a small number of devices that induce negligible changes on the existing flow, but are not pertinent for large TPPs that are able to cause large flow changes.
In the present study we go beyond the above limitations by using the time-dependent currents produced by a realistic, state-of-theart, global eddy-resolving ocean circulation model. Our strategy is threefold. First, following the approach used by, for example, Chang et al. (2015) 8 , we use a multi-parameter criterion to search for the most favourable sites. Further, we implement virtual TPPs, all identical, at the selected sites located in the WBCs of the world ocean. Because the constraint of the TPP on the current is the same at every location, the sensitivity of the various sites can be compared. The formulation of the effect of the turbines is comparable with that followed by Haas et al. (2014) 19   the energy potential of tidal streams [20][21][22] although these studies follow different paradigms. Finally, we perform simulations of the modifications to global ocean currents caused by large TPPs, thereby providing valuable insights into the assessment of the potential impact of TPPs on circulation patterns and the accuracy of power-potential estimates.

Great ocean currents in eddy-resolving simulations
We use here the last five years of a global 33-year-long (1979-2011) simulation (the control experiment) performed with the ORCA12 (refs. 23,24 ) eddy-resolving global ocean circulation model. ORCA12 has a horizontal grid resolution of 1/12° (9.25 km at the equator and 7 km at mid-latitudes). See Methods for a detailed model description. The model is routinely used by the Copernicus Marine Environment Monitoring Service 25 (CMEMS) for operational forecasts 26,27 and by the research community for climate-oriented long-term simulations [28][29][30][31] and process studies 32,33 . This model has been shown to accurately reproduce most of the observed characteristics of the major currents and mesoscale turbulence of the upper ocean [34][35][36] .
ORCA12 currents are assessed against observations using the mean current climatology (1/2° resolution) of the National Oceanic and Atmospheric Administration (NOAA) Surface Velocity Program (SVP) drifter data 37 and the multi-observation threehourly analysis (1/4° resolution) of CMEMS 38 (see Methods). The agreement between modelled and observed currents is quite good for magnitude, width and path ( Fig. 1) although slightly better between model and surface drifter currents in terms of amplitude than for CMEMS. The only noticeable discrepancy is the meander of the Kuroshio south of Kyushu Island being excessively amplified by the model. Thanks to the high resolution, the model shows sharper WBCs and stronger-than-observed flows in small straits like in the Cozumel Channel at the tip of the Yucatan Peninsula, or in the Tsugaru Strait separating the islands of Hokkaido and Honshu, where we installed virtual TPPs (Fig. 2). The agreement is also good between CMEMS and ORCA12 on the instantaneous currents (see Methods). A comparison with local current-meter mooring observations in the Florida Current (See Supplementary Note 1) shows that the grid-sized model currents reproduce the vertical shear of the current rather well, but have a slightly lower magnitude than the locally observed current, as is true for most global eddying ocean circulation models with an equivalent resolution 4,15,16,39 .
Summarizing, the modelled currents show surface intensification, time-mean velocity and variability closely comparable with observations. The surface intensification (a fundamental feature of great ocean currents) suggests that TPPs designed to extract energy should be placed at a depth close to the surface. Thus, we consider a depth range from 20-46 m, consistent with the 20 to 60 m operating depth range suggested in similar studies 3,4,13,15,16 .

Sites of high power potential
Five years (2007-2011) of the control experiment were used to estimate the theoretically available power 40  grid point. The TAP represents the power that would be available if the bulk kinetic energy of the current flowing through the model grid cell were converted into electrical energy. It is calculated using the instantaneous current velocities of the control experiment (see Methods). The five-year mean global TAP shown in Fig. 2 is a first global scale representation of the power that would be available to a 10-km-wide TPP operating in that depth range if the energy extraction had no impact on the flow. WBCs are characterized by TAP values of 100 MW to over 1,200 MW. The TAP of 379 MW in the Florida current (location g7) corresponds to a power density of 1,457 W m −2 , equivalent to ~20 GW for a 100-km-wide and 150-m-deep cross-section, numbers that compare well with the estimates provided by previous studies 4,11,12 . A great value added by our study is the global ocean view of TAP with a spatial resolution of 10 km in a realistic depth range.
The global TAP map derived from the control experiment (Fig. 2) provided a pre-selection of 42 sites of high energy potential where virtual TPPs could be implemented. The final selection of the sites utilized five-year-long time series of the modelled currents because the interannual variability was found to be large in certain areas. Following the approach proposed by Chang et al.  Table 1) were located in the Pacific Ocean along the path of the Kuroshio from the Philippines to Japan and in the Atlantic Ocean along the Yucatan Peninsula and Florida. For example, location gs10 in the Florida Current (see Supplementary Table 1) was considered to be favourable, with a ~35 km distance from the coast (see below in the discussion of the Florida Current), a long-term annual mean current velocity of 1.31 ± 0.14 m s −1 , a minimum velocity in the 20-46 m layer never decreasing below 0.71 m s −1 and a high directional steadiness (the current direction is within the angular sector of 0° to 11°). At the same time, some of the 42 selected sites were not favourable enough: despite a strong annual mean current (above 0.8 m s −1 ), they showed alternating periods of several months with strong declines in current velocity (see Supplementary Note 2).
TPPs, all assumed to be identical, were virtually implemented at 16 of the 42 'most favourable' pre-selected sites (Fig. 2), each of which had a TAP above ~120 MW, and their impacts were simulated in a turbine experiment. The methodology used here to represent TPPs is generally similar to that used for the representation of farms of turbines in ocean models [19][20][21][22] . Specifically, a quadratic drag force added in the momentum equations acts only at the modelled grid points corresponding to the exact locations of the sites. At the same time, because of the much coarser resolution of our global  model compared with strongly localized tidal stream applications, our representation of the drag force is estimated to roughly mimic the integral (that is, over a grid cell) effect of a large pool (several hundred) of turbines on the flow (see Methods). The drag coefficient, C D , which sets the strength of the drag force, depends on the local grid-cell size, with typical values ranging between 2.6 × 10 -3 and 3.2 × 10 −3 . The turbine experiment was performed for year 2008, and the solution was compared with the control experiment for the same year, thereby allowing an assessment of the effect of the TPPs. This effect is quantified by the harnessable power (HP), defined as the maximum resource available for energy extraction after accounting for the effects of the TPPs, and estimated from the power dissipated by the drag force (see Methods). The HP and TAP are related by the reduction coefficient Cp (%) that measures how much the initial expectation of energy production is reduced (see Methods). The effects of virtual TPPs on the ocean currents strongly vary from site to site. The HP represents a reduction of the TAP by 29% to 89% (Fig. 2, insets). Importantly, the large scatter of the reduction coefficients (Fig. 3a) suggests a lack of a robust empirical relationship between the TAP and the HP, implying that the HP cannot be accurately projected from the TAP. The sites of the greatest reduction of TAP are located north of Luzon Island in the Philippines (87%), east of Taiwan (75%) and south of Kyushu Island in Japan (89%). In the Gulf-Stream, the greatest TAP reduction is observed in the northern part of the Florida Current (53%) and west of Cozumel Island in the Caribbean Sea (53% and 64%). Remarkably, the sites where the TAP hints at a high energy potential may be not efficient with respect to the HP because of the considerable reduction in the current speed after implementing the virtual TPPs. A set of sensitivity experiments with a smaller drag force (C D divided by a factor of two, thus varying from 1.3 × 10 −3 to 1.6 × 10 −3 ; Supplementary Note 6) or with a finer grid resolution (1/36°; see Methods) confirmed the robustness of these results (Fig. 3b).

regional effects of TPPs on ocean currents
To understand the wide spread of variation in the reduction factor Cp we considered four sites in the Florida Current (where the turbine experiment shows a significant reduction of the mean currents) and the two sites north of Luzon Island (where the turbine experiment identified the highest reduction in the mean TAP).
In the Florida Current (Fig. 4), the volume transport estimates at the sections located upstream and downstream of the sites are practically unaffected by the virtual TPPs (see Supplementary Table 3), implying that the four TPPs together have a negligible impact on the total transport of the Florida Current. However, the main current stream has shifted offshore upstream and is significantly reduced in the wake of the TPPs implemented at sites gs10, gs15 and gs17 (Fig. 4a-c).
Thus, these TPPs also contribute to the reduction of the TAP at site gs7 located nearly 100 km downstream, which partially explains the strongest reduction of the TAP (53%) at this location. Initially, site gs7, which has an annual mean TAP of 379 MW (Fig. 2), appears to be very favourable for implementing a TPP. Being located ~15 km offshore of Palm Beach (local ocean depth is 208 m), it is well suited for the moorage of subsurface mechanical devices.
The currents in the control experiment (Fig. 5) are quite persistent (1.2 m s −1 to 1.7 m s −1 89% of the time) with a weak seasonal cycle and a remarkably steady direction (0° to 5° 91% of the time). The local drag mimicking a TPP at this site contributes to slowing down the current. However, additional effects are caused by the offshore current shift because of the impact of the three TTPs located upstream. As a result, the mean current velocity decreases by 22% from 1.39 m s −1 to 1.08 m s −1 , the persistence of the current weakens and the directional steadiness decreases (the frequency of occurence of currents in the angular sector 0° to 5° has dropped to 60%). Thus, a proper estimation of the HP in the turbine experiment shows that this site remains favourable for energy extraction, although the efficiency is notably lower than the initial estimates based on the large value of the TAP.   Note that the impacts described above are also seen in simulations with the TPP-induced drag force reduced by half (see Supplementary Note 6) with only a slight reduction of amplitude, and are amplified when the model resolution is increased (from 1/12° to 1/36°, see Methods).
The passage between Luzon and Camiguin in the northern Philippines is another location of great potential and presents the largest TAP among all the selected sites (Fig. 3). In the control experiment (Fig. 6a), a branch of the Kuroshio Current is identified westward of Camiguin. The large velocities (above 1.5 m s −1 ) at the southern tip of the island produce a mean TAP of 855 MW at site ks16 and 1,060 MW at site ks17 (Fig. 2), which is equivalent to the power of a standard nuclear plant (~1,000 MW). Implementation of the two virtual TPPs near Camiguin Island in the turbine experiment has a dramatic effect, shifting the current eastward by ~30 km, where the current joins the main path of the Kuroshio on the eastern side of the island (Fig. 6b,c). This shift results in a considerable reduction of the current velocity in the channel (from 1.70 to 0.89 m s −1 at site ks16; Fig. 7a and Supplementary Table 2b) and a strong decrease in the available power (87%, HP = 110 MW, TAP = 855 MW; Fig. 2 and Fig. 7a). This decrease is one of the strongest decreases in energy efficiency caused by the implementation of virtual TPPs of all the sites studied. A similar estimate for site ks17 leads to a 56% reduction (Fig. 2). Altogether, these two sites provide a HP of 575 MW compared with the 1,915 MW expected from the TAP estimation. A similar current shift leading to a drastic reduction of the TAP is also observed near Cozumel Island in the Caribbean Sea (sites gs1 and gs4; see Supplementary Note 3) and to the south of Kyushu Island (site ks10), where the greatest reduction in the TAP is observed (89%). As in the case of the Gulf-Stream, the impacts described above are qualitatively insensitive to a reduction of the drag force (see Supplementary Note 6) or to the increase in model resolution (see Methods).
Finally, variations caused by eddies or meanders generated by local dynamical instabilities could also affect the available power, thereby decreasing the ability to effectively manage the energy delivered by a TPP. An analysis of this effect is proposed in Supplementary Note 4.

Discussion
We demonstrated that global ocean model simulations resolving ocean mesoscale turbulence can help to pre-select sites for potential energy extraction by immersed large TPPs and also to improve the assessment of the true potential of each site by explicitly simulating the circulation changes caused by power plant feedback on the ocean currents. While the large-scale volume transport was not seriously affected by the TPPs, for certain locations this feedback could drastically reduce (by more than 80%) the available power by locally changing the path of the ocean currents by a few tens of kilometres. Because we simulated only a small number of virtual power plants, this effect will likely be amplified under intense current exploitation scenarios. Importantly, the effects of TPPs on ocean currents could also be non-local, influencing the efficiency of other TPPs located up to 100 km downstream. Each site appears to be peculiar regarding this feedback, which seriously limits the possibility to accurately project power estimates exclusively from observations of currents, and thus implying the importance of developing future modelling approaches.
In this respect, the use of high-resolution models covering domains of several hundred kilometres is challenging and critical for the proper design of TPPs, because such a framework is best suited to account for the non-local effects typically missed in local simulations. Moreover, the design requirements should be dependent on the arrangement of pre-existing or future-planned TPP implementations. Advances in ocean-flow energy extraction technologies expected in the near future will require improvement in the accuracy of TPP efficiency estimates. The emphasis here should be on the representation of the drag induced by a large TPP, for which the use of individual hydrodynamic models of the ocean turbines 6 is a challenge. The large investments required for mastering the key technologies, the optimization of construction, operational and maintenance costs, and the need for new regulations suggests that power plants should have a size of several hundred energy production devices. However, models of large power plants will demand quite high resolution to account for the interactions between individual devices as well as for the compound effect of large power plants on ocean-current systems. Incorporation of such high-resolution models into regional circulation models will definitely require a careful optimization. Large-scale model studies like that presented here will be critical for providing guidance on where to place the domain lateral boundaries and will allow researchers to accurately account for the circulation changes upstream and downstream of the location initially selected for the power plant.

Modelling summary.
We developed a methodology based on model simulations of ocean currents that provides guidance to select global ocean sites where the properties of large wind-driven ocean currents are favourable to the implementation of power plants constructed to include submerged rotors driven by the motion of the water. The methodology is also able to estimate the flow changes induced by the implementation of a large power plant in the current and  to assess the renewable energy that could be recovered and the potential impacts to the environment.
Ocean circulation model. The ORCA12 model used in this study represents an implementation of the NEMO 42 ocean and sea-ice general ocean circulation model on a quasi-isotropic grid of 1/12° horizontal resolution for the world oceans. ORCA12 is being used for global operational forecasts 25 at Mercator-Ocean (https://www.mercator-ocean.fr/), and these operational products are available through the portal of CMEMS. ORCA12 is also used for long-term climatic hindcasts by the DRAKKAR consortium 24 , which conducted the simulations used here. The ORCA12 grid has a horizontal grid resolution that gets finer towards high latitudes (9.25 km at the equator and 7 km at Cape Hatteras, 35° N). The representation of the coastline at that resolution is far from perfect, and caution is required when assessing the accuracy of the available power. The vertical discretization uses 46 vertical levels, and the spacing of these levels increases with depth, similar to many DRAKKAR model configurations 43 44 . Moreover, tides are not considered in this global implementation of NEMO, although they have been explicitly implemented in regional applications of this model 45 .
As shown in Fig. 1, modelled time-mean currents compare rather well with those of the two global observation-based datasets that we chose for reference. The assessment of the modelled currents is complemented by a comparison of the distribution of the speed of instantaneous modelled currents with those of the CMEMS multi-observation dataset over the whole year of 2008 ( Supplementary  Fig. 1). In the Gulf-Stream sector, the model indicates a greater occurrence of weak (<0.1 m s −1 ) and strong (>1.2 m s −1 ) currents, but the model and CMEMS agree very well for speeds greater than 0.8 m s −1 . In the Kuroshio, the consistency between the model and CMEMS is also very good for the high velocities, although the model shows systematically smaller currents than the observation-based analyses for speeds greater than 1.4 m s −1 .
Multifactor criteria for site selection. Chang et al. (2015) 8 defined an index for the selection of potential sites for ocean-current power generation using observations of surface currents in the Kuroshio. This index is based on an optimal combination of four weighted factors (distance to the coast, water depth, mean flow speed and maximum flow speed). Twelve sites were selected as most favourable near Vietnam, Japan, Taiwan and the Philippines, with power densities in the range of 500 to 1,400 W m −2 . Following the same approach but adding currents derived from 20-yr satellite altimetry and from local radar measurements, Tseng et al (2018) 10 established more advanced selection criteria accounting for the current stability.
Our method of searching for the grid points at which to implement virtual TPPs, while conceptually close to the method referred to above 8,10 , follows a more ad hoc application (as we do not weight the various criteria used for selection). To select the most appropriate locations that present a significant ocean-current energy potential, we use five years (2007-2011) of the ORCA12 model outputs of the control experiment and identify regions of the world ocean where wind-driven ocean currents are greater on average than 0.8 m s −1 . Regions that present the highest nearshore currents include the Gulf-Stream and the North Brazil Current in the Atlantic, the Kuroshio Current in the Pacific and the Agulhas Current in the Indian Ocean. The present study focused on the Gulf-Stream and the Kuroshio. The criteria for selecting the model grid points of interest were that the (five-year) mean current velocity should be greater than 0.80 m s −1 in a depth range from 20-46 m, and the distance from the shore should be preferably less than 20 km but always less than 50 km; strict requirements for the steadiness of the current speed and direction were not included, but were qualitatively assessed, and sites with large variability were discarded.
Over 42 locations (longitude, latitude and depth) were selected along the paths of the Gulf-Stream and the Kuroshio (Supplementary Table 1), where the relevant properties of the currents were analysed (amplitude, direction, vertical shear and steadiness) and found appropriate for implementing the virtual TPPs. (Examples of this analysis are shown in Supplementary Note 2).
Control experiment. The initial simulation, referred to as the 'control experiment' , was run for three decades (from 1979 to 2011) by the DRAKKAR consortium. It was driven by three-hourly analyses of atmospheric surface conditions provided by the ERA interim atmospheric reanalysis 46 . Five-day means of all prognostic model variables (including currents) were saved at every grid point every five days over the entire period. At a resolution of 1/12°, the model represents large ocean currents and the associated mesoscale turbulence well. The control experiment, in which the currents were not perturbed by parameterizations mimicking the presence of TPPs in the flow, was used to calculate the TAP (definition below) and identify regions where the available power could be large. The last five years of the simulation (2007-2011) were used to identify these regions because the interannual variations of currents can be large (Supplementary Note 2). The TAP calculated at every model grid point of this control experiment in the 20-46 m depth range is shown in Fig. 2. The WBCs are clearly regions where 100 MW to >1,000 MW is theoretically available through 10-km-wide sections in the depth range from 20-46 m. From that map and basic statistics of the currents over the five-year period, the 16 most favourable locations were selected (circles in the insets of Fig. 2) for the implementation of virtual TPPs.
Turbine experiment. The strategy of our second simulation, referred to as the 'turbine experiment' , was to include a parameterization that mimics the effect of an array of turbines at the 16 previously identified locations. The effect of the array of turbines (that is, the virtual TPP) is parameterized by the introduction of a drag force in the momentum equation applied at the grid points and depths corresponding to the 16 locations. The formula and strength of the drag force are described hereafter. This simulation was run for one year and covered the year of 2008 (its initial conditions were that of 31 December 2007 in the control experiment). The usual storage strategy (five-day means stored every five days) at every model grid point (horizontal and vertical) was applied such that five-day mean estimates of the HP (definition below) could be provided. The comparison of data from the year 2008 in the two simulations represents a method of assessing the impact of the array of turbines on the flow. Comparable simulations have been performed using regional models for limited areas [19][20][21][22] (with more advanced representations of turbines permitted by much higher grid resolution), but the present simulation was carried out on a global model. In addition, because we proceeded with twin simulations, changes seen between the control and turbine experiments remain indicative of changes that could be expected in the real ocean, despite the model flaws identified in a comparison with data from observations. Drag force mimicking TPPs. The mathematical formulation of the drag force that represents the virtual TPPs somewhat follows that used to model large turbine arrays [19][20][21][22] . A TPP is a pool of turbines distributed within a model grid cell. Its impact on the flow is parameterized by a drag force acting on the full grid cell and over the depth that corresponds to the total area swept by the rotors of the turbines. At every selected grid point (longitude, latitude and depth) where the effect of a power plant should be represented, a drag force is added to the model momentum equations with the following form: C D is a (dimensionless) drag coefficient to be determined, ρ is the density of seawater, U is the current velocity and L 2 is the horizontal surface of the model grid cell. The challenge is to estimate the drag coefficient C D so that it realistically mimics the integral effect of turbines on the flow. This complex process is described in detail in Supplementary Note 8. We provide here a comprehensive summary.
Our estimation of C D is based on the one-dimensional theory for an ideal wind turbine 40 , which states that the thrust of a single rotor is equal to the drag force exerted by the turbine on the flow in the stream-wise direction. The thrust for a single rotor T i is calculated with the following formula: [19][20][21][22] where U i is the current speed in the turbine, C T is the thrust coefficient and a is the area swept by the rotor of the turbine. We consider that the overall thrust of a pool containing a large number, N, of identical turbines operating in the same flow regime (U i = U) is the sum of the thrusts of all rotors. The mathematical expression of the overall thrust is as follows: In comparable studies [19][20][21][22] the drag force F (equation (1)) is often set equal to the overall thrust T (equation (3)). The underlying hypothesis is a hypothesis of perfect extracting devices: devices exert no drag on the flow other than the thrust of the rotors. In practice, extracting devices will be imperfect and the total drag force will be the addition of the thrust force exerted by the rotors (given by equation (3)) and an additional drag force resulting from all the other effects resisting the flow (for example, floatability, frames carrying rotors, wings, cables and so on). We assume here that the additional drag force is quadratic and can be expressed as a proportion (α) of the overall thrust of the rotors. The drag force and the thrust are thus related by F = (1 + α)T: Equation (4) yields the following for the drag coefficient: Equation (5) establishes the relationship that links the drag coefficient and the number of rotors operating in the grid cell. Because the turbines are distributed within a full grid cell, we assume that the total area represented by the N rotors is equal to the face of the grid cell, that is, where Δz is the vertical extent of the grid-cells where the drag force is applied. For simplicity we use C T = 1 (which means that all the energy from the flow is transferred to the rotors and we do not carry this coefficient in the following). This yields, for the value of the drag coefficient C D used in this study, the following: The turbine experiment described in the paper uses α = 1. It corresponds to a case of imperfect extracting devices, and the drag coefficient is C D = Δz/L. In this experiment, it is hypothesized that only half of the dissipated energy contributes to the HP (see also Supplementary Note 8). Typical values of C D range between 2.6 × 10 −3 and 3.4 × 10 −3 at the latitudes (0° to 35° N) where most TPPs have been implemented.
A sensitivity experiment using α = 0 has been carried out. It corresponds to the case of perfect turbines and a drag coefficient C D = ½Δz/L (see also Supplementary Note 8). The results of this experiment are used to assess the robustness of our results to the uncertainty related to the arbitrariness in the choice of α. These results are reported in Fig. 3b and Supplementary Note 6. Since energy extracting devices are designed to reduce the additional drag, these two cases can be seen as the upper and lower limits for the drag coefficient, with realistic cases lying in between.
Considering that rotors have a diameter of Δz (to cover the relevant depth range) and thus sweep an area a = π/4(Δz) 2 , the number of turbines can be estimated using equation (6). Using a reference value L = 10 4 m, which is close to the model grid cell at the equator, and a value of Δz = 26 m, which represents the slice of ocean in the depth range from 20-46 m where the virtual turbines are operating, yields a farm of nearly 500 rotors. A conceptual TPP consistent with the calculations above, and thus consistent with our model set-up, could be achieved with approximately 315 C-plane turbines 7,13 (each with two rotors of 23 m in diameter) as described in Supplementary Note 5. Such a number of devices is quite consistent with previous conceptual studies 6,13,47 . Note that the area occupied by 315 C-plane turbines represents less than 3% of the grid-cell area, which allows us to assume there are no interactions between individual turbines within the pool, in the calculation of the HP.
Definition of TAP and HP. We define TAP as follows. Great ocean currents are several tens of kilometres wide, and assessments of their potential as sources of renewable energy are generally based on an estimation of the total hydrokinetic energy flux through wide sections of the current 15 . For an ocean current flowing with a velocity V 0 through a vertical section of area A, this energy flux, also called the TAP, is calculated as follows: 2,15,40,41 where ρ is the density of seawater. With this definition, the TAP represents the maximum power that would be available if the totality of the kinetic energy of the current flowing through the section was converted to electrical energy. Because the presence of turbines will significantly modify the flow, the TAP, which is calculated using the undisturbed flow, systematically overestimates the available power. The HP (see Supplementary Note 8 for a detailed calculation) is calculated as the power that could be extracted by the rotors in the case where the rotors could convert 100% of the incident energy into electrical power under the hypothesis of a maximum thrust (that is, the rotors have 100% efficiency). Since such a level of efficiency is possible only in theory, the HP is referred to as the maximum power theoretically available for energy extraction. With this definition, it can be considered as an analogue to the TAP, but for the case that accounts for TPPs in the flow (both are energy fluxes through the same area). Calculating the HP requires quantitative knowledge of the turbine-induced flow changes, and such knowledge is obtained here via the turbine experiment. In the present experimental setting, the HP is directly calculated from the power dissipated by the drag force, P DISS = FU: Consistent with the value of the drag coefficient given in equation (5), another mathematical expression of the HP can be directly obtained from the thrust (equation (3)) of the N rotors (each of area a): An estimate of the extracted power can be obtained from the HP by applying a coefficient of efficiency (for example 30 to 40% for a marine turbine 19 ).

Calculations of TAP and HP.
In the present study, we calculated the TAP at every model grid point for every five days of the period 2007-2011 using equation (7) with the following parameters: seawater density of ρ = 1,024 kg m −3 ; V 0 that is the five-day mean current speed in the control experiment; and A that is the area of a 10-km-long cross-section centred at the model grid points in the depth range 20-46 m (that is, A = L 0 Δz with L 0 = 10 4 m and Δz = 26 m). The reason for choosing the same value for L 0 at all grid points (rather than using, for example, the varying model grid size L) was to have TAP estimates that could be compared independently of the location. Figure 2 shows the five-year mean (2007-2011) TAP for the global ocean. We then calculated the HP at identical 10-km-wide sections centred at the 16 TPP sites for every five days of the year 2008. The HP was calculated using equation (9) with L 2 = LL 0 for consistency of the use of L 0 in the calculation of the TAP. The above calculations allow for a strict comparison between the TAP and the HP at the 16 TPP sites for the year 2008, with differences being quantified via the reduction coefficient, Cp (%) = (1 -HP/TAP) × 100. A large value of Cp will indicate that the power available for extraction is significantly less than initially expected. The HP is much more relevant than the TAP for assessing renewable energy sources represented by great ocean currents. However, this parameter still overestimates the power that could be delivered by a pool of turbines because it does not account for turbine efficiency. Only a portion of this resource can be harnessed practically.
Sensitivity to model grid resolution. We investigated the robustness of our results via an increase of the grid resolution from 1/12° to 1/36° (increasing the number of grid points by a factor of nine). This investigation was performed in a regional model of the western North Atlantic that includes the Florida Current and a large part of the Caribbean Sea (Supplementary Note 7). A control experiment (no turbine) and a turbine experiment with virtual TPPs were implemented in the Florida Current and the Cozumel Channel at the same locations and the exact same size (that is, nine grid points in the 1/36° model) as in the 1/12° experiment (gs1 and 4 in the Cozumel Channel and gs7, 10, 15 and 17 in the Florida Current; Fig. 2). The drag coefficient was as in the 1/12° experiment. As expected, the currents were stronger at 1/36° resolution (Supplementary Note 7). Changes induced by the TPPs (Supplementary Fig. 2) exhibited very similar patterns to those seen in the 1/12° experiments, with a significant reduction of the current at and in the wake of each individual TPP, and a general reduction/increase in the nearshore/offshore part of the current (indicating a current shift) that was seen far upstream from the location of the TPPs. Because the currents are stronger in the higher resolution model, the amplitude of the changes are also greater (consistent with the quadratic form of the drag force). The scatter plot of the HP as a function of TAP (Fig. 3b) shows that if the TAP is generally greater in the 1/36° experiment due to the greater amplitude of the currents, the reduction of the HP is always greater, although it remains proportionally comparable (similar Cp values), which is consistent with the pattern of the changes in the currents.

Data availability
The current velocity dataset that supports the findings of this study is available at https://doi.org/10.5281/zenodo.3631372. This dataset comprises the five-day zonal and meridional velocity components, from the surface down to 100 m depth, for the regions of the Gulf-Stream and Kuroshio that were considered in this study, for the 1/12° control and turbine simulations used in this study. These components are available in a NetCDF format that contains the appropriate metadata. Documentation that indicates the model coordinates of the points where the drag force was implemented is also provided. Source Data is provided for Figs. 1-7. The five years (2017-2011) of the 1/12° global model outputs of the control experiment and of the 1/36° control and turbine simulations are available only on request because of the large amount of storage space required. Contact: jean-marc. molines@univ-grenoble-alpes.fr.