Recent developments in the application of finite element methods to nonlinear problems

FINITE ELEMENTS IN ANALYSIS A N D DESIGN

z; '. " F'

ELSEVIER

Finite Elements in Analysis and Design 18 (1994) 1 15

Recent developments in the application of finite element methods to nonlinear problems D.R.J. O w e n * , D j o r d j e Peri6 Department of Civil Engineering, University College of Swansea, Swansea SA2 8PP, Wales, UK

Abstract

Computational strategies relevant to finite strain elasto-plastic problems are discussed and the necessity for both a sound theoretical framework and numerically consistent linearisation procedures emphasised. Strain localisation phenomena are briefly outlined and a material modulus for large inelastic deformations at finite strain is generalised to account for the effects of curvature of the yield surface, presence of ductile damage and rate dependence. The treatment of contact/friction problems within the framework of the theory of elasto-plasticity is also discussed. Numerical results are presented for some representative examples.

1. Introduction The last five years or so have witnessed rapid advances in the development of numerical solution procedures for nonlinear problems which can be attributed to progress on two fronts. Firstly, establishment of a theoretical framework for finite deformation mechanics and associated numerical implementation has reached a sufficient level of maturity and, secondly, advances in workstation technology have made the necessary computational power available. The formulation of finite strain plasticity problems has been the subject of vehement debate over the last decade with differing opinions on issues such as the notions of intermediate configuration, material frame invariance in conjuction with deformation gradient decomposition, decomposition of strain rates and related aspects, etc. Only recently has some consensus been reached on an appropriate constitutive theory based on tensorial state variables to provide a theoretical framework for the macroscopic description of a general elasto-plastic material at finite strains. In computational circles, effort has been directed at the formulation of algorithms for integration of the constitutive equations relying on operator split methodology. In particular the concept of

*Corresponding author. 0168-874X/94/$07.00 © 1994 Elsevier Science B.V. All rights reserved SSDI 0 1 6 8 - 7 8 4 X ( 9 4 ) 0 0 0 0 9 - 5

2

D.R.J. Owen, D. Peric'/Finite Elements in Analysis and Desiyn 18 (1994) 1 15

consistent linearisation has been introduced to provide quadratically convergent solution procedures. By employing logarithmic stretches as strain measures a particularly simple model for large inelastic deformations at finite strains is recovered. In particular, the effects of finite strains appear only at the kinematic level and the integration algorithms and corresponding consistent tangent operators for small strain situations can be directly employed. One immediate consequence of the ability to numerically simulate problems involving finite strain elasto-plastic deformations is the manifestation of strain localisation conditions for certain classes of problem. This localisation of deformation is a physical phenomenon characterised by the development of large strains inside a narrow band of the material, leading almost inevitably to a failure mechanism. The finite element modelling of strain localisation behaviour is currently the subject of considerable research, since the solution exhibits spurious mesh size dependence unless regularisation procedures are adopted. One way of avoiding such difficulties is to employ a rate dependent formulation through use of an elasto-viscoplastic model in place of the rate independent elasto-plastic approach. The formulation of solution procedures for finite strain elasto-viscoplastic problems follows that for the rate independent case and the use of logarithmic strains as a kinematic variable offers similar benefits. Practical applications involving both finite strains and inelastic material deformation often require the numerical modelling of other physical phenomena. For example, since the development of finite strains is frequently a precursor to failure, the modelling of failure mechanisms such as ductile damaging may be necessary. The inclusion of damage effects within a finite strain elastoplastic setting can again be accomplished through an operator split implementation based on consistent linearisation principles. A further class of nonlinear problem for which considerable advances in numerical modelling have been made in recent years is that of contact-friction behaviour. Contact-friction phenomena arise in many important practical areas and numerical treatment in the past has, of necessity, relied on temperamental and poorly convergent algorithms. This situation has changed markedly with recognition of the complete analogy that exists between contact-friction behaviour and the classical theory of elasto-plasticity. Hence, the operator split algorithms and consistent linearisation procedures developed for the latter case translate directly to contact-friction models to provide robust and rapidly convergent numerical solutions. The remainder of this paper discusses some of the above issues in more detail and provides numerical examples illustrating the recent advances made in the finite element analysis of nonlinear problems.

