Considerations of Blood Properties, Outlet Boundary Conditions and Energy Loss Approaches in Computational Fluid Dynamics Modeling

Ji Young Moon; Dae Chul Suh; Yong Sang Lee; Young Woo Kim; Joon Sang Lee

doi:10.5469/neuroint.2014.9.1.1

Abstract

Despite recent development of computational fluid dynamics (CFD) research, analysis of computational fluid dynamics of cerebral vessels has several limitations. Although blood is a non-Newtonian fluid, velocity and pressure fields were computed under the assumptions of incompressible, laminar, steady-state flows and Newtonian fluid dynamics. The pulsatile nature of blood flow is not properly applied in inlet and outlet boundaries. Therefore, we present these technical limitations and discuss the possible solution by comparing the theoretical and computational studies.

Go to :

Modeling blood vessels is one of the main interests in computational simulation of hemodynamics. By materializing the geometry of certain blood vessel cases and running simulations, we can obtain a tremendous amount of data regarding certain afflictions such as aneurysm and stenosis. However, accurate modeling of the circulatory system is still a challenging objective. There are many technical difficulties to overcome in simulating blood flow, including the non-Newtonian properties of blood, the complicated shapes and elasticity of vessels and the pulsatile blood flow [1-4].

Researchers have dedicated effort to developing accurate methods for circulatory system modeling by simplifying or approximating properties that are difficult to precisely imitate [5-8]. One of these simplified methods includes an assumption which regarded blood as a Newtonian fluid, ignoring its characteristics and treating it as a water-like liquid. Also the usage of simple boundary conditions is inevitable, since the complicated vessel shape of the entire vascular system can neither be measured nor materialized completely.

Based on the recent articles published in Stroke [9-12], we identified three major weak points which are overlooked in the computational modeling and simulation processes of the circulatory system. First, the non-Newtonian property is scarcely applied to blood flow. As blood is not a pure liquid but rather a compound of tissue particles such as red blood cells, its viscosity differs greatly with change in shear rate (Fig. 1A). The second consideration is the inappropriate outlet boundary conditions of simulation models. The effect from outside of simulation domain, including backward pressure by wave reflection, must be considered at the outlet boundary for precise modeling (Fig. 1B). Third, when energy loss term is used in result analysis, the effect by vessel's geometry itself is not reflected properly. Blood flow properties are affected not only by abnormal vessel conditions such as aneurysm, but also by the cross-sectional area or directional change of blood vessels. These factors must be considered for accurate calculation of energy loss (Fig. 1C).

Fig. 1

The three main considerations in circulatory system modeling: A. non-Newtonian properties of blood, B. outlet boundary conditions of the model, and C. additional energy loss due to vessel geometry.

These technical issues have emerged to the surface partly due to recent progress of computational simulation methodology. In this paper, the three main modeling considerations are pointed out in detail. The importance of these three categories is explained thoroughly by numerical explanation and comparison with the works of other researchers.

Non-Newtonian Properties

Blood is a non-Newtonian fluid or more precisely, a pseudo-plastic fluid [13]. This type of fluid is less like water and more like toothpaste or ketchup, flowing more smoothly as the velocity of flow increases. This smooth flowing phenomenon at higher velocity is related with decreased viscosity and increased shear rate of the blood. Such phenomenon is called as the shear thinning behavior of the blood as a non-Newtonian fluid (Fig. 1A).

Blood flow is pulsatile; therefore, properties such as pressure, velocity, flow rate, and viscosity change as a function of time. Wall shear rate is also a function of time (Figs. 2A, B), indicating that if wall shear rate changes, viscosity might also change. This phenomenon cannot be explained with Newtonian assumptions. However, the Carreau-Yasuda model, one of the non-Newtonian models, demonstrates the change in viscosity according to time (cardiac cycle). This is an important reason why non-Newtonian analysis must be used in the pulsatile cardiac conditions.

Fig. 2

