A Photovoltaic Power Curtailment Method for Operation on Both Sides of the Power-Voltage Curve

Abstract: Massive integration of non-dispatchable energy into electric power systems is a challenging task. Electric power systems are becoming increasingly vulnerable in terms of frequency stability, as renewable energy displaces conventional synchronous generation from the energy mix. For this reason, grid codes are starting to demand different ancillary services from renewable generators, such as frequency control. In contrast to wind generators, which can deliver to the grid part of the kinetic energy stored in their rotating mass, photovoltaic generators must provide this service using batteries or power curtailment methods. The latter approach is preferable regarding the initial investment and its implementation cost, and several methods have been presented in the literature for this purpose. However, there is no consensus in which is the most appropriate side for operating the photovoltaic system in the curtailed mode. As both possible options have advantages and drawbacks, this paper proposes a novel photovoltaic power curtailment strategy that allows operation on both sides of the power-voltage curve depending on the needs. Moreover, in order to estimate the output characteristic of the photovoltaic system, a real-time nonlinear least squares curve fitting is applied. The proposed methodology has been tested in a simulation environment and the results show that this strategy achieves the requested active power reserves, regardless of the operation side.


Introduction
Global warming is a real problem. The emission of greenhouse gases into the atmosphere, mainly derived from industrial activity, has caused the Earth's temperature to rise 1 • C above the pre-industrial level [1]. In this context, the design and control of future electric power systems play a fundamental role. As a matter of fact, several authors have discussed the feasibility and viability of reaching 100% renewable-electricity systems [2,3], as well as energy policies that minimize dependence on fossil fuels [4]. It is well-known that Renewable Energy Sources (RES) are displacing conventional power plants from the generation mix [5]. However, in contrast to conventional generation, RES are considered non-dispatchable sources of energy, as they depend on meteorological conditions to produce their output power. This means that RES power forecasting is volatile and, to a certain extent, unpredictable [6]. Apart from the energy allocation problem, the displacement of synchronous generators by renewable ones involves the deterioration of some ancillary services such as frequency control. RES used to be interfaced with the utility grid via power electronics, which have no rotational parts and therefore do not contribute to the total inertia of the power system. Inertia is a fundamental parameter that directly affects the frequency nadir and the rate of change of frequency following a power unbalance between generation and demand, both key issues in the operation of electric power Let us suppose that the operating point of the PV module is point 1, and therefore the generated power is equal to P MPP . In these conditions, if the system operator imposes a power reserve level of 50%, there are two possible operating points with P RES -point 2 and point 3. While operating at point 2, if there is a sudden decrease in irradiance, the operation point falls to 4, which is a stable point. However, when there is a sudden drop of irradiance, and if point 3 is chosen to operate at P RES , the operating point falls to 5, with no power production. Thus, in some circumstances, it is preferable to operate on the right part of the P-V curve, and in some others, on the left. The method presented in this paper allows left and right side operation. Section 2 describes the PV model used. Section 3 analyzes the curve fitting problem. Section 4 details the design of the voltage control system. Section 5 presents simulation results showing the method performance. Finally, Section 6 draws conclusions from the work.

