{'fraction_non_alphanumeric': 0.03740267593922443, 'fraction_numerical': 0.012699372590520825, 'mean_word_length': 4.356385132399384, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Maximum power operation of interacting molecular motors

We study the mechanical and thermodynamic properties of different traffic models for kinesin which are relevant in biological and experimental contexts. We find that motor-motor interactions play a fundamental role by enhancing the thermodynamic efficiency at maximum power of the motors, as compared to the non-interacting system, in a wide range of biologically compatible scenarios. We furthermore consider the case where the motor-motor interaction directly affects the internal chemical cycle and investigate the effect on the system dynamics and thermodynamics.


I. INTRODUCTION
Molecular motors are biological machines that harness chemical energy and convert it into motion or useful mechanical work. These molecular machines are responsible for performing tasks as diverse as DNA replication and repair, RNA transcription, protein synthesis and intracellular transport [1]. The development of sophisticated single-molecule experimental techniques and the emergence of several theoretical frameworks for modelling single motors has over the past two decades led to a great amount of accumulated knowledge on mechanics and thermodynamics of single molecular machines [2][3][4][5]. However, many motors, such as, e.g., kinesin or dynein motors involved in cellular cargo transport, move on crowded filamenteous tracks known as microtubules where they can encounter other motors. Such encounters give rise to non-negligible motor-motor interactions and affect the resulting motion of the motors. Molecular motor traffic is therefore an important and widely studied phenomenon, which is typically modelled by using exclusion processes on lattices [6][7][8][9][10][11][12][13][14][15][16].
In this paper we extend and study in detail two different traffic models for kinesin, which were first used in [17] to study the efficiency of kinesin operating under external mechanical load force in the maximum power regime. Model I represents the simplest possible description of molecular motor traffic which is obtained by neglecting the internal conformational states of the molecular motor and modelling the system as a standard asymmetric simple exclusion process (ASEP) [16]. In this framework, the motor dynamics is represented as a continuous-time, stochastic jumping process of particles on a discrete lattice. Kinesin is powered by the ATP hydrolysis reaction, in which an adenosine triphosphate (ATP) molecule is hydrolyzed into an adenosine diphosphate (ADP) and a phosphate (P i ) molecule. The forward stepping is thus associated with ATP hydrolysis, while the backward jumps must proceed through ATP synthesis in this minimal model in order to obtain a thermodynamically consistent description. The particle stepping is subject to an exclusion rule; if the particle attempts to step either forward or backward to a neighbouring site that is already occupied by another motor, the step will be rejected. The ASEP and extensions hereof have been used in several studies of molecular motor traffic [8,9,15].
However, all these models do not take into account that the stepping of a molecular motor is a complex process consisting of a series of transitions between different internal motor states. Especially, it is known that backward steps can occur as a consequence of ATP hydrolysis as well as ATP synthesis [18,19]. It is therefore crucial to incorporate detailed kinetic models for the internal conformations into models of molecular traffic, as stressed in several works [10][11][12][13]. This goal is achieved within model II that is an extension of the ASEP and combines the thermodynamically consistent descriptions of kinesin's stepping presented in [18,19] with the standard exclusion process formalism. To our best knowledge, our model II, as introduced in ref. [17], is the first of its kind to incorporate a kinetic motor model with several mechanochemical cycles into a description of molecular motor traffic.
An obvious quantity to consider when dealing with motors of any kind is their thermodynamic efficiency. For traditional heat engines operating between two thermal baths, the efficiency is constrained by the well-known Carnot's law. Molecular machines, on the other hand, operate in environments at constant temperature, and their efficiency is thus bounded by the thermodynamic limit which is equal to 1. However, both in the case of heat engines and isothermal machines, this lossless limit can only be achieved for reversible processes, or infinitely close to thermodynamic equilibrium. The corresponding power output is therefore zero and thus of limited practical interest. As a consequence, the concept of efficiency at maximum power (EMP) in the context of microscopic engines has recently received considerable attention in the literature since it provides a quantitative measure of the trade-off between power and efficiency [20][21][22][23]. While these works were concerned with the EMP of single motors, the EMP in systems of interacting molecular motors has only been considered in ref. [17]. We exploit model I and II to investigate the effect of mutual motor interactions on the kinetics and thermodynamics in the maximum power regime, and in particular on the EMP. Furthermore, for both models we consider two different types of boundary conditions. The case of open boundary conditions where the system primarily exchanges particles with the reservoir at the ends of the filament appears to be relevant for intracellular transport, where the cell products have to be transported over relatively long distances between different cellular regions. On the other hand, bulk-dominated binding and unbinding, a scenario known as Langmuir kinetics [12,15], is relevant for experimental in vitro studies of molecular motor traffic. Moreover, we consider the case where, besides the mechanical transitions, also the chemical transitions are affected by the ASEP exclusion rule, and discuss the consequences of this, so-called, strong exclusion rule mimicking the possibility that molecular machines can shut down their motor in the presence of high traffic. The idea that non-steric collective effects can suppress dissipation in molecular motor traffic and thereby increase efficiency has been proposed previously in [7,24].
We solve the appropriate equations of motion for the various combinations of models, boundary conditions and exclusion types to obtain dynamical quantities such as velocity and hydrolysis rate for fixed chemical input and externally applied mechanical load force. This enables us to characterize the maximal power regime by considering the variation of the power output with the load force for a given value of the chemical driving. Once the optimal load force is known, all the relevant quantities such as EMP or velocity at maximum power can be obtained straightforwardly from the previously calculated dynamics.
We find the remarkable effect that the EMP in manymotor systems is enhanced, as compared to the single motor case, by the mutual exclusion interactions under a wide range of biologically applicable conditions. The enhancement is caused by an altered characteristic forcevelocity relation as a consequence of steric motor-motor interactions. Furthermore, in the case of open boundary conditions we observe for both models a trade-off behaviour between efficiency and velocity in the maximum power regime when the force dependence of the unbinding mechanism is altered. This observation has interesting prospects in terms of switching between fast and efficient transport in artificial many-motor systems. In this respect, the so-called 'molecular spiders' are promising examples of highly tunable artificial molecular motors [25].
The paper is organized as follows. Sections II and III are concerned with model I. In section II we introduce model I and describe the motor mechanics in the meanfield limit under different boundary conditions, while sec. III provides a discussion of thermodynamic efficiency, velocity and motor density in the maximum power regime. The following two sections deal with the more detailed model II. In subsections IV A and IV B we present the two variations of the model with internal states. We obtain an analytic solution for the open-boundary problem in the mean-field approximation by employing the maxi- The motor is thus represented as a particle that can jump forward or backward with rates p and q, respectively, if the target side is accessible. The lattice constant corresponds to the motor step size and is denoted by a. a) Periodic boundary conditions (PBC). b) PBC and Langmuir kinetics: If a lattice site is empty, a motor can attach to the site with rate α. Conversely, a motor bound to the track at a given site can detach from the corresponding site with rate β. c) Open boundary conditions (OBC): Motors can bind to the filament at the left end with rate α, if the first lattice site is empty, and detach from the right end with rate β. mal current principle (MCP) [26] as discussed in section IV A 2. In section V we investigate the EMP of interacting kinesin motors calculated within model II and the corresponding velocity and motor density at maximum power. Finally, section VI discusses the biological implications of our findings and provides some concluding remarks.

