A Comparative Study of Smart THD-Based Fault Protection Techniques for Distribution Networks

The integration of Distributed Generators (DGs) into distribution systems (DSs) leads to more reliable and efficient power delivery for customers. However, the possibility of bi-directional power flow creates new technical problems for protection schemes. This poses a threat to conventional strategies because the relay settings have to be adjusted depending on the network topology and operational mode. As a solution, it is important to develop novel fault protection techniques to ensure reliable protection and avoid unnecessary tripping. In this regard, Total Harmonic Distortion (THD) can be used as a key parameter for evaluating the grid’s waveform quality during fault events. This paper presents a comparison between two DS protection strategies that employ THD levels, estimated amplitude voltages, and zero-sequence components as instantaneous indicators during the faults that function as a kind of fault sensor to detect, identify, and isolate faults. The first method uses a Multiple Second Order Generalized Integrator (MSOGI) to obtain the estimated variables, whereas the second method uses a single SOGI for the same purpose (SOGI-THD). Both methods rely on communication lines between protective devices (PDs) to facilitate coordinated protection. The effectiveness of these methods is assessed by using simulations in MATLAB/Simulink considering various factors such as different types of faults and DG penetrations, different fault resistances and fault locations in the proposed network. Moreover, the performance of these methods is compared with conventional overcurrent and differential protections. The results show that the SOGI-THD method is highly effective in detecting and isolating faults with a time interval of 6–8.5 ms using only three SOGIs while requiring only 447 processor cycles for execution. In comparison to other protection methods, the SOGI-THD method exhibits a faster response time and a lower computational burden. Furthermore, the SOGI-THD method is robust to harmonic distortion, as it considers pre-existing harmonic content before the fault and avoids interference with the fault detection process.


Introduction
In recent years, the integration of Distributed Energy Resources (DERs) into Distribution Systems (DSs) has gained considerable attention due to the significant benefits it offers to the overall energy system. DERs, including solar and wind power systems, can contribute to improving energy efficiency and reliability in the distribution grid. They can reduce the dependence on traditional fossil fuel sources, mitigate environmental impacts, and suppsort the integration of new energy technologies, leading to a more sustainable and resilient energy industry [1,2].
However, integrating DERs into DSs requires careful management due to the significant impact it has on the power flow direction. In conventional DSs, the power flow is one-way, but with DER integration, it can be bi-directional, presenting a challenge for the network operator, particularly regarding fault protection [3,4]. Fault protection is a complex issue in DSs, especially with Distributed Generators (DGs), as power flow changes distinguish between phase-to-phase and phase-to-phase-to-ground faults and exh slower tripping responses compared to other methods.
To overcome the previous issues, two THD-based protection techniques for dete and isolating different types of faults (symmetrical and unsymmetrical) in MV DSs presented in [36,37]. The methods are based on the measurement of the THD levels o grid voltages, defined as , the estimated amplitude voltages ̃, and the sequence components 0 . Figure 1 depicts the conceptual diagram of the methods first approach in [36] calculates these variables using an MSOGI approach, conside only the fundamental and the triple-n harmonics of the grid voltages: the third, sixth ninth harmonics. The second approach in [37] uses a SOGI with a few additional m operations to obtain the variables, which was implemented into a DSP from Texas In ments with a low computational burden. The purpose of this paper is to assess the effectiveness of these techniques tog with traditional overcurrent and differential protections across a range of scenarios. cifically, the objectives of this study can be summarized as follows:

MSOGI-THD
- The MSOGI-THD and SOGI-THD protection methods are explained and mod using Matlab/Simulink to automatically detect symmetrical and unsymmetrical events in DSs; -A robust fault detection algorithm for the SOGI-THD is proposed. This approac creases the system's reliability and accuracy and speeds up the protection syst response in the event of a fault, regardless of the harmonics present in the grid b the fault. This approach is also computationally efficient; -A finite state machine has been used to locate and isolate faults in different locat detect permanent faults, and ignore temporary faults; -A comparison study between the two THD-based approaches and the convent methods is proposed under various conditions, such as changing the fault types penetration, fault resistances, and fault locations in the network. Moreover, the s assesses the robustness of these methods against communication delays and the ence of harmonics before fault events. Additionally, the paper evaluates the comp tional burden of the different THD methods when implemented on a digital proce The rest of the paper is structured as follows: Section 2 presents the MSOGI-THD SOGI-THD approaches. Section 3 presents the proposed fault classification algorithm The purpose of this paper is to assess the effectiveness of these techniques together with traditional overcurrent and differential protections across a range of scenarios. Specifically, the objectives of this study can be summarized as follows: − The MSOGI-THD and SOGI-THD protection methods are explained and modeled using Matlab/Simulink to automatically detect symmetrical and unsymmetrical fault events in DSs; − A robust fault detection algorithm for the SOGI-THD is proposed. This approach increases the system's reliability and accuracy and speeds up the protection system's response in the event of a fault, regardless of the harmonics present in the grid before the fault. This approach is also computationally efficient; − A finite state machine has been used to locate and isolate faults in different locations, detect permanent faults, and ignore temporary faults; − A comparison study between the two THD-based approaches and the conventional methods is proposed under various conditions, such as changing the fault types, DG penetration, fault resistances, and fault locations in the network. Moreover, the study assesses the robustness of these methods against communication delays and the presence of harmonics before fault events. Additionally, the paper evaluates the computational burden of the different THD methods when implemented on a digital processor. The rest of the paper is structured as follows: Section 2 presents the MSOGI-THD and SOGI-THD approaches. Section 3 presents the proposed fault classification algorithm and the FSM in detail. Simulations and comparisons with OCR and DR systems are carried out to validate the performance of each method in Section 4. Finally, the conclusion is provided in Section 5.

THD Measurement Methods
The THD is a crucial indicator for evaluating the quality and performance of the power grid. This study utilizes THD as an indicator of system faults and develops fast and sensitive fault protection methods that are both cost-effective and minimize the need for complex equipment. By relying on THD to detect faults, the efficiency and effectiveness of the protection system could be improved while minimizing the risk of costly and timeconsuming equipment failures. Overall, this study offers valuable insights into the potential of THD as a fault detection mechanism.
To calculate THD, the square root of the sum of the harmonic components of a signal, squared and divided by the fundamental component, is used. This is per the standard definition outlined in references [38,39] and given by: where h and A h are the harmonic order and the amplitude of the h-th-harmonic component, respectively, for h = 1, and A 1 is the amplitude of the fundamental component. In this study, the MSOGI-THD and SOGI-THD methods are employed to acquire the harmonics required for determining the THD.

