Computational Modelling of Bone Pathologies

Fracture healing is a complex regenerative process that gradually restores the functional and mechanical properties of bone. However, at a rate of about 5%-10% of the millions of fractures that occur annually, the healing process is not successfully completed within a few months, but shows a significant delay called delayed union or stops at an intermediate stage called non-union.

Nowadays, ultrasound is considered as an innovative tool for the evaluation of the bone status and health. The research interest has been mainly focused on the assessment of osteoporosis, while recently the application of ultrasound to the monitoring of the healing process has attracted many research groups worldwide.

Two and three dimensional computational models of healing long bones have been developed, and the performance of experimental measurements on sheep tibial models and the time-frequency representation of propagating guided waves towards quantitatively charactering the progress of healing has been studied.

Bone modelling was initially approached based on simplified geometries and using the classical theory of linear elasticity. The Finite Difference Method (FDM) was applied for the numerical solution of the wave propagation problem. Moreover, the use of guided waves was incorporated based on the theory that describes guided waves in plates (Lamb Theory). In more recent studies, the main objective has been the realistic modelling of the geometry and microstructure of the bone and callus tissues, the incorporation of enhanced elasticity theories and the use of Boundary Element methods for solving the wave propagation problem. In this direction, we have analytically derived dispersion curves for guided wave propagating in bones by incorporating higher order gradient theories and especially on the general theory proposed by Mindlin (Form II). Among the most recent achievements is the realistic reconstruction of healing bones that can reflect accurately bone geometry and callus heterogeneity by using Scanning Acoustic Microscopy (SAM) images.

SAM is a relatively new imaging technique that offers not only a detailed representation of the bone and callus microstructure, but also calibrated maps of the material properties (acoustic impedance). Therefore, by exploiting the acoustic impedance and measuring (or estimating) the density, the calibrated maps of the elastic coefficient c33 can be determined. Towards this, significant collaboration with the German research group of the Julius Wolff Institute (Berlin Brandenburg School for Regenerative Therapies, Charite-Universitatsmedizin Berlin) which has presented an extensive scientific work related to SAM in bones has been established. SAM images are employed, and their calibrated maps of elastic properties as obtained from experimental studies on tibial osteotomies, in order to study how callus inhomogeneity, porosity and gradual consolidation affect wave propagation.

_{a)SAM image from a sheep tibia at the 9th postoperative week. b)Simulation of wave propagation in models reconstructed from SAM images.}

Trying to elucidate the impact of each individual bone-and callus-related parameter on wave propagation characteristics, wave dispersion and attenuation parameters in specific segments of healing bones are estimated. Porous media in bone and callus are approximated using an iterative methodology, called iterative effective medium approximation (IEMA). Numerical solution is achieved by using the BEM method and the results are also compared with analytical predictions.

_{Segment extracted from a callus region as obtained from SAM; 2) Binary where bone and fluid is represented, 3) IEMA effectrive medium and 4) corresponding model of image 1 where callus was described by 15 materials.}

The ongoing research on the ultrasonic assessment of bone fracture healing has resulted in more realistic, inhomogeneous and anisotropic computational bone models and the development of enhanced elasticity theories. All achievements and findings can be also found useful for the study of osteoporosis as well as be extended to other fields where ultrasound is used to assess the properties of composite materials.

Bone fracture healing is a complex regenerative process which advances in stages of callus (i.e. the soft tissue that develops at the fracture site) formation and consolidation. During each of the healing stages the material and geometrical properties of callus are bone gradually modified until the bone restores its structural integrity and mechanical function. Monitoring of fracture healing refers to the evaluation of the status of the healing bone, the detection of complications in the process, and the accurate assessment of the endpoint of healing. In clinical practice, evaluation of fracture healing is performed by serial clinical and radiographic examinations, both of which depend on the orthopaedic surgeon’s expertise and clinical judgment. Therefore, the development of more objective and quantitative means of monitoring the healing process will significantly assist in the early diagnosis of possible complications and the timely identification of the healing endpoint.

Several non-invasive quantitative methods for monitoring fracture healing have been reported in the literature including bone densitometry, vibrational analysis and the attachment of strain gauges to external fixation devices. Animal and clinical studies have also reported on the use of axial-transmission quantitative ultrasound (Fig.1) for the examination of fractures in long bones, such as tibia and radius [1-3,6]. A review of the ultrasonic bone monitoring literature and a presentation of the current status of knowledge can be found in [6].

_{Fig. 1: The axial transmission technique used for the ultrasonic evaluation of fracture healing in long bones}

