Content area

Abstract

Guaranteeing the effective coordination of directional overcurrent relays (DOCRs) within microgrids (MGs) is a complex nonlinear problem due to bidirectional power flows, varying fault current levels, and the need for adaptive operation across multiple grid configurations. To address this challenge, this paper proposes a hybrid matheuristic approach combining a Biased Random-Key Genetic Algorithm (BRKGA) with Mixed-Integer Linear Programming (MILP). This formulation treats the selection of relay characteristic curves as a decision variable, allowing for simultaneous optimization of time multiplier settings (TMS), plug setting multipliers (PSM), and curve types. The BRKGA handles the global search, while the embedded MILP decoder performs exact optimization under fixed conditions. The proposed BRKGA–MILP method was tested on the IEC benchmark microgrid under multiple operating modes. Compared with conventional MILP-based coordination, it achieved up to 18.31% reduction in total relay operating times (11.81% on average) while maintaining proper coordination time intervals (CTI). Relative to previous heuristic and hybrid approaches, the method improved protection speed by up to 14.87%. These results indicate that the proposed framework effectively enhances coordination performance in adaptive microgrid protection, particularly under bidirectional power flows and varying fault current levels.

Full text

Turn on search term navigation

1. Introduction

1.1. Motivation

The transition toward sustainable energy has emerged as a top priority for governments and companies worldwide, driven by the need to reduce carbon emissions and ensure a sustainable energy supply. However, decreasing dependence on fossil fuels is a complex challenge that requires significant changes in the composition of the energy mix [1,2]. This transition demands not only major investments in renewable energy but also a restructuring of generation and consumption systems to ensure an efficient and equitable shift [3,4].

Microgrids (MGs) play a crucial role in the success of the energy transition as they help integrate multiple renewable and conventional generation sources while optimizing demand management. By combining various generation and storage technologies, MGs facilitate greater penetration of renewable energy and improve efficiency in electricity distribution, reducing losses. To guarantee the secure and dependable operation of MGs, protection coordination serves a fundamental role. Protecting these MG systems has been a challenge because traditional distribution networks operate radially, without power backflows, whereas MGs can experience multi-directional power flows. Additionally, the presence of multiple generation sources increases short-circuit currents, further complicating the design of effective protection schemes. Furthermore, MGs can operate in different modes, either connected to the main grid or in islanded mode, making the coordination of overcurrent relays even more complex [5].

Recent studies on protection coordination in MGs have focused on seeking strategies to enhance their security and stability [6]. In this regard, this paper proposes a hybrid method that combines exact techniques with heuristic methods, developing a matheuristic for the efficient coordination of directional overcurrent relays (DOCRs) within MGs. The proposed approach seeks to improve the reliability of the microgrid, ensuring an efficient response to faults under various operational scenarios.

1.2. Literature Review

DOCRs have emerged as a reliable and widely used solution for microgrid protection due to their ability to detect and isolate faults effectively [7]. The study by [7] presents an analysis of DOCR coordination in distribution systems with distributed generation, emphasizing the importance of adaptive settings to maintain selectivity under varying fault levels. The authors validated their approach using a medium-voltage test feeder, showing that adaptive coordination significantly reduces total operating time during islanded conditions.

However, ensuring proper integration and optimal coordination within MGs remains an ongoing research challenge. This difficulty stems from the complex operating conditions of MGs, which involve varying fault currents, bidirectional power flows, and transitions between grid-connected and islanded modes. Consequently, researchers continue to explore advanced techniques to enhance DOCRs coordination, including adaptive protection schemes, intelligent optimization algorithms, and real-time control strategies.

Given that protection coordination in MGs is inherently nonlinear and non-convex, many approaches rely on metaheuristic strategies. Metaheuristics are optimization algorithms, often inspired by natural processes, that enable the approximation of solutions to complex problems within reasonable computation times [8]. Ref. [8] provides a comprehensive overview of metaheuristic optimization in electrical engineering, describing how algorithms such as PSO and GA are applied to multi-objective problems in protection coordination, power flow, and unit commitment. The paper highlights their efficiency in finding near-optimal solutions within limited computational time, especially in non-convex search spaces.

Ref. [9], a Particle Swarm Optimization (PSO) method is introduced to achieve optimal coordination of DOCRs in systems with high penetration of distributed generation, focusing on improving performance and minimizing the total time dial setting. In ref. [9] applied PSO to a 9-bus distribution network with multiple DG units, optimizing the TMS values and curve selection of each relay. Their results demonstrated a reduction in total operating time of up to 20% compared to traditional linear programming, proving the suitability of PSO for nonlinear relay coordination problems.

Other approaches for optimal protection coordination include Genetic Algorithms (GAs) [10], Ant Colony Optimization (ACO) [11], Firefly Algorithm (FA) [12], Honey Bee Algorithm (HBA) [13], Bat Algorithm (BA) [14], Gravitational Search Algorithm (GSA) [15], Differential Evolution (DE) [16], and Harmony Search (HS) [17], among others. Each of these algorithms introduced distinct search and adaptation mechanisms. For example, ref.  [11] validated ACO on a 14-bus IEEE system showing superior convergence speed, while [12] used FA for time–current coordination achieving better global exploration. Ref. [13] demonstrated that HBA improved reliability indices in networks with multiple DGs. These studies collectively highlight the adaptability of swarm-based algorithms in complex protection systems.

Although metaheuristics provide robust global search capabilities and are highly effective in navigating and exploring complex solution spaces, they often suffer from slow convergence rates and risk getting trapped in local optima. To address these limitations, hybrid optimization techniques have emerged as a powerful alternative. By combining the complementary strengths of different algorithms, hybrid methods enhance solution quality, improve convergence reliability, and increase computational efficiency.

In [18], a GA was used to initialize values for the Time Multiplier Setting (TMS) and Plug Setting (PS), which were then refined using Nonlinear Programming (NLP) to yield globally optimal solutions. In ref. [18] implemented their GA–NLP hybrid on IEEE 8- and 14-bus systems, reporting faster convergence and lower total operating times than standalone GA or NLP methods. Their work demonstrated the benefit of combining stochastic global exploration with deterministic refinement for DOCR coordination.

Similarly, a hybrid GA and an Efficient Heuristic Algorithm (EHA) were proposed in [19], where GA selected PS and curve types, while the EHA efficiently solved the remaining linear optimization task. In ref. [19] tested the GA–EHA hybrid on a 30-bus distribution system. The hybrid method reduced the computational burden by nearly 40% while maintaining high coordination accuracy, confirming that heuristic refinement can enhance traditional GA performance.

Other hybrid GA implementations are reported in [20,21]. Ref. [20] introduced a multi-objective GA framework balancing coordination time and relay sensitivity, validated on a microgrid with inverter-based DGs. Their model outperformed conventional GA in maintaining coordination selectivity under variable fault scenarios.

Multiple hybrid approaches have been developed to improve the performance of PSO when solving the optimal protection coordination problem. In [22,23], a Nelder-Mead-enhanced PSO was introduced, where PSO performed global exploration while the simplex method accelerated local convergence. Refs. [22,23] validated their hybrid PSO–Nelder–Mead approach on standard IEEE test systems, showing that the inclusion of simplex-based refinement reduced convergence iterations by 30% and improved the precision of TMS values.

Further enhancements include PSO combined with DE [24], PSO-GA hybrids [25], and PSO with Local Search Algorithms to ensure feasibility and constraint enforcement [26,27]. These hybrid models combine the fast global search of PSO with the fine-tuning capabilities of evolutionary or local search methods. For instance, ref. [26] tested PSO-LSA on a 33-bus distribution network and achieved better constraint satisfaction, while [25] demonstrated enhanced solution diversity with a PSO-GA model.

Simulated Annealing (SA) is effective at escaping local optima, making it a strong candidate for hybridization with swarm-based algorithms. In [28], SA was embedded within PSO to enhance global search capabilities. The author of [28] implemented PSO–SA for relay coordination in a DG-integrated system, confirming through MATLAB simulations that the hybrid method reached stable convergence and reduced the number of infeasible solutions by 25%.

Similarly, FA and Whale Optimization Algorithm (WOA) have been hybridized with GAs and SA to enhance convergence and solution quality [29,30]. In [31], a hybrid Harris Hawks Optimization (HHO) with Sequential Quadratic Programming (SQP) demonstrated improved performance for constrained DOCR coordination problems. HS, when combined with SA as in [32], benefited from improved exploration and escape from local optimal solutions. The author of [31] evaluated the HHO–SQP model on a 6-bus microgrid, showing that it achieved 15% faster convergence compared to standalone HHO, while Korashy et al. (2019) demonstrated that HS–SA hybrids provided better selectivity margins under varying DG outputs.

Other notable hybridizations include the combination of Water Cycle Algorithm (WCA) with Moth-Flame Optimization (MFO) to leverage both exploration and exploitation properties [33]. Gradient-Based Optimizer (GBO), while efficient in local search, suffers from premature convergence. This limitation was addressed by combining it with Success-History based Adaptive Differential Evolution with Linear Population Size Reduction (LSHADE), enhancing diversity and global exploration [34]. Additionally, the hybridization of WOA and Grey Wolf Optimizer (GWO), termed HWGO, provided a synergistic balance between exploration and exploitation phases [35]. The Invasive Weed Optimization (IWO) algorithm was also improved by incorporating Sequential Quadratic Programming (SQP) for better convergence and computational efficiency on standard IEEE test systems [33]. These recent hybrid approaches, such as HWGO and LSHADE–GBO, have demonstrated outstanding robustness when applied to IEEE 14- and 33-bus test systems, consistently outperforming single metaheuristics in both convergence rate and CTI satisfaction.

Building upon the concept of hybridization between metaheuristics, a more structured approach emerges with matheuristics. Matheuristics represent a sophisticated class of hybrid algorithms that merge exact mathematical programming—such as Linear Programming (LP) or Mixed-Integer Linear Programming (MILP)—with metaheuristics like GAs, PSO, or FA. These strategies aim to harness the guaranteed convergence of mathematical models while benefiting from the flexible exploration capabilities of metaheuristics. For instance, ref. [36] proposed a GA-LP hybrid where GA coded the PS values, and LP calculated TMS for various topologies. In refs. [37,38], PSO determined PS values while LP computed TMS settings. Likewise, hybridizations of LP with ABC [30], Biogeography-Based Optimization (BBO) [39], FA [40], and SA [41] demonstrated superior convergence speed and accuracy compared to standalone metaheuristics. Literature reviews on the topic of protection coordination methodologies in MGs can be consulted in [8,42]. These matheuristic frameworks bridge theoretical exactness and algorithmic adaptability, offering optimized DOCR settings validated in IEEE standard test systems and actual microgrid case studies. Their growing use signals a shift toward mixed-integer formulations enhanced by AI-based optimization.

As the literature review shows, many metaheuristics have been applied to solve DOCR coordination problems; however, very few methods employ matheuristic algorithms that combine mathematical formulations and heuristics to improve solution quality. Advancing beyond existing metaheuristic approaches, this study introduces a framework that integrates MILP with the Biased Random Key Genetic Algorithm (BRKGA). As pointed out in [43], in BRKGAs, mathematical programming is frequently applied to the decoding process. Matheuristics leverage on the mutual benefits of mathematical models and metaheuristic techniques improving both solution quality and convergence rate, albeit at the cost of higher computational effort [44]. Within this hybrid structure, the BRKGA performs global exploration across non-convex solution spaces, while MILP refines solutions with a exact solution decoder. A key innovation lies in treating protection curve selection as a decision variable, enabling adaptive curve assignment within the optimization process—an enhancement over traditional methods with fixed curve types.

1.3. Contributions and Paper Organization

As evidenced in the literature review, numerous approaches have been proposed for the optimal coordination of DOCRs in microgrids. Nonetheless, hybrid matheuristic algorithms have typically relied on a single standard curve as a fixed parameter. In contrast, the hybrid algorithm proposed in this study introduces curve selection as a decision variable, which enhances the performance of protection coordination. Additionally, the use of BRKGA has not been reported among the metaheuristics applied to this topic—let alone within hybrid optimization frameworks. In this context, the main features and contributions of this paper are:

A non-linear programming model for optimal protection coordination in MGs is proposed and evaluated, enhancing modeling accuracy and flexibility compared to previously reported approaches.

