# Modelling and Analysis of Complex Viscous Fluid in Thin Elastic Tubes

Authors: Yufang Gao; Zongguo Zhang
Journal: Complexity
 (2020)

---

## Abstract

Cardiovascular disease is a major threat to human health. The study on the pathogenesis and prevention of cardiovascular disease has received special attention. In this paper, we have contributed to the derivation of a mathematical model for the nonlinear waves in an artery. From the Navier–Stokes equations and continuity equation, the vorticity equation satisfied by the blood flow is established. And based on the multiscale analysis and perturbation method, a new model of the Boussinesq equation with viscous term is derived to describe the propagation of a viscous fluid through a thin tube. In order to be more consistent with the flow of the fluid, the time-fractional Boussinesq equation with viscous term is deduced by employing the semi-inverse method and the fractional variational principle. Moreover, the approximate analytical solution of the fractional equation is obtained, and the effect of viscosity on the amplitude and width of the wave is studied. Finally, the effects of the fractional order parameters and vessel radius on blood flow volume are discussed and analyzed.


---

## Body

## 1. Introduction

In recent years, in the field of biological rheology, the rheology which is related to the blood, blood vessels, and heart that constitute the human blood circulation has been developed rapidly. Blood circulation is a complex system that can cause great damage to the whole body when a certain organ is diseased. In order to better understand the physiological and pathological behavior of the human cardiovascular system, it is necessary to deeply study the dynamics of blood flow in the arteries and the mechanical factors of blood flow.Soliton phenomena exist in many fields [1–3]. Organism is a completely nonlinear complex medium, and both the blood composition and the structure of blood vessels show obvious nonlinear characteristics. The nonlinearity of blood flow has long been discovered by Womersley [4, 5] and McDonald [6, 7], which also provides a new direction and way for people to understand the law of life movement. Since then, many scholars have begun to study the field and made great progress. Ravindran et al. [8] derived the nonlinear Schödinger equation (the terms of pseudodifferential operators) governing the modulation of periodic waves. The KdV equation describing blood flow is obtained by Sigeo [9], to explain the steepness of pressure waves during propagation. Hashizume [10] analyzed the propagation of pressure waves from a theoretical perspective. Liu [11] combined arterial flow as a balanced flow with periodic small pulsatile flow and analyzed the effect of vascular elasticity on blood flow. Demiray [12, 13] considered propagation of wave through a viscous incompressible fluid contained in a prestressed thin elastic tube. Choy [14] deduced the mathematical model of nonlinear wave modulation of artery with stenosis.Fractional derivative theory and methods [15–18] are widely used in the study of nonequilibrium systems of various intermediate processes and critical phenomena in physics and mechanics, especially in nonlinear science [19–22]. Fractional differential equations are transformed in a standard differential equation by replacing the time derivative or the space derivative with the fractional derivative. Compared with integer-order differential equations, the most important advantage of fractional derivative equations is that it can better fit some natural physical processes and dynamic system processes. Moreover, the study of solving partial differential equations also has a new exploration, such as the extended tanh method [23], the exp-function method [24], the variational-iteration method [15, 25], the Hirota bilinear method [19], and the G′/G-expansion method [26].The structure of the full article is as follows: in Section2, the Boussinesq equation with viscous term is derived by the multiscale analysis and perturbation method and used for the first time to describe blood flow. Based on the new model, we obtain the time-fractional Boussinesq equation with viscous term in Section 3. In Section 4, the approximate analytical solution of the above fraction equation and the viscous effect is discussed. Finally, the effects of fractional order, vascular radius, and blood flow velocity on stroke volume are analyzed and studied.


## 2. Derivation of the Boussinesq Equation

