A Rogeau

¦¦ R Girard
email: robin.girard@mines-paristech.fr

G Kariniotakis
email: georges.kariniotakis@mines-paristech.fr

A generic GIS-based method for small Pumped Hydro Energy Storage (PHES) potential evaluation at large scale

Keywords: Pumped hydro energy storage, GIS, Potential, Sobol analysis, Renewable energies, France

The increasing share of weather-dependent renewable energies in power systems creates a need for energy storage technologies to reduce the impacts of variable production. The most mature technology to store energy on the grid remains Pumped Hydro Energy Storage (PHES). The potential of high-energy sites has already been assessed in Europe by the EU JRC, considering mostly dams and reservoirs from global European databases which include only massive water bodies. This paper focuses on estimating the potential for small-PHES, proven to have lower environmental impact and an positive impact on grid balance and reliability. A generic method is designed, able to evaluate a global PHES storage capacity at large scale. It considers both existing lakes and natural depressions suitable to be filled for PHES purposes. The volume of filled lakes is estimated using the surrounding topography. The method is organized so that the "heavy" calculations, i.e. sink detection, volume evaluation, constraints verification etc. are run only once. Consequently, the actual potential estimation phase only includes fast calculations and can be integrated in a loop for carrying out a sensitivity analysis. The proposed method is then applied considering France as a test case. Suitable environmental, land-use and structural constraints are applied to eliminate irrelevant sites. The analysis leads to an estimated value of the small-PHES potential in France, which ranges from 14 GWh when only existing lakes are considered to 33 GWh when lakes and depressions are considered. These estimations represent respectively 8% and 18% of the current hydro storage capacity in France. Thanks to a global sensitivity analysis, factors like the maximum distance between lakes, the maximum altitude of the sites, and the distance to the electrical grid are shown to have the most influence on the global evaluated potential, which is further sensitized. Lastly, another application is suggested that makes it possible to select the connections to be built first within a restricted area, based on a cost-per-energy-like approach. It uses the connections between reservoirs detected at large geographical scale.

Introduction

Over the last few decades, renewable energy sources (RES) have continuously increased their share in the world energy market. In fact, worldwide RES installed capacity went from 800 GW in 2004 -mostly from hydropower sources (715 GW) -to almost 1850 GW by the end of 2015 [START_REF]Renewables 2016 Global Status Report[END_REF] -1064 GW from hydropower. During this period, the capacities of wind power (WP) and solar photovoltaic (PV) installations have multiplied by thirteen, going from 50 to 660 GW [START_REF]Renewables 2016 Global Status Report[END_REF]. Both PV and WP production are highly variable and not entirely dispatchable due to their inherent dependence on weather conditions. To handle RES production fluctuations, different solutions exist, such as demand management, the combination of RES plants with other means of conventional generation, and also energy storage. Several storage technologies exist [START_REF] Zach | Report summarizing the current Status, Role and Costs of Energy Storage Technologies, Tech. Rep., stoRE -Facilitating energy storage to allow high penetration of intermittent renewable energy[END_REF], including electrochemical devices such as batteries or fuel cells and mechanical solutions such as pumped hydro energy storage (PHES).

Water lifting devices have been invented, used and improved by humans for thousands of years [START_REF] Yannopoulos | Evolution of Water Lifting Devices (Pumps) over the Centuries Worldwide[END_REF]. Following the evolution of these technologies, PHES is now one of the most mature and cost-effective solutions for absorbing RES power fluctuations. Rehman et al. [START_REF] Rehman | Pumped hydro energy storage system: A technological review[END_REF] and Uniyal [START_REF]Old concept, new opportunities. Pump-storage[END_REF] estimated the global pumped storage installed capacity as slightly more than 100 GW in 2010.

Large PHES facilities are usually employed as grid-scale solutions to ensure that electricity production matches demand at any time, especially thanks to their fast response. PHES offers considerable storage capacities at affordable prices [START_REF] Manwaring | Challenges and Opportunities For New Pumped Storage Development[END_REF].

In Europe, high-potential sites for large PHES have been assessed by an EU JRC study [START_REF] Gimeno-Gutirrez | Assessment of the European potential for pumped hydropower energy storage based on two existing reservoirs[END_REF], which considered artificial reservoirs larger than 100,000 m 3 separated by at least 150 m of head. This resulted in a total European feasible potential varying from 36 to 4000 GWh for reservoirs at a distance of between 3 and 20 km. For France, this potential reaches 4 GWh if the distance between lakes is limited to 3 km and rises to 500 GWh at 20 km. Large PHES schemes were initially considered as the only viable type of pumped storage due to elevated inversion costs and scale effects [START_REF] Manwaring | Challenges and Opportunities For New Pumped Storage Development[END_REF]. However, small pumped hydro energy storage (small-PHES) has recently been studied more deeply. Small-PHES has a lower environmental impact in terms of CO 2 equivalent [START_REF] Zhang | Review on the externalities of hydropower: A comparison between large and small hydropower projects in Tibet based on the CO2 equivalent[END_REF]. Large-PHES has a national and continental role to play, while the role of small-PHES is local and regional. The use of local small-PHES utilization can reduce peak demand and thus avoid or defer electrical grid reinforcement [START_REF] Nourai | Installation of the First Distributed Energy Storage System (DESS) at American Electric Power (AEP)[END_REF][START_REF] Denholm | The Role of Energy Storage with Renewable Electricity Generation[END_REF]. Small-PHES can also mitigate voltage issues in medium and low voltage lines, due to the injection of decentralized production [START_REF] Gabathuler | Le pompageturbinage a petite echelle pour le stockage local d energie[END_REF]. In brief, small-PHES can contribute to the decentralized integration of intermittent RES. In this perspective, PHES systems have been miniaturized and are even starting to be used in cities, installed on top of buildings [START_REF] De Oliveira E Silva | Pumped hydro energy storage in buildings[END_REF]. The present paper thus focuses on evaluating the potential of sites connecting close reservoirs of limited volume that are suitable for installing small-PHES plants. This method does not aim to reduce the PHES costs in itself. It aims to find spots requiring few civil work and thus exploitable with little expenditure. For this purpose, a new computationally efficient methodology applicable at large scale (i.e. national) is proposed. The method employs Geographic Information Systems (GIS) to detect reservoirs, associate those that could host a small-PHES plant, and finally apply the different constraints to derive a feasible potential. One disadvantage of small pumped hydro energy storage is the investment cost, given the low storage capacity. Small PHES indeed suffers from scale effects [START_REF] Crettenand | The Facilitation of Mini and Small Hydropower in Switzerland[END_REF].The method developed considers natural opportunities, such as existing lakes and natural terrain depressions, which could be operated at a limited cost. The proposed methodology is applied to the case of France to evaluate small-PHES potential at national scale.

