SSDP Model with Inflow Clustering for Hydropower System Operation

Sampling stochastic dynamic programming (SSDP), which considers the uncertainty of streamflow, is a popular and useful method for solving release decisions of reservoirs. It is easy to implement the long-term operation for cascaded hydropower system with poor inflow prediction ability. Furthermore, SSDP describes the temporal and spatial structure of the stochastic streamflow processes implicitly through inflow scenarios instead of representing the multivariate distribution of inflow by conditional probabilities in stochastic dynamic programming (SDP). However, computation time of SSDP procedure will increase exponentially with the growth inflow scenarios. Thus, clustering algorithm is employed to reduce the number of inflow scenarios in order to improve the efficiency and operability of SSDP in practical applications. The calculation results of SSDP with inflow clustering are analyzed with different cluster numbers. The principle of how to find the least inflow scenarios to represent all inflow sequences has also been proposed. The least inflow scenarios and relevant probabilities found by clustering algorithm can approximate the empirical distribution of all streamflow scenarios used in this study without obviously decreasing the energy and exacerbating the shortage of firm power.


Introduction
Hydropower plays an essential role in the operation of the power grid because of its largescale installed capacity and flexible performance in electric load adjustment. Since an unprecedented hydropower boom in 2000 (Wu et al. 2020), China has been the biggest hydropower producer in the world, with a total capacity of 391 GW in 2021. With the consistent development of hydropower in China, the cascaded hydropower stations will account for more shares in the power grid.
How to operate these cascaded hydropower systems to achieve the established goals becomes a significant challenge for operators under the uncertain inflows. In the past decades, various methods have been reported to optimize the operation of hydropower stations, such as Heuristic algorithm (Arunkumar and Jothiprakash 2013;Afshar et al. 2015;Rahimi et al. 2020), Mixed-Integer programming (Catalão et al. 2010;Tong et al. 2013;Doganis and Sarimveis 2014) and dynamic programming (DP) evolved methods (Mariño and Mohammadi 1983;Lei et al. 2017;He et al. 2019). Among these types of methods, DP-evolved solutions are the most accepted algorithms because the water release policies of hydropower stations can be resolved into some continuous decision processes.
Considering the uncertainty of natural inflow, stochastic dynamic programming (SDP) was proffered by Loucks et al. (1981), and many SDP-based techniques were applied to guide the release decisions of hydropower stations (Carpentier et al. 2014;Xu et al. 2014). In the recursive equation of SDP, it captures the information about streamflow by setting extra hydrologic state variables. The recursive equation of SDP can't evaluate all possible values of state variables, instead, it must calculate a set of discretization state variables. However, as the number of cascaded reservoirs increases and higher requirement for calculation accuracy, SDP suffers from the curse of dimensionality, result in exponential increasing in time consumption. SDP also can't describe the cross correlation among inflows in different reservoirs and autocorrelation of inflow at each reservoir concisely. To solve this problem, Kelman et al. (1990) put forward sampling SDP (SSDP) to combine SDP with a scenariobased procedure. In this study, each specific inflow scenario is made of the local inflow of every reservoir in a year representing a stochastic process with 12-month. SSDP uses the inflow sequences to take into account the uncertainty of inflow implicitly that aims to maintain the temporal correlation or long-term persistence of the inflow concisely. SSDP have been studied by some researchers and presents an acceptable performance (Zhao et al. 2017;Haguma and Leconte 2018;Mohanavelu et al. 2022). Nevertheless, the recursive equation of SSDP will become much tougher to solve due to the increase of streamflow scenarios. To overcome this drawback, several researchers have introduced scenario reduction techniques to decrease the number of inflow scenarios. Clustering is a machine learning technique that aims to divide unlabeled samples into multiple clusters. Anvari et al. (2014) had illustrated the improvement in time consumption by using K-means to cluster the ensemble streamflow prediction traces from 225 to 3 including dry, normal and wet conditions. In addition, Côté and Arsenault (2019) employed a clustering algorithm to reduce the available inflow scenarios from 63 to 7 through a parameter of normalized distance matrix. To make SSDP model more efficient, it's an appropriate solution to cut down the inflow scenarios by clustering all inflow sequences into several groups. Subsequently, the most representative inflow sequence in each group and relevant probability distribution can be obtained to approximate the empirical distribution of all inflow scenarios. K-Means is a suitable method to separate unmarked data samples into several sections. Nevertheless, K-Means calculates the average value of all samples in a cluster as the class center in every group, it breaks the natural temporal and spatial structure of the 12-month streamflow sequence. Furthermore, previous studies specified a fixed cluster number, without analyzing the rule of how to get an appropriate cluster number of inflow scenarios.
In this study, SSDP is adopted to make a long-term operation schedule considering the uncertainty of natural inflow. To alleviate the curse of dimensionality and retain the spatial persistence of the streamflow scenarios (Wu et al. 2022), K-medoids (another clustering method which is similar to K-means) is employed to clustering the inflow scenarios into several groups. Additionally, this paper utilizes disparate parameters to represent inflow scenarios and qualitatively ascertains the rule of how to obtain an appropriate number of clusters. To test the performance of proposed model, it is implemented in a multi-reservoir hydropower system of Wu River, China, which has six hydropower stations in series and two mainly reservoirs. This study can also provide valuable references for the operation of cascaded hydropower system.

