*2.4. Problem Set Up*

We focus our attention onto a displacement controlled problem under plane strain assumptions in order to highlight the capability of the proposed approach. In the first stage, a displacement linearly increasing with time is applied along the direction normal to the finite layer, up to a given final value of ∆*z*,0 = 2*g*0, reached at time *t* = *t*0, which is then held constant. At this point, a tangential displacement with a constant horizontal velocity is applied to the indenter, which starts sliding. The indenter profile is analytically expressed by:

$$\frac{h(\mathbf{x},t)}{\mathbf{g}\_0} = 1 - \cos\left[\frac{2\pi}{\lambda\_0}(\mathbf{x} - vt)\right] \tag{12}$$

While the velocity of the application of normal load is the same for all the simulations, and assumed to be quasi-static, for what concerns the horizontal load different sliding velocities have been considered in the range *v<sup>i</sup>* = 10 (*i*−10)/3 [m/s], *i* = [1, . . . , 10], with their numerical value being summarised in Table 2.

**Table 2.** Range of horizontal velocities employed.


A regularized Coulomb frictional law [8] is considered, with *f* = 0.2 being the friction coefficient. Figure 4 lists the remaining geometric parameters that describe the problem set, together with the rheological model that is employed for modelling viscoelasticity, which has already been thoroughly discussed in Section 2.3: three different simulations are performed, each of them characterised by one, two, or three terms of a Prony series used for modelling a linear viscoelastic material. The model geometry and applied velocities are the same in all of the cases considered. Finally, periodic boundary conditions have been introduced in correspondence of the two vertical sides of the domain, in order to simulate a semi-indefinite contact in the horizontal direction. The simulations have been performed using the Finite Element Analysis Program FEAP [10], where the MPJR formulation has been implemented as a user element routine. The validation of the proposed computational method is provided in Appendix A.

**Figure 4.** Sketch of the model, *b* = 1, *λ*<sup>0</sup> = *b*, *g*<sup>0</sup> = 5 × 10 <sup>−</sup>4*λ*0.
