October 17, 2017, Tuesday
University of Colorado at Boulder Search A to Z Campus Map University of Colorado at Boulder CU 
Search Links


MBW:Womersley Arterial Flow

From MathBio

Jump to: navigation, search

Author: Tracey Morland

This article summarizes the results and approach introduced by Womersley in his famous paper on pulsatile flow in arteries [1]. His approach uses concepts from fluid mechanics, including Poiseuille Flow, to model the pressure gradient and flow velocity in an arterial pulse. The model considers the flow of blood in a rigid tube. Womersley's number is also defined, and the usefulness of Womersley flow is discussed.


Overview

This is a model of blood flow as a Newtonian fluid through a rigid tube propagated only by a pressure gradient. The Newtonian fluid is described by the Navier-Stokes second-order PDE's, and simplified by assuming a Poiseuille flow. A Fourier series accounts for the periodicity of the pressure gradient, and the coefficients are determined by the Fast Fourier Transform (FFT). Simplifications and substitutions allow for the system of PDE's to be approximated by Ordinary Differential Equations with solutions in the form of Bessel functions.

Biological Context

The human heart beats over 2.5 billion times in an average lifetime and about 100,000 times per day, and in one day your blood travels 12,000 miles, [2] or roughly the distance of traveling from Denver to Tokyo and back again. The heart is composed of four chambers. The sinoatrial nodes (SA) nodes in the right atrium (RA) initiate the electric pulse and cause the right ventricle (RV) to fill with blood. The action potential is propagated through the atria via the atrial cells [3]. The RV then contracts, sending blood to the pulmonary artery where it is then sent to the lungs. The blood, now fully oxygenated, now returns to the heart and fills the left atrium (LA). The LA contracts sending blood into the left ventricle (LV). Blood is then pumped from the LV into the aortic artery which sends the blood to the rest of the body. For an animation of the contraction of the heart valves, see NOVA: Map of the Human Heart.

Figure 1: The Human Heart, Nucleus Medical Art, Inc./Getty Images [2]

How fast does the blood flow through the arteries with each beat? Is there any phase-lag between the shock wave (the pulse you feel in your wrist or neck) and the flow of blood? These are questions addressed in Womersley's paper reviewed in this article. The more we can understand the nature of how the arterial pulse works and perform accurate calculations of arterial blood flow, the more we can detect heart disease and defects, such as in patients with diabetes and atherosclerosis [3].

History

Most of the information on the nature of the arterial pulse has evolved from the study of fluid mechanics. Various aspects of the pulse can be included in a model, including the elasticity of the artery, fluid viscosity, pressure gradient, and the presence of arterial bifurcations (branches) to name a few. In 1808, Thomas Young, a British scientist and physician, connected the elastic nature of the arteries to pulse wave velocity [4]. Then, in 1878, Moens and Korteweg independently derived a mathematical model relating arterial elasticity, or stiffness, to pulse wave velocity. Today it is called the Moens-Korteweg equation and it is dervied from Newton's second law of motion F=ma.

The simplest model for pulsatile flow is the [Windkessel model], developed by Otto Frank in 1899 [5]. This approach is not accurate enough to be used for quantitative analysis, but it provides a simple foundation on which to build more complicated models. In the model, blood storage is simplified to a single chamber, called a Windkessel, and the pressure of this chamber varies periodically over time [6]. Inflow into the chamber is from the heart, and outflow is to the outer arteries, veins and capillaries and is represented as a simple resistance vessel [3]. There are two parts to this model: compliance (representing the elastic nature of arteries) and resistance. The resistance in the system causes blood to enter the arteries at a higher rate than it flows out. Thus, there is storage of blood in the arteries [7].

Figure 2: Windkessel Model [8]

In 1970, a three-part Windkessel model was created by Westerhof [9]. Westerhof's model incorporates impedence which, combined with Frank's model, includes aspects of wave propagation. Westerhof's work was inspired by the previous work on pulsatile flow established by Womersley in his 1955 paper, as well as McDonald in his classic book [10] [11].

