Research on the Human–Robot Collaborative Disassembly Line Balancing of Spent Lithium Batteries with a Human Factor Load

Abstract

The disassembly of spent lithium batteries is a prerequisite for efficient product recycling, the first link in remanufacturing, and its operational form has gradually changed from traditional manual disassembly to robot-assisted human–robot cooperative disassembly. Robots exhibit robust load-bearing capacity and perform stable repetitive tasks, while humans possess subjective experiences and tacit knowledge. It makes the disassembly activity more adaptable and ergonomic. However, existing human–robot collaborative disassembly studies have neglected to account for time-varying human conditions, such as safety, cognitive behavior, workload, and human pose shifts. Firstly, in order to overcome the limitations of existing research, we propose a model for balancing human–robot collaborative disassembly lines that take into consideration the load factor related to human involvement. This entails the development of a multi-objective mathematical model aimed at minimizing both the cycle time of the disassembly line and its associated costs while also aiming to reduce the integrated smoothing exponent. Secondly, we propose a modified multi-objective fruit fly optimization algorithm. The proposed algorithm combines chaos theory and the global cooperation mechanism to improve the performance of the algorithm. We add Gaussian mutation and crowding distance to efficiently solve the discrete optimization problem. Finally, we demonstrate the effectiveness and sensitivity of the improved multi-objective fruit fly optimization algorithm by solving and analyzing an example of Mercedes battery pack disassembly.

1. Introduction

With the burgeoning consciousness of environmentally friendly, low-carbon initiatives, coupled with the rapid and substantial advancements in lithium battery technology, as well as the resolute implementation of pertinent policies, the electric vehicle market is currently in the midst of an unprecedented upsurge in growth [1,2]. According to the IEA 2023, as shown in Figure 1, the global ownership of new energy vehicles exceeded 26 million in 2022. China leads the global market with a staggering 13.8 million electric vehicles, accounting for half of the world’s total. Commonly, lithium batteries of electric vehicles have a service life of 5–8 years [3]. Once the power battery’s capacity diminishes to 80% of its original level, it is deemed as a spent lithium battery (SLIB) [4]. According to China’s Ministry of Public Security, it is estimated that in 2021, China had a staggering 15 GWh of retired waste power batteries, totaling over 70,000 tons. By 2030, China’s retirement SLIB of electric vehicles will reach 200–500 million tons, with an annual retirement rate of 15% [5,6]. Although SLIB may no longer be capable of powering electric vehicles, they still possess inherent utility. With the direct disposal rather than recycling of these batteries, there is undeniably a substantial squandering of valuable resources [7]. SLIB from electric vehicles can be used for power grids, non-interruptible power supplies, and additional energy storage devices. When the battery capacity falls further, the cell or module is re-supplied. Therefore, recycling waste power cells is considered an essential means to solve the problem of SLIB.

Disassembly sequence planning has a relatively complete theoretical system and research methods [8]. In recent years, as a result of the rapid advancements in power cell technology, several researchers have undertaken studies on the disassembly sequence planning of power cell packs, modules, and monomers and have successfully yielded significant results. Alfaro et al. 2020 proposed a computational model for SLIB disassembly [9]. The model features a cost assessment of four scenarios for the reuse, remanufacturing, recycling, and disposal of the target parts during the disassembly process. Wegener et al. 2014 obtained the battery disassembly sequence based on the priority relation matrix, but the evaluation of the disassembly sequence was relatively simple [10]. Harper et al. 2019 summarized the research status and technical difficulties of battery disassembly planning in SLIB recycling for electric vehicles [11]. Schwarz et al. 2018 used a predetermined practice method to conduct the virtual disassembly of power battery packs, helping decision makers integrate information such as time, risk, complexity, and operating tools in the disassembly process [12]. Cong et al. 2023 proposed a multi-objective optimization mathematical model for SLIB disassembly that considered economy, time, safety, and energy consumption [13]. They designed a hybrid gene-firework algorithm based on priority graphs to solve the planning problem. To obtain the optimal disassembly sequence, Xiao et al. 2022 presented a Bayesian network method for dynamic disassembly based on the SLIB disassembly diagram model [14]. Baazouzi et al. 2023 establish a multi-method self-configuration simulation model considering the battery design and the configuration and layout of the disassembly station [15]. Xiao et al. 2023 discussed the possibility of echelon utilization in future EV decommissioning battery recycling from two aspects: disassembly optimization and human–robot collaboration [16]. Tan et al. 2021 proposed a battery disassembly framework with improved disassembly efficiency as the core, using improved automated robotic arms and special tools to shorten the disassembly time [17]. Wu et al. 2023 proposed a method of power battery disassembly sequence planning based on knowledge graph representation [18]. However, most of these studies used one-sided and sequential disassembly methods without incorporating parallel disassembly techniques. In terms of characteristics, the parallel disassembly method is more in line with the requirements of the SLIB disassembly line. Therefore, it is necessary to perform parallel disassembly line planning for decommissioned power cell packs.

