HIL-Assessed Fast and Accurate Single-Phase Power Calculation Algorithm for Voltage Source Inverters Supplying to High Total Demand Distortion Nonlinear Loads

: The dynamic performance of the local control of single-phase voltage source inverters (VSIs) can be degraded when supplying to nonlinear loads (NLLs) in microgrids. When this control is based on the droop principles, a proper calculation of the active and reactive averaged powers (P–Q) is essential for a proﬁcient dynamic response against abrupt NLL changes. In this work, a VSI supplying to an NLL was studied, focusing the attention on the P–Q calculation stage. This stage ﬁrst generated the direct and in-quadrature signals from the measured load current through a second-order generalized integrator (SOGI). Then, the instantaneous power quantities were obtained by multiplying each ﬁltered current by the output voltage, and ﬁltered later by utilizing a SOGI to acquire the averaged P–Q parameters. The proposed algorithm was compared with previous proposals, while keeping the active power steady-state ripple constant, which resulted in a faster calculation of the averaged active power. In this case, the steady-state averaged reactive power presented less ripple than the best proposal to which it was compared. When reducing the velocity of the proposed algorithm for the active power, it also showed a reduction in its steady-state ripple. Simulations, hardware-in-the-loop, and experimental tests were carried out to verify the e ﬀ ectiveness of the proposal.


Introduction
The Smart Grid, as a concept, can be defined in terms of its outcomes as an electrical system that operates in an efficient manner, that provides a reliable energy supply and a power quality for the needs of a digital economy, that demonstrates a resilient performance against uncertainties or grid faults, and that integrates a large variety of distributed energy resources (DERs) in the conventional electrical grid, especially renewable energy sources (RESs), according to the Department of Energy of the United States [1] and the European Commission [2]. In addition, following the United Nations Sustainable Development Goals [3], the energy supply must be affordable, reliable, sustainable, and accessible for all users, especially when involving DERs based on RESs, as can be seen in References [4,5]. These DERs can be composed of a mixture of energy production units, energy storage systems (ESSs), and loads that operate jointly in clusters that are connected or not to the main electrical grid infrastructure [6]. This set The calculation of the P-Q parameters is fundamental for these droop-operated VSIs. When functioning in islanded mode, single-phase MGs are weaker in terms of stability, as compared with three-phase VSI-based systems, especially when sharing NLLs. In these conditions, accurate and fast calculations of P-Q become crucial [34,35]. Different solutions for the calculation of single-phase P-Q can be found in the literature. It is usual to see that, for the P-Q calculations, a voltage signal delayed 90 • to the measured one in the point of common coupling (PCC) is employed. After that, the instantaneous active and reactive powers are calculated and later conditioned to obtain the P-Q parameters. For the in-quadrature voltage, one method consists of applying a transport delay by digital means without filtering the signal in magnitude [36,37]. Other approaches are based on the p-q theory, which implies a dq-Synchronous Reference Frame technique [38,39]. In this case, it is also acquired the in-quadrature signal of the PCC current. Moreover, second-order generalized integrators (SOGIs) were used thanks to their low-pass filter capability (SOGI-LPF), as in References [40,41]. Coming up next, the calculations of the instantaneous power quantities signals are conditioned for the extraction of the averaged active power (P), and the averaged reactive power (Q). One method consists of the application of SOGI-LPF [42]. Thus, a better strategy for improving the transient response against abrupt load changes is performed by adding a final stage for removing the double frequency power components, either with a final low-pass filter (LPF) stage in Reference [41] or without final LPF [35]. In the same conditions, in Reference [43] a different approach based on the discrete Fourier Transformation was utilized for directly obtaining P and Q. All of these strategies have in common that were only designed for droop-operated single-phase VSI sharing linear loads, introducing a significant delay that constraints the velocity of the response in front of abrupt load changes. The common goal of all those techniques is to obtain accurate and fast P-Q calculations to enhance the robustness of the parallelization of VSIs when sharing loads. Then, in Reference [35], it was proposed a method that introduced a pre-filtering of the measured current before the calculation of the instantaneous powers that resulted in being faster in front of NLL. This pre-filtering was achieved by the band-pass filter (BPF) capability of a double SOGI approach (DSOGI), followed by the remotion of the double frequency power components like Reference [34], but without the LPF final stage. The fastest solution from those benchmarked was [35], but at the cost of increased complexity in the calculation scheme and a worsening of the reactive averaged power calculation. The algorithm was designed based on the high distortion of the NLL, increasing the signal conditioning blocks with respect to previous proposals. Thus, the NLL was presented as a highly polluted in harmonics current, but was not quantitatively characterized.
The algorithm proposed in this work is designed by firstly characterizing the NLL accordingly to the well-known standard IEEE std 519-2014 [44]. Consequently, the topology of the calculation structure is modified, and the conditioning signal blocks are reduced with respect to Reference [35]. The objective is to obtain a faster and more accurate calculation of the P-Q parameters. Figure 1 represents a basic scheme of a single-phase MG topology containing two VSI sharing an NLL. The block "P-Q Power Calculation" is the research object of this paper. The voltage v o (t) and the current i o (t) are measured at the PCC when the switches S 0 and S 1 are open. The averaged active and reactive power, P and Q, are respectively calculated in the "P-Q Power Calculation" block from these measurements. Then a voltage reference v ref (t) is generated. Later, the inner control loops based on [45] employ this reference and finally generate a pulse-width modulation (PWM) for the switching of the H-Bridge.

