Elucidating the Onset of Plasticity in Sliding Contacts Using Differential EBSD Tomography

Depending on the mechanical and thermal energy introduced to a dry sliding interface, the near-surface regions of the mated bodies may undergo plastic deformation. In this work, we use “diﬀerential EBSD tomography” to highlight changes to the microstructure near tribological fcc alloy surfaces, allowing us to detect subtle diﬀerences in lattice orientation and small distances in grain boundary migration. The analysis approach compares electron backscatter diﬀraction images generated from large-scale molecular dynamics simulations with their undeformed counterparts via a simple image analysis ﬁlter. We use our visualization method to discuss the acting microstructural mechanisms in a load- and time-resolved fashion, focusing on sliding conditions that lead to twinning, partial lattice rotation, and grain boundary dominated processes. Extracting and laterally averaging the color saturation value of the generated tomographs allows us to produce quantitative time- and depth-resolved maps that give a good overview of the progress and severity of near-surface deformation. Corresponding maps of the lateral standard deviation in the color saturation show evidence of homogenization processes occurring in the tribologically loaded microstructure, frequently leading to the formation of a well-deﬁned separation between deformed and undeformed regions. When integrated into a computational materials engineering framework, our approach could help optimize material design for tribological and other deformation problems. Abstract


Introduction
Metal interfaces in relative motion are ubiquitous of many relevant machine components, such as bearings, piston rings and gears. The mated metallic surfaces are usually only in contact at the level of individual roughness features, assuming a focus on dry or boundary-lubricated situations [1,2]. The sliding surface is therefore repeatedly exposed to non-homogeneous pressure distributions, with fluctuating stress conditions applying near the contact regions as the asperities of the counterbody pass [3]. The applied stresses can lead to irreversible processes that alter the material's microstructure in the near-surface zone significantly, ultimately resulting in a reduced or increased grain size depending on the loading conditions and the initial microstructure [4,5,6]. These microstructural changes have as a consequence that the mechanical properties of these near-surface zones such as toughness or hardness are altered, which can change the material's reaction to external loads and thus directly affect its friction and wear performance [4,7,8,9]. Tailoring materials and their interacting surfaces to meet specific requirements during their service life is becoming increasingly important to improve the performance or extend the lifetime of machine parts. A comprehensive description of the microstructural changes that occur in the near-surface of metals cannot be described using macroscopic bulk material properties, although their local properties will greatly influence the friction behavior between two surfaces [10,11,12]. In general, for conventional coarse-grained metals, sliding leads to grain refinement and the formation of a well-defined nanocrystalline layer [13,14,15]. Grain coarsening has also been experimentally observed but in connection with ultrafine-grained metals [10]. The reason for the clear and distinct separation between deformed and undeformed regions in sliding metals is still the subject of ongoing research. Single-pass experiments performed on pure copper have shown the formation of a sharp and narrow band of dislocations, dubbed dislocation trace line [16]. The presence of dislocation trace lines has been recently observed independently of the applied load and size of the counterbody [17]. For high applied contact pressures, grain refinement is initiated in the region between the dislocation trace line and the sliding interface. This observation supports the thesis that the formation of an array of dislocations might be the first step that determines the extent of the nanocrystalline layer formed after a higher number of sliding cycles. Surprisingly, the authors found that the use of identical contact pressures applied via the combination of different loads and diameter of the spherical indenters led to radically different microstructural features. The latter is crucial, as contact models relying on conventional continuum approaches may be unable to capture the complexity of the dislocation phenomena occurring beneath a sliding asperity.
The analysis through modeling of the highly dynamic processes occurring in the near-surface zones of a tribologically loaded system requires a simulation approach that can handle the occurring non-homogeneous loads to describe the deformation behavior in the first couple of nanometers from the sliding interface. Progress in high performance computing has made molecular dynamics (MD) simulations and other meshless simulation methods viable tools for studying processes occurring in sliding or abrasive contacts [18,19]. Mesoscopic approaches such as smooth particle hydrodynamics or the material point method are computationally cheap and have already been successfully used to study tribological aspects of scratching [20], milling [21], and grinding [22]. However, their strengths lie mainly in the stable treatment of large deformations, but they fall short of being able to resolve microstructural developments near the mated surfaces. Conversely, it has already been demonstrated in the early 2000s that plastic deformation processes like dislocation generation in nanocrystals can be modeled using polycrystalline MD [23,24,25]. One important advantage of the atomistic approach using accurate interatomic potentials is that it does not require parameters such as activation energies for structural changes like grain growth, or the assumption of constitutive material laws as in continuum mechanics [26]. Polycrystalline MD models featuring tens of millions of atoms have now matured to the level of being able to make predictions about the outcome of nanoscopic sliding processes in terms of the microstructural evolution during plastic deformation [27,28]. These predictions can be qualitatively translated to the micro scale because the simulated grains are sufficiently large to correctly reproduce a realistic material response, thus generating results of surface engineering relevance [29,30].
Changes of the microstructure, especially subtle ones, already occur during running-in and sometimes even after a single sliding pass of the counterbody [17,31,32]. Since it is known that running-in greatly determines the subsequent tribological behavior of the material, it is highly important to understand the microstructural evolution happening during this stage. In particular, for subtle changes in the microstructure it is beneficial to look at the changes in crystallographic orientation/microstructure with respect to the initial state rather than the typical electron backscatter diffraction (EBSD) visualization, which shows the current state of crystallographic orientation. Prasad et al. have impressively demonstrated such a technique to visualize microstructural changes happening in single crystal nickel subjected to sliding [33]. However, this type of visualization is usually almost impossible to carry out experimentally, as the exact initial microstructure cannot be known for a polycrystalline material after the tribological experiment. Using computer simulations, however, the same initial microstructure can be used for varying conditions and material compositions, making it easy to compare the current microstructure of a tribologically loaded surface to its initial state via exact superposition. Additionally, experimental studies at this length scale are costly and time-consuming, which is why computer simulations have become an important means to study and quantify complex phenomena as well as to identify the associated mechanisms [34,35]. Finally, the time span observable with MD is perfectly reasonable to study the phenomena happening during the early stages of the sliding process, which greatly influences its later stages.
[Note to Editor: this paragraph might be moved to a different section of the special issue, if inappropriate here. As this article and the entire special issue are dedicated to the memory of Mark Robbins, this may be the right place to relate a brief appropriate story. I (Stefan J. Eder) was discussing some of the groundwork for this article [32] with him standing in front of my poster at the Gordon Research Conference 2018 in Lewiston, ME, both with a cold New England pale ale in our hands. As his poster was right next to mine, we had each just spent 15 minutes "walking" the other through our work. He had many helpful suggestions, several of which made it into my subsequent work. He mentioned David Rigney, who he thought might get a kick out of this kind of simulation and visualization approach, and suggested that I should just look him up, write him a note and attach the above paper (which had just been accepted). Upon my return to Vienna I wrote to Dave, who answered right away although it was the 4th of July and he was about to turn 80 that year. Two days later, he sent me an annotated list of all his publications with little notes below the papers that I might find useful, which absolutely made my day. Although Mark was a busy man, he never kept me waiting for an answer. Two days before he passed away, I invited him to contribute to a special issue that I am guest editing for an open-access journal. Needless to say, he was the first person to reply, kindly declining because he said he preferred to publish in more classical journals, such as Tribology Letters.] In this contribution, we will first provide a concise overview of our computational model and simulation setup, followed by a description of how to produce, read, and quantify differential EBSD ("diffEBSD") tomographs. We then showcase the capabilities of the analysis approach at the example of 8 representative CuNi systems, varying the alloy composition, the contact pressure, as well as the temperature to sample the relevant regions of deformation mechanism map introduced in [30] and extended in [36]. The last part of the paper focuses on the interpretation of the time-and depth-resolved maps of the diffEBSD color saturation and its standard deviation, allowing us to quantify plastic deformation and its lateral inhomogeneities.

