Power System Voltage Correction Scheme Based on Adaptive Horizon Model Predictive Control

Model predictive control (MPC) is commonly used to compensate for modeling inaccuracies and measurement noise in voltage control problems. The length of the prediction horizon and control horizon of a MPC-based method has significant impact on the control performances. In existing relevant works, those horizon parameters are determined off-line based on experience or enumeration, and keeps constant during the entire receding-horizon optimization process. This paper presents a system voltage correction scheme based on adaptive horizon model predictive control (AH-MPC). The reactive power compensation and voltage regulation devices are coordinated to maintain the system voltages within a desired range. An evaluation index is proposed to determine the horizon parameters, which reflects the maximum voltage regulation ability with the current parameter configuration. Within each sampling interval, the horizon parameters are updated according to the evaluation index and real-time measurements periodically, which comprehensively considers the system uncertainties and voltage recovery speed, and the computational effort is remarkably reduced. The validation and effectiveness of the proposed method is verified by the simulation analysis on the test system.


Introduction
The power system voltages should be maintained within predefined limits according to the operation requirements. Nowadays, the increasing penetration of distributed renewable energy has had significant impacts on system behavior due to its variability and uncertainty. Besides, the power systems are operated more closely to their physical limits. The risk of temporary voltage problems or even voltage collapse has increased [1]. When abnormal or unstable voltage profiles are detected, appropriate countermeasures should be taken to enhance voltage quality, power system security and operating economy.
A reliable real-time optimal corrective voltage control scheme requires valid power system prediction models and feedback signals. Due to the inevitable modeling inaccuracies and measurement noise, it is impossible to build an accurate mathematical model to describe all the dynamic features of power systems, thus non-optimal or even inappropriate control decisions may be made. Model predictive control (MPC) is an appealing approach for dealing with the above issues [2], which is especially suitable for complicated systems with changeable dynamics and uncertain factors. The prediction model of a MPC-based method traces the operating state of practical system, and the control sequences are updated and corrected by real-time measurements.
In order to reduce prediction complexity for real-time application, linear prediction models are commonly applied in a MPC-based voltage control scheme. In [3,4], the coordinated voltage control is realized using tree search optimization method based on MPC. The authors of [5,6] apply trajectory sensitivity and MPC to solve the optimal voltage control problem. In [7], the prediction model of a MPC-based real time voltage control scheme is simplified using wide area measurement information in receding horizons. In [8,9], the search solution space of the MPC-based voltage controller is significantly reduced with the prepared control knowledge. A distributed model predictive control scheme for interconnected power systems is presented in [10]. In [11], MPC is applied in a large power system with most effective controls identified.
In recent years, the increasing penetration of renewable power can introduce high variability into transmission and distribution systems, which can be well handled by MPC-based methods. A hierarchical controller based on MPC is proposed to regulate voltages in system with high penetration of wind power in [12]. Model predictive voltage control methods of wind farms are introduced in [13][14][15]. Reference [16] proposes a centralized controller based on MPC to mitigate voltage violations in active distribution networks (ADN) by coordinating power outputs of distributed generation units and the voltage set-point of the transformer load tap changer. The impact of the model predictive control of the ADN voltages to the transmission grid is discussed in [17]. A two-level real-time voltage correction scheme is proposed in [18]. The local level provides the fast response to a disturbance, and the central level coordinates the various distributed generation units relying on the concept of MPC. In [19], the performance of a MPC-based voltage controller is verified by the field test in a real distribution system.
The length of the control horizon and prediction horizon of the existing MPC-based voltage control method has significant impact on the voltage control performances. A larger value of control horizon leads to smoother control but increased computation time. A smaller value of control horizon leads to less control steps and more aggressive control decisions, but it weakens the ability to deal with the modeling inaccuracies and measurement noise problem. The length of the prediction horizon should be equal to the length of the control horizon unless the controller is requested to consider changes happening beyond the control horizon [16]. References [3,4,[6][7][8][9][10][11][12][13][14][15]18,19] apply constant horizon parameters, but how those parameters are determined is not discussed. The horizon parameters in [5] is chosen to be the time in which the predicted voltage drops by a predefined percentage at the initial control instant, and the horizon parameters are kept constant during the entire receding-horizon optimization process. In [16,17], the controller uses constant horizon parameters unless the prediction horizon is increased to anticipate load tap changer actions. The control methods discussed so far determines the control horizon parameters off-line based on experience or enumeration without quantitative analysis.
In this paper, a coordinated voltage correction scheme based on adaptive horizon model predictive control (AH-MPC) is proposed to maintain a satisfactory voltage profile while minimizing the control cost. Generator voltage references, shunt capacitors and on-load tap changers are considered as candidate control actions. The proposed scheme is based on the principles of MPC. However, it differs from the standard MPC-based method by way the horizon parameters are chosen and embedded in the optimization procedure. The main contribution of this paper is developing a principle of determine the horizon parameters, which is a key feature missing in previous control schemes. An index is formulated to evaluate and ensure the applicability of the horizon parameter configurations, and barely consume extra computation time. Within each sampling interval, the length of the horizon parameters are updated according to the proposed evaluation index and real-time measurements, aiming at better control performances compared with the traditional MPC-based methods.
This paper is organized as follows: Section 2 introduces the standard MPC-based voltage control model. Section 3 presents the determination of the adaptive horizon parameters. The voltage correction scheme based on AH-MPC is presented in Section 4. The simulation results and conclusions are presented in Sections 5 and 6 respectively.

Voltage Prediction
MPC finds a sequence of control actions applied within the control horizon in order to achieve the control objectives within the prediction horizon considering equality and inequality constraints. At each sampling point, only the first step of the so-computed sequence is applied. The principle of MPC is shown in Figure 1. The length of the control horizon and prediction horizon are t c and t p respectively. t s is the sampling interval. N c and N p are the number of sampling intervals within the control horizon and prediction horizon respectively, t c = N c t s , t p = N p t s . correction scheme based on AH-MPC is presented in Section 4. The simulation results and conclusions are presented in Section 5 and Section 6 respectively.

