1 Introduction

Assembly lines have been widely utilized in modern industries to assemble products, and such systems help to improve the productivity in the production floor Battaïa and Dolgui (2013). The traditional assembly line consists of a set of connected stations and one worker is allocated to each station to perform the assembly tasks. On the assembly line, the assembly tasks are divided into several task sets and the workers on stations operate the corresponding task set from the first station to the last station in a sequence. The assignment of the tasks determines the line balance and based on this assembly line balancing problem (ALBP) is formulated to optimize the task assignment with one or several optimization criteria. There are two important constraints in the ALBP: precedence constraint and cycle time constraint. Precedence constraint demands that the predecessors of one task must be completed before operating this task; cycle time constraint demands that all the tasks on stations must be completed within the pre-determined cycle time (CT).

U-shaped assembly line balancing problem (UALBP) has received increasing attentions since it offers higher flexibility than a normal assembly line. U-shaped assembly line could be divided into two sub-lines: the sub-line on the entrance side and the sub-line on the exit side. On the traditional U-shaped assembly line, human workers are allocated to the stations to operate all the tasks. The worker on one station first operates the tasks on the entrance side of this station, later moves to the exit side and operates the tasks on the exit side of this station, finally comes back to the entrance side. Faced with increased labor cost and concerns of efficiency, one of the trends that is being adopted in industries is the replacement of the human workers with robots. U-shaped assembly line with robots to complete all the tasks is referred to as robotic U-shaped assembly line (Li et al. 2019a). This type of line will help to reduced labor cost, higher accuracy and higher stability. The robotic U-shaped assembly line might be criticized for being expensive and lacks flexibility, another trend that is seen recently is to utilize collaborative robots, referred to as cobots, to operate some tasks solely or assist the human workers to operate some tasks (Weckenborg et al. 2020). The utilization of cobots achieves the human–robot collaboration, which combines the advantages of the human workers and cobots to help and improve the productivity of the system. Specifically, human–robot collaboration inherits the higher flexibility, higher adaptability and high decision-making skills of human workers and the high strength, endurance and accuracy of cobots (Weckenborg et al. 2020). When designing the U-shaped assembly line with cobots, U-shaped assembly line balancing problem with collaborative robots (CRUALBP) should be studied to optimize the task assignment and worker and robot allocation.

Figure 1 depicts the U-shaped assembly line with workers and cobots. In the U-shaped assembly line with workers, each station is allocated with one worker. The worker on one station first operates the tasks on the entrance side, later operates the tasks on the exit side and finally comes back to the entrance side to operate the corresponding tasks. In the U-shaped assembly line with robots, robots replace the workers completely to assemble the products and all the tasks are operated by robots. In the U-shaped assembly line with workers and cobots, one station might be allocated with one worker, one cobot or a pair of worker and cobot. On the station with one worker and one cobot, cobot can operate one task separately or assist the worker to operate the common tasks collaboratively.

Fig. 1
figure 1

U-shaped assembly line with workers and cobots

Although there are some reported literature on robotic U-shaped assembly line balancing problem (RUALBP), there is no study tackling the CRUALBP. Hence, this study presents the first attempt to solve the considered CRUALBP, where different types of cobots with different cost should be selected and allocated to the stations within the maximum budget allowed. This study presents three main contributions as follows. 1) This study considers the CRUALBP where several cobots are selected within the budget constraint and this done for the first time in the literature. 2) Three mixed-integer linear programming (MILP) models are formulated to optimize the cycle time with the given station number. The developed MILP models can solve the small-sized instances optimally utilizing the CPLEX solver. 3) Improved artificial bee colony algorithm (IABC) and improved migrating bird optimization algorithm (IMBO) are developed to solve the large-sized instances within the acceptable computation time. The implemented algorithms utilize new encoding scheme and decoding procedure, and propose several improvements to achieve the proper balance between exploitation and exploration. The comprehensive comparative study is carried out to evaluate the performances of these algorithms, and the proposed algorithms obtain competing performance in comparison with other 12 implemented algorithms.

The structure of this paper is presented as follows. Section 2 provides the literature review and Sect. 3 presents the problem description and formulation. Subsequently, Sect. 4 describes the implemented metaheuristic methodologies and Sect. 5 illustrates a small-sized example. Section 6 carries out the comparative study, where the formulated model and algorithms are evaluated. Finally, Sect. 7 provides the conclusions and directions for future research.

2 Literature review

Since the first research paper published in 1955, several variants of ALBPs have been studied by many researchers (Battaïa and Dolgui 2013; Eghtesadifard et al. 2020). Based on the assembly layout, the assembly lines can be divided into straight line, U-shaped line, two-sided line, multi-manned line and parallel line. Based on the process alternative, there are three types: human worker, robot (or cobot) and worker and robot in collaboration. As there is no study considering all the characteristics of the problem in this study, this section mainly reviews the recent studies on the UALBP with human workers, RUALBP and assembly line balancing problem with collaborative robots (CRALBP).