Simulation, Visualization, and Analysis
All our MD simulations were performed using the open-source code LAMMPS [37], the de-facto standard for meshless simulations in academic tribological research. The polycrystalline model of the CuNi base body (see Fig. 1) measures 85 × 85 × 40 nm 3 , contains approximately 25 million atoms, and was prepared as detailed further in Refs. [32,38]. The fcc CuNiX samples, where X is the atomic percentage of Ni (0%, 5%, 25%, 60%, and 100%), interact according to an embedded atom potential with parameters taken from [39]. The system was subjected to heat treatment and cooling to 300 K/600 K as described in [30,36]. Rigid boundary constraints were applied to the lower 3Å of the model, and periodic boundary conditions were applied in the lateral (x and y) directions. The counterbody is a rough bcc Fe monocrystal with a Gaussian root-mean square surface roughness of 0.5 nm, a fractal dimension of 2.186, and a characteristic lateral asperity extent of 33 nm, values that are similar to what has previously been applied [40], in an attempt to keep asperity flank angles low. We reduced the counterbody thickness to several monolayers and kept it rigid to have most of the computational resources available for the polycrystalline CuNi sample, implying that the counterbody is much harder than the sample. Lennard-Jones potentials controlled the interactions between the counterbody and the sample, which implies a third body, with the global energy parameter ε = 0.095 eV obtained as described in [41], while σ Fe−Cu = 0.224 nm and σ Fe−Ni = 0.221 nm were calculated from interatomic first-neighbor distances using the Lorentz-Berthelot mixing rules.
The counterbody was moved across the surface of the sample at a sliding velocity of v x =80 m/s and a small orthogonal component of v y = 9 m/s to prevent roughness features from coming into repeated contact with their own sliding grooves. During the simulation runs, the normal pressure σ z on the sample was kept constant during sliding up to a total sliding distance of 560 nm. The normal pressure was adapted to the composition, the system temperature, and the desired deformation mechanism regime as discussed in depth in [30,36], focusing on parameter sets where (partial) lattice rotation, twinning, and grain refinement dominate.
A Langevin thermostat with a time constant of 0.5 ps acted on all the non-rigid atoms to drain frictional heat from the system. This is equivalent to the phonons in the system being coupled to the electrons, which act as a heat bath to mimic the electronic contribution to the thermal conductivity in a metal [42] in an attempt to approximate the macroscopic heat conductivity. The thermostat acts only perpendicular to the directions of sliding and normal pressure in order not to interfere with these external constraints. In his foundational studies on shearing liquids [43], Mark Robbins observed that if the thermostat acted on degrees of freedom in the direction of shear, artifacts such as shear localization were enhanced. These artifacts are minimized when the thermostat operates on degrees of freedom not directly coupled to the shear force. Mark referred to this as a profile unbiased thermostat.
Many of the computational tomographs of the CuNi alloy system discussed throughout this work are colored according to grain orientations as in electron backscatter diffraction (EBSD), using the inverse pole figure (IPF) coloring standard. The orientations were calculated using polyhedral template matching [44] as implemented in OVITO [45], and the color rendition was carried out using the MTEX toolbox [46,47] for Matlab, with the orientations projected onto the axis parallel to the normal loading direction (z).
In this work, we discuss some of the more subtle microstructural changes in the system that occur during the onset of plasticity, especially migrating grain boundaries and slight lattice rotations. These may be visualized and analyzed using differential EBSD tomography (referred to as "diffEBSD" henceforth), which compares the structure of a given time step with the initial one. In Ref. [48], this approach was introduced in a rather pedestrian fashion using an image processor, which allowed the analysis of only a handful of frames.
Here, we produced some code that can quickly produce a time-resolved set of diffEBSD images from a set of EBSD images, bringing out the full benefit of diffEBSD analysis by allowing the production of movies. We initially filtered out the gray representation of the rigid counterbody from the initial time step which otherwise would lead to confusing artifacts. Then, the difference between the initial and the current time step are calculated as for all three color channels, where A 0 , A in , and A out are the initial, the input, and the output RGB values (0-255), respectively. As can be seen, this filter is symmetrical with respect to A 0 and A in . Care must be taken that the data types of A 0 and A in allow negative values (int16), as most software will by default import RGB data as an unsigned integer (uint8), which leads to saturation artifacts in the difference term.
To be able to export the diffEBSD image, A out has to be reconverted to uint8. The resulting representation of the system is nearly white in all regions where no changes occurred and colorful in regions that changed by a large amount, with pastel hues where changes were slight, see Fig. 2. It is therefore perfectly suited for identifying small deviations in the courses of grain boundaries and especially slight lattice rotations. Once noteworthy changes have been identified using diffEBSD, the respective EBSD images should be revisited for a thorough discussion of the results, as the differential image may obscure which feature belongs to which configuration, and the diffEBSD color scheme is not unambiguous per definition. This can be seen in the bottom right of Fig. 2, where the resulting diffEBSD colors are shown as a trapezoidal matrix as a function of the initial and final EBSD colors aligned with the matrix axes. To extend the benefit from the diffEBSD analysis approach and furnish it with some quantification, we produced time-and depth-resolved maps representing the saturation value of the HSV (hue-saturation-value) color model, S HSV , of the diffEBSD images. This quantity is defined as and S HSV = 0 otherwise. S HSV was laterally averaged over 20 evenly spaced tomographic sections, and the standard deviation of the value distribution was also recorded to retain some information about its lateral variance. The saturation value itself was then analyzed to study the time dependent extent of subtle and stronger microstructural changes into the depth of the sample, while its standard deviation was used to quantify lateral equilibration processes within the evolving microstructure.