Smart manufacturing has emerged as a prominent focus of research within the industrial sector, and the advancement of sophisticated remanufacturing technologies is intricately linked to the concept of smart manufacturing. As a result, several researchers have applied automated disassembly techniques in the field of SLIB recycling and conducted related studies. Herrmann et al. 2012 took the implementability and necessity of the battery pack disassembly process as indicators and evaluated the automation potential of each step through the method of weight factors [19]. Hellmuth et al. 2021 proposed an automatic disassembly evaluation method and verified it with two disassembly examples of power batteries [20]. Chou et al. 2021 developed an automatic battery pack disassembly task planning system to realize the identification and positioning of SLIB components [21]. Zhou et al. 2020 presented a battery recycling path of disassembly for detection-grouping [22]. Compared to manual operation, the proposed method can significantly save operation time. Wegener et al. 2015 designed a robot-assisted disassembly system [23]. This system was able to efficiently perform simple repetitive tasks to reduce the operator’s workload. Kay et al. 2022 adopted a customized robotic arm to realize cutting and clamping operations in the battery disassembly process and realized the automatic/semi-automatic disassembly of the battery pack [24]. Li et al. 2019 designed an automatic disassembly device for flake batteries [25]. Schäfer et al. 2020 developed an automated system that could replace individual cells in a battery module while maintaining integrity [26]. To realize the automatic disassembly task planning of SLIB, Yu et al. 2022 proposed a disassembly task planning method based on ontology and partial failure rules [27]. Aimed at the uncertainties in the disassembly process of lithium batteries, Qu et al. 2023 developed a human–robot collaborative disassembly method based on digital twins [28]. It significantly increased the flexibility of the disassembly operation. Wang et al. 2023 proposed a SLIB disassembly planning method based on semantic ontology structure and knowledge graph [29]. Zhang et al. 2023 designed a knowledge-driven flexible human–robot hybrid disassembly line [30]. It could split the SLIB disassembly task layer by layer into knowledge-based primitive-level subtasks. Considering the risk and complex characteristics of battery components, Wu et al. 2022 proposed a human–robot collaborative disassembly planning model for SLIB [31].

The disassembly line balancing of waste power cells is a key link to maximize the combined benefits of disassembled SLIB. To a certain extent, it also reflects the level of human–robot collaborative disassembly techniques. The introduction of the human–robot cooperative mode and the consideration of workload brings a new approach to the problem of disassembly line balancing. The human–robot collaborative disassembly process and the corresponding disassembly safety are undetermined due to uncertainties in the internal and external characteristics factors of the used power cells. Disassembly planning relies heavily on the personnel’s work experience and related experience manuals, which makes it difficult to support the large-scale development of disassembled SLIB. Therefore, this work proposes an optimization method for human–robot collaborative disassembly line balancing, builds an optimization model for the route of the disassembly process, taking into account the workload, and uses a modified heuristic algorithm to find the optimal solution. It effectively optimizes the safe disassembly process of SLIB and comprehensively optimizes the economic benefits and smoothness of disassembly. In comparison with existing studies, the following contributions are made.

(1)

Establishing a human–robot collaborative disassembly line balancing model for SLIB considering time smoothness and workload smoothness.

(2)

Drawing an improved multi-objective optimization algorithm with Gaussian mutation and crowding distance.

(3)

Validating the effectiveness and sensitivity of models and algorithms by combining an SLIB disassembly practical case.

The following is the structure of this work. Section 2 constructs a human–robot collaborative disassembly model. In Section 3, an improved multi-target fruit fly algorithm is proposed. The case application’s results and analysis are provided in Section 4. Finally, the conclusion and future research directions are summarized.

2. Human–Robot Collaborative Disassembly Modeling

2.1. Problem Description

Human–robot collaborative disassembly (HRCD) means that the cooperative robot and the worker are assigned to the same station and complete some disassembly tasks together with the worker. Depending on the execution logic of the task in the workstation, the HRCD has different distribution patterns. Human–robot collaborative asynchronous parallel disassembly means that the employee and the cooperative robot are assigned to the same station, and the two can cooperate to complete the same task, or they can perform different tasks in parallel. When the task is executed in the human–robot cooperative mode, both workers and robots must be idle at the same time. When the task is completed, the worker and the robot perform other tasks in parallel. To further illustrate the disassembly mode allocation, Figure 2 shows a schematic of the human–machine hybrid parallel disassembly mode for six tasks. Tasks 3, 4, and 5 are assigned to the human–robot collaboration station. Task 3 involves a disassembly process that is primarily carried out by robots. In contrast, Task 4 encompasses a collaborative, human-centered approach to dismantling. Task 5 indicates that both humans and robots have decision-making abilities when executing tasks.

