TOUGHENING AND FAILURE MECHANISM OF COMPOSITE MATERIALS FOR IMPACT
AND DAMAGE TOLERANCE
By
Wei Zhang

A DISSERTATION
Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of
Mechanical Engineering-Doctor of Philosophy
2015

ABSTRACT
TOUGHENING AND FAILURE MECHANISM OF COMPOSITE MATERIALS FOR IMPACT
AND DAMAGE TOLERANCE
By
Wei Zhang
Impact resistance is a material’s ability to resist the force of a sudden impact. Compared to static
loading, materials are more susceptible to a shock force applied over a short duration. Avoiding the
catastrophic failure during impact events is the major challenge in the design of impact resistant
material. Various approaches have been applied to increase the material’s impact resistance. One
approach is to combine two or more materials with distinguished physical or chemical properties
to produce a new material with better performance without changing the characteristics of the constituents. The proposed study focuses on this approach, specifically, characterizing the toughening
and failure mechanism of fiber-reinforced materials or laminate panels for impact and damage tolerance. Two types of materials are chosen for the research work. One is the transparent laminated
glassy panel and the other is the natural fiber reinforced composite material, cortical bone. The material properties and failure mechanism were characterized through uniaxial tensile testing (both
quasi-static and high strain rate), Charpy impact testing, and drop-weight impact testing. Postfailure examinations were carried out using a scanning electron microscope and bright field microscope. Numerical approaches including traditional finite element analysis method, multiscale
modeling technique based on the classical computational homogenization scheme, and XFEM
technique for discontinuities were applied to investigate the correlation between impact resistance
and damage tolerance of the materials to their microstructure, hierarchical structures, and physical
and mechanical properties of each constituent.

ACKNOWLEDGMENTS

The author would like to express her sincere appreciation to her advisor Dr. Srinivasan Arjun
Tekalur for his persistent support, guidance, and advice throughout the duration of this research
work. The author would also like to express her deepest appreciation to Dr. Alejandro Diaz, Dr.
Seungik Baek, and Dr. Laura McCabe for all their valuable suggestions, help, and encouragement
in this research.
The author would also like to thank Andy VanderKlok, Abhishek Dutta, Aiswarya Venkadachalam, David Gonzalez, and Oishik Sen, for providing an amicable environment in the laboratory.
The author would like to thank the support of Gail Berry, Mike McClean, Adam, and Todd for their
help throughout the course of the research. The author would like to acknowledge Ziwei Zhong,
Brinn Cochrane, Brooke Peruski, and Steven Utz for helping run the experiments and machine
the specimens. The authors also would like to thank Bellingar Packing (Ashley, Michigan) for
providing the bovine femurs.
Last the author would like to express her gratitude to her parents, sister, and friends (Miao
Wang, Shutian Yan, Lisi Pei, Danghe Shi, Wu Zhou, Fang Hou, and Christopher Cater) for their
support, and thanks as well to everyone else who helped complete this research.

iii

TABLE OF CONTENTS

LIST OF TABLES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii
LIST OF FIGURES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viii
Chapter 1
Introduction . . . . . . . . . . . . . . . . . . . .
1.1 Objective and Expected Significance . . . . . . . . . . .
1.1.1 Specific Aim 1 . . . . . . . . . . . . . . . . . .
1.1.2 Specific Aim 2 . . . . . . . . . . . . . . . . . .
1.2 Background Rational and Review of Relevant Literature .
1.2.1 Transparent Laminated Panel . . . . . . . . . . .
1.2.2 Damage Accumulation in Bone . . . . . . . . . .
REFERENCES . . . . . . . . . . . . . .
. . . . .

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

. 1
. 1
. 1
. 3
. 5
. 5
. 9
. 13

Chapter 2
Impact Behavior and Dynamic Failure of PMMA and PC Plates
2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2 Experimental Procedures . . . . . . . . . . . . . . . . . . . . . . . . . .
2.3 Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.3.1 Monolithic PMMA Plates . . . . . . . . . . . . . . . . . . . . . .
2.3.2 Layered Plates . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.3.2.1 Bi-layered Structures of PMMA on the Top of PC . . .
2.3.2.2 Bi-layered Structures of PC on the Top of PMMA . . .
2.3.2.3 Bi-layered Structures of PMMA on the Top of PMMA .
2.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
REFERENCES . . . . . . . . . . . . . .
. . . . . . . . . . . . . .

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.

18
18
19
19
19
24
25
27
27
28
33
35

Failure and Toughening Mechanism of PMMA and PC Laminated Glassy
Panel Subjected to Impact Loading . . . . . . . . . . . . . . . . . . . . .
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.3.1 Verification of Material Model for Impact Testing Simulation . . . . . . . .
3.3.2 FEM of Low Velocity Impact Test . . . . . . . . . . . . . . . . . . . . . .
3.3.3 FEM of High Velocity Impact Test . . . . . . . . . . . . . . . . . . . . . .
3.3.4 Optimal Design of Laminate Panel Subjected to Impact Test . . . . . . . .
Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

37
37
41
45
45
46
48
50
52
54

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

Chapter 3
3.1
3.2
3.3

3.4
3.5

iv

REFERENCES . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . 56

Chapter 4

The Effects of Damage Accumulation on Tensile Strength and Toughness of Compact Bovine Bone subjected to Static Tensile Loading, Charpy
Impact Loading, and SHPB High Strain Rate Tensile Loading . . . . . .
4.1 The Effects of Damage Accumulation on Tensile Strength and the Toughness of
Compact Bovine Bone subjected to Static Tensile Loading . . . . . . . . . . . . .
4.1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.1.2 Material and Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.1.3 Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.1.3.1 Imaging Longitudinal Specimens Loaded in the Same Strain Rate
of 0.0001 s−1 . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.1.3.2 Imaging Strain Rate Effects . . . . . . . . . . . . . . . . . . . .
4.1.3.3 Imaging Orientation Effects at a Loading Rate of 0.0001 s−1 . . .
4.1.3.4 SEM Post-failure Examination . . . . . . . . . . . . . . . . . . .
4.1.4 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2 Microcracking Morphology of Compact Bovine Bone Subjected to Impact Loading
4.2.1 Non-notched Beam Specimens of Compact Bovine Femur Subjected to
Charpy Impact Loading . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2.1.1 Material and Methods . . . . . . . . . . . . . . . . . . . . . . .
4.2.1.2 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . .
4.2.2 High Strain Rate Tensile Tests of Compact Bovine Bone by Split Hopkinson Bar(SHPB) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2.2.1 Material and Methods . . . . . . . . . . . . . . . . . . . . . . .
4.2.2.2 Results and Discussion . . . . . . . . . . . . . . . . . . . . . . .
4.2.3 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
REFERENCES . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . .

61
61
64
67
67
69
70
71
72
79
79
79
83
87
87
89
94
97

Chapter 5
5.1
5.2

5.3

5.4

Investigating the Fracture and Toughening Mechanism of Cortical Bone
Using XFEM and Multiscale Modeling Methods . . . . . . . . . . . . . .
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5.2.1 Multi-Scale Modeling of Uniaxial Tensile Testing . . . . . . . . . . . . . .
5.2.2 Characterizing the Fracture Toughness of Cortical Bone with Various Microstructures by XFEM . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5.2.3 Semi-Concurrent Multi-scale Modeling Approach for Modeling Impact
Testing on Cortical Bone by ABAQUS/EXPLICIT . . . . . . . . . . . . .
Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5.3.1 Multiscale Modeling of Bovine Bone under Uniaxial Tensile Loading . . .
5.3.2 Characterizing the Fracture Toughness of Cortical Bone with Various Microstructures by XFEM . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5.3.3 Semi-Concurrent Multi-scale Modeling of Impact Testing on Cortical Bone
Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

61

v

102
102
104
104
109
112
117
117
118
122
123

REFERENCES . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . . 129

Chapter 6
Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 133
6.1 New and Original Contributions Made . . . . . . . . . . . . . . . . . . . . . . . . 133
6.2 Future Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 135

vi

LIST OF TABLES

Table 2.1

Mechanical and physical properties of PMMA and PC plate. . . . . . . . 19

Table 2.2

Threshold energy of mono-PMMA plates under constant impact mass. . . 23

Table 2.3

Threshold energy of mono-PMMA plates under constant impact energy
12 J. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24

Table 2.4

Summary of impact scenarios of bi-layered structures. . . . . . . . . . . . 25

Table 2.5

Threshold energy, energy dissipated, and time to crack of bi-layered PMMA
and PC plates with weak bonding under constant impact mass 6.61 kg. . . 29

Table 2.6

Peak force and time to crack of bi-layered PMMA and PC plates with
strong bonding under constant impact mass 6.61 kg. . . . . . . . . . . . . 29

Table 3.1

Material properties for PC and PMMA. . . . . . . . . . . . . . . . . . . . 42

Table 3.2

Various configurations of the laminate plate panel simulated. . . . . . . . 43

Table 4.1

Comparison of mechanical tensile properties among different investigators. 74

Table 4.2

Summaries of all the testing results. . . . . . . . . . . . . . . . . . . . . 83

Table 4.3

Comparison of mechanical properties of cortical bone subjected to high
strain rate loading among different investigators. . . . . . . . . . . . . . . 92

Table 5.1

Cortical bone with various microstructures simulated in the uniaxial tensile testing. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108

Table 5.2

Elastic properties of osteon and matrix. . . . . . . . . . . . . . . . . . . . 109

Table 5.3

The geometric parameters and elastic properties of RVE for all the cases
studied. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117

Table 5.4

Summaries of calculated Young’s moduli and their corresponding microstructures. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119

Table 5.5

Summary of fracture toughness for all the cases studied. . . . . . . . . . . 120

vii

LIST OF FIGURES

Figure 2.1

Instron Dynatup 9250 drop weight tower: (a) before clamping; (b) after
clamping. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

Figure 2.2

Post-failure pictures of mono-PMMA plates under constant impact mass:
velocities are 0.7 m/s, 1.0 m/s, 2.0 m/s, and 3.0 m/s respectively from left
to right. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

Figure 2.3

Reaction force versus time of mono-PMMA plates under constant impact
mass of 6.61 kg. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21

Figure 2.4

Reaction force versus displacement of mono-PMMA plates under constant impact mass. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22

Figure 2.5

Peak force and time to crack of mono-PMMA plates under constant impact mass. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23

Figure 2.6

Peak force and time to crack of mono-PMMA plates under constant impact velocities (0.7 m/s and 1.0 m/s). . . . . . . . . . . . . . . . . . . . . 23

Figure 2.7

Peak force and time to crack of mono-PMMA plates under constant impact energy. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24

Figure 2.8

Post-failure pictures of bi-layered structures with PMMA on the top of
PC plates under constant impact mass: the velocities are 0.7 m/s, 1.0 m/s,
2.0 m/s, 3.0 m/s, and 5 m/s respectively from left to right; the first row is
the top PMMA plates and the second row is the bottom PC plates. . . . . 26

Figure 2.9

a)Reaction force history of bi-layered structures with PMMA on the top
of PC with weak bonding under constant impact mass 6.61 kg;b)Impact
energy history of bi-layered structures with PMMA on the top of PC with
weak bonding under constant impact mass 6.61 kg . . . . . . . . . . . . . 26

Figure 2.10

a)Reaction force history of bi-layered structures with PMMA on the top
of PC with weak bonding under constant impact mass 6.61 kg; b)Impact
Energy history of bi-layered structures with PMMA on the top of PC with
weak bonding under constant impact mass 6.61 kg; . . . . . . . . . . . . 27

viii

Figure 2.11

Post-failure pictures of bi-layered structures with PC on the top of PMMA
plates under constant impact mass: the velocities are 0.7 m/s, 1.0 m/s, 2.0
m/s, 3.0 m/s, and 5 m/s respectively from left to right; the first row is the
top PC plates and the second row is the bottom PMMA plates. . . . . . . 28

Figure 2.12

Pos-failure pictures of bi-layered PMMA plates under constant impact
mass: the velocities are 0.7 m/s, 1.0 m/s, 2.0 m/s and 3.0 m/s respectively
from left to right; the first row is the top PMMA plates and the second
row is the bottom PMMA plates. . . . . . . . . . . . . . . . . . . . . . . 28

Figure 2.13

Comparison of the reaction force of the mono-PMMA plates to the bilayered PMMA and PC structures with weak bonding and strong bonding
under constant impact mass. . . . . . . . . . . . . . . . . . . . . . . . . 30

Figure 2.14

Comparison of the reaction force of the mono-PMMA plates to the bilayered PMMA/PMMA structures with weak bonding and strong bonding under constant impact mass: (a) Time to crack; (b) Peak force; (c)
Threshold energy. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31

Figure 2.15

Post-failure pictures of the bonded bi-layered structures under constant
impact mass: velocities are 1.0 m/s, 2.0 m/s and 3.0 m/s respectively from
left to right; the first row is the bonded PMMA/PMMA plates, the second row is the bonded PC/PMMA plates and the third row is the bonded
PMMA/PC plates. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32

Figure 3.1

Schematic representation of amorphous polymers under compression. . . 41

Figure 3.2

Geometry and mesh of FE simulations. . . . . . . . . . . . . . . . . . . . 43

Figure 3.3

Schematic representation of the optimization process . . . . . . . . . . . 44

Figure 3.4

Comparisons of the numerically predicted uniaxial stress-strain curves
(compression) of PC to the experimental measurements. . . . . . . . . . . 45

Figure 3.5

Comparisons of FEM low velocity impact simulation results to experimental results: a) Reaction force vs. time for monolithic PMMA subjected to impact velocity of 2.0 m/s; b) Impact energy vs. time for monolithic PMMA subjected to impact velocity of 2.0 m/s; c) Reaction force
vs. time for monolithic PC subjected to impact velocity of 2.5 m/s; d)
Impact energy vs. time for monolithic PC subjected to impact velocity of
2.5 m/s. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46

ix

Figure 3.6

FEM results of layered PMMA/PC laminate panel and monolithic PC
subjected to impact velocity of 5 m/s by a spherical projectile with radius of 6.35 mm: a) Kinetic energy vs. time; b) Distribution of various
energies dissipated calculated by ABAQUS/EXPLICIT. Note: ALLSE,
strain energy; ALLFD, fiction dissipation; ALLPD, plastic dissipation;
ALLDMD, damage dissipation. . . . . . . . . . . . . . . . . . . . . . . . 47

Figure 3.7

Stress contours of 8-layered PMMA/PC laminate, monolithic PC, and
monolithic PMMA subjected to impact velocity of 5 m/s by a spherical
projectile with radius of 6.35 mm at different time frames. . . . . . . . . . 48

Figure 3.8

FEM results of layered PMMA/PC laminate panel and monolithic PC
subjected to impact velocity of 300 m/s by a blunt projectile with radius of
2.73 mm : a) Kinetic energy vs. time; b) Distribution of various energies
dissipated calculated by ABAQUS/Explicit. . . . . . . . . . . . . . . . . 49

Figure 3.9

Stress contours of 8-layered PMMA/PC laminate, monolithic PC, and
monolithic PMMA subjected to an impact velocity of 300 m/s by a blunt
projectile with radius of 2.73 at different time frames. . . . . . . . . . . . 50

Figure 3.10

Fraction of kinetic energy of 8 layers of PMMA and PC laminate panel
subjected to impact velocity of 300 m/s by a blunt projectile with radius
of 1.73-5 mm and impact velocity of 5 m/s by a spherical projectile with
radius of 6.35-1.73 mm. . . . . . . . . . . . . . . . . . . . . . . . . . . . 51

Figure 3.11

Non-linear optimization of residual kinetic energy of 8-layered PMMA/PC
laminate subjected to impact velocity of 300 m/s by a blunt projectile with
radius of 2.73 mm: a) multiplier of yielding stress of PC as variable; b)
multiplier of stress at failure of PMMA as variable. . . . . . . . . . . . . 52

Figure 3.12

Non-linear optimization of residual kinetic energy of 8-layered PMMA/PC
laminate subjected to impact velocity of 300 m/s by a blunt projectile with
radius of 2.73 mm in terms of two variables: multiplier of yielding stress
of PC and multiplier of stress at failure of PMMA. . . . . . . . . . . . . . 52

Figure 4.1

Geometry and dimensions of dog-bone specimens. (a), specimen orientation diagram (1, longitudinal specimens, loaded along the osteon direction; 2, transverse specimens, loaded perpendicular to the osteon direction); (b), the dimension of the longitudinal specimens (the length of the
transverse specimens is half of the longitudinal specimens) . . . . . . . . 65

Figure 4.2

Representative images showing real time optical effects during tensile
stretching of bovine bone sample. The frame rate of this is 1000 fps. . . . 66

x

Figure 4.3

Correlations of high speed images with stress-strain curves for two typical longitudinal compact bovine femur specimens tested under a constant
strain rate of 0.0001 s−1 : (a) stress-strain curve and the corresponding
high speed images for a typical specimen with a higher failure strain; (b)
stress-strain curve and the corresponding high speed images for a typical
specimen with a lower failure strain; (c) the corresponding light intensity
and the stress relation for the specimen shown in (a); (d) the corresponding light intensity and the stress relation for the specimen shown in (b). . . 68

Figure 4.4

Strain rate effects on the mechanical behavior of the longitudinal compact
bovine femur specimens: (a) stress-strain curve with the corresponding
high speed images for a typical specimen under a strain rate of 0.0012s−1 ;
(b) the corresponding relation of stress and microcracking illumination
for the specimen; (c) the variation of failure strain and ultimate stress with
strain rate, represents Failure strain and represents Ultimate stress;
(d) the variation of Young’s modulus and energy absorption with strain
rate, represents Young’s modulus and represents Energy absorption
capacity. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70

Figure 4.5

The orientation effects on the mechanical properties of the compact bovine
femur specimens under a constant strain rate of 0.0001 s−1 : (a) the stressstrain curves of the longitudinal and transverse specimens; (b) the comparison of the Young’s modulus and ultimate stress of the transverse specimens to the longitudinal specimens; (c) plot of stress and light intensity
changing with time for a typical transverse specimen, represents stress
and solid line represents light intensity; (d) plot of stress and light intensity changing with time for a typical longitudinal specimen, represents
stress and solid line represents light intensity . . . . . . . . . . . . . . . . 71

Figure 4.6

Microstructure effects on the mechanical properties of the longitudinal
compact bovine femur specimens subjected to a constant strain rate of
0.0001 s−1 : (a) SEM fracture surface of a typical specimen with a failure
strain 2.4%; (b) SEM fracture surface of a typical specimen with a failure
strain of 0.9%; (c) plots of failure strain and energy absorption capacity to
the percent area of osteons of the specimens.
represents failure strain,
♦ represents energy absorption capacity. Dot line is the linear trend line
for failure strain; solid line is the linear trend line for energy absorption
capacity. Osteons characterized by cylindrical lamellae surrounding the
Haversian system, which are circled in a, b. Primary bone characterized
by laminar structures (see arrows in b). The area percentages of osteon in
the whole fracture surface were calculated following Saha’s methods [13]. 73

xi

Figure 4.7

Representative SEM post-failure flat surface images (close to the fracture
edge) of four typical longitudinal specimens under a constant strain rate
of 0.0001 s−1 . Example features are labeled: (mac) macrocrack with a
crack length over 500 µm, (mic) microcrack with a crack length smaller
than 500µm. Long linear macrocracks and many microcracks were observed on the specimen with a large failure strain (a), but not on the specimens with small failure strain (c, d). . . . . . . . . . . . . . . . . . . . . 74

Figure 4.8

Representative SEM post-failure images of two typical specimens (under
a strain rate of 0.0001 s−1 ) for orientation evaluation: (a) the flat surface
of the transverse specimen; (b) the flat surface of the longitudinal specimen. Many linear microcracks were observed on the transverse specimens compared to the longitudinal specimen with fewer microcracks;
(c) the corresponding fracture surface of the transverse specimen; (d) the
corresponding fracture surface of the longitudinal specimen. Smoother
tearing surfaces (indicated by arrows) were observed on the transverse
specimen in compared to the rough fracture surface of the longitudinal
specimen. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75

Figure 4.9

Representative SEM post-failure fracture surface images of two typical
longitudinal specimens for the strain rate effects evaluation: (a), (b) osteon pullout of the specimen under a strain rate of 0.0003s−1 ; (c), (d)
osteon-pull out of the specimen under a strain rate of 0.0012 s−1 . Osteon
is comprised of the Haversian canal (pointed out by arrows) surrounded
by cylindrical lamellae. The close view images show that lower strain
rate specimen exhibit rougher fracture surface (the ridge rings pointed
out by arrows (b)) as compared to the higher strain rate specimen with a
smoother tearing surface (see arrow in (d)). . . . . . . . . . . . . . . . . 76

Figure 4.10

Bone specimens: a, A diagram shows the orientation of bone (osteon
grows along bone axis); b, Dimensions of bone specimen; c, A typical
specimen before testing; d, A typical specimen after testing. . . . . . . . . 79

Figure 4.11

Tabletop Charpy impact system. . . . . . . . . . . . . . . . . . . . . . . 80

Figure 4.12

Schematic representative of the table-on Charpy impact system. . . . . . . 81

Figure 4.13

The locations and orders of the viewpoints examined and a sample of
microscope image. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82

Figure 4.14

Impact force history of bovine cortical bone subjected to Charpy impact
loading at various impact energy levels. . . . . . . . . . . . . . . . . . . 84

Figure 4.15

Impact energy history of bovine cortical bone subjected to Charpy impact
loading at various impact energy levels. . . . . . . . . . . . . . . . . . . 84
xii

Figure 4.16

Representative post failure microscopic images at each impact energy
level. The images are labeled with the run case number (E2 Run2 etc),
specimen status after impact (broken or no broken), and the area percentage occupied by the secondary osteon (32%, etc.). Features are labeled
as Mic: microcrack, Osteon: secondary osteon, and Laminar: laminar
structure. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85

Figure 4.17

a, Correlation of mean microcracking numbers to the area percentage occupied by secondary osteon; b, Correlations of microcracking to the impact energy applied. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86

Figure 4.18

SHPB setup for high strain rate testing (300 900 s−1 ). . . . . . . . . . . . 88

Figure 4.19

A typical specimen fixture used in the SHPB testing. . . . . . . . . . . . 88

Figure 4.20

A typical incident and transmitted waves measured in the SHPB testing
of bovine bone specimens. . . . . . . . . . . . . . . . . . . . . . . . . . 89

Figure 4.21

Schematic representation of interfaces between incident bar with sample
and sample with transmission bar. . . . . . . . . . . . . . . . . . . . . . 90

Figure 4.22

Stress-strain relationships of bovine cortical bone subjected to high strain
rate tensile loading by SHPB setup. . . . . . . . . . . . . . . . . . . . . . 91

Figure 4.23

a) Correlation of ultimate stress and strain at ultimate stress to the strain
rate; b)The polynomial regression model of ultimate stresses at different
strain rates. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92

Figure 4.24

Typical stress-strain curves and their corresponding high speed images. . . 93

Figure 4.25

The time history of strains of the specimen showed in 4.24. . . . . . . . . 94

Figure 5.1

Schematic picture of the corner nodes numbering for a RVE in 3D view. . 106

Figure 5.2

a) Microstructure of cortical bone: 1) Osteon, 2) Cement line, 3) Interstitial bone, b) RVE. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108

Figure 5.3

Meshes and dimensions of the C(T) specimen: a) Schematic drawing of
C(T) specimen (Part 2 is the concerned zone with detailed microstructures of osteon bone; Part 1 is represented with bulk bone matrix without
showing detailed microstructures): a=70 mm, W=140 mm, B=0.5 mm;
b) Meshes of the C(T) specimen; c) Schematic drawing of the concerned
zone (Part 2): 1: osteon, 2: cement line, 3: interstitial bone, d) Meshes of
Part 2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110

xiii

Figure 5.4

Computational framework of the concurrent multi-scale modeling scheme. 113

Figure 5.5

Schematic picture of RVE with both local and global coordinate systems. . 114

Figure 5.6

Numerical implementation flow chart of the concurrent multiscale modeling scheme in ABAQUS/EXPLICIT . . . . . . . . . . . . . . . . . . . 116

Figure 5.7

Effective Young’s modulus versus osteon size for cortical bone with 35%
area occupied by secondary osteon: (a) Longitudinal direction (b) Transverse direction. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118

Figure 5.8

Stress contours of a macro-model and the unit cell corresponding to one
of the integration points in macro-model: a) Longitudinal direction; b)
Transverse direction. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119

Figure 5.9

a) Force-displacement curve; b) Calculated fracture toughness. . . . . . . 120

Figure 5.10

Stress contour at the last increment of the four of the cases (Case 1, 3, 4,
and 5) studied in the CT fracture toughening simulations . . . . . . . . . 121

Figure 5.11

Fracture toughness versus cracking extension. . . . . . . . . . . . . . . . 122

Figure 5.12

Geometry and mesh of macro-model of impact testing simulated. . . . . . 123

Figure 5.13

Stress contour of the two micro-models studied at impact time of 0.6 ms. . 124

xiv

Chapter 1
Introduction

1.1

Objective and Expected Significance

Laminated transparent glassy panels have a wide application in both military structures such as
transparent armour and gun shields, as well as civilian structures such as residential and commercial aquariums, automobile windshields, aircraft windows, etc. Implementing an optimal design
to maximize its impact resistance, as well as keeping it lightweight and transparent is of great engineering interest in the design and manufacture of laminated glassy panel. As a fiber reinforced
natural composite material, bone is lightweight, rigid, and tough. The unique biological microstructure of bone provides its special toughening mechanism to achieve its functional integrity
including rigid support and metabolic function. Study of the relationship of its mechanical behavior, loading, and structure helps to develop understanding of the toughening mechanism of bone,
improve its clinical treatment, and inspire new materials that mimic artificial bone.

1.1.1

Specific Aim 1

Transparency, light weight, and high impact resistance are three major factors needed to be considered when designing armour plates. The amorphous polymer, polycarbonate (PC), has been
widely used as armour materials due to abilities to undergo large plastic deformation before failure, its light weight, its good transparency, and its ease of manufacturing. However, due to its

1

ductility, PC is susceptible to penetration with localized plastic deformation in projectile testing.
To increase energy dissipation and avoid localized deformation, commercial glass panels are typically designed as laminated structure built with two materials with contrasting material properties:
ductility and brittleness. The reason is that brittle material reinforces the glassy panel by dissipating energy through cracking, and ductile material strengthens the brittle material by avoiding
shattering into pieces. Glass is one of the typical brittle materials chosen for the laminated glassy
panel due to its transparency. Recently, it was found that poly polymethyl methacrylate (PMMA),
a transparent thermoplastic material, is a good alternative to glass due to its lighter weight and
higher shatter-resistance. The density of PMMA, 1190 Kg/m3 , is the same as PC and is half of
glass, 2400-2800 Kg/m3 . Additionally, the impact resistance of PMMA is almost double that of
glass. The relatively rigid PMMA compensates the large compliance of PC and reinforces PC by
increasing the area involved in damage. Furthermore, PMMA and PC are of almost the same index of refraction, which improves transparency of the glass panel. The PMMA and PC layers are
bonded together with a soft interlayer, impact resistant adhesive.
The mechanical behavior of amorphous polymers is pressure sensitive. Thus, the yield stress
under compression is different from the yield stress under tension. In tension, PMMA and PC
exhibit a significantly contrasting failure mechanism. PC shows higher impact resistance and larger
plastic deformation, and PMMA instead exhibits brittle failure and almost zero plastic deformation
but a higher tensile strength. In compression, these two materials follow a same trend: large plastic
deformation with strain relaxation after yielding, followed by a strain hardening. PC and PMMA
are highly strain rate dependent materials. That means their yielding stress changes with strain
rate. Due to the large plastic deformation, classical continuum mechanics are not accurate enough
to predict the mechanical response of amorphous polymers. Many constitutive models based on
large deformation theory have been developed for amorphous materials. Among the constitutive
2

models, Arruda and Boyce’s [1] three dimensional model based on a rubber elasticity spring system
was the most popular used in modeling the mechanical response of PMMA and PC subjected to
various strain rate.
Impact resistance of a glassy laminated panel depends on many factors, including the hierarchical configuration, thickness of each layer, mechanical properties of each layer, weight and
geometry of the projectile, and impact velocity. The objective of the research work was to experimentally and numerically investigate the toughening mechanism of laminate PMMA and PC plates
subjected to impact testing. First, low-velocity drop weight tests were conducted to evaluate the
threshold condition and failure mechanism of PMMA and PC subjected to impact loading. The experimental methods and results are presented in Chapter 2. Then, the material model based on the
existing Arruda and Boyce’s physical model for amorphous polymers was coded and implemented
in the commercial finite element modeling (FEM) software ABAQUS/EXPLICIT to predict the
strain rate and pressure dependent behavior of PC. The failure behavior of monolithic PMMA and
PC and layered PMMA/PC laminate subjected to low and high speed impact loading were then
numerically investigated by ABAQUS. Lastly, an optimization approach was implemented using
MATLAB and ABAQUS. The detailed methods, results, and finding of the numerical investigation
are presented in Chapter 3.

1.1.2

Specific Aim 2

The second type of material studied is the fiber-reinforced cortical bone. As the major constitution
of the skeleton of vertebrates, bone not only provides rigid support and protection to the soft tissues
and organs but also participates in the fundamental metabolic activity such as the regulation of
calcium. According to their shapes, bone can be classified as long bone, flat bone, short bone, etc.
The femur and tibia are two typical long bones, which are the major support of vertebrate. The
3

diaphysis of long bone is made of the dense cortical bone; and the two ends, called epiphysis, are
made of the phony trabecular bone.
On the microstructural level (from 10-500 mum) [2], mature cortical bone (also called osteon
bone) can be considered a fiber-reinforced composite with fibers, cylindrical osteons, distributed in
the interstitial lamellae matrix. Osteon, the major functional unit of mature cortical bones, grows
along the loading axis (longitudinal direction) and comprises a center Haversian system surrounded
by concentric lamellae. The longitudinal Haversian system surrounds blood vessels and nerves,
which provide nutrition supplies and communicate with osteocyte contained in lacunae through
small transverse canals, canaliculi. Osteons are connected to each other by oblique channels called
Volkmann’s canals. In the nanostructural level (from a few hundred nanometers to 1 mum) [2],
lamellae are formed by type I collagen fibrils and embedded mineral. The boundary between
osteon and surrounding bone tissue, called the cement line, exhibits ductile properties due to its
low mineral composition. Immature bone (laminar bone) is characterized as laminar structure.
During its life span, bone consistently undergoes targeted and non-targeted remodeling in order to
assure its functional integrity. Non-targeted remodeling is a bone adaption to meet the metabolic
needs, and targeted remodeling is triggered by microdamage induced in daily fatigue activities.
Deterioration of bone is a process of damage accumulation in the form of microcracking spanning
a wide range of dimensional scales.
Osteon bone is an orthotropic composite material with a higher Young’s modulus and strength
along the longitudinal direction, which is the growing direction of ostoen or the loading axis, compared to the transverse direction, which is perpendicular to loading axis. Osteon with a central
hole, cement line, and interstitial bone matrix are the three major components of osteon structure. The dimension, properties, micro-arrangement, and density of osteon are predicted to have
a strong correlation to bone’s overall mechanical performance, including microcracking initiation
4

and propagation and remodelling of bone.
The objective of this research work was to investigate the failure and toughening mechanism
of the unique osteon structure. Both FEM modeling and experimental methods were used to conduct the study. First, damage accumulation and effects on the strength and toughness of compact
bovine bone was studied through uniaxial tensile testing, Charpy impact testing, and tensile split
hopkinson bar (SHPB) testing. The experimental procedures and findings are presented in Chapter
4. Second, numerical studies were conducted to investigate the significant role of cement lines,
density and size of osteons, and the properties of each constituent on the fracture mechanism of
cortical bone through multi-scale modeling and XFEM technique in ABAQUS. The corresponding
methods and results are presented in Chapter 5.