II. MODEL I
We start out by considering a simplified motor traffic model, in which we neglect the internal conformational states of a molecular motor and represent the stepping kinetics as a Poissonian process, see fig. 1. The motor thus moves on the lattice track by performing forward or backward steps of length a with jumping rates p or q, respectively. However, the presence of multiple motors on the microtubule gives rise to motor-motor interactions, which we take to be steric in our description. Hence, the motor dynamics is modified by an exclusion rule implying that the stepping can only proceed if the target site in question is unoccupied.
The chemical free energy of the ATP hydrolysis reaction driving the motor can be written in terms of the reactant concentrations as [27] where K eq is the equilibrium constant of the hydrolysis reaction, and T = 4.1 pN nm is the room temperature (here and in the following, we take k B = 1). The released chemical energy is used to perform mechanical work against a constant external load force f < 0, which represents the effect of the cargo particle on the motor. Since the model motor is tightly coupled, i.e. one ATP molecule is hydrolyzed (synthesized) for every forward (backward) mechanical step, the input work and the output work when completing a forward step are w in = ∆µ and w out = −f a, respectively. Thermodynamic consistency thus requires that the jumping rates fulfill the local detailed balance condition [28], This relation allows for a parametrization of the rates which for constant ADP and phosphate concentrations reads where θ is the load distribution factor representing the coupling of the mechanical force to the kinetic model parameters.
In the mean-field description the spatial correlations between the different sites are neglected, and the equation of motion governing the bulk dynamics on the lattice readsρ where ρ j is the motor density at lattice site j = 1, . . . , N , and J i+1/2 is the mean-field probability current expressed as We proceed by considering different types of boundary conditions for the system and their effect on the resulting motor dynamics. We focus on periodic boundary conditions (PBC), periodic boundary conditions with Langmuir kinetics (PBC-LK) describing attachment and detachment dynamics in the bulk of the filament, and open boundary conditions (OBC) where motor binding and unbinding occurs at the filament ends. We do not consider the most general, yet technically more involved, case of open boundaries and Langmuir kinetics, since in relevant biological and experimental setups the system can be well approximated by one of the above-mentioned special cases as we argue in sec. VI.

A. Periodic boundaries (PBC)
The simplest description is obtained when neglecting the exchange of molecular motors with a reservoir, i.e. the fact that the motors can attach and detach from the filamenteous track they are moving on. Such a situation corresponds to employing periodic boundary conditions (PBC), see fig. 1a. In this scenario all the lattice sites become equivalent, and the steady-state motor density is thus independent of the lattice site, i.e. ρ j ≡ ρ for all j. The master equation (4) is hence trivially satisfied in the steady-state.
The probability current takes the form, cf. eq. (5), where the density ρ serves as a free parameter. Here, we have introduced the single motor velocity in the absence of interactions and particle exchange with reservoirs, v 0 = a(p − q). The motor velocity as a function of ρ is hence given by Fig. 2 illustrates J PBC and v PBC for ∆µ = 20T and for three different values of the load force f . In the absence of internal motor states, the molecular traffic model exhibits particle-hole symmetry, i.e. a symmetry upon interchanging ρ and (1 − ρ), and J PBC is thus symmetric around ρ = 1/2 as can be seen in fig. 2a.

B. Periodic boundaries and Langmuir kinetics (PBC-LK)
Next, we consider the effect of binding and unbinding from the filament on the periodic system by introducing Langmuir kinetics (LK), see fig. 1b. In addition to the on-lattice dynamics, a motor can now bind (unbind) from the filament with rate α (β). The master equation in the steady-state thus assumes the form As a result, the motor current remains unchanged, i.e. J LK PBC = J PBC , while the motor density attains the Langmuir equilibrium value, ρ LK = α/(α + β). We note that the single motor velocity in a non-interacting system with LK is equal to v 0 , since the interaction with the reservoir only affects the motor run length and not the propagation along the track. The motor thus progresses with the velocity along the filament. Furthermore, we note that the motor velocity is a linear function of the Langmuir density, see also fig. 2b, which is not the case for the more detailed model IIb as discussed in sec. IV B. In general, the detachment rate β increases with opposing loads [29,30], and we take the force dependence to be [18,31] where β 0 is the zero-force unbinding rate, and the parameter φ quantifies the effect of the force on the detachment mechanism. As a result, the motor density is a decreasing function of the applied load, and the motor dynamics in the presence of interactions differs from that of noninteracting motors.