MSOGI-THD
The method described in [36] is used to obtain the harmonic components of the grid in order to calculate the THD. In [36], an MSOGI is employed together with a Frequency Locked-Loop (MSOGI-FLL) to obtain the fundamental and the 3rd, 6th, and 9th harmonics of the grid voltages. The MSOGI uses multiple SOGIs operating in parallel and a crossfeedback cancellation network to remove the unnecessary components from the input signal, ensuring that each SOGI receives only the necessary component. Figure 2 illustrates this configuration, where the input voltage signal is v in , the error is e, and the estimated angular frequency is ω.
The MSOGI system uses an FLL to track the operative frequency of the grid. The FLL is linked to the first SOGI, which provides the fundamental component. The FLL delivers ω to the paralleled SOGIs. The MSOGI can be affected by faults, which distort the estimated ω and cause further distortions to the system. Then, a saturation block is applied to the FLL output to limit the distortion in ω to be restricted within ±1 Hz of the grid's nominal frequency. In the MSOGI-FLL, the parameter ξ is set to 1/ √ 2 to achieve an optimal relationship between rejection to harmonic distortion and transient response speed [40].
In [36], the THD is calculated using the triple-n harmonics (i.e., the 3rd, 6th, and 9th) since these harmonics are typically present only in the power inverter's AC-side neutral point [41]. This ensures that the THD is not influenced by other harmonics that might exist in the grid prior to the fault. When a fault occurs in the grid, there is a sudden drop in the voltage phases, exciting almost all harmonic components. Consequently, the triple-n harmonics are also excited, leading to a sharp pulse in the THD. This allows fast detection despite the presence of other harmonics in the grid before the fault. Therefore, during normal operation without faults, the presence of harmonics in the grid does not affect the THD calculation because they are not considered in the computation.  The MSOGI-FLL is used for each phase of the grid voltage and used to determin fundamental and 3rd, 6th, and 9th amplitudes of the harmonic components by using in-phase and quadrature-phase outputs of each SOGI as follows: where ℎ is the index of the harmonic component. To calculate the THD for each pha the grid, the formula in (1) is used, which involves basic mathematical operations suc summing the squared harmonic components, taking the square root, and dividin shown in Figure 3. The resulting THD is denoted as � . It should be noted that m plication by 100 is only necessary to convert the value to a percentage scale, while sa tion is applied to prevent a potential division by 0 during the initial stages of the sys Additionally, a 2nd-order low-pass filter is designed to achieve a balance between transient response speed and the levels of distortion in the THD signal. Figure 3. Block diagram of the MSOGI-THD measurement method.

SOGI-THD
The method for measuring THD outlined in [37] is employed in this paper, in w The MSOGI-FLL is used for each phase of the grid voltage and used to determine the fundamental and 3rd, 6th, and 9th amplitudes of the harmonic components by using the in-phase and quadrature-phase outputs of each SOGI as follows: where h is the index of the harmonic component. To calculate the THD for each phase of the grid, the formula in (1) is used, which involves basic mathematical operations such as summing the squared harmonic components, taking the square root, and dividing, as shown in Figure 3. The resulting THD is denoted as v THD . It should be noted that multiplication by 100 is only necessary to convert the value to a percentage scale, while saturation is applied to prevent a potential division by 0 during the initial stages of the system. Additionally, a 2nd-order low-pass filter is designed to achieve a balance between the transient response speed and the levels of distortion in the THD signal.
Sensors 2023, 23, x FOR PEER REVIEW 5 The MSOGI-FLL is used for each phase of the grid voltage and used to determine fundamental and 3rd, 6th, and 9th amplitudes of the harmonic components by using in-phase and quadrature-phase outputs of each SOGI as follows: where ℎ is the index of the harmonic component. To calculate the THD for each phas the grid, the formula in (1) is used, which involves basic mathematical operations suc summing the squared harmonic components, taking the square root, and dividing shown in Figure 3. The resulting THD is denoted as � . It should be noted that m plication by 100 is only necessary to convert the value to a percentage scale, while sat tion is applied to prevent a potential division by 0 during the initial stages of the sys Additionally, a 2nd-order low-pass filter is designed to achieve a balance between transient response speed and the levels of distortion in the THD signal.

SOGI-THD
The method for measuring THD outlined in [37] is employed in this paper, in w the SOGI approach is used to obtain the harmonic components required for calcula THD [42,43]. The technique described in [37] uses the definition of THD outline

SOGI-THD
The method for measuring THD outlined in [37] is employed in this paper, in which the SOGI approach is used to obtain the harmonic components required for calculating THD [42,43]. The technique described in [37] uses the definition of THD outlined in [38,39], which is calculated using Equation (1). A SOGI-FLL is employed for 1 of the 3 phases of the grid voltage, while only a SOGI is utilized for the remaining two phases, as ω is provided to them by the FLL stage of the 1st SOGI.
Consider a grid voltage composed of a fundamental and harmonic component as follows: where ω i and ϕ h represents the grid frequency and the phase angle of the h-th-harmonic component, respectively. A SOGI filter is used for each phase to extract its fundamental component, A 1 , by using its in-phase and quadrature-phase outputs and (2), while the rest of the harmonic components are given through its error signal, e(t). Thus, by squaring e(t) and applying trigonometric identities, it is found that the dc component is equivalent to half the sum of the square of the amplitude of the harmonic components. The dc component can be extracted by applying a low-pass filter with an appropriate cut-off frequency, resulting in the average value of e(t) as follows: Now, the THD is obtained, by multiplying (4) by 2 to remove the 1/2 gain and by taking the square root, as: where µ = e 2 /A 2 . Figure 4 depicts the THD's block diagram. The THD output is denoted as v THD . Note that multiplication by 100 is used to show the value on a percentage scale, and saturation is employed to prevent an eventual division by 0. The low-pass filter is designed to achieve a balance between speed transient response and distortion levels in the THD signal.
Sensors 2023, 23, x FOR PEER REVIEW phases of the grid voltage, while only a SOGI is utilized for the remaining tw � is provided to them by the FLL stage of the 1st SOGI. Consider a grid voltage composed of a fundamental and harmonic com follows: where and ℎ represents the grid frequency and the phase angle of the h-t component, respectively. A SOGI filter is used for each phase to extract its fu component, 1 , by using its in-phase and quadrature-phase outputs and (2), w of the harmonic components are given through its error signal, ( ). Thus, b ( ) and applying trigonometric identities, it is found that the dc component is to half the sum of the square of the amplitude of the harmonic components. T ponent can be extracted by applying a low-pass filter with an appropriate quency, resulting in the average value of ( ) as follows: Now, the THD is obtained, by multiplying (4) by 2 to remove the 1/2 g taking the square root, as: Figure 4 depicts the THD's block diagram. The THD outpu as � . Note that multiplication by 100 is used to show the value on a perce and saturation is employed to prevent an eventual division by 0. The low-p designed to achieve a balance between speed transient response and distorti the THD signal. The zero-sequence component is calculated as follows: And, Figure 5 illustrates the block diagram used to obtain the 0-sequence Note that in the SOGI-THD method, only the in-phase voltages of the SOGI u phase of the voltage grid are necessary, i.e., 0 , 0 , and 0 . For the case of t THD, only the 1st SOGI of the scheme in Figure 2 is used to achieve the 0-seq ponent.
The measurement of THD using the MSOGI-THD and SOGI-THD meth able to detect faults since they suppose a sudden sharp change in grid voltages all harmonic components. Therefore, the THD obtained using these methods for a fast of faults, which is the core of the proposed detection algorithms.
As the behavior of , ̃ , and 0 will differ depending on The zero-sequence component is calculated as follows: And, Figure 5 illustrates the block diagram used to obtain the 0-sequence component. Note that in the SOGI-THD method, only the in-phase voltages of the SOGI used in each phase of the voltage grid are necessary, i.e., v a0 , v b0 , and v c0 . For the case of the MSOGI-THD, only the 1st SOGI of the scheme in Figure 2 is used to achieve the 0-sequence component.
The measurement of THD using the MSOGI-THD and SOGI-THD methods will be able to detect faults since they suppose a sudden sharp change in grid voltages that excites all harmonic components. Therefore, the THD obtained using these methods can be used for a fast of faults, which is the core of the proposed detection algorithms. Zero sequences measurement

