Christoph
Klein
ab,
Christopher R.
Iacovella
ab,
Clare
McCabe
abc and
Peter T.
Cummings
*ab
aDepartment of Chemical and Biomolecular Engineering, Vanderbilt University, Nashville, Tennessee 37235, USA. E-mail: peter.cummings@vanderbilt.edu
bVanderbilt Facility for Multiscale Modeling and Simulation (MuMS), Vanderbilt University, Nashville, Tennessee 37235, USA
cDepartment of Chemistry, Vanderbilt University, Nashville, Tennessee 37235, USA
First published on 9th March 2015
The tribology of 2-methacryloyloxyethyl phosphorylcholine monolayers in water is studied using molecular dynamics simulations. Our results show two distinct shear regimes where the first is dominated by hydration lubrication, exhibiting near zero friction coefficients, and the second by chain–chain interactions, resembling monomer-supported lubrication. These results provide insight into the hydration lubrication mechanism – a phenomena thought to underlie the extremely efficient lubrication provided by surfaces functionalized with polyzwitterionic polymer brushes and the mammalian synovial joint.
Previous studies have demonstrated the utility of molecular simulation in elucidating tribological mechanisms19–23 but have not specifically focused on hydration lubrication. Other simulation studies have investigated both neutral24 and zwitterionic25–27 monolayer systems in aqueous environments. However, these simulation studies considered densely populated monolayers that do not compare to the environment existing within the experimentally studied brush systems thought to exhibit hydration lubrication. Moreover the foci of the studies considering zwitterionic monolayers were not centered directly on the behavior within zwitterionic brush films but rather on dense phospholipid and carboxybetaine coated surfaces and their interactions with solvated ions.
In this work, atomistic molecular dynamics simulations are used to study the frictional properties of sparse, silane-functionalized MPC (Fig. 1A) monolayers in an aqueous environment. More specifically, the aim is not to directly model the behavior of pMPC brush films but rather gain a deeper understanding of the molecular level environment surrounding individual monomers found in the brushes. Additionally, the monolayers studied herein differ from previously studied monolayers in two primary ways. First, the zwitterionic nature of the monomers clearly distinguishes them from neutral monolayers such as alkylsilanes. Second, in contrast to previously studied densely packed monolayers (both neutral and zwitterionic), our sparse monolayers are designed to mimic the aqueous environment found in experimentally studied pMPC brushes where water can interpenetrate the polymer film. Clearly the experimental brush systems are more complex due to the network of chain–chain interactions but the goal of this study is to clearly capture the MPC–water interactions, largely eliminating the chain–chain interactions to establish a baseline behavior.
Two systems were prepared with average chain surface densities of 1.02 and 2.03 chains per nm2 which corresponds to monomers attached to 20 and 40% of the under-coordinated surface oxygens, respectively. The chain arrangement for both systems is shown in the ESI (Fig. S1†) where opposing substrates were mirrored. The sparser chain density (1.02 chains per nm2) was chosen because monomers in the experimentally studied polymer brush systems are likely to possess higher mobility than would be granted by a densely packed monolayer on a flat substrate. The denser system (2.03 chains per nm2) compares to chain densities found in previous simulation studies of zwitterionic monolayers (1.98 chains per nm2)25 and thus provides a reference point. The previously studied monolayers, however, consisted of longer, more linear monomers than MPC and were formed in regular lattice arrangements that did not allow for water to penetrate into the monolayer. Additionally, observing similar results across two chain densities allows us to rule out anomalous behavior arising from one particular system. Both systems were then solvated in TIP3P water to create a density of ∼1.0 g cm−3 in the interstitial water layer post-compression, at the largest separation distances where friction simulations were performed. The separation distance, D, between opposing monolayers is defined as the distance between the outermost layers of silicon atoms in each substrate. A schematic overview of the system setup is given in Fig. 1B. The number of water molecules remains constant for all separation distances thus representing the edge case where compression happens fast enough to prevent significant squeeze out of water molecules during the elapsed simulation time. This model makes sense for a closed system where squeeze out occurs on a much longer timescale than the simulations performed here; e.g., we are considering a nanoscale segment of what would typically be a milli- to centimeter scale substrate where water may not be able to rapidly escape from the highly absorbent zwitterionic films. The other edge case, not considered in this work, would consider complete equilibration with bulk water outside of the pore at each separation distance.
Simulations were performed in LAMMPS29 using OPLS-aa parameters for both the MPC monomers30 and the substrate.31 To ensure charge neutrality of the system, partial charges were adjusted from their values for free chains, which by construction are essentially charge neutral; these adjusted charges were localized to the region around the attachment site to avoid significantly altering more relevant components of the original model. The partial charges used in this work are summarized in Table S1.† The rRESPA multi-timestep integrator32 was applied to the MPC chains and water molecules with timesteps of 0.25 fs for bonds, 0.5 fs for angles/dihedrals and 1 fs for non-bonded interactions. Periodic boundary conditions were applied in the x- and y-directions to simulate an infinite surface. Electrostatic interactions were modeled using the 2D particle–particle particle–mesh algorithm33 where the real space component is calculated within a 1 nm radius. Lennard–Jones interactions were calculated with a 1 nm cutoff radius. The temperature of the water and chain molecules was maintained at 300 K using the Nosé–Hoover thermostat.34,35
Two distinct sets of simulations were performed. First, compression runs were carried out by fixing the bottom substrate in space and displacing the top substrate down at 0.001 nm ps−1. Second, using configurations obtained from the compression runs as input, systems with fixed separation were examined under shear. Shearing simulations were performed at separation distances of 4.1–2.8 nm and 3.8–2.5 nm in increments of 0.05 nm for chain densities of 2.03 and 1.02 chains per nm2, respectively, for a total of 54 shearing simulations. Shear stresses were induced by displacing the top substrate at 10 m s−1 in the x-direction. From each independent shearing simulation, 15–30 ns of post-equilibration data were collected; systems with higher separation distances were sampled for longer periods of time due to larger fluctuations in shear forces. Normal and shear forces were calculated by summing the forces experienced by the bottom substrate plus chains in the z- and x-directions, respectively. These forces were sampled every 1 ps and block averaged across each simulation to yield the data points and standard errors shown in Fig. 4.
Although shearing speeds of 10 m s−1 are greater than those typically studied experimentally, they are frequently encountered in practical applications. For example, micro- and nanoelectromechanical devices routinely operate at speeds on the order of 10–20 m s−1 and natural joints can experience shearing on the order of 1 m s−1, e.g., when a kicking motion is made where angular velocities can reach over 20 rad s−1.36 Considering artificial joints or cartilage scaffolds as real world applications for materials lubricated via the mechanisms studied here, speeds experienced in knee joints are of obvious interest. Additionally, both experimental5 and simulation studies28 of shearing, functionalized surfaces show little dependence on shearing velocities across several orders of magnitude. This is not to suggest that there is no dependence but that the effects are minimal within the, respectively, accessible ranges of shearing speeds for experimental and simulation studies.
![]() | ||
Fig. 2 Normal force responses during compression of MPC systems. In main plot, shown forces are acting on substrate only in order to highlight the effects of chain density. Summing the forces acting on both the substrate and chains shifts the curves down the y-axis, yielding a base line normal force (at high separation) of approximately zero; additionally, shifting the separation distance by 0.22 nm for the denser system shows almost complete overlap in trend behavior for the two systems (shown in inset). The diagonal solid line approximately denotes the two regimes (a and b) discussed in text and highlighted in Fig. 3. The upper and lower horizontal dashed lines correspond to ranges of separations across which we performed shearing simulations for the respective chain densities. |
Two clear regimes are observed with respect to compression (hereinafter referred to as a and b). In systems of 2.03 and 1.02 chains per nm2, the transitions between these two regimes respectively occur at approximately 3.4 nm & 3.1 nm. This transition is consistent between the two surface densities, although, as noted, shifted to lower separations for the less dense monolayer. In regime a, as highlighted in Fig. 3, an interstitial water layer separates the opposing monolayers. The transition, from a to b, approximately corresponds to the first contact of terminal groups of chains attached to opposing substrates, rather than just a pure interstitial water layer (compare Fig. 3a and b). Throughout the transition and beyond, a large spike in normal force occurs. Examination of the total amount of water within the monolayer region as a function of separation reveals a similar trend to Fig. 2; relatively constant values within regime a and an increase as we cross over into regime b (see Fig. S2†). These phenomena suggests a high resistance to compression marked by swelling of the monolayer with water, which is qualitatively consistent with experimental studies of polymeric brushes built from the same monomers.4–8
![]() | ||
Fig. 3 Snapshots of systems with 1.02 (left) and 2.03 (right) chains per nm2 at various separation distances. Regime (a) interstitial water layer separates surfaces. Regime (b) interstitial water layer no longer exists and chains from opposing substrates contact one another. Substrates are shown as black van der Waals spheres, chains as red van der Waals spheres and water is depicted as blue van der Waals spheres in the left image for every configuration and a surface map in the right image for the same configuration. All visualizations were performed using VMD.37 |
![]() | ||
Fig. 4 Block averaged shear stresses as a function of normal force per area. Colors correspond to Fig. 2. Error bars represent standard errors of the means; x-axis error bars are plotted but smaller than data markers for all points. Friction coefficients, μ, correspond to the slope of Fsversus Fn as determined by linear least squares fitting. R-values are shown in parentheses. The friction coefficients and fitted slopes near a and b respectively correspond to the regimes discussed in text and highlighted in Fig. 3. |
As normal loads trend towards zero, shown in the top left inset of Fig. 4, the average shear stresses fluctuate about zero. Portions of this regime correspond to the normal force fluctuations observed during compression beginning around 3.8 nm for the sparser system and 4.1 nm for the denser system (Fig. 2). We suspect that some of this behavior is caused by the transition from the attractive regime, where the two substrates are not yet in compressive contact, to the repulsive regime.
At normal forces exceeding ∼35 nN (1.8 GPa), which corresponds to regime b, shear stresses monotonically increase in a linear trend. In regime b, we find friction coefficients of 0.072 for the sparser system and 0.077 for the denser system. This magnitude of friction coefficients is comparable to those reported in previous experimental and simulation studies of neutral, densely packed monolayers without water38–42 but again at higher normal loads. Due to consistent chain–chain interactions between opposing monolayers, the shear response in this regime is likely more similar to that of systems without water which gives rise to the comparable friction coefficients. Just as in other studies without water, chains from opposing monolayers provide the primary normal load and shear stress response as opposed to the water molecules exchanging via the hydration lubrication mechanism. The chains are far less mobile than water molecules and are thus unable to provide a comparably fluid shear response yet they are able to provide a strong resistance to compression as shown in Fig. 2.
We also note that additional shearing simulations conducted at 1 m s−1 demonstrate the same trends; these results are summarized in Fig. S3.†
![]() | ||
Fig. 5 Velocity profiles of system components in the shear direction. The dashed lines denote the monolayer film thickness, wherein 90% of the monolayer atoms remain over the course of a simulation. Film thickness calculation is detailed in Fig. S4† and trends as a function of separation are shown in Fig. 6. Note that in the bottom two panels, interpenetration of chains from opposing monolayers occurs for portions outside the defined film thickness cutoff thus displacing the interstitial water layer. |
Fig. 7 further details the distribution of water across the pores. The bulk of the water remains in the interstitial area but spreads into the outer parts of the monolayers (approximately denoted by the dashed lines) where it appears to bind near where the choline headgroups are located. As the systems are compressed and they transition from regime a to b the interstitial water does densify but water is also pushed into the monolayers as is particularly evident in the sparser system.
![]() | ||
Fig. 7 Water density across xy-planes as a function of normal force. The y-axis denotes the normalized position in the pore from the bottom substrate at 0 to the top substrate at 1. Solid vertical lines denote approximate transitions between regimes a and b as discussed in text and highlighted in Fig. 3. The dashed lines denote the monolayer film thickness, wherein 90% of the monolayer atoms remain over the course of a simulation and which the choline headgroups frequently protrude from. Note that these are local water densities which do not include the monomers. |
![]() | ||
Fig. 8 Water flux across xy-planes as a function of normal force. The y-axis denotes the normalized position in the pore from the bottom substrate at 0 to the top substrate at 1. Solid vertical lines denote approximate transitions between regimes a and b as discussed in text and highlighted in Fig. 3. The dashed lines denote the monolayer film thickness, wherein 90% of the monolayer atoms remain over the course of a simulation and which the choline headgroups frequently protrude from. |
Footnote |
† Electronic supplementary information (ESI) available: Details of system charges and chain arrangements as well as film thickness results. See DOI: 10.1039/c4sm02883j |
This journal is © The Royal Society of Chemistry 2015 |