Polymer dynamics in time-dependent periodic potentials

Dynamics of a discrete polymer in time-dependent external potentials is studied with the master equation approach. We consider both stochastic and deterministic switching mechanisms for the potential states and give the essential equations for computing the stationary state properties of molecules with internal structure in time-dependent periodic potentials on a lattice. As an example, we consider standard and modified Rubinstein-Duke polymers and calculate their mean drift and effective diffusion coefficient in the two-state non-symmetric flashing potential and symmetric traveling potential. Rich non-linear behavior of these observables is found. By varying the polymer length, we find current inversions caused by the rebound effect that is only present for molecules with internal structure. These results depend strongly on the polymer type. We also notice increased transport coherence for longer polymers.


I. INTRODUCTION
There has been considerable progress in the research of Brownian motors during the last decade (see e.g. [1,2,3]). Starting with the simple pointlike Brownian particles with timedependent driving forces, research has expanded towards the more complex objects such as interacting Brownian particles (e.g. [4,5,6,7,8,9]) and polymers [10,11,12]. In this paper, we study polymer motion with discrete lattice models, which allows us to consider different kinds of microscopic polymer dynamics in detail. Aside from being a purely theoretical branch of study, analysis of simplified discrete non-equilibrium particle models has became an important tool for studying biologically inspired Brownian motor systems (e.g. [13]).
Discrete models have been applied widely to single particle ratchet problems (e.g. [14,15,16,17]). We expand this picture by considering a generalized Rubinstein-Duke model (RD model) [18,19] for polymer motion in discrete time-dependent potentials. An interesting question is, what kind of dynamics lies beyond simple pointlike particles and how one can calculate its properties such as the effective diffusion coefficient and the drift. Although there are plenty of studies concerning the behavior of the RD polymer in zero (or uniform) field (e.g. [20,21]), only recently a ratchet mechanism (tilting ratchet) has been considered in this context [22].
Especially because of the high complexity of Brownian motors with internal structure, most studies of these systems have applied the Monte Carlo method. However, since the ratchet systems are quite sensitive to the values of parameters, and drifts generated by pure ratchet mechanism are usually very small, Monte Carlo simulations tend to be very timeconsuming and inaccurate. In this paper we study these systems with the master equation approach. The results obtained this way are accurate enough to reveal the details of the dynamics.
The purpose of this paper is to give a hands-on example of how one applies the master equation method to systems involving time-dependent periodic potentials and complex molecules by using a modified RD model polymer as a prototype of such molecules. We perform calculations for short linear polymers in a non-equilibrium environment generated by flashing and traveling ratchets. To test the significance of the polymer type, reptating or not, we compare the motion of the RD polymer with the dynamics of a modified version of the RD polymer with less constrained microscopic movement.
The paper is organized as follows: In Section II we expand the modified RD polymer model to periodic time-dependent potentials and give equations for the calculation of the drift and diffusion coefficients, in Section III we present results of the calculations for short polymers and finally in Section IV we give our conclusions and discuss the implications to applications.

