Reversible to irreversible transitions for cyclically driven particles on periodic obstacle arrays

The transition from reversible to irreversible dynamics was first studied in cyclically sheared dilute colloidal systems where the particles undergo only contact interactions with each other.¹ In these studies, a dilute assembly of colloidal particles in a very viscous medium was sheared periodically, causing collisions among the particles. Over time, if the shear amplitude or density was small enough, the colloidal particles would organize into an arrangement where collisions were avoided and would then return to the same positions after each shear cycle. For larger shear amplitude or particle density, collisions were unavoidable, the motion was irreversible, and the colloidal particles exhibited a long time motion similar to Brownian diffusion both parallel and perpendicular to the shearing direction even though, due to the high viscosity, thermal effects were negligible. Additional studies demonstrated that such systems are generally always in an irreversible state initially but can organize to a reversible state after a number of cycles, which diverges near the critical shear amplitude or density, suggesting that the reversible (R) to irreversible (IR) transition is an example of a nonequilibrium phase transition.² For dilute particles, the transition to the reversible state is called random organization since the system forms a disordered configuration in which collisions between particles do not occur. In additional studies of periodically driven dilute colloidal particles, it was argued that these randomly organized states exhibit hyperuniformity.^3–6 Transitions to irreversible states have also been studied in other periodically driven systems, including granular matter^7,8 and amorphous solids,^9–15 where the particles are always in contact and the transition occurs for a critical shear amplitude. Such systems can also display reversibility spanning multiple cycles^9,15–17 and a variety of memory effects.^10,18–20

In another class of systems of collectively interacting particles that are cyclically driven over quenched disorder,²¹ such as vortices in type-II superconductors,^22–26 magnetic skyrmions,^27,28 and colloidal particles,²⁹ the drive is applied uniformly to all the particles, and in the absence of quenched disorder, the assembly moves back and forth uniformly in a reversible manner. When quenched disorder is present, however, plastic deformations can occur that permit particles to move relative to one another from cycle to cycle.²¹ In vortex systems, studies of R–IR transitions in simulation²² and experiments^23–26 show that diverging time scales appear near critical drive amplitudes²⁵ and critical densities²⁶ with critical exponents similar to those observed in the dilute colloidal systems. Stoop and Tierno²⁹ applied backward and forward pulse driving to hard sphere colloids moving over random obstacle arrays and obtained a variety of different dynamical phases as a function of obstacle density.

R–IR transitions can also occur for particles moving over a periodic array of obstacles or a periodic substrate.²¹ The dynamics of particles coupled to periodic substrates have been explored for superconducting vortices,^26,30–32 magnetic skyrmions,^33,34 and colloidal systems.^35–40 A key aspect of systems with periodic substrates is that the dynamics depend strongly on the direction θ of the drive relative to a substrate symmetry direction. For example, in a square obstacle array, at θ = 0°, the particles can flow easily between the obstacles without collisions. Similarly, other drive angles such as θ = 45° and θ = 90° are also aligned with easy flow directions.^35,41 At incommensurate angles, particles cannot easily travel in a straight line without encountering an obstacle and the flow is more disordered. As the direction of a dc drive is changed relative to the substrate, a series of directional or symmetry locking effects appear, in which the particle motion becomes locked to certain symmetry directions of the substrate even when the drive is not aligned precisely along those directions.^{33,35,39,40,42} This directional locking effect has been studied as a method for particle separation^36,43–45 and in the context of transitions from ordered to disordered flow.^40,41,46 In previous work on dc driven disks moving through square obstacle arrays, it was shown that the system is susceptible to jamming for flow along certain non-symmetry angles.^46,47

Particles cyclically driven over a periodic substrate array provide a convenient system to study reversible to irreversible transitions since the effective frustration of the array can be tuned simply by changing the orientation of the drive relative to the symmetry directions of the array. Frustrated flow occurs for driving angles at which the particles strike the centers of the obstacles so that there is equal probability for the particle to divert to the right or the left of the obstacle. This can lead to highly degenerate flow. In this work, we examine a monodisperse assembly of disks interacting with a square obstacle array under periodic square wave driving. When the drive is applied along θ = 0° or θ = 45°, the system readily organizes to a reversible pattern forming state in which the particles return to the same positions after each drive cycle. For driving at incommensurate angles with a fixed drive amplitude, we find that there is a critical disk density above which an irreversible state forms that exhibits diffusive dynamics. Below this density, the number of cycles required to reach a reversible state varies with density as a power law with the same exponents found for periodically sheared colloidal particles² and periodically driven superconducting vortex systems.^23–26 In general, we find that the disks form a disordered or fluid-like state under irreversible flow, while a pattern forming or ordered configuration appears when the system reaches a reversible state. The critical density for the R–IR transition is non-monotonic as a function of θ, reaching maximum values for commensurate driving angles and showing a global minimum near the arctangent of the inverse of the golden ratio. Our results could be tested in a variety of systems such as colloidal particles, superconducting vortices, or magnetic skyrmions under periodic driving coupled to a periodic array of obstacles or a periodic substrate.

