2. Left ventricle

The left ventricle is one of the four heart chambers, and (in the normal circulation) one of two pumps that push blood on its long journey to ensure the nourishment and oxygenation of organs and tissues. It is located in the heart's lower left portion and is roughly conical in shape, with a thick muscular wall that is essential for generating the necessary force to pump oxygen-rich blood throughout the body. During one cardiac cycle, the left ventricle undergoes roughly four phases: passive filling (diastole), isovolumetric contraction, ejection (systole) and isovolumetric relaxation. In the diastolic phase, the ventricular muscle is relaxed, with oxygenated blood flowing in from the left atrium. The mitral valve, located between the left atrium and the left ventricle, is open, while the aortic valve, which connects the left ventricle and the aorta, is sealed due to the higher pressure in the aorta. As the ventricular muscles contract, the pressure increases sharply from 3–12 to 100–140 mmHg within a short time span of 200 ms. The volume of blood within the chamber remains constant since two valves are closed in the isovolumetric contraction phase. Then, as the left ventricular pressure exceeds that in the aorta, the aortic valve opens and blood is pumped into the aorta. A subsequent muscle relaxation lowers the pressure, and after isovolumetric relaxation, the mitral valve reopens, refilling the left ventricle. This cycle produces what is known as the ‘pressure-volume loop’, with the enclosed area quantifying the work done during each heart beat.

Due to the substantial change in one cardiac cycle, the flow inside the left ventricle is strongly unsteady. For example, during the rapid passive filling phase, the mitral jet, emanating from the mitral valve (diameter ${\approx }2\ \text {cm}$), generally has a peak velocity of 0.8 $\text {m s}^{-1}$ and attains a Reynolds number of 4000. For a classical jet, this would entail a turbulent dynamics. An annular vortex is formed downstream of the mitral orifice resulting from the roll-up of the shear layer emanating from the leaflets of the mitral valve. This three-dimensional vortex ring develops into a large-scale, asymmetric, clockwise rotating vortex at the centre of the left ventricle by the end of diastole (Reference Le, Sotiropoulos, Coffey and KeefeLe et al. 2012; Reference Le and SotiropoulosLe & Sotiropoulos 2012; Reference Fortini, Querzoli, Espa and CenedeseFortini et al. 2013). This facilitates filling by reducing convective losses (Reference Charonko, Kumar, Stewart, Little and VlachosCharonko et al. 2013). Its rotation also supports the ejection of blood flow toward the aortic outflow track during systole with minimal flow energy dissipation (Reference Pedrizzetti and DomenichiniPedrizzetti & Domenichini 2005).

Considering this intricate haemodynamics, both in vivo imaging as well as computational modelling offer promising ways to non-invasively quantify left ventricular haemodynamics and assess both healthy and pathological behaviour. The imaging of such cardiac flows is traditionally performed using Doppler echocardiography and, more recently, using ‘four-dimensional (4-D) flow MRI’. The latter is a relatively recent development in cardiac magnetic resonance imaging (MRI) that allows the imaging of time-resolved flows and anatomies in three dimensions (Reference Markl, Frydrychowicz, Kozerke, Hope and WiebenMarkl et al. 2012; Reference Stankovic, Allen, Garcia, Jarvis and MarklStankovic et al. 2014; Reference Azarine, Garçon, Stansal, Canepa, Angelopoulos, Silvera, Sidi, Marteau and ZinsAzarine et al. 2019; Reference Garcia, Barker and MarklGarcia, Barker & Markl 2019). Four-dimensional flow MRI has been used extensively in the imaging of cardiac flows as well as flows within several large blood vessels to identify aberrant flow patterns and their relationship to cardiac, vascular and valvular diseases. While current clinical use remains confined to large academic centres, new imaging protocols and analysis techniques are poised to move this imaging into the clinical standard of care. However, a challenge with non-invasive in vivo flow imaging is its limited resolution, which does not allow fine-grained imaging of all flow features as well as flows in smaller blood vessels. Computational modelling is therefore a viable alternative to quantify left ventricular haemodynamics accurately and non-invasively.

However, performing realistic and accurate simulations still poses considerable challenges; accurate simulations of left ventricular flow require patient-specific geometries, physiological wall motion generated by active tissue mechanics and fully coupled fluid–solid interaction with the inclusion of the mitral and aortic valves at the inlet and outlet. Significant progress has been made by using various simplifying assumptions to reduce the complexity of these models given specific research interests. For example, the flow can be simplified using lumped-parameter models (see § 7) to efficiently obtain pressure and volume curves, but these models are inherently incapable of providing the full flow field inside the chamber (Reference Hirschvogel, Bassilious, Jagschies, Wildhirt and GeeHirschvogel et al. 2017; Reference Augustin, Gsell, Karabelas, Willemen, Prinzen, Lumens, Vigmond and PlankAugustin et al. 2021). The geometry of the left ventricle can be idealized (Reference Cheng, Oertel and SchenkelCheng, Oertel & Schenkel 2005; Reference Domenichini, Pedrizzetti and BaccaniDomenichini, Pedrizzetti & Baccani 2005; Reference Meschini, De Tullio, Querzoli and VerziccoMeschini et al. 2018; Reference Dedè, Menghini and QuarteroniDedè, Menghini & Quarteroni 2021) or patient specific (Reference Le and SotiropoulosLe & Sotiropoulos 2012; Reference Vedula, Seo, Lardo and MittalVedula et al. 2016; Reference Hirschvogel, Bassilious, Jagschies, Wildhirt and GeeHirschvogel et al. 2017). In this context, it is worth noting that the interior surface of the cardiac chambers is not smooth. It features papillary muscles, trabeculae of different sizes and false tendons, presenting a sponge-like structure extending into the chamber (see § 11). From simulations that have incorporated these detailed endocardial structures from in vivo imaging, they are thought to play a role in ensuring complete closure of the mitral valve, efficient filling and ejection of blood as well as reduced stress in the myocardium wall (Reference Madu and D'CruzMadu & D'Cruz 1997; Reference Kulp, Gao, Zhang, Qian, Voros, Metaxas and AxelKulp et al. 2011; Reference Chnafa, Mendez and NicoudChnafa, Mendez & Nicoud 2016; Reference Lantz, Henriksson, Persson, Karlsson and EbbersLantz et al. 2016; Reference Vedula, Seo, Lardo and MittalVedula et al. 2016; Reference Sacco, Paun, Lehmkuhl, Iles, Iaizzo, Houzeaux, Vázquez, Butakoff and Aguado-SierraSacco et al. 2018b). Another simplification in studies that have focused on ventricular filling is the absence of the mitral and aortic valves (Reference Cheng, Oertel and SchenkelCheng et al. 2005), or simplified valves incorporated as resistive lumped-parameter elements (Reference Dedè, Menghini and QuarteroniDedè et al. 2021). More recently, 3-D valves, simulated using Arbitrary Lagrangian–Eulerian or immersed boundary methods, have been used to capture realistic flow dynamics in the healthy or pathological left ventricle (Reference Vedula, Seo, Lardo and MittalVedula et al. 2016; Reference Meschini, De Tullio, Querzoli and VerziccoMeschini et al. 2018; Reference Meschini, Viola and VerziccoMeschini, Viola & Verzicco 2019). Finally, the ventricular wall motion can also be simplified as prescribed motion on moving meshes or boundaries constructed from imaging data (Reference Domenichini, Pedrizzetti and BaccaniDomenichini et al. 2005; Reference Chnafa, Mendez and NicoudChnafa, Mendez & Nicoud 2014). However, to model the contraction and relaxation of tissues accurately, fluid–solid interaction simulation of the left ventricle can also be expanded to include electrophysiology and active contraction (Reference Vigmond, Clements, McQueen and PeskinVigmond et al. 2008; Reference Quarteroni, Lassila, Rossi and Ruiz-BaierQuarteroni et al. 2017; Reference Viola, Meschini and VerziccoViola, Meschini & Verzicco 2020). We recommend the review of Reference VerziccoVerzicco (2022), which describes the physiology and numerical methods coupling electrophysiology, solid and fluids of the whole heart in detail.

The high Reynolds number in ventricular flow raises questions about the validity of the laminar flow assumption. Research, including both in vivo studies (Reference Zajac, Eriksson, Dyverfeldt, Bolger, Ebbers and CarlhällZajac et al. 2015) and numerical methods (Reference Chnafa, Mendez and NicoudChnafa et al. 2014; Reference Dedè, Menghini and QuarteroniDedè et al. 2021), indicates that turbulence occurs during late diastole and is initiated by the breakup of large vortex structures. Turbulence models, such as large eddy simulation and Reynolds-averaged Navier–Stokes, have shown improved results in left ventricular haemodynamic analyses compared with laminar flow models and experimental data, due to their ability to capture transient or turbulent states involving small-scale phenomena that cannot be retrieved under the laminar hypothesis (Reference Chnafa, Mendez and NicoudChnafa et al. 2016; Reference Jahanzamin, Fatouraee and Nasiraei-MoghaddamJahanzamin, Fatouraee & Nasiraei-Moghaddam 2019). In fact, a recent study by Reference Saqr, Tupin, Rashad, Endo, Niizuma, Tominaga and OhtaSaqr et al. (2020) demonstrated that turbulence is inherent to cardiovascular flows, not only in the left ventricle but throughout the body. However, it is important to note that fully developed turbulence is not typically observed under physiological conditions (Reference Chnafa, Mendez and NicoudChnafa et al. 2014; Reference Dedè, Menghini and QuarteroniDedè et al. 2021). Large eddy simulation approaches are generally better suited for transitional flow regimes than Reynolds-averaged Navier–Stokes approaches, since the latter assume that turbulence is fully developed. In the context of finite-element modelling, the variational multiscale method (Reference Hughes, Feijóo, Mazzei and QuincyHughes et al. 1998; Reference Hughes, Mazzei and JansenHughes, Mazzei & Jansen 2000; Reference Bazilevs, Calo, Cottrell, Hughes, Reali and ScovazziBazilevs et al. 2007), which is a stabilization technique to control numerical oscillations associated with high Reynolds numbers, models turbulence in a similar manner to large eddy simulation and has been widely applied in both cardiac flow (Reference Quarteroni, Lassila, Rossi and Ruiz-BaierQuarteroni et al. 2017; Reference Dedè, Menghini and QuarteroniDedè et al. 2021; Reference Zingaro, Dede, Menghini and QuarteroniZingaro et al. 2021) and cardiovascular flow simulations (Reference Hughes, Takizawa, Bazilevs, Tezduyar and HsuHughes et al. 2020; Reference Liu, Yang, Lan and MarsdenLiu et al. 2020).

Left ventricular haemodynamics is significantly affected by the functioning of the mitral valve, which ensures the correct direction of blood flow from the left atrium to the left ventricle. Mitral valve disease is the most common valvular heart disorder, particularly in ageing populations, with a prevalence of more than 10 % in people over 75 years (Reference Nkomo, Gardin, Skelton, Gottdiener, Scott and Enriquez-SaranoNkomo et al. 2006). The two main types of mitral valve disease are mitral valve stenosis and mitral valve regurgitation. Mitral valve stenosis occurs due to the narrowing or stiffening of the mitral valve, which prevents it from fully opening. As a result, flow from the left atrium to the left ventricle is obstructed and, in turn, the heart pumps less efficiently. Mitral valve stenosis causes a pressure increase in the left atrium, which induces atrial wall thickening and often atrial dilatation. The augmented atrial pressure may also produce blood congestion in the pulmonary veins (Reference Maeder, Weber, Buser, Gerhard, Haager, Maisano and RickliMaeder et al. 2018). From CFD simulations, the mitral jet, generated during diastole, is seen to shrink in size and strengthen when the stenosis gets more severe. As a consequence, the kinetic energy of the flow, the shear stresses acting on the tissues and the transvalvular pressure drop increase (Reference Meschini, Viola and VerziccoMeschini et al. 2019). This also increases mitral valve regurgitation, where the incomplete closure of the mitral valve allows blood to flow backward into the atrium. As mitral valve regurgitation worsens, the heart must work harder to pump blood to the body, leading to a progressive increase in left ventricle volume with normalization of wall stress. In terms of the fluid dynamics, large-scale flow recirculation, smaller vortices scattered throughout the left ventricle and the presence of a regurgitant jet which promotes increased turbulence in left ventricle have been observed in patients with mitral regurgitation (Reference Al-Wakeel, Fernandes, Amiri, Siniawski, Goubergrits, Berger and KuehneAl-Wakeel et al. 2015; Reference Bennati, Giambruno, Renzi, Di Nicola, Maffeis, Puppini, Luciani and VergaraBennati et al. 2023).

Common therapies for valvular diseases include valve repair or replacement with prosthetic heart valves. Present-day prosthetic heart valves can be classified into two broad categories, rigid-leaflet or mechanical valves, and valves with flexible leaflets made either of synthetic materials or biological tissue (the latter are called bioprosthetic valves). However, all known heart valve prostheses have been associated with various clinical complications due to the non-physiological flow patterns they induce upon implantation. Experimental techniques such as laser Doppler velocimetry and particle image velocimetry have been extensively employed to study prosthetic heart valve flows in vitro, showing that the transmitral vortex formation could be influenced by the saddle annulus dynamics and also by the angle of valve opening (Reference Kheradvar and FalahatpishehKheradvar & Falahatpisheh 2012). Fluid–solid interaction simulations have also revealed drastically different flow fields created by different prosthetic heart valves, and identified non-physiologic flow patterns causing blood cell damage and platelet activation (Reference BorazjaniBorazjani 2013). Moreover, studies have shown that, with impaired left ventricular function, the flow disturbances introduced by the bi-leaflet mechanical valve reduce the penetration capability of the mitral jet and weaken the recirculation in the ventricular apex (Reference Meschini, De Tullio, Querzoli and VerziccoMeschini et al. 2018). These effects are discussed in detailed prior reviews on heart valve fluid mechanics (Reference PeskinPeskin 1982; Reference Yoganathan, Chandran and SotiropoulosYoganathan, Chandran & Sotiropoulos 2005; Reference Sotiropoulos, Le and GilmanovSotiropoulos et al. 2016).

