The roadway was modeled as a three-dimensional rigid body with R3D3 elements. The geometry of the roadway is based on a laser measurement of the asphalt surface. The measured raw data was first filtered and smoothed before the surface was meshed with an edge length of 1

$\mathrm{m}$$\mathrm{m}$. The pavement is fixed firmly in space. As pattern geometries the patterns BB (big block), BBs (big block with sipes) and SB (small block) shown in

Figure 1 were modeled with C3D8R elements. Due to the underintegrated element types the ’hourglass control’ recommended by ABAQUS is used by default. The blocks have a height of

$10.5$ $\mathrm{m}$$\mathrm{m}$, which corresponds to the tread depth of the test tires of

$7.5$ $\mathrm{mm}$ plus the intermediate structure above the steel belt with a height of 3

$\mathrm{mm}$. The tread block was meshed with an edge length of

$0.5$ $\mathrm{mm}$. The convergence of the model was investigated in [

29].

Figure 8 shows the FE model consisting of roadway and block. The concentrated load

${F}_{z}\left(t\right)$ is converted to a surface load and acts on the surface

${A}_{Top}$. The constant slip velocity

${v}_{S}$ is specified for a node in the center of

${A}_{Top}$ as a start and boundary condition. All other nodes of the

${A}_{Top}$ surface are rigidly coupled to it, which corresponds to an infinitesimal flat, rigid and massless slab on

${A}_{Top}$. Both the contact between rubber and road surface and the self-contact of the rubber were implemented with the KINEMATIC formulation, because the use of the PENALTY formulation in combination with the subroutine

vfric led to erroneous results. In this case, the specification of a constant coefficient of friction within the subroutine did not lead to the same result as the specification of a constant coefficient of friction without subroutine. The KINEMATIC formulation leads to the same result with a constant coefficient of friction with and without subroutine. This is again almost identical to the result of the PENALTY formulation without subroutine. Due to the computing time and the strongly non-linear behavior, caused by the large local deformations, an explicit integration method was used. The friction law presented in

Section 3.1 was implemented in the subroutine

vfric. A detailed description of this implementation can be found in

Section 3.2.2.

#### 3.2.2. Local Friction Law

The friction law described in

Section 3.1 was implemented with the Abaqus subroutine

vfric using the programming language

Fortran 77 [

33]. The roadway was defined in contact as master (surface) and the tread block as slave (surface). Based on Equation (

18), the friction value

${\mu}_{{s}_{j}}\left({t}_{i+1}\right)$ is calculated for each node

${n}_{{s}_{j}}$ of the slave surface at each time step

${t}_{i+1}$ (see

Figure 8). The index

j numbers all slave nodes in contact with the road surface. The local surface pressure

${\overline{p}}_{j}$ at the corresponding node is available from the simulation. The time increment

${}_{\Delta}t$ is determined by the step size of the solver. The parameters

$\lambda $,

${\mu}_{0}$ and

${\mu}_{\infty}$ are passed to the program as constants before the simulation starts.

Table 4 shows all vectors used in the implementation in Abaqus. For the node

${n}_{{s}_{j}}$ the calculation rule for the friction coefficient is

Here

${\mu}_{{s}_{j}}\left({t}_{i}\right)$ describes term

$I.b$ from Equation (

18) and thus the part of the coefficient of friction assigned to the rubber surface, since the node

${n}_{{s}_{j}}$ has moved in the current time increment

${}_{\Delta}t$ with the local slip velocity

${v}_{{s}_{j}}\left({t}_{i}\right)$ by the distance

${}_{\Delta}x$ to the position of the node

${n}_{{s}_{j}}$ at time

${t}_{i+1}$. The friction value

${\overline{\mu}}_{{m}_{j}}\left({t}_{i}\right)$ is the friction value on the road at time

${t}_{i}$ at the position of the node

${n}_{{s}_{j}}$ at the end of the current time increment

${t}_{i+1}$. The coefficient of friction on the master surface represents the stationary portion of the fluid film on the road surface. Therefore, the friction coefficient of the master surface is used at the position where the slave node is located at time

${t}_{i+1}$, and thus at the end of the current time increment (see term

$I.a$ in Equation (

18)). By formulating the surface as R3D3 elements, the node

${n}_{{s}_{j}}$ is located on a triangular surface, which is defined by the nodes of the master surface

${n}_{{m}_{k}}\forall k\in {M}_{j}$ as shown in

Figure 10.

${M}_{j}$ describes the set of all master nodes defining the area where the slave node

${n}_{{s}_{j}}$ is located at time

${t}_{i+1}$. The coefficient of friction

${\overline{\mu}}_{{m}_{j}}\left({t}_{i}\right)$ as a description for term

$I.a$ from Equation (

18) can be approximated as the mean value of the coefficients of friction at the three adjacent master nodes and is therefore

