Multitask Emergency Logistics Planning under Multimodal Transportation

Abstract

Multitask emergency logistics planning is a complex optimization problem in practice. When a disaster occurs, relief materials or rescue teams should be dispatched to destinations as soon as possible. In a nutshell, the problem can be described as an optimization of multipoint-to-multipoint transportation delivery problem in a given multimodal traffic network. In this study, a multimodal traffic network is considered for emergency logistics transportation planning, and a mixed-integer programming (MIP) formulation is proposed to model the problem. In order to solve this model, we propose a two-layer solution method. The inner layer is to manage the single-task route recommendation, for which we develop a shortest-path algorithm with the multimodal traffic network. Here, the optimal substructure of the algorithm and its time complexity are presented. With the route of each task calculated by the single-task solver, a general optimization algorithm based on improved particle swarm optimization (PSO) is proposed at the outer layer to coordinate the execution of each task constrained by the limited transportation capacity, so as to derive solutions for multi-commodity emergency logistics planning. Extensive computational results show that the proposed method can find solutions of good quality in reasonable time. Meanwhile, through the sensitivity analysis of the algorithm, we find the appropriate parameters for general optimization algorithm to solve the problem proposed in this paper. The proposed approach is effective and practical for solving multitask emergency logistics planning problem under multimodal transportation, which can find a satisfactory solution in an acceptable time.

1. Introduction

Considering the severe impact of emergencies (e.g., natural hazards and outbreak of epidemics), the timeliness and rationality of emergency logistics planning are critical to reduce the loss caused by emergencies. Currently, emergency logistics is an important field of research in emergency management, as the practice of prevention and rescue during all major natural disasters are accompanied by many logistics activities, such as acquisition, storage, deployment, transportation, distribution and recovery of emergency supplies. In a broad sense, emergency logistics is the process of planning, organizing and controlling the delivery of aid to victims. It also addresses the integration of social sciences and analytical research, incorporating mathematical characterization and modelling of diverse aspects of relief efforts, together with an intense field work [1,2].

Unfortunately, experience with natural disasters globally shows that the vast majority of disaster management is not well planned in terms of logistics. Meanwhile, traffic networks can be seriously damaged after disasters, but transportation route planning is typically neglected.

One of the early works related to emergency logistics management was promoted by Sheu [3]. Based on this, a recent survey of emergency logistics indicates that most of the relevant research focuses on the transportation of supplies and wounded people, applying models such as Dynamic Network Flow (DNF), Static Network Flow (SNF), Route Enumeration (RE), Location-Routing Problem (LRP), Vehicle Routing Problem (VRP), Multimodal Transportation System (MTTS) and many more [4]. However, to the best of our knowledge, there is a lack of research on multitask emergency logistics planning under multimodal transportation. The optimization of emergency logistics as a post-disaster decision-making activity also has a lot of value in practice and thus deserves further investigation.

In essence, the optimization problem of multi-commodity emergency logistics planning can be abstracted as a multipoint-to-multipoint transportation planning problem in a given multimodal traffic network. The characteristic of multimodal transportation is to transport supplies through two or more different modes of transportation, such as highway, railway, aviation, inland waterways, and short-distance or long-distance maritime transportation. An illustration of the emergency logistics planning problem is shown in Figure 1. With an increasing attention to emergency logistics, emergency logistics under multimodal transportation has also been a heated research field, extending to many real-world applications [5,6,7].

The focus of this study is emergency logistics planning for multi-commodity delivering under multimodal transportation. Many studies regard emergency logistics as a complex nonlinear network planning problem involving multiple constraints and multiple objectives [8,9,10,11,12,13]. Based on the existing literature, we take into account several characteristics that are crucial in multitask emergency logistics planning as follows:

• Limited available infrastructure capacity: emergency logistics rely on transportation routes with limited capacity, and the capacity of loading and unloading nodes is typically limited, which can easily lead to traffic network congestion in mass transportation;
• Insufficient available transport capacity: when performing emergency logistics tasks, we are usually faced with a relatively fragile traffic network and cannot undertake large-scale air transportation missions. For example, the design of the load and size of the bridges, culverts, railway beds and highway beds restrict the passage of some heavy rescue equipment;
• Many tasks involved: emergency logistics typically involve multiple emergency tasks, and these tasks are scattered with different destinations, which incur great difficulties to model and solve the emergency logistics planning problem.

Based on the above realistic difficulties, this paper conducts research on multitask emergency logistics planning under multimodal transportation. In this study, we first develop a mixed-integer programming (MIP) model for the emergency logistics planning problem. Considering the complexity of this model, we decompose the solution procedure into two layers. The inner layer is the route planning for a single task, and a shortest-path algorithm under a multimodal traffic network is devised. After the route of each task is calculated, the alternative route set is obtained. In the outer layer, multi-task emergency logistics planning should be optimized according to the information of each task and its multiple alternative routes, and the limited capacity in the multimodal transportation network. Based on the mathematical model, a multi-modal optimization algorithm based on the particle swarm optimization (PSO) algorithm is proposed. Through this work, good quality solutions can be found and the multitask emergency logistics planning problem can be solved effectively.

The remainder of this paper is organized as follows. In Section 2, we review the related literature. In Section 3, we present a MIP model for the multitask emergency logistics problem. In Section 4, a shortest-path algorithm under a multimodal traffic network is proposed to solve the single-task model. We also construct the real transport time calculation algorithm to recalculate the shortest-path with the minimum weighted time obtained by the previous algorithm. In Section 5, we propose a general optimization algorithm of multitask emergency logistics planning to coordinate the execution of each task constrained by the limited transportation capacity. Computational results are presented and discussed in Section 6. Finally, in Section 7, we conclude the work and point out some directions for future research.

2. Literature Review

Since its introduction in the early 1980s [14], emergency logistics has drawn a lot of attention and many mathematical models have been developed to solve related problems. Many extensions have also been made to incorporate various characteristics of real-world problems, and the primary research directions have been optimization modelling, algorithm design and model evaluation.

Based on the MTTS, the mathematical model that we present for multitask emergency logistics has its roots in the work that first established a multicommodity and multimodal emergency supply distribution model based on a time-space network [15]. In his work, the large-scale disaster relief transportation problem was described as a multicommodity and multimodal network flow model with a single objective function. Since then, the relevant research had focused on the mixed integer programming model and the programming problem under uncertainty.

Ransikarbum and Mason [16] presented a multiobjective optimization model in an integrated network for making strategic decisions in the supply distribution and network restoration phases of humanitarian logistics operations. To use the scarce resources efficiently and achieve the best possible emergency relief, Najafi et al. [17] proposed a multi-objective, multi-mode, multi-commodity, and multi-period stochastic model for managing the logistics of both commodities and injured people in earthquake response. Yoon and Albert [18] formulated a Markov decision process model that dispatches multiple types of ambulances to multiple patient priorities in emergency medical services, which determines the type of ambulances to be deployed for patients in real-time. Barbarosoǧlu and Arda [19] proposed a scenario-based two-stage stochastic programming model to plan disaster relief transportation, and developed a multi-commodity, multi-mode network flow formula to describe logistics on urban transportation networks. Gao et al. [20] proposed a bi-objective stochastic optimization model to rebalance and transport commodities with the multi-modal transportation system, aiming at solving the multi-commodity rebalancing issue for humanitarian logistics.

At the same time, some scholars also studied the problem of planning and scheduling in emergency logistics. Ding et al. [21] put forward an emergency material scheduling model with multiple logistics supply points to multiple demand points based on the grey interval number, aiming at addressing the problem of post-disaster emergency material dispatching from multiple supply points to multiple demand points. Song et al. [22] designed a multi-stage emergency logistics scheduling model for multi-supply points and multi-demands points, aiming at minimizing total dispatching cost and emergency points’ satisfaction.

It can be seen that problems have been investigated from different perspectives. However, there is a lack of general description and mathematical modeling for multitask emergency logistics planning under multimodal transportation.

In this study, we investigate the optimization problem of emergency logistics. Solution algorithms of such problems can be categorized as either exact and heuristic. The multitask emergency logistics planning optimization problem investigated in this study is a typical NP-hard problem [23,24]. Due to the complexity of mathematical models of such practical problems, the exact algorithm cannot find the optimal solution in a limited time. Then many scholars had focused on designing and improving heuristic optimization algorithms to find corresponding problems’ satisfactory solutions. Yi and Kumar [25] first proposed an ant colony optimization (ACO) heuristic algorithm for disaster relief operations. Since then, many scholars have tried to use the optimization algorithm combining PSO and ACO to solve the multimodal transportation scheme, and achieved good results [26,27,28]. Li and Su [29] solved the route selection problem in multimodal transportation mainly by using improved genetic algorithm. Fazayeli et al. [30] gave a presentation of a two-part genetic algorithm for addressing the distribution problem in multimodal transport network. Das et al. [31] modified genetic algorithm with an in-vitro-fertilization-based crossover as well as a generation-dependent mutation, and the numerical results obtained by this algorithm are superior to other heuristic algorithms. To deal with the heterogeneous fleet vehicle routing problem, Subramanian et al. [32] proposed a hybrid algorithm composed by an Iterated Local Search based heuristic and a Set Partitioning formulation.