We consider a two dimensional system of size L × L containing a square array of N_obs = 81 circular obstacles of lattice spacing a = 4.0 and radius r_obs. We impose periodic boundary conditions in the x- and y-directions and place N_d monodisperse repulsive disks in the sample. The dynamics of disk i are governed by the following overdamped equation of motion:

This model was previously employed to study locking and clogging effects for dc driven disks.^42,46 The velocity of the disk at position r_i is v_i = dr_i/dt, and we set the damping constant α_d to unity. The first term on the right is the disk–disk interaction force $Fidd$ represented by a short-range harmonic repulsive potential $Fidd=∑j≠iNdk(2rd−rij)Θ(2rd−rij)r̂ij$⁠, where r_d is the disk radius, r_ij = |r_i − r_j|, and $r̂ij=(ri−rj)/rij$⁠. In this work, we fix r_d = 0.55. The disk–obstacle force $Fiobs$ is given by a similar repulsive harmonic interaction, $Fiobs=∑k≠iNobsk(rd+robs−rik)Θ(rd+robs−rik)r̂ik$⁠. We set the harmonic spring constant k = 300, which is the same for both disks and obstacles, large enough to prevent the disk–disk or disk–obstacle overlap from becoming larger than one percent of the disk radius for the densities and driving forces we consider. The density ϕ is defined to be the area covered by the obstacles and mobile disks, $ϕ=Nobsπrobs2/L2+Ndπrd2/L2$⁠.

The square wave driving force remains constant for a fixed period of time in the forward direction prior to reversing and has the form $FD=Acos(θ)x̂+Asin(θ)ŷ$⁠, where A is the drive amplitude and θ is the direction of the drive relative to the x-axis of the periodic array. We fix the duration of the drive to T/2 = 2 × 10⁵ simulation time steps spent on each half cycle, and we vary the density ϕ and the drive amplitude A. Each simulation time step advances time by δt = 0.002 in dimensionless units. The distance d_T that an isolated particle in the absence of obstacles would move during half of a drive cycle is d_T = ATδt/2. Throughout most of this work, we use a drive and period such that in units of the obstacle lattice constant a, d_Ta = 3.07. We can characterize the system by measuring the change in the position of the disks from one cycle to the next, $R(n)=Nd−1∑iNd|ri(t0+nT)−ri(t0+(n−1)T)|2$⁠, where t₀ is an initial reference time. If the motion is reversible, R_n = 0. We also measure the total net displacement d(n) as a function of cycle n, $d(n)=Nd−1∑iNd|ri(t0+nT)−ri(t0)|2$⁠. In an irreversible state, d(n) grows continuously, while for a reversible state, d(n) saturates to a finite value.

In Fig. 1(a), we show a snapshot of the obstacle locations and mobile disks for a system with r_obs = 1.0, A = 0.031 623, and a driving angle of θ = 18.435° for N_d = 269, giving an overall system density of ϕ = 0.3962. Even after 1500 ac drive cycles, the system remains in an irreversible or fluctuating state. When the number of mobile disks is reduced to N_d = 239, giving ϕ = 0.3716, the system organizes into a reversible state with an ordered structure, as illustrated in Fig. 1(b).

In Fig. 2(a), we plot R_n vs cycle number n for the system in Fig. 1 for increasing total densities of ϕ = 0.335, 0.3496, 0.3569, 0.364 27, 0.368 67, 0.3716, 0.3789, 0.3862, and 0.3962. For ϕ < 0.3789, R_n goes to zero, indicating that after an initial transient of some length, the disks return to the same positions after every driving cycle and the system behaves reversibly. As ϕ increases, the number of cycles τ required to reach a reversible state also increases. For example, at ϕ = 0.3716, it takes τ = 1360 cycles to reach the reversible state illustrated in Fig. 1(b), while at lower densities such as ϕ = 0.335, τ = 12. Figure 2(b) shows the corresponding d vs n. For ϕ ≥ 0.3789, d continues to grow as a function of time, while below this density, it saturates to a finite value. The dashed line indicates a fit to d ∝ n^1/2. Since n also corresponds to an elapsed time, this implies that the displacements are growing as t^1/2 and thus have Brownian characteristics, similar to the behavior of the displacement found in the irreversible states of sheared colloidal systems.¹ 