In the basic UALBP with human workers, all the workers are considered to have the same skillset and only task assignment needs to be optimized. Miltenburg and Wijngaard (1994) presented the pioneering work on the UALBP and they develop the dynamic programming procedure to minimize the station number. Since then, many exact methods, heuristics and metaheuristics are developed. Regarding the exact methods, Scholl and Klein (1999) presented a branch and bound algorithm known as ULINO. Fattahi and Turkay (2015) formulated a corrected model to handle the precedence relation, and later Li et al. (2017) presented a new MILP model to solve the UALBP optimally. Li et al. (2018b) developed the branch, bound and remember algorithm to minimize the station number, which outperforms all the exact methods and is the best exact method. For the heuristic and metaheuristic methodologies, they include simulated annealing algorithm (Erel et al. 2001), genetic algorithm (Hwang et al. 2008), ant colony algorithms (Baykasoglu and Dereli 2009; Li et al. 2017; Sabuncuoglu et al. 2009), hybrid improvement heuristic approach (Özcan and Toklu 2009), critical path method (Avikal et al. 2013) and others, to cite just a few. Recently, Li et al. (2020) proposed an enhanced beam search algorithm, which outperforms the published algorithms and updates the upper bounds for many instance. Among these methodologies, branch, bound and remember algorithm is the best exact method (Li et al. 2018b) and the enhance beam search is state-of-the-art method for the UALBP (Li et al. 2020). There are also some studies considering the human workers with different skills, Oksuz et al. (2017) formulated the U-shaped assembly line worker assignment and balancing problem and they also develop the artificial bee colony algorithm and genetic algorithm to tackle it in short computational time. Afterward, Zhang et al. (2019a) presented an enhanced migrating birds optimization algorithm, which outperforms several other algorithms. Zhang et al. (2020) developed restarted iterated Pareto greedy algorithm to optimize the ergonomic risk and cycle time simultaneously, and this algorithm outperforms several other multi-objective algorithms in the comparative study.

For the studies on the RUALBP, Nilakantan and Ponnambalam (2016) formulated this problem to minimize the cycle time for the first time and develop the particle swarm optimization to solve this problem. Li et al. (2019a) presented two new MILP models with the objective of minimizing the cycle time and develop an improved migrating birds optimization algorithm, which outperforms seven other implemented algorithms. Zhang et al. (2019b) considered the energy consumption and the cycle time in the RUALBP. They build the multi-objective model and develop the Pareto artificial bee colony algorithm to achieve a set of non-dominated solutions. Zhang et al. (2019c) tackled the multi-objective RUALBP to minimize the carbon emission, noise emission and cycle time. They formulated this problem with a multi-objective model and developed the hybrid Pareto grey wolf optimization to solve this problem. The proposed algorithm outperforms five multi-objective algorithms in the comparative study.

In the literature, there are several different types of human–robot collaboration in the assembly lines and they are summarized as follows. Çil et al. (2018) formulated the semi-robotic ALBP and present a case study. Weckenborg and Spengler (2019) formulated the mathematical model to optimize the cost in the assembly line with collaborative robots taking into account ergonomics. Samouei and Ashayeri (2019) considered the semi-automated operations in the mixed-model ALBP and they formulate two mathematical models to solve this problem. Dalle Mura and Dini (2019) presented a genetic algorithm to tackle the CRALBP. Weckenborg et al. (2020) formulates the CRALBP to minimize the cycle time, which is capable of solving the small-size instance optimally. They also develop a hybrid algorithm to tackle both the small-size and large-size instances. Yaphiar et al. (2020) presented a new MILP model for the mixed-model CRALBP. More recently, Çil et al. (2020) formulated the mixed-model ALBP with physical human–robot collaboration, which solves the small-size instance optimally. They also develop the improved bee algorithm and artificial bee colony algorithm to tackle the large-size instance efficiently. Rabbani et al. (2020) considered the human–robot collaboration in mixed-model four-sided ALBP. They formulated the MILP model and developed the particle swarm optimization algorithm to solve the considered problem. More recently, Li et al. (2021) formulated the model to minimize the cycle time and cobot cost. They also developed the Pareto migrating bird optimization algorithm to optimize two objectives simultaneously. Koltai et al. (2021) analyzed the task assignment and cycle times when robots are added with mathematical models. Nourmohammadi et al. (2022) optimized the cycle time and the number of cobots and workers in the CRALBP, and they provided the mathematicl model and simulated annealing algorithm to solve the CRALBP. Weckenborg et al. (2022) considered the ergonomic and economic objectives with collaborative robots and exoskeletons utilizing one mixed-integer programming model. The literature on the CRALBP are summarized in Table 1.

Table 1 Summary of literature on the CRALBP

From the literature review presented, the problem setting is different from that in the related and published studies and there is no study considering the CRUALBP. To fill this gap, this study formulates this problem for the first time and proposes two recent algorithms to solve this problem in reasonable computation time.

3 Problem description and formulation

This section first provides the problem description in Sect. 3.1 and later presents the MILP model of this problem in Sect. 3.2.

3.1 Problem description

For the human–robot collaboration, there are five scenarios in the assembly lines (Hashemi-Petroodi et al. 2020; Krüger et al. 2009) as follows. And all of them have real applications in the diverse industrial context. 1) Shared workspace, separate cooperation and sequential operations: A worker and a cobot are allocated to the same station but not the common workspace (e.g., the left side and the right side of one station), and they operate the tasks separately and they can only operate the tasks in sequence. 2) Shared workspace, separate cooperation and simultaneous operations: A worker and a cobot are allocated to the same station but not the common workspace, and they operate the tasks separately and they can operate the different tasks in parallel simultaneously. 3) Shared workspace, collaborative cooperation and sequential operations: A worker and a cobot are allocated to the same station but not the common workspace, and they can perform a common task in collaboration and they can only operate the tasks in sequence. 4) Shared workspace, collaborative cooperation and simultaneous operations: A worker and a cobot are allocated to the same station but not the common workspace, and they can perform a common task in collaboration and they can operate the different tasks in parallel simultaneously. 5) Common workspace, collaborative cooperation and sequential operations: A worker and a cobot are allocated to the same station and the same workspace, and they can perform a common task in collaboration and they can only operate the tasks in sequence.