Description of the System under Test: VSI Supplying to a Nonlinear Load
The NLL is a load which voltage-current characteristic is not linear. Different types of NLLs are found in the literature, which may be classified as attending different criteria. They are either characterized only considering the V-I characteristic [46], or based on power quality parameters associated with the active, reactive, and distorted power [47] or on the measured current [48]. In References [49,50], the classification is focused on the harmonic components only present in the Electronics 2020, 9, 1643 4 of 24 current through the NLL. Similarly, Reference [51] studies an NLL accordingly to its current harmonic components and its total harmonic distortion in the current. The NLL study in Reference [52] is based on the frequency spectrum and the total harmonic distortion in voltage. The particular NLL load of this work (Z NL ) consists of an unbalanced diode-bridge rectifier (DBR) that supplies power to a resistive-capacitive (R-C) load, shown in Figure 2. The R-C load parameters are listed in Table 1 and are characterized accordingly to power quality parameters listed in Table 2 power to a resistive-capacitive (R-C) load, shown in Figure 2. The R-C load parameters are listed in Table 1 and are characterized accordingly to power quality parameters listed in Table 2 Figure 2. Nonlinear load ZNL, consisting of an unbalanced diode-bridge rectifier that supplies to an R-C load, based on Reference [35].
After reaching a steady-state, the switch S2 closes at t = ton, triggering an abrupt change in the value of P and Q. This load is the same as the one employed in Reference [35], and the voltage measured in the PCC follows Expression (1). The current in the PCC can be described by (2), and draws a highly distorted waveform: Electronics 2020, 9, x FOR PEER REVIEW 4 of 24 power to a resistive-capacitive (R-C) load, shown in Figure 2. The R-C load parameters are listed in Table 1 and are characterized accordingly to power quality parameters listed in Table 2 Figure 2. Nonlinear load ZNL, consisting of an unbalanced diode-bridge rectifier that supplies to an R-C load, based on Reference [35].
After reaching a steady-state, the switch S2 closes at t = ton, triggering an abrupt change in the value of P and Q. This load is the same as the one employed in Reference [35], and the voltage measured in the PCC follows Expression (1). The current in the PCC can be described by (2), and draws a highly distorted waveform:

