Faculty Bibliography

In this work, we extend our previous esophageal transport model using an immersed boundary (IB) method with discrete fiber-based structural model, to one using a continuum mechanics-based model that is approximated based on finite elements (IB-FE). To deal with the leakage of flow when the Lagrangian mesh becomes coarser than the fluid mesh, we employ adaptive interaction quadrature points to deal with Lagrangian-Eulerian interaction equations based on a previous work (Griffith and Luo [1]). In particular, we introduce a new anisotropic adaptive interaction quadrature rule. The new rule permits us to vary the interaction quadrature points not only at each time-step and element but also at different orientations per element. This helps to avoid the leakage issue without sacrificing the computational efficiency and accuracy in dealing with the interaction equations. For the material model, we extend our previous fiber-based model to a continuum-based model. We present formulations for general fiber-reinforced material models in the IB-FE framework. The new material model can handle non-linear elasticity and fiber-matrix interactions, and thus permits us to consider more realistic material behavior of biological tissues. To validate our method, we first study a case in which a three-dimensional short tube is dilated. Results on the pressure-displacement relationship and the stress distribution matches very well with those obtained from the implicit FE method. We remark that in our IB-FE case, the three-dimensional tube undergoes a very large deformation and the Lagrangian mesh-size becomes about 6 times of Eulerian mesh-size in the circumferential orientation. To validate the performance of the method in handling fiber-matrix material models, we perform a second study on dilating a long fiber-reinforced tube. Errors are small when we compare numerical solutions with analytical solutions. The technique is then applied to the problem of esophageal transport. We use two fiber-reinforced models for the esophageal tissue: a bi-linear model and an exponential model. We present three cases on esophageal transport that differ in the material model and the muscle fiber architecture. The overall transport features are consistent with those observed from the previous model. We remark that the continuum-based model can handle more realistic and complicated material behavior. This is demonstrated in our third case where a spatially varying fiber architecture is included based on experimental study. We find that this unique muscle fiber architecture could generate a so-called pressure transition zone, which is a luminal pressure pattern that is of clinical interest. This suggests an important role of muscle fiber architecture in esophageal transport.

In standard models of cardiac electrophysiology, including the bidomain and monodomain models, local perturbations can propagate at infinite speed. We address this unrealistic property by developing a hyperbolic bidomain model that is based on a generalization of Ohm's law with a Cattaneo-type model for the fluxes. Further, we obtain a hyperbolic monodomain model in the case that the intracellular and extracellular conductivity tensors have the same anisotropy ratio. In one spatial dimension, the hyperbolic monodomain model is equivalent to a cable model that includes axial inductances, and the relaxation times of the Cattaneo fluxes are strictly related to these inductances. A purely linear analysis shows that the inductances are negligible, but models of cardiac electrophysiology are highly nonlinear, and linear predictions may not capture the fully nonlinear dynamics. In fact, contrary to the linear analysis, we show that for simple nonlinear ionic models, an increase in conduction velocity is obtained for small and moderate values of the relaxation time. A similar behavior is also demonstrated with biophysically detailed ionic models. Using the Fenton-Karma model along with a low-order finite element spatial discretization, we numerically analyze differences between the standard monodomain model and the hyperbolic monodomain model. In a simple benchmark test, we show that the propagation of the action potential is strongly influenced by the alignment of the fibers with respect to the mesh in both the parabolic and hyperbolic models when using relatively coarse spatial discretizations. Accurate predictions of the conduction velocity require computational mesh spacings on the order of a single cardiac cell. We also compare the two formulations in the case of spiral break up and atrial fibrillation in an anatomically detailed model of the left atrium, and we examine the effect of intracellular and extracellular inductances on the virtual electrode phenomenon.

