J. agric. Engng Res. (1996) 63, 121 – 128
An Endochronic Constitutive Model for Grain En Masse S. Xu; Q. Zhang; M. G. Britton Department of Agricultural Engineering, University of Manitoba, Winnipeg, Manitoba, Canada R3T 5V6 (Receiy ed 16 September 1994; accepted in rey ised form 18 September 1995)
An endochronic model has been develeloped to predict the stress – strain responses of grain en masse. The model uses intrinsic time to trace the effect of deformation history on the constitutive behaviour of grain en masse. A stress-based free energy density representation has been proposed to relate strains to stresses. Predicted stress – strain behaviour compared well with triaxial test data reported in the literature.
and Peters7 developed an endochronic model for granular materials. The model was shown to be adequate for predicting stress – strain behaviour of wheat en masse.8
÷ 1996 Silsoe Research Institute
1. Introduction
The finite element method has been shown to have great potential in predicting structural and functional behaviour of bulk grain storage systems.1 One of the most important steps in finite element predictions of behaviour of grain storage systems is constitutive modelling of grain en masse. A constitutive model describes the responses of the material to deformations or loads. Various non – linear, elastoplastic and elasto-viscoplastic theories have been applied to grain en masse.2–4 Endochronic theory, founded on the irreversibility principle of thermodynamics, provides a unified approach to describe the plastic behaviour of materials. The word endochronic means internal (endo) time (chronic). The theory uses a monotonically increasing time-like variable, termed ‘‘intrinsic time’’, to describe the loading history-dependent behaviour of materials. Intrinsic time is different from the real time. An endochronic model predicts timeindependent behaviour of materials unless the real time is included as a component in the definition of intrinsic time. Intrinsic time was first defined by Valanis5 as an equivalent length of the total strain path. He later modified the definition of intrinsic time as an equivalent length of the plastic strain path.6 Based on this definition of endochronic time, Valanis 0021-8634 / 96 / 020121 1 07 $18.00 / 0
121
Notation A Ad Ah A(r) C (r) C0 C1 D (r) d» h d» ph deij de pij dz dzd dzh E (ijr ) e pij Fd Fh G H (r) J J1 K k (r) L n
free energy density deviatoric component of free energy density hydrostatic component of free energy energy stored in spring element r material constant material parameter material parameter material constant hydrostatic strain increment plastic hydrostatic strain increment deviatoric strain increment tensor plastic deviatoric strain increment tensor intrinsic time increment deviatoric component of intrinsic time increment hydrostatic component of intrinsic time increment internal deviatoric strain plastic deviatoric strain tensor deviatoric hardening function hydrostatic hardening function elastic shear modulus internal hydrostatic strain hydrostatic kernal function simplified hydrostatic kernel function (constant) elastic bulk modulus spring constant of spring element r Laplace transformation operator number of terms in kernel function or number of elements in mechanical model ÷ 1996 Silsoe Research Institute
122
S. XU ; Q. ZHANG ; M. G. BRITTON
p q (r) r Q (r) Q (ijr ) s0 sij s` z9 a
b (r) b1 G G1 » ph » pij k l (r) s s (r) s0 sh s ij s` Fd u u i i
Laplace transformation parameter frictional limit of frictional element r , or internal variable running index hydrostatic internal variable deviatoric internal variable deviatoric stress at yield under triaxial loading deviatoric stress tensor deviatoric stress at failure under triaxial loading dummy variable in integration material parameter (slope of (s` 2 s ) vs. zd curve plotted in a semi-log scale) intermediate variable material parameter coupling kernel function simplified coupling function (constant) plastic hydrostatic strain plastic strain tensor coupling parameter intermediate variable stress stress on spring element r confining pressure in triaxial test (hydrostatic stress at zh 5 0). hydrostatic stress stress tensor hydrostatic stress at zh 5 ` deviatoric kernel function absolute value norm of tensor
The Valanis and Peters model7 describes the stress responses of materials to deformation, i.e. it is a strain-based theory. An alternative way of characterizing constitutive behaviour of materials is to relate strains to stresses (stress-based theories). If the response of a material is reversible, meaning no energy is dissipated during deformation or loading, the result from a stress-based model should be identical to that from a strain-based model. For dissipative materials, strain- and stress-based models are complementary, i.e. a strain-based model offers an upper bound to the solution, whereas a stress-based model gives a lower bound in limit analysis.9 Energy dissipation occurs in the grain mass because of interparticle friction. Therefore, predictions by a strain-based theory might not be identical to those of a stress-based theory for grain en masse, depending on loading conditions. Finite element algorithms based on strain-based constitutive laws use displacements as primary variables and forces are calculated from displacements as
secondary variables. Most existing finite element models fall in this category, partly owing to the lack of stress-based constitutive laws for grain en masse. In grain storage systems, however, forces (pressures) may be of primary interest in many cases. For example, pressures and forces dictate structural design of grain bins. To date, no stress-based endochronic models have been developed for grain en masse. The object of this study was to develop and validate a stress-based endochronic formulation for grain en masse and compare it with the strain-based theory of Valanis and Peters.7 This stress-based constitutive model will provide a basis for formulating force-based finite element algorithms for predicting structural and functional behaviour of bulk grain storage systems.
2. Model development
2.1 . Intrinsic time
For granular materials which dilate during shearing, Valanis and Peters7 proposed the following definition of intrinsic time based on experimental observations dz 2 5 de pijde pij 1 k 2(d» ph)2
(1)
where: dz is intrinsic time increment, de pij is plastic deviatoric strain increment tensor, d» ph is plastic hydrostatic strain increment, k is a coupling parameter. The two plastic strain components are determined from the total strain increments as follows de pij 5 deij 2
dsij ds h and d» ph 5 d» h 2 2G K
(2)
where: G is elastic shear modulus, K is elastic bulk modulus, deij is deviatoric strain increment tensor, d» h is hydrostatic strain increment, dsij is deviatoric stress increment tensor. Corresponding to the two strain components, intrinsic time is decomposed into a deviatoric component and a hydrostatic component7 dzd 5
dz dz and dzh 5 Fd k Fh
(3)
where: dzd is deviatoric component of intrinsic time increment, dzh is a hydrostatic component of intrinsic time increment, Fh is the hydrostatic hardening function, Fd is the deviatoric hardening function.
123
ENDOCHRONIC CONSTITUTIVE MODEL FOR GRAIN
2.2 . Stress – strain relationship By defining an appropriate form of free energy stored in the material, strains may be related to stresses as follows
» pij 5
A s ij
(4)
1 1 1 A(r) 5 » (r)s (r) 5 (r) [s 2 q (r)]2 2 2k
(5)
The total free energy stored in the system is the sum of energy stored in all spring elements
OA n
Variables q
r 51
(0)
O k1 [s 2 q n
(r)
???q
5 1–2
(n)
r 51
(r)
(r) 2
]
(6)
q(0) q(n) k(1)
and
(r)
i sij 2 Q (ijr ) i 2
(r)
us h 2 Q
r 51
(7)
n
Ah 5
1 – 2
r 51
k(2)
e pij 5
Ad Ah and » ph 5 sij s h
σ
(8)
where: e pij is plastic deviatoric strain tensor, e ph is plastic hydrostatic strain. By comparing Eqns (6) and (7), it can be seen that Q (ijr ) and Q (r) represent ‘‘internal stresses’’ which govern energy dissipation of the material. The corresponding internal strains can be determined by differentiating the free energy densities [Eqn (7)] with respect to internal stresses E (ijr ) 5 2
Ad Ah and H (r ) 5 2 Q (ijr ) Q (r)
SdEdz D . 0
Q (ijr )
q(2) q(1)
(9)
(r ) ij
SdHdz D . 0 (r)
and Q (r)
(10)
h
Substituting Eqn (9) into (10) yields 2Q (ijr )
k(n)
q(n)
u
where: Ad and Ah are deviatoric and hydrostatic components of free energy, respectively, C (r) and D (r) are material constants, sij is deviatoric stress tensor, s h is hydrostatic stress, Q (ijr ) and Q (r) are internal variables, u u denotes an absolute value, i i denotes the norm (‘‘length’’) of tensor. Following Eqn (4), the two strain components are calculated from the derivatives of free energy with respect to the corresponding stresses
d
q(2)
(r) 2
where: E (ijr ) is internal deviatoric strain, H (r) is internal hydrostatic strain. To satisfy the condition of positive internal dissipation for an irreversible process,7 the dissipation rate, defined as the product of the internal stress and the rate of internal strain, must be positive, i.e.
σ
q(1)
Ad 5 1–2
are internal variables, which,
σ
q(0)
OC OD n
where: » pij is plastic strain tensor, s ij is stress tensor, A is free energy density. Grain en masse is elastoplastic in general. A simple mechanical model (Fig. 1 ) consisting of plastic (friction) and elastic (spring) elements is introduced to illustrate the calculation of free energy for grain en masse under uniaxial loading conditions.10 The maximum stresses that frictional elements can sustain (frictional limits) are q (0) , q (1) , q (2) , . . . , q (n) , and spring constants of spring elements are k (1) , k (2) , . . . , k (n). When the system is loaded by a stress s , energy is dissipated by friction elements and stored in spring elements. For an individual spring element r , the stress acting on it s (r) is equal to the total stress s less the stress on the friction element q (r) , i.e., s (r) 5 (s 2 q (r)) (Fig. 1 ) , and the corresponding strain of the spring element is » (r) 5 s (r) / k (r). Therefore, the energy stored in a single spring element A(r) is determined as
A5
together with stress s , define the state of the material. When the applied stress reaches the frictional limit q (0) , plastic flow occurs (the failure state is reached). Corresponding to the two strain components, the total free energy is assumed to be composed of a deviatoric component and a hydrostatic component, each is determined by generalizing Eqn (6)
ε
Fig. 1. A mechanical model for elastoplastic granular materials
2Q (r)
S D d A S D.0 dz Q
d Ad .0 dzd Q (ijr ) h (r)
(11) (12)
h
The above inequalities can be satisfied mathemati-
124
S. XU ; Q. ZHANG ; M. G. BRITTON
cally in many possible ways. A simple way is to express Q (ijr ) and Q (r) as positive linear functions of internal strain rates Q (ijr ) 1 a (11r ) Q (r) 1 a (22r )
S D d A S D50 dz Q
d Ad 50 dzd Q (ijr )
(13)
h (r)
(14)
h
where a (11r ) and a (22r ) are positive intermediate variables. Equations (13) and (14) are the evolution equations for deviatoric and hydrostatic internal variables, respectively. These equations assume that deviatoric and hydrostatic responses are independent of each other. For grain en masse, however, shearing may cause volume changes (dilatancy). Therefore, Eqn (14) is modified to account for shear-volume coupling by including shear work, which is calculated as the product of the shear strain e pij and the shear stress rate (dsij / dz ) integrated over the range of hydrostatic intrinsic time7 Q (r) 1 a (22r )
S D
d Ah 5 2a (21r) dzh Q (r)
E e Sdsdz D dz9 (15) zh
ij
p ij
where a is a positive intermediate variable. Substitution of Eqn (7) into Eqns (8) and (13), and performing the Laplace transformation on the resulting equations gives
#e pij 5
OC
# Q
(r ) ij
# (ijr )] [#sij 2 Q
(r)
r 51
(16)
(r)
(r ) ij
(17)
where p is a transformation parameter; and an overhead bar (x# ) indicates that the variable is in the # (ijr ) and transformed domain. Solving Eqn. (17) for Q substituting it into Eqn. (16) yields
O C F1 1 a 1 C pG p1 ps# n
#e pij 5
(r)
(r ) 11
r 51
ij
(r)
(18)
Using the convolution theorem of Laplace transformation expressed by Eqn (19)11
HE
zd
L
J
f1(zd 2 z 9)f2(z 9) dz 9 5 Lh f1(zd )jLh f2(zd )j (19)
0
the inverse Laplace transformation of Eqn (18) gives e pij 5
E
zd
Fd (zd 2 z 9)
0
OC n
Fd 5
dsij dz 9 dz 9
[1 2 e 2b
(r)
r 51
b (r) 5
1 a C (r) (r ) 11
(r)z d
]
(23)
where: C 0 , C 1 and b 1 are material parameters. Under one-dimensional conditions, Eqn (23) is a mathematical representation of a mechanical model consisting of a friction element and a friction-spring unit, as shown in Fig. 1 when n 5 1 . Here q (1) and q (0) define the initial yield point and the failure point, respectively, and zd represents a deviatoric plastic strain scale. Following a similar procedure, a constitutive equation similar to Eqn. (20) is obtained for the hydrostatic response
» ph 5
E
zh
ds h dz 9 1 dz 9
J (zh 2 z 9)
0
E
zh
G(zh 2 z 9)e pij
0
O D [1 2 e a G5 O [1 2 e a kF n
J5
2l (r)zh
(r)
dsij dz 9 dz 9 (24) (25)
]
r 51
n
# )50 2 a C p (#sij 2 Q (r ) 11
Fd < C 0 zd 1 C 1(1 2 e 2b 1zd)
0
(r ) 21
n
where: z 9 is a dummy variable, Fd is the deviatoric kernel function, b (r) is an intermediate variable. Equation (20) may be viewed as the summation of an incremental constitutive equation de pij 5 Fd dsij over the entire intrinsic time history, where Fd is equivalent to a tangent compliance. For materials with an initial yield point and failure point, the deviatoric kernel function expressed by Eqn (21) may be approximated by one linear term and one exponential term
r 51
(r ) 21 (r ) 22 h
l (r) 5
2l (r)zh
]
1 a D (r ) (r ) 22
(26) (27)
where: J is the hydrostatic kernel function, G is the coupling kernel function, l (r) is an intermediate variable. For granular materials which are subject to large plastic deformations and negligible elastic deformations, the above equations may be simplified12 as J < J1 zh 21
G < G1(k Fh ) zh
(28) (29)
(20)
where: J1 is a simplified hydrostatic kernel function (constant), G1 is a simplified coupling function (constant).
(21)
3. Determination of model parameters
(22)
Model parameters for the stress-based theory were evaluated using triaxial test data reported by Zhang et
125
ENDOCHRONIC CONSTITUTIVE MODEL FOR GRAIN
al.4 for wheat at a moisture content of 8?1% w.b. and a bulk density of 817 kg / m3. Their data, collected from conventional cylindrical triaxial tests, include four sets of data from axial compression tests at confining pressures of 20?7, 34?5, 48?3 and 62?1 kPa, and one set from hydrostatic compression tests. Axial compression data for confining pressure of 48?3 kPa and hydrostatic compression data were used in determining model parameters. Because the definition of intrinsic time is the same as that used by Valanis and Peters7 in their strain-based theory, all parameters associated with the intrinsic time remain identical to those determined by Xu et al.8 for the grain (Table 1). Detailed procedures for determining these parameters can also be found in Ref. 8. Determination of parameters associated with the kernel functions are discussed in the following sections. Under triaxial loading conditions, the deviatoric stress may be expressed as7 s 5 s` 1 (s0 2 s`)e
2a zd
(30)
Fd 5 s 0 1 bs
de p de p 5 Fd < Fd dzd dz
(31)
Table 1 Model parameters determined for wheat at a moisture content of 8?1% w.b. and a bulk density of 817 kg / m3 2067† 2750† 1?13† 0?7† 8?8 3 1025† 73?2† 27 0?154 2?92 3164 8?0 3 1024
† Reported by Xu et al.8 Detailed procedures for determining these parameters can also be found in Ref. 8.
(32)
where b has a constant value of 0?408. Combining Eqns (31) and (32) yields 7
de p 5 (s 0 1 bs )dzd
(33)
Substituting Eqn (30) into (33) and integrating the resulting equation gives e p 5 (s 0 1 bs `)zd 2
b (s` 2 s0)(1 2 e 2azd) a
(34)
The deviatoric strain can also be calculated from Eqns (20) and (23) as follows e p 5 (C 0 zd 1 C 1)(s 2 s0) 2 a (s` 2 s0)(I1 1 I2) (35)
E z9e I 5E e zd
I1 5
az 9
dz 9 5
0
zd
where: s is the deviatoric stress under triaxial loading, s` is the deviatoric stress at failure (at zd 5 ` ) under triaxial loading, s0 is the deviatoric stress at yield (at zd 5 0) under triaxial loading, a is a material constant. When the material is subjected to axial compression as in the conventional triaxial tests, the plastic volumetric strain increment d» p… is negligible compared with the plastic deviatoric strain de p , therefore de p / dz < 1 as indicated by Eqn (1) (de pij 5 de p under axial compression, thus de pijde pij 5 (de p )2). It follows from Eqn (3) that
Elastic bulk modulus (K ) , kPa Elastic shear modulus (G ) , kPa Coupling parameter (k ) Hydrostatic hardening function (Fh ) F1 a c Simplified coupling function (G1) , kPa21 Hydrostatic kernel function (J1) , kPa21 Deviatoric kernel function (Fd ) C0 , kPa21 b1 C1 , kPa21
Under axial compression, the deviatoric hardening function Fd may be expressed as7
b 1zd 1(b 12a )z 9
2
zd 2azd 1 2azd e 2 2 (e 2 1) a a
(36)
dz 9
0
3
E
1 (e 2azd 2 e 2b 1zd) b1 2 a zd e
if b 1 ? a
2a zd
(37)
if b 1 5 a
If b 1 ? a , comparison of Eqn (35) with (34) results in e 2azd 5 e 2b 1zd
(38)
which indicates that b 1 must be equal to a . Eqn (30) suggests that a is the slope of (s` 2 s ) vs. zd curve plotted on a semi-log scale. The a value was determined to be 3164 by Xu13 for the wheat used in this study. Substituting Eqn (30) into (35) and rearranging the equation yields
S
e p 5 C 0 zd (s` 2 s0) 2 C1 2
S
2 (s` 2 s0) C 1 2
D
C0 (s` 2 s0) a
D
C0 1 C 1 zd e 2azd a
(39)
Differentiating Eqn (39) with respect to zd , and then letting zd approach infinity, Eqn (40) is obtained C0 5 lim
zd 5`
Fs 21 s dzde G p ij
`
0
(40)
d
Following a similar procedure, J1 in the hydrostatic kernel function is determined as J1 5 lim
zh 5`
Fs 21 s ddz» G p h
`
0
h
(41)
126
S. XU ; Q. ZHANG ; M. G. BRITTON
de p dzd
U
5 C 1(a 2 1)(s` 2 s0)
(42)
zd 50
Another expression for strain derivative may be obtained from Eqn (33) de p dzd
U
5 a 0 1 bs 0
(43)
zd 50
Equating Eqns (42) and (43) yields C1 5
C0 2 b a 21
(44)
The coefficient C1 is calculated to be 8?0 3 1024 using Eqn (44). Parameter G1 was calculated to be 27 kPa21 using Eqn (24) at » ph 5 0 .
0·04 0·03 Hydrostatic strain
where: s ` is hydrostatic stress at zh 5 ` , s 0 is confining pressure in triaxial test (at zh 5 0). The failure and yield stresses and strain rates at failure determined from data reported by Zhang et al.4 are summarized in Table 2. The strain rate (d» ph / dzh ) at zh 5 ` was determined from the slope of the curve for hydrostatic strain versus hydrostratic intrinsic time at a large zh value (Fig. 2 ) (the slope remained almost constant for zh . 0?3). The strain rate (de pij / dzd ) at zd 5 ` was determined in a similar fashion from plot of deviatoric strain versus deviatoric intrinsic time. The hydrostatic stress s ` was calculated as (s 0 1 0?408s`).7 From these data C 0 and J1 were determined to be 2?94 kPa21 and 0?154 kPa21, respectively. Differentiating Eqn (39) with respect to zd , and then letting zd approach zero yields
Confining pressure (s 0) , kPa Deviatoric stress failure (s `) , kPa Hydrostatic stress at failure (s `) , kPa Deviatoric stress at yield (s0) , kPa Slope (a ) Deviatoric strain rate at failure (de pij / dzd ) Hydrostatic strain rate at failure (d» ph / dzh ) † Reported by Xu et al.8
48?3 57?8† 71?9 33?4† 3164 71?6 3?62
1
–0·02 0
0·5
1
1·5
2
Hydrostatic intrinsic time
Fig. 2. Variation of hydrostatic plastic strain (» ph) with hydrostatic intrinsic time (zh)
tions, shown in Figs 3 and 4 , compared favourably with the experimental data reported by Zhang et al.4 The predicted deviatoric stresses and volumetric strains closely followed the experimental data with average differences of 0?83 kPa and 1?08 3 1023, respectively. The proposed stress-based theory is in good agreement with the strain-based theory of Valanis and Peters7 (Figs 3 and 4 ). It is interesting to note that predictions of volumetric strains averaged from stress theory and strain theory are in closer agreement with the experimental data than either theory alone (Fig. 4 ).
5. Conclusions 1. The stress-based endochronic theory adequately predicts the stress – strain behaviour of grain en masse under axial compression. Predictions of
30 25 Deviation stress, kPa
Table 2 Failure stresses and strain rates derived from triaxial test data for wheat at a moisture content of 8?1% w.b. and bulk density of 817 kg / m3 (Ref 9)
0·01
–0·01
4. Model verification The model parameters listed in Table 1 were used to predict the stress – strain behaviour of wheat en masse under axial compression with a constant confining pressure of 20?7 kPa. The integrations of Eqns (20) and (24) were carried out numerically. Model predic-
0·02
20 15 10 5 0
0
0·02
0·04
0·06
0·08
0·1
Deviatoric strain
Fig. 3. Comparison of predicted dey iatoric stresses with data reported by Zhang et al.4 for triaxial loading at a confining pressure of 20?7 kPa. Predictions by the strain-based theory were reported by Xu et al.8 $, Measured; ——, stress theory; - - - - -, strain theory
ENDOCHRONIC CONSTITUTIVE MODEL FOR GRAIN
0·005 2
Volumetric strain
0·003 3
0·001 –0·001
4
–0·003 5
–0·005 0
0·02
0·04
0·06
0·08
0·1
Deviatoric strain
Fig. 4. Comparison of predicted y olumetric strains with data reported by Zhang et al.4 for triaxial loading at a confining pressure of 20?7 kPa. Predictions by the strain-based theory were reported by Xu et al.8 $, measured; —— , stress theory; - - - - , strain theory , n , ay erage
6
7
8
deviatoric behaviour is closer to experimental data than predictions for volumetric behaviour. 2. Predictions by the stress-based endochronic theory are consistent with those of the strainbased theory of Valanis for Peters7 for wheat en masse under axial compression loading.
9
10
Acknowledgements We thank the Natural Science and Engineering Research Council of Canada for financial assistance through a Research Grant.
11
12
References 13 1
Puri V M; Manbeck H B Potential of finite element method in modeling load responses of particulate
127
materials. ASAE Paper 91 – 4076 1991. St. Joseph: American Society of Agricultural Engineers Smith D L O; Lohnes R A Bulk strength and stress – strain behaviour of four agricultural grains. Journal of Powder and Bulk Solids Technology 1983, 7(4): 1 – 7 Li Y; Puri V M; Manbeck H B Elastic – viscoplastic cyclic constitutive model parameter determination and evaluation for wheat en masse. Transactions of the ASAE 1990, 33(6): 1984 – 1995 Zhang Q; Puri V M; Manbeck, H B Determination of elastoplastic constitutive parameters for wheat en masse. Transactions of the ASAE 1986, 29(6): 1739 – 1746 Valanis K C A theory of viscoplasticity without a yield surface, Part I: General theory. Archives of Mechanics 1971, 23(4): 517 – 533 Valanis K C Fundamental consequences of a new intrinsic time measure, plasticity as a limit of the endochronic theory. Archives of Mechanics 1980, 32(2): 171 – 191 Valanis K C; Peters J F An endochronic plasticity theory with shear-volumetric coupling. International Journal for Numerical and Analytical Methods in Geomechanics 1991, 15: 77 – 102 Xu S; Zhang Q; Britton M G An endochronic finite element model for predicting loads in grain storage structures. Transactions of the ASAE 1993, 36(4): 1191 – 1199 Gao Y On the complementary bounding theorems for limit analysis. International Journal of Solid Sand Structures 1988, 24(6): 545 – 556 Xu S; Zhang Q; Britton M G Physical interpretation of endochronic constitutive models for granular materials. Paper No. 94-4028. 1994, ASAE, St. Joseph, USA Hirsch K A Handbook of mathematics. New York: Van Nostrand Reihold, 1985 Valanis K C; Fan J; Nashiro T Endochronic plasticity: theory and comparison with experiment. In Computer Modelling of Fabrication Processes and Constitutive Behaviour of Metals. Too, J. J. Ed. CANMET Research and Technology Centre, 1986 Xu S Endochronic finite element analysis of loads in grain storage bins. M.Sc. Thesis. 1992, University of Manitoba, Winnipeg, Manitoba, Canada