1. Introducción
En la actualidad el interés sobre la simulación numérica del comportamiento mecánico de la arteria aorta ascendente radica en que proporciona información relevante para la clínica.
La onda de presión es un indicador del estado de salud del sistema circulatorio en general, la salud de las arterias, el corazón, el riñón y la sangre misma pueden monitorearse a través de la onda de presión.
El potencial del modelado numérico radica no sólo en la comprensión del fenómeno en casi su totalidad, ya que a partir de la observación del fenómeno se proponen modelos que explican el mecanismo de funcionamiento del proceso, sino que también surgen variables poco evidentes que pueden ser nuevos indicadores del estado de salud de las arterias, como son la verticidad y la turbulencia en el flujo sanguíneo, así como el trabajo cardiaco en el corazón y los diagramas fase en el monitoreo de la presión arterial. Numerosos estudios [1, 2, 3, 4], han intentado reproducir esta onda de presión utilizando diferentes hipótesis sobre el comportamiento del flujo y de la pared arterial, se han implementado también una gran cantidad de métodos numéricos para la solución tanto del flujo sanguíneo, del medio continuo y su interacción. Diferentes modelos también se han usado para representar la naturaleza hiperelástica de la pared arterial [5, 6].
Este trabajo reproduce la onda de presión proponiendo que la mayor influencia sobre su forma es precisamente el comportamiento mecánico de la arteria aorta, se utiliza un fluido viscoso y se emplea el método de elementos finitos para la solución de ambos modelos junto con su interacción. Las principales propuestas para el estudio del flujo sanguíneo y su interacción con la pared arterial son las siguientes:
- Para el análisis de fuerzas se utiliza el modelo de cadena generalizado propuesto por Quarteroni [7], con un material de Ogden hiperelástico isótropo.
- Los campos de presión que forzan el modelo de la pared arterial se obtienen de la solución numérica de las ecuaciones de Navier-Stokes bajo el esquema de Chorin-Teman.
- La interacción fluido-estructura se resuelve mediante el método ALE (Arbitrary Lagrangian Eulerian).
Las condiciones iniciales y de frontera, así como las propiedades de la arteria y de la sangre que se usan en los modelos a resolver son las que se reportan para una persona sana en promedio, obteniendo con esto que las soluciones son también las reportadas para una persona sana en promedio. Así pues, en este trabajo se logra reproducir, con las consideraciones mencionadas, la curva de presión sanguínea de la arteria aorta ascendente.
2. Metodología
2.1. La pared arterial
Los materiales isótropos han mostrado ser efectivos para modelar el tejido blando, y en la industria se usan de manera frecuente precisamente para obtener modelos de diferentes materiales como hules y cauchos [8,9], sin embargo, dentro de la gama de materiales isótropos existentes se encuentra el modelo de Ogden, el cual se utiliza en este trabajo con un reajuste en sus parámetros para hacer rígido el material a pequeñas deformaciones, haciendo con lo anterior que el modelo se ajuste de manera satisfactoria al material de la pared arterial.
Los materiales que son capaces de recuperar su forma original después de sufrir una deformación debido a la acción de una fuerza son llamados elásticos. Algunos de estos materiales son del tipo hiperelásticos. Los materiales hiperelásticos, se encuentran caracterizados por la expresión de su función de energía libre de deformación, W, la cual describe como se almacena la energía en el cuerpo. Los modelos teóricos que describen el comportamiento mecánico de los materiales elásticos se basan en el estudio de sus geometrías y de acciones de fuerza básicos como, compresión y tracción uniaxial, tracción cortante, etc. En el 2004 Ogden [5], propone un modelo basado en las deformaciones principales del material y se expresa de la siguiente manera:
donde αi, µi y λi son las constantes de ajuste del material y el segundo término que incluye al determinante del tensor de deformación J = det(F) representa la relación de volúmenes entre la condición deformada e inicial.
2.1.1. El modelo de la pared arterial
La geometría para el estudio de la arteria aorta ascendente se asumirá como un cilindro recto y se considerarán las siguientes suposiciones:
- Espesor delgado y esfuerzo plano,
- Geometría cilíndrica y desplazamiento radial,
- Gradientes de pequeñas deformaciones,
- Incompresibilidad de la pared arterial.
La configuración en el tiempo t de la superficie del cilindro está dada por:
Se indicará con n el vector normal a la superficie de la arteria y se tomó un diferencial de superficie dS, para el análisis de fuerzas sobre la arteria, las cuales son:
Fuerzas envolventes de la arteria. Sobre toda la pared arterial se ejerce una presión externa constante Pext, donde la fuerza resultante actuando sobre dS es simplemente:
Fuerzas del fluido. Las fuerzas que ejerce el fluido sobre la pared arterial son representadas por el esfuerzo de Cauchy sobre la pared. Si se indica con Tf el tensor de esfuerzo para el fluido, tenemos:
Considerando los esfuerzos radiales σθ, con un término de amortiguamiento, los esfuerzos longitudinales σz, con un término de viscoelasticidad y haciendo el análisis de esfuerzos obtenemos el modelo a resolver, que recibe el nombre de modelo de cadena generalizado [10], por tomar en cuenta un pre-esfuerzo en el eje longitudinal de la arteria.
2.2. El flujo sanguíneo
La sangre es una solución compuesta por un 45% de elementos sólidos, (glóbulos rojos, glóbulos blancos y plaquetas), y del 55% de agua, lo que hace que la sangre sea cinco veces más viscosa que el agua, y el hecho de que haya mayor cantidad de agua en este compuesto hace que la sangre se considere un flujo incompresible, y es necesario considerar estas características a la hora de realizar el modelado numérico del flujo sanguíneo [11].
Las principales cantidades que describen el flujo sanguíneo son la velocidad u y la presión. Para el cómputo numérico unas de las principales dificultades que se presentan son las siguientes:
El flujo es transitorio. - El flujo sanguíneo es pulsátil, esto es depende del tiempo, y sabiendo entonces que es pulsátil, se aproxima el comportamiento como un fenómeno armónico. El ciclo de este movimiento armónico se compone de dos partes:
- La apertura de las valvas aórticas, con la correspondiente eyección de fluido, que corresponde a la sístole y es donde se encuentra el máximo del flujo.
- La diástole que corresponde a la contracción de las valvas aórticas con el mínimo del flujo.
Interacción del fluido con la pared arterial. - En las arterias mayores este aspecto es importante. Por ejemplo, en la aorta la fuerza ejercida por el flujo sobre la pared arterial produce un desplazamiento de entre un 5% y un 10% en la diástole y la sístole respectivamente [12]. Este es un desplazamiento significativo, ya que éste afecta al fluido. La interacción permite también la propagación de la onda del pulso de presión. Una observación obvia pero importante es que la carencia en el desplazamiento de la pared arterial indicaría:
2.2.1. El modelo del flujo sanguíneo
Una cuestión que debe ser considerada a la hora de pretender hacer un modelo es el incremento en el radio de la arteria, el cual puede variar hasta en un 10 %, debido a las fuerzas ejercidas por el flujo sanguíneo. Cuando el flujo interactúa con una porción de arteria se calcula la solución del flujo en un dominio computacional, Ωt, que varía con el tiempo, como se ilustra en la figura 1, en la cual, Las paredes Γtw están moviéndose, mientras que la sección correspondiente a la entrada y la salida, Γtin y Γtout, permanecen fijas.
La frontera Ωt se divide en dos partes. La primera parte coincide con las fronteras físicas del fluido, i.e., la pared arterial. Que se denotan por Γtw y que son las que se mueven bajo la acción del fluido. La otra parte denotada por ∂Ωt corresponde a fronteras ficticias o artificiales, las cuales delimitan la región de interés.
En el caso del presente trabajo, estas entradas ficticias son la entrada y salida del fluido de la sección de la arteria a estudiar, y están indicadas por Γtin y Γtout.
El problema queda planteado entonces:
Sea el dominio del fluido Ω ∈ ℜ3. Dados ud(t), tn(t) y u0(x), encontrar un campo vectorial u(x, t):Ω × [0, T] → ℜ3 y un campo escalar p(x, t):Ω × [0, T] → ℜ tales que: .
Existen numerosos esquemas utilizados para la solución numérica de las ecuaciones de Navier-Stokes, Váldes [4], por ejemplo, propone el método α − generalizado de pasos fraccionarios para la convergencia de las ecuaciones; por otro lado, Quarteroni [10], implementa el esquema de Yosida para el mismo fin. En este trabajo se usa el esquema de Chorin-Temam.
Para solucionar la interacción entre la estructura y el fluido se usa el método ALE, el cual consiste en agregar un dominio de referencia extra al dominio material, propio de la formulación lagrangiana, y al dominio espacial, propio de la formulación euleriana.
2.3. Construcción del modelo
2.3.1. La pared arterial
La geometría.- Una sección de arteria sana cuenta con las características mostradas en la tabla 1, según la American Heart Association, y que son las que se usaron en el modelo implementado en presente trabajo, con una longitud de 5 cm.
Pared Arterial | Unidades | Dominio 1 |
---|---|---|
Densidad (ρ) | [Kg/m3] | 960 |
Diámetro(d) | cm | 2.5 |
Espesor | cm | .27 |
Modelo del Material | Hiperelástico | |
Esfuerzo Cortante Inicial (µ) | Pa | 6204106 |
Esfuerzo Total (κ) | Pa | 20*6204106 |
Ecuaciones dictaminantes.- En el caso de la pared arterial, la conservación de momentum genera el sistema de ecuaciones dictaminantes que resuelven el problema desde un punto de vista mecánico:
donde la primer integral representa el trabajo interno y δEij es la variación del tensor de deformación de Green-Lagrange, Sij es el segundo tensor de esfuerzos de Piola-Kirchhoff y Ω0 es el dominio de análisis de la arteria. La segunda integral representa el trabajo cinético o inercial donde δui es el vector de los desplazamientos virtuales y üi es la aceleración. El lado derecho de la ecuación 4 se denomina trabajo externo donde el término ti − 0 hace referencia a las fuerzas superficiales que actúan sobre la estructura.
En la tabla 2, se muestran las características utilizadas para la generación del modelo simulado en elemento finito.
Número de grados de libertad | 12948 |
Número de puntos de malla | 4316 |
Número de elementos | 18447 |
Elementos tetraédricos | 18447 |
Elementos tipo prisma | 0 |
Número de elementos de contorno | 3924 |
Elementos triangulares | 3924 |
Elementos tipo cuadrilátero | 0 |
En la figura 2 se muestra el mallado del modelo generado.
En la figura 3 se muestran los campos de velocidades y esfuerzos de Von Mises.
2.3.2. El flujo sanguíneo
Las propiedades físicas de la sangre que se usaron en el modelo se muestran en la tabla 3.
Flujo Sanguíneo | Unidades | Dominio 2 |
---|---|---|
Densidad (ρf) | [Kg/m3] | 1060 |
Modelo del Material | Flujo laminar incompresible | |
Viscosidad | [Ns/m2] | .005 |
Presión de entrada | [Pa] | 11000 |
Presión de salida | [Pa] | 0 |
Ecuaciones dictaminantes.- El sistema a resolver surge tanto de la ecuación de conservación de momentum como de continuidad, dadas en la forma mostrada en la ecuación 5:
A partir de la cual se desarrolla la formulación euleriana usada en la dinámica de fluidos. Para poder resolver el sistema anterior se "debilita" multiplicándola por una función de prueba e integrándolo, quedando de la manera presentada en las Ec. 7 y 8.
Cabe mencionar que son las ecuaciones de Navier-Stokes para flujo incompresible y en donde se utilizan las condiciones homogéneas de contorno de Dirichlet:
2.4. Resultados
La función de impulso es una función senoidal, obtenida de Vasava [13], que toma la siguiente forma:
Las presiones alcanzadas se visualizan en 3D en la figura 4.
Las presiones alcanzan un máximo de hasta 16,500 pascales esto es 123.76 mmHg y con un mínimo de 10,500 pascales equivalentes a 78.75 mmHg que son las que se reportan en la clínica para una persona sana en promedio, La curva de presión entonces se muestra en la figura 5.
En la figura 6 se muestra el desplazamiento de la pared arterial obtenida como resultado del modelo implementado. Se puede apreciar la similaridad cualitativa con la curva de presión reportada en la clínica, se distingue la onda dicrota, para la cual hasta el momento no hay consenso de su origen, pero que es evidente en los experimentos que se trata de un rebote de la elasticidad de la arteria y no un reflujo por vacío o un rebote por el choque con la bifurcación de la arteria.
3. Conclusiones
El flujo sanguíneo se considera un flujo laminar, incompresible y que cumple con las condiciones de no deslizamiento en las paredes de la arteria, encontrando con esto las presiones reportadas en la clínica y además la forma de la curva de presión es cualitativamente similar a la reportada en la literatura médica.
Se utiliza el método de los elementos finitos para resolver las ecuaciones de Navier-Stokes mediante el esquema de Chorin-Temam y con las condiciones iniciales y de frontera propios del flujo sanguíneo para una persona en promedio sana, para la pared arterial se propone el modelo de Ogden con las propiedades reportadas en la literatura especializada, para una persona sana en promedio.
La interacción entre el fluido y la pared arterial se resuelve mediante el método ALE, obteniendo con todo lo anterior una curva de presión cualitativa y cuantitativamente similar a la encontrada en la clínica, esto es, las presiones calculadas corresponden a las presiones reportadas en la clínica y se observa claramente en la forma de la onda depresión dícrota que caracteriza a la onda de presión y que se relaciona con la rigidez de la arteria, es decir, la magnitud y la forma de la incisura dícrota indica de manera directa el estado de salud de la arteria.
Los resultados obtenidos proporcionan visualizaciones de diferentes escenarios clínicos, proporcionando así una herramienta de comprensión del fenómeno para cada caso. Mediante esta simulación se resalta entonces la importancia de la hipótesis sobre que el origen de la onda dícrota radica en la naturaleza elástica de la pared arterial. Un modelo de interacción unidimensional resuelto de forma monolítica tiene un bajo coste computacional y se puede usar para acoplarlo a los modelos tridimensionales de las arterias en el extremo final.