For both the conventional an pumped hydropower systems, the location of the power plant is important, as the extracted energy and power depend, among other factors, on the gross head, water flow and volume of water available. Selecting optimal spots that respect these parameters thus appears to be crucial. For this purpose, automatic numeric methods have been developed to select the best spots, taking advantage of the emergence of GIS. Thereby, Larentis et al. [START_REF] Larentis | Gis-based procedures for hydropower potential spotting[END_REF] present a method to detect the optimal location in a river basin to install a run-of-river hydropower plant. They use GIS to detect where to install a dam and build the powerhouse. Kusre et al. [START_REF] Kusre | Assessment of hydropower potential using {GIS} and hydrological modeling technique in Kopili River basin in Assam (India)[END_REF] evaluate suitable places in an Indian valley, using hydrological tools and soil, land use and weather data. Both focus on power generation, without storage and optimize the location within restricted areas, mostly valleys.

Other studies concern the PHES scheme, and are carried out in limited areas to evaluate either the potential of a zone, or a particular spot to install a PHES plant. Ahmadi and Shamsai [START_REF] Ahmadi | Preliminary Site Selection of Pumped Storage Hydropower Plants -A GIS-based approach[END_REF] focus on a single massive lake, looking for the best location to construct the corresponding reservoir using criteria ranking through raster analysis. This is a localized analysis, with an emphasis on very detailed geological or environmental data, which makes application of this method at a larger scale difficult. Mailler et al. [START_REF] Mailler | Pumped-storage with small hydropower -an option?[END_REF] focus on small hydropower schemes, taking advantage of existing reservoirs in Switzerland, such as artificial snow reserves. GIS data are not used and the study is conducted in a very restricted area. Kucukali [START_REF] Kucukali | Finding the most suitable existing hydropower reservoirs for the development of pumped-storage schemes: An integrated approach[END_REF] estimates the best location for conversion of existing hydropower plant to PHES schemes comparing a few sites in Turkey, with ranking based on geological, social and environmental criteria. Garcia [START_REF] Garcia | Testing a GIS-based methodology for optimal location of Pumped Storage power plants in Norway[END_REF] focuses on one very restricted region, looking for PHES opportunities within a group of many lakes. Several potential connections were identified and the selection was then optimized depending on criteria such as stored energy potential, and the cost of stored energy potential. A common characteristic of these four studies is that they are conducted on a small scale and based on very detailed datasets (geological composition, ground slopes, urban planning, etc.), which are difficult to obtain for a large area.

Other works are based on the detection of high power (and energy) capacity PHES sites. Connolly et al. [START_REF] Connolly | Development of a computer program to locate potential sites for pumped hydroelectric energy storage[END_REF] proposes a method that exploits a Digital Elevation Model (DEM) to identify restructurable hills and build reservoirs in the south of Ireland. This method implies considerable investment costs, which are only profitable if they correspond to high-power projects. Four major works study energy storage possibilities at a larger scale [START_REF] Douglas | Assessment of Opportunities for New United States Pumped Storage Hydroelectric Plants Using Existing Water Features as Auxiliary Reservoirs[END_REF][START_REF] Fitzgerald | A GIS-based model to calculate the potential for transforming conventional hydropower schemes and non-hydro reservoirs to pumped hydropower schemes[END_REF][START_REF] Gimeno-Gutirrez | Assessment of the European potential for pumped hydropower energy storage based on two existing reservoirs[END_REF][START_REF] Lu | A GIS-based assessment of Tibet's potential for pumped hydropower energy storage[END_REF]. The first one is focused on taking advantage of existing dams, but also includes a section on matching existing lakes in the United States. It is oriented to large PHES schemes, but only focuses on specific sites and does not offer a global potential evaluation in terms of energy or power. Similarly, the second study only searches for high-power locations throughout Turkey, with high head or volume characteristics (150 m and more than 10 6 m 3 ). It also includes a search for flat locations to dig new reservoirs. Likewise, the third document is a global study of potential PHES sites. However, it includes a European assessment of energy available per country, depending on a maximum authorized distance between reservoirs of 3 to 20 km. Only heads greater than 150 m are considered. Lu and Wang [START_REF] Lu | A GIS-based assessment of Tibet's potential for pumped hydropower energy storage[END_REF] present an evaluation of Tibet's potential for PHES. It considers connection with very high heads (greater than 500 m) and large reservoir areas (greater than 60,000 m 2 ). All the above works exclude connections with low head and low potential, which might be beneficial for small-pumped storage due to restricted civil works, as we aim to evaluate here.

Pauwels et al. [START_REF] Pauwels | Assessment of the potential of Pumped Hydro Energy Storage[END_REF] conducted a study at the scale of a French county (about 4000 km 2 ). It matches natural terrain depressions within a wide constrained zone. Existing lakes are not taken into account. Moreover, no optimization through interconnected lakes is performed, thus leading to an overestimation of energy potential. The evaluation by Gimeno-Gutirrez and Lacal-Arntegui [START_REF] Gimeno-Gutirrez | Assessment of the European potential for pumped hydropower energy storage based on two existing reservoirs[END_REF] presents a substantial shortfall in the potential reservoir dataset, by considering only existing artificial reservoirs.

Considering only a share of potential reservoirs can result in a considerable shortfall in evaluations. Moreover, in order to avoid matching one lake with several others, it is necessary to select the most effective connections, which results in more precise estimations of the overall potential. In Fitzgerald et al. [START_REF] Fitzgerald | A GIS-based model to calculate the potential for transforming conventional hydropower schemes and non-hydro reservoirs to pumped hydropower schemes[END_REF], a sensitivity analysis is carried out featuring the buffer size, minimum head, and maximum slope of flat terrains compatible with PHES installations. Nonetheless, each criterion is treated separately, which makes it difficult to understand the global behaviour of the output. Similarly, Garcia [START_REF] Garcia | Testing a GIS-based methodology for optimal location of Pumped Storage power plants in Norway[END_REF] presents a sensitivity analysis of the evaluated potential, drawing an abacus representing the evolution of the detected energy depending on storage duration and distance between reservoirs. It clearly identifies the latter factor as more decisive. However, the analysis is only carried out for two input parameters. These sensitivity analyses remain focused on a limited number of parameters and are not globally applied to all of the parameters of the corresponding evaluation method.