In Fig. 3, we plot the trajectories of the disks from the system in Fig. 2 to illustrate more clearly the difference between the irreversible and reversible dynamics. Figure 3(a) shows that the trajectories in the irreversible state at ϕ = 0.3962 fill space, and the disks are translating in both the x and y directions. In Fig. 3(b), the reversible state at ϕ = 0.3496 contains much more ordered trajectories and the motion is always confined between the rows of obstacles with no hopping from row to row. If the trajectory plot in the reversible state is extended over a larger number of cycles, exactly the same trajectory pattern appears.

In Fig. 4, we plot the number of cycles τ required to reach the reversible state as a function of ϕ − ϕ_c, where we assume a critical density of ϕ_c = 0.3726. The line indicates a fit to $τ=(ϕ−ϕc)−ν$ with ν = 1.36. The previous work on R–IR transitions in two-dimensional (2D) sheared colloidal systems showed a similar divergence in the time to reach the reversible state with ν = 1.33,² while studies of superconducting vortices driven over random disorder gave exponents of ν = 1.38 for critical drive amplitudes and ν = 1.32 for critical densities.²⁵ These exponents are close to those expected for 2D directed percolation, where ν = 1.295.⁴⁸ 

The images in Fig. 3 clarify why the R–IR transition is connected to percolation. In a reversible state, the trajectories do not simultaneously percolate in both the x and y directions, while in an irreversible state, the trajectories are mixing. The percolation transition could be considered to occur at the point where the trajectories just begin to overlap in both the x and y directions.

We next consider the R–IR transition for driving at commensurate angles. In Fig. 5, we plot R_n vs n for the same system in Fig. 2 but at a driving angle of θ = 45° for total densities of ϕ = 0.5124 to 0.562 25. When ϕ < 0.5446, the system organizes to a reversible state in a time τ that grows with increasing ϕ. In Fig. 4(b), we plot τ vs ϕ − ϕ_c for a critical density of ϕ_c = 0.546 12 and a line indicating a fit to τ = |ϕ − ϕ_c|^−ν with ν = 1.38. This result indicates that the R–IR transition shown in Fig. 4 persists for driving along θ = 45°; however, the critical density ϕ_c is higher. In Fig. 6, we illustrate the disk configurations above and below the critical R–IR transition density. Figure 6(a) shows a reversible state at ϕ = 0.5036 where the disks form ordered one-dimensional (1D) chains aligned with the drive along θ = 45°. In Fig. 6(b), the same system at ϕ = 0.577 is in an irreversible state where the disk positions are disordered.

By conducting a series of simulations for fixed driving amplitude A = 0.031 623 and varied θ, we explore the dependence of ϕ_c on the driving angle θ. We plot the reversible and irreversible regions as a function of ϕ vs θ in Fig. 7. Note that for symmetry reasons, the pattern shown in Fig. 7 repeats in an inverted fashion over the range θ = 45° to θ = 90°. At θ = 0°, the system remains in a reversible state up to the largest values of ϕ we consider, ϕ = 0.61. For larger densities, jamming effects become important and we would need to switch to a different disk initialization algorithm. It may be possible that additional R–IR transitions occur at higher disk densities when jammed states begin to appear; however, this is beyond the scope of the present work. For 0 < θ < 7.5°, we find that the R–IR transition occurs near a critical density ϕ_c = 0.575. There is a peak in ϕ_c near θ = 26.565°, which corresponds to a commensurate angle of θ = arctan (1/2). When the driving angle is close to lattice symmetry directions such as 0/1, 1/2, or 1/1, which correspond to θ = 0°, 26.5°, and 45°, respectively, ϕ_c reaches its highest values. Under these commensurate angles, the disks can move easily along straight lines while avoiding collisions with the obstacles. There is no noticeable peak in ϕ_c at θ = arctan (1/3) or θ = arctan (2/3), and the disk dynamics for these driving angles are similar to what is found at incommensurate driving angles. A minimum in ϕ_c occurs near θ = 31° or close to θ = arctan (0.618), where 0.618 is the inverse of the golden ratio from the Fibonacci sequence. The incommensuration is maximized at the inverse golden ratio where the driven disk collides with the largest possible number of obstacles while moving through the system. The variations in the extent of the reversible regions should also depend on the radius r_d of the mobile disks. If smaller disks were used, other possible commensuration effects could appear depending on how many rows of mobile particles can fit along 45° or other commensurate angles. The line dividing reversible and irreversible states in Fig. 7 is also affected by the driving amplitude A and would shift to lower ϕ with increasing A and to higher ϕ with decreasing A. For a sufficiently small driving amplitude, the behavior would be reversible over the entire range of ϕ shown in the figure.

