Stability and dynamic analysis of a grid-connected environmentally friendly photovoltaic energy system

Photovoltaic (PV) system is the cleanest form of electricity generation, and it is the only form with no effect on the environment at all. However, some environmental challenges persist, which must be overcome before solar energy may be used to represent a source of truly clean energy. This paper aims to study the stability and dynamic behavior of a grid-connected environmentally friendly photovoltaic energy system using the bifurcation theory. This theory introduces a systematic method for stability analysis of dynamic systems, under changes in the system parameters. To produce bifurcation diagrams based on the bifurcation theory, a parameter is constantly changed in each step, using MATLAB and AUTO, and eigenvalues are monitored simultaneously. Considering how the eigenvalues approach the system’s imaginary axis in accordance with the changes in the targeted parameter, the occurred saddle-node and Hopf bifurcations of the grid-connected PV system are extracted. Using the obtained bifurcations, the system’s dynamic stability limits against changes in controlled (controller coefficients) and systematic parameters (such as the Thevenin impedance network) are found.


Introduction
In recent years, a large number of studies have been conducted on renewable energies such as solar, wind, sea wave, tidal, and geothermal. Not only are they renewable, but they are also environmentally friendly and make no or little, bearable pollution. Today, for the above reasons, the implementation of such resources in power systems is more than it has ever been (Khokhar et al. 2015). This has resulted in less number of casualties and lower costs of upgrading supply and transmission lines while increasing the reliability and security of electrical distribution and, in turn, decreasing the usage of fossil fuels (Jahannoosh et al. 2020).
Photovoltaic systems (PV) are one of the technologies employing solar energy, by converting the received energy into electricity (Parizad and Hatziadoniu 2020). Over the past few years, photovoltaic panels have had the highest growth rate of installation capacity compared to other technologies (Naderipour et al. 2019). The lifetime of about 20 years; ease of installation, setting up, and maintenance; tolerance toward all weather conditions; and low noise pollution are only some of the reasons for the preference of photovoltaic panels over other renewable resources (Naderipour et al. 2017).
Germany, Belgium, and Spain are leading countries in encouraging electricity consumers to use the photovoltaic systems . It is predicted, however, that in the future, with the development of technology and lowered costs, the installation rate of photovoltaic panels in sunlightreceiving countries is increased (Marra et al. 2013).
Responsible Editor: Philippe Garrigues * Amirreza Naderipour namirreza@utm.my 1 Therefore, due to the rising need for multi-megawatt photovoltaic systems, and their undoubted positive effect on the producer and the consumer, researching the control, dynamic specification, and transient operation of PV systems of a big scale seems inevitable (Jafar-Nowdeh et al. 2020).
The broadest technical term related to PV systems is maximum power point tracking (MPPT). There are 19 different methods of MPPT reviewed in Esram and Chapman (2007). General investigations have been carried out on three-phase and single-phase converter circuits for the application of PV in Haeberlin (2001). Stability of scattered and single-phase PV systems is studied in Wang and Lin (2000). The effect of network impedance changes on the stability of a singlephase closed-loop PV system is investigated in Liserre et al. (2006). Moreover, a control scheme is suggested to ensure the generation of sufficient attenuation. Liserre et al. (2006) proposes a reduced-order model for a PV system, in which the model and controllers are based on the voltage mode control. This model can be applied in studies of power system simulation in the time domain. This reference, however, does not mention the stability analysis and method of designing controllers. A novel compensation procedure is explained in Liserre et al. (2006) for DC interface voltage control loop to remove nonlinear effects of PV panels on the stability of the closed-loop. This proposed compensator allows the DC interface voltage controller to be designed, ignoring the working point of the PV system. It can also calculate the eigenvalues of a closed-loop PV system without sensitivity to sunlight and DC interface voltage levels.
A new issue facing today's developed power systems is the occurrence of load instability through voltage instability. Voltage stability is a subcategory of power systems' stability and their dynamic behavior following a disorder. A power system never remains in stable conditions for a long time. A disorder in a power system is a sudden change or a sequence of changes in a parameter or working point of the system. When the reason for voltage instability is the gradual increase of load or changes in the parameters of the system, the bifurcation theory shows how the voltage has collapsed. The power system is a nonlinear dynamic system. Once the system's load gradually increases, the equilibrium point of the system and, subsequently, the eigenvalues of the system's Jacobian matrix change. If the bifurcation parameter passes the bifurcation amount, it results in voltage instability, and in some cases, this instability is in the form of a rapid and steady reduction of the voltages.
An approach taken in the stability analysis of dynamic systems under parameter change is known as the continuation method which is discussed in bifurcation theory. The detailed explanation of the method can be found in Kuznetsov (2013) and Wiggins (2003).
Local borders of saddle-node and Hopf bifurcation of a modeled power system are extracted by a set of parameter-dependent differential-algebraic equations in Guoyun et al. (2005). The differences of voltage collapse in a simple double-engine power system, for each of the two modelsone due to the occurrence of saddle-node and the other for the Hopf bifurcation-are shown in Gu et al. (2007). The effect of static VAR compensator (SVC) on bifurcations, chaos, and voltage collapse of a power system containing nonlinear load is investigated in Subramanian et al. (2011). Moreover, the effect of adding some flexible AC transmission system (FACTS) equipment on the occurrence of Hopf bifurcation in a sample power system is discussed in Srivastava and Srivastava (1998). The results indicate that the addition of such equipment as a controllable series capacitor, phase angle regulator, and SVC to the studied power system can remove Hopf bifurcation from the system. Mahmud et al. (2019) presented to depict a comprehensive dynamic phasor model for DC and fundamental components of the SSSP gridconnected PV device. The active condition during the capacitance, as opposed to current grid-connected VSI networks, is another unique feature of this article. The eigenvalues and stability analysis of the dynamic phasor model are also examined, which are the essential evaluation for any dynamic structure. The evaluation was carried out with a variety of switching frequencies in order to provide a greater understanding of the accuracy of the VSI's switching frequency.
In El Aroudi et al. (2019), the nonlinear dynamics of a PV-fed high-voltage gain single-switch quadratic boost converter loaded by a grid-interlinked DC-AC inverter was studied. The operation of the input port of the converter is configured to maintain stability on a slow time scale, utilizing a resistive control method. However, at the quick time scale, the device will exhibit undesired subharmonic instabilities. In Aly (2016), a PV generator, a battery storage system (BSS), a PV inverter, a battery inverter/charger, and an AC control load are used with this system. In order to coordinate the electricity between the PV generator and the load, the PV generator is an intermittent source of energy, so this method needs to be coupled with the BSS (Bawazir and Cetin 2020). El Aroudi et al. (2019), the nonlinear dynamics of a PV-fed high-voltage gain singleswitch quadratic boost converter loaded by a gridinterlinked DC-AC inverter was studied.
The zero dynamic architecture method of feedback linearization is used to regulate the grid current and DC link voltage, which partly linearizes the device and allows controller design for reduced-order PV systems. Wang et al. (2020) used a modeling approach focused on matrix variables for a distributed PV grid-connected device. The main principle of the modeling approach is to transform the complex model comprising several PV-DC optimizer generation units into an average model consisting of just two normal submodules by constructing the variables generated by the block matrix. In Negi et al. (2020), it addresses an enhanced model guide adaptive control methodology for the design of voltage source converter control parameters and increases PV generation stability in various grid capabilities. In Huang et al. (2017), the integrated mathematical model of a single-phase two-phase grid-connected photovoltaic system is developed, in which both the DC-DC converter and the DC-AC converter are also considered to be the characteristic of the PV array; (2) an observer-pattern modeling approach is used to remove timevarying variables; and (3) the system's stability is analyzed using sensitivity and self-value. However, the stability of a grid-connected PV power system using bifurcation theory has, never before, been studied comprehensively.
This paper provides a systematic method for evaluating the stability of complex processes under adjustments in the parameters of the structure. A parameter is continuously modified in each stage to generate bifurcation diagrams based on the bifurcation principle, and eigenvalues are controlled simultaneously. The saddle-node and Hopf bifurcations of the grid-connected PV device are derived in light of how the eigenvalues approach the imaginary axis of the system in conjunction with the adjustments in the intended parameter.
This paper is organized as follows. The "Mathematical model of the PV system" section is discussed in the second section. In this section, details of the mathematical model, including the voltage-current characteristic in solar and power-voltage characteristic in solar, are explained. The third section presents the "Modeling of power grid-connected photovoltaic system" section. The "Dynamic modeling of a gridconnected photovoltaic power system" and "The PV system control" sections are presented in fourth and fifth sections, respectively. The sixth section presents the "Analyzing the proposed method in PV system" section. In this section, the dynamic behavior of a grid-connected photovoltaic power system is analyzed and elaborated. Simulation results are presented in the "Result and discussion" section. Finally, the "Conclusions" section is presented in the eight section.