3. Problem Description and Mathematical Formulation

To minimize the impact of a disaster, a critical issue in emergency logistics is to provide rescue materials and perform rescue operations correctly and in a timely manner. Wang et al. [33] mentioned that the occurrence of natural disasters or accidents, such as COVID-19, can obstruct or interrupt traffic connectivity and may affect the transportation of essential materials. Therefore, the problem to be solved by multitask emergency logistics planning optimization can be described as follows: using limited traffic resources, on the premise of meeting all types of real situations and superior requirements, choose the appropriate route and time, and complete the delivery of each task from its departure point to the destination point in the shortest time possible.

Most importantly, it is necessary to clarify all types of key information required in the optimization process:

• Traffic network (i.e., the topology of the emergency logistics transportation system): the traffic network is composed of nodes and arcs, where nodes are typically cities, stations and transfer stations; and arcs are typically transportation lines between two nodes;
• Loading and unloading capacity (i.e., the maximum periodic loading and unloading capacity of a node for a specified transportation mode);
• Transport capacity: the transport capacity refers to the periodic maximum number of batches passing arcs between two nodes, which is jointly determined by the type of route, logistics support capacity and other factors;
• Transportation routes (i.e., the collection of routes that connect the departure point and destination point of each task);
• Periodic delivery quantity (i.e., the number of batches delivered periodically for each task): if the number of task batches is greater than the maximum single-day transport capacity of its transportation route, it will likely take several periods to complete the delivery. The maximum single-day transport capacity is determined by the loading and unloading capacity of the departure point, destination point and transfer node as well as the transport capacity of its transportation route.

Constraints that must also be considered in the optimization process should be clarified. The problem investigated in this study is essentially a type of complex optimization problem that must find the optimal solution according to the optimization objective under various constraints. The main constraints considered are as follows:

• Traffic resource constraints: these constraints primarily include the maximum capacity of each route, the maximum loading and unloading capacity of each node, and the requirements of each task on the type of transportation mode;
• Emergency principle constraints: these constraints primarily include the last arrival time of each task, arrival order of the batches in different tasks, batch continuity in time and space during transportation;
• Natural environment Constraints (e.g., logistic requirements of high-altitude areas and the transport time limitations of severe weather);
• Transfer condition constraints: these constraints primarily include the maximum capacity of each node and the priority requirements for different modes of transport. In theory, transfer is possible between any two modes of transportation. However, in the real transport process, we consider that emergency logistics should take the time benefit as the first goal. The transport time is faster with the mode of transportation with higher priority, and the transport distance should be as far as possible with a high priority without transfer. Therefore, we assume that priority exists between different traffic networks and stipulate that transfer can only occur from a high-priority to a low-priority traffic network.

The traffic capacity, transport route selection and schedule of tasks must be considered to avoid congestion at key routes or nodes. To model the optimization problem of multitask emergency logistics planning, we first must clarify the task list information and transport capacity data of the traffic network. This list is a general emergency logistics planning list that is used to transport the given task list from the designated departure point to the designated destination point in the existing traffic network, with the shortest total completion time as the goal.

To describe the model clearly,

Table 1. Notations of related parameters used in the models.

Notations	Definition
F	The set of transportation modes
f	The index of transportation modes and $f\in F$
T	The set of planning horizon
t	The index of time and $t\in T$
$N_f$	Subnetwork of transportation mode f
$l_i^f$	Periodic loading capacity of node i under transportation mode f
$u{l}_i^f$	Periodic unloading capacity of node i under transportation mode f
$r_{ij}^f$	Periodic traffic capacity of arc $(i,j)$ under transportation mode f
$t{r}_{ij}^f$	Transportation time of arc $(i,j)$ under transportation mode f
M	The set of tasks
m	The index of tasks and $m\in M$
$u_m$	The number of batches in task m
${\underline{u}}_m$	Minimum number of batches in task m transported in a period
$s{w}_m^{f_1,f_2}$	Whether task m can transfer from transportation mode $f_1$ to $f_2$ ($f_1$ has a higher priority than $f_2$)
$s_m$	The departure point of task m
$q_m$	The destination point of task m
$e{s}_m$	The earliest departure time of task m
$l{f}_m$	The last arrival time of task m
U	A large constant

shows the notations of the related parameters used in the models, and

Table 2. Notations of decision variable used in the model.

Notations	Definition
$Z_t$	Binary variable indicating whether there is any task being executed on period t and $Z_t\in \left\{\begin{array}{c}0,1\end{array}\right\}$
$st{a}_{mt}^f$	Binary variable indicating whether task m begins delivering from departure point $s_m$ on period t under the transportation mode f and $st{a}_{mt}^f\in \left\{\begin{array}{c}0,1\end{array}\right\}$
$fi{n}_{mt}^f$	Binary variable indicating whether task m finishes delivering from departure point $s_m$ on period t under the transportation mode f and $fi{n}_{mt}^f\in \left\{\begin{array}{c}0,1\end{array}\right\}$
$mas{k}_{mij}^f$	Binary variable indicating whether task m chooses the arc $(i,j)$ to transport under the transportation mode f and $mas{k}_{mij}^f\in \left\{\begin{array}{c}0,1\end{array}\right\}$
$x_{mijt}^f$	Binary variable indicating whether task m chooses the arc $(i,j)$ to transport on period t under the transportation mode f and $x_{mijt}^f\in \left\{\begin{array}{c}0,1\end{array}\right\}$
$y_{mijt}^f$	The number of batches in task m transport in arc $(i,j)$ on period t under the transportation mode f
$swi{t}_{mi}^{f_1,f_2}$	Whether task m transfer from transportation mode $f_1$ to $f_2$ at node i ($f_1$ has a higher priority than $f_2$)
$aux{l}_{mit}^f$	Task m’s auxiliary variable of loading capacity at node i on period t under the transportation mode f
$auxu{l}_{mit}^f$	Task m’s auxiliary variable of unloading capacity at node i on period t under the transportation mode f
$X_{mt}$	Binary variable indicating whether there are any batches in task m deliver from departure point $s_m$ on period t and $X_t\in \left\{\begin{array}{c}0,1\end{array}\right\}$

shows the notations of the decision variables used in the models.

The topology of an emergency logistics transportation system can be described by a directed network consisting of the set of vertices and the set of arcs. The vertices in the network may denote terminals and the junctions between two arcs, and the arcs represent the pathways connecting two vertices. We consider different transportation modes, which are often jointly used in emergency logistics. The transportation network can be denoted as $N=(V,E)$, and the subnetwork of transportation mode f can be abstracted as $N_f=(V_f,E_f)$, where the index of node is $i\in V_f$ and $(i,j)\in E_f$ represents the arc from i to j under transportation mode f. The node set is $V=\left\{\begin{array}{c}{V}_1\cup V_2\cup \dots \cup V_{\left|F\right|}\end{array}\right\}$, and the arc set is $E=\left\{\begin{array}{c}{E}_1\cup E_2\cup \dots \cup E_{\left|F\right|}\end{array}\right\}$.

According to the problem description, a general multitask emergency logistics planning optimization model is established in this section. Several relevant assumptions should be considered as follows:

• Each task to be delivered can be split into batches. Each batch in each task cannot be split during transportation, and a batch is the smallest unit;
• The transportation of each task must be continuous, in terms of both time and space;
• Transfer can be made instantaneous and can only be carried out from a transportation mode of a higher priority to that of a lower priority.

The problem can be modelled as follows:

$$\begin{array}{cc}\hfill \mathrm{minimise}& \sum _{t\in T}{Z}_t\hfill \end{array}$$

(1)

$$\begin{array}{cc}\hfill \mathrm{subject}\mathrm{to}& \sum _{t\in T}\sum _{f\in F}\sum _{(i,q_m)\in E_f}{y}_{m,i,q_m,t}^f=u_m,\forall m\in M\hfill \end{array}$$

(2)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{j:(s_m,j)\in E_f}mas{k}_{m,s_m,j}^f=1,\forall m\in M\hfill \end{array}$$

(3)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{i:(i,q_m)\in E_f}mas{k}_{m,i,q_m}^f=1,\forall m\in M\hfill \end{array}$$

