ABSTRACT

A short review of some reference solutions for the magnetohydrodynamic flow of blood is proposed in this paper. We present in details the solutions of Hartmann (1937), of Vardanyan (1973) and of Sud et al. (1974). In each case, a comparison is provided with the corresponding solution for the flow without any external magnetic field, namely Poiseuille (plane or cylindrical) and Womersley. We also present a synopsis of some other solutions for people who would like to go further in this topic. The interest in MHD flow of blood may be motivated by many reasons, such as Magnetic Resonance Imaging (MRI), Pulse Wave Velocity measurement, magnetic drug targeting, tissue engineering, mechanotransduction studies, and blood pulse energy harvesting… These fundamental solutions should also be used as particular limiting cases to validate any proposed more elaborated solutions or to validate computer codes.

Keywords:

Magnetohydrodynamic Flow of Blood, Magnetic Field in Biomedical Applications, Hartmann, Womersley

1. Introduction

The motivation to study the flow of blood in the presence of an external magnetic field has changed over the years: in the years 1960; the aim was, for example, to study the influence of a magnetic field on people who worked in factories, or to try to use the Lorentz force to slow down the speed of blood in case of haemorrhage. The development of the Magnetic Resonance Imaging (MRI) technology has induced a new interest for MHD flow of blood [1] [2] [3] . A survey of the recent literature [4] shows that such studies now have applications also in Pulse Wave Velocity measurement [5] [6] , magnetic drug targeting [7] [8] [9] [10] [11] , tissue engineering [12] [13] [14] [15] , mechanotransduction studies [16] , blood pulse energy harvesting [17] , and others… [18] .

The aim of this paper is to focus on some reference solutions for the magnetohydrodynamic flow of blood, although a short synopsis of some other solutions is proposed at the end of the paper for people who would like to go further in this topic. The fundamental solutions reviewed here should be used as particular limiting cases to validate any proposed more elaborated solutions. We present in details the calculations of Hartmann (1937) [19] [20] (unidirectional flow between two parallel plates under the action of a constant pressure gradient), of Vardanyan (1973) [21] (unidirectional flow in a rigid cylindrical vessel under the action of a constant pressure gradient) and of Sud et al. (1974) [22] (unidirectional flow in a rigid cylindrical vessel under the action of a pulsatile pressure gradient). In these three calculations, the walls of the duct are considered electrically insulating, and the induced electric and magnetic fields are neglected. This is the “low magnetic Reynolds number approximation”. In each case, a comparison with the corresponding solution for the flow without any external magnetic field, namely Poiseuille (plane or cylindrical) and Womersley [23] , is also provided.

2. Solution of Hartmann (1937)

The fluid flows in the direction (Oz), between two parallel non conducting plates located at x = −a and x = +a, under the action of a constant pressure gradient (∂P/∂z). The plates have an infinite width (in the direction (Oy)). An external transverse magnetic field B₀ is applied along the direction (Ox), and an external electric field E₀ along the direction (Oy). The components of the velocity u are u(0, 0, u).

Since blood is an incompressible fluid and the flow is unidirectional, the continuity equation (divu = 0) simply reduces to ∂u/∂z = 0. With the hypothesis of “infinite width along y”, we also have: ∂u/∂y = 0. The flow is stationary, so that ∂u/∂t = 0.

The low magnetic Reynolds number approximation enables the original MHD equations to be reduced to the Navier-Stokes equations including the Lorentz force j^B:

$\overrightarrow{O}=-\overrightarrow{\nabla }P+\overrightarrow{j}\wedge \overrightarrow{B}+{\mu }_f\Delta \overrightarrow{u}$ (1)

where j, the electric current density (A/m²), is obtained from Ohm’s law:

$\overrightarrow{j}=\sigma \left(\overrightarrow{E}+\stackrel{\rightharpoonup }{u}\wedge \stackrel{\rightharpoonup }{B}\right)$ (2)

σ is the electric conductivity of the fluid (S/m) and μ_f is its dynamic viscosity (Pa.s).

Due to the particular geometry of the problem studied, j has only one component not equal to zero, j_y, given by:

$j_y=\sigma \left(-E_0+u{B}_0\right)$ (3)

The “z” projection of Equation (1) is:

$\frac{\partial P}{\partial z}=-j_y{B}_0+{\mu }_f\frac{{\partial }^{2}u}{\partial x^2}$ (4)

Combining Equations (3) and (4), one obtains the differential equation to solve for u(x):

$\frac{\partial P}{\partial z}-\sigma E_0{B}_0={\mu }_f\frac{{\partial }^{2}u}{\partial x^2}-\sigma B_0^{2}u\left(x\right)$ (5)

Introducing the Hartmann number H_a,

$H_a=B_{0}a\sqrt{\frac{\sigma }{{\mu }_f}}$ (6)

we get:

$\frac{1}{{\mu }_f}\frac{\partial P}{\partial z}-\frac{\sigma E_0{B}_0}{{\mu }_f}=\frac{{\partial }^{2}u}{\partial x^2}-\frac{H_a^2}{a^2}u\left(x\right)$ (7)

With fixed, rigid plates, the no-slip condition at the walls yields: $u\left(x=-a\right)=u\left(x=a\right)=0$ , andthe solution of Equation (7) can be written as:

$u\left(x\right)=\frac{1}{H_a^2}\left(\frac{H_a^2{E}_0}{B_0}-\frac{a^2}{{\mu }_f}\frac{\partial P}{\partial z}\right)\left(1-\frac{ch\left(H_a\frac{x}{a}\right)}{ch{H}_a}\right)$ (8)

If u₀ is defined as u(x = 0), we have:

$\frac{u\left(x\right)}{u_0}=\frac{\left(ch\left(H_a\right)-ch\left(H_a\frac{x}{a}\right)\right)}{\left(ch\left(H_a\right)-1\right)}$ (9)

This solution is illustrated in Figure 1, where the non-dimensional velocity u(x)/u₀ is plotted versus x/a, for different values of the Hartmann number. The

Hartmann number measures the extent to which the magnetic forces prevail over viscous forces. The velocity profiles are more and more flattened as H_a increases, and they are stretched parallel to the direction of B₀. Since the curves are presented in non-dimensional form, the retardation effect of the Lorentz force cannot be seen (u₀ also changes with H_a, so that the ratio u(x)/u₀ is always equal to 1 at x = 0). Near the walls, there are thin boundary layers, called “Hartmann boundary layers”, where viscous drag drives the flow to zero and where shear frictions are increased.

If we consider the case of H_a going to zero, ch(H_a) may be approximated as:

$ch\left(H_a\right)=1+\frac{H_a^2}{2}+\frac{H_a^4}{24}+\cdots$ (10)

and u(x)/u₀ tends towards (1 − x²/a²) (Poiseuille solution between parallel plates).

The flow rate per unit width (dy = 1) (in the case E₀ = 0) is calculated as:

$Q_{Hart}={\displaystyle {\int }_{x=-a}^{x=+a}u\left(x\right)\text{d}x}$ (11)

where u(x) is given by Equation (8). This yields:

$Q_{Hart}=\frac{2a^3}{H_a^2{\mu }_f}\left(-\frac{\partial P}{\partial z}\right)\left\{1-\frac{sh\left(H_a\right)}{H_{a}ch\left(H_a\right)}\right\}$ (12)

Considering that:

$sh\left(H_a\right)=H_a+\frac{H_a^3}{6}+\frac{H_a^5}{120}+\cdots$ (13)

it is easy to show that:

$Q_{Hart}\to \frac{2a^3}{3{\mu }_f}\left(-\frac{\partial P}{\partial z}\right)=Q_{Pois.}^{plane},\text{when}{H}_a\to 0$ (14)

In the general case, we have:

$\frac{Q_{Hart.}}{Q_{Pois.}^{plane}}=\frac{3}{H_a^2}\left\{1-\frac{sh\left(H_a\right)}{H_{a}ch\left(H_a\right)}\right\}$ (15)

This quantity tends towards 1 when H_a tends towards zero.

