I mCHIGALi‘sBTATE ES
. (.1 UNIVERSITY
,1 ) 1 ,» EAST LANSING. MICH 48824-1048

This is to certify that the
thesis entitled

MULTl-LAYERED EXPERIMENTAL AND NUMERICAL
ANALSYSIS OF STAMP THERMOFORMING PROCESSING
OF NATURAL FIBER REINFORCED POLYPROPYLENE
SHEETS

presented by
BARBARA NICOLE RODGERS

has been accepted towards fulfillment
of the requirements for the

M. S. degree in Mechanical E ineering

 

CQWM

Major Professor’s Signature
t1 / 15 / 0 ‘1‘

Date

MSU is an Affirmative Action/Equal Opportunity Institution

 

PLACE IN RETURN BOX to remove this checkout from your record.
TO AVOID FINES return on or before date due.
MAY BE RECALLED with earlier due date if requested.

 

DATE DUE

DATE DUE

DATE DUE

 

 

 

 

 

 

 

 

 

 

 

 

 

 

6/01 c-JCIRCIDatoDuepes-pts

 

MULTI-LAYERED EXPERIMENTAL AND NUMERICAL ANALSYSIS OF
STAMP THERMOFORMING PROCESSING OF NATURAL FIBER
REINFORCED POLYPROPYLENE SHEETS

BY

Barbara Nicole Rodgers

A THESIS

Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of

MASTER OF SCIENCE

Department of Mechanical Engineering

2004

ABSTRACT
MULTl-LAYERED EXPERIMENTAL AND NUMERICAL ANALSYSIS OF
STAMP THERMOFORMING PROCESSING OF NATURAL FIBER
REINFORCED POLYPROPYLENE SHEETS
BY
Barbara Nicole Rodgers

Due to increasing ecological concerns, much attention has been given to
the study of natural fiber composites to replace non-biodegradable glass and
woven fabric composites. This study focuses predominantly on kenaf natural
fibers with a polypropylene matrix. The first objective was to produce a fabrication
process that yielded composites with comparable properties to existing synthetic
fiber composites. A detailed outline of this optimized process for compression
moulding is presented. In addition, material characterization was done which
included squeeze flow, tensile, and flexural testing. The results show that this
newly developed fabrication process produced composites with competitive
flexural and tensile strength to other natural fiber composites both using
compression moulding and other fabrication processes.

The fabricated kenaf/polypropylene composites were then formed using a
stamp thermo-fon'ning press. The material showed good formability. In addition,
numerical analysis was done which included single and double layered modeling.
An updated material law constitutive relationship for multiple preferred fiber
orientations was used to model the laminae, while two interface models were
used: linear sOftening and plasticity. Additional research was done on other

possible material and interface models to be used in the future.

DEDICATION

“I can do all things through Christ which strengtheneth me. ” -Philipians 4: 13

ACKNOWLEDGMENTS
Thank You!
Giving all honor and glory to God, who is truly the head of my life.

To my beloved parents, I will always be grateful for your years of
unconditional love and support.

To my loving husband, who is always right by my side. Without his love
this definitely would have been a lot harder to achieve.

Special thanks to my advisor Dr. Farhang Pourboghrat for his leadership
and guidance throughout my graduate experience. In addition, thanks to Dr.
Michael Zampaloni for being willing to share his expertise in order to help me.
Furthermore, to Stacey Yankovich, with whom I worked closely with during this
entire process.

This work would also have not been possible without assistance from the
Composite Materials and Structures Center at Michigan State University. In
particular, Dr. Drzal, Dr. Misra, and Dr. Huda, who were all very helpful in providing
information or materials needed to make this research a successful endeavor.

Sincerely,
Barbara Nicole Rodgers
December 2004

East Lansing, Michigan

iv

TABLE OF CONTENTS

LIST OF TABLES ......................................................................................... vi
LIST OF FIGURES ....................................................................................... vii
ABBREVIATIONS ........................................................................................ viii
CHAPTER 1 INTRODUCTION AND LITERATURE REVIEW ..................... 1
1.1 Introduction to Biocomposites ................................................. 1
1.2 History of Composite Modeling ............................................... 5
1.2.1 Cohesive Zone Model .................................................. 6
1.2.2 Continuum Damage Model .......................................... 9
1.2.3 Rate Dependent Models ................................................. 12
1.2.4...|nterface Models ............................................................ 16
CHAPTER 2 BIOCOMPOSITE .................................................................... 23
2.1 Fabrication of Kenaf/Polypropylene Composite ...................... 23
2.2 Material Characterization — Directionality ............................... 29
2.3 Material Characterization - Properties .................................... 33
2.4 Compression Moulding Process and Composite Comparison.38
CHAPTER 3 EXPERIMENTAL WORK ....................................................... 42
3.1 Introduction to Stamp Thermoforming .................................... 42
3.2 Experimental Stamping Results .............................................. 44
CHAPTER 4 NUMERICAL ANALYSIS ........................................................ 47
4.1 Material Model ........................................................................ 47
4.2 Numerical Stamping Results ................................................... 50
4.2.1 Two Preferred Fiber Orientations ................................. 50
4.2.2 Three Preferred Fiber Orientations .............................. 57
4.2.3 Four Preferred Fiber Orientations ................................ 59
4.3 Discussion of Results ............................................................. 62
CHAPTER 5 CONCLUSIONS ..................................................................... 69
CHAPTER 6 FUTURE WORK ..................................................................... 70
BIBLIOGRAPHY .............................................................................................. 71

Table 1.
Table 2.
Table 3.
Table 4.
Table 5.

Table 6.

LIST OF TABLES

Plasticity model properties of the interface .................................... 20
Tensile testing results for all specimens — Young’s Modulus ......... 35
Tensile testing results for all specimens —- Poisson’s Ratio ........... 36
Flex testing results for all specimens — Modulus of Elasticity .......... 37
Material properties for simulation ..................................................... 37
Sensitivity analysis of Young’s Modulus ........................................... 61

LIST OF FIGURES

Figure 1. Linear function of pressure-overclosure relationship ............................ 18
Figure 2. (a) Carver Press (b) Tetrahedron Press ............................................. 24
Figure 3. Long Kenaf fibers sandwiched between two PP sheets .................... 25
Figure 4. Chopped Kenaf fibers sandwiched between two PP sheets ............. 25
Figure 5. Premixed chopped Kenaf fibers and PP ............................................ 25
Figure 6. Multilayered composite fabrication process ....................................... 26
Figure 7. Multilayered fabricated composite ......................................................... 27
Figure 8. Squeeze flow test sample before testing .............................................. 29
Figure 9. Schematic of squeeze flow test experimental setup [27]. .................... 30

Figure 10.
Figure 11.
Figure 12.
Figure 13.

Figure 14.

Squeeze flow test sample after testing .............................................. 30
Normalized length versus angle for all samples ................................ 31
Normalized length versus angle for best samples ........................... 32
Kenaf/PP and PP tensile specimens ............................................... 33

Stress strain curve for 40° fiber orientation specimens at room
temperature ..................................................................................... 34

Figure 15. Stress strain curve for 170° fiber orientation specimens at room

Figure 16.

Figure 17.
Figure 18.

Figure 19.

temperature ..................................................................................... 34

Stress strain curve for Polypropylene/Epolene specimens at room
temperature ........................................................................................ 35

Tensile strength comparison ............................................................... 39
Flexural strength comparison .............................................................. 39

Flexural strength comparison for various manufacturing processes
............................................................................................................ 4O

Figure 20.
Figure 21.
Figure 22.
Figure 23.
Figure 24.
Figure 25.
Figure 26.
Figure 27.
Figure 28.
Figure 29.
Figure 30.

Figure 31.

Figure 32.
Figure 33.

Figure 34.

Figure 35.
Figure 36.

Figure 37.

Figure 38.

Figure 39.

Modulus per cost comparison ............................................................ 41

Specific modulus comparison .......................................................... 41
Stamp thermo-hydroforming press ................................................... 42
Stamp thermoforming schematic ...................................................... 43
Stamp thermoformed part with tearing - top view ............................... 44
Stamp thermoformed part A, (a) top view, (b) side view ...................... 45
Stamp thermoformed part B, (a) top view, (b) side view ...................... 45
Stamp thermoformed part C, (a) top view, (b) side view ...................... 46

Deformed blank for v=0.2 (a)Vf=O.218 (b) Vg=0.384 (c) VFO.6 ......... 53
Deformed blank for v=0.3 (a)Vf=O.218 (b) Vf=0.384 (c) VFO.6 ......... 54
Deformed blank for v=0.4 (a)Vf==O.218 (b) VF—‘O.384 (c) Vf=0.6 ........ 55

Deformed blank for v=O. 3 and VFO. 384 with a plasticity model

interface .......................................................................................... 56
Deformed blank for 3 preferred fiber orientations VUMAT .............. 57
Stress contour plot for 3 preferred fiber orientations VUMAT .......... 57

Internal (a) and Kinetic (b) energy plots for 3 preferred fiber orientations

VUMAT ............................................................................................ 58
Deformed blank for 4 preferred fiber orientations VUMAT .............. 59
Stress contour plot for 4 preferred fiber orientations VUMAT ........... 59

Internal (a) and Kinetic (b) energy plots for 4 preferred fiber orientations

VUMAT ............................................................................................ 60
Front angled views of formed part and simulations ......................... 62
Top views of experimental and simulation results .............................. 63

viii

ABBREVIATIONS

Pi. - Force

Kij. — Diagonal Matrix of Penalty Stiffnesses
Dj. - Displacement

F. — Load Surface

f(P). — Normalized measure of P

P. — Pressure

‘L'i. — Interlaminar Direct and Shear Stresses
kio. — Penalty Stiffnesses of the Interface

5,. - Displacements

kio — Initial Stiffness

k, —105~107mm"

Tic. — Interlaminar Tensile and Shear Strengths

fs(r,) — Failure Function
fg(G,) — Failure Criterion

G, - Energy Release Rates

G, G", G... - Critical Energy Release Rates in Mode I, II, and Ill

HEIGI) — Damage Surface

(p — Function to Control Shrinkage Rate
t — Traction Vector

F - Surface

[u] - Vector of Displacement Jumps
K — Elastic Interface Stiffnesses
f(t,}t) — Viscoplastic Potential

h — Hardening/Softening Coefficient
t — Time

N — Model Parameter

7A4 — Material Parameter

A — Viscoplastic Multiplier

D, - Damage Variables

Yi— Conjugate to Damage Variables
E11 — Longitudinal Young’s Modulus
E22 — Transverse Young’s Modulus
En — Young’s Modulus for ith Fiber Direction
E, - Young’s Modulus of the Fiber
Em — Young’s Modulus of the Matrix
v12 , v21 — Poisson’s Ratio

v, - Poisson’s Ratio of the Fiber

vm — Poisson’s Ratio of the Matrix
G12 — Shear Modulus

G13, G23 — Transverse shear Moduli,
Gf — Shear Modulus of the Fiber

Gm — Shear Modulus of the Matrix

N - Number of Plies

[Q] — Stiffness Matrix

[S] - Compliance Matrix

[T (9)] — Transformation Matrix

