Recurrence dynamics of particulate transport with reversible blockage: from a single channel to a bundle of coupled channels

TL;DR

This paper models particulate flow through channels with reversible blockage, analyzing how system parameters affect throughput and open probability, and compares coupled versus uncoupled channel efficiencies.

Contribution

It provides analytic solutions for small capacities and numerical analysis for larger systems, revealing complex behaviors and efficiency differences in coupled versus uncoupled channels.

Findings

Exiting flux can increase monotonically or have a maximum depending on deblocking time.

Large N systems show abrupt transition from few blockages to permanent blockage.

Coupled channels are more efficient for N=1, with complex behavior for N>1.

Abstract

We model a particulate flow of constant velocity through confined geometries, ranging from a single channel to a bundle of $N_c$ identical coupled channels, under conditions of reversible blockage. Quantities of interest include the exiting particle flux (or throughput) and the probability that the bundle is open. For a constant entering flux, the bundle evolves through a transient regime to a steady state. We present analytic solutions for the stationary properties of a single channel with capacity $N\le 3$ and for a bundle of channels each of capacity $N = 1$. For larger values of $N$ and $N_c$, the system's steady state behavior is explored by numerical simulation. Depending on the deblocking time, the exiting flux either increases monotonically with intensity or displays a maximum at a finite intensity. For large $N$ we observe an abrupt change from a state with few blockages to one in which the bundle is permanently blocked and the exiting flux is due entirely to the release of blocked particles. We also compare the relative efficiency of coupled and uncoupled bundles. For $N=1$ the coupled system is always more efficient, but for $N>1$ the behavior is more complex.

Figures

Equations

$$p_o(\lambda )=\frac{⟨t⟩}{⟨t⟩+{\tau }_b},$$

$$j(\lambda )=\frac{⟨m⟩+N}{⟨t⟩+{\tau }_b}.$$

$$⟨m⟩+N=\lambda ⟨t⟩.$$

$$j(\lambda )=\lambda p_o(\lambda ).$$

$$p_o(\lambda )\simeq \frac{1}{1+\frac{(\lambda \tau {)}^N{\tau }_b}{(N-1!)\tau }}.$$

$$p_o(\lambda )\simeq \frac{N}{N+\lambda {\tau }_b},$$

$$j(\lambda )\simeq \frac{N\lambda }{N+\lambda {\tau }_b}.$$

$$\frac{d{p}_o(t)}{dt}=-\lambda p_o(t)+\lambda p_o(t-{\tau }_b).$$

$$j(t)=\lambda p_o(t-{\tau }_b),$$

$${\tilde{p}}_o(u)=\frac{1}{\lambda +u-\lambda e^{-u{\tau }_b}},$$

$$\tilde{j}(u)=\frac{\lambda e^{-u{\tau }_b}}{\lambda +u-\lambda e^{-u{\tau }_b}}.$$

$$p_o(t)=\sum _{n=0}^{\infty }\frac{[\lambda (t-n{\tau }_b){]}^n}{n!}{e}^{-\lambda (t-n{\tau }_b)}\theta (t-n{\tau }_b),$$

$$p_o(\lambda )=\frac{1}{1+\lambda {\tau }_b},$$

$$j(\lambda )=\frac{\lambda }{1+\lambda {\tau }_b}.$$

$$\frac{d{p}_0}{dt}$$

$$\frac{d{p}_1}{dt}$$

$$\frac{d{p}_2}{dt}$$

$$p_0(0)=1,p_1(0)=p_2(0)=0.$$

$$\left[u+\lambda (1-e^{-\tau (u+\lambda )})\right]{\tilde{p}}_0(u)-\lambda e^{-u{\tau }_b}{\tilde{p}}_1(u)$$

$$-\lambda (1-e^{-\tau (u+\lambda )}){\tilde{p}}_0(u)+(\lambda +u){\tilde{p}}_1(u)$$

$$-\lambda (1-e^{-u{\tau }_b}){\tilde{p}}_1(u)+u{\tilde{p}}_2(u)$$

$${\tilde{p}}_0(u)$$

$${\tilde{p}}_1(u)$$

$${\tilde{p}}_2(u)$$

$$\Delta =(u+\lambda {)}^2-\lambda (\lambda +u)e^{-\tau (u+\lambda )}-{\lambda }^2{e}^{-u{\tau }_b}(1-e^{-\tau (u+\lambda )}).$$

$$j(t)=\lambda e^{-\lambda \tau }{p}_0(t-\tau )+2\lambda p_1(t-{\tau }_b).$$

$$\tilde{j}(u)=\frac{\lambda }{\Delta }[(u+\lambda )e^{-(u+\lambda )\tau }+2\lambda e^{-u{\tau }_b}(1-e^{-\tau (u+\lambda )})].$$

$$p_o(\lambda )=\frac{2-e^{-\lambda \tau }}{2+\lambda {\tau }_b-(1+\lambda {\tau }_b)e^{-\lambda \tau }},$$

