Simulation of flow in an artery under pathological hemodynamic conditions: The use of a diagnostic disease descriptor.

A numerical model for simulating and predicting blood flow dynamics in diseased arterial vessels has been developed. The time-dependent one-dimensional hyperbolic system of quasilinear partial differential equations which incorporates a diagnostic disease descriptor ( k D ) was used to simulate transient flow distribution for idealized healthy and diseased states. Blood flow simulations in the iliac arteries over about 125% of a cardiac cycle were generated and calibrated using the k D values from 0 to 3 representing hypothetical diseased states. Early results indicate that disease conditions induce abnormal flow in the artery, generating disorder and increased amplitude of blood pressure, flow and distensibility with increasing numerical values of the disease factor k D . More so, the prospective use of the k D -approach with documentation of in vivo adverse flow visualizations for diagnostic purposes was decisively discussed.

Introduction

Arterial diseases are renowned precursors to killer cardiac conditions like strokes and heart attacks. With a great number of human mortalities in the developed world attributed to abnormal blood flow in arteries (Liu, 2000), investigative studies on the pathological flows in arterial vessels are being increasingly pursued by researchers in recent times (Cornet, 2018; Gamilov et al., 2019; Bunonyo et al., 2020; Bertaglia et al., 2020; Lopes et al., 2020; Leguy, 2019). Besides investigations which contribute to elucidating the complex physiologic processes involved in the genesis and progress of these deadly vascular diseases, others that offer potential clinical assistance for their easy diagnosis and treatment are especially important.

While the pathological flow of blood is now widely studied, hydrodynamic investigations of the flow of blood under normal physiologic conditions provided directions for descriptive studies on diseased flow systems; this is an adventure with a long history. Quarteroni (2001a) chronicled the important studies and mathematical developments in cardiovascular flow research of several centuries that led to the evolution of its present-day sophisticated methods and models. It was noted that the now popular Bernoulli and Poiseuille equations in fluid mechanics owe credits to the 18th and 19th century inquiries into blood flow of these great scientists respectively. With progressive research, the quantitative study of vascular flows in the Sixties was approached from the electrical analogy perspective (Formaggia et al., 1999) enabling investigators to build complex network models to simulate macroscopic flow in the circulatory system (see, for example, Westerhof et al., 1969). However, animal experiments held out in the Seventies as the main mode of cardiovascular investigations (Quarteroni, 2001a).

Now, numerical and computational techniques have the day, increasingly being employed to expound vascular hemodynamics. With the advent of sophisticated computing technologies is the unfolding of the finest models for adverse flow representation in diseased vascular regions. Several models of different scales and complexities coupled to depict flows in large parts or the whole sections of the cardiovascular system are now being simulated via the aid of computational tools – a procedure known as multiscale modelling. Geometrical multiscale modelling, as reviewed by Team REO (2004), is a numerical modelling concept advanced by the need for accurate description of flow in sensible pathological vascular regions as well as the critical requirement of boundary conditions for such regions. While Formaggia et al. (1999) provided a preliminary analysis of the multiscale model; Quarteroni and Veneziani (2003) did explore some of its applications.

Besides aiding in the understanding of poorly understood disease processes of flow, computational studies of blood flow are also being used to predict the outcome of vascular surgeries and in the development of treatment plans. In their paper, Taylor et al. (1999) promulgated the new paradigm of “predictive medicine”, a simulation-based medical planning procedure through which alternative treatment plans for a patient are tested in a virtual environment. Quantitative studies which furnish mathematical models and numerical simulations of disease-induced transient arterial flows are geared towards this end; this goal is further progressed in this work. A numerical blood flow model incorporating a diagnostic disease descriptor is introduced to calibrate diseased states in arteries for clinical diagnosis procedures and treatment. This new concept has been demonstrated using a simulation software utilized with published physiologic data on the iliac arteries.

Model development and numerical analysis

Although the complex three-dimensional Navier-Stokes model equations can provide full details of an arterial flow field, a major drawback is that they are computationally expensive (Quarteroni, 2001b). Three-dimensional numerical models can deliver full flow-field simulations of blood flow across and along the three spatial directions but this comes at a very high computational cost of simulation and power.