Another pathology of the left ventricle is cardiomyopathy (figure 2a), which refers to a condition where the heart muscle becomes enlarged, thickened or stiff, consequently diminishing the heart's ability to pump blood and potentially leading to heart failure. The most prevalent type is dilated cardiomyopathy, which is characterized by progressive enlargement and deterioration of the left ventricle's contraction capacity. The reduced cardiac function is reflected in a weaker diastolic mitral jet which penetrates less into the wider left ventricle, and a smaller vortex in the centre of the cavity that does not support flow redirection toward the outflow, resulting in increased residual volume during pumping (Reference Mangual, Kraigher-Krainer, De Luca, Toncelli, Shah, Solomon, Galanti, Domenichini and PedrizzettiMangual et al. 2013; Reference Eriksson, Bolger, Ebbers and CarlhällEriksson et al. 2016; Reference RivaRiva et al. 2022).

Figure 2.

(a) Average flow and coagulation patterns in the left ventricle for two patients with cardiomyopathy, from patient-specific simulations of flow coupled with biochemical coagulation reactions. Note the different flow fields and coagulation patterns, highlighting the importance of patient-specific risk. Courtesy of J.-H. Seo and R. Mittal, Johns Hopkins University. (b) Simulated flow through normal and bicuspid aortic valves. Note the asymmetric jets through bicuspid valves. The middle panel shows the more common phenotype (fusion of the left and right coronary cusp leaflets) and the right panel shows the less common phenotype (fusion of the right and non-coronary cusp leaflets) (Reference Kaiser, Shad, Schiavone, Hiesinger and MarsdenKaiser et al. 2022).

Left ventricular assist devices (LVADs) have become a key element for the treatment of patients suffering from terminal heart failure. These devices serve different purposes, such as a bridge to transplant (Reference MillerMiller et al. 2007; Reference AaronsonAaronson et al. 2012), a bridge to recovery or as long-term therapy (Reference SlaughterSlaughter et al. 2009). The LVADs pump blood from the left ventricle to the aortic root, providing a mechanical counterpart that works in conjunction with the patient's own left ventricle. Modern LVADs employ continuous flow physiology, which means they continuously pump blood throughout the cardiac cycle. The LVADs provide significant unloading of the left ventricle, manifested by decreased end-diastolic pressure, end-diastolic volume and pressure-volume area (Reference UrielUriel et al. 2016). Simultaneously, a decrease in the left ventricular volume is evident almost immediately in patients on LVAD support (Reference McCarthy, Nakatani, Vargo, Kottke-Marchant, Harasaki, James, Savage and ThomasMcCarthy et al. 1995). Besides the left ventricle, device implantation substantially modifies the haemodynamic conditions within the aorta. Due to the high velocity introduced by the pump, the aortic valve may remain persistently closed, especially when native cardiac ejection is low. Moreover, the adverse event profile of contemporary left ventricular assist devices includes suction events (where devices unload too much so that the ventricular chamber collapses) (Reference Vollkron, Voitl, Ta, Wieselthaler and SchimaVollkron et al. 2007; Reference Griffin and KatzGriffin & Katz 2014), thrombosis within the device (Reference Eckman and JohnEckman & John 2012), right ventricular failure, device-related infections, stroke, dysrhythmia and aortic incompetence (Reference ParkPark et al. 2012; Reference Kirklin, Naftel, Pagani, Kormos, Stevenson, Blume, Myers, Miller, Baldwin and YoungKirklin et al. 2015). Computational fluid dynamics has been employed to characterize and optimize numerous different designs for LVAD functional performance towards improved long-term outcomes (Reference Fraser, Taskin, Griffith and WuFraser et al. 2011). Computational fluid dynamics results identify that, with device support, the aortic root is largely a region of vascular stasis, with significantly higher shear stress on the leaflet tips in patients with aortic incompetence (Reference ShadShad et al. 2021). Strong correlations between device outflow graft orientation and wall shear stress, pressure and turbulent energy dissipation in the ascending aorta are found. Based on these flow patterns, a more distal insertion has been recommended to produce at least occasional parallel flow through the aortic valve and ensure intermittent aortic valve opening (Reference Osorio, Osorio, Ceballos, Tran, Clark, Divo, Argueta-Morales, Kassab and DeCampliOsorio et al. 2013; Reference Karmonik, Partovi, Loebe, Schmack, Weymann, Lumsden, Karck and RuhparwarKarmonik et al. 2014).

3. Aortic flow

Oxygenated blood pumped by the left ventricle to the rest of the circulation first passes through the aorta, the largest blood vessel in the human body at roughly 3 cm diameter. Continuing our Lagrangian journey, from the outlet of the left ventricle, the aorta ascends towards the thoracic region of the body, where it arches around and then descends into the abdomen. Along the way, several arteries branch off to carry blood to the heart, head, torso, kidneys and other organs. The aortic root features a valve which consists of three passive leaflets that are housed within three corresponding pouches in the aortic wall. These pouches are called aortic sinuses. The aortic valve, aortic sinuses and the curvature of the aortic arch interact with the flow emanating from the left ventricle to exhibit an array of interesting fluid mechanics that has fascinated scientists for centuries. This is famously evidenced by the illustrations of da Vinci in 1513 showing regions of recirculation around the aortic valve (Reference BoonBoon 2009).

Of particular interest has been the fluid mechanics of the opening and closure of the aortic valve and its relation to the geometry of the aortic sinus (Reference Yoganathan, He and JonesYoganathan, He & Jones 2004; Reference Sotiropoulos, Le and GilmanovSotiropoulos et al. 2016). The three passive leaflets of the aortic valve respond to the haemodynamics within the upstream left ventricle and the downstream aorta to provide minimal flow resistance while open during ventricular contraction at systole, and to minimize backflow during ventricular filling at diastole. The opening is driven by the pressure difference between the contracting left ventricle and the aorta, i.e. the leaflets open when the pressure within the left ventricle exceeds the aortic pressure. Interestingly, the closure of the aortic valve begins before the end of this ejection phase – with an estimated 75 % of the valve closure occurring while blood continues to flow from the left ventricle to the aorta (Reference Bellhouse and TalbotBellhouse & Talbot 1969). This ensures that the valve is completely closed with minimal reversed flow at the end of systole.

The hydrodynamic mechanism driving this has intrigued numerous researchers (Reference Bellhouse and BellhouseBellhouse & Bellhouse 1968; Reference Peskin and WolfePeskin & Wolfe 1978), starting with the first modern study of this dynamics by Reference Henderson and JohnsonHenderson & Johnson (1912). The seminal works of Bellhouse and colleagues highlighted the importance of the sinus vortex, which results from flow over the top of the aortic valves into the aortic sinuses, in initiating aortic valve closure using laboratory experiments accompanied by elegant theoretical calculations of the pressure induced by these vortices on the leaflet surfaces (Reference BellhouseBellhouse 1969; Reference Bellhouse and TalbotBellhouse & Talbot 1969). They also hypothesized the role of these vortices in maintaining the valve open position through a balance between sinus flow and the flow through the valve. While the physiological significance of the sinus vortex is still a topic of investigation (Reference Moore and DasiMoore & Dasi 2014, Reference Moore and Dasi2015), it is remarkable that the contemporary results exhibit a striking parallel with the illustrations of Da Vinci over 400 years earlier (Reference Gharib, Kremers, Koochesfahani and KempGharib et al. 2002).

The interaction of the flow with the valve leaflets and the aortic sinus is also important to the pathophysiology of aortic valve disease. Pathological conditions relating to aortic root diseases and congenital defects, such as bicuspid aortic valves where the valve forms with two leaflets instead of three, lead to abnormal flow patterns. In these cases, studies have reported asymmetric jets emanating from the valve opening, which impinge on the aortic wall to cause unphysiological shear stresses (Reference Girdauskas, Borger, Secknus, Girdauskas and KuntzeGirdauskas et al. 2011; Reference Barker, Markl, Bürk, Lorenz, Bock, Bauer, Schulz-Menger and Von Knobelsdorff-BrenkenhoffBarker et al. 2012; Reference Verma and SiuVerma & Siu 2014; Reference Saikrishnan, Mirabella and YoganathanSaikrishnan, Mirabella & Yoganathan 2015; Reference Kaiser, Shad, Schiavone, Hiesinger and MarsdenKaiser et al. 2022). This is shown in figure 2(b). The replacement of diseased aortic valves with mechanical and prosthetic valves also alters the haemodynamics. This is most commonly done using mechanical valves that consist of flat pivoting metal plates to mimic leaflets. The interaction of flow with the rigid leaflets gives rise to highly complex unphysiological flow within the ascending aorta, such as interacting vortex rings and shear layers emanating from the tips of the flat-plate leaflets (Reference Dasi, Simon, Sucosky and YoganathanDasi et al. 2009; Reference Sotiropoulos, Le and GilmanovSotiropoulos et al. 2016). Moreover, studies have reported increased shear regions within the housing of prosthetic valves that can damage blood cells (Reference Dasi, Simon, Sucosky and YoganathanDasi et al. 2009). Prosthetic valves with flexible leaflets also exhibit the unphysiological flow patterns described above, albeit to a lesser degree depending on the properties of the leaflet and housing. In addition, geometrical alterations to the aortic root through surgical procedures to treat aortic valve diseases have been shown to result in secondary flows that exhibit turbulent transition and associated stress fluctuations (Reference Stein and SabbahStein & Sabbah 1976; Reference PeacockPeacock 1990; Reference Moore and DasiMoore & Dasi 2014). The interaction of these unphysiological flow features with the valve leaflets is believed to exacerbate the stiffening or calcification of aortic valve leaflets as well as complications at the aortic wall. While the exact mechanism leading to leaflet calcification is not clear, there are strong correlations between flow-induced mechanical stresses, such as oscillatory shear stresses from flow recirculation, and the inflammation of endothelial cells on the leaflets (Reference Ge and SotiropoulosGe & Sotiropoulos 2010; Reference Gould, Srigunapalan, Simmons and AnsethGould et al. 2013). The flow through asymmetrically calcified or stenosed valves has been shown to develop into a jet that impinges on the aortic wall, where the curvature of the aortic arch introduced helical Dean vortex flows, intermittent turbulent stresses and audible sounds called murmurs (Reference Zhu, Seo and MittalZhu, Seo & Mittal 2018, Reference Zhu, Seo and Mittal2019; Reference Bailoor, Seo, Dasi, Schena and MittalBailoor et al. 2021).

The rich haemodynamics associated with aortic flows has led to efforts to non-invasively quantify these flow in vivo using 4-D flow MRI in search of haemodynamic biomarkers of disease (Reference Von Knobelsdorff-Brenkenhoff, Trauzeddel, Barker, Gruettner, Markl and Schulz-MengerVon Knobelsdorff-Brenkenhoff et al. 2014; Reference Burris and HopeBurris & Hope 2015; Reference Garcia, Barker and MarklGarcia et al. 2019; Reference Zimmermann, Bäumler, Loecher, Cork, Marsden, Ennis and FleischmannZimmermann et al. 2023). In this context, 4-D flow MRI is a promising method to non-invasively evaluate haemodynamic patterns that can aid in the diagnosis of valvular diseases as well as in the evaluation of prosthetics. However, the challenge with 4-D flow MRI is accurately resolving the flow dynamics and associated wall interactions, which has led to the emergence of computational methods to simulate these flows. The first method aimed at simulating the fluid dynamics of cardiac valves was introduced by Peskin and colleagues just over 50 years ago (Reference PeskinPeskin 1972; Reference Peskin and McQueenPeskin & McQueen 1980). This technique, which is now referred to as the immersed boundary method, continues to be the mainstay of computational models for cardiac valves (Reference Mittal and IaccarinoMittal & Iaccarino 2005; Reference Griffith and PatankarGriffith & Patankar 2020). Since its introduction, the immersed boundary method has been coupled with finite-element formulations of the mechanics and material properties of valve leaflets (Reference Boffi, Gastaldi, Heltai and PeskinBoffi et al. 2008; Reference Griffith and LuoGriffith & Luo 2017; Reference Kaiser, Shad, Hiesinger and MarsdenKaiser et al. 2021). Sharp-interface immersed boundary methods have been developed, which can resolve the surface of the leaflets more accurately and therefore provide improved surface force estimates (Reference Gilmanov and SotiropoulosGilmanov & Sotiropoulos 2016; Reference Seo, Zhu, Resar and MittalSeo et al. 2020; Reference Bailoor, Seo, Dasi, Schena and MittalBailoor et al. 2021). In the finite-element context, these immersed boundary fluid–structure interaction methods are broadly referred to as immersed finite-element or immersogeometric methods (Reference Hsu, Kamensky, Bazilevs, Sacks and HughesHsu et al. 2014, Reference Hsu, Kamensky, Xu, Kiendl, Wang, Wu, Mineroff, Reali, Bazilevs and Sacks2015; Reference Kamensky, Hsu, Schillinger, Evans, Aggarwal, Bazilevs, Sacks and HughesKamensky et al. 2015).

These methods have been extensively used to study the haemodynamics of native as well as prosthetic valves in healthy and diseased states (see Reference Votta, Le, Stevanella, Fusini, Caiani, Redaelli and SotiropoulosVotta et al. (2013) and Reference Sun, Mao and GriffithSun, Mao & Griffith (2019) for comprehensive reviews of applications). In fact, one such simulation-based product has recently received United States Food and Drug Administration approval to support clinical decision making for transcatheter aortic valve replacements (www.dasisim.com). However, an ongoing challenge for more widespread clinical acceptance of these simulation tools is the lack of validation with in vivo flow fields. The majority of validation studies have focused on the haemodynamics of mechanical heart valves measured in laboratory experiments (Reference Dasi, Ge, Simon, Sotiropoulos and YoganathanDasi et al. 2007; Reference Ge, Dasi, Sotiropoulos and YoganathanGe et al. 2008; Reference Nobili, Morbiducci, Ponzini, Del Gaudio, Balducci, Grigioni, Maria Montevecchi and RedaelliNobili et al. 2008; Reference Guivier-Curien, Deplano and BertrandGuivier-Curien, Deplano & Bertrand 2009; Reference Jun, Saikrishnan, Arjunon, Yun and YoganathanJun et al. 2014). In comparison, fewer studies have validated simulation models of heart valves with flexible leaflets (Reference WuWu et al. 2016; Reference Sigüenza, Pott, Mendez, Sonntag, Kaufmann, Steinseifer and NicoudSigüenza et al. 2018; Reference Lee, Rygg, Kolahdouz, Rossi, Retta, Duraiswamy, Scotten, Craven and GriffithLee et al. 2020), especially with patient-specific flow fields (Reference Emendi, Sturla, Ghosh, Bianchi, Piatti, Pluchinotta, Giese, Lombardi, Redaelli and BluesteinEmendi et al. 2021; Reference Kaiser, Schiavone, Elkins, McElhinney, Eaton and MarsdenKaiser et al. 2023).

While a major focus of studies on aortic haemodynamics has been on the function of the aortic valve, the influence of aortic flow on the wall of the aorta has important health consequences. In particular, shear stresses induced on the aortic wall due to the jet emanating from the aortic valve as well as the resultant turbulent flow within the aorta can damage the epithelial cell lining of the aortic wall. The peak Reynolds number in the aorta can be as high as ${\sim }4000$, which is significantly higher than the transitional Reynolds number in pipes (Reference KuKu 1997). However, interesting work has shown that the nature of cardiac pulsatility, including the specific waveform of aortic flow, is close to the optimal behaviour required to combine sufficient flow rates with low wall shear stresses. This is due to the combination of a short ejection phase followed be a prolonged diastolic phase that inhibits turbulence (Reference Scarselli, Lopez, Varshney and HofScarselli et al. 2023).

Studies have also reported that 60 %–80 % of patients with bicuspid aortic valves develop dilation of the aortic wall (Reference Verma and SiuVerma & Siu 2014). While the exact mechanism leading to aortic dilation is unknown, altered haemodynamics associated with the impinging jet arising from bicuspid valves is thought to play a role (Reference Girdauskas, Borger, Secknus, Girdauskas and KuntzeGirdauskas et al. 2011; Reference Barker, Markl, Bürk, Lorenz, Bock, Bauer, Schulz-Menger and Von Knobelsdorff-BrenkenhoffBarker et al. 2012; Reference Rodríguez-PalomaresRodríguez-Palomares et al. 2018; Reference Dux-SantoyDux-Santoy et al. 2020). Furthermore, aortic dilation increases the risk of aortic dissection, which is a tear in the inner wall of the aorta that introduces a parallel flow channel within the aortic wall, referred to as the false lumen (in contrast to the native ‘true’ lumen), and poses the risk of rupture of the aortic wall (Reference Verma and SiuVerma & Siu 2014). Although the clinical management of aortic dissections is primarily based on anatomical imaging (Reference Spinelli, Benedetto, Donato, Piffaretti, Marrocco-Trischitta, Patel, Eagle and TrimarchiSpinelli et al. 2018), studies have shown that haemodynamics plays a significant role in the prognosis (Reference Bäumler, Vedula, Sailer, Seo, Chiu, Mistelbauer, Chan, Fischbein, Marsden and FleischmannBäumler et al. 2020). For example, while the development of clots to occlude the false lumen has been shown to improve patient outcomes, the structure of these clots can cause elevated pressure within the false lumen and increase rupture risk (Reference Tsai, Schlicht, Khanafer, Bull, Valassis, Williams, Berguer and EagleTsai et al. 2008; Reference Berguer, Parodi, Schlicht and KhanaferBerguer et al. 2015).

Another pathological condition that poses the risk of aortic wall rupture is the occurrence of dilations in the abdominal aorta, referred to as abdominal aortic aneurysms (Reference Sakalihasan, Limet and DefaweSakalihasan, Limet & Defawe 2005). For patients with aortic dissections and aortic aneurysms, the stratification of rupture risk depends on thrombus formation due to flow stagnation as well as stress distributions induced on the aortic wall through its interaction with the haemodynamics (Reference Bluestein, Niu, Schoephoerster and DewanjeeBluestein et al. 1997; Reference Vorp, Lee, Wang, Makaroun, Nemoto, Ogawa and WebsterVorp et al. 2001; Reference FrancisFrancis 2002; Reference Sho, Sho, Hoshina, Kimura, Nakahashi and DalmanSho et al. 2004). A challenge in this context is the quantification of aortic wall mechanics, which is often heterogeneous and difficult to measure in vivo. Conversely, there are also haemodynamic risks associated with the constriction of the aortic wall. This narrowing of the descending aorta, which is referred to as aortic coarctation and accounts for 6 %–8 % of congenital heart defects, can lead to aortic rupture, hypertension and stroke (Reference Singh, Hakim, Sharma, Roy, Panse, Chandrasekaran, Alegria and MookadamSingh et al. 2015). These are therefore active avenues for the future application of patient-specific computational models to aid in clinical decision making (for example, Reference Bluestein, Niu, Schoephoerster and DewanjeeBluestein et al. 1997; Reference Salsac, Sparks and LasherasSalsac, Sparks & Lasheras 2004; Reference Les, Shadden, Figueroa, Park, Tedesco, Herfkens, Dalman and TaylorLes et al. 2010; Reference Alimohammadi, Sherwood, Karimpour, Agu, Balabani and Díaz-ZuccariniAlimohammadi et al. 2015; Reference Chen, Peelukhana, Berwick, Kratzberg, Krieger, Roeder, Chambers and KassabChen et al. 2016; Reference Qiao, Zeng, Ding, Fan, Luo and ZhuQiao et al. 2019; Reference Bäumler, Vedula, Sailer, Seo, Chiu, Mistelbauer, Chan, Fischbein, Marsden and FleischmannBäumler et al. 2020, among many others). Reference Mourato, Valente, Xavier, Brito, Avril, de Sá, Tomás and FragataMourato et al. (2022) provide a comprehensive review of fluid–structure interaction computational models pertinent to aortic wall diseases.

4. Coronary haemodynamics

While the pumping action of the heart ensures the delivery of blood to the far reaches of our circulatory system, the coronary arteries bear the unenviable responsibility of delivering nutrition and oxygen-laden blood to the heart muscle itself. The importance of the coronary circulation is underscored by the fact that the heart has the highest oxygen consumption per gram of any organ; it extracts 70 %–80 % of the oxygen delivered, compared with 30 %–40 % utilized by skeletal muscles (Reference Goodwill, Dick, Kiel and TuneGoodwill et al. 2017). It is no wonder that coronary artery disease, which inhibits the transport of blood to the myocardium due to the narrowing or occlusion of coronary arteries, is the leading cause of death worldwide (Reference Ralapanawa and SivakanesanRalapanawa & Sivakanesan 2021; Reference TsaoTsao et al. 2022).

The left and right coronary arteries arise on either side of the aortic sinus and carry 3 %–5 % of the volume of blood pumped by the heart into the myocardial tissue (Reference Guyton and HallGuyton & Hall 2006). They have diameters of ${\approx }3$ mm at their origin on the aortic sinus, and proceed through a series of bifurcations that leads to small arteries and arterioles with diameters ${\approx }100\ \mathrm {\mu }$m, and eventually capillaries with diameters ${\approx }10\ \mathrm {\mu }$m (Reference Lee and SmithLee & Smith 2012). Although the peak Reynolds number in coronary arteries is typically ${<}500$, coronary flow can exhibit higher Reynolds number jets downstream of constrictions that occur due to disease. This has consequences for the progression of coronary artery disease as well as the stability of calcified plaques that cause these constrictions (Reference Falk, Shah and FusterFalk, Shah & Fuster 1995; Reference Samady, Eshtehardi, McDaniel, Suo, Dhawan, Maynard, Timmins, Quyyumi and GiddensSamady et al. 2011). A peculiar aspect of the coronary circulation is that, while maximum flow through most vascular networks in the body occurs during systole, a majority of the flow through coronary arteries occurs during diastole. This is due to the so-called ‘cross-talk’ between the contracting myocardial tissue and the coronary circulation from large vessels to arterioles and capillaries embedded within the myocardial tissue (Reference Westerhof, Boer, Lamberts and SipkemaWesterhof et al. 2006). The disparate length scales involved in the coronary circulation and its interaction with the beating (and moving) heart, combined with the importance of coronary flows to cardiovascular function, present a host of challenges in understanding, imaging, as well as modelling coronary artery blood flows.

Most studies on coronary artery haemodynamics have focused on the flow through the large arteries that can be clinically imaged, either invasively using coronary angiography or non-invasively using computed tomography. Efforts to measure and model coronary artery flow are motivated by several randomized clinical trials that have firmly established the importance of considering haemodynamics in clinical care. Coronary artery disease causes the narrowing and stiffening of coronary artery walls, which limits the delivery of blood to the myocardium. This has been observed to occur in regions of the vasculature that experience low and oscillating wall shear stress in laboratory experiments as well as animal models (Reference Giddens, Zarins and GlagovGiddens et al. 1993). Coronary artery disease is usually treated either by inserting stents within affected arteries to relieve the narrowing, or in cases with more severe occlusion, surgically bypassing the location of narrowing in arteries by using grafts (Reference YusufYusuf et al. 1994; Reference WindeckerWindecker et al. 2014). In both cases, the long-term outcomes for patients receiving these treatments has been shown to significantly improve when treatments are guided by haemodynamic measurements of pressure loss within the affected artery, rather than based on anatomical measures of the narrowing in the artery (Reference MeijboomMeijboom et al. 2008; Reference ToninoTonino et al. 2009; Reference PijlsPijls et al. 2010). There is also evidence that short as well as long-term treatment outcomes are dependent on local blood flow patterns at the site of stenting or grafting, such as flow reversal and competitive flow between the native diseased artery and bypass graft (Reference Pagni, Storey, Ballen, Montgomery, Chiang, Etoch and SpencePagni et al. 1997; Reference Sabik, Lytle, Blackstone, Khan, Houghtaling, Cosgrove, Shumway, Tatoulis and BednarskiSabik et al. 2003; Reference Tokuda, Song, Oshima, Usui and UedaTokuda et al. 2008), as well as regions of stagnation or abnormal wall shear stress (Reference Malek, Alper and IzumoMalek, Alper & Izumo 1999; Reference ChatzizisisChatzizisis et al. 2008; Reference HumphreyHumphrey 2008).

While the major focus of studies on coronary artery blood flows has been on the flows within the coronary arteries themselves, the location of coronary arteries and their proximity to the aortic valve also raises interesting questions. Owing to their origin within the aortic sinus, researchers have long hypothesized about the influence of recirculating flow within the aortic sinus (due to interaction with the aortic valve) on coronary circulation, and vice versa. Reference Bellhouse, Bellhouse and ReidBellhouse, Bellhouse & Reid (1968) used a laboratory model to relate flow stagnation on the ridge of the aortic sinus to high pressure that would be beneficial to driving coronary flow. They showed a significant drop in coronary flow when the coronary arteries originate outside the aortic sinus. There is also a difference in the fluid dynamics between the two sinus lobes that contain the origin of coronary arteries, or coronary sinuses, and the lobe that does not. Reference Moore and DasiMoore & Dasi (2015) showed that the presence of coronary arteries originating from the coronary sinuses pulls the sinus vortex closer to the walls. This leads to higher flow velocities, shear stresses and washout in the coronary sinuses compared with the non-coronary sinus. In addition, the aortic leaflets corresponding to coronary sinuses open farther into the sinus and exhibit better mechanics. These haemodynamic differences have physiological significance because low recirculation velocities lead to less washout in the non-coronary sinus, and low and oscillatory shear stresses have been linked to calcification (Reference Chandra, Rajamannan and SucoskyChandra, Rajamannan & Sucosky 2012). Interestingly, these observations also present a possible mechanism for the observation that the non-coronary sinus is often the first to calcify – eventually leading to aortic valve calcification (Reference Freeman and OttoFreeman & Otto 2005).

The importance of fluid mechanics in these critical clinical questions has spurred development of patient-specific computational models of coronary artery blood flows for personalized risk assessment and treatment planning. Current state-of-the-art methods use non-invasive clinical imaging to build patient-specific coronary artery anatomical models (Reference Updegrove, Wilson, Merkow, Lan, Marsden and ShaddenUpdegrove et al. 2017). The flow within these anatomies is simulated using either high-fidelity 3-D flow simulation methods, or reduced-order methods that evaluate the flow field either only along the centreline of the arteries or at specifically chosen points in the anatomy, such as bifurcations, by assuming Poiseuille flow throughout (see § 7). Owing to the complex geometries of vascular trees, 3-D simulations are usually performed using stabilized finite-element methods with unstructured meshes (Reference Taylor, Hughes and ZarinsTaylor et al. 1998; Reference Marsden and Esmaily-MoghadamMarsden & Esmaily-Moghadam 2015). A significant challenge is that these image-based models are usually restricted to the large vessels in the coronary tree that are visible from volumetric clinical imaging. Therefore, the effect of the smaller downstream arteries and their interaction with the contracting myocardium needs to be modelled appropriately. This is usually accomplished by coupling the flow within the large arteries included in these anatomical models with lumped-parameter boundary conditions at their outlets (Reference Kim, Vignon-Clementel, Coogan, Figueroa, Jansen and TaylorKim et al. 2010). These boundary conditions are designed to account for the resistance as well as dilation of downstream blood vessels not included in the model. They also recapitulate the key features of coronary physiology in which coronary flows peak at diastole due to the ‘intramyocardial pressure’ experienced by the coronary circulation via its interaction with myocardial contractions. In addition, the compliance of the coronary artery walls as well as their movement with the contracting myocardium pose additional challenges in the modelling of coronary haemodynamics. Although the former is commonly included in models using simplified as well as high-fidelity fluid–structure interaction techniques, such as the coupled momentum method (Reference Figueroa, Vignon-Clementel, Jansen, Hughes and TaylorFigueroa et al. 2006) and arbitrary Lagrangian–Eulerian (Reference Seo, Schiavazzi, Kahn and MarsdenSeo et al. 2020), studies have shown that the effect of wall compliance can be small in certain situations (Reference ToriiTorii et al. 2009; Reference Eslami, Tran, Jin, Karady, Sotoodeh, Lu, Hoffmann and MarsdenEslami et al. 2020). However, although the movement of coronary arteries due to myocardial contraction can have significant effects on the haemodynamics, this is rarely included in simulations (Reference Santamarina, Weydahl, Siegel and MooreSantamarina et al. 1998).

As expected, the modelling techniques and boundary conditions used to simulate coronary haemodynamics include several parameters relating to the coupling between the microcirculation and large vessels that cannot be measured on a patient-specific basis. Importantly, they determine how the flow entering the coronary artery tree is split amongst all the vessels in the tree. While there are physiological and empirical measures relating coronary artery diameters to the flow splits amongst the arteries (Reference MurrayMurray 1926; Reference Zhou, Kassab and MolloiZhou, Kassab & Molloi 1999), these are not patient specific and do not account for vascular disease. It is not surprising that parameters affecting the amount of flow through each coronary artery have a measurable impact on our ability to predict pressure losses resulting from narrowing or aneurysms in these arteries (Reference Wellnhofer, Osman, Kertzscher, Affeld, Fleck and GoubergritsWellnhofer et al. 2010; Reference Sankaran, Kim, Choi and TaylorSankaran et al. 2016; Reference Vardhan, Gounley, Chen, Kahn, Leopold and RandlesVardhan et al. 2019). Advances in imaging methods, combined with patient-specific modelling, are now aiming to directly measure myocardial perfusion to address this need (Reference Nieman and BallaNieman & Balla 2020; Reference Menon, Khan, Sexton, Richter, Nieman and MarsdenMenon et al. 2024).

In spite of the challenges outlined above, computational modelling of coronary artery haemodynamics is a powerful tool for non-invasive clinical risk assessment and treatment planning. The computation of the pressure loss resulting from diseased coronary arteries using patient-specific computational models has successfully predicted disease severity and reduced unnecessary invasive procedures in clinical trials (Reference KooKoo et al. 2011; Reference MinMin et al. 2012; Reference Taylor, Fonte and MinTaylor et al. 2013; Reference Mathew, Gottbrecht and SalernoMathew, Gottbrecht & Salerno 2018). This has led to the approval of the first simulation-based risk assessment tool, $FFR_{{CT}}$, for clinical use by the United States Food and Drug Administration (www.heartflow.com). In addition, computational modelling also provides insights not accessible in clinical imaging. For example, the multi-directionality and topological skeleton of wall shear stress, i.e. regions of expansion and constriction in the wall shear stress field, have emerged as possible predictor of atherosclerosis (Reference Arzani and ShaddenArzani & Shadden 2018; Reference Kok, Molony, Timmins, Ko, Boersma, Eshtehardi, Wentzel and SamadyKok et al. 2019; Reference HoogendoornHoogendoorn et al. 2021; Reference Mazzi, De Nisco, Hoogendoorn, Calò, Chiastra, Gallo, Steinman, Wentzel and MorbiducciMazzi et al. 2021). In addition, computational modelling has highlighted the role of coronary flow helicity in promoting healthy wall shear stress (Reference Morbiducci, Ponzini, Grigioni and RedaelliMorbiducci et al. 2007; Reference De Nisco, Kok, Chiastra, Gallo, Hoogendoorn, Migliavacca, Wentzel and MorbiducciDe Nisco et al. 2019), as well as other insights into wall shear stress patterns that are related to the progression of coronary artery disease and coronary artery bypass graft failure (Reference Nordgaard, Swillens, Nordhaug, Kirkeby-Garstad, Van Loo, Vitale, Segers, Haaverstad and LovstakkenNordgaard et al. 2010; Reference Ramachandra, Kahn and MarsdenRamachandra, Kahn & Marsden 2016; Reference Khan, Tran, Zhu, Boyd, Packard, Karlsberg, Kahn and MarsdenKhan et al. 2020; Reference CandrevaCandreva et al. 2022).

The importance of haemodynamics in risk assessment and treatment has also been demonstrated in other afflictions of the coronary arteries. For example, Kawasaki disease is a paediatric inflammatory condition that leads to coronary artery aneurysms in some untreated patients (Reference McCrindleMcCrindle et al. 2017). A risk associated with these aneurysms is clot formation as a result of flow stagnation and recirculation. While current clinical guidelines recommend treatment based purely on aneurysm size, recent studies have demonstrated the use of computational modelling to reveal blood flow patterns within aneurysms that might provide a better indicator of patients’ risk (Reference Sengupta, Kahn, Burns, Sankaran, Shadden and MarsdenSengupta et al. 2012, Reference Sengupta, Kahn, Kung, Esmaily Moghadam, Shirinsky, Lyskina, Burns and Marsden2014; Reference Grande Gutierrez, Mathew, McCrindle, Tran, Kahn, Burns and MarsdenGrande Gutierrez et al. 2019; Reference Menon, Seo, Fukazawa, Ogawa, Kahn, Burns and MarsdenMenon et al. 2023). The use of computational modelling for non-invasive clinical risk assessment has also been explored in patients affected by anomalous aortic origin of coronary arteries, a congenital defect that causes the coronaries to arise from outside their normal aortic sinus location (Reference Jiang, Khan, Ghobrial, Rogers, Pettersson, Blackstone and MarsdenJiang et al. 2022). In addition to clinical risk stratification, there have also been efforts to employ patient-specific computational models of coronary haemodynamics for virtual planning of treatments and surgical procedures (Reference Balossino, Gervaso, Migliavacca and DubiniBalossino et al. 2008; Reference Sankaran, Moghadam, Kahn, Tseng, Guccione and MarsdenSankaran et al. 2012; Reference Morlacchi and MigliavaccaMorlacchi & Migliavacca 2013; Reference Ballarin, Faggiano, Manzoni, Quarteroni, Rozza, Ippolito, Antona and ScrofaniBallarin et al. 2017; Reference Seo, Ramachandra, Boyd, Marsden and KahnSeo et al. 2021). Due to the importance of haemodynamics in the clinical diagnoses, treatment planning and disease progression of such a central part of our circulatory system, these efforts to develop personalized and predictive models of coronary artery haemodynamics, as well as developments in imaging methods to measure coronary flows, have the potential to make a significant impact on cardiovascular health.

5. Cerebral haemodynamics

Following the branching off of coronary arteries at the aortic sinuses, the next set of arteries to arise from the aorta is located downstream at the aortic arch. Of these, the most important are the left and right common carotid arteries, the latter bifurcating from the brachiocephalic artery. The carotid arteries are responsible for supplying the face, neck and, crucially, the brain. The left and right internal carotid arteries, which are branches of the left and right common carotid arteries, travel up the neck into the base of the skull, after which they bifurcate into a network of cerebral arteries that perfuse the brain. In humans, the brain accounts for only 2 %–3 % of body mass, but receives approximately 15 %–20 % of cardiac output and consumes 20 % of the available oxygen in the body under normal conditions (Reference Claassen, Thijssen, Panerai and FaraciClaassen et al. 2021). The brain's high metabolic rate and the limited ability for anaerobic metabolism in neurons underscore the importance of maintaining healthy cerebral blood flow. The disruption of cerebral blood flow therefore has severe consequences that can lead to life-threatening strokes, permanent brain damage, as well as other maladies related to impaired removal of metabolic wastes from the brain (Reference Tarasoff-ConwayTarasoff-Conway et al. 2015).

One of the most widely investigated topics relating to cerebral circulation, is the growth, rupture and treatment of cerebral aneurysms. These aneurysms typically form in the Circle of Willis, which is a traffic circle-like junction of several arteries at the base of the brain. Research has shown that the most common locations of aneurysm development coincide with regions of elevated flow and wall shear stress, usually occurring due to pathological vessel narrowing or at bifurcations between arteries (Reference Sforza, Putman and CebralSforza et al. 2009). Haemodynamics also plays a role in the growth and rupture of cerebral aneurysms due to the interaction between flow-induced stresses and cellular function on blood vessel walls. In fact, computational modelling has revealed that a significant proportion of cerebral aneurysms exhibit flow instabilities that lead to high-frequency oscillations even at modest Reynolds numbers of ${\sim }300$ (Reference Khan, Toro Arana, Najafi, MacDonald, Natarajan, Valen-Sendstad and SteinmanKhan et al. 2021). Interestingly, both low and high wall shear stress have been implicated in unfavourable aneurysm progression, with each pathological state being associated with a different vascular response (Reference Shojima, Oshima, Takagi, Torii, Hayakawa, Katada, Morita and KirinoShojima et al. 2004; Reference Jou, Lee, Morsi and MawadJou et al. 2008; Reference Cebral, Mut, Weir and PutmanCebral et al. 2011; Reference Xiang, Natarajan, Tremmel, Ma, Mocco, Hopkins, Siddiqui, Levy and MengXiang et al. 2011). While there is general agreement on the mechanism of aneurysm initiation – triggered by cell damage caused by high flow and wall shear stress states that have been observed in animal models (Reference StehbensStehbens 1989) – there is disagreement on the mechanisms associating abnormal wall shear stress with aneurysm progression and rupture (Reference Xiang, Tutino, Snyder and MengXiang et al. 2014). Nevertheless, these observed links between haemodynamics and aneurysm progression, and the inability of in vivo imaging to measure critical haemodynamic stimuli, have accelerated the use of computational models to identify haemodynamic predictors of aneurysm progression (Reference Chung and CebralChung & Cebral 2015). Studies have introduced several haemodynamic metrics to quantify aneurysm rupture risk, such as high oscillatory shear index and particle residence time (Reference Gonzalez, Cho, Ortega and MoretGonzalez et al. 1992; Reference Xiang, Natarajan, Tremmel, Ma, Mocco, Hopkins, Siddiqui, Levy and MengXiang et al. 2011). The risk of aneurysm rupture has also been correlated with the presence of complex vorticity patterns (Reference Xiang, Natarajan, Tremmel, Ma, Mocco, Hopkins, Siddiqui, Levy and MengXiang et al. 2011) as well as the size, strength and location of the jet that emanates from the neck of the aneurysm and impinges on the aneurysm wall (Reference Cebral, Castro, Burgess, Pergolizzi, Sheridan and PutmanCebral et al. 2005).

Haemodynamics has also played a central role in the treatment of cerebral aneurysms to prevent rupture. Increasingly, clinical interventions have moved away from surgery, towards the implantation of endovascular devices that promote the formation of stabilizing clots within aneurysms (Reference Brisman, Song and NewellBrisman, Song & Newell 2006; Reference Brinjikji, Rabinstein, Nasr, Lanzino, Kallmes and CloftBrinjikji et al. 2011). The development of such clots, either forming spontaneously or promoted by implanted devices, has been shown to reduce the incidence of rupture (Reference Etminan and RinkelEtminan & Rinkel 2016). The most commonly implanted devices are coils and flow diverters, both of which impede flow into the aneurysm and lead to flow stagnation, low wall shear stress and, consequently, clot formation (Reference Ngoepe, Frangi, Byrne and VentikosNgoepe et al. 2018). Although the biochemistry of clot formation is currently still under investigation, computational haemodynamics models have proven to be remarkably effective at predicting the potential for clot formation and its location based on flow stagnation and residence time (Reference Sadasivan, Lieber, Gounis, Lopes and HopkinsSadasivan et al. 2002; Reference Rayz, Boussel, Ge, Leach, Martin, Lawton, McCulloch and SalonerRayz et al. 2010; Reference De SousaDe Sousa et al. 2016), as well as simulating the process of clot formation by coupling fluid and chemical transport with models of thrombogenesis (Reference Bedekar, Pant, Ventikos and SundaramBedekar et al. 2005; Reference Ngoepe and VentikosNgoepe & Ventikos 2016).

However, broader clinical acceptance of computational models for patient-specific risk stratification and thrombosis prediction has been slow. One reason for this is the lack of clear mechanobiology that relates the observed haemodynamic patterns to pathophysiology, which makes conflicting findings regarding high and low wall shear stress difficult to interpret (Reference Xiang, Tutino, Snyder and MengXiang et al. 2014). In the prediction of thrombosis, an interesting challenge is the wide range of time scales involved – the fluid transport occurs at the scale of seconds, clot formation occurs in minutes, while clot maturation and stabilization occurs over months (Reference Szikora, Turányi and MarosfoiSzikora, Turányi & Marosfoi 2015). Moreover, these interacting time scales influence the properties of blood itself, which changes during the process of clotting. Cerebral haemodynamics models also contend with the familiar challenges of patient specificity in terms of model parameters and boundary conditions (Reference Levitt, McGah, Aliseda, Mourad, Nerva, Vaidya, Morton, Ghodke and KimLevitt et al. 2014; Reference Steinman and PereiraSteinman & Pereira 2019). An added difficulty is that the biochemical coagulation profile of each person is different, and this needs to be accounted for in patient-specific simulations of clot formation (Reference Ngoepe, Frangi, Byrne and VentikosNgoepe et al. 2018). Finally, the lack of many large patient cohort studies has hampered their translation to clinical use (Reference Sforza, Putman and CebralSforza et al. 2009).

A unique aspect of the cerebral circulation that also benefits from patient-specific haemodynamic measurements and simulations is the presence of collateral arteries due to the anatomy of the Circle of Willis. Collaterals provide alternate flow paths that allow cerebral circulation to be robust to occlusions in specific parts of the vascular anatomy. However, they also complicate the clinical management of cerebral stenoses, because anatomical measurements of vessel narrowing do not correlate well with their effect on cerebral blood flow due to the presence of alternate flow paths via collaterals (Reference PowersPowers 1991). Furthermore, cerebral circulation exhibits a high degree of autoregulation, i.e. the blood vessels adapt to changes in upstream pressure and flow to preserve homeostatic cerebral blood flow volumes (Reference Claassen, Thijssen, Panerai and FaraciClaassen et al. 2021). The study of the haemodynamics associated with cerebral and carotid stenosis and their treatments at the patient-specific level has therefore emerged as another area for the application of personalized imaging and computational simulations based on imaging of patient-specific cerebral artery anatomies (Reference Birchall, Zaman, Hacker, Davies and MendelowBirchall et al. 2006; Reference Liang, Fukasaku, Liu and TakagiLiang et al. 2011).