Fault Classification Algorithm
An algorithm has been developed for identifying and detecting faults a cation. This algorithm uses measurements of the , ̃, and 0 phase voltage signals to determine the type of fault. For identifying phase-to-0 will be used. Figure 6 illustrates a single-line diagram of the DS used for testing the m the system parameters are provided in Table 1. The diagram consists of a m multiple DGs connected to various buses. The high-voltage (HV) grid has a r of 66 kV and a rated power of 25 MVA, and it is connected to three Distri (DLs; DL1, 2, and 3) through a step-down HV/MV transformer. This transfo figured in a star/delta (YNd11) setup, which means its ground connection can impact the short-circuit current in the network and, therefore, the behavior o tion system. To overcome this issue, a zig-zag transformer is used to creat neutral in the power system. This common configuration provides a balance system even in the presence of ground faults. In this work, the zig-zag transfo nected to the delta side of the voltage transformer, which is grounded throu Bus 2, as shown in Figure 6 [44]. Each DL has two PDs located on either sid with an FSM intended for fault isolation. The PDs contain a fault detection r cuit breaker that trip and isolate the line based on a message from the FSM occurs. Different local loads (L1 to L3) are connected to the end of each DL. Tw and 2) are connected, through an MV/LV transformer, to different locations o Loads and DGs are connected to the low-voltage side using a delta/star groun configuration [44]. Communication channels are employed between the PDs and between the FSMs themselves, for transmitting trip signals to coordinat faulted locations. The proposed approach's behavior will be tested under va ios in different locations, defined as F1 and F2, in Figure 6.  As the behavior of THD abc , A abc , and V abc0 will differ depending on the type of fault, it is essential to have a fault classification algorithm to accurately identify the type of fault.

Fault Classification Algorithm
An algorithm has been developed for identifying and detecting faults at each PD location. This algorithm uses measurements of the THD abc , A abc , and V abc0 of the threephase voltage signals to determine the type of fault. For identifying phase-to-phase faults, V abc0 will be used. Figure 6 illustrates a single-line diagram of the DS used for testing the methods, and the system parameters are provided in Table 1. The diagram consists of a main grid and multiple DGs connected to various buses. The high-voltage (HV) grid has a rated voltage of 66 kV and a rated power of 25 MVA, and it is connected to three Distribution Lines (DLs; DL1, 2, and 3) through a step-down HV/MV transformer. This transformer is configured in a star/delta (YNd11) setup, which means its ground connection can significantly impact the short-circuit current in the network and, therefore, the behavior of the protection system. To overcome this issue, a zig-zag transformer is used to create an isolated neutral in the power system. This common configuration provides a balanced and stable system even in the presence of ground faults. In this work, the zig-zag transformer is connected to the delta side of the voltage transformer, which is grounded through zig-zag at Bus 2, as shown in Figure 6 [44]. Each DL has two PDs located on either side of the line, with an FSM intended for fault isolation. The PDs contain a fault detection relay and circuit breaker that trip and isolate the line based on a message from the FSM when a fault occurs. Different local loads (L1 to L3) are connected to the end of each DL. Two DGs (DG1 and 2) are connected, through an MV/LV transformer, to different locations of the system. Loads and DGs are connected to the low-voltage side using a delta/star grounded (Dyn11) configuration [44]. Communication channels are employed between the PDs and the FSM, and between the FSMs themselves, for transmitting trip signals to coordinate and isolate faulted locations. The proposed approach's behavior will be tested under various scenarios in different locations, defined as F1 and F2, in Figure 6.

Fault Classification Algorithm Stages
In this work, the algorithm uses three stages to perform a rapid and secure decision to detect and identify types of faults that may occur in the grid.

Fault Classification Algorithm Stages
In this work, the algorithm uses three stages to perform a rapid and secure decision to detect and identify types of faults that may occur in the grid.

Pre-Processing Stage
This stage involves measuring the three-phase voltages, , at each PD in the time domain. From these measurements, the , ̃, and 0 are calculated using the MSOGI-THD and SOGI-THD approaches, as previously explained.

MSOGI-THD
In this case, the fault is detected using a threshold, . The transient response induced by the fault in is compared with , and the detection is activated when > . is set to 5%, which is the recommended level for voltage harmonic distortion for DS in the IEEE standard 519-2014 [45].

SOGI-THD
In this case, the fault is detected using the transient response and low-pass filtering of THD, noted as , which allows the detection of faults even when there is harmonic distortion in the grid. The design of this proposal is illustrated in Figure 7.

Pre-Processing Stage
This stage involves measuring the three-phase voltages, v abc , at each PD in the time domain. From these measurements, the THD abc , A abc , and V abc0 are calculated using the MSOGI-THD and SOGI-THD approaches, as previously explained.

MSOGI-THD
In this case, the fault is detected using a threshold, α MSOGI . The transient response induced by the fault in THD abc is compared with α MSOGI , and the detection is activated when THD abc > α MSOGI . α MSOGI is set to 5%, which is the recommended level for voltage harmonic distortion for DS in the IEEE standard 519-2014 [45].

SOGI-THD
In this case, the fault is detected using the THD abc transient response and low-pass filtering of THD, noted as σ LPF abc , which allows the detection of faults even when there is harmonic distortion in the grid. The design of this proposal is illustrated in Figure 7.
The detection is done based on the difference between THD abc and σ LPF abc . at the precise instant a fault occurs. This difference is known as α abc : where α abc represents the difference calculated for each phase. The presence of a fault is detected when α abc exceeds a predefined threshold, α abc ≥ α o . In this case, α o was set to 25% to ensure a sufficient margin for avoiding false triggering caused by unexpected transients. The σ LPF abc was obtained by applying a low-pass filter with an appropriate cutoff frequency to THD abc , thus obtaining its average.

MSOGI-THD
In this case, the fault is detected using a threshold, . The transient respo induced by the fault in is compared with , and the detection is activa when > . is set to 5%, which is the recommended level for volt harmonic distortion for DS in the IEEE standard 519-2014 [45].

