A novel Fibonacci-inspired enhancement of the Bees Algorithm: application to robotic disassembly sequence planning

Abstract

The Bees Algorithm (BA) is an intelligent, nature-inspired metaheuristic algorithm first introduced in 2005 and based on the foraging activity of honeybees. Enhancements to the BA have been made continually since its introduction, with some aimed at reducing the number of user-determined parameters. The goal is to achieve optimal results without extensive efforts required for parameter tuning. This article presents an enhanced version of the BA. The enhancement uses the Fibonacci sequence, which can be found in the family tree of the male bee (drone), to give the numbers of bees assigned to conduct local searches. The proposed BA is simple, having only four parameters that must be set by the user, as opposed to five in the basic combinatorial version. To evaluate the performance of the new algorithm in robotic disassembly sequence planning, a case study is conducted on gear pump disassembly, a commonly used benchmark problem in robotic disassembly. A comparison is made between the proposed Fibonacci Bees Algorithm and the Basic Bees Algorithm (BBA). The results are assessed using a statistical performance metric and a performance evaluation index. The findings demonstrate that the proposed algorithm is better than the BBA for more complex problems while exhibiting a similar performance on simpler tasks. Future research will explore broader applications for both continuous and combinatorial problems, demonstrating the versatility of the algorithm across different domains.

Keywords:

• Bees Algorithm
• disassembly
• Fibonacci
• optimisation
• remanufacturing

REVIEWING EDITOR:

• Zude Zhou, Wuhan University of Technology, China

SUBJECTS:

• Operation Research
• Optimization Techniques
• Computational Techniques
• Design of Experiments General

1. Introduction

The Bees Algorithm (BA) has demonstrated success in solving complex real-world problems (Castellani & Pham, 2023; De la Torre Gutiérrez & Pham, 2018). The BA, developed in 2005 (Pham et al., 2005, 2006), draws inspiration from the foraging behaviour of honeybees. The BA employs a population of artificial ‘bees’ to navigate and exploit the search space, thereby emulating the foraging behaviour exhibited by bees in their quest for nectar sources. Since its inception, the BA has garnered considerable attention from researchers and has been applied in numerous domains (Castellani & Pham, 2023).

However, a problem with most versions of the BA is that they have many parameters (more than five) which the user has to set, making the algorithm difficult for the uninitiated to apply. This article presents a new version of the BA with just four user-selected parameters. In the proposed algorithm, the numbers of bees assigned to conduct local searches are determined with the help of the Fibonacci sequence which can be observed in the family tree of a male bee (drone). This natural inspiration has not only expanded the landscape of optimisation techniques but also furnished researchers and practitioners with a versatile tool adaptable to a wide range of complex problems across many domains. This study introduces a novel approach that significantly contributes to the parameter configuration and enhancement of the BA.

The article is structured as follows: Section 2 provides an overview of the BA and its variants, highlighting applications in the context of robotic disassembly planning. Section 3 focusses on the proposed enhancement of the BA, drawing inspiration from the Fibonacci sequence to yield the Fibonacci Bees Algorithm or BA_F. A review of other metaheuristics inspired by Fibonacci is provided to highlight the novelty of BA_F. In Section 4, the chosen case study is presented along with a detailed justification for its selection. The evaluation of the performance of BA_F is elaborated upon in Section 5. Section 6 presents the experiments conducted using both the BA and the enhanced version in RDSP. Section 7 discusses the results obtained. Finally, the study concludes by summarising the key findings and insights obtained from the research and suggesting areas for further investigation.

2. Bees Algorithm

The BA, originally developed for continuous problems, has six user-determined parameters (Pham et al., 2005) and is commonly referred to as the Basic Bees Algorithm (BBA). The pseudocode of the BBA (Pham et al., 2006) is given in Algorithm 0.

Algorithm 0:

Basic Bees Algorithm Pseudocode

Input: n: number of scout bees, m: number of selected sites, e: number of elite sites, nsp: recruited bees for other selected sites, nep: recruited bees for elite sites, ngh: initial size of neighbourhood

1 Function BBA(n, m, e, nsp, nep, ngh):

2 InitialisePopulation with random solutions

3 EvaluateFitness of the population

4 while stopping criterion not met do

5  Forming new population

6  SelectSites for neighbourhood search

7  RecruitBees for selected sites (more bees for the best e sites) and evaluate fitness

8  SelectFittestBees from each patch

9  AssignRandomBees to search randomly and evaluate their fitness

10 end

11 return BestBee

The neighbourhood (ngh) for combinatorial problems depends on the design of the local search, which may involve operators like swap, insert, reverse, mutation, 2-OPT and 3-OPT. In the combinatorial version of BBA (Ismail et al., 2020; Liu et al., 2018a, b; Zeybek & Koç, 2015; Zeybek et al., 2021), where the swap, insert, and reverse operators are used, the neighbourhood size is considered equal to the sequence length. Therefore, for the combinatorial version, only five parameters need to be set: n, m, e, nsp and nep.

The introduction of site abandonment and neighbourhood shrinking to BBA, resulting in two additional parameters that need to be configured (Castellani & Pham, 2023; Pham & Castellani, 2009, 2015), represents a noteworthy enhancement and gives rise to a variant of the algorithm known as the Standard Bees Algorithm (SBA). SBA requires 7–8 parameters to be set by the user. Alternative methods for recruitment, neighbourhood modification and site abandonment have been extensively explored in the literature (Castellani & Pham, 2023). For a comprehensive understanding of the variants of the BA developed until 2017, interested readers are encouraged to refer to the survey paper by Hussein, Sahran and Sheikh Abdullah, which provides an in-depth analysis of the various modifications and enhancements in the algorithm (Hussein et al., 2017).

The utilisation of metaheuristics has been the subject of critique owing to their dependence on parameter values. This reliance can result in a lengthy procedure, as the optimal settings for each problem must be determined while considering their unique attributes. Previous researchers opted for the application of Taguchi and Design of Experiments (DoE) methods to discern the optimal parameter setting. However, these techniques are time-consuming (Hartono, Ismail et al., 2023). Several efforts have been made to reduce the parameter tuning in the BA, such as the application of fuzzy selection for self-tuning (Pham & Darwish, 2008). The Ternary Bees Algorithm, introduced by Laili et al. (2019), uses only three bees and incorporates a single parameter setting for the site abandonment threshold. Another reduced-parameter version of BA is BA₂, which has two parameters. BA₂ is inspired by the traplining foraging behaviour of honeybees (Ismail, 2021; Ismail & Pham, 2023). The achievement of an adaptive parameter-free version and the potential for further parameter removal remain uncertain in the pursuit of enhanced versions of the BA (Castellani & Pham, 2023). However, it is apparent that reducing the number of parameters significantly decreases the effort required for parameter tuning.

BA₂, initially developed for continuous problems and later extended to address combinatorial problems by incorporating specialised local operators tailored for the vehicle routing problem (VRP), has been previously investigated by Ismail (2021). However, the primary focus of this study is on application to robotic disassembly problems. Robotic Disassembly Sequence Planning (RDSP) entails planning the sequence and actions required to disassemble objects, considering factors such as interdependencies, sizes, shapes, connections, accessibility, tools and multiple dimensions. These inherent complexities differentiate robotic disassembly from the VRP, which typically deals with fewer dimensions. Therefore, direct comparison with both the continuous version of BA₂ or its combinatorial version of VRP falls outside the scope of this research. While earlier research by Ismail (2021) drew inspiration from the traplining behaviour of bees to reduce the parameter count of the BA, this study takes a different approach by drawing inspiration from the Fibonacci sequence-based family tree pattern of honeybees. Building upon this unique source of inspiration, the proposed enhancement to the algorithm introduces a new configuration with four parameters. By exploring this alternative parameter setting, the research aims to investigate potential benefits and performance improvements. This novel approach contributes to the field by offering a fresh perspective on parameter configuration in the BA for robotic disassembly problems.

Numerous methodologies and algorithms have been developed to address the challenges posed by the problem of optimally planning disassembly (Gao et al., 2020). A limited number of previous studies have utilised exact methods to address either simple or prototypical problems (Ming et al., 2019; Yin et al., 2022). While mathematical and exact methods can provide optimal solutions for small-scale instances, their applicability to larger problems is limited by the nondeterministic polynomial (NP)-hard nature of the problem (Cevikcan et al., 2020; Gao et al., 2020; Guo et al., 2022; McGovern & Gupta, 2011). The number of possible subassemblies (N_n), the number of complete disassembly sequences (c_n) and the total number of disassembly sequences (P_n) in a disassembly problem can be determined theoretically using Eqs. (1)–(3) (Lambert, 1997). For instance, with four parts, there are 15 subassemblies, 15 complete disassembly sequences and a total of 41 disassembly sequences. In the case of 10 parts, the number of subassemblies increases to 1023, the complete disassembly sequences amount to 34,459,425, and the total number of disassembly sequences is 314,726,297. These theoretical calculations demonstrate the exponential growth in the number of disassembly sequences as the number of products increases. Such insights shed light on the combinatorial nature of disassembly problems and underscore the challenges associated with exploring all possible disassembly sequences. (1) $$N_n=\left(2.N_{n-1}\right)+1$$(1) (2) $$c_n=\left(2n-3\right).c_{n-1}$$(2) (3) $$\begin{array}{l}{P}_1=1\\ P_2=\frac{1}{2}\mathrm{ }\left(\begin{array}{c}1\\ 2\end{array}\right).P_1+1\\ \begin{array}{l}{P}_3=\left(\begin{array}{c}1\\ 3\end{array}\right).P_2+1\\ ⋮\\ P_{10}=\left(\begin{array}{c}1\\ 10\end{array}\right).P_9+\left(\begin{array}{c}2\\ 10\end{array}\right).P_2.P_8+\left(\begin{array}{c}3\\ 10\end{array}\right).P_3.P_7\\ +\left(\begin{array}{c}4\\ 10\end{array}\right).P_4.P_6+\frac{1}{2}.\left(\begin{array}{c}5\\ 10\end{array}\right).P_5.P_5+1\end{array}\end{array}$$(3)

Another theoretical calculation suggests that in the context of planning the optimal sequence for a product with n parts, the exploration of the search tree involves examining n! nodes, and the collision checking operation requires k x n! computations, where k represents the average number of directions tested for component removal (Costa et al., 2018). These theoretical calculations, coupled with the findings from previous research, provide compelling evidence to support the assertion that the disassembly problem is inherently complex and exhibits exponential growth in complexity. Metaheuristic algorithms have garnered significant attention from researchers due to their effectiveness in navigating the expansive search space of disassembly problems and achieving near-optimal solutions (Cevikcan et al., 2020; Gao et al., 2020; Ma et al., 2011; McGovern & Gupta, 2011).

