Expectation–maximisation model for stochastic distribution network planning considering network loss and voltage deviation
2018; Institution of Engineering and Technology; Volume: 13; Issue: 2 Linguagem: Inglês
10.1049/iet-gtd.2018.5813
ISSN1751-8695
AutoresBaoxia Qi, Jiajia Chen, Yanlei Zhao, Pihua Jiao,
Tópico(s)Optimal Power Flow Distribution
ResumoIET Generation, Transmission & DistributionVolume 13, Issue 2 p. 248-257 Research ArticleFree Access Expectation–maximisation model for stochastic distribution network planning considering network loss and voltage deviation Baoxia Qi, Baoxia Qi School of Electrical and Electronic Engineering, Shandong University of Technology, Zibo, 255000 People's Republic of ChinaSearch for more papers by this authorJiajia Chen, Corresponding Author Jiajia Chen jjchen@sdut.edu.cn orcid.org/0000-0001-6027-7712 School of Electrical and Electronic Engineering, Shandong University of Technology, Zibo, 255000 People's Republic of ChinaSearch for more papers by this authorYanlei Zhao, Yanlei Zhao School of Electrical and Electronic Engineering, Shandong University of Technology, Zibo, 255000 People's Republic of ChinaSearch for more papers by this authorPihua Jiao, Pihua Jiao School of Electrical and Electronic Engineering, Shandong University of Technology, Zibo, 255000 People's Republic of ChinaSearch for more papers by this author Baoxia Qi, Baoxia Qi School of Electrical and Electronic Engineering, Shandong University of Technology, Zibo, 255000 People's Republic of ChinaSearch for more papers by this authorJiajia Chen, Corresponding Author Jiajia Chen jjchen@sdut.edu.cn orcid.org/0000-0001-6027-7712 School of Electrical and Electronic Engineering, Shandong University of Technology, Zibo, 255000 People's Republic of ChinaSearch for more papers by this authorYanlei Zhao, Yanlei Zhao School of Electrical and Electronic Engineering, Shandong University of Technology, Zibo, 255000 People's Republic of ChinaSearch for more papers by this authorPihua Jiao, Pihua Jiao School of Electrical and Electronic Engineering, Shandong University of Technology, Zibo, 255000 People's Republic of ChinaSearch for more papers by this author First published: 18 December 2018 https://doi.org/10.1049/iet-gtd.2018.5813Citations: 3AboutSectionsPDF ToolsRequest permissionExport citationAdd to favoritesTrack citation ShareShare Give accessShare full text accessShare full-text accessPlease review our Terms and Conditions of Use and check box below to share full-text version of article.I have read and accept the Wiley Online Library Terms and Conditions of UseShareable LinkUse the link below to share a full-text version of this article with your friends and colleagues. Learn more.Copy URL Share a linkShare onFacebookTwitterLinkedInRedditWechat Abstract The penetration of wind power is increasing in distribution network for reducing reliance on fossil fuels and covering continuously increasing demand for energy. However, it is argued that the forecasting error of wind power cannot be avoided even using the best forecasting approach. Therefore, in this study, a mean-variance-skewness based expectation maximisation (MVSbEM) model has been proposed by maximisation of the mean and skewness while simultaneously minimisation of the variance to obtain the optimal trade-off relationship between the profit and risk of distribution network planning (DNP) considering uncertain wind power integrated. In the MVSbEM model, the indexes of network loss, voltage deviation, and investment cost are concurrently taken into account under several kinds of actual operation constraints. In addition, the authors have made a full investigation on the MVSbEM by considering different forecasting errors, power factors of wind power, the different forecasting wind speeds, the number of wind turbines as well as the lines and substations upgrading. Furthermore, in order to reduce the computational burden, the Latin hypercube sampling method is used to sample uncertain wind speed. The feasibility and effectiveness of the MVSbEM model have been comprehensively evaluated on a modified IEEE 33-bus system. 1 Introduction Distribution network planning (DNP) is one of the most important issues in power systems, which aims to minimise the total generation cost while satisfying various operational constraints. In recent years, the wind power has been increasingly utilised and integrated into distribution network because of its sustainable, environment-friendly and lower marginal operation cost, which definitely contributes to DNP, i.e. for deducing generation cost and improving the flexibility of DNP [1–3]. However, the wind power still exhibits some obvious drawbacks that hinder its development such as uncertain and intermittent. In addition, the forecasting error of wind power cannot be avoided even using the best forecasting approach due to the chaotic nature of the weather system and climate change [4, 5]. As a consequence, in the distribution network, the wind power leads to the difficulty in determining the optimal dispatching solution for stochastic DNP in an uncertain environment. In the past few years, the methods solving this kind of stochastic DNP problem are classified into two groups, which are the fuzzy and probabilistic methods. In the fuzzy method, the wind power is deemed as the fuzzy variable by the fuzzy set theory [6, 7]. In this method, the dispatcher's attitude can be well reflected to choose the dispatching solution [8]. However, this method depends mainly on the selection of fuzzy functions, which will lead to the solution subjected to strong subjectiveness in some cases and could not adjust well to the actual needs [9, 10]. As for the probabilistic method, the uncertainty of wind power, wind speed or their corresponding forecasting errors are usually assumed to obey probability distribution. For instance, the authors of [11, 12] argued that the Weibull distribution is more effective for describing the probabilistic characteristic of wind speed during a long period (e.g. a year). Therefore, the Weibull distribution is not suited to be used for DNP. Moreover, the authors of [13, 14] have proven that the wind speed forecasting error follows the Gaussian distribution during a short period. Since the DNP is a non-convex and stochastic optimisation problem with complex inequality and equality constraints, it is necessary to consider both economic and reliability issues simultaneously in distribution network operation. In [15], a nonlinear interval number programming has been developed for assessing the uncertain optimisation problem by using a weighting method. In [14], the authors proposed a weighted mean variance model to measure the uncertain wind power and distribution loads. However, it is argued that the weighting method is inappropriate for a non-convex and stochastic power system optimisation problem [16]. In addition, the weighting coefficient applied to balance the return and risk tends to be very subjective and is usually assigned more than less arbitrarily by the analyst [17–19]. In addition, the Monte Carlo method used in [17, 20] does not consider the probability information of stochastic variables, and it needs a large scale repeated samplings. In view of this, Li et al. [21] proposed the mean-variance (MV) model and Latin hypercube sampling (LHS) method to assess the uncertain wind power from the perspective of the profit and risk. In this way, the stochastic method can address the stochastic dispatching problem by providing the optimal dispatching solution considering both economic risk and economic return of uncertain wind power. However, the simulation results proved that the MV model is used well to assess the wind power integrated problem only if the returns follow a normal distribution. Therefore, the MV model cannot be used directly for assessing the uncertain returns unless there are reasons to believe that the asset returns are symmetrically distributed around the mean [16]. Thanks to the MV-skewness (MVS) model [16], which was proposed to overcome the shortcoming of the MV model by considering both the profit and the risk as for the uncertain environment. Inspired by this point of view, in this study, an MVS based expectation maximisation (MVSbEM) model has been proposed for DNP with respect to the uncertain wind power penetration. The MVSbEM model can well reflect the relationship between the profit and risk of the DNP with uncertain wind power penetration, by maximising the expected value and skewness while minimising the variance of the objective function of the DNP. In addition, in the MVSbEM model, the indexes of network loss, the voltage deviation together with the investment cost are considered simultaneously. Moreover, the LHS method is adopted to sample uncertain wind speeds to reduce the computational burden of repeated samples. Furthermore, the simulated annealing particle swarm optimisation (SAPSO) algorithm is used to search for the optimal dispatching solution of the MVSbEM model. Compared to the previous research work, the main contributions of the study can be expressed as follows: (i) in order to cope with the uncertainty of the DNP with wind power integrated, an MVSbEM model is proposed, which aims to maximise both the expected profit and skewness while simultaneously minimising the risk. In the proposed model, the lines upgrading, the network loss, and voltage deviation are considered to reduce the network loss and improve the stability of the distribution network. (ii) To reduce the computational burden, the Fisher z transformation-based LHS method is derived in consideration that the cumulative density function of wind speed is irreversible and the detailed steps for simulating the uncertain wind speed is provided. (iii) To obtain the optimal dispatching solution of the DNP with wind power integrated, the SAPSO algorithm is developed to balance the local exploitation and global exploration of particle swarm optimisation (PSO). The remainder of this paper is organised as follows. Section 2 introduces the MVSbEM model considering the profit, risk, and skewness for assessing DNP. The optimisation algorithm of the SAPSO is discussed in Section 3. Section 4 carries out experiments and discusses simulation results for MVSbEM, MV based expectation maximisation (MVbEM), and conditional-value-at-risk based expectation maximisation (CVaRbEM) and consider the impact of distribution network upgrading (DNU). In the end, a conclusion is given. 2 MVSbEM model 2.1 DNP model The DNP problem with uncertain wind power penetration aims to minimise the overall cost, including the cost of network loss and voltage deviation, and the investment cost of wind power. Besides, it should also satisfy a set of operational constraints, which is presented in the following: (1) (2) (3)where (4)where F, , , , and indicate the overall cost, the cost of network loss, the cost of voltage deviation, the investment cost of wind power, and the save electricity cost, respectively, and are the total generation time and utilisation time of the wind turbine (WT), respectively, is the power factor of wind power, r and m are the utilisation rate and the service life of WT, is the maximum allowable deviation, denotes the operation and maintenance cost of WT, , , , and represent the sale price, the price index of voltage deviation, the integrated capacity cost of wind power and feed-in tariffs, respectively, h and g represent equality constraints and inequality constraints, respectively, X and U indicate state variables and decision variables, respectively, is the uncertain wind power. The total network loss of the distribution network is (5)where L denotes the total number of system branches, is the resistance of the i th branch, , , and are active power, reactive power, and bus voltage of the i th branch, respectively. The voltage deviation of all load buses can be mathematically expressed as follows: (6)where N represents the total number of load buses and is the reference voltage of the i th load bus, which is taken as 1 p.u. and is the voltage of the i th load bus. 2.2 MVSbEM model The MVS model proposed by Tarja Joro and Paul Na [22] is an extension of the traditional MV model in dealing with portfolio optimisation problems. In addition, the feasibility and effectiveness of the MVS model have been approved in the electricity market portfolio for dealing with investment returns of non-normal distribution [16, 23]. Moreover, the skewness is an important index in investors’ decision-making processes because of the skewness on behalf of the risk level [16]. Therefore, this study proposes a multi-objective MVSbEM model to optimise problem of the DNP from the perspective of the profit, risk, and skewness. In this study, the profit M is defined as the decrease of the overall cost of the DNP when the wind power is integrated into a distribution network, which is written as [16] (7)where denotes the return with respect to the i th sample of wind power, indicates the overall cost corresponding to (1) without considering wind power penetration, indicates the overall cost corresponding to (1) with respect to the i th sample of wind power, and represents the total number of samples. The risk is manifested by profits brought by wind power integrated, which means the extent of fluctuations between the actual profits and their expectation under different scenarios of uncertain wind power. The risk of the DNP with wind power penetration can be calculated as follows [16]: (8)In addition, studies by Arditti show that investors prefer positive skewness in the uncertain environment [24]. Harvey and Siddique further argue that the skewness preference is consistent with the Arrow–Pratt notion of risk aversion [25]. Therefore, in this study, the skewness is also formulated to assess the profit (9)In this subsection, we also give a necessary analysis for the MVSbEM model. First, let denotes the utility function and denotes the utility value of the i th wind power integrated. Then, consider Taylor's series expansion of around the expected profit M can be represented as follows [16, 22]: (10)Then taking the expectation of (10) yields (11)where represents the k th central moment of the i th return . Therefore, the skewness is important in investors’ decision making due not only to empirical but also theoretical [16, 25]. This conclusion is also consistent with the notion of behavioural theory that people attach different values to (monetarily) equal losses and gains when making decisions under risk [26]. Thus, it is necessary to consider the skewness in DNP with uncertain wind power. In this study, let us only consider the third moment approximation of the expected utility of (11), because higher moments are not very meaningful for power systems planning [16]. Therefore, the objective function of the DNP based on the MVSbEM is formulated as follows: (12) (13) (14)where the utility function is defined as [16]. Moreover, in order to fairly compare the models, the objective functions of the MVbEM model and CVaRbEM model are, respectively, given by (15) and (16). In the MVbEM model, the optimisation objective consists of the profit and the risk , therefore the objective function of the DNP based on MVbEM can be expressed as (15)As for the CVaRbEM model, the objective function can be calculated as follows [27]: (16)where y represents the predefined expectation profit and is the confidence level. 2.3 Wind power characterisation 2.3.1 Latin hypercube sampling Li et al. [14] argued that Gaussian distribution is more effective than the Weibull distribution for describing the probabilistic characteristic of wind speed during a short period, and has been widely used in DNP. Therefore, this study uses the Gaussian distribution to describe the uncertain wind speed v (17)where and represent the forecasting value and forecasting error of wind speed, respectively. The LHS method is an equal probability sampling method based on the inverse function of the cumulative density function of randomness variables [28, 29]. Compared with the Monte Carlo sampling method, the LHS method does not require a large number of samples. However, the cumulative density function of Gaussian distribution is not integrable, and the inverse function cannot be calculated directly, therefore the wind speed cannot be sampled equally by LHS. Thanks to the Fisher z transformation proposed by Fisher [30], it can be used to approximate the cumulative density function of Gaussian distribution. The wind speed sampling procedure can be divided into the following four steps: firstly, calculate the cumulative density function of wind speed using (18); secondly, the approximate function of (18) is received by the Fisher z transformation as shown in (19); thirdly, by calculating the inverse function of (19), the expression of v is obtained shown by (20); finally, the wind speed can be obtained by equal probability sampling. The cumulative density function of v can be expressed as (18)By using the Fisher z transformation, (18) can be calculated as follows: (19)where and . From (19), we can obtain the inverse function of (18) (20)where indicates the cumulative probability. Assuming that the sample size is , y is divided into equivalent lengths, where each length is . The value of the sampling probability is intermediate value in the interval, therefore the wind speed v can be calculated by (20). 2.3.2 Wind power output mathematical model The mathematical model between active wind power and wind speed v is given as follows [31]: (21)where is the rated active power of the WT, , , and v are the rated wind speed, the cut-out wind speed, the cut-in wind speed, and the actual wind speed, respectively. 2.4 System constraints The power flow calculation of the distribution network is an important part for power system analysis, which contributes to analysing and studying the distribution loss and voltage deviation. In this study, the forward push-back method is used to calculate the distribution network power flow considering WT integration. Since the convergence of the forward push-back method is not affected by the ratio of the high resistance and the reactance in the distribution network. In addition, it can take full advantage of the feature that the distribution network has a unique path from any given bus to the source bus to the correct voltage and current. Therefore, the forward push-back method is widely used in distribution network power flow calculation [32]. 2.4.1 Bus voltage constraints The voltage drop in a line is given as follows: (22)where and are the voltage of the i th and j th load bus, respectively, , , and are the resistance, reactance, and impedance of a line between bus i and bus j, respectively; are the upper and lower voltage limits of the i th load bus, respectively. 2.4.2 Current constraints (23)where is the current of the i th branch, is the upper limit current of the i th branch. 2.4.3 Power balance constraints The requirements of active and reactive power balance at each bus are as follows: (24) (25)where and L is the total number of buses, and are active and reactive load demand at bus i, respectively, and are active and reactive power of WT of the i th bus, respectively, Bij and Gij are transfer susceptance and conductance, respectively, between buses i and j. 2.4.4 Wind power constraints (26)where and are the upper and lower active power limits of wind power, respectively, and are the upper and lower reactive power limits of wind power, respectively. 3 SAPSO algorithm and implementation steps In recent years, PSO algorithm has aroused intense interest in power system planning due to its flexibility, robustness, and versatility in seeking a globally optimal solution without requirements on the optimisation problems. However, it is argued that the PSO algorithm still faces the challenge in premature convergence for non-convex and stochastic power system optimisation problem [33]. In view of this, a hybrid algorithm named simulated annealing PSO (SAPSO) has been developed for DNP. In the search process, the SAPSO algorithm accepts not only the better but also the worse neighbouring solutions with a certain probability. The probability of accepting a worse solution is larger at higher initial temperature. As the temperature decreases, the probability of accepting worse solutions gradually approaches zero. This feature means that the SAPSO algorithm makes it possible to jump out of a local optimum to search for the global optimum [34]. The SAPSO implementation steps for DNP are given as follows: Step 1: Set up the basic parameters for each particle, e.g. population size and simulated annealing constant . Set , the maximum number of generations . Step 2: Randomly initialise the particle and velocity in the population. Step 3: Calculate the fitness values of each particle according to (12). Step 4: Store the current particle as the initial solutions as , and then store a globally optimal solution from the current particle as . Step 5: Determine the initial temperature T according to (12) and (27) (27) Step 6: According to (28), determine the fitness value TF of each particle at the current temperature (28)where T denotes the current temperature. Step 7: Utilise the roulette strategy to determine the global optimal replacement for do if ; = ; end if end for Step 8: Update each particle according to (29), recalculate the fitness value according to (12) and update and (29)where , , and are learning factors, and are evenly distributed random numbers between 0 and 1. Step 9: . Step 10: , if , then go to Step 11. otherwise, go to Step 6. Step 11: Stop the procedure. 4 Numerical results In this study, in order to validate the effectiveness of the proposed MVSbEM model, simulations are conducted on the modified IEEE 33-bus test system considering uncertain wind power integrated. This test system is a 12.66 kV radial distribution system, which consists of 33 buses and 32 branches with a total load of 5.08 MW and 2.547 MVar. The load data and line data are given in [35]. The rated active power of wind power is set to 2 MW. The rated wind speed, the cut-in wind speed, and the cut-out wind speed are set to 12.5, 4, and 20 m/s, respectively. The number of sampling size of wind power is 100. Additionally, the basic parameters of SAPSO are set as follows: , the group size N is 100, annealing constant is 0.5 and is set to 30. The WT operation parameters are shown in Table 1. Table 1. WT parameters Parameters maximum annual power generation time 8760 h maximum annual utilisation hours 3000 h power factor 0.95 utilisation rate 10% service life 75,000 h annual operation costs $31.30/kW sales price $0.12/kW voltage deviation index price $72374.24 [36] loss capacity cost $986.15/kW feed-in price $0.08/kW line upgrade fixed cost $31,100 [37] line capacity cost $1.50×105/MVA [38] substation upgrade fixed cost $30,153 [39] 4.1 Comparison of MVSbEM, MVbEM, and CVaRbEM To test the effectiveness of the proposed MVSbEM, the optimal number of WTs has been investigated by MVSbEM. The optimal expectation , profit , risk , skewness , network loss and total voltage deviation obtained on the different WT numbers are shown in Table 2, respectively. Table 2. Numerical results obtained by MVSbEM model under a different number of WTs penetration E M V S 0 — — — — 0.3730 0.3332 1 0.4751 0.7217 0.0157 −1.0849 0.1723 0.0410 2 0.7302 1.2156 0.0091 −2.1621 0.1695 0.0158 3 0.7042 1.1461 0.0100 −1.7352 0.2934 0.0203 It can be seen from Table 2, the WT integration can reduce the network loss of the distribution network, more accurately, network loss reducing 0.0201, 0.2035, 0.0796 MW than that without wind power penetration, respectively. When two WTs integrated into the distribution network, it can obtain the optimal expectation, which the expectation increases 53.69 and 3.69% than that of one WT and three WTs integrated, respectively. Therefore, a high penetration level of wind power is not beneficial for system operation due to the index of risk. In addition, when two WTs integrated, the voltage deviation has reduced 61.46 and 22.17% than that of one WT and three WTs integrated, respectively. Thus, in the following discussion, we only consider two WTs integrated DNPs. From Table 2, we can also find that the skewness does exist, which means that the MVSbEM is more appropriate for DNP. The simulation results also demonstrate that the MVSbEM can obtain a trade-off solution between the risk and profit. Moreover, the voltage profile is provided by Fig. 1. It is apparent from the figure that the voltage of each load bus satisfies the voltage constraints. Fig. 1Open in figure viewerPowerPoint Voltage profile with a different number of WTs integration Moreover, to validate the effectiveness of the SAPSO algorithm, the SAPSO algorithm is compared with the PSO algorithm. The comparisons of the convergence speed between these two algorithms are provided by Fig. 2. It is quite apparent from the figure that the proposed SAPSO performs better in global searching in comparison with PSO. Additionally, the optimal expectation, profit, risk, and skewness obtained based on the SAPSO and PSO are shown in Table 3. It is obvious that the expectation and the profit obtained by SAPSO are higher than those obtained by PSO, increasing by 1.11% expectation value and 1.66% profit. In addition, the risk obtained by SAPSO has decreased by 37.67% than that of the PSO. Therefore, the SAPSO is more efficient than the PSO in dealing with DNP. Fig. 2Open in figure viewerPowerPoint Comparison of convergence rates between SAPSO and PSO Table 3. Simulation results obtained by MVSbEM model between SAPSO and PSO Algorithm E M V S SAPSO 0.7302 1.2156 0.0091 −2.1621 0.1695 0.0158 PSO 0.7222 1.1958 0.0146 −2.0939 0.1698 0.0159 Furthermore, in order to make a comprehensive comparison of the proposed MVSbEM, in this subsection, the optimal solution of MVSbEM is compared with that of the MVbEM model and CVaRbEM model under two WTs integrated. It should be mentioned that in the MVbEM model the objective function is similar to MVSbEM by ignoring the skewness. The optimal solutions of MVSbEM, MVbEM, and CVaRbEM are listed in Table 4. It is evident that MVSbEM provides the highest expectation than MVbEM and CVaRbEM. In addition, the profit obtained by MVSbEM is larger than that of MVbEM and CVaRbEM, which increases by 0.32 and 1.34%, respectively. Moreover, the wind power penetration level of MVSbEM increases by 0.66 and 107.20% than that of CVaRbEM and MVbEM, respectively. This demonstrates that the CVaRbEM seems to be too conservative in dealing with uncertain wind power. It can also be seen from the table that the skewness exists, which means the returns do not follow a normal distribution and the MVbEM model is inappropriate for measuring DNP. Therefore, the proposed MVSbEM model is effectiveness in a deal with the problem of the DNP. Table 4. Numerical results obtained among MVSbEM, MVbEM, and CVaRbEM models Model E M V S MVSbEM 0.7302 1.2156 0.0091 −2.1621 0.9990 1.4580 MVbEM 0.7292 1.2117 0.0092 −2.1273 0.9829 1.4580 CVaRbEM 0.7230 1.1995 0.0096 −2.1016 0.9531 0.2327 4.2 Comparison of MVSbEM, CVaRbEM, and MVbEM under different forecasting wind speeds To further validate the effectiveness of MVSbEM, the hourly optimal wind power, expectation and risk under different forecasting wind speeds are provided in Table 5. Here, the deviations of each hour are set to . From the table, it is clear that in most cases the MVSbEM outperforms MVbEM and CVaRbEM in finding the trade-off solutions between the expectation and risk. This is due to the fact that in MVSbEM, mean , variance and skewness are considered in the objective function to maximise M and S while simultaneously minimise V. Therefore, in most cases, the MVSbEM can obtain higher values of M and S, and obtain lower values of V. Note that in MVSbEM, the variance is developed to assess the risk of DNP considering uncertain wind power penetration. To further investigate the proposed model, the risk measure proposed in CVaR is used to assess the risk of the optimal solution obtained by MVSbEM. In Table 5, it should be mentioned that the CVaRbEM is simulated on the risk level of 0.9. The risk level of the optimal solution obtained by MVSbEM is 0.88, which means that the risk gained by MVSbEM is higher than that of CVaEbEM by using the risk measure of CVaRbEM. Therefore, the effectiveness of different risk measure models depends on the uncertainty nature and available data about the uncertain parameters of the model, which cannot be easily compared with each other. In this study, the proposed MVSbEM model does not need a predefined risk level, and can automatically balance the risk and profit in the optimisation process. Table 5. Numerical results obtained by MVSbEM, MVbEM, and CVaRbEM under different forecasting wind speeds Hour MVSbEM MVbEM CVaRbEM E V E V E V 1 7.5 1.0175 7.2 1.0182 0.6709 0.0155 0.6527 0.0172 0.6574 0.0167 2 6.6 0.8000 7.0 0.9400 0.6284 0.0190 0.6241 0.0193 0.6224 0.0194 3 7.2 0.9018 7.8 1.2600 0.6939 0.0127 0.7073 0.0119 0.6892 0.0129 4 6.4 0.7300 6.7 0.8292 0.5953 0.0217 0.5913 0.0219 0.5933 0.0214 5 6.8 0.8700 6.6 0.8000 0.6081 0.0208 0.6077 0.0208 0.6081 0.0208 6 6.8 0.8700 6.9 0.9100 0.6318 0.0188 0.6318 0.0188 0.6291 0.0190 7 7.3 1.0600 7.1 0.9800 0.6691 0.0157 0.6624 0.0163 0
Referência(s)