Anales AFA - XVI Meeting on Recent Advances of Physics of Fluids and its Applications 57-61
ESTUDIO NUMÉRICO DE LUBRICACIÓN VISCOELASTOHIDRODINÁMICA EN
PRÓTESIS DE RODILLA CON FLUIDO NO NEWTONIANO
NUMERICAL ANALYSIS OF VISCOELASTOHYDRODYNAMIC LUBRICATION IN
KNEE PROSTHESIS WITH NON-NEWTONIAN FLUID
L. E. Robledo Blasco*1, B. A. Weiss1, M. E. Berli1y J. Di Paolo1
1Grupo Biomecánica Computacional (GBC), Facultad de Ingeniería, Universidad Nacional de Entre Ríos,
Ruta 11 km 10, Oro Verde, Entre Ríos, Argentina.
Recibido: 28/10/2021; Aceptado: 17/01/2022
En este trabajo se diseñó un modelo computacional de lubricación de película delgada aplicada a una prótesis de rodilla
con el fin de realizar una comparación entre tres modelos constitutivos, es decir, la ley de potencias, el modelo de
Carreau-Yasuda y el modelo de Cross. El modelo equivalente de la prótesis de rodilla consiste en un cilindro rígido sobre
un plano deformable. El componente deformable que representa la base tibial se asumió como un sólido estándar lineal
(SEL) viscoelástico. Las ecuaciones gobernantes se resolvieron simultáneamente con la determinación de una frontera
libre móvil utilizando el software COMSOL Multiphysics. Los resultados obtenidos mostraron que el modelo Cross
presenta el mayor valor de tasa de corte, el menor espesor de película y la viscosidad dinámica con menor variación a
lo largo del canal de lubricación, alcanzando un valor mínimo de viscosidad de 0,02 Pa.s. El modelo Carreau-Yasuda
presenta el valor más alto de coeficiente de fricción, siendo un 21,5% superior al modelo de Ley de Potencias y un
3,67% superior al modelo de Cross.
Palabras Clave: sólido estándar lineal, modelo unidimensional, modelado computacional.
In this work, a computational model of thin film lubrication applied to a knee prostheses was designed in order to make
a comparison between three constitutive models, i.e. a power law, the Carreau-Yasuda model and the Cross model.
The equivalent model of the knee prosthesis was modeled as a rigid cylinder on a deformable plane. The mechanical
behavior of the deformable component representing the tibial base was assumed as a viscoelastic Standard Linear Solid
(SLS). The governing equations were solved simultaneously with the determination of a free moving boundary by
implementing it in COMSOL Multiphysics software. The results obtained showed that the Cross model presents the
highest shear rate value, the lowest film thickness and the dynamic viscosity with less variation along the lubrication
channel, reaching a minimum viscosity value of 0.02 Pa.s. The Carreau-Yasuda model presents the highest value of
friction coefficient, being 21.5% higher than for the Power Law model and 3.67 % higher than the Cross model.
Keywords: Linear Standard Solid, one-dimensional model, computational modelling.
https://doi.org/10.31527/analesafa.2022.fluidos.57 ISSN 1850-1168 (online)
I. INTRODUCCIÓN
La rodilla es la articulación sinovial más extensa del cuer-
po humano. Cuando el daño sobre las superficies articulares
es severo, la única solución posible es reemplazar las super-
ficies articulares por un conjunto protésico artificial.
Una prótesis de rodilla se compone de un elemento fe-
moral metálico, un elemento tibial metálico y un elemento
tibial polimérico que se apoya sobre el elemento tibial me-
tálico (Fig. 1).
Si bien las prótesis de rodilla ofrecen una mejora en la
calidad de vida de los pacientes, su vida útil no supera los
15 años. Esto, sumado al incremento en el promedio de vi-
da de la población conduciría a un incremento en el núme-
ro de procedimientos de reemplazo de prótesis de rodilla
(artroplastia de revisión). Uno de los factores que definen
la durabilidad de estos implantes reside en la capacidad de
lubricación de las superficies articulares. Los procesos de
lubricación y deformación de la articulación natural y arti-
ficial resultan de gran interés para la comunidad médica e
* lerblasco@gmail.com
ingenieril.
Gran parte de los estudios sobre el comportamiento de
la articulación artificial se basan en modelos computacio-
nales. La mayoría de estos trabajos consideran al fluido si-
novial como un fluido Newtoniano [1], a pesar de que su
comportamiento no Newtoniano está ampliamente probado
tanto para las articulaciones sanas como para las patológi-
cas y periprostéticas [2,3]. Asimismo, la mayor parte de la
bibliografía de referencia adopta modelos simplificados de
deformación elástica para el material polimérico deforma-
ble, a pesar de que el mismo se comporta como un material
viscoelástico [4,5].
El modelado del fluido sinovial y del material poliméri-
co puede resultar en diferencias notables en los resultados
obtenidos [6-9].
Por ello, en este trabajo se aplica un modelo computacio-
nal de lubricación viscoelasto-hidrodinámica en prótesis de
rodilla y se compara el comportamiento de tres modelos
constitutivos de fluido sinovial no Newtoniano.
Es usual la utilización de condiciones de frontera de
Sommerfeld en la resolución del problema de lubricación,
©2022 Anales AFA 57
FIG. 1: Esquema del conjunto articular protésico.
esto conduce a presiones subambientales que no se obser-
van en la realidad. La eliminación de estos resultados anó-
malos se resuelve utilizando condiciones de contorno de
Reynolds, las cuales especifican que la presión alcanza el
valor ambiental con un gradiente nulo sobre la frontera del
dominio, la cual es a priori desconocida y se clasifica como
libre. Es usual la utilización de métodos aproximados para
determinar la localización de esta frontera libre, como así
también es habitual el forzado de las presiones subambien-
tales a cero cuando se aplican condiciones de frontera de
Sommerfeld [6]. En este trabajo se aplican condiciones de
frontera de Reynolds determinando la frontera libre a partir
de un mapeo lineal.
II. MÉTODOS
Modelo Geométrico
Se construyó un modelo unidimensional estacionario de
la articulación artificial considerando al componente femo-
ral metálico como indeformable, mientras que el compo-
nente tibial se lo consideró un material deformable.
El conjunto osteoartromuscular y ligamentario le permite
a la articulación de la rodilla tener principalmente un único
grado de libertad, por lo que en este caso es válido con-
siderar un modelo geométrico 1D, de tipo cilindro que se
desplaza sobre un plano (véase Fig. 2). El radio del cilindro
equivalente R, se obtiene a partir de las curvaturas de las
superficies femoral y tibial, de radios R1yR2, respectiva-
mente:
R=(R1+R2)
R1R2
(1)
La ecuación que define la altura del canal se obtiene al restar
las parábolas que mejor se ajustan a las curvaturas de las
superficies articulares:
h=h0+r2
2R+d(r)(2)
Donde h0es el entrecruzamiento entre los materiales, Res
el radio equivalente y d(r)es el desplazamiento vertical que
sufre el material deformable en cada punto de la superficie.
FIG. 2: Modelo geométrico.
El componente tibial se consideró como un material de-
formable del tipo Sólido Estándar Lineal (SEL) que obede-
ce al comportamiento definido por la expresión 3:
du(r)
dr +E′′
2
η′′Vu(r) = E′′
1+E′′
2
E′′
1η′′Vp(r) + 1
E′′
1
d p(r)
dr (3)
Donde u(r)es la deformación del material polimérico,
p(r)es la presión aplicada sobre el material, Ves la velo-
cidad tangencial media, η′′ = (4/3)ηyE′′
1yE′′
2son los
módulos equivalentes del material definidos como:
E′′
1=E1(1v)
(1+v)(12v)(4)
E′′
2=E2(1v)
(1+v)(12v)(5)
Estos módulos de elasticidad equivalentes corresponden
al modelo de deformación de columna en el que se consi-
dera al componente tibial como un material de baja rigidez
considerando una deformación plana [5] en la que –median-
te el coeficiente de Poisson v- se incluyen las deformaciones
en sentido perpendicular a la dirección de aplicación de la
fuerza de compresión o descompresión.
A partir de estos módulos elásticos se define un módulo
equivalente Eeq que representa el comportamiento elástico
total del modelo. Este módulo equivalente se define como:
Eeq =E′′
1+E′′
2
E′′
1E′′
2
(6)
Para cuantificar la influencia del comportamiento viscoso
frente al comportamiento elástico del material, se define un
número adimensional a partir de parámetros constructivos y
operativos del modelo:
R. Blasco et al. / Anales AFA - XVI Meeting on Recent Advances of Physics of Fluids and its Applications 57-61 58
ND =η′′V
E′′
2p(2RL)(7)
siendo Lel espesor del material deformable.
Modelos constitutivos de fluido no Newtoniano
En este trabajo el fluido sinovial se modela como fluido
newtoniano generalizado, a partir de tres modelos constitu-
tivos. En todos ellos, la viscosidad aparente (µ) es función
de la tasa de deformación ( ˙
γ).
El primero de ellos es la Ley de Potencias, de acuerdo a
la cual, la viscosidad aparente (µ
p) el fluido se define como:
µ
p=m(˙
γ)n1(8)
siendo mel índice de consistencia del fluido y nel índice
exponencial. Dado que el fluido sinovial es pseudoplástico,
ntomará valores menores a 1.
El segundo de ellos es el modelo de Cross, siendo la vis-
cosidad aparente:
µ
Cr =µin f +µ0µin f
1+α˙
γβ(9)
donde µ0es la viscosidad a tasas de corte bajas, µin f es
la viscosidad correspondiente a elevados valores de tasa de
corte, αyβson constates de ajuste del modelo.
El tercer modelo considerado es el de Carreau-Yasuda,
donde la viscosidad aparente se define como:
µ
Ca =µin f + (µ0µin f )(1+ (λ˙
γ)a)(n1)/a(10)
donde µ0yµin f corresponden a los valores extremos de
viscosidad aparente a bajas y altas tasas de corte, λes una
constante temporal, nyason valores de ajuste del compor-
tamiento exponencial del modelo.
Los parámetros constitutivos de los modelos de Cross y
de Carreau se obtuvieron de estudios del comportamien-
to reológico de muestras de fluido sinovial patológico [3],
mientras que los valores de los parámetros constitutivos pa-
ra la ley de Potencias se obtuvieron a partir de los resultados
de los otros dos modelos.
TABLA 1: Valores de los parámetros de diseño.
Descripción Valor
RRadio equivalente del cilindro 0.035 m
E1yE2
Módulos elástico de los ele-
mentos en serie y paralelo del
material deformable
1 GPa
vMódulo de Poisson 0.4
η′′ Viscosidad del amortiguador 2.93 ×107Pa.s
µ0
Viscosidad aparente del fluido
a bajas tasas de corte 0.01 Pa.s
µin f
Viscosidad aparente del fluido
a tasas de corte elevadas 2.00 ×102Pa.s
LEspesor del material deforma-
ble 1.00 ×103m
VVelocidad tangencial prome-
dio 1.91 ×102m/s
WxCarga lineal maáxima 7.35 ×104N/m
yAncho condiíleo 1.5 ×102m
Ecuación de lubricación
Los modelos de fluido no Newtoniano considerados for-
man parte de los fluidos Newtonianos Generalizados, en los
cuales se pueden aplicar las ecuaciones de Navier-Stokes.
Se adoptaron las siguientes hipótesis simplificatorias para
la deducción de la ecuación de lubricación: estado estacio-
nario, frontera de volumen de control fija, fluido isotérmico
e incompresible.
A partir de estas consideraciones, se obtuvo la ecuación
de lubricación de Reynolds para fluidos Newtonianos Ge-
neralizados:
d
dr h3d p
dr =12µVd
dr (h)(11)
Condiciones operativas de la rodilla
Durante la marcha se observa un desplazamiento relati-
vo entre las superficies articulares, donde el movimiento de
flexión-extensión (FE) en la dirección antero-posterior es el
más importante. La velocidad tangencial media (V) se defi-
ne como el promedio de la suma de las velocidades de las
superficies articulares (V1yV2). En este trabajo se adopta
una velocidad tangencial promedio V=0.0191 m/s.
La carga (W) es la fuerza a la que está sometida la ar-
ticulación. En este trabajo, la carga se obtiene como post-
procesado, a partir de la integración de la presión en el do-
minio:
W=YZ
p(r)dr =YWx(12)
donde la integral corresponde a la carga lineal WxyYes la
profundidad cada cóndilo en el eje sagital (profundidad del
cilindro).
Durante la fase de apoyo simple del ciclo de marcha, la
carga aplicada alcanza a ser tres veces el peso corporal de la
persona. Considerando una persona de 75 kg de masa (pe-
so de aproximadamente 735 N), la carga máxima a la que
estará sometida la articulación total será de 2205 N. Consi-
derando una distribución uniforme de la carga entre los dos
contactos cóndilo-tibiales y adoptando un ancho condíleo
Y=0.015 m, la carga lineal (W=Rp(r)dr) aplicada en
cada cóndilo es de 7.35 ×104N/m.
El coeficiente de fricción (φ) se define como el cociente
entre la fuerza de fricción ejercida sobre la superficie y la
carga. Dado que se considera un modelo unidimensional,
en este trabajo se adopta la carga lineal (Wx) y la fuerza de
fricción ejercida sobre la línea de contacto correspondiente
a la superficie deformable. Esta fuerza se obtiene al integrar
los esfuerzos cortantes. El coeficiente de fricción se define,
entonces, como:
φ=R(µ2V
h+d ph
dr2)dr
Wx
.(13)
Técnica de resolución
Para resolver las ecuaciones gobernantes se utilizó el
software comercial COMSOL Multiphysics, ya que el alto
grado de acoplamiento entre las ecuaciones y la no lineali-
dad de las mismas requiere de métodos numéricos para su
resolución en forma simultánea. Debido a esto, se discretizó
R. Blasco et al. / Anales AFA - XVI Meeting on Recent Advances of Physics of Fluids and its Applications 57-61 59
el dominio por medio de elementos cuadráticos de Lagran-
ge.
Para utilizar este software, las ecuaciones gobernantes se
adimensionalizaron y debilitaron. Además, debido a la con-
dición de frontera libre impuesta en el sistema, fue necesa-
rio implementar una estrategia de mapeo de las ecuaciones a
un dominio fijo. Para este trabajo se realizó un mapeo lineal
obedeciendo la Ec. (14). Las ecuaciones fueron resueltas en
el nuevo dominio (de frontera fija y coordenada x, Fig. 3)
y las soluciones fueron luego mapeadas hacia el dominio
físico (de frontera libre y coordenada r)
d2r
dx2=0.(14)
El alto grado de acoplamiento entre las ecuaciones y la no
linealidad de las mismas requieren la realización de una
continuación paramétrica de orden cero en h0, desde cargas
nulas o despreciables (correspondientes a elevados h0) has-
ta alcanzar la carga lineal deseada (Wx=7.35 ×104N/m).
FIG. 3: Conformación de los dominios físico y de referencia utili-
zado durante el proceso de mapeado.
III. RESULTADOS
La altura del canal de lubricación para los casos en es-
tudio se muestra en la Fig. 4. Se observa una diferencia de
altura mínima de 2% entre el modelo de Ley de Potencias
y el de Carreau mientras que la máxima diferencia en este
sentido observada entre el modelo de Cross y de Carreau es
de 9,4%.
La viscosidad a lo largo del canal de lubricación para los
tres modelos se presenta en la Fig. 6. El modelo de Ley de
Potencias presenta la mínima viscosidad aparente.
El modelo de Cross presenta los valores de viscosidad
aparente con menor variación a lo largo del canal de lubri-
cación, alcanzando un valor de viscosidad mínima de 2.00
×102Pa.s.
El modelo de Carreau presenta valores de viscosidad mí-
nima cercanos a 2.3 ×102Pa.s, mientras que el modelo de
Ley de Potencias presenta un valor de viscosidad mínima de
1.95 ×102Pa.s. El descenso de este último valor por de-
bajo de la viscosidad esperada para tasas de corte elevadas
(de 0.02 Pa.s) se debe principalmente a las características
constructivas de este modelo de fluido no Newtoniano, el
cual desciende indefinidamente a medida que la tasa de cor-
te aumenta.
La Fig. 6muestra la tasa de corte a lo largo del canal
de lubricación. Los modelos de Cross y de Carreau poseen
los valores más extremos, mientras que el modelo de Poten-
FIG. 4: Altura del canal de lubricación para los modelos de fluido
no Newtoniano analizados.
FIG. 5: Viscosidad aparente para los modelos de fluido no New-
toniano analizados.
FIG. 6: Tasa de corte a lo largo del canal para los modelos de
fluido no Newtoniano analizados.
cias posee valores de tasa de corte intermedios. Al presentar
mayores tasas de corte, la viscosidad disminuirá debido al
comportamiento no Newtoniano del fluido, tal como se pue-
de ver en la Fig. 5.
Los tres modelos muestran valores de presión máxima
muy similares con una diferencia menor a 1%. Esto se de-
be a que la carga impuesta en los tres casos de estudio es
R. Blasco et al. / Anales AFA - XVI Meeting on Recent Advances of Physics of Fluids and its Applications 57-61 60