A hybrid matheuristic based on a BRKGA is developed, embedding an MILP model within the decoder. This hybridization allows the MILP to perform an exact evaluation under fixed conditions. Such an approach—previously unexplored in BRKGA applications for power systems—enables a more effective solution process.

The proposed hybrid algorithm incorporates an exact optimization model. Specifically, the proposed MILP formulation includes curve selection as a decision variable, enabling the optimal joint selection of curve type, plug setting multiplier (PSM), and time multiplier setting (TMS). Unlike previous studies, where the curve type is fixed and only the TMS is optimized, this approach provides greater flexibility in protection coordination.

The performance of the proposed approach was evaluated against existing methods using standard benchmark systems from the literature, demonstrating improved results in terms of tripping time.

The rest of the paper is organized as follows. Section 2 provides a methodological framework of coordination studies. Section 3 presents the mathematical modeling of the protection coordination problem. Section 4 describes the proposed BRKGA matheuristic approach. Section 5 details the IEC benchmark microgrid used for testing. Section 6 presents and analyzes the computational results. Lastly, Section 7 presents the conclusions and outlines avenues for future work.

2. Methodological Framework of Overcurrent Coordination Studies

In compliance with IEEE 242-2001 [45], clasical overcurrent coordination is performed through an iterative trial-and-error methodology, wherein the time-current characteristic (TCC) curves of protective devices connected in series are graphically represented on log-log scales to achieve selective coordination. The determination of device settings involves a trade-off between two competing objectives: the provision of optimal equipment protection and the assurance of maximum service continuity. Consequently, complete selective coordination cannot always be guaranteed across all system configurations. This procedure is applicable only to radial systems without distributed generation along the feeder. The following procedure presents a systematic approach for constructing time-current coordination plots and establishing the settings of time-overcurrent protective devices within the electrical system.

Figure 1 illustrates the procedural workflow for conducting an overcurrent coordination study, as structured in accordance with IEEE Std 242-2001 (Buff Book) [45], Section 15.7. The process begins by defining the system configuration, which specifies the arrangement and operating mode of the electrical network, followed by the collection of essential system data, encompassing quantitative parameters such as impedance values, voltage levels, load profiles, and transformer characteristics. Once the data are compiled, suitable electrical power analysis software is selected to support the coordination study.

Subsequently, a short-circuit analysis is performed to determine the available fault currents at key buses and terminals throughout the system. Based on these results, protective devices—such as circuit breakers, relays, and fuses—are selected or verified, ensuring appropriate ratings and types. TCC curves for all overcurrent protective devices are then plotted sequentially.

Once plotted, the TCC curves are analyzed to assess time selectivity and coordination between devices. If selectivity is not achieved, iterative adjustments to the settings or types of protective device are performed, and the coordination analysis is re-evaluated. Upon successful coordination, final protective device settings are established and documented, along with the corresponding TCC plots, analysis rationale, and recommendations. The process concludes with the implementation of the approved settings in the field and the establishment of a monitoring plan to accommodate future system modifications.

Figure 2 presents a comprehensive visualization of the protection coordination methodology applied to a power system. The illustration is divided into three distinct parts (labeled a, b, and c), each corresponding to a specific aspect of the protection analysis and design process. It should be emphasized that this methodology assumes a radial system configuration without distributed generation along the feeder, where fault current decreases with distance from the source. In systems with distributed generation, bidirectional fault currents, and variable short-circuit levels can significantly complicate coordination, requiring more advanced protection schemes.

Part (a) displays a one-line diagram of a power system, where buses are interconnected through transmission lines and ultimately connected to a central power source or substation. Square boxes in the diagram represent DOCRs assigned to protect each line. A fault scenario is indicated by a lightning bolt symbol, annotated with the label Iif, denoting the fault current observed at the relay locations. This current serves as a reference for evaluating the response of both primary and backup relays.

Part (b) presents a plot illustrating the relationship between fault current magnitude and distance from the source. The vertical axis denotes the fault current, whereas the horizontal axis represents the electrical distance. The plot reveals a decreasing trend in fault current as distance increases. This inverse relationship is attributed to the rising equivalent impedance seen by the source as the fault location moves farther away, leading to a decrease in the resulting fault current.

Part (c) shows an inverse TCC curve diagram, which constitutes the core of the coordination analysis. The horizontal axis represents current in kiloamperes, while the vertical axis denotes operating time in seconds. The curve labeled Iif corresponds to the simulated fault current for the given scenario. Two inverse-time protection curves are shown in Figure 2c: the red curve represents the primary relay (Rmain), and the blue curve represents the backup relay (Rbackup). The intersection of the fault current (green line) with the red curve determines the operating time tif for the main relay, while its intersection with the blue curve yields the operating time tjf of the backup relay. Here, the subscript i refers to the primary relay, and j to the backup relay.

The difference between these two time values defines the Coordination Time Interval (CTI), a critical metric for ensuring selective operation of protection layers. An adequate CTI allows the primary relay to clear the fault before the backup relay operates, which maintains coordination and prevents unnecessary disconnections.

Minimizing the fault clearing time is a key objective of protection coordination. Faster fault clearance reduces the amount of energy released during the disturbance, thereby mitigating thermal and mechanical stress on equipment and limiting the risk of damage to critical infrastructure. Additionally, shorter clearing times contribute to improved system stability and enhance the overall safety of both personnel and equipment. Therefore, the protection coordination strategy must aim to achieve the fastest possible operation of protective devices, particularly the primary relay, without compromising coordination requirements with backup devices.

Figure 2 thus provides a visual summary of the methodology outlined in the preceding section, illustrating the process from system modeling and fault simulation to relay coordination. As stated in IEEE Standard 242 [45], this procedure reflects the iterative nature of protection coordination in practice, aimed at achieving both sensitivity and selectivity under diverse fault conditions.

The following overview will discuss several advanced functionalities that improve the coordination and reliability of protection systems in contemporary power networks. Among these are adaptive overcurrent protection, which dynamically modifies relay settings based on real-time system changes, as well as centralized and decentralized control schemes, each offering distinct advantages in decision-making and fault management; and fundamental protection principles, such as sensitivity and selectivity, which are essential for ensuring prompt and reliable fault clearance.

This study explicitly considers adaptive and directional protection schemes to improve the reliability and selectivity of fault detection in microgrids. Adaptive protection is incorporated through the definition of operating modes (OMs) that model changes in the MG configuration, such as the disconnection or activation of DGs and the transition between islanded and grid-connected operation. Adaptive overcurrent protection enhances system resilience by dynamically adjusting relay settings—such as pickup current and time delay—in response to real-time variations in load, network topology, or DG output. Given that most commercial relays have four sets of settings, this study considers four OMs. Directional protection introduces an additional layer of selectivity by analyzing the phase relationship between voltage and current to accurately determine the direction of fault currents. This capability is especially crucial in MG environments characterized by bidirectional power flows and the possibility of islanded operation, where traditional nondirectional overcurrent protection may fail to correctly identify the fault source or coordinate with other devices.

Protection control architectures can be centralized, decentralized or hybrid. In a centralized scheme, a supervisory unit collects system-wide data and issues coordinated commands to relays, achieving optimal decision-making at the cost of increased communication dependency. Decentralized schemes delegate decision logic to local controllers embedded in each device, improving response speed and reducing vulnerability to communication failures but requiring careful design to maintain global selectivity. Hybrid approaches seek to combine centralized oversight with local autonomy, balancing coordination quality and operational robustness.

Underpinning these features and architectures are four key coordination principles: selectivity, sensitivity, speed and reliability. Selectivity guarantees that only the relay nearest to a fault operates first, thereby isolating the faulted section without interrupting healthy portions of the network. Sensitivity ensures that relays detect faults of minimum magnitude across their designated zones. Speed aims to minimize fault clearing time, reducing the energy released during faults and limiting thermal or mechanical stress on equipment. Reliability reflects the ability of relays and the coordination scheme to perform correctly under all defined operating and fault scenarios.

Traditional protection coordination schemes, designed for radial distribution systems with unidirectional power flow and predictable fault currents, are inadequate for MGs due to their unique operational characteristics. MGs often operate in both grid-connected and islanded modes, leading to bidirectional power flows and varying fault current levels. The integration of Distributed Energy Resources (DERs), such as inverter-based generators, results in lower fault current contributions compared to conventional systems, challenging the sensitivity and selectivity of overcurrent protection devices [46]. Additionally, the dynamic topology and intermittent nature of DERs can cause coordination issues among protective devices, leading to potential misoperations or failure to isolate faults effectively [47]. These complexities necessitate the development of adaptive and intelligent protection schemes tailored to the specific needs of MGs. As a result, there has been a growing research interest in advanced protection coordination strategies, including optimization-based methods, adaptive relaying, and communication-assisted schemes, to enhance protection reliability and ensure secure microgrid operation under diverse scenarios [48].

These advanced features are incorporated into a mathematical optimization framework to reduce the overall operating times of DOCRs while ensuring selectivity between primary and backup relays through the enforcement of CTI constraints. Additionally, the model includes constraints that ensure sufficient sensitivity to specified fault currents and validate the protection system’s performance across a range of fault scenarios, thereby guaranteeing the reliability of the coordination scheme. By integrating adaptive and directional protection capabilities, OMs adjust protection parameters to accommodate the changing conditions of MG, enhancing both the security and operational resilience of the system. The following section will provide a detailed mathematical formulation of these concepts, outlining the optimization model and its application in the coordination of protection within MGs.

3. Mathematical Modeling of the Protection Coordination

This section presents the mathematical formulation of the protection coordination problem addressed in this study. The proposed model is designed to determine optimal settings for DOCRs in an MG environment by minimizing total relay operating times, subject to selectivity and operational constraints. The formulation begins by defining the sets and parameters, and then presents the key equations that characterize the behavior of the protection system.

The main goal of the model is to minimize the total operating time of the relays through the coordinated selection of their TMS and associated standard inverse-time characteristic curves (SCCs). Selectivity between primary and backup relays is enforced using a fixed CTI of 0.3 s, consistent with values commonly adopted in the specialized literature [49,50,51]. This fixed CTI allows for benchmarking and facilitates comparison with existing methodologies.

3.1. Sets

The model is structured using three main sets. Set R contains all the DOCRs deployed within the system. Set F comprises all relevant fault conditions detected by the relays, including variations due to different network configurations or OMs. Set C includes the available standard time-current characteristic curves, as defined by IEC and IEEE standards, from which the optimization model can select appropriate profiles for each relay.

3.2. Parameters

In the proposed model, the Plug Setting Multiplier (PSMif), defined in Equation (1), quantifies the intensity of a fault f as perceived by relay i. It is calculated as the proportion of Iif to the relay’s pickup current (Ipickupi), plus a small design constant ηi that allows for a reduction in relay operating times.

(1)PSMif=IifIpickupi+ηiiR,fF

To linearize the relay operating time function for its use in a mixed-integer optimization framework, the parameter βfic is introduced, as shown in Equation (2). This expression is derived from the standard time-current characteristic equation, where Ac and Bc are curve-specific coefficients for each cC. These constants are established by international standards, specifically IEC 60255-3 and IEEE C37.112 [52], which include some parameters for standard inverse-time overcurrent characteristics used in relay coordination studies:

(2)βfic=AcPSMifBc1iR,fF,cC

To enforce selectivity between relays, the binary parameter BUij is defined in Equation (3) to model backup relationships. If relay j is designated as the backup for relay i, then BUij=1; otherwise, BUij=0. This relationship is asymmetric by nature, consistent with the selectivity principle, where BUijBUji:

(3)BUij=1,ifrelayjisthebackupforrelayi0,otherwisei,jR

4. Proposed BRKGA Matheuristic

Equation (4) defines the operating time tif of relay i in response to fault f, based on its time multiplier setting (TMSi) for each relay i, the selected standard characteristic curve parameters (Ac, Bc), and the fault current magnitude Iif. The term IifIpickupi+ηi represents the plug setting multiplier (PSMif), which accounts for the ratio of the fault current to the pickup current of the relay, with an added small design constant ηi to enhance responsiveness and reduce relay operating time. The expression corresponds to the general form of an inverse-time overcurrent protection curve defined by international standards (IEC 60255-3 and IEEE C37.112 [52]). The numerator contains the curve constant Ac scaled by TMSi, while the denominator captures the nonlinearity of the inverse-time characteristic through the exponent Bc. The right-hand side of the equation introduces the linearized equivalent βfic (see Equation (2)), which simplifies the model representation in the optimization framework.

