On the Implementation of a Non-linear Viscoelastic Model into Coupled

Blood Flow-biochemistry Model

Tom´aˇs Bodn´ar

and Ad´elia Sequeira

Department of Technical Mathematics, Faculty of Mechanical Engineering, Czech Technical University in Prague,

Karlovo N´amˇest´ı 13, 121 35 Prague 2, Czech Republic

Department of Mathematics and CEMAT, Instituto Superior T´ecnico, Technical University of Lisbon,

Av. Rovisco Pais, 1049–001 Lisbon, Portugal

Keywords:

Blood, Viscoelastic Model, Coagulation, Finite-volume, Runge-Kutta.

Abstract:

This paper presents selected numerical results obtained using a macroscopic blood coagulation model coupled

with a non-linear viscoelastic model for blood ﬂow. The governing system is solved using a central ﬁnite-

volume scheme employing an explicit Runge-Kutta time-integration. An artiﬁcial compressibility method is

used to resolve the pressure. A non-linear TVD ﬁlter is applied for stabilization. A simple test case of blood

ﬂow over a clotting surface in a straight 3D vessel is solved. This work merges and signiﬁcantly extends our

previous studies (Bodn´ar and Sequeira, 2008) and (Bodn´ar et al., 2011a).

1 INTRODUCTION

Mathematical and numerical modeling of blood re-

lated phenomena is a very challenging problem.

The ﬂow of blood is difﬁcult to solve mainly

due to its complex rheological behavior which can

be highly non-Newtonian in certain ﬂow regimes

(Robertson et al., 2008; Robertson et al., 2007; Galdi

et al., 2008). Phenomena like shear-thinning vis-

cosity, viscoelasticity or thixotropy can be observed.

Corresponding models for blood have been developed

including several of these features. To date, there is

no single, generally accepted rheological model for

blood. The modeling is done on case to case basis de-

pending on ﬂow conditions and actual needs in pre-

dicting the blood behavior. The work presented here

follows on previous studies (Bodn´ar and Sequeira,

2010; Bodn´ar et al., 2009; Bodn´ar et al., 2011b) de-

scribing the shear-thinning and viscoelastic behavior

of blood in simple geometries. Even if the blood

ﬂow problem is considered as a purely mechanical

phenomena, its mathematical modeling and numer-

ical simulation is a subject of several serious chal-

lenges as demonstrated e.g. in (Janela et al., 2010;

Gambaruto et al., 2011) or (Pirkl et al., 2011).

The rheology of blood as well as its ﬂow is heav-

ily affected by the underlying microstructure and

is closely related to the biochemistry of reacting

blood constituents. Blood coagulation has been rec-

ognized for a long time as one of the most com-

plex problems in biology as described recently e.g.

in (Fasano et al., 2012). There has been several

blood coagulation models developed in the past based

on different modeling strategies (Mann et al., 2006;

Ataullakhanov et al., 2002; Zarnitsina et al., 1996;

Butenas and Mann, 2002) or (Kuharsky and Fogelson,

2001; Anand et al., 2003). One of the major problems

is the high complexity of the chemical system which

crucially depends on the supply of chemicals by the

ﬂow. This leads to an important dependence of the

coagulation process on the blood ﬂow. On the other

hand the ﬂow is determined by the extent of the clot,

forming an obstacle to blood ﬂow.

This is why any future successful model of blood

coagulation necessarily has to be coupled with the

blood ﬂow. The blood itself exhibits a very com-

plicated rheological behavior including phenomena

such as viscoelasticity and shear-thinning. The aim

of this paper is to present a successful way of cou-

pling two of the most complex macroscopic contin-

uum based models of blood ﬂow and biochemistry

proposed in (Anand et al., 2003; Anand et al., 2008).

The work presented herein merges and signiﬁcantly

extends our previous studies (Bodn´ar and Sequeira,

2008) and (Bodn´ar et al., 2011a). Due to the lack

of space here, we refer the reader to these papers for

many details concerning the model development and

implementation as these are necessary to fully under-

652

Bodnár T. and Sequeira A..

On the Implementation of a Non-linear Viscoelastic Model into Coupled Blood Flow-biochemistry Model.