C. Open boundaries (OBC)
Finally, we turn our attention to the problem of open boundary conditions (OBC). In this case we assume that the system only exchanges particles with the particle reservoirs at the lattice boundaries as depicted in fig. 1c. Motors can bind to the first lattice site with rate α, if the site is unoccupied, and unbind from the last lattice site with rate β. The master equations for the boundary sites are thusρ where J j±1/2 is defined in eq. (5). In the steady-state limit, current conservation entails that the probability current is independent of the lattice site, i.e. J j±1 ≡ J.
In the following we consider the phase diagram for the (mean-field) probability current J in the thermodynamic limit. Since the typical microtubule length, L = 10 µm [27], is long compared to the kinesin motor step size a = 8 nm, the thermodynamic limit defined by N = L/a → ∞ constitutes a good approximation to the finitesize dynamics. The phase diagram can be obtained by several mean-field techniques such as recursion relations, hydrodynamic equations or the maximal current hypothesis [26,[32][33][34] and happens to coincide with the phase diagram obtained in the thermodynamic limit from the exact analytical solution of the ASEP [16]. The diagram consists of three regions termed the lowdensity (LD) phase, the high-density (HD) phase and the maximal current (MC) phase, respectively. In the LD phase the bulk density is determined by the attachment dynamics and equals ρ LD = α/(p − q). The density in the HD phase is dictated by the detachment dynamics and takes the value ρ HD = 1 − β/(p − q). Finally, in the MC phase the density is ρ MC = 1/2 and is independent of α and β. The LD-MC and the HD-MC boundaries are characterized by the critical values α c and β c , respectively, and the particle-hole symmetry implies that α c = β c = (p − q)/2. At the LD-HD boundary, given by α = β and α, β ≤ (p − q)/2, the low-density and the high-density regions coexist on the lattice and are separated by a fluctuating domain wall. Hence, the location of the phase transitions varies with the chemical input ∆µ and the load force f , since the hopping rates p and q are functions of ∆µ and f , cf. eq. (3). In fig. 3 the phase diagram is depicted for ∆µ = 20T and for different values of f . We note that the MC phase grows with increasing values of the load force. The motor density thus varies with f even when the detachment dynamics is load independent, i.e. φ = 0, cf. eq. (10), as opposed to the system with PBC-LK. The probability current and the motor velocity are obtained from the phase diagram as J = J PBC (ρ i ) and v = v PBC (ρ i ), respectively, where i denotes the appropriate phase, i = LD, HD, MC. In summary, using eq. (7) for v PBC the single particle velocity for the interacting system reads  It is worth noting that v is always smaller than v 0 and attains its maximum in the LD phase as expected, since the jamming is minimal in this phase. Moreover, in the HD phase where it reaches the smallest possible value, the velocity only depends on the detachment rate β.

III. MODEL I: MAXIMUM POWER REGIME
The goal of this section is to investigate the operation of our model machine in the maximum power regime when the particle exchange with the reservoirs occurs either at the ends of the filament (OBC) or in the bulk of the system (PBC-LK) as described above. We start out by considering the efficiency at maximum power (EMP) in order to characterize the balance between the output power and thermodynamic efficiency in the presence of motor-motor interactions. Next, we study the corresponding motor velocity and density.

A. Efficiency at maximum power
Since the chemical and mechanical currents are tightly coupled in model I, the hydrolysis rate r describing the number of consumed ATP molecules per time unit is proportional to the motor velocity, i.e. r = v/a. The input power P in = ∆µr is thus proportional to the delivered output power P out = −f v, and the efficiency η of the Pout (pN nm/s) motor is independent of the mechanical and chemical currents, Here, f s = −∆µ/a is the stalling force of the motor for which the forward and backward rates are equal, p = q, and the velocity vanishes. It is clear from eq. (13) that the maximum efficiency, η = 1, is obtained under stalling conditions corresponding to equilibrium where the power output vanishes. A more relevant quantity to consider is therefore the efficiency at maximum power as discussed in the Introduction. For fixed ∆µ we thus calculate the motor efficiency when the power is optimized with respect to the load force. The equation for the maximizing force f * thus reads and the EMP is simply obtained as η * = f * /f s . The above procedure is then repeated for increasing values of the ATP concentration, and hence ∆µ. It is worth noting that the concept of EMP is different from the maximal possible efficiency at fixed ∆µ, as given by the condition ∂η/∂f = 0 and studied in, e.g., ref. [19]. The results obtained for the EMP for model I are reported in fig. 4 for different values of the load sharing parameter θ and for different load dependencies of the detachment rate characterized by the parameter φ. The EMP, η * LK , for the system with PBC-LK is shown in fig. 4a together with the EMP for the corresponding non-interacting system, η * 0 . We observe an increase in η * LK with respect to η * 0 due to a change in the characteristic force-velocity relation v(f ), and hence P out (f ), as illustrated in fig. 5. Since both v 0 and ρ LK are decreasing functions of the load, it follows from eq. (9) that the velocity decreases more slowly with the force in the presence of steric interactions. This, in turn, leads to a higher optimal force f * , as compared to the non-interacting system, and hence to a higher EMP. The enhancement of the non-int

EMP occurs for all values of θ and all non-zero values of
Hence, larger values of φ result in a greater effect of motor-motor interactions on the EMP. Furthermore, the increase in the EMP is more pronounced for small values of θ, as has been discussed previously in [17]. In the linear regime given by ∆µ/T → 0 we recover the well-known result η * → 1/2 for tightly coupled systems [20][21][22]. Moreover, η * LK → η * 0 for ∆µ/T → ∞, since the single-motor velocity v 0 is the dominating contribution to v LK , eq. (9), for large ∆µ. Fig. 4b shows the EMP, η * OBC , for the system interacting with the environment at the boundaries together with the corresponding EMP in the absence of interactions, η * 0 . The EMP exhibits different behaviour with ∆µ depending on the value of φ. If the detachment process depends weakly on the load force (or in the limit φ = 0), the EMP displays two different regimes. When the chemical free energy is smaller than the critical value ∆µ c,1 (see App. B for a discussion of the critical values), the maximum power is reached in the MC phase. Since the velocity in this phase is proportional to v 0 , cf. eq. (12), we obtain η * OBC = η * 0 in this regime, see fig. 6I-IIa. For ∆µ > ∆µ c,1 the maximum of P out is no longer achieved in the MC phase, but rather at the dynamical phase transition between the HD and the MC phase, as depicted in figs. 6I-IIb. As a consequence, the optimal force f * OBC is larger than the corresponding force f * 0 for the noninteracting system, which in turn leads to a higher EMP, in analogy to the observation made for the system with PBC-LK. In particular, for ∆µ/T → ∞, the HD-MC phase transition occurs infinitely close to stall, which implies that η * OBC → 1 in this limit. We note, however, that this limit is caused by the symmetry of the phase diagram and is unphysical. We shall see in the following that under more realistic conditions, such as dissipation by futile hydrolysis, as considered in model II, or nonnegligible load dependence of β, the EMP decays to zero as ∆µ approaches infinity.
Besides the two regimes characterized above, the EMP exhibits a third, qualitatively different, regime whenever φ = 0, see figs. 6IIIa-c. When ∆µ is larger than the critical value ∆µ c,2 with ∆µ c,2 > ∆µ c,1 (see App. B for a discussion of the critical value ∆µ c,2 ), the motors operate at maximum power at the discontinuous phase transition between the HD and the LD phases as shown in fig. 6IIIc. Hence, the maximizing force only depends on φ, and is (numerically) larger than f * 0 . Since the stalling force increases with ∆µ and f * OBC is constant, the EMP thus goes to zero for ∆µ → ∞ as illustrated in fig. 4b. Furthermore, the EMP decreases with increasing φ. It is worth noting that the dependence on φ is contrary to the one obtained for PBC+LK, where the EMP, and hence the enhancement in the EMP, increased when increasing φ. Moreover, it follows from eq. (14) that the EMP η * OBC in this regime is independent of θ, while the enhancement in the EMP represented by the ratio η * OBC /η * 0 increases with increasing θ, which is again at variance with the behaviour observed for PBC+LK. Thus, the detachment dynamics at the boundaries and in the bulk, respectively, affects the thermodynamics of molecular machines differently.
In conclusion, for ∆µ > ∆µ c,1 , we observe, as for the case of PBC-LK, that motor-motor interactions lead to an increase in the EMP, as compared to the noninteracting system, due to a change in the characteristic response of the velocity to the external driving.

B. Velocity and density at maximum power
We proceed by discussing the velocity at maximum power (VMP) and the density at maximum power (DMP) depicted in figs. 7 and 8, respectively, for the two cases of PBC-LK and OBC. When the particle exchange occurs in the bulk, fig. 7a shows that the VMP, v * LK , is smaller than the corresponding VMP for the non-  1 − ρ * LK , see eq. (9), and is hence roughly proportional to v * 0 . Moreover, we note that the DMP is rather insensitive to changes in the load dependence parameters θ and φ when these parameters attain realistic values, see fig. 8a. This in turn leads to a negligible dependence of the VMP on θ and φ.
The VMP in the case of OBC is depicted in fig. 7b and exhibits three different regimes in analogy to the EMP, see fig. 4b. In the first regime, characterized by ∆µ < ∆µ c,1 , maximum power is achieved in the MC phase for f * = f * 0 , and the VMP is thus simply v * OBC = v * 0 /2, cf. eq. (12). The corresponding DMP is ρ * OBC = 1/2 as shown in fig. 8b. For ∆µ c,1 < ∆µ < ∆µ c,2 (or ∆µ c,1 < ∆µ for φ = 0), the system operates at maximum power at the MC-HD boundary. The DMP is thus again and is considerably smaller than the velocity v * 0 of noninteracting motors in the maximum power regime. In the third regime defined by ∆µ c,2 < ∆µ, the power output is optimal at the LD-HD phase transition. Here, the output power and the motor velocity as functions of f are discontinuous due to a discontinuity in the density. However, keeping in mind that the maximum power regime corresponds to a specific value of the force, and a molecular motor might operate at a slightly different force, one may relax the requirement down to, say, 90% of maximum power. This places the system in the low density regime, and in figs. 7b and 8b we thus report for ∆µ c,2 < ∆µ the low-density values for the VMP, v * , respectively. The transition to this regime is accompanied by a drastic decrease in the DMP. It is also worth noting that the VMP in this region increases with decreasing values of θ. Furthermore, it is interesting to observe that the EMP, as well as the enhancement in the EMP, increases significantly when decreasing φ, while the VMP decreases notably when φ is decreased. Thus, in artificial many-motor systems, altering the force dependence of the detachment dynamics can serve as a control mechanism for the trade-off between efficient and fast transport.

IV. MODEL II
We will now generalize model I to include several biochemical motor states. In order to achieve this goal we consider two different network models for kinesin's mechanochemical cycles and incorporate them into the standard ASEP. Model IIa is a six-state model introduced in [18] to describe a single kinesin motor operating under a constant external load force, while model IIb represents an extension of the previous model [18] to seven states [19]. As shown in [19], the two single motor models provide similar predictions for the thermodynamics of kinesin for ∆µ 20T , which includes the biologically relevant range for the free energy, ∆µ/T ∼ 20 − 25 [35]. Since model IIa, due to the presence of only one mechanical stepping transition, admits an analytic solution for OBC by using the maximal current principle, we will primarily use model IIa to investigate the effect of interactions on the open-boundary dynamics. When the Langmuir kinetics dominates, however, a detailed description of the detachment dynamics is necessary, which requires the presence of the seventh state, and model IIb is therefore the appropriate to consider in this case. Furthermore, for both model IIa and IIb we explore the possibility that steric interactions affect the internal conformational motor dynamics leading to mechanical stepping, as well as the stepping transitions themselves. In the following we refer to this type of exclusion rule as strong exclusion.
Below we describe the two models IIa and IIb in detail and discuss the mean-field solutions of the models.

A. Model IIa
This model is based on the six-state network of biochemical states introduced in [18] and depicted in fig. 9a. Kinesin is a two-headed motor, and each of the two motor heads goes through an ATP hydrolysis cycle. Since the rate of phosphate release from a motor head is limited by the preceding hydrolysis reaction, i.e. P i is released immediately after ATP hydrolysis, these two subsequent transitions can be combined into a single one. Each head can therefore either be empty, contain bound ATP or bound ADP. The discrete state space of kinesin is thus composed of 3 2 = 9 states that differ in chemical composition. However, in Ref. [18] Liepelt and Lipowsky argue that the chemical processes in the heads are coordinated, and hence that only states where the two heads have different chemical composition should be connected in the network. They thus arrive at the reduced state space consisting of six chemical states. The motor hydrolyzes an ATP molecule and performs a forward step through the F + = |12341 dicycle, while the B + = |45614 dicycle represents ATP hydrolysis leading to backward step-ping, see fig. 9a. Furthermore, the network includes the futile hydrolysis dicycle D + = |1234561 , where two ATP molecules are consumed while no stepping, and hence output work, is performed by the motor. In order to satisfy the constraints dictated by the local detailed balance requirement, the network model must contain the reverse dicycles F − = |14321 , B − = |41654 and D − = |1654321 , in which ATP synthesis takes place.
The mechanical steps where the two heads switch positions on the track correspond to the transitions between states 1 and 4 as shown in fig. 9. The standard ASEP exclusion rule modifies these mechanical stepping transitions, see fig. 9b, while all the other transitions representing internal conformational changes remain unaffected. Furthermore, as for the standard ASEP we allow the motors to bind to the microtubule at the left end with rate α and detach at the right end with rate β. In the mean-field approximation the evolution of the system in the bulk is described by the master equations, with the obvious notation k + 1 = 1 if k = 6 and k − 1 = 6 if k = 1. Here, ω ij is the transition rate for going from state i to state j. The density of motors in the chemical state i at lattice site j is denoted by ρ i j , while ρ j = 6 i=1 ρ i j is the total motor density at site j [36]. Note that the equations for ρ 1 j and ρ 4 j are modified as compared to the non-interacting system due to mutual exclusion. Open boundary conditions entail that eqs. (15) should be modified for the first and the last lattice sites, analogously to the case of model I discussed in section II C. However, we will show in the following that only the total incoming (outgoing) probability current at the left (right) boundary enters the mean-field solution. Hence, the detailed dynamics of the motors at the boundary sites is unimportant as long as the flux of motors entering or leaving the microtubule remains unchanged. We solve the model introduced above under steadystate conditions and in the thermodynamic limit using a mean-field approach known as the maximum current principle (MCP) [26]. According to the MCP, the boundaries of the lattice are substituted by reservoirs of particles, and the dynamics between the reservoirs and the lattice is assumed to be identical to that in the bulk. The (constant) densities of the left and the right reservoirs are denoted by ρ l and ρ r , respectively. The mechanical probability current through the open system in the thermodynamic limit, J m , is then predicted by the MCP to be [26] where J m PBC (ρ) is the steady-state mechanical probability current through a corresponding homogeneous system with periodic boundary conditions (PBC) and density ρ. Furthermore, the MCP implicitly states that the bulk density of the system takes on the value that maximizes (or minimizes) the mechanical current through the lattice as specified by eq. (16). Hence, we proceed by considering the problem of motor traffic on a periodic lattice.

Periodic boundaries
For the periodic system all the lattice sites are equivalent, and the densities thus become independent of the position on the track, ρ i j ≡ ρ i . In this case, the governing equations of motion, eq. (15), admit an analytical solution for ρ i in the steady-state limit, i.e. forρ i = 0. Since the master equations are linearly dependent, the homogeneous density ρ = 6 i=1 ρ i serves as an independent parameter under PBC. Each mechanochemical cycle C then carries a unique probability current, where C is either the forward stepping cycle F , the backward stepping cycle B or the futile hydrolysis cycle D.
The mechanical probability current on the lattice can hence be expressed as while the chemical current accounting for fuel consumption becomes where the second equalities follow from one specific representation of J(F ) and J(B), which is, however, not unique due to current conservation. Figure 10a illustrates J m,c PBC (ρ) for ∆µ = 20T and for different values of the load force. We note that J m PBC is not symmetric around the value ρ = 1/2, as it is in the case of the standard ASEP [16], since the particle-hole symmetry is broken by the presence of internal conformational states, as also pointed out in [12,13].
The motor velocity, v PBC , and the ATP hydrolysis rate, r PBC , can be obtained from the probability currents as  [18] for the Carter and Cross experiment [37]. and are plotted in fig. 10b. Equations (18)- (20) underline that kinesin is a loosely coupled motor since the average number of steps per unit time, v PBC /a, is in general smaller than the hydrolysis rate r PBC as depicted in the inset of fig. 10b. The motor velocity vanishes, i.e. v PBC = 0, under stalling conditions characterized by J(F ) = J(B), see eq. (18). For fixed chemical input ∆µ this condition determines the stall force f s , which is identical to the one obtained for non-interacting motors. Futile hydrolysis implies, however, that r is non-vanishing at stall, and the coupling between the chemistry and the mechanical motion is low in this limit, see inset of fig. 10b.

Open boundaries: Maximal current principle
We now return to the original problem with open boundary conditions (OBC). The MCP prediction, eq. (16), can be thought of as a variational statement for the bulk density that together with the expression for J m PBC , eq. (18), leads to three qualitatively different phases depending on the relative values of the reservoir densities ρ l and ρ r , as depicted in fig. 11. Figure 11a illustrates the so-called low-density (LD) phase that occurs when ρ l is smaller than the densityρ = argmax ρ J m PBC (ρ) that maximizes J m PBC , and when ρ r is smaller than ρ max r , where ρ max r fulfills J m PBC (ρ l ) = J m PBC (ρ max r ). The bulk density attains the value ρ LD = ρ l in this case, since the current takes the value J m = J m PBC (ρ l ), as predicted by eq. (16). On the other hand, when ρ r is larger thanρ, and , the system is in the high-density (HD) phase, and the bulk density is ρ HD = ρ r , since the current takes the value J m = J m PBC (ρ r ) according to eq. (16), see fig. 11b. Finally, for relatively low right reservoir densities, ρ r <ρ, and relatively high left reservoir densities, ρ l >ρ, the system bulk dynamics becomes independent of the boundary conditions, and the bulk density assumes the maximal current value, ρ MC =ρ.
The reservoir densities can now be calculated in terms of the transition rates by employing current conservation. Equating the bulk current in the LD phase with the current at the left boundary yields which can be solved analytically for ρ l in terms of the attachment rate α and the single-motor jumping rates ω ij . Similarly, current conservation in the HD phase determines the right reservoir density ρ r as a function of ω ij and the detachment rate β through where we assume that track detachment in model IIa can proceed from any mechanochemical state. The procedure sketched in figure 11 together with eqs. (21)-(22) allow us to compute analytically the complete mean-field phase diagram for any choice of the parameter set. The phase diagram projected onto the attachment rate α and the detachment rate β is shown in fig. 12 for ∆µ = 20T and for three different values of the load force. The transition between the LD and the MC phase takes place at α = α c , where α c fulfills ρ l (α c ) =ρ, cf. fig. 11a and 11c. The phase boundary between the HD and the MC phases occurs at β = β c , where β c is given by ρ r (β c ) =ρ, see fig. 11b-c. Finally, the LD-HD boundary is characterized by J m PBC (ρ l (α)) = J m PBC (ρ r (β)) as indicated in fig. 11a-b. We note that the MC phase grows with increased load in agreement with the phase diagram obtained for model I, see fig. 3. Furthermore, the phase diagram is generally asymmetric in α and β due to the lack of particle-hole symmetry, and the LD-HD boundary is highly non-linear for small to intermediate load forces. The non-linearity arises from the interplay of the ASEP exclusion rule, imposed on the mechanical transitions, with the internal chemical transitions that are not affected by steric interactions. The velocity and the hydrolysis rate are simply determined from the phase diagram as v OBC = v PBC (ρ i ) and r OBC = r PBC (ρ i ), respectively, where i is either LD, HD or MC. Finally, it is interesting to note that by using the MCP together with the same arguments discussed in this section, one finds that for any exclusion process with internal dynamics whose current J m PBC in PBC has a single maximum as the one depicted in fig. 11, the phase diagram is qualitatively identical to the one presented in fig.  12 with three different phases.

Strong exclusion
In the following we consider the consequences of strong exclusion for many-motor kinetics. In this scenario the ASEP exclusion rule affects all the transitions between internal states, and not only the stepping transitions where spatial displacement occurs, as illustrated in fig.  9c. The mean-field steady-state master equations for PBC thus become  Figure 13. Model IIb: a) The extended seven-state network model for single kinesin with two mechanical stepping transitions that is presented in [19] and used to model the experimental data of Carter and Cross [37]. b) The standard ASEP exclusion rule prevents the motor from performing the two stepping transitions between states 1 and 4 and 3 and 6, respectively, if the neighbouring site is occupied.
with k = 1, . . . , 6, and can easily be solved to obtain for the mechanical probability current. Here, v 0 is the single motor velocity of non-interacting motors as obtained from the six-state model [18]. Similarly, the chemical current can be written as where r 0 is the ATP hydrolysis rate of non-interacting motors. By comparing eq. (24) with the corresponding eq. (6) for model I it is tempting to believe that strong exclusion gives rise to qualitatively similar thermodynamics of interacting motors as the tightly coupled model I. However, there is an important difference between the two cases. For model I, the coupling ratio v/ra is by construction equal to the corresponding ratio v 0 /r 0 a = 1 for non-interacting motors. While eqs. (24)- (25) imply v/ar = v 0 /ar 0 in the presence of strong exclusion, a single motor, as described by model IIa, is in general not tightly coupled, i.e. v 0 /r 0 a < 1. Hence, when chemical transitions are directly affected by steric exclusion, only futile hydrolysis events due to traffic jams are prevented, and the coupling ratio becomes independent of motor density. It is also worth noting that strong exclusion restores the particle-hole symmetry in model IIa. As a consequence, the density in the MC phase is ρ = 1/2, and the phase diagram becomes symmetric, with α c = β c , as already discussed for model I.

