1 Introduction

Computer graphics applied to medicine is a field of research with many open lines and continues to rise due to the advancement of new technologies. These advances have allowed, in the field of traumatology and bone fractures, the creation of computer assisted techniques that allow the reduction of intervention time and minimise the risk of error of specialists through the study of fractures.

The study of bone fractures is a challenging field due to the complexity of the hierarchical structure of bones. It is not possible to monitor a fracture in real time. Moreover, the issue of obtaining bone models for simulations means that most work focuses on bone at macroscopic level without taking into consideration the smaller scale levels. However, several studies of bone structures from a microscopic point of view suggest that they have a direct impact on all the phenomena that occur in bone [1,2,3,4,5].

The aim of this work is to propose a statistical approach for fracturing bone models from a microscopic point of view starting from a model at macroscale. The previous literature has focused on the study of the cross section and the use of 2D sketches to represent cortical tissue without longitudinal axis information. Therefore, to validate the results, our simulations have also been performed on the cross section of the bones using a 2D representation of the cortical bone. For this purpose, the process initiates with a 3D bone model and the selection of an area of the bone. In this selected area, a representation of the cortical region is generated to facilitate the simulation of the fracture. This step allows us to reproduce the clinical conditions of different patients. The approach for fracturing followed is similar to the one taken by other authors on the propagation of light on biological tissue proposed by Kumar et al. [6].

This work solves the main problem of simulations performed using the traditional method, the finite element method (FEM). In FEM simulations, a minimal portion of the cortical area of the bone is discretised due to the high computational resources required to carry out the fracture simulations. In these simulations, the mesh representing the cortical tissue at the microscopic level represents less than 5% of the cortical area of the bone in the cross section. Our microscopic representation of cortical tissue, used for fracture simulation, represents a complete section of the bone. Previous studies are based on a very specific area of a cross section and require a manual segmentation step. Moreover, the generated models enable simulations under various clinical conditions by replicating conditions of a real patient with different morphological characteristics in a matter of seconds. Additionally, the algorithm described in this work allow the realistic simulation of a fracture in the whole cortical area representation, which has been statistically validated and approved by a group of experts. It is possible to study how a fracture line interacts as it advances through the cortical tissue and how bone microstructures influence fracture propagation. The results obtained allow us to study the different situations that can occur during the growth of a fracture along the entire bone as well as generate examples for future studies. Moreover, deep learning-based models could be trained both to classify the different fracture patterns and to improve the accuracy of the fracture lines obtained. Moreover, deep learning based models could be trained with the result obtained both to classify the different fracture patterns and to improve the accuracy of the fracture lines obtained.

The structure of the article is as follows: firstly, a review of previous literature has been carried out to analyse the latest scientific advances in bone simulation and fracturing at microscale. The method followed for the simulation, as well as the parameters and values used, are described in the next section. Section 4 show the results obtained and they are discussed in the following. Finally, the conclusions and future work are presented.

2 Related Work

Sabet et al. [7] described the several levels of scale in which bones can be divided as well as the harmony they work to perform a variety of mechanical, biological and chemical functions. The hierarchical structure of bones provides high stiffness, strength, toughness and energy absorption. The structure of bones ranges from the macroscale to the nanoscale. The macroscale level corresponds to the representation of the whole bone. At mesoscale level, two types of bones can be distinguished: cortical and trabecular bone. The cortical bone is composed of a dense and hard shell that contributes mainly to its rigidity and toughness. However, the trabecular part of the bone fills the interior space of the bones built up by a porous network and it is the main responsible for absorbing energy and distributing the loads they bear through their porous structure. One level below corresponds to the microscale level which distinguishes small bone structures that make up both cortical and trabecular bone tissue. The principal bone structures that make up these tissues are osteons and trabeculae respectively. The remaining levels, down to the nanoscale, allow us to identify everything from collagen fibres to their assembly into larger structures known as lamellae but these structures are too small to be taken into consideration in the study of the fracture of bone models, although some authors have explored these aspects further [7].