1.2

Background Rational and Review of Relevant Literature

In literature, some aspects of the previously mentioned aims have been investigated. This section
provides an overview on current understandings relevant to the problems.

1.2.1

Transparent Laminated Panel

Many impact testing methods have been used to evaluate the impact resistance of transparent armor
materials. The testing methods include low velocity drop weight impact tests [3, 4, 5, 6, 7, 8, 9,
10, 11], high velocity impact tests by use of an air gun [12, 13, 14, 15], and a shock wave impact
by use of explosives [10, 11].
Among them, the low velocity drop weight impact test is popularly used to evaluate the threshold conditions indicating the onset of damage and to investigate the failure mechanism of the tested
material. Liu [3] conducted a low velocity drop weight impact test on a monolithic PMMA plate
5

using an instrumented drop weight impact tester. The variations of force with time and strain energy release rate were examined. The brittle-ductile transition temperature of PMMA was also
evaluated to be in the range of 85◦ C to 96◦ C [3]. A low velocity impact and quasi-static failure of PMMA, based on global and local strain measurements, were conducted by Pearson and
his colleagues [4]. In their work, Fiber Bragg Grating (FBG) sensors, were used to measure local strains, which were combined with the measurements from quasi-static indentation and lowvelocity impact test to examine the failure mechanisms of PMMA. Multiple strikes were applied
to the PMMA sample until failure. Stenzler [16] investigated the impact mechanics of transparent
multi-layered polymer composites using an instrumented intermediate velocity impact test facility. The multi-layered structures with various bonding adhesives were investigated to evaluate the
impact performance [16].
The ballistic impact response of PMMA plate structures attracts high attention due to its wide
application in military defense [13, 17]. Sarva [17] studied the impact of a PMMA /PC (polycarbonate) hierarchical assembly by using a projectile at the impact velocities ranging from 300 to 550
m/s. It was shown that a hierarchical assembly exhibited higher penetration resistance compared to
uniform distributed PMMA discs structure [17]. Hsieh [14] measured the ballistic impact response
of coextruded PC/PMMA multilayered composites; the thickness of the PMMA plate was proved
to play an important role in improving the ballistic performance of the same.
With the development of computational techniques, finite element modeling has been widely
used in the study of the mechanical response of materials subjected to dynamic loading [3, 15, 19,
20, 21]. Finite element models were developed to predict the dynamic response of PMMA plates
[3] and fibrous composite panels [19] subjected to a low velocity drop weight. A computational
model [15] was used to study the ballistic impact response of PMMA plate structures. All the work
shows that finite element simulation is a powerful tool to study the failure and fracture mechanism
6

of materials.
Some study has been done on the optimization of layered PMMA/PC plate panels. The general
procedures they followed includes thoroughly understanding the impact event, quantitatively finding the analytical solutions of the non-stationary strain problem, constructing the corresponding
objective function according to the purpose of the optimization, and solving the objective function
to achieve the desired optimization. The challenge of the optimization problems for the PMMA/PC
layered plates is the difficulty in configuring the analytical solutions due to the unstable and instantaneous nature of brittle cracking. The classical laminate theory uses effective and reasonable
assumptions to reduce the three-dimensional elasticity problems to solvable two-dimensional mechanics of deformable body problems [18]. Thus, the classical laminate theory is widely used to
map the stress states of laminated plate subjected to stationary loading. However, to solve the
layered plates subjected to impact loading, attempts have to be made to refine the classical theory
so that it can be used to solve the dynamic loading problems [22, 23, 24, 25]. Smetankina [22]
used a refined laminate theory on the basis of the hypothesis of the refined theory of Timoshenko
to investigate the strain-stress state of multilayer plates under impulse and impact loading. In their
work the transverse shear was taken into account and the displacement and the external load are
expanded into Fourier series functions satisfying the problem boundary conditions. The numerical
results had a good agreement with the experimental results. Shupikov [23] used Smetankina’s refined laminate theory to investigate the optimal design of multilayer plates. In their work attempts
were made to minimize the stress in layers and the total mass of the multilayered plates. Shokuhfar
[24] investigated the optimal design of smart hybrid composite plates subjected to low-velocity impact using the plate theory developed by Whitney [25]. The interaction between the projectile and
the plate was modeled by a system having two-degrees of freedom, consisting of spring-masses.
While the classical laminate theory are widely used to study the impact loading problem, many
7

publications [26, 27, 28] deal with the dynamic loading problem from the aspect of wave propagation by application of wave theory. Chen [28] developed an analytical solution to the problem of one-dimensional high amplitude wave propagation in layered heterogeneous systems. The
approach of his solution was to convert the initial velocity boundary value problem to a time
dependent stress boundary value problem and the stress time history was obtained by means of
superposition. Luo [27] investigated the stress wave propagation for a layered structure under
impulsive loads. It was shown that the reduced amplitude and elongated pulse duration could be
obtained with proper selection of materials and dimensions of layers when an incident stress pulse
passed through a layered structure. Taha [26] developed a multi-objective optimization approach
for design of blast-resistant composite laminates using carbon nanotubes. A simplified approach
for wave propagation in bounded elastic media was adopted.
The optimal design of layered structure could also be implemented based on the concept
[29, 30] of energy balance, during an impact event. Ray [30] optimized the energy dissipation
characteristics of active constrained layer damping treatments of plates using a rational design approach. In his work, the governing equation was given based on the system energy balance, which
includes the potential energy, kinetic energy, the work done by the piezo-electric control force, and
the work dissipated in the viscoelastic core. Kam [29] investigated the optimal arrangements of
laminated composite plates with maximum stiffness to the side constraints by applying the laminate theory as developed by Whitney. In Kam’s work, the objective function was formulated in
terms of the energy balance. The effects of aspect ratio, length to thickness ratio, and number of
layers on the optimal lamination configuration were studied.
The literature review shows that there has been extensive research focused on the experimental
and FEM study of laminated composite material, as well as optimal design of laminated composite
material based on analytical analysis in the past decades. However, the unique characteristic and
8

failure mechanism of PMMA and PC subjected to impact loading has yet not been fully understood.

1.2.2

Damage Accumulation in Bone

As a tough and lightweight natural composite, bone works as a scaffold to provide support and protect to internal organs. Understanding the failure mechanisms in bone is of great medical interest
for preventing bone fracture failure.
The mechanical properties of compact bone have been widely studied in the past 50 years.
Mcelhane [31] investigated the dynamic response of both fresh human and bovine femur subjected
to dynamic compressive loading at strain rate up to 1000 s−1 . It was found that strain rates between
0.1-1 s−1 were the critical strain rates at which the energy absorption capacity and the maximum
strain to failure exhibit large variation. The Young’s modulus and ultimate stress increased with
the strain rates. Saha and Hayes [32, 33]measured tensile impact properties of human and bovine
compact bone using a pendulum-type instrumented tensile impact tester at a strain rate of 133
s−1 . Tensile strength, energy absorption capacity, and modulus of elasticity were found to be
negatively correlated with the percentage area of secondary osteons. Crowninshiled and Pope [34]
examined the response of compact bovine bone in tension at various strain rates (from 0.01 to 200
s−1 ). Ultimate tensile stress both in longitudinal and transverse directions increased with the strain
rate. Energy capacity reached maximum at a strain rate of 0.1 s−1 . Lewis and Goldsmith [35]
designed a biaxial SHPB apparatus to measure the dynamic material properties of bovine femur in
compression, tension, torsion, and combined torsion and compression. It was found that the bone
exhibited viscoelastic properties in compression, and the fracture stress increased with strain rates.
Currey [36] measured the material properties of compact bovine bone in tension at various strain
rates (from 0.00013 to 0.16 s−1 ). It was found that the strain rate effect on Young’s modulus could
be considered simply as the result of bone being a viscoelastic material, but the strain rate effect
9

on yield strength and breaking strength could not be predicted directly from the viscoelasticity
of bone with using any further facts. Wright and Hayes [45] conducted uniaxial tensile tests of
cortical bovine bone at seven strain rates (from 5.3 x 10−4 to 237 s−1 ). The ultimate strength,
modulus of elasticity, and energy absorption capacity were found to increase with strain rates.
Robertson and Smith [37] investigated the compressive strength of the porcine mandible and its
dependence upon microstructure and strain rate. Similar conclusions as Mcelhane were obtained.
Ferreira [38] investigates the strain rate effects on the Young’s modulus and ultimate strength of
bovine femur bone subjected to compressive loading by SHPB. Ultimate strength was found to
increase with the increase of strain rate. However, Young’s modulus was found to decrease with the
increase of strain rate. The study also showed that experimental results are scattering due to bone’s
highly heterogeneous structure. Pithious et al [39] measured the stress-strain relationship in both
quasi-static and dynamic compressive loading (1000 s−1 ). It was found that the bovine compact
bone was three to four times more brittle under a dynamic load than under a quasi-static one.
Adharapurapu et al. [40] used a two-bar three point bending apparatus, modified from SHPB, to
measure the dynamic fracture toughness of bovine bone at high strain rates. The dynamic fracture
toughness was significantly lower than the static one. The cleavage-like fracture was observed on
the fracture surface at very high loading rates of 2 x 105 s−1 . Hansen et al [41] examined the strain
rate effects on the material properties of human cortical bone in compression with strain rates from
0.14 to 29.1 s−1 , and in tension with strain rates from 0.08 to 17 s−1 . There was a significant
difference found between the trend of the ultimate stress and failure strain in compression and
tension. The ultimate stress in compression increased with strain rate, but the ultimate stress in
tension decreased with strain rate. Stress and strain at yield decreased with strain rate for both
compression and tension. However, the effects of strain rates on the post yield properties of bone
not following the simple linear relationship implied that strain rate has a stronger effect on post
10

yield deformation than on initiation of yielding.
The morphology and composition of the osteonal cement line in human bone were studied by
Schaffler [42]. His results showed that collagen is absent in cement lines, which contain significantly less calcium and phosphorus but more sulfur than surrounding bone matrix. The intrinsic
and extrinsic mechanisms were widely accepted to be the major toughening mechanisms for bone
subjected to quasi-static loading. Intrinsic mechanisms [43] are microstructural damage mechanisms that operate ahead of the crack tip, and extrinsic mechanisms [43] act to shield the crack
from the applied driving force and operate principally behind the crack tip in the crack wake. In
detail, the toughening mechanism can be classified into crack deflection by osteons, crack bridging
by collagen fibers, uncracked-ligament bridging, diffuse microcracking, and void. Vanishth et al.
[44, 45] revealed a rising R-curve behaviour in bone. Rising R-curve behaviour is a direct result
of the extrinsic toughening mechanisms. Kulin et al. [46] studied the effects of aging and loading
rate on the compressive failure, fracture toughness, and resistance of equine cortical bone by use
of split-Hopkinson pressure bar (SHPB). The post fracture surface analysis by use of both of the
confocal microscopy and scanning electronic microscopy (SEM) showed a transition of ductile to
brittle. It was attributed to the viscoelastic behavior of collagen which comprised 50 percent of
bone structure. However, the failure mechanism of bone in high strain rate is still not fully uncovered due to the challenge of recording the high speed fracture process during high strain rate
loading. It is pertinent to further investigate the toughening mechanism of bone subjected to high
strain rate loading.
Although the existing research has shown that the mechanical properties of bone have some
correlation to its microstructure, the role of cement lines and osteons on the overall mechanical
properties of cortical bone and the initiation and propagation of microcracks has not been discovered yet. Thus, the proposed work was to conduct a series of research through experimental
11

investigation and FEM analysis to address the above issues.

12

REFERENCES

13

REFERENCES
[1] Ellen M Arruda and Mary C Boyce. Evolution of plastic anisotropy in amorphous polymers
during finite straining. International Journal of Plasticity, 9(6):697–720, 1993.
[2] Jae-Young Rho, Liisa Kuhn-Spearing, and Peter Zioupos. Mechanical properties and the
hierarchical structure of bone. Medical engineering & physics, 20(2):92–102, 1998.
[3] Yanxiong Liu and Benjamin Liaw. Drop-weight impact tests and finite element modeling of
cast acrylic/aluminum plates. Polymer Testing, 28(8):808–823, 2009.
[4] JD Pearson, MA Zikry, M Prabhugoud, and K Peters. Measurement of low velocity and
quasi-static failure modes in pmma. Polymer composites, 28(3):381–391, 2007.
[5] Herzl Chai and Guruswami Ravichandran. On the mechanics of fracture in monoliths and
multilayers from low-velocity impact by sharp or blunt-tip projectiles. International journal
of impact engineering, 36(3):375–385, 2009.
[6] Md Akil Hazizan and WJ Cantwell. The low velocity impact response of foam-based sandwich structures. Composites Part B: Engineering, 33(3):193–204, 2002.
[7] Ik Hyeon Choi and Cheol Ho Lim. Low-velocity impact analysis of composite laminates
using linearized contact law. Composite structures, 66(1):125–132, 2004.
[8] Md Akil Hazizan and WJ Cantwell. The low velocity impact response of an aluminium
honeycomb sandwich structure. Composites Part B: Engineering, 34(8):679–687, 2003.
[9] Michael Hebert, Carl-Ernst Rousseau, and Arun Shukla. Shock loading and drop weight
impact response of glass reinforced polymer composites. Composite Structures, 84(3):199–
208, 2008.
[10] Srinivasan Arjun Tekalur, Arun Shukla, and Kunigal Shivakumar. Blast resistance of polyurea
based layered composite materials. Composite Structures, 84(3):271–281, 2008.
[11] Shinji Ogihara, Tomoyuki Ishigure, and Akira Kobayashi. Study on impact perforation fracture mechanism in pmma. Journal of materials science letters, 17(8):691–692, 1998.
[12] S Sahraoui and JL Lataillade. Deformation and fracture of pmma at high rates of loading.
Journal of applied polymer science, 51(9):1527–1532, 1994.

14

[13] L Roy Xu and Ares J Rosakis. An experimental study of impact-induced failure events in
homogeneous layered materials using dynamic photoelasticity and high-speed photography.
Optics and Lasers in Engineering, 40(4):263–288, 2003.
[14] Alex J Hsieh and John W Song. Measurements of ballistic impact response of novel coextruded pc/pmma multilayered-composites. Journal of reinforced plastics and composites,
20(3):239–254, 2001.
[15] Sai Sarva, Adam D Mulliken, Mary C Boyce, and Alex J Hsieh. Mechanics of transparent
polymeric material assemblies under projectile impact: simulations and experiments. Technical report, DTIC Document, 2004.
[16] Joshua S Stenzler and NC Goulbourne. Impact mechanics of transparent multi-layered polymer composites. In Society for experimental mechanics-SEM annual conference and exposition on experimental and applied mechanics, page V3, 2009.
[17] S Sarva, AD Mulliken, and MC Boyce. The mechanics of large-strain inhomogeneous deformation of polymeric materials under dynamic loading conditions. In Journal de Physique IV
(Proceedings), volume 134, pages 95–101. EDP sciences, 2006.
[18] YM Gupta and JL Ding. Impact load spreading in layered materials and structures: concept and quantitative measure. International Journal of Impact Engineering, 27(3):277–291,
2002.
[19] U Faroop and Karl Gregory. Computational modeling and simulation of low velocity impact
on fibrous composite panels drop weight un partitioned model. Journal of Engineering and
Applied Sciences, pages 1819–6608, 2009.
[20] Myungsoo Park, Jeonghoon Yoo, and Dong-Teak Chung. An optimization of a multi-layered
plate under ballistic impact. International journal of solids and structures, 42(1):123–137,
2005.
[21] Pizhong Qiao, Mijia Yang, and Florin Bobaru. Impact mechanics and high-energy absorbing
materials: review. Journal of Aerospace Engineering, 21(4):235–248, 2008.
[22] NV Smetankina, S Yu Sotrikhin, and AN Shupikov. Theoretical and experimental investigation of vibration of multilayer plates under the action of impulse and impact loads. International journal of solids and structures, 32(8):1247–1258, 1995.
[23] AN Shupikov, SV Ugrimov, AV Kolodiazhny, and VG Yareschenko. High-order theory
of multilayer plates. the impact problem. International journal of solids and structures,
35(25):3391–3403, 1998.
15

[24] A Shokuhfar, SMR Khalili, F Ashenai Ghasemi, K Malekzadeh, and S Raissi. Analysis and
optimization of smart hybrid composite plates subjected to low-velocity impact using the
response surface methodology (rsm). Thin-Walled Structures, 46(11):1204–1212, 2008.
[25] JM Whitney and NJ Pagano. Shear deformation in heterogeneous anisotropic plates. Journal
of Applied Mechanics, 37(4):1031–1036, 1970.
[26] Mahmoud M Reda Taha, Arife B Colak-Altunc, and Marwan Al-Haik. A multi-objective optimization approach for design of blast-resistant composite laminates using carbon nanotubes.
Composites Part B: Engineering, 40(6):522–529, 2009.
[27] Xiaobo Luo, Amjad J Aref, and Gary F Dargush. Analysis and optimal design of layered
structures subjected to impulsive loading. Computers & Structures, 87(9):543–551, 2009.
[28] X Chen, N Chandra, and AM Rajendran. Analytical solution to the plate impact problem
of layered heterogeneous material systems. International journal of solids and structures,
41(16):4635–4659, 2004.
[29] TY Kam and FM Lai. Maximum stiffness design of laminated composite plates via a constrained global optimization approach. Composite structures, 32(1):391–398, 1995.
[30] MC Ray and A Baz. Optimization of energy dissipation of active constrained layer damping
treatments of plates. Journal of Sound and Vibration, 208(3):391–406, 1997.
[31] J H McElhaney. Dynamic response of bone and muscle tissue. Journal of Applied Physiology,
21(4):1231–1236, 1966.
[32] S Saha and WC Hayes. Tensile impact properties of human compact bone. Journal of biomechanics, 9(4):243–251, 1976.
[33] S Saha and WC Hayes. Relations between tensile impact properties and microstructure of
compact bone. Calcified tissue research, 24(1):65–72, 1977.
[34] RD Crowninshield and MH Pope. The response of compact bone in tension at various strain
rates. Annals of Biomedical Engineering, 2(2):217–225, 1974.
[35] JL Lewis and W Goldsmith. The dynamic fracture and prefracture response of compact bone
by split hopkinson bar methods. Journal of biomechanics, 8(1):27–40, 1975.
[36] JD Currey. The effects of strain rate, reconstruction and mineral content on some mechanical
properties of bovine bone. Journal of Biomechanics, 8(1):81–86, 1975.
16

[37] Diane Margel Robertson and Dennis C Smith. Compressive strength of mandibular bone as a
function of microstructure and strain rate. Journal of biomechanics, 11(10):455–471, 1978.
[38] F Ferreira, MA Vaz, and JA Simoes. Mechanical properties of bovine cortical bone at high
strain rate. Materials characterization, 57(2):71–79, 2006.
[39] Martine Pithioux, Damien Subit, and Patrick Chabrand. Comparison of compact bone failure
under two different loading rates: experimental and modelling approaches. Medical engineering & physics, 26(8):647–653, 2004.
[40] Raghavendra R Adharapurapu, Fengchun Jiang, and Kenneth S Vecchio. Dynamic fracture
of bovine bone. Materials Science and Engineering: C, 26(8):1325–1332, 2006.
[41] Ulrich Hansen, Peter Zioupos, Rebecca Simpson, John D Currey, and David Hynd. The effect
of strain rate on the mechanical properties of human cortical bone. Journal of biomechanical
engineering, 130(1):011011, 2008.
[42] Mitchell B Schaffler, David B Burr, and Richard G Frederickson. Morphology of the osteonal
cement line in human bone. The Anatomical Record, 217(3):223–228, 1987.
[43] RK Nalla, JS St¨olken, JH Kinney, and RO Ritchie. Fracture in human cortical bone: local
fracture criteria and toughening mechanisms. Journal of biomechanics, 38(7):1517–1525,
2005.
[44] D Vashishth, KE Tanner, and W Bonfield. Contribution, development and morphology of microcracking in cortical bone during crack propagation. Journal of Biomechanics, 33(9):1169–
1174, 2000.
[45] D Vashishth, KE Tanner, and W Bonfield. Experimental validation of a microcracking-based
toughening mechanism for cortical bone. Journal of Biomechanics, 36(1):121–124, 2003.
[46] Robb M Kulin, Fengchun Jiang, and Kenneth S Vecchio. Loading rate effects on the r-curve
behavior of cortical bone. Acta biomaterialia, 7(2):724–732, 2011.

17

Chapter 2
Impact Behavior and Dynamic Failure of
PMMA and PC Plates

2.1

Introduction

In this chapter, low velocity drop weight tests were conducted on the monolithic PMMA and bilayered PMMA and PC at various impact conditions to evaluate the threshold conditions indicating
the onset of damage and failure mechanism of PMMA and PC plates. The experimental results
will also serve as the validation for the numerical study in Chapter 3. The failure mechanism
of monolithic PMMA was first studied through drop-weight impact tests with impact conditions
of constant impact mass with various impact velocities, constant energy with various impact velocities and mass, and constant velocities with various impact mass. The threshold conditions,
energy dissipation, and crack pattern were evaluated. Then, the bonded and unbonded bi-layered
PMMA and PC structures was tested to investigate the effects of the interfacial bonding on the
overall impact performance of layered PMMA plates. The failure mechanism was transversely
compared through the monolithic PMMA plate, the bonded PMMA layered structures to the unbonded PMMA layered structures.

18

Table 2.1
and physical
properties
of PMMA
Table 1. Mechanical
andMechanical
physical properties
of PMMA
and PC plate
Specimen

Density
( kg/mm3)

PC
PMMA

1.2x10-6
1.19x10-6

Modulus of
Elasticity
( MPa )
2300
3100

Poisson’s
Ratio

Tensile
strength
( MPa )
65.5
74.4

0.3
0.33

and PC plate.

Shear
Strength
( MPa )
41.4
68.9

Table 2 . Threshold energy of mono-PMMA plates under constant impact mass
Impact Velocities (m/s)
0.7
1.0
2.0
Average Threshold Energy (J)
2 .0 r 0
3.1 r 0.88
2.8 r 0.59
Confidence (%)
95
95
95

2.2

Experimental Procedures

Compression
Strength
( MPa )
86
120.7

3.0

2.1 r 0.09
95

The drop-weight impact tests were carried out through an Instron Dynatup 9250 drop weight tower
Table 3 . Threshold energy of mono-PMMA plates under constant impact energy 12 J
Impact Velocities (m/s)
1.88
1.35
1.11
0.94

(Fig. 2.1), which was equipped with a load cell above the tup to record the force history during
Drop Weight (Kg)

6.61

12.91

18.57

25.57

2.4 r 0.25on the tup
2.4 rinitial
0.29 velocity
Average
Threshold detector
Energy (J)was mounted
2.6 to
r 0measure
.59
2.5 r 0of
.49the impact.
the test.
A velocity
the
Confidence (%)

95

95

95

95

During the test, the 12.7 mm diameter hemispherical tup was raised to a specific height above
Table 4 . Summary of impact scenarios of bi-layered structures
(m/s) MassThe
(Kg)
the Models
specimen and was allowed to drop the weight ontoVelocity
the specimen.
PC/PMMA_NA_CM_1

0.70

7.25

Energy (J) range
drop-weight
1.78

PC/PMMA_NA_CM_2
1.00
7.25commercial
3.55 PMMA and
was 2.7-14.7 kg for the low weight
drop. The tested specimens
were
PC/PMMA_NA_CM_3
PC/PMMA_NA_CM_4
which were square plates (100 mm x
PC/PMMA_AD_CM_1
PC/PMMA_AD_CM_2
plate and 3.2 mm forPC/PMMA_AD_CM_3
each of the plates

PC on the top of PMMA

PC plates,

monolithic

2.00
7.25
14.50
3.00
7.25
32.63
100 mm) with a thickness of 6.4 mm for the
1.00
7.25
3.55
2.00
7.25
14.50
in a layered
configuration.
The 32.63
mechanical and
3.00
7.25

PMMA/PC_NA_CM_1

0.70

7.25

1.78

PMMA/PC_NA_CM_3

2.00

7.25

14.50

physical properties are summarized
in Table 2.1. When the
PMMA/PC_NA_CM_2
1.00impact tower
7.25was armed,
3.55the specimen
wasPMMA
fixed on
bythe
a hydraulic
with a diameter of 76.2 3.00
mm. To avoid
double impact,
top of PC clampPMMA/PC_NA_CM_4
7.25
32.63 a rebound
PMMA/PC_AD_CM_1

1.00

7.25

3.55

PMMA/PC_AD_CM_2
7.25were performed
14.50
brake cylinder caught the tup after
the first impact. All 2.00
the experiments
at room
PMMA/PC_AD_CM_3

3.00

7.25

32.63

PMMA/PMMA_NA_CM_1
0.70 results were
7.25 obtained
1.78
temperature. The same tests were
repeated twice, and similar
repeatedly.

PMMA on the top of PMMA

PMMA/PMMA_NA_CM_2
PMMA/PMMA_NA _CM_3
PMMA/PMMA_NA _CM_4
PMMA/PMMA_AD_CM_1
PMMA/PMMA_AD_CM_2
PMMA/PMMA_AD_CM_3

2.3

Experimental Results

2.3.1

Monolithic PMMA Plates

1.00
2.00
3.00
0.70
1.00
2.00

7.25
7.25
7.25
7.25
7.25
7.25

3.55
14.50
32.63
1.78
3.55
14.50

To examine the failure mechanism of the mono-PMMA plates, multiple tests were conducted. The
combinations included constant impact mass (6.61 kg) with various impact velocities, constant
19

(a)

(b)

Figure
1. Instron
Dynatup
9250(a)
drop
weight
tower: (b) after clamping.
Figure 2.1 Instron Dynatup
9250
drop weight
tower:
before
clamping;
(a) before clamping; (b) after clamping.

energy (12 J) with various impact velocities and mass, and constant velocities (0.7 m/s and 1.0
m/s) with various impact mass.

Figure 2. Post-failure pictures of mono-PMMA plates under constant impact mass. The
velocities are 0.7 m/s, 1.0 m/s, 2.0 m/s and 3.0 m/s respectively from left to right.

Figure 2.2 Post-failure 4000
pictures of mono-PMMA plates under constant impact mass: velocities are
0.7 m/s, 1.0 m/s, 2.0 m/s,
and 3.0 m/s respectively from left to right.
3500
Reaction Force (N)

CM_V=0.7m/s

The post-failure

3000

CM_V=2.0m/s

2500
photographs

of the monolithic PMMA plates under constant drop weights are

2000

shown in Fig. 2.2. It was observed that no visible crack was induced on the specimen under the
1500

impact velocity of 0.7 1000
m/s, and 5 radial cracks propagated on the specimen under the velocities
500

of 1.0 m/s and 2.0 m/s respectively. And 6 radial cracks were induced on the specimen under the
0

0

2

4

6

8

10

velocities of 3.0 m/s. When the impact velocity reached 5.0 m/s, the specimen was penetrated.
Time (ms)

The curves of force versus time and force versus displacement of the specimen under the impact
Figure 3. Reaction force versus time of mono-PMMA plates under constant impact mass

velocities of 0.7 m/s and 2.0 m/s are shown in Fig. 2.3 and Fig. 2.4 respectively. A smooth force
against time curve was obtained when the velocity of 0.7 m/s was applied. However, when the
20

specimen was subjected to the impact velocity of 2.0 m/s, the force curve first gradually reached
a peak point and then underwent a sudden drop followed by oscillations and eventually rested on
zero. The first sudden drop of the force curve indicated the onset of damage. Fig. 2.4 shows that the
(a)

(b)

unloading process of the specimen under an impact velocity of 0.7 m/s did not fully follow the same
Figure 1. Instron Dynatup 9250 drop weight tower:
(a) before clamping; (b) after clamping.

track of the loading process. The phenomenon indicated that the PMMA plate exhibited elastic
hysteresis behavior under impact loading. The corresponding hysteresis energy, the work done
during the loading and unloading cycles, was calculated to be around 0.7 J. When the specimen was
broken, the area under the force-displacement curve before the cracking point was the threshold
energy. The threshold energy of the specimen subjected to an impact velocity of 2.0 m/s was
Figure 2. Post-failure pictures of mono-PMMA plates under constant impact mass. The

calculated
to are
be 2.8
(Fig.
velocities
0.7 Jm/s,
1.02.4).
m/s, 2.0 m/s and 3.0 m/s respectively from left to right.

4000

Reaction Force (N)

3500

CM_V=0.7m/s

3000

CM_V=2.0m/s

2500
2000
1500
1000
500
0
0

2

4

6

8

10

Time (ms)

Figure 3. Reaction force versus time of mono-PMMA plates under constant impact mass

Figure 2.3 Reaction force versus time of mono-PMMA plates under constant impact mass of 6.61
kg.

The peak force and time to crack of the monolithic PMMA plates under constant impact mass
are shown in Fig. 2.5. The plots showed that higher velocities induced higher reaction force and
needed less time to initiate crack. The time to crack was around 3.0 ms, 2.1 ms, 0.76 ms, 0.73,
and 0.26 ms respectively when the specimens were subjected to the impact velocities of 0.7 m/s,
21

4000
3500

CM_V=0.7m/s

3000

CM_V=2.0m/s

2500
2000
1500
1000
500
0
0

2

4

6

8

10

12

Displacement (mm)

Figure 2.4 Reaction force versus displacement of mono-PMMA plates under
constant impact mass.
4000
3500

1.0 m/s, 2.0 m/s, 3.0 m/s, and 5 m/s respectively. The corresponding threshold energies examined
3000

2500
were in the range of 1.9 - 2.9 J, as summarized in Table 2.2. When the
impact energy reached
2000

the threshold value, cracks initiated on the bottom surface of the PMMA plates and propagated
1500

to the top due to the maximum tensile stress induced on the bottom caused
1000by bending. Fig. 2.6
500

showed the peak force and time to crack of the monolithic specimen under constant velocities, 0.7
0

m/s and 1.0 m/s. The curves also show that higher velocities induced higher 0reaction force
2 and
need less time to initiate crack. The peak force and time to crack of the monolithic PMMA plate
under constant energy, 12 J, are shown in Fig. 2.7. The curves showed that higher drop-weight
and lower impact velocity needed a longer time to fail compared to lower drop-weight and higher
impact velocity. Thus, it can be concluded that the dynamic failure of PMMA is more sensitive
to impact velocity than impact mass. The threshold energy of the monolithic PMMA plate under
constant impact energy was kept in the range of 2.4 to 2.6 J, as summarized in Table 2.3.

22

4