As we all know, the cardiovascular system is a complete closed conduit system. Thus, the blood can be considered as an incompressible non-Newtonian fluid. When dynamic equilibrium of the blood is disturbed by a pressure pulse generated by the motion of the heart, a harmonic wave type of motion will be developed in the blood. Although previous researchers have done some research about this question, the viscosity of fluid is often ignored. We consider blood vessel as cylindrical shape and adopt cylindrical coordinate system to describe the motion of blood. In this paper, we setx as the central axis of a blood vessel and axial coordinate and set r as the radial coordinates. Therefore, based on Womersley theory [4], basic equations of pulsatile flow in arteries are described as follows:(1)ρ∂u∂t+u∂u∂x+v∂u∂r−η∂2u∂t2+1r∂u∂r+∂2u∂x2=−∂P∂x,(2)ρ∂v∂t+u∂v∂x+v∂v∂r−η∂2v∂t2+1r∂v∂r+∂2v∂x2−vr2=−∂P∂r,(3)∂u∂x+1r∂rv∂r=0,where u and v are the axial and radial velocities of blood, respectively; p is the atmospheric pressure; ρ is the blood density; and η represents the viscosity coefficient of blood. From (3), we can define the flow function of blood which satisfies(4)u=1r∂Ψ∂r,v=−1r∂Ψ∂x.Substituting (4) into (1) and (2), we can gain the vorticity equation of blood flow with viscosity term:(5)ρ∂∂t+1r∂Ψ∂r∂∂x−1r∂Ψ∂r+2r2∂Ψ∂xD2Ψ−η∂2∂r2−1r∂∂r+∂2∂x2D2Ψ=0,where D2≡∂2/∂r2−1/r∂/∂r+∂2/∂x2. Through analysis, we assume that the basic law of blood flow can be expressed by the Poiseuille flow. Therefore, the speed of blood can be expressed as(6)u¯r=umax1−rR2,where umax=ΔpR2/4μΔL is the maximum velocity of central axis of the blood vessel. ΔL represents the length of the blood fluid, and μ represents viscosity coefficient of blood. Δp is the pressure difference between two ends of blood, and R is the radius of the blood vessel. Next, when considering the viscosity of the vascular wall, we can get the revised form of the Poiseuille flow:(7)u¯r=umax2r2R2lnrR+1−rR2.We set basic flow function of blood asψ¯. Therefore, it meets(8)u¯=1r∂ψ¯∂r,where we define disturbed flow function φr,x,t which represents disturbance of blood flow. Therefore, the flow function can be expressed as(9)Ψr,x,t=ψ¯r+φr,x,t.And then, from (8) and (4), we can obtain(10)u=u¯+1r∂φ∂r,v=−1r∂φ∂x.Substituting (9) and (10) into (5) and considering that there are other effects of dissipation which can offset 2aη/ρ, we obtain(11)∂∂t+1r∂φr∂∂x−1r∂φ∂x∂∂r+2r2∂φ∂x+u¯∂∂xD2φ−a∂φ∂x=ηρD2D2φ,where a=4umax/R2 and limr⟶0u¯r=umax. The above equation is a vorticity equation of disturbed flow function φ which contains information about the effects on the heart when physiological or psychological condition changes dramatically.In the previous studies, researchers have obtained the KdV equation from the Navier–Stokes equation and continuity equation of blood flow. However, with the development of nonlinear theory increasing, more equations which have stronger nonlinear properties are needed to describe problems in reality. Therefore, the Boussinesq equation model will be deduced in this paper. In order to obtain this equation, we firstly define coordinate transformations with long wave approximation and adopt the space-time transformation as follows:(12)X=εx−ct,T=ε2t,r=r,and we set(13)ηρ=ε3Q,where ε is a small parameter. Substituting (12) into (11) yields(14)∂∂Tε2+1r∂φ∂r∂∂Xε−1r∂φ∂X∂∂rε+2r2∂φ∂Xε+u¯−c∂∂Xε∂2∂r2−1r∂∂r+∂2∂X2ε2φ−a∂φ∂Xε=Qε3∂2∂r2−1r∂∂r+∂2∂X2ε2α∂2∂r2−1r∂∂r+∂2∂X2ε2φ.Then, we extend perturbed flow function aboutε as follows:(15)φ=ε2φ0+ε3φ1+ε4φ2+⋯.Due to the fact that the velocity of wave propagation in blood is much faster than that of blood flow on the axis of the blood vessel, we can obtainu¯≪c. Thus, substituting (15) into (14), we can yield all levels of approximation equations about ε:(16)ε3:Γφ0=0,ε4:Γφ1=−∂∂T∂2∂r2−1r∂∂rφ0,ε5:Γφ2=−∂∂T∂2∂r2−1r∂∂rφ1+c∂3∂X3φ0−1r∂φ0∂r∂∂X−1r∂φ0∂X∂∂r+2r2∂φ0∂X∂2∂r2−1r∂∂rφ0+Q∂2∂r2−1r∂∂r∂2∂r2−1r∂∂rφ0,where Γ≡−c∂/∂X∂2/∂r2−1/r∂/∂r−a∂/∂X.Forε3, it is easy to see that φ0 has a solution, and the form is φ0=AX,Tϕ0r; then, we have(17)−c∂2∂r2−1r∂∂rϕ0−aϕ0=0.Forε4, by analysis, we assume ∂φ1/∂X=∂A/∂Tϕ1r, where ϕ1r satisfies the following equation:(18)−c∂2∂r2−1r∂∂rϕ1−aϕ1=acϕ0.Let us consider the equation of orderε5. By substituting φ0=AX,Tϕ0r, ∂φ1/∂X=∂A/∂Tϕ1r, (17), and (18) into the equation and multiplying the both sides of the order ε5 equation by ϕ0/r∂/∂X, we have(19)∂2∂X2ϕ0r∂2∂r2−1r∂∂rφ2+a∂2φ2∂X2ϕ0r1c=−∂2A∂T2ϕ0rac2ϕ1+ϕ0c−∂4A∂X4ϕ02r−ac2ϕ02r3∂2A2∂X2−Q∂A∂X1cϕ0r∂2∂r2−1r∂∂r2ϕ0.And then integrating it with respect tor from 0 to R, we can have the following equation:(20)∂2A∂T2+a1∂2A2∂X2+a2∂4A∂X4+a3∂A∂X=0,where(21)a0=∫0Rϕ0rac2ϕ1+ϕ0cdr,a1=∫0Ra/c2ϕ02/r3dra0,a2=∫0Rϕ02/rdra0,a3=Q∫0R1/cϕ0/r∂2/∂r2−1/r∂/∂r2ϕ0dra0.Equation (20) is a new Boussinesq equation with viscosity effect which can be used to describe the state of the disturbance flow of blood with viscous properties. Comparing with the previous KdV model, the above model has stronger nonlinear properties. Besides, it can be concluded that there is disturbance flow in the form of solitary waves in the arteries.


