Quantum computing for several AGV scheduling models

TL;DR

This paper explores using quantum computing to improve scheduling of automated guided vehicles, achieving significant time savings.

Contribution

The paper introduces quantum computing for AGV scheduling and proposes two new QUBO models solved with an optical quantum computer.

Findings

Quantum computing reduces AGV scheduling computation time by 92% on average.

Two QUBO models are developed for different scheduling objectives.

Optical coherent Ising machine effectively solves AGV scheduling problems.

Abstract

Due to the high degree of automation, automated guided vehicles (AGVs) have been widely used in many scenarios for transportation, and traditional computing power is stretched in large-scale AGV scheduling. In recent years, quantum computing has shown incomparable performance advantages in solving specific problems, especially Combinatorial optimization problem. In this paper, quantum computing technology is introduced into the study of the AGV scheduling problem. Additionally two types of quadratic unconstrained binary optimisation (QUBO) models suitable for different scheduling objectives are constructed, and the scheduling scheme is coded into the ground state of Hamiltonian operator, and the problem is solved by using optical coherent Ising machine (CIM). The experimental results show that compared with the traditional calculation method, the optical quantum computer can save 92% computation time on average. It has great application potential.

Figures

The AGV scheduling problem and a feasible solution. All AGVs start from a fixed start node, perform transportation tasks, and reach the end node after performing all tasks. ’S’ represents the starting point of a transportation task, and ’E’ represents the end point of a transportation task. Different colors represent different AGV’s mission routes.

Schematic diagram of coherent Ising mechanism construction and principle.

Computation time versus number of tasks for MIP model. (a) Represents the change of MIP model computation time with the number of tasks under the goal of minimizing the total travel time. (b) Represents the change of MIP model computation time with the number of tasks under the goal of minimizing the makespan.

Evolution diagram of Hamiltonian with time under the objective functions of minimizing the total travel time in node model. (a) Represents the evolution diagram of Hamiltonian with time under the example of 4 tasks. (b) Represents the evolution diagram of Hamiltonian with time under the example of 5 tasks. (c) Represents the evolution diagram of Hamiltonian with time under the example of 6 tasks. (d) Represents the evolution diagram of Hamiltonian with time under the example of 7 tasks.Figure 5Evolution diagram of Hamiltonian with time under two objective functions of arc mode. (a) represents

Schematic diagram of quantum computing solutions for node model. (a) represents solution under 4 tasks. (b) represents solution under 5 tasks. (c) represents solution under 6 tasks. (d) represents solution under 7 tasks.Figure 7Schematic diagram of quantum computing solutions for arc model. (a) represents the solution of 4 task under the objective function of minimizing the total time. (b) represents the solution of 4 task under the objective function of minimizing the makespan.

Full text

Introduction

The scale of the logistics industry has maintained a considerable growth rate, and many human and material resources have been invested in it, thus creating the labor-intensive industry. Improving the automatic and intelligent level of the logistics industry has become an important issue for industry and academia. In recent years, some of the industry’s leading enterprises have already carried out technological reforms. For example, the retail giant, Amazon, has established a huge logistics center, that uses intelligent sorting technology, delivery drones, automated guided vehicles (AGVs), etc. China’s e-commerce giant, Jingdong, has also set up an ‘Asia One’ warehouse, in which more than 100 AGVs are used for transportation operations at the same time^1^. In addition, technical innovation also occurred in ports. In 1993, the world’s first automated wharf was built in the Amsterdam Port in the Netherlands, and dozens of AGVs were used for container transshipment. With the introduction of technology and the accumulation of operation experience, the construction of automatic terminals has been tried all over the world, and the usage of AGVs has gradually increased. Currently, the use of AGVs has penetrated all aspects of logistics, transportation and production, which has greatly promoted the level of industrial automation and intelligence and improved efficiency.

To meet the needs of application scenarios, the amount of parallel work of AGVs is increasing, which brings great difficulty to the AGV scheduling. AGV scheduling problem is a difficult combinatorial optimization problem, and a large number of researchers have devoted themselves to this field and made some contributions. Singh et al.^2^ considered AGV scheduling with battery constraints, developed a mixed-integer programming model with the objective of minimizing the combined task delay cost and AGV transportation cost, and designed a customized adaptive large-neighborhood search algorithm to solve the model. Zhang et al.^3^ studied an AGV scheduling problem in matrix manufacturing plants, and proposed a mixed integer programming model to minimize the generalized transportation cost, based on which an improved iterative greedy algorithm was designed and compared with six other algorithms to show its superior solution performance. For the scheduling problem of AGVs in smart factories, Zhang et al.^4^ proposed a self-organized dynamic scheduling method, that groups multiple AGVs to perform tasks among themselves and uses improved gene expression programming to learn dynamic scheduling rules. The numerical experimental results show that the method can considerably reduce system costs. Wang and Zeng^5^ studied the port AGV scheduling and path planning problem under conflict-free paths, established a mixed integer model with the objective of minimizing task completion time, and proposed a customized branch and bound algorithm combined with a heuristic algorithm to solve the small-scale problem, and further developed a two-stage greedy heuristic algorithm to quickly obtain a satisfactory solution for the large-scale problem. Sagar and Jerald^6^ proposed a real-time scheduling strategy for AGVs based on deep reinforcement learning technology, established a Markov decision model for real-time scheduling, and developed a Q-learning algorithm. The superiority of the method is shown through numerical experiments. Considering the scheduling and path planning model of shop floor AGVs, Saidi et al.^7^ developed a discrete-time model and proposed a two-stage ant colony algorithm to solve the model. From the above literature, the research on scheduling problems of AGVs covers several scenarios such as workshops and terminals. Researchers have built mixed integer programming models, integer programming models, Markov decision process models, etc. The methods used are scheduling rules, exact algorithms, heuristic algorithms, reinforcement learning algorithms, etc. From the results, it is observed that the exact algorithm can generate optimal scheduling solutions, however, its computational time is prohibitively slow, rendering it impractical for large-scale problems. Inexact algorithms exhibit favorable efficiency but often converge to local optima. The provision of high-quality scheduling solutions within a short timeframe poses a significant challenge.