Most of the literature to date has focused on the study of fractures at the macroscale level. These studies are usually aimed at healing the bone fracture through fracture reduction [8,9,10] and the analysis of bone health problems such as biodegradation or the corrosion rate [11]. Hambli et al. [12] designed a system to study the process of bone remodelling and the behaviour of a fracture through the use of continuum damage mechanics to predict the fracture pattern using the whole bone. These articles usually ignore the architecture of the bones at microscopic levels although it is scientifically documented that bone microstructures are involved in phenomena such as mechanical loads and fracture propagation. For example, fracture reduction simulations do not consist only of the joining of two separate pieces at macroscale level, it involves an internal process in which the microstructures that make up the bone and are visible in these lower hierarchical levels are regenerated through a process called remodelling [13, 14]. In addition, the bone microstructural elements are responsible for determining the path that the fracture follows during its propagation [7, 15,16,17,18]. Thus, the lack of analysis at the microscopic level of the bone makes the studies incomplete and potential for improvement, for example, by studying how to accelerate the remodelling process to speed up recovery after a fracture reduction.

Fig. 1
figure 1

Scheme for obtaining the microscopic model carried out by Demirtas et al. [2, 3] through manual segmentation. a Represents a sample transverse microscopy image of cortical bone, b is a detailed 2D sketch based on image a and c represents the 3D finite element model obtained by extruding the 2D sketch of b

The solutions presented in the literature for the fracturing of cortical bone tissue at microscopic levels are mostly based on FEM for the high reliability and accuracy approach in biomechanical applications [19]. These types of fractures can be classified as physical based fractures. Andreaus and Colloca [20] have studied the structural behaviour under physiological loadings and boundary conditions in intact and implanted human femurs obtaining as result important information about the stress that the bones receive in the most common movements. Abdel-Wahab et al. [17, 21] has analysed the propagation of a fracture through bone tissue under the influence of tissue heterogeneity as well as the influence of osteons. Some studies have developed a FEM-based approach which demonstrated the importance of osteons in determining fracture growth as well as the effect of microscale properties in assessing the risk of a whole bone fracture [22,23,24]. Gao et al. [25] focus on the behaviour of the cortical bone under dynamic loading. Budyn et al. [26] tried to implement a multiscale approach to experiment with fracture forces, from the microscale level, and study their influence on fracture growth. Strength and elastic modulus were found to be altered due to orientation variation in lamellae and existence of osteocyte lacunae [27]. Tang et al. [1] analysed the shear deformation and fracture of human cortical bone using a fluorescence microscopy and a device to damage the cortical tissue sample until it fractures. As a result, they conclude how the Haversian lamellae play an important role in the fracture growth process. Demirtas et al. [2, 3] focused on the study of the distribution and density of bone microstructures and how this affects fracturing. Other authors have studied the effect of features as the age or the gender on fracture propagation and bone modelling [4]. Microscale studies are necessary because from a macroscopic point of view it is not possible to determine the reason for varying the growth of a fracture line [5].

The previous literature reviewed in this section, based on FEM, uses images of segmented bone structures from microscope images as a starting point. The main problem is that the microscope images represent a very small region of the cortical tissue. Therefore, they are not very representative. The three-dimensional representations that can be found in the literature are extrusions of the two-dimensional models as shown in Fig. 1 and are not very common due to the lack of information on the longitudinal axis of the bone tissue. Therefore, this way of representing the tissue is inaccurate and difficult to verify as the information on the longitudinal axis is omitted. For that reason, our proposal is also based on the use of a two-dimensional representation. This work provides a solution to the problems that arise by allowing the fracturing of models that represent the entire cortical area realistically regardless of the number of bone microstructures that make up the tissue, as well as establishing a set of rules to determine the path that the fracture follows over the bone tissue. The parametric representation of cortical tissue used for fracturing has been validated through statistical analysis and visual assessment [28]. The nomenclature, abbreviations, units and parameters follow the recommendations of the American Society for Bone and Mineral Research (ASBMR) [29]. This recommendations also defines the meaning of each parameter. The parameters established are configurable so that the fracture and the simulations can be adapted to specific cases. Moreover, the results obtained with this method could be exported and used to generate data sets and allow a more accurate assessment and study of the progression of a fracture along the bone tissue.

3 Methods

The method proposed in this paper is based on the microscopic representation of cortical tissue proposed by Perez-Cano et al. [28]. In contrast to the referred literature shown in the previous section, the model used represents a complete cross section of a human bone. In this way, the fracture simulation does not focus on a specific area of the bone, it allows the fracture of the entire cortical bone. The method described in the following subsections describes a process and establishes a set of rules that determine the growth of a fracture. A flowchart that summarises the following steps can be seen in Fig. 2.

Fig. 2
figure 2

Flowchart of the complete process

3.1 Model