Mathematical Background

As mentioned above, Womersely's model makes use of Poiseuille flow, and is a simplification of the Navier-Stokes equations. Therefore, it is important to have some background on these approaches. The Navier-Stokes equations can be used to completely model the motion of incompressible, Newtonian fluids. However, these equations are very difficult to analyze since they are non-linear, second order partial differential equations, and only in a few special cases can their exact solutions be found [12]. The equations, simplified using the continuity equation, for the x,\;y,\;z directions are listed below [12]. In the right-hand side of Equation (1), and similarly for Equations (2) and (3), the term {\frac  {\partial \omega }{\partial x}} represents the pressure force, the term \rho g_{{x}} represents the weight of the fluid, and the second-order partials in parentheses represents the viscous forces.

Navier-Stokes Equations:

              {\begin{aligned}\rho ({\frac  {\partial \omega }{\partial t}}+\omega {\frac  {\partial \omega }{\partial x}}+v{\frac  {\partial \omega }{\partial y}}+w{\frac  {\partial \omega }{\partial z}})=\rho g_{{x}}-{\frac  {\partial P}{\partial x}}+\mu ({\frac  {\partial ^{{2}}\omega }{\partial x^{{2}}}}+{\frac  {\partial ^{{2}}\omega }{\partial y^{{2}}}}+{\frac  {\partial ^{{2}}\omega }{\partial z^{{2}}}})(1)\\\rho ({\frac  {\partial \omega }{\partial t}}+\omega {\frac  {\partial \omega }{\partial x}}+v{\frac  {\partial \omega }{\partial y}}+w{\frac  {\partial \omega }{\partial z}})=\rho g_{{y}}-{\frac  {\partial P}{\partial y}}+\mu ({\frac  {\partial ^{{2}}\omega }{\partial x^{{2}}}}+{\frac  {\partial ^{{2}}\omega }{\partial y^{{2}}}}+{\frac  {\partial ^{{2}}\omega }{\partial z^{{2}}}})(2)\\\rho ({\frac  {\partial \omega }{\partial t}}+\omega {\frac  {\partial \omega }{\partial x}}+v{\frac  {\partial \omega }{\partial y}}+w{\frac  {\partial \omega }{\partial z}})=\rho g_{{z}}-{\frac  {\partial P}{\partial z}}+\mu ({\frac  {\partial ^{{2}}\omega }{\partial x^{{2}}}}+{\frac  {\partial ^{{2}}\omega }{\partial y^{{2}}}}+{\frac  {\partial ^{{2}}\omega }{\partial z^{{2}}}})(3)\end{aligned}}

These equations are frequently written in cylindrical form. For example, Equation (1) can be written as:

   \rho ({\frac  {\partial \omega }{\partial t}}+\omega {\frac  {\partial \omega }{\partial x}}+v_{{r}}{\frac  {\partial \omega }{\partial r}}+{\frac  {v_{{\theta }}}{r}}{\frac  {\partial \omega }{\partial \theta }})=-{\frac  {\partial p}{\partial x}}+\rho g_{{x}}+\mu ({\frac  {\partial ^{{2}}\omega }{\partial x^{{2}}}}+{\frac  {\partial ^{{2}}\omega }{\partial r^{{2}}}}+{\frac  {1}{r}}{\frac  {\partial \omega }{\partial r}}+{\frac  {1}{r^{{2}}}}{\frac  {\partial ^{{2}}\omega }{\partial \theta ^{{2}}}})   (4)

A simplification of the Navier-Stokes equations can be made assuming a Poiseuille flow in which the velocity of the fluid is described by the following equation.

                                                               \omega ={\frac  {1}{4\mu }}{\frac  {dP}{dx}}[r^{{2}}-R^{{2}}]    (5)

