A bi-level math-heuristic algorithm for integrated optimization of task scheduling and storage location assignment in automated storage and retrieval systems

Abstract

The operation of high-bay warehouses for semi-finished cut tobacco is crucial in cigarette manufacturing and is primarily responsible for receiving inbound shipments, aging storage, and on-demand outbound delivery. Many local cigarette manufacturers have implemented automated storage and retrieval systems (AS/RSs) to enhance manufacturing efficiency. Nevertheless, the AS/RS systems encounter task scheduling and storage allocation issues. These challenges include the inability to dynamically adjust slotting promptly in response to real-time inventory status and ineffective task scheduling when meeting order deadlines, which limits the operational efficiency and responsiveness of the system. To address these limitations, this study designs a classification-based storage strategy and constructs a mathematical model that combines dynamic programming with integer linear programming to minimize task completion time and latency. A two-layer mathematical programming heuristic algorithm based on an improved dung beetle optimization (IDBO) algorithm and an assignment model is introduced. In the outer layer of this model, inbound tasks strictly follow the first-come, first-served principle. In addition, the IDBO algorithm is employed to optimize the execution sequence of outbound tasks. Moreover, to enhance solution quality, this study combines an improved circle chaotic map, an adaptive t-distribution perturbation strategy, and a competition mechanism. In the inner layer of the proposed model, the movement cost of a stacker crane is modeled as edge weights in a bipartite graph between the tasks and storage locations due to a fixed task sequence. Furthermore, the assignment problem is efficiently solved by the Hungarian algorithm, yielding an optimized task-to-location matching scheme. The simulation results show that the proposed mathematical programming heuristic algorithm outperforms other optimization methods across different order scales. The proposed algorithm provides satisfactory solutions, thereby significantly improving the operational efficiency and order fulfillment capacity of the AS/RS system. This study performs a collaborative optimization of task scheduling and storage location assignment in high-bay warehouses using a bi-level solution strategy. This provides strong support for meeting the production requirements of the cigarette manufacturing industry and other complex warehousing systems, highlighting the significant engineering application value of the proposed method.

Introduction

The complexity of global supply chains and advancements in information and communication technologies, automation, and Industry 4.0 have been continuously increasing. Accordingly, warehousing systems have evolved toward higher levels of automation, intelligence, and lean operations, and have been widely adopted in the manufacturing sector [1]. In cigarette production, the semi-finished tobacco warehouse represents a critical factor between the production and packaging stages. The core functions of a warehouse include tobacco storage, aging, and dynamic allocation. The operational efficiency of a warehouse affects the continuity of production and directly influences logistics costs and supply chain stability [2]. In traditional cigarette factories, warehouse operations have largely relied on manual labor, resulting in inefficiencies in location assignment and task scheduling, which further causes underutilization of storage space and slow responsiveness for inbound and outbound tasks. In recent years, automated storage and retrieval systems (AS/RS) have been gradually introduced into the tobacco industry to improve overall warehouse efficiency, emerging as a key technology for optimizing storage space and enhancing material handling efficiency [3, 4].

Automated high-bay warehouses for semi-finished products typically comprise high-rise racking systems, stacker cranes, conveyor mechanisms, and a warehouse management system. These systems often use rail-guided vehicles (RGVs) to enable highly efficient automated storage operations [5]. Compared to traditional manual storage practices, AS/RSs offer greater adaptability, which makes them well-suited for complex warehousing scenarios with diverse product types, various batches, and frequent turnover. For instance, industries with storage and handling requirements, such as electronics and automotive manufacturing, can benefit from the flexibility offered by AS/RSs. Although AS/RSs have been successfully applied to various industries, they still face three main challenges in practical operations:

1. In practical AS/RS operations, task assignment often adheres to a first-come-first-served (FCFS) policy, without adequately considering the impact of storage zone constraints on the operational efficiency of a stacker crane. However, inefficient task sequencing may cause the stacker crane to frequently travel across different storage zones, thereby significantly increasing the operational time;

2. Storage location assignment lacks a real-time feedback mechanism. As a result, adjusting the task scheduling dynamically based on current inventory status and operational demands can be challenging, leading to inefficient and delayed resource utilization;

3. In industrial practice, inbound and outbound warehouse operations are closely linked to order fulfillment processes. However, most existing optimization approaches primarily aim to minimize the travel time of stacker cranes, overlooking the deadline constraints associated with the outbound tasks. This omission may result in production cycle delays and inadequate responsiveness to customer orders.

Integrated optimization of task scheduling and storage location assignment in AS/RSs primarily employs a two-stage optimization framework, where the problem is decomposed into two loosely coupled subtasks, solved in isolation: storage slotting and task scheduling. However, most recent studies have focused on minimizing order completion time while omitting consideration of task urgency or time constraints. In addition, existing methods often rely on predefined, static slotting policies that fail to meet fluctuating and time-sensitive storage requirements in dynamic operational environments. To address this problem, this study proposes an integrated optimization framework for task scheduling and real-time slotting decisions in AS/RSs based on a category-based storage strategy. In addition, a bi-level optimization model integrating dynamic programming (DP) and integer linear programming (ILP) is developed to enable the coordinated optimization of task sequencing and storage location assignment. Moreover, the sequencing process of inbound and outbound tasks is optimized by introducing task timeliness constraints, thereby improving the operational efficiency of the warehouse and enhancing order fulfillment punctuality.

The main contributions of this study can be summarized as follows:

1. Joint optimization modeling: A unified mathematical model combining the ILP and DP approaches is constructed to jointly optimize the scheduling sequences of the inbound and outbound tasks and the rules for storage location assignment. The total task completion time and delay penalties are minimized while explicitly incorporating the zone-based storage constraints and task urgency requirements. This formulation enhances the scheduling efficiency and real-time responsiveness of AS/RS operations under dynamic storage requirements;

2. Bi-level mathematical programming and heuristic collaboration: A mathematical programming heuristic (MPH) method is developed. In the outer layer of the proposed method, the IDBO algorithm is employed to optimize outbound task sequencing. Moreover, an enhanced circle chaotic mapping is introduced to improve population initialization, and an adaptive t-distribution perturbation strategy and a competition mechanism are incorporated to avoid premature convergence to local optima. In the inner layer, the Hungarian algorithm is employed to efficiently solve the task–location assignment problem, thereby enabling effective resource allocation under dynamic conditions;

3. Performance validation and adaptability assessment: Extensive simulation experiments are conducted across varying order scales to evaluate the optimization performance of the proposed method. The results indicated that compared to the existing approaches, the proposed method can significantly reduce both task completion time and delay penalties, demonstrating strong system adaptability and practical engineering applicability.

The remainder of this paper is organized as follows. Recent studies on warehouse scheduling are reviewed in Sect. 2. The warehouse scheduling problem is described in Sect. 3. The mathematical model of the integrated task scheduling and slot allocation problem is formulated in Sect. 4. A two-level metaheuristic algorithm combining an improved dung beetle optimization algorithm with an assignment model is introduced in Sect. 5. In Sect. 6, the proposed model and algorithm are numerically verified. Lastly, the study is concluded, and future research directions are discussed in Sect. 7.

Literature review

The warehouse scheduling problem addressed in this study involves the interaction between the task scheduling and slot allocation processes. Unlike the approaches that consider only a single scheduling dimension, a method of the mutual influence between task dispatching and slotting decisions can significantly improve the overall operational efficiency of a system.

The warehouse scheduling problems can be roughly categorized into three classes based on the nature of decision variables:

1. Task scheduling optimization problems,

2. Storage location assignment problems, and

3. Integrated optimization problems with both scheduling and slotting decisions.

The first category includes the task scheduling problems where the storage locations for inbound and outbound tasks are predetermined. The primary objective of this category is to sequence outbound tasks to minimize the total travel distance of the stacker crane. In practice, many enterprises adopt simple dispatching rules to perform task scheduling, such as the FCFS and shortest-distance-first policies. However, these rule-based methods often omit considering zone-related constraints imposed by category-based storage strategies. As a result, the outbound tasks may be restricted to specific storage zones, leading to poor alignment between the task sequencing and slotting processes and degrading the overall system performance.

The second category concerns storage location assignment problems, where the sequence of inbound and outbound tasks is predetermined. The objective of this category is to select optimal storage locations for inbound items while simultaneously minimizing retrieval and storage time and improving space utilization. Most existing slotting strategies focus on static allocation, where storage locations are assigned in advance before task arrival. However, these approaches lack the ability to dynamically adjust storage decisions based on real-time inventory fluctuations and evolve order requirements.

The third category addresses the integrated optimization of task scheduling and storage location assignment based on a list of outbound tasks and available storage locations. The goal of this category is to determine an efficient execution sequence for outbound tasks while selecting suitable storage slots for inbound tasks. Most existing studies on this category have employed a sequential or two-stage approach, either optimizing task sequencing first and then assigning storage locations, or vice versa. However, such decoupled strategies often fail to achieve global optimality due to the inherent interdependencies between the two decision components. Moreover, they generally pay limited attention to the time-sensitivity of outbound tasks, thereby making meeting service-level requirements in real-world warehouses, where certain orders have strict delivery deadlines, challenging.