A. Wall shear rate predicted using the flow rate and radii common carotid, brachial, and femoral arteries reported by Stroev et al.[30] (Reprinted with permission) B. Viscosity calculated using shear rate from Carreau-Yasuda model (

= (1+(λη)^a)^(n-1)/a, µ₀=0.0022 Paty c a=0.644, n=0.392, λ=0.110s) in common carotid, brachial, and femoral arteries compared to Newtonian viscosity[31].

Several studies have emphasized the importance of this non-Newtonian behavior of blood by comparing Newtonian and non-Newtonian modeling research. Gijsen et al. [14] showed that applying non-Newtonian fluid flow at arteries with bifurcations greatly affects the flow features compared to Newtonain fluid flow. Molla and Paul [15] also compared Newtonian and non-Newtonian models in their research of arterial stenosis. Application of non-Newtonian properties was emphasized at specific diseases such as aneurysm. Wall shear stress, which is one of the most important factors for diagnosis of aneurysm, is the function of viscosity and shear rate. Assuming blood as Newtonian fluid means that the effect of viscosity term is ignored. Fisher and Rossmann [16] compared analysis of wall shear stress (WSS) distribution at aneurysm with both non-Newtonian (Casson, Carreau, Power law) and Newtonian models. The result showed that WSS of non-Newtonian models is more influential in aneurysms than that of Newtonian model. Also, Kim et al. [17] applied non-Newtonian properties to the analysis of stents in aneurysms, proposing that viscosity and velocity profiles in the aneurysm are largely determined by non-Newtonian properties of blood flow.

Go to :

Outlet Boundary Condition

Boundary conditions are the given conditions at the domain boundaries and which account for outer domain properties, are one of the most important factors for a modeling process. As the simulation domain is restricted to only a small proportion of the total circulatory system, the conditions outside of the domain cannot be directly applied to the simulation. Directly measured properties from a specimen may be applied at the boundaries; however, it is difficult to fully demonstrate the entire vascular system due to the problems of universality, simplicity, and technical deficiency. Simplified numerical boundary conditions are instead applied in most circulatory system modeling. There are several methods for selecting the outlet boundary condition. The' zero pressure outlet' or 'traction-free'boundary condition is used in various models due to its simplicity [18]. However, when the purpose of modeling is to analyze blood flow, the zero pressure outlet should not be used due to its inaccurate simulation results. The zero pressure boundary is as akin to assuming that the vessel is cut and exposed to atmospheric air conditions. Therefore, it neglects the change in pressure and flow rate by' wave reflection' from downstream vessels, resulting in considerable difference from actual blood flow.

Wave reflection is an influx of blood due to the reflected wave from downstream arteries. This is mainly due to the pulsatile flow of blood and elasticity of blood vessels. As pulsatile antegrade flow passes through the elastic vessel, this vessel repeatedly shrinks and expands, causing retrograde flow through the whole vascular system. This wave reflection results in considerable changes in blood pressure and flow rate due to overlap with the forward flow of blood. This wave reflection can be used in diagnosis of cardiovascular disease. As wave reflections and arterial stiffness are associated with an increased risk of cardiovascular disease, the measurement of central arterial pressure can provide a valid aid in cardiovascular risk assessment [19, 20].

To demonstrate this wave reflection, the term 'impedance'is applied to outlet boundary conditions. The well-known outlet conditions in theoretical hemodynamics include resistance, Windkessel, and impedance boundary models [18].

Resistance model substitutes the outer domain vessels into one single variable, while Windkessel model uses three variables [21]. These two models can develop wave reflection phenomenon, yet still inappropriate because it cannot show the arrival time of wave reflection correctly. However, impedance boundary model reflects the entire vessel geometry and property. By constructing the 1D vessel system that represent the real vessels, the variable called 'impedance' is calculated from the peripheral blood vessel to the boundary of simulation domain. Unlike the other models, this impedance variable includes both the data of magnitude and arrival time of wave reflection [18]. Because it considers wave reflection from bifurcations and geometry of vessels, the impedance model successfully demonstrates the pulsatile aspects of real blood flow. Comparison of zero pressure boundary condition and impedance boundary condition is shown in Fig. 3.

