The influences of including solar photovoltaic system on distribution network operation using one of Bahir Dar radial distribution network as a case study based on selected parameters

Abstract

There are considerable power losses as well as under voltage issues in a distribution network. Utilities can improve this by including solar PV installations into their distribution networks. However, the integration of solar PV systems could have an influence on network performance. The BATA, one of Bahir Dar’s 34-bus radial distribution feeders, is utilized as a case study. This region possesses ample solar resources suitable for the deployment of solar energy systems. The PSO is employed to ascertain the precise location and scale of the PV system. The MATLAB software is utilized to apply the backward/forward load flow algorithm for calculating the voltage profile and power losses at each bus. Bus 34 and 744 kW are identified to be the optimal PV site and size respectively. The study evaluates the performance of the distribution network with and without the integration of solar systems, considering parameters such as voltage profile, power losses, and harmonic distortion. After integrating the optimal solar PV system, the minimum voltage profile rises from 0.8392 p.u. to 0.954 p.u, active and reactive power losses decrease by 63.27% and 66.723%, respectively and compliance with IEEE standards is upheld for both total voltage harmonic distortion and total current demand distortion. On the contrary, excessive-integration of solar PV systems led to abundance harmonic distortion, overvoltage, and high-power losses. The integration of PV is economically viable, with a payback period of 3.85 years. Overall, Optimizing the integration of the solar PV system is crucial for a substantial improvement in network performance.

Keywords:

• Power loss
• voltage profile
• harmonic distortion
• economic analysis

REVIEWING EDITOR:

• Qingsong Ai, Wuhan University of Technology, China

SUBJECTS:

• Electrical & Electronic Engineering
• Power & Energy
• Technology

1. Introduction

An electrical power system consists of a generation system, a transmission system, and a distribution system. A power generation station transforms fuel energy into electricity, the transmission system links the generation stations and distribution substations, and the distribution system delivers power to consumers. Radial structures are utilized in distribution systems to achieve operational simplicity. Through an interconnected transmission network, the primary distribution substation receives electricity from the generating stations (Kaur & Singh, 2017). The Radial Distribution System (RDS) network operates passively and conveys power from the substation to consumers. Consequently, in RDS, the power flow occurs in a single direction.

The current distribution system is incapable of meeting the desired power demand in terms of both quality and reliability. This challenge can be minimized by using non-traditional renewable energy supplies such as solar and wind, which are clean and safe for the environment, while employing efficient technologies (Kadir et al., 2010; Sahoo & Kulkarni, 2019). The optimal integration of renewable energy resources into the distribution network can reduce power loss and improve voltage profile (Eftekharnejad et al., 2013; Kadir et al., 2010; Sahoo & Kulkarni, 2019; Singh et al., 2019). Because the sun is a powerful unlimited resource that can readily generate electric energy via the photoelectric effect, photovoltaic (PV) generation is the best and has a high potential among alternative renewable generators (Han et al., 2017). The integration of a solar PV system into the distribution network increases network performance (Abdel-Salam et al., 2016; Almeida et al., 2019; Devi & Geethanjali, 2014; Elnozahy, 2013).

The primary elements that determine the performance of solar PV systems are the solar insolation level, ambient temperature, and sunshine hour (Pawar & History, 2019). A large number of PV cells are connected in series and parallel to produce a high power output (Xiong et al., 2020). The modules are also linked together to provide the necessary level of power (Nagaraja & Manohar, 2016). The normal test parameters are 1000 W/m², an AM 1.5 sun spectrum, and a module temperature of 25 °C (Chandrasiri, 2019; Nagaraja & Manohar, 2016).

The extent of penetration of solar PV systems influences network performance. When the solar PV system exceeds the load power demand, the excess power is fed back into the grid (Alboaouh & Mohagheghi, 2020; Ibraheem et al., 2014). When the solar power produced is insufficient to meet load demand, the grid supplies power to the load. The voltage profile improves and the power loss decreases up to a certain level of penetration, but raising the PV production beyond the maximum penetration capacity causes poor system functioning and large system power losses (Al Momani et al., 2017). Harmonics are formed when a PV system is connected to a distribution network due to the inverter’s conversion of DC to AC. Harmonic distortion has a negative impact on network performance, causing various electrical devices to fail (Xiong et al., 2020). Equation (1)(1) $$\text{PV Penetration Level}\left(\mathrm{\%}\right)=\frac{\text{Total generated PV power}\left(\mathrm{\text{MW}}\right)}{\text{Total system power}\left(\mathrm{\text{MW}}\right)}*100\%$$(1) can be used to calculate the PV penetration level (Eftekharnejad et al., 2013; Sharew et al., 2021). (1) $$\text{PV Penetration Level}\left(\mathrm{\%}\right)=\frac{\text{Total generated PV power}\left(\mathrm{\text{MW}}\right)}{\text{Total system power}\left(\mathrm{\text{MW}}\right)}*100\%$$(1)

Several scientists debated various types of literature about the effects of solar photovoltaic system integration on distribution networks a few decades ago.

Solanki et al. (2012) proposed utilising OpenDSS to analyse high penetration solar PV on a distribution system at various PV penetration levels (0%, 10%, 30%, and 50%). It has been discovered that the impact of PV integration on the distribution system in terms of profile and losses is proportional to the amount of solar PV generation. Mulenga (2015) investigated the effects of solar PV grid connection at various PV penetration levels. The research was carried out without PV and with various amounts of PV penetration (0%, 30%, 60%, and 90%) utilizing NEPLAN, Paladin Design base, and PVSyst software. After integrating a solar system, the voltage profile rises. Al Momani et al. (2017) proposed photovoltaic performance in terms of power loss. DIGSILENT power factory software assesses the influence of various penetration level PV systems on radial feeders in terms of power loss. The power loss of the system is reduced; however, this is dependent on the extent of PV penetration. Balamurugan et al. (2012) proposed utilizing power loss as a measure to assess the performance of distributed solar PV generators connected to a distribution network. It has been discovered that raising the penetration level of solar PV systems reduces active and reactive power losses up to the maximum achievable capacity of PV integration. Farhoodnea et al. (2012) highlighted the influence of highly penetrated solar PV in distribution systems. The infusion of high active power from the PV causes the voltage profile to rise above the standard limit.

Shafiullah (2015) addressed the harmonic impact of solar PV integration with the distribution network. PV values of 0%, 50%, and 100% are used. The harmonic distortion grows when solar PV injects more power into the distribution system. Sahoo and Kulkarni (2019) argued that PV integration has an influence on harmonics of the distribution networks. The percentage of photovoltaics (10–60%) is considered. The harmonic analysis also revealed that the harmonics were within the acceptable operating range up to 38.2% PV penetration.

The study examined the impact of various levels of solar PV integration on the distribution system and determined that the integration of PV systems has either a positive or negative impact depending on the level of solar PV penetration. For this article, the various levels of photovoltaic system are integrated with the selected Bahir Dar distribution feeder called BATA to examine its impact based on the selected parameters. In addition, solar energy resource potential assessment and economic feasibility of solar photovoltaic system integration is also included in this study.

2. The viability study for solar energy potential of the area