of 24
After reaching a steady-state, the switch S2 closes at t = t on , triggering an abrupt change in the value of P and Q.
This load is the same as the one employed in Reference [35], and the voltage measured in the PCC follows Expression (1). The current in the PCC can be described by (2), and draws a highly distorted waveform: v where V o and I o are the voltage and current amplitudes, respectively; ω o is the fundamental frequency of the system (100π rad/s); ϕ o is the phase-shift between the fundamental components of the voltage and the current; h is the harmonic index; I h is the amplitude of the harmonic components of the intensity; and the phase-shift ϕ h corresponds to each current harmonic component. The term I DC corresponds to a DC offset present in the load current. Figure 3a detailed the local control structure for a single-phase VSI supplying the NLL of Figure 2. This model is based on References [45,53], and this work only focuses on the block named "POWER CALCULATION BLOCK," that generates the voltage reference for the voltage and current inner loops. Those consist of a proportional-integral (PI) block for the voltage and a proportional-resonant (PR) loop for the current. A final stage for the PWM for the switching of an H-bridge is shown. where Vo and Io are the voltage and current amplitudes, respectively; is the fundamental frequency of the system (100π rad/s); is the phase-shift between the fundamental components of the voltage and the current; h is the harmonic index; Ih is the amplitude of the harmonic components of the intensity; and the phase-shift corresponds to each current harmonic component. The term corresponds to a DC offset present in the load current. Figure 3a detailed the local control structure for a single-phase VSI supplying the NLL of Figure  2. This model is based on References [45,53], and this work only focuses on the block named "POWER CALCULATION BLOCK," that generates the voltage reference for the voltage and current inner loops. Those consist of a proportional-integral (PI) block for the voltage and a proportional-resonant (PR) loop for the current. A final stage for the PWM for the switching of an H-bridge is shown. Figure 3b displays the harmonic distribution of the measured current. Regarding the recommended limitations in [44], the individual harmonic distortion in the current (THDi) and the total distortion demand (TDD) are out of limits. Moreover, a DC component in (2) and pictured in Figure 3b is not allowed by Reference [44]. Finally, the individual harmonic distortion in voltage, THDv, and its total harmonic distortion, THD, are within limits.   (b) nonlinear current harmonic distribution for Z NL without filtering, after S 0 is turned on.
Electronics 2020, 9, 1643 6 of 24 Figure 3b displays the harmonic distribution of the measured current. Regarding the recommended limitations in [44], the individual harmonic distortion in the current (THDi) and the total distortion demand (TDD) are out of limits. Moreover, a DC component in (2) and pictured in Figure 3b is not allowed by Reference [44]. Finally, the individual harmonic distortion in voltage, THDv, and its total harmonic distortion, THD, are within limits. Table 2 summarizes this characterization of Z NL :

P-Q Calculation Algorithms
The droop equations for a mainly inductive system are as follows: where ω * is the calculated angular frequency of the system, V * the calculated voltage amplitude, V n is the rated value for the voltage, ω n is the frequency rated value, and m and n are the droop coefficients. With these parameters, it is generated the sinusoidal voltage reference v re f (t) necessary for the inverter inner control loops of Figure 3a: The calculation of P-Q is done as follows, in the time domain: where (6) is the calculated active instantaneous power [38], P is the averaged active power in (7), and a double frequency pulsating component, p(t), is represented by (8). For the instantaneous reactive power calculation, an in-quadrature signal for the voltage is employed, represented by (9). Then an instantaneous reactive power quantity is calculated: where q(t) is the instantaneous reactive power, (10), Q is the averaged reactive power in (11), and q(t) is a double frequency component. Only in Reference [35], the highly distorted load current was pre-filtered prior to the instantaneous power calculations in (6) and (10), using a DSOGI approach and its BPF capability. Thus, the LPF of a SOGI was employed as a quadrature signal generator (QSG) delaying in π 2 rad the voltage signal. The BPF and LPF transfer functions for a SOGI are described in (13) and (14), respectively: Electronics 2020, 9, 1643 7 of 24 H q (s) = 2ξω 2 s 2 + 2ξωs + ω 2 (14) where ξ is the damping factor, and ω is the center frequency of the system. Figure 4 shows the structure of a SOGI, with its BPF and LPF magnitude and phase Bode plots.
Electronics 2020, 9, x FOR PEER REVIEW 7 of 24 In Figure 4b, the SOGI-BPF is more selective while reducing the damping factor. The same occurs with its LPF capability. The attenuation presents a rate of −30dB/decade for the BPF and −60dB/decade for the LPF. However, the LPF transfer function presents a drawback in attenuation for frequencies below the fundamental.
The following Figure 5 shows three averaged active and reactive power calculation: The obtained averaged active powers are named as PavC, PavT, and PavM, for Figure 5a-c, respectively. In the same manner, the averaged reactive powers are named as QavC, QacT, and QavM.   . Structure of second-order generalized integrator (SOGI) and bode diagrams with various damping factors: (a) structure of a SOGI with its damping factor, ξ, and center frequency, ω, where V in is the input signal, V d the direct output signal, filtered by a band-pass filter (BPF); V q the in-quadrature output signal, filtered by a low-pass filter (LPF). (b) Magnitude and phase Bode plots for BPF in (13), varying ξ from 0.1 to 0.9. (c) The magnitude and phase Bode plots for LPF in (14) varying ξ from 0.1 to 0.9.
In Figure 4b, the SOGI-BPF is more selective while reducing the damping factor. The same occurs with its LPF capability. The attenuation presents a rate of −30 dB/decade for the BPF and −60 dB/decade Electronics 2020, 9, 1643 8 of 24 for the LPF. However, the LPF transfer function presents a drawback in attenuation for frequencies below the fundamental.
The following Figure 5 shows three averaged active and reactive power calculation: In Figure 5a, Pavc and Qavc are obtained after an LPF stage. On Figure 5b, the oscillatory double frequency components are extracted by a SOGI-BPF tuned at 2 and then removed from the calculated instantaneous powers. Lastly, an LPF similar to Figure 5a is employed to reduce the steady-state ripple. The latest algorithm in Figure 5c also removes the double frequency components but eludes the final LPF stage.
Besides, Figure 5 shows the increasing complexity of the algorithms, especially Figure 5c, when supplying to the NLL in Reference [35]. Therefore, the next section proposes a novel algorithm that calculates in a faster manner the averaged active power and, in a more accurate manner, the averaged reactive power, while reducing the complexity of the calculation scheme.