B. Model IIb (PBC-LK)
We now turn to the seven-state model for single kinesin molecules operating under external load that is introduced by Liepelt and Lipowsky in [19]. This model extends the previously considered six-state model [18] by introducing an additional motor state, state 7, and an additional stepping transition, see fig. 13a. Since both motor heads are loosely bound to the filament in the newly introduced state, the motor is most likely to detach from the track from state 7. Hence, the seven-state model admits a precise description of the detachment pathway, which is necessary for considering Langmuir-type kinetics. The master equations for model IIb with Langmuir kinetics are discussed in App. C. In order to solve these equations under steady state conditions, one has to resort to numerical techniques. Furthermore, we consider the effect of strong exclusion (SE) on the dynamics within this model. In analogy to strong exclusion dynamics in model IIa, cf. eqs. (24)- (25), the mechanical and chemical probability currents are given by where the Langmuir density ρ LK is obtained as Here, α is the rate with which free motors bind to an empty lattice site, while β is the bulk detachment rate characterizing the unbinding of motors from an occupied lattice site. The quantity ρ 7 0 is the probability to find a non-interacting motor in the detachment state 7, while v 0 and r 0 denote, respectively, the motor velocity and hydrolysis rate of non-interacting motors calculated from the seven-state model. The motor velocity and hydrolysis rate of single motors in the presence of interactions are thus linear functions of the motor density.
In fig. 14a we report the probability currents J m LK and J c LK as a function of the Langmuir density ρ LK obtained by varying the ratio α/β 0 of the binding rate to the (zeroforce) unbinding rate. The curves correspond to the free energy ∆µ = 20T and to two different values of the external load force, f = 0 and f = −3 pN. Furthermore, we plot the numeric results obtained for the mechanical exclusion (ME) rule as dictated by the standard ASEP as well as the curves for the strong exclusion (SE) rule, eq. (26). We note that the motor stall force for ME is a function of the rates α and β, and as such varies with ρ LK . The density ρ LK , which is obtained numerically for ME and is given by eq. (27) for SE, is in general a fast increasing function of α/β 0 as shown in fig. 14b. Interestingly, ρ LK is essentially the same for standard and strong exclusion, even though the probability currents are largely affected when adding exclusion on chemical transitions. Finally, fig. 14c shows the resulting motor velocity v LK = J m LK a/ρ LK and hydrolysis Langmuir density ρLK as a function of the ratio α/β0. c) The velocity vLK and the hydrolysis rate rLK as a function of ρLK obtained from a) and b). Parameter values are those obtained in [19] for the Carter and Cross experiment [37]. rate r LK = J c LK /ρ LK . For ME (curves with empty symbols in the figure), the velocity is essentially constant and equal to the velocity of non-interacting motors, i.e. v LK ≃ v 0 , for a wide range of densities, 0 < ρ LK 0.9, which is in agreement with the experimental observations of Seitz and Surrey on the kinetics of conventional kinesin (kinesin-1) motor traffic at zero force [38]. It is also worth noting that our result for the velocity in the presence of Langmuir kinetics is at variance with that of Klumpp, Chai and Lipowsky [12]: by using a simplified two state model, they conclude that the internal cycle has a small effect on the linear decrease of velocity with increasing density, such a decrease being expected for systems with no internal dynamics, cf. eq. (9). However, it is interesting to note that the linearly decreasing density-velocity relationship, eq. (28), is recovered for SE, i.e. when the exclusion rule affects both chemical and mechanical transitions. As opposed to Seitz and Surrey [38], in a recent experiment Leduc et al. [39] observe a linear relationship between motor density and velocity for kinesin-8 motors moving on crowded filaments. Hence, our modelling suggests that the contradictory experimental results for v(ρ) could be attributed to different extents of strong exclusion that different motor types might exhibit, in the different experiments.  [18] for the Carter and Cross experiment [37], second row (b): same parameters as in (I) but with the chemical load factors χ1 = 0.3 and χ2 = 0.4 taken from the seven-state model [19] (see text). The detachment rate force sensitivity is φ = 0 (I) and φ = 0.1 (II), respectively. The value of the zero-force unbinding rate is β0 = 3 s −1 . The biologically relevant regime ∆µ/T ∼ 20 − 25 [27] is marked by vertical dashed lines.

