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 transversely isotropic, transversely homogeneous microstructural-statistical model of articular cartilage

Articular cartilage is a multi-phasic, composite, fibre-reinforced material. Therefore, its mechanical properties are determined by the tissue microstructure. The presence of cells (chondrocytes) and collagen fibres within the proteoglycan matrix influences, at a local and a global level, the material symmetries. The volumetric concentration and shape of chondrocytes, and the volumetric concentration and spatial arrangement of collagen fibres have been observed to change as a function of depth in articular cartilage. In particular, collagen fibres are perpendicular to the bone-cartilage interface in the deep zone, their orientation is almost random in the middle zone, and they are parallel to the surface in the superficial zone. The aim of this work is to develop a model of elastic properties of articular cartilage based on its microstructure. In previous work, we addressed this problem based on Piola's notation for fourth-order tensors. Here, mathematical tools initially developed for transversely isotropic composite materials comprised of a statistical orientation of spheroidal inclusions are extended to articular cartilage, while taking into account the dependence of the elastic properties on cartilage depth. The resulting model is transversely isotropic and transversely homogeneous (TITH), the transverse plane being parallel to the bone-cartilage interface and the articular surface. Our results demonstrate that the axial elastic modulus decreases from the deep zone to the articular surface, a result that is in good agreement with experimental findings. Finite element simulations were carried out, in order to explore the TITH model's behaviour in articular cartilage compression tests. The force response, fluid flow and displacement fields obtained with the TITH model were compared with the classical linear elastic, isotropic, homogeneous (IH) model, showing that the IH model is unable to predict the non-uniform behaviour of the tissue. Based on considerations that the mechanical stability of the tissue depends on its topological and microstructural properties, our long-term goal is to clearly understand the stability conditions in topological terms, and the relationship with the growth and remodelling mechanisms in the healthy and diseased tissue.     

Finite element modeling of impulsive excitation and shear wave propagation in an incompressible,transversely isotropic medium

Elastic properties of materials can be measured by observing shear wave propagation following localized, impulsive excitations and relating the propagation velocity to a model of the material. However, characterization of anisotropic materials is difficult because of the number of elasticity constants in the material model and the complex dependence of propagation velocity relative to the excitation axis, material symmetries, and propagation directions. In this study, we develop a model of wave propagation following impulsive excitation in an incompressible, transversely isotropic (TI) material such as muscle. Wave motion is described in terms of three propagation modes identified by their polarization relative to the material symmetry axis and propagation direction. Phase velocities for these propagation modes are expressed in terms of five elasticity constants needed to describe a general TI material, and also in terms of three constants after the application of two constraints that hold in the limit of an incompressible material. Group propagation velocities are derived from the phase velocities to describe the propagation of wave packets away from the excitation region following localized excitation. The theoretical model is compared to the results of finite element (FE) simulations performed using a nearly incompressible material model with the five elasticity constants chosen to preserve the essential properties of the material in the incompressible limit. Propagation velocities calculated from the FE displacement data show complex structure that agrees quantitatively with the theoretical model and demonstrates the possibility of measuring all three elasticity constants needed to characterize an incompressible, TI material.     

Spherical indentation method for determining the constitutive parameters of hyperelastic soft materials

A comprehensive study on the spherical indentation of hyperelastic soft materials is carried out through combined theoretical, computational, and experimental efforts. Four widely used hyperelastic constitutive models are studied, including neo-Hookean, Mooney–Rivlin, Fung, and Arruda–Boyce models. Through dimensional analysis and finite element simulations, we establish the explicit relations between the indentation loads at given indentation depths and the constitutive parameters of materials. Based on the obtained results, the applicability of Hertzian solution to the measurement of the initial shear modulus of hyperelastic materials is examined. Furthermore, from the viewpoint of inverse problems, the possibility to measure some other properties of a hyperelastic material using spherical indentation tests, e.g., locking stretch, is addressed by considering the existence, uniqueness, and stability of the solution. Experiments have been performed on polydimethylsiloxane to validate the conclusions drawn from our theoretical analysis. The results reported in this study should help identify the extent to which the mechanical properties of hyperelastic materials could be measured from spherical indentation tests.     

Mechanical characterization of anisotropic planar biological soft tissues using finite indentation: experimental feasibility