Time (

4500

3.5

4000

3

3500
2.5

3000

2

2500
2000

Peak force

1500

Time to crack

1.5
1

1000
0.5

500
0

0
0

1

2
3
Velocity (m/s)

4

5

6

4000

4

3500

3.5

3000

3

2500

2.5

2000

2

1500

1.5

Peak force V=0.7 m/s
Peak force_V=1.0 m/s
Time to crack V=0.7 m/s
Time to crack_V=1.0 m/s

1000
500

Time (ms)

Force (N)

Figure 2.5 Peak force and time to crack of mono-PMMA plates under constant impact mass.

1
0.5

0

0
0

5

10
Mass (Kg)

15

20

25

Figure
Peak
force
andtime
timetotocrack
crack of
of mono-PMMA
under
constant
impact
velocities
Figure
2.6 6.
Peak
force
and
mono-PMMAplates
plates
under
constant
impact
velocities
m/s)
(0.7 (0.7
m/s m/s
andand
1.01.0
m/s).
2.5

5000

Table 2.2 Threshold energy of mono-PMMA plates under constant impact mass.
4500

2

4000

0.7
N/A

1.0
2.0
3.0
5.0
2.9 ± 0.21 2.8 ± 0.28 2.1 ± 0.07 1.9 ± 0.02
1.5

3000
2500

Impact
Velocities (m/s)
2000
Drop Weight (kg)
1500
Average Threshold Energy (J)
1000
500

Impact Velocities
(m/s)
0
Drop Weight
(kg)
0
5
Average Threshold Energy (J)

1.88
1.35
1.11
0.94 1
6.61
12.91
18.57
25.57
Peak 2.6
force
2.4 ± 0.18
± 0.42 2.4 ± 0.21 2.5 ± 0.35
Time to crack

23

Time (ms)

Force (N)

Impact Velocities (m/s)
3500
Average Threshold Energy (J)

0.5

0.7
0.7
0.7
1.0 0
1.0
1.0
6.61
12.91
18.57
6.64
12.91
18.57
10
15
20
25
30
N/A (Kg)2.0 ± 0.16 2.0 ± 0.11 2.8 ± 1.34 2.5 ± 0.07 2.4 ± 0.11
Mass

0

5

10
Mass (Kg)

15

20

25

Figure 6. Peak force and time to crack of mono-PMMA plates under constant impact velocities
(0.7 m/s and 1.0 m/s)
2.5

5000
4500

2

4000

1.5

3000
2500
2000

1

1500

Time (ms)

Force (N)

3500

Peak force

1000

0.5

Time to crack

500
0

0
0

5

10
15
Mass (Kg)

20

25

30

Figure
7. Peak
force
and
time
crackof
of mono-PMMA
mono-PMMA
plates
constant
impact
Figure
2.7 Peak
force
and
time
totocrack
platesunder
under
constant
impact
energy.
Impact
Velocities
(m/s)
0.7
1.0
2.0
3.0 energy
5.0

Average Threshold Energy (J)

N/A

2.9 ± 0.21 2.8 ± 0.28 2.1 ± 0.07 1.9 ± 0.02

Table 2.3 Threshold energy of mono-PMMA plates under constant impact energy 12 J.
Impact Velocities (m/s)
Drop Weight (kg)
Average Threshold Energy (J)

2.3.2

1.88
1.35
1.11
0.94
6.61
12.91
18.57
25.57
2.4 ± 0.18 2.6 ± 0.42 2.4 ± 0.21 2.5 ± 0.35

Layered
ImpactPlates
Velocities (m/s)

0.7
0.7
0.7
1.0
1.0
1.0
Drop Weight (kg)
6.61
12.91
18.57
6.64
12.91
18.57
Average
Threshold
Energy
(J)
N/A
2.0
±
0.16
2.0
±
0.11
2.8
±
1.34
2.5
±
0.07
2.4
± 0.11
The mechanical performance of the composite materials can be changed or improved by increasing

or decreasing the comprising layers, or by using different materials in each layers. In the current
study, we seek to examine how the ductile PC affects on the impact performance of brittle PMMA
plate, and furthermore how the interstitial layer affects on the overall impact behavior of the layered
structures. To achieve these objectives, various configurations were designed: PMMA on the top
of PC plates with weak and strong bonding, PC on the top of PMMA plates with weak and strong
bonding, and PMMA on the top of PMMA plates with weak and strong bonding. The weak
bonding was created by placing two plates together without adhesive. The strong bonding was
achieved by using Loctite E-30CL glass bonder epoxy between layers. The tensile strength of
the E-30CL epoxy is 55 MPa. When cured for 5 days at 25◦ C, the shear strength of the epoxy
24

PC
PMMA

1.2x10-6
1.19x10-6

2344
3100

0.38
0.38

65.5
74.4

41.4
68.9

86
120.66

Table 2 . Threshold energy of mono-PMMA plates under constant impact mass
Impact Velocities (m/s)

0.7

1.0

2.0

3.0

Confidence (%)

95

95

95

95

was Average
13.4 MPa
with the
substrate
3.1MPa
 0.88with the2.substrate
2of
. 0 PC,
 0 and 2.1
8  0.59 of PMMA.
Threshold
Energy
(J)
2.1  0.The
09 adhesive
was dispensed using a manual dispensing gun and was cured at room temperature. The layered
Table 3 . Threshold energy of mono-PMMA plates under constant impact energy 12 J

combinations included PMMA on the top of PC plates, PC on the top of PMMA plates, and biImpact Velocities (m/s)

1.88

1.35

1.11

0.94

layered
plate structures. The6.61
impact scenarios
are summarized
DropPMMA
Weight (Kg)
12.91of the bi-layered
18.57 structures
25.57
Average Threshold Energy (J)

in Table
2.4. (%)
Confidence

2.4  0.25

2.6  0.59

2.4  0.29

2.5  0.49

95

95

95

95

Table 2.4 Summary of impact scenarios of bi-layered structures.

Table 4 . Summary of impact scenarios of bi-layered structures
Models

Velocity (m/s)

Mass (Kg)

Energy (J)

PC on the top of PMMA

PC/PMMA_NA_CM_1
PC/PMMA_NA_CM_2
PC/PMMA_NA_CM_3
PC/PMMA_NA_CM_4
PC/PMMA_AD_CM_1
PC/PMMA_AD_CM_2
PC/PMMA_AD_CM_3

0.70
1.00
2.00
3.00
1.00
2.00
3.00

7.25
7.25
7.25
7.25
7.25
7.25
7.25

1.78
3.55
14.50
32.63
3.55
14.50
32.63

PMMA on the top of PC

PMMA/PC_NA_CM_1
PMMA/PC_NA_CM_2
PMMA/PC_NA_CM_3
PMMA/PC_NA_CM_4
PMMA/PC_AD_CM_1
PMMA/PC_AD_CM_2
PMMA/PC_AD_CM_3

0.70
1.00
2.00
3.00
1.00
2.00
3.00

7.25
7.25
7.25
7.25
7.25
7.25
7.25

1.78
3.55
14.50
32.63
3.55
14.50
32.63

PMMA on the top of PMMA

PMMA/PMMA_NA_CM_1
PMMA/PMMA_NA_CM_2
PMMA/PMMA_NA _CM_3
PMMA/PMMA_NA _CM_4
PMMA/PMMA_AD_CM_1
PMMA/PMMA_AD_CM_2
PMMA/PMMA_AD_CM_3

0.70
1.00
2.00
3.00
0.70
1.00
2.00

7.25
7.25
7.25
7.25
7.25
7.25
7.25

1.78
3.55
14.50
32.63
1.78
3.55
14.50

2.3.2.1

Bi-layered Structures of PMMA on the Top of PC

The post-failure photographs of the bi-layered structures of PMMA on the top of PC with weak
bonding are shown in Fig. 2.8. The results showed that no visible crack was observed on both PC
and PMMA plates when the specimens were subjected to the impact velocities of 0.7 m/s and 1.0
m/s respectively. However, when higher velocities (2.0 m/s, 3.0 m/s, and 5.0 m/s), were applied,
25

4 to 5 radial cracks were observed on the PMMA plates. The PC plates were kept unbroken for
all the velocities. Furthermore, no ductile deformation was observed on the PMMA plates. The
corresponding forces against time curves and impact energy against time curves are shown in Fig.
2.9. The unbroken specimens were indicated by the smooth curves.

Figure 2.8 Post-failure pictures of bi-layered structures with PMMA on the top of PC plates under
constant impact mass: the velocities are 0.7 m/s, 1.0 m/s, 2.0 m/s, 3.0 m/s, and 5 m/s respectively
from left to right; the first row is the top PMMA plates and the second row is the bottom PC plates.

a)

b)
30

PMMA/PC_0.7 m/s
PMMA/PC_1 m/s
PMMA/PC_2 m/s
PMMA/PC_3 m/s
PMMA/PC_5 m/s

2500
2000
1500

PMMA/PC_0.7 m/s
PMMA/PC_1 m/s
PMMA/PC_2 m/s
PMMA/PC_3 m/s
PMMA/PC_5 m/s

25
20
15

1000

10

500

5

0

0
0

2

4

6
8
Time (ms)

10

12

0

2

4

6
Time (ms)

8

10

12

Figure 2.9 a)Reaction force history of bi-layered structures with PMMA on the top of PC with weak
bonding under constant impact mass 6.61 kg;b)Impact energy history of bi-layered structures with
PMMA on the top of PC with weak bonding under constant impact mass 6.61 kg

26

2.3.2.2

Bi-layered Structures of PC on the Top of PMMA

Fig. 2.10 shows the reaction force and time and impact energy and time relationships of the bilayered structures of PC on the top of PMMA with weak bonding under constant impact mass 6.61
Kg. Compared to the bi-layered structures of PMMA on the top of PC, only one smooth curve was
obtained when placing PC on the top of PMMA which means cracks occurred on all specimens
except one under the impact velocity of 0.7 m/s. The phenomenon also can be observed on the
post failure pictures, Fig. 2.11 . When the velocities were larger than 0.7 m/s, radial cracks were
induced on all the PMMA plates. No crack was observed on the PC plates.

a)

b)
PC/PMMA_0.7 m/s
PC/PMMA_1 m/s
PC/PMMA_2 m/s
PC/PMMA_3 m/s
PC/PMMA_5 m/s

1600
1400
1200
1000
800
600
400
200
0
0

2

4

6
8
Time (ms)

10

PC/PMMA_0.7 m/s
PC/PMMA_1 m/s
PC/PMMA_2 m/s
PC/PMMA_3 m/s
PC/PMMA_5 m/s

20
18
16
14
12
10
8
6
4
2
0
0

12

2

4

6
8
Time (ms)

10

12

Figure 2.10 a)Reaction force history of bi-layered structures with PMMA on the top of PC with
weak bonding under constant impact mass 6.61 kg; b)Impact Energy history of bi-layered structures with PMMA on the top of PC with weak bonding under constant impact mass 6.61 kg;

2.3.2.3

Bi-layered Structures of PMMA on the Top of PMMA

The post-impact pictures of the bi-layered PMMA/PMMA structures with weak bonding are shown
in the Fig. 2.12. Brittle failure, revealed by radial cracking pattern, was observed on each PMMA
plate of the bi-layer PMMA/PMMA structures with weak bonding when the impact velocity was
larger than 0.7 m/s.

27

Figure 2.11 Post-failure pictures of bi-layered structures with PC on the top of PMMA plates under
constant impact mass: the velocities are 0.7 m/s, 1.0 m/s, 2.0 m/s, 3.0 m/s, and 5 m/s respectively
from left to right; the first row is the top PC plates and the second row is the bottom PMMA plates.

Figure
Pos-failurepictures
picturesof
ofbi-layered
bi-layered PMMA
PMMA plates
Figure
2.1212.Pos-failure
platesunder
underconstant
constantimpact
impactmass:
mass:thethe vevelocities
are
0.7
m/s,
1.0
m/s,
2.0
m/s
and
3.0
m/s
respectively
from
left
to
right;
the row
first is the
locities are 0.7 m/s, 1.0 m/s, 2.0 m/s and 3.0 m/s respectively from left to right; the first
row is the top PMMA plates and the second row is the bottom PMMA plates.
top PMMA plates and the second row is the bottom PMMA plates.

2.4

Discussion

The peak force and the time to crack of bi-layered PMMA and PC plates with weak bonding are
summarized in Table 2.5. The results show that in the same impact conditions PMMA on the

28

Table 2.5 Threshold energy, energy dissipated, and time to crack of bi-layered PMMA and PC
plates with weak bonding under constant impact mass 6.61 kg.

Specimen
PMMA/PC
PC/PMMA
PMMA/PMMA_1
PMMA/PMMA_2
PMMA/PC
PMMA/PMMA_3
PC/PMMA
PMMA/PMMA_4
PMMA/PMMA_5
PMMA/PC
PC/PMMA
PMMA/PC
PC/PMMA
PMMA/PC
PC/PMMA

Velocity (m/s) Threshold Energy (J)
0.7
N/A
0.7
0.7
2400 1.5
1
2300 N/A
1.0
1.39
2500
1.0
2
2230 1.6
3
3170
2.0
2.0
2.0
1.7
3.0
1.9
3.0
1.6
5.0
0.8
5.0
0.7

Time to Crack (ms)
N/A
3.1 4.0
2.2 N/A
0.6
0.55 2.8
0.4 2.3
1.4
1.7
1.3
0.3
0.3

Energy Dissipated (J)
2.1
2.1
2.4
3.1
2.3
2.2
3.5
2.7
2
3.2
2.1
1.1
3.2
0.4
6.9
6.5
12.4
9.6
24.1
17.8

Table 2.6 Peak force and time to crack of bi-layered PMMA and PC plates with strong bonding
under constant impact mass 6.61 kg.
Specimens
PMMA/PC_AD
PC/PMMA_AD
PMMA/PC_AD
PC/PMMA_AD
PMMA/PC_AD
PC/PMMA_AD

Velocity (m/s)
1.0
1.0
2.0
2.0
3.0
3.0

Peak Force (N)
2500
2200
3500
2400
3800
3800

Time to Peak Force (ms)
4.5
2
1.7
1
1.1
0.6

Threshold Energy (J) Energy Dissipated (J)
N/A
2.5
2.1
3.9
5.8
9.8
2.1
8.8
6.8
14.3
2
12.4

top PC needed longer time to failure compared to the bi-layered structures with PC on the top
of PMMA. Compared to the post-failure pictures of the bi-layered structures with PC on the top
of PMMA with weak bonding shown in Fig. 2.5, less damage was observed on the bi-layered
structures with PMMA on top of PC shown in Fig. 2.8. The same results were observed on the
strong bonding cases, the peak force, time to crack, and threshold energy which are summarized
in Table 2.6. Thus, it can be concluded that PMMA on the top of PC can provide a better impact
resistance compared to PC on the top of PMMA.
The comparisons were also made through the bi-layered PMMA and PC structures with strong
bonding, the bi-layered PMMA and PC structures with weak bonding, and the mono-PMMA plates

29

2400
2100
2500
2600
2300

PC/PMMA_AD
PC/PMMA_NA
Mono_PMMA

2.5
2

Peak Force (N)

Time to crack (ms)

3

1.5
1
0.5
0
0

1

2
Velocity (m/s)

3

4000
3500
3000
2500
2000
1500
1000
500
0

PC/PMMA_AD
PC/PMMA_NA
Mono_PMMA
0

4

(a)

Peak Force (N)

Time to crack (ms)

PMMA/PC_AD
PMMA/PC_NA
Mono_PMMA

3
2
1
0
0

1
2
Velocity (m/s)

2
Velocity (m/s)

3

4

(b)

5
4

1

3

4

4000
3500
3000
2500
2000
1500
1000
500
0

PMMA/PC_AD
PMMA/PC_NA
Mono_PMMA

0

1
2
Velocity (m/s)

3

4

(d)

(c)

Figure
14. Comparison
of reaction
the reaction
force
mono-PMMA plates
plates to
Figure 2.13
Comparison
of the
force
of of
thethemono-PMMA
to the
thebi-layered
bi-layeredPMMA
PMMA
and
PC
structures
with
weak
bonding
and
strong
bonding
under
constant
impact
mass
and PC structures with weak bonding and strong bonding under constant impact mass.

under constant impact mass, as shown in Fig. 2.13. It was found that compared to the monolithic
PMMA plates and the strong bonded bi-layered structures, the bi-layered plates with weak bonding
took a longer time to failure, and induced lower reaction force on the specimen. It is therefore the
bi-layered PMMA and PC structures with weak bonding that showed advantage over both the
monolithic PMMA and the strong bonded bi-layered PMMA and PC structures during impact
events.
Fig. 2.14 shows the comparisons of the reaction forces, time to crack, and threshold energy
of the bi-layered PMMA/PMMA plates to the mono-PMMA plates under constant impact mass.
The results showed that the bi-layered structures led to lower reaction forces compared to the
monolithic PMMA plates. The peak force of the bi-layered structure seemed to be half that of the

30

4000

Mono_PMMA
PMMA/PMMA_AD
PMMA/PMMA_NA

2.5
2

Peak Force (N)

Time to crack (ms)

3

1.5
1
0.5

3500
3000
2500
2000
1500

Mono_PMMA
PMMA/PMMA_AD
PMMA/PMMA_NA

1000
500
0

0
0

1

2
Velocity (m/s)

3

4

0

1

2
Velocity (m/s)

3

4

(b)

(a)

Threshold Energy (J)

3.5
3
2.5
2
1.5

Mono_PMMA
PMMA/PMMA_AD
PMMA/PMMA_NA

1
0.5
0
0

1
2
Velocity (m/s)

3

4

(c)

Figure 2.14 Comparison of the reaction force of the mono-PMMA plates to the bi-layered PMMA/PMMA structures with weak bonding and
strong
bonding under constant impact mass: (a)
Impact
direction
Time to crack; (b) Peak force; (c) Threshold energy.
mono-layered PMMA plates. The bi-layered structure took a longer time to initiate cracks and
attain peak force. It is therefore concluded that layered-structures are beneficial to improve the
impact resistance of the PMMA plates. For this particular adhesive, no benefits were observed
on bi-layered PMMA/PMMA structures with strong bonding either. The average time to reach
the peak force was around 2.2 ms for the bi-layered structure, and 1.7 ms for the mono-PMMA
plate. The bi-layered structures seemed to delay the crack initiation by 29% of the mono-PMMA
Figure 16. Geometry and mesh of FE simulations

plate. Compared to the average threshold energy (2.5 J) of the monolithic PMMA plates, that (2.0
J) of the bi-layered structures (except the bi-layered structure with PMMA on the top of PC with

31

strong bonding) was lower than the monolithic PMMA plates. The lower threshold energy of the
bi-layered structures is because the threshold energy only depended on the PMMA plates since the
PC plates were not broken during the impact events. The high threshold energy of the bi-layered
structures with PMMA on the top of PC with strong bonding meant that cracking was induced not
only on the PMMA plates but also on the PC plates

Figure 13. Post-failure pictures of the bonded bi-layered structures under constant impact

Figure
2.15 Post-failure pictures of the bonded bi-layered structures under constant impact mass:
mass: the velocities are 1.0 m/s, 2.0 m/s and 3.0 m/s respectively from left to right; the first
velocities
m/s, PMMA/PMMA
2.0 m/s and 3.0
m/sthe
respectively
from
leftbonded
to right;PC/PMMA
the first plates
row is the
row isare
the 1.0
bonded
plates,
second row
is the
and PMMA/PMMA
the third row is the
bonded
plates.
bonded
plates,
the PMMA/PC
second row
is the bonded PC/PMMA plates and the third row
is the bonded PMMA/PC plates.

From previous results, we can see that although both the PMMA plates were broken when
the PMMA/PMMA bi-layered structures with weak bonding were subjected to impact velocity
larger than 1.0 m/s, as opposed to the monolithic PMMA plates, it seemed that the natural gap
32

between the two PMMA plates reduced the impact damage by some scale. The interlayer, the
natural gap, therefore, has some effects on the impact performance of the tested specimen. Thus,
to further evaluate how the soft interlayers affect the impact resistances of the bi-layered structures,
the strong bonded bi-layer structures were tested. The post-failure results are shown in Fig. 2.15.
In the mixed mode, brittle cracking indicated by the radial crack and ductile deformation indicated
by the circular crack were observed on some specimens. Cracks were induced on all the PMMA
plates except for the configuration with PMMA on top PC subjected to the impact velocity of 1
m/s. Localized delimitations were observed.

2.5

Conclusion

The low velocity impact response of the monolithic PMMA plate and a multi-layered PMMA/PC
and PMMA/PMMA plates with weak and strong interface bonding were studied by both experimental and finite element methods. The threshold impact energies of the monolithic PMMA plates
were found to be in the range of 2.0 to 3.1 J. It was found that the time to crack depended on
the impact velocity. Higher velocity took less time to initiate crack. The PMMA plate exhibits
elastic hysteresis behavior during drop weight impact test. The failure of PMMA is more sensitive
to impact velocity than impact mass. For the same impact energy, higher impact mass and lower
velocity needed a longer time to failure compared to lower impact mass and higher impact velocity.
The findings agreed well with literature [1, 2, 3, 4].
The bi-layered structures exhibited lighter deformation compared to the monolithic PMMA
plates. The layered structures seemed to delay the time to failure and create a lower reaction force.
Under the same impact conditions, the average time to failure of the bi-layer structures with the
PMMA on top of PC was 170% of the monolithic PMMA plates. When designing a bi-layered

33

PMMA and PC structure, a higher impact resistant can be obtained by putting the PMMA plate
facing the impactor because the maximum tensile stress is induced in the bottom PC layer, which
is able to undergo more deflection before cracking. For the particular Loctite E-30CL glass bonder
Epoxy adhesive, bi-layer configuration with weak bonding has an advantage over the strong bonding case during an impact event. Results of this work clearly demonstrate that layered-structures
are beneficial to improve the impact resistance of the PMMA plates and interfacial bonding plays
a role in determining the extent of the same.

34

REFERENCES

35

REFERENCES
[1] Yanxiong Liu and Benjamin Liaw. Drop-weight impact tests and finite element modeling of
cast acrylic plates. Polymer Testing, 28(6):599–611, 2009.
[2] JD Pearson, MA Zikry, M Prabhugoud, and K Peters. Measurement of low velocity and quasistatic failure modes in pmma. Polymer composites, 28(3):381–391, 2007.
[3] Michael Hebert, Carl-Ernst Rousseau, and Arun Shukla. Shock loading and drop weight impact response of glass reinforced polymer composites. Composite Structures, 84(3):199–208,
2008.
[4] Shinji Ogihara, Tomoyuki Ishigure, and Akira Kobayashi. Study on impact perforation fracture
mechanism in pmma. Journal of materials science letters, 17(8):691–692, 1998.

36

Chapter 3
Failure and Toughening Mechanism of
PMMA and PC Laminated Glassy Panel
Subjected to Impact Loading

3.1

Introduction

Transparent laminate composite is usually made of two materials with contrast material properties, ductile and brittle, to enhance the overall fracture toughness of the laminate [1, 2, 3, 4]. The
transparent glassy polymer, polymethylmethacrylate (PMMA), brittle in nature, and polycarbonate
(PC), ductile in nature, are the two typical materials used in laminate composite due to their light
weight, transparency, low cost, and machinability. PMMA dissipates energy through cracking and
usually shatters into pieces when subjected to impact loading. On the contrary, PC absorbs energy
by plastic deformation and usually is penetrated with a local deformation with four steps of failure
behavior: the elastic dishing, petalling, deep penetration, and cone cracking and plugging [5]. It
is hypothesized that the glassy laminate composite manufactured by alternatively gluing the thin
layers of PMMA and PC together has a higher energy absorption capacity than the monolithic PC
or PMMA when subjected to impact loading. Although this phenomenon has been experimentally
observed and proved [6, 7, 8], the toughening mechanism including the interaction between im-

37

pact resistance and impact velocities, geometry and aspect ratio of projectile, and thickness and
mechanical properties of each layers has not been yet fully understood.
As amorphous thermoplastic, PMMA and PC are pressure, loading rate, and temperature dependent. The pressure dependent behavior is reflected by the different stress-strain response under
tension and compression [9, 10, 11, 12, 13, 14]. The deformation of amorphous polymer is governed by the state of orientation of the chain network of the microstructure. During the process
of deformation, the chain network changes from random orientation to a highly orientated state,
which is the preferred state in strengthening the polymer [15, 9, 12]. The pressure dependent
behavior of amorphous polymer is the result of different orientation processes in tension and compression: planar orientation process in compression and uniaxial orientation process in tension
[9] respectively. PMMA and PC are highly strain rate dependent, exhibiting an increase in yield
strength with an increase in strain rate [16].
The impact behavior of PMMA and PC subject to various strain rates has been studied extensively through use of various impact testing methods including low velocity drop-weight, impact
tests [17, 18, 19, 20, 21], high velocity impact by use of an air gun [2, 22, 23, 4, 24, 25], and shock
wave impact by use of explosives [26, 27]. It was found that multiple radial cracks are usually
induced in the PMMA plate subjected to low to intermediate impact load, and the energy release
rate of PMMA is found to be in the range of 252-728 J/m [19]. The thickness or volume fraction
of the PMMA plate was observed to play an important role in improving the ballistic performance
of multilayered composites subject to high velocity impact ranging from 300 to 550 m/s [1, 2].
The failure of a PMMA and PC laminate plate with a higher volume fraction of PMMA is dominated by brittle cracking. On the contrary, the failure is dominated by ductile deformation. The
experimental observations showed that the impact resistance of PMMA and PC laminate depends
on many variables. To avoid a case-by-case or trail and error study, an optimal design procedure
38

for maximizing the impact-resistance of glassy panels should be conducted.
The challenge of the optimization methods for the PMMA/PC layered plates is the difficulty
in configuring the analytical solutions for dynamic problems due to the unstable and instantaneous
nature of brittle cracking subjected to impact loading. The classical laminate theory, wave theory,
and energy method are the three major approaches in constructing the analytical solutions for
optimizing the impact resistance of laminate structures. Although the classical laminate theory is
widely used to map the stress states of laminated plate subjected to stationary loading, solving the
layered plates subjected to impact loading is still challenging. Attempts have to be made to refine
the classical theory so that it can be used to investigate the strain-stress state of multilayer plates
under impulse and impact loading [28, 29, 30]. Smetankina [30]used a refined laminate theory on
the basis of the hypothesis of the refined theory of Timoshenko to investigate the strain-stress state
of multilayer plates under impulse and impact loading. In his work, the transverse shear was taken
into account and the displacement and the external load are expanded into Fourier series functions
satisfying the problem boundary conditions. The numerical results had a good agreement with
the experimental results. Shupikov [29] used Smetankinas refined laminate theory to investigate
the optimal design of multilayer plates. In his work, attempts were made to minimize the stress
in layers and the total mass of the multilayered plates. The approach of wave propagation in
layered heterogeneous systems was to convert the initial velocity boundary value problem to a time
dependent stress boundary value problem, and the stress time history was obtained by means of
superposition [31, 19, 32]. Luo [19] investigated the stress wave propagation for a layered structure
under impulsive loads. It was shown that the reduced amplitude and elongated pulse duration could
be obtained with proper selection of materials and dimensions of layers when an incident stress
pulse passed through a layered structure. The optimal design of layered structures could also be
implemented based on the concept energy balance [33, 34] during an impact event. Ray [34]
39

optimized the energy dissipation characteristics of active constrained layer damping treatments of
plates using a rational design approach. In his work, the governing equation was given based on
the system energy balance, which includes the potential energy, kinetic energy, and the work done
by the piezo-electric control force, and the work dissipated in the viscoelastic core.
The previous literature shows that there has been extensive research focused on the experimental and FEM study of laminated composite material, as well as optimal design of laminated
composite material based on analytical analysis in the past decades. However, the unique characteristic and failure mechanism of PMMA and PC subjected to impact loading has yet not been fully
understood. Thus, in this research of work, a material model for amorphous materials with large
deformation was first implemented in ABAQUS/EXPLICIT based on the existing Boyce and Arruda’s [9] physical model to predict the pressure and strain rate dependent stress-strain relationship
of PC. Second, the impact behavior of monolithic PMMA and PC and layered PMMA/PC laminate
at low and high impact loading was numerically investigated by use of ABAQUS/EXPLICIT and
the failure mechanism, crack pattern, and energy dissipation were transversely compared. Lastly,
to further evaluate material characteristics such as strain-rate sensitivity and pressure dependent
behavior of the monolithic components on the ballistic performance of these PMMA/PC laminated systems, an optimal study by incorporating MATLAB optimization tool box and ABAQUS
FEM simulations were conducted. The role of aspect ratio, yielding stress of PC, and failure
stress of PMMA on the impact behavior were studied. This work serves as a design tool for effectively choosing materials and structures for optimizing ballistic performance of ceramic-plastic
laminates-based armor.

40

2. Methods

3.2

Methods
As amorphous polymer, both PMMA and PC are strain rate, temperature, and pressure dependent

As amorphous polymers, both PMMA and PC are strain rate, temperature, and pressure dependent
materials. A typical stress-strain curve of PMMA and PC under compression is shown in Figure 1.The

materials. A typical stress-strain curve of PMMA and PC under compression is shown in Fig.
deformation includes three parts: linear elastic deformation, strain softening, and strain hardening.

3.1. The deformation includes three parts: linear elastic deformation, strain softening, and strain
From the view of microscopic structure of PMMA and PC, the linear elastic region is where the

hardening. From the point view of microscopic structure of PMMA and PC, the linear elastic
applied force has not passed the intermolecular interactive barriers and the molecular structures rotate
region
is where the applied force has not passed the intermolecular interactive barriers and the
or slide structures
along each other.
After breaking
the intermolecular
segments changed
to new
molecular
are rotating
or sliding
along eachforce,
other.molecular
After breaking
the intermolecular

force,
molecular
to new positions
causes
thetostrain
softening. ofThe
positions
whichsegments
caused thechange
strain softening.
The strainwhich
hardening
is due
the contribution
thestrain
hardening
is due
to the which
contribution
the segments
reordering,
which tends
to orientmolecular
parallel with
segments
reordering
tends tooforient
parallel with
loading direction
to maximize
loading direction to maximize molecular force.
force.

Strain
hardening
Yielding

Strain
softening

Figure 1 Schematic representation of amorphous polymers under compression.

Figure 3.1 Schematic representation of amorphous polymers under compression.
The default ABAQUS model fails to describe the strain rate and pressure dependent behavior. Thus,

In the past decades, many constitutive equations have been proposed to predict the lineara user subroutine to implement the constitutive behavior of PC is needed. Many constitutive equations

softening-hardening behavior of amorphous glassy polymers. Among them, Arruda and Boyce’s
have been proposed to predict the linear-softening-hardening behavior of amorphous glassy polymers.

three-dimensional physical-based model [9] is one of the most adopted. Boyce’s model was de8

veloped based on Harvard and Thackray [35] and Argon’s physical model [15], which divided

41

resistance inside the amorphous materials into two parts: one is the molecular resistance; another
is the entropic resistance. The entropic resistance contributes the strain hardening of the materials,
and the molecular resistance accounts for the strain softening. Thus, a material subroutine based on
Arruda and Boyces constitutive model for amorphous polymers was developed and implemented
in ABAQUS to predict the strain rate and pressure dependent behavior of PC subjected to uniaxial
loading at different strain rate. The elastic material properties used in the constitutive model are
shown in Table 3.1.
Table 3.1 Material properties for PC and PMMA.

Material

Density
( kg/mm3 )

Modulus of
Elasticity ( MPa )

Poisson's
Ratio

PC

1.2x10-6

2300

0.3

PMMA

1.19x10-6

3100

0.33

The PMMA was simulated as a brittle material using the brittle cracking material model of
ABAQUS. When elements are under compression, the stress-strain response is linear. The failure
of an element is caused by tensile failure. Thus, when the maximum stress at the element is larger
than the post yielding stress, the crack will initiate until failure of the element. Failed elements
will be removed from computation. The elastic material properties of PMMA are shown in Table
3.1.
The geometry and mesh of the FEM analysis are shown in Fig. 3.2. The overall thickness
of the laminate panel simulated is 6.4 mm. The projectile was modeled as a rigid body, and an
8-node C3D8 linear brick element was used for modeling the specimen. Two types of projectile
were simulated. One is spherical projectile with a radius of 6.35 mm (Type 2), which was used for
the low velocity impact test. The other was the blunt projectile with a radius of 2.73 mm (Type
1) which was used for the high velocity impact test. The peripheries of the laminate plates were
42

Impact force, energy, and stress and strain contour were obtained. The various configurations of the
laminate structure simulated are summarized in Table 2.

linear brick
element was used for modeling the specimen. Two types of tup were simulated. One is
Impact

Direction
spherical tup with a radius of 6.35 mm (type 2) which is for the low velocity impact test. Another is
the blunt projectile with a radius of 2.73 mm (type 1) which is for the high velocity impact test. The

A impact direction. The mesh density of the specimen
A
peripheries of the laminate plates are fixed in the
Type 1 Projectile
was chosen based on the low cost of computation without compromising
the accuracy of the results.
Impact force, energy, and stress and strain contour were obtained. The various configurations of the
laminate structure simulated are summarized in Table 2.