Mathematical model of the PV system
A building block of PV arrays is the solar cell which is a semiconductor p-n junction. Solar cell directly converts the Sun's energy to electricity. Its circuit is shown in Fig. 1.
In this figure, the source of I ph current indicates the cell's photon current. This current is the short-circuit current of a set of PV panels. D is the sign for nonlinear resistance of p-n junction. Respectively, R s and R sh are symbols for intrinsic parallel and series resistance of the cell. Generally, the amount of R sh is great and R s is inconsiderable, and so both are ignored for easier analysis. Larger units of PV cells are categorized and create panels. A PV array can be one single panel or a set of panels, connected in a series or parallel mode-forming a large-scale PV system. To simulate this array, modeling is required. The below equation clarifies the mathematical PV model based on its current-voltage characteristic (Hussein et al. 1995): In this equation, I rs is the cell's reverse saturation current in AMPs, T is the p-n junction temperature in Kelvins, q is the electric charge, K is the Boltzmann constant, and A is the p-n junction's ideal index. The photon current I ph is a function in the form of temperature T and solar radiation S, as follows: where T r stands for reference cell temperature, I scr is a PV cell's short-circuit current in reference temperature and radiation, and K t is the temperature coefficient of the short-circuit current. According to the following equation, the cell's reverse saturation current I rs changes with temperature (Vachtsevanos and Kalaitzakis 1987).
where I rr is reverse saturation current at T r temperature, E G is the cell's semiconductor bandgap energy, and T c is the cell's temperature changes due to solar radiation in Celsius. Using the following equation, the generative power of the PV array is computed: Figures 2, 3, 4, and 5 show current-voltage and powervoltage characteristics of the PV array, as well as the effects of temperature and radiation parameters. As shown, the power curve in terms of voltage has one maximum point for each radiation level and a certain temperature. Known as MPPT, the goal, through controlling theV dc , is to match the operating point of the system with the above maximum point in different conditions.  (Mutoh et al. 2006) Modeling of power grid-connected photovoltaic system Figure 6 shows a single-line photovoltaic system diagram connected to a power grid at point of common coupling (PCC). The main PV system parts include PV panels array, a voltage-source converter (VSC) and a three-phase interface reactor. To ensure the large-enough production of power, PV panels array is constructed of a parallel connection of a set of n p panels. Moreover, to make sure of an adequately great DC voltage, each set includes n s PV panels, connected in a series format. PV array is connected in parallel to an interface capacitor (C) and the DC output of the voltage source converter. VSC is controlled based on the pulse width modeling (PWM) method. Interface reactor connects the voltage source converter's AC output to its corresponding PCC phase. Respectively, resistance and inductance of the interface reactor, combined with the resistance and inductance Thevenin network, are shown by symbols R and L. R also includes the VSC power key resistance. P s and Q s show the real powers and reactive output of the PV system, which are later delivered to the power network. PV system is connected to the power network, presented here as voltage source V s .