A Poiseuille flow assumes that the flow is steady, uniform (over a cross-section), laminar, and axially symmetric within a cylindrical tube. Under these assumptions, the term {\frac  {\partial \omega }{\partial t}}=0 (since there is no change in velocity over time). Furthermore, the terms {\frac  {\partial \omega }{\partial x}} and {\frac  {\partial \omega }{\partial \theta }} equal zero.

Figure 3: Poiseuille Flow in a Rigid Tube [13]

(For another project using Poiseuille flow, see MBW:Optimum Design of Blood Vessel Bifurcation)

Notation

Before examining Womersley's model, it is important to define some notation that will be used.


\rho liquid density
\mu viscosity
v={\frac  {\mu }{\rho }} kinematic viscosity
p_{{1}},\;p_{{2}} pressures at ends of pipe
R radius of pipe
r distance from center of pipe
l length of pipe
\omega longitudinal velocity of liquid
f={\frac  {n}{2\pi }} fundamental frequency (typically the heart rate in radians)
T_{{0}} period of the pressure gradient wave

Womersley Flow Defined

We now take an in-depth look at Womersly's paper and model for pulsatile flow. To provide the basis for Womersley's model, we begin with a more complete derivation of Poiseuille's formula for steady flow and also include the pressure gradient. A constant pressure gradient throughout a pipe of length l is defined as:

                                                                          {\frac  {p_{{1}}-p_{{2}}}{l}}   (6)

After accounting for the simplifications of the Navier-Stokes equations based on Poiseuille flow, the equation of motion is then:

                                                      {\frac  {d^{{2}}\omega }{dr^{{2}}}}+{\frac  {1}{r}}{\frac  {d\omega }{dr}}+{\frac  {p_{{1}}-p_{{2}}}{\mu l}}=0   (7)

This has Equation (5) as its solution.

Now, Womersley expresses the pressure gradient as a periodic function of time with frequency f=n/2\pi to represent the arterial pulse. The representation of this function is done using Fourier series.

Using Fourier Series to Represent Pressure Gradient

Since the change in pressure gradient is periodic, it can be expressed using the following function:

                                                              {\frac  {\partial P}{\partial x}}=Re[\sum _{{n=0}}^{{\infty }}a_{{n}}e^{{ifnt}}]   (8)

where a_{{n}}=A_{{n}}+iB_{{n}}.

Using Euler's formula, this can be expressed as a Fourier series with Fourier coefficients A_{{0}},\;A_{{1}},...,\;A_{{n}},\;B_{{0}},\;B_{{1}},...,\;B_{{n}}.

                                       {\frac  {\partial P}{\partial x}}=A_{{0}}+\sum _{{n=1}}^{{N}}A_{{n}}cos(fnt)+\sum _{{n=1}}^{{N}}B_{{n}}sin(fnt)   (9)

The coefficients are calculated as follows:

                                              A_{{n}}={\frac  {2}{T_{{0}}}}\int _{{0}}^{{T_{{0}}}}{\frac  {dP}{dx}}cos(fnt)dt   (10)

                                              B_{{n}}={\frac  {2}{T_{{0}}}}\int _{{0}}^{{T_{{0}}}}{\frac  {dP}{dx}}sin(fnt)dt   (11)

These coefficients can be easily calculated using Matlab or Mathematica. A Fast Fourier Transform (FFT) algorithm can also be used instead to calculate these coefficients. See [14] page 196 for more information and sample Matlab code.

Solving for Flow Velocity (u)

Since we are assuming a Poiseuille flow that changes over time, we have that {\frac  {\partial \omega }{\partial x}}=0 and {\frac  {\partial \omega }{\partial \theta }}=0. After some algebra we get that for a single harmonic, n:

                                              {\frac  {\partial ^{{2}}w}{\partial r^{{2}}}}+{\frac  {1}{r}}{\frac  {\partial w}{\partial r}}-{\frac  {1}{v}}{\frac  {\partial w}{\partial t}}=-{\frac  {a_{{n}}e^{{ifnt}}}{\mu }}   (12)