Impact
Type 2 Projectile
Direction

A-A

Figure
Geometry and
and mesh
Figure
3.22Geometry
meshofofFE
FEsimulations
simulations.
Table 2 Various configurations of the laminate plate panel simulated
Impact

Mass of A

Radius of

Numbers

A Type

Sequences

fixed in the impact direction. The mesh density of the specimen was chosenbased
on the low cost
PC
of
of
Projectile
of
Case Velocity Projectile
 PMMA
Type 1 Projectile
Plates
(Kg)
(m)
(m/s)
Interface
No.
of computation
without compromising
the accuracyLayer
of the results.
Impact force, energy,
and stress
1

5

0.00273

7.1

8

MCMCMCMC

0.96

bounded

2
5
8
MCMCMCMC
0.96 structure
bounded
and strain
contour
were0.00635
obtained.7.1The various configurations
of the laminate
simulated
3

300

0.00273

0.0032

8

MCMCMCMC

0.36

bounded

are summarized
3.2. The0.0032
σ pc and σPMMA
stress 0.36
of PC and
failure stress
4
300 in Table
0.005
8 are the yielding
MCMCMCMC
bounded
of PMMA respectively.Type 2 Projectile

A-A
10

Figure
2 Geometry and
of FE simulations
Table 3.2 Various
configurations
ofmesh
the laminate
plate panel simulated.
Table 2 Various configurations of the laminate plate panel simulated
Numbers

Sequences

Impact

Radius of

Mass of

Case

Velocity

Projectile

Projectile

of

of

No.

(m/s)

(m)

(Kg)

Layer

Plates

 PC
 PMMA

Type
of
Interface

1

5

0.00273

7.1

8

MCMCMCMC

0.96

bounded

2

5

0.00635

7.1

8

MCMCMCMC

0.96

bounded

3

300

0.00273

0.0032

8

MCMCMCMC

0.36

bounded

4

300

0.005

0.0032

8

MCMCMCMC

0.36

bounded

10

The impact resistance depends on thickness of each layer, material properties of each layer,
sequences of layers, aspect ratio of plates to projectile, and impact velocities, as well as the mass

43

of projectile. Thus, optimization is conducted based on some of the above design variables. The
objective The
function
to maximize
energy absorption,
whichproperties
was equivalent
the minimum
impactwas
resistance
depends the
on thickness
of layers, material
of layers,tosequences
of
residual
velocity.
during
animpact
impact
event isasshown
layers,
aspect The
ratioenergy
of platesbalance
to projectile,
and
velocities,
well asasmass of projectile. Thus,
optimization is conducted based on the some of above design variables. The objective function is to

1
Edissipated = Kinitial − Kresidual = m(v2initial − v2residual )
(3.1)
maximize the energy absorption which is equivalent to2the minimum residual velocity. The energy
wherebalance
Kdissipated
, the
, initial
during
an dissipated
impact eventenergy;K
is shown initial
in Equation
(1).kinetic energy; Kresidual , residual energy; m,
1
the mass of projectile; vinitial , the initial
2velocity
2 of the projectile; and vresidual , the residual velocity
(1)
Edissipated  K initial  K residual  m(vinitial
 vresidual
)
2
of the projectile.
The optimal study was performed by incorporating MATLAB and FEM simulations by ABAQUS.

The optimal study was performed by incorporating MATLAB and FEM simulations by ABAQUS.
The MATLAB
fmincon
functionwas
wasused
used as
as the
the searching
to find
the minimum
of objective
The MATLAB
fmincon
function
searchingmethod
method
to find
the minimum
of objecfunction. The
The flow
of of
thethe
optimization
study study
is shown
in Figurein3.Fig. 3.3.
tive function.
flowchart
chart
optimization
is shown

Go back to the first step
Repeat the iteration

Run MATLAB fmincon script to pass design variables, Xi, to
ABAQUS/Explicit

Run ABAQUS/Explicit model with the assigned
design variables, Xi

Access to ABAQUS output database to get the
objective function: f =-energy dissipation

Check constraints of Xi and gradient of f
to judge if satisfied
No
f*=f, X*=Xi
Yes
Problem Solved
Figure 3 Schematic representation of the optimization process of constrained optimization
Figure 3.3 Schematic representation of the optimization process
problem.

44

11

3.33. Results
Results
3.1. Verification of material model for impact testing simulation

3.3.1

Verification of Material Model for Impact Testing Simulation

The material model for PC was tested in a single element subjected to uniaxial compressive load at

The material model for PC was tested in a single element subjected to uniaxial compressive load
low and high strain rates. The stress-strain relationships from the simulations were compared to the

at low and high strain rates. The stress-strain relationships from the simulations were compared to
experimental results in the literature (Mulliken and Boyce, 2006).The plots are shown in Figure 4.

the experimental results in the literature [16]. The plots are shown in Fig. 3.4. Due to the lack of
Due to lack of experimental data for the strain rate of -5/s, the experimental result at strain rate of -1/s

experimental data for the strain rate of -5 s−1 , the experimental result at the strain rate of -1 s−1 was
was used as reference data to FEM result at -5/s. The comparisons show that the material model could

used as reference data for the FEM result at -5 s−1 . The comparisons show that the material model
predict
the pressure
and strain
dependent
behavior ofbehavior
PC with acceptable
accuracy.
could
predict
the pressure
andrate
strain
rate dependent
of PC with
acceptable accuracy.
180

True Stress (MPa)

160
140
120
100
80
60

Exp 5/s
FEM 5/s
EXP 5050/s
FEM 5050/s

40
20
0
0

0.2

0.4

0.6

0.8

1

True Strain

4 Simulations and experimental curves of the uniaxial compression responses of PC at
FigureFigure
3.4 Comparisons
of the numerically predicted uniaxial stress-strain curves (compression)
low and high strain rates.
of PC to the experimental measurements.
The stress-strain relationship predicted by the material mode at different strain rates will be used as

The
stress-strain relationship predicted by the material mode at different strain rates was used
the input of PC in the impact testing simulations. The simulation results of monolithic PC and PMMA
as the material model input for PC in the impact testing simulations. The simulation results of
subjected to low velocities are compared to the ones from the low velocity impact testing (W. Zhang,

monolithic PC and PMMA subjected to low velocities are compared to the ones from the low
2011) (Figure 5). The FEM simulations of low velocity impact testing have a good agreement with the

velocity impact testing [36] (Fig. 3.5). The FEM simulations of low velocity impact testing have
experimental results. There are no force and energy curves for high velocity impact testing available in

a good agreement with the experimental results. There are no force and energy curves for high
the existing literatures. Sarva (Sarva et al., 2006b) measured the fraction of residual and initial kinetic

velocity impact testing available in the existing literature. Sarva [24] measured the fraction of
residual and initial kinetic energy of monolithic PC plate and 6-layered laminate PC and PMMA
12

45

energy of monolithic PC plate and 6-layered laminate PC and PMMA plate with a dimension of
100x100mm and thickness of 6.35 mm subjected to impact velocities of 300-500 m/s. Thus, the FEM

plate with a dimension of 100x100 mm and thickness of 6.35 mm subjected to impact velocities
simulations refereed to the exrpiemtal results from Sarva as baseline results. The FEM simulation

of 300-500 m/s. Thus, the FEM simulations referred to the experimental results from Sarva as
results
and discussion
are simulation
showing in section
FEMdiscussion
of high velocity
impact testing.
baseline
results.
The FEM
resultsofand
are shown
in section 3.3.3.
b)

3500
3000
2500
2000
1500
1000
500
0

EXP/PMMA
FEM/PMMA

10
8

2.0 m/s

6

EXP/PMMA
FEM/PMMA

4

2.0 m/s

2

0

c)

5

10
Time(ms)

15

0

20

0

d)

6000

EXP/PC
FEM/PC

5000

5

10
Time (ms)

2000

25

10
EXP/PC
FEM/PC

5

1000

2.5 m/s

0

0
5 Time(ms) 10

20

15

3000

0

15

20

2.5 m/s

4000

E(J)

Force (N)

12

E(J)

Force (N)

a) 4000

15

0

5

10

15

Time (ms)

Figure 5 Comparisons of FEM low velocity impact simulation results to experimental results: a)

Figure 3.5 Comparisons of FEM low velocity impact simulation results to experimental results: a)
Reaction force vs. time for monolithic PMMA subjected to impact velocity of 2.0 m/s; b) Impact
Reaction
force vs. time for monolithic PMMA subjected to impact velocity of 2.0 m/s; b) Impact
energy vs. time for monolithic PMMA subjected to impact velocity of 2.0 m/s; c) Reaction force
energy vs. time for monolithic PMMA subjected to impact velocity of 2.0 m/s; c) Reaction force
vs. time
monolithic
subjectedtotoimpact
impact velocity
velocity ofof2.5
d) d)
Impact
energy
vs. time
for for
vs. time
for for
monolithic
PCPC
subjected
2.5m/s;
m/s;
Impact
energy
vs. time
monolithic
subjected
impactvelocity
velocity of
of 2.5
monolithic
PC PC
subjected
to to
impact
2.5 m/s.
m/s.
3.2. FEM of low velocity impact test
The kinetic energy history and various types of energy dissipation of layered PMMA and PC

3.3.2

FEM of Low Velocity Impact Test

laminate panel and monolithic PC subjected to impact velocity of 5m/s by a spherical projectile with

The kinetic energy history and various types of energy dissipation of layered PMMA and PC
radius of 6.35mm are shown in Figure 6. The results showed that monoclinic PC is beneficial to the

laminate panel and monolithic PC subjected to impact velocity of 5m/s by a spherical projectile
with radius of 6.35 mm are shown in Fig. 3.6. The results showed that monolithic PC was beneficial
13

to the low velocity impact testing with an aspect ratio (a ratio of the diameter of projectile to the
diameter of the plate) of 0.167. The energy dissipated by a monolithic PC was about 1.58 times
46

low velocity impact testing with an aspect ratio (a ratio of the diameter of projectile to the diameter of
the plate) of 0.167. The energy dissipated by a monolithic PC is about 1.58 times higher than
8-layered PMM/PC laminate (Figure 6). The impact event of monolithic PC is dominated by plastic

higher than 8-layered PMM/PC laminate (Fig. 3.6). The impact event of monolithic PC was
dissipation. However, strain energy is the major form of energy dissipation for brittle PMMA. The

dominated by plastic dissipation. However, strain energy was the major form of energy dissipation
results showed that although the strain energy dissipated by the layered laminate increased from 0.4 J

for brittle PMMA. The results showed that although the strain energy dissipated by the layered
to 5.4
J compared
monolithic
the total energy
dissipated actually
becausedissipated
of the
laminate
increased
fromto0.4
J to 5.4PC,
J compared
to monolithic
PC, thedecreased
total energy
energy because
dissipatedthe
deceased
42.1 Jdissipated
to 19.9 J. decreased from 42.1 J to 19.9 J.
actuallyplastic
decreased
plasticfrom
energy
a)

b)
Energy Dissipated (J)

Kentic Energy (J)

110
2L
4L
8L
Mono_PC

90
70
50
30

50
45
40
35
30
25
20
15
10
5
0

ALLDMD
ALLPD
ALLFD
ALLSE

8

0

2

Time (ms)

4

6

4

2

1

Laminate Configurations

Figure Figure
3.6 FEM
results
of layered
PMMA/PC
laminate
panel
andand
monolithic
PCPCsubjected
6 FEM
results
of layered
PMMA/PC
laminate
panel
monolithic
subjectedtotoimpact velocity
5 m/s of
by5m/s
a spherical
projectile
with radius
of 6.35
a) Kinetic
energy
vs.VS.
time;
impact of
velocity
by a spherical
projectile
with radius
of mm:
6.35mm:
a) Kinetic
energy
b) Distribution
of
various
energies
dissipated
calculated
by
ABAQUS/EXPLICIT.
Note:
ALLSE,
time; b) Distribution of various energies dissipated calculated by ABAQUS/Explicit. Note:
strain energy;
fictionALLFD,
dissipation;
ALLPD,
plastic
dissipation;
dissiALLSE,ALLFD,
strain energy;
fiction
dissipation;
ALLPD,
plasticALLDMD,
dissipation;damage
ALLDMD,
pation. damage dissipation.
The stress contours of half plates of the 8-layered PMMA/PC laminate, monolithic PC, and

The stress contours of half plates of the 8-layered PMMA/PC laminate, monolithic PC, and
monolithic PMMA at different time frames are shown in Figure 7. It was shown that crack of

monolithic PMMA at different time frames are shown in Fig. 3.7. It was shown that the crack in
monolithic PMMA first initiated on the back plate and then prorogated toward to radial direction

the monolithic PMMA first initiated on the back of the plate and then propagated radially. In total,
of the plate. Totally 4 radial cracks were formed and this has a good agreement to the

4 radial cracks were formed and this has a good agreement with the experimental observations [36].
experimental observations(W. Zhang, 2011). As expected, the deformation of monolithic PC was
As expected,
the deformation of monolithic PC was dominated by local plastic deformation. The
dominated
by local
plasticmitigated
deformation.
8-layered
PMMA/PC
laminate
locked
8-layered
PMMA/PC
laminate
the The
radial
cracking
in the PMMA
layers,
butthethebrittle
plastic

deformation in the PC layers was not magnified.

47

14

cracking in a local area under the projectile, but the plastic deformation is still constrained in the
local area.
Mono_PMMA

Mono_PC

8L

0.8ms

1.6ms

3.2ms

7 Stress
contoursofof8-layered
8-layered PMMA/PC
PMMA/PC laminate,
PC,PC,
andand
monolithic
FigureFigure
3.7 Stress
contours
laminate,monolithic
monolithic
monolithic
PMMA
subjected
to impact
velocityofof55m/s
projectile
withwith
radius
of 6.35mm
PMMA
subjected
to impact
velocity
m/sby
bya aspherical
spherical
projectile
radius
of 6.35at mm
different
time
frames.
at different
time
frames.
3.3. FEM of high velocity impact test

3.3.3

FEM of High Velocity Impact Test
The kinetic energy history and various types of energy dissipation of layered PMMA and PC

The kinetic energy history and various types of energy dissipation of layered PMMA and PC lamlaminate panel and monolithic PC subjected to impact velocity of 300 m/s by a blunt projectile with a

inate panel and monolithic PC subjected to impact velocity of 300 m/s by a blunt projectile with
radius of 2.73 mm are shown in Figure 8. The results showed that the energy dissipated by a 8-layered

a radius of 2.73 mm are shown in Fig. 3.8. The results showed that the energy dissipated by a
PMMA/PC laminate increased 5 percent than the monolithic PC. The amount of energy dissipated in

8-layered PMMA/PC laminate increased 5 percent more than the monolithic PC. The amount of
different
forms in
including
strain
energy,
plastic dissipation,
frictionplastic
dissipation,
and damage
dissipation
energy
dissipated
different
forms
including
strain energy,
dissipation,
friction
dissipacalculated
by ABAQUS/EXPLICIT
are by
shown
in Figure 8b. It can be
that although
the volume
tion, and
damage
dissipation calculated
ABAQUS/EXPLICIT
areseen
shown
in Fig. 3.8b.
It can be

seen that although the volume fraction of PC deceased from 1 to 0.5 in 8-layered PMMA/PC lami15

nate compared to the monolithic PC, the plastic dissipation only deceased about 24 percent. These
results indicated that the PMMA layer enhanced the PC layer by amplifying the local plastic defor48

fraction of PC deceased from 1 to 0.5 in 8-layered PMMA/PC laminate compared to the monolithic
PC, the plastic dissipation only deceased about 24 percent. These results indicate that the PMMA layer
enhanced the PC layer by amplifying the local plastic deformation of PC. Compared to monolithic PC,

mation of PC. Compared to monolithic PC, the strain energy, friction dissipation in the 8-layered
the strain energy, friction dissipation in the 8-layered PMMA/PC laminate increased from 1.7 J to 10.4

PMMA/PC laminate increased from 1.7 J to 10.4 J and 3.9 J to 5.6 J respectively. Compared to the
J and 3.9 J to 5.6 J respectively. Compared to the experimental results shown by Sarva (Sarva et al.,

experimental results shown by Sarva [24], although the total energy dissipation predicted by the
2006b), although the total energy dissipation predicted by the simulations do not match (lower than

simulations did not match (lower than experimental measurements), the computational trend that
experimentallaminate
measurements)
are lower,
the computational
that layered
PMMA/PC
laminate
layered PMMA/PC
improving
the energy
dissipationtrend
appeared
to match
the experimental
improving the energy dissipation appeared to follow the experimental trend.

Kentic Energy (J)

a)160

b)
2L
4L
8L
Mono_PC

140

Energy Dissipated (J)

trend.

120

100

80
0

0.01

0.02
0.03
Time (ms)

45
40
35
30
25
20
15
10
5
0

ALLDMD
ALLPD
ALLFD
ALLSE

8

0.04

4

2

1

Laminate Configurations

Figure 3.8 FEM results of layered PMMA/PC laminate panel and monolithic PC subjected to
Figure 8 FEM results of layered PMMA/PC laminate panel and monolithic PC subjected to
impact velocity of 300 m/s by a blunt projectile with radius of 2.73 mm : a) Kinetic energy vs.
subjected to impact velocity of 300m/s by a blunt projectile with radius of 2.73 mm : a) Kinetic
time; b) Distribution of various energies dissipated calculated by ABAQUS/Explicit.
energy VS. time; b) Distribution of various energies dissipated calculated by ABAQUS/Explicit.

The stress contours of half plates of the 8-layered PMMA/PC laminate, monolithic PC, and

The stress contours of the half plates of the 8-layered PMMA/PC laminate, monolithic PC,
monolithic PMMA at different time frames are shown in Figure 9. The simulation contour showed the

and monolithic PMMA at different time frames are shown in Fig. 3.9. The simulation contour
failure of monolithic PMMA is dominated by brittle cracking. Many spalls were observed on the back

showed that the failure of monolithic PMMA was dominated by brittle cracking. Many spalls were
surface of the immediate impact region. The failure mechanism of PMMA has a good agreement with

observed on the back surface of the immediate impact region. The failure mechanism of PMMA
The deformation
monolithic
the agreement
experimentalwith
observation
shown by Sarva
(Sarva et al.,
2006b).
has a good
the experimental
observation
shown
by Sarva
[24]. Theofdeformation
PC was
by localby
plastic
which occurred
in occurred
the immediate
zone
of monolithic
PCdominated
was dominated
localdeformation
plastic deformation,
which
in theimpact
immediate

impact zone under the projectile. Stress contour of the 8-layered PMMA/PC laminate showed that
16

the laminate structures magnified the plastic deformation zone and decreased the brittle damage
area.

49

under the projectile. The 8-layered PMMA/PC laminate magnified the plastic deformation zone
and decreased the brittle damage area.
Mono_PMMA

Mono_PC

8L

0.006ms

0.02ms

0.04ms

Figure 9 Stress contours of 8-layered PMMA/PC laminate, monolithic PC, and monolithic
Figure
3.9 Stress contours of 8-layered PMMA/PC laminate, monolithic PC, and monolithic
PMMA
subjected
impactvelocity
velocity of
of 300
300m/s
projectile
withwith
radius
of 2.73
at
PMMA
subjected
to to
an an
impact
m/sbybya ablunt
blunt
projectile
radius
of 2.73
at
different
time
frames.
different time frames.
3.4. Optimal design of laminate panel subjected to impact test

3.3.4

Optimal Design of Laminate Panel Subjected to Impact Test
The study in the previous two sections showed that the 8-layer laminate PMMA and PC panel

The study in the previous two sections showed that the 8-layer laminate PMMA and PC panel
subjected to high velocity impact through a smaller projectile has more benefits than the low velocity

subjected to high velocity impact with a smaller projectile has more benefits than the low velocity
impact with a larger projectile. The fraction of kinetic energy vs aspect ratio of the diameter of the

impact with a larger projectile. Thus, to evaluate the effect of aspect ratios on the impact behavior
plates and diameter of the projectile at the impact speeds of 5 m/s and 300 m/s of the monolithic PC

of the PMMA and PC laminate more simulations with different aspect ratios were conducted. The
and 8-layer laminate are plotted Figure 9. For both of the monolithic PC and PMMA, the fraction of

fraction of kinetic energy versus various aspect ratios of the monolithic PC and 8-layer laminate
kinetic energy decreased with aspect ratio. When the aspect ratio decreased to about 0.07, both the
subjected
to the impact velocities of 5 m/s and 300 m/s are plotted in Fig. 3.10. For both of the

monolithic PC and PMMA, the fraction of kinetic energy decreased with aspect ratio. When the
17

aspect ratio decreased to about 0.07, both the high and low velocity showed that the 8-layered
laminate panel started to dissipate more energy than the monolithic PC.

50

high and low velocity showed that the 8-layered laminate panel started to dissipate more energy than
the monolithic PC.

Fraction-Kentic Energy

1.0
0.8
0.6
0.4

Mono_PC_5 m/s
8L_5 m/s
Mono_PC_300 m/s
8L_300m/s

0.2
0.0
0

0.05

0.1
Aspect Ratio

0.15

0.2

Figure 10 Fraction of kinetic energy of 8 layers of PMMA and PC laminate panel subjected to

Figure 3.10 Fraction of kinetic energy of 8 layers of PMMA and PC laminate panel subjected to
impact velocity of 300m/s by a blunt projectile with radius of 1.73–5 mm and impact velocity of
impact
velocity of 300 m/s by a blunt projectile with radius of 1.73-5 mm and impact velocity of 5
5 m/s by a spherical projectile with radius of 6.35-1.73 mm.
m/s by
a spherical projectile with radius of 6.35-1.73 mm.
Non-linear optimizations of energy dissipation capacity of 8-layered PMMA/PC laminate were

Non-linear optimization of the energy dissipation capacity of 8-layered PMMA/PC laminate
conducted by use of MatLab optimization tool box, fmincon. The residual kinetic energy of the

were conducted by use of MATLAB optimization tool box, fmincon. The residual kinetic energy
laminate system was minimized in terms of two variables, respectively: one is the multiplier of

of the laminate system was minimized in terms of two variables: one is the multiplier of postpost-yielding stresses of PC at impact speed of 300 m/s within the range of 0.5 to 2.0 (case 1); another

yielding stresses of PC at impact speed of 300 m/s within the range of 0.5 to 2.0 (Case 1); another
is multiplier of stress at failure of PMMA at impact speed of 300 m/s within the range of 0.6 to 1.4

is multiplier of stress at failure of PMMA at impact speed of 300 m/s within the range of 0.6 to 1.4
iterationprocesses
processes of
casescases
are shown
in Figure
11. For
bothFor
of the
(Case(case
2). 2).
TheThe
iteration
ofthe
thetwo
twooptimization
optimization
are shown
in Fig.
3.11.
both of
cases, the
energy
has ahas
negative
correlation
with thewith
multipliers.
For case 1, For
it took
6 1, it
the cases,
the residual
residualkinetic
kinetic
energy
a negative
correlation
the multipliers.
case

took iterations
6 iterations
to seek
the minimum
kineticwith
energy
with
value
J at aofmultiplier
to seek
the minimum
residualresidual
kinetic energy
a value
of a84.6
J atof
a 84.6
multiplier
1.99
of 1.99,
which
to ayielding
post yielding
of PC
MPa.
2, it took
which
yieldsyields
to a post
stress of stress
PC 182.5
MPa.182.5
For case
2, itFor
tookCase
16 iterations
to 16
seekiterations
the
to seek the minimum residual kinetic energy with a value of 98.2 J at a multiplier of 1.38, which
minimum residual kinetic energy with a value of 98.2 J at a multiplier of 1.38 which yields to a failure

yields to a failure stress of PMMA 351 MPa.
stress of PMMA 351 MPa.

To evaluate the joint effects of the above two variables on the energy dissipation capacity of
18
8-layered PMMA/PC laminate, an optimization was also performed by use of both two multipliers

as variables to minimize the residual kinetic energy of the laminate system subjected to impact
velocity of 300 m/s. The corresponding optimization iteration process is shown in Fig. 3.12. The
51

70

Residual En

Residual E

80

Min. Point

60
50

96

Min. Point

94
92
90

0

1

2

3
4
Iterations

5

6

0 1 2 3 4 5 6 7 8 9 1011121314151617
Iterations

7

Figure 11 Non-linear optimization of residual kinetic energy of 8-layered PMMA/PC laminate
subjected to impact velocity of 300m/s by a blunt projectile with radius of 2.73 mm: a) multiplier
of yielding stress of PC as variable; b) multiplier of stress at failure of PMMA as variable.
To evaluate the joint effects of the above two variables on the energy dissipation
capacity of
Min. Point
Min. Point
8-layered PMMA/PC laminate, an optimization was also performed by use of both two multipliers as
variables to minimize the residual kinetic energy of the laminate system subjected to impact velocity
Figure 11 Non-linear optimization of residual kinetic energy of 8-layered PMMA/PC laminate

Figure
3.11
Non-linear
optimization
of residual process
kineticisenergy
ofFigure
8-layered
PMMA/PC
laminate
of 300subjected
m/s. The
corresponding
optimization
shown
12.2.73
The
optimization
to
impact velocity
of 300m/siteration
by a blunt projectile
withinradius
of
mm:
a) multiplier
subjected to impact velocity of 300 m/s by a blunt projectile with radius of 2.73 mm: a) multiplier
of yielding stress of PC as variable; b) multiplier of stress at failure of PMMA as variable.
ofshowed
yielding
of PC asfunction
variable;
multiplier
of stress
at Jfailure
of PMMA
asmultipliers
variable.
thatstress
the objective
wasb)oscillating
at value
of 84
with both
of the two
To evaluate the joint effects of the above two variables on the energy dissipation capacity of

oscillating at showed
their upper
bound
at the end of
the iteration.
Considering at
theathree
cases,
it seemed
optimization
that
the objective
function
was oscillating
value
of 84
J with that
both of the
8-layered PMMA/PC laminate, an optimization was also performed by use of both two multipliers as

two
oscillating
at their
upper
bound
at therate
end
of the iteration.
themultipliers
objective function
is relatively
more
sensitive
to strain
hardening
of PC than PMMA.
variables to minimize the residual kinetic energy of the laminate system subjected to impact velocity
105
Residual Energy (J)

100 optimization iteration process is shown in Figure 12. The optimization
of 300 m/s. The corresponding
95
90

showed that the objective function was oscillating at value of 84 J with both of the two multipliers
85
80

oscillating at their upper bound
75 at the end of the iteration. Considering the three cases, it seemed that
70
0
4
8 12 16 20 24 28 32
the objective function is relatively
more
sensitive
to strain
rate hardening of PC than PMMA.
Iteration
s

Figure 12 Non-linear optimization of residual kinetic energy of 8-layered PMMA/PC laminate

Figure 3.12 Non-linear optimization of residual kinetic energy of 8-layered PMMA/PC laminate
subjected to impact velocity of 300m/s by a blunt projectile with radius of 2.73 mm in terms of
subjected
to impact velocity of 300 m/s by a blunt projectile with radius of 2.73 mm in terms of
two
variables:
multiplier of
of PC
andand
multiplier
of stress
at failure
of PMMA.
two variables: multiplier
ofyielding
yieldingstress
stress
of PC
multiplier
of stress
at failure
of PMMA.
4. Discussion
At low strain rate (<1s-1), the stress-strain relationships of PMMA and PC follow a similar trend

3.4

Discussion

19

Figure 12 Non-linear optimization of residual kinetic energy of 8-layered PMMA/PC laminate
subjected
to impact
velocity
of 300m/s by relationships
a blunt projectile
with radius
mm ainsimilar
terms of
At a low
strain rate
(≤1 s−1
), the stress-strain
of PMMA
and of
PC2.73
follow
trend
two variables: multiplier of yielding stress of PC and multiplier of stress at failure of PMMA.

(Figure 1): linear elastic, strain rate dependent yielding, strain softening, and strain hardening.
4. Discussion

With the increase of strain rate, PMMA starts to show its brittle nature with a different stress-strain
At low strain rate (<1s-1), the stress-strain relationships of PMMA and PC follow a similar trend

response: linearly increasing to a high peak stress, a steep softening, and failure. Due to these
19

52

different material characteristics, PMMA and PC exhibit different failure behavior under impact
loading: the failure of PMMA is governed by brittle cracking and the failure of PC is governed
by plastic deformation. Thus, plastic and brittle material models were used to simulate the impact
performance of PMMA and PC respectively by ABAQUS/EXPLICIT. In conjunction with damage
initiation and evolution material failure model, the rate and pressure dependent stress-strain curves
of PC predicted by Arruda and Boyces three-dimensional physical-based models [9] was used as
the material inputs of ABAQUS/EXPLICIT to predict the impact behavior of PC plates subjected
to low (5 m/s) and high speed (300 m/s) impact loading. The material model of brittle cracking
of ABAQUS/EXPLICIT was used to simulate the brittle failure behavior of PMMA subjected to
impact loading.
It was found that at the impact conditions of the aspect ratio of 0.17, impact speed of 5 m/s,
and impact mass of 7.1 kg, 8-layered PMMA/PC did not show benefits over the monolithic PC
in terms of energy dissipation. However, in the same impact conditions, with the decrease of the
aspect ratio to a critical value, 0.07, the 8-layered PMMA/PC started to dissipate more energy than
the monolithic PC. The same phenomena was observed in the high speed impact simulations with
a impact mass of 0.0032 kg and impact speed of 300 m/s. When the aspect ratio was small, the
deformation of the PC plate was governed with local plastic deformation. However, in the layered
PMMA/PC laminate, the hard PMMA layers severed as back plates to support the ductile PC
layers, and the ductile PC layers served as barriers to prevent the large crack of PMMA coalescing
to lead to a catastrophic failure. Thus, at small aspect ratios, the layered PMMA/PC laminate
has an improved impact resistance compared to the monolithic PC. With the increase of aspect
ratios, the benefits of the laminate structure tended to disappear. The reason was that the failure
mechanism of PMMA and PC is different when it was impacted with small and large aspect ratios.
As shown in Fig. 3.9, the failure of the PMMA plates was governed with brittle spall-induced
53

damage at small aspect ratios. However, at large aspect ratios (Fig. 3.7), a large scale of radial
damage was formed at the bottom surface of the plate due to the high bending force induced by the
big projectile. For the PC plate, with the increase of aspect ratio, more areas (the immediate impact
zone under the projectile) were involved in the plastic deformation. Thus, the role of PMMA plates
in magnifying the plastic damage zone was not significant.
For the purpose of developing a hybrid ductile/brittle transparent laminate that could stand a
wide range of impact loading, an optimization of 8-layered PMMA/PC laminate was conducted.
The residual kinetic energy was shown to have a negative correlation with the post yielding stress
of PMMA and PC. With an increase in post yielding stress of PMMA from 254 MPa to 351 MPa,
the energy dissipated only increased 2.8 J. Relatively, with an increase in the post yielding stress of
PC from 91.7 MPa to 182.3 MPa, the energy dissipated increased 17.2 J. Thus, it seemed that the
impact resistance of 8-layered PMMA/PC laminate is more sensitive to the strain rate hardening
of PC.

3.5

Conclusion

The failure behavior of monolithic PMMA and PC and layered PMMA/PC laminate subjected
to low and high speed impact loading were numerically investigated. Boyces three-dimensional
physical-based models [9] were implemented in ABAQUS to predict the strain rate and pressure
dependent behavior of PC. The failure mechanism, fragmentation pattern, and energy dissipation
were compared through the simulations. The failure of monolithic PMMA was dominated by brittle cracking. Many spalls were observed on the back surface of the immediate impact region at
small aspect ratio. The deformation of monolithic PC was dominated by local plastic deformation
which occurred in the immediate impact zone under the projectile. For both of the monolithic PC

54