Dynamic modeling of a grid-connected photovoltaic power system
The dynamic equation of DC interface voltage based on power balance is stated as follows: where P dc indicates power delivery of voltage source converter on its DC side ( Figure 1). By ignoring VSC losses, P dc   can be estimated from P t which is the real power component of voltage source converter on its AC side. Also P t is the sum of P s and the real power absorbed by the interface reactor. Thus: Phasor from a three-phase quantityx abc (t) is defined as follows: By replacing P dc from Eq. (6) in Eq. (5), we have: Equation (8) specifies that V dc and consequentlyP pv can be controlled by the AC current of the voltage source converter. The current dynamic of the AC side of the voltage source converter is defined by the space-phasor equation as follows: The output voltage signal of voltage source converter's AC side, V ! t , is controlled by: where m ! indicates the corresponding space-phasor to PWM modulating signal, which is normalized toward the triangle carrier signal peak. By replacing V ! t from Eq. (10) in Eq. (9), the following equation is obtained: Equations (8) and (10) show the PV system's state-space representation, where V 2 dc and i ! are the state variables, m ! is the controlled input, and V ! s and S are the external inputs. It should be noted that V ! s is the Thevenin bus voltage of the power grid load.

The PV system control
In Figure 6, the main objective of the PV system is to control the DC interface voltage in order to control or maximize the capability of the PV array. Figure 6 shows that (1) voltage source converter via pulse width modulation and photovoltaic control system through a phase-locked loop (PLL) is synthesized with the grid voltage. Thus three-phase AC signals will be transferred to the appropriate corresponding dq reference. As a result, instead of the main signals with sinusoidal variation, their DC equivalent is entered and processed in controllers.
(2) The error between DC interface voltage square and corresponding reference value has been processed by k v (s) compensator, and subsequently its output is added to a leading Fig. 6 Single-line schematic diagram of a grid-connected single-stage photovoltaic (PV) system signal in order to generate i dref current command signal. The leading compensator has neutralized PV array nonlinear and unstable properties and improved PV system's stability. DC interface source voltage is generally obtained from one of the MPPT methods and variant in an allowed range. Restriction on DC interface voltage guarantees voltage source converter appropriate and safe operation. (3) i dref command signal is entered in a control-current plan in DQ source device which forces i d to trace i dref . In the next section, it is discussed that Q s is proportional to i q . If i dref equals to zero, thus Q s equal to zero in stable condition and power factor of the power grid will become unite in the PV system's connection point. Also, operation with unite power factor leads to minimum domain of VSC line current in sample real power dispatch.