Heart valve tissue engineering offers a promising alternative for current treatment and replacement strategies, e.g., synthetic or bioprosthetic heart valves. In vitro mechanical conditioning is an important tool for engineering strong, implantable heart valves. Detailed knowledge of the mechanical properties of the native tissue as well as the developing tissue construct is vital for a better understanding and control of the remodeling processes induced by mechanical conditioning. The nonlinear, anisotropic and inhomogeneous mechanical behavior of heart valve tissue necessitates a mechanical characterization method that is capable of dealing with these complexities. In a recent computational study we showed that one single indentation test, combining force and deformation gradient data, provides sufficient information for local characterization of nonlinear soft anisotropic tissue properties. In the current study this approach is validated in two steps. First, indentation tests with varying indenter sizes are performed on linear elastic PDMS rubbers and compared to tensile tests on the same specimen. For the second step, tissue constructs are engineered using uniaxial or equibiaxial static constrained culture conditions. Digital image correlation (DIC) is used to quantify the anisotropy in the tissue constructs. For both validation steps, material parameters are estimated by inverse fitting of a computational model to the experimental results.     

Robust strategies for automated AFM force curve analysis-II: adhesion-influenced indentation of soft, elastic materials

In the first of this two-part discourse on the extraction of elastic properties from atomic force microscopy (AFM) data, a scheme for automating the analysis of force-distance curves was introduced and experimentally validated for the Hertzian (i.e., linearly elastic and noninteractive probe-sample pairs) indentation of soft, inhomogeneous materials. In the presence of probe-sample adhesive interactions, which are common especially during retraction of the rigid tip from soft materials, the Hertzian models are no longer adequate. A number of theories (e.g., Johnson-Kendall-Roberts and Derjaguin-Muller-Toporov), covering the full range of sample compliance relative to adhesive force and tip radius, are available for analysis of such data. We incorporated Pietrement and Troyon's approximation (2000, "General Equations Describing Elastic Indentation Depth and Normal Contact Stiffness Versus Load," J. Colloid Interface Sci., 226(1), pp. 166-171) of the Maugis-Dugdale model into the automated procedure. The scheme developed for the processing of Hertzian data was extended to allow for adhesive contact by applying the Pietrement-Troyon equation. Retraction force-displacement data from the indentation of polyvinyl alcohol gels were processed using the customized software. Many of the retraction curves exhibited strong adhesive interactions that were absent in extension. We compared the values of Young's modulus extracted from the retraction data to the values obtained from the extension data and from macroscopic uniaxial compression tests. Application of adhesive contact models and the automated scheme to the retraction curves yielded average values of Young's modulus close to those obtained with Hertzian models for the extension curves. The Pietrement-Troyon equation provided a good fit to the data as indicated by small values of the mean-square error. The Maugis-Dugdale theory is capable of accurately modeling adhesive contact between a rigid spherical indenter and a soft, elastic sample. Pietrement and Troyon's empirical equation greatly simplifies the theory and renders it compatible with the general automation strategies that we developed for Hertzian analysis. Our comprehensive algorithm for automated extraction of Young's moduli from AFM indentation data has been expanded to recognize the presence of either adhesive or Hertzian behavior and apply the appropriate contact model.     

Robust strategies for automated AFM force curve analysis--I. Non-adhesive indentation of soft, inhomogeneous materials

The atomic force microscope (AFM) has found wide applicability as a nanoindentation tool to measure local elastic properties of soft materials. An automated approach to the processing of AFM indentation data, namely, the extraction of Young's modulus, is essential to realizing the high-throughput potential of the instrument as an elasticity probe for typical soft materials that exhibit inhomogeneity at microscopic scales. This paper focuses on Hertzian analysis techniques, which are applicable to linear elastic indentation. We compiled a series of synergistic strategies into an algorithm that overcomes many of the complications that have previously impeded efforts to automate the fitting of contact mechanics models to indentation data. AFM raster data sets containing up to 1024 individual force-displacement curves and macroscopic compression data were obtained from testing polyvinyl alcohol gels of known composition. Local elastic properties of tissue-engineered cartilage were also measured by the AFM. All AFM data sets were processed using customized software based on the algorithm, and the extracted values of Young's modulus were compared to those obtained by macroscopic testing. Accuracy of the technique was verified by the good agreement between values of Young's modulus obtained by AFM and by direct compression of the synthetic gels. Validation of robustness was achieved by successfully fitting the vastly different types of force curves generated from the indentation of tissue-engineered cartilage. For AFM indentation data that are amenable to Hertzian analysis, the method presented here minimizes subjectivity in preprocessing and allows for improved consistency and minimized user intervention. Automated, large-scale analysis of indentation data holds tremendous potential in bioengineering applications, such as high-resolution elasticity mapping of natural and artificial tissues.     

Non-ideal effects in indentation testing of soft tissues