Power products
Filtering Filtering io(t)  (b) advanced P-Q calculation based on the schemes and algorithms reported in Reference [41]; and (c) P-Q calculation based on Reference [35].
The obtained averaged active powers are named as PavC, PavT, and PavM, for Figure 5a-c, respectively. In the same manner, the averaged reactive powers are named as QavC, QacT, and QavM.
In Figure 5a, v oqC (t) is the π 2 delayed voltage signal, as described in (9), achieved by using a time-delay block [47]. In Figure 5b,c, a SOGI-LPF was employed as QSG for the obtention of v oqT (t) and v oqM (t), as well as BPF for the direct component of the voltage.
The next common stage of the three calculation algorithms consists of obtaining the instantaneous active and reactive powers, as described in (6) and (10). Note that the measured current is directly applied, except in Figure 5c, where it utilizes a BPF for the current [35].
In Figure 5a, Pavc and Qavc are obtained after an LPF stage. On Figure 5b, the oscillatory double frequency components are extracted by a SOGI-BPF tuned at 2ω 0 and then removed from the calculated instantaneous powers. Lastly, an LPF similar to Figure 5a is employed to reduce the steady-state ripple. The latest algorithm in Figure 5c also removes the double frequency components but eludes the final LPF stage. Besides, Figure 5 shows the increasing complexity of the algorithms, especially Figure 5c, when supplying to the NLL in Reference [35]. Therefore, the next section proposes a novel algorithm that calculates in a faster manner the averaged active power and, in a more accurate manner, the averaged reactive power, while reducing the complexity of the calculation scheme.