II. THE MODEL AND METHODS
The RD model was originally developed to model the random motion of a flexible polymer in a confined medium with static obstacles (e.g. pores in gel) that the polymer must bypass, therefore causing the polymer reptation. By assuming that the network of obstacles can be modeled (on average) by a lattice-like structure, that the correlation length between the polymer segments is smaller than the distance between the obstacles, and that only the polymer heads are able to move into previously unoccupied cells (lattice sites), the problem can be discretized to a simple particle hopping model [18]. Soon after the original model was expanded [19] to be suitable for external potentials (e.g. static field), pure theoretical research of the model started to flourish such as in Refs. [20,23]. Technically the RD model is a spin-1 chain with special kind of nearest neighbor interactions between the particles (reptons). By assuming that the reptons experience random "pushes" by the environment modelled with a continuous time Markov process with exponentially distributed waiting times, we can construct the stochastic generator of the system.
In order to study the effect of the intrinsic transition rules of the polymer in timedependent periodic potentials on long-time dynamics, we will compare the results of the RD model to the results of a non-reptating polymer which allow the breaking of the reptation tube. In this paper we call this extended model the free-motion model (FM model). In Fig. 1 there is an illustration of an example configuration of a six repton polymer with arrows indicating all allowed moves for both RD and FM models (see also Ref. [24]). All repton transitions are between nearest neighbor lattice sites only. Similar extensions have been studied previously in a different context in Refs. [24,25,26].
As an environment for the polymers, we assume a discrete periodic potential V (x) such that V (x + L) = V (x). To make contact with Kramers rate theory (see e.g. review [27]) and the previous work related to discrete ratchets [16], we define the transition rate from state i to j by and choose α = β = 1 to define the time and energy scales along with the lattice constant 1 to define the spatial length scale. We shall next define the time-evolution operators for the RD and FM models (readers not interested in the formal development may skip the rest of this section).
The mathematical model for the polymer, which contains the RD model as a special case but also allows breaking of the reptation tube if wanted, is constructed as follows (see e.g. [28]). Within the most compact, the inner coordinate representation, every bond between reptons can be in three states; up (state A), down (state B) or flat (state 0). In An N-repton polymer has N − 1 bonds. The state corresponding to polymer configuration y is thus given by a 3 N −1 -dimensional state vector |Ψ y .
The non-zero elements of the local creation and annihilation operators defining the dynamics of the bonds are The operators a and b produce changes in the local bond configuration as indicated in Fig. 1.
To extend the model to include a periodic potential V , we must add an additional state.
One repton is chosen as a marker repton that keeps track of the polymer position within the potential. The transition rates of a single repton now depend on the position of the marker repton and all other bonds separating it from the marker. Either of the head reptons is the most convenient choice for the marker repton, hence we choose here the repton labeled 1 (see Fig. 1). The dimension of the marker state is L, so the dimension of the total system of equations becomes L × 3 N −1 .
By denoting where the operator A applies to bond 1 and the marker repton, M applies to bulk reptons and B applies to bond N − 1. The explicit forms of these operators are given in Appendix A.

A. Time-dependent potentials
The time-dependence of the environment can break the detailed balance and may result in a directed drift. We assume that the switching between the distinct environments is independent of the polymer state in the potential i.e. there is no feedback from the polymer.
The switching mechanism between the potentials can be either stochastic or deterministic.  [29]).
First assume Markovian switching between the potentials. We must include an additional state that keeps track of the current potential where |Ψ s is the state vector of the potential with dimension S i.e. the number of different potentials. Since there is no feedback mechanism, adding this new state is straightforward.
The non-zero state and transition operator elements for the potential state are where 1 ≤ s ≤ S. By defining the operatorĥ like this, we consider only cyclic transitions between the potentials (i.e. 1 → 2 → · · · → S → 1 → . . . ) to preserve the analogy with the deterministically switching potentials. The stochastic generator becomes With deterministic switching, the stochastic generator is given by where H s is the operator of the type (1) in the potential s and T = S i=1 T i is the timeperiod and H(t + T ) = H(t). In this case there exists a T -periodic stationary solution.
Once H is given, the time-evolution of the system is governed by the master equation is the probability vector. Since H s 's are generally nonsymmetric, q(t) usually has an oscillating behavior.

B. Drift and diffusion
We are interested in the drift and diffusion of the center of mass of the polymer. The velocity and the diffusion coefficient can be defined as where x CM is the center of mass of the polymer. Here v and D eff could also be defined for single reptons instead of the center of mass and this local approach naturally leads to the same longtime values.
From the previous we define the Peclet number where we choose the length scale ℓ = 1. Since our polymer is simply composed of several neighbor-hopping random walkers, we can generalize the formalism of Ref. [17] (which generalizes the ideas of Ref. [29]). First define where p n,l,y is the probability to find the marker-repton in the position l + nL with the polymer inner configuration y ′ and the re-defined state y includes both the marker-repton position (l) and the inner configuration (y ′ ) within the L-periodic potential. Assume that the stochastic generator H of the total system is defined by the rates Γ i,j from state i to j.
It can be shown by using the definitions above, by taking the time-derivatives and using the master equation (see e.g. [30] for a similar calculation) that with arrows indicating the direction (right or left) of those repton transitions that lead from the state y to states i, neglecting all the rest. Since this expression is for the center of mass, Similarly we get The evolution equations for q y (t) and s y (t) can be found by differentiating (2) and (3) in time and using the master equation once more. We arrive at Note that all transitions are assumed to be between nearest neighbor lattice sites only.
Otherwise transitions of certain length should be collected in their own sets according to their hopping distances (≤ L), which would appear as coefficients of additional sum-terms in the equations. In the matrix form where H sign has the structure This operator is easily built while building the stochastic generator itself. Since v(t) in Eq. (5) is governed by Eq. (4), these systems must be solved simultaneously. See also Refs. [31,32] where similar approach has been applied to find the drift and the effective diffusion coefficient for complex molecules.