$$j(\lambda )=\lambda \frac{2-e^{-\lambda \tau }}{2+\lambda {\tau }_b-(1+\lambda {\tau }_b)e^{-\lambda \tau }}.$$

$$\lambda ⟨t⟩=\frac{2e^{\nu }\sinh (g\nu )+g{e}^{\lambda \tau }}{-g-2\sinh (g\nu )e^{-\nu }+e^{\nu }\left(\sinh (g\nu )+g\cosh (g\nu )\right)}+1$$

$$\lambda ⟨t⟩=\frac{2e^{\nu }\sin (g\nu )+g{e}^{\lambda \tau }}{-g-2\sin (g\nu )e^{-\nu }+e^{\nu }\left(\sin (g\nu )+g\cos (g\nu )\right)}+1$$

$$j(\lambda )=\frac{\lambda ⟨t⟩}{⟨t⟩+1/{\mu }^{\ast }}.$$

$$⟨t⟩=\frac{1}{\lambda }\sum _{j=0}^{N-1}\frac{N!}{(j+1)(N-j-1)!}{\left(\frac{\mu }{\lambda }\right)}^j.$$

Full text

I Introduction

Many physical and biological processes feature particle flow in confined geometries. Examples include vehicular and pedestrian traffic Helbing (2003); Appert-Rolland et al. (2010), filtration of particle suspensions, and the flow of macromolecules or ions through micro- or nanochannels Finkelstein and Andersen (1981); Kelkar and Chattopadhyay (2007); Daiguji et al. (2004). In biology, an important example is flux regulation mediated by pore-forming membrane proteins. The transport of ions and water through these channels is primarily a single-file process, i.e. cations and water molecules cannot pass each other within the channel Daiguji et al. (2004). Moreover, these flux regulatory channels can be clogged by toxins or medicines, with significant consequences. Recent studies of tracer diffusion of biased (or active) particles in a crowded, narrow channel revealed a non-trivial relationship between the exerted force and the mean velocity, as well asymmetric density profiles of the environment Bénichou et al. (2018a); Bertrand et al. (2018); Krapivsky et al. (2014); Bénichou et al. (2018b).

Some information systems, such as telecommunication Senderovich et al. (2015) and computing networks Tadayon and Aissa (2014); Ezaki et al. (2015), as well as trunked mobile radio systems and air traffic Janssen and Van Leeuwaarden (2008); Barcelo et al. (1996); Masrub (2014) are also amenable to the channel description.

A blockage may be caused by either ‘extrinsic’ or ‘intrinsic’ mechanisms. The former refers to the situation where the number of particles present somehow exceeds the channel carrying capacity, and will be the focus of this article. The latter mechanism arises from collective effects such as encountered in filtration processes. In this case, while isolated particles can pass through a mesh hole, clogging occurs when two or more particles arrive in near concurrence, causing one to impede the other. This effect, due to the delicate interplay between the spatio-temporal closeness of the particles and the confining geometry, could be seen as setting the capacity of a channel to greater than one, which establishes a connection between both types of blocking mechanisms. A model based on this phenomenology successfully accounted for experimental data Roussel et al. (2007); Redner and Datta (2000). Various approaches, including the Totally Asymmetric Simple Exclusion Process (TASEP)Mallick (2011, 2015) have been applied to model these phenomena. Gabrielli et al. introduced a class of stochastic models in which blockages occur when the carrying capacity of a channel is exceeded Gabrielli et al. (2013, 2015). For these stochastic models, the particle velocity within the channels is identical, and the mean particle density is low enough to prevent exclusion effects. The blockage is triggered when the number of particles within the channel at a given time exceeds the channel capacity. The original model considered one channel with capacity $N=2$, i.e., two particles must be simultaneously present in the channel to block the system. Particles enter at random times according to a Poisson process of intensity $\lambda$ and exit, if no blockage occurs, after a fixed transit time $\tau$. Subsequently, several generalizations were studied, including a higher blocking threshold ($N>2$) Barré et al. (2015), an inhomogeneous entering flux Barré and Talbot (2015a), and multiple channels Barré and Talbot (2015b, 2016). When the blockage is reversible, the system is reactivated after a constant waiting time, $\tau_{b}$. This mechanism gives rise to a transient regime leading to a steady state Barré et al. (2013). In this article we associate the last two generalizations.

Queuing theory Adan and Resing (2002); Medhi (1991) provides an alternative description of reversible blocking phenomena. This framework is typically used to describe customers arriving at a server, according to one given distribution, and receiving service according to another distribution. If all the elementary steps are Markov processes, the time evolution of the state probabilities can be described by systems of linear differential equations. In accordance with this approach, we recently introduced Markovian models of blockage Page et al. (2018); Barré et al. (2018) for which exact solutions can be obtained.

