In order to extract a significant amount of power from solar panels, considerable space is required (Tripathi et al., 2021), and hence replacing a conventional generator with a single solar farm may not be feasible. As a consequence, solar energy is captured at numerous locations in a given region that receives an equivalent amount of solar radiation.

The modeling of intermittency associated with distributed power generation is a major challenge in analyzing the impact on power system reliability. In contrast to the conventional generating units, which are usually represented by a two-probability state (up state and down state) model, solar energy is associated with multiple probability states (Ostovar et al., 2021) due to varying levels of irradiance at different intervals of time. The multiple probability states of solar energy are approximated using standard probability density functions (PDFs). However, the pattern of solar irradiance varies with changes in season, and hence it is difficult to determine the PDF to accurately represent them (Sreenivasulu et al., 2021). Therefore, an accurate generalized probabilistic model to represent the uncertainty associated with solar power needs a careful approach. The interaction of the probabilistic solar generation model with the discrete stochastic load model (DSLM) is another important aspect for the computation of the reliability of the electric grids (Vardi, Zahavi, and Avi-Itzhak,1977). Conventionally, the interaction of the generation model with the load model is represented by the convolution process, which becomes complex and time-consuming with an increase in the number of generating units and multistate load model. The computations can be drastically reduced if the interactions of the models are analyzed in the frequency domain (Proakis and Manolakis,2021).

The motivation of this research is to model the uncertainty and intermittency associated with solar irradiance, which can be easily combined with load uncertainty and outage of conventional generating units to obtain the reliability of the system with reduced computations. The impact of the inclusion of solar distributed generation (DG) on carbon emissions is also investigated to achieve carbon neutrality.

Many research studies have been proposed in the recent past for reliability modeling of solar irradiance and integration of the solar system with the conventional systems. The focus of this work was to track the maximum irradiance for obtaining maximum solar power (Yan et al., 2018; Li et al., 2021; de Vries, Loonen, and Hensen 2021), assess reliability based on power loss due to variable insolation (Pandit et al., 2021), determine the reliability of electric vehicles (Hashemi-Dezaki 2019; Talukdar et al., Deka, and Goswami 2021), and optimize the cost for reliability evaluation (Ram, Bhandari, and Kumar 2021). However, the major challenge for the present and upcoming grid-connected solar DG system is to adequately supply the load demand (Okundamiya 2021), which needs to be addressed. A significant amount of work has been reported for adequacy assessment of wind farm integrated power system (Ak et al., 2018; Zhu and Zhang, 2018); however, limited work has been reported on adequacy of solar DG considering the intermittency of solar irradiation. The power generated by a PV cell depends on the intensity of solar irradiation, which further depends on the location of the solar unit. The uncertainty of solar generation is represented by probability density functions, which vary with the change in seasons of the year (Vale 2015). Normal distribution function (Kim et al., 2020) and Weibull distribution (Afzaal et al., 2020) function have been used to represent the global solar irradiation data. Hence, a common simplified model needs to be developed to represent solar data to integrate the solar DG with the conventional reliability assessment schemes. The adequacy assessment schemes also involve multiple iterations, which need to be further addressed.

Reliability analysis of solar photovoltaic systems has been reported in considering the components of the solar PV systems (Pradeep Kumar and Fernandes 2017; Gautam and Kaushika 2002). The methods take into account the solar PV arrays and the grid-connected converter, but the solar irradiation intermittency is not taken into account. Although intermittency of solar irradiation has been taken into account by Yin, Molini, and Porporato (2020), the research does not consider a generalized model for it. The topological variations for solar irradiation have been taken into account for reliability studies in Carpio (2021); Almaktar, Elbreki, and Shaaban (2021), but the methods for reliability evaluation are complex. The consideration of irradiation uncertainty having multiple states in generation makes the reliability evaluation process complex. The FFT-based algorithm used in Lakshmi et al. (1995) is imperative to make the computations simpler. Also, carbon neutrality with the inclusion of solar DG needs to be taken into account (Kumar, Mohanta, and Reddy 2015; Wyrwa et al., 2022).

Solar PV has been globally deployed; hence, the solar irradiation uncertainty modeling in a generalized way can facilitate significantly toward impact assessment on reliability and carbon neutrality. The representation of solar irradiation using different PDFs at different locations is cumbersome, and hence a generalized model for the representation of the probabilities of solar irradiation has been developed. The model is also consistent with topological variations so that interactions of a multistate solar power model with the stochastic load model as well as with the random outage probabilities of conventional generating units are blended suitably using the frequency domain approach for computation. The algorithm used for generation planning considering carbon neutrality is computationally efficient with lesser memory requirements. In a nutshell, the unique contributions of the present research as compared to the existing literature are as follows:

1) Development of a generalized probabilistic model for the representation of solar intermittency taking topological variations into account.

The remainder of the study is organized as follows: Section 2 formulates a probabilistic model for solar power uncertainty; the section proposes a generalized model for solar irradiation based on Weibull distribution. Section 3 discusses the discrete production simulation considering uncertainty in solar DG for reliability analysis of smart grids. Section 4 analyzes the impact of solar generation considered in the present work on carbon neutrality and Section 5 provides a discussion on the results of the proposed methodology.

It is vital to get the data in the time-series format in order to describe the solar data using the Weibull distribution. The present study utilizes the data from various locations in California, the United States, because that state receives the maximum solar irradiation all year round. Two locations, namely, farm 1 and farm 2 were chosen to represent a distributed generation. Hourly data of each location are considered for a duration of 2 years (2018–2019) from 7:00 a.m. to 5:00 p.m. The data points can be represented by Eq. 1. Here ,   I corresponds to irradiance data, ( i = 1,2 , … , n ) corresponds to the number of days, and ( j = 1,2 , … , m ) corresponds to the hour taken into consideration. The sample irradiation values of farm 1 are presented in Table 1 for a duration of 2 years, where the duration of (7:00 a.m. to 8:00 a.m.) is represented by first hour, duration (8:00 a.m. to 9:00 a.m.) is represented by second hour, …., and (4:00 p.m. to 5:00 p.m.) is represented by 10th hour. I = [ d 11 d 21 ⋯ d m 1 d 12 d 22 ⋯ d m 2 ⋮ d i j ⋱ ⋮ d 1 n ⋯ ⋯ d m n ] (1) f ( I ) = α β ( I β ) α − 1 ⁡ exp ( − I β ) α (2) x r = β [ ln ( 1 1 − r ) ] 1 α (3)

The data given by Eq. 1 are a time-series data for a duration of 10 h each day, and hence the data of the duration considered (1 year) can be approximated to Weibull distribution f ( I ) as given by Eq. 2 (Phan and McCool 2009). Here, α = ( s t d / I ¯ ) corresponds to the shape parameter of the distribution, β   = I ¯ / Γ ( 1 + 1 α ) corresponds to the scale parameter of the distribution, I ¯ corresponds to the mean of the solar irradiation data, s t d is the standard deviation of the irradiation data, and Γ is the gamma function. The PDF plots for farm 1 and farm 2 considering 2 years of data are shown in Figure 1. The r t h quantile of Weibull distribution in terms of the shape and scale parameters of Weibull distribution is given by Eq. 3. Here, r t h quantile represents an irradiation value equal to or less than x r , and the x r values for different quantiles are presented in Table 2.

The Weibull parameters of different farms have been tabulated in Table 4. It can be observed from Table 2 that the quantile which includes the maximum solar irradiation value for farm 1 is 93.02%, and the quantile which includes the maximum solar irradiation value for farm 2 is 89.13%. Hence, it can be inferred that the solar data considered spans about 90% of the distribution, thus making Weibull distribution suitable for solar irradiation. The Weibull parameter fluctuates with change in location, even within the same region, making it challenging to determine an appropriate model to represent solar irradiation. A generalized model for solar irradiance would be useful for analyzing solar data at various locations in practical applications. The approach is useful in places where there is a shortage of historical data. The purpose of the framework is to accurately determine the output power of solar DGs installed in various geographic regions. The PDF plot of two solar farms is shown in Figure 1 in terms of the two-parameter Weibull distribution. This section determines the probability distribution of solar irradiance in terms of variance, which is obtained from the expected values of the PDF plot. The combined PDF plot for farm 1 and farm 2 is shown in Figure 2, which is obtained by taking the average of the probabilities of the two wind farms. The standard deviation S ( X ) =     E ( X i − E ( X ) ) 2 of the PDF can be obtained from the expected value E ( X ) =   ∑ P i X i of the combined PDF. Here, P i corresponds to the probability of X i (irradiance value) in the combined PDF plot of the solar irradiation at different sites.