If H_a = 1, Q_Hart./Q_Pois. = 0.715. This means that, under a same pressure gradient (-∂P/∂z), the flow rate in the case of Hartmann flow is reduced when compared to the classical Poiseuille-plane flow. This is due to the Lorentz force that slows the flow down. Other typical values are illustrated in Table 1.

3. Solution of Vardanyan (1973)

The situation studied is illustrated in Figure 2 [24] [25] . The cylindrical vessel is

assumed to be rigid. The flow is stationary (the pressure gradient, ∂P/∂z, is constant). The magnetic field B₀, is oriented along (Ox). The current density j is obtained from Ohm’s law (Equation (2), with E = 0):

In the cylindrical frame defined in Figure 2, we have:

$u\wedge B=\left|\begin{array}{c}0\\ 0\\ u\end{array}\right|\wedge \left|\begin{array}{c}{B}_0\cos \left(\theta \right)\\ -B_0\sin \left(\theta \right)\\ 0\end{array}\right|,\text{and}\text{consequently}j\wedge B=\left|\begin{array}{c}0\\ 0\\ -\sigma u{B}_0^2\end{array}\right|$ (16)

The longitudinal projection of the momentum equation (Navier Stokes) is reduced to:

${\mu }_f\left(\frac{{\partial }^{2}u}{\partial r^2}+\frac{1}{r}\frac{\partial u}{\partial r}\right)-\sigma B_0^{2}u=G,\text{with}G=\frac{\partial P}{\partial z}$ (17)

because we consider a stationary flow (∂u/∂t = 0) and because ∂u/∂z = 0 (this comes from the continuity equation).

At the wall, the boundary condition is: $u\left(r=a\right)=0$ .

Introducing the Hartmann number H_a (as defined in Equation (6)), Equation (17) may be written as:

$\frac{a^2}{H_a^2}\left(\frac{{\partial }^{2}u}{\partial r^2}+\frac{1}{r}\frac{\partial u}{\partial r}\right)-u\left(r\right)-\frac{G{a}^2}{{\mu }_f{H}_a^2}=0$ (18)

The solution of Equation (18) is:

$u\left(r\right)=-\frac{G{a}^2}{{\mu }_f{H}_a^2}\left[1-\frac{I_0\left(H_a\frac{r}{a}\right)}{I_0\left(H_a\right)}\right]$ (19)

where I₀ is the modified Bessel function of the first kind.

Using the power-series expansion of I₀:

$I_0\left(r\right)=1+{\displaystyle \sum _{k=1}^{\infty }\frac{{\left(\frac{r}{2}\right)}^{2k}}{{\left(k!\right)}^2}}$ (20)

it is easy to see that:

$\text{when}{H}_a\to 0,u\left(r\right)\to \frac{1}{4{\mu }_f}\left(-\frac{\partial P}{\partial z}\right)\left(a^2-r^2\right)$ (21)

which is the classical Poiseuille flow in a cylindrical tube.

If $u_0=u\left(r=0\right)$ , Equation (19) can be also expressed as:

$\frac{u\left(r\right)}{u_0}=\frac{\left(I_0\left(H_a\right)-I_0\left(H_a\frac{r}{a}\right)\right)}{\left(I_0\left(H_a\right)-1\right)}$ (22)

Equation (22) is illustrated in Figure 3.

As expected, the influence of B₀ is the same as for Hartmann’s solution (Figure 1). (The physical conditions are the same; only the geometry of the problem differs). As shown on Figure 2, the induced currents are captured inside the vessel, due to the fact that the wall is electrically insulating. This creates current loops inside the vessel and, on the return path, a compensating Lorentz force, that limits the flow retardation. Neglecting the induced fields (and, by the way, the current loops), the flow retardation is overestimated.

One may also be interested in calculating the flow rate from Equation (19):

$Q_{Vard.}={\displaystyle {\int }_0^{a}2\text{π}ru\left(r\right)\text{d}r}=-2\text{π}\frac{G{a}^2}{{\mu }_f{H}_a^2}{\displaystyle {\int }_0^a\left[r-\frac{r{I}_0\left(H_a\frac{r}{a}\right)}{I_0\left(H_a\right)}\right]\text{d}r}$ (23)