The fracture studies reviewed in the previous section are based on manual segmentation of bone structures through microscope images (Fig. 3). The issue with the use of these models is that the microscope images only contain the representation of small sections of the cortical tissue, so the simulations carried out on them are not accurate (Fig. 4). In this work, a representation of the entire cortical bone has been utilized. A representation of all the bone microstructures comprising the bone in a cross section has been proposed, starting from a macroscopic representation of the bone model, based in a previous study conducted by Pérez-Cano et al. [28].

Fig. 3
figure 3

Manual segmentation of cortical tissue in which osteons are delimited

Fig. 4
figure 4

Representation of the cortical bone. a Represents a cut in the diaphysis area of a femur to be analysed using a microscope. b Shows the cortical bone tissue view from an optical microscope while c shows the cortical bone after conducting a series of analyses using different laser intensities and microscope colours

The procedure to model the cortical tissue is composed of two stages. The first phase consists of obtaining a 3D representation of the bone order to determine the area in which to carry out the representation of the microstructures. Not all techniques are valid to represent the bone in 3D as it is necessary to represent the information of the internal part to delimit the cortical area of the bone. Therefore, the best option is the segmentation of CT images of real bones. A combination of the algorithm developed by Paulano et al. [30] and algorithms based on threshold and rank growth from the platform 3DSlicer [31] has been used to segment the bones. Simpleware Software, Invivo5, Romexis, Materialise Mimics or OsiriX [32] are other platforms that allow segmentation of bone models through medical imaging. However, 3D Slicer obtains better dimensional and geometric accuracy in the 3D reconstruction of bone models [33].

After obtaining a 3D representation of the bone model, the second step involves selecting a bone section through a plane. This enables the delimitation of the inner and outer boundary of the cortical area in which the microstructures of the bone tissue will be procedurally generated using the alpha shapes algorithm. The inner and outer boundary of the bone allows us to calculate the tissue volume parameter, through Gauss area formula [34], which calculates the area of irregular polygons. This parameter represents a section of cortical bone to generate the bone microstructures. Therefore, the outer and inner areas of the cortical area of the bone are calculated and subtracted to determine the area between both limits. In this area, the bone microstructures are generated using a set of estimated and input parameters. The osteonal density is the single input parameter. This parameter establishes a relationship between the mass and volume of all osteons in an area selected for the generation. Regarding the estimated parameters, these parameters have been determined from studies which analyse the microstructural characteristics of the bones. The model generator sets constant parameters for the size and shape of the microstructures, based on the average results of previous studies [18, 35,36,37,38,39]. The resulting values of the generator are random and determined by the mean value and a normal distribution over a range of values extracted from the literature. The generator distributes the osteons randomly following the analysis carried out by Gao et al. [25]. In the generating process of the osteons, the Haversian canals that osteons have inside are also generated. In addition to these parameters, Table 1 shows other parameters that have been calculated from the input parameters and estimates. These parameters allow the models generated to be validated and compared with analyses carried out by other authors. The work conducted by Pérez-Cano et al. [28] details all aspects related to the generation of the microstructures in the selected area. Figure 5 represents a bone section with all of its constituent microstructures.

Table 1 Calculated parameters and formulas used in the study
Fig. 5
figure 5

Representation of a cross section of a fibula with all the bone microstructures that make up the tissue

Fig. 6
figure 6

Grids that manage the distribution of bone microstructures

One of the key points to consider for fracturing is the way in which the information of bone microstructures is structured. Two grids are generated, adjusted to the boundaries of the area in which the bone microstructures are generated (Fig. 6). The first one is a large grid that is responsible for storing the information of all the bone microstructures that are found within each of the quadrants. Moreover, it stores the information concerning all osteons whose centre is within the grid cell. This grid is mainly involved in the distribution of osteons and the prevention of unwanted overlaps. In addition, the characteristics of each osteon such as the area, the position of Haversian canal or the points that delimit the border of the osteons (cement line) are also stored.

The second level grid has a larger number of divisions and is used to distribute the osteons on the cortical area. These grid cells are divided into boundary cells, free cells or occupied cells. Thus, the distribution of the bone microstructures can be managed with precision. The grid is also responsible for refining the representation of the bone structures to avoid that the limits of the cortical area are overstepped. In addition, this grid checks the overlaps to discard those that are undesirable. The information from both grids is used to minimise the computation time in the fracture simulation.

3.2 Fracture Initiation