Indentation has several advantages as a loading mode for determining constitutive behavior of soft, biological tissues. However, indentation induces a complex, spatially heterogeneous deformation field that creates analytical challenges for the calculation of constitutive parameters. As a result, investigators commonly assume small indentation depths and large sample thicknesses to simplify analysis and then restrict indentation depth and sample geometry to satisfy these assumptions. These restrictions limit experimental resolution in some fields, such as brain biomechanics. However, recent experimental evidence suggests that conventionally applied limits are in fact excessively conservative. We conducted a parametric study of indentation loading with various indenter geometries, surface interface conditions, sample compressibility, sample geometry and indentation depth to quantitatively describe the deviation from previous treatments that results from violation of the assumptions of small indentation depth and large sample thickness. We found that the classical solution was surprisingly robust to violation of the assumption of small strain but highly sensitive to violation of the assumption of large sample thickness, particularly if the indenter was cylindrical. The ramifications of these findings for design of indentation experiments are discussed and correction factors are presented to allow future investigators to account for these effects without recreating our finite element models.     

Mechanical characterization of anisotropic planar biological soft tissues using large indentation: a computational feasibility study

Traditionally, the complex mechanical behavior of planar soft biological tissues is characterized by (multi)axial tensile testing. While uniaxial tests do not provide sufficient information for a full characterization of the material anisotropy, biaxial tensile tests are difficult to perform and tethering effects limit the analyses to a small central portion of the test sample. In both cases, determination of local mechanical properties is not trivial. Local mechanical characterization may be performed by indentation testing. Conventional indentation tests, however, often assume linear elastic and isotropic material properties, and therefore these tests are of limited use in characterizing the nonlinear, anisotropic material behavior typical for planar soft biological tissues. In this study, a spherical indentation experiment assuming large deformations is proposed. A finite element model of the aortic valve leaflet demonstrates that combining force and deformation gradient data, one single indentation test provides sufficient information to characterize the local material behavior. Parameter estimation is used to fit the computational model to simulated experimental data. The aortic valve leaflet is chosen as a typical example. However, the proposed method is expected to apply for the mechanical characterization of planar soft biological materials in general.     

Elastic response, buckling, and instability of microtubules under radial indentation

We tested the mechanical properties of single microtubules by lateral indentation with the tip of an atomic force microscope. Indentations up to approximately 3.6 nm, i.e., 15% of the microtubule diameter, resulted in an approximately linear elastic response, and indentations were reversible without hysteresis. At an indentation force of around 0.3 nN we observed an instability corresponding to an approximately 1-nm indentation step in the taxol-stabilized microtubules, which could be due to partial or complete rupture of a relatively small number of lateral or axial tubulin-tubulin bonds. These indentations were reversible with hysteresis when the tip was retracted and no trace of damage was observed in subsequent high-resolution images. Higher forces caused substantial damage to the microtubules, which either led to depolymerization or, occasionally, to slowly reannealing holes in the microtubule wall. We modeled the experimental results using finite-element methods and find that the simple assumption of a homogeneous isotropic material, albeit structured with the characteristic protofilament corrugations, is sufficient to explain the linear elastic response of microtubules.     

Versatile and inexpensive Hall-Effect force sensor for mechanical characterization of soft biological materials

Mismatch of hierarchical structure and mechanical properties between tissue-engineered implants and native tissue may result in signal cues that negatively impact repair and remodeling. With bottom-up tissue engineering approaches, designing tissue components with proper microscale mechanical properties is crucial to achieve necessary macroscale properties in the final implant. However, characterizing microscale mechanical properties is challenging, and current methods do not provide the versatility and sensitivity required to measure these fragile, soft biological materials. Here, we developed a novel, highly sensitive Hall-Effect based force sensor that is capable of measuring mechanical properties of biological materials over wide force ranges (μN to N), allowing its use at all steps in layer-by-layer fabrication of engineered tissues. The force sensor design can be easily customized to measure specific force ranges, while remaining easy to fabricate using inexpensive, commercial materials. Although we used the force sensor to characterize mechanics of single-layer cell sheets and silk fibers, the design can be easily adapted for different applications spanning larger force ranges (>N). This platform is thus a novel, versatile, and practical tool for mechanically characterizing biological and biomimetic materials.     

Contact analysis of biphasic transversely isotropic cartilage layers and correlations with tissue failure.