This paper is organized as follows: We first consider the single channel model in Sec. II, for which exact solutions are available for $N=1,2,3$, and simulation results for greater values of $N$. We also compare our results with the previously introduced Markovian models. In Sec. III a channel bundle model consisting of $N_{c}$ identical channels each of capacity $N$ and sharing an incoming flux is presented. Exact results are obtained for $N=1$, as well as numerical simulation results for larger values of $N$. These allow us to observe and explain some results in the limit of large $N$ or $N_{c}$. In Sec. IV, we compare the efficiency of different configurations of channels in conveying a particulate flux of given intensity. Finally, in Sec. V we summarize our results.

II Single channel Model

Particles with identical constant velocities, are injected into a channel of length $L$, according to a Poisson distribution of mean intensity $\lambda$. Given no blockage occurs, the particle transit time is $\tau=L/v$. An instantaneous blockage occurs if $N$ particles are simultaneously present, and lasts for time $\tau_{b}>\tau$. In the limit $\tau_{b}\rightarrow\infty$, there is no steady state and the exiting flux falls to zero Barré et al. (2015). Here, we focus instead on reversible blockages, during which, the $N$ particles are retained, and no more may enter the channel. After the deblocking time, the channel instantaneously releases all $N$ particles, resetting to the empty state, thereby allowing new ones to enter. The dynamics is therefore a recurring cycle of alternating open and closed states, that ultimately leads to a stationary state.

In an average steady state recurrence cycle, the channel is open for an average time $\langle t\rangle$ and blocked for a fixed time $\tau_{b}$. The stationary probability that the system is open is therefore

[TABLE]

where the denominator represents the total mean time of a recurrence. The stationary output flux is then given by the ratio of the mean number of particles released during one cycle to the cycle period,

[TABLE]

$\langle m\rangle$ is the mean number of output particles between two successive blockages. By equating the number of entering particles in one period to those exiting we obtain the following ‘number balance’

[TABLE]

Finally from the above three equations we deduce that

[TABLE]

The latter relation does not depend on the existence of a cycle, as it is the result of number conservation. The output flux is equal to the entering one minus the part that is rejected when the channel is in the closed state.

By taking the limit $\lambda\tau\ll 1$, the mean blockage time and the mean number of exiting particles between blockages behave asymptotically as $\langle t\rangle\gg\tau_{b}$ and $\langle m\rangle\gg N$, respectively. Therefore, for a given $\tau_{b}$, the probability that the channel is open is close to unity and the flux $j(\lambda)\simeq\lambda$. Blockages rarely occur at low $\lambda$. In this limit, the mean blockage time can be estimated by noting that a blockage occurs when a batch of particles enters in a finite duration $\tau$, leading to $\langle t\rangle=\tau\frac{(N-1)!}{(\lambda\tau)^{N}}$ Barré et al. (2015). Expanding $p_{o}(\lambda)$ to first order gives

[TABLE]

When $\lambda\tau\gg 1$, blockages are very frequent, and both the mean number of exiting particles between blockages, $\langle m\rangle$, and mean time between blockages, $\langle t\rangle$, approach zero. The resulting flux consists entirely of successive releases of the blocked particles, $j(\lambda)=N/\tau_{b}$ and $p_{o}(\lambda)\simeq 0$. In this limit, $\langle t\rangle$ corresponds to the time necessary for $N$ particles to enter an empty channel, $N/\lambda$. The open probability and the flux in this high intensity limit are therefore

[TABLE]

and

[TABLE]

II.1 Solvable models: $N\leq 3$

For small capacities, $N\in[1,2,3]$, the time evolution of the process can be expressed by analytically tractable differential or integro-differential equations Barré et al. (2013). For larger values of $N$, the time evolution cannot be solved by any known means.

We first consider $N=1$, which corresponds to a stochastic switch. The transit time $\tau$ is an irrelevant variable because no particle can exit the channel without having already blocked it. For $N>1$, it is possible for particles to pass through the channel without causing a blockage. Let $p_{o}(t)$ denote the probability that the channel is open at time $t$. Its time evolution obeys

[TABLE]

The loss term corresponds to the entrance of a particle at time $t$, while the channel is open, causing the channel to block. The gain term corresponds to the exit of a particle that became blocked at time $t-\tau_{b}$, with the subsequent reopening of the channel at time $t$.

The mean output flux at time $t$ is given by

[TABLE]

which corresponds to the release of a blocked particle that entered at $t-\tau_{b}$. Applying the time Laplace transform, $\tilde{f}(u)=\int_{0}^{\infty}dte^{-ut}f(t)$, to Eqs.(8) and (9) gives,

[TABLE]

and

[TABLE]

Expanding the denominator of Eq.(10) in terms of $\lambda e^{-u\tau_{b}}/(\lambda+u)$, allows one to easily invert the Laplace transform, term by term, giving,

[TABLE]

where $\theta(t)$ is the Heaviside function. The stationary open probability, $p_{o}(\lambda)$, and flux, $j(\lambda)$, can be obtained from Eqs.(10) and (11) by using $\tilde{f}(u)\simeq\frac{f(\lambda)}{u}$,