Understanding the interaction between the valves and walls of the heart is important in assessing and subsequently treating heart dysfunction. This study presents an integrated model of the mitral valve (MV) coupled to the left ventricle (LV), with the geometry derived from in vivo clinical magnetic resonance images. Numerical simulations using this coupled MV-LV model are developed using an immersed boundary/finite element method. The model incorporates detailed valvular features, left ventricular contraction, nonlinear soft tissue mechanics, and fluid-mediated interactions between the MV and LV wall. We use the model to simulate cardiac function from diastole to systole. Numerically predicted LV pump function agrees well with in vivo data of the imaged healthy volunteer, including the peak aortic flow rate, the systolic ejection duration, and the LV ejection fraction. In vivo MV dynamics are qualitatively captured. We further demonstrate that the diastolic filling pressure increases significantly with impaired myocardial active relaxation to maintain a normal cardiac output. This is consistent with clinical observations. The coupled model has the potential to advance our fundamental knowledge of mechanisms underlying MV-LV interaction, and help in risk stratification and optimisation of therapies for heart diseases.

Each year, approximately 300,000 heart valve repair or replacement procedures are performed worldwide, including approximately 70,000 aortic valve replacement surgeries in the United States alone. Computational platforms for simulating cardiovascular devices such as prosthetic heart valves promise to improve device design and assist in treatment planning, including patient-specific device selection. This paper describes progress in constructing anatomically and physiologically realistic immersed boundary (IB) models of the dynamics of the aortic root and ascending aorta. This work builds on earlier IB models of fluid-structure interaction (FSI) in the aortic root, which previously achieved realistic hemodynamics over multiple cardiac cycles, but which also were limited to simplified aortic geometries and idealized descriptions of the biomechanics of the aortic valve cusps. By contrast, the model described herein uses an anatomical geometry reconstructed from patient-specific computed tomography angiography (CTA) data, and employs a description of the elasticity of the aortic valve leaflets based on a fiber-reinforced constitutive model fit to experimental tensile test data. The resulting model generates physiological pressures in both systole and diastole, and yields realistic cardiac output and stroke volume at physiological Reynolds numbers. Contact between the valve leaflets during diastole is handled automatically by the IB method, yielding a fully competent valve model that supports a physiological diastolic pressure load without regurgitation. Numerical tests show that the model is able to resolve the leaflet biomechanics in diastole and early systole at practical grid spacings. The model is also used to examine differences in the mechanics and fluid dynamics yielded by fresh valve leaflets and glutaraldehyde-fixed leaflets similar to those used in bioprosthetic heart valves. Although there are large differences in the leaflet deformations during diastole, the differences in the open configurations of the valve models are relatively small, and nearly identical hemodynamics are obtained in all cases considered.

The smallest flying insects commonly possess wings with long bristles. Little quantitative information is available on the morphology of these bristles, and their functional importance remains a mystery. In this study, we (1) collected morphological data on the bristles of 23 species of Mymaridae by analyzing high-resolution photographs and (2) used the immersed boundary method to determine via numerical simulation whether bristled wings reduced the force required to fling the wings apart while still maintaining lift. The effects of Reynolds number, angle of attack, bristle spacing and wing-wing interactions were investigated. In the morphological study, we found that as the body length of Mymaridae decreases, the diameter and gap between bristles decreases and the percentage of the wing area covered by bristles increases. In the numerical study, we found that a bristled wing experiences less force than a solid wing. The decrease in force with increasing gap to diameter ratio is greater at higher angles of attack than at lower angles of attack, suggesting that bristled wings may act more like solid wings at lower angles of attack than they do at higher angles of attack. In wing-wing interactions, bristled wings significantly decrease the drag required to fling two wings apart compared with solid wings, especially at lower Reynolds numbers. These results support the idea that bristles may offer an aerodynamic benefit during clap and fling in tiny insects.

Rupture of aneurysms and acute dissection of the thoracic aorta are life-threatening events which affect tens of thousands of people per year. The underlying mechanisms remain unclear and the aortic wall is known to lose its structural integrity, which in turn affects its mechanical response to the loading conditions. Hence, research on such aortic diseases is an important area in biomechanics. The present study investigates the mechanical properties of aneurysmatic and dissected human thoracic aortas via triaxial shear and uniaxial tensile testing with a focus on the former. In particular, ultimate stress values from triaxial shear tests in different orientations regarding the aortas orthotropic microstructure, and from uniaxial tensile tests in radial, circumferential and longitudinal directions were determined. In total, 16 human thoracic aortas were investigated from which it is evident that the aortic media has much stronger resistance to rupture under 'out-of-plane' than under 'in-plane' shear loadings. Under different shear loadings the aortic tissues revealed anisotropic failure properties with higher ultimate shear stresses and amounts of shear in the longitudinal than in the circumferential direction. Furthermore, the aortic media decreased its tensile strength as follows: circumferential direction >longitudinaldirection> radial direction. Anisotropic and nonlinear tissue properties are apparent from the experimental data. The results clearly showed interspecimen differences influenced by the anamnesis of the donors such as aortic diseases or connective tissue disorders, e.g., dissected specimens exhibited on average a markedly lower mechanical strength than aneurysmatic specimens. The rupture data based on the combination of triaxial shear and uniaxial extension testing are unique and build a good basis for developing a 3D failure criterion of diseased human thoracic aortic media. This is a step forward to more realistic modeling of mechanically induced tissue failure i.e. rupture of aneurysms or progression of aortic dissections.

It has long been recognized that aortic root elasticity helps to ensure efficient aortic valve closure, but our understanding of the functional importance of the elasticity and geometry of the aortic root continues to e-volve as increasingly detailed in vivo imaging data become available. Herein, we describe a fluid-structure interaction model of the aortic root, including the aortic valve leaflets, the sinsuses of Valsalva, the aortic annulus, and the sinotubular junction, that employs a version of Peskin's immersed boundary (IB) method with a finite element (FE) description of the structural elasticity. As in earlier work, we use a fiber-based model of the valve leaflets, but this study extends earlier IB models of the aortic root by employing an incompressible hyperelastic model of the mechanics of the sinuses and ascending aorta using a constitutive law fit to experimental data from human aortic root tissue. In vivo pressure loading is accounted for by a backward displacement method that determines the unloaded configurations of the root model. Our model yields realistic cardiac output at physiological pressures, with low transvalvular pressure differences during forward flow, minimal regurgitation during valve closure, and realistic pressure loads when the valve is closed during diastole. Further, results from high-resolution computations indicate that although the detailed leaflet and root kinematics show some grid sensitivity, our IB model of the aortic root nonetheless produces essentially grid-converged flow rates and pressures at practical grid spacings for the high-Reynolds number flows of the aortic root. These results thereby clarify minimum grid resolutions required by such models when used as stand-alone models of the aortic valve as well as when used to provide models of the outflow valves in models of left ventricular fluid dynamics.

We develop a rigid multiblob method for numerically solving the mobility problem for suspensions of passive and active rigid particles of complex shape in Stokes flow in unconfined, partially confined, and fully confined geometries. As in a number of existing methods, we discretize rigid bodies using a collection of minimally resolved spherical blobs constrained to move as a rigid body, to arrive at a potentially large linear system of equations for the unknown Lagrange multipliers and rigid-body motions. Here we develop a block-diagonal preconditioner for this linear system and show that a standard Krylov solver converges in a modest number of iterations that is essentially independent of the number of particles. Key to the efficiency of the method is a technique for fast computation of the product of the blob-blob mobility matrix and a vector. For unbounded suspensions, we rely on existing analytical expressions for the Rotne-Prager-Yamakawa tensor combined with a fast multipole method (FMM) to obtain linear scaling in the number of particles. For suspensions sedimented against a single no-slip boundary, we use a direct summation on a graphical processing unit (GPU), which gives quadratic asymptotic scaling with the number of particles. For fully confined domains, such as periodic suspensions or suspensions confined in slit and square channels, we extend a recently developed rigid-body immersed boundary method by B. Kallemov, A. P. S. Bhalla, B. E. Griffith, and A. Donev (Commun. Appl. Math. Comput. Sci. 11 (2016), no. 1, 79-141) to suspensions of freely moving passive or active rigid particles at zero Reynolds number. We demonstrate that the iterative solver for the coupled fluid and rigid-body equations converges in a bounded number of iterations regardless of the system size. In our approach, each iteration only requires a few cycles of a geometric multigrid solver for the Poisson equation, and an application of the block-diagonal preconditioner, leading to linear scaling with the number of particles. We optimize a number of parameters in the iterative solvers and apply our method to a variety of benchmark problems to carefully assess the accuracy of the rigid multiblob approach as a function of the resolution. We also model the dynamics of colloidal particles studied in recent experiments, such as passive boomerangs in a slit channel, as well as a pair of non-Brownian active nanorods sedimented against a wall.