For the considered CRUALBP problem, a worker and a cobot can be allocated to the same station and they can complete the tasks solely or in collaboration, but they cannot operate the different tasks in parallel simultaneously. Namely, this study takes the first scenario, the third scenario and the fifth scenario into account. The applications of second scenario and the fourth scenario might refer to Weckenborg et al. (2020). In the CRUALBP, one station might be allocated with a worker, a cobot and a pair of worker and cobot. If one worker and robot are allocated to the same station, one task on this station is operated by one of the three process alternatives: operated by a human worker, operated by cobot and operated by the worker and cobot in collaboration. Hence, the CRUALBP consists of three interrelated sub-problems: 1) task assignment to station; 3) the allocation of the workers and cobots; 3) the selection of the process alternative. And there are three main constraints needed to be considered when solving the CRUALBP: precedence constraint, cycle time constraint and budget constraint. Precedence constraint demands that one task can be assigned when all its predecessors or successors have been assigned. Cycle time constraint demands that all the task must be completed within the cycle time. Budget constraint requires that the costs of buying cobots must be equal to or less than the maximum budget provided in advance.

Figure 2 provides an example task assignment and worker and cobot allocation on U-shaped assembly line. As you can see, there are four stations utilized and 11 tasks are allocated to the entrance side and the exit side of this line. And one worker, one worker and one cobot, one worker and one worker are allocated to station 1, station 2, station 3 and station 4, respectively. For the task assignment on station 2, it is observed that task 2 is allocated to the entrance side and is operated by the worker and cobot in collaboration, task 3 is allocated to the entrance side and is operated by the worker and cobot in collaboration, task 7 is allocated to the exit side and is operated by the worker and finally task 9 is allocated to the exit side and is operated by the worker and cobot in collaboration.

Fig. 2
figure 2

Task assignment and worker and cobot allocation on U-shaped assembly line

The main assumptions of the considered CRUALBP are provided at first based on Weckenborg et al. (2020).

  1. (1)

    One type of product is produced on the U-shaped assembly line.

  2. (2)

    The precedence relations between tasks and the operation times of tasks by process alternatives are known and fixed.

  3. (3)

    The operation times of tasks are determined by the process alternative and one task might need different time to be completed when being operated by different process alternatives

  4. (4)

    Only one type of worker is considered (the workers are assumed to have the same skills).

  5. (5)

    There are different types of cobots available and they might have different skills and different purchasing cost.

  6. (6)

    The total cost of cobots must be within the maximum budget allowed to buy the cobots.

  7. (7)

    It is not considered that a worker and cobot operate different tasks at the same time on the same station (the second scenario and the fourth scenario are not considered).

  8. (8)

    The setup times, loading and unloading time and maintenance operation are not considered.

3.2 Model formulation

Based on Li et al. (2017), this section formulates three mathematical models: Model 1, Model 2 and Model 3 for the considered problem. For clarification, the indices, parameters and decision variables in these models are first explained as follows. Among the indices and parameters, the \(p\) is utilized to indicate the process alternative. For one task, there are \(2\cdot nr+1\) process alternatives available. Here, this study sets that \(p=1\) when one human worker operates one task, \(p=1+r\) when cobot \(r\) operates one task solely, and \(p=\mathrm{nr}+1+r\) when a worker and cobot \(r\) operate one task collaboratively.

Indices and parameters:

\(nt\)

Number of tasks.

\(ns\)

Number of available stations.

\(nr\)

Types of available cobots.

i, j

Task index; \(i,j\in \left\{\mathrm{1,2},\ldots ,nt\right\}.\)

k, l

Station index; \(k,l\in \left\{\mathrm{1,2},\ldots ,ns\right\}\).

r

Cobot index; \(r\in \left\{\mathrm{1,2},\ldots ,nr\right\}\).

\(p\)

Process alternative index, \(p\in \left\{\mathrm{1,2},\ldots ,nr+1,nr+2,\ldots ,2\cdot nr+1\right\}\).

\(I\)

Set of tasks, \(I=\left\{\mathrm{1,2},\ldots ,nt\right\}\).

\(K\)

Set of stations, \(K=\left\{\mathrm{1,2},\ldots ,ns\right\}\).

\(R\)

Set of cobot types, \(R=\left\{\mathrm{1,2},\ldots ,nr\right\}\).

\(P\)

Set of process alternatives, \(P=\left\{\mathrm{1,2},\ldots ,nr+1,nr+2,\ldots ,2\cdot nr+1\right\}\).

\({t}_{ip}\)

Operation time of task \(i\) by processing alternative \(p\).

\( \wp \)

Set of pairs of tasks (\(\wp =\{\cdots ,r,\cdots \}\)), where \(r=\left(i,j\right)\) is one pair of tasks in which task \(i\) is the immediate predecessor of task \(j.\)

\({c}_{r}\)

Purchasing cost of cobot \(r\).

\(\psi \)

A very large positive number.

\(MB\)

Maximum budget allowed to buy the cobots.

\(MW\)

Maximum number of available workers.

Decision variables:

\(CT\)

Cycle time.

\({a}_{ikp}\)

1, if task \(i\) is operated by process alternative \(p\) on station \(k\); 0, otherwise. (Only utilized in Model 1)

\({u}_{i}\)

1, if task \(i\) is allocated to the entrance side of one station; 0, otherwise. (Only utilized in Model 1)

\({x}_{ikp}\)

1, if task \(i\) is operated by process alternative \(p\) on the entrance side of station \(k\); 0, otherwise. (Only utilized in Model 2 and Model 3)

\({y}_{ikp}\)

1, if task \(i\) is operated by process alternative \(p\) on the exist side of station \(k\); 0, otherwise. (Only utilized in Model 2 and Model 3)

\({w}_{rk}\)

1, if cobot \(r\) is allocated to station \(k\); 0, otherwise.

\({v}_{k}\)

1, if a human worker is assigned to station \(k\); 0, otherwise.

The main difference between the CRALBP and CRUALBP is the precedence constraint. On the straight line, one task is assignable when all its predecessors have been allocated; on the U-shaped line, one task is assignable to the entrance side when all its predecessors have been allocated or the exist side when all its successors have been allocated. On the basis of the model in Fattahi and Turkay (2015), the first mathematical model, referred to as Model 1, is formulated utilizing expression (1) to expression(11), where the decision variable \({a}_{ikp}\) is utilized to determine the allocated station that one task is allocated to and decision variable \({u}_{i}\) is proposed to determine the side (entrance side or exit side) that one task is allocated to.