V. MODEL II: MAXIMUM POWER REGIME
In the following we investigate the operation of loosely coupled molecular machines in the maximum power regime. Model IIa with six internal states is used to investigate systems where the particle exchange predominantly occurs at the lattice boundaries, while the sevenstate model IIb allows us to consider motors interacting with a bulk reservoir.

A. Efficiency at maximum power
The optimal force f * maximizing the output power for model II is calculated as described in sec. III for model I. Since the mechanical and chemical currents are not tightly coupled, the EMP for a fixed value of ∆µ is then obtained as where r * = r(f * ) is the hydrolysis rate at maximum power. We find that the EMP is only susceptible to variations of a few parameters, namely the binding rate α, the unbinding force dependence parameter φ, and the chemical load factors χ 1 and χ 2 that describe the load dependence of the transition rates associated to chemical transitions [40]. There is in general no consensus in literature on the biologically applicable values for these parameters since they are highly model dependent and/or experiment dependent as will be discussed below. In fig.  15 we plot the EMP for a few combinations of parameter values for α, φ and χ i . The values are chosen such as to show the typical qualitative behaviour for the EMP when the parameters are varied within a range compatible with experimental results reported in literature. We therefore believe that the results for the EMP shown in fig. 15 are representative for the rather broad range of biologically relevant parameter values. The solid lines in fig. 15(Ia) show the EMP predicted for model IIa with mechanical exclusion calculated using the parameter values obtained in [18] by fitting the Carter and Cross' experimental data [37] and for force independent detachment dynamics, i.e. φ = 0. For both α = 5 s −1 (green crosses) and α = 50 s −1 (magenta squares) the single-motor EMP exhibits an enhancement as compared to noninteracting motor systems under the same conditions when the (dimensionless) free energy ∆µ/T lies in the range 18 − 22. Thus, the region of EMP enhancement partly overlaps with the biologically compatible range ∆µ/T ∼ 20 − 25 [27]. As for model I, the EMP increase is induced by the presence of phase transitions in the collective dynamics of interacting motors as depicted in fig. 16 for a specific set of parameters. Here, we plot the output power P out as a function of the load force f for various values of ∆µ. For small values of ∆µ the maximum power output is reached in the MC phase, and the EMP η * is essentially equal to the EMP η * 0 of noninteracting motors, see fig. 16a. For 15.7 ∆µ/T 18.6 the motors operate at maximum power in the LD phase, see fig. 16b, and η * < η * 0 since the velocity v in the LD phase decreases more rapidly with f than the velocity v 0 of non-interacting motors. When ∆µ/T increases beyond the value 18.6, the maximum of P out lies at the boundary between the HD and LD phase as illustrated in fig. 16c. The optimal force f * is thus greater than f * 0 , and the EMP η * exceeds the noninteracting value η * 0 as shown in fig. 17a. At ∆µ/T ∼ 22 the force f * becomes so large that η * = η * 0 . For even larger values of ∆µ, the maximizing force f * approaches the stall force, and η * → 0 due to dissipation. Furthermore, for ∆µ/T 22 the maximum power occurs at the MC-HD boundary, see fig. 16d.
For a larger value of the binding rate, α = 50 s −1 , the EMP also exhibits an enhancement, which is, however, smaller than for α = 5 s −1 . We also note that the maximum of the output power lies at the MC-HD boundary in the region of enhancement as shown in 17b. We stress, nevertheless, that the physical mechanism leading to the EMP increase is identical for the two situations, namely that motor-motor interactions affect the characteristic response of the velocity to external mechanical driving. For α = 5 s −1 the EMP boost disappears rapidly with increasing values of φ, see fig. 15(IIa), because the maximum power is obtained in the LD regime as pictured in fig. 17c. The EMP curves for α = 50 s −1 are more robust to changes in φ, since the motors operate deeper in the HD phase for f = 0, see fig. 12. This effect is enhanced by the fact that the HD-LD boundary is highly non-linear as discussed in sec. IV A 2.
Together with the EMP results obtained for the ME rule and discussed above, we plot in fig. 15 the EMP predicted for the SE case (dot-dashed lines). We observe that the EMP increase is present and pronounced for both values of α when φ = 0, fig. 15(Ia). However, even for SE, the effect disappears for α = 5 s −1 when increasing φ for the same reasons as for ME. On the other hand, for α = 50 s −1 the EMP increase is robust towards variations in φ, because the maximum of the output power lies at the MC-HD transition as for ME. Hence, interestingly, for small α the dependence of the EMP on φ is determined by the intrinsic dissipation present in the single-motor model rather than by traffic jam related futile hydrolysis.
In figs. 15(Ib)-(IIb) we report the EMP where the chemical parameters take on the values obtained by fitting the seven-state model [19] to the experimental data of Carter and Cross, while all the other model parameters are kept constant. For this parameter set we observe an approximately twofold increase in the EMP as compared to the non-interacting case within a wide range of free energies, ∆µ/T ∼ 18 − 40, and for all combinations of parameter values and exclusion types except one. However, the absence of EMP enhancement for the case α = 5 s −1 , φ = 0.1, with mechanical exclusion seems to be due to the fact that we have considered the simplified six-state model: indeed, by considering the same parameter values for the seven-state model, we find an increase in the EMP, see App. D. It is worth to note that, for negligible values of φ, the EMP for interacting molecular machines is essentially independent of the binding rate and type of exclusion considered.
A discussion of the parameter values used in the present study is now in order.
By comparing the values for the chemical load parameters for the six-state model, χ 1 = 0.15 and χ 2 = 0.25, with the corresponding values for the seven-state model, χ 1 = 0.3 and χ 2 = 0.4, we conclude that they depend strongly on the underlying model used to fit the experimental data. Furthermore, the parameters χ i assume different values when the six-state model is used to fit two different sets of experimental data [18]. Thus, the values for χ i are also sensitive to experimental conditions. Lastly, it is worth noting that the parameters for model II have been obtained in refs. [18,19] by fitting experimental data for single motor movement, and thus it is also possible that the parameters χ i assume different values in the presence of multiple motors on the track. This conjecture is based on the fact that the description of the dynamics of a single motor as a stochastic hopping process between discrete mechanochemical states is merely an approximation to the microscopic picture where the motor movement is represented as Brownian diffusion on   fig. 15 (row I). We consider the ME rule here.
a continuous free energy landscape [20,41,42]. From this point of view, the states of the discrete model correspond to the potential minima of the energy landscape, and the transition rates of the discrete model are obtained from the minimal energy paths connecting the potential minima. Since motor-motor interactions alter the energy landscape, it is thus likely that the parameters representing the coupling of the external force to the mechanochemistry will change as a consequence. In ref. [20], for example, we address the question of how the load factor θ introduced in model I changes with the external force, in a system with continuous phase space.
To our best knowledge, there are very few precise measurements of φ available in literature. The value φ = 0.1 for the force dependence of the unbinding rate is obtained in ref. [18] by fitting a simplified seven-state model with only one mechanical transition to the Visscher et al. experiment [43]. Since the parameters χ i describing the effect of the external load force on chemical transitions are both model dependent and experiment dependent, it is hence likely that φ would be so as well and thus attain a value different from 0.1 when the full seven-state model is fitted to the experimental data of Carter and Cross.
The value of the binding rate α is widely assumed to depend linearly on the (local) motor concentration [10,12]. Since intracellular kinesin concentrations vary over several orders of magnitude [27], there is therefore reason to believe that α attains a wide range of values in biological systems. For example, the values used for α in the present paper lie within the biological range 0.8 − 80s −1 obtained from the binding constant 8.27 × 10 7 m −1 s −1 [44] and intracellular kinesin concentrations (10 − 1000) · 10 −9 m [10].
In summary, considering the results for the EMP presented in this section and the above discussion on parameter uncertainty, we believe that the EMP enhancement observed takes place for biological systems under a wide range of biological conditions. Finally, we consider the EMP for model IIb which accounts for the behaviour of the efficiency when particle exchange predominantly occurs in the bulk. As for model IIa, we investigate the dependence of the EMP on α, φ and exclusion type and plot typical curves in fig. 18. We note, however, that, as opposed to the system with OBC, the EMP calculated within model IIb changes in a continuous manner when the parameters are varied since no dynamical phase transitions are present in the case of Langmuir kinetics. For φ = 0, see fig. 18a, the EMP in the presence of motor-motor interactions is smaller or approximately equal to the EMP of non-interacting motors for a wide range of α-values and for both mechanical and strong exclusion. The EMP is generally decreasing with increasing α. The EMP for φ = 0.1 is shown in fig. 18b. In the presence of exclusion on both mechanical and chemical transitions (SE), the EMP exhibits a considerable increase as compared to η * 0 in the region ∆µ/T ∼ 5 − 22. Furthermore, the EMP increases with α for these values of the free energy as opposed to the observed general tendency. Interestingly, as opposed to model IIa, the motor density variations induced by imposing steric exclusion on mechanical stepping only are not sufficient to observe EMP enhancement in the case of LK.