## 3. Derivation of the Time-Fractional Boussinesq Equation

In this section, we seek for the time-fractional Boussinesq equation by using the semi-inverse method and the fractional variational principle. Firstly, we use a potential functionVXX,T, where AX,T=VXX,T gives the potential equation of (20) in the following form:(22)VXTT+a1A2XX+a2VXXXXX+a3AX=0.The functional of the potential equation (22) can be represented by(23)JVX=∫X˜dX∫T˜dTVXc0VXTT+c1a1A2XX+c2a2VXXXXX+c3a3AX,where c0, c1, c2, and c3 are the Lagrangian multipliers which can be obtained later by taking the variation of (23). AX and A2XX are considered as the fixed functions.Making use of integrating (23) by parts and assuming VTT˜=VXX˜=VXXX˜=VXXXX˜=0, we can gain(24)JVX=∫X˜dX∫T˜dT−c0VXT2+c1a1VXA2XX+c2a2VXXX2+c3a3VXAX.Applying the variation optimum condition, we can derive the Euler equation as follows:(25)2c0VXTT+c1a1A2XX+2c2a2VXXXXX+c3a3AX=0.Comparing (25) with (22), we can obtain the following Lagrangian multipliers:(26)c0=12,c1=1,c2=12,c3=1.Therefore, the Lagrangian form of (20) can be given by(27)LVXT,VX,VXX,VXXX=−12VXT2+a1VXA2XX+12a2VXXX2+a3VXAX.At this time, the Lagrangian form of the time-fractional Boussinesq equation can be represented as(28)FD0TβVX,VX,VXX,VXXX=−12D0TβVX2+a1VXA2XX+12a2VXXX2+a3VXAX,where based on the left Riemann–Liouville compression derivative, the fractional derivative D0TβVX is defined as(29)D0TβfT=1Γ1−βddT∫aTdτT−τ−βfτ,0≤β≤1,T∈a,b.Therefore, the time-fractional form of (20) can be expressed by(30)JVX=∫X˜dX∫T˜dTFDTβ0VX,VX,VXX,VXXX.Thus, the variation of functional (30) can be obtained:(31)δJVX=∫X˜dX∫T˜dT∂F∂D0TβVXδD0TβVX+∂F∂VXδVX+∂F∂VXXδVXX+∂F∂VXXXδVXXX.As we all know, the fractional integration by parts is obtained by the following rule:(32)∫badTfTDaTβgT=∫abdTgTDTbβfT,fT,gT∈a,b,based on the right Riemann–Liouville fractional derivative, D0TβfT is obtained by(33)DTTβfT=−1Γ1−βddT∫Tbdττ−T−βfτ,0≤β≤1,T∈a,b.Using (31) and (32), we get(34)δJVX=∫X˜dX∫T˜dTDTbβ∂F∂D0TβVX+∂F∂VX−∂∂X∂F∂VXX+∂2∂X2∂F∂VXXXδVX.By using variation principle, whenδJV=0, we obtain the optimization of the variation of functional. And then, we yield the Euler–Lagrange equation of (20) as the following equation:(35)−D0Tβ∂F∂D0TβVX+∂F∂VX−∂∂X∂F∂VXX+∂2∂X2∂F∂VXXX=0.Substituting (28) into (35) and letting D0TβVX=AX,T, (35) is turned into the following form:(36)D0T2βAX,T+a1A2XX+a2AXXXX+a3AX=0.Equation (36) is a new model, which is gained by the time-fractional method. Therefore, it is named as time-fractional Boussinesq equation.