Simultaneously, the process of disassembly typically necessitates the expenditure of both physical and mental exertion on the part of the worker, coupled with a manual workload. Physical load refers to the amount of physical work that the human body undergoes per unit of time. A large workload increases the physical load. The mental load is a term that corresponds to the physical load and is used to describe the mental stress or information processing ability of a person while performing a task. Utilizing an ensemble model approach and computer simulation method, the physical energy consumption of the task is assessed based on both the static pose hazard level and dynamic energy consumption level. The mental load is assessed using the task rigidity method.

2.2. Symbolic Description

The symbolic variables used in the construction model are as follows:

i is the task number, $i\in \left(\mathrm{1,2},...,I\right)$.

j is the workstation number, $j\in \left(\mathrm{1,2},...,J\right)$.

k is the station direction number, $k\in \left(\mathrm{1,2}\right)$; when the direction is left, k = 1, to the right, k = 2.

r refers to the types of robots, $r\in \left(1,2,...,R\right)$.

p represents the decision variables selected for the disassembly mode, $p\in \left(1,2,3\right)$; when a task is performed by a human, p = 1; when a task is performed by a human, p = 2; otherwise, p = 3.

$t_i$ is the disassembly completion time of task i.

$t_i^s$ is the disassembly start time of task i.

$t_i^p$ is the time during which task i is executed by mode p.

$l_{pi}$ is the worker’s workload when task i is executed by mode p.

${\alpha }_1,{\alpha }_2$ are weighting coefficients.

$∆$ is a very large number.

${ST}_{jk}$ is the effective working time of the k side of workstation j.

${WL}_{jk}$ is the human workload on the k side of workstation j.

CT is the working beat.

$T_{max}$ is the maximum operating time on all workstations.

$L_{max}$ is the maximum workload of humans at all stations.

$C_r$ is the unit equipment wear cost of robot r.

$C_h$ represents the unit operating expenses of humans.

$C_{max}^R$ is the maximum budget for buying robots.

$N_{max}^h$ is the maximum number of humans involved in the operation.

$L_{max}^h$ is the maximum workload that a human can bear.

$x_{ijkp}$ if task i is assigned to the k side of workstation j and disassembly mode p is executed, $x_{ijkp}=1$; otherwise, $x_{ijkp}=0$.

$y_{rjk}$ if robot r is assigned to the k side of workstation j, $y_{rjk}=1$; otherwise, $y_{rjk}=0$.

$z_{jk}$ if a human is assigned to the k side of workstation j, $z_{jk}=1$; otherwise, $z_{jk}=0$.

$x_{i,o}$ if task i and task o are assigned to the same workstation, and i takes precedence, $x_{i,o}=1$; otherwise, $x_{i,o}=0$.

2.3. Mathematical Model Construction

Taking the total cost and the comprehensive smoothing index as the optimization objectives, this paper establishes the human–robot collaborative parallel disassembly line balancing model considering the workload. It is worth noting that the comprehensive smoothness index refers to the weighted sum of the worktime smoothness index and the workload smoothness index. To facilitate the construction of the mathematical model, we propose some hypotheses.

(1)

The structural information of the product to be dismantled is known and can be completely disassembled, and each task can be performed by an operator or robot.

(2)

The robot tool conversion time is ignored, and there is no fault in the execution of the equipment.

(3)

The disassembly time for the manual completion of the same task is the same, and the unit labor cost is the same.

(4)

The total disassembly time is determined by humans, robots, or both.

(5)

Under manual disassembly and the human–robot collaborative disassembly mode, the workload consumed by the task is fixed.

(6)

Each workstation can only be assigned a maximum of one employee and one robot.

(7)

Only one type of battery is disassembled.

The human–robot parallel collaborative disassembly line balancing optimization model considering the workload can be described as follows. In the human–robot parallel collaborative mode, the selection, and allocation of employees and robots are determined under the conditions of satisfying practical constraints, such as priority relationship constraints, time constraints, and workload threshold constraints, so as to achieve the optimization of the beat, cost, and comprehensive smoothness index.

$$\mathrm{m}\mathrm{i}\mathrm{n} f_1=CT$$

(1)

$$\mathrm{m}\mathrm{a}\mathrm{x} f_2={\sum }_{r=1}^R{\sum }_{j=1}^J{\sum }_{k=1}^2{C}_r\cdot y_{rjk}+{\sum }_{j=1}^J{\sum }_{k=1}^2{C}_h\cdot z_{jk}$$

(2)

$$\mathrm{m}\mathrm{i}\mathrm{n} f_3={\alpha }_1\cdot \sqrt{{\sum }_{j=1}^J{\sum }_{k=1}^2{\left(T_{max}-{ST}_{jk}\right)}^2/2\cdot J}+{\alpha }_2 \sqrt{{\sum }_{j=1}^J{\sum }_{k=1}^2{z_{jk}\left(L_{max}-{WL}_{jk}\right)}^2/{\sum }_{j=1}^J{\sum }_{k=1}^K{z}_{jk}}$$

(3)

where $f_1$ can minimize the disassembly line time; $f_2$ can minimize the total cost of the robot and human; and $f_3$ is the minimum comprehensive smoothness index. The first half represents the working time smoothness index, and the second half represents the workload smoothness index.