Phase-locked loop space-phase variable
Space-phase variables in the model of the PV system are imaged on the dq source frame to analyze and control the system. This is done by the equation of dq source frame as follows: If x ! shows a state variable, then d x ! dt is obtained as: where ω is the angular velocity of dq source frame.
Transmission from space-phasor to dq source frame is given as: In the system shown in Figure 6, three-phase signals vary with ω 0 grid frequency. Thus, variables of dq source frame will be time-independent (at steady state), if angular velocity of dq source frame (ω) is been equaled toω 0 . This aim will be achieved by a PLL shown in Figure 7. First, V s is divided into d and q axis components from (15), then V sq gets into H(s) compensator to obtain ω.
At steady state, when ω becomes equal toω 0 , V sq should become equal to zero. Thus, H(s) should have one integrator leastwise (one axis onS = 0). H(s) is considered as a proportional-integral (phase-locked loop) series compensator with a low-pass first-degree transfer function, as follows: where Ω(s) is the Laplace transform of ω(t) and V sq (s) is the Laplace transform of V sq (t). If two state variables are defined as ζ 1 (t) = L −1 {Z 1 (s)} and ζ 2 (t) = L −1 {Z 2 (s)}, the equations are given as: According to the Eqs. (17) and (18): and by calculating the inverse Laplace from Eqs. (17) and (19): According to Eqs. (17) and (18), Eq. (16) is obtained as By calculating the inverse Laplace from the equations above: Thus, according to Eqs. (14), (20), (21), and (23), the model of PLL state space is summarized as follows: Equation (24) shows PLL as a dynamic system, where V sd and V sq are its inputs; ζ 1 , ζ 2 , and ρ are its state variables; and ω and ρ are its output.
The apparent output power of the PV system is as follows: According to Eq. (12), in the dq source machine, (25) is obtained as: Equalizing V sq to zero will make the real output power to be simplified as follows: Equation (27) indicates that P s is proportional to i d and it can be controlled by that. As indicated in (8), P s is controlled to adjust DC interface voltage and control or maximize power taken from PV array. Thus, as it can be seen in Figure 6, controlling the P s leads to controlling the i d . By equalizing V sq to zero, and according to Eq. (26), the reactive output power of the PV system is similarly simplified as: Thus, Q s can be controlled by i q to adjust the power factor of the PV system in the PCC.