Expression (1) is objective function to minimize the cycle time. Constraint (2) is the task occurrence constraint and it demands that each task must be allocated to one station and operated by one process alternative available. Constraint (3) and constraint (4) deal with the precedence constraint. Constraint (3) indicates that the predecessor \(i\) of task \(j\) must be allocated to the former or same station when both of them are allocated to the entrance side; Constraint (4) indicates that the predecessor \(i\) of task \(j\) must be allocated to the latter or same station when both of them are allocated to the exit side. Constraint (5) deals with the cycle time constraint and it demands that the total operation time on one station must be less than or equal to the cycle time. Constraint (6) and constraint (7) determine the cobot allocation. Constraint (6) requires that cobot r must be allocated to station k when there is one task on this station is operated by cobot r solely and a worker and cobot \(r\) in collaboration. Constraint (7) demands that at most one cobot is allocated to one station and it is not allowed for one station to have more than one cobot. Constraint (8) is the budget constraint, indicating the total cost of buying the cobots must be less than or equal to the maximum budget. Constraint (9) determines the worker allocation and it requires that a worker must be allocated to station k when there is one task on this station is operated by a worker solely and a worker and one cobot in collaboration. Constraint (10) restricts that the number of workers is less than or equal to the maximum number of available workers. Constraint (11) restricts that \({a}_{ikp},{u}_{i},{w}_{rk}\) and \({v}_{k}\) are the 0–1 variables.

$$\mathrm{Min} \,\,\mathrm{CT}$$
(1)
$$\sum_{k\in K}\sum_{p\in P}{a}_{ikp}=1 \, \, \,\,\forall i\in I$$
(2)
$$\sum_{k\in K}\sum_{p\in P}{k\cdot a}_{ikp}-\sum_{l\in K}\sum_{p\in P}{l\cdot a}_{jlp}\le ns\cdot \left(1+{u}_{i}-2{u}_{j}\right) \, \,\,\forall \left(i,j\right)\in \wp$$
(3)
$$\sum_{k\in K}\sum_{p\in P}{k\cdot a}_{jkp}-\sum_{l\in K}\sum_{p\in P}{l\cdot a}_{ilp}\le ns\cdot {u}_{i} \, \,\,\forall \left(i,j\right)\in \wp$$
(4)
$$\sum_{i\in I}\sum_{p\in P}{{t}_{ip}\cdot a}_{ikp}\le CT \, \, \,\,\forall \, k\in K$$
(5)
$$\psi \cdot {w}_{\mathrm{rk}}\ge \sum_{i\in I}\left({a}_{ik,1+r}+{a}_{ik,nr+1+r}\right) \,\,\forall \, k\in K,r\in R$$
(6)
$$\sum_{r\in R}{w}_{\mathrm{rk}}\le 1 \,\,\forall \, k\in K$$
(7)
$$\sum_{r\in R}\left({c}_{r}\cdot \sum_{k\in K}{w}_{\mathrm{rk}}\right)\le MB$$
(8)
$$\psi \cdot {v}_{k}\ge \sum_{i\in I}{a}_{ik,1}+\sum_{i\in I}\sum_{p\in \{nr+2,\ldots ,2\cdot nr+1\}}{a}_{ikp} \,\,\forall \, k\in K$$
(9)
$$\sum_{k\in K}{v}_{k}\le \mathrm{MW}$$
(10)
$${a}_{ikp},{u}_{i},{w}_{rk},{v}_{k}\in \left\{\mathrm{0,1}\right\}$$
(11)

The second model, referred to as Model 2, is built based on Urban and Chiang (2006). Model 2 utilizes the \({x}_{ikp}\) and \({y}_{ikp}\) to determine the task assignment and is formulated with expression (12) to expression (22). The objective in expression (12) optimizes the cycle time. Constraint (13) is the occurrence constraint, indicating that one task must be allocated to the entrance side or exit side of one station and be operated by one process alternative available. Constraint (14) and constraint (15) tackle the precedence constraint. Constraint (14) demands that the predecessor \(i\) of task \(j\) must be allocated to the former or same station when both of them are allocated to the entrance side; Constraint (15) demands that the predecessor \(i\) of task \(j\) must be allocated to the latter or same station when both of them are allocated to the exit side. Constraint (16) is the cycle time constraint. Constraint (17) and constraint (18) deal with cobot allocation and constraint (19) is the budget constraint. Constraint (20) deals with the worker assignment and constraint (21) restricts the number of human workers.

$$\mathrm{Min}\,\, CT$$
(12)
$$\sum_{k\in K}\sum_{p\in P}\left({x}_{ikp}+{y}_{ikp}\right)=1 \, \, \,\,\forall i\in I$$
(13)
$$\sum_{k\in K}\sum_{p\in P}\left(\mathrm{ns}-k+1\right)\cdot \left({x}_{ikp}-{x}_{jkp}\right)\ge 0 \, \,\,\forall \left(i,j\right)\in \wp$$
(14)
$$\sum_{k\in K}\sum_{p\in P}\left(\mathrm{ns}-k+1\right)\cdot \left({y}_{jkp}-{y}_{ikp}\right)\ge 0 \, \,\,\forall \left(i,j\right)\in \wp$$
(15)
$$\sum_{i\in I}\sum_{p\in P}{t}_{ip}\cdot \left({x}_{ikp}+{y}_{ikp}\right)\le CT \, \, \,\,\forall \, k\in K$$
(16)
$$\psi \cdot {w}_{\mathrm{rk}}\ge \sum_{i\in I}\left({x}_{ik,1+r}+{y}_{ik,1+r}+{x}_{ik,nr+1+r}+{y}_{ik,nr+1+r}\right) \,\,\forall \, k\in K,r\in R$$
(17)
$$\sum_{r\in R}{w}_{\mathrm{rk}}\le 1 \,\,\forall \, k\in K$$
(18)
$$\sum_{r\in R}\left({c}_{\mathrm{r}}\cdot \sum_{k\in K}{w}_{\mathrm{rk}}\right)\le MB$$
(19)
$$\psi \cdot {v}_{k}\ge \sum_{i\in I}\left({x}_{ik,1}+{y}_{ik,1}\right)+\sum_{i\in I}\sum_{p\in \{nr+2,\ldots ,2\cdot nr+1\}}\left({x}_{ikp}+{y}_{ikp}\right) \,\,\forall \, k\in K$$
(20)
$$\sum_{k\in K}{v}_{k}\le \mathrm{MW}$$
(21)
$${x}_{ikp},{y}_{ikp},{w}_{\mathrm{rk}},{v}_{\mathrm{k}}\in \left\{\mathrm{0,1}\right\}$$
(22)