A comparative review of recent literature is presented in Table 1, where the differences in research focus and methodological characteristics across the three categories are clearly illustrated. The comparison summarizes and contrasts representative studies in several key aspects: task-slotting coordination, objective function formulation, storage strategy configuration, and solution methodologies. Most existing studies have modeled the task scheduling and slot allocation processes as independent problems in separate phases [6, 7, 8, 9, 10, 11, 12, 13–14]. In addition, recent studies have aimed to minimize the task completion time, with limited attention given to the task time sensitivity and delivery deadlines [3, 6, 7, 8, 9–10, 15, 16–17]. However, static or random slotting methods have been predominantly adopted for storage strategies in recent work [8, 16, 18], while category-based and dynamic slotting have remained underexplored. This misalignment with the practical requirements of modern warehousing systems (e.g., storage efficiency and operational flexibility) highlights a significant gap in current research. Furthermore, regarding solution methodology, most existing studies have employed heuristic algorithms to reduce the combinatorial complexity of inbound and outbound task scheduling in warehousing systems [3, 6, 7, 12, 13]. However, fewer studies have used mathematical programming techniques to enhance the search efficiency and feasibility of the obtained solutions.

Table 1. Key characteristics of related studies

References	Decision variables		Objective function(s)		Storage strategy		Solution approach
	Task scheduling	Slotting	Min. completion time	Min. tardiness	Categorical storage	Dynamic slotting	Math programming	Heuristic
Xu et al.(2020)	–	√	√	–	–	√	–	√
Silva et al.(2020)	–	√	√	–	√	–	–	√
Bolanos et al.(2020)	–	√	√	–	√	–(SFS)	√	–
Winkelmann et al.(2024)	–	√	√	–	√	–	–	√
Mokarrari et al.(2025)	√	–	√	–	–	–	√	–
Balogh et al. (2022)	√	–	–	√	–	–	√	–
Yao et al.(2024)	√	–	√	–	–	–	–	√
Falq et al.(2021)	√	–	√	√	–	–	–	√
Chen et al.(2023)	√	–	√	–	–	–	√	–
Zhang et al.(2025)	√	√	√	–	–	√	–	√
Tang et al.(2020)	√	√	√	–	–	√	√	–
Liu et al.(2025)	√	√	√	–	–	–(SFS)	–	√–
He et al.(2022)	√	√	√	–	√	√	√	–
Rizqi et al.(2024)	√	√	√	–	–	–	–	√
Ours	√	√	√	√	√	√	√	√

SFS, Static-fixed storage; √, Optimization addressed; –, Optimization not addressed

Warehouse scheduling problems can be formulated using three main types of models: (1) integer programming (IP) models [19, 20]; (2) queueing network (QN) models [21, 22]; and (3) DP models [23, 24]. The IP models have been widely used to solve the task scheduling and slotting problems under deterministic conditions. However, when the task requirements of AS/RSs exhibit dynamic fluctuations, the IP models often cannot adapt promptly to varying inputs, resulting in suboptimal responsiveness and reduced real-time applicability of the optimization outcomes. The QN models have been commonly employed to analyze the throughput, task waiting times, and resource utilization of AS/RSs, which makes them well-suited for evaluating the overall operational efficiency of warehouses. However, the QN models do not apply to task path optimization and dynamic scheduling tasks, as they lack the structural capacity to represent real-time decision-making and spatial constraints. The DP models are well-suited for solving the scheduling optimization problems that involve multi-stage decision-making because they allow for the dynamic adjustment of task execution sequences and storage location choices. However, with the increase in the problem scale and the number of decision variables, the state space expands exponentially, leading to the so-called state explosion problem, which makes the application of the DP models to large-scale AS/RS systems challenging. The empty storage locations created by outbound operations can be immediately used by subsequent inbound tasks because of the dynamic nature of slotting. As a result, the status of storage locations on the racks continuously evolves with task scheduling decisions. Consequently, the decision-making at each stage of the DC operations becomes interdependent through these state transitions. Therefore, this study develops an approach that integrates the DP models with ILP to optimize task sequencing and storage location assignment.

Once the mathematical model is established, an efficient solution method should be properly selected to obtain a feasible solution that satisfies all constraints in a limited computation time. Heuristic algorithms have been widely employed in warehouse scheduling to reduce the complexity of task sequencing problems. These algorithms can be roughly categorized into rule-based heuristics and intelligent optimization heuristics. Rule-based heuristic algorithms are grounded in predefined dispatching rules and heuristic criteria. They offer fast computation and are well-suited for small-scale and real-time scheduling problems. These algorithms typically rely on simple decision rules, such as priority-based scheduling, greedy strategies, and local search, to rapidly generate feasible solutions [25]. However, due to their limited global search capability, rule-based heuristic algorithms are prone to being trapped in local optima and struggle to satisfy complex scheduling requirements.

Intelligent optimization heuristics explore larger solution spaces using advanced search mechanisms to identify superior solutions, which makes them suitable for solving complex and large-scale scheduling problems. These algorithms are inspired by biological behaviors, physical phenomena, and mathematical models and employ adaptive search strategies to optimize task scheduling decisions. Representative intelligent optimization heuristics include the genetic algorithm (GA), ant colony optimization, and simulated annealing [26]. However, the number of decision variables in task scheduling and storage allocation increases exponentially with the scale of warehouse systems. In such high-dimensional search spaces, heuristic algorithms can exhibit slow convergence and struggle to identify high-quality solutions in a limited time interval. Hence, this study proposes an MPH that integrates the IDBO algorithm and an assignment model to solve the warehouse scheduling problem by optimizing the coordinated decision-making of task scheduling and storage location assignment. Unlike conventional heuristic methods, the proposed method employs the IDBO algorithm to sequence inbound and outbound tasks. In addition, the IDBO algorithm is coupled with an IP model to identify high-quality storage locations, thereby enhancing solution stability and global optimality.

Although the previous work has demonstrated significant advancement in the field, several research gaps remain:

1. The existing studies have not sufficiently addressed the coordinated optimization of task scheduling and storage location assignment. Most proposed approaches adopt a sequential optimization strategy, where task scheduling is optimized first, followed by storage assignment, or vice versa. However, these decoupled strategies overlook the mutual influence between the task sequencing and storage assignment processes. As a result, poor alignment may occur between the task execution order and storage location allocation, which impacts the overall operational efficiency of the system. To address this shortcoming, this study develops a mathematical programming heuristic algorithm that integrates an improved dung beetle optimization (IDBO) algorithm with an assignment model. The main goal is to achieve a synchronized optimization of the task scheduling and storage location assignment processes, thereby enhancing the overall operational efficiency of the system;

2. The application of dynamic storage location assignment in complex warehousing environments remains limited. Most recent studies have employed static slotting strategies, where storage locations are pre-assigned before task arrival, lacking effective mechanisms for adaptive adjustment. However, outbound operations continuously release empty storage locations in real-world AS/RS environments. In addition, fixed slotting strategies fail to accommodate dynamic task requirements and may reduce the utilization level of storage resources. By contrast, the proposed method can dynamically adjust storage allocation during task scheduling optimization, thereby significantly improving space utilization efficiency.

Problem description

This section provides a detailed explanation of the integrated optimization problem of task scheduling and storage location assignment in a category-based AS/RS. First, the system structure, stacker crane operations, and the principles of category-based storage are introduced. Then, the optimization objectives and modeling assumptions are explained, laying the groundwork for subsequent mathematical modeling and algorithm design. The intralogistics layout of a high-bay warehouse for semi-finished tobacco products is illustrated in Fig. 1.

[See PDF for image]

Fig. 1

Intralogistics layout of the high-bay warehouse for semi-finished tobacco products

High-bay warehouse system for semi-finished tobacco products

The high-bay warehouse system used in this study adopts a category-based storage strategy. In this strategy, tobacco cases are classified according to their historical turnover frequency and stored in designated storage zones. As illustrated in Fig. 1, the warehouse is physically structured into multiple parallel aisles, which are equipped with a stacker crane that operates independently along a dedicated rail, performing storage and retrieval operations on both sides of the aisle. Each stacker crane handles only a single tobacco case per operation. The starting and ending points of task execution are located at a unified input/output (I/O) port at the aisle's front end. In addition, the storage racks are spatially divided into three functional zones to accommodate the scheduling requirements of tobacco cases with varying turnover frequencies, as illustrated in Fig. 2. The category-specific zones are arranged in order of increasing distance from the I/O port, corresponding to high, medium, and low-frequency tasks, which are denoted as Classes A, B, and C, respectively, to minimize travel costs for high-frequency operations.

[See PDF for image]

Fig. 2

The schematic diagram of the storage racks

The inbound and outbound tasks adopt a DC strategy, where each stacker crane run performs a storage task followed immediately by a retrieval task. This approach enhances equipment utilization while reducing the empty travel cost. The specific operational procedure is as follows. The incoming tobacco cartons are transferred by the RGV and unloaded onto the conveyor at the aisle’s front end. Then, cartons are transported by the conveyor to the aisle’s I/O point, where they are received and delivered to the designated storage location $S1$ by the stacker crane. Subsequently, the stacker crane moves to the storage location $R1$ associated with the outbound task, retrieves the goods, and returns to the I/O point, thereby completing a full DC operation. All inbound tobacco cases are sequentially queued on the conveyor belt at the front end of each aisle upon arrival, waiting to be retrieved by the stacker crane for storage. The sequence of cases cannot be altered once they enter the belt because of the structural constraints of the conveyor. As a result, inbound tasks strictly follow the FCFS policy in actual operations.

Optimization objectives and model assumptions

The problem addressed in this paper involves the integrated optimization of task sequencing and storage location assignment, where the decision variables are highly interdependent. Therefore, a joint optimization method is required to achieve optimal system-level performance. The main optimization objective is to minimize the total task completion time and the tardiness of outbound tasks. The task completion time reflects the throughput capacity of a system, while task tardiness captures the time-sensitivity requirements of customers.

The following assumptions are made to simplify the modeling process and ensure computational feasibility:

1. A warehouse adopts a category-based storage strategy;

2. No inventory shortages occur in a warehouse, and the stacker cranes can perform both single- and dual-command cycles;

3. The horizontal and vertical speeds of all stacker cranes are constant;