Voltage Prediction
MPC finds a sequence of control actions applied within the control horizon in order to achieve the control objectives within the prediction horizon considering equality and inequality constraints. At each sampling point, only the first step of the so-computed sequence is applied. The principle of MPC is shown in Figure 1. The length of the control horizon and prediction horizon are tc and tp respectively. ts is the sampling interval. Nc and Np are the number of sampling intervals within the control horizon and prediction horizon respectively, tc = Ncts, tp = Npts.  Real-time optimization requires reliable and quick decisions. A linear model based on trajectory analysis is applied to predict the trajectory of the system with different control sequences. The discrepancies arising from the model mismatch problems and system uncertainties can be resolved by the closed-loop feature of the centralized MPC-based controller [2]. The trajectory sensitivity analysis can provide valuable insights evaluating the effects of slight control variations on the system evolution. The detailed theory and analytical basis of trajectory sensitivity is discussed in [20].
The system model can be expressed by the differential-algebraic Equations (DAE) (1) and (2) subject to some constraints, the detailed introduction of the DAE model can be found in [1].
where x is the vector of state variables, which is related to the generator and its regulation system, the load recovery behavior et al., y is the vector of algebraic variables such as active and reactive powers P and Q, voltage magnitudes V and phase angles θ, y = [P, Q, V, θ]. u is the vector of control variables. In this paper, the set points of automatic voltage regulators (AVR), the switching of shunt capacitors and the position of on-load tap changers (OLTCs) are considered as candidate control actions. The sensitivity Equations (3) and (4) are obtained by differentiating (1) and (2) with respect to u. Real-time optimization requires reliable and quick decisions. A linear model based on trajectory analysis is applied to predict the trajectory of the system with different control sequences. The discrepancies arising from the model mismatch problems and system uncertainties can be resolved by the closed-loop feature of the centralized MPC-based controller [2]. The trajectory sensitivity analysis can provide valuable insights evaluating the effects of slight control variations on the system evolution. The detailed theory and analytical basis of trajectory sensitivity is discussed in [20].
The system model can be expressed by the differential-algebraic Equations (DAE) (1) and (2) subject to some constraints, the detailed introduction of the DAE model can be found in [1].
where x is the vector of state variables, which is related to the generator and its regulation system, the load recovery behavior et al., y is the vector of algebraic variables such as active and reactive powers P and Q, voltage magnitudes V and phase angles θ, y = [P, Q, V, θ]. u is the vector of control variables. In this paper, the set points of automatic voltage regulators (AVR), the switching of shunt capacitors and the position of on-load tap changers (OLTCs) are considered as candidate control actions. The sensitivity Equations (3) and (4) are obtained by differentiating (1) and (2) with respect to u.
4 of 18 0 = g x (t)x u (t) + g y (t)y u (t) + g u (t) (4) where y u (t) is composed of the trajectory sensitivity of all the algebraic variables with respect to the candidate control variables at a time instant t, which can be obtained by solving (3) and (4).
. When the time domain simulation of the system (1) and (2) is based on the trapezoidal numerical integration, the Jacobian obtained during the numerical integration process is the same as the coefficient matrix of solving the sensitivity Equations (3) and (4). When the result of time domain simulation is obtained, the calculation of trajectory sensitivity involves little extra computation [20]. When control action ∆u is applied, the corresponding voltage variation at a sampling point k can be approximately expressed as ∆V k = V u (t n + kt s )∆u, V u (t n + kt s ) is the voltage trajectory sensitivity at sampling point k with respect to control u, and also one of the elements that constitute y u (t).
The voltage prediction based on trajectory sensitivity analysis is graphically illustrated in Figure 2. According to the proportion and superposition characteristics of linear systems, the predicted voltage at a sampling point k is composed of two parts: the voltage without controlV * k during the prediction horizon and the voltage variation ∆V k cause by the control actions applied at each sampling points. At the initial control time,V * k is solved based on trapezoidal numerical integration of (1) and (2). ∆V k is solved based on the Jacobian obtained during the numerical integration process. The vector consist of all the predicted voltages at each sampling point can be expressed as: where yu(t) is composed of the trajectory sensitivity of all the algebraic variables with respect to the candidate control variables at a time instant t, which can be obtained by solving (3)

and (4). yu(t) = [Pu(t), Qu(t), Vu(t), θu(t)].
When the time domain simulation of the system (1) and (2) is based on the trapezoidal numerical integration, the Jacobian obtained during the numerical integration process is the same as the coefficient matrix of solving the sensitivity Equations (3) and (4). When the result of time domain simulation is obtained, the calculation of trajectory sensitivity involves little extra computation [20]. When control action u  is applied, the corresponding voltage variation at a sampling point k can be approximately expressed as is the voltage trajectory sensitivity at sampling point k with respect to control u, and also one of the elements that constitute yu(t).
The voltage prediction based on trajectory sensitivity analysis is graphically illustrated in Figure 2. According to the proportion and superposition characteristics of linear systems, the predicted voltage at a sampling point k is composed of two parts: the voltage without control * ˆk V is solved based on the Jacobian obtained during the numerical integration process. The vector consist of all the predicted voltages at each sampling point can be expressed as: Voltage magnitude time t n +t s t n +kt s t n +2t s … … … …

Optimization Model
The proposed scheme aims at minimizing the system voltage deviation and the corresponding control cost. The optimization problem can be formulated as follows: Subject to:

Optimization Model
The proposed scheme aims at minimizing the system voltage deviation and the corresponding control cost. The optimization problem can be formulated as follows: Subject to: V min ≤V N p ≤ V max (8) u min ≤ u k ≤ u max (9) |∆u k | ≤ ∆u max (10) where V r is the vector consists of the control references.V k is the predicted voltage at sampling point k. ∆u k is the vector consists of the control variable changes at sampling point k. Q and R are both diagonal matrices which define the weight coefficients of the controlled voltages and the cost of the control actions respectively. The voltage constraints in the prediction horizon is shown in Equation (8),V N p represents the predicted voltage at the sampling point N p . The objective function sets a target of minimizing the voltage deviations at every sampling point under the constraints of control input (9) and (10), and the voltage constraints are only enforced at the end of the prediction horizon. The inequality constraints (9) stipulates the upper and lower bounds of the control inputs.
When undesirable bus voltage magnitudes are detected, appropriate voltage control actions should be applied as soon as possible to reduce the negative influence to the consumers and avoid the deterioration of system voltages. The single-step optimization approach or large control variations in one control step may suffers from deficiencies such as low optimization efficiency or undesirable control results [21]. Thus the constraint of control variation in one control step (10) is included in the model.

Adaptive Horizon Parameters
MPC is a receding-horizon multi-step optimization method. The computation time and control performance are closely related to the horizon parameter selection, which determines when and how the system variables will reach their target states. This section presents the determination of adaptive horizon parameters, which are periodically updated according to the real-time measurements and system trajectory prediction during the entire optimization process of MPC.

Principles for the Horizon Parameter Settings
At each sampling point, the MPC-based method obtains an optimal control sequence [u 1 , u 2 , . . . , u Nc ] distributed within the control horizon. The control variables of the control sequence keeps constant between sampling intervals, and only changes at sampling points within the control horizon. The basic principles for the horizon parameters are listed as follows:

•
Sampling interval t s . The length of sampling interval t s should fulfill the requirements of measurement collection, data processing and optimization, which change little during the entire receding horizon optimization process. A constant value t s = 10 s is applied in the proposed scheme, which also matches the value used in the coordinated secondary voltage control in operation in practical systems [22]. • Control horizon t c . With a constant sampling interval, the optimal length of the control horizon is related to the voltage evolution and response feature, which varies over the receding horizon optimization process. In this paper, the system voltage is predicted using a linear model as is introduced in Section 2. The sensitivity feature changes little if no disturbances occur during the optimization process [21], and the solved control variation is proportional to the detected voltage deviation. When severe voltage deviation is detected, the length of the control horizons should be short enough to meet the requirement of the voltage recovery speed and online computation time, and should be long enough to deal with modeling inaccuracies and measurement noise. As the receding horizon optimization progressed, the controlled voltages are gradually closer to its desired range, and smaller control variations are needed to accomplish the control objective.
In this situation, shorter control horizon parameters could reduce settling times of the optimization process and computational effort.

•
Prediction horizon t p . The system response caused by the control actions are predicted within the prediction horizon, thus the length of the prediction horizon should be greater or equal to the length of the control horizon. When the controlled voltages are gradually stabilized to the steady state values, the necessity to consider changes happening beyond the control horizon is reduced from a computational viewpoint. Therefore, at the beginning of the MPC implementation, the length of the prediction horizon is set longer than the length of control horizon. As the receding horizon optimization progressed, the length of the prediction horizon gradually approaches the length of the control horizon, which could also save unnecessary computation time

Voltage Regulation Ability Evaluation
Under the condition of fulfilling the control variation constraint (10), large control variations are forced to be applied in multi-steps. The limited candidate control actions in one control step have theoretical maximum voltage regulation ability for every bus. The procedure to estimate the ability is summarized as follows: (1) At t n , the voltage evolution trajectories are predicted within a predetermined prediction horizon [t n , t n + t pn ] based on the trapezoidal numerical simulation. t pn is the length of prediction horizon applied in the n-th optimization. If the predicted voltage deviation of load bus i satisfies the inequality constraint (11) at the initial time instant t 1 , the voltage magnitude of bus i is determined as the evaluation target.
is the voltage magnitude of bus i without control at the prediction time instant (t 1 + t p1 ) obtained by numerical simulation. DB is the control deadband. V max and V min initially presented in Equation (8) are subject to: The trajectory sensitivity matrix of the evaluation targets at time instant t n + t pn with respect to the control action of each step can be expressed as: where S k(i, j) represents the trajectory sensitivity of the i-th evaluation target with respect to the j-th candidate control in the k-th control step. N i and N j are the total number of evaluation targets and available control actions respectively. (3) Based on the voltage prediction and trajectory sensitivity analysis, the maximum control target variation caused by controls applied at sampling point k is evaluated as follows.
) ≥ 0, the optimization tends to increase the j-th control in order to raise (when S k(i,j) ≥ 0 and V r ≥V * (i) (t n + t pn )) or reduce (when S k(i,j) ≤ 0 and V r ≤V * (i) (t n + t pn )) the voltage magnitude of bus i towards its desired value. Under the inequality constraints (9) and (10), the predicted voltage magnitude variation of bus i reaches its maximum value when the variation of control j equals its maximum positive value: where ∆V k(i,j)max (t n + t pn ) is the predicted maximum variation of the evaluation target i at a future time instant t n + t pn , which is caused by the j-th control action applied at the sampling point k. ∆u (j)max is the j-th element in vector ∆u max , representing the maximum value of the corresponding control variation. u (j)max is the j-th element in vector u max .
Similarly, when S k(i,j) · (V r −V * (i) (t n + t pn )) < 0, the optimization tends to reduce the j-th control in order to reduce (when S k(i,j) > 0 and V r <V * (i) (t n + t pn )) or raise (when S k(i,j) < 0 and V r >V * (i) (t n + t pn )) the voltage magnitude of bus i towards its desired value. The maximum voltage variation value of bus i: where u (j)min is the j-th element in vector u min .
The maximum voltage variation of the i-th evaluation target considering all candidate control actions applied at sampling point k can be expressed as: (4) Considering multi control steps applied at the first m sampling points, The maximum voltage regulation ability at time instant t n + t pn caused by the control sequence (∆u 1 ,∆u 2 , . . . , ∆u m ) is evaluated as follows: where V a(i)m is the evaluated voltage regulation ability of the i-th evaluation target when the length of control horizon is mt s .

