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.