Based on the standard deviation S ( X ) obtained, the generalized model considers solar irradiation values up to ( 5 σ ) , which accounts for very high irradiance despite their low probability. The distribution is divided into N D number of steps with a step size of ( 5 σ / N D ) , and the midpoint of each step M P S i ( i = 1,2,3 … ,   N D ) is given by (4). For example, if we consider a 50-step model, then the size of each step is 0.1 σ , and the midpoint of steps M P S i are given as: 0.05 σ ,   0.15 σ ,   0.25 σ , … ,   4.95 σ . If the number of 24-h solar irradiation values obtained for a given time interval is N y s , and N b s i is the number of values in step i , then the probability P b i =   N b s i / N y s of step i . M P S i = 5 σ 2 N D + 5 σ N D ( i − 1 ) (4)

The generalized solar irradiation model can be used to obtain the probability of a particular site if the mean ( I ¯ ) and standard deviation ( s t d ) of the irradiance of the data at the site are known.

The conventional generating unit is represented by a two-state model, with the states representing zero power generation or “down state” and rated capacity generation or “up state”. However, due to solar power intermittency caused by variation of solar irradiation due to movement of the sun, it is pragmatic to represent the solar power with a multistate model. The multistate modeling helps to incorporate the random variation of solar energy generation correlating with random irradiance. Furthermore, the generation reliability analysis requires the computation of power generated at a particular site with a solar farm so as to check whether it meets the load demand. The amount of power generated from a solar farm location varies with the amount of solar irradiation. The solar power curve, which is a correlation between solar irradiation and generated power as shown in Figure 3, can be used to calculate the power generated at a location. The power varies non-linearly when the solar irradiation is less than a predefined set radiation point R s p (usually considered 150 W/m²). The power varies linearly when the irradiation lies between R s p and solar irradiation at a predefined standard value in a given environment R p s e (usually taken as 1000 W/m²). The generated power is the rated value P r beyond R p s e . This has been explained in Section 2.2.

Figure 4 depicts the process of constructing a standard solar power model. The model requires solar data from a place with similar mean and standard deviation solar irradiation statistics. To develop a site-specific solar irradiation model, the annual mean and standard deviation of solar irradiation for that location are needed. Finally, solar power generated at a specific site for various irradiation values is calculated using solar power curve data.