6 — Angle between Material Frame and Fiber Directions
[51, - Transformed Reduced Stiffness Matrix

t — Ply Thickness

Vfo" Initial Fiber Volume Fraction

V; — Fiber Volume Fraction

[0].. - Stress Tensor

[8].. - Strain Tensor

p — Density of the Composite Material

fvd — Damping Forces

p0 — Critical Damping Fraction

m - Nodal Mass

k9. — Nodal Contact Stiffness

Eimeme - Young’s Modulus of the Interface
Gimme — Shear Modulus of the Interface
timenace — Thickness of the Interface

de — Nominal Strain Increment

ot — Elastic Trial Stress

At - Relative Displacement

011—— Normal Stress

oy— Yield Stress

Eh — Hardening Modulus of the Interface
dam - Plastic Strain Increment

k — Thermal Conductivity

Tmaste, - Temperature of Master Node
Tsjave — Temperature of Slave Node

A - Nodal Area

At — Stable Time Increment

L° - Characteristic Element Length

0,, — Dilation Wave Speed of the Material

Chapter 1

INTRODUCTION AND LITERATURE REVIEW

1.1 — Introduction to Biocomposites
In the past few decades, research and engineering has shifted from

monolithic materials to fiber-reinforced polymer materials. Although glass fibers
are the most widely used reinforcements of plastics due to their low cost and
good mechanical properties, they have quite a few drawbacks. These include
high density the fact that they are not renewable, recyclable, biodegradable, or
002 neutral. These fibers also are abrasive to machines and present a health
risk when inhaled. Therefore, natural fibers have come to the forefront in
research, and unlike glass fibers, do not share these same shortcomings [1].

In recent years, there have been many discussions about preserving the
world’s natural resources and recycling. This has led to an interest in the use of
biodegradable materials that come from renewable raw materials in composites.
A biocomposite consists of a biodegradable polymer as a matrix and a biofibre as
reinforcement. In most applications of composites, it is most beneficial to have a
matrix that is resilient to the environment, thus, for now, composite research
mostly focuses on natural fibers within a non-biodegradable synthetic polymer
matrix [2]. While, natural fibers traditionally were used to reinforce thermosets,
greater attention is now being placed on using natural fibers to reinforce

thermoplastics, such as polypropylene, due to its recyclability [3].

The automobile industry, in particular, is constantly looking for ways to
decrease weight and cost of interior and engine components. It is estimated that
in the near future polymer and polymer composites will make up 15% of each
automobile [4]. Some applications are in the door and instrument panels, glove
boxes, package trays, arm rests, and seat backs [5]. These natural fiber
composites are an alternative to the wood or glass reinforced composites
currently used. The main advantages include low cost, low density, good specific
strength, enhanced energy recovery, reduced tool wear and dermal/respiratory
irritation, and biodegradability [6]. The disadvantages include a low processing
temperature and its general hydrophilic (water absorbent) nature, which is not
compatible with hydrophobic (non-water absorbent) polymeric matrices [2]. This
problem can be remedied by surface treatment of the fiber or by adding a
coupling agent, which will be discussed later in detail. In addition, the use of
surface treatments or a coupling agent can enhance the mechanical properties of
natural fiber-reinforced composites by ensuring good adhesion between the fiber
and the matrix [7].

For this study, the natural fiber chosen was Kenaf with a Polypropylene
(PP) matrix. Kenaf fiber (Hibiscus cannabinus) is native to Africa and is also
grown in the United States in tropical and subtropical regions. The approximate
cost of Kenaf bast fibers are 20 - 25 9: per pound (40 — 45 a: per kg) versus 90¢
per pound ($2 per kg) for E-glass fiber [7]. Natural fibers characteristics can vary
considerably depending on source, age, extraction period and technique,

weather variability, and the quality of the soil and climate of the particular region

it is grown in [8]. Yet, studies have shown that fibers from the same batch that
contain similar histories grown in moderate conditions result in composites that
have low variation in properties. In these conditions, Kenaf fiber can grow as
much as 9.8 ft (3 m) in 3 months. Although research has only recently begun to
take a serious interest in Kenaf as reinforcement, mankind has used Kenaf even
in ancient times as canvas, rope, and sacking [9].

The main constituents of the natural fiber are cellulose, hemicellulose, and
Iignin. Cellulose is a hydrophilic glucan polymer consisting of linear chains of
anhydroglucose units that contain alcoholic hydroxyl groups. This is the reason
all natural fibers are hydrophilic: their moisture content reaches 8 — 12.6%. The
chemical structure of cellulose is the same for different natural fibers although the
degree of polymerization varies, which is related to the mechanical properties of
the fiber [6].

The Iignin is a biochemical polymer that structurally supports the plant. It
is generally resistant to microbial degradation, but surface treatment of the fibers
can tend to make it susceptible to the cellulose enzyme. When the plant cell
walls are synthesized, cellulose and hemicellulose are present first, and then
Iignin fills the spaces cementing them together. This process is known as
lignification and causes a stiffening of the cell walls in the fiber [2].

Polypropylene (PP) belongs to a class of materials known as polyolefins.
These polymers account for approximately 60 % of the total polymer produced.
They are so widely used because they have good chemical and physical

properties including low cost, good recyclability, and good processibility. Out of

all the polyolefins, polypropylene is known to be one of the most recyclable.
Polyolefins are fabricated by various methods including injection moulding,
extrusion, compression moulding, injection blow molding, and calendaring into
films, fiber, pipe, and sheet products of various shapes and sizes [10].
Polypropylene is a thermoplastic material that is produced by a
polymerization process that aligns propylene molecules, or monomer units, into
long molecules or chains. This process produces a semicrystalline solid that is
characterized by high rigidity, heat resistance, and melting point (~ 1570), low
density, and relatively good impact resistance. In the solid state, the properties
are dependent upon the amount of crystalline and amorphous regions produced
during forming, while, in the melting state, they are determined by the average
length of the polymer chains and their distribution in a product. Therefore, the
properties can be altered simply by changing the chain lengths, distribution, or by
adding a comonomer such as ethylene into the polymer chains or impact modifier
into the resin mixture [11]. These properties make polypropylene ideal for use as

a matrix in the stamp thermo-hydroforming process.

1.2 — History of Composite Modeling
The finite element analysis of failure in composites has been addressed

using many different types of interface elements and interface laws. The most
typical mode of interfacial damage in laminated composites is delamination. It
occurs between the constituent layers in a composite due to their weak
interlaminar strengths. This can happen under various types of loading including
transverse concentrated loads in the form of low velocity impacts, which is the
case in stamp thermoforming. The initiation of delamination is commonly
predicted using stress-based methods, while the propagation of delamination is
described using fracture mechanics [16].

When it comes to simulating delamination, or the separation of lamina in
composites, there are two model groups — fracture mechanics or damage
mechanics/softening plasticity models. If there are no non-linearities in the
material behavior, linear elastic fracture mechanics (LEFM) provides an accurate
method for prediction of delamination, provided that the initial location of
delamination is known. Some examples are the virtual crack closure, virtual
crack extension, J-integral, and stiffness derivative methods. The idea is that a
critical energy release rate value must be reached before delamination can
propagate [17].

These methods are often straight forward for two dimensional analyses
since the crack forms and propagates in one direction, which is not the case in
three dimensions. In addition, the initial location of crack formation, or flow, is

known. Therefore, great effort has been put into developing models that use a

damage mechanics and/or softening plasticity model combined with fracture
mechanics to not only capture the initiation of the delamination, but its
propagation. This work began with a fictitious crack model (FCM) developed in
1976 by Hillerborg et al. FCM uses a model where stresses act across a crack
while it is narrowly opened and introduced fracture mechanics into finite element
analysis [18]. This model was the basis for many other models. The models
discussed in this study are the cohesive zone model (CZM), the continuum
damage modeVpenalty technique, and rate dependent models.

1.2. 1 - Cohesive Zone Model
The history of the cohesive zone model can be traced as far back as 1960

when Dugdale developed a method to trace the spread of plasticity from the
center of stress concentration as loads increase. The yielding at the end of a slit
in a sheet was investigated to find a relationship between the external load
applied and the amount of plastic yielding. In order to verify the accuracy of the
model, panels containing both internal and external slits were loaded in tension
so that the lengths of the plastic zones could be measured [19]. This model
eventually took many different forms. In 2000, Borg et al simulated delamination
using the explicit FE code LS-DYNA and a Cohesive Crack Model [20]. In this
model, crack propagation was restricted to a known surface and the presumed
path was modeled as an adhesive interface initially tying both surfaces (laminae)
together. Coincident nodes were used along the interface and the problem was
simply reduced to determining a force-displacement relationship between nodes

as follows.

Pi = KijD j (1)
In this equation, P is the force, D is displacement, and K is a diagonal matrix of 3
orthogonal springs representing the penalty formulation. This penalty stiffness is
valid until a certain combination of these forces is reached then before
delamination will propagate. The amount of dissipated energy will shrink the size
of the maximum load surface which will in turn reduce the stiffness making it
possible for an elastic unloading. The initial maximum load surface is shown
below.

F=f(P)—1=0 (2)
In this equation P is the vector of adhesive forces, and f is the normalized measure
of P, force. The function f is determined using an experimental program to assure
it approaches unity as the load state reaches its limiting value. The work is
validated by using a basic delamination model of two steel cubes, and doing
delamination simulations for all 3 modes: double cantilever beam (Model I), end
notch flexural (Mode 2), and mixed mode bending (Mode 3) [20].

Another approach used a cohesive zone model that allows for displacement
jumps across a crack represented by extra degrees of freedom at the coincident
interface nodes. This model allows for a path of discontinuity that is independent
of the mesh structure. There are no restrictions on the type of finite element used
and the jumps are continuous across the element boundaries making it more
versatile for different meshes [21].

Cohesive zone models can also be used in conjunction with other existing

models to simulate different material behaviors. In 2000, Espinosa et al used a

cohesive zone model in conjunction with an anisotropic Viscoplastic model to
simulate delamination in glass fiber reinforced plastic (GRP) composites. The
large deformations in the material are described with the Viscoplastic model in
Lagrangian coordinates with its coefficients determined experimentally. The
interface consisted of zero thickness interface elements that were embedded
between the laminae. The tensile and shear fractions of these elements were
calculated from interface cohesive laws that described these tractions in terms of
normal and tangential displacement jumps. Once the displacements reach a
specified critical value, the interface elements fail, and delamination takes place
[22].

In recent years, a two-parameter cohesive zone model has emerged. In
2003, this type of model was coupled with finite element analysis to study fracture
of fiber composites and adhesively bonded joints. This model was based on
fracture models and criteria stating that two parameters exist that describe
fracture, Le. a critical limiting value of maximum stress, omax, and fracture energy,
Gc. The aim was to determine if the maximum stress values was unique for a
given problem or if there was any physical significance to this parameter. Both
analytical and finite element analysis was done to validate this model [23].

In 2004, Borg et al used a cohesive zone based delamination model and
applied it to a low velocity transverse impact simulation. The purpose was to
model delamination initiation and growth with this specific type of deformation.
This model was simulated in LS-DYNA and applied to a glass fiber/epoxy three

layer cross-ply laminate. The composite was subjected to a non-penetrating

transverse impact using a spherical impactor. In addition, a simulation was done
on an eight layer laminate of Carbon fiber. This study improved upon existing
models by making it possible to determine the local fiber orientations above and
below the interface relative to the propagation direction of the delamination. In
addition, an interpolation scheme is used to determine arbitrary energy release
rates from a finite number of experiments. Although this was a 2D model, it has
very promising results for use in a 3D application [24].

Overall, the cohesive zone model has been shown to be in good agreement
with experimental data in all cases. This model is an easier model to implement,
and does not require the addition of interface elements since it can also be used
with coincident nodes between laminae. It also has shown the ability to be
combined with different material models. Disadvantages include the determination
of the specific properties used in this model for a specific laminate, and overall
implementation in ABAQUS/Explicit, which would require the development of a
user subroutine.

1.2.2 - Continuum Damage Model/Penalty Technique
In the continuum damage model a composite laminate is made up of

stacked homogeneous isotropic layers. The damage surface is based on fracture
mechanics - assuming initial defects/cracks and a stress based failure criterion -
which predicts initiation of delamination. A zero thickness layer is used between
each lamina that unifies initiation and propagation of the delamination. The
interlaminar direct and shear stresses are related to the displacements 8 by the

constraints of the penalty stiffness of the interface k as follows [16].

7i = k?6. <3)
The stiffness is initially very large and displacements are zero to simulate a

perfect bond, which degrades as damage increases. The following equation was

used to estimate the initial stiffness values.
kio = kiTic (4)

In this equation 11c, 12c, and 1:9,c are the interlaminar tensile and shear strengths and
k, is approximately 105 ~ 107 mm". The maximum displacement before
delamination is very small and usually below 10'5 ~ 10‘7 mm. As the loading
increases the delamination damage initiates and slowly develops at the interface.

In order to represent micro-defects in the context of continuum damage
mechanics, a damage parameter was required. The effective properties of the
composite, i.e. stiffness, can then be represented using this damage parameter.
The parameter, I», is 0 when the composite is in an undamaged state and 1 when