The constraint condition is as follows:

$${\sum }_{j=1}^J{\sum }_{k=1}^2{\sum }_{p=1}^3{x}_{ijkp}=1$$

(4)

$${\sum }_{j=1}^J{\sum }_{k=1}^2{\sum }_{p=1}^{3}J\cdot x_{r,o}\le {\sum }_{i=1}^I{\sum }_{k=1}^2{\sum }_{p=1}^{3}j\cdot x_{ijkp}$$

(5)

$$t_i^s+\mathsf{\Delta }\left(1-{\sum }_{p=1}^3{x}_{ijkp}\right)+\mathsf{\Delta }\cdot \left(1-{\sum }_{k=1}^2{\sum }_{p=1}^3{x}_{ojkp}\right)\ge t_o^s+t_o$$

(6)

$$t_o^s+\mathsf{\Delta }\left(1-{\sum }_{p=1}^3{x}_{ijkp}\right)+\mathsf{\Delta }\cdot \left(1-{\sum }_{p=1}^3{x}_{ojkp}\right)+\mathsf{\Delta }\left(1-v_{io}\right)\ge t_i^s+t_i$$

(7)

$$t_i^s+\mathsf{\Delta }\left(1-{\sum }_{p=1}^3{x}_{ijkp}\right)+\mathsf{\Delta }\cdot \left(1-{\sum }_{p=1}^3{x}_{ojkp}\right)+\mathsf{\Delta }\cdot v_{io}\ge t_o^s+t_o$$

(8)

$$t_i+\mathsf{\Delta }\left(1-x_{ijkp}\right)\ge t_{pi}$$

(9)

$$t_i^s+t_i\le CT, t_i^s\ge 0$$

(10)

$${\sum }_{r=1}^R{y}_{rjk}\le 1$$

(11)

$$\mathsf{\Delta }\cdot z_{jk}\ge {\sum }_{i=1}^I{x}_{ijk}+{\sum }_{i=1}^I{\sum }_{p=1}^3{x}_{ijkp}$$

(12)

$${\sum }_{j=1}^J{\sum }_{k=1}^2{z}_{jk}\le N_{max}^h$$

(13)

$${\sum }_{r=1}^R{\sum }_{j=1}^J{\sum }_{k=1}^2{C}_r\cdot y_{rjk}\le C_{max}^R$$

(14)

$${WL}_{jk}={\sum }_{i=1}^I{x}_{ijkp}\cdot l_{pi}$$

(15)

$$T_{max}\le L_{max}^h$$

(16)

Constraint (4) comprises task assignment constraints, and each task must be assigned to a certain side of the station and executed by a certain assembly mode. Constraint (5) requires that the priority of the task must be assigned to the front. Constraint (6) requires that the start time of task i must be greater than or equal to the end time of task o. Constraint (7) requires that when task i precedes task o at the same station, the start time of task o is no less than the end time of task i. Constraint (8) requires that when task o precedes task i at the same station, the start time of task i is no less than the end time of task o. Constraint (9) specifies the processing time of task i. Constraint (10) is a beat constraint, requiring that the completion time of all tasks must be within a given beat. Constraint (11) requires that only one robot can be assigned to each station. Constraint (12) determines the assignment of employees, and if task i performs disassembly mode p, one employee must be assigned to station j. Constraint (13) is the constraint on the number of employees. Constraint (14) is the purchase cost where the robot is required to be constrained. Constraint (15) is used to obtain the workload at each station. Constraint (16) is a constraint on the human workload threshold.

3. Improved Multi-Objective Fruit Fly Optimization Algorithm

The fruit fly optimization algorithm (FOA), as a heuristic algorithm for solving continuous nonlinear problems, has the advantages of fewer parameters, simple calculation, and a strong search ability. The core idea of FOA is to imitate the foraging process of fruit flies. Initially, fruit flies are attracted by the olfactory signals of high food concentration, utilizing their keen sense of smell to pinpoint the area. They then employ their remarkable visual acuity to locate the source of nourishment and proceed to navigate toward it with precision. Subsequently, they continue their quest for sustenance upon arrival at the designated location. Compared with traditional intelligent algorithms, FOA has an excellent global optimization mechanism. FOA tends to fall into local optimality during optimization, especially in the process of solving complex functions. Aiming at the defects of FOA, an improved multi-objective fruit fly optimization algorithm (IMFOA) is proposed by introducing chaos initialization and inertia weight updating strategies. The detailed steps are described as follows.

3.1. Encoding and Decoding