and layered PMMA/PC laminate, the fraction of kinetic energy decreased with aspect ratio. There
is found to exist a critical aspect ratio, 0.07 for both the high velocity impact (300 m/s) and low
velocity impact (5 m/s). When the aspect ratio decreased below this critical value, the 8-layered
laminate started to dissipate more energy than the monolithic PC. The 8-layered PMMA/PC laminate magnified the plastic deformation zone and decreased the brittle damage area. An optimal
study by incorporating MATLAB optimization tool box and ABAQUS FEM simulations was conducted to maximize the impact resistance of the transparent laminated panel. It was found that
the impact resistance or the overall ballistic performance of these PMMA/PC laminated systems is
relatively more sensitive to strain rate hardening of PC than PMMA.

55

REFERENCES

56

REFERENCES
[1] Alex J Hsieh, Daniel DeSchepper, Paul Moy, Peter G Dehmer, and John W Song. The effects
of pmma on ballistic impact performance of hybrid hard/ductile all-plastic-and glass-plasticbased composites. Technical report, DTIC Document, 2004.
[2] Alex J Hsieh and John W Song. Measurements of ballistic impact response of novel coextruded pc/pmma multilayered-composites. Journal of reinforced plastics and composites,
20(3):239–254, 2001.
[3] Parimal J Patel, Gary A Gilde, Peter G Dehmer, and James W McCauley. Transparent ceramics for armor and em window applications. In International Symposium on Optical Science
and Technology, pages 1–14. International Society for Optics and Photonics, 2000.
[4] S Sahraoui and JL Lataillade. Deformation and fracture of pmma at high rates of loading.
Journal of applied polymer science, 51(9):1527–1532, 1994.
[5] SC Wright, NA Fleck, and WJ Stronge. Ballistic impact of polycarbonatean experimental
investigation. International journal of impact engineering, 13(1):1–20, 1993.
[6] Ravi Shankar Kalluri. Failure of transparent polymer composite laminated glass panels under
impact loading. PhD thesis, University of Missouri–Columbia, 2007.
[7] Qasim H Shah. Impact resistance of a rectangular polycarbonate armor plate subjected to
single and multiple impacts. International Journal of Impact Engineering, 36(9):1128–1135,
2009.
[8] Xin Sun, Kevin C Lai, Tara Gorsich, and Douglas W Templeton. Optimizing transparent armor design subject to projectile impact conditions. Advances in Ceramic Armor IV: Ceramic
Engineering and Science Proceedings, Volume 29, Issue 6, pages 15–22, 2008.
[9] Ellen M Arruda and Mary C Boyce. Evolution of plastic anisotropy in amorphous polymers
during finite straining. International Journal of Plasticity, 9(6):697–720, 1993.
[10] C Bauwens-Crowet. The compression yield behaviour of polymethyl methacrylate over a
wide range of temperatures and strain-rates. Journal of Materials Science, 8(7):968–979,
1973.
[11] C Bauwens-Crowet, J-C Bauwens, and Georges Homes. The temperature dependence of
yield of polycarbonate in uniaxial compression and tensile tests. Journal of Materials Science, 7(2):176–183, 1972.
57

[12] MC Boyce and EM Arruda. An experimental and anaiytical investigation of the large strain
compressive and tensile response of glassy polymers. Polymer Engineering & Science,
30(20):1288–1298, 1990.
[13] Mary C Boyce, Ellen M Arruda, and R Jayachandran. The large strain compression, tension,
and simple shear of polycarbonate. Polymer Engineering & Science, 34(9):716–725, 1994.
[14] F Rietsch and B Bouette. The compression yield behaviour of polycarbonate over a wide
range of strain rates and temperatures. European Polymer Journal, 26(10):1071–1075, 1990.
[15] AS Argon. A theory for the low-temperature plastic deformation of glassy polymers. Philosophical Magazine, 28(4):839–865, 1973.
[16] AD Mulliken and MC Boyce. Mechanics of the rate-dependent elastic–plastic deformation of
glassy polymers from low to high strain rates. International journal of solids and structures,
43(5):1331–1356, 2006.
[17] Herzl Chai and Guruswami Ravichandran. On the mechanics of fracture in monoliths and
multilayers from low-velocity impact by sharp or blunt-tip projectiles. International journal
of impact engineering, 36(3):375–385, 2009.
[18] Md Akil Hazizan and WJ Cantwell. The low velocity impact response of foam-based sandwich structures. Composites Part B: Engineering, 33(3):193–204, 2002.
[19] Yanxiong Liu and Benjamin Liaw. Drop-weight impact tests and finite element modeling of
cast acrylic/aluminum plates. Polymer Testing, 28(8):808–823, 2009.
[20] JD Pearson, MA Zikry, M Prabhugoud, and K Peters. Measurement of low velocity and
quasi-static failure modes in pmma. Polymer composites, 28(3):381–391, 2007.
[21] Joshua S Stenzler and NC Goulbourne. The effect of polyacrylate microstructure on the
impact response of pmma/pc multi-laminates. International Journal of Impact Engineering,
38(7):567–576, 2011.
[22] JW McCauley, E Strassburger, P Patel, B Paliwal, and KT Ramesh. Experimental observations on dynamic response of selected transparent armor materials. Experimental Mechanics,
53(1):3–29, 2013.
[23] Shinji Ogihara, Tomoyuki Ishigure, and Akira Kobayashi. Study on impact perforation fracture mechanism in pmma. Journal of materials science letters, 17(8):691–692, 1998.

58

[24] S Sarva, AD Mulliken, and MC Boyce. The mechanics of large-strain inhomogeneous deformation of polymeric materials under dynamic loading conditions. In Journal de Physique IV
(Proceedings), volume 134, pages 95–101. EDP sciences, 2006.
[25] L Roy Xu and Ares J Rosakis. An experimental study of impact-induced failure events in
homogeneous layered materials using dynamic photoelasticity and high-speed photography.
Optics and Lasers in Engineering, 40(4):263–288, 2003.
[26] Michael Hebert, Carl-Ernst Rousseau, and Arun Shukla. Shock loading and drop weight
impact response of glass reinforced polymer composites. Composite Structures, 84(3):199–
208, 2008.
[27] Srinivasan Arjun Tekalur, Arun Shukla, and Kunigal Shivakumar. Blast resistance of polyurea
based layered composite materials. Composite Structures, 84(3):271–281, 2008.
[28] A Shokuhfar, SMR Khalili, F Ashenai Ghasemi, K Malekzadeh, and S Raissi. Analysis and
optimization of smart hybrid composite plates subjected to low-velocity impact using the
response surface methodology (rsm). Thin-Walled Structures, 46(11):1204–1212, 2008.
[29] AN Shupikov, SV Ugrimov, AV Kolodiazhny, and VG Yareschenko. High-order theory
of multilayer plates. the impact problem. International journal of solids and structures,
35(25):3391–3403, 1998.
[30] NV Smetankina, S Yu Sotrikhin, and AN Shupikov. Theoretical and experimental investigation of vibration of multilayer plates under the action of impulse and impact loads. International journal of solids and structures, 32(8):1247–1258, 1995.
[31] X Chen, N Chandra, and AM Rajendran. Analytical solution to the plate impact problem
of layered heterogeneous material systems. International journal of solids and structures,
41(16):4635–4659, 2004.
[32] Mahmoud M Reda Taha, Arife B Colak-Altunc, and Marwan Al-Haik. A multi-objective optimization approach for design of blast-resistant composite laminates using carbon nanotubes.
Composites Part B: Engineering, 40(6):522–529, 2009.
[33] TY Kam and FM Lai. Maximum stiffness design of laminated composite plates via a constrained global optimization approach. Composite structures, 32(1):391–398, 1995.
[34] MC Ray and A Baz. Optimization of energy dissipation of active constrained layer damping
treatments of plates. Journal of Sound and Vibration, 208(3):391–406, 1997.

59

[35] RN Haward and G Thackray. The use of a mathematical model to describe isothermal stressstrain curves in glassy thermoplastics. Proceedings of the Royal Society of London. Series A.
Mathematical and Physical Sciences, 302(1471):453–472, 1968.
[36] Wei Zhang, Srinivasan Arjun Tekalur, and L Huynh. Impact behavior and dynamic failure of
pmma and pc plates. In Dynamic Behavior of Materials, Volume 1, pages 93–104. Springer,
2011.

60

Chapter 4
The Effects of Damage Accumulation on
Tensile Strength and Toughness of Compact
Bovine Bone subjected to Static Tensile
Loading, Charpy Impact Loading, and
SHPB High Strain Rate Tensile Loading

4.1

The Effects of Damage Accumulation on Tensile Strength
and the Toughness of Compact Bovine Bone subjected to
Static Tensile Loading

4.1.1

Introduction

As a tough and light-weight natural composite, bone works as a scaffold to provide support and
protection to internal organs. Understanding the failure mechanisms in bone is of great medical
interest for preventing bone fracture failure and improving bone fracture treatment. Bone fractures
are a result of tensile failure subsequent to bending or torsion because bone is more susceptible to
61

tensile failure compared to compressive failure[1, 2, 3, 4]. Bone is a strain rate dependent material due to its viscoelastic properties where increases in loading rate have resulted in increases in
Young’s modulus and ultimate stress [5, 6, 7, 8, 9, 10, 11]. The mechanical properties of compact
bone are also related to its hierarchical microstructure. The strength, energy absorption, and modulus of elasticity decrease with the percentage area of osteons [11, 12, 8, 13, 14]. Interestingly,
little emphasis has been placed on the role of microcracking on the tensile failure behavior of bone
although, like quasi-brittle materials, the failure of bone is characterized as the result of damage
accumulation in the form of microcracks.
Microcrack accumulation has been shown to yield a loss of stiffness in compact bone under
cyclic fatigue testing [15, 16]. Microcracking was also found to be a primary toughening mechanism that resulted in increased R-curve behavior [17, 18, 19, 20] applied a series of fluorescent
chelating agents to monitor the crack growth during cyclic compression fatigue finding that microcracking initiated in the interstitial bone tissue in the early stage of the fatigue life and that
longer microcracks developed in the longitudinal direction as opposed to the transverse direction.
In recent efforts, Zioupos [21] examined the plastic strain of compact canine femurs, subjected to
tensile loading at various strain rates, to the number of cracks (imaged by staining), concluding
that the plastic strain had a positive correlation with the number of cracks. However, each of the
above studies was limited by the resolution limits imposed by Frost’s [22] bulk fuchsin staining
technique, which requires the tested specimens to be cut into microsections. Sectioning can introduce new microcracks during cutting; and, furthermore, the small field of view compromises the
accuracy of the crack length measurement. Several attempts have been made to avoid these limitations. Zioupus [23] used the acoustic emission method to monitor crack growth in bone and antler
specimens under uniaxial tensile loading where microcracks were found to initiate upon yielding
in both antler and bovine tibia specimens. However, the shortcoming of the technique is that the
62

growth and distribution of microcracks cannot be directly visualized. Using the laser illumination
technique, Reilly [24] found that the microcracks developed on the tensile side of bone specimens
subjected to four-point bending. These tensile microcracks first appeared at a strain of 0.004 and
the peak growth of the microcracks occurred when the strain reached 0.008. However, due to the
step loading and examination, the stress-strain curve was not continuously obtained throughout the
test.
In the present dissertation, an innovative optical technique was used to monitor the damage
accumulation during the process of tensile failure of compact bovine femur. A high speed camera and a LED panel were set up to record the light illumination of the microcracks initiated and
propagating on the specimen surface when stretched under a multipurpose tensile testing system
(MTS). By using the strong light scattering effects of crack edges, microcracks are visible as features with higher light intensity. To evaluate the correlation between the developing microcracks
and the bone microstructure, specimens from several bovine femurs of different ages were examined. Strain rate effects were also investigated by examining 10 specimens from one piece of
bovine femur with strain rates ranging from 0.0001 to 0.012 s−1 . As a natural composite material,
the mechanical properties of bone depend on the fiber or osteon direction, which runs along the
loading axis of bone. It has been shown that osteons potentially act as barriers to prevent coalescence of the microcracks initiated in the interstitial bone tissue [25, 26]. Thus, the failure pattern
of bone loaded in the transverse (perpendicular to the osteon direction) and the longitudinal direction (along osteon direction) are expected to be different. To evaluate such orientation effects,
the compact bovine specimens machined from one piece of bovine femur bone were loaded in the
longitudinal and transverse directions respectively. The stress-strain curves were correlated to the
high-speed images to examine the microcrack growth in the specimens. Post-failure examinations
were performed with a scanning electron microscope (SEM). Both the front flat surface and the
63

fracture surface were examined.

4.1.2

Material and Methods

The tensile specimens were machined from fresh compact femurs (12 to 30 months old) provided
by a local beef abattoir. The mid-diaphysis of the bovine femur was first cut into rough rectangular beams using a band saw. Then, the beams were machined into 60x8x2 mm (longitude) or
30x8x2 mm (transverse) rectangular shapes using a milling machine (ProCut RL80RF311) at a
low feed speed. Subsequently, a Dremel hand tool was used to shape the rectangular beams into
the dog-bone shape specimens with the exact same geometry. Finally, the specimens were manually polished with 80, 100, 600, 1000 grit. During the entire machining process, the specimens
were kept wet with the phosphate buffered saline (PBS) solution. The geometry and dimensions
are shown in Fig. 4.1. After machining, the specimens were covered with gauze dipped in PBS
solution and preserved in a -20◦ C freezer until testing. Specimens were thawed just prior to testing.
Overall, 28 specimens were machined, including 12 specimens from one femur (30-months old)
for strain rate effects, 10 specimens from several femurs with various ages (12 to 30 months old)
for microstructure effects, and 6 specimens from one femur (30-months old) to examine orientation
effects.
Uniaxial tensile testing was performed using a MTS Insight tensile testing system with a 10
kN load cell. A high speed camera (Phantom V12) was set facing the front side of the specimen to
record the fracture process. In order to visualize the microcracks, the tests were performed in a dark
room and a strong cold light source was set at the side of the camera to illuminate the specimen. A
laser extensometer (EIR Model LE-05) was set up facing the back side of the specimen to measure
the extension. Both the high speed camera and the laser extensometer were triggered at the same
time as the MTS started to apply the load. The rationale of the research approach is not directed
64

1

Figure 4.1 Geometry and dimensions of dog-bone specimens. (a), specimen orientation diagram
2
(1, longitudinal specimens,
loaded along the osteon direction; 2, transverse specimens, loaded
perpendicular to the osteon direction); (b), the dimension of the longitudinal specimens (the length
of the transverse specimens is half of the longitudinal specimens)
to quantify the exact amount of microcracking in bone but rather indirectly provide a ”real-time”
relationship between microcracking and strain development. It is proposed that there exist critical
strain rate and strain levels to initiate and propagate microcracking in bone. It has never been
clearly shown what these levels are. More importantly, the orientation of bone and the osteon size
affects initiation and propagation of microcracking. The experimental design is geared to bring out
these effects. Since we indirectly correlated microcracking by the diffuse light pattern, questions
might arise on the validity of using the light intensity as a measure of microcracking. It is true that
many materials exhibit a change in color or appearance while subjected to loading (e.g., crazing
in polymers). But it is different in the case of bone where the change in light intensity is directly
related to microcracking. Validation for the current is provided by Zioupos et al. [21]and in our
experimental measurements themselves where the light diffusion gradually increases as the load
is applied and reaches peak value at peak load (Fig. 4.2). Exact quantification of light intensity
to density of microcracking could be undertaken but again is not the primary focus of the current
work.
The specimens were covered with sterile PBS soaked gauze before testing, where the duration

65

120

14
12
10
8
6

80

100
Stress (MPa)

Stress (MPa)

20
18
16

0.0049 ε

4
2
0

60
40

0.0203 ε

20
0

0

0 .0 0 2

0 .0 0 4

0

0 .0 0 6

100

100

80

80

60

0.0243 ε

Stress (MPa)

120

Stress (MPa)

120

40
20

0 .0 2

0 .0 3

60

0.0298 ε

40
20

0
0

0 .0 1

0 .0 2

0

0 .0 3

0

S tr ai n (m m / m m )
120

120

100

100

80

80

60

0.0310 ε

40

0 .0 2

0 .0 4

S tr ain (m m / m m )

Stress (MPa)

Stress (MPa)

0 .0 1

S tr ai n (m m / m m )

S tr ain (mm/mm)

20

60

0.0315 ε

40
20

0
0

0 .0 2

0

0 .0 4

0

S tra i n (m m /m m )

0 .0 2

0 .0 4

S tr ain (m m / m m )

Fig. 2 Representative
images showing
real showing
time opticalreal
effects
during
tensile effects
stretchingduring
of bovine
bone sample.
The frame
Figure
4.2 Representative
images
time
optical
tensile
stretching
of rate
of this isbone
1000 fps
.This technique
is validated
Zioupos
(Zioupos
bovine
sample.
The frame
rate ofbythis
is 1000
fps.et al., 2008).

of the test lasted about one minute. Thus, the specimens were considered wet. For the strain
rate evaluation, the cross-head speeds ranged from 0.5 mm/min to 10 mm/min corresponding to
strain rates of 0.0001 to 0.0012 s−1 . All other tests were performed at a cross-head speed of 0.5
mm/min. The frame rates of the high speed camera were 200 s−1 , 800 s−1 , and 1500 s−1 for tests
with cross-head speeds of 0.5 mm/min, 1 to 2 mm/min, and larger than 3 mm/min, respectively.

66

19

Images from the camera were used to calculate the light intensity in each frame by a MATLAB
code. The light intensity was represented by the normalized gray scale value in each frame, which
was calculated by summarizing the gray scale value of each pixel in each frame after subtracting
the corresponding value in the reference frame and then normalizing by the numbers of pixel in
the frame.
The post-failure fracture damage of the tested specimens was examined in a JEOL JSM-6400V
(lanthanum hexaboride electron emitter) scanning electron microscope. Digital images were acquired using analySIS Pro software version 3.2 (Olympus Soft Imaging Solution Corp., Munster,
Germany). Before the examination, samples were coated with gold ( 20 nm thickness) in an Emscope Sputter Coater model SC 500 (Ashford, Kent, England) purged with argon gas. Samples
were mounted on aluminum stubs using high vacuum carbon tab.

4.1.3

Experimental Results

4.1.3.1

Imaging Longitudinal Specimens Loaded in the Same Strain Rate of 0.0001 s−1

A higher light intensity was expected when the microcracks initial 0.0001 s−1 ted on the specimen
surface due to the light scattering effect of the crack edges. High speed grayscale images extracted
from the high speed videos were correlated to the stress-strain curves. The microcrack illumination
of two typical specimens, one with dispersed microcracks and one with coalesced microcracks, are
shown in Fig. 4.3a, b, respectively. It was found that microcrack illumination began to occur at
yielding and the specimens with dispersed microcracks developed a larger plastic strain as compared to the specimen with local cracks, which lead to quick failure after yielding. The yielding
point was determined by Currey’s [7] method, which is the intersection of the tangents of the linear portion and the post-yield portion of the stress-strain curve. The changing in the light intensity

67

with time was calculated for each specimen. The plots of light intensity and stress over time for
the two typical specimens ( Fig. 4.3a, b) are shown in Fig. 4.3c, d. The calculated results agreed
with the visual observation from the images that, for both of these two specimens, the significant
increases in light intensity occurred at onset of yielding. After initiation, light intensity gradually
increased with time until failure. The specimen with a larger plastic strain (typically 2.0%) developed more microcracks leading to a higher microcracking illumination compared to the specimen
with a smaller plastic strain (typically 0.5%).

Figure 4.3 Correlations of high speed images with stress-strain curves for two typical longitudinal
compact bovine femur specimens tested under a constant strain rate of 0.0001 s−1 : (a) stress-strain
curve and the corresponding high speed images for a typical specimen with a higher failure strain;
(b) stress-strain curve and the corresponding high speed images for a typical specimen with a lower
failure strain; (c) the corresponding light intensity and the stress relation for the specimen shown
in (a); (d) the corresponding light intensity and the stress relation for the specimen shown in (b).

68

4.1.3.2

Imaging Strain Rate Effects

Twelve specimens were tested under strain rates of 0.0001 to 0.0012 s−1 to evaluate strain rate
effects on the stress-strain behavior of compact bovine femur bone. The experimental results
showed that the stress-strain curves fell into two groups: one turned as lower strain rate group
with strain rates lower than 0.0007 s−1 ; and the other turned as higher strain rate group with the
strain rates equal to or greater than 0.0007 s−1 . In the high strain rate group, specimens showed
viscoelastic behavior with non-linear stress-strain curves, even in the elastic region. The Young’s
modulus for this group was calculated as dynamic Young’s modulus by dividing the ultimate stress
by failure strain. The average Young’s modulus of the higher strain rate group was 41.9 GPa,
which was about twice that of the lower strain rate group (21.9 GPa). However, the failure strain
of the higher strain rate group (0.3%), was about half of the failure strain of the lower strain rate
group (0.7%). The difference between the average ultimate stress for the two groups (which was
123.4 MPa for the higher strain rate group and 105.2 for the lower one) was not significant. The
plots of Young’s modulus, energy absorption capacity, ultimate stress, and the failure strain as
a function of strain rate are shown in Fig. 4.4c, d. Young’s modulus and ultimate stress have
a positive correlation with strain rate, and the energy absorption capacity and the failure strain
have a negative correlation with strain rate. The correlation of the stress-strain curve and the high
speed images for a typical specimen subjected to a higher strain rate is shown in Fig. 4.4a. The
corresponding plots of stress-microcracking illumination and strain-microcracking illumination are
shown in Fig. 4.4b. It was shown that the microcracking illumination of the specimen initiated at
the onset of yielding, then rapidly increased until failure.

69

b)
140

4.5
Stre ss

100

120

4

Light Inte nsity

80

3.5

100

60

Stress (MPa)

Stress (MPa)

120

40
20
0
0

0.001
0.002
S train (m m /m m )

0.003

3

80

2.5

60

2
1.5

40

Light Intensity

a) 140

1
20

0.5

0

0
0

2
Time (s)

3

4

d)
0.012

160

60

0.01
0.008

100

0.006

80
60

0.004

40
0.002

20

0
0

0.0005
0.001
Strain rate (s -1 )

0
0.0015

Stress (MPa)

120

Young's modulus (GPa)

140

0.8
0.7

50

0.6

40

0.5

30

0.4
0.3

20

0.2

10

0.1

0
0

0.0005
0.001
Strain Rate(s -1 )

Energy absorption
capacity(J/mm 3 )

c)

Strain

1

0
0.0015

Figure 4.4 Strain rate effects on the mechanical behavior of the longitudinal compact bovine femur
Fig. 4 Strain rate effects on the mechanical behavior of the longitudinal compact bovine femur specimens:
specimens:
(a) stress-strain curve with the corresponding high speed images for a typical speci(a) stress-strain curve with the corresponding high speed images for a typical specimen under a strain rate of 0.0012s-1
men
under a strain rate of 0.0012s−1 ; (b) the corresponding relation of stress and microcracking
(b) the corresponding relation of stress and microcracking illumination for the specimen
illumination
the specimen;
(c) the stress
variation
of failure
andFailure
ultimate
stress
strainUltimate
rate,
(c) the variationfor
of failure
strain and ultimate
with strain
rate, ∆strain
represents
strain
and □with
represents
represents Failure strain and represents Ultimate stress; (d) the variation of Young’s modulus
stress
and
energy
absorption
with strain
rate,absorption
represents
Young’s
modulusYoung’s
and modulus
represents
(d) the
variation
of Young’s modulus
and energy
with strain
rate, ∆ represents
and □ Energy
represents
absorption
capacity.
Energy absorption
capacity
4.1.3.3

Imaging Orientation Effects at a Loading Rate of 0.0001 s−1

The stress-strain curves of the transverse and longitudinal specimens are shown in Fig. 4.5a. It
can be seen that the transverse specimens exhibited a strain softening behavior unlike the strain
hardening that was typically observed on the longitudinal specimens. The ultimate stress and the
21

Young’s modulus of the longitudinal specimens were about twice of the transverse specimens (Fig.
4.5b). The plots of the stress-microcracking illumination and strain- microcracking illumination

70

of a typical transverse specimen are shown in Fig. 4.5c, d, respectively. Unlike the smooth growth
of the post-yield light intensity observed on the longitudinal specimens, it was found that the light
intensity of the transverse specimens had a significant jump after yielding.

a)

b)

Longitudinal

100

80

80

60

60

40

Transverse

36.2

40

20

21.8

20

0
0

0.002

0.004

0.006

0.008

10.9

0

0.01

Young's modulus
(GPa)

Strain(mm/mm)

d)

50

Stress
Light Inte nsity

20
10

20
Time (s)

30

50
40

80

30

60
20

40

10

20

0
10

Ultimate stre ss
(MPa)

Stre ss
Light Intensity

100

40
30

0

120

Stress (MPa)

45
40
35
30
25
20
15
10
5
0

Light Intensity

Stress (MPa)

c)

93.5

Longitudinal
Transve rse

0

40

Light Intensity

Stress(MPa)

100

0
0

10

20

30

40

50

Time (s)

28
29
30
31
32
33
34
35
36
37
38

Fig.
3 The
orientaiton
effects
on the
mechanical
properties
of the
compact
bonvine
under a
Figure
4.5 The
orientation
effects
on the
mechanical
properties
of the
compact
bovinefemur
femurspecimens
specconstant
of
0. longitudinal 0.0001s-1:
imens under a constant strain
strain rate of 0.0001 rate
s−1 : (a) the stress-strain
curves of the
and transverse
specimens;
comparison
modulus and and
ultimatetransverse
stress of the specimens
(a)
the
stress-strain(b) the
curves
of of the
theYoung’s
longitudinal
transverse
specimens
longitudinal
specimens;
(c) plot
of of
stress
light intensity
changing
(b)
the comparison
of to
thethe
Young's
modulus
and ultimate
stress
the and
transverse
specimens
to the longitudinal
with
time
for
a
typical
transverse
specimen,
represents
stress
and
solid
line
represents
light
intenspecimens
sity;
(d)ofplot
of stress
and intensity
light intensity
changing
withfor
time
for a typical
longitudinal
(c)
plot
stress
and light
changing
with time
a typical
transverse
specimen,specimen,
□ represents stress and
represents stress and solid
solid
lineline represents light intensity
represents
light
intensity
(d) plot of stress and light intensity changing with time for a typical longitudinal specimen, □ represents stress and
solid line represents light intensity

4.1.3.4

SEM Post-failure Examination

The post-failure photographs of the tested specimens were examined by the SEM. For specimens
subjected to the same loading rate, the examination focused on the histological structures and the
qualities aspects of microcrack lengths and densities. Fig. 4.6a, b shows the fracture surfaces of
71

the specimens at a failure strain of 2.4% (a) and 0.9% (b). The SEM images show the hierarchical structure of the compact bone, which includes planar lamellae (primary bone), and osteons
characterized by cylindrical lamellae surrounding the Haversian channel. Six typical specimens
were selected to calculate the osteon area percentage for the entire fracture surface by following
Saha’s method [13]. The results (Fig. 4.6c) agreed well with the literature showing that the energy
absorption capacity and failure strain are negatively correlated with the area percentage of the osteons [7, 8]. Four typical specimens subjected to the same loading rate (0.0001 s−1 ) were chosen to
examine the flat surfaces close to the fracture surface (Fig. 4.7). Higher micro-crack densities and
long linear macrocracks over five hundred microns (Fig. 4.7a) were observed on the specimens
with higher failure strains.
The typical post-failure SEM images of the longitudinal and transverse specimens are shown
in Fig. 4.8. Compared to the microcracks that developed on the longitudinal specimens (Fig.
4.8b), the crack length and density of the transverse specimens were much higher (Fig. 4.8a).
The corresponding fracture surfaces of the two typical specimens were shown in Fig. 4.8 c, d,
respectively. It can be seen that the fracture surface of the transverse specimen (c) was much
smoother in compared to the longitudinal specimen (d). The fracture surfaces of the specimens for
the strain rate effects were examined. The SEM images for two typical specimens are shown in
Fig. 4.9. Osteon pullout was observed on both of the specimens. But the specimen subjected to
the strain rate of 0.0012 s-1 seemed to show more sharp tearing fracture surface (Fig. 4.9c, d) as
compared to the specimen under lower strain rate (Fig. 4.9a, b).

4.1.4

Discussion

Innovative use of a relatively simple optical technique was used to study the accumulation of damage in the form of microcracking in compact bovine bone specimens subjected to tensile loading.
72

a)

b)

2.5

3

Energy absorption (J/mm )

c) 2.5

2

1.5

1.5

1

1

0.5

Strain%

2

0.5

0

0
0

20

40
60
80
Area occupied by osteon%

100

3

Energy absorption (J/mm )

Figure 4.6 Microstructure effects on the mechanical properties of the longitudinal compact bovine
femur specimens subjected to a constant strain rate of 0.0001 s−1 : (a) SEM fracture surface of a
typical specimen with a failure strain 2.4%; (b) SEM fracture surface of a typical specimen with a
failure strain of 0.9%; (c) plots of failure strain and energy absorption capacity to the percent area of
osteons of the specimens. represents failure strain, ♦ represents energy absorption capacity. Dot
line is the linear trend line for failure strain; solid line is the linear trend line for energy absorption
capacity. Osteons characterized by cylindrical lamellae surrounding the Haversian system, which
2.5
are circled in a, b. Primary bone characterized by laminar structures (see arrows
in b). The area
percentages of osteon
in the
whole fracture surface
were
following
Saha’s
methods [13].
Fig.
5 Microstructure
effects on
thecalculated
mehcanical
properties
2 of the longitudinal compact
bovine femur speicmens subjected to a constant strain rate of 0.0001s-1: (a) SEM fracture
surface
a
typcial
specimen and with
a1.5to the
failure
strain
2.4%.
This powerful approach
relates of
the orientation,
microstructure,
strain rate
resulting
mi(b) SEM fracture surface of a typical specimen with a failure strain of 0.9%.
1
(c) plots
of failure
and energy fracture
absorption
to the were
percent
area of secondary
cro damage in femoral
bovine
bone. strain
The post-failure
andcapacity
flat surfaces
examined
0.5 energy absorpiton capacity.
osteons of the specimens. ∆ represents failure strain,◊ represents
by SEM. The tensileDot
properties
are
in
good
agreement
with
the
literatures
(Table
4.1). trendline
The ini- for energy
line is the linear trendline for failure strian, solid line is the linear
0
absorption capacity. Secondary osteons (SO) characterized by cylindrical lamellae
0
20
40
60
80
tiation and propagation
of the microcracks in the specimens is indicated by its light
illumination
surrounding the Haversian system, which are circled in a, b. Primary osteons (PO)
Area of osteon%
characterized by laminar structures, which are indicated by white arrows in b.The area
visualized by a high speed camera. The correlation of the high speed images to the stress-strain
percentages of the secondary osteon in the whole fracure surface were calculated following
Saha’s methods
[39]. Magificatin
ofafter
the SEM
images:
100x.
curves showed that microcracks
initiated
immediately
yielding
(Fig.
4.2a, b, and Fig. 4.3a). It

was found that the specimens loaded along the osteon direction showed different failure patterns as
compared to the specimens loaded along the transverse osteon direction. For the longitudinal spec73

5

mac

mic
mic

mic

Failure strain: 2.4%

Failure strain: 1.9%

Failure strain: 0.9%

Failure strain: 0.5%

Fig. 6 Representative
SEM post-failure
flat surface flat
images
(close to
the fracture
of four
typical edge)
longitudinal
Figure
4.7 Representative
SEM post-failure
surface
images
(closeedge)
to the
fracture
of
−1
specimens
under
a
constant
strain
rate
of
0.0001s-1.
Example
features
are
labeled:
(mac)
macrocrack
with
a crack
four typical longitudinal specimens under a constant strain rate of 0.0001 s . Example features
length
over 500
μm, macrocrack
(mic) microcrack
a crack
lengthover
smaller
μm. microcrack
Long linear macrocracks
are
labeled:
(mac)
withwith
a crack
length
500than
µm,500
(mic)
with a crackand
many microcracks were observed on the specimen with a large failure strain (a), but not on the specimens with
length smaller than 500µm. Long linear macrocracks and many microcracks were observed on the
small failure strain (c, d). Magnification: 200x.
specimen with a large failure strain (a), but not on the specimens with small failure strain (c, d).