4. The time for inbound and outbound operations (i.e., handling and unloading) is negligible, and only travel time for movement is considered;

5. Each DC processes a single pair of the inbound and outbound tasks;

6. Inbound tasks adopt the FCFS policy, and the sequence of outbound tasks is adjustable.

DP model construction

This study proposes a framework that combines DP integrated with ILP formulation to achieve efficient scheduling of inbound and outbound tasks and rational allocation of storage resources in an automated warehousing system. The DP model is employed to characterize the sequencing process of outbound tasks, as well as the dynamic evolution of storage states across decision stages. In addition, the ILP model is used to determine the storage location assignment for inbound and outbound tasks for a given static storage configuration. As illustrated in Fig. 3, the DP model defines each decision stage based on a stacker crane's completion time of a DC operation. The ILP model is solved at each stage based on the current storage states to minimize the task completion time and outbound task tardiness. The ILP model assigns specific storage locations to the pending inbound and outbound tasks, computing their completion times and tardiness values. After executing each stage task, the storage state is updated and used as an input for the subsequent stage. This process is continued until all tasks are scheduled.

[See PDF for image]

Fig. 3

The flowchart of the DP model

Symbol definitions and path calculation

The symbols in this study, as well as their definitions, are listed in Table 2. If the decision variable $x_{{\text{bhk}}_1{j}_1{k}_2{j}_2}^o$ equals one, the inbound cargo is assigned to the storage location $(k_1,j_1)$ on the $h$ th rack in the $b$ th scheduling phase. As a result, the outbound cargo is retrieved from a storage location $(k_2,j_2)$ on the same rack. These two actions complete the oth outbound task, performing a dual-command (DC) operation. The completion time for the path between the inbound and outbound storage locations is denoted as $T_{\text{total}}$, and it is calculated using the Chebyshev distance method. This value represents the travel time of a stacker crane when executing the DC task. The total time $T_{\text{total}}$ is calculated by:

1

$$T_{\text{total}}=\underset{I/O\to \text{inbound}}{\underbrace{max\left(\frac{k_1}{v_x},\frac{j_1}{v_y}\right)}}+\underset{\text{inbound}\to \text{outbound}}{\underbrace{max\left(\frac{|k_1-k_2|}{v_x},\frac{|j_1-j_2|}{v_y}\right)}}+\underset{\text{outbound}\to I/O}{\underbrace{max\left(\frac{k_2}{v_x},\frac{j_2}{v_y}\right)}},$$where $(k_1,j_1)$ denotes the inbound storage location, $(k_2,j_2)$ indicates the outbound storage location, and $v_x$ and $v_y$ represent the horizontal and vertical speeds of the stacker crane, respectively.

Table 2. Model symbols and definitions

Type	Symbol	Meaning
Task set	${\mathcal{T}}^{\text{out}}$	Outbound task set
	$o\in {\mathcal{T}}^{\text{out}}$	Outbound task index, $o=1,2,...,n$
	${\mathcal{T}}^{\text{in}}$	Inbound task set
	$\theta \in {\mathcal{T}}^{\text{in}}$	Inbound task index, $\theta =1,2,...,n$
	${\chi }_o$	Item type associated with an outbound task o
	$D_o$	The deadline of an outbound task o
	$L_o$	The tardiness of an outbound task o
Shelf structure parameters	$k$	Horizontal storage index, $k=1,2,...,K$
	$j$	Vertical storage index, $j=1,2,...,J$
	$h$	The hth rack in the aisle, $h=1,2$
Stacker crane operational parameters	$v_x$	The horizontal travel speed of the stacker crane
	$v_y$	The vertical travel speed of the stacker crane
	$p_b$	The execution time of the DC tasks by the stacker crane in a stage b
	$T_b$	The completion time of the DC tasks by the stacker crane in a stage b
Others	$i$	Types of goods, $i=1,2,3$
	$Z_i$	The set of categorized storage zones assigned to inbound and outbound goods, where i denotes three different categories of storage zones, and $i=1,2,3$
	$Z_{(k,j)}$	Each storage location ($k,j$) is mapped to the corresponding categorized storage zone
	${\tau }_i$	The aging time requirement of an item type $i$
	${\tau }_i^a$	The elapsed aging time of an item type $i$
	$M$	A sufficiently large positive constant used in the constraint
	$w$	The weighting coefficient of the objective function,$w\in \left[0,1\right]$
	${\text{SE}}^{\text{fb}}$	Given the scheduling task sequence f, this denotes the set of available empty storage slots in the rack area assigned to the stacker crane before stage b begins
	${\text{SR}}_i^{\text{fb}}$	Given the scheduling task sequence f, this denotes the set of storage slots containing an item type iii in the rack area assigned to the stacker crane before stage b begins

Stage division and state transition of DP model

In this study, n DC operations are treated as stages to capture dynamic changes in storage states during the task execution process, where each stage corresponds to the execution of a single DC operation. The state variables are denoted by ${\text{SE}}^{f(b+1)}$ and ${\text{SR}}_i^{f(b+1)}$, and they are, respectively, expressed by Eqs. (2) and (3).

2

$$\begin{array}{c}\text{SE}_{}^{f(b+1)}=\text{SE}_{}^{\text{fb}}-\{k_1{j}_1|x_{{\text{bhk}}_1{j}_1{k}_2{j}_2}^o=1,k_1{j}_1\in \text{SE}_{}^{\text{fb}},k_2{j}_2\in {\text{SR}}_i^{\text{fb}}\}\\ +\{k_2{j}_2|x_{{\text{bhk}}_1{j}_1{k}_2{j}_2}^o=1,k_1{j}_1\in \text{SE}_{}^{\text{fb}},k_2{j}_2\in {\text{SR}}_i^{\text{fb}}\}\end{array}$$

3

$${\text{SR}}_i^{f(b+1)}={\text{SR}}_i^{\text{fb}}-\{k_2{j}_2|x_{{\text{bhk}}_1{j}_1{k}_2{j}_2}^o=1,\underset{\phantom{a}}{\phantom{a}}{k}_1{j}_1\in {\text{SE}}_{\phantom{a}}^{\text{fb}},\underset{\phantom{a}}{\phantom{a}}{k}_2{j}_2\in {\text{SR}}_i^{\text{fb}}\},$$

Expression in Eq. (2) represents the removal of storage locations $k_1{j}_1$ used for inbound tasks. Once the empty storage locations $k_2{j}_2$ generated after the outbound task are added, the set of empty storage locations in a warehouse can be updated. Similarly, Eq. (3) indicates the removal of storage locations $k_2{j}_2$ released after task completion, updating the set of occupied storage locations in the warehouse.

The optimal value function of DP is defined by:

4

$$f_b(o,R)=\underset{d\in R}{min}\left\{f_{b-1}(d,R\backslash \{d\})+{\text{Obj}}_o(S{E}^{\text{fb}},S{R}_i^{\text{fb}})\right\}$$

The initial conditions are given by:

5

$$f_b(o,\mathrm{\varnothing })=\underset{d\in R}{min}\left\{{\text{Obj}}_o({\text{SE}}^0,{\text{SR}}_i^0)\right\}$$

In Eq. (4), $f_b(o,R)$ is the minimum cumulative operation time of the system when the oth outbound task is scheduled after completing a task set R; $f_{b-1}(d,R\backslash \{d\})$ denotes the minimum operation time at stage (b – 1) after removing task d from a set $R\backslash \{d\}$, where $r=1,2,\dots ,n$, and $b=1,2,\dots ,n$. In Eq. (5), $f_b(o,R)$ denotes the minimum objective value when executing an outbound task o with an initially empty task set; ${\text{Obj}}_o$ is the weighted sum of the completion time and tardiness for task o in system states $SE_{}^{\text{fb}}$ and ${\text{SR}}_i^{\text{fb}}$; $Ob{j}_o$ is the parameter obtained by solving the integer linear programming model.

Intra-stage integer linear programming model

The DP framework incorporates an integer linear programming sub-model within each stage to derive the optimal scheduling solution at each decision point. The objective at each stage is to minimize the weighted sum of the completion time and tardiness of the oth task, which can be expressed as follows:

6

$$min{\text{Obj}}_o=w·T_b^{}+(1-w)·L_o$$

The constraints are defined as follows:

1. Unique assignment constraint: A single inbound location and a single outbound location are selected to execute the dual-command (DC) operation as follows:

   7

   $$\sum _{h=1}^2\sum _{k_1=1}^K\sum _{k_2=1}^K\sum _{j_1=1}^J\sum _{j_2=1}^J{x}_{{\text{bhk}}_1{j}_1{k}_2{j}_2}^o=1$$

1. Category zone consistency constraint: Inbound goods must be assigned to empty storage locations within the corresponding category zones. However, outbound goods must be retrieved from occupied storage locations in the same category zones, which can be formulated as follows:

   8

   $$x_{bh{k}_1{j}_1{k}_2{j}_2}^o\le 1,k_1{j}_1\in \text{SE}_{}^{\text{fb}};Z_{(k_1,j_1)}\in Z_i;k_2{j}_2\in {\text{SR}}_i^{\text{fb}};Z_{(k_2,j_2)}=Z_i$$

1. Aging time constraint: Selecting storage locations for outbound tasks is restricted to the locations where the aging duration meets the preset requirement for a threshold, ensuring that the tobacco undergoes sufficient aging before outbound operations, which can be expressed by:

   9

   $${\tau }_i^a\ge {\tau }_o-M\times \left(1-x_{{\text{bhk}}_1{j}_1{k}_2{j}_2}^o\right)$$