Study Area
The Wu River originates from the eastern edge of Wumeng Mountains, with two sources in the South and North. The Sancha River in the south is 322 km long, which is the main source of the Wu River. The Liuchong River in the north is 210 km long. When the two sources meet, it is called Wu River. The Wu River is the largest tributary on the south of the upper reaches of the Yangtze River, with a total drainage area of 87900 km 2 , a total length of 1050 km from the origin to the estuary, and a drop of 1787.46 m. The annual average precipitation in the western part of Wu River Basin is about 1000 mm, 1000-1200 mm in the central part and about 1400 mm in the southeast part. The rainy season mostly occurs from May to October. Wu River is one of the main power sources of "West to East power transmission" projects of China. The studied cascaded system on Wu River is on the main stem in Guizhou Province, including Hongjiadu, Dongfeng, Suofengying, Wujiangdu, Goupitan and Silin stations seen from Fig. 1, with a total installed capacity of 7195 MW. Inflow data from year 1953 to 2008 are used for solving and simulating. There are 12 × 56 networks to be trained in each round of SSDP program. Table 1 shows the beneficial storages of Hongjiadu and Goupitan are much larger than the other reservoirs, and they are main long-term regulation reservoirs of this system. The regulation abilities of Dongfeng, Suofengying, Wujiangdu and Silin are neglected in this study, the reservoir levels are set to be fixed values of 969 m, 836 m, 759 m and 438 m correspondingly.

