Numerical simulation of stress amplification induced by crack interaction in human

femur bone
Noor Alia, Ruslizam Daud, Mohammad Fadzli Ramli, Wan Zuki Azman, Ahmad Faizal, and Siti Aisyah

Citation: AIP Conference Proceedings 1660, 070008 (2015); doi: 10.1063/1.4915726

View online: http://dx.doi.org/10.1063/1.4915726
View Table of Contents: http://aip.scitation.org/toc/apc/1660/1
Published by the American Institute of Physics
 Numerical Simulation of Stress Amplification Induced by
Crack Interaction in Human Femur Bone
Noor Aliaa, Ruslizam Daudb, Mohammad Fadzli Ramlia,
Wan Zuki Azmana, Ahmad Faizalb and Siti Aisyaha
a
Institute of Engineering Mathematics, Pauh Putra Campus, Universiti Malaysia Perlis,
02600 Arau, Perlis, Malaysia
b
Fracture and Damage Mechanic Research Group, School of Mechatronic Engineering,
Pauh Putra Campus, Universiti Malaysia Perlis, 06200 Arau, Perlis, Malaysia

Abstract. This research is about numerical simulation using a computational method which study on stress amplification
induced by crack interaction in human femur bone. Cracks in human femur bone usually occur because of large load or
stress applied on it. Usually, the fracture takes longer time to heal itself. At present, the crack interaction is still not well
understood due to bone complexity. Thus, brittle fracture behavior of bone may be underestimated and inaccurate. This
study aims to investigate the geometrical effect of double co-planar edge cracks on stress intensity factor ( ) in femur
bone. This research focuses to analyze the amplification effect on the fracture behavior of double co-planar edge cracks,
where numerical model is developed using computational method. The concept of fracture mechanics and finite element
method (FEM) are used to solve the interacting cracks problems using linear elastic fracture mechanics (LEFM) theory.
As a result, this study has shown the identification of the crack interaction limit (CIL) and crack unification limit (CUL)
exist in the human femur bone model developed. In future research, several improvements will be made such as varying
the load, applying thickness on the model and also use different theory or method in calculating the stress intensity factor
( ).
Keywords: Stress amplification; stress intensity factor; crack interaction limit; crack unification limit.
PACS: 62.20.-x, 62.20.M-, 62.20.mt

INTRODUCTION
Bone is the main component in the skeletal system in our body. Bone plays important roles, especially to support
our body. There are 206 bones inside an adult skeleton system (Marieb, 2009). The bone that is responsible to carry
our body weight is located at the lower part of our body which is called “femur” bone. The femur bone is located in
the thigh region. The femur bone is the strongest bone in the human body. If there is a fracture or crack happening
within the bone, there must be large force acting on it. Fracture in a bone can cause impairment to it, where the bone
cannot bear the load of the person’s weight. Fracture happens when there is an initial crack occurring when certain
amount of load and stress acting on the bone repeatedly. The stress at the tip of a sharp crack has the highest stress,
which can lead to failure of the material. This explains why small cracks can contribute to the failure of the entire
structure. There are a few parameters of fracture that are widely used which include the stress intensity factor (K),
Elastic energy release rate ( ), J-integral ( ), and crack tip opening displacement (CTOD)
(Rudraraju, 2004).
In this research, the main focus is to simulate stress amplification induced by crack interaction. It is based on the
stress intensity factor ( ) calculation using the stress singularity. This simulation is based on the changes of stress
amplification in different crack intervals towards the human femur bone. Stress amplification in a solid, with the
existence of crack is the amount of load that exhibits onto the solid perpendicularly. Other than that, crack
unification limit (CUL), crack interaction limit (CIL) and stress interaction limit (SIL) also can be calculated.
There are several biomechanics properties of the bone that need to be further understood. One of them is the
stress amplification interaction within the bone. This property is important because the effect of stress amplification
on the crack may provide accurate estimation of healing processes in bone fracture. Stress amplification can be
evaluated using numerical simulation. Numerical simulation is important because it provides tools that can analyze
the crack problem in a given boundary condition. This project examined the maximum limit of crack interaction
based on several parameters. This is vital because in real life situations, human femur bone will be completely

International Conference on Mathematics, Engineering and Industrial Applications 2014 (ICoMEIA 2014)
AIP Conf. Proc. 1660, 070008-1–070008-8; doi: 10.1063/1.4915726
© 2015 AIP Publishing LLC 978-0-7354-1304-7/$30.00

070008-1
damaged if the parameters exceed certain values such as the fact that the load is exerted too much on the bone or
that the crack length is too large.

STRESS AMPLIFICATION SIMULATION

There are three parameters which are used to characterize the stress which cause a component to fracture. These
include the stress intensity factor ( ), elastic energy release rate ( ) and J-Integral ( ). The stress intensity factor
( ) is the stress field that surrounds the crack tip of an isotropic linear elastic material. The stress can be formulated
into this equation:
(1)

Where is the stress intensity factor and is only valid near the crack tip. The represents the stress field around
the crack. Thus the value of determining the level of stress magnitude around the crack tip which has been done
by previous researcher (Rudraraju, 2004). The value of should not exceed the value of the critical crack intensity
factor ( ) because if the value of is bigger than , the material will fracture
(Huang, Z. ,2011), as expressed in Eq.(2)
. (2)

Every material has its own specific Poisson’s ratio and Young’s Modulus value. Poisson’s ratio is the ratio of
the transverse direction strain which is normal to the acting load with the direction of the acting strain or load.
Young modulus is also known as elastic modulus. It measures the stiffness of a material. If the value is higher, it
means the material is harder and stiffer. Thus, it is hard to stretch the material. The mechanical properties of the
femur bone (polyethylene) are: Young’s modulus = 1.3 GPA, density = 950 kg/m3, Poisson’s ratio = 0.42
(Fedida et. al., 2005). The crack will propagate and start to initiate from the crack tip in Mode I, II and III. In this
project, the crack failure is composed of mixed mode situations which combine both Mode I and Mode II crack
deformation. It can be expressed as:
(3)

where i and j represent the direction of crack propagation (x and y directions) and I and II represent the mode of
crack deformation. Based on the theory of linear elasticity, when the stress acts on the crack tip, the crack introduces
discontinuity in the elastic (Pithioux et. al., 2004). All the calculations of are done near the crack tip. FIGURE 1
illustrates the crack tip coordinate based from Barsoum singularity derivation (Barsoum, 1975).

FIGURE 1. Geometrical definition of multiple edge cracks with Mode I loading,

(b) Eight node quadrilateral element (c) Complete singular finite element at crack tip.

070008-2
 The stress acted on the crack tip can be derived as the interior asymptotic expansion. The derivation of the stress
at the crack tip can be expressed as

(4)

Based on Eq. (4), for the value which intends to propagate towards the 0 value where represents the finite
stress at the crack tip. Hence, the corresponding value of the stress field is indicated by

(5)
(6)

Analytical data are important to validate the system whether the computational method is accurately the same
with the analytical method. This means that the computational method can be used in real world phenomena. Based
on previous study (Mohanty and Verma, 2010), crack propagation is characterized by changes of crack length
values and number of cycle load acting on the crack tip. The analytical equation used to determine the value of the
stress intensity factor ( ) is based on Brown and Srawley (1966)

(7)

where is the load applied to specimen and is the ratio of crack length and width of the specimen.
For this research, single and double crack tip design has been constructed. The geometrical analysis and design
phase were then carried out using ANSYS software. Once the design was verified according to ASME standard for
crack interaction, meshing process with 8 node quadrilateral element was developed. By using Barsoum singular
element and using Displacement Extrapolation Method, stress intensity factor ( )can be calculated (Henshell and
Shaw, 1975; Barsoum, 1976 and Barsoum 1977). Several parameters were implemented into the system. Then, the
analytical results were compared with the computational results for verification and validation. The results were then
documented for further improvements on the design. The value of is determined by computational method using
ANSYS software where assumption of the model is under the LEFM theory. It uses behaviour in strain and
stress in the element of crack tip. For and modes, the stress intensity factor value is solved using linear
displacement extrapolation methods (DEM) (Kuang et. al., 1993). FIGURE 2 shows the implementation of Barsoum
singular element for double edge crack interaction.