As such, the use of one-dimensional models which are of intermediate level of complexity to simplify calculations and save computing costs are on the increase (see, for example, Canic et al., 2006; Formaggia et al., 2003; Matthys et al., 2007; Olufsen et al., 2000; Perdikaris and Karniadakis, 2014; Reymond et al., 2009; Smith et al., 2002; Canic et al., 2004; Formaggia et al., 1999). Likewise, dynamic flow in an artery is here modelled by one-dimensional continuity and flow motion equations obtained by applying basic mass and force balances to a short differential segment of blood. Detailed mathematical analysis leading to the derivation of these is presented elsewhere (Oshin, 2006).

Although there are other one-dimensional models for arterial flow in literature, none of the models proposed or incorporated a numerical index which can help to indicate the progress and influence of disease on the flow.

The nomenclature and list of symbols used in the mathematical analysis that follows in this section is given in Table 1.

Table 1

List of mathematical symbols.

Symbol	Meaning
a	Constant defined in equation (3)
A	Artery cross-sectional area
A0	Cross-sectional area of artery in the unstretched state
B	Constant in equation (49)
c	Variable defined in equation (21)
C	Characteristic curve
D	Artery diameter
D0	Diameter of artery in the unstretched state
e	Arterial wall thickness
E	Young's modulus
F	Variables defined in equations (32) and (34)
G	Variables defined in equations (33) and (35)
k0	Constant defined by equation (19)
kD	Disease factor
K	Blood bulk modulus of elasticity
K⁎	Combined modulus of elasticity
L	Length of artery
M	Matrix in equation (9)
M′	Matrix in equation (13)
N	Matrix in equation (10)
N′	Matrix in equation (14)
p	Internal blood pressure
p0	Extravascular pressure
pss	Steady state flow pressure
q	Flow rate
qss	Steady state flow rate
t	Time
Δt	Time increment
T	Cardiac period
ν	Variable defined by equation (20)
x	Longitudinal dimension
Δx	Increment in length
Y	Variable vector in equation (8)
Greek symbols
λ	Variable defined in equations (22) and (23)
μ	Blood viscosity
ρ	Blood density
σ	Poisson ratio
τ0	Wall shear stress
Superscripts
+	Pertaining to slope of characteristic curve
−	Pertaining to slope of characteristic curve
Subscripts
i	Referring to grids along vessel length

Model development

The continuity model for the flow is obtainable as equation (1): The internal blood pressure p and the flow rate q are both dependent on time t and the distance x along the artery measured from the upstream end; and . K is the bulk modulus of elasticity of blood. The cross-sectional area of the artery, A is related to the net pressure on the arterial wall by with representing the external pressure on the wall and the cross-sectional area at zero net pressure. The zero-order equation (2), obtained by considering the static equilibrium of the arterial wall in the radial direction, is one in a family of structural models used to depict the elastic behaviour of vascular walls. Other models with higher complexity and order have been used by researchers; for instance, the Voigt model which takes into account viscoelastic behaviour was employed by Formaggia et al. (1999).

The Poisson ratio σ, Young's modulus E and wall thickness e of the artery are absorbed by the constant a in (2) above defined by equation (3) thus: Now, an equivalent bulk modulus which expresses the combined elasticity of the blood and the vessel wall, indicated by the first and second terms of equation (4) respectively, is defined as: so that the continuity equation takes the final form Furthermore, the flow motion model is presented as equation (6) where ρ and μ retain their usual identities as density and viscosity respectively.

The relevant simplifying assumptions made in arriving at (6) include that

The artery is horizontal, parallel to the chosen horizontal datum.

The frictional resistance for unsteady flow and hence the shear stress at the vascular wall can be substituted for by that for laminar steady flow obtained from the Poiseuille equation.