Fig. 3

Comparison between the zero pressure boundary condition and the impedance boundary condition. Zero pressure boundary condition assumes that pressure at the outlet boundary is zero, as if bleeding were occurring, without considering any influences from the vessels outside of the simulation domain. However, in impedance boundary case, as shown at right, the term Z is adopted to demonstrate the wave reflection caused by the outer domain. By regarding the blood flow of the whole vascular system as if it were electrical current running through a circuit, and so adopting the impedance term to express resistance of vessel to blood flow, the wave reflection is successfully accounted for at the outlet of the simulation domain.

The difference of those two boundary conditions compared in a vessel mentioned in Fig. 3, was clearly demonstrated by Vignon-Vlementel et al [18]., as shown in Fig. 4. When there is 75% area reduction stenosis in a branch of the bifurcation vessel model, mean flow of constant pressure boundary condition is divided into 70% and 30% each, because the fluid resistance by the geometry of computational domain is only considered. For the impedance boundary condition, on the other hand, the flow is also affected by entire circulatory system so that the difference between those two vessels is smaller than that of the constant pressure boundary condition. These results agree with representative clinical observations for iliac artery stenosis with less than an 85% area reduction [18]. For these reasons, the wave reflection should be considered to perform an accurate simulation for the purpose of analysis.

Fig. 4

The simulation of mean flow rate at the stenotic blood vessel with two different boundary conditions: the constant pressure and the impedance boundary condition [18]. (Reprinted with permission) A. Geometry of a bifurcated blood vessel with stenosis that reduces the cross sectional area by 75%. B. Flow distribution between normal and stenotic iliac vessels. Flow rate is calculated at the bottom end of each bifurcation. For the case of constant pressure boundary, the left branch of bifurcation with stenosis has a drastically different flow rate compared to the right branch, which is completely unrealistic. In the case of impedance boundary, however, the difference between the two branches is much smaller as the pressure from the outer domain has been applied properly.

Go to :

Energy Loss Due to the Vessel's Pathway

Interpreting computational results is as much important as performing accurate simulation. There are various ways to interpret simulation result data, such as plotting pressure and flow rate or calculating wall shear stress [16]. One of the possible methods is energy loss analysis. As the loss of energy is equal to the amount of energy applied to the vessel wall, it is directly related the risk of conditions such as the rupture of aneurysm [12, 22-24]. However, there are many factors, which are mostly related to the geometry of vessel, to be considered when using the energy loss equation. It might cause a considerable error on result if those factors are not reflected properly.

To predict clinical treatment effects, Bernoulli's energy loss concept is used [25]. If P_i, V_i are pressure and velocity at the initial cross sectional area domain respectively and P_f, V_f the pressure and velocity at the final cross sectional area domain, then the energy loss equation can be written as:

(1)

This Bernoulli's energy loss is purely due to mechanical energy loss. The energy loss consists of a limited number of ' loss factors', so some papers have already claimed that it cannot be applied in all clinical pathologies such as aneurysms or stenosis [25-27]. As energy loss can change by depending on its geometrical conditions, the following issues should be considered when using the loss term.

The main properties of energy loss by the pathway of blood are the change in cross sectional area and the direction of flow axis, based on Bernoulli's equation and momentum conservation [26, 28]. From Fig. 5A, if ΔE_geometry is the energy loss from a change in geometry, such as cross-sectional area or directional change, then the mechanical energy balance (Bernoulli's equation) is

(2)

Fig. 5

A. Schematic representation of a stream of fluid traversing a vessel from an initial position i to a final position f. Uniform flow experiences a sudden change in its cross-sectional area A and direction θ. The symbols A_i and A_f each represent the cross-sectional area, and V_i and V_f represent the corresponding average velocities of flow. B. The relationship between geometrical energy loss and cross-sectional area ratio. Each line represents the degree of change in flow angle. Each dot represents the point where the loss becomes minimum at the listed angles (0°, 30°, 45°, 60°, 90°).[26, 28] (Reprinted with permission)