Determination of the Horizon Parameters
Based on the voltage regulation ability evaluation, an applicability judgment index for control horizon parameter settings is proposed: where I s(i)m is the applicability judgment index of the i-th evaluation target when the control horizon is set to mt s . It should be noted that the proposed scheme will not be triggered if V r =V * (i) (t n + t pn ). When the predicted voltage without control is lower than its reference value (V r >V * (i) (t n + t pn )), Equation (13) can be written as: Similarly, when V r <V * (i) (t n + t pn ), Equation (13) can be written as: If I s(i)m > 0, it indicates that the controller is able to bring the voltage magnitude of bus i back to its desired range with the current horizon parameter settings. Let I sm be a vector that consists the applicability judgment index of all the evaluation targets, I sm = (I s(1)m , I s(2)m , . . . . . . , I s(Ni)m ). The initial value of m is set to 1. If ∀I s(i)m ≤ 0 in vector I sm , the proposed algorithm will keep increasing m until all the elements in I sm meet I s(i)m > 0, then mt s is determined as the control horizon of the n-th optimization and the prediction horizon of the (n + 1)-th optimization, mt s = t cn = t p(n+1) .
If it is evaluated that the candidate control variables cannot achieve ∀I s(i)m > 0 according to Equations (13)- (16), m will be reset to 1, and control actions with high cost such as load shedding will be added to the candidate control actions. In this paper, it is assumed that the power system is equipped with adequate reactive power compensation and voltage regulation devices for a voltage correction scenario.
As is presented in Figure 3, the length of the control horizon and prediction horizon are adaptively adjusted during the receding-horizon optimization process.
of m is set to 1. If ∀ Is(i)m ≤ 0 in vector Ism, the proposed algorithm will keep increasing m until all the elements in Ism meet Is(i)m > 0, then mts is determined as the control horizon of the n-th optimization and the prediction horizon of the (n + 1)-th optimization, mts = tcn = tp(n+1).
If it is evaluated that the candidate control variables cannot achieve ∀ Is(i)m > 0 according to Equations (13)- (16), m will be reset to 1, and control actions with high cost such as load shedding will be added to the candidate control actions. In this paper, it is assumed that the power system is equipped with adequate reactive power compensation and voltage regulation devices for a voltage correction scenario.
As is presented in Figure 3, the length of the control horizon and prediction horizon are adaptively adjusted during the receding-horizon optimization process.

Voltage Correction Scheme Based on AH-MPC
The structure of the proposed scheme is shown in Figure 4. The system measurements are collected and processed by a state estimator in the control center. The results are then used to initialize the state variable of the system model and predict the future system trajectory. Compared with a standard MPC controller, the byproduct of the voltage prediction results are used to update the horizon parameter settings. With the numerical simulation results of * V and the trajectory sensitivity information already available, the calculation of the horizon parameters based on Equations (13)-(17) requires solving a set of linear algebraic equations, which involves little extra computation.
With adaptive horizons applied, the objective function (6) and constraint (8) in the standard MPC-based optimization model are rewritten as: where Npn and Ncn are the number of prediction steps and control steps applied in the n-th optimization. Figure 5 shows the flow chart for the voltage correction scheme based on AH-MPC. The procedure of the proposed scheme in the n-th optimization include the following.


At a sampling point tn, trapezoidal numerical simulation is used to obtain the predicted system trajectory within the prediction horizon [tn, tn + tpn].

Voltage Correction Scheme Based on AH-MPC
The structure of the proposed scheme is shown in Figure 4. The system measurements are collected and processed by a state estimator in the control center. The results are then used to initialize the state variable of the system model and predict the future system trajectory. Compared with a standard MPC controller, the byproduct of the voltage prediction results are used to update the horizon parameter settings. With the numerical simulation results ofV * and the trajectory sensitivity information already available, the calculation of the horizon parameters based on Equations (13)- (17) requires solving a set of linear algebraic equations, which involves little extra computation. With adaptive horizons applied, the objective function (6) and constraint (8) in the standard MPC-based optimization model are rewritten as: N p(n+1) = N cn = t cn /t s (22) where N pn and N cn are the number of prediction steps and control steps applied in the n-th optimization. Figure 5 shows the flow chart for the voltage correction scheme based on AH-MPC. The procedure of the proposed scheme in the n-th optimization include the following.

•
At a sampling point t n , trapezoidal numerical simulation is used to obtain the predicted system trajectory within the prediction horizon [t n , t n + t pn ].

•
The trajectory sensitivity is obtained as a by-product of the first step.
• Based on the predicted system trajectory and the trajectory sensitivity results, the adaptive horizon parameters are determined as is introduced in Section 3.

•
The optimal control model is formulated according to Equations (7), (9), (10), (20), (21) and (22). • By solving the optimal control model, the first step of the so computed control sequence is applied in the system at t n + t s .

•
The above procedure is repeated at the initial time of the (n + 1)-th optimization.
If disturbances are detected by phase measurement unit (PMU) measurements and monitoring systems, all the parameters are initialized at the next sampling point.


The trajectory sensitivity is obtained as a by-product of the first step.  Based on the predicted system trajectory and the trajectory sensitivity results, the adaptive horizon parameters are determined as is introduced in Section 3.