it is in a fully delaminated state. The tractions are then represented as follows.
_ 0
Ti — k, (1- (0)5. (5)

This gives the constitutive law for an interface with delamination damage. As the
effective stiffness decreases damage grows. The growth is driven by interlaminar
stresses according to the appropriate damage law described next.
For the damage surface, conventional stress-based failure criterion is used
to predict the initiation of the damage.
1:s(1'i)‘1=0 (6)

Here 13 is a failure function which was chosen using appropriate criterion.

10

2
f5 =l+é+fi ift1>0 (7)

2
71c 72c 75c
Once the delamination exists, fracture mechanics predicts propagation using the

failure criterion fg shown below.

rg(G,)—1 =0 (8)

Here, Gr are the energy release rates corresponding to each fracture mode (I, II,

III). The function i9 is chosen to have the following form.

a r3 r
_ i E _§LII_
fg—[Grc] +[GIIC] +[GIIIIC] (9)

Where, 6., G... and G... are the critical energy release rates in Mode I, II, and Ill,
respectively. This completes a softening interfacial constitutive relationship in
which the shrinkage rate can be controlled by a function (p that is a term in the
damage surface equation.
Flrl.Gi)=fs(ri)-l1—¢(fg)l=o (10)

The damage surface is in the Ti and G space. The function (p is monotonically
increasing and satisfies the conditions <p(0) = 0 and cp(1) = 1. This allows the
damage surface to shrink in the stress space with increasing damage. The
damage is irreversible and non-linear. Penetration is also taken into account and if
it does occur, the stiffness is returned to the initial value and the interface will not
degrade further [16].

This method was coupled with the Lagrange multiplier method in 2003 by

Bruno et al. The purpose was to enforce interface displacement continuity

11

between the plate elements and to determine an accurate and direct approach to
evaluate the direct energy release rates [25]. As before, the initial adhesion of
both layers before delamination was simulated by treating interface stiffnesses as
penalty parameters.

This created an accurate laminate model using multilayer shear deformable
plate modeling coupled with interface elements. These elements were used to
simulate mixed mode delamination in composites.

1.2.3 - Rate Dependent Models
There were two different rate dependent interface models found in the

literature used for composites: a Viscoplastic model and a time-dependent elastic
damage model. They were used to simulate interface degradation in polymer
matrix composites, which includes debonding and rate-dependent delamination.
The basic interface model for both consisted of a relationship between the
traction vector t relevant to a surface I‘, which divided a solid into two parts, 9+
and 52-. The vector of displacement jumps between the two parts was
represented by [u]. The laminate is represented by a stacking sequence of
layers and interfaces. The interface is a zero-thickness layer which transfers
stress from one layer to the next [26].

The Viscoplastic interface model wasgoverned by the following equations.
[U] = [Ule + IUI"p (11)
t = K[u]° where K = diag(Ki), i=1 ,2,3 (12)

N af(t,/l)

[GIVP = nut/1)), a, (13)

12

 

f(t,,1) = Jajt12+ aztg + a3(t3)+2 — 1 + M (14)

T
,1: L‘flmr’l’] [[u]Vp]dr' (15)

(.)+ = .11. > 0(a)+ = ort- s o (16)

 

Equation 11 shows the addition of the elastic and Viscoplastic discontinuities.
The elastic portion is governed by equation 12 where K is the elastic interface
stiffnesses that relate the traction, t, to the elastic discontinuities. Equation 13
governs the evolution of the Viscoplastic displacement discontinuities through a
Perzyna kind law, where f (LA) is a Viscoplastic potential defined by equation 14.
Here, the a, terms are directly related to the values of tractions at the beginning
of inelastic behavior, while h relates to the hardening/softening in the model.
When h is positive, a softening effect takes place, and the interface degrades.
Equation 14 defines the Viscoplastic multiplier and contains 1' which is time. In
both equations, N is a model parameter, A is a Viscoplastic multiplier and t; are
tractions. The final equation, 16, simply shows the meaning of the symbol used
for the positive part of what it contains. In equation 14, the positive part of the
traction normal to the interface is used to keep Viscoplastic displacement
discontinuities from developing in the normal direction due to interface
compression.

The main characteristics of the Viscoplastic model are that it has
irreversible Viscoplastic displacement discontinuities, an evolution law of the

Perzyna kind, a unilateral effect in the normal direction, and a softening behavior

13

to simulate interface degradation. The model depends on nine parameters as

follows where F, L, and T denote force, length, and time, respectively.

Ki[F/L3], 1:123 (17)
a,[L4 /F2], i=1 ,2,3 (18)
1111/14 (19)
rlF/ LT] (20)
N (21)

The time-dependent damage model was an anisotropic time-dependent
elastic-damage model with 3 damage variables and is governed by the following

equafions.

t1 = (1'01 )K1[U]1. t2 = (1'02 )K2IU]2 (22)
t3 = (1'03)K§ (Ma) + + K§< [U13 )_ (23)

1 1 1 2
Y1 = 2 K1IUI12 . Y2 = 2 Kzlulg. Y3 = 2 K3+([U]3)+ (24)

 

A A A A A

f(Y,A)= a1Y1+a2Y2+a3Y3—1+h)\ (25)
N

O A A

D, = yi<f(Y,A)> , i=1 ,2,3 (26)
+

o 3 '2

,1 = Di (27)

i=1

Equations 22 and 23 show the elastic-damage behavior of the interface where K)

are the initial interface stiffnesses and D are the three damage variables for each

14

mode. The differences in equation 23 for the normal direction are to introduce
the elastic behavior in the model depending on the sign of the normal
displacement discontinuities. The conjugate to the damage variables, Y,, are
shown in equation 24. Equations 25 — 26 show the evolution of the three
damage variables. The first equation, 25, shows the potential function or the
amount of accumulated damage, and is dependent on Y and it. Equation 26

defines the damage variable rates as a function of the positive part of the

A
potential function and the material parameters N and 71- Therefore, a softening

effect is introduced if and only if the damage variables always increase.

The damage interface model also depends on nine parameters as shown

below.
Ki[F/L3], i=1 ,2,3 (28)
QirL/F], i=1 ,2,3 (29)
G (30)
;[1/T] (31)
N I (32)

In summary, the first interface model is Viscoplastic with irreversible
Viscoplastic displacement discontinuities, an evolution law of the Perzyna kind, a
unilateral effect in the normal direction, and a softening behavior to simulate the
degradation of an interface. The second model is a time-dependent elastic

damage model, anisotropic, contains three damage variables, and also a

15

unilateral effect in the normal direction. This model simulates interface
degradation by using an evolution law of damage variables that is governed by a
potential function.

The main advantage to using one of these rate-dependent interface
models is reduced computing time due to improved time step integration and
element formulation. The disadvantages, however, are the extensive number of
parameters needed for both models. The problem is in parameter identification
for a particular interface. In addition, extensive experimentation must be done to
show if the delamination occurring with a particular material is indeed time
dependent. Therefore, it has good potential, but it is also hard to implement.

Overall, the models discussed in this literature review were very complex
and required many input parameters. In some cases, there was no clear method
on how to determine these inputs. Therefore it was established that the use of a
less complex model would be beneficial in this study. The models used have
input parameters that are common material properties, and were either
determined experimentally in this study, or can be in future work.

1.2.4 Interface Models
There were two interface models used in this study, partially because of

the simplicity to implement. The first model was a linear softening interface. The
linear softening was a built in function in ABAQUS and is defined by a contact
pressure-overclosure relationship [28]. There are two main types of contact
pressure-overclosure relationships. The first is “hard” contact that minimizes

penetration of slave nodes into the master surface and blocks the transfer of

16

tensile stress across the interface. The second “softened” contact relationship
used a contact pressure that is either a linear function, exponential function, or a
piecewise linear (tabular) function of the clearance between the surfaces. In
addition, a viscous damping relationship can be defined to further affect the
interface behavior.

The contact pressure-overclosure relationship chosen was of the
“softened” form and was prescribed by a linear law. This “softened” contact
pressure-overclosure relationship can be used to model a soft, thin layer on
either one or both surfaces. In ABAQUS/Explicit the contact can be enforced
using a kinematic or penalty constraint method. The kinematic predictor/correct
contact algorithm doesn’t allow any penetrations. Therefore, the penalty contact
algorithm is often chosen. Although the enforcement of contact constraints is
weaker, it allows for a more general treatment of different types of contact. In
addition, with the use of the penalty constraint method, as the contact stiffness
increases the time increment decreases. With penalty enforcement, which was
used here, the contact collisions are elastic except in the case a contact damping
is specified. If this is the case, some energy will be absorbed by the impact.
This energy tends to increase as the contact stiffness increases. The damping is
defined by a critical damping fraction, which is a unitless damping coefficient in
terms of the fraction of critical damping associated with a particular contact

stiffness. The damping forces are then calculated as follows.

rvd = 110,14mkevg, (33)

17