[TABLE]

and

[TABLE]

These results can be easily inferred from Eqs. (1) and (2) by setting $\langle m\rangle=0$ and $\langle t\rangle=1/\lambda$. The exiting particle flux is controlled by the incoming flux $\lambda$ and the time of blockage $\tau_{b}$.

We now consider the $N=2$ model, i.e. blockage occurs when two particles are simultaneously in the channel, for which exact results have already been obtained Barré et al. (2013). Here we propose an alternative, simpler derivation using the state probabilities of the channel. Let $p_{0}(t)$, $p_{1}(t)$ denote the probability that an open channel contains zero or one particle respectively, and $p_{2}(t)$ be the probability that it contains two particles and is therefore blocked. The time evolution of the process is given by

[TABLE]

with the following initial conditions:

[TABLE]

In Eq.(15), the loss term corresponds to the entrance of a particle in the empty channel at time $t$. The two gain terms $\lambda e^{-\lambda\tau}p_{0}(t-\tau)$ and $\lambda p_{1}(t-\tau_{b})$ correspond to a particle exiting the channel at time $t$ and a channel release (with a blockage occurring at time $t-\tau_{b}$), respectively. In Eq.(16), the two loss terms describe either a particle exiting the occupied channel at time $t$ or a particle entering the occupied channel. The gain term corresponds to a particle entering a free channel. In Eq.(17), the loss term corresponds to a channel release and the gain term to a particle entering a channel with one particle already inside. Summing the three equations verifies that the total probability is conserved: $p_{0}(t)+p_{1}(t)+p_{2}(t)=1$.

Taking the Laplace transform of Eqs.(15-17) gives

[TABLE]

These simultaneous equations may be solved to give

[TABLE]

where

[TABLE]

The mean exiting flux $j(t)$ is the sum of two contributions: the exit of a particle from an open channel and the release of of two particles from a closed channel. $j(t)$ is therefore given by

[TABLE]

By using Eqs.(22),(23), the Laplace transform of the output flux $\tilde{j}(u)$ is

[TABLE]

As expected, we recover the results of Ref.Barré et al. (2013) and the time-dependent mean flux can be obtained by a Laplace inversion of Eq.(27).

We here focus on the key quantities, namely the stationary probability $p_{o}(\lambda)$ that the system is open and the mean flux $j(\lambda)$. $p_{o}(\lambda)$ is the sum of the two stationary probabilities $p_{0}$ and $p_{1}$, each obtained by evaluating $\lim_{u\rightarrow 0}u\tilde{p}_{i}(u)$ with $i=0,1$:

[TABLE]

and

[TABLE]

Figure 1 displays $j(\lambda)$ versus $\lambda\tau$ for different integer values of $\tau_{b}$ from $3$ to $12$. The dashed lines correspond to the asymptotic values of the exiting flux $\lim_{\lambda->\infty}j(\lambda)=\frac{2}{\tau_{b}}$. One first observes that the stationary flux reaches the asymptotic values more rapidly as $\tau_{b}$ increases. Moreover, $j(\lambda)$ displays a maximum when $\tau_{b}$ is larger than $7$. It is possible to obtain the exact value of $\tau_{b}$ for which the flux $j(\lambda)$ displays a maximum at a finite value of $\lambda$Page et al. (2018) by solving $\frac{\partial j(\lambda)}{\partial\lambda}=0$. A real solution for $\lambda$ exists if $\tau_{b}>6.2$. Note that for $N=1$ the flux is always a monotonically increasing function of $\lambda$.

For $N=3$, the complete kinetic description of the model is cumbersome so we restrict our attention to the stationary quantities for which analytical expressions have been obtained Barré et al. (2015). In particular, the mean time to blockage starting from an empty channel is given by

[TABLE]

for $\lambda\tau>2\ln(2)$ and

[TABLE]

for $\lambda\tau<2\ln(2)$, where $g=\sqrt{|1-4e^{-\lambda\tau}|}$ and $\nu=\frac{\lambda\tau}{2}$ (note that these correct the expressions given in Barré et al. (2015)).

The two stationary quantities $p_{o}(\lambda)$ and $j(\lambda)$ are obtained by inserting this result in Eqs.(1) and(2). Fig. 2 displays $j(\lambda)$ as a function of $\lambda\tau$ for different values of $\tau_{b}$. There are several differences compared with the $N=2$ model. First, a maximum exiting flux occurs if the blockage time $\tau_{b}>3.6$, which is significantly smaller for $N=2$ ($\tau_{b}>6.2$). Second, the asymptotic values are reached at a lower value of $\lambda$, and finally, the intensity $\lambda$ at which $j(\lambda)$ is maximum is also shifted towards larger intensity.

II.2 Simulation results: $N>3$