The optimal control model is formulated according to Equations (7), (9), (10), (20), (21) and (22).  By solving the optimal control model, the first step of the so computed control sequence is applied in the system at tn + ts.  The above procedure is repeated at the initial time of the (n + 1)-th optimization.
If disturbances are detected by phase measurement unit (PMU) measurements and monitoring systems, all the parameters are initialized at the next sampling point.   The trajectory sensitivity is obtained as a by-product of the first step.  Based on the predicted system trajectory and the trajectory sensitivity results, the adaptive horizon parameters are determined as is introduced in Section 3.


The optimal control model is formulated according to Equations (7), (9), (10), (20), (21) and (22).  By solving the optimal control model, the first step of the so computed control sequence is applied in the system at tn + ts.  The above procedure is repeated at the initial time of the (n + 1)-th optimization.
If disturbances are detected by phase measurement unit (PMU) measurements and monitoring systems, all the parameters are initialized at the next sampling point.

Simulations
The proposed method is applied to the nine-bus system and the New England test system. Dynamic recovery load models [23,24] are used in the simulation, with α s = β s = 0, α t = β t = 2 and T p = T q = 60. The voltage reference V r = 1, and the voltage control deadband is 0.1 p.u., which indicates that the voltage magnitude should be kept within [0.95 p.u., 1.05 p.u.] with control. The sampling interval t s is 10 s. The initial prediction horizon t p1 is 60 seconds. The diagonal elements of Q and R are set to be 1 because the correction of each bus voltage is equally important and no expensive control actions such as load shedding are considered in the simulation cases. The trajectory prediction and adaptive horizon parameter determination is realized based on the quasi-steady state model [1], and the actual control performance is illustrated by the full-time domain simulation on a detailed system model.
The model prediction was realized with the Matlab/Simulink-based tool. The optimization problem was solved by GAMS [25]. The calculations were completed using a windows machine with Intel Core 2.50 GHz CPU and 8 GB of RAM.

Case A: Nine-Bus Test System
Simulation on a simple test system is presented in this section to demonstrate the basic idea and implementation procedure of the proposed scheme. The topology of the test system is shown in Figure 6. The available candidate control actions include the shunt capacitors at bus6, the AVR set points of generator at bus2 and bus3, the OLTCs of the transformer between bus9 and bus3. The constraint of control variation in one control step for the shunt capacitors, AVR set points and the OLTCs are 0.1 p.u., 0.04 p.u. and 1.67% respectively.

Appl. Sci. 2018, 8, x FOR PEER REVIEW 10 of 18
The proposed method is applied to the nine-bus system and the New England test system. Dynamic recovery load models [23,24] are used in the simulation, with αs = βs = 0, αt = βt = 2 and Tp = Tq = 60. The voltage reference Vr = 1, and the voltage control deadband is 0.1 p.u., which indicates that the voltage magnitude should be kept within [0.95 p.u., 1.05 p.u.] with control. The sampling interval ts is 10 s. The initial prediction horizon tp1 is 60 seconds. The diagonal elements of Q and R are set to be 1 because the correction of each bus voltage is equally important and no expensive control actions such as load shedding are considered in the simulation cases. The trajectory prediction and adaptive horizon parameter determination is realized based on the quasi-steady state model [1], and the actual control performance is illustrated by the full-time domain simulation on a detailed system model.
The model prediction was realized with the Matlab/Simulink-based tool. The optimization problem was solved by GAMS [25]. The calculations were completed using a windows machine with Intel Core 2.50 GHz CPU and 8 GB of RAM.