The fracture starts with the selection of any point that constitutes the outer boundary of the cortical area representation. The simulator also needs a direction in which the fracture will start to propagate as well as the input of parameters relating to the crack conditions that occur (Fig. 7). Before starting propagation, a ray is projected to a position far away that provides the number of collisions with the borders of the cortical area, so that the fracture only begins the propagation if a valid direction is chosen, that is, a direction towards the inner part of the cortical area is selected.

Fig. 7
figure 7

Selection of a fracture direction. The blue dots determine the possible directions that can be selected

Fig. 8
figure 8

Interaction of a fracture with an osteon during its growth from the study conducted by Tang et al. [1]. a Represents a crack deflected by the cement line, b shows a crack deflected by the lamellae that surround the Haversian canal and c a crack through the Haversian canal. d Shows a crack deflected along the cement line

As mentioned above, the objective of this method is to simulate a realistic-looking fracture of the entire cortical area at microscale using a statistical model. The parameters that are taken into consideration to set the fracture initiation are the energy applied to the model as well as the maximum energy it is able to resist before a fracture. These values are used to delimit the length of the fracture line and will be discussed in the following sections. Factors such as modulus of elasticity, toughness and other mechanical properties have not been taken into consideration since some are unknown and others have been previously studied in the literature using the FEM method. However, the complexity and computational cost of these simulations imply that fracturing a complete cortical section requires significant resources and time. In addition, manual segmentation of bone microstructures is necessary before the simulation. Thus, the method we propose focuses on replicating real clinical cases of a complete bone section to conduct a fracture simulation than emulate the physical conditions of the bone tissue. This fracture simulation studies the influence of bone microstructures on the growth of a fracture. Once the values described have been provided, the algorithm starts to propagate the fracture in the indicated direction.

3.3 Propagation Criteria

In terms of fracture propagation, many studies report the influence of osteons and the heterogeneous characteristics of the tissue during fracture growth [2,3,4, 17, 21, 23, 24, 26, 40]. Tang et al. [1] performed a microscopic study about the mechanical properties of the cortical tissue. In this study they applied forces on cortical tissue samples analysing the fractures produced at both macroscopic and microscopic levels. They report that three situations can occur when a fracture line reaches an osteon:

  1. 1.

    Be deflected by the cement line and went around the osteon (Fig. 8a and d).

  2. 2.

    Propagate into the osteon and then being deflected by the lamellar structure without arriving at the Haversian canal (Fig. 8b).

  3. 3.

    Pass through the Haversian canal without being deflected (Fig. 8c).

The method proposed considers these three cases when propagating the fracture in the indicated direction. The fracture does not propagate completely in a straight line, it has a margin of ±10 degrees with respect to the previous direction. This improves the realism of the results obtained in the simulation since the fracture lines present small perturbations in their progress through the cortical tissue. Once fracturing is initiated following an initial direction, the fracture line still propagating until a stop condition occurs. The percentage indicated above has been established after several tests and visual evaluations of the results, although this could be changed.

Fig. 9
figure 9

It represents the different types of cells that a fracture may encounter in its path. The green cells are the empty ones, the yellow ones are the ones that are part of the cement line of the osteons and the grey part represents the interior of an osteon

The data structure employed to store information about cortical tissue consists of two grids. Since osteons act as barriers to fracture growth, and various scenarios may arise involving them, the smaller grid is primarily used to manage the collisions between the fracture line and the osteons. Figure 9 represents main types of cells that make up the grid. The green cells represent the empty cells so the fracture line can be propagated through them without any additional calculation. The yellow areas represent the cement lines of the osteons so it is necessary to determine how the fracture growth continues after the collision with an osteon. In each cycle the fracture advances from a grid cell to an adjacent cell until a stop condition occurs. If this cell is free of osteons, the fracture continues growing. In case the grid cell is occupied or partially occupied, the osteons that are in the area are identified and the collision is being analysed. The cement line of an osteon and the growth of the fracture line are both determined by polylines. Therefore, an intersection between two segments is performed to determine whether a collision with the cement line of the osteon occurs as well as the point in which it occurs. In addition, the angle at which the collision is produced is also analysed. Without any collision with the osteon, the fracture line will continue its path as it implies that it has passed very close to the border of an osteon but there is no collision.