Before proceeding with photovoltaic power installation for both grid-connected and independent PV systems, it is critical to determine whether the solar resource is practical in the specific site under consideration (Bawazir & Cetin, 2020). Photovoltaic power generation is a method of generating electricity by converting radiant light energy into electricity using semiconductor materials. A solar photovoltaic (PV) system generates electricity by utilizing the hours of sunlight. The sunlight hour is the duration of time that the sun shines in a specific region at a specific time. Most design issues are regarded to be worth 4–6 h each day (Ibrahim et al., 2016). It is mostly presented as an average of several years. Table 1 shows the sunshine hour statistics for Bahir Dar City for ten consecutive years, as recorded by the metrology department.

Table 1. Monthly average sunshine hour data.

	Region; Gojjam Station: Bahir Dar
	Element: Monthly sunshine hour
Month	Jan	Feb	Mar	Apr	May	Jun	July	Aug	Sept	Oct	Nov	Dec
2010	9.2	10.0	8.7	8.6	8.7	7.1	4.8	5.1	6.2	8.8	10	10.3
2011	9.8	10.0	8.8	8.8	7.7	7.0	4.7	4.1	5.7	8.8	9.6	9.6
2012	9.2	8.6	9.6	8.9	8.2	6.6	4.2	4.5	5.8	8.8	9.3	10.3
2013	9.0	9.8	10.2	10.7	8.0	8.0	5.1	4.2	5.8	8.5	9.7	9.4
2014	8.0	9.0	9.2	9.0	7.2	7.6	4.1	5.1	7.7	8.2	9.0	8.5
2015	8.7	8.3	8.6	8.8	8.4	7.2	4.5	4.8	6.4	9.0	9.5	9.3
2016	9.2	10.5	8.0	9.9	8.2	7.1	5.2	4.3	6.1	9.1	9.3	9.9
2017	10.1	10.1	9.2	9.4	8.2	5.6	4.5	4.0	5.9	9.5	8.6	9.9
2018	9.7	10.1	9.3	9.8	8.6	6.6	4.2	4.3	6.8	8.3	9.7	9.9
2019	9.8	9.8	9.1	7.9	7.7	6.9	4.9	4.7	6.0	8.6	9.7	9.7

According to Table 1, the monthly average sunlight hour of Bahir Dar City meets the minimum sunshine hour requirement for developing a solar power system. It indicated that the area has an abundance of untapped solar resources. The radiation of the area is then estimated using the recorded average duration of the sunshine hour (Sahoo & Kulkarni, 2019; Sharew et al., 2021). The monthly average daily solar radiation of the site is calculated using the Angstrom – Prescott estimation model (Podder et al., 2014). (2) $${\mathrm{ }H}_o=\hspace{0.17em}\frac{24}{\pi }{G}_{\mathit{\text{SC}}}\left(1+\mathrm{\text{cos}}\mathrm{ }\frac{360\mathit{\text{nd}}}{365}\right)\mathrm{*}(\mathrm{\text{cos}}\phi \mathrm{ }\mathrm{\text{cos}}\delta \mathrm{ }\mathrm{\text{sin}}{\omega }_s+\frac{{\pi \omega }_s}{{180}^0}\mathrm{ }\mathrm{\text{sin}}\phi \mathrm{ }\mathrm{\text{sin}}\delta )$$(2) Where, $$\delta ={23.45}^o\hspace{0.17em}\mathrm{\text{sin}}\hspace{0.17em}\left({360}^o\frac{284+n_d}{365}\right)$$

Cos ${\omega }_s$ = −tan $\phi$ tanδ $$H=H_{0(a+\frac{\mathit{\text{bn}}}{N})}$$ where $G_{\mathit{\text{sc}}}$ is the solar constant = 1.376 kW/m²; $\phi$ is the latitude in degree; $\delta$ is the solar declination in degree; ${\omega }_s$ is the sunset hour angle in degree; $\mathit{\text{nd}}$ is a day of the year from January 1 to December 31 taking January 1st as 1; $N$ is the monthly average of the maximum possible hours of sunshine; n is the monthly average daily hours of sunshine; $H_o$ is the monthly average daily extraterrestrial solar radiation; H is the monthly average daily global solar radiation; a and b are empirical coefficients.

The computed monthly average radiation of the site is displayed in Table 2. As shown in Table 2, the city’s solar radiation is found to be within permissible limits (4–6 kWh/m²/day), implying that the area is rich in solar resources.

Table 2. The average solar irradiance of Bahir Dar City.

Months	Jan	Feb	Mar	Apr	May	Jun	Jul	Aug	Sep	Oct	Nov	Dec
Solar irradiance(kWh/m^2/day)	4.76	5.2	5.5	5.6	5.26	4.74	4.02	4.04	4.67	4.87	5.01	4.67

3. Algorithm of particle swarm optimization

Optimization is defined as the process of either minimizing or maximizing tasks. As minimizing any function is equivalent to maximizing its additive inverse, the terms minimization and optimization can be applied interchangeably. Due to this rationale, optimization has become highly significant in numerous fields (Akula, 2013). Kennedy and Eberhart introduced a solution to the nonlinear and intricate optimization problem by studying the behavior of a flock of birds. They formulated the idea of optimizing a function by employing a swarm of particles. To address optimization challenges, the PSO algorithm draws inspiration from the behavior of animals. In the context of PSO, ‘swarm’ refers to the population, and each ‘particle’ represents an individual member of that population. Every particle explores the entire space by moving randomly in various directions. It retains the best solutions encountered thus far and keeps track of the positions of neighboring particles. In this manner, all particles in the swarm strive to advance towards improved positions until the swarm attains an optimal solution (Akula, 2013). This research employs Particle Swarm Optimization (PSO) to identify the optimal location and size of the PV system. PSO is a method for acquiring optimal values of essential parameters in a system by utilizing swarms of insects. This approach begins by initializing particles as a proposed solution to a specified problem and trains them to explore for an optimal location (Akula, 2013; Thani et al., 2017). PSO has the capability to yield results within a brief execution time and iterations, rendering it highly favorable (Akula, 2013; Devi & Geethanjali, 2014; Sharew et al., 2021). The positions of the particles are modified, as illustrated in the diagram depicted in Figure 1 (Alrashidi & Alhajri, 2010; Kefale et al., 2021; Shafiullah, 2015).

Figure 1. Graphical expression of PSO (Sharew et al., 2021).

X_i-The current position of the given particle

V_i-The current velocity of the given particle

The distance between the personal best (P_best) and current position of the given particle

The distance between the global best (P_Gbest) and current position of the given particle (3) $$V_{\mathit{\text{ij}}}\left(t+1\right)=w\mathrm{*}{v}_{\mathit{\text{ij}}}\left(t\right)+r_1{c}_1\left(p_{\mathit{\text{ij}}}\right(t)-x_{\mathit{\text{ij}}}(t\left)\right)+r_2{c}_2\left(g\right(t)-x_{\mathit{\text{ij}}}(t\left)\right)$$(3) (4) $$x_{\mathit{\text{ij}}}\left(t+1\right)=v_{\mathit{\text{ij}}}\left(t+1\right)+x_{\mathit{\text{ij}}}\left(t\right)$$(4) where w is inertia weight; $v_{\mathit{\text{ij}}}\left(t\right)$ is the velocity of agent at time t; $x_{\mathit{\text{ij}}}\left(t\right)$ is the position of agent at time t; $p_{\mathit{\text{ij}}}\left(t\right)$ is the personal best position of agent starting from initialization through time t; $g\left(t\right)$ is the global best position of agent starting from initialization through time t; c1 and c2 are positive acceleration constants utilized to ascertain the contribution levels of the cognitive and social components, respectively; $r_1$ and $r_2$ are random numbers between 0 and 1.