(4)tif=TMSi·AcIifIpickupi+ηiBc1=TMSi·βfic

As noted in [53], the multiplicative term TMSi·βfic in the relay operating time equation is naturally suited to a MILP formulation. However, studies in the specialized literature [51,54] indicate that allowing targeted variation of the PSM can further reduce operating times. When the PSM is treated as a decision variable, the resulting model becomes non-linear. To address this issue, this study proposes a hybrid solution methodology. In this case, the non-linearities are effectively handled using BRKGA, which generates candidate parameter vectors ηi. Each sampled vector ηi is then fixed within the MILP formulation, which subsequently solves for the optimal TMS and curve type assignments under the given βfic values. By combining the global search capability of the BRKGA with the exact optimization of MILP, the proposed approach yields competitive solutions for the coordination of protection relays.

Figure 3 illustrates a TCC diagram for DOCRs. It highlights the impact of varying ηi, a decision variable in the BRKGA, which is processed by the decoder. The vertical axis represents time in seconds, while the horizontal axis displays current in kiloamperes. Protection coordination curves are typically shown on a log-log scale to clearly demonstrate the coordination between relay curves. The figure depicts the characteristic curves of the main relay (Rmain) and its corresponding backup relay (Rbackup). The fault current (Iif) is the current measured by both the main and backup relays during a fault event. The relay trip time corresponds to the intersection between the main relay curve and the fault current. DOCRs are configured with a pickup current (Ipickup), which determines the threshold above which the overcurrent protection detects faults. The decision variable ηi modifies the main relay curve by shifting it to the right. As shown in Figure 3, the dashed line represents the original relay curve with η=0 (no effect), while the solid blue line represents the modified curve with η>0. This adjustment impacts both the main and backup relay curves. The CTI is a predefined time delay that ensures the backup relay only operates if the main relay fails to clear the fault. Regardless of the ηi value, the system must maintain a defined CTI to uphold the selectivity criterion.

In general, the BRKGA algorithm influences the protection zones of the relays but enhances trip times by optimizing coordination settings. By adjusting η, the relays can tolerate higher admissible overloads, reducing the trip times of both main and backup relays. This optimization improves overall system protection by adapting the pickup current and prioritizing fault-prone areas, thereby reducing fault clearance time and enhancing performance metrics.

4.1. BRKGA Overview

The operational stages of the Biased Random-Key Genetic Algorithm (BRKGA) are outlined in Figure 4. As a variant of genetic algorithms, the BRKGA framework is designed for efficiently solving combinatorial optimization problems. Its problem-independent design utilizes random keys (real-valued vectors) for solution encoding and incorporates biased selection to enhance convergence toward high-quality solutions. The decoder is a key component that translates the random-key vector into a feasible solution for the problem [55,56]. This design allows BRKGAs to be flexible and easily adaptable to different optimization problems [43,57]. As a result, they have been successfully applied to electric distribution network configuration problems with renewable sources [58] and network reconfiguration [59], as well as to multi-objective distribution network reconfiguration for annual energy loss reduction [60], and for the multi-objective optimization of energy storage for transmission loss reduction in meshed power distribution networks [61].

As illustrated in Figure 4, the BRKGA evolves through a sequence of well-defined stages, which are iteratively executed for a predetermined number of generations or until a convergence criterion is satisfied. These stages are described as follows:

Initialization: A population of individuals is randomly generated, where each individual is represented by a real-valued vector whose elements lie within the interval [0,1].

Evaluation: The fitness of each individual is assessed using a decoder function, which maps the random keys to feasible solutions in the problem space and assigns a corresponding objective value.

Selection: A subset of elite individuals with the highest fitness values is selected to guide the search process in subsequent generations.

Biased Crossover: New offspring are produced by performing biased crossover operations between elite and non-elite individuals. The crossover mechanism favors the inheritance of genes from elite parents with a predefined probability.

Replacement: The next generation is constructed by combining elite individuals, newly generated offspring, and a fraction of randomly introduced individuals to preserve genetic diversity and avoid premature convergence.

In ref. [57] the authors extended the standard BRKGA by introducing the MP-BRKGA-IPR framework, which incorporates several enhancements to improve evolutionary dynamics and solution quality. First, the multiple-parent mechanism (MP-BRKGA) strengthens the standard algorithm’s double elitism by allowing offspring to inherit genes from several elite individuals, thereby increasing genetic diversity and exploration capability. Second, the algorithm integrates a local search component through Implicit Path-Relinking (IPR), which intensifies the search around high-quality solutions to accelerate convergence. Additionally, MP-BRKGA-IPR adopts an island model, where a set of independent populations evolve separately and periodically exchange their best individuals, helping maintain diversity and prevent premature convergence. Together, these mechanisms make MP-BRKGA-IPR a more robust and efficient variant than the standard BRKGA.

4.2. Solution Representation

In the BRKGA algorithm, candidate solutions are encoded as random-key vectors. Each individual is represented by a vector (η1,,ηi,,η|R|) of real values within the interval [0,1], where the vector length corresponds to the number of relays iR in the system. Each component (value) ηi encodes a decision variable associated with relay i as stated in Equation (2).

Once the values of ηi are assigned by the BRKGA for all relays, the MILP model is employed to decode the chromosome into a feasible solution and to evaluate the corresponding relay tripping times. This decoding process is illustrated in Figure 5, which represents a low-level matheuristic cooperation strategy, wherein a component of the metaheuristic is replaced by an exact optimization model [62,63]. In the proposed framework, the BRKGA explores the space of random keys, while the MILP-based decoder determines the optimal microgrid protection coordination settings, as detailed in the following section.

4.3. MILP Decoder for BRKGA-Based Protection Coordination

This section provides a comprehensive explanation of the MILP model and its constraints. The objective function is to minimize the total tripping time (tif) of DOCRs in the MG for each fault scenario. The model incorporates decision variables that represent the offline protection settings for each relay, enabling coordination across various Operative Modes (OMs). A centralized controller determines the current OM based on system conditions—such as DG availability or grid connection status—and communicates the corresponding relay settings to minimize operating times. Each OM defines a set of settings, supporting the implementation of adaptive protection by allowing the model to adjust coordination settings dynamically. The decision variables are the time multiplier settings (TMSic) and the selected standard coordination curve (SCC), indexed by cC. The model is formulated as an assignment-like optimization problem. Consequently, a value of TMSic is defined for each relay–curves pair, but its activation is governed by an auxiliary binary variable xic, which ensures that exactly one curve is selected for each relay. This formulation allows the application of upper and lower bounds to the selected TMSic values, while maintaining logical consistency in the curve assignment process. The MILP model is given by Equations (5)–(13) as indicated below.

4.3.1. Objective Function

Equation (5) defines the objective function of the proposed optimization model, which seeks to minimize the cumulative tripping times of all relays across all possible fault scenarios. This ensures that the protection system responds as quickly as possible while maintaining selectivity between primary and backup relays.

(5)miniRfFtif

Here, tif represents the tripping or operating time of relay i when fault f occurs. This variable is computed as a function of the relay’s time multiplier setting (TMSic), the inverse-time curve parameters (Cc), and the coefficient βfic, as explained below.

4.3.2. Relay Operating Time Calculation

The constraint given by Equation (6) defines the relationship between the operating time tif of relay i under fault f and the decision variables that control its performance. The value of tif is obtained by summing the contributions of all possible curve types cC, where βfic acts as a sensitivity coefficient that depends on the PSM and represents the inverse-time characteristic of the relay. The binary variable xic determines whether a given curve c is active for relay i.

(6)tif=cC(βficTMSic+xicCc)iR,fF,βfic0

This expression ensures that the operating time is linearly related to the time multiplier setting TMSic and curve constant Cc. Since tif appears in the objective function, minimizing this variable effectively improves the speed of fault clearance across the network.

4.3.3. Coordination Time Interval Constraint

Coordination between relays is essential to guarantee that primary relays operate before their corresponding backup relays during fault conditions. Equation (7) enforces this principle by maintaining a minimum CTI between each primary relay i and its designated backup relay j. The CTI ensures that backup relays only trip if the primary relay fails, thus preserving selectivity.

(7)cC(βfjcTMSjc+xjcCc)tifCTIi,jR:ij,fF,βfic0,βfjc0,BUij=1

In this case, BUij is a binary parameter that indicates whether relay j acts as the backup of relay i. The constraint is enforced only when both βfic and βfjc are positive—meaning that the corresponding fault current exceeds the relay’s pickup current.

4.3.4. Curve Selection Constraint

Each relay must be assigned exactly one inverse-time characteristic curve from the available set C. The constraint given by Equation (8) enforces this exclusive selection through the binary variable xic, which takes a value of 1 if curve c is assigned to relay i and 0 otherwise:

(8)cCxic=1iR

This constraint ensures that unselected curves are effectively deactivated in the model. As a result, each relay’s operating time tif depends only on the parameters and TMS corresponding to the selected curve.

4.3.5. Upper and Lower Bounds on the TMS

To ensure feasible relay settings, the TMS values are constrained within predefined upper and lower limits. Equations (9) and (10) impose these bounds and link them to the curve-selection variable xic. When a curve is not selected (xic=0), the corresponding TMSic is forced to zero, preventing it from influencing the relay operation.

(9)TMSicTMSubxiciR,cC

(10)TMSicTMSlbxiciR,cC

Typical values for the TMS bounds are TMSlb=0.001 and TMSub=1.0, which define a realistic operational range consistent with relay manufacturer recommendations.

4.3.6. Domain of Decision Variables

The domains and nature of the decision variables are defined in Equations (11)–(13). These expressions restrict the variables to meaningful and physically interpretable values.

(11)tif0iR,fF

(12)TMSic0iR,cC

(13)xic{0,1}iR,cC

In this case, tif is a continuous non-negative variable representing the relay i’s tripping time in seconds for fault f. The variable TMSic is also continuous and adjusts the operating time depending on the selected curve type c. Finally, xic is a binary decision variable ensuring that exactly one curve is assigned to each relay.

As previously mentioned, the MILP is used as a decoder to evaluate the objective function, which minimizes the total relay tripping times for a specific OM. The parameter ηi, generated by the BRKGA for each relay, determines the value of βfic in the MILP formulation according to Equation (2). Although ηi is not directly visible within the MILP, it implicitly guides the optimization by shaping the coefficients of the objective function and constraints. Thus, the interaction between the BRKGA and the MILP occurs directly. The BRKGA explores the solution space by proposing values of ηi, which are then embedded as parameters in the MILP through βfic, enabling the MILP to perform exact evaluation under the given conditions.

The BRKGA was implemented using the Python library for MP-BRKGA-IPR presented in [57]. Likewise, the MILP decoder was modeled in pyomo [64] and solved using Highs [65]. The BRKGA parameters were initially set following the standard configuration proposed by [57]. After conducting preliminary experiments, these values were adjusted to favor faster convergence and shorter running times. The final configuration was as follows: a population size of 20 individuals per population, an elite percentage of 30%, a mutants percentage of 15%, 2 elite parents and 3 total parents per mating, a logarithmic inverse bias function to guide the selection of elite individuals, and 3 independent populations evolving in parallel to preserve diversity. No path-relinking parameters are needed because the Python implementation of [57] does not include this feature.

It is important to note that the proposed BRKGA–MILP algorithm functions as an off-line coordination tool rather than a real-time adaptive protection scheme. The algorithm is designed to perform comprehensive optimization of relay settings under multiple operating modes prior to system operation, thus ensuring optimal coordination in all expected scenarios. While this off-line process involves higher computational times, it provides robust settings that can later be implemented in real-time systems using pre-computed lookup tables or reduced search spaces. Future work will explore strategies such as parallelization, surrogate modeling, and incremental learning to further reduce computational effort and enable faster adaptive updates.

5. IEC Benchmark Microgrid Description