SOGI-THD
In this case, the fault is detected using the transient response and low-p filtering of THD, noted as , which allows the detection of faults even when ther harmonic distortion in the grid. The design of this proposal is illustrated in Figure 7.    Figure 8 depicts the behavior of the THD for phase "b" during a fault at 0.2 s, and simultaneously, there is a fifth harmonic with a 5% amplitude distortion in the grid voltage. Prior to the fault, the harmonic had been there for a sufficient duration to make the output of the LPF equal to the input, so α b = 0. When a fault occurs, a sudden, sharp peak-pulse waveform is produced in THD b that disappears over time, following an exponential decay pattern. This transient is filtered by the low-pass filter (σ LPF b ) which has a much slower response than the THD b and takes a lot of time to be observed at its output. However, at the moment of the fault, σ LPF b . remains close to the THD b level that existed prior to the fault. Thus, the difference, α b , closely corresponds to the difference between the actual THD b and the past THD b level that the grid had at the moment of the fault. Therefore, the utilization of the low-pass filter in this manner allows for the detection of sudden changes in THD b caused by faults, regardless of the prior harmonic distortion levels that might be present in the grid voltage. Upon detection, the algorithm reads and records the value of σ LPF b , referred to as σ LPFm b , which will be used later for the identification process. The detection is done based on the difference between and . at the precise instant a fault occurs. This difference is known as : represents the difference calculated for each phase. The presence of a fault is detected when exceeds a predefined threshold, ≥ . In this case, was set to 25% to ensure a sufficient margin for avoiding false triggering caused by unexpected tran sients. The was obtained by applying a low-pass filter with an appropriate cutoff frequency to , thus obtaining its average. Figure 8 depicts the behavior of the THD for phase "b" during a fault at 0.2 s, and simultaneously, there is a fifth harmonic with a 5% amplitude distortion in the grid volt age. Prior to the fault, the harmonic had been there for a sufficient duration to make the output of the LPF equal to the input, so = 0. When a fault occurs, a sudden, sharp peak-pulse waveform is produced in that disappears over time, following an expo nential decay pattern. This transient is filtered by the low-pass filter ( ) which has a much slower response than the and takes a lot of time to be observed at its output However, at the moment of the fault, . remains close to the level that existed prior to the fault. Thus, the difference, , closely corresponds to the difference between the actual and the past level that the grid had at the moment of the fault Therefore, the utilization of the low-pass filter in this manner allows for the detection o sudden changes in caused by faults, regardless of the prior harmonic distortion levels that might be present in the grid voltage. Upon detection, the algorithm reads and records the value of , referred to as , which will be used later for the identifi cation process.

Fault Identification Stage
In this stage, the identification process, in case of using any of the two THD ap ̃

Fault Identification Stage
In this stage, the identification process, in case of using any of the two THD approaches, is done based on the behavior of A abc and V abc0 , depending on the fault case.
The MSOGI-THD approach uses a threshold named ∆v for comparing with A abc . ∆v is set to 7.5%, which is the acceptable range for grid voltage drop at the medium voltage level as recommended by the technical requirements in the Spanish grid code for reliable energy integration [46]. In contrast, the SOGI-THD approach employs two thresholds for identifying faults, σ LPFm abc for comparing with THD abc , which should be verified as THD abc > σ LPFm abc upon detection, and ∆v, for comparing with A abc . In this case, ∆v is also set to 7.5% [46].
In both THD approaches, during a fault event, a sudden drop in A abc is produced due to the fall in the voltages of the faulted phases, which makes any of the estimated phases be below the threshold, i.e., A abc < (1 − ∆v) pu. Then, this condition is considered a fault event. At this point, V abc0 is used by both methods (THD approaches) to differentiate between a phase-to-phase (2PH) and phase-to-phase-ground (2PH-G) fault, as the initial conditions for both types of faults are the same. In the event of a 2PH fault, the zerosequence voltages are zero, while they are non-zero in the case of a 2PH-G fault [47]. Therefore, V abc0 is utilized to determine the type of fault and adopt a decision. The faults were grouped into 11 categories, numbered from zero to 10, see Figure 1. Figure 9 depicts the fault classification algorithm structure. as > upon detection, and Δ , for comparing with . In this case, Δ is also set to 7.5% [46].
In both THD approaches, during a fault event, a sudden drop in is produced due to the fall in the voltages of the faulted phases, which makes any of the estimated phases be below the threshold, i.e., < (1 − Δ )pu. Then, this condition is considered a fault event. At this point, is used by both methods (THD approaches) to differentiate between a phase-to-phase (2PH) and phase-to-phase-ground (2PH-G) fault, as the initial conditions for both types of faults are the same. In the event of a 2PH fault, the zerosequence voltages are zero, while they are non-zero in the case of a 2PH-G fault [47]. Therefore, is utilized to determine the type of fault and adopt a decision. The faults were grouped into 11 categories, numbered from zero to 10, see Figure 1. Figure 9 depicts the fault classification algorithm structure.

Finite State Machine
In the grid-line diagram of Figure 6, three FSMs have been defined to locate and isolate faults in different locations of the DS. Fault isolation refers to the process of disconnecting a faulted part of a system, in this case, DL, while keeping the rest of the system operational by maintaining power to the rest of the DLs and the corresponding loads they serve [47]. This allows for repairs to be made to the faulted section without disrupting power to the entire system or loads. Fault isolation is an important aspect of power system protection and helps to minimize the impact of faults on the system and the customers.
In addition, the FSMs are designed to detect permanent faults and ignore temporary faults, thus avoiding unnecessary tripping of protection devices, which can save time and resources. Note that a permanent fault is expected to persist, while a temporary fault is expected to clear on its own without any protection action. The duration of a temporary fault is usually less than 100 ms [48,49].
In this system, each FSM is allocated to a DL with its two PDs and comprises six distinct states, as depicted in Figure 10. The figure provides a description of the FSM for

Finite State Machine
In the grid-line diagram of Figure 6, three FSMs have been defined to locate and isolate faults in different locations of the DS. Fault isolation refers to the process of disconnecting a faulted part of a system, in this case, DL, while keeping the rest of the system operational by maintaining power to the rest of the DLs and the corresponding loads they serve [47]. This allows for repairs to be made to the faulted section without disrupting power to the entire system or loads. Fault isolation is an important aspect of power system protection and helps to minimize the impact of faults on the system and the customers.
In addition, the FSMs are designed to detect permanent faults and ignore temporary faults, thus avoiding unnecessary tripping of protection devices, which can save time and resources. Note that a permanent fault is expected to persist, while a temporary fault is expected to clear on its own without any protection action. The duration of a temporary fault is usually less than 100 ms [48,49].
In this system, each FSM is allocated to a DL with its two PDs and comprises six distinct states, as depicted in Figure 10. The figure provides a description of the FSM for the two THD approaches, where the blue color denotes the MSOGI-THD, the red color represents the SOGI-THD, and the green color signifies the common conditions of both approaches.

State S1. Normal Operation
The FSM typically begins in state S1 and stays there, while the system operates within a boundary of 1 ± ∆ pu around nominal values and < for the MSOGI THD and | | < for the SOGI-THD. Therefore, at S1, the FSM is waiting for a fault.

State S2. Fault Detection
When a fault occurs, the PD in charge of a DL detects it by observing > and ̃ < (1 − Δ ) pu for the MSOGI-THD or ≥ for the SOGI-THD and ̃ < (1 − Δ )pu. At this point, a Wi-Fi protocol communication [50] is initiated between the PD and the FSM of a faulted DL to send a fault message to the FSM. Upon receiving the message, the FSM transitions to state S2.
It is important to mention that the detection time of a fault by the algorithm inside the PD of a DL is influenced by the distance between the PD and the fault location. Specifically the PD closest to the fault is expected to detect the fault in the shortest time. On the othe hand, the Wi-Fi communication delay for a fault message to be transmitted will be smal enough, less than 1 ms [51]. This gives the FSM a chance to have a high level of responsive ness to quickly locate the fault and adopt the decision to start the isolation process.
In the event of a communication failure, the PD located at the other end of the faulted DL is responsible for transmitting the fault message to the FSM. This is because it is the second-fastest PD to detect the fault. As a result, the FSM of the faulted DL that is closes to the fault will generally locate the fault more quickly than the other FSMs. This feature