The prevalence of metaheuristic approaches is evident in robotic disassembly research, where many researchers have adopted them (Gao et al., 2020; Hartono, 2023). Notably, metaheuristic algorithms offer the advantages of relatively short computation times, meeting practical requirements (Gao et al., 2020) and facilitating efficient decision-making in real-world scenarios (Laili et al., 2022b). In a real-world context, the efficient discovery of near-optimal solutions holds greater significance than achieving exact solutions. Moreover, real products exhibit complexity that further amplifies the challenge. During the initial stages of research, many studies employed toy problems to explore the mathematical and conceptual aspects, which is understandable considering the inherent complexity of the disassembly problem.

The BA has demonstrated its robustness and effectiveness in a wide range of remanufacturing applications (Caterino et al., 2022, 2023; Kerin et al., 2023; Zeybek, 2023), establishing it as a reliable optimisation approach within this domain. In particular, in the field of robotic disassembly, the BA is recognised as one of the most popular metaheuristic algorithms (Chen, Xu, et al., 2020; Cui et al., 2023; Hartono et al., 2022a, 2022b, 2023; Laili et al., 2019, 2022a, 2022b; Liu, Zhou, Pham, Xu, Ji, et al., 2018a; Liu, Zhou, Pham, Xu, Yan, et al., 2018b; Liu et al., 2020; Liu, Liu, et al., 2023; Liu, Xu, et al., 2023; Yang et al., 2022; Ye et al., 2022). Consistently outperforming other algorithms, the BA is highly regarded for its exceptional performance and capabilities. Thus, it remains a compelling choice for addressing optimisation problems in robotic disassembly. However, so far, only one robotic disassembly study (Laili et al., 2019) has been conducted on parameter reduction.

To ensure a meaningful and relevant evaluation, this research appropriately employs a comparison with Enhanced Discrete Bees Algorithm (EDBA), which is the most commonly adopted BA for RDSP. This choice enables an assessment of the effectiveness of the proposed enhancements in addressing the specific challenges posed by robotic disassembly problems. By aligning with the objectives and scope of this study, this approach facilitates a comprehensive examination of the proposed enhancements and their impact on addressing the complexities associated with robotic disassembly.

As previously mentioned, the BA was initially developed for the continuous domain. However, when the algorithm is applied to combinatorial problems such as the Travelling Salesman Problem (TSP), the local search component needs to be modified accordingly. In the case of the TSP, the BBA incorporates local search operators such as swap, insert and reverse (Ismail et al., 2020), which are specifically designed to manipulate the sequence of cities. Similarly, in the context of robotic disassembly problems, the BBA’s local search is adapted to include swap, insert and mutation operators (Liu, Zhou, Pham, Xu, Ji, et al., 2018a; Liu, Zhou, Pham, Xu, Yan, et al., 2018b) that are tailored to address the specific requirements and constraints of these combinatorial problems. In this research, the BA and the proposed BA utilise the swap, insert and mutation operators in their local search, as described by Liu, Zhou, Pham, Xu, Ji, et al. (2018).

3. Fibonacci Bees Algorithm (BA_F)

Leonardo Fibonacci, a mathematician from the 13th century, is widely recognised for his significant contributions to the field of mathematics. One of his most notable achievements was the introduction of the Fibonacci sequence (Scott & Marketos, 2014). This sequence is a series of numbers, where each number is the sum of the two preceding numbers. It begins with 0 and 1, and the subsequent numbers are generated by adding the previous two numbers. The sequence progresses as follows: 0, 1, 1, 2, 3, 5, 8, 13, 21, 34 and subsequent numbers. The pattern extends indefinitely, revealing a sequence that possesses unique mathematical properties and occurs in diverse natural phenomena. The pattern can be found in a variety of natural phenomena, including leaf arrangements, spirals in seashells, tree branches and flower petals (Omotehinwa & Ramon, 2013; Scott & Marketos, 2014). The bee lineage parallels the Fibonacci sequence, establishing a link between bee ancestry and the numerical pattern seen in the Fibonacci sequence. The family tree of a male bee or drone exhibits a fascinating relationship with the Fibonacci sequence, as the number of ancestors in each successive generation corresponds to the Fibonacci numbers (Omotehinwa & Ramon, 2013; Scott & Marketos, 2014). Drones originate from unfertilised eggs through a reproductive process known as parthenogenesis. Parthenogenesis involves the development of embryos without fertilisation, resulting in male offspring that possess a mother but do not have a father. In contrast, female bees, which include both worker bees and queen bees, originate from fertilised eggs. Figure 1 illustrates the family tree of a drone in a simple ancestry model (Koshy, 2017; Scott & Marketos, 2014). Interestingly, upon counting the total number of bees in each generation, it is revealed that they conform to the Fibonacci sequence. This correlation emphasises the prevalence of Fibonacci patterns in honeybee population and the drone’s mathematical structure in its genealogical lineage.

Figure 1. Fibonacci sequence in the family tree of a drone (adapted from Koshy, 2017).

To investigate whether previous enhancements of the BA have incorporated the inspiration from the Fibonacci sequence, given its relevance to bees, the following steps were conducted. An extensive search was conducted using the Scopus database, employing the keywords “bee” and “Fibonacci”. This search yielded a total of nine articles; however, none of them were found to be directly relevant to the BA. Additionally, a search was performed using the Scite database, but no pertinent results were found regarding the utilisation of the Fibonacci sequence in the BA. These searches across multiple databases indicate that, as of the present, there is no scholarly evidence to suggest that prior versions of the BA have integrated the Fibonacci sequence as a means of enhancement.

Further exploration was conducted using the Scopus database, employing the keywords “Fibonacci” and “heuristic”, with the aim of identifying algorithms that draw inspiration from Fibonacci. The initial search yielded a total of 62 articles. After refining the search based on titles and language criteria (English only) and availability of content, 43 articles were retained for further analysis. Upon examining the contents, nine articles were deemed irrelevant and were excluded from the analysis. Additionally, 18 articles that were identified as improvements to the original algorithm were also removed. Of the remaining articles, 15 focussed on metaheuristic algorithms, many of which were improved or hybridised with other algorithms, such as the Genetic Algorithm (GA) (Gudla & Ganguli, 2005; Sha & Pan, 2018; Yeh, 2006), Salp Swarm Algorithm (SSA) (Tran-Ngoc et al., 2023), Wolf Pack Algorithm (WPA) (Zhao et al., 2020), Tabu Search (Aras & Aksen, 2008), Particle Swarm Optimisation (PSO) (Khosla & Verma, 2022) and Grey Wolf Optimiser (GWO) (Khosla & Verma, 2022). These articles covered a range of inspirations related to Fibonacci, including the Fibonacci search, Fibonacci heaps, Fibonacci trees, the golden ratio and the application of Fibonacci indicators in the stock market.

Table 1 summarises the Fibonacci inspirations of these algorithms and their applications. The Automated Hybrid Genetic Algorithm (AHGA) is utilised for benchmark function optimisation, wherein a combination of a GA and a local optima avoidance mechanism is employed (Gudla & Ganguli, 2005). Fibonacci numbers are used to select the number of GA banks in AHGA. The Memetic Algorithm (MA) is utilised for solving the multi-stage supply chain network problem, with the option of Fibonacci number generation to select the best offspring regardless of its origin (Yeh, 2006). The Fibonacci Tree Optimisation (FTO) algorithm is designed specifically for power grid optimisation, utilising a data storage structure known as a Fibonacci tree (Da et al., 2018). The Fibonacci indicator algorithm (FIA) finds its application in benchmark function optimisation, incorporating Fibonacci retracement and time zone techniques commonly utilised in stock market trading (Etminaniesfahani et al., 2018). For wellbore trajectory optimisation in the oil and gas industry, the Fibonacci sequence-based quantum GA (FSQGA) is employed, reflecting the rotation angle step in the quantum GA (Sha & Pan, 2018). The Golden Ratio Optimisation Method (GROM) is applied to benchmark function optimisation, with the movement direction of the algorithm following the formulation of the golden ratio (Nematollahi et al., 2020). Similarly, the Improved Wolf Pack Algorithm (IWPA) is employed in benchmark function optimisation, where the step length of scout wolves is adjusted based on the Fibonacci sequence (Zhao et al., 2020). The Child Drawing Development Optimisation (CDDO) algorithm uses the golden ratio to calculate two solution points in benchmark function optimisation (Abdulhameed & Rashid, 2022). The Opposition-based Learning PSO GWO (Opp-PSOGWO) algorithm is applied to benchmark function optimisation, generating an opposite population for searching using the Fibonacci sequence (Khosla & Verma, 2022). In the area of bridge structural health monitoring, a Fibonacci sequence-based optimisation algorithm (FSOA) rearranges the population of the algorithm using the golden ratio (Tran-Ngoc et al., 2023). Additionally, the table includes algorithms employed for solving the collection centre location problem inspired by the Fibonacci search (Aras & Aksen, 2008), project scheduling problem using the Multi-Objective Fibonacci-Based Algorithm (MOFA) inspired by the Fibonacci sequence (Hosseinian & Baradaran, 2020), and VRP using Brand Local Search Algorithm (BLSA), Adaptive Variable Neighbourhood Search (AVNS) and Evolutionary Local Search (ELS) that inspired by Fibonacci heaps (Wei et al., 2014; Zachariadis & Kiranoudis, 2010; Zhang et al., 2015).

Table 1. Metaheuristics inspired by Fibonacci.

Author(s)	Name	Fibonacci	Application	Hybrid with
Gudla and Ganguli (2005)	AHGA	Sequence	Benchmark function	GA
Yeh (2006)	MA	Sequence	Multistage supply chain network	GA, greedy
Aras and Aksen (2008)	Not named	Search	Collection centre location problem	Tabu search
Zachariadis and Kiranoudis (2010)	BLSA	Heaps	VRP
Wei et al. (2014)	AVNS	Heaps	VRP
Zhang et al. (2015)	ELS	Heaps	VRP
Da et al. (2018)	FTO	Tree	Power grid
Etminaniesfahani et al. (2018)	FIA	Retracement and time zone	Benchmark function
Sha and Pan (2018)	FSQGA	Sequence	Oil and gas (wellbore trajectory)	GA
Nematollahi et al. (2020)	GROM	Golden ratio	Benchmark function
Zhao et al. (2020)	IWPA	Sequence	Benchmark function	WPA
Hosseinian and Baradaran (2020)	MOFA	Sequence	Project scheduling problem
Abdulhameed and Rashid (2022)	CDDO	Golden ratio	Benchmark function
Khosla and Verma (2022)	Opp-PSOGWO	Sequence	Benchmark function	PSO, GWO
Tran-Ngoc et al. (2023)	FSOA	Golden ratio	Bridge structural health monitoring	SSA

The findings presented in Table 1 provide additional support for the distinctiveness of the proposed Fibonacci Bees Algorithm (BA_F) compared to previous research. While previous algorithms utilised Fibonacci numbers in various ways, such as selecting the number of banks, reflecting rotation angles, adjusting step lengths and generating populations, BA_F introduces a novel framework specifically inspired by the Fibonacci sequence. In the process of improving BA, various ideas were experimented with and tested to identify improvements. The chosen and most effective idea is outlined below.