Results and Discussion
Out of the large amount of available data, we have selected 8 representative systems with parameter sets that showcase the influence of composition, normal pressure, and temperature on the microstructural development. These influences correspond to the variation of the stacking fault energy, the mechanical driving force, as well as a diffusive component acting on the system, respectively. The selected systems cover the regimes of incipient twinning and lattice rotation as well as grain refinement (as detailed in [30,36]), with one system exhibiting almost no plastic deformation (pure Ni) and two systems bordering the regime that features strong plastic deformation with the formation of a shear layer (highest load / Cu-content). In Fig. 3, we have arranged representative diffEBSD tomographic slices in such a way that time progresses from left to right, and the variation of an individual parameter is denoted by grouping within a colored frame. The system CuNi60 at 300 K and 0.5 GPa, framed in black, is chosen as a base system as it features incipient plastic deformation and can therefore be compared to all other systems. All systems of other alloy compositions, but otherwise identical loading parameters, are grouped together in a green frame. Those where only the normal pressure is varied are surrounded by a blue frame, and a variation of temperature is denoted by a pink frame.
At first sight, the diffEBSD images show the influence of changing composition (green frame), where both the depth and the severity of microstructural change decreases with increasing Ni content. The increase in Ni results in a higher stacking fault energy, which increases the resistance of the polycrystal to plastic deformation, ultimately leading to a negligible plastic deformation as observed for pure Ni at 300 K and 0.5 GPa [49]. The rise in stacking fault energy changes the deformation mechanisms from shear layer formation, via grain refinement and incipient twinning to purely elastic deformation with increasing Ni content [30,36]. For the reference system CuNi60, the diffEBSD images clearly show incipient twinning, partial lattice rotation in the grain to the left, where the upper half of the grain is sheared and rotated, as well as slight grain boundary migration (stronger hues near the grain boundaries). In particular, grain boundary migration would not be apparent in the normal visualization with EBSD-IPF, demonstrating the strength of the diffESBD visualization technique. When going to higher Cu-content or higher loads and, therefore, stronger plastic deformation, the diffEBSD visualization comes to its limits since the images become very busy, and in particular the onset of shearing leads to a highly incongruous superposition of micrographs. This becomes obvious when looking at the composition with the lowest Ni-content (CuNi5) or CuNi60 at the highest load of 1.1 GPa, where it seems that a new grain is developing in the middle until 2-3 ns. However, when looking at the conventional EBSD color scheme in Fig. 4, whose arrangement is identical to that in Fig. 3, it becomes obvious that the upper layer starts shearing to the right, therefore the grains tilt slightly in sliding direction. In the diffEBSD this grain boundary migration and shearing make it seem as if a new grain is developing. Similar observations can be made for CuNi60 at the higher temperature, where a small grain seems to develop in the middle of the slice. Apart from this, twinning, partial lattice rotation, and dislocation activity can be very nicely seen in the upper part of the slice with a clear separation between deformed and undeformed zone [16]. In the diffEBSD it thus becomes clear that merely the part above the separating line rotates while the bottom part maintains its original orientation. Therefore, it is important to stress that diffEBSD is a tool to make subtle microstructural changes visible, but their analysis should always be performed with the aid of the original EBSD-IPF visualization.
Additionally, as was mentioned earlier, the coloring scheme for diffEBSD cannot be unique by definition, so that it is sensible to review the times and points of interest in the unambiguous EBSD-IPF color scheme known to every materials scientist, see Fig. 4. An explicit display of the initial configuration was not necessary, as the pure Ni system in row 4 exhibits such little microstructural change that it can double as a base for comparison. Figure 5 visualizes the trajectories of the orientations of the four marked grains shown in the EBSD panel over time under selected loading conditions. As before, the system CuNi60 at 300 K and 0.5 GPa is central in this comparison, as it is subjected to variations in composition (green axis "nickel content", top), normal pressure (blue axis, bottom), and temperature (red axis, right). Note that only those orientations that feature sufficient population densities were considered in this visualization to keep it as simple as possible.
The systems with varying composition show that lattice rotation is mostly constrained to the three grains marked with white, green, and blue pentagons in the figure being closest to the surface, which has similarly been shown in earlier studies [50]. Thereby, the two smaller grains marked with white and green pentagons are more prone to lattice rotation than the larger grain (blue pentagon), which can be traced back to the dependence of the rotation rate on the grain diameter [51]. The sub-surface grain marked with a red pentagon only shows minor lattice changes for CuNi25 and CuNi60, while it remains entirely unaltered in the pure Ni system. In general, pure Ni shows hardly any microstructural evolution, as already mentioned in the section discussing the micrographs. A decrease in Ni content (CuNi60) leads to a higher activity in the grains marked with white and green symbols, which rotate noticeably, particularly during the first 2-3 ns of the simulation. The surface grain (blue) and the sub-surface one (red) start deforming plastically for this composition, albeit weakly. A further reduction of the Ni content (CuNi25) results in a more pronounced lattice rotation of the two grains marked in white and green when compared to higher Ni contents, with the respective trajectories following more elaborate paths. The major feature of this composition in terms of texture is that the grain marked in blue that remains unchanged for higher Ni contents, shows in this case considerable twinning. However, as observed for higher Ni contents, the grain marked in red hardly deforms as it has no direct contact to the sliding interface. At the highest Cu content (CuNi5), we observe lattice rotations that populate large portions of the orientation triangle, including bifurcations that might be associated with grain refinement (potentially resulting from partial rotation of  grains) as well as intermittent orientations that may arise from constant lattice rotation within a shear layer. In this case, also the sub-surface grain marked with red symbols rotates towards a higher-index orientation. If we allow ourselves to fix the composition (CuNi60) and move along the pressure axis, we see similar trends as those observed for the composition variation. An increase of the contact pressure to 0.9 GPa leads to lattice rotations with trajectories similar to those observed for CuNi25 at 0.5 GPa, with apparent twinning and clear lattice rotations in the grains marked with green and white. Increasing the pressure further to 1.1 GPa leads to a complex set of trajectories reminiscent of those for CuNi5 at 0.5 GPa, with the surface grain marked with white symbols rotating first in the direction of [1 0 0] and then on to [1 1 1] as a shear layer starts evolving.
Finally, when increasing the temperature to 600 K, the map again looks reminiscent of that representing CuNi25 at 300 K, featuring twinning and the following of elaborate orientation trajectories mentioned earlier, with some additional twinning in the grain marked with white symbols. The latter seems to be more particular to this diffusive loading of the system, leading to more dislocation activity inside the grain. Despite similar trajectories shown here independent of load, composition, or temperature is varied, a preferential grain orientation due to mechanical loading cannot be observed.
In Fig. 6, we compare the time-and depth resolved, laterally averaged HSV saturation value of the diffEBSD tomographs for our representative selection of systems. This quantity gives a good overview over when, to which depth into the sample, and how strongly microstructural changes (including lattice rotation and twinning) occur. Note that the pure Ni system was omitted in this comparison, as it did not feature much contrast regardless of time and depth. The lower and upper bounds of the logarithmic color axis were chosen ex post to yield the highest contrast without the need to re-scale between systems. Although the HSV color scheme is periodic (it begins and ends at red), it is clear which occurrences of red denote an HSV saturation of 0.08 or below (lower left) and which ones denote 0.6 or above (top right). As before, the base system for this comparison is CuNi60 at 300 K and 0.5 GPa, which exhibits the weakest microstructural response in this representation. As either the copper content, the normal pressure, or the temperature are increased, we observe that the onset of plasticity occurs earlier. At the highest pressures or Cu contents plasticity can be observed only about 200 ps after sliding starts. Furthermore, the rate at which the plastic deformation propagates into the depth of the sample increases, as does the equilibrium depth of the plastically deformed zone, which may or may not have been reached after 5 ns depending on the composition and boundary conditions. As with the trajectories in orientation space discussed in Fig. 5, pairwise similarities exist between some systems. In order of increasing deformation intensity, CuNi60 at 300 K and 0.7 GPa behaves like CuNi60 at 600 K and 0.5 GPa; CuNi60 at 300 K and 0.9 GPa corresponds to CuNi25 at 300 K and 0.5 GPa; and CuNi60 at 300 K and 1.1 GPa shares similarities with CuNi5 at 300 K and 0.5 GPa. That said, there seems to exist a consistent difference between the onset of plasticity induced by higher pressure or temperature compared to an increased Cu content. In the former two cases, we observe very early slight increases in the HSV saturation value that permeate large portions of the sample, whereas the latter case keeps deeper lying regions undeformed for longer periods of time, but plastic deformation then sets in much faster and to a deeper extent.
After studying the time-and depth development of the mean saturation value, it is also highly instructive to compare the corresponding maps of the lateral standard deviation of the HSV saturation, see Fig. 7. This quantity reflects how homogeneously the plastic deformation is distributed laterally at a given time and depth. As before, pure Ni was excluded from this figure, as it does not exhibit enough plastic deformation to produce any contrast. Again, the color range was adjusted for maximum contrast and is the same for all the maps. While the saturation maps in Fig. 6 generally showed a monotonic increase of the depth to which the system is affected by plastic deformation, possibly saturating to some constant value, these maps feature interesting fluctuations that correspond to dynamic equilibration processes. All maps except those for CuNi60 at 0.5 GPa reflect such behavior, which is basically a lateral homogenization of the plastically deformed zone. This homogenization process is particularly illustrative for CuNi60 at 0.7 GPa. In this case, the HSV saturation value (Fig. 6) exhibits a monotonic growth of the deformation depth with sliding distance/time. However, this microstructural refinement does not proceed homogeneously, as highlighted by the lateral standard deviation of the HSV saturation. The latter features a decrease in value after 3.5 ns, suggesting a homogenization of the deformed region.
To discuss this effect with an example, the HSV saturation for pure Cu at 0.4 GPa is shown in Fig. 8, along with its standard deviation and corresponding representative diffEBSD tomographs. The tomographs are particularly suitable for visualizing the homogenization of the microstructure as after 5 ns sliding, all cross-sections show a well-defined deformed region comprising the 10 nm right below the surface. However, previous time steps reveal that the development of this discontinuity is inhomogeneous. According to the tomographic images, the inhomogeneity is caused by two main factors. One contribution is given by the diverse deformation stages that can be observed in the slices of the polycrystalline aggregate as a consequence of the varying deformation imparted by the rough counterbody. Due to its surface topography, the imparted deformation has a stochastic component that tends to homogenize for sufficiently large sliding distances. However, an inhomogeneous evolution of the deformed region can be also observed within a single tomograph, even though the amount of strain is identical in this case (or similar, if the asperity of the counterbody is not sliding fully parallel to the slice). Let us consider the third slice (at y = 40 nm) as an example. Initially, after 1 ns, only the upper part of one grain is deformed, while the neighboring regions start following after 2 ns. A similar observation can be made in the first slice (y = 6 nm), where until 4 ns, the entire plastic deformation is confined to a single grain. The reason for the inhomogeneous deformation within one slice seems to be related to the different lattice orientations of the grains, as certain orientations will be oriented more favorably for undergoing slip by having slip systems with a higher Schmid factor. As deformation progresses, the previously deformed grains undergo strain hardening, making them progressively more resistant against plastic deformation and allowing the neighboring undeformed grains to catch up with the deformation, which ultimately leads to a well-defined discontinuity with increasing sliding   time. Even though this homogenization of the deformed sub-surface region was observed in a number of testing conditions, it seems to be more evident at intermediate loads with respect to the varying resistance to plastic deformation for the respective compositions (region of grain boundary dominated processes). For instance, the most exemplary cases were Cu and CuNi5 at 0.4 GPa, CuNi25 at 0.5 GPa, CuNi60 at 0.7 GPa, and Ni at 1.5 GPa. At lower normal pressures, where plasticity was confined to the immediate vicinity of the interface, or at high loads, where plastic deformation immediately spread throughout the polycrystalline aggregate, such a homogenization process could not be observed. This irregular formation of the deformed layer seems to contradict the experimental observations made in Ref. [17]. However, the authors found that the dislocation trace line was wavy and rather difficult to identify in some cases. These cases corresponded to those obtained with the highest contact pressure (1953 MPa), but achieved using a combination of low load and small indenter diameter. For higher loads, maintaining the contact pressure with a correspondingly larger indenter diameter or smaller contact pressures, the dislocation trace line was again present and defined below the surface. This could be an indication that particular loading conditions lead to an initial inhomogeneous formation of the discontinuity between the deformed and the undeformed sub-surface region, even though further research would be required to verify this statement.