In Fig. 8(a), we illustrate the disk positions in a reversible state at θ = 0° and ϕ = 0.577, where the disks form two nearly filled rows moving in the x direction. We note that not all of the reversible states are associated with ordered disk arrangements. For example, at θ = arctan (1/2), where a peak in ϕ_c appears in Fig. 7, the system forms the disordered but repeatable pattern shown in Fig. 8(b) for ϕ = 0.4302 in the reversible state.

Up to this point, we have concentrated on samples with r_obs = 1.0, but there can be a reentrant R–IR transition as r_obs is varied. In a system with θ = 7.5°, Fig. 7 indicates that an irreversible state appears only for large values of ϕ. In this case, a transition occurs from 1D reversible motion of disks along the x direction to 2D irreversible motion. If we reduce r_obs, the transition to irreversible motion shifts to lower ϕ because the disks can more readily move in two dimensions instead of remaining locked in a 1D channel. If, however, r_obs is reduced even further, collisions with the obstacles become less frequent and the system can once again organize into a reversible state. This is illustrated in Fig. 9 where we plot the locations of the reversible and irreversible regimes as a function of r_obs vs A for a system with N_d = 279 at θ = 7.5°. Here, when r_obs = 2.0, ϕ = 0.401. For large r_obs > 0.9, the system is always in a reversible state regardless of the value of A, and the disks form 1D chains. Furthermore, for A < 0.04, the system is always in a reversible state since the disks do not move far enough during a single drive cycle to collide with the obstacles. When A > 0.04 and 0.4 < r_obs < 0.9, irreversible behavior appears, while for r_obs < 0.4, there is a reversible state in which the disks are moving. The result is the appearance of a reentrant R–IR transition as a function of r_obs for all but the smallest values of A. Similar reentrant transitions should occur near commensurate driving angles such as θ = 45°. In contrast, for incommensurate angles, the system will remain in an irreversible state down to much smaller r_obs since even relatively small moving disks continue to collide with the obstacles due to the driving direction. The existence of reentrance will also depend on the mobile disk density since for low mobile disk densities, the system will generally be able to organize into a reversible state.

We can also observe a R–IR transition at fixed ϕ under increasing A, as illustrated in Fig. 10 for a system with ϕ = 0.549 and θ = 45°. When A = 0.031 62, this system is in an irreversible state. In Fig. 10(a), we plot R_n vs n for A = 0.028 46, 0.027 67, 0.0268, 0.025 247, 0.0253, 0.024 35, and 0.022 136. The system organizes to a reversible state when A < 0.028 46, and the number of cycles τ needed to reach the reversible state decreases with decreasing A. In Fig. 10(b), we plot τ vs A − A_c, where A_c = 0.0285. The line is a fit to τ = |A − A_c|^−ν with ν = 1.4, an exponent slightly larger than what we observe when varying the total density. The exponents are not accurate enough to determine whether the two transitions are in different universality classes. If we focus only on the points closest to A_c, a lower value of ν can be fit, suggesting that both transitions in fact fall in the same universality class.

We have numerically examined the reversible to irreversible transition for periodically driven disks moving through a two-dimensional square periodic obstacle array. For fixed ac drive amplitude, we find that there is a critical density at which the system is able to organize into a reversible state instead of remaining in an irreversible state. The number of cycles required to reach the reversible state diverges as a power law with an exponent ν ≈ 1.36. This is close to the value of ν observed for periodically sheared colloidal particles and periodically driven superconducting vortices, suggesting that the reversible–irreversible transitions of all of these systems fall into the same universality class. The critical density at which the transition occurs is non-monotonic as a function of the angle between the applied drive and a symmetry direction of the obstacle array. The highest critical densities appear for commensurate driving angles such as θ = 0° and θ = 45°. We find the same power law exponents for both incommensurate and commensurate angles. We obtain the lowest critical density for θ = arctan (0.618), which is the inverse of the golden ratio. This degenerate driving direction produces the highest frequency of collisions between disks and obstacles.

We acknowledge the support of the U.S. Department of Energy through the LANL/LDRD program for this work. This work was carried out under the auspices of the NNSA of the U.S. DoE at LANL under Contract No. DE-AC52-06NA25396 and through the LANL/LDRD program.

There are no conflicts of interest to disclose.

The data that support the findings of this study are available from the corresponding author upon reasonable request.