1. Task completion time constraint: The completion time of the current-stage task, denoted by $T_b^{}$, must be equal to or longer than the completion time of the previous-stage task plus the stacker crane's operation time, ensuring sequential task execution, which can be expressed as follows:

   10

   $$T_b^{}\ge T_{b-1}^{}+p_{b-1}$$

1. Task tardiness definition constraint: Tardiness $L_o$ of an outbound task o is defined as a difference between its completion time $T_{\text{b}}$ and the corresponding deadline $D_{\tau }$ as follows:

   11

   $$L_o\ominus T_b-D_{\tau }-M\times (1-x_{{\text{bhk}}_1{j}_1{k}_2{j}_2}^o),$$

   12

   $$L_o\ge 0.$$

Algorithm design

As previously explained, the integrated optimization of the task scheduling and storage location assignment tasks represents a combinatorial optimization problem characterized by multi-constraint coupling and multi-stage interdependence. This problem is subject to multiple interrelated constraints, including task-specific storage zones, dynamic evolution of storage states, and the reachability of stacker crane paths. These factors introduce various challenges, including strong coupling between task sequencing and storage state transitions, rapid growth in decision variable dimensions, and historical dependence in decision-making across stages. Although this problem can be theoretically solved by exact algorithms, the computational complexity increases exponentially with the task volume and the number of storage locations. As a result, traditional methods often fail to obtain high-quality solutions within an acceptable time interval in practical applications. Therefore, in this paper, a mathematical programming heuristic algorithm that integrates an assignment model with an IDBO algorithm is developed to solve the high-dimensional and complex scheduling problem. The mathematical programming-based heuristic algorithm combines the global search capability of metaheuristics and the local refinement performance of mathematical models, which makes it suitable for the real-time optimization of large-scale scheduling systems. As shown in Fig. 4, the proposed algorithm consists of two main layers, which are as follows:

1. The IDBO layer encodes the execution sequence of inbound and outbound tasks in the outer-layer module. The IDBO algorithm guides the search process toward a globally optimal task sequencing solution by simulating four social behaviors of dung beetles: rolling, breeding, foraging, and stealing behaviors [27]. Simultaneously, constraints of the first-come and FCFS for inbound tasks are satisfied;

2. The bipartite graph layer is constructed in the inner-layer module based on the task sequence obtained from the outer layer. In this graph, the edge weights represent the travel cost of a stacker crane between the tasks and storage locations. The matching process is then formulated as a minimum-cost assignment problem. Finally, this problem is solved by the Hungarian algorithm to obtain an optimal location assignment for a given task sequence.

[See PDF for image]

Fig. 4

The flowchart of the mathematical programming heuristic algorithm

Encoding scheme

The warehouse task scheduling problem represents a discrete combinatorial optimization problem, whereas the DBO algorithm was originally developed for continuous domains. Therefore, this study adopts the largest ranking value (LRV) rule to map continuous priority values into a discrete task sequence by adapting the DBO algorithm to this discrete setting. Priority values are used to represent candidate solutions in the DBO population. As illustrated in Fig. 5, the LRV rule transforms a set of continuous values representing an individual in the population into a task sequence by sorting them in descending order. This approach facilitates the application of the DBO algorithm to the warehouse task scheduling problem, thereby improving solution efficiency and accuracy.

[See PDF for image]

Fig. 5

Illustration of the LRV rule encoding process

This study proposes a two-layer encoding structure to comprehensively describe a feasible warehouse operation scheduling scheme, representing the scheduling sequence of inbound and outbound tasks while maintaining operability for subsequent optimization procedures. Specifically, the first layer of the encoding structure represents the scheduling sequence of inbound tasks, where a continuous priority value is assigned to each task. The inbound task sequence is fixed according to the task arrival order, ensuring FCFS compliance.

The second layer of the encoding structure represents the scheduling sequence of outbound tasks, where each task is assigned a continuous priority value. The outbound task sequence is generated by sorting these values in descending order. It should be noted that each individual in the population comprises two sets of continuous values, having a length n, thereby forming a complete scheduling solution for the inbound and outbound tasks.

Population initialization based on improved circle chaotic map

Individuals may cluster excessively in local regions in random population initialization, causing uneven distribution and reduced global search capability. To address this problem, this study introduces a chaotic map to enhance the initial population diversity. Using chaotic maps can provide a uniform distribution during initialization, which improves the algorithm’s global exploration capability. Commonly used chaotic mapping methods include the logistic, cubic, and circle maps. Compared to the other chaotic maps, circle maps exhibit a more uniform distribution in non-boundary regions of the population and achieve excellent performance in mitigating the edge aggregation phenomenon [28]. The circle chaotic map can be mathematically expressed as follows:

13

$$x_{n+1}=mod(x_n+0.2-\frac{0.5}{2\pi }sin(2\pi \times x_n),1).$$

However, standard circle maps suffer from uneven distribution of chaotic values across the search space, resulting in local concentration of individuals in certain regions. Therefore, this study proposes an improved version of the circle chaotic map to optimize the spatial distribution of individuals, enhance population diversity, and improve global search performance. The modified formulation is defined as follows:

14

$$x_{n+1}=mod(2x_n+0.3-\frac{0.7}{2\pi }sin(2\pi \times x_n),1).$$

In this study, a dimensionality of $n=3000$ is selected. Compared to the original circle chaotic map, the improved map has a significantly more uniform distribution of chaotic values, as illustrated in Fig. 6. This prevents the excessive clustering of individuals in the search space, ensuring population diversity and enhancing global exploration capability during the optimization process.

[See PDF for image]

Fig. 6

Comparison of population distributions of different chaotic maps

Dynamic slotting strategy under category-based storage

An effective slotting strategy based on the task sequences generated by the IDBO algorithm and its encoding scheme is employed in this study to assign storage locations for inbound and outbound tasks. Due to the category-based zoning constraints, each task is assigned to an appropriate location within its designated storage zone. Specifically, outbound tasks are assigned to the slot closest to the I/O point. In contrast, inbound tasks are allocated to the slot, which results in the shortest DC operation. This strategy enables reassigning outbound locations as target locations for subsequent inbound tasks, thereby avoiding underutilization of high-quality storage locations, which could negatively affect operational efficiency.

In addition, a bipartite graph matching framework is introduced to model the mapping relationships between the tasks and storage locations, as well as to address the assignment problem between the inbound and outbound tasks under these constraints. Specifically, a bipartite graph $G=\left(U,V,E\right)$ is defined, where a left node set, denoted by U, represents the set of pending outbound or inbound tasks, and the right node set, labeled as V, represents the set of available outbound or inbound storage locations. The cardinalities of U and V are equal, that is, $\left|U\right|=\left|V\right|$, which conforms to the standard structure of an assignment problem. In addition, dummy tasks and locations are introduced to balance the two sets in the case of a mismatch in the node number. Each edge $(j,k)\in E$ represents an assignment of a task $j\in U$ to a storage location $k\in V$, and it is associated with a weight $c_{\mathit{jk}}$, which indicates the operational distance incurred by the assignment. The assignment cost for edges involving dummy tasks and locations is set to ensure that virtual assignments do not affect the overall objective value. Hence, an IP model is designed to achieve optimal matching that minimizes the total operational cost, which can be formulated as follows:

15

$$min\sum _{j\in RR}\sum _{k\in L}{c}_{\mathit{jk}}{y}_{\mathit{jk}},$$

The constraints are defined as follows:

1. Each outbound task is assigned to at most one storage location, which can be expressed by:

   16

   $$\sum _{k\in V}{y}_{\mathit{jk}}=1,\forall j\in U;$$

2. Each outbound storage location is assigned to at most one outbound task as follows:

   17

   $$\sum _{j\in U}{y}_{\mathit{jk}}=1,\forall k\in V;$$

3. The decision variable is binary, representing whether the allocation is made, and it is defined by:

   18

   $$y_{\mathit{jk}}\in \left\{0,\left(1\right\}\right).$$

Fitness evaluation

The main objective of this study is to minimize the weighted sum of the completion time for inbound and outbound tasks and the tardiness of outbound tasks. Therefore, the objective function is defined as a fitness function in the proposed heuristic algorithm as follows:

19

$${\text{fitness}}_p=w\times \sum _{b=1}^n{T}_b+(1-w)\times \sum _{o=1}^n{L}_o,\underset{\phantom{a}}{\phantom{a}}p=1,2,...,N_{\text{pop}},$$where ${\text{fitness}}_p$ is the fitness value of a dung beetle individual $p$; $N_{\text{pop}}$ is the population size; $T_b$ is the total completion time of all tasks executed by a stacker crane; $L_o$ is the tardiness of the oth outbound task; $w$ is the weight of the objective for minimizing task completion time, and it was set to 0.5 in the experiments.

Four types of individuals in dung beetle optimization algorithm

1. Rolling dung beetles

This behavior simulates dung beetles rolling a dung ball through the search space along a given direction. The rolling process comprises two modes: obstacle-free and obstacle-present modes. Without obstacles, rolling beetles rely on celestial cues such as sunlight intensity to guide their movement. The position update in the obstacle-free mode is defined by:

20

$$P_i^{t+1}=P_i^t+\alpha ·k·P_i^{t-1}+b·\left|P_i^t-P_{\text{worst}}^t\right|,$$where $P_i^t$ is the position of the ith dung beetle in the population during the tth iteration of the rolling operation; $\alpha$ is a natural coefficient assigned to the value of − 1 or one, with $\alpha =1$ indicating a forward movement toward the target, and $\alpha =-1$ showing a movement away from the target; $k\in \left(0,0.2\right]$ is the deflection coefficient that adjusts the angular deviation; $b\in \left(0,1\right)$ is a constant used to simulate variations in light intensity; $P_{\mathit{worst}}^t$ is the worst-performing individual in the current population; $\left|P_i^t-P_{\mathit{worst}}^t\right|$ is a term used to modulate the brightness-based perturbation strength.

