Direct evidence of plastic events and dynamic heterogeneities in soft-glasses
Abstract
By using fluid-kinetic simulations of confined and concentrated emulsion droplets, we investigate the nature of space non-homogeneity in soft-glassy dynamics and provide quantitative measurements of the statistical features of plastic events in the proximity of the yield-stress threshold. Above the yield stress, our results show the existence of a finite stress correlation scale, which can be mapped directly onto the cooperativity scale, recently introduced in the literature to capture non-local effects in the soft-glassy dynamics. In this regime, the emergence of a separate boundary (wall) rheology with higher fluidity than the bulk, is highlighted in terms of near-wall spontaneous segregation of plastic events. Near the yield stress, where the cooperative scale cannot be estimated with sufficient accuracy, the system shows a clear increase of the stress correlation scale, whereas plastic events exhibit intermittent clustering in time, with no preferential spatial location. A quantitative measurement of the space-time correlation associated with the motion of the interface of the droplets is key to spot the long-range amorphous order at the yield stress threshold.
I Introduction
Soft amorphous materials, such as emulsions, foams, microgels and colloidal suspensions, display complex flow properties, intermediate between the solid and the liquid state of matter: they are solid at rest and able to store energy via elastic deformation, whereas they flow whenever the applied stress exceeds a critical yield threshold. The yielding behavior makes such systems as interesting for applications as challenging from the fundamental point of view of out-of-equilibrium statistical mechanics Larson (1999). Some of these systems, referred as simple yield stress fluids (including nonadhesive emulsions and microgels), were shown to flow via a sequence of reversible elastic deformations and local irreversible plastic rearrangements, associated with a microscopic yield stress. These physical ingredients lie at the core of mesoscopic models for soft-glassy dynamics Sollich et al. (1997); Sollich (1998); Fielding et al. (2000, 2009); Hebraud and Lequeux (1998); Picard et al. (2002); Bocquet et al. (2009); Pouliquen and Forterre (2009); Mansard et al. (2013). A challenging question concerns the emergence of features that are non-homogeneous in space (like, for example, shear bandings), where the global rheology is unable to properly capture the complex space-time behavior of the system. One needs to properly bridge between local and global rheology of the soft-glasses, an issue that has been recently addressed in several papers Goyon et al. (2008); Katgert et al. (2010); Sbragaglia et al. (2012). In Goyon et al. (2008, 2010); Sbragaglia et al. (2012); Bocquet et al. (2009) it was suggested that such a bridge can be established by introducing a cooperativity scale which determines correlations (non-local effects) in the flow rheology. The underlying idea is that correlations among plastic events exhibit a complex spatio-temporal scenario: they are correlated at the microscopic level with a corresponding cooperativity flow behavior at the macroscopic level. It is the aim of this paper to study the nature of space non-homogeneity in soft-glassy dynamics and to understand the link with correlations emerging from the dynamics of plastic events. More precisely, we investigate the above issues by using a mesoscopic approach based on the Lattice Boltzmann method Benzi et al. (2009, 2010, 2013), which allows the simulation of emulsion droplets and their interface motion under different load conditions. The simulations provide access to a broad spectrum of scales of motion at a very competitive computational cost, a fact that is instrumental for large-scale simulations of yielding materials, where the dynamics of a collection of a substantial number of droplets needs to be accounted for. The peculiar features of plastic events are investigated below and above the yield stress threshold. Above the yield stress, the “fluidity” model recently introduced by Goyon et al.Goyon et al. (2008, 2010); Jop et al. (2012); Mansard et al. (2013) captures the essential features of the flow: fluidity changes near the boundaries on a scale which is close to the stress correlation scale and to the characteristic scale of plastic events. Near the yield stress, however, the cooperativity scale can not be estimated with enough accuracy, whereas the stress correlation scale shows a clear increase. In this regime, plastic events do not show any preferential location and the system starts to behave as an elastic medium, characterized by near zero fluidity (i.e. large viscosity) and with a long-range amorphous order. Our findings echo some recent results on slowly driven thermal glassesChikkadi et al. (2011) and on driven athermal amorphous materials Maloney and Robbins (2009); Olsson and Teitel (2007).
Ii Dynamic rheological model
We resort to a lattice kinetic model that has already been described in several previous papers Benzi et al. (2009, 2010). Here, we just recall its basic features. We start from a mesoscopic lattice Boltzmann model for non ideal binary fluids, which combines a small positive surface tension, promoting highly complex interfaces, with a positive disjoining pressure, inhibiting interface coalescence. The mesoscopic kinetic model considers two fluids and , each described by a discrete kinetic distribution function , measuring the probability of finding a particle of fluid at position and time , with discrete velocity , where the index runs over the nearest and next-to-nearest neighbors of in a regular two-dimensional lattice Benzi et al. (2009, 2009). In other words, the mesoscale particle represents all molecules contained in a unit cell of the lattice. The distribution functions evolve in time under the effect of free-streaming and local two-body collisions, described, for both fluids (), by a relaxation towards a local equilibrium () with a characteristic time scale : {dmath} f_ζi(r+c_i,c_i;t+1) -f_ζi(r,c_i;t) = -1τLB (f_ζi-f_ζi^(eq) )(r,c_i;t)+F_ζi(r,c_i;t). The equilibrium distribution is given by
(1) |
with a set of weights known a priori through the choice of the quadrature Sbragaglia et al. (2007); Shan (2006). Coarse grained hydrodynamical densities are defined for both species as well as a global momentum for the whole binary mixture , with . The term is just the -th projection of the total internal force which includes a variety of interparticle forces. First, a repulsive () force with strength parameter between the two fluids
(2) |
is responsible for phase separation Benzi et al. (2009). Furthermore, both fluids are also subject to competing interactions whose role is to provide a mechanism for frustration () for phase separation Seul and Andelman (1995). In particular, we model short range (nearest neighbor, NN) self-attraction, controlled by strength parameters , ), and “long-range” (next to nearest neighbor, NNN) self-repulsion, governed by strength parameters , ) {dmath} F^(F)_ζ(r)=-G_ζζ,1 ψ_ζ(r) ∑_i ∈NN w_i ψ_ζ(r+c_i)c_i -G_ζζ,2 ψ_ζ(r) ∑_i ∈NNN w_i ψ_ζ(r+c_i)c_i with a suitable pseudo-potential function Shan and Chen (1993); Sbragaglia and Shan (2011). Despite their inherent microscopic simplicity, the above dynamic rules are able to promote a host of non-trivial collective effects Benzi et al. (2009, 2010). By a proper tuning of the phase separating interactions (2) and the competing interactions (2), the model simultaneously achieves small positive surface tension and positive disjoining pressure . This allows th simulations of droplets of one dispersed phase into the other (see left panel of figure 2) which are stabilized against coalescence. Once the droplets are stabilized, different packing fractions and polydispersity of the dispersed phase can be achieved. In the numerical simulations presented in this paper, the packing fraction of the dispersed phase in the continuum phase is kept the same and approximately equal to . The model gives direct access to the hydrodynamical variables, i.e., density and velocity fields, as well as the local (in time and space) stress tensor in the system, the latter characterized by both the viscous (fluid) as well as the elastic (solid) contributions. Thus, it is extremely useful to properly characterize the relationship between the droplets dynamics, their plastic rearrangements, and the stress fluctuations Sbragaglia et al. (2012).
Iii Numerical evidence of plastic events
To place our results within the proper perspective, we first analyze the global rheological properties of the system under investigation. The computational domain is a rectangular box of size ( is the stream-flow direction) covered by lattice sites. The simulations, performed on latest generation Graphics Processing Units (GPU) Bernaschi et al. (2009), require a few GPU hours for one million time steps, the typical time span of a single run. All quantities will be given in lattice Boltzmann units (lbu) and brackets will be used to indicate averages, either in time (), in space (()), or both. Two different boundary conditions are considered: (a) planar Couette Flow with steady velocity at the boundaries ; (b) Oscillating Strain conditions with strain and boundary velocity . In figure 2 we show a zoom of the configurations resembling the initial conditions (the same for both boundary conditions). For the Oscillatory Strain boundary conditions, the frequency is chosen to guarantee that the stress, , and the strain, , are homogeneous in for very small . We then write , where denotes the maximum value of . In figure 1 we show the resulting shear-stress relation following the definition of the global shear and stress discussed in Benzi et al. (2013). Our simulations provide a yield stress value of about lbu independently of the two load conditions used. The stress is compatible with a Herschel-Bulkley (HB) relation Larson (1999)
(3) |
with . Thus, the material in point shows a non-trivial rheology. The bottom panel of figure 1 reports the normalized velocity profiles, , as a function of the reduced position in a confined steady Couette Flow for different values of the nominal shear . In absence of non-local effects, one would expect the reduced velocity profiles to be a straight line. This is however not the case: the normalized profiles collapse on the same master curve independent of the applied shear, and emphasize that the shear rate is greater at the wall than in the center of the channel. This non-local effect has been discussed in terms of plastic rearrangements of the flow Goyon et al. (2008, 2010); Mansard et al. (2013); Sbragaglia et al. (2012). It is therefore of great interest to provide direct dynamic evidence of such plastic events.
To develop a systematic analysis of plastic events, we perform a Voronoi tessellation
Next, we analyze plastic events and their space distribution, to characterize the transition at the yield stress. In particular, we consider the Couette Flow and compute the number of plastic events which occur, for any , at location . Results are displayed in figure 4 for different values of the average stress . The clear feature emerging from figure 4, is that below the yield stress (), plastic events are distributed almost uniformly in , whereas for there exists a preferential location near the boundary, with a characteristic thickness of the order of . It turns out that such thickness is also close to . Thus, two main messages are conveyed by figures 3 and 4: most plastic events show the same characteristic scale , while their number increases by increasing ; above the yield stress, a preferential concentration of plastic events occurs near the boundaries in a layer of thickness . Although it is not surprising that most of the plastic events concentrate near the boundaries, the fact that for this does not occur, appears to be non-trivial.
Iv Connection with Fluidity Model
Hereafter, we consider the results of section III and establish a connection between the plastic events in droplets rearrangements and the corresponding cooperative flow behavior at the hydrodynamic scale. A step towards this goal has been taken in recent worksGoyon et al. (2008, 2010); Bocquet et al. (2009), where the rate of plastic events is connected to the “fluidity” field, defined as the ratio between the shear rate and the stress, . By using a kinetic model for the elasto-plastic dynamics of the stress distribution function, the local fluidity is shown to obey (in the steady state) the following equation
(4) |
where the scale is a measure of the non-locality of the cooperativity within the flow. The quantity is the bulk fluidity, i.e. the value of the fluidity in the absence of spatial heterogeneities. The bulk fluidity depends upon the local shear rate only, whereas depends upon the position in space. Its value is equal to when the stress and the shear rate are constant in an unbounded geometry, i.e. without the perturbing effects of the boundaries. The fluidity model has been tested with considerable success both in experiments Jop et al. (2012); Goyon et al. (2008, 2010) and in molecular dynamics simulations Mansard et al. (2013). Under the hypothesis of low cooperativity, the model predicts proportionality between the fluidity and the rate of plastic events Jop et al. (2012); Mansard et al. (2013). This feature is strikingly robust, as also evidenced by the work of Nicolas & Barrat, based on a different mesoscopic model of interacting elasto-plastic blocks Nicolas and Barrat (2013). Thus, an increase of the number of plastic events near the boundary should be correlated to a corresponding increase of the fluidity. Also, we may argue that the cooperativity scale should be of the order of , a statement that echoes the results presented by Mansard et al. Mansard et al. (2013), where molecular dynamics simulations with the “bubble model” of DurianDurian (1995) were compared with the fluidity model. A good agreement was found by using a value of of the order of 5 bubbles radii. To check the validity of this interpretation, we investigate the behavior of the Couette Flow above the yield stress . In this case, the mean shear stress is spatially homogeneous, which considerably simplifies the solution of equation (4). We first consider the fluidity averaged in time and in the stream-flow direction. For such a case, the fluidity is predicted to obey a non-local equation of the form Goyon et al. (2008, 2010)
(5) |
where (and hence ) is a constant in the stationary Couette flow. The solution of the fluidity equation requires boundary conditions, i.e. one has to prescribe the value of the fluidity close to the boundaries. When the boundary condition is the same, , the expression of the shear rate reduces to:
(6) |
Several remarks are in order. First, the predictions for the velocity profiles exhibit features in qualitative agreement with the simulation results reported in the bottom panel of figure 1: even though the shear stress is homogeneous, velocity profiles are not straight lines from wall to wall. Deviations from a linear profile close to the wall, extend over a characteristic distance fixed by . The fluidity profiles for different values of the nominal shear rate are reported in the left panel of figure 5. All the numerical simulations are performed above the yield stress . Starting from the wall region, the fluidity field decays towards the bulk value , which can also be deduced from the rheological flow curve reported in figure 1. As for the wall fluidity, , we directly measure it at the distance evidenced by the vertical dots in the left panel of figure 5 and compare it with the bulk fluidity, , in the middle panel of figure 5. The existence of a specific wall rheology is clear: the wall fluidity is significantly larger than the bulk fluidity Goyon et al. (2010). To double-check the quantitative consistency with equation (6), we rescaled all profiles with respect to the wall fluidity, , and studied the quantity . The profiles of the rescaled fluidity collapse on the same curve, consistently with equation (6) and a value of . In line with the notion of cooperativity Goyon et al. (2008, 2010), describing the characteristic scale for non-local effects in the soft-glassy dynamics, we find that and the characteristic scale of plastic events, , are close to each other.
V Connection with Stress Correlation
We can gain further insights by studying the correlation of the stress in the material and exploring its connection with the results presented in sections III and IV. In particular, we measure the stress correlation scale in the system. Let be the stress correlation function defined as:
(7) |
where is the average of the stress along the mainstream direction, i.e. , and where the subscript denotes the connected correlation function. We estimate as the distance away from the location where the correlation function is . In the bottom panel of figure 6, we plot . As one can appreciate, above the yield stress (), the stress correlation scale and the cooperativity scale are basically the same (up to a scale factor, close to ). However, very close to the yield stress, shows a fast growth at decreasing shear. Such an increase of can actually be explained by resorting to a very simple scalar rheological model for the stress field Picard et al. (2002); Takeshi and Sekimoto (2005); Janiaud et al. (2006):
(8) | |||||
(9) |
where is the average mainstream flow speed, the shear, the molecular dynamic viscosity and the elastic modulus. Equation (8) is the momentum conservation relation, while equation (9) is a phenomenological model for the evolution of the stress. Finally, is a relaxation time, diverging close to the yield stress Benzi et al. (2013). Such kind of models have been known for long in the literature Picard et al. (2002); Takeshi and Sekimoto (2005); Janiaud et al. (2006): equation (8) is usually considered in the stationary state, on account of inertia being totally negligible Picard et al. (2002). In the stationary state, with an average stress above the yield stress , equation (9) is consistent with the Herschel-Bulkley global flow curve, equation (3) with
(10) |
Equations (8) and (9) can also be written as:
(11) | |||||
(12) |
where . Finally, ignoring the inertial term represented by (11), we can refer to (12) in order to understand the behavior of the stress correlation functions. Let us remark that, in most cases, such as those presented here, the molecular viscosity is much smaller than the “solid” contribution and we can estimate . This ensures a negligible difference between and . Equation (12) shows that the stress correlation scale should be of the order , which diverges close to in agreement with our findings. In particular, by using equation (3), we can predict . In the inset of the bottom panel figure 6, we plot versus in log-log scale, which shows that our prediction is consistent with numerical data with . Close to the yield stress, we cannot measure by using equation (6), since the fluidity near the wall becomes close to the bulk fluidity and both tend to zero. We were able to obtain accurate measurements only down to , where the cooperative length does not show any substantial variation (see triangles in the bottom panel of figure 6). The computation of , instead, does not result from any best fit procedure and it is an independent measure of space correlations. We note that the increase of near is also consistent with the results shown in figure 4, suggesting that below the yield stress, the system behaves as an elastic medium with long-range order, where plastic events occur without any preferential location.
Other key signatures of the physics below the yield stress are provided in figures 7 and 8. In figure 7, we monitor the space-time distribution of the stress-field, , as well as its average along the mainstream direction, , in a Couette Flow for . The diagonal stripes in the top panel of figure 7 provide a neat signature of propagating stress-waves, which become apparent in close connection with the drop of the average stress, as shown in the bottom panel. This shows that the dynamics of the system supports propagation of stress-waves, in connection with the occurrence of stress-releasing plastic events. Plastic events also show an intermittent clustering in time, as evidenced in figure 8, where we report the area related to plastic events in a Couette Flow simulation at . As we can see, a substantial number of plastic events occur in quite short time intervals, after which quiescent periods are observed. It is tempting to speculate that there are “avalanches” of plastic events. It seems that the stress-waves generated by the first plastic event, trigger a number of other events, each generating stress-waves and, eventually, triggering further plastic events Picard et al. (2004). Such intermittent clustering in time can be indeed quantified by looking at the probability density distribution , being the time interval between two successive plastic events. In figure 9, we show for , and . A striking feature emerges from figure 9: the clustering properties of plastic events are peculiar of the pre-yield condition. Only for , we observe a long tail in which is a clear signature of time intermittency or clustering in the plastic events. Actually, we observe that the tail increases in the course of the numerical simulations; i.e. the system shows aging.
Vi Interface Correlations and Elastic Stress
Although the simple model reported in equation (12) seems to explain the behavior of near the yield stress, it does not reveal much of the underlying physics. Basically, the statement or is a shorthand to characterize the yield stress transition. A more interesting question concerns the physical mechanism characterizing the transition, i.e. the reason why the stress correlation scale increases and/or the relaxation time increases. In this section, we provide a quantitative answer to this question by further exploring the space-time correlations of the elastic stress of the system. We concentrate on the space-time correlations of the motion of the interface which, we argue, are responsible for the increase in the stress correlation scale discussed in the previous section. As shortly outlined earlier on, the picture we have in mind is the following: very close to the yield stress, plastic events take place in an otherwise elastic material Picard et al. (2004, 2005). During a plastic event, the whole interface moves and changes the local stress fluctuations, as well as the interface configuration. Because of the effect of the stress waves, which propagate after the end of the plastic event, the interface may become locally unstable and there is a relatively high probability to trigger further plastic events Krysac and Maynard (1998). Since the motion of the interface induces a large change in the stress fluctuations, the stress correlation scale is large (order the system size). Our qualitative description highlights the link between an increase in the relaxation time of the system () and the increase of space correlation.
To define a quantitative measure of the correlations associated with the motion of the interface, we introduce the phase field . Next, we define the overlap as:
(13) |
where . The physical meaning of is the following: for constant , let us indicate by the space-time average of ; then provides a quantitative measure of how much two field configurations, separated by a time , are on average correlated. Thus, to compute space-time correlations, we need to evaluate the space correlation of
where . Since the change in the configuration of the phase field is due to the interface motion, is a quantitative measure of the space-time correlations of the interface dynamics. By using the Voronoi construction, we have identified plastic events as changes in the topological configuration of the interface. Such changes are actually instantaneous. However, a careful inspection of the dynamics shows that the interface motion associated to a local plastic event, takes a finite time . The value of is not fixed, although it does not show large variations among different plastic events. The characteristic time can be estimated of the order of , where is the stress-wave velocity: for a time scale much longer than it is unlikely that any locally confined source of energy does not radiate out the region where the plastic event occurs. A few numbers may help elucidating the picture. Using and in lbu, we obtain lbu. Note that the time for a stress-wave to propagate from one boundary to the other is , whereas the time scale induced by the external driving is in the range across the yield stress transition, where the longer time refers to the value at . We then consider (the connected correlation function of ) for (figure 10, left panel) and (figure 10, right panel) and for . For , based on the qualitative picture previously described, a clear correlation is expected. Note the long tail in the correlation function for large at : this quantitative measure indicates that the whole interface is spatially correlated on time scales smaller than and comparable to the time scale of the plastic event . However, for , the correlation increases with time due to the propagation of stress waves, whereas for the long tail in the correlation length disappears. Moreover, when the system starts to flow at , stress waves no longer propagate and consequently one cannot observe long-range correlations in the interface motions. Figure 10, therefore, supports our view and indicates the interface motion as the source of the large scale correlation in the stress fluctuations.
Vii Summary and outlook
We have presented quantitative measurements of the statistics and correlations of plastic events, as they arise in the proximity of the yield-stress threshold, obtained by using simulations of concentrated emulsion droplets under soft-glassy conditions. We provide two basic results. First, above the yield stress, the typical spatial scale of the plastic events, , is in a good quantitative match with the cooperativity scale, , introduced by previous authors Goyon et al. (2008, 2010). Both scales are close to the correlation scale of the fluctuating stress within the material, . Among others, a notable result emerging from the above findings is the spontaneous segregation of the plastic events within a near-wall layer of thickness . Second, below the yield stress, shows a clear increase and plastic events exhibit intermittent clustering in time, while showing no preferential locations. This is understood in terms of the long-range amorphous order emerging at the yield stress threshold, where one cannot purport the system as an assembly of mesoscopic elements: the whole interface configuration comes into play during plastic events and the “energy landscape” should be classified in terms of interface configurations with large space-time correlations.
Another important aspect emerging from our analysis is the key role of stress-waves. Usually, having slow flows of soft-glassy materials in mind, one neglects inertial effects in developing mesoscopic models for elasto-plastic materials Picard et al. (2002, 2004). In this work, inertia is not invoked to explain the non-linear rheology of the system, but to allow the propagation of sound waves in the solid, which proves key to sustain long-range dynamic correlations. At low shear rates, experiments are performed to ensure a uniform strain in the system and a nearly constant stress. This is certainly the case when one considers linear rheology in a Couette Flow configuration at very low frequency. Also, the computations performed with the Oscillatory Strain display a clear uniform rate strain and uniform stress for small . However, close to the yield stress, space fluctuations of the stress and the interfaces are crucial to correctly describe the dynamics of the system. As we have seen, stress-waves are able to trigger plastic events and produce an avalanche. Stress-waves can exist only by assuming the active presence of inertial terms. As a matter of fact, mesoscopic models which describe the deformation of elastic solids, do make use of inertia terms Dahmen et al. (2009); Takeshi and Sekimoto (2005). A recent study by Salerno & RobbinsSalerno and Robbins (2013) shows indeed that inertia can strongly influence activity bursts and avalanches in sheared disordered solids.
Overall, all the simulation results presented in this paper refer to a situation where , with the diffusive time associated with molecular viscosity (see equations (8)-(9)), the elastic time for a stress-wave to propagate from one boundary to the other (see section VI) and , where is the frequency at which the storage modulus and the loss modulus cross each other, i.e. . In our case, and is found to be of the order of the characteristic time of plastic events (see section VI), with . A close look at some experimental dataGoyon et al. (2008, 2010); Mason (1999), reveals that is in the range and in the range , thus suggesting that the adopted ordering of time scales is reasonable.
Finally, we wish to highlight the importance of “randomness” and disorder in the initial condition Katgert et al. (2013), which provides a nontrivial feedback to the dynamics. All the simulations presented here have been performed with a small but not negligible polydispersity in the initial configuration. For an ordered hexagonal packing of monodisperse droplets, the yield stress and strain follow from Princen theory Princen (1983); Kraynik (1988). Even a small polidispersity changes the yield strain and opens the way to a much richer and complex dynamics. However, the role of polydispersity or space randomness in the system is still not clearly understood. In particular, preliminary results suggest that an increase in the polidispersity is equivalent to increase the level of “noise” in the system and change the space-time correlations. Although most of the above discussions are rather speculative, we argue that our work may enhance the interest in discussing space-time correlation near the yield stress transition and provide some insights to develop a complete theory of soft-glassy rheology.
We are particularly grateful to A. Scagliarini for his careful reading of the manuscript. M. Sbragaglia & R. Benzi kindly acknowledges funding from the European Research Council under the Europeans Communityâs Seventh Framework Programme (FP7/2007-2013)/ERC Grant Agreement No. 279004. P. Perlekar & F. Toschi acknowledge partial support from the Foundation for Fundamental Research on Matter (FOM), which is part of the Netherlands Organisation for Scientific Research (NWO).
Footnotes
- The Voronoi tessellation has been performed by using the Voro++ library. https://math.lbl.gov/voro++
References
- R. G. Larson, The Structure and Rheology of Complex Fluids, Oxford University Press, 1999.
- P. Sollich, F. Lequeux, P. Hebraud and M. E. Cates, Phys. Rev. Lett., 1997, 78, 2020.
- P. Sollich, Phys. Rev. E, 1998, 58, 738.
- S. M. Fielding, P. Sollich and M. E. Cates, J. Rheol., 2000, 44, 323.
- S. M. Fielding, M. E. Cates and P. Sollich, Soft Matter, 2009, 5, 2378.
- P. Hebraud and F. Lequeux, Phys. Rev. Lett., 1998, 81, 2934.
- G. Picard, A. Ajdari, L. Bocquet and F. Lequeux, Phys. Rev. E, 2002, 66, 051501.
- L. Bocquet, A. Colin and A. Ajdari, Phys. Rev. Lett., 2009, 103, 036001.
- O. Pouliquen and Y. Forterre, Phil. Trans. R. Soc. London A, 2009, 367, 5091.
- V. Mansard, A. Colin, P. Chauduri and L. Bocquet, Soft matter, 2013, 9, 7489–7500.
- J. Goyon, A. Colin, G. Ovarlez, A. Ajdari and L. Bocquet, Nature, 2008, 454, 84–87.
- G. Katgert, B. Tighe, M. Mobius and M. V. Hecke, Europhys. Lett., 2010, 90, 54002.
- M. Sbragaglia, R. Benzi, M. Bernaschi and S. Succi, Soft Matter, 2012, 8, 10773–10782.
- J. Goyon, A. Colin and L. Bocquet, Soft Matter, 2010, 6, 2668–2678.
- R. Benzi, M. Sbragaglia, S. Succi, M. Bernaschi and S. Chibbaro, J. Chem. Phys., 2009, 131, 104903.
- R. Benzi, M. Bernaschi, M. Sbragaglia and S. Succi, Europhys. Lett., 2010, 91, 14003.
- R. Benzi, M. Bernaschi, M. Sbragaglia and S. Succi, Europhys. Lett., 2013, 104, 48006.
- P. Jop, V. Mansard, P. Chaudhuri, L. Bocquet and A. Colin, Phys. Rev. Lett., 2012, 108, 148301.
- V. Chikkadi, G. Wegdam, D. Bonn, B. Nienhuis and P. Schall, Phys. Rev. Lett., 2011, 107, 198303.
- C. E. Maloney and M. O. Robbins, Phys. Rev. Lett., 2009, 102, 225502.
- P. Olsson and S. Teitel, Phys. Rev. Lett., 2007, 99, 178001.
- R. Benzi, S. Chibbaro and S. Succi, Phys. Rev. Lett., 2009, 102, 026002.
- M. Sbragaglia, R. Benzi, L. Biferale, S. Succi, K. Sugiyama and F. Toschi, Phys. Rev. E, 2007, 75, 026702.
- X. Shan, Phys. Rev. E, 2006, 73, 047701.
- M. Seul and D. Andelman, Science, 1995, 267, 476–483.
- X. Shan and H. Chen, Phys. Rev. E, 1993, 47, 1815–1819.
- M. Sbragaglia and X. Shan, Phys. Rev. E, 2011, 84, 036703.
- M. Bernaschi, L. Rossi, R. Benzi, M. Sbragaglia and S. Succi, Phys. Rev. E, 2009, 80, 066707.
- A. Nicolas and J.-L. Barrat, Phys. Rev. Lett., 2013, 110, 138304.
- D. Durian, Phys. Rev. Lett., 1995, 75, 4780.
- O. Takeshi and K. Sekimoto, Phys. Rev. Lett., 2005, 95, 108301.
- E. Janiaud, D. Weaire and S. Hutzler, Phys. Rev. Lett., 2006, 97, 038302.
- G. Picard, A. Ajdari, F. Lequeux and L. Bocquet, Eur. Phys. J. E, 2004, 15, 371–381.
- G. Picard, A. Ajdari, F. Lequeux and L. Bocquet, Phys. Rev. E, 2005, 71, 010501.
- L. C. Krysac and J. D. Maynard, Phys. Rev. Lett., 1998, 20, 4428–4431.
- K. A. Dahmen, Y. Ben-Zion and J. T. Uhl, Phys. Rev. Lett., 2009, 102, 175501.
- K. M. Salerno and M. O. Robbins, arXiv:1309.1872.v1, 2013, 1–14.
- T. Mason, Current Opinion in Colloid & Interface Science, 1999, 4, 231–238.
- G. Katgert, B. P. Tighe and M. V. Hecke, Soft Matter, 2013, 9, 9739.
- H. M. Princen, J. Colloid Interface Sci., 1983, 91, 160â175.
- A. M. Kraynik, Annu. Rev. Fluid Mech., 1988, 20, 325â357.