In recent years, significant advancements have been made in both theoretical understanding and practical applications of quantum computing. The fundamental distinction between quantum computers and traditional computers lies in their reliance on quantum mechanical principles. Quantum computers utilize quantum bits (qubits) as the fundamental units of information storage^8^, which can exist in superposition states of both 1 and 0, enabling them to hold exponentially more information compared to traditional computers. It is well recognized that quantum computers offer substantial advantages, particularly in addressing specific problems such as combinatorial optimization, often described as the superiority of quantum computing. Many combinatorial optimization problems are NP-hard, presenting significant challenges for traditional computers to solve. Combinatorial optimization problem can be mapped to the ground state search problem of Ising model. Hardware systems can be built in many different ways to simulate the process of Hamiltonian reduction, such as adiabatic quantum computing (AQC), quantum annealing (QA), etc. However, it is always a difficult problem to improve the connection density between qubits, which will affect the efficiency of problem solving^9,10^. Coherent Ising Machine (CIM) is a quantum computer developed according to the optical principle^11–17^, which can work at room temperature and deal with large-scale problems, such as compression sensor problems^18^ and polyhedron problems^19^. CIM uses laser pulses in optical fiber as qubits for quantum calculation. The early prototype of CIM is injection synchronous laser Ising machine. The number of coupled lasers in this scheme is proportional to the square of qubits, which is quite difficult. On this basis, optical delay linear CIM and measurement feedback CIM using nonlinear optical crystal instead of laser are developed. The latter uses measurement feedback to avoid the challenge that the former needs to control a large number of optical delay lines accurately^14^. The machine used in this study is measurement feedback CIM.

AGV scheduling problem can be understood as a kind of routing problem. Most traditional solutions to routing problems require sacrificing large amounts of computational resources and Osaba et al.^20^ indicated that quantum computing techniques have great potential in the area of solving routing and optimization problems. In the early days, Goswami et al.^21^ developed a phase estimation technique to solve the traveling salesman problem (TSP), using IBM’s quantum simulator to provide results for four city cases. Then, researchers tried to solve more complex problems with quantum computing. Feld et al.^22^ presented a quadratic unconstrained binary optimization (QUBO) formulation for solving the vehicle routing problem with capacity constraints, evaluated the solution quality and computation time and compared it with classical solution methods. Bao et al.^23^ proposed a two-stage QUBO formulation of the vehicle routing problem with balanced pickup, mapping the first stage to a clustering problem and describing the second stage as a TSP problem, and evaluated it against traditional methods in terms of numerical experimental results. Harwood et al.^24^ tried to establish a qubo model to describe the vehicle routing problem by using the modeling idea of node and arc, and evaluated the model by using analog quantum devices. Geitz et al.^25^ built a QUBO model to solve the job-shop scheduling problem, using quantum computers or simulators, constrained programming and tabu search. The calculation results proved the effectiveness of quantum computing in small-scale situations. And the established QUBO model can be extended to AGV scheduling problem. Ohzeki et al.^26^ formulated an Ising model for the collision-free scheduling problem of AGVs within a factory setting. They utilized a quantum annealing machine to solve the model, with results demonstrating the potential application of quantum annealing machines in addressing real-world industrial challenges. Based on the above cases, it can be seen that some researchers have begun to use quantum computing to solve practical problems in the field of optimization. However, the research on quantum computing related to AGV scheduling has just started, and many researchers used simulators to solve them, because the current physical real machine resources are scarce, and the scale of solving problems is still relatively small, and it is easy to make mistakes and lacks the running data of physical real machines.

Our contributions

In a word, most of the existing AGV scheduling research adopts traditional models and methods, which can not effectively meet the actual needs of large-scale scheduling. Quantum computing has great application potential in solving specific problems that traditional computers cannot solve, and researchers have tried in optimization fields. However, as far as we know, there are few literatures about AGV scheduling using quantum computing technology. Based on these facts, the idea of carrying out this research came into being. The main contributions of this paper are summarized as follows.

1. In traditional research on the AGV scheduling problem, the computation time increases greatly with an increase in the number of AGVs and tasks. We introduce quantum computing technology into the research of the AGV scheduling problem and construct new QUBO models of AGV scheduling. In real scenarios, dispatchers often set different scheduling objectives according to the nature of the work, among which minimizing the total AGV travel time and minimizing the task completion time (makespan) are the two most common objectives. According to the different objectives, we have deduced different QUBO models, and given the model solutions and related theoretical basis under two different objectives.
2. We use traditional computer and CIM to carry out numerical experiments on the traditional model and QUBO model proposed by us respectively. The experimental results show that the computation speed of CIM is much faster than that of traditional computer, and the average calculation time is saved by 92%, which proves that CIM has great application potential in solving AGV scheduling problem and similar combinatorial optimization problems.

AGV scheduling model

AGV scheduling problems have many classifications according to different scenarios and considerations. For example, consider the time window of the task, joint optimization of scheduling and path, cooperation with other devices, charging strategy and so on. Due to the limitation of quantum bits of CIM, it is impossible to solve the AGV scheduling problem in complex scenes^27–29^. Therefore, we simplify the problem and keep the essence of AGV scheduling problem. On this basis, we construct the AGV scheduling model. In this section, we present the classical AGV scheduling model based on mixed integer programming (MIP), and propose two new models, which we call the node and arc models.

Problem description

Figure 1. The AGV scheduling problem and a feasible solution. All AGVs start from a fixed start node, perform transportation tasks, and reach the end node after performing all tasks. ’S’ represents the starting point of a transportation task, and ’E’ represents the end point of a transportation task. Different colors represent different AGV’s mission routes.

We consider an AGV scheduling problem (to make the problem more general, we do not set up a working scenario), as shown in Fig. 1. Given an AGV set, all AGVs have a unified starting node and ending node, and all AGVs start to accept tasks from the unified starting node until all transportation tasks are completed, and then return to the unified ending node. In the AGV scheduling problem, we are given a set of transportation tasks, each with a starting point and an ending point. For an AGV, the process of completing the transportation task can be described as first arriving at the starting point of the task to load the transported goods, then transporting them to the ending point of the task, and then driving to the starting point of the next task to perform the next task. The travel time of the AGV between any two task points is known. We consider two optimization objectives, the first is to minimize the total AGV travel time and the second is to minimize the maximum task completion time (makespan). The first objective is generally used when the task is not urgent to achieve a reduction in total system energy consumption, while the second objective is set to complete transportation tasks quickly.