2. Finite strain elasto-plasticity

2.1. Constitutive equations for large strain elasto-plasticity It is only recently that the constitutive theory for finite strain plasticity, commonly based on the multiplicative decomposition of the deformation gradient, and the associated numerical integration algorithms, have received a complete rational basis. Notable contributions in this area are presented in [-1-9]. It is generally accepted that a sufficiently general constitutive model with possible wide application, including most of metal plasticity, together with considerable

D.R.J. Owen, D. Peric'/ Finite Elements in Analysis and Design 18 (1994) 1 15

3

simplifications in the numerical treatment of finite strain elasto-plasticity, may be defined by employing a quadratic strain energy function in the form of a scalar symmetric function of its principal stretches 21~),(i = 1, 2, 3) W(,~,~I), ,~,(2),e/]'~3)) = #((ln ,,~1)) 2 ~- (In ,~2)) 2 -~- (In 2e(3)!~2).._~ 12 (ln je)2

(2.1)

where /~, 2 are non-zero material constants and je = Z(1)Z(2)2(3) 0e "e e is the Jacobian. By denoting a logarithmic strain as E:= In U and its work conjugate stress as T:= RTTR the constitutive equations for large deformations of elasto-plastic solids at finite strains, with attention restricted to isotropic elasto-plasticity, can be represented in a form similar to the standard small strain plasticity equations as ( see [1, 2, 5 9] for details) F = FeF p

(2.2a)

T = c3e~W = 2 # E e + tc tr [Ee ] 1

(2.2b)

sym[ FV FP- a] = 9C~r Cb(T ' K )

(2.2c)

/¢ = 9H(E p, K)

(2.2d)

Here ~ ( T , K ) is the yield surface, K is the scalar variable describing isotropic hardening while function H(Er, K ) defines the hardening law. Finally, loading/unloading conditions may be formulated in the standard Kuhn Tucker form q> ~< 0,

9 >~ 0,

9~ = 0.

(2.3)

2.2. Numerical integration o f the large strain elasto-plastic constitutive model

Procedures for numerical integration of the rate type constitutive equations which are particularly suitable for implementation are connected with general operator split methodology, where the original problem of evolution is solved through composition, applying first the elastic and then the plastic algorithm. As a crucial point in the derivation of the algorithm, an exponential approximation is used in the discretisation of the plastic flow rule in the plastic corrector stage first employed in its original form in the computational literature by Weber and Anand [9] (see also [1,2,7,8]). It leads to the following incremental evolution equation FP+ 1 = exp [A])~T(~) n + 1] FPn

(2.4)

It can be easily shown that the above approximation results in Ee+

1 =

L]tTe ' n + ltrial

__ A ~ r ~ . +

1

(2.5)

which is valid whenever ~.-~n+ trial "--'--~Jt n+l/]~'e trial~T t r i--n+l al~"e l! and 0 r ~ , + l commute. In conjunction with the quadratic strain energy function expressed in terms of logarithmic stretches (2.1), Eqs. (2.2-2.5) lead to a particularly simple format of the integration algorithms for multiplicative finite strain plasticity which is equivalent to the small strain plasticity algorithms. For convenience the integration algorithm for the Jz version of multiplicative finite strain plasticity is briefly reviewed in the following.

4

D.R.J. Owen, D. Peric'/Finite Elements in Analysis and Design 18 (1994) 1 15

2.2.1. Elastic predictor

From the known configuration z.(X) assuming that the incremental displacement u(x,) is given, the updated configuration and incremental deformation gradient can be evaluated through )f,+l=X,+u

and

F,=I

+Grad[u].

(2.6)