Using X = H_ar/a, and the following relation

${\displaystyle {\int }_0^{H_a}X{I}_0\left(X\right)\text{d}X}={\left[X{I}_1\left(X\right)\right]}_0^{H_a}=H_a{I}_1\left(H_a\right)$ (24)

the flow rate Q_Vard is found to be:

$Q_{Vard.}=-\text{π}\frac{G{a}^4}{{\mu }_f{H}_a^2}\left(1-\frac{2}{H_a}\frac{I_1\left(H_a\right)}{I_0\left(H_a\right)}\right)$ (25)

Using the expansions in ascending powers of H_a for I₀(H_a) and I₁(H_a):

$I_0\left(H_a\right)=1+\frac{H_a^2}{4}+\frac{H_a^4}{64}+\cdots \text{and}{I}_1\left(H_a\right)=\frac{H_a}{2}+\frac{H_a^3}{16}+\frac{1}{6}\frac{H_a^5}{2^6}+\cdots$ (26)

it may be shown that:

$Q_{Vard}\to \frac{\text{π}{a}^4}{8{\mu }_f}\left(-\frac{\partial P}{\partial z}\right)=Q_{Pois.},\text{when}{H}_a\to 0$ (27)

Consequently,

$\frac{Q_{Vard.}}{Q_{Pois.}}=\frac{8}{H_a^2}\left\{1-\frac{2I_1\left(H_a\right)}{H_a{I}_0\left(H_a\right)}\right\}$ (28)

For H_a = 1, under the same pressure gradient, the flow rate is reduced: Q_Vard./Q_Pois = 0.86; conversely, if one wants the same flow rate in the presence of the magnetic field than in its absence, it will be necessary to increase the pressure gradient: G_Vard/G_Pois = 1.16. Other typical values are illustrated in Table 2, and the comment is the same as for the data of Table 1.

4. Solution of Womersley (1955) (No Magnetic Field)

The situation considered is an unidirectional flow (in the z direction): u(r, t), in a rigid cylindrical tube with radius a; the pressure gradient is now pulsatile (harmonic):

$\frac{\partial P}{\partial z}=G{\text{e}}^{i\omega t}$ (29)

where ω is the angular frequency and i² = −1.

The continuity equation indicates that ∂u/∂z = 0. Consequently, the equation to solve for the longitudinal component of the velocity, u, is

$\rho \frac{\partial u}{\partial t}=-\frac{\partial P}{\partial z}+{\mu }_f\left(\frac{{\partial }^{2}u}{\partial r^2}+\frac{1}{r}\frac{\partial u}{\partial r}\right)$ (30)

associated with the boundary condition at the wall: $u\left(r=a\right)=0$ , any time.

The fluid density is denoted ρ.

The solution is searched in the form:

$u\left(r,t\right)=u^{*}\left(r\right){\text{e}}^{i\omega t}$ (31)

where u^*(r) is the solution of:

$0=-\frac{G}{\rho }-i\omega u^{*}\left(r\right)+\frac{{\mu }_f}{\rho }\left(\frac{{\partial }^2{u}^{*}}{\partial r^2}+\frac{1}{r}\frac{\partial u^{*}}{\partial r}\right)$ (32)

The general solution of Equation (32) is:

$u^{*}\left(r\right)=\frac{iG}{\rho \omega }\left(1-\frac{J_0\left(i^{3/2}\alpha \frac{r}{a}\right)}{J_0\left(i^{3/2}\alpha \right)}\right)$ (33)

where α is known as the “Womersley number”:

$\alpha =a\sqrt{\frac{\rho \omega }{{\mu }_f}}$ (34)

and J₀(X) is the Bessel function of the first kind and of zero order.

Noting that

${\displaystyle {\int }_0^{\alpha i^{3/2}}X{J}_0\left(X\right)\text{d}X}={\left[X{J}_1\left(X\right)\right]}_{X=0}^{X=\alpha i^{3/2}}$ (35)

