Investigation of role of cartilage surface polymer brush border in lubrication of biological joints

Although experimental evidence has suggested that the polymer brush border (PBB) on the cartilage surface is important in regulating fluid permeability in the contact gap, the current theoretical understanding of joint lubrication is still limited. To address this research gap, a multiscale cartilage contact model that includes PBB, in particular its effect on the fluid permeability of the contact gap, is developed in this study. Microscale modeling is employed to estimate the permeability of the contact gap. This permeability is classified into two categories: For a gap size > 1 µm, the flow resistance is assumed to be dominated by the cartilage roughness; for gap size < 1 µm, flow resistance is assumed to be dominated by the surface polymers extending beyond the collagen network of the articular cartilage. For gap sizes of less than 1 µm, the gap permeability decreases exponentially with increasing aggrecan concentration, whereas the aggrecan concentration varies inversely with the gap size. Subsequently, the gap permeability is employed in a macroscale cartilage contact model, in which both the contact gap space and articular cartilage are modeled as two interacting poroelastic systems. The fluid exchange between these two media is achieved by imposing pressure and normal flux continuity boundary conditions. The model results suggest that PBB can substantially enhance cartilage lubrication by increasing the gap fluid load support (e.g., by 26 times after a 20-min indentation compared with the test model without a PBB). Additionally, the fluid flow resistance of PBB sustains the cartilage interstitial fluid pressure for a relatively long period, and hence reduces the vertical deformation of the tissue. Furthermore, it can be inferred that a reduction in the PBB thickness impairs cartilage lubrication ability.


Introduction
Articular cartilage is a biological tissue located in the diarthrodial joints of vertebrate animals. It encompasses the ends of long bones in a synovial fluid-filled lined cavity. Although articular cartilage is only 2−4 mm thick, it can sustain extreme biomechanical conditions. For example, in knee joints, cartilage must withstand a vertical load up to three times the body weight during walking [1] while having a remarkably low initial coefficient of friction, on the order of 0.01 [2]. For comparison, even the best-manufactured bearing (e.g., Teflon) can only achieve a coefficient of friction 0.05−0.08 under a 3.4 MPa static load [2].
In addition to earlier lubrication theories (e.g., weeping [3] and boosted lubrication [4]), the concept www.Springer.com/journal/40544 | Friction of "hydration lubrication" [5,6] is delved to the effects on cartilage tribological performance of the "surface amorphous layer" [7], which includes polymer brushes tethered to the cartilage surface. As shown in Fig. 1, cartilage comprises chondrocytes and an extracellular matrix (proteoglycan and type II collagen), exhibiting a zonal structure throughout its thickness [8,9]. Cartilage surface roughness is formed by bundles of collagen fibrils within the superficial zone [9]. The reported roughness heights depend on the length scale at which they are measured. For example, in the normal human knee cartilage, for small length scales of 1−2 μm, asperities are small and measured in tens of nanometers, whereas, for the typical contact measurement length scale exceeding 500 μm, the reported roughness heights are relatively consistent, i.e., 5−10 μm [9]. Most importantly, transmission electron microscopy images revealed that an acellular, non-collagenous amorphous layer appeared on top of collagen fibrils [10]. The tethered layer within the amorphous layer was formed by polymers embedded in the cartilage surface extending beyond the collagen fibril defined surface. The polymers included molecules such as hyaluronan (HA), aggrecan (GAG), lubricin, phospholipids, and various other proteins [11], which formed a "polymer brush border (PBB)" on the cartilage surface. The thickness of this PBB was in the range of 200 nm (approximately the height of a single lubricin molecule [12]) to a few microns, varying with species, the joint type, or the age [13]. It is known that the negative fixed charge density of GAG molecules provides resistance to fluid flow in articular cartilage [14]; therefore, the authors postulated that PBB tethered to the surface of articular cartilage might reduce the permeability in the contact gap to lateral fluid flow. This reduced permeability might affect cartilage lubrication. It is hypothesized that the negatively charged polymers on the cartilage surface may support large contact stresses without being salvaged out of the gap (unlike the remainder of the amorphous layer) because these polymers are tethered to the cartilage surface [15,16].
Evidence supporting the idea that PBB may be important in joint lubrication has been reported in several recent experimental studies. For example, it was observed that the selective digestion of HA and GAG increased the friction force on cartilage samples by 2 and 10 times, respectively [17]. In addition, the initial friction coefficients for cartilage samples with PBB removed were higher than those of intact samples [18]. However, our current theoretical understanding regarding the roles of PBB in cartilage lubrication is still limited.
Most cartilage lubrication models are typically formulated based on the assumption of a perfectly smooth surface [19]. This assumption disregards both the effect of surface roughness and the "contact gap space", which is created between opposing surfaces as surface asperities begin to form contact with one another (i.e., when surface asperities are in contact initially, the gap size in the normal human knee cartilage is h = 10−20 μm, i.e., approximately twice the roughness height for cartilage-on-cartilage contact [9]).  A recent coupled contact model developed by the authors was used to investigate cartilage lubrication in the mixed-mode regime (i.e., considering the surface roughness and contact gap space). This model revealed that interstitial fluid exuded from cartilage tissue into the contact gap by asperity contact significantly extended the mixed-mode duration [20]. However, the viscosity of the interstitial fluid exuded from the underlying cartilage into the gap was relatively low (approximately that of water, i.e., ~0.001 Pa·s). Therefore, the viscosity of the support fluid would presumably be decrease as synovial fluid is diluted by the exudate from the articular cartilage. When the viscosities of synovial fluid are in the range of 0.01−0.1 Pa·s (corresponding to the shear rates of physiological activities of 10 2 −10 4 s −1 [21]), the simulation results of the coupled contact model [20] suggested that the gap fluid pressure could be sustained for only a relatively short time compared with the experimentally measured times for cartilage consolidation. These modeling results suggest that our initial contact model is incomplete; hence, we performed further investigations.
The focus of this study was to investigate the possible role of PBB in cartilage lubrication. Specifically, we hypothesized that for narrow contact gap sizes, PBB could potentially provide sufficient resistance to the exudate fluid flow to maintain the fluid pressure in the contact gap for a duration that represented a significant fraction of the consolidation time of articular cartilage (on the order of 1 h).
Although a previous study [22] attempted to model surface polymers as a second softer biphasic tissue on top of cartilage, the permeability of the soft layer was not assessed (it was simply assumed to be the same as that of the cartilage tissue). In this study, more sophisticated models were developed to accurately evaluate the permeability at the contact interface. The modeling involves dividing the contact gap into two layers and establishing two sets of nonlinear relationships between gap permeability and gap size in their respective layers as the contact persists. By assuming that the thickness of PBB is approximately 1 μm [23,24], we hypothesized that for gap sizes > 1 μm, fluid flow was primarily resisted by the surface roughness obstruction effect with a viscous synovial fluid and that the contact gap permeability could be estimated using a microscale computational fluid dynamics (CFD) model. Meanwhile, for gap sizes < 1 μm, we assumed that the contact gap permeability was dominated by the charged polymers tethered to the cartilage surface, the gap permeability was decreased exponentially with GAG concentration, and the GAG concentration varied inversely with the gap size. With the contact gap permeability estimated for all gap sizes, this permeability was then employed in a macroscale cartilage contact model, in which both the contact gap and cartilage tissue were modeled as two interacting poroelastic systems. The fluid exchange between these two systems was achieved by imposing pressure and flux continuity boundary continuity conditions. The numerical procedure adopted in this study is depicted schematically in Fig. 2. Next, we described the methods in more detail, specifically the development of our multiscale mathematical model and its numerical solution. Subsequently, we compared the model predictions of important gaps and tissue parameters in the presence and absence of PBB. In a series of parametric studies, we evaluated the effects of the initial gap size, viscosity of synovial fluid, and thickness of PBB on cartilage lubrication. Based on these studies, we concluded that PBB is crucial to the lubrication of normal synovial joints.