From the logarithmic stretches, the elastic trial stress is evaluated as Ttrial n + l z 2/~dev [Ee+l] + ~c(tr [Ee+l"])l.

(2.7)

If the consistency condition d~trial "~n+l

= IId e v [ T ,trial + l ] l I -,~33K(ff?), <~ 0

(2.8)

is not satisfied than a plastic corrector phase is performed. 2.2.2. Plastic corrector

The plastic corrector phase is particularly simple for the plane strain case of J2-plasticity and consists of the radial mapping of the trial elastic stress onto the yield surface. This procedure may be summarized as ]Vn+ 1 ~-

dev r T trial q t_-,+ lJ trial Ildev[T.+,] II

(/~P)n-t- 1 = (EP)n "~

~32A~{'

dev[T,+ 1] = dev rTtrialqt-.+lJ- 2/~A 7N,+ 1.

(2.9a) (2.9b) (2.9c)

For success of the proposed scheme a proper linearization of the incremental equations is of crucial importance, which leads to the so-called consistent tangent elasto-plastic modulus. Due to the simple format of the above large strain elasto-plastic model, the results from small strain theories are directly applicable. Discussion on the recent advances and details of implementation of the small strain algorithms are given by Crisfield [10], Dutko et al. [11], Mitchell and Owen [12], Peri6 [4] and Simo and Hughes [13].

3. Contact-friction modelling Following standard formalism of the theory of elasto-plasticity, additive decomposition of the tangential velocity at the contact interface is adopted, i.e. gT- = ,~fr + g~

(3.1)

where gT = (I - N ® N ) " u is the increment of tangential displacement. Furthermore, a perfect friction law is assumed ( in the sense as introduced by Curnier 1-14] ) by stating that the friction force is proportional to the normal force and is independent of the other state variables, which leads to the slip criterion

4~ = I/pr II + vF lien I / - r

(3.2)

D.R.J. Owen, D. PeriC/ Finite Elements in Analysis and Design 18 (1994) 1 15

5

in which PN and Pr denote the normal and tangential force components, respectively, vF is the coefficient of friction and r characterizes the adhesion. With these assumptions introduced, the constitutive equations for frictional contact take the simple format of classical elastoplasticity gv = g) + gP

(3.3a)

PN = ON ue

(3.3b)

PT = D r ue

(3.3c)

ger = )~t?O(P, r)

(3.3d)

i" = )~h(ger, r).

(3.3e)

Here, D r = - k r ( I - N ® N ) and DN = - - k N N ® N are the tangential and normal parts of the elastic modulus tensor, 0(p, r) is the slip potential and the function h(g~-, r) defines the hardening (softening) law. Finally, loading/unloading conditions may be formulated in the standard Kuhn Tucker form 4) ~< O,

2 ~ O,

).q5 = O,

(3.4)

The effectiveness of the numerical implementation again relies on the use of a consistent tangent modulus as illustrated by Peri6 and Owen [15].

4. Further aspects of the deformation of inelastic solids at finite strains

The localisation of deformation is a physical phenomenon characterised by the development of large strains inside a narrow band of the material, without significant contribution to the overall deformation of the body. Although the onset of localisation is associated with the loss of ellipticity in a local material region ( see [16,17] and references therein ), the solution of the complete boundary value problem in a large scale numerical computation becomes a necessity. The physical mechanism of localisation is highly problem-dependent, influencing conversely the material description. For example, shear band formation in a thin sheet under uniaxial tension may be well described utilizing simple J2 flOW theory of plasticity as shown in [18], but, under biaxial tension prediction the onset of localised necking requires relaxation of the normality condition through vertex formation on the yield surface, or increase of the curvature of the yield surface as will be discussed below. In the presence of small strains, the loss of ellipticity and localisation of deformation must be initiated either by the lack of normality of the plastic flow, vertex formulation on the yield surface or through the explicit inclusion of strain softening. In the presence of large strains the situation is more complex due to the geometrical effects. In this section certain aspects of phenomenologically based refinements in constitutive modelling of inelastic solids at finite strains are discussed. In particular, the influence of the curvature of the yield surface is examined utilizing the anisotropic yield criterion recently proposed by Barlat and

6

D.R.J. Owen, D. PeriC/Finite Elements in Analysis and Desifn 18 (1994) 1 15

Lian 1-19]. Additionally, by introducing damage parameters the effect of continuous deterioration of material properties on the localisation of deformation is illustrated. Finally, the rate dependence of material is introduced within the viscoplastic constitutive model. Refined material characterisation is important in the accurate modelling of material behaviour and, particularly, it is considered essential for failure description.

4.1. Yield sur[ace representation For thin sheet metal forming operations, due to the processing of the material, the plastic behaviour of textured polycrystalline sheet is predominantly anisotropic. Use of the Hill anisotropic criterion, which contains no shear stress, is restricted to a planar isotropy or for the cases where the principal stress axes coincide with the anisotropy axes. Full planar anisotropy is described by the yield function introduced by Barlat and Lian [19] which for the plane stress state is of the superquadric form

f(~r, 6) = aIK1 + K21M + aJK1 - K21M + (2 - a)12K2l M - 2~ Ivt

(4.1)

in which

K1 - o,,~, +2 hayy and

K2

=

X/( axx--2

hffyY)2 .~_p2Oxy2

(4.2)

where a, h, p and M are material constants and 6 is the yield stress from a uniaxial tension test. For a given value of the exponent M material constants a, h, p can be evaluated using R values, i.e. plastic strain ratios of the in-plane strains to the thickness strain obtained from uniaxial tension tests in three different directions. The role of the material constant M is illustrated in Fig. 1 obtained by plotting the function (4.1) for various M values in the normalized oxx and ayy plane for the isotropic case (i.e. for a = h = p = 1) and taking ffxy = 0. The resulting set of functions span the set of yield surfaces which include the standard von Mises and Tresca yield surfaces for M = 2 and M --, oo, respectively. An important property of the yield function described by equations (4.1) and (4.2) is its convexity when constants a, h, p are positive and M > 1, as proved by Barlat and Lian [19].

4.2. Elasto-plastic damaging solids As experimentally verified for many materials (see [20, 21]) the nucleation and growth of voids and microcracks which accompany large plastic flow causes a reduction of the elastic modulus as well as material softening and can be strongly influenced by the triaxiality of the state of stress. In this case, the prediction of rupture with reasonable accuracy requires consideration of the coupling between elastoplasticity and damage. With attention restricted to isotropic solids, a scalar parameter D is introduced, which represents the density of micro-cracks and micro-cavities in the body. The elasticity law is then given by o" = C :Ee(1 -- D) where C is the elastic modulus which is assumed constant and isotropic.

(4.3)

D.R.J. Owen, D. Peric'/ Finite Elements in Analysis and Design 18 (1994) 1-15

7

2.0

:::

1.5

I /

M-2'

. . . . . . . .

i

. L.......................

1.0

0.5

0,0

-0.5

........!

-1.0

-1.5 - - ~ -1.5

-1.0

-0.5

0.0

0.5

1.0

1.5

0"1/~

Fig. 1. Isotropic yield surfaces for several values of material constant M and for shear stress axv = O.

Stresses are assumed to satisfy the von Mises yield criterion in the form f(o', 5, D) - ~

1-D

5

(4.4)

where a(Y) is the uniaxial yield stress. The evolution laws for internal variables are obtained in a standard way employing the second law of thermodynamics in the form of the Clausius-Duhem inequality and adopting the normal dissipation assumption (see [-21, 22]), which results in ~P

1 x/~ dev [a] =Y(1--D) 2 x/-)~ 1

(_ y~o

b= YZ--b\ So ]

(4.5a) (4.5b)

where So and So are material constants, and Y is the damage energy release rate expressible in terms of the stresses and accumulated damage, D. A computational framework for the fully coupled elasto-plastic damage model at finite strains is provided in [8].

8

D.R.J. Owen, D. Perid/Finite Elements in Analysis and Design 18 (1994) 1-15

4.3. Elasto-viscoplastic solids The constitutive model for elasto-viscoplasticity based on a v o n Mises yield criterion and power law isotropic hardening is represented as f(~, #) = F(~r) - #

(4.6a)

F(o-) = ~

(4.6b)

~Y= h0(e,o + gyp),

(4.6c)

~vP = 7(5u)O~f(~r, ~)

(4.6d)

g~'P = J H -1 = )'(7 ~)

(4.6e)

where ho, eo and n are material constants, g~v is the equivalent viscoplastic strain and H = d~(g~'P)/dgvp is the standard hardening modulus. Furthermore, the viscoplastic flow potential is introduced as 7'(o-,

~) =

-

1

(4.7)

where N is the material parameter representing the rate-sensitivity of the material. The evolution problem described by (4.6) and (4.7) is widely accepted as a description of rate-dependent deformations of solids and has a firm experimental basis. Its detailed study is provided in [-4].

5. Numerical examples In this section some illustrative examples are briefly presented and full details of these problems can be found in the references cited. A full Newton Raphson numerical solution procedure is employed in each case and a quadratic rate of asymptotic convergence is exhibited whenever the solution is within the radius of convergence. Example 5.1 (Failure analysis of Scordelis-Lo roof [23]). The geometry, material characteristics and loading conditions are shown in Fig. 2. Computation is performed employing standard-arc length control [10] and the localised mode of failure is illustrated in Fig. 3. A plastic hinge line first forms along the Xz axis, which is followed by a second hinge in the X1 direction across the full mid-section of the shell. Example 5.2 (Stretching of a circular thin sheet by a hemispherical punch. elasto-plastic material). The geometry and material characteristics for this example are shown in Fig. 4. The material is assumed to follow Barlat's yield criterion (4.1 and 4.2) where parameter M is varied. The solution for the standard von Mises material represented by M = 2 is also provided. To solve this problem the blank is discretised by 736 constant strain triangular finite elements, and 2145 and 612 triangular flat elements are employed to discretize the surfaces of the punch and die, respectively.

D.R.J. Owen, D. Peril~Finite Elements in Analysis and Design 18 (1994) 1-15

Geometry:

9

-•2L.•

Lo = 7.6m R = 7.6m

2B.

t = 0.076m 0 = 40 ° Material properties: E = 2.1 × 1 0 4 N / r a m 2

u=O.O ay = 4.2N/mrn ~ Loading conditions: Self weight, qo = 4 . 0 k N / m 2

Fig. 2. Failure analysis of Scordelis-Loroof: Geometry, material characteristics and loading conditions

Results are obtained for a coefficient of friction between tools and blank of vv = 0.30. Fig. 5 gives the deformed mesh at punch displacement Dp = 40.0 mm. The punch force versus punch displacement diagram, presented in Fig. 6(a), gives comparison between results obtained for various M-values. The maximum punch force decreases with increase of the M-value, which indicates a strong influence of the curvature of the yield surface on the initiation of strain localisation and its development. The distribution of true strain in the radial direction is shown in Fig. 6(b) for various M-values and for punch displacement Dp = 35 mm. A typical localisation behaviour may be observed where strain accumulates in a narrow zone, reaching high levels and leading to failure. The appearance of localisation and associated failure is less pronounced with decrease of the curvature on the yield surface, specified by decrease of M-value.

Example 5.3 (Stretching of a circular thin sheet by a hemispherical punch: elasto-plastic damaging material). The geometry and material characteristics for this example correspond to those shown in Fig. 4, with a continuous damage described by material constants So = 1.0 and So varied between 250.0 and 1000.0. The punch force versus punch displacement diagram, presented in Fig. 7(a), shows that the maximum punch force decreases with reducing So-values, with the falling branch becoming progressively steeper. At the same time, the maximum punch force is attained at significantly lower values of the punch displacement, indicating that higher values of damage energy release rate promote localisation and significantly decrease the bearing capacity of the material.

10

D.R.J. Owen, D. Perid/ Finite Elements in Analysis and Design 18 (1994) I 15

a

b

d Fig. 3. Failure analysis of Scordelis-Lo roof: Deformed finite element mesh (16 x 16) at various stages of loading. Each quadrilateral consists of two triangles. (a) q/qo = 1.335, (b) q/qo = 0.505, (c) q/qo - 0.280, (d) q/qo = 0.351.

The distribution of true strain in the radial direction is shown in Fig. 7(b) for several values of material p a r a m e t e r So, a n d for punch displacement Dp = 3 0 m m . The underlying elasto-plastic solution with no d a m a g e in the material matrix is also presented, indicating the significant effect of d a m a g e energy release rate on the localisation. E x a m p l e 5.4 (Stretchin 9 o f a circular thin sheet by a hemispherical punch: elasto-viscoplastic material [16]). The g e o m e t r y a n d material characteristics for this example correspond to those s h o w n in

D.R.J. Owen, D. Peric'/ Finite Elements in Analysis and Design 18 (1994) 1-15

~

~

/

I ~

,-~ blank

xl

to ~

~__ d

die i

,

11

Geometry: Ro= 59.18mm to =-1.0mm Rp-- 50.80turn Rd=6"35mm Materialproperties: E-- 69004N/mm2 u=0.3 -=c, 589. (1.0.10 -4 + ~P)°'216N/mm2

2R.

Fig. 4. Stretching of a circular thin sheet by a hemispherical punch: Geometry and material characteristics.

• Dp = 40.0 mm

Fig. 5. Stretching of a circular thin sheet by a hemispherical punch: Deformed finite element mesh at punch displacement Dp = 40.0 mm.

Fig. 4, with elasto-viscoplastic material behaviour described by the viscosity parameter 7 = 0.002 s-1. The punch force versus punch displacement diagram, presented in Fig. 8(a)-(b), gives comparison between results obtained by changing the loading program and the rate sensitivity of the material. The loading program was altered by varying the punch speed between Vp = 0.5-500.0 mm/s, increasing the speed by a factor of ten for every new test. Tests were performed for the rate sensitivity exponents N = 10, 20, 50 and 100. For comparison, the underlying rate independent solution is also presented. The major influence of the rate of loading, demonstrated by the effect on the maximum value of the punch force, is observed for the highly rate-sensitive material represented by N = 10. The distribution of true strain in the radial direction is shown in Fig. 8(c). for the punch speed Vp = 0.5 mm/s and the rate sensitivity exponents N = 10, 20, 50 and 100 for punch displacement Dp = 40 mm. For the rate independent solution a typical localisation behaviour is observed where strain accumulates in a narrow zone, reaching high levels and leading to failure. This trend is closely followed by the low rate-sensitive material represented by N = 100, while the appearance of localisation is completely averted in the higher rate-sensitive material characterized by N = 10.

D.R.J. Owen, D. Peric'/ Finite Elements in Analysis and Desiqn 18 (1994) 1 15

12 80

;(.)

70

'°q(b) __

LEGEND: =7o, M:

•

.

.

.

.

1

.

...... :: M=4

0.8

60 Z LL

i

50

.~

t

:

~

Exponent

it

. . . . . . . . .

M:/

t

:::... . . . .

0.6

N

O 40

g

IT 0.4

g 3o

O-

20

0.2 .~,~

10

0.0 ~ 10

20

30

40

0.0

0.2

0.4

Punch travel, Us (mm)

0.6

0.8

1.0

Initial p o s i t i o n , R / R o

Fig. 6. Stretching of a circular thin sheet by a hemispherical punch: (a) Punch force versus punch displacement curves for various M-values. (b) Distribution of true strain in the radial direction plotted over the initial configuration for various M-values at punch displacement Dp = 35.0 mm.

0.8

70

--

t

6o :I

Z LL

parameter:

I . . . . . . s=1000.

4

I ...... s= soo.

...s _2so.

L

I .....

j

¢:

50

o

/

40

i

•

-

._= 0.5

'

" .:

•.~ 0.4

%... :

IT

0.3

20

0.2

10

0.1

"%

•

'

i

I0

. . . .

i

. . . .

20 Punch travel, U3(rnm)

parameter:

...... ... S=1000 o o o ... S=500 - - , - ... S=250

-%

g 3o a.

0

Damage

Pi @i !i !i !i Ji !i !i !L

0.6

/ , .......

LEGEND: No damage

.I

I

cff o

0.7

i ....... Damage

(b)

i

30

. . . .

I

40

0.0

. . . .

0.0

i

0.2

. . . .

i

0.4

. . . .

i

. . . .

0.6

i

0.8

. . . .

i

1.0

Initial p o s i t i o n , FI/Ro

Fig. 7. Stretching of a circular thin sheet by a hemispherical punch: (a) Punch force versus punch displacement curves for several values of material parameter So. (b) Distribution of true strain in the radial direction plotted over the initial configuration for several values of material parameter So for punch displacement Dp = 30.0 mm.

D.R.J. Owen, D. Peric'/Finite Elements in Analysis and Design 18 (1994) 1 15 90

=+]ia)-

-

•'1 l

":J I

2

o Punch spe~
•Im> Ii ............ j --........... s.om}'°mm'° Ioo - - ..... o.~