The main contribution of this paper is to propose a generic method able to evaluate the potential of small-PHES over a large geographical zone (several hundreds or thousands of square kilometres. To ensure that a significant share of potential reservoirs are considered, the method includes the dataset of both existing water bodies and natural terrain depressions . It can be applied at a large scale thanks to the global datasets selected. Connections compatible with small-PHES purposes (e.g. not too massive) are evaluated, thus opening up a new perspective in comparison with former studies. A pairing of potential reservoirs is carried out authorizing a maximum of one connection per reservoir, to increase the accuracy of the evaluation. Moreover, the method includes an optimization step to select the most appropriate connections when a lake is concerned by various matchings. A global Sensitivity Analysis is incorporated into the method and is feasible thanks to the low execution time of the evaluation process. As a second contribution, the developed method is applied to the whole of France (about 675,000 km 2 ). This paper presents the overall estimated potential but also the spatial distribution of detected connections. The proposed global SA is applied to the same case study. Additionally, a potential avenue for further development is to use the obtained connections database as an alternative way of making a best site selection analysis.

Method description

The proposed methodology consists in identifying existing reservoirs -either natural water bodies (referred to as lakes) or sinks in the topography (referred to as depressions) -and searching for another compatibly located potential reservoir nearby in order to reference them as potential small-PHES sites. The method is divided into four main parts.

Main model assumptions

The proposed approach considers only closed-loop systems, excluding rivers and the sea as a lower reservoir. Only connections between two inland water-bodies are considered.

Reservoir detection

In a simple closed-loop hydro storage scheme, two water bodies are used to pump water up and down, in order to store and generate electricity. The primary type of reservoir suitable for this purpose is existing lakes. Given that they are already filled with water and hermetic, they constitute the easiest resource for installing a PHES plant.

Existing lakes

Databases containing water body information on a large scale -such as the Global Reservoir and Dam (GRanD), the European catchments and Rivers network system (Ecrins) or the DAM Positioning (DAMPOS) databases -focus on very large reservoirs. For smaller water bodies, national institutes are better suited, such as NGI in Belgium, SNIG in Portugal and IGN in France, all of which provide appropriate national data-bases.

These databases, which are often supported by GIS, generally include the water bodies, their altitude, regime and nature, but no information on volume. However, this information is essential for hydro storage because it defines the energy available in a connection. Some studies have been carried out to estimate reservoir volume using surface area as the only input. Liebe et al. [START_REF] Liebe | Estimation of small reservoir storage capacities in a semi-arid environment: A case study in the Upper East Region of Ghana[END_REF] and Sawunyama et al. [START_REF] Sawunyama | Estimation of small reservoir storage capacities in Limpopo River Basin using geographical information systems (GIS) and remotely sensed surface areas: Case of Mzingwane catchment[END_REF] developed methods estimating the volume of a reservoir as a function of its surface area in Ghana and Zimbabwe. Nonetheless, both works concerned small reservoirs and required manual categorization of lakes. In this work, we evaluate the volume of lakes using the method proposed in Hollister et al. [START_REF] Hollister | Predicting Maximum Lake Depth from Surrounding Topography[END_REF], by which we can estimate a natural WB's volume from its surrounding topography. This approach is computationally efficient for the large-scale study we are aiming at, compared to alternatives that would give an exact estimation of the volume. Altitude information can be found using global Digital Elevation Models (DEM) and especially the SRTM model, available worldwide and free of charge. National institutes also often possess their own national DEM. Typically, the spatial resolution of such models is 25 m.

Natural depressions

Beyond existing lakes, the second objective of this step is to identify potential reservoirs by evaluating natural terrain depressions that might be flooded. The same DEM as in 2.1.1 can be used to carry out this step.

First, all depressions in the DEM must be detected. The methodology described by Planchon and Darboux [START_REF] Planchon | A fast, simple and versatile algorithm to fill the depressions of digital elevation models[END_REF] allows us to obtain all terrain depressions, of any size and on any territory. To ensure that very large valleys are not included in the research, we can set a maximum surface area. A threshold of 1 square kilometre seems reasonable for such small-PHES schemes. Depressions exceeding this limit are then scanned downwards to identify smaller internal depressions. Volume calculation is straightforward. The smaller sinks can be eliminated using a minimum energy criterion (20 kWh/m for example, corresponding to Volume 8000 m 3 ).

Connection creation 2.2.1. Reservoir matching

One of the objectives of the proposed method is its computational efficiency. To achieve a fast simulation, the matching detection is achieved following two major aspects:

First, a database is created, containing all potential small-PHES connections. To do so, two settings are defined, regarding the minimum head and the maximum distance between lakes (See Table 1). These are chosen at not very restrictive values so that a maximum amount of pairs are included in the dataset. Concretely, a reservoir will be paired with every other reservoir presenting a head of minimum 10 m within a representative diameter of 3 500 m. For each pair of lakes that satisfies the above criteria, the following parameters are estimated using the available terrain databases and then included

Maximum distance between reservoirs

Minimum head between reservoirs 3 500 m 10 m in the database with all of the potential small-PHES connections:

• Distance between reservoirs (m)

• Gross head between reservoirs (m)

• Energy of the connection (kWh)

• Volume of the smallest reservoir -available to transfer (m 3 )

• Depth of the shallowest reservoir (m)

• Altitude of the upper reservoir (m)

• Shortest distance to the electrical grid (m)

The energy of the connection is calculated as follows:

E ηρgHV (1)
where η is the efficiency of the pumped storage plant (chosen as 0.7 in our case, the typical efficiency of small-PHES [START_REF] Papaefthimiou | Operating policies for wind-pumped storage hybrid power stations in island grids[END_REF]), ρ is the water density (kg/m 3 ) and g the gravity acceleration (m/s 2 ). H is the altitude difference between the bottom of the upper reservoir and the surface of the lower reservoir, called the gross head of the connection (m). V is the available water volume to transfer (m 3 ), calculated as 90% of the volume of the smaller of the two involved reservoir.

The shortest distance to the electrical grid is the distance from the lower reservoir to the high-or medium-voltage electrical grid. This electric network database, generally produced by the national geographic institute or the grid owner itself, has to be imported for the study zone. Installing a power plant very far from the existing network would make it difficult to carry out the project, as it adds costs and involves civil works.

Attaching acceptability constraints information

Some of the potential reservoirs detected, i.e. existing lakes and natural terrain depressions, might be ineligible for PHES application. Many zones are preserved from human presence and constructions to protect fauna and flora. Many kinds of preservation zones exist but they do not all ban energyrelated constructions. The World Database Protection Area (WDPA) 1 provides an inventory of these protected zones worldwide. These zones can be regulated by regional, national or international instances, and involve different restrictions.

Additional constraints that can lead to the exclusion of some locations relate to land use. These may include urban areas, leisure areas (i.e. navigation, swimming, etc.) or professional activity zones (i.e. flooding arable land might make a project difficult to complete). Available land-use datasets differ depending on the study zones, but a good option for Europe is the Corine Land Cover 2 . The 2015 publication relating to 2012 land cover is available for the whole of Europe and describes land use, dividing activities into several sections. Whatever the dataset, the following five categories of land use can be considered as relevant:

• Urban zones: Most of the land is covered by a waterproofed surface. Most of it features buildings or roads.

• Leisure activities: The land is used for camping, sports and other leisure activities.

• Arable land: Land capable of being ploughed and used to grow crops.

• Permanent farming: Permanent farming land, including crops such as orchards, vineyards etc.

• Pasture: Mostly land used to feed livestock.

Moreover, one of the main reasons for eliminating a natural depression from potential PHES sites is when it features human constructions. The detection of any building within natural depressions must be verified using a building database (obtained from national institutes for example). Each of the constraints cited above is applied to the whole set of small-PHES compatible connections, in order to verify whether the connected reservoirs are present in a sensitive area. If so, the connection involving these reservoirs could be discarded in an upcoming process. Further evaluation of the storage energy available across the study zone is then carried out, depending on which of these acceptability constraints are considered to be incapacitating. Secondly, for each detected connection, the inclusion of the future pipe, represented by a line between the two matched reservoirs, is verified in line with the constraints concerning urban zones. Consequently, the three following pieces of acceptability information are attached to the detected connection pairs:

• Inclusion of matched reservoirs in land-use zones, environmental protection areas (vector with binary values 0 or 1)

• Presence of buildings in any reservoir (0,1 vector)

• Inclusion of the pipe in a urban zone (0,1 vector)

Finally, the detected connection pairs are divided into three main groups, defined as shown in Table 2.

2 http://land.copernicus.eu/pan-european/ corine-land-cover/clc-2012/view

Constraints application

Among the numerous connection pairings detected, some have to be discarded using the parametrized constraint data generated in 2.2.2 and before the next step of optimization. Initially, a number of pairings is eliminated after being checked against the land-use constraints. As a default configuration, restrictions on such constraints are suggested in Table 3. Likely No-go Pasture Suitable environmental preservation rankings must be defined depending on the preservation area dataset used for the study zone. Generally, international preservation zones and natural reserves can be considered as No-go zones.

Connection pairings with configurations that do not correspond to a small-PHES scheme also have to be eliminated. For this purpose, nine main criteria are selected; their respective thresholds are defined in Table 4. We choose a minimum

Optimization of inter-connected pairings

Among the selected connection pairings, some may include reservoirs featuring in more than one pairing. We consider that the small-PHES scheme can involve only two lakes (one reservoir can not be linked to several others). An optimization is done as explained below to obtain non-interconnected pairings.

The problem to solve can be seen as a multidimensional knapsack problem [START_REF] Garey | Computers and Intractability; A Guide to the Theory of NP-Completeness[END_REF], defined as follows:

max Z n j1 c j x j (2) subject to n j1 a ij x j ¤ 1, di v1, mw (3)
x j t0, 1u, dj v1, nw (4)
where j is the connection number and i is the lake number.

x j is the binary representing whether connection j must be built or not and c j is the "benefit" of the connection. Equation (3) represents the fact that one lake i can be paired only once, where a ij is a matrix describing whether connection j concerns lake i. Equation ( 4) describes x j as a binary. As a result, connections are selected so that each reservoir is connected no more than once. In some cases, more than one connection can remain within a group of interconnected reservoirs (if there are more than 3 reservoirs).

Benefit definition

Different benefit vectors can be defined before running the optimization tool:

• Energy, as calculated in [START_REF]Renewables 2016 Global Status Report[END_REF]. Maximizing the energy would provide the biggest storage capacity in the study zone.

• Power, calculated as P ρgHQη, where Q is the flow rate in the pipe (m 3 /s). This strategy implies defining the maximum water flow for each connection and favours the pairs with high head, which provide bigger power.

• Cost-energy ratio, calculated as:

CoE I E ( 5
)
where I is the investment costs and E the energy available in the connection.

This cost, expressed in e/KWh, allows us to choose more expensive connections if they offer a considerable energy advantage. This strategy also creates the possibility of eliminating the least effective connections (below a defined threshold).

• Other functions, such as levelized cost of energy, calculated as follows:

LCOE °n t1 It Mt p1 rq t °n t1 Et p1 rq t (6)
where I t is the investments costs at year t, M t is the operation and maintenance costs at year t, E t is the benefit (here energy or power) in the year t, and r is the discount rate.

• Other benefit functions can be used, such as cost of power.

In case of using energy or power as a benefit vector, the above optimization problem can be adapted to give priority to solutions including existing lakes, which are generally less expensive to fill than depressions. This translates into selecting T1 connections first, then T2, and ultimately the remaining T3 connections.

2.5. Sensitivity analysis 2.5.1. Impact of the defined exclusion zones In section 2.3, connection pairings are eliminated from the evaluation because of their inclusion in constraint zones defined by the user as No-go. The binary state "No-go/Suitable" is generally not so clear in reality, and the impact of this uncertainty can be significant on a large study area. For this reason, it is recommended to identify the sensitive constraint areas that are most likely to be between the No-go and suitable status, and evaluate the variation of the energy storage capacity whether these specific constraints are accepted or not. A quantification of the output variation around a default state should be achieved, constraint by constraint.

Sensitivity analysis method

As described in Padey et al. [START_REF] Padey | From LCAs to Simplified Models: A Generic Methodology Applied to Wind Power Electricity[END_REF], a sensitivity analysis (SA) is an approach for investigating results variability as a function of the input parameters. Most SAs only consider variation around a nominal value of one parameter at a time, or from a worst case to best case scenario. This results in misestimations of the influence of each parameter on the function output. The Global Sensitivity Analysis (GSA) proposed by Sobol [START_REF] Sobol | Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates[END_REF] overcomes these limitations. Moreover, this model does not require linear modelling assumptions.

The Sobol indices reflect the shares of variance in the output that are explained by variations in the input parameters. The higher the value of a Sobol index, the higher the impact of the corresponding parameter on the model output. The overall methodology is derived from Sobol [START_REF] Sobol | Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates[END_REF] and is divided into 4 sub-steps:

• Stochastic generation of n scenarios,

• Computation of the impacts, applying the method with the n scenarios as inputs,

• Estimation of the Sobol indices for each parameter,

• Selection of the key parameters explaining most of the variability.

Sobol indices definition

To quantify the sensitivity of the resource estimation with respect to the various inputs considered in our model, we choose to use Sobol indices [START_REF] Sobol | Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates[END_REF]. These are global sensitivity indices designed to estimate the relative share of output variance determined by the considered inputs or their interactions. Sobol indices are very popular because their computation does not rely on any (e.g. linear) hypothesis about the underlying model function f that links input X i , i φ : t1, . . . , pu to output Y f pX φ q. However, they assume that inputs are uncorrelated. Then, the following generic decomposition also has uncorrelated terms: Various summary statistics based on Sobol indices may then be proposed. In Homma and Saltelli [START_REF] Homma | Importance measures in global sensitivity analysis of nonlinear models[END_REF], the total sensitivity of a system S Ti to a given parameter X i is defined by adding up all indices that involve parameter X i , i.e. S Ti °I,iI S I Extending that idea, it may be interesting to evaluate the total sensitivity of a system S T L with some general aspect (e.g. technological characteristics, relief constraints, etc.) whose definition potentially involves several parameters L φ, i.e. V ). This can be done through Monte Carlo simulations using two independent, equally sized, random samples of inputs: pX j,1 φ q and pX j,2 φ q, where j 1, . . . , n [START_REF] Homma | Importance measures in global sensitivity analysis of nonlinear models[END_REF]:

Y EpY q Iφ f I pX I q with f I pX I q ErY |X I s ¡ I 1 I f I 1 pX I 1 q (7)
f0 ÊpY q 1 n n j1 f pX j,p1,2q φ q (11) V 1 n n j1 f 2 pX j,p1,2q φ q ¡ f 2 0 (12) v arpErY |X I sq 1 n n j1 f pX j,1 φ qfpX j,1 I , X j,2 Ī q ¡ f 2 0 ( 13
)
where Ī φzI. Estimating Sobol indices requires pK 1qn evaluations of f , K being the number of considered indices.

Application of the method to France

The method described in Section 2 is applied to a whole country, i.e. France (Area: 675,000 km 2 ). This section presents the software employed to implement the method, the datasets used for this purpose and the results of the evaluation. A very local study case is presented to demonstrate the efficiency of the optimization tool.

Dataset choice and connection detection

Existing water body data is provided by the French National Geographic Institute (IGN) TOPO database. This contains all WBs over 20 metres long detected by satellite. It represents about 820,000 lakes, rivers and ponds. Rivers are excluded.

The DEM data, used to evaluate the volume of lakes and detect sinks, comes from the IGN ALTI database, offering 25 m x 25 m resolution of the whole of France. Due to the large size of the territory to process, the whole dataset (WB and DEM) was re-organized into a 15 m x 15 m-cell grid. This represents 138 land rasters measuring about 70 km down one side. To be certain that all depressions and further connections are taken into account, even between raster-cells, the original rasters were extended with 5-kilometre buffers, on the right and bottom sides, as shown in Figure 2. To confirm the ef-Figure 2: One raster as built in the developed method, with 5 kilometres buffer ficiency of the Hollister et al. [START_REF] Hollister | Predicting Maximum Lake Depth from Surrounding Topography[END_REF] method, we applied it to 37 French lakes with previously estimated volumes 3 , multiplying the estimation by 0.65, in accordance with the method described in Hollister et al. [START_REF] Hollister | Predicting Maximum Lake Depth from Surrounding Topography[END_REF]. Results are presented on Figure 3 and the correlation is considered sufficient. This method eliminates lakes with an area of less than about 8000 m 2 because of the DEM resolution. A reservoir of this area can be considered as too small for pumped storage application and does not represent a major loss of information in the evaluation process. The environmental protection zones considered for this study are defined as follows, and their constraint ranking is presented in Table 5.