Proposed P-Q Calculation Algorithm
The proposed algorithm for the calculation of P-Q is presented in Figure 6. This new algorithm is intended to ease the dynamic performance of the system by reducing the settling-time during abrupt load changes.
The proposed algorithm for the calculation of P-Q is presented in Figure 6. This new algorithm is intended to ease the dynamic performance of the system by reducing the settling-time during abrupt load changes. Figure 6a shows that the measured current, ( ), is conditioned through SOGI-0, obtaining a BPF filtered component, ( ) , and an LPF filtered component, ( ) . Then, each current component is directly multiplied by the voltage, computing the instantaneous active power, ( ), and the instantaneous reactive power, ( ). Those quantities are expected to contain the averaged active and reactive powers Pav and Qav, respectively, plus a certain amount of undesired harmonic components, similarly to Reference [35] and following Equations (6) and (10). Therefore, a final LPF stage is applied to each instantaneous power to acquire active and reactive quantities with the lesser possible steady-state harmonics. Those last active and reactive power quantities are, respectively, ( ) and ( ) , and contain the desired Pav and Qav. Comparing this scheme to those pictured in Figure 5, note that the voltage signal is not conditioned due to the specific NLL characterized according to Table 2. The proposed scheme also presents fewer signal conditioners and control parameters than in Reference [35], in Figure 5c, showing a simplified calculation structure.    Figure 6a shows that the measured current, i o (t), is conditioned through SOGI-0, obtaining a BPF filtered component, i od (t), and an LPF filtered component, i oq (t). Then, each current component is directly multiplied by the voltage, computing the instantaneous active power, p i (t), and the instantaneous reactive power, q i (t). Those quantities are expected to contain the averaged active and reactive powers Pav and Qav, respectively, plus a certain amount of undesired harmonic components, similarly to Reference [35] and following Equations (6) and (10). Therefore, a final LPF stage is applied to each instantaneous power to acquire active and reactive quantities with the lesser possible steady-state harmonics. Those last active and reactive power quantities are, respectively, p F (t) and q F (t), and contain the desired Pav and Qav. Comparing this scheme to those pictured in Figure 5, note that the voltage signal is not conditioned due to the specific NLL characterized according to Table 2. The proposed scheme also presents fewer signal conditioners and control parameters than in Reference [35], in Figure 5c, showing a simplified calculation structure. Figure 6b is the frequency-domain analytical representation of the calculated quantities and the transfer functions for each SOGI employed in Figure 6a. The signals reported in this scheme are the Laplace Transform of those indicated in Figure 6a: Analytically, the LPF and BPF transfer functions of SOGI0 are represented, respectively, by H d0 (s) and H q0 (s): Those are tuned at ω 0 , and their selectivities are controlled through the damping factor, ξ i . Hence, the direct and the in-quadrature filtered currents, i od (s) and i oq (s), respectively, are as follows: and the time-domain expressions are: where h 1 is a harmonic index. For (19), I dh is the harmonic amplitude, and ϕ dh is its phase-shift. In (20), I qh is the harmonic amplitude and ϕ qh , is its phase-shift. Note that I dh I qh and ϕ dh ϕ qh . Later, following the scheme in Figure 6a, the instantaneous active and reactive powers, p i (t) and q i (t), are as follows: Pav and Qav are the averaged active and reactive power outputs, respectively. P i (s) and Q i are the domain frequency expression of (21) and (22), related in Figure 6b.
Therefore, for SOGI1 and SOGI2, they utilize their LPF capability as follows: where H P (s) and H Q (s) are the transfer functions for the LPF capability of SOGI1 and SOGI2, respectively. Note that both transfer functions are essentially the same, only differentiated by the h 1 and h 2 coefficients. Those coefficients (h 1 , h 2 ) [0.05, 0.5] are employed for the attenuation of subharmonics reducing the SOGI1 and SOGI2 LPF cutoff frequencies. Next, Figure 7 shows their magnitude and phase plots. For the sake of simplicity, h 1 = h 2 = h i : Finally, back to the time domain, the following can be found: The expressions ( ) and ( ) contain Pav and Qav. The parameter k is a harmonic index similar to h in (21), being ( ) and ( ) the undesired oscillatory components. For the attenuation of those, SOGI1 is tuned at a frequency = ℎ and SOGI2 at = ℎ . Low values for h1 and h2 will lead to a substantial reduction of these components. The next section includes simulations to study the values of , , h1, and h2 for a more accurate and faster calculation of the P-Q parameters.

Simulation Results
The proposed algorithm is simulated to compare its dynamic performance against the structures shown in Figure 5, with an abrupt load change after closing S1. A similar steady-state ripple for the active power calculation is set as a reference for the analysis. The parameters of the simulations are listed in Table 3. Then, P F (s) and Q F (s) are the result of the filtering of (23) and (24) by their respective LPF-SOGI: Finally, back to the time domain, the following can be found: The expressions p F (t) and q F (t) contain Pav and Qav. The parameter k is a harmonic index similar to h in (21), being p k (t) and q k (t) the undesired oscillatory components. For the attenuation of those, SOGI1 is tuned at a frequency ω = h 1 ω 0 and SOGI2 at ω = h 2 ω 0 . Low values for h 1 and h 2 will lead to a substantial reduction of these components. The next section includes simulations to study the values of ξ P , ξ i , h 1 , and h 2 for a more accurate and faster calculation of the P-Q parameters.