1
3 2
5 4
b
b

a1
V-Axis

c1
a1 a2 bne b

c1 c2 H U-Axis s H
a2

c2

1
3 2
5 4
h

m=8

FIGURE 2. Barsoum singular element for double crack interaction (a) (b) .

070008-3
 Based on FIGURE 2, the arrow represents the stress or load acting perpendicularly on the model 100 MPa,
is the crack length, is the crack interval, s is the crack distance, is crack tip 1, is crack tip 2 while the
numbers 1-5 represent s the nodes on the crack tip.

(8)

(9)

where is the Young Modulus, since plain strain are used, is the length of element and are
both the displacement indisplacement in the internal Cartesian coordinate system while is the Poisson’s ration.

RESULTS AND DISCUSSIONS

The results which based on the graph also show that there is a pattern of crack interaction limit (CIL) `and crack
unification limit (CUL) in the model constructed. The numerical results were then compared with the analytical
result using Brown & Srawley (1966) theory for the purpose of validation. The numerical values obtained from
computational method are compared with analytical values which based on previous experimental research. Figure 3
illustrates the differences in stress intensity factor, value for single edge crack, double edge crack and analytical
result. Based on all the graph constructed, the x-axis represents the ratio of while the y-axis represent the shape
correction factor, or where . Figure 3 displays the differences of value based on
several cases. For analytical results, Brown & Srawley (1966)’s theory where the value of stress or load acted on the
model which is 100 MPa. Based on the Figure 3, it can be seen that if the value of the crack length increase, the
value of b decreased. It can be observed that stress intensity factor ( ) for single edge crack is about equal with the
analytical Brown & Srawley. This is important in validation of developing the algorithm for stress amplification in
bone. Thus, the changes in crack length values are considered for calculating the value of . These two equations
below, as used in ANSYS software, seek to calculate the value of and which used five nodes of path
network based on FIGURE 1. Generally it also observed that at lower crack width ratio, the intersection of CUL is
identified while at higher ratio (0.42-0.45), the intersection of CIL is clearly shown.

070008-4
 6
Single edge crack
Double edge crack (crack tip 1)
Double edge crack (crack tip 2)
5 Analytical result (Brown & Srawley)

Shape correction factor, Y

4

CIL region
3

2 CUL region

0
0.02
0.03
0.05
0.07
0.08

0.12
0.13
0.15
0.17
0.18

0.22
0.23
0.25
0.27
0.28

0.31
0.33
0.35
0.36
0.38

0.42
0.43
0.45
0.47
0.48
0.4
0.1

0.2

0.3
Crack to width ratio, a/w

FIGURE 3. Graph of differences in stress intensity, K value for single edge crack, double edge crack (crack tip 1), double edge
crack (crack tip 1), and analytical result.

The graph in FIGURE 4 focuses on the similarities between values in the single edge crack and the analytical
result or theoretical results based on Brown and Srawley’s (1966) theory. This has shown that the numerical results
obtained from the finite element simulation have been verified. Comparison with the analytical data prove that the
data has been verified and error value is acquired by MATLAB software where the error is very small; 0.5665%.

FIGURE 4. Graph of comparison in K values between numerical data and theoretical data.

FIGURE 5 until 7 shows the intersection point of single edge crack data with analytical data and double crack
data for both crack tip 1, c1 and crack tip 2, c2. All the three intersection points in these figures give the value of

070008-5
crack interaction limit (CIL) range of the developed model. Figure 5 represents the intersect point between single
crack and with the analytical data. The intersection point from the results occur at while .

FIGURE 5. Intersection point between single crack and analytical data.

The intersection point in FIGURE 6 shows the value while where the intersection occur
intersection point between single crack and double crack at crack tip 1.

FIGURE 6. Intersection point between single crack and double crack (crack tip 1).

From the intersection point in FIGURE 7, it can be seen that the point intersects and . Hence, the
range of crack interaction limit (CIL) is and .

070008-6
 FIGURE 7. Intersection point between single crack and double crack (crack tip 2).

CONCLUSION
This project has used various types of parameter in investigating the fracture behavior effects on the stress
intensity factor (K) in human femur bone. Biomechanical properties of the bone have been implied in the
computational method which includes Young’s modulus and Poisson’s ratio. The results show that values are
directly proportional to values. However, if the difference value of the crack distance is used, it proves that
value will be oppositely proportional to the values. The higher the value of , the lower the value of obtained.
On the other hand, the results have also shown that at a higher rate of , the double edge crack shows that the
model becomes more unstable compared to the single edge crack. The results obtained also prove that there are
crack unification limit (CUL) also crack interaction limit (CIL) in the model developed. However, the range value of
the CIL can only be acquired because the graph intersects with each other. There is no exact value for the CUL
obtained in this project but it can still be seen that the graph converges with each other, which clearly indicates that
there is a CUL in the model. Lastly, the corresponding computational results of crack propagation have been
obtained and verified using theoretical results using Brown & Srawley’s theory.

ACKNOWLEDGMENTS
First of and foremost, praise to Allah who had given me the strength, courage and good health to complete this
research. There are several people were involved in the completion of this study, those who have provided me with
invaluable information’s, experimental feedback, views, comments and support. Without them, I would never have
been able to finish this project.

REFERENCES
1. A. Carpinteri, R. Brighenti and S. Vantadori, Engineering Fracture Mechanics 71, 485-499 (2004).
2. B. Gross, J. E. Srawley and W. F. Brown Jr, “Stress-intensity Factors for a Single-edge-notch Tension Specimen by
Boundary Collocation of a Stress Function,” (1964).
3. Huang, B. W., Huang, M. Y., Tseng, J. G., Chang, C. H., Wang, F. S., Lin, A. D. and Tsai, Y. C., Dynamic Characteristics of
a Hollow Femur 9, 723-726 (2012).

070008-7
4. J. H. Kuang and L. S. Chen, Engineering Fracture Mechanics 46, 735-741 (1993).
5. J. J. Kruzic, J. A. Scott, R. K. Nalla, R. O. Ritchie, Abdullah and Shahrum, Journal of Biomechanics 39, 968–972 (2006).
6. M. Doblare, J. M. Garcia and M. J. Gomez, Engineering Fracture Mechanics 71, 1809–1840 (2004).
7. M. R. Ayatollahi and R. Hashemi, Composite Structures 78, 602-609 (2005).
8. M. S. Hasan, in Biomedical Engineering Conference, Proceedings of the 1996 Fifteenth Southern IEEE Xplore, 494 – 496
(1996).
9. M. B. Shaha, J. L. Ferracaneb and J. J. Kruzica, Dental Materials 25, 760-770 (2009).
10. M. Bessho, I. Ohnishi, T. Matsumoto, J. M. S. Ohashi, K. Tobita, M. Kaneko and K. Nakamura, Bone 45, 226-231 (2009).
11. P. Hutar, L. Náhlík and Z. Knésl, Procedia Engineering 45, 653–657 (2007).
12. R. Daud, A. K. Ariffin, S. Abdullah and A. E. Ismail, Meccanica 47, 1141-1156 (2012).
13. R. Fedida, Z. Yosibash, C. Milgrom and L. Joskowicz, Proceedings of ICCB05-II International Conference on
Computational Bioengineering 10, 85-96 (2005).
14. R. Hambli and E. Lespessailles, Journal of the Mechanical Behavior of Biomedical Materials 17, 89-106 (2013).
15. R. S. Barsoum, International Journal of Fracture 10, 603-605 (1975).
16. S. Benbarek, B. A. B. Bouiadjra, B. M. E. Mokhtar, T. Achour and B. Serier, Design and Computation of Modern
Engineering Materials 54, 73-87 (2014).
17. T. Hao and T. Yifei, Procedia Environmental Sciences 10, 640 – 646 (2011).
18. X. F. Hu and W. A. Yao, Mechanics Research Communications 48, 247–256 (2012).
19. Huang, Z., Himes, J. H. and McGovem, P. G., American Journal of Epidemiology 144, 124-134 (1996).

070008-8