DOI: 10.5220/0004621306520657

In Proceedings of the 3rd International Conference on Simulation and Modeling Methodologies, Technologies and Applications (BIOMED-2013), pages

652-657

ISBN: 978-989-8565-69-3

 2013 SCITEPRESS (Science and Technology Publications, Lda.)

stand the work presented below.

2 MATHEMATICAL MODELS

The mathematical model consists of several parts de-

scribing the ﬂow of blood and a biochemical cascade

leading to clot formation.

2.1 Blood Flow Model

The ﬂow is described by a non-linear shear-thinning

viscoelastic model following the thermodynamic

framework established in (Rajagopal and Srinivasa,

2000) and extended for blood ﬂow in (Anand and

Rajagopal, 2004). The set of governing equations

is based on the conservation of mass (reduced to

divergence-free constraint) and conservation of linear

momentum for an incompressible ﬂuid.

divu

u = 0. (1)

= divT

T (2)

The stress tensor T

T is split as follows:

T = − p1

1+ µB

(t)

+ η

D (3)

where D

D denotes the symmetric part of the velocity

gradient tensor and 1

1 stands for identity tensor. The

upper-convected time-derivative of the elastic stretch

tensor B

(t)

is given by

▽

(t)

= −2K



tr(B

(t)

) − 3λ



(t)

− λ1

(4)

where the coefﬁcient λ depends on the trace of the

inverse of the tensor B

(t)

according to

λ =



−1

p(t)



. (5)

The remaining model coefﬁcients for blood are taken

exactly from (Anand and Rajagopal, 2004):

= 0.01 Pa· s µ = 0.1611 N/m

n = 0.5859 K = 58.0725 s

−1

More details and explanation concerning this model

and its implementation can be found e.g. in our pre-

vious paper (Bodn´ar et al., 2011a).

The subscript κ

(t) is used to emphasize that the stretch

is expressed withrespect to natural (time dependent) conﬁg-

uration κ

(t). This notation follows exactly the original pa-

pers (Rajagopal and Srinivasa, 2000; Anand and Rajagopal,

2004) and (Bodn´ar et al., 2011a) where the model has been

introduced and used.

2.2 Biochemistry Model

The biochemistry model is based on a coupled set of

advection-diffusion-reaction (ADR) equations. It has

been originally developed in (Anand et al., 2003) and

further extended in (Anand et al., 2008). It describes

the spatio-temporal evolution of concentrations C

23 chemical constituents (enzymes, zymogens, pro-

teins, etc.).

= div(K

gradC

) + R

(6)

The non-linear chemical reaction terms R

are mainly

based on second order or Michaelis-Menten kinetics.

As an example let’s mention the reaction term R

the equation for ﬁbrin

[IIa][I]

+ [I]

−

[PLA][Ia]

+ [Ia]

(7)

The concentrations of thrombin (denoted by [IIa]),

ﬁbrinogen (denoted by [I]), ﬁbrin (denoted by [Ia])

and plasminogen (denoted by [PLA]) are used to eval-

uate the reaction term R

. The chemical kinetics rates

, h

and constants K

, H

are known (taken from

(Anand et al., 2008)). The values of the diffusion pa-

rameters K

and the exact form of the reaction terms

are given in (Bodn´ar and Sequeira, 2008), where

the model has been for the ﬁrst time implemented and

used in 3D simulations.

2.3 Coupling Strategy

The coupling between blood ﬂow and biochemistry is

based on the ﬁbrin concentration. Fibrin is one of the

major constituents of clots (Fasano et al., 2012) and

thus it can be used as an indicator of the clot forma-

tion. The main idea is to make the material properties

of blood/clot dependent on ﬁbrin concentration. For

low ﬁbrin concentration the ﬂuid behaves like blood,

while for high ﬁbrin concentration it changes its be-

havior to a clot-like medium. In our model the ﬂuid

viscosity is multiplied by a factor

that locally de-

pends (linearly, up to a certain saturation value η

∗

) on

ﬁbrin concentration [Ia]. The viscosity η

is multi-

plied by a non-dimensional factor

= min



∗

− 1

clot

[Ia],η

∗



(8)

In our study we have used η

