SIMULIA LS Med Devices Ebook
2015-10-19
: Ensight Simulia-Ls-Meddevices-Ebook SIMULIA-LS-MedDevices-ebook AbaqusRUM_2015 CEI Houston
Open the PDF directly: View PDF .
Page Count: 34
Download | |
Open PDF In Browser | View PDF |
Optimize Medical Device Reliability through Multiphysics Simulation Life Sciences Industrial Applications of Multiphysics Simulation to Improve Device Performance, Reduce Time & Costs, and Support Regulatory Approval Content Simulation Now Recognized by FDA as Essential to Medical Device Evaluation 3 SIMULIA Community News The Living Heart Project: A robust and integrative simulator for human heart function 6 Modeling Hemodynamics with Abaqus/CFD Steady State Solver: FDA Benchmark Nozzle Model 7 Brian Baillargeon, Nuno Rebelo, and David D. Fox (Dassault Systèmes) Robert L. Taylor (University of California at Berkeley, Berkeley, California) Ellen Kuhl (Stanford University, Stanford, California) European Journal of Mechanics A/Solids Santanu Chandra, PhD and Richard Swift, PhD, PE (MED Institute Inc.) Ramji Kamakoti, PhD (Dassault Systèmes) SIMULIA Community Conference, Providence, Rhode Island Novo Nordisk Makes Designing Injection Pens a Snap (Fit) 16 Optimization of Surgical Positioning in Total Hip Replacement 18 Computational Study of Cortical Bone Screw Pullout using the eXtended Finite Element Method (XFEM) 19 SIMULIA Spearheads The Living Heart Project 25 Abdominal Aorta Blood Flow Analysis with Abaqus/CFD 26 Fluid-Structure Interaction (FSI) of a Prosthetic Aortic Valve using Abaqus/Explicit Smoothed Particle Hydrodynamics 31 Leading innovators in diabetes care use realistic simulation to improve product integrity from design to manufacture SIMULIA Community News Jacob Elkins, John Callaghan, Douglas Pedersen, and Thomas Brown (University of Iowa, Iowa City, Iowa) SIMULIA Community News Emer M. Feerick and Patrick McGarry (National University of Ireland, Galway, Ireland) SIMULIA Community Conference, Providence, Rhode Island Bringing the Medical Community Together for Improved Patient Care SIMULIA Community News Sharat Prasad (Dassault Systèmes) Abaqus Technology Brief Cheryl Liu (Dassault Systèmes) Abaqus Technology Brief www.3ds.com/simulia 2 Compliance Simulation Now Recognized by FDA as Essential to Medical Device Evaluation FDA sees increasing numbers of applications that include simulation As Life Sciences engineers embrace simulation, they are achieving increasingly accurate levels of precision when evaluating device function, including the ability to evaluate aspects of device performance not possible with bench tests alone. As a result, the Food and Drug Administration’s (FDA) Center for Devices and Radiological Health (CDRH) is seeing a growing number of submissions for medical devices that include a simulation-data component. Cheryl Liu, Senior Technical Marketing Specialist, Life Sciences Industry, SIMULIA One of the toughest design engineering challenges is making a medical device that works flawlessly with the human body. The unique anatomy and physiology of every patient create physical complexities, and ever-shifting functional parameters, that must be thoroughly accounted for when producing a therapeutic product that may need to last a lifetime. The CDRH is responsible for regulating firms that manufacture, repackage, relabel, and/or import medical devices sold in the U.S. The submissions for these therapeutic devices typically contain data from four types of evaluation models—animal, bench, computational, and human—to demonstrate a reasonable assurance of safety and effectiveness. When a company submits simulation metrics that supplement bench testing, this can help promote approval by demonstrating both the integrity of the proposed device and the required realistic device failure analysis. As the ultimate safety-and-effectiveness regulatory body between medical device manufacturers and patients, the FDA recognizes the value of such advancing technologies—and its own need to stay abreast of them—and has now begun actively encouraging the use of simulation in device evaluation. Domestic inpatient procedures involving medical devices— stents, heart valves, dental implants, spine and joint implants, surgical tools, blood pumps, endovascular grafts, drug-eluting devices, and more—totaled 46 million in the U.S. alone in 2006, according to the Centers for Disease Control and Prevention (CDC). It’s a global market that is growing along with aging populations everywhere. Computer simulation, already widely accepted in many industries, is increasingly being viewed as an important tool by medical device companies and their designers. It helps them visualize what they cannot see, explore the design space more fully, refine their ideas faster and more accurately—and reduce expensive prototyping and testing. However, the FDA has also put the industry on notice that verification and validation (V&V) must go hand-in-hand with the use of simulation in applications. The CDRH is looking to quantify when a computational model is credible enough, and whether its intended purpose is appropriate for a regulatory submission. Unclear reporting standards, insufficient data about geometries and boundary conditions, lack of validation metrics, incomplete understanding of physiological loads in the body, and variations in patient populations—any and all of these uncertainties can impact the relevance of simulation outputs. Solid mechanics simulations can help determine proper implant size, evaluate manufacturing tolerances, compare design geometries, or consider next-gen devices. Fluid dynamics can be employed to identify high-shear stresses on blood vessels, regions of low flow, and potential for blood damage. And simulation-based product development processes can be linked in automated workflows, optimizing huge quantities of design data to provide exquisitely fine-tuned results that are of particular value for creating patient-specific medical devices. a) b) c) Figure 1. This example of modeling and simulation of a medical device shows an aortic valve geometry (a), a model of the effect of blood flow on the valve in a blood vessel (b), and an Abaqus finite element analysis (FEA) of the stress on the valve leaflets during the diastolic phase (c). This work was performed by SIMULIA in conjunction with the FDA’s Center for Devices and Radiological Health (CDRH). www.3ds.com/simulia Previous 3 Next Contents Compliance SIMULIA contributes to advancement of knowledge The Virtual Physiological Patient Noticing that a significant proportion of the applications they have seen in recent years have included simulations with Abaqus finite element analysis (FEA), the CDRH reached out to us in 2010 for support in developing their own internal framework, and in-house expertise, for validating and regulating industry-submitted simulations. DBS Massachusetts General Hospital ODEC SIMULIA presented at the FDA’s 3rd workshop on Computational Modeling of Medical Devices the same year. We continue to deliver on-site training courses to FDA reviewers about best practices in modeling and simulation and to partner with the FDA on aortic valve model development (see Figure 1). The FDA has also presented at our SIMULIA Community Conference and Regional User Meetings. Realizing the importance of model V&V in early 2011, ASME and FDA launched the V&V 40 subcommittee to develop V&V guidelines for the medical device industry specifically; we are actively participating, along with others in the industry and software communities. Johns Hopkins University CABAS IMCAFS Georgia Institute of Technology Penn Medicine As one outcome of these efforts, the FDA will publish a guidance document titled “Reporting Computational Modeling Studies in Medical Device Regulatory Submissions” in 2013. Appendices will cover fluid and mass transport, solid mechanics, electromagnetism, control loops, thermal transport, and ultrasound. Publication date updates can be found on the CDRH website at www.fda.gov/MedicalDevices/ DeviceRegulationandGuidance/GuidanceDocuments/ ucm371016.htm BioEVAL University of Pittsburg CHAP Dartmouth Hitchcock M2S SENTARA Penn Medicine HIPS Cleveland Clinic Massachusetts General Hospital BARD Boston Scientific Abbott ndc SIMULIA As knowledge about the importance of simulation grows, another priority for the FDA is the creation of a publicly available ‘Virtual Physiological Patient’ of human body computer models in different disease states (see Figure 2). This is not intended to be a single model encompassing every function and disease at once. Rather, the project will comprise a library of verified and validated submodels and data based on the combined expertise of those groups in the relevant disciplines, i.e., cardiology, orthopedics, software, and so forth. VFAM IT’IS Foundation The goal of the Virtual Physiological Patient project is a shared point of reference that will improve understanding of model attributes and limitations, and provide discrete models, data, and simulations validated for regulatory evaluation. Peer review by experts in academia, government, and industry will ensure robust V&V and provide periodic assessment. SIMULIA is contributing expertise to a group that is developing a computational model for the evaluation of a diseased femoral artery for stent evaluation. Previous University of Colorado School of Medicine The University of Iowa ASPECT The ‘Virtual Patient’ idea is born www.3ds.com/simulia CPI/CFD Cleveland Clinic Mississippi State University Penn State RIT University of Pittsburg Figure 2. Broad cross-industry collaboration between medical device manufacturers, academia, and software companies is being harnessed for the FDA’s Virtual Physiological Patient project. 4 Next Contents Compliance Newly launched public-private partnership benefits all parties the function of the living body is critical to every medical-device developer, and sharing the data that lie at the core of that understanding can be accomplished without infringing on any one company’s patents. Concurrent with the development of the Virtual Physiological Patient concept, the FDA is reaching outward to device manufacturers, software providers, and medical professionals to form a Regulatory Science Public-Private Partnership. Launched in December of 2012, the partnership is called the Medical Device Innovation Consortium (MDIC). The FDA views modeling and simulation as incentives to innovation that can reduce the time and cost of device design, assessment, and manufacturing. It is in all our interests—the medical device industry, the regulatory agency, and software companies—to collaborate to ensure that the power of simulation is increasingly utilized to solve the wide range of challenges in medical device development. We can all agree that the ultimate goal is the safety and effectiveness of medical devices for every physician who uses them, and every patient who needs them. The idea is to create an opportunity for information gathering in a pre-competitive state, i.e., not device-specific, but diseasespecific. For example, if the heart valve community were interested in a comprehensive evaluation of the structure and function of heart valves, costs could be minimized through nonprofit group funding and participation in the development of tools and resources for modeling and simulating of a range of valves. All results would be shared. End-stage renal disease is another area recently identified by the FDA as a priority. Industry forums on this topic are already underway. Special thanks to Dr. Tina Morrison, Dr. Nandini Duraiswamy, and Dr. Donna Lochner of the FDA for their assistance in preparing this article. Read more about how the FDA is promoting innovation in “High Stakes Balancing Act” in Compass magazine— www.compassmag.3ds.com. The medical device industry can only benefit from such endeavors. Individual device design copyrights certainly need to be protected, but the tradition of publishing evidencebased research results in order to move the entire body of medical knowledge forward has resonated in the life sciences throughout the history of medicine. A deep understanding of www.3ds.com/simulia Previous For More Information www.fda.gov www.deviceconsortium.org 5 Next Contents Cardiovascular The Living Heart Project: A robust and integrative simulator for human heart function Brian Baillargeon a, Nuno Rebelo a, David D. Fox b, Robert L. Taylor c, Ellen Kuhl d,* Dassault Systèmes Simulia Corporation, Fremont, CA 94538, USA Dassault Systèmes Simulia Corporation, Providence, RI 02909, USA c Department of Civil and Environmental Engineering, University of California at Berkeley, Berkeley, CA 94720, USA d Departments of Mechanical Engineering, Bioengineering, and Cardiothoracic Surgery, Stanford University, Stanford, CA 94305, USA a b ABSTRACT The heart is not only our most vital, but also our most complex organ: Precisely controlled by the interplay of electrical and mechanical fields, it consists of four chambers and four valves, which act in concert to regulate its filling, ejection, and overall pump function. While numerous computational models exist to study either the electrical or the mechanical response of its individual chambers, the integrative electro-mechanical response of the whole heart remains poorly understood. Here we present a proofof-concept simulator for a four-chamber human heart model created from computer topography and magnetic resonance images. We illustrate the governing equations of excitation-contraction coupling and discretize them using a single, unified finite element environment. To illustrate the basic features of our model, we visualize the electrical potential and the mechanical deformation across the human heart throughout its cardiac cycle. To compare our simulation against common metrics of cardiac function, we extract the pressure-volume relationship and show that it agrees well with clinical observations. Our prototype model allows us to explore and understand the key features, physics, and technologies to create an integrative, predictive model of the living human heart. Ultimately, our simulator will open opportunities to probe landscapes of clinical parameters, and guide device design and treatment planning in cardiac diseases such as stenosis, regurgitation, or prolapse of the aortic, pulmonary, tricuspid, or mitral valve. Read the entire paper at http://dx.doi.org/10.1016/j.euromechsol.2014.04.001. Mechanical displacement © 2014 The Authors. Published by Elsevier Masson SAS. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/3.0/). * Corresponding author. E-mail address: ekuhl@stanford.edu (E. Kuhl). www.3ds.com/simulia Previous 6 Next Contents Cardiovascular Modeling Hemodynamics with Abaqus/CFD Steady State Solver: FDA Benchmark Nozzle Model Santanu Chandra, PhD and Richard Swift, PhD, PE (MED Institute Inc.) Ramji Kamakoti, PhD (Dassault Systèmes) Abstract: Understanding the blood flow dynamics (hemodynamics) and the fluid forces exerted on the blood by implantable medical devices and predicting blood damage is an intricate part of interventional medical device design. Computational techniques such as Computational Fluid Dynamics (CFD) are increasingly being used as a tool for describing complex hemodynamics and calculating the fluid forces such as pressure and shear stress. However, this technique is challenged by the lack of standardized methods for validation and verification of the results. In an effort to standardize the process, the FDA made an initiative in 2008 to engage the Medical Device/CFD Community worldwide to participate in an open challenge that involves CFD computation of blood flow and its experimental validation with Particle Image Velocimetry (PIV) on a benchmark nozzle model. The goal of the CFD phase of the study was to understand the limitations of CFD, understand the variability that arises due to choice of software or solver, and the influence of user expertise level and diverse modeling and meshing techniques. This information would then be employed to identify best practices and define the critical techniques necessary to achieve credible results. Though a variety of software packages have been used by 28 different groups from 6 countries, Abaqus/CFD has not been publicly applied in this challenge. In this study we aim to assess the performance of Abaqus/CFD in modeling hemodynamics using the FDA Benchmark Nozzle model. The Nozzle model consists of a tube with a straight section, followed by a conical section, a section with reduced tube diameter and a section with a sudden expansion in tube diameter. The device was designed to include accelerating flow, decelerating flow, variations in shear stress and velocities, and recirculating flow, all of which may be present in a medical device and relate to blood damage. Five different flow rates corresponding to fully laminar, transitional and fully turbulent flow were modeled with laminar and turbulent solvers. In this study we have used pure hex mesh as well as the steady state solver introduced recently to solve both laminar and turbulent flow problems and compared the results with published experimental data. The hex-dominant mesher available in the 3DExperience platform was also tested with this model (results not presented here). For solving the turbulent flow at higher flow rates we tested three models, namely Spalart-Allmaras, k-ε and k-ω models. Results of interest were the axial velocity, pressure and fluid shear stress along the centerline, axial velocity along a radius of cross section at selected planes, and wall shear stress. Overall, Abaqus/CFD results were found to match very well with experiment results and faired competently with other software results. The k-ω model was found to perform best among the three turbulent models. This detailed study provides valuable insight into effective strategies for modeling hemodynamics using Abaqus/CFD. Keywords: CFD, hemodynamics, FDA Nozzle, Sudden expansion, steady state solver, turbulent models, Spalart-Allmaras, RNG k-ε model, SST k-ω model 1. Introduction Computational Fluid Dynamics (CFD) techniques are intensely used in the design and analysis of components in heavy industries, such as the aerospace and automotive industries. Following the same trend this technology is being used in the medical device industry to design, develop and analyze devices like ventricular assist devices, prosthetic heart valves, stents, blood filters and hemodialysis catheters. However this technique can not only help in the design process of a device but also can predict the change in the hemodynamic environment (i.e., the change blood flow dynamics and/or the fluid forces exerted on the device). This information is crucial as blood flowing through medical devices is subject to hemolysis or thrombosis. Understanding the effect of hemodynamic forces like shear stress and exposure time on hemolysis has been an area of active research focus for the past decade. Though crucial understanding has been achieved and several power law based empirical formulas have been modeled, the use of CFD and blood damage models in predicting medical device safety has not been adequately proven yet. To successfully implement these techniques for medical device design and/or for regulatory reviews a better understanding of the computational model and experimental validation must be achieved and steps towards “standardized” practices should be taken. www.3ds.com/simulia Previous 7 Next Contents Cardiovascular With an aim to utilize this technology to its fullest and standardize it, FDA initiated a critical path initiative named as” Computational Fluid Dynamics (CFD)/Blood Damage Project”. The initiative consist of two projects, the first one is CFD/PIV analysis of a Nozzle model that represents an idealized medical device geometry and the second one is a CFD/PIV and blood damage experiments on a heart pump model. The goal of the first project, which has been evaluated in this study, was to assess the current state of the art in CFD analysis in modeling hemodynamics. The experimental PIV analysis was performed in 4 labs, 3 from University research groups and one from the FDA. The averaged data was used for validation purposes. Twenty eight groups from 8 countries participated in the CFD challenge and submitted their simulation results. This project was initiated in 2008 and the experimental results were published in 2010 (Hariharan, 2010). The computational results with experimental validation were published in 2011 (Stewart, 2011). Abaqus/CFD, being a relatively recent addition to the SIMULIA brand has not been documented in this challenge. Therefore the need to validate the standard hemodynamics process and assess the performance of Abaqus/CFD was the prime motivation for this study. The easy accessibility of the experimental results from the FDA repository is an additional reason for choosing this benchmark model as a validation standard. The FDA Nozzle model consists of a straight tube, a conical section that connects to a straight throat region of smaller diameter, and a sudden expansion region where the diameter suddenly increases back to the original tube diameter. The device was designed to include accelerating flow, decelerating flow, sudden expansion/gradual diffusion, variations in shear stress and velocities, and recirculating flow, all of which may be related to blood damage in medical devices. For this geometry, flow in one direction is a sudden expansion problem, whereas flow in the opposite direction is a conical diffusor problem. In this study we concentrated on the sudden expansion problem. The model was simulated with 5 different flow rates that cover the range of fully laminar to fully turbulent flow. Through this study we present the capabilities of Abaqus/CFD as a comprehensive tool offered by SIMULIA, capable of solving for laminar, transitional and turbulent flow with accuracy. Abaqus/CFD can be used within the Abaqus/CAE environment in all three stages of analysis - modeling/meshing, solving and post-processing. In version 6.13, Abaqus/CFD has introduced a steady state solver that is robust and uses a 2nd order accurate SIMPLE based algorithm. Abaqus/CFD has also introduced three turbulent flow models, the first one is a one equation Spalart-Allmaras model, and the other two are two equation RNG k-ε and SST k-ω models. The first two models are available in CAE and the last one is available through insertion in the input deck. Details of the modeling technique followed are discussed in the method section, and the simulation results are presented in the Results and Discussion section. 2. Methods Nozzle Geometry The dimensions of the nozzle geometry were obtained from the published literature (Stewart, 2011). The FDA Nozzle model consists of a straight tube of diameter (D) 12 mm, a conical section that reduces the tube diameter to 4 mm in the throat region, and a small straight throat section of 40 mm length which is followed by a sudden increase in tube diameter (D) to 12 mm. The inlet length and the outlet length were left for the analyst to choose. We have chosen the inlet length to be 10 times D. We have used a parabolic inlet velocity profile therefore larger inlet length was unnecessary. For assessing the outlet length, simulations were performed with different outlet length and the velocity profile at the outlet observed. Our objective was to obtain a parabolic velocity profile at outlet. Results with shorter outlet length resulted in a velocity profile where the peak velocity is eccentric. As the length was increased, the outlet velocity was fully developed and peak velocity was concentric. Depending on this criterion we chose the outlet length to be 60D. www.3ds.com/simulia Previous 8 Next Contents Cardiovascular The 3D model of the CFD domain as shown in Figure 1a was developed using the solid modeling features in Abaqus/CAE . Figure 1. (a) FDA Nozzle model geometry with cross sectional surfaces identified by numbers b) final hex mesh of the model c) view of the radial section of the mesh at inlet Meshing Linear Hex elements were used to mesh nozzle model in this study. The model volume was subdivided along the axial direction and also along the radial direction as shown in in Figure 1b. There were 13 cross sectional surfaces of interest identified for post processing of results that were used for volume subdivision along the axis. The subdivision in the radial direction helped in prescribing adequate seeds near the walls. A mesh refinement was performed by changing the global seed from 0.001 to 0.004 along with adequate seeding of radial and axial edges. The final mesh consists of 1.39 million elements. The final mesh pattern is shown in Figures 1b and 1c. Another meshing technique was also tested and used for this model. It is the hex-dominant mesher available as part of the R2014 3DExperience platform. It is capable of handling complex geometries in a robust manner, and creating high-quality Hex-dominant 3D meshes. The Hex-dominant 3D mesher includes surface wrapping technology to automatically extract the fluid domain by creating a clean closed ‘wetted’ surface for subsequent volume meshing. A snapshot of the mesh created for the nozzle geometry is shown in Figure 2 (a). A closer look of the mesh, shown in Figure 2 (b) shows the smooth transition from a finer mesh (near the wall) to a coarser mesh away from the wall. The robustness and the time it takes to create high quality meshes, makes this mesher useful for creating CFD meshes for both simple and highly complex geometries. The results obtained by using this mesh are not included in this paper, but the results were very comparable to the more traditional hexahedral element mesh created using Abaqus/CAE. www.3ds.com/simulia Previous 9 Next Contents Cardiovascular (a) (b) Figure 2: Nozzle mesh using R2014 3DExperience Mesher Materials Blood is modeled as incompressible fluid with a density of 1056 kg/m3. It is assumed be a Newtonian fluid with dynamic viscosity of 0.0035 N.s/m2 as per the guidelines of the challenge (Stewart, 2011) Boundary Conditions Simulations were performed with 5 different inlet velocities ranging from inlet Reynolds number 167 to 2167 as documented in (Stewart, 2011).The flow rates were labeled as Q1, Q2, Q3, Q4 and Q5.The average inlet velocity, flow rate, average velocity at throat and throat Reynolds number are presented in Table 1. Though the inlet Reynolds number indicates the flow to be in a fully laminar zone at Q1 to a transitional zone at Q5, the throat Reynolds number indicates that the flow type ranged from fully laminar at Q1 to fully turbulent at Q5 in the sudden expansion region. In this study we have chosen to represent the flow rate Q1 as fully laminar flow, flow rate Q3 as representative of transitional flow and flow rate Q5 as fully turbulent flow. A parabolic inlet velocity profile Vz was used to model the inlet flow distribution following Equation 1 1. where Vmax is 2 times the average inlet velocity, r is the radial coordinate , R is the inlet radius and a = Vmax and b = Vmax/R2. The outlet boundary condition was set to zero pressure and a no-slip boundary condition was applied to the wall surface of the model Table 1. Inlet flow rate, velocity, Reynolds number used in the simulation Inlet Re Average Inlet Velocity Flow Rate Peak Velocity (Parabolic) Rei m/s m3/s m/s a 167 0.05 5.22e-6 0.09 Q2 667 0.18 2.08e-5 Q3 1167 0.32 3.65e-5 Q4 1667 0.46 Q5 2167 0.60 Q1 www.3ds.com/simulia Velocity at throat Throat Re b m/s Ret 0.09 2563.52 0.42 500 0.37 0.37 10234.73 1.66 2000 0.64 0.64 17906.93 2.90 3500 5.21e-5 0.92 0.92 25579.14 4.14 5000 6.77e-5 1.20 1.20 33251.35 5.39 6500 Previous 10 Parabolic Profile Next Contents Cardiovascular Analysis Settings A steady state incompressible fluid flow analysis was performed in Abaqus/CFD. The steady state method is a new addition to Abaqus/CFD version 6.13 and is based on a second order accurate SIMPLE (Semi-implicit Method for Pressure Linked Equation) based algorithm. Depending on the problem, a proper choice of an under-relaxation factor can significantly increase the convergence rate of the solution. For this study we chose 0.7 for the Momentum equation and 0.3 for the pressure equations and obtained a converged solution in 1000-1500 iterations. For the lowest flow rate Q1 the laminar flow model was used. For all other higher flow rates, we have used the laminar flow model along with three turbulent flow models. Abaqus/CFD offers a Spallart-Allmaras model and RNG k-ε model in CAE and the SST k-ω model can be accessed through the input deck. The Spallart-Allmaras model is a one equation Reynolds Averaged Navier Stokes (RANS) turbulent model that solves for the Reynolds’s stress when the Spalart-Allmaras variable kinematic turbulent viscosity is known. Abaqus/CFD recommends calculating the kinematic turbulent viscosity (ν∼) as 3 to 5 times of the kinematic viscosity (ν) of the fluid. We used a value of 1.3258e-5 Stoke, which is 4 times of kinematic viscosity. RNG k-ε model and SST k-ω models are two-equation RANS models that solve for Reynold’s stress by including two transport equations to represent the turbulent properties of the flow, like turbulent kinetic energy(k) and turbulent dissipation (ε) and /or specific dissipation(ω). The RNG k-ε model is a variant of k-ε models where the Re-Normalization Group (RNG) method is used to re-normalize the Navier stokes equation. The parameters k and ε can be calculated using Equation 2 if the inlet velocity (u0), turbulent intensity (I) and turbulent length scale (l) is known. 2. Cµ is an empirical constant whose value is 0.085. The SST k-ω model is a two equation eddy viscosity model where a Shear Stress Transport (SST) formulation is combined with a standard k-ω model and represents the “best of both” of the models. The variables for the transport equation are available in Abaqus/CFD whereas the parameters k, and ω can be calculated if the inlet velocity, turbulent length scale and turbulent intensity are known using Equation 3. 3. For this study the turbulent intensity and length scales were not available so we assumed the values to be 0.057 and 0.84 mm, given that turbulence intensity (I ) at core of a fully developed pipe flow can be approximated by and the turbulent length scale at inlet can be approximated by 0.07*L, where L is the characteristic length, in this case it is the inlet diameter of the model. Table 2 summarizes the parameter values used for 4 different flow rates (Q2, Q3, Q4 and Q5). Table 2. Turbulence parameters used for RNG k-ε model and SST k-ω model, for flow rates corresponding to throat Reynolds number u (inlet) k eps (ε) omega(ω) Turbulence Intensity length scale Q2:Re = 2000 0.1839 1.64E-04 3.96E-04 28.31 0.057 8.40E-04 Q3:Re = 3500 0.3218 5.04E-04 2.12E-03 49.53 0.057 8.40E-04 Q4:Re = 5000 0.4606 1.03E-03 6.23E-03 70.89 0.057 8.40E-04 Q5:Re = 6500 0.5986 1.74E-03 1.36E-02 92.13 0.057 8.40E-04 0 Simulations were performed on a desktop computer equipped with an Intel CoreTM i7 CPU and 16 GB of RAM. Simulation times were between 5 hours to 12 hours. www.3ds.com/simulia Previous 11 Next Contents Cardiovascular 3. Results & Discussion 3.1 Centerline Velocity and Pressure Results were post-processed in the Abaqus/CFD visualization module. Snapshots of the velocity magnitude plot from three different flow rates are shown in Figure 3. The snapshot shows the jet length in the sudden expansion zone for laminar flow model with flowrate Q1 and turbulent flow model with flowrates Q3 and Q5. It should be noted here that there are some instabilities seen for the Laminar case near the outlet. This is an outcome of the steady state plot shown at an instant prior to fully converged solution. Figure 3 : Velocity magnitude plot obtained from laminar model for flow rate Q1 and turbulent model (SST k-ω model) are plotted for flow rates Q3 and Q5. A direct comparison of the Abaqus/CFD results with the experimental data obtained from FDA repository is presented in Figure 4. Figure 4 : Centerline velocity in the axial direction, obtained from laminar and turbulent simulations, are plotted for flow rates (a) Q1 , (b) Q3 and (c) Q5. Similarly, wall pressure, represented by the centerline pressure in the axial direction is plotted for flow rates (d) Q1 (e) Q3 and (f) Q5. Flow rates Q1, Q3 and Q5 are representative of fully laminar flow, transitional flow and fully turbulent flow. www.3ds.com/simulia Previous 12 Next Contents Cardiovascular Centerline velocity in the axial direction and pressure data were recorded using the ‘path’ feature in the visualization module and plotted in Figure 4. The data obtained were for three different flow rates that are representative of fully laminar (Q1), transitional (Q3) and fully turbulent (Q5) flow types. Due to very little radial pressure gradient, wall pressure is represented by the centerline pressure. The experimental data presented in Figure 4 are an average of 3 different experimental data sets obtained from the FDA repository [3]. The experimental and computational pressure values are referenced with respect to the pressure at the location of the sudden expansion (z=0.0). Our simulation results matched very closely with the experimental data for laminar, transitional and turbulent flow. As the flow through the nozzle before the sudden expansion (z=0) region is truly laminar, the centerline axial velocity characterized by the laminar flow model was in excellent agreement with the experimental data for all three flow rates. Flow beyond the sudden expansion zone is complicated and none of the models were an exact match with the experimental data. Similar behavior was reported by the other contributors to the challenge and this was comprehensively discussed in the FDA publication (Stewart, 2011). At the sudden expansion zone the laminar flow model did well for lower flow rates and demonstrated fluctuations for transitional and fully turbulent flow rates. Here, the SST k-ω model and Spallart-Allmaras model achieved the best agreement with the experimental data. Figure 4(d-f) represents the pressure plots. Clearly pressure from laminar flow model for flow rate Q1 (Figure 3d) did not match with the experimental data. As reported (Stewart, 2011) results from other computational studies did not match with the experiment either. The errors are attributed to experimental error in pressure measurement, and were partly due to wrong reference scaling and were partly for not using differential pressure transducers and were discussed in the literature (Stewart, 2011). For higher flow rates (Q2 & Q3), the largest discrepancy in pressure was observed for the RNG k-ε model. The discrepancy in the pressure characteristics by the RNG k-ε model can be attributed to its inability to characterize turbulence in multiple different length scales as the eddy viscosity is determined from a single selected turbulence length scale. It is also well known that the RNG k-ε does not perform well for transitional Reynolds number, which we believe is the case in cases Q2 and Q3. In addition, the impact of wall mesh resolution requirement (y+ near the wall) could also be a contributing factor to the deviation from experimental results. Future work will involve using the newly introduced realizable k-ε turbulence model, which uses a hybrid wall functions approach thereby alleviating the need for creating meshes with specific near wall resolution (any y+ value near the wall) as opposed to the RNG k-epsilon, which had a stringent near wall mesh requirement (y+>20 near the wall). Similar behavior of the k-ε model was observed by users of other software as well (Stewart, 2011). Though none of the turbulent flow models were an exact match with the experimental data at each location, the Spalart-Allmaras and SST k-ω models did a superior job in characterizing the axial velocity and pressure before and beyond the sudden expansion zone, whereas the laminar model results deviated slightly in the sudden expansion region. The next section describes the axial velocity profile in the radial direction at specific cross sectional surfaces along the axial length of the model. 3.1. Axial velocity profile at selected cross sections Figure 5 : a) Axial velocity in radial direction at specific cross sectional surfaces (Surf-1,Surf-3,Suf-6,Surf-8 and Surf-13) are obtained from simulation and plotted along each row for three different flow rates (b) Q1 (c) Q3 and (d) Q5. www.3ds.com/simulia Previous 13 Next Contents Cardiovascular The nozzle model was subdivided in such a way that 13 cross sectional surfaces were created as shown in Figure 1a. Out of the 13 surfaces, results from only the following surfaces of interest were presented in this study: Surf-1 is closest to inlet, Surf-3 is at the conical section, Surf-6 is at the throat region, Surf-8 is in the sudden expansion zone and Surf 13 is the farthest away from the sudden expansion region. Axial velocities along the radial direction at these cross-sections were recorded for three strategic flow rates (Q1,Q3 and Q5) and presented in Figure 5 (a, b and c) respectively. The experimental data presented here is from one set of experiments obtained from the FDA repository [3] that best matches the average experimental data. As expected, the velocity profile from the laminar flow model matched very well before the sudden expansion zone for all three flow rates. As flow rate increases, the turbulent flow models did a better job in matching the experiment results at and beyond the sudden expansion region, as demonstrated clearly by the Surf-13 results described in Figure 5b and Figure 5c. 3.2. Radial shear rate profile at selected cross sections Figure 6 : Shear Rate in radial direction at specific cross sectional surfaces (Surf-1, Surf-3, Suf-6, Surf-8 and Surf-13), obtained from simulation and plotted along each row for three different flow rates (b) Q1 (c) Q3 and (d) Q5. We have further post processed hemodynamic parameters such as shear rate at selected cross sectional surfaces (Figure 6). Shear Rate is the ratio of velocity and distance and fluid shear stress can be calculated by multiplying shear rate and dynamic viscosity. Results from the laminar flow models matched the shear rate characteristics very well for lower flow rate Q1 and transitional flow rate Q3. With flow rate Q3, a larger discrepancy in shear rate was observed at Surface 10 and none of the models were effective in capturing the true behavior of the flow observed through experiment. However, as the flow rate increased to Q5, the turbulent flow models were successful in capturing the flow behavior at this location. As flow rate Q5 is still in a transitional zone in the inlet tube, the laminar flow model did a better job in capturing the flow behavior at this location. We have further post-processed the wall shear stress (data not presented here) and found very good correlation with the experimental data. Overall, the simulation results are in very close agreement with the experimental data and were satisfactory, particularly in light of the wide variability of data presented in the challenge. 4. Conclusion Abaqus/CFD was used to model, simulate and visualize the flow through a benchmark FDA nozzle model and the simulation results were validated by experimental data obtained from the FDA repository and presented in this technology brief. Abaqus/CFD now has the steady state feature, three turbulent models (Spalart-Allmaras model, RNG k-ε and SST k-ω) and a laminar flow model. All www.3ds.com/simulia Previous 14 Next Contents Cardiovascular of these features were put on test for modeling hemodynamics. The results illustrate the challenging nature of this benchmark and highlight areas where particular user attention is required. This study also revealed that it can be challenging to predict which model type is appropriate for any general medical device geometry for a given flow rate, as laminar, transitional and turbulent flow regions may co-exist in different locations at the same time. This also emphasizes the fact that CFD is not a “push button” technology, and user skill is an important factor in achieving accurate results, while the necessary details in modeling, meshing and analysis set-up play a crucial role in achieving higher levels of accuracy. 5. References 1. Abaqus User’s Manual, Version 6.10-1, Dassault Systèmes Simulia Corp., Providence, RI. 2. FDA Critical Path Initiative-Computational Round Robin #1 (https://fdacfd.nci.nih.gov/interlab_study_1_nozzle) 3. Hariharan, P. et al., Multi-laboratory Particle Image Velocimetry Analysis of the FDA Benchmark Nozzle Model to Support Validation of Computational Fluid Dynamics Simulations, Journal of Biomechanical Engineering, April 2011, Vol. 133 / 041002-1 4. Stewart, S. et al., Assessment of CFD Performance in Simulations of an Idealized Medical Device: Results of FDA’s First Computational Inter-laboratory Study, Cardiovascular Engineering and Technology, Vol. 3, No. 2, June 2012 ( 2012) pp. 139–160 www.3ds.com/simulia Previous 15 Next Contents Drug Delivery Novo Nordisk Makes Designing Injection Pens a Snap (Fit) Injection typically involves twisting a short needle onto the pen, turning a dial to the required dose, and pushing a button to deliver the medication under the skin. After a given number of doses are injected, the cartridge is exchanged for a new one (with a durable device) or discarded (with disposable pens). Audible clicks that occur at key stages of this procedure reassure the patient that they are engaging the device correctly at each step. It looks pretty easy. A one-minute video of a woman checking her blood sugar and then using a NovoPen to inject insulin is available here: http://www.novonordisk.com/press/ broadcastroom/default.asp. But every one of those reassuring clicks represents a challenge that has been overcome by the engineers who created the pens. So do the clicks the patient never hears: those that occur as the pen parts are assembled in the factory before use. Leading innovators in diabetes care use realistic simulation to improve product integrity from design to manufacture “Parts that click into place with ‘snap fit’ instead of screw connectors are very efficient to assemble within mass production,” says Torben Strøm Hansen, principal scientist in the Device R&D division of Novo Nordisk, near Copenhagen, Denmark. “Snap fit is the commonly used way to connect parts in our device mechanisms, and it also signals reliability when the internal components have optimal connections that don’t rattle. It’s very efficient when designed correctly.” Getting those designs correct from the start is the task that Hansen and his Mechanical Analysis team focus on in close collaboration with Novo Nordisk’s mechanical designers. “Even though an injection pen is not that big, there are a lot of fine details in its design,” he says. Whatever the configuration of device, the plastic-polymer components must withstand the rigors of both manufacturing and patient use, performing as required at different temperatures and loads. Many medical conditions can be treated with tablets, but others require injections under the skin in order for therapeutic drugs to reach the bloodstream. In the case of insulin administration for diabetes treatment, patients need to self-inject the drug daily. Making those injections easy and safe is of prime importance for Novo Nordisk, the Danish company that has been a world leader in the production of insulin ever since it was discovered by Canadian scientists in the 1920s. The company innovated beyond standard syringe technology to produce the world’s first patient-friendly self-injection system, the NovoPen, some 25 years ago. With more than 350 million diabetics worldwide—8.3 percent of the global population and growing, according to the International Diabetes Federation—demand for insulin pens will likely remain strong into the foreseeable future. Since effective control of the disease is dependent on consistent use of the drug, these delivery systems need to be portable, easy-to-use, reliable, and even resistant to minor misuse by patients. Small medical device, big design task An insulin pen may be small, but it is a precision instrument with a number of complex parts that must work in perfect tandem. Some pens are durable, containing a replaceable drug cartridge, while other disposable ones come pre-filled with the drug. www.3ds.com/simulia Previous (Top) CAD image of diabetes pen components. Grey and red cylindrical parts are snapped onto the green. (Bottom) An injectionmolded ratchet component from a medical device used by Novo Nordisk for a benchmark study. 16 Next Contents Drug Delivery Teaming up to model polymer behavior To ensure the integrity of their designs, Hansen and his team in the Device Simulation department rely on computer simulation with Abaqus finite element analysis (FEA). “More than a decade ago, my colleagues and I explored a number of commercial software codes,” he says. “We chose Abaqus because it was a well-integrated solution that provided both implicit and explicit capabilities, and could model the nonlinear behavior of the fine details in our designs correctly, including the high number of interfaces in contact.” FEA analysis demonstrates the ultimate straining of the snap fit on the device components during assembly. Over time, the group’s device models have become more refined, sophisticated, and computationally demanding. The SIMULIA Polymer Customer Review Team (PCRT) has worked closely with Novo Nordisk all along to provide updated enhancements in Abaqus that enable the company to model and predict the complexities of polymer behavior with increasing accuracy and efficiency. (The PCRT includes members from the automotive, high-tech, life sciences, and consumer goods industries.) Final snap deformation model (left) in Abaqus FEA (using a parallel framework to capture the changing behavior of the polymer material under cyclic loading) shows a more accurate plastic deformation of 0.12mm. The earlier elastic-plastic model (right) over-predicted that deformation would be 0.66mm. The snap-fit challenge A recent focus on snap fits in insulin pens demonstrates the challenges the team has faced when modeling polymers. “We concentrated on snap fits because they demonstrate almost ideal cyclic loading, with parts repeatedly loading and unloading from single to multiple cycles,” says Hansen. framework makes use of an arbitrary number of viscoelastic networks and an elastic equilibrium network to create a specific nonlinear viscoelastic model that is used to predict and track changes in the internal structural networks of a polymer as the material responds to repeated cyclic loads during snap fit. The material parameters in the FEA model are updated at each time step to reflect the new, altered state of the polymer. Since every type of polymer shows a different response to temperature, load, etc., the team continues to explore ways to identify the material characteristics of different polymer networks. During such cycles, the viscous nature of the thermoplastic material determines how the bouncing back to ‘normal’ occurs. Prediction through analysis of such time-dependent behavior is key to the device development process. Since devices can be subjected to different environments, including elevated temperatures, the function of the device must be as unaffected as possible by such changes and always comply with the specification even though the material properties of the components vary. Even just sitting on a pharmacy shelf or in the medicine cabinet, polymer materials are prone to creep and relaxation over time at rates that can vary with the temperature. Some polymers are also more complex than others: those used in durable devices may contain carbon or glass fillers that show anisotrophic behavior, which can be hard to predict. Not only are such advanced models useful to designers finetuning the latest pen configuration, the data can help inform manufacturing processes in the factory. “We have a processsimulating capability, through Moldflow, for which Abaqus has an interface. This allows us to input the stress fields that result from the molding process right into our models,” says Hansen. “As a result, we have greater insight into our manufacturing process and are more able to design parts that have a very low level of residual stresses in critical regions. “SIMULIA is working closely with us to provide capabilities we need,” says Hansen. “Having material models incorporating time-dependent viscous behavior is very important for our work and we’re now able to simulate both creep and relaxation with Abaqus. We are investigating how well the model will adapt to different kinds of thermoplastics, which may require different networks. Calibration will be key going forward.” “Modeling these diverse material characteristics as well as the behavior of the polymer as the load induces larger strains closer to yield is difficult,” says Hansen. “In order to predict such viscoelasticity precisely, we needed a more refined model that goes beyond a mainstream elastic-plastic approach.” Modeling material that’s constantly changing For More Information The team is now using the ‘parallel rheological framework’ methodology available in Abaqus to model polymer nonlinear viscoelasticity with greater accuracy than ever before. The www.3ds.com/simulia Previous www.novonordisk.com www.3ds.com/SCN-February2013 17 Next Contents Orthopedics Optimization of Surgical Positioning in Total Hip Replacement c) d) Jacob Elkins, John Callaghan, Douglas Pedersen, and Thomas Brown (Department of Orthopaedics and Rehabilitation, University of Iowa, Iowa City, Iowa) a) When joint pain and loss of mobility occur as a result of end-stage osteoarthritis or other severe hip pathologies, over 250,000 people choose to have total hip replacement (THR) surgery. Even though THRs are one of the most successful surgical inventions in medical history, they do fail. THR failures are often grouped as “early” or “late,” with early failure usually due to dislocation of the head from the cup, and late failure frequently due to adverse biologic reaction to wear debris generated at the bearing surface. Despite nearly six decades of investigation, the ideal surgical orientation of THR components remains unclear. Positioning of total hip bearings involves significant tradeoffs, as cup orientations most favorable in terms of stability are not necessarily ideal in terms of reduction of contact stress and wear potential. Previous studies and models have not addressed these potentially competing considerations for optimal THA function. Additionally, it is currently unknown whether the ideal orientation varies on implant parameters, such as variations in femoral head size. We, therefore, investigated optimal surgical cup orientation with a previously generated and physically validated finite element (FE) model of metal-on-metal THR. Figure 1. The FE model consisted of bony anatomy (a) and the hip soft tissues (b, anterior region of capsule rendered transparent for clarity). Four values of femoral anteversion were considered (c) as were five distinct femoral head sizes (d). Fitted ellipse a) f) c) g) d) Isosurface h) Figure 2. For every combination of femoral head size and femoral anteversion (20 such combinations total), the Stability Score (a) and Wear Score (b) are combined to determine the Performance Score (c). The optimal orientation is determined as the center of an ellipse fitted to an isosurface of scores > 90 (d). When considering all 20 combinations, regressions could be performed demonstrating optimal surgical orientation (e-h). Discussion The “landing zone” of ideal cup orientation did not increase with increased head size, challenging the presumption that larger heads are more forgiving in terms of stability and durability. Additionally, ideal cup positioning was considerably more sensitive to cup anteversion than to inclination. Finally, the current investigation is the first to quantitatively suggest that ideal cup positioning varies with both femoral anteversion and femoral head size. The FE model consisted of bony anatomy and the hip soft tissues (see Figure 1). Five dislocation-prone motions as well as gait were considered, as were permutations of femoral anteversion (0° to 30°), femoral head diameter (32 mm to 48 mm), cup inclination (25° to 75°), and cup anteversion (0° to 50°), resulting in 4,320 distinct FE simulations. A novel metric (“Performance Score”) was developed to delineate optimized cup orientation by considering both surface wear and component stability (see Figure 2 A-D). Positioning THR bearings involves significant tradeoffs with regard to stability and long-term bearing wear. The computational analysis identified optimal orientations to balance these considerations. These tradeoffs help explain the alarming rates of adverse local tissue response reported for large head metalon-metal THR devices that have demonstrated an improvement in joint stability. The conclusions from this study can readily be translated to other hard bearing surfaces—including ceramics and highly cross-linked polyethylene—suggesting careful consideration of the choices and compromises in THA design are required for all bearing couples. All FE simulations were performed using Abaqus/Explicit. Results Ideal cup position was substantially more sensitive to cup anteversion than to inclination. Regressions demonstrated strong correlations between optimal cup inclination vs. head diameter (Pearson’s r = -0.88), between optimal cup inclination vs. femoral anteversion (r = 0.96), between optimal cup anteversion vs. head diameter (r = 0.99) and between cup anteversion and femoral anteversion (r = -0.98) (see Figure 2 E-H). Previous b) e) Method www.3ds.com/simulia b) For More Information www.uiortho.com 18 Next Contents Orthopedics Computational Study of Cortical Bone Screw Pullout using the eXtended Finite Element Method (XFEM) Emer M. Feerick, Patrick McGarry (National University of Ireland, Galway, Ireland) Abstract: A study of screw pullout from cortical bone has been conducted using the UDMGINI subroutine with the eXtended finite element method (XFEM). XFEM alleviates mesh dependency when computationally modeling crack initiation and propagation. Cortical bone is a naturally occurring composite with a distinctive aligned microstructure that leads to anisotropic material behavior and damage. In this study, using a UDMGINI subroutine, stress components relative to a predefined osteon orientation were computed and an anisotropic damage criterion was used to determine damage initiation and to predict crack propagation. A 2D model of a single screw embedded in cortical bone was generated. A displacement boundary condition was applied to the top surface of the screw. During pull-out, contact interactions were implemented between the newly formed surfaces during crack propagation, eliminating unphysical over closure of elements. The model predicted the correct patterns of damage and crack propagation for both longitudinal (osteons aligned parallel to screw axis) and transverse (osteons aligned perpendicular to screw axis) pullout tests compared to experimental observations. In both cases crack propagation was predicted in the direction of osteon alignment. Keywords: Screw Pullout, Abaqus/Standard, Extended Finite Element Model (XFEM), Crack Propagation, Damage, Failure. 1. Introduction 1.1 The Extended Finite Element Method The eXtended finite element method (XFEM) is a method developed for computationally modeling crack initiation and propagation while alleviating mesh dependency. It can be implemented to simulate initiation and propagation of a discrete crack along an arbitrary, solution dependent path. Additionally contact interactions can be applied to the newly formed surfaces exposed as a result of crack propagation. Thus the effects of friction between newly exposed surfaces can be accounted for. Previously, work has been conducted and microstructure models developed to model fracture of cortical bone. These models incorporated a damage process that initiated damage within an entire element and when the properties of the element degraded such that they could no longer carry load they were deleted from the mesh nucleating voids within the material that replicated experimentally observed fracture mechanisms (Feerick et al. 2011). The possibility to alleviate mesh dependency within these materials offers an advance on previous models as the density of the meshes becomes less influential. Also a material model that is capable of predicting the correct crack patterns without incorporating a detailed geometry of the microstructure will also facilitate the development of larger macro scale models requiring significantly lower computational power. A limitation of some current cortical bone failure models is that they would be too computationally expensive to incorporate in macro scale models of whole bones. Previous studies that generated models of cortical bone using XFEM incorporated detailed microstructures (Budyn et al. 2010; Abdel-Wahab et al. 2012). If these previously developed models were applied to a macro scale simulation of whole bones the computational demand would not be viable. However in the present study we investigate a new feature of Abaqus 6.11 which now facilitates the use of a user subroutine (UDMGINI) for an anisotropic damage initiation criterion without the need for complex microstructure geometry. 1.2 Cortical Bone Cortical Bone is a naturally occurring composite, consisting of several constituents that dictate the overall response of the material during loading. In the present study we simplify the complex microstructure of bone to osteons (the fiber) embedded in an interstitial matrix. Osteons consist of a series of concentric cylinders of lamellae with a central vascular canal. Osteons are typically 200 µm and are 1-2mm in length and run parallel to the long axis of the bone. The osteons are surrounded by an interstitial matrix consisting of hydroxyapatite. (Rho et al. 1998; Cowin 2001). 1.3 Screw Pullout Screws are used for orthopedic applications throughout the human body ranging from fracture fixation plates to spinal fusion rods. Maximizing the screw insertion depth in cortical bone is regarded as a way of significantly increasing the pullout strength of the screw (Pollard et al. 2010). A computational model capable of predicting the failure modes and loads that occur during screw pullout would provide a design tool for the evaluation of future designs as well as design optimization. It is important to understand www.3ds.com/simulia Previous 19 Next Contents Orthopedics the failure mechanisms that occur during screw pullout as this may determine the angle of insertion and geometry of screw threads to maximize the holding power of the screw. 2. Materials & Methods 2.1 UDMGINI Failure Criteria A new feature of Abaqus 6.11 is the ability to incorporate a user subroutine UDMGINI for user defined damage initiation criterion. The user may define as many failure indexes as required. In the present study we define two failure indexes for damage initiation in the x and y directions. Failure index one defines the fiber failure index for initiation. The UTS of the fiber as well as the ultimate shear failure of the fiber is defined by the user. Once the value of σf =1 the damage is initiated. Where: σf σ11 σ12 σf f σfτ f is the fiber damage initiation criterion is the current stress in the local x direction is the current stress in the local x-y direction is the UTS of the fiber is the shear failure strength of the fiber Failure index 2 defines the matrix failure index for initiation. The UTS of the matrix as well as the ultimate shear failure of the fiber is defined by the user. Once the value of σm=1 damage is initiated. Where: σm σ22 σ12 σ mf σ mτ f is the fiber damage initiation criterion is the current stress in the local y direction is the current stress in the local x-y direction is the UTS of the matrix is the shear failure strength of the matrix As each of the failure indexes are calculated during an iteration of the simulation an array referred to as FNormal is compiled. This array contains the normal direction to the fracture line for each failure mechanism. This means that as damage is initiated in a particular failure index, the direction in which the crack will initiate is returned by the subroutine. The crack direction is defined in terms of the local orientation assigned to the material. For the present study if failure index 1 (fiber failure) is initiated the crack will propagate in the local y direction, perpendicular to the fiber direction. However, if failure index 2 (matrix failure) is initiated the crack will propagate in the local x direction, parallel to the fiber direction. 2.2 Single Element Tests 2D single element tests were conducted to examine crack propagation based upon osteon alignment. A unit cell geometry was created. The material properties assigned to the model are summarized in Table 1. An orientation was applied to the material so that the local x direction represents the direction of the fibers. A single element was held fixed as shown in Figure 1. The block was held fixed in the x and y directions at the bottom left corner. Also the top left and right hand corners were held fixed in the x direction. A 1mm displacement boundary condition was applied in the vertical direction. Longitudinal (fibers orientated parallel to the direction of loading) and transverse (fibers orientated perpendicular to the direction of loading) simulations were conducted. www.3ds.com/simulia Previous 20 Next Contents Orthopedics All simulations were conducted using Abaqus/Standard 6.11 with a CPS4 element. For simulation of crack growth a process of damage evolution (DE) was applied to the model. Energy dissipation was used to determine crack growth. The value selected for cortical bone energy release was based upon those reported in the literatrure (Abdel-Wahab et al. 2012). Table 1. Cortical bone material properties. E1 17100 MPa σff 233 MPa E2 10100 MPa σ fτ f 85 MPa E12 3300 MPa σ mf 51 MPa ν 0.3 σ mτ f 34 MPa DEm 0.3 N/mm DEf 0.3 N/mm Figure 1. Single element boundary conditions. 2.3 Screw Pullout A 2D model was developed for screw pullout to illustrate the potential application of XFEM for longitudinal and transverse pullout simulations. A limitation of XFEM is that axisymmetric elements cannot be used with the enrichment feature required for XFEM. Thus a plane stress model is presented here. The boundary conditions applied to the model are shown in Figure 2. A displacement boundary condition of 5mm was applied to the top of the screw. The edge of the screw was held fixed in the global x direction. The edge of the cortical bone was held fixed in both the global x and y directions. Local orientations were applied to assign fiber directions for each simulation. For a longitudinal simulation the fiber direction (local x) was orientated parallel to the axis of the screw. For a transverse simulation the fiber direction (local x) was orientated perpendicular to the axis of the screw. 35,000 CPS4 elements were used with element size of 0.05mm. www.3ds.com/simulia Previous 21 Next Contents Orthopedics Figure 2. Screw pullout boundary conditions for a longitudinal and transverse simulation. 3. Results 3.1 Single Element Tests The results for single element simulations are shown in Figure 3. A longitudinal simulation predicts vertical crack formation. A transverse simulation predicts horizontal crack formation. In both a longitudinal and transverse simulation failure index two reaches the value of 1 first. Thus cracks will propagate parallel to the direction of the fiber. Hence, cracks grow in the vertical direction for a longitudinal simulation, while cracks grow in the horizontal direction for a transverse simulation. Figure 3. Single element results for longitudinal and transverse simulations. 3.2 2D Screw Pullout Simulation The results for longitudinal and transverse screw pullout simulations are summarized in Figure 4. A longitudinal pullout simulation predicts localized deformation at the tips of the screw threads. Crack formation is predicted vertically upwards from the screw thread tips with the material between the screw threads removed with the screw leaving no thread definition on the fracture www.3ds.com/simulia Previous 22 Next Contents Orthopedics surface. A transverse pullout simulation predicts horizontal crack formation with deformation extending much further into the surrounding bone. In both longitudinal and transverse simulations failure index 2 reaches the value 1 first. For both simulations the cracks propagate parallel to the fiber direction results to significantly different failure modes as the screw is removed from the cortical bone. Figure 4. Crack patterns for longitudinal and transverse screw pullout simulations. 4. Discussion Results highlight the ability of the models to predict significantly different failure mechanisms for cortical bone. This is achieved by simply altering the orientation assigned to the material. It has previously been shown that fracture patterns during screw pullout from cortical bone are dependent upon osteon alignment (Feerick et al. 2011). In the past complex geometries of cortical bone with XFEM have been used to model the alternate fracture patterns of bone depending on osteon alignment (Budyn et al. 2010; AbdelWahab et al. 2012). However, in the present study, we have shown that without generating complex microstructure geometry and using the UDMGINI sub routine it is possible to define crack directions based upon the fiber (osteon) orientation. 5. Conclusion Screw pullout from cortical bone is one example of an application for the UDMGINI subroutine with XFEM. However the same method could be applied to cortical bone simulations for a wide range of orthopedic applications. These applications could range from orthopedic device design evaluation to modeling fracture incidence in whole bones. The use of the model is not limited to applications that consider cortical bone but rather any composite material containing a known fiber alignment. 6. References 1. Abdel-Wahab, A. A., Maligno, A. R. and Silberschmidt, V. V., 2012. “Micro-scale modelling of bovine cortical bone fracture: Analysis of crack propagation and microstructure using X-FEM”. Computational Materials Science 52, 128-135. 2. Budyn, É. and Hoc, T., 2010. “Analysis of micro fracture in human Haversian cortical bone under transverse tension using extended physical imaging”. International Journal for Numerical Methods in Engineering 82, 940-965. 3. Cowin, S. C., 2001. Bone mechanics handbook, CRC press USA. 4. Feerick, E.M.; Mullett, H; FitzPatrick, D; McGarry, J.P. “Computational Investigation of Cortical Bone Failure Mechanisms During Screw Pullout”, 19th Annual Symposium on Computational Methods in Orthopaedic Biomechanics (Pre-ORS), Long Beach, California, January 12, 2011. www.3ds.com/simulia Previous 23 Next Contents Orthopedics 5. Pollard, J. D., Deyhim, A., Rigby, R. B., Dau, N., King, C., Fallat, L. M. and Bir, C., 2010. “Comparison of Pullout Strength between 3.5-mm Fully Threaded, Bicortical Screws and 4.0-mm Partially Threaded, Cancellous Screws in the Fixation of Medial Malleolar Fractures”. The Journal of foot and ankle surgery : official publication of the American College of Foot and Ankle Surgeons 49, 248-252. 6. Rho, J.-Y., Kuhn-Spearing, L. and Zioupos, P., 1998. “Mechanical properties and the hierarchical structure of bone”. Medical Engineering & Physics 20, 92-102. 7. Acknowledgements The authors would like to acknowledge the assistance of Cheryl Liu at Dassault Systemes Simulia Corp for her assistance in the development of the failure criteria for the UDMGINI subroutine for the applications presented in this study. www.3ds.com/simulia Previous 24 Next Contents Cardiovascular SIMULIA Spearheads The Living Heart Project Bringing the medical community together for improved patient care The unveiling of Dassault Systèmes’ Living Heart Project at this year’s SIMULIA Community Conference (SCC) is certainly fitting. It represents the culmination of many months of collaboration, planning, and software development that were launched into action after the same event in 2013. While reviewing presentations in the medical special-interest group at the SCC in Vienna last year, Steven Levine, Senior Director, SIMULIA Portfolio Management, observed that simulation technology was approaching a tipping point in replicating the physics of living systems. Aware of the impact of cardiovascular disease, currently the #1 cause of death worldwide, he and his team had already performed introductory studies to explore the possibility of simulating a human heart. However, to make the translation from proof-of-concept to clinical practice, it would require an effort beyond what SIMULIA could provide alone. watching the project with great interest, as part of a Regulatory Science initiative. Its aim is to inspire leadership from the community—and SIMULIA and collaborators have stepped up through The Living Heart Project.” The plan was simple, yet visionary: use SIMULIA’s technology and The 3DEXPERIENCE® platform to challenge the scientific community to create a realistic 3D simulation of a beating, human heart—one that could be employed as a foundation to build personalized models based on a patient’s measured data to gain insight into potential cardiovascular diseases and methods of treatment. Levine approached Bernard Charlès, Dassault Systèmes CEO, with the idea—and The Living Heart Project was born. The project has high aspirations. Simulating an entire human heart that beats realistically is a highly coupled, multiscale, multi-physics problem layered on top of a complex materials engineering problem. Initially, the focus will be electromechanical applications for device design, such as pacemaker leads, stents, and artificial valves. However, project collaborators are already exploring its applicability to study treatment for heart disease or personalized ventricular assist devices. Future versions will include blood flow and thrombosisrelated applications, with increasingly detailed models of electrical pathways eventually reaching the cellular level. The knowledge developed under this project could one day pave the way for fully 3D-printed bioficial hearts. While simulation has gained broad adoption in many industries, healthcare remains the last frontier. Strides have certainly been made in medical device development, but the extreme complexity of the human body continues to hamper the quest for lifelike accuracy. What’s more, without such virtual models, regulators often have no choice but to accept laboratory testing as validation for new medical devices, with failures and recalls grabbing news headlines with disappointing regularity. Without question, simulation has great potential to reduce testing risks and costs while further benefiting patient safety and health. “This project began with a challenge to realistically simulate the physics of the heart,” says Levine. “As it evolves, we already see it beginning to transform how people think about what is possible and how this could serve as a platform to translate science into meaningful medical practice. We can imagine that this will lead to a new paradigm for data delivered directly into the hands of physicians; not simply 2D patient scan data, but fully analyzed 3D models where abnormal behavior is clearly identified and treatment options evaluated. The 3DEXPERIENCE platform could connect an ecosystem of providers to allow any physician access to these state-of-theart diagnostics.” In launching the project, SIMULIA reached out to researchers, device developers, and cardiac physicians, emphasizing the validity—and potential—of the effort. It was quickly recognized that, if successful, they would have an important platform for medical innovation. “There are millions of patients out there who really need this technology today, including my own daughter,” Levine says. “Participation has gained momentum beyond our expectations. Before we release it, the first model will be tested by the community in their respective disciplines. In fact, the FDA is www.3ds.com/simulia Previous For More Information www.3ds.com/heart 25 Next Contents Cardiovascular Abdominal Aorta Blood Flow Analysis with Abaqus/CFD Dassault Systèmes SIMULIA would like to thank Simpleware for providing the geometric model of the abdominal aorta. The progression of disease in arteries is affected by local blood flow characteristics. Arterial structural features, such as branches, bifurcations, and irregular shape and curvature changes, introduce complexities to the blood flow field. Further insight on how such mechanical factors affect arterial health can be gained by employing computational fluid dynamics tools. With this approach, accurate and detailed quantitative data on hemodynamic factors can be obtained. In this Technology Brief, we demonstrate the use of Abaqus/ CFD to simulate complex blood flow in an abdominal aorta and its branches. Blood velocity profiles and wall shear stresses are computed, and flow recirculation is visualized. It will be shown that blood flow analyses using Abaqus/CFD can provide comprehensive data that would be difficult to obtain from experimentation. Figure 1: Abdominal aorta incompressible viscous flow solver based on an advanced second-order projection method that uses a node-centered finite-element discretization for the pressure. This hybrid approach guarantees accurate solutions and preserves the local conservation properties associated with traditional finite volume methods. Abaqus/CFD can be used with unstructured grids and with various element types such as hexahedral, tetrahedral, wedge, and pyramid. The Abaqus/CFD solution methodology incorporates iterative Krylov solvers with algebraic multi-grid (AMG) preconditioning for solving the pressurePoisson equation. Background Many research studies have highlighted the importance of hemodynamic factors in vascular disease progression. The prevalence of disease in specific arterial locations is often correlated with flow features such as wall shear stress, stagnation, and recirculation zones. For example, atherosclerotic lesions are more likely to be found along the outer wall of the carotid sinus region of the carotid artery where the wall shear stress is low, while they are less likely to occur along the inner wall of carotid sinus where the wall shear stress is high [1]. Similarly, atherosclerotic disease is more likely to develop in the abdominal aorta below the diaphragm [2]. Dilatation of the aortic wall in the abdomen as seen in aortic aneurysms is also affected by blood flow behavior; the growth and rupture of such aneurysms are affected by the hemodynamics [3]. Renal artery stenosis may also alter flow characteristics; narrowing of the artery can cause the blood flow regime to transition from primarily laminar to turbulent. Changes in the shear stress distribution affect the blood pressure as well as the evolution of atherosclerosis [4]. Surgical procedures such as coronary artery bypass, where grafts reroute blood flow across a blocked artery, can cause changes in flow characteristics that initiate a biological response such as local growth and remodeling of surrounding tissues. Key Abaqus/CFD Features and Benefits • Transient incompressible viscous fluid flow analysis • Non-Newtonian viscosity models • Luminal time dependent pressure waveform boundary conditions • Pulsatile flow velocity boundary condition, derived from Womersley theory, implemented in a velocity user subroutine • Third-party boundary layer mesh import in Abaqus/CAE Geometry and Model Patient-specific abdominal aorta geometry is obtained using the ScanIP software from Simpleware. ScanIP reads CT scan files and creates a high quality triangulation for STL export. The STL file of the abdominal aorta is imported into CATIA V6R2013. A surface reconstruction is performed to obtain the CAD model of the abdominal aorta. The geometry is then trimmed at the boundaries to make it suitable for CFD analysis. Figure 1 shows the final CAD model of the abdominal aorta and branches. In this Technology Brief, we show that Abaqus/CFD can be used to model the pulsatile blood flow in the abdominal aorta and its various branches that supply blood to the organs in the abdomen. Abaqus/CFD uses a time accurate transient www.3ds.com/simulia Previous 26 Next Contents Cardiovascular Mesh Table 1: Non-Newtonian fluid properties for the Carreau-Yasuda model The CAD model is imported into ANSA v13.1.4 and meshed with a boundary layer meshing technique. A total of 6 wedge layers are generated along the wall. The mesh consists of 192375 nodes and 568769 elements. Of these, 314873 are tetrahedral elements, 249660 are wedge elements and 4236 are hexahedral elements. The mesh is exported in the Abaqus input file format and then imported into Abaqus/CAE where the CFD analysis attributes are defined. A close-up view of the boundary layer mesh is shown in Figure 2. Shear viscosity at low shear rate 0.025 Shear viscosity at high shear rates 0.0035 Time constant 25 Flow behavior index .025 Material constant 2.0 Figure 3: Viscosity v. shear rate for non-Newtonian blood properties Figure 2: Boundary layer mesh of the abdominal aorta the major branches of the abdominal aorta: the superior mesenteric artery, the left and right renal arteries, the inferior mesenteric artery, and the left and right iliac arteries. Material Blood is modeled as an incompressible fluid with a density r = 1050 kg/m3 and a non-Newtonian viscosity described by the Carreau-Yasuda model. The model properties are listed in Table 1. The Carreau-Yasuda model describes the shear thinning behavior of blood. It is often a reasonable approximation to treat blood as a Newtonian fluid in blood vessels greater than 0.5 mm in diameter. In vessels of such large diameter, the viscosity is relatively constant due to high rates of shear. However, non-Newtonian effects need to be accounted for in smaller vessels. The dependence of viscosity on the shear rate, on a log-log scale, is shown in Figure 3. A no-slip/no-penetration wall boundary condition is applied on the wall surface of the abdominal aorta. A time-dependent luminal pressure waveform as shown in Figure 4 is applied at the left and right iliac artery outlets. The pressure waveform is a triphasic pulse appropriate for normal hemodynamic conditions in the abdominal aorta below the renal artery ([7], [8]). The volumetric flow rates at the abdominal aorta inlet and superior mesenteric, left and right renal, and inferior mesenteric artery outlets are shown in Figure 5 ([5], [6]). These flow rate waveforms are represented as complex Fourier series and their coefficients are utilized for evaluating the velocity boundary condition at the inlet and outlets (except for the left and right iliac artery): Analysis Procedure A transient incompressible laminar fluid flow analysis is performed in Abaqus/CFD. Abaqus/CFD uses a projection method that enables segregation of pressure and velocity fields for efficient solution of the incompressible Navier-Stokes equations. The advective and diffusive fluxes are both treated implicitly. The solution method uses a second-order accurate least-squares gradient estimation. Simulation is performed for 4 cardiac cycles, with each cycle lasting 0.9 sec. Boundary Conditions Solving the transient Navier-Stokes and continuity equations requires specification of appropriate initial and boundary conditions. Blood enters the abdominal aorta at a level below the celiac artery. The model has six outlet sections representing www.3ds.com/simulia Previous 27 Next Contents Cardiovascular Womersley theory [9] evaluates the velocity profile at any crosssection for the unsteady, laminar flow of an incompressible fluid through a pipe of constant radius, when a time-varying pressure gradient is applied. The profile is defined as Here Qk represents the complex Fourier coefficients of the flow waveforms. The velocity boundary condition is specified at the inlet and outlets through user subroutine SMACfdUserVelocityBC. The user subroutine calculates a spatially- and time-varying velocity boundary condition based on Womersley theory and the volume flow rates shown in Figure 5. In this equation, αn is the non-dimensional Womersley number, where R denotes the radius, ρ the density, μ the viscosity and ω the circular frequency. Real(*) denotes the real part of a complex number and J0 and J1 are Bessel functions of the first kind of order 0 and 1, respectively. The choice of the Womersley velocity profile is a significant improvement over spatially-constant or parabolic velocity profiles since it captures the flow reversal at artery outlets and provides a more realistic boundary condition for pulsatile flows. Womersley theory, however, is limited to circular crosssections and hence, a mapping method ([10]) is used to map the evaluated velocity profile to non-circular cross-sections. Figure 4: Luminal pulsatile pressure waveform Velocity User Subroutine Implementation User subroutine SMACfdUserVelocityBC is written in C and provides access to the necessary model and solution quantities. Access is also provided to the MPI communicator and MPI routines for parallel implementation. The following steps are performed in the user subroutine to calculate velocities on the inlet and outlets: 1. Find the facets of elements lying on the edge of the surface. This is accomplished using a convex hull algorithm for finding edge points amongst a point cloud. 2. Transform all surface points to a local coordinate system specified by the user and with its origin lying at the surface centroid. 3. Express the points in a polar coordinate system with its origin located at the surface centroid. 4. Evaluate a normalized radius [10]: a. For each internal point (r, θ), find the two edge nodes i and j with θ coordinates such that θi ≤ θ ≤ θj. Figure 5: Flow rates for abdominal aorta, superior mesenteric artery, left and right renal artery, and inferior mesenteric artery www.3ds.com/simulia Previous b. Evaluate a normalized radius R by linearly interpolating the radius of edge nodes, Ri and Rj. The normalized radius serves as the radius for all points lying on the line connecting an edge point to the centroid. 28 Next Contents Cardiovascular 5. Evaluate the Womersley velocity using the radius r, normalized radius R, current time value, and Fourier coefficients for the flow waveform. 6. Transform the velocity to the global system for boundary condition specification. The above procedure can be skipped for circular surfaces. For these cases, the radius of the circular surface serves as the normalized radius for all internal points. Results and Discussion The kinetic energy of the fluid domain is plotted in Figure 6. It can be seen that the analysis reaches a steady periodic condition after 1 cardiac cycle. Figure 7 shows the velocity vectors along the midsagittal plane in the infrarenal segment of the abdominal aorta at various times during the cardiac cycle. A mildly recirculating flow pattern can be seen near the anterior wall; this vortex is formed during late systole and remains until peak systole when the effect of strong systolic acceleration results in a completely attached flow along the anterior wall. The particular location of the vortex is due to the curvature of the infrarenal segment which forces the flow to expand suddenly about the corner. The anterior portion of the abdominal aorta displays lower velocities than the core and posterior portion. Figure 6: Total kinetic energy of the flow Wall shear stresses during the cardiac cycle are contoured in Figure 8. The shear stress magnitude can be used to estimate the possibility of rupture of aneurysms. Figure 9 shows the surface traction vectors on the infrarenal segment of the abdominal aortic walls. Relatively low tractions are found compared to the branches where the velocities are higher. The primary movement of surface traction vectors is from the top of the vessel to the bottom, following the movement of the vortex during the cardiac cycle. Figure 7: Velocity vector field (m/sec) along the midsagittal plane in the infrarenal segment of abdominal aorta at (a) the beginning of the cycle, (b) peak systole, (c) end systole, and (d) mid diastole. Midsagittal plane definition (far right) and schematic representation of various time points during cardiac pressure cycle In Figure 10, streamlines are shown for the end systole. The streamlines show the formation of the vortex in the infrarenal segment and significant mixing of the blood. Figure 11 plots the percentage error in the mass balance, defined as the ratio of the sum of the inlet and outlet volume flow rates to the inlet volume flow rate. Summary In this technology brief, we demonstrate that Abaqus/CFD can be used to model pulsatile blood flow in the abdominal aorta and the various branches that supply blood to the organs in the abdomen. Newtonian and non-Newtonian blood properties can be modeled and specialized boundary conditions can be implemented with a user subroutine. This study can be easily extended to perform CFD analyses of diseased arteries. Parametric and quantitative studies on various hemodynamic www.3ds.com/simulia Previous Figure 8: Wall shear stress (N/m2) at (a) the beginning of the cycle, (b) peak systole, (c) end systole, and (d) mid diastole 29 Next Contents Cardiovascular 3. David A. Vorp, “Biomechanics of abdominal aortic aneurysm,” Journal of Biomechanics, Vol. 40, pp. 18871902, 2007. 4. F. Liang, Ryuhei Yamaguchi, H. Liu, “Fluid dynamics in normal and stenosed human renal arteries: an experimental and computational study,” Journal of Biomechanical Science and Engineering, Vol. 1, pp. 171-182, 2006. 5. D. N. Ku, S. Glagov, J. E. Moore, and C. K. Zarins, “Flow patterns in the abdominal aorta under simulated postprandial and exercise conditions: an experimental study,” Journal of Vascular Surgery, Vol. 9, pp. 309–316, 1989. 6. Charles A. Taylor, Thomas J.R. Hughes, and Christopher K. Zarins, “Finite Element Modeling of Three-Dimensional Pulsatile Flow in the Abdominal Aorta: Relevance to Atherosclerosis,” Annals of Biomedical Engineering, Vol. 26, pp. 975-987, 1998. Figure 9: Surface traction (N/m2) vector plot in the infrarenal segment at (a) the beginning of the cycle, (b) peak systole, (c) end systole, and (d) mid diastole 7. C. J. Mills, I. T. Gabe, J. H. Gault, D. T. Mason, J. Ross Jr., E. Braunwald, and J. P. Shillingford, “Pressure-flow relationships and vascular impedance in man,” Cardiovascular Research, Vol. 4, pp. 405-417, 1970. 8. C. M. Scotti, E. A. Finol, “Complaint biomechanics of abdominal aortic aneurysms: A fluid-structure interaction study,” Computers and Structures, Vol. 85, pp. 1097-1113, 2007. 9. J. R. Womersley, “Method for the calculation of velocity, rate of flow and viscous drag in arteries when the pressure gradient is known,” Journal of Physiology, Vol. 127, pp. 553-563, 1955. 10. J. P. Mynard, P. Nithiarasu, “A 1D arterial blood flow model incorporating ventricular pressure, aortic valve and regional coronary flow using the locally conservative Galerkin (LCG) method,” Communications in Numerical Methods in Engineering, Vol. 24, pp. 367-417, 2008. 11. Abaqus 6.13 Analysis User’s Manual, Dassault Systèmes Simulia Corp, 2013. Figures 10 & 11: End systole streamlines (left), mass balance percentage error (right) factors could be performed as well. Such studies can offer insight into the temporal and spatial variations of velocity and pressure fields as well as variations of wall shear stresses in vascular geometries. References 1. C. K. Zarins, D. P. Giddens, B. K. Bharadvaj, V. S. Sottiurai, R. F. Mabon, and S. Glagov, “Carotid bifurcation atherosclerosis: Quantitative correlation of plaque localization with flow velocity profiles and wall shear stress,” Circulation Research, Vol. 53, pp. 502–514, 1983. 2. J. C. Roberts, C. Moses, and R. H. Wilkins, “Autopsy studies in atherosclerosis: I. Distribution and severity of atherosclerosis in patients dying without morphologic evidence of atherosclerotic catastrophe,” Circulation, Vol. 20, pp.511–519, 1959. www.3ds.com/simulia Previous 30 Next Contents Cardiovascular Fluid-Structure Interaction Analysis of a Prosthetic Aortic Valve using Abaqus/Explicit Smoothed Particle Hydrodynamics Dassault Systèmes SIMULIA would like to thank Dr. Nandini Duraiswamy of the U.S. Department of Health and Human Services FDA Center for Devices and Radiological Health. Durability is a key measurement of prosthetic heart valve function. Assessment of fatigue life requires accurate estimates of the stresses induced during the cardiac cycle. Finite element (FE) studies have been used to estimate peak stresses in valves [1], and computational fluid dynamics (CFD) studies have been used to model blood flow around valves [2]. Fluid-structure interaction (FSI) studies are less common, in part because the closure of the valve creates CFD domain pinching. CFD can model the behavior of the blood, and a coupled fluidstructure interaction (FSI) analysis can capture the effect of the blood on the valve during the cardiac cycle. There is a final condition during the cycle that presents a challenge to coupled FSI: the fluid domain pinches during valve closure, which is a condition most CFD packages cannot handle. The Smoothed Particle Hydrodynamics (SPH) analysis technique, available in Abaqus/Explicit 6.11-1, addresses this challenge and makes modeling heart valves for the entire cardiac cycle possible, thus increasing the accuracy of prosthetic valve stress analysis. The smoothed particle hydrodynamic (SPH) analysis method in Abaqus/Explicit overcomes this difficulty. In this Technology Brief, the SPH technique will be used to determine the FSI response of a generic prosthetic heart valve. Background Key Abaqus Features and Benefits There are two principal modes of aortic valve disease: aortic stenosis, in which the valve no longer fully opens, and aortic regurgitation, in which the valve no longer fully closes. Either condition can eventually require the implantation of a prosthetic valve to replace the underperforming original. • Abaqus/Explicit Smoothed Particle Hydrodynamics capability for analyses involving extreme deformation • Robust hyperelastic material modeling • General contact capability for simplified definition of contact interactions Surgically implanted or transcatheter-delivered bioprosthetic aortic valve leaflets undergo dynamic cyclic loading and large deformation during the cardiac cycle. This can cause fatigue failure of the leaflets, compromising valve function and potentially affecting the patient. Accurate stress analysis of the valve during operation is therefore essential for designing durable aortic valves and improving patient outcomes. Analysis Approach Smoothed Particle Hydrodynamics SPH offers several advantages over CFD and coupled EulerianLagrangian methods in tracking free surface boundaries, handling small material-to-void ratios, and modeling extreme deformation with fragmentation. The latter capability makes it ideal for simulating the behavior of blood during valve closure and pressure changes. The operating conditions of the aortic valve are complex. The pressure on the aorta side of the leaflets is lower than that on the ventricular side when the ventricle is pumping oxygenated blood into the aorta, and the pressure on both sides varies depending on the stage of the cardiac cycle. This can be modeled by applying dynamic pressure loads (corresponding to loads measured in the aorta and left ventricle) directly onto the leaflets, which is an improvement in accuracy compared to previous analyses that used only static load conditions. SPH is part of a larger family of meshless numerical methods that define a body by a collection of points, instead of using nodes and elements. The SPH method implemented in Abaqus 6.111 uses a cubic spline kernel for interpolation, applying either a fixed or a variable “smoothing” length to particles. Internally, particle connectivity is determined based on smoothing length. The particles can contact Lagrangian bodies (in this case, the Even this method, however, does not account for the inertial and viscous effects of blood contacting the leaflets during flow. www.3ds.com/simulia Previous 31 Next Contents Cardiovascular valve leaflets) through the Abaqus/Explicit general contact feature. In addition, particles can be “glued” to Lagrangian bodies through *TIE constraints. SPH supports an extensive library of solid and fluid materials, including user materials. For this particular simulation, a finite volume of blood near the aortic valve was modeled with one-node PC3D elements. All particles had the same volume initially. There were 4956 particles, each with a radius of 1 mm. Material Modeling A generic aortic valve was meshed with shell (S4) elements. The valve had a diameter of 26mm and a thickness of 0.5mm. The junction between the aorta and the left ventricle was represented with a rigid tube, and two rigid plates were used to apply pressure on either side of the fluid particles (Figure 1). Figure 2: Leaflet material test data and Marlow model representation Figure 1: Aortic valve model with SPH particles The material for the valve was modeled with the Marlow isotropic hyperelastic representation, the general first-invariant hyperelastic material model in Abaqus. This model can exactly duplicate physical test data from one of several standard modes of loading (uniaxial, biaxial, or planar). It works well in situations where extensive data for one of the test modes is available. For the present analysis, uniaxial tensile test was used. (Figure 2). Figure 3: Pressure load profiles [1] Boundary and Loading Conditions Results and Conclusions Translational degrees of freedom were fixed for the valve edges. Left ventricle and aorta pressure profiles were applied to the end plates (Figure 3) [1]. The pressure profiles start from the point at which the pressure inside the left ventricle and the aorta are the same since the initial condition of the valve was stress-free. The same pressure was applied to the fluid as an initial condition. The end plates were not allowed to rotate, and because the finite volume of fluid is incompressible, the two rigid plates were constrained to have the same displacement along the axial direction using an equation constraint. Peak stress in the valve leaflets occurs during the diastolic phase, when the valve leaflets are closed. Higher stresses are observed in the FSI analysis using SPH than the reference model (Figure 4). In addition, the distribution of stresses is also different. Stress hot spots are observed in the middle of the leaflets as well as near the corners where two leaflets meet. This shows that, in addition to the pressure loads, the inertia effect of the fluid also influences the stress analysis results. The present SPH simulation capability is an important step toward providing prosthetic valve designers with increased simulation accuracy and the data needed to design more durable valves. As a reference model, a second analysis was run with the same (uniform) pressure profiles directly applied on the valve leaflets without the fluid. All other conditions were the same as the FSI model. www.3ds.com/simulia Previous 32 Next Contents Cardiovascular Figure 4: Mises stress on the leaflets during diastolic phase, SPH model (left) and reference model (right) References 1. H. Kim, J. Lu, M. S. Sacks and K. B. Chandran, “Dynamic Simulation of Bioprosthetic Heart Valves Using a Stress Resultant Shell Model”, Annals of Biomedical Engineering, Vol. 36, No. 2, 2008 2. W. Sun, K. Li, and E. Sirois, ”Simulated Elliptical Bioprosthetic Valve Deformation: Implications for Asymmetric Transcatheter Valve Deployment”, Journal of Biomechanics 43(16), 2010, 3. Abaqus References 4. For additional information on the Abaqus capabilities referred to in this brief, please see the following Abaqus 6.13 documentation reference: 5. Analysis User’s Guide: hydrodynamic analyses www.3ds.com/simulia 15.2.1 Smoothed particle Previous 33 Contents Dassault Systèmes, the 3DEXPERIENCE Company, provides business and people with virtual universes to imagine sustainable innovations. Its world-leading solutions transform the way products are designed, produced, and supported. Dassault Systèmes’ collaborative solutions foster social innovation, expanding possibilities for the virtual world to improve the real world. The group brings value to over 170,000 customers of all sizes in all industries in more than 140 countries. For more information, visit www.3ds.com. Americas Dassault Systèmes 175 Wyman Street Waltham, Massachusetts 02451-1223 USA Europe/Middle East/Africa Dassault Systèmes 10, rue Marcel Dassault CS 40501 78946 Vélizy-Villacoublay Cedex France Asia-Pacific Dassault Systèmes K.K. ThinkPark Tower 2-1-1 Osaki, Shinagawa-ku, Tokyo 141-6020 Japan ©2014 Dassault Systèmes. All rights reserved. 3DEXPERIENCE, the Compass icon and the 3DS logo, CATIA, SOLIDWORKS, ENOVIA, DELMIA, SIMULIA, GEOVIA, EXALEAD, 3D VIA, BIOVIA, NETVIBES, and 3DXCITE are commercial trademarks or registered trademarks of Dassault Systèmes or its subsidiaries in the U.S. and/or other countries. All other trademarks are owned by their respective owners. Use of any Dassault Systèmes or its subsidiaries trademarks is subject to their express written approval. Our 3DEXPERIENCE Platform powers our brand applications, serving 12 industries, and provides a rich portfolio of industry solution experiences.
Source Exif Data:
File Type : PDF File Type Extension : pdf MIME Type : application/pdf PDF Version : 1.7 Linearized : No Create Date : 2014:07:02 16:00:53-04:00 Creator : Adobe InDesign CS6 (Windows) Modify Date : 2014:07:14 15:11:07-04:00 Has XFA : No XMP Toolkit : Adobe XMP Core 5.2-c001 63.139439, 2010/09/27-13:37:26 Metadata Date : 2014:07:14 15:11:07-04:00 Creator Tool : Adobe InDesign CS6 (Windows) Instance ID : uuid:1ee0c640-d6d0-4942-a0b6-f0ce3a6d04ef Original Document ID : xmp.did:40F58A5B04DAE311A105EF713449498F Document ID : xmp.id:40E42C262302E411B948A0887374B344 Rendition Class : proof:pdf History Action : converted History Parameters : from application/x-indesign to application/pdf History Software Agent : Adobe InDesign CS6 (Windows) History Changed : / History When : 2014:07:02 16:00:53-04:00 Derived From Instance ID : xmp.iid:3FE42C262302E411B948A0887374B344 Derived From Document ID : xmp.did:40F58A5B04DAE311A105EF713449498F Derived From Original Document ID: xmp.did:40F58A5B04DAE311A105EF713449498F Derived From Rendition Class : default Format : application/pdf Producer : Adobe PDF Library 10.0.1 Trapped : False Page Layout : SinglePage Page Count : 34EXIF Metadata provided by EXIF.tools