Case A: Nine-Bus Test System
Simulation on a simple test system is presented in this section to demonstrate the basic idea and implementation procedure of the proposed scheme. The topology of the test system is shown in Figure 6. The available candidate control actions include the shunt capacitors at bus6, the AVR set points of generator at bus2 and bus3, the OLTCs of the transformer between bus9 and bus3. The constraint of control variation in one control step for the shunt capacitors, AVR set points and the OLTCs are 0.1 p.u., 0.04 p.u. and 1.67% respectively. The voltage evolution with the proposed control scheme is shown in Figure 8. The step-like voltage variations are caused by the discrete control feature of the reactive power compensation and voltage regulation devices, which is accordant with the voltage response feature in real-life applications. At the initial time instant t1 = 10 s, the proposed controller predicts the system trajectory over a predefined prediction horizon tp1 = 60 s. The voltage magnitude of load bus 6 and bus 8 fall out of the lower limit 0.95 p.u at the end of the prediction horizon t1 + tp1 = 70 s, thus the voltage magnitudes of these two buses are determined as the evaluation targets. The horizon parameters are updated at the initial time of each optimization according to the applicability judgment index Ism. The value of Ism and the applied horizon parameters in this case is listed in Table 1, where n represents the n-th optimization. The voltage evolution with the proposed control scheme is shown in Figure 8. The step-like voltage variations are caused by the discrete control feature of the reactive power compensation and voltage regulation devices, which is accordant with the voltage response feature in real-life applications. At the initial time instant t 1 = 10 s, the proposed controller predicts the system trajectory over a predefined prediction horizon t p1 = 60 s. The voltage magnitude of load bus 6 and bus 8 fall out of the lower limit 0.95 p.u at the end of the prediction horizon t 1 + t p1 = 70 s, thus the voltage magnitudes of these two buses are determined as the evaluation targets. The horizon parameters are updated at the initial time of each optimization according to the applicability judgment index I sm . The value of I sm and the applied horizon parameters in this case is listed in Table 1, where n represents the n-th optimization.    10 20 When n = 1 and m = 1, the applicability judgment index Is1 is obtain according to Equation (17). As is shown in Table 1, the value of the element Is(1)1 = −0.00195 < 0, thus the proposed controller increases m and calculate Is2. Since all the elements of Is2 are positive, the control horizon of the current optimization tc1 and the prediction horizon of the next optimization tp2 are both determined as 20 s. The control model of the current optimization is then formulated and solved. Considering the time consuming of communication and calculation, the above procedure is started at 10 s, and the first step of the so-computed control sequence is applied a sampling interval later. When n = 2 and m = 1, all the elements of Is1 are positive, thus tc2 = tp3 = 10 s. According to the system trajectory prediction implemented during the third optimization, the voltage magnitudes of all the system buses has restored to the desired range and stops declining, the receding-horizon optimization is then completed. The adaptive horizons applied in the simulation is graphically illustrated in Figure 9. The control actions are listed in Table 2, where Uref,Gi represents the AVR set point of generator i, Ci is the applied shunt capacitor at bus i, nt(i − j) is the ratio of the transformer connecting bus i and j.    10 20 When n = 1 and m = 1, the applicability judgment index Is1 is obtain according to Equation (17). As is shown in Table 1, the value of the element Is(1)1 = −0.00195 < 0, thus the proposed controller increases m and calculate Is2. Since all the elements of Is2 are positive, the control horizon of the current optimization tc1 and the prediction horizon of the next optimization tp2 are both determined as 20 s. The control model of the current optimization is then formulated and solved. Considering the time consuming of communication and calculation, the above procedure is started at 10 s, and the first step of the so-computed control sequence is applied a sampling interval later. When n = 2 and m = 1, all the elements of Is1 are positive, thus tc2 = tp3 = 10 s. According to the system trajectory prediction implemented during the third optimization, the voltage magnitudes of all the system buses has restored to the desired range and stops declining, the receding-horizon optimization is then completed. The adaptive horizons applied in the simulation is graphically illustrated in Figure 9. The control actions are listed in Table 2, where Uref,Gi represents the AVR set point of generator i, Ci is the applied shunt capacitor at bus i, nt(i − j) is the ratio of the transformer connecting bus i and j.  When n = 1 and m = 1, the applicability judgment index I s1 is obtain according to Equation (17). As is shown in Table 1, the value of the element I s(1)1 = −0.00195 < 0, thus the proposed controller increases m and calculate I s2 . Since all the elements of I s2 are positive, the control horizon of the current optimization t c1 and the prediction horizon of the next optimization t p2 are both determined as 20 s. The control model of the current optimization is then formulated and solved. Considering the time consuming of communication and calculation, the above procedure is started at 10 s, and the first step of the so-computed control sequence is applied a sampling interval later. When n = 2 and m = 1, all the elements of I s1 are positive, thus t c2 = t p3 = 10 s. According to the system trajectory prediction implemented during the third optimization, the voltage magnitudes of all the system buses has restored to the desired range and stops declining, the receding-horizon optimization is then completed. The adaptive horizons applied in the simulation is graphically illustrated in Figure 9.
The control actions are listed in Table 2, where U ref,Gi represents the AVR set point of generator i, C i is the applied shunt capacitor at bus i, n t(i − j) is the ratio of the transformer connecting bus i and j.

Case B: New England Test System
In this case, the validation and effectiveness of the proposed method further verified by the simulation on the New England 39-bus test system. Figure 10 shows the topology of the test system. The candidate control actions include the shunt capacitors at bus 8, the AVR set points of all the 10 generators, the OLTCs of the transformers connecting bus12 and 11, 12 and 13, 19 and 20. The control parameter settings are the same as case A.
At t = 10 s, generator 3 are tripped. As is shown in Figure 11.The voltage magnitude of load bus 4, 7, 8 and 12 drops below 0.95 p.u. after the fault occurs.

Case B: New England Test System
In this case, the validation and effectiveness of the proposed method further verified by the simulation on the New England 39-bus test system. Figure 10 shows the topology of the test system. The candidate control actions include the shunt capacitors at bus 8, the AVR set points of all the 10 generators, the OLTCs of the transformers connecting bus12 and 11, 12 and 13, 19 and 20. The control parameter settings are the same as case A.
At t = 10 s, generator 3 are tripped. As is shown in Figure 11.The voltage magnitude of load bus 4, 7, 8 and 12 drops below 0.95 p.u. after the fault occurs.