A weighting of the values of

${\mu}_{{m}_{k}}\forall k\in {M}_{j}$ corresponding to their distance to

${n}_{j}$ was not done to reduce the complexity of the model.

${\mu}_{ini{t}_{j}}\left({t}_{i}\right)$ thus describes the new coefficient of friction by transporting a part of the viscous residual water film with the slip velocity

${\overrightarrow{v}}_{s}$ and leaving the other part of the residual water film on the road surface. This coefficient of friction now changes due to further squeezing of the water film, which is described by term

$II$ from Equation (

24). For the surface pressure

${\overline{p}}_{j}\left({t}_{i+1}\right)$ at the slave node

${n}_{{s}_{j}}$${p}_{k}={\overline{p}}_{j}\left({t}_{i+1}\right)$ applies and thus for the characteristic time at node

${n}_{{s}_{j}}$
applies. Due to the predefined implementation of the subroutine

vfric by

Abaqus, the surface pressure

${p}_{j}\left({t}_{i+1}\right)$ at the end of the current time increment is used to calculate the friction coefficient gradient using the characteristic diagram

${\mu}_{k}({p}_{k},\tau )$. Physically, this means that first the fluid film is divided between the two surfaces, then the incremental displacement is performed due to the local slip, and finally the squeezing is calculated based on the water film averaged by the incremental displacement. A different order would also be conceivable, but since the model is freely parameterized anyway, this is not of further importance. Thus, the friction coefficient

${\mu}_{{s}_{j}}\left({t}_{i+1}\right)$ is completely defined and is used for the calculation of the friction force at the corresponding slave node. This new friction value will be stored for all slave nodes

${n}_{{s}_{j}}$ until the next time step and will be used again as

${\mu}_{{s}_{j}}\left({t}_{i}\right)$ in the next time step. This data can be stored within the subroutine in the

statev variable provided by

Abaqus. Furthermore, in the next time step a value for

${\mu}_{{m}_{k}}\left({t}_{i}\right)$ is needed again (see Equation (

24)), which is why the calculated value of

${\mu}_{{s}_{j}}\left({t}_{i}\right)$ is also stored proportionally in

${\mu}_{{m}_{k}}\left({t}_{i+1}\right)$. The formulaic description is

where

${S}_{k}$ describes the set of all slave nodes that are located at time

${t}_{i+1}$ on one of the areas defined by

${n}_{{m}_{k}}$. This relationship is shown in

Figure 11. Since there are several slave nodes on the areas defined by

${n}_{{m}_{k}}$, the new friction value

${\mu}_{{m}_{k}}\left({t}_{i+1}\right)$ for a master node is calculated as the average of the friction values

${\mu}_{{s}_{j}}\left({t}_{i+1}\right)\forall j\in {S}_{k}$ of all slave nodes directly adjacent at time

${t}_{i+1}$. This means that at the end of each time increment the friction coefficient vectors

${\overrightarrow{\mu}}_{s}$ and

${\overrightarrow{\mu}}_{m}$ approximate the same course of a real friction coefficient

$\mu (x,y,{t}_{i+1})$ in space. A weighting of the values of

${\mu}_{{s}_{j}}\forall j\in {S}_{k}$ corresponding to their distance to

${m}_{k}$ was omitted analogous to Equation (

25). Alternatively, it would also be possible to use only the four slave nodes that form the segment of the slave surface with which the master node

${n}_{{m}_{k}}$ is in contact. Which of these four nodes are in contact cannot be determined directly from the subroutine. Therefore, the coarser approach of the mean value of all Slave nodes located on an area defined by

${n}_{{m}_{k}}$ is chosen here. The storage of

${\mu}_{{m}_{k}}\left({t}_{i+1}\right)$ is not provided in

Abaqus and is therefore carried out in

Fortran 77 via a storage location for global variables (Common block).

Figure 12 shows the program flow of the subroutine

vfric. First, Equation (

25) is used to calculate the friction value

${\overrightarrow{\overline{\mu}}}_{m}$ on the master interface. The result is used in Equation (

24), part I together with the friction value on the slave surface

${\overrightarrow{\mu}}_{s}\left({t}_{i}\right)$ to calculate the initial friction value

${\overrightarrow{\mu}}_{init}$. Thus, according to Equation (

26) using

${\overrightarrow{\overline{p}}}_{s}$, the characteristic time

$\tau $ can be calculated, which in turn is used in Equation (

7) to read out the gradients

$\partial {\mu}_{k}\phantom{\rule{0.166667em}{0ex}}\xb7\phantom{\rule{0.166667em}{0ex}}\partial {\tau}^{-1}$. This allows the friction value

${\overrightarrow{\mu}}_{s}\left({t}_{i+1}\right)$ to be determined in Equation (

24), which is used by