Simulation Results
The proposed algorithm is simulated to compare its dynamic performance against the structures shown in Figure 5, with an abrupt load change after closing S 1 . A similar steady-state ripple for the active power calculation is set as a reference for the analysis. The parameters of the simulations are listed in Table 3.
In Reference [35], it was demonstrated that the more suitable algorithm in the presence of an NLL was the calculation of PavM. The following figure shows a family of Pav plots after varying its control parameters and comparing it against PavM. Figure 8a shows a family of P av plots varying, 0.1 ≤ ξ i ≤ 0.7075, while keeping constant h 1 = 0.25 and ξ p = 0.7075. It can be seen that, when increasing the damping factor, the transient response is faster. However, there is an undesired overshoot when ξ i > 0.2. Therefore, the fastest configuration avoiding the overshoot is the one with ξ i = 0.2.   Figure 9b is the detail of the steady-state ripple, where it can be seen that the fastest algorithm, Pav, shows a similar ripple to the other algorithms. In these conditions, the calculation of Pav results faster than PavM. Therefore, these parameters were chosen for the study of the employed for this last simulation.    Figure 9a shows the calculation of PavC, PavT, and PavM, compared with the proposed P av with ξ i = 0.2. and h 1 = 0.25. Figure 9b is the detail of the steady-state ripple, where it can be seen that the fastest algorithm, Pav, shows a similar ripple to the other algorithms. In these conditions, the calculation of Pav results faster than PavM. Therefore, these parameters were chosen for the study of the employed for this last simulation. In Figure 10, the spectrum and the steady-state THDs of Pav and PavM are compared. From the observation of Figures 9 and 10, it is deduced that, for 0.15 ℎ1 < 0.25, there is a family of Pav calculations that result to be faster and more accurate than PavM.
The study of the reactive power quantities is done with h2 = 0.1, = 0.7075 and with = 0.2. Figure 11 shows the calculation of Qav, QavM, QavT, and QavC when the abrupt load change occurs. In Figure 10, the spectrum and the steady-state THDs of Pav and PavM are compared. From the observation of Figures 9 and 10, it is deduced that, for 0.15 ≤ h 1 < 0.25, there is a family of Pav calculations that result to be faster and more accurate than PavM.  The study of the reactive power quantities is done with h 2 = 0.1, ξ p = 0.7075 and with ξ i = 0.2. Figure 11 shows the calculation of Qav, QavM, QavT, and QavC when the abrupt load change occurs. Figure 11a shows that QavM is the worst option in terms of steady-state ripple when supplying the NLL. Moreover, Qav presents a similar time response than QavM and the lowest steady-state ripple (see Figure 11b). Thus, similarly to the active power calculation, if the transient time needs to be reduced, that will succeed at the cost of a higher ripple. Nevertheless, although the reactive power calculation is necessary for the droop control, a variation of less than 10 VAr barely influences the droop reference generation due to the nature of the NLL. Considering the calculated mean value of Pav and Qav, it can be extracted a power factor (PF) equals to 0.9976. If the value of Qav increases up to 35 VAr, then PF = 0.9928. Although the load has not been characterized according to the PF, it indicates that the reactive power variations in mean value barely influence its value.
(c)  Then, it is deduced that the proposed calculation method for Pav can be faster than PavM for a range of values of ξ i and h 1 , keeping ξ p = 0.7075. Moreover, when increasing the SOGI1-LPF capability, it maintains a better or similar settling-time while reducing the steady-state ripple. Thus, when Pav reaches a similar settling-time of that of PavM by reducing the h 1 parameter, the THD falls from 1.32% to 0.59% ( Figure 10). The commented results are shown in Table 4: Table 4. THD, settling-time, and time-delay for the simulation of Pav. As it can be seen in Table 4, Pav settling-time is 37.5% shorter than PavM when ξ i = 0.2 and h 1 = 0.25 with a similar ripple (Figure 11c). The time-delay for both calculations is almost the same, 38 ms for P avM and 40 ms for P av . Then, when h 1 = 0.15, the THD falls drastically down to 0.59% while keeping an 18% shorter settling-time. However, in this last case, the time-delay is higher in 20% for the Pav calculation concerning PavM. For this final reason, the chosen set of parameters for comparing Pav against PavM is are

Calculated THD Settling-Time (ms)/% Reduction with Respect to PAVM
On the other hand, the simulation results for the reactive power calculation algorithms are compared in Table 5, in terms of settling-time and comparing its steady-state ripple through a THD, with respect to the DC component analysis.  Table 5 shows how the conventional droop method, QavC, is the best option for reducing steady-state ripple in reactive power. However, its settling-time is the worst with a value of 780 ms. On the other hand, the best time is achieved by the QavM algorithm but at the cost of a higher THD = 7.85%. The conclusion is that the calculation of Q is more accurate through the proposed algorithm, with similar transient velocity. Therefore, the chosen set of parameters for the proposed algorithm is as follows: The next section pretends to assess the obtained results from the carried out simulations.

Hardware in the Loop Assessment
For the assessment of the proposed algorithm, HIL tests are carried out. Those tests compare the calculation of PavM and QavM against Pav and Qav with the chosen parameters from the simulations.
For this purpose, a real-time interface platform based on dSPACE 1006© (dSPACE Inc.50131 Pontiac Trail Wixom, MI, USA 48393-2020) digital platform is operated. The control structure presented in Figure 2 is firstly discretized in Matlab/Simulink/SimPowerSystems© (The MathWorks Inc., Natick, MA, USA) and then compiled in C code for its download in the dSPACE. Moreover, this RTI platform supports the model libraries of physical/electrical plants from Matlab/Simulink/SimPowerSystems© (The MathWorks Inc., Natick, MA, USA). Those libraries correspond to the modeled H-bridge, the LCL filter, and the NLL under test. The electronic central unit (ECU) of the dSPACE compiles the control algorithms on its multiprocessor core. The control-desk software permits the configuration and control of the tests, acting as a human-machine interface. The switching frequency of 10 kHz for the VSI is emulated, setting the sample time at Ts = 100 µs. Note that the discretization of the integrators employed in the SOGI has been achieved through a third-order method: A first HIL Test-1 is then carried out to compare the proposed algorithm against PavM, PavT, and PavC. The parameters for Pav are listed in Table 6.  Figure 12 shows the load current plots; the active and reactive power after an abrupt load change is done manually.