We develop an immersed boundary (IB) method for modeling flows around fixed or moving rigid bodies that is suitable for a broad range of Reynolds numbers, including steady Stokes flow. The spatio-temporal discretization of the fluid equations is based on a standard staggered-grid approach. Fluid-body interaction is handled using Peskin's IB method; however, unlike existing IB approaches to such problems, we do not rely on penalty or fractional-step formulations. Instead, we use an unsplit scheme that ensures the no-slip constraint is enforced exactly in terms of the Lagrangian velocity field evaluated at the IB markers. Fractionalstep approaches, by contrast, can impose such constraints only approximately, which can lead to penetration of the flow into the body, and are inconsistent for steady Stokes flow. Imposing no-slip constraints exactly requires the solution of a large linear system that includes the fluid velocity and pressure as well as Lagrange multiplier forces that impose the motion of the body. The principal contribution of this paper is that it develops an efficient preconditioner for this exactly constrained IB formulation which is based on an analytical approximation to the Schur complement. This approach is enabled by the near translational and rotational invariance of Peskin's IB method. We demonstrate that only a few cycles of a geometric multigrid method for the fluid equations are required in each application of the preconditioner, and we demonstrate robust convergence of the overall Krylov solver despite the approximations made in the preconditioner. We empirically observe that to control the condition number of the coupled linear system while also keeping the rigid structure impermeable to fluid, we need to place the immersed boundary markers at a distance of about two grid spacings, which is significantly larger from what has been recommended in the literature for elastic bodies. We demonstrate the advantage of our monolithic solver over split solvers by computing the steady state flow through a two-dimensional nozzle at several Reynolds numbers. We apply the method to a number of benchmark problems at zero and finite Reynolds numbers, and we demonstrate first-order convergence of the method to several analytical solutions and benchmark computations.

BACKGROUND AND PURPOSE: Texture analysis has been applied to medical images to assist in tumor tissue classification and characterization. In this study, we obtained textural features from parametric (relative CBV) maps of dynamic susceptibility contrast-enhanced MR images in glioblastoma and assessed their relationship with patient survival. MATERIALS AND METHODS: MR perfusion data of 24 patients with glioblastoma from The Cancer Genome Atlas were analyzed in this study. One- and 2D texture feature ratios and kinetic textural features based on relative CBV values in the contrast-enhancing and nonenhancing lesions of the tumor were obtained. Receiver operating characteristic, Kaplan-Meier, and multivariate Cox proportional hazards regression analyses were used to assess the relationship between texture feature ratios and overall survival. RESULTS: Several feature ratios are capable of stratifying survival in a statistically significant manner. These feature ratios correspond to homogeneity (P = .008, based on the log-rank test), angular second moment (P = .003), inverse difference moment (P = .013), and entropy (P = .008). Multivariate Cox proportional hazards regression analysis showed that homogeneity, angular second moment, inverse difference moment, and entropy from the contrast-enhancing lesion were significantly associated with overall survival. For the nonenhancing lesion, skewness and variance ratios of relative CBV texture were associated with overall survival in a statistically significant manner. For the kinetic texture analysis, the Haralick correlation feature showed a P value close to .05. CONCLUSIONS: Our study revealed that texture feature ratios from contrast-enhancing and nonenhancing lesions and kinetic texture analysis obtained from perfusion parametric maps provide useful information for predicting survival in patients with glioblastoma.