VSC current control
In the previous section, it was determined that i d and i q must be controlled. In this part, to make sure that i d and i q follow i dref and i qref command reference values, a current-control layout is suggested. Besides, current-control strategy upgrades VSC protection against overload and external errors if it puts limitation upon i dref and i qref . Current-control system is designed based on Eqs. (29) and (30). These equations are formulated in the resource system d g according to Eqs. (11) and also (12) and (13): i d and i q are the state variables, and m y and m q are the control inputs (according to Eqs. (29) and (30)). i d and i q dynamic are interdependent and nonlinear. In order to make they are dynamic isolated and linear, m s and m q are determined based on general law followed blew: where u d and u q are two new controlled inputs. By replacing m s and m q in Eqs. (29) and (30), we will have: Equations (33) and (34) show an independent linear firstorder quadratic system in which i s and i q can be controlled, respectively, by u s and u q . Figure 8 shows the block diagram of the current-control layout in d g resource system (realization of Eqs. (31) and (32)). This figure represents u s as the output number one of compensator k d (s) which processes the error signal of e d equal to i dref − i d . It should be noted that the factor of 2/V dc is used as a leading signal for producing m s and m q to separate the dynamic of i d and i q from the dynamic of VSC. The PWM modulator signals are generated by changing m s and m q into m a , m b , and m c . These signals generate pulses anent VSC gates.
Since the controlled systems of (33) and (34) are similar, k d (s) and k q (s) can be designed similarly as follows: where k p and k i are, respectively, the proportional and integral gains of the PI controller. If k p andk i be selected as follows: Then, the transfer function of d and q axis current controller closed loop is first-order in which time constant T i must be small to reply the fast current-control as well as being selected big enough; so that 1/τ i (which is the width of current-control closed-loop band) becomes significantly smaller than VSC frequency switch (e.g., 10 times).
According to Figure 8 and Eqs. (33)-(37), G i (s) is obtained as: In this article, T i is selected equal to 0.5ms. Notably, the proposed controlled way separates the dynamic i d and i q from ω,V dc ,V sd ,and V sq .

