Manipulating Shear-induced Non-equilibrium Transitions in Colloidal Films by Feedback Control

Using Brownian Dynamics (BD) simulations we investigate non-equilibrium transitions of sheared colloidal films under controlled shear stress s xz. In our approach the shear rate _ g is a dynamical variable, which relaxes on a time scale s c such that the instantaneous, configuration-dependent stress s xz (t) approaches a pre-imposed value. Investigating the dynamics under this " feedback-control " scheme we find unique behavior in regions where the flow curve s xz (_ g) of the uncontrolled system is monotonic. However, in non-monotonic regions our method allows to select between dynamical states characterized by different in-plane structure and viscosities. Indeed, the final state strongly depends on s c relative to an intrinsic relaxation time of the uncontrolled system. The critical values of s c are estimated on the basis of a simple model.


Introduction
So matter under shear ow can display rich nonlinear behaviour involving transitions between different dynamical states, 1,2 spontaneous spatial symmetry-breaking, 3 shear-banding 4-7 (for recent reviews, see ref. [8][9][10][11][12], rheochaos [13][14][15] and heterogeneous local dynamics. 16,17 These intriguing phenomena oen strongly affect the rheological properties of the system. Understanding shear-induced effects in, e.g., complex surfactant solutions 18 or liquid crystals, 1,2 colloids, 19-21 so glasses, 16,17,22 and active uids, 23 is thus a major current topic connecting non-equilibrium physics and so material science. A quantity of particular interest is the ow (or constitutive) curve, 7,12 that is, the shear stress s as function of the shear rate _ g, both of which can serve as experimental control parameters.
In many examples, the curve s( _ g) behaves not only nonlinear (indicating shear-thinning, 24,25 shear-thickening, 25,26 sometimes connected irregular (chaotic) rheological response 15,27 ), but becomes also multivalued, i.e. different shear rates lead to the same stress. In complex uids of e.g. wormlike micelles, this multivalued property is in fact a universal indicator of a shearbanding instability, specically, gradient banding, where the (formerly homogeneous) system separates in gradient direction into coexisting bands characterized by a smaller and a larger local shear rate 12 (note that this is different from the more exotic vorticity banding, i.e., the formation of bands along the vorticity direction as discussed e.g., in ref. 11, 12 and 28). In so (colloidal) glasses, multivalued functions s( _ g) occur as transient phenomena aer a sudden switch-on of shear stress (Bauschinger effect), 29 or in the vicinity of the so-called yield stress; 30 in these systems one observes strong dynamical heterogeneities. 17 A further intriguing feature is that close to such nonmonotonicities, the system's behaviour strongly depends on whether one uses the shear stress or the shear rate as a control parameter. 14,31 In fact, both choices are common in rheological experiments. 22,32,33 Here we present BD computer simulation results of yet another system with multivalued stress-shear curve, that is, a thin colloidal lm sheared from an (equilibrium) state within in-plane crystalline order in a planar Couette geometry. As shown in previous experimental 34 and simulation 21,35 studies, such lms can display successive non-equilibrium transitions from square crystalline over molten into hexagonal states. Here we demonstrate that the structural transitions lead again to non-monotonic ow curves, with a shape reminding that of materials which perform a solid-to-liquid transition beyond a critical (yield) stress. 36 Based on this nonlinear behaviour, we then investigate the lms in presence of controlled shear stress. In fact, so far most simulations of sheared colloids have been done under xed shear rate _ g, exceptions being e.g. ref. 17, 30 and 37, where constant s has been realized by xing the force acting on the atoms forming the walls. Here we introduce an alternative, easy-to-apply scheme to control s which has been previously used by us in continuum approaches of sheared liquid crystals. 31 In that scheme _ g (which directly enters the BD equations of motion) becomes a dynamical variable whose time dependence is governed by a relaxation equation involving a time scale s c . The relaxation is based on the difference between the instantaneous (conguration-dependent) stress s(t) and a preimposed value s 0 . The idea of adapting the shear rate _ g is inspired by experiments of stresscontrolled systems. 32,38 Our scheme differs from earlier schemes 17,30,37 where the desired value s 0 is imposed instantaneously. Moreover, due to the coupling to the particle positions our method corresponds to a "feedback" (closedloop) control scheme, which is similar in spirit to e.g. a Berendsen thermostat for temperature control. 39 However, here the choice for s c is found to be crucial for the observed dynamical behaviour. In particular we demonstrate that, if our scheme is applied within the multivalued regime of s( _ g), the nal state strongly depends on the magnitude of s c relative to important intrinsic time scales of the system. Thus, the stress-control can be used to deliberately select a desired dynamical state.

Model and simulation details
We consider a model colloidal suspension where two particles with distance r ij interact via a combination of a repulsive Yukawa potential, u Yuk (r ij ) ¼ W exp[Àkr ij ]/r ij and a repulsive sosphere potential u SS (r ij ) ¼ 43 SS (d/r ij ) 1/2 involving the particle diameter d. 21 The interaction parameters are set in accordance to a real suspension of charged silica macroions [for details, see ref. 40 and 41]. Specically, at the density considered (see below), the interaction strength W/(k B Td) z 123 (where k B T the thermal energy) and the inverse screening length k z 3.2d. Spatial connement is modeled by two plane parallel, smooth, uncharged surfaces separated by a distance L z along the z direction and of innite extent in the x-y plane. We employ a purely repulsive uid-wall decaying as z À9 [see ref . 21]. This is motivated by previous investigations of the equilibrium layer structure, where we found a very good agreement with experiments. 40,41 Our investigations are based on standard BD simulations in three dimensions, where the position of particle i is advanced according to 42 where F i is the total conservative force acting on particle i. Within the framework of BD, the inuence of the solvent on each colloidal particle is mimicked by a friction constant and a random Gaussian displacement dW i . The friction constant is set to (D 0 /k B T) À1 , where D 0 is the short-time diffusion coefficient. The value dW i has zero mean and variance 2D 0 dt for each Cartesian component. The time scale of the system was set to s ¼ s 2 /D 0 , which denes the so-called Brownian time. We impose a linear shear prole [see last term in eqn (1)] representing ow in x-and gradient in z-direction, characterized by uniform shear rate _ g. We note that, despite the application of a linear shear prole, the real, steady-state ow prole can be nonlinear. 43 This approximation has also been employed in other recent simulation studies of sheared colloids; 19,44,45 the same holds for the fact that we neglect hydrodynamic interactions.
One quantity of prime interest in our study is the x-z component of the stress tensor, Thus, we consider only the conguration-dependent contribution to s xz ; the kinematic contribution (which involves the velocity components in x-and z-direction) is negligible under the highly conned conditions here. We note that, apart from the kinematic contribution, eqn (2) also neglects higher-order contributions involving gradients. 46,47 In continuum approaches, one typically includes non-local terms which are essential for the description of interfaces between shearbands. 7,14 The importance of such terms in our highly conned system, which is characterized by pronounced layer formation, remains to be investigated.
Based on the shear stress, we introduce a feedback scheme as follows. Starting from an initial value for _ g we calculate, in each time step, the conguration-dependent stress s xz from eqn (2) and adjust _ g via the relaxation equation where s 0 is a pre-imposed value of s xz , and s c determines the time scale of relaxation. Also, h 0 is the shear viscosity obtained for _ g / 0 (Newtonian regime). This control scheme is inspired by experiments under xed stress [see, e.g. ref. 32], where the adaptation of the shear rate to a new stress value always takes a nite time.
From a more formal point of view, we note that through eqn (3), _ g becomes an additional dynamical variable. Therefore, and since s xy (t) depends on the instantaneous conguration {r i (t)} of the particles, simultaneous solution of the N + 1 equations of motion (1) and (3) forms a closed feedback loop with global coupling. Interestingly, this feedback scheme is in accordance with the common view that, in a stable system, the shear stress s xz should increase with the applied shear rate. This can be shown (at least for a homogeneous system) by a linear stability analysis of eqn (3) as outlined in Appendix A.
In our numerical calculations, we focus on systems at high density, specically rd 3 ¼ 0.85, and strong connement, L z ¼ 2.2d. The corresponding equilibrium system forms a colloidal bilayer with crystalline in-plane structure. 41 We also show some data with L z ¼ 3.2d, yielding three layers. The values L z ¼ 2.2 and 3.2 have been chosen because the equilibrium lattice structure is square 41 (other values would lead to hexagonal equilibrium structures which do not show the shear-induced transitions discussed here). The number of particles was set to N ¼ 1058 and 1587, the width of the shear cell to L z 23.8d and 24.2d for L z ¼ 2.2d and 3.2d, respectively. Periodic boundary conditions were applied in ow (x) and vorticity (y) direction. The time step was set to dt ¼ 10 À5 s where s ¼ d 2 /D 0 the time unit.
The system was equilibrated for 10 Â 10 6 steps (i.e. 100s). Then the shear force was switched on. Aer the shearing was started the simulation was carried out for an additional period of 100s. During this time the steady state was reached. Only aer these preparations we started to analyze the considered systems.

Shear-induced transitions
We start by considering ow curves for systems at constant _ g. The functions s xz ( _ g) for both, two-and three layer systems are plotted in Fig. 1, where we have included data for the viscosities h ¼ s xz / _ g. Note that we have dened h via the externally applied rate _ g rather than via the effective shear rate within the system (characterized by the average velocity of the layers 21 ), which can show signicant deviations from _ g. As seen from Fig. 1, both the bi-and the trilayer systems are characterized by a nonmonotonic ow curve, accompanied by pronounced shearthinning. At small shear rates, the systems display linear stress, related to an Newtonian response of the square in-plane lattice structure. An exemplary simulation snapshot for the three-layer system is shown in Fig. 2a [see ref. 21 for corresponding results for the bilayer]. In fact, within the Newtonian regime (square state) the layer velocities are zero, i.e., the particles are "locked". 21 In this state the lattice structure persists and the system reacts to the displacement of the particles in elastical manner. Further increase of the shear rate then destroys the square order. In the bilayer, the system then enters a "shearmolten" state characterized by the absence of translational order within the layers (as indicated, e.g., by a non-zero longtime diffusion constant in lateral directions 21 ). At the same time, the entire system starts to ow, that is, the layer velocities increase from zero to non-zero values. 21 From the function s xz ( _ g) plotted in Fig. 1a, the appearance of the shear-molten state is indicated by a sudden decrease of s xz , implying the onset of shear-thinning. In fact, with the shear-molten regime the slope of the ow curve is negative. For bulk systems, such a behavior implies that the system is mechanically unstable. 48,49 Here we are considering a strongly conned system, where the macroscopic arguments cannot be immediately applied. Nevertheless, it is an interesting question to which extent the ow curve pertains to a true steady state in the regime where the shear rate has values corresponding to shear-molten congurations. Investigating the shear stress as function of strain (see Appendix B) we nd that s xz assumes indeed a constant value on the time scale of our simulations; however, the relaxation time is extremely long (see next section). We also mention that our observation of shear-molten (long-time) congurations in the regime, where s xz decreases with _ g, is consistent with ndings in an earlier theoretical study of a colloidal suspension under shear. 50 Somewhat different behaviour is found in the trilayer system which displays, before melting, an intermediate state [see Fig. 2b]: this state is characterized by a non-zero layer velocity. In addition, the middle layer separates into two sublayers, in which the particles are ordered in "lanes" [see inset of Fig. 2b] and move with the velocity of the corresponding outer layer (a more detailed discussion of this "laned" state will be given elsewhere 51 ). Only further increase of _ g then yields a shearmolten state characterized by a decreasing ow curve (in analogy to the bilayer). Finally, both systems transform into a state with in-plane hexagonal ordering [see Fig. 2c] and low viscosity. In this hexagonal state the layer velocity is non-zero, 21 that is, the systems ows. As demonstrated earlier by us 21 the mechanism of relative motion involves collective oscillations of the particles around lattice sites, consistent with recent experiments of 3D sheared colloidal crystals. 20 Regarding the stress, we see that the hexagonal regime is (in both systems) characterized by a slight increase of s xz with _ g. As a consequence, there is a parameter range (indicated as region "II" in Fig. 1) where the ow curve is multivalued, that is, different _ g lead to the same s xz . In many contexts (such as for worm-like micelles), multivalued ow curves are associated with the phenomenon of shear-banding, that is, the separation of the system into spatial regions with different shear rates. In our case, where the system consists of two or three layers such a separation can not occur. Instead, the non-monotonic stress curve is a consequence of the structural transitions of the system induced by the shear.

Intrinsic time scales
Before exploring the impact of shear-stress control, which involves a time scale itself through the parameter s c [see eqn (3)], we take a closer look at the intrinsic time scales  In the main figures, the three colours correspond to particles in the three layers. In the insets, the two colours indicate the particles of the middle layer which are positioned at z > 0 (red) and z < 0 (yellow), respectively (where z ¼ 0 is the middle plane).
characterizing the uncontrolled systems. We focus on the bilayer (the same ndings apply qualitatively on the trilayer) and consider the response of the unsheared equilibrium system, which is in a square state, to a sudden switch-on at time s on of shear with rate _ g new . The resulting time dependence of the instantaneous stress is plotted in Fig. 3.
If _ g new s has a value pertaining to the square state, the shear stress jumps at s on to a non-zero value but then settles quickly to its steady-state value [see Fig. 1]. At shear rates corresponding to the shear-molten state we can also observe a relaxation at some non-zero value. However, it should be emphasized that this value is transient in character. The true, steady state value is only achieved at much longer simulation times (see the stress-strain relations presented in Appendix B). Finally, for shear rates related to the hexagonal steady state ( _ g new s > _ g hex s z 257, [see Fig. 1]), we observe a well-pronounced stress overshoot, similar to what is observed e.g. in so glassy systems, 29 wormlike micelles 52 and polymer melts. 53 Closer inspection shows that the actual value of s 1 as well as the functional behavior of s xz (t) strongly depends on the distance between _ g new s and the threshold between shear-molten and hexagonal state, _ g hex s: the smaller this distance is, the larger becomes s 1 , and the more sensitive it is against small changes of the shear rate. Moreover, a sudden quench deep into the hexagonal state leads to an oscillatory relaxation of the stress s xz (t) [see curves _ g new s ¼ 400, 500], with s 1 (which now corresponds to the relaxation time of the envelope) being still quite large. Taken together, for _ g new s > _ g hex s, s 1 can be tted according to (see inset in Fig. 3) where we nd a 1 /s ¼ 0.21, b 1 ¼ 0.52 (setting _ g 1 s ¼ _ g hex s). The oscillations occurring at large _ g new induce yet a different time scale s os , which is smaller than s 1 . Specically, we nd s os /s ¼ O (10 À3 ). The physical reason for these oscillations is the "zig-zag" motion of particles in adjacent layers. 21 This motion is accompanied by periodic variations of nearest-neighbor distances, and thus, pair forces, which eventually leads to oscillations of s xz (t).
Furthermore it is interesting to relate the relaxation times emerging from Fig. 3 to the structural transition from square to hexagonal state in Fig. 1. To this end we consider the Peclet number Pe^_ gs Pe where s Pe is a "typical" relaxation time. 54 Identifying s Pe with s 1 and considering shear rates _ g close to the transition from the quadratic into the shear-molten state, we nd Pe ¼ O (10 0 ). In other words, the shear-induced structural transitions happen at Pe S 1, consistent with our expectations. 54 For comparison we have also investigated the reverse situation, where the system is initially in a hexagonal steady state at shear rate _ g init ¼ 400. We then suddenly change the shear rate towards a much smaller value and consider the relaxation towards the square equilibrium state. The corresponding behaviour of s xz (t) is shown in Fig. 4. Again we nd that, the smaller the difference _ g new s À _ g 2 s is, the larger s 2 becomes (and the more pronounced is the sensitivity to small changes in _ g new ). The resulting relaxation times can be also tted via eqn (4) with a 2 /s ¼ 22.58, b 2 ¼ 1.73 and _ g 2 s ¼ 215, whereby _ g 2 s ¼ _ g sq s denotes the threshold between the square and the molten states. The result for this t is visualized in the inset of the Fig. 4.

Impact of feedback control
We now discuss the impact of our shear stress control scheme dened in eqn (3). The latter involves the zero-shear viscosity, h 0 , which is estimated from Fig. 1 to h 0 ¼ 0.086/ds and 0.323/ds for the bilayer and trilayer, respectively.
The overall dynamical behaviour under feedback control strongly depends on the value of s 0 (imposed shear stress) relative to the ow curve of the original system [see Fig. 1]. We can differentiate between regimes I, II, and III, which are indicated in Fig. 1.  For a s 0 chosen in region I, the response of the system is unique, that is, the nal state is independent of the control time scale s c , as well as of the initial shear rate _ g init and the initial microstructure. Thus, when starting from a square state, with a corresponding initial shear rate _ g init , the system immediately settles at this state. As a more critical test of the injectivity of the ow curve in region I, we plot in Fig. 5a and b the functions _ g(t) and s xz (t) for the bilayer system at s 0 ds 2 ¼ 6 and various s c , starting from a hexagonal conguration (and _ g init s ¼ 400). In all cases, the shear rate decreases towards the value _ gs z 70 and the structure relaxes into the square state pertaining to the value s xz ds 2 ¼ 6 in the uncontrolled system. This indicates that the square state in region I is indeed the only xed point of the dynamics. We also see from Fig. 5a that the relaxation time into this steady state increases with s c . Fig. 5b additionally shows that s xz (t) displays a pronounced peak. The peak indicates the time window in which the initial hexagonal ordering transforms into a square one. In fact, the high values of s xz at the peak reect the large friction characterizing the intermediate molten state. Similar behaviour occurs in region I of the trilayer system [see Fig. 5c and d] where, however, uctuations of s xz (t) are generally larger.
We now choose s 0 within region II of the ow curve, where there are three different shear rates (and thus, three xed points) pertaining to the same stress [see Fig. 1]. We focus on systems which are initially in a square congurations, whereas the initial shear rate _ g init has a value pertaining to the hexagonal state (other initial conditions will be discussed below). The impact of s c on the time dependence of _ g(t) and s xz (t) is shown in Fig. 6. For small values of the control time scale the systems stays in the initial lattice conguration, i.e., _ g relaxes towards the value pertaining to the square state ( _ gs z 90). Different behaviour occurs at larger values of s c /s: although the initial structure is square, the nal state is hexagonal, and the shear rate essentially remains at its high initial value. We stress that these ndings crucially depend on the choice of _ g init . In particular, the dependency of the long-time behaviour on s c /s only arises for large values of _ g init ; for small values the system remains in the square state irrespective of s c . An overview of the nal dynamical states in the feedback-controlled bilayer at s 0 ds 2 ¼ 8 and various combinations of _ g init and s c /s (assuming a square initial structure) is given in Fig. 7. The colour code indicates the ratio of local bond-order parameters hJ 6 /J 4 i [for a denition of the J n see e.g., ref . 21]. The restriction to values hJ 6 /J 4 i # 6 is related to the actual values observed in the simulations. From Fig. 7 one clearly sees that for initial shear rates _ g init s > _ g hex s z 257, the nal state of the feedbackcontrolled system depends on s c /s. This is in contrast to the uncontrolled system which becomes hexagonal for all _ g init > _ g hex . For a hexagonal initial conguration the diagram (not shown here) looks similar from a qualitative point of view; however, the range of control times where the system retains a hexagonal state despite of _ g init < _ g quad (with _ g quad being the threshold between square/molten states) is much smaller.   We conclude that, by varying s c and the initial structure, we can "switch" between the two stable, steady-state congurations arising in the multivalued region of the uncontrolled system. That these states are stable also under feedback (stress) control is supported by the linear stability analysis presented in Appendix A. Indeed, the dynamics under feedback control never evolves towards the intermediate, shear-molten states, consistent with the view that these states are mechanically unstable. This holds also in region III of the ow curve of the uncontrolled system, e.g., for s 0 ds 2 ¼ 16(5) for the bilayer (trilayer): here, a small value of s c yields relaxation towards the square state, whereas for large s c , the system evolves into a hexagonal state. Finally, we note that completely analogous behaviour is found in the trilayer system [see Fig. 6b and c] for a s 0 pertaining to the regime where square, molten and hexagonal states exist.

Transition line
The most signicant observation from Fig. 7 is that at high values of _ g init , the feedback-controlled system can achieve either the hexagonal or the square conguration, provided that we start from a square conguration and choose s c /s accordingly. We now propose a simple model which allows us to estimate the transition values of the control time, s trans c .
The physical idea behind our model is that, with the initial conditions described above, relaxation into the hexagonal state only occurs if the reorganization time s reorg required by the system to transform from a square into a hexagonal conguration, is smaller than the time s decay in which _ g decays to a value pertaining to the square state. We can estimate s decay from eqn (3) if we assume, for simplicity, a linear relationship s xz (t) ¼ m _ g(t) (note that such a relationship is indeed nearly fullled within the square and hexagonal states, see Fig. 1). Under this assumption eqn (3) can be easily solved, yielding _ g(t) ¼ m À1 e Àmt/h0sc (m _ g init À s 0 + e mt/h0sc s 0 ).
From eqn (5) we nd that the decay time of _ g to the threshold value _ g hex s z 257 (below which the hexagonal state of the uncontrolled system is unstable) is given by To estimate the reorganization time s reorg (from the initial square into a hexagonal conguration), we assume that its dependence on _ g init is analogous to that of the relaxation time s 1 introduced for the uncontrolled system [see eqn (4)]. Specically, we make the ansatz As stated above, a crucial assumption of our model is that the system can only reach the hexagonal state if s reorg does not exceed s decay . Note that the latter involves (in fact, is proportional to) the time s c . By equating expressions (6) and (7) for s decay and s reorg , respectively, we can therefore nd an expression for the minimal control time, s trans c , above which the system reaches the hexagonal state, that is Due to the square initial conguration, we set m ¼ h 0 and s xz (t) ¼ m _ g(t) as dened in our ansatz. The remaining parameters a 0 and b 0 are determined by tting the numerical results for s c /s at the boundary [see Fig. 7] to expression (8), yielding a 0 /s ¼ 54.127 and b 0 ¼ 1.503. The resulting line s trans c ( _ g init ) is included in Fig. 7, showing that our estimate describes the transition between square and hexagonal states very well.
Similar considerations are possible, when we use a hexagonal initial lattice structure. Choosing then a small value of g init s we nd that we can switch between hexagonal and square state. This is illustrated in Fig. 8. To obtain the corresponding transition values of s c , we use the same strategy as before, but take a different ansatz for the stress. Specially, we set s xz (t) ¼ n + m _ g(t) which approximately describe the ow curve in the hexagonal state of the uncontrolled system. From the results plotted in Fig. 1 we nd n ¼ 7.0477/ds 2 and m ¼ 0.0025. The analog of the eqn (8) then reads with a 0 /s ¼ 0.012 and b 0 ¼ 0.237 [see Fig. 8]. The result is visualized in Fig. 8. Comparing the typical control time scales at the transition with those seen in Fig. 7 we nd that s c /s which is necessary to switch from square into the hexagonal state [see Fig. 7] is about two orders of magnitude larger than switching from hexagonal into the square state. We suspect that this difference results from the differences of the slope of the shear stress in the square and hexagonal regimes. All simulations were started from the hexagonal initial structure and the imposed shear stress was set to s 0 ds 2 ¼ 8.

Conclusions
Using numerical simulation we have studied the complex dynamical behaviour of sheared colloidal lms under a specic type of shear-stress control. Our approach involves relaxation of the shear rate in a nite relaxation time s c , until the instantaneous stress matches its desired value. This approach is inspired by rheological experiments 32 where the instantaneous shear rate as function of time can be measured. Focusing on systems with multivalued ow curves (resulting from successive non-equilibrium transitions) we have found that, by tuning s c and the initial conditions, it is possible to select a specic dynamical state. In the present system these are either a state with square in-plane ordering and high viscosity, or a hexagonal state with low viscosity. Therefore, our study suggests a way to stabilize states with desired rheological properties, particularly shear viscosities. Moreover, we have proposed a model which relates the transition values of s c to relevant intrinsic relaxation times under sudden change of _ g. Although most of our results pertain to a colloidal bilayer, the fact that we found analogous results for trilayers suggests that the proposed technique can also be applied for systems with larger number of layers. In fact, we think that this method, aer some minor adaptations such as consideration of the kinematic (and, possibly, also the non-local) contributions in eqn (2), should also be applicable and fruitful in bulk systems. Indeed, we expect the method to allow for state selection in any shear-driven system with multivalued ow curve. For example, in an earlier study we have used an analogous approach (based, however, on continuum equations) to select states and even suppress chaos in shear-driven nematic liquid crystals. 31 It therefore seems safe to assume that the capabilities of the present scheme are quite wide. For colloidal layers one may envision, e.g., stabilization of time-dependent structures such as oscillatory density excitation, which may have profound implications for lubrication properties. 55 Finally, our ndings are quite intriguing in the broader context of manipulating nonlinear systems by feedback control. In our case, the feedback character stems from the fact that the stress control involves the conguration-dependent instantaneous stress. Mathematically, this scheme can be viewed as feedback control with exponentially distributed time-delay 56 (as can be seen by formally integrating eqn (3) and inserting it into eqn (1)). Similar schemes are used to stabilize dynamical patterns in laser networks, 57 neural systems, 58 and more generally, coupled oscillator systems. 59 The implications of these connections are yet to be explored. Á s c h 0 À vs xz À g : ; t Á s c h 0 vg : g : 0 À g : À g : For long times we expect the rst term on the right side of eqn (10) to vanish, since s xz ( _ g 0 ,t) / s 0 . To linear order, eqn (10) then reduces to d dt g : z À 1 s c h 0 vs xz À g : ; t Á vg : Noting that the values of s c and h 0 are both positive, we can follow that the feedback controlled shear rate approaches a steady-state value only if vs xz vg c . 0: This corresponds to the usual criterion of mechanical stability. 48

B. Strain-stress relation under constant shear rate
In this Appendix we present results for stress-strain relations at different xed values of _ g. These can be obtained from the data shown in Fig. 3 by rescaling the time axis with the applied shear rate. Numerical results are shown in Fig. 9. Similar to the stresstime relations shown in Fig. 3, one observes simple, monotonic behavior for the case _ g new s ¼ 200 (quadratic regime), while the curves for larger shear rates display pronounced stress overshoots. The width of these overshoots is largest at _ g new s ¼ 240, where the system is in the molten state. This is consistent with the appearance of a particularly large intrinsic relaxation time as discussed in Section 4.