Let's consider the following simple solution to Equation (12) since the velocity also changes periodically over time:

                                                                      \omega =u(r)e^{{ifnt}}   (13)

Next we substitute this equation for \omega into Equation (12) and divide both sides by e^{{ifnt}} to get the following ordinary differential equation that is not dependent on time:

                                                      {\frac  {d^{{2}}u}{dr^{{2}}}}+{\frac  {1}{r}}{\frac  {du}{dr}}-{\frac  {ifn}{v}}u={\frac  {-a_{{n}}}{\mu }}   (14)

This equation can be rewritten to be a Bessel zero-order differential equation using the fact that -i=i^{{3}}:


                                                      {\frac  {d^{{2}}u}{dr^{{2}}}}+{\frac  {1}{r}}{\frac  {du}{dr}}+{\frac  {i^{{3}}fn}{v}}u={\frac  {-a_{{n}}}{\mu }}   (15)

Solving the Bessel Differential Equation

The general form for a zero-order Bessel differential equation is:


                                                                {\frac  {d^{{2}}u}{dr^{{2}}}}+{\frac  {1}{r}}{\frac  {du}{dr}}+\lambda ^{{2}}u=0   (16)

In our case note that \lambda ^{{2}}={\frac  {i^{{3}}fn}{v}}. So the homogeneous version of Equation (15) fits the above form and its solution is of the general form:

                                             u_{{n}}(r)=c_{{1}}J_{{0}}(ri^{{{\frac  {3}{2}}}}{\sqrt  {{\frac  {fn}{v}}}})+c_{{2}}Y_{{0}}(ri^{{{\frac  {3}{2}}}}{\sqrt  {{\frac  {fn}{v}}}})   (17)

The term The term Y_{{0}} must be discarded since u(r) has the requirement that it must be finite at the center (origin) of the pipe. Therefore, c_{{2}}=0. Since Equation (15) is non-homogeneous, we use the technique described in [14] and try the simple solution u(r)=c_{{3}}. Then we have that the terms involving the derivative of u are zero, and thus we have that {\frac  {i^{{3}}fn}{v}}u={\frac  {-a_{{n}}}{\mu }}. Therefore, using the fact that v={\frac  {\mu }{\rho }}, we have that:

                                                                          c_{{3}}={\frac  {a_{{n}}}{i\rho fn}}   (18)

Then the solution looks like:

                                                     u_{{n}}(r)=c_{{1}}J_{{0}}(ri^{{{\frac  {3}{2}}}}{\sqrt  {{\frac  {fn}{v}}}})+{\frac  {a_{{n}}}{i\rho fn}}   (19)

We can now solve for c_{{1}} by using the no-slip condition that u=0 at the boundary r=R. Therefore we have:

                                                  0=c_{{1}}J_{{0}}(ri^{{{\frac  {3}{2}}}}{\sqrt  {{\frac  {fn}{v}}}})+{\frac  {a_{{n}}}{i\rho fn}} (20)

                                                    c_{{1}}={\frac  {a_{{n}}}{i\rho fnJ_{{0}}(Ri^{{{\frac  {3}{2}}}}{\sqrt  {{\frac  {fn}{v}}}})}} (21)

Then finally we get that:

                                             u_{{n}}(r)={\frac  {a_{{n}}}{i\rho fn}}[{\frac  {J_{{0}}(\lambda r)}{J_{{0}}(\lambda R)}}-1] (22)


The last step is to add the steady flow velocity term u_{{0}} and thus we have:

                                                              u(r)=u_{{0}}+\sum _{{n=1}}^{{\infty }}u_{{n}}(r)   (23)

Womersley Number

\alpha =r{\sqrt  {{\frac  {\omega }{v}}}}

The quantity r{\sqrt  {{\frac  {f}{v}}}} in Equation (19) is called the Womersley Number, \alpha . It is a dimensionless parameter that represents the ratio of transient forces, originating from the pulse wave, to the viscous force, or shear force. To get a feel for the magnitude of \alpha , an example problem is presented.