BA_F is based on the observed Fibonacci sequence-based family tree pattern in drones, where the number of drone ancestors follows the Fibonacci sequence, as mentioned earlier. However, BA_F introduces a different approach. Instead of employing the Fibonacci sequence to count ancestors, BA_F uses it to determine the number of bees sent to the flower patches. This decision is motivated by the inherent growth pattern represented by the Fibonacci sequence, which aligns with the objective of maximising the foraging efficiency of the algorithm. By employing a ranking-based approach that recruits bees according to the Fibonacci sequence for targeting the most promising patches, BA_F aims to exploit the patches more effectively, potentially yielding improved results. This departure from the conventional interpretation of the Fibonacci sequence within the context of drone ancestry showcases an innovative adaptation of the concept to address the specific requirements of the proposed algorithm.

In addition, unlike the original BA, which employs “elite” and “other selected” sites, BA_F focusses solely on the selected sites, eliminating the need to distinguish between “elite” and “other selected” categories in the selected sites. Within these selected sites (sites with the best fitness values), bees conduct their search, and the allocation of bees to exploit the selected sites (nr) is determined by a ranking system that follows the Fibonacci sequence. The top-ranked bee receives the maximum number from the specified Fibonacci sequence while lower-ranking bees are assigned decreasing numbers from the sequence. For example, with three patches, m = 3 and maxnr = 8. The 1st patch attracts 8 recruited bees (nr₁ = maxnr = 8), the 2nd patch attracts 5 recruited bees (nr₂ = 5) and the 3rd patch attracts 3 recruited bees (nr₃ = 3). Although it is not certain if the allocation of foragers in nature follows the Fibonacci sequence, this differential allocation does reflect the fact that more foragers are allocated to higher-quality flower patches (Shackleton et al., 2023; Visscher & Seeley, 1982). After a specified number of tries, if the fittest bee in a patch remains unchanged, nr for that patch becomes 0 and a fresh set of patches is initialised.

Both BBA and EDBA share the same parameters and structure, with the only difference being the incorporation of a check for disassembly sequence feasibility after random initialisation. By disregarding that minor point, comparing EDBA and BA_F becomes equivalent to comparing BBA and BA_F. This article will compare the proposed BA_F and EDBA in the context of RDSP. The problem is the minimisation of disassembly time. The comparison will utilise the Statistical Performance Metric (SPM) and the Performance Evaluation Index (PEI) presented in Hartono (2023). The objective function (disassembly time) and number of function of evaluations (NFE) will be recorded for comparison purposes. In this work, the input data is a disassembly information matrix based on the modified interference matrix proposed by Liu, Zhou, Pham, Xu, Ji, et al. (2018a). This matrix ensures that the disassembly process takes into consideration the presence of fasteners while being feasible by adhering to precedence constraints. Using the matrix thus produces feasible disassembly sequences. In Algorithm 1 (EDBA) and Algorithm 2 (BA_F), the input is the robotic disassembly information matrix, denoted as dis_m.

Algorithm 1:

The pseudo-code of EDBA for RDSP

Input: n: number of scout bees, m: number of selected sites, e: number of elite sites, nep: recruited bees for elite sites, nsp: recruited bees for other selected sites, dis_m: robotic disassembly information matrix

Output: RDSP(sequence, direction, tool, f:cost)

1 Function EDBA(n, m, e, nsp, nep):

2  Start

3  initialRDSP ← GlobalFDS(dis_m: sequence, direction)  // Generate initial population of feasible disassembly sequences (FDS)

4  while stopping criterion not met do

5    Evaluate population fitness

6    f ← FVALUE(initialRDSP)

7    Sort population according to f

8    Select m sites for local search

     // Generate local sites with waggle dance

9   for EliteSite(1 to e) do

     // Assign best elite site bee

10    BestEliteSiteBee ← the scout bee that found the elite site

11    for RecruitedEliteSiteBee(1 to nep) do

      // Do feasibility check

12     while feasibility not met do

13       RecruitedEliteSiteBee ←WaggleDance(dis_m : sequence, direction)

      //Mutate the disassembly direction

14     RecruitedEliteSiteBee ←Mutation(dis_m : direction)

15     Evaluate fitness of RecruitedEliteSiteBee

16     if RecruitedEliteSiteBee is better than BestEliteSiteBee then

        //Update BestEliteSiteBee

17       BestLocalBee = RecruitedEliteSiteBee

18   for OtherSelectedSite(1 to (m-e)) do

     // Assign best other selected site bee

19    BestOtherSelectedSiteBee ← the scout bee that found the other selected site

20    for RecruitedOtherSelectedSiteBee (1 to nsp) do

      // Do feasibility check

21     while feasiblity not met do

22        RecruitedOtherSelectedSiteBee ← WaggleDance(dis_m : sequence, direction)

       // Mutate the disassembly direction

23     RecruitedOtherSelectedSiteBee ← Mutation(dis_m : direction)

24     Evaluate fitness of RecruitedOtherSelectedSiteBee

25     if RecruitedOtherSelectedSiteBee is better than

        BestOtherSelectedSiteBee then

         // Update BestOtherSelectedSiteBee

26        BestLocalBee = RecruitedOtherSelectedSiteBee

    // Assign remaining scout bees for global search

27   for RemainingScoutBee(1 to (n-m)) do

28     RemainingScoutBee ← GlobalFDS(dis_m : sequence, direction)

29   Evaluate fitness of the new population

30   Sort population according to f

      // Store the best RDSP with minimum cost f

31   Best RDSP = Bee with minimum cost f

32  return Best RDSP (Bee with minimum cost f)

Algorithm 2:

The pseudo-code of BA_F for RDSP

Input: n: number of scout bees, m: number of selected sites, nr: number of bees recruited for selected sites using ranking based recruitment, max_rv: maximum number of re-visits before the nr is set to zero, dis_m: robotic disassembly information matrix

Output: RDSP(sequence, direction, tool, f:cost)

1 Function BA_F(n, m, nr, max_rv):

2  Start

3  initialRDSP ← GlobalFDS(dis_m: sequence, direction)  // Generate initial population of feasible disassembly sequences (FDS)

4  revisit = 0

5  while stopping criterion not met do

6    Evaluate population fitness

7    f ← FVALUE(initialRDSP)

8    Sort population according to f

9    Select m sites for local search

    // Generate local sites with waggle dance

10   for SelectedSite(1 to m) do

      // Assign best local bee

11     BestLocalBee ← the scout bee that found the site

      // Allocation of bees to the selected sites

12     nr ← number from Fibonacci sequence

13     for RecruitedBee (1 to nr) do

        // Do feasibility check

14       while feasiblity not met do

15         RecruitedBee ← WaggleDance(dis_m : sequence, direction)

        // Mutate the disassembly direction

16       RecruitedBee ← Mutation(dis_m : direction)

17       Evaluate fitness of RecruitedBee

18       if RecruitedBee is better than BestLocalBee then

          // Update BestLocalBee

19          BestLocalBee = RecruitedBee

20       else

21         revisit ← revisit + 1

22         if revisit ≥ max_rv then

            // Abandon site and randomly generated a new site

23           NewSite ← FDS(dis_m : sequence, direction)

    // Assign remaining scout bees for global search

24   for RemainingScoutBee(1 to (n-m)) do

25      RemainingScoutBee ← GlobalFDS(dis_m : sequence, direction)

26   Evaluate fitness of the new population

27   Sort population according to f

    // Store the best RDSP with minimum cost f

28   Best RDSP = Bee with minimum cost f

29   return Best RDSP (Bee with minimum cost f)

4. Case study description

The case study chosen for this research is focussed on industrial gear pumps, which have been frequently used in previous studies on robotic disassembly of EoL products (Belhadj et al., 2022; Chen, Xu, et al., 2020; Chen, Yao, et al., 2020; Frizziero et al., 2022; Hartono et al., 2022a, 2022b, 2023; Laili et al., 2022c; Liu et al., 2018a, 2018b, 2020; Liu, Liu, et al., 2023; Liu, Xu, et al., 2023; Liu & Wang, 2017; Ramírez et al., 2020). Furthermore, the broad industrial application and low wear of gear pumps make them an ideal candidate for demonstrating the optimisation results and studying complete disassembly without destruction, highlighting their suitability as a test subject in this research. The results of this study on the robotic disassembly of gear pumps contribute to the broader context of robotic disassembly research, as they demonstrate the potential of using optimisation algorithms to disassemble EoL products efficiently and sustainably.

The two gear pumps in this research have different flow rates. Gear Pump A has a flow rate of 7.5 l/min, while Gear Pump B has a flow rate of 10 l/min. The data utilised in this study were derived from academic sources, and the 3D models were acquired from Grabcad (Grabcad Community, 2020; Liu, Zhou, Pham, Xu, Ji, et al., 2018a; Ramírez et al., 2020). The assumption for this study, taken from Liu, Zhou, Pham, Xu, Ji, et al. (2018a), is as follows: (1) the disassembly product is the same with geometric information and components that are known and applicable for repetitive disassembly; (2) All parts can be disassembled completely. The properties of the gear pumps can be found in Table 2.

Table 2. Properties of Gear Pumps A and B.

Item	Name	Disassembly point			Disassembly tool	tb(xi)
		X	Y	Z
Gear Pump A
1	Bolt A	49.4	105.5	−12.6	Spanner-I	3
2	Bolt B	74.4	81.0	−12.6	Spanner-I	3
3	Bolt C	74.4	45.0	−12.6	Spanner-I	3
4	Bolt D	49.4	20.5	−12.6	Spanner-I	3
5	Bolt E	24.4	45.0	−12.6	Spanner-I	3
6	Bolt F	24.4	81.0	−12.6	Spanner-I	3
7	Cover	49.4	63.0	−20.6	Gripper-II	4
8	Gasket	49.4	105.5	1.4	Gripper-I	3
9	Gear A	49.4	81.0	3.4	Gripper-I	6
10	Gear B	49.4	45.0	3.4	Gripper-I	6
11	Driven Shaft A	49.4	81.0	−7.6	Gripper-I	4
12	Base	49.4	81.0	49.4	Gripper-II	8
13	Driven Shaft B	49.4	45.0	152.4	Gripper-I	4
14	Packing Gland	49.4	45.0	91.4	Gripper-I	2
15	Gland Nut	49.4	45.0	96.4	Spanner-II	3
Gear Pump B
1	Bolt A	59.1	114.0	−48.4	Spanner-I	4
2	Bolt B	90.3	89.0	−48.4	Spanner-I	4
3	Bolt C	90.3	33.0	−48.4	Spanner-I	4
4	Bolt D	59.1	8.0	−48.4	Spanner-I	4
5	Bolt E	27.9	33.0	−48.4	Spanner-I	4
6	Bolt F	27.9	89.0	−48.4	Spanner-I	4
7	Cover	59.1	82.0	−64.6	Gripper-II	5
8	Gasket	59.1	114.0	−31.4	Gripper-I	4
9	Gear A	59.1	82.0	−30.9	Gripper-I	6
10	Gear B	59.1	40.0	−30.9	Gripper-I	6
11	Shaft A	59.1	40.0	−48.9	Gripper-I	4
12	Base	59.1	114.0	7.1	Gripper-II	4
13	Shaft B	59.1	82.0	136.1	Gripper-I	8
14	Gland A	59.1	94.8	34.1	Gripper-I	3
15	Gland B	59.1	94.8	41.1	Gripper-I	3
16	Gland C	59.1	94.8	48.1	Gripper-I	3
17	Gland D	59.1	94.8	55.1	Gripper-I	3
18	Gland E	59.1	82.0	79.1	Gripper-I	3
19	Bolt stud A	35.1	82.0	89.1	Spanner-II	3
20	Bolt stud B	83.1	82.0	89.1	Spanner-II	3
21	Nut A	35.1	82.0	84.1	Spanner-III	4
22	Nut B	83.1	82.0	84.1	Spanner-III	4
24	Nut C	35.1	82.0	87.1	Spanner-III	4
25	Nut D	83.1	82.0	87.1	Spanner-III	4