State S1: Normal Operation
The FSM typically begins in state S1 and stays there, while the system operates within a boundary of 1 ± ∆v pu around nominal values and THD abc < α MSOGI for the MSOGI-THD and |α abc | < α o for the SOGI-THD. Therefore, at S1, the FSM is waiting for a fault.

State S2: Fault Detection
When a fault occurs, the PD in charge of a DL detects it by observing THD abc > α MSOGI and A abc < (1 − ∆v) pu for the MSOGI-THD or α abc ≥ α o for the SOGI-THD, and A abc < (1 − ∆v) pu. At this point, a Wi-Fi protocol communication [50] is initiated between the PD and the FSM of a faulted DL to send a fault message to the FSM. Upon receiving the message, the FSM transitions to state S2.
It is important to mention that the detection time of a fault by the algorithm inside the PD of a DL is influenced by the distance between the PD and the fault location. Specifically, the PD closest to the fault is expected to detect the fault in the shortest time. On the other hand, the Wi-Fi communication delay for a fault message to be transmitted will be small enough, less than 1 ms [51]. This gives the FSM a chance to have a high level of responsiveness to quickly locate the fault and adopt the decision to start the isolation process.
In the event of a communication failure, the PD located at the other end of the faulted DL is responsible for transmitting the fault message to the FSM. This is because it is the second-fastest PD to detect the fault. As a result, the FSM of the faulted DL that is closest to the fault will generally locate the fault more quickly than the other FSMs. This feature is utilized to achieve faster and more efficient fault isolation, thereby minimizing the impact of the fault. Now, once the faulted DL's FSM enters S2, after receiving a fault detection message from the PD in charge, a timer named T c is started inside to measure the fault duration. Additionally, a new fault signal is sent to the other FSMs via communication links, advising them of the event and allowing them to stop their processes. For example, if a fault occurs at DL1, FSM1 will detect the fault faster than the other FSMs.
The variability in the transmission time of a fault message between FSMs can be attributed to the communication technology employed and additional delays, referred as T d . This system utilizes the IEC 61850 protocol, which is estimated to have a delay of 10 ms [52]. As a result, this protection system is reliant on a communication delay between the FSMs to coordinate and operate effectively. However, delays between 10 ms and 100 ms can occur, which may cause an FSM to receive a fault message during its process. In cases where the delay exceeds 100 ms, the FSM may make an independent decision, potentially resulting in the tripping of breakers.
Meanwhile, in S2, if the fault disappears, as in a temporary fault, then THD abc < α MSOGI for the MSOGI-THD, or |α abc |< α o for the SOGI-THD, and A abc will return to its normal value. The FSM will then return to its normal operation state (S1). Additionally, if the FSM receives a fault message from the other FSMs, it will transition to state S6, which is used to hold non-faulted FSMs.
The THDs exhibit a spike when a fault occurs, which gradually decreases over time and eventually disappears. Therefore, whenever THD abc < α MSOGI for the MSOGI-THD, or |α abc |< α o for the SOGI-THD, while A abc < (1 − ∆v) pu still holds, the FSM enters state S3 to check for the permanence of the fault. This peak behavior that exponential decay can be observed in the figures presented in the results section of the paper.

State S3: Waiting for the Fault to Be Permanent
In this state, the timer T c continues to count up, and once it reaches 100 ms, the FSM transitions to state S4 for fault isolation. However, if the fault disappears, i.e., A abc returns to its normal value before T c reaches 100 ms, the FSM returns to S1 (normal operation) without tripping the circuit breakers. In the event that the FSM receives a fault message during state S3, and T c is less than 100 ms, it transitions to state S6, indicating that the fault is not located within the FSM's protection area.

State S4: Fault Isolation
The FSM enters this state only when T c reaches 100 ms, indicating that the fault is permanent. Therefore, the FSM sends an active trip signal to the nearest circuit breaker to isolate the fault. Furthermore, if any of the non-faulted FSMs also enter state S4, this implies that they did not receive a fault message from the faulted FSM due to a communication loss. After the tripping, T c is reset to 0, and the FSM goes to state S5.

State S5: Grid Reconnection
In this state, the FSM waits for the grid operator or maintenance personnel to reconnect the line and return the FSM to its normal operation state (S1).

State S6: Holding Non-Faulted FSMs
Finally, the FSM transitions to this state upon receiving a fault signal from a faulted FSM in any of the S1, S2, or S3 states. At S6, the FSM keeps waiting until the fault either disappears or is isolated. This waiting action is accomplished using a count-up timer called T w . When T w reaches 100 ms, the FSM returns to its normal operating state (S1).

Results
In this section, a comparison is conducted first between the MSOGI-THD and SOGI-THD protection methods and then with other conventional protection methods. The goal is to evaluate the performance and effectiveness of each method in detecting and isolating faults in the proposed DS (see Figure 6). The methods are analyzed using MAT-LAB/Simulink software under various scenarios, including changes in DG penetration, fault type, fault resistance, and location, as well as in case of pre-existing harmonics and additional communication delay. The parameters of the network in Figure 6 are shown in Table 1.