Test
; HIL Test-1 0.2 0.7075 0.25;0.10 HIL Test-2 0.2 0.7075 0.15;0.10 Figure 12 shows the load current plots; the active and reactive power after an abrupt load change is done manually.  Then, in Figure 12b, Pav is compared against PavM, PavT, and PavC. The shortest transient response corresponds to Pav, while the steady-state ripple is kept constant. Figure 12c shows the transient response for the calculated reactive powers. There, it can be  Then, in Figure 12b, Pav is compared against PavM, PavT, and PavC. The shortest transient response corresponds to Pav, while the steady-state ripple is kept constant. Figure 12c shows the transient response for the calculated reactive powers. There, it can be appreciated that the response rapidity is similar between Qav and QavM. Moreover, the lowest steady-state ripple corresponds to Qav. The relevant comparison here is between the proposed algorithm and that based on Reference [35]. The results are exposed in Table 7. The settling-time for Pav is shorter in 35.7% with respect to PavM, similar to the simulation results with the same control parameters. The time-delay was found to be similar for Pav and PavM. However, in the reactive power, this time-delay is larger for Qav than for QavM, even when the Qav settling-time is a 7.1% minor than QavM.
A second HIL Test-2 is carried out, reducing the h 1 to 0.15 for Pav, to assess the simulation tests achieved in this sense. Figure 13 shows the active power responses during an abrupt load change.
Electronics 2020, 9, x FOR PEER REVIEW 19 of 24 Figure 13. HIL active averaged power calculation, when an abrupt load change occurs, through HIL emulation, using a dSPACE-RTI setup at Aalborg Microgrid Laboratory: active power calculation, Pav (green), PavM (red), and its steady-state ripple detail.

Experimental Results
An experimental test was carried out, employing the load described in Table 2. The experimental setup is prepared to evaluate the model simulated and HIL-assessed for the active power calculation in the presence of a measured TDD = 124.9% in current at the PCC. The experimental setup is shown in Figure 14. It is composed of a VSI Danfoss © FC302, 2.2 kW rated, interfaced to a real-time dSPACE 1006 platform, for the switching signals for the H-bridge and the measured parameters. The current of the NLL was monitored by using a Fluke 435-II Power Quality and Energy Analyzer. The power calculation algorithms tested correspond to PavM and Pav. The results obtained are displayed in Figure 15 and Table 8.  Figure 13 shows that a low h 1 coefficient in Pav allows a smaller steady-state ripple while maintaining a similar settling-time both for active and reactive powers, Figures 13 and 12c, respectively. However, the time-delay for Pav is 50 ms, superior to the reported 40 ms in Table 7 for PavM.

