Photovoltaic Cell and Module I-V Characteristic Approximation Using Bézier Curves

The aim of this work was to introduce new ways to model the I-V characteristic of a photovoltaic (PV) cell or PV module using straight lines and Bezier curves. This is a complete novel approach, Bezier curves being previously used mainly for computer graphics. The I-V characteristic is divided into three sections, modeled with lines and a quadratic Bezier curve in the first case and with three cubic Bezier curves in the second case. The result proves to be accurate and relies on the fundamental points usually present in the PV cell datasheets: V o c (the open circuit voltage), I s c (the short circuit current), V m p (the maximum power corresponding voltage) and I m p (the maximum power corresponding current), and the parasitic resistances R s h 0 (shunt resistance at I s c ) and R s 0 (series resistance at V o c ). The proposed algorithm completely defines all the implied control points and the error is analyzed. The temperature and irradiance influence is also analyzed. The model is also compared using the least squares fitting method. The final validation shows how to use Bezier cubic curves to accurately represent the I-V curves of an extensive range of PV cells and arrays.


Introduction
The forecast of the total photovoltaic (PV) installations, offered by Bloomberg New Energy Finance (BNEF) predicts an optimistic growth at 111 GW in 2018, rising to 121 GW in 2019, along with a polysilicon factory growth boom and module prices drop to US$ 0.30/W [1,2]. This robust growth explains the high interest in PV research, modeling, and simulation-along with design and development of PV equipment.
The electrical characteristics of the PV cell and PV modules have been of interest for several decades, and different models have been proposed. Phang and Chan [3] were among the first to propose a solution for PV cell parameter extraction. Garrido-Alzar [4] uses a double exponential model to extract the PV cell parameters using the experimental I-V curve. Villalva et al. [5] developed an algorithm to find the parameters defining the I-V characteristic for the single diode model of a PV cell, using the Newton-Raphson method and imposing a minimum error threshold for the maximum power point. Babu and Gurjar [6] introduced a simplified two-diode model for the PV module, while Cubas et al. [7] used the Lambert W-Function for finding the solar panel equivalent circuit parameters and also proposed an LTSpice model according to their findings. Temperature influence has been studied by Chander et al. [8]. Ishaque and Salam [9] use differential evolution to find the PV modules parameters. Franzitta et al. provided extensive studies of the most widely used models both in single diode case [10] and for the double-diode version [11], introducing a criterion for rating the accuracy and usability of the analyzed models. In our previous paper [12], we proposed a complete SPICE model including all the parameters variation and selfheating. In a recent work, Cuce et al. [13] claim a good accuracy for their electrical model for a PV module and they also discuss energy and exergy efficiency as a reliable substitute for the fill factor. All the aforementioned works use an electrical model to describe the behavior of the circuit and rely on a specific circuit to generate the I-V characteristic of the PV cell or module.
This paper introduces a new approach. This time, the cell or module are not involved at the electrical level, being defined by just the specific points P sc , P mp , P oc and by the parallel and series resistances R sh0 and R s0 , specified at I sc and V oc , respectively. Careful inspection of the typical I-V characteristic of the PV module or PV cell (Figure 1a,b) shows a similar pattern in all curves. Our aim was to find a way to model it using smooth curves and the datasheet information currently available. Bernstein polynomials have been studied since the beginning of the 20 th century and they form the foundation for Bézier curves [14]. The core applications for graphics came first in 1959 when the French mathematician Paul de Casteljau developed an algorithm able to evaluate a family of specific curves at Citroën. In 1962 the French engineer Pierre Bézier also used them to design automobile bodies at Renault and afterwards they achieved wide acceptance [15].
Bézier curves are largely used in computer graphics [16,17] and also in the time domain for smoothing the trajectory of the robotic arms, for an accurate gluing or welding path or for trajectory generation [18]. Further development for shape representation was proposed by Jalba et al. [19].
The current proposal analyses the use of Bézier curves [20] in order to accurately represent the I-V characteristic of a PV cell or module. Although Bézier curves are commonly used as 2D curves (which are usually not functions) [21], our approach uses such curves to approximate an I-V characteristic (a function). A complete mathematical solution is provided, separately validated for a PV cell and a PV module and the error is analyzed. The results are also studied for different temperature and irradiances and finally compared with the ones offered by the least squares fitting method. In depth analysis of Bézier cubic curves fitting for 18 PV arrays and cells (from various manufacturers and different technologies) is also performed.
There are other possible approaches, for example spline fitting [22,23], well known Hermite interpolation, or spiral curve fitting [24]. Our aim was solely the Bézier curves.
The remainder of this paper is organized as follows: Section 2 briefly analyzes the definition of the quadratic and cubic Bézier curves and their equations, focusing on the basic knowledge needed in the subsequent paragraphs. Section 3 deals with the information usually provided by the PV cell or module manufacturers in their datasheets. The proposed models are covered in Section 4, with Section 4.1. introducing the approximation with two segments and one quadratic Bézier curve, while Section 4.2. deals with the better approximation based on three Bézier cubic curves. The models are rated at the reference temperature, 25 • C. The results are also compared with the least squares fitting method in Section 4.3. In Section 4.4. the influence of the external parameters is analyzed. The proposal is verified in Section 4.5. against a large range of PV cells and modules and the results show a good fit. Discussion and conclusions are provided in the next Sections.