• Conservatory: Lands acquired by the national coastal and lakeside conservation authority, established to protect coastal areas and sensitive green spaces.

• Biotope protection: Natural zones harbouring some wild/protected animals or plants.

• Natura 2000: Protected because it shelters a significant share of the population of a vulnerable species. Different status level, not accessible in the database.

• National and regional natural parks: No fauna and flora preservation rules, but observation and control of activities in the area.

• National and regional reserves: Territory very strictly preserved, generally ruling out human activity.

• International protected site: Highly sensitive or remarkable areas protected by international institutions, like Ramsar or Unesco. [START_REF] Reidy | Development of Bulk Energy Storage and Natura[END_REF].

In order to control the presence of buildings in natural depressions, the BATI layer from the IGN TOPO database is used. This contains all constructions higher than 20 m 2 , including houses, public authorities, religious buildings and sports facilities. The intersection is checked between each detected depression and building. The electric grid is also imported from the IGN TOPO database. It contains high voltage lines, from 400kV to 63kV. Employing the databases mentioned above, about 850,000 connections are detected throughout France.

For the optimization, the knapsack problem could be solved using a Mixed Integer Linear Programming solver (MILP) such as CPLEX [START_REF]IBM ILOG CPLEX Optimizer[END_REF] for example, and the solution can be estimated using an alternative solution. We encoded the Primal Effective Capacity Heuristic (PECH) algorithm of Akay et al. [START_REF] Akay | A Greedy Algorithm for the General Multidimensional Knapsack Problem[END_REF]. This second alternative is selected for two main reasons:

• It gives similar results compared to the exact calculation, as described in Akay et al. [START_REF] Akay | A Greedy Algorithm for the General Multidimensional Knapsack Problem[END_REF]. For the expected size of the inter-connected groups (hundreds of connections between tens of lakes maximum), this algorithm has been proven to be over 90% accurate, and possibly 95 to 99%.