## 4. Solution of the Time-Fractional Boussinesq Equation

In the previous section, a new time-fractional Boussinesq equation had been gained. It can be used to describe the flow characteristics of blood with viscous properties. Due to this reason, we decide to yield the solution of the time-fractional Boussinesq equation to further characterize the flow state of blood.Let we define the fractional complex transform as follows:(37)τ=TαΓ1+α,where α is the unknown constant. And then, we can yield the fractional derivatives into classical derivatives with (37):(38)∂αA∂Tα=∂2A∂τ2.It is obvious that (36) can be deduced to(39)∂2A∂τ2+a1A2XX+a2AXXXX+a3AX=0.Next, we consider to gain the solution of (39). Due to the fact that a3 describes viscous properties of blood, we set μ=a3 as a small parameter. Assuming that A0=A0X,μ1/2τ and considering to get the properties of (39), two time scales are defined as follows:(40)ζ=τ,η˜=μ1/2τ,and A is expanded as(41)AX,τ=A1X,ζ,η˜+μA2X,ζ,η˜.Thus, we obtain all levels of approximation equations aboutμ:(42)μ0:∂2A1∂ζ2+a1∂2A12∂X2+a2∂4A1∂X4=0,μ1:∂2A2∂ζ2+a1∂2A1A2∂X2+a2∂4A2∂X4=−∂2A1∂η˜2−∂A1∂X.Forμ0, it is obvious that the solution of the equation can be expressed as(43)A1X,ζ,η˜=A0sech2a1A06a21/2X−νζ+c˜,where ν=−2a1A0/31/2 and A0 represents the maximum amplitude at initial moment. Then, for μ1, we use the equation to determine the form of A0. Assuming that(44)A2=BI,I=X+ζ.By substituting (44) into the order one of μ equation from (42), the following equation is obtained:(45)∂2B∂I2+a1∂2∂I2A1B+a2∂4B∂I4=MA1,where MA1=−∂2A1/∂η˜2−∂A1/∂X. In order to describe the structure of A2, we consider to choose proper function G which should be orthogonal to MA1 as(46)∫−∞+∞GIMA1dI=0.Multiplying both sides of (45) by G, integrating them by parts and using BI=±∞=0 and dB/dII=±∞=0, we can see that GI satisfies(47)d2GdI2+a1A1d2GdI2+a2d4GdI4=0.By assumingG±∞=0, the solution of (47) can be obtained as(48)G=A0sech2a1A06a21/2I+c˜.Substituting (48) into (46) and using (45), we obtain(49)A0=e−μτ−X.Therefore, we obtain the approximate solution of (39) as follows:(50)AX,ζ,η˜=A0sech2a1A06a21/2X−νζ+c˜,where A0X,μ1/2τ=e−μτ−X. The above equation is the solution of (39) and can be used to describe the flow state of blood more specifically.Finally, the solution of the time-fractional Boussinesq equation is(51)AX,T=exp−μTαΓ1+α−Xsech2a16a2exp−μTαΓ1+α−X1/2X−νTαΓ1+α+c˜,where ν=−2a1/3exp−μTα/Γ1+α−X1/2.In order to study the influence of blood viscosity and fractional order on the evolution of blood flow, we change the value ofμ, α, X, and T to get related results. From Figure 1, we can see that when viscosity parameter μ=0, the amplitude of solitary waves is largest; at the same time, the width of solitary waves is the smallest. With the decrease of the amplitude, the width increases gradually. Figure 2 depicts the evolution of solitary waves in blood under different viscosity coefficient μ. It is evident that amplitude of the solitary waves becomes smaller with time T from Figure 2(a), while the width of the solitary waves becomes larger with time T. And we can know that the viscosity effect leads to the amplitude of the solitary waves to decrease and the width of the solitary waves to increase. In addition, in contrast with Figure 2(b), we find with the increasing viscosity coefficient μ, the rate of amplitude reduction has been accelerated and the direction of the wave changes from the positive direction of X to the negative direction of it. Figure 3 shows the solitary waves under different fractional order α. We find the peak value of amplitude is not linear with α, and there is a minimum of the peak with the corresponding critical value of α. According to the critical value, the peak value of amplitude increases firstly and then decreases.Figure 1
Variation of amplitude and width according to the viscosity coefficientμα=1,T=5,X=0.005,a2=−1,a1=−20,c˜=0.Figure 2
The evolution of solitary waves in blood with viscosity: (a)α=1,a2=−1,a1=−20,c˜=0,μ=0.005; (b) α=1,a2=−1,a1=−20,c˜=0,μ=0.45.
(a)(b)Figure 3
The evolution of solitary waves with different numbers of fractional orderαT=1,a2=−1,a1=−20,c˜=0,μ=0.005.