Example 1: A 20 kg dog has an average heart rate of 90 beats per minute (bpm) and a 70 kg human has an average of 70 bpm [15].The density of blood is 1060 kg/m^{{3}}. Using the following calculation ([15] Equation 3.7, p. 27) we calculate the diameter of the aorta to be 2.03 cm for humans, and 1.33 cm for dogs.

                                                                    D=.48W^{{.34}}

where W is the weight of the animal. From [16] we get that the blood viscosity for a human is 0.006 Ns/m^{{2}} and for a dog is 0.0056 Ns/m^{{2}}.

Now we calculate Womersley's number: r{\sqrt  {{\frac  {f}{v}}}}=r{\sqrt  {{\frac  {f\rho }{\mu }}}}. Note that we need to convert bpm to rad/s so we multiply the bpm by {\frac  {2\pi }{60}}={\frac  {\pi }{30}}.


For a human, we have: \alpha ={\frac  {2.03}{100}}m\cdot {\sqrt  {{\frac  {70({\frac  {\pi }{30}})rad/s\cdot 1060kg/m^{{3}}}{0.006Ns/m^{{2}}}}}}=23.1


For a dog, we have:\alpha ={\frac  {1.33}{100}}m\cdot {\sqrt  {{\frac  {90({\frac  {\pi }{30}})rad/s\cdot 1060kg/m^{{3}}}{0.0056Ns/m^{{2}}}}}}=17.76


This tells us that the oscillatory inertial forces become more important than the viscous force as the size of the animal increases, or particulary as the size of the blood vessels increases. (Note that this example is adapted from a similar example found in [14], p. 30).

Flow Rate (Q)

In the next part of the paper, Womersley derives the flow rate (Q) of the fluid passing through a cross-sectional area of the pipe, or artery in this case. This is accomplished by the integrating the velocity over a differential area:

                                                                      Q=2\pi \int _{{0}}^{{R}}u(r)rdr   (21) 

For steady flow, recall that u={\frac  {p_{{1}}-p_{{2}}}{\mu l}}(r^{{2}}-R^{{2}}). After integrating, we get that:

                                                                     Q={\frac  {p_{{1}}-p_{{2}}}{8\mu l}}\pi R^{{4}}   (22)

which is Poisueille's formula. Now, we can substitute in the formula we found for u(r), and then use the fact that \int xJ_{{0}}(x)dx=xJ_{{1}}(x), then after some calculus and algebraic simplification we get that:

                                Q_{{n}}={\frac  {\pi R^{{2}}}{\mu }}Re[\,{\frac  {a_{{n}}}{ifn}}\{1-{\frac  {2\alpha i^{{3/2}}}{i^{{3}}\alpha ^{{2}}}}{\frac  {J_{{1}}(\alpha i^{{3/2}})}{J_{{0}}(\alpha i^{{3/2}})}}\}e^{{ifnt}}]   (23)

Recall that the Q_{{n}} need to be summed for each n=1,2,3..., and added to the average flow rate given by the constant term in the Fourier series Q_{{0}}=a_{{0}}{\frac  {8\mu }{\pi R^{{4}}}} [14] to get:

                                                                    Q=Q_{{0}}+\sum _{{n=1}}^{{\infty }}Q_{{n}}   (24)

Womersley now derives a way to calculate Q without Bessel functions, perhaps because back in 1955 these functions were numerically more challenging to calculate. He uses modulus and phase functions [1] to accomplish this. The details of these functions are out of the scope of this review, but the following relations are used:

                                                                  J_{{0}}(xi^{{3/2}})=M_{{0}}(x)e^{{i\theta _{{0}}}}   (25)

                                                                  J_{{1}}(\alpha i^{{3/2}})=M_{{1}}e^{{i\theta _{{1}}}}   (26)