The third model, referred to as Model 3, is formulated based on Model 2 and Li et al. (2017). Model 3 has the same decision variables as Model 2, whereas Model 3 utilizes constraint (23) to replace constraint (14) and constraint (15). Namely Model 3 consist of expression (12), constraint (13), constraint (23) and constraints (1622). Here, constraint (23) is the precedence constraint and it demands that, for a pair of task \(i\) of task \(j\) where task \(i\) is the predecessor of task \(j\), 1) task \(i\) must be allocated to the former or same station when both of them are alloated to the entrance side, 2) task \(i\) must be allocated to the latter or same station when both of them are allocated to the exit side, 3) task \(i\) must be allocated to the entrance side when they are allocated to entrance side and exit side, respectively.

$$\sum_{k\in K}\sum_{p\in P}k\cdot \left({x}_{ikp}-{x}_{jkp}\right)+\sum_{k\in K}\sum_{p\in P}\left(2\cdot ns-k\right)\cdot \left({y}_{ikp}-{y}_{jkp}\right)\le 0 \, \,\,\forall \left(i,j\right)\in \wp$$
(23)

All the three models are MILP models and they are capable of tackling the small-size instances optimally utilizing the CPLEX solver. This study will further evaluate and compare these models in Sect. 6.2.

4 Implemented metaheuristic methodologies

As the considered problem is one NP-hard problem in nature, the formulated models cannot achieve satisfying solutions when solving large-size instance. Hence, this section develops two recent algorithms to tackle the large-size instances effectively. These algorithms utilize new encoding scheme and decoding procedure to achieve feasible solutions. The encoding and decoding, applied algorithms and the utilized neighbor structures are presented in the following subsections in detail.

4.1 Encoding and decoding

As the CRUALBP involves task assignment and the worker and cobot allocation, this study utilizes two vectors for encoding: permutation vector and process alternative selection vector. Task permutation vector determines the task assignment; process alternative selection vector determine the worker and cobot allocation. Figure 3 illustrates one example solution presentation with 11 tasks, four stations and four types of cobots. In the decoding scheme, task permutation is one sequence of all the tasks and the tasks in the former positions are supposed to have the higher priorities. For instance, task 1 and task 11 are in the first two positions and they are allocated to the first station in the solution presentation. Process alternative selection vector is one set of process alternatives on stations. If the code in the \(k\) th position is \(p\), it donates that process alternative \(p\) should be allocated to station \(k\). For instance, one worker with the \(p\) of 1 is allocated to the first station and one worker and cobot 3 with the \(p\) of 8 are allocated to the second station.

Fig. 3
figure 3

Solution presentation

As the encoding scheme is not a feasible solution, this study develops the one decoding procedure in Algorithm 1, where the \({CT}_{Init}\) is the initial cycle time. As you can see, the worker and cobot allocation is first determined utilizing the process alternative selection vector. Subsequently, the tasks are assigned to the entrance side and exit side of stations in sequence. During task assignment, one task is assignable when it satisfies the precedence constraint and the cycle time constraint. Precedence constraint demands that the predecessors or successors have been assigned; cycle time constraint demands that this task can be completed within the initial cycle time \({CT}_{Init}\) by at least one process alternative when \(k<ns\). It is to be noted that, this decoding procedure ignores the cycle time for the last station to obtain the feasible solution for each solution presentation.

figure a

How to determine the proper initial cycle time \({CT}_{Init}\) is one important factor when solving the type II CRUALBP. One possible method is determining the best cycle time with increments. This method starts testing \({\mathrm{CT}}_{LB}\) as the initial cycle time \({CT}_{Init}\), where \({CT}_{LB}=\left(\sum_{i}\underset{p}{\mathrm{min}}{t}_{ip}\right)/\mathrm{ns}\). Subsequently, \({CT}_{Init}\) is updated with \({\mathrm{CT}}_{\mathrm{Init}}\leftarrow {\mathrm{CT}}_{\mathrm{Init}}+1\) until the tasks of the last station can be completed within \({CT}_{Init}\). One possible drawback of this method is that the decoding procedure is conducted for many times to obtain the proper initial cycle time for one encoding scheme. Hence, on the basis of Li et al. (2019a), Li et al. (2019b) and several others, this study proposes the iterative mechanism to update the initial cycle time during the algorithm’ evolution. The proposed iterative mechanism is presented in Algorithm 2, where \({CT}_{Best}\) is the best cycle time during the evolution. It is clear that all the solutions utilize the same initial cycle time and the cycle time is reduced gradually utilizing the \({\mathrm{CT}}_{\mathrm{Init}}\leftarrow {\mathrm{CT}}_{\mathrm{Best}}-1\). The preliminary experiments show that iterative mechanism produces superior performance and hence it is utilized in this study.

figure b

4.2 Improved artificial bee colony algorithm

