A Layer-Based Model for Frictional Sliding of Pillar Arrays
Jasreen Kaur, Xuemei Xiao, Preetika Karnal, Chung Yuen Hui, Anand Jagota

TL;DR
This paper introduces a new model to simulate and understand friction between interdigitated micropillar surfaces, inspired by biological systems.
Contribution
A novel layer-based physical model is developed to predict friction and deformation mechanisms in micropillar arrays.
Findings
The model accurately predicts shear friction force in sliding experiments.
Deformation mechanisms vary with misorientation and height overlap of the pillars.
Interpillar coupling significantly influences the overall frictional behavior.
Abstract
Bioinspired micropatterned surfaces have been studied as a way to enhance or control interfacial mechanical properties. Here, we study friction for the relative sliding of two interdigitated pillar surfaces. The overall frictional force originates from individual pillar–pillar interactions across the interface as well as the interpillar coupling within each substrate. In this study, we develop a layer-based physical model to simulate the contact and sliding behavior of these pillar arrays. The system is modeled as comprising four layers of nodes, with uniform shear displacement applied to the top layer, while the bottom one is held fixed. Nodes in the inner two layers represent the joints between the pillars and the substrate. The model predicts shear friction force and reveals underlying deformation mechanisms at various misorientations and height overlaps, in good agreement with…
Genes, proteins, chemicals, diseases, species, mutations and cell lines named across the full text — each resolved to its canonical identifier and authoritative record.
Click any figure to enlarge with its caption.
1
2
3
4
5
6
7
8
9
10
11
12- —Division of Civil, Mechanical and Manufacturing Innovation10.13039/100000147
Peer Reviews
No public reviews on file for this paper yet. If you reviewed it on a platform where reviews are public (OpenReview, ICLR, NeurIPS, ICML), you can paste yours below so the community can read it here.
Videos
No videos yet. Explain this paper in a talk, walkthrough, or lecture? Add one.
Taxonomy
TopicsAdhesion, Friction, and Surface Interactions · Force Microscopy Techniques and Applications · Advanced MEMS and NEMS Technologies
Introduction
1
Bioinspired structures have been recognized for their ability to enhance surface mechanical properties. ?−? ? Drawing inspiration from natural designs like the hierarchical features of a gecko’s foot ?,? and the head-arresting systems in? dragonflies, these structures provide switchable adhesion and friction properties. Such mimicry has led to surfaces that significantly improve interfacial mechanical properties, ?−? ? enabling advancements across various fields including soft robotics? and medical devices.?
We use bioinspired pillar structures to study how friction properties can be modulated by altering the design of near-surface architecture at the scale of a few microns. Most previous studies of biomimetic surfaces focus on single-sided patterning where only one surface is structured. ?,? Studies involving patterning on both interfacing surfaces are limited despite their potential to mimic more realistic biological interfaces. Bioinspired pillar-array structures allow us to study the effect of patterns on both sides and how friction can be modulated on such an interface. Shape complementary interfaces, if well-aligned, show strong enhancement of adhesion due to crack-trapping and frictional losses.? Shape complementary surfaces can provide significant increase in adhesion as shown by Vajpayee et al.? for rippled surfaces and Singh et al.? for ridge-channel surfaces.
Guduru et al.? showed that surface waviness in elastic contacts causes unstable detachment and increased energy dissipation, enhancing pull-off force and offering insight into tough adhesion mechanisms. Chen et al.? also showed increase in adhesion in shape complementary pillar surfaces. Velcro? serves as a classic example of how structural complementarity can be harnessed to enhance adhesion.
When two regularly patterned surfaces are brought in contact, these shape complementary structures show spontaneous eruption of Moiré patterns on the interface due to the presence of an orientation or a lattice mismatch.? Moiré patterns arising from an orientation mismatch generate arrays of screw dislocations, whereas those caused by a lattice mismatch produce arrays of edge dislocations. ?−? ? When both lattice and orientation mismatches are present, the resulting Moiré patterns consist of arrays of mixed dislocations. These dislocations, acting as interfacial defects, store elastic energy that is released upon the opening of the interface. As a consequence, Singh et al.? showed that the increase in the density of these dislocations in ridge/channel geometries decreases adhesion. The sliding friction of ridge/channel shape complementary interfaces is likewise strongly enhanced, but this enhancement is strongly attenuated by dislocations.? The relative sliding of the interface is accommodated by the glide of dislocations. For a related but different analysis of adhesion of interdigitated pillar–pillar interfaces, see ref ?
Previous experiments on sliding friction of soft pillar interfaces have demonstrated that friction depends on the normal load, and overall friction arises due to the effect of individual pillar–pillar interaction. ?,? In the present study, we develop a layer-based model that simulates two arrays of pillars coming in contact and sliding relative to each other. This model is based on individual pillar deformation, allowing for pillar–pillar interactions across the interface, as well as on the same side of the interface. The model is in quantitatively good agreement with the overall friction measurements. It also provides insight into how local behavior scales to collective frictional performance, an aspect difficult to isolate experimentally. We also performed finite element (FE) simulations to corroborate our experimental results.
Experimental Materials and Methods
2
Sample Fabrication
2.1
Pillar samples are fabricated by molding the polydimethylsiloxane (PDMS) elastomer into an etched silicon master, where the pillar geometry is defined by photolithography. The PDMS precursor (silicone elastomer base) is combined with the cross-linker (curing agent, Sylgard 184 Silicone Elastomer kit, Dow Corning) in the ratio 10:1 by weight. The mixture is then poured on an etched silicon master and cured at room temperature for 2 days for PDMS to solidify. We also cure samples at different temperatures and curing times that result in slightly different periodic spacings due to thermal shrinkage. The samples used in this manuscript are, in particular, cured at 60 or 110 °C for 2 h. The cured PDMS is then carefully peeled off the silicon master. A complementary sample is fabricated in the same manner. More details on the procedure can be found in previous studies. ?,? Flat controls were also fabricated by using a similar procedure. Figure shows a 3D view of the pillar interface obtained by scanning electron microscopy.
3D SEM image of micropillars arranged in a square array with a minimum interpillar spacing of 20 μm. The pillars are 10 μm in diameter and ∼16 μm in height.
Friction Measurement
2.2
Friction experiments are performed using a custom-built apparatus shown in Figurea. We mounted our samples on top and bottom stages and performed the sliding friction experiments. The typical size of a fabricated sample is 30 mm in length and a width of 10 mm on top and a smaller sample of size 4 × 4 mm for the bottom stage. More details on the experiments can be found here.? The experimental setup comprises two platforms, where the test samples are secured, each equipped with a load sensor to measure forces: one for detecting lateral (shear) forces and another for vertical (normal) forces. The movement along these axes is regulated by independent vertical and horizontal motors, while an additional rotational motor adjusts the lower platform’s angle to introduce misorientation, θ. These motors are integrated with a motion control system, and the entire operation is managed by using custom-developed LabVIEW software. A high-resolution camera captures the real-time behavior of the pillars at the contact interface as sliding occurs.
Experimental setup and definition of misorientation and lattice mismatch. (a) Custom-built tribometer for measuring sliding friction between pillar arrays. The tribometer consists of top and bottom stages where pillar samples are mounted to perform the sliding friction experiment. The bottom sample is kept smaller in comparison to the top to ensure full contact while sliding. The camera on top helps us visualize the interface. (b) Definition of misorientation, θ, which is defined as angular orientation of one sample with respect to the other. (c) Definition of lattice mismatch, λ=a1a2 , which is defined as the ratio of lattice spacing of the two samples in contact.
Frictional sliding tests were performed at five distinct misorientation angles (θ = 0°, 5°, 15°, 30°, and 45°) and for five different normal loads, ranging from 0.075 to 0.4 N. The definitions of the types of mismatches, orientation mismatch, θ, and a lattice mismatch, λ, are illustrated in Figureb,c, respectively. Both patterned and flat control surfaces were tested, with each condition being repeated three times. The system can operate in two modes: maintaining a fixed normal load or controlling normal displacement. In the load-controlled mode, the applied normal force is regulated through normal displacement, which is continuously controlled by a PID controller within the LabVIEW software.
Results and Discussion
3
Moiré Pattern and Dislocations
3.1
When two micropillar samples with complementary shapes are initially brought into contact, Moiré patterns naturally emerge at the interface, as observed in prior studies. ?,? The periodic square arrangement of the pillars results in a structured pattern where light and dark repeating square regions appear in the interface. The size of these regions is influenced by both the angular misalignment, θ, and differences in lattice spacing, λ, which is the ratio of lattice spacing on two sides of the interface.?
Figure shows the interfacial Moiré pattern between two pillar samples in contact with the lattice mismatch, λ = 1.006, and misorientation from θ = 0° to θ = 45°. A rotational misalignment (θ > 0° when λ = 1) leads to the formation of screw dislocations, whereas a mismatch in lattice spacing (λ ≠ 1 when θ = 0°) results in arrays of edge dislocations. When both misalignment and lattice mismatch are present (θ > 0° and λ > 1), dislocations appear with mixed screw and edge character. As illustrated in Figurea, when the two samples are perfectly aligned (θ = 0°, λ ≠ 1), the Moiré pattern only appears because of a lattice mismatch but are not clearly visible because the lattice mismatch is very small (λ = 1.006), leading to a large Moiré period. However, as the angular mismatch increases, the density of these patterns increases systematically, as shown in Figureb,c, and patterns are clearly discernible. Additional cases for λ = 1 and 1.023 are shown in Figures S1 and S3. Figure shows the corresponding images from the coordinates before the sliding simulation starts in the layer-based model for λ = 1.006 and θ = 0°, 5°, and 15°. (model discussed later in Section).
Formation of Moiré patterns and dislocations for λ = 1.00 at the interface of pillar samples on the micron scale. Appearance of the Moiré pattern on pillar structures at several misorientations, θ (a) θ = 0°, (b) θ = 5°, (c) θ = 15°, (d) θ = 30°, and (e) θ = 45°. The black lines in (c) show dislocation lines. As θ increases, the density and periodicity of the Moiré pattern become more pronounced and finer, indicating increased local mismatch and higher density of dislocation-like features.
Representation of Moiré patterns from the layer-based model at angles, θ = 0° to 15° and λ = 1.006°. (a) θ = 0°, (b) θ = 5°, and (c) θ = 15°. (Compare with Figure a–c.)
These Moiré patterns correspond to dislocation arrays. Therefore, an increase in pattern density corresponds to a higher number of dislocations at the interface, as shown in eqs and ? in Section S1.2. We have previously demonstrated? that the overall frictional force in micropillar arrays originates from the interactions between individual pillars. There, we considered the limit where pillar compliance is much larger than substrate compliance so that pillars deform independently of each other, and the total friction force can simply be calculated as the sum of forces between individual interacting pillar pairs. More generally, pillars interact with one another due to substrate elasticity. This interaction can change the friction force and displacement patterns and is included in this work.
Friction
between Pillar Surfaces
3.2
We measured the friction between two pillar samples using a custom-built tribometer. Sliding friction experiments were conducted under different normal forces ranging from 0.075 to 0.4 N and at different misorientations, θ = 0°, 5°, 15°, 30°, and 45°. Forces thus measured were converted to stress by dividing by area of the sample ∼16 mm^2^. Figurea shows shear stress versus displacement for θ = 0°, λ = 1.006, at several normal stresses. It can be seen that shear stress increases with normal stress. Figureb shows typical friction stress versus normal stress for different misorientations in comparison with a control. The shear stress plotted in Figureb is the average value of the shear stress plotted in Figurea during steady sliding, here taken to be the region from 1 to 3 mm sliding distance. The rest of the cases from experiments for λ = 1, 1.023 are shown in Figures S2 and S4.
Sliding friction stress and its average value. (a) Shear stress vs displacement curves at a lattice mismatch of λ = 1.006 and misorientation angle θ = 0°, shown for multiple applied normal stresses. These results highlight how increasing normal load influences the onset and magnitude of shear stress. (b) Friction stress as a function of normal stress at the same mismatch λ = 1.006, plotted for various misorientation angles (θ = 0°, 5°, 15°, 30°, and 45°). Results are compared with those from a flat control sample showing a significant increase due to surface patterning.
In previous work, we performed single pillar-pair experiments on mm scale pillars to visualize the motion of pillars during sliding.? The single pillar, fabricated at the millimeter scale, has a height of 4.8 mm and a diameter of 3 mm, maintaining the same aspect ratio as that of the micropillars. These experiments showed that overall friction force arises from individual pillar–pillar interaction. We measured the progression of sliding and bending of contacting pillar-pairs as they slid apart. (See ref ? which shows the shear and normal force curves for single pillar-pair experiments.) These experiments were performed at several diametric overlaps ( 1, 0.75, 0.5, 0.25, and 0) and various heights of contact (h c) ranging from 4.8 mm to 0.8 mm. Here, is normalized lateral offset and d O refers to the lateral offset (as shown in Figuree) between two cylindrical pillars, with 1 indicating full alignment (maximum overlap) and corresponding to incipient contact with no overlap. The height of contact, h c (Figuref), represents the vertical height of the overlap between the interacting pillars. We performed quartic fits to shear force versus displacement experiments for different heights of contact. We also showed corresponding normal forces for different heights of contact with quartic fits. The parameters from these fits were then used in our simulations to model pillar interaction forces. These interpillar interaction measurements are combined with the layer-based model discussed below to calculate friction stress. This approach was necessary because the pillars are quite stocky and deformations are large. This made classical beam bending theory to model pillar–pillar interactions inapplicable to our experiments.
Schematic drawings of the samples and definition of variables. (a) Schematic of the cross-sectional view of pillar arrays showing pillars in contact. Defining geometry of pillar arrays in contact consisting of layers I, II, III, and IV. Nodes in the top layer, I, are given incremental displacements in the “y” direction. Nodes in layers II and III move because of unbalanced forces, which stops when equilibrium is achieved. The procedure is repeated incrementally. (b) Top view of layers I and II, given a uniform sliding displacement, Δy, in layer I. (c) 3D rendering of the upper part of the specimen (left) showing labels for node location and a pair of pillars in contact (right) with x–y locations described at each node. (d) Pillar–pillar contact defining shear displacement in the deformed state (left). Pillar contact defining force observed on the bending pillar (right). (e) Depiction of the diametric overlap (δ x ) between two pillars in contact and d refers to the distance between the centers of interacting pillars. δ y refers to the relative sliding displacement in the shear direction, y. (f) Description of height overlap (h c) for interacting pillars and 2L denotes the gap between surfaces of interacting pillars.
Layer-Based Model
3.3
Our previous work showed that when pillar–pillar contact compliance dominates substrate compliance, the total friction force arises from summing the contributions of individual pillar interactions during sliding. To model more complex scenarios, we have developed a layer-based model, as described below. This model allows for greater flexibility and accounts for substrate compliance, and selecting variables such as pillar height, diameter, density, and substrate elasticityparameters that are challenging to manipulate in experimental settings.
The layer-based model represents the pillar interface as a four-layer system of nodes (Figurea). The four layers are defined as follows: layer I corresponds to the top surface of the upper substrate, layer II represents the base of the pillars attached to the upper substrate; layer III represents the base of the pillars on the lower substrate, and layer IV corresponds to the bottom surface of the lower substrate, which is fixed.
The interface is parallel to the x–y plane. The interface bears and transmits shear load in this x–y plane and a normal load in the “z” direction. Figurea shows a cross-sectional view of the interface with a uniform shear displacement (in the x–y plane) applied to the top layer (layer I). Figureb provides a top view of the interface, showing nodes in layers I and II. Each node in layer II (smaller dark blue circles) lies at the center of the pillar, where it joins its substrate. For each node in layer II, there is a corresponding node in layer I (larger light blue circles). In the absence of any shear deformation, the “x–y” coordinates of a node in layer II and the corresponding node in layer I are identical. The zoomed-in region of Figureb highlights the relative y displacement between layers I and II, denoted as Δu _ iy _ ^II^, that arises when a motion, say Δy, is applied to layer I. Figurec (left) is a 3D rendering of the top part of the specimen showing layers I and II. The top part interacts with a similar lower part (layers III and IV). Figurec (right) shows a configuration of two interacting pillars, one from each part. The x–y location of node “i” in layer I of the configuration is denoted by (x _ i _ ^I^,y _ i _ ^I^) and so-on for the other layers. Stiffness coefficients of the substrate are denoted as “k s” with additional subscripts as needed. Figured (left) represents bending of pillars in contact, showing deformation in upper and lower substrates as Δu _ ix _ ^II^, Δu _ iy _ ^II^, Δu _ jx _ ^III^, and Δu _ jy _ ^III^. Figured (right) shows bending of pillars in contact and substrate force and beam force as observed during bending.
Starting with an equilibrium state, we apply a uniform incremental displacement to the top layer I in the sliding or “y” direction. This sets up unbalanced forces on the contact layers II and III, which initiates motion of their nodes, and the system relaxes to a new equilibrium state. The net force on nodes in layers II and III has two contributions: a substrate force, F s, due to substrate deformation and a beam force, F b, due to the bending of pillars (Figured). The substrate force, F s, is the restoring force exerted on nodes in contact layers (layers II and III) due to their elastic coupling to layers I and IV, respectively, and to other nodes within their own layer. The beam force arises due to pillar–pillar contact and deformation as detailed later. After equilibrium, this process is repeated with an incremental increase in displacement applied to layer I. To achieve equilibrium, for each increment of displacement of layer I, we solve the following equation of motion until the right-hand side of eq is smaller than a preset value
where (F) is the total force vector, (F s) is a vector of forces from the substrate, (F b) is a vector of forces due to bending, both functions of (u), which is the vector of nodal displacements, η represents a viscous damping coefficient, which modulates the rate of change of displacement, and describes the velocity of nodes. Note that the force and displacement vectors have dimensions of 2(N + M) × 1, where N and M are the number of nodes in layers 2 and 3, respectively (i.e., we balance forces in the “x” and “y” directions on each node: see eqs S23 and S24 for details).
Equation is a numerical tool to iteratively relax the system toward force balance. The left-hand side represents a damping term that gradually diminishes as the system evolves. At each increment, the solution is updated to reduce the net residual force, and equilibrium is considered achieved when the average of forces on the right-hand side reduces below a predefined fraction of a characteristic force value.
This approach is particularly advantageous in handling discontinuities that arise from sudden changes in contact conditions, such as when pillar–pillar contacts break or form, where traditional quasi-static solvers may fail or exhibit numerical instability.
We next relate forces to the nodal displacements. Since the substrate is much thicker (∼1 mm) than the characteristic in-plane distances between pillars (∼20 μm), it is appropriate to approximate both the top and bottom substrates as elastic half-spaces connected to their respective pillars, each pillar being represented by a node. Consider node “i” in layer II (Figurea). Let the shear displacements of this node relative to layer I be Δu _ ix _ ^II^ and Δu _ iy _ ^II^, where the superscript denotes the layer number and the subscript indicates the node number and direction (x or y). Displacements in x and y directions on node i and between layers I and II are defined as
These shear displacements are due to the force on node i and a sum of influences on node i from all other nodes j on substrate II. This relationship is expressed using stiffness coefficients derived from elasticity theory for an elastic half-space (details in Section S2)
Δu _ ix _ ^II^ is the total x-direction displacement on layer II at node i, where Δu _ si,x _ ^II^ is the x-component of the substrate force at pillar i of layer II, F _ sj,x _ ^II^ is the substrate force at pillar j in the x-direction, F _ sj,y _ ^II^ is the force in the y-direction at pillar j of layer II, and ** N ** is the total number of nodes in layer II. Similarly, Δu _ iy _ ^II^ is the total y-direction displacement on layer II at node i.
Substrate shear stiffness (k _ sii _) relates force on pillar i to its contribution to displacement and is defined from ref ?
where G is the shear modulus, R is the pillar radius, and v is Poisson’s ratio. The stiffnesses k _ sx,ij _ and k _ sy,ij _ are likewise derived from ref ?
k _ sx,ij _ captures how much x-force pillar i experiences due to a x-displacement at neighboring pillar j, where Δx _ ij _ ^II^ = x _ i _ ^II^ – x _ j _ ^II^ and , where Δr _ ij _ ^II^ is the distance between centers of pillars i and j. In the results presented later, Poisson’s ratio, for the elastomer used in experiments, ν = 0.5 and shear modulus G = 0.67 MPa. Similarly,
A third substrate stiffness, k _ fxy,ij _, quantifies x–y coupling
These stiffnesses encode how forces on the pillar base relate to its elastic deformation. The same process is followed for the lower part of the specimen with layer III playing the role of layer II and layer IV playing the role of layer I; in this case, layer IV is held fixed. Details of this derivation are described in Section S2.
These expressions can be compactly written in matrix form to describe the displacement of all of the pillars in the system. Letting (Δu) be the vector of displacements and (F s) the vector of substrate forces, we write
[C]s is the substrate compliance matrix, consisting of self-compliances (diagonal) and interaction compliances (off-diagonal) (see eq S24). The equation above applies to layer II. A similar equation applies to any pillar j on layer III to calculate substrate forces on layer III defined as
where “k” refers to a neighboring pillar on layer III and M is the number of nodes on layer III. Sliding experiments and simulations are generally conducted with a smaller upper part sliding on a larger lower part and so M need not equal N.
Next, we move on to the calculation of pillar interaction forces, which occur due to the interaction between pillars on layers II and III. These forces are caused by the relative sliding motion and contact between individual pillars during bending. The magnitude of beam force depends primarily on the overlaps, δ_ x , δ y , and height of contact (h c) between each interacting pillar pair (as shown in Figuree,f). These two geometric parameters determine the extent of pillar deformation, with larger overlaps resulting in higher restoring pillar interaction forces, as described previously.? To identify interacting pairs, the model calculates the center-to-center distance, d, between each pillar on layer II and all pillars on layer III. The pillar interaction is considered active when d ≤ 2R, where R is the pillar radius. In Figuree, blue circles represent the undeformed cross section of pillar i in the top layer (layer II), and the orange circle represents the undeformed cross-section of pillar j on the bottom layer (layer III). The blue circle moves in the “y” direction, and the lowest one represents incipient contact between pillars i and j. The next higher blue circle represents the case where pillar–pillar overlap is maximal, which we define as the diametral offset. (The two pillars will never make contact if the diametral offset exceeds the pillar diameter.) The variable δ y _ is the shear displacement between the two interacting pillars in the sliding direction. This is the quantity that is fed into the empirical quartic force–displacement relation derived from single-fiber shear tests. This relation describes how much lateral force a pillar–pillar interaction generates as they slide past each other, and δ_ y _ changes accordingly.
The model uses empirical quartic polynomial fits with coefficients (p_1_, p_2_, p_3_, and p_4_) derived from single-pillar-pair shear experiments to relate pillar interaction forces to displacements, δ_ y _ and δ_ x _;? see Figuree. The pillar interaction force in the shear direction, y, on pillar i is calculated as
where . This force is turned off when δ_ x _ > R, the condition under which two contacting pillars slip (see FE results in the next section).
The “x” force on pillar i is calculated as
where δ_ x _ = 2R – |x _ i _ ^II^ – x _ j _ ^III^| is the overlap in the “x” direction.
These fits were originally obtained from millimeter-scale single-pillar-pair experiments. To apply them at the microscale of the simulation, we rescale the displacements and forces under the assumption that stress remains scale-invariant. Specifically, if the geometric scale factor between mm-scale experiments and multipillar micron-scale experiments is α, we scale down our displacements by α and forces by α^2^. The pillars used in single-pillar experiments have a diameter of 3 mm and a height of 4.8 mm, while the microscale pillars have a diameter of 10 μm and a height of 16 μm. The geometric scaling factor α = 300. This scaling allows us to conduct simulations at the mm scale and to scale forces and displacements to the micron scale. By the same token, we can state that for a given geometry, stresses increase with shear modulus. This is tested by direct FE modeling; see Section S4.4.
The total force in eq, (F), is column vector for “x” and “y” components at each node and can be written as a concatenation of forces on nodes in layers II and III. (Nodes in layers I and IV are fixed.) For layer II, we write
and for layer III
Once forces at nodes in layers II and III are known, for each increment, we solve eq using simple Euler forward integration in time until equilibrium is achieved. We validate our numerical approach by comparing with independently calculated solutions for some simple cases.
Simulation
Results
3.4
To compute predicted friction forces, we chose a system of ∼1700 pillars sliding over ∼3200 pillars, a system large enough to capture interfacial dislocation arrays and result in steady-state overall shear stress while sliding. This lattice size was selected because the shear stress curves converge beyond this scale, indicating that larger systems do not significantly alter the results shown in Figurea–c.
Friction simulation results for λ = 1.006 at various misorientation angles θ = 0°, 5°, and 15° and five different normal loads. (a) Shear stress vs shear displacement at θ = 0°. (b) Shear stress vs shear displacement at θ = 5°. (c) Shear stress vs shear displacement at θ = 15°.
Figure presents the simulated shear stress versus displacement curves for several lattice misorientation angles. As shown in Figurea (0°), the shear stress exhibits a periodic oscillation arising from the registry between the two aligned lattices. However, this oscillatory response progressively diminishes with increasing misorientation: for 5° and 15° (Figureb,c), the curves become increasingly smooth. In all cases, the shear stress is computed by summing the lateral force on each top-layer pillar at every time step and dividing it by the total area of the top lattice.
We orient one of the lattices at a misorientation, θ, with respect to center to generate Moiré patterns or dislocations on the interface as shown in Figurea–c. Full simulation videos corresponding to cases in Figurea–c are shown in the Supporting Information (Section S6), namely: SV1, SV2, and SV3. The rest of the figures for cases of λ = 1, 1.023 are shown in Figures S5 and S7. Figure shows the spatial maps of displacement between interacting pillars across top and bottom layers before sliding at various rotational offsets (0°, 5°, and 15°). Videos SV1, SV2, and SV3 show simulations corresponding to the three misorientation angles. When sliding begins, the lattice mismatch or an orientation mismatch gives rise to interfacial features resembling Moiré patterns. As the rotation angle increases, the pattern transitions from a uniform grid (0°) to Moiré-like interference structures (5°–45°) as also shown in Figure, reflecting how relative alignment, θ, affects the interface. As shown in Figures and ?, the density of Moiré patterns increases with an increase in θ. Additional simulation videos for a smaller lattice size (∼200 pillars sliding over 700 pillars) for θ = 0° to 45° are available as Videos SV4, S5, S6, S7, S8, and S9 in the Supporting Information.
In our previous modeling framework,? pillar interactions were treated as isolated pairwise contacts without accounting for the underlying substrate’s compliance. This simplification assumed that each pillar deformed independently, neglecting the mechanical coupling introduced by the shared substrate. As a result, shear stress versus displacement curves for θ = 0° exhibited a prominent sawtooth curve-like behavior, reflecting abrupt transitions in contact states between individual pillars. In contrast, the current simulation model incorporates substrate compliance, recognizing that all of the pillars are elastically coupled through a continuous substrate. This allows long-range stress redistribution and coordinated deformation across the pillar array. Consequently, for the same θ = 0° alignment case, the characteristic sawtooth profile observed previously transitions to a smoother force response as shown in Figurea, reflecting the collective influence of pillar and substrate-mediated interactions. Figureb,c refers to simulation results for θ = 5° and θ = 15° respectively. As shown, the shear stress becomes smoother as the misorientation increases, reflecting similar data from experiments. Note that the simulation (Figure) does not capture the static friction peak associated with the initiation of sliding (Figurea). The static friction peak is associated with different processes than sliding or dynamic friction; it often involves fracture or adhesion.? Since our focus in this work is on sliding friction, we have not included adhesion in our models and consequently do not capture the static friction peak. We only show some of the cases here, and the rest of the cases of θ = 30°, 45°, and λ = 1, 1.023 for all θs are shown in Figures S6, S8, and S10.
Normal coupling between pillars is neglected based on the assumption that vertical compliance is localized beneath each pillar and that normal forces can be determined independently via local force balance. Since the model focuses on shear deformation with a one-way coupling to normal force, accounting for cross-pillar normal interactions is not required. In contrast, shear deformation involves a distributed substrate strain and must include coupling between pillars.
In a fully coupled model, both shear and normal displacements interact across and within layers. In our approach, we primarily model shear deformation and compute normal forces independently for each pillar through local force balance. While normal force coupling between pillars on the same layer is not explicitly included, it is indirectly captured because the normal force depends on the shear forces, which are coupled through the substrate.
We proceeded to calculate the normal force. Until this point, our analysis focused on shear forces and displacements in the x–y plane. In a full model, one would introduce additional kinematic variables such as “z” displacement and additional forces such as the “z” component of force. Here, we simplify the model to be in the x–y plane and provide one-way coupling via vertical force balance. In particular, this approach obviates the need to account for the z-displacement coupling between pillars.
We begin by examining a representative interaction between two pillars during sliding, as illustrated in Figure. Upon initial contact, the pillars start bending to accommodate the imposed sliding displacement. This bending persists until reaching a mechanical limit, beyond which the pillars can no longer sustain the deformation solely through bending. At this threshold, the pillars initiate sliding apart. During the bending process, both normal and shear forces arise concurrently. To understand their interplay, we employ a force balance approach that couples these two forces. As sliding initiates, assuming point contact, the pillars bend by an angle α, defined by the relationship: , where δ_ y _ represents the shear displacement and 2L denotes the effective contact height between the pillars. The shear force (F) is defined as the overall shear force observed on a single pillar in the x and y directions. Here, we assume Coulombic friction, and the total shear force F is resisted by the tangential component of the contact normal force (C N) and normal component of τ, where τ = μC N, and μ is the coefficient of friction between pillars. The normal contact force (C N) decomposes into components as C _N_cos α and C _N_sin α in horizontal and vertical directions, respectively. τ also decomposes into τ sin θ and τ cos θ. The shear force can be written as F = C N cos α + μC N sin α. This gives, F = C N(cos α + μ sin α), and thus, . The total normal force on a single pillar (N) can be defined as a function of shear force. From the force balance, we get normal force as N = C N(sin α – μ cos α). Substituting C N to get normal force from given shear force as
Free body diagram illustrating the interaction of two contacting pillars undergoing lateral sliding. Upon contact, bending induces both a shear force and a normal force. The angle of bending, α, reflects the geometry via tanα=δy2l , where δ y is the shear displacement and 2L is the effective contact height.
This final expression clearly illustrates how the normal force acting on pillars systematically depends on the known shear force measurements, geometry (α, δ_ y _, l), and frictional parameter (μ). It underscores the interdependence of normal and shear forces.
Figurea,b shows a snapshot of the interface during sliding friction experiment and simulation at θ = 5°. Figurea corresponds to Moiré patterns as seen in micropillar experiments, and Figureb corresponds to Moiré patterns as they appear in our layer-based model simulation. The color scale represents the amount of relative shear displacement between pillar pairs: black indicates maximum displacement in sliding direction, while lighter shades of red indicate lesser relative displacement between the two interacting pillars. Sliding friction simulation is performed for 1700 pillars sliding over a lattice of 3200 pillars. The normal stress is computed as the sum of normal loads on all pillars on the upper part, divided by the interfacial area. Figurec shows the friction stress versus normal stress from simulations at θ = 0° to 45° at λ = 1.006 in comparison to stress data from the experiments. As shown in Figurec, our model results compare well with those of our experiments.
Comparison of experimental and simulated images of the interface and friction stress for varying normal stress. (a) Interfacial dislocations as seen in experiments during sliding. (b) Simulated sliding of a lattice containing 1700 pillars over 3200 pillars. The color bar refers to the local value of pillar deflection. (c) Friction stress vs normal stress from simulation in comparison to our stress data from the experiments.
FEA Simulation of Pillar
Layer Sliding
3.5
We perform FE analysis of pillar sliding to capture the overall sliding and forces and to validate the layer-based model using ABAQUS 2023 (Dassault Systèmes, Providence, RI, USA). As shown in Figure, the model consists of two parts: a bottom part containing 25 pillars (yellow) arranged in a square array with a center-to-center spacing of a = 4R = 6 mm along two orthogonal in-plane directions and a top part containing a single pillar (green) located at the center. All pillars have a diameter of D = 2R = 3 mm and a height of h = 4.8 mm. The bottom part includes a substrate measuring 30 × 30 mm, while the substrate of the top part is smaller, 12 × 12 mm, to reduce computational cost. Both parts have a uniform thickness of T = 12 mm in the vertical direction. Figurea shows the top view of the model geometry and Figureb shows the front view. Interaction between the top and bottom parts is defined using a friction coefficient of 0.4. To prevent interpenetration between pillars on the same part and between pillars and the substrate surface, we apply hard contact.
Geometry of the simulation model. (a) Top view and (b) front view. The model consists of a top part with one central pillar (green) and a bottom part with 25 periodically arranged pillars (yellow). All pillars have diameter D = 3 mm and spacing a = 6 mm. Both substrates have a thickness of T = 12 mm. The bottom part is fixed at its bottom; the top part moves diagonally while being constrained in height.
At the start of the simulation, the center of the top part is positioned at the edge of the bottom part with a height overlap between the two sets of pillars defined as h c (Figureb). The pillars are initially separated and not in contact. A displacement-controlled motion is applied to the top surface of the top part, moving it diagonally in the plane at a velocity of , while keeping its height fixed. The bottom part is fully fixed at its bottom surface. Both parts are modeled as incompressible neo-Hookean materials, with a shear modulus of μ_0_ = 1 MPa and a mass density of ρ = 965 kg/m^3^.? The simulations are performed using the dynamic implicit, quasi-static solver in ABAQUS, meaning that inertial effects are neglected; the density is included primarily to regularize the problem and to enable the simulation of unstable events, which cannot be captured using a purely static solver. Incompressibility is enforced by setting bulk modulus K = 2000 MPa. Following our previous work and to match the experimental results, the shear modulus is scaled to μ = 0.65 MPa. Accordingly, the forces used for comparison with the layer-based model results are obtained by multiplying the original FE output by a scaling factor of 0.65.
All elements are of type C3D8RH. The mesh is shown in Figure, where Figurea presents the mesh for the bottom part and Figureb presents the mesh for the top part. Note that the scale is not consistent between the two views; Figureb is enlarged to show mesh details. For h c = 2.4 mm and h c = 3.6 mm, the mesh size on the pillars (top and bottom parts) is uniformly 0.2 mm along both radial and axial directions. In the in-plane substrate regions, the mesh size increases from 0.2 mm at the pillar edge to 0.8 mm at the substrate edge for the top part and to 0.4 mm at the substrate edge for the bottom part. For the top part, we added an additional square partition centered on the part, with a side height of 4R and applied a mesh size of 0.4 mm on this square to facilitate a smoother transition from the fine mesh near the pillar to the coarser mesh toward the edge of the substrate. Along the thickness direction of the substrates, the mesh is refined to 0.2 mm within the first 0.6 mm from the interface with pillars, then gradually coarsens to 0.8 mm over the next 1.8 mm, and remains at 0.8 mm through the remainder of the domain. For the case of h c = 1.2 mm, a finer mesh is used for the single pillar in the top part: 0.1 mm along the axial direction and radially refined from 0.2 mm at the center to 0.1 mm near the edge. All other mesh settings remain the same as in the h c = 2.4 mm and h c = 3.6 mm cases. A mesh convergence test was performed by using a smaller system consisting of only two diagonally positioned pillars on the bottom part (with a correspondingly smaller substrate size). Further details are provided in the Supporting Information, Figures S17–S19.
FE mesh of the simulation model. (a) Mesh of the bottom part with 25 pillars; edges of pillars are highlighted in yellow circles. (b) Mesh of the top part with one central pillar highlighted in green. Note that the scale differs between (a,b); the top part is enlarged to show mesh details.
The simulation results of shear and normal reaction forces applied to the top surface of the top part are shown in Figurea,b, respectively, as a comparison with results from our layer-based simulation for a pillar sliding over the 5 × 5 matrix of bottom pillars. Each plot includes data for three cases: h c = 1.2 mm, h c = 2.4 mm, and h c = 3.6 mm. All reaction forces are scaled by a factor of 0.65 to account for the adjusted shear modulus. The validation linear scaling is confirmed in Figure S19. Spikes appear in the force curves when the top pillar transitions from contact with one bottom pillar to the next, likely due to the high contact stiffness introduced by the penalty method. The raw simulation data, including these spikes, is provided in Figures S13–S15 of the Supporting Information. The strong agreement in normal force results validates the reliability of our layer-based model in capturing the dependence of normal force on shear deformation, supporting its use in predictive studies.
Comparison of results from FEA and the layer-based model. (a) Shear force vs shear displacement for a configuration with one top pillar and 25 bottom pillars, shown for three different contact heights. (b) Normal force vs shear displacement for the same configuration and contact height conditions.
Summary and Conclusion
4
This study focuses on understanding sliding friction in microscale pillar arrays using a layer-based model. This model simulates interactions between top and bottom pillar arrays by calculating forces based on their proximity. The model captures the origin of shear friction as arising from both pillar–pillar interactions at the sliding interface and elastic coupling between pillars within the substrate. The system is represented as four layers of nodes, where displacement is applied to the top layer and deformation propagates through substrate and beam-like interactions across the interface. By accounting for misorientation and height overlap between pillars, the model successfully predicts friction forces and deformation patterns that align well with experimental observations from multifiber sliding tests. Our model successfully captures frictional forces for an array of sliding pillars and helps us visualize the Moiré patterns while sliding.
Supplementary Material
The reference list from the paper itself. Each links out to its DOI / PubMed record.
- 1Scherge, M. ; Gorb, S. Biological Micro- and Nanotribology: Nature’s Solutions (Nanoscience and Technology); Springer, Berlin; New York, 2001; p 304.
- 2Arzt E.Gorb S.Spolenak R.″From micro to nano contacts in biological attachment devices,″Proc. Natl. Acad. Sci. U.S.A.200310019106031060610.1073/pnas.153470110012960386 PMC 196850 · doi ↗ · pubmed ↗
- 3Bhushan B.″Biomimetics: lessons from nature-an overview,″Philos. Trans. R. Soc., A 200936718931445148610.1098/rsta.2009.001119324719 · doi ↗ · pubmed ↗
- 4Autumn K.″Gecko Adhesion: Structure, Function, and Applications,″MRS Bull.200732647347810.1557/mrs 2007.80 · doi ↗
- 5Autumn K.″Evidence for van der Waals adhesion in gecko setae,″Proc. Natl. Acad. Sci. U.S.A.20029919122521225610.1073/pnas.19225279912198184 PMC 129431 · doi ↗ · pubmed ↗
- 6Gorb S. N.″Evolution of the dragonfly head-arresting system,″Proc. R. Soc. London, Ser. B 1999266141852553510.1098/rspb.1999.0668 · doi ↗
- 7Arzt E.″Biological and artificial attachment devices: Lessons for materials scientists from flies and geckos,″Mater. Sci. Eng., C 20062681245125010.1016/j.msec.2005.08.033 · doi ↗
- 8Glassmaker N. J.Jagota A.Hui C. Y.Kim J.″Design of biomimetic fibrillar interfaces: 1. Making contact,″J. R. Soc. Interface 200411233310.1098/rsif.2004.000416849150 PMC 1618938 · doi ↗ · pubmed ↗