Studied Hydropower Operation Model
Developing a long-term operation strategy for cascaded hydropower stations is challenging due to the uncertainty of inflow and the poor accuracy of inflow forecasting in long-term. According to the monthly operation of the cascaded hydropower system, resolution on the last period has an immense impact on future available water volume and reservoir level (Cheng et al. 2014;Shen et al. 2021). Different current release decision of upstream reservoir will consequently lead to dissimilar ending storage for upstream reservoir and diverse expected future benefits for all reservoirs (Kişi 2009). Thus, making a release decision of a leading reservoir should consider conditions in more aspects.
For the step-to-step decision processes, DP is appropriate to deal with this problem. The deterministic DP is described as: Equations (1) and (2) are the recursive equation of DP and the mass balance equation respectively. Where f t (S t ) and f t+1 (S t+1 ) are the objective functions of all reservoirs for reservoir storage state vector S t ; is the benefit function at period t when S t , Q t and R t are given; e t (S t , S t+1 ) is the evaporation loss at period t.
Using Eq. (1) can work out a release result for the cascaded reservoirs efficiently when the inflow of every reservoir is known. However, inflow can't be forecasted exactly with existing techniques (Alfieri et al. 2020). Thus, stochastic DP is a natural extension to be employed to make a strategy that maximizes the expected future benefits. The classic recursive equation of SDP can be concluded as follows (Loucks et al. 1981;Stedinger et al. 1984): , R i * t is the target release of reservoir i at period t; and the actual release R t can be solved by Eq. (4). Where E Q t |X t stands for the conditional probability of Q t when X t is given. One more hydrologic state variable X t is supplemented to make up the recursive equation of SDP.
In general, the enhanced operation rules can be reaped by adding hydrologic state variables. They are supplemented to describe the persistence of streamflow with a Markov process which is usually represented by the previous month's inflow. Many improved SDP methods have been reported accord with actual situations, such as Gal (1979) adopted two preceding inflows as the hydrologic state variable. Stedinger et al. (1984) employed snowpack to forecast the snowmelt season's runoff then let the forecast of remainder seasonal runoff as hydrologic state variable. In the recursive equation of SDP seen in Eq. (3), the newly added hydrologic state variable can capture the characteristic of inflow with an explicit probabilistic description through classifying the stochastic variables into discrete cases with representative values.
Allowing for the complex temporal and spatial structure in the streamflow processes, SSDP is an effective method to improve the solution. In the recursive equation of SDP, the correlation or long-term persistence of the streamflow also can't be hold concisely. Therefore, Kelman et al. (1990) presented SSDP to reflect the various characteristics of inflow by incorporating the inflow scenario into the recursive equation directly instead of classifying the inflow and the transition probability into discrete cases. The formulation of SSDP can be described as follows: Where m represents various streamflow scenarios (m = 1, 2, …, M). In the traditional recursive equation of SDP, the optimal release and optimal future value function (FVF) are both determined by Eq. (3). Different from SDP, in the SSDP model, the optimal target release for each stage is calculated by Eq. (5), however, the FVF for each period is updated separately for each scenario by Eq. (6). The actual release is also computed by Eq. (4). ( The essential advantage of SSDP compared to SDP is that SSDP uses multiple streamflow sequences to describe the potential empirical distribution of inflow in the recursive equation. Consequently, the time correlation or long-term persistence of the inflow which may exist in history inflow scenarios can be conveyed concisely. The more realistic operation policies can be obtained spontaneously without explicit probabilistic description for the stochastic streamflow processes. For the transition probabilities among inflow scenarios of P[j|i] (where i and j refer to different inflow scenarios) in this study, simplistically assuming that there are few year-toyear annual correlations between different streamflow scenarios, hence the probability of each inflow scenario is identical which can be expressed as 1/M. For simplification, the inflow forecast ability for the cascaded reservoirs is ignored due to the poor accuracy of long-term forecasting. The recursive equations of classic SSDP model can be rewritten as Eqs. (7) and (8)