The continuity and flow motion equations: (5) and (6), constitute the model describing the transient one-dimensional flow of blood in an artery both of which can be jointly written in vector-matrix form of equation (7) (Oshin, 2006): where Y, and are as defined in equations (8), (9) and (10) respectively: The system (7) is sufficient to capture normal haemodynamic changes; however, the forcing effect of progressing vascular disease on the flow has not been reflected. Not only is it necessary to provide a model modification for sufficiently accurate simulation of disease-induced adverse blood dynamics, it is also required that such a model simulation be able to deliver an acceptable measure of disease quantification.

In clinical routines, many indices used for diagnosis often give false estimations, there are anticipations that numerical simulations may give indications to define new indices, simple enough to be used in clinical examinations, but more precise than those currently used (Team REO, 2004). In this work therefore, a new index is introduced for the calibration, diagnosis and numerical quantification of vascular disease. This is a dimensionless diagnostic disease descriptor which is used to modify the one-dimensional dynamic model of flow in an artery presented above in (7). This approach was employed by Abhulimen and Susu (2004, 2007) in modelling leak in liquid pipelines.

If there is a disturbance causing instability of blood flow in the artery, incorporating an arterial disease descriptor to the artery flow system (7) results in a modified form of the model equations re-presented as Dividing equation (11) by results in equation (12) where and are as defined by equations (13) and (14): These modified one-dimensional quasilinear partial differential equations describe blood flow in a diseased artery.

Numerical analysis

The method of characteristics solution of transient flow equations has been derived and well-documented in standard literatures and book monograms (see Douglas et al., 1995; Streeter et al., 1998). When this is applied to the hyperbolic flow model (12) assuming ideal incompressibility of blood so that the elasticity moduli become, the result is these pairs of ordinary differential equations: subject to And subject to This analysis generated new variables, hybrids of the key ones, worth defining for their physical significance to the understanding of the flow phenomenon. They are as defined in equations (19)–(23). The variables ν and c have comparable features to the flow and sonic velocities respectively. They have been so defined for ease of analysis. The variable pair and with pseudo-eigenvalue identities satisfy the linear combination of the flow motion and continuity equations in a manner defined thus by equation (24): where equation of motion and is given by equation (25) and equation of continuity which is given by equation (26) It is essential to know that equations (15) and (17) are only valid along the characteristics curves described by (16) and (18) respectively; these are labelled C+ and C− in Fig. 1. The numerical technique is such that the sought transient flow simulation in the diseased arterial vessel in the x-t plane is achieved by dividing it into an even number of N sections, each of length (Fig. 2). It is projected that flow measurements between two numerical nodes in the diseased artery will give flow signatures different from those received from healthy arteries so that the characteristic distance-time curve is used to simulate the sickly conditions of flow at regular time steps Δt defined by equation (27). The general equations therefore for computing flow rate, q and pressure, p at each space-time nodes from 2 to N derived from the discretization via a first-order finite difference approximation of equations (15) to (18) are where and are defined by equations (30) and (31) respectively The terms in [] refer to conditions at time t+Δt and subscript i refers to the artery section number. For compactness, and combine known upstream variables at a previous time node as given by equations (32) and (33): While and similarly bring together flow variables at an adjacent downstream node one time-step earlier according to equations (34) and (35): Now, equations (28) and (29) become respectively thus: With , , , , and known, the solution of equations (36) and (37) is given by equation (38): The flow simulation consists of finding p and q for alternate grid points along , then proceeding to , and so on, until the desired time duration has been covered. End or boundary conditions of the vessel are introduced every other time step after the initial conditions.

Figure 1

Characteristic curve in the x-t plane.

Figure 2

Characteristic solution in the x-t grid for an arterial vessel.

At the end of the artery only one of the compatibility equations is available in the two variables. At the upstream end (Fig. 3(a)), equation (37) holds along the C− characteristic, and for the downstream end (Fig. 3(b)), equation (36) is valid along the C+ characteristic. These are linear equations in and . Each conveys to its respective boundary the complete behaviour and response of blood in the artery during the transient. An auxiliary equation is needed in each case, which specifies p, q or some relation between them. Such supplementary equation is the boundary condition. Although boundary condition specification has been noted to be one of the most challenging issues in modelling blood flow (Taylor, 2000), end conditions were carefully formulated for the sample simulation in the iliac artery that follows.