(4)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{j:(i,j)\in E_f}mas{k}_{mij}^f=\sum _{f\in F}\sum _{j:(j,i)\in E_f}mas{k}_{mji}^f,\forall m\in M,\forall i\in V\setminus \left\{\begin{array}{c}{s}_m,q_m\end{array}\right\}\hfill \end{array}$$

(5)

$$\begin{array}{cc}\hfill & x_{mijt}^f\le mas{k}_{mij}^f,\forall m\in M,\forall f\in F,\forall (i,j)\in E,\forall t\in T\hfill \end{array}$$

(6)

$$\begin{array}{cc}\hfill & x_{mijt}^f·{\underline{u}}_m\le y_{m,i,j,t}^f\le x_{mijt}^f·u_m,\forall m\in M,\forall f\in F,\forall (i,j)\in E,\forall t\in T\hfill \end{array}$$

(7)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{j:(i,j)\in E_f}{y}_{m,i,j,t-t{r}_{ji}^f}^f=\sum _{f\in F}\sum _{k:(j,k)\in E_f}{y}_{m,j,k,t}^f,\forall m\in M,\forall (i,j)\in E,\forall t\in T\hfill \end{array}$$

(8)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{t\in T}st{a}_{mt}^f=1,\forall m\in M\hfill \end{array}$$

(9)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{t\in T}fi{n}_{mt}^f=1,\forall m\in M\hfill \end{array}$$

(10)

$$\begin{array}{cc}\hfill & \sum _{t\in T}st{a}_{mt}^f=\sum _{t\in T}fi{n}_{mt}^f,\forall m\in M,\forall f\in F\hfill \end{array}$$

(11)

$$\begin{array}{cc}\hfill & \sum _{k=1}^{t}st{a}_{mk}^f\ge \sum _{k=1}^{t}fi{n}_{mk}^f,\forall m\in M,\forall f\in F,\forall t\in T\hfill \end{array}$$

(12)

$$\begin{array}{cc}\hfill & \sum _{k=1}^{t}st{a}_{mk}^f-\sum _{k=1}^{t}fi{n}_{mk}^f=\sum _{(s_m,j)\in E_f}{x}_{m,s_m,j,t}^f,\forall m\in M,\forall f\in F,\forall t\in T\hfill \end{array}$$

(13)

$$\begin{array}{cc}\hfill & (\sum _{k=1}^{t}st{a}_{mk}^f-\sum _{k=1}^{t}fi{n}_{mk}^f)·{\underline{u}}_m\le \sum _{(s_m,j)\in E_f}{y}_{m,s_m,j,t}^f,\forall m\in M,\forall f\in F,\forall t\in T\hfill \end{array}$$

(14)

$$\begin{array}{cc}\hfill & \sum _{(s_m,j)\in E_f}{y}_{m,s_m,j,t}^f\le (\sum _{k=1}^{t}st{a}_{mk}^f-\sum _{k=1}^{t}fi{n}_{mk}^f)·u_m,\forall m\in M,\forall f\in F,\forall t\in T\hfill \end{array}$$

(15)

$$\begin{array}{cc}\hfill & \sum _{t=1}^{e{s}_m-1}\sum _{f\in F}\sum _{(s_m,j)\in E_f}{y}_{m,s_m,j,t}^f=0,\forall m\in M\hfill \end{array}$$

(16)

$$\begin{array}{cc}\hfill & \sum _{t=l{f}_m+1}^{\left|T\right|}\sum _{f\in F}\sum _{(i,q_m)\in E_f}{y}_{m,i,q_m,t}^f=0,\forall m\in M\hfill \end{array}$$

(17)

$$\begin{array}{cc}\hfill & swi{t}_{mi}^{f_1,f_2}\le s{w}_m^{f_1,f_2},\forall i\in V\setminus \left\{\begin{array}{c}{s}_m,q_m\end{array}\right\},\forall f_1,f_2\in F,\forall m\in M\hfill \end{array}$$

(18)

$$\begin{array}{cc}\hfill & \sum _{j:(j,i)\in E_{f_1}}mas{k}_{mji}^{f_1}\ge swi{t}_{mi}^{f_1,f_2},\forall i\in V\setminus \left\{\begin{array}{c}{s}_m,q_m\end{array}\right\},\forall m\in M,\forall f_1,f_2\in F\hfill \end{array}$$

(19)

$$\begin{array}{cc}\hfill & \sum _{j:(i,j)\in E_{f_2}}mas{k}_{mij}^{f_2}\ge swi{t}_{mi}^{f_1,f_2},\forall i\in V\setminus \left\{\begin{array}{c}{s}_m,q_m\end{array}\right\},\forall m\in M,\forall f_1,f_2\in F\hfill \end{array}$$

(20)

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill & \sum _{j:(j,i)\in E_{f_1}}mas{k}_{mji}^{f_1}=\sum _{f_2\in F\setminus f_1}swi{t}_{mi}^{f_1,f_2}+\sum _{k:(i,k)\in E_{f_1}}mas{k}_{mik}^{f_1},\forall i\in V\hfill \\ \hfill & \forall f_1,f_2\in F,\forall m\in M\hfill \end{array}\end{array}$$

(21)

$$\begin{array}{cc}\hfill & \sum _{m\in M}{y}_{m,i,j,t}^f\le r_{ij}^f,\forall f\in F,\forall (i,j)\in E,\forall t\in T\hfill \end{array}$$

(22)

$$\begin{array}{cc}\hfill & x_{mijt}^f\le Z_t,\forall m\in M,\forall f\in F,\forall (i,j)\in E,\forall t\in T\hfill \end{array}$$

(23)

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill & aux{l}_{mit}^{f_2}+U·(1-\sum _{f_1\in F\setminus f_2}swi{t}_{mi}^{f_1,f_2})\ge \sum _{j:(i,j)\in E_{f_2}}{y}_{m,i,j,t}^{f_2},\forall m\in M,\hfill \\ \hfill & \forall f_1,f_2\in F,\forall i\in V,\forall t\in T\hfill \end{array}\end{array}$$

(24)

$$\begin{array}{cc}\hfill & aux{l}_{mit}^{f_2}+U·\sum _{f_1,f_2\in F}swi{t}_{mi}^{f_1,f_2}\ge 0,\forall m\in M,\forall f_1,f_2\in F,\forall i\in V,\forall t\in T\hfill \end{array}$$

(25)

$$\begin{array}{cc}\hfill & \sum _{m\in M}aux{l}_{mit}^f+\sum _{m:s_m=i}\sum _{j:(j,i)\in E_f}{y}_{m,i,j,t}^f\le l_i^f,\forall f\in F,\forall i\in V,\forall t\in T\hfill \end{array}$$

(26)

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill & auxu{l}_{mit}^{f_1}+U·(1-\sum _{f_2\in F\setminus f_1}swi{t}_{mi}^{f_1,f_2})\ge \sum _{j:(j,i)\in E_{f_1}}{y}_{m,j,i,t}^{f_1},\forall m\in M,\hfill \\ \hfill & \forall f_1,f_2\in F,\forall i\in V,\forall t\in T\hfill \end{array}\end{array}$$

(27)

$$\begin{array}{cc}\hfill & auxu{l}_{mit}^{f_1}+U·\sum _{f_1,f_2\in F}swi{t}_{mi}^{f_1,f_2}\ge 0,\forall m\in M,\forall f_1,f_2\in F,\forall i\in V,\forall t\in T\hfill \end{array}$$

(28)

$$\begin{array}{cc}\hfill & \sum _{m\in M}auxu{l}_{mit}^f+\sum _{m:q_m=i}\sum _{j:(j,i)\in E_f}{y}_{m,j,i,t}^f\le u{l}_i^f,\forall f\in F,\forall i\in V,\forall t\in T\hfill \end{array}$$

(29)

$$\begin{array}{cc}\hfill & Z_t+1\le Z_t,\forall t\in T\hfill \end{array}$$

(30)

Objective function (1) minimizes the total completion time of multitask emergency logistics, i.e., the arrival time of the last task. Constraints (2) state that all tasks have been transported. Constraints (3)–(5) guarantee that the transportation route of each task is physically continuous. Constraints (6)–(7) limit the number of batches of each task that can be delivered on period t, and Constraints (8) define the continuity of the transportation in the traffic network. Constraints (9)–(12) state that each task should be delivered in a continuous period of time. Constraints (13) stipulate that delivery can only be implemented within the selected period. Constraints (14) and (15) define the upper and lower bounds of the number of batches transported from departure point $s_m$ on period t.