∗

= 100 and C

clot

1000 nM.

The clot is thus characterized as a highly viscous

ﬂuid. As a result, the region occupied by this simu-

lated clot represents an obstacle to the ﬂow of blood,

The subscript Ia refers to the chemical notation for ﬁb-

rin.

OntheImplementationofaNon-linearViscoelasticModelintoCoupledBloodFlow-biochemistryModel

653

with much lower viscosity. This effect is even sig-

niﬁcantly magniﬁed due to the shear-thinning non-

Newtonian behavior of blood, leading to a further in-

crease of ﬂuid viscosity in regions of low shear.

Changes in viscosity modify the local ﬂow ﬁeld

which consequently affects the concentration ﬁeld

that leads to further changes in viscosity. In this way

the two-way biochemistry-ﬂow coupling is enforced.

More details on this simple coupling technique

can be found in (Bodn´ar and Sequeira, 2008) where it

has been used together with a generalized Newtonian

model for blood, i.e. neglecting the viscoelasticity.

In the present study the technique has been extended

to the new non-linear viscoelastic model as suggested

in (Anand and Rajagopal, 2002; Anand et al., 2005;

Anand et al., 2008).

3 NUMERICAL METHODS

The system of governing equations is rather complex

and highly non-linear. This is why the numerical dis-

cretization has been chosen to be as simple and pre-

dictable as possible. We do not claim this choice is

optimal, it only serves as the ﬁrst step that allows us to

evaluate the underlying mathematical model and test

its applicability in simple conﬁgurations. The semi-

discretization approach is adopted to ﬁrst discretize

the PDEs in space and then integrating the resulting

system of ODEs in time. The same discretization is

employed for ﬂow variables, viscoelastic stress ten-

sor and concentrations in the biochemistry model.

The space discretization is based on a simple cen-

tral ﬁnite-volume discretization on a structured grid

with hexahedral cells. A multiblock grid topology

with wall-ﬁtted cells was used. The viscous ﬂuxes are

also discretized by the ﬁnite-volume technique over a

diamond-shaped cells adjoint to primary control vol-

umes faces. This approach was used in our previous

papers (Bodn´ar and Pˇr´ıhoda, 2006; Bodn´ar and Se-

queira, 2010; Bodn´ar et al., 2011b) or in (Keslerov´a

and Kozel, 2011; Keslerov´a, 2013).

The time integration is performed using a Runge-

Kutta (RK) multistage scheme. A speciﬁc advection-

diffusion optimized RK method has been used to re-

duce the computational cost. The basic idea behind

this subclass of RK methods is to split the space-

discretization operator into its inviscid and viscous

part. The inviscid part is evaluated in every stage of

RK method while the viscous ﬂuxes are only evalu-

ated in few stages. This corresponds to an operator

splitting technique with different RK methods (coef-

ﬁcients) used for the advection step and another for

the diffusion step. This approach allows to save sev-

eral (very expensive) evaluations of diffusive ﬂuxes

per time-step, while retaining the rather large stability

region of the RK method. For details see (Jameson

et al., 1981; Jameson, 1991) or (Bodn´ar and Sequeira,

2010; Bodn´ar et al., 2011b).

Along with these two basic components of the

numerical solver a speciﬁc stabilization technique is

used to avoid non-physical numerical oscillations due

to the central discretization. The non-linear TVD ﬁl-

ter (Engquist et al., 1989; Shyy et al., 1992) was

used to smooth the concentration ﬁelds, as reported

recently in (Bodn´ar, 2012).

4 NUMERICAL RESULTS

The numerical test case follows almost exactly the se-

tups used in (Bodn´ar and Sequeira, 2008) and (Bodn´ar

et al., 2011a) where we refer the complete parame-

ter set for this simulation. The geometry represents

a straight section of a blood vessel with diameter

6.2mm and length 31mm with the grid shown in the

ﬁgure 1. The clotting surface is simulated in a region

Figure 1: Grid structure.

that is formed by the intersection of a sphere with the

blood vessel wall. The evolution in space and time of

the clot is tracked down.

The model is very complex and thus the simu-

lations generate large amounts of data to be visual-