5. Performance evaluation

The literature review reveals that prior research has not fully utilised the potential of statistics, as they were only applied to the final outcomes. In addition, contradictory results have been observed across performance metrics (Laili et al., 2021). To bridge these identified gaps, two approaches are used in this research. One approach involves the application of statistical tests, while the other utilises a straightforward yet versatile metric. These proposed methods aim to address shortcomings and provide valuable insights into the performance analysis of the algorithm. The specifics of the first approach are elaborated upon in Section 5.1, while the details of the second approach are provided in Section 5.2.

5.1. Statistical performance metric

The SPM, introduced by Hartono et al., employs a statistical test not only for evaluating algorithm performance but also for identifying optimal parameter settings (Hartono et al., 2022a). This method follows a well-defined selection process for the appropriate statistical test. The advantage of this approach lies in its ability to statistically analyse observed differences and determine the optimal parameter settings for the chosen algorithm. The methodology is outlined as follows: Initially, the results are visualised using descriptive statistics. Subsequently, the assumption checklist is performed. If the number of experiments conducted exceeds 30, it is necessary to test whether the data adhere to certain assumptions, such as normality and homogeneity. In the case of any assumption violations, a nonparametric test is carried out. Conversely, if the number of experiments is below 30, a nonparametric test is employed. Finally, if the results indicated statistical significance, a post hoc test was conducted. The decision process can be observed in Algorithm 2. The choice of software for conducting the statistical analysis is dependent on the researcher’s preference and may involve the use of either commercial or freely available statistical software. In this research, statistical tests were performed using IBM SPSS Statistics 29 and MATLAB 2020b.

By employing this method, researchers can rigorously analyse the performance of an algorithm and identify the optimal parameter settings in a statistically sound manner. This systematic approach adds credibility to the experimental results and contributes to the advancement of algorithm optimisation techniques.

Algorithm 3:

Statistical Performance Metric

Require: data

1  Function Perform Descriptive Statistics And Box plot (data):

2  Calculate  descriptive statistics for the data;

3  Generate a boxplot to visually represent the results;

4  Function PerformAssumptionChecklist(data):

5   Check  the data size;

6   if the data size is below the threshold (eg 30), then

7    Perform  a nonparametric test;

8    Choose an appropriate nonparametric test;

9    return;

10   end

11   Test  for normality;

12   Test  for homogeneity of variances;

13  Function DecisionStepOne(data):

14   assumptionFail  ← false;

15   PerformAssumptionChecklist(data);

16   if any assumption fails (data size, normality, or homogeneity) then

17    Perform a nonparametric test based on the failed assumption;

18    assumptionFail  ← true;

19   end

20   if all assumptions are met and assumptionFail is false then

21    Perform  a parametric test;

22   end

23  Function DecisionStepTwo (result):

24   if result of the test is statistically significant then

25    Conduct post hoc tests to determine specific group differences;

26    return;

27   end

28     return;

29  Main Algorithm;

30  PerformDescriptiveStatisticsAndBoxplot(data);

31  DecisionStepOne(data);

32  DecisionStepTwo (result);

5.2. Performance Evaluation Index

As previously mentioned, relying solely on diagrams or figures to illustrate these indicators is insufficient, as it neglects the possibility of conflicting indicators, an aspect that has not received sufficient attention in scholarly investigations. As previously mentioned, there is a possibility of conflicting results obtained with different indicators (Laili et al., 2021). This research uses the PEI, a new index for evaluating the performance of algorithms introduced by Hartono et al. (2023). The PEI can be used for MO as well as SO optimisation. The concept behind the proposed PEI is derived from formulas commonly used in the Multiple-Criteria Decision-Making (MCDM) literature. The index provides decision makers with a valuable tool for expediently evaluating multiple criteria and facilitating prompt decision-making. To calculate the index, two approaches are discussed in the MCDM literature. The first approach involves assigning weights to different criteria and summing the scores to obtain an overall score. However, this approach has certain limitations, including the need for additional computations such as normalisation and the inadequate accommodation of conflicting objectives related to specific criteria maximisation or minimisation. To address these limitations, Tofallis (2014) introduced the second approach, known as the multiplicative approach. This approach resolves the limitations of the first approach by considering the multiplicative combination of scores, resulting in a simplified computation process and improved accommodation of conflicting objectives. In contrast to the additive method, the multiplicative approach does not require the normalisation or re-scaling of criteria, as the final outcome is unaffected by these operations (Tofallis, 2014). The multiplicative approach is followed by the PEI methodology, whereby the PEI is derived through the multiplication of indicators that are desired to have higher values and the division of indicators that are preferred to have lower values. Equal weights (ω) were assigned to all the functions, as all the indicators were deemed to be of equal significance. It should be noted that the value of ω is subject to the discretion of decision makers, who may set it based on their individual preferences.

As previously mentioned, the PEI can be used for both SO and MO problems. The example provided here pertains to a MO problem, while the next section will show the use of PEI in an SO application. In this example, a common MO nondominated performance metric is considered: the Hypervolume Indicator (HI), the Pareto Optimal Solutions (POSs) and NFE. A higher HI is desirable because it signifies a wider range of POSs. Having a higher number of POSs is also desirable. On the other hand, NFE serves as a reliable measure of computational complexity and is independent of the computer system. In this case, a lower NFE is preferred. The mathematical expression is represented by Equation (4)(4) $$\mathrm{\text{PEI}}=\left[\mathrm{H}{\mathrm{I}}^{{\omega }_1}\mathrm{P}\mathrm{\text{OS}}{\mathrm{s}}^{{\omega }_2}\right]/\mathrm{\text{NF}}{\mathrm{E}}^{{\omega }_3}$$(4) as follows: (4) $$\mathrm{\text{PEI}}=\left[\mathrm{H}{\mathrm{I}}^{{\omega }_1}\mathrm{P}\mathrm{\text{OS}}{\mathrm{s}}^{{\omega }_2}\right]/\mathrm{\text{NF}}{\mathrm{E}}^{{\omega }_3}$$(4)

The PEI is a versatile metric that can be tailored by researchers to align with their preferred evaluation criteria through the modification of equations. The addition of supplementary metrics can be incorporated into the equations based on whether higher or lower values are desired. As previously mentioned, metrics with higher desired outcomes are included in the numerator, while those with lower desired outcomes are placed in the denominator. By consolidating multiple performance indices into a single metric, this method streamlines the evaluation process and offers valuable insights for decision-making in complex optimisation scenarios.

6. Experiments

The selected test problem for investigation is the RDSP problem with an SO approach. As previously explained, the Enhanced Discrete Bees Algorithm (EDBA) proposed by Liu, Zhou, Pham, Xu, Ji, et al. (2018a) serves as the BBA employed to address the challenges of the robotic disassembly problem. The objective in this RDSP problem is to minimise the disassembly time, as expressed in Equation (5)(5) $$Z=\sum _{i=0}^{N_p-1}{t}_b\left(x_i\right)+\sum _{i=0}^{N_p-2}{t}_z\left(x_i,x_{i+1}\right)+\sum _{i=0}^{N_p-2}{t}_t\left(x_i,x_{i+1}\right)+\sum _{i=0}^{N_p-2}{m}_t\left(x_i,x_{i+1}\right)$$(5) . (5) $$Z=\sum _{i=0}^{N_p-1}{t}_b\left(x_i\right)+\sum _{i=0}^{N_p-2}{t}_z\left(x_i,x_{i+1}\right)+\sum _{i=0}^{N_p-2}{t}_t\left(x_i,x_{i+1}\right)+\sum _{i=0}^{N_p-2}{m}_t\left(x_i,x_{i+1}\right)$$(5) where

• Z is the total disassembly time.

• N_p is the number of total parts.

• t_b(x_i) is the basic time for disassembling part x_i.

• t_z(x_i, x_i+1) is the penalty time for disassembly direction changes between part x_i and x_i+1.

• t_t(x_i, x_i+1) is the penalty time for disassembly tool changes between part x_i and x_i + 1.

• m_t(x_i, x_i+1) is the moving time between part x_i and x_i+1.

6.1. Experimental setup and metrics

There are two commonly employed stopping criteria in optimisation algorithms: iterations or NFE. In this research, the stopping criteria for the algorithm are based on the number of iterations. The parameter settings used for EDBA are based on the recommended settings from Liu, Zhou, Pham, Xu, Ji, et al. (2018a). In their work, the authors applied population sizes of 10, 20, 30, 40 and 50 for Gear Pump A for 100, 200, 300 and 400 iterations. Similarly, for Gear Pump B, the same population numbers were employed between iterations 100 and 600. The values for the numbers of selected sites (m), elite sites (e), recruited bees for other selected sites (nsp) and recruited bees for elite sites (nep) were assigned as follows: 4, 1, 1 and 2, respectively. According to their findings, the best results were obtained for Gear Pump A in 300 iterations using a population of 20 and for Gear Pump B in 500 iterations using a population of 40 (Liu, Zhou, Pham, Xu, Ji, et al., 2018a). In this research, to ensure a fair comparison, the iteration and population sizes for all algorithms have been set to the same values.