Next, we elaborate on the symbolic settings in the problem as follows.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$V=\{1,\ldots ,k,\ldots ,K \}$$\end{document} set of AGVs,

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ R=\{1,\ldots ,r,{{r}^{'}},{{r}^{''}},\ldots ,n-1 \}$$\end{document} set of actual tasks,

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{R}_{1}}=\{0,1,\ldots ,r,{{r}^{'}},{{r}^{''}},\ldots ,n-1\}$$\end{document} set of actual tasks and virtual start task,

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{R}_{2}}=\{1,\ldots ,r,{{r}^{'}},{{r}^{''}},\ldots ,n\}$$\end{document} set of actual tasks and virtual end task,

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{R}^{'}}=\{0,1,\ldots ,r,{{r}^{'}},{{r}^{''}},\ldots ,n\}$$\end{document}  set of all tasks,

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$A=\{(r,{{r}^{'}})\mid r,{{r}^{'}}\in {{R}^{'}}\}$$\end{document} arc set that consists of all valid task pairs that can be conducted adjacently,

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\pi =\{1,\ldots ,t,\ldots ,N \}$$\end{document} set of the sequence of tasks performed by an AGV

a        a single task arc

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _{r}^{+}$$\end{document}      an arc with task r as the left node

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\theta _{r}^{-}$$\end{document}      an arc with task r as the right node

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$r_{\text{s}}$$\end{document}       starting point of task r,

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$r_{\text{d}}$$\end{document}       ending point of task r.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{e}_{{{r}_\text{s}}{{r}_\text{d}}}}$$\end{document}     indicates the travel time of the AGV from two points \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{r}_\text{s}}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{r}_\text{d}}$$\end{document} ,

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{e}_{{{r}_\text{d}}r_\text{s}^{'}}}$$\end{document}     indicates the travel time of the AGV from two points \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{r}_\text{d}}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$r_\text{s}^{'}$$\end{document} ,

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{c}_{r{{r}^{'}}}}$$\end{document}      contains two parts of time, the first part is the time from the start of task r to the end of task r, and the second part is the time from the end of task r to the start of task \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${r}^{'}$$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{c}_{r{{r}^{'}}}}={{e}_{{{r}_{s}}{{r}_{d}}}}+{{e}_{{{r}_{d}}r_{s}^{'}}}$$\end{document} .

Mixed integer programming model

In this subsection, we introduce the classical model of AGV scheduling. The classical model is formed as a mixed integer programming, and it has the following variables.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{y}_{r,{{r}^{'}},k}}$$\end{document}       binary variable, equal to 1 if task r is performed directly prior to \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{r}^{'}}$$\end{document} by AGV k, 0 otherwise;

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$f^\text{s}_{rk}$$\end{document}          arrival time of AGV k at the start node of task r;

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$f^\text{d}_{rk}$$\end{document}          arrival time of AGV k at the end node of task r;

T           the makespan for AGVs to perform transportation tasks.

The first optimization objective of the MIP model is to minimize the total AGV travel time, and its model is presented as follows.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \text {(P1)}\qquad&\text {min}{} & {} \sum _{k\in V} f_{nk}^\text{d} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}&\text {s.t.}{} & {} \sum _{r'\in R_2} y_{0,r',k}=1&\forall k\in V \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}{} & {} {}&\sum _{r\in R_1}y_{r,n,k}=1&\forall k\in V \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}{} & {} {}&\sum _{r'\in R_2}\sum _{k\in V} y_{r,r',k}=1&\forall r\in R \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}{} & {} {}&\sum _{r'\in R_1}y_{r',r,k}-\sum _{r'\in R_2} y_{r,r',k}=0&\forall r\in R,\forall k\in V \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}{} & {} {}&f_{0k}^\text{s}=f_{0k}^\text{d}=0&\forall k\in V \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}{} & {} {}&f_{rk}^\text{s}+e_{r_\text{s}r_\text{d}}=f_{rk}^\text{d}&\forall k\in V, \forall r\in R' \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}{} & {} {}&f_{rk}^\text{d}+e_{r_\text{d}r'_\text{s}}\le f_{r'k}^\text{s}+M(1-y_{r,r',k})&\forall r,r'\in R,\forall k\in V \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}{} & {} {}&y_{r,r,k}=0&\forall r\in R,\forall k\in V \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}{} & {} {}&{{y}_{r,{{r}^{'}},k}}\in \{0,1\}\text { }&\forall {r,}{{{r}}^{'}}\in {{R}^{'}},\forall k\in V \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}{} & {} {}&f_{rk}^{s},f_{rk}^{d}\ge 0\text { }&\forall {r}\in {R,}\forall {k}\in {V}. \end{aligned}$$\end{document}

The objective function (1) is to minimize the total travel time of the AGV. Constraints (2) and (3) ensure that all AGVs need to complete the virtual start task and the virtual end task, and constraints (4) guarantee that all actual tasks are uniquely assigned to a particular AGV. Constraints (5) ensure that each AGV completes its task satisfying the flow balance. Then, constraints (6) guarantee that the virtual start task starts and ends at moment 0. Constraints (7) states that the time to reach the end of a task is equal to the time to reach the start of that task plus the transport time from the start to the end. Constraints (8) indicates that the time to reach the start of a task is later than the time to reach the end of the previous task plus the transportation time required to reach the start of that task from the end of the previous task, where M is a sufficiently large value. The last constraints (9) eliminates task self-citation. Constraints (10) and (11) represent range limits of variables.

A slight modification of the above model can be used as a model for minimizing the makespan for AGVs to perform transportation tasks, which is given as follows, Eqs. (2)–(11).

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \text {(P2)}\qquad&\text {min}{} & {} T \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}&\text {s.t.}{} & {} T\ge f_{nk}^{d}&\forall k\in V. \end{aligned}$$\end{document}