Abaqus to calculate the friction forces. Within the subroutine, the friction value

${\overrightarrow{\mu}}_{m}\left({t}_{i+1}\right)$ is then stored to the corresponding master node using Equation (

27).

Figure 13 shows the local coefficient of friction on the road surface at the end of simulation for pattern BBs. Where no contact is established, the viscous water film remains undisturbed and the coefficient of friction is equal to

${\mu}_{0}$. On the roughness peaks the water film has been squeezed out and the coefficient of friction

${\mu}_{\infty}$ is present. On flatter elevations, where lower contact pressures have prevailed, the coefficient of friction is between

${\mu}_{0}$ and

${\mu}_{\infty}$. At 50

$\mathrm{k}$$\mathrm{m}$ ${\mathrm{h}}^{-1}$ the tread block was longer in the tire contact area and therefore more time was available to squeeze out the water. This is reflected in the slightly higher overall friction coefficient level. It should be noted that the sliding path covered by the tread block during the simulation is independent of the driving speed, since the sliding speed is directly proportional to the driving speed, which in turn is inversely proportional to the contact time.

In some cases, areas can be seen where the coefficient of friction at 140

$\mathrm{k}$$\mathrm{m}$ ${\mathrm{h}}^{-1}$ is higher than at 50

$\mathrm{k}$$\mathrm{m}$ ${\mathrm{h}}^{-1}$. On the one hand, the overall greater dynamics of the tread block movement at 140

$\mathrm{k}$$\mathrm{m}$ ${\mathrm{h}}^{-1}$ can lead to contact being established briefly in some areas where there would otherwise have been no contact at all between the rubber and the road surface. This can be seen in the top right of

Figure 13b. On the other hand, the material damping in the rubber can lead to locally higher contact pressures and thus to a faster increase in the coefficient of friction. Both effects are negligible for the global friction force.

#### 3.2.3. Wiping Edge

As an extension of the friction law, a so-called wiping edge effect was modeled. The model idea is that an edge of the tread block does not float on the viscous residual water film in front of it, but wipes it away. Since a high coefficient of friction is quickly established under an edge due to the high surface pressure, this is carried along when the edge slips over the road surface. This is illustrated in

Figure 14.

At the time

${t}_{1}$ a certain coefficient of friction is present at the points

${Q}_{1}$ and

${P}_{1}$. At point

${P}_{2}$ there is the initial friction value

${\mu}_{0}$. Up to time

${t}_{2}$${Q}_{1}$ now slides to the right by the distance

${}_{\Delta}x$. Instead of forming the mean value of the friction values at

${Q}_{1}$ and

${P}_{2}$ according to Equation (

24), the following now applies for term I of Equation (

24)

So that we ignore the low coefficient of friction or high viscous residual water film at

${P}_{2}$ and use the higher coefficient of friction of

${Q}_{1}$ instead. Whether it is really an edge is checked in the simulation using two conditions. First, the friction value at the slave node at the beginning of the time step must be higher than on the master surface. The condition for this is

Physically, this means that the edge slides from a dry area of the roadway back onto a wet part of the roadway. Accordingly, the coefficient of friction on the master surface would be lower than on the slave node. However, Equation (

28) adjusts this to the friction coefficient at the slave node and the water film in front of the edge is thus ’wiped away’. The threshold value

${\alpha}_{w}$ ensures that minor irregularities in the friction coefficient curve are not mistakenly classified as an edge. In addition, the condition

must be fulfilled for at least one

$k\in {M}_{j}$. This means that at least one node of the road segment on which the master node

${n}_{{s}_{j}}$ is currently located must have the initial friction value

${\mu}_{0}$, which corresponds to an undisturbed viscous residual water film at this node. If this condition is not fulfilled, all neighboring nodes of the master were already in contact at an earlier time and there is no undisturbed viscous residual water film that could be wiped away from the tread block edge. A restriction of the nodes that could be considered to be edges to certain areas was not carried out due to the sometimes very large tread block deformations. In principle, every node of the tread block can be identified as an edge, as long as the conditions from Equations (

28) and (

30) are fulfilled.

In

Figure 15 all nodes of the roadway are shown on which the wiping edge effect was detected until the end of the simulation. The color scale marks the friction coefficient change caused by the wiping edge effect according to Equation (

28). Apart from a few false detections, the 4 edges of pattern BBs (b), or individual edges of pattern BB (a) are clearly visible. In pattern BBs, the front edge is much more pronounced than the three following edges. This is due to the fact that from a certain point in time the rear edges slip over areas of the carriageway on which the water film has already been squeezed out by the preceding block segment. Thus, Equation (

30) is not fulfilled and no wiping edge effect is detected. Nevertheless, there is a high coefficient of friction in this area, as can be seen in

Figure 13b.