1st International Conference on Mathematical and Computational Biomedical Engineering - CMBE2009
June 29 - July 1, 2009, Swansea, UK
P.Nithiarasu and R.L?ohner (Eds.)
EXACT SOLUTION TO A REFINED CONTACT PROBLEM
FOR BIPHASIC CARTILAGE LAYERS
G. Mishuris
Wales Institute of Mathematical and Computational Sciences, Institute of Mathematics and Physics,
Aberystwyth University, Ceredigion, SY23 3BZ, Wales, UK, ggm@aber.ac.uk
I. Argatov
Laboratory of Friction and Wear, Research Institute of Mechanical Engineering Problems, V.O.,
Bolshoy pr., 61, 199178 St. Petersburg, Russia, argatov@home.ru
ABSTRACT
We revisit the axisymmetric contact problem for a biphasic cartilage layer and consider a refined for-
mulation taking into account the both normal and tangential displacements at the contact interface. The
obtained analytical solution is valid for arbitrary time and increasing loading conditions. We compare
it with the classic result and indicate cases where the difference could be pronounced.
Key Words: contact problem, biphasic cartilage layer, tangential displacement, analytical solution.
1 INTRODUCTION
Solution to biomechanical contact problems for biological joints are of great importance for medical ap-
plications, especially in analysing different stages of osteoarthitis. As the cartilage layers are extremely
thin in comparison with the characteristic size of the bones while the physics of the problem is rather
complicated, direct computations can meet with serious numerical difficulties [1]. This explains those
several attempts which have been made to take advantage from existence of the small dimensionless
parameter in the problem by employing asymptotic analysis [2, 3]. Few solutions are available [2?4]
for monotonic loading where the latter [4] has been extended for the case of an arbitrary time dependent
load in [5]. All these solutions are based on the simplified kinematic relationship in contact zone when
the tangential displacement can be negligible in the analysis.
In real joints, the structural and biomechanical properties of the cartilage layer may change dramatically
in the pathologic stage [5]. Therefore, the refined kinetic relationship which takes into account the tan-
gential displacements of the boundary points of a biphasic cartilage layer should be rather imposed for
correct predictions of biomechanical parameters in the contact region of the articular cartilage instead
of the classic kinetic relationship.
We revisit the original asymptotic formulation of the contact problem and improve it by taking into ac-
count the tangential displacements at the contact interface. Namely, using the well-known asymptotic
relationships for biomechanical parameters valid inside the biphasic cartilage layer [2], we construct a
new more accurate solution which is also capable to represent the limit situation when the biomechan-
ical properties of the contacting cartilage layers differ drastically. We restrict ourself a monotonic load
and essentially use in our analysis the approach proposed earlier in [6] for the case of the classic static
mechanical contact problem. For the details of the original problem formulation, we refer prospec-
tive reader to the paper [2] where, from our best knowledge, the asymptotical approach to the contact
problem for the biphasic cartilage layer has been developed.
2 PROBLEM FORMULATION AND MAIN RESULTS
We consider a thin linear biphasic cartilage layer indented by a spherical punch. The refined linearized
kinetic relationship which takes into account the both normal w(r;t) and tangential u(r;t) displace-
ments of the boundary points of the cartilage layer can be written as follows [6] (see Fig. 1a):
?w(r;t)+ rR
0
u(r;t) = ?0(t)? r
2
2R0; r ? a(t): (1)
Here, ?0(t) is the punch displacement, R0 is the radius of the punch surface, a(t) is the contact radius.
bone
cartilage
impermeable spherical punch
Figure 1: a) - geometry of the problem; b) - support of the solution of Eq. (2).
a) b) t
ra0 = a(0)
r = a(t)
?+??
R0
h ?0(t)
a(t)
Expressing the displacements w(r;t) and u(r;t) in terms of the contact pressure P(r;t) by means of
formulas [2], we reduce Eq. (1) to the following integro-differential equation (we assume that r ? a(t)):
1
r
@
@r
r@P(r;t)@r
?
+?
Z t
0
1
r
@
@r
r@P(r;?)@r
?
d? +?r@P(r;t)@r = m(Cr2 ??0(t)): (2)
Here we used the notation ? = 3?sk=h2, ? = 3=(2hR0), m = 3?s=h3, C = 1=2R0.
The radius of the contact area a(t) is determined from the boundary conditions [2, 4]:
P(a(t);t) = 0; @rP(r;t)
flfl
flr=a(t) = 0; @rP(r;t)
flfl
flr=0 = 0: (3)
Denoting a non-decreasing external load by F(t), we will have the following equilibrium equation:
2?
Z a(t)
0
P(?;t)?d? = F(t): (4)
Notice that because the load F(t) is non-decreasing, the contact radius a(t) increases monotonously.
Notice also that in case ? = 0 the contact problem (2) ? (4) coincides with that studied in [2, 4].
The following exact equation connects the the unknown punch displacement ?0(t) and the contact load
F(t):
?0(t) = Ca
2(t)
2 +
2?
m
F(t)
?a2(t); (5)
while the radius a(t) of the contact area and the contact load F(t) satisfy the equation:
?m
48 Ca
6(t) = F(t)+?
Z t
0
F(?)d? +?
(
a2(t)
4 F(t)+
?
2
Z a(t)
0
?4@?P(?;t)d?
)
: (6)
The contact pressure, P(r;t), will be considered separately in the domains ?? and ?+ (see Figure 1b)).
If we denote by a?1(r) the unknown inverse to the function a(t) then
b(r) = 0; r ? a0; b(r) = a?1(r); a0 ? r: (7)
In the domain ??, integrating Eq. (2) leads to the following result:
?
m
@P?(r;t)
@r = Cr +C0M(t;r)e
??r2=2 + e??r
2=2
r
Z t
0
M(t??;r)(2C???1 +?00(?))d?: (8)
Here we introduced the notation C0 = 2C??1 +?0(0), ? is an arbitrary value 0 < ? < ? and
M(t;r) = 12?i
Z i1??
?i1??
1
s
exp
? ??r2
2(s+?)
?1
?
estds: (9)
In the domain ?+, we have
?
m
@P+(r;t)
@r =
C
r a
2
0 +
1
re
??r22
?
C0M(t;r)+ 1r
Z t
0
M(t??;r)
h2C?
? +?
0
0(?)
i
d?
? 2C?r e??
r2?a20
2
n
M(t;
q
r2 ?a20)+?
Z t
0
M(t??;
q
r2 ?a20)d?
o
? 2Cr
Z r
a0
e??r
2??2
2 ?M(t?b(?);
q
r2 ??2)d?
+ 1r
Z r
a0
e??r
2??2
2 @tM(t?b(?);
q
r2 ??2)b0(?)(C?2 ??0(b(?)))d?: (10)
3 NUMERICAL EXAMPLES AND CONCLUSIONS
In order to illustrate the constructed solution of the refined biomechanical contact problem and compare
it with that from [4], we will adopt for the computations the same material properties and geometrical
parameters of a typical human cartilage reported in [2,4]. Namely we assume that ?s = 0:25 MPa,
k = 2?10?3 mm4N?1s?1, h = 1 mm, R0 = 400 mm. In the sequel we also consider some variation
of this set to analyze the sensitivity of the solution to the variation of material characteristics. In all the
numerical tests, a constant contact load F(t) = F0 = 100 N is applied. The respective results for the
main contact parameters are presented in Fig. 2.0 50 100 150 20010.8
11
11.2
11.4
a(t) [mm]
0 50 100 150 2000.14
0.15
0.16
0.17
t [s]
2? 0(t) [mm]
classic modelnew model
classic modelnew model
0 2 4 6 8 10 120
0.20.4
0.60.8
P(r,t) [MPa]
0 2 4 6 8 10 12
?100
?50
0
r [mm] dP(r,t)/dr [MPa/m]
classic model new model
classic model new model
t=0 s
t=200 s t=100 s
t=0 s
t=100 s t=200 s
Figure 2: Comparison of the radius of the contact zone, a(t), the indentation parameter ?0(t), the contact
pressure P(r;t) and its derivative with respect to the radial coordinate r for times t = 0;100;200 [s]
computed in accordance with the model from [4] and the present one for the set of parameters from [2]
The difference between the obtained solution and the classic one [4] is very small and increases with
time though. Thus, for the contact zone radius, a(t), the maximal discrepancy between the solutions is
of 0:074% for t = 200 s. The punch indentation, ?0(t), differs for these two models on 3:5% and does
not practically change in time. (We depict the doubled punch indentation, 2?0(t), for easy comparison
with the results from [4]). Finally, the maximal pressure difference between the models appears at point
r = 0 and increases from 0:85% and 3:57% with time.0 50 100 150 2007.5
8
8.5
9
a(t) [mm]
0 50 100 150 2000.07
0.08
0.09
0.1
t [s]
2? 0(t) [mm]
classic modelnew model
classic modelnew model
0 2 4 6 80
0.51
1.5
P(r,t) [MPa]
0 2 4 6 8?300
?200?100
0
r [mm] dP(r,t)/dr [MPa/m]
classic model new model
classic model new model
t=0 s
t=200 s
t=100 s
t=0 s
t=100 s t=200 s
Figure 3: Comparison of the contact zone radius, a(t), the indentation parameter ?0(t), the contact
pressure P(r;t) and its derivative, @rP(r;t), for moments t = 0;100;200 [s] computed in accordance
with the model [4] and the present one for h0 = h=2 while others parameters are the same as above.
The second set of the biomechanical parameters chosen for the numerical simulation defers from the
previous one only by the thickness of the cartilage layer h0 = h=2, while all others stay the same.
Now the aforementioned discrepancy between the solution becomes more pronounced (Fig. 3). Thus,
the difference between the contact zone takes its maximal value of 0:6% for t = 200 s. The punch
indentation, ?0(t), differs even further on to reach its maximum 9% for t = 0. Finally, the maximum
pressure difference between the models lies in range between 0:8% and 5% increasing with time.
Summarizing, for the set of parameters discussed above, the difference between the classic solution
and the refined one is not essential, but for the punch indentation and the contact pressure should be
important in some cases though. These cases may include inverse problems of determining biomechan-
ical and the geometrical parameters of the cartilage layer from experiments. Note also, that the refined
contact model always yields higher contact pressure maximum with the longest contact zone and, as a
direct consequence, to a higher pressure gradient.
REFERENCES
[1] W. Wilson, C.C. van Donkelaar, R. van Rietberger, R. Huiskes, The role of computational mod-
els in the search for the mechanical behaviour and damage mechanisms of articular cartilage.
Medical Engng & Physics, 27, 810-826, 2005.
[2] G.A. Ateshian, W.M. Lai, W.B. Zhu, V.C. Mow, An asymptotic solution for the contact of two
biphasic cartilage layers. J. Biomechanics, 27, 1347-1360, 1994.
[3] J.Z. Wu, W. Herzog, J. Ronsky, Modeling axi-symmetrical joint contact with biphasic cartilage
layers ? An asymptotic solution. J. Biomechanics, 29, 1263-1281, 1996.
[4] J.Z. Wu, W. Herzog, M. Epstein, An improved solution for the contact of two biphasic cartilage
layers. J. Biomechanics, 30, 371-375, 1997.
[5] J.Z. Wu, W. Herzog, M. Epstein, Joint contact mechanics in the early stages of osteoarthitis.
Medical Engng & Physics, 22, 1-12, 2000.
[6] I.I. Argatov, Approximate solution of an axisymmetric contact problem with allowance for
tangential displacements on the contact surface. J. Appl. Mech. and Techn. Physics, 45,
118-123, 2004.