Of particular interest is the haemodynamics of the growth and rupture of atherosclerotic plaque in carotid stenosis due to the risk of dislodged plaque travelling downstream, occluding other parts of the circulation and causing stroke (Reference Carr, Nemoto, Schwartz and ShaddenCarr et al. 2013; Reference Mukherjee, Jani, Selvaganesan, Weng and ShaddenMukherjee et al. 2016). The location of these plaques has been shown to most often correlate with regions of low wall shear stress and flow disruption (Reference Zarins, Giddens, Bharadvaj, Sottiurai, Mabon and GlagovZarins et al. 1983; Reference Ku, Giddens, Zarins and GlagovKu et al. 1985). The development of carotid re-stenosis following treatment has also been correlated with wall shear stress patterns that exhibit high temporal variation and strong expansion/contraction through the analysis of the wall shear stress topological skeleton (Reference Morbiducci, Mazzi, Domanin, De Nisco, Vergara, Steinman and GalloMorbiducci et al. 2020). In addition, the narrowing of these arteries due to the presence of plaque leads to accelerated and highly transient downstream flows that exhibit elevated shear stress and stress gradients as well as turbulent fluctuations, all of which are thought to cause cell damage and plaque rupture (Reference Davies, Remuzzi, Gordon, Dewey and GimbroneDavies et al. 1986; Reference Gertz and RobertsGertz & Roberts 1990; Reference Loree, Kamm, Atkinson and LeeLoree et al. 1991; Reference DePaola, Gimbrone, Davies and DeweyDePaola et al. 1992; Reference Lee, Lee, Fischer, Bassiouny and LothLee et al. 2008).

Finally, a fascinating and relatively unexplored haemodynamic phenomenon associated with cerebral circulation is the transport and clearance of brain interstitial fluid, which occupies interstitial regions of the brain surrounding blood vessels and neurons and carries metabolic wastes and toxins (Reference Shetty and ZaniratiShetty & Zanirati 2020). In particular, the impaired clearance of amyloid betas is a hallmark of Alzheimer's disease (Reference Tarasoff-ConwayTarasoff-Conway et al. 2015). While the bulk flow of interstitial fluid has recently been identified as one of several simultaneous mechanisms for the clearance of such wastes, the mechanisms driving this flow are currently unknown. Due to the fact that interstitial fluid envelopes cerebral vasculature, studies have suggested that arterial wave propagation as a result of pulsatile blood flow are a major driver of interstitial fluid transport (Reference Iliff, Wang, Zeppenfeld, Venkataraman, Plog, Liao, Deane and NedergaardIliff et al. 2013; Reference Ladrón-de Guevara, Shang, Nedergaard and KelleyLadrón-de Guevara et al. 2022). However, this has been debated due to the strength of these pulsations, with larger amplitude and lower-frequency motions caused by vasomotion proposed as an alternate mechanism (Reference Diem, Sharp, Gatherer, Bressloff, Carare and RichardsonDiem et al. 2017; Reference van Veluw, Hou, Calvo-Rodriguez, Arbel-Ornath, Snyder, Frosch, Greenberg and Bacskaivan Veluw et al. 2020). In spite of these disagreements, the role of circulatory phenomena in the clearance of toxins, as well as the aforementioned haemodynamic phenomena associated with aneurysms and strokes, have important consequences for brain health and are promising avenues for ongoing and future research.

6. Peripheral vascular flows

As we travel beyond the main arteries discussed so far, the circulatory system begins to branch off into the peripheral vascular system, which consists of the blood vessels outside of the heart and brain. The peripheral vascular system plays a critical role in transporting blood to and from the body's tissues, and it can be affected by various conditions. The most notable is peripheral artery disease, which affects over 8 million individuals aged 40 and older in the United States (Reference Aday and MatsushitaAday & Matsushita 2021). This disease is a slow and progressive disorder characterized by narrowing, blockage or spasms in the arteries, commonly caused by atherosclerosis – a condition marked by plaque accumulation on the arterial walls. This buildup of plaque narrows the arteries, impeding blood flow and often leading to complete occlusion. Since peripheral artery disease is a subtype of atherosclerotic disease, it is not surprising that people with peripheral artery disease have a higher risk of other atherosclerotic diseases such as coronary artery disease and stroke (Reference Aday and MatsushitaAday & Matsushita 2021).

Similar to the cerebral flow, wall shear stress has been emphasized to play an important role in peripheral artery disease: low wall shear stress may promote the build-up of plaque while high wall shear stress could increase the vulnerability of plaques (Reference Samady, Eshtehardi, McDaniel, Suo, Dhawan, Maynard, Timmins, Quyyumi and GiddensSamady et al. 2011; Reference Casa, Deaton and KuCasa, Deaton & Ku 2015). Simulations of local flow through thrombosis have quantified abnormal mean and oscillatory wall shear stress, which, when combined, may serve as risk factors for rupture (Reference Xu, Liu, Song, Chen, Wang, Zhang, Yan, Liu and HuangXu et al. 2016; Reference Ferrarini, Finotello, Salsano, Auricchio, Palombo, Spinella, Pane and ContiFerrarini et al. 2021; Reference Wang, Serracino-Inglott and FengWang, Serracino-Inglott & Feng 2021). Additionally, simulations have highlighted the importance of patient specificity of peripheral vascular anatomies. Reference Wood, Zhao, Zambanini, Jackson, Gedroyc, Thom, Hughes and XuWood et al. (2006) investigated the differences in flow patterns, curvature and tortuosity in the superficial femoral artery between young men and women. Their MRI-based computational modelling revealed that men exhibited greater tortuosity and curvature in the superficial femoral artery, along with lower mean wall shear stress, compared with women, which might contribute to the understanding of the two- to five-fold higher incidence of peripheral vascular disease in men.

Finally, a unique challenge of modelling peripheral artery blood flow is that lower limb movement significantly affects the haemodynamics in the lower peripheral vessel. To study this, Reference Colombo, Luraghi, Cestariolo, Ravasi, Airoldi, Chiastra and PennatiColombo et al. (2020) mimicked hip rotation, knee flexion and complete movement of walking with idealized femoropopliteal arteries and moving-boundary models. Interestingly, they showed that leg movements induce higher time-averaged wall shear stress, with flow rate amplitude and walking period being the most influential parameters that affect wall shear stress.

7. Microvascular and capillary flows

From centimetre-scale blood vessels arising from the heart, to micron-scale capillaries perfusing the various organs in the body, the interacting scales in the circulatory system span roughly six orders of magnitude. For example, distal capillary beds, consisting of millions of micron-sized blood vessels, respond to biomechanical stimuli and metabolic demand to change the flow within upstream millimetre-sized blood vessels, which can be as much as four-fold in the coronary circulation (Reference Wilson, Wyche, Christensen, Zimmer and LaxsonWilson et al. 1990), via autoregulatory mechanisms (Reference Goodwill, Dick, Kiel and TuneGoodwill et al. 2017). The flow of blood, which enters the aorta and its branches exhibiting Newtonian-like flow in the large vessels, ceases to behave as a continuum fluid as it traverses the $\mathrm {\mu }$m-sized arterioles and capillaries. This leads to a nonlinear as well as non-monotonic relationship between fluid viscosity and vessel diameter (Reference SecombSecomb 2017). The study of such multi-scale interaction necessitates multi-scale and multi-physics modelling approaches to tackle the prohibitive computational cost associated with simulating blood flow in these disparate scales.

The earliest physiologically motivated models of flow in the circulatory system considered the system as a so-called ‘lumped-parameter network’, drawing analogies between electric and hydraulic circuits (Reference Westerhof, Bosman, De Vries and NoordergraafWesterhof et al. 1969). These highly simplified models, also known as ‘zero-dimensional models’, have basic building blocks comprising resistors (R) to model viscous resistance within blood vessels derived from Poiseuille flow, capacitors (C) to model vessel wall distensibility and inductors (L) to model fluid inertia. The simplified equations that encapsulate these physical phenomena to relate flow ($Q$) and pressure ($P$), which are analogous to electric current and voltage, respectively, of blood with density $\rho $ and viscosity $\mu $ within a blood vessel with length $l$, radius $r$, wall thickness $h$ and Young's modulus $E$ are given by

(7.1) \begin{equation} \left.\begin{array}{c} \text{Resistance: } {\rm \Delta} P = RQ, \quad \text{where } R = \dfrac{8\mu l}{{\rm \pi} r^4},\\ \text{Capacitance: } {\rm \Delta} \dot{P} = \dfrac{Q}{C}, \quad \text{where } C = \dfrac{3{\rm \pi} l r^3}{2Eh},\\ \text{Inductance: } {\rm \Delta} P = L\dot{Q}, \quad \text{where } L = \dfrac{\rho l}{{\rm \pi} r^2}. \end{array}\right\} \end{equation}

Lumped-parameter models use suitable arrangements of the above resistors, capacitors and inductors to provide macroscopic descriptions of blood flow rates and pressures at specific locations in hydraulic circuits via a set of ordinary differential equations satisfying global mass conservation but only partial momentum conservation (Reference Pfaller, Pham, Verma, Pegolotti, Wilson, Parker, Yang and MarsdenPfaller et al. 2022). Figure 3 shows two examples of such zero-dimensional models. Although the simplicity of these models does not conventionally allow accurate modelling of complex flow topologies, such as vessel bifurcations and curvature, several studies have proposed improvements and corrections to account for flow nonlinearities while preserving the simplicity of these models (Reference Mynard and Valen-SendstadMynard & Valen-Sendstad 2015; Reference Chnafa, Valen-Sendstad, Brina, Pereira and SteinmanChnafa et al. 2017; Reference Mirramezani and ShaddenMirramezani & Shadden 2020). The minimal computational expense of zero-dimensional models has led to their wide use in low-order modelling of multi-scale circulatory phenomena by effectively ‘lumping’ together the resistance and capacitance of entire micro-vascular/capillary beds into single electric circuit-like components (Reference Vignon-Clementel, Alberto Figueroa, Jansen and TaylorVignon-Clementel et al. 2006; Reference Kim, Vignon-Clementel, Coogan, Figueroa, Jansen and TaylorKim et al. 2010).

Figure 3.

Examples of multiscale computational models. (a) A model of Fontan surgery that allows the quantification of local as well as whole-body haemodynamics (Reference SchwarzSchwarz et al. 2021). (b) A coronary flow model coupling 3-D flow in the aorta and large coronary vessels, 0-D flow in smaller myocardial vessels (shown in blue) and systemic circulation and Darcy flow in the myocardial microvasculature (Reference Menon, Seo, Fukazawa, Ogawa, Kahn, Burns and MarsdenMenon et al. 2024).

Another class of simplified vascular flow models that has enabled the modelling in large, multi-scale circulatory networks at reasonable computational expense are so-called ‘one-dimensional’ flow models (Reference Formaggia, Nobile, Quarteroni and VenezianiFormaggia et al. 1999, Reference Formaggia, Gerbeau, Nobile and Quarteroni2001). The description of flow within blood vessels is obtained by averaging the 3-D Navier–Stokes equations across the vessel cross-section and projecting this on the 1-D curve tracing the centre of blood vessels ($z$). This yields partial differential equations describing the flow rate, $Q(t,z)$, and pressure, $P(t,z)$, within vessels of cross-sectional area $S(t,z)$

(7.2)\begin{equation} \left.\begin{array}{c} \displaystyle \dfrac{\partial Q}{\partial t} + \dfrac{4}{3}\dfrac{\partial }{\partial z}\left(\dfrac{Q^2}{S}\right) + \dfrac{S}{\rho}\dfrac{\partial p}{\partial z} = Sf - \dfrac{8\mu}{\rho}\dfrac{Q}{S} + \dfrac{\mu}{\rho} \dfrac{\partial^2 Q}{\partial z^2},\\ \displaystyle \dfrac{\partial S}{\partial t} + \dfrac{\partial Q}{\partial z} = 0. \end{array}\right\} \end{equation}

The above equations are closed using constitutive models for the wall deformation that relate $P(t,z)$ to $S(t,z)$. These models provide higher-fidelity descriptions of the fluid dynamics than lumped-parameter models, particularly in vasculature with more complicated topology, while still remaining significantly less computationally expensive than full 3-D flow simulations.

These 0-D and 1-D descriptions of the circulation have been widely used in multi-scale flow simulations of extensive vascular networks (Reference Westerhof, Bosman, De Vries and NoordergraafWesterhof et al. 1969; Reference AvolioAvolio 1980; Reference OlufsenOlufsen 1999). They are most often used in conjunction with high-fidelity 3-D models of blood flow in a multi-scale manner, where 3-D blood flow simulations are performed in limited vascular segments of interest and 0-D or 1-D models approximate downstream/upstream circulatory behaviour that is not included in 3-D blood flow simulations (Reference Migliavacca, Balossino, Pennati, Dubini, Hsia, De Leval and BoveMigliavacca et al. 2006; Reference Vignon-Clementel, Alberto Figueroa, Jansen and TaylorVignon-Clementel et al. 2006; Reference Kim, Vignon-Clementel, Coogan, Figueroa, Jansen and TaylorKim et al. 2010; Reference Sankaran, Kim, Choi and TaylorSankaran et al. 2016; Reference Grande Gutierrez, Mathew, McCrindle, Tran, Kahn, Burns and MarsdenGrande Gutierrez et al. 2019; Reference Grande Gutiérrez, Sinno and DiamondGrande Gutiérrez, Sinno & Diamond 2022; Reference Menon, Seo, Fukazawa, Ogawa, Kahn, Burns and MarsdenMenon et al. 2024). This is shown in figure 3. Such multi-scale models have allowed researchers to investigate the relationship between widely disparate circulatory scales, such as the effect of exercise on changes in cardiac function and pulmonary circulation (Reference Kung, Pennati, Migliavacca, Hsia, Figliola, Marsden and GiardiniKung et al. 2014), haemodynamic changes in distal/proximal circulation to compensate for local abnormalities within specific vascular regions (Reference Balossino, Pennati, Migliavacca, Formaggia, Veneziani, Tuveri and DubiniBalossino et al. 2009; Reference SchwarzSchwarz et al. 2021) and the effects of local surgical interventions (e.g. shunt sizing) on global physiologic parameters such as oxygen delivery and cardiac workload in patients with congenital heart disease (Reference Migliavacca, Balossino, Pennati, Dubini, Hsia, De Leval and BoveMigliavacca et al. 2006; Reference Moghadam, Migliavacca, Vignon-Clementel, Hsia and MarsdenMoghadam et al. 2012; Reference Esmaily-Moghadam, Hsia and MarsdenEsmaily-Moghadam, Hsia & Marsden 2015).

Although their minimal computational expense has allowed the study of physiological interactions amongst different scales, the primary challenge associated with such multi-scale models is the relationship between their parametrization and the associated physiology they aim to capture. For example, in models consisting of 3-D descriptions of vascular haemodynamics coupled with reduced-order representations of the beating heart, it is not immediately obvious how the parameters of the heart model, usually consisting of resistors, time-varying capacitors and valves, relate to patients’ physiological cardiac volume, ejection fraction and other clinically relevant metrics. The need to associate these low-order flow models with clinical measurements for personalized modelling has led to development of parameter estimation and automated tuning schemes, which ideally should account for uncertainties (Reference Tran, Schiavazzi, Ramachandra, Kahn and MarsdenTran et al. 2017; Reference Fleeter, Geraci, Schiavazzi, Kahn and MarsdenFleeter et al. 2020; Reference Seo, Schiavazzi, Kahn and MarsdenSeo et al. 2020; Reference Menon, Khan, Sexton, Richter, Nieman and MarsdenMenon et al. 2024).

While the above low-order modelling strategies can be used to provide macroscopic descriptions of flow and pressure through lumped-parameter representations of microvasculature, they are computationally intractable when modelling flow that is resolved at the scale of each blood vessel in dense capillary networks within tissues. This is due to the complexity and computational expense of modelling millions of interconnected vessels, as well as the network-like structure of capillary beds where simplified Poiseuille flow, without accounting for viscous effects at junctions between blood vessels, may not be valid. A spatially averaged/coarse-grained description of flow is therefore required in such microvascular networks. In fact, dimensional analysis as well as microfluidic measurements have shown a linear relationship between the pressure gradient and spatially averaged flow within microvasculature (Reference FungFung 1997; Reference Stauber, Waisman, Korin and SznitmanStauber et al. 2017) – akin to Darcy's law governing flow through porous materials. Recent work has therefore led to the development of computational models for microvascular/capillary tissue perfusion that represent dense capillary networks and interstitial regions of tissues as porous materials.

These models have been applied to study capillary flows, oxygen or nutrient exchange, as well as the transport of therapeutics in various organs such as in alveoli within the lungs (Reference Erbertseder, Reichold, Flemisch, Jenny and HelmigErbertseder et al. 2012; Reference Zurita and HurtadoZurita & Hurtado 2022), within tumours (Reference Secomb, Hsu, Park and DewhirstSecomb et al. 2004; Reference Stylianopoulos and JainStylianopoulos & Jain 2013; Reference d'Espositod'Esposito et al. 2018), in the brain (Reference Sweeney, Walker-Samuel and ShipleySweeney, Walker-Samuel & Shipley 2018; Reference Kim, Rundfeldt, Lee and LeeKim et al. 2023) and in the heart (Reference Chapelle, Gerbeau, Sainte-Marie and Vignon-ClementelChapelle et al. 2010; Reference Cookson, Lee, Michler, Chabiniok, Hyde, Nordsletten, Sinclair, Siebes and SmithCookson et al. 2012; Reference MichlerMichler et al. 2013; Reference Di Gregorio, Fedele, Pontone, Corno, Zunino, Vergara and QuarteroniDi Gregorio et al. 2021, Reference Di Gregorio2022; Reference PapamanolisPapamanolis et al. 2021; Reference Kim, Rundfeldt, Lee and LeeKim et al. 2023). While the development of tissue perfusion models is relatively nascent compared with the 3-D and lumped-parameter models discussed above, in the context of flow through cardiac tissue, porous medium flow models have been coupled with high-fidelity and reduced-order models of flow in upstream coronary arteries (Reference HydeHyde et al. 2014; Reference Di Gregorio, Fedele, Pontone, Corno, Zunino, Vergara and QuarteroniDi Gregorio et al. 2021; Reference PapamanolisPapamanolis et al. 2021; Reference Kim, Rundfeldt, Lee and LeeKim et al. 2023; Reference Menon, Khan, Sexton, Richter, Nieman and MarsdenMenon et al. 2024) as well as the elasto-mechanics describing cardiac contraction and its effect on the pressure driving coronary flow within the cardiac tissue (Reference Chapelle, Gerbeau, Sainte-Marie and Vignon-ClementelChapelle et al. 2010; Reference Cookson, Lee, Michler, Chabiniok, Hyde, Nordsletten, Sinclair, Siebes and SmithCookson et al. 2012; Reference Di Gregorio, Fedele, Pontone, Corno, Zunino, Vergara and QuarteroniDi Gregorio et al. 2021). Figure 3(b) shows an example of a Darcy myocardial flow model coupled to upstream 3-D flow in coronary arteries and 0-D distal circulation models. However, as with the low-order micro-circulation models discussed above, a challenge with these approaches is the correspondence between modelling parameters, such as Darcy permeability and pressure, with quantities that can be measured in vivo using current clinical techniques. Recent advances in imaging methods hold promise for more accurate prescription of these parameters.

Delving even deeper into smallest scales of the circulatory system, i.e. the cellular level, fluid dynamics continues to play a prominent role. Typical Reynolds numbers in capillaries are $O(10^{-3})$, where Stokes flow assumptions are valid (Reference SecombSecomb 2017). However, the challenge at these scales is that the continuum description of blood breaks down. It is more accurately characterized as a complex fluid or suspension of predominantly red blood cells, white blood cells and platelets within plasma. These constituents exhibit a complex dynamics, such as migration away from vessel walls (creating a cell-free layer), intricate orientations and motions, as well as non-uniform concentrations within vessels as well as between different branching vessels (Reference Chien, Usami, Dellenback and GregersenChien et al. 1967a,Reference Chien, Usami, Dellenback, Gregersen, Nanninga and Guestb; Reference ChienChien 1970; Reference Fedosov, Noguchi and GompperFedosov, Noguchi & Gompper 2014; Reference Takeishi, Rosti, Imai, Wada and BrandtTakeishi et al. 2019). This occurs in response to various factors including their concentration, which varies with scale from ${\sim }50\,\%$ in most of the body to ${\sim }20\,\%$ in small vessels with $30\mathrm {\mu }$m diameters, flow-induced shear forces within vessels and their different material properties (Reference Zhao, Shaqfeh and NarsimhanZhao, Shaqfeh & Narsimhan 2012; Reference Qi and ShaqfehQi & Shaqfeh 2017). The cell-scale dynamics within small vessels has several implications, including for clot formation, drug delivery and, most importantly, blood viscosity (Reference ChienChien 1970). Interestingly, the dynamics and deformability of blood's constituents lead to the effective viscosity of blood being only 4–5 times that of plasma, compared with an equivalent suspension of rigid spheres which would be 10–100 times as viscous as its solvent (Reference VerziccoVerzicco 2022). Shear-induced cell migration and deformation within vessels also leads to blood acting as a shear-thinning fluid and produces a non-monotonic trend in blood viscosity with respect to vessel diameter (Reference Pries, Neuhaus and GaehtgensPries, Neuhaus & Gaehtgens 1992). Effective blood viscosity first decreases with increasing diameter due to the cell-free near vessel walls. This is referred to as the Fåhræus–Lindqvist effect (Reference Fåhræus and LindqvistFåhræus & Lindqvist 1931). It reaches a minimum in vessels with diameters of the order of red blood cell diameters, where red blood cells line up in a single-file-like motion, and then increases in vessels with smaller diameters due to the deformation of red blood cells. A review of these effects is available in Reference SecombSecomb (2017), and efforts to simulate them are discussed in Reference FreundFreund (2014).

8. Venous flow

After perfusing tissues throughout the body, blood begins its return journey to the heart via the venous system. The venous circulation consists of veins that carry blood depleted of oxygen and nutrients, and stowed with metabolic wastes, from deep within organs, tissues and cells to the great veins that return to the heart. Compared with arteries, veins operate under lower pressure with low pulsatility, and therefore have thinner walls. For the large veins close to the heart, such as the vena cava, there is some mild pulsatility due to right atrial contraction and respiratory effects.

A unique feature of venous haemodynamics is the presence of one-way valves within most veins in the body, adding complexity to venous flow simulation. These valves are essential for maintaining unidirectional venous flow, especially given the low pressure under which they operate. Their dysfunction, which impairs the ability of the venous network to return blood to the heart, is a primary cause of chronic venous inefficiency and venous hypertension (Reference Bergan, Schmid-Schönbein, Smith, Nicolaides, Boisseau and EklofBergan et al. 2006). Although these valves operate in very low-speed flows, the formation of an axial and helical jet, the presence of flow recirculation within the sinus, and the pulsatility of flow resulting from the opening and closure of valves have been observed and are hypothesized to play roles in healthy venous function (Reference Hamer, Malone and SilverHamer, Malone & Silver 1981; Reference Lurie, Kistner, Eklof and KesslerLurie et al. 2003; Reference Lurie and KistnerLurie & Kistner 2013). In this context, another interesting feature is that when multiple valves are present at venous junctions, the valve openings are consistently observed to be at an angle to each other and this angle correlates with the distance between the valves (Reference Lurie and KistnerLurie & Kistner 2012). This asymmetry naturally gives rise to helical flows that have been shown to improve efficiency and reduce flow stagnation – thus preventing thrombosis/clotting (Reference Lurie and KistnerLurie & Kistner 2013; Reference Chen, Diaz, Lurie, Chambers and KassabChen et al. 2018).

The formation of clots within veins – or deep venous thrombosis – presents a significant health concern, with an annual incidence of approximately 1 per 1000 adults in the United States (Reference CushmanCushman 2007). Physical immobility is one of the most significant risk factor for deep venous thrombosis because the absence of muscular activity leads to a loss of oscillatory blood flow around valve areas in veins, thereby promoting thrombus or clot formation (Reference WelshWelsh et al. 2019). Subsequently, these clots can travel to the lungs and get lodged in the pulmonary arteries, leading to a potentially fatal pulmonary embolism. In its prevention, alongside anticoagulation therapy, placing a filter in the inferior vena cava serves as an alternative strategy. This mechanical device is designed to filter thrombi carried by the blood flow within the vein. Using in vitro experiments and numerical simulation, the design of these filters has been characterized and optimized with a focus on maximizing thrombus capture efficiency, minimizing areas of low wall shear stress, and ensuring minimal disruption to blood flow (Reference Singer, Wang and DiachinSinger, Wang & Diachin 2010; Reference Feng, Li and FengFeng, Li & Feng 2022). We recommend the recent review by Reference Watson, Saaid, Vedula, Cardenas, Henke, Nicoud, Xu, Hunt and ManningWatson et al. (2024) for a discussion of modelling efforts in this domain to understand flow-mediated mechanisms driving deep vein thrombosis.

Beyond the pulmonary circulation, abnormalities in venous flow can have far-reaching health effects because the entire circulatory system is intricately linked with other bodily fluid systems. For instance, global multi-scale models of human circulation with its interactions with central nervous system fluids have been developed to explore how venous irregularities might influence the development and progression of neurological diseases (Reference Müller and ToroMüller & Toro 2014; Reference Toro, Celant, Zhang, Contarino, Agarwal, Linninger and MüllerToro et al. 2022).

A theme of central focus, with regard to the haemodynamics of the venous system, is palliative approaches for patients with single ventricle physiology. This is a group of congenital heart defects in which one ventricle (either the left or right) is too underdeveloped to function normally. Typically, these patients undergo a series of three surgeries between birth and 4 years of age to reconfigure their systemic and pulmonary circulations sequentially so that they are driven by the single functioning ventricle. These surgeries are primarily centred around the great veins that return blood to the heart, with the aim of maintaining passive flow in these and the pulmonary vessels due to the absence of the second functioning ventricle. While these surgeries are life saving, the associated complications have spurred the use of CFD for improving surgical designs virtually. In these applications, 3-D computational models are often integrated with lumped-parameter networks, thus offering insights into the local haemodynamics while also capturing the interactions with the broader circulatory system (Reference Corsini, Cosentino, Pennati, Dubini, Hsia and MigliavaccaCorsini et al. 2011; Reference KungKung et al. 2020; Reference SchwarzSchwarz et al. 2021).

During the stage one, Norwood surgery, the pulmonary circulation is connected to the aorta with a shunt to provide blood flow to the lungs. Early computational studies proposed and examined several variations in the position and size of the shunt, and found the larger shunts result in poorer oxygen delivery while leading to increased pulmonary flow (Reference Migliavacca, Pennati, Dubini, Fumero, Pietrabissa, Urcelay, Bove, Hsia and de LevalMigliavacca et al. 2001; Reference Hsia, Cosentino, Corsini, Pennati, Dubini and MigliavaccaHsia et al. 2011). The effect of various degrees of pulmonary banding in the hybrid procedure has been simulated to determine the optimal systemic to pulmonary flow ratio and oxygen delivery (Reference Corsini, Cosentino, Pennati, Dubini, Hsia and MigliavaccaCorsini et al. 2011). Moreover, modelling provides valuable quantitative insight into the optimization of the neo-aortic reconstruction. With patient specific geometries, a smooth aortic arch and a larger anastomosis have been shown to result in higher energy efficiency in the blood flow (Reference Qian, Liu, Itatani, Miyaji and UmezuQian et al. 2010; Reference Itatani, Miyaji, Qian, Liu, Miyakoshi, Murakami, Ono and UmezuItatani et al. 2012).

In the second stage, the Glenn procedure, the superior vena cava is directly connected to the pulmonary arteries. By allowing passive blood flow to the lungs from the upper body, it reduces the workload on the single ventricle. The connection can be achieved through either the hemi-Fontan or bidirectional Glenn (Reference Douglas, Goldberg, Mosca, Law and BoveDouglas et al. 1999; Reference Kogon, Plattner, Leong, Simsic, Kirshbom and KanterKogon et al. 2008), but questions remain as to which results in the best haemodynamic results. Patient-specific simulation initially indicated that the Glenn connection is haemodynamically more efficient compared with hemi-Fontan, due to its direct end-to-side anastomosis (Reference Pekkan, Dasi, de Zélicourt, Sundareswaran, Fogel, Kanter and YoganathanPekkan et al. 2009). Despite what appears to be significant local power loss differences, follow up studies found no notable distinctions in haemodynamic or physiological parameters between the two models (Reference KungKung et al. 2020).

After the third surgery, known as the Fontan procedure, the total cavopulmonary connection is established with inferior and superior vena cava connected to the pulmonary arteries. This distinct circulation causes elevated central venous pressure and reduced pulmonary pulsatility, which are associated with morbidities such as liver fibrosis, exercise intolerance and heart failure (Reference Goldberg, Shaddy, Ravishankar and RychikGoldberg et al. 2011). Besides assessing the cavopulmonary connection at rest, computational frameworks have also offered valuable insights into the impacts of respiration and exercise in Fontan patients (Reference Marsden, Vignon-Clementel, Chan, Feinstein and TaylorMarsden et al. 2007; Reference Kung, Pennati, Migliavacca, Hsia, Figliola, Marsden and GiardiniKung et al. 2014). They reveal that the elevated resistance at connection and nonlinearly increased energy loss significantly restrict the Fontan circulation's capacity to augment cardiac output during exercise, demonstrating a limited ability to accommodate exercise conditions (Reference Whitehead, Pekkan, Kitajima, Paridon, Yoganathan and FogelWhitehead et al. 2007; Reference Sundareswaran, Pekkan, Dasi, Whitehead, Sharma, Kanter, Fogel and YoganathanSundareswaran et al. 2008). The connection was then proposed to be offset to prevent the head-on collision of flow (Reference De Leval, Dubini, Migliavacca, Jalali, Camporini, Redington and PietrabissaDe Leval et al. 1996), but this design introduces uneven hepatic flow distribution and therefore uneven growth in pulmonary vascular tree. Subsequently, the bifurcated Y-graft was introduced and improved energy loss and hepatic flow distribution compared with the offset design in computational studies that were later validated clinically (Reference Marsden, Bernstein, Reddy, Shadden, Spilker, Chan, Taylor and FeinsteinMarsden et al. 2009; Reference Yang, Vignon-Clementel, Troianowski, Reddy, Feinstein and MarsdenYang et al. 2012; Reference Martin, Feinstein, Chan, Marsden, Yang and ReddyMartin et al. 2015; Reference Trusty, Wei, Slesnick, Kanter, Spray, Fogel and YoganathanTrusty et al. 2019, Reference Trusty, Wei, Sales, Kanter, Fogel, Yoganathan and Slesnick2020). Figure 4(a) shows examples of flow fields resulting from these offset and Y-graft configurations. It is worth mentioning that, although power loss is a major focus of these computational studies, Fontan haemodynamics from 100 patients showed that Fontan power loss varies from patient to patient, and elevated levels are correlated with lower systemic flow and cardiac index (Reference Haggerty, Restrepo, Tang, de Zélicourt, Sundareswaran, Mirabella, Bethel, Whitehead, Fogel and YoganathanHaggerty et al. 2014). As their predictive power has grown, computational simulations continue to be used in investigations of new palliative techniques, such as models of tissue-engineered vascular grafts in the Fontan pathway (Reference SchwarzSchwarz et al. 2021).

Figure 4.

(a) Haemodynamics associated with standard offset grafts and novel Y-grafts in the Fontan procedure to treat single ventricle physiology (Reference Marsden, Reddy, Shadden, Chan, Taylor and FeinsteinMarsden et al. 2010; Reference Yang, Chan, Reddy, Marsden and FeinsteinYang et al. 2015). The highly disturbed flow and the split between left and right pulmonary arteries are ongoing challenges. (b) Patient-specific models for treating pulmonary stenoses. Comparing exclusively proximal stenting vs more extensive stenting showed improved pressure alleviation with extensive stenting (Reference Lan, Yang, Feinstein, Kreutzer, Collins, Ma, Adamson and MarsdenLan et al. 2022).

9. Right heart

After its long journey from large blood vessels through the depths of tissues and eventually through the venous system, de-oxygenated blood returns to the right side of the heart, which consists of the right atrium and the right ventricle. The superior and inferior vena cava first direct de-oxygenated blood into the right atrium. From there, it flows into the right ventricle, which is the secondary pumping chamber of the heart, responsible for driving flow to the lungs to be re-oxygenated. Blood moving from the right atrium to the right ventricle flows through the tricuspid valve that separates these chambers, and is initially driven by its own inertia. Eventually, an active right atrial contraction completes the ventricular filling. Subsequently, the right ventricle contracts to increase the blood pressure therein, thus closing the tricuspid valve and opening the pulmonary valve to pump blood into the lungs. The left and right sides of the heart, working in parallel, operate synchronously to ensure equal volumes of blood flow, maintaining equilibrium in the circulatory system. The atrial and ventricular septums separate the right and left atria and ventricles, respectively.

Despite similarities to the left heart, the right atrium possess unique anatomical and functional traits. The right atrium connects to the right ventricle through the tricuspid valve, while the left atrium links to the left ventricle via the bicuspid mitral valve. Both atria feature trabeculated appendages that may enhance stagnant flow and promote clot formation (or thrombosis). However, the right atrial appendage has a larger junction than its left counterpart, which likely accounts for reduced thrombosis in the right atrium appendage compared with the left (see § 11) (Reference Lang, Cameli, Sade, Faletra, Fortuni, Rossi and Soulat-DufourLang et al. 2022). Flow patterns within the right atrium are significantly influenced by the inflows from the superior and inferior vena cava. Studies have observed that these incoming flows engage in a circular motion rather than colliding directly, which results in predominantly rotational energy that tends to form a single vortex in the right atrium (Reference Kilner, Yang, Wilkes, Mohiaddin, Firmin and YacoubKilner et al. 2000). The kinetic energy is better preserved compared with configurations with multiple vortices or less prevalent helical flow, and conversely, atria without this single vortex flow tend to rely more on energy from atrial contractions (Reference Dewhurst, Coats, Parikh, Hollingsworth and GanDewhurst et al. 2020).

Often termed the ‘forgotten ventricle’, the right ventricle also exhibits distinct characteristics from the left ventricle. It is a thin-walled crescent-shaped chamber that wraps around the left ventricle and pumps against the low-impedance pulmonary circuit at a lower pressure of approximately 5–20 mmHg (667–2650 Pa). As such, its mass is roughly 1/6th that of the left ventricle, and its stroke work is approximately a quarter that of the left ventricle. Interestingly, although both right and left ventricular contractions are predominantly driven by shortening of the long-axis, left ventricular contractions significantly affect right ventricular function. The bulging of the ventricular septum due to high pressure in left ventricle reduces the right ventricle volume, and about 20 % to 40 % of the right ventricular systolic pressure and volume outflow result from left ventricular contraction (Reference Santamore and Dell'ItaliaSantamore & Dell'Italia 1998; Reference Gaine, Naeije and PeacockGaine et al. 2014). The right ventricle is also more afterload sensitive. When afterload rises, the right ventricle experiences a significant surge in energy use and a decrease in efficiency compared with the left ventricle, causing the right ventricle's stroke volume to decrease more rapidly. This sensitivity could explain why systemic arterial hypertension is better tolerated than pulmonary arterial hypertension (Reference Bravo, Navarro, Dhaliwal, Khorsandi, Keenan, Mudigonda, O'Brien and MahrBravo et al. 2022).

The flow field in the right ventricle differs from the left ventricle largely due to its different geometry. A number of MRI and CFD studies have shown that the general flow path through the right ventricle is dominated by the smooth curvature of the chamber with a less acute change in direction between the inflow and the outflow tract (Reference Kilner, Yang, Wilkes, Mohiaddin, Firmin and YacoubKilner et al. 2000; Reference Fredriksson, Zajac, Eriksson, Dyverfeldt, Bolger, Ebbers and CarlhällFredriksson et al. 2011). A ring-like vortex forms behind the tricuspid valve at the early diastole, which interacts with the narrower boundaries of the right ventricle, particularly on the septal side that is flatter, quickly breaking into a weakly turbulent flow pattern (Reference Kilner, Yang, Wilkes, Mohiaddin, Firmin and YacoubKilner et al. 2000; Reference Markl, Kilner and EbbersMarkl, Kilner & Ebbers 2011). Thus, the flow in the right ventricle features a higher dissipative state, whereas the left ventricle preserves more circulation and kinetic energy throughout diastole. Additionally, its shape helps streamline the flow, promoting the formation of a helical flow towards the pulmonary orifice (Reference Mangual, Domenichini and PedrizzettiMangual, Domenichini & Pedrizzetti 2012). During systolic contraction, the residual vorticity extends along the outflow tract, adding a slightly helical pattern to the largely irrotational velocity field, and assists in mixing of the blood (Reference Collia, Zovatto, Tonti and PedrizzettiCollia et al. 2021). A number of CFD simulations have also characterized the effect of endocardial structures, such as trabeculae and papillary muscles, on disrupting the main vortex ring and dissipating energy in intra-ventricular flows (Reference Sacco, Paun, Lehmkuhl, Iles, Iaizzo, Houzeaux, Vázquez, Butakoff and Aguado-SierraSacco et al. 2018a). However, as with the left ventricle, this is usually ignored in both simulations and imaging due to the challenge of simultaneously resolving both the macro-scale flow features and the geometric details of the endocardial surface.

Flow modelling techniques have also been expanded to describe changes linked to pathological conditions, serving as markers for impaired right ventricular function. For example, pulmonary arterial hypertension patients, who have increased pressure and diminished ejection fraction, exhibit a significant increase in right ventricle kinetic energy work density and pronounced pulmonary artery energy loss across the cardiac cycle, which correlates with right ventricle ejection fraction (Reference Han, Witschey, Fang-Yen, Arkles, Barker, Forfia and HanHan et al. 2015). Moreover, haemodynamic analysis has revealed an increase in both right atrium and right ventricle vorticity at late diastole in pulmonary artery hypertension patients, suggesting its potential interplay between right ventricle fluid dynamics and tissue mechanics leading to pulmonary artery hypertension (Reference Fenster, Browning, Schroeder, Schafer, Podgorski, Smyser, Silveira, Buckner and HertzbergFenster et al. 2015).

10. Pulmonary circulation

The pulmonary circulation carries oxygen-depleted blood from the right side of the heart to the lungs for oxygen enrichment, and then channels the oxygen-saturated blood back to the left side of the heart. Starting from the right ventricular outflow tract, blood is pumped into the main pulmonary artery, which splits into the left and right pulmonary arteries. As these arteries penetrate the lungs, they further subdivide into smaller vessels, following a pattern analogous to the bronchial tree in the respiratory system. At their finest level, these vessels sprawl out across the alveoli, the tiny air sacs where the critical process of gas exchange occurs. Here, carbon dioxide is expelled from the blood and replaced with fresh oxygen. This newly oxygenated blood is then collected by the pulmonary veins and directed back to the left atrium of the heart, ready to be pumped to the rest of the body.