Regarding the constraint description for the transport time window, Constraints (16) and (17) enforce the transportation of each task within the time window. Constraints (18) state that transfer occurs only on node i with transfer capability. Then, Constraints (19)–(21) ensure the connectivity between transfer node i and the transport route of task m.

To describe the traffic capacity, Constraints (22) stipulate that the number of batches transported in the route should be no larger than the periodic traffic capacity of each arc of task m’s transport line, and Constraints (23) determine whether there is any task being executed on period t. Constraints (24)–(26) also limit the loading capacity of transport node i, and Constraints (27)–(29) limit the unloading capacity of transport node i. Constraints (30) ensure the early delivery of each task.

4. Solution Method for Single-Task Route Planning

Single-task route planning, as a simple variant of the problem investigated, is also a very common problem in practice; thus, it is necessary to devise an efficient solution method for this basic problem. On the other hand, the multitask emergency logistics planning problem is highly complex, thus we base its solution algorithm on the solution of the single-task planning problem, which is discussed in detail in Section 5. In this section, we first describe the single-task path planning problem and establish a route planning model in Section 4.1. Then, we propose a shortest-path algorithm under a multimodal traffic network in Section 4.2 to solve the route planning model. Extensions of the algorithm to incorporate more timing aspects are discussed in Section 4.3, and an algorithm to retrieve the real transport time is proposed.

4.1. Single-Task Route Planning Model

For an emergency logistics planning problem with one transportation task, the key problem is to optimize the transportation route from its departure point to the destination point. In practice, factors such as the batches of tasks, route transport time, transport capacity, transport continuity of tasks and transfer time should be considered comprehensively to recommend the optimal transport route to arrive and perform in the fastest time.

We now develop a single-task route planning model. For simplicity, the index m in the model has been neglected. For the single task, the problem can be modelled as:

$$\begin{array}{cc}\hfill \mathrm{minimise}& \sum _{f\in F}\sum _{(i,j)\in E_f}t{r}_{ij}^f·mas{k}_{ij}^f+(\sum _{t\in T}{X}_t-1)\hfill \end{array}$$

(31)

$$\begin{array}{cc}\hfill \mathrm{subject}\mathrm{to}& \sum _{f\in F}\sum _{(s,j)\in E_f}mas{k}_{s,j}^f-\sum _{f\in F}\sum _{(j,s)\in E_f}mas{k}_{j,s}^f=1\hfill \end{array}$$

(32)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{(q,j)\in E_f}mas{k}_{q,j}^f-\sum _{(j,q)\in E_f}\sum _{f\in F}mas{k}_{j,q}^f=-1\hfill \end{array}$$

(33)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{(i,j)\in E_f}mas{k}_{ij}^f-\sum _{f\in F}\sum _{(j,i)\in E_f}mas{k}_{ji}^f=0,\forall i\in V\setminus \left\{\begin{array}{c}s,q\end{array}\right\}\hfill \end{array}$$

(34)

$$\begin{array}{cc}\hfill & \sum _{j:(j,i)\in E_{f_1}}mas{k}_{ji}^{f_1}\ge swi{t}_i^{f_1,f_2},\forall i\in V\setminus \left\{\begin{array}{c}s,q\end{array}\right\},\forall f_1,f_2\in F\hfill \end{array}$$

(35)

$$\begin{array}{cc}\hfill & \sum _{j:(i,j)\in E_{f_2}}mas{k}_{ij}^{f_2}\ge swi{t}_i^{f_1,f_2},\forall i\in V\setminus \left\{\begin{array}{c}s,q\end{array}\right\},\forall f_1,f_2\in F\hfill \\ \hfill & \sum _{f_2\in F\setminus f_1}swi{t}_i^{f_1,f_2}+\sum _{j:(j,i)\in E_{f_2}}mas{k}_{ji}^{f_2}=\sum _{f_3\in F\setminus \{f_1,f_2\}}swi{t}_i^{f_1,f_2}+\sum _{k:(i,k)\in E_{f_3}}mas{k}_{ik}^{f_3},\hfill \end{array}$$

(36)

$$\begin{array}{cc}\hfill & \forall i\in V,\forall f_1,f_2,f_3\in F\hfill \end{array}$$

(37)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{(j,i)\in E_f}{y}_{j,i,t-t{r}_{ji}^f}^f=\sum _{f\in F}\sum _{(i,k)\in E_f}{y}_{i,k,t}^f,\forall t\in T,\forall i\in V\setminus \left\{\begin{array}{c}s,q\end{array}\right\}\hfill \end{array}$$

(38)

$$\begin{array}{cc}\hfill & y_{ijt}^f\le mas{k}_{ij}^f·r_{ij}^f,\forall (i,j)\in E_f,\forall f\in F,\forall t\in T\hfill \end{array}$$

(39)

$$\begin{array}{cc}\hfill & \sum _{(s,j)\in E_f}{y}_{s,j,t}^f\le l_s^f,\forall t\in T,\forall f\in F\hfill \end{array}$$

(40)

$$\begin{array}{cc}\hfill & \sum _{(i,q)\in E_f}{y}_{i,q,t}^f\le u{l}_q^f,\forall t\in T,\forall f\in F\hfill \end{array}$$

(41)

$$\begin{array}{cc}\hfill & \sum _{j:(i,j)\in E_{f_1}}{y}_{ijt}^{f_1}\le swi{t}_j^{f_1,f_2}·u{l}_j^{f_1},\forall t\in T,\forall f_1,f_2\in F\hfill \end{array}$$

(42)

$$\begin{array}{cc}\hfill & \sum _{j:(i,j)\in E_{f_1}}{y}_{ijt}^{f_1}\le swi{t}_j^{f_1,f_2}·l_j^{f_2},\forall t\in T,\forall f_1,f_2\in F\hfill \end{array}$$

(43)

$$\begin{array}{cc}\hfill & \sum _{t\in T}\sum _{f\in F}\sum _{(s,j)\in E_f}{y}_{s,j,t}^f=\sum _{t\in T}\sum _{f\in F}\sum _{(i,q)\in E_f}{y}_{i,q,t}^f=u\hfill \end{array}$$

(44)

$$\begin{array}{cc}\hfill & swi{t}_i^{f_1,f_2}\le s{w}^{f_1,f_2},\forall i\in V,\forall f_1,f_2\in F\hfill \end{array}$$

(45)

$$\begin{array}{cc}\hfill & \sum _{j:(i,j)\in E_f}{y}_{i,j,t-tf}^{f_1}·swi{t}_j^{f_1,f_2}=\sum _{j:(j,k)\in E_f}{y}_{j,k,t}^{f_2}·swi{t}_j^{f_1,f_2},\forall t\in T,\forall f_1,f_2\in F\hfill \end{array}$$

(46)

$$\begin{array}{cc}\hfill & \sum _{f\in F}\sum _{(s,j)\in E_f}{y}_{s,j,t}^f\ge \sum _{f\in F}\sum _{(s,j)\in E_f}{y}_{s,j,t+1}^f,\forall j\in V\hfill \end{array}$$

(47)

$$\begin{array}{cc}\hfill & X_t·U\ge \sum _{f\in F}\sum _{(s,j)\in E_f}{y}_{s,j,t}^f,\forall t\in T\hfill \end{array}$$

(48)

The objective function (31) minimizes the total transport time of a single task consisting of transport time on the route and the number of periods required to deliver all the batches. Constraints (32)–(34) state that the selected route is physically continuous from departure point s to destination point q. Specifically, when $mas{k}_{ij}^f=1$, the corresponding arc $(i,j)$ under transportation mode f is selected as the transport route. Otherwise, the arc $(i,j)$ will not appear on the selected route. Constraints (35)–(37) ensure the connectivity between transfer node i and the selected route. When $swi{t}_i^{f_1,f_2}=1$, the transfer from transportation mode $f_1$ to $f_2$ is completed at node i. Then, Constraints (38) ensure that the number of batches transported on the selected routes is continuous and consistent. The number of batches transported per period is limited by Constraints (39)–(41), which must be less than the periodic traffic capacity of each selected arc, the periodic loading capacity of departure point s and the unloading capacity of destination point q. Concurrently, Constraints (42) and (43) limit the number of batches transported in route less than the periodic loading and unloading capacity of the transfer node i. Constraint (44) guarantees that all batches in task m have been transported. Constraints (45) and (46) ensure that transfer occurs only when the conditions are met and the number of batches is consistent before and after transfer. Finally, Constraints (47) state that the task should be finished as early as possible, and Constraints (48) are used to calculate the number of periods needed to complete delivering all batches.

4.2. Shortest-Path Algorithm under Multimodal Traffic Networks

