Dynamic Optimal Power Flow (DOPF) is a typical complex multi-constrained non-convex non-linear programming problem when considering the valve-point effect of conventional generators. Moreover, it is mixed integer when considering the discreteness of reactive compensation devices for power system with Fixed Speed Wind Generators (FSWGs). DOPF model in FSWGs integrated power system is established in this study and then a hybrid Taguchi-Particle Swarm Optimization (TPSO) algorithm is proposed for solving the established DOPF model. This hybrid algorithm combines the well-known Particle Swarm Optimization (PSO) algorithm with the established Taguchi method which has been a important tool for robust design. This study clearly presents the improvements obtained despite the simplicity of the hybridization process. The Taguchi method is run only once in every iteration and therefore does not give significant impact in terms of computational cost. The method creates a more diversified population, which also contributes to the success of avoiding premature convergence. IEEE 39-bus system is used to illustrate the effectiveness of the proposed method compared with those obtained from PSO algorithm. The test results show that the proposed method is effective and has a certain practicality.
PDF Abstract XML References Citation
How to cite this article
Wind energy is the worlds fastest growing renewable energy source. With the increasing levels of wind generator penetration in modern power systems, one of major challenges in the present and coming years is the optimization control, such as optimal power flow including wind farms (Li and Chen, 2008).
The well-known advantages of FSWGs (fixed speed wind generators) are it is robust, easy and relatively cheap for mass production, there are still widely use in power system with FSWGs. Since FSWGs always draw reactive power from the grid, in most cases, capacitors are connected in parallel to the generator to compensate for the reactive power consumption, the P-bus power flow model of FSWGs bus is proposed in this study.
In this study, the problems of DOPF (Dynamic Optimal Power Flow) including FSWGs are researched. The expectation model of wind generators' active power outputs is adopted. DOPF is a method to schedule the online generator outputs with the predicted load demands over a certain period of time so as to operate an electric power system most economically. Normally, it is solved by dividing the entire dispatch period into a number of small time intervals, then a static economic dispatch has been employed to solve the problem in each interval (Pothiya et al., 2008; Uturbey and Costa, 2007). In this study, DOPF model, which takes all conventional units cost minimum as the objective function and takes the whole time and the inherent relations of different intervals into account in FSWGs integrated system, is established.
DOPF is a typical complex multi-constrained non- convex non-linear programming problem when considering the valve-point effect of conventional generators. Moreover, it is mixed integer when considering the discreteness of FSWGs reactive compensation devices. However, both lambda-iterative and gradient technique methods in conventional approaches to the problems are calculus-based techniques and require a smooth and convex cost function and strict continuity of the search space.
On intelligent algorithms, easy to handle discrete variables, no special requirements for the objective function, no problems of dimension disaster, strong global search capability and the greater probability in finding out the global optimum solution are the prominent advantages of these algorithms. However, the computational accuracy and convergence rate are not high for most intelligent algorithm.
In current study, we see a myriad of improved algorithms. Many of these algorithms incorporate some improvement strategies in the algorithms itself while some improved algorithms are the result of hybridizing two algorithms (He et al., 2008; Yuan et al., 2009). The Taguchi method is an established approach for robust design, applying the idea from statistical experiment design for evaluating and improvements in products, processes and equipment (Taguchi et al., 2000, 2005; Ross, 1989). The key to this concept is to improve the quality of a product by minimizing the effect of the causes of variation without the elimination of the relevant causes. The two major tools used in the Taguchi method are: (1) Signal-to-Noise Ratio (SNR) which measures quality and (2) orthogonal arrays which are used to study many design parameters simultaneously (Taguchi et al., 2000, 2005; Ross, 1989). This study adopts the concept of Taguchi method into PSO algorithm to improve its convergence which is so called Taguchi-PSO (TPSO). By applying the algorithm to dealing with IEEE 39-bus system, compared with PSO algorithm, the experimental results show that the algorithm is an effective way to solving DOPF problems in FSWGs integrated power system because of higher computational accuracy, better robustness and convergence performance.
IEEE 39-bus system is used to illustrate the effectiveness of the proposed method compared with those obtained from PSO algorithm. The test results show that the proposed method is effective and has a certain practicality.
DOPF MODEL IN FSWGS INTEGRATED POWER SYSTEM
Due to the random variation of the wind velocities and load demands, it is difficult to research the DOPF in the power system including wind farms. For simplifying this problem, the dividing-stage strategy is adopted in this study. According to the wind velocity forecasting curves and the load forecasting curves in the planning horizon, the expectations of wind generators' outputs and the load demands at dispatch interval can be calculated.
Constraints: Constraints include equality and inequality constraints. The equation constraint is the power flow formulation constraint while inequality constraints including generator power output, ramp rate and bus voltage are as in Eq. 1-4. For the FSWGs bus, absorbed reactive power from the power system is bounded as in Eq. 4. The constraints of real power generation limit and the ramp rate are taken into account as in Eq. 1.
where, Pi,min and Pi,max are the maximum and minimum limits of the power generation of unit i, Pit is the real power output of unit i at the tth interval, Pit-1 is the real power output of unit i at the t-1th interval; URi is the up-ramp limit of the ith generator (in units of MW/time-period) and DRi is the down-ramp limit of the ith generator (in units of MW/time-period); ΔT is time interval, Ngen is the number of conventional generating units and N is the number of system buses (excluding slack bus); Vit is the voltage magnitude output of bus i at the tth interval; QtGi is the reactive power output of conventional generating unit i at the tth interval; Qabt is the absorbed reactive power of wind farm bus at the tth interval; max is the maximum value of the variable, min is the minimum value of the variable.
After calculating the power flow, the state variables, power loss and real power output of the slack bus generator corresponding to the current control variables are available. The real power output of the slack bus generator will be set to the limit if it violates the limit. After handling overlimit of the real power output of the slack bus generator, the system power balance constraints as in Eq. 5 must meet, otherwise adding Eq. 5 as penalty terms to the objective function to form a generalized objective function. Details of the generalized objective function used in this study are given in Eq. 8.
where, ΔPt is the unbalance of the real power at the tth interval, Ngen-1 represents the number of conventional generating units excluding the slack bus, Pslt is the real power output of the slack bus generator after handling its overlimit at the tth interval, Plst is the total power loss at the tth interval, Pldt is the total load expectation at the tth interval, Ptw,av is the expectation of wind generators' real power outputs at the tth interval.
Objective function: Due to the fact that wind generation does not consume the fuel, the utility must purchase all the energy produced by wind generating units. Consequently, the objective is to minimize the following total incremental fuel cost function F associated to Ngen dispatchable units for T intervals in the given time horizon, subject to the above-mentioned equality and inequality constraints.
The inclusion of valve-point loading effects makes the modeling of the fuel cost function of the unit more practical. This increases the non-linearity and local optima in the solution space. Also, the solution procedure can easily trap in the local optima in the vicinity of optimal value. The fuel cost function of the ith unit F (Pit) with valve-point loadings are represented as follows:
where ai, bi and ci are cost coefficients and ei, fi are constants from the valve-point effect of the ith generating unit.
Evaluation function: We must define the evaluation function for evaluating the fitness of each individual in the population. In the most of the nonlinear optimization problems, the constraints are considered by generalizing the objective function using penalty terms.
FSWGs always draw reactive power from the grid. According to the principle of the nearest reactive power compensation, the required reactive power of wind farm doesnt absorb from the system as much as possible, but mainly from the wind farm reactive compensation devices. Because reactive power compensation devices are group switched capacitors in wind farms, the maximum of absorbed reactive power of wind farm from the grid is set to one set of group switched capacitors which is considered by generalizing the objective function using penalty terms. Considering efficiency maximization of wind farm reactive power compensation devices, the injected reactive power of wind farm bus is not constrained.
To sum up, the above problems are generalized as follows:
where, KV, KQ, KC and KD are variable overlimit penalty coefficients, Vit is the voltage magnitude of bus i at the tth interval (excluding the slack bus and PV bus ); QtGi is the reactive power output of generator i at the tth interval; Qabt is the absorbed reactive power of wind farm bus at the tth interval; Vilim and QGilim denote the violated upper or lower limits; Qablim denotes the violated lower limit that is defined as:
where, Qabmax is the upper limit value of the absorbed reactive power of wind farm bus which is set to 1Mvar in the study.
In this study, KV, KQ and KC are set to 1,1,1 respectively. Because the unbalance of the real power ΔP is hard to meet, an adaptive penalty function to handle penalty coefficient KD is adopted, where, k is the algorithms current iteration number; β is a relative violated value of the constraints, γ is a multi-stage assignment value, a is the power of the penalty value.
Meanwhile, several experiments have been done in order to obtain the penalty parameters. In this study, if β≤1 then α = 1, otherwise α = 2. Furthermore, if β≤0.001, then γ = 1, else, if β≤0.01 then γ = 10, else, if β≤0.1 then γ = 30, else, if β≤1 then γ = 100, otherwise γ = 300.
THE CONCEPT OF TAGUCHI METHOD
Taguchi methods were developed by Genichi Taguchi as an aid to improving the implementation of total quality control in Japanese industry. He described his methods as 'quality engineering and, in fact, the methods take an engineering approach to the understanding of process information. The methods are based on the design of experiments to provide near optimal quality characteristics for a specific objective. Taguchi defines quality in a negative manner as the loss imparted to society from the time the product is shipped (Taguchi et al., 2000).
The method is based on several statistical concepts which have proven to be valuable tool in the subject of quality improvement. Many Japanese manufacturers have applied this approach for improving product and process quality with unprecedented success. Taguchi essentially utilizes the conventional statistical tools, which has been simplified by identifying a set of stringent guidelines for experiment layout and the analysis of results. Recently Western industries have begun to recognize Taguchis method as simple but highly effective approach in improving product and process quality. Taguchi method has a wide range of applications in other areas (Cornejo-Mazon et al., 2008; Nava-Arenas et al., 2008). Two major tools used in the Taguchi method are the orthogonal array and the SNR. This study incorporates the two-level orthogonal arrays (Table 1), whereby the two levels are denoted by number 1 and 2 in the Table 1. The basic concept of the structure and usage of two-level orthogonal arrays are briefly described here. This is adequate in the context of describing the development of the hybrid TPSO algorithm.
|Table 1:||L8 (27) orthogonal array|
The general form of two-level standard orthogonal arrays can be represented by the following mathematical term:
where, n = 2k, n is the number of experimental runs, k is a positive integer which is greater than 1, 2 is the number of levels for each factor and n-1is the number of columns in the orthogonal array. The size of the orthogonal array to be used depends on the problem instance. At a time only two particles (two factors) are executed at two different levels for n number of experiments. The two-level standard orthogonal arrays most often used in practice are L4(23), L8(27), L16(215), L32(231), L64(263), L128(2127) and so on.
Orthogonal arrays are readily composed and are available from many texts (Ross, 1989). The way that they are constructed is to have (1) each level of every factor appear the same number of times in every column of the array and (2) each combination of factors between any two columns appears the same number of times. For any given values of (number of levels) and (number of factors), (number of tests) is determined as the smallest number that satisfies (1) and (2) above.
In this study, a two-level orthogonal array is used. There are m factors, where m is the number of design factors (variables) and each factor has two levels. If, m < n-1, only the first columns are used, while the other n-1- m columns are ignored. For example, there are 6 factors with two levels for each factor. We only need 6 columns to allocate these factors and L8(27) is enough for this purpose because it has 7 columns. Table 1 shows an orthogonal array L8(27). The number on the left of each row is called the run number or experiment number and runs from 1 to 8. In fact, each row represents a particle with only first 6 columns. In our methodology, we randomly choose two particles in every iteration. For this example, perform 8 experiments based on Table 1 whereby the number (either 1 or 2) of each column is the dimension from first and second particle, respectively.
Many designed experiments use matrices called orthogonal arrays for each experimental run and for analyzing the data. The array is called orthogonal because all columns can be evaluated independently of one another, as shown in Table 1. As one column represents a particle, therefore each column has its own fitness value (f1, f2, f8). The primary goal in conducting this matrix experiment is to determine the best or the optimal level for each factor. In taguchi based PSO, two tools of taguchi method i.e. orthogonal array and SNR are used. The basic idea is to identify the one best particle in an orthogonal array experiment, which consists of a set of experiments conducted at two different levels.
Besides two-level orthogonal arrays, there exists three-level and four-level orthogonal arrays which can provide better precision to the specified process. The detailed description of how these orthogonal arrays are formed can be found in the books written by Ross (1989), Roy (1990) and Taguchi et al. (2000).
OVERVIEW OF PARTICLE SWARM OPTIMIZATION
Kennedy and Eberhart (1995) first introduced the Particle Swarm Optimization (PSO) method. It is one of the optimization method categorized in the family of evolutionary computation. The method has been found to be robust in solving real-world problems featuring non-differentiability, high dimension, multiple optima and non-linearity which is widely used in various fields (Heng et al., 2006; Sutha and Kamaraj, 2008; Lai-Jun et al., 2009; Yang et al., 2009; Zakermoshfegh et al., 2008). PSO algorithm is a model that mimics the movement of individuals (fishes, birds, or insects) within a group (school, flock and swarm). Similar to GA, a PSO consists of a population refining its knowledge of the given search space. PSO is inspired by particles moving around in the search space. The individuals in a PSO thus have their own positions and velocities. Each particle moves in the search space with velocity which is dynamically adjusted and balanced based on its own best movement (pBest) and the best movement of the group (gBest).
Instead of using evolutionary operators such as selection, mutation and crossover, each particle in the population moves in the search space with velocity which is dynamically adjusted and all particles are assumed to be of no volume. Each particle keeps track of its coordinates in the search space, which are associated with the best solution it has achieved so far. This value is known as pBest. Another best value that is tracked by the global version of the particle swarm optimizer is the overall best value or the best solution in the population is called gBest.
The PSO concept consists of, at each time step, changing the velocity of each particle toward its pBest and gBest solutions. The movement is weighted by a random term, with separate random numbers being generated toward pBest and gBest values. For example, the ith particle consisting d dimensions is represented as Xi = (Xi,1, Xi,2, Xi,3, , Xi,d). The same notation applied to the velocity, Vi = (Vi,1, Vi,2, Vi,3, , Vi,d). The best previous position of the ith particle is recorded and represented as pBesti = (pBesti,1, pBesti,2, pBesti,3, pBesti,d). In the case of minimization that we consider in this study, the value of pBesti with lowest fitness is known as gBest. The modification of velocity and position can be calculated using the current velocity and the distance from pBesti to gBest as shown in the following formulas
where, i ∈ 1...N, j ∈ 1...d, t ∈ 1...T with N is the number of population size, d is the number of dimension and T is the number of maximum generation.
The position, X of each particle is updated for every dimension for all particles in each iteration. This is done by adding the velocity vector to the position vector, as described in Eq. 12 above. In Eq. 11, w is known as the inertia weight. This parameter was introduced by Shi and Eberhart (1998) to accelerate the convergence of PSO. Suitable selection of w provides a balance between global and local explorations, thus requiring less iteration on average to find sufficiently optimal solution. Low values result in particles moving in the region far from the optimum value before being tugged back. On the other hand, high values result in abrupt movement toward target regions.
The parameters ρ1 and ρ2 are set to constant values, which are normally given as 2.0 whereas r1 and r2 are two random values, uniformly distributed in [0,1]. The constants, ρ1 and ρ2 represent the weighting of the stochastic acceleration terms that pull each particle toward pBest and gBest positions.
TAGUCHI-PARTICLE SWARM OPTIMIZATION
Here, we illustrate with a simple example the concept of Taguchi-PSO. The Taguchi method is performed after the entire population is updated via PSO equations. Following this Taguchi method, further elaborations are given in a list of steps below. All the pseudocodes are written in programming style for L8(27),other two-level standard orthogonal arrays such as L4(23), L16(215), L32(231), L64(263) and L128(2127) etc., can be analogized.
|•||Step 1: Choose two particles randomly from the current population, noted as P1 and P2|
|•||Step 2: Create 8 new particles; noted as T1, T2 T8 from P1 and P2. (Each particle represents one experiment as shown in Table 1)|
|•||Step 3: Assign all the relevant dimensions (columns A-G) to the 8 new particles T1, T2 T8 from P1 and P2 based on these columns in orthogonal array given in Table 1. Number 1 means the relevant dimension of new particle is copied from P1 whereas number 2 means the relevant dimension is from P2|
|•||Step 4: Evaluate the fitness of T1, T2 T8 (refer to last column of Table 1)|
|•||Step 5: Calculate the SNR of various factors (A, B, , G) defined as:|
For instance, for the case of factor = B and level = 1, we have:
The same calculation is applied to the similar factor with level = 2:
|•||Step 6: Assign the relevant dimension to the optimal particle Po based on the following rule, taking factor =A in this case|
|•||Step 7: Replace the worse particle among the two randomly chosen if the optimal particle is better than both. Also, update the gBest as follows|
|•||Step 8: Continue with the next iteration until maximum generation is reached|
Step 6 above is important as a better dimension value is assigned to optimal particle Po from two randomly chosen particles. This technique provides an interesting philosophy of vertical optimization which is very useful for many real world optimization problems. For better performance in regards to convergence, the two particles in step 1 of Taguchi method are randomly chosen from pBest vector. As of such, the Taguchi method will therefore improve the dimensions value among pBest, an effective way leading to faster convergence and better accuracy. We employ the similar methodology in this work for better performance.
For DOPF problem including FSWGs, there are T dispatches by Ngen-1 conventional generating units and many wind farms (for simplicity, one wind farm is adopted, many wind farms may be analogized). A particle array of control variable vectors is:
where, P is particle vector, g is the number of population particles, Pnt is the real power output of nth generating unit at the tth interval, Qct is the number of group switched capacitors in wind farm at the tth interval:
|•||Step 1. Initialization: For the complete g population particles, the candidate solution of each individual particle is randomly initialized within the feasible range in such a way that it should satisfy the constraint given by Eq. 1|
|•||Step 2. Power flow and fitness calculation: Through the power flow calculation including wind farms, the state variables, power loss and real power output of the slack bus generator corresponding to the current control variables have been able to get. The real power output of the slack bus generator will be set to the limit if it violates the limit. After handling overlimit of the real power output of the slack bus generator, the system power balance constraints as in Eq. 5 must meet, otherwise adding Eq. 5 as penalty terms to the objective function to form a generalized objective function. In this study, Eq. 8 is used as the fitness or evaluation function. This is a generalized fitness function used to evaluate the fitness of the candidate solution of each particle individual. Also, record the particles position with the global best fitness as gBest, record the current position of each particle as its current pBest. Set the iteration count k = 0|
|•||Step 3. TPSO method: Solve DOPF in FSWGs integrated power system using TPSO method|
|•||Step 4. Update velocity and position: If k =kmax, iteration stop; otherwise, update wk,set k= k+1 and update velocity and position by Eq. 11-12|
|•||Step 5. Update pBest and gBest: After power flow and fitness calculation, update pBest and gBest. Go to step 3|
To verify the effectiveness and efficiency of the adopted TPSO for DOPF problems including FSWGs, IEEE 39-bus power system is used as the test systems. The procedure has been implemented in MATLAB 7.0 programming language and numerical tests are carried on a Pentium 4 2.4 G computer. The wind farm including 60 wind generators with the same type, the rating power of which reaches 36 MW, is connected to the system at the bus 14. For the single FSWG, its rating voltage is 690 V, the excitation reactance xm is 2.2059O, the sum of the stator reactance x1 and rotor reactance x2 are 0.1998O. For simplifying the analysis, the load size is considered invariable in the planning horizon. The planning horizon is divided into 12 intervals and every interval is 0.5 h. The wind generators' outputs are shown in Table 2. IEEE 39-bus system data are given by Zimmierman et al. (2009) and Victoire and Jeyakumar (2005). The parameters of the conventional generating units are shown in Table 3 and 4.
The maximum capacity of the reactive power compensation is 50 Mvar which constitutes of 50 set of group switched capacitors with per unit 1 Mvar. So, the upper limit value of the absorbed reactive power of wind farm bus is set to1Mvar in the study.
To demonstrate the superiority of the proposed approach for DOPF problems, simulation results have been compared with PSO method. Owing to the randomness in intelligent algorithms, two algorithms are executed 20 times when applied to the test system.
|Table 2:||The wind farm data in different periods|
|Table 3:||The parameters of conventional generating units|
|Table 4:||The parameters and limits of conventional generating units|
|Table 5:||Best solution obtained using TPSO method|
|Total production cost : 457039.8613 $/h|
|Table 6:||Best solution obtained using PSO method|
|Total production cost: 458852.7267 $/h|
For DOPF problem including FSWGs, Table 5 and 6 list the best control variables found by TPSO and PSO algorithm for IEEE 39-bus system, respectively. It is clearly shown that, by using TPSO, the total production cost savings of 1812.8654$/h is obtained compared with PSO algorithm. Hence, it is justified that TPSO approach gives the exact minimum dispatch solution. From Table 7, the best, worst and average cost values are 457039.8613 $/h, 458177.5291 $/h, 457518.5863 $/h and 458852.7267 $/h, 460137.6551 $/h, 459394.2806 $/h, respectively with TPSO and PSO after 20 independent trials. From the results, the superiority of TPSO strategies over PSO can be noticed. The difference between the best and worst solutions are 1137.6678 $/h with TPSO. At the same time, the difference between the best and worst solutions is 1284.9284$/h with PSO. Moreover, the best and worst solutions obtained by TPSO are very close to the average value, which proves that TPSO is more robust and consistent. In conclusion, it is clearly shown that TPSO is the most accurate and gives the exact minimum dispatch solution.
|Table 7:||Comparison of best, worst and average cost values in the IEEE 39-bus system|
|Table 8:||Average execution time comparison in the IEEE 39-bus system|
The average execution time taken to complete the fixed number of iterations (Tfix) and the average execution time taken to converge into the lower solution range (Tlow) for 20 trials are shown in Table 8 for IEEE 39-bus system.
For IEEE 39-bus system, PSO takes an average execution time of 4235.83 sec to complete 150 iterations. PSO converges faster than TPSO by reason of the small sub-memeplex generation number of TPSO. In comparison to PSO, TPSO has additional components, i.e., the Taguchi method. This extra burdens increase the execution time of TPSO. TPSO takes 4317.79 sec more than PSO to complete 150 iterations. Nevertheless, TPSO takes only 2140.45 sec to converge into the lower solution range (457039-457107 $/h), PSO are not able to converge into the lower solution range.
Considering the valve-point effect and ramp rate limits of conventional generators, as well as the discreteness of the reactive compensation devices including FSWGs, DOPF model, which takes the all conventional units cost minimum as the objective function and takes the whole time and the inherent relations of different stages into account in wind power integrated system, is established. The P-bus model of FSWG bus is adopted in power flow calculation in this study. A novel TPSO is proposed for solving the established DOPF model and the detailed methods of the algorithm are given. For DOPF problem including FSWGs, according to the principle of the nearest reactive power compensation, the required reactive power of wind farm doesnt absorb from the system as much as possible, but mainly from the wind farm reactive compensation devices. IEEE 39-bus system is used to illustrate the effectiveness of the proposed method compared with those obtained from PSO algorithm. The test results show that TPSO gives the exacter minimum dispatch solution and achieves better economic efficiency compared with PSO algorithm.
The authors would like to thank anonymous reviewers. This work is supported by the Science Research Program of Education Bureau of Hubei Province under Grant No. D20092906.
- He, D., F. Wang and Z. Mao, 2008. A hybrid genetic algorithm approach based on differential evolution for economic dispatch with valve-point effect. Int. J. Electr. Power Energy Syst., 30: 31-38.
- Li, H. and Z. Chen, 2008. Overview of different wind generator systems and their comparisons. IET Renew. Power Gener., 2: 123-138.
- Pothiya, S., I. Ngamroo and W. Kongprawechnon, 2008. Application of multiple tabu search algorithm to solve dynamic economic dispatch considering generator constraints. Energy Convers. Manage., 49: 506-516.
- Uturbey, W. and A.S. Costa, 2007. Dynamic optimal power flow approach to account for consumer response in short term hydrothermal coordination studies. IET Gener. Transm. Distrib., 1: 414-421.
- Victoire, T.A.A. and A.E. Jeyakumar, 2005. A modified hybrid EP-SQP approach for dynamic dispatch with valve-point effect. Int. J. Electr. Power Energy Syst., 27: 594-601.
- Yuan, X., L. Wang, Y.C. Zhang and Y.B. Yuan, 2009. A hybrid differential evolution method for dynamic economic dispatch with valve-point effects. Expert Syst. Applied, 36: 4042-4048.