REDUCED ORDER MODEL FOR COMPUTING BLOOD FLOW DYNAMICS
20220237864 · 2022-07-28
Inventors
Cpc classification
G16H50/20
PHYSICS
A61B2034/105
HUMAN NECESSITIES
A61B34/10
HUMAN NECESSITIES
International classification
A61B34/10
HUMAN NECESSITIES
Abstract
A computer-implemented method can include generating centerlines of a patient's cardiovascular network, determining geometric features of the cardiovascular network based on the centerlines and a three-dimensional (3D) computer model of the cardiovascular network, constructing a lumped parameter network (LPN) of resistors corresponding to the cardiovascular network, and solving a system of equations corresponding to flow and pressure for the LPN model.
Claims
1. A computer-implemented method for non-invasive assessment of a patient's hemodynamics information, comprising: generating a plurality of centerlines corresponding to at least one portion of a cardiovascular system of the patient; determining a plurality of geometric features corresponding to the at least one portion of the cardiovascular system based on the plurality of centerlines and a three-dimensional (3D) computer model of the at least one portion of the cardiovascular system; constructing a lumped parameter network (LPN) of resistors corresponding to the at least one portion of the cardiovascular system; and solving a system of equations corresponding to flow and pressure for the LPN.
2. The computer-implemented method of claim 1, further comprising providing a 3D representation corresponding to a result of the solving.
3. The computer-implemented method of claim 1, further comprising constructing the 3D computer model of the at least one portion of the cardiovascular system of the patient.
4. The computer-implemented method of claim 3, wherein the 3D computer model of the at least one portion of the cardiovascular system is an aortic model.
5. The computer-implemented method of claim 3, wherein the 3D computer model of the at least one portion of the cardiovascular system is an aorto-femoral model.
6. The computer-implemented method of claim 3, wherein the 3D computer model of the at least one portion of the cardiovascular system is a coronary model.
7. The computer-implemented method of claim 3, wherein the 3D computer model of the at least one portion of the cardiovascular system is a cerebrovascular model.
8. The computer-implemented method of claim 3, wherein the 3D computer model of the at least one portion of the cardiovascular system is a pulmonary model.
9. The computer-implemented method of claim 3, wherein the 3D computer model of the at least one portion of the cardiovascular system is a congenital heart disease model.
10. The computer-implemented method of claim 1, wherein the plurality of geometric features includes a vascular obstruction.
11. The computer-implemented method of claim 1, wherein the plurality of geometric features includes a vessel curvature.
12. The computer-implemented method of claim 1, wherein the plurality of geometric features includes a vessel bifurcation.
13. The computer-implemented method of claim 1, wherein the LPN describes a mathematical relationship between blood flow, blood pressure, and their derivatives.
14. The computer-implemented method of claim 1, wherein blood flow rate is linearly proportional to pressure drop across the network in the LPN.
15. One or more tangible, non-transitory computer-readable media storing instructions that, when executed by a processor, cause the processor to perform the computer-implemented method of claim 1.
Description
BRIEF DESCRIPTION OF THE DRAWINGS
[0009]
[0010]
[0011]
[0012]
[0013]
[0014]
[0015]
[0016]
[0017]
[0018]
[0019]
[0020]
[0021]
[0022]
[0023]
[0024]
[0025]
[0026]
[0027]
[0028]
[0029]
[0030]
[0031]
[0032]
[0033]
[0034]
[0035]
[0036]
[0037]
[0038]
[0039]
[0040]
[0041]
[0042]
[0043]
[0044]
[0045]
[0046]
[0047]
[0048]
[0049]
[0050]
[0051]
[0052]
DETAILED DESCRIPTION
[0053] Implementations of the disclosed technology are generally directed to a distributed lumped parameter (DLP) modeling framework arranged to rapidly compute blood flow and pressure in vascular domains. This can be achieved by developing analytical expressions describing expected energy losses along vascular segments, including from viscous dissipation, unsteadiness, flow separation, vessel curvature and vessel bifurcations. This methodology can be applied to solve for unsteady blood flow and pressure in a variety of complex 3D image-based vascular geometries, which are typically approached using time-dependent computational fluid dynamics (CFD) simulations.
[0054] Implementations of the DLP framework disclosed herein generally demonstrate consistent agreement with CFD simulations for a range of flow conditions and vascular anatomies, with mean errors in flow rate and pressure less than 7% over the range models considered. The computational cost of the DLP framework is orders of magnitude lower than the computational cost of CFD, which opens new possibilities to use hemodynamics modeling for timely decision support in clinical settings, or to improve the ability to perform data assimilation, optimization, parameter tuning, uncertainty quantification and sensitivity analysis for cardiovascular modeling applications.
[0055] Implementations of the disclosed technology are generally directed to a distributed lumped parameter (DLP) framework to compute temporal flow and pressure in different vasculature models by taking into account various sources of energy dissipation. This approach results in ODEs that are significantly easier to solve than the aforementioned 1D and 2D PDE-based models. This approach can be even more accurate depending on the types of energy losses considered.
[0056] Most ROM approaches described above are reasonable in smaller vessels where the flow becomes 1D and fully developed. However, most ROMs produce erroneous results for flow in larger arteries, as typically derived from medical imaging, where flow conditions are known to be complex and unsteady.
[0057] Therefore, to evaluate the accuracy of the proposed DLP methodology, this framework is applied to a diverse range of healthy and diseased patient-specific cardiovascular anatomies including aortic, aorto-femoral, cerebrovascular, coronary, pulmonary and congenital heart disease models, and have compared the DLP results to the “gold-standard” from 3D time-dependent (“3Dt”) CFD simulations. The ability of this framework to be automated provides many advantage for both exploratory and translational applications.
[0058] Similar to most fluid mechanic models, a relationship can be developed between flow and pressure from conservation of mass and balance of momentum principles. Namely, the formulation can be in terms of cross-sectionally averaged pressure P(x) and volumetric flow rate Q(x), where x denotes the axial coordinate along a given vessel. The cross-sectionally averaged Navier-Stokes equations, under suitable assumptions, reduce upon linearization to
[0059] where A.sub.0 is the reference cross-sectional area, ρ is fluid density and f characterizes viscous dissipation. These equations have been widely used to describe flow and pressure in vascular domains, but are generally inaccurate in realistic geometries that contain energy loss due to curvature, flow separation, sudden expansions, etc. Therefore, as a basic adaptation, the dissipation term is redefined on the right hand side of Eq. (2) to account for losses neglected in deriving this equation. Namely, for a vascular segment of length L and (spatially varying) radius R, a momentum balance is considered of the form
[0060] An appropriate vascular resistance can be defined by considering expected sources of energy dissipation, that when used in Eq. (3) and along with Eq. (1) and appropriate boundary conditions, provides accurate estimation of hemodynamics in realistic vascular geometries, particularly those derived from medical imaging and normally approached by CFD simulation.
[0061] Viscous Dissipation
[0062] Assuming Poiseuille flow, the pressure drop across a cylindrical vessel of length L is given by ΔP=R.sub.vQ, where R.sub.v=8 μL/πR.sup.4 is the hydraulic resistance, with μ being the blood viscosity and R being the vessel radius. An integral form of this equation is considered to better account for potential variations of radius along the length of a vessel as
[0063] For vessels of non-circular cross-section, R(x)=√{square root over (A(x)/π)} denotes the effective radius (henceforth radius) at local axial coordinate x∈(0, L).
[0064] The Poiseuille law, ΔP=8 μLQ/πR.sup.4, can be viewed a special case of the Darcy-Weisbach equation
[0065] where V=Q/πR.sup.2 is the average fluid velocity, and λ.sub.s is the viscous friction factor. The viscous friction factor of a straight vessel is λ.sub.s=64/Re, which yields Poiseuille's law, where Re denotes the Reynolds number. This viewpoint can be used next to incorporate pressure loss due to vessel curvature.
[0066] Curvature
[0067] The presence of centrifugal forces in a curved vessel results in a skewed velocity profile, which generally increases viscous dissipation as compared to flow through a straight vessel. Assuming steady and fully developed flow, the viscous friction factor of a curved vessel (λ.sub.c) can be related to the viscous friction factor of a straight vessel (λ.sub.s) with similar length and radius
[0068] where K=Re√{square root over (R/a)} is the Dean number and a is the vessel curvature. Therefore, Eq. (4) can be modified to
[0069] Note that R, a and K generally vary as a function of x.
[0070] Expansions
[0071] The energy dissipation due to sudden expansions can be taken into account by using a semi-empirical model
[0072] where A.sub.s and A.sub.0 are cross-sectional areas describing the expansion as described below, and |Q| denotes the absolute value of flow rate (in case of reverse flow). Eq. (8) with K.sub.t=1 can be derived from a control volume energy balance of inviscid flow through a sudden expansion. K.sub.t is added an empirical correction factor to account for losses from flow separation. A previous study showed that K.sub.t=1.52 produced consistently accurate prediction of pressure drop compared to 3D simulations in realistic vascular (coronary) stenoses.
[0073] In this study, the local minima and maxima of R(x) were computed along each vessel segment. Then A.sub.s was computed at each local minimum. For each local minimum, A.sub.0 was computed from the mean value of the radius at the local maxima immediately proximal and distal to the corresponding local minimum. For vascular segments with multiple sudden expansions, these losses were added in series
[0074] where n is the number of sudden expansion locations over the given vascular segment. Aggregate losses due to sudden expansions were added in series to the viscous pressure loss derived above.
[0075] Bifurcations
[0076] Energy losses at vascular junctions were determined from geometric (branch angle) and hydraulic (flow split) considerations. Namely, the amount of pressure drop due to a vascular bifurcation was calculated as
[0077] where Q.sub.dat and A.sub.dat are the flow rate and cross-sectional area of the datum supplier vessel, respectively. λ.sub.j=Q.sub.j/Q.sub.dat defines the volumetric flow split, ψ.sub.j=A.sub.dat/A.sub.j defines the area split and φ.sub.j=3(π−θ.sub.j)/4 is defined from the angle θ.sub.j between a datum supplier and child branch. Hence the nonlinear resistance
[0078] can be added in series to the resistances due to viscous and sudden expansion effects in the child branch.
[0079] Unsteadiness
[0080] Flow unsteadiness affects the momentum balance in Eq. (3) in two main ways. First, it explicitly leads to a change in fluid inertia as captured in the first term. Second, it implicitly changes the cross-sectional profile of the flow, which in turn changes the shear rate and viscous dissipation captured in the second term. To account for the later, the viscous resistance may be modified according to the Womersley number.
[0081] The Womersley number, α=R√{square root over (ρω/μ)}, quantifies the relative importance between pulsatile inertial effects and viscous effects. The heart rate is generally the dominant frequency and thus used to define ω. For image-based vascular models, ω, ρ and μ vary minimally and thus intra- and inter-model variations of a are due to variations in vessel radius R.
[0082] As the Womersley number increases the velocity profile changes from parabolic to flat, increasing viscous losses. Based on the well-known Womersley's solution for pulsatile flow in a straight, rigid vessel, wall shear stress can be computed and then used to compute the (viscous) resistance due to Womersley flow. This can be used to define a function ζ(α) describing the viscous resistance ratio between Womersley and Poiseuille flow. Hence, the (nominally Poiseuille) viscous resistance introduced in Eq. (4) can be modified to
[0083] where ζ is a function of x because R and hence α is a function of x. In this study, the maximum between curvature effect and pulsatile effect ζ was used in each artery to modify the Poiseuille resistance.
[0084] Implementation
[0085] Based on the above considerations, the resistance term in Eq. (3) can be given by
[0086] and Eq. (3) along with conservation of mass can be used to solve for flow rate through each vascular segment and pressure at vascular junctions.
[0087] Certain embodiments include an in-house Python framework to completely automate the DLP modeling procedure based on image-based geometry. Namely, the only input to the framework might be the geometric model (e.g. STL surface mesh) and specification of boundary conditions. Automated procedures can arranged to compute vessel lengths, bifurcation angles, and the variation of cross-sectional area and curvature along each vascular segment. A variety of boundary conditions can be implemented typical to current state-of-the-art image-based modeling, including specification of a time-varying flow or pressure waveform, or coupled LPN models at the inlet or outlets.
[0088] The end result is generally a system of nonlinear AEs and ODEs that describe the basic conservation laws and coupled boundary conditions. An implicit Euler scheme may be used to handle time-stepping. This results in a set of nonlinear algebraic equations to be solved each time-step, which is handled by iterative linearization. The solution of this system includes flow rate at all vascular segments and pressure values at all vascular junctions and boundaries.
[0089] In certain implementations, the DLP framework may be applied to a variety of image-based vascular models and the results were compared to those derived from CFD. An example included constructing a 3D model of the arteries from corresponding image data, employing and tuning boundary conditions representative of downstream and upstream physiology, finite element simulation of the 3D time-dependent Navier-Stokes equations, and post-processing to compute the flow rate and pressure at locations for comparison with the DLP model. All CFD simulations were run assuming rigid vessel walls. The domain was discretized using linear tetrahedral elements that included boundary layer meshing for all simulations. Mesh independence studies were undertaken for all 3D simulations to ensure that the results in the final meshes were not significantly affected by inadequate mesh resolution. Solutions were run until the pressure fields at the inlets and outlets did not change more than 1.0% compared to the solutions at the same time point in the previous cardiac cycle. Fluid properties (blood density ρ=1.06 g/cm.sup.3 and blood viscosity ρ=0.04 dyn/(cm.sup.3 s)) were assumed common among all models.
[0090] The DLP model is not tuned to the 3D CFD results. To ensure a consistent comparison, the same inlet and outlet boundary conditions were used. To quantify the comparison, the relative errors were calculated between the mean values of the temporal flow rate and pressure from the DLP model against those derived from CFD simulations at the inlets and outlets of each model. Values at the outlets were chosen because this evaluates the ability of the DLP model to serve as a surrogate for the entire 3D domain in applications such as UQ, parameter tuning, etc., and because these locations are generally where integrated differences become compounded (i.e., error is expected to accumulate as greater segments of vasculature is traversed). Nonetheless, the temporal variations in the flow rate and pressure computed from the DLP model were also compared with those from the 3D CFD simulations at several locations in the interior of the domain. For perspective, there was a comparison against results from a basic DLP model for which the Poiseuille resistance was assumed. This model can be viewed as the baseline model (cf. Eq. (2)) for which the corrections to the dissipation term proposed herein are neglected.
[0091] In the example, the DLP framework, on a single CPU core, took about 334±49 seconds to complete steps 102-110 of
[0092] Framework
[0093]
Aortic Models
[0094] Four aortic models 202, 302, 402, and 502 were considered including one normal, two coarctation (moderate and severe), and one aneurysmal, as illustrated by
[0095] For all cases, phase-contrast (PC) MRI was used to measure flow in the ascending aorta. The respective volumetric flow waveforms from PC-MRI for each patient was mapped to the inlet of the CFD models using a time varying parabolic flow profile as the inflow boundary condition. Three-element “RCR” Windkessel models were coupled at all outlets. Regional mesh refinement was used for cases with aortic coarctation to resolve the complex flow features in the stenotic regions.
[0096] Aorto Femoral Models
[0097] Four aorto-femoral models 902, 1002, 1102, and 1202 were considered as illustrated by
[0098]
[0099] For Model A (902), an aortic flow waveform was adapted to have a mean cardiac output of 4.6 L/min (female). For Model B (1002) the supraceliac aorta blood flow waveform was derived from PC-MRI data. For Models C (1102) and D (1202) individualized inflow boundary conditions were determined based on the Baker equation, relating body surface area to cardiac output, and assuming that 70% of the cardiac output is distributed to the supraceliac aorta. The resulting mean flows were used to generate inflow waveforms by scaling a gender-matched representative supraceliac aortic flow waveform. Three element “RCR” Windkessel models were applied at each outlet. The values of the three parameters for each outlet are determined based on flow distributions to the outlets obtained from clinical PC-MRI measurements for Model B (1002), or from literature data for Models A (902), C (1102), and D (1202).
[0100] Coronary Models
[0101] Four patient-specific models 1602, 1702, 1802, and 1902 of the proximal aorta and major coronary arteries were considered as illustrated by
[0102] In all cases, aortic flow was prescribed at the inlet, an “RCR” Windkessel of the systemic circulation was coupled at the aortic outlet, and coronary-specific LPNs that consider the time-dependent intramyocardial pressure were coupled at the coronary outlets (each outlet having a separate such LPN). The effect of intramyocardial pressure is modeled by scaling the typical left and right ventricle pressures to recover realistic coronary flow waveforms. Namely, the left coronary waveform is described by low systolic flow and high diastolic flow, while the right coronary flow demonstrates two peaks one in systole and one in diastole. The LPN parameters of the systemic and coronary outlets were tuned to match target pressure and flow splits to the aorta and systemic and coronary outlets.
[0103] Cerebrovascular Models
[0104] Four cerebrovascular models 2402, 2502, 2602, and 2702 were considered as illustrated by
[0105] For all models, a characteristic vertebral blood flow waveform was scaled to match time-averaged PC-MRI measurements in the vertebral arteries and used to prescribe an inflow boundary condition to the computational models assuming a parabolic profile. Resistance boundary conditions were used at the outlets. A total resistance was calculated and distributed amongst the outlets by assuming all outlets act in parallel with values inversely proportional to the outlet area.
[0106] Pulmonary Models
[0107] Four patient-specific pulmonary models 3102, 3202, 3302, and 3402 extended from the main pulmonary artery to various levels of branching in the left and right pulmonary arteries (LPA, RPA) were considered as illustrated by
[0108] Pulmonary blood flow waveforms 3104, 3204, 3304, and 3404 from PC-MRI are applied to the inlet of each model illustrated by
[0109] Congenital Heart Disease Models
[0110] Models 3802, 3902, 4002, and 4102 from four pediatric patients with congenital heart disease were considered as illustrated by
[0111] PC-MRI data was used to prescribe an inflow waveform to the inlets of computational domains. Inflow waveforms prescribed at the inferior and superior vena cava (IVC, SVC), internal jugular vein (IW) and brachiocephalic vein (BrS) are shown in
DISCUSSION
[0112] Implementations of the disclosed technology are generally directed to a DLP framework to predict temporal flow rate and pressure waveforms distributions in 3D image-based vascular models. The proposed DLP framework can provide consistent prediction with CFD simulations with mean errors <7% in a range of cardiovascular modeling applications. This framework can be fully automated based on an input model geometry (and specification of boundary conditions), and generally requires less than 1/1000 of the computational cost of corresponding CFD simulations.
[0113] The instant disclosure includes the most comprehensive comparisons of a ROM to the current state-of-the-art image-based CFD modeling. A broad range of vascular anatomies and variations in geometric features were considered. Namely, the two key non-dimensional parameters of Reynolds (Re) and Womersley (α) numbers varied from ˜100-3,000 and ˜1-20, respectively. Notably, these ranges in Re, α, and vessel diameter span the spectrum of expected values encountered in most image-based computational modeling of blood flow. In addition to the dynamic and geometric difference explored, the number of vascular junctions across the models varied from 4 to 96. Interestingly, the junctions include both converging and diverging flows, where converging flows sometimes result in a negative loss coefficient, and in presence of back flow, some bifurcations can have both diverging and converging flows at different time points of the cardiac cycle. Additionally, vascular expansions, as quantified by the area reduction ratio A.sub.s/A.sub.0, ranged up to ˜90%, and the mean curvature ratio δ=a/R ranged up to −0.5 in the considered models, which again cover a broad range of expected values.
[0114] To provide more general insight into the contribution of each source of energy dissipation, the modification due to flow separation downstream of expansions, when present, may have the highest contribution to accurately predict flow and pressure in the corresponding model. Viscous dissipation introduced by curvature effects generally has the second highest contribution in energy losses. This can be confirmed by comparing the second and third panels of error plots for coronary models with negligible unsteadiness effect, particularly for the healthy coronary Model C or Models B and D with mild stenoses.
[0115] Moreover, for cerebrovascular models the results are improved mostly by considering the curvature effects, since the Womersley number is ≈2 in these models and there is no sudden changes in cross-sectional area. The Womersley effect has the next highest impact on modifying Poiseuille resistance model, where viscous resistance can increase by a factor of 5 in for example aortic branches. The Womersley effect can be seen for example in the aorto-femoral models where the errors approximately tripled without considering this effect. Finally, the energetic losses due to bifurcations was found to be relatively negligible, except in cases where numerous junctions exist, such as pulmonary models.
[0116] Given the accuracy of the DLP approach, implementations of the disclosed technology can potentially replace CFD in applications where time is of the essence, thus opening the door to broader use of image-based modeling in, for example, clinical settings. A major impact of this approach is to applications where numerous CFD simulations would be desirable, such as for model parameter or boundary condition tuning, uncertainty quantification, sensitivity analysis, surgery/device design, data assimilation, machine learning training, etc. In such scenarios, it could be the de facto model, or used in conjunction with more expensive CFD simulations to reduce the number of full CFD model evaluations. This could be highly significant as the above types of engineering analyses are playing an increasing role in image-based modeling research and translational applications.
[0117] The following discussion is intended to provide a brief, general description of a suitable machine in which embodiments of the disclosed technology can be implemented. As used herein, the term “machine” is intended to broadly encompass a single machine or a system of communicatively coupled machines or devices operating together. Exemplary machines may include computing devices such as personal computers, workstations, servers, portable computers, handheld devices, tablet devices, and the like.
[0118] Typically, a machine includes a system bus to which processors, memory such as random access memory (RAM), read-only memory (ROM), and other state-preserving medium, storage devices, a video interface, and input/output interface ports can be attached. The machine may also include embedded controllers such as programmable or non-programmable logic devices or arrays, Application Specific Integrated Circuits (ASICs), embedded computers, smart cards, and the like. The machine may be controlled, at least in part, by input from conventional input devices such as keyboards and mice, as well as by directives received from another machine, interaction with a virtual reality (VR) environment, biometric feedback, or other pertinent input.
[0119] The machine may utilize one or more connections to one or more remote machines, such as through a network interface, modem, or other communicative coupling. Machines can be interconnected by way of a physical and/or logical network, such as an intranet, the Internet, local area networks, wide area networks, etc. One having ordinary skill in the art will appreciate that network communication may utilize various wired and/or wireless short range or long range carriers and protocols, including radio frequency (RF), satellite, microwave, Institute of Electrical and Electronics Engineers (IEEE) 545.11, Bluetooth, optical, infrared, cable, laser, etc.
[0120] Embodiments of the disclosed technology may be described by reference to or in conjunction with associated data including functions, procedures, data structures, application programs, instructions, etc. that, when accessed by a machine, may result in the machine performing tasks or defining abstract data types or low-level hardware contexts. Associated data may be stored in, for example, volatile and/or non-volatile memory, such as RAM and ROM, or in other storage devices and their associated storage media, which can include hard-drives, floppy-disks, optical storage, tapes, flash memory, memory sticks, digital video disks, biological storage, and other non-transitory, physical storage media.
[0121] Associated data may be delivered over transmission environments, including the physical and/or logical network, in the form of packets, serial data, parallel data, etc., and may be used in a compressed or encrypted format. Associated data may be used in a distributed environment, and stored locally and/or remotely for machine access.
[0122] Having described and illustrated the principles of the invention with reference to illustrated embodiments, it will be recognized that the illustrated embodiments may be modified in arrangement and detail without departing from such principles, and may be combined in any desired manner And although the foregoing discussion has focused on particular embodiments, other configurations are contemplated.
[0123] Consequently, in view of the wide variety of permutations to the embodiments that are described herein, this detailed description and accompanying material is intended to be illustrative only, and should not be taken as limiting the scope of the invention. What is claimed as the invention, therefore, is all such modifications as may come within the scope and spirit of the following claims and equivalents thereto.