The study carried out by Tang et al. [1] analysed the frequency at which the situations occurred when a fracture line collided with an osteon. In 50% of the collision cases, the fractures propagated into the osteon and were deflected by the concentric layers that surround the Haversian canal, but did not reach into it (Fig 8b). Moreover, only 12% of the cases reported the first situation in which the cement line deflected the fracture line (Fig. 8a and d) The remaining percentage (38%) correspond to the third situation in which the fracture line penetrated the osteons and reached the Haversian canals (Fig. 8c).

These previous studies also show that the first osteons encountered by a fracture line are traversed regardless of the angle at which it collides. That is due to the fact that they are close to the point of impact and the energy applied to the model is sufficient to allow the growth of the fracture line. Additionally, it is observable that the further the line advances, the greater the resistance of osteons to fracture due to the energy loss caused by the absorption of energy from cortical tissue as the fracture line progresses and collides with osteons. The direction in which the fracture advances is the most important factor to determine the fracture path and length since the cement lines of the osteons are the main barrier and reduce the potential for the energy to advance through the bone tissue.

The algorithm developed considers the features shown in the previous paragraphs. The fracture begins at the boundaries of the bone and propagates until a stop condition occurs. These conditions are detailed in the next subsection and are based on the loss of energy as the fracture line grows along the bone tissue, similar to the propagation of the light in biological tissue described by Kumar et al. [6]. The fracture line advances regardless of the osteons it encounters up to 50% of the maximum distance expected for a fracture line, as previously indicated. The microstructural barrier effect is dependent on the crack length at the time of encountering an osteon by the large amount of energy with which it propagates at the beginning of the fracture [41]. In this situation, an osteon may actually act as a weakness in the bone and facilitate fracture propagation. There will be cases in which it will penetrate the Haversian canals and cases in which it will not, according to the direction of the fracture. Once 50% of the maximum length of the fracture line is covered, collisions and the angle at which they occur are considered. If the angle of the collision formed by the fracture line with the polyline representing the cement line is between 60 and 120 degrees, the line penetrates the osteon, giving rise to situations 2 and 3 (Fig. 10a). However, if the angle is above or below this interval, the cement line will act as a barrier and the fracture will advance along its edge (Fig. 10b and c). The percentages and degrees have been adjusted by running several tests and setting the parameters until obtaining a percentage of cases of collisions similar to the analysis conducted by Tang et al. [1]. The evaluation of these parameters has been analysed in Sect. 4. As mentioned above, all these parameters and values can be customised to adapt them to the fracturing simulation.

Fig. 10
figure 10

Interaction of fracture growth upon collision with an osteon. a Represents a collision at an angle close to 90 degrees causing the line to penetrate the osteon (angle between 60 and 120 degrees). b and c represent collisions outside the range between 60 and 120 degrees causing the cement line to act as a barrier and deflect the fracture line

3.4 Stop Criteria

The propagation of a fracture line at the microscope level can end in two ways: fracturing the model or not. A model is considered to be fractured when the fracture propagation line crosses the bone boundaries twice, while a fracture line that only crosses the bone boundaries once is called a crack. It is not necessary that the line crosses the outer and inner limit to refer to a fracture; it could cross the outer limit twice, causing a detachment of a part of the cortical zone without reaching the centre of the cortical area. Likewise, it would not be representative of the simulation of a fracture starting at the inner edge of the bone, so in no case will it cross it twice. The fracture line will advance until the model is fractured or until enough energy is dissipated, causing the fracture to cease propagation after colliding with bone structures.

The main elements involved in determining the growth arrest of a fracture line are the energy applied to the model, the energy required to fracture the entire cortical area, the covered distance by the fracture line and the number of osteons found in its path. Once the user select the initial point of the fracture line and the direction, the distance between the starting point and the nearest collision point with the boundary of the cortical tissue representation is calculated in the chosen direction. This distance will be the equivalent to the distance that the fracture line can be traversed if the maximum energy that the model can resist is applied as input. In addition, the fracture line also considers the osteons that are found in its path once the line has advanced following the indications of the literature [1]. Studies have shown how osteons act as barriers to fracture growth and their ability to absorb part of the energy. The fracture line only considers the encounters with osteons after growth the 50% of the maximum path the fracture can travel with the energy established. Until it reaches this length, the line passes through the osteons regardless of the angle with which it collides. Subsequently, each time the fracture line collides with an osteon, its maximum propagation length is reduced a percentage. This simulates the absorption of energy by collisions with bone microstructures. All these values can be customised to adapt them to the circumstances of the fracture to be simulated. These values have been established after several tests and following the analysis of the microscopic images of bone tissue reviewed in the previous work.

