An Optimal Control Scheme for Load Bus Voltage Regulation and Reactive Power-Sharing in an Islanded Microgrid

: This paper presents an optimal control scheme for an islanded microgrid (MG), which performs reactive power-sharing and voltage regulation. Two-fold objectives are achieved, i.e., the load estimation strategy, firstly, approximates the MG’s impedance and transmits this information through a communication link. Based on approximated impedance information, an optimal regulator is then constructed to send optimal control commands to respective local power controllers of each distributed generation unit. An optimal regulator is a constraints optimized problem, mainly responsible to restore the buses’ voltage magnitudes and realize power-sharing proportionally. The important aspect of this control approach is that the voltage magnitude information is only required to be transferred to each inverter’s controller. In parallel, a secondary control layer for frequency restoration is implemented to minimize the system frequency deviations. The MATLAB/Simulink and experimental results obtained under load disturbances show the effectiveness for optimizing the voltage and power. Modeling and analysis are also verified through stability analysis using system-wide mathematical small-signal models.


Introduction
One of the important factors of autonomous MG is to drive the electrical network with high security and reliability at various operating conditions. Control of voltage and reactive power is a salient step in order to obtain a reliable smart network [1][2][3]. Such issues have been investigated by many researchers as it is considered one of the leading cause for voltage collapse, which may drive the power system blackout [4]. Voltage vulnerability can be defined as the failure of bus voltage to return at a nominal value under disturbance [5,6]. Most of the power shutdown came about from voltage uncertainty, which predominantly came up due to failure of the control system to draw out adequate reactive power to reinforce the voltage at fault-finding grid buses [7]. The issue principally emerged from depending on an on-off control facility, without using an automatic control system. The design of voltage control and the power-sharing layer was first set up by European countries and later propagated by the United States, Asia, and Africa, respectively [4]. Each country employed the control layers in accordance with their own methodology, which is sufficient to remove the voltage fluctuations and well address the reactive power sharing.
The real power-frequency p-f and reactive power-voltage Q-V droop layers are employed to resolve the power control issues [6,[32][33][34][35]. Authors presented an improved power sharing control strategy to regulate the load bus voltage but at the cost of inverter terminal voltage deviation [27]. A "Q-V dot droop" scheme was developed but the reactive power sharing was not evident when the local loads were added [36]. The droop control schemes based on virtual impedance are observed as an evident strategy to deal with reactive power-sharing problems [30,31,37,38]. The power sharing can be improved by virtual impedance techniques, but voltage droop and virtual impedance deteriorate the inverter terminal voltage quality [39,40]. Voltage control methods are presented in [35] in order to reduce the trade-off between bus voltage and reactive power. Such schemes develop additional voltage offset for compensation by depending on the communication link, which may lead to voltage instability due to time delay. A distributed voltage control scheme is developed to regulate system voltage without a communication link, but the reactive power-sharing performance has not. Recently, a Kalman filter-based state estimator scheme was used for reactive power sharing and bus voltage regulation, however, this strategy requires a high bandwidth data rate and thus increases complexity when generation units are located at long distances [41,42].
Power systems are predominantly working under stressed operating conditions due to speedy growing load demand. The voltage stability has set off vital concerns in presentday electrical power system operations which depends on the proportional realization of reactive power. Even the simplest version of the distribution system faces challenges such as load bus voltage restoration and at the same time minimizing the circulating reactive power in grid-forming nodes. Therefore, in this paper, an effort has been made to resolve such issues faced by an islanded MG distribution network. To the best of our knowledge, the simple quadratic optimized cost function-based techniques for optimizing the load bus voltage and power are not fully investigated and further research is required to enhance the existing studies. In some articles, various methods have been sought to achieve these aims by acquiring real-time signals at each load bus, however, it overburdens communication by taking up high communication bandwidth. On that account, this paper presents an optimal control method, an extension of previous work [34]. The presented method does not require accurate grid impedance information, eventually reduces the data bandwidth. We summarize the major contribution of this work as follows:


This study is conducted in a radial feeder system accompanied by three load buses connected with two grid forming nodes.  The approximated impedance information of all load-buses is estimated through local load agents, which reduces the bandwidth data requirement.


An optimal regulator, an optimized cost function, is constructed to send optimal control commands to the inner layer of each DG unit. Moreover, the secondary control layer is responsible to restore the frequency deviations.  The MATLAB/Simulink (R2018a, product registered with TJU, China) and experiment results show the effectiveness of the proposed methodology. For the stability analysis, separate Simulink and a linear analysis tool are considered to analyze the complex system by perturbing dynamical equations of the same.
The remainder of this paper is organized as follows. Section 2 presents the design procedure with details of the system used in the proposed study. Section 3 outlines the operation principle of the proposed optimal control strategy. Section 4 provides the derivations of the small-signal model for this MGs system with the employed controls. Section 5 elaborates the stability analysis of the system under the discussed controls. Further, it also presents the simulation and experimental results in detail with and without the proposed control scheme, and finally, Section 6 concludes the paper.

Design Procedure
This section details the general procedural design flow of the proposed methodology, which is primarily composed of load estimation and an optimal regulator, illustrated by Figure 1. In the load estimation strategy, the load impedance of i th (i = 1, 2, 3) bus is calculated through measured voltage and current at each bus node. Based on such impedance information, the grid impedance is approximated in terms of one node, as discussed in detail in Section 3.1. An optimal regulator sends the control commands uci (uc1, uc2) to inner control layers, according to the received approximated grid impedance information. To adequately describe the problem and presented technique, we now systematically go about describing the mechanism of optimizing power and voltage, power flow control, and power steady-state equations with the necessary mathematical representations.

Proposed Optimal Control Strategy
This section describes the operation principle of the proposed methodology and details the associated MG setup employed in this research. For simulation studies, a radial type of feeder is used as shown in Figure 2, where two DG units are interfaced using three phase, three-wire power electronics inverters, feeding three load buses Vbus1, Vbus2, and Vbus3 through feeder impedances with connected loads, i.e., Zload1, Zload2, and Zload3. The working principle of load estimation and an optimal regulator is discussed in Sections 3.1 and 3.2, respectively. System control parameters for the stability analysis, simulation, and experiment are given in Tables 1 and 2

Load Estimation: The Approximation of Grid Impedance
The load estimation strategy does not require the accurate grid impedance information of each load bus, as the computational complexity increases proportional to the number of load buses. In this presented study, three load buses connected in a radial configuration are considered, as illustrated in Figure 3. The technique, firstly sense the voltage Vbus and current Ibus at each local load node. Sinus-based park transformation is then employed from measured three-phase signals to a dq0 rotating reference frame. The angular position of the rotating frame is given by input ϑi in rad/s, which vary from 0 to 2π. Deduce the dq0 components from three-phase bus voltage signal by performing Clarke transformation, yields vs = (vd + j.vq) = (vα + j.vβ).e -j (ϑ i t-π/2) and vdqo = Kptv.vabc, where, vdq0 = [vd vq v0] t and vab c= [va vb vc] t . The term Kpv denotes the park transformation matrix as given in Appendix A. The dq0 components are expressed as: Similarly, deduce the dq0 components from every i th (i = 1, 2, 3) load current, results is = (id + j.iq) = (iα + j.iβ).e -j (ϑ i t-π/2) and idqo = Kpti.iabc, where, idq0 = [id iq i0] t and iabc = [ia ib ic] t . The dq0 components of idq0 is written as: The load impedance information of the i th bus (i = 1, 2, 3) is estimated by the ratio of Equation (1) to Equation (2), yielding: Once the load impedance of each bus is estimated through local node agents, then the impedance approximation procedure is adopted, expressed by Equations (4)- (7), and illustrated in Figure 3. Wherein, the two optimal control commands uci (uc1, uc2) are the voltage inputs to the inner control layers of DG1 and DG2, as discussed in detail in Section 3.2. For the simulation studies, uc1 of DG1 has been considered as variable frequency while uc2 of DG2 only voltage magnitude is used. Knowing the parameters, e.g., Si, load bus impedance Zloadi (i = 1,2,3), and optimal control commands uci (i = 1,2), we can now approximate the three load buses impedances in term of 1 load node as a mathematical representation depicted in Figure 3.
Where, the ZaLine,i (ZaLine,1, ZaLine,2) and ZaL,i (I = 1) are the approximated transmission lines and load impedances, respectively. Once the grid impedances are approximated then the optimal regulator, discussed in Section 3.2, is constructed to compute the control commands based on received impedance information.

Optimal Regulator and Algorith Flowchart
The proposed optimal regulator, an equality optimized cost function, computes control commands according to MG's approximated load impedance information. To minimize cost of function, different optimization techniques are being employed, e.g., linear, quadratic, and higher-order optimization methods. Solvability of linear optimization is easier, but it only works with linear variables and its problem formulation is peculiar. The solutions of higher-order cost functions are inconvenient. Therefore, in this study, a constraints quadrative optimization problem is considered which is convenient for problem formulation and solvability. Some cost functions have constraints reactive power, i.e., Q1 = Q2 results in a tradeoff between real power and voltage magnitude. Since the reactive power requirement is not stringent, therefore, we considered real power equality con-

S1 S2
straints, i.e., P1 = P2. The desired control commands are acquired by computing the optimization cost function that minimizes the reactive power sharing error ∆Qi and bus voltage error ∆Vi, expressed by Equation (8). 22 1 1 where ωbj and ωQ are the weights of bus voltage and reactive power error, respectively. The term Vjref = 300 V (j = 1, 2, 3) is the nominal bus voltage, while nb = 3 is the number of load buses that has been chosen for study, as shown in Figure 2. Figure 4 shows the algorithm flow chart to find the minimum of optimized cost function and is expressed by Equation (9): The minima of the augmented function ∧ (v1r, v1q, v2r, λ) are located where all of the partial derivatives, i.e., (∂∧/∂v1r = 0, ∂∧/∂v1q = 0, ∂∧/∂v2r = 0 and ∂∧/∂λ = 0) are equal to zero. As the minima "fvalue" of function and control inputs uci, corresponding to different disturbance loads, are illustrated by Figure 4b,c. This algorithm numerically approximates partial derivatives as opposed to analytical evaluation. The important aspect of the numerical derivative strategy is that it can work for multiple variables without even changing the code. Computing the algorithm, the optimal regulator sends the required control commands uci (uc1, uc2) to respective power controllers of the DG units, and the system works the same until the upcoming sampling update. The term uci compensates bus voltage deviation and shares reactive power proportionally. Periodically, the optimal regulator revises its original control plan due to measurements feedback, in line with the latest approximated grid impedance information.

MG Power Flow Control
The complex power delivered is expressed in terms of network voltages and admittances. Single-phase power relations are utilized as all the quantities are assumed to be in per unit. Therefore, the complex power delivered to the ith bus is written as: where the term Vi is a phasor term, and the magnitude and angle are expressed by Vi = |Vi|∠φi. The term Yij is the admittance matrix element, define Gij and Bij as the real and imaginary parts, i.e., Yi = Gij + jBij. By applying the algebraic multiplication to Equation (A6), and then collecting the real part Pi and imaginary part Qi, results in: Power angle φ in medium voltage lines are small, and we can assume sinφ = φ and cosφ = 1. Therefore, Equations (12) and (13) are re-expressed as: where φij = XPi/ViVj and Vi-Vj = XQi/Vi. Equations (14) and (15) illustrate a direct relationship among the real power Pi and power angle φ as well as reactive power Qi and voltage magnitude (Vi-Vj).

Mechanism of Phasor Implementation and Reactive Power Sharing
The output (uc1, uc2) of the proposed optimal regulator is a set of voltage phasors. The direct implementation of voltage phasers of i th inverter in MG requires very fast and reliable communication infrastructure. Nevertheless, in the proposed methodology, the phases have been adjusted through active power and frequency (p-w) droop control loops, only voltage magnitudes are updated in accordance with optimal regulator output. This combination will realize real power proportionally, which is also considered as equality constraint (Pi -Pi + 1 = 0) in an optimized cost function, thus the accurate phase angle will be automatically modified by the system.
The real power and frequency (p-ω) droop expression is illustrated by: where ωref, and mPi are the nominal frequency and droop coefficients, respectively. The slop of frequency droop characteristic is controlled by mPi parameters. Output frequency of inverter connected with dc power source is adjusted by real power and inner controller. Every ith DG unit is composed in its d-q frame, depending on their angle. By transformation equations, shown in Equation (18), each inverter is transferred to the d-q frame. The proposed power controller is shown in Figure 5a, where the real power and frequency (p-ω) droop layer is responsible for sending the angle to the ith inverter as given by Equation (19).
The instantaneous real power (Pi) can be expressed in d-q frame, yields by Pi = (ωc/(s + ωc)) and pi = void. iodi + voqi.ioqi. Here, the term  c is the cutoff frequency of low pass filters.
The vodqi and iodqi are the measured voltage and current in the d-q frame. Block diagrams of voltage and current control layers are illustrated in Figure 5b. The differential algebraic equations of voltage loop are expressed by: where the term dqi  is the auxiliary state variable and wb is the nominal angular frequency. The differential algebraic equations of the current layer are written as: where the dqi  is the auxiliary state variable and ildqi are the quadrature components of ili.

Secondary Layer: The Frequency Restoration Methodology
The aim of the secondary control layer is to observe and minimize the system frequency deviations. The frequency restoration technique is shown in Figure 1 and expressed by, where the terms ωref and ωavg are the nominal and measured frequencies, respectively. The frequency correction term δωi is responsible for regulating the system frequency at scheduled values, i.e., 50 Hz, while Kpf and Ki are the proportional and integral gains, respectively, for controllers.

Small Signal Analysis of the MG System
This section details the stability analysis of the system by varying P-w droop gain mPi and communication delay Td. The intermittent latencies and failures of component communication links may result in power imbalances, deviations in voltages, and frequency. For this purpose, we separately considered the MATLAB/Simulink and linear analysis tool to analyze the complex system by perturbing dynamical equations of the same. The small-signal modeling strategy is based on three important sub-modules, i.e., inverter, network, and loads. State equations of the network and the connected loads with ith DG inverter are presented in the below sections.

Power Controller
In the proposed methodology, the conventional active power and frequency (P-w) droop layer is used, where the instantaneous power p = vodi.iodi + voqi.ioqi and q = vodi.ioqi − voqi.iodi is obtained from the measured voltage and current signals. The small-signal model of active power can be obtained by linearization, yielding: where vodqi = [vodi voqi] T and iodqi = [iodi ioqi] T . Since, the traditional reactive power and voltage magnitude (Q-V) droop layer is replaced with "optimal regulator", as discussed in Section 3.2. Therefore, the state equations for the modified layer of the power controller are described below, accordingly. The small-signal modeling of optimized cost function approximates complicated multivariable function J near vector inputs V*. This task can be achieved more efficiently by quadratic approximation using the information provided by the second partial derivatives as compared to local linearization. The approximation of multi-dimensional input of function J, Equation (9), is expressed by: The vector V* is the optimal values ** 12 [ ] T vv while vector V are the variable inputs which have the same matrix size of vector V*. Finally, perturbation in ∆v1 and ∆v2 with respective to approximated load impedance variation ∆ZL can be achieved, as expressed by Equations (28)- (34), at optimal values of V* and respective load impedance ZLi.
Further, the mathematical modeling of gird side filter, complete ith inverter, combined model of N inverters, and network load model are given in Appendix C, Appendix D, Appendix E, and Appendix F, respectively.

Complete MG Model
The node voltages are used as inputs to every subsystem. In the assurance of a welldefined node voltage, a virtual resistor rN is supposed among every m node and ground; its symbolic form is defined by Diagonal elements of matrix RN are equal to rN with a matrix size of 2 m * 2 m. The mapping matrix MINV and Mload are the size of 2 m *2 m and 2 m*2 p, respectively. M NET matrix maps the connecting transmission lines onto the network nodes having a size of 2 m*2 n. Finally, the complete MG small-signal model is given below.

Stability, Simulation, and Experimental Verification
This section details the stability analysis, simulation, and experimental results. The simulations are conducted on MATLAB/Simulink on the circuit configuration shown in Figure 2. Wherein the two paralleled connected DG1 and DG2 units are feeding three load buses. In contrast, the experiments are conducted for three phases, 50 Hz scaled islanded MG prototype as depicted in Figure 6, where the system consists of two identical paralleled DG units with a maximum rating of 3 A, 30 V. Ethernet module USR TCP 232 is used for serial communication transmission of data amongst the data terminal equipment (DTE) and data communication equipment (DCE). The whole platform of islanded MGs is controlled by a desktop LabVIEW control system. The system and controller parameters that have been used for the proposed methodology are given in Table 2.

Stability Analysis Results
The complete model of the test system under the presented control scheme is used to analyze stability under varying communication latencies and control gains. MATLAB/Simulink and linear analysis tools have been used in order to obtain modeling results by analyzing the complex system through perturbing dynamical equations. Figure  7 shows the stability plot for the proposed control framework by varying real power droop gain mP to trace the network trajectory. The control values of the droop gain mP where system poles appear to be in the vicinity of the unit origin are considered to be maximum allowable limits. Therefore, using the poles-zero evolutions, the control gain sensitivity, and system stability are predicted. The maximum and minimum mP gain values used for the proposed control scheme are 1x9.5 −09 and 1x9.5 −05 , respectively.  is captured by starting with time delay td = 0 and increasing it step by step to maximum time delay of td = 3 s. From Figure 8, poles are observed in the stable region for the time delay values 0 > td < 1 s and outside the unit circle for the td > 1 s which makes the MG system divergent and unstable.

Case Simulation Studies
The circuit and control parameters for simulation studies are shown in Table 1. The verification is composed of two cases to validate the effectiveness of the proposed control scheme. In case 1, disturbance load Ld2/Rd2 has been added at bus 2, and later in case 2, disturbance loads (Ld1/Rd1, Ld2/Rd2) have been added at bus 1 and bus 2, simultaneously.

Case 1: Power Sharing and Bus 2 Voltage Control
The primary objective of this case study is to observe the power-sharing and load bus 2 voltage performance. The voltage magnitude of bus 2 should be at its original state, i.e., Vref = 300 V, in the presence of load disturbance. Initially, the voltage error weights wb1, wb2, and wb3 for three buses are set at a default value of 1. A disturbance load Ld2/Rd2 with values of Ld2 = 69mH and Rd2 = 9.9 ohms is now added at bus 2 at 0.3 s. The behavior of voltage regulation, active, and reactive power sharing under disturbance load is shown in Figure  9. An obvious active, reactive power sharing error is observed once the Ld2/Rd2 is added as illustrated in Figure 9a,b, which is later compensated to almost zero once the proposed control scheme is activated at 0.3 s as depicted in Figure 9c,d. All three load buses voltage curve behavior is shown in Figure 9f. A rough average of 7 V voltage deviations are observed in Figure 9g which is compensated to almost zero once the proposed strategy is activated. The 6V voltage deviation occurs at bus 2 and even the proposed control scheme is activated, shown in Figure 9e by the red cure. This deviation is removed and bus 2 voltage is regulated by varying the respective voltage error weight wb2, set from 1 to 300, as depicted in Figure 9g.

Load Variation at Bus 1 and Bus 2 Simultaneously
To further investigate the superior control performance of the proposed optimal control strategy two disturbance loads (Ld1 = 69 mH, Rd1 = 9.9 ohm, Ld3= 69 mH and Rd3 = 9.9 ohm) are exerted at bus 1 and bus 3, simultaneously. It is observed that the disturbance leads to a heavy power-sharing error and respective bus voltage drop, shown in Figure  10. Active and reactive power response with conventional control scheme is demonstrated by Figure 10a average of 5 V volts deviation at bus1 and bus 3, and even the control strategy was activated. This deviation is improved and restored to a nominal voltage by varying both buses voltage error weights, set from default value 1 to 300 (ωb1 = 300 and ωb3 = 300) as demonstrated in Figure 10g. Distributed secondary control strategy is used to regulate the system frequency. In Figure 11, it is noticed that activation of the proposed strategy at time 0.02 s, yields frequency restoration within an acceptable range of ±0.5 Hz. Figure 11. System frequency regulation curve.

Experimental Results
An accompanying experiment is conducted to validate simulation results, wherein the disturbance load Ld2/Rd2 is exerted at bus 2 at time 0.02 s. Power-sharing error and voltage deviations caused by disturbance are depicted in Figure 12. Real and reactive power sharing curve response with conventional scheme is demonstrated by Figure 12a,b, respectively. From Figure 12c,d, the real and reactive power sharing errors are compensated almost to zero once the proposed optimal control methodology is activated. All three load buses voltage response curve are collectively shown by Figure 12e, where the initial conditions of parameters are same for both conventional and proposed control strategies. Figure 12f illustrates the average buses voltage response, where it is observed the average bus voltage deviation is restored to its original value (black cure) once the proposed control strategy is activated. Figure 13 shows the variation of optimal control commands uc1 (28.25-0.3j V) and uc2 (27.30 V), before and after the disturbance load exertion. Good agreement is observed by comparing the simulation and experimental results. The effectiveness of the proposed optimal control scheme is thus confirmed.

Conclusions
In this paper, an optimal control technique is presented for islanded MG. This study performs two tasks: load estimation and an optimal regulator. The load estimation strategy approximated the grid impedance and transmitted the information through a communication channel. An equality constraints optimized problem, called an optimal regulator, sends optimal control commands to respective inner voltage control layers of each DG unit. In parallel, a secondary layer is implemented to regulate system frequency. To observe the effectiveness of the proposed methodology, we have exerted the disturbance loads at bus 2 as well as bus 1 and 2, simultaneously. In accordance with the proposed optimized cost function, the MG is at the optimum operating point as verified through the case simulation and hardware studies. Modeling and analysis of three-phase inverterbased MG was also verified through stability analysis using system-wide mathematical small-signal models.

Informed Consent Statement: Not applicable
Acknowledgments: The authors would like to thank the editor in chief, the associate editors, and the reviewers for their worthy suggestions, time, and effort in improving and finalizing this paper. The authors would further like to thank Sajid A. Khan (sajidkhan@cuilahore.edu.pk) for providing the support and resources.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix C. Grid Side Filter Model
The small-signal model of LC output filter and coupling inductance by assuming that voltage provided by the inverter is the same as the demand voltage is expressed as: . .

Appendix D. Complete Model of an ith Inverter
The output currents are the output variables of i th inverters which can be expressed in vector form ∆iodqi. ∆ioDQi, which is a small-signal output current, expressed as:

Appendix E. Combined Model of N Inverters
In MG, the N number of DG inverters may be operating as a source at variable distances among each other. In this section, our approach is to discuss the possible sub-models form of all individual i th to k th DG inverters and combine them to the existing corresponding network. The combined small-signal model of "N" number of DG inverters can be represented as:

Appendix F. Network and Load Model
If an i th feeder line is connected between node j and k, then algebraic equations can be presented as the following: