Constitutive Model for Full-Range Elasto-Plastic Behavior of Structural Steels With Yield Plateau Formulation (Hu, Et Al. 2016)

Engineering Structures xxx (2016) xxx–xxx

Contents lists available at ScienceDirect

Engineering Structures journal homepage: www.elsevier.com/locate/engstruct

Constitutive model for full-range elasto-plastic behavior of structural steels with yield plateau: Formulation and implementation Fangxin Hu, Gang Shi ⇑, Yongjiu Shi Key Laboratory of Civil Engineering Safety and Durability of China Education Ministry, Department of Civil Engineering, Tsinghua University, Beijing 100084, PR China

a r t i c l e

i n f o

Article history: Received 15 October 2015 Revised 16 February 2016 Accepted 22 February 2016 Available online xxxx Keywords: Constitutive model Plastic behavior Structural steel Yield plateau Formulation Implementation

a b s t r a c t Performance-based engineering methodologies allow for the design of more reliable seismic resistant structures. Nonetheless, to implement this technique, an accurate constitutive model to predict the elasto-plastic behavior of structural steel components or systems under various loadings is needed to properly evaluate their strength, deformation and energy absorption capacities in case of severe earthquakes. Such a model should also be relatively simple to use for practical purposes in engineering. With these objectives in mind, a new constitutive model is formulated to describe the elasto-plastic behavior of structural steels with yield plateau. This model uses nonlinear kinematic hardening to trace well the significant Bauschinger effect in full-range cyclic loadings, and couples nonlinear isotropic hardening with a memory surface in the plastic strain space to account for the stabilization phenomena of cyclic softening and hardening. An impermanent bounding surface in the stress space is employed to correctly describe the yield plateau response. The consistency condition is investigated in detail, which results in implications that will greatly facilitate the calibration of material dependent parameters. The implementation technique is also presented for three-dimensional and two-dimensional problems respectively. Using the resulting integration algorithms, the proposed constitutive model is successfully incorporated into ABAQUS/Standard by the UMAT subroutine feature. Ó 2016 Elsevier Ltd. All rights reserved.

1. Introduction 1.1. Background and motivation A review of the lessons learned from severe earthquakes that have occurred [1–4] has led to the conclusion that civil engineering structures are exposed to an increasing seismic risk in urban areas and their safety are always challenged due to the intrinsic extreme uncertainties of earthquakes. One of the most effective efforts to reverse this situation in future severe earthquakes is through the development of performance-based seismic design as a more reliable method to resist the effects of earthquakes [5,6]. Implementation of such a design concept results in the well-known multi-level seismic design criteria to achieve different performance levels subjected to different earthquake ground motions, and also gives rise to the capacity design rules allowing specially designed regions within a structure to dissipate energy under earthquake loading, in current seismic standards and code provisions [7–10]. To facilitate a reasonable performance-based design, the detailed and accurate knowledge of deformation and energy absorption capacities of

⇑ Corresponding author.

structural components, assemblages and systems is required in addition to the strength capacity. While experimental methods can be used to determine their seismic performance, numerical methods making use of more and more advanced computational capability for the analysis of structures are becoming more and more popular, especially taking into account the high cost of large-scale tests of structural systems. Regarding the numerical analysis, both the geometrical nonlinearity and material nonlinearity should be well simulated to improve its accuracy. The geometrical nonlinearity has already been well incorporated in many commercial general-purpose finite-element programs; however, the material nonlinearity is not so easily captured, especially under repeated or cyclic loading conditions. It’s well recognized that accurate description of the microscopic material behavior involving cyclic plasticity is essential to predict the macroscopic structural behavior [11,12], and is also the basis before proposing appropriate damage models to characterize the material failure under uncertain loadings [13–16]. Constitutive models for cyclic plasticity of structural steels have developed remarkably in the past decades, and a comprehensive review of the field can be found in the literature of Chaboche [17,18] and Ohno [19]. In the classical single-surface theories of plasticity, it is necessary to give a yield function, a flow rule that

http://dx.doi.org/10.1016/j.engstruct.2016.02.037 0141-0296/Ó 2016 Elsevier Ltd. All rights reserved.

Please cite this article in press as: Hu F et al. Constitutive model for full-range elasto-plastic behavior of structural steels with yield plateau: Formulation and implementation. Eng Struct (2016), http://dx.doi.org/10.1016/j.engstruct.2016.02.037

2

F. Hu et al. / Engineering Structures xxx (2016) xxx–xxx

defines the evolution law of plastic strain and a hardening rule that specifies the post-yielding behavior [20]. The von Mises yield criterion and associated flow rule are commonly accepted for structural steels, while isotropic hardening and/or kinematic hardening can be assumed, which indicate the change in size and translation of the yield surface, respectively. The simplest kinematic hardening is the linear one proposed by Prager [21] and Ziegler [22], which was then generalized to a mixed hardening rule by Hodge [23], but they are usually inaccurate in predicting the actual hardening behavior. Armstrong and Frederick [24] developed a nonlinear kinematic hardening rule, which was then greatly improved by Chaboche and Rousselier [25,26] in that multiple kinematic hardening components evolve independently according to the Armstrong–Frederick rule. This evolution law has been modified again to better describe the ratchetting effects due to its large overestimations [27,28]. The isotropic hardening can be linear or nonlinear too, with respect to the equivalent plastic strain, but based on the extensive experimental results it is generally assumed in a nonlinear form and to saturate eventually, as proposed by Zaverl and Lee [29] and used by Chaboche [30]. Moreover, the saturation of isotropic hardening has been found to depend on the plastic strain amplitude. Thus, Chaboche et al. [31] introduced the concept of memory surface in plastic strain space; Ohno [32] then generalized it to measure the cyclic loading amplitude and applied it successfully in cyclic plasticity models [33]. To describe the cyclic behavior in the yield plateau of mild steels, Rodzik [34] extended the integral description of the back stress evolution model and proposed a new concept of the yield plateau region; Yoshida [35] proposed a model to capture the phenomena of sharp yield point and the subsequent abrupt yield drop; Ucak and Tsopelas [36,37] assumed a threshold controlled by both the equivalent plastic strain and the loading amplitude which is measured by the memory surface, to distinguish between the plateau and hardening regions. The memory surface can also be defined in the stress space, such as by Yoshida and Uemori [38] to better describe the work hardening stagnation effect appearing in the context of sheet metal forming [39]. Jia and Kuwamura [40] then modified the original Yoshida–Uemori model to simulate behavior of structural steels with yield plateau at large cyclic plastic strain. Another modeling approach is the multisurface model proposed independently by Mróz [41] and Iwan [42] based on the underlying idea of the multilayer model by Besseling [43]. The multisurface model utilizes a family of surfaces in stress space to trace the nonlinear plasticity. Due to its complexity in implementation and large storage requirement in numerical analysis, the popular twosurface models were proposed successively by Dafalias and Popov [44,45], Krieg [46] and Tseng and Lee [47], which employ only two surfaces in stress space, i.e., a bounding surface [48–50] (or limit surface) enclosing a yield surface (or loading surface). The plastic modulus is then defined as a function of the distance between the surfaces at the loading point. Petersson and Popov [51,52] extended the Dafalias–Popov model by using a number of intermediate surfaces between the yield and bounding surfaces to remedy the spurious prediction if unloading and reloading are involved without substantial plastic flow, and the cyclic softening, cyclic hardening and mean-relaxation phenomena exhibited by structural steels can be represented [53] as well as the cyclic behavior of structural members [54]. Minagawa et al. [55–57] modified the Petersson–Popov model further by introducing the accumulated effective plastic strain and the effective plastic strain increment as the state variables. To trace the cyclic behavior within the yield plateau, Shen et al. [58–60] developed another modified uniaxial two-surface model and generalized it to a multiaxial one [61–63]; Mahan et al. [64] developed a conceptually different model where a fixed nonhardening bounding surface is assumed to capture the plateau response. Goto et al. [65] refined the

Dafalias–Popov model to formulate a three-surface model, which is characterized by a discontinuity surface inserted between the yield and bounding surfaces, thus representing the transition behavior from the yield plateau to strain hardening region. The three-surface model was later modified also by Goto et al. [66] who introduced a new internal variable to consider the behavior under large plastic strain. In addition, much more complicated behavior including the non-proportional and distortional cyclic hardening has been extensively investigated, such as by Ortiz and Popov [67], McDowell [68,69], Chang and Lee [70,71] and Jiang [72]. Since the fiber-model approach for nonlinear analysis of slender beams and columns is more and more popular, emphasis has been also given to the formulation of one-dimensional uniaxial models, such as those proposed by Cofie and Krawinkler [73], Santhanam [74], Hays [75], Mahan et al. [64] and Wang et al. [76–78]. Although there are still other concepts of constitutive model such as the endochronic theory initially developed by Valanis [79,80] and extensively studied by other researchers [81,82], the classical theory of yield surface based on the superposition of kinematic hardening and isotropic hardening is still popular for modeling material plasticity in structural analysis due to its relatively simple implementation. Thus, it has been incorporated into most commercial finite-element programs, e.g. the so-called Chaboche model in ABAQUS [83]. However, Huang and Mahin [84] have illustrated that although the nonlinearity can be described to some extent, a significant drawback of the Chaboche model is, the calibration may result in widely different sets of kinematic hardening parameters depending upon the strain range of interest. The response in the yield plateau of mild structural steels cannot be obtained either. As for the material parameters, most constitutive models introduced above are not easy to calibrate, and substantial cyclic coupon tests are generally required especially for those multisurface models, e.g., three types of cyclic loading are necessary to calibrate 18 parameters in the two-surface model by Shen et al. [63], or the parameters even have to be determined using the test results of structural members whose cyclic behavior is essentially the one to predict in the three-surface model by Goto et al. [65]. Therefore, from the view point of numerical analysis in engineering design, if we intend to predict the full-range elasto-plastic behavior and stress–strain distributions in steel structures typically under severe earthquakes, a constitutive model that can properly describe the yield plateau phenomena, cyclic hardening and softening behavior and is relatively simple for a straightforward calibration of material dependent parameters, is in urgent need. 1.2. Scope A constitutive model in the framework of classical theory using a yield surface and isotropic and kinematic hardening to describe the nonlinear behavior is explored in this paper. Cyclic softening and hardening and the Bauschinger effect can be described by a carefully devised combination of isotropic and kinematic hardening rules, and a memory surface with a bounding surface is used to correctly represent the plateau behavior and hardening stabilization. It must be mentioned that, several implications for the material dependent parameters can be derived from the formulation, which will greatly facilitate the calibration. The nonproportional loading effects on cyclic hardening is neglected in this work to further simplify the calibration, which is conceived to be reasonable in practical use for structural analysis. In this view, this paper is organized as follows: Initially, Section 1 introduces the background, motivation and scope of this study, whilst the critical observations and assumptions relative to the material cyclic behavior of structural steels with yield plateau are summarized in Section 2. Successively, the complete mathematical formulation of the proposed constitutive model is



illustrated in Section 3. Furthermore, Section 4 presents the implementation algorithm of the proposed model in both threedimensional and two-dimensional numerical problems by employing the UMAT subroutine feature provided in ABAQUS/Standard [83]. Finally, in Section 5, main conclusions are drawn. The detailed calibration procedure of the proposed constitutive model is presented in a companion paper as well as its validation and application at the level of materials, members, connections and structures.

2. Observations on material behavior During monotonic loading, the uniaxial behavior of structural steels with yield plateau is essentially linear elastic from its virgin state up to the sharp point of yield followed by an abrupt yield drop and the yield plateau. The plastic deformation along the yield plateau is caused by the Lüders band propagation [35,85]. From a macroscopic point of view, such a behavior can be treated as perfectly plastic, and the length of the yield plateau can be quantified by the plastic strain at the end of the plateau, as schematically illustrated in Fig. 1. At the end of yield plateau, the material usually starts hardening nonlinearly up to the point corresponding to the ultimate stress when necking initiates in a typical tension coupon test. After ultimate stress is reached, the true stress–strain relationship is almost linear, which indicates the hardening modulus is approximately a constant [86,87]. During cyclic loading, the material behavior of structural steels with yield plateau is much more complex than monotonic loading and is closely related to the strain amplitude of cycling. Generally, cyclic loops converge to a stabilized saturation dependent only on the loading amplitudes. The curve joining the tips of the stabilized cycles as shown in Fig. 2(a) is known as the cyclic backbone curve, which is an important ingredient in uniaxial constitutive models [88,89]. There also exists one cyclic strain amplitude epst : for cyclic amplitudes not greater than epst , the stabilized stress amplitude will not exceed virgin yield stress, which characterizes cyclic softening behavior; otherwise the observed behavior will be cyclic hardening and the stabilized stress amplitude will exceed virgin yield stress [36,53], as schematically illustrated in Fig. 2(b) and (c). The shape of cyclic curves shows significant Bauschinger effect, even for very small strain amplitude [90,91]. This suggests a contraction of elastic range under small strain amplitude and then an expansion under large strain amplitude, which is also verified by the experimental studies [92–94]. In addition, under constant strain amplitude cycling, the stabilization of softening or hardening behavior accompanied by the mean-stress relaxation will be observed [40,53], as shown in Fig. 2(d). Once the strain exceeds current constant amplitude, the stabilized cycle then exhibits softening or hardening again. Based on above observations, a reasonably accurate constitutive model for full-range elasto-plasticity of structural steels with yield

True stress σu

linear

3

plateau must represent the monotonic, small and large strain amplitude cyclic behavior, and must also be capable of predicting the cyclic softening or hardening stabilization and mean-stress relaxation under random cyclic loading. Therefore, several assumptions can be made in order to formulate such a model: (1) There exists a plateau region followed by a hardening region for both monotonic and cyclic loading conditions, and the transition from the plateau region to the hardening region depends on both loading amplitude and the accumulated plastic strain. Such an assumption has also been employed by Ucak and Tsopelas [36]. (2) Within the plateau region, only cyclic softening behavior can take place which characterizes the contraction of elastic range. (3) Within the hardening region, only cyclic hardening behavior can take place which characterizes the expansion of elastic range. (4) No matter in which region, the material is capable of memorizing its current loading amplitude and then determines the status of cyclic softening or hardening behavior. The change of elastic range or radius of yield surface takes place only if the plastic strain state evolves outside current memory. With this rule the stress can stabilize under cyclic loadings after several constant amplitudes. 3. Mathematical formulation According to the observations and assumptions in Section 2, a constitutive model for full-range cyclic behavior of structural steels with yield plateau is proposed here based on the classical theory of rate-independent plasticity. It includes general equations to define the yield surface and flow rule in the stress space, and defines the memory surface in the plastic strain space to memorize the plateau region and loading history. Especially, the kinematic hardening and isotropic hardening equations in plateau and hardening regions are carefully devised to express those characteristics in cyclic plasticity with less material parameters and more consistent internal variables. 3.1. General equations Generally, to describe the metal elasto-plasticity, it is necessary to decomposite the total strain and then defines the elastic and plastic behavior independently. 3.1.1. Strain rate decomposition In the development of stress–strain relationship of the incremental type, time and temperature effects are ignored herein and the total strain rate tensor is assumed to be the sum of the elastic and plastic strain rate components as

e_ ¼ e_ e þ e_ p

ð1Þ

where bold face letters indicate second-order tensors.

nonlinear

σy

εst

εu

True strain plateau hardening region region

Fig. 1. Monotonic behavior of structural steels with yield plateau.

3.1.2. Elastic behavior The elastic behavior is modeled as linear isotropic and the elastic strain rate is defined by the incremental Hooke’s law with two material parameters such as the bulk modulus K and the shear modulus G, i.e.

r_ ¼ ce : e_ e

ð2Þ

e

where c represents the fourth-order elasticity tensor

1 ce ¼ K1 1 þ 2G i 1 1 3

ð3Þ


4


True stress Cyclic backbone curve

(a)

True stress

(b)

Stabilized cycle ε st

p

p

− ε st

True strain

True strain

Plateau region

(c)

True stress

True stress

(d)

σy

p Hardening − ε st region

Plateau region

ε stp

True strain

True strain

Mean-stress at each half cycle

Hardening stabilization

Stabilized cycle

Fig. 2. Cyclic behavior of structural steels with yield plateau: (a) steady-state response; (b) softening under small amplitudes; (c) hardening under large amplitudes and (d) hardening stabilization and mean-stress relaxation under constant amplitudes.

where 1 is the second-order identity tenor, and i is the fourth-order symmetric identity tensor. They can be expressed as

_ Note that substituting Eq. (8) into Eq. (10) results in e_ p ¼ k.

1 ¼ dij ei ej

3.2. Memorization effects

i¼

ð4Þ

1 dik djl þ dil djk ei ej ek el 2

ð5Þ

where dij is Kronecker delta and ei is unit vector. 3.1.3. Plastic behavior The plastic behavior is modeled as hydrostatic pressure independent. The generally accepted von Mises yield criterion is adopted to define the yield surface