Also, the force-momentum relationship (momentum conservation) is

(3)

(P_i - P_f) = ρQ (V_f - V_i cosθ)/A_f = ρQ² (A_i - A_f cosθ)/A_iA_f²

where ρ is density of blood, θ is angle change in direction, Q is flow rate (m³/s) and A is cross sectional area. Using Equations 2 and 3, we can derive the following equation:

(4)

From Equation 4, energy loss can occur due to change in cross-sectional area or flow angle. Fig. 5B shows this relationship between energy loss and change of cross sectional area or angle. The velocity difference at initial and final points is related to the momentum conservation equation, Δ(mv)=Ft, where m is mass, v is velocity, F is force and t is time. A change in 'mv' means that force is applied and the velocity and direction of blood flow are changed by the force of collisions with the vessel wall. The altered velocity and direction result in a change in fluid properties, which are directly related to energy loss.

From Equation 4,

(5)

Fig. 5B shows that energy loss is also affected by changes in cross sectional area and flow angle. The larger the flow angle changes, the larger the energy loss becomes. In the case of vessel contraction, the larger the change in cross-sectional area, the larger energy loss is. The dots presented at 0-60 degrees represent the cross-sectional area where energy loss is at a minimum. Energy loss tends to decrease until its minimum value, then increase again. As, energy loss approaches

the cross-sectional area change is larger at each change of angle. In the case of a change in 90°, energy loss decreases to

without a minimum value.

The level of change in cross-sectional area and direction of flow also affects energy loss. This means that if the rate of change in cross-sectional area or the angle of path is larger, then the value of energy loss also increases [26, 29].

A. Change in cross-sectional area but no alteration in direction of flow

The following Equations 6 and 7 show that variation in vessel geometry affects energy loss:

(6)

where P₀=P_i and P_N = P_f , A₀ = A_i and A_N = A_f

(7)

We can rewrite Equation 7, ΔE_geometry (θ=0) = function of β·B_i,f , B_i,f =

, which is term of Bernoulli's equation.

In Fig. 6A, it is assumed that the vessel having some gradient of cross sectional area is divided by' N' number of steps. Then the momentum conservation between steps (n-1) and (n) can be written as Equation 8.When θ=0, from Equations 2 and 6, Equation 7 can be derived to find the value of ΔE_geometry by using 'β =

' which represent cross-sectional area ratio.

is from the Bernoulli's equation and

is the loss factor based on ratio of change in cross sectional area. As observed in Fig. 6A,

gets smaller with increasing N, meaning that the greater the cross-sectional area changes, the less the energy loss is.

Fig. 6