When obstacles are encountered, dung beetles dance to adjust their direction and recalculate a new movement trajectory. The position update in the obstacle-present scenario is defined by:

21

$$P_i^{t+1}=P_i^t+tan(\theta )·\left|P_i^t-P_i^{t-1}\right|,$$where $\theta \in \left[0,\pi \right]$ represents the rolling direction; when $\theta =0,\pi /2,or\pi ,$ the position of the dung beetle remains constant, indicating a navigation failure or hesitation in movement;

1. Reproductive dung beetles

For the reproduction behavior, a dynamic boundary selection strategy is adopted to simulate the safe oviposition area of the dung beetle, which is expressed by:

22

$$\left\{\begin{array}{c}L_b^{\ast }=max\{P_{\text{lbest}}·(1-R),L_b\}\\ U_b^{\ast }=min\{P_{\text{lbest}}·(1+R),U_b\}\end{array}\right),$$where $L_b^{\ast }$ and $U_b^{\ast }$ denote the lower and upper bounds of the spawning region, respectively; $P_{\text{lbest}}$ is the current local best position.

The dynamic contraction coefficient is defined as $R=1-t\left./\phantom{1-t{T}_{max}}\right)T_{max}$, where $T_{max}$ is the maximum number of generations, and t is the current generation index. Parameters $L_b$ and $U_b$ in Eq. (22) indicate the lower and upper bounds of the entire population’s activity range, respectively. Once the spawning region is determined, the initial position of the offspring beetle is dynamically calculated during the evolutionary process as follows:

23

$$P_i^{t+1}=P_{\mathit{lbest}}+b_1·(P_i^t-L_{}^{\ast })+b_2·(P_i^t-U_{}^{\ast }),$$where $b_1$ and $b_2$ are 1 × D-dimensional independent random vectors generated to introduce diversity in the offspring generation process;

1. Foraging dung beetles

The foraging region is also dynamically adjusted using a boundary selection strategy as follows:

24

$$\left\{\begin{array}{c}L_b^{\ast }=max\{P_{\text{gbest}}·(1-R),L_b\}\\ U_b^{\ast }=min\{P_{\text{gbest}}·(1+R),U_b\}\end{array}\right),$$where $P_{\text{gbest}}$ represents the current global best position, and $L_b$ and $U_b$ denote the lower and upper bounds that define the local foraging region, respectively.

Foraging dung beetles perform localized search behavior in this bounded region. The position update is expressed by:

25

$$P_i^{t+1}=P_i^t+c_1·(P_i^t-L_{_b^{}}^{\ast })+c_2·(P_i^t-U_b^{\ast }),$$where $c_1$ is a D-dimensional random vector sampled from a normal distribution, and $c_2$ is a random scalar from the (0,1) interval;

1. Stealing dung beetles

Stealing dung beetles attempt to snatch food near the best food source, which corresponds to the global best position in the population. The position update for stealing dung beetles is defined by:

26

$$P_i^{t+1}=P_{\text{gbest}}^t+S·g·\left(\left|P_i^t-P_{\text{lbest}}\right|+\left|P_i^t-P_{\text{gbest}}\right|\right),$$where $S$ is a constant scaling factor, and $g$ is a 1 × D-dimensional random vector sampled from a normal distribution.

Competition mechanism

The imperialist competitive algorithm (ICA) is a metaheuristic inspired by the socio-political process of imperialistic competition and colonization. Due to its distinctive competitive mechanism, the ICA has demonstrated promising performance in solving many optimization problems [29]. Inspired by the core idea of the ICA, this study integrates its competitive mechanism into the dung beetle optimization algorithm to enhance global search capability. Specifically, the competition is manifested in the interaction between stealing and rolling dung beetles, simulating a dynamic rivalry for resource occupation and dominance.

In each iteration, a stealing dung beetle randomly selects a rolling dung beetle and moves toward its position to initiate a theft preparation phase. This process is analogous to the assimilation mechanism in the ICA, where colonies gradually converge toward their imperialist. In addition, during the update process, a stealing dung beetle may perform a random directional jump with a certain probability. This jump mimics the revolution mechanism in the ICA, thereby helping the algorithm escape from local optima. Once a stealing beetle has successfully seized the dung ball from a rolling beetle, its role is swapped. This process corresponds to the inversion mechanism in the ICA, which enhances the algorithm’s adaptability and solution diversity.

The two mechanisms are as follows:

1. Assimilation mechanism

In each iteration, a stealing dung beetle randomly selects a rolling dung beetle and moves for a certain distance toward its position to initiate the theft preparation phase. This process resembles the assimilation mechanism in the ICA, where the stealing beetle uses the positional information on the rolling beetle to guide its movement. The stealing beetle continuously adjusts its position and incrementally absorbs the superior traits of the target via assimilation. The specific position update rule is defined by:

27

$$P_i^{\text{tnew}}=P_i^t+p\times q\times \frac{P_i^t-R_{\text{rand}}^t}{\left.∥P_i^t\right)-\left(R_{\text{rand}}^t\right.∥},$$where $P_i^t$ is the position of the stealing beetle in the tth generation; $R_{\mathit{rand}}^t$ is the randomly selected position from the rolling beetle population; $p$ is a random scalar that controls the step size of the stealing beetle’s movement; $q$ is a directional perturbation coefficient that introduces additional randomness into the update process;

1. Revolution mechanism

The revolution mechanism simulates the revolution process in the ICA. A stealing beetle, with a certain probability during its movement, may be attracted by another rolling beetle, causing a directional jump and a deviation from its original trajectory.

The movement amplitude $p$ is calculated by:

28