Inflow Clustering Method
A lot of improved methods have been proposed to alleviate the curse of dimensionality in SSDP. The precise and realistic policies can be figured out by classifying the storage of reservoir or reservoir level with more discretization state variables. However, it will exponentially increase time consumption, especially for a multi-reservoir system. The curse of dimensionality will be also aggravated by employing a large number of runoff scenarios in SSDP. To mitigate the curse of dimensionality, it's an appropriate solution to decrease the inflow scenarios by clustering all inflow sequences into several groups then select a representative inflow scenario from each group.
Clustering is a machine learning technique, which aims to divide unlabeled samples into multiple clusters. The samples in the same cluster are similar, while samples in different groups have high degree of difference. Similarity and dissimilarity are usually measured by the Euclidean distance among diverse samples. Clustering is used in this study to divide the streamflow scenarios into several clusters, then choose the cluster center in each group to be sampled in the SSDP procedure. K-means algorithm is a typical clustering algorithm in data mining and can be described as: assuming an integer K and a dataset S with n data samples, the objective function is to find K cluster centers to minimize the total sum of squares errors (TSSE). The formula of TSSE is: Where n is the number of samples, x i stands for each data sample, c(i) is the cluster center of sample i. The implementation of K-means algorithm is easy and concrete steps can be concluded from Arora and Varshney (2016). Since K-means algorithm is not effective when used with global cluster (Arora and Varshney 2016) and calculates the average value of all samples in a group as the cluster center which will be used for sampling in SSDP, so that K-medoids algorithm which is based on an object representative technique (Park and Jun 2009) is employed for clustering in this study. The implementation steps of K-medoids are similar to K-means. It chooses a sample as the cluster center rather than calculating all samples' average value when determining the centroid in a cluster. Furthermore, the persistence of the central streamflow scenarios will be maintained by applying a natural streamflow scenario and corresponding probability distribution to represent one inflow group directly. This method caters to the advantages of SSDP which tersely captures the potential temporal and spatial structure in the stochastic streamflow processes. In cascaded hydropower system, the reservoir operation policies will be impacted greatly by the water volume of inflow and the distribution of inflow in time. In this study, the similarities are appraised by water volume and distribution over time (Faber and Stedinger 2001), different with Côté and Leconte (2016), assuming the Wu River possesses few correlations between different years and strong temporal correlation within year. Therefore, the annual inflow is divided into flood season (May-October) inflow and dry season (November-April) inflow for every reservoir to represent an inflow scenario as follows: Where i is the label of inflow scenario. Q i flood and Q i dry are the corresponding total inflows in flood season and dry season. i is the converted data sample with 12 dimensions. There are six hydropower plants in Wu River in series as shown in Fig. 1, for ease of expression, adopting hjd, df, sfy, wjd, gpt and sl as abbreviations in Eqs. (11), (12) and (13). For a cascaded hydropower system, the release policies of leading station have enormous affect for whole system when it possesses long-term regulated reservoir. In the Wu River cascaded hydropower system, the upstream Hongjiadu reservoir has carryover year storage, and the installed capacity is much smaller than the Goupitan reservoir. So that the cascaded firm power is satisfied through the drawing of Hongjiadu reservoir in most periods, unless the reservoir level is very low. While the reservoir level of Goupitan reservoir is kept high, unless the cascade firm power can't be satisfied or there is spill risk in future periods. Whereupon the annual inflow of Hongjiadu and annual inflows of the other five reservoirs are applied to express a streamflow sequence for Wu River as below: Where i denotes the transformed sample with 2 dimensions. Taking into account the different expressions of a streamflow scenario, two indexes are employed: the first one is dividing the annual inflow of each reservoir into two parts, total inflows in flood season and total inflows in dry season, with a total of 12 dimensions to represent a streamflow scenario expressed as 12-D; another is taking the annual streamflow of Hongjiadu that is the leading hydropower station in Wu River with long-term regulated reservoir and the sum of annual streamflow of other five stations including Dongfeng, Suofengying, Wujiangdu, Goupitan and Silin with two dimensions to describe a streamflow sequence concluded as 2-D. Different from existing researches, this study focuses on how to get the appropriate cluster number by evaluating the TSSE index and clustering result instead of specifying wet, normal and dry conditions or giving a fixed cluster number simply.

