Introduction: Fractures of the tibial shaft are relatively common injuries. There are indications that tibial shaft fractures share characteristics in terms of site, type and local fracture mechanisms. In this study, we aimed to set up a mathematical, computer-based model using finite elem-ent analysis of the bones of the lower leg to examine if such a model is adequate for prediction of fracture locations and patterns. In future studies, we aim to use these biomechanical results to examine fracture prevention, among others, and to simulate different types of osteosynthesis and the process of bony healing. The biomechanical results are the basis for fracture healing, biomechanical fall analysis and stability analysis of osteosynthesis.
Material and methods: A finite element model of the bony part of the lower leg was generated on the basis of computed tomography data from the Visible Human Project. The data consisted of 21,219 3D elements with a cortical shell and a trabecular core. Three types of load of torsion, a direct lateral load and axial compression were applied.
Results: The finite element linear static analysis resulted in relevant fracture localizations and indicated relevant fracture patterns.
Conclusion: In the present study, we have successfully simulated fracture mechanisms, obtained adequate fracture locations and achieved an indication of the fracture morphology. The method of fracture simulation employed showed good agreement with known clinical data and data from prior mechanical testing. This substantiates the validity of fracture simulation for future studies examining tibial fractures, fracture healing and prevention.
Fractures of the tibial shaft are relatively common injuries, when e.g. front seat occupants are involved in frontal crashes . Although the fracture mechanism is predominantly caused by axial compression of the leg, a substantial bending moment develops in the tibial shaft due to its natural curvature and to transverse loading when the tibial shaft collides with the dashboard . This can lead to tibial shaft fractures . In general, experimental crash tests have shown that tibial fractures may be expected when the sum of the normal stresses due to axial compression and the resultant bending moment exceed the strength of the tibial cortical bone (tibia index) . In a review study of tibial shaft fractures, 88.9% were located in middle third or distal third part of the shaft . These fractures were caused by a variety of occurrences including falls, various sports, pedestrian crashes, etc. All the studies referred to above gave an indication that tibial shaft fractures share a number of features in terms of site, type and local fracture mechanisms.
In this study, we aimed to set up an advanced, mathematical computer-based model using finite element (FE) analysis of the bones of the lower leg to examine whether an adequate computer model for predicting lower leg fracture initiation could be generated. We wanted to substantiate the validity of the model by predicting the fracture location to get an indication of the morphology of tibial shaft fractures for various common mechanical fracture loads (axial compression, torsion and bending).
MATERIAL AND METHODS
Finite element model
FE analyses were performed using the finite element code COSMOS/M. The three-dimensional "lower leg" FE model (LL FE model) used included the tibia and the fibula. The bony structures were generated by segmentation of a data set of computed tomographies from the Visible Human Project . Segmentation is a data reduction tool used to simplify the geometrical shape of the bony structures of the structure without losing accuracy. Segmentation was necessary for the FE code to handle the large amounts of data involved. The cortical bone was simulated using tetrahedral elements (TETRA 4). The LL FE model comprised a total of 21,219 volumetric elements. The material property for trabecular bone was a uniform distribution of Young’s modulus of 100 MPa and a Poisson’s ratio of 0.2. Correspondingly, Young’s modulus for cortical bone was 12,000 MPa and a Poisson’s ratio of 0.3 . The thicknesses of cortical and trabecular bone were taken from the Visible Human . The interosseos membrane and the tibiofibular joints were simulated using binding of the proximal and distal contact areas of the bones.
The fracture locations for this particular FE model were determined as the superficial areas of maximal Von Mises stress (σ Von Mises) in the cortical bone. Moreover, the σVon Mises were compared with a ratio of failure (RF) between σVon Mises and ultimate stress (σlim ) as in Scheilo et al . The RF had to exceed 1 (RF > 1) for failure. The RF was determined as:
The σlim was determined in a linear relationship with the ash density in accordance with Keyak et al :
Two different ash densities were determined from two different density-elasticity relationships between the apparent density by Carter and Hayes  and Scheilo et al . Ash density was obtained from a ratio of 0.6 between ash and the apparent density, within the range 0.55-0.63 identified in the literature [8, 10]. The two density-elasticity relationships were specific for the hip region, but the results for both density-elasticity relationships with regard to the ultimate stress were similar (106 and 110 MPa, respectively). These results coincide with experimental results for ultimate stress for the tibial cortical bone, which is approximately 110 MPa for the tibia cortical bone with a yield stress amounting to 100 MPa .
Loading and boundary conditions
Three types of load were applied to the LL FE model:
1. Axial compression load. A total load of 800 N was distributed on the upper surface of the lower leg FE model. The compression load was applied as 60% to the medial side and 40% to the lateral side , as the knee contact force is not uniformly shared between the condyles of the tibia . The lower end of the lower leg FE model was completely immobilized in all directions.
2. Torsional load. Torsional loads of 1,000 N in external rotation were applied to the upper and lower ends of the LL FE model. We were not able to obtain data with values for the torsional loads by which fractures to the tibia are known to occur, so the loads were assumed to be an adequate rotational moment (for the sake of convenience). The middle part of the tibial shaft was immobilized in the axial direction for analysis stability.
3. Transverse load.. Bending load was applied as nodal loads to the middle parts of the tibial and fibular shafts. A 5,000 N load was applied to the tibia and 2,500 N was applied to the fibula. The upper and lower ends of the LL FE model were completely immobilized in all directions.
Static analyses were performed for all three loading cases. Fracture locations and patterns were estimated by evaluating the maximal nodal Von Mises stresses as described previously. The LL FE model is illustrated in Figure 1.
Convergence test for the FE model was performed to ensure that the FE model had an appropriate number of elements. During this test, the number of elements in the FE model was increased until the point at which the calculated results converge to one exact solution, thus yielding the appropriate number of elements for in the FE model.
The distribution of Von Mises stress of the torsional load is illustrated in Figure 2.
Maximal Von Mises stress was located in the distal third of the lower leg, and Von Mises stresses were distributed in a spiral pattern. The RF was 0.46.
The distribution of Von Mises stress of the axial load is illustrated in Figure 3.
Maximal Von Mises stress was located in the middle third of the lower leg, and Von Mises stresses were distributed in an oblique pattern. The RF was 0.26.
The distribution of Von Mises stress of the lateral load is illustrated in Figure 4.
Maximal Von Mises stress was located in the middle third of the lower leg, and Von Mises stresses were distributed in a transverse pattern. The RF was 0.46.
The purpose of the present study was to evaluate the use of the maximal Von Mises stress calculated by finite element analysis as a measure for examination of the fracture mechanisms of the lower leg. We included experimental material data and geometric data of the bones of the lower leg and were able to simulate adequate fracture localizations and to obtain an indication of fracture morphology of the most common lower leg fractures . We applied adequate fracture loads with regard to direction and were able to generate fractures located in the middle and distal third of the lower leg, as expected. Moreover, the fracture patterns generated by the Von Mises stress distribution showed adequate spiral, transverse and oblique type patterns for torsional, lateral direct and axial loads, respectively. While static analysis may not be the best way to obtain true predictions of fracture patterns, the von Mises stress distribution method, on the other hand, yields clinically correct locations of fracture initiations prior to initial failure. After a fracture, stress redistribution at the fracture site would be expected, which would influence subsequent fracture patterns. Static FE analysis gives an estimate of the location of the fracture onset and the level of fracture risk. In this study, RF values for the various loads were somewhat lower than expected, which is presumably so because the estimated loads were too small. To the best of our knowledge, no in vivo data are available on this subject. However, certain in vitro mechanical test data on the whole tibia bone indicate that we have probably underestimated the loads – at least for axial compression (which had the smallest RF) – by a factor 4 . Furthermore, studies have indicated that other mechanical signals, such as maximal principal strain or distorsion energy [7, 14], may be better suited for fracture analysis. This may also explain some of the discrepancy in RF values. However, the generated FE results were substantiated by known clinical and experimental studies with regard to locations and patterns [2, 3], thereby indicating that simulation of tibial shaft fractures is a valid path for further research, i.e. in fracture prevention, simulation of different types of osteosynthesis and the process of bony healing. We acknowledge the need for further improvement of the LL FE model as outlined above. To summarize, we have performed a general evaluation of the FE model and the computed results of a reaction force balance and the performed convergence tests, but have not correlated the computed results with an experimental cadaver test. The present LL FE model was based on geometrical data from one test subject (from The Visible Human Project) and we acknowledge that the computed stress results, and thus the fracture morphology and site, may be influenced by anatomical variation, i.e. the curvature of the tibia.
In this study, a linear analysis was chosen since the aim was not to investigate the whole fracture process, but to examine if we could reproduce mechanisms of fracture at onset. The use of static analysis is justifiable, because bone has been observed to have a quasi-brittle failure mode when tested experimentally . Earlier studies have also indicated that even if bone is considered to be of anisotropic or orthotropic nature, it is justified to use isotropic material properties for maximal Von Mises stress analysis for long bones and also for uniaxial yield strain for the lower leg [15, 16]. However, the assumption of isotropic behaviour may have influenced the computed maximal Von Mises stress distribution, and thereby the fracture estimation. All in all, when considering the relative crudity of the LL FE model and the estimated loads, the use of isotropic material properties instead of anisotropy seems justifiable taking into account the lack of overall knowledge of anisotropic fracture properties of bone in FE models .
In this study, we excluded the effects of lower leg muscles. It has been demonstrated that the calf muscles influence the forward bowing (bending) of the tibia, thus affecting fracture development . However, the muscles were not included, since – in the context of a fracture – the effects caused by muscles are currently unknown. The tibiofibular joints and the interosseos membrane were simulated by binding the tibiofibular joints. Again, appropriate mechanical properties have not been reported, hence justifying the crude simulation of these parts. In conclusion, we acknowledge that we would need to address the above mentioned considerations in future studies. However, in the present study we have been able to simulate fracture mechanisms, obtain adequate fracture locations and get an indication of fracture morphology. This indicates that the LL FE model may be used to examine e.g. fracture prevention and to simulate different types of osteosynthesis and the process of bony healing in future studies.
Correspondence: Christian Wong, Department of Orthopaedics, University Hospital of Hvidovre, 2650 Hvidovre, Denmark. E-mail: email@example.com
Accepted: 24 February 2010
Conflicts of interest: None
Taylor A, Morris A, Thomas P et al. Mechanisms of lower extremity injuries to front seat car occupants—an in depth accident analysis. Proceedings of the International IRCOBI Conference on the Biomechanics of Impact 1991:53-72.
Funk JR, Rudd RW, Kerrigan JR et al. The effect of tibial curvature and fibular loading on the tibia index. Traffic Inj Prev 2004;5:164-72.
Ivarsson BJ, Manaswi A, Genovese D et al. Site, type, and local mechanism of tibial shaft fracture in drivers in frontal automobile crashes. Forensic Sci Int 2008;175:186-92.
Mertz HJ. Anthropometric test devices. I: Nahum AM, Melvin JW, eds. Accidental injury: Biomechanics and prevention. New York: Springer-Verlage, 1993:292-310
Ackerman MJ, Spitzer VM, Scherzinger AL et al. The Visible Human data set: an image resource for anatomical visualization. Medinfo 1995;8 Pt 2:1195-8.
Shirazi-Adl SA, Shrivastava SC, Ahmed AM. Stress analysis of the lumbar disc-body unit in compression. A three-dimensional nonlinear finite element study. Spine 1984;9:120-34.
Schileo E, Taddei F, Cristofolini L et al. Subject-specific finite element models implementing a maximum principal strain criterion are able to estimate failure risk and fracture location on human femurs tested in vitro. J Biomech 2008;41:356-67.
Keyak JH, Lee IY, Skinner HB. Correlations between orthogonal mechanical properties and density of trabecular bone: use of different densitometric measures. J Biomed Mater Res 1994;28:1329-36.
Carter DR, Hayes WC. The compressive behavior of bone as a two-phase porous structure. J Bone Joint Surg Am 1977;59:954-62.
Goulet RW, Goldstein SA, Ciarelli MJ et al. The relationship between the structural and orthogonal compressive properties of trabecular bone. J Biomech 1994;27:375-89.
Ebacher V, Tang C, McKay H et al. Strain redistribution and cracking behavior of human bone during bending. Bone 2007;40:1265-75.
Duda GN, Mandruzzato F, Heller M et al. Mechanical boundary conditions of fracture healing: borderline indications in the treatment of unreamed tibial nailing. J Biomech 2001;34:639-50.
Morrison JB. The mechanics of the knee joint in relation to normal walking. J Biomech 1970;3:51-61.
Keyak JH, Rossi SA, Jones KA et al. Prediction of femoral fracture load using automated finite element modeling. J Biomech 1998;31:125-33.
Peng L, Bai J, Zeng X et al. Comparison of isotropic and orthotropic material property assignments on femoral finite element models under two loading conditions. Med Eng Phys 2006;28:227-33.
Chang WC, Christensen TM, Pinilla TP et al. Uniaxial yield strains for bovine trabecular bone are isotropic and asymmetric. J Orthop Res 1999;17: 582-5.
Devas MB. Stress fractures of the tibia in athletes or shin soreness. J Bone Joint Surg Br 1958;40-B:227-39.