~-,+o

.....

;

- - - - .....

ram, m~,m,d.~:

I

:l

.......

T. . . . T LEGEND:

80

/"

/

(b)

.....+/ ~ . ~ - - --.--. - - - - + ...... -- ....

-

-

=i

I

i

P,m , 3 ~ - , , ~

no~:

..

|

. +~..,c~mx,a:J{mm.~}: I

......... mo.o

~

...... 5.o A60 Z v LL -50

/---:;r--?

13

i-

+

+

;

+ /

i

i ~ - ~

/..+;:.2.2 i /..'+,..~/..:-.~ + +

i

i

,2

~ 100

j:: 40 30

+

+

i

20

........ _ .......... + - - + + . + . . . , . . ,

+

....

+ ...........+ , +

T ~ ~ _ ....

o 0

u .... 5

, .... 10

..,,,=,o , ....

15

i .... 2O

i .... 25

i ....

I ....

3O

.... !...;

,

35

5

40

10

15

(c)

!

• Dp

0.9

=

40.0mm

r ~ ....

~

30

35

40

uEaE.o:

i J Rate s e r l ~ y : ........ N-10 ! / ....... N-20

..... ~-'[

~ ~] ...... ......

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

0.8

.... + 25

Punch travel, Us (ram)

Punch travel, Us (ram) 1.0