The problem of parallel disassembly sequence planning involves not only the constraint relationship between disassembly tasks but also the disassembly mode and the assignment of disassembly operators. Therefore, this work adopts multi-segment structure coding. In multi-segment structure coding, each drosophila individual is composed of a disassembly task segment, disassembly mode segment, and operator segment. For the disassembly task layer $v_1=\left\{x_1,x_2,..., x_n\right\}$, n indicates the total number of disassembly tasks. Each element represents a task to be disassembled. For the type of disassembly layer $v_2=\left\{y_1,y_2,..., y_n\right\}$, each element represents the disassembly mode corresponding to the disassembly task in $v_1$, that is, human disassembly, robot disassembly, and human–robot collaborative disassembly. For the operator layer $v_3=\left\{z_1,z_2,..., z_n\right\}$, each element represents the operator performing the corresponding disassembly task in segment $v_1$. When the corresponding disassembly mode is worker disassembly or robot disassembly, each element represents an operator, that is, a worker or a robot. When the corresponding disassembly mode is human–robot collaborative disassembly, each element represents two operators, that is, a worker and a robot. The coding process is shown in Figure 3, where three operator types perform the disassembly operation. For the type layer, type 1 pertains to the exclusive execution of tasks by humans. On the contrary, type 2 denotes solitary task performance by the robot. Type 3 signifies a collaborative effort between humans and robots in executing tasks. As for the operator layer, tasks {4, 7, 5} are dismantled by humans. Tasks {2, 6} are undertaken in the robot disassembly mode, with disassembly by a robot. Tasks {8, 1, 3} are in the human–robot collaborative disassembly mode, which is disassembly by humans and robots. It is worth noting that Task 8 is a human-centered HRCD with operation codes (1, 3). Task 1 is an HRCD centered on robots, with operation codes (2, 3). In addition, Task 3 is an equivalent HRCD with operation codes (3, 3).

According to the beat and workload threshold, the decoding adopts the mode of task collection. The specific decoding steps are as follows.

Step 1: Select a sequence based on the disassembly mode of the combined vector to determine the disassembly mode assigned to both sides of the current workstation.

Step 2: For unassigned tasks, select the disassembly mode with the shortest processing time to execute the task. If the disassembly time is the same, the disassembly mode with the least workload is selected. If the workload is the same, the manual disassembly mode is selected.

Step 3: Obtain the task set of the workstation according to the beat constraint and workload threshold constraint.

Step 4: If the workstation’s set of assignable tasks is empty and all tasks have been assigned, terminate the process. If the task set is empty but the task is not assigned, return to step 1. If the set of assignable tasks is not empty, proceed to Step 5.

Step 5: If none of the assignable task sets are empty, select the set with the largest remaining capacity and randomly select a set when the remaining capacity is equal.

Step 6: Determine whether the assignable task is executed at the earliest start time on the selected station.

Step 7: Select the task in front of the task sorting from the assignable task set, assign the task to the selected set, and update the task set capacity; proceed to Step 2.

3.2. Population Initialization

To avoid FOA falling into local optimal, the chaos mapping theory is introduced. The chaotic variable search has better overall coordination performance than random search and has stronger ergodic, regularity, and randomness to the initial value. Logistic mapping, as a typical chaotic system, has been widely used in the optimization of different fields. Logistic mapping is calculated as follows.

$$a_{i+1}=\lambda \times a_i\times \left(1-a_i\right)$$

(17)

where $a_i$ is the chaos variable, $a_i\in \left[0,1\right]$; $\lambda$ is the control parameter, $\lambda \in \left[0,4\right]$.

When $\lambda =4$,${ 0\le a}_i\le 1$, logistic is in a completely chaotic state. For the problem variable $x_i$, if it satisfies $x_i\in \left[{\delta }_i, {\phi }_i\right]$, the mapping between the chaos variables can be expressed as follows.

$$a_i=\left(x_i-{\delta }_i\right)/\left({\phi }_i-{\delta }_i\right)$$

(18)

$$x_i={\delta }_i+a_i·\left({\phi }_i-{\delta }_i\right)$$

(19)

where ${\delta }_i$ and ${\phi }_i$ are the lower and upper limits of the variable, respectively.

According to the above analysis, the initial position of individual flies can be obtained by introducing Logistic mapping into the initialization of the drosophila population $\left(x\_axis,y\_axis\right)$.