• It allows for much faster computing of the optimization, and this is its biggest advantage. Indeed, the total number of potential connections can be very high: avoiding having to calculate the benefit of every combination saves considerable time.

Finally, the benefit function for the optimization process has been defined as the Energy. Indeed, we do not possess robust information to deduce an optimization of costs. The costs of making a natural depression watertight can represent a major share of the total expenditure if the soil is unsuitable, while it can be almost zero otherwise. Unfortunately, no such cost estimation has been considered possible at a large scale.

A cost study lacking this information does not allow us to compare type 1, 2 and 3 connections and would imply an erroneous selection of pairings.

Overall potential results

After selecting the potential connections using the default restriction values (Table 4), the optimization tool gives us the results presented in Table 6. A France-wide evaluation obtains a result of about 33 GWh. The optimization step, created to ensure that every lake is related to a single connection, divides the potential by 2, confirming the capacity of the algorithm to achieve an accurate potential evaluation. The variant of the method proposed to favour the selection of pairs including existing lakes (building Type 1 and 2 first) indeed increases their proportion in the total number of connections but slightly reduces the potential. Moreover, when counting only connections involving at least one existing lake -with lower financial costs -total energy is about 22 GWh. Finally, only linking existing lakes represents a potential of about 14 GWh .

According to a report by Eurelectric [START_REF] Eurelectric | Hydro in Europe: Powering Renewables[END_REF], France possesses 184 GWh of pumped storage capacity, employing 9 plants. No 4 presents the geographical distribution of the optimized connections across France, in terms of energy. It also gives the proportions of this energy per type of connection, as defined in Table 2. Two main observations can be made regarding this map:

• Most of the storage capacity is present in four mountainous regions: The Pyrenees, the Central Mountain Range, the Jura and the Alps.

• As shown in the pie charts, some areas offer more T1 pairing than others. In the Pyrenees (south-west France), a significant share of the total potential is provided by Type 1 connections. On the contrary, most of the energy available in the Central Mountain Range (central France) and the Jura (north-east France) is required to flood natural depressions.

Since type 1 connections are easier to build, technically and financially, we show that significant potential is concentrated in south-west France. The effects of the different steps of the evaluation process on the number of reservoirs, connections and total energy is presented in Figure 5. The cumulative dis- ). 80% of the potential is reached by building only 425 connections (37% of the total number).

Sensitivity analysis results

Five zones have been considered sensitive as defined in 2.5.1, and the variation of the estimated potential around the default acceptability constraints rankings (see Tables 3 and5) is presented in Figure 7. Several connections are included in Natura 2000 areas and national and regional parks, which are restricted human activity zones. Some of the detected spots might face environmental challenges. We undertook our Global Sensitivity Analysis using 2 matrices of 400 samples of the 9 inputs described in Table 4, corresponding to 4000 evaluations of the potential. The Sobol indices estimated for two cases (including all connections and only T1) are displayed in Table 7. Estimation errors are detected, such as negative values and totals differing from 1. Both can be explained by the fact that first, the number of simulation is limited due to the time-consuming evaluation process (see 3.5), and second, we calculated only the first order of the Sobol indices. In doing so, we do not consider the potential impact of the inputs interactions in the output variability. Homma and Saltelli [START_REF] Homma | Importance measures in global sensitivity analysis of nonlinear models[END_REF] introduced new indices expressing the total effect of a variable X I , including the fraction of variance accounted for by X I alone and the fraction accounted for by its combination with any other variable. These calculation of these total indices has been expressed by Jacques et al. [START_REF] Jacques | Sensitivity analysis in presence of model uncertainty and correlated inputs[END_REF], and the results for our study case are presented in Table 8. Information can be extracted from the analysis.For the overall potential -including T1, T2 and T3 connections -two parameters present a high influence on the global potential: maximum distance between lakes and minimum head-distance ratio. They have a direct, and genuine influence on the evaluated potential, although their impact increases when combined with other variables. Maximum altitude has significant direct impact on the output, but this is not much greater when com- Concerning the second case focusing on T1 connections, maximum altitude and maximum distance between lakes are the two most impacting parameters both directly and combined. Maximum volume seems to have a relative pact if combined with other inputs.

Reasonable variations of the evaluation output

Considering the most influential factors detected above, three abacuses have been drawn representing the variation of the potential depending on these parameters, shown in Figures 8 to 10. We make them vary within a range of reasonable values, while less impacting parameters are set to their default value (Table 4). Similarly, Figure 11 represents the abacus for T1 connections. Maximum altitude and distance to the grid vary in reasonable ranges, while other parameters are at their default state (Table 4). The vertical fracture visible on Figures 9 and 11 at about 2500 m altitude reveals the fast increase of the potential when accepting connection above 2 300 m. Several natural lakes must be present at these altitudes.

Oppositely, the impact is limited above 2 500 m.

Computational performance

The method developed was encoded and run under RStudio [START_REF]RStudio: Integrated Development Environment for R[END_REF] software, using R [START_REF] Core | R: A Language and Environment for Statistical Computing, R Foundation for Statistical Computing[END_REF]. The simulations were performed on a 3.50 GHz desktop computer with 16 GB of RAM. This work aims to permit rapid evaluation of a large zone of energy capacity. It was thus built so that all the cumbersome and time-consuming tasks are run in one go. The actual potential evaluation steps, which require less calculation power, are run in the evaluation loop. The time-processes of the different steps are detailed in Table 9. The selection of the compatible lakes, applying environmental, land-use and technical constraints, only takes a few seconds (less than 2 seconds) while most of the time required is for the optimization process. In the simulation loop, the calculation time varies depending on the restrictiveness of the inputs, which in turn result in discarding more or less connections. The observed simulation times varied between 30 and 3,500 seconds, and a correlation curve can be drawn as on Figure 12.

Best site selection approach

Ahmadi and Shamsai [START_REF] Ahmadi | Preliminary Site Selection of Pumped Storage Hydropower Plants -A GIS-based approach[END_REF] presented a method based on GIS data layers to rank different potential location for building an upper reservoir. Criteria are defined and each location receives a grade (ranging from 0 to 100). This makes it possible to choose the most favourable spots to be evaluated in more precise manner. This section presents a similar approach, ranking connections to select the best sites for building a small-PHES.

Part of the dataset of detected connections are considered and ranked to find which connections should be built first. A French county is considered as a demonstration case, its characteristics are detailed in Table 10.