Failure of articular cartilage has been investigated experimentally and theoretically, but there is only partial agreement between observed failure and predicted regions of peak stresses. Since trauma and repetitive stress are implicated in the etiopathogenesis of osteoarthritis, it is important to develop cartilage models which correctly predict sites of high stresses. Cartilage is anisotropic and inhomogeneous, though it has been difficult to incorporate these complexities into engineering analyses. The objectives of this study are to demonstrate that a transversely isotropic, biphasic model of cartilage can provide agreement between predicted regions of high stresses and observed regions of cartilage failure and that with transverse isotropy cartilage stresses are more sensitive to convexity and concavity of the surfaces than with isotropy. These objectives are achieved by solving problems of diarthrodial joint contact by the finite-element method. Results demonstrate that transversely isotropic models predict peak stresses at the cartilage surface and the cartilage-bone interface, in agreement with sites of fissures following impact loading; isotropic models predict peak stresses only at the cartilage-bone interface. Also, when convex cartilage layers contacted concave layers in this study, the highest tensile stresses occur in the convex layer for transversely isotropic models; no such differences are found with isotropic models. The significance of this study is that it establishes a threshold of modeling complexity for articular cartilage that provides good agreement with experimental observations under impact loading and that surface curvatures significantly affect stress and strain within cartilage when using a biphasic transversely isotropic model.     

A transversely isotropic biphasic model for unconfined compression of growth plate and chondroepiphysis.

Using the biphasic theory for hydrated soft tissues (Mow et al., 1980) and a transversely isotropic elastic model for the solid matrix, an analytical solution is presented for the unconfined compression of cylindrical disks of growth plate tissues compressed between two rigid platens with a frictionless interface. The axisymmetric case where the plane of transverse isotropy is perpendicular to the cylindrical axis is studied, and the stress-relaxation response to imposed step and ramp displacements is solved. This solution is then used to analyze experimental data from unconfined compression stress-relaxation tests performed on specimens from bovine distal ulnar growth plate and chondroepiphysis to determine the biphasic material parameters. The transversely isotropic biphasic model provides an excellent agreement between theory and experimental results, better than was previously achieved with an isotropic model, and can explain the observed experimental behavior in unconfined compression of these tissues.     

A transversely isotropic hyperelastic constitutive model of the PDL. Analytical and computational aspects

This study describes the development of a constitutive law for the modelling of the periodontal ligament (PDL) and its practical implementation into a commercial finite element code. The constitutive equations encompass the essential mechanical features of this biological soft tissue: non-linear behaviour, large deformations, anisotropy, distinct behaviour in tension and compression and the fibrous characteristics. The approach is based on the theory of continuum fibre-reinforced composites at finite strain where a compressible transversely isotropic hyperelastic strain energy function is defined. This strain energy density function is further split into volumetric and deviatoric contributions separating the bulk and shear responses of the material. Explicit expressions of the stress tensors in the material and spatial configurations are first established followed by original expressions of the elasticity tensors in the material and spatial configurations. As a simple application of the constitutive model, two finite element analyses simulating the mechanical behaviour of the PDL are performed. The results highlight the significance of integrating the fibrous architecture of the PDL as this feature is shown to be responsible for the complex strain distribution observed.     

The dependence of transversely isotropic elasticity of human femoral cortical bone on porosity

The objective of this study was to examine the dependence of the elastic properties of cortical bone as a transversely isotropic material on its porosity. The longitudinal Young's modulus, transverse Young's modulus, longitudinal shear modulus, transverse shear modulus, and longitudinal Poisson's ratio of cortical bone were determined from eighteen groups of longitudinal and transverse specimens using tensile and torsional tests on a servo-hydraulic material testing system. These cylindrical waisted specimens of cortical bone were harvested from the middle diaphysis of three pairs of human femora. The porosity of these specimens was assessed by means of histology. Our study demonstrated that the longitudinal Young's and shear moduli of human femoral cortical bone were significantly (p<0.01) negatively correlated with the porosity of cortical bone. Conversely, the elastic properties in the transverse direction did not have statistically significant correlations with the porosity of cortical bone. As a result, the transverse elastic properties of cortical bone were less sensitive to changes in porosity than those in the longitudinal direction. Additionally, the anisotropic ratios of cortical bone elasticity were found to be significantly (p<0.01) negatively correlated with its porosity, indicating that cortical bone tended to become more isotropic when its porosity increased. These results may help a number of researchers develop more accurate micromechanics models of cortical bone.     

Static indentation of anisotropic biomaterials using axially asymmetric indenters--a computational study