The Dijkstra algorithm was first proposed by [34], which has been a classic algorithm for solving the shortest-path problem. The algorithm can effectively solve the single-modal shortest-path problem and it takes the starting point as the centre point to extend outward until the end point has been reached, then the shortest-path between two points is obtained.

So far, the literature regarding solution methods for finding the shortest-path on a multimodal traffic network is sparse. Based on the traditional Dijkstra algorithm, we propose a shortest-path path algorithm under multimodal traffic networks according to the characteristics of emergency logistics. The traffic networks are based on the descriptions in Section 3.

In our shortest-path algorithm under a multimodal traffic network, network splitting is used to depict various traffic networks. An illustration of network splitting in a double modal (rail-highway) network is shown in Figure 2. The connection lines refer to the transfer between different transportation modes. The nodes in the high-priority network are connected only in one direction to the nodes in the low-priority network, and the reverse connection time is infinite, which prohibits the transfer.

When we only need to consider the first part of objective (31) in Section 4.1 (i.e., the transport time on the route), an exact algorithm can be achieved using the Dijkstra algorithm based on the aforementioned transformation of the multimodal network. The primary idea of the minimum time path algorithm is described as follows. Firstly, we need to determine the highest priority subnetwork $N_f$ where departure point s is located and mark the index as i, then the departure point goes to $s_i$. We thus introduce two sets P and H, where P records the nodes where the minimum time has been found, and H records the nodes where the minimum time has not been found. Considering the split networks, the nodes are assigned with indices related to their traffic subnetworks. Then we traverse the node with the minimum time from H, add this node to P, and update the paths of nodes in H. We loop through this process until the destination point $q_j(j\ge i)$ is reached (destination point may exist in different subnetworks), and then compare the times of different destination points $q_j(j\ge i)$. Finally, the shortest-path with minimum time can be identified.

The pseudocode of the proposed shortest-path algorithm under a multimodal traffic network is shown in Algorithm 1.

Algorithm 1: Shortest-Path Algorithm under Multimodal Traffic Network

We also now demonstrate the correctness of the algorithm by proving the optimal substructure of the shortest-path.

Proposition 1.

Given that $SP(s,q,f^s,f^q)$ is the shortest-path from s to q using transportation modes with priority levels from $f^s$ to $f^q$, with $f^i$ denoting the transportation mode selected to reach node i. Assuming that $SP(s,q,f^s,f^q)=\left\{\begin{array}{c}{v}_s^{f^s},\dots ,v_k^{f^k},\dots ,v_j^{f^j},\dots v_q^{f^q}\end{array}\right\}$, and $\left\{\begin{array}{c}{v}_k^{f^k},\dots ,v_j^{f^j}\end{array}\right\}$ is an arbitrary segment along the $SP(s,q,f^s,f^q)$, there must be $SP(k,j,f^k,f^j)=\left\{\begin{array}{c}{v}_k^{f^k},\dots ,v_j^{f^j}\end{array}\right\}$, which is shortest-path from node k to j using transportation modes with priority levels from $f^k$ to $f^j$. So the shortest-path has the property of an optimal substructure.

 Proof of Proposition 1.

We already know that the shortest-path from node s to q is $SP(s,q,f^s,f^q)$, so $SP(s,q,f^s,f^q)=SP(s,k,f^s,f^k)+SP(k,j,f^k,f^j)+SP(j,q,f^j,f^q)$. If $SP(k,j,f^k,f^j)$ is not the shortest path from k to j using transportation modes with priority levels from $f^k$ to $f^j$, there must exist another shortest-path named $\overline{SP}(k,j,f^k,f^j)$. In this case, $\overline{SP}(s,q,f^s,f^q)=SP(s,k,f^s,f^k)+\overline{SP}(k,j,f^k,f^j)+SP(j,q,f^j,f^q)$ is the new shortest-path from node s to q with a lower transportation cost, which leads to a contradiction.

Therefore, the shortest-path algorithm under the multimodal traffic network proposed in this study must be able to provide the shortest-path with the minimum time from the departure point to the destination point if the destination point is reachable.    □

For the shortest-path algorithm under a multimodal traffic network, we assume that the number of nodes is $\left|V_1\right|,\left|V_2\right|\dots ,\left|V_{\left|F\right|}\right|$ and the number of arcs is $\left|E_1\right|,\left|E_2\right|\dots ,\left|E_{\left|F\right|}\right|$ in subnetwork $N_1,N_2\dots ,N_{\left|F\right|}$, respectively. We also let the average number of arcs per node be $k={\sum }_{i=1}^{\left|F\right|}{E}_i/{\sum }_{i=1}^{\left|F\right|}{V}_i$. The calculation formula of the time complexity based on the small root binary heap operation is given as follows:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill & TimeComplexity\hfill \\ \hfill & =(\sum _{i=1}^{\left|F\right|}{V}_i-1)·(T_{Extract-min}+T_{Delete}+T_{Decrease-key}·k)\hfill \\ \hfill & =(\sum _{i=1}^{\left|F\right|}{V}_i-1)·(2·log\left(\sum _{i=1}^{\left|F\right|}{V}_i\right)+k·log\left(\sum _{i=1}^{\left|F\right|}{V}_i\right))\hfill \\ \hfill & =\sum _{i=1}^{\left|F\right|}{V}_i·(2+k)·log\left(\sum _{i=1}^{\left|F\right|}{V}_i\right)\hfill \\ \hfill & =(2·\sum _{i=1}^{\left|F\right|}{V}_i+\sum _{i=1}^{\left|F\right|}{E}_i)·log\left(\sum _{i=1}^{\left|F\right|}{V}_i\right)\hfill \\ \hfill & =(\sum _{i=1}^{\left|F\right|}{V}_i+\sum _{i=1}^{\left|F\right|}{E}_i)·log\left(\sum _{i=1}^{\left|F\right|}{V}_i\right)\hfill \end{array}\end{array}$$

(49)

4.3. Extensions of the Algorithm

In Section 4.2, the transport time is used as the objective function to obtain the optimal transport path for the single task. However, the time that we consider in the objective can still be extended to incorporate different timing aspects in practice, from the departure point to the destination point.

In practical route planning for a single task, the transfer time between different transportation modes should be taken into account. Given the transfer time $T{F}^{f_1,f_2}$ ($f_1$ has a higher priority than $f_2$), the algorithm can be directly modified as follows: line 25 in Algorithm 1 must be modified to $C_{kv}=T{F}^{f_o,f_e}$. Therefore, even with transfer time considered, the complexity and optimality of Algorithm 1 do not change.

In most cases, the available transportation resources in emergency logistics are not enough to complete the transportation of each task at once. Therefore, the number of periods required to deliver all the batches should also be considered in addition to transport time on the route. The problem coincides with the arc-dependent shortest-path problem that has been proved to be NP-complete [35,36].

To approximately solve this problem, the minimum weighted time path algorithm under multimodal traffic networks is proposed based on the shortest-path algorithm. The algorithm should be modified as follows: line 33 in Algorithm 1 must be modified to $C_{kv}=t{r}_{kv}^{f_o}+R/r_{kv}^{f_o}$. As a weighted amount of transport time, $R/r_{kv}^{f_o}$ states that if its value is smaller (i.e., $r_{kv}^{f_o}$ is larger), the probability of selecting this route to transport is higher. Then, we can obtain the optimal path with this minimum weighted time.

Therefore, after obtaining the minimum weighted time path, we must recalculate the transportation time according to the traffic time on the road, traffic capacity and the number of batches. To this end, an algorithm to calculate the real transport time is proposed.

The primary idea of the algorithm is described as follows. After the minimum weighted time path is obtained, the threshold value on the path (i.e., the minimum number of batches that can pass per period) is calculated and is considered to be the number of batches that can be delivered per period for the task. The continuous transport periods and total route transport time are also calculated, and the sum of these two parts is the real transport time.

The pseudocode of the real transport time calculation algorithm is shown in Algorithm 2.

Algorithm 2: Real Transport Time Calculation Algorithm

5. Solution Method for the Optimization Model of Multitask Emergency Logistics Planning

In Section 4, we design a solution method to effectively solve the single-task route planning problem. In practice, multitask emergency logistics planning is more general and complex problem. The general multitask emergency logistics planning optimization model is proposed in Section 3. To address this problem in this section, we design a general optimization algorithm for multitask emergency logistics planning aiming at the shortest task completion time, based on routes recommended for each task using algorithm discussed in Section 4. Based on the general framework of PSO discussed in detail in Section 5.1, we improve the model with a two-dimensional array encoding method and introduce the crossover operation of the genetic algorithm (GA) in Section 5.2. Finally, the steps of the general optimization algorithm for multitask emergency logistics planning are presented in Section 5.3.