where T denotes the makespan. The objective function is to minimize T, and constraints (13) denotes that T must be no less than the time required by the last AGV to complete the task.

QUBO model

QUBO and Ising model

This section mainly describes the concepts of QUBO model and Ising model and their relationship. QUBO is an expression of optimization problem, and its goal is to find binary variables that minimize quadratic polynomials. Ising model was first put forward and applied in statistical physics. It describes a system composed of interacting units, in which each spin particle must have two possible random states (such as + 1 and − 1), and then it was introduced into the field of mathematics as a model to describe a series of optimization problems. Many combinatorial optimization problems can be expressed in the form of quadratic unconstrained binary optimization or ising model, and they can be transformed into each other^30–32^. The general expression of QUBO model is shown in Eq. (14).

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} H(x)={{x}^\texttt{T}}Qx+{{c}^\texttt{T}}x,x\in {{\{0,1\}}^{z}} \end{aligned}$$\end{document}

where x is a z dimensional vector of binary variables, Q is the quadratic coefficient matrix, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{c}^\texttt{T}}$$\end{document} is the coefficient matrix of the primary term.

The above model in the form of QUBO can be easily transformed into an Ising model, and the variable range of the Ising model is \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\left\{ -1,1 \right\} $$\end{document} . Specifically, it can be realized by variable substitution \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\sigma }_{i}}=2{{x}_{i}}-1$$\end{document} . Then, the optimization function can be expressed in the following form.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} H(\sigma )=-\sum \limits _{i,j}{{{J}_{ij}}}{{\sigma }_{i}}{{\sigma }_{j}}-\sum \limits _{i}{{{h}_{i}}}{{\sigma }_{i}} \end{aligned}$$\end{document}

where the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\sigma }_{i}}$$\end{document} is spin variable, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{J}_{ij}}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{h}_{i}}$$\end{document} are the quadratic and linear coefficients. The solution of Ising problem is to find the ground state of Hamiltonian. CIM solves the Ising problem according to the principle of minimum gain, and can find the ground state or low energy state of Ising Hamiltonian. The method is to map the QUBO problem into a fully connected Ising Hamiltonian with programmable parameters, and obtain the solution of the problem through controllable quantum phase transition^33,34^.

Node model

In this section, we will describe the node model. The core idea of node model is to regard tasks as nodes and the order of task execution as the order of vehicles passing through nodes.The node model has the QUBO form and is suitable for quantum computing, the variables of the model are described as below.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{x}_{r,t,k}}$$\end{document} binary variable, equal to 1 if task r is assigned to AGV k as it is t-th task, 0 otherwise.