The proposed IABC is adopted from that in (Çil et al. 2020), where IABC is the best performer among the 11 compared algorithms. The main procedure of the IABC is provided in Algorithm 3. Here, parameter PS is the population size, parameter limit is the number of times before abandoning the individual remained unchanged, and parameter \(r\) within \(\left[\mathrm{0,1}\right]\) is the probability to conduct the local search. IABC starts with initializing the swarm, and afterward conducts the modified employed bee phase, modified onlooker phase, modified scout phase and local search until meeting the termination criterion.

IABC has made following modifications on the employed bee phase, onlooker phase and scout phase, and also optionally utilizes the local search. In this algorithm, modified onlooker phase selects the best and non-duplicated PS solutions from the incumbent swam and the new swam to obtain a high-quality and diverse swarm. Modified scout phase replaces the abandoned individuals with the high-quality neighbor solution. Meanwhile, local search aims at emphasizing the exploitation capacity. In short, the modifications enhance the exploitation and exploration capacity of the proposed IABC method. As you will see in Sect. 6.3, the proposed IABC outperforms the original artificial bee colony algorithm by a significant margin.

figure c

4.3 Improved migrating bird optimization algorithm

The proposed IMBO is adopted from (Janardhanan et al. 2019), where IMBO outperforms five other algorithms. The main procedure of the IMBO is illustrated in Algorithm 4. The proposed IMBO include four original parameters and three new parameters. The four original parameters include the population size (PS), the consecutive times before the leader replacement (m), the number of neighbor solutions for the leader (k) and the number of shared neighbors (x). The three new parameters include the initial temperature (\({T}_{\mathrm{init}}\)), the cooling rate (\(\alpha \)) and parameter rt where T is re-initialized with \(T={T}_{\mathrm{init}}\) when the best cycle time has not been updated for rt times. The proposed IMBO starts with initializing the population, and afterward conducts the modified leader improvement, modified block improvement and leader replacement repeatedly until the termination criterion is met. In the main loop, leader replacement is conducted after executing the modified leader improvement and modified block improvement for m consecutive times.

The modified leader improvement and modified block improvement propose several modifications to enhance the exploitation and exploration capacity as follows. 1) The incumbent individual is replaced with the neighbor solution once the neighbor solution achieves the same or better performance. This modification speeds up the evolution process by replacing the poor-quality individual. 2) New acceptance criterion based on the simulated annealing algorithm is developed to accept the worse solution with a certain probability. Especially, if the best cycle time has not been updated for rt times, the temperature is re-initialized to increase the possibility of accepting the worse solution. This utilization of new acceptance helps the algorithm to escape from being trapped into local optima. 3) The neighbor solution of one individual is set to a very large number when it shares the same fitness value of the original individual. In other words, the neighbor solutions with the same fitness value of the original individual cannot be shared by the individual in the back. This modification aims at preventing the premature of the algorithm by preserving the diversity of the population. As you will see in Sect. 6.3, the proposed IMBO outperforms the original migrating bird optimization algorithm statistically.

figure d

4.4 Utilized neighbor structure

The neighbor structures have important effect on the algorithms’ performance. As this encoding scheme has two vectors, it is necessary to design effective neighbor operators to modify the two vectors to obtain high-quality neighbor solutions. As seen in Fig. 4, this study utilizes the swap operator and insert operator to modify the task permutation vector. Swap operator selects two tasks in different positions randomly and exchange the positions of the two selected tasks. Inset operator selects one task randomly and insert it in a new position different from its original position. As seen in Fig. 4, to modify the process alternative selection vector, this study utilizes both swap operator and mutation operator. Swap operator selects two process alternatives in two different stations randomly and exchanges the selected two process alternatives. Mutation operator selects one process alternative on one station randomly and replaces it with a different process alternative. As two vectors are involved and each vector has two neighbor operators, for one individual, this study first selects one vector randomly to be modified and later selects one of the two corresponding neighbor operators to modify the selected vector.

Fig. 4
figure 4

Utilized neighbor operators

5 An illustrated example

This section provides an example to highlight the features of the considered problem. This instance has 11 tasks, four types of cobots (\(2\cdot 4+1\) process alternatives) and four stations. Table 2 illustrates the precedence relation and the operation times of task by process alternatives (\(P=\left\{\mathrm{1,2},\ldots ,4+\mathrm{1,4}+2,\ldots ,2\cdot 4+1\right\}\)). The cost of these four cobots are 10.11, 12.79, 18.55 and 20.83 units.

Table 2 Precedence relation between tasks and the operation times of tasks

Table 3 presents the achieved cycle times under different budgets on the straight assembly line and U-shaped assembly line. As you can see, the cycle time could be reduced when the budget increase. Specifically, on the U-shaped assembly line, the cycle time is 12 when the budget is 0.0 (no cobot is utilized). And the cycle time is reduced to 8 when the budget increase to 70. Comparing the cycle times on the straight assembly line and U-shaped line, it is clear that for some situations U-shaped assembly line obtains the smaller cycle time. For instance, when the budget is 20.0, straight assembly line obtains a cycle time of 11 whereas U-shaped assembly line obtains a cycle time of 10. Figure 5 depicts the detailed task assignment worker and robot allocation on straight assembly line (Fig. 5a) and U-shaped assembly line (Fig. 5b) within the budget of 20.0. As it can be seen, the workers and cobot on the straight line operate the tasks from the first station to the last station in sequence. The U-shaped assembly line, on the contrary, is divided into entrance side and exit side. The process alternative on station 1 first operates task 1 on the entrance side and later operates task 11 on the exit side, the process alternative on station 2 first operates task 2 and task 3 on the entrance side and later operates task 7 and task 9 on the exit side, the process alternative on station 3 operates the task 4, task 6 and task 10 on the exit side in sequence, and the process alternative on station 4 first operates task 5 on the entrance side and later operates task 8 on the exist side. One clear advantage of the U-shaped assembly line is that it has lager flexibility and it might have the smaller cycle time than the straight assembly line.