The hourly solar power generated P s h i ( i = 1,2,3 … , N D ) at a particular site corresponding to solar irradiation value M P S i ( i = 1,2,3 … , N D ) is given by Eq. 4. Table 3 lists the values of solar irradiation and equivalent output power for a solar DG considering a rated power of 1500 MW and σ = 229.6488 for a 30-step model. P s h i =   { P r × M P S i 2 R p s e × R s p ;   0 ≤ M P S i < R s p P r × M P S i R p s e ; R s p ≤ M P S i < R p s e P r ; M P S i ≥ R p s e (5)

The power generated is zero for the duration for which solar irradiation is not present, which has not been taken into account for the power calculation state. The probability P p r z for zero power output is obtained from Eq. 6, where N z is the number of data points for zero solar irradiation. Similarly, rated power output is obtained beyond R p s e . Therefore, the probability of step bands beyond it can be combined to obtain the probability for rated power output P p r r given by Eq. 7. P p r z = N z N s y (6) P p r r = ∑ P b i ,   f o r   M P S i ≥ R p s e (7)

The outage of the conventional unit is represented by a two-state model. If the outage of each generating unit G 1 , with an outage capacity U 1 , is treated as a single event, then the probability of failure state is represented by forced outage rate (FOR) of the generator having the value " q " . The probability of the generator being in the normal state is by " p = ( 1 − q ) " . The DSLM is constructed by convolving all the generating units one at a time, with the order of convolution based on the priority of operation. The DSLM after consideration of outage of a single unit is pictorially represented in Figure 5; here, D S L M ( 0 ) ( l ) is the original load duration curve, and D S L M ( 1 ) ( l ) is the DSLM after consideration of outage of one generating unit G 1 , with an outage capacity of U 1 . The outage of the unit can be mathematically represented by Eq. 9, which involves the convolution of DSLM with the outage probability of the unit. Here, p 1 is the probability of G 1 being in a normal state, and q 1 is the probability of outage of G 1 . In general, the outage of unit G i ( i = 1,2,3 , … ,   N G ) with capacities U i ( i = 1,2,3 , … , N G ) can be represented by (9), which is a convolution of DSLM and outage of a generator. Here, p i is the probability of G i being in a normal state, and q i is the probability of outage of G i . If the maximum capacity of the generating units is U G T , then the final DSLM has a maximum load of ( l m + U G T ) . The present work considers DSLM and probabilistic DSLM synonymously. D S L M i ( l ) = D S L M i − 1 ( l ) ⊗ G i = p i × D S L M i − 1 ( l ) + q i × D S L M i − 1 ( l − U i ) (9)

The uncertainty of solar power generation is represented by multiple power states, necessitating the use of the multistate algorithm. If the uncertainty of the solar power is represented by N κ states, then the probability of P s h i ( i = 1,2,3 , … , N κ ) solar DG states are p r s i ( i = 1,2,3 , … , N κ ) . The power states of the solar DG, where P s h ν is the rated capacity of the solar DG, and the probabilities of all the power states of the solar DG must satisfy ∑ p r s i = 1

If the outage of N G conventional units has been taken into account by convolving the initial D L S M with all the generating units, then the final D S L M after convolution is D S L M N G . This obtained D S L M needs to be convolved with the stochastic power generating states of the solar unit. If the solar unit S G j is added to the system, then for a generating state P s h i having probability p r s i is mathematically represented by the shifting of DSLM by P s h i to the right. The final DSLM is obtained by convolving all the generating states of solar DG with the previously obtained DSLM and is given by Eq. 10. D S L M j ( l ) = ∑ i = 1 N κ [ p r s i . D S L M j − 1 ( l − P s h i ) ] (10)

Equation 10 reduces to Eq. 9 if N κ = 2 , and for j = 1 , D L S M j − 1 is taken as the final DSLM obtained by convolving all the conventional generating units with initial DSLM. The pictorial representation of convolution of all the conventional generating units to obtain D S L M N G and its convolution with probabilistic solar DG power model is shown in Figure 6.

Convolution of generating units with initial DSLM in the time domain necessitates several iterations, and the number of data points after convolution rises. This increases the memory requirements of the processor, and massive computations make the processor slow. If the time domain signals are converted to the frequency domain, the large number of addition and multiplications reduces to point by point data multiplication (Proakis and Manolakis 2021). The data are converted from the time domain to the frequency domain using a fast Fourier transform (FFT) algorithm. FFT requires the initial DSLM data and the data related to the uncertainty of the generating units to be sampled, and the number of sampling points must satisfy N s a m p = 2 c , where C is an integer. A discrete time-domain signal s ( z ) can be represented by " N s a m p " impulses scaled by a factor " a m ( m = 1,2,3 , … , N s a m p ) " and uniformly shifted by a factor " Δ z " , given by (11). The signal is converted to the frequency domain using (12), where W = e j 2 π / N s a m p , n = ( 1,2,3 , … , N s a m p ) , and m = ( 1,2,3 , … , N s a m p ) . The signal is converted back to the frequency domain using inverse FFT (IFFT) given by Eq. 13. The number of sampling points is chosen based on N s a m p ≥ T 0 T    (21) (Lakshmi et al., 1995), where T 0 = l m + ( U N G + P s h ν ) to consider sampling based on maximum loading, and the sampling interval T is chosen as 1 MW. When the number of generating units and states increases, the computational time and memory requirements for the frequency domain drop significantly. s ( z ) = ∑ k = 0 N s a m p − 1 a k × δ ( z − k . Δ z ) (11) S ( n ) = ∑ k = 0 N s a m p − 1 a k W − n k (12) y k = ∑ n = 0 N s a m p − 1 S ( n ) W n k s ( n ) = ∑ k = 0 N s a m p − 1 y k × δ ( z − Δ z ) } (13)

The DSLM obtained after all the units have been convolved with the initial DSLM shifts the initial DSLM based on the maximum generating capacity of the system. If the system has " λ " number of generating units with a total generating capacity ( U t ) , then the final DSLM after convolving all the generating units is D S L M λ ( l ) . The peak load in the final DSLM is equal to ( l m a x + U t ) ; the reliability indices, namely, loss of load probability (LOLP) and expected energy not served (EENS), which can be obtained from the final obtained DSLM, are given by Eqs 14, 15. The proposed methodology first requires solar data to formulate a generalized solar generation model so that a multistate generation model for any site can be obtained using mean and standard deviation. The obtained data are convolved with the DSLM obtained by convolving the conventional generating units with the initial DSLM. Finally, the reliability indices are obtained to determine the reliability of the system. The proposed methodology is shown in Figure 7. L O L P = D S L M λ ( l max + U t ) (14) E E N S = ∑ U t l max + U t D S L M λ ( l ) (15)