.... 20

~/--

~

N-50 N-100 ~

0.7 +---~

0.6

-~ 0.5

0.4 0.3

0.2

0.1 +

0.o

! 0.0

0.1

0.2

0.3

0A 0.5 0.6 Initial position, R / R o

0.7

0.8

0.9

"1.0

Fig. 8. S t r e t c h i n g of a circular t h i n sheet by a h e m i s p h e r i c a l p u n c h : ( a ) - ( b ) P u n c h force versus p u n c h d i s p l a c e m e n t curves for v a r i o u s p u n c h speeds a n d several values of the rate sensitivity e x p o n e n t . (c) D i s t r i b u t i o n of true strain in the radial d i r e c t i o n p l o t t e d o v e r the initial c o n f i g u r a t i o n for p u n c h speed Vp = 0.5 m m / s a n d several values of the rate sensitivity exponent.

14

D.R.J. Owen, D. Peric'/ Finite Elements in Analysis and Design 18 (1994) 1-15

6. Concluding remarks Some recent advances in the finite element analysis of nonlinear problems have been reviewed, indicating the progress that has been made both in the theoretical understanding of inelastic material behaviour and associated numerical implementation. It should be, however, emphasised that significant developments are also occurring in other areas of computational mechanics. For example, the use of parallel processing techniques and associated hardware for the solution of finite element problems is currently an active research topic [24]. Parallel processing offers a natural approach to improving computational power and the next generation of workstations will increasingly employ multi-processor architectures to achieve significantly advanced performance. Another area of fruitful research involves the development of discrete element techniques in which the individual deformable elements are in frictional contact only [25]. The use of combined finite/discrete element approaches is potentially attractive for the simulation of fracturing solids; particularly under dynamic conditions.