4 Results

Fig. 11
figure 11

a and b represents the first 10 complete fracture simulation in the same model with different osteonal density. a Shows the first model with an osteonal density of 15.78 while b has an osteonal density of 16.04

Table 2 Percentages of each type of collision managed on each cortical tissue model compared with the study conducted by Tang et al. [1]

Three different models have been generated using one fibula. Each generated model represented a different cross section of the bone with different geometry and area. For the generation of the bone microstructures that make up the cortical tissue, an osteonal density of 12.69\(mm^2\), 15.78\(mm^2\) and 16.04\(mm^2\) was established. These density values has been extracted from clinical cases analysed by Yeni et al. [42]. These values corresponds to the bone density of different cadavers of women and men over 50 years of age. With this osteonal density, the osteons represent between 40% and 50% of the bone volume. Fifty fracture lines have been launched with sufficient energy settings to ensure complete bone fracture in each model. Figure 11 only shows the first 10 simulated fracture lines in the first model with varying osteonal density. Each fracture simulation take less than a second to manage several fractures simultaneously to analyse how the bone microstructures affect the fracture propagation. Table 2 shows the way in which the collisions were managed according to the weights established in the previous section. As can be seen, the observed percentages closely align with those previously described by Tang et al. [1].

Figure 13 shows the fracturing of the entire model. Figure 13a and b shows the advance of a fracture line through the model using the same input configuration (starting point, direction and energy). Thus, it is possible to check the effect of the random variable and the encounters with the bone structures in the growth of the fracture line. Moreover, the model used as a starting point is the same in all cases. These figures show the three situations identified in the study conducted by Tang [1]. As in their study, collisions predominate in which penetration occurs in the internal part of the osteon, and as the distance travelled increases, the possibility of these acting as barriers and deflecting the fracture increases. Figure 12 illustrates how the energy and the number of osteons encountered by a fracture in its path determine the fracture length, resulting in a crack or a complete fracture of the bone. This is the result of the energy absorption from the collisions. All simulations use the same input configuration (starting point, direction and energy). Figure 14a and b highlights the different situations that can occur when a collision occurs with an osteon.

Fig. 12
figure 12

This figure shows several fracture simulation in the same model with the same initial configuration. The distance travelled by each fracture line depends on the number of osteons it encounters

Fig. 13
figure 13

a and b represents the fracture of a model. Both figures use the same inputs

Fig. 14
figure 14

a and b shows a more detailed view of the area in which a fracture line collides with an osteon

In addition to check the way in which collisions are handled in previous studies, a second validation process was carried out by consulting three experts in the field. Two of them had expertise in bone fracturing with expertise in the different hierarchical levels of the bones and the other expert has experience in simulation using the FEM method. To evaluate the results, they have been provided with a set of 15 images of bone fractures and asked to evaluate the results obtained through a questionnaire using a 5-point Likert scale with 5 being the highest score.

In the questionnaire, they were asked about the variability of situations shown in fracture cases at the microscopic level and the realism of the simulation. Regarding the assessment of the variability of situations in fracture growth, the three experts provided an average score of 4.66 out of 5. They found that the shape of the fractures had a realistic appearance as they advanced along the surface of the osseous section. However, the average score obtained for the degree of realism of the simulation was 3.33. According to expert feedback, the lack of textures and an inability to discern the inherent roughness of the material, resulting in a perceived decrease in realism. In addition, when they zoomed in too much on the image, they could identify some pronounced changes in the direction of the fracture as it progressed.

5 Discussion

The algorithm developed focuses on cortical bone fracturing at the microscopic level starting from a bone model at macroscopic level in 3D. For this purpose, the first step was to select an area of the bone and proceed to generate a 2D representation of the microstructures that make up the bone tissue in the selected section. The main objective is to obtain a realistic-looking fracture that is determined by the microstructures of the bone tissue in which it is produced. Thus, several parameters and the criteria has been proposed to control the propagation of the fracture line on the generated model. The algorithm generates fractures along the model until it reaches a stop condition. Two scenarios are considered: either a complete fracture or a crack. The fracturing of the model occurs when the fracture line crosses the boundaries of the cortical area twice. Additionally, a crack is a fracture without a complete model breakage. This is due to insufficient energy initially applied. Fracture propagation typically follows a specific direction, but there is some randomness involved, leading to slight deviations in the propagation path. Thus, the resulting fractures are prevented from being linear. In addition, a number of rules that determine how fracture lines behave when they encounter an osteon in their path have been indicated above. In contrast to previous work, the results show how the complete fracture of the entire cortical section can be simulated in a matter of seconds as opposed to FEM-based simulations where small areas of the cortical tissue are fractured and need to be manually segmented. In addition, the experts determined that the synthetic representation generated with the proposed algorithm is faithful to the nature of real fractures and received a high score in the evaluation process. While the visual appearance of the bone tissue can be improved, the algorithm proposed fulfills its intended purpose efficiently following an statistical approach.