A. Gradual change of cross-sectional area and magnitude the function of β for increasing values of N (# of changes) at

= 1/10. B. Change in angle where there is no change in cross-sectional area. The graph at the bottom shows that a larger number of transitions causes the angle change to become more gradual resulting in lower 'loss factor'. Both cases have directional change of 60°[26]. (Reprinted with permission)

B. No change in cross-sectional area but an alteration in direction of flow

When the angle changes by 'N' number of steps, the energy loss between n-1 and n steps is demonstrated by Equation 8 which is derived from Equation 4:

(8)

The sum of the ΔE_i,f (θ) values in all N steps can be expressed by Equation 9:

(9)

ΔE_i,f (θ) is

multiplied by 'energy loss factor'. This factor decreases as N increases (Fig. 6B). To sum up, the energy loss decreases as the angle change becomes more gradual

Equation 9 can be rewritten as:

(10)

It has been shown that, even without an aneurysm there is energy loss due to the geometry of the vessel itself. Hence, energy loss near an aneurysm consists of loss from the aneurysm itself and the loss from the geometry of the vessel, expressed as ΔE=ΔE_aneurysm+ΔE_geometry. These equations show that the geometry of vessel must be considered when using the energy loss in diagnosis.

Go to :

SUMMARY

Three important modeling considerations that are rarely used in circulatory system modeling are explained in this paper. The non-Newtonian properties of blood are found to be critical for the accuracy of modeling. Choosing the appropriate outlet boundary condition also makes an important difference in the results of a simulation. Finally, the energy loss related with blood vessel geometry must be included when the energy term is used in numerical calculations.

There are many researchers developing better modeling techniques and numerical methods, yet these new methods are not routinely applied in circulatory system modeling. Establishing accurate and trustworthy modeling of the circulatory system requires the newest techniques and implementation.

Go to :

Acknowledgment

This work (NRF-2013R1A2A2A01015333) was supported by the Mid-Career Researcher Program through and National Research Foundation of Korea (NRF) grant funded by the Ministry of Science, ICT and Future Planning (MSIP).

Go to :

References

1. Chen J, Lu XY. Numerical investigation of the non-Newtonian blood flow in a bifurcation model with a non-planar branch. J Biomech. 2004; 37:1899–1911. PMID: 15519598.

2. Olufsen MS, Peskin CS, Kim WY, Pedersen EM, Nadim A, Larsen J. Numerical simulation and experimental validation of blood flow in arteries with structured-tree outflow conditions. Ann Biomed Eng. 2000; 28:1281–1299. PMID: 11212947.

3. Darne B, Girerd X, Safar M, Cambien F, Guize L. Pulsatile versus steady component of blood pressure: a cross-sectional analysis and a prospective analysis on cardiovascular mortality. Hypertension. 1989; 13:392–400. PMID: 2522417.

4. Johnson DA, Spaeth JR, Rose WC, Naik UP, Beris AN. An impedance model for blood flow in the human arterial system. Part i: Model development and matlab implementation. Comput Chem Eng. 2011; 35:1304–1316.

5. Suh DC, Park ST, Oh TS, Park SO, Lim OK, Park S, et al. High shear stress at the surface of enhancing plaque in the systolic phase is related to the symptom presentation of severe M1 stenosis. Korean J Radiol. 2011; 12:515–518. PMID: 21852914.

6. Suh DC, Ko YB, Park ST, Yoon K, Lim OK, Oh JS, et al. Computational flow dynamics of the severe m1 stenosis before and after stenting. Neurointervention. 2011; 6:13–16. PMID: 22125742.

7. Park ST, Kim JK, Yoon KH, Park SO, Park SW, Kim JS, et al. Atherosclerotic carotid stenoses of apical versus body lesions in high-risk carotid stenting patients. AJNR Am J Neuroradiol. 2010; 31:1106–1112. PMID: 20093309.

8. Oh TS, Ko YB, Park ST, Yoon KH, Lee SW, Park JW, et al. Computational flow dynamics study in severe carotid bulb stenosis with ulceration. Neurointervention. 2010; 5:97–102.

9. Shojima M, Oshima M, Takagi K, Torii R, Hayakawa M, Katada K, et al. Magnitude and role of wall shear stress on cerebral aneurysm: computational fluid dynamic study of 20 middle cerebral artery aneurysms. Stroke. 2004; 35:2500–2505. PMID: 15514200.

10. Xiang J, Natarajan SK, Tremmel M, Ma D, Mocco J, Hopkins LN, et al. Hemodynamic-morphologic discriminants for intracranial aneurysm rupture. Stroke. 2011; 42:144–152. PMID: 21106956.

11. Luo B, Yang X, Wang S, Li H, Chen J, Yu H, et al. High shear stress and flow velocity in partially occluded aneurysms prone to recanalization. Stroke. 2011; 42:745–753. PMID: 21233477.

12. Takao H, Murayama Y, Otsuka S, Qian Y, Mohamed A, Masuda S, et al. Hemodynamic differences between unruptured and ruptured intracranial aneurysms during observation. Stroke. 2012; 43:1436–1439. PMID: 22363053.

13. Johnston BM, Johnston PR, Corney S, Kilpatrick D. Non-Newtonian blood flow in human right coronary arteries: steady state simulations. J Biomech. 2004; 37:709–720. PMID: 15047000.

14. Gijsen FJ, van de Vosse FN, Janssen JD. The influence of the non-Newtonian properties of blood on the flow in large arteries: steady flow in a carotid bifurcation model. J Biomech. 1999; 32:601–608. PMID: 10332624.

15. Molla MM, Paul MC. LES of non-Newtonian physiological blood flow in a model of arterial stenosis. Med Eng Phys. 2012; 34:1079–1087. PMID: 22153320.

16. Fisher C, Rossmann JS. Effect of non-Newtonian behavior on hemodynamics of cerebral aneurysms. J Biomech Eng. 2009; 131:091004. PMID: 19725693.

17. Kim YH, Xu X, Lee JS. The effect of stent porosity and strut shape on saccular aneurysm and its numerical analysis with lattice Boltzmann method. Ann Biomed Eng. 2010; 38:2274–2292. PMID: 20300847.

18. Vignon-Clementel IE, Figueroa CA, Jansen KE, Taylor CA. Outflow boundary conditions for 3D simulations of non-periodic blood flow and pressure fields in deformable arteries. Comput Methods Biomech Biomed Engin. 2010; 13:625–640. PMID: 20140798.

19. Salvi P, Bellasi A, Di Iorio B. Does it make sense to measure only the brachial blood pressure? Blood Purification. 2013; 36:21–22. PMID: 23735431.

20. Weber T, Auer J, O'Rourke MF, Kvas E, Lassnig E, Berent R, et al. Arterial stiffness, wave reflections, and the risk of coronary artery disease. Circulation. 2004; 109:184–189. PMID: 14662706.

21. Nichols W, O'Rourke MF. Mcdonald's blood flow in arteries: Theoretical, experimental, and clinical principles. London: Edward Arnold;1998.

22. Passerini T, Sangalli L, Vantini S, Piccinelli M, Bacigaluppi S, Antiga L, et al. An integrated statistical investigation of internal carotid arteries of patients affected by cerebral aneurysms. Cardiovasc Eng Technol. 3:26–40.

23. Qian Y, Takao H, Umezu M, Murayama Y. Risk analysis of unruptured aneurysms using computational fluid dynamics technology: preliminary results. AJNR Am J Neuroradiol. 2011; 32:1948–1955. PMID: 21903914.

24. Qian Y, Liu J, Itatani K, Miyaji K, Umezu M. Computational hemodynamic analysis in congenital heart disease: simulation of the norwood procedure. Ann Biomed Eng. 2010; 38:2302–2313. PMID: 20195758.

25. Akins CW, Travis B, Yoganathan AP. Energy loss for evaluating heart valve performance. J Thorac Cardiovasc Surg. 2008; 136:820–833. PMID: 18954618.

26. Ascuitto RJ, Kydon DW, Ross-Ascuitto NT. Streamlining fluid pathways lessens flow energy dissipation: relevance to atriocavopulmonary connections. Pediatr Cardiol. 2003; 24:249–258. PMID: 12522658.

27. Guillot M, Ross-Ascuitto N, Ascuitto RJ. Fluid-flow energetics for curved or angulated pathways associated with staged operations for the modified fontan procedure. Congenit Cardiol Today. 2012; 10:1–13.

28. Ascuitto RJ, Kydon DW, Ross-Ascuitto NT. Pressure loss from flow energy dissipation: relevance to Fontan-type modifications. Pediatr Cardiol. 2001; 22:110–115. PMID: 11178663.

29. Ensley AE, Lynch P, Chatzimavroudis GP, Lucas C, Sharma S, Yoganathan AP. Toward designing the optimal total cavopulmonary connection: an in vitro study. Ann Thorac Surg. 1999; 68:1384–1390. PMID: 10543511.

30. Stroev PV, Hoskins PR, Easson WJ. Distribution of wall shear rate throughout the arterial tree: a case study. Atherosclerosis. 2007; 191:276–280. PMID: 16828101.

31. Moore S, David T, Chase JG, Arnold J, Fink J. 3D models of blood flow in the cerebral vasculature. J Biomech. 2006; 39:1454–1463. PMID: 15953607.

Go to :