The update velocity has three parts as shown in Equation (3)(3) $$V_{\mathit{\text{ij}}}\left(t+1\right)=w\mathrm{*}{v}_{\mathit{\text{ij}}}\left(t\right)+r_1{c}_1\left(p_{\mathit{\text{ij}}}\right(t)-x_{\mathit{\text{ij}}}(t\left)\right)+r_2{c}_2\left(g\right(t)-x_{\mathit{\text{ij}}}(t\left)\right)$$(3) . These are

• The inertia part: $w*v_{\mathit{\text{ij}}}\left(t\right)$

• The cognitive part: $r_1{c}_1(p_{\mathit{\text{ij}}}(t)-x_{\mathit{\text{ij}}}(t))$

• The social part: $r_2{c}_2(g(t)-x_{\mathit{\text{ij}}}(t))$

The previous particle position is modified by incorporating the new velocity vector. As illustrated in Figure 2, the iterative process persists until a desired point is achieved.

Figure 2. PSO algorithm flowchart.

Benefits of the Particle Swarm Optimization (PSO) algorithm (Akula, 2013):

• This algorithm finds application in both scientific research and engineering challenges due to its straightforward implementation

• In contrast to other optimization techniques, this algorithm has a minimal impact on the parameters affecting the optimal solution, primarily because it involves fewer parameters.

• It requires simple analysis

• The optimal value can be quickly acquired within a short period

• In comparison to other optimization methods, this algorithm relies less on the selection of initial values

• A straightforward idea is at play here

4. Modelling of the distribution system

Bahir Dar serves as the capital city of the Amhara region and is situated in the northern hemisphere at a latitude of 11.5936403 and a longitude of 37.39077. The urban area is experiencing swift growth in its industrial and commercial domains. The power demand is met through two substations situated within the city of Bahir Dar. The power demands are primarily met and linked through the Tiss Abay I, Tiss Abay II, Beles, and Fincha generation stations, with outgoing 230 kV connections to the Alamata and Gonder-Metema substations. There is a 400/230 kV substation (substation II), two 230/132/15 kV and 230/66/15 kV substations (substation II), and one 66/45/15 kV substation (substation I) transformer, all of which provide power to the town. The system comprises eleven radial feeders, with seven (Adet, Tis Abay, Ghion, Papyrus, Industry, Bata, and Airforce) deriving from substation II, and the remaining four (Bete-mengist, Gambi, Hamusit, and Boiler) from substation I. For this case study, the Bata feeder is chosen among the seven 15 kV outgoing feeders from Bahir Dar substation II due to significant power interruptions in the area. It is necessary to create a model of the distribution network configuration for the purpose of performing load flow analysis and determining the necessary parameter values. The Viso software is employed for drawing the selected radial distribution system in Bahir Dar. Consequently, the radial distribution feeder structure named Bata, as illustrated in Figure 3, was selected for the study. It consists of 34 buses, where Bus-1 functions as the reference bus, four act as common connection buses, and the remaining 29 serve as load buses.

Figure 3. Single-line schematic depiction of the Bata feeder (Sharew et al., 2021).

The overhead distribution line in the case study area spans a length of 18.77 kilometers. Every one of the 34 nodes maintains a voltage level of 15 kV. The feeder is supplied with 1.86 megawatts of active power and 1.25 megavolt-amperes of reactive power (Sharew et al., 2021). The impedance of the overhead transmission line can be computed based on the manufacturer’s resistance, specified length, and total line length. Distribution feeders make use of stranded conductors. The Bata distribution feeder utilizes stranded conductors of AAC-25, AAC-50, and AAC-95 types. The impedance for the positive sequence of line segments is determined by multiplying the impedance per kilometer by the segment length. In the case of short transmission lines, the shunt admittance of overhead wires is neglected. The line data for the network is calculated based on the information provided in Table 3. The resistance per kilometer for each conductor is extracted from the manufacturer’s data sheet.

Table 3. Overhead medium voltage conductor size.

Conductor type	Nominal area (mm^2)	Resistance (R0) (Ω/km)	Reactance (X0) (Ω/km)
AAC	25	1.1810	0.8345
AAC	50	0.5787	0.4362
AAC	95	0.3085	0.2147