Here, the damping force, fvd , is a function of the critical damping fraction, (to,
nodal mass, m, nodal contact stiffness, kg, and the relative rate of motion

between the two surfaces, er1-

In using the linear pressure-overclosure relationship, the user must define

the slope, k, of the relationship. A general graph representing this relationship is

shown below.

 

 

Pressure, p
X

 

 

 

fl

Overclosure, h

 

 

 

Figure 1. Linear function of pressure-overclosure relationship.

An approximation for the slope of the pressure-overclosure relationship is given
in equation 34.

k : EInterface (34)
tint erface

In this equation, Eimedace is the Young’s Modulus of the matrix material, and

tintedace is the thickness of the interface.

18

An additional element to this linear softening model is friction. Friction is
used in this case to prevent slipping regardless of contact pressure. The friction
is characterized as rough and once two surfaces contact and undergo rough
friction, they should remain closed. This is primarily used with the surface
behavior definitions with the no separation option to restrict motions normal to the
surfaces. The shear traction slope for this option must also be defined as

follows.

ShearTractionSlope = M (35)
tint erface

This parameter is used to define a tangential softening in ABAQUS/Explicit.

The second interface model used in this study was a plasticity model.
Unlike the linear softening model, it is not a built in function in ABAQUS. Use of
this model requires a user subroutine to define the interfacial constitutive
relationship. In ABAQUS/Explicit this is known as a VUINTER. The particular
user subroutine modeled a uniform thickness interface bonded to both surfaces.
This interface material is characterized by uniaxial plasticity in the normal
direction with linear hardening, while the shear behavior continues to be elastic.
It is also important to note that membrane straining of the interface does not
affect the stress being transmitted across the interface. In addition, heat transfer
can be incorporated in the form of conductance that is independent of the
pressure or gap distance.

The plasticity model requires 7 parameters as shown in Table 1.

19

Table 1. Plasticity model properties of the interface.

 

 

 

 

 

 

 

 

 

 

Gap Cut Off 0.000001 m
Young’s Modulus 1.797 GPa
Poisson’s Ratio 0.42
Initial Yield Stress 16.4 GPa
Hardening Modulus 0.638 GPa
Thickness 0.00001 m
Conductivity 0.1W/m-K

 

The parameters determined through experimentation were the Young’s Modulus
and Poisson’s Ratio. The remaining parameters such as initial yield stress,
hardening modulus, and conductivity were taken from material databases which
are a collection of average material properties for various polypropylene types.
The data was interpolated to find approximate values for the particular
polypropylene used in this study.

The subroutine begins by zeroing and reading in the needed values. Then
it is determined whether any two nodes directly across the interface have passed
the gap cut off limit. The shear modulus is then computed using the Young’s
modulus, E, and the Poisson’s ratio, v. This value could just be read in from the
material properties, but for simplicity, the program was left as originally written.

_ E
_2(1+v)

 

(36)

If the limit has been. crossed, the interface at that point is considered separated
or debonded, and the stress and heat flux are set to zero, therefore no additional
stress or heat transfer can occur. If the two nodes are still bonded, the nominal
strain increment, OE, and elastic trial stress, at, are computed as follows. In these

equations t is the thickness of the interface and E is the Young’s modulus of the

20

interface and At is the relative displacement between the two surfaces. The
indices 1, 2, and 3 correspond to 3 dimensions.

At1 - A12

d5“ = T191912 =—-t--.d813 = —— M3

t (37)

“11:011+E*d€11
012 = 012 +(3"‘(11‘312 (38)
013 = 013 +G*d€13

The next step is to determine how much to scale the normal stress, 011 due to
any yielding of the interface, If the absolute value of the elastic trial stress for the
normal component is greater than the yield stress, (IV of the interface then the
plastic strain increment is calculated and is shown in equation 39 where E, is the
hardening modulus of the interface.

39
E+Eh ( )

01=l011|1 d5P1 =

 

The change in yield stress is then calculated using this plastic strain increment.
0y = 0y + Eh * d€p1 (40)
In addition, the corrected normal stress is calculated as follows.
011=Si90(0y.011) (41)
If an initial temperature of the material is specified in the simulation, along

with a thermal conductivity, k, the heat fluxes are calculated using the

temperature of the master node, Tmaster, and the temperature of the slave node,
Tslave-

flux = kA Easter t‘ Tslave (42)

21

To conclude the subroutine, if the temperature the particular node being
analyzed is equal to or less than zero, the flux is calculated as follows where the

master node temperature is taken as that of the first master node.

flux : kA M (43)

22

Chmwv2

BACKGROUND

In this study collaboration was done with the Composite Materials and
Structures Center at Michigan State University. Researchers in that group have
extensive experience fabricating and testing many types of natural fiber
composites and expressed an interest in testing the formability of these materials.
The material recommended was a Kenaf-Polypropylene composite.

2.1 Fabrication of Kenaf/Polypropylene Composite

In order to fabricate the composites needed for use in the stamping
machine, many approaches were used. To begin, the fiber is baked in a vacuum
oven for approximately 6 - 8 hours. The fiber is either chopped into pieces 0.197
— 0.315 inches (5 — 8 mm) in length or left long depending on the application.
There were two pressing machines used in this study: the Carver Laboratory Press
Model SP-F6030 and the Tetrahedron manufactured in San Diego, CA. The
Carver press was used for fabrication of smaller composites, while the
Tetrahedron had a larger platen and computerized, making it easier to fabricate

the larger composites.

23

 

(a) (b)

Figure 2. (a) Carver Press (b) Tetrahedron Press.

The first approach consisted of fabricating 2 — 0.022 lb (10 9) PP sheets
and sandwiching 0.019 lb (8.5 g) of long and chopped Kenaf fibers between these
sheets. A 5” x 7” x 0.039” inch (127 x 177.8 x 1 mm) stencil was used to ensure
proper sizing of the sample.

The second approach was to premix chopped Kenaf fibers and PP powder
in a kitchen mixer, prior to pressing. It was observed during the mixing process
that most of the PP powder would not distribute evenly through the fiber, but tend
to fall to the bottom of the bowl. An attempt was made to do hand mixing in order
to more evenly distribute the fiber with little improvement. In order to preserve the
30% fiber composition of the composite, the same amount of fiber and matrix
materials was used. The resulting composites from both approaches are shown in

Figures 3 - 5.

24

 

.llIltmi

Figure 4. Chopped Kenaf fibers sandwiched between two PP sheets.

 

Figure 5. Premixed chopped Kenaf fibers and PP.

25

As shown, these methods created many voided areas. In the final attempt
to produce a more even fiber distribution the PP was changed to the Microfine
Polyolefin Powder and a multilayered sprinkling and sifting method was used.

The new process is outlined in detail in the following flow chart.

 

Kenaf/Polypropylene Composite
Fabrication Process

Bake kenaf fiber>
in vacuum oven,
6- 8' hours
Weigh/Chop Weigh/Mix Preheat
Kenaf Fiber PP and Press to
(5 — 8 mm) Coupling 1900 (375F)
Agent
Top bottom steel
plate with non-stick
paper, then stencil

Qverinq Process
1) Sprinkle layer of PP mixture to coat bottom
2) Sprinkle half of chopped fiber randomly
3) Sprinkle layer of PP mixture in center
4) Sprinkle half of chopped fiber randomly
5) Sprinkle layer of PP Imixture to coat top

Cover with non-
stick paper,
then steel plate

Pressing Process
1) 15 minutes, 0 lbs, 1900
2) 10 minutes, 10 kip, 1900
3) 5 minutes, 24 kip, 1900
4) Cool to 100 0 and remove

I
Remove composite
from stencil

Figure 6. Multilayered composite fabrication process.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

26

The number of layers in the final composites fabricated for material
characterization and stamping was contingent upon how many layers were
needed in the completed laminate. For example, two layers would consist of
three layers of PP, and two layers of fiber. The composite was then pressed

using the same process and the resulting composites are shown below.

lllwwltillrr

.1 I

I

 

IIII‘ ..

Figure 7. Multilayered fabricated composite.

It is important to note that the preceding samples were all done with
unbaked raw fiber with no surface treatment. In order to make useable samples
for stamp thermo-forming process, a coupling agent was introduced in the form
of a powder added to the PP powder before fabrication. The coupling agent was
Epolene Wax G-3015P from the Eastman Chemical Company in Kingsport,
Tennessee. The purpose of the coupling agent is to make the fiber more
compatible with the matrix and to provide better adhesion of the matrix to the
fiber.

Another improvement to the composites was seen when the fiber was pre-

baked in a vacuum oven at 86°F (30°C) for a minimum of 4 hours. It was found

27

that after baking the fiber would lose approximately 30% of its mass by weight.
The extraction of this moisture reduced some of the kink in the fiber, and made it
straighter for better chopping. In addition, it allowed for better adhesion between
the fiber and the matrix, which can be damaged if water is present. After baking,
the fibers, which were originally chopped using scissors, were chopped using a
paper cutter. This allowed for better chopping consistency that produced a more
even length of fiber.

For the stamping process, a 6.5-7” (165.1 — 177.8 mm) diameter sample
was created by using the 8” (203.2 mm) square stencil and using a band saw to
cut the circles, followed by sanding the edges with a belt sander to ensure that
the edges did not get caught in the die during forming. The individual mass
amounts calculated for the 8” inch (203.2 mm) square stencil of 0.079 inch (2
mm) thickness was calculated to be 0.154 lb (69.83 g) of PP, 0.069 lb (31.27 g)
of fiber, and 0.007 lb (3.13 g) of coupling agent. For the 12” square stencil of 3
mm thickness, which was used for cutting all material samples for testing, the
masses were calculated to be 0.505 lb (229.14 9) of PP and 0.231 lb (104.7 g) of

fiber, with 0.023 lb (10.5 g) of the coupling agent.

28

2.2 Material Characterization - Directionality

To begin characterization of the composite, the first test conducted was a
squeeze flow test [27]. The purpose of this test was to determine preferred fiber
orientations (in angles) to be used in the tensile testing and simulation. This was
done by constructing 3 large 12” double layered composites, which were then cut
into 3” circles using a drill press. The zero direction was marked on each sample
and it corresponded to the zero direction on the stencil when each composite
layer was fabricated. The center was also clearly marked for identification. A

sample is shown below.

 

Figure 8. Squeeze flow test sample before testing.
A schematic of the squeeze flow experimental setup is shown in Figure 9.
Prior to testing, the presses were heated to 375 °F (190 °C). The squeeze flow
samples were placed between metal plates and then placed into the presses
under nominal pressure and allowed to heat for forty minutes. After that time, the
platens were compressed to a pressure of 12,000 lbs (53.4 kN) for sixty seconds.

The samples were then removed and allowed to cool in the metal plates.

29

.UUUII.

"111.1111...‘ 111mm mm . 1.111113 "111111.811 11111115111111-..11111..1111t.111V
Il' ..II.I'1" (IIIII' .III'III‘AIIII" TIIIITZI I'uII". IIIII’TTIIIWIHIIIVFIIIII‘TEmit—11W" ITI'SIT TI"

 

 

     

    

Metal

Heated
Plates

Platens

 

 

 

 

 

 

 

 

 

Fiber Reinforced Sample

Figure 9. Schematic of squeeze flow test experimental setup [27].
After all 24 samples were pressed; they were labeled in increments of 10
degrees going completely around the sample starting with the zero direction as
shown in Figure 10. This was done so that measurements could be taken at

each increment to determine in which directions the fibers tended to spread the

most.

 

Figure 10. Squeeze flow test sample after testing.

30

Each sample was then analyzed in two different ways. First, at each
angle increment a measurement in millimeters was taken from the center to the
edge. Second, the same measurement was taken from the center but excluding
the edge areas with large PP concentrations and no fiber. The preferred fiber
orientations correspond to the peaks or valleys in the graphs of the
measurements. The initial four angles chosen were 40°, 170°, 250°, and 340°
degrees. In order to confirm the accuracy of these results; a few of the samples
(i.e. 2, 3, 9, 10, 17, and 22) were taken out due to the fiber clumping during testing.
The normalized charts for all samples and only the best samples are shown below
with the preferred fiber orientations represented by the vertical lines on each

graph.

Normalized Squeeze Flow Data - All Samples

2.45

 

340

2.4 1

 

2.35

 

2.3 3 ,.

RIRo

2.25 J '-
2.2 f
2.15 7

 

 

 

 

2.1

 

 

Angle (Degrees)

Figure 11. Normalized length versus angle for all samples.

31

Normalized Squeeze Flow Data - Best
Samples

 

 

 

 

 

 

 

 

Angle (Degrees)

Figure 12. Normalized length versus angle for best samples.

32

2.3 Material Characterization - Properties
In this study, material characterization has two primary parts —

characterization of the composite and characterization of the matrix/interface
material — polypropylene. To begin, testing was done at room temperature to
determine the tensile properties, i.e. Young’s Modulus and Poisson’s ratio of the
material. The tests were performed on the UTS Machine, Model SFM 20 of
United Calibration Corporation, with the twin screw loaded frame.

The standard used was ASTM Standard D638. The required minimum
length and width of each specimen was 6.5” (165.1 mm) and 0.5” (127 mm),
respectively. Two of the preferred fiber orientation angles were used in testing —
40° and 170° degrees, since the two preferred fiber orientation VUMAT was
initially chosen to be used in the simulation. Test specimens were cut in each of
these directions for testing. In addition, a sheet of PP with Epolene (coupling
agent) was made at the same volume ratio as the composite, i.e. 3% coupling

agent. A sample of each is shown in Figure 13.

 

Figure 13. Kenaf/PP and PP tensile specimens.

33

The testing consisted of two tests for each type of sample. The first
procedure exerted a tensile force on the specimens until breakage. This was
used to determine the Young’s Modulus and other important values. The
resulting stress strain curves are shown in Figures 14 - 16. In addition, the

average values for the Young’s Modulus and other values are shown in Table 2.

Stress-Strain Curve for Kenaf/PP Composites - 40
Degrees

 

7000
6000

g 5000 /

9: 4000

§ 3000 //

2000

1000 [X

If

0 . .
0 0.5 1 1 .5
Extension (%)

Figure 14. Stress strain curve for 40° fiber orientation specimens at room
temperature.

 

 

 

 

Str

 

 

 

 

 

Stress-Strain Curve for Kenaf/PP Composites - 170
Degrees

9000
8000 -1--- ~-— ——— .—— e— a m _...-__._. ._ 2 _ EL _

 

 

 

 

 

Extension (%)

Figure 15. Stress strain curve for 170° fiber orientation specimens at room
temperature.

34

Stress-Strain Curve for Polypropylene/Epolene

4500

 

4000-l

3500..m _

Stress (psi)
5

30001 -
2500.._._,-

1500- -~w ~

1000-.- —

500 -
O

 

 

 

 

1.5

Extension (%)

2 2.5

Figure 16. Stress strain curve for Polypropylene/Epolene specimens at room
temperature.

It is important to note that from these curves, the behavior of the material is

not linear and appears to have some plasticity or viscous properties. For this

study, the material model was predetermined. In the future the material behavior

could play an important part in the future material model chosen for simulation. In

addition, the Young’s Modulus was calculated using the extension range of

approximately 0.25 — 0.85 %.

Table 2. Tensile testing results for all specimens — Young’s Modulus.

 

 

 

 

 

 

 

 

 

 

Spfigge" T3223” $132331 ”5133111? Elanrsztlion 211332113...
(psi) (1981) (%) CA)
409 519.9325 8595.15 592134.4 1.435 1.42
1709 512.88 7263.87 701258.25 1.594 1.57
PP 223.08 4017.468 2630366 2.752 2.748

 

In the second procedure, the specimens were only allowed to extend 0.5

in. (12.7 mm).

35