The IEC benchmark MG, described in detail in [53], serves as a standardized platform for evaluating and comparing protection, control, and optimization strategies in MGs under both grid-connected and islanded conditions. Originally introduced by Ustun et al. [66] and subsequently extended by Kar et al. [67], this benchmark is recognized for its capacity to emulate realistic operational environments and foster reproducibility in simulation-based research.

The benchmark includes a utility connection modeled with a rated short-circuit power of 1000 MVA, operating at 120 kV and 60 Hz, it provides a robust grid for fault and stability studies, featuring four DG units:

DG1 and DG2: synchronous generators rated at 9 MVA, 2.4 kV, with an inertia constant of 1.07 s.

DG3: an inverter-based wind farm of 6 MVA (three turbines of 2 MVA each), operating at 575 V, with an inertia constant of 0.62 s.

DG4: a DFIG-based wind farm with a total capacity of 9 MVA (six turbines of 1.5 MVA), also operating at 575 V.

Voltage transformation is managed through four transformers:

TX1 (15 MVA, 120/25 kV) connects the utility to the microgrid.

TX2 and TX3 (12 MVA each, 2.4/25 kV) link synchronous DGs.

TX4 (10 MVA, 575/25 kV) links inverter-based DGs.

All transformers operate at 60 Hz and share common per-unit impedance values, such as a winding resistance R1=0.00375 pu and a leakage reactance X1=0.01 pu.

The benchmark also includes five 30 km radial distribution lines modeled with detailed sequence parameters:

R1=0.413Ω/km,L1=3.32×103H/km,C1=5.01×109F/km,R0=0.1153Ω/km,L0=1.05×103H/km,C0=11.33×109F/km

Protection components include 15 overcurrent relays (R1 to R15), each associated with pickup current thresholds and current transformer ratios. These elements support testing of coordination schemes.

The overall system load is modeled as 22 MW and 10 MVAR, distributed across six loads (L1 to L6), connected via six buses (B1 to B6). The configuration allows for the simulation of faults, load variations, and operational transitions, which are critical for evaluating protection coordination, adaptive control, and resilience strategies.

5.1. IEC Benchmark Microgrid Layout

Figure 6 illustrates the structure of the IEC benchmark microgrid, designed to emulate a typical radial medium-voltage distribution system with multiple DERs. The system connects to the main utility grid through a transformer (TX1) at the Point of Common Coupling (PCC), labeled as Bus 1 (B1), and is protected by the main circuit breaker (CB_MAIN).

From B1, the system branches into several feeders through distribution lines DL1 to DL5, each protected by overcurrent relays (R1 to R15) and sectionalizing breakers (CB_LOOP1 and CB_LOOP2). Each feeder is associated with a set of loads (L1 to L6), buses (B2 to B6), and distributed generators (DG1 to DG4), connected via step-up or step-down transformers (TX2 to TX5) according to their voltage levels.

This layout enables simulation of various fault scenarios, protection schemes, and islanding conditions. The hierarchical configuration of relays and breakers supports studies in coordination, selectivity, and stady-state response under both normal and disturbed operations.

5.2. Three-Phase Fault Characterization

Table 1 presents the fault scenarios considered in this study, labeled F1 through F5. Each fault corresponds to a solidly grounded three-phase short circuit occurring at the midpoint of a designated distribution line, DL-1 to DL-5. The spatial layout and connectivity of these lines within the benchmark IEC microgrid are depicted in Figure 6. The table lists the short-circuit currents calculated under four distinct operative modes (OM1 to OM4), expressed in kiloamperes. These values were derived using the IEC 60909 standard [68] and obtained through simulation in DIgSILENT PowerFactory.

The variation in short-circuit levels across OMs reflects differences in topology, generation availability, and grid interconnection status. Notably, OM4 represents a fully islanded configuration, resulting in significantly reduced fault current magnitudes due to the absence of utility grid contribution. These short-circuit current levels serve as inputs to the optimization model, as they directly affect the calculation of pickup ratios and, consequently, the sensitivity and selectivity of the overcurrent protection settings.

5.3. Adaptive System Configurations

This study defines four representative OMs to reflect the variability of microgrid operating conditions, as summarized in Table 2. These OMs are configured based on the availability of DGs and the status of grid connectivity. The structural arrangement of DG units and their interconnection with the system buses are illustrated in Figure 6.

In Table 2, the designation “enabled” indicates that the associated generation unit is operational and injecting power into the network, whereas “disabled” denotes that the unit is inactive. Throughout all OMs, switches CB_LOOP1 and CB_LOOP2 are assumed to remain open, maintaining a radial configuration. Operating modes OM1 through OM3 represent grid-connected conditions, each differing in the level of DG participation. In contrast, OM4 corresponds to an islanded scenario, where the microgrid functions independently from the main grid. These four configurations were chosen to reflect typical relay capabilities, as most commercial protection relays support up to four distinct settings groups. Furthermore, the last column of Table 2 describes the characteristics of each operative mode in terms of DG availability and Grid Connection (GC) status. This description supports the interpretation of the corresponding tripping times and protection behavior under each configuration. The comparison of tripping times with those reported in other studies was restricted to cases employing similar test systems, since differences in network topology and configuration limit the possibility of direct quantitative comparisons.

6. Computational Results

Table 3 presents a comparative analysis of the total tripping time across different operative modes (OM1–OM4) for various protection coordination strategies. Four operative modes were considered for consistency with commercial relays, which typically support four adaptive setting groups. The studies reviewed span a progression of methodological approaches, ranging from traditional metaheuristic-based coordination to more advanced hybrid optimization frameworks. Early works, such as [49,50], applied a GA to optimize the settings of overcurrent relays. These methods demonstrated the effectiveness of metaheuristics in handling the nonlinear and discrete nature of the coordination problem. Later, refs. [51,54] further explored adaptive coordination by incorporating non-standard characteristics and multiple metaheuristic techniques, including GA, PSO, TLBO, and Shuffled Frog-Leaping Algorithm (SFLA), thus enhancing flexibility in microgrid protection. More recently, ref. [69] introduced a ROMI-PSO algorithm, which explicitly handles both continuous and discrete decision variables, such as time dial settings and characteristic curve types. In contrast to these metaheuristic-based strategies, ref. [53] proposed a MILP model, offering an exact optimization formulation capable of guaranteeing global optima under linearized assumptions. Building on this evolution, the hybrid BRKGA-MILP approach introduced in this work combines the global search capabilities of the BRKGA with the precision of MILP solvers.This hybrid matheuristic framework leverages both the flexibility of metaheuristics and the structural rigor of mathematical programming. Table 3 shows variations in total tripping time across operative modes, highlighting the influence of each optimization methodology on overall protection system performance.

Although numerous studies addressing DOCR coordination have been reported in the specialized literature, not all of them employ the same benchmark system, relay characteristics, or fault scenarios. Consequently, direct numerical comparison of tripping times with all cited references is not feasible. The analysis presented in Table 3 focuses on approaches implemented under equivalent test conditions, ensuring consistency in the comparative evaluation of total operating times.

Table 4 compares the performance of the proposed BRKGA–MILP framework against the MILP baseline and the best-performing method reported in prior studies. For each operative mode (OM1–OM4), absolute tripping times (in seconds) are presented alongside percentage reductions relative to both references. The “MILP (%)” column quantifies the improvement of BRKGA–MILP over the exact MILP formulation, while the “BestPrev (%)” column highlights additional gains with respect to the strongest previous metaheuristic approach. Results show consistent reductions in total tripping time, with improvements of up to 18.31% compared to MILP and up to 14.87% compared to the best previous method, confirming the advantage of the hybrid strategy in microgrid protection coordination.

6.1. Operative Mode 1

Table 5 presents the relay settings for OM1 obtained using both the MILP [53] and BRKGA-MILP approaches. Table 5 includes the assigned values for the variable ηi, the selected standard characteristic curve (SCCi), and the TMS for each relay. It is important to note that in the MILP approach, etai is fixed to zero as it serves as a parameter and is not optimized. In contrast, the BRKGA-MILP approach optimizes ηi, resulting in non-zero values for this variable. Cases where no changes were applied are indicated with a dash. In this context, the symbol “–” denotes relays that are inactive under the analyzed operating mode, meaning they are not assigned any settings or are located on feeders that remain de-energized. This notation explicitly distinguishes inactive relays from those actively participating in the coordination process. The results show that relays R2, R4, R6, R7, R10, and R12 exhibit different ηi values in the hybrid approach, leading to variations in their corresponding settings.

Each OM involves a distinct set of active relays and fault conditions. For this reason, odd-numbered tables present the optimized decision variables (ηi, SCCi, and TMSi) for each relay, while even-numbered tables summarize the resulting tripping times (tif) and total clearing times used to evaluate the objective function. Given that relay activation patterns differ across modes, tabular reporting ensures a consistent and technically accurate comparison framework for both the MILP and BRKGA-MILP approaches.

Table 6 presents the tripping times (tif) obtained for OM1 using both the MILP and BRKGA-MILP approaches. The table lists the primary and back-up relays operating for each fault scenario, along with their respective tripping times. In both approaches, some relays exhibit identical response times, while others show variations due to the optimization of the parameter ηi in the BRKGA-MILP method. Notably, for faults F2, F3, and F4, the hybrid approach results in reduced tripping times compared to the MILP method. Conversely, for fault F5, a slight increase in tripping time is observed. The total tripping time across all fault scenarios is 3.43 s for the MILP approach while it is 2.94 s for the proposed BRKGA-MILP.

6.2. Operative Mode 2

Table 7 presents the relay settings for OM2 using both the MILP and BRKGA-MILP approaches. As in the previous case, the table includes the assigned values for the parameter ηi, which is optimized in the BRKGA-MILP approach. In contrast to the MILP approach, where ηi is set to zero for all relays, the hybrid approach results in varying non-zero values for ηi. These values influence the selection of the SCC and the TMS. For example, relays R1, R2, R4, and R6 show notable changes in ηi, leading to distinct choices for their corresponding SCCi and TMSi settings. The results indicate that multiple relays exhibit changes in ηi, particularly R1, R2, R4, R5, and R11, suggesting a redistribution of coordination parameters. Additionally, modifications in the selected characteristic curves are observed, such as the transition from IEEE-LTI to IEEE-VI for R7 and from IEEE-STI to IEEE-LTEI for R14. These variations reflect the adaptive behavior introduced by the BRKGA-MILP. Note that the EI curve is the dominant curve in the hybrid selection.

Table 8 presents the tripping times obtained for OM2 using both the MILP and BRKGA-MILP approaches. The table lists the operating relays and their corresponding tripping times (tif) for each fault scenario. The results indicate that while both approaches yield similar relay operation sequences, the BRKGA-MILP approach achieves lower total tripping times, reducing the overall system response from 8.54 s to 7.27 s. Notably, in faults F2, F3, and F4, several relays exhibit reduced tripping times under the hybrid approach, particularly R4, R6, R7, and R8. These differences highlight the impact of optimizing ηi through BRKGA, leading to adjustments in relay coordination that enhance the efficiency of fault isolation in the microgrid.

6.3. Operative Mode 3

Table 9 provides the assigned values for ηi, the selected SCCi, and the TMSi for each relay. As observed in previous cases, the MILP approach treats ηi as a fixed parameter set to zero, meaning that only SCCi and TMSi are optimized in this framework. In contrast, the BRKGA-MILP approach integrates a BRKGA to explore a larger solution space by optimizing ηi, allowing for a more flexible coordination of relays.

Furthermore, the modifications in SCCi and TMSi reflect the adaptability of BRKGA-MILP in optimizing fault-clearing times while ensuring selectivity constraints are met. The shifts observed in some relays indicate that BRKGA-MILP prioritizes configurations that reduce tripping times and enhance fault isolation without compromising system reliability. Notably, some relays remain unchanged across both approaches, suggesting that for these cases, the MILP solution was already optimal or had minimal room for improvement.

Table 10 presents the tripping times for OM3 under both the MILP and BRKGA-MILP approaches. The table provides the clearing times (tif) for each relay in response to different fault scenarios. A comparison of the total operating times between the two methodologies reveals that the BRKGA-MILP approach achieves a faster overall response, reducing the total tripping time from 7.21 s (MILP) to 5.89 s (BRKGA-MILP).