DG Penetration Protection Test
In this section, the two protection methods are evaluated in the event of a fault occurring at 0.2 s with varying levels of DG penetration. The conditions used in the evaluation are given in Table 2.  12 depict the behavior of the MSOGI-THD and SOGI-THD protection methods, respectively, during a 3PH-G with 2 DGs at F1; see Figure 6. Notice in both figures how THD abc increases abruptly at 0.2 s while A abc drops toward 0 pu due to the low fault resistance (r = 0.001 Ω). In Figure 11, the fault is detected at 0.2065 s when THD abc > 5%. At 0.209 s, the condition A abc < 0.925 is met; therefore, the fault is identified as an (ABC) fault.
In Figure 12, the fault is detected at 0.204 s when α abc > 25 and σ LPF abc is read and stored in memory (σ LPFm abc ). At 0.206 s, the conditions A abc < 0.925 and THD abc > σ LPFm abc are met; therefore, the fault is identified as an (ABC) fault. In both methods, once PD2 detects and identifies the fault, a fault message is sent to FSM1 using wireless communication channels. It is worth noting that if PD2 fails for any reason, or in the case of communication problems, PD1 will be the second fastest to be in charge of detecting the fault.
When FSM1 receives the fault message from PD2, it immediately moves to S2 and sends a fault message to other non-faulted FSMs while starting a count-up on its own T c timer. In this fault-case scenario, an extra delay of 25 ms has been intentionally added to the FSM's communication links to show its effect on the protection system. As a result, the non-faulted FSMs detect the fault and move to S2 before receiving the message. For instance, in Figure 11, FSM2 and 3 detect the fault (i.e., move to S2) after 23 ms due to the distance between the fault location and the PDs of DL2 and 3.
Once the FSMs receive the message, 25 ms later due to the delay, they transition from S2 to S6 (hold state), and T w starts counting up to wait for the fault to be isolated. The THD abc 's behavior in FSM1 has a peak that decays exponentially to zero after a short time. Then, when THD abc returns to be below 5% in Figure 11 or when the absolute value of α abc is less than the threshold α o (|α abc | < α o ) in Figure 12, while A abc < 0.925 pu still holds, FSM1 transitions from S2 to S3, and T c continues counting up. Upon reaching T c = 100 ms, the fault is declared permanent, and FSM1 moves from S3 to S4, causing the circuit breakers in DL1 to trip and isolate the fault. Additionally, T c is reset to 0, and FSM1 transitions from S4 to S5 to wait for actions to reconnect DL1 to the grid. Note that in FSM1, a momentary disruption in the grid occurs when transitioning from S3 to S4 due to DL1 disconnection, causing THD abc to spike again. However, this does not affect the FSM since it is in state S4 and does not consider THD. In non-faulted FSMs, once T w reaches 100 ms, they return to S1 (normal operation).
ures how increases abruptly at 0.2 s while drops toward 0 pu due to the low fault resistance ( = 0.001 Ω). In Figure 11, the fault is detected at 0.2065 s when > 5%. At 0.209 s, the condition ̃< 0.925 is met; therefore, the fault is identified as an (ABC) fault.
In Figure 12, the fault is detected at 0.204 s when > 25 and is read and stored in memory ( ). At 0.206 s, the conditions ̃< 0.925 and > are met; therefore, the fault is identified as an (ABC) fault. In both methods, once PD2 detects and identifies the fault, a fault message is sent to FSM1 using wireless communication channels. It is worth noting that if PD2 fails for any reason, or in the case of communication problems, PD1 will be the second fastest to be in charge of detecting the fault.    When FSM1 receives the fault message from PD2, it immediately moves to S2 and sends a fault message to other non-faulted FSMs while starting a count-up on its own timer. In this fault-case scenario, an extra delay of 25 ms has been intentionally added to the FSM's communication links to show its effect on the protection system. As a result, the non-faulted FSMs detect the fault and move to S2 before receiving the message. For in-  Figure 13a,b shows the grid currents (i abc ) during a 3PH-G fault with 2DGs at DL1 for the MSOGI-THD and SOGI-THD protection systems, respectively. Note in Figure 13a that, due to the tripping signal sent by FSM1, both PD1 and PD2 of DL1 are tripped at 0.309 s. In Figure 13b, note that the PDs trip at 0.306 s.   Figure 6. During this test, a 5th harmonic distortion was introduced in the grid before the fault. Figure 14 indicates that remains unaffected by the 5th harmonic distortion because the MSOGI-THD design does not consider this harmonic in his scheme. In contrast, Figure 15 shows that the SOGI-THD detected the 5th harmonic distortion. However, since < 25, it will not be considered a fault, and the algorithm will not trigger any protective action.
It is worth noting that, at the time of the fault (0.2 s) in both Figures, only increases abruptly, whereas remains at 0%. Whereas ̃ drops towards 0.5 pu due to the absence of the ground connection and fault resistance ( = 0.1 Ω), while ̃ remains unaffected (1 pu). In Figure 14, the fault is detected at 0.2055 s when > 5%. At 0.2085 s, the condition ̃< 0.925 is met and 0 is checked to be zero, leading to the identification of the fault as a (BC) fault. Figure 15 shows that the fault is detected at 0.204 s when > 25 and is read and stored in memory ( ). At 0.2065 s, the algorithm identifies the fault as a (BC) type since the ̃< 0.925 , > , and 0 = 0 conditions are met. It is worth mentioning that, in both protection methods, after PD6 detects the fault, a fault message is transmitted wirelessly to FSM3.  Figure 6. During this test, a 5th harmonic distortion was introduced in the grid before the fault. Figure 14 indicates that THD abc remains unaffected by the 5th harmonic distortion because the MSOGI-THD design does not consider this harmonic in his scheme. In contrast, Figure 15 shows that the SOGI-THD detected the 5th harmonic distortion. However, since α abc < 25, it will not be considered a fault, and the algorithm will not trigger any protective action.
It is worth noting that, at the time of the fault (0.2 s) in both Figures, only THD bc increases abruptly, whereas THD a remains at 0%. Whereas A bc drops towards 0.5 pu due to the absence of the ground connection and fault resistance (r = 0.1 Ω), while A a remains unaffected (1 pu). In Figure 14, the fault is detected at 0.2055 s when THD bc > 5%. At 0.2085 s, the condition A bc < 0.925 is met and V abc0 is checked to be zero, leading to the identification of the fault as a (BC) fault. Figure 15 shows that the fault is detected at 0.204 s when α bc > 25 and σ LPF bc is read and stored in memory (σ LPFm bc ). At 0.2065 s, the algorithm identifies the fault as a (BC) type since the A bc < 0.925, THD bc > σ LPFm bc , and V abc0 = 0 conditions are met. It is worth mentioning that, in both protection methods, after PD6 detects the fault, a fault message is transmitted wirelessly to FSM3. Upon receiving the fault message from PD6, FSM3 enters state S2 and sends a notification to the other non-faulted FSMs while initiating the countdown on its timer. These FSMs receive the message after 10 ms due to the normal delay of the communication links Upon receiving the fault message from PD6, FSM3 enters state S2 and sends a notification to the other non-faulted FSMs while initiating the countdown on its timer. These FSMs receive the message after 10 ms due to the normal delay of the communication links Upon receiving the fault message from PD6, FSM3 enters state S2 and sends a notification to the other non-faulted FSMs while initiating the countdown on its T c timer. These FSMs receive the message after 10 ms due to the normal delay of the communication links and move from S1 to S6 (hold state), and T w starts counting up. The THD abc 's behavior of FSM3 suffers from a peak that decays exponentially to zero in a short time. Then, when THD bc returns to THD bc < 5% in Figure 14 or when |α bc | < α o in Figure 15 while A bc < 0.925 pu still holds, FSM3 goes from S2 to S3, and T c continues counting up. Whenever T c reaches 100 ms, the fault is considered permanent, and FSM3 moves to the fault isolation state (S4), sending a signal to trip the circuit breakers in DL3. Then, T c is reset, and FSM3 moves to S5, where it waits for reconnection actions. While FSM3 is in state S4, all THD signals experience a spike resulting from DL3 disconnection, which causes an abrupt change in the grid voltages. However, since FSM3 is in S4, this does not affect its operation. At the non-faulted FSMs, the timer T w continues counting up, and once it reaches 100 ms, they return to the normal operation state (S1). Figure 16a,b shows the grid currents (i abc ) during a 2PH fault without DGs at DL3 for the MSOGI-THD and SOGI-THD protection systems, respectively. It can be observed that a fifth harmonic of 5% is present before the fault. Note in Figure 16a that, due to the tripping signal sent by FSM3, both PD5 and PD6 of DL3 are tripped at 0.3085 s. In Figure 16b, note that the PDs trip at 0.3065 s. and move from S1 to S6 (hold state), and starts counting up. The 's behavior of FSM3 suffers from a peak that decays exponentially to zero in a short time. Then, when returns to < 5% in Figure 14 or when | | < in Figure 15 while ̃ < 0.925 pu still holds, FSM3 goes from S2 to S3, and continues counting up. Whenever reaches 100 ms, the fault is considered permanent, and FSM3 moves to the fault isolation state (S4), sending a signal to trip the circuit breakers in DL3. Then, is reset, and FSM3 moves to S5, where it waits for reconnection actions. While FSM3 is in state S4, all THD signals experience a spike resulting from DL3 disconnection, which causes an abrupt change in the grid voltages. However, since FSM3 is in S4, this does not affect its operation. At the non-faulted FSMs, the timer continues counting up, and once it reaches 100 ms, they return to the normal operation state (S1). Figure 16a,b shows the grid currents ( ) during a 2PH fault without DGs at DL3 for the MSOGI-THD and SOGI-THD protection systems, respectively. It can be observed that a fifth harmonic of 5% is present before the fault. Note in Figure 16a that, due to the tripping signal sent by FSM3, both PD5 and PD6 of DL3 are tripped at 0.3085 s. In Figure  16b, note that the PDs trip at 0.3065 s.