ized, analyzed, and understood. In this paper only a

few snapshots of results are presented to demonstrate

some of the most important types of model outputs.

The evolution in time of some of the coagulation

factors can be seen in the ﬁgures 3, 2 and 4. The con-

centration is visualized in a single point located in the

center of the clotting surface on the vessel wall. Only

the initial 600 seconds of clotting are shown. These

graphs show the nature of the coagulation process,

initially very fast, with rather slow long term evolu-

tion in the later phase.

The spatial extent of the clot can be shown using

SIMULTECH2013-3rdInternationalConferenceonSimulationandModelingMethodologies,Technologiesand

Applications

654

0 100 200 300 400 500 600

Concentration [nM]

Time [s]

Figure 2: Thrombin (Factor IIa) concentration.

500

1000

1500

2000

2500

3000

3500

4000

4500

5000

0 100 200 300 400 500 600

Concentration [nM]

Time [s]

Figure 3: Fibrin (Factor Ia) concentration.

0 100 200 300 400 500 600

Concentration [nM]

Time [s]

Figure 4: Tissue Plasminogen Activator (tPA) concentra-

tion.

the contours of ﬁbrin. Figure 5 shows the contours of

ﬁbrin concentration on the surface of the blood ves-

sel (in the range 50–5000nM in 15 levels with ex-

ponential distribution) and their variation in time. It

is evident that the core of the clot is rather stable in

time, while the downstream concentration ﬁeld varies

in time as a consequence of the interaction with the

blood ﬂow.

The ﬂow ﬁeld is modiﬁed in the clot region as de-

Figure 5: Fibrin concentration evolution

(Time=200s,400s,600s).

Figure 6: Axial velocity in near-wall layer.

picted in ﬁgure 6, showing the contours of the axial

velocity in the near-wall layer (the ﬁrst grid cell). The

red color means high velocity, while blue color is used

in low speed regions. The ﬂow is not only deceler-

ated in the clotting area, but it is also deﬂected in the

tangential direction as it is shown in ﬁgure 7. The

green color is used in unaffected regions, while the

red/blue is used for positive/negative values. This is

only shown to demonstrate the qualitatively correct

behavior of the ﬂow, that is stopped inside the clot-

ting area and forced to ﬂow around this region. The

exact contour values are not given as they depend on

the distance from the wall which was chosen arbitrar-

ily just for illustration.

OntheImplementationofaNon-linearViscoelasticModelintoCoupledBloodFlow-biochemistryModel

655

Figure 7: Tangential velocity in near-wall layer.

5 CONCLUSIONS

The numerical study presented in this paper has

demonstrated the successful implementation of both,

the blood ﬂow and the biochemistry model. The

model is now more complex in comparison with the

one presented in our previous work on blood coagu-

lation (Bodn´ar and Sequeira, 2008). The viscoelas-

tic extension of the model should allow to extend the

range of applicability of the model to critical ﬂow

regimes. The price to pay for this non-linear vis-

coelastic model extension is an important increase of

computational cost. The original, generalized New-

tonian model with shear-thinning viscosity, contained

4+ 23 = 27 PDEs to solve in 3D. The new model has

6 more equations for the components of the viscoelas-

tic stress tensor. This means that we have to solve now

4+ 6+ 23 = 33 equations for the coupled ﬂow + rhe-

ology + biochemistry model.

Future research will focus on performance and ro-

bustness improvements of the model and numerical

solvers. The stability issues raised in (Sequeira et al.,

2011) should be also addressed in the context of this

new model. Both of these topics will be important in

future applications of the model requiring long time

clot evolution simulations.

ACKNOWLEDGEMENTS

The ﬁnancial support for the present project was

partly provided by the Czech Science Foundation un-

der the Grant No.201/09/0917 and by the Portuguese

Science Foundation under the Project EXCL/MAT-

NAN/0114/2012.

REFERENCES

Ambrosi, D., Quarteroni, A., and Rozza, G., editors (2012).

Modeling of Physiological Flows, volume 5 of Mod-

eling, Simulation & Applications. Springer.

Anand, M. and Rajagopal, K. (2004). A shear-thinning

viscoelastic ﬂuid model for describing the ﬂow of