To ensure a meaningful comparison with the best results obtained from EDBA, it was decided to use 100, 200, 300, 400 and 500 iterations for both gear pumps, along with population sizes of 21, 31, 41 and 51. It is important to note that in the context of the EDBA article, the term “population” refers to the number of scout bees (n). Therefore, when calculating the total population using the EDBA parameter settings and applying Equation (6)(6) $${\mathrm{\text{BA}}}_{\mathrm{\text{pop}}}=(e*\mathit{\text{nep}})+\left(\right(\mathrm{m }‐\mathrm{ }e\left)\mathit{\text{* nsp}}\right)+(\mathrm{n }‐\mathrm{ }m)$$(6) , the corresponding values for EDBA are 21, 31, 41 and 51 for the population sizes. (6) $${\mathrm{\text{BA}}}_{\mathrm{\text{pop}}}=(e*\mathit{\text{nep}})+\left(\right(\mathrm{m }‐\mathrm{ }e\left)\mathit{\text{* nsp}}\right)+(\mathrm{n }‐\mathrm{ }m)$$(6)

6.2. Experimental parameter setting

In this section, the determination of the experimental parameters employed in the study is presented. The objective was to find the optimal parameter settings for BA_F and assess its performance in comparison to EDBA. The steps taken to determine the optimal parameter settings for BA_F are as follows:

1. Experimental Design: Define the experimental parameters. The numbers of selected sites (m), the maximum number of recruited bees around the selected sites (maxnr) and the maximum number of re-visits (max_rv) are configured as depicted in Table 3.

2. Initial Runs: Perform 10 runs for each scenario, E1 to E100.

3. Top Results Identification: Identify the best results from the initial runs.

4. Iterations and Performance Validation: Evaluate the performance of the selected best results from the initial runs across 50 runs using 100 iterations.

5. Descriptive Statistics: Analyse the descriptive statistics of the results.

6. Best Result Selection: Identify and select the best result.

Table 3. Experimental design for BA_F.

		max_rv
m	maxnr	0	1	5	10	15
2	3	E1	E2	E3	E4	E5
	5	E6	E7	E8	E9	E10
	8	E11	E12	E13	E14	E15
	13	E16	E17	E18	E19	E20
	21	E21	E22	E23	E24	E25
3	3	E26	E27	E28	E29	E30
	5	E31	E32	E33	E34	E35
	8	E36	E37	E38	E39	E40
	13	E41	E42	E43	E44	E45
	21	E46	E47	E48	E49	E50
4	3	E51	E52	E53	E54	E55
	5	E56	E57	E58	E59	E60
	8	E61	E62	E63	E64	E65
	13	E66	E67	E68	E69	E70
	21	E71	E72	E73	E74	E75
5	3	E76	E77	E78	E79	E80
	5	E81	E82	E83	E84	E85
	8	E86	E87	E88	E89	E90
	13	E91	E92	E93	E94	E95
	21	E96	E97	E98	E99	E100

Since BA_F uses a ranking-based mechanism to give the number of bees sent to the selected sites, for comparison purposes, it is necessary to calculate the population, as shown in Equation (7)(7) $${\mathrm{\text{BA}}}_{\mathrm{F}}\mathrm{\text{population}}=\left(n-m\right)+\sum _{i=1}^m\left({\mathit{\text{nr}}}_i\right)$$(7) , to obtain comparable population sizes for BA_F and EDBA. Table 4 presents values calculated using m = 5 and various maxnr with a target maximum population of 51 bees. The maximum value of nr is nr1. (7) $${\mathrm{\text{BA}}}_{\mathrm{F}}\mathrm{\text{population}}=\left(n-m\right)+\sum _{i=1}^m\left({\mathit{\text{nr}}}_i\right)$$(7)

Table 4. Example of BA_F for a population of 51.

Experiments	n	maxnr = nr1	nr2	nr3	nr4	nr5	BAF population
E76	48	3	2	1	1	1	51
E81	44	5	3	2	1	1	51
E86	37	8	5	3	2	1	51
E91	25	13	8	5	3	2	51
E96	6	21	13	8	5	3	51

The experiment was then conducted 10 times on Gear Pump B for a maximum of 100 iterations each. Gear Pump B was selected as it was more complex than Gear Pump A. The parameter settings and average results of BA_F are presented in Table 5. Using the best parameter settings from Liu, Zhou, Pham, Xu, Ji, et al. (2018a), EDBA achieved a disassembly time of 150.113 with 5550 function evaluations.

Table 5. BA_F results of the initial runs (10 independent runs).

Experiments	m	maxnr	n	NFE	Disassembly time
					max_rv = 0	max_rv = 1	max_rv = 5	max_rv = 10	max_rv = 15
E1 to E5	2	3	48	5348	173.8317	171.5567	160.185	145.065	151.3192
E6 to E10	2	5	45	5345	164.3175	165.0625	151.9225	145.7392	151.4817
E11 to E15	2	8	40	5340	160.9383	154.3392	147.3933	148.6567	143.74
E16 to E20	2	13	32	5332	150.0267	151.0992	148.7533	142.4158	142.9308
E21 to E25	2	21	19	5319	148.09	146.1292	146.4492	144.4858	144.3008
E26 to E30	3	3	48	5448	172.5958	168.9592	155.4217	146.0192	155.5033
E31 to E35	3	5	44	5444	160.6658	167.5758	149.49	143.3275	146.0117
E36 to E40	3	8	38	5438	161.3008	156.2958	146.8142	146.0458	144.9183
E41 to E45	3	13	28	5428	153.625	148.8067	143.9108	144.1133	143.245
E46 to E50	3	21	12	5412	151.1475	145.6175	141.955	142.8275	143.6292
E51 to E55	4	3	48	5548	173.1675	173.7192	153.4317	149.655	147.2775
E56 to E60	4	5	44	5544	165.3075	169.1067	148.6542	145.0658	145.8825
E61 to E65	4	8	37	5537	152.553	160.473	146.28	145.689	144.548
E66 to E70	4	13	26	5526	149.4717	150.3658	143.9317	142.383	142.0917
E71 to E75	4	21	8	5508	148.7025	147.88	144.2342	143.1633	143.0742
E76 to E80	5	3	48	5648	173.423	177.248	154.938	144.268	147.38
E81 to E85	5	5	44	5644	162.987	165.988	149.543	144.179	148.528
E86 to E90	5	8	37	5637	156.544	157.107	144.883	143.512	143.217
E91 to E95	5	13	25	5625	149.994	152.092	144.622	143.033	141.015
E96 to E100	5	21	6	5606	151.319	150.675	143.292	142.158	147.656

Bold values indicate the consistent lowest value.

The preliminary findings indicate that all scenarios require fewer NFEs compared to the benchmark algorithm, EDBA. Furthermore, the majority of results across all scenarios demonstrate lower average disassembly times compared to EDBA. These results suggest that BA_F offers higher accuracy and efficiency. However, further experiments were necessary before definitive conclusions could be drawn. Based on the average results and NFE, experiments E48, 49, E50, E68, E69 and E70 were selected for further investigations through 50 independent runs. The statistical results obtained are shown in Table 6.

Table 6. Descriptive statistics of E48, E49, E50, E68, E69 and E70 (50 runs).

	N Statistic	Minimum statistic	Maximum statistic	Mean		Std. deviation statistic
				Statistic	Std. error
E48	50	135.32	154.15	142.3295	0.54402	3.84683
E49	50	136.68	158.32	142.1598	0.61608	4.35638
E50	50	135.32	154.04	142.1537	0.54545	3.85691
E68	50	137.25	151.38	143.0255	0.51073	3.61138
E69	50	135.32	148.13	141.6850	0.39180	2.77047
E70	50	135.32	153.49	141.6220	0.52304	3.69846
Valid N (listwise)	50

Table 6 shows that scenario E69 yielded a lower standard error of the mean (SEM). This lower SEM signifies a more accurate estimate of the population mean. Moreover, the standard deviation is also smaller, indicating reduced data variability. This characteristic is desirable regardless of whether the objective is to minimise or maximise a certain parameter. The preference is for results to cluster closely around the mean, as it signifies greater precision in achieving the intended outcome. Therefore, scenario E69 with m = 4, maxnr = 13, n = 26, max_rv = 10 emerges as the optimal choice. Fifty independent runs were subsequently conducted with the above parameter settings to compare BA_F and EDBA. The best parameter settings for BA_F specified a minimum population size of 31. Consequently, the comparative assessment with a population size 21 was eliminated.

6.3. Experimental results

Experiments were run 50 times with the maximum number of iterations set at 100, 200, 300, 400 and 500, while the population sizes used were 31, 41 and 51. For ease of reference, a grouping code was employed, where “100 31” indicates the experiment with a maximum of 100 iterations and a population size of 31. This coding scheme helps to conveniently identify the specific parameter settings used in each experiment. The best disassembly results using EDBA and BA_F are presented in Tables 7 and 8 for Gear Pump A and Gear Pump B, respectively. These tables display the sequence, direction, tool and total time of disassembly. The best result aligns with the best-known values for the case study, as reported in Liu, Zhou, Pham, Xu, Ji, et al. (2018a), who conducted experiments using 1000 runs and obtained values of 87.5731 and 135.3167 for Gear Pumps A and B, respectively. Tables 9 and 10 display the average values of the disassembly time, SEM, average NFE (Ave_NFE) and delta. Delta is the percentage difference between the average objective value (disassembly time) to the best-known value. The boxplots for Gear Pumps A and B are presented in Figures 2 and 3, respectively.

Figure 2. EDBA and BA_F results (Gear Pump A).

Figure 3. EDBA and BA_F results (Gear Pump B).

Table 7. Gear Pump A best results (EDBA and BA_F).

	Best disassembly output
Sequence	2-1-6-5-4-3-7-10-9-11-8-13-15-14-12
Direction	2-2-2-2-2-2-2-2-2-2-2-2-1-1-1
Tool	1-1-1-1-1-1-4-3-3-3-3-3-2-3-4
Time	89.5731
Sequence	2-1-6-5-4-3-7-10-11-9-8-13-15-14-12
Direction	2-2-2-2-2-2-2-2-2-2-2-2-1-1-1
Tool	1-1-1-1-1-1-4-3-3-3-3-3-2-3-4
Time	89.5731
Sequence	3-4-5-6-1-2-7-10-11-9-8-13-15-14-12
Direction	2-2-2-2-2-2-2-2-2-2-2-2-1-1-1
Tool	1-1-1-1-1-1-4-3-3-3-3-3-2-3-4
Time	89.5731
Sequence	3-4-5-6-1-2-7-10-9-11-8-13-15-14-12
Direction	2-2-2-2-2-2-2-2-2-2-2-2-1-1-1
Tool	1-1-1-1-1-1-4-3-3-3-3-3-2-3-4
Time	89.5731
Sequence	5-4-3-2-1-6-7-10-11-9-8-13-15-14-12
Direction	2-2-2-2-2-2-2-2-2-2-2-2-1-1-1
Tool	1-1-1-1-1-1-4-3-3-3-3-3-2-3-4
Time	89.5731
Sequence	5-4-3-2-1-6-7-10-9-11-8-13-15-14-12
Direction	2-2-2-2-2-2-2-2-2-2-2-2-1-1-1
Tool	1-1-1-1-1-1-4-3-3-3-3-3-2-3-4
Time	89.5731
Sequence	6-1-2-3-4-5-7-10-11-9-8-13-15-14-12
Direction	2-2-2-2-2-2-2-2-2-2-2-2-1-1-1
Tool	1-1-1-1-1-1-4-3-3-3-3-3-2-3-4
Time	89.5731
Sequence	6-1-2-3-4-5-7-10-9-11-8-13-15-14-12
Direction	2-2-2-2-2-2-2-2-2-2-2-2-1-1-1
Tool	1-1-1-1-1-1-4-3-3-3-3-3-2-3-4
Time	89.5731