Each connection is ranked from 0 to 100 in line with 3 criteria defined as follows, using a linear model. • Surface of ground to be sealed (R Surf ): If a natural depression is to be filled, the cost is assumed as proportional to the volume. It is doubled if two natural depressions are paired.

• Length of the conveyance (R Dist ): The greater the distance between reservoirs, the more expensive the construction of the penstock will be. The length is calculated as follows:

L a

H 2 D 2 (14)
• Distance to electrical grid (R Dist grid ): Likewise, the greater the distance between the connection and the electrical grid, the more expensive the connection to the network will be.

Each connection receives a corresponding grade for costs (R C ), defined as follows: R C aR Surf αp1 ¡ aqR Dist p1 ¡ αqp1 ¡ aqR Dist grid [START_REF] Kusre | Assessment of hydropower potential using {GIS} and hydrological modeling technique in Kopili River basin in Assam (India)[END_REF] 0 ¤ a 1: weight of the depression sealing in the cost-like function The a factor represents the cost of sealing a natural depression, and clearly depends on the nature of the terrain. In a restricted area, this cost can be considered as constant for all of the connections, but if precise and more accurate data exist concerning the terrains nature and associated sealing costs, a can be used to weight the cost-like function.

Likewise, the connection energy is ranked from 1 to 100, with 100 being attributed to the highest energy value (R E ). Finally, the "cost of energy"-like value is obtained for each connection.

COE R

C R E (16)
The results are displayed in Table 11. These results lead to two T1 -T1 -T1 323 -99 -48 major remarks. First, Type 1 connections are widely present in the results. This can be explained by their lower cost-like ranking (no surface to be sealed) and their generally higher energy. Nonetheless, in cases where the ground does not need to be sealed, Type 2 connections can be competitive. Secondly, seven connections seemed more worthwhile than others. These connections are represented on Figure 12. Thanks to this analysis, we have been able to detect which connections might be economically advantageous. This can be of great help for a local authority keen to start implementing small-PHES on its territory.

Conclusions and perspectives

Conclusions

The aim of this work was to build a method able to evaluate the small-PHES potential in a massive study zone. The proposed algorithm makes it possible to assess the energy that could be stored using existing lakes or small depressions on a terrain that could be filled. Both types of reservoirs play an important role in the evaluation given that their respective shares, taken separately, mount to less than half of the total estimated potential. Moreover, the proposed method permits fast calculation on the scale of the whole country, thanks to a global dataset generation that only needs running once. Default restrictions are suggested (e.g.. volume, depth, head, etc.) to run this analysis on any area. Applying the method to an entire country, in this case France, allowed us to affirm that the potential increases to about 33 GWh, which is almost 20% of the current national energy storage capacity. Connections are geographically divided between the Pyrenees, Central Mountain Range, Alps and Jura. The proportions of the three pairing types considered (linking existing lakes or natural terrain depressions) vary significantly within each region; composed mostly of existing lakes in the Pyrenees for example, while the potential in the Jura is based on natural depressions. Considering existing lakes only, the potential is estimated at about 14 GWh (about 8% of French current hydro storage capacity), and should be easy to mobilize. These orders of magnitude are rather high, and so we conclude that this avenue merits further study in the future. This paper also presents a sensitivity analysis, in the form of a Sobol study. The most influential inputs are the maximum distance between lakes, their maximum altitude, and the maximum distance to the electrical grid. The corresponding abacuses are drawn to represent the output variations. Global potential varies between 16 and 54 GWh, depending on the restrictiveness of the selected inputs.

However, the different limitations -especially the lack of economic studies and concrete sensitive zone restrictions -imply that this potential may vary in a reasonable range.

This methodology is particularly suitable for large-scale studies. It is important to note that the same process can be carried out for the case of an other country or region. The minimum necessary datasets are a DEM (SRTM is provided by NASA worldwide) and a natural waterbody database. Estimation is improved when data on environmental protection and land use are available.

The resulting connection database could also be used by hydro-power plant developers or local decision-makers to prospect for new sites to install facilities. With this mind, the last part of the study presents a different perspective. The best locations in a French county are sought, using a method similar to a "cost of energy" estimation. Some connections clearly stand out and present favourable conditions for small-PHES implementation. Moreover, the visualisation tends to confirm the effectiveness of the method, which produces coherent pairings.

Possible future ameliorations

Further work could be carried out on this topic. Indeed, even though the value obtained for the global French potential is probably close to reality, some details could further increase accuracy.

Connections chosen by the optimization tool within a group of inter-connected lakes are subject to discussion. The "energy" criterion is used to select the pairs, without any cost dimension. A precise cost function would eliminate the weakest connections and identify the most profitable combinations among connection groups, independently of their type. The prioritization tool would thus be obsolete. Lake usage is not defined in the database used for this study. Some lakes used for leisure purposes were excluded using the CLC database, but it is probable that a number of existing lakes are used for other purposes, or are private, and therefore cannot be considered for PHES purposes. If a database of this kind were to be developed, a share of non-feasible reservoirs could be excluded. Moreover, due to the method used to estimate the volume of natural lakes (i.e. from surrounding topography), ponds are incorrectly evaluated. The only mean to ascertain their volume is to use real-life measurements.

The study does not include roads and transport infrastructure. Verifying the natural depression intersection with such constraints would discard other potential reservoirs. This data was too heavy to consider in this analysis, due to the very high amount of roads on the territory.