$$p(t)=\left(2-\frac{t}{T_{\text{max}}}\right)·\left.∥P_i^t\right)-\left(R_{\mathit{rand}}^t\right.∥.$$

The directional perturbation coefficient $q$ is defined by:

29

$$q=\left\{\begin{array}{cc}tan(\phi ),\hfill & 0\le \xi \le 0.5\\ 0,\hfill & 0.5<\xi <1\end{array}\right),$$where $\phi$ is a random number uniformly distributed in the $[-\frac{\pi }{4},\frac{\pi }{4}]$ interval, and $\xi$ is a random number from the $[0,1][0,1]$ range;

1. Inversion mechanism

Once a stealing dung beetle has successfully seized the dung ball from a rolling beetle, its role is immediately swapped. This mechanism is analogous to the inversion strategy in the ICA, where dominance shifts based on the competitive outcomes. The corresponding update rule is defined by:

30

$$\text{switch}(P_i^{\text{tnew}},R_{\text{rand}}^t)=\left\{\begin{array}{cc}{P}_i^{\text{tnew}},\hfill & iff(P_i^{\text{tnew}})<f(R_{\text{rand}}^t)\\ \mathrm{\varnothing },\hfill & \text{otherwise}\end{array}\right).$$

Adaptive student's t-distribution perturbation strategy

In the original DBO, position updates in the foraging phase are achieved through simple local random perturbations. However, this strategy lacks sufficient global exploration capability due to its limited perturbation scope, especially during the early stages of the iteration process. Therefore, an adaptive perturbation strategy based on the Student’s t-distribution is introduced into the foraging phase. This strategy adaptively adjusts the number of degrees of freedom of the t-distribution to dynamically balance wide-ranging exploration in the early iterations with fine convergence in the later stages. As a result, the algorithm’s ability to escape local optima is enhanced. The update formula for the perturbed new position is defined as follows:

31

$$P_i^{t+1}=P_i^t+P_i^t·\delta ,$$where $\delta$ is a perturbation term sampled from a Student’s t-distribution whose number of degrees of freedom,$v$, is dynamically adjusted according to the current iteration number as follows:

32

$$v=v_{min}+\left(v_{max}-v_{min}\right)·\frac{t}{T_{max}},$$where $v_{max}$ and $v_{min}$ denote the maximum and initial degrees of freedom, which are typically set to 30 and 2, respectively; t is the current iteration number; $T_{max}$ is the total number of iterations.

This design enables the algorithm to maintain a broader search range in the early iteration stages while gradually converging in the later stages of the iteration process.

Case study and experimental analysis

Warehouse configuration and task parameter settings

A simulation platform was developed based on the real-world configuration of a semi-finished tobacco automated high-bay warehouse in a cigarette manufacturing plant to evaluate the effectiveness and generalizability of the proposed integrated optimization algorithm. The simulation environment was implemented in Python 3.9 and executed on a pc running Windows 10, equipped with an Intel(R) Core (TM) i7-14700HX 2.10 GHz CPU and 8 GB RAM.

The specific simulation setup was as follows:

1. Warehouse configuration parameters

The target system was a single-aisle warehouse with racks on both sides, each divided into three independent storage zones. The system comprised 29 columns and 11 levels per side, yielding 319 storage locations. The average storage utilization rate was set to 70%, and the initial warehouse state was generated randomly. The stacker crane operated at a horizontal speed of 3 m/s and a vertical speed of 0.7 m/s;

1. Generation parameters of inbound and outbound tasks

The total number of inbound and outbound tasks, N, was selected from the set {20, 40, 60, 80}. The outbound tasks were assigned to storage zones according to the following distribution: 30% to Zone A, 40% to Zone B, and 30% to Zone C. The inbound tasks followed a different allocation: 40% to Zone A, 30% to Zone B, and 30% to Zone C. The numbers of inbound and outbound tasks were kept symmetric to construct dual-command cycles;

1. Due-time modeling for outbound tasks

The due times of the outbound tasks were modeled following the approach proposed in [30], and the due-date assignment was formulated as follows:

33

$$D_o\sim u\left[min\{p_o|o\in {\mathcal{T}}^{\text{out}}\},\frac{2\times (1-\gamma )\times \sum _{o\in {\mathcal{T}}^{\mathit{out}}}{p}_o+min\{p_o|o\in {\mathcal{T}}^{\text{out}}\}}{z}\right],$$where $p_o$ denotes the DC operation time of a task o, which is defined by the original task execution sequence generated based on the slotting strategy described in Sect. 5.3; $\gamma \in \{0.6,0.7,0.8\}$ is the due-date tightness factor, where a smaller γ value corresponds to a looser deadline $U(·)$, which denotes a uniformly distributed random variable; $z$ is the number of stacker cranes, whose value is obtained by calculating the slot allocations using the initial sequence of the inbound and outbound tasks.

Algorithm parameter selection

To further enhance the solution performance of the improved dung beetle algorithm, this study used a full factorial experimental design to optimize the combination of the key control parameters. Compared to the orthogonal design, the full factorial design required more experimental runs. Still, it enabled a simultaneous estimation of the main effects and the relationships between the parameters, thereby overcoming the limitation of neglecting inter-factor influences. The selected core parameters included the population size, and the stability and convergence quality of the algorithm under different parameter settings were evaluated using the signal-to-noise ratio (S/N) value. The S/N ratio was defined as follows:

34

$$S/N=-10\times {log}_{10}(\frac{1}{n}\sum _{i=1}^n{y}_i^2),$$where n is the number of the algorithm’s runs under each parameter setting, and $y_i$ is the objective value obtained in the ith experiment.

Based on relevant literature [27], this study selected four key control parameters of the dung beetle algorithm for a full factorial experiment. Each parameter was configured at three levels (i.e., low, medium, and high levels) to capture the influence of different parameter gradients on the algorithm’s performance. Population size was set to (30, 50, 100), deflection coefficient k was set to (0.05, 0.1, 0.15), constant b was set to (0.1, 0.3,0.7), and constant S was set to (0.1, 0.5, 1.5). In the dung beetle optimization algorithm, the population comprised 20% rolling beetles, 20% reproducing beetles, 25% foraging beetles, and 35% stealing beetles. The due-date tightness factor $\gamma$ was tested using the values of 0.6, 0.7, and 0.8. By using the aforementioned fixed parameters and adjustable control factors, this study constructed a full factorial design comprising 3⁴(81) parameter combinations. However, due to the large number of combinations, the detailed table is not provided in the manuscript. The order size was set to 40, and each combination was independently executed ten times, with the average value used to analyze the impact of different parameter levels on the algorithm’s performance, as well as to determine the optimal parameter configuration for the subsequent experimental analysis (Table 3).

Table 3. Results of the full factorial experiments

No	Population size	Deflection coefficient k	Constant b	Constant S	Objective (Mean)
1	30	0.05	0.1	0.1	44.95
2	30	0.05	0.1	0.5	41.85
3	30	0.05	0.1	1.5	46.82
4	30	0.05	0.3	0.1	49.31
5	30	0.05	0.3	0.5	44.75
6	30	0.05	0.3	1.5	50.75
7	30	0.05	0.7	0.1	46.58
8	30	0.05	0.7	0.5	48.47
9	30	0.05	0.7	1.5	49.94

The experimental results of the signal-to-noise ratio are illustrated in Fig. 7. The S/N ratio reached the maximum for the population size of 50, k = 0.05, b = 0.1, and S = 0.5, indicating the best algorithmic performance under this configuration. Therefore, these parameter settings were adopted for all subsequent comparative experiments.

[See PDF for image]

Fig. 7

The results of S/N analysis

Comparison with heuristic algorithms

Before conducting the comparative experiments, the key parameters of the GA and the original DBO were tuned using a full factorial design to ensure a fair comparison. For the GA, following relevant literature [31], the crossover probability was set to $p_c\in \{0.6,0.7,0.8\}$, and the mutation probability was $p_m\in \{0.01,0.02,0.05,0.10,0.20\}$, yielding 15 parameter combinations. In addition, each combination was independently executed ten times, performing the same number of function evaluations, and the average value was recorded. For the DBO, the same parameter ranges were used as those in the IDBO, thereby ensuring that all comparative algorithms obtained reasonably optimized parameter configurations within a valid range. The S/N ratio results are shown in Fig. 8.

[See PDF for image]

Fig. 8

The main effects of the GA and DBO’s parameters on the S/N ratio

As shown in Fig. 8, optimal parameters of the GA were $p_c=0.7$ and $p_m=0.1$, and those of the DBO were k = 0.1, b = 0.3, and S = 0.5. Moreover, using the same experimental settings (e.g., the identical computational platform, a population size of 50, and 200 algorithm iterations), comparative experiments were conducted on the IDBO, GA, and DBO algorithms. Each algorithm was independently executed ten times, and the best solution was recorded. In the subsequent analysis, N denoted the task size, CT indicated the task completion time, Delay showed the task delay time, and Obj referred to the objective function value.

As presented in Table 4, the three algorithms exhibited different performances in terms of the CT, Delay, and Obj metrics, with the IDBO algorithm consistently achieving the best results across all scales. With the increase in the task size, the CT value increased for all the algorithms, but the IDBO algorithm consistently had a lower CT value than the other algorithms. For instance, for N = 160, the IDBO algorithm reaches 650.50, outperforming the GA (660.33) and DBO (668.55) algorithms. In terms of the Delay value, the IDBO algorithm performed significantly better than the other two algorithms, with an average of 144.87; the Delay values of the GA and DBO algorithms were 169.43 and 180.46, respectively. Moreover, the advantage of the IDBO algorithm became even more pronounced as the scale increased. Regarding the Obj metric, the IDBO algorithm achieved an average value of 237.94, thereby outperforming the GA (253.62) and DBO (261.83) algorithms, achieving relative improvements of 6.18% and 9.12%, respectively.

Table 4. Comparison results of different heuristic algorithms

N	$\gamma$	IDBO			GA			DBO
		CT	Delay	Obj	CT	Delay	Obj	CT	Delay	Obj
40	0.6	78.46	12.52	45.49	79.14	20.46	49.80	80.89	21.46	51.18
70	0.7	172.04	56.68	114.36	173.57	70.53	122.05	176.38	78.39	127.39
100	0.8	286.07	175.96	231.02	302.71	205.71	254.21	310.65	227.52	269.09
130	0.7	467.98	205.26	336.62	473.25	239.89	356.57	479.43	253.21	366.32
160	0.6	650.50	273.93	462.22	660.33	310.56	485.45	668.55	321.74	495.15
Average		331.01	144.87	237.94	337.80	169.43	253.62	343.18	180.46	261.83

The task size of $N\in \{40,60,80\}$ was selected for the test cases to verify the convergence performance of the proposed algorithm; the due-date tightness parameter $\gamma$ was set to 0.6, and the number of iterations was fixed at 200. The convergence behavior of the IDBO, GA, and DBO algorithms was comparatively analyzed under these conditions.

As shown in Fig. 9, the IDBO algorithm has significantly superior convergence speed and accuracy compared to the GA and DBO algorithms, reflecting its strong global search capability and good robustness. This advantage primarily originated from the balance between diversity and convergence maintained in the IDBO algorithm during the initialization and search phases. The improved Circle chaotic mapping enhanced the uniform distribution of the initial population in the search space, thus effectively avoiding the local initial-solution problem in conventional algorithms. Moreover, it improved the quality of the initial solutions and reduced convergence barriers in the early iterations. In addition, the introduction of the adaptive t-distribution perturbation strategy and the imperial competition mechanism increased the flexibility of the search process, allowing it to escape local optima and thereby ensured that the algorithm maintained a convergent trend while exploring new solution spaces. Particularly as the task scale increased, the IDBO algorithm consistently preserves a stable convergence trend, whereas the GA and DBO algorithms exhibited a marked decline in both the convergence speed and accuracy. This further validated the adaptability and optimization advantages of the IDBO algorithm in addressing complex problems.

[See PDF for image]

Fig. 9

Convergence curves of different algorithms for varying task sizes

Comparison with different strategy variants

Two representative baseline strategies reflecting the rule- and model-based scheduling characteristics in practical applications were selected to comprehensively evaluate the scheduling performance of the proposed method in scenarios with varying complexity. Specifically, the FCFS strategy was adopted as a rule-driven strategy, and the sequential slotting and picking (SSP) strategy served as a model-driven counterpart.

The FCFS technique is a commonly adopted warehouse storage policy in enterprise operations, where the heuristic rule allocates tasks strictly based on their arrival sequence, without introducing task prioritization or relying on a scheduling optimization model. The SSP approach is a heuristic scheduling framework for single-aisle AS/RSs under a random storage policy. This method constructs an IP model to minimize the DC operation time and achieve optimal matching between storage locations and inbound or outbound tasks. The system updated the shelf status after each iteration, reconstructed the model accordingly, and repeated the optimization cycle until all tasks were completed. This study used the task-to-zone constraints and employed a dynamic slotting mechanism based on the iterative updating mechanism of the original SSP framework to enhance the SSP adaptability. A comparative evaluation of the proposed algorithmic enhancements was conducted, and the experimental results are illustrated in Fig. 10.

[See PDF for image]

Fig. 10

Comparison results of the MHP, SSP, and FCFS scheduling strategies in terms of the task completion time, delay time, and objective function value: (a)–(c) $\gamma =0.6$; (d)–(f)$\gamma =0.7$; (g)–(i) $\gamma =0.8$

This study introduced the GAP metric to quantify the deviation of each strategy from the proposed method and improve the performance evaluation between different scheduling strategies. This metric was defined as follows:

35

$${\text{GAP}}_x=\frac{f_x-f_{\text{MHP}}}{f_{\text{MHP}}},x\in \{\text{FCFS},\text{SSP}\},$$where $f_x$ is the optimal solution obtained using a scheduling strategy $x$, and $f_{\text{MHP}}$ is the optimal solution achieved by the proposed method (Table 5).

Table 5. Comparison results of the proposed method and the SSP and FCFS algorithms

Metric	Average GAP (%) of SSP relative to the proposed method	Average GAP (%) of FCFS relative to the proposed method
Task completion time	− 12.33	380.14
Task delay time	82.73	304.60
Objective value	13.78	359.34

In the comparative experiments with two typical scheduling algorithms (i.e., the SSP and FCFS scheduling algorithms), the proposed MHP algorithm demonstrated significant advantages across most performance indicators. Specifically, compared to the SSP algorithm, the MHP algorithm achieved average improvements of 65.8% and 9.96% in the task delay time and objective function value, respectively. Although the MHP algorithm performed slightly worse in task completion time (with an average difference of − 7.27%) than the other algorithms, it exhibited greater overall scheduling quality and computational efficiency advantages. Compared to the traditional FCFS strategy, the MHP algorithm achieved significantly improved performance and stronger practical applicability. The GAP values of the FCFS algorithm for the task completion time, delay time, and objective function value were 359.69%, 259.17%, and 316.05%, respectively. These values indicated that the lack of effective task sequencing and resource coordination mechanisms severely limited the algorithm’s ability to optimize scheduling under complex constraints. These results verified the effectiveness and practical value of the proposed algorithm in optimizing warehouse task scheduling and storage slot allocation.

Next, comparative experiments were conducted across four task scales (i.e., 30, 60, 90, and 120) to validate the optimization performance of the proposed dynamic allocation strategy, setting parameter $\gamma$ to a fixed value of 0.6. The proposed approach was compared with the static-fixed and dynamic-nearest allocation strategies, as presented in Table 6.

Table 6. Comparison results of the objective function value of different allocation strategies

Task number	Objective function value (MHP)	Relative GAP (%)—static-fixed allocation	Relative GAP (%)—dynamic-nearest allocation
30	34.55	19.29	11.73
60	73.08	97.01	15.02
90	133.90	159.17	13.77
120	258.80	156.46	13.73

As shown in Table 6, the proposed MHP strategy consistently achieved the lowest objective function value for all task scales. Regarding the relative GAP, the difference between the static-fixed allocation strategy and the MHP strategy increased significantly with the task number, from 19.29 for 30 tasks to 159.17% for 90 tasks. The warehouse operation efficiency tended to stabilize after a task scale of 90, indicating a saturation point in system performance improvement. In contrast, the dynamic-nearest allocation strategy achieved relatively lower GAP values across the task scales, ranging from 11.73 to 15.02%, and then slightly decreased from 13.77 to 13.73%. Nevertheless, a noticeable performance gap remained compared to the proposed MHP strategy. The MHP strategy demonstrated superior optimization performance and scalability across the task scales, which highlighted its practical significance for real-world warehouse scheduling systems.

Model validation using visual components simulation platform

To verify the applicability of the proposed algorithm, this study constructed a simulation model of inbound and outbound operations in a high-bay house using the Visual Components platform. The operational performance of the proposed algorithm was comparatively analyzed using different strategies to validate the feasibility of the integrated scheduling model, as well as the applicability of the proposed improved algorithm.

Simulation model construction

This study selected an automated high-bay warehouse of a cigarette manufacturing plant as a case study and constructed a three-dimensional simulation model using the Visual Components platform. The warehouse layout consisted of two rows of racks and a single stacker crane, where each row of racks was divided into three independent storage zones. The number of inbound and outbound tasks was set to 40. Meanwhile, to simulate the uncertainty of on-site order arrivals, the task inter-arrival times were assumed to follow a truncated normal distribution with a mean of 9 s, a standard deviation of 2.5 s, and a truncation interval of [1, 23] s. This ensured a more realistic capture of the fluctuation characteristics of real operations. The specific simulation parameter settings are summarized in Table 7.

Table 7. Key parameter settings of the constructed Visual Components simulation model

Simulation parameter	Value	Unit
Cargo length, width, and height	1.2	$\text{m}$
The number of inbound tasks (A/B/C categories)	8/6/6	Boxes
The number of outbound tasks (A/B/C categories)	6/8/6	Boxes
The number of rack rows	29	Columns
The number of rack levels	11	Levels
The number of rack bays	2	Bays
Stacker crane horizontal/vertical running speed	3/1	$\text{m}·{\text{s}}^{-1}$
Stacker crane horizontal/vertical acceleration	2/1	$\text{m}·{\text{s}}^{-2}$
Stacker crane fork extension/retraction time	5	$\text{s}$

In the constructed Visual Components simulation model, the order tasks, optimized by the MPH, were imported into the order input module. The order information is presented in Table 8. The simulation runs were conducted for both the original order sequence (before optimization) and the MPH-optimized sequence, and the arrival time of the first task was defined as the simulation start time.

Table 8. Inbound and outbound task schedule

Task No.	Cargo category	Target location (bay, level, column)	Task type	Task No.	Cargo category	Target location (bay, level, column,)	Task type
1	A	(2, 1, 1)	Inbound	21	A	(2, 4, 1)	Inbound
2	A	(1, 1, 1)	Outbound	22	B	(1, 7, 1)	Outbound
3	A	(1, 2, 1)	Inbound	23	B	(2, 8, 2)	Inbound
4	A	(2, 2, 1)	Outbound	24	C	(2, 10, 2)	Outbound
5	A	(2, 5, 2)	Inbound	25	C	(2, 10, 2)	Inbound
6	C	(1, 10, 1)	Outbound	26	B	(1, 7, 2)	Outbound
7	B	(1, 9, 2)	Inbound	27	A	(1, 3, 1)	Inbound
8	C	(1, 10, 2)	Outbound	28	A	(1, 4, 1)	Outbound
9	C	(1, 10, 1)	Inbound	29	B	(1, 9, 3)	Inbound
10	C	(2, 10, 1)	Outbound	30	C	(2, 10, 3)	Outbound
11	B	(2, 7, 1)	Inbound	31	A	(1, 4, 1)	Inbound
12	B	(1, 6, 1)	Outbound	32	A	(1, 5, 1)	Outbound
13	C	(1, 10, 2)	Inbound	33	C	(2, 10, 1)	Inbound
14	A	(1, 3, 1)	Outbound	34	C	(1, 11, 1)	Outbound
15	C	(1, 10, 3)	Inbound	35	A	(1, 5, 1)	Inbound
16	B	(1, 6, 2)	Outbound	36	A	(2, 5, 1)	Outbound
17	B	(1, 6, 1)	Inbound	37	C	(1, 11, 1)	Inbound
18	B	(2, 6, 1)	Outbound	38	B	(2, 7, 2)	Outbound
19	A	(1, 4, 2)	Inbound	39	B	(1, 7, 1)	Inbound
20	B	(2, 6, 2)	Outbound	40	B	(1, 8, 2)	Outbound

Simulation results analysis

The total operation time of different strategies in the experimental simulation of 40 tasks is shown in Fig. 11, where it can be seen that the optimization strategy demonstrated significant advantages in all comparisons. The total operation times of the four baseline strategies, namely the nearest inbound-random outbound, random inbound-random outbound, random inbound-nearest outbound, and nearest inbound-nearest outbound strategies, were all within the range of 1934–1960 s, with differences of less than 1.5%. These results indicated that relying only on the heuristic rules based on the nearest or random selection approach provided limited improvement in the overall efficiency and might cause the problems of local congestion and redundant routing. By contrast, the proposed strategy effectively reduced the total operation time to 1601 s, achieving an average reduction of approximately 345 s compared to the baseline methods, with an optimization margin of 17–18%. This verified the effectiveness and applicability of the MPH-based dynamic scheduling method in improving warehouse operational efficiency and system stability.

[See PDF for image]

Fig. 11

Comparison results of the task completion time of different strategies

Conclusion and future work

In this paper, task scheduling and slotting in semi-finished product high-bay automated warehouses within the cigarette manufacturing process are conducted. The problems of low operational efficiency and delayed responsiveness in traditional warehouse scheduling strategies under complex constraints are effectively addressed. In addition, a joint optimization model is developed combining the ILP and DP approaches to minimize the weighted sum of task completion time and order delay. Moreover, a bi-level collaborative optimization method that combines mathematical programming and heuristic algorithms is designed. In the proposed method, the outer layer performs global task sequencing, and the inner layer conducts local storage location optimization. The simulation results indicate that, across different task scales, the proposed algorithm consistently outperforms the GA and DBO algorithms in terms of performance. Based on the results, the proposed method can improve the objective function value by 13.78% and 359.34% compared to the SSP and FCFS scheduling strategies, respectively, ensuring optimal performance across all task scales. Finally, the proposed model is validated by tests on the Visual Components simulation platform, and the results demonstrate that the proposed method can reduce the total operation time by 17–18% compared with the four baseline strategies.

The AS/RS operation involves task scheduling, slotting, and routing processes of RGVs. In future work, the configuration and path planning of multiple RGVs could be considered, particularly in inbound and outbound staging areas. In addition, advanced models incorporating these factors could be developed to enhance the overall operational efficiency of warehouse systems. Furthermore, for multiple RGVs deployed on a single track, potential conflicts could be considered to ensure safe and efficient operations. Addressing the aforementioned challenges would help enhance the robustness of scheduling systems in complex and dynamic warehouse environments.

However, although the proposed heuristic algorithm demonstrates robust stability and high-quality solutions across a range of scheduling scenarios, it is dependent on manually predefined rules and heuristic logic. Therefore, it lacks a mechanism for adaptively responding to dynamic changes in the environment. Thus, in practical environments characterized by continuous changes in the task scale, job characteristics, and resource availability, the proposed algorithm can exhibit a relatively static adjustment strategy. The proposed algorithm struggles to effectively leverage informative patterns embedded in historical scheduling processes to guide future decision-making. Accordingly, future work could consider introducing reinforcement learning and other advanced technologies into online learning and policy-optimization capabilities, enabling the proposed algorithm to continuously adjust its decision-making strategies based on actual operational conditions. This would further enhance the proposed algorithm’s adaptability and intelligence in dynamic scheduling environments.

Acknowledgements

This research was supported by the Ningbo University “double world-class project” cooperation special directional entrusted scientific and technological cooperation projects (HX2025000080, HX2024000404), and the Ningbo Municipal National Natural Science Foundation (2023J113).

Funding

Ningbo University “double world-class project” cooperation special directional entrusted scientific and technological cooperation projects, HX2025000080, Rui Wang, HX2024000404, Rui Wang, Natural Science Foundation of Ningbo Municipality, 2023J113, Rui Wang.

Data availability

The data are available from the corresponding author on reasonable request.

Declarations

Conflict of interests

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Emilio C. N. Silva

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

References

1. Guchhait, R; Sarkar, B. A decision-making problem for product outsourcing with flexible production under a global supply chain management. Int J Prod Econ; 2024; 272, 109230. [DOI: https://dx.doi.org/10.1016/j.ijpe.2024.109230]

2. Shi D, Dong G, Chen E, et al (2024) Optimization of storage paths for finished cigarette logistics distribution based on improved GA-A. Informatica. 48(18)

3. Zhang, W; Deng, Z; Zhang, C et al. Integrated optimization of storage space allocation and crane scheduling in automated storage and retrieval systems. Robot Comp Integr Manuf; 2025; 93, 102918. [DOI: https://dx.doi.org/10.1016/j.rcim.2024.102918]

4. Duan Q, Jin C, Zhang S et al (2025) Intelligent inventory model for cigarette automated storage systems based on 5G IoT technology/2025 Asia-Europe conference on cybersecurity, internet of things and soft computing (CITSC). IEEE. pp 267–271

5. Gao, J; Xie, W; Shi, D et al. Synchronized optimization of the logistics system of a tobacco high-bay warehouse under production task fluctuations. Eng Res Express; 2024; 6, 4 045580. [DOI: https://dx.doi.org/10.1088/2631-8695/ad9f24]

6. Xu, X; Ren, C. Research on dynamic storage location assignment of picker-to-parts picking systems under traversing routing method. Complexity; 2020; 2020, 1 1621828.

7. Silva, A; Coelho, LC; Darvish, M et al. Integrating storage location and order picking problems in warehouse planning. Transp Res Part E: Log Transp Rev; 2020; 140, 102003. [DOI: https://dx.doi.org/10.1016/j.tre.2020.102003]

8. Bolanos Zuniga, J; Saucedo Martinez, JA; Salais Fierro, TE et al. Optimization of the storage location assignment and the picker-routing problem by using mathematical programming. Appl Sci; 2020; 10, 2 534. [DOI: https://dx.doi.org/10.3390/app10020534]

9. Winkelmann, D; Tolkmitt, F; Ulrich, M et al. Integrated storage assignment for an E-grocery fulfilment centre: accounting for day-of-week demand patterns. Flex Serv Manuf J; 2024; [DOI: https://dx.doi.org/10.1007/s10696-024-09549-7]

10. Mokarrari, KR; Sowlati, T; English, J et al. Optimization of warehouse picking to maximize the picked orders considering practical aspects. Appl Math Model; 2025; 137, 115585.4801413 [DOI: https://dx.doi.org/10.1016/j.apm.2024.06.037]

11. Balogh, A; Garraffa, M; O’Sullivan, B et al. MILP-based local search procedures for minimizing total tardiness in the no-idle permutation flowshop problem. Comput Oper Res; 2022; 146, 105862.4436369 [DOI: https://dx.doi.org/10.1016/j.cor.2022.105862]

12. Yao, S; Guo, Y; Yang, B et al. Single-objective flexible job-shop scheduling problem based on improved dung beetle optimization. STEM Educ; 2024; 4, 3 pp. 299-327. [DOI: https://dx.doi.org/10.3934/steme.2024018]

13. Falq, AE; Fouilhoux, P; Kedad-Sidhoum, S. Mixed integer formulations using natural variables for single machine scheduling around a common due date. Discrete Appl Math; 2021; 290, pp. 36-59.4186609 [DOI: https://dx.doi.org/10.1016/j.dam.2020.08.033]

14. Lei, C; Jiangfeng, C; Yonghuai, Z et al. Multi-task textile stereoscopic warehouse location allocation and task sequence optimization method. Comput Integr Manuf Syst; 2023; 29, 4 pp. 1371-1385.

15. Tang, HT; Yan, WJ; Chen, QF et al. Integrated optimization of location assignment and job scheduling in automated storage and retrieval system. Comp Sci; 2020; 47, 5 pp. 204-211.

16. Liu, Z; Lu, J; Ren, C et al. Joint optimization of storage assignment and order batching in robotic mobile fulfillment system with dynamic storage depth and surplus items. Computers Ind Eng; 2025; 200, 110767. [DOI: https://dx.doi.org/10.1016/j.cie.2024.110767]

17. He, L; Tao, Y; Luo, J et al. Job integrated optimization of automated storage/retrieval systems based on two-stage wolf pack algorithm. Chin Mech Eng; 2022; 33, 21 2538.

18. Rizqi, ZU; Chou, SY; Khairunisa, A. Energy harvesting for automated storage and retrieval system with sustainable configuration of storage assignment and input/output point. Transp Res Part E: Log Transp Rev; 2024; 192, 103781. [DOI: https://dx.doi.org/10.1016/j.tre.2024.103781]

19. Ariyanti, FD; Paramaputra, BE EDP Sciences. Reduce overtime of distribution centre by re-layout and employee shift scheduling use class based storage and integer linear programming. E3S Web Conf; 2023; 426, 01060. [DOI: https://dx.doi.org/10.1051/e3sconf/202342601060]

20. Wang, IL; Wang, TH. Efficient routing in robotic movable fulfillment systems with integer programming: a rolling horizon and heuristic approach. Robot Comp Integr Manuf; 2025; 91, 102849. [DOI: https://dx.doi.org/10.1016/j.rcim.2024.102849]

21. Ekren, BY; Akpunar, A. An open queuing network-based tool for performance estimations in a shuttle-based storage and retrieval system. Appl Math Model; 2021; 89, pp. 1678-1695.4148240 [DOI: https://dx.doi.org/10.1016/j.apm.2020.07.055]

22. Wang, W; Wu, Y; Zheng, J et al. A comprehensive framework for the design of modular robotic mobile fulfillment systems. IEEE Access; 2020; 8, pp. 13259-13269. [DOI: https://dx.doi.org/10.1109/ACCESS.2020.2966403]

23. Annear, LM; Akhavan-Tabatabaei, R; Schmid, V. Dynamic assignment of a multi-skilled workforce in job shops: an approximate dynamic programming approach. Eur J Oper Res; 2023; 306, 3 pp. 1109-1125.4531503 [DOI: https://dx.doi.org/10.1016/j.ejor.2022.08.049]

24. Koulamas, C; Kyparisis, GJ. A classification of dynamic programming formulations for offline deterministic single-machine scheduling problems. Eur J Oper Res; 2023; 305, 3 pp. 999-1017.4514749 [DOI: https://dx.doi.org/10.1016/j.ejor.2022.03.043]

25. Murad, SA; Muzahid, AJM; Azmi, ZRM et al. A review on job scheduling technique in cloud computing and priority rule based intelligent framework. J King Saud Univ Comp Inf Sci; 2022; 34, 6 pp. 2309-2331. [DOI: https://dx.doi.org/10.1016/j.jksuci.2022.03.027]

26. Gülmez, E; Koruca, HI; Aydin, ME et al. Heuristic and swarm intelligence algorithms for work-life balance problem. Comput Ind Eng; 2024; 187, 109857. [DOI: https://dx.doi.org/10.1016/j.cie.2023.109857]

27. Xue, J; Shen, B. Dung beetle optimizer: a new meta-heuristic algorithm for global optimization. J Supercomput; 2023; 79, 7 pp. 7305-7336. [DOI: https://dx.doi.org/10.1007/s11227-022-04959-6]

28. Lucyshyn, T; Des Enffans d’Avernas, LV; Holzer, C. Influence of the mold material on the injection molding cycle time and warpage depending on the polymer processed. Polymers; 2021; 13, 18 3196. [DOI: https://dx.doi.org/10.3390/polym13183196]

29. Atashpaz-Gargari E, Lucas C (2007) Imperialist competitive algorithm: an algorithm for optimization inspired by imperialistic competition[C]//2007 IEEE congress on evolutionary computation. Ieee. pp 4661–4667

30. Scholz, A; Schubert, D; Wäscher, G. Order picking with multiple pickers and due dates–simultaneous solution of order batching, batch assignment and sequencing, and picker routing problems. Eur J Oper Res; 2017; 263, 2 pp. 461-478.3669045 [DOI: https://dx.doi.org/10.1016/j.ejor.2017.04.038]

31. Samsuria, E; Mahmud, MSA; Wahab, NA et al. An improved adaptive fuzzy-genetic algorithm based on local search for integrated production and mobile robot scheduling in job-shop flexible manufacturing system. Computers Ind Eng; 2025; 204, 111093. [DOI: https://dx.doi.org/10.1016/j.cie.2025.111093]