Modeling of Photovoltaic Cell, Module and Array
In this work, the Single-Diode Model (SDM) of the photovoltaic cell is adopted. Figure 2 presents the equivalent circuit of the SDM where I ph,cell is the photocurrent, I s,cell is the saturation current of the diode, n is the diode ideality factor, and R s,cell and R sh,cell are the series and shunt resistance, respectively [17]. The electrical variables of interest are the PV cell terminal voltage (V pv,cell ) and current (I pv,cell ). Applying Kirchhoff's first law, it is possible to relate both magnitudes through the implicit Equation (1) where V T,cell is the thermal voltage of the diode, calculated as (2): where k is the Boltzmann constant, q is the electron charge and T cell is the cell temperature expressed in Kelvin. To obtain the PV module model, it'd be necessary to take into account the number of series-connected cells (n cell ) in the PV module. Both resistances and the thermal voltage are multiplied by this value: R s,module = n cell R s,cell R sh,module = n cell R sh,cell .
The PV array parameters can be calculated by introducing the number of series-connected modules (n s ) and parallel-connected strings (n p ) to form the PV array: I s,array = n p I s,module V T,array = n s V T,module R s,array = n s n p R s,module R sh,array = n s n p R sh,module .
From now on, the parameters I ph , I s , n, R s and R sh will refer to parameters of the PV array model, without the need of using the subscript. The SDM model is defined for specific values of the incident irradiance and the temperature of the PV array. At Standard Test Conditions (T re f = 25 • C and G re f = 1000 W/m 2 ), the SDM parameters are known as reference parameters and denoted as: I ph0 , I s0 , n 0 , R s0 and R sh0 . De Soto et al. [18] developed the translation equations of the reference parameters to different operating conditions as follows: where α Isc is the short-cirucit temperature coefficient, where Eg is the energy bandgap of the p-n junction's material, Equation (1) can be converted into an explicit equation using Lambert W function, as explained in Reference [19]. Introducing the dependence of the SDM parameters with irradiance and temperature leads to the expression of the photovoltaic current determined by the environmental conditions and the photovoltaic voltage as follows: Figures 3 and 4 show the accuracy of the proposed model compared with the data available in the Matlab Electrical Simscape model [20]. Figure 3 shows the I-V curves for three different irradiance levels (600 W/m 2 , 1000 W/m 2 and 1400 W/m 2 ) with the same temperature (25 • C) and Figure 4 shows the I-V curves for three different temperatures (25 • C, 50 • C and 75 • C) and the same irradiance (1000 W/m 2 ). The array is formed by 8-series connected 1Soltech 1STH-250-WH solar panels. For the implementation of the Lambert W function in this work, the built-in lambertw Matlab function has been used. It is worth noting that, for real-time implementation, approximations of the Lambert W function as References [21][22][23][24] are more suitable in order to reduce the computation time.

Non-linear Least Squares Curve Fitting
Traditionally, PV systems have been operated in order to extract the maximum available power employing MPPT algorithms. However, in some applications such as power smoothing or power reserve control, PV systems need to be operated in a deloaded mode. In these cases, it is important to estimate the PV system operating conditions, just to know the available power (MPP) or the entire I-V curve. One approach to solve this issue is to apply nonlinear least squares curve fitting [25]. In this work, the curve fitting problem has been solved with the Levenberg-Marquardt method [26,27] by means of the Matlab embedded function lsqnonlin, based on Reference [28]. In nonlinear least squares problems, a set of m measurements (V pv ,I pv ) is fitted to a nonlinear model by minimizing the vertical distance between the measurements and the model as: where θ = {I ph , I s , n, R s , R sh } and f (V pv,i , θ(G, T cell )) is: From (20), it can be seen that the objective function depends on measurements, the incident irradiance and the cell temperature. Therefore, for a specified set of (V pv ,I pv ) measurements, it is possible to study the change in the objective function caused by variations in irradiance and temperature. A detailed analysis has been accomplished considering three different measurement datasets along the I-V curve, as depicted in Figure 5.   As can be seen in Figure 6, there are multiple local minima when the operation point lays on the left part of the I-V curve. In addition, these minima have similar values and are very differentiated in terms of temperature, which makes escaping from the local minimum of the initial solution complicated. A different outcome can be achieved when measurements are close to the MPP. Figure 7 illustrates this case. In contrast to the previous set of measurements, when the PV array is operating in the MPP region, a global minimum can be found and the optimization problem can be solved satisfactorily. In this last scenario, when the PV system is operating on the right part of the I-V characteristic, the behavior of the objective function is depicted in Figure 8. Again, there are multiple local minima with similar objective function values. Nevertheless, the temperature range and the objective function jumps between local minima have been noticeably reduced, which facilitates the search.
It is worth noting that, in the three cases analyzed, the measurements belong to a unique curve. But, in fluctuating weather conditions, irradiance, temperature, or both, vary between measurements. This case is illustrated in Figure 9, where irradiance is increasing and cell temperature is decreasing from the first measurement (purple) to the last measurement in the window (yellow). By computing the objective function for different combinations of irradiance and temperature, one can find that there are up to ten local minima, as shown in Figure 10. This study reveals that it is not possible to estimate the cell temperature while operating on the left part of the I-V curve. Nonetheless, and as the cell temperature dynamics are slower than irradiance dynamics, operation on the left-hand side is feasible, at least for short periods.