This was done to protect the strain gauges being used to

 

determine the Poisson’s Ratio of the composites and polypropylene. The results

are shown in Table 3.

Table 3. Tensile testing results for all specimens - Poisson’s Ratio.

 

 

 

 

 

 

 

 

 

 

Specimen Poisson’s Tensile Tensile Modulus 0f Break
T 9 Ratio (lbs) Strength Elasticity Elongation
VP (lbs) (psi) (%)
409 0.35 496.552 6825.21 773300 0.538
170° 0.4 507.29 5676.35 659594 0.528
PP/Epolene 0.42 223.132 3961 .54 31 3886 0.538

 

 

Testing showed a large Poisson’s ratio for the PP/Epolene mixture. The
normal value for polypropylene without a coupling agent would be around 0.1.
This further shows that adhesive strength of the coupling agent. It not only binds
the fiber and matrix in the composite, it also adds strength to the matrix material
alone.

In addition to tensile testing, a flex test was done using ASTM standard
D790 to determine flexural properties and other properties associated with
composite performance. This test utilizes a three-point loading system applied to
a simply supported beam that subjects the composites to a mixture of tension,
compression, and shear forces. The specimens used for these tests were very
small with minimum dimensions of 0.5 inches (12.7 mm) and 2.5 inches (63.5
mm) for width and length, respectively. The minimum thickness of each sample

was approximately 0.14 inches (3.5 mm). The results are shown in Table 4.

36

Table 4. Flex testing results for all specimens — Modulus of Elasticity.

 

Bend

 

 

 

 

 

 

 

 

 

Specimen Width Thickness Max Force MOE
. . Strength .
ID (In) (In) (lbs) (psi) (ps1)
409 0.536 0.1424 25.8 8073.4 5268148
1 709 0.5356 0.147 30 8823.6 591 062.6

 

Completed testing gave accurate values of the Young’s Modulus and
Poisson’s Ratio’s for two of four preferred fiber orientations, 40° and 170°. The
values for the final two preferred fiber orientations, 250° and 340°, were taken as
an average of those determined experimentally. Shear modulus values for the
fiber were approximated as follows since proper testing methods were not
available.

6: E
2(1+v)

 

(44)

In addition, the Poisson’s ratio and shear modulus for the polypropylene were

taken from previous work in literature. Values for the transverse shear moduli,

G13 and 023, were also taken to be 65 kPa [27].

Table 5. Material properties for simulation.

 

 

 

 

 

 

 

 

 

 

 

 

 

Name Kenaf/PP Name Kenaf/PP
Er, 40° 4.153 GPa V1, 340° 0.375
Ef1170° 4.835 GPa Vm, a_ll_ 0.42

E1, 250° 4.495 GPa 611,400 1.538 GPa
Er, 340° 4.495 GPa Gt, 170° 1.727 GPa
Emfl 1.797 GPa (5., 250. 1.63 GPa
Vf, 40° 0.35 Gf1340° 1.63 GPa
Vt, 170° 0.40 GM 0.2 MP8
Vf, 250° 0.375 Vf 0.384

 

 

37

Table 4. Flex testing results for all specimens — Modulus of Elasticity.

 

 

 

Specimen Width Thickness Max Force 83:23,} MOE
lD (In) (In) (lbs) (psi) (pSI)
40° 0.536 0.1424 25.8 8073.4 5268148
170° 0.5356 0.147 30 8823.6 591 062.6

 

 

 

 

 

 

 

 

Completed testing gave accurate values of the Young’s Modulus and
Poisson’s Ratio’s for two of four preferred fiber orientations, 40° and 170°. The
values for the final two preferred fiber orientations, 250° and 340°, were taken as
an average of those determined experimentally. Shear modulus values for the
fiber were approximated as follows since proper testing methods were not
available.

6: E
2(1+v)

 

(44)

In addition, the Poisson’s ratio and shear modulus for the polypropylene were
taken from previous work in literature. Values for the transverse shear moduli,

G13 and 023, were also taken to be 65 kPa [27].

Table 5. Material properties for simulation.

 

 

 

 

 

 

 

 

 

 

 

 

 

Name Kenaf/PP Name Kenaf/PP
E1, 40o 4.153 GPa v1, 340° 0.375
E14709 4.835 GPa Vm, 3L 0.42

Er, 250° 4.495 GPa Gt 40° 1.538 GPa
Ef, 340° 4.495 GPa GI. 170° 1.727 GPa
Em, a_|| 1.797 GPa Gt, 250° 1.63 GPa
v1, «,0 0.35 6., 3400 1.63 GPa
Vf, 170° 0.40 Gm, 21' 0.2 MP3
Vf, 250° 0.375 V1 0.384

 

 

37

2.4 Compression Moulding Process and Composite Comparison
In order to evaluate the effectiveness of the Kenaf/PP composites and the

fabrication process previously outlined, a comparison was done to other natural
fiber composites. To begin the Kenaf/PP composites made in this study were
compared to other compression moulded natural fiber-polypropylene composites
with 40% wt. fiber in Figure 17. These composites were fabricated using natural
fibers spread between polypropylene films [2]. The Kenaf/PP composites
fabricated using the fabrication process used in this study outperformed all fibers
with the exclusion of hemp. The properties of the composites improved by
introducing the optimized fabrication technique and the addition of the coupling
agent, which proves that better matrix fiber adhesion was achieved allowing the
composite to reach its full capabilities. In addition, the Kenaf/PP composites
fabricated in this study showed a tensile strength of within 725.10 psi (5 MPa) of
the Hemp/PP composites found in literature. This was a good achievement
considering the fact that hemp fiber alone has a tensile strength of 79.77 —
130.53 ksi (550-900 MPa) while Kenaf ranges from only 41.19 — 116.03 ksi (284-

800 MPa) [2,12].

38

Tensile Strength Comparison of Compression
Moulded Composites

0')
O

 

woos-01
COCO

Tensile Strength (MPa)
8

 

O

Kenaf/PP Kenaf/PP Coir/PP Hemp/PPSlsaI/PP
(30%) (40%) (40%) (40%) (40%)

Figure 17. Tensile strength comparison.

The flexural strength of the Kenaf/PP composites in this study were also
compared with previous results from other compression moulding studies [2].
The new process has produced nearly double the flexural strength of some other
natural fiber composites, including other Kenaf/PP composites. The exception to
this gain was once again the Hemp/PP composites. Yet, it is believed that
Hemp/PP composites produced with the new technique would outperform these

layered hemp composites [2].

Flexural Strength Comparison of
Compression Moulded Composites

 

 

Kenaf/P Kenaf/PP Coir/PP Hemp/PP Sisal/PP
P (30%) (40%) (40%) (40%) (40%)

Figure 18. Flexural strength comparison.

39

A comparison was also done between the compression moulding process
and other manufacturing processes that are currently in use to fabricate
composites. Figure 19 shows a comparison of the new compression moulded
Kenaf/PP composites to natural fiber composites using resin transfer moulding
[13], and an extrusion and injection moulding process [14]. Although fiber types
and weights differ between the studies, it can be concluded that the new
optimized compression moulding process is a very competitive manufacturing

p I'OCGSS .

Flexural Strength Comparison of Various
Manufacturing Processes

 

 

0103\103
0000
l
l
l
l

Flexural Strength (MPa)
.5.
O

 

 

 

4ND»)
GOOD
I
l

 

 