As a result of strong time correlations between the transiting particles, it is not possible to obtain analytic solutions for $N>3$. We therefore used numerical simulations to investigate these cases. In order to benchmark our code, we compared the simulation results for the stationary flux for $N=3$ with the exact expressions for three different values of $\tau_{b}$. In Fig. 2 we observe perfect agreement between the analytical expressions and the simulation results.

Figure 3 displays the stationary exiting flux as a function of the intensity for different values of $N$ and $\tau_{b}=4$. At low intensity, the flux increases linearly, and at high intensity the asymptotic behavior of the simulation results is well-desribed by Eq. (7). The behavior in the intermediate region is due to complex dynamics that alternates between blockages and sequences of uninterrupted transport. For $N>3$, the stationary flux $j(\lambda)$ may display a maximum at finite $\lambda\tau$, whose amplitude increases with $N$. The stationary flux also exhibits a minimum which is always smaller than the asymptotic value, $N/\tau_{b}$.

We performed a systematic study of the behavior of the exiting flux as a function of $\lambda,\tau_{b}$ and $N$. The flux always displays a maximum when $\tau_{b}$ exceeds a threshold value. Figure 4 shows that the critical value of $\tau_{b}$ decreases rapidly with $N$, showing that the feature observed in Fig. 1 is very general and occurs for smaller values of $\tau_{b}$ when $N$ increases. For $\tau_{b}$ below the critical value, the stationary flux is a monotonically increasing function of $\lambda\tau$. The right panel of Fig. 4 shows the values of $\lambda\tau$ corresponding to the critical values of $\tau_{b}$.

The behavior of the open probability $p_{o}(\lambda)$, shown in Fig. 5, is consistent with Eq. (4). In particular one observes the appearance of a plateau whose length increases with $N$ (roughly as $N/2$). This corresponds to the situation where blockage events are rare and the output flux is close to $\lambda$.

II.3 Markovian versus non-Markovian models

The physical assumption of constant transit and deblocking times $\tau$ and $\tau_{b}$, respectively, is responsible for strong memory effects which prevent analytical solutions for general $N$ from being obtained. We therefore recently introduced Markovian models Page et al. (2018); Barré et al. (2018), where the average transit and deblocking times are stochastic variables given by exponential distributions of intensity $\mu$ and $\mu^{*}$, respectively. The kinetic description of the Markovian model is given by a set of differential equations for the time evolution of the state probabilities $P(i,t)$ with $i\in[0\cdots N]$ giving the number of particles in the channel. Unlike the non-Markovian model, analytic solutions for the steady state properties can be obtained for arbitrary $N$ (some generalizations of the Markovian models for which time-dependent solutions can be obtained and could be investigated in the future Leonenko (2009); Escobar et al. (2002)).

The channel is open for an mean time $\langle t\rangle$ and blocked for a mean time $1/\mu^{*}$. The stationary flux is obtained using the previously employed recurrence arguments, giving

[TABLE]

The average time for which the Markovian system is open in a recurrence cycle is Page et al. (2018)

[TABLE]

To compare the two models, $\mu$ and $\mu^{*}$ must be related to $\tau$, $\tau_{b}$ and $\lambda$. Eq. (32) with Eq. (2) shows that $\mu^{*}$ must equal $1/\tau_{b}$. To obtain an expression for $\mu$, we consider the system’s behavior at low and high intensity. When $\lambda\tau<<1$, the non-Markovian transit time is equal to $\tau$. The mean transit time is $1/\mu$ in the the Markovian model. A first approach is to therefore set $\mu=1/\tau$. When $\lambda\tau>>1$, we expect $\mu$ to decrease to zero. Figure 6 shows that the stationary flux of the Markovian model is always larger than that of the non-Markovian model. Even though the two models behave similarly at small and large input intensity, they increasingly deviate for intermediate intensities with increasing $N$.

To obtain an exact mapping (in the steady state) we equate the mean blocking time of the two models. For $N=2$ we equate $\langle t\rangle$ given by Eq. (33) with the result for the non-Markovian model Gabrielli et al. (2013), $\langle t\rangle=(2-e^{-\lambda\tau})/(\lambda(1-e^{-\lambda\tau}))$. The expressions are identical when

[TABLE]

With this mapping, we recover the aforementioned expected limiting behaviour for both extremes of entering flux intensity. We emphasize that the transient regimes of the two models are different (See the Appendix A for a similar model where time-dependent analytic solutions are obtained).

The same procedure can be carried out for $N=3$ using Eqs. (30) and (31), but the resulting expression for $\mu$ is considerably more complex. For general $N$ we therefore propose the following ansatz, taking a similar form as the mapping for $N=2$:

[TABLE]

which behaves as $1/\tau$ at low intensity and approaches zero exponentially at large intensity. Substituting Eq.(35) into Eq.(32), produces a lower bound of the stationary flux (full curves). Furthermore, the maximum of the flux is underestimated and shifted to a smaller intensity than in the non-Markovian model. For $N=3$ and $N=4$, the curves are very close to the results of the non-Markovian model. For $N>4$, the ansatz leads to a significant underestimation of the exiting flux for small $\lambda\tau$.

