Dependence of the transportation time on the sequence in which particles with different hopping probabilities enter a lattice

TL;DR

This paper studies how the order of particle entry affects transportation time in a modified TASEP model, revealing that sequence optimization can reduce time and that sometimes random sequences are preferable when considering sorting costs.

Contribution

It introduces a new multi-species TASEP model considering entry sequence effects and analyzes how sequence optimization impacts transportation efficiency.

Findings

Group sequences can minimize transportation time compared to random sequences.

Including sorting costs can make random sequences more advantageous.

Theoretical analysis supports simulation results for specific cases.

Abstract

Smooth transportation has drawn the attention of many researchers and practitioners in several fields. In the present paper, we propose a modified model of a totally asymmetric simple exclusion process (TASEP), which includes multiple species of particles and takes into account the sequence in which the particles enter a lattice. We investigate the dependence of the transportation time on this `entering sequence' and show that for a given collection of particles group sequence in some cases minimizes the transportation time better than a random sequence. We also introduce the `sorting cost' necessary to transform a random sequence into a group sequence and show that when this is included a random sequence can become advantageous in some conditions. We obtain these results not only from numerical simulations but also by theoretical analyses that generalize the simulation results for some special cases.

Tables

Table 1. Table 1: Classification of Regions

Region No.	Range
1
2
3
4	$1-\sqrt{1-p_S}\le \alpha$

Table 2. Table 2: Classification of subregions in Region 2

Subregion No.	Range
2–$1$
2–$2$
…	…
2–$v$
…	…
2–$\left(u-1\right)$

Table 3. Table 3: Classification of subregions in Region 3

Subregion No.	Range
3–$u$
3–$\left(u+1\right)$
…	…
3–$v$
…	…
3–$S$

Table 4. Table 4: Sign of f ​ ( α ) 𝑓 𝛼 f(\alpha) and d ​ f ​ ( α ) / d ​ α 𝑑 𝑓 𝛼 𝑑 𝛼 df(\alpha)/d\alpha

Region No.	$f\left(\alpha \right)$	$df\left(\alpha \right)/d\alpha$
1	$-$	U
2	U	U
3	U	$+$
4	$+$	0

Equations

$${\Delta }_T=\frac{T_{\mathrm{G}}-T_{\mathrm{R}}}{T_{\mathrm{R}}}.$$

$$\sum _{i=0}^2{P}_{i0}+\sum _{i=0}^2{P}_{i\ast }=1.$$

$$\begin{array}{rl}{Q}_2(p_1,p_2;r)& =p_1{P}_{10}+p_2{P}_{20}\\ & =\frac{p_1{p}_{2}A}{\{(1-r)p_1+r{p}_2\}A+p_1{p}_{2}B},\end{array}$$

$$A=\alpha \beta (\alpha +\beta -\alpha \beta )$$

$$B={\alpha }^2+{\beta }^2-{\alpha }^2\beta -\alpha {\beta }^2+\alpha \beta .$$

$$Q_1(p)=\frac{pA}{pB+A}.$$

$$p_{\mathrm{h}}=\frac{p_1{p}_2}{(1-r)p_1+r{p}_2}.$$

$$Q_1(p=p_{\min })<Q_S<Q_1(p=p_{\mathrm{h}}).$$