B. Velocity and density at maximum power
We proceed by discussing the VMP and DMP of interacting kinesin motors as obtained within model II. Fig.  19 shows the VMP v * (a) and DMP ρ * (b) for the pa- . Parameter values are as obtained in [19], and β0 = 3 s −1 . rameters used to calculate the EMP in fig. 15b. The corresponding results for the second set of parameters used in fig. 15a are similar and therefore not presented here. First of all, we note that, contrary to the behaviour predicted by the tightly-coupled model I, see fig. 7b, v * at most differs from the VMP of non-interacting motors v * 0 by one order of magnitude. This effect is due to the saturation of the velocity at biological values of the chemical input as discussed in App. A. Secondly, the VMP is generally lowered when the exclusion rule affects both the mechanical and chemical transitions. Furthermore, in analogy to model I, the VMP and EMP exhibit  fig. 18b. a trade-off behaviour; the EMP decreases for increasing values of φ, while the VMP increases. As for model I, the DMP never exceeds 0.5 for SE, since maximum power is never achieved in the high-density regime. For the ME rule, however, the DMP reaches ρ * ∼ 0.85 for α = 50 s −1 since the density in the MC phase is shifted to higher values, see fig. 11. Moreover, the DMP is sensitive to the value of α.
Finally, in fig. 20 we report the VMP (a) and the DMP (b) for model IIb. Interestingly, for ME, the VMP can exceed the corresponding VMP v * 0 in non-interacting systems for high values of ∆µ. However, the high VMP comes at the expense of the EMP η * as can be seen in fig. 18b. For strong mechanical and chemical exclusion, where an EMP enhancement is observed for a wide range of ∆µ-values, the velocity is reduced in analogy to model IIa. In analogy to model I, the DMP does not vary significantly with ∆µ. Contrary to model IIa, the DMP is essentially unaffected by chemical exclusion and is above 0.5 for all values of ∆µ.