Acknowledgement The authors are indebted to M. Dutko, G.P. Mitchell and E.A de Souza Neto for their involvement in various phases of this work.

References [1 ] A. Cuitifio and M. Ortiz, "A material-independent method for extending stress update algorithms from small-strain plasticity to finite plasticity with multiplicative kinematics", Eng. Comput. 9, pp. 437~451 (1992). [2] A.L. Eterovi6 and K.-J. Bathe, "A hyperelastic-based large strain elasto-plastic constitutive formulation with combined isotropic-kinematic hardening using logarithmic stress and strain measures', Int. J. Numer. Methods Eng. 30, pp. 1099-1114, 1990. [3] R.A. Eve, B.D. Reddy and R.T. Rockafellar, "An internal variable theory ofelastoplasticity based on the maximum plastic work inequality", Quart. Appl. Math. 48, pp. 59 83, 1990. [4] D. Peri6, "On a class of constitutive equations in viscoplasticity: Formulation and computational issues", Int. J. Numer. Methods Eng. 36, pp. 1365-1393, 1993. [5] D. Peri6, D.R.J. Owen and M.E. Honnor, "A model for large strain elastoplasticity based on logarithmic strain: computational issues", Comput. Methods Appl. Mech. Eng. 94, pp. 35-61, 1992 [6] D. Peri6 and D.R.J. Owen, "A model for large deformations of elasto-viscoplastic solids at finite strains: computational issues", in: E. Stein and D. Besdo (eds.), Finite Inelastic' Deformations Theory and Applications, Springer, Berlin, 1992. [7] J.C. Simo, "Algorithms for static and dynamic multiplicative plasticity that preserve the classical return mapping schemes of the infinitesimal theory", Comput. Methods Appl. Mech. Eng. 99, pp. 61-112, 1992. [8] E.A. de Souza Neto, D. Perib and D.R.J. Owen, "A computational model for ductile damage at finite strains: Algorithmic issues and applications", Eng. Comput. ll, pp. 257-281, 1994. [9] G. Weber and L. Anand, "Finite deformation constitutive equations and a time integration procedure for isotropic, hyperelastic-viscoplastic solids", Comput. Methods Appl. Mech. Eng., 79, pp. 173 202, 1990. [10] M.A. Crisfield, Non-linear Finite Element Analysis of Solids and Structures. Volume 1: Essentials, Wiley, Chichester, 1991.