Case B: New England Test System
In this case, the validation and effectiveness of the proposed method further verified by the simulation on the New England 39-bus test system. Figure 10 shows the topology of the test system. The candidate control actions include the shunt capacitors at bus 8, the AVR set points of all the 10 generators, the OLTCs of the transformers connecting bus12 and 11, 12 and 13, 19 and 20. The control parameter settings are the same as case A.
At t = 10 s, generator 3 are tripped. As is shown in Figure 11.The voltage magnitude of load bus 4, 7, 8 and 12 drops below 0.95 p.u. after the fault occurs.  According to Equation (11), the voltage magnitudes of load bus 4, 7, 8 and 12 are determined as the evaluation targets. Figure 12 presents the corrected voltage trajectory with the proposed scheme started at t = 10 s. The applied control actions are listed in Table 3. As the receding horizon optimization progresses, the system voltages are gradually brought back inside the limits at t = 40 s. The applicability judgment index Ism is listed in Table 4. When large voltage deviation of the control targets is detected at the initial control instant t = 10 s, the proposed scheme applies three-step optimization according to the value of Is3. As the voltages gradually approaching the desired range, the declining voltage trajectories are stabilized, and smaller control variations are needed to complete the optimization process. In the third optimization, the control variation needed is small enough to be applied in one step without causing negative effect listed in [21], and shorter prediction length beyond the control horizon is needed due to the steadier voltage evolution trend, which reduces the computation time and settling times of the entire optimization process. The applied decreasing prediction horizon and control horizon are graphically illustrated in Figure 13.  According to Equation (11), the voltage magnitudes of load bus 4, 7, 8 and 12 are determined as the evaluation targets. Figure 12 presents the corrected voltage trajectory with the proposed scheme started at t = 10 s. The applied control actions are listed in Table 3. As the receding horizon optimization progresses, the system voltages are gradually brought back inside the limits at t = 40 s. The applicability judgment index I sm is listed in Table 4. When large voltage deviation of the control targets is detected at the initial control instant t = 10 s, the proposed scheme applies three-step optimization according to the value of I s3 . As the voltages gradually approaching the desired range, the declining voltage trajectories are stabilized, and smaller control variations are needed to complete the optimization process. In the third optimization, the control variation needed is small enough to be applied in one step without causing negative effect listed in [21], and shorter prediction length beyond the control horizon is needed due to the steadier voltage evolution trend, which reduces the computation time and settling times of the entire optimization process. The applied decreasing prediction horizon and control horizon are graphically illustrated in Figure 13. According to Equation (11), the voltage magnitudes of load bus 4, 7, 8 and 12 are determined as the evaluation targets. Figure 12 presents the corrected voltage trajectory with the proposed scheme started at t = 10 s. The applied control actions are listed in Table 3. As the receding horizon optimization progresses, the system voltages are gradually brought back inside the limits at t = 40 s. The applicability judgment index Ism is listed in Table 4. When large voltage deviation of the control targets is detected at the initial control instant t = 10 s, the proposed scheme applies three-step optimization according to the value of Is3. As the voltages gradually approaching the desired range, the declining voltage trajectories are stabilized, and smaller control variations are needed to complete the optimization process. In the third optimization, the control variation needed is small enough to be applied in one step without causing negative effect listed in [21], and shorter prediction length beyond the control horizon is needed due to the steadier voltage evolution trend, which reduces the computation time and settling times of the entire optimization process. The applied decreasing prediction horizon and control horizon are graphically illustrated in Figure 13.      Figure 13. Adaptive horizons applied in case B.

Simulating Modeling Inaccuracies
In this section, the modeling inaccuracies and disturbance are considered in the course of controlling the system. The fault scenario, system parameter and other controller parameter settings are the same as case B. At t = 25 s, A disturbance is imposed to the load at bus 3, bus 4, bus 6 and bus7 with 5% increase of the active and reactive power.
If the disturbance can be measured and sent to the control center, the voltage behavior based on AH-MPC is shown in Figure 14. The system voltage drops dramatically due to the disturbance at t = 25. At the next sampling point t = 30 s, the optimization is carried out based on the updated system measurements.
V (p.u.) Figure 14. Voltage evolution using AH-MPC with disturbance occurred during the control process.
If the disturbance is ignored by the control center, there will be modeling inaccuracy problems at all future control steps after t = 25 s. The voltage behavior based on AH-MPC is shown in Figure  15. Although the prediction model is not completely updated after the disturbance, the controller will

Simulating Modeling Inaccuracies
In this section, the modeling inaccuracies and disturbance are considered in the course of controlling the system. The fault scenario, system parameter and other controller parameter settings are the same as case B. At t = 25 s, A disturbance is imposed to the load at bus 3, bus 4, bus 6 and bus7 with 5% increase of the active and reactive power.
If the disturbance can be measured and sent to the control center, the voltage behavior based on AH-MPC is shown in Figure 14. The system voltage drops dramatically due to the disturbance at t = 25 s. At the next sampling point t = 30 s, the optimization is carried out based on the updated system measurements.   Figure 13. Adaptive horizons applied in case B.

Simulating Modeling Inaccuracies
In this section, the modeling inaccuracies and disturbance are considered in the course of controlling the system. The fault scenario, system parameter and other controller parameter settings are the same as case B. At t = 25 s, A disturbance is imposed to the load at bus 3, bus 4, bus 6 and bus7 with 5% increase of the active and reactive power.
If the disturbance can be measured and sent to the control center, the voltage behavior based on AH-MPC is shown in Figure 14. The system voltage drops dramatically due to the disturbance at t = 25. At the next sampling point t = 30 s, the optimization is carried out based on the updated system measurements.
V (p.u.) Figure 14. Voltage evolution using AH-MPC with disturbance occurred during the control process.
If the disturbance is ignored by the control center, there will be modeling inaccuracy problems at all future control steps after t = 25 s. The voltage behavior based on AH-MPC is shown in Figure  15. Although the prediction model is not completely updated after the disturbance, the controller will If the disturbance is ignored by the control center, there will be modeling inaccuracy problems at all future control steps after t = 25 s. The voltage behavior based on AH-MPC is shown in Figure 15. Although the prediction model is not completely updated after the disturbance, the controller will get get the actual voltage responses according to measured feedback signals at each future sampling point, and the structure of receding-horizon multi-step optimization offers chance to correct the consequences of modeling errors.
From Figures 14 and 15, we can see that the proposed control scheme is still valid. This is mainly because the disturbance and modeling inaccuracy problems are compensated by the closed-loop feature of the centralized controller. V (p.u.) Figure 15. Voltage evolution using AH-MPC with both disturbance and model mismatch occurred during the control process.

Control Performance Comparison
The control performance of the proposed method is compared with the ordinary MPC-based method using constant horizon parameters (t c = 30 s, t p = 60 s). The fault scenario, system parameter and other controller parameter settings are the same as case B. The control result is shown in Figure  16. All the voltages meet the constraint requirements after 4 optimizations are implemented. The performance comparison is listed in Table 5. An index of average voltage deviation of all the load buses is introduced as:  From Figures 14 and 15, we can see that the proposed control scheme is still valid. This is mainly because the disturbance and modeling inaccuracy problems are compensated by the closed-loop feature of the centralized controller.