VI. DISCUSSION AND CONCLUSIONS
In the present paper we have considered two different models for kinesin motors interacting through steric exclusion on a filamentous track. Furthermore, we have investigated the effect of different binding and unbinding kinetics on the overall dynamics of the interacting system as well as on the operation in the maximum power regime. In the following we discuss the relevance of our findings to intracellular transport and in the context of state-of-the-art multiple-motor experiments.
In the present work we have for simplicity been concerned with single-motor-single-cargo systems. There is indeed experimental evidence that justifies this assumption. It has been found in several in vivo studies on different motor types (see [45] and references therein) as well as in a recent in vitro experiment on kinesin-1 [46], that only a single motor is actively engaged at any one time in carrying its cargo. Furthermore, another experimental study [47] performed with multiple kinesin-1 motors concludes that cooperative effects for this motor are small, albeit they increase with efficiency. Nonetheless, it is widely believed, primarily based on results from in vitro experiments, that multiple motors are involved in cargo transport in the cell giving rise to non-steric coupling between the motors. A recent state-of-the-art in vivo study of intracellular transport of lipid droplets in Drosophila embryos [30] corroborates this conclusion. The typical run length of a single kinesin motor is known to be around 1 µm as reported in, e.g., [12,18,38]. Since microtubule lengths are at least one order of magnitude larger than the single-motor run length [27], one would therefore expect the Langmuir dynamics to be non-negligible when the cargo is transported by a single motor. However, several experimental as well as theoretical works report that the cargo run length increases by up to several orders of magnitude when multiple motors carry a single cargo [8,[48][49][50]. Even under the assumption that only one motor is active at any time, the presence of other motors on the cargo would likely lead to an enhancement in the run length. Hence, under physiological conditions one would expect the cargo run lengths to be comparable to the microtubule lengths. Furthermore, the density of the cargo is not homogeneous across the cell, so one would rather expect that the cargo concentration gradient led to the cargo binding primarily occurring at the filament end closer to the production site. Based on the above arguments, we thus believe that intracellular traffic is best described by a model with open boundary conditions.
In that context we find that the EMP of kinesin motors is increased, as compared to the non-interacting case, due to a change in the characteristic force-velocity relation, v(f ), as discussed in sec. V A. Remarkably, for the two kinesin models studied, it occurs in the biologically relevant portion of the parameter space and for different types of exclusion rules. Still, the results have been obtained for the single-motor-single-cargo dynamics, and future developments of the present research would be concerned with intracellular transport involving several motors on a cargo. However, we believe that our conclusions can be carried over to the many-motor-single-cargo system, as discussed below. When a cargo is transported by multiple motors, a question arises regarding the sharing of the load force amongst the motors, as also discussed in [51]. Leidel et al. [30] arrive at the conclusion that the load force is shared equally by the motors, and that only one type of motors is involved in pulling the droplets at any one time, i.e. there is no tug-of-war mechanism as has been suggested by several authors, see, e.g., ref. [9] and references therein. It is certainly true that the presence of multiple motors on a cargo would modify v(f ), as also discussed in [8,11,[51][52][53]. However, based on the assumption of homogeneous load-sharing and identical motors discussed above, we believe that v(f ) would only change quantitatively, and not qualitatively, see, e.g., [8], showing thus a monotonic behaviour as for the case we consider here. Hence, our single-motor-single-cargo model system can be thought of as a renormalized version of the multiplemotor-single-cargo system. One may therefore expect the multiple-motor-single-cargo systems to undergo a dynamic phase transition from a high-density phase to a phase with a lower density, resulting in an enhancement in the EMP, similar to the one we find for our models.
Despite the recent tremendous progress achieved within experimental tracking and manipulating techniques, it is as yet not possible to study under controlled conditions the dynamics of multiple motors operating under externally applied load forces. One class of stateof-the-art motor traffic experiments deals with traffic of single motors under zero-load conditions [38,39]. In this case, the density of the free unbound motors surrounding the filament is expected to be homogeneous across the system, and the run length is short compared to the filament length, as discussed above. Under such conditions we therefore expect the Langmuir dynamics to dominate. Furthermore, it has been proposed that the use of particular tracking techniques such as quantum-dot labelling can influence motor-motor interactions [54]. Hence, the extent of exclusion can well depend on the experimental setup, which the different density-velocity relations, v(ρ), could be a sign of, see sec. IV B. We expect that future motor traffic experiments would provide a better estimation for the model parameters which have only been extracted from experimental data on single motors so far.
Finally, we note that the enhancement in EMP is present for each variant of model kinesin we considered, for a wide choice of parameters: this shows that this effect does not originate from the details of the microscopic dynamics, but is rather a result of the cooperative nature of the system. We therefore believe that the interaction-caused enhancement of the efficiency due to a change in the force-velocity response is a generic feature for nanomachines optimized for transport in the maximum power regime.  [18] for the Carter and Cross experiment [37] (magenta, squares) and for the Visscher et al. experiment [43] (green, crosses), respectively. The corresponding values of the microscopic rate ω0 that yield the best fit of the velocity are ω0 = 1.00 × 10 −7 s −1 and ω0 = 1.33 × 10 −7 s −1 , respectively, for the two experiments.  [37], we obtain ω 0 = 1.00 × 10 −7 s −1 , while using the parameters obtained by Visscher et al.'s data [43] leads to ω 0 = 1.33 × 10 −7 s −1 .
It is worth noting, however, that energy dissipation due to loose coupling results in saturation of the velocity at ∆µ 20k B T , see fig. 21, which is not present in the tightly coupled model I.
For the attachment rate α and the bare detachment rate β 0 we take 5 s −1 [8] and 3 s −1 [18], respectively. For the force dependence of the detachment rate Liepelt and Lipowsky find φ = 0.1 [18]. However, since the EMP strongly depends on φ, cf. fig. 4, we consider different values of φ in sec. III. The mechanical load sharing parameter θ is sensitive to the choice of model and to experimental conditions [18,19]. In model I we use the values θ = 0.65 and θ = 0.3 as obtained in [18] for the Carter and Cross and for the Visscher et al. experiment, respectively. place from different internal states, the chemical potential differences involved in the two processes are different. Hence, thermodynamic consistency requires that an additional term due to particle exchange with the reservoir is introduced in the input power. Since experimental values for these chemical potential differences are not available at present, we have chosen to minimize the number of unknown parameters by taking state 7 as the binding state.
Appendix D: Model IIb. OBC The mean-field equations for the seven-state model IIb with OBC can be solved numerically with the help of the MCP. The solution is, however, difficult to obtain due to numerical errors arising as a consequence of the exponential dependencies of the transition rates on the free energy ∆µ and load force f . However, we have calcu-lated the EMP with mechanical exclusion for α = 5 s −1 , φ = 0.1 and for a limited range of ∆µ values, see fig.  22. As opposed to what we find for the six-state model, fig. 15(IIb), the resulting EMP exhibits a significant increase, as compared to the non-interacting system, for ∆µ/T ∼ 19 − 23.