Definition of the Bézier Curves
A quadratic Bézier curve ( Figure 2) can be specified by three control points [15]: the curve goes through the ends P 0 and P 2 and approximates P 1 .

Figure 2.
A quadratic Bézier curve representation. P 0 and P 2 are the end points, the control point P 1 is approximated and the curve is tangent to P 0 P 1 and P 2 P 1 segments at P 0 and P 2 respectively. The curve equation is as follows [20]: where t varies between 0 and 1. Equation (1) can be expressed for the x and y coordinates: The derivative of (1) is: At the end points, (3) becomes (4): A cubic Bézier curve ( Figure 3) can be specified by four control points [15]: the curve goes through the ends P 0 and P 3 and approximates P 1 and P 2 . The analytical expression of the curve is a cubic polynomial. The curve is tangent at P 0 to P 0 P 1 and at P 3 to P 3 P 2 .
The equation for the Bézier cubic curve is [20]: The previous equation can be expressed for the x and y coordinates: The derivative of (5) is: At the end points, (7) becomes (8): Figure 3. A cubic Bézier curve representation. P 0 and P 3 are the end points, the control points P 1 and P 2 are approximated, and the curve is tangent to P 0 P 1 and P 2 P 1 segments at P 0 and P 3 respectively.

Materials and Methods
The first PV cell used in our work is a high efficiency, silicon monocrystalline 156 × 156 mm 2 cell [25] and has the main characteristics summarized in Table 1. Series resistance at V oc 3.8 mΩ The MSMD290AS-36_EU Monocrystalline PV module formerly manufactured by MünchenSolar Germany [26] has been well documented and studied by Cubas et al. [6]. Its main electrical data is listed in Table 2 and this information is used in Section 4.4. to evaluate the influence of the temperature and irradiance to our Bézier curves based model. For studying and representing Bézier curves, an interesting application which allows draggable control points was developed by Mugnaini [27]. For computing the coordinates on the curves we used the Kronecker tensor product found as in [28]. An example for the Bézier least square fitting method is given in [29].

Results
It must be stressed that all the physical actual values involved in Section 4.1., Section 4.2., and Section 4.3. are specified at 25 • C, being reference values. The irradiance is also standard (1000 W/m 2 ). As it is demonstrated in Section 4.4., the same method is suitable for different temperatures and irradiances.

