Model-predictive control of power supply for particle accelerators

Model-predictive control of power supply for particle accelerators 1 corresp first-author Corresponding author, qianxp@lzu.edu.cn Corresponding author, qianxp@lzu.edu.cn qianxp@lzu.edu.cn 1 1 School of Nuclear Science and Technology, Lanzhou University, Lanzhou730000, China School of Nuclear Science and Technology, Lanzhou University, Lanzhou730000, China In this paper, model-predictive control (MPC) is proposed for controlling power source of accelerators. The system state equation is employed as the predictive model. With MPC, the difference between possible output and the ideal output is forecasted and decreased, so that the system can trace the ideal trail as closely and quickly as possible. The results of simulations and experiments show that this method can reduce influence of low frequency noise. Model predictive control Current power source State equation I. INTRODUCTION 1 Particle accelerators are widely used in physics research, clinic treatments and environmental protection. Performance of the power source is important to accelerator operation as it determines how well the magnetic field accelerates the particles, hence the high performance requirement on DC power supply of accelerators [ 1 1 , 2 2 ]. If the stability requirement is just d-3, a normal negative feedback control will be good enough. When a power stability of d-5 ∼ d-6 is required, however, the changes in amplification, power grid noise, and many other factors, begin to affect the power stability. Against such a requirement on power stability, the PID (proportion, integration and differentiation) method is often used, because it is easy and effective. However, this method cannot satisfy requirements all at once on accuracy, quickness and controlling simplification [ 3 3 - 7 7 3-7 ]. Model-predictive control (MPC) is one of the modern control theories developed based on the rapid progresses of computer technology. It grew quickly in 1970’s and has been widely used in industries. One of its advantages is that MPC does not need high precision model demand for the controlled system. By using the rolling optimization, the error of actual output from the ideal output is minimized. The system’s next prospective output will be justified in advance by the feedback of the error. Therefore, influence of the model uncertainty, undesirable noise and digital time delay can be reduced, and the system robustness is improved [ 8 8 - 13 13 8-13 ]. In this paper, together with the noise, the given input is added to the power source after it is modified with MPC. The control system is a close-loop structure. The results of a Matlab simulation show that influence of the noise at low frequency can be reduced effectively. A real model is set up and FPGA in VHDL language is devised so that the control method can actualize on it. Test results show that the method is effective. II. EQUIVALENT MODEL OF CURRENT POWER SOURCE 1 Generally, a current power source can be simplified as in Fig. 1 Fig. 1 . L 3 , R 3 , i L 3 and V O are the inductor, resistor, and output current and voltage, respectively, of the magnetic field coil; C 1 , C 2 , L 1 , R 1 and R 2 form the filter section. When T 1 and T 4 is conducting, V 1 = E ; and when T 2 and T 3 is conducting, V 1 = - E ; so the following equivalent equations can be derived: <math display="block" id="M1"><mrow><msub><mi>V</mi><mn>1</mn></msub><mo>=</mo><msub><mi>R</mi><mn>1</mn></msub><msub><mi>i</mi><mrow><mi>L</mi><mn>1</mn></mrow></msub><mo>+</mo><msub><mi>L</mi><mn>1</mn></msub><mtext>d</mtext><msub><mi>i</mi><mrow><mi>L</mi><mn>1</mn></mrow></msub><mo>/</mo><mtext>d</mtext><mi>t</mi><mn>,</mn></mrow></math> <math display="block" id="M2"><mrow><msub><mi>i</mi><mrow><mi>L</mi><mn>1</mn></mrow></msub><mo>=</mo><msub><mi>C</mi><mn>1</mn></msub><mtext>d</mtext><msub><mi>V</mi><mtext>O</mtext></msub><mo>/</mo><mtext>d</mtext><mi>t</mi><mo>+</mo><msub><mi>i</mi><mrow><mi>L</mi><mn>3</mn></mrow></msub><mo>+</mo><msub><mi>C</mi><mn>2</mn></msub><mtext>d</mtext><msub><mi>V</mi><mrow><msub><mi>C</mi><mn>2</mn></msub></mrow></msub><mo>/</mo><mtext>d</mtext><mi>t</mi><mn>,</mn></mrow></math> <math display="block" id="M3"><mrow><msub><mi>V</mi><mtext>O</mtext></msub><mo>=</mo><msub><mi>i</mi><mn>3</mn></msub><msub><mi>R</mi><mn>3</mn></msub><mo>+</mo><msub><mi>L</mi><mn>3</mn></msub><mtext>d</mtext><msub><mi>i</mi><mn>3</mn></msub><mo>/</mo><mtext>d</mtext><mi>t</mi><mn>,</mn></mrow></math> <math display="block" id="M4"><mrow><msub><mi>V</mi><mtext>O</mtext></msub><mo>=</mo><msub><mi>V</mi><mrow><msub><mi>C</mi><mn>2</mn></msub></mrow></msub><mo>+</mo><msub><mi>R</mi><mn>2</mn></msub><msub><mi>C</mi><mn>2</mn></msub><mtext>d</mtext><msub><mi>V</mi><mrow><msub><mi>C</mi><mn>2</mn></msub></mrow></msub><mo>/</mo><mtext>d</mtext><mi>t</mi><mn>.</mn></mrow></math> From Eqs. (1)– (4), the state equation can be obtained: <math display="block" id="M5"><mrow><mrow><mo>[</mo><mrow><mtable><mtr><mtd><mrow><msub><mover accent="true"><mi>i</mi><mo>˙</mo></mover><mrow><msub><mi>L</mi><mn>1</mn></msub></mrow></msub></mrow></mtd></mtr><mtr><mtd><mrow><msub><mover accent="true"><mi>i</mi><mo>˙</mo></mover><mrow><msub><mi>L</mi><mn>3</mn></msub></mrow></msub></mrow></mtd></mtr><mtr><mtd><mrow><msub><mover accent="true"><mi>V</mi><mo>˙</mo></mover><mrow><msub><mi>C</mi><mn>2</mn></msub></mrow></msub></mrow></mtd></mtr><mtr><mtd><mrow><msub><mover accent="true"><mi>V</mi><mo>˙</mo></mover><mtext>O</mtext></msub></mrow></mtd></mtr><mtr><mtd><mrow/></mtd></mtr></mtable></mrow><mo>]</mo></mrow><mo>=</mo><mrow><mo>[</mo><mrow><mtable><mtr><mtd><mrow><mfrac><mrow><msub><mi>R</mi><mn>1</mn></msub></mrow><mrow><msub><mi>L</mi><mn>1</mn></msub></mrow></mfrac></mrow></mtd><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd><mtd><mrow><mfrac><mn>1</mn><mrow><msub><mi>L</mi><mn>1</mn></msub></mrow></mfrac></mrow></mtd></mtr><mtr><mtd><mn>0</mn></mtd><mtd><mrow><mfrac><mrow><msub><mi>R</mi><mn>3</mn></msub></mrow><mrow><msub><mi>L</mi><mn>3</mn></msub></mrow></mfrac></mrow></mtd><mtd><mn>0</mn></mtd><mtd><mrow><mfrac><mn>1</mn><mrow><msub><mi>L</mi><mn>3</mn></msub></mrow></mfrac></mrow></mtd></mtr><mtr><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd><mtd><mrow><mfrac><mn>1</mn><mrow><msub><mi>R</mi><mn>2</mn></msub><msub><mi>C</mi><mn>2</mn></msub></mrow></mfrac></mrow></mtd><mtd><mrow><mfrac><mn>1</mn><mrow><msub><mi>L</mi><mn>3</mn></msub></mrow></mfrac></mrow></mtd></mtr><mtr><mtd><mrow><mfrac><mn>1</mn><mrow><msub><mi>C</mi><mn>1</mn></msub></mrow></mfrac></mrow></mtd><mtd><mrow><mfrac><mn>1</mn><mrow><msub><mi>C</mi><mn>1</mn></msub></mrow></mfrac></mrow></mtd><mtd><mrow><mfrac><mn>1</mn><mrow><msub><mi>R</mi><mn>2</mn></msub><msub><mi>C</mi><mn>1</mn></msub></mrow></mfrac></mrow></mtd><mtd><mrow><mfrac><mn>1</mn><mrow><msub><mi>R</mi><mn>2</mn></msub><msub><mi>C</mi><mn>1</mn></msub></mrow></mfrac></mrow></mtd></mtr><mtr><mtd><mrow/></mtd><mtd><mrow/></mtd><mtd><mrow/></mtd><mtd><mrow/></mtd></mtr></mtable></mrow><mo>]</mo></mrow><mrow><mo>[</mo><mrow><mtable><mtr><mtd><mrow><msub><mi>i</mi><mrow><msub><mi>L</mi><mn>1</mn></msub></mrow></msub></mrow></mtd></mtr><mtr><mtd><mrow><msub><mi>i</mi><mrow><msub><mi>L</mi><mn>3</mn></msub></mrow></msub></mrow></mtd></mtr><mtr><mtd><mrow><msub><mi>V</mi><mrow><msub><mi>C</mi><mn>2</mn></msub></mrow></msub></mrow></mtd></mtr><mtr><mtd><mrow><msub><mi>V</mi><mtext>O</mtext></msub></mrow></mtd></mtr><mtr><mtd><mrow/></mtd></mtr></mtable></mrow><mo>]</mo></mrow><mo>+</mo><mrow><mo>[</mo><mrow><mtable><mtr><mtd><mrow><mfrac><mn>1</mn><mrow><msub><mi>L</mi><mn>1</mn></msub></mrow></mfrac></mrow></mtd></mtr><mtr><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mn>0</mn></mtd></mtr><mtr><mtd><mrow/></mtd></mtr></mtable></mrow><mo>]</mo></mrow><msub><mi>V</mi><mn>1</mn></msub><mo>,</mo></mrow></math> <math display="block" id="M6"><mrow><mtable columnalign="left"><mtr columnalign="left"><mtd columnalign="left"><mrow><mi>Y</mi><mo>=</mo><mo stretchy="false">(</mo><mtable><mtr><mtd><mn>0</mn></mtd><mtd><mn>1</mn></mtd><mtd><mn>0</mn></mtd><mtd><mn>0</mn></mtd></mtr></mtable><mo stretchy="false">)</mo><mrow><mo>[</mo><mrow><mtable><mtr><mtd><mrow><msub><mi>i</mi><mrow><msub><mi>L</mi><mn>1</mn></msub></mrow></msub></mrow></mtd></mtr><mtr><mtd><mrow><msub><mi>i</mi><mrow><msub><mi>L</mi><mn>3</mn></msub></mrow></msub></mrow></mtd></mtr><mtr><mtd><mrow><msub><mi>V</mi><mrow><msub><mi>C</mi><mn>2</mn></msub></mrow></msub></mrow></mtd></mtr><mtr><mtd><mrow><msub><mi>V</mi><mtext>O</mtext></msub></mrow></mtd></mtr><mtr><mtd><mrow/></mtd></mtr></mtable></mrow><mo>]</mo></mrow><mo>+</mo><mo stretchy="false">(</mo><mn>0</mn><mo stretchy="false">)</mo><msub><mi>V</mi><mn>1</mn></msub><mn>,</mn></mrow></mtd></mtr></mtable></mrow></math> where, V 2 = E ( t 2 - t 1 )/( t 2 + t 1 ), t 1 and t 2 are conducting time of T 1 and T 4 , and T 2 and T 3 , respectively, within one time period. III. ESTABLISHMENT OF PREDICTIVE CONTROL MODEL 1 The overall MPC control structure is shown in Fig. 2 Fig. 2 . A predictive model is designed to forecast a prospective output y p ( n +1). As the predictive model is not the real model, so the output y p ( n +1) is modified to get the final prospective output y pr ( n +1), from the real model output y ( n ). Then an optimization is performed to minimize a cost function of the tracking error and the finale input u ( n ) is obtained. At every new sampling period, the measured real output modifies the forecasted output on the basis of the predictive model. Then, a new round of optimization follows, and this forms a close-loop control [ 13 13 - 20 20 13-20 ]. A. The predictive control model 2 The predictive control model forecasts the subsequent output of the system on basis of the actual output just measured and the future time inputted. Any kind of model having some features of the system can be used as predictive model. In this paper, the predictive model is state equation of the power source. The discrete state equations can be described as: <math display="block" id="M7"><mrow><mi>x</mi><mo stretchy="false">(</mo><mi>n</mi><mo>+</mo><mn>1</mn><mo stretchy="false">)</mo><mo>=</mo><mi>G</mi><mi>x</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mo>+</mo><mi>H</mi><mi>u</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mn>,</mn></mrow></math> <math display="block" id="M8"><mrow><mi>y</mi><mo stretchy="false">(</mo><mi>n</mi><mo>+</mo><mn>1</mn><mo stretchy="false">)</mo><mo>=</mo><mi>C</mi><mi>x</mi><mo stretchy="false">(</mo><mi>n</mi><mo>+</mo><mn>1</mn><mo stretchy="false">)</mo><mn>,</mn></mrow></math> where G = I + AT , and H = BT . I is the identity matrix, T is the sample period. A , B , C is state matrix, input matrix and output matrix of the system’s state equation, respectively. Then, the predictive output in time n +1 can be obtained: <math display="block" id="M9"><mrow><msub><mi>y</mi><mi>p</mi></msub><mo stretchy="false">(</mo><mi>n</mi><mo>+</mo><mn>1</mn><mo stretchy="false">)</mo><mo>=</mo><mi>C</mi><mi>G</mi><mi>x</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mo>+</mo><mi>C</mi><mi>H</mi><mi>u</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mn>,</mn></mrow></math> where x ( n ) and u ( n ) are state variables and input at n time. B. Feedback adjustment and rolling optimization 2 In the MPC, feedback adjustment and rolling optimization are required. Since the predictive model is not the real model, so its output y pr ( n +1) need to be justified by the error e ( n ) of the actual output from the last predictive output <math display="block" id="M10"><mrow><msub><mi>y</mi><mrow><mtext>pr</mtext></mrow></msub><mo stretchy="false">(</mo><mi>n</mi><mo>+</mo><mn>1</mn><mo stretchy="false">)</mo><mo>=</mo><msub><mi>y</mi><mtext>p</mtext></msub><mo stretchy="false">(</mo><mi>n</mi><mo>+</mo><mn>1</mn><mo stretchy="false">)</mo><mo>+</mo><mi>h</mi><mi>e</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mn>,</mn></mrow></math> <math display="block" id="M11"><mrow><mi>e</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mo>=</mo><mi>y</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mo>−</mo><msub><mi>y</mi><mtext>p</mtext></msub><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mn>,</mn></mrow></math> where y ( n ) and y p ( n ) are the actual output and predictive output of n time, respectively; and h is error correction coefficient ( h =1). This feedback adjustment makes the predictive model more assimilate to the real model. C. The rolling optimization 2 The purpose of predictive control is to make the system output trace the ideal output as quickly and closely as possible. So the control objective is to minimize the sum of square errors of the prospective value from the given value. The input, however, should not be too large. Thus, the quadratic performance expression is: <math display="block" id="M12"><mrow><mi>J</mi><mo>=</mo><mi>p</mi><msup><mrow><mo stretchy="false">[</mo><msub><mi>y</mi><mrow><mtext>pr</mtext></mrow></msub><mo stretchy="false">(</mo><mi>n</mi><mo>+</mo><mn>1</mn><mo stretchy="false">)</mo><mo>−</mo><msub><mi>y</mi><mi>l</mi></msub><mo stretchy="false">(</mo><mi>n</mi><mo>+</mo><mn>1</mn><mo stretchy="false">)</mo><mo stretchy="false">]</mo></mrow><mn>2</mn></msup><mo>+</mo><mi>q</mi><mi>u</mi><msup><mrow><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo></mrow><mn>2</mn></msup><mn>,</mn></mrow></math> where p and q are weighting coefficients. It shows that as the process goes on, the optimum varies on line. So this optimization relates only to dynamic performance of the system. To minimize the objective function, the derivative of J should be equal to zero and then the control value is: <math display="block" id="M13"><mrow><mi>u</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mo>=</mo><mi>z</mi><mi>p</mi><mi>C</mi><mi>H</mi><mo>/</mo><mo stretchy="false">[</mo><mi>p</mi><msup><mrow><mo stretchy="false">(</mo><mi>C</mi><mi>H</mi><mo stretchy="false">)</mo></mrow><mn>2</mn></msup><mo>+</mo><mi>q</mi><mo stretchy="false">]</mo><mn>,</mn></mrow></math> <math display="block" id="M14"><mrow><mi>z</mi><mo>=</mo><msub><mi>y</mi><mi>l</mi></msub><mo stretchy="false">(</mo><mi>n</mi><mo>+</mo><mn>1</mn><mo stretchy="false">)</mo><mo>−</mo><mi>C</mi><mi>G</mi><mi>x</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mo>−</mo><mi>e</mi><mo stretchy="false">(</mo><mi>n</mi><mo stretchy="false">)</mo><mn>.</mn></mrow></math> IV. SIMULATION AND EXPERIMENTAL RESULTS 1 Based on the model in Fig. 1 Fig. 1 , a simulation was performed with Matlab. Specifications of that model are: L 1 = 0.3 mH, R 1 = 0.01 Ω, L 3 = 91.4 mH, R 3 = 0.0796 Ω, C 1 = 10 μF, C 2 = 47 μF, and R 2 = 1 Ω. The sampling frequency is 10 kHz, and the total simulating time is 10 s. As shown in Fig. 2 Fig. 2 , noise ( n ) of different frequencies was added to the power source to check their influence on the controlled object. The results are given in Figs. 3 Figs. 3 and 4 4 Six simulations were performed with two kinds of input: ramp input ( Fig. 3 Fig. 3 ) and sine input ( Fig. 4 Fig. 4 ). The noise added to the system was in sine wave. In order to find how different noises affect the output, two kinds of frequency were included: high and low. The specific input and noise for each simulation are given in Table 1 Table 1 , where n low stands for low frequency noise and n high , high. 1 1 Figs. 1 1 Inputs 1 1 n low 1 1 n high 1 1 3(a) 1 1 u = t 1 1 1 * (2 π * 2 * t ) 1 1 0 1 1 3(b) 1 1 u = t 1 1 1 * (2 π * 2 * t ) 1 1 0.1 * sin(2 π * 30 * t ) 1 1 3(c) 1 1 u = t 1 1 1 * (2 π * 2 * t ) 1 1 0.1 * sin(2 π * 100 * t ) 1 1 4(a) 1 1 u =10 * ( t ) 1 1 3 * (2 π * 2 * t ) 1 1 0 1 1 4(b) 1 1 u =10 * ( t ) 1 1 3 * (2 π * 2 * t ) 1 1 0.1 * sin(2 π * 30 * t ) 1 1 4(c) 1 1 u =10 * ( t ) 1 1 3 * (2 π * 2 * t ) 1 1 0.1 * sin(2 π * 100 * t ) In Figs. 3 Figs. 3 and 4 4 , the red curves are system output controlled by MPC, and the black curves are without the MPC. The X -axis is the time and the Y -axis is the system’s output of corresponding input. In Fig. 3(a) Fig. 3(a) , the input is in ramp wave. Compared to the red curve, the black one reflects the influence of the noise on the output. In fact the red curve is also the ideal output that has no disturbance on the system. So it is clear the MPC diminished the noise influence. In Fig. 3(b) Fig. 3(b) , in addition to the same input and noise as in Fig. 3(a) Fig. 3(a) , a 30 Hz noise was added to the system. We can see the red one is much more close to the ideal curve than the black one. This is because that the influence of low frequency noise was greatly reduced by MPC. However, comparing to the ideal output in Fig. 3(a) Fig. 3(a) , there are ripple waves in the red curve, caused by the 30 Hz noise. In Fig. 3(c) Fig. 3(c) , the added is noise which is even higher (100 Hz), and bigger ripple waves are seen. So, the MPC decreases influence of low frequency noise greatly but it can do little against higher frequency noises. In Fig. 4 Fig. 4 , the input is sine wave. The red curve in Fig. 4(a) Fig. 4(a) is the ideal output without the disturbance in the input. From Figs. 4(b) Figs. 4(b) and Fig. 4(c) Fig. 4(c) the red curves are much more similar to the ideal one than the black ones, but there are ripples, too, at noises of higher frequencies. The MPC works against low frequency noise of sine wave. A real model with specifications shown in Fig. 5 Fig. 5 was made. Fig. 5 Fig. 5 differs from Fig. 1 Fig. 1 in that there is resistance for the capacity and inductance in the real model. Hardware implementation that actualized the proposed method was achieved by using GX-SOC/SOPC-CIDE simulator box. The system input was sine waveform with some regular noise of different frequencies. Fig. 6 Fig. 6 shows the experimental output waveform without and with the MPC. One sees that the curve in Fig. 6(b) Fig. 6(b) is more similar to the standard sine wave than that is in Fig. 6(a) Fig. 6(a) , because of the noise reduction by the MPC. Therefore, the proposed control strategy is validated. V. CONCLUSION 1 A real object to be controlled is much more complicated than the theoretical model. Because of unexpected disturbance, inaccuracy of components and computational delay, the real output differs with the ideal output. This problem must be taken into account in design of the controller algorithm to improve its performance and robustness. In this paper, we present an MPC method to reduce errors of the real output from the ideal output, and obtain a quick and accuracy respond for the power source that used in the particle accelerators. Unlike conventional control theories, MPC can forecast the output according to information of both the future reference input and the last output. By optimally determining the future input on every sampling period, the proposed method can adjust the input in advance to diminish influence of the unexpected factors as much as possible, so that the real output can trace the ideal output closely. The control strategy employs the state equation as the predictive model, and noises of different frequency are added to the system to search the impact of the MPC on the system. The proposed control scheme is simple and computationally efficient, requiring just a few operations. The results of simulations and experiments show improvement on reduction of the noise influence. References Chen Y X , Yan H B , Huang Y Z . High Pow Laser Part Beams , 2010 , 9 : 2138 - 2142 . (in Chinese) High Pow Laser Part Beams Xia J W , Zhan W L , Wei B W . Atom Energ Sci Technol , 2009 , 43 : 150 - 158 . (in Chinese) Atom Energ Sci Technol Shi C F , Huang Y Z , Chen Y X . Atom Energ Sci Technol , 2013 , 47 : 1652 - 1655 . (in Chinese) Atom Energ Sci Technol Wang B Y and Ma X L . Pow Electron , 2010 , 2 : 67 - 69 . (in Chinese) Pow Electron Wang R K , Huang Y Z , Chen Y X . Nucl Electron Detec Tech , 2013 , 33 : 558 - 561 . (in Chinese) Nucl Electron Detec Tech Huang Y Z , Chen Y X , Zhou Z Z . Nucl Phys Rev , 2011 , 28 : 296 - 299 . (in Chinese) Nucl Phys Rev Wang B Y . Ph. D. Thesis , Tianjin University , 2006 . (in Chinese) Ph. D. Thesis Wang Z P , Xie Y X , Wang Y . Auto Inf Eng , 2013 , 34 : 34 - 40 . (in Chinese) Auto Inf Eng Liu Y , Wang H Q , Li P . Cont Theory App , 2009 , 26 : 107 - 110 . (in Chinese) Cont Theory App Morel F , Lin-Shi X , Rétif J M , et al . Electric Pow Syst Res , 2008 , 78 : 1437 - 1447 . Electric Pow Syst Res Grimble C J and Ordys A W . Annu Rev Control , 2001 , 25 : 13 - 24 . Annu Rev Control Bleris L G , Garcia J , Kothare M V , et al . J Process Contr , 2006 , 16 : 255 - 264 . J Process Contr Wang X G , Xie Y X , Shuai D X . J South China Univ Tech , 2009 , 11 : 56 - 63 . (in Chinese) J South China Univ Tech Xu Z L and Luo C F . Pow Syst Protect contl , 2009 , 37 : 83 - 86 . (in Chinese) Pow Syst Protect contl Errouissi R and Ouhrouche M . Math Comput Simulat , 2010 , 81 : 394 - 406 . Math Comput Simulat Bleris L G , Garcia J , Kothare M V . J Process Contr , 2006 , 16 : 255 - 264 . J Process Contr Kasserm A M . Electrical Power and Energy Systems , 2012 , 41 : 124 - 132 . Electrical Power and Energy Systems Blanco E , Prada C , Cristea S , et al . Control Eng Pract , 2009 , 17 : 1136 - 1147 . Control Eng Pract Wang S B , Li C W , Gong W . Contr Instrum Chem Indust , 2010 , 31 : 85 - 88 . (in Chinese) Contr Instrum Chem Indust Bouafia A , Gaubert J , Krim F . Energ Convers Manage , 2010 , 51 : 2473 - 2481 . Energ Convers Manage 10.13538/j.1001-8042/nst.25.050203 Supported by National Natural Science Foundation of China (No. 11027508) 2013-11-15 2014-06-27 2014-09-25 2014-10-20 2021-04-11 2014 Nuclear Science and Techniques 5 25