Anales AFA - XVI Meeting on Recent Advances of Physics of Fluids and its Applications 36-40
HACIA PRÓTESIS DE RODILLA MÁS DURADERAS: EL MODELO DE PHAN-THIEN Y
TANNER Y LA LUBRICACIÓN ELASTOHIDRODINÁMICA
TOWARDS MORE DURABLE KNEE PROSTHESES: THE PHAN-THIEN AND TANNER
MODEL AND THE ELASTO-HYDRODYNAMIC LUBRICATION
B. A. Weiss*1, S. Ubal1y J. Di Paolo1
1Grupo Biomecánica Computacional (GBC) Facultad de Ingeniería Universidad Nacional De Entre Ríos,
Ruta Prov. 11 km 10 (3100) Oro Verde Entre Ríos Argentina.
Recibido: 9/11/21; Aceptado: 17/11/21
El fluido sinovial que lubrica las prótesis de rodilla usualmente se asume como newtoniano a pesar de ser viscoelástico.
En este trabajo se presenta un modelo basado en la ecuación de Reynolds y en la ley reológica de Phan-Thien y Tanner
(PTT). El modelo, de tipo elipsoide-sobre-plano, contempla la curvatura de las superficies articulares y asume que el
componente tibial, típicamente fabricado en polietileno de ultra alto peso molecular (UHMWPE) se deforma de acuerdo
al modelo de columna. Utilizando el software comercial COMSOL Multiphysics versión 5.2 se resolvieron simultánea-
mente las ecuaciones gobernantes del flujo de lubricación no newtoniano para encontrar la componente del tensor de
tensiones totales para el modelo PTT, en la dirección perpendicular al flujo, y la forma del canal de flujo. A partir de
estos resultados se calcularon las fuerzas de fricción y los coeficientes de fricción sobre las superficies articulares. Se
realizaron simulaciones para distintas condiciones operativas correspondientes a la lubricación elastohidrodinámica,
mecanismo de lubricación predominante en la fase de apoyo simple del ciclo de marcha. Los resultados manifiestan
la importancia de un adecuado modelo para obtener información cuantitativa, considerar al fluido como newtoniano
conduciría a sobreestimaciones del espesor de película lubricante. Los resultados sugieren que, mientras el régimen
de lubricación sea de película delgada, la viscosuplementación y consiguiente aumento de la constante de relajación
λconduciría a menores coeficientes de fricción. Sin embargo, para que el régimen de lubricación corresponda al de
película delgada, la constante de relajación del fluido debería ser inferior a 2.5×105s.
Palabras Clave: fluido no newtoniano, elasticidad, lubricación.
The synovial fluid is usually assumed to be Newtonian despite its viscoelastic behavior. In this work a model based on
the Reynolds equation and the Phan-Thien and Tanner rheological law (PTT) is presented. The model, of the ellipsoid-
on-plane type, contemplates the curvature of the articular surfaces and assumes the tibial component, typically made
of ultra-high molecular weight polyethylene (UHMWPE), is deformed according to the column model. Using com-
mercial software COMSOL Multiphysics version 5.2, the governing equations of non-Newtonian lubrication flow were
simultaneously solved to find the component of the total stress tensor for the PTT model, in the direction perpendi-
cular to the flow, and the lubricating film thickness. From these results the friction forces and friction coefficients on
the articular surfaces were calculated. Simulations were carried out for different operating conditions corresponding
to elasto-hydrodynamic lubrication, the predominant lubrication mechanism in the simple support phase of the gait
cycle. The results show the importance of an adequate model to obtain quantitative information, considering the fluid
as Newtonian would lead to overestimations of the lubricant film thickness. The results suggest that, as long as thin film
lubrication is the lubrication regime, the viscosupplementation and consequent increase in relaxation constant λwould
lead to lower coefficients of friction. However, for the lubrication regime to be thin film lubrication, the relaxation
constant of the fluid should be less than 2.5×105s.
Keywords: non-newtonian fluid, elasticity, lubrication.
https://doi.org/10.31527/analesafa.2022.fluidos.36 ISSN 1850-1168 (online)
I. INTRODUCCIÓN
Los implantes metal-polietileno de ultra alto peso mole-
cular (UHMWPE) son los más populares para el reemplazo
de articulaciones de rodilla. Las consecuencias del desgaste
de las superficies protésicas constituyen el principal factor
limitante de su vida útil [1], la cual resulta insatisfactoria
para pacientes menores de 60 años, que deben someterse
a nuevas intervenciones quirúrgicas (llamadas “revisiones”)
para reemplazar dichos implantes [1]. Las revisiones son ci-
rugías complicadas, riesgosas y costosas; además, los nue-
* brenda.weiss@uner.edu.ar
vos implantes suelen ser menos exitosos que los primeros
[2].
Al igual que las articulaciones naturales, el implante está
lubricado con fluido sinovial, siendo esto un factor funda-
mental para reducir la fricción y el desgaste de las superfi-
cies en contacto.
Este trabajo tiene como objetivo investigar el fenómeno
de contacto lubricado entre las superficies de los compo-
nentes de prótesis de rodilla, con el fin de predecir com-
binaciones de parámetros geométricos, físicos y operativos
que aseguren la existencia del contacto lubricado, minimi-
cen la fricción evitando el contacto sólido-sólido entre las
©2022 Anales AFA 36
partes y de esta forma permitan incrementar la vida útil del
implante.
II. MÉTODOS
La prótesis de rodilla posee dos cóndilos convexos articu-
lando, cada uno de ellos, con un cóndilo cóncavo presente
en el platillo tibial. En este trabajo se considera el modelo
de elipsoide sobre plano [3] para modelar un cóndilo femo-
ral articulando con el correspondiente cóndilo tibial (Fig.
1). A partir de los radios de curvatura de estas superficies
(Rf x yRtx) se obtienen los radios equivalente del elipsoide
en cada una de las direcciones principales: Rxen la direc-
ción anteroposterior o dirección de la marcha (Ec. (1)), y Ry
en dirección medial-lateral (Ec. (2) [4]):
Rx=Rf xRtx
Rtx Rf x
,(1)
Ry=2.844Rx.(2)
FIG. 1: Modelo de elipsoide sobre plano, componentes del vector
velocidad (u,v,w) sistema de coordenadas y dominio (celeste). e es
el espesor del sustrato deformable, Rxy Ryson los radios princi-
pales del elipsoide.
En un estado indeformado (Fig. 2), la altura del canal de
lubricación h será igual a la distancia instantánea entre el
sustrato (fijo) y la superficie del elipsoide (móvil), la cual,
está dada por h(x,y) = h00 +x2/(2Rx) + y2/(2Ry), donde
h00 es la menor distancia entre las superficies indeformadas
(en x=y=0). Adicionando el desplazamiento δzen cada
punto del sustrato deformable debido a la deformación del
mismo, se obtiene la ecuación de altura del canal:
h=h00 +x2
2Rx
+y2
2Ry
+δz.(3)
FIG. 2: Vista del plano y=0, L=2l es la longitud del dominio en
dirección x, l es la longitud del dominio en dirección y (debido a
la condición de simetría), h es el espesor de la película lubricante.
El elipsoide se considera rígido y el sustrato deformable
adopta las características mecánicas del UHMWPE (Tabla
TABLA 1: Parámetros adoptados
Descripción Valor
lSemi-longitud en dirección Rxy longitud en dirección Ry0.02 m
RxRadio equivalente en dirección x 123.7 mm
RyRadio equivalente en dirección y 352 mm
EMódulo de elasticidad del UHMWPE 1 GPa
vCoeficiente de Poisson del UHMWPE 0.4
eEspesor del sustrato deformable 10 mm
ξConstante de deslizamiento del modelo PTT 0.4
η0Viscosidad dinámica del solvente del fluido sinovial 0.1 Pa s
1). En este trabajo se adopta el modelo de columna (modelo
simplificado de elasticidad lineal).
En la fase de apoyo simple del ciclo de marcha se evi-
dencian las máximas cargas en la articulación de rodilla, y
el mecanismo de lubricación predominante es la lubricación
elastohidrodinámica (EHL).
Si bien el líquido sinovial es un fluido viscoelástico y
pseudoplástico, usualmente se lo considera como un fluido
newtoniano o newtoniano generalizado [4-7]. En este traba-
jo se lo modela como un fluido viscoelástico de Phan Thien
y Tanner (PTT). En consecuencia, en estado estacionario,
se obtiene una ecuación de Reynolds modificada, que con-
templa a este tipo de fluidos (más detalles en [8,9]):
xh3
12η πzz
xb1+
yh3
12η πzz
yb2=ue
h
x
(4)
ηη0
1+λ2ξ(2ξ)( b
˙
γ2
x+b
˙
γ2
y)(5)
b1=λ
ξ η b
˙
γ2
y+2ηb
˙
γ2
x
x+ (2ξ)λ η d
˙
γx˙
γy
y(6)
b2= (2ξ)λ η d
˙
γx˙
γy
x+λ
ξ η b
˙
γ2
x+2ηb
˙
γ2
y
y(7)
donde πzz es la componente normal del tensor total de
tensiones, λes la constante de relajación del fluido, ξes
un parámetro del modelo PTT, ηes la viscosidad aparente
del fluido. La velocidad de arrastre ue= (ut+uf)/2 siendo
utla velocidad del componente tibial y ufla velocidad del
componente femoral, ambas en dirección x. Los valores de
uepueden determinarse a partir del ángulo de flexión y el
desplazamiento antero-posterior en el plano sagital. En este
trabajo se utiliza la velocidad de arrastre obtenida a partir
de la norma ISO 14243-3 [4,7], la cual describe un ensayo
de desgaste en prótesis de rodilla.
Dada la no linealidad de las componentes del tensor de
tensiones viscosas, se asume que el flujo viscoelástico con-
siste de un flujo base Newtoniano que sufre una perturba-
ción no-newtoniana permitiendo estimar los cuadrados de
los gradientes de velocidad a través del espesor de la pelí-
cula lubricante a partir de los valores medios de los cuadra-
dos de los gradientes de velocidad newtoniana: b
˙
γ2
x,b
˙
γ2
y,d
˙
γ2
x˙
γ2
y.
Se requiere entonces resolver el problema análogo newto-
Weiss et al. / Anales AFA - XVI Meeting on Recent Advances of Physics of Fluids and its Applications 36-40 37
TABLA 2: Casos de estudio [4].
Tiempo (s) Velocidad de arrastre ue(mm/s) Carga soportada por cóndilo medial (N)
0,03 16.8 1192
0,08 40.9 812
0,13 19.6 1659
0,45 14.9 1558
0,50 57.6 1192
niano (con idénticas condiciones geométricas, operativas y
mecánicas pero considerando al fluido sinovial como New-
toniano).
La capacidad de carga (W) de la película lubricante se
obtiene al integrar la componente normal del tensor de ten-
siones (πzz) en el dominio ():
W=Z
πzzd(8)
Dada la simetría del problema respecto al eje x, las ecua-
ciones gobernantes (Ecs. (3), (4)y(8)) se resolvieron en la
mitad del dominio, utilizando el software comercial COM-
SOL Multiphysics (basado en el método de elementos fini-
tos) y adoptando las condiciones de frontera que se mues-
tran en la Fig. 3. Para facilitar la convergencia, la obtención
de resultados se realiza para cargas gradualmente crecien-
tes, utilizándose como inicialización en cada cálculo los re-
sultados obtenidos para una carga inferior.
FIG. 3: Dominio (, en gris) simétrico respecto a eje x, y con-
diciones en las fronteras: πzz =0en las fronteras Γ1,Γ2yΓ4,
πzz/y=0 (simetria) en Γ3.
Dado que la carga y la velocidad de arrastre varían a lo
largo de un ciclo de marcha, se seleccionaron como casos
de estudio a los instantes correspondientes a máximos de
carga o velocidad, a lo largo de la fase de apoyo simple (ver
Tabla 2).
III. RESULTADOS
La carga aplicada es la misma para el caso PTT y el new-
toniano, por lo que la distribución de presiones y de la com-
ponente normal del tensor de tensiones es idéntica, eviden-
ciando un máximo en el centro del canal. La Fig. 4muestra
la máxima componente normal del tensor de tensiones ob-
tenidos en este trabajo. Se observa una buena concordancia
entre estos resultados y los obtenidos por Su et al [4]. El
efecto de considerar la viscoelasticidad del fluido se apre-
cia, principalmente, en el espesor de la película lubricante.
La Tabla 3presenta el mínimo espesor de película lubri-
cante de fluido PTT (hmin) y de fluido newtoniano (hN
min) pa-
ra los distintos casos de estudio, con λ= 106s, 2.5×105s
FIG. 4: Máxima componente normal del tensor de tensiones
(πzz,max, símbolos llenos) vs tiempo (t), considerando fluido PTT.
En línea llena se presenta la máxima presión (p00 = pmax obteni-
da por Su et al. [4], considerando fluido newtoniano e idénticos
parámetros geométricos, físicos y operativos (Rx, Ry, e, E, v, ue,
W , η0) que en este trabajo.
y 3×105s. En esta tabla el número adimensional de Weis-
senberg (We), número que representa la relación entre las
fuerzas elásticas y viscosas, se calcula a partir del espesor
mínimo obtenido en cada caso, y su correspondiente cons-
tante de relajación (λ) y velocidad de arrastre (ue). En flui-
dos newtonianos We=0.
En la Tabla 3se aprecia que, para todos los tiempos en
estudio, considerar λ= 106s conduce a números de Weis-
senberg del orden de 102y a espesores de película lubri-
cante que difieren menos del 0.3% respecto al correspon-
diente caso newtoniano. Además, la viscosidad aparente al
considerar λ= 106s es superior al 99% de la viscosidad a
tasa de corte cero η0(considerada para el caso newtoniano).
La Fig. 5presenta el espesor de película lubricante (h) a
lo largo de la línea de simetría (y=0) para el caso corres-
pondiente al tiempo t=0.13 s, momento en que tiene la
máxima carga sobre la articulación.
FIG. 5: Espesor de película lubricante (h) vs coordenada x, en
la línea de simetría (y=0), considerando distintas constantes de
relajación del fluido (λ, leyenda) para el caso correspondiente al
tiempo t=0.13 s.
Weiss et al. / Anales AFA - XVI Meeting on Recent Advances of Physics of Fluids and its Applications 36-40 38
TABLA 3: Número de Weissenberg (We) y mínimo espesor de película lubricante de fluido PPT (hmin) y de fluido newtoniano (hN
min)
considerando λ=106,2.5×105y3×105s con η0=0.1Pa.s.
t(ue)λ(s) hmin (µm) hN
min (µm) (hminhN
min)˙
100%
hN
min
We(λue/hmin)ηmin
η0
0.03
1060.49
0.49
-0.09 0.03 0.997
2.5×1050.32 -34.83 1.33 0.323
3×1050.28 -43.34 1.83 0.249
0.08
106 0.79
0.80
-0.19 0.05 0.993
2.5×1050.44 -44.52 2.32 0.183
3×105 0.36 -54.28 3.38 0.134
0,13
106 0.51
0.51
-0.11 0.03 0.996
2.5×1050.33 -35.54 1.27 0.271
3×1050.29 -42.09 1.69 0.205
0.45
1060.45
0.45
-0.07 0.03 0.997
2.5×1050.31 -31.95 1.22 0.332
3×1050.28 -38.57 1.62 0.256
0.5
1060.90
0.91
-0.24 0.06 0.989
2.5×1050.40 -55.95 3.60 0.126
3×1050.37 -58.80 4.62 0.091
En la Fig. 5y en la Tabla 3se aprecia que el aumento de
la constante de relajación del fluido λconduce a la dismi-
nución del espesor de la película lubricante, observándose
un mínimo en la desembocadura del canal. Por lo tanto, asu-
mir al fluido sinovial como fluido Newtoniano conduciría a
estimaciones de espesores de película lubricante mayores a
los que en realidad se producirían.
En la Fig. 6se presentan la altura mínima y central ob-
tenidas con distintas constantes de relajación λ, y los resul-
tados obtenidos por Su et al. [4] al simular la lubricación
elastohidrodinámica dependiente del tiempo. Las grandes
diferencias observadas podrían atribuirse al fenómeno de
aplastamiento de la película lubricante, no simulado en este
trabajo.
FIG. 6: Espesor mínimo (hmin, símbolos vacíos) y central (h00,
símbolos llenos) de película lubricante vs tiempo (t), para distintas
constantes de relajación (λ). Resultados considerando fluido PTT
con η0=0,1 Pa.s. Se presentan el espesor central (línea llena) y
mínimo (línea punteada) de película lubricante obtenidos por Su
et al. [4], considerando fluido newtoniano e idénticos parámetros
geométricos, físicos y operativos (Rx, Ry, e, E, v, ue, W , η0) que
en este trabajo.
El régimen de lubricación depende del espesor de la pe-
lícula lubricante y de la rugosidad equivalente de las super-
ficies articulares. Contemplando la situación más favorable,
con 0.1 µm de rugosidad del UHMWPE, hdebería ser su-
perior a 0.3 µm para que el régimen de lubricación sea de
película delgada. Los resultados presentados en la Tabla 3
sugieren que si λ>2.5 ×105s el régimen de lubricación
correspondería a la lubricación mixta (0.1 µm<hmin <0.3
µm ).
FIG. 7: Coeficiente de fricción a) en dirección x (fx) y b) en di-
rección y (fy) vs tiempo (t), para distintas constantes de relajación
(λ).
Weiss et al. / Anales AFA - XVI Meeting on Recent Advances of Physics of Fluids and its Applications 36-40 39
El aumento de λconduce a una disminución en la míni-
ma viscosidad aparente, alcanzando viscosidades entre 0.01
y 0.03 Pa.s (0.1ηmin/η00.3, Tabla 3) acordes a resul-
tados experimentales [10].
Las Figs. 7y8presentan los coeficientes de fricción ( fx
yfy) y los coeficientes de fricción relativos al caso análogo
newtoniano ( fx/fN
xyfy/fN
y). Cuando λ=106s, se tienen
coeficientes de fricción relativos aproximadamente 1 para
todos los instantes del ciclo de marcha simulado. El aumen-
to de λconduce a la reducción del coeficiente de fricción y
coeficiente de fricción relativo. Estos resultados indicarían
la importancia del comportamiento viscoelástico del fluido
sinovial en la buena performance de las articulaciones si-
noviales jóvenes y sanas (las cuales presentan las mayores
constantes de relajación λ[11]).
FIG. 8: Coeficiente de fricción relativo a) en dirección x ( fx/fN
x)
y b) en dirección y ( fy/fN
y) vs tiempo (t), para distintas constantes
de relajación (λ).
IV. CONCLUSIONES
Asumir al fluido sinovial como fluido Newtoniano podría
conducir a estimaciones de espesores de película lubricante
mayores a los que en realidad se producirían. De aquí la
importancia de un apropiado modelado de este fluido, tanto
en experimentos numéricos como en ensayos de prótesis.
Los resultados sugieren que, mientras el régimen de lu-
bricación sea de película delgada (lo cual sucede con λ
2.5×105s), la viscosuplementación y consiguiente au-
mento de λconduciría a menores coeficientes de fricción.
Las diferencias observadas al comparar los espesores
central y mínimo con los resultados de Su et al. [4] podrían
atribuirse a la lubricación por aplastamiento de la película
lubricante, no simulada en este trabajo. La resolución del
modelo completo, incluyendo el aplastamiento y la lubrica-
ción elastohidrodinámica, será explorada en futuros traba-
jos para corroborar esta hipótesis.
REFERENCIAS
[1] P. F. Sharkey, W. J. Hozack, R. H. Rothman, S. Shastri
y S. M. Jacoby. Why Are Total Knee Arthroplasties Failing
Today? Clinical Orthopaedics and Related Research 404,
7-13 (2002).
[2] D. J. Paulo. Biotribology (Wiley, London, 2013).
[3] Z. M. Jin, D. Dowson y J. Fisher. Contact Pressure Pre-
diction in Total Knee Joint Replacements Part 1: General
Elasticity Solution for Elliptical Layered Contacts. Procee-
dings of the Institution of Mechanical Engineers, Part H:
Journal of Engineering in Medicine 209, 1-8 (1995).
[4] Y. Su, P. Yang, Z. Fu, Z. Jin y C. Wang. Time-dependent
elastohydrodynamic lubrication analysis of total knee re-
placement under walking conditions. Computer Methods
in Biomechanics and Biomedical Engineering 14, 539-548
(2011).
[5] J. D. Paolo y M. E. Berli. Numerical analysis of the effects
of material parameters on the lubrication mechanism for
knee prosthesis. Computer Methods in Biomechanics and
Biomedical Engineering 9, 79-89 (2006).
[6] M. Mongkolwongrojn, K. Wongseedakaew y F. E. Ken-
nedy. Transient elastohydrodynamic lubrication in artificial
knee joint with non-Newtonian fluids. Tribology Interna-
tional 43, 1017-1026 (2010).
[7] ISO 14243-3. Implants for Surgery - Wear of Total Knee-
Joint Prostheses - Part 3: Loading and Displacement Para-
meters for Wear-Testing Machines with Displacement Con-
trol and Corresponding Environmental Conditions for Test.
[8] B. A. Weiss, B. Bou-Said, S. Ubal y J. D. Paolo. The
Phan-Thien and Tanner Model Applied to the Lubrication
of Knee Prostheses. Journal of Biomechanical Engineering
141 (2019).
[9] B. A. Weiss, S. Ubal, M. E. Berli y D. P. José. Investigacio-
nes teórico-computacionales para el desarrollo de prótesis
articulares de mayor duración para miembro inferior. Su-
plemento Ciencia, Docencia y Tecnología 10 (2020).
[10] B. Petcharatana. Rheology of synovial fluid with and
without viscosupplements in patients with osteoarthritis:
a pilot study Tesis doct. (University of British Columbia,
2009).
[11] V. O. Ribitsch. Viscoelastic behaviour of synovial fluids and
artificial replacements 1.aed. (Springer - Verlag New York,
1990).
Weiss et al. / Anales AFA - XVI Meeting on Recent Advances of Physics of Fluids and its Applications 36-40 40