III Bundle model

We now consider a bundle of $N_{c}$ identical channels. Each channel has the same properties as the single channel model above, i.e. blockage occurs when $N$ particles are present in a channel at the same time. In the following we assume that the total intensity, $\Lambda=\lambda N_{c}$, is constant and is equally distributed over the open channels. Thus, after $k$ blockages the intensity on each of these open $N_{c}-k$ channels is

[TABLE]

Since a blocked channel releases all particles after finite time $\tau_{b}$, the system’s mean output flux evolves towards a non-zero stationary value. The bundle has two states: open in which at least one of the constituent channels is open and closed if all the constituent channels are blocked. If a particle arrives while the bundle is in the latter state, it is rejected.

Equations (1,2,3) cannot be applied to the channel bundle in the steady state, as it does not cycle between closed and empty states for finite intensity $\Lambda$. In the limit of very large intensities, however, we have

[TABLE]

and

[TABLE]

where $\langle M\rangle$ is the mean number of exiting particles that are not the due to blockage releases, and

[TABLE]

In this limit the intensity is so high that all channels block instantaneously and simultaneously and the blocked particles are released after a time $\tau_{b}$. The exiting flux is entirely the result of these releases.

The analogue of Eq. (4),

[TABLE]

is valid for arbitrary intensity since, as for the single channel case, it is result of the conservation of particle number.

III.1 Exact solution: $N=1$

When $N>1$ a particle may traverse the channel in a time $\tau$ without causing a blockage. In comparison, the $N=1$ model is singular as no unimpeded transit is possible: each entering particle causes a blockage that lasts for a fixed time, $\tau_{b}$. The variable $\tau$ is thus absent in this model.

Despite the relative simplicity of the $N=1$ model, its dynamics cannot be written as a system of differential equations for the state probabilities $P(i,t)$, where $i$ denotes the number of blocked channels at time $t$ (in Appendix A the full time dependent solution for $N_{c}=2$ is presented). However, in the stationary state, by applying detailed balance (known as the “rate up - rate down” principle in queuing theory), one has

[TABLE]

Solving the difference equation and applying conservation of the total probability leads to

[TABLE]

The stationary exiting flux is given by Eq. (40) with

[TABLE]

The result can be written in the form

[TABLE]

where $\Gamma(n,x)$ is the incomplete gamma function. The asymptotic behavior at small intensity $\Lambda$ is

[TABLE]

whereas at large intensity the flux behaves as,

[TABLE]

whose leading term is in accordance with Eq. (39). For all values of $N_{c}$, $J(\Lambda)$ is always a monotonically increasing function of $\Lambda$.

We note that the expression for $P(N_{c})$ is Erlang’s first formula Takacs (1969); Medhi (1991) for a stochastic queuing process with $N_{c}$ servers with exponential entry and service time distributions under the condition that when all servers are busy an arrival is rejected. Both models are ‘birth-and-death’ processes that have the same stationary solution. Their transient regimes, however, are significantly different. See Appendix A.

III.2 Simulation results: $N>1$

For the multichannel models, no exact solution can be obtained for $N>1$. Therefore we have performed numerical simulations to obtain the stationary exiting flux, $J$ and the stationary probability that at least one channel is open, $P_{o}$ for bundles composed of different numbers of channels with increasing capacity $N$ and for $\tau_{b}=4\tau$. All quantities were investigated as a function of the mean incoming flux $\Lambda$. As discussed in the previous section, the stationary flux rapidly displays a maximum at a finite value of $\lambda\tau$ when $N>1$.

Figure 7 shows $J$ as a function of $\Lambda\tau$ for $N_{c}=2,10$ and for $N=2,3,4,6,10$. When $\Lambda\tau\ll 1$, the rate of incoming particles is very small and the finite capacity of the channel is rarely reached, meaning that blockage events are scarce. The stationary exiting flux is therefore equal to the input flux, $J\simeq\Lambda$. This behavior is observed for a larger range of $\Lambda\tau$ for larger values of $N$ and $N_{c}$ .

Figure 8 shows the rescaled flux $J/N_{c}$ versus the rescaled intensity $\Lambda/N_{c}$ for different values of $N_{c}$. In the low intensity regime $J$ is equal to $\Lambda$ and displays a finite discontinuity at $\Lambda_{c}$. At high intensities the curve evolves towards an asymptote and is quite well described by

[TABLE]

We observe an abrupt change of kinetic behavior: Below the critical value $\Lambda_{c}$, almost all particles cross the bundle without triggering a significant number of blockages, whereas for larger $\Lambda$, all channels are closed and the stationary flux is essentially given by the release of blocked particles.

Fig. 9 shows the rescaled flux $J/N$ versus the rescaled intensity $\Lambda/N$ for different values of $N$ for a given $N_{c}=10$. At a low input intensity, the exiting flux $J$ is equal to $\Lambda$ until it reaches a maximum close to a critical value that closely follows a logarithmic law, as shown in the inset of Fig. 9. For higher values of input intensity, the rescaled exiting flux, for all values of $N$, rapidly collapses to a single curve, whose best fit is again given by Eq. (47).