Table 3 Achieved cycle times under different budgets
Fig. 5
figure 5

Task assignment worker and robot allocation on straight assembly line and U-shaped assembly line

6 Computational study

This section first presents the experimental design in Sect. 6.1, including the tested instances and running environments. Later, Sect. 6.2 evaluates the formulated model, where the formulated model is tested on the small-size instances. At last, Sect. 6.3 evaluates the developed algorithms and 12 other implemented algorithms statistically.

6.1 Experimental design

This study utilizes the instances consisting of 93 instances. These instances have 22 precedence diagrams, and each precedence diagram corresponds to several station numbers. The smallest-size instance has only 7 tasks and the largest-size instance has 297 tasks. For each instance, there are four types of cobots available and there are nine possible process alternatives for each task. For the operation times of tasks by process alternatives, the operation time by worker is set as the original operation times in the literature, the operation time of robots are set to be larger than that by worker, and the operation times by the collaboration are set to be smaller than that by worker.

To evaluate the performance of the proposed algorithms, they are compared with 12 algorithms. They include the late acceptance hill-climbing algorithm (LAHC)(Li et al. 2019b), simulated annealing algorithm (SA) (Li et al. 2019b), genetic algorithm (GA) (Li et al. 2019b), discrete particle swarm optimization algorithm (DPSO) (Li et al. 2019b), original cuckoo search algorithm (OCS), discrete cuckoo search algorithm (DCS)(Li et al. 2018a), original artificial bee colony algorithm(OABC), two improved artificial bee colony algorithms(ABC1 and ABC2) (Li et al. 2019b), original bee algorithm (OBA), improved bee algorithm (Çil et al. 2020) and original migrating bird optimization algorithm (OMBO).

All the implemented algorithms utilize the same encoding and decoding in Sect. 4.1 and neighbor structures presented in Sect. 4.4. And they terminate when the computation time reaches \(nt\cdot nt\cdot \tau \) milliseconds (ms). And in this study the value of \(\tau \) is set to 10, 20, 30, 40, 50 and 60, respectively, to observe the algorithms’ performance under different running times. To have enough data for statistical analysis, all the algorithms solve the tested instances in ten repetitions and the achieved cycle times are recorded. After completing all the experiments, this study utilizes the relative percentage deviation (RPD) to transfer the achieved results with \(\mathrm{RPD}=100\cdot \left({\mathrm{CT}}_{\mathrm{some}}-{\mathrm{CT}}_{\mathrm{Best}}\right)/{\mathrm{CT}}_{\mathrm{Best}}\), where \({CT}_{some}\) is the cycle time by one algorithm in one run and \({CT}_{Best}\) is the best cycle by all algorithms in ten runs. Here, the smaller value of RPD indicates the better performance.

The mathematical model is solved utilizing the CPLEX solver on the platform of the General Algebraic Modeling System 23.0. As the model cannot solve the large-size instance in acceptable time, only the instances with less than or equal to 70 tasks are solved by the model. And the model terminates once the optimal solution is obtained and the optimality of the achieved solution is proved or the computational time reaches 3600 s (s). The implemented algorithms are programmed with the C +  + programming language of the Microsoft Visual Studio 2015. And the algorithms terminate when the given termination criterion of an elapsed computation time is satisfied. The real experiments are conducted on a set of virtual computers of a tower type of server. Each virtual computer has one virtual processor and 2 GB RAM memory and the tower type of server is equipped with two Intel Xeon E5-2680 v2 processors and 64 GB RAM memory.

6.2 Evaluating the model

This section evaluates the model by solving the small-size instances with budget of 20.0 in Table 4. In this table, LB-All in the third column is the maximum value of the lower bound or the optimal solution by the three models. UB, LB and Time are the upper bound of the cycle time, the lower bound of the cycle time and the consumed computation time, respectively. The last three columns present the best results by IABC, the results by IMBO in ten repetitions and the utilized computation time (\(nt\cdot nt\cdot 60\) ms). Here, one model is capable of solving one instance optimally within the given computation time when UB is equal to the LB. As you can see, the models can solve the small-size instances optimally with less than or equal to 28 tasks, whereas they might not obtain the optimal solution for the remained instances. The proposed algorithms are also achieving the optimal solutions for the small-size instances with less than or equal to 28 tasks, and they might achieve clear superiority over the mathematical models when solving the large-size instances. For instances, IABC and IMBO outperform Model 1, Model 2 and Model 3 when solving the P70 with 13, 18 and 23 stations.

Table 4 Results by the formulated models