Controlling front voltage DC
Front voltage DC is controlled i d based on relation (8). This relation in dq tracking device is equal to: In dq tracking device, first and second phrases show productive real power, respectively, by PV service and PV system. Third and fourth phrases show absorbed power, respectively, by the inductance and front reactor resistance VSC. Disregarding the third and the fourth phrases in relation (39) and because V sq is tuned in value zero using PLL, front voltage dynamic DC decreases in one model order and converts to the following: This equation shows a controlled system. The time constant τ i in relation (38) is reassumed to be at least small until the relation i d ≈ i dref is true and its input will be equal to i dref . In relation (40), V dc 2 are both a state variable and an output. The V sd is an external input (disorder). The controlled system in Eq. (40) is nonlinear because (1) P pv = f(V dc , S, T) is a nonlinear function of V dc . (2) V dc is a nonlinear function of state variable V dc 2 on its own. In order to decrease the nonlinear impacts of the controlled system, the following controlled law is considered for i dref : where u v is a new controlled input and γ is again that can be unit or null. The relation P pv = f(V dc , S, T) is a product of V dc and i pv . Substituting i dref of (41) for (40) results in the following: Relation (42) shows that if γ = 1, PV service characteristic impact on front voltage DC disappears, and an effective controlled system will be an integrator. Although the product V sd u v introduces a nonlinear property, because the V sd processing is a relatively constant variable, its impact is inconsiderable. The output u v is a compensator as the following: In which the error signal is e v ¼ V 2 dcref −V 2 dc . The above relation shows a compensator PI which has been series by low-pass filter of the first order. α 1 and α 2 are, respectively, proportional and integral gains of compensator PI. α 3 is a pole of the transfer function of low-pass filter. Similar to (16), the mode space model of (43) will be as the following: Figure 9 shows the diagram block of front voltage DC control plan which was introduced in the former part. It is noticeable that if γ = 1, the effective control system has been made up of u v to V 2 dc process by a series integrator with a firstorder transfer function G i (s).
Bifurcation is the emergence of phases in a non-equivalent topological format under changes in parameters (Kuznetsov 2013;Wiggins 2003). Therefore, bifurcation is changing the topologic scheme of the system, which is the result of parameter changes and passing a (critical) amount of bifurcation.

Analyzing the proposed method in PV system
In this section, with the help of given models in previous sections, the dynamic behavior of a grid-connected photovoltaic power system is analyzed and elaborated. System simulation in the time domain is programmed in the environment of MATLAB software, and the proposed method is applied to the system. By altering the monitorial and non-monitorial parameters which affect the stability of the quiescent point and using the linearization of the system around it, the Jacobean matrix of the system is extracted, and by analyzing how specific amounts of the system change and also using the proposed method, the type of stability points of the system, place of occurrence, and the type of occurred bifurcations are checked. This simulation is done in AUTO software. System information including parameters of the PV system is given in Table 1 (Yazdani and Dash 2009).
The PV system which is displayed in Fig. 6 with suggested monitorial technique comprises 11 variables of state and the differential-algebraic equations of which is given compactly below: Finally, computer simulation of the PV system can be done by using (45) to (66) relations.

Simulation of photovoltaic system in time domain
In this section, the simulation of the system is conducted with the assumption, V dcref = 1126V, y=1 and s=100. The system is  Parameters of the PV system n p = 176, n s = 1500, I rs = 1.2 × 10 −7 A, q = 1.602 × 10 −19 C, k = 1.38 × 10 −23 J/K, A = 1.92, I scr = 8.03A, k t = 0.0017A/K, T = 2980K, T r = 300K, C = 5000μF, R = 1mΩ, L = 100μH, ω 0 = 2π × 60rad/s, θ 0 = 0 ∘ , b V s ¼ 391V Controller parameters β 1 = 307.3(Vs 2 ) −1 , β 2 = 202.9(Vs 3 ) −1 , β 3 = 600s −1 , k p = 0.2Ω, k i = 0.6Ωs −1 , α 1 = − 0.77A(V 2 s) −1 , α 2 = − 328.2A(Vs) −2 , α 3 = 909s −1 run on MATLAB in the given conditions, and after reaching the stability point at t=4s, the step-by-step changes of the radiation levels 100 to 150 is considered as a disorder. The results of the dynamic simulation of the system are demonstrated in Figures 10, 11, and 12. Figure 10 demonstrates waves of the state variables ζ 1 ,ζ 2 , P, X 1 , and X 2 and frequency of the AC system related to the PLL (f) output. As it is seen, the PLL can follow the network frequency well. In addition, the mentioned disorder does not affect the performance of the PLL. The six remaining variables including current components dq, y 1 , y 2 , y 3 , and the DC interface capacitor voltage of the PV system (V dc ) are demonstrated in Fig. 11. As expected, since the source current of the q axis is assumed zero (i qef =0) the controller has been able to follow the real current of the q axis (i q ). Q is also zero in a stable state and the PV system delivers the electrical energy to the power network in unitypower-factor conditions. In Figure 12, V sq and V sd waves, error signals u v , e v , PV (P pv ) output power array, and source current of the d axis (i dref ) are given. It is noticeable that the system's dynamic behaviors are satisfactory while confronting external disorders.