I-V Characteristic Approximation with Two Segments and a Quadratic Bézier Curve
The first approximation implies five control points: P sc (0, I sc ), P a (x a , y a ), P b (x b , y b ), P c (x c , y c ), and P oc (V oc , 0) and is made of two segment lines I sc P a and P b V oc and one quadratic Bézier curve defined by the endpoints P a and P b and the control point P c (Figure 4). It has already been proven [3] that the slopes of the lines can be written as (9) and (10): Thus, the equation for first line is (11): By choosing x a in the linear region (eg 0.6V oc ), one can find y a from the above equation, so P a is completely defined.
For the second line, the next equation is valid (12): The P c (V c , I c ) control point has therefore the coordinates defined by (13): For x b , it must be emphasized that its position is on the end of the curve, a realistic value being 0.9V oc . The maximum power point is positioned on the second curve, so solving (14) gives t mp : Replacing the positive solution for t mp in (2) yields y b as in (15): Now all the control points of the plot are completely defined. The results are summarized in Table 3. The application code written for the coordinate finding can be retrieved from [30]. The final plot is represented in Figure 5, where one can observe an excellent correspondence between the actual PV cell I-V characteristic, represented with black dots and the P sc P a segment (blue line), a fair correlation for the second range, approximated by the Bézier quadratic curve (red line) and some modest results in the third region (magenta line). The same conclusion arises from Figure 6, where the relative error has been plotted. It is worth mentioning that although the relative error is quite high above 0.64 V (0.92V oc ), the absolute error is in fact less than 0.7 A in a region where the cell normally should not operate. Looking for a more accurate model is the reason we came up with the second scenario, where the I-V characteristic is entirely modeled with cubic Bézier curves.

I-V Characteristic Approximation with Three Cubic Bézier Curves
In order to have a general solution, we analyzed the case where all three regions are covered with cubic Bézier curves. This implies 12 control points (Figure 7), i.e. 24 coordinates to be found. The first curve, represented in Figure 8, is described by the control points P 00 (0, I sc ), P 01 , P 02, and P 03 . It turns out that the linear approximation of the first region of the I-V curve has an error below 0.5% if P 03x = V oc 2 . During various simulations we also discovered that all P j,k points can be evenly arranged, with j = {1, 2}; k = {1, 2, 3}. This leads to P 01x = P 03x 3 and P 02x = 2P 03x 3 . For the y coordinates, P 0ky = I sc − P 0kx

R sh0
, with k = {1, 2, 3}. Now the first curve is completely defined. The second curve (Figure 7) is described by the control points P 10 = P 03 , P 11 , P 12 , and P 13 . We observed that P 13x = 0.75V oc offers a very good fit of the curve for this type of PV cell. With the same evenly arrangement for the x coordinates, P 1kx = P 10x + k(P 13x −P 10x ) 3 , with k = {1, 2}. P 11 is also located in the linear region of the I-V curve, so P 11y = P 03y − P 11x R sh0 . This leaves P 12y and P 13y as unknowns at this stage.
The third curve (Figure 7) is described by the control points P 20 = P 13 , P 21 , P 22, and P 23 (V oc , 0). Using the same assumptions as for the second curve, P 2kx = P 20x + k(P 23x −P 20x ) 3 , with k = {1, 2}. It is obvious that P 20x = P 13x and P 23x = V oc . The segment P 22 P 23 is tangent to the curve at the point P 23 , . This leaves P 21y as an additional unknown at this step.
For continuity reasons, P 12 P 13 and P 20 P 21 segments belong to the same line. This implies that the derivatives of the second curve at P 13 and of the third curve at P 20 are equal (16): (16) Which means that: The control point P 11 is placed on the second curve, so (18) can be written: Solving the previous equation and keeping only the real solution for t 11 , (19) is also valid: Finally, the graph also goes through the maximum power point P mp V mp , I mp , yielding Equation (20): Keeping only the real solution for t mp , (21) is also valid: The linear system made of equations (17), (19), and (21) give the last three unknown coordinates P 12y , P 13y , and P 21y . The results are summarized in Table 4. The application code written for coordinate finding can be retrieved from [30].  Figure 9 shows the location of the control points with respect to the I-V characteristic of the PV cell. The control points P 00 , P 01 , P 02 , P 03 = P 10 and P 11 are collinear and with P mp , are all placed on the I-V characteristic.  Figure 10 shows the modeled characteristic (red, green, and blue lines) overlapping in most areas with the practical I-V characteristic (black markers). The application code can also be retrieved from [30]. The relative error of the Bézier modeled I-V characteristic against the actual data taken from [12] is shown in Figure 11. It must be emphasized that in the 0-0.94 V oc range, the relative error is below 1%. Above 0.94 V oc the absolute error is less than 72 mA, while the reference I sc = 9.207A.

