Discretization error when using finite element models: Analysis and evaluation of an underestimated problem

Mesh convergence tests are often insufficiently performed in finite element analyses. There are many parameters which may have an effect on the mesh convergence behavior. The aim of this study was to identify the influence of different parameters on the mesh convergence behavior.For this purpose we used a simplified axis-symmetrical model of a single pedicle screw flank with surrounding bone to simulate a pull-out test. In parameter studies, the flank radii and the contact conditions at the bone–screw interface were varied. These parameter studies were carried out using an implicit and explicit solver. Thereby, the convergence criteria and the number of the substeps for the implicit nonlinear iteration process as well as the velocity and the material density for the explicit approach were considered.The mesh convergence behavior was influenced by varying the flank radii and the contact conditions. The implicit calculations led to a reaction force, which converged rapidly to a certain value with increasing mesh density, whereas the maximum von-Mises stress showed substantial convergence problems. The number of substeps and the convergence criteria of the iteration process strongly influenced the implicit solutions. In contrast, the maximum von-Mises stresses resulting from explicit calculations converged to a certain value after only a few refinement steps. Different pull-out velocities substantially affected the mesh convergence behavior, while the material density showed only a negligible influence.The results indicated the need to perform an appropriate mesh convergence test when using finite element methods. We were able to show that different parameters strongly influence the mesh convergence behavior and we demonstrated that convergence tests do not always lead to a satisfactory or acceptable solution.     

An anisotropic elastic-viscoplastic damage model for bone tissue

A new anisotropic elastic-viscoplastic damage constitutive model for bone is proposed using an eccentric elliptical yield criterion and nonlinear isotropic hardening. A micromechanics-based multiscale homogenization scheme proposed by Reisinger et al. is used to obtain the effective elastic properties of lamellar bone. The dissipative process in bone is modeled as viscoplastic deformation coupled to damage. The model is based on an orthotropic ecuntric elliptical criterion in stress space. In order to simplify material identification, an eccentric elliptical isotropic yield surface was defined in strain space, which is transformed to a stress-based criterion by means of the damaged compliance tensor. Viscoplasticity is implemented by means of the continuous Perzyna formulation. Damage is modeled by a scalar function of the accumulated plastic strain ${D(\kappa)}$ , reducing all element s of the stiffness matrix. A polynomial flow rule is proposed in order to capture the rate-dependent post-yield behavior of lamellar bone. A numerical algorithm to perform the back projection on the rate-dependent yield surface has been developed and implemented in the commercial finite element solver Abaqus/Standard as a user subroutine UMAT. A consistent tangent operator has been derived and implemented in order to ensure quadratic convergence. Correct implementation of the algorithm, convergence, and accuracy of the tangent operator was tested by means of strain- and stress-based single element tests. A finite element simulation of nano- indentation in lamellar bone was finally performed in order to show the abilities of the newly developed constitutive model.     

Automated three-dimensional finite element modelling of bone: a new method

Three-dimensional finite element stress analysis of bone is a key to understanding bone remodelling, assessing fracture risk, and designing prostheses; however, the cost and complexity of predicting the stress field in bone with accuracy has precluded the routine use of this method. A new, automated method of generating patient-specific three-dimensional finite element models of bone is presented — it uses digital computed tomographic (CT) scan data to derive the geometry of the bone and to estimate its inhomogeneous material properties. Cubic elements of a user-specified size are automatically defined and then individually assigned the CT scan-derived material properties. The method is demonstrated by predicting the stress, strain, and strain energy in a human proximal femur in vivo. Three-dimensional loading conditions corresponding to the stance phase of gait were taken from the literature and applied to the model. Maximum principal compressive stresses of 8–23 MPa were computed for the medial femoral neck. Automated generation of additional finite element models with larger numbers of elements was used to verify convergence in strain energy.     

Computer aided stress analysis of long bones utilizing computed tomography

