Friction on Mesoscopic Scale

Collective effects at frictional interface  

Equations (21−24) form the complete set of equations which describes evolution of the tribological system; in a general case it has to be solved numerically. However, a qualitative picture may be obtained analytically. The interface dynamics depends on whether or not the λ-blocks undergo the elastic instability, i.e., on the ratio of the stiffness of the λ-block  K_λ≈ (2c_l²+ 3c_t²)ρλ_c  [as follows from the discretized version of Eq.(23)] and the effective critical stiffness parameter of the interface  K_λ^*_eff = β₀K_λ^*,  where  K_λ^*~ K_λsx_c/Δx_s  and  K_λs= N_λk_c. If the elastic instability does not appear, then a local perturbation at the interface relaxes, spreading over an area of size λ_s − the screening length considered in next section. In the opposite case, when the elastic instability does emerge (locally), in may propagate through the interface. Below we consider a simplified one-dimensional version, which allows us to get some analytical results (the using of the 1D model is also supported by the fact that the largest forces near the broken contact are just ahead/behind it − see Fig.I1). Recall that the interaction between the λ-blocks is weaker than in the short-range zone, it follows the law δf ∝ r^–3 which determines, e.g., the block-block interaction strength κ_λ in Eq.(26) below (although the interaction is power-law, we may consider the interaction of nearest neighbors only, because excitations at the interface, such as "kinks", are localized excitations, and the role of long-range character of the interaction reduces to modification of their parameters).

Elastic screening length

Let us introduce the dimensionless variable ε_n= x_cΣ_m≠n Γ_m(X_m) Π_mn. The shift of forces in the nth block due to broken contacts in the neighboring blocks may be accounted by a renormalization of the rate:

Pn(u) → Pn [ (1+εn) u ] .

(27)

Indeed, when contacts in the neighboring blocks break, then the forces in the given block increase, ε_n> 0, and the contacts in the given block should start to break earlier, i.e., their threshold distribution effectively shifts to lower values.