Data Fitting Using the Least Squares Method
Data fitting using the least squares method is a standard approach in data analysis [31,32]. A good overview of curve fitting using Bézier cubic curves in image processing is given by Shao et al. in [33], while Zhao et al. [34] extend this method using a genetic algorithm for parameter optimization for Bézier curve fitting. In Section 4.2. we have shown that for the studied PV cell, the best results arise when the x coordinates of the middle end points are set at 0.5V oc and 0.75V oc respectively. A similar conclusion arises if the least squares method is used for the same cell modeling. Running the least squares method for the MSMD290AS-36_EU Monocrystalline PV module proved that the minimum error occurs when the control end points are set again at 0.5V oc and 0.75V oc respectively. Table 5 summarizes the data fitting results for the same PV cell used in Section 4.1. and Section 4.2., where the results from the two approaches are very close. The graphical representation of the date fitting is given in Figure 12, where only the endpoints are represented.  Figure 13 depicts the relative error of the modeled I-V characteristic compared with the actual data taken from our previous work [12]. In the 0-0.96V oc range, the relative error is below 2%. Furthermore, above 0.96V oc , the absolute error is less than 66 mA, while the reference short circuit current is 9.207 A.

Parameter Variation
In order to further validate the proposed method, in this section we analyze the temperature and irradiance influence for the MSMD290AS-36_EU Monocrystalline PV module. An extensive study of the parameter influence over the PV cell can be found in [12]. It is important to notice that the Bézier approximation is not related to any of the parameter variation, just to the specified points P sc , P mp , P oc and the parasitic resistances R sh0 and R s0 as already stated. The challenge becomes in this case the finding of the new position for the control points and the new values for the parasitic resistances.
Villalva et al. [5] accurately describe the short circuit current variation as in (22): Ishaque and Salam [9] propose for the V oc,cell the following variation (23): Equation (23) proved to be too conservative in this case, as larger V oc,cell variations were observed. A better approximation is the empirical law (24): A possible way for defining R sh behavior is suggested in [12], as in (25) with k Rsh estimated as 8 for the best fit.
For R s , a linear variation law (26) is given in [12] with α Rs = −0.01K −1 , again for the best fit: Figure 14 shows the irradiance influence for the I-V module characteristic, where the approximated data using our proposed method is plotted with solid lines and the experimental data is represented with markers. I sc , V oc , R sh , and R s were computed using (22), (24)(25)(26) respectively.

Figure 14.
Bézier approximation of the I-V irradiance dependent characteristics for the MSMD290AS-36_EU monocrystalline PV module. The lines represent the computed curves, whereas the markers represent the actual data.
The temperature dependent Bézier curves resulted from our algorithm compared with the actual data are introduced in Figure 15. Once again, the results show a very good correlation between the modeled data and the actual data.