A computer aided analysis method has been developed which utilizes computed tomography (CT) and a finite element (FE) computer program to determine the stress-displacement pattern in a long bone section. The CT data file provides the geometry, the apparent density and the elastic properties for the three-dimensional FE model. A developed pre-processor generates the FE model of a human diaphyseal tibia section which is then analyzed by the SAP IV finite element program. The results obtained are sorted and displayed by a developed post-processor and compared with stresses and deformations from the literature. The model generation method was verified by applying it to a model of simple geometry and boundary conditions, then comparing the results with the analytical solution of the same problem. The convergence behavior of nodal displacements was tested as a function of mesh refinement. This method provides an automatic, versatile, non-invasive and accurate tool of long bone modeling for finite element stress analysis.     

The use of strain energy as a convergence criterion in the finite element modelling of bone and the effect of model geometry on stress convergence

Keyak and Skinner¹ state that the adequacy of a finite element mesh may be verified by confirming that convergence of strain energy has been achieved. The paper concludes by stating that convergence of strain energy does not ensure that a particular mesh is adequate for producing accurate stress/strain results, and that only qualitative results can be achieved from 3D models of bone. This conclusion is challenged and the effect of element size on stress convergence is clarified. The influence of element size on stress convergence in regions of surface irregularity is studied and current methods are discussed that facilitate convergence.     

Sensitivity Analysis and Optimal Shape Design for Bone-Prosthesis Interfaces in a Femoral Head Surface Replacement

A numerical optimization procedure has been applied for the shape optimal design of a femoral head surface replacement. The failure modes of the prosthesis that were considered in the formulation of the objective functions concerned the interface stress magnitude and the bone remodelling activity beneath the implant. In order to find a compromising solution between different requirements demanded by the two objective functions, a two step optimization procedure has been developed. Through step 1 the minimization of interface stress was achieved, through step 2 the minimization of bone remodelling was achieved with constraints on interface stresses. The results obtained provided an optimal design that generates limited bone remodelling activity with controlled interface stress distribution. The computational procedure was based on the application of the finite element method, linked to a mathematical programming package and a design sensitivity analysis package.     

Implementation of controlling strategy in a biomechanical lower limb model with active muscles for coupling multibody dynamics and finite element analysis

Computational biomechanics for human body modeling has generally been categorized into two separated domains: finite element analysis and multibody dynamics. Combining the advantages of both domains is necessary when tissue stress and physical body motion are both of interest. However, the method for this topic is still in exploration. The aim of this study is to implement unique controlling strategies in finite element model for simultaneously simulating musculoskeletal body dynamics and in vivo stress inside human tissues. A finite element lower limb model with 3D active muscles was selected for the implementation of controlling strategies, which was further validated against in-vivo human motion experiments. A unique feedback control strategy that couples together a basic Proportion-Integration-Differentiation (PID) controller and generic active signals from Computed Muscle Control (CMC) method of the musculoskeletal model or normalized EMG singles was proposed and applied in the present model. The results show that the new proposed controlling strategy show a good correlation with experimental test data of the normal gait considering joint kinematics, while stress distribution of local lower limb tissue can be also detected in real-time with lower limb motion. In summary, the present work is the first step for the application of active controlling strategy in the finite element model for concurrent simulation of both body dynamics and tissue stress. In the future, the present method can be further developed to apply it in various fields for human biomechanical analysis to monitor local stress and strain distribution by simultaneously simulating human locomotion.     

Heterogeneous linear elastic trabecular bone modelling using micro-CT attenuation data and experimentally measured heterogeneous tissue properties

High-resolution voxel-based finite element software, such as FEEBE developed at the NCBES, is widely used for studying trabecular bone at the micro-scale. A new approach to determine heterogeneous bone tissue material properties for computational models was proposed in this study. The specimen-specific range of tissue moduli across strut width was determined from nanoindentation testing. This range was mapped directly using linear interpolation to that specimen's micro-computed tomography (microCT) grey value range as input material properties for finite element analysis. The method was applied to cuboid trabecular bone samples taken from eight, 4-year-old (skeletally mature) ovine L5 vertebrae. Before undergoing experimental uniaxial compression tests, the samples were microCT scanned and 30 microm resolution finite element models were generated. The linear elastic finite element models were compressed to 1% strain. This material property assignment method for computational models accurately reproduced the experimentally determined apparent modulus and concentrations of stress at locations of failure.     