The first optimization objective of the node model is to minimize the total AGV travel time, and the model is shown below.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \text {(P3)}\qquad {{H}_{1}}&={{\partial }_{1}}{{H}_{A}}+{{\partial }_{2}}{{H}_{B}}+{{\partial }_{3}}{{H}_{C}}+{{\partial }_{4}}{{H}_{D}}+{{\partial }_{5}}{{H}_{E}}+{{\partial }_{6}}{{H}_{F}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{A}}&=\sum \limits _{k\in V}{\sum \limits _{(r,{{r}^{'}})\in A}{\sum \limits _{t=1}^{N-1}{{{c}_{r{{r}^{'}}}}}}}{{x}_{r,t,k}}{{x}_{{{r}^{'}},t+1,k}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{B}}&={{\sum \limits _{k\in V}{(1-{{x}_{0,1,k}})}}^{2}}\ \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{C}}&=\sum \limits _{k\in V}{\sum \limits _{r\in R}{\sum \limits _{t=2}^{N-1}{{{x}_{n,t,k}}}}}{{x}_{r,t+1,k}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{D}}&=\sum \limits _{r\in R}{{{\left( 1-\sum \limits _{k\in V}{\sum \limits _{t=2}^{N-1}{{{x}_{r,t,k}}}}\right) }^{2}}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{E}}&=\sum \limits _{k\in V}{\sum \limits _{\begin{array}{c} r,{{r}^{'}}\in {{R}^{'}} \\ r\ne {{r}^{'}} \end{array}}{\sum \limits _{t=1}^{N}{{{x}_{r,t,k}}{{x}_{{{r}^{'}},t,k}}}}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{F}}&=\sum \limits _{k\in V}{\sum \limits _{t=1}^{N-1}{{{\left( \sum \limits _{r\in {{R}_{1}}}{{{x}_{r,t,k}}}-\sum \limits _{{{r}^{'}}\in {{R}_{2}}}{{{x}_{{{r}^{'}},t+1,k}}}\right) }^{2}}}}.\ \end{aligned}$$\end{document}

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\partial }_{i}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(i=1,\ldots ,6)$$\end{document} are weights correspond to each objective function. The objective function (17) is to minimize the total travel time of all the AGVs. The minimization function (18) and (19) ensures that for each AGV, the virtual start task and the virtual end task must be the first task and the last task, respectively. For each actual task, we want it to be assigned exactly to one AGV, so we add the minimization function (20). We also consider the minimization function (21) in order to make each AGV perform at most one task on each task sequence. For a particular AGV, the order of its tasks must be continuous, based on which we set the minimization function (22).

The model described below is modified from the above model to accommodate the goal of minimizing the makespan.

The least time-consuming task model requires finding the AGV with the longest task execution time and minimizing its execution task time, which leads to inequality constraints, as follows.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}&T\ge \sum \limits _{(r,{{r}^{'}})\in A}{\sum \limits _{t=1}^{N-1}{{{c}_{r{{r}^{'}}}}}}{{x}_{r,t,k}}{{{x}_{{{r}^{'}},t+1,k}}}\quad \forall k\in V \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}&{{H}_{Y}}=T \end{aligned}$$\end{document}

The objective function (24) is to minimize the makespan T.

Then we add slack variables to transform the above inequality constraint into an equation, as follows.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}&{{H}_{Z}}={{\sum \limits _{k\in V}{\left( T-\sum \limits _{(r,{{r}^{'}})\in A}{\sum \limits _{t=1}^{N-1}{{{c}_{r{{r}^{'}}}}{{x}_{r,t,k}}{{x}_{{{r}^{'}},t+1,k}}-{{T}_{k}}}}\right) }}^{2}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}&{T=}{{\delta }_{1}}\text {+2}{{\delta }_{2}}\text {+}\ldots \text { +}{{\text {2}}^{m}}{{\delta }_{m}}\ \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned}&{{{T}}_{k}}{=}{{\delta }_{1k}}\text {+2}{{\delta }_{2k}}\text {+}\ldots \text { +}{{\text {2}}^{{{m}^{'}}}}{{\delta }_{{{m}^{'}}k}} \end{aligned}$$\end{document}

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{T}_{k}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(k\in V)$$\end{document} is the slack variable, and both T and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{T}_{k}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(k\in V)$$\end{document} need to be represented by binary variables. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\delta }_{i}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$ (i=1,2,\ldots ,m)$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\delta }_{ik}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(i=1,2,\ldots ,{m}^{'}, {k}\in {V})$$\end{document} are the discretized auxiliary variables we introduce, whose number is related to the size of the arithmetic case and needs to be estimated. A large number of auxiliary variables will make the difficulty of solving soar. In general, to represent an integer between 0 and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\varsigma $$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\left[ {{\log }_{2}} \varsigma \right] +1$$\end{document} discretized auxiliary variables need to be introduced, where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$[\varsigma ]$$\end{document} denotes the largest integer that does not exceed \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\varsigma $$\end{document} . Of course, if there are non-integer values introduced in the calculation example, then it is necessary to introduce an approximate representation. First, we introduce the precision matrix as follows:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} p=\left[ {{\left( \frac{1}{2} \right) }^{0}},{{\left( \frac{1}{2} \right) }^{1}},{{\left( \frac{1}{2} \right) }^{2}},\ldots \right] . \end{aligned}$$\end{document}

Then the real numbers \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\varpi }_{i}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(i=1,\ldots ,L)$$\end{document} in some interval can be approximated as follows.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{\varpi }_{i}}\approx {{p}^\texttt{T}}{{b}_{i}}\text {=}\left[ {{\left( \frac{1}{2} \right) }^{0}},{{\left( \frac{1}{2} \right) }^{1}},{{\left( \frac{1}{2} \right) }^{2}},\ldots \right] ^\texttt{T}{{b}_{i}} \end{aligned}$$\end{document}

where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{b}_{i}}\in {{\{0,1\}}^{L}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(i=1,\ldots ,L)$$\end{document} . If the error is expressed in terms of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\phi $$\end{document} , the error satisfies the following relation:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \phi \le {{\left( \frac{1}{2}\right) }^{L}}. \end{aligned}$$\end{document}

In this way, we can rewrite (24) and (25) as follows.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$T={{\delta}_{1}}+2{{\delta}_{2}}+\ldots+{{2}^{m}}{{\delta }_{m}}+{{\left(\frac{1}{2}\right)}^{0}}{{\sigma }_{1}}+{{\left(\frac{1}{2}\right)}^{1}}{{\sigma}_{2}}+\ldots+{{\left(\frac{1}{2}\right)}^{L-1}}{{\sigma }_{L}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{T}_{k}}={{\delta}_{1k}}+2{{\delta}_{2k}}+\ldots+{{2}^{{{m}^{'}}}}{{\delta }_{{{m}^{'}}k}}+{{\left(\frac{1}{2}\right)}^{0}}{{\sigma }_{1k}}+{{\left(\frac{1}{2}\right)}^{1}}{{\sigma}_{2k}}+\ldots+{{\left(\frac{1}{2}\right)}^{{{L}^{'}}-1}}{{\sigma }_{{{L}^{'}}k}}$$\end{document}

The positive real numbers can be approximated using integer auxiliary variables. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\sigma }_{j}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(j=1,2\ldots ,L)$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\sigma }_{jk}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(j=1,2,\ldots ,{{L}^{'}},k\in V)$$\end{document} are used to approximate the decimal part. The number of binary variables used is related to the required approximate accuracy, as shown in Formula (28).

Thus, the model under this objective can be represented as follows, Eqs. (18)–(22), (24)–(25), (31)–(32).

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text{(P4)}\qquad{{H}_{2}}={{\varepsilon}_{1}}{{H}_{Y}}+{{\varepsilon}_{2}}{{H}_{B}}+{{\varepsilon}_{3}}{{H}_{C}}+{{\varepsilon}_{4}}{{H}_{D}}+{{\varepsilon}_{5}}{{H}_{E}}+{{\varepsilon}_{6}}{{H}_{F}}+{{\varepsilon}_{7}}{{H}_{Z}}$$\end{document}

In (33), \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\varepsilon }_{i}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(i=1\ldots ,7)$$\end{document} are weights for each objectives. The model does not satisfy the QUBO form, because \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{H}_{Z}}$$\end{document} is a quadrinomial binary polynomial, which needs to be degenerated. Next, we provide some analysis of the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{H}_{Z}}$$\end{document} . First, to make it easier to show our results, we perform a variable substitution, as follows:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \eta \text { =}\sum \limits _{(r,{{r}^{'}})\in A}{\sum \limits _{t=1}^{N-1}{{{c}_{r{{r}^{'}}}}}}{{x}_{r,t,k}}{{x}_{r,t+1,k}}. \end{aligned}$$\end{document}

The total number of tasks is N, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\eta $$\end{document} contains \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$({{N}^{2}}-3N+3)(N-1)$$\end{document} monomials. Then \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{H}_{Z}}$$\end{document} can be expanded as follows.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{Z}}=\sum \limits _{k\in V}{({{\eta }^{2}}-2\eta T+2\eta {{T}_{k}}+{{T}^{2}}+T_{k}^{2}-2T{{T}_{k}})} \end{aligned}$$\end{document}

In the \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{H}_{Z}}$$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\eta }^{2}}$$\end{document} is a quadrinomial binary polynomial, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$-2\eta T$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$2\eta {{T}_{k}}$$\end{document} are cubic binary polynomials, and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{T}^{2}}$$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$T_{k}^{2}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$-2T{{T}_{k}}$$\end{document} are all quadratic binary polynomials, so we need to reduce \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\eta }^{2}}$$\end{document} , \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$-2\eta T$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$2\eta {{T}_{k}}$$\end{document} . The number of quardrinomial binary monomials in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\eta }^{2}}$$\end{document} is represented by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\tau }_{4}}$$\end{document} , and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\tau }_{4}}$$\end{document} is as follows:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{\tau }_{4}}={{(({{N}^{2}}-3N+3)(N-1))}^{2}}|V|. \end{aligned}$$\end{document}