Control Performance Comparison
The control performance of the proposed method is compared with the ordinary MPC-based method using constant horizon parameters (t c = 30 s, t p = 60 s). The fault scenario, system parameter and other controller parameter settings are the same as case B. The control result is shown in Figure 16. All the voltages meet the constraint requirements after 4 optimizations are implemented. The performance comparison is listed in Table 5. An index of average voltage deviation of all the load buses is introduced as: get the actual voltage responses according to measured feedback signals at each future sampling point, and the structure of receding-horizon multi-step optimization offers chance to correct the consequences of modeling errors. From Figures 14 and 15, we can see that the proposed control scheme is still valid. This is mainly because the disturbance and modeling inaccuracy problems are compensated by the closed-loop feature of the centralized controller. V (p.u.) Figure 15. Voltage evolution using AH-MPC with both disturbance and model mismatch occurred during the control process.

Control Performance Comparison
The control performance of the proposed method is compared with the ordinary MPC-based method using constant horizon parameters (t c = 30 s, t p = 60 s). The fault scenario, system parameter and other controller parameter settings are the same as case B. The control result is shown in Figure  16. All the voltages meet the constraint requirements after 4 optimizations are implemented. The performance comparison is listed in Table 5. An index of average voltage deviation of all the load buses is introduced as:   The integration time ∆t = 300 s. Smaller value of ∆V offset indicates faster voltage recovery speed and better control performance.
The comparison results of Figures 12 and 16 and Table 5 show that both methods are able to stabilize the system voltage within a desire range, and the AH-MPC-based method outperforms MPC-based scheme with constant horizon in computation effort and voltage recovery speed.
The proposed scheme applies decreasing prediction parameters and control horizon. When n = 1, both methods applies the same horizon parameters and cost almost the same computation time.
With the increase of n, the time consuming per optimization of the AH-MPC-based scheme is remarkably reduced. This is mainly because the numerical simulation needed in the model prediction is implemented within a shorter prediction horizon, and the control sequence consists of less control steps is computed. Moreover, the proposed scheme takes less settling times, thus there are less optimal control model to be formulated and solved, which also leads to decreased computational effort.

Discussion
A principle of determine the horizon parameters based on measurements and prediction is developed in this paper, which is integrated with the standard MPC-based method as is graphically illustrated in Figure 4. The determination of horizon parameters costs little computation effort, and the performance of the controller can be improved, which has been verified by the simulations. It should be noted that the basic idea of adaptive horizon parameter settings is not restricted to the method presented in Section 2, but can be combined with any MPC-based method using linear prediction model.
The time-consuming nature of the optimization is affected by many aspects, such as hardware configuration, system scale and prediction method. Simulation results of a load increase scenario in IEEE 118 bus test system [26] is presented to further verify the practicability of the proposed method. The active and reactive powers of loads at buses 82, 83, 84, 85 and 86 are linearly increased at a rate of 10% per second until t = 20 s. The voltage trajectories of bus 83 with and without control are shown in Figure 17. With the actions of the AVR set points of generators at bus 77, 85, 87, 89 and 90, the system voltages are kept within a desired range. The maximum computation time for one optimization (when n = 1) is 5.18 s. The comparison results of Figures 12 and 16 and Table 5 show that both methods are able to stabilize the system voltage within a desire range, and the AH-MPC-based method outperforms MPC-based scheme with constant horizon in computation effort and voltage recovery speed.
The proposed scheme applies decreasing prediction parameters and control horizon. When n = 1, both methods applies the same horizon parameters and cost almost the same computation time. With the increase of n, the time consuming per optimization of the AH-MPC-based scheme is remarkably reduced. This is mainly because the numerical simulation needed in the model prediction is implemented within a shorter prediction horizon, and the control sequence consists of less control steps is computed. Moreover, the proposed scheme takes less settling times, thus there are less optimal control model to be formulated and solved, which also leads to decreased computational effort.

Discussion
A principle of determine the horizon parameters based on measurements and prediction is developed in this paper, which is integrated with the standard MPC-based method as is graphically illustrated in Figure 4. The determination of horizon parameters costs little computation effort, and the performance of the controller can be improved, which has been verified by the simulations. It should be noted that the basic idea of adaptive horizon parameter settings is not restricted to the method presented in Section 2, but can be combined with any MPC-based method using linear prediction model.
The time-consuming nature of the optimization is affected by many aspects, such as hardware configuration, system scale and prediction method. Simulation results of a load increase scenario in IEEE 118 bus test system [26] is presented to further verify the practicability of the proposed method. The active and reactive powers of loads at buses 82, 83, 84, 85 and 86 are linearly increased at a rate of 10% per second until t = 20 s. The voltage trajectories of bus 83 with and without control are shown in Figure 17. With the actions of the AVR set points of generators at bus 77, 85, 87, 89 and 90, the system voltages are kept within a desired range. The maximum computation time for one optimization (when n = 1) is 5.18 s. The measurement collection and optimization starts at each sampling point, and the control decisions are applied to the practical system a sampling interval (t s = 10 s) later, which offers sufficient time for the telecommunication and computation. In a real-life application, it is confidently expected that the computation time can be further reduced:

•
The simulations were carried out on a personal computer with the Matlab/Simulink-based tool. Both the hardware and software configuration can be significantly improved in a real-life application [27].

•
The optimization problem can be solved concentrating on areas close to the disturbance with a simple representation used for other areas [21].

Conclusions
This paper proposes a voltage correction scheme based on AH-MPC. By coordinating the control actions of the AVR set points, OLTCs and shunt capacitors, the system voltages are kept within a desired range after the disturbance occurs. The inevitable prediction model mismatch problem is compensated by the closed-loop feature of the centralized controller. This paper also presents a real-time determination method for horizon parameters which is commonly ignored in the related research. The horizon parameters are updated at each sampling point according to the real-time feedback signals and system prediction results, which effectively improves the control performance compared with the traditional MPC-based method. The information needed to determine the horizon parameters can be obtained as by-products of the prediction process in a standard MPC-based method, thus little extra computation time is needed.