$$\displaystyle\left\{\begin{array}[]{ll}\alpha=1-\sqrt{1-p}\ \land\ 1-\sqrt{1-p}<\beta<1\\ \beta=1-\sqrt{1-p}\ \land\ 1-\sqrt{1-p}<\alpha<1\end{array}\right..$$

$${\Delta }_Q=\frac{Q_2-Q_1(p=p_{\mathrm{h}})}{Q_1(p=p_{\mathrm{h}})},$$

$$T_{\mathrm{R}}\approx \frac{N}{Q_S}$$

$$T_{\mathrm{G}}\approx \sum _{s=1}^S\frac{r_{s}N}{Q_1(p=p_s)}.$$

$$T_{\mathrm{H}}\approx \frac{N}{Q_1(p=p_{\mathrm{h}})}$$

$$T_{\mathrm{M}}\approx \frac{N}{Q_1(p=p_{\min })},$$

$$T_{\mathrm{H}}<T_{\mathrm{R}}<T_{\mathrm{M}}.$$

$$f(\alpha ,\beta ;p_1,......,p_S;r_1,......,r_S)=T_{\mathrm{H}}-T_{\mathrm{G}}.$$

$$\displaystyle Q_{1}(p)=\left\{\begin{array}[]{ll}\displaystyle\alpha\frac{p-\alpha}{p-\alpha^{2}}&{\rm for\ LD\ phase.}\\ \displaystyle\frac{1}{2}(1-\sqrt{1-p})&{\rm for\ MC\ phase.}\end{array}\right.$$

$$T_{\mathrm{G}}\approx \sum _{s=1}^S\frac{r_{s}N(p_s-{\alpha }^2)}{\alpha (p_s-\alpha )}$$

$$T_{\mathrm{H}}\approx \frac{N(p_{\mathrm{h}}-{\alpha }^2)}{\alpha (p_{\mathrm{h}}-\alpha )}.$$

$$f(\alpha )=\frac{N(p_{\mathrm{h}}-{\alpha }^2)}{\alpha (p_{\mathrm{h}}-\alpha )}-\sum _{s=1}^S\frac{r_{s}N(p_s-{\alpha }^2)}{\alpha (p_s-\alpha )}.$$

$$f(\alpha )<0;$$

$$\displaystyle Q_{1}(p)=\left\{\begin{array}[]{ll}\displaystyle\frac{1}{2}(1-\sqrt{1-p})&{\rm for}\ p=p_{1},......,p_{v-1},p_{\rm h}.\\ \displaystyle\alpha\frac{p-\alpha}{p-\alpha^{2}}&{\rm for}\ p=p_{v+1},......,p_{S}.\end{array}\right.$$

$$T_{\mathrm{G}}\approx \sum _{s=1}^v\frac{2r_{s}N}{1-\sqrt{1-p_s}}+\sum _{s=v+1}^S\frac{r_{s}N(p_s-{\alpha }^2)}{\alpha (p_s-\alpha )}$$

$$T_{\mathrm{H}}\approx \frac{N(p_{\mathrm{h}}-{\alpha }^2)}{\alpha (p_{\mathrm{h}}-\alpha )}.$$

$$\begin{array}{rl}f(\alpha )& \approx \frac{N(p_{\mathrm{h}}-{\alpha }^2)}{\alpha (p_{\mathrm{h}}-\alpha )}\\ & -\sum _{s=1}^v\frac{2r_{s}N}{1-\sqrt{1-p_s}}-\sum _{s=v+1}^S\frac{r_{s}N(p_s-{\alpha }^2)}{\alpha (p_s-\alpha )}.\end{array}$$

$$\underset{\alpha \to q_h-0}{\lim }\frac{df(\alpha )}{d\alpha }>0,$$

$$\displaystyle Q_{1}(p)=\left\{\begin{array}[]{ll}\displaystyle\frac{1}{2}(1-\sqrt{1-p})&{\rm for}\ p=p_{1},......,p_{v-1},p_{\rm h}.\\ \displaystyle\alpha\frac{p-\alpha}{p-\alpha^{2}}&{\rm for}\ p=p_{v},......,p_{S}.\end{array}\right.$$

$$T_{\mathrm{G}}\approx \sum _{s=1}^{v-1}\frac{2r_{s}N}{1-\sqrt{1-p_s}}+\sum _{s=v}^S\frac{r_{s}N(p_s-{\alpha }^2)}{\alpha (p_s-\alpha )}$$

$$T_{\mathrm{H}}\approx \frac{2N}{1-\sqrt{1-p_{\mathrm{h}}}}.$$

$$\begin{array}{rl}f(\alpha )& \approx \frac{2N}{1-\sqrt{1-p_{\mathrm{h}}}}\\ & -\sum _{s=1}^{v-1}\frac{2r_{s}N}{1-\sqrt{1-p_s}}-\sum _{s=v}^S\frac{r_{s}N(p_s-{\alpha }^2)}{\alpha (p_s-\alpha )}.\end{array}$$

Full text

††preprint: APS/123-QED

I INTRODUCTION

Transportation systems are key topics in social or biological systems Appert-Rolland et al. (2017). In social systems, researchers have sought to obtain smooth transportation in various situations, such as production flow Croom et al. (2000); Ezaki et al. (2015), vehicular traffic Taniguchi et al. (2015); Woelki (2013); Imai and Nishinari (2015); Yamamoto et al. (2017), and pedestrian evacuation Christensen and Sasaki (2008); Pelechano and Malkawi (2008); Yanagisawa et al. (2009); Ezaki et al. (2012). On the other hand, for biological systems, intracellular transportation along microtubules has been vigorously investigated Rothman (1994); Müller et al. (2008); Bressloff and Newby (2013).

Among various transportation models, the asymmetric simple exclusion process (ASEP), pioneered by MacDonald and Gibbs MacDonald et al. (1968); MacDonald and Gibbs (1969), has attracted much attention. It is a stochastic process on a one-dimensional lattice in which particles move asymmetrically. A derivative of ASEP, in which particles are allowed to hop unidirectionally (left to right in the present study) is called a totally asymmetric simple exclusion process (TASEP). In the field of nonequilibrium statistical mechanics, researchers have applied TASEP to various transportation problems, such as molecular-motor traffic Parmeggiani et al. (2003); Chowdhury et al. (2005); Chou et al. (2011); Appert-Rolland et al. (2015), vehicular traffic Helbing (2001); Imai and Nishinari (2015); Ito and Nishinari (2014); Tsuzuki et al. (2018); Woelki (2013); Yamamoto et al. (2017), and the exclusive-queuing process Yanagisawa et al. (2010); Arita and Yanagisawa (2010); Arita and Schadschneider (2015), especially since the TASEP with open boundary conditions has been solved exactly de Gier and Nienhuis (1999); Blythe and Evans (2007); Derrida et al. (1993).

In practice, researchers struggle to achieve smooth operation for various tasks, smooth logistics for various products, and an effective evacuation method for pedestrians in various situations, such as exit plans from sports stadiums and concert venues. To attain smooth flow in such situations, we often consider the sequence in which we perform tasks and pedestrians move because this sequence may affect the total performance of the systems. For example, slow pedestrians may block fast ones at the back of a narrow street, which worsens pedestrian flow. To investigate how the abovementioned sequences affect pedestrian flow in various systems, herein, we propose a modified TASEP comprising a finite number of multi-species particles, in which the entering sequences of the particles are considered.

In the proposed model, we consider the number of particles to be finite and study the transportation times of those particles. Note that we do not consider the steady state of the system itself. Minemura et.al Minemura et al. (2010) investigated the transportation time for a hopping probability that depends upon the lot size, using the single-species TASEP with a finite number of particles. Other related works Adams et al. (2008); Cook et al. (2009); Cook and Zia (2009, 2012) also adopted a finite number of particles. In those models, however, particles circulated through a system comprising a lattice and a particle pool while the input or output rate was varied.

Additionally, the concept of multiple particles has already been extensively studied  Angel (2006); Ferrari et al. (1994); Arita (2006); Mallick (1996); Derrida and Evans (1999); Ayyer et al. (2009); Privman (2005); Ferrari et al. (2007); Cantini (2008); Uchiyama (2008); Crampe et al. (2015); Aneva (2003); Arita et al. (2009); Duchi and Schaeffer (2005); Evans (1997); Krug and Ferrari (1996); Mallick et al. (1999); Crampe et al. (2016); Arita and Mallick (2013); Prolhac et al. (2009); Ayyer and Roy (2017); Arndt et al. (1998); Evans (1996); Sasamoto (2000). For example, second class particles were introduced in Refs. Angel (2006); Ferrari et al. (1994); Arita (2006); Mallick (1996); Derrida and Evans (1999); Ayyer et al. (2009); Privman (2005); Ferrari et al. (2007); Cantini (2008); Uchiyama (2008); Crampe et al. (2015) and more than two-species particles were introduced in Refs. Aneva (2003); Arita et al. (2009); Duchi and Schaeffer (2005); Evans (1997); Krug and Ferrari (1996); Mallick et al. (1999); Crampe et al. (2016); Arita and Mallick (2013); Prolhac et al. (2009); Ayyer and Roy (2017). However, most of these studies focused on mathematically exact solutions to the systems under consideration by using such as Matrix Product Ansatz and did not much consider the application of the studied model to real-world situations. Furthermore, owing to their simplicity, periodic-boundary conditions have been adopted in many studies Angel (2006); Arita et al. (2009); Derrida and Evans (1999); Arndt et al. (1998); Evans (1996, 1997); Krug and Ferrari (1996); Sasamoto (2000); Mallick et al. (1999); Prolhac et al. (2009); Arita and Mallick (2013). Studies on multi-species ASEP with open boundaries and random updating were undertaken only recently Arita (2006); Ayyer et al. (2009); Privman (2005); Uchiyama (2008); Crampe et al. (2016, 2015); Ayyer and Roy (2017); for example, Ref. Ayyer and Roy (2017) obtained the exact phase diagram for a multi-species (more than two-species) ASEP. The present investigation primarily focuses on the problem of minimizing the transportation time, adopting open-boundary conditions and parallel updating. With the same boundary conditions and updating rules as the present study, Ref. Bengrine et al. (1999) adopted particles with disorder, whereas jumping particles were introduced in Ref. Duchi and Schaeffer (2008). Note that majority related works considering multi-species particles assume that swapping between different types of particles, i.e., bidirectional particle hopping, can occur, whereas our model prohibits swapping 111 To be precise, if we regard vacant sites as particle 0, swapping between particle 0 and other kinds of particles (1, 2, 3,……) occurs. However, swapping between particles 1, 2, 3,…… does not occur. .

In contrast, to the best of our knowledge, no TASEP investigations that focus on the entering sequence of the particles (the key highlight of our model) have been reported thus far. Herein, we have considered two special types of sequences in particular: ‘random sequences’ and ‘group sequences,’ and we have compared the transportation times for these two types of sequences. In association with the entering sequence, we have introduced the sorting cost in our model. Without sorting, particles are typically transported at random, i.e., in a random sequence. Therefore, considering the cost of sorting particles from a random sequence into a group sequence is useful. In the present study, we define this sorting cost and compare the results obtained with and without sorting.

We have determined the dependence of the transportation time on the entering sequence of the particles from numerical simulations based on our model. Moreover, we find that the optimal sequence can vary, depending upon choice of parameter set, when the sorting cost is considered. In addition, we have succeeded in obtaining mathematical proofs of the simulation results for some special cases.

The remainder of the present study is organized as follows. Section II describes the details of our proposed model and some important parameters, modifying the original TASEP. In Sec. III, we present and discuss the results of numerical simulations using the modified TASEP. Section IV presents theoretical analyses of the simulation results for some special cases. The paper concludes in Sec. V.

II MODEL DESCRIPTION

II.1 Original (single-species) TASEP with open-boundary conditions

The original TASEP with open-boundary conditions is defined as a one-dimensional lattice of $L$ sites, labeled from left to right $i=0,1,......,L-1$ (see Fig. 1). Each site can be either empty or occupied by a single particle. In the present study, we adopt discrete time steps and parallel updating. In parallel updating, the states of all the particles on the lattice are determined simultaneously in the next time step. Notably, we can use random updating, which is usually adopted in the ASEP; however, we intentionally adopt parallel updating in the present study (see the specific reasons in 222This system can be simulated also with random updating; however, we adopt parallel updating for the following three reasons. First, we consider that the results will be mainly applied to real systems, where particles and agents move simultaneously. It is generally natural to adopt parallel updating in the contexts of pedestrian dynamics or traffic flows Nagel and Schreckenberg (1992); Burstedde et al. (2001). Second, it is more natural to consider that time for swapping two particles is compatible with a one-time step with parallel updating than with random-sequential updating. Finally, the simulation times can be suppressed drastically with parallel updating compared to with random updating.). Particles enter the lattice from the left boundary with probability $\alpha$, and leave the lattice from the right boundary with probability $\beta$. In the bulk of the lattice, if the right-neighboring site is empty, a particle hops to that site with probability $p$; otherwise it remains at its present site. Our modified TASEP differs from this original one in the following four ways.

II.2 Difference 1: Finite number of particles

First, the number of particles $N$ is finite, as illustrated in Fig. 2. The system evolves until the $N$th particle leaves the lattice. We define the transportation time $T$ as the time gap between the start of the simulation and the time when the $N$th particle leaves the lattice.

II.3 Difference 2: Multi-species particles

Second, our model adopts multi-species particles, i.e., particles with different hopping probabilities, as illustrated in Fig. 3. Specifically, each of the $N$ particles is allocated to one of $S$ species, where $1\leq S\leq N$. Particles that belong to each species $s$ ($s=1,2,......,S$) all have the same hopping probability $p=p_{s}$ ($0<p_{s}\leq 1$). Note that with $S=1$ our model reduces to the single-species TASEP, whereas with $S=N$ all particles have different hopping probabilities. The fraction of all the $N$ particles allocated to each species $s$ is defined as $r_{s}$, obviously satisfying $\sum_{s=1}^{S}{r_{s}}=1$.

II.4 Difference 3: Consideration of entering sequence of particles

Third, we consider the sequence in which the particles enter the lattice (i.e., the ‘entering sequence’), which is the most important feature in our model. Specifically, particles form a queue before the left boundary and enter the lattice according to the sequence, as illustrated in Fig. 4. In the present study, we investigate two types of sequences: ‘random sequences’ and ‘group sequences,’ as illustrated in Fig. 5.

In a random sequence, particles line up randomly regardless of their hopping probabilities. A random sequence thus has $N!/\prod_{s=1}^{S}(r_{s}N)!$ patterns. Note that in real situations without any controls, random sequences can be assumed to occur spontaneously.

On the other hand, in a group sequence, particles form groups of the same species and line up group by group. There are $S!$ possible patterns of group sequences, which are clearly among the random sequences.

For the case $S=N$, where all hopping probabilities are different, we bunch the particles with similar hopping probabilities close together with each other as much as possible, imaginarily considering them as ‘continuous groups.’ Consistent with this idea, we define a group sequence with $S=N$ as either an ascending or a descending sequence. Note that we define such a sequence by considering the rightmost particle to be the first particle in the sequence.

We define the transportation times for the random and group sequences to be $T_{\rm R}$ and $T_{\rm G}$, respectively.

II.5 Difference 4: Introduction of the sorting cost

Finally, we introduce the cost of sorting the particles and investigate the effect of the sorting cost on the transportation time. Here, we define the sorting cost as the minimal number of exchanges $K(\tau_{\rm a},\tau_{\rm b})$ necessary to sort the particles form sequence $\tau_{\rm b}$ to sequence $\tau_{\rm a}$, where $\tau_{\rm a}$ and $\tau_{\rm b}$ represent the sequence after sorting and before sorting, respectively. Note that the arguments of $K(\tau_{\rm a},\tau_{\rm b})$ will be abbreviated in obvious cases.

In the present study, $\tau_{\rm a}$ ($\tau_{\rm b}$) correspond to $\tau_{\rm G}$ ($\tau_{\rm R}$), where $\tau_{\rm R}$ and $\tau_{\rm G}$ represent a random sequence and a group sequence, respectively. The sequence $\tau_{\rm G}$ can differ depending upon $\tau_{\rm R}$; that is, $\tau_{\rm G}$ is determined so that the number of exchanges is minimized for each $\tau_{\rm R}$. Figure 6 shows two examples for which $K(\tau_{\rm a},\tau_{\rm b})=2$ when $N=6$ and $S=3$. Note that we do not consider the distance between the exchanged particles.

We define the number of time steps necessary to sort the particles to be $\lambda K$, where the parameter $\lambda$ is the ratio of the sorting cost to number of TASEP time steps.

III SIMULATION RESULTS

In this section, we use numerical simulations to investigate the dependence of the transportation time on the entering sequence of the particles.

In all the simulations below, we set $L=200$ and $N=10,000$; we validate this selection of $L$ an $N$ in Appendix A. We determine the value of $T$ for each parameter, and average $T$ over 100 trials for Fig. 7 and over 10 trials for Figs. LABEL:fig:deltaT, LABEL:fig:lambda, and LABEL:fig:withsort).

III.1 Without sorting cost ($\lambda=0$)

In this subsection, we set $\lambda=0$, i.e., we do not include the sorting cost.

In Fig. 7 we plot the simulation values of the number of particles that have not yet exited the lattice at time $t$ for $S\in\{2,3,N\}$. We fix $(\alpha,\beta)=(1,1)$ for (a)–(c) and $(\alpha,\beta)=(0.1,0.2)$ for (d)–(f). In the figures, we refer to the number of particles that have not yet exited the lattice at time $t$ simply as the ‘remaining particles.’ The simulation starts at $t=0$, and the number of particles becomes 0, i.e., the $N$th particle exits the lattice, at $t=T$.

We note two important phenomena in Figs. 7 (a)–(c). First, surprisingly, $T_{\rm G}$ is smaller than $T_{\rm R}$ for all three values of $S$ when $\alpha=\beta=1$. This result implies that the group sequences yield smoother transportation than the random ones for the cases $(\alpha,\beta)=(1,1)$. Second, $T_{\rm G}$ seems not to depend upon the order of each group in the group sequence, which can take $S!$ possible patterns.

On the other hand, in Figs. 7 (d)–(f), unlike the cases in Figs. 7 (a)–(c), the difference between $T_{\rm R}$ and $T_{\rm G}$ seems to vanish.

In order to compare the difference between $T_{\rm R}$ and $T_{\rm G}$ for various $(\alpha,\beta)$, we define $\Delta_{T}$ as the ratio of the change from $T_{\rm R}$ to $T_{\rm G}$; that is,

[TABLE]

From this definition of $\Delta_{T}$, $\Delta_{T}<0$ $(\Delta_{T}>0)$ indicates that group (random) sequences are preferable for smooth transportation. Note that in the following, to calculate $\Delta_{T}$ we assume that each group in a group sequence is arranged in ascending order in terms of species number $s$.

The simulation values of $\Delta_{T}$ for various $(\alpha,\beta)$ with (a) $S=2$, (b) $S=3$, and (c) $S=N$ are plotted in Fig. LABEL:fig:deltaT. Note that the black lines represent the boundaries between the low-density/high-density (LD/HD) and the maximal current (MC) phases of the single-species TASEP with hopping probability $p_{1}$ (boundary A) in Fig. LABEL:fig:deltaT (a), and $p_{1}$ (boundary B1) and $p_{3}$ (boundary B3) in Fig. LABEL:fig:deltaT (b), respectively.

Figure LABEL:fig:deltaT shows that for all three values of $S$, $\Delta_{T}$ is small in the region where ${\rm min}(\alpha,\beta)$ is relatively large. [In Fig. LABEL:fig:deltaT (b), $\Delta_{T}$ finally yields to a constant value in the upper-right region beyond boundary B3.] On the other hand, $\Delta_{T}$ is small in the region where ${\rm min}(\alpha,\beta)$ is relatively small. [In Fig. LABEL:fig:deltaT (a) and (b), $\Delta_{T}$ is almost 0, especially in the lower-left region beyond the boundary A or B1.] Here, we term the region with $\Delta_{T}<0$ as the ‘group-advantageous region’ ($T_{\rm R}>T_{\rm G}$), whereas we designate the region with $\Delta_{T}\approx 0$ as a ‘neutral region’ ($T_{\rm R}\approx T_{\rm G}$), if it exists.

These results indicate that group sequences can make transportation smoother than random sequences when the system is mainly governed by the bulk region of the lattice, but the dependence on the type of sequences vanishes (or decreases) when the system is mainly governed by the boundaries.

III.2 With sorting cost ($\lambda>0$)

In this subsection, we consider the sorting cost by varying $\lambda$ for the same parameter sets in the previous subsection. Appendix B presents specific schemes for obtaining the minimal number of exchanges necessary to sort the particles in the simulations.

Figure LABEL:fig:lambda plots $\Delta_{T}$ for (a) $S=2$, (b) $S=3$, and (c) $S=N$ as functions of $\lambda$ for various $(\alpha,\beta)\in\{(0.1,0.2),(0.6,0.6),(1,1)\}$, which are plotted as black crosses in Fig. LABEL:fig:deltaT. We emphasize again that in the region with $\Delta_{T}<0$ group sequences are preferred, even if when the sorting cost is considered, whereas in the region with $\Delta_{T}>0$ random sequences are preferred. Note that the cases with $\lambda=0$ correspond to those obtained without considering the sorting cost.

As discussed in the previous subsection, we note that $\Delta_{T}\leq 0$ for almost all $(\alpha,\beta)$ when $\lambda=0$, indicating that sorting is almost always beneficial for smooth transportation. However, once the sorting cost is considered, the sign of $\Delta_{T}$ can become positive, especially in the region where ${\rm min}(\alpha,\beta)$ is relatively small, indicating that sorting is not always beneficial. Note that the curves of $(\alpha,\beta)=$(0.6, 0.6) and (1, 1) are observed to overlap each other in Fig. LABEL:fig:lambda (b) and (c), unlike Fig. LABEL:fig:lambda (a). This happens because there is no difference in $\Delta_{T}$ at these two points when $\lambda=0$, as we can see in Figs. LABEL:fig:deltaT (b) and (c).

Figure LABEL:fig:withsort plots $\Delta_{T}$ for (a) $S=2$, (b) $S=3$, and (c) $S=N$ for various $(\alpha,\beta)$ with $\lambda=1$. In this figure, we note the existence of a new region in which $\Delta_{T}>0$, which we term a ‘random-advantageous region’ ($T_{\rm R}<T_{\rm G}$). This new region widens as $\lambda$ increases, finally resulting in the complete disappearance of the group-advantageous region for large enough $\lambda$. Note that Fig. LABEL:fig:withsort (c) exhibits only a random-advantageous region.

IV THEORETICAL ANALYSES

In this section, we show that the simulation results can be theoretically reproduced in some special cases. Specifically, we have succeeded in obtaining a mathematical proof of the appearance of the group-advantageous region for any group number $S(>1)$ when $\lambda=0$.

IV.1 Approximate flow of a multi-species TASEP

In this subsection, before calculating $T$, we briefly discuss the steady-state flow $Q_{S}$ of the multi-species TASEP that corresponds to a random sequence. We write $Q_{S}=Q_{S}(p_{1},......,p_{S};r_{1},......,r_{S})$, with the arguments abbreviated in obvious cases. When the flow $Q$ is simulated for each parameter set, we first evolve the system for $10^{5}$ time steps and then average over the next $10^{6}$ time steps.

IV.1.1 $L=2$

This subsection presents the derivation of an approximate $Q_{S}$ based on a Markov chain model. Due to the difficulty of considering general values of $L$ (the length of the lattice) and $S$ (the number of particle species), we consider the simplest case—with $L=2$ and $S=2$. As two species of particles exist—that is, particles with hopping probability $p_{1}$ and particles with $p_{2}$—each site may have three states: ‘unoccupied (state 0),’ ‘occupied by a particle 1 (state 1),’ and ’occupied by a particle 2 (state 2).’ This results in 9 possible states; however, noting that it is not necessary to distinguish the particle at site 1 because it always leaves the lattice with probability $\beta$, the number of possible states can be reduced to 6. Here, we define the probability distribution $P_{ij}$ ($i=0,1,2,j=0,\ast$), where $i$ and $j$ represent the state number of site 0 and 1, respectively. Note that state $\ast$ indicates either of state 1 or 2.

The master equations for the steady state are summarized in Eq. (2), using the relation $r_{1}(=r)+r_{2}=1$. Note that $r_{1}$ and $r_{2}$ are replaced with $r$ and $1-r$, respectively, for the case $S=2$. In addition, $P_{ij}$ must satisfy the normalization condition

[TABLE]

From Eqs. (2) and (3), the flow of the system can be written as a function of $p_{1}$, $p_{2}$, and $r$; that is, $Q_{2}(p_{1},p_{2};r)$, is given by the following expression:

[TABLE]

where

[TABLE]

and

[TABLE]

The specific forms of the probability distributions are summarized in Appendix C.

For $r=1$ and $p_{1}=p$, the system reduces to the single-species TASEP with the flow $Q_{1}(p)$, where

[TABLE]

Note that the flow of the single-species TASEP for general $L$ is exactly solved in Ref. Evans et al. (1999). Therefore, assuming that the value $p=p_{\rm h}$ satisfies the condition $Q_{2}(p_{1},p_{2};r)=Q_{1}(p)$, we can derive

[TABLE]

The quantity $p_{\rm h}$ is termed the harmonic mean of $p_{1}$ and $p_{2}$. Accordingly, for $L=2$, $Q_{2}$ is equivalent to $Q_{1}(p=p_{\rm h})$. This relation holds for any species number $S(>2)$, as we show in Appendix D.

Figure LABEL:fig:L2 compares the simulation and theoretical curves for various $\beta\in\{0.2,0.6,1\}$ with (a) $S=2$, (b) $S=3$ and (c) $S=N(=10,000)$. In all the figures, the simulations show very good agreement with our exact analyses.

IV.1.2 General $L(>2)$

For general $L$ and $S$, it is complicated to solve the master equations. Therefore, in this subsection, we instead assume an inequality, based on the results in the previous subsection and the qualitative discussions, and confirm the validity of the inequality by the simulations.

First, for general $L$ and $S$, $Q_{S}$ is clearly larger than $Q_{1}(p=p_{\rm min})$, where $p_{\rm min}={\rm min}\{p_{1},p_{2},......,p_{S}\}$.

In addition, for $L>2$, a platoon can be observed in the bulk of the lattice, in which a slower particle behaves as a bottleneck, and faster particles behind it cannot hop with a probability larger than that of the smaller one, i.e., less than $p_{\rm h}$, as shown in Fig. 12. This phenomenon suppresses the flow, implying that $Q_{S}$ is smaller than $Q_{1}(p=p_{\rm h})$.

Consequently, $Q_{S}$ satisfies the following inequality;

[TABLE]

In this subsection, we hereafter consider the case $S=2$.

Figure LABEL:fig:phasediagram shows the phase diagrams obtained by plotting the simulation values for (a) the single-species TASEP with $p=p_{\rm min}$, (b) the two-species TASEP, and (c) the single-species TASEP with $p=p_{\rm h}$, respectively. Note that $Q_{2}$ ($Q_{1}$) are the simulation (theoretical) values (and similarly hereafter).

Comparing these three figures shows that Eq. (9) obviously holds. In addition, as in Figs. LABEL:fig:phasediagram (a) and (c), we find that three different phases—HD, LD, and MC—also exist in Fig. LABEL:fig:phasediagram (b). Due to Eq. (9), the boundaries between the LD (HD) and MC phases of Fig. LABEL:fig:phasediagram (b) lie between those of Figs. LABEL:fig:phasediagram (a) and LABEL:fig:phasediagram (c). Note that the black lines in Figs. LABEL:fig:phasediagram (a) and (c) are theoretical boundaries, based on the fact that the boundary between the LD and MC phases of the single-species TASEP with hopping probability $p$ de Gier and Nienhuis (1999) can be written as

[TABLE]

Here, as for $\Delta_{T}$, we define $\Delta_{Q}$ as the ratio of the change from $Q_{1}(p=p_{\rm h})$ to $Q_{2}$; that is,

[TABLE]

and we note that $\Delta_{Q}=0$ when $Q_{1}(p=p_{\rm h})=0$.

Figure 14 shows $\Delta_{Q}$ for various ($\alpha,\beta$), for the fixed parameter set $(p_{1},p_{2};r)=(0.5,1;0.5)$. The black lines represent the boundaries between the LD/HD and MC phases of the single-species TASEP with hopping probability $p_{1}$ (boundary 1) and the single-species TASEP with $p=p_{\rm h}$ (boundary 2). Therefore, the lower-left (upper-right) region beyond boundary 1 (boundary 2) corresponds to the LD/HD (MC) phases both for the two-species TASEP and for the single-species TASEP with $p=p_{\rm h}$. This figure confirms that $\Delta_{Q}$ starts from 0 in the LD/HD phase, decreases, and finally yields to a constant value in the MC phase as $(\alpha,\beta)$ approaches the upper right.

Figure LABEL:fig:Q2 plots $\Delta_{Q}$ as a function ($p_{1},p_{2}$) for various $(\alpha,\beta)\in\{(0.1,0.2),(0.2,0.1),(1,1)\}$, fixing $r=0.5$. Note that both the single-species TASEP with hopping probability $p_{\rm h}$ and the two-species TASEP exhibit the LD, HD, and MC phases with $(\alpha,\beta)=$(0.1, 0.2), (0.2, 0.1), and (1, 1), respectively. This is because $\forall(p_{1},p_{2})$ $(0.2\leq p_{1},p_{2}\leq 1)$  $0.1<1-\sqrt{1-p_{\rm min}}\leq 1-\sqrt{1-p_{\rm h}}$. For example, Fig. LABEL:fig:phasediagram confirms that those three points exist within each corresponding phase for $(p_{1},p_{2})=(0.5,1)$.

In Figs. LABEL:fig:Q2 (a) and (b), we find that $\Delta_{Q}$ is approximately 0, whereas $\Delta_{Q}$ deviates from 0 in Fig. LABEL:fig:Q2 (c), as is also observed in Fig. 14. These phenomena can be explained as follows.

First, in the LD/HD phase, $Q_{2}$ is mainly governed by the input/output probability, leading to $\Delta_{Q}\to 0$, i.e., $Q_{2}$ approaches $Q_{1}(p=p_{\rm h})$. This is because $Q_{2}$ deviates from $Q_{1}(p=p_{\rm h})$ mainly due to the existence of platoons, which do not influence the flow much in this phase. Note that $\Delta_{Q}$ decreases as $\alpha$ or $\beta$ approaches 0.2, because the influence of platoons increases, approaching the MC phase of the two-species TASEP.

On the other hand, in the MC phase, $Q_{2}$ is mainly governed by the bulk region of the lattice. Therefore, the existence of platoons has a more critical influence on $Q_{2}$, causing $\Delta_{Q}$ to deviate from 0; i.e., $Q_{2}<Q_{1}(p=p_{\rm h})$. Especially as $|p_{1}-p_{2}|$ increases, the extent of the deviation also increases. This is because the effect of platoons increases when there is a large gap between $p_{1}$ and $p_{2}$.

Figure LABEL:fig:QLdepend plots (a) $Q_{2}/Q_{1}(p=p_{\rm h})$, (b) $Q_{3}/Q_{1}(p=p_{\rm h})$, and (c) $Q_{N}/Q_{1}(p=p_{\rm h})$ for various $L\in\{2,3,4,10,100,200\}$. Both of $Q_{S}$ ($S=2,3,10,000$) and $Q_{1}(p=p_{\rm h})$ are obtained by the simulations. In all the figures, we observe that for $L>2$, $Q_{S}/Q_{1}(p=p_{\rm h})$ generally becomes less than 1, i.e., $Q_{S}<Q_{1}(p=p_{\rm h})$, especially when $\alpha$ increases, i.e., the system approaches and exhibits the MC phase. Note that for $S=10,000$, the difference among $L=2,3,4$ is unclear; however, it becomes obvious when $L\geq 10$. Those results imply that Eq. (9) and its qualitative discussions can be applicable for general $L$ and $S$.

IV.2 Relation between $T_{\rm R}$ and $T_{\rm G}$

without the sorting cost

Hereafter, we assume $p_{1}<p_{2}<......<p_{S}$, $S>1$, $\forall s\ r_{s}>0$, and $\alpha>0$.

In this subsection, we fix $\lambda=0$, i.e., we do not consider the sorting cost. If for any number of particle species $s$, $r_{s}N$ is large enough for $T_{\rm G}$ and $T_{\rm R}$ to be determined by the steady-state flow (see Appendix E), we obtain

[TABLE]

and

[TABLE]

Note that this approximation immediately implies the independence of $T_{\rm G}$ from the order of the group sequence. Strictly speaking, $T_{\rm G}$ can differ depending on the order of each group in the group sequence. However, that difference can be ignored for large $N$ (see Appendix E).

In addition, we define the transportation times of the particles with the same hopping probabilities $p_{\rm h}$ and $p_{\rm min}$ as

[TABLE]

and

[TABLE]

respectively. From Eqs. (9), (14), (16) and (17), we immediately obtain the inequality

[TABLE]

In the following, we show that a general relation between $T_{\rm R}$ and $T_{\rm G}$ can be obtained mathematically for general $S(>1)$. We emphasize that this relation can be proven by comparing $T_{\rm H}$ and $T_{\rm G}$, and not by comparing $T_{\rm R}$ and $T_{\rm G}$ directly. Here, we introduce the new function $f(\alpha,\beta;p_{1},......,p_{S};r_{1},......,r_{S})$, which is defined as follows:

[TABLE]

Because we can assume $\alpha<\beta$ without loss of generality, we adopt this assumption in the following discussion, writing in abbreviated form $f(\alpha,\beta;p_{1},......,p_{S};r_{1},......,r_{S})=f(\alpha)$. Note that for cases with $\alpha\geq\beta$, the theoretical results can be obtained simply by replacing $\alpha$ (LD) with $\beta$ (HD).

A contour map of $f(\alpha)$ in the $(\alpha,\beta)$ plane exhibits four large regions, which are summarized in Tab. 1.

In the following subsections, we examine the behavior of $f(\alpha)$ according to this classification.

IV.2.1 Region 1: $\alpha<1-\sqrt{1-p_{1}}$

In this region, all the steady-state phases of the single-species TASEP for any $p_{s}$ exhibit the LD phase. Here, the steady-state flow for the single-species TASEP with parallel updating de Gier and Nienhuis (1999) is given by

[TABLE]

Therefore, we obtain

[TABLE]

and

[TABLE]

From Eqs.(23) and (24), we obtain $f(\alpha)$ as

[TABLE]

After some calculations, we obtain

[TABLE]

the detailed derivation is given in Appendix F.

IV.2.2 Region 2: $1-\sqrt{1-p_{1}}\leq\alpha<1-\sqrt{1-p_{\rm h}}$

This region is further divided into ($u-1$) subregions, as summarized in Tab. 2.

In Subregion 2–$v$ ($v=1,2,......,u-1$), the single-species TASEP with $p=p_{1},......,p_{v},p_{\rm h}$ exhibits the MC phase, whereas that with $p=p_{v+1},......,p_{S}$ displays the LD phase. Therefore, using Eq. (22), we obtain

[TABLE]

From Eqs. (15), (16), and (29), $T_{\rm G}$ and $T_{\rm H}$ can be written as follows;

[TABLE]

and

[TABLE]

From Eqs. (30) and (31), we thus obtain $f(\alpha)$ in the form

[TABLE]

For $1-\sqrt{1-p_{v}}\leq\alpha<1-\sqrt{1-p_{v+1}}$, due to $\alpha<p_{v+1}<......<p_{S}$ and $\alpha<p_{\rm h}$, we obtain $p_{s}-\alpha>0$ for $s=v+1,......,S$ and $p_{\rm h}-\alpha>0$. The function $f(\alpha)$ is continuous and differentiable with respect to $\alpha$, including at each boundary (see Appendix G). However, the signs of $f(\alpha)$ and $df(\alpha)/d\alpha$ are not specified. Note that the following condition

[TABLE]

where $q_{h}=1-\sqrt{1-p_{\rm h}}$ indicates that $f(\alpha)$ increases monotonically at least near the boundary between Subregion 2–($u-1$) and Region 3. This is discussed in Appendix H.

IV.2.3 Region 3: $1-\sqrt{1-p_{\rm h}}\leq\alpha<1-\sqrt{1-p_{S}}$

Similarly to Region 2, this region is further divided into ($S-u+1$) subregions, as summarized in Tab. 3. Note that Subregion 3–1 vanishes in the case $p_{\rm h}=p_{u}$, resulting in ($S-u$) subregions.

In Region 3–$v$ ($v=u,u+1,......,S$), the single-species TASEP with $p=p_{1},......,p_{v-1},p_{\rm h}$ exhibits the MC phase, whereas that with $p=p_{v},......,p_{S}$ displays the LD phase. Therefore, using Eq. (22), we obtain $Q_{1}(p)$

[TABLE]

From Eqs. (15), (16), and (36), $T_{\rm G}$ and $T_{\rm H}$ can be written as follows;

[TABLE]

and

[TABLE]

From Eqs. (37) and (38), $f(\alpha)$ becomes

[TABLE]

For $1-\sqrt{1-p_{v-1}}\leq\alpha<1-\sqrt{1-p_{v}}$, due to $\alpha<p_{v}<......<p_{S}$, we obtain $p_{s}-\alpha>0$ for $s=v,......,S$. Similarly to Region 2, $f(\alpha)$ is continuous and differentiable with respect to $\alpha$ including at each boundary.

After some calculations,we find that $df(\alpha)/d\alpha$ satisfies

[TABLE]

in each subregion, as discussed in detail in Appendix I. Eq. (40) indicates that $f(\alpha)$ is a monotonically increasing function of $\alpha$ throughout Region 3.

IV.2.4 Region 4: $1-\sqrt{1-p_{S}}\leq\alpha$

In this region, all the steady-state phases of the single-species TASEP for any hopping probabilities $p_{s}$ and $p_{\rm h}$ exhibit the MC region. Note that this region vanishes in the case $p_{S}=1$ because $1-\sqrt{1-p_{S}}=1$.

Therefore, we obtain

[TABLE]

and

[TABLE]

both of which are independent of $\alpha$. From Eqs. (41) and (42),we thus obtain $f(\alpha)$ as

[TABLE]

After some calculations, we obtain

[TABLE]

the detailed derivation of which is given in Appendix J.

IV.2.5 Relation between $T_{\rm R}$ and $T_{\rm G}$

With the results of Subsec. IV.2.1–IV.2.4, we can obtain a general relation between $T_{\rm R}$ and $T_{\rm G}$ for some special cases. Table 4 summarizes the signs of $f(\alpha)$ and $df(\alpha)/d\alpha$ in each region. Note that ‘U’ indicates that the sign is unclear.

Considering Tab. 4 and the continuity of $f(\alpha)$ including at each boundary (see Appendix G), we find from the intermediate value theorem that $\exists\alpha_{\rm cr}$ such that $f$ satisfies

[TABLE]

in Region 2 or 3. The specific conditions that $\alpha_{\rm cr}$ must satisfy are given in Appendix K.

Defining $\alpha_{\rm cr,max}$ as the largest value among the quantities $\alpha_{\rm cr}$, we obtain $f(\alpha)>0$—i.e., $T_{\rm H}>T_{\rm G}$—in the region where $\alpha>\alpha_{\rm cr,m}$. This is because $f(\alpha)$ is continuous and increases monotonically from a point in Region 2 (and through Region 3), to yield $f(\alpha)>0$ in Region 4.

Considering Eq. (18), we finally obtain

[TABLE]

in the region $\alpha>\alpha_{\rm cr,m}$. Eq. (46) means that $\Delta_{T}<0$, reproducing the simulation results in the region where ${\rm min}(\alpha,\beta)$ is relatively large. This result indicates that the group-advantageous region must appear even in a case with $p_{S}=1$, for which Region 4 vanishes.

In analogy with the discussion above, we can also predict that a region with $\Delta_{T}<0$ must appear in the case $S=N$.

IV.3 Relation between $T_{\rm R}$ and $T_{\rm G}$ with sorting cost

In this subsection, we discuss the change in the relation between $T_{\rm R}$ and $T_{\rm G}$ when $\lambda>0$, i.e., when the sorting cost is included. In the following, we first obtain a general formula for the sorting cost and then evaluate upper and lower limits to $\lambda$.

IV.3.1 General formula for the sorting cost

First, we calculate mathematically the averaged minimal number of exchanges necessary to sorting the particles from random to group sequences.

We here define $\overline{K}$ as the averaged value of $K$, using the fact that $\tau_{\rm R}$ can take $N!/\prod_{s=1}^{S}(r_{s}N)!$ patterns with equal probability. We thus have

[TABLE]

If $K^{\prime}(\tau_{\rm G},\tau_{\rm R})$ is the minimal number of exchanges necessary to sort the particles from a random sequence $\tau_{\rm R}$ to a given fixed group sequence $\tau_{\rm G}$, then $K^{\prime}(\tau_{\rm G},\tau_{\rm R})$ satisfies

[TABLE]

Note that the number of elements of $\{K(\tau_{\rm G},\tau_{\rm R})\}$ is equal to that of $\{\tau_{\rm G}\}$ from the definition. Eqs. (47) and (48) indicate that the best group sequence $\tau_{\rm G}$ can vary depending on the particular random sequence $\tau_{\rm R}$.

Due to the difficulty of a general calculation of $\overline{K}$, we instead calculate $\overline{K}^{\prime}$, which is defined as follows:

[TABLE]

where $\tau_{\rm G}$ is a fixed sequence out of the set {$\tau_{\rm G}$} for all possible $\tau_{\rm R}$.

For $S=2$ and $S=N$, $\overline{K}^{\prime}$ can be generally calculated as

[TABLE]

the detailed derivations of which are discussed in Appendix L.

Figure 17 shows the ratio $\overline{K}/\overline{K}^{\prime}$ for (a) $S=2$ and (b) $S=N$. Both figures show that $\overline{K}/\overline{K}^{\prime}\approx 1$, i.e., $\overline{K}\approx\overline{K}^{\prime}$, indicating that there is no problem in substituting $\overline{K}^{\prime}$ for $\overline{K}$ for large enough $N$.

In the following calculations, we therefore use $\overline{K}^{\prime}$ instead of $\overline{K}$ because $\overline{K}^{\prime}$ can be represented by a general formula, whereas $\overline{K}$ cannot.

IV.3.2 Upper and lower limits to $\lambda_{\rm cr}$

We first define $\lambda=\lambda_{\rm cr}\geq 0$, as the value for which $T_{\rm R}=T_{\rm G}$. Note that $\lambda_{\rm cr}$ is defined to be equal to 0 if $T_{\rm R}\leq T_{\rm G}$ when $\lambda=0$. From the definition of $\lambda_{\rm cr}$, the random-advantageous region appears when $\lambda>\lambda_{\rm cr}$. Based on the discussions in Subsec. IV.2 and IV.3.1, we here evaluate $\lambda_{\rm cr}$ for $S=2$ and $\alpha<\beta$.

To take into account the sorting cost, we add the term $\lambda K\ (\lambda>0)$ to $T_{\rm G}$; that is,

[TABLE]

Conversely, we do not add that term to $T_{\rm R}$ because a random sequence means a sequence without sorting.

We also define $\lambda_{\rm H}$ and $\lambda_{\rm M}$ as the values of $\lambda$ for which $T_{\rm H}=T_{\rm G}$ and $T_{\rm M}=T_{\rm G}$, respectively. Note that $\lambda_{\rm H}$ can have negative values, because $T_{\rm H}$ can be less than $T_{\rm G}$.

From Eq. (18) and Subsec. IV.2, when $\lambda=0$ the relations among $T_{\rm R}$, $T_{\rm G}$, $T_{\rm H}$, and $T_{\rm M}$ must satisfy one of the following three inequalities:

[TABLE]

[TABLE]

or

[TABLE]

Therefore, the relations of $\lambda_{\rm cr}$, $\lambda_{\rm H}$, and $\lambda_{\rm M}$ can be written as follows:

[TABLE]

where we note that by definition $\lambda_{\rm cr}\geq 0$, whereas $\lambda_{\rm H}$ can be either negative or positive, while $\lambda_{\rm M}$ must be positive.

In Region 1, i.e., $\alpha<1-\sqrt{1-p_{1}}$, where $T_{\rm H}<T_{\rm G}$ with $\lambda=0$, $\lambda_{\rm H}$ must be negative, while $\lambda_{\rm M}$ must be positive, and satisfy

[TABLE]

In Region 2, i.e., $1-\sqrt{1-p_{1}}\leq\alpha<1-\sqrt{1-p_{\rm h}}$, $\lambda_{\rm H}$ can be either negative or positive, whereas $\lambda_{\rm M}$ must be positive. The quantities $\lambda_{\rm H}$ and $\lambda_{\rm M}$ satisfy

[TABLE]

and

[TABLE]

respectively.

In Region 3, i.e., $1-\sqrt{1-p_{\rm h}}\leq\alpha<1-\sqrt{1-p_{2}}$, $\lambda_{\rm H}$ can be either negative or positive, and $\lambda_{\rm M}$ must be positive. Thus, $\lambda_{\rm H}$ and $\lambda_{\rm M}$ satisfy

[TABLE]

and

[TABLE]

respectively.

In Region 4, i.e., $1-\sqrt{1-p_{2}}\leq\alpha$, due to $T_{\rm G}<T_{\rm H}<T_{\rm M}$ when $\lambda=0$, $\lambda_{\rm H}$ and $\lambda_{\rm M}$ must both be positive. The quantities $\lambda_{\rm H}$ and $\lambda_{\rm M}$ therefore satisfy

[TABLE]

and

[TABLE]

respectively.

Table VI summarizes the upper and lower limits to $\lambda_{\rm cr}$ in each region. Furthermore, Fig. 17 shows the simulation values (black circles) and the theoretical existence range of $\lambda_{\rm cr}$ (yellow region) as functions of $\alpha$. Note that the we calculated the simulation values using with 10-trial-averaged values of $T_{\rm R}$, $T_{\rm G}$, and $K$.

We can interpret Fig. 17 as demonstrating that a group (random) sequence is preferable in the region below (above) the black line. In the blue (green) region, a group (random) sequence is in fact theoretically verified to be preferable. Comparing Figs. 17 (a) and (b), the simulation values approach the lower limit—i.e., the accuracy of approximating $T_{\rm R}$ by $T_{\rm H}$ increases—as $|p_{1}-p_{2}|$ decreases. We admit that the yellow region is extensive, especially when $|p_{1}-p_{2}|$ is relatively large; however, we emphasize that the simulation values always exist within the expected region and that the region can be limited easily without numeric calculations, which is convenient for applications to actual situations.

V CONCLUSION

In the present study, we have used a modified TASEP to analyze the dependence of the transportation time on the entering sequences of particles, using both the numerical simulations and theoretical analyses.

Here, we summarize a number of important results. In Sec. III, we discovered that there exists an important ‘group-advantageous region’ where $T_{\rm R}>T_{\rm G}$ when ${\rm min}(\alpha,\beta)$ is relatively large and the sorting costs are neglected. When sorting costs are introduced, a new region called a ‘random-advantageous region’ appears with $T_{\rm R}<T_{\rm G}$. In addition, the group-advantageous region shrinks and finally disappears as $\lambda$ increases. We explored these phenomena for various $S\in\{2,3,N\}$.

In Sec. IV, we analyzed the simulation results by employing mathematical approaches for certain special cases. Using some approximations, we have shown theoretically that without the sorting cost the group-advantageous region must appear for any parameter sets ($S,p_{s},r_{s}$). Moreover, we have succeeded in deriving the upper and lower limits to the value of $\lambda_{\rm cr}$ where $T_{\rm R}=T_{\rm G}$ by obtaining a general formula for the sorting cost.

Our findings can be applied to real-world situations, such as providing efficient operation for various tasks and smooth logistics for various products and yielding an effective evacuation method for pedestrians. Specifically, for smooth operation, we can determine whether we should begin tasks without considering the operation sequence or otherwise. Similarly, for smooth logistics, we can select whether the products should be bunched with nearly equal sizes. Furthermore, for ensuring effective evacuation of pedestrians, we can determine whether the bunching of pedestrians having nearly equal velocities should be conducted before transportation. The criteria for these judgments depend on the magnitude of the consideration or bunching cost ($\lambda$). Note that these magnitudes significantly differ from each other, i.e., considering only the sequence of tasks is typically deemed cheaper (have a smaller $\lambda$) than sorting various pedestrians and products.

ACKNOWLEDGMENTS

This work was partially supported by JST-Mirai Program Grant Number JPMJMI17D4, Japan, JSPS KAKENHI Grant Number JP15K17583, and MEXT as ’Post-K Computer Exploratory Challenges’ (Exploratory Challenge 2: Construction of Models for Interaction Among Multiple Socioeconomic Phenomena, Model Development and its Applications for Enabling Robust and Optimized Social Transportation Systems) (Project ID: hp180188).

Appendix A Validity of our selection of $L$ and $N$

In this Appendix, we briefly discuss the validity of selecting $L=200$ and $N=10,000$.

As finite-size effects may occur for small $L$, we compare the simulation values of $Q_{2}$ for $L=200$ and $L=1,000$. Figure 18 shows the ratio $Q_{2}/Q_{2}^{\prime}$, where $Q_{2}$ and $Q_{2}^{\prime}$ represent the flow of the multi-species TASEP with $L=200$ and $L=1,000$, respectively, as functions of $\alpha$ for various $\beta\in\{0.2,0.6,1\}$. The result that $Q_{2}/Q_{2}^{\prime}\approx 1$ indicates that the effect can be ignored for $L=200$. Thus, we choose $L=200$ to decrease the simulation time.

On the other hand, the assumption that $T$ is determined by a steady-state flow may be inappropriate for small $N$. Therefore, we have compared the results for $N=10,000$ and $N=20,000$, in both cases for $L=200$. Figure 19 shows the ratio $T/T^{\prime}$, where $T$ and $T^{\prime}$ represent the transportation times for $N=10,000$ and $N=20,000$, respectively, as functions of $\alpha$ for various $\beta\in\{0.2,0.6,1\}$. The result that $T/T^{\prime}\approx 0.5$, i.e., that $T$ is proportional to $N$, indicates that the assumption can be regarded as valid for $N=10,000$. Thus, we choose $N=10,000$ similarly to decrease the simulation time.

Appendix B Simulation schemes for obtaining the minimal number of necessary exchanges

In this Appendix, we briefly describe the specific simulation schemes we used to obtain $K$. We emphasize that the cost of counting or comparing particles and the distances between exchanged particles are both ignored in the following.

First, for $S=2$, $\tau_{\rm a}=\tau_{\rm G}$ can have only one of two patterns. Once $\tau_{\rm G}$ is fixed to be either of these two sequences, we can immediately obtain the number of particles placed at the wrong areas in sequence $\tau_{\rm b}=\tau_{\rm R}$, which is twice as large as the number of necessary exchanges (see also Appendix L). Consequently, comparing the results for the two $\tau_{\rm G}$ gives the smaller number as $K$.

Second, for $S=3$, $\tau_{\rm a}=\tau_{\rm G}$ can have six patterns. Once $\tau_{\rm G}$ is fixed at one of these six sequences, we can immediately obtain the number of particles placed at the wrong areas in any sequence $\tau_{\rm b}=\tau_{\rm R}$. After selecting one species, which we first replace at the correct location, we exchange all particles of that species that are placed in the wrong areas in sequence $\tau_{\rm b}=\tau_{\rm R}$. The subsequent procedure is similar to the case for $S=2$. Consequently, comparing the six results for each $\tau_{\rm G}$ again gives the smallest number as $K$. Note that we can similarly calculate the numbers for general $S>4$.

Finally, for $S=N$, $\tau_{\rm a}=\tau_{\rm G}$ can have one of two patterns: either an ascending or a descending sequence. One exchange is needed for each particle in $\tau_{\rm b}=\tau_{\rm R}$ for which there exists a particle with a smaller (larger) hopping probability than the noted particle. This is termed a ‘selection sort.’ This procedure starts from the leading particle. Consequently, by comparing the results for the two $\tau_{\rm G}$, the smaller number is again selected as the minimal number of necessary exchanges.

Appendix C Probability distribution with $L=2$ and $S=2$

Here, we summarize the probability distributions with $L=2$ and $S=2$, which can be obtained from Eqs. (2) and (3). The specific forms are described as follows:

[TABLE]

where

[TABLE]

[TABLE]

and

[TABLE]

Appendix D $Q_{S}$ for general $S$ with $L=2$

In this appendix, we prove that for general $S$ with $L=2$, $Q_{S}$ is equal to $Q_{1}(p=p_{\rm h})$.

From the results with $L=2$ and $S=2$ (see Appendix C), we can conjecture the probability distributions for general $S$ with $L=2$ as

[TABLE]

where $s=1,2,......,S$.

On the other hand, the master equations of the steady state are summarized as $2(S+1)$ equations:

[TABLE]

where $s=1,2,......,S$. In addition, $P_{ij}$ must satisfy the normalization condition

[TABLE]

We can confirm that Eqs. (69) satisfy Eqs. (70) and (71). With Penron-Frobenius theorem regarding stochastic matrix, this indicates that Eqs. (69) are unique solutions for Eqs. (70) and (71).

From Eqs. (69), the flow of the system is given by the following expression:

[TABLE]

Appendix E Validity of the approximation for $T$

In this Appendix, we briefly demonstrate the validity of Eq. (15).

Figure 20 (a) shows the ratio $T_{\rm G,sim}/T_{\rm G,theo}$ as a function of $\alpha$ for various $\beta\in\{0.2,0.6,1\}$ with $S=2$ and Fig. 20 (b) shows the same ratio for $S=3$. Note that $T_{\rm G,sim}$ and $T_{\rm G,theo}$ represent the values of $T_{\rm G}$ from the simulations and that given by Eq. (15), respectively. Both figures show that $T_{\rm G,sim}/T_{\rm G,theo}\approx 1$, indicating that Eq. (15) provides a good approximation for $T_{\rm G}$.

Strictly speaking, $T_{\rm G,sim}/T_{\rm G,theo}$ must be larger than 1 on average. This is mainly due to the fact that $T_{\rm G,sim}$ includes $T_{1}$, which is the time required for the first particle to reach the right-hand boundary, whereas $T_{\rm G,theo}$ ignores that time. This also indicates that $T_{\rm G}$ can differ depending on the order of each group in the group sequence (i.e., the hopping probability of the leading group). However, this difference has little influences on the theoretical results, as explained below.

First, $T_{1}$ can be estimated as

[TABLE]

where $s=1,2,......,S$ and the time steps before the first particle enters the lattice are assumed to be small enough to be ignored. The quantities $T_{1}$ and $T_{G}$ without $T_{1}$ satisfy

[TABLE]

and

[TABLE]

respectively. Therefore, $T_{1}/T_{\rm G}$ reduces to

[TABLE]

Under the proposition that $N$ is large enough, we can assume $L/2N<<1$ ($L/2N=0.01$ in the present study). In fact, observing the time series of the flows (199-steps central moving average) in Fig. 21, we find that nearly the entire duration during transportation can be regarded to be in the steady state for large enough $r_{s}N$. Note that we calculate the flows at time $t$ by averaging number of moving particles per bond between $t-1$ and $t$. Moreover, all the transportation times $T$ ($T_{\rm R}$, $T_{\rm H}$, and $T_{\rm M}$) originally include $T_{1}$, so that this term disappears when they are subtracted from each other. Consequently, $T_{1}$ (and therefore, the dependence of $T_{\rm G}$ on the order of each group in the group sequence) can be assumed to be ignorable.

Appendix F Discussion of the sign of $f(\alpha)$

in Region 1

In this Appendix, we give a detailed derivation of Eq. (26) for Region 1, where $\alpha<1-\sqrt{1-p_{1}}$.

Eq. (25) gives $f(\alpha)$

[TABLE]

where

[TABLE]

The quantity $C$ is calculated as follows:

[TABLE]

where

[TABLE]

By regarding the sum of the term with $(s,t)=(x,y)$ and that with $(s,t)=(y,x)$ as a new term for $\exists(x,y)\ (x,y=1,2,......,S,x<y)$, we can rewrite Eq. (79) as follows:

[TABLE]

Because $p_{s}-\alpha>0$ ($\forall s$) and $p_{\rm h}-\alpha>0$, we obtain $C>0$.

Considering $\alpha-1<0$ and $C>0$, we finally obtain

[TABLE]

Appendix G Continuity and differentiability of $f(\alpha)$ at each boundary

In this Appendix, we briefly discuss the continuity and differentiability of $f(\alpha)$ at each boundary.

Defining $g(x)$ for $0<x\leq 1$ as

[TABLE]

where $0<p\leq 1$, the following equations hold:

[TABLE]

and

[TABLE]

where $q=1-\sqrt{1-p}$. Therefore, $g(x)$ is continuous and differentiable at $x=q=1-\sqrt{1-p}$, resulting in the continuity and differentiability of $g(x)$ for $0<x\leq 1$.

As a result, because $f(\alpha)$ is represented as a linear sum of terms $g(\alpha)$, where $p$ is substituted for $p_{s}$ or $p_{\rm h}$ ($0<p_{s},p_{\rm h}\leq 0$), $f(\alpha)$ is clearly continuous and differentiable at each boundary.

Appendix H Discussion of the sign of $df(\alpha)/d\alpha$

in Subregion 2–$v$

In this Appendix, we discuss the sign of $df(\alpha)/d\alpha$ in Subregion 2–$v$, i.e., $1-\sqrt{1-p_{v}}\leq\alpha<1-\sqrt{1-p_{v+1}}$

From Eq. (32), $df(\alpha)/d\alpha$ can be calculated as follows:

[TABLE]

For $1-\sqrt{1-p_{v}}\leq\alpha<1-\sqrt{1-p_{v+1}}$, the following two inequalities hold:

[TABLE]

and

[TABLE]

We cannot specify the sign of $df(\alpha)/d\alpha$ in this subregion from Eqs. (88), (89), and (90). However, near the boundary between Subregion 2–$(u-1)$ and 3–$u$, we obtain the following conditions:

[TABLE]

and

[TABLE]

where $q_{h}=1-\sqrt{1-p_{\rm h}}$. Therefore, noting the obvious continuity of $df(\alpha)/d\alpha$ for $1-\sqrt{1-p_{v}}\leq\alpha<1-\sqrt{1-p_{v+1}}$, the region of $df(\alpha)/d\alpha>0$ must exist at least in Subregion 2–$(u-1)$.

Appendix I Discussion of the sign of $df(\alpha)/d\alpha$

in Subregion 3–$v$

In this Appendix, we give a proof on Eq. (40) in Subregion 3–$v$, i.e., $1-\sqrt{1-p_{v-1}}\leq\alpha<1-\sqrt{1-p_{v}}$.

From Eq. (39), $df(\alpha)/d\alpha$ can be calculated as follows:

[TABLE]

For $s=v,......,S$, the quantity $p_{s}-2\alpha+\alpha^{2}$ satisfies

[TABLE]

From Eqs. (93) and (94), we finally obtain

[TABLE]

Appendix J Discussion of the sign of $f(\alpha)$

in Region 4

In this Appendix, we give a detailed derivation of Eq. (44), where $1-\sqrt{1-p_{S}}\leq\alpha$.

From Eqs. (41) and (42), $f(\alpha)$ can be represented as follows:

[TABLE]

and

[TABLE]

From Eqs. (43), (96) and (97), $f(\alpha)$ is given by

[TABLE]

where

[TABLE]

and

[TABLE]

From Eqs. (99) and (100), $E^{2}-F^{2}$ becomes

[TABLE]

Here, regarding the sum of the term with $(s,t)=(x,y)$ and that with $(s,t)=(y,x)$ as a new term for $\exists(x,y)\ (x,y=1,2,......,S,x<y)$, Eq. (101) can be rewritten as

[TABLE]

Due to Eq. (102) and the non-negativity of both $E$ and $F$, we have $E>F$, thereby resulting in

[TABLE]

Appendix K Specific conditions on $\alpha_{\rm cr,max}$

Here, we discuss the specific conditions on $\alpha_{\rm cr,max}$.

Table VI summarizes explicit expressions for $f(\alpha=\alpha_{\rm cr})=0$, where the upper (lower) expression holds in Region 2 (Region 3). For $S=2$, the lower expression becomes a quadratic equation in $\alpha_{\rm cr}$. However, the upper expression becomes a quartic equation that is too difficult to solve analytically those conditions. Note that for $S>2$, both equations become more than quartic.

From its definition of $\alpha_{\rm cr,max}$, $\alpha_{\rm cr,max}$ can be written as

[TABLE]

where $\{\alpha_{\rm cr}\}$ represents the set of $\alpha_{\rm cr}$.

Appendix L Derivation of $\overline{K}^{\prime}$

In this Appendix, we derive the approximate averaged minimal number of exchanges $\overline{K}^{\prime}$ necessary to sort the particles for two special cases: $S=2$ and $S=N$.

L.1 $S=2$

First, for a general calculation of $\overline{K}^{\prime}$, $\tau_{\rm G}$ has to be fixed to be either of the two possible patterns. Once $\tau_{\rm G}$ is fixed, $K(\tau_{\rm G},\tau_{\rm R})$ can be determined uniquely for all possible $\tau_{\rm R}$. Without loss of generality, we can assume $rN\leq(1-r)N$ and $\tau_{\rm G}$ can be fixed as illustrated in the lower panel of Fig. 22.

Suppose that for $\tau_{\rm R}$, $k\ (0\leq k\leq rN)$ particles of species 1 are located in the Area 2, ($k$ particles of species 2 are located in the Area 1, conversely) as described in the lower of Fig. 22. Under this supposition, $k$-time exchanges are necessary for sorting particles from $\tau_{\rm R}$ to $\tau_{\rm G}$. Considering that $\tau_{\rm R}$ satisfying this supposition possibly has $\binom{rN}{k}\times\binom{(1-r)N}{k}$ sequences, $a_{N}=\sum_{\forall{\tau_{\rm R}}}K^{\prime}(\tau_{\rm G},\tau_{\rm R})$ can be written as follows;

[TABLE]

Using the Vandermonde convolution formula, Eq. (105) can be rewritten as follows:

[TABLE]

Because the sequence $\tau_{\rm R}$ can take any of $N!/\{(rN)!((1-r)N)!\}$ possible patterns with equal probability, we can finally reduce $\overline{K}^{\prime}$ to

[TABLE]

Figure 23 compares the simulation (circles) and theoretical (curves) values for various $N\in\{1,000({\rm red}),5,000({\rm green}),10,000({\rm blue})\}$ for $S=2$. The simulations show very good agreement with our exact analysis.

L.2 $S=N$

When $S=N$, $\tau_{\rm G}$ also has to be fixed as either of the two possible patterns—an ascending or a descending sequence—for a general calculation of $\overline{K}^{\prime}$, as illustrated in the upper panel of Fig. 24. Once $\tau_{\rm G}$ is fixed, $K^{\prime}(\tau_{\rm G},\tau_{\rm R})$ can be determined uniquely for all possible $\tau_{\rm R}$.

If we regard the entire sequence as consisting of two parts—the first (blue) particle and other ($N-1$) particles, as described in the lower panel of Fig. 24—the sorting procedure can also be divided into two parts: sorting ($N-1$) particles plus the last exchange for the first particle. If the first particle corresponds to the particle with hopping probability $p_{l}$ $(l=1,2,......,N)$, and noting that the sequence for the remaining ($N-1$) particles has $(N-1)!$ possible patterns, we can calculate the quantity $b_{N,l}=\sum_{\forall{\tau_{{\rm R},l}^{\prime}}}K(\tau_{\rm G},\tau_{{\rm R},l}^{\prime})$ as follows:

[TABLE]

where $N>1$ and $\tau_{{\rm R},l}^{\prime}$ represents the sequence for which the first particle is the particle with hopping probability $p_{l}$. Note that the last sort is not necessary in the case where $l=1$.

Therefore, for $N>1$, we can write $a_{N}=\sum_{\forall{\tau_{\rm R}}}K(\tau_{\rm G},\tau_{\rm R})$:

[TABLE]

Dividing both sides of Eq. (109) by $N!$, we obtain

[TABLE]

where $c_{N}=a_{N}/N!$ and $N>2$. With the initial condition $c_{1}=a_{1}=0$, $a_{N}$ is finally reduced to

[TABLE]

which we note holds for the case $N=1$.

The sequence $\tau_{\rm R}$ can take $N!$ patterns with equal probability, and therefore, $\overline{K}^{\prime}$ is finally reduced to

[TABLE]

Figure 25 compares the simulation (circles) and theoretical (line) values for $S=N$. The simulations again show a very good agreement with our exact analysis.