imens, the microcracking increased gradually after yielding until failure and microcracks tended
to toughen the specimens by showing strain hardening behavior. The specimens that exhibited
more microcracks were shown to absorb more plastic energy before failure and the microcracking
tended to be more localized with
increasing strain
rate. However, Young's
for the transverse
specimens, the
Ultimate Failure
Loading Specimen
-1
Specimen

Strain rate s

Modulus
(GP
)
10.5-36.8

Stress
(MP
)
125-348

9.3-19.4
10-21.7
18-23.3
17.7-40.4
11-12
12.4-48.3

223-310
33-50
NR
99-271
110-250
133-400

strain%
0.3-1.2
1.5-3.2
0.2-0.4
NR
NR
3-6
0.7-2.9

Young's
Modulus (GPa)
18.3-54.8
11.3-17.5
12.6
11-12
22.8
17.7-40.4

Ultimate
Stress (MPa)
83.0-135.8
105-193
121.3
110-250
107.4
99.2-271.4

Failure
strain%
0.2-2.4
0.9-1.8
1.3
3-6
NR
NR

Energy-absortion
6
3
capacity (10 J/m )
0.2-1.8
NR
NR
3.8-6.0
NR
0.4-1.3

B26_T2_2

B26_B2_1_I B16_M1_4B15_M1_3

Type
condition
Tensile
PartiallyWet
CompressivDry
Tensile
Dry
Tensile
Dry
Tensile
Wet
Tensile
Wet
CompressivWet

microcracking instantly increased to its peak value remaining constant till failure. Strain softening
Present results
F. Ferreira (2006)
M.Pithious (2004)
F. Katsamanis (1990)
T.M. Wright (1976)
R.D. Crowninshield (1974)
J.H. McElhaney(1966)

290-820
368-795
1000 m/s
100
0.00053-237
0.00167-250
0.001-1500

Bovine Femur
Bovine Femur
Bovine Shaft
Human Femur
Bovine Femur
Bovine Tibia
Bovine Femur

behavior was observed on the transverse specimens.

Table 4.1 Comparison of mechanical tensile properties among different investigators.

Specimen
condition
Present results
0.0001-0.0012
Bovine Femur
Wet
M.Pithious (2004)
0.5-500 mm/min Bovine Shaft
Dry
S.Saha (1973)
133
BovineTibia+Femur Wet
R.D. Crowninshield (1974) 0.00167-250
Bovine Tibia
Wet
J.D. Currey (1975)
0.00013-0.16
Bovine Femur
Wet
T.M. Wright (1976)
0.00053-237
Bovine Femur
Wet
Strain rate s

-1

Specimen

7

Cross head speed
(mm/min)

Microcracking in bone has been shown to contribute to the rising R-curve behavior [27, 25, 17,
-1

Strain Rate s

0.5

0.5

0.5

0.5

0.0001

0.0001

0.0001

0.0001

18, 19] because its toughening mechanism is similar to brittle materials such as ceramics [28, 29].
Young's modulus
(GPa)
Ultimate stress(MPa)
Ultimate strain
Energy density(N/mm^2)
Failure strain%

20.57

21.94

19.61

16.68

94.42
112.93
117.13
It was also found that microcracking tended to impair the mechanical
properties
of compact114.46
bone

74

0.0102
0.7335
1.02

0.0054
0.3627
0.54

0.0192
1.8244
1.92

0.0237
2.3603
2.37

A
Longitudinal

Transverse

A

b)

Flat surfaces

a)

d)

Fracture surfaces

c)

Fig. 6 Representative SEM post-failure images of two typical specimens (under a strain rate of 0.0001s-1) for
Figure 4.8 Representative SEM post-failure images of two typical specimens (under a strain rate
orientation evaluation:
of 0.0001 s−1 ) for orientation evaluation: (a) the flat surface of the transverse specimen; (b) the flat
(a) the flat surface of the transverse specimen, Magnification: 19x.
surface of the longitudinal specimen. Many linear microcracks were observed on the transverse
(b) the flat surface of the longitudinal specimen, Magnification: 19x. Many linear microcracks were observed on
specimens compared to the longitudinal specimen with fewer microcracks; (c) the corresponding
the transverse specimens compared to the longitudinal specimen with fewer microcracks.
fracture surface of the transverse specimen; (d) the corresponding fracture surface of the longitu(c) the corresponding fracture surface of the transverse specimen, Magnification: 100x.
dinal specimen. Smoother tearing surfaces (indicated by arrows) were observed on the transverse
(d) the corresponding fracture surface of the longitudinal specimen, Magnification: 100x.
specimen in compared to the rough fracture surface of the longitudinal specimen.
Smoother tearing surfaces were observed on the transverse specimen indicated by the white arrow in compared to
the rough fracture surface of the longitudinal specimen.

during cyclic loading [15, 30, 24, 16]. The present study focused on the effects of microcracking
on the mechanical behavior of compact bovine specimens subjected to uniaxial tensile loading.
For longitudinal specimens, our results agreed well with the literature [23, 21]. That is, specimen
post-yield behavior was affected by the extent of microcracking and the specimen subjected to the
low strain rate loading showed more microcracking when compared to the specimens subjected to

75

8

b)

c)

d)

0.0012s-1

0.0003s-1

a)

