Towards a Mathematical Model for the Solidification and Rupture of Blood in Stenosed Arteries

by   Fatima Abbas, et al.
Lebanese University

In this paper, we present a mathematical and numerical model for blood solidification and its rupture in stenosed arteries. The interaction between the blood flow and an existing stenosis in the arterial wall is modeled as a three dimensional fluid-structure interaction problem. The blood is assumed to be a non-Newtonian incompressible fluid with a time-dependent viscosity that obeys a modified Carreau's model and the flow dynamics is described by the Navier-Stokes equations. Whereas, the arterial wall is considered a hyperelastic material whose displacement satisfies the quasi-static equilibrium equations. Numerical simulations are performed using FreeFem++ on a two dimensional domain. We investigate the behavior of the viscosity of blood, its speed and the maximum shear stress. From the numerical results, blood recirculation zones have been identified. Moreover, a zone of the blood of high viscosity and low speed has been observed directly after the stenosis in the flow direction. This zone may correspond to a blood accumulation and then solidification zone that is subjected to shear stress by the blood flow and to forces exerted by the artery wall deformation. Therefore, this zone is thought to break and then to release a blood clot that leads to the occlusion of small arterioles.



There are no comments yet.


page 11

page 12

page 19

page 20


Non-Newtonian and poroelastic effects in simulations of arterial flows

In this paper, we focus on investigating the influence on hydrodynamic f...

A new two-dimensional blood flow model with arbitrary cross sections

A new two-dimensional model for blood flows in arteries with arbitrary c...

Proper tangential vector field for the wall shear stress evaluation in numerical blood flow simulations in human carotis tree

In numerical simulations of the cardiovascular system the realistic depi...

HEMELB Acceleration and Visualization for Cerebral Aneurysms

A weakness in the wall of a cerebral artery causing a dilation or balloo...

Periodic-corrected data driven coupling of blood flow and vessel wall for virtual surgery

Fast and realistic coupling of blood flow and vessel wall is of great im...

A direct Eulerian GRP scheme for a blood flow model in arteries

In this paper, we propose a direct Eulerian generalized Riemann problem ...
This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1. Introduction

Cardiovascular diseases, mainly due to atherosclerosis, are causes with the highest percentage leading to death worldwide. Curiosity of finding cures for these diseases has enthused mathematicians to study them from their mathematical viewpoint. Consequently, many computational techniques and models have been developed. These tools aim at describing the blood flow and are effictive in studying the response of the arterial wall under certain conditions, the characteristics of the blood components in addition to those of the heart [3, 16, 17, 21, 22, 24]. Indeed, mathematical models with numerical analysis and simulations play an important role in providing knowledge and insights that are unnoticeable clinically. The rheological behavior of blood is captured by deriving constitutive models that constitute a constructive tool in the diagnosis of the pathologies, investigating appropriate remedies and proposing preventive therapies [1, 14]. Recently, the lumen-wall modeling has been adopted using fluid-structure interaction (FSI) model. In the FSI, the behavior of the blood and the arterial wall are taken into consideration, so that one is capable of representing them by their appropriate dynamics and models. A most commonly used method when dealing with FSI systems is the Arbitrary Lagrangian-Eulerain (ALE) method [8] that is effective when combining the fluid formulation in the Eulerian description and the structure formulation in the Lagrangian description. An overview of FSI in biomedical applications has been considered in the book [2]. In particular, modeling of cardiovascular diseases has been highlighted in [9, 18, 25]. In lumen-wall modeling some difficulties encounter due to the complexity of the arterial wall formed of several layers, each with its own unique mechanics and thickness. Assuming that the arterial wall is negligibly thin, or the ratio of the arterial wall thickness to the aretry raduis is small, reduced shell or membrane models have been employed [4, 13]. Introduction of computational model with FSI in order to investigate the wall shear stresses, blood flow field and recirculation zones in stenosed arteries have been studied in [3] where the blood is considered to be an incompressible Newtonian fluid, whereas in the case of a compressible non-Newtonian fluid these factors have been analyzed in [6].