IV Flux optimization

Here we compare the stationary flux of a bundle of coupled channels with that of a bundle of uncoupled channels and one high capacity (HC) channel. The transport efficiency is measured by the difference in output flux, $\Delta J$. The systems are chosen so that in the limits of low and high input flux intensity, $\Delta J=0$. In the low intensity limit, since blocking events are rare, the exiting flux is equal the input flux $\Lambda$, irrespective of the configuration. In the high intensity limit, Eqs.(6) and (39) demonstrate that the exiting fluxes of the single high capacity, bundled uncoupled or coupled channels are also be equal. Since the stationary flux of a bundle of coupled channels displays non-trivial behavior with increasing $N_{c}$ and $N$, we therefore expect non-trivial behavior of $\Delta J$.

IV.1 Coupled versus uncoupled channels

We first compare a channel bundle composed of $N_{c}$ channels, each of capacity $N$. The entering flux, $\Lambda$, is equally distributed over the coupled open channels. In contrast, the $N_{c}$ independent channels, each of capacity $N$, each receive an incoming flux of intensity $\Lambda/N_{c}$. Fig. 10 illustrates the configurations compared. The difference in the output flux of the two configurations is defined as:

[TABLE]

Simulation results comparing the output flux of each configuration for $N_{c}=2,5$, $\tau_{b}=4$, are shown in Fig. 11. The differences are shown in the inset of the figure. At low intensity, for all configurations, the output flux is approximately equal for each setup, and linearly increases with $\Lambda$ until a critical value, which itself is a monotonically increasing function of $N$.

The behavior for $N=1$ can be understood quantitatively using the results of Sec. III.1 and Eq. (14). At low intensity, the flux difference is

[TABLE]

and at high density we find

[TABLE]

and one can confirm that $\Delta J>0$ for $0<\Lambda<\infty$. We conclude that the coupled channels are always more efficient than the uncoupled ones. The difference is maximized for a finite value of $\Lambda$.

For all $N>1$, we note the appearance of two maxima in the flux difference with an intervening minimum. The increased complexity is due to the presence of two characteristic times, the transit time $\tau$ and the blockage time $\tau_{b}$ (while the $N=1$ system has only the latter).

Fig. 12 shows the flux difference between the bundle configurations at $Nc=2$, $N=2$, as a function of intensity of entering flux, $\Lambda$, for different values of $\tau_{b}$. For $\tau_{b}>2$ the behavior is more complex after the first maximum, with the appearance of a minimum followed by a second maximum before tending towards zero.

IV.2 Single HC channel verses coupled LC channels

The transport efficiency of a single ’high capacity’ (HC) channel, with a bundle of several coupled channels, of a proportionately reduced capacity, is now compared. The two configurations are illustrated in Fig. 13. Figure 14 shows the difference in the output stationary state flux between a single HC channel and a bundle of $N_{c}=2,5$ coupled LC channels, with different capacities. The flux difference in this case is

[TABLE]

For all $N$, the flux difference displays a minimum, followed by a maximum, before tending towards zero for increasing intensity. The amplitudes of the maxima are always greater than those of the minima.

V Discussion

We have presented a model of blockage in channel bundles that may be relevant for various physical processes. The bundle consists of $N_{c}$ constituent channels, each with a capacity of $N$. A particle transits through an open channel in time $\tau$, but if $N$ particles are simultaneously present in a channel, it is blocked for a time $\tau_{b}$, before being emptied. While blocked, the entering flux is redistributed over the remaining open channels. A bundle of channels is open if at least one of its constituent channels is not blocked. If the entering stream is of constant intensity the bundle evolves to a stationary state with a steady exiting flux, or throughput, that depends on the intensity, $\tau_{b}$ and $N$. In the steady state the exiting flux is simply related to the probability that the bundle is open. For a single channel with capacity $N>1$ the exiting flux displays a maximum value at finite intensity if $\tau_{b}$ is sufficiently large. If not, the exiting flux increases monotonically with the intensity. A Markovian model with stochastic transit and blockage times, inspired by queuing theory, can be made to display the same steady state behavior with an appropriate mapping between the two models’ parameters. The transient behavior is, however, quite different. This suggests that, in the steady state, the details of the transport mechanisms and the triggered releases are not important. For large $N$, the models display an abrupt change from a state with few blockages to one in which the bundle is permanently blocked and the output flux is entirely due to the release of blocked particles. This behavior raises new questions about whether more general relationships describing the abrupt transitions in dynamics may be obtained for general $N$ and $N_{c}$. The transport efficiency of a bundle in which the entering flux is equally distributed over the open channels was also compared with a bundle composed of independent channels. For $N=1$, the coupled channels always have a higher throughput, but for larger values of $N$ the behavior is more complex.