`
Fig. 7 Representative SEM post-failure fracture surface images of two typical longitudinal specimens for the strain
Figure
4.9evaluation:
Representative SEM post-failure fracture surface images of two typical longitudinal
rate
effects
specimens
for
the strain
effectsunder
evaluation:
(a), of
(b)0.0003
osteon
(a), (b) osteon pull-out
of therate
specimen
a strain rate
s-1.pullout of the specimen under a
−1
strain
of 0.0003s
; (c),
(d) osteon-pull
outrate
of of
the0.0012
specimen
(c),
(d) rate
osteon-pull
out of the
specimen
under a strain
s-1. under a strain rate of 0.0012
−1
s .osteon
Osteon
is comprised
of the
Haversian
(pointed
arrows)
The
is comprised
of the
Haversian
canalcanal
(pointed
out byout
theby
white
arrowsurrounded
surroundedby
bycylindrical
cylindrical lamellae.
lamellae.
The
close
view
images
show
that
lower
strain
rate
specimen
exhibit
rougher
The close view images show that lower strain rate specimen exhibit rougher fracture surfacefracture
(the ridge rings
surface out
(theby
ridge
out by
arrows (b))toasthecompared
to therate
higher
strainwith
rate aspecimen
pointed
the rings
white pointed
arrows (b))
as compared
higher strain
specimen
smoother tearing
with a smoother tearing surface (see arrow in (d)).
surface (arrow in (d)).

high strain rate loading. The light illumination of microcracking showed that the microcracking
tended to be more coalesced in specimens subjected to a high strain rate. The steeper slope of
the microcracking illumination suggested a higher propagation speed for the high strain rate specimens. Type I collagen, which is the major organic composition of bone lamellae, contributes to the
viscoelastic behavior of cortical bone. With an increase in the loading rate, the collagen tends to
be more brittle [27], which is hypothesized to cause the changing in the elastic behavior of bone,
9

76

leading to a decrease in the energy absorption capacity. Microcracking has been shown to only
affect the post-yield behavior of bone. And with an increase in the loading rate, the point-yield
strain decreased. Thus, the effects of the post yield behavior on the overall mechanical response
of bone tended to decrease. This suggests that when tested under a high strain rate, the effects
of type I collagen on the elastic part of the stress-strain curve seemed to be the critical factor to
determine bone’s failure behavior when compared to the microcracking effects of the mineralized
bone tissue.
Although it has been found that diffuse microcracks may develop at the submicron level [31],
the hierarchical structure plays a role in determining the bone behavior [32]. At the micro structural level, compact bone can be considered as a bone matrix with osteon fiber distributed along the
longitudinal direction. Osteon resulting from bone remodeling is comprised of Haversian channels
surrounded by concentric lamellae. Osteon has less mineralization and a relatively weaker connection with the surrounding bone tissue when compared to the primary bone deposited as the bone
grows. The post-failure microstructure examination by SEM indicated that the specimens with a
higher area percentage of primary bone exhibited the ability to develop more dispersed microcracking than the bone with a higher area percentage of osteons. These results again agreed with
similar trends in the literature [7, 8].
Bone is a natural composite where osteons, growing along the loading axis, serve as fibers to
enhance the mechanical properties of bone [33, 34]. The orientation effects on the compression
mechanical properties of compact bone have been studied by some researchers [5, 9, 27, 35, 36].
It was found that compact bone was stronger in the longitudinal direction than in the transverse
direction and that the ultimate strength and the Young’s modulus of the longitudinally loading specimens were higher than the specimens loaded in the transverse direction. Reilly [37] conducted
tensile, compression, and torsion tests on compact human and bovine bone, where the Young’s
77

modulus and ultimate stress were 22 GPa and 133 MPa for the longitudinal specimen, and 11.3
GPa and 51 MPa for the transverse specimen. The corresponding measurements from the present
study were 21.8 GPa and 93.5 MPa, and 10.9 GPa and 36.2 MPa, respectively. The correlation
of the stress-strain curves and the light illumination of the microcracking showed that microcracking tended to toughen the longitudinal specimens and soften the transverse specimens. The
instant jump of the light intensity to the peak plateau (Fig. 4.4c, d) after yielding indicates that
many microcracks were induced on the transverse specimens suddenly. The post-failure images by
SEM further confirmed the microcracking illumination results by showing that comparatively large
numbers of microcracks with lengths up to hundreds of microns were developed on the transverse
specimens. The post-failure fracture surfaces of the two-directional specimens clearly showed the
different failure pattern of the two directions. The transverse specimens showed much smoother
fracture surfaces when compared to the longitudinal specimens, which were characterized by osteon pullout. This can be explained by the pullout along the osteon direction (the fiber direction
for a composite material) and would be harder than pulling by the cement line. Transverse specimens showed that many microcracks instantly initiated in the cement lines, which resulted in the
stiffness loss in the specimen. The present results show that, to some extent, microcracking serves
as a toughening mechanism; but if it is over a critical value, microcracking weakens the specimen
leading to failure. It is pertinent to further investigate the toughening mechanism of microcracking
in terms of density and length from the aspect of fracture mechanics.

78

loaded in the transverse direction. The current technique has shown potential in relating damage

4.2
Microcracking Morphology of Compact Bovine Bone Subaccumulation real time in bone samples subjected to tensile loading conditions. This information
will bejected
extremely
in relating
the role of micro-damage accumulation in initiating failure
tohelpful
Impact
Loading
and remodeling in bone. The work [61] has been published in the Journal of Biomechanics.

4.2.1

Non-notched Beam Specimens of Compact Bovine Femur Subjected

3.2.1.2. Damage Accumulation and Effects on the Strength and Toughness of Compact

to Charpy Impact Loading

Bovine Bone Subjected to Dynamic Loading
4.2.1.1
Material
and Methods morphology of compact bone subjected to impact loading, nonTo examine
the microcracking

notched
beam specimens
of compact
femur were
loaded
with
instrumented
Charpy
All
the specimens
were machined
out ofbovine
fresh compact
bovine
femurs
(60anmonths
old) provided
under various
impact energy
levels.
Thefemur
middle
sections
the tested specimens
byimpact
a localsystem
beef abattoir.
The mid-diaphysis
of the
bovine
was
cut intoofrectangular
samples.
were the
examined
the by using
basic
technique[62]
to buffered
investigate
the
During
entire machining
process,the
the Frost’s
specimens
werefusion
kept wet
with phosphate
saline
microcracking
morphology.
Thethe
damage
mechanism
and energy
were
(PBS)
solution. After
machining,
specimens
were covered
withdissipation
gauze dipped
inalso
PBSanalyzed.
solution
and preserved in a -20◦C freezer until testing. Specimens were thawed prior to testing.
3.2.1.2.1. Material and Method

c

b
a

Thickness
6mm

d

Figure
Bone
specimens:
A diagramshows
showsthe
theorientation
orientation of
of bone
Fig. 84.10
Bone
specimens:
a, aa, diagram
bone (osteon
(osteon grows
growsalong
along
bone
axis);
b,
Dimensions
of
bone
specimen;
c,
A
typical
specimen
before
testing;
d,
A
typical
bone axis); b, dimensions of bone specimen; c, a typical specimen before testing; d, a typical
specimen
after
testing.
specimen
after
testing.

AllThe
thedimensions
specimens were
machined
out of
compact
(60 months
old) provided
of bone
specimens
arefresh
shown
in Fig.bovine
4.10. femurs
To examine
the microcracking
by a local subjected
beef abattoir.
Theloading,
mid-diaphysis
of the
bovine
femur ofwas
cut into
rectangular
morphology
to impact
non-notched
beam
specimens
compact
bovine
femur
were loaded with an instrumented Charpy impact system under different impact energy levels.

25
79

A tabletop Charpy impact system with light impact mass (0.6 kg) was designed to perform the
low energy (ranging from 0.1 to 2.3 J) impact tests. The Charpy impact system is shown in Fig.
4.11, 4.12. The load cell, projectile, and hammer weighed 0.6 kg. When the hammer system was
dropped at different heights, various impact energy was imposed on the tested specimens. The
projectile was perpendicular to the specimen at the moment of impact.

Swing bar

Hammer

Load cell Projectile Upright 3-Point bending fixture

vinitial

Oscilloscope

Figure 4.11 Tabletop Charpy impact system.
= 2 gL cos θ

The length of the swing bar, L, was 0.39 mm. The initial velocity was calculated as
Load cell:
Mass of Hammer + Load Cell :
1.004Voltage = 125 µε
m = 0.60 Kg
vinitial = 2gL cos θ Force = 7.55ε
Legth of Bar:L = 0.39m
a=F /m
CalculateVinitial = 2 gL cos θ
v = vbar
, s =upright.
si + fvdt
i + and
∫ adtthe
where θθ was
theinitial
initialangle
angle between the center lines of swing
is the

(4.1)

⇒ E = ∫ Fds
The load cell was calibrated before the testing and it yields relationships
between voltage, strain

(ε), force (F), and strain (ε) as

Vinitial = 2 gL cos θ
1.004Voltage = 125µε
80

(4.2)

Table-on Charpy Impact Apparatus

Three Point Bending Fixture

Table-on Charpy Impact Apparatus

Three Point Bending Fixture

Figure 4.12 Schematic representative of the table-on Charpy impact system.

F = 7.55ε

(4.3)

Thus, the velocity (v) and displacement (s) of the projectile at any time during the impact event
can be calculated as

a = F/m

(4.4)

v = vinitial +

adt

(4.5)

s = sinitial +

vdt

(4.6)

where m is the total mass (0.6 kg) of the projectile, hammer, and load cell. Thus, impact energy
(E) can be calculated
E=

Fds
81

(4.7)

2

1

4

9

8

5

6

7

6mm

3

8mm

Figure 4.13 The locations and orders of the viewpoints examined and a sample of microscope
image.
The impact energy
set up+to
be Cell
0.14,: 0.31, 0.41, 0.53, and 0.67 J respectively. Four
Mass ofwas
Hammer
Load

m = 0.60 Kg
Loadspecimens
cell:
specimens were tested in each energy level. After the tests, the tested
were soaked
Legth of Bar:L = 0.39m
1.004Voltage = 125 µε
in the basic fusion
solution
according
basic fusion technique [22]. After bulk stain,
Calculate
Vinitial
= 2 gLto
cosFrost’s
θ
Force = 7.55ε

θ is the initial angle
a =the
F /fracturing
m
a thin section (thickness, 3 mm) was cut off from the end closest to
surface by

v = vi + ∫ adt , s = si + fvdt
a band saw. The surface interested was the cross section of the beam specimen
perpendicular
⇒ E = ∫ Fds
to the osteon growing direction. The thin sections were first polished by a Leco wafering saw
with a precision of 0.01 mm and then followed by final polishing by a Leco lapcloth polishing

Vinitial = 2 gL cos θ
machine with precision
of 1 µm. The final polished ground sections were 8x6x0.3 mm. Then
microcracks and microstructures of each of the ground sections were visualized and photographed
by an Olympus BHS microscope with a Spot digital camera (Diagnostics Inc). Microscopic images
82

Table 4.2 Summaries of all the testing results.

Run No.
E2_Run1
E2_Run2
E2_Run3
E2_Run4
E3_Run1
E3_Run2
E3_Run3
E3_Run4
E4_Run1
E4_Run2
E4_Run3
E4_Run4
E5_Run1
E5_Run2
E5_Run3
E5_Run4
E6_Run1
E6_Run2
E6_Run3

Impact Energy
Energy
Applied (J) Dissipated (J)
0.14
0.082
0.14
0.090
0.14
0.109
0.14
0.088
0.31
0.284
0.31
0.287
0.31
0.201
0.31
0.306
0.41
0.393
0.41
0.248
0.41
0.402
0.41
0.390
0.53
0.487
0.53
0.489
0.53
0.517
0.53
0.511
0.67
0.460
0.67
0.624
0.67
0.507

Area Percentage Occupied Cracks
Sample
by Secondary Osteon %
No.
Status After Testing
100
0
No. Broken
31.85
0
No. Broken
26.78
0
No. Broken
100
0
No. Broken
95
0
No. Broken
90
0
No. Broken
20
1.36
Broken
100
0
No. Broken
15
0
No. Broken
5
0
Broken
100
0
No. Broken
100
0
No. Broken
8
0
Broken
20
0
No. Broken
100
0
No. Broken
5
0
No. Broken
5
1.1
Broken
25
0
No. Broken
5
0.77
Broken

at different magnifications (the maximum field of view was 2.25 x 1.6 mm at magnification x5)
0.67
0.6
Broken
were taken for
each ground0.120942067
section. The80locations and orders of the
viewpoints
examined and a

sample of microscope image are shown in Fig. 4.13.

4.2.1.2

Results and Discussion

The testing results are summarized in Table 4.2. The relationships of impact energy with time and
impact force with time are shown in Fig. 4.14, 4.15 respectively. It was found that no specimens
fractured at the impact energy of 0.14 J; one specimen fractured at each of the impact energies of
0.31, 0.41, and 0.53 J; and three specimens fractured at the impact energy of 0.67 J. One of the
specimens which fractured at the impact energy of 0.67 J possessed pre-existing defects and was
therefore discarded and omitted from Fig. 4.14, 4.15, and Table 4.2).
Two representative post-failure microscopic images at each impact energy level are shown in
83

0.8
E2_Run1
E2_Run2
E2_Run3
E2_Run4
E3_Run1
E3_Run2
E3_Run3
E3_Run4
E4_Run1
E4_Run2
E4_Run3
E4_Run4
E5_Run1
E5_Run2
E5_Run3
E5_Run4
E6_Run1
E6_Run2
E6_Run4

0.7

Energy (J)

1

0.6
0.5
0.4
0.3
0.2
0.1
0
0

50

100

150

200

Time (ms)

Figure 4.14 Impact force history of bovine cortical bone subjected to Charpy impact loading at
various impact energy levels.

1400
E2_Run1
E2_Run2
E2_Run3
E2_Run4
E3_Run1
E3_Run2
E3_Run3
E3_Run4
E4_Run1
E4_Run2
E4_Run3
E4_Run4
E5_Run1
E5_Run2
E5_Run3
E5_Run4
E6_Run1
E6_Run2
E6_Run4

Impact Force (N)

1

1200
1000
800
600
400
200
0
0

50

100
Time (ms)

150

200

Figure 4.15 Impact energy history of bovine cortical bone subjected to Charpy impact loading at
various impact energy levels.
Mass of Hammer + Load Cell :
m = 0.60 Kg
Legth of Bar:L = 0.39m

84
Load cell:
1.004Voltage = 125 µε

E2_Run2
No Broken
32 %

E2_Run1
No Broken
100 %

Mic
E3_Run3
Broken
20 %

Laminar

E3_Run4
No Broken
100 %

Osteon

E4_Run2
Broken
5%

E4_Run3
No Broken
100 %

E5_Run1
Broken
8%

E5_Run3
No Broken
100 %

E6_Run1
Broken
5%

E6_Run2
No Broken
25 %

Mic

Figure 4.16 Representative post failure microscopic images at each impact energy level. The images are labeled with the run case number (E2 Run2 etc), specimen status after impact (broken or
no broken), and the area percentage occupied by the secondary osteon (32%, etc.). Features are
labeled as Mic: microcrack, Osteon: secondary osteon, and Laminar: laminar structure.
85

Fig. 4.16. The correlations of the microcracking numbers to the area percentage occupied by secondary osteon and the impact energy applied are shown in Fig. 4.17. The post-failure examinations
showed that no microcracks were visible when the specimen was not broken. The microcracks in
the fractured specimens were linear and were usually larger than hundreds of micrometers in length
(Fig. 4.16). Fig. 4.17a showed that no specimens fractured or no microcracks were observed when
the area occupied by secondary osteon is larger than 25%. The results showed that the threshold
energy was lager than 0.14 J. When the applied energy is between 0.31 to 0.53 J, the fracture possibility of specimens was random (Fig. 4.17b) and may rely on the microstructures of the sample.
The testing results seemed to conclude that osteon bone can sustain higher impact energy compared
to the laminar bone (Fig. 4.16).
1.6
1.4
1.2
1
0.8
0.6
0.4
0.2
0
0

1.6
1.4
1.2
1
0.8
0.6
0.4
0.2
0

50
100
150
Area Percentage Occupied by
Secondary Osteon %

0

0.2

0.4

0.6

0.8

Applied Impact Energy (J)

Figure 4.17 a, Correlation of mean microcracking numbers to the area percentage occupied by
secondary osteon; b, Correlations of microcracking to the impact energy applied.

86

4.2.2

High Strain Rate Tensile Tests of Compact Bovine Bone by Split Hopkinson Bar(SHPB)

4.2.2.1

Material and Methods

A designed tensile Split Hopkinson Bar (SHPB) was used to conduct the high strain rate tensile
testing (Fig. 4.18). It includes an incident bar (a solid aluminum bar with a big end at the interface
with the momentum trap), a transmitted bar (a hollow aluminum bar), a striker (a hollow aluminum
cylinder sliding along the incident bar), and a striker stopper also called momentum trap (a solid
aluminum bar). The tensile specimen is placed between the incident bar and the transmitted bar
with specially designed threaded steel fixture (Fig.4.19).
The tensile specimens were machined from the fresh, compact, bovine femur following the
same method as the tensile specimens used for static tensile testing in section 4.1 and are of the
same dimensions as the specimen shown in Fig.4.1. The specimens were first glued to the fixture
with Loctite E-30CL glass bounder Epoxy adhesive and then were further fixed with four screws
in each end (Fig. 4.19). Although the specimens were kept moist with gauze soaked with PBS
solution during the preparation of the specimens for the test, the two end of the dog-bone specimen
were in dry conditions for the purpose of curing the epoxy. Thus, the conditions of the specimen
tested were considered as partially wet.
The incident bar was struck by the striker, which was powered by compressed gas. A tensile
wave was created and then propagated through the incident bar toward the specimen. Upon reaching the specimen, part of the pulse was reflected back to the incident bar, and the rest of the pulse
was transferred through the specimen to the transmitted bar. The pulse transferring through the
specimen caused the deformation of the specimen. A full-bridge strain-gauge configuration was
installed on both incident and transmitted bars to capture the wave pulses traveled inside the bar.
87

Figure 4.18 SHPB setup for high strain rate testing (300 900 s−1 ).
Figure 2. Microscope images of
compact bovine bone.

I
nci
dentbarend

Tr
ansmi
ssi
onbarend

Figure 4.19 A typical specimen fixture used in the SHPB testing.
The pulse waves
in bar
the incident
and transmitted
Incident
end
Transmitted
bar bars
end were measured as voltages through an amplifier and an oscilloscope and then were converted to strain. A card-stock paper with a thickness of
0.5 mm was used for pulse shaping to obtain a clean pulse. It has been shown that, by using strong
light scattering effects of crack edges, light illuminations of the microcracks on the specimen surIncident bar end
Transmitted bar end
face were observed on the specimens subjected to quasi-static tensile testing (see section 4.1) [38].
Thus, in order to visualize the fracture behavior of bovine cortical bone subjected to high strain
rate loading, a high speed camera (Phantom v12.1) with a LED panel was also set up to record the
fracture process during the SHPB testing.

88

4.2.2.2

Results and Discussion

A typical incident and transmitted waves measured are shown in Fig. 4.20. It can be seen that only
a small portion of the incident pulse was transmitted to the transmission bar. This was due to the
low impedance in the bone sample and the area mismatch between the incident and transmitted bar
with the dog-bone shaped bone specimens.

Incident

Transmitted

Relected

Figure 4.20 A typical incident and transmitted waves measured in the SHPB testing of bovine bone
dot{\varepsilon }=\frac{{{C}_{0}}}{L}({{\varepsilon }_{I}}(t)-{{\varepsilon
specimens.
}_{R}}(t)-{{\varepsilon }_{T}}(t))
The stress and strain fields are derived from the classic wave theory [39]. Stress in the specimen

σs =

A1 E(εI (t) + εR (t)) + A2 EεT (t)
2As

(4.8)

Strain in the specimen
ε=

C0
L

t
0

[εI (t) − εR (t) − εT (t)]dt

(4.9)

Strain rate
ε˙ =

C0
(εI (t) − εR (t) − εT (t))
L

(4.10)

where E is the Young’s modulus of the incident and the transmitted bar (aluminum), 70 GPa;

89

Incident bar end

Transmitted bar end

εI , εR , εT is the incident, reflected, and transmitted strain field respectively; A1 , A2 , As is cross-

Incidentbar,
bartransmitted
end
Transmitted
bar endrespectively; C0 is the wave speed in the
section area of incident
bar, and sample
incident and transmitted bar, 5092 m/s; and L is the gauge length of the bone sample. A schematic
representation of the interfaces between incident bar with sample and sample with transmission
bar is shown in Fig. 4.21.

Figure 4.21 Schematic representation of interfaces between incident bar with sample and sample
with transmission bar.

The stress and strain relationships obtained from the high stain rate tensile testing are shown
in Fig. 4.22. The correlation of the ultimate stress and the strain at ultimate stress with the corresponding strain rate are shown in Fig. 4.23. The results showed that ultimate stress increased with
the increase of strain rate before a critical strain rate, 600 s−1 . After passing the critical point, the
ultimate stress decreased. The average Young’s modulus measured was 20.3 GPa. The results have
a good agreement with the existing literature [11, 2, 14, 9, 40]. Crowninshiled and Pope [2] examined the response of compact bovine bone in tension at various strain rates (from 0.01 to 200 s−1 ).
Ultimate tensile stress both in longitudinal and transverse increased with the strain rate. Wright
and Hayes [14] conducted uniaxial tensile tests of cortical bovine bone at seven strain rates (from
5.3 x 10−4 to 237 s−1 ). The ultimate strength, modulus of elasticity, and energy absorption capacity were found to increase with strain rates. McElhaney [11] investigated the dynamic response of
both fresh human and bovine femur subjected to dynamic compressive loading at strain rates up
to 1000 s−1 . The Young’s modulus and ultimate stress increased with the strain rates. Katsmanis
90

[40] used the SHPB to measure the Young’s modulus and Poisson’s ratio of human femur cortical
bone in both dynamic compressive and tensile loading. Young’s modulus in dynamic loading was
found to be higher than in static loading, but no difference was found between the Young’s modulus in tension and compression. The studies also showed that scattering of results is significant
due to its highly heterogeneous structure. The results from the literature and the current study are
summarized in Table 4.3.
400
350
Max_strainrate_400 /s

Stress(MPa)

300

Max_straintate_510 /s
Max_strainrate_460 /s

250

Max_strainrate_820 /s
200

Max_strainrate_380 /s
Max_strainrate_290 /s

150
100
50
0
0

0.5

1

1.5

2

2.5

strain%

Figure 4.22 Stress-strain relationships of bovine cortical bone subjected to high strain rate tensile
loading by SHPB setup.

The maximum frame rate of the Phantom high camera is 1,000,000 FPS with a resolution of
128 x 8 pixels . The wave speed can be roughly calculated as

E/ρ. The average density of

cortical bone is about 1900 kg/m3 and the Young’s modulus is about 21.5 GPa. Thus, the wave
speed in cortical bone is calculated to be 3364 m/s. The length of the dog-bone shaped tensile
specimen is 60 mm. Thus, it takes about 17.8 µs for the wave to pass through the specimen. It
has been demonstrated that when the raising time of the incident pulse is long enough to allow the
91

Max Strain rate

400

510

460

380

1.5

1

200
0.5

100

Max stress
Strain% at Max Stress

0
0

200

400

600

Strain rate

800

Ultimate stress(MPa)

300

290

476.67

400
Strain% at Max Stress

Ultimate stress(MPa)

400

820

300
200
100

0

y = -0.0027x^2+3.1628x-601.11
R^2 = 0.6674

0

1000

0

( s 1 )

200

400

600

Strain rate

800

1000

1

(s )

Figure 4.23 a) Correlation of ultimate stress and strain at ultimate stress to the strain rate; b)The
polynomial regression model of ultimate stresses at different strain rates.
Table 4.3 Comparison of mechanical properties of cortical bone subjected to high strain rate loading among different investigators.
Strain rate s

Specimen

Loading Type

Specimen
condition

Young's
Ultimate
Failure Energy-absortion
Modulus (GPa) Stress (MPa) strain% capacity (106J/m3)

290-820
368-795
1000 m/s
100
0.00053-237
0.00167-250
0.001-1500

Bovine Femur
Bovine Femur
Bovine Shaft
Human Femur
Bovine Femur
Bovine Tibia
Bovine Femur

Tensile
Compressive
Tensile
Tensile
Tensile
Tensile
Compressive

PartiallyWet
Dry
Dry
Dry
Wet
Wet
Wet

10.5-36.8
9.3-19.4
10-21.7
18-23.3
17.7-40.4
11-12
12.4-48.3

-1

Present results
F. Ferreira (2006)
M.Pithious (2004)
F. Katsamanis (1990)
T.M. Wright (1976)
R.D. Crowninshield (1974)
J.H. McElhaney(1966)

125-348
223-310
33-50
NR
99-271
110-250
133-400

0.3-1.2
1.5-3.2
0.2-0.4
NR
NR
3-6
0.7-2.9

NR
NR
NR
NR
0.4-1.3
3.8-6.0
2.2-4.8

wave to reflect back and forth three to four times inside the sample, force equilibrium could be
achieved within the initial time [41, 42, 43, 44]. Thus, in order to capture the deformation process
of the bone sample during the high strain rate tensile testing, the camera has to be able to record
a video with a length of about 50-80 µs (about three to four times the time required for a wave to
pass through the specimen) and with a reasonable frame rate. For the present testing, the horizontal
resolution of the high speed image was set to be 500-256 pixels and vertical resolution was 100-60
pixels. This setup yielded a frame rate of about 240,000 FPS. Typical stress-strain curves and their
corresponding high speed images are shown in Fig. 4.24.
As shown in Fig. 4.25, it took about 50 µs to fracture the sample(from 0 to 0.003ε). The time
interval was 8 µs from 0.002ε to 0.003ε. With the current high speed camera used, there was only
92

0.002ε

0.003ε

0.004ε

Figure 4.24 Typical stress-strain curves and their corresponding high speed images.

93

Strain (mm/mm)

0.005
0.004
0.003
0.002
0.001
0

0

10

20

30

40

50

60

Time (μm)

Figure 4.25 The time history of strains of the specimen showed in 4.24.
one image captured at each of these two strains (Fig. 4.24)and no image captured between the two
points. Macro-cracking was visualized in the high speed image at 0.003ε and there is no microor macro-cracking observed in the high speed image at 0.002ε. There is no micro-cracking being
visualized for all the tested samples and the failure behavior followed the similar trend as the one
shown in Fig. 4.24, 4.25. It implied that the microdamage or microcracking once initiated would
coalesce and extend rapidly (less than several microseconds) to cause the catastrophic failure of the
sample. Thus, compared to the ductile failure behavior exhibited during quasi-static tensile tensile,
the failure of the bovine cortical bone subjected to high strain rate tensile testing is governed by
brittle failure. The observations agree with the existing literature [2, 14, 9].

4.2.3

Conclusion

Damage accumulation in compact bovine femur subjected to uniaxial tensile loading was examined
by strong light illumination effects of microcracking. Imaging was done using a high-speed camera
capturing images at 200 to 1500 FPS. The tensile tests were performed in a multipurpose tensile
testing system with cross-head speeds ranging from 0.5 to 10 mm/min which leads to strain rates

94

of 0.0001 to 0.0012 s−1 (physiologically relevant to walking and running [10]). The post-failure
images were examined in a scanning electron microscope (SEM); and the effects of microstructure,
strain rate, and orientation were evaluated. Correlation of the high-speed images with stress-strain
curves indicated that optically visible microcracks were most likely initiated at yielding, and the
specimens with dispersed microcracks exhibited a higher energy absorption capacity compared
to the specimens with coalesced local cracks. It was found that damage accumulation negatively
correlates to strain rate and that transverse specimens exhibited a different failure pattern compared to the longitudinal specimens. Strain hardening and softening were found in the longitudinal
and transverse specimens respectively. The microcracking in the transverse specimens instantly
increased to a peak after yielding compared to the gradual growth until failure in the longitudinal specimens. The average Young’s modulus (21.5 GPa) and ultimate stress (93.5 MPa) of the
specimens loaded in the longitudinal direction were more than twice that of the specimens (10.9
GPa and 36.2 MPa respectively) loaded in the transverse direction. The current technique has
shown potential in relating real time damage accumulation in bone samples subjected to the static
tensile loading condition. This information will be helpful in relating the role of micro-damage
accumulation in initiating failure and remodeling in bone.
The microcracking morphology of compact bovine femur subjected to impact loading was investigated by an instrumented Charpy impact system under different impact energy levels. The
post failure examinations were performed by a bright-field microscope. The microcracks were
counted for each of the samples and were correlated to its microstrutures. The results showed that
no samples were broken under the impact energy of 0.14 to 0.67 J when the area percentage occupied by the secondary osteon was higher than 25 %. The microstructures of the broken specimens
featured as laminar structures. The threshold energy for the non-notched rectangular specimens
was larger than 0.14 J. The experiments showed that osteon bone can sustain higher impact energy
95

compared to the laminar bone.
The microcracking morphology and high strain rate failure behavior of compact bovine femur
subjected to high strain rate tensile loading were studied through a designed tensile SHPB system.
The ultimate stress and strain at ultimate stress were measured to be in ranges of 125-348 MPa
and 0.3-1.2 mm/mm respectively with a strain rate of 290-820 s−1 . The Young’s modulus was in
the range of 10.5 to 36.8 GPa. The ultimate stress was found to increase with strain rate. The
fracture processes were recorded through a high speed camera with a frame rate of about 240,000
FPS. The high-speed images were synchronized to the stress-strain curves. Unlike the static tensile
testing, no optically visible microcracks were captured during the high strain rate tensile testing.
The results indicted that the failure of the cortical bovine bone was governed by brittle failure at
a high strain rate. The microdamage initiated and coalesced in a short period of time (less than
several microseconds) to cause the catastrophic failure of the sample.
In conclusion, the experiments showed that the microcracking toughening mechanism is the
major toughening mechanism of cortical bone in low strain rate loading. In an impact or a high
strain rate, the failure behavior depends more on its micro-structure. In the next chapter, the correlation between the microstructures and the macro performance of osteon bone is studied through
computational methods.

96

REFERENCES

97

REFERENCES
[1] Dennis R Carter, William E Caler, Dan M Spengler, and Victor H Frankel. Fatigue behavior
of adult cortical bone: the influence of mean strain and strain range. Acta Orthopaedica,
52(5):481–490, 1981.
[2] RD Crowninshield and MH Pope. The response of compact bone in tension at various strain
rates. Annals of Biomedical Engineering, 2(2):217–225, 1974.
[3] JD Currey and K Brear. Tensile yield in bone. Calcified tissue research, 15(1):173–179,
1974.
[4] Charles H Turner. Bone strength: current concepts. Annals of the New York Academy of
Sciences, 1068(1):429–446, 2006.
[5] Raghavendra R Adharapurapu, Fengchun Jiang, and Kenneth S Vecchio. Dynamic fracture
of bovine bone. Materials Science and Engineering: C, 26(8):1325–1332, 2006.
[6] Martine Pithioux, Damien Subit, and Patrick Chabrand. Comparison of compact bone failure
under two different loading rates: experimental and modelling approaches. Medical engineering & physics, 26(8):647–653, 2004.
[7] JD Currey. The effects of strain rate, reconstruction and mineral content on some mechanical
properties of bovine bone. Journal of Biomechanics, 8(1):81–86, 1975.
[8] S Saha and WC Hayes. Tensile impact properties of human compact bone. Journal of biomechanics, 9(4):243–251, 1976.
[9] F Ferreira, MA Vaz, and JA Simoes. Mechanical properties of bovine cortical bone at high
strain rate. Materials characterization, 57(2):71–79, 2006.
[10] Ulrich Hansen, Peter Zioupos, Rebecca Simpson, John D Currey, and David Hynd. The effect
of strain rate on the mechanical properties of human cortical bone. Journal of biomechanical
engineering, 130(1):011011, 2008.
[11] James H McElhaney et al. Dynamic response of bone and muscle tissue. J. appl. Physiol,
21(4):1231–1236, 1966.
[12] Diane Margel Robertson and Dennis C Smith. Compressive strength of mandibular bone as a
function of microstructure and strain rate. Journal of biomechanics, 11(10):455–471, 1978.
98

[13] S Saha and WC Hayes. Relations between tensile impact properties and microstructure of
compact bone. Calcified tissue research, 24(1):65–72, 1977.
[14] TM Wright and WC Hayes. Tensile testing of bone over a wide range of strain rates: effects of
strain rate, microstructure and density. Medical and biological engineering, 14(6):671–680,
1976.
[15] David B Burr, Mark R Forwood, David P Fyhrie, R Bruce Martin, Mitchell B Schaffler,
and Charles H Turner. Bone microdamage and skeletal fragility in osteoporotic and stress
fractures. Journal of Bone and Mineral Research, 12(1):6–15, 1997.
[16] MB Schaffler, EL Radin, and DB Burr. Long-term fatigue behavior of compact bone at low
strain magnitude and rate. Bone, 11(5):321–326, 1990.
[17] D Vashishth, JC Behiri, and W Bonfield. Crack growth resistance in cortical bone: concept
of microcrack toughening. Journal of Biomechanics, 30(8):763–769, 1997.
[18] D Vashishth, KE Tanner, and W Bonfield. Contribution, development and morphology of microcracking in cortical bone during crack propagation. Journal of Biomechanics, 33(9):1169–
1174, 2000.
[19] D Vashishth, KE Tanner, and W Bonfield. Experimental validation of a microcracking-based
toughening mechanism for cortical bone. Journal of Biomechanics, 36(1):121–124, 2003.
[20] Fergal J OBrien, David Taylor, and T Clive Lee. Microcrack accumulation at different intervals during fatigue testing of compact bone. Journal of Biomechanics, 36(7):973–980,
2003.
[21] Peter Zioupos, Ulrich Hansen, and John D Currey. Microcracking damage and the fracture
process in relation to strain rate in human cortical bone tensile failure. Journal of biomechanics, 41(14):2932–2939, 2008.
[22] HM Frost. Presence of microscopic cracks in vivo in bone. Henry Ford Hosp Med Bull,
8(2):35, 1960.
[23] P Zioupos, JD Currey, and AJ Sedman. An examination of the micromechanics of failure of
bone and antler by acoustic emission tests and laser scanning confocal microscopy. Medical
engineering & physics, 16(3):203–212, 1994.
[24] Gwendolen C Reilly and John D Currey. The development of microcracking and failure in
bone depends on the loading mode to which it is adapted. Journal of experimental biology,
202(5):543–552, 1999.
99

[25] Ravi K Nalla, Jamie J Kruzic, and Robert O Ritchie. On the origin of the toughness of
mineralized tissue: microcracking or crack bridging? Bone, 34(5):790–798, 2004.
[26] RK Nalla, JS St¨olken, JH Kinney, and RO Ritchie. Fracture in human cortical bone: local
fracture criteria and toughening mechanisms. Journal of biomechanics, 38(7):1517–1525,
2005.
[27] Robb M Kulin, Fengchun Jiang, and Kenneth S Vecchio. Loading rate effects on the¡ i¿
r¡/i¿-curve behavior of cortical bone. Acta biomaterialia, 7(2):724–732, 2011.
[28] AG Evans and KT Faber. Crack-growth resistance of microcracking brittle materials. Journal
of the American Ceramic Society, 67(4):255–260, 1984.
[29] Anthony G Evans and Robert M McMeeking. On the toughening of ceramics by strong
reinforcements. Acta Metallurgica, 34(12):2435–2441, 1986.
[30] David B Burr, Charles H Turner, Pratap Naick, Mark R Forwood, Walter Ambrosius, M Sayeed Hasan, and Ramana Pidaparti. Does microdamage accumulation affect the mechanical
properties of bone? Journal of Biomechanics, 31(4):337–345, 1998.
[31] Gagik P Parsamian and Timothy L Norman. Diffuse damage accumulation in the fracture process zone of human cortical bone specimens and its influence on fracture toughness. Journal
of Materials Science: Materials in Medicine, 12(9):779–783, 2001.
[32] JW Ager, G Balooch, and RO Ritchie. Fracture, aging, and disease in bone. Journal of
materials research, 21(08):1878–1892, 2006.
[33] JC Behiri and W Bonfield. Orientation dependence of the fracture mechanics of cortical bone.
Journal of biomechanics, 22(8):863–872, 1989.
[34] Jae-Young Rho, Liisa Kuhn-Spearing, and Peter Zioupos. Mechanical properties and the
hierarchical structure of bone. Medical engineering & physics, 20(2):92–102, 1998.
[35] R Shahar, P Zaslansky, M Barak, AA Friesem, JD Currey, and S Weiner. Anisotropic poisson’s ratio and compression modulus of cortical bone determined by speckle interferometry.
Journal of biomechanics, 40(2):252–264, 2007.
[36] Y Tanabe, S Tanaka, M Sakamoto, T Hara, H Takahashi, and Y Koga. Influence of loading
rate on anisotropy of compact bone. Le Journal de Physique IV, 1(C3):C3–305, 1991.

100

[37] Donald T Reilly and Albert H Burstein. The elastic and ultimate properties of compact bone
tissue. Journal of biomechanics, 8(6):393–405, 1975.
[38] Wei Zhang, Srinivasan Arjun Tekalur, Melissa Baumann, and Laura R McCabe. The effects of damage accumulation on the tensile strength and toughness of compact bovine bone.
Journal of biomechanics, 46(5):964–972, 2013.
[39] Karl F Graff. Wave motion in elastic solids. Courier Corporation, 2012.
[40] Fotios Katsamanis and Demetrios D Raftopoulos. Determination of mechanical properties of
human femoral cortical bone by the hopkinson bar stress technique. Journal of biomechanics,
23(11):1173–1184, 1990.
[41] Kenneth S Vecchio and Fengchun Jiang. Improved pulse shaping to achieve constant strain
rate and stress equilibrium in split-hopkinson pressure bar testing. Metallurgical and materials transactions A, 38(11):2655–2665, 2007.
[42] W Chen, B Zhang, and MJ Forrestal. A split hopkinson bar technique for low-impedance
materials. Experimental Mechanics, 39(2):81–85, 1999.
[43] Guruswami Ravichandran and Ghatuparthi Subhash. Critical appraisal of limiting strain rates
for compression testing of ceramics in a split hopkinson pressure bar. Journal of the American
Ceramic Society, 77(1):263–267, 1994.
[44] Lal Ninan, J Tsai, and CT Sun. Use of split hopkinson pressure bar for testing off-axis
composites. International Journal of Impact Engineering, 25(3):291–313, 2001.

101

Chapter 5
Investigating the Fracture and Toughening
Mechanism of Cortical Bone Using XFEM
and Multiscale Modeling Methods

5.1

Introduction

Osteon (also called secondary osteon) is the major bone unit of cortical bone in the microstructure level (from 10-500 µm) [1]. The primary constituent of osteon bone is osteon, which grows
along the loading axis and comprises a center Haversian system surrounded by cylindrical lamellae formed by collagen fibrils and embedded mineral in nanostructure level (from a few hundred
nanometers to 1 µm). Osteon bone resembles the fiber-reinforced composite with a higher Young’s
modulus and tensile strength along the osteonal direction or the loading axis, and lower ones in the
transverse direction which is perpendicular to the loading axis. The boundary between osteons and
surrounding bone tissue (interstitial bone), called the cement line, exhibits different morphology
and composition compared to its surrounding osteon [2] and is hypothesized to serve as a barrier
or ductile interface between osteons and interstitial bone. Although extensive studies have been
focusing on the fracture and toughening mechanism of cortical bone, the correlations between microstructure and microcracking, fracture toughness, and failure mechanism of bone has not yet

102

been fully understood.
It has been shown that crack deflection, microcracking, and uncracked ligament bridging are
the major toughening mechanisms of cortical bone [3, 4, 5, 6, 7]. Micocracking formed ahead of
the crack tip increasing the inherent fracture resistance of bone is the intrinsic toughening mechanism; the uncracked ligaments operated in the wake of the crack tip preventing the two newly
formed crack edges from falling apart is the extrinsic toughening mechanism. The raising R-curve
behavior observed [8, 6, 7] is believed to be due to contributions of the above two toughening
mechanisms. The deflection of cracks is potentially due to the interference of cement lines, which
prevent cracks from penetrating through osteons to lead to a catastrophic failure. A cement line is
shown to be of lower mineralization compared to its surrounded osteon [9]. Thus, this soft interface
separates osteons from interstitial bone to prevent microcracks from propagating into osteons. The
above hypothesis supports the claim that osteon bone has a higher fracture toughness compared to
primary bone, although its Young’s modulus and strength is lower than the primary bone [10].
Limited in the ability to capture moving microcracks, differentiate the roles of bone constituents, and bridge macro- and micro-scale mechanics, experimental methods alone seem unable
to provide full insight into the fracture and toughening mechanism of bone. Thus, finite element
modeling (FEM) has been widely used in bone fracture study. Cohesive element with linear fracture mechanics method (LFEM) [11] is one of the major techniques used to study the interaction of
microcracking and microstructure. Although the cohesive element method is capable of capturing
interface delaminating behaviors, the cracking path is specified and is mesh dependent. Thus, it
may not be applicable in the case of path independent microcracking.
A newly developed FEM technique, XFEM, solved the singular problems at the crack tip by
adding a discontinuous jump function across the crack surfaces to the traditional shape function
[12, 13, 14]. XFEM does not require conforming mesh to the geometry of discontinuities and
103

cracks can propagate along an arbitrary path. Thus, XFEM is suitable for modeling fracture of
cortical bone.
In the present dissertation, the XFEM technique for discontinuities is used to investigate the
interaction between microcracking, fracture, microstructure, and its mechanical properties. The
fracture toughness of a compact tensile specimen (CTS) of cortical bone is simulated. The simulations studied the significant role of cement lines in cortical bone fracture mechanism by calculating
fracture toughness at various modulus ratios, various fracture energy ratios, and initial damage ratios of cement line and osteon. The osteon size effect on the fracture toughness is also investigated.
To study the role of density and size of osteons and the properties of each constituent on the tensile
strength of cortical bone, a multiple-scale modeling approach [15] based on the classical homogenization scheme was used to investigate the tensile behavior of cortical bone subjected to uniaxial
tensile loading. To evaluate the orientation effects, the simulated compact bovine was loaded in the
longitudinal and transverse directions, respectively. Lastly, a semi-concurrent modeling scheme for
impact testing simulations was developed and implemented in ABAQUS/EXPLICIT 6.11 to study
the role of cement lines on the impact behavior of bone.

5.2
5.2.1

Methods
Multi-Scale Modeling of Uniaxial Tensile Testing

The multiscale modeling approach proposed by Kouznetsova [15] bridges the microstructure to
the macroscale behavior by the classic computational homogenization scheme [16, 17]. The deformation gradient from the macro-model is imposed on the microstructural representative volume
element (RVE) to form a microscopic boundary value problem. The results obtained by solving
the microscopic boundary problem are the variable input to compute the macroscopic output.
104

The scheme of macro-micro-macro transition relies on three hypotheses: volume averaging
theorem, separation of microscale and macroscale, and the global periodicity. The macro-micro
couplings including deformation gradient coupling, stress coupling, and internal energy coupling
are implemented by the volume averaging theorem. The deformation gradient coupling [15] can
be determined from the following equations:
The stress couplings

FM =

PM =

1
V0

1
V0

Fm dV0

(5.1)

Pm dV0

(5.2)

V0

V0

1
PM • FMT
det(FM )

(5.3)

Pm : δ FmT dV0 = PM : δ FmT

(5.4)

σM =
The internal work coupling

1
V0

where V0 is the initial RVE volume; FM is the macroscopic deformation gradient tensor; Fm is
the microscopic deformation gradient tensor; PM is the macroscopic first piola-Kirchhoff stress
tensor; Pm is the microscopic first piola-Kirchhoff stress tensor; σM is the macroscopic Cauchy
stress tensor; and σm is the microscopic Cauchy stress tensor.
The separation of two scales requires that the dimension of RVE must be much smaller than the
dimension of the macro-model so that the gradient is constant for RVE. Based on the continuity
condition of composite materials, a periodic boundary condition is assumed [15].
105

continuity condition of composite materials, a periodic boundary condition is assumed (Gougeon
et al., 2003b; Marechal et al., 2004).

Figure
Schematic picture
thethe
corner
nodes
numbering
for afor
RVE
in 3D in
view
Figure
5.11 Schematic
pictureofof
corner
nodes
numbering
a RVE
3D view.

Thecomputational
computational homogenization
multiscale
problem
can can
be solved
numerically
using FE
The
homogenization
multiscale
problem
be solved
numerically
using the

FEM
method.
Fig.
5.1 shows
a schematic
3D with
RVEthe
with
the numbers
of the 8 corner
nodes. The
method.
Figure
1 shows
a schematic
3D RVE
8 corner
nodes numbering.
The corner
corner
the master
to preserve
the periodicity
nodesnodes
of 1, of
2, 1,
4, 2,
5, 4,
are5,theare
master
nodes. nodes.
In orderInto order
preserve
the periodicity
of RVE at of
theRVE at
thedeformed
deformedstate,
state,
the displacements of the master corner nodes and the surface pairs (top and
the displacements of the master corner nodes and the surface pairs (top and
bottom, right and left, and front and back) should be correlated as follows [15]:

bottom, right and left, and front and back) should be correlated as following (Marechal et al.,
7

utop = ubottom + u5 − u1

(5.5)

uright = ule f t + u2 − u1

(5.6)

u f ront = uback + u1 − u4

(5.7)

where u1 , u2 , u4 , u5 are the initial displacement vectors of the four master corner nodes; the utop
and ubottom , uright and ule f t , u f ront and uback are the initial displacement vectors for the perspective
pair of nodes on the three pair of surfaces.
Thus, the above applied constraints implied that the total deformation of the RVE can be imposed by the displacement of the four corner nodes. PM and σ M can be computed by [15]

106

PM =

1
V0

f1e x01 + f2e x02 + f4e x04 + f5e x05 =

1
V

f1e x1 + f2e x2 + f4e x4 + f5e x5 =

σM =

1
V

1
V0

∑

f pe x p

(5.8)

p=1,2,4,5

∑

f pe x p

(5.9)

p=1,2,4,5

where f p is the reaction external forces acted on four corner nodes, 1,2,4,5. Tangent stiffness tensor
can be computed as [15]

4 τ
CM

=

LT
1
ij
( x (i) KM x ( j) ) : 4 I S + 4 I • τM
∑
∑
V0 i j

1
τ
δ τM = 4CM
: Dδ M τM = det(FM )σM Dδ M = (Lδ M + LδTM )
2

(5.10)
(5.11)

ij

where KM is the component of RVE reduced stiffness matrix,KM ,4 I S , and 4 I are the symmetric
fourth-order unit tensor and fourth-order unit tensor, respectively.
In the microscopic level, osteon with a central canal (Harvension system), cement line, and
interstitial bone (matrix) are the three major continents of osteon bone. Fig. 5.2a shows a typical
microscopic image of a bovine cortical bone. A microstructural RVE for osteon bone chosen to
represent the microstructure of cortical bone in the multiscale modeling is shown in Fig. 5.2b.
The RVE was cubic with a dimension of 500x500x500 mum. The diameters of the central canals
and width of the cement lines were kept as constant values, 35 mum and 5 mum respectively. The
diameters of osteons were selected in the range of 100 to 200 mum [14]. The various microstructure
configurations for the RVE are summarized in Table 5.1. To study the osteon size effects, the area
percentage occupied by the secondary osteon was kept constant and the diameters of the osteons
were varied (Cases 1 to 4 (longitudinal direction) and Cases 5 to 8 (transverse direction) in Table
107

Areathe diameter of osteons
5.1). To study the effects ofYoung's
the area percentage occupied by osteons,
Diameter of
Loading
Modulus of Numbers of
percentage
Osteons
Osteon (mm)
Cement
occupied by direction
were kept constant and the area percentage was varied (Case 4&9, longitudinal direction and Case
Line
Osteon
12.66/22.78
5
0.14927
35
Transverse
5
0.14927
Transverse
8&10, transverse direction). All the microstructure
configurations
were35
generated
by a MATLAB
35
Transverse
5
0.14927

code with random distributions of osteons in the matrix. The dimensions and loading conditions of
the macro-model resembled the 2uniaxial tensile
testing by Zhang [18]: a dog bone specimen with
3
11.51/20.71

5.755

NA

12.66/22.78
12.66/22.78
dimensions of 60x8x2 mm was
loaded in
a constant speed of 0.1 mm/s.
11.51
0.17/0.18
0.153/0.162
0.17
5
0.14927
35
Transverse

NA
NA
0.153/0.162
NA
0
NA
0
Transverse

5.755

0.
2mm

Figure 5.2 a) Microstructure of cortical bone: 1) Osteon, 2) Cement line, 3) Interstitial bone, b)
RVE.
Table 5.1 Cortical bone with various microstructures simulated in the uniaxial tensile testing.
Case
No.
1
2
3
4
5
6
7
8
9
10

Numbers of Diameter of
Osteons Osteon (mm)
3
5
7
10
3
5
7
10
5
5

0.19271
0.14927
0.12616
0.10555
0.19271
0.14927
0.12616
0.10555
0.10555
0.10555

Area Percentage
Occupied by
Osteon %
35
35
35
35
35
35
35
35
17.5
17.5

Loading
Direction
Longitudinal
Longitudinal
Longitudinal
Longitudinal
Transverse
Transverse
Transverse
Transverse
Longitudinal
Transverse

1 20.01 3 0.19271
2 19.72 5 0.14927
3 19.44 7 0.12616
4 19.00 0.10555 3
5 10.92 3 0.19271
6 10.60 5 0.14927
7 10.30 7 0.12616
8 9.81 10 0.10555
9 20.25 5 0.10555
10 10.81

The simulations were implemented by the commercial FEA software ABAQUS/Standard 6.11
with two user subroutines coded by FORTRAN, which was compiled by Intel (R) Visual Fortran
108

Compiler Version 11.1.072. The materials of the RVE components were defined by using the
ABAQUS built-in material models. Osteons and interstitial bone (matrix) were treated as transversely isotropic and cement lines were treated as isotropic. The material properties of osteons and
interstitial bone (matrix) used in the simulations are shown in Table 5.2.
Table 5.2 Elastic properties of osteon and matrix.

5.2.2

Elastic Properties

Osteon

Matrix

ET (GPa)
EL (GPa)
νT
νL
GTL (GPa)
GTT (GPa)

11.51
20.71
0.17
0.18
8.4
3.57

12.66
22.78
0.153
0.162
8.4
3.93

Characterizing the Fracture Toughness of Cortical Bone with Various
Microstructures by XFEM

The dimensions and meshes of the C(T) specimen [19] simulated are shown in Fig. 5.3. It included
two parts: Part 2, the area within the crack front, was the concerned zone; Part 1, the area away
from the crack front, was less concerned. Thus, Part 2 was modeled with detailed microstructures
of cortical bone including osteon, cement line, and interstitial bone, while part 1 was modeled
as a bulk bone matrix. Boundary surfaces between Part 1 and Part 2 were tied together. The
microstructures were created using the same method described in the previous section, as well
as the dimensions of the three constituents also followed the same rules in the previous sections.
The loading conditions complied to the ASTM E399 [19] by fixing the three translation degrees
of freedom of one loading hole and moving another loading hole in z direction (Fig. 5.3) with a
constant low speed (1 mm/s).
109

a)

b)
Part1
Part2

Part1

c)
Part2

3

2

d)

1

Part2

3 Meshes
dimensions
C(T)specimen:
specimen: a)
a) Schematic
Schematic drawing
FigureFigure
5.3 Meshes
andand
dimensions
ofof
thetheC(T)
drawingofofC(T)
C(T)specimen
specimen
(Part 2(Part
is the
zonezone
with
detailed
microstructures
bone;Part
Part11isisrepresented
represented
2 isconcerned
the concerned
with
detailed
microstructuresofofosteon
osteon bone;
with bulkwith
bonebulk
matrix
without
showing
detailed
microstructures):
a=70
mm,
W=140
mm, B=0.5
bone matrix without showing detailed microstructures ): a=70mm, W=140mm,
mm; b)
Meshes b)
ofMeshes
the C(T)
specimen;
c) Schematic
drawing
of theof concerned
zone
(Part
2):2):1:
B=0.5mm;
of the
C(T) specimen;
c) Schematic
drawing
the concerned
zone
(Part
osteon, 2: cement line,
3: interstitial
bone,
of Part
1, osteon;
2, cement
line;d)3,Meshes
interstitial
bone;2.d) Meshes of part 2.

The fracture toughness, KQ(MPa·m1/21/2
), is calculated as(Miehe and Apel, 2004):

The fracture toughness, KQ (MPa.m

(

), is calculated as [19]

)

K Q = PQ / BW 1 / 2 • f (a / W )
where:

(7)

KQ = PQ /BW 1/2 • f (a/W )

(2 + a / W )(0.886 + 4.64a / W − 13.32a 2 / W 2 + 14.72a 3 / W 3 − 5.6a 4 / W 4 )
where f (a / w) =
(1 − a / W ) 3 / 2

(5.12)
(8)

PQ is the load (KN) determined from the load-displacement curve following the rule of ASTM
2 /W 2 +(CM)
3 −in
4 /W 3;
4 ) W is the
E399(Miehe
Apel,
2004);
is the specimen
thickness
as3 /W
shown
Figure
(2 +and
a/W
)(0.886
+B
4.64a/W
− 13.32a
14.72a
5.6a

f (a/w) =

(1 − a/W )3/2

(5.13)

specimen width (CM) as shown in Figure 3; a is the crack length (CM) as shown in Figure 3.

110

12

PQ is the load (kN) determined from the load-displacement curve following the rule of ASTM
E399 [19]; B is the specimen thickness (cm) as shown in Fig. 5.3; W is the specimen width (cm)
as shown in Fig. 5.3; a is the crack length (cm) as shown in Fig. 5.3.
The simulations were implemented by the commercial FEA software, ABAQUS/Standard. As
quasi-brittle material, the fracture of cortical bone has been shown to be strain controlled [4]. Thus,
the maximum principal strain failure criterion was chosen for the criterion of damage initiation
and a mixed-mode, energy-based damage evolution law was chosen for the criterion of damage
0
propagation. Initiation occurs when the maximum principal strain reaches a critical value (εmax
)

f=

(εn )
0
εmax

(5.14)

The relevant material data for the damage initiation and propagation of cortical bovine bone
0 = 40MPa, G = 1374 N/m, G
are as follows [14, 11]: σmax
IC
IIC = 4016 N/m, GIIIC = 4016 N/m.

Where GIC is the normal mode fracture energy, GIIC is shear mode fracture energy for failure in
the first shear direction, and GIIC is the shear mode fracture energy for failure in the second shear
direction. Since there is a lack of respective material data for individual constituents, it is assumed
the osteon has the same properties as above, and the matrix takes a 10 percent increase of the
matrix principal stress and a 10 percent decrease in Young’s modulus. By changing the cement
line properties, various configurations are studied to evaluate the interface effects on the fracture
toughness of osteon bone.

111

5.2.3

Semi-Concurrent Multi-scale Modeling Approach for Modeling Impact Testing on Cortical Bone by ABAQUS/EXPLICIT

To investigate the role of cement lines on the impact behavior of bone, a semi-concurrent modeling
scheme for impact testing simulations was developed and implemented in ABAQUS/EXPLICIT
6.11. The stress coupling in this approach is also based on the averaging theorem used in the
computational homogenization scheme [15, 17, 20, 21]. It provides correct results only when the
localization in RVE is absent. Thus, this method cannot predict failure. In this work, a method
was developed to perform a halfway coupling from the macro-model to the micro-model. The
boundary conditions of macro-model were passed to the micro-model through shape functions.
Thus, deformation of the micro-model was concurrent to the macro-model. The analysis result of
the micro-model provided an insightful view of how cement lines and osteons reacted to impact
loading.
ABAQUS/EXPLICIT uses the explicit method [22] to solve dynamic equilibrium equations at
each time increment, which does not enforce equilibrium of the internal forces with the externally
applied loads. It updates the stiffness matrix based on geometry changes and material changes at
the end of each time increment and then uses this stiffness matrix to calculate the next increment of
load (or displacement) applied to the system. Thus, to avoid the drifting of results from the correct
solution, a small time increment should be used in ABAQUS EXPLICIT. ABAQUS/STANDARD
enforces equilibrium at each time increment by Newton-Raphson iterations. Thus, large time increment could be used in ABAQUS/STANDARD. However, it would require a reconstruction and
a tangent stiffness matrix inversion at each iteration. The construction of the stiffness matrix can be
computationally costly. Choosing the correct solver based on the analysis type is important for improving the accuracy and time efficiency of a FEM analysis. ABAQUS/EXPLICIT works better for

112

impact problems and is chosen as the solver for the semi-concurrent approach.The computational
framework of the semi-concurrent multi-scale modeling scheme is shown in Fig. 5.4.

Macro-model

Projectile

UM

σM

Micro-model
Boundary value problem
RVE

RVE

Figure 5.4 Computational framework of the concurrent multi-scale modeling scheme.

The plate simulated in the macro-model was a quarter of a plate with dimensions of 25x25
mm. The two surfaces of the corner under the projectile (Fig. 5.4) were simulated as symmetric
along the x and y axes. One of the remaining two surfaces of the plate was set free, and another
was fixed in the z direction. The plate was meshed with the element size equaling the size of
RVE, 500x500x500 µm. A linear brick element with reduced integration point was chosen as
the element type for both of the macro-model and micro-model. At each time increment, the
boundary conditions of each element in the macro-model were passed to the 8 corner nodes of
a corresponding RVE through ABAQUS subroutines. Then, in the micro-model, the boundary
conditions of all the surface nodes were mapped through the nodal displacements of the 8 corner
nodes of the RVE by the first order brick shape function
113

Macro-model

Projectile

u = (1/8)(1 − g)(1 − h)(1 − r)u1 + (1/8)(1 + g)(1 − h)(1 − r)u2
+(1/8)(1 + g)(1 + h)(1 − r)u3 + (1/8)(1 − g)(1 + h)(1 − r)u4

σM

UM

+(1/8)(1 − g)(1 − h)(1 + r)u5 + (1/8)(1 + g)(1 − h)(1 + r)u6
+(1/8)(1 + g)(1 + h)(1 + r)u7 + (1/8)(1 − g)(1 − h)(1 + r)u8

(5.15)

Micro-model
where u1 , u2 , u3 , u4 , u5 , u6 , u7 , u8 are the displacements
the 8problem
corner nodes passed from the
Boundaryofvalue
RVE
macro-model; g, h, r are the local coordinates inside of RVE (Fig. 5.5).
RVE

Figure 5.5 Schematic picture of RVE with both local and global coordinate systems.
Through the nodal mapping of the isoparametric interpolation, the boundary conditions of each
element in the macro-model were copied to RVE. The advantage of this method is that the periodicity and continuity conditions are persevered without applying the prescribed linearly periodic
boundary conditions for RVE in the traditional computational homogenization scheme [15, 17, 20].
The average stress theorem was applied to the first Piola-Kirchhoff stress tensor as well as the
Cauchy stress tensor. For the macro first Piola-Kirchhoff stress, PM , is given in terms of micro first
Piola-Kirchhoff stress Pm as
114

PM =

1
V0

V0

(5.16)

Pm dV0

where V0 is the initial volume of RVE. After applying the divergence theorem, the PM can be
calculated as [15, 17]
PM =

1
V0

(5.17)

px0 dΓ0
Γ0

where p = n0 · PmT is the first Piola-Kirchhoff stress vector; x0 is the position vector; Γ0 is the
surface of RVE.
After applying the boundary conditions, PM , is further simplified as

PM =

1
V0

f1e x01 + f2e x02 + f3e x03 + · · · · · · · · fne x0n =

1
V0

∑

f pe x p

(5.18)

p=1,2,3.....n

where f pe is the reaction (resulting) external forces acted on the surface nodes of RVE [15, 17]. The
Cauchy stress tensor σM then can be calculated as

σM =

1
PM · FMT
det(FM )

(5.19)

where FM is the deformation gradient at the integration point of element in the macro-model.
The concurrent multiscale framework was implemented through ABAQUS/EXPLICIT through
subroutine Vumat, Vdisp, and ABQMAIN. The flow chart showing the numerical implementation
of the framework is shown in Fig. 5.6.
When Vumat was called at each integration point at every time increment, the nodal displacements (obtained through user subroutine ’Post macro.for’) of the element containing the integration point was passed to the subroutine Vdisp. In the meantime, the deformation gradient of the
integration point was passed to ’Post micro.for’ for calculating Cauchy stress. Then the Cauchy
115

ABAQUS
Macro-model

ABAQUS
Miro-model

macro_impact.inp

micro_structure.inp

*Heading
…..

*Heading
…..

*user material
*boundry, user
Vumat.for
Subroutine Vumat
Vdisp.for

DFGRD1(I,J)=FM
Do K=1, nblock
……..
Execute Post_macro.for
……..
Execute micro-structure
analysis job
……...
Execute Post_macro.for
………
STRESS(nblock,
ndir+nshr)= σM
……....
END DO
RETURN

disp.txt

Subroutine Vdisp
Do K=1, nblock

UM

rval(j, k)=…….

fm.txt

END DO
RETURN

FM
Post_micro.for
Subroutine ABQMAIN
Do K=1, nblock

Post_macro.for
Subroutine ABQMAIN
Do K=1, nblock
CALL
INITPF(FNAME,NRU,L
RUNIT,LOUTF)
…………………
CALL
DBFILE(LOP,ARRAY,J
RCD)
………….
Obtain nodal
displacements for the
element containing the
integration point
…………………
END DO
RETURN

σM

CALL
INITPF(FNAME,NRU,L
RUNIT,LOUTF)
…………………
CALL
DBFILE(LOP,ARRAY,J
RCD)
……………..
Obtain initial coordinate
sand reaction force for
surface nodes of RVE
…………………
Calculate Cauchy Stress

σM
…………..

END DO
RETURN

Figure 5.6 Numerical implementation flow chart of the concurrent multiscale modeling scheme in
ABAQUS/EXPLICIT
116

Osteon=5

Transverse
Time
0
1

S22(MPa)
0
31.5174

Table 5.3 The geometric parameters and elastic properties of RVE for all the cases studied.
Case No.
Young's Modulus of Osteon
ET/EL(GPa)
Young's Modulus of Matrix
ET/EL(GPa)
Young's Modulus of Cement Line
E (GPa)
Poisson's Ratio of Osteon
vT/νL
Poisson's Ratio of Matrix
vT/νL
Poisson's Ratio of Cement Line
v
Shear Modulus of Osteon
GTL/GTT(GPa)
Shear Modulus of Matrix
GTL/GTT(GPa)
Shear Modulus of Cement Line
GTL/GTT(GPa)
Numbers of Osteons
Diameter of Osteon (mm)
Area Percentage Occupied
by Osteon %
Loading Direction

1

2

11.51/20.71

11.51/20.71

12.66/22.78

12.66/22.78

11.51/20.71

5.755

0.17/0.18

0.17/0.18

0.153/0.162

0.153/0.162

0.17/0.18

0.4

8.4/8.4

8.4/8.4

3.57/3.93

3.573.93

8.4/8.4

N/A

5
0.14927

5
0.14927

35

35

Transverse

Transverse

stress tenser was sent back to Vumat to complete one loop of analysis.
The geometric parameters and elastic properties are shown in Table 5.3. The geometry and
mesh of the macro-model and micro-model are shown in Fig. 5.4. The mass of the projectile is 3.2
g and the impact speed is 300 m/s.

5.3
5.3.1

Results
Multiscale Modeling of Bovine Bone under Uniaxial Tensile Loading

The effective Young’s moduli (both longitudinal and transverse directions) of osteon bone were
calculated by the classical homogenization scheme to study the effects of osteon size and area
117

Time
0
1

osteon bone were calculated by the classical homogenization scheme to study the effects of
osteon size and area percentage occupied by secondary osteons. The simulations resemble the real
life uniaxial
tensile by
testing
with bottom
of the
specimen
fixed and
top end
percentage
occupied
secondary
osteons
ondog-bone
the macro
performance
of moving
corticalthe
bone.
The simwith aresembled
constant speed
of 0.1
Thetensile
reaction
force with
of thefixing
top end
used to
ulations
the real
life mm/s.
uniaxial
testing
the isbottom
of calculate
the dog-bone

specimen
and moving
end with aofconstant
speed
of 0.1
force
at the top
engineering
stress andthe
thetop
displacement
the top end
is used
to mm/s.
calculateThe
thereaction
engineering
strain.
endFigure
is used4 to
calculate
stress and
displacement
the top end
is used
to the
calculate
showed
bothengineering
of the longitudinal
andthe
transverse
Young’sat modulus
decrease
with
the engineering strain. Fig. 5.7 shows both of the longitudinal and transverse Young’s moduli
increase of numbers of osteons (decrease of osteon diameter). Stress Contours of macromodels

decrease with an increase in osteon quantity (decrease in osteon diameter).The stress contours of
and the unit cell corresponding to one of the integration point in macromodels are shown in

a macro-model and the unit cell corresponding to one integration point are shown in Fig. 5.8 to
Figure 5 to illustrate how the unit cell interacts with the macromodel. To study the effects of the

illustrate how the unit cell interacts with the macro-model. To study the effects of the area occuarea occupied by second osteon, another two cases with area percentage occupied by second

pied by secondary osteon, another two cases with area percentage occupied by secondary osteon
osteon of 17.5% are studied. All calculated Young’s moduli for various cases are summarized in

of 17.5% were studied. All calculated Young’s moduli for various cases are summarized in Table
5.4.Table 3.
b)

12
Young's Modulus (GPa)

21
Young's Modulus (GPa)

a)

20
19
18
17
16
15

11.5
11
10.5
10
9.5
9
8.5
8

0

2

4

6

8

10

12

0

Numbers of Osteon

2

4
6
8
Numbers of Osteon

10

12

Figure 4 Effective Young's modulus vs. osteon size for cortical bone with 35% area occupied by

Figure 5.7 Effective Young’s modulus versus osteon size for cortical bone with 35% area occupied
secondary osteon: (a) Longitudinal direction (b) Transverse direction.
by secondary osteon: (a) Longitudinal direction (b) Transverse direction.
14

5.3.2

Characterizing the Fracture Toughness of Cortical Bone with Various
Microstructures by XFEM

118

Young's
Modulus of Numbers of
Osteons
Cement
Line
12.66/22.78
5
5
5

a)

2
11.51/20.71
12.66/22.78
11.51
0.17/0.18
0.153/0.162
0.17
5
0.14927
35
Transverse

b)

Diameter of
Osteon (mm)
0.14927
0.14927
0.14927

3
NA
12.66/22.78
NA
NA
0.153/0.162
NA
0
NA
0
Transverse

Area
Loading
percentage
occupied by direction
Osteon
35
Transverse
35
Transverse
35
Transverse

5.755

5.755

Figure 5.8 Stress contours of a macro-model and the unit cell corresponding to one of the integraFigure 5 Stress Contours of macromodels and the unit cell corresponding to one of the
tion points
in macro-model: a) Longitudinal direction; b) Transverse direction.
integration point in macromodels: a)Longidinal direction; b) Transverse direction.
Table5.4
3 Summaries
calculated
Young's
modules
and and
its corresponding
microstructures
Table
Summariesofof
calculated
Young’s
moduli
their corresponding
microstructures.
Case
No.

Numbers

8

Case
Numbers of
Modulus Osteons
No.(GPa)
Osteons
(GPa)
20.011
20.01 3
3
19.72 5
5
19.722
19.44 7
7
19.443
4
19.00
10
19.005
10.92 10
3
10.60 3
5
10.926
7
10.30
7
10.60
5
8
9.81
10
10.309
20.25 7
5
10.81 10
5
9.8110

Diameter ofArea Percentage
Area percentage
Loading
Diameter of
Loading
Occupied
by
Osteon
(mm)
occupied by Osteon
Osteon
(mm)
Direction direction
Osteon %
0.19271
35
0.19271
35
LongitudinalLongitudinal
0.14927
35
LongitudinalLongitudinal
0.14927
35
0.12616
35
LongitudinalLongitudinal
0.12616
35
0.10555
35
Longitudinal
0.10555
35
0.19271
35
TransverseLongitudinal
0.14927
35
TransverseTransverse
0.19271
35
0.12616
35
Transverse
0.14927
35
Transverse
0.10555
35
Transverse
0.12616
35
0.10555
17.5
LongitudinalTransverse
0.10555
17.5
TransverseTransverse
0.10555
35

9

20.25

5

0.10555

17.5

Longitudinal

10

10.81

5

0.10555

17.5

Transverse

1
2
3
4
5
6
7

Young's modulus
Young's

of

It has been shown in the literature[10, 3, 4, 5]that the microstructure of cortical bone has strong
effects on its fracture toughness. The cement line seems to serve as a barrier to stop the microcrack
from propagating through osteons. Although the osteon bone shows a lower Young’s modulus
15

119

a (Cm)
0.7

W (cm)
1.4

a (Cm)
0.70759093

W (cm)
1.4

a (Cm)

W (cm)

a/w
0.5

f(a/W)
9.659079

B(mm)
0.05

a/w
f(a/W)
0.505422093 9.822543

B(mm)
0.05

a/w

f(a/W)

B(mm)

0.70759
1.4 bones, the matrix-cement line-osteon
0.505421429
9.822522
0.05
and tensile strength than
primary
microstructure
enhances
a (Cm)

W (cm)

a/w

f(a/W)

B(mm)

its resistance to fracture.
Thus,1.4
in this section, the fracture toughness0.505421429
of various microstructure
0.70759
9.822522
0.05

configurations was
evaluated W
to (cm)
study the interaction of osteon-cement line-interstitial
a (Cm)
a/w
f(a/W)bone with
B(mm)
0.7

1.4

0.5

9.659079

0.05

cracking. The force-displacement curves and the calculated fracture toughness for various configurations are shown in Fig. 5.9. The results for all the simulated configurations are shown in Table
5.5.
45
40
35
30
25
20
15
10
5
0

7

1
2
3
Case No.
4
5
6

0

6

GC/GO
εC/εO
EC/EO
1
0.5
1
2
2
0.5
1
1
3
1
1
1
4
0.5
0.5
1
5
2
1
1
2
1
6
0.5
0.1
0.2
7
0.5
1 0.3
0.5
Displacement
6 (mm)
0.5
1
2
7
0.5
1
1
8
1
1
1
Figure 5.9 a) Force-displacement curve; b)
9
0.5
2
1
0.5
2
0.5
10

5

Orientation

4 Transverse
3 Transverse

K
4.05
3.18
3.03
3.19
3.32

Force

Transverse
Transverse
55.02
1 Transverse
57.31
0 Transverse
Transverse
1
2
3
4
5
Transverse
4.05No.
Case
Transverse
3.18
Transverse
Calculated fracture 3.03
toughness.
Transverse
3.01
Transverse
2

46
45
43_new
50
51
48
6 52

Table 5.5 Summary of fracture toughness for all the cases studied.
Case No.

EC/EO

GC/GO

εC/εO

Orientation

K (MPa.m1/2)

1
2
3
4
5
6

0.5
0.5
1
0.5
2
0.5

1
1
1
0.5
1
1

2
1
1
1
1
0.5

Transverse
Transverse
Transverse
Transverse
Transverse
Transverse

6.6
5.7
5.51
5.71
5.67
3.9

46
45
43_new
50
51
52

C3D8R_50_1

The maximum
principal stress
contours at the last increment
of four studied
cases are
shownf(a/w)
in
Time
Force Da(mm) Mesh
a (cm)
w (cm)
a/w
Fig. 5.10. The top four contours show the crack extension on the CT specimen surface, the middle
four contours show the crack extension in the osteon fibers, and the last four contours show the
crack extension in the cement line. As shown in Table 5.5, the Young’s modulus of the cement line
is half that of the Osteon for Case 1, which yielded to a highest fracture toughness of all the cases
120

Case1

Case3

Case4

Case5

Case1
Osteon

Case3
Osteon

Case4
Osteon

Case5
Osteon

Case1
Cement
line

Case3
Cement
line

Case4
Cement
line

Case5
Cement
line

Figure 5.10 Stress contour at the last increment of the four of the cases (Case 1, 3, 4, and 5) studied
in the CT fracture toughening simulations
121

studied. The second highest fracture toughness is Case 4 in which the Young’s modulus is also half
that of the osteon. The contour plots in Fig. 5.10 give an insightful view of how cracks propagated
through the osteons and cement lines inside the bone matrix. It shows that cracks only propagated
through the first two layers of osteons in Cases 1&4, and cracks penetrate three layers of osteons in
Cases 3&5. The contours also showed that cracks tended to deflect along the longitudinal direction
of osteons in Cases 1&4 instead of extending across osteons.

0.05 The 0.1
0.3configuration at different crack extension was also
fracture 0.15
toughness0.2
of each0.25
simulated
Displacement (mm)
calculated. The results are shown in Fig. 5.11.
9
8
7
6
5
4
3
2
1

1

2

3

4

5

6

0
0

0.1

0.2

0.3
∆a (mm)

0.4

0.5

0.6

Figure 5.11 Fracture toughness versus cracking extension.

5.3.3

Semi-Concurrent Multi-scale Modeling of Impact Testing on Cortical
Bone

52
43
The boundary
conditions applied to the micro-model were passed from the element 25 of the

51
45
macro-model
(Fig. 5.12) at step time 0.6 ms. The results of the two cases (Table 5.5) studied
4 are shown
0.5 in Fig.
0.65.13. It can be seen that the maximum tensile stresses in osteons 3&4 were
122

constrained to the osteons (see the stresses contour in the top figures in Fig. 5.13). However, in
Case 1 the stresses in the two osteons extended to the surrounding interstitial bone. Furthermore,
tensile stress zones were observed in Osteons 3&5 in Case 1, but not in Case 2. The results
indicated that soft interface introduced by cement lines interfered with the wave traveling inside
the bone. It worked as a barrier to separate osteons from the surrounding bone. Therefore, the
cracks and failure were constrained in the local area, which prevented the macroscopic failure of
cortical bone.

Element 25

Macro-Model
Deformed Contour at 0.6 ms

Figure 5.12 Geometry and mesh of macro-model of impact testing simulated.

5.4

Discussion

The multiscale modeling using classical homogenization method for evaluating the material properties of heterogeneous materials [23, 24, 25, 17, 20, 16] was used to evaluate the macroscopic
123

Case 1
Step Time: 0.6 ms

Osteon 1

Osteon 5

Osteon 3

Case 2
Step Time: 0.6 ms

Osteon 2

Osteons

Osteon 4

Figure 5.13 Stress contour of the two micro-models studied at impact time of 0.6 ms.

124

Young’s modulus of cortical bone with various microstructure configurations. The Young’s modulus in both longitudinal and transverse directions was found to decrease with the increase of
numbers of osteon or decrease of the diameter of osteon with a constant area percentage occupied
by secondary osteon. The reason is that increasing the numbers of osteon increases the porosity of
the cortical bone because each osteon comprises a central canal, Harvention system.
The results also showed that macroscopic Young’s modulus in the longitudinal direction decreased with the area percentage occupied by the secondary osteon with a constant osteon size.
The results agree with the experimental observations in the literature. Saha and Hayes [10, 26]
measured the tensile strength and the modulus of elasticity, which was found to be negatively
correlated with the area percentage of secondary osteons.
However, the central canal and cement effects on the Young’s modulus in the transverse direction are slightly different than in the longitudinal direction. The Young’s modulus is 10.81 GPa
with an area percentage of 17.5 and 5 osteons (Case 10). This value is slightly lower than the value,
10.98 GPa, with an area percentage of 35 and 3 osteons(Case 6) but larger than all the other cases.
Thus, it seems that there exists an optimal diameter of osteons for the Young’s modulus in the
transverse direction. This difference between the effects of area occupied by secondary osteons on
the Young’s modulus of longitudinal and transverse directions is likely due to the effect of cement
lines. Case 6, 7, and 8 include more osteons or more soft interfaces than case 10. Thus, lower
Young’s moduli were observed in Case 6, 7, and 8 than Case 10. However, compared to Case 10,
Case 5 comprised less cement lines. Therefore, Case 5 exhibited a higher Young’s modulus than
Case 10. The experimental results of the uniaxial tensile testing of cortical bone in Chapter 4 also
showed that cement lines weaken the cortical bone in transverse direction with that many microcracks instantly initiated in the cement lines resulting in a quick stiffness loss in the transverse
specimen.
125

The fracture toughness of cortical bone was evaluated by XFEM method to investigate the role
of cement lines in the fracture process of cortical bone. The fracture toughness calculated for all
the cases studied are in the range of 3.9-6.6 MPa.m1/2 (Table 5.5). The results agree with the
literature [4, 3, 7]very well. The results also showed that a case with a Young’s modulus of cement
lines half of its surrounded osteon and a cracking initiation stain double of the surrounded osteon
exhibited the highest fracture toughness among all the cases studied. However, a brittle cement
line interface seemed to decrease the fracture toughness. The results agree with the experimental
observation that the cement line acts as a barrier to stop or deflect the crack propagation. Crack
deflection, microcracking, and uncracked ligament bridging are the major toughening mechanism
of cortical bone [4, 3, 7, 6, 8, 5]. The simulations support the hypnosis that the osteon-cement
line-interstitial bone has benefits in the bone fracture process and properties of cement lines play a
critical role in it.
The effects of fracture energy including normal modes and two shear models of cement lines
on the fracture toughness of the C(T) specimen were also evaluated. In Case 4 and Case 2, all the
other factors are kept the same besides the values of fracture energy. However, very close results
were obtained for Case 2 and 4. Compared to the fracture initiation criterion, fracture evolution
seemed to have less effect on the results.
The rising fracture resistance with the crack extension was also studied. The fracture toughness
versus crack extension is plotted in Fig. 5.7. Linear regression lines are added to each case.
The fracture toughness at zero crack extension for all the configurations studied was the same,
2.72MPa.m1/2 . This is due to the first layer at the precrack front was the interstitial bone, which
has the same material properties for all the cases. After crack initiation, rising R-curve behavior
was observed all the cases, which agrees with the literature very well [4]. Each case exhibited a
different rising slope. Case 1 with a soft and ductile interface showed the steepest rising R-curve
126

and case 6 with brittle interface seemed to exhibit a most flat rising R-curve among all the cases.
The role of cement lines on the impact behavior of cortical bone was also evaluated through the
use of semi-concurrent modeling approach. The results showed that a RVE with a soft interface
introduced by cement lines exhibited less damage in both interstitial bone and osteons than the
RVE without cement line. Thus, the soft interface also benefits impact behavior of cortical bone,
which may explain the higher impact resistance of osteon bone observed in the Charpy impact
testing in Chapter 4.
Osteon bone comprising three constituents (secondary osteon, cement line, and interstitial
bone) in the microscopic level is hypothesized to exhibit fiber-reinforced composite material behavior. With different morphology and composition compared to its surrounded osteon, the cement
line is assumed to act as a barrier or soft interface to prevent microcrack penetrating through osteon
to lead to a catastrophic failure. Understanding the micro interactions between cement, osteon, interstitial bone, and micorcracking helps disclose the fracture mechanism of cortical bone. In the
present paper, the correlation between the micro components and the macro performance of osteon bone was studied using multi-scale modeling methods and XFEM technique. By using the
classical computational homogenization method, the tensile moduli of cortical bone with various
microstructures were evaluated. The Young?s moduli in both longitudinal and transverse directions
were found to have a positive correlation with the diameter of osteons when the area percentage
occupied by secondary osteons was kept constant. And Young’s modulus in longitudinal direction
was found to decrease with the increase of the area percentage occupied by secondary osteons.
However, it seems that there exists an optimal diameter of osteons benefiting the Young’s modulus
in the transverse direction due to effects of the weaker interface introduced by cement lines. The
fracture toughness of cortical bone with various micro configurations was characterized using the
ABAQUS/XFEM technique. The highest fracture toughness was the case with a soft and duc127

tile cement line interface. Raising R-curve behavior was observed. A semi-concurrent modeling
scheme for impact loading simulations was developed and implemented in ABAQUS/EXPLICIT
6.11 to investigate the role of cement lines on the impact behavior of cortical bone. The study
showed that soft interface introduced by cement lines interfering with the wave traveling inside the
RVE (representative volume element) resulted in a mitigating in the micro damage compared to
the RVE without cement lines.

128

REFERENCES

129

REFERENCES
[1] Jae-Young Rho, Liisa Kuhn-Spearing, and Peter Zioupos. Mechanical properties and the
hierarchical structure of bone. Medical engineering & physics, 20(2):92–102, 1998.
[2] Mitchell B Schaffler, David B Burr, and Richard G Frederickson. Morphology of the osteonal
cement line in human bone. The Anatomical Record, 217(3):223–228, 1987.
[3] Ravi K Nalla, Jamie J Kruzic, and Robert O Ritchie. On the origin of the toughness of
mineralized tissue: microcracking or crack bridging? Bone, 34(5):790–798, 2004.
[4] RK Nalla, JS St¨olken, JH Kinney, and RO Ritchie. Fracture in human cortical bone: local
fracture criteria and toughening mechanisms. Journal of biomechanics, 38(7):1517–1525,
2005.
[5] Fergal J O’Brien, David Taylor, and T Clive Lee. Microcrack accumulation at different
intervals during fatigue testing of compact bone. Journal of Biomechanics, 36(7):973–980,
2003.
[6] D Vashishth, KE Tanner, and W Bonfield. Contribution, development and morphology of microcracking in cortical bone during crack propagation. Journal of Biomechanics, 33(9):1169–
1174, 2000.
[7] D Vashishth, KE Tanner, and W Bonfield. Experimental validation of a microcracking-based
toughening mechanism for cortical bone. Journal of Biomechanics, 36(1):121–124, 2003.
[8] D Vashishth, JC Behiri, and W Bonfield. Crack growth resistance in cortical bone: concept
of microcrack toughening. Journal of Biomechanics, 30(8):763–769, 1997.
[9] David B Burr, Mitchell B Schaffler, and Richard G Frederickson. Composition of the cement
line and its possible mechanical role as a local interface in human compact bone. Journal of
Biomechanics, 21(11):939–945, 1988.
[10] S Saha and WC Hayes. Relations between tensile impact properties and microstructure of
compact bone. Calcified tissue research, 24(1):65–72, 1977.
[11] Susan Mischinski and Ani Ural. Finite element modeling of microcrack growth in cortical
bone. Journal of Applied Mechanics, 78(4):041016, 2011.
[12] Ted Belytschko and Tom Black. Elastic crack growth in finite elements with minimal remeshing. International journal for numerical methods in engineering, 45(5):601–620, 1999.
130

[13] Ted Belytschko, Nicolas Mo¨es, Shuji Usui, and Chandu Parimi. Arbitrary discontinuities in
finite elements. International Journal for Numerical Methods in Engineering, 50(4):993–
1013, 2001.
´
[14] Elisa
Budyn and Thierry Hoc. Multiple scale modeling for cortical bone fracture in tension using x-fem. European Journal of Computational Mechanics/Revue Europ´eenne de
M´ecanique Num´erique, 16(2):213–236, 2007.
[15] Varvara Kouznetsova, Wouter Wilson, and Marc Geers. Multi-scale modeling of textile composites thiam wai chua mt 10.19. 2011.
[16] V Kouznetsova, WAM Brekelmans, and FPT Baaijens. An approach to micro-macro modeling of heterogeneous materials. Computational Mechanics, 27(1):37–48, 2001.
[17] Varvara Kouznetsova, Marc GD Geers, and WA Marcel Brekelmans. Multi-scale constitutive
modelling of heterogeneous materials with a gradient-enhanced computational homogenization scheme. International Journal for Numerical Methods in Engineering, 54(8):1235–1260,
2002.
[18] Wei Zhang, Srinivasan Arjun Tekalur, Melissa Baumann, and Laura R McCabe. The effects of damage accumulation on the tensile strength and toughness of compact bovine bone.
Journal of biomechanics, 46(5):964–972, 2013.
[19] ASTM Designation. E399-90. standard test method for plane-strain fracture toughness of
metallic materials. 1991 Annual Book of ASTM Standards, 3:485–51, 1997.
[20] VG Kouznetsova, MGD Geers, and WAM Brekelmans. Multi-scale second-order computational homogenization of multi-phase materials: a nested finite element solution strategy.
Computer Methods in Applied Mechanics and Engineering, 193(48):5525–5550, 2004.
[21] Stefan Loehnert, L Krstulovic-Opara, Matej Vesenjak, and Dana Mueller-Hoeppe. Homogenization principle based multi-scale modeling of cellular structures. Journal of the Serbian
Society for Computational Mechanics/Vol, 4(2):97–109, 2010.
[22] ABAQUS Version. 6.11 documentation. Dassault Systemes Simulia Corp., Providence, RI,
USA, 2011.
[23] Christian Miehe, J¨org Schr¨oder, and Martin Becker. Computational homogenization analysis in finite elasticity: material and structural instabilities on the micro-and macro-scales
of periodic composites and their interaction. Computer Methods in Applied Mechanics and
Engineering, 191(44):4971–5005, 2002.
131

[24] Christian Miehe and Nikolas Apel. Anisotropic elastic–plastic analysis of shells at large
strains. a comparison of multiplicative and additive approaches to enhanced finite element
design and constitutive modelling. International journal for numerical methods in engineering, 61(12):2067–2113, 2004.
[25] PE McHugh, RJ Asaro, and CF Shih. Computational modeling of metal matrix composite
materials. Modeling the deformation of crystalline solids, 1991.
[26] S Saha and WC Hayes. Tensile impact properties of human compact bone. Journal of biomechanics, 9(4):243–251, 1976.

132

Chapter 6
Conclusion

6.1

New and Original Contributions Made

Low velocity drop weight impact tests were conducted to evaluate the threshold conditions indicating the onset of damage and failure mechanism of monolithic and bi-layered PMMA and
PC plates. A material model based on the Arruda and Boyces physical model for amorphous
polymers was coded and implemented in ABAQUS/EXPLICIT to predict the strain rate and pressure dependent behavior of polymers. The impact behavior of monolithic PMMA and PC and
layered PMMA/PC laminate at low and high impact loading were numerically investigated by
use of ABAQUS/EXPLICIT. The failure mechanism, crack pattern, and energy dissipation were
transversely compared. The FEM simulation provided a comprehensive insight into the energy
dissipation and toughening mechanism of the impact events. An optimal design was performed
interactively using the optimization tool in MATLAB and commercial finite element modeling
(FEM) software ABAQUS/EXPLICIT to study the effects of aspect ratio, material properties of
PMMA and PC, etc., on the toughening mechanism of the laminate panel. This work serves as a
design tool for effectively choosing materials and structures for optimizing ballistic performance
of ceramic-plastic laminate-based armor.
The prior literature on bone failure and fracture mechanism showed that the microcracking
toughening mechanism is one of the major toughening mechanisms of cortical bone by use of
the post-processing microcracking detection technique. In this research work, the microcracking
133

toughening mechanism was visualized concurrently with the fracture process of cortical bones
subjected to uniaxial tensile loading. Strain hardening and softening were observed in the longitudinal and transverse specimens respectively. The microcracking in the transverse specimens were
found to instantly increase to peak after yielding compared to the gradual growth until failure in
the longitudinal specimens. Transverse specimens showed that many microcracks instantly initiated in the cement lines, which resulted in the stiffness loss in the specimen. The present results
showed that, to some extent, microcracking serves as a toughening mechanism; but if it is over a
critical value, microcracking weakens the specimen instead. The Charpy impact testing and the
SHPB high strain rate uniaxial tensile testing were conducted to study the failure behavior of the
cortical bone under high strain rate loading. The results showed that microcracking seems to not
be the toughening mechanism of cortical bone at high strain rate failure. Rather, it is the unique
osteon-cement line-interstitial bone microstructure, which plays an important role in toughening
the failure behavior of the cortical bone at high strain rate. A multiscale modeling based on the
classic first order homogenization scheme was implemented in ABAQUS/STANDARD to study
the constitutive behavior of the cortical bone in terms of the osteon size, orientation, and area percentage occupied by osteon. The computational results verified the experimental observations that
tensile modulus of cortical bone has a negative correlation with the area percentage occupied by
osteon. XFEM element with linear fracture mechanics method was used to evaluate the role of cement line, osteon, interstitial bone, and their corresponding mechanical properties on the fracture
toughness of bone. The results numerically proved that the experimental observation that cement
lines, the soft interface between osteon and its surrounding bone matrix, act as barriers to stop or
deflect the crack propagation. A semi-concurrent modeling scheme for impact loading simulations
was proposed and implemented in ABAQUS/EXPLICIT to investigate the role of cement lines on
the impact behavior of the cortical bone. The study showed that the soft interface introduced by
134

cement lines interferes with the wave traveling inside the RVE, resulting in a greater mitigation in
the micro damage compared to an RVE without cement lines. The research work provides inspiration for bone replacement material and impact resistance material design and synthesis, as well as
improving the understanding of bone fracture.

6.2

Future Work

As for future work, three research directions extended from this work would be interesting: methods bridging failures of RVE to the failure of macro-model for both ABAQUS/STANDARD and
ABAQUS/EXPLICIT, an optimization scheme for multiscale modeling, and mesh independent
cracking method. Although multiscale modeling methods have been studied and developed extensively in the past decades, bridging the micro and macro failure is still challenging nowadays.
Similar difficulties exists in the mesh independent cracking method. Although ABAQUS/XFEM is
a crack path and remeshing independent, it fails to converge when multi-cracks are initiated in one
element. Advance in the computational failure and cracking method is the key for resembling the
real life failure behavior of materials in simulations. Optimization method for multiscale modeling
is also of great engineering interest for material optimal design.

135