An isogeometric analysis framework for ventricular cardiac mechanics

The finite element method (FEM) is commonly used in computational cardiac simulations. For this method, a mesh is constructed to represent the geometry and, subsequently, to approximate the solution. To accurately capture curved geometrical features many elements may be required, possibly leading to unnecessarily large computation costs. Without loss of accuracy, a reduction in computation cost can be achieved by integrating geometry representation and solution approximation into a single framework using the Isogeometric Analysis (IGA) paradigm. In this study, we propose an IGA framework suitable for echocardiogram data of cardiac mechanics, where we show the advantageous properties of smooth splines through the development of a multi-patch anatomical model. A nonlinear cardiac model is discretized following the IGA paradigm, meaning that the spline geometry parametrization is directly used for the discretization of the physical fields. The IGA model is benchmarked with a state-of-the-art biomechanics model based on traditional FEM. For this benchmark, the hemodynamic response predicted by the high-fidelity FEM model is accurately captured by an IGA model with only 320 elements and 4,700 degrees of freedom. The study is concluded by a brief anatomy-variation analysis, which illustrates the geometric flexibility of the framework. The IGA framework can be used as a first step toward an efficient workflow for an improved understanding of, and clinical decision support for, the treatment of cardiac diseases like heart rhythm disorders

New hybrid quadrature schemes for weakly singular kernels applied to isogeometric boundary elements for 3D Stokes flow

This work proposes four novel hybrid quadrature schemes for the efficient and accurate evaluation of weakly singular boundary integrals (1/r kernel) on arbitrary smooth surfaces. Such integrals appear in boundary element analysis for several partial differential equations including the Stokes equation for viscous flow and the Helmholtz equation for acoustics. The proposed quadrature schemes apply a Duffy transform-based quadrature rule to surface elements containing the singularity and classical Gaussian quadrature to the remaining elements. Two of the four schemes additionally consider a special treatment for elements near to the singularity, where refined Gaussian quadrature and a new moment-fitting quadrature rule are used. The hybrid quadrature schemes are systematically studied on flat B-spline patches and on NURBS spheres considering two different sphere discretizations: An exact single-patch sphere with degenerate control points at the poles and an approximate discretization that consist of six patches with regular elements. The efficiency of the quadrature schemes is further demonstrated in boundary element analysis for Stokes flow, where steady problems with rotating and translating curved objects are investigated in convergence studies for both, mesh and quadrature refinement. Much higher convergence rates are observed for the proposed new schemes in comparison to classical schemes

Isogeometric Analysis for High Order Geometric Partial Differential Equations with Applications

In this thesis, we consider the numerical approximation of high order geometric Partial Differential Equations (PDEs). We first consider high order PDEs defined on surfaces in the 3D space that are represented by single-patch tensor product NURBS. Then, we spatially discretize the PDEs by means of NURBS-based Isogeometric Analysis (IGA) in the framework of the Galerkin method. With this aim, we consider the construction of periodic NURBS function spaces with high degree of global continuity, even on closed surfaces. As benchmark problems for the proposed discretization, we propose Laplace-Beltrami problems of the fourth and sixth orders, as well as the corresponding eigenvalue problems, and we analyze the impact of the continuity of the basis functions on the accuracy as well as on computational costs. The numerical solution of two high order phase field problems on both open and closed surfaces is also considered: the fourth order Cahn-Hilliard equation and the sixth order crystal equation, both discretized in time with the generalized-alpha method. We then consider the numerical approximation of geometric PDEs, derived, in particular, from the minimization of shape energy functionals by L^2-gradient flows. We analyze the mean curvature and the Willmore gradient flows, leading to second and fourth order PDEs, respectively. These nonlinear geometric PDEs are discretized in time with Backward Differentiation Formulas (BDF), with a semi-implicit formulation based on an extrapolation of the geometry, leading to a linear problem to be solved at each time step. Results about the numerical approximation of the two geometric flows on several geometries are analyzed. Then, we study how the proposed mathematical framework can be employed to numerically approximate the equilibrium shapes of lipid bilayer biomembranes, or vesicles, governed by the Canham-Helfrich curvature model. We propose two numerical schemes for enforcing the conservation of the area and volume of the vesicles, and report results on benchmark problems. Then, the approximation of the equilibrium shapes of biomembranes with different values of reduced volume is presented. Finally, we consider the dynamics of a vesicle, e.g. a red blood cell, immersed in a fluid, e.g. the plasma. In particular, we couple the curvature-driven model for the lipid membrane with the incompressible Navier-Stokes equations governing the fluid. We consider a segregated approach, with a formulation based on the Resistive Immersed Surface method applied to NURBS geometries. After analyzing benchmark fluid simulations with immersed NURBS objects, we report numerical results for the investigation of the dynamics of a vesicle under different flow conditions