Experimental Results
An experimental test was carried out, employing the load described in Table 2. The experimental setup is prepared to evaluate the model simulated and HIL-assessed for the active power calculation in the presence of a measured TDD = 124.9% in current at the PCC. The experimental setup is shown in Figure 14. It is composed of a VSI Danfoss© FC302, 2.2 kW rated, interfaced to a real-time dSPACE 1006 platform, for the switching signals for the H-bridge and the measured parameters. The current of the NLL was monitored by using a Fluke 435-II Power Quality and Energy Analyzer. The power calculation algorithms tested correspond to PavM and Pav. The results obtained are displayed in Figure 15 and Table 8.
An experimental test was carried out, employing the load described in Table 2. The experimental setup is prepared to evaluate the model simulated and HIL-assessed for the active power calculation in the presence of a measured TDD = 124.9% in current at the PCC. The experimental setup is shown in Figure 14. It is composed of a VSI Danfoss © FC302, 2.2 kW rated, interfaced to a real-time dSPACE 1006 platform, for the switching signals for the H-bridge and the measured parameters. The current of the NLL was monitored by using a Fluke 435-II Power Quality and Energy Analyzer. The power calculation algorithms tested correspond to PavM and Pav. The results obtained are displayed in Figure 15 and Table 8. First of all, the load draws an asymmetrical current, as shown in Figure 15a, with +2.19A/−1.46A peak values. The difference between the measured and the simulated and HIL tested current is due to power losses in the whole system. Figure 15b shows the harmonic spectrum measured by employing a Fluke 435-II Power Quality and Energy Analyzer, which yields a TDD = 124.9%, compatible with the simulated in Figure 3b.  First of all, the load draws an asymmetrical current, as shown in Figure 15a, with +2.19A/−1.46A peak values. The difference between the measured and the simulated and HIL tested current is due to power losses in the whole system. Figure 15b shows the harmonic spectrum measured by employing a Fluke 435-II Power Quality and Energy Analyzer, which yields a TDD = 124.9%, compatible with the simulated in Figure 3b. Figure 15c presents the PavM calculation, with the detail of the steady-state ripple. In the same manner, Figure 15d shows the proposed Pav calculation, with the detail of its steady-state ripple. Note that the transient for Pav results faster than for PavM, as expected. The measured settling-times are exposed in Table 8, along with the simulation and HIL results for active power.
The Simulation, the HIL test and experimental results for the settling-time in active power calculation are resumed in Table 8. The results compared here correspond to the Pav and PavM keeping the same steady-state ripple. The proposed Pav results to be 37.5% faster than PavM in the Simulation, 35.7% in the HIL tests, and 30% in the experimental test. The settling-time measured in the experimental test is a 40% higher than the obtained in the simulation, for Pav. In contrast, the measured settling-time for the HIL test is 16.67% higher than the simulated for Pav. Those differences are attributed to the 3rd order integrator employed both in HIL and the experimental setup. For the experimental setup, it may be considered latencies due to the internal communications and data acquisition boards, as well as nominal values biases of components and energy losses. Electronics 2020, 9, x FOR PEER REVIEW 20 of 24 Figure 15c presents the PavM calculation, with the detail of the steady-state ripple. In the same manner, Figure 15d shows the proposed Pav calculation, with the detail of its steady-state ripple. Note that the transient for Pav results faster than for PavM, as expected. The measured settling-times are exposed in Table 8, along with the simulation and HIL results for active power.

Conclusions
The proposed method enhances the dynamical performance in terms of rapidity and accuracy of the droop-based local control, which degrades in the presence of NLLs like the employed in this work, which was characterized considering IEEE std 519-2014. Only PavM demonstrated its suitability in the presence of NLL in Reference [35], focusing the calculation effort on the obtention of the fundamental component of the current and avoiding a final LPF stage. However, this previous work did not differentiate types of NLL. Oppositely, the proposed method characterizes the NLL, and then the algorithm architecture is decided. Thus, the implemented algorithm results less complex than those compared with, when supplying an NLL.
The proposed algorithm was compared with the previously studied algorithms and assessed through Matlab/Simulink simulation and a HIL test. Finally, an experimental test for the active power Pav and PavM was carried out to evaluate the proposed model. The main conclusions after analyzing the simulation, HIL, and experimental results are summarized as follows, in terms of transient response velocity and steady-state accuracy: Velocity: • Reduction of settling-time in 30% for the calculation of Pav with respect to Pav M [35] while keeping a similar steady-state ripple (Table 8). Accuracy: • Active Power: Reduction in 47.78% of the steady-state calculated THD with respect to DC in the simulations for Pav, when the settling-time is similar (Table 4).

•
Reactive Power: Reduction in 68.66% of the steady-state calculated THD with respect to DC in the simulations for Qav, when the settling-time is similar ( Table 5).
As expected, the settling-time for the Pav calculation during an abrupt load change resulted in being smaller than the other compared methods. Moreover, from Table 8, it can be seen how the relative reduction of settling-time was preserved in all the scenarios, i.e., between 30% and 37.5%. That leads to a faster operation in the droop controlled VSI in the presence of high TDD NLLs, which points to an increase in single-phase MG stability when sharing NLLs.
Concerning the enhancement of the accuracy, it is noteworthy that the calculation of Qav results to be more accurate than the other methods, see Table 5. Regarding the active power calculation, the steady-state ripple can be smoother by reducing the settling-time for Pav, see Figures 8b and 13.
Future investigations are intended to study the same issues when other NLL types are present, considering a well-known standard as IEEE std 519-2014. Further studies are also planned for the parallelization of single-phase VSI against different types of NLL to study its dynamic performance and control stability.