Currently, numerous computational models are simple when describing the cardiovascular diseases and the related processes such as inflammations, coagulation, plaque growth, clot formation, etc.. Indeed, they are managed to capture only some essential features of the processes that take place in the cardiovascular system. Further, these models neglect some of the blood components with their characteristics as well as the arterial wall layers and their own mechanics. Consequently, more suitable models are needed through which the physiological parameters associated to the blood, atherosclerosis and clots must be investigated clinically. In addition, the arterial wall must be considered as a multi-component structure taking into account the effects of the plaque growth and the clot formation on their mechanism. Further, timescales of the biological phenomena, the pulse duration and the time of the plaque growth must be analyzed.

In the present work, the interaction between the blood modeled by a modified Carreau’s model and the hyperelastic incompressible arterial wall has been considered. In a first step, we introduce the FSI system which is composed of the incompressible Navier-Stokes equations representing the blood flow dynamics and the quasi-static equilibrium equations describing the elastic large deformation of the arterial wall. In addition, coupling conditions that ensure a global energy balance of the FSI system have been imposed on the common interface. Variational formulation has been presented and its discrete formulation has been derived. In fact, the Navier-Stokes equations have been semi-discretized in time, while, the nonlinear material equilibrium equations have been solved using the Newton-Raphson method. Numerical simulations have been performed using FreeFem++ [12]. A deep analysis has been made for better understanding of the behavior of the blood flow, the blood viscosity, the maximum shear stress and the recirculation zones. Based on the numerical results, a location where the blood is thought to accumulate and solidify has been identified. Finally, the factors affecting this zone have been investigated, in particular, this zone is subjected to forces exerted by the artery wall and the blood flow. Consequently, numerical simulations for the deformation of this zone have been performed in order to understand its rupture and thus the release of a clot that will lead to the occlusion of small arterioles.

2. The Fluid-Structure Interaction Problem

The total domain representing the artery in the actual configuration at time is composed of two sub-domains and representing the lumen of the artery and the arterial wall, respectively.

We denote by the reference configuration of the structure of density . Its deformation is described by the displacement field that satisfies the quasi-static incompressible equilibrium equations. The evolution of the structure domain can also be given by the deformation map defined in terms of the displacement as . Its deformation gradient

which is a second order tensor is given by

. Its associated Jacobian is .

On the other hand, we describe the blood flow dynamics by the incompressible Navier-Stokes equations on the sub-domain . We denote by

the velocity of the blood and its pressure, respectively. Further, the blood is assumed to be an incompressible fluid with a constant density .
The sub-domain of moving boundaries evolves from some reference configuration according to an ALE map given by


that is, .
The ALE map is considered to be an extension of the displacement of the interface , that is to say


The operator stands for an extension of the displacement of the boundary . Possible extensions can be found in [19, Section 5.3], [5, Chapter 2] (harmonic, biharmonic, wislow, etc.). In particular we consider the harmonic extension as we will see in Subsection 2.1. The deformation gradient associated to is defined by where the symbol indicates the gradient with respect to the variable . Its Jacobian is .

Here and throughout the context, denotes the displacement of the domain which we set to be . Formulating the Navier-Stokes equations in the ALE frame results a new variable that describes the velocity of the domain . It is related to the displacement by the relation . It is worth to point out that . One must distinguish between the physical velocity of the particles and the velocity of the fluid domain .

In what follows, we refer to the elements in the reference configuration by ” ”. In fact the velocity and the pressure of the blood are given on the reference configuration by


The Cauchy stress tensor is expressed in terms of the strain tensor as


where represents the blood viscosity that will be detailed in the sequel. In the reference configuration , the stress tensor is given by

The arterial wall is assumed to be a hyperelastic material then it is characterized by the existence of an energy density function such that the first Piola-Kirchhoff stress tensor . Further, due to the incompressible behavior of the material its Piola-Kirchhoff stress tensor is modified to the form

The variable , called the hydrostatic pressure, plays the role of the Lagrange multiplier associated to the incompressibility condition .

On the fluid domain , a volumetric force is applied. Moreover, a velocity is enforced on the inlet of the artery . On the contrary, a free-exit condition given by is enforced on the outlet .
On the other hand, a volumetric force is applied on the structure domain which is assumed to be fixed on the boundary , that is to say, on .
On the interface , surface forces and are exerted from the fluid domain and the structure domain, respectively.