A Spline-Based Partial Element Equivalent Circuit Method for Electrostatics

This contribution investigates the connection between Isogeometric Analysis (IgA) and the Partial Element Equivalent Circuit (PEEC) method for electrostatic problems. We demonstrate that using the spline-based geometry concepts from IgA allows for extracting circuit elements without an explicit meshing step. Moreover, the proposed IgA-PEEC method converges for complex geometries up to three times faster than the conventional PEEC approach and, in turn, it requires a significantly lower number of degrees of freedom to solve a problem with comparable accuracy. The resulting method is closely related to the isogeometric boundary element method. However, it uses lowest-order basis functions to allow for straightforward physical and circuit interpretations. The findings are validated by an analytical example with complex geometry, i.e., significant curvature, and by a realistic model of a surge arrester

Doctor of Philosophy

dissertationThe medial axis of an object is a shape descriptor that intuitively presents the morphology or structure of the object as well as intrinsic geometric properties of the object’s shape. These properties have made the medial axis a vital ingredient for shape analysis applications, and therefore the computation of which is a fundamental problem in computational geometry. This dissertation presents new methods for accurately computing the 2D medial axis of planar objects bounded by B-spline curves, and the 3D medial axis of objects bounded by B-spline surfaces. The proposed methods for the 3D case are the first techniques that automatically compute the complete medial axis along with its topological structure directly from smooth boundary representations. Our approach is based on the eikonal (grassfire) flow where the boundary is offset along the inward normal direction. As the boundary deforms, different regions start intersecting with each other to create the medial axis. In the generic situation, the (self-) intersection set is born at certain creation-type transition points, then grows and undergoes intermediate transitions at special isolated points, and finally ends at annihilation-type transition points. The intersection set evolves smoothly in between transition points. Our approach first computes and classifies all types of transition points. The medial axis is then computed as a time trace of the evolving intersection set of the boundary using theoretically derived evolution vector fields. This dynamic approach enables accurate tracking of elements of the medial axis as they evolve and thus also enables computation of topological structure of the solution. Accurate computation of geometry and topology of 3D medial axes enables a new graph-theoretic method for shape analysis of objects represented with B-spline surfaces. Structural components are computed via the cycle basis of the graph representing the 1-complex of a 3D medial axis. This enables medial axis based surface segmentation, and structure based surface region selection and modification. We also present a new approach for structural analysis of 3D objects based on scalar functions defined on their surfaces. This approach is enabled by accurate computation of geometry and structure of 2D medial axes of level sets of the scalar functions. Edge curves of the 3D medial axis correspond to a subset of ridges on the bounding surfaces. Ridges are extremal curves of principal curvatures on a surface indicating salient intrinsic features of its shape, and hence are of particular interest as tools for shape analysis. This dissertation presents a new algorithm for accurately extracting all ridges directly from B-spline surfaces. The proposed technique is also extended to accurately extract ridges from isosurfaces of volumetric data using smooth implicit B-spline representations. Accurate ridge curves enable new higher-order methods for surface analysis. We present a new definition of salient regions in order to capture geometrically significant surface regions in the neighborhood of ridges as well as to identify salient segments of ridges

Multiscale Modeling of Granular Materials

Granular materials have a “discrete” nature whose global mechanical behaviors are originated from the grain scale micromechanical mechanisms. The intriguing properties and non-trivial behaviors of a granular material pose formidable challenges to the multiscale modeling of these materials. Some of the key challenges include upscaling of coarse-scale continuum equation form fine-scale governing equations, calibrating material parameters at different scales, alleviating pathological mesh dependency in continuum models, and generating unit cells with versatile morphological details. This dissertation aims to addressing the aforementioned challenges and to investigate the mechanical behavior of granular materials through multiscale modeling. Firstly, a three-dimensional nonlocal multiscale discrete-continuum model is presented for modeling the mechanical behavior of granular materials. We establish an information-passing coupling scheme between DEM that explicitly replicates granular motion of individual particles and a finite element continuum model, which captures nonlocal overall response of the granular assemblies. Secondly, a new staggered multilevel material identification procedure is developed for phenomenological critical state plasticity models. The emphasis is placed on cases in which available experimental data and constraints are insufficient for calibration. The key idea is to create a secondary virtual experimental database from high-fidelity models, such as discrete element simulations, then merge both the actual experimental data and secondary database as an extended digital database to determine material parameters for the phenomenological macroscopic critical state plasticity model. This expansion of database provides additional constraints necessary for calibration of the phenomenological critical state plasticity models. Thirdly, a regularized phenomenological multiscale model is investigated, in which elastic properties are computed using direct homogenization and subsequently evolved using a simple three-parameter orthotropic continuum damage model. The salient feature of the model is a unified regularization framework based on the concept of effective softening strain. The unified regularization scheme is employed in the context of constitutive law rescaling and the staggered nonlocal approach to alleviate pathological mesh dependency. Lastly, a robust parametric model is presented for generating unit cells with randomly distributed inclusions. The proposed model is computationally efficient using a hierarchy of algorithms with increasing computational complexity, and is able to generate unit cells with different inclusion shapes

An immersed methodology for fluid-structure interaction using NURBS and T-splines: theory, algorithms, validation, and application to blood flow at small scales

[Abstract] Mesh-based immersed approaches shine in a variety of fluid-structure interaction (FSI) applications such as, e.g., simulations where the solid undergoes large displacements or rotations, particulate flow problems, and scenarios where the topology of the region occupied by the fluid varies in time. In this thesis, a new mesh-based immersed approach is proposed which is based on the use of di erent types of splines as basis functions. This approach is put forth for modeling and simulating di erent types of biological cells in blood flow at small scales. The specific contributions of this thesis are outlined as follows. Firstly, a hybrid variational-collocation immersed technique using nonuniform rational B-splines (NURBS) is presented. Newtonian viscous incompressible fluids and nonlinear hyperelastic incompressible solids are considered. Our formulation boils down to three coupled equations which are the linear momentum balance equation, the mass conservation equation, and the kinematic equation that relates the Lagrangian displacement with the Eulerian velocity. The latter is discretized in strong form using isogeometric collocation and the other two equations are discretized using the variational multiscale (VMS) paradigm. As usual in immersed FSI approaches, we define a background mesh on the whole computational domain and a Lagrangian mesh tailored to the region occupied by each solid. Besides of using NURBS for creating these meshes, the data transfer between the background mesh and the Lagrangian meshes is carried out using NURBS functions in such a way that no interpolation or projection is needed, thus avoiding the errors associated with these procedures. Regarding the time discretization, the generalized- method is used which leads to a fully-implicit and second-order accurate method. The methodology is validated in two- and three-dimensional settings comparing the terminal velocity of free-falling bulky solids obtained in our simulations with its theoretical value. Secondly, we extend our algorithms in order to use analysis-suitable T-splines (ASTS) as basis functions instead of NURBS. This required to develop isogeometric collocation methods for ASTS which was an open problem. The data transfer between meshes changes significantly from NURBS to ASTS due to the fact that their geometrical mappings are local to patches and elements, respectively. ASTS possess two main advantages with respect to NURBS: (1) ASTS support local h-refinement and (2) ASTS are unstructured. The ASTSbased method is validated solving again the aforementioned benchmark problems and showing the potential of ASTS to decrease the amount of elements needed, thus enhancing the e ciency of the method. Thirdly, capsules, modeled as solid-shell NURBS elements, are proposed as numerical proxies for representing red blood cells (RBCs). The dynamics of capsules are able to reproduce the main motions and shapes observed in experiments with RBCs in both shear and parabolic flows. Hemorheological properties as the Fåhræus and Fåhræus-Lindqvist e ects are captured in our simulations. In order to obtain the aforementioned results, it is essential to adequately satisfy the incompressibility constraint close to the fluid-solid interface, which is an arduous task in immersed approaches for fluid-structure interaction. Finally, compound capsules are presented as numerical proxies for cells with nucleus such as, e.g., white blood cells (WBCs) and circulating tumor cells (CTCs). The dynamics of hyperelastic compound capsules in shear flow are studied in both two- and three-dimensional settings. Moreover, we simulate how CTCs manage to pass through channel narrowings, which is an interesting characteristic of CTCs since it is used in experiments to sort CTCs from blood samples

Methods based on B-splines for model representation, numerical analysis and image registration