Overview of this study
A theoretical model was developed in this study to investigate the role of PBB in cartilage lubrication. As a case study to test the hypothesis that the loading of PBB significantly affects the gap fluid pressure and its duration, and in vitro indentation on a large cartilage disc was simulated computationally. As shown in Fig. 3(a), an unconfined compression experiment was simulated in the model. The model geometry was simplified to be axisymmetric, representing an explant obtained from one of the tibial and femoral condyles. During the numerical experiment, the cartilage disc was immersed in synovial fluid and vertically compressed for 1 h by a rigid, impermeable, perfectly smooth indenter. On the indenter, a uniformly distributed quasistatic load at  t = 1 MPa (i.e., 314 N) was applied. The simulation began when the indenter established contact with the highest asperity of the cartilage surface (i.e., at the onset of the mixed-mode lubrication). As depicted by the microscopic view shown in Fig. 3(b), once contact was initiated, an interconnected pore space (termed the "contact gap") was formed. The initial gap size  | https://mc03.manuscriptcentral.com/friction h 0 was equal to the peak asperity height R p of the surface roughness (5−10 μm [9]). Under a persistent and constant load, the gap reduced gradually as the consolidation of the cartilage tissue progressed. This gap closure is described by the gap height h, The fluid flow in the contact gap was governed by Darcy's law, and the permeability was dependent on the gap size. When the gap height h was greater than 1 μm, the asperities provided resistance to radial flow in the gap, with flow resistance primarily originating from the roughness obstruction owing to the viscous drag of synovial fluid flowing around the asperities. However, when h was less than 1 μm, the surface-tethered polymer brushes occupied most of the contact gap, and PBB was assumed to contribute primarily to the radial flow resistance.
The permeability of PBB was dependent on the GAG concentration in the gap [25]. As the GAGs protruded into gaps or were bound to HA protruding from the cartilage surface [11], the actual spatial variation of the GAG concentration was expected to vary with distance from the cartilage surface into the gap space. To estimate the gap permeability due to PBB, we identified the constraints that bounded its magnitude. Hence, we first assumed that the switch in the primary source of permeability that occurred at h = 1 μm was continuous. Next, we assumed that as h→0, the permeability in the contact gap approached the permeability of the underlying cartilage tissue. Finally, between these two bounding permeabilities, we assumed that the logarithm of the gap permeability was decreased linearly with the gap size. Based on these assumptions, we can define the permeability in PBB at all times when h is less than 1 μm.
Some additional key assumptions employed in the model were as follows: For simplicity, in this analysis, it is assumed that the viscosity of synovial fluid remains constant at 0.01 Pa·s during indentation. However, the viscosity of synovial fluid is, in fact, shear rate dependent (0.01 Pa·s corresponds to a shear rate > 1,000 s −1 [26]); As described in our previous study [20], the model assumes an exponential constitutive equation exists that can describe the relationship between gap closure and contact stress.