The values of RL (calculated resistance) and XL (calculated reactance) can be found by utilizing R0 (ohm/km) and X0 (ohm/km) for each type of distribution line branch in the following manner: (5) $${\mathrm{R}}_{\mathrm{L}}=\left({\mathrm{R}}_0\times \text{length of the line branch}\right)$$(5) (6) $${\mathrm{X}}_{\mathrm{L}}=({\mathrm{X}}_0\times \text{length of the line branch}$$(6)

The per unit value of the parameters should be determined. Per unit is the ratio of actual value of any quantity to the base or reference value of the same quantity with identical unit. It is important for analysis of distribution system by converting the parameters, like voltage, current, power, and impedance with the chosen base values for ease of calculation. The per unit values are dimensionless. (7) $$\mathrm{\text{Quantity}}\left(\mathrm{\text{pu}}\right)=\frac{\text{Actual quantity}\left(\text{normal unit}\right)}{\text{Base quantity}\left(\text{normal unit}\right)}$$(7)

The per unit values of active and reactive power at each bus can be calculated as follow using the selected base values: (8) $$P_{\mathit{\text{pu}}}=\frac{P}{S_{\mathit{\text{base}}}}$$(8) (9) $$Q_{\mathit{\text{pu}}}=\frac{Q}{S_{\mathit{\text{base}}}}$$(9)

The base values for both apparent power and voltage are chosen as shown below.

Sbase = 100 MVA (common power)

Vbase = 15 KV (system voltage as common)

4.2. Load profile under base case conditions

The computed values for RL, XL, and the load data of the BATA feeder are presented in Table 4 (Alrashidi & Alhajri, 2010).

Table 4. The feeder line and load data (Sharew et al., 2021).

Sending node (i)	Receiving node (j)	Conductor type	Length (km)	Resistance RL (Ω)	Reactance XL (Ω)	Receiving end load (kW)	Receiving End load (kVAr)
1	2	AAC95	1.8599	0.57378	0.39932	24.6659	12.4079
2	3	AAC95	0.2911	0.08980	0.06250	9.959	6.689
3	4	AAC95	0.3809	0.11751	0.08178	69.335	62.259
4	5	AAC95	0.6789	0.20944	0.14576	25.769	10.589
5	6	AAC95	0.2899	0.08943	0.06224	0.000	0.000
6	7	AAC95	1.1061	0.34123	0.23748	25.9649	15.85
7	8	AAC95	0.5549	0.17119	0.11914	0.000	0.00
8	9	AAC95	0.4741	0.14626	0.10179	47.99	25.799
3	10	AAC50	0.9031	0.52262	0.39393	81.699	71.109
10	11	AAC50	0.2731	0.15804	0.11913	136.739	63.839
6	12	AAC25	0.4989	0.58920	0.41633	75.448	54.219
12	13	AAC95	0.7111	0.21937	0.15267	0.000	0.00
13	14	AAC50	0.457	0.26447	0.19934	63.899	44.648
13	15	AAC50	1.2759	0.73836	0.55655	80.399	52.797
15	16	AAC25	0.1031	0.12176	0.08604	41.979	21.798
16	17	AAC25	0.7061	0.83390	0.58924	88.279	66.929
17	18	AAC50	0.3471	0.20087	0.15141	5.079	3.059
8	19	AAC50	0.2451	0.14184	0.10691	112.99	46.819
19	20	AAC95	1.1231	0.34648	0.24113	85.069	62.195
20	21	AAC25	0.5181	0.61188	0.43236	55.599	38.989
21	22	AAC50	0.5110	0.29572	0.22291	0.000	0.00
22	23	AAC50	0.4241	0.24543	0.18499	27.219	12.019
22	24	AAC50	0.2810	0.16262	0.12257	28.039	12.589
24	25	AAC50	2.1210	1.22742	0.92518	55.809	37.959
25	26	AAC25	0.1880	0.22203	0.15689	74.599	55.229
26	27	AAC25	0.2751	0.32489	0.22957	46.519	24.539
22	28	AAC25	0.6310	0.74521	0.52657	48.395	93.79
28	29	AAC50	0.3279	0.18976	0.14303	9.299	5.499
29	30	AAC50	0.3181	0.18409	0.13876	34.449	17.368
30	31	AAC50	0.1971	0.11406	0.08598	25.779	48.739
29	32	AAC50	0.2419	0.13999	0.10552	78.499	47.7968
32	33	AAC50	0.3790	0.21933	0.165320	98.00	66.00
33	34	AAC95	0.0910	0.028074	0.0195377	93.499	86.00

5. Analysis of the power flow in a Radial Distribution System

The research focuses on a distribution network of radial type, which possesses distinctive characteristics. The fundamental attributes of this feeder include its radial configuration, utilization of multi-phase conductors, operation with unbalanced loads, and a high R/X ratio (Kaur & Singh, 2017). The features of radial feeders result in conventional power flow methods such as Conventional Gauss-Seidel, Newton-Raphson, and Fast Decoupled Load Flow algorithms exhibiting suboptimal convergence properties, making them inadequate for meeting the requirements of the distribution system (Tahir et al., 2018). Therefore, it is imperative to devise a suitable power flow technique for conducting load flow calculations in radial systems. The backward/forward sweep load flow analysis stands out as the appropriate and favored technique. It is necessary to set the magnitude of voltage and the source phase angle, while providing the load demand at each bus (Kaur & Singh, 2017). During the backward sweep, Kirchhoff’s current law (KCL) and Kirchhoff’s voltage law (KVL) are employed to compute line part currents and voltages, starting from the far end of the feeder and progressing back to the source. The forward sweep is conducted to update node voltages along the feeder whenever the computed voltage deviates from the initially specified voltage beyond an acceptable threshold. The process iterates between the forward and backward sweeps until the voltage error at the source is brought within the acceptable range (Devi & Geethanjali, 2014; Kaur & Singh, 2017; Tahir et al., 2018). In this study, the load flow technique involves the application of either the forward sweep, backward sweep, or a combination of both.

5.1. Load flow analysis using forward/backward sweep

This method is formulated through the injection of bus values into the branch current matrix, conversion of branch current to the bus voltage matrix, and incorporation of the equivalent current injection. The load at each bus, denoted as Si, is expressed as (Devi & Geethanjali, 2014; Sharew et al., 2021): (10) $$S_{\mathit{\text{LI}}}=P_{\mathit{\text{Li}}}+j{Q}_{\mathit{\text{Li}}}i=1\dots N$$(10)

Step-1: Algorithm for Load Flow Analysis Using Backward Sweep

The current injection at the ith bus during the kth iteration is expressed as: (11) $$I_i{}^k=I_i^r(V_i^k)+j{I}_i^i(V_i^k)={\left(\frac{S_i}{V_i{}^k}\right)}^{\ast }={\left(\frac{P_i+j{Q}_i}{V_i{}^k}\right)}^{\ast }$$(11) where $V_i{}^k$ – The bus voltage at the ith node during the kth iteration; $I_i^r$ – The actual current injection at the ith node during the kth iteration; $I_i^i$ – The imaginary current injection of ith node at kth iteration (Sharew et al., 2021).

The matrix relationship is established utilizing the distribution model exemplified in the diagram presented in Figure 4.

Figure 4. Illustrative distribution system.

Equation (11) enables the determination of the injected currents. Applying Kirchhoff’s current law (KCL) to the distribution network allows the calculation of branch currents. Subsequently, branch currents can be represented as functions of the equivalent current injection. The formulation of branch currents B1, B2, B3, B4, and B5 is expressed as (Devi & Geethanjali, 2014). (12) $${\mathrm{B}}_1={\mathrm{I}}_{\mathrm{L}2}+{\text{I}}_{\mathrm{L}3}+{\mathrm{I}}_{\mathrm{L}4}+{\text{I}}_{\mathrm{L}5}+{\mathrm{I}}_{\mathrm{L}6}$$(12) (13) $${\mathrm{B}}_2={\mathrm{I}}_{\mathrm{L}3}+{\mathrm{I}}_{\mathrm{L}4}+{\text{I}}_{\mathrm{L}5}+{\mathrm{I}}_{\mathrm{L}6}$$(13) (14) $${\mathrm{B}}_3={\mathrm{I}}_{\mathrm{L}4}+{\text{I}}_{\mathrm{L}5}$$(14) (15) $${\mathrm{B}}_4={\text{I}}_{\mathrm{L}5}$$(15) (16) $${\mathrm{B}}_5={\mathrm{I}}_{\mathrm{L}6}$$(16)

Therefore, the relationship between the bus current injection and the branch currents can be denoted as, (17) $$\left[\begin{array}{c}{\mathrm{B}}_1\\ {\mathrm{B}}_2\\ {\mathrm{B}}_3\\ {\mathrm{B}}_4\\ {\mathrm{B}}_5\end{array}\right]=\left[\begin{array}{c}1 1 1 1 1\\ 0 1 1 1 1 \\ 0 0 1 1 0\\ 0 0 0 1 0\\ 0 0 0 0 1\end{array}\right]\left[\begin{array}{c}{\mathrm{I}}_{\mathrm{L}2}\\ {\mathrm{I}}_{\mathrm{L}3}\\ {\mathrm{I}}_{\mathrm{L}4}\\ {\mathrm{I}}_{\mathrm{L}5}\\ {\mathrm{I}}_{\mathrm{L}6}\end{array}\right]$$(17)

Equation (17)(17) $$\left[\begin{array}{c}{\mathrm{B}}_1\\ {\mathrm{B}}_2\\ {\mathrm{B}}_3\\ {\mathrm{B}}_4\\ {\mathrm{B}}_5\end{array}\right]=\left[\begin{array}{c}1 1 1 1 1\\ 0 1 1 1 1 \\ 0 0 1 1 0\\ 0 0 0 1 0\\ 0 0 0 0 1\end{array}\right]\left[\begin{array}{c}{\mathrm{I}}_{\mathrm{L}2}\\ {\mathrm{I}}_{\mathrm{L}3}\\ {\mathrm{I}}_{\mathrm{L}4}\\ {\mathrm{I}}_{\mathrm{L}5}\\ {\mathrm{I}}_{\mathrm{L}6}\end{array}\right]$$(17) can be written in the following manner: (18) $$\left[\mathrm{B}\right]=\left[\mathrm{\text{BIBC}}\right]\left[\mathrm{I}\right]$$(18) Where: BIBC- BIBC is a matrix that transforms bus injection into branch current, and it is an upper triangular matrix composed solely of zeros and ones. Step-2: Algorithm for Load Flow Analysis Using Forward Sweep

The nodal voltage vector V is iteratively updated from the source to the load, adhering to Kirchhoff’s Voltage Law (KVL). This process involves the use of the previously calculated branch current vector B and the transformation from branch current to bus voltage (BCBV). The relationship between branch current and bus voltage, as illustrated in Figure 4, can be expressed as (Devi & Geethanjali, 2014): (19) $$\begin{array}{l}{V}_2=V_1-B_1{Z}_{12}\\ V_3=V_2-B_2{Z}_{23}=V_1-B_2{Z}_{23}\\ V_4=V_3-B_3{Z}_{34}=V_1=B_1{Z}_{12}-B_2{Z}_{23}-B_3{Z}_{34}\\ V_5=V_4-B_4{Z}_{45}=V_1=B_1{Z}_{12}-B_2{Z}_{23}-B_3{Z}_{34}-B_4{Z}_{45}\\ V_6=V_3-B_5{Z}_{56}=V_1=B_1{Z}_{12}-B_2{Z}_{23}-B_5{Z}_{56}\end{array}$$(19)

Where V_i – Bus i‘s voltage, and Z_ij – The line impedance between buses i and bus j.

As per Equation (19)(19) $$\begin{array}{l}{V}_2=V_1-B_1{Z}_{12}\\ V_3=V_2-B_2{Z}_{23}=V_1-B_2{Z}_{23}\\ V_4=V_3-B_3{Z}_{34}=V_1=B_1{Z}_{12}-B_2{Z}_{23}-B_3{Z}_{34}\\ V_5=V_4-B_4{Z}_{45}=V_1=B_1{Z}_{12}-B_2{Z}_{23}-B_3{Z}_{34}-B_4{Z}_{45}\\ V_6=V_3-B_5{Z}_{56}=V_1=B_1{Z}_{12}-B_2{Z}_{23}-B_5{Z}_{56}\end{array}$$(19) , it is clear that the substation voltage, branch currents, and line characteristics can collectively represent the bus voltage. With additional buses undergoing analogous operations, the relationship between branch currents and bus voltages can be expressed as follows: (20) $$\left[\begin{array}{c}{V}_2\\ V_3\\ V_4\\ \begin{array}{c}{V}_5\\ V_6\end{array}\end{array}\right]=\left[\begin{array}{c}{V}_1\\ V_1\\ V_1\\ \begin{array}{c}{V}_1\\ V_1\end{array}\end{array}\right]-\left[\begin{array}{c}\begin{array}{cccc}{Z}_{12}& 0& 0& 0\\ Z_{12}& Z_{23}& 0& 0\\ Z_{12}& Z_{23}& Z_{34}& 0\\ Z_{12}& Z_{23}& Z_{34}& Z_{45}\end{array}\begin{array}{c}0\\ 0\\ 0\\ 0\end{array}\\ \begin{array}{cccc}{Z}_{12}& Z_{23}& 0& 0\end{array}\begin{array}{c}{Z}_{56}\hfill \end{array}\end{array}\right]\left[\begin{array}{c}{B}_1\\ B_2\\ B_3\\ \begin{array}{c}{B}_4\\ B_5\end{array}\end{array}\right]$$(20)

Where

The sample network furnishes the Branch Current to Bus Voltage (BCBV). (21) $$\begin{array}{l}\text{[BCBV ]=}\left[\begin{array}{l}\begin{array}{llll}{Z}_{12}& 0& 0& \begin{array}{ll}0& 0\end{array}\end{array}\\ \begin{array}{llll}{Z}_{12}& Z_{23}& 0& \begin{array}{ll}0& 0\end{array}\end{array}\\ \begin{array}{llll}{Z}_{12}& Z_{23}& Z_{23}& \begin{array}{ll}0& 0\end{array}\end{array}\\ \begin{array}{llll}{Z}_{12}& Z_{23}& Z_{23}& \begin{array}{ll}{Z}_{23}& 0\end{array}\end{array}\\ \begin{array}{llll}{Z}_{12}& Z_{23}& 0& \begin{array}{ll}0& Z_{23}\end{array}\end{array}\end{array}\right]\\ \end{array}$$(21)

The bus voltage at the (k + 1)th iteration can be expressed in a general form as: (22) $$\left[V^{k+1}\right]=\left[V_1\right]-\left[\mathit{\text{BCBV}}\right]\left[B\right]$$(22)

Typically, the node is represented by the letter 'i', and the iteration number is denoted by ‘k’. (23) $$\begin{array}{l}{\mathrm{I}}_{i=1,i}^k={\mathrm{I}}_i^k+{\mathrm{I}}_{i,i+1}^k\\ {\mathrm{V}}_i^k={\mathrm{V}}_{i-1}^k-Z_{i-1,i}*{\mathrm{I}}_{i-1,i}^{k-1}\end{array}$$(23)

The total active and reactive power losses in the distribution system can be mathematically represented using the following formulas: (24) $$P_{\mathit{\text{loss}}}=\sum _{i=1}^{\mathit{\text{nb}}}{\left|I\left(i\right)\right|}^{2}R\left(i\right)$$(24) (25) $$Q_{\mathit{\text{loss}}}=\sum _{i=1}^{\mathit{\text{nb}}}{\left|I\left(i\right)\right|}^{2}X\left(i\right)$$(25) Where nb is the total number of branches in the distribution system.

6. Problem formulation

The active and reactive power equations are created between the buses in the network using a single line diagram depicted in Figure 5. Mathematical expressions for the total active and reactive power losses can be derived for both scenarios, before and after the integration of solar PV.

Figure 5. A main feeder’s single-line schematic.

The following formulas can be used to determine the active and reactive power between buses i and i + 1 prior to PV system integration. (26) $$\begin{array}{l}{P}_{i+1}=P_i-P_{\mathit{\text{loss}},i}-P_{\mathit{\text{loss}},i+1}\\ =P_i-\frac{R_i}{V_i{}^2}\left(P_i{}^2+Q_i{}^2\right)-P_{\mathit{\text{loss}},i+1}\end{array}$$(26) (27) $$\begin{array}{l}{Q}_{i+1}=Q_i-Q_{\mathit{\text{loss}},i}-Q_{\mathit{\text{loss}},i+1}\\ =Q_i-\frac{X_i}{V_i{}^2}\left(P_i{}^2+Q_i{}^2\right)-Q_{\mathit{\text{loss}},i+1}\end{array}$$(27)

After that, the total power losses of the feeder, Ploss, and Qloss, can be calculated by adding the losses of each feeder line section, which is provided as: (28) $${\mathrm{ }P}_{T,\mathit{\text{loss}}}={\sum }_{i=1}^N\frac{R_i}{{V_i}^2}({P_i}^2+{Q_i}^2)$$(28) (29) $${\mathrm{ }Q}_{T,\mathit{\text{loss}}}={\sum }_{i=1}^N\frac{X_i}{{V_i}^2}\left({P_i}^2+{Q_i}^2\right)$$(29)

After PV integration, let the apparent power flow be $P_{\mathit{\text{pv}},\mathit{\text{Loss}}(i,i+1)}$ and ${ Q}_{\mathit{\text{pv}},\mathit{\text{Loss}}(i,i+1)}.$ Using the same methods as previously demonstrated, the power losses between these buses are calculated (Sharew et al., 2021). (30) $${\mathrm{ }P}_{\mathit{\text{pv}},\mathit{\text{Loss}}(i,i+1)}=\frac{{\left(P_{i-P_{\mathit{\text{pv}}}}\right)}^2+{Q_i}^2}{{V_i}^2}{R}_i\mathrm{ }$$(30) (31) $${\mathrm{ }Q}_{\mathit{\text{pv}},\mathit{\text{Loss}}(i,i+1)}=\frac{{\left(P_{i-P_{\mathit{\text{pv}}}}\right)}^2+{Q_i}^2}{{V_i}^2}{X}_i$$(31)

During solar PV system integration to the distribution system, the network must adhere to the operational range (limiting values) listed below (Sharew et al., 2021).

6.1. Voltage constraint

According to the IEEE agreement, bus voltage levels (pu) shall not exceed the limit of 0.95 pu to 1.05 pu.

6.2. Current constraint

The current in each branch should not exceed its maximum capacity. This constraint can be described as: $$I_{\mathit{\text{ij}}}\le {I_{\mathit{\text{ij}}}}^{\mathrm{\text{max}}}$$

Where ${I_{\mathit{\text{ij}}}}^{\mathrm{\text{max}}}$ and $I_i$ are maximum permissible current and current of a branch between nodes i and j, respectively

6.3. Harmonic constraint

The total and individual voltage harmonic distortion levels are restricted by the IEEE's permitted operating range. $${\mathit{\text{THD}}}_v\le {\mathit{\text{THD}}}_v^{\mathrm{\text{max}}}$$ $${\mathit{\text{IHD}}}_v\le {\mathit{\text{IHD}}}_v^{\mathrm{\text{max}}}$$

Optimization can be described as a process that seeks to identify the maximum or minimum value within a given function or procedure. The term ‘objective function’ refers to the function that we seek to either minimize or maximize. The objective function is adjusted in the manner shown in Equation (32)(32) $$\mathrm{F}=\mathrm{W}1\text{*Δ}{P}_{\mathit{\text{loss}}}^{\mathit{\text{PV}}}+W2\mathrm{*}\Delta Q_{\mathit{\text{loss}}}^{\mathit{\text{PV}}}+\mathrm{W}3\mathrm{*}\sum _{i=1}^n{(V_i-1)}^2$$(32) (Sharew et al., 2021):

Maximize (32) $$\mathrm{F}=\mathrm{W}1\text{*Δ}{P}_{\mathit{\text{loss}}}^{\mathit{\text{PV}}}+W2\mathrm{*}\Delta Q_{\mathit{\text{loss}}}^{\mathit{\text{PV}}}+\mathrm{W}3\mathrm{*}\sum _{i=1}^n{(V_i-1)}^2$$(32)

Subjected to: $$0.95\le V_i\le 1.05$$ $$P_V\le {\mathit{\text{Pv}}}^{\mathrm{\text{max}}}$$ $$Q_V\le {\mathit{\text{Qv}}}^{\mathrm{\text{max}}}$$ $$\mathrm{W}1+\mathrm{W}2+\mathrm{W}3=1$$

Where, W1, W2 and W3 are the inertia weights

Figure 6 shows how to model a distribution network with a connected solar PV installation using ETAP software.

Figure 6. PV integrated BATA feeder single-line schematic for harmonic analysis (Sharew et al., 2021).

7. Harmonic analysis

This research examines the effects of voltage harmonic distortion on the radial distribution network in conjunction with different levels of solar PV system integration, employing the ETAP software. The study explores the incorporation of solar photovoltaic systems into the distribution system, leading to the generation of harmonics that impact network performance and power quality. The generated harmonics are a direct result of the fundamental frequency. An abundance of harmonic generation within a network has the potential to induce failures in both static and rotating electrical equipment (Pawar & History, 2019; Rupa & Ganesh, 2014; Sharew et al., 2021).

Equation (33)(33) $$h=\left(n\mathrm{*}p\right)\pm 1$$(33) defines the harmonic number formula (Ciprés, 2011; Sharew et al., 2021). (33) $$h=\left(n\mathrm{*}p\right)\pm 1$$(33) where h denotes for the harmonic number; n can be any number, and p denotes the circuit’s pulse count.

Total harmonic distortion is computed by dividing the collective RMS value of the harmonic components in a signal by the fundamental value of that signal. As indicated below, Equations (34)(34) $${\mathit{\text{IHD}}}_v\left(\mathrm{\%}\right)=\frac{V_{h,\mathit{\text{rms}}}}{{V1,}_{\mathit{\text{rms}}}}\mathrm{*}100$$(34) and (35)(35) $${\mathrm{ }\mathit{\text{THD}}}_V\left(\mathrm{\%}\right)=\raisebox{1ex}{$\sqrt{\sum _{\mathrm{n}=2}^{\mathrm{\infty }}{{(\mathit{\text{Vn}},}_{\mathit{\text{rms}}})}^2}$}\!\left/ \!\raisebox{-1ex}{$\left({V1,}_{\mathit{\text{rms}}}\right)$}\right.\mathrm{*}100$$(35) can be employed to specify the individual and cumulative voltage harmonic distortions at each bus (Kadir et al., 2010; Sharew et al., 2021). (34) $${\mathit{\text{IHD}}}_v\left(\mathrm{\%}\right)=\frac{V_{h,\mathit{\text{rms}}}}{{V1,}_{\mathit{\text{rms}}}}\mathrm{*}100$$(34) (35) $${\mathrm{ }\mathit{\text{THD}}}_V\left(\mathrm{\%}\right)=\raisebox{1ex}{$\sqrt{\sum _{\mathrm{n}=2}^{\mathrm{\infty }}{{(\mathit{\text{Vn}},}_{\mathit{\text{rms}}})}^2}$}\!\left/ \!\raisebox{-1ex}{$\left({V1,}_{\mathit{\text{rms}}}\right)$}\right.\mathrm{*}100$$(35)

Equation (35)(35) $${\mathrm{ }\mathit{\text{THD}}}_V\left(\mathrm{\%}\right)=\raisebox{1ex}{$\sqrt{\sum _{\mathrm{n}=2}^{\mathrm{\infty }}{{(\mathit{\text{Vn}},}_{\mathit{\text{rms}}})}^2}$}\!\left/ \!\raisebox{-1ex}{$\left({V1,}_{\mathit{\text{rms}}}\right)$}\right.\mathrm{*}100$$(35) gives us the overall harmonic distortion of voltage; $$V_{\mathit{\text{rms}}}=V_{\mathit{\text{dc}}}$$ $${V1,}_{\mathit{\text{rms}}}=\frac{{4V}_{\mathit{\text{dc}}}}{\sqrt{2\Pi }}$$

Where, Vh,rms is the voltage at h order’s rms value Vn,rms is the voltage’s rms value V1,rms is the voltage’s root mean square (rms) at fundamental

The IEEE standards define specific operational ranges for harmonic distortion levels. Table 5 illustrates the allowable limits for voltage harmonic distortion at the Point of Common Coupling (PCC) (Minh Quan Duong, 2019).

Table 5. Distortion ranges for voltage harmonics.

Bus voltage V at PCC	Individual harmonic distortion IHD (%)	Total harmonic distortion THD (%)
V ≤ 1.0 kV	5.0	8.0
1 kV < V ≤ 69 kV	3.0	5.0
69 kV < V ≤ 161 kV	1.5	2.5
161 kV < V	1.0	1.5

Table 6 provides a summary of the permissible operational ranges for current harmonic distortion at the Point of Common Coupling (PCC) (Minh Quan Duong, 2019).

Table 6. Current harmonics distortion limits.

The maximum harmonic current distortion, presented as a percentage of the total current (IL)
Order of individual harmonics (odd harmonics)
ISC/IL	3 ≤ h < 11	11 ≤ h < 17	17 ≤ h < 23	23 ≤ h < 35	35 ≤ h < 50	TDD
<20	4.0	2.0	1.5	0.6	0.3	5.0
20 < 50	7.0	3.5	2.5	1.0	0.5	8.0
50 < 100	10.0	4.5	4.0	1.5	0.7	12.0
100 < 1000	12.0	5.5	5.0	2.0	1.0	15.0
>1000	15.0	7.0	6.0	2.5	1.4	20.0

ISC = maximum short-circuit current at PCC; IL = maximum demand load current.

The total demand distortion of current at the PCC under normal load operation circumstances can be computed as (Nagaraja & Manohar, 2016): (36) $${\mathit{\text{TDD}}}_i=\raisebox{1ex}{$\sqrt{\sum _{\mathrm{n}=2}^{\mathrm{\infty }}{I_n}^2}$}\!\left/ \!\raisebox{-1ex}{$I_L$}\right.$$(36)

The relation between total harmonic distortion of current and total demand distortion of current can be calculated as: (37) $${\mathit{\text{THD}}}_i\mathrm{*}{I}_1={\mathit{\text{TDD}}}_i\mathrm{*}{I}_L$$(37)

8. Simulation results and discussion

The MATLAB software is utilized to conduct a backward-forward sweep load flow study, aiming to ascertain the voltage profile and power losses in the base case. The influence of integrating solar PV systems with the distribution network is examined at various penetration levels, considering measures such as power loss, voltage profile, and harmonic distortion. As system input data, base MVA = 100 and base KV = 15 are used to change to per-unit values. The simulation is executed in five distinct scenarios (cases): the base case, 20% PV integration, 40% PV integration (at the optimal level), 60% PV integration, and 80% PV integration.

8.1. Base case simulation results

To ascertain the voltage profile at individual buses and calculate the real and reactive power losses in each line, a load flow analysis has been conducted using the backward/forward sweep algorithm in the base case. The majority of buses have voltage profiles that are below the 0.95 pu minimal threshold. This implies that the selected distribution network is experiencing an under-voltage problem. Bus 34 displays a minimum voltage profile of 0.8393 per unit, accompanied by a voltage deviation that accounts for 16.07% of the nominal voltage. In the base case, the system experiences active power losses of 179.17 kilowatts and reactive power losses of 111.80 kilovolt-amperes. The Particle Swarm Optimization (PSO) algorithm identifies the optimal location for the solar PV installation at bus 34.

8.2. Results of the simulation for each scenario

The voltage profiles for each case are depicted in Figure 7.

Figure 7. Comparison of network voltage profiles for all scenarios.

In Figure 7, both Scenario I and Scenario II exhibit minimum voltage profiles below the specified threshold value of 0.95 per unit. Scenarios IV and V exceeded the acceptable operating limit due to certain bus voltages surpassing the maximum threshold value of 1.05 per unit. The voltage profiles adhere to the acceptable IEEE standard limit only in the third scenario, which is the optimal level.

The power losses of the network in each scenario are illustrated in Figure 8.

Figure 8. Power loss results for all scenarios.

The effect of solar PV system integration on power loss reductions is visualized in Figure 8. Scenario III, corresponding to the optimal PV integration level at 40%, results in the minimum power losses. For scenarios IV (60%) and V (80%), the losses exceed those observed at the optimal integration level. The findings indicated that the incorporation of solar photovoltaic (PV) technology decreases power losses up to the point of optimal integration. However, integrating solar photovoltaic (PV) technology beyond its optimal size led to increased power losses.

Utilizing Particle Swarm Optimization (PSO), it has been identified that bus 34 is the optimal location for integrating PV systems. Therefore, this bus is considered a common coupling point. The extent of PV integration with the distribution network should be such that the harmonic distortion level remains within the permissible range. Figure 9 illustrates the voltage harmonic distortion across different penetration levels.

Figure 9. THDv in various PCC scenarios.

Table 5 displays the suggested total voltage distortion range at the point of common coupling. This case study uses a bus voltage of 15 kV. Therefore, 5% is the highest allowable total bus voltage harmonic distortion in this scenario at PCC.

The value of 0.75% is established as the base case total voltage harmonic distortion. Figure 9 illustrates that cases 1 and 2 fall within the permitted ranges. But, the remaining cases are beyond the standard ranges. Similar to how PV power introduced into the system increases, voltage harmonic distortion also rises. The overall voltage harmonic distortion level decreased to 4.07%, which is within IEEE standards, when the 40% optimal size PV system (744 kW) was integrated.

Table 7 illustrates the total current harmonic distortion (THDi) at the Point of Common Coupling (PCC) for the chosen feeder in each specific case. The calculation of total current demand distortion involves the utilization of the total harmonic distortion, as indicated in the provided Equation (370. A three-phase fault is induced at the Point of Common Coupling (PCC) to assess the short-circuit current in each scenario.

Table 7. THDi and TDDi at PCC for each cases.

Cases	THDi	I1(A)	IL(A)	Isc(A)	Isc/IL	TDDi	Max.TDDi
Case0	0.56	4.77	73.33	15970	217.78	0.04	15.0
Case2	37.72	11.30	73.33	16090	219.42	5.81	15.0
case4	33.33	25.15	73.33	16200	220.92	11.43	15.0
Case6	31.73	38.99	73.33	16320	222.56	15.57	15.0
Case8	30.88	52.73	73.33	16430	224.06	22.21	15.0

Figures 10 and 11 display the current total harmonic distortion and demand distortion of the distribution network, respectively.

Figure 10. THDi under different cases at PCC.

Figure 11. TDDi under different cases at PCC.

Figure 10 illustrates the total current harmonic distortion in the network under various solar PV system penetration levels at the point of common coupling. Based on these findings, the total current demand distortion is subsequently determined.

Table 6 displays the recommended upper limit for total current demand distortion at the point of common coupling. The highest achievable current demand distortion in this scenario at the Point of Common Coupling (PCC) is 15%.

The total current demand distortion for the chosen distribution feeder at the Point of Common Coupling (PCC) in the base case is determined to be 0.04%. As depicted in Figure 11, both Case 1 and Case 2 fall within the limits specified by the IEEE standard. However, the remaining cases surpass the prescribed standard limit. The rise in solar PV power injected into the system correlates with an increase in the total current demand distortion.

9. Economic analysis of integrating PV system

The economic viability of integrating a PV system into the current distribution network is assessed based on the cost savings resulting from reduced losses after PV integration and the cost associated with unsold power due to electric power failures. (38) $$\text{Cost of Energy due to power loss}=\text{Power loss}*\text{time}\ast \text{tariff for electricity}$$(38)

In the base case, the total power loss is 179.17 kW, and this decreases to 65.82 kW following the integration of the optimal-sized PV system. Taking into account the electricity price from EEU at 0.6943 Birr per kilowatt-hour, the annual cost of energy due to power loss can be computed as (IEEE Power & Energy Society, 2014): $$\begin{array}{l}\text{Cost of Energy due to power loss}\\ =179.\text{17 kW}*8760\text{hr}/\text{yr}*0.\text{6943Birr}/\text{Kwh}\\ =1,089,724.\text{124 Birr}\end{array}$$

Following the integration of a properly sized PV system at a suitable location, the power loss decreased to 65.82 kW. The cost of energy attributable to power loss decreases to 400,321.716 birr after the integration of the PV system. Furthermore, the calculated value for the cost of unsold power due to electric power failures amounts to 4.525 million birrs annually. The annual cost savings from the saved energy amount to 5.214 million birrs.

The installation cost of solar-based PV varies based on the technology and manufacturers. However, using the USA as a benchmark due to its technological advancement, the calculated investment cost of PV is (http://admin.theiguides.org/Media/Documents/Electricity%20Tariffs.pdf): $$\begin{array}{l}\text{The installation cost of PV}\\ =\text{power generated by PV}*\text{cost in dollar per watt}\\ ={\mathrm{P}}_{\text{PV}}*\text{USD}/\text{watt}\end{array}$$

The total installation price of solar panels also varies from state to state. Based on benchmark data of the USA for 2018, the average installed cost of solar power is 0.7$/W at which 20 to 25 years lifetime and efficiency of 14–16% {Formatting Citation}. The investment cost of PV is 744 kW*0.7USD/W*35 birr/USD = 18.228 million birr. The installation cost and maintenance cost are taken to be 10% of the total cost which yields 1.8228 million birrs for the installation and the maintenance cost. Therefore, the total capital cost is 20.051 million birrs.

Calculating the payback period of the finance invested for PV installation also very essential to assess economic feasibility. The payback period is the period or length of time required to cover the cost expended in an investment. It is the ratio of the investment cost of the PV system to the annual energy saved. The payback period is the basic determinant factor to undertake a project because projects having longer payback period is not a feasible investment. The payback period is the ratio of the capital cost to the saving cost as shown in Equation (39)(39) $$\begin{array}{l}\text{Payback period}\left(\mathrm{\text{yr}}\right)=\frac{\text{Capital cost}\left(\mathrm{\text{birr}}\right)}{\text{Saving cost}\left(\mathrm{\text{birr}}/\mathrm{\text{year}}\right)}\mathrm{ }\\ =3.85\mathrm{ }\mathrm{\text{years}}\end{array}$$(39) . (39) $$\begin{array}{l}\text{Payback period}\left(\mathrm{\text{yr}}\right)=\frac{\text{Capital cost}\left(\mathrm{\text{birr}}\right)}{\text{Saving cost}\left(\mathrm{\text{birr}}/\mathrm{\text{year}}\right)}\mathrm{ }\\ =3.85\mathrm{ }\mathrm{\text{years}}\end{array}$$(39)

10. Conclusion

This article examines the performance of incorporating solar photovoltaic (PV) systems into the distribution system, focusing on power losses, voltage profile, and harmonic distortion. The study utilizes the BATA radial distribution feeder, one of the 34-bus Bahir Dar feeders. The MATLAB software employs the backward-forward sweep method to calculate the bus voltage profile and power losses in branches. Through the application of particle swarm optimization, Bus 34 has been identified as the optimal location for the integration of a solar power system. The ETAP software is utilized to evaluate the extent of harmonic distortion. Different levels of solar photovoltaic (PV) system integration into the feeder are considered to assess the impact of solar PV on the distribution system, specifically in terms of power losses, voltage profile, and harmonic distortion. The investigation of solar resources confirms that solar energy is plentiful and feasible in the region under study. Upon evaluating power losses and voltage profiles, it is concluded that the optimal size for the PV system is 744 kW, equivalent to 40%. Integrating the PV system into the selected feeder beyond its optimal size resulted in challenges such as overvoltage issues and significant power losses. Integrating the PV system optimally enhances the minimum voltage profile, increasing it from 0.8392pu to 0.954pu. Additionally, active and reactive power losses are reduced by 63.27% and 66.723%, respectively. Moreover, the total harmonic distortion is determined to be 4.07%, a value comfortably below the IEEE working limit. Incorporating non-optimal levels of solar photovoltaic (PV) into the distribution network led to significant power losses, overvoltage, undervoltage, and the emergence of excessive harmonic distortion. It was mentioned that incorporating solar PV systems alone would not enhance network performance and could potentially exacerbate the situation. The analysis of cost-effectiveness and payback period indicates that a savings of 5.24 million birrs annually can be achieved four years after the integration of the solar PV system. Overall, the integration of an appropriately sized PV system at an optimal location improves the performance of the distribution system.

Acknowledgments

I would like to express my deepest appreciation to Mr. Elias Mandefro Getie, a researcher and University industry Linkage coordinator of Bahir Dar institute of technology who has supported me by editing some contents during this paper preparation. I would like also to express my deepest thanks to Mr. Habtemariam Aberie Kefale, a researcher and Vice dean for Faculty of Electrical and Computer Engineering who has supported me by giving his invaluable encouragement.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Data availability statement

Data can be accessed via author request.

Additional information

Notes on contributors

Estifanos Abeje Sharew

Estifanos Abeje Sharew has MSc degree in Power Systems Engineering. He is a senior Electrical Engineering full time lecturer @ Bahir Dar Institute of Technology, Bahir Dar University, Bahir Dar, Ethiopia and also a researcher and Assistant registrar for Faculty of Electrical and Computer Engineering.