Title Page
 Table of Contents
 Classical models for powder...
 A time dependent model for alumina...
 Realization of cristescu model...
 Biographical sketch

Title: Numerical modeling of compaction of particulate systems
Full Citation
Permanent Link: http://ufdc.ufl.edu/UF00100673/00001
 Material Information
Title: Numerical modeling of compaction of particulate systems
Physical Description: Book
Language: English
Creator: Wang, Weijian, 1971-
Publisher: State University System of Florida
Place of Publication: <Florida>
Publication Date: 1999
Copyright Date: 1999
Subject: Aerospace Engineering, Mechanics, and Engineering Science thesis, M.S   ( lcsh )
Dissertations, Academic -- Aerospace Engineering, Mechanics, and Engineering Science -- UF   ( lcsh )
Genre: government publication (state, provincial, terriorial, dependent)   ( marcgt )
bibliography   ( marcgt )
theses   ( marcgt )
non-fiction   ( marcgt )
Summary: ABSTRACT: The goal of this dissertation is to combine some constitutive models for particulate systems with finite element methods (FEM) in order to understand the mechanical behavior of cohesive powders during the compaction process and make some proper prediction for manufacture control and tools design. One classical constitutive model -- the Drucker-Prager model -- is implemented with advanced FEM software ABAQUS. Also, a more complicated model -- the Cristescu model -- will be introduced. This model will hopefully be used with ABAQUS. Some algorithms will be introduced in this thesis.
Summary: KEYWORDS: powder
Thesis: Thesis (M.S.)--University of Florida, 1999.
Bibliography: Includes bibliographical references (p. 62-63).
System Details: System requirements: World Wide Web browser and PDF reader.
System Details: Mode of access: World Wide Web.
Statement of Responsibility: by Weijian Wang.
General Note: Title from first page of PDF file.
General Note: Document formatted into pages; contains v, 64 p.; also contains graphics.
General Note: Vita.
 Record Information
Bibliographic ID: UF00100673
Volume ID: VID00001
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: oclc - 45265664
alephbibnum - 002484126
notis - AMJ9740


This item has the following downloads:

wang ( PDF )

Table of Contents
    Title Page
        Page i
        Page ii
    Table of Contents
        Page iii
        Page iv
        Page v
        Page 1
        Page 2
        Page 3
    Classical models for powder materials
        Page 4
        Page 5
        Page 6
        Page 7
        Page 8
        Page 9
        Page 10
        Page 11
        Page 12
        Page 13
        Page 14
        Page 15
        Page 16
        Page 17
        Page 18
        Page 19
        Page 20
        Page 21
        Page 22
        Page 23
        Page 24
        Page 25
        Page 26
        Page 27
        Page 28
    A time dependent model for alumina powder - cristescu model
        Page 29
        Page 30
        Page 31
        Page 32
        Page 33
        Page 34
        Page 35
        Page 36
        Page 37
        Page 38
        Page 39
        Page 40
        Page 41
        Page 42
        Page 43
        Page 44
        Page 45
        Page 46
        Page 47
        Page 48
        Page 49
        Page 50
        Page 51
        Page 52
        Page 53
        Page 54
    Realization of cristescu model with FEM
        Page 55
        Page 56
        Page 57
        Page 58
        Page 59
        Page 60
        Page 61
        Page 62
        Page 63
    Biographical sketch
        Page 64
Full Text







I would like to express my deepest gratitude to Dr. Nicolae D. Cristescu, my

advisor, for his guidance, tremendous help, and encouragement for the whole two years

passed. Dr. Bhavani V. Sankar, Dr. Peter G. Ifju and Dr. Edward K. Walsh are so kind to

serve on the committee. I am so grateful to their support and help. I also appreciate Dr.

Oana Cazacu for her guidance and discussions with me last year. Many other thanks

should go to my friends. I would have given up without their encouragement.

Finally, the financial support of the Engineering Research Center (ERC) at

University of Florida is gratefully acknowledged.


ACKNOWLEDGMENTS ...................................................... ii

ABSTRACT ...................................................................................... v


1 IN T R O D U C T IO N ................................................................................ 1

2 CLASSICAL MODELS FOR POWDER MATERIALS ............................ 4

2 .1 In tro d u ctio n ........................................................................... 4
2.2 Critical State Model .......... ..................................... 5
2.3 Drucker-Prager/Cap M odel .................................. ............. 7
2 .4 F riction ................................................. ........... .. . ........ 10
2.5 Compaction Cases Study ................................ ......... 11
2.5.1 B background ......................... ......................... ......... 11
2.5.2 Computational Model ................................. ......... 12
2.5.3 FEM Model and Mesh .................... 14
2.5.4 Material Parameters ............................................. 14
2.6 R esult and A analysis .................................. .............. .............. 18
2 .7 C onclu sion .................................................. ........... ......... 27

M ODEL ........................................... ......... 29

3.1 Introduction .......... ...... .......................... ........... 29
3.2 Description of Cristescu Model .................................... .......... 30
3.2.1 General Form .......... .......................... ............. 30
3.2.2 M material P aram eters ...................................................... 31
3.2.3 Determination of the Yield Function H(c) ....................... 33
3.2.4 Determination of the Viscoplastic Strain Rate
Orientation Tensor N(c) ......................... 34
3.3 Simulation Algorithm for General Loading Test ........................... 34
3.3.1 Description ..................................... ............. ......... 34
3.3.2 Result and A analysis ........................................ ......... 37
3.4 Transient Creep Computation ................................................ 39
3.4.1 Description ..................................... ............. ......... 39
3.4.2 Result and Analysis ...................................................... 41
3.5 Analytical Solution for Transient Creep ........................................ 51

3.5.1 D description .......................................................... 51
3.5.2 Result and Analysis ..................................................... 52
3 .6 C o n clu sio n ........................... ........ ..... ............. ............. 54
4 REALIZATION OF CRISTESCU MODEL WITH FEM ............................ 55
4 .1 G e n e ra l ........... ... ..... ............................................. .......... 5 5
4.2 Algorithm for Stress Control Test ............................................ 55
4.3 Algorithm for Strain Control Test ....................................... 57
4.4 Algorithm for ABAQUS ........................ .......... 58

REFERENCES ........................................ ........... ......................... 62

BIO GRA PH ICAL SK ETCH ........................................................................ 64

Abstract of Thesis Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Master of Science



Weijian Wang

August 1999

Chairman: Graduate Research Professor Nicolae D. Cristescu
Major Department: Aerospace Engineering, Mechanics, and Engineering Science

The goal of this dissertation is to combine some constitutive models for particulate

systems with finite element methods(FEM) in order to understand the mechanical behavior

of cohesive powders during the compaction process and make some proper prediction for

manufacture control and tools design. One classical constitutive model the Drucker-

Prager model is implemented with advanced FEM software ABAQUS. Also, a more

complicated model the Cristescu model will be introduced. This model will

hopefully be used with ABAQUS. Some algorithms will be introduced in this thesis.


Powder compaction is important in many aspects of industry, including the flow

and storage of agricultural products, the forming of ceramic components, powder

metallurgy parts, and pharmaceutical tablets, road and dam construction, subway and

tunnel excavation, mine and oil-well locating, and rock slides and avalanches. For

example, it is estimated that more than 80% of all medication doses are administered as

tablets, that is, in unit dosage forms prepared by compacting powders in dies.

The object of powder pressing is to quickly and reproducibly form a defect-free,

homogeneously dense powder to net shape. The pressing pressure should be high enough

to impart sufficient strength for subsequent handling but low enough to avoid excessive

wear on the press and tooling. Pharmaceutical tablets must have sufficient strength to

survive packaging and transportation but must be weak enough to disintegrate after

administration. Ceramic compacts may require higher strength to withstand the high

stresses during green compact machining.

In recent decades, a substantial amount of research and development has been

conducted in different areas and significant technological advances have been made.

Scientists have been trying to find appropriate models that can describe the characteristics

of different kinds of powders and predict the behaviors of the powder precisely during the

manufacturing process for the industry. It is of great significance because this will bring

more economic effects to the industry and engineers can control the particle compaction.

In modeling of powder compaction, constitutive relation and numerical methods

are developed for simulating the response of powder during compaction. A good model

should be able to predict the state of stress, strain and the displacements everywhere in a

powder sample subjected to some external and internal forces. Traditional soil mechanics

and rock mechanics concepts and theories are introduced and adopted into the constitutive

models for the powders. At the same time, finite element methods (FEM) have become

the most important numerical methods because they can accurately predict the density

variation in a compaction procedure and can provide information of defect detection.

More and more scientists are engaged in the work of combining classical or promising

constitutive models with FEM and compare with the results from the experiments and

improve them. This shows a great potential in the research work of powder materials.

Many theoretical models have been issued for understanding the constitutive

relations of powders. Among them the Cam-Clay model21 and the Drucker-Prager

model13 are undoubtedly the most popular and most widely used. They will be introduced

briefly in chapter two and some computation examples will be introduced at the same

time. These examples combine the Drucker-Prager model with some advanced FEM

software (ABAQUS)14. They not only show interesting results but also tell us how

significant it is to implement the constitutive model into the modem numerical methods.

The advantages of the Cam-Clay model and the Drucker-Prager model are their easy-to-

use property and comparatively accurate agreement with experimental results under

various conditions. Unavoidably, the disadvantages limit the use of these models when

more experimental conditions are to be considered. For example, both of these two

models do not consider the time effects while it is very important in the manufacturing

process. This requires some more sophisticated models. Cristescu model is one of these.

Dr. Nicolae D. Cristescu has been working in the area of geomaterials and powder

materials for dozens of years.4,6-10,16-18 His models for rock materials have been widely

used in scientific research and industry. One of his models for ceramic powder is

introduced in chapter three. Efforts will be made to implement this model with ABAQUS

software since ABAQUS provides some useful interfaces for users. Although no

satisfactory result can be presented in this dissertation, I would still like to introduce and

discuss the methods and algorithms briefly in chapter four.


2.1 Introduction

The powder behaves quite differently than the same material in solid phase and in a

more complex way, but under essentially monotonic loading conditions rather simple

constitutive models provide useful design information. The body composed of a large

amount of particles is thought of as a continuous media as are other solid materials, so

that the constitutive model of similar form to that used for other solid materials is also

used for powder materials.

The strain rate is decomposed into two parts, sE is the elastic strain rate and sF is

the irreversible (plastic) strain rate.

S= sE sI (2.1)

The elastic part is defined by a differential Hooke's Law, which is

