Abstract
An original nonsequential Monte Carlo simulation tool is developed. It permits to compute the optimal dispatch of classical thermal generation in order to minimize pollutants emissions in presence of wind power and under operating constraints.
1. Introduction
Nowadays, most of the conventional electrical parks are still using fossil resources like coal or oil. Those primary resources involve the emission of gaseous pollutants like carbon oxides (COx) or oxides of nitrogen (NOx). Recently, following the Kyoto agreements, a great research effort has been made in order to reduce those emissions worldwide. In this context, one of the most promising alternative resources is certainly wind power. Given the fluctuating behaviour of the wind and due to several operating constraints (cost, reliability, etc.) related to electrical systems, it is important to adequately dispatch conventional thermal generation in order to guarantee a minimization of the gaseous pollutants emissions and to simultaneously face the requirements of modern networks. In that way, [1, 2] have proposed analytical models in order to minimize emissions under cost and load covering constraints. However, even if they represent a first approach, those models are limited by several hypotheses. For example, it can be quoted that the hourly load fluctuations cannot be analytically taken into account and that the computations are thus made for a single load value in [1, 2]. Moreover, investigated correlation scenarios between wind parks are relatively limited with the analytical approach and wind generation is generally supposed to be entirely correlated over the considered territory. Finally, as only one single load value is considered, it is not possible to compute reliability indices like the Loss of Load Probability (LOLP) or the Loss of Load Expectation (LOLE) [3, 4] with the existing analytical tools. In order to avoid the drawbacks of the analytical approach, it is therefore proposed in this paper to develop a nonsequential Monte Carlo simulation tool [5, 6] that permits to consider not only a fluctuating behaviour for the load, but also several correlation scenarios between wind parks. Thanks to the developed tool, an optimal dispatch of the classical thermal generation can be realized for each generated system state, reliability indices can be evaluated, and wind generation long-term impact on gaseous pollutants emissions can be studied. Concerning this last point, it is important to note that the simulated states are independent steady-state ones and that, therefore, short-term transitions between system states cannot be studied here. By consequence, the eventual increase of CO2 emissions due to the fast starting of coal thermal units in case of unexpected lack of wind is beyond the scope of the present study. Nevertheless, in order to prove the utility of the developed simulation tool for the system operators from the long-term electrical dispatch point of view, it has been applied to a slightly modified version of the Roy Billinton Test System (RBTS) [7].
Finally, the present paper is organized as follows. In a first section, the general algorithm of the developed Monte Carlo simulation tool is presented. Then, an optimization model for minimizing the Nox emissions of thermal units is described. In a third part, some simulation results and applications of the developed tool are proposed in order to show the impact of wind power on pollutants emissions. Finally, a conclusion pointing out the main interests of the present work is proposed.
2. General Algorithm of the Developed Monte Carlo Simulation Tool
Monte Carlo simulations are generally used to simulate the actual process and random behaviour of a given electrical system. The pursued objectives are, for example, the computation of reliability indices, the search for profitable investments scenarios, and so on [4, 6]. Theoretically, there are two basic techniques used when Monte Carlo methods are considered for power system applications, these methods being known as the sequential and nonsequential techniques [6, 8]. In that way, in the present study, a non sequential Monte Carlo algorithm (Figure 1) has been developed under Matlab to evaluate reliability indices of interest, gaseous pollutants mean emissions, adequate dispatch of classical units, and so forth. This algorithm is only limited to Hierarchical Level (HLI) [4]. That means that the total system generation is examined with the total system load requirement on a pooled basis. Moreover, the constraints (ampacity) on the transmission system are ignored.
The developed Monte Carlo simulation could theoretically incorporate any number of system parameters and states but it has been assumed in our calculation that a classical generation unit was only able to reside in one of the following two states: fully available and unavailable. Moreover, in the established non sequential simulation, only hourly uncorrelated states are considered as it is supposed that a generation unit outage state does not condition or influence its state during the next or previous hours of simulation. Consequently, at the start of each hour, a uniformly distributed random number (u) on the interval is sampled for each classical (thermal, hydroelectric,…) generation unit in order to decide its operation state, using the following process:(i)If u ≤ Forced Outage Rate (FOR) [4], the classical unit is considered unavailable.(ii)If u > FOR, the classical unit is decided to be fully available.Concerning the hourly load of the system, its determination in a nonsequential algorithm can be practically based on a random sampling over its cumulative distribution function [5] or, more precisely, established by the use of modulation diagrams of the annual peak load value [9].(i)Diagram of weekly modulation of the annual peak load: this last one permits to calculate the peak load of the current week on the basis of the annual peak load value. This diagram contains thus 52 modulation rates of the annual peak load value.(ii)Diagram of the hourly modulation of the weekly peak load: it permits to calculate the hourly load for each hour of the week. This diagram contains thus 24 modulation rates of the weekly peak load value.Using this last methodology, no random sampling is needed in order to generate the hourly load of the system. More easily, the program just considers, in the weekly modulation diagram, the rate corresponding to the current week during the simulation process (Figure 1). Then, it associates to the generated weekly peak load the rate of the hourly modulation diagram corresponding to the investigated hour of the day. Also note that, as the consumption during one week can change from one day to the other (days of the week, Saturday or Sunday), several diagrams of hourly modulation can be defined during one week. Moreover, seasonal aspects can also be taken into account by defining periods during the year and by changing the set of hourly modulation diagrams associated to the load from one period to the other.
From the wind generation point of view, several modeling methods like the ones based on autoregressive moving average models have already been proposed in the literature and require increased computing capacities [10]. Here, as we consider a nonsequential approach, it is not necessary to establish such a complex chronological model for wind speeds but it will be more efficient to introduce a random sampling based on a defined statistical distribution [5, 9] associated to each wind park. Practically, the nonsequential wind speed sampling is generally based on a direct inversion of the cumulative distribution function [5, 9]. More precisely, the process can be described as follows:(i)sampling of a uniformly distributed number “u” on the interval ;(ii)application of that sampled random number “u” to the cumulative distribution function in order to determine the associated wind speed.Given the expected wind correlation scenario, the same random number “u” can be sampled for all the wind parks. In that case, they are all supposed to be entirely correlated. At the opposite, the parks can be entirely independent by sampling an independent random number “u” for each of them. Otherwise, more accurate correlation scenarios can be reached by combining for each wind park an adequately computed noise with a single mean distribution [11]. Finally, note that a power curve is associated to each wind park and permits to convert the sampled wind speed into generated power.
At the end of the system states generation process, reliability indices like the LOLP and the LOLE [3, 4] can be computed based on the simulated states during which the hourly load cannot be fully covered with the available (classical + wind) generation.
Moreover, during the Monte Carlo simulation process, each generated system state must be analyzed in order to minimize (under cost and load covering constraints) the gaseous pollutants emissions by adequately dispatching the available (classical + wind) generation. In that way, the implemented optimization model is presented in Section 3.
3. Optimization Model Analyzed during Each Generated System State
Before presenting the optimization model developed in this section, it is important to note that only coal and oil thermal units, hydroelectric generation, and wind parks are considered in this paper. Nevertheless, the proposed model could be easily extended to other kinds of production means like nuclear units
The optimization model presented in this section tends to minimize NOx emissions under cost and load covering constraints. In that way, the pursued objective is based on the definition of an Environmental Impact Index (EII) expressed in tons per hour. The contribution of classical thermal (coal or oil) units to this index can be defined as follows [1, 2]: where is the generated power (expressed in per unit) of the considered classical thermal (oil, coal,…) unit and are empirical coefficients related to this unit. Note that one of the main characteristics of the EII index comes from its U-shape in function of the generated thermal power. Finally, note that hydroelectric and wind generation units can be approximated as zero emission units [12] and do therefore not have any contribution to the EII.
The cost is generally a strictly growing evolution of the generated power [13]. In that way, the cost associated to classical thermal units ($/h) is generally characterized by a fuel cost square function of the generated power [13] and a starting investment cost quite limited in comparison with large hydroelectric units [14]: where considered are empirical coefficients related to the unit.
Concerning the hydroelectric parks cost ($/h), as already mentioned, the starting investment cost is generally quite high and the operating cost is quite reduced (as the fuel cost is zero). The latter is mainly due to maintenance activities and is practically modelled by a linear evolution depending of the generated power [14]: where is the generated power (expressed in per unit) of the hydroelectric unit and are empirical coefficients related to this unit.
Due to the multiple encouraging governmental policies in many countries, wind power has to be considered as a must run zero cost () generation unit in the proposed tool. Indeed, thanks to this hypothesis, wind power will be considered before classical units in order to cover the load and will be in agreement with the reality [15].
Practically, in order to express the cost constraint, a maximal hourly cost is defined [1, 2] and the constraint is defined as where, , and , are, respectively, the number of installed wind, thermal and hydroelectric units.
The load constraint associated to the optimization problem consists in exactly covering the load with the optimally dispatched classical and wind generation units. This constraint can thus be written as where , , and are, respectively, the hourly optimal generation of wind, thermal, and hydroelectric units.
Using expressions (1), (4), and (5), the minimization model of NOx emissions under cost and load covering constraints can be expressed as For efficiency and dimensioning reasons, each classical unit must also be constrained by lower () and upper () boundary values: where is the total number of classical units .
The developed optimization model is thus characterized by an objective function and constraints derivable and continuous in each point. Consequently, the optimal solution associated to each generated system state is calculated by the use of the fmincon solver from the Matlab Optimization Toolbox. This optimizer provides local optimization based on a gradient descent algorithm and can thus be limited in the case of problems with multiple local optima. However, the objective function (6) is always characterized by a single optimum and is thus perfectly suited for the application of the fmincon solver. Finally, note that this solver can also handle several types of constraints (linear or nonlinear, equality, inequality, etc.).
4. Calculated Indices and Useful Results Provided by the Monte Carlo Simulation Tool
4.1. Calculated Reliability Indices
In the Monte Carlo Simulation tool, the power system is thus modeled by specifying a set of “events,” where an event is a random occurrence that changes the system state. In the present study, the events recognized by the established program are the changes in load, the possible failure of a generating unit, or the variation of the power produced by a wind park. Each simulated system state is then defined in terms of available margin, which is the difference between the available (classical + wind) generation capacity and the load value. In order to proceed to this step, the total available system capacity must be superimposed on the load during each simulated hour to define the following.(i)Healthy state: the total available capacity is greater than the corresponding hourly load.(ii)risky state: the total available capacity is less than the corresponding hourly load value.Two complementary indices among which the LOLP [3, 4] are then defined as where and are, respectively, the total number of hourly healthy and risky simulated states; being the total amount of simulated years.
The number of hours per year during which the available total (classical + wind) generation cannot meet the load is defined as the Loss of Load Expectation(LOLE) index and is obtained by multiplying the LOLP index with the annual number of hours.
4.2. Optimization Results
Risky states cannot be optimized as the available generation is lower than the load. Consequently, the load constraint (8) cannot be met during those states. An indicator of risky states has thus been defined and set to 1 when such a state is simulated. Moreover, during a risky state, all the available generation has to operate at its maximal point in order to face the lack of production.
During healthy states, the previous indicator is set to 0 and it is worth to solve the optimization problem. In that way, in some healthy cases, the cost constraint cannot be met with the available generation and the load can thus not be fully covered with an adequate cost. Such cases are also characterized by another indicator set to 1 when the maximal cost constraint is violated. Healthy states presenting an optimized solution that verifies both load and cost constraints have both indicators set to 0. Moreover, thanks to the Monte Carlo environment, a focus can be made on the optimal configuration of each of those healthy states. This property of the developed tool could be useful for the electrical system operator in order to adequately dispatch the generation park when gaseous pollutants emissions are taken into account. Finally, at the end of the simulation, a mean EII index can be computed and permits, for example, to quantify the impact of wind power on the long-term reduction of gaseous pollutants emissions.
5. Simulation Results on the Slightly Modified RBTS Test System
In order to test the developed Monte Carlo simulation tool, a slightly modified version of the academic RBTS test system has been considered. In its initial version [7], the RBTS was consisting of 7 hydroelectric units and 4 classical thermal parks. Unfortunately, wind generation was not taken into account in that version of the RBTS. In order to introduce that renewable energy into the Roy Billinton Test System, reference [16] has consequently proposed the addition of wind power based on Weibull statistical distributions in the context of HLI adequacy studies. In the present paper, two wind parks with an installed capacity of 2 p.u., each, are considered. Both are subject to Weibull distributions with scale (A) and shape (B) parameters being, respectively, (; ) and (; ) [17]. Wind power is finally obtained by use of a variable speed classical power curve for each park [18]. In order to keep a sufficient number of simulated risky states and to test our simulation tool in every condition, the annual peak load value of the RBTS has been increased by 3.5 p.u. (annual peak load value = 22 p.u.) jointly to the introduction of wind generation. This new load value is afterwards modulated by use of the Belgian hourly and weekly modulation rates [9]. Finally, the empirical emission coefficients () associated to each thermal unit “” are extracted from [1, 2] and are listed in Table 1.
In the same way, the empirical cost coefficients associated to thermal and hydroelectric parks are taken from [1, 2, 14] and are, respectively, summarized in Tables 2 and 3.
Finally, note that the installed classical (thermal + hydroelectric) generation capacity is reaching 24 p.u. in our development and that the number of simulated system states is defined such as 10−4 accuracy can be reached on the computed indices.
5.1. Results Collected with the Developed Simulation Tool for the Modified RBTS with Independent Wind Parks
In this subsection, both wind parks are supposed to be entirely independent (results with entirely correlated wind parks are proposed in Section 5.2) and the obtained results are shown for one week in Figures 2(a), 2(b), and 3.
(a)
(b)
It can clearly be observed in Figure 2 that the optimal solution from the gaseous pollutants emissions point of view adequatly follows the load. However, in some cases like the one that can be pointed out during hour 14, the available “classical + wind” generation is greater than the load but the optimally dispatched “classical + wind” generation stays below the load. This result can be explained by looking at Figure 3. Indeed, it is due to the fact that the cost constraint ( $/h) cannot be faced with the available generation and the load that has to be covered.
The “life” of each unit of the considered system can be analyzed with the proposed program. Indeed, for a given unit, the optimal generation level calculated during each system state by the gaseous pollutants minimization model (Section 3) can be displayed. In that way, the optimal generation calculated for thermal unit 1 and the one computed for hydroelectric unit 1 are, respectively, plotted in Figures 4(a) and 4(b). On the same figure, indicators of outage state for the considered unit, LOLP system state, maximal cost constraint violation, and the available capacity of the unit are also proposed. In order to compare the behavior of classical units with the one of wind parks, Figure 4(c) presents, for wind park 1 (Weibull parameters: and ), the same quantities as the ones plotted in Figures 4(a) and 4(b). Several observations can be made on basis of Figures 4(a), 4(b) and 4(c).(i)First of all, it can be seen that the optimal wind generation is always equal to the available one. This can be explained by the fact that wind power is considered as a “zero cost zero emission must run” generation mean and that the installed wind generation (4 p.u. in total) is very limited in comparison with the considered load characteristics (minimal load value of 5.6 p.u.).(ii)Secondly, it can be observed that the optimal generation for hydroelectric unit 1 is nearly always equal to the available generation whereas the optimal generation of thermal unit 1 is often reduced in comparison with the available generation of that unit. This observation is explained by the fact that hydroelectric units are “zero emission” units. Moreover, hydroelectric unit 1 is quite cheap in comparison with the other classical parks (see Tables 2 and 3) and is thus very interesting when it comes to the minimization of pollutants emissions under cost constraints.(iii)Thirdly, some states during which the optimal generation is greater than the available one are computed. Those solutions must be rejected and have no meaning as the optimal generation cannot practically be greater than the available one. In fact, such “nonrealistic” solutions are reached when the maximal cost constraint could not be faced by the simulator and that this last one could consequently not find an adequate solution to the optimization problem.(iv)As expected, the available wind generation of wind park 1 is much more fluctuating than the one of both considered classical units. Indeed, whereas a classical unit is always (excepted during some outage states) able to operate at its installed capacity, wind generation is depending on the available wind speed.Note that the graphs of Figure 4 can be extended to each unit of the slightly modified RBTS test system and the same analysis can be made for each of them.
(a)
(b)
(c)
At the end of the simulation, the LOLP, LOLE (in hours/years), annual number of maximal cost constraint violation states, and mean EII (in tons/hour) are computed by the software and are displayed on the Matlab workspace.
5.2. Impact of Wind Generation on the Collected Reliability (LOLP, LOLE), Cost and Emission (Mean EII) Indices
Thanks to the developed simulation tool, it is possible to quantify the impact of wind generation on the electrical system adequacy and on gaseous pollutants emissions. In that way, the same wind parks as the ones defined in the introduction of Section 5 have been considered and the installed wind capacity has been increased from 4 p.u. to 8 p.u. (4 p.u. at each wind park). Moreover, in order to also quantify the impact of wind correlation on the obtained simulation results, two extreme scenarios have been computed: on the one side, both wind parks have been supposed to be entirely correlated and on the other side, the same wind parks have been considered as entirely independent. Finally, note that all the simulations have been realized with the same unchanged load and classical generation characteristics as the ones detailed for the modified RBTS test system.
From Table 4, it can be easily observed that the increase of wind penetration (entirely correlated or not) is improving both emissions and reliability indices. This result is quite logical as the load and classical generation characteristics are not changed while the installed wind capacity is increased. Therefore, as wind generation is considered as a zero cost nonpolluting generation mean, it will be preferentially used to cover the load and will consequently decrease the gaseous pollutants emissions. Moreover, as the load is not changed, the increase of wind penetration also increases the installed total (wind + classical) generation capacity and tends thus to decrease the risk of facing states of load nonrecovering.
On the other side, when the correlation between wind parks is analyzed, it can be seen that the power smoothing implied by the independence scenario between wind parks (an important production for one wind park does not necessarily involve a consequent production for the other as the parks are supposed to be independent) tends to reduce the risk of load nonrecovering and leads to reduced gaseous pollutants emissions in comparison with the “entire correlation” scenario. Indeed, when wind parks are entirely correlated, the simulated global wind power is much more fluctuant. Therefore, in some states during which the available wind generation is very low (for all the wind parks due to the entire correlation between them), the risk of not being able to cover the load is increased and much more classical thermal generation is needed with, as a consequence, an increase of the gaseous pollutants emissions.
6. Conclusions
In this paper, an original Monte Carlo simulation tool has been developed in order to optimally dispatch, from the gaseous pollutants emissions point of view, classical generation under load covering, and cost constraints. With the developed tool, the random behaviour of wind generation and of the hourly load can also be taken into account and represents therefore a consequent improvement in comparison with the quite limited existing analytical approach. Moreover, thanks to the proposed software, the long-term impact of wind generation on the electrical system adequacy and on gaseous pollutants emissions can be evaluated in different wind correlation scenarios and allows thus to better quantify the real influence of wind generation. Finally, it is thought that the developed tool could be a great help for electrical system operators in order to adequately dispatch the available classical generation in presence of wind power and under reliability, cost, and minimized gaseous pollutants emissions constraints.