Note that

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{x}^{n}}=x (n\in {{\mathbb {Z}}^{+}}),x\in (0,1) \end{aligned}$$\end{document}

where \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\mathbb {Z}}^{+}}$$\end{document} represents the set of all positive integers. Equation (37) indicates that any power of the binary variable itself is equal to itself, and thus, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau _{4}^{'}$$\end{document} quardrinomial monomial can be directly reduced, and its specific number is as follows.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \tau _{4}^{'}=(({{N}^{2}}-3N+3)(N-1))|V| \end{aligned}$$\end{document}

Therefore, the number of quadrinomial polynomials that truly need to be reduced is \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\tau }_{4}}-\tau _{4}^{'}$$\end{document} terms. Looking at the part of cubic polynomial, the number of cubic monomials in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{H}_{Z}}$$\end{document} is represented by \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\tau }_{3}}$$\end{document} , and then \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\tau }_{3}}$$\end{document} is as follows:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{\tau }_{3}}=(m+L+{{m}^{'}}+{{L}^{'}})({{N}^{2}}-3N+3)(N-1)|V|. \end{aligned}$$\end{document}

Then, the number of all polynomials in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{H}_{Z}}$$\end{document} that need to be descended \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau $$\end{document} is as follows:

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \tau ={{\tau }_{3}}+{{\tau }_{4}}-\tau _{4}^{'}. \end{aligned}$$\end{document}

At least \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\tau $$\end{document} binary auxiliary variables must be introduced to complete the descending order according to a paper^35^. Due to the number of auxiliary variables introduced, the processing is more complex and it is difficult to calculate using existing quantum computers. So this model will not be introduced so far in this paper.

Arc model

In this section we will describe the arc model. The core idea of arc model is that the sequence of tasks before and after execution is regarded as an arc connected between nodes, and building the model with arc as the basic unit can reduce the dimension. The arc model also has a QUBO form with the same parameter settings as the node model, and the decision variables are shown below.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\upsilon }_{a,t,k}}$$\end{document} binary variable, equal to 1 if arc a is assigned to AGV k in sequence t, 0 otherwise.

As with the node model, we first explore the model with the optimization objective of minimizing the total travel time.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\text{(P5)}\qquad{{H}_{3}}={{\beta }_{1}}{{H}_{G}}+{{\beta }_{2}}{{H}_{J}}+{{\beta }_{3}}{{H}_{L}}+{{\beta }_{4}}{{H}_{M}}+{{\beta }_{5}}{{H}_{O}}+{{\beta }_{6}}{{H}_{Q}}+{{\beta }_{7}}{{H}_{S}}+{{\beta }_{8}}{{H}_{U}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{G}}&=\sum \limits _{k\in V}{\sum \limits _{a\in A}{\sum \limits _{t=1}^{N-1}{{{c}_{a}}}}}{{\upsilon }_{a,t,k}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{J}}&={{\sum \limits _{k\in V}{\left( 1-\sum \limits _{a\in \theta _{0}^{+}}{{{\upsilon }_{a,1,k}}}\right) }}^{2}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{L}}&=\sum \limits _{k\in V}{\sum \limits _{{\begin{matrix} a\ne {{a}_{1}} \\ a\in \theta _{n}^{-} \\ {{a}_{1}}\in A \end{matrix}}}{\sum \limits _{t=1}^{N-2}{{{\upsilon }_{a,t,k}}{{\upsilon }_{{{a}_{1}},t+1,k}}}}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{M}}&=\sum \limits _{k\in V}{\sum \limits _{{\begin{matrix} a\ne {{a}_{1}} \\ a\in A \\ {{a}_{1}}\in A \end{matrix}}}{\sum \limits _{t=1}^{N-1}{{{\upsilon }_{a,t,k}}}}}{{\upsilon }_{{{a}_{1}},t,k}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} H_O&= \sum _{r\in R}\left( 1-\sum _{k\in V}\sum _{a\in \theta _r^-}\sum _{t=1}^{N-2}v_{a,t,k}\right) ^2 \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} H_Q&= \sum _{k\in V}\left( 1-\sum _{a\in \theta _0^+}\sum _{t=1}^{N-2}v_{a,t,k}\right) ^2 \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} H_S&= \sum _{k\in V}\left( 1-\sum _{a\in \theta _n^-}\sum _{t=2}^{N-1}v_{a,t,k}\right) ^2 \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{U}}&=\sum \limits _{k\in V}{\sum \limits _{r\in R}{\sum \limits _{t=1}^{N-2}{{{\left( \sum \limits _{a\in \theta _{r}^{-}}{{{\upsilon }_{a,t,k}}}-\sum \limits _{{{a}_{1}}\in \theta _{r}^{+}}{{{\upsilon }_{{{a}_{1}},t+1,k}}}\right) }^{2}}}}}. \end{aligned}$$\end{document}

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${{\beta }_{i}}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(i=1,\ldots ,8)$$\end{document} are weights for the corresponding objective. The objective function (42) is to minimize the total travel time of the AGV. We want all AGVs to be executed in the first order a task arc that starts with a virtual start task, so we add the minimization function (43). The minimization function (44) means that the last task completed by each AGV must be a virtual end task. We want to complete at most one task arc in a certain order of an AGV, so we add the minimization function (45). The minimization function (46) indicates that for each actual task, we want a certain task arc with it as the starting node to be assigned to an AGV in a certain order of completion, and we want the virtual start task to be completed only once for each AGV, so we add the minimization function (47). Similarly, for each AGV, its virtual end task can only be completed once, so we add the minimization function (48). The minimization function (49) ensures that for each AGV, it is must to satisfy the flow balance when performing the task arc.