He also uses the fact that the real part of a_{{n}}e^{{ifnt}} is M_{{n}}\;cos(fnt+\phi ) where M_{{n}}={\sqrt  {A^{{2}}+B^{{2}}}}=|a_{{n}}| since a_{{n}}=A_{{n}}+iB_{{n}} [14]. Then, after some simplifications Womersley gets the following, simplified formula for Q (not including the steady flow term Q_{{0}}):

                                            Q={\frac  {\pi MR^{{4}}}{\mu }}{\frac  {M'_{{10}}}{\alpha ^{{2}}}}\;sin(fmt+\phi +\epsilon _{{10}})   (27)

where m=1,2,...,6 and the values for \alpha ^{{2}}, M, M'_{{10}}, \phi , and\epsilon _{{10}} are all given in a table in Womersley's paper. Let's look at an example problem presented in the paper to get more of an idea of the values for Q.


Example 2: The radius of an artery is 0.15 cm, the viscosity is \mu =.04P, and the pulse frequency is 180 bpm. The value for \alpha ^{{2}}, M, M'_{{10}}, \phi , and \epsilon _{{10}} are given in the following two tables:

Table1 morland.png
Table2 morland.png

The expression for Q is then then sum of the six terms:

{\begin{aligned}Q=3.56\;sin(t+71^{{\circ }}13')+2.71\;sin(2t-49{}^{{\circ }}13')\\-1.2\;sin(3t+50{}^{{\circ }}33')-.28\;sin(4t-20{}^{{\circ }}17')\\+.22\;sin(5t+99{}^{{\circ }}28')-.13\;sin(6t+5{}^{{\circ }}47')\end{aligned}}


This result for Q is plotted against the pressure gradient below:

Q.png

Analysis and Conclusions

Womersley's arterial flow model gives formulas to calculate the flow rate of a viscous fluid through a rigid tube under a periodic pressure gradient, described by a Fourier series. The solution is then extended to calculate the flow rate over a cross-sectional area of the tube.

The graph in Example 2 depicts the pressure gradient (P) juxtaposed with the flow rate (Q), and clearly shows the phase-lag between the two curves. This implies that the pulse wave is first sent through the body (seen as an increase in the pressure gradient) and then the blood flow follows. Typically the phase-lag is about 90 degrees, except at the boundary layer where it is about 45 degrees [17]. Also noteworthy is that this phase-lag is only present in large arteries; in smaller arteries less than 4mm in diameter, the flow behavior closely follows Poiseuille's formula [17]. Interestingly, the pulse wave travels about 5 times the maximum blood velocity [17].

Also evident from the graph is that the direction of flow velocity is actually reversed as seen when the velocity becomes negative. However, this is the topic of another paper written by Womersley et. al in 1955 [18].

Womersley's number has provided both fluid mechanics and biological sciences with a means to measure the inertial forces versus the viscous forces. It is as significant in analyzing unsteady flow as the Reynolds number is in measuring steady flow [19].

Womersley's model has been foundational to many models of arterial blood flow, but itself is limited to modeling laminar flow through a rigid, cylindrical tube. In reality, there are elastic effects in the blood vessel and the pressure gradient may depend on other factors besides time. In addition, Womersley's model considers only longitudinal velocity, whereas there may also be a radial part.

Recent Extension

In 2011, the Womersley article was cited in a study on Effects of vessel wall elasticity and non-Newtonian rheology on blood flow regime and hemodynamic parameters distribution by Foad Kabinejadian and Dhanjoo N. Ghista. This paper is a follow up of an earlier study, wherein they computationally simulated blood flow using Womersley's model. They were focused on vessel intersections and junctions. The current paper revisits the computational simulation with compliant walls and non-Newtonian fluid. They compared the results of their two studies and found that the compliant and non-Newtonian model was consistent with the results of earlier studies and observational data. They conclude that the inclusion of wall compliance and non-Newtonian rheology in flow simulation of blood vessels can be essential in quantitative and comparative analyses [20]. For more information on human blood flow rheology, read Viscoelastic Versus Newtonian Behavior from the Physical Review Letters.

References

  1. Womersley, J.R. “Method for the calculation of velocity, rate of flow and viscous drag in arteries when the pressure gradient is known,” J. Physiol. 127 (1955) 553-563. http://www.ncbi.nlm.nih.gov/pmc/articles/PMC1365740/.
  2. 2.0 2.1 NOVA. “Amazing Heart Facts.” http://www.pbs.org/wgbh/nova/heart/heartfacts.html.
  3. 3.0 3.1 3.2 Keener, J. and J. Sneyd. Mathematical Physiology II: Systems Physiology. Second Edition. Springer, New York, 2009.
  4. Young, Thomas. “The Croonian Lecture: On the Functions of the Heart and Arteries.” Phil. Trans. R. Soc. Lond. 1809 99, 1-31. http://rstl.royalsocietypublishing.org/content/99/1.2.full.pdf+html.]
  5. Frank, Otto. “Die Grundform Des Arteriellen Pulses Erste Abhandlung: Mathematische Analyse.” Zeitschrift f ¨ ur Biologie , 37 (1899). 483–526.
  6. Ethier, C. Ross , Craig A. Simmons. “Introductory Biomechanics: From Cells to Organisms.” Cambridge University Press. 2007.
  7. Wikipedia. “Windkessel Effect.” http://en.wikipedia.org/wiki/Windkessel_effect.
  8. Saouti, N., et al. “The arterial load in pulmonary hypertension.” Eur Respir Rev 2010; 19: 117, 197–203. http://err.ersjournals.com/content/19/117/197.full.pdf+html.
  9. Westerhof, Nico., Noodergraaf, A. “Errors in measurement of hydraulic input impedence.” J Biomech. 1970, 3: 351-356.
  10. McDonald, DA. “Regional Pulse Wave Velocity in the Arterial Tree.” Journal of Applied Physiology. 1968, 24:73-78.
  11. Gross, D.R., “Animal Models in Cardiovascular Research.” 3rd Edition. Springer Science + Business Media, LLC. New York: 2009.
  12. 12.0 12.1 Young, Donald F., Munson, Bruce R., and Theodore H. Okiishi. “A Brief Introduction to Fluid Mechanics.” Second Edition. John Wiley & Sons, Inc. New York, 2001.
  13. La Barck Jr., Anthony Joseph. “A Computational Study of Blood Flow and Oxygen Transport during Reperfusion and Post-conditioning.” Thesis. Rowan University, 2010.
  14. 14.0 14.1 14.2 14.3 14.4 Waite, Lee, and Jerry Fine. “Applied Biofluid Mechanics.” McGraw-Hill Professional. Chicago. 2007.
  15. 15.0 15.1 Li, John K-J. “Comparative Cardiovascular Dynamics of Mammals.” CRC Press. Boca Raton, 1996.
  16. Windberger, U. et. al. “Whole blood viscosity, plasma viscosity and erythrocyte aggregation in nine mammalian species: reference values and comparison of data.” The Physiological Society. Experimental Physiology (2003) 88.3, 431–440.
  17. 17.0 17.1 17.2 Lighthill, James (Sir). “Mathematical Biolfluiddynamics.” Society for Industrial Mathematics, Philidelphia: 1975.
  18. Womersley, J.R., Hale, J.F., and McDonald, D.A. “Velocity profiles of oscillating arterial flow, with some calculations of viscous drag and the reynolds number.” J. Physiol, 128:629-640. 1955.
  19. Chandran, K. B., et al. “Biofluid Mechanics: the Human Circulation.” CRC: Taylor and Francis, Boca Raton: 2007.
  20. Kabinejadian, Foad and Ghista, Dhanjoo. “Compliant model of a coupled sequential coronary arterial bypass graft: Effects of vessel wall elasticity and non-Newtonian rheology on blood flow regime and hemodynamic parameters distribution.” Elsevier Ltd, 2011. http://www.sciencedirect.com/science/article/pii/S1350453311002633.