blood. Intern. Journal of Cardiovasc. Medicine and

Sci., 4(2):59–68.

Anand, M., Rajagopal, K., and Rajagopal, K. R. (2003). A

model incorporating some of the mechanical and bio-

chemical factors underlying clot formation and dis-

solution in ﬂowing blood. Journal of Theoretical

Medicine, 5(3–4):183–218.

Anand, M., Rajagopal, K., and Rajagopal, K. R. (2005).

A model for the formation and lysis of blood clots.

Pathophysiol. Haemost. Thromb., 34:109–120.

Anand, M., Rajagopal, K., and Rajagopal, K. R. (2008).

A model for the formation, growth and lysis of clots

in quiescent plasma. A comparison between the ef-

fects of antithrombin III deﬁciency and protein C de-

ﬁciency. J. Theor. Biol., 253:725–738.

Anand, M. and Rajagopal, K. R. (2002). A mathematical

model to describe the change in the constitutive char-

acter of blood due to platelet activation. Comptes Ren-

dues Mechaniques, 330:557–562.

Ataullakhanov, F., Zarnitsina, V., Pokhilko, A., Lobanov,

A., and Morozova, O. (2002). Spatio-temporal dy-

namics of blood coagulation and pattern formation. A

theoretical approach. International Journal of Bifur-

cation and Chaos, 12(9):1985–2002.

Bodn´ar, T. (2012). On the use of non-linear TVD ﬁlters

in ﬁnite-volume simulations. In Algoritmy 2012 Pro-

ceedings of Contributed Papers and Posters., pages

190–199, Bratislava. Slovak University of Technol-

ogy, Faculty of Civil Engineering.

Bodn´ar, T. and Pˇr´ıhoda, J. (2006). Numerical simulation of

turbulent free-surface ﬂow in curved channel. Journal

of Flow, Turbulence and Combustion, 76(4):429–442.

Bodn´ar, T., Rajagopal, K., and Sequeira, A. (2011a). Sim-

ulation of the three-dimensional ﬂow of blood using a

shear-thinning viscoelastic ﬂuid model. Mathematical

Modelling of Natural Phenomena, 6(5):1–24.

Bodn´ar, T. and Sequeira, A. (2008). Numerical simulation

of the coagulation dynamics of blood. Computational

and Mathematical Methods in Medicine, 9(2):83–104.

Bodn´ar, T. and Sequeira, A. (2010). Numerical study of the

signiﬁcance of the non-Newtonian nature of blood in

steady ﬂow through a stenosed vessel. In Rannacher,

R. and Sequeira, A., editors, Advances in Mathemati-

cal Fluid Mechanics, pages 83–104. Springer Verlag.

Bodn´ar, T., Sequeira, A., and Pirkl, L. (2009). Numeri-

cal simulations of blood ﬂow in a stenosed vessel un-

der different ﬂow rates using a generalized Oldroyd-B

model. In Numerical Analysis and Applied Mathemat-

ics, volume 2, Melville, New York. American Institute

of Physics.

Bodn´ar, T., Sequeira, A., and Prosi, M. (2011b). On the

shear-thinning and viscoelastic effects of blood ﬂow

under various ﬂow rates. Applied Mathematics and

Computation, 217(11):5055–5067.

Butenas, S. and Mann, K. G. (2002). Blood coagulation.

Biochemistry (Moscow), 67(1):3–12. Translated from

Biokhimiya, Vol. 67, No. 1, 2002, pp. 5–15.

Cangiani, A., Davidchack, R., Georgoulis, E., Levesley,

A. G. J., and Tretyakov, M., editors (2013). Numer-

SIMULTECH2013-3rdInternationalConferenceonSimulationandModelingMethodologies,Technologiesand

Applications

656

ical Mathematics and Advanced Applications 2011.

Springer. Proceedings of ENUMATH 2011, the 9th

European Conference on Numerical Mathematics and

Advanced Applications, Leicester, September 2011.

Engquist, B., L¨otstedt, P., and Sj¨ogreen, B. (1989). Nonlin-

ear ﬁlters for efﬁcient shock computation. Mathemat-

ics of Computation, 52(186):509–537.