Next, we show a model with the goal of minimizing the makespan as follows, Eqs. (42)–(49).

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} \text {(P6)}\qquad {{H}_{4}}&={{\lambda }_{1}}{{H}_{V}}+{{\lambda }_{2}}{{H}_{J}}+{{\lambda }_{3}}{{H}_{L}}+{{\lambda }_{4}}{{H}_{M}}+\nonumber \\&{{\lambda }_{5}}{{H}_{O}}+{{\lambda }_{6}}{{H}_{Q}}+{{\lambda }_{7}}{{H}_{S}}+{{\lambda }_{8}}{{H}_{U}}+{{\lambda }_{9}}{{H}_{W}} \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{V}}&=T \end{aligned}$$\end{document} \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} {{H}_{W}}&=\sum \limits _{k\in V}{{{\left( T-\sum \limits _{a\in A}{\sum \limits _{t=1}^{N-1}{{{c}_{a}}{{\upsilon }_{a,t,k}}}}-{{T}_{k}}\right) }^{2}}}. \end{aligned}$$\end{document}

Numerical experiments

In this section, we provide rich numerical experimental results, and research data can be obtained on public databases^36^. In the first subsection we use Gurobi^37^ solver to solve the MIP model proposed above on a traditional computer, and show its computing performance under different problem scales. In the second subsection we use optical quantum computer to solve the problem cases of node model and arc model at different scales. And the computation performance is compared with that of traditional computers.

The CIM we used is provided by Beijing Qboson Quantum Technology Co.Ltd, and its structure and principle diagram are shown in Fig. 2. The components of this CIM are mainly composed of optical parts and electrical parts. The optical part of the machine is composed of pulsed laser, erbium-doped fiber amplifier(EDFA), fiber rings and periodically poled lithium niobate (PPLN) crystals, while the electrical part is mainly composed of optical balanced homodyne detectors (BHD), analog-to-digital/digital-to-analog (AD/DA) converter and field-programmable gate array (FPGA). The laser emits laser with a repetition frequency of 100mhz, which is amplified by EDFA, and then the amplified laser frequency is doubled by PPLN crystal to generate 780 nm laser, which is used as the pump source to synchronously pump the phase sensitive amplifier, forming degenerate optical parametric oscillation(DOPO). There are 211 oscillation pulses in the fiber ring, and the time interval between adjacent pulses is 10 ns, so the transmission time of optical pulses in the ring is 2.11 µs. Then, the laser output in the fiber ring and the laser with the fundamental frequency of 1560 nm are determined by BHD, and the FPGA obtains the feedback signal of the next round trip according to the interaction intensity between spins in Ising Hamiltonian, which is used as the control signal of the intensity modulator (IM), and its sign defines the phase shift (0 or \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\pi $$\end{document} ) of the phase modulator (PM)^9,14,34,38^.

To compare the performance between CIM and traditional computer, we also run our experiments with Gurobi 9.5.1 on a Mechrevo computer with 2.8 GHz Intel Core i7 CPU and 8GB memory, using up to four threads. The task points used in this experiment are randomly selected on the two-dimensional axis, ranging from 10 to 90, and then Euclidean distance is used as the length between two points, and we design the speed of each AGV to be constant, the time passing through the unit distance is the unit time.Figure 2. Schematic diagram of coherent Ising mechanism construction and principle.

Computing on a traditional PC

In this section, we use Gurobi to solve the mixed integer programming model of AGV scheduling for two optimization objectives. In “Number of tasks” we show experiments on the variation in computation time with the number of tasks, while in “Number of AGVs” we show experiments on the variation in computation time with the number of AGVs. We set a time limit to 1800 seconds for each run.

Number of tasks

In general, an increase in the number of tasks leads to a slower generation of AGV scheduling solutions. In this subsection, we investigate the effect of task number variation on the computational speed of the three models proposed in this paper. To achieve this goal, we generate instances of 4 tasks to 12 tasks with a fixed number of AGVs of 2 and obtain the computational time graphs shown in Fig. 3, where the left figure takes the minimum total travel time as the objective function, and the right figure takes the minimum makespan as the objective function. The legend section represents the model number, which corresponds to the previous section number.Figure 3. Computation time versus number of tasks for MIP model. (a) Represents the change of MIP model computation time with the number of tasks under the goal of minimizing the total travel time. (b) Represents the change of MIP model computation time with the number of tasks under the goal of minimizing the makespan.

In Fig. 3, we find that the computing speed of mixed integer programming model gradually slows down with the increase of the number of AGV tasks, and the computation time increases sharply when the number of tasks reaches a certain critical value, which is a common property reflected by two different objective functions. Especially when the number of tasks increases to 12, the computing time has exceeded 1800s, which reflects the weakness of traditional models in the face of large-scale problems.

Number of AGVs

In this subsection, we hope to explore the influence of the change of AGV number on the computing time of mixed integer programming model, so we fixed the number of tasks as 10 and 11, and set the number of AGVs in the range from 2 to 8. We show the results in Tables 1 and 2. Notation “–” in the table implies that the corresponding model failed to obtain an optimal solution within the time limit. An instance with a tasks and b AGVs are denoted by “a–b”.Table 1. Computation time (s) under objective of total travel time.instancesMIP10–234.7610–3108.3010–4150.4010–5255.1010–6349.0010–7398.0010–8–11–2340.2011–31262.0011–4–11–5–11–6–11–7–11–8–Table 2. Computation time (s) under objective of makespan.instancesMIP10–217.7510–328.8510–431.1710–561.0110–658.8810–779.9610–850.7911–2338.5011–3393.0011–4186.7011–5204.1011–6188.5011–7503.2011–8405.90

From the results obtained in Tables 1 and 2, we can conclude that there is no strict correlation between the number of AGV and the computing time. Table 1 shows the computational performance of the three models under the objective of minimize the total travel time, and we can observe that the computational performance of the MIP model is very poor, and many groups of experiments failed to obtain an optimal solution within the limited time. Table  2 shows the computational performance of the MIP model under the objective of minimize the makespan, and the model performs much better, the optimal solution is obtained in the limited time in all groups of experiments. However, its computation time is still at a great disadvantage.

Computational experiment on CIM