Owing to the central role of blood oxygenation and its subsequent transport, the focus of pulmonary haemodynamics studies has generally revolved around two important pathologies – pulmonary arterial hypertension and pulmonary stenosis. Pulmonary arterial hypertension is a chronic disease characterized by progressive elevation in pulmonary artery pressure due to increased pulmonary vascular resistance. The progression of pulmonary artery hypertension is driven by the remodelling of vascular tissue in response to elevated haemodynamic pressure. The wall thickness and vascular compliance decrease during this remodelling, which leads to increased right ventricle afterload and ventricle–vascular decoupling (Reference Kheyfets, O'Dell, Smith, Reilly and FinolKheyfets et al. 2013). A primary clinical diagnostic indicator of pulmonary arterial hypertension is the pulmonary vascular resistance. It is calculated as the ratio of the pressure drop through the pulmonary circulation to pulmonary blood flow and signifies the structural changes at the pulmonary arteriole level. Interestingly, this widely used clinical metric does not reflect the fact that blood flow is pulsatile with reflective waves. Therefore, several studies have proposed pulmonary impedance to better predict patient outcomes by taking into account not only resistance but also elasticity and inertia to characterize the steady and oscillatory behaviours of blood flow and dynamic interactions between the heart and the lungs, especially under varying physiological conditions (Reference Hunter, Lee, Lanning, Ivy, Kirby, Claussen, Chan and ShandasHunter et al. 2008; Reference NaeijeNaeije 2012; Reference KheyfetsKheyfets et al. 2015).

Pulmonary stenosis, or narrowing of the pulmonary vessel lumen, is another significant cause of pulmonary flow disruption. Peripheral pulmonary artery stenosis is commonly associated with congenital heart diseases, including Williams syndrome, Alagille syndrome and tetralogy of Fallot. Severe peripheral pulmonary artery stenosis results in right ventricular hypertension, right ventricular hypertrophy and eventually right ventricular failure. Additionally, if pulmonary blood flow is not promptly restored, it can lead to irreversible damage to the pulmonary artery (Reference Mainwaring and HanleyMainwaring & Hanley 2016). Interventional procedures alleviate the stenosis with catheter-based angioplasty or surgical homograft patches. Nonetheless, such surgeries are challenging since they require long hours of cardiopulmonary bypass and address only the most proximal pulmonary arteries. Moreover, post-operative complications include persistent right ventricle hypertension, tears in the vessel wall, aneurysm formation, in-stent restenosis and even death (Reference Cunningham, McElhinney, Gauvreau, Bergersen, Lacro, Marshall, Smoot and LockCunningham et al. 2013; Reference Felmly, Mainwaring, Collins, Lechich, Martin, Ma and HanleyFelmly et al. 2023).

In light of these debilitating pathologies, both CFD and 4-D flow MRI have been used to investigate blood flow patterns that provide insights into the physiological mechanisms behind disease progression and clinical outcomes. For example, MRI data have revealed that pulmonary artery hypertension patients exhibit reduced flow, peak velocities, wall shear stress, helicity, vorticity and strain (Reference Tang, Pickard, Chan, Tsao, Taylor and FeinsteinTang et al. 2012; Reference Barker, Roldán-Alzate, Entezari, Shah, Chesler, Wieben, Markl and FrançoisBarker et al. 2015; Reference Schäfer, Barker, Kheyfets, Stenmark, Crapo, Yeager, Truong, Buckner, Fenster and HunterSchäfer et al. 2017; Reference Sieren, Berlin, Oechtering, Hunold, Drömann, Barkhausen and FrydrychowiczSieren et al. 2019). Increased recirculation and large vortex formation along the inferior curvature of the main pulmonary artery are also reported in pulmonary artery hypertension patients, which can serve as marker for pulmonary artery hypertension risk stratification. This formation of large-scale vortices stems from vessel dilation and diminished right ventricular function, creating shear layers with dramatically different velocities resulting in flow instabilities. This also explains the decreased wall shear stress at the proximal artery wall, which leads to endothelial dysfunction and maladaptive mechanobiological responses (Reference Reiter, Reiter, Kovacs, Kainz, Schmidt, Maier, Olschewski and RienmuellerReiter et al. 2008, Reference Reiter, Reiter, Kovacs, Stalder, Gulsun, Greiser, Olschewski and Fuchsjäger2013).

While the resolution limitations of MRI hinder its effectiveness in distal branches where disease is thought to initiate, computational simulations can provide an in-depth assessment of haemodynamics in these farther reaches of the pulmonary arterial tree. Automatic algorithms have been developed to segment the pulmonary vasculature below 1 mm resolution (Reference Doel, Gavaghan and GrauDoel, Gavaghan & Grau 2015). With more detailed anatomical models, blood flow simulations have reported decreased wall shear stress in the proximal arterial tree in severe pulmonary artery hypertension, as well as increased oscillatory shear index and relative residence time in the distal arteries, which indicate more stagnant blood flow (Reference Tang, Pickard, Chan, Tsao, Taylor and FeinsteinTang et al. 2012; Reference KheyfetsKheyfets et al. 2015; Reference Yang, Dong, Rabinovitch, Chan, Marsden and FeinsteinYang et al. 2019; Reference Pillalamarri, Piskin, Patnaik, Murali and FinolPillalamarri et al. 2021). Computational models have allowed the exploration of pulmonary artery hypertension linked to congenital heart defects, and found increased wall shear stress with increasing ventricular septal defect size in paediatric pulmonary artery trees (Reference Dong, Lan, Yang, Rabinovitch, Feinstein and MarsdenDong et al. 2021). Moreover, recent work has attempted to couple models of pulmonary haemodynamics with vascular growth and remodelling to simulate mechanobiological metrics of interest and potential treatment strategies for pulmonary artery hypertension patients (Reference Szafron, Yang, Feinstein, Rabinovitch and MarsdenSzafron et al. 2023). Image-based CFD can also serve as a virtual platform to evaluate interventions on a patient-specific basis to optimize surgical procedures. Reference Yang, Feinstein and Vignon-ClementelYang, Feinstein & Vignon-Clementel (2016) and Reference Yang, Hanley, Chan, Marsden, Vignon-Clementel and FeinsteinYang et al. (2018) performed virtual pulmonary artery reconstruction for patients with Alagille syndrome and achieved accurate predictions of postoperative pulmonary artery flows via adaptive outflow boundary conditions. Recently, Reference Lan, Yang, Feinstein, Kreutzer, Collins, Ma, Adamson and MarsdenLan et al. (2022) performed virtual intervention in Williams and Alagille syndrome patients and showed that stenting only proximal lesions is often insufficient to reduce overall vascular resistance and thus right ventricular pressures (figure 4b).

A challenge in such computational models is that of prescribing the boundary conditions required to accurately simulate the haemodynamics both before and after clinical interventions. For example, outflow boundary conditions substantially impact the flow distribution in the distal vessels, which is a significant predictor of pulmonary hypertension severity (Reference KheyfetsKheyfets et al. 2015). Models also need to consider flow autoregulation in distal pulmonary vasculature to correctly predict haemodynamics post-treatments (Reference Yang, Feinstein and Vignon-ClementelYang et al. 2016, Reference Yang, Hanley, Chan, Marsden, Vignon-Clementel and Feinstein2018). Another challenge is assessing blood flow in distal pulmonary arteries with diameters between 50 and 500 $\mathrm {\mu }$m, which are beyond typical modelling resolution. To address this, 3-D models have been combined with 0-D morphometric trees that consist of resistances in parallel and series based on existing pulmonary artery morphometric data of vessel segment diameters, lengths and branching patterns (Reference Huang, Yen, McLaurine and BledsoeHuang et al. 1996). Results from models with morphometric trees have uncovered increased wall shear stress in the distal small vessels compared with controls (Reference Yang, Dong, Rabinovitch, Chan, Marsden and FeinsteinYang et al. 2019), as well as elevated pulmonary vascular resistance and hypertension due to rarefaction in the pulmonary circulation (Reference Olufsen, Hill, Vaughan, Sainsbury and JohnsonOlufsen et al. 2012). They have also reported decreased shear stress along with increased cyclic stretch due to pulmonary vascular resistance at the microvascular level (Reference Bartolo, Qureshi, Colebank, Chesler and OlufsenBartolo et al. 2022). While it is challenging to accurately recapitulate the wide spatial scales and physiological responses at these scales in computational models, continued progress is important for clinical translation, especially due to the difficulty in accurately resolving these scales using clinical imaging.

11. Left atrium

At the end of our tour, we return to the left side of the heart where oxygen-rich blood returning from the lungs via the pulmonary veins collects in the left atrium. From here, blood is directed back into the left ventricle, via contraction of the left atrium, to restart its journey through the circulatory system.

Haemodynamic studies have revealed intricate flow patterns within the left atrium. The inflow from pulmonary veins are asymmetric, with left pulmonary veins forming a rotational loop while right pulmonary veins are more aligned with the wall (Reference Fyrenius, Wigström, Ebbers, Karlsson, Engvall and BolgerFyrenius et al. 2001). Such asymmetric filling generates asymmetric atrial flow pathways and tends to preserve blood momentum during the cardiac cycle (Reference Kilner, Yang, Wilkes, Mohiaddin, Firmin and YacoubKilner et al. 2000). During systole and mid-diastole, discrete left atrial vortices have been observed both from clinical (Reference Fyrenius, Wigström, Ebbers, Karlsson, Engvall and BolgerFyrenius et al. 2001; Reference ParkPark et al. 2013) and computational data (Reference Vedula, George, Younes and MittalVedula et al. 2015; Reference Zingaro, Dede, Menghini and QuarteroniZingaro et al. 2021). The formation of atrial vortices is conjectured to produce a better ‘washout’ effect in the left atrium and avoid intra-atrial flow stasis (Reference Fyrenius, Wigström, Ebbers, Karlsson, Engvall and BolgerFyrenius et al. 2001).

In terms of the haemodynamics, the most interesting feature of the left atrium is the left atrial appendage. Anatomically, the left atrium has a smooth-walled interior, except for the left atrial appendage protruding from its upper left corner. The left atrial appendage is a small and ear-shaped pouch, featuring a ridged interior due to the presence of pectinate muscles. It is more distensible than the rest of the left atrium and is hypothesized to act as a decompression chamber during left ventricular systole and periods of high atrial pressure (Reference Al-Saady, Obel and CammAl-Saady, Obel & Camm 1999). However, its structure has a propensity to cause increased flow stasis – which can lead to the formation of blood clots and the risk of stroke. This is particularly important in the context of atrial fibrillation, which is the most prevalent and clinically important cardiac arrhythmia. In patients with this condition, the heart's electrical signals are chaotic, leading to the inefficient atrial pumping. This can result in an irregular and frequently rapid heart rate, along with the risk of the blood pooling and potential clotting. If untreated, atrial fibrillation increases the risk of heart-related complications, such as stroke or heart failure. Notably, over 90 % of left atrial thrombi (or blood clots) originate in the left atrial appendage (Reference Al-Saady, Obel and CammAl-Saady et al. 1999). Flow patterns within the left atrial appendage have long been suggested as potential predictors for these events (Reference Agmon, Khandheria, Gentile and SewardAgmon et al. 1999; Reference Kamp, Verhorst, Welling and VisserKamp et al. 1999).

Therefore, the integration of computational methods with medical imaging has been used to provide a comprehensive view of the haemodynamics within the left atrial appendage. Decreased mean and peak velocities in left atrial appendage have been reported for atrial fibrillation patients, thus making the washout field incomplete and more prone to fluid stagnation. A large vortex has been observed in the centre of the left atrium for atrial fibrillation patients, instead of the multiple discrete vortices in healthy left atrium (Reference ParkPark et al. 2013). Four typical left atrial appendage morphologies and their impacts on the haemodynamics have been investigated using CFD, showing the close relation between morphology and thrombosis, in agreement with clinical reports. By modelling flow in distinct left atrial appendage morphologies, simulation results have underscored the critical importance of patient-specific left atrial appendage geometric characteristics in influencing thromboembolic risks (Reference Bosi, Cook, Rai, Menezes, Schievano, Torii and BurriesciBosi et al. 2018; Reference Masci, Barone, Dedè, Fedele, Tomasi, Quarteroni and CorsiMasci et al. 2019). While various CFD models have been employed for the left atrium to address potential shifts from transitional to turbulent flows, it is noteworthy that, in atrial fibrillation contexts, models without turbulence remain viable (Reference Dueñas-Pamplona, García, Sierra-Pallares, Ferrera, Agujetas and Lopez-MinguezDueñas-Pamplona et al. 2021). For more details about computational modelling of the left atrium and the left atrial appendage, we refer to the review by Reference Valvez, Oliveira-Santos, Piedade, Gonçalves and AmaroValvez et al. (2023).

12. Future outlook

As we conclude our journey through the human circulatory system, we hope to have impressed upon the reader the central role that fluid dynamics plays in a wide range of circulatory phenomena. The decades of research that have uncovered these insights have been primarily motivated by one goal – to use physics and engineering insights in advancing clinical care. To that end, several ongoing challenges remain – especially in the context of widespread clinical translation of computational models.

Chief amongst these is the pressing need for more comprehensive verification and validation of realistic patient-specific computational models and their associated insights in large-scale clinical studies (Reference SteinmanSteinman et al. 2012). This is a theme that was discussed several times in preceding sections of this review, which highlights its importance in all areas of the circulation. While there have been efforts to validate computational models with in vitro flows (Reference Emendi, Sturla, Ghosh, Bianchi, Piatti, Pluchinotta, Giese, Lombardi, Redaelli and BluesteinEmendi et al. 2021; Reference Kaiser, Schiavone, Elkins, McElhinney, Eaton and MarsdenKaiser et al. 2023; Reference Zimmermann, Bäumler, Loecher, Cork, Marsden, Ennis and FleischmannZimmermann et al. 2023), there is an urgent need for in vivo validation using patient-specific clinical data.

A related ongoing challenge is the need to develop methods for more personalized computational models that are informed by a variety of patient-specific clinical measurements – not just anatomical imaging, as is currently the norm. This requires the combination of accurate models with data-driven methods to tune model parameters based on routine as well as novel clinical measurements/imaging (Reference Tran, Schiavazzi, Ramachandra, Kahn and MarsdenTran et al. 2017; Reference Menon, Seo, Fukazawa, Ogawa, Kahn, Burns and MarsdenMenon et al. 2024). The translation of such models to clinical use also critically depends on the use of rigorous uncertainty quantification techniques (Reference Fleeter, Geraci, Schiavazzi, Kahn and MarsdenFleeter et al. 2020).

We should also highlight emerging tools to process and interpret the fluid dynamics of cardiovascular function. This includes methods that enable new and higher quality insights from clinical imaging, such as super-resolution and de-noising (Reference Fathi, Perez-Raya, Baghaie, Berg, Janiga, Arzani and D'SouzaFathi et al. 2020). In this context, there have also been several recent efforts to enhance the utility of low-resolution clinical imaging by augmenting such measurements with physics-based and data-enabled techniques to extract the clinically relevant haemodynamics that is otherwise under-resolved (Reference Bakhshinejad, Baghaie, Vali, Saloner, Rayz and D'SouzaBakhshinejad et al. 2017; Reference Habibi, D'Souza, Dawson and ArzaniHabibi et al. 2021; Reference Zhang, Brindise, Rothenberger, Markl, Rayz and VlachosZhang et al. 2022a,Reference Zhang, Rothenberger, Brindise, Markl, Rayz and Vlachosb). Finally, new analysis techniques focused on topological features of the flow field itself, including Lagrangian coherent structures and the topology of the wall shear stress field, have recently been shown to provide promising clinical indicators for clinical risk and disease progression (Reference Arzani, Gambaruto, Chen and ShaddenArzani et al. 2017; Reference Arzani and ShaddenArzani & Shadden 2018; Reference Morbiducci, Mazzi, Domanin, De Nisco, Vergara, Steinman and GalloMorbiducci et al. 2020; Reference Mazzi, De Nisco, Hoogendoorn, Calò, Chiastra, Gallo, Steinman, Wentzel and MorbiducciMazzi et al. 2021).

As important as fluid dynamics is to cardiovascular health, future treatments will hinge on more comprehensive predictive models that incorporate new findings from a range of different disciplines, including mechanobiology, genetics, and cell signalling (Reference Schwarz, Pegolotti, Pfaller and MarsdenSchwarz et al. 2023a). This has led to the development of multi-physics models, such as whole-heart models, which combine haemodynamics with the spatio-temporal behaviour of the heart at molecular, cellular and tissue scale – from ionic activation in electrophysiology to the tissue mechanics leading to cardiac contraction – to provide a more detailed description of cardiac function (Reference Sugiura, Okada, Washio and HisadaSugiura et al. 2022; Reference VerziccoVerzicco 2022; Reference Bucelli, Zingaro, Africa, Fumagalli, Dede’ and QuarteroniBucelli et al. 2023). Models combining haemodynamics with mechanobiology are also being developed to predict vascular growth and remodelling, which occur due to the interplay between haemodynamics and biomechanical stimuli experienced by vessel walls, in realistic patient-specific models (Reference Schwarz, Pfaller, Szafron, Latorre, Lindsey, Breuer, Humphrey and MarsdenSchwarz et al. 2023b). Finally, while computational modelling has provided many of the crucial insights discussed in this review, such modelling is often not possible in time frames that allow their widespread clinical uptake. To address this, machine learning methods are being developed to automatically construct patient-specific anatomical models as well as extract new insights from routine imaging (Reference Fathi, Perez-Raya, Baghaie, Berg, Janiga, Arzani and D'SouzaFathi et al. 2020; Reference Kong, Wilson and ShaddenKong, Wilson & Shadden 2021; Reference Arzani, Wang, Sacks and ShaddenArzani et al. 2022). There is also a need for real-time cardiovascular simulations, through reduced-order modelling and machine learning methods (Reference Pegolotti, Pfaller, Rubio, Ding, Brufau, Darve and MarsdenPegolotti et al. 2023).

Efforts to address these challenges will yield validated, accurate and robust digital twins that will accelerate the promise of personalized and non-invasive in silico medicine.