Comparison with Different Protection Methods
The system presented in Figure 6 was used to test the effectiveness of the conventional overcurrent and differential protections [53] and compare them with the MSOGI-THD and SOGI-THD methods. Recent research shows that integrating DGs into the power grid poses a threat to conventional protections, particularly in complex systems with multiple protection devices. As a result, coordinating the Overcurrent Relays (OCRs) can be difficult, leading to false tripping and unnecessary service interruptions [54].
The Differential Relay (DR) is a commonly employed protection that depends on the differential current value to address the issue of bidirectional power flow in a DS. Nevertheless, if the current transformers are saturated or not properly configured, the DR may

Comparison with Different Protection Methods
The system presented in Figure 6 was used to test the effectiveness of the conventional overcurrent and differential protections [53] and compare them with the MSOGI-THD and SOGI-THD methods. Recent research shows that integrating DGs into the power grid poses a threat to conventional protections, particularly in complex systems with multiple protection devices. As a result, coordinating the Overcurrent Relays (OCRs) can be difficult, leading to false tripping and unnecessary service interruptions [54].
The Differential Relay (DR) is a commonly employed protection that depends on the differential current value to address the issue of bidirectional power flow in a DS. Nevertheless, if the current transformers are saturated or not properly configured, the DR may experience tripping problems. Moreover, the settings of both the OCR and DR must be updated due to the grid's changing circumstances [55][56][57].

Fault Detection Performance Analysis for the Protection Approaches
A 3PH-G (ABC) fault with 2DGs located at DL1 in F1 is considered in this section at 0.2 s with a fault resistance of 0.001Ω to compare the performance of the OCR and DR with the results of the proposed THD methods presented in Section 5.1.1. Tables 3 and 4 show the settings used for the OCR and DR, respectively. It's worth noting that the OCR was designed with the IEEE extremely inverse curve, which ensures the fastest disconnection in comparison with the other curves [58]. Table 3. OCR Settings.

Parameter Value
Differential current (pu) 1.08 Biased characteristic (K) 0.5 Current transformer ratio (CT) 500:1 The performance of the OCR and DR protection methods is shown in Figure 17. The OCR and DR had been designed to take action after a delay of 100 ms of the fault detection for simplicity of comparison. experience tripping problems. Moreover, the settings of both the OCR and DR must be updated due to the grid's changing circumstances [55][56][57].

Fault Detection Performance Analysis for the Protection Approaches
A 3PH-G (ABC) fault with 2DGs located at DL1 in F1 is considered in this section at 0.2 s with a fault resistance of 0.001Ω to compare the performance of the OCR and DR with the results of the proposed THD methods presented in Section 5.1.1. Tables 3 and 4 show the settings used for the OCR and DR, respectively. It's worth noting that the OCR was designed with the IEEE extremely inverse curve, which ensures the fastest disconnection in comparison with the other curves [58].  The performance of the OCR and DR protection methods is shown in Figure 17. The OCR and DR had been designed to take action after a delay of 100 ms of the fault detection for simplicity of comparison.  Table 5 shows the fault detection and clearing times for each method. Note that the DR and OCR of DL1 had cleared the fault at 0.315 s and 0.4750 s, respectively, which, considering the 100 ms, means that the DR takes 15 ms to detect the fault, whereas the OCR takes 175 ms. Figure 17 indicates that utilizing the OCR protection method is not advisable, as it allows the fault to continue for several cycles, posing a risk of equipment  Table 5 shows the fault detection and clearing times for each method. Note that the DR and OCR of DL1 had cleared the fault at 0.315 s and 0.4750 s, respectively, which, considering the 100 ms, means that the DR takes 15 ms to detect the fault, whereas the OCR takes 175 ms. Figure 17 indicates that utilizing the OCR protection method is not advisable, as it allows the fault to continue for several cycles, posing a risk of equipment damage. It can be inferred that the DR, MSOGI-THD, and SOGI-THD methods can rapidly isolate the faulted DL in the system. Nevertheless, the SOGI-THD protection approach demonstrated the fastest response, clearing the fault at 0.306 s, meaning that it takes only 6 ms for detection, as shown in Figure 13b.

Different Fault Resistance Protection Tests
Changing the fault impedance can pose a challenge in detecting faults in distribution networks for many protection techniques. Therefore, this section explores the behavior of the protection methods when an unsymmetrical 1PH-G fault with the 2DGs occurs at the F2 location in DL3 at 0.2 s; see Figure 6. The fault resistance is changed to r = 6 Ω. Figure 18 shows the performance of the DR, MSOGI-THD, and SOGI-THD methods during the fault. It can be observed that the SOGI-THD is the fastest at clearing the fault. The SOGI-THD, MSOGI-THD, and DR clear the fault at 0.3068 s, 0.3092 s, and 0.316 s, respectively. damage. It can be inferred that the DR, MSOGI-THD, and SOGI-THD methods can rapidly isolate the faulted DL in the system. Nevertheless, the SOGI-THD protection approach demonstrated the fastest response, clearing the fault at 0.306 s, meaning that it takes only 6 ms for detection, as shown in Figure 13b.

Different Fault Resistance Protection Tests
Changing the fault impedance can pose a challenge in detecting faults in distribution networks for many protection techniques. Therefore, this section explores the behavior of the protection methods when an unsymmetrical 1PH-G fault with the 2DGs occurs at the F2 location in DL3 at 0.2 s; see Figure 6. The fault resistance is changed to = 6 Ω. Figure  18 shows the performance of the DR, MSOGI-THD, and SOGI-THD methods during the fault. It can be observed that the SOGI-THD is the fastest at clearing the fault. The SOGI-THD, MSOGI-THD, and DR clear the fault at 0.3068 s, 0.3092 s, and 0.316 s, respectively.  Figure 18. currents for the grid during a 1PH-G fault with 2DGs at DL3 using (a) SOGI-THD, (b) MSOGI-THD, and (c) DR protection methods.  The protection methods were tested under various fault resistances, types, locations, and DG penetrations. The results consistently showed that in all scenarios, the SOGI-THD method cleared the fault faster than the other ones, as seen in Figure 19. The protection methods were tested under various fault resistances, types, locations, and DG penetrations. The results consistently showed that in all scenarios, the SOGI-THD method cleared the fault faster than the other ones, as seen in Figure 19.