The best site selection approach may be improved, if reliable geology data are available. Each connection sealing cost could then be weighted depending on the nature of the terrain. It is likely that the global cost grade could be parametrized with cost dependence to the slope of the pipe, the distance to a road, etc. Storage using the sea or a river as the lower reservoir was not considered. Using the sea would involve a different technology which is still in development and thus more expensive; river storage would entail difficulties in processing the river data we possess, such as river size and water flow, to ensure that they can be used as a source of water. Their inclusion in the evaluation could increase the global potential.

  which allows us to decompose the variance V of Y in the same way. Normalizing by V finally generates Sobol indices S I : 1 Iφ S I , with S I V I {V and V I varpf I pX I qq (8)

4 .

 4 Sobol indices estimationFirst, note that thanks to the definition of f I pX I q in Equation 7 we can recursively estimate Sobol indices, i.e. S I varpErY |X I sq V estimate S I , given that we are only considering first order indices, we only need to estimate varpErY |X I sq (and

Figure 3 :

 3 Figure 3: Results of the estimation of lake volumes on a sample of lakes with known volumes.

Figure 5 :

 5 Figure 5: Evolution of number of connections and energy available in steps 3 and 4 (see Fig. 1). Cs: Connections; Rs: Reservoirs

Figure 6 :

 6 Figure 6: Cumulative distribution of energy in connections.

Figure 7 :

 7 Figure 7: Variation of the evaluated storage capacity depending on some exclusion zones' status.

Figure 4 :

 4 Figure 4: Distribution map of feasible connections in France, and proportion of types per county (Source of background map: http://www.cartesfrance.fr/)

Figure 8 :

 8 Figure 8: Energy detected (GWh) depending on the maximum distance between lakes and the maximum altitude.

Figure 9 :

 9 Figure 9: Energy detected (GWh) depending on the minimum head-distance ratio and the maximum altitude.

Figure 10 :

 10 Figure 10: Energy detected (GWh) depending on maximum distance between lakes and minimum head-distance ratio.

Figure 11 :

 11 Figure 11: Energy detected depending on maximum distance between lakes and maximum altitude. Only Type 1 connections are considered.

Figure 12 :

 12 Figure 12: Evolution of process-time depending of the number of selected lakes (i.e. not restricted by constraints)

Table 1 :

 1 Basic settings for the connections dataset building.

Table 2 :

 2 Denomination of connections depending on the nature of connected reservoirs

	Name	Initials Description
	Type 1	T1	Pairs composed of two existing reser-
			voirs
	Type 2	T2	Pairs composed of one existing lake
			and one terrain depression detected
			by the DEM method
	Type 3	T3	Pairs composed of two terrain de-
			pressions detected by the DEM
			method

Table 3 :

 3 Ranking of the different constraints, environmental and land-use

	Restriction	Ranking
	Urban zones	No-go
	Leisure activities	No-go
	Arable land	No-go
	Permanent farming

Table 4 :

 4 Limitation for different connection parameters

	Criterion	Limit	Default value
	Head	Minimum	40 m
	Distance	Maximum	2 500 m
	Volume	Minimum	0 m 3
	Volume	Maximum	750,000 m 3
	Depth	Minimum	0 m
	Grid distance	Maximum	10,000 m
	Altitude	Maximum	2 500 m
	Energy	Minimum	5000 kWh
	Head/Distance ratio	Minimum	0.1

energy threshold so that the minimum sizing accepted correspond to a 500 kW plant active for 10 hours. The minimum Head/Distance ratio corresponds to a 10% slope

[START_REF] Berezinskii | Design of penstocks of pumped-storage stations[END_REF]

[START_REF] Lu | A GIS-based assessment of Tibet's potential for pumped hydropower energy storage[END_REF]

.

Table 5 :

 5 Ranking of the different environmental and land-use constraints

	Restriction	Ranking
	Conservatory	No-go
	Biotope protection	No-go
	Natura 2000	Suitable
	Nat. & Reg. parks	Suitable
	Nat. & Reg. reserves	No-go
	International protected site	No-go
	Natura 2000 is not considered as No-go despite its environ-

mental protection status. Such areas do not in themselves exclude PHES plants and a guide to promote PHES project development within Natura 2000 zones has even been published by stoRE

Table 6 :

 6 Results of the different optimization configurations for the default restrictions. For proportions: Pink: Type 1, Gold: Type 2, Green: Type 3 has been carried out since this report. The assessed potential represents an energy increase of about 18% considering all connection types, and 8% using only T1. Figure

	Connection	Prioritization	Energy before optim.	Final Energy	Proportions of built connections
	Nature		(GWh)	(GWh)
	All	No	78.42	33.56
	All	Yes	78.42	32.93
	T1 and T2	Yes	51.18	21.76
	T1	/	30.35	14.35
	new major project

Table 7 :

 7 Sobol indices estimation for different evaluation scenarios: Case 1: Any kind of connection, with prioritization ; Case 2: Only existing lakes

	Parameter	Case 1	Case 2
	Min. Head	-0.003	-0.006
	Max. Distance	0.285	0.119
	Min. Volume	-0.012	-0.067
	Max. Volume	0.015	-0.013
	Min. Depth	0.023	-0.067
	Max. Distance/grid	0.005	0.012
	Max. Altitude	0.326	0.544
	Min. Energy	0.017	-0.059
	Min. H/D ratio	0.514	0.083
	Total	1.17	0.423

Table 8 :

 8 Sobol total indices estimation for different evaluation scenarios: Case 1: Any kind of connection, with prioritization ; Case 2: Only existing lakes

	Parameter	Case 1	Case 2
	Min. Head	0.391	0.111
	Max. Distance	0.678	0.375
	Min. Volume	0.394	0.114
	Max. Volume	0.465	0.321
	Min. Depth	0.411	0.118
	Max. Distance/grid	0.418	0.115
	Max. Altitude	0.467	0.784
	Min. Energy	0.362	0.104
	Min. H/D ratio	0.647	0.114
	bined with other parameters.

Table 9 :

 9 Process times of the different steps of the method.

	Step	Process	Principal difficulties
		time
	Reservoirs	10-15 days Data import, Volume
	detection		evaluation, Sink detec-
	Connections	3-5 days	tion Matching of reservoirs,
	detection		Fusion of all rasters,
			Verification of accept-
			ability constraints
	Constraints	2 seconds	Binary vectors selection
	application
	Optimization	30-3 500	Optimization	within
		seconds	each	inter-connected
			lakes group

Table 10 :

 10 Characteristics of the demonstration case county.

	Name	Hautes-Pyrenees
	Area	4 446 km 2
	Number of connections detected	59
	Number of T1 connection	35
	Number of T2 connections	16
	Total evaluated energy	4.15 GWh

0 ¤ α 1: balance between pipe and electrical connections weights.

Table 11 :

 11 Results

of the "Cost of energy"-like approach.

Table 12 :

 12 Visualization of the connections selected by the COE method. Source of background maps: Google Earth, accessed on Mar 03, 2017

	ID Visualisation	ID Visualisation
	1	5
	2	6
	3	7
	4

http://blog.protectedplanet.net/

http://adour-garonne.eaufrance.fr/

Acknowledgement

The research work in this paper was carried out in the frame of the research and demonstration project "Micro-step de BERRIEN" funded in part by the French Energy Agency ADEME under the program "Energies Decarbonnees" from "Investissements d'Avenir".