Results of Inflow Clustering
The cluster number in the range of 2 to 49 was computed by K-medoids in order to get appropriate cluster numbers through analyzing the variation trend of the TSSE curve as shown in Fig. 2 and the result of clustering in Fig. 3. The clustering effect of 12-D is not as good as that of 2-D because the data samples of 12-D have more dimensions, more discrete the samples will become. The aggregation of 12-D is also not high enough, causes poor clustering effect. Thus, subsequent analysis is mainly conducted for 2-D. The inflection point of the TSSE value is applied as the criterion to determine a preliminary number of clusters (Park and Jun 2009). According to the definition of TSSE, when cluster number increases constantly, TSSE value will gradually decrease to 0. Finally, the number of clusters will equal to the number of inflow scenarios. When the K value approaches the optimal number of clusters, the TSSE value will descend slowly as K increases rather than decline steeply. The cluster number at this turning point is the preliminary value of K that needed. It can be known from Fig. 2 that the inflection point appears at K = 5. To improve the operability and practicability of SSDP with clustering algorithm, it is inappropriate to determine a fixed value of K. Accordingly, an expanded range of cluster numbers based on the initial value will be applied to analyze the effect of different cluster numbers on the computation efficiency and the rule of how to decide a suitable cluster number, the introductory cluster number 5 is extended to the range of 3 to 8. The clustering results of 2-D with 3-8 classes achieved by K-medoids are shown in Fig. 3. Where the red dots stand for the central inflow scenarios for each class. Subsequently, the rule of how to derive the best cluster number from inflow clustering can be observed. For ease of explanation, some hydrologic conditions will be defined. The probability distributions of all inflow scenarios are fitted by the Peason-III curve. The separator lines are marked in Fig. 3 corresponding to the probability distributions of 30% and 70% respectively. Every subgraph is divided into three sections and individually corresponding to wet, normal and dry conditions in hydrology from top to bottom. For the history inflow scenarios used in this study, they mainly distribute around the separator lines of wet and dry conditions. When there are few classes, clustering algorithm tends to divide the inflow scenarios into several contrasting zones. For example, when the cluster number is 3 as shown in Fig. 3a, the divided areas can be described as three blocks which corresponding to wet, dry and especially dry conditions. Subsequently, with the raising of cluster number, K-medoids will start to separate some extreme scenarios into a separate category. Giving an example, when cluster number grows to 4 and 5 based on 3, the scenarios in upper part of cluster 3 are divided to additional clusters seen from Fig. 3b, c; as it expands to 6, Fig. 3c, d show that the especial dry inflow scenarios are isolated at this time. After classifying these rare inflow scenarios, the most concentrated area of inflow data will be sundered. When the cluster number is going to rise from 6, the inflow scenarios around separator lines begin to be split into several parts that can be concluded from the division shown in Fig. 3e, f based on cluster 2 and cluster 3 in Fig. 3d.
Clustering the inflow scenarios into several pieces, it aims to obtain the similar empirical distribution of all inflow samples by least clusters. In the inflow series of a watershed, the most concentrated sequences that need not to be divided into several sections, but inflow scenarios in the especially wet and dry conditions needed. Therefore, the least number of clusters acquired by clustering algorithm should be the one after separating special cases (such as cluster 1, cluster 4, cluster 5 and cluster 6 in Fig. 3d), in this paper, 6 is considered as the best cluster number for all inflow scenarios.