Notes: Direction: 1 = Y+ direction, 2 = Y– direction. Tool: 1 = Spanner-I, 2 = Spanner-II, 3 = Gripper-I, 4 = Gripper-II.

Table 8. Gear Pump B best result (EDBA and BA_F).

	Best disassembly output
Sequence	24-22-20-23-21-19-1-2-3-4-5-6-7-13-10-9-8-12-14-15-16-17-18-11
Direction	1-1-1-1-1-1-2-2-2-2-2-2-2-2-2-2-2-2-2-2-2-2-1-1
Tool	3-3-2-3-3-2-1-1-1-1-1-1-5-4-4-4-4-5-4-4-4-4-4-4
Time	135.3167
Sequence	24-22-20-23-21-19-1-6-5-4-3-2-7-13-10-9-8-12-14-15-16-17-18-11
Direction	1-1-1-1-1-1-2-2-2-2-2-2-2-2-2-2-2-2-2-2-2-2-1-1
Tool	3-3-2-3-3-2-1-1-1-1-1-1-5-4-4-4-4-5-4-4-4-4-4-4
Time	135.3167

Notes: Direction: 1 = Y+ direction, 2 = Y– direction. Tool: 1 = Spanner-I, 2 = Spanner-II, 3 = Spanner-III, 4 = Gripper-I, 5 = Gripper-II.

Table 9. EDBA and BA_F average results for Gear Pump A.

Iteration_ population	EDBA				BAF
	Mean	SEM	Delta (%)	Ave_NFE	Mean	SEM	Delta (%)	Ave_NFE
100_31	90.5237	0.33	0.95	3530	90.0960	0.213	0.52	3506
200_31	90.0357	0.27	0.46	7030	89.7208	0.064	0.15	7006
300_31	89.8065	0.23	0.23	10,530	89.5731	0.000	0.00	10,506
400_31	89.6047	0.03	0.03	14,030	89.5731	0.000	0.00	14,006
500_31	89.7387	0.14	0.17	17,530	89.5731	0.000	0.00	17,506
100_41	90.6014	0.32	1.03	4540	89.7404	0.137	0.17	4516
200_41	89.8637	0.21	0.29	9040	89.5731	0.000	0.00	9016
300_41	89.6364	0.04	0.06	13,540	89.5731	0.000	0.00	13,516
400_41	89.5747	0.00	0.00	18,040	89.5731	0.000	0.00	18,016
500_41	89.5731	0.00	0.00	22,540	89.5731	0.000	0.00	22,516
100_51	89.7875	0.08	0.21	5550	89.5731	0.000	0.00	5526
200_51	89.7066	0.08	0.13	11,050	89.5731	0.000	0.00	11,026
300_51	89.6047	0.03	0.03	16,550	89.5731	0.000	0.00	16,526
400_51	89.6047	0.03	0.03	22,050	89.5731	0.000	0.00	22,026
500_51	89.6047	0.03	0.03	27,550	89.5731	0.000	0.00	27,526

Bold values indicate the consistent lowest value.

Table 10. EDBA and BA_F results for Gear Pump B.

Iteration_ population	EDBA				BAF
	Mean	SEM	Delta (%)	Ave_NFE	Mean	SEM	Delta (%)	Ave_NFE
100_31	150.9088	1.36	15.59	3530	146.6653	1.08	11.35	3506
200_31	146.918	0.90	11.60	7030	140.9998	0.36	5.68	7006
300_31	145.5225	0.91	10.21	10,530	140.5303	0.42	5.21	10,506
400_31	142.217	0.70	6.90	14,030	139.7332	0.35	4.42	14,006
500_31	142.0165	0.45	6.70	17,530	139.1772	0.27	3.86	17,506
100_41	148.5027	0.97	13.19	4540	142.6655	0.63	7.35	4516
200_41	143.5897	0.56	8.27	9040	140.5633	0.34	5.25	9016
300_41	142.4695	0.63	7.15	13,540	139.8480	0.19	4.53	13,516
400_41	142.0715	0.81	6.75	18,040	139.0933	0.23	3.78	18,016
500_41	141.9507	0.54	6.63	22,540	138.6493	0.21	3.33	22,516
100_51	148.5675	1.08	13.25	5550	141.6850	0.39	6.37	5526
200_51	143.7058	0.59	8.39	11,050	139.9833	0.29	4.67	11,026
300_51	142.9958	0.69	7.68	16,550	139.2005	0.15	3.88	16,526
400_51	142.0362	0.60	6.72	22,050	139.2242	0.15	3.91	22,026
500_51	142.0932	0.64	6.78	27,550	138.8473	0.14	3.53	27,526

6.4. Statistical performance metric results

The statistical assumptions, as discussed in Section 5.1, were checked prior to conducting the subsequent statistical test. From the boxplots, it is evident that the data do not follow a normal distribution. The statistical tests for normality and homogeneity, as presented in Appendix A, confirm that the results deviate from both normality and homogeneity assumptions. Therefore, as the data violate the statistical assumptions required for conducting parametric tests, nonparametric tests are employed. The significance of differences between the mean ranks of various iterations and population sizes on EDBA and BA_F results is examined using the Kruskal–Wallis ANOVA test, as depicted in Tables 11 and 12. The test results indicate the rejection of the null hypothesis (p < 0.05), providing strong evidence of a difference between at least one pair of iterations and population sizes. Consequently, a post hoc test, namely, the Dunn–Sidak test, is conducted. The results of the Dunn–Sidak tests are presented in Figures 4 and 5.

Figure 4. Dunn–Sidak test results (Gear Pump A).

Group 1 = 100_31 EDBA, Group 2 = 100_31 BA_F, Group 3 = 200_31 EDBA, Group 4 = 200_31 BA_F, …, Group 29 = 500_51 EDBA, Group 30 = 500_51 BA_F.

Figure 5. Dunn–Sidak test results (Gear Pump B).

Group 1 = 100_31 EDBA, Group 2 = 100_31 BA_F, Group 3 = 200_31 EDBA, Group 4 = 200_31 BA_F, …, Group 29 = 500_51 EDBA, Group 30 = 500_51 BA_F.

Table 11. Kruskal–Wallis test results (Gear Pump A).

Kruskal–Wallis ANOVA table
Source	SS	df	MS	Chi-sq	Prob > Chi-sq
Columns	8.18798e + 06	29	282,344	263.11	1.4359e-39
Error	3.84605e + 07	1470	26,163.6
Total	4.66485e + 07	1499

Table 12. Kruskal–Wallis test results (Gear Pump B).

Kruskal–Wallis ANOVA table
Source	SS	df	MS	Chi-sq	Prob > Chi-sq
Columns	9.89515e + 07	29	3,412,119.9	533.75	3.24619e-94
Error	1.78946e + 08	1470	121,731.9
Total	2.77897e + 08	1499

The Dunn–Sidak results offer detailed insights into the groups that display variation in performance based on iteration within each group size. Specifically, for Gear Pump A, the results show significant differences of groups 1 and 11 (100_31 EDBA and 100_41 EDBA) to the other groups. It is noteworthy that BA_F yielded consistent outcomes across all parameters, as evidenced by the groups with even-numbered labels. However, in the case of Gear Pump B, which has a larger dataset in comparison to Gear Pump A, both EDBA and BA_F displayed more pronounced disparities across all groups.

These findings emphasise the importance of selecting optimal parameter settings for disassembly time optimisation, taking into account both iteration size and population size. The observed variability in performance highlights the need for careful parameter selection to achieve optimal results.

In the subsequent step, groups with significant differences from the previous statistical test are eliminated, and all iterations and population sizes are combined to focus on the groups that demonstrate similar performance with the lowest mean objective values. This step ensures that the final test results only include the groups that yield the best results. The outcomes of the final statistical tests can be found in Figure 6 and Tables 13 and 14.

Figure 6. EDBA and BA_F boxplot final results.

Table 13. Kruskal–Wallis final results (Gear Pump A).

Kruskal–Wallis ANOVA table
Source	SS	df	MS	Chi-sq	Prob > Chi-sq
Columns	103,511.1	21	4929.1	27.21	0.1641
Error	4,077,771.9	1078	3782.72
Total	4,181,283	1099

Table 14. Kruskal–Wallis final results (Gear Pump B).

Kruskal–Wallis ANOVA table
Source	SS	df	MS	Chi-sq	Prob > Chi-sq
Columns	111,867.9	6	18,644.7	11.44	0.0756
Error	3,299,830.1	343	9620.5
Total	3,411,698	349

6.5. PEI results

As previously highlighted, PEI is designed to serve as a versatile index. In this experiment, the objective is to minimise the disassembly time while also aiming for the lowest value of the NFE. To achieve this, both metrics are placed in the denominator during the PEI calculation, as shown in Equation (8)(8) $$\mathrm{\text{PEI}}=1/\left[\mathrm{\text{Ob}}{\mathrm{j}}^{{\omega }_1}\mathrm{NF}{\mathrm{E}}^{{\omega }_2}\right]$$(8) , with equal weight assigned to all the metrics. It is important to note that since none of the metrics are desirable to have the highest value, a value of 1 is assigned to the numerator. (8) $$\mathrm{\text{PEI}}=1/\left[\mathrm{\text{Ob}}{\mathrm{j}}^{{\omega }_1}\mathrm{NF}{\mathrm{E}}^{{\omega }_2}\right]$$(8)

To ensure a meaningful comparison, only the groups that exhibited no significant differences in means were considered for the calculation. Figures 7(a) and 7(b) and Appendix B show the PEI calculation results, scaled to units of 10⁻⁶. The figures provide a visual representation of the PEI values for the selected groups, facilitating comparison of their performance. Note that, although this work employs both SPM and PEI, they can be used on their own as they are independent performance measures. SPM is more suitable for a single-objective (SO) problem, while PEI can be utilised for both SO and multi-objective (MO) optimisation.

Figure 7. PEI (histogram) and Average Disassembly Time (dot): Higher PEI and Lower Time are Better. Higher PEI represents highest performance and Lower time represents best disassembly time (desirable).

Note: The line shows the highest PEI; EDBA is represented by blue and BA_F by purple colour

7. Discussion