The BRKGA-MILP approach dynamically adjusts the coordination of relays, prioritizing configurations that minimize fault-clearing times while maintaining selectivity constraints. The reduction in individual relay tripping times, especially for faults F2, F3, and F4, suggests that BRKGA-MILP successfully identifies relay settings that improve the speed of protection while ensuring proper fault isolation. Additionally, the BRKGA-MILP solution exhibits a more adaptive behavior by modifying the sequences of relay responses for certain faults. For example, relays R5 and R7 exhibit notable differences in tripping times between the two approaches.

6.4. Operative Mode 4

Table 11 presents the relay settings for both the exact model (MILP) and the hybrid model (BRKGA-MILP). Several differences can be observed in the parameter values between the two approaches.

The main difference is the TMSi, since the values in the hybrid model tend to be lower than those in the exact model. For instance, R1’s TMSi decreases from 0.1975 to 0.0291, and R4’s TMSi decreases from 0.4893 to 0.0641. This pattern is also observed in other relays, such as R5, R8, and R11. The selected SCCi also differs between the two models. The exact model predominantly uses conventional curve types such as EI, STI, and VI. The hybrid model, however, includes additional curve types like IEEE-LTVI (R11), IEEE-LTI (R13), and IEEE-LTEI (R14).

Table 12 presents the tripping times for OM4. The total system response time is reduced in the BRKGA-MILP case (7.94 s vs. 8.08 s), indicating an improvement in fault clearance. Most relays exhibit similar tripping times across both methods, with variations typically within 0.01 to 0.05 s. Notable reductions occur in R8 and R11 for F4, where BRKGA-MILP achieves faster operation (0.6489 s vs. 0.7013 s for R8 and 0.9489 s vs. 1.0013 s for R11).

6.5. Performance Analysis of MILP and BRKGA-MILP

Table 13 compares the computational performance of MILP and BRKGA-MILP across different OMs in terms of the objective function (OF) value and computational time in seconds. While BRKGA-MILP consistently achieves lower OF values compared to MILP, indicating improved optimization performance, this improvement comes at a significantly higher computational cost. For instance, in OM1, BRKGA-MILP reduces the OF from 3.43 to 2.93 but increases the computational time from 0.2280 s to 64.81 s. Similarly, in OM3, the OF is reduced from 7.21 to 5.89, yet the computation time rises from 0.3341 s to 178.64 s. This trend highlights the trade-off between solution quality and computational efficiency, with MILP providing faster results while BRKGA-MILP prioritizes optimality at the expense of longer processing times.

6.6. Convergence Analysis of the BRKGA–MILP Matheuristic

Figure 7 displays the convergence trajectories obtained for the four operative modes (OM1–OM4). For each mode, five independent runs are plotted over 60 generations, where the horizontal axis represents the generation index and the vertical axis denotes the best fitness found up to that point. Across all modes, the curves rapidly descend during the initial generations and then stabilize, indicating that the BRKGA–MILP algorithm reliably converges to similar near-optimal solutions with limited evolutionary effort.

It is important to note that the vertical axis in each subplot operates within a very narrow numerical range. For example, in OM1 the best fitness varies only between approximately 2.941 and 2.938. Because each subplot spans such small magnitudes, the visual scale is significantly zoomed in, which amplifies the apparent differences among the five runs. Despite this magnification, the numerical deviations are minimal across all operative modes, and the best solutions achieved by the runs are extremely close to each other. This behavior indicates that the BRKGA–MILP algorithm consistently converges to statistically similar near-optimal solutions, reflecting a robust and stable performance across independent executions. Thus, although the curves may appear visually separated due to axis scaling, the underlying values confirm a high degree of reproducibility and convergence reliability.

Figure 8 illustrates the evolution of the mean fitness of the three independent populations of the BRKGA component across 60 generations for each operative mode (OM1–OM4). The horizontal axis represents the generation index, while the vertical axis corresponds to the average fitness value of each population at every generation. Each color (green, blue, and orange) denotes one of the three independent populations maintained by the BRKGA, and each population is evaluated through five independent runs, represented by the five trajectories of the same color. In addition, the solid black curve depicts the instantaneous average fitness computed across all three populations, providing a global measure of the algorithm’s progress.

At the beginning of the evolutionary process, the mean fitness values of the three populations display noticeable dispersion, which reflects the exploratory behavior intentionally promoted to diversify the search space and avoid premature convergence to local optima. However, as generations progress, the trajectories of the populations progressively move closer together, and the global average (black line) stabilizes. This alignment indicates that the algorithm identifies increasingly similar incumbent solutions across populations. When the three populations converge toward nearly identical mean fitness values, it suggests that they are evolving toward high-quality regions of the solution space, effectively capturing a common near-optimal solution. This collective convergence provides additional evidence of robustness, since the independent evolutionary paths reach consistently high quality solutions, despite starting from different initial conditions.

6.7. Computational Evaluation and Solution Robustness

To assess the computational performance and solution robustness of the proposed BRKGA–MILP matheuristic. A total of 30 independent execution runs were performed for each operative mode (OM1–OM4). The box plots in Figure 9 summarize the distribution of the Objective Function Value obtained across these runs. Regarding the objective function value, the results demonstrate minimal variation across the 30 runs for all operative modes. The median objective function value was nearly identical across all modes, with interquartile ranges (IQRs) that are tightly clustered around the median. This high degree of consistency is a direct consequence of the matheuristic’s structure: the BRKGA component determines the fixed input configuration, and the subsequent MILP subproblem guarantees a global optimal solution for that specific configuration. Consequently, the matheuristic consistently converges to configurations yielding solutions with nearly the same objective value.

Table 14 summarizes the results of 30 independent executions of the BRKGA–MILP model for the four operative modes (OM1–OM4), in 10 generations. For each mode, the table reports the objective value obtained in every run, followed by descriptive statistics (maximum, minimum, mean, median, and standard deviation) that quantify the variability of the stochastic search process.The results show that the BRKGA–MILP model exhibits very high numerical stability, with minimal standard deviations across all operative modes (ranging from 4×104 in OM1 to 1.7×102 in OM3). The closeness between the mean and median confirms the robustness and reproducibility of the algorithm, as all runs converge to nearly identical solutions.

Because the objective function is minimized, the comparison with the MILP benchmark is performed using the worst BRKGA–MILP outcome for each operative mode, i.e., the maximum value among the 30 runs. This provides a conservative evaluation of the hybrid algorithm. Even under this strict criterion, the BRKGA–MILP approach achieves a noticeable reduction in total operating time for OM1–OM3, while OM4 shows a smaller but still positive improvement. These results demonstrate that the BRKGA layer effectively enhances the search process and consistently yields solutions that outperform the standalone MILP formulation, even when compared against the worst-case BRKGA-MILP performance.

7. Conclusions

This study developed and validated an innovative hybrid matheuristic approach for optimal directional overcurrent relay coordination in MGs. The proposed BRKGA-MILP framework demonstrates superior performance compared to conventional methods, achieving up to 18.31% reduction in total relay operating times relative to the MILP reference (with an overall reduction of 11.81%) and up to 14.87% compared with the best previous approaches, while strictly maintaining all coordination time interval constraints across various OMs. The key innovation of the algorithm lies in its simultaneous optimization of protection curve selection, time multiplier settings, and plug setting multipliers, overcoming limitations of traditional approaches that fix curve types as predetermined parameters. This integrated optimization capability proves particularly valuable for microgrids with dynamic operating conditions, where conventional protection schemes often struggle to maintain proper coordination.

The computational experiments conducted on the IEC benchmark microgrid validate the effectiveness of the method in both grid-connected and islanded configurations. While the hybrid approach requires greater computational effort than standalone MILP (typically 3–5 min per OM), the significant improvement in solution quality justifies this investment for protection system design. The successful implementation demonstrates practical applicability for modern power systems with high renewable penetration and variable fault current profiles.

Future research should focus on integrating communication-assisted protection schemes, extending the methodology to meshed microgrid topologies, and developing hardware-in-the-loop validation platforms. Additional opportunities exist in combining the proposed approach with machine learning techniques for real-time setting adaptation. This work provides power engineers with an advanced tool for addressing the growing challenges of microgrid protection in the era of energy transition and decentralized power systems. The framework’s ability to handle complex coordination scenarios while improving protection speed positions it as a valuable contribution to the field, particularly for systems requiring adaptive protection schemes to accommodate renewable energy integration and operational flexibility requirements.

Although the proposed BRKGA–MILP algorithm has been validated through detailed numerical simulations using the IEC benchmark microgrid, it has not yet been tested in a real or hardware-in-the-loop (HIL) platform. Future work will address this limitation by developing a HIL-based implementation of the coordination model, which will enable the evaluation of its robustness under realistic conditions such as parameter uncertainties, measurement noise, and communication delays. This extension will provide a deeper understanding of the algorithm’s practical applicability and support its future integration into experimental microgrid protection systems.

In addition to the improvements observed in relay operating times, the convergence and statistical analyses performed across 30 independent BRKGA–MILP executions provide further evidence of the robustness and repeatability of the proposed approach. The results show that all runs converge to statistically indistinguishable objective values, with very small variances across operative modes. This behavior confirms the algorithm’s ability to consistently identify high-quality solutions despite the stochastic nature of BRKGA. Furthermore, the comparison of the worst-case BRKGA–MILP result with the MILP baseline shows that the matheuristic outperforms the exact model in all operating modes. These findings reinforce the relevance of integrating a metaheuristic layer into the coordination problem, particularly for complex microgrid scenarios where the MILP alone may face combinatorial limitations. The results suggest promising future research directions involving the integration of the MILP model with random-key, local-search-based metaheuristics that do not rely on an evolutionary paradigm [70].

Author Contributions

Conceptualization, L.F.S.-M., S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Data curation, L.F.S.-M. and S.D.S.-Z.; Formal analysis, L.F.S.-M., S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Funding acquisition, J.M.L.-L., N.M.-G. and J.G.V.; Investigation, L.F.S.-M., S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Methodology, L.F.S.-M., S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Project administration, S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Resources, L.F.S.-M., S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Software, L.F.S.-M., and S.D.S.-Z.; Supervision, S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Validation, L.F.S.-M., S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Visualization, L.F.S.-M., S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Writing—original draft, L.F.S.-M., S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V.; Writing—review and editing, L.F.S.-M., S.D.S.-Z., J.M.L.-L., N.M.-G. and J.G.V. All authors have read and agreed to the published version of the manuscript.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The original contributions presented in this study are included in the article. Further inquiries can be directed to the corresponding author.

Acknowledgments

The authors gratefully acknowledge the financial support provided by the Colombian Ministry of Science, Technology, and Innovation “MinCiencias” through “Patrimonio Autónomo Fondo Nacional de Financiamiento para la Ciencia, la Tecnología y la Innovación, Francisco José de Caldas” (Perseo Alliance, Contract No. 112721-392-2023) and the Institución Universitaria Pascual Bravo, Proyecto de Investigación “Regulación de frecuencia eléctrica empleando resortes eléctricos en sistemas con alta penetración de energías renovables”, código PCT00036.

Conflicts of Interest

The authors declare no conflicts of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1 Flowchart summarizing the procedural steps for performing an overcurrent coordination study in accordance with IEEE Std 242-2001.

View Image -

Figure 2 General overview of the protection coordination methodology. (a) One-line diagram of the power system showing the location of the main (Rmain) and backup (Rbackup) relays relative to the fault (Iif). (b) Fault current profile illustrating the reduction of current magnitude as the distance from the source increases. (c) Inverse Time-Current Characteristic (TCC) curves demonstrating selective relay coordination, including the operating times (tif, tjf) and the resulting Coordination Time Interval (CTI).

View Image -

Figure 3 Effect of η variation on the inverse time coordination of DOCRs.

View Image -

Figure 4 Biased random key overview. Source: [56].

View Image -

Figure 5 Overview of the hybridization between the BRKGA and the MILP model for relay coordination.

View Image -

Figure 6 IEC benchmark microgrid topology with feeders, distributed generators, transformers, protection devices, and sectionalizing breakers.

View Image -

Figure 7 Convergence behavior of the BRKGA–MILP matheuristic for OM1–OM4. Each subplot shows the best fitness over 60 generations across five independent runs.

View Image -

Figure 8 Mean-fitness evolution of the three independent BRKGA populations for OM1–OM4 over 60 generations. Each color represents one population across five runs.