Photovoltaic Voltage Control
For the implementation of the proposed methodology, a detailed model of the photovoltaic and the voltage control system has been developed. As can be seen in Figure 11, the PV system is composed of a PV array and a boost converter that increases the array output voltage to meet the voltage requirements of the DC-Link. For the purpose of this study, it has been assumed that the voltage at the DC-Link is constant and it has been modeled by a DC voltage source. The voltage at the PV array terminals and, therefore, the PV array operation side are controlled by the voltage control system, as explained below using a block diagram. The voltage control system receives voltage and current measurements from PV array terminals continuously. These measurements are discretized by means of two Zero Order Hold (ZOH) blocks, with a sample frequency of 10 kHz. Subsequently, measurements are stored in 100-sample buffer blocks, before feeding a triggered Matlab function. The function also receives the active power reserve command and the desired operation side as arguments. Every 10 ms, the Matlab function executes the lsqnonlin function, and solves the nonlinear least squares curve fitting problem. As a result, the actual P-V curve is estimated and the MPP can be easily identified. Furthermore, the reference voltage is obtained taking into account the reserve command and the operation side. This reference voltage is then compared to the actual PV voltage, and the error is fed to a PID controller, which produces the duty cycle signal. The duty cycle, whose value is bounded within the interval [0,1], is the input for a Pulse Width Modulation (PWM) generator, which generates the control signal that opens and closes the switch of the boost converter. Thus, the control loop is closed and the active power reserve and the operation side can be controlled. It is well known that both the PV array and the boost converter have nonlinear characteristics, which hinders the PID controller tunning. To overcome this issue, a small-signal model of the system is developed as in Reference [29]: whereṽ pv andĩ L are the state-space variables,d is the control variable and L, C in and V o are parameters of the system shown in Figure 11. R pv is the dynamic resistance of the PV array, defined as the quotient between a small increment in PV voltage (ṽ pv ) and a small increment in PV current (ĩ pv ): From the state-space Equation (22), it is possible to obtain the transfer function of the linearized model and express it in the canonical form of a second order transfer function: where the damping factor ξ, the undamped natural frequency ω n and the static gain K o of the plant can be expressed as a function of the model parameters as follows: It can be seen from (24)-(26) that the plant model G(s) and the damping factor ξ are affected by the dynamic resistance and, therefore, by the operating point of the PV system. To guarantee a robust control of the PV system, it is necessary to contemplate all the possible operating conditions. For illustration, Figure 12 represents the evolution of the dynamic resistance along the I-V curve for five combinations of irradiance and temperature. As depicted, R pv can take values in a wide range from −3000 to near 0 Ω. In order to design a robust regulator, affine parameterization is applied.

Affine Parameterization
Affine parameterization [30] is a straightforward method to tune a compensator in terms of robustness and stability of the closed-loop system. The design process starts by defining the expected closed-loop characteristics of the system: the closed-loop damping factor ξ cl and the closed-loop natural frequency ω cl . Table 2 shows the percentage of overshoot for different closed-loop damping factor values. Once the closed-loop damping factor is selected, the closed-loop natural frequency can be derived from the settling time expression (29): The compensator's transfer function can be derived directly as follows: After the calculation of C(s), the closed-loop robustness and stability must be verified. To do this, the Bode plot of the open loop transfer function C(s)G(s) can be used. Figure 13 shows Bode plots of the open-loop transfer function for five different operating conditions. The gain and phase margins guarantee a closed-loop stable operation for the five possible scenarios presented.

Simulation Results
The effectiveness of the proposed control strategy has been tested in three simulation case studies in MATLAB/Simulink.