To have a better observation of the models’ performance, Table 5 presents the numbers the instances which are solved optimally (#OPT), the average value of the RPD values (RPD-Avg), the average computation time of all the instances (Time-Avg). As the lower bounds can be achieved by the models, here RPD is calculated with \(\mathrm{RPD}=100\cdot \left({\mathrm{CT}}_{\mathrm{some}}-{\mathrm{CT}}_{\mathrm{LB}}\right)/{\mathrm{CT}}_{\mathrm{LB}}\), where \({CT}_{some}\) is the cycle time by one method and \({CT}_{LB}\) is the maximum value of the lower bounds by three models. As you can see, among the three models, Model 1 is the best performer in terms of the #OPT and RPD-Avg; Model 3 is the best performer in terms of the Time-Avg. For the comparison between the models and algorithms, it is clear that IABC and IMBO obtains the smaller values of RPD-Avg within the short computation time, demonstrating that the algorithms outperform the MILP models when solving the large-size instances in term of the solution quality.

Table 5 Summarized results by models

To further analyze the reasons leading to the different performances of the models, Table 6 provides the model statistics when solving the illustrated example in Sect. 5 as an example. As you can see, Model 1 has the smallest number of single variables, the smallest number of nonzero elements and the smallest number of discrete variables, and hence Model 1 achieves the best performance in terms of the #OPT and RPD-Avg. Model 3 has the smallest number of blocks of equations and the smallest number of single equations, and hence Model 3 achieves the best performance in terms of the Time-Avg.

Table 6 Model statistics when solving the illustrated example in Sect. 5

6.3 Evaluating the implemented algorithms

As the parameter values have great impact on the algorithms’ performance, the parameters of each algorithm are calibrated before running the final experiments utilizing the Design of Experiments approach and multi-factor analysis of variance (ANOVA) technique. Firstly, two or three values are selected for each parameter with other parameter fixed. Secondly, all the combinations of the parameter values solve a set of 10 instances with different sizes for 10 times with an elapsed computation time of \(nt\cdot nt\cdot 60\) ms. Finally, ANOVA test is conducted to select the best parameter values where the parameters are selected as controlled factors and the average RPD of all the instances in one run is selected as response variable. This parameter calibration method has been widely in literature (see (Li et al. 2019b), (Çil et al. 2020) and many others) and hence the detailed information are not provided here for space reasons. After parameter calibration, the selected parameter values of the IABC are as follows: population size is 5, parameter limit is 500, and parameter \(r\) is 0.1. And in the modified scout phase, this algorithm obtains 10 neighbor solutions of one abandoned solution by conducting the neighbor operators for two times. The selected parameter values of the IMBO are as follows: population size is 5, parameter k is 11, parameter x is 5, parameter m is 20, parameter \({T}_{\mathrm{init}}\) is 0.5, parameter \(\alpha \) is 0.95 and parameter rt is 500.

After conducting all the experiments, the achieved cycle times are transferred into RPD values. As 10 runs and six computation times are conducted, there are a number of \(93\cdot 10\cdot 6\) data to evaluate one algorithm statistically. Table 7 provides the average RPD by implemented algorithms when τ is equal to 40 and 60. In this table, Instance donates the size of the instances, where the number is the number of tasks. Num. presents the number of cases with different station numbers. The detailed results of the cycle times in each run are not exhibited due to space limit, and they are available upon request. As you can see, under \(\tau =40\), IMBO is the best performer with the overall average RPD of 0.34, IABC is the second-best performer with the overall average RPD of 0.37, and SA is the third-best performer with the overall average RPD of 0.52. When \(\tau =60\), IMBO is also the best performer, IABC is the second-best performer and SA is the third-best performer again.

Table 7 Average RPD values by implemented algorithms when \(\tau =40, 60\)

To have a direct observation of the algorithms’ performance under different computation times, Table 8 provides the overall average RPD values under τ = 10, 20, 30, 40, 50 and 60, respectively. It is observed that all the algorithms produce improvements when the computation time increases. Meanwhile, the proposed IMBO and IABC are the two best performers under all the termination criteria, demonstrating the superiority of the IMBO and IABC over the remained algorithms.

Table 8 Overall average RPD values under different computation times

This study also conducts the statistical analysis to make sure that the observed differences are statistically significant. The ANOVA test is again utilized to analyze the results here, where algorithms and termination criteria are the controlled factors and the average RPD in one run is the response variable. However, the initial ANOVA test shows that normality of the residuals is not satisfied. In fact, this situation lies behind that algorithms have quite different performances, and the normality of the residuals can be satisfied if we only consider the best three algorithms. Meanwhile, the nonparametric Friedman test is also applied to ascertain the findings by the ANOVA test. As the Friedman test can only consider one controlled factor, the Friedman test is conducted for six times utilizing the different results under six termination criteria. The ANOVA test shows that there is statistically significant difference between the algorithms, the termination criteria and the interactions of the two controlled factors. Friedman test shows that there is statistically significant difference between the algorithms under all the computation times. Figure 6 illustrates the Means plot and 95% Tukey HSD confidence intervals for the interactions. Clearly, IABC and IMBO are the two best performers under all the computation time and they outperform the remained algorithms statistically. All in all, the comparative study and the statistical analysis verify the superiority of the IMBO and IABC over the compared algorithms, and, as a consequence, they could be regarded as the effective and powerful methods in solving the considered problem.

Fig. 6
figure 6

Means plot and 95% Tukey HSD confidence intervals for the interactions between algorithms and termination criteria

7 Conclusions and future research

Usage of collaborative robot (cobot) gains more and more applications in recent days to assist the human workers in the assembly line. This research studies the U-shaped assembly line balancing problem with cobots to minimize the cycle time, where several cobots of different types are selected under the budget constraint. Three mixed-integer linear programming (MILP) models are developed to formulate this new problem, and these models are capable of solving the small-size instance optimally. In addition, this study develops two algorithms, improved artificial bee colony algorithm and improved migrating bird optimization algorithm to tackle the large-size instance effectively. The two algorithms utilize the new encoding scheme and decoding procedure to obtain feasible solution and propose several improvements to enhance the exploitation and exploration capacity.

A set of experiments are conducted to evaluate the formulated models and developed algorithms. The experimental results demonstrate that the utilization of collaborative robots helps reduce the cycle time. And the comparative study between the models shows that the models can solve the small-size instances optimally, whereas they cannot obtain satisfying solutions for the large-size instances within the acceptable computation time. The statistical analysis between the algorithms demonstrates that proposed algorithms are the two best performers and they outperform the other 12 implemented methodologies.

The findings of this study can be utilized in the U-shaped assembly line to obtain the smaller cycle time and the developed algorithm might be integrated into the decision support system to achieve high-quality solution to assist the line manager to design the U-shaped assembly lines. Future research stems from developing the branch, bound and remember algorithm to solve the large-size instance optimally to extending the considered problems. This considered problem might be extended by considering the cobot breakdown, the mixed-model U-shaped assembly line, and the multiple objectives utilizing the multi-objective algorithms.