Time-independent stationary states
When H is time-independent, we can take the limit t → ∞ and define the steady-state parameters as By using these we get well-defined stationary values for the velocity and the effective diffusion coefficient. Now Q y 's and S y 's are found by solving the equations So far equations like these have been solved exactly only for a single particle on a periodic lattice. The first solution was given in Ref. [29] for the nearest neighbor hopping particle with arbitrary transition rates. This has been later extended e.g. for parallel one-dimensional lattices in Ref. [30]. However, for more complex systems (like RD polymers), solutions cannot be found by exact methods and numerics must be applied. The structure of H also raises some issues. Since the determinant of H is always zero, mathematically there is no unique solution for the non-homogeneous linear set of equations in (7). This can be easily seen by using the fact that an ergodic stochastic system always has a non-trivial stationary state, Eq. (6) is a generalization of the result derived in Ref. [20]. This can be seen by considering the case L = 1 without external potentials (i.e. y's are simply inner configurations, v = 0 and Q := Q y = 3 −N +1 ), so that with a y := −2S y /hQ we have (for lattice constant 1) Here we used the fact that, in this case, every state y has a weight 1/3 N −1 and r y /l y 's can be interpreted as the "number of arrows" for right/left transitions out from the state y. In Ref. [20], symmetry properties (reflections) of polymer configurations were used to find a unique solution for a y 's, but this is not possible when external potentials are present and the problem is non-symmetric. However, numerical linear algebra tools can be used to find the solution.

Time-dependent stationary states
When H is time-dependent, we must integrate equations (4) and (5) In practice, these are calculated by integrating in time long enough so that results have converged.

C. Fast and slow switching regimes
When the switching times of the potential are close to the characteristic timescales of the system (i.e. relaxation times), the behavior depends heavily on the switching type and lifetimes of the states. However, when the potential changes very rarely or extremely fast, the system becomes independent of the switching type and even of the relative life-times of the states.
First assume that the total mean switching period T → 0 such that T i > 0 for all mean lifetimes of the potentials (1 ≤ i ≤ S). In this case particles experience an effective average potential ("mean-field" [36]) and the transition rates become where Γ k i,j are the transition rates of the stochastic generator of type (1) in the potential k and x k = T k /T 's are weight factors determined by the mean life-times of the potentials.
This leads to a mean-field stochastic generator with dimension L × 3 N −1 . This approach was used in Ref. [33] to solve exactly the single particle dynamics in two arbitrary alternating periodic potentials. Although this mean-field limit is mathematically well defined, from the physical point of view it's artificial since real-world systems have inertia, and changing the potential state takes some finite time (e.g. charge re-distribution to build up an electric field). So the velocity always goes to zero in the fast switching limit. Now assume that T i ≫ τ i for all 1 ≤ i ≤ S where τ i is the longest relaxation time of the system in the potential i. This means that the system always converges close to the stationary state in the current potential before the potential is switched to the next one. By the model assumptions, drift is always zero at the stationary state in all potentials. Let d j|i denote the mean travel distance of the molecule center of mass within the potential j using the stationary state of the potential i as an initial state and then letting the system fully relax [37]. Summing over all d j|i 's gives the total expected distance within one time-period T , and by assuming cyclic switching of the states, we define The sign of d determines the drift direction in the large T limit and the asymptotic drift That internal molecular states may have strong influence on the dynamics can be already seen in the slow switching regime. Letting the molecule first find its equilibrium in some non-flat potential and then turning the potential off may indeed result in directed motion of the molecule after the switching, due to internal relaxation, whereas a single particle would be immobile in the mean. These rebounds might be dominating and define the sign of d.

A. Choice of the potentials
We have numerically analyzed RD and FM models with the polymer length of N = 1...11 reptons and with two potentials (S = 2) and stochastic switching. All calculations were done with MATLAB. We used the standard Runge-Kutta 4 method to integrate (4) in time and a trapezoid method to calculate the resulting integral in (9). The Arnoldi and BiConjugate gradient stabilized methods were used to solve homogeneous and non-homogeneous systems in (7).
We consider two basic potential types: flashing and traveling ratchets. The first type is the most general non-symmetric potential that has been extensively used in studies of Brownian motors and the latter one is a generic example of asymmetrically placed symmetric potentials and has been recently used with single particle models [17,34]. We consider the simplest case L = 3, which is the smallest possible length that can form both of these potentials with the ratchet effect. See Fig. 2  (N ≫ 11) than we can efficiently handle. We require that the polymer must be able to cover several potential periods when fully extended. We also set V max = 1/2, which we where T on/off are the corresponding mean life-times of the potentials (see Fig. 2). In previous studies (e.g. [15,16]), only symmetric flashing x = 1/2 was considered. This results in zero drift for T → 0, which also happens in all real systems (for all x). However, with x = 1/2, this does not happen for the models considered here. The drift changes its sign as a function is present in Fig. 4 (c).
The behavior of the Peclet number is very clear and similar in every case in Fig. 4: the larger the polymer, the larger the Peclet number. Thus the transport of longer polymers is more coherent than of shorter ones. Similar behavior was found in a continuum model consisting of elastically coupled Brownian particles [4]. By comparing the values of the Peclet number between polymer types, we see no significant differences between the curves.
There is a slight difference for large values of T , where the Peclet number remains larger for FM polymers. This holds with every choice of parameters, excluding the possible current inversion points (e.g. the interval ln (T ) = −1...0 in Fig. 4 (c)).
Next we take a closer look at the asymptotic behavior at T → ∞. In Fig. 5 we have plotted the mean travel distance d defined in the Eq. (9). For N = 1, 2 there are no bulkreptons so the mean travel distances of RD and FM polymers may differ for N ≥ 3 only.
The calculation reveals that for long RD polymers (N > 5, a 'critical length') the rebound effect wins (i.e. d > 0) and the polymer starts traveling backwards while the single particle and FM polymers are traveling to the expected negative direction. The rebound effect is also present in FM polymers, but it is not strong enough to reverse the drift direction. For RD polymers with L > 3 with feasible polymer lengths our Monte Carlo test simulations do not display this kind of an anomalous current inversion, suggesting that it may be related to spatial discretization and that longer-range interactions (e.g. stiffness) between reptons need to be introduced to see such inversions for L > 3.  We also note that a similar effect of multiple current inversions with tightly connected Brownian particles (rods) was reported in Ref. [12]. In that work, however, current inversions were not found for objects able to vary their length (rotating rods) in the ratchet direction, whereas the polymers in our work are able to vary their length between 1...N and still have drift inversion.
The reason for the stronger rebound effect of the RD polymer is caused by the strong tendency to enter (possibly deformed) U-shaped configurations because of the strict reptation rule. After the potential is turned off, this shape unwinds and causes the drift. This also happens with time-dependent fields [22]. Since FM polymers lack the reptation rule, there is not as much variation in their shape as RD polymers have, thus resulting in a weaker rebound effect.

C. Traveling ratchet potential
Let now T = T 1 + T 2 for the mean life-times T 1 and T 2 of the potentials depicted in the right column of Fig. 2 and define the symmetry parameter x = T 1 /T . A Similar drift and diffusion behavior as previously reported in Ref. [17,34] for a single particle is expected.
In T → 0, the traveling potential creates a similar effective rate structure as the non-symmetric flashing ratchet.
Next we fix x = 1/4 and examine the T dependence in detail. In Fig. 7, we have plotted v and the Peclet number for ln (T ) = −4...7.5. As N > 2, drift inversions can be seen around ln (T ) ≈ 2 for both polymer types. As before, the single particle remains the fastest for small T , but eventually the drift curves begin to intersect as T gets larger and the single particle is not always the fastest (see e.g. the N = 3 FM polymer in Fig. 7 (a), right column). The (≪ 0.1) for both potential types, indicating very low coherence of transport.

IV. CONCLUSIONS
We studied the ratchet effect with discrete polymer models in time-dependent potentials using the master equation approach. We gave general equations for calculating the effective diffusion coefficient and drift in time-dependent periodic systems. Using these equations, we performed calculations in the flashing and traveling ratchet potentials for short discrete polymers with the Rubinstein-Duke model and a relaxed version of this model allowing tube breaking. We found complex dynamics that results from the non-pointlike structure of the polymers by the coupling between the potential and polymer internal states. By varying the potential switching rates, we found non-trivial inversions of the polymer drift direction, which cannot occur with simple pointlike non-interacting particles. We also found that the Peclet number grows as the polymer gets longer and is largely independent of the polymer type thus allowing more coherent transport for longer polymers. The overall polymer dy-namics in ratchet potentials was found to be very model specific. The discretization of the problem in this work may be far from many real-world applications but, nevertheless, since our model catches the essential characteristics of the Brownian motor system, we expect that similar properties could be found in the nano-scale objects that can be described with discrete states instead, such as molecular motors with internal structure. Drift inversions are especially interesting since they facilitate more efficient separation methods of molecules.