View Image -

Figure 9 Distribution of the Objective Function Value for the BRKGA–MILP matheuristic across 30 independent runs in the four operative modes (OM1–OM4).

View Image -

Short-circuit current levels (kA) for each fault and operative mode.

DL Fault OM1 OM2 OM3 OM4
DL-1 F4 5.130 6.989 6.428 3.188
DL-2 F3 8.375 11.891 10.116 3.558
DL-3 F5 3.695 5.904 4.223 3.590
DL-4 F2 5.130 8.725 6.363 4.158
DL-5 F1 3.695 6.296 4.293 3.941

Definition of operative modes based on generation availability.

OMs GD1 GD2 GD3 GD4 Utility Description
OM1 disabled disabled disabled disabled enabled Only grid
OM2 enabled enabled enabled enabled enabled GD+GC
OM3 enabled enabled disabled disabled enabled 50% GD+GC
OM4 enabled enabled enabled enabled disabled Islanded

Comparative analysis of total tripping time.

OM [50] [49] [54] [51] [69] MILP [53] BRKGA-MILP
OM1 7.53 6.64 4.99 3.86 3.09 3.43 2.94
OM2 19.18 17.48 13.66 8.58 8.77 8.54 7.27
OM3 14.04 12.67 10.71 8.39 6.91 7.21 5.89
OM4 15.56 15.56 12.63 8.83 9.08 8.08 7.94

Comparative performance of BRKGA–MILP with MILP baseline and best previous method.

OM MILP BestPrev BRKGA–MILP MILP (%) BestPrev (%)
OM1 3.43 3.09 2.94 14.29 4.85
OM2 8.54 8.54 7.27 14.87 14.87
OM3 7.21 6.91 5.89 18.31 14.76
OM4 8.08 8.08 7.94 1.73 1.73
Total 27.26 26.62 24.04 11.81 9.69

Relay settings for OM1.

MILP BRKGA-MILP
Relay η i SCC i TMS i Relay η i SCC i TMS i
R1 - - - R1 - - -
R2 0 EI 0.0001 R2 0.2960 STI 0.0001
R3 - - - R3 - - -
R4 0 VI 0.3884 R4 0.9751 EI 0.1433
R5 - - - R5 - - -
R6 0 LTI 0.1053 R6 0.5168 EI 0.8967
R7 0 SI 0.0000 R7 0.1254 IEEE-VI 0.0001
R8 - - - R8 - - -
R9 - - - R9 - - -
R10 0 EI 0.0001 R10 0.2165 STI 0.0001
R11 - - - R11 - - -
R12 0 EI 0.0001 R12 0.2983 EI 0.0001
R13 - - - R13 - - -
R14 - - - R14 - - -
R15 - - - R15 - - -

Tripping times for OM1.

MILP BRKGA-MILP
Fault Relay tif [s] Fault Relay tif [s]
F1 R2 0.0000 F1 R2 0.0000
F1 R4 0.3000 F1 R4 0.3000
F2 R4 0.2127 F2 R4 0.1537
F2 R6 0.5127 F2 R6 0.4537
F3 R6 0.3092 F3 R6 0.1696
F3 R7 0.6092 F3 R7 0.4911
F4 R7 0.7594 F4 R7 0.4911
F4 R12 0.0000 F4 R12 0.0000
F5 R6 0.7232 F5 R6 0.8797
F5 R10 0.0000 F5 R10 0.0000
Total time 3.43 Total time 2.94

Relay settings for OM2.

MILP BRKGA-MILP
Relay η i SCC i TMS i Relay η i SCC i TMS i
R1 0 EI 0.1975 R1 0.8644 EI 0.0233
R2 0 EI 0.0001 R2 0.7615 STI 0.0001
R3 0 STI 0.0001 R3 0.0737 STI 0.0001
R4 0 LTI 0.0556 R4 0.9281 EI 0.2446
R5 0 EI 0.1904 R5 0.7297 EI 0.0283
R6 0 LTI 0.1068 R6 0.5942 EI 0.8103
R7 0 IEEE-LTI 0.4723 R7 0.4047 IEEE-VI 0.0001
R8 0 EI 0.1480 R8 0.6536 EI 0.0167
R9 0 STI 0.0001 R9 0.4306 STI 0.0001
R10 0 EI 0.0001 R10 0.2165 EI 0.0001
R11 0 IEEE-LTVI 0.0429 R11 0.8937 EI 0.0132
R12 0 EI 0.0001 R12 0.3510 EI 0.0001
R13 0 IEEE-STI 0.0781 R13 0.9870 IEEE-VI 0.0076
R14 0 IEEE-STI 0.0294 R14 0.3645 IEEE-LTEI 0.0039
R15 0 IEEE-VI 0.0001 R15 0.9330 IEEE-VI 0.0001

Tripping times for OM2.

MILP BRKGA-MILP
Fault Relay tif [s] Fault Relay tif [s]
F1 R1 0.2362 F1 R1 0.2294
F1 R2 0.0000 F1 R2 0.0001
F1 R4 0.3000 F1 R4 0.3001
F1 R13 0.5362 F1 R13 0.5294
F2 R1 0.3001 F2 R1 0.3001
F2 R3 0.0001 F2 R3 0.0001
F2 R4 0.1890 F2 R4 0.1219
F2 R6 0.4890 F2 R6 0.4219
F2 R15 0.4911 F2 R15 0.4924
F3 R5 0.0883 F3 R5 0.0819
F3 R6 0.2831 F3 R6 0.1453
F3 R7 0.5831 F3 R7 0.4911
F3 R8 0.5831 F3 R8 0.4453
F3 R15 0.4912 F3 R15 0.4946
F4 R5 0.3000 F4 R5 0.3000
F4 R7 0.7722 F4 R7 0.4913
F4 R8 0.5026 F4 R8 0.3699
F4 R11 0.8026 F4 R11 0.6699
F4 R12 0.0000 F4 R12 0.0000
F5 R6 0.7972 F5 R6 1.0819
F5 R9 0.0001 F5 R9 0.0001
F5 R10 0.0000 F5 R10 0.0000
F5 R14 0.3001 F5 R14 0.3001
F5 R15 0.4913 F5 R15 0.0000
Total time 8.54 Total time 7.27

Relay settings for OM3.

MILP BRKGA-MILP
Relay η i SCC i TMS i Relay η i SCC i TMS i
R1 - - - R1 - - -
R2 0 EI 0.0001 R2 0.5558 STI 0.0001
R3 - - - R3 - - -
R4 0 VI 0.4548 R4 0.7232 SI 0.0000
R5 0 EI 0.0267 R5 0.9473 STI 0.0001
R6 0 VI 0.9695 R6 0.6642 EI 0.7385
R7 0 IEEE-LTI 0.4793 R7 0.9826 IEEE-VI 0.0001
R8 0 EI 0.1495 R8 0.6889 EI 0.0156
R9 - - - R9 - - -
R10 0 EI 0.0001 R10 0.6760 STI 0.0001
R11 0 IEEE-EI 0.3291 R11 0.2178 STI 0.5764
R12 0 EI 0.0001 R12 0.5205 EI 0.0001
R13 - - - R13 - - -
R14 - - - R14 - - -
R15 0 IEEE-VI 0.0069 R15 0.0620 IEEE-LTI 0.1994

Tripping times for OM3.

IEC and IEEE Curves Only IEEE Curves
Fault Relay tif [s] Fault Relay tif [s]
F1 R2 0.0000 F1 R2 0.0001
F1 R4 0.3000 F1 R4 0.3001
F2 R4 0.1992 F2 R4 0.1356
F2 R6 0.4992 F2 R6 0.4356
F2 R15 0.4992 F2 R15 0.4356
F3 R5 0.1222 F3 R5 0.0003
F3 R6 0.2890 F3 R6 0.1499
F3 R7 0.5890 F3 R7 0.4912
F3 R8 0.5890 F3 R8 0.4499
F3 R15 0.5005 F3 R15 0.4495
F4 R5 0.3000 F4 R5 0.0000
F4 R7 0.7478 F4 R7 0.4916
F4 R8 0.5077 F4 R8 0.3724
F4 R11 0.8077 F4 R11 0.6724
F4 R12 0.0000 F4 R12 0.0000
F5 R6 0.7537 F5 R6 0.9651
F5 R10 0.0000 F5 R10 0.0001
F5 R15 0.5098 F5 R15 0.5388
Total time 7.21 Total time 5.89

Relay settings for OM4.

MILP BRKGA-MILP
Relay η i SCC i TMS i Relay η i SCC i TMS i
R1 0 EI 0.1975 R1 0.7385 EI 0.0291
R2 0 STI 0.0001 R2 0.9246 STI 0.0001
R3 0 STI 0.0001 R3 0.6674 STI 0.0001
R4 0 EI 0.4893 R4 0.8478 EI 0.0641
R5 0 EI 0.4488 R5 0.8198 EI 0.0612
R6 0 IEEE-VI 0.0235 R6 0.4524 IEEE-VI 0.0049
R7 - - - R7 - - -
R8 0 EI 0.2065 R8 0.9362 EI 0.0160
R9 0 STI 0.0001 R9 0.4783 STI 0.0001
R10 0 STI 0.0001 R10 0.6623 STI 0.0001
R11 0 SI 0.1940 R11 0.9890 IEEE-LTVI 0.0107
R12 0 STI 0.0001 R12 0.3047 STI 0.0001
R13 0 VI 0.1462 R13 0.9805 IEEE-LTI 0.1631
R14 0 STI 0.3300 R14 0.3004 IEEE-LTEI 0.0045
R15 0 IEEE-VI 0.0219 R15 0.0629 IEEE-VI 0.0154

Tripping times for OM4.

MILP BRKGA-MILP
Fault Relay tif [s] Fault Relay tif [s]
F1 R1 0.2362 F1 R1 0.2309
F1 R2 0.0000 F1 R2 0.0001
F1 R4 0.3000 F1 R4 0.3001
F1 R13 0.5362 F1 R13 0.5309
F2 R1 0.3001 F2 R1 0.3001
F2 R3 0.0001 F2 R3 0.0001
F2 R4 0.2171 F2 R4 0.2141
F2 R6 0.5171 F2 R6 0.5141
F2 R15 0.5171 F2 R15 0.5141
F3 R5 0.2081 F3 R5 0.2048
F3 R6 0.5137 F3 R6 0.5107
F3 R8 0.8137 F3 R8 0.8107
F3 R15 0.5331 F3 R15 0.5286
F4 R5 0.3000 F4 R5 0.3001
F4 R8 0.7013 F4 R8 0.6489
F4 R11 1.0013 F4 R11 0.9489
F4 R12 0.0000 F4 R12 0.0001
F5 R6 0.5370 F5 R6 0.5380
F5 R9 0.0001 F5 R9 0.0001
F5 R10 0.0000 F5 R10 0.0001
F5 R14 0.3001 F5 R14 0.3001
F5 R15 0.5505 F5 R15 0.5448
Total time 8.08 Total time 7.94

Comparison of Objective Function (OF) Values and Computational Time for MILP and BRKGA-MILP Across Different OMs.

MILP BRKGA-MILP
OM OF [s] Time [s] OF [s] Time [s]
OM1 3.43 0.2280 2.93 64.81
OM2 8.54 0.4689 7.27 195.88
OM3 7.21 0.3341 5.89 178.64
OM4 8.08 0.2114 7.94 149.35

Statistical summary of 30 independent BRKGA–MILP runs for each operative mode (OM1–OM4).