The FSI model describing the blood-wall interaction is obtained by the coupling between the incompressible Navier-Stokes equations which are formulated in the ALE frame and the quasi-static incompressible elasticity equations formulated in the Lagrangian frame on the reference configuration . The FSI system is




such that


where and are given by (3). Further, is the transformation of to the reference configuration. Figure 1 illustrates a 3D model of an artery in the actual configuration including the boundaries.

Figure 1. A 3D model of an artery in the actual configuration.
Remark 2.1.

From Expression (2) we get that the ALE map and the structure deformation coincide on the interface , that is to say,

Remark 2.2.

Due to the incompressibility condition the -th component of is

where is the Kronecker symbol. The shear stress components are and , whereas and are the normal stress components.

In a two dimensional space the maximum shear stress- an effective parameter in studying the forces exerted on a fluid- is given by the expression [24]


The variational formulation associated to System (5) is

, , ,

, ,

such that




for all and .
The coupling conditions on the interface are given in the strong form as


The spaces and are respectively

2.1. The Discrete Variational Formulation of the FSI Problem

The variational formulation (7)-(10) of the FSI stands for the incompressible homogeneous Navier-Stokes equations coupled with the quasi-static incompressible equilibrium equations. We assume that no external forces are exerted on neither the fluid domain nor the structure domain, i.e, and . Consider a time step and finite element partitions and for the fluid and the solid sub-domains respectively, of a maximum diameter denoted by . Our aim is to approximate the solution at time , for , in the finite element spaces. The approximation of the solution at time is denoted by .

Semi-Discretization in Time of the Fluid Sub-Problem

In order to guarantee the existence and uniqueness of the solution of the discrete fluid sub-problem when performing the numerical simulations, we use the penalty method [11]. This method consists of replacing the natural weak formulation by a regular one by adding a term multiplied by a sufficiently small parameter . Indeed, writing the modified formulation in a matrix form results a positive definite matrix, which assures the existence and the uniqueness of the solution of the discrete sub-problem. The weak formulation associated to the Navier-Stokes equations obtained upon adding a negligible parameter is then semi-discretized in time, that is, the convective term and the viscosity are considered at the instant , whereas other terms are considered at time . The discrete formulation reads


where the non-linear convective term can be approximated by [12, Section 9.5, p. 267]


Notice that, since the strain rate tensor is symmetric, then we have which gives

Newton-Raphson Method for the Structure Sub-Problem

Regarding the structure sub-problem, at the time iteration we will solve the non-linear problem (9) using Newton-Raphson method. The variational formulation corresponding to the structure sub-problem at the iteration is


The method depends on linearizing the structure sub-problem (13) with respect to the unknowns and . We start initialization by considering a suitable choice of the initial values (). In particular, we link the iterations of the Newton-Raphson method with the time iteration by considering and . Then, we solve iteratively the obtained system corresponding to the Newton-Raphson method until its solution converges to a solution of the non-linear System (9). To ensure the existence of the solution of the structure problem (13) we use the penalty method by modifying System (13) through adding the penalized term with . For simplicity of notation, in what follows we omit the subscript of the deformation , that is, we write .

We proceed to derive the formulation of the structure sub-problem corresponding to the Newton-Raphson method. Let us define the following space

Given , a tolerance and

we construct iteratively the two sequences and by solving for the following system:

Set . Repeat: for , while , find (,) in satisfying


for all .

Set and .

When the condition is fulfilled then convergence of the Newton-Raphson method is achieved. Thus, the solution of the structure sub-problem (13) is given by , for the last value of for which the Newton-Raphson method converges.

Space Discretization of the FSI Problem

Space discretization of the variational formulation is carried out using the finite element method (FEM) [10]. We consider the two finite element spaces associated to the fluid weak formulation

and those associated to the structure weak formulation

where , , and are finite dimensional subspaces. The functional spaces associated to the velocity and displacement fields are considered to be , whereas those associated to the pressures (fluid and hydrostatic) are considered to be . In what follows, all terms are discretized in space as mentioned above, so that the approximation of solution in finite element spaces is verifying (11) and (14).

Finally, the discrete variational formulation reads:

Given and a tolerance , find such that


where the non-linear convective term is approximated by the Expression (12).
The coupling conditions on are


Fix . Set and . Repeat: for , find (,) in satisfying