In this subsection, we use the CIM to solve the QUBO model and compare its performance with the MIP solver on traditional PC. Since the maximum number of Quantum bits of the CIM used in this research is 100, all the comparative examples in this section limit the number of variables to 100. Based on this, we completed six groups of computation experiments with a quantum computer. In all the experiments, the number of AGV was limited to two. For the numerical experiment of node model, we set the number of tasks to 4 to 7, while for the numerical experiment of arc model, we fixed the number of tasks to 4.Figure 4. Evolution diagram of Hamiltonian with time under the objective functions of minimizing the total travel time in node model. (a) Represents the evolution diagram of Hamiltonian with time under the example of 4 tasks. (b) Represents the evolution diagram of Hamiltonian with time under the example of 5 tasks. (c) Represents the evolution diagram of Hamiltonian with time under the example of 6 tasks. (d) Represents the evolution diagram of Hamiltonian with time under the example of 7 tasks.Figure 5. Evolution diagram of Hamiltonian with time under two objective functions of arc mode. (a) represents the evolution diagram of Hamiltonian with time under the objective function of minimizing the total travel time. (b) Represents the evolution diagram of Hamiltonian with time under the objective function of minimizing the makespan.

In Figs. 4 and  5, we plot the evolution of Ising Hamiltonian with time under node model and arc model, respectively. According to the above explanation of CIM principle construction, we can know that the time interval between every two adjacent data points in the figure is 2.11 microseconds. The Hamiltonian decreases with the passage of time, and the phase transition occurs as the power of the pump light gradually increases to the oscillation threshold. The solution obtained when reaching the lowest energy state is the result of CIM solution, and the corresponding time at this time is the computation time.Figure 6. Schematic diagram of quantum computing solutions for node model. (a) represents solution under 4 tasks. (b) represents solution under 5 tasks. (c) represents solution under 6 tasks. (d) represents solution under 7 tasks.Figure 7. Schematic diagram of quantum computing solutions for arc model. (a) represents the solution of 4 task under the objective function of minimizing the total time. (b) represents the solution of 4 task under the objective function of minimizing the makespan.

Figures 6 and 7 show the schematic diagrams of CIM’s solution under node model and arc model. Among them, Fig. 6 contains four parts, which respectively represent the schematic diagrams under 4-7 tasks, while Fig. 7 contains two parts, which respectively represent the schematic diagrams of solving two objective functions under 4 tasks. In fact, Ising model can be transformed into the corresponding representation of maximum cut problem^15^. The maximum cut problem is usually used as a measure and demonstration basis for the complexity of quantum computing problems and the distribution of solutions. The figure shows the solution of our problem expressed by the representation method of maximum cut problem. Points with different colors indicate that they are in different groups, and the connecting lines of points in the same group are gray, while those of points in different groups are red. In these figures, we can clearly perceive the complexity of each model at different scales.

Here, we compare the performance of node model and arc model on quantum computer with that of mixed integer programming model on traditional computer, and the comparison results are shown in Table 3. An instance with a tasks and b AGVs are denoted by “a–b”.Table 3. Computation time (ms) of traditional computer and CIM.InstancesObjectiveMip (solved by traditional computer)Node (solved by CIM)Arc (solved by CIM)TimeObjTimeObjTimeObj4–2Total travel time15.652.822.0052.822.2552.825–2Total travel time31.257.771.0057.77––6–2Total travel time93.760.110.7560.11––7–2Total travel time234.365.350.7365.35––4–2Makespan15.644––2.3044

Due to the limitation of hardware, the comparison of large-scale examples cannot be carried out. However, from Table 3, we can see that the solutions obtained by CIM are all optimal solutions. And the CIM is much faster than the traditional computer in small-scale examples. We can observe that CIM has obvious performance advantages over traditional computers in small-scale examples. In particular, when the scale increases, the time required for CIM does not increase significantly as that of traditional computers. This shows that CIM has great development and application potential. In addition, there is little difference in computing performance between node model and arc model on quantum computer. Node model is slightly faster than arc model, but arc model is more universal than node model. In order to measure the improvement of CIM’s computing efficiency compared with the traditional computer in the given example, we propose the following computation formula.

\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\begin{aligned} IMP=\frac{Q_{TRA}-Q_{CIM}}{Q_{TRA}} \end{aligned}$$\end{document}

IMP represents the calculation speed improvement rate, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q_{TRA}$$\end{document} and \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$Q_{CIM}$$\end{document} respectively represent the computation time of traditional computer and CIM on the same example, and both node model and arc model participate in the comparison. After calculating the IMP of all examples, we find that the computation efficiency of CIM is \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$92\%$$\end{document} faster than that of traditional computers.

Conclusion and future research

We applied quantum computing technology to the research on AGV scheduling, and proposed QUBO models that adapts to solve the problem under two different criteria, minimizing total AGV travel time and makespan. Compared with the traditional MIP model, numerical experiments were carried out on traditional computers and CIM. The experimental results proved the superiority and great potential of quantum computing in this field. Of course, due to the limitation of hardware, there are still some shortcomings in this study, which can not show the advantages of quantum computing in large-scale situations. It is believed that with the continuous development of quantum computing technology, the outstanding performance of quantum computing will be demonstrated in solving large-scale problems in the future. In addition, we also summarized the situation that this study can expand, as follows.

First, the model can be improved and expanded. The model we considered above applies to AGVs with a uniform start node and a uniform end node. Realistic scenarios exist where AGVs have different start and end nodes, and our model is easily expandable for these types of cases. The solution is to involve the number of virtual start and end tasks based on the number of AGVs. Due to the difficulty of mapping large-scale optimization problems to the QUBO form with the small number of bits currently available for quantum computing, our proposed model is a pure scheduling problem that does not consider path optimization. Of course, subsequent researchers can consider this possibility in the case of mature technology. We believe that a two-layer planning model can be built on the basis of the existing model, where the scheduling and path planning problems are computationally solved alternately in two sub-models, which reduces both the modeling difficulty and the number of bits used in the solution of a single model.

Second, quantum computer and traditional computer can be combined to study AGV scheduling problem, and their respective characteristics can be better used to improve the efficiency of solving this problem.

Supplementary Information

Supplementary Information 1.Supplementary Information 2.Supplementary Information 3.