Conclusion and Outlook
The early time evolution of near-surface microstructure in metal sliding contacts play a significant role in the long-time tribological behavior of bearings, piston rings, and gears. For example, near-surface microstructures that harden to the point of brittleness can spall off and result in enhanced wear. Here, we addressed the question of early time microstructure evolution in sliding contacts through a series of largescale MD simulations in nanoscale polycrystalline fcc metals. Atoms in the sliding counterbody are held fixed, mimicking a much harder slider. In order to understand the influence of process conditions (temperature and normal pressure) and alloy composition, we fixed the initial microstructure for all simulations and correlated the evolving microstructure to this initial microstructure using a novel diffEBSD visual analysis. While EBSD has become a standard tool to quantify a microstructure, diffEBSD compares microstructures at two different times to expose changes in microstructure, such as grain boundary motion, lattice rotation, twin formation, and dislocation plasticity.
What can be learned from this? To illustrate the power of the analysis, we performed simulations of copper/nickel alloys ranging from pure copper to pure nickel. Copper readily displays plasticity and grain reorientation, while adding nickel increases the stacking fault energy, and we observe increasing resistance to plastic deformation. Pure nickel is dominated by elastic deformation. Increasing both temperature and pressure enhances plastic deformation. The extent of the plastically deformed layer is revealed by analyzing the time-and depth-resolved HSV color saturation value of the laterally averaged diffEBSD images. Our analyses serve to provide insight into the initial development stages of the discontinuity that arises between deformed and undeformed regions.
The two-time correlation revealed by the diffEBSD offers an interesting prospect. While both temperature and pressure are process parameters, both alloy composition and initial microstructure are material design parameters. The diffEBSD analysis, when integrated into a computational materials engineering framework, can enable optimal material design for tribological and other deformation problems. Figure 1 Exploded 3D view of the MD system consisting of a plastically deformed fcc CuNi base body and a rigid rough counterbody. The colors represent different grain orientations, and grain and twin boundaries have been emphasized using an edge lter. On how to produce and read a diffEBSD image. The left column shows the initial EBSD0 tomograph and a subsequent one at time t together with the diffEBSD tomograph at time t. The right side shows the color legend, which gives a sense of the non-uniqueness of the diffEBSD color scheme. The diffEBSD colors in the trapezoidal matrix result from applying Eq. (1) to the color bars aligned with the axes. The horizontal color bar corresponds to the path followed by the white arrow in the EBSD-IPF orientation map, spiraling inwards from low-index to higher index orientations. Time-resolved comparison of diffEBSD tomographs. The system CuNi60 at 300 K and 0.5 GPa (black frame) can be compared to all other systems in terms of composition (green frame), pressure (blue frame), or temperature (pink frame), with all other parameters kept constant. Time-resolved comparison of EBSD tomographs. The system CuNi60 at 300 K and 0.5 GPa (black frame) can be compared to all other systems in terms of composition (green frame), pressure (blue frame), or temperature (pink frame), with all other parameters kept constant.  Time-and depth-resolved HSV saturation value of laterally averaged diffEBSD images for a set of representative systems. Note the logarithmic color scale.  Lateral dynamic equilibration of plastic deformation in pure Cu sliding at 300 K and 0.4 GPa. At the top, the timeand depth-resolved values of the diffEBSD HSV saturation value and its standard deviation are mapped. For illustration of the lateral inhomogeneity of the deformation process, ve diffEBSD tomographs, equally spaced along the y-direction perpendicular to normal pressure and sliding, are shown at 1 ns intervals.

Supplementary Files
This is a list of supplementary les associated with this preprint. Click to download.