Indentation has historically been used by biomechanicians to extract the small strain elastic or viscoelastic properties of biological tissues. Because of the axisymmetry of indenters used in these studies however, analysis of the results requires the assumption of material isotropy and often yields an "effective" elastic modulus. Since most biological tissues such as bone and myocardium are known to be anisotropic, the use of conventional indentation techniques for estimating material properties is therefore limited. The feasibility of using an axially asymmetric indenter to determine material directions and in-plane material properties for anisotropic tissue is explored here using finite element analysis. The load versus displacement curves as would be measured by an indenter depend on the orientation of the indenter cross section relative to the in-plane material axes, thus suggesting a method for determining the underlying material directions. Additionally, the stiffness of the tissue response to indentation is sensitive to the values of the in-plane anisotropic material properties and prestretches, and thus test results can be used to back out relevant constitutive parameters.     

A transversely isotropic constitutive model of excised guinea pig spinal cord white matter

Narrowing of the spinal canal generates an amalgamation of stresses within the spinal cord parenchyma. The tissue’s stress state cannot be quantified experimentally; it must be described using computational methods, such as finite element analysis. The objective of this research was to propose a compressible, transversely isotropic constitutive model, an augmentation of the isotropic Mooney–Rivlin hyperelastic strain energy function, to describe the guinea pig spinal cord white matter. Model parameters were derived from a combination of inverse finite element analysis on transverse compression experiments and least squared error analysis applied to quasi-static longitudinal tensile tests. A comparison of the residual errors between the predicted response and the experimental measurements indicated that the transversely isotropic constitutive law that incorporates an offset stretch reduced the error by a factor of four when compared to other commonly used models.     

Confined and unconfined stress relaxation of cartilage: appropriateness of a transversely isotropic analysis.

Previous studies have shown that stress relaxation behavior of calf ulnar growth plate and chondroepiphysis cartilage can be described by a linear transverse isotropic biphasic model. The model provides a good fit to the observed unconfined compression transients when the out-of-plane Poisson's ratio is set to zero. This assumption is based on the observation that the equilibrium stress in the axial direction (deltaz) is the same in confined and unconfined compression, which implies that the radial stress deltar = 0 in confined compression. In our study, we further investigated the ability of the transversely isotropic model to describe confined and unconfined stress relaxation behavior of calf cartilage. A series of confined and unconfined stress relaxation tests were performed on calf articular cartilage (4.5 mm diameter, approximately 3.3 mm height) in a displacement-controlled compression apparatus capable of measuring delta(z) and delta(r). In equilibrium, delta(r) > 0 and delta(z) in confined compression was greater than in unconfined compression. Transient data at each strain were fitted by the linear transversely isotropic biphasic model and the material parameters were estimated. Although the model could provide good fits to the unconfined transients, the estimated parameters overpredicted the measured delta(r). Conversely, if the model was constrained to match equilibrium delta(r), the fits were poor. These findings suggest that the linear transversely isotropic biphasic model could not simultaneously describe the observed stress relaxation and equilibrium behavior of calf cartilage.     

Comparison of single-phase isotropic elastic and fibril-reinforced poroelastic models for indentation of rabbit articular cartilage

Classically, single-phase isotropic elastic (IE) model has been used for in situ or in vivo indentation analysis of articular cartilage. The model significantly simplifies cartilage structure and properties. In this study, we apply a fibril-reinforced poroelastic (FRPE) model for indentation to extract more detailed information on cartilage properties. Specifically, we compare the information from short-term (instantaneous) and long-term (equilibrium) indentations, as described here by IE and FRPE models. Femoral and tibial cartilage from rabbit (age 0–18 months) knees (n=14) were tested using a plane-ended indenter (diameter=0.544 mm). Stepwise creep tests were conducted to equilibrium. Single-phase IE solution for indentation was used to derive instantaneous modulus and equilibrium (Young's) modulus for the samples. The classical and modified Hayes’ solutions were used to derive values for the indentation moduli. In the FRPE model, the indentation behavior was sample-specifically described with three material parameters, i.e. fibril network modulus, non-fibrillar matrix modulus and permeability. The instantaneous and fibril network modulus, and the equilibrium Young's modulus and non-fibrillar matrix modulus showed significant (p<0.01) linear correlations of R²=0.516 and 0.940, respectively (Hayes’ solution) and R²=0.531 and 0.960, respectively (the modified Hayes’ solution). No significant correlations were found between the non-fibrillar matrix modulus and instantaneous moduli or between the fibril network modulus and the equilibrium moduli. These results indicate that the instantaneous indentation modulus (IE model) provides information on tensile stiffness of collagen fibrils in cartilage while the equilibrium modulus (IE model) is a significant measure for stiffness of PG matrix. Thereby, this study highlights the feasibility of a simple indentation analysis.