No.
PV Type Tech n s V oc (V) Using the Villalva algorithm [5], the main parameters were computed and are listed in Table 7. For the last three cases, due to different technology, interesting values for the diode ideality factor a occur.  Table 8 lists the control points computed for Bézier curve fitting, where in all cases the control end points are set at 0.5V oc and 0.75V oc respectively. Only two x coordinates are presented, as the others are evenly spaced and can be easily computed. Selected plots of the PV devices are presented in Figure 16. In the second column of Table 9, the average of the current (I) relative error is displayed in order to evaluate the fit quality. Maximum current error (absolute and relative, listed as positive values) is associated in each case with the coordinates where it appears. The maximum power point is also investigated as an absolute and relative error and finally the computed value is listed. The current (I) relative error is below 1.18% in all cases, and the P mp relative error is even better (lower than 1%).    An additional plot of the relative error is introduced in Figure 17a,b for two of the analyzed PV devices. The same pattern arises as in Figure 11, with small errors up to 0.94V oc . The curvature plots from Figure 17c,d were computed using the Herman and Klette method [50], with 20 forward and backward points, using an adapted script from [51].

Discussion
In all studied cases, the x coordinates can be evenly spaced. Both for the PV cells and for the PV modules, the first Bézier cubic curve was very close to a straight line and ended at x 03 = 0.5V oc , while the middle curve ended at x 13 = 0.75V oc . The relative error is less than or equal to 1.18% for all the studied PV devices, while neglecting the values in the 0.94V oc − V oc region, where the absolute error is tens of mA and the relative error is misleading, (Figures 11 and 17).
The differences between the results obtained with the proposed method and the least squares method are negligible (less than 1% for all coordinates). R s and R sh can be easily derived from the manufacturer's datasheet, using for example the method proposed by Vilallva et al. [5]. In most cases R sh is close to R sho (Table 7). Larger differences occur for the series parasitic resistance: for the same PV cell. For different irradiances and temperatures, Section 4.4 provides all the necessary formulae.
In most cases V mp > P 13x . The method was also valid for Onyx Ref 10, Onyx Ref 30, and 6.5 Wp L Cell, where the previous relation was not satisfied.
In relation to the Maximum Power Point Tracking (MPPT), it must be emphasized that the Bézier curves are inherently smooth. This reduces the risk for the algorithm of becoming stuck in some false area/minimum of the curve. Furthermore, as the generation of the Bézier curves is so easy, the MPPT simulation can be further simplified.

Conclusions
A novel method for modeling a PV cell or a PV module I-V characteristic was introduced. To the best of our knowledge, Bézier curves have not been used to model the I-V characteristic of PV devices before. The method proved high accuracy and was validated both in the case of a single PV cell and in the case of a whole PV module, for different technologies and manufacturers. The method was also used in the case of varying irradiance and temperature. The proposed method can be used for implementing hardware solar array simulators, for teaching or remote study. It is far easier to use the proposed method to find the I-V characteristic of a PV cell or module when compared with solving the exponential equations associated with the single or double diode model largely used today. A common microcontroller can compute the points on the I-V curve with a minimum of resources, inherently increasing the computing speed and the response of the system.
The advantage of our method relies on the ease of I-V characteristic generation: if we exclude V oc and I sc , only 10 different values (P 23x = V oc , P 00y = I sc , P 01y , P 02y , P 03y , P 11y , P 12y , P 13y , P 21y , P 22y ) have to be stored-SAS manufacturers usually use 1024 or more double points to accurately define the I-V characteristic. Alternate use is for any graphical plot of the I-V (and subsequently P-V) curves. Furthermore, little knowledge of the device itself is required, as only common data from the datasheet is needed.
As future work, fitting spiral segments will be investigated and compared with our current method.
Author Contributions: Both authors contributed to this research. Aurel Gontean conceived and designed the study and carried out the simulations. Roland Szabo analyzed the data. Aurel Gontean wrote the paper and reviewed the manuscript. Both authors read and approved the manuscript.
Acknowledgments: This work was supported by both the Romanian National Authority for Scientific Research and Innovation, CNCS/CCCDI-UEFISCDI within PNCDI III, project number PN-III-P2-2.1-PED-2016-0074 and by the Politehnica University Timisoara, according to the Administration Board research policy.

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