Figure 3

Boundary conditions. (a) Upstream end. (b) Downstream end.

Computer simulation of flow in a diseased iliac artery

The authors chose to apply the developed model for the simulation of blood flow in the iliac artery due to overwhelming information which established it as a pathological region and the available data on its anatomic and physiologic parameters. Clinically, it is observed that atherosclerotic disease develops first in the abdominal aorta, and is much more common in this section of the aorta below the diaphragm than the segment of the aorta above the diaphragm, the thoracic aorta (Taylor, 2000). As an extension of the abdominal aorta, plaque deposition would normally occur also in either of the iliac arteries supplying blood to the legs, which is a cause of lower extremity occlusive disease as shown in Fig. 4 (Society for Vascular Surgery, 2007).

Figure 4

Iliac artery narrowed by plaque (Society for Vascular Surgery, 2007).

To achieve the simulation, a total of four boundary conditions are required; a possible list includes the following indicated by equations (39) to (42): (a) The initial conditions of pressure and flow rate in each of the vessels can be obtained from the steady-state analysis. Under steady state condition, the transient equations of continuity and motion respectively become equations (43) and (44): Integrating equation (43) through the entire length of the vessel informs equation (45) where the superscript “ss” signifies steady state. To obtain this constant rate of flow at steady state, information about the cardiac output is required. Cardiac rate and stroke volume can be taken as 72 beats/min and 70 ml/beat respectively (Ganong, 2003), so that the cardiac output is 5040 ml/min or 8.4×10−5 m3/s which is the amount of blood delivered by the left ventricle and conveyed by the aorta. It would be assumed that about half of this enters the abdominal aorta. Taylor (2000) also reported that under resting conditions, nearly 30% of the blood that enters the abdominal aorta flows down the infrarenal segment through the bifurcation into the legs, the iliac arteries. So the steady state flow through the vessel is taken as 1.26×10−5 m3/s.

Initial conditions

Boundary conditions

The initial condition of pressure is given as equation (46) When this is integrated using exit condition, the following nonlinear relation (equation (47)) results: where B is given by equation (48) (b) Pulsatile flow boundary condition, derived from Womersley theory, is prescribed for the inflow boundary by equation (49): where T is the time period in seconds per heart beat (≈0.8s).

Zero net pressure boundary condition is imposed at the outflow boundary as found in (Taylor, 2000) according to equation (50): The anatomic dimensions of the iliac arteries and the properties of the blood through them (Table 2) were obtained from (Canic et al., 2004).

Table 2

Parameters for the iliac arteries (Canic et al., 2004).

PARAMETERS	VALUES
Characteristic diameter, D_0	0.005 m
Characteristic length, L	0.065 m
Wall thickness, e	0.002 m
Young's modulus, E	105
Poisson ratio, σ	0.5
Dynamic viscosity, μ	3.5x10−3 kg/m-s
Blood density, ρ	1050 kg/m3

The extra-vascular pressure, may be assumed to be essentially equal to the atmospheric environmental pressure (Sandquist et al., 1982). So Pa.

After the careful resolution of all the required parameters and boundary conditions as spelt out above, the computer program “ARTSIM” depicted by the flowchart of Fig. 5 was utilised to simulate the transient blood flow in the iliac artery over a period of about 125% of a cardiac cycle at zero and integer values of the disease quantification factor to reveal remarkable results discussed below.

Figure 5

Simulation Program ARTSIM.

Results and discussion

Pressure and flow variation with disease progression

The manner in which arterial disease, here quantified by , affect flow of blood in arteries was typified by the simulation results. At each node position considered, the pressure wave at (which is equivalent to a no-disease normal flow situation) has an initial period of irregular propagation which smoothens out eventually to give a smooth wave. With increasing , as can be observed from Figure 6a, Figure 6b, Figure 6c, Figure 6d, irregularity in flow escalates depicting increasing turbulence. The turbulence in flow with increasing is also accompanied by increasing pressure wave amplitude. This reveals that with progress of disease, arterial blood pressure increases, a common effect readily observable in patients with vascular diseases such as the narrowing down of arterial lumen by plaques (atherosclerosis).