Cartilage tissue model
The cartilage tissue model adopted in this study was established within the poroelastic framework [27,28]. Three primary features of the extracellular matrix (i.e., GAG-dependent permeability, GAGdependent compressive modulus, and tensioncompression nonlinearity) were incorporated in the model. Zhang et al. [29] validated this cartilage tissue model against experimental measurements [30]. This cartilage model is summarized below.
It is assumed that the cartilage tissue spatially overlaps a combination of a solid matrix and a fluid phase [31]. Based on poroelasticity theory, the fluid flow in the cartilage tissue is governed by Darcy's law as follows: where and  s = 0.2 is the solid volume fraction [33].
In this study, it was assumed that the cartilage tissue experienced a small deformation; therefore, a constant fluid and a solid volume fraction were adopted in the formulation. The conservation of momentum equation for the cartilage tissue is expressed as where t σ is the total applied stress tensor, which is the sum of the solid matrix stress s σ and fluid stress f σ in the tissue.
where σ s E is the incremental effective stress due to the deformation of the solid phase, and I is the identity tensor.
In this study, it could be reasonably assumed that the cartilage tissue experiences negligible rotations www.Springer.com/journal/40544 | Friction as the cartilage was deformed by a rigid indenter. In this case, an infinitesimal strain formulation could be used to simulate the mechanical behavior of the cartilage. A nonlinear elastic material model was developed to simulate the stress-stiffening behavior of cartilage under tension and compression.
Cartilage is anisotropic and inhomogeneous, exhibiting tension-compression nonlinearity. The compressive stiffness in the model was governed by a nonlinear deformation-dependent GAG concentration, whereas, the tensile stiffness was regulated by the collagen volume fraction and direction. The details of the constitutive model are described below.
Experimental results suggested that the cartilage tissue permeability K c was dependent on the GAG concentration [25]. In this study, K c was assumed to be isotropic and can be obtained by calibration with experimental observations as follows [25]: where n = 5.4 × 10 −22 m 2 and m = −2.37 are empirical parameters obtained from a previous study [33], c  is the viscosity of the interstitial fluid (0.0007 Pa·s) and c agg c is the actual GAG concentration in the cartilage tissue.
It is noteworthy that the actual GAG concentration (i.e., milligrams of GAG/mL of extrafibrillar volume) was higher than the apparent GAG concentration (i.e., milligrams of GAG/mL of cartilage tissue). Miramini et al. [33] explained this difference. The relationship is expressed as follows: where c agg, 0 c is the initial apparent GAG concentration, ξ is the volume fraction of the collagen network, which is approximately equivalent to the solid volume fraction of the tissue (i.e., 45% in superficial zone, 30% in middle zone, 25% in deep zone [33]), s ( ) J t is the volumetric change of the solid phase, which is equal to the Jacobian determinant of the deformation gradient of the solid phase F S (i.e., J s (t) = det(F s ). Furthermore, the GAG concentration was inhomogeneously distributed throughout the cartilage depth. The measured apparent GAG concentration was the lowest in the superficial zone (approximately 25 mg/mL at the tissue surface) and increased linearly with the depth to approximately 120 mg/mL in the deep zone [34]. Therefore, in this study, c agg c = 25 and 120 mg/mL were adopted at z = 0 and −3 mm, respectively, and the values in between were obtained by linear interpolation. The permeability computed at the tissue surface K c (z = 0 mm) was 5 × 10 −15 m 2 /(Pa·s), which is within the range of typical values measured from healthy cartilage samples (0.5 × 10 −15 −8 × 10 −15 m 2 /(Pa·s) [35]).
The aggregate modulus (H −A ) was dependent on the GAG concentration. It has been experimentally demonstrated that the aggregate modulus at the equilibrium state was increased with GAG content and decreased with increasing water content [36]. A quadratic equation has been proposed to capture the relationship between the actual GAG concentration ( c agg c ) and the aggregate modulus [29].

 
Furthermore, the elastic compressive modulus of cartilage tissue E c can be computed using Eq. (8) [9]: where 1  = 0.25 MPa and 2  = 0.0155 MPa are empirical constants obtained from a previous study [33], and  is the Poisson's ratio of the (drained) GAG matrix, typically set as zero [35]. The tensile and shear resistance of the cartilage were provided by the collagen network, and the moduli were dependent on the collagen volume fraction, which varied with the cartilage zones. The cartilage tissue was partitioned into three zones along with the depth: the "superficial zone" (SZ, 5% of cartilage thickness), the "middle zone" (MZ, 45% of the cartilage thickness), and the remaining "deep zone" (DZ). For simplicity, constant tensile and shear moduli were adopted in this study. Following Miramini et al. [33], the shear moduli are 3, 3, and 2 MPa for SZ, MZ, and DZ respectively. The tensile moduli are different in two directions: (1) in z-axis, they are 25, 10, and 15 MPa for SZ, MZ, and DZ, respectively; (2) in r-axis, they are 100, 30, and 10 MPa for SZ, MZ, and DZ, respectively [33].

Contact gap model
A contact gap model was proposed by Liao et al. [20]. In the current study, we extended this model by incorporating PBB into the contact gap model. The contact gap model was formulated as a poroelastic system comprising three sets of equations: (1) Darcy's law governing the fluid flow in the gap; (2) the mass balance equation and fluid pressure continuity regulating the fluid exchange between the contact gap and cartilage tissue; (3) the momentum balance equation stating the stress equilibrium state in the contact gap, in which an exponential constitutive equation is assumed for the asperity local deformation under effective contact stresses.

Contact gap flow and gap permeability
Owing to the large length scale difference between the contact gap (5−10 μm) and cartilage thickness (2−4 mm), the fluid flow in the contact gap was approximated as a one-dimensional problem and modeled based on Darcy's law as follows: where g d v is the Darcy velocity of the gap flow, g p is the gap fluid pressure, and g K is the gap permeability. In Darcy's law, permeability can be regarded as the ratio of the intrinsic permeability of the pore network (in this case, the gap space) to the fluid viscosity. The intrinsic permeability is governed by the pore size, shape, and connectivity. The gap permeability for h > 1 μm can be numerically computed by upscaling a microscale gap flow model through a homogenization process [37,38]. The methodology and simulation process are outlined in Fig. 4.
To accurately simulate the gap flow, two cartilage surface roughness profiles measured from the bovine lateral femur (LF) and medial tibia (MT) using a Dektak stylus profilometer from previous studies [37,39] were adopted in this study. The samples (n = 3 for each condyle) were stored at -20 ℃; prior to testing, they were thawed in phosphate-buffered saline (PBS) to room temperature. During the measurement, the moisture of the samples was maintained throughout the imaging process. The samples were scanned at a speed of 100 μm/s with a resolution of 0.33 μm/pt on the y-axis and 200 scans were conducted on the x-axis, and the resolution of the z-axis was 8 nm [37,39]. Both roughness profiles were for an area measuring 1,000 μm × 1,000 μm. These profiles were imported to a CFD model as representative elementary volumes for the microscale gap flow [37]. An isothermal, laminar, and incompressible fluid flow with constant viscosity in the contact gap was assumed for the microscale CFD model, which was modeled based on the Navier-Stokes equations without the body force.
where  g = 1,225 kg/m 3 , which is the density of synovial fluid; u is the fluid velocity vector;  g = 0.01 Pa·s, which is the viscosity of synovial fluid in the gap. of the CFD mode [38], and the gap permeability is expressed as As the permeability was dependent on the gap size, separated CFD models were developed for gap sizes larger than 1 μm, and the results were approximated using a trendline.
Quantitative information regarding the GAG concentration within PBB is limited. However, the GAG concentration on the cartilage tissue surface has been determined previously via magnetic resonance imaging (MRI) to be approximately 25 mg/mL [34]. The gap permeability at h = 1 μm computed by the CFD method above was on the order of 10 -11 -10 -10 m 2 /(Pa·s). If we assume that the gap permeability is continuous over h = 1 μm, then by extrapolating the experimental curve of the hydraulic conductivity of the proteoglycan solution [25], the magnitude of permeability estimated by the CFD method is approximately equivalent to an average GAG concentration of 0.65 mg/mL in the contact gap. This concentration is extremely small, i.e., less than the minimum value of the typical concentration range tested in the experiment [25], and its effect on the gap permeability is minimal. Therefore, it is reasonable to disregard the effect of the GAG concentration on the gap permeability at h = 1 μm. Furthermore, at h = 0, it is reasonable to assume that the GAG concentration in the contact gap approaches the GAG concentration on the cartilage tissue surface ( c agg c , set as 25 mg/mL in this study [34]). For gap size between 0 and 1 μm, similar to the GAG profile in the cartilage tissue [34], it is assumed that the "effective GAG concentration" in the contact gap varies inversely with the gap size; therefore, its depth-dependent concentration profile during deformation can be simplified to where PBB t = 1 μm is the thickness of PBB.
According to previous experimental findings [25], the gap permeability in PBB (h < 1 μm) was assumed to decrease exponentially with the effective GAG concentration, approximated as follows: where a = [ln(Kc|z = 0 mm/Kg|h = 1 μm)]/ PBB t and is an empirical constant. By substituting Eq. (12) into Eq. (13), the variation of the gap permeability in PBB during deformation can be defined, as shown in Fig. 5. It exhibits a linear variation on a semi-log plot between the two bounding permeabilities.

Fluid exchange
Two flow paths exist for the fluid in the contact gap. One path is along with the lateral gap space and is modeled by Darcy's law, as shown in Eq. (9). The where  g v is the volumetric strain of the contact gap and s is the fluid exchange between the gap and tissue per unit volume per unit time.
An integration of s over the gap space reveals the fluid flow rate into or out of the contact gap space; it is associated with the vertical component (z-axis) of the Darcy velocity of cartilage tissues at the contact surface, which is detailed in Section 2.2.1. It is noteworthy that s may be into or out of the tissue; however, our previous study showed that s is a "source" term (i.e., s > 0), meaning that the interstitial fluid in the cartilage tissue "weeps" (or exudes) into the gap space from the cartilage [20]. However, the fluid exudate from the cartilage had much lower viscosity compared with the synovial fluid, which exhibited an effect that increased the gap permeability, thereby reducing the duration of elevated gap fluid pressure and accelerating gap closure. Hence, in this study, we investigated a model that considered PBB.

Constitutive relationship for gap space
A constitutive equation is required to describe asperity compression during the closure of the gap space. First, the effective stress principle must be defined based on porous media theory. The principle of effective stress states that the total stress is supported by the solid phase stress (  c ) and fluid phase stress ( g p ) within the contact gap, as follows: where  c is the effective asperity contact stress. It is noteworthy that Eq. (15) represents the vertical stress equilibrium state along the z-axis. The volumetric strain of the gap is primarily related to the reduction in gap size h under the asperity contact stress c  . Because the stress-strain relationship of the cartilage tissue is exponential, as per experimental observations [40], it is reasonable to assume an exponential constitutive equation for the asperity deformation (i.e., gap closure) under contact stress c  , as follows: where 0 h is the initial (i.e., the first asperity contact) gap size, which is equal to the peak roughness height in our case, i.e., 0 h = R p . As shown in Fig. 4, the values of 0 h are approximately 5 and 9 μm for the LF and MT surface, respectively. β is the stiffness of the cartilage asperity, which was set as 1/5 of the cartilage tissue aggregate modulus H -A , as reported by Graindorge et al. [22]. The cartilage tissue aggregate modulus H -A is detailed in Section 2.2.1.

Boundary and initial conditions
To couple the governing equations of the cartilage tissue and contact gap, a few boundary conditions [41] must be employed at the contact interface (z = 0 mm), as follows: It is noteworthy that the unit normal to the contact interface is denoted by n. Equations (17a) and (17b) ensure the continuity of the fluid flux (i.e., the Darcy velocity) normal to the cartilage-gap boundary and the fluid pressure across the cartilage-gap boundary, respectively. In addition, the total surface traction  t normal to the contact gap and cartilage tissue is continuous.
For both the cartilage tissue and contact gap, the fluid pressure at the perimeter edge is equal to the reference ambient fluid pressure, typically set as zero [42].
The osteochondral junction was assumed to be a fixed and impermeable surface; hence, these boundary conditions were applied in the model [33].
Furthermore, an initial condition was required for the contact gap model. This is expressed by the stress equilibrium state shown in Eq. (15). The analysis starts when a contact is established with the highest asperity. At this instant, the contact stress www.Springer.com/journal/40544 | Friction c  = 0; therefore, the total applied load is assumed to be solely resisted by the fluid pressure in the gap space (i.e., g   p t = ), and the gap size is at its maximum extent (i.e., h = 0 h ). Mathematically, When the gap begins to close, the surface asperities are deformed and the gap fluid is squeezed out; consequently,  c increases while g p decreases.

Computational modelling
Computational models in both the microscale (CFD model) and macroscale (cartilage contact model) were conducted using the commercial software package COMSOL Multiphysics (version 5.3, COMSOL, Inc.).
A typical microscale CFD model and its boundary conditions are shown in Fig. 4. For an approximation, fluid exchange and fluid-structure interaction were not considered in the microscale model. The input pressure was 100 kPa, which resulted in the same initial fluid pressure gradient as that in the macroscale model. The gap closure was modeled by slicing up the asperities by the upper wall at different gap sizes, as the Poisson's ratio of the cartilage extracellular matrix was approximately zero [35]. A mesh sensitivity test was performed in advance with "coarse", "normal" and "fine" mesh options in COMSOL [43]. The results showed that the differences in u vol.ave were less than 1%. To balance simulation time and accuracy, the "normal" mesh was adopted; 296,147 free tetrahedral elements were used, including two boundary layers (131,873 elements) at the upper and lower walls and 1,738 corner elements. The parallel sparse direct solver was selected, and the relative tolerance was set to 0.001.
The geometry and dimensions of the macroscale integrated contact model are shown in Fig. 3. The dimensions of the model were obtained from MRI images [44]. The cartilage thickness remained relatively constant around the center (3 mm thick) and was then gradually decreased toward the edge of the tibial plateau (1 mm thick), from a distance of 2/3 the cartilage disc radius (10 mm).The model was meshed using 1896 free tetrahedral elements, in which the average element size was 0.1 mm. The numerical analysis was halted after 1 h of simulation for a 1-MPa indentation. The model was solved by the time-dependent implicit solver using the backward differentiation formula time stepping method. The relative tolerance was set to 10 -3 .

Results and discussion
In this section, the effects of PBB on cartilage lubrication are first analyzed. Subsequently, the effects of the roughness vertical height, synovial fluid viscosity, and PBB thickness on cartilage lubrication are discussed.

Effects of PBB
To investigate the effects of PBB on cartilage lubrication, we considered models that either included or disregarded the presence of a PBB. Both models were used to simulate an LF surface as an example surface. The gap permeability curves are shown in Fig. 5. For gap sizes of less than 1 μm, it was observed that the gap permeability for the model without PBB was decreased gradually as the gap began closing. This occurred because only the gap flow resistance arising from roughness obstructions and the viscous synovial fluid were considered. Meanwhile, the gap permeability with PBB decreased more rapidly with the closing gap size, compared with the case without PBB. It is noteworthy that for the part where h > 1 μm, the permeability curves were independent of the presence or absence of PBB.
A few metrics were used to understand the interaction between the cartilage tissue and the contact gap, as well as their synergistic effects on cartilage lubrication. Specifically, the metrics were the gap fluid load support fraction, cartilage interstitial fluid pressure, and cartilage vertical strain along the z-axis.
The fluid pressure distribution in the contact gap during the first 30 min of contact is shown in Fig. 6(a), in which both cases with (solid lines) and without PBB (dashed lines) are plotted together for comparison. As shown, the gap fluid pressure decayed gradually toward the contact center, and it declined more rapidly without PBB. For example, for a 10-min contact, the gap fluid pressure without PBB was decreased to 1/10 of the applied load, whereas with PBB, the gap fluid pressure was approximately five-fold greater, 120 Friction 10(1): 110-127 (2022) | https://mc03.manuscriptcentral.com/friction particularly at the area near the contact center (r < 5 mm); nevertheless, the gap fluid pressure can still be maintained at approximately 65% of the applied load. It is more meaningful to analyze the gap fluid load support f W , which is obtained by integrating the gap fluid pressure over the cartilage surface. It is a key parameter in evaluating cartilage lubrication performance, equivalent to the monitoring of hydrodynamic lubrication. The coefficient of friction is directly proportional to the normal load supported by the solid phase at asperity contacts  Fig. 6(a). In both cases, the gap fluid support fraction decreased with loading time, but the support decreased the most rapidly without PBB. For example, at 1,200 s, the gap fluid load support without PBB was almost exhausted (~2% of the total load). By contrast, the gap fluid load support fraction still remained at approximately 40% of the total load with PBB (26 higher than that without PBB). The gap fluid load support fraction decreased gradually and was less than 20% after a 1-h indention. The results above indicate that at 1,200 s when considering PBB, the asperity solid-to-solid load support was approximately 40% less than that of the model without PBB, implying a 40% smaller coefficient of friction.
To further investigate the effect of PBB, the contour plots of interstitial fluid pressure (i.e., the fluid pressure within the cartilage) are shown in Fig. 6(b). Without PBB, the interstitial fluid pressure decreased rapidly. For example, the fluid pressure near the center of the tissue (r = 0, z = -1.35 mm) decreased to 0.04 MPa after a 30-min indentation, whereas in the same condition, the fluid pressure with PBB was more than 10 times greater at 0.52 MPa. The fluid pressurization in the gap was due to the interstitial cartilage fluid exuding into the contact interface, as clearly indicated by the streamlines intersecting the cartilage surface. During the early contact stage (at t = 60 s), with and without PBB, the interstitial cartilage fluid wept into the contact www.Springer.com/journal/40544 | Friction gap space. Without PBB, this weeping process was continuous and occurred over the entire contact surface. However, owing to the fluid flow resistance provided by the surface polymers in the gap space, the weeping process decelerated significantly and occurred primarily in the region close to the contact center, as indicated by the streamlines at t = 900 and 1,800 s. Furthermore, more fluid had to be exuded from the "side outlet", thereby involving much longer drainage paths. Therefore, the rate of interstitial fluid pressure drop was decreased significantly by PBB. The slowdown of the weeping process can be further quantified by comparing the fluid exudation volume of the two cases after a 1-h indentation, as shown in Table 1. Without PBB, the fluid exudation volume to the contact interface after a 1-h indentation was 0.117 mL (88% of its total exudation volume), which was almost seven times that of its counterpart with PBB (0.017 mL).
The results suggest the critical role of PBB in prolonging the load support by maintaining the interstitial fluid pressure in the tissue, extending the weeping lubrication period, and extending the duration of the hydrodynamic mode of lubrication. Without PBB, the cartilage interstitial fluid can be rapidly squeezed out under high contact loading, thereby increasing solid-to-solid contracts at the interface as well as increasing the associated frictional wear. In other words, the fluid permeability in the contact gap with PBB will be lowered, rendering it more difficult for the fluid to exude from the articular cartilage and to be squeezed out of the contact gap space. The surface polymers fixed to the cartilage surface cannot be squeezed out of the gap space, unlike other components of the amorphous layer on the cartilage surface. Therefore, weeping exudation from the articular cartilage in the presence of PBB serves to extend the duration of the hydrodynamic lubrication mode, thereby reducing friction between contacting cartilage surfaces. To further investigate the role of PBB, we next consider its effect on the cartilage biomechanical performance. A comparison of the average cartilage vertical strains over the contact interface of the two cases is shown in Fig. 6(c). When PBB was included in the model, the average vertical strain after 30 min of loading was 11% (compared with 16% without PBB), a prediction that matched reasonably well with in vivo measurements under similar loading conditions. Halonen et al. [46] utilized computed tomography arthrography to measure the cartilage strain. The test subject was standing on one leg supporting approximately half of the bodyweight with the aid of harnesses (386 N). The total knee joint reaction force was reported to be approximately the full body weight (107%) [46]. If we assume that the load is approximately equally shared by both joint condyles, then the total force on one condyle is approximately half the body weight (386 N), which is comparable to the loading condition of our modeling in this study (314 N). The strains after 30 min of contact were obtained by comparing computed tomography (CT) images, and they were 12% and 10% for the lateral and medial tibia, respectively [46]. This comparison verifies the model predictions performed in this study for the model with PBB. This suggests that it is essential to include PBB in cartilage contact modeling for an accurate simulation.
After 1 h of indentation, the cartilage strain in the model without PBB reached an equilibrium state at 16% strain, which was approximately 19% higher than the strain with PBB (at 13%). This was due to the additional fluid exudation that occurred at the contact interface in the model without PBB. As shown in Table 1, the total fluid exudation without PBB (0.133 mL) was exactly 19% higher than that of its counterpart model with PBB (0.112 mL).
In summary, the study of the two models with and without PBB demonstrates that PBB can provide significant additional resistance to the exudate fluid flow along the contact gap, offering two benefits. First, the flow resistance in the gap space limits the rate of fluid exudation to the contact gap, thereby maintaining the interstitial fluid pressure inside the cartilage tissue and reducing tissue strain. Second, the fluid load support in the contact gap can be | https://mc03.manuscriptcentral.com/friction maintained for a much longer period that is comparable to the consolidation time of the cartilage, as the consolidation process results in exudate flowing into the contact gap space. This behavior increases the fraction of hydrodynamic lubrication at the contact interface, and hence reduces the contact friction and surface wear.

Effect of synovial fluid viscosity
The viscosity of the healthy synovial fluid can vary by several orders of magnitude owing to its shear rate dependence. Using the MT surface (Fig. 4) and its gap permeability (Fig. 5), a parametric study was performed to assess the effect of the synovial fluid viscosity on cartilage lubrication. Three constant viscosity values for synovial fluid with 10-fold differences (i.e., 1, 0.1, and 0.01 Pa·s) were compared, as well as that of water (0.001 Pa·s), whereas all the other model parameters were fixed. The results are shown in Fig. 7(a). The gap fluid support fraction decreased as the viscosity decreased. For example, at t = 1,800 s, if we consider a gap fluid support of 1 Pa·s as the reference point, every 10-fold decrease in the viscosity magnitude (to 0.1, 0.01 and 0.001 Pa·s) will result in a reduction in the gap fluid support by 23%, 40%, and 58%, respectively. In other words, for a viscosity reduction of 1,000 times (from 1 to 0.001 Pa·s), the gap fluid load support declines by 58%.
This viscosity result suggests that synovial fluid can enhance the fluid support fraction in the gap space, thereby reducing the friction coefficient. Next, we compare experimental findings with our model predictions. Forster and Fisher [47] measured the initial friction coefficient of cartilage on metal contact, in which Ringer's solution or synovial fluid was used as the lubricant. The most significant differences in friction coefficient that they recorded were at loading times of 20 and 45 min, where the coefficient of frictions using synovial fluid (μ 0 = 0.18 and 0.26) were only 21% and 16% less than those using Ringer's solution (μ 0 = 0.22 and 0.31, respectively). The experimental results matched reasonably well with our computational predictions, particularly in terms of the percentage of difference plotted in Fig. 7(b).
It can be reasonably assumed that the viscosity cases of 0.01 and 0.001 Pa·s in Fig. 7(a) correspond to the synovial fluid and Ringer's solution case of Forster and Fisher [47]. By regarding the contact gap as a porous medium, the effective coefficient of friction eff  can be computed based on biphasic lubrication theory as follows [41]: where eq  is the coefficient of friction in the equilibrium state. As shown in Fig. 7(b), by assuming eq  = 0.3 [45], the predicted coefficients of friction of η = 0.01 Pa·s at 20 and 45 min (i.e., 0.21 and 0.25) were 20% and 18% less than those of η = 0.001 Pa·s at 20 and 45 min (i.e., 0.26 and 0.30), respectively. The percentage of differences was similar to the experimental measurements [47]. This study provides a reasonable and possible theoretical explanation for the experimental observations. However, it is noteworthy that the different amounts of lubricin in the experiments may affect the friction measurements www.Springer.com/journal/40544 | Friction (particularly at 45 min when the cartilage was near equilibrium [47], i.e., mainly in the boundary lubrication regime), and the effect of lubricin on cartilage lubrication is beyond the scope of this study. Nevertheless, as suggested by Forster and Fisher [47], the friction coefficient of articular cartilage is primarily affected by fluid load support.

Effect of PBB thickness
The effect of PBB thickness on cartilage lubrication is considered in this section. Experimental studies have shown that PBB thickness can range from 200 nm to a few microns [13]. In this study, three thickness values (200 nm, 500 nm, and 1 μm) were investigated using the LF surface, and the viscosity was maintained at 0.01 Pa·s.
The results of PBB thickness are shown in Fig. 8(a). In general, the cartilage lubrication varied inversely with the PBB thickness. Considering t = 1,800 s as an example, by reducing the PBB thickness by 50% (i.e., from 1 μm to 500 nm), the gap fluid load support decreased by 33%, whereas another 60% reduction in thickness (i.e., from 500 to 200 nm) resulted in a further 84% decrease in the gap fluid load support. Figure 8(b) shows the time-dependent coefficients of friction as the PBB thickness was varied between 200 nm and 1 μm, and water (0.001 Pa·s) was used as the lubricant. The predicted time-dependent variation of the coefficient of friction for the PBB thickness of 400 nm (using Eq. (21), in which eq  is 0.33, as per the experimental measurements) can yield a reasonably close approximation to (particularly after 300 s) the experimental measurement reported by Accardi et al. [48] (the test was performed in PBS at a contact pressure of 1.2−1.8 MPa).

Limitations
This study has some limitations. First, to simplify the model complexity, the shear rate dependent viscosity of synovial fluid was not considered in the model. Our previous study showed that with reductions in gap pressure gradients with gap closure, the viscosity may increase [37]. The increase in viscosity might decrease the gap permeability and hence further prolong the gap fluid load support. Second, the interaction deformation between the asperity and the bulk tissue was not considered in the simulation. In future studies, a relationship between the local asperity deformation and tissue bulk consolidation must be established through experimental observations (e.g., measuring the surface roughness at different cartilage strains). Third, the model can be improved by the availability of experimental data regarding GAG content or the fixed charge density in PBB. This can be achieved using high-strength MRI scanners with sufficiently high resolution or the tracer cation method using 22 Na [49]. Furthermore, it may be more beneficial to directly model the interactions of the surface polymers and solvent molecules in normal joint motions [27].

Conclusions
In this study, PBB on the cartilage surface was integrated using a coupled contact model. The effect 1) PBB can substantially enhance cartilage lubrication by increasing the gap fluid load support fraction and hence improve the hydrodynamic mode of lubrication. Based on the case study and using the specified parameters, PBB increased the fluid support by 26 times at a 20-min indentation compared with the model without PBB.
2) Weeping (fluid exudation) and hydrodynamic lubrication reduced friction synergistically. The exudation of interstitial fluid from the articular cartilage into the contact gap space prolonged the hydrodynamic mode of lubrication. Owing to the resistance of PBB to the lateral transport of the exudate along the contact gap space, less fluid was required to sweep into the contact gap while the gap fluid pressure was maintained. Hence, the interstitial fluid pressure within the articular cartilage tissue can be maintained for a longer period, and cartilage deformation can be reduced compared with similar load durations.
3) Synovial fluid improved fluid support in the gap space relative to saline water, as reducing the viscosity magnitude by 1,000 times (from 1 to 0.001 Pa·s) reduced the gap fluid support at 30 min of indentation by 58%.
4) The PBB thickness significantly affected the cartilage lubrication performance; a 60% reduction in the PBB thickness (from 500 to 200 nm) resulted in an 84% decrease in the gap fluid load support at 30 min of indentation.