Fasano, A., Santos, R., and Sequeira, A. (2012). Blood

coagulation: a puzzle for biologists, a maze for math-

ematicians. In (Ambrosi et al., 2012), chapter 3, pages

41–75.

Galdi, G., Rannacher, R., Robertson, A., and Turek, S., edi-

tors (2008). Hemodynamical Flows - Modeling, Anal-

ysis and Simulation, volume 37 of Oberwolfach Sem-

inars. Birk¨auser.

Gambaruto, A., Janela, J., Moura, A., and Sequeira, A.

(2011). Sensitivity of hemodynamics in a patient

speciﬁc cerebral aneurysm to vascular geometry and

blood rheology. Mathematical Biosciences and Engi-

neering, 8(2):409–423.

Jameson, A. (1991). Time dependent calculations using

multigrid, with applications to unsteady ﬂows past air-

foils and wings. In AIAA 10th Computational Fluid

Dynamics Conference, Honolulu. AIAA Paper 91-

1596.

Jameson, A., Schmidt, W., and Turkel, E. (1981). Numeri-

cal solutions of the Euler equations by ﬁnite volume

methods using Runge-Kutta time-stepping schemes.

In AIAA 14th Fluid and Plasma Dynamic Conference,

Palo Alto. AIAA paper 81-1259.

Janela, J., Moura, A., and Sequeira, A. (2010). Absorbing

boundary conditions for a 3D non-Newtonian ﬂuid-

structure interaction model for blood ﬂow in arter-

ies. International Journal of Engineering Science,

48(11):1332–1349.

Keslerov´a, R. (2013). Numerical study of effect of stress

tensor for viscous and viscoelastic ﬂuids ﬂow. In

(Cangiani et al., 2013), pages 529–538. Proceedings

of ENUMATH 2011, the 9th European Conference on

Numerical Mathematics and Advanced Applications,

Leicester, September 2011.

Keslerov´a, R. and Kozel, K. (2011). Numerical solution

of laminar incompressible generalized newtonian ﬂu-

ids ﬂow. Applied Mathematics and Computation,

217(11):5125–5133.

Kuharsky, A. L. and Fogelson, A. L. (2001). Surface-

mediated control of blood coagulation: the role of

binding site densities and platelet deposition. Biophys-

ical Journal, 80:1050–1074.

Mann, K., Brummel-Ziedins, K., Orfeo, T., and Butenas, S.

(2006). Models of blood coagulation. Blood Cells,

Molecules, and Diseases, 36:108–117.

Pirkl, L., Bodn´ar, T., and T˚uma, K. (2011). Viscoelastic

ﬂuid ﬂows at moderate weissenberg numbers using ol-

droyd type model. In AIP Conference Proceedings,

volume 1389, pages 102–105. American Institute of

Physics.

Rajagopal, K. and Srinivasa, A. (2000). A thermodynamic

frame work for rate type ﬂuid models. Journal of Non-

Newtonian Fluid Mechanics, 80:207–227.

Robertson, A., Sequeira, A., and Kameneva, M. (2007).

Hemodynamical Flows: Modelling, Analysis and Sim-

ulation, chapter Hemorheology. Birkh¨auser. in press.

Robertson, A., Sequeira, A., and Owens, R. (2008). Car-

diovascular Mathematics, chapter Rheological mod-

els for blood. Springer-Verlag. to appear.

Sequeira, A., Santos, R. F., and Bodn´ar, T. (2011). Blood

coagulation dynamics: Mathematical modeling and

stability results. Mathematical Biosciences and En-

gineering, 8(2):425–443.

Shyy, W., Chen, M.-H., Mittal, R., and Udaykumar, H.

(1992). On the suppression of numerical oscillations

using a non-linear ﬁlter. Journal of Computational

Physics, 102:49–62.

Zarnitsina, V. I., Pokhilko, A. V., and Ataullakhanov,

F. I. (1996). A mathematical model for the spatio–

temporal dynamics of intrinsic pathway of blood co-

agulation - I. The model description. Thrombosis Re-

search, 84(4):225–236.

OntheImplementationofaNon-linearViscoelasticModelintoCoupledBloodFlow-biochemistryModel

657