Different THD Methods Computational Burden Assessment
The THD protection method explained in [24] is compared with the SOGI-THD and MSOGI-THD methods to validate and emphasize the advantages of the proposed approaches. In [24], a Fast Fourier transform (FFT) was used to achieve the required harmonics and THD for each phase of the voltage signal and define the protection algorithm.
According to [36,43], the SOGI-FLL and FFT methods were implemented into a DSP, and the computational burden was computed in the number of processor cycles (c) used to compute the SOGI, the FLL, and the FFT blocks. In [43], it was reported that the SOGI needs 149c, the FLL 49c, and the FFT 7019c per phase. Table 6 shows the number of SOGIs and the processor cycles required for executing each method. It can be seen from Table 6 that the SOGI-THD method requires the lowest number of SOGIs and processor cycles when it is executed in the DSP. This justifies the fact that it is much more computationally affordable without affecting the THD calculation or the operation of the protection system. Table 6. Number of SOGIs and cycles required to compute the THD of the three-phase system in each method.

Method
Number of SOGIs Number of Cycles (c) SOGI-THD 3 447 MSOGI-THD  12  1788  FFT-THD  Not applicable  21,057 Additionally, the THD protection techniques proposed in this study are compared alongside other methods in similar conditions. Table 7 compares the methods based on their response speed, accuracy, cost of implementation, the use of inverter-based systems, communication needs, and the grid configuration used for testing. Table 8 presents a summary of the tripping times of these methods and their respective advantages and drawbacks. It can be concluded that the SOGI-THD approach may be a feasible solution for ensuring reliable and fast protection under various conditions when using communication lines.

Different THD Methods Computational Burden Assessment
The THD protection method explained in [24] is compared with the SOGI-THD and MSOGI-THD methods to validate and emphasize the advantages of the proposed approaches. In [24], a Fast Fourier transform (FFT) was used to achieve the required harmonics and THD for each phase of the voltage signal and define the protection algorithm.
According to [36,43], the SOGI-FLL and FFT methods were implemented into a DSP, and the computational burden was computed in the number of processor cycles (c) used to compute the SOGI, the FLL, and the FFT blocks. In [43], it was reported that the SOGI needs 149c, the FLL 49c, and the FFT 7019c per phase. Table 6 shows the number of SOGIs and the processor cycles required for executing each method. It can be seen from Table 6 that the SOGI-THD method requires the lowest number of SOGIs and processor cycles when it is executed in the DSP. This justifies the fact that it is much more computationally affordable without affecting the THD calculation or the operation of the protection system.
Additionally, the THD protection techniques proposed in this study are compared alongside other methods in similar conditions. Table 7 compares the methods based on their response speed, accuracy, cost of implementation, the use of inverter-based systems, communication needs, and the grid configuration used for testing. Table 8 presents a summary of the tripping times of these methods and their respective advantages and drawbacks. It can be concluded that the SOGI-THD approach may be a feasible solution for ensuring reliable and fast protection under various conditions when using communication lines.  Possibility of communication failures [37] SOGI-THD 6-8.5 ms The same merits of the MSOGI-THD. In addition, it is faster and more reliable with a higher THD threshold, requires fewer SOGIs, and has a lower computational burden.
Possibility of communication failures.

Discussion
This study aims to compare the effectiveness of two THD-based protection methods, MSOGI-THD and SOGI-THD, against traditional overcurrent and directional protection methods in a radial distribution system (DS) that includes distributed generation (DG) connections. The analysis takes into account various fault scenarios, types of faults, fault resistances, DG penetrations, and fault locations. The impact of communication delays and the presence of distorting harmonics on the effectiveness of these protection methods have also been examined.
The results indicate that both MSOGI-THD and SOGI-THD protection methods are highly effective in detecting faults in distribution systems with varying levels of DG penetration and fault resistance. In addition, these methods are capable of differentiating between actual faults and harmonic distortions that may exist in the system prior to a fault, which helps prevent unnecessary protective actions. Moreover, the results showed that communication delays had little impact on the performance of the protection process, which points out that the THD-based protection methods may be well suited for use in large-scale power grids where communication delays are common.
In addition, the comparison found that the THD-based protection methods performed better than traditional OCR and DR protection methods in terms of response time and efficiency in isolating faulted sections of the distribution system. The OCR protection method takes longer to detect faults, which can delay the clearing of the fault and increase the risk of equipment damage.
Finally, the work also tried to compare the performance of an FFT-THD-based protection method with the MSOGI-THD and SOGI-THD methods. An FFT-THD method is found to be more computationally expensive and requires a significant number of processor cycles. In contrast, the SOGI-THD method is found to significantly reduce fault detection time and provide accurate fault detection while requiring a low computational burden, making it a possible option to be used for improving the reliability and safety of power grids with DGs.

Conclusions
A comparative study of two THD protection techniques for DS has been proposed in this paper. These techniques are based on the THD of the grid voltages, the estimated amplitude voltages, and the zero-sequence component to define an algorithm to detect, identify, and isolate faults in the grid. The first method employs an MSOGI-THD to obtain the estimated variables, whereas the second method employs a SOGI-THD for the same purpose.
The simulation results in Section 5 showed that, when compared to MSOGI-THD and traditional OCR and DR protection methods, the SOGI-THD technique is the fastest in detecting and isolating faults. The method proved to be effective under various circumstances, including different fault types, DG penetration levels, fault resistances, and fault locations within the studied network.
The study findings demonstrate that the SOGI-THD protection method achieves a faster fault detection response in the protection system. The time response of the SOGI-THD method has been measured to be between 6 and 8.5 ms across all cases examined. In contrast, the MSOGI-THD method showed a fault detection time response that was between 7 and 10 ms. The DR approach exhibited a longer fault detection time, exceeding 15 ms. The OCR approach was found to be undesirable, given its fault detection time of more than 150 ms.
In terms of fault detection accuracy, the SOGI-THD approach outperforms OCRs and DRs. The SOGI-THD has a high detection threshold and the ability to operate regardless of the harmonic distortion of the grid at the moment of the fault. Conversely, OCRs have limited accuracy and may not accurately detect faults under certain conditions, such as high impedance faults. While DRs are more accurate than OCRs, they may require more complex algorithms and higher computational resources.
Compared to other methods, the SOGI-THD approach requires the least computational burden from digital processors. It employs only three SOGIs and 447c for execution, involving few mathematical operations to obtain the estimated variables. In contrast, the MSOGI-THD requires 12 SOGIs and 1788c for execution, while the FFT-THD needs 21057c. This indicates that SOGI-THD is the most computationally efficient method.
Furthermore, the practical advantages of the SOGI-THD method were demonstrated in terms of the trip time response speed of the PD when compared to other protection methods that operate under similar conditions.
In future work, the method's findings can be verified experimentally to ensure their reliability. Moreover, the feasibility of implementing the methods in more complex DSs, such as ring DS, that feature varying high voltage levels can be investigated. Furthermore, further research could be directed toward minimizing the communication requirements within the system.

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