Results of SSDP with Clustered Inflow Scenarios
The central inflow scenarios in each cluster and corresponding probabilities are employed in the clustering SSDP model, then conduct a comparative analysis between the clustering SSDP model and the original model. In the long-term operation of a cascaded hydropower system, both firm power and energy generation are the most critical indicators. The total square sum of lacks between power and firm power (TSSLP) for every model is calculated by: where P i,t is the power of hydropower station i at period t, P firm is the firm power of whole system, T1 is the period set when ∑ n i=1 P i,t < P firm . In this study, the firm power of Wu River cascaded hydropower system is set to be 1900 MW. The objective function is set to minimize the value of E − ⋅ TSSLP , where E is energy generation of whole system and is the penalty coefficient for TSSLP. Due to the probability of each inflow scenario is identical as 1/M, a/M will be used in Eq. (7) for each inflow group after clustering, where a is the number of inflow scenarios in each cluster. In the 56 years simulated, the TSSLP and energy generation of different models are shown in Fig. 4. From the perspective of firm power, all computation results of 2-D are better than that of 12-D due to the data samples of 2-D fit the main factors of inflow scenarios in Wu River, while 12-D don't. When the cluster number is 6, the TSSLP of 2-D is very close to the original model which is only 19% higher. For the highest TSSLP value in 2-D (the clusters number is 3), it is 53.5% larger than the original one which is much less than the best one of 12-D. In terms of power generation, the calculation results of 12-D are more impressive, even the minimum energy is also higher than the without clustering one, and the maximum rise in power generation reaches 0.33%. As a contrast, the power generation is larger than the without Fig. 4 Energy generation and TSSLP clustering case only when the number of clusters is 7 and 8 in 2-D. However, for the objective function which includes energy and TSSLP indexes, the results of 2-D are much better than that of 12-D. It can be concluded that it is proper to clustering inflow sequences with 2-D. In the results of 2-D, the objective function is increasing all the time with the additive cluster number from 3 to 7. Nevertheless, when the cluster number is greater than 6, the clustering algorithm begins to divide the samples in the most concentrated area, and the promotion for objective function is very lightly, hence 6 is the least cluster number for the data set used in this study. As shown in Fig. 4, for the model 2-5 (5 is the inflection point of TSSE curve) and model 2-6 (which achieves good simulation effect for all inflow scenarios with the least clusters), more information is analyzed in the following.
Two inflow clustering consequences with cluster numbers of 5 and 6 with SSDP expressed as model 2-5 and model 2-6 are computed for comparative analysis with the original model. Table 2 shows the energy in wet season and dry season and the increase/ decrease percentages calculated by model 2-5 and model 2-6 contrast to the original model. Wu River cascaded hydropower system reduces the energy in wet season and improves the power generation in dry season for both model 2-5 and model 2-6, the electricity transfers from wet season to dry season within a year can be realized more positively. For the model 2-5, the energy of Wu River cuts down 0.78% in wet season, compared to model 2-6, it reduces by 0.75%. In dry season, the power generation of model 2-5 increased by 1.13%, which is smaller than that of model 2-6 increased by  Figure 5 shows the alteration tendency of the average power in each month for each station calculated by model 2-5, model 2-6 and the original model. The temporal distributions of power include Hongjiadu and Goupitan reservoirs that both hold carryover year storage changed significantly compared with the other four runoff hydropower stations for model 2-5 and model 2-6. It can be known from Fig. 5 that in wet season of model 2-6, the average power curves of the four runoff hydropower stations are almost the same with the original model, the difference is that Hongjiadu generates more energies while Goupitan reduces. A conclusion can be also summarized that the reservoir which has regulation storage is more important for the cascaded hydropower system on the temporal distribution of power because the average power processes of Dongfeng, Suofengying, Wujiangdu are nearly the same with their upstream reservoir, Hongjiadu. The average power curve of Silin is also similar to Goupitan shown in Fig. 5.
Finally, the time consumption of multifarious models are exhibited in Fig. 6. For the model 2-3, it barely spends 3.63% of time of the original model in computing. Even for the most time-consuming one corresponding to model 2-8, a substantial reduction can be also achieved contrast to the primal model. This study also stresses the fact that the proposed SSDP with inflow clustering model is effective in mitigating the disaster of dimensionality.

Conclusions
In this paper, sampling stochastic dynamic programming (SSDP) and clustering algorithm are combined to make long-term release decisions for the six hydropower stations on Wu River with two regulation reservoirs. For the classical SSDP model, it is trapped in the disaster of dimensionality caused by a large number of inflow scenarios and abundant hydrologic state vectors. Clustering inflow scenarios to reduce the number of sequences can effectively alleviate this trouble. K-medoids algorithm can perfectly fit the advantage of SSDP that concisely considers the temporal and spatial correlation in stochastic streamflow process. In addition, in order to enhance the operability in the practical application, the best cluster number of inflow scenarios is not determined by a specified value but decided on the effect of clustering which can reflect the empirical distribution of historical inflow scenarios better. Furthermore, the principle of how to find the least and suitable cluster number is reported qualitatively. The least number of clusters required should be the one after separating special cases. The computation results also show that the proposed method is practical and effective in mitigating the disaster of dimensionality compared with original SSDP model without obvious distinctions in decreasing the energy and exacerbating the shortage of firm power.