The thesis consists of inter-connected parts for modeling and analysis using newly developed isogeometric methods. The main parts are reproducing kernel triangular B-splines, extended isogeometric analysis for solving weakly discontinuous problems, collocation methods using superconvergent points, and B-spline basis in image registration applications. Each topic is oriented towards application of isogeometric analysis basis functions to ease the process of integrating the modeling and analysis phases of simulation. First, we develop reproducing a kernel triangular B-spline-based FEM for solving PDEs. We review the triangular B-splines and their properties. By definition, the triangular basis function is very flexible in modeling complicated domains. However, instability results when it is applied for analysis. We modify the triangular B-spline by a reproducing kernel technique, calculating a correction term for the triangular kernel function from the chosen surrounding basis. The improved triangular basis is capable to obtain the results with higher accuracy and almost optimal convergence rates. Second, we propose an extended isogeometric analysis for dealing with weakly discontinuous problems such as material interfaces. The original IGA is combined with XFEM-like enrichments which are continuous functions themselves but with discontinuous derivatives. Consequently, the resulting solution space can approximate solutions with weak discontinuities. The method is also applied to curved material interfaces, where the inverse mapping and the curved triangular elements are considered. Third, we develop an IGA collocation method using superconvergent points. The collocation methods are efficient because no numerical integration is needed. In particular when higher polynomial basis applied, the method has a lower computational cost than Galerkin methods. However, the positions of the collocation points are crucial for the accuracy of the method, as they affect the convergent rate significantly. The proposed IGA collocation method uses superconvergent points instead of the traditional Greville abscissae points. The numerical results show the proposed method can have better accuracy and optimal convergence rates, while the traditional IGA collocation has optimal convergence only for even polynomial degrees. Lastly, we propose a novel dynamic multilevel technique for handling image registration. It is application of the B-spline functions in image processing. The procedure considered aims to align a target image from a reference image by a spatial transformation. The method starts with an energy function which is the same as a FEM-based image registration. However, we simplify the solving procedure, working on the energy function directly. We dynamically solve for control points which are coefficients of B-spline basis functions. The new approach is more simple and fast. Moreover, it is also enhanced by a multilevel technique in order to prevent instabilities. The numerical testing consists of two artificial images, four real bio-medical MRI brain and CT heart images, and they show our registration method is accurate, fast and efficient, especially for large deformation problems

Discrete Element Modeling of the Grading- and Shape-Dependent Behavior of Granular Materials

Granular materials, such as sand, biomass particles, and pharmaceutical pills, are widespread in nature, industrial systems, and our daily life. Fundamentally, the bulk mechanical behavior of such materials is governed by the physical and morphological features of and the interactions among constituent particles at the microscopic scale. From a modeling standpoint, the particle-based discrete element method (DEM) has emerged as the most prevalent numerical tool to model and study the behavior of granular materials and the systems they form. A critical step towards an accurate and predictive DEM model is to incorporate those physical and morphological features (e.g., particle size, shape, and deformability) pertaining to the constituent particles. The main objective of this dissertation is to approach an accurate characterization and modeling of the grading- and shape-dependent behavior of granular materials by developing DEM models that incorporate realistic physical and morphological features of granular particles. Revolving around this objective, three studies are presented: image-based particle reconstruction and morphology characterization, grading and shape-dependent shearing behavior of rigid-particle systems, and granular flow of deformable irregular particles. The first study presents a machine learning and level-set based framework to re- construct granular particles and to characterize particle morphology from X-ray computed tomography (X-ray CT) imaging of realistic granular materials. Images containing detailed microstructure information of a granular material are obtained using the X-ray CT tech- nique. Approaches such as the watershed method in two dimensions (2D) and the combined machine learning and level set method in three dimensions (3D) are then utilized and implemented to segment X-ray CT images and to numerically reconstruct individual particles in the granular material. Based on the realistic particle shapes, particle morphology is characterized by descriptors including aspect ratio, roundness, circularity (2D) or sphericity (3D). The particle shapes or morphology provide important constraints to develop DEM models with particle physical and morphological features conforming to the specific granular material of interest. In the second study, DEM models incorporated with realistic particle sizes and shapes are developed and applied to study the shearing behavior of sandy soils. The particle sizes and shapes are obtained from realistic samples of JSC-1A Martian regolith simulant. Irregular-shape particles are represented by rigid clumps based on the domain overlapping filling method. The effects of particle shape irregularity on the shearing behavior of granular materials are investigated through direct shear tests, along with the comparisons from spherical particles with or without rolling resistance. The micro-mechanisms of shape irregularity contributing to the shear resistance are identified. The last study investigates the effects of particle deformability (e.g., compression, deflection or torsion), together with particle sizes and shapes, on the granular flow of flexible granular materials. A bonded-sphere DEM model is implemented with the capability of embodying various particle sizes and irregular shapes, as well as capturing particle deformability. This approach is then applied to simulate and study the behavior of flexible granular materials in cyclic compression and hopper flow tests. The effects of particle size, shape and deformability on the bulk mechanical behavior are investigated on the basis of the DEM simulation results. The importance of particle deformability to the DEM simulations of flexible granular materials is demonstrated