The experimentation process for determining the optimal parameters setting for BA_F was conducted systematically in six steps, including statistical tests for result analysis. This approach was implemented to ensure that the best parameters are obtained for the subsequent comparison with EDBA. The best parameters for the RDSP problem were found to be m = 4, maxnr = 13, max_rv = 10. The study conducted by Liu, Zhou, Pham, Xu, Ji, et al. (2018a) used a population range of 21–51. The number of scout bees (n) in BA_F was chosen in the range 6–26 to yield a similar population size of 31–51. As previously mentioned, the population size 21 was eliminated from the final experiments. After determining the optimal parameters settings, 50 independent runs were conducted for both EDBA and BA_F.

The minimum disassembly time can be achieved using both EDBA and BA_F, as presented in Tables 6 and 7 for both gear pumps. These tables give the sequence, direction, tools and total time required to disassemble the gear pumps. However, the average results reveal differences between EDBA and BA_F. The average results of EDBA and BA_F, presented in Tables 8 and 9, indicate that BA_F outperforms EDBA. Examination of the outcomes for Gear Pump A shows that BA_F consistently exhibits superior performance in comparison to EDBA across all 15 parameter configurations. BA_F achieves a perfect accuracy rate in 12 instances, as evidenced by a delta value of 0%, indicating that the accuracy is 100%. In contrast, EDBA only achieves 100% accuracy once with the same settings. The SEM and delta exhibit a decreasing trend as larger parameter values are used. EDBA performed best in 500 iterations and with a population size of 41, giving an average NFE of 22,540. In contrast, BA_F demonstrated the best performance in almost all cases. The average NFE was 5526 that was the lowest NFE required to achieve optimal outcomes in 100 iterations and with a population size of 51. Similarly, on Gear Pump B, it is evident that BA_F outperforms EDBA for all population sizes and iterations. Figure 2 displays the boxplot results for Gear Pump A, indicating that BA_F achieves better results. The algorithm demonstrates the ability to find the minimum disassembly time with minimal data spread, particularly at higher iterations for all tested population sizes. In Figure 3, the boxplot results for Gear Pump B exhibit improved data visualisation, allowing for better observation of data spread. In this boxplot, it becomes even clearer that BA_F outperforms EDBA in terms of mean, minimum disassembly time and data spread across all tested iterations and population sizes.

As mentioned earlier, research by Liu, Zhou, Pham, Xu, Ji, et al. (2018a) indicates that EDBA performs best on Gear Pump A in 300 iterations and with a population of 21 (NFE = 7520), while Gear Pump B achieves optimal performance in 500 iterations with a population of 41 (NFE = 22,540). However, it is important to emphasise that the conclusions drawn in this research are not solely based on statistically descriptive results. To ensure the reliability and validity of the findings, additional statistical tests are imperative to validate the parameter settings that yield the best performance. Relying solely on the statistically descriptive results may not provide the conclusive evidence required to draw robust conclusions. Since the experiments involved 50 independent runs, it is necessary to assess whether parametric statistical tests can be employed, as discussed in Section 5.1. The statistical assumptions for parametric tests include having a sample size of more than 30, data following a normal distribution, and data belonging to the same populations. However, the boxplot results indicate that the data do not follow a normal distribution, and thus, normality and homogeneity tests were conducted. The results confirm that the data violate the assumptions of normality and homogeneity, necessitating the use of nonparametric tests despite the sample size exceeding 30. As previously mentioned, the Kruskal–Wallis test indicated a significant difference among the iterations within the same population size groups (see Tables 11 and 12). To further investigate these differences, the Dunn–Sidak post hoc test was employed to identify specific groups that exhibited statistically significant variations in mean values.

Figures 4 and 5 clearly demonstrate noticeable differences, particularly in the outcome obtained from EDBA. Subsequently, the groups with significant differences are removed, and the statistical tests are repeated. The boxplot results, as shown in Figure 6, visually indicate that the data distributions in each group are similar. The Kruskal–Wallis test (see Tables 13 and 14) confirmed that there was no statistically significant difference among the groups (p > 0.05). Therefore, post hoc tests are not needed, as the Kruskal–Wallis test has already confirmed the similarity among the groups.

The final step involves the utilisation of the proposed index, as detailed in Section 5.2. In this study, the index is calculated based on the disassembly time and NFE. The calculation for the index has been previously explained in the previous section, and the obtained results, depicted in Figure 7, exhibit the PEI values in conjunction with the average disassembly time. By examining Figure 7, the parameter settings that yield the highest PEI can be identified, while the average disassembly time has been graphically presented for analytical purposes. As mentioned previously, the PEI is designed to be interpreted such that a higher index indicates better performance. The PEI results reveal that the best PEI of 9.87 × 10⁻⁶ for Gear Pumps A was obtained in 100 iterations and with a population size of 41 for BA_F while for EDBA the optimum PEI of 4.23 × 10⁻⁶ was found in 300 iterations and with a population size of 31. The results indicate that BA_F with a smaller NFE can produce statistically similar outcomes to those with a higher NFE. Figure 7(a) displays the average disassembly time as a dot, indicating that BA_F consistently outperforms EDBA in all instances. Directly comparing similar parameter settings for the two algorithms reveals that BA_F consistently yields higher PEI values – a highly desirable outcome for both gear pumps. In terms of average disassembly time, BA_F yields lower values, aligning with the desired objective of minimising the disassembly time. The results obtained from the PEI analysis of Gear Pump B (Figure 7(b)) clearly indicate that BA_F outperforms EDBA. The optimal results were achieved after 400 iterations and with a population size of 31, yielding a PEI value of 2.04 × 10⁻⁶.

These findings highlight the robustness of the chosen parameter settings, which cannot be easily discerned based on statistical description alone. It is evident from the results that the average disassembly time exhibits a pattern of decreasing values with higher population size and iteration. However, when considering the average disassembly time and the PEI together, it becomes apparent that the parameter settings yielding the lowest average disassembly time may not always be the optimal ones. This is because the NFE also plays a crucial role, as a higher NFE indicates that the algorithm requires more evaluations to achieve the best results. Therefore, a balance between the average disassembly time and the NFE needs to be considered to determine the best parameter settings.

By utilising SPM and PEI, it becomes possible to identify the optimal parameter settings and evaluate the performance of the algorithms. From a statistical standpoint, it is evident that BA_F outperforms EDBA with BA_F demonstrating a faster convergence (smaller NFE). This observation is further supported by the average results presented in Tables 9 and 10, where it can be seen that BA_F consistently achieves lower disassembly times across various iterations and population sizes. BA_F consistently outperforms EDBA for both gear pumps, particularly in larger datasets, such as Gear Pump B, demonstrating that it is capable of handling more complex problems. Moreover, the statistical tests and PEI calculations consistently indicate that BA_F outperforms EDBA. These findings provide valuable insights into the comparative performance of the algorithms and highlight the advantages of employing BA_F in terms of achieving lower disassembly times and higher PEI values.

Furthermore, the statistical analysis has revealed interesting findings, particularly evident in the statistical results presented in Table 9 and the normality test detailed in Appendix Table A.1, focussing on the case study of Gear Pump A. It is worth noting that the best-known solution for Gear Pump A stands at 89.5731. As illustrated in Table 9, BAF demonstrates remarkable accuracy across various combinations of iterations and population sizes, consistently yielding a mean value of 89.5731, indicating a 100% accuracy rate. Upon closer examination in Appendix Table A.1, it becomes evident that the data cannot conform to a normal distribution due to its precise alignment with the best-known solution. This non-normality aligns with our analysis, which suggests that in optimisation problems focussed on minimising or maximising objectives, deviation from the normal distribution is to be expected. In this context, it is crucial to emphasise that in problems where the objectives are centred on minimisation or maximisation, truly accurate results align precisely with a single point, inherently deviating from a normal distribution. Similarly, when results are near-accurate, they tend to cluster closely around that optimal point, leading to a skewed distribution pattern. As a result, we recommend that future researchers incorporate normality tests and examine normal graphs. It is anticipated that such analyses will confirm the departure from a normal distribution, particularly when results are highly accurate. In such instances, nonparametric tests should be preferred over parametric tests due to the observed accuracy of the results.

8. Conclusion

The enhancement of the BA was achieved by reducing the number of parameters to 4 and by simplifying the algorithm steps, removing elite sites and incorporating the Fibonacci sequence to guide the recruitment of bees for exploring the local search space. The SPM and PEI analyses conducted in this study indicate that BA_F exhibits better performance than the benchmark EDBA in this RDSP problem. EDBA has been compared against other algorithms, including EDBA variants (EDBA without mutation operator (EDBA-WMO)), GA with Precedence Preserving Crossover (GA-PPX), and Self-Adaptive Simplified Swarm Optimisation (SASSO), and has demonstrated better performance. Therefore, the fact that BA_F outperforms EDBA in this study means that BA_F also outperforms these other algorithms. This further strengthens the evidence supporting the effectiveness of BA_F as a robust and efficient algorithm in the context of RDSP. This conclusion is supported by several key factors, including the lower NFE values and the ability to achieve a minimum disassembly time. BA_F consistently outperforms EDBA on both gear pumps. These findings highlight the effectiveness of BA_F as an optimal choice for achieving efficient and effective RDSP. The SPM and PEI are both valuable tools in the field of algorithm analysis. These tools serve as effective means to identify optimal parameter settings and evaluate the performance of algorithms. SPM allows for statistical comparison and analysis of algorithm performance based on specific metrics, providing insights into the differences among various parameter settings. The PEI is a versatile metric that combines multiple metrics, such as disassembly time and NFE, into a single index to provide a comprehensive evaluation of algorithm performance. By combining SPM and PEI, researchers and practitioners can make informed decisions regarding parameter settings and algorithm selection to optimise their outcomes in various domains. Moreover, the validation of the proposed approach through a case study provides evidence of the algorithm’s robustness and effectiveness in successfully solving the robotic disassembly planning problem. Although this research specifically focusses on addressing the challenges of the RDSP problem, it is important to consider the wider applicability of the proposed algorithm in future research. Future studies will explore and compare the proposed enhancements with the BA2 approach in a broader context, encompassing both continuous and combinatorial problems. The parameters used in this study are smaller compared to those in TSP and VRP problems, which typically employ population sizes of 200–400 or even larger. This suggests the need for further research on the potential applications of BA_F in other contexts. This would enable a comprehensive evaluation and understanding of the algorithm’s performance across different problem domains. Furthermore, investigating analogies with other features of honeybees, such as their learning process, olfactory system and hive organisation, could lead to more improvements to the BA. The algorithm can also be utilised for parameter tuning in machine learning including deep learning, thereby increasing its usefulness as an advanced artificial intelligence tool. This extension and exploration of the BA in different contexts will contribute to its versatility and applicability in various optimisation domains. In addition, the findings shows that BA_F demonstrated exceptional accuracy in optimising Gear Pump A, leading to a non-normal distribution due to precise results. Therefore, for such an accurate and near-accurate result, the null hypothesis that the distribution is normal is expected to be rejected, and non-parametric tests should be conducted. Notably, none of the previous research in this area has observed this phenomenon, underscores the contribution to existing knowledge that results from this research. The MATLAB code of BAF for RDSP is available in github. The github address is given in Appendix C.