## 5. Analysis and Discussion

In this section, we consider to calculate the value which can represent the state of blood by using the results from the previous sections. First, we search for the solution of (17) and (18). Taking ϕ0r=ϕ0c/az=zwz and substituting them into (17) yield(52)z2d2wdz2+zdwdz+z2−1=0.And the solution of (52) is expressed as(53)wz=C1J1z+C2N1z,where C1andC2 are the constants and J1z and N1z are expressed as(54)J1z=∑m=1∞−1mm!Γm+2z22m+1,N1z=J1zcosπ−J−1zsinπ.Thus, the solution of (17) can be described as(55)ϕ0r=C1racJ1rac.Similarly, by substituting (55) into (18), we obtain(56)ϕ1r=C3racJ1rac,where C3 is a constant. We take the physiological parameters of the human body as follows:(57)R=1.25×10−2m,umax=0.4m/s,L=6.3×10−2m,c=5.5m/s.And setQ=−0.01, C1=0.1, and C3=−500. Therefore, by using (21), we obtain the coefficients of (36) as follows:(58)a0=−8.52776×10−1,a1=−2.5917,a2=−5.9084×10−7,a3=3.7×10−3.Due to that the basic flow in a human artery can be assumed to be the Poiseuille flow, the stream function can be expressed as(59)ψ¯=umaxr22−r44R2.So the stream function can be obtained as follows:(60)Ψ=umaxr22−r44R2+ε2exp−μTαΓ1+α−Xsech2a16a2exp−μTαΓ1+α−X1/2X−νTαΓ1+α+c˜C1racJ1rac,where X=εx−ct and T=ε2t following (12).Because in the cylindrical coordinate system, the difference of the flow function values of any two points in the plane passing through the symmetry axis multiplied by2π is equal to the value, which is the flow rate of the rotating surface formed by the rotation of any curve passing through these two points around the symmetric axis, we decide to calculate blood flow by using the values of the outlet and entrance of blood vessels. We take the flow function of the central position of the vessel at the entrance of the aorta as Ψ1, where X=0 and r=0, and take the flow function at the junction of heart and aorta as Ψ1, where x=0 and r=R. Thus, we obtain(61)Ψ1=0,Ψ2=14umaxR2+ε2exp−με2tαΓ1+α+εctC1racJ1racsech2a16a2exp−με2tαΓ1+α+εct1/2−εct−−2a13exp−με2tαΓ1+α+εct1/2ε2tαΓ1+α+c˜.Thus, the blood flow at the outlet of the heart can be expressed as(62)Qheart=12πumaxR2+2πε2exp−με2tαΓ1+α+εctC1RacJ1Racsech2a16a2exp−με2tαΓ1+α+εct1/2−εct+2a13exp−με2tαΓ1+α+εct1/2ε2tαΓ1+α+c˜.In a cardiac cycle, the time-dependent image of the blood flow volume of the cardiac outlet is shown in Figure4.Figure 4
The blood flow volume of cardiac outlet at the different values of time andαc˜=0.From Figure4, we can see that, when α is smaller, this model can be used to describe people with less blood flow or faster blood flow conduction. And when α is large, this model can be used to describe people with slow blood flow in the body, or blood output from the heart is not easy to circulate throughout the body. These people may have lesions in some parts of the blood vessels, causing obstruction to blood flow. Or they are in good physical condition, but they are at high risk of cardiovascular and cerebrovascular diseases. Furthermore, we can find the appropriate value of α to describe different populations by combining the actual monitoring data.As we all know, stroke volume of healthy people is basically stable. When stroke volume is not stable for a long time, we can infer that there may be something wrong with that person’s cardiovascular system. If the amount of stroke volume is smaller, it will lead to insufficient oxygen supply to other tissues and organs. On the contrary, when stroke volume is too large, it will increase the load of human blood vessels, thus causing some cardiovascular and cerebrovascular diseases. According to the above introduction, we can find that stroke volume in each cycle of the heart is an important index which can reflect the health status of human cardiovascular system. Thus, we give the description of the stroke volume as follows:(63)V=∫0T˜Qheartdt,where T˜ represents the pulsatile cycle of the heart. By calculating, with umax=0.4 and R=1.25×10−2, we obtain the stroke volume in each cycle of the heart as 78.54ml which is close to the theoretical calculation value 78.5ml. Therefore, this model can well reflect the actual blood flow of the human body.Considering that there is always a certain degree of difference in genetic characteristics of different people, it is easy to misjudge the disease condition if we use a unified health standard to measure the population with different characteristics. Thus, we consider to yield several different standards which can be used to measure the health of different groups of people by adjusting the value of fractional orderα. Figure 5 shows that the stroke volume increases with the increase of α. V has a larger growth rate in the range of α value from 0.4 to 0.8, while it grows slowly in other value ranges. Thus, we can consider selecting the appropriate value of α from 0.4 to 0.8 with the actual test data.Figure 5
Relationship between fractional order and stroke volumeR=1.25×10−2,umax=0.4,c˜=0.Because the radius of the blood vessel can influence stroke volume, we consider changing the value ofR to get several results as shown in Figure 6.Figure 6
Relationship between blood vessel radius and stroke volumeα=0.5,c˜=0.In Figure6, the values of α from bottom to top are 0.2, 0.4, 0.6, and 0.8. We can find that under the condition of the same vessel radius, stroke volume increases with the increase of the maximum velocity of blood vessel center. In addition, under the same umax, the value of V is also related to the size of the vascular radius. As the vessel radius increases, stroke volume increases as well.


## 6. Conclusion

In this paper, we derive time-fractional Boussinesq equation with viscous term to describe blood flow and discuss the effects of fractional order and viscosity on blood flow volume. This indicates that local hemodynamic factors play an important role in the occurrence and development of cardiovascular disease. Therefore, it is of great significance to study the changes of blood flow in blood vessels to reveal the pathogenesis of major vascular diseases.



---

*Source: 9256845-2020-09-15.xml*