Case Study #1
The first study consists of a scenario with constant irradiance and temperature conditions. In this context, while operating at MPP, the command reserve is set to 50% at t = 1 s. Then, at t = 1.6 s, the operation side is shifted, so the PV system passes from left to right side operation and vice versa. The simulation compares the performance of the control system when the PV system operates on the left and on the right sides of the P-V curve. Figure 14a illustrates the evolution of photovoltaic voltage for both cases (MPP-left-right and MPP-right-left) and the maximum power point voltage. As evidenced, in the MPP-left-right case, the photovoltaic voltage is below maximum power point voltage from t = 1 s to t = 1.6 s and above it from t = 1.6 s to t = 2 s. Figure 14b depicts the level of reserves achieved in both cases. According to the results, the control system is able to provide the required amount of reserves while operating on both sides of the power-voltage characteristic. This is confirmed in Figure 14c, where the actual P-V curve at t = 1.5 s is represented together with the MPP and both measurement windows.

Case Study #2
In this scenario, the irradiance and the required reserves vary while the temperature is kept constant. The PV system operates on the left part of the P-V curve. Figure 15a shows the available power and the photovoltaic power over the simulation period. As can be seen in Figure 15b, the requested reserve is satisfied at every moment of the simulation. Operation on the left part of the P-V curve is guaranteed, as shown in Figure 15c, because the photovoltaic voltage is always lower or equal to the maximum power point voltage. Finally, the irradiance evolution and estimation are represented in Figure 15d. Despite the discharge level, the proposed estimation method shows an efficient performance.

Case Study #3
Case Study #3 compares the proposed method with the one presented in Reference [13]. The temperature, irradiance and required reserve profiles have been replicated from the aforementioned work. The objective of this simulation is to get the PV system to provide the requested amount of reserves while operating on the right-hand side of the P-V curve. Figure 16a,b illustrate the estimated temperature and irradiance by both methods and their real values. Figure 16c shows the evolution of the ripple for both strategies. As can be seen, while the Batzelis' approach uses a ripple control to follow the ripple reference of 3%, the proposed methodology takes advantage of the inherited ripple. It should be pointed out that the inherited ripple is greater than the reference just in the case of a power reserve level of 90%, from t = 12 s to t = 19 s. However, when the PV system operates at its maximum available power point, from t = 0 s to t = 5 s and from t = 27 s to t = 30 s, the inherited ripple is very low, reducing power and voltage fluctuations. Figure 16d depicts the power evolution along the simulation. As evidenced, the Batzelis' and proposed methods have the same behavior. Figure 16e illustrates the requested power reserves and the one achieved by both methods. As it is clear, both strategies are able to follow the power reserve command. Finally, Figure 16f shows the squared 2-norm of the residual as a result of the curve fitting process. According to this result, it seems that the objective function proposed in this paper is more suitable than the one of Batzelis' proposal.

Conclusions
In this work, a novel photovoltaic power curtailment method has been developed and implemented.
Unlike the existing methods in the literature, the proposed strategy allows PV systems to operate on the right and on the left side of the MPP, depending on the needs. It has been shown that left-hand side operation is preferable in case of fluctuating irradiance conditions. However, this operation side suffers from a limited power regulation interval. On the other hand, right-hand side operation has a wider power regulation range, but it may be unstable in a variable irradiance scenario. The proposed strategy avoids these limitations. In the curtailed mode, a curve fitting algorithm estimates the actual P-V curve and the maximum available power. The solution of the curve fitting problem has been analyzed, and the conclusion was that it is impossible to estimate the temperature when the PV system operates on the left-hand side of the MPP. Thus, our proposal is to operate photovoltaic systems on the left part of the P-V curve when there are sudden changes in irradiance and operate them on the right when curtailment is long enough that temperature may change. To evaluate the developed strategy, three case studies have been thoroughly designed. In the first case study, it has been demonstrated that with this methodology it is possible to shift the operation side without altering the available power reserve. In the second one, a scenario with a trapezoidal irradiance profile is presented. When power curtailment is required, the PV system is operated on the left part of the MPP. In these conditions, the incident irradiance and the MPP were estimated with a high level of accuracy. The last scenario compared the proposed strategy with the one suggested by Batzelis [13]. The irradiance, temperature and power reserve profiles have been replicated from the aforementioned work. With these conditions, the proposed methodology showed a better performance in terms of the curve fitting process. In the last scenario, it can be also seen that the inherited ripple may suffice in order to obtain an adequate window measurement.
Finally, for future research, the authors suggest an experimental validation of the proposed methodology and the realization of studies on frequency stability in weak electric power systems with high photovoltaic penetration.