Kenaf/PP (30%) Natural Fiber Hemp/PP (30%
Compression Composites (20.6%) Extrusion and
Moulding Resin Transfer Injection

Moulding Moulding

Figure 19. Flexural strength comparison for various manufacturing processes.
Using the elastic modulus data from the material testing, it is possible to
compare the benefits of using this Kenaf composite over other natural fibers as
well as E-glass. In Figures 20 and 21, it is shown that the Kenaf/PP composites
in this study have a higher Modulus/Cost and Specific Modulus than sisal, coir,
and E—glass fibers. Therefore, not only has Kenaf been proven to be a viable

replacement for more expensive and non-biodegradable glass fibers, the new

optimized fabrication process has been proven to produce competitive

composites to those already in use.

Comparison of Modulus Per Cost for Various
Fibers

 

Modulus/Cost
(E-Modulusl($lkg))

 

 

Kenaf Sisal Coir E-Glass

Figure 20. Modulus per cost comparison.

Comparison of Specific Modulus for Various
Fibers

 

Specific Modulus
(E-Modulus/Denslty)

 

Kenaf Hemp Sisal Coir E-Glass

Figure 21. Specific modulus comparison.

41

Chapter 3

EXPERIMENTAL WORK

3.1 Introduction to Stamp Thermoforming
Experimental stamping simulations were used in evaluating the accuracy of

the numerical results. The machine used was a stamp thermo-hydroforming press

pictured below.

 

Figure 22. Stamp thermo-hydroforming press.

The press consisted of a 4 in. (101.6 mm) hemispherical punch to form the parts
with a female die. Nothing was used to constrict movement of the sheets as they
were formed in order to see the wrinkling behavior. A gap of 0.6 in. (15.24 mm)

was maintained between the blank holder and the die. The samples were drawn

42

to a depth of 1.47 in (37.338 mm) and allowed to partially cool before being

removed. A schematic is also shown to illustrate each part of the press.

Blank
Holder

  

 

Thermoplastic
Sheet

Rigid
Punch

Figure 23. Stamp thermoforming schematic.

  

    

 

 

    

 

 

0.6 inch
Gap

 

 

 

   

Initially the bottom die portion was heated to approximately 350°F (177°C), which
was measured using a thermocouple. The samples were also preheated in an
oven to the same temperature for approximately 20 minutes. This allowed for

better forming and melting of the polypropylene without burning the natural fibers.

43

3.2 Experimental Stamping Results
In order to form parts properly numerous samples were fabricated and

tested. The parameters for stamping previously described were the optimized
values found from trial and error. Initially, it was found that during fabrication, if
any resin rich areas were in the initial samples, a tearing would occur during
forming. In addition, some of the first samples tested were not preheated, which

also led to tearing. A sample of this is shown below.

 

Figure 24. Stamp thermoformed part with tearing - top view.
Therefore, care was taken in fabrication so resin rich areas did not occur. The

following parts were successfully formed using the stamping press.

"'v ‘ "WI‘II‘I‘W “111'"
111““? 1‘1““

   

(a) (b)

Figure 25. Stamp thermoformed part A, (a) top view, (b) side view.

 

"1il 11111111111

(a) (b)

Figure 26. Stamp thermoformed part B, (a) top view, (D) side view.

45

[[ljjllilw,1... I
I : HI'II'J‘ ‘

11111-11 ,
I‘ll

.. ‘,- 111". (.1111.- “Iii"
.' M :1“ WI

11' ”mi

1
l .

       
    
 

. 1“”

1 (.11.
‘IllI 11]

Mill“

 

13“ (Illli

11..I1

...1.1l:‘l‘

(a) (b)

Figure 27. Stamp thermoformed part 0, (a) top view, (D) side view.

It is important to note the similarities amongst all these formed samples. To
begin, the top view of each part shows that the material forming is asymmetric. In
addition, each part consistently formed 7 wrinkles. The wrinkles were not very
round (or wave like), but rectangular from the side view, which can be seen in the
rectangular outline in Figure 27(b). In addition, there are parts in the wrinkling that
exhibit a pinching effect and form triangular corners as shown in the circled outline
in Figure 27(a). Also, no separation of the layers was seen. These behavioral

similarities are what will be used in the comparison to the numerical simulations.

Chapter 4

NUMERICAL ANALYSIS

The complete numerical analysis of the stamp thermoforming process was
done using ABAQUS/Explicit. The material layers were modeled using a multiple
preferred fiber orientations updated material law through a user subroutine,
VUMAT. The interface layer was modeled using two methods as outlined in
Chapter 1: an existing linear softening algorithm, and a plasticity model, both

used in ABAQUS.

4.1 Material Model
The constitutive relationship used in this study was developed by Mike

Zampaloni et al and is based on preferred fiber orientations [27]. These
orientations can be non-orthogonal, and are determined using squeeze flow
testing as outlined in Chapter 2. The approach is to model the fiber-reinforced
material as unidirectional layers. Each layer is considered independently and a
stiffness matrix is created with the summation of all these layers. This model
used an updated material law which took into account the fact a material’s
constitutive relationship changes with deformation.

The stiffness matrix [Q] is determined based on its principal geometrical

axes. This is done for each preferred fiber direction.

Q11 Q12 Q13
[Qij]= Q21 Q22 Q23 (45)
O31 O32 Q33

47

In this equation, the subscript 6 was changed to 3 for simplicity. Since each layer
is very thin, a plane stress assumption was used. The components of the

stiffness matrix can be written in terms of the engineering constants as shown in

equations 46 through 51.
E
Q11 = —1—1—— (45)
(1" 0,2021)
522
= —— (47)
22 (1‘ ”12021)
E
Q12 = _U_12_11_ (48)
(1 - 012021)
033 = G12 (49)
013 = 031 = Q23 = Q32 = 0 (50)
021 = Q12 (51)

In these equations, E11 represents the Young’s Modulus in the 1-direction and
E22 is the Young’s Modulus normal to the 1-direction. The Poisson’s Ratio and
the Shear Modulus in the 1-direction are v12 and G12, respectively. For any other
preferred fiber orientation, the same directions hold. Therefore, E11 is the
Young’s Modulus in the fiber direction, E22 is the Young’s Modulus in the matrix
normal to the fiber direction, and v12 and G12 are the Poisson’s Ratio and Shear
Modulus in the fiber direction.

Since the material properties are measured along the fiber directions, they
must be transformed into the material frame to assure continuity. The

transformation matrix is shown in equation 41, where G is the angle between the

48

material frame and the fiber directions. The m is the cosine of 6 and n is the sine

of 0. The final constitutive relation for a single fiber orientation is shown in

equation 53.
m2 n2 2mn
T(6)= n2 m2 —2mn (52)
2 2

[011: Q21 Q22 Q23 [Elk (53)

 

 

For multiple preferred fiber orientations, a combined stiffness matrix is

created to take into account the behavior in all the directions.
[Gisheet = [511 + [612 + +[61n (54)
[0] = [aisheet [5] (55)
This constitutive relationship was implemented through a user subroutine,

VUMAT, in an ABAQUS/Explicit analysis. Three versions of the VUMAT were

used in this study: 2, 3, and 4 preferred fiber orientations.

49

4.2 Numerical Stamping Results
4.2. 1 — Two Preferred Fiber Orientations
Initially, the two preferred fiber orientation VUMAT was used in

conjunction with both interface models. Many factors were found to affect the
quality and validity of the simulation including mesh size, punch velocity, and
total simulation time. The mesh was reduced from 12800 TRIA elements for
each layer (2 total) to 5712 QUAD elements, and finally to 2723 TRIA elements.
This drastically reduced computation time and file size for each simulation. In
addition to the change in mesh, the velocity curve was changed to a triangular
shape. Since the simulation is a quasi-static process, there were a few things
that needed to be considered.

Applying ABAQUS/Explicit to a quasi-static problem required special
considerations. By definition a static solution is over a long period of time, but
that can be impractical when it comes to simulation. The long time period would
require an excessive number of small time increments. In order to counteract
this dilemma, it is imperative to model the process in the shortest time period in
which the inertial forces remain insignificant [28].

This is achieved by changing loading rates, using smooth amplitude
curves, and even in some instances using mass scaling. The loading rate can be
increased so the event occurs in less time as long as dynamic effects remain
insignificant. A smooth amplitude curve may also be used. If there are sudden
movements they can induce noise which will lead to an inaccurate solution. It

was found in this study that the best amplitude curve was triangular, which

50

allowed for the velocity to steadily increase and then decrease to the finish.
Another important factor is mass scaling. ABAQUS has a built in mass scaling
function that can be used either as *FIXED MASS SCALING or *VARIABLE
MASS SCALING. The same result can be accomplished by artificially increasing
the material density, although care must be taken to not increase the density too
much or an inaccurate solution will also occur. The minimum stable time

increment of all elements can be expressed as follows.

L6

Cd

At (56)

In this equation, L9 is the characteristic element length and cd is the dilation wave
speed of the material, which is given by the following equation where E is the

Young’s Modulus and p is the material density.

’E
= _ 57
Cd ,0 ( )

Therefore, if the density is increased by a factor ’of f2 then the wave speed
decreases by a factor off, and the stable time increment increases by a factor of
f. This leads to fewer increments required to perform the analysis, thus
decreasing the time of the simulation.

In conclusion, under the assumption that a quasi-static analysis in real
time would be virtually the same as a completely static solution, it is often
necessary to change simulation parameters such as loading rate, mass, and
amplitude in order to achieve an accurate simulation. In this study, the velocity

was modeled to ensure that the internal energy leveled off, and the kinetic

51

energy gradually increased then decreased to the end. The total time also
played a role in the simulation, not only for computational time, but to achieve the
best model of the quasi-static process.

Initially, the initial volume fraction was under scrutiny. Using the density of
87.40 Ib/ft3 (1.4g/cc), the volume fraction was calculated to be 0.218. Upon
visual observation of the composite, this number seemed to be very low
compared to the large volume of fiber put into each composite sheet. Therefore,
a density test was done on the fiber. The fiber was put in a glass of water to see
if it floated. Since the fiber did float, it was determined that perhaps the density
was actually less than 62.43 lb/ft3 (1 g/cc), the density of water. Therefore, after
a literature search, it was found that the density of cubed Kenaf fiber is 55.56
Ib/ft3 (.89 g/cc) which gave an initial volume fraction of 0.384 [15]. In order to test
the sensitivity of the initial volume fraction, simulations were done with 0.218,
0.384, and 0.6. In addition, the Poisson’s ratio was taken to be 0.2, 0.3, and 0.4
as another sensitivity test. For this sensitivity analysis, only the linear softening
model was used, which only used two of the material parameters. The results
are shown below in Figures 28 to 30. It is important to note the differences in
each simulation. Very little changes were noticed when varying the Poisson’s

ratio or volume fraction as shown.

52

1
.1‘1‘

. 3311., ..
“9.1111111” ‘V’ ‘

741$ .1

z x
ar-

‘1‘ (lit. 3 31.15;.“
4“)“ - 1

4“"

. - up,

.. 4 II ‘

' “ ‘ I‘m“ a:

_ m'm'n‘r‘ Vfib

i':':':'4'-w- 4 7 "$5519. .6

.‘ :"mu'A .11‘ 4 A': 2.0‘: < .
l A'A'A7‘15') 4 ‘v‘ - ”u 5

Manama; . ”y ‘

an»... ., L." .

xv.»

V I2.

Av. 45"; ‘

r .~.~ 45"“ ‘ v.V-L

. n,~$‘\vflAv‘vA:¢‘unu s
‘ SVAvAVAVA'.. —

 

 
   
 
   

(PW v‘

> VAT ”11‘
, figflyfik 454“ 11‘“

  
      

T8883:
. 4 ., 101:“ ..
.~ 4110,49 a; 1 5"
""5141933114947 7 Al ' ' .:
-

I‘ i”
t." '
415111.44
_ 51‘»‘§‘u.:v .
‘ :5 0'55 ‘ ’ “AV

0

,.

1
"'1" "I
£11111“
1 I ‘
. ‘ ‘g A A
:f’,‘ ‘ I
-=

- <
#3.:
"' '«y‘ <u
Jug.»

‘ V

a .
74h
'4

AV
V \

4» .
”-st S‘VA::

AVQV

    

,..v'.( ,43. j .
= 2:5»: , kfiwfij'fizfil) .

elk
"
., .. 15.2.45
' 4’

    
  

  
 
   

  

  
     

  

  

 

G.
: m;
V) v

V
A
315:.
- 4»
v31¢dfi$
ng‘
. v.vg‘

Figure 28. Deformed blank for v=0.2 (a) Vf=0.218 (b) Vf=0.384 (c) Vf=0.6.
53

  
 
 
 
 
  
   

    
  
  
  
 
  
    

Afllpk‘m ' ‘33.:
‘ ‘,Iii‘dwz§" ' 1
wé‘lfifii‘fi; '“5‘41 12;. .-
01") 1" «9.1555: n

d 4’ " \‘ifinfiu mu
53$: ’ gala! ““':'5:1:=2i}§§5§‘
v1.1

_.-
e A
. «f
A}.
-2“§V‘§)‘:v
‘ \ :v: >

e
F

< u.“
- x
.13»
"

   

 
 

 
   

‘ A a v 'I. ,
A'A'AV‘eA:."9

'9." “:67: 9
Q

' v
TA:

I"

 
 

A

'AVA'

A
"A

  

b

‘1

H

mm

‘ I
1

' A
0'
‘ 1
0
'I
A 1 1‘4"'
I'A'A'A‘

3.
a
.

V
VI"!
1“,
VI'A'I
I.
1:1 v
’17..

u
i
0'
.1
I‘ll!“
I'I'I‘I
y .
1 1‘ .‘

4m:
9

A"

5
1
v‘

u

.1 g
s ..
v “agha',

:V‘VAV‘V)
“9.999

A
.1

1

1
q
’4

("1"
\
‘\

.~’1 I 111311.13,”
11115111118 VK‘V. .,

'4?) (1'1

,111‘72» 4A

12'.- . . ..
'A Vgls ~ ‘ ‘e
”0595“”

1.
. 1., ‘II'
« . .‘ “"5
‘ AV‘V‘ I. .
'n‘v‘v." v j»
, A .v‘ 1 I)
v ¢'$"'A'Anvnl¥‘
"4 :‘ag'gnvanAw‘ 4
a," . ..’.. . h“
an

5 "<5

An- 1' > wp“

.- egg; “<95,“

m r . ..“v.v.v.uv:,,uvg. a
7 - Av‘v.v‘vp , - -

2 ‘A. T
- .3.“
T.

. -. .12).
, _ 4. .,
IL 0’4 981-,
' ¥ ‘S‘wm.
”teem
‘Vflb‘éyfiif-

A"
4 '4 5‘ ‘ ,
nr‘E‘ . v
Evan's': 1»
A'LVA'VAVA 5

 

Figure 29. Deformed blank for v=0.3 (a) Vf=0.218 (b) Vf=0.384 (c) Vf=0.6.

54

  
 
  
 
  
 
  

'1' ' ~' ‘
‘(‘§‘4115:L" .
. , \‘ ‘l‘k‘ln'h RV
5.4.13»; 11, AVA? ‘ “4316835555
v P‘V‘nf‘v 11.1. c)‘ (A ‘v‘v.v_ ‘55 a
A :v‘vflg'hgnl’ 5 “"5 >— 5‘: g
A!" ‘-‘u.(
- ' g. A
. 3:525:44

  
 
  

’ ’14:; .
k 1) we! '
u AV“ (¥ 0
:1 “stems.
A AVAV“ v1- . «v
'e"":‘ 4
$15‘
v‘-
v

   

 
   
 

1
AVA“
9' w
41‘ ‘

V).
p 0 ? '5-1515.’
‘ 4. 0,0» '0'51
‘5». 45:0. '3-
A=E<> ’4» 00 "4 ‘ '
v4»<.< p 4%, 4’4,
9% «.1 "0 ‘74? ‘
<35“ fizga.‘ 4'4
25:15 11:: ’<
«ti-*1 ‘ 1
> v; ’ . »
1.564%“-

". . :‘
(3:4 ,‘€:¢:T“"".:;
A A'A'Ar. .'
I VA 41‘ L“ .gr. .. ,
I 0..
11 v “'M ‘1‘“ ~ ‘ 1v
. 1 y .
3:11:31 x." 1"” ' V I ‘9“ 1911351.. “‘3;
"""“" Q 9:41:15: .
v ‘ v.‘ 4 ":p‘)‘ “

1". 5‘ C v Q
ass-17:35:55.
gVAV“ AVAV,VA V‘.,v2nn

shun" — — 7 7 '

u

33% ‘
11I A» «111‘
. ”my, 1 M
4:33:53 4'

11.8118:
::f~ AI.

I
A

1
1%
V

 

'
‘ V V
"VL v Av L A
3’ ‘69e5':1:.:p".
v"
’ 1
' £5415
‘3 .45
5.,
5%:
>15
5":
.

Figure 30. Deformed blank for v=0.4 (a) Vf=0.218 (b) V1=0.384 (c) V1=0.6.

§

55

From initial observation, 15 wrinkles were associated with each blank,
which is double an actual part. Unlike actual results, the wrinkles appeared
evenly spaced around the part. The only Poisson’s ratio that showed a slight
difference was 0.3 with a volume fraction of 0.384. In Figure 29 (b) the deformed
part has 15 wrinkles, but less material is found on one side of the formed
hemisphere than the other, similar to the behavior when forming. Therefore, it
was determined that a Poisson’s ratio of 0.3 was a key value, which is
corroborated by the tensile testing results that show Poisson’s ratios of 0.35 and
0.4. The volume fraction chosen was 0.384. In addition, the plasticity interface
model was used in a simulation in order to see if there were any changes. The
result is shown in Figure 31. Similar to the linear softening interface, this model
was unable to capture the behavior found during experimentation. Therefore, it
was determined that the two preferred fiber orientation was insufficient, and work

began on the three and four preferred fiber orientation VUMAT’s.

 
    

 
  
    

“11,-; "mil“: “

,‘I‘M’Av “W: “1 Y
“igwwiv‘vfiz$€$\w .
«itmtgA Aw...

4A3”$\§‘ ..

‘ 'Ngfix§“ " .1-
fl§flb\“e '
‘ v;

      

 
        
       

  

‘ A‘A”

v2
1
VI"%

 
 

V:

. A ”a
y u v
,:.u"
u"

 
 

u
{‘0

 

r

   

     
  

      

Figure 31. Deformed blank for v=0.3 and Vf=0.384 with a plasticity model
interface.

56

4.2.2 — Three Preferred Fiber Orientations

The program for three preferred fiber orientations used the same
constitutive model as before, but with an extension for the third fiber direction.
Since the sensitivity analysis was already done, the Poisson’s ratio was set to
0.35, 0.4, and 0.375 for the 40°, 170°, and 250° fiber directions, respectively. The
initial volume fraction was set to .384. Difficulty arose when using the 250°
preferred fiber orientation, therefore 70° was used instead since the results for the
squeeze flow show symmetry from 0° — 180°and 180° — 360°. The resulting

deformation is shown below.

   
 
  

   
  
   

  

AW) 1
V ‘

V
74"4'

   
  
    
 
 
 
 
  
  

 
  
   
     

   
  
   
  
   
 
 
 

 

7’
a
‘-
4'4
:51? f

a
.17

'. u
an! V‘v‘n‘
A74 navy")

AV
E

B

    
 

      
  
    
  

 

      
  

n

  
   

  

L

is
'47
v‘

    
  
  
  

  

     

‘. «if

wfi,b
fifl:y ¢
4m

47"
7474

A
v
5

i!
7'
74

 

   
    
   

  

‘4
II

      

31..
fit

      

 

E
I
9'
e

A
‘V

W
A
in
V‘
A”
V"
' A

M

i

n
\

V
\

 

S, Mises

SNEG, [fraction = —l.0)
lAve. Crit.: 75%)
+8.8142+08
+8.0959+08
+7.377e+08
+6.658e+08
+5.939e+08
+5.220e+08
+4.502e+08
+3.783e+08
+3.064e+08
+2.345e+08
+1.627e+08
+9.078e+07
+1.890e+07

 

 

 

 

 

Figure 33. Stress contour plot for 3 preferred fiber orientations VUMAT.

57

[x109]

- 1mm 1 I l I I I l I r

 

p
8
l
l

4.00 L- -—

2.00— -—

ALLIE In ELSET BLANK

 

 

o.w l l 1 l L L l l
0.00 020 0.40 0.00 can mo
Time

(a)

 

 

2.40—

200—

1.” —-

1.20—

 

 

—— —
o.”

0.40—- —l

ALLKE In ELSET BLANK

 

 

 

o.w l l I l l l
0.00 0.20 0.40 0.” 0.80 1.00

Time

(b)

Figure 34. Internal (a) and Kinetic (b) energy plots for 3 preferred fiber orientations
VUMAT.

An important observation of the deformed part in Figure 32 is that there are

13 wrinkles. This is fewer than with the 2 preferred fiber orientations, but still

58

almost double the amount found in the actual formed part. The internal and kinetic
energy both show good quasi-static behavior.

4.2.3 — Four Preferred Fiber Orientations
The program for four preferred fiber orientations again used the same

constitutive model as before, but with an extension for the third and fourth fiber
directions. The Poisson’s ratio was set to 0.375 for the new direction and the initial

volume fraction to .384. The resulting deformation is shown below.

    
 

    
    

“a
u. “an“.
a» . 4
' =53

  

A" ‘ 4‘,
Emma

 
 
     

 

S. Mises
SNEG. {fraction = -1.0)
(Ave. Crit.: 75%]

+1.220e+09
+1.120e+09
+1.020e+09
+9.203e+08
+8.203e+08
+7.204e+08
+6.204e+08
+5.204e+08
+4.204e+08
+ .20

 

+2.054e+07

 

 

 

 

Figure 36. Stress contour plot for 4 preferred fiber orientations VUMAT.

59

 

12.00 — l l ' l ' r l ' a
10.00 I —
8.00 — _
0.00 — —

4.00 — —

ALLIE in ELSET BLANK

2.00 — -

 

 

0.00 L I 1 I 1 I 4
0.“) 0.20 0.40 0.60 0.80 1 .00

Time

 

(a)

 

4.00 '
3.50
3.00
2.50
2.00
1.50

1.00

ALLKE in ELSET BLANK

0.50

I I I I I I l I I I I I I I I I

 

 

 

0.00 .
0.00 0.20 0.40 0.00 0.80 1.00
Time

(b)

 

Figure 37. Internal (a) and Kinetic (b) energy plots for 4 preferred fiber orientations
VUMAT.

60

This simulation again showed 15 wrinkles as with the two preferred fiber
orientations. The internal and kinetic energy are also acceptable. Since
previously for the two preferred fiber orientations a sensitivity analysis was done
for the volume fraction and Poisson’s ratio, one was also conducted on the three
and four preferred fiber orientation models. This time the analyzed parameter was
the Young’s Modulus. The following table shows the results of this analysis.

Table 6. Sensitivity analysis of Young’s Modulus.

 

 

 

 

 

 

 

 

 

 

 

 

Young’s Modulus of Fiber % Complete # of Wrinkles
Directions — 3PFO
E=E*107 85 12+
[5:9108 85 12+
E=E*109 100 15
E=E*10‘” 100 23+
E=E*10" 100 cc
Young’s Modulus of Fiber % Complete # of Wrinkles
Directions - 4PFO
E=E*107 100 13
E=E*103 100 13-14
E=E*109 100 15
E=E*1010 90 27+
E=E*"1011 100 00

 

 

 

 

 

The results show that as the Young’s Modulus is decreased, the number of
wrinkles decrease, and vice versa. Yet, they do not decrease enough to match the
actual formed part. Therefore, it was determined that for a Young’s modulus
between an order of the 7th to 9th power, there is not a lot of change in the

simulation. For higher powers, the simulation becomes unstable.

61

4.3 Discussion of Results

The results of the simulations were somewhat different than the

experimental results. To begin this analysis a picture of different view of the

simulations versus the actual experiments are shown below.

 
      

111‘," WM”? 1 '
WVflVfiVQ ,
yLVAV5§L At t
, zigzag“
@1339"
Wm?
. “"1“!“ I

 
 
 
 

   

   

 
 
  
  

        
 

ll
Experimentally Two Preferred Fiber
Formed Part Orientations

0'1 ‘1
“'47 4
‘W'Vfififl

\V.v.”':,

Mr
4; Mai
.4471qu

 

Three Preferred Four Preferred
Fiber Orientations Fiber Orientations

Figure 38. Front angled views of formed part and simulations.

The results do not show a good correspondence between the wrinkling

behavior of the experimentally formed parts and the simulations. Where the

experimental part only has 7 wrinkles, the simulated parts have 15 wrinkles,

nearly double. The two preferred fiber orientation model used with the linear

62

softening showed the most unusual behavior. It was highly asymmetric, and had
one half of the formed part with more material on that side than the other, which
is usually how the actual part formed. Yet, the 3 and 4 preferred fiber
orientations didn’t show this behavior. In addition, the wrinkles tended to be

flatter and wider than the 2 preferred fiber orientations model. This is more of a

resemblance to the actual part.

     

eunnn
oran-
, r4-

 
 
   
  
 
 
  

A
«v
v, V

       
   
  

v
\

  
   
   
  
 

.t . ~
"3'“
.. “m
“53"",
qua “n.
‘

.r.‘ f“.
107,1'2‘

V

     

 

 

Experimentally Formed Two Preferred Fiber
Part Orientations

    
 
   
 
 
 
 

 

0

  
   

'1
"7‘57”.
000

 

    
   

0’ ,I. ~
I\ 4
- .v;
y ‘ : in. .
4' .
‘vl'l'l'A‘ ' “1‘“
"In“?
ram“
‘

kg Aflf‘ '

‘4 ‘
V I; ‘ V
'51:. ‘ A :: .
2‘ i ‘ V’ > 1
‘= ”$10. M'AWAVAVIVIM .
‘ l

 

       
   
   
  

 

      
  

     
       

~ -- . r V, i
«new: , i.e.tumwmvm.,
‘ v, ‘ v;V.p-;V: '0': 'I "AMN$V‘EF7 "V"!'I:;VA
4A» 59"“‘5'55 .e'» "0 "“‘Rfi’ V)" 90""
‘ f: .‘ r w '1'
sv'ffihhif‘wjfifégg‘ =.£=;=1=‘.:i'~;1 E‘s?! ’1’“
.' "'H :1”, if: 3: v . l

.1
5' w ‘
..
..

n
'u
,y

g1 x

1“
. <9

           

v
.. I
u

 

Three Preferred Fiber Four Preferred Fiber
Orientations Orientations

Figure 39. Top views of experimental and simulation results.

63

The top views further show the asymmetry of the parts. Yet, it is not an
exact replica of what was found in experiment. One reason for this could have
been the fact that the updated material law model uses a constitutive relationship
with linear elastic stress-strain relations. By definition, an elastic material is one
that returns to its original (unloaded) shape once the applied forces are applied.
This type of elastic behavior is before permanent defamation. In this case a more
suitable model may be a viscoelastic material. In a viscoelastic material the state
of stress is a function of the strain and the time rates of change of the strain as
well. Therefore it uses a combination of elastic (spring-like) and viscous (dashpot-
like) elements to form a viscous-elastic model. Many materials such as glass,
ceramics, plastics, synthetic rubbers, and even biomaterials are considered to be
linear viscoelastic materials [29]. Therefore, this may be a more suitable model to
use.

In 2000, Klasztomy et al derived a viscoelastic model for unidirectional
fibrous polymeric composites. The composite was modeled as a viscoelastic
isotropic polymer matrix and elastic monotropic fibers. The Mittag-Leftler fractional
exponential functions were used in order to model the shear/bulk creep in the
matrix. Therefore, the viscoelastic model for the matrix was described with 2
elastic and 6 viscoelastic constants and the elastic fibers were described with 5
elastic constants. Collectively both models were used to derive the coupled
constitutive equations and the composite together was modeled as a
homogeneous monotropic material with 5 elastic constants, and 27 viscoelastic

constants (i.e. 9 long lasting compliance ratios, 9 retardation times, and 9 fractions

64

defining the order of the fractional exponential functions). This study also provides
a way to theoretically predict the viscoelastic constants [30].

In addition, Holzapfel et al used a viscoelastic model for a fiber-reinforced
composite material that exhibits direction-dependent properties and sustain finite
strains without a significant volume change. The composites used were made of a
soft matrix material reinforced by two different types of fibers, or two fiber
directions. Like the model used in this study, a global response is given by a
summation of all individual responses of the material. This work developed a
closed-form expression for the fourth-order elasticity tensor. Constitutive models
were presented for special cases including orthotropic, transversely isotropic, and
isotropic hyperelastic materials at finite strains with and without dissipation. The
model was validated using ZD and 3D simulation results and comparing those to
experimental results of a pressurized laminated circular tube, which has a strong
anisotropic response [31]. This is considered to be on a macro mechanical or
continuum approach.

Additional works on viscoelastic constitutive models on a microscopic level
were done by R. M. Haj-Ali and A. Muliana in 2004. The approach was to idealize
each unidirectional lamina using the Aboudi four-cell micro model. This is coupled
with incremental formulation in interfaces of the average stress and strain in each
sub cell. The matrix sub cells are described by the Scharpery non-linear
viscoelastic model and the fiber is considered transversely isotropic and linear
elastic. The sub cell constitutive relationships are embedded in a numerical

stress-update algorithm. This framework can easily include temperature, moisture,

65

and physical aging effects. It was implemented within a shell-based non-linear
finite element analysis by assuming plane stress. The formulation was validated
using several experimental off-axis specimen creep tests and was applied to a
laminated panel and a composite ring, both showing agreement [33].

The time dependent response of polymeric composite systems was studied
using classical homogenization methods in 2002 by Sejnoha and Zeman. The
study focused on random, non-periodic material systems with loading that
promotes the viscoelastic deformation of the material. Two modeling approaches
were used. In the first, the material was represented using volume elements with a
small number of particles to statistically represent the microstructure of the
composite. These elements are periodically dispersed and a finite element
analysis can be carried out. The second approach is based on the Hashin-
Shtrikman variational principles. The randomness of the fiber is incorporated using
statistical descriptors. This work was only applicable to microstructures that can
be described by the two-point probability function [34].

Overall, it is possible that any of these methods could be used in place of
the current, elastic representation. In addition, the processing temperature of the
composite could have an effect on its behavior. It has been shown in a recent
study of ramie fibers, another natural fiber, that processing temperatures between
180 — 200 °C for a given period of time can cause degradation of the mechanical
properties of exposed fibers [6]. Since the melting temperature of the
polypropylene is approximately 160°C and the composites are fabricated and

formed at 190°C it will melt during forming and possibly expose fibers to the

66

temperature degradation. This can be remedied by surface modification of the
natural fibers, which also improves the mechanical properties [6].

Another important issue to consider is the moisture absorption of the fiber.
Natural fibers, such as Kenaf, are hydrophilic which means they absorb water. An
average range of moisture absorption is 5 — 10% although in some cases it has
been reported to be as high as 20% (~30% was found in this study). This moisture
can affect the final properties of the composite once it is processed by creating
voids between the fiber and the matrix which can interfere with the fiber/matrix
adhesion [13].

In addition, the preferred fiber orientations can be taken into account. The
method used selected at most four preferred fiber orientations, and measurements
were only taken at every 10 degrees. A more accurate representation may be
found by taking measurements every 2 - 5 degrees and by using an Orientation
Distribution Function (ODF). This would allow for a multi-orientation approach and
could prove to capture the behavior of this material more accurately.

Furthermore, the Young’s modulus of the fiber is on the same order as the
matrix. The updated material law model was originally built based on representing
the matrix as an elastic material, which was adequate because of the higher
strength of the glass fibers. Therefore, the fibers controlled the behavior and
material properties of the composite. In the Kenaf/polypropylene composites used
in this study, the Young’s modulus of the fiber is the same order of magnitude as

the matrix. Therefore, a more accurate model for the matrix material, which is not

67

elastic, should be used. This is another reason a different material model such as
viscoelastic or Viscoplastic should show better results.

In conclusion, any of these affects could have led to poor representation
during simulation. A viscoelastic model could prove to be a better representation
of the material behavior. Additionally, fiber treatment could be done that would not
only increase the overall strength and properties of the fibers, but will also keep
them from absorbing moisture and experiencing degradation during processing. In
addition, the use of an ODF may also more accurately capture the material

behavior.

68

Chapter 5

CONCLUSIONS

In the present work, an optimized process for the fabrication of Kenaf
natural fiber and polypropylene composites has been presented. This process
proved to provide good adhesion between the fiber and matrix. In addition, an
evenly dispersed fiber distribution was achieved. The composites also showed
comparable tensile and flexural strength with other natural fiber and polypropylene
composites fabricated by compression moulding and other fabrication processes.

In addition, simulations were performed using ABAQUS/Explicit with a
VUMAT (user defined material constitutive relationship) coupled with a VUINTEFI
(user defined interface model) in an attempt to represent experimental forming of
the composites. The numerical simulations did not show good resemblance.
Possible reasons for the discrepancies include the use of an elastic constitutive
relationship and the addition of a coupling agent into the polypropylene powder
instead of fiber treatment. The stress-strain curves from tensile testing clearly
show non-linear behavior. Therefore, a viscoelastic constitutive model could prove
to be more appropriate. Additionally, fiber treatment would not only help keep the
fibers from being damaged during heat processing, it would also decrease the
moisture absorption of the fibers, further improving the consistency and
performance of the composites. Furthermore, the use of an ODF would allow for a

multi-orientation representation of the material.

69

Chapter 6

FUTURE WORK

There are many options for future work with the kenaf and polypropylene
composites presented in this study. These include but are not limited to the
varying of the fiber volume fraction and using fiber treatment instead of adding the
coupling agent to the polypropylene powder. In addition, more testing could be
completed including tensile testing at elevated temperatures and Dynamic
Mechanical Analysis (DMA) testing for further characterization.

Numerical work could include the use of different material model — possibly
viscoelastic, or one that is found to more closely resemble the actual behavior of
the composite material. In addition, different interface models can be incorporated.
The linear softening and plasticity models could be used with a different material
model or any of the models presented in literature could be substituted. Also, the
use of an ODF could improve the simulations representation of the actual material
behavior.

Additional work could also use the sandwiching of the kenaf and

polypropylene sheets with a metal such as aluminum.

7O

BIBLIOGRAPHY

COMPOSITE MATERIAL CHARACTERIZATION

I.

Mohanty, A.K.., Drzal, L.T., and Misra, M. “Novel hybrid coupling agent as
an adhesion promoter in natural fiber reinforced powder polypropylene
composites.” Journal of mrials Science Letters 21 (2002): 1885-1888.

Wambua, P., lvens, J., and Verpoest, l. “Natural fibres: can they replace
glass in fiber reinforced plastics?” Composites Science and Technology
63(2003): 1259 — 1264.

Feng, D., Caulfield, DE, and Sanadi, A. R. “Effect of Compatibilizer on
the Structure-Property Relationships of Kenaf-Fiber/Polypropylene
Composites.” Polymer Composites 22(2001): 506 — 517.

Joshi, S. V., Drzal, L. T., Mohanty, A. K., and Arora, S. ”Are natural fiber
composites environmentally superior to glass fiber reinforced
composites?” Composites: Part A 35 (2004): 371 — 376.

. Mohammad, P., and Sain, M. M. ’Sheet-Molded Polyolefin Natural Fiber

Composites for Automotive Applications.” Macromolecular Materials and
Engineering 288(2003): 553 — 557

Mohanty, A.K.., Misra, M., and Hunrichsen, G. “Biofibres, biodegradable
polymers and biocomposites: An overview.” Macromolecular Materials and
Engineering 276/277 (2000): 1-24.

Mohanty, A.K.., Drzal, L.T., and Misra, M. “Engineered natural fiber
reinforced polypropylene composites: influence of surface modifications and
novel power impregnation processing.” J. lthesion Sci. Technol., 16
(2002): 999 - 1015.

Biagiotti, J., Fiori, 8., Torre, L., Lopez-Manchado, M. A., and Kenny, J. M.
“Mechanical Properties of Polypropylene Matrix Composites Reinforced
With Natural Fibers: A Statistical Approach.” Polymer Composites,
25(2004): 26 — 36.

Nishino, T., Hirao, K., Kotera, M., Nakamae, K., Inagaki, H. “Kenaf
reinforced biodegradable composites.” Composites Science and
Technology 63 (2003): 1281 - 1286.

10. Karian, Harutun G. Handbook of Polvproovlene and Polypropylene

Composites. New York: Marcel Dekker, Inc., 1999.

71

11.Chung, T. C. Mike. Functionalization of Polyolefins. New York: Academic
Press, 2002.

12.Mohanty A.K., Misra M., and Drzal L.T... “Sustainable Bio-Composites
from Renewable Resources: Opportunities and Challenges in the Green
Materials World.” Journal of Polvmers and the Environment 10(2002): 19
- 26.

 

13. Rouison, D., Sain, M., and Couturier. M. “Resin transfer molding of natural
fiber reinforced composites: cure simulation.” Composites Science fl!
Technology 64(2004): 629 - 644.

14.Mohanty, A.K., Wibowo, A., Misra, M., and Drzal, L.T. “Effect of process
engineering on the performance of natural fiber reinforced cellulose
acetate biocomposites.” Composites Pa¥rt A: Applied Science and
Manufacturing 35(2004): 363 - 370.

15.Webber Ill, Charles L., Bledsoe, V.K., and Bledsoe, R. E. “Kenaf
Harvesting and Processing.” Trenps in New Crops and New Uses (2002):
340-347.

DELAMINATION

16.Zou, 2., Reid, SR, and Li, S. “A continuum damage model for
delaminations in laminated composites.” Journal of m Mechanics and
Physics of Solid_s 51 (2003): 333-356.

17.Alfano, G., and Crisfield, M.A. “Finite element interface models for the
delamination analysis of laminated composites: mechanical and
computational issues.” International Journal for Ntflertcal Methods in
Engineering 50 (2001 ): 1701 -1 736.

18. Hillerborg, A., Modeer, M., and Petersson, P.E. “Analysis of crack formation
and crack growth in concrete by means of fracture mechanics and finite
elements.” Cement and Concrete Researih 6(1976): 773 - 782.

19.Dugdale, D.S. ”Yielding of steel sheets containing slits.” Journal of the
Mecha_nics and Physics of Solids 8 (1960): 100 - 104.

20. Borg, R., Nilsson, L., and Simonsson, K. “Simulation of delamination in fiber

composites with a discrete cohesive failure model.” Compositesjciencp
and Technology 61 (2001): 667-677.

72

21.Wells, ON, and Sluys, L.J. “A new method for modeling cohesive cracks
using finite elements.” lntemational Journal for Ninerical Methods in
Engineering 50(2001): 2667 - 2682.

22.Espinosa, H.D., Dwivedi, S., and Lu, H.C. “Modeling impact induced
delamination of woven fiber reinforced composites with contact/cohesive

laws.” Computer Methods in Applied Mechanics and Engineering
183(2000): 259 - 290.

 

23.Blackman, B. R. K., Hadavinia, H., Kinloch, A.J., and Williams, J.G. “The
use of a cohesive zone model to study the fracture of fibre composites and
adhesively-bonded joints.” International Journal of Fracture 119(2003): 25 -
46.

24. Borg, R., Nilsson, L., and Simonsson, K. ”Simulation of low velocity impact
on fiber laminates using a cohesive zone based delamination model.”
Composites Science and Technology 64(2004): 279 - 288.

25. Bruno, D., Greco, F., and Lonetti, P. “A coupled interface-multilayer
approach for mixed mode delamination and contact analysis in laminated
composites.” International Journal of Solids and Structures 40(2003): 7245 -
7268.

26.Corigliano, A. and Ricci, M. “Rate-dependent interface models: formulation
and numerical applications.” lntemational Journal of Solids and Structures
38(2001): 547 — 576.

27.2ampaloni, MA, A multi-preferred fiber orientation constit_u_t_ive model for
mer mat reinforcegthermoplastics with a random orientation applied to the
stamp thermo-hygroforming process. Michigan State University, 2003.

28. Hibbit, Karlsson & Sorensen Inc. ABAQUS Manuals. 2004.

29.Ugural, Ansel C. and Fenster, Saul K. Advanced Strength and Applied
Elasticity Thirg_§dition. Upper saddle River, NJ: Prentice Hall, 1995.

30.Klasztorny, M. and Wilczynski, A. P. ”Constitutive Equations of
Viscoelasticity and Estimation of Viscoelastic Parameters of Unidirectional
Fibrous Polymeric Composites.” Journal of Composite Materia_|3_34(2000):
1624 — 1639.

31.Holzapfel, G. A., and Gasser, T. C. “A viscoelastic model for fiber-
reinforced composites at finite strains: Continuum basis, computation
aspects and applications.” Computer Methods in Applied Mecha_nicsa_nd
Engineering 190(2000): 4379 — 4403.

73

AAAAAAAAA

lIIIIIIIIlllllIllll llillllliliii

2732