Figure 6a

Flow Pressure comparison for different k values at x = 0.00000 m.

Figure 6b

Flow Pressure comparison for different k values at x = 0.01625 m.

Figure 6c

Flow Pressure comparison for different k values at x = 0.03250 m.

Figure 6d

Flow Pressure comparison for different k values at x = 0.04875 m.

From Figure 6a, Figure 6b, Figure 6c, Figure 6d, it can also be observed that at a given hypothetical disease state (indicated by a given numerical value of ), the pressure decreases along the length of the artery from upstream to downstream. This is typical of arterial flows where pressure change drives blood flow from the upstream end to the downstream end.

The proposed index () for the quantification of disease in arterial flows holds prospect as a useful index for a scaled calibration of extent of disease in arterial flow. It is a simple “grading scale” of arterial flow abnormality proposed such that very low, moderate and high values respectively depict mild, medial and maximally-developed diseased conditions of flow. It is projected that the can be utilized clinically for diagnosis using it with already established methods and techniques of performing in vitro visualizations of in vivo vascular flows. Since abnormal flow visualizations and measurements in diseased arteries may be electronically documented; such rich clinical observatory databank correlated and calibrated with a simple index as the disease factor proposed here holds great prospect for easy diagnosis of arterial disease.

The modifying factor is introduced into the model equation as a diagnostic disease descriptor to capture the changes induced into the flow phenomenon by progressing disease. While the changes may be specific or generic in nature, it is only captured hypothetically or ideally and not specific to a named disease or specific patient. As such, the work is rather a rough initial trial of a novel method of studying or quantifying the effect of disease on arterial blood flow. Future work will be devoted to a more expansive study on the quantification of the effects of disease on key components or physical terms of the model equations using the same methodology of introducing hypothetical modifiers or weighting parameters to capture the pathology-induced flow in an artery.

Variation of distensibility with disease progression

Distensibility is a measure of the stiffness or rigidity of an artery. Its influence on the flow of blood through the artery lumen is important especially in disease conditions such as atherosclerosis. Mathematically, distensibility is defined as Using the structural model depicting the elastic behaviour of the arterial wall, equation (2), distensibility can be written in terms of flow pressure as The numerical estimates of the distensibility is as shown in Figure 7, Figure 8. With increasing kD values, the Figures show that distensibility increases confirming that with the progression of disease, the artery becomes more stiff or rigid. Fig. 9 depicts the relationship between distensibility and the disease factor kD. The curve fits well to a quadratic relationship with R2 value of 0.9936 showing that arterial stiffness or rigidity is amplified by progressing disease.

Figure 7

Plot of Numerical estimates of distensibility over numerical range of time for k = 0 and k = 1.0.

Figure 8

Plot of Numerical estimates of distensibility for different time node points for k = 2.0, 3.0 and 5.0.

Figure 9

Correlation of distensibility and disease factor k.

Conclusion

This work has established that arterial disease can be quantified in scaled degrees, generating simulations for the gradual deterioration of flow. A numerical index was introduced to quantitatively indicate the presence and progression of arterial disease and its interference with normal hemodynamic flow. It has also put forth a direction of research for which a simple “grading scale” of arterial flow abnormality is sought to be designed so that very low, moderate and high values respectively depict mild, medial and maximally-developed diseased conditions of flow. With established methods and techniques for performing in vitro visualizations of in vivo vascular flows already existent, most accurate of which include the photochromic dye tracer technique (Ojha et al., 1988; Couch et al., 1996), abnormal flow visualizations and measurements in diseased arteries may be electronically documented. Such rich clinical observatory databank correlated and calibrated with a simple index as the disease factor proposed here holds great prospect for easy diagnosis of arterial disease.

Declarations

Author contribution statement

Oshin, Temitope A.; Abhulimen, Kingsley E.: Conceived and designed the experiments; Performed the experiments; Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data; Wrote the paper.

Funding statement

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

Data availability statement

No data was used for the research described in the article.

Declaration of interest's statement

The authors declare no conflict of interest.

Additional information

No additional information is available for this paper.