F 1 1 (
sE- [- +( ) l (2.2)
_J 2G 3K 2G 1(

in which K and G are elastic moduli.

For metals, sometimes the elastic part is enough for many applications. Some solid

materials accept a linear stress-strain relationship under a considerably large load. While

for powders the situation is quite different. The yield stress of the powder is often much

lower than that of metals. The post-yield behavior cannot be described by the elastic

theory, and the constitutive equation for powder is of great concern for the industry. One

comparison between the implementation of elastic theory and plastic theory will be shown

in this chapter.

There are many kinds of plastic models that can be used to describe the plastic

behavior of the powder. We cannot use the perfect plasticity models because the powder

exhibits work hardening. Fortunately ABAQUS provides several options for non-metal

materials. The plastic part is defined by the flow rule:

k g (2.3)

where g is a "plastic potential" and X is the flow parameter. (Mises, 1928)

2.2 Critical state model

The critical state model (Cam-clay model) was developed by Roscoe and Schofield

et al., 1968 and Parry, 1972, at Cambridge University. This model is used often for

modeling non-cohesive materials. It contains mainly two parts: a yield surface, and a

critical state surface that is the locus of effective stress.

The equation of the yield surface is

f(p,q) (P- 1)2 +( q 2 -1= 0 (2.4)
p a Ma

The equation of the critical state line is

q=Mp (2.5)

p and q are the mean stress and the equivalent stress respectively. They are defined

P =- (, + CY 2 +3)
S:S (2.6)
q= 2:S

S is the stress deviator tensor.

S = T +pI (2.7)

T is the stress tensor.

M and P are user-specified material parameters and a is the evolution parameter.

Figure 2-1 shows how the yield surface of the Cam-clay model looks like.

Critical State Line

" Yield surface

Figure 2-1. Typical yield surface for Cam-clay model

Cam-clay model is often used for soil and sand materials. It doesn't consider the

cohesion of the material. Unfortunately for most of the metallic powders and

geomaterials, cohesion cannot be neglected. Thus another model is introduced here, it is

the modified Drucker-Prager/cap model.

2.3 Drucker-Prager/Cap model

A typical yield surface of this model is shown in Figure 2-2. This surface is shown in

the q-p space in which p and q are the mean stress and the equivalent stress and S is the

stress deviator tensor. Their definitions are the same as that defined in (2.6) and (2.7).

Cap, F,

R(d+P tanp3)

Figure 2-2. Typical yield surface for Drucker-Prager/Cap model

Similar to the Cam-clay model, the yield surface in the modified Drucker-Prager/cap

model is expressed by two segments: a shear failure surface and a cap that intersects the

equivalent pressure axis. There is also a transition region between these segments, which

provides a smooth connection. The cap surface serves two main purposes. It bounds the

" --------------------------

yield surface in pure hydrostatic compression. It also controls the volume dilatancy when

the material yields in shear.

The equation of the shear failure surface is written as,

q -ptan3- d = 0 (2.8)

where P3 is the internal friction angle and d is the cohesion of the material.

The equation of the shear failure surface is written as,

/ ,2

(p-pa)2 + -R(d+pa tan3)= 0 (2.9)

l^ cosp]

where R and ca are material parameters that control the shape of the cap. Pa is an

evolution parameter and defined as

Pb Rd
Pa = b R (2.10)
1+Rtan p

The parameter ca is a small number, typically 0.01 ~ 0.05. It is used to provide a

smooth connection between the shear failure surface and the cap surface.

The cap surface has an elliptical shape with a constant eccentricity in the meridian q-

p plane.

The function involved in the transition segment is

(P-pa)2+ q-(1- )(d+patanp) -u(d+patanp) =0 (2.11)
1 cosp

Plastic flow is non-associated to the shear failure surface and the transition segment,

but is associated to the cap surface. As shown in Figure 2-3, the flow potential surface is

made up of two parts -

shear surface respectively.

The functions of Go

G, = (p -Pa)2 +

Gs = [(p-Pa)tanf


Go and Gs. They are associated with the cap surface and the

and Gs are written as,

Rq )2
cos P3 (2.12)

2 +( q )2


G, (Shear failure)

- \ ..- --- - \. I ap) (l4Ua-a Seci


Figure 2-3. Flow surface for Drucker-Prager/Cap model

While in the manufacturing of powder materials, the variation of the density is of

great concern and additional derivation about the density variation is given here.

We use the formula

v V (2.13)
Vo Vo
no" '

in which ev represents the volumetric strain, Vo and Vi are the initial and final

volumes respectively.


m m

8 = PI o = PO 1 (2.14)
m p1


p1 =Po (2.15)

Equation 2.16 will be used to compute the density variation during the compaction

described in the following examples.

2.4. Friction

Usually friction exists between the powder sample and the container wall. Thus

friction cannot be neglected in the computational model. The standard friction model is

the classical Coulomb friction, which is accepted and used in ABAQUS. [11][12] This model

states that at the contact surfaces the two bodies do not slide over each other as long as

the shear stress magnitude is less than P3, the friction coefficient times the pressure stress

between them. For alumina powder we assume P3 is a constant. As a matter of fact, we

do not consider the deformation of the container wall. So the container wall is always

modeled as a rigid surface in a FEM model by ABAQUS.

2.5 Compaction Cases study

2.5.1 Background

In the industry of manufacture or transport and handling, the shapes for products are

versatile and the containers for storing powders have various shapes. See Figure 2-4 and

Figure 2-5.

a. Solid cylinder

b. Cylindrical pipe

c. T-shape d. Cup-shape

Figure 2-4. Different shape powder products

r to I

a. Pyramid-shape container b. Cone-shape container
Figure 2-5. Different shape containers

2.5.2 Computational model

In this chapter, two cases are used as typical samples for computation. See Figure

2-6 and Figure 2-7. In case I, the powder is surrounded by the die with a cylindrical

shape. Correspondingly, a cylindrical shaped powder sample forms after the compaction.

There are two punches acting on the top and bottom of the powder sample. One can

assert a downward load on the top punch or an upward load on the bottom punch.

The second case is a little more complicated. Different from the last case, there are

two punches at the bottom of the sample. Thus two forces can be loaded from the

bottom. Finally a cup shaped sample forms.

In the uniaxial test, the lateral surfaces of the sample are regarded as fixed because

of the rigid die wall of the container. In the computation, we have 2 =3 = 0 (in

Cartesian coordinates) or sr = -e = 0 (for cylindrical coordinates).







Figure 2-6. Case I- cylindrical sample apparatus for compaction


7 Punch

Die Die



P3 P3

Figure 2-7. Case II Cup shape sample apparatus for compaction

2.5.3 FEM model and mesh

In both cases 8-node rectangular element is used. Due to the symmetry, only half of

the cross section is considered. Figure2-8, Figure 2-9 and Figure 2-10 describe the model

size, boundary conditions and load conditions.

In the first case, there are altogether 49 elements. The left side is the axis of the

sample. Because of the symmetry, only the horizontal displacement is prescribed. The

bottom is assumed to be fixed. The load acts on the top. Two kinds of loads are

considered in this case. One is a uniform displacement load, that is, the top surface will

keep parallel to the horizontal plane. The other one is the uniform distributed load.

2.5.4 Material Parameters

Table 2-1 lists the parameters used in the computation. The values of these

parameters are obtained from the experiments.

Bulk modulus, K 6777.35MPa

Shear modulus, G 3531.86MPa

Elastic modulus, E 9027.43MPa

Poisson's ratio, v 0.278

Internal friction angle, f3 16.50

Cohesion, d 5.5Mpa

Cap shape parameter, R 0.558

Transition surface parameter, ac 0.03

Table 2-1 Table of Material parameters




-"- 1424 1824 2224

1222 1622 2C22
1422 1822 2222

#1^, 1220 1620 2C20
14:0 1820 2220

89) r (B1@ (B (B@ B
, 1218 1618 2018
1418 18 8 2218

116 1216 1616 2 16
- 1416 18 6 2216

S114 1614 2(14
1414 1814 2:14

#^, 1212 1612 2 12
1412 18 2 2212

1010l 1410 1810 2210

0.01 m --

Figure 2-9. FEM mesh and boundary conditions for cylindrical shape sample

0.01 m








0.385 m



@ 0 @ @@@@
0 @ @@@@
@ @ @ @ @@(@
@ @ @ @ @@@@
1 3 5 7 9 11

13 --

15 --

19 --


- -H- 8013

- -H- 8015



-- 8021

Figure 2-10. FEM mesh for cup shape sample

801 803 805


1 18 48 1 1 88 1 1

i | 8019

I~p1' 8 8j.1~


21 8021

Figure 2-11. Boundary conditions for cup shape sample



The cap hardening curve shown in Figure 2-11 is also from experimental data. It

determines the corresponding parameter option in the ABAQUS input file.

0 05 0 1 0 15 02 0 25 03 0 35
plastic volumetric strain

Figure 2-11. Hardening curve

2.6 Result and analysis

2.6.1 Case 1.1

Cylindrical shape sample, distributated load acting on the top

This is a load control case. Note that the top surface is not flat any more after the

compaction. Due to the friction existing between the sample and the wall, the

displacement at the edge is smaller than at the center. Thus the initial horizontal lines form

into the curves with the concave downward shape, as shown in Figure 2-12.

Figure 2-13, Figure 2-14 and Figure 2-15 show the distribution of the axial stress,

the volumetric strain and the density after the compaction respectively. All of them

distribute stratiformly because of the existence of friction. If there is no friction, their

distribution should be uniform over the cross-section.


The volumetric strain and the density actually show the same result because they

differ only by a constant (see equation 2.11). Density is of great concern in the industry

because people want to know how the powder mass is distributed after the compaction,

and this distribution determines the quality and the characteristic of the products.

LFigure 2-12. Deformation after loading

Figure 2-12. Deformation after loading

T. 49L-01
- 97L-01
-6 45F-1>o
-.5 DIX-ot
-S 42k-01
-4, DOz~O1
-A 311,01
-3 671.Ot
-3 35LOt
- 2. BJE+o1
-3 1z.Oi
- 9 OIO1
I-31 O

Figure 2-13. Axial stress distribution after loading



-- 3. 3; i-01
-- 7. 1z-O
3. 111-01

-B3 2I-01

-2 Mlr-01
-2 r51-01

-J 96L-01i
-1. 93i-Oi
--] TO -01


Figure 2-14. Volumetric strain distribution after loading

-i. AE, QD
91. BE OD
. ,IE UD

2 2DE D0
S24E, 00D
,2 DE 3 OD

3E tE UD


Figure 2-15. Density distribution after loading

2.6.2 Case 1.2 Cylindrical shape sample, uniform displacement on the top

Figure 2-16, Figure 2-17, Figure 2-18 and Figure 2-19 show the results under the

uniform displacement load. They are very similar to the results under the uniform

distributed load. The volumetric strain and the density distribute more evenly than in the

former case.


Figure 2-16. Deformation after loading


.4 *~5*qI

4 ~4E.!l
4 ~1L.!1
~. OJt' II
-.k ....I
-: T~t'II
4 4E'*1
-i h.t
-1 1~54*i


Figure 2-17. Axial stress distribution after loading

TraHl rJitli
-.9 aoE-1

s i- .
-5, OT 1- I 1
-i *;|.i -I

-t j11 .11

Figure 2-18. Volumetric strain distribution after loading

ti M-00
tZ, iiEtDO

+2, 222E4O

+2.3 +00
+2. 56O +00
*2. 56i+SO
+2. 67 +00
+2. 1t OQ

3, DE+tDO0
t6, 5EtaO

Figure 2-19. Density distribution after loading

The difference between the elastic model and plastic model is apparent, as shown

in Figure 2-20, especially for powder materials. Thus, only the elastic model is not

enough for the study of the powders and geomaterials.








0.20 0.25 0.30 0.35

Figure 2-20. Comparison of the results from the elastic and the plastic model

Friction does not play a big role in this case. In the elastic stage, there is almost no

effect from friction. However, only the basic Coulomb friction is implemented in this

chapter. If another more complicated friction model is used, we may see how it works in

more detail.


without friction

with friction

0.05 0.10 0.15
axial strain

0.20 0.25 0.30 0.35

Figure 2-21. Comparison of the cases with and without friction

0.05 0.10 0.15






J 30

" 20

0 -

Figure 2-22 and Figure 2-23 show the density variation of such a point lying at the

center of the half cross section under two different load conditions. They exhibit the same

tendency except at the two ends, that is, at the beginning and at the end of loading that

could cause some unstable results.



E 1.50




0.000 0.001 0.002 0.003
displacement (m)

0.004 0.005 0.006

Figure 2-22. Under displacement control, the density variation at point 1818





0 10 20 30
stress (MPa)

40 50 60

Figure 2-23. Under stress control, the density variation at point 1818

Figure 2-24 shows the general trend of the stress-strain curve. It conforms to the

tendency exhibited in the hardening curve of this material. (See Figure 2-11)







0.00 0.05 0.10 0.15 0.20 0.25 0.30
axis strain

Figure 2-24. Stress-strain curve at point 1818 under displacement control

2.6.3 Case 2 Cup shape sample, uniform displacement on the top

It is available to handle those samples with a little more sophisticated shapes using

Drucker-Prager/Cap model. Figure 2-25, Figure 2-26 and Figure 2-27 show the

distribution of the axial stress, the volumetric strain and the density respectively. Tran et

al did the computations for the models with cylinder shape and T-shape.161 This is

interesting and challenging because either the industrial products or containers have

versatile shapes, sometimes some of them are even weired. Sophisticated shape brings

difficulty to the modeling work. One needs to consider the complicated boundary

conditions and some negative effects caused by the complicated geometry such as stress

concentration. In the example above, we should be aware that around the interior corner

of the cup singularities may exist. We need to use more advanced theory if we want to

explain some uncommon pheonomena. Here we consider only the horizontal and vertical

section. Along the vertical section the density distributes in layers while along the

horizontal section the density is almost homogeneous everywhere. This is reasonable and

was proved by the experiments.[61

-4 96tOoQ
-4. 02 E+00
-). 08E+00
-2 .L4+DD
-J. 2.E+00

+1. L E+00

Figure 2-25. Axial stress distribution after loading

1" 3L-V

Figure 2-26. Volumetric distribution after loading

iti SsBOs-0

t1], U EtDO
4], I Et00
4 1, l 6 E4-0i0
.a, 78s0Q
1. 00 E+00
J] .2Et0- -


Figure 2-27. Density distribution after loading

2.7 Conclusion

Drucker-Prager/cap model can be used to describe the material behaviors during

loading and unloading. It can be used for the models with different shapes shown in

Figure 2-4 and Figure 2-5 or even more sophisticated geometry shapes. Some

experimental data show excellent agreement with the computational results obtained using

this model. However, the plots of the contours of the density show that the distribution of

the density is not uniform. Only the density values in the middle are close to the

experimental data. At the same time the Drucker-Prager/cap model is relatively simple, it

assumes that the material parameters, such as the elastic modulus and friction angle, are

constant. As a matter of fact they are not invariable throughout the loading and unloading

procedure. And it doesn't take into account the time effects. We know that for powder


materials the time effect is very important. Also, it is believed that it is not appropriate for

some kinds of powder like consolidated clays and loose sands. Therefore more advanced

models are expected to be used to simulate more complicated experimental conditions and



3.1 Introduction

For the geomaterials and metallic powders, it is necessary to incorporate the

irreversible volumetric changes and time effects into the constitutive equation. The

uniaxial compression tests are revealing that the behavior of the powder material is always

nonlinear. The volume is first compressible and afterwards dilatant. Just before failure the

volume increases very much. Most of the volume change is irreversible and should be

considered into the constitutive model. Since before failure the volume is dilatant, it is a

property that may be used to predict failure. Time effect is another important factor. It is

apparent in a uniaxial test that the faster the loading rate is, the higher the stress-strain

curve is. That is, for the same deformation different stress levels can be reached.

During the decades of 1980s and 1990s, Dr.N. Cristescu and his colleagues studied

extensively various geomaterials.4'6-10'16-18 They found, first, that the irreversible

volumetric changes must be introduced into the constitutive equations. Secondly, in a

uniaxial compression test the volumetric deformation is non-linear. The volume shows

compressibility at the beginning and dilatancy afterwards. Most of the volume change is

irreversible. The strain rate or loading history has a significant influence on dilatancy,

failure and in a lesser degree on compressibility. Secondly, time effects are also very

important. When loading rate is increased the whole stress-strain curve raises starting

from the lower-value stress and strain. Thirdly, they introduced a very important concept

- the compressibility/dilatancy boundary. Most powder materials are compressible for

small octahedral shear stress, but dilatant for higher values. This boundary indicates

where the material is compressible and where it is dilatant.7'9

3.2 Description of Cristescu Model

To completely describe the constitutive relations for an elasto-plastic material, four

different assumptions are required:

1. A yield function, which generalizes the concept of the yield stress to 2-D and 3-D

stress status.

2. A relationship between the directions of the principal plastic strain increments and

the principal stresses.

3. A flow rule that specifies the relative magnitudes of the incremental plastic strains

when the material is yielding.

4. A hardening law, a relationship between the equivalent stress and the equivalent

plastic strain or the stress works that is done on the material when it is yielding.

Let us see how Cristescu model works.

3.2.1 General Form

The general form of Cristescu model9 is

1 1 1 W(t) F (3.1)
+(- )ol+k 1- (3.1)
2G 3K 2G H(cy) /Bc

in which,

s is the total strain rate.

G and K are elastic moduli.

k is the viscosity parameter.


H(c) is the yield function.

F(Cy,T) is the flow potential function. For alumina powder, F can be replaced with a tensor

evaluated function N(c) which defines the orientation of the viscoplastic strain rate.4

3.2.2 Material Parameters

In this model, the material parameters such as elastic modulus are not constant any

more. Some fitting formulas are obtained from the hydrostatic and triaxial experimental


K(cy) = K paa exp b -b
L Pab (3.2)

E(c) = E pa Pexp(-dG)

where K", E" are the constant limiting values that K and E approach as the pressure

increases. c, P3, b, and d are dimensionless material constants.

Figure 3-1 and Figure 3-2 show how K and E vary with the mean stress. K and E

are both asymptotically approaching to some constant values when the stress level



10000 -

9000 -


0- 7000 -


0 5000-

= 4000-

3000 -

2000 -

1000 -

0 20 40 60 80 100
mean stress(MPa)

Figure 3-1. Bulk modulus variation




2 500 -





0 100 200 300 400 500 600
mean stress(MPa)

Figure 3-2 Elastic modulus variation

3.2.3 Determination of the Yield Function H(y)

For the isotropic materials the yield function H(y) depends on stress invariants only.

It can be determined from the data obtained in drained triaxial compression tests in which

the first stage is hydrostatic and the second is deviatoric. Thus we assume that H(c) has

two parts:

H(c) = Hh (G) + Hd (Y, T) (3.3)

In the hydrostatic stage, Hd(o, c) = 0.

Hh ()= ah ()2 +bh () (3.4)
Pa Pa

Hd (,T) = a(- +a2)-m()9 +B( )2 (3.5)
Pa Pz Pa

where ah, bh, a1, a2, m and B are fitting constants.

3.2.4 Determination of the Viscoplastic Strain Rate Orientation Tensor N(G)

N(c) is a tensor function to describe the orientation of the irreversible strain

increments.4 For isotropic material, N(c) is assumed of the form

N(y) = N1I+N2- (3.6)

where cy' is the stress deviator.


; =k I- (NI+N2 -) (3.7)
H )T

kNi and kN2 can be calculated by the following formulas,

kN, = q(C) + sgn X(o, Tc)Y(c, -c)

cy a-- (3.8)
kN2 =m(_)14 m2(- ) +n2(J)2exp n2(- +r
Pa Pa Pa Pa

in which

p(c) = ak C2 + bkC
X(c,'c) = 0.555 -c (3.9)
S2Y(oT)= z -y )-y exp(zb, )
WQ(,'c)=Zc -2-) 2c --

Pa Pa Pa

where ak, bk, Zc, ya, yc, Zb are material constants.

3.3 Simulation Algorithm for General Loading Test

3.3.1 Description

In a general loading test, one can get the stress-strain curves from different test

conditions such as loading time, loading speed, ambient temperature. The influences of

these conditions are important for evaluating the mechanical behaviors of the material.

Here, the test procedure simulated with the algorithms introduced in this chapter always

starts from a hydrostatic loading. For example, if a hydrostatic loading is carried out up to

C71 = 02 = C3 = 294kPa, afterwards only the axial load is increased while the other two

stress components are kept constant. There are three kinds of algorithms for three

different types of test, the universal continuous loading test, and the creep test. There are

two kinds of creep tests, one considers the transient stage only and the test stops as the

stabilization is reached. In another algorithm a fixed time length is considered at each load

step. In the general test, we don't care about the particular load path. After reaching the

yielding, the irreversible strain rate is

S(1 W)(kNI + kN2 -) (3.10)

Next, we make the inner product with stress tensor a on both sides. Recall the

definition of W, W(t) = f : s'dr On the left side,

W = (1 W)(kN I + kN2 ) (3.11)

in which,

l l = 3C (3.12)

= j (cj C61) = cj 3C2 (3.13)


dW = [3kNkNN +C( -3Z2 3)]dt (3.14)


Integrate it in the time interval [tn, tn+i],

W'n+ ft [3kN, N2 ( 3C2)]dt (3.15)


H ln(1- ) Wn+1 = tn [3kNC + (kN2 3C2)]dt (3.16)

Let FI = 3kNY + kN2 ( o 3C2)

1 Wn+

exp(- 1n FIdt)
H tn



Wn,+ = H+, [1-(1- ) exp(- Fdt)]-
n+l H t

This is the iteration scheme for calculating W at every step. Numerically, the

integration in the formula above can be calculated using the trapezoidal rule.

jtn1 F~dt


- (FI + FIn +)

Table 2-1 lists the parameter values used in the computation.

Table 2-1 Table of material parameters
Confining pressure 294 kPa
Axial stress increment 20 kPa
Time increment 0.01 second
ah 5.833x10-6
bh 1. Ox 10-6
a, 0.27
a2 37
m 8.5
B 4.6x10-5
ze 1.127x10-15
ya 0.842
ye -56.06
Zb -0.01
ak -8.0x10-11
bk 1.641 x 10-6
ml 2.38x1033
m2 -0.038
ni 4.0x10-7
n2 -0.005
r 8.0x 10-4
P3 6.95x105
d 0.002
E __7.0x 105

3.3.2 Result and Analysis

Confining pressure is an important factor in the triaxial test. As shown in Figure 3-

3, if the confining pressure is higher, then higher stress level is needed to get to the same

deformation. At the same time, more volumetric change can be obtained under the higher

confining pressure.

392kPa 700 (kPa)


-0.01 -0.008 -0.006 -0.004 -0.002 0 0.002 0.004

Figure 3-3. Stress-strain curves for different confining pressures in a general test

Figure 3-4 and Figure 3-5 show the influences of some factors, load increment and

time increment. Actually, these two figures express the same thing, the role of the loading

rate Ao/At in the general test. For any time dependent constitutive model like this, we

should show clearly how the time effects are described. In the triaxial test the loading rate

is often adjusted to different value levels as required in order to get different results.

Ac=20kPa Ac=30kPa


200 400 600

Figure 3-4 Irreversible stress work for different load increments, same time increments

At=0.01s At=0.001s


t (kPa)

Figure 3-5.Irreversible stress work for different time increments, same load increments







Viscosity is another important factor. Figure3-6 exhibits its influence on stress work

W. If the viscosity is higher, which means the powder is more viscous, more difficult to

"flow", then more stress work should be done under the same stress level.

100 200 300 400

Figure 3-6.Irreversible stress work under different values of k

3.4 Transient Creep Computation

3.4.1 Description

Creep test is one of the most important tests to study the constitutive relation of

the materials, especially for the time dependent models. In the creep test, the load (stress)

level is kept constant at each step while the deformation still continues. There may be

many loading steps in a creep test. The mechanism of creep is discussed in detail in the

books by Boresi, Malvern and Cristescu.3'10'19

In this section, only the transient stage is considered. The calculation stops when

the deformation becomes stable. When the deformation approaches stabilization, the

irreversible strain rate approaches to zero. This rule is used as a condition to know about

the current state of the material. In order to compute the irreversible strain rate, the

irreversible stress work is required. The difficulty in the algorithm for the transient creep

is the iteration scheme the irreversible strain rate and the irreversible stress work.

A brief computing flow chart is introduced in this section. See Figure 3-7.

In a deviatoric test, the confining pressure keeps constant in the whole test

procedure while the hydrostatic stress and the deviatoric stress keep changing. We need to

know the influence from the hydrostatic stress and the deviatoric stress separately. Thus,

the irreversible stress work can be decomposed into two parts, the volumetric part and the

deviatoric part.

W(t) = J; dT

= t (76 1+ )(1' 6 + E, )d' (3.20)
to j d

Let o (3.21)
Wd =J t dT


W = Wv+ Wd (3.22)


Calculate the stress at the
current step

Calculate a, a', and c for the
current step

Calculate material parameters
for the current step
E(c), K(c), G(c)

Calculate H and W for the
current step

(s =( Nextl)( ad tep

SCalculate totalversible tra instrain

4Yes current step
w- = ft a:FrNext load stepd

Figure 3-7 The algorithm for the transient stage in the creep test

3.4.2 Result and analysis

In this test, after each load increment, the stress will be kept constant until the

material reaches the stabilization. The time to get to the stabilization is different at

different stress levels. Finally, if we connect the successive steps, we can not get a linear

response, as Figure 3-8 shown.


450 -

400- 03 = 294 kPa

350- Acy = 20 kPa

300 --
250 -

200 -

150 -

100 -


0 10 20 30 40 50 60 70 80 90

Figure 3-8 Octahedral shear stress in triaxial test

Figure 3-9 is the typical stress-strain curve that can be seen in many experimental

results. Axial strain si and lateral strain 82 expand in different directions. Volumetric

strain Ev is also presented in this figure. It shows clearly the transition between the

compressibility and dilatancy.

Figure 3-10, Figure 3-11 and Figure 3-12 show how the irreversible stress work

and its components as they vary with the increasing load. Although the total irreversible

stress work increases monotonically as the load increases, its volumetric part shows a

complete different tendency. This illustrates the concept of compressibility and dilatancy.

When the volume of the sample begins to increase, the sign of the volumetric irreversible

stress work becomes negative.

than that of the deviatoric part.

500 -





r (kPa)
200 -





However, the magnitude of volumetric part is much less

-0.05 0 0.05


0.1 0.15 0.2 0.25

Figure 3-9 Stress-strain curves





Wv (kPa) -4 -



-10 -

100 200 300
r (kPa)

400 500

Figure 3-10 Volumetric irreversible stress work as a function of c



100 -

Wd (kPa)


0 ---------- ------------

0 100 200 300 400 500
c (kPa)

Figure 3-11 Deviatoric irreversible stress work as a function of c


120 -

100 -

80 -

W (kPa)
60 -

40 -


0 100 200 300 400 500
c (kPa)

Figure 3-12 Total irreversible stress work as a function of c

In Figure 3-13, we observe clearly the transition from the compressibility to

dilatancy. Also, the confining pressure is one factor that can influence the result. The

higher the confining pressure is, the higher the load is needed to get to the same

deformation. At the compressibility/dilatancy transition points that locate on the rightmost

of the curves respectively, the volumetric change is larger if the confining pressure is


Figure 3-14 shows the deviatoric irreversible stress work for the different confining

pressures. The deviatoric part dominates the magnitude of the irreversible stress work.

As already mentioned, the higher the confining pressure is, the higher load level is needed

to get the same stress work. Moreover, under the same stress level, the results of the

deviatoric part of the irreversible stress work are much different for different confining


392kPa 600(kPa)

294kPa 500




-0.01 -0.008 -0.006 -0.004 -0.002 0 0.002 0.004 0.006 0.008

Figure 3-13 Influence of confining pressure on the volumetric strain

In Figure 3-15, we can clearly observe that the volume firstly decreases, which

means the compressibility. Afterwards it increases, which means dilatancy.










294kPa 392kPa

100 200 300 400

Figure 3-14 Variation of the irreversible stress work as a function of C for different
confining pressure

x 10-3

0 10 20 30 40 50 60 70 80 90

Figure 3-15 Volumetric strain during the loading procedure

Figure 3-16 shows the results under different load increments. The confining

pressure is kept constant (294kPa) while the axial stress increment takes three different

values from low to high, lOkPa, 20kPa, and 30kPa. If the load increment is lower, the

value of ev is larger at the C/D transition point, which means higher compressibility is

obtained. At the same time, the time to get to the C/D boundary is longer if the load

increment is lower. The constitutive relation here is only related to the load increment,

instead of the loading rate. To make this point more straightforward, let us see what will

happen if the load changes dramatically, see Figure 3-17.

X 10-3

3.5 -







20 40 60 80
Time (s)

100 120

Figure 3-16 Volumetric strain for different loading increment and same confining pressure

The curves in Figure 3-17 exhibit the same tendency as that shown in Figure 3-16.

When the load increment is larger, the number of the step is fewer.

X 10-3
3 i--------------------------
2- CD e f

1 B d
A c E

-1 failure


-33 = 294 kPa


0 5 10 15 20 25
Time (s)

Figure 3-17 Volumetric strain for larger loading increment, same confining pressure

The viscousity parameter k may vary due to some reasons(temperature, etc). As

shown in Figure 3-18, if k is smaller, the material is more flowable and longer time is

needed to get to the stabilization.

Figure 3-19 shows the Compressibility/Dilatancy boundary, failure surface, the

curve OH/Wc = 0 and several yield surfaces.

The equation of the Compressibility/Dilatancy boundary is

0.555 c = 0 (3.23)

The equation of the failure surface is

ar + b, c = 0 (3.24)

For this kind of alumina powder, ar = 0.679, br = 9.308.

From = 0, we can derive out

2ah bh am 2 m F2 90
( ) + ( + a2) () = 0 (3.25)
Pa Pa Pa Pa Pa Pa

As a matter of fact, this curve passes in the 2-D plane through the points of the

maxima on the yield surfaces H = constant.

X 10-3

3.5 -

k=10.0 k=.0

Sv 2-

1.5 -


0.5 -

0 10 20 30 40 50 60 70 80
Time (s)
Figure 3-18. Volumetric strain curve for different value of k

The points shown in Figure 3-20 correspond to those in Figure 3-17. It shows the

load path in the whole test procedure, from where the deviatoric test begins, how it passes

through different yield surfaces, the C/D boundary and reach the failure surface. The

small letters correspond to the case with the axial stress increment of Aci=98kPa while the

capital letters correspond to the case with Acy=196kPa. We can see that these two stress

states coincide every second point. Connecting these points we can get a line which has a

slope of 1.4142. It is just the equation of the deviatoric loading path, which is

S- = const (3.26)

Expand (3.26),

( --2= ( +2 2)2 2(Y -2) =2 (3.27)

Equation (3.27) means the confining pressure doesn't change along this loading

path, that is, in the deviatoric test, although cy and -c keep increasing. It is correct for

Karman test in which 02 = 3.

800 C/D bo ndary

600 Failure


H=8.03kPa H=19.22kP

0 500 1000 1500 2000
a (kPa)

Figure 3-19 Compressibility/Dilatancy boundary, failure surface 8H/8y = 0, and several
yield surfaces


Failure C/D boundary

k F
& 400 1

200- e0C Arctan(1.4142
c, B
0 A
0 500 1000 1500 2000
cy (kPa)

Figure 3-20 Deviatoric test procedure in T-c- plane

3.5 Analytical Solution for the Creep Test with Fixed Time Interval

3.5.1 Description

In this type of creep test, we assume that at time to a load starts. The stresses

increase quite fast from the previous stress level and afterwards they are kept constant for

a fixed time interval. More loading steps like this repeat with the increase of the stress

level if required. A new formula can be obtained.6

1 W =(t) 1 (- exp -(kN 1 kN2 ).(to- t) (3.28)
H(y,T-) [ H(cy, -) 1 1 TI

where o(t) = o(to) = constant and W(to) is the initial value of W at time to. If the

initial state is on a yield surface then W(to) = H(co). The formula for calculating strain at

arbitrary time is

I L W(to kN + kN y(t0
s(t)= sE + exp (kN,1+kN, (to-t)
l(kN + kN2 ). L


with the initial conditions

= 0, W(t0) Wo

I 3K 2G 2G

This time the time interval during which the load at each step is constant. It can be

set up at different magnitudes.

3.5.2 Result and Analysis

Figure 3-21 exhibits the same story as Figure 3-17 does. For bigger value of k, the

time is enough for the material to reach the stabilization while it is not enough for smaller

k. Thus we see two different curves, one with the apparent creep stage while the other

only with the transient stage.

Figure 3-22 shows the results for different step time intervals. If the time interval

is relatively small, we can get the transient creep only.

X 10-3 ..


2 0k=.1

1 -

Sv 0


At = 10 sec
-2 -

0 10 20 30 40 50 60 70 80

Figure 3-21 Volumtric strain change for different k values

x 10-3

3 At=3sec


0 10 20

40 50 60 70 80 90

Figure 3-22 Volumetric strain change for different step time intervals


3.6 Conclusion

Cristescu model is in reasonable qualitative agreement with the experiment

results,4'6'10 although a quantitative comparison is not shown here. It describes the time

effects in the compaction and dilatant process. A more complete discussion could be done

if it could be combined with FEM methods. This will be introduced in next chapter.


4.1 General

Cristescu model is a time-dependent semi-linear model. Usually, the variable

stiffness method is the most acceptable method to solve the non-linear problem.24 In this

method an effective iteration scheme needs to be erected to describe the constitutive

relation. Any computation process should proceed in an incremental manner considering

suitable small time intervals. In each time interval the stress level should be known so that

the increments of strain can be found.

Owen and Hinton present the formulation of elastic/viscoplasticity theory in

discrete time form in the book. They provide the detailed algorithm and an example

program for solving some simple problems.20 Jin presented his algorithm in his dissertation

for alumina powder combining with a large amount of experimental work in his research.17

These are all good reference for designing the FEM algorithm for powder materials.

4.2 Algorithm for Stress Control Test

In stress control test, the loading rate or the stress increment is known. The

unknown is the strain rate or the strain increment.

Let us consider an arbitrary load step and define it as the n+lth step between the

state n and state n+1. According to the basic assumption,

A+= A n+ As n+(4.1)

The elastic part can be calculated using Hooke's Law,

Asn+ (D n+) 1A~n+1 (4.2)

in which Dn+1 is the elastic matrix at the n+lth step. In Cristescu model, the

elastic moduli are not constant. D is the function of 7.

Implement the implicit Euler method to get the irreversible strain increment,

AsI = At (l- 0)sn+s0n+ (4.3)

0 is a constant and 0 < 0 < 1. To get the stable result, 0 should satisfy 0 < 0 < 1.

Plug (4.2) and (4.3) into (4.1),

n+1 Atn+i = D -1A + At(1 0)sn+s (4.4)

Equation (4.4) forms the iteration scheme of two successive steps.

Define two state variables p and q.

qn = D-'Ac + At(1-6)sn+sn (4.5)

pn = n At0sn (4.6)

Then equation (4.4) changes into

Pn+l = qn (4.7)

If the state variables at the nth step are given, that is, tn,sn,Cn,sn,Sn,pn are


Conduct another loop from the nth step to the n+lth step. Now consider the ith step

in this loop.

q = p, + P 6s, (4.8)

From equation (4.6),

S= 1 (4.9)

So, equation (4.8) changes into,

6s, = qn p, (4.10)


1+1 = + -6 ,
P +1 = 81 AtOg1+1

The loops stops when 6Ei < tolerance and then continue the n+2th step.

4.3 Algorithm for Strain Control Test

In strain control test, the strain rate or the strain increment is known. The

unknown is the loading rate or the stress increment.

Let us consider an arbitrary load step and define it as the n+lth step between the

state n and state n+1. From equations (4.1), (4.2), (4.3)

D lo +1 + AtOn+ = A At(1 0)+D l c (4.12)

Equation (4.12) forms the iteration scheme of two successive steps.

Define two state variables p and q.

q. = As At(1 0)n + D 'A (4.13)

p. = D c +AtO n (4.14)

Then equation (4.12) changes into

Pn.+ = qn (4.15)

If the state variables at the nth step are given, that is, tn,sn,CY,sn,Sn,pn are


Conduct another loop from the nth step to the n+lth step. Now consider the ith step

in this loop.

qn = p, + O o, (4.16)

From equation (4.14),

p D 1 + AtO (4.17)

So, equation (4.16) changes into,

P' = q n P, (4.18)


Pl+1 = ,1+1 AtOsg+,

The loops stops when 60i < tolerance and then continue the n+2th step.

4.4 Algorithm for ABAQUS

Since ABAQUS is hopefully to be implemented with Cristescu model, we need to

pay more attention to the characteristics provided by ABAQUS. Besides those internal

constitutive models, ABAQUS provides user subroutine as the interface with the users.

User subroutine is a powerful and flexible tool for analysis because users can define the

constitutive relations by themselves. UMAT is one of these subroutines that could be

used. First, let us look at the interface defined in this subroutine.


Figure 4-1. The interface of UMAT subroutine

In Figure 4-1, we notice that there are many parameters defined for this subroutine.

Among them, the most important parameter is DDSDDE that is elasticity (Jacobian)

matrix of the constitutive model. DDSDDT is defined in an incremental form,

OAA/0As. It has to be updated every iteration step. Figure 4-2 is my planning algorithm.

This method actually uses the predictor-corrector method. In detail, first we need the

elasticity matrix from the last step to calculate the new stress level, that is, to get a pseudo

new stress. Then, calculate values the modulus and construct the new elasticity matrix.

Recalculate the new stress with this new elasticity matrix. Next, consider the constitutive


equation to find the strain rate. Implementing the classical numerical methods such as the

trapezoidal method, we can get the strain tensor in the current step.

The algorithms introduced in Section 4.2 and Section 4.3 cannot be used with

ABAQUS. Another algorithm is designed, as shown in Figure 4-2. Basically, it uses the

predictor-corrector method to calculate the elasticity matrix.

No results from any of these algorithms can be presented in this dissertation.

Figure 4-2. Flow chart of the algorithm for Cristescu model combining with ABAQUS


[1] Aydin I., Briscoe, B.J., Sanliturk, K.Y., "The internal form of compacted ceramic
components: a comparison of a finite element modeling with experiment," Technology and
Medicine, Powder Technology, 1996, Vol.89.

[2] Bickford W.B., Finite Element Method (2nd Edition), Richard D.Irwin, Inc., Boston,
1990 and 1994.

[3] Boresi, A.P., Schmidt, R.J., Sidebottom, O.M., Advanced Mechanics of Materials (5th
Edition), John Wiley & Sons, Inc, New York, 1988.

[4] Cazacu, 0., Jin J., Cristescu, N. D., "A New Constitutive Model for Alumina Powder
Compaction," KONA No. 15, 1997.

[5] Chotourou, H., Gakwaya, A., Guillot, M., and Hrairi, M., "Implementing a cap
material model for the simulation of metal powder compaction," AMD-Vol.216, Net
Shape Processing of Powder Materials, ASME 1995.

[6] Cristescu, N.D., "A procedure for determining the constitutive equations for materials
exhibiting both time-dependent and time-independent plasticity," Int. J. Solids Structures,
1972, Vol.8.

[7] Cristescu, N.D., "A procedure to determine nonassociated constitutive equations for
geomaterials," Int. J. Plasticity, 1994, Vol.10.

[8] Cristescu, N.D., "Plasticity of Porous and Particulate Materials," Journal of
Engineering Materials and Technology, April 1996, Vol.118.

[9] Cristescu, N.D., Cazacu, 0., and Jin, J., "Constitutive equation for compaction of
ceramic powders," IUTAM Symposium on Mechanics of Granular and Porous Materials,
University of Cambridge, July 15-17, 1996.

[10] Cristescu, N.D., Hunsche, U., Time Effects in Rock Mechanics, John Wiley & Sons,
Chichester, 1998.

[ 11] Desai, C.S., Siriwardane, H.J., Constitutive Laws for Engineering Materials i ilh
Emphasis on Geologic Materials, Prentice-Hall, Inc., Englewood Cliffs, New Jersey,

[12] DiMaggio, F.L., Sandler, I.S., "Material model for granular soils," Journal of the
Engineering Mechanics Division, 1971, Vol.97.

[13] Drucker, D.C., and Prager, W., "Soil mechanics and plastic analysis or limit design,"
Quarterly of Applied Mechanics, Vol.10, 1952.

[14] Hibbitt, Karlsson and Sorensen, Inc., ABAQUS Theory Manual Version 5.4, 1994.

[15] Hibbitt, Karlsson and Sorensen, Inc., ABAQUS User's Manual Version 5.4, 1994.

[16] Jin, J., "Elastic/Viscoplastic Models for Geomaterials and Powder-like Materials with
Applications," Ph.D. Dissertation, University of Florida, 1996.

[17] Jin, J., Cazacu, 0., and Cristescu, N.D., "An elastic/viscoplastic consolidation model
for alumina powder," Proc. 5th World Congress of Chemical Engineering, Vol.5, July
1996, San Diego.

[18] Jin, J., Cristescu, N.D., "Experimental study on the consolidation of alumina
powder," Preprint ASME Mechanics and Materials Conference. The Johns Hopkins
University, 1996.

[19] Malvern, L.E., Introduction to the Mechanics of a Continuous Medium, Prentice-
Hall, Inc., Englewood Cliffs, New Jersey, 1969.

[20] Owen, D.R.J. and Hinton, E., Finite Elements in Plasticity: Theory and Practice,
Pineridge Press, Swansea, U.K., 1980.

[21] Schofield, A., Wroth, P., Critical State Soil Mechanics, McGraw-Hill, London,

[22] Tran, D.V., Lewis, R.W., Gethin, D.T., and Ariffin, A.K., "Numerical Modeling of
Powder Compaction Processes: Displacement Based Finite Element Method," Powder
Metallurgy, 1993, Vol.36, No.4.

[23] Watson, S.A., Adams, M.J., Rough, S.L., Briscoe, B.J., and Papathanasiou, T., "The
modeling of the influence of wall friction on the characteristics of pressed ceramic parts,"
IUTAM Symposium on Mechanics of Granular and Porous Materials, Kluwer Academic
Publishers, Netherlands, 1997.

[24] Zienkiewicz, O.C., The Finite Element Method in Engineering Science (The second
edition, expanded and revised, edition of The Finite Element Method in Structural and
Continuum Mechanics), McGraw-Hill, London, 1971.


Weijian Wang was born on May 10, 1971, in Hebei Province, People's Republic of

China. He received his B.E. degree in the Department of Engineering Mechanics,

Tsinghua University in July 1993. Afterwards he worked as an assistant engineer in China

Aerospace Agency. In May 1997, he was admitted into the Department of Aerospace

Engineering, Mechanics, and Engineering Science, University of Florida. Two years later,

in August 1999, he received his M.S. degree from University of Florida.

University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - - mvs