Designing Superlubricious Hydrogels from Spontaneous Peroxidation Gradients

Hydrogels are hydrated three-dimensional networks of hydrophilic polymers that are commonly used in the biomedical industry due to their mechanical and structural tunability, biocompatibility, and similar water content to biological tissues. The surface structure of hydrogels polymerized through free-radical polymerization can be modified by controlling environmental oxygen concentrations, leading to the formation of a polymer concentration gradient. In this work, 17.5 wt % polyacrylamide hydrogels are polymerized in low (0.01 mol % O2) and high (20 mol % O2) oxygen environments, and their mechanical and tribological properties are characterized through microindentation, nanoindentation, and tribological sliding experiments. Without significantly reducing the elastic modulus of the hydrogel (E* ≈ 200 kPa), we demonstrate an order of magnitude reduction in friction coefficient (from μ = 0.021 ± 0.006 to μ = 0.002 ± 0.001) by adjusting polymerization conditions (e.g., oxygen concentration). A quantitative analytical model based on polyacrylamide chemistry and kinetics was developed to estimate the thickness and structure of the monomer conversion gradient, termed the “surface gel layer”. We find that polymerizing hydrogels at high oxygen concentrations leads to the formation of a preswollen surface gel layer that is approximately five times thicker (t ≈ 50 μm) and four times less concentrated (≈ 6% monomer conversion) at the surface prior to swelling compared to low oxygen environments (t ≈ 10 μm, ≈ 20% monomer conversion). Our model could be readily modified to predict the preswollen concentration profile of the polyacrylamide gel surface layer for any reaction conditions—monomer and initiator concentration, oxygen concentration, reaction time, and reaction media depth—or used to select conditions that correspond to a certain desired surface gel layer profile.


Minimizing Fitting Error
The lsqcurvefit Matlab function was used to fit the approach curves for both the microindentation and nanoindentation measurements.This function uses the squared 2-norm of the residual (Eqn.S1), and minimizes this value to find the desired variable, x.

Small Angle Approximation
The small angle approximation for the trigonometric function cos is as follows: The indentation depth into the sample can be predicted via geometry (Figure S2).Eqn.S2 can be substituted for cos(θ) so: To stay within 1% error for the small angle approximation, θ = 38°.With a probe radius of R = 2.5 µm for the nanoindentations, d = 0.55 µm.

Theoretical Noise Floor for Friction Coefficient Measurements
The noise floor, or minimum detectable friction coefficient, for the tribometer was estimated with the following equation: where K f is the tangential spring constant of the double-leaf cantilever, x is the minimum detectable displacement by the capacitance probes, and F n is the applied normal force.For our experiments, K f = 100 µN/µm.F n = 1 mN, and the capacitance probes (Lion Precision, C5R-0.80-2.0)used to measure the displacement have a 5 nm resolution (x = 5 nm).

Contact Area and Pressure Estimations
The contact radius, a, was estimated based on Hertzian contact mechanics using the following equation:  =  (S5) where R is the probe radius and d is the indentation depth.The maximum contact radius was determined based on the maximum indentation depth reached at a normal force of F n = 1 mN.
The maximum contact pressure was estimated using the following equation:  The vast structural dispersity arising from polyacrylamide propagation mechanisms, including MBAm crosslinker incorporation and reaction pathways, is explored in this section and summarized in Figure S3.All individual radicals (highlighted in teal) are assumed to have equivalent reactivity in all of these reactions in the model, in order to avoid accounting for the vast dispersity of local polyacrylamide chain radical structures.These radicals would all have an equivalent-but different-reactivity with oxygen, as denoted by the difference in k p and k s in Table S5.While Figure S3 shows a number of possible reaction pathways after incorporation of a crosslinker into a growing oligomer, note that MBAm could also be the first monomer incorporated, and that the monomers can be incorporated in any sequence.
One of these reaction pathways is intramolecular cyclization, wherein the other end of an MBAm monomer already incorporated into the polymer polymerizes into the chain.While a ring formed between radicals and crosslinks within the same MBAm crosslinker was shown here, note that this process can occur for much larger cycles.Cycles are more likely to form at lower monomer concentration and are considered defects because they render the crosslinker elastically ineffective.
When the other end of an incorporated MBAm monomer does not experience cyclization, it is likely that growth of a second polyacrylamide chain will be initiated at that double-bond, establishing MBAm an elastically-effective crosslink.Note that the depiction in Figure S3 shows the example of this initiation prior to the continued growth of the existing chain, but due to the short radical lifetimes it is likely that this initiation event occurs at a later time.Aside from cyclization, the incorporation of each end of MBAm into the gel matrix are decoupled.