Sensitivity Analysis and Optimal Shape Design for Bone-Prosthesis Interfaces in a Femoral Head Surface Replacement

Abstract

A numerical optimization procedure has been applied for the shape optimal design of a femoral head surface replacement. The failure modes of the prosthesis that were considered in the formulation of the objective functions concerned the interface stress magnitude and the bone remodelling activity beneath the implant. In order to find a compromising solution between different requirements demanded by the two objective functions, a two step optimization procedure has been developed. Through step I the minimization of interface stress was achieved, through step 2 the minimization of bone remodelling was achieved with constraints on interface stresses.

The results obtained provided an optimal design that generates limited bone remodelling activity with controlled interface stress distribution.

The computational procedure was based on the application of the finite element method, linked to a mathematical programming package and a design sensitivity analysis package.     

Evaluation of bone plate with low-stiffness material in terms of stress distribution

The aim of this study is to evaluate a newly developed bone plate with low-stiffness material in terms of stress distribution. In this numerical study, 3D finite element models of the bone plate with low-stiffness material and traditional bone plates made of stainless steel and Ti alloy have been developed by using the ANSYS software. Stress analyses have been carried out for all three models under the same loading and boundary conditions. Compressive stresses occurring in the intact portion of the bone (tibia) and at the fractured interface at different stages of bone healing have been investigated for all three types of bone-plate systems. The results obtained have been compared and presented in graphs. It has been seen that the bone plate with low-stiffness material offers less stress-shielding to the bone, providing a higher compressive stress at the fractured interface to induce accelerated healing in comparison with Ti alloy and stainless-steel bone plate. In addition, the effects of low-stiffness materials with different Young's modulus on stress distribution at the fractured interface have been investigated in the newly developed bone-plate system. The results showed that when a certain value of Young's modulus of low-stiffness material is exceeded, increase in stiffness of the bone plate does not occur to a large extent and stress distributions and micro-motions at the fractured interface do not change considerably.     

Finite element implementation of a new model of slight compressibility for transversely isotropic materials

