Finite element models of human trabecular bone from microCT images in the Ciclope open-source package: convergence and performance of a tetrahedral mesh compared to the standard voxel-based mesh

Finite element (FE) models derived from microcomputed tomography (μCT) images are often used to non-invasively estimate the mechanical properties of trabecular bone. Traditionally, image 3D voxels are directly converted into hexahedral elements, generating voxel-FE models that are reliable but computationally very demanding. Tetrahedral FE models (tetra-FE) offer greater computational efficiency but are less frequently used, likely due to complexities in mesh parameter definition and risk of element distortion. Ciclope is an open-source package, written in Python, that facilitates the preprocessing of μCT data to generate voxel-FE and tetra-FE. This study aimed to verify the robustness and accuracy of the computationally efficient tetra-FE approach, through a convergence test on the main mesh parameters and a comparison of results with the standard voxel-FE approach. Forty-two cylindrical human trabecular bone samples (diameter = 5 mm, height = 10 mm, bone volume fraction range 8%-38%) were extracted from larger samples imaged using μCT (19.5 μm voxel size). Voxel- and tetra-FE models were generated using the Ciclope workflow. Common preprocessing steps included inversion of image greyscale, bone segmentation via a predetermined threshold, and removal of unconnected voxel clusters. Voxel-FE models (8-node hexahedra) were then generated and solved. Four-node tetra-FE models were generated, performing convergence tests that sequentially assessed the impact of element size and distortion. Bone was modelled as a linear elastic, homogeneous material (E=18000 MPa, ν=0.3), simulating static compression through vertical displacement (0.4% shortening) of the upper base and full constraint of the lower base. The target outcome was the Apparent Young's Modulus (Eapp), calculated as: Eapp = (F/A)/ε where F is the total reaction force, A is the cross-sectional area, and ε is the nominal strain. Both Relative (RPV, against the next finer mesh) and absolute (APV, against the finest mesh) percentage variations of Eapp were adopted as convergence metrics. Correlation and discrepancies in results between tetra-FE and voxel-FE were used to analyse the comparability between the two models. Convergence tests indicated that a mesh size factor 1.5 times the voxel size and a maximum element circumradius (chosen as element distortion index) of 3 times the voxel size are needed to produce APV and RPV variations lower than 2%. Tetra-FE models configured for 2% convergence significantly reduced computational cost (1/30 of degrees of freedom and 1/6 of elements) and exhibited an almost ideal correlation with the reference voxel-FE models (Eapp-tetra = 1.07*Eapp-voxel + 13.47, R² = 0.998) but slightly overestimated Eapp (11±3%). This overestimation could not be explained by a volumetric offset of the tetrahedral meshing, which resulted negligible. Despite the rather fine meshing adopted, the known stiffening effect of tetrahedral elements was instead identified as the main responsible for Eapp overestimation, as preliminary analyses conducted with quadratic (10-node) tetra-FE reduced Eapp discrepancy with voxel-FE to less than 2%. In summary, this work demonstrated that tetra-FE models of trabecular bone from μCT images can be a valid alternative to the standard voxel-FE models for the estimate of bone apparent elastic modulus. Future efforts will be devoted to the development and full characterisation of 10-node tetra-FE models to avoid the slight Eapp overestimation seen in tetra-FE models.