Figure S3
Structures formed during initiation and propagation of polyacrylamide; radicals are highlighted in teal for easy visualization.On initiation, the chain grows with AAm until incorporation of a crosslinker, denoted by example oligomer Species A. Once incorporated into a chain, the resulting radical structure can either react with AAm or MBAm monomer or react with the double bond at the other end of the MBAm crosslinker to form a cyclic defect.Additionally, the other double bond may come in contact with another radical initiator instead, resulting in polyacrylamide chains continuing to grow simultaneously and independently at both ends of the crosslinker, and thus the establishment of a crosslink.

Calculated with 17.5 wt.% gel protocol
Using the values listed in Table S5, the following constants were calculated with Eq.S7-S11: , and dimensionless time as , we obtain the following = Similarly, the dimensionless consumption of TEMED follows: Next, nondimensionalize the reaction-diffusion equation for oxygen.Here, we introduce ] and as dimensionless depth from the surface of the reacting medium: The concentration of reactive radicals [R*] is given by: or in terms of dimensionless variables, Plug this expression into the partial differential equation for and simplify to obtain: ) By defining the constants and (see main text), the final expression for oxygen diffusion is: Next, the reaction-diffusion equation for monomer follows .This expression , the dimensionless monomer concentration at a given  =

[M]
[]  depth in the reaction medium.Note that the conversion of monomer, , which is reported in the  main text and results, satisfies .
Plug in [R*] to obtain: ) Using the same definition of and defining a new Damköhler number, , the final differential equation for monomer conversion is obtained: ) Literature Comparison Table S7.Comparison of the elastic modulus, friction coefficient, and estimations of the surface gel layer thickness of polyacrylamide hydrogels cast with different monomer (acrylamide, AAm) and crosslinker (N,N'-methylenebisacrylamide, MBAm) concentrations and polymerization conditions to our data (last row).The listed friction coefficient values were determined by digitizing the data from literature at the listed testing parameters (sliding speed, v, and estimated contact pressure, P). [10][11][12][13][14][15][16][17][18][19][20] Unless otherwise stated, the tribological testing parameters are as follows: v = 0.5 mm/s and P ≈ 6-11 kPa.Gemini indicates a self-mated gel-on-gel configuration during tribological experiments (e.g., the probe material was composed of the same hydrogel material).*Testing parameters are as follows: v = 100 µm/s and P ≈ 1 kPa ** Poisson's ratio of 0.5 was used to estimate the elastic modulus for our data (17.5 wt.% polyacrylamide hydrogels).10 Friction Coefficient

Minimum Film Thickness Calculations
The Sommerfeld number, S, is a dimensionless parameter that is proportional to the viscosity of the lubricating fluid, the sliding velocity, v, and the applied normal force, F n , as . 21The      Stribeck curve plots the friction coefficient as a function of the Sommerfeld number for two sliding interfaces and partitions it into four main lubrication regimes: boundary, mixed, elastohydrodynamic lubrication (EHL), and hydrodynamic lubrication.Boundary lubrication occurs when the two sliding interfaces are in direct contact with each other whereas EHL and hydrodynamic lubrication occurs when the fluid film is thick enough to separate the two sliding interfaces.Due to this lubricating fluid layer, the friction coefficients measured within this regime are lower than those in boundary lubrication.The lubrication mode is highly dependent on applied force and sliding speed, with boundary lubrication typically occurring at lower sliding speeds and higher applied forces.Using soft-elastohydrodynamic lubrication theory developed by Hamrock and Dowson, 22 the minimum fluid film thickness, h min , for soft-EHL lubrication was estimated between 33 -36 nm using Eqn.S23.
where R is the probe radius of curvature, is the viscosity of the fluid, is the sliding velocity,  o  is the applied normal load, and where is the reduced elastic modulus of the  n  ′ = 2 *  * sample.The minimum fluid film thickness was calculated by using the viscosity of water ( = 8.9 o x 10 -4 Pa s), = 0.5 mm/s, = 1 mN, and = 200 and 240 kPa (the minimum and maximum •   n  * measured for the polyacrylamide hydrogels).If h min was significantly greater than the surface  * roughness of the hydrogel, then the measured superlubricity would be more likely due to a thick fluid film layer rather than the hydrogel surface.Since the hydrogels were cast in open air, it is difficult to estimate the surface roughness.However, using the scaling relationship µ ξ -1 ~ developed by Urueña for self-mated gel-on-gel contacts 20 and following a Monte Carlo analysis outlined by Pitenis et al., 23 the estimated mesh sizes at the surfaces of the hydrogels at 0.01 mol% and 20 mol% O 2 are ξ = 2.2 ± 0.4 nm and ξ = 31.8± 14.9 nm, respectively.Since h min is on the same order of magnitude as the estimated mesh size of the 20 mol% O 2 gel, it is likely that there is still contact between the gel surface and probe, indicating that the measured superlubricity is not solely due to a fluid film layer.
At an applied normal force of 1 mN, the maximum pressure during sliding ranged between 9 -11 kPa, depending on E* of the polyacrylamide hydrogel.Since E* scales with osmotic pressure and the applied contact pressures are much lower than E*, it is unlikely that fluid flow or draining occurred, suggesting that the friction coefficients measured are not due to fluid flow. 24e track length was chosen so the distance was at least 4 times the contact area diameter to ensure that the probe moved out of its initial contact area zone when sliding.The minimum track lengths necessary for an applied normal force of 1 mN for hydrogels with reduced elastic moduli of 200 kPa and 240 kPa are 1.8 mm and 1.7 mm, respectively.Hertzian contact mechanics was used to estimate the contact area diameter, d (Eqn.S24).Previous measurements of 17.5 wt.% polyacrylamide hydrogels cast against polystyrene in ambient conditions (20 mol% O 2 ) and tested in self-mated contact have reported friction coefficients of µ = 0.037 at an applied normal force of F n = 2 mN and sliding speed of v = 0.5 mm/s (Table S7). 20For self-mated gel-on-gel sliding configurations, friction coefficients generally decrease with increasing applied normal force since the contact radius does not scale with normal load. 25For the 17.5 wt.% hydrogels herein, the friction coefficient slightly increased from µ = 0.021 ± 0.006 to µ = 0.026 ± 0.003 at F n = 2 mN.

Extracting Surface Layer Thickness and Gradient
A series of studies by Baselga et al. in the late 1980s sought to determine the process of network formation for AAm/MBAm gelation with APS and TEMED initiator under various conditions. 29ne such study determined that a sample with initial [AAm] init = 1.28 mol/L and [MBAm] init = 0.074 mol/L reached the gel point at 6% monomer conversion based on digitization of graphical data provided in that paper. 29In a separate study, Baselga et al. found that doubling the concentration of monomer and crosslinker, while halving the weight percent of crosslinker in the initial feedstock-similar to our experimental system with [AAm] init = 2.46 mol/L and [MBAm] init = 0.05 mol/L-approximately halved the conversion (by weight) at the gel point. 30As a safe estimate, we assumed gelation occurred at 5% monomer conversion (rather than 3%) in this work.Thus, according to our model, the entirety of this surface gel layer reached gelation.
The first point where the monomer conversion exceeds 5% is taken to be the upper surface of the gel.To minimize error of fit while capturing surface shape, a depth of gel is chosen for fitting such that the surface gel layer represents about 10%.In other words, for the 0.01 mol% O 2 case, the top 100 µm of gel is fit to extract surface gel layer thickness and gradient, whereas the top 300 µm is fit for the 20 mol% O 2 case.Fitting is conducted on this cropped data at each time point using the MATLAB shape language modeling (SLM) engine.A degree 1 fit was conducted using 3 knots and free interior knots; i.e., two line segments were fit to the data.The thickness of the surface gel layer was taken to be the depth corresponding to the intersection of these two line segments.The slope of the gradient layer was taken to be the slope of the line at the surface of the gel.

Figure S1
Figure S1 Representative microindentation approach curves with the Hertz contact mechanics model overlayed in red fitting from the point of contact to F n = 1 mN for the gels cast at (a) 0.01 mol% O 2 and (b) 22 mol% O 2 .

Figure S2
Figure S2 Schematic of probe with radius of curvature, R, indenting a sample with indentation depth, d.

Figure S4
Figure S4 Comparison of friction coefficients from literature of polyacrylamide hydrogels with varying monomer and crosslinker concentrations molded against different materials (polystyrene, glass, PDMS, polyethylene, PTFE, and air-liquid interface ("free")).All samples are cast at 20 mol% O 2 unless labeled otherwise.Sliding experiments were conducted with a self-mated gel-ongel configuration unless otherwise stated.Tribological testing parameters are as follows: v = 0.5 mm/s and P ≈ 6-11 kPa.

𝑑 = 2 (Figure S5
Figure S5 Friction coefficient as a function of sliding distance over 30 cycles for the (a) 0.01 mol% O 2 and (b) 20 mol% O 2 .The friction does not significantly change over time.Note the difference in scale between the two plots.

Figure S6 11 . 2 Friction
Figure S6 Friction coefficient as a function of applied normal force, F n , for the 17.5 wt.% polyacrylamide hydrogels at 0.01 mol% O 2 and 20 mol% O 2 .(a) Friction coefficient slightly increased with increasing normal force at 0.01 mol% O 2 .(b) Friction coefficient slightly decreased with increasing normal force at 20 mol% O 2 .

Figure S10
Figure S10Representative log-log plots of the derivative (dF n /dd) of the nanoindentation approach curves as a function of force following the methods outlined by Garcia et al.26The dashed red line represents the Hertzian contact model while the solid red line represents the Fredrickson high-penetration brush model.11,27,28(a) At 0.01 mol% O 2 , the Hertz contact model fits the data better.(b) At 20 mol% O 2 , the Fredrickson high-penetration brush contact model fits the data better.

Figure S11
Figure S11 Fits of monomer conversion for the top 100 µm of gel for the hydrogels polymerized at 0.01 mol% O 2 .The gray lines show the exact line segment fit of monomer conversion, with the intersection of the line segments denoting the transition from surface gel layer to the bulk gel.Fits are shown at various reaction times.SLM algorithm fits well at all time points.

Figure S12
Figure S12 Fits of monomer conversion for the top (a) 300 µm and (b) 100 µm of gel for the hydrogels polymerized at 20 mol% O 2 .The gray lines show the exact line segment fit of monomer conversion, with the intersection of the line segments denoting the transition from surface gel layer to the bulk gel.Fits are shown at various reaction times.The SLM algorithm fits well at all time points except 3 minutes, where the root mean square error (RMSE) value is above the cutoff; clearly the surface and bulk layers are not yet distinct at 3 minutes.

Figure S13
Figure S13Root mean square error (RMSE) of the SLM fit for 0.01 mol% O 2 (dark blue) and 20 mol% O 2 (light blue).Horizontal lines indicate the RMSE cutoff value, one standard deviation above the average of the non-zero RMSE values at each O 2 concentration.Above the cutoff, the fits do not accurately reflect the shape of the surface layer, corresponding to an ill-defined distinction between the surface gel layer and bulk gel as the gel forms at early times.The red dots represent the reaction times in which the corresponding thickness and gradient of the surface gel layer were excluded from the main text.

Figure S14 ( a )
Figure S14 (a) Surface gel layer thickness and (b) surface gel layer thickness gradient as a function of reaction time with insets highlighting the decrease in surface gel layer thickness and gradient with increasing time.The red dots represent the reaction times in which the RMSE of the SLM fit was greater than the cutoff value and excluded from the main text.

Figure S15
Figure S15Unlike the 0.01 mol% O 2 case, in the 20 mol% O 2 polymerization shown here, oxygen inhibits the polymerization at the surface of the reaction so that even after 15 minutes, no significant monomer conversion or gelation occurs until about 300 µm into the reacting medium, resulting in a thinner gel.

Figure S17
Figure S17 The following are plotted at varying oxygen concentrations corresponding to the various profiles plotted in Figure S16.(a) Conversion of monomer at the surface of the gel.(b) Thickness of the surface gel layer from SLM fit of the profiles.(c) Monomer conversion gradient in the surface gel layer, also extracted from the SLM fit of the profiles.(d) Thickness of the unpolymerized liquid layer at the upper surface of the reacting medium.(e) Monomer conversion in the gel bulk, as measured by the conversion deep into the gel, which is affected by oxygen diffusion only at early times during polymerization.(f) Root mean square error (RMSE) of fit corresponding to (b) and (c).Note that the RMSE cutoff for the 0.01 mol% and 20 mol% O 2 data reported in the text is around 0.04, so all values reported here fit well.

Table S2 .
Volume and mass of each component in the pre-polymerized hydrogel solution, along with the corresponding wt.% and mol.%.

Table S4 .
Estimated contact radius, contact area, and maximum contact pressure based on Hertzian contact mechanics.

7 Reaction-Diffusion Constants and Component ConcentrationsTable S5 .
Kinetic constants from literature and experimental reaction conditions used in our model.
We note that in simulation, the dimensionless time is significantly larger than the equivalent  numerical value of real time t.When t = 1 sec, this corresponds to dimensionless time  =   O 2  2 = 2.96 x 10 -5 .Therefore, t = 15 min corresponds to = 0.0266272. 8 Reaction-Diffusion Model Derivation Consumption of APS and TEMED follow , where C represents [APS] or

Table S8 .
Compilation of the mechanics, friction coefficient, and surface gel layer thickness for the 17.5 wt.% polyacrylamide hydrogels at 0.01 mol% and 20 mol% O 2 .