the flow rate Q_Wom. may be obtained as follows:

$Q_{Wom.}={\text{e}}^{i\omega t}{\displaystyle {\int }_0^{a}2\text{π}r{u}^{*}\left(r\right)\text{d}r}=\frac{i\left(-G\right)\text{π}{a}^2}{\rho \omega }\left[F_{10}\left(\alpha \right)-1\right]{\text{e}}^{i\omega t}$ (36)

with:

$F_{10}\left(\alpha \right)=\frac{2J_1\left(i^{3/2}\alpha \right)}{i^{3/2}\alpha J_0\left(i^{3/2}\alpha \right)}$ (37)

Using the power expansions of J₀(X) and J₁(X),

$J_0\left(X\right)=1-\frac{X^2}{4}+\frac{X^4}{2^6}+\cdots \text{and}{J}_1\left(X\right)=\frac{X}{2}-\frac{X^3}{16}+\frac{1}{3}\frac{X^5}{2^7}+\cdots$ (38)

it is possible to come back to Poiseuille law when ω à 0.

Since $Q_{Pois}=\left(-G\right)\text{π}{a}^4/8{\mu }_f$ , the ratio Q_Wom/Q_Pois is given by:

$\frac{Q_{Wom}}{Q_{Pois}}=\frac{8i}{{\alpha }^2}\left[F_{10}\left(\alpha \right)-1\right]{\text{e}}^{i\omega t}$ (39)

We chose to illustrate Equations (31) and (33) in non-dimensional form, as follows:

$\frac{u\left(r,t\right)}{u\left(0,0\right)}=\frac{\left(J_0\left(i^{3/2}\alpha \right)-J_0\left(i^{3/2}\alpha \frac{r}{a}\right)\right)}{\left(J_0\left(i^{3/2}\alpha \right)-1\right)}{\text{e}}^{i\omega t}$ (40)

The results are presented in Figure 4, for different values of the Womersley number, α, and at different times (to show the pulsatility effect). The value α = 14.358 is obtained with physiological data: a = 0.01 m, ρ_f = 1050 kg/m³, f = 75 b.p.m., μ_f = 4 ´ 10⁻³ Pa.s.

When the Womersley number is large, the effect of the viscosity of the fluid does not propagate very far from the wall. In the central portion of the vessel, the transient flow is determined by the balance of the inertial forces and pressure forces as if the fluid was non viscous, and consequently, the profile is relatively blunt (in contrast to the parabolic profile of the Poiseuillean flow, which is determined by the balance of viscous and pressure forces) [26] [27] .

5. Solution of Sud et al. (1974)

The situation studied is the same as the one solved by Womersley [23] (pulsatile pressure gradient, rigid cylindrical tube), except that, as in the study of Vardanyan [21] , a transverse magnetic field B₀ is present. The longitudinal projection of momentum equation turns out to be:

$\rho \frac{\partial u}{\partial t}=-\frac{\partial P}{\partial z}+{\mu }_f\left(\frac{{\partial }^{2}u}{\partial r^2}+\frac{1}{r}\frac{\partial u}{\partial r}\right)-\sigma B_0^{2}u$ (41)

with the pressure gradient ∂P/∂z given by Equation (29).

At time t = 0, the solution for the velocity u is the stationary solution of Vardanyan [21] .

The solution of Equation (41) may be written as:

$u\left(r,t\right)=\frac{-G{a}^2}{{\mu }_f\left(H_a^2+i{\alpha }^2\right)}\left[1-\frac{J_0\left(\sqrt{-H_a^2-i{\alpha }^2}\frac{r}{a}\right)}{J_0\left(\sqrt{-H_a^2-i{\alpha }^2}\right)}\right]{\text{e}}^{i\omega t}$ (42)

where J₀ is the Bessel function of the 1 srt kind and of zero order, α is the Womersley number (same as in Equation (34)), and H_a is the Hartmann number (same as in Equation (6)).

It is easy to see that, when H_a = 0, the solution of Sud et al. reduces to the classical solution of Womersley. And, in the limiting case where ω à 0 (α à 0), it reduces to the solution of Vardanyan, because:

$u\left(r,t\right)=\frac{-G{a}^2}{{\mu }_f{H}_a^2}\left(1+i\omega t+\frac{{\left(i\omega t\right)}^2}{2}+\cdots \right)\left[1-\frac{J_0\left(i{H}_a\frac{r}{a}\right)}{J_0\left(i{H}_a\right)}\right]$ (43)

and consequently:

$u\left(r,t\right)\to \frac{-G{a}^2}{{\mu }_f{H}_a^2}\left[1-\frac{J_0\left(i{H}_a\frac{r}{a}\right)}{J_0\left(i{H}_a\right)}\right]$ (44)

Remembering the properties of Bessel functions $\left(J_0\left(iy\right)=I_0\left(y\right)\right)$ , we recognize the solution of Vardanyan.

Then we plot: $u\left(r,t\right)/u\left(0,0\right)$ versus r/a:

$\frac{u\left(r,t\right)}{u\left(0,0\right)}=\left[\frac{\left[J_0\left(\sqrt{-H_a^2-i{\alpha }^2}\right)-J_0\left(\sqrt{-H_a^2-i{\alpha }^2}\frac{r}{a}\right)\right]}{\left[J_0\left(\sqrt{-H_a^2-i{\alpha }^2}\right)-1\right]}\right]{\text{e}}^{i\omega t}$ (45)

The results are presented in Figure 5 for α = 14.358 and in Figure 6, for α = 5, in order to allow a comparison with the graphs of Figure 4. The curves present similar shapes, but, as shown in the calculations of Hartmann and of Vardanyan, the Lorentz force induces a flattening of the profiles. Magnetic forces also play a

role in the damping of the oscillations. This is more evident for the lower value of α (α = 5).

6. Discussion and Conclusion

An overview of some other analytical solutions is given in Table 3, according to the assumptions made in each of them for the fluid problem and for the electrical problem. All the solutions cited in Table 3 assume that the vessel wall is not conducting and that blood is newtonian.

In the case of a vessel wall with an electrical conductivity (σ_w ≠ 0), an analytical solution is unavailable and thus numerical solutions must be resorted to. Kinouchi et al. [30] were probably the first group to propose such a solution for the flow of blood. Due to the conductivity of the vessel wall, the induced currents are no more captured inside the vessel; instead, they are transmitted from the blood into the extravascular tissue. As explained in [24] , in that case, the retardation effect of the Lorentz force on the flow is higher than in the case of an insulating wall.

Finally, the results given in this paper may be discussed in terms of four non-dimensional numbers: the Reynolds number (R_e = ρu₀a/μ_f), the magnetic Reynolds number (R_em = au₀σμ_m, where μ_m is the fluid magnetic permeability), the Womersley number (defined in Equation (34)), and the Hartmann number (defined in Equation (6)). In the case of blood flow, typical values may be: a = 0.01 m, u₀ = 0.4 m/s, ρ = 1050 kg/m³, μ_f = 4 ´ 10⁻³ Pa.s, ω = 7.854 rd/s, B₀ = 1.5 T, σ = 0.5 S/m, and μ_m = 4π ´ 10⁻⁷ H/m. This yields: R_e = 1050, α = 14.36, H_a = 0.1677, R_em = 2.51 ´ 10⁻⁹. It is clear that the magnetic Reynolds number is quite small and that the “low magnetic number approximation” is totally justified. The calculations also demonstrate that, with a value of 0.17 for H_a, the influence of the magnetic field on the flow and pressure of blood is quite negligible. However, as shown in [24] , even if the induced magnetic field is negligible, the induced electric potentials may not be so small, since they are related to the magnetic field by a factor (1/μ_m).

Declarations

Competing interests: none.

Funding: none.

Ethical approval: not required.

Cite this paper

Drochon, A., Beuque, M. and Rodriguez, D.A.-A. (2018) A Review of Some Reference Analytic Solutions for the Magnetohydrodynamic Flow of Blood. Applied Mathematics, 9, 1179-1192. https://doi.org/10.4236/am.2018.910078

References