FEM method uses all the elements surrounding the fracture in question to determine its growth. In these studies, the representations used are based on small sections of the cortical area in which only a few microstructures are visible. Therefore, the simulations that are carried out are inaccurate because many of the structures that make up the cortical tissue and should affect its growth are ignored. The high computational cost makes it very difficult to carry out a fracture simulation on a model with all the structures that make up the cortical tissue using the traditional method. In a cut of the diaphysis area of a femur there may be more than 1000 osteons while FEM studies usually show no more than 10 structures in their simulations. In addition, the bone model used to perform these simulations need to be segmented by hand. So, these studies require a process to obtain images from a microscope and manual segmentation before performing the physical simulation.

There is no limitation with the configuration of the proposed algorithm, so the user can adjust it to his needs, although a series of values based on previous literature are proposed. The fracture result obtained after the simulation is determined by the configuration established by the user. Moreover, each fracture simulation take less than a second so several fracture lines can be simulated at the same time to analyse and compare it.

The rules established allow the fracture growth to be precisely oriented. As seen above, the results obtained with the default values have been positively evaluated by the experts. The parameters established allow to reproduce the conditions of different clinical. The evaluation of the results demonstrates the variability of the cases in fracture and how the proposed algorithm handles collisions between the fracture line and the osteons in a manner similar to that observed in previous works. In addition, it is possible to replicate the conditions of different clinical cases to carry out the simulation as well as to launch multiple fractures due to the speed of the proposed simulation method.

Considering the results and the evaluation realised by the experts, the results can be considered promising, although the algorithm should to improve by increasing the number of physical parameters taken into consideration during the simulation. In this work, the fracture is not approached from a physical point of view, but as a statistical simulation of the entire bone. FEM simulations are more interesting for studying the physical and mechanical properties of bone microstructures. However, they cannot be carried out on a complete section of the bone and require a greater amount of time due to the high computational cost, the need to obtain a microscope image to segment it, and the manual generation of the model for simulation. In the visual aspect, some improvements could be made to improve the perception of the experts about the method. Textures obtained from microscope images could be used to increase the degree of realism as well as to add some roughness to the borders of the osteons so that the propagation through them, when they act as a barrier, is not linear.

6 Conclusion and Future Work

In this paper, a method for fracturing bone tissue representations at the microscopic level has been developed using as starting point a 3d model of the bone at the macroscale level. To move from the macroscale to the microscale model, a specific area of the bone was selected and a 2D model was generated where all the bone microstructures that make up a section of the bone are represented.

The proposed algorithm allows to fracture a complete cortical bone section considering all the microstructures that compose the cortical tissue, without any limitation in the representation and in a few seconds. The traditional method for this process is unfeasible due to the significant amount of time required to prepare the simulation and its high computational cost as described in the previous section. The algorithm establishes a set of rules that determine the growth of the fracture using an statistical approach as well as allow to replicate the conditions of different clinical cases to carry out the simulation. Several simulations can be launched at the same time can be due to the speed of the proposed method. To mitigate some of the limitations of this model, an attempt could be made to combine the proposed model with the results obtained from a simplified FEM model. Advances in studies on the longitudinal axis of the bone may also allow progress to be made in obtaining a 3D representation of the bone as well as its fracture, obtaining more accurate results.

Finally, the results obtained with this method could be exported and used to generate data sets for later use by the scientific community. In addition, they allow a more accurate assessment and study of the progression of a fracture along the bone tissue. Moreover, deep learning-based models could be trained both to classify the different fracture patterns and to improve the accuracy of the fracture lines obtained. So far, the simulations available focus on small areas of the bone tissue in which there are large differences between the bone microstructures. Moreover, the settings can be adjusted with data obtained from other simulations for the entire cortical section in seconds and without having to manually segment microscope images.