In the first computational study of ultrasound propagation in a healing long bone, a 2-D model of an elastic isotropic plate was developed. A 2-mm fracture gap was modeled at middle of the plate’s length and the consolidation of callus was simulated by a simple 7-stage process. The callus tissue was assumed to be homogeneous and isotropic with properties evolving throughout the stages. Axial transmission of ultrasound was simulated by two longitudinal transducers placed in direct contact with the plate’s upper surface. The bone plate was assumed to be in vacuum (free plate) neglecting thus the presence of soft tissues.

This issue was addressed in a subsequent study in which the model of the healing bone was enhanced by assuming the callus tissue to be an inhomogeneous material consisting of six different ossification regions and by investigating different cases of boundary conditions to simulate the presence of overlying soft tissues (Fig. 2). The healing course was simulated by a three-stage process in which the properties of each region evolved corresponding to various soft tissue types that participate in the healing process.

_{Fig. 2: 2-D computational model of the healing bone and the transmitter-receiver configuration. (a) Model in which the bone is immersed in blood occupying the semi-infinite spaces, (b) model which has two soft layers, one corresponding to blood, whereas the other to blood bone marrow.}

In a subsequent computational study, we further extended our work by considering the 3-D geometry and anisotropy of the bone and the callus tissue. First, a hollow cylinder was examined for two cases of elastic symmetry: isotropy and transverse isotropy. Next, the real geometry of an intact tibia was modeled (Fig. 3). The callus tissue was again considered as an inhomogeneous material consisting of six ossification regions. The bone was considered free of tractions on the inner and outer surfaces. Solution to the 3-D elastic wave propagation problem was performed using the explicit elastodynamics finite element analysis.

_{Fig. 3: The model of the diaphyseal segment of cortical bone incorporating the fracture callus (sagittal section). The transmitter-receiver configuration is also illustrated.}

We performed traditional velocity measurements to investigate whether the nature of the first arriving signal (FAS) wave throughout the healing stages can provide quantitative assessment of the healing process. However our findings reveal that FAS velocity measurements are generally sub-optimal in terms of their sensitivity to alterations that occur in deep cortical layers. To this end, we further investigate the use of guided waves as an alternative means of monitoring bone fracture healing. A time-frequency (t-f) methodology was followed for the representation of the propagating guided modes. The method was first applied to the signals obtained from the free intact plate model. Mode identification was performed using group velocity dispersion curves predicted by the Lamb theory. It was demonstrated that the characteristics of the guided waves, such as the dispersion of velocity, are strongly affected by the irregular geometry and anisotropy of the cortical bone and callus, and by the boundary conditions induced by the soft tissues.

In all the works dealing with wave propagation in bones, the bone is mimicked as a linear elastic and homogenized plate or hollow cylinder. However, considering human bones as a linear elastic material with microstructure, its dynamic mechanical behaviour cannot be described adequately by the classical theory of linear elasticity, which is associated with concepts of homogeneity and locality of stresses. When the material exhibits a non homogeneous behaviour and microstructure is comparable to the macrostructure, microstructural effects become important and the state of stress has to be defined in a non-local manner. These microstructural effects can be successfully modelled in a macroscopic framework by employing enhanced elastic theories such as the couple stresses theory proposed by Cosserat brothers and generalized later by Eringen as micropolar elastic theory, the general higher order gradient elastic theory proposed by Mindlin and the non-local theory of elasticity of Eringen.

Higher-order gradient theories can be considered as generalizations of linear theory of elasticity, utilizing the displacement vector to describe the deformation of the continuum and introducing in both potential and kinetic energy higher order terms associated with internal length scale parameters that correlate microstructural effects with the macrostructural behaviour of the considered material. In our current study we adopt the simplest form of gradient theory (Mindlin FormII) to theoretically determine the group velocity dispersion curves of the guided modes propagating in isotropic bone-mimicking plates with microstructural effects. Two additional terms are now included in the constitutive equations representing the characteristic length of bone: a) the gradient coefficient g, introduced in the strain energy to balance the dimensions of strains and strain gradients and b) the micro-inertia term h, introduced in the kinetic energy. The plate was assumed free of stresses, as in the classical Lamb problem, and free of double stresses. The characteristic length h were assumed at the order of the osteons size. For each of the two cases for h, three sub-cases for were assumed, resulting totally in six different combinations between g and h. In the first two subcases g was assumed higher and smaller than h respectively, whereas in the third subcase the value of gwas assumed equal to that of g. The group velocity dispersion curves of guided waves were numerically obtained and compared with the corresponding ones of the Lamb modes. The results clearly indicate that when the two elastic constants have different values microstructure has a significant impact on mode dispersion by inducing both material and geometrical dispersion. The gradient theory of elasticity can provide supplementary information to better understand guided wave propagation in real bones.