R¼ f ¼r

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 ðs aÞ : ðs aÞ R ¼ 0 2

ð6Þ

is the effective von Mises stress and R is the radius of yield where r surface; a is the back stress tensor and s is deviatoric stress tensor which is defined as

1 s ¼ r trðrÞ1 3

ð7Þ

where ‘‘tr” indicates the trace operator and r is the stress tensor. Associated flow rule is assumed for which the plastic strain rate tensor can be expressed as

e_ p ¼ k_

@f 3 ¼ nk_ @r 2

ð8Þ

where k_ is a positive scalar determined using the consistency condition df ¼ 0, and n is the flow direction normal to the yield surface

As illustrated in Section 2, a memory surface is assumed in the plastic strain space to determine the transition from the plateau region to the hardening region, and also to determine the status of isotropic hardening. Such a memory surface has been proposed by Chaboche et al. [31] and then generalized by Ohno [32] as

g¼

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 p ðe nÞ : ðep nÞ r ¼ 0 3

ð11Þ

where n and r is the center and radius of memory surface respectively, whose evolution are defined as

2 n_ ¼ ð1 cÞHðgÞ m : e_ p m 3

ð12Þ

r_ ¼ cHðgÞhn : mie_ p

ð13Þ

where c is a scalar material parameter which determines the rate of expansion of the memory surface; HðgÞ is Heaviside function, i.e. HðgÞ ¼ 1 if g = 0 and HðgÞ ¼ 0 if g < 0; h i denotes the Macaulay brackets, i.e. hai ¼ a if a > 0 and hai ¼ 0 if a < 0; m is the direction pffiffiffiffiffiffiffiffi normal to the memory surface with the norm kmk ¼ 3=2

m¼

ep n r

ð14Þ

Note that for the sake of implementation simplicity, n is not a unit pffiffiffiffiffiffiffiffi vector and its norm is knk ¼ 2=3. The equivalent plastic strain rate is defined as

It is clear that the material parameter c controls the relative rate of the isotropic (increase of radius defined by Eq. (13)) and kinematic components (translation defined by Eq. (12)) in the evolution of memory surface, and it is usually taken as 0 < c 6 1=2 [32]. To distinguish between the plateau and hardening region, the following criterion proposed by Ucak and Tsopelas [36] is assumed

e_ p ¼

if r 6 epst or ep 6 epst ! plateau region if r > epst and ep > epst ! hardening region

n¼

sa r

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 p p e_ : e_ 3

ð9Þ

ð10Þ

ð15Þ


5


where epst is the plastic strain at the end of yield plateau in monotonic stress–strain curve (see Fig. 1), epst is a threshold amplitude of plastic strain (see Fig. 2(b) and (c)) and they are both material dependent parameters. This criterion ensures that the material will not experience any cyclic hardening for fully reversed loading amplitudes smaller than epst ; while for amplitudes larger than epst , the cyclic hardening will not take place until the accumulated plastic strain reaches epst . 3.3. Plateau region According to Mahan et al. [64] and Ucak and Tsopelas [36], in the plateau region a bounding surface can be incorporated into the deviatoric stress space which is identical to the yield surface of the virgin material but is not allowed to translate or change size, as shown in Fig. 3(a). This bounding surface is defined as

l¼

rffiffiffiffiffiffiffiffiffiffiffiffi 3 s : s ry ¼ 0 2

ð16Þ

where ry is the initial yield stress of the virgin material. Such a bounding surface will ensure the equivalent stress in the plateau region cannot exceed ry , which is consistent with the observations in Section 2. Contraction of the yield surface is then described with the following nonlinear isotropic softening equation s R_ ¼ b

ry þ Q s R e_ p

ð17Þ

s