D.R.J. Owen, D. Peric'/ Finite Elements in Analysis and Design 18 (1994) 1-15

15

[11] M. Dutko, D. Peri6 and D.R.J. Owen, "Universal anisotropic yield criterion based on superquadric functional representation. Part 1: Algorithmic issues and accuracy analysis", Comput. Methods Appl. Mech. Eng. 109, pp. 73-93, 1993. 1-12] G.P. Mitchell and D.R.J. Owen, "Numerical solution for elastic-plastic problems", Eng. Compup. 5, pp. 274~284, 1988. 1-13] J.C. Simo and T.J.R. Hughes, "General return mapping algorithms for rate-independent plasticity", in: C.S. Desai et al., (eds.), Constitutive Laws for Engineerin 9 Materials: Theory and Applications, Elsevier Science Publ., New York, pp. 221-231, 1987. [14] A. Curnier, "A theory of friction", Int. J. Solids Struct. 20, pp. 637-647, 1984. 1,15] D. Peri6 and D.R.J. Owen, "Computational model for 3-D contact problems with friction based on the penalty method", Int. J. Numer. Methods Eng. 35, pp. 1289-1309, 1992. [16] A. Needleman and V. Tvergaard, "Finite element analysis of localization in plasticity", in: J.T. Oden and G.F. Carey (eds.), Finite Elements V. Special Problems in Solid Mechanics, Prentice-Hall, Englewood Cliff, N J, pp. 94-157, 1984. [17] J.R. Rice, "The localization of plastic deformation", in: W.T. Koiter (ed.), Theoretical and Applied Mechanics, North-Holland, Amsterdam, pp. 207 220, 1977. 1-18] D. Peri6, G.P. Mitchell and D.R.J. Owen, "Some computational issues in strain localisation problems for rate independent solids", Int. J. Numer. Methods Eng. 32, pp. 79 101, 1991. [19] F. Barlat and J. Lian, "Plastic behavior and stretchability of sheet metals. Part I: A yield function for orthotropic sheets under plane stress condition", Int. J. Plasticity, 5, pp. 51-66, 1989. [20] J. Lemaitre, "A continuous damage mechanics model for ductile fracture", J. Eng. Mater. Techn. 107, pp. 83 89, 1985. 1,21] J. Lemaitre and J.-L. Chaboche, Mechanics of Solid Materials, Cambridge University Press, Cambridge, 1990. 1,22] J. Lubliner, Plasticity Theory, Macmillan, New York, 1990. [23] D. Peri6 and D.R.J. Owen, "The Morley shell element for large deformation problems: simplicity versus sophistication", in: N. Bi6ani6 et al. (eds.), Nonlinear Engineering Computations, Pineridge Press, Swansea, 121-142, 1991. 1,24] G.F. Paw and D.R.J. Owen, "Domain decomposition procedures for finite element analysis on transputer based parallel computers", Compup. System Eng. 2, pp. 451460, 1991. [25] D.R.J. Owen, A. Munjiza and N. Bi6ani6, "A finite element-discrete element approach to the simulation of rock blasting problems", in: FEMSA92 - 11th Syrup. Finite Elements in South Africa, Cape Town, 1992.