Making the transition from the discrete sliding blocks to a continuum sliding interface, Π_mn→ Π(r'−r)  and  ε_n→ ε(r), we obtain a master equation of the form:

∂Q[u; X(r); r] ∂X(r) + ∂Q[u; X(r); r] ∂u + P([1+ε(r)] u) Q[u; X(r); r] = δ(u) Γ[X(r); r] ,

(28)

where we again assumed that the contacts are reborn with zero stretchings, R(u) = δ(u),

ε(r) = xcλc–2 ⌠ ⌡ |r'-r| ≥ λc d2r' Γ[X(r'); r'] Π(r'−r)

(29)

and

Γ[X(r); r] = ⌠ ⌡ dξ P( [1+ε(r)] ξ ) Q[ξ; X(r); r] .

(30)

In the long-wave limit, when |dε(r)/dr| << ε(r)/λ_c, we may assume that the interface is locally equilibrated, i.e., the distribution of forces on contacts is close to the steady-state solution of the master equation, Q[u; X(r); r] ≈ Q_s(u; r), which depends parametrically on the coordinate r through the function ε(r) entered into the expression for the rate P( [(1+ε(r)] u ). The stationary solution of the ME is known analytically, thus we may find the function (30), Γ(r) = [1 + ε(r)]/x_c. Together with Eq.(29) this gives a self-consistent equation on the function ε(r):

ε(r) = λc–2 ⌠ ⌡ |r'-r| ≥ λc d2r' [(1+ε(r')] Π(r'−r) .

(31)

Taking into account the interaction of nearest neighboring λ-blocks only and expanding ε(r) in Taylor series, we obtain the equation

ε(r) = Π0 [ 1 + ε(r) + 1/2 λc2ε"(r) ] ,

(32)

where Π₀= νΠ(λ_c) = νN_λκ_λ~ νκλ_c/a  and  ν = 2 − 4 is the number of nearest neighbors. Writing ε(r) = ε₀+ Δε(r), where ε₀= Π₀/(1−Π₀),  Eq.(32) may be rewritten as

λs2 Δε"(r) = Δε(r) ,

(33)

where λ_s= λ_c(ε₀/2)^1/2 is the characteristic screening length in the sliding interface.

From the known analytical solution of the ME for the steady state, we may predict the dependence of screening length on temperature and the sliding velocity v. In particular, if T > 0, then λ_s∝ v^–1/2→ ∞  when  v → 0.

Self-healing frictional crack as a solitary wave

Figure 0: Left: Schematic picture showing the propagation of the friction fronts generated at weak points in the center of the system. The thick (red) lines represent λ-contacts which have been stretched by the shear. The (blue) dotted lines represent λ-contacts which have relaxed by sliding. As time grows the relaxed region extends and its boundaries make two kinks separating domains where the λ-contacts are in different states. Right: domino effect.

When one of the "collective contacts" (the λ-contact) breaks, it may initiate a chain reaction, with contacts breaking domino-like one after another. This scenario may be described accurately by reducing the system of contacts to a FK-like model. Recall that the FK model describes a chain of harmonically interacting atoms subjected to the external periodic potential V_sub(x) of the substrate. If the atoms are additionally driven by an external force  f,  then the equations of motion for the atomic coordinates u_n take the form

m d2un/dt2 + mη dun/dt − g (un+1+ un-1− 2un) + V 'sub(un) = f ,

where m is the atomic mass, g is the strength of elastic interaction between the atoms, and η is an effective damping coefficient which describes dissipation phenomena such as the excitation of phonons in the substrate. The main advantage of using the FK model is that its dynamics is well known. Mass transport along the chain is carried by kinks (antikinks) − local compressions (extensions) of otherwise commensurate structure. The kink is a well-defined topologically stable excitation (quasiparticle) characterized by an effective mass m_k which depends on the kink velocity v_k,  m_k= m_k0(1−v_k²/c²)^–1/2 (the relativistic Lorentz contraction of the kink width when its velocity approaches the sound speed c). Therefore, the maximal kink velocity is  v_{k max}= c. In the discrete chain, kinks move in the so-called Peierls-Nabarro (PN) potential, whose amplitude is much lower than that of the primary potential V_sub(x). Therefore, the kink motion is activated over these barriers, and its average minimal velocity  v_{k min} is nonzero. The steady-state kink motion is determined by the energy balance: the incoming energy (because of action of the external driving force  f) should go to creation of new "surfaces" (determined by the amplitude of the substrate potential) plus excitation of phonons by the moving kink (described by the phenomenological damping coefficient η), so that v_k(f) = f /(m_kη).

FK-like model. Thus, let us consider a chain of λ-contacts ("atoms" of mass m = ρλ_c³), coupled harmonically with an elastic constant g, driven externally through a spring of elastic constant K with the spring's end moving with a velocity v. Using the discretized version of Eq.(23), the elastic constants may be estimated as  g ≈ 2Eλ_c(1−ν)/[(1−2ν)(1+ν)]  and  K ≈ λ_cρc_t². The λ-contacts are coupled "frictionally" with the bottom substrate; the latter is described by the nonlinear force F_s(u). The equation of motion of the discrete chain is

Let us define the function F(u) = F_s(u) + Ku − f.  The degenerated ground states of the chain are determined by the equation F(u) = 0. Let the right-hand side (n → ∞) of the chain be unrelaxed, k_cu_R + Ku_R = f,  or

uR = f / (kc+ K) ,

(36)

while the left-hand side (n → −∞) already undergone relaxation, k_c(u_L− u_c) + Ku_L = f,  or

uL = (f + kcuc) / (kc+ K) .

(37)

Thus, the FK-like model of friction (the FK-ME model) is described by Eqs.(34) and (35) with the boundary conditions given by Eqs.(36) and (37).

Figure 6a: The effective FK-like "potential" Veff(u) = ∫du F(u)  for uc= 1, k = 1, K = 0.4 and different values of the driving force:  f = 0 (dashed, no driving),  f = fm= Kuc (red dotted, when the second local minimum appears),  f = fmin= (0.5k+K )uc (solid magenta, the two local minima are characterized by the same energy), and  f = fmax= (k+K)uc (blue das-dotted, the left minimum disappears). Inset shows graphical solution of the equation F(u) = 0, or F(u) = f − Ku.

Continuum-limit approximation. Let the system be overdamped (d²u/dt²=0); later on we shall remove this restriction. In the continuum-limit approximation, n → x=na  (a=1), the motion equation takes the form

mηut − a2guxx+ F(u) = 0,      F(u)|x → ±∞ = 0 .

(38)

We look for a solution in the form of a wave of stationary profile (the solitary wave), u(x,t) = u(x−v_kt), so that u_t= −v_ku' and u_xx= u". In this case Eq.(38) takes the form

mηvku' + a2gu" = F(u) ,

(39)

which may be solved analytically by standard methods.

A solution of Eq.(39) with these boundary conditions exists only for a certain value of the kink velocity v_k, defined by the equation

(mηvk)2 = ga2(kc+ K)(2−β)2/(β−1) ,

(40)

where β = k_c/(k_*−K) and k_* = f/u_c. The solitary-wave solution exists for forces  f_min< f < f_max only. The minimal force which supports the kink motion − the Griffith threshold − is given by

fmin= (1/2 kc+ K) uc .

(41)

The maximal force, for which a kink may exist, is given by

fmax = (kc+ K) uc ;

(42)

at higher forces, the barriers of F_s(u) are degraded, the stationary ground states disappear, and the whole chain must switch to the sliding state.

From Eq.(40) we can find the kink velocity as a function of the driving force. At low velocities

vk ≈ (f − fmin) /mkη ,

(43)

where we introduced the effective kink (crack) mass

mk = m /   4a uc   √ g kc   ( 1 +  K kc  )  ,

(44)

while at  f → f_max the velocity tends to infinity,

mηvk ≈ √ gkc(kc+ K) a2uc (fmax− f)  .

(45)

The latter limit should be corrected by taking into account inertia effects. The term md²u/dt² in Eq.(34) gives mv_k²u" for the solitary-wave solution, so it can be incorporated if we substitute in the above equations g → g_eff = g (1 − v_k²/c₀²),  where c₀= (ga²/m)^1/2 is the sound speed along the chain. The high-velocity limit now takes the form

vk ≈ c0 /   √ 1 +   mη2 (fmax− f) kc (kc+ K) uc  .

(46)

Simulation. The continuum-limit approximate is accurate for the case of strong interaction between the contacts, g >> 1; in the opposite limit one has to resort to computer simulation. We solved Eq.(34) by the Runge-Kutta method. As the initial state, we took the chain of length N  (typically N = 3×10³ or 3×10⁴) with periodic boundary conditions and all contacts relaxed, but the threshold breaking value for two central contacts was taken lower than for the other contacts. Then the driving force increases because of stage motion, two central contacts break first and initiate two solitary waves of subsequent contact breaking which propagate in the opposite directions through the chain. The value k_c' of the lower threshold of the central contacts determines the driving force and therefore the kink velocity; the lower this threshold, the lower the threshold force for the motion to start. As soon as the kink motion is initiated, the k_c values of the central contacts are restored to the same value as for other contacts (otherwise these contacts will act as a source of creation of new and new pairs of kinks), and we change the stage motion to in the opposite direction, v > 0 → v_b< 0, so that the driving force linearly decreases with time (see Fig.7b) and the average chain velocity  <du_i/dt> = N^–1Σ_idu_i/dt  decreases as well (Fig.7c) until the motion stops (Fig.7a).

Contrary to the continuum-limit approximation, in the discrete chain of contacts the kink oscillates during motion (see Fig.7d) − the well-known discreteness effect of the FK model due to existence of the PN barriers. The stronger the elastic interaction between the contacts, the larger the kink "width" and the smaller the kink oscillations (compare Figs.8a and 8b). The amplitude of oscillations also depends on the shape of the "substrate potential" − it is larger for a sawtooth potential F_s(u), but smaller for a smoother shapes. Recall that the λ-contacts are characterized by a smooth dependence F_s(u) as follows from the solution of the ME (for example, see ). The PN oscillations determine the lowest average kink velocity. Therefore, the lowest velocity allowed for the frictional crack propagation, v_{k min}, is determined by the parameters g and λ_c − the larger are g and λ_c, the smaller is v_{k min}.

Discussion. The FK-ME model used here is rather close to the well-known 1D Burridge-Knopoff (BK) model of earthquakes with a velocity-weakening friction law. The difference is in the interface force F_s(u): we use the function derived from the ME-EQ model (with well-defined parameters which may be extracted from experiments or calculated from first principles), whereas the BK model adopts a phenomenological velocity-dependent function for F_s. Nevertheless, the qualitative behavior of the two models is similar, the BK model also exhibits solitary-wave dynamics as was demonstrated numerically by Schmittbuhl et.al. [Schmittbuhl J., Vilotte J.-P. and Roux S.: Propagative Macrodislocation Modes in an Earthquake Fault Model. Europhys. Lett. 21, 375-380 (1993)]. In our case, however, by reducing the model to the FK-ME one, we can describe the solitary waves analytically.

In the simulation we started from the well-defined initial configuration, when all contacts are relaxed except the one or two where kink's motion is initiated. If one starts from a random initial configuration, we expect that kinks will emerge at random places and several kinks may propagate through the system simultaneously.

Also we assumed that all λ-contacts are characterized by the same F_s(u) dependence and thus have the same threshold values F_th. This is correct if the number of original contacts within a single λ-contact, N_λ= (λ_c/a_c)², is infinite. Otherwise, different λ-contacts will have different threshold values, although the distribution of the thresholds is narrower that the distribution of thresholds of single asperities by a factor √N_λ. A narrow distribution of thresholds will nevertheless have a qualitative effect because rupture fronts may stop when they meet a λ-contact with a threshold above the driving force.

Our approach may also incorporate the existence of disorder and defects always present in real materials. On the one hand, defects may nucleate kinks (cracks); on the other hand, the kink propagation may be slowed down up to its complete arrest due to pinning by the defects.

Thus, reducing the EQ-ME model of friction to the FK-ME one, we described analytically avalanche-like dynamics of the frictional interface − the solitary wave of contacts breaking. The analogy with the FK model may be extended even further: (i) Effects of nonzero temperature may be considered; at T > 0 the sliding kinks will experience an additional damping, while the immobile (arrested) kinks will slowly move (creep) due to thermally activated jumps; (ii) As was shown, a fast driven kink begins to oscillate due to excitation of its shape mode, and then, with the further increase of driving, the kink is destroyed; this effect is similar to what is observed in fracture mechanics, where cracks begin to oscillate and then branch; (iii) One may suppose that the damping coefficient η in the equation of motion (34) depends on the kink velocity, η(v); in fracture mechanics, this coefficient defines the rate at which the energy is removed from the crack edge, thus it plays a crucial role.

Published in: O.M. Braun and M. Peyrard, Phys. Rev. E 85 (2012) 026111 "Crack in the frictional interface as a solitary wave"

Onset of sliding: propagation length of the self-healing crack

Recent experiments allowed visualization of the onset of sliding in a variety of experimental situations. These experiments have revealed that the onset of sliding is always mediated by the propagation of micro-slip (rupture) fronts, called precursors, which separate a stuck region and a slipping region along the contact interface. In side-driven plane-on-plane contacts, the precursors were nucleated near the trailing edge of the contact, where the loading is applied. The first precursor propagates through the shear and normal stress fields produced by the initial loading of the system only. After the precursor has stopped, the mechanical state of the interface is modified all along the ruptured path, so that the next fronts propagate through the stress field prepared by the previous ones.

Several models have been proposed to investigate the series of precursors to sliding. Braun et al. (2009) used dynamic simulation in which the slider was modeled as a chain of blocks coupled by springs, while the interface was treated as a set of "frictional" springs with random breaking thresholds. They reproduced both a series of precursors and a modification of the stress field by the precursors as observed at experiments. These results stimulated the study of different simplified models. In these models the interface was described by the phenomenological Amontons-Coulomb law of friction: the interface is pinned until the local ratio of shear to normal stress reaches the static friction coefficient μ_s, and then the interface is locally slipping and characterized by a kinetic friction coefficient μ_k< μ_s. All these models produce a series of precursors nucleated at the trailing edge. Maegawa et al. (2010) further showed that an asymmetric normal loading influences the precursor length. Scheibert and Dysthe (2010) showed that a similar effect is obtained for a symmetric normal loading, as soon as the friction force is not applied exactly in the plane of the contact interface. Amundsen et al. (2012) showed that, to be realistic, 1D models should involve an intrinsic length scale for the interfacial stress variations. Trømborg et al. (2011) showed that, in 2D models, such a length scale naturally arises when introducing the system's boundary conditions and bulk elasticity. In all the above simplified models, each broken contact point repins when its slipping velocity goes back to zero. As a consequence, the precursor front separates a stick region from a region which is slipping behind it.

Above, however, we developed the multiscale model of the frictional interface that does not use the phenomenological Amontons-Coulomb friction law: the friction between two solids results from multiple individual micro-junctions which break under stress and form again elsewhere. As a consequence, the rupture mode of the interface corresponds to the self-healing crack, i.e. the slipping part of the interface is confined between the rupture front and a repinning front; such cracks may be treated as solitary aves. In the previous section we used a too simplified 1D model of the self-healing crack propagation: first, we postulated the elastic interaction between the nearest neighboring contacts, and second, we assumed the uniform constant driving force acting on the contacts. Both these issues, however, are produced by the elastic body of the slider. A full treatment of the problem requires a solution of Eqs.(21−24) which can be done numerically only. Below we propose a more involved 1D model, where the slider elastic body is modeled  as an additional 1D chain of blocks. This model, being still oversimplified, nevertheless allows a rigorous analytical treatment, describes the kinematics of the onset of sliding and predicts the propagation length Λ of the first precursor as well as the shear stress field associated with it.

1D model of elastic slider. We consider that the macroscopic contact is made of a large number of micro-junctions in parallel, with an average distance a_c. Let an individual junction (micro-contact, bridge, solid island, etc.) have an average elastic constant k_c. Elastic theory introduces the elastic correlation length  λ_c below which the frictional interface behaves rigidly. The set of N_c= λ_c²/a_c² contacts within the area λ_c² forms an effective macro-contact (the λ-contact introduced above) with the elastic constant k = N_ck_c. The λ-contacts are elastically coupled through the deformation of the slider's bulk. If we denote by u_i the displacement of the point on the slider's bottom to which the ith λ-contact is attached, then the elastic energy stored between two nearest neighboring λ-contacts i and j in a non-uniformly deformed slider may formally be written as (1/2)K (u_i− u_j)², where K is the slider rigidity defined below.

The system can be described by three variables: u_t(x) describes the displacement field of the TB, u(x) is the displacement field of the blocks in the IL, and u_b(x) is the displacement of attachment points of the IL to the BB. The stress in the IL is then given by

σc(x) = k [ u(x) − ub(x) ] /λc2 ,

(J4)

while the "driving" stress is equal to

σd(x) = KT [ ut(x) − u(x) ] /λc2 .

(J5)

Self-healing crack-like rupture. The frictional behaviour of each junction connecting the IL to the BB is the following: it acts as a spring of elastic constant k as long as its stretching remains below a critical value u_s= f_s/k but breaks and then repines with zero stretching when the local shear stress reaches σ_s= f_s/λ_c².

When a loading is applied to the slider, it is transmitted to each λ-contact due to elasticity of the slider. The leftmost λ-contact is the first to reach its breaking threshold, slide and locally relax the interface. This causes an extra stress on the neighboring contacts, which tend to slide too, so that sliding events propagate as a kink, extending the initial relaxed domain. Such a scenario may be described as a propagation of the solitary wave, i.e., the rupture mode corresponds to the self-healing crack with crack's width equal to the lateral size of one contact λ_c.

In order to simplify the analysis, we suppose the existence of a hierarchy of times: we assume that the pushing rate is adiabatically slow, v_d→ 0, while the sound speed (which determines the rate of propagation of the elastic stress) is very large, v_R→ ∞. Therefore, when a front propagates, it is accompanied by a quasi-static stress field defined by mechanical equilibrium of the forces acting on the interface layer (thus, inertia effects are ignored here). In the following, we show that there exist a typical length scale Λ which characterizes the stress accumulation and thus the precursor propagation length.

Equations. In the discrete model we have x = iλ_c, u_t(x) → u_t,i , u(x) → u_i  and u_b(x) → u_bi. In equilibrium, the displacements should satisfy

A key trick of our analytical approach is that the general solution of the equation  u''(x) = κ²[u(x) − w(x)]  is  u(x) = De^±κx − κ ∫^xdξ w(ξ) sinh[κ(x − ξ)].

As we impose the position of the left-hand side (trailing edge) of the TB we have

ut(0) = Ut ,

(J11)

while the left-hand side of the IL is free so that

K λcu'(0) + KT [Ut− u(0)] = k [u(0) − ub(0)] .

(J12)

To simplify analytical consideration, we ignore the right-hand side boundary conditions assuming that the chain is infinite, L→ ∞, so that u_t(∞) = u(∞) = u_b(∞) = 0.

Let initially, at t=0, the system be completely relaxed: u_t(x) = u(x) = 0 and u_b(x) = u_b⁽⁰⁾(x) = 0. We then start to push slowly the trailing edge of the TB, so that U_t=v_dt  with v_d→ 0. The solution of Eqs.(J8−J12) is:

ut(x) = D10e–κ1x + D20e–κ2x ,

(J13)

u(x) = α1D10e–κ1x + α2D20e–κ2x ,

(J14)

where

κ1,22 =  1 2   [ (κ2 + κT2) ± √   D   ] ,

(J15)

α1,2 =  1 2   [ 1 − ( κ2 ± √ D ) / κT2 ] ,

(J16)

D = (κ2 − κT2)2 + 4βκ2κT2 ,

(J17)

D10 = (α2Ut − uc) / (α2− α1) ,

(J18)

D20 = (uc− α1Ut) / (α2− α1) ,

(J19)

and

uc = Ut   λc β κ2 (κ1− κ2) + (KT /K) √ D λcκT2 (α2κ2− α1κ1) + λc2κ2 √ D  .

(J20)

The positions U_t and u_c grow together until the trailing edge of the TB reaches a position U_t0 where the left end of the IL achieves the threshold value u_s at some time t₁= U_t0/v_d. At t = t₁ the left-most λ-contact breaks and attaches again with zero stretching, u_b(0) = u_b⁽¹⁾(0) = u_s, and the whole system relaxes, adjusting itself to a new configuration with u_t(0) = U_t0, u_b(0) = u_s and u_b(x ≥ λ_c) = 0. After relaxation, all contacts have shifted to the right, so that the IL stress σ_c(x) for x ≥ λ_c increases. As a consequence, the stress on the second contact has grown above the threshold σ_s so that it must break too. This domino-like process will continue until the stress at some (s₀th) contact will remain below the threshold. This first passage distance defines a characteristic length Λ = s₀λ_c which controls the kinematics of the onset of sliding.

We emphasize that the driving stress σ_d(x), Eq.(J5), changes self-consistently during crack propagation, as u(x) changes. The moving self-healing crack leaves behind itself relaxed contacts. When the crack reaches a position s, the current shapes of the displacement fields u_t(x; s) and u(x; s) are given by solution of Eqs.(J8−J10). The current displacement of the bottom surface of the IL, u_b(x; s), is determined by the current propagation length of the crack:

ub(x; s) = u(x+0; x)   for   x < s ,

(J21)

while ahead the moving crack the contacts are still loaded,

ub(x; s) = ub(0)(x)   for   x > s .

(J22)

With these conditions, the functions u_t(x; s) and u(x; s) take the following form. For the crack tail (x < s):

ut(x; s) = [ D11(s) sinh(κ1x) + D12(s) cosh(κ1x) ] + [ D21(s) sinh(κ2x) + D22(s) cosh(κ2x) ] + ut(0; x) ,

(J23)

u(x; s) = α1[ D11(s) sinh(κ1x) + D12(s) cosh(κ1x) ] + α2[ D21(s) sinh(κ2x) + D22(s) cosh(κ2x) ] + u(0; x) ,

(J24)

where

ut(x0; x) = ⌠ ⌡ x x0   dξ ub(ξ) { b1sinh[κ1(x − ξ)] + b2sinh[κ2(x − ξ)] } ,

(J25)

u(x0; x) = ⌠ ⌡ x x0   dξ ub(ξ) { α1b1sinh[κ1(x − ξ)] + α2b2sinh[κ2(x − ξ)] } ,

(J26)

with b₁= (1−β) κ²/κ₁(α₂− α₁) and b₂= − b₁κ₁/κ₂. Ahead of the crack (x > s):

ut(x; s) = D1(s) e–κ1(x–s) + D2(s) e–κ2(x–s) + ut(s; x) ,

(J27)

u(x; s) = α1D1(s) e–κ1(x–s) + α2D2(s) e–κ2(x–s) + u(s; x) .

(J28)

The coefficients D₁(s), D₂(s), D₁₁(s), D₁₂(s), D₂₁(s) and D₂₂(s) in Eqs.(J23−J28) are functionals of the function u_b⁽⁰⁾(x) and should be determined by the two left-hand-side boundary conditions (J11) and (J12) and the following four continuity conditions at x = s ± 0:

ut(s − 0; s) = ut(s + 0; s) ,

(J29)

ut'(s − 0; s) = ut'(s + 0; s) ,

(J30)

u(s − 0; s) = u(s + 0; s) ,

(J31)

u'(s − 0; s) = u'(s + 0; s) ,

(J32)

which finally lead to a set of integral equations that has to be solved self-consistently.

Numerics. A typical solution of the equations for the first crack passage, for a completely relaxed IL initially (u_b⁽⁰⁾(x) = 0), is shown in Fig.J2.

When the slider's height is much larger than λ_c, the displacement field in the TB is approximately unchanged upon propagation of the crack, even if the displacement in the IL almost double upon crack propagation (Fig.J2a). Because the λ-contacts repin immediately after breaking, they also immediately start to load again as the next λ-contacts relax and the crack propagates. As a consequence, the stress behind the crack is increasing with the distance to the crack (Fig.J2b). Ahead of the crack, the relaxed λ-contacts induce an extra loading of the unbroken contacts ahead of the crack, characterized by a stress peak at the crack location (Fig.J2b). This extra stress Δσ_c is enough to bring initially less stressed contacts up to their threshold. However, because the initial stress is a decaying function of x, the extra stress is sufficient to bring the contacts above their threshold only over a finite length Λ (Fig.J2c).

Analytics. For the first precursor arising from a completely relaxed IL (u_b⁽⁰⁾(x) = 0), we found some analytical results. The full proof is available in SI; here we will only provide the main result and the assumptions made to get it.

We assume that the elastic correlation length is much smaller than the slider's thickness (λ_c<< H). Since λ_c is typically in the micrometer range, this assumption is valid for a large range of tribological systems. We also assume that the interfacial stiffness is much smaller than the internal stiffness of the slider (k << K). For rough interfaces, this assumption is also generally true. Based on the results of Fig.J2a, we assume that the displacement field in the TB is given by u_t(x) ≈ U_t0e^–κ₂x and does not change during crack propagation. Under these assumptions, we find an analytical expression for the length Λ of the first precursor (see Eq.(41) in SI). This expression approximates, for Λ >> κ^–1 (i.e., for large H/λ_c and/or for large k/K), to

Λ ≈ κ2–1 ln{ 2 / [(1 + β − 2κ2/κ) Ψ1] } ,

(J33)

where Ψ₁= 1 + λ_cκ₂ /(1 + λ_cκ) .

Figure J3: Dependence of the characteristic length Λ on model parameters. Symbols are for numerics, solid lines for analytics (see SI), dotted lines for approximate analytics, Eq.(J33). (a) Λ vs H/λc for different values of the interface stiffness k/K = 0.03 (diamonds, red), 0.01 (squares, blue) and 0.003 (circles, black). (b) Λ vs k/K for fixed H/λc= 10 (bottom, red), 25 (middle, blue, dots), 50 (top, black). ν = 0.3,  L/λc= 100.

The analytical results together with the numerical ones are presented in Fig.J3.  Λ is found to increase quasi-linearly with the height of the slider (Fig.J3a) and to increase quasi-logarithmically with the stiffness of the interface (Fig.J3b). Comparison with the numerical results shows a quantitative agreement. We also checked that molecular dynamics simulations of the model agree with the analytical results.

Discussion. We have studied the onset of sliding due to propagation of the self-healing precursor cracks. The characteristic length Λ of the first precursor is controlled by two parameters, H/λ_c and k/K, determined by the slider and interface properties. Importantly, Λ does not depend on the threshold value σ_s. This results from the fact that we use linear elasticity and immediate repinning of λ-contacts in a completely relaxed state. As a consequence, all curves Fig.J2c are unmodified if σ_s is changed, and thus the length Λ is also unchanged.

We have considered a constant value of σ_s along the interface. If the normal load is non-uniform, e.g., σ_n(x) = σ_n0+ εx as in top-driven systems with friction-induced torque, then the thresholds will also depend on x, σ_s(x) = σ_s0+ μ_sεx, and  Λ will depend on the asymmetry parameter ε roughly as Λ(ε) ≈ Λ(0)(1−ε).

We have assumed that broken contacts are immediately restored with zero stretching. Instead, one may assume that the contacts repin after some delay time τ_d. Assuming a constant τ_d, the width of the self-healing crack will increase from λ_c to  λ_c+ vτ_d, with v the (minimal) propagation velocity of the crack, and the length Λ will also increase to Λ + vτ_d.

In the presence of inertia, which has been neglected in the present quasi-static model, the TB blocks will be able to continue increasing the load on their right-hand neighbours even after repinning of their individual λ-contacts. This extra loading is expected to help breaking further contacts and thus increase the precursor length. In contrast, the crack nucleation will be unaffected because it originates from a previously prepared static state of the system.

We have assumed that all λ-contacts behave as classical springs with a sawtooth-shape dependence u(σ). In real systems, they are composed of N_c micro-junctions with a distribution of thresholds P_c(f_s). If P_c(f_s) is wide enough, the elastic instability will disappear, and the λ-contacts will smoothly slide with the pushing velocity v_d after reaching the stretching ~u_s, and no discrete precursor will thus be observed.

Here we have considered the first precursor only. The system is left in the stress state shown in solid line in Fig.J2b. If the system is further loaded, the left-most λ-contact will again be the first to reach its breaking threshold, and a second precursor will nucleate. It will propagate through the stress field left by the previous event, itself modified by the additional loading, and this scenario will repeat itself. We have run numerical simulations for the series of precursors and found that all precursors propagate roughly over the same length Λ. This is not in agreement with experiments as well as with previous models based on phenomenological Amontons-Coulomb law, in which successive precursors propagate over increasing lengths. This discrepancy is resolved, as we checked numerically, if aging of the λ-contacts − newborn contacts have to be characterized by lower breaking thresholds σ_s which then grow with time, e.g., see or − is taken into account. The previously broken region of the interface is characterized by "younger" contacts and therefore by smaller thresholds, than the "old" unbroken part of the system. Hence, the next crack will pass this region more easily and propagate deeper into the system. This process repeats itself until a precursor would span the whole interface and trigger macroscopic sliding of the interface. Note that in the case of aging, it appears a competition between the aging, loading and crack propagation time scales.

Finally, we considered the 1D model of the slider. In this case the stress falls off according to the exponential law, therefore the propagation length (J33) depends logarithmically on the model parameters. In a 3D model instead the stress decays according to a power law; this will modify the expression (J33). Unfortunately, we have no analytics for the 3D model.

Next: Conclusion

Back to main page or tribology page

Last updated on May 5, 2014 by O.Braun.  Copyright © by O.Braun.  Translated from L^AT_EX by T_TH