where b and Q s are both short-range material parameters representing the rate of softening and the saturated change (ry < Q s < 0Þ of the size of yield surface, respectively. With initial values of R and ep equal to ry and zero in the virgin state, integration of Eq. (17) yields s

s R ¼ ry þ Q 1 exp b ep

ð18Þ

which indicates that the saturated size of the yield surface is ry þ Q s . Based on the assumptions in Section 2, cyclic softening of the yield surface defined by Eq. (17) or (18) can only take place when the loading is directed outwards of current memory surface shown in Fig. 3(b) and defined in Eqs. (11)–(14) with initial values of zero and short-range scalar material parameter cs , i.e. s R_ ¼ b

ry þ Q s R e_ p when g ¼ 0 and n : m > 0

ð19Þ

R_ ¼ 0; otherwise

To simulate simultaneously the plateau behavior during monotonic loading and softening behavior during cyclic loading in the plateau region, it is imperative to deal with the initial monotonic plastic loading and subsequent unloading or reloading independently.

(a)

3.3.1. Initial loading During the monotonic plastic loading on yield plateau, the bounding surface is treated as the loading surface and consistency is guaranteed by assuming the outward normal direction of yield surface in Eq. (6) and the bounding surface in Eq. (16) to be identical [36], i.e. @f =@ r ¼ @l=@ r, which yields the following evolution of back stress on the yield plateau

a¼ 1

R

ry

s

ð20Þ

It is easy to recognize that during the initial plastic loading, although the yield surface contracts and translates according to Eqs. (19) and (20), it is always in contact with the bounding surface at the loading point. 3.3.2. Unloading and reloading Once unloading–reloading occurs on the yield plateau, the back stress no longer evolves as in Eq. (20). Now, the back stress in Eq. (20) is assumed to be decomposed into several short-range components as proposed by Chaboche and Rousselier [25,26], and for each component, the evolution law initially proposed by Armstrong and Frederick [24] is used to describe the nonlinear kinematic hardening in the form

2 3

a_ sj ¼ C sj e_ p csj asj e_ p

ð21Þ

where C sj and csj are short-range material dependent parameters. For the simulation in engineering application, it is generally recognized that superposition of two short-range back stresses is sufficient to produce acceptable results in describing the cyclic behavior at small plastic strain, i.e.

a¼

2 X

asj

ð22Þ

j¼1

Eq. (21) can be manipulated to obtain

0 1 rffiffiffi rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi! rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s n : a 3 s s 3 3 j C s qffiffiffiffiffiffiffiffiffiffiffiffiffi csj d a : aj ¼ B as : asj C @ Adep 2 j 2 j as : as 2 j j

6

j

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi! 3 s s C sj csj a : aj dep 2 j

ð23Þ

where the inequality n : asj 6 knkasj is used and asj denotes the qffiffiffiffiffiffiffiffiffiffiffiffiffi norm of asj , i.e. asj ¼ asj : asj . Integration of Eq. (23) yields

h i rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s p s 3 s s C j exp cj ðe þ const Þ a : aj 6 s 2 j cj csj

ð24Þ

(b) Strain space

Stress space R s

s 1

O Bounding surface

σy

+

p

s 2

Yield surface

O

ξ

r

Memory surface

Fig. 3. Multiaxial stress and strain space representation in plateau region: (a) yield surface and bounding surface and (b) memory surface.


6


where const represents a scalar constant determined by the initial conditions. With asj and ep denoting the uniaxial back stress and plastic strain, Eqs. (21) and (24) can be rewritten in a given monotonic load path (e.g. a half-cycle) as

a_ sj ¼ C sj e_ p csj asj je_ p j a ¼ s j

C sj

csj

þ a

s j;0

C sj

csj

!

ð26Þ p 0

to Fig. 2(b), the equivalent stress always locates inside the yield plateau upon the stabilization of cyclic softening in the plateau region, which makes it reasonable to assume the following condition adopted by Ucak and Tsopelas [36]

csj

ð32Þ

C sj 3 3 asj : n P s knkasj P 0 c 2 cj 2

C sj s j

ð33Þ

by making use of Eq. (24). On the other hand, the following relation can be obtained by using Eqs. (6) and (9)

3 3 a : ðs aÞ 1 3 3 ¼ s : s a : a R2 a:n¼ 2 2 R 2R 2 2

ð34Þ

With reference to Fig. 3(a), the following triangle inequality holds

rffiffiffiffiffiffiffiffiffiffiffiffiffiffi rffiffiffiffiffiffiffiffiffiffiffiffi 3 3 a:a s : s R 6 2 2

ð35Þ

Substituting Eq. (35) into Eq. (34) yields

¼ Q

s

ð27Þ

Eq. (27) means the saturated absolute value of total back stress equals to the saturated contraction of the yield surface. Thus, it ensures that the yield surface approaches to a limiting surface which is identical to the initial yield surface or bounding surface. The consistency condition of the yield surface leads to

@f @f @f _ f_ ¼ R¼0 : r_ þ : a_ þ @r @a @R

ð28Þ

Assuming the isotropic softening is activated and substituting Eqs. (8), (17) and (21), (22) into Eq. (28) while making use of @f =@ a ¼ @f =@ r and @f =@R ¼ 1, the following relation can be obtained

1 @f k_ ¼ : r_ H @r H¼

!

where the inequality holds since

exp csj ep ep0

where indicates the direction of loading, a and e are initial values of the back stress and plastic strain respectively. Eqs. (24) and (26) are asymptotic functions which approach to C sj =csj regardless for uniaxial and of initial conditions ( asj 6 C sj =csj qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s s s s 3=2aj : aj 6 C j =cj for multiaxial loading). Hence, with reference

j¼1

2

X C sj 3 a:n P min csj j cs 2 j¼1 j

ð25Þ

s j;0

2 X C sj

! X 2 2 X C sj 3 s @f s s s s ¼ C j cj aj : cj s aj : n @r cj 2 j¼1 j¼1

2 X

@f s þ b ry þ Q s R C sj csj asj : @r j¼1

ð29Þ

ð30Þ

where H is hardening modulus. Note that k_ is a positive scalar, which requires the hardening modulus H to be positive, since for plastic loading the term ð@f =@ rÞ : r_ in Eq. (29) is always positive [20]. Thus, the material dependent parameters that appear in the kinematic hardening and isotropic softening have to be carefully calibrated such that the following equation holds for any given load path

2 X

@f s P b R Q s ry C sj csj asj : @ r j¼1

ð31Þ

To satisfy the internal consistency in Eq. (31) regardless of the loading amplitude, Ucak and Tsopelas [36] adopted a crucial assumption s that the material parameter b in Eq. (17) is large enough to achieve close to instantaneous saturation of the size of yield surface in the initial plastic loading, such that the right term in Eq. (26) can be treated as nearly zero. However, such an assumption will definitely overestimate the cyclic softening in plateau region, since experimental results have shown that the yield surface does contract gradually with the equivalent plastic strain within the yield plateau [62,93]. To further elaborate on the restrictions associated with material parameters, substituting Eq. (8) into the term on the left hand side of Eq. (31) leads to

3 a:n6 2

rffiffiffiffiffiffiffiffiffiffiffiffi 3 s : s R 6 ry R 2

ð36Þ

by making use of the bounding condition in Eq. (16). Combining Eqs. (27), (32) and (36) leads to

2

X

@f P min csj R Q s ry C sj csj asj : j @r j¼1

ð37Þ

At this instant, it is obvious that the consistency condition presented in Eq. (31) is always satisfied as long as

s min csj P b j

ð38Þ

which implies that the rate of saturation of the size of yield surface should be less than that of any short-range back stress component. Eqs. (27) and (38) will be used as guidance for calibration of material parameters. Another issue is the determination of initial values of back stress components at the first unloading from the yield plateau. Ucak and Tsopelas [36] have proposed a convenient approach to ensure the internal consistency, which reads

asj ¼

s 1 Cj s s a Q cj

ð39Þ

where a is determined according to Eq. (20). So far, a complete framework for cyclic behavior in the plateau region is elaborated as illustrated in Fig. 3. The material dependent parameters should satisfy Eqs. (27) and (38) such that the consistency condition always holds. Nonetheless, the bounding condition in Eq. (16) could be violated in cycles with very small plastic strain amplitudes, since the evolution of back stress according to Eq. (21) could reach some magnitude larger than the contraction of yield surface before sufficient saturation. Hence, Eq. (16) is monitored at any time and activated with Eq. (20) when violated. Once unloading–reloading from the bounding surface takes place again, the back stress evolution is switched from Eq. (20)–(22) with the initial condition Eq. (39). 3.4. Hardening region Once the criterion in Eq. (15) is satisfied, hardening region is activated. At this instant, the bounding surface defined by Eq. (16) vanishes as shown in Fig. 4(a), and the existing memory


7


(a)

(b)

Stress space

Strain space

R s 1

+

s 2

+

l 1

r

Yield surface

s +

l 2

O

σy

Memory surface Memory is reset

p

Bounding surface vanishes

O

εp0

≤ε

p st

Existing memory is erased

Fig. 4. Multiaxial stress and strain space representation in hardening region: (a) yield surface and vanished bounding surface and (b) reset memory surface.

of the material defined by Eqs. (12)–(14) is erased and reset with the following initial conditions [36]

n ¼ ep0 and r ¼ 0

ð40Þ

where e is the plastic strain at the onset of hardening region, as shown in Fig. 4(b). The evolution law of the memory surface is defined in the same Eqs. (12)–(14) but with initial conditions in Eq. (40) and long-range scalar material parameter cl . Cyclic hardening behavior is described with the following nonlinear isotropic hardening equation

memory surface and cyclic hardening instead of softening. The consistency condition can be given as well, and it is easy to find that the hardening modulus is definitely positive due to progressive hardening of the yield surface.

p 0

l R_ ¼ b

ry þ Q s þ Q l R e_ p

ð41Þ

l

where b and Q l are both long-range positive material parameters representing the rate of hardening and the saturated change of the size of yield surface, respectively. For convenience, the size of yield surface is assumed to have been saturated from contraction when hardening region is activated, which indicates the parameter s b controlling the rate of saturation should be large enough to achieve close to full saturation for an equivalent plastic strain at least being epst . Thus, with initial values of R and ep equal to ry þ Q s and ep0 (ep0 P epst Þ at the onset of hardening region, integration of Eq. (41) yields

h

i l l R ¼ ry þ Q s þ Q l 1 exp b ep0 b ep

ð42Þ

which indicates that the saturated size of the yield surface is

ry þ Q s þ Q l . Based on the assumptions in Section 2, cyclic hardening of the yield surface defined by Eq. (41) or (42) can only take place when the loading is directed outwards of current memory surface to simulate the stabilization of hardening, i.e. l R_ ¼ b

ry þ Q s þ Q l R e_ p when g ¼ 0 and n : m > 0

R_ ¼ 0; otherwise

ð43Þ

In the hardening region, two long-range back stress components are activated to simulate better cyclic behavior at large plastic strain and added to the existing short-range components, while they also evolve in the Armstrong–Frederick hardening rule [24], i.e., Eqs. (21) and (22) are now modified as

2 3

a_ lj ¼ C lj e_ p clj alj e_ p a¼

2 X

2 X

j¼1

j¼1

asj þ

alj

4. Numerical implementation Although the mathematical equations formulated in Section 3 are generally given in rate form, in numerical analysis the state variables of a material constitutive model such as stress, elastic and plastic strains are updated from one discrete time to the next, known as ‘‘time integration”. The integration procedure is conceived as a strain-driven algorithm, where the object is to calculate the current values of state variables based on the preceding values for a given strain increment. It is also necessary to compute the algorithmic consistent tangent modulus which is used to carry out the iterative solution of the global equilibrium equations [95]. Therefore, the corresponding discretization and update procedure of these state variables and the consistent tangent modulus, i.e., implementation of the proposed constitutive model is presented in this section. The general procedure for implementation in three-dimensional (or plane-strain) problems is described first, followed by the special treatment in two-dimensional (e.g. planestress) cases.

4.1. Three-dimensional problem For a three-dimensional problem, all six components of the state variables such as stress and strain tensors are evaluated. First, to discretize the state variables, an implicit backward Euler difference scheme is employed. Next, the state variables are updated using a two-step algorithm: an elastic trial predictor followed by a plastic corrector that performs projection of the trial states onto the yield surface, as illustrated by Fig. 5. The elastic predictor

trial

σ i +1

ð44Þ

ð45Þ

where C lj and clj are long-range material dependent parameters. So far, a concise framework for cyclic behavior in the hardening region is also established as illustrated in Fig. 4. The basic idea behind is similar to that in the plateau region, except for the reset

σ i +1

∂Eσ

Eσ (elastic domain) σi Fig. 5. Geometric illustration of closest point projection for implementation.


8


obtains the trial elastic states by freezing plastic flow during the increment (i.e. Dk ¼ 0Þ as

ep;trial ¼ epi iþ1

ð46Þ

e rtrial iþ1 ¼ ri þ c : De

ð47Þ

strial iþ1

¼r

trial iþ1

1 tr rtrial iþ1 1 3

ð48Þ

ee;trial ¼ eei þ De iþ1

ð49Þ

e

ð50Þ

p;trial iþ1

¼e

p i

De ¼ 0

Thus, upon exploiting kniþ1 k ¼ as

Dk ¼

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 trial s : strial iþ1 ry 2 iþ1

Once Dk is determined, the state variables can be updated accordingly. To update the memory surface, the following discretization is given

niþ1 ¼ ni þ ð1 cs ÞCDkmiþ1

ð63Þ

riþ1 ¼ r i þ cs CDk

ð64Þ

Then the plastic return-mapping corrector updates the states as

miþ1 ¼

epiþ1 ¼ ep;trial þ Dk iþ1

where

riþ1 ¼ r

trial iþ1

ð52Þ

c : De e

ð53Þ

p

e

De

p

¼e

ð55Þ

þ De

p

3 ¼ niþ1 Dk 2

Riþ1 siþ1 ¼ 1

ry

ð67Þ

r i þ C Dk

pffiffiffiffiffiffiffiffi 3=2, another scalar equation is

ð68Þ

which indicates that the memory surface should be updated if the updated plastic strain locates outside it. To update the size of yield surface, Eq. (19) is discretized as

Riþ1 ¼ Ri þ b

s

ry þ Q s Riþ1 CDk

ð69Þ

which gives

Riþ1 ¼

Ri þ

ry þ Q s bs CDk

ð70Þ

s

1 þ b C Dk

The consistent tangent modulus necessary for a rapid convergence of the solution of global equilibrium equations requires the derivative of the stress riþ1 with respect to the total strain eiþ1 . According to Eqs. (47), (48), (53), (57), (61) and (62), the following partial derivatives can be obtained @strial

iþ1 strial strial : 2G i 13 1 1 @ Dk 1 iþ1 : @ eiþ1 ffi ¼ niþ1 : i 1 1 ¼ niþ1 ffi ¼ iþ1 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi @ eiþ1 2G 3 strial : strial 3 3 trial 2G s : strial 2 iþ1

iþ1

2 iþ1

iþ1

ð71Þ

3Gniþ1 Dk

siþ1 aiþ1 ¼ Riþ1 niþ1

ð59Þ

@niþ1 ¼ @ eiþ1 ¼

ð60Þ

Combining Eqs. (59) and (60) leads to

strial iþ1 ry þ 3GDk

epiþ1 ni

ð57Þ

The yield function gives

niþ1 ¼

miþ1 ¼

*rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi +

p

2 p CDk ¼ e ni : eiþ1 ni ri 3 iþ1

resulting in

siþ1 aiþ1 ¼

ð66Þ

ð58Þ

ry

strial iþ1

ð56Þ

4.1.1. Plateau region, on the bounding surface In the plateau region, when plastic flow occurs on the bounding surface such as in the initial monotonic loading or in occasions for which the bounding condition Eq. (16) is violated in subsequent unloading and reloading, the back stress in Eq. (20) is discretized as

Riþ1

C ¼ H g iþ1 hniþ1 : miþ1 i

and making use of kmiþ1 k ¼ obtained as

where subscript i represents the increment number, Dk is equivalent plastic strain increment to be determined, and De and Dep are strain increment and plastic strain increment, respectively. The plastic corrector can be implemented based on the theory of radial return algorithm originally proposed by Wilkins [96], which was then generalized as the closest point projection method by Simo and Hughes [97]. A similar fully implicit integration technique for nonlinear isotropic and kinematic hardening plasticity has been proposed by Doghri [98], while Mahnken [99], Hartmann and Haupt [100] have shown that such a computation problem can be reduced to a one-dimensional scalar equation, in which only the plastic multiplier (i.e. Dk) appears as the unknown. The development of integration algorithms for the proposed constitutive model in Section 3 is shown below.

aiþ1

ð65Þ

r iþ1

ð54Þ

p eeiþ1 ¼ ee;trial iþ1 De p;trial iþ1

epiþ1 niþ1

Substituting Eqs. (63) and (64) into Eq. (65) yields

p siþ1 ¼ strial iþ1 2GDe

p iþ1

ð62Þ

3G

ð51Þ

p

pffiffiffiffiffiffiffiffi 2=3, a scalar equation is obtained

ð61Þ

@strial iþ1 @ eiþ1

@ Dk ry þ 3GDk 3Gstrial iþ1 @ eiþ1

ry þ 3GDk

2

2G i 13 1 1 3Gniþ1 niþ1 ry þ 3GDk

@ riþ1 @ Dk @niþ1 ¼ ce 3Gniþ1 3GDk @ eiþ1 @ eiþ1 @ eiþ1

ð72Þ

ð73Þ

Thus, substituting Eqs. (3), (71) and (72) into Eq. (73), the consistent tangent modulus is obtained as


9


cep iþ1 ¼

@ riþ1 @ eiþ1

¼ K1 1 þ

The state variables are updated after convergence of Eqs. (79)–(84). Because eiþ1 is fixed in the plastic corrector

1 2G i 1 1 3 ry þ 3GDk

ry

0¼

ry 3Gniþ1 niþ1 ry þ 3GDk

ð74Þ

4.1.2. Plateau region, inside the bounding surface In the plateau region, when plastic flow occurs inside the bounding surface in subsequent unloading and reloading, the back stress should be now discretized using Eq. (21) as

asj

iþ1

¼ asj þ C sj niþ1 csj asj i

resulting in

"

siþ1 aiþ1 ¼ fiþ1 3G þ

2 X j¼1

where

fiþ1 ¼

strial iþ1

Dk

iþ1

ð75Þ

# C sj Dkniþ1 1 þ csj Dk

ð76Þ

asj i s 1 þ c D k j j¼1 2 X

ð77Þ

The yield function gives Eq. (60). Hence, the direction niþ1 is determined exclusively in termsrof ffiffiffi fiþ1

niþ1

f ~ iþ1 ¼ ¼ iþ1 where r r~ iþ1

3 kf k 2 iþ1

ð78Þ

Therefore, combining Eqs. (60), (76)–(78) renders the following nonlinear scalar equation " #

~ iþ1 3G þ 0 ¼ hðDkÞ ¼ r

2 X j¼1

C sj Dk Riþ1 1 þ csj Dk

ð79Þ

¼

1 @h ðkÞ Dk h DkðkÞ @ Dk

@ Dk 3G @h niþ1 where k ¼ ¼ ðDkÞ @ eiþ1 k @ Dk

@h ðkÞ 3 ðkÞ Dk ¼ niþ1 @ Dk 2 2 3

2 2 csj asj X C sj 6X 7 i :4

2 5 3G

2 j¼1 1 þ cs DkðkÞ j¼1 1 þ cs DkðkÞ j j @Riþ1 ðkÞ Dk @ Dk

ð81Þ

ð82Þ

8 when CDkðkÞ ¼ 0 @ CDk ðkÞ < 0 @n

Dk ¼ DkðkÞ þ DkðkÞ @Diþ1k ðDkðkÞ ÞþnðkÞ :ðep ni Þ i iþ1 : @ Dk when CDkðkÞ > 0 CDkðkÞ þr i

ð83Þ where according to Eqs. (77) and (78)

P2 @niþ1 ðkÞ Dk ¼ @ Dk

j¼1

csj asj

i

1þcsj DkðkÞ

2

3

6 ðkÞ P2 2 432 niþ1 : j¼1

csj asj

i

1þcsj DkðkÞ

ð86Þ

@niþ1 @niþ1 @fiþ1 i 32 niþ1 niþ1 ¼ : ¼ @ eiþ1 @fiþ1 @ eiþ1 r~ iþ1 2 3

s s X 2 c a j j 1 @ D k 6 7 i : 42G i 1 1 þ 5

2 3 @ eiþ1 j¼1 1 þ cs Dk j

ð87Þ

Thus, substituting Eqs. (3), (86) and (87) into Eq. (73), the consistent tangent modulus is obtained as

cep iþ1 ¼

@ riþ1 1 ¼ K1 1 þ 2Ghiþ1 i 1 1 3 @ eiþ1 ~ c iþ1 3Ghiþ1 niþ1 niþ1 e

ð88Þ

where

hiþ1 ¼ 1

3GDk r~ iþ1

ð89Þ

~hiþ1 ¼ 3G ð1 hiþ1 Þ k

ð90Þ

2 3

s s X 2 c a j j 9G Dk 3 6 7 i ¼ i niþ1 niþ1 : 4

2 niþ1 5 ~ iþ1 kr 2 s j¼1 1 þ c Dk j 2

ð91Þ

4.1.3. Hardening region In the hardening region, the back stress definition is switched from Eqs. (20)–(22), (44), (45). A similar derivation to Section 4.1.2 can be carried out. The nonlinear scalar equation to solve Dk is obtained as

where the index k refers to the iteration steps. The update procedure for memory surface and the size of yield surface is the same with that in Section 4.1.1. Thus, from Eqs. (50), (56), (57), (68) and (70)

@Riþ1 ðkÞ ry þ Q s Ri @ CDk ðkÞ s Dk Dk ¼b 2 @ Dk @ Dk s 1 þ b CDkðkÞ

ð85Þ

and

e c iþ1 ð80Þ

@h @ Dk ðDkÞ þ 3Gniþ1 @ Dk @ eiþ1

So

which is solved with the following Newton method

Dkðkþ1Þ ¼ DkðkÞ

@s trial @h @h @ Dk @h ¼ ðDkÞ þ : iþ1 trial @ eiþ1 @ Dk @ eiþ1 @siþ1 @ eiþ1

7 ðkÞ 2 5niþ1

0 ¼ hðDkÞ

"

~ iþ1 3G þ ¼r

2 X j¼1

!# C sj C lj Dk Riþ1 þ 1 þ csj Dk 1 þ clj Dk

ð92Þ

where

rffiffiffi 3 kf k and fiþ1 2 iþ1

1 0 2 asj alj X i i A @ ¼ strial þ iþ1 s 1 þ cj Dk 1 þ clj Dk j¼1

r~ iþ1 ¼

ð93Þ

To update the memory surface and the size of yield surface, Eqs. (63), (64) and (70) are now modified as

niþ1 ¼ ni þ 1 cl CDkmiþ1

ð94Þ

riþ1 ¼ ri þ cl CDk

ð95Þ

r~ ðkÞ iþ1 ð84Þ

Riþ1 ¼

l Ri þ ry þ Q s þ Q l b CDk l

1 þ b C Dk

ð96Þ


10


Eq. (92) is solved with Newton method as shown in Eq. (80), but now

@h ðkÞ 3 ðkÞ ¼ niþ1 Dk @ Dk 2 0 1

s s l l 2 c a c a X j j j j B C i i : @ 2 þ 2 A 3G ðkÞ ðkÞ s l j¼1 1 þ cj Dk 1 þ cj Dk 0 1 s l 2 X C C B C j j @ 2 þ 2 A ðkÞ ðkÞ s l j¼1 1 þ cj Dk 1 þ cj Dk

@Riþ1 DkðkÞ @ Dk

Beginning of Analysis Define Initial Conditions Start of Step Start of Increment Calculate Integration Point Field Variables from Nodal Values Start of Iteration Calculate Δε

ð97Þ Call UMAT

where

NDI=3

@Riþ1 DkðkÞ @ Dk

l r y þ Q þ Q R i @ C Dk DkðkÞ ¼b 2 @ Dk l 1 þ b CDkðkÞ s

l

ð98Þ

P2 B j¼1 @

csj asj

i

2 þ

1þcsj DkðkÞ

@niþ1 DkðkÞ ¼ @ Dk

2

clj alj

i

1þclj DkðkÞ

r~ ðkÞ iþ1 0

p

csj asj

i

l=

2 þ

clj alj

i

1þclj DkðkÞ

13 C7 ðkÞ 2 A5niþ1

ðkÞ r~ iþ1

3 2

No

Inside the Bounding Surface

On the Bounding Surface

Update σ, SDVs, ∂σ/∂ε

After convergence, the consistent tangent modulus is obtained as Eq. (88) but in which

9G2 Dk 3 e i niþ1 niþ1 c iþ1 ¼ ~ iþ1 kr 2 2 0 3

1 s s l l 2 c a c a X j j j j 6 B C 7 i i : 4 @ 2 þ 2 A niþ1 5 s l j¼1 1 þ cj Dk 1 þ cj Dk

s: s − σ y < 0

Yes

ð99Þ

Obtain Element Stiffness Matrix [Kel] Definition of System Stiffness Matrix [K] Definition of Loads {F} Solve [K]{U}={F}

ð100Þ

Converged? Yes Write Output

4.2. Two-dimensional problem No

For a two-dimensional problem such as the plane-stress case, the normal and transverse shear stresses are assumed to be zero, i.e., r3i 0 for i ¼ 1; 2; 3. This form of the equations make a simple radial return algorithm no longer applicable, but as proposed by Simo and Taylor [101], it’s viable to perform it directly in the constrained plane-stress subspace, where the general equations in tensor form in Section 3 can be rephrased in vector form. Denote the vector space of symmetric second order tensors by S, then the plane-stress subspace is denoted by

ð101Þ

Similarly, the subspace of deviatoric symmetric second-order tensors is denoted by

SD :¼ fs 2 Sjs13 ¼ s23 ¼ 0; tr½s ¼ skk ¼ 0g

ð102Þ

Thus, Eqs. (1), (2), (6) and (8) are rewritten as

e ¼ ee þ ep

ð103Þ

r ¼ Ce

ð104Þ

e

No

Check the Bounding Condition

SP :¼ fr 2 Sjr13 ¼ r23 ¼ r33 ¼ 0g

p

Hardening Region

Plateau Region

C 2 A

1þcsj DkðkÞ

p

Yes

1

63 ðkÞ P2 B 42 niþ1 : j¼1 @

Plane Stress Elements

r ≤ ε st or ε ≤ ε st ?

and @ CDk DkðkÞ =@ Dk is defined in Eq. (83) but in which

0

NDI=2

3D Solid Elements

End of Step?

No

Yes

Fig. 6. Flow chart of UMAT in ABAQUS/Standard.

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 T R¼ f ¼r v Pv R ¼ 0 where v ¼ r b 2

e_ p ¼ k_

@f 3 Pv _ k ¼ @r 2 r

ð105Þ

ð106Þ

where r ¼ ½r11 r22 r12 T , e ¼ ½e11 e22 e12 T , b ¼ ½b11 b22 b12 T is the back stress vector hypothetically defined in SP , and the elastic stiffness matrix C and projection matrix P are defined as

2 1 m E 6 m 1 C¼ 4 1 m2 0 0

2 3 2 1 0 1 6 7 0 7 5 and P ¼ 4 1 2 0 5 3 1m 0 0 6 2 0

3

ð107Þ

where E is elastic modulus and m is the Poisson ratio. The back stress evolution in Eq. (20), (21) and (44) are now rewritten as, while making use of Eq. (106)



R b¼ 1 r

ry

b_ sj ¼

! C sj s s _ v cj bj k R

b_ lj ¼

! C lj v clj blj k_ R

ð108Þ

ð109Þ

ð110Þ

Following a similar integration algorithm to that in Section 4.1, the update procedure for the state variables and consistent tangent matrix can be established as well. For the sake of simplicity, the detailed development of these algorithms is omitted here. 4.3. Implementation using UMAT subroutine in ABAQUS/Standard The proposed constitutive model is implemented in the general-purpose finite-element package program ABAQUS/Standard [83] by the UMAT subroutine feature. UMAT is called for each material point at each iteration of every increment. It is provided with the material state at the start of the increment (stresses and solution-dependent state variables, i.e. SDVs including elastic and plastic strains, back stresses, the center and radius of memory surface, equivalent plastic strain and a predefined variable to indicate whether the material is in plateau region or hardening region) and with the strain increment, the stresses and those solutiondependent state variables will be updated at the end of the increment. In addition, UMAT also provides the material Jacobian matrix, i.e., the consistent tangent modulus. The flow chart of UMAT in ABAQUS/Standard is shown in Fig. 6. 5. Conclusions A new constitutive model is developed and aimed at providing an easy to calibrate, yet complete and rigorous, approach to describe the structural steel response with yield plateau under various loading conditions, monotonic or cyclic, periodic or random. Such an ability to address all these various loadings is based on the fundamental principles of classical theory of plasticity. The model consists of a yield surface with isotropic and kinematic hardening, a memory surface to measure the dependency of response on the loading amplitude, and a bounding surface specifically to deal with the behavior within yield plateau. While there are many models for plasticity of structural steels, what distinguishes the present model is the particular handling of restrictions on material dependent parameters from the consistency condition, and correct description of cyclic softening and hardening behavior using unified kinematic hardening and isotropic hardening coupled with the memorization effects. This essential characteristic is necessary for an easy implementation and application for practical purposes in structural engineering. Acknowledgement The authors gratefully acknowledge the support for this work, which was sponsored by the National Natural Science Foundation of China (No. 51478244), and the Excellent Young Scientist Programme of the National Natural Science Foundation of China (No. 51522806). References [1] Miller DK. Lessons learned from the Northridge earthquake. Eng Struct 1998;20(4–6):249–60. [2] Mahin SA. Lessons from damage to steel buildings during the Northridge earthquake. Eng Struct 1998;20(4–6):261–70.

11

[3] Nakashima M, Inoue K, Tada M. Classification of damage to steel buildings observed in the 1995 Hyogoken-Nanbu earthquake. Eng Struct 1998;20(4– 6):271–81. [4] Watanabe E, Sugiura K, Nagata K, Kitane Y. Performances and damages to steel structures during the 1995 Hyogoken-Nanbu earthquake. Eng Struct 1998;20(4–6):282–90. [5] Bertero RD, Bertero VV. Performance-based seismic engineering: the need for a reliable conceptual comprehensive approach. Earthquake Eng Struct Dynam 2002;31(3):627–52. [6] Krawinkler H, Zareian F, Medina RA, Ibarra LF. Decision support for conceptual performance-based design. Earthquake Eng Struct Dynam 2006;35(1):115–33. [7] ANSI/AISC 341-10. Seismic provisions for structural steel buildings. Chicago (IL, US): American Institute of Steel Construction (AISC); 2010. [8] Eurocode 8 (EN1998-1-2004). Design of structures for earthquake resistance – Part 1: General rules, seismic actions and rules for buildings. Brussels (Belgium): European Committee for Standardization (CEN); 2004. [9] AIJ. Guidelines for limit state design of steel structures. Tokyo (Japan): Architectural Institute of Japan; 1998 [in Japanese]. [10] GB 50011-2010. Code for seismic design of buildings. Beijing (China): China Architecture & Building Press; 2010 [in Chinese]. [11] Franco JM, Cahís X, Gracia L, López F. Experimental testing of a new antiseismic dissipator energy device based on the plasticity of metals. Eng Struct 2010;32(9):2672–82. [12] Al-Bermani FGA, Zhu K. Nonlinear elastoplastic analysis of spatial structures under dynamic loading using kinematic hardening models. Eng Struct 1996;18(8):568–76. [13] Ayhan B, Jehel P, Brancherie D, Ibrahimbegovic A. Coupled damage–plasticity model for cyclic loading: theoretical formulation and numerical implementation. Eng Struct 2013;50:30–42. [14] Kim D, Dargush GF, Basaran C. A cyclic two-surface thermoplastic damage model with application to metallic plate dampers. Eng Struct 2013;52:608–20. [15] Mendes LAM, Castro LMSS. A simplified reinforcing steel model suitable for cyclic loading including ultra-low-cycle fatigue effects. Eng Struct 2014;68:155–64. [16] Pereira JCR, de Jesus AMP, Xavier J, Fernandes AA. Ultra low-cycle fatigue behaviour of a structural steel. Eng Struct 2014;60:214–22. [17] Chaboche JL. A review of some plasticity and viscoplasticity constitutive theories. Int J Plast 2008;24(10):1642–93. [18] Chaboche JL. Constitutive equations for cyclic plasticity and cyclic viscoplasticity. Int J Plast 1989;5(3):247–302. [19] Ohno N. Recent topics in constitutive modeling of cyclic plasticity and viscoplasticity. Appl Mech Rev 1990;43(11):283–95. [20] Chen WF, Saleeb AF. Elasticity and plasticity. Beijing (China): China Architecture & Building Press; 2005. [21] Prager W. A new method of analyzing stresses and strains in work-hardening plastic solids. J Appl Mech 1956;23:493–6. [22] Ziegler H. A modification of Prager’s hardening rule. Quart Appl Mech 1959;17:55–65. [23] Hodge PG. Discussion of Prager (1956). J Appl Mech 1957;24:482–4. [24] Armstrong PJ, Frederick CO. A mathematical representation of the multiaxial Bauschinger effect. CEGB Report RD/B/N731. Berkeley (UK): Central Electricity Generating Board; 1966. [25] Chaboche JL, Rousselier G. On the plastic and viscoplastic constitutive equations – Part I: Rules developed with internal variable concept. J Pressure Vessel Technol 1983;105(2):153–8. [26] Chaboche JL, Rousselier G. On the plastic and viscoplastic constitutive equations – Part II: Application of internal variable concepts to the 316 stainless steel. J Pressure Vessel Technol 1983;105(2):159–64. [27] Chaboche JL. On some modifications of kinematic hardening to improve the description of ratchetting effects. Int J Plast 1991;7(7):661–78. [28] Dafalias YF, Kourousis KI, Saridis GJ. Multiplicative AF kinematic hardening in plasticity. Int J Solids Struct 2008;45(10):2861–80. [29] Zaverl J, Lee D. Constitutive relations for nuclear reactor core materials. J Nucl Mater 1978;75(1):14–9. [30] Chaboche JL. Time-independent constitutive theories for cyclic plasticity. Int J Plast 1986;2(2):149–88. [31] Chaboche JL, Dang-Van K, Cordier G. Modelization of the strain memory effect on the cyclic hardening of 316 stainless steel. In: Proceedings of international conference on structural mechanics in reactor technology (SMIRT5). Amsterdam, Netherlands: North-Holland Publishing Co.; 1979. 9–21 August, Paper No. L11/3. [32] Ohno N. A constitutive model of cyclic plasticity with a nonhardening strain region. J Appl Mech 1982;49(4):721–7. [33] Ohno N, Kachi Y. A constitutive model of cyclic plasticity for nonlinear hardening materials. J Appl Mech 1986;53(2):395–403. [34] Rodzik P. Cyclic hardening rule for structural steels with yield plateau. J Eng Mech 1999;125(12):1331–43. [35] Yoshida F. A constitutive model of cyclic plasticity. Int J Plast 2000;16(3– 4):359–80. [36] Ucak A, Tsopelas P. Constitutive model for cyclic response of structural steels with yield plateau. J Struct Eng (ASCE) 2011;137(2):195–206. [37] Ucak A, Tsopelas P. Accurate modeling of the cyclic response of structural components constructed of steel with yield plateau. Eng Struct 2012;35:272–80.


12


[38] Yoshida F, Uemori T. A model of large-strain cyclic plasticity describing the Bauschinger effect and workhardening stagnation. Int J Plast 2002;18(56):661–86. [39] Yoshida F, Uemori T, Fujiwara K. Elastic-plastic behavior of steel sheets under in-plane cyclic tension-compression at large strain. Int J Plast 2002;18(5– 6):633–59. [40] Jia LJ, Kuwamura H. Prediction of cyclic behaviors of mild steel at large plastic strain using coupon test results. J Struct Eng (ASCE) 2014;140(2):441–54. [41] Mróz Z. On the description of anisotropic workhardening. J Mech Phys Solids 1967;15(3):163–75. [42] Iwan WD. On a class of models for the yielding behavior of continuous and composite systems. J Appl Mech 1967;34(3):612–7. [43] Besseling JF. A theory of elastic, plastic and creep deformations of an initially isotropic material showing anisotropic strain-hardening, creep recovery and secondary creep. J Appl Mech 1958;80:529. [44] Dafalias YF, Popov EP. A model of nonlinearly hardening materials for complex loading. Acta Mech 1975;21(3):173–92. [45] Dafalias YF, Popov EP. Plastic internal variables formalism of cyclic plasticity. J Appl Mech 1976;43(4):645–51. [46] Krieg RD. A practical two surface plasticity theory. J Appl Mech 1975;42 (3):641–6. [47] Tseng NT, Lee GC. Simple plasticity model of two-surface type. J Eng Mech 1983;109(3):795–810. [48] Dafalias YF. Bounding surface plasticity. I: Mathematical foundation and hypoplasticity. J Eng Mech 1986;112(9):966–87. [49] Dafalias YF, Herrmann LR. Bounding surface plasticity. II: Application to isotropic cohesive soils. J Eng Mech 1986;112(12):1263–91. [50] Anandarajah A, Dafalias YF. Bounding surface plasticity. III: Application to anisotropic cohesive soils. J Eng Mech 1986;112(12):1292–318. [51] Petersson H, Popov EP. Constitutive relations for generalized loadings. J Eng Mech Div 1977;103(EM4):611–27. [52] Popov EP, Petersson H. Cyclic metal plasticity: experiments and theory. J Eng Mech Div 1978;104(EM6):1371–88. [53] Elnashai AS, Izzuddin BA. Modeling of material non-linearities in steel structures subjected to transient dynamic loading. Earthquake Eng Struct Dynam 1993;22(6):509–32. [54] Mizuno E, Kato M, Fukumoto Y. Multi-surface model application to beamcolumns subjected to cyclic loads. J Constr Steel Res 1987;7(4):253–77. [55] Minagawa M, Nishiwaki T, Masuda N. Prediction of hysteretic stress–strain relations of tension-compression steel members by a multi surface plasticity model. J Struct Eng (JSCE) 1986;32A:193–206 [in Japanese]. [56] Minagawa M, Nishiwaki T, Masuda N. Modeling cyclic plasticity of structural steels. Struct Eng/Earthquake Eng 1987;4(2):361s–70s. [57] Minagawa M, Nishiwaki T, Masuda N. Prediction of uni-axial cyclic plasticity behaviours of structural steel in plastic flow region. J Struct Eng (JSCE) 1989;35A:53–65 [in Japanese]. [58] Shen C, Tanaka Y, Mizuno E, Usami T. A two-surface model for steels with yield plateau. Struct Eng/Earthquake Eng 1992;8(4):179s–88s. [59] Tanaka Y, Mizuno E, Shen C, Usami T. A cyclic two-surface model for steel with yield plateau. J Struct Eng (JSCE) 1991;37A:1–14 [in Japanese]. [60] Shen C, Mizuno E, Usami T. Further study on two-surface model for structural steels under uniaxial cyclic loading. Struct Eng/Earthquake Eng 1993;9 (4):257s–60s. [61] Shen C, Mizuno E, Usami T. A generalized two-surface model for structural steels under cyclic loading. Struct Eng/Earthquake Eng 1993;10(2):59s–69s. [62] Mamaghani IHP, Shen C, Mizuno E, Usami T. Cyclic behavior of structural steels. I: Experiments. J Eng Mech 1995;121(11):1158–64. [63] Shen C, Mamaghani IHP, Mizuno E, Usami T. Cyclic behavior of structural steels. II: Theory. J Eng Mech 1995;121(11):1165–72. [64] Mahan M, Dafalias Y, Taiebat M, Heo YA, Kunnath SK. SANISTEEL: simple anisotropic steel plasticity model. J Struct Eng (ASCE) 2011;137(2):185–94. [65] Goto Y, Wang Q, Obata M. FEM analysis for hysteretic behavior of thin-walled columns. J Struct Eng (ASCE) 1998;124(11):1290–301. [66] Goto Y, Jiang K, Obata M. Stability and ductility of thin-walled circular steel columns under cyclic bidirectional loading. J Struct Eng (ASCE) 2006;132 (10):1621–31. [67] Ortiz M, Popov EP. Distortional hardening rules for metal plasticity. J Eng Mech 1983;109(4):1042–57. [68] McDowell DL. A two surface model for transient nonproportional cyclic plasticity. Part 1: Development of appropriate equations. J Appl Mech 1985;52(2):298–302. [69] McDowell DL. A two surface model for transient nonproportional cyclic plasticity. Part 2: Comparison of theory with experiments. J Appl Mech 1985;52(2):303–8. [70] Chang KC, Lee GC. Biaxial properties of structural steel under nonproportional loading. J Eng Mech 1986;112(8):792–805. [71] Chang KC, Lee GC. Constitutive relations of structural steel under nonproportional loading. J Eng Mech 1986;112(8):806–20.

[72] Jiang W. Study of two-surface plasticity theory. J Eng Mech 1994;120 (10):2179–200. [73] Cofie NG, Krawinkler H. Uniaxial cyclic stress-strain behavior of structural steel. J Eng Mech 1985;111(9):1105–20. [74] Santhanam TK. Model for mild steel in inelastic frame analysis. J Struct Div 1979;105(ST1):199–220. [75] Hays CO. Inelastic material models in earthquake response. J Struct Div 1981;107(ST1):13–27. [76] Shi YJ, Wang M, Wang YQ. Experimental and constitutive model study of structural steel under cyclic loading. J Constr Steel Res 2011;67:1185–97. [77] Wang M, Shi YJ, Wang YQ. Equivalent constitutive model of steel with cumulative degradation and damage. J Constr Steel Res 2012;79:101–14. [78] Shi G, Wang M, Bai Y, Wang F, Shi YJ, Wang YQ. Experimental and modeling study of high-strength structural steel under cyclic loading. Eng Struct 2012;37(4):1–13. [79] Valanis KC. A theory of viscoplasticity without a yield surface. Part I: General theory. Arch Mech 1971;23(4):517–33. [80] Valanis KC. A theory of viscoplasticity without a yield surface. Part II: Application to mechanical behavior of metals. Arch Mech 1971;23 (4):535–51. [81] Sugiura K, Lee GC, Chang KC. Endochronic theory for structural steel under nonproportional loading. J Eng Mech 1987;113(12):1901–17. [82] Hsu SY, Jain SK, Griffin OH. Verification of endochronic theory for nonproportional loading paths. J Eng Mech 1991;117(1):110–31. [83] ABAQUS. ABAQUS Analysis User’s Guide (Version 6.13). Dassault Systèmes Simulia Corp.: Providence, RI, U.S.; 2013. [84] Huang YL, Mahin SA. Simulating the inelastic seismic behavior of steel braced frames including the effects of low-cycle fatigue. PEER Report 2010/104. Berkeley (CA, U.S.): Pacific Earthquake Engineering Research Center. University of California; 2010. [85] Hall EO. Yield point phenomena in metals and alloys. New York (U. S.): Plenum Press; 1970. [86] Bridgman PW. Studies in large plastic flow and fracture: with special emphasis on the effects of hydrostatic pressure. New York (U.S.): McGrawHill Book Company, Inc.; 1952. [87] Jia LJ, Kuwamura H. Ductile fracture simulation of structural steels under monotonic tension. J Struct Eng (ASCE) 2014;140(5):472–82. [88] Yamada S, Imaeda T, Okada K. Simple hysteresis model of structural steel considering the Bauschinger effect. J Struct Constr Eng 2002(559):225–32 [in Japanese]. [89] Matsumoto Y, Chung K, Yamada S. Experimental study on the hysteresis behavior of structural steel under multi-axial cyclic loading. J Struct Constr Eng 2005(588):181–8 [in Japanese]. [90] Lehmann T, Raniecki B, Trampczynski W. The Bauschinger effect in cyclic plasticity. Arch Mech 1985;37(6):643–59. [91] Nakajima N, Yamada M. Hysteresis loop characteristics of structural carbon steel SS400 in very large plastic zones. J Struct Constr Eng 2000(536):17–22 [in Japanese]. [92] Fujimoto M, Hashimoto A, Nakagomi T, Yamada T. Study on fracture of welded connections in steel structures under cyclic loads based on nonlinear fracture mechanics. Part 1: Formulation of multi-axial stress-strain relations of structural steel for cyclic loads. J Struct Constr Eng 1985(356):93–102 [in Japanese]. [93] Nishimura N, Ono K, Ikeuchi T, Shinke T. Experimental investigation on hysteretic behavior of structural steels in plastic range. Steel Constr Eng 1994;1(1):173–82 [in Japanese]. [94] Nishimura N, Ono K, Ikeuchi T. A constitutive equation for structural steels based on a monotonic loading curve under cyclic loading. Doboku Gakkai Ronbunshu 1995(513):27–38. I-31 [in Japanese]. [95] Simo JC, Taylor RL. Consistent tangent operators for rate-independent elastoplasticity. Comput Methods Appl Mech Eng 1985;48(1):101–18. [96] Wilkins ML. Calculation of elastic-plastic flow. In: Alder B et al., editors. Methods of Computational Physics 33. New York (U.S.): Academic Press; 1964. [97] Simo JC, Hughes TJR. Computational Inelasticity. New York (U.S.): SpringerVerlag New York, Inc.; 1998. [98] Doghri I. Fully implicit integration and consistent tangent modulus in elastoplasticity. Int J Numer Meth Eng 1993;36(22):3915–32. [99] Mahnken R. Improved implementation of an algorithm for non-linear isotropic/kinematic hardening in elastoplasticity. Commun Numer Methods Eng 1999;15(10):745–54. [100] Hartmann S, Haupt P. Stress computation and consistent tangent operator using non-linear kinematic hardening models. Int J Numer Meth Eng 1993;36 (22):3801–14. [101] Simo JC, Taylor RL. A return mapping algorithm for plane stress elastoplasticity. Int J Numer Meth Eng 1986;22(3):649–70.


Constitutive Model for Full-Range Elasto-Plastic Behavior of Structural Steels With Yield Plateau Formulation (Hu, Et Al. 2016)

Recommend Documents