5.1. Particle Swarm Optimization

PSO was first proposed by [37] and is an evolutionary computing or bionic optimization algorithm that is inspired by the foraging of birds. The basic idea of PSO is as follows:

• for the solution space of the optimization problem, a group of particles are randomly initialized as the feasible solution of the objective function, and the fitness function is defined to describe the advantages and disadvantages of each solution;
• then, each particle moves randomly in the solution space, while the particle’s position is determined by the particle’s velocity. The incumbent solution found by the particle itself and found by the particle swarm will affect the particle velocity, respectively;
• all particles are searched step-by-step in this manner to find a high-quality solution in the solution space.

We let the vector $p_i=(p_1,p_2\dots ,p_n)$ represent the position and the vector $v_i=(v_1,v_2\dots ,v_n)$ represent the velocity of particle i in space. Each particle will have an adaptive value determined by the objective function, and the particle itself will know its current position as well as the incumbent solution that it has found to date. Each particle also knows where the other particles in the swarm have found the incumbent best position, and the particle compares its own experience with the best experience of the other particles to determine the next direction of movement.

The position of the incumbent solution found by the particle so far is:

$$Pbes{t}_i=(Pbes{t}_{i1},Pbes{t}_{i2},\dots ,Pbes{t}_{in})$$

(50)

The location of the incumbent solution of the entire particle swarm found up to now is denoted as:

$$Gbest=(Gbes{t}_1,Gbes{t}_2,\dots ,Gbes{t}_n)$$

(51)

When iterated to the $k+1$ generation, the position of each particle is updated by the following formulae:

$$\begin{array}{c}{v}_{id}^{k+1}=\omega ·v_{id}^k+c_1·r_1·(Pbes{t}_{id}^k-p_{id}^k)+c_2·r_2·(Gbes{t}_d^k-p_{id}^k)\end{array}$$

(52)

$$\begin{array}{c}{p}_{id}^{k+1}=p_{id}^k+v_{id}^{k+1}\end{array}$$

(53)

where $c_1,c_2$ are learning factors, $\omega$ is the inertia factor, and $r_1,r_2$ are random numbers between $[0,1]$.

5.2. Two-Dimensional Array Encoding and Cross Operation

The crossover operation is a common strategy in GA. Because the optimization of multitask emergency logistics planning is a complex combinatorial optimization problem, it is necessary to introduce an exchange operator into PSO to further enhance the search capability of the algorithm.

We assume that emergency tasks on the task list have a starting transport order, and each task has n alternative transport routes; thus, we construct the solution as a two-dimensional array $Sol(order,route)$, where $Sol\left(order\right)$ denotes the starting transport order for each task (i.e., $tas{k}_i$ begins transport in $ith$ among all tasks); and $Sol\left(route\right)$ represents the route number chosen for $tas{k}_i$ and $rout{e}_i\in \left\{\begin{array}{c}1,2\dots ,n\end{array}\right\}$. The array of solutions can be expressed as:

$$\begin{array}{c}\hfill \begin{array}{c}\hfill Sol=\left\{\begin{array}{cc}\hfill & Sol\left(order\right)=(tas{k}_1,tas{k}_2,\dots ,tas{k}_{\left|M\right|})\hfill \\ \hfill & Sol\left(route\right)=(rout{e}_1,rout{e}_2,\dots ,rout{e}_{\left|M\right|})\hfill \end{array}\right.\end{array}\end{array}$$

(54)

Now, we introduce the concepts of commutators and commutator series. The commutator is defined as $s=Swap(i,j)$, indicating that the commutator s acts on the sequence p to interchange elements at position i with position j. In PSO, the commutator can be viewed as the velocity of the particle, which can change the position of p. The commutator series $ss$ is defined as a sequential set of commutators ($ss=(Swa{p}_1,Swa{p}_2,\dots )$), which represents multiple operations on the sequence p. We also introduce three different operators ⊖, ⊕ and ⊗ as follows:

• operator ⊖ retrieves the commutator series $ss$ required to transform the latter sequence into former sequence. For instance, given that sequence $a=(a_1,a_2,\dots ,a_{n-1},a_n)$ and $a^{\prime }=(a_1,a_{n-1},\dots ,a_n,a_2)$, then $ss=(Swa{p}_1,Swa{p}_2,\dots )=a\ominus a^{\prime }$ means that an obtained commutator series $ss$ can transfer sequence $a^{\prime }$ to sequence a;
• operator ⊕ merges two commutator series into one commutator series, e.g., $(Swa{p}_1,Swa{p}_2)\oplus \left(Swa{p}_3\right)=(Swa{p}_1,Swa{p}_2,Swa{p}_3)$;
• operator ⊗ operates on a sequence according to the commutator series $ss$.

According to the above definition, the particle velocity and the update formula are modified as follows:

$$\begin{array}{c}{v}_{id}^{k+1}=c_1·r_1·(Pbes{t}_{id}^k\ominus p_{id}^k)\oplus c_2·r_2·(Gbes{t}_d^k\ominus p_{id}^k)\end{array}$$

(55)

$$\begin{array}{c}{p}_{id}^{k+1}=p_{id}^k\otimes v_{id}^{k+1}\end{array}$$

(56)

where the velocity update formula (55) indicates that there is a commutator series $ss$ that transforms $p_{id}^k$ into $Pbes{t}_{id}^k$ after $ss$. The position update formula (56) means that $p_{id}^k$ obtains a new sequence $v_{id}^{k+1}$ after the action of the commutator series $p_{id}^{k+1}$. $r_1$ and $r_2$ represent the probability of performing a commutator, equivalent to the operation of randomly discarding a portion of the commutators. Because $Sol(order,route)$ is a two-dimensional array, the crossover operation must operate on the starting transport order array and the route number array.

5.3. General Optimization Algorithm for Multitask Emergency Logistics Planning

The primary goal of the general optimization algorithm for multitask emergency logistics planning is to initialize the population parameters and positions, including the size of the population and the velocity and position of each particle (i.e., randomly generate the starting transport order array and route number array of tasks). Then, the particle fitness value is calculated. Next, the velocity vector and position vector of each particle are updated. Finally, the local position vector and global position vector are updated. We loop through this process until a certain number of iterations are reached, and then we terminate the algorithm.

The pseudocode of the general optimization algorithm of multitask emergency logistics planning is shown in Algorithm 3.

Algorithm 3: General Optimization Algorithm of Multitask Emergency Logistics Planning

6. Computational Results

To verify the efficacy of the model and algorithm, the computational results of the instances are reported in this section. We implemented the single-task route planning approach and multitask emergency logistics planning approach in Python 3.8. All computations were performed on a computer with an Intel(R) Core(TM) i7-1065G7 running at 1.50 GHz and 16.0 GB RAM.

6.1. Test Set

Due to the novelty of the problem and the complexity of the problem setting, there are no benchmark instances available in the literature, and we generate instance sets randomly to test the algorithms. In each instance of the test set, the following elements are included:

• city node information contains the periodic loading and unloading capacity in different transportation modes of each city node. We consider the different numbers of nodes in the test set, ranging from 20 to 65. The periodic loading and unloading capacity range from $[1,10]$, $[5,15]$, $[5,10]$ in air, rail and highway, respectively;
• primary arcs data describe the structure of different traffic networks; the periodic transport capacity; time consumption; and distance of each arc. We primarily consider air, rail and highway three modes of transportation networks and generate these three modes of transportation network topologies based on the number of nodes $\left|V\right|$ in each instance. The number of arcs $\left|E\right|$ in the air, rail and highway modes of transportation network ranges from $[1,0.5|V\left|\right]$, $\left[0.6\right|V|,0.9|V\left|\right]$, $\left[\right|V|,1.5|V\left|\right]$, respectively. The distance of each arc ranges from $[200\mathrm{km},400\mathrm{km}]$, $[150\mathrm{km},300\mathrm{km}]$, and $[50\mathrm{km},200\mathrm{km}]$ by air, rail and highway, respectively. the time consumption of each arc is the ratio of distance to the speed of its mode of transportation. The periodic transport capacity of each arc ranges from $[1,5]$, $[5,20]$, and $[4,15]$ by air, rail and highway, respectively;
• task lists include the specific information of each task, including its designated transportation mode, departure point, destination point, higher priority task and alternative paths. The three alternative routes for each task are obtained under different circumstances using the solution method presented in Section 4. Some tasks, such as specifying air transport, are point-to-point transport and have only one path; thus, multiple alternative paths cannot be recommended. $\left|M\right|$ ranges from $\left[\right|V|-2,|V|+10]$ in the first 20 instances, and from $\left[2\right|V|,3|V\left|\right]$ in the last 10 instances.

The results include the transportation route, departure time, arrival time, transfer node, periodic transport batches, the objective value Z and CPU runtime. The objective value Z refers to the transportation completion time of the last task.

6.2. Results and Comparison

6.2.1. Results of Instances with Small $\left|M\right|$

In the first experiment, we set the population size parameter to $size=100$, the number of iterations is $ite{r}_{num}=50$, and the learning factors are $c_1=c_2=1$ and $r_1=0.7$, $r_2=0.8$. We set the loading time and unloading time of each node to 6 hours. To achieve the smallest objective value, the calculation results for instances 1 to 20 are shown in

Table 3. Computational results of instances 1 to 20.

Instance	$\left|\mathit{V}\right|$	$\left|\mathit{E}\right|$	$\left|\mathit{M}\right|$	Z (Period)	CPU Runtime (s)
instance1	20	110	25	10	36.335
instance2	20	122	24	13	48.913
instance3	25	130	30	8	52.650
instance4	25	138	24	6	51.160
instance5	30	162	32	8	65.500
instance6	30	182	29	7	48.742
instance7	35	210	36	7	78.947
instance8	35	198	41	10	91.720
instance9	40	228	52	11	103.285
instance10	40	292	52	11	157.455
instance11	45	306	61	15	158.226
instance12	45	336	63	17	163.906
instance13	50	328	50	9	163.127
instance14	50	320	59	23	183.715
instance15	55	360	59	8	183.043
instance16	55	366	57	10	155.145
instance17	60	376	73	12	192.868
instance18	60	378	72	13	192.188
instance19	65	478	72	13	227.981
instance20	65	452	84	10	227.794

. $\left|M\right|$ in instance 1 to 20 ranges from 34 to 84, which is small compared to the number of nodes $\left|V\right|$. The detailed result of each instance is an emergency logistic plan and includes a lot of details; thus, we do not list the results of each instance. Instead, we show instance 1’s results in

Table 4. Emergency logistic plan of instance1.

Task Name	Departure Point	Arrival Point	Departure Day	Arrival Day	Batch Number	Periodic Batch Arrangement	Transport Route
task1	Node2	Node1	3	3	2	[2]	Node2–C0(air)–Node1
task2	Node20	Node17	0	3	15	[4, 4, 4, 3]	Node20–A0(air)–Node17
task3	Node14	Node13	0	10	10	[1, 1, 1, 1, 1, 1, 1, 1, 1, 1]	Node14–A279(air)–Node13
task4	Node13	Node14	3	5	2	[1, 1]	Node13–A279(air)–Node14
task5	Node16	Node11	3	4	4	[3, 1]	Node16–A7(air)–Node11
task6	Node16	Node6	0	4	23	[6, 6, 6, 5]	Node16–K545(rail)–Node19–G664(rail)–Node10–Y06(road)–Node6
task7	Node17	Node1	3	10	29	[5, 5, 3, 3, 3, 5, 5]	Node17–K42(rail)–Node11–G3(rail)–Node3–G53(rail)–Node20–X411(road)–Node9–X8(road)–Node1
task8	Node11	Node16	3	8	30	[1, 4, 8, 8, 8, 1]	Node11–G3(rail)–Node3–G53(rail)–Node20–K344(rail)–Node15–G01(rail)–Node13–K4(rail)–Node19–K545(rail)–Node16
task9	Node9	Node6	3	4	12	[8, 4]	Node9–G6(road)–Node6
task10	Node13	Node2	0	4	21	[5, 5, 5, 5, 1]	Node13–K4(rail)–Node19–G664(rail)–Node10–G749(rail)–Node2
task11	Node9	Node8	0	4	18	[3, 3, 3, 5, 4]	Node9–K1(rail)–Node20–K344(rail)–Node15–G01(rail)–Node13–K4(rail)–Node19–K545(rail)–Node16–K29(rail)–Node8
task12	Node10	Node11	2	6	20	[2, 7, 7, 4]	Node10–G664(rail)–Node19–K4(rail)–Node13–G01(rail)–Node15–Y225(road)–Node20–G6(road)–Node3–Z03(road)–Node11
task13	Node11	Node5	0	3	21	[7, 7, 7]	Node11–G3(rail)–Node3–G53(rail)–Node20–K344(rail)–Node15–G01(rail)–Node13–S1(road)–Node5
task14	Node8	Node5	3	8	22	[4, 4, 4, 4, 4, 2]	Node8–S674(road)–Node12–Y566(road)–Node15–S7(road)–Node13–S1(road)–Node5
task15	Node11	Node6	0	4	17	[5, 5, 5, 2]	Node11–G3(rail)–Node3–G53(rail)–Node20–K344(rail)–Node15–G01(rail)–Node13–K4(rail)–Node19–K545(rail)–Node16–K29(rail)–Node8–G993(rail)–Node18
task16	Node7	Node16	0	0	5	[5]	Node7–D60(rail)–Node19–K545(rail)–Node16
task17	Node8	Node9	3	7	24	[5, 5, 5, 5, 4]	Node8–K29(rail)–Node16–K545(rail)–Node19–K4(rail)–Node13–G01(rail)–Node15–K344(rail)–Node20–K1(rail)–Node9
task18	Node11	Node4	0	2	16	[6, 6, 4]	Node11–S880(road)–Node17–Z786(road)–Node4
task19	Node2	Node11	3	5	27	[11, 11, 5]	Node2–K2(rail)–Node1–D62(rail)–Node9–K1(rail)–Node20–G53(rail)–Node3–G3(rail)–Node11
task20	Node20	Node17	0	3	10	[2, 2, 4, 2]	Node20–G6(road)–Node3–Z03(road)–Node11–S880(road)–Node17
task21	Node4	Node1	0	5	30	[6, 6, 6, 6, 6]	Node4–K15(rail)–Node7–D60(rail)–Node19–G664(rail)–Node10–Z601(road)–Node2–S0(road)–Node1
task22	Node7	Node9	3	5	9	[6, 3]	Node7–G466(rail)–Node13–G01(rail)–Node15–Y225(road)–Node20–X411(road)–Node9
task23	Node3	Node1	3	5	14	[6, 6, 2]	Node3–G6(road)–Node20–X411(road)–Node9–X8(road)–Node1
task24	Node15	Node3	0	3	16	[5, 5, 3, 3]	Node15–Y225(road)–Node20–G6(road)–Node3
task25	Node11	Node13	0	4	15	[2, 2, 4, 4, 3]	Node11–G1(road)–Node18–Z5(road)–Node19–S08(road)–Node13

as an example.

Table 3 shows that the CPU runtime goes up with the increasing number of nodes $\left|V\right|$, arcs $\left|E\right|$ and $\left|M\right|$. Due to the complexity of the problem and algorithm, the CPU runtime is large. Moreover, the objective value Z is not directly related to the number of nodes $\left|V\right|$, arcs $\left|E\right|$ and $\left|M\right|$.

6.2.2. Results of Instances with Large $\left|M\right|$

In this experiment, the computational results for instances 21 to 30 are shown in

Table 5. Computational results for instances 21 to 30.

Instance	$\left|\mathit{V}\right|$	$\left|\mathit{E}\right|$	$\left|\mathit{M}\right|$	Z (Period)	CPU Runtime (s)
instance21	20	120	58	18	71.536
instance22	25	140	52	28	84.455
instance23	30	178	73	18	130.854
instance24	35	222	103	26	204.818
instance25	40	312	104	18	237.317
instance26	45	258	121	25	246.684
instance27	50	292	140	28	311.662
instance28	55	348	137	32	298.705
instance29	60	386	158	30	336.258
instance30	65	358	150	56	397.488

. The difference from the last experiment is that $\left|M\right|$ in instances 21 to 30 ranges from 52 to 158, which is relatively large compared to the number of nodes $\left|V\right|$. We also set the population size parameter as $size=100$, the number of iterations is $ite{r}_{num}=50$, and the learning factors are $c_1=c_2=1$ and $r_1=0.7$, $r_2=0.8$.

Comparing Table 3 and Table 5 shows that the CPU runtime is closely related to $\left|M\right|$ with the same number of nodes $\left|V\right|$ and arcs $\left|E\right|$. Concurrently, an increase in $\left|M\right|$ leads to more conflicts between different tasks, and the objective value Z also increases.

6.2.3. Sensitivity Analysis

In this section, we discuss the impact of the different parameters on the performance of the presented algorithm.

First, to improve computational efficiency, we also optimize the parameters of the PSO. The parameters that have the greatest influence on CPU runtime of the algorithm are population size and iteration times. With a fixed population size, we try to find a value of $ite{r}_{num}$ that leads to a satisfying solution without sacrificing the overall run time. We set the population size parameter as $size=100$ and the learning factors $c_1=c_2=1$ and $r_1=0.7$, $r_2=0.8$. The CPU runtime and the objective value of each instance under different numbers of iterations are shown in

Table 6. CPU runtime and objective value under different iteration times.

	${\mathit{iter}}_{\mathit{num}}=50$		${\mathit{iter}}_{\mathit{num}}=30$		${\mathit{iter}}_{\mathit{num}}=20$		${\mathit{iter}}_{\mathit{num}}=10$		${\mathit{iter}}_{\mathit{num}}=5$
Instance	CPU Runtime (s)	Z (Period)	CPU Runtime (s)	Z (Period)	CPU Runtime (s)	Z (Period)	CPU Runtime (s)	Z (Period)	CPU Runtime (s)	Z (Period)
instance1	36.335	10	21.398	10	20.265	10	13.247	10	7.631	10
instance3	52.650	8	27.013	8	19.951	8	13.618	9	11.956	9
instance5	65.500	8	53.363	8	30.317	8	22.275	10	17.18	10
instance7	78.947	7	48.231	7	33.194	7	21.05	9	15.035	9
instance9	103.285	11	74.073	11	55.746	11	31.584	12	22.711	13
instance11	158.226	15	95.339	15	66.237	15	41.694	15	29.536	18
instance13	163.127	9	98.162	9	70.314	9	43.009	9	29.239	12
instance15	183.043	8	116.62	8	76.799	10	45.784	10	30.998	9
instance17	192.868	12	115.921	12	85.789	13	50.901	13	35.105	14
instance19	227.981	13	145.831	13	102.033	14	60.927	14	45.779	14
instance21	71.536	18	58.663	18	32.347	18	21.355	18	14.512	20
instance23	130.854	18	83.617	18	53.516	18	34.96	19	23.624	20
instance25	237.317	18	164.489	21	115.582	22	66.999	21	45.646	22
instance27	311.662	28	196.318	34	136.812	35	86.469	35	59.745	38
instance29	336.258	30	228.203	31	176.163	33	107.372	33	76.182	37

. From the results in Table 6, it is clear that the solutions improve, and computation times increase, when $ite{r}_{num}$ increases. With the population size fixed at $size=100$, we find that for some simple multimodal traffic networks and fewer tasks, lowering $ite{r}_{num}$ to 20 can still lead to solutions with the same objective value. However, when we must manage a complex multimodal traffic network with many transportation tasks, $ite{r}_{num}$ need to be set over 50 to ensure convergence to near the optimal value.

To demonstrate the impact of the learning factors $c_1$ and $c_2$, we consider instance 10 as an example. We set the population size parameter as $size=100$, and the number of iterations $ite{r}_{num}=20$, $r_1=0.7$ and $r_2=0.8$.

Table 7. CPU runtime and objective value of instance10 with different $c_1$, $c_2$.

${\mathit{c}}_1$, ${\mathit{c}}_2$	Z (Period)	CPU Runtime (s)
$c_1=c_2=0.5$	11	60.203
$c_1=c_2=1$	11	56.518
$c_1=c_2=1.5$	12	57.916
$c_1=c_2=2$	15	59.642
$c_1=c_2=2.5$	17	58.53
$c_1=c_2=3$	16	59.592
$c_1=c_2=3.5$	18	57.509
$c_1=c_2=4$	26	67.2

shows the different results of instance 10 with different learning factors. We find that different learning factors have a strong impact on the calculation results. When $c_1$, $c_2$ are larger, the objective value Z of the satisfactory solution is larger. The reason for this situation is that when the learning factor is small, the particles will wander outside the target region. If the learning factor is large, the particles will cross the target region and form a large adaptive value fluctuation. However, for computing efficiency, the learning factor does not affect the CPU runtime.

6.3. Discussion

By considering multimodal traffic networks with different numbers of nodes V, arcs E and different numbers of tasks, we demonstrate the efficacy of the proposed method in solving multitask emergency logistics planning problems. However, there may be a different plan for multitask emergency logistics with the same objective value Z. Even the optimal solution may also not be unique, which also provides multiple references for the plan-maker when making plans.

From the multitask emergency logistics planning experiment, we find that there are three main factors affecting the computational efficiency and accuracy of calculation results. The first factor is the arrival order of transportation. Considering the special needs of some tasks in emergency logistics, and due to the different functions of each task, there is often a certain arrangement. Therefore, in the task list, we consider the arrival order of different tasks; thus, the final calculation result of objective value Z will also be affected by this issue.

The second factor that affects the computational results is the formulation of alternative routes. In the formulation of alternative routes for each task, we use the solution method for single-task route planning in Section 4 to provide three alternative routes considering different conditions. Because the transportation routes given for individual task are optimal, the diversity of route planning may be lost. If the alternative route is given artificially, although it is not optimal with respect to single-task transportation, it may achieve a better result in multitask emergency logistics planning.

The last factor is the different parameters of the proposed algorithm. Considering the scale of the problem tackled in this study and the complexity of the proposed algorithm, we reduce the number of iterations of the PSO algorithm to improve computational efficiency. Although it is possible that the computational results will not converge to the optimal value, it is still acceptable that the results will converge closer to the optimal value. Different learning factors also strongly affect the results. To find a satisfactory solution in a reasonable run time, the value of the learning factors should be $c_1=c_2\le 1$.

In related research, Ergün et al. [38] also equally committed to emergency logistics planning in natural disasters. They solved the game theoretical model by cooperative game theory. The experimental results show that it can improve the execution efficiency of emergency logistics management after natural disasters. However, because it essentially solves the problem of single-point to single-point multimodal transportation, its practicability and generality is not ideal. In comparison to the work by Wang et al. [39], they proposed a heuristic algorithm based on nested partitions to solve the deterministic dynamic emergency logistics planning model. By solving linear programming relaxation subproblems and iteratively fixing integer variables, a high-quality feasible solution is finally obtained. However, the complexity of the model is relatively not high and the descriptions of emergency logistics constraints are not sufficient. Sabouhi et al. [40] studied the multi-point to multi-point emergency logistics planning problem. Compared with their work, the advantages of ours are that our model considers multimodal transportation and solves such problem through proposed algorithms. Although a lot of work has been done in related research, the actual constraints of the model are different due to the different problems considered. Moreover the case data used are also disparate, so it is tough to compare the calculation results at different level.

Overall, the MIP model of multitask emergency logistics planning under multimodal proposed in this paper is more general than the existing research, since it satisfies more practical constraints such as transportation mode priority constraints, transport time and space continuity constraints for each task, emergency principle constraints, etc. Through the two-layer solution process is obviously easy to operate and can effectively solve practical multitask emergency logistics planning problems as well. Although it may not be possible to find the optimal solution in a limited time, we can find the satisfactory solution within acceptable time to the decision maker, which is more in line with the actual needs of emergency logistics.

7. Conclusions

Multitask emergency logistics planning under multimodal transportation is important in emergency management. In this paper, we consider multitask emergency logistics planning under different modes of traffic networks. A MIP model for the multitask emergency logistics planning problem is established. Considering the complexity of this model, we decompose it into two layers for optimization. To manage the single-task route recommendation as the inner layer, we develop a shortest-path algorithm with multimodal traffic networks. We present the optimal substructure of the algorithm and its time complexity. Then, with the routes of each task recommended by the single-task solver, a general optimization algorithm based on improved PSO is proposed to coordinate the execution of each task constrained by the limited transportation capacity at the outer layer, which solves the multitask emergency logistics planning problem.

Computational results demonstrate that the proposed approach is effective and practical for solving multitask emergency logistics planning problems, without costing too much time and computing effort to struggle for the optimal solution. The two-layer solution process we proposed can effectively solve multitask emergency logistics planning problems under multimodal transportation. The result of this research assists emergency logistics management organizations, such as natural disaster management organization for a better tactical and operational reaction in a critical situation. Furthermore we summarize three main factors that affect the computational efficiency and accuracy of calculation results, which can be helpful for practitioners when the proposed method is needed in this paper to solve the emergency logistics planning in the future.

We believe that this study opens up several directions for future research. First of all, an effective exact algorithm should be studied and developed, as the solution method designed in this work is a heuristic procedure. On the other side, we consider a fixed multimodal traffic network structure in emergency logistics during our work. However, it can incur more difficulty if the network structure is dynamic or uncertain, which is common and should be considered. So for future works, emergency logistics planning under uncertainty such as interval, polyhedral, fuzzy or stochastic, can be further study in these fields [40,41]. Refer to Li et al. [42] and Zhou et al. [43], the potential research areas on the topic of COVID-19-induced emergency air transportation can also be considered.