Abbreviations
AHGA	=	Automated Hybrid Genetic Algorithm
AVNS	=	Adaptive Variable Neighbourhood Search
BA	=	Bees Algorithm
BAF	=	Fibonacci Bees Algorithm
BA2	=	Bees Algorithm with 2 parameters
BBA	=	Basic Bees Algorithm
BLSA	=	Brand Local Search Algorithm
CDDO	=	Child Drawing Development Optimisation
DoE	=	Design of Experiments
EDBA	=	Enhanced Discrete Bees Algorithm
ELS	=	Evolutionary Local Search
EoL	=	End-of-Life
FIA	=	Fibonacci Indicator Algorithm
FSOA	=	Fibonacci Sequence-based Optimisation Algorithm
FSQGA	=	Fibonacci Sequence-based Quantum GA
FTO	=	Fibonacci Tree Optimisation
GA	=	Genetic Algorithm
GA-PPX	=	GA with Precedence Preserving Crossover
GROM	=	Golden Ratio Optimisation Method
GWO	=	Grey Wolf Optimiser
HI	=	Hypervolume Indicator
IWPA	=	Improved Wolf Pack Algorithm
MA	=	Memetic Algorithm
MCDM	=	Multiple-Criteria Decision-Making
MO	=	Multi-Objective
MOFA	=	Multi-Objective Fibonacci-Based Algorithm
NFE	=	Number of Function Evaluations
NP	=	Non-Deterministic Polynomial
Opp-PSOGWO	=	Opposition-based Learning PSO GWO
PEI	=	Performance Evaluation Index
POSs	=	Pareto Optimal Solutions
PSO	=	Particle Swarm Optimisation
RDSP	=	Robotic Disassembly Sequence Planning
SASSO	=	Self-Adaptive Simplified Swarm Optimisation
SEM	=	Standard Error of the Mean
SSA	=	Salp Swarm Algorithm
SBA	=	Standard Bees Algorithm
SO	=	Single-Objective
SPM	=	Statistical Performance Metric
TSP	=	Travelling Salesman Problem
VRP	=	Vehicle Routing Problem
WPA	=	Wolf Pack Algorithm

Acknowledgements

The authors would like to extend their sincere gratitude to Dr Jiayi Liu for generously providing the EDBA code for our research. The author would also like to thank the editor and the reviewers for their constructive feedback, which helped to improve the manuscript.

Disclosure statement

The authors declare no competing interests.

Additional information

Funding

This work was partly funded by the Engineering and Physical Sciences Research Council (EPSRC), UK, grant no. EP/N018524/1. Natalia Hartono’s doctoral study is supported by the Indonesian Endowment Fund for Education (LPDP), Ministry of Finance, Republic of Indonesia, grant no. 201908222215155 under BUDI LN scheme.

Notes on contributors

Natalia Hartono

Natalia Hartono earned her doctoral degree in Mechanical Engineering from the University of Birmingham, UK, in 2023, specialising in intelligent optimisation, the bees algorithm, sustainability, robotic disassembly, and remanufacturing. She has made contributions to these fields, with publications in reputable journals, research poster and thesis presentation competitions, and presentations at international conferences. She is part of the Bees Algorithm Research Group and Autonomous Remanufacturing Group at the University of Birmingham. Dr Hartono is well-known as the co-chair of the International Workshop Series on the Bees Algorithm and Its Applications and the co-editor of the book “Intelligent Production and Manufacturing Optimisation – The Bees Algorithm Approach” published by Springer.

Appendix A.

Statistical test results

 

Table A1. Gear Pump A normality test.

Test of normality
		Kolmogorov–Smirnov^a			Shapiro–Wilk
	Iter_pop	Statistic	df	Sig.	Statistic	df	Sig.
GearA_EDBA	100_31	0.364	50	<0.001	0.461	50	<0.001
	100_41	0.423	50	<0.001	0.522	50	<0.001
	100_51	0.482	50	<0.001	0.406	50	<0.001
	200_31	0.480	50	<0.001	0.269	50	<0.001
	200_41	0.478	50	<0.001	0.197	50	<0.001
	200_51	0.531	50	<0.001	0.248	50	<0.001
	300_31	0.517	50	<0.001	0.127	50	<0.001
	300_41	0.540	50	<0.001	0.198	50	<0.001
	300_51	0.536	50	<0.001	0.125	50	<0.001
	400_31	0.536	50	<0.001	0.125	50	<0.001
	400_41	0.536	50	<0.001	0.125	50	<0.001
	400_51	0.536	50	<0.001	0.125	50	<0.001
	500_31	0.528	50	<0.001	0.167	50	<0.001
	500_41	–	50	–	–	50	–
	500_51	0.536	50	<0.001	0.125	50	<0.001
GearA_BAF	100_31	0.456	50	<0.001	0.395	50	<0.001
	100_41	0.509	50	<0.001	0.169	50	<0.001
	100_51	–	50	–	–	50	–
	200_31	0.527	50	<0.001	0.349	50	<0.001
	200_41	-	50	–	–	50	–
	200_51	–	50	–	–	50	–
	300_31	–	50	–	–	50	–
	300_41	–	50	–	–	50	–
	300_51	–	50	–	–	50	–
	400_31	–	50	–	–	50	–
	400_41	–	50	–	–	50	–
	400_51	–	50	–	–	50	–
	500_31	–	50	–	–	50	–
	500_41	–	50	–	–	50	–
	500_51	–	50	–	–	50	–

^aLilliefors Significance Correction.

Table A2. Gear Pump B normality test.

Test of normality
		Kolmogorov–Smirnov^a			Shapiro–Wilk
	Iter_pop	Statistic	df	Sig.	Statistic	df	Sig.
GearB_EDBA	100_31	0.129	50	0.036	0.889	50	<0.001
	100_41	0.123	50	0.058	0.928	50	0.005
	100_51	0.144	50	0.011	0.915	50	0.002
	200_31	0.130	50	0.035	0.898	50	<0.001
	200_41	0.165	50	0.002	0.866	50	<0.001
	200_51	0.130	50	0.033	0.906	50	<0.001
	300_31	0.191	50	<0.001	0.825	50	<0.001
	300_41	0.159	50	0.003	0.841	50	<0.001
	300_51	0.183	50	<0.001	0.814	50	<0.001
	400_31	0.219	50	<0.001	0.735	50	<0.001
	400_41	0.278	50	<0.001	0.521	50	<0.001
	400_51	0.200	50	<0.001	0.838	50	<0.001
	500_31	0.178	50	<0.001	0.777	50	<0.001
	500_41	0.183	50	<0.001	0.850	50	<0.001
	500_51	0.229	50	<0.001	0.717	50	<0.001
GearB_BAF	100_31	0.146	50	0.009	0.900	50	<0.001
	100_41	0.171	50	<0.001	0.866	50	<0.001
	100_51	0.145	50	0.010	0.936	50	0.010
	200_31	0.110	50	0.179	0.971	50	0.255
	200_41	0.225	50	<0.001	0.838	50	<0.001
	200_51	0.240	50	<0.001	0.870	50	<0.001
	300_31	0.175	50	<0.001	0.911	50	0.001
	300_41	0.264	50	<0.001	0.865	50	<0.001
	300_51	0.335	50	<0.001	0.777	50	<0.001
	400_31	0.114	50	0.123	0.963	50	0.116
	400_41	0.268	50	<0.001	0.842	50	<0.001
	400_51	0.350	50	<0.001	0.770	50	<0.001
	500_31	0.167	50	0.001	0.948	50	0.029
	500_41	0.282	50	<0.001	0.858	50	<0.001
	500_51	0.395	50	<0.001	0.708	50	<0.001

^aLilliefors Significance Correction.

Table A3. Gear Pumps A and B homogeneity test.

Test of homogeneity of variance
		Levene statistic	df1	df2	Sig.
GearA_EDBA	Based on mean	11.653	14	735	<0.001
	Based on median	3.872	14	735	<0.001
	Based on median and with adjusted df	3.872	14	262.360	<0.001
	Based on trimmed mean	7.551	14	735	<0.001
GearA_BAF	Based on mean	17.344	14	735	<0.001
	Based on median	4.337	14	735	<0.001
	Based on median and with adjusted df	4.337	14	94.091	<0.001
	Based on trimmed mean	8.133	14	735	<0.001
GearB_EDBA	Based on mean	7.254	14	735	<0.001
	Based on median	5.442	14	735	<0.001
	Based on median and with adjusted df	5.442	14	539.180	<0.001
	Based on trimmed mean	6.670	14	735	<0.001
GearB_BAF	Based on Mean	29.396	14	735	<0.001
	Based on median	22.333	14	735	<0.001
	Based on median and with adjusted df	22.333	14	232.077	<0.001
	Based on trimmed mean	26.876	14	735	<0.001

Appendix B.

Final PEI values

Table B1. PEI values – Gear Pump A.

Iteration_population Algorithm	Average disassembly time	Average NFE	PEI (x 0.000001)
300_31 EDBA	89.8065	10,530	4.23
300_31 BAF	89.5731	10,506	4.25
400_31 EDBA	89.6047	14,030	3.18
400_31 BAF	89.5731	14,006	3.19
500_31 EDBA	89.7387	17,530	2.54
500_31 BAF	89.5731	17,506	2.55
100_41 BAF	89.7404	4516	9.87
200_41 BAF	89.5731	9016	4.95
300_41 EDBA	89.6364	13,540	3.30
300_41 BAF	89.5731	13,516	3.30
400_41 EDBA	89.5747	18,040	2.48
400_41 BAF	89.5731	18,016	2.48
500_41 EDBA	89.5731	22,540	1.98
500_41 BAF	89.5731	22,516	1.98
100_51 BAF	89.5731	5526	8.08
200_51 BAF	89.5731	11,026	4.05
300_51 EDBA	89.6047	16,550	2.70
300_51 BAF	89.5731	16,526	2.70
400_51 EDBA	89.6047	22,050	2.02
400_51 BAF	89.5731	22,026	2.03
500_51 EDBA	89.6047	27,550	1.62
500_51 BAF	89.5731	27,526	1.62

Table B2. PEI values – Gear pump B.

Iteration_population Algorithm	Average disassembly time	Average NFE	PEI (x 0.000001)
400_31 BAF	139.7332	14,006	2.04
500_31 BAF	139.1772	17,506	1.64
400_41 BAF	139.0933	18,016	1.60
500_41 BAF	138.6493	22,516	1.28
300_51 BAF	139.2005	16,526	1.74
400_51 BAF	139.2242	22,026	1.30
500_51 BAF	138.8473	27,526	1.05

Appendix C.

MATLAB Code

https://github.com/NataliaHartonoFung/Fibonacci-Bees-Algorithm-disassembly-sequence.