Published by SUT, Shanghai, China
Applied Mathematics and Mechanics (English Edition. Vol. II. No. I, Jan. 1990)
ON U N I L A T E R A L
CONTACT LARGE D E F O R M A T I O N
FRICTION ( I I ) - - N O N L I N E A R
P R O B L E M WITH
FINITE E L E M E N T TECHNIQUE
A N D ITS A P P L I C A T I O N
Shang Yong (l:':J ~ ) Chcq Zhi-da ([~i~ff.i~) (Beijiug Graduate School, Chimt University oJ'Mining, Beij'iug)
A~bstraet Based tm the v(a'iatit2aal cquatiou derived in re.['. [1]. a nonlinear #wremental F.E. equation is formulated .[br unilateral eoutaet ehlstic atul plastic large deJormation problems. A m'w technique- co.movhtg coardinate.[bffte element method is #ttroduced, attd tt practical mathematical modelJbr large &J'ormation contact problem is ok,scribed. To show the ~Jfectiveness of the method, problems of contact htrge dt:[brmation of caatilever beam, circularplate, as well as metalring are computed. Compared with e.~7~erimcats, the results show good a,ffre6'meals.
I.
I n t r o d u e t l o n of CMC N o n l i n e a r FEM
It is certain that the reference configuration is very important for defining most of the variables ia mech'mics of large deformatiotl, l-or large deformation FEM, two desc;iption techniques are commonly used. One is theTotal-Lagrangian (T.L.) technique, by which initial configuration is chosen as reference bases, the other is Updated-Lagrangian (U.L.) techrique, which considers ( i - I)th incremental step.as the reference configuration of ith in~_remental step. U.L. method is in nature of T.L. method by each incremental step, and its main difference from T.L. lies on reference configurations updated for each step. Co-moving coordinate (CMC) systeln is composed of two coordinate systems, i.e. fixed coordinate system and dragging-along (or intrinsic) system. The latter is embedded in the deformable body and varies with body deformation. It is more reasonable to describe nonhomogeneous deformation field by virtue of local base vector :n co-moving coordinate systeJn. As large displacement is progressing, most. of the elements will undergo large translation, large rotation a'ad large deformation. Therefore, how to increase the computing accuracy and efficiency becomes one of the important problems concerned in computational mechanics. It is clear that a key to the increase ofcomputingspeed and accuracy is to apply the reasonable finite deformation theory and correct descriptiol.l.tcchnique in nonlinear F.E. procedure. The proposition of CMC nonlinear FEM is just for the sake above-mentioned. In our research, the successful results have been shown in the computation of large deformation problems for metal forming, rock deformation and structtu'al elements. before deriving the incremental finite clemeut equations, it is worthwhile to summarize shortly the main features of c M C nonlinear I-EM: a. Taking advantage of the fact that co-moving coordinates are intrinsic wu'iables in dcfor!ning process, direct interpolate method for discrete elements can be applied with respect to current configuration. Actually three coordinate systems, i.e. spatial fixed, co-moving and element local I
2
Shang Yong and Chen Zhi-da
coordinate systems, are used i,1 the technique, the relation between these can be shown as:
X-_QXx_-_QXQO~ where QKand QO mean ',,:'aematic and geometric translbrmations, and X,x mid /5 are spaiial fixed, co-moving and ele,n:nt local coordimdes rcspectively. b. Let g~ o r 9~tlenote covariant and conlr;wariant vcctors in co-moving coordinate system respectively. Nodal displacemclat and vclocily vcctors can be cxprcssed as: U ~f|J
yj
V ~ u j gl
Besides taking uj(o,'~j), displacement (or velocity) COtllpOllents in Cartesian coordinate systenl as fundamental interpolating quantities like T.L. and U,L. techniques, contravariant components of nodal displacem'ent (or velocity) vector in torrent co-moving coordinate can also be taken as interpolating quantities. In this way, computing process may be simplified for large deformation problems. c. In the new technique, the reference is the end configuration of a current step. Based on that, displacement component can be written as: (I}$ II
t(-[~
u
-" - - r - ' l + a .
(d)/
)'mii'l"
xx~
v~
a,
.
T
g* fo-t
Vl
~.
3
Fig. I
Fig. 2
where (I)
(t '-:-I)
(d)
F ] is Lhe dcfornlation gradient in current (i)th step. d. During deforming process, the effect ofconfi'guration on 10ad can be updated step by step, and the shortcoming of 'dead load' appearing in most of the program is overcome. e. For frictioll contact boundary, it is convenient to define tangential slip displacement and friction force along basic vectors of the coordinate system.
lI.
Nonlinear I n c r e m e n t a l F.E. Equations
a. Nonlinear geometric m a t r i x In current conliguration, il)terpolation forms of velocity, co-moving coordinate and spatial fixed coordinate in an element arc expressed as:
On l.!,filato'al Contact Large Dcl'tmuatioll i)roblcm wilh Friction (11)
m,
v ' = }"~N,(~J')o.l
_
3
,tt,~' [
|-I
m
(2.1) &-I
.~:,= } ,,N,(~j).,:I A " = \--' Nt C,U)A" ;
i
where N(~') is intc,'polation function, ,*' means clement nodal functional v.'due.. According to S-R theorem, the physical component of finite strain ,'ate is,
Si
(5' IIJ-- fi'll ',')
(2.2)
and its matrix form i~ (the upper symbol ' h ' means physical component of tensor and will be omitted in the following sections for simplification), t,~}--{S'I, S'], S'], 2,.~'1, 2,-~'~, 22~}T Note that the definition of
(2.3)
v'll., is based on current configuration in co-moving coordinate system .v'llj=v ; j-I- l"~,tv"
(2.4)
where F ' j is Christoffcl symbol. By eqs. (2.1)-(2.4). it iseasy to get the expl'cssion for t3"} in the form of Iv}: which is a physical component of an element nodal velocity along the current co-moving coordinate
IS} =
[133 {v} = ([B]L + [ B ] # ) Iv}
(2.5)
where [B]zis the linear geometric matrix we are familar with, and [Bills the nonlincar geometric matrix derived in this paper. Its ,rob-matrix can be shown as
[B,]#=
FI, FI2
l'~r I'l~
F;, 1"12
l-'lz+[':l
F'~z-i-l"ll FI3"FI'Iz FI3.FF~I
FI2+FI, 1"13+F12. FI3+F|t
Fl~+I'~2 FI3+F~, It is not difficult to cstimate I'~j
9N,
(2.6)
by mnnerical solution, tile general expression is ~kr Ogk,
Ogu
0.,r k
(2.7)
f o r / = i ( o r / = j ) ; the above equation will be simplified into
1 g.
F~=-~where
O,\'* U~J----- 8X~
is matrix tensor, according to cq. (2.1)
Og,t
O.~:*
OX z &~ O.,r
(2.~)
4
Shang Yong t!nd Chert Zhi-da
ox'
N,
),
i
in
l'!
Og,j
___.) O ' X ~
OX ~ ,
It should be mcntioncd that the following transformation must be considcrcd in computation: {Nz({]'),j}={1/9,,f~-~(jj~}~'[J]-'{Nl(es ' ) , $ , }
(2.9)
where I f ' l = [Ox~/O6 J] is geometric Jacobian matrix. The nonllinear geometric matrix given in the previous part follows the.fact ofcontravariants of velocity vector as interpolation quaqtities. Altermttively, Letting velocity be decomposed along the patial fixed coordimlte system, following the same procedure, another form of noalitmar geometric matrix caa be obtained, too. By relation
we can prove that
o'tl~=, u;'-' ~ Substituting it intoEq. (2.2), we have
{S}=[B]{~}
(2.10)
a sub-matrix of ,mnlinear gcomctric matrix [BJ is cxprcsscd as -]Vt t l r l i
I'B,] =
N~, zrl2 N, ,.at, 3 N,,~rll+N,,lrlz Nl,arn+N,,zrl.~ -N/,3rtl-bN,,lr,t
N i ,p lF:,l
.Nil ! 11"31
N~, zr2, JV~,3r23 N,,zr2~+N,,lr2z N,,arz,q-Nl,zr,a N,,.~rzlWN,,lr2a
g i , zr3z N,3, rz.~ N,,zrstWN~,lr3, N. ,,3rs2+N,,zr~3 N',,,~rs,-t-N,,tr3z
(2.11)
whcrc r,j is an element of [F]-",
IA-X-I Ox .1
[F]=L III.
Elasto-Plastic Constitutive
Matrixes
Two forms of clasto-plastic constitutiv~ equations are used ill this paper ill order to treat deformation in different case. Ill this section, we only show the matrix expression for computation in FEM. a. The e l a s t o - p l a s t i c c o n s t i t u t i v e e q u a t i o n w i t h single i n t e r a a l v a r i a b l e for finite deformation proposcd by Cheu Zhi-da (1964):
On Unilateral Contact Large Deformation Problem with Friclion (11)
5
r.i = 2 G ( 1 - - f l ) i ' l
e)=
(3.1)
I --2v
where r ; is objective rate of stress deviator, and v are elastic constants. Eq. (3. I) gall be rewritten as
i'}
)
,|k
is the rate of finite strain deviator, and G,E
";
Vi =C,r5 [
(3.2)
where CI~---- 1 E --v
i ] ( t _ l , ) , , ; , , : + ( , ~9 ,,
_k/_~_)6;,~l)] ,
1 f i = l --" 3G
(3.3)
dcr de
The matrix form of (3.2) is
{V}=[C](S'} Letting /3=0 , [C] will become an elastic constitutive matrix. b. F i n i t e d e f o r m a t i o n .I., f l o w t h e o r y la) This theory is due to Hutchinson, the constitutive law may be expressed as ,
Vi=
~."
l--v
,5, t- I '-z) ,a'cS'-"--'-<'1le' q - , - i j,,, t
( 5 , ! _ ~
i,
~
The first two terms just correspond to the case
(3.,1)
f l = 0 .inEq. (3.3).For matrix expression, we have,
{V}= (EC]o- EC],) {.~}
(3.5)
where [C]o is the elastic part, [C]v is tim plastic part, L' [C]~,= 1-t-v 9
a {r}.{r}. r q
h
(3.6)
9
(2= (. L--I-l,) ~---I- -~-~o~ 2
W=
~,
~-,E,=at
0
in c!astic region
1
in plastic region
On the otl.,:r hand. we have tile formula of objective rate of stress 171 9
,
"t.
,
V ~ = o'~-t-cri6 j - - a k S
r',k
~
(3.7)
It shows the time rate of stress only, and the other terms show the influence of the change of configuration on stress rate. By changing cq. (3.7) into a concise alternative form:
6
Shang Yong and Chert Zhi-da
*;=v;-I-AI~S'~
(a.Ta)
where
-],j" = ~ . ('~"'~I+ L'ri' -o;,5~ - 6 ' a l ) z
we usually use A to represent its matrix form. Fitmlly, Ihe matrix expression for ./clasto-pktstic large deformatiol~ constitutive equation can b e sIiOWll .'is
~a,t = ( [ c ] , IV.
[c],-t- [A]) iS' I- = [C] IS'I-
(3 .a)
Nu,nerical Model for Large Deformation Contact Proble,n
Let a COlll:l.Cl.system be composed of three objects. Body A is a rigid contactor, del'ornmble body B is tile target, and C is the rigid basis. For the convenic,lce of dcscribillg the motion. Ihree COol'di.ate systems are used. A rigid coordinate system {y,~} is fixed in A, co-movi,g COol'dilmte system IXI is embedded in B, aud {X} is selected as a fixed spatial rcfci'cnce system.
G,
i'ig. 3
Fig. 4
With respect to [ Y.,}, tile surface equation of contact botmdary of rigid body A is expressed as
Transforming it imo [,V}-system, we have:
XJ=X,(u,,
,f , r
yz), t)
In giving a displacement {A u j }of b o d y / I , then at ith step, we have I| ~
I'|-
I'
tl ~
A" j =-- ,\" ~-I- c\ u ~
(,1.1)
(1)
Similarly, a corresponding displacement A u ~. tm the boundary ol" Ii will occur simultaneously; the fixed coordinate value at a poim on B b~,tmdary will be (I)
~I-
l ~
II)
.\'~--- x ~ + A ,, '.
(4:2)
In order to make the study easily, toni:tel I~otmdal'y may be devMcd into two parts. One is called contacted region, .'rod the olher prc-conlact region. In I~l'e-contact region, an important problem is to set up the contact criterioa. In our discussion, coustraining condition for contact deformation can be defined as:
On Unilatcral Contact Large Deformation Problem with Friction (11)
njAu~ -- n.sAu,h --g~O
7
(4.3)
As shown in Fig. 4,. g means tile dista,lce between two corresponding points (m and e)on two contacting boundaries and can be determined by the following method: taking a point m on precontact boundary of body/L let L bc the commot~ tangc,~t I-Phmcat the criterion point C'. The length ol'normal line between two comact.ing bodies is tlcfincd as thc gap distance at m, i.e. 5 , = m a , which is a l'tmclion of spatial coordilmtC lind time. .-:h,
N
ma I ,n6". l
is the direction cosiae of normal. In r ?oints is
region, contact criterion for two corresl'Jonding
g=O
(4.a)
When a polar on body B comes intc, cc,ntach its motion will be subjected to the constraint of boundary curve ofrigid body A or basis C. and ils dimension will be changcd fi'om independcnt to dependent. As I\-~r tim treatment Ibr large deformation contact boundary, most skills will be included ia the program, and we will not have a more detailed discussion here. V.
Nonlinear Incremental
F.E. Equation
In Part l Ift, we have obtaincd the variatiomd inequality f'or large deformatio,i contact problem with friction:
a(v, v ' - v ) - ~ ( v ' ) - f ~ ( v ) +
(K,v~-O), v',-v.
< ( / , u'-o)+(b(v,), v'--v)
(Vv'EU)
(,5.1)
where ~k ( v ' ) rci~rescnting friction power, is usually a nonconvex and tmn-diffcrential functional. For convcl~ie,~ce ol'comfmtation, it is necessary Io change the inequality i,~to variational equation in common sense. The following approximate treatment is ad.optcd. Let fi'ictio ~ between two contact surfaces obey Coulomb's law (although it is vcry approximate). Then we have
~)k(V')= I 9
"[--It.cr,,IrIv', IclFc, I'cz
(d-l)
Fui'thcr. if we use cy. . which is defined on the previous step. instead of crn in the above equation. ~be will become diffcrcntial. Besides. the influance ~rf varying contact boundary on variation can bc considered in the proccss of program design instead of appearing in vari.'ltional equation. So Eq. (5.1) can be rewritten as: ,(v,
,'iv I - I - ( -1- ( l e , , ' \ , ' - - a ) , , ~ ' , , ) = ( f ,
&,)-I- Jl .
I ', "r ," 3 l v , l d F , . ~
r
According to Part I of this paperl'l, its full CXl',rcssion is
9
l'c+~
:
.
"t
(5.2)
8
Shang Yong and Chen Zhi-da =J
TJ6v~dFI+
~p f s ' 3 ~
u o" [~'~.~'-116S~dD
+ [ ua.~ I v, l dF,., 3 / "c1 where / - ' o = F ~ , + I ' o , , pemdty fi~ctor,
(5.a)
Fo~ is pre-contact region and
~'1)
fl)
II-I)
/-'~z is contacted region,
e
is a
OX k ax *
In co-moving coordinate system, the rclalion between load components in current increment and that in initial configuratio,i is cl~
tl ~
(I)
CIl)
,pfj= .pf~l,'~
,T j = .7'J;'~,
(a)
where 'i" and 'o' represent current and initial increments respcclively. It should bc mentioncd that the effect on load owing to configuration change is considered by introduci,lg dclbrmation gradient 1;'~ , Substituting (a) iuto (5.3). it is not difficult to rearrange it as matrix form : /'
C
t'j
=
J'l.] { . T(''k r [ F ] r [ N J d l ' / q --
-
l'fd
{ . p fTM }r[l;'jr[NJd~2
J'.~"-"}~ tr [ / 7 ] r [ B ] d ~ - l - +I
j
+ 'l
Ya
g{q}[~,; ]dl ' ~
#{ ( ic- Or . } r l r } [ N. ~l d F o ~
where IN1 is the shape function matrix on contact boundary, {q} is' the vector of normal direction and {r} is the vector of tangent d:rcction on Fo~ and Foz.Transporting in both sides of the equation, we will obtain the incremental finite clement equation for elastic and plastic contact large deformation with friction:
([K]+
rE~.)~a,:r = ~/er - ~or - U e r . - U e ~ ,
(5 .a)
where [K] is the stilTness matrix contributed by the coupling of geometric and material nonlinearity, [_/:E1, isthel~enalty stiffness malrix.{ Ir {A~}o a.dl 1~}/are the cquivalent nodal force vectors contributed by external load, slress in ( i - I) slops, pelmlty term and friction respectively. VI.
Numerical
Examples
Based on tile theories and the leclmklues hr tile previous sections, a nonlinea, finite element p,'ogra,n has been complied, and some typical elastic and plastic large deformation contact problem~ are numerically solved. a. E l a s t i e c o n t a c t ; l a r g e d e f o r m a t i o n
of cantilever beam
Geometric sizes arc shown in Fig. 6, material conslants: E= 1000kg/cm, v = 0 . 3 . In this example, we use plane stress clement and i0crcmcntal wu'ying stilTness algorithm under given load. The whole deformation process is devidcd into 20 iucrctuental steps, and 4 steps are used from the
On Unilateral Contact Large .Dclbrmation I'roblcm with Friction (11)
9
begi,ming to the initiation of free end contact. The curve of contacted boundary length with load and geometric configurations during deformatio,1 process are show. in Fig. 5 and Fig. 6. q '
~(kg/cm) 6
'l~!:l
t.J
~',
.I
I '."
Icm
I L.-I:! _ir/.,a
5
8
(a)
3
1
Fig. 5Curve of load-contact boundary length
i
,~.~ocn,
1--
-
9
-
(b)
I
';='~. 8kg/cm s
"
9
l
.
-'/iXItSE I-I-FM-EJSI-FccI-I-~%_r~,~ .l N
,
I
"" ~ " " " . . . .
I. . . . . . .
'.
~
2 I-'~
I "j
1~' "e)
I
(a) ii .
~
Fig. 6 Conligurations in different deformatio, stages ofcantilever I~ca,n(friction coefficient between beam and basis 0.1)'
.
I : ~ . . ,..,.__,_~_.~_~_ . . . . . . . I
b)
/_
./l/-d-,"
/
/
/
/
/
7
_/ I
S
L trig. 7 Coxnputing configurations of thick circular pl:tte contact large defortnatio,i
Fig. 8 Metal ring in compression
b. C o n t a c t l a r g e d e f o r m a t i o n of t h i c k c i r c u l a r plate In this example, axisymmetrical continuum clement is used. The plate whh well-distributed loadis assumed with sliding support in horizontal direction. Up to final configuration shown in Fig. 7(a), 20 increment steps are used by incremental varying stiffness alogrithm with given load. ( as shown in Fig. 7) c. Radial c o m p r e s s i o n e l a s t o - p l a s t i c large d e f o r m a t i o n o f m e t a l ring This is a typical problem of both tmilatc,'i:d varying boundary contact and coupling of nonlinear geometric, material and boundary condition. Yclla ;.Ind Reid (1980)1-'1studied instability of metal ring with radial compression by expc,'iment. Talay and Bcckerl31 have computed this p,'oblem by using rate-dependent plastic constitutive equation. Comparcd wRhRef. [3] and the
10
Shang Yong and Chen Zhi-da
exl)crilrienl in ref. [2], the good results are confirmed by our con'iputhlg progranl. In our computation, the nunlber of elements is only half of ref. [3], 108 increment steps are used for the whole deformatioll process, the maxinmm iterative time in each increment is 3 steps and CPU time Jess lhan 40 minutes in UNIVAC compttier. It is seen from the restllt that the iuteraction between rigkl plate :lnd ring boundary follows the process: contact -~ friction ~slip-~,separation-+. filihlre. This is a very complicated free boundary acting probleln. We have coluputcd the problenl orrubber ring compression hirge defornlatioil in 198,11"i.Separation ol'colll;Ict boundary for rubber ring is not evident as metal one. Geomciric size: External diluiler OD= l l.43ciu, wall Ihickncss t=0.866cnl, thickness in axial direction h = Icnl. Material constant: E = 1890000kg/cnl-' (270(}t)ksi), r = 0.3 t s - r curve is shown in l'ig. 9. mid coiulmting con|]gur~ttions for dilTcrent increnleilt ;ire shown in Fig. 10.
....
,i
. '"...,
l]
5.'.3
~.. 0.V7
Fig. 9
0.04
ut-e
U.06
experinlr
-I 2
~ l t / ~
cuivl~
"',.
"~
r~l
---"
EXl'ohnenl
----('illllllllllll ----it
\,
"
_A
II1
iOII IJl
i
('OIIIllUlalioII {hi Ilii~ i,;ll~cr)
I )iSlqlicl21ill: nl
Pig. l0 Load-radial displaconenl l:iii'vo
Vll.
tl 1.l
0,68
pCki~/i,,)
4
J
t--
to)--
Fig. I I Colilpul;iliOllal COllligUl'alions o[ inctal fillg iii radial COilipressioli
Conclusion
As cOlnl)arcd with TolaI-Lagrangiall or Ul~thitcd-Lagrangian techniques and Green's I]nile Sll'aill del]nition used ill most o f tile largo tlcforniation [ : E M , CO-lllovhlg coordinate tr and
On Unilateral Colll;.Ict l_arge Delbrmalion Problem with Friction (ll)
II
S--R decomposition theorenl ;ire more reasonable and el'fcctive, l-or S-R theorem gives a consistent definition of finite strain and finite rotation, while tile classical definition of Green's strain is not compatiblc. Based on tile new theory and technique, a new tOtal'rating program or nonlinear finite element med~od has beel-~established. Through the typical problems, the rcsl.ltl.',; show highcr accuracy and computing speed may be achived by the new method Ibr solving elasto-i~lastic large deibrmation contact (with fi'iction) problems. Since the computing results are sensitive to the constitutive hLw oJ" the material in the state of large deformation and large rotation, it is still important to examine the constitutive law for large plastic delbrmalion, constitutive relation betweca~ friction and slip, and also the algorithm for variational inequality for frictional contact problems. References
[ I ] Shang Yong and Chen Zhi-da. On unilateral contact large deformation problem with friction, I-incremental variational equations, Appl. Muth. aml Mech.. 10, 12 (1989), 1107-1117. [2 ] Reddy, T. Yella and S. R. Reid, Phenomeml associated with the crushing of met~tl tubes between rigid plates, bit. J. Solids and Struct., 16 (1980), 545- 562. [ 3 ] Taylm, L. M. and E. B. Backer, Some computatiom!l aspects of large deformatimr, ratedep~ Jldct~t plasticity problcms, Conq~. Moths. Appl. Mech. Eng., 41 (1983), 251-278. [ 4 ] Odell, J. T. and G. l:. Carey, Finite Element, VoJ. V, Prentice-Hall, Inc. (1984). [ 5 ] Cheu Zhi-da, Rational Mechanics, China Institute of Mining and Technology Press, Xuzhou, (1988). (in Chinese) [ 6 ] Shang Yong, Theoretical and finite element analysis for elastic and plastic large deformation contact problems with fi'iction, Ph. D. Dissertation, China institute of Mining (1987). (in Chin,.-se) [ 7 ] Sh,'mg Yong and Chen Zhi-da, On Objective stress rate in co-moving coordinate, Appl. Math. aml Mech. 10, 2 (1989), 103- 112, [ 8 ] Neale, K. W., Phenomenological constitutive laws in finite plasticity, SM Archives, 6 (1981), 79 - 129. [ 9 ] Shang Yong and Chen Zhi-da, Analysis of large elastic-plastic deformation by finite element method using co-moving coordinate, Prec. ICNM, Shanghai (Ed by Chien Wei-zang), Science Press (1985), 1294- [301.