$$\left\{\begin{array}{c}x\_axis=\mathrm{r}\mathrm{a}\mathrm{n}\mathrm{d}\left(\mathrm{d}\mathrm{o}\mathrm{m}\mathrm{a}\mathrm{i}\mathrm{n} \mathrm{o}\mathrm{f} \mathrm{d}\mathrm{e}\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{i}\mathrm{t}\mathrm{i}\mathrm{o}\mathrm{n}\right)\\ y\_axis=\mathrm{r}\mathrm{a}\mathrm{n}\mathrm{d}\left(\mathrm{d}\mathrm{o}\mathrm{m}\mathrm{a}\mathrm{i}\mathrm{n} \mathrm{o}\mathrm{f} \mathrm{d}\mathrm{e}\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{i}\mathrm{t}\mathrm{i}\mathrm{o}\mathrm{n}\right)\end{array}\right.$$

(20)

$$\left\{\begin{array}{c}{x}_i=\lambda \times x\_axis\times \left(1-x\_axis\right)\\ y_i=\lambda \times y\_axis\times \left(1-y\_axis\right)\end{array}\right.$$

(21)

where the $\mathrm{r}\mathrm{a}\mathrm{n}\mathrm{d}\left(\mathrm{d}\mathrm{o}\mathrm{m}\mathrm{a}\mathrm{i}\mathrm{n} \mathrm{o}\mathrm{f} \mathrm{d}\mathrm{e}\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{i}\mathrm{t}\mathrm{i}\mathrm{o}\mathrm{n}\right)$ is a random number in the domain$\left[0,1\right]$.

3.3. Olfactory Search

In the process of olfactory search, to improve the diversity of solutions and the convergence speed of the algorithm, the Gaussian variation operation is introduced to generate multiple random directions and random search distances for individual positions. The variation factors are adjusted according to the relevance of each fruit fly population variation and the weight of the optimization goal until the new variable meets its defined range. The Gaussian variation process is as follows.

$$y_k^{\prime }=y_k+s_g·N_k\left(\mathrm{0,1}\right)$$

(22)

$$s_g=random\left(+,-\right)\sqrt{2\mathit{ln}\left({\omega }_g\sqrt{2\pi }\right)}$$

(23)

where $y_k$ and $y_k^{\prime }$ represent the values of the first decision variable before and after variation, respectively. $s_g$ is the variation length; $g$ is the number of Gaussian variations, $1\le g\le TG$; ${\omega }_g$ is a random number in (0, fgau(0)); and $N_k\left(\mathrm{0,1}\right)$ represents a Gaussian random number with a mean of 0 and a variance of 1 when the kth decision variable is variable.

3.4. Visual Operation

In order to avoid falling into local optimality in visual operation, the Pareto non-dominated sorting method is adopted. The drosophila melanoglypha constantly updates its position according to the optimal domain, and the neighborhood solutions are Pareto screened to obtain SN optimal solutions. We introduce the dissimilarity to represent the distance between the Pareto solution and its parent. The optimal solution replaces the parent individual, updates the individual position, and finally determines the optimal solution. The formula for calculating heterogeneity is as follows.

$$D\left(X_i,X_j\right)={\sum }_{k=1}^{n}sgn\left|x_{ik}-x_{jk}\right|$$

(24)

where $n$ represents the spatial coordinate dimension; $x_{ik}$ is the dimensional component of the space coordinate; and $x_{jk}$ represents the k-dimensional component of the spatial coordinate $X_j$.

3.5. External File Update

The non-inferior solutions obtained by the drosophila optimization algorithm are constantly retained and updated to the external file. The solution with the larger crowding distance is used as the new individual of the initial population in the next iteration; that is, the larger the crowding distance, the more uniform the distribution of the solution. The crowding distance of the boundary individual is defined as $l_1^k=l_D^k=1$, and the formula for calculating the non-edge crowding distance of the multi-objective non-inferior solution $X$ is shown as follows.

$$L\left(X\right)={\sum }_{k=1}^m\left(\frac{f_k^{i+1}-f_k^{i-1}}{f_k^{max}-f_k^{min}}\right)$$

(25)

where $i$ is the number of solutions, $i=2,...,D-1$; $m$ is the number of objective functions; and $f_k^{max}$ and $f_k^{min}$ represent the maximum and minimum values of the first objective function, respectively.

4. Computational Experiment and Analysis

4.1. Experimental Design

As an important part of resource recovery for new energy vehicles, battery packs have the characteristics of high disassembly and reuse rates and great economic value. Therefore, an example of Mercedes-Benz EQS580 battery pack disassembly is used in this study to verify the model and algorithm, and the battery pack structure is shown in Figure 4.

Table 1. Relevant information of disassembly experiment.

No.	Name of Part	Direction	Priority Task	Operate Time (s)	Quantity
1	Low-voltage harness	+Z	1	13	1
2	Housing cover	+Z	2, 3	62	1
3	Upper cover seal	−Z	1, 2	8	1
4	High-voltage harness	−X	1, 5, 6	16	1
5	High-voltage connectors	+X	4, 6, 7	36	1
6	Electronics	−Z	5, 7	28	4
7	Top plate	+Z	4, 6, 9	84	1
8	Cells	+Z	8	425	216
9	Cell models	−Z	7, 10	246	12
10	Housing	−Z	9, 7, 12	45	1
11	Battery cooling	+X	10, 12	22	1
12	Underbody cover	+Y	10	18	1
13	Lower cover seal	+Z	11, 12	10	1

shows the relevant information on the disassembly experiment of the EQS580 battery pack, including the disassembly direction, the immediate task, the task completion time, and the number of components. The disassembly time adopts stopwatch timing and the Three Sigma criterion. The average of multiple disassembly is used as the time for task completion. To determine the disassembly direction, the geometric center of the battery pack is taken as the coordinate origin, and the three-dimensional coordinate system is used to disassemble the battery pack. The algorithm running environment is MATLAB 2016, Windows 10, Intel Core i5, 2.8 GHz CPU, 8 GB RAM. FOA algorithm parameters are set with the maximum number of iterations at 100, a population size of 100, and all experiments are run independently 10 times.

4.2. Result Analysis

Table 2. Objective function value corresponding to Pareto solution.

Solution	$\mathbf{Disassembly} \mathbf{Sequence}$	${\mathit{f}}_1$ (s)	${\mathit{f}}_2$ (CNY)	${\mathit{f}}_2$
1	Human: 2→3→5→7→9 Robot: 1→4→6 Human–robot: 8→10→13→12→11	1058	33	124
2	Human: 2→3→5→7 Robot: 1→4→6→11 Human–robot: 8→9→10→13→12	1074	32	116
3	Human: 1→2→3→4→12 Robot: 5→6→11 Human–robot: 7→8→9→10→13	1182	42	118
4	Human: 3→4→8→9 Robot: 1→2→5→6→11 Human–robot: 7→10→13→12	1013	31	132
5	Human: 1→2→3→4→8→9 Robot: 5→6→11 Human–robot: 7→10→13→12	1006	34	107
6	Human: 3→4→8→9→13 Robot: 1→2→5→6 Human–robot: 7→10→11→12	1122	38	103
7	Human: 3→4→7→8→10 Robot: 1→2→5→6 Human–robot: 9→11→13→12	1046	41	123
8	Human: 2→3→10→13→12 Robot: 1→4→6 Human–robot: 5→7→8→9→11	1025	38	126

shows a group of Pareto solutions obtained by IMFOA. In the table, the first column is the number of the disassembly scheme, and the second column is the disassembly sequence of the scheme. The remaining three columns are the three objective function values of disassembly time, disassembly cost, and comprehensive smoothness. As can be seen from Table 2, the value of $f_1$ ranges from 1006 to 1182 s. The maximum value of $f_2$ is scheme 3 (42 Chinese Yuan), and the minimum value is scheme 4 (31 Chinese Yuan). $f_3$ ranges from 103 to 132. For the time index, scheme 5 is the best, and the corresponding disassembly sequence is human (1, 2, 3, 4, 8, 9), robot (5, 6, 11), and human–robot collaboration (7, 10, 13, 12). In terms of economy, scheme 4 is the best, with the corresponding disassembly sequence of human (3, 4, 8, 9), robot (1, 2, 5, 6, 11), and human–robot collaboration (7, 10, 13, 12). For the smoothness index, the optimal scheme is scheme 6, with the corresponding disassembly sequence of human (3, 4, 8, 9, 13), robot (1, 2, 5, 6), and human–robot collaboration (7, 10, 11, 12).

4.3. Algorithm Comparative Analysis

In order to better illustrate the robustness of the IMFOA algorithm, it is necessary to compare the proposed algorithm with traditional FOA and NSGA-II. The framework of NSGA-II and FOA is very similar to IMFOA. In addition, NSGA-II is widely used as a benchmark algorithm in many studies, providing benchmark results when evaluating new multi-objective optimization algorithms. All parameter settings are the same for all algorithms to allow for reasonable comparisons. To further demonstrate the superior performance of the improved algorithm, NSGA-II and FOA are executed under the same conditions. These two algorithms adopted the same coding and decoding rules, modeled by objective function and the screened Pareto disassembly scheme. To reflect the performance of different algorithms, three performance indicators were used to compare the running results of the algorithms. Since the algorithm terminates with a fixed number of iterations, the shorter the execution time (ET) on each indicator, the better the performance and efficiency of the algorithm. In addition, two evaluation indexes, including the convergence measure (CM) and hypervolume (HV) of the Pareto optimal solution set, are introduced to provide the combined information on the convergence and diversity of the obtained solution set. For solution sets with the same convergence, CM preferentially selects individuals that are evenly distributed along Pareto fronts. The smaller the value of CM, the better the convergence of the solution. HV reflects the volume of the hypercube surrounded by all reference points in the Pareto solution set in the target space. The larger the index value, the more uniform the Pareto solution distribution.

Table 3. Calculation results of IMFOA, FOA, and NSGA-II.

PopSize	Gmax	Pareto Solutions
		NSGA-II			FOA			IMFOA
		f1	f2	f3	f1	f2	f3	f1	f2	f3
30	30	1093 1134 1145 1179	36 42 44 46	138 134 129 135	1175 1167 1057 1088 1231	42 39 38 43 44	136 138 125 117 129	1058 1038 1054 1078 1124 1175	34 39 47 41 37 36	128 113 125 124 112 128
50	50	1138 1146 1178 1096 1112	36 39 40 39 41	129 133 127 130 126	1142 1137 1074 1073 1182 1164	41 40 36 39 42 39	118 124 113 128 125 119	1039 1064 1063 1143 1034 1155 1073 1140	39 42 38 38 44 37 36 40	133 118 124 117 121 120 119 125
100	100	1156 1168 1175 1203 1197	744 721 729 694 713	128 119 124 121 122	1123 1168 1185 1184 1075 1174 1173	41 34 35 36 42 39 38	127 132 114 121 126 129 127	1058 1074 1082 1013 1006 1022 1046 1025	33 32 42 31 34 38 41 38	124 116 118 132 107 103 123 126

shows the results of the experiment, and it can be seen that non-dominated solutions obtained from IMFOA are more numerous and more evenly distributed than NSGA-II and FOA. Generally speaking, as the population size increases, there is a progressive proliferation of Pareto solutions. This phenomenon indicates a positive correlation between population size and the abundance of Pareto solutions. In terms of the number of solutions, the results of IMFOA are superior to those of NSGA-II and MDFOA. The above comparison results show that IMFOA proposed in this study is superior to FOA and NSGA-II in solving the constructed mathematical model. At the same time, IMFOA is more efficient and robust in solving the constructed model.

For multi-objective optimization problems, it is necessary to evaluate some comprehensive indexes based on the Pareto algorithm to further test the performance of MSFOA. In this paper, six calculation examples in the literature [8,32] were selected for comparative testing, and the results are shown in

Table 4. CM, ET, and HV results of the three algorithms under different examples.

Case	Algorithm	ET(s)	CM (Mean ± Std Dev)	HV (Mean ± Std Dev)
P1	NSGA-II	4.2	4.32 ± 4.62 × 10−4	2.36 ± 4.37 × 10−5
	FOA	3.4	1.32 ± 2.21 × 10−4	1.93 ± 3.46 × 10−4
	IMFOA	2.8	1.36 ± 1.51 × 10−4	2.82 ± 3.14 × 10−4
P2	NSGA-II	10.2	1.25 ± 3.34 × 10−3	1.76 ± 6.53 × 10−4
	FOA	6.3	2.57 ± 1.84 × 10−3	2.25 ± 3.44 × 10−4
	IMFOA	3.5	1.26 ± 4.45 × 10−4	2.32 ± 4.16 × 10−3
P3	NSGA-II	13.6	2.36 ± 1.37 × 10−2	1.72 ± 3.15 × 10−4
	FOA	9.5	2.53 ± 3.34 × 10−2	2.36 ± 1.34 × 10−4
	IMFOA	5.7	1.55 ± 3.38 × 10−3	1.52 ± 2.75 × 10−4
P4	NSGA-II	20.1	3.24 ± 2.51 × 10−3	1.21 ± 3.32 × 10−4
	FOA	12.3	1.63 ± 2.44 × 10−3	3.53 ± 2.35 × 10−5
	IMFOA	8.5	1.23 ± 3.68 × 10−4	1.48 ± 3.25 × 10−4
P5	NSGA-II	18.3	3.43 ± 1.83 × 10−3	1.33 ± 2.04 × 10−5
	FOA	10.4	2.21 ± 3.86 × 10−3	2.02 ± 4.11 × 10−5
	IMFOA	7.6	1.55 ± 2.17 × 10−4	2.53 ± 3.22 × 10−4
P6	NSGA-II	22.7	4.01 ± 2.35 × 10−3	1.75 ± 4.30 × 10−4
	FOA	13.6	3.39 ± 1.27 × 10−3	2.53 ± 1.14 × 10−4
	IMFOA	8.5	2.66 ± 3.26 × 10−4	2.28 ± 3.41 × 10−3

. As can be seen from Table 4, the ET value of IMFOA is significantly lower than that of the other two algorithms. The more complex the example, the more clear the advantage of the IMFOA algorithm in running time. For the CM index, NSGA-II and FOA have little difference, but they are larger than IMFOA. FOA and IMFOA have higher HV values than NSGA-II, but IMFOA performs better overall.

5. Conclusions

With the acceleration of the intelligent process and the deepening of the concept of green manufacturing, human–robot collaborative disassembly has become an important means to realize the efficient operation of the disassembly line. Some advancements have been achieved in the exploration of collaborative human–robot disassembly line balancing. However, researchers have paid little attention to the effect of uncertain disassembly on the workload. In this paper, a mathematical model of the human–robot collaborative disassembly line balancing model is established. To obtain a more optimized algorithm, an improved multi-objective fruit fly optimization algorithm is proposed. By comparing some traditional classical algorithms and newly developed meta-heuristic algorithms, the effectiveness of the proposed algorithm is proved. This study provides a new perspective for intelligent disassembly, improves the compatibility of human–robot collaborative work and human factor safety, and promotes the development of human-centered intelligent manufacturing. In this work, the balancing problem of the human–robot collaborative disassembly line considering human-centered time variations is studied. However, there are many directions to expand in the future. For example, the robot avoids the operator in the process of task grasping. The integration of technologies, such as digital twins and blockchain, is also a direct direction that can be explored in depth. For the model presented in this paper, we encourage researchers to try to develop precise algorithms to solve it.