Metric OM1 OM2 OM3 OM4
2.9395 7.2366 5.8602 7.9231
2.9400 7.2446 5.8533 7.9202
2.9389 7.2520 5.8772 7.9282
2.9391 7.2390 5.8882 7.9346
2.9396 7.2448 5.8704 7.9368
2.9393 7.2481 5.8915 7.9358
2.9390 7.2510 5.8450 7.9305
2.9391 7.2335 5.8527 7.9093
2.9392 7.2591 5.8667 7.9253
2.9390 7.2424 5.8757 7.9331
2.9396 7.2279 5.8883 7.9226
2.9387 7.2505 5.8510 7.9078
2.9398 7.2346 5.8654 7.9402
2.9388 7.2366 5.8584 7.9218
2.9388 7.2444 5.8680 7.9173
2.9387 7.2447 5.8718 7.9200
2.9397 7.2406 5.9158 7.9272
2.9393 7.2407 5.8766 7.9159
2.9389 7.2348 5.8984 7.9272
2.9387 7.2338 5.8622 7.9123
2.9391 7.2402 5.8526 7.9239
2.9392 7.2370 5.8613 7.9180
2.9388 7.2403 5.8653 7.9172
2.9389 7.2372 5.8527 7.9206
2.9393 7.2592 5.8697 7.9339
2.9387 7.2411 5.8668 7.9218
2.9387 7.2334 5.8550 7.9202
2.9389 7.2351 5.8515 7.9243
2.9399 7.2550 5.9014 7.9223
2.9393 7.2427 5.8623 7.9180
Max 2.9400 7.2592 5.9158 7.9402
Min 2.9387 7.2279 5.8450 7.9078
Mean 2.9391 7.2420 5.8692 7.9236
Median 2.9391 7.2406 5.8660 7.9224
Std. Dev. 0.0004 0.0078 0.0170 0.0080
MILP 3.43 8.54 7.21 8.08
Difference 0.4900 1.2808 1.2942 0.1398

References