Modelling transversely isotropic materials in finite strain problems is a complex task in biomechanics, and is usually addressed by using finite element (FE) simulations. The standard method developed to account for the quasi-incompressible nature of soft tissues is to decompose the strain energy function (SEF) into volumetric and deviatoric parts. However, this decomposition is only valid for fully incompressible materials, and its use for slightly compressible materials yields an unphysical response during the simulation of hydrostatic tension/compression of a transversely isotropic material. This paper presents the FE implementation as subroutines of a new volumetric model solving this deficiency in two FE codes: Abaqus and FEBio. This model also has the specificity of restoring the compatibility with small strain theory. The stress and elasticity tensors are first derived for a general SEF. This is followed by a successful convergence check using a particular SEF and a suite of single-element tests showing that this new model does not only correct the hydrostatic deficiency but may also affect stresses during shear tests (Poynting effect) and lateral stretches during uniaxial tests (Poisson's effect). These FE subroutines have numerous applications including the modelling of tendons, ligaments, heart tissue, etc. The biomechanics community should be aware of specificities of the standard model, and the new model should be used when accurate FE results are desired in the case of compressible materials.     

A 3D skeletal muscle model coupled with active contraction of muscle fibres and hyperelastic behaviour

This paper presents a three-dimensional finite element model of skeletal muscle which was developed to simulate active and passive non-linear mechanical behaviours of the muscle during lengthening or shortening under either quasi-static or dynamic condition. Constitutive relation of the muscle was determined by using a strain energy approach, while active contraction behaviour of the muscle fibre was simulated by establishing a numerical algorithm based on the concept of the Hill's three-element muscle model. The proposed numerical algorithm could be used to predict concentric, eccentric, isometric and isotonic contraction behaviours of the muscle. The proposed numerical algorithm and constitutive model for the muscle were derived and implemented into a non-linear large deformation finite element programme ABAQUS by using user-defined material subroutines. A number of scenarios have been used to demonstrate capability of the model for simulating both quasi-static and dynamic response of the muscle. Validation of the proposed model has been performed by comparing the simulated results with the experimental ones of frog gastrocenemius muscle deformation. The effects of the fusiform muscle geometry and fibre orientation on the stress and fibre stretch distributions of frog muscle during isotonic contraction have also been investigated by using the proposed model. The predictability of the present model for dynamic response of the muscle has been demonstrated by simulating the extension of a squid tentacle during a strike to catch prey.     

A model of anterior cruciate ligament reconstructive surgery: A validation construct and computational insights

This study sought to develop a computational framework that emulates the anterior cruciate ligament reconstruction surgery using transtibial portal technique. The proposed model included the tibia–femoral and patella–femoral joints, articular cartilage and menisci. Key surgical parameters were incorporated including bone-patellar-tendon-bone graft excision and pre-tensioning, tunnel morphology, bone plugs, and bone plug fixation. Several simulation steps were parameterized to reflect the clinically reported surgical procedure. Our focus was to explore the intra-operative effects of variations in tunnel directions on the selected metrics of joint mechanics during Lachman and Anterior Drawer tests. A mathematical construct capable of transforming the limited and heterogeneous experimental and surgical data to evidence-based validation was adopted to ensure the viability of the finite element models. We found that the proposed models, subject to a variation in tunnel directions, resulted in simulation outputs that favor the reported experimental data of Lachman and Anterior Drawer tests under uncertainty. Simulation results for a population of three-dimensional tunnel orientations provided insights into the graft–tunnel contact mechanics and the spatial stress distribution in the graft, insights that have been anecdotally observed in prior experimental studies. The intraarticular graft tension was found to be higher than the estimated in tunnel graft force, and larger differences were found for the least inclined tunnels exhibiting higher contact pressures, transverse bending and twisting of the graft and Von-Mises stress at the graft–femoral tunnel interface. Conversely, tunnels with high inclination angles exhibited higher intraarticular graft tension and Von-Mises stress at the graft–tibial bone plug interface.     

A method for measuring strains in soft tissue

A finite element based method has been developed for measuring strains in soft tissue. An array of markers is placed on the tissue surface and treated as nodes of a four node isoparametric element. The displacements of the marker centroids are directly measured using a high sensitivity television camera. Finite element method mathematics are then used to calculate the plane strain tensor at any point inside the element. The method has been implemented using non-rectangular elements that are approximately 2 mm on each side.     

Design of a custom angled abutment for dental implants using computer-aided design and nonlinear finite element analysis

Computer-aided design/computer-aided manufacturing (CAD/CAM) custom abutments have been attracting more and more attention due to their advantages of accuracy fit and esthetic emergence profile. However, the CAD key technology for custom abutments has been seldom studied as well as their biomechanical behavior. This paper explored a novel method to design a CAD/CAM custom angled abutment, evaluated the biomechanical performance of the whole system and compared the difference between the custom and the conventional abutment through 3D nonlinear finite element analysis (FEA). Firstly, the digital data of the dental casts at the healing abutment level was acquired by optical scanner. Thus the position of the healing abutment and the implant can be determined by CAD technology. The custom angled abutment was then designed according to the need of restoration and esthetics with CAD software. The described system can eliminate wax and cast, create an esthetic anatomical emergence profile and provide a satisfactory angle correction. Simulation results indicate that there was no distinct difference in the stress distribution and magnitude of implant-bone interface and screw using the custom or the conventional angled abutment.     

A finite shell element for heart mitral valve leaflet mechanics, with large deformations and 3D constitutive material model

This paper presents a shell finite element formulation appropriate for simulating the heart valve leaflet mechanics, including three-dimensional (3D) stress and strain effects. A 4-node mixed-interpolation shell is formulated in convected coordinates. This shell model is made capable of handling arbitrary 3D material models by use of an algorithm that satisfies the shell stress assumption at every element integration point. A method for tracking the fiber direction is incorporated. The resulting shell element operates under the same conditions as a standard 4-node shell element with 5 degrees of freedom per node, but extends the modeling capabilities to handle large-deformation and anisotropic behavior.     

Finite element estimates of interface stress in the trans-tibial prosthesis using gap elements are different from those using automated contact

When compared with automated contact methods of finite element (FE) analyses, gap elements have certain inherent disadvantages in simulating large slip of compliant materials on stiff surfaces. However, automated contact has found limited use in the biomechanical literature. A non-linear, three-dimensional, geometrically accurate, FE analysis of the trans-tibial limb-socket prosthetic system was used to compare an automated contact interface model with a gap element model, and to evaluate the sensitivity of automated contact to interfacial coefficient of friction (COF). Peak normal stresses and resultant shear stresses were higher in the gap element model than in the automated contact model, while the maximum axial slip was less. Under proximally directed load, compared with automated contact, gap elements predicted larger areas of stress concentration that were located more distally. Gap elements did not predict any relative slip at the distal end, and also transmitted a larger proportion of axial load as shear stress. Both models demonstrated non -linear sensitivity to COF, with larger variation at lower magnitudes of COF. By imposing physical connections between interface surfaces, gap elements distort the interface stress distributions under large slip. Automated contact methods offer an attractive alternative in applications such as prosthetic FE modeling, where the initial position of the limb in the socket is not known, where local geometric features have high design significance, and where large slip occurs under load.     

On the importance of considering porosity when simulating the fatigue of bone cement

Fatigue cracking in the cement mantle of total hip replacement has been identified as a possible cause of implant loosening. Retrieval studies and in vitro tests have found porosity in the cement may facilitate fatigue cracking of the mantle. The fatigue process has been simulated computationally using a finite element/continuum damage mechanics (FE/CDM) method and used as a preclinical testing tool, but has not considered the effects of porosity. In this study, experimental tensile and four-point bend fatigue tests were performed. The tensile fatigue S-N data were used to drive the computational simulation (FE/CDM) of fatigue in finite element models of the tensile and four-point bend specimens. Porosity was simulated in the finite element models according to the theory of elasticity and using Monte Carlo methods. The computational fatigue simulations generated variability in the fatigue life at any given stress level, due to each model having a unique porosity distribution. The fracture site also varied between specimens. Experimental validation was achieved for four-point bend loading, but only when porosity was included. This demonstrates that the computational simulation of fatigue, driven by uniaxial S-N data can be used to simulate nonuniaxial loadcases. Further simulations of bone cement fatigue should include porosity to better represent the realities of experimental models.     

Analysis of effects of friction on the deformation behavior of soft tissues in unconfined compression tests

Frictionless specimen/platen contact in unconfined compression tests has traditionally been assumed in determining material properties of soft tissues via an analytical solution. In the present study, the suitability of this assumption was examined using a finite element method. The effect of the specimen/platen friction on the mechanical characteristics of soft tissues in unconfined compression was analyzed based on the published experimental data of three different materials (pigskin, pig brain, and human calcaneal fat). The soft tissues were considered to be nonlinear and viscoelastic; the friction coefficient at the contact interface between the specimens and platens was assumed to vary from 0.0 to 0.5. Our numerical simulations show that the tissue specimens are, due to the specimen/platen friction, not compressed in a uniform stress/strain state, as has been traditionally assumed in analytical analysis. The stress of the specimens obtained with the specimen/platen friction can be greater than those with the frictionless specimen/platen contact by more than 50%, even in well-controlled test conditions.