Investigating the bifurcation phenomenon in photovoltaic system
In this section, the dynamic performance of the photovoltaic system is investigated under controlled and uncontrolled parameter changes such as proportional integral (PI) coefficients or the parameters of power system like the Thevenin impedance network and the referenced amounts in the control system. Furthermore, the effect of their changes on the system's bifurcation and finally its stability is analyzed. In other words, by selecting different parameters of the photovoltaic system as the bifurcation parameter, the possibility of the occurrence of different kinds of bifurcations is evaluated in the system. Investigated parameters include all photovoltaic variable parameters of the array such as temperature and radiant level, resistance and inductance of reactor, and q (I qref ) axis current reference input (which consists of resistance and inductance of the Thevenin network). It means that R and L are the domains of voltage network ( b V s ) and the controlled parameters of k p , k i , α 1 ,α 2 , α 3 , β 1 , β 2 ,andβ 3 . After comprehensive investigations on the mentioned parameters in the photovoltaic system in Figure 6 and the given model in Eqs. 45-66 in AUTO software, a supercritical Hopf bifurcation in the parameter of α 3 = 873.0017163779 is seen. In this condition, a double paired-pole with amounts of s = ± j0.285 is located on the imaginary axis and the system's working point loses its asymptotic stability. As α 3 decreases, the stable operating point of the system is lost and it takes a swinging state. This is the previously mentioned limit cycle in the supercritical Hopf bifurcation which was explained previously. The dynamic behavior of some state variables α 3 = 873 is shown in Figure 13. Some of these variables are not swinging. As the bifurcation parameter of α 3 decreases, fluctuation domain of these variables increases, and they may exceed their expected amount and system's relays order them to stop. In order to show a limit cycle, two state variables of γ 1 and γ 2 are drawn based on each other in Figure 14.

Conclusions
In the present paper, the analytic model of a gridconnected environmentally friendly photovoltaic energy system is shown in detail based on the differentialalgebraic equations and the bifurcation theory is discussed. This theory introduces a systematic method to analyze the stability of dynamic systems under changes in the system's parameters. Two different bifurcations happen in the system known as saddle-node and Hopf, respectively. The way of system collapse is completely different after the occurrence of these two bifurcations. In the first one, the system's state moves suddenly and steeply toward instability, whereas in the second one, the instability happens in a swinging pattern. The first parameter, which results in the presence Fig. 11 State variables i d , i q , γ 1 , γ 2 , and γ 3 and DC linked voltage (V dc ) Fig. 12 Signals of V sq and V sd , error signals (u v , e v ), and output power of PV arrays (P pv ) and i dref of one or two eigenvalues on the imaginary axis, is registered and considered as the bifurcation parameter. Finally, after system simulation in the time domain and extraction of operating point, the bifurcation theory is applied to the photovoltaic system. The effects of changing the system's parameters, both controlled and uncontrolled, is then investigated during the occurrence of different types of bifurcations in the system. Specifying the stability limits of controlled parameters can prove to be important in the design and execution stages.