Cortical Bone Model with a Microcrack under Tensile Loading

: The fracture mechanics of cortical bone has received much attention in biomedical engineering. It is a fundamental question how the material constants and the geometric parameters of the cortical bone affect the fracture behavior of the cortical bone. In this work, the plane problem for cortical bone with a mi‐ crocrack located in the interstitial tissue under tensile loading was considered. Using the solution for the continuously distributed edge dislocations as Green 􀆳 s functions, the problem was formu‐ lated as singular integral equations with Cauchy kernels. The nu‐ merical results suggest that a soft osteon promotes microcrack propagation, while a stiff osteon repels it, but the interaction effect between the microcrack and the osteon is limited near the osteon. This study not only sheds light on the fracture mechanics behavior of cortical bone but also offers inspiration for the design of bioin‐ spired materials in biomedical engineering.


Introduction
Through million years of evolution, biological materials such as bones and teeth have evolved with outstanding performance and efficiency [1,2] .Bone can achieve diverse mechanical, biological, chemical and biomedical functions through its natural hierarchical structures [2][3][4] .Bone consists of cortical and trabecular bone.Cortical bone may contain defects, including inclusions, holes and cracks, due to fatigue and cyclic loading [1,5] .Understanding the fracture behavior of the cortical bone is necessary to prevent the fracture of cortical bone and can offer inspiration for the design of bioinspired materials in biomedical engineering in the near future.
Much effort has been devoted to studying the fracture mechanics of cortical bone.It has been observed experimentally that the fracture properties of cortical bone are substantially anisotropic [6,7] , and the fracture toughening mechanism of cortical bone has been evaluated [8] .Many researchers have considered different testing methods and test geometries to understand the fracture behavior of cortical bone [9] .Numerical simulations have been developed to understand the fracture characteristics of cortical bone.For example, the finite element method is used for microfractures in cortical bone [10,11] , and the results indicate that the material properties and morphological parameters of the microstructure greatly influence the fracture behavior of bone.By using the extended finite element method, anisotropic damage initiation criteria are developed to predict the crack initiation and propagation in cortical bone [12] .The fracture characteristics of cortical bone at the micro level are studied by extended isogeometric analysis, and the effect of uncertainties in the osteon Youngs modulus, cement line thickness and porosity percentage on the stress intensity factor values and crack path trajectories is considered [13] .Maghami et al [14] investigated the role of material heterogeneity on crack growth trajectory in human cortical bone using a phase field fracture framework, and the material heterogeneity of microstructural and crack-microstructure interactions play important roles in bone fragility.
Theoretical models have been developed to explore the interaction effect between osteons and microcracks.The microstructure and mechanical behavior of cortical bone share several similarities with fiber-ceramic matrix composites: osteon is analogous to fiber, and interstitial tissue is analogous to matrix [10,15] .The hierarchical structure of cortical bone can be compared to a fiberreinforced composite material [10,15] .Based on the fiberceramic matrix composite model of cortical bone, the integral equation method [16,17] and complex variable method [18] are used in the theoretical modeling of the fracture mechanics of cortical bone, and the interaction effect between microcracks and an osteon has been studied.The result shows that a softer osteon promotes microcrack toward the osteon, while a stiffer osteon repels the microcrack from the osteon, but this interaction effect is limited to the vicinity of the osteon [16][17][18] .However, the specific form of the kernel of the integral equation satisfying the problem has not yet been given, which is very important in the theory of integral equations [19,20] .
From the literature, we learned that the fracture mechanics of cortical bone in the transverse direction are greater than those in the longitudinal direction [8] .The bone is the weakest in tension and the strongest in compression and fails mainly due to tensile forces [21] .Failure occurs not by the accumulation and interaction of many small microcracks but by one microcrack which grows large enough to cause failure [22] .Thus, the objective of this study is to analyze the interaction effect between a microcrack and an osteon under tensile loading and to study the effect of the material constants and geometric parameters of the cortical bone on the fracture mechanics behavior of the cortical bone.
In this paper, by using the Kolosov-Muskhelishvili complex potential [23][24][25] , the elastic fields of the problem of an osteon embedded in the interstitial tissue without microcrack under tensile loading are solved.Using the solution for the continuously distributed edge dislocations as Green  s functions, the problem for the cortical bone with a microcrack located in the interstitial tissue is formulated as singular integral equations with Cauchy kernels.Material constants and geometric parameters, including the shear modulus ratio between the osteon and the interstitial tissue, the radius of the osteon, the distance of the left tip of the microcrack from the cement line and the length of the microcrack, are involved in the numerical analysis of the stress intensity factors of the microcrack tip.The numerical results suggest that a soft osteon promotes microcrack propagation, while a stiff osteon repels it near the osteon.These results correspond well with the results given in the literature [11,[15][16][17] .Additional predicted numerical results need to be confirmed with experiments.

Description of the Model
The model of cortical bone consists of a solid circular osteon embedded in the interstitial tissue, as shown in Fig. 1.The shear modulus of the interstitial tissue and the osteon are G 1 and G 2 , respectively.K 1 = 3 -4v 1 and K 2 = 3 -4ν 2 are material constants, where v 1 and v 2 are Poissons ratios.The radius of the osteon is R. A microcrack L lies along the horizontal axis of the interstitial tissue and is at a distance d from the cement line.The model is under a uniform remote uniaxial tensile loading σ 0 in the y direction.

Solution of the Model
The solution of the model may be obtained through the superposition of two basic problems.In the first problem, an osteon is embedded in the interstitial tissue without the microcrack under uniform remote uniaxial tensile loading σ 0 .In the second problem, only the stress disturbance due to the existence of the microcrack in the interstitial tissue is considered.The only external loads are the microcrack surface tractions, which are equal in magnitude and opposite in sign to the stress obtained in the first problem along the presumed location of the microcrack.

Complex Variable Formulation of the Stress Fields of the First Problem
Let z = x + iy = re iθ , where r is measured from the center of the osteon and θ is measured from the direction of loading σ 0 .
The basic equations are as follows [23][24][25] : where Φ (z) and Ψ (z) are two holomorphic functions that need to be solved.
The boundary conditions may be expressed as By using the Kolosov-Muskhelishvili complex potential [23][24][25] , the two holomorphic functions are obtained: where Then, the stress fields may be formulated as

The Integral Equations for Solving the Second Problem
Let the interstitial tissue contain two edge dislocations at point (x = c 0 y = 0) with Burgers vectors b x and b y .Using the solution for the edge dislocations as Greens functions [26,27] , the stresses at a point P ( xy ) in the inter- stitial tissue may be expressed as π ( ) (11)   where h xx1 , h xx2 , h yy1 , h yy2 , h xy1 , and h xy2 are the horizontal and vertical components of σ xx , σ yy and σ xy , respectively [26,27] .
If we assume that b x and b y are continuously distributed functions of t 0 on L, the second problem may be formulated as (13) where A + B ( ) ( ) where l is the half length of the microcrack, x c is the midpoint of the microcrack, and k 1 ( tt 0 ) and k 2 ( tt 0 ) are bounded.
σ n (t) and σ t (t) are the normal and tangential stresses on L in the first problem.
The continuity of displacement requires that After separating the singularity of the dislocation density, b y ( t 0 ) and b x ( t 0 ) can be expressed as b y ( t 0 ) = F n ( ) where F n ( t 0 ) and F t ( t 0 ) are two nonsingular functions on Substituting equations ( 16) and ( 17) into equations ( 12)-( 15), then ∫ -1 Integral equations ( 18)- (21) give the solution to the problem of a microcrack located along the horizontal axis of the interstitial tissue under tensile loading.

The Stress Intensity Factors of the Microcrack Tips
The stress intensity factors of the microcrack tips are given as follows [28] .The microcrack opening displacement near the tip of the microcrack is dg n ( ) At the right hand microcrack tip, dg n ( ) The dislocation density b y (t) can be expressed as b y (t) = F n ( ) Hence, at the end of the microcrack t = +1, the microcrack tip stress intensity factor is related to lim r ® 0 r dg n ( ) From equations ( 23) and ( 26), we have A similar argument may be used to obtain the stress intensity factor of the other end of the microcrack, Similarly, we have

Results and Discussion
Stress intensity factors are presented for a variety of microcrack-osteon geometries.For the soft osteon, the shear modulus are G 1 = 8.08 GPa and G 2 = 7.31 GPa.For the stiff one, the shear modulus are G 1 = 6.15GPa and G 2 = 7.31 GPa.The Poisson ratio of both osteon and interstitial tissue is v 1 = v 2 = 0.4.A uniform remote uniaxial tensile loading σ 0 = 10 MPa.By using the Gauss-Chebyshev quadrature method [20,28,29] , the integral equations ( 18)-( 21) are numerically solved.
Here, microcrack tips (a) and (b) correspond to the left tip of the microcrack and the right tip of the microcrack, respectively.The normalized stress intensity factor of microcrack tips (a) and (b) versus the shear modulus ratio m = G 2 G 1 near the osteon is shown in Fig. 2.
We assume that the distance d = 10 μm, the fixed microcrack length l = 100 μm and the radius of the osteon R = 100 μm. Figure 2 indicates that for the case of m = 1, the normalized stress intensity factor of microcrack tips (a) and (b) is equal to 1, which means that there is no interaction effect between the microcrack and the osteon.In the case of 0 < m < 1, the enhancement effect of the osteon on the microcrack occurs; in the case of m > 1, the shielding effect of the osteon on the microcrack occurs in front of the microcrack tips (a) and (b).These results are in accordance with the results that osteons act as barriers to microcrack propagation [30][31][32][33] .Moreover, when m changes from 1 to 0, the effective normalized stress intensity factor increases, which means that the enhancement effect becomes increasingly significant; when m changes from 1 to 5, the effective normalized stress intensity factor decreases, which means that the shielding effect becomes increasingly significant.For a fixed value of m, the effective normalized stress intensity factor of the microcrack tip (a) is always larger than that of the tip (b) because the microcrack tip (a) is nearer to the osteon than the tip (b).Furthermore, as clearly shown in Fig. 2, when m changes from 0 to 5 continuously, the effective normalized stress intensity factor of microcrack tips (a) and (b) decreases continuously.The validity of the results still needs to be confirmed.
The normalized stress intensity factor of microcrack tips (a) and (b) versus the normalized osteon radius R l is shown in Fig. 3.We assume that the distance d = 10 μm and the fixed microcrack length l = 100 μm.For the soft osteon (m = 7.31 8.08), the stress intensity factor of microcrack tips (a) and (b) increases as the radius of the osteon increases.For a fixed value of R l, the normalized stress intensity factor value of the microcrack tip (a) is larger than that of the tip (b).For the stiff osteon (m = 7.31 6.15), the opposite tendency is found for the same reason.Figure 3 indicates that when the radius of the osteon is very small, the interaction effect between the microcrack and the osteon can be ignored.When the radius of the osteon increases, the interaction effect between the microcrack and the osteon also increases in front of the microcrack tips near the osteon.This is expected, as the interaction effect between the microcrack and the osteon will be more significant when the radius of the osteon is larger.
We assume that the radius of the osteon is R = 100 μm, and the fixed microcrack length is l = 100 μm.The normalized stress intensity factor of microcrack tips (a) and (b) versus the normalized distance d R is shown in Fig. 4. For the soft osteon (m = 7.31 8.08), the stress intensity factor of the microcrack tip (a) increases sharply as the microcrack approaches the osteon, and the stress intensity factor of the microcrack tip (b) increases gradually as the microcrack approaches the osteon.For the homogeneous osteon (m = 1), the stress intensity factor of the microcrack tips (a) and (b) is equal to 1, which indicates that the interaction effect between the microcrack and the osteon can be ignored.For the stiff osteon (m = 7.31 6.15), the stress intensity factor of the microcrack tip (a) decreases sharply as the microcrack approaches the osteon, and the stress intensity factor of the microcrack tip (b) decreases gradually as the microcrack approaches the osteon.The numerical results indicate that the interaction effect between the osteon and microcrack is limited to the neighborhood of the osteon.The above results are in accordance with the literature [15,16] .
Since the normalized stress intensity factor of the microcrack tip (b) approaches the constant 1, which can be seen in Fig. 4, we should focus on the normalized stress intensity factor of the microcrack tip (a).The normalized stress intensity factor of the microcrack tip (a) is plotted against the dimensionless distance d R for various shear modulus ratios m = 0.50.70.80.91.11.2 1.31.52.0, as shown in Fig. 5.As the shear modulus mismatch increases, the interaction effect between the microcrack tip (a) and osteon will be enhanced.However, the interaction effect is limited near the osteon.The numerical results indicated in Fig. 5 are consistent with the results given in the Ref. [15].
We assume that the radius of the osteon is R = 100 μm.With a fixed microcrack tip (b), the normalized stress intensity factor of microcrack tips (a) and (b) versus the normalized microcrack length l R is shown in Fig. 6.For the soft osteon (m = 7.31 8.08), the stress intensity factor of the microcrack tips (a) and (b) increases with the microcrack length, which is consistent with the results obtained in Refs.[33-35].Initially, the stress intensity factors of the microcrack tips (a) and (b) are almost the same.When the microcrack tip approaches the osteon, there is a sharp enhancement in the stress intensity factor of the microcrack tip (a), which indicates that when the microcrack approaches the osteon, the soft osteon promotes microcrack propagation.For the homogeneous osteon (m = 1), the stress intensity factor of the microcrack tips (a) and (b) increases with the microcrack length, and the stress intensity factor of the microcrack tips (a) is equal to that of the microcrack tip (b), which indicates that the homogeneous osteon has no effect on the stress intensity factor of the microcrack tips (a) and (b).For the stiff osteon (m = 7.31 6.15), the normalized stress intensity factor of the microcrack tips (a) and (b) increases with the microcrack length.When the microcrack tip approaches the osteon, there is a sharp shielding in the stress intensity factor of the microcrack tip (a).This indicates that when the microcrack approaches the osteon, the stiff osteon repels microcrack propagation, which is in accordance with the results that the microcrack is likely to be deflected and often stops growing when encountering the cement line [33] , and that the cement line acts as a barrier for crack growth [36] .The numerical results indicate that a soft osteon promotes microcrack propagation, while a stiff osteon repels microcrack propagation, which is in agreement with the literature [15,37] .
The normalized stress intensity factor of the microcrack tip (a) is more influenced than that of the tip (b), which can be seen in Fig. 6.Thus, we should focus on the normalized stress intensity factor of the microcrack tip (a).The normalized stress intensity factor of the microcrack tip (a) is plotted against the normalized microcrack length l R for various shear modulus ratios m = 0.50.70.80.91.11.21.31.52.0, as shown in Fig. 7.When the microcrack tip (a) approaches the osteon, the more shear modulus ratio difference will appear, and the more effect of the osteon on the stress intensity factor of the microcrack tip (a) will have.The result indicates that the soft osteon represents the enhancement effect of the stress intensity factor, while the stiff osteon represents the shielding effect, which is in accordance with the conclusions given in the Ref. [15].

Conclusion
A cortical bone model with a microcrack located along the horizontal axis of the interstitial tissue under tensile loading was developed, and singular integral equations with Cauchy kernels were obtained in this work.The numerical result indicates that the stress intensity factor of the microcrack is dominated by the material constants and geometric parameters of the cortical bone.On the one hand, some of the numerical results are in accordance with the corresponding results given in the studies [15][16][17]33] ; on the other hand, additional predicted numerical results have been obtained, which need to be confirmed via future experiments or computer simulations. Thenumerical results suggest that a soft osteon promotes microcrack propagation, while a stiff osteon repels it; perhaps the explanation for the result is that the microstructure of the cortical bone and the mineralization of the osteon determine the fracture mechanics of cortical bone.However, this effect is limited near the osteon.Further study may consider the fracture mechanics of periodic osteons and microcracks, which share similarities with cortical bone fracture mechanics.

Fig. 1
Fig. 1 Cortical bone model with a microcrack under tensile loading

Fig. 2 Fig. 3 Fig. 4
Fig. 2 The normalized stress intensity factor of microcrack tips (a) and (b) versus the shear modulus ratio m = G 2 G 1 near the osteon

Fig. 5 Fig. 6
Fig. 5 The normalized stress intensity factor of the microcrack tip (a) is plotted against the normalized distance d R for various shear modulus ratios

Fig. 7
Fig. 7 With a fixed microcrack tip (b), the normalized stress intensity factor of the microcrack tip (a) is plotted against the normalized microcrack length l R for various shear modulus ratios