Acknowledgements.

We thank Jacques Resing, Karim Guerouate et Valentin Wiener for useful discussions.

Appendix A Full Solution of the Channel Bundle with $N_{c}=2$, $N=1$

This model can be solved exactly by using an approach similar to that used to solve the single channel model with $N=3$ Barré et al. (2015). This requires introducing the partial probabilities $Q(n,t)$ that at time $t$, $n$ particles have entered such that at least one channel is open. The open probability of the bundle is then given by

[TABLE]

Since blockage is not possible for $n=0,1$ the first two probabilities are trivial:

[TABLE]

For $n\geq 2$, the probability $Q(n,t)$ can be written as:

[TABLE]

The $n$-fold integral corresponds to all events of $n$ incoming articles at the entry of the channel with the associated constraint ($\delta$ function) which occurs between [math] and $t$. The first product of Heaviside functions expresses the fact that a new particle can enter if both channels are not blocked, which imposes the condition that the time since the entry of the second last particle is larger than $\tau_{b}$; The last product corresponds to both complementary situations: either the last entered particle leads to a complete blockage of two channels (The exponential factor expresses that no particle can enter in a duration $\tau-t_{k}$) or only one channel is blocked by the last entered particle. The last Heaviside function imposes the requirement that at least one channel be open at time $t$.

Taking the Laplace transform of $Q(n,t)$, the integral over $t_{0}$ is trivial and the integral over $t^{\prime}$ must be split into two parts, which finally gives

[TABLE]

with $t_{n-1}\equiv t$. The function $r(n-1,t,u)$ is defined as

[TABLE]

By using Eq.(57), one infers the recurrence relation between $r(n,u,t)$ and $r(n-1,u,t)$

[TABLE]

with the initial condition $r(1,t,u)=1$

We define the generating function $G(z,t,u)$

[TABLE]

Multiplying Eq. (58) by $z^{n}$ and summing over $n$, gives the generating function $G(z,t,u)$ obeys to a similar equation

[TABLE]

For $t\geq\tau_{b}$, the right-hand side of Eq.(A) is independent of $t$, which gives that $G(z,t,u)$ is then independent of $t$. For $t<\tau_{b}$, by taking the two first derivatives of Eq.A with respect to $t$, one can rewrite the integral equation, Eq.(A) as an ordinary differential equation

[TABLE]

whose solution is given by

[TABLE]

where $r_{1}$ and $r_{2}$ are the solutions of the characteristic equation

[TABLE]

The functions $A(u)$ and $B(u)$ are determined by using the boundary conditions. At $t=0$, by using Eq. (A), $G(z,0,u)$ is expressed as

[TABLE]

which gives

[TABLE]

Now, by using the first derivative of Eq. (A) at $t=\tau_{b}$, one obtains the second boundary equation

[TABLE]

and by using Eq.(62), one obtains

[TABLE]

In the Laplace space, the open probability $\tilde{P}_{o}(u)$ is given by

[TABLE]

By using Eq.(A) and Eq.(62), this can be expressed as

[TABLE]

The functions $A(u)$ and $B(u)$ can be determined by using the boundary conditions, but the lengthy expressions are not displayed here. The output flux $J(t)$ is entirely due to the release of blocked particles and is therefore given by

[TABLE]

The stationary open probability is $P_{o}(\infty)=\lim_{u\rightarrow 0}u\tilde{P}_{o}(u)$. By using Eq.(69), the stationary open probability $P_{o}$ is then given by

[TABLE]

After a tedious but straightforward calculation, one finally obtains

[TABLE]

which corresponds to the Erlang formula, Eq. (43).

Let us now compare with the stochastic model with $N_{c}=2$ and $N=1$ with parameters $\Lambda$ and $\mu^{*}$. The state of the system is defined by the probabilities $P_{0}(t)$, $P_{2}(t)$, $P_{2}(t)$ of having zero, one and two particles in the bundle at time $t$, respectively. They evolve according to the coupled differential equations

[TABLE]

Taking the Laplace transform, and solving the linear system of algebraic equations, one obtains the Laplace transform of the open probability, $\tilde{P}_{o}(u)$

[TABLE]

Finally, calculating the inverse Laplace transform and substituting $\mu^{*}=1/\tau_{b}$, one obtains

[TABLE]

with $\tau_{b}^{\prime}=\frac{\tau_{b}}{\sqrt{1+4\Lambda\tau_{b}}}$.

Figure 15 shows the time evolution of the open probability, $P_{o}(t)$, of the non-Markovian model and of the stochastic model for $\Lambda=2$ and $\tau_{b}=4$. While both models converge to the same stationary state (with the same open probability and output fluxes), the relaxation is significantly different, except at short time (typically for $t<1/\Lambda$). In particular, for the non-stochastic model, the deterministic blockage release mechanism leads to slow, oscillatory relaxation towards the stationary state, whereas the relaxation is roughly exponential and faster for the stochastic model.