1. Muhire, F.; Turyareeba, D.; Adaramola, M.S.; Nantongo, M.; Atukunda, R.; Olyanga, A.M. Drivers of green energy transition: A review. Green Energy Resour.; 2024; 2, 100105. [DOI: https://dx.doi.org/10.1016/j.gerr.2024.100105]

2. Wang, Y.; Yuan, Z.; Luo, H.; Zeng, H.; Huang, J.; Li, Y. Promoting low-carbon energy transition through green finance: New evidence from a demand-supply perspective. Energy Policy; 2024; 195, 114376. [DOI: https://dx.doi.org/10.1016/j.enpol.2024.114376]

3. Li, J.; Khan, I. Transition to a sustainable energy balance and conversion factors in economic management: Basics for transitioning to a low-carbon economy. Energy; 2024; 313, 133724. [DOI: https://dx.doi.org/10.1016/j.energy.2024.133724]

4. Igliński, B.; Kiełkowska, U.; Mazurek, K.; Drużyński, S.; Pietrzak, M.B.; Kumar, G.; Veeramuthu, A.; Skrzatek, M.; Zinecker, M.; Piechota, G. Renewable energy transition in Europe in the context of renewable energy transition processes in the world. A review. Heliyon; 2024; 10, e40997. [DOI: https://dx.doi.org/10.1016/j.heliyon.2024.e40997]

5. Saldarriaga-Zuluaga, S.D.; Lopez-Lezama, J.M.; Muñoz-Galeano, N. Protection Coordination in Microgrids: Current Weaknesses, Available Solutions and Future Challenges. IEEE Lat. Am. Trans.; 2020; 18, pp. 1715-1723. [DOI: https://dx.doi.org/10.1109/TLA.2020.9387642]

6. Welankiwar, A.; Sharma, R.; Kumar, B. Challenges and mitigation techniques in adaptive protection for microgrids: Comprehensive review. Renew. Energy Focus; 2024; 51, 100652. [DOI: https://dx.doi.org/10.1016/j.ref.2024.100652]

7. Jaramillo Serna, J.d.J.; López-Lezama, J.M. Alternative Methodology to Calculate the Directional Characteristic Settings of Directional Overcurrent Relays in Transmission and Distribution Networks. Energies; 2019; 12, 3779. [DOI: https://dx.doi.org/10.3390/en12193779]

8. Foqha, T.; Alsadi, S.; Omari, O.; Refaat, S.S. Optimization Techniques for Directional Overcurrent Relay Coordination: A Comprehensive Review. IEEE Access; 2024; 12, pp. 1952-2006. [DOI: https://dx.doi.org/10.1109/ACCESS.2023.3347393]

9. Ghanbari, M.; Gandomkar, M.; Nikoukar, J. Protection Coordination of Bidirectional Overcurrent Relays Using Developed Particle Swarm Optimization Approach Considering Distribution Generation Penetration and Fault Current Limiter Placement. IEEE Can. J. Electr. Comput. Eng.; 2021; 44, pp. 143-155. [DOI: https://dx.doi.org/10.1109/ICJECE.2020.3018876]

10. Rezaei, N.; Uddin, M.N.; Amin, I.K.; Othman, M.L.; Marsadek, M. Genetic Algorithm-Based Optimization of Overcurrent Relay Coordination for Improved Protection of DFIG Operated Wind Farms. IEEE Trans. Ind. Appl.; 2019; 55, pp. 5727-5736. [DOI: https://dx.doi.org/10.1109/TIA.2019.2939244]

11. Labrador Rivas, A.E.; Gallego Pareja, L.A. Coordination of directional overcurrent relays that uses an ant colony optimization algorithm for mixed-variable optimization problems. Proceedings of the 2017 IEEE International Conference on Environment and Electrical Engineering and 2017 IEEE Industrial and Commercial Power Systems Europe (EEEIC/I&CPS Europe); Milan, Italy, 6–9 June 2017; pp. 1-6. [DOI: https://dx.doi.org/10.1109/EEEIC.2017.7977750]

12. Tharakan, K.I.; Swathika, O.V.G. Optimum coordination of using overcurrent relay using firefly and ant colony optimization algorithm. Proceedings of the 2017 International Conference on Computing Methodologies and Communication (ICCMC); Erode, India, 18–19 July 2017; pp. 617-621. [DOI: https://dx.doi.org/10.1109/ICCMC.2017.8282541]

13. Rashtchi, V.; Gholinezhad, J.; Farhang, P. Optimal coordination of overcurrent relays using Honey Bee algorithm. Proceedings of the International Congress on Ultra Modern Telecommunications and Control Systems; Moscow, Russia, 18–20 October 2010; pp. 401-405. [DOI: https://dx.doi.org/10.1109/ICUMT.2010.5676606]

14. Sampaio, F.C.; Belmino, L.M.; Leão, R.P.S.; Sampaio, R.F.; Melo, L.S.; Barroso, G.C. Optimal Directional Overcurrent Relays Coordination using Directional Bat Algorithm. Proceedings of the 2019 IEEE PES Innovative Smart Grid Technologies Conference - Latin America (ISGT Latin America); Gramado, Brazil, 15–18 September 2019; pp. 1-6. [DOI: https://dx.doi.org/10.1109/ISGT-LA.2019.8894935]

15. Srivastava, A.; Tripathi, J.M.; Krishan, R.; Parida, S.K. Optimal Coordination of Overcurrent Relays Using Gravitational Search Algorithm With DG Penetration. IEEE Trans. Ind. Appl.; 2018; 54, pp. 1155-1165. [DOI: https://dx.doi.org/10.1109/TIA.2017.2773018]

16. Moirangthem, J.; Krishnanand, K.R.; Saranjit, N. Optimal coordination of overcurrent relay using an enhanced discrete differential evolution algorithm in a distribution system with DG. Proceedings of the 2011 International Conference on Energy, Automation and Signal; Bhubaneswar, India, 28–30 December 2011; pp. 1-6. [DOI: https://dx.doi.org/10.1109/ICEAS.2011.6147114]

17. Barzegari, M.; Bathaee, S.; Alizadeh, M. Optimal coordination of directional overcurrent relays using harmony search algorithm. Proceedings of the 2010 9th International Conference on Environment and Electrical Engineering; Prague, Czech Republic, 16–19 May 2010; pp. 321-324. [DOI: https://dx.doi.org/10.1109/EEEIC.2010.5489935]

18. Bedekar, P.P.; Bhide, S.R. Optimum Coordination of Directional Overcurrent Relays Using the Hybrid GA-NLP Approach. IEEE Trans. Power Deliv.; 2011; 26, pp. 109-119. [DOI: https://dx.doi.org/10.1109/TPWRD.2010.2080289]

19. Kida, A.A.; Gallego, L.A. A high-performance hybrid algorithm to solve the optimal coordination of overcurrent relays in radial distribution networks considering several curve shapes. Electr. Power Syst. Res.; 2016; 140, pp. 464-472. [DOI: https://dx.doi.org/10.1016/j.epsr.2016.05.029]

20. Sadeghi, S.; Hashemi-Dezaki, H.; Entekhabi-Nooshabadi, A.M. Optimal Protection Scheme of Micro-Grids Considering N-1 Contingency By A New Hybrid GA-PSO-LP Optimization Algorithm. Proceedings of the 2021 11th Smart Grid Conference (SGC); Tabriz, Iran, 7–9 December 2021; pp. 1-5. [DOI: https://dx.doi.org/10.1109/SGC54087.2021.9664076]

21. Srinivas, S.T.P.; Swarup, K.S. A hybrid GA — Interval linear programming implementation for microgrid relay coordination considering different fault locations. Proceedings of the 2017 7th International Conference on Power Systems (ICPS); Pune, India, 21–23 December 2017; pp. 785-790. [DOI: https://dx.doi.org/10.1109/ICPES.2017.8387396]

22. Liu, A.; Yang, M.T. A New Hybrid Nelder-Mead Particle Swarm Optimization for Coordination Optimization of Directional Overcurrent Relays. Math. Probl. Eng.; 2012; 2012, 456047. [DOI: https://dx.doi.org/10.1155/2012/456047]

23. Liu, A.; Yang, M.T. Optimal Coordination of Directional Overcurrent Relays Using NM-PSO Technique. Proceedings of the 2012 International Symposium on Computer, Consumer and Control; Taichung, Taiwan, 4–6 June 2012; pp. 678-681. [DOI: https://dx.doi.org/10.1109/IS3C.2012.176]

24. Mousavi Motlagh, S.H.; Mazlumi, K. Optimal Overcurrent Relay Coordination Using Optimized Objective Function. ISRN Power Eng.; 2014; 2014, pp. 1-10. [DOI: https://dx.doi.org/10.1155/2014/869617]

25. Khurshaid, T.; Wadood, A.; Farkoush, S.G.; Kim, C.H.; Cho, N.; Rhee, S.B. Modified Particle Swarm Optimizer as Optimization of Time Dial Settings for Coordination of Directional Overcurrent Relay. J. Electr. Eng. Technol.; 2019; 14, pp. 55-68. [DOI: https://dx.doi.org/10.1007/s42835-018-00039-z]

26. Wadood, A.; Kim, C.H.; Khurshiad, T.; Farkoush, S.; Rhee, S.B. Application of a Continuous Particle Swarm Optimization (CPSO) for the Optimal Coordination of Overcurrent Relays Considering a Penalty Method. Energies; 2018; 11, 869. [DOI: https://dx.doi.org/10.3390/en11040869]

27. Irfan, M.; Wadood, A.; Khurshaid, T.; Khan, B.M.; Kim, K.C.; Oh, S.R.; Rhee, S.B. An Optimized Adaptive Protection Scheme for Numerical and Directional Overcurrent Relay Coordination Using Harris Hawk Optimization. Energies; 2021; 14, 5603. [DOI: https://dx.doi.org/10.3390/en14185603]

28. Habib, K.; Lai, X.; Wadood, A.; Khan, S.; Wang, Y.; Xu, S. Hybridization of PSO for the Optimal Coordination of Directional Overcurrent Protection Relays. Electronics; 2022; 11, 180. [DOI: https://dx.doi.org/10.3390/electronics11020180]

29. Kida, A.A.; Labrador Rivas, A.E.; Gallego, L.A. An improved simulated annealing–linear programming hybrid algorithm applied to the optimal coordination of directional overcurrent relays. Electr. Power Syst. Res.; 2020; 181, 106197. [DOI: https://dx.doi.org/10.1016/j.epsr.2020.106197]

30. Khurshaid, T.; Wadood, A.; Gholami Farkoush, S.; Yu, J.; Kim, C.H.; Rhee, S.B. An Improved Optimal Solution for the Directional Overcurrent Relays Coordination Using Hybridized Whale Optimization Algorithm in Complex Power Systems. IEEE Access; 2019; 7, pp. 90418-90435. [DOI: https://dx.doi.org/10.1109/ACCESS.2019.2925822]

31. ElSayed, S.K.; Elattar, E.E. Hybrid Harris hawks optimization with sequential quadratic programming for optimal coordination of directional overcurrent relays incorporating distributed generation. Alex. Eng. J.; 2021; 60, pp. 2421-2433. [DOI: https://dx.doi.org/10.1016/j.aej.2020.12.028]

32. Korashy, A.; Kamel, S.; Jurado, F.; Youssef, A.R. Hybrid Whale Optimization Algorithm and Grey Wolf Optimizer Algorithm for Optimal Coordination of Direction Overcurrent Relays. Electr. Power Components Syst.; 2019; 47, pp. 644-658. [DOI: https://dx.doi.org/10.1080/15325008.2019.1602687]

33. Srinivas, S.T.P.; Shanti Swarup, K. Application of improved invasive weed optimization technique for optimally setting directional overcurrent relays in power systems. Appl. Soft Comput.; 2019; 79, pp. 1-13. [DOI: https://dx.doi.org/10.1016/j.asoc.2019.03.045]

34. Rizk-Allah, R.M.; El-Fergany, A.A. Effective coordination settings for directional overcurrent relay using hybrid Gradient-based optimizer. Appl. Soft Comput.; 2021; 112, 107748. [DOI: https://dx.doi.org/10.1016/j.asoc.2021.107748]

35. Usama, M.; Moghavvemi, M.; Mokhlis, H.; Mansor, N.N.; Farooq, H.; Pourdaryaei, A. Optimal Protection Coordination Scheme for Radial Distribution Network Considering ON/OFF-Grid. IEEE Access; 2021; 9, pp. 34921-34937. [DOI: https://dx.doi.org/10.1109/ACCESS.2020.3048940]

36. Noghabi, A.S.; Sadeh, J.; Mashhadi, H.R. Considering Different Network Topologies in Optimal Overcurrent Relay Coordination Using a Hybrid GA. IEEE Trans. Power Deliv.; 2009; 24, pp. 1857-1863. [DOI: https://dx.doi.org/10.1109/TPWRD.2009.2029057]

37. Bashir, M.; Taghizadeh, M.; Sadeh, J.; Mashhadi, H.R. A new hybrid particle swarm optimization for optimal coordination of over current relay. Proceedings of the 2010 International Conference on Power System Technology; Hangzhou, China, 24–28 October 2010; pp. 1-6. [DOI: https://dx.doi.org/10.1109/POWERCON.2010.5666600]

38. Damchi, Y.; Mashhadi, H.R.; Sadeh, J.; Bashir, M. Optimal coordination of directional overcurrent relays in a microgrid system using a hybrid particle swarm optimization. Proceedings of the 2011 International Conference on Advanced Power System Automation and Protection; Beijing, China, 16–20 October 2011; pp. 1135-1138. [DOI: https://dx.doi.org/10.1109/APAP.2011.6180976]

39. Khurshaid, T.; Wadood, A.; Gholami Farkoush, S.; Kim, C.H.; Yu, J.; Rhee, S.B. Improved Firefly Algorithm for the Optimal Coordination of Directional Overcurrent Relays. IEEE Access; 2019; 7, pp. 78503-78514. [DOI: https://dx.doi.org/10.1109/ACCESS.2019.2922426]

40. Ramli, S.P.; Mokhlis, H.; Wong, W.R.; Muhammad, M.A.; Mansor, N.N. Optimal coordination of directional overcurrent relay based on combination of Firefly Algorithm and Linear Programming. Ain Shams Eng. J.; 2022; 13, 101777. [DOI: https://dx.doi.org/10.1016/j.asej.2022.101777]

41. Misaghi, M.; Yaghoobi, M. Improved invasive weed optimization algorithm (IWO) based on chaos theory for optimal design of PID controller. J. Comput. Des. Eng.; 2019; 6, pp. 284-295. [DOI: https://dx.doi.org/10.1016/j.jcde.2019.01.001]

42. Kumar, P.; Singh Rana, A. Review of optimization techniques for relay coordination in consideration with adaptive schemes of Microgrid. Electr. Power Syst. Res.; 2024; 230, 110240. [DOI: https://dx.doi.org/10.1016/j.epsr.2024.110240]

43. Londe, M.A.; Pessoa, L.S.; Andrade, C.E.; Resende, M.G. Biased random-key genetic algorithms: A review. Eur. J. Oper. Res.; 2024; 321, [DOI: https://dx.doi.org/10.1016/j.ejor.2024.03.030]

44. Maniezzo, V.; Stützle, T.; Voß, S. Matheuristics; Springer: Berlin/Heidelberg, Germany, 2021.

45. IEEE. Recommended Practice for Protection and Coordination of Industrial and Commercial Power Systems (IEEE Buff Book). IEEE Std 242-2001 (Revision of IEEE Std 242-1986) [IEEE Buff Book]; IEEE: Piscataway, NJ, USA, 2001; pp. 1-710. [DOI: https://dx.doi.org/10.1109/IEEESTD.2001.93369]

46. Waqar, A.; Hussain, B.; Ahmad, S.; Yahya, T.; Sarwar, M. A Communication-less Protection Strategy to Ensure Protection Coordination of Distribution Networks with Embedded DG. arXiv; 2019; [DOI: https://dx.doi.org/10.48550/arXiv.1906.07571] arXiv: 1906.07571

47. Mohammadreza Alizadeh, M.R.H. The Power System and Microgrid Protection—A Review. Appl. Sci.; 2020; 10, 8271. [DOI: https://dx.doi.org/10.3390/app10228271]

48. Fang, X.; Chen, J.; Li, F.; Tang, Y.; Yang, Y. A review of protection strategies and coordination techniques for microgrids. Electr. Power Syst. Res.; 2022; 208, 107925. [DOI: https://dx.doi.org/10.1016/j.epsr.2022.107925]

49. El-Naily, N.; Saad, S.M.; Hussein, T.; Mohamed, F.A. A novel constraint and non-standard characteristics for optimal over-current relays coordination to enhance microgrid protection scheme. IET Gener. Transm. Distrib.; 2019; 13, pp. 780-793. [DOI: https://dx.doi.org/10.1049/iet-gtd.2018.5021]

50. Saad, S.M.; El-Naily, N.; Mohamed, F.A. A new constraint considering maximum PSM of industrial over-current relays to enhance the performance of the optimization techniques for microgrid protection schemes. Sustain. Cities Soc.; 2019; 44, pp. 445-457. [DOI: https://dx.doi.org/10.1016/j.scs.2018.09.030]

51. Saldarriaga-Zuluaga, S.D.; López-Lezama, J.M.; Muñoz-Galeano, N. Adaptive protection coordination scheme in microgrids using directional over-current relays with non-standard characteristics. Heliyon; 2021; 7, e06665. [DOI: https://dx.doi.org/10.1016/j.heliyon.2021.e06665]

52.IEEE C37.112-1996 IEEE Standard Inverse-Time Characteristic Equations for Overcurrent Relays; IEEE Standard Institute of Electrical and Electronics Engineers: New York, NY, USA, 1996; Available online: https://ieeexplore.ieee.org/document/9767803 (accessed on 25 November 2025).

53. Serna-Montoya, L.F.; Saldarriaga-Zuluaga, S.D.; López-Lezama, J.M.; Muñoz-Galeano, N. Optimal Microgrid Protection Coordination for Directional Overcurrent Relays Through Mixed-Integer Linear Optimization. Energies; 2025; 18, 2035. [DOI: https://dx.doi.org/10.3390/en18082035]

54. Saldarriaga-Zuluaga, S.D.; López-Lezama, J.M.; Muñoz-Galeano, N. Optimal Coordination of Overcurrent Relays in Microgrids Considering a Non-Standard Characteristic. Energies; 2020; 13, 922. [DOI: https://dx.doi.org/10.3390/en13040922]

55. Bean, J.C. Genetic algorithms and random keys for sequencing and optimization. ORSA J. Comput.; 1994; 6, pp. 154-160. [DOI: https://dx.doi.org/10.1287/ijoc.6.2.154]

56. Gonçalves, J.F.; Resende, M.G. Biased random-key genetic algorithms for combinatorial optimization. J. Heuristics; 2011; 17, pp. 487-525. [DOI: https://dx.doi.org/10.1007/s10732-010-9143-1]

57. Andrade, C.E.; Toso, R.F.; Gonçalves, J.F.; Resende, M.G. The multi-parent biased random-key genetic algorithm with implicit path-relinking and its real-world applications. Eur. J. Oper. Res.; 2021; 289, pp. 17-30. [DOI: https://dx.doi.org/10.1016/j.ejor.2019.11.037]

58. Cavalheiro, E.M.; Vergilio, A.H.; Lyra, C. Optimal configuration of power distribution networks with variable renewable energy resources. Comput. Oper. Res.; 2018; 96, pp. 272-280. [DOI: https://dx.doi.org/10.1016/j.cor.2017.09.021]

59. De Faria, H.; Resende, M.G.; Ernst, D. A biased random key genetic algorithm applied to the electric distribution network reconfiguration problem. J. Heuristics; 2017; 23, pp. 533-550. [DOI: https://dx.doi.org/10.1007/s10732-017-9355-8]

60. Raposo, A.A.; Rodrigues, A.B.; da Silva, M.d.G. Robust meter placement for state estimation considering distribution network reconfiguration for annual energy loss reduction. Electr. Power Syst. Res.; 2020; 182, 106233. [DOI: https://dx.doi.org/10.1016/j.epsr.2020.106233]

61. Mikulski, S.; Tomczewski, A. Use of energy storage to reduce transmission losses in meshed power distribution networks. Energies; 2021; 14, 7304. [DOI: https://dx.doi.org/10.3390/en14217304]

62. Jourdan, L.; Basseur, M.; Talbi, E.G. Hybridizing exact methods and metaheuristics: A taxonomy. Eur. J. Oper. Res.; 2009; 199, pp. 620-629. [DOI: https://dx.doi.org/10.1016/j.ejor.2007.07.035]

63. Talbi, E.G. Combining metaheuristics with mathematical programming, constraint programming and machine learning. Ann. Oper. Res.; 2016; 240, pp. 171-215. [DOI: https://dx.doi.org/10.1007/s10479-015-2034-y]

64. Bynum, M.L.; Hackebeil, G.A.; Hart, W.E.; Laird, C.D.; Nicholson, B.L.; Siirola, J.D.; Watson, J.P.; Woodruff, D.L. Pyomo-Optimization Modeling in Python; Springer: Berlin/Heidelberg, Germany, 2021; Volume 67.

65. Huangfu, Q.; Hall, J.J. Parallelizing the dual revised simplex method. Math. Program. Comput.; 2018; 10, pp. 119-142. [DOI: https://dx.doi.org/10.1007/s12532-017-0130-5]

66. Ustun, T.S.; Ozansoy, C.; Zayegh, A. Modeling of a Centralized Microgrid Protection System and Distributed Energy Resources According to IEC 61850-7-420. IEEE Trans. Power Syst.; 2012; 27, pp. 1560-1567. [DOI: https://dx.doi.org/10.1109/TPWRS.2012.2185072]

67. Kar, S.; Samantaray, S.R.; Zadeh, M.D. Data-Mining Model Based Intelligent Differential Microgrid Protection Scheme. IEEE Syst. J.; 2017; 11, pp. 1161-1169. [DOI: https://dx.doi.org/10.1109/JSYST.2014.2380432]

68.IEC 60909 Short-Circuit Currents in Three-Phase a.c. Systems; IEC Standard International Electrotechnical Commission: Geneva, Switzerland, 2001.

69. Plongkrathok, C.; Kaiyawong, K.; Chayakulkheeree, K. Adaptive Directional Overcurrent Relay Coordination in Microgrid Considering Optimal Full Setting Parameters. Electr. Power Components Syst.; 2024; 52, pp. 985-1000. [DOI: https://dx.doi.org/10.1080/15325008.2023.2237973]

70. Chaves, A.A.; Resende, M.G.; Schuetz, M.J.; Brubaker, J.K.; Katzgraber, H.G.; de Arruda, E.F.; Silva, R.M. A Random-Key optimizer for combinatorial optimization: AA Chaves et al. J. Heuristics; 2025; 31, 32. [DOI: https://dx.doi.org/10.1007/s10732-025-09568-z]

© 2025 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/). Notwithstanding the ProQuest Terms and Conditions, you may use this content in accordance with the terms of the License.