Treatment of frictionless contact boundaries by direct minimization

TREATMENT

OF FRICTIONLESS CONTACT BY DIRECT MINIMIZATION

BOUNDARIES

A\TU\E F. S.+YECH+and FK;\>‘;K K. Tsot C. F. Braun and Company. Alhambra. C;\. U.S..\.

.Abstract-The class of problems considered involves elastic structures with small deflections. but vrith nonlinear boundary conditions. The boundary conditions are dependent on the loading and. in general. are not known (I priori. Examples of these conditions are resting and lifting,supports and gap boundaries where a boundary point is permitted to displace within prescribed limtts. In this paper. the supports are assumed to be frictionless. Friction will be dealt with in a separate paper. A generalized form of the principle of minimum complementary potential. with a force formulation. is used. It is shown that this direct minimization leads to a quadratic programming problem which is solved using an algorithm described in the paper. Even though the problem area is nonlinear. it is shown that the solution is independent of the loading path. Thus. it is not necessary to incrementally follow the loading history and the solution can instead be obtained in one shot for any prescribed loading. Examples are given to illustrate the nonlinearity of the problem, loading path independence and the solution techniaue. This solution techniaue has been imolemented in Braun’s proprietary computer programs for the analysis of piping systems.

I. ISTRODUCTION

imization of the total complementary potential energy. It has been shown[2] that the principles of total potential and total complementary potential of the theory of elasticity can be extended to cover these nonlinear boundary conditions. First, both minimum principles will be reviewed. Next, the force approach will be used to formulate the problem using the principle of minimum complementary potential. An algorithm to solve the resulting quadratic programming problem is described next. Finally, numerical examples are given to show the nonlinearity and loading path independence of this class of problems. This solution method has been successfully implemented in C. F. Braun’s computer program for the analysis of piping systems. This program handles a wide variety of boundary conditions including gaps and friction. It contains finite elements as well as pipes and elbows.

In this paper we are dealing with structures made of linear. elastic materials, with small deflections, but with special boundary conditions. Looking at Fig. I. the supports at points 3 and 4 are regular fixed supports. However, the supports at points I, 2, and 5 are different. At point I, the support is such that the displacement ([I CC6,. The prescribed gap 6, can be positive, negative, or zero. If it turns out that. if for the given loading 1/I < 6,, then the support exerts no reaction force. Nonzero reaction is possible only when M, = 6,. and when this occurs, it would be acting downwards. At point 2. the displacement ~12 8:. A negative value for S2is shown. but it can be +. -, or zero. Thus the support at Z can exert only positive (up) reaction when 11: = 8:. or be force-free uhen 11: > 6:. An example of this type of support is the resting support where the structure rests on a “block” but is not restrained from moving upwards. Nonzero 62 can represent a settlement of the support at 2. 4t point 5. we have a double acting gap such that 6, < lls G 6,. The reaction force can now be + , -, or zero. Friction at the supports is not considered in this paper. but will be the subject ofa subsequent paper. These types of supports occur in the analysis of piping systems and, up to now. most piping analysis computer programs have used a trial and error approach to iteratively solve this problem[ I]. Besides being costly, the trial and error approach does not always converge[Z]. The approach used in this paper is the direct min-

2. hllNI\ttN Genern/i:ed

principle

PRINCtPLES[3. JJ

of minimum

fotol

potenticd

This principle states that among all kinematically admissible displacement fields II,! (i = I, 2, 3) which do not violate the boundary conditions, the exact solution lli is the one which reduces the total potential energy fl; to an absolute minimum. For an elastic body with volume V. with SF being the portion of the surface on which the forces F; (i = I. 2. 3) are prescribed, the total potential corresponding to II,! is defined as

t Principal Engineer, Computer Systems Department. $ Section Leader, Autodrafting Computer Systems Department.

and Stress Section,

where u,;, E,;. are the stresses and strains derived 39

at point 2. the contribution 6:. anJ Fz 2 U. at point

is F;

Fig.

--+

The

principle

tential from

rc;“. And

pressed

the

minimum

principle

can

be ex-

the force lem

in

terms

potential

I(,. To be kinematically field

admissible.

II,* has to satisfy

nonlinear ferring

of the exact

boundaries

certain

solution

a displacement

requirements

are present.

when

For esample,

re-

I,

to Fig.

FOR\lLL~\TIOY

is used

of

order

stress

field.

structure

to

I will

Step I: Remove to form

structure

supports the

In

removal

This principle stress

boundary

states that among fields

I$

conditions,

one which

reduces

the total

For

an elastic the portion

placements complementary fined

the external

body

complementary

with

volume

C’. with

corresponding

E$ are the strains

corresponding

the corresponding

surface

stable

boundaries. is shown

The

in Fig. 2(a). boundaries,

be

SL:

the disthe total

to ui\ is de-

the minimum

to a:‘, and P,;‘,’are

as shown

in Fig.

Step 3: Apply

2(a).

the unknown

Fj. one at a time, structure

of course.

principle

[(I].

corresponding I,

at

can be expressed

number

them

that

systems.

is the total

complementary

o,,.

Nonzero

cluded

in the second

integral

if they

were

displacements.

imposed

be statically

admissible.

isfy certain

requirements referring

potential

gaps should over

In addition.

to Fig.

I, the contribution

coefficients. in a reduced

flexibility

matrix

since

only few boundaries the strain

to the

conditions. even

This

for large

are usually

is

piping

nonlinear.

energy

U = I/” + ; {F)T[~~] {f)

+ {o}“’ {F}.

(8)

to

boundaries. I.

where

u? is the strain

under

the external

have

x 6, and F, s 0,

to the nonlin-

c$ has to sat-

for nonlinear again

values

of

be in-

SL. in eqn t-l) as

a stress field

unit

Z(b)-(d).

2, and 5 due to unit loads at namely (l,, ti.j = I. 2. 5) are.

boundary

property,

Srep 4: Compute

F2. the

the size of [cl] corresponds

of nonlinear

an important

in Figs.

F,,

Since

forces

vre can apply

as shown

the flexibility

Notice

reaction

on the basic structure.

is linear.

Let us group

tractions.

as

is F,

problems.

the structure

at the nonlinear

the same boundaries.

at point

practical

nonlinear

loads,

The displacements

For example,

most

the apthe basic

poten-

on which

I. 3, 3) are specified.

potential

solution

that

minimum.

ear boundaries

n,

under

assumes

the

basic

the exact

structure

in our example

of the unknowns,

where

conditions

o,, is the

as

Thus.

The

ad-

violate

solution

of the surface

II, (i =

all statically

do not

the exact

tial 112 to an absolute being

which

boundary

to leave

of

Let the displacements

missible

admissible

are followed.

be used as an example.

This

exist

basic structure under

steps

the basic loads.

is stable.

enough

in for

structure.

Srep 1: Analyze external

than

advantages

a statically

the nonlinear

the basic

the prob-

friction.

construct

in Fig.

Thus

a statically

rather

has definite

the following

po-

paper.

Casting

forces.

the case of nonzero

In

after

field.

unknown

of displacements.

plied

complementary

to construct

stress (or force)

terms

to

5;.

as the basis of this

method

treating ITA is the total

x

of minimum

is selected

admissible

as

where

is Fi

3. FORCE

u5 ‘F5

I.

j. 17,

.G 6; and for the range F, s 0.

the contribution

!e

is fZ

5. for the range Fc 2 0. the contribution

the integral

I

to the integrsl

to the integral

over

energy

loads.

to be computed.

of the basic structure

which

actually

does

and {F) is the vector

SC,

known

reactions

(6)

Srep 5:

Add the complementary

at the nonlinear

not

of un-

boundaries.

potential

for non-

Treatment

of frictionless

contact

boundaricz

(0)

(cl

(b)

(d)

Fig. 2. (a) Basic structure;

zero gaps to form the total complementary potential rs:. . At point I the complementary potential is -F, 6,. At point 2 the complementary potential is - F2 6:. At point 5, it is - F>S; when F5 5 0 and is -F5 ST;when Fs 2 0.

(b); (c); (d).

Srep 6: Formulate

the constraints.

I

F, s 0.

at point Z

Fz a 0.

at point 5

F5 2 0 1 F5 s 0

at point

For example,

(IO) when 6, is used. when S,L is used.

Thus. in general

n!; = I/” + ; {F}‘[u]{F} =

U= + ({u)”

-

{6})‘{F)

+ {u}“‘{F}

- {S}‘{F}

+ ?{F}‘[o]{F}.

(9)

where {S} is the appropriate vector of gaps corresponding to {F}. We see that the first term on the right-hand side of (9) is a constant. the second term is linear in (F}. and the third term is the quadratic term. In addition. because of the properties of the flexibility matrix [a]. we see that II!; is a positive definite quadratic form of the unknowns {F}.

We see that the constraints are linear in the unknowns {F}. Srep 7: Perform the minimization of HI{, eqn (9). subject to the constraints t IO). Since the function IIS is quadratic and the inequality constraints are linear. this minimization problem is a quadratic programming problem. An algorithm to solve this problem will be described later. The solution {F} is the one corresponding to the absolute minimum of II; without violating the constraints. Step 8: After the unknowns {F) are found, we apply a simple superposition of the results of the basic structure under external loads and the basic structure under the forces {F}.

In our example.

the final

result\

get

will

correspond

we

f-i

< Fig.

The algorithm

to

Fig.

‘1~11 -

F,

A Fig.

Irb,

F:

-

x Fig. 2ic)

-

2tdt.

The

in our

.l’orc, The

analysis

the external

of the basic

Braun‘s gram

reactions.

general

has finite

purpose

elements

is based

Thus.

while

force

lowed

of course.

boundaries blocks

case

Thus.

by a

the

each

step.

system

(frame)

including

mod-

programming

problem

This

is

the

method.

where

,=I

N;,_K,.K,

c

+

j=l

that violate

inequality

J’

:‘; reaches

Thus

zero.

1.~ is com-

or (16)

to

AXJV_k

increment

any AX,.

=

Gauss-Seidel

vve have one extra

the

in state

-).I,

iterative requirement,

given

of the constraints vve see from

variables

by

(14).

eqn

(13)

(16)

Correthat

the

is k=

C:,l..Ui,

I.2 ,...(

M.

(17)

minthe increment variable

(I 1)

hrxr ,

Axi=

AX,. which

equal

to zero

_-%ZL, (‘i,

I=,

~,i and hi are constants

N variables

derivative

to increase.

in

by 11,.

will

make

is given

the Xth

by

.\

.v

.L

c c

is reduced

incremented

\vell-known

namely.

state ;

incremented,

the function

Y, to reach

However.

doesn’t

Thus,

=

are

AX, = -I-! i’,. Cl,i

as formu-

as follows:

imize

?’

start

by requiring

change

1. QllADR.ATIC PROCR.ASI~IISG

lated in this paper can be expressed

that

( I-l).

condition

v,ariables

as

are fol-

continua

tl, in eqn (13) are zero. satisfies

the decision

of the solu-

elements.

quadratic

until

then it will

puted

sponding

The

such

point.

.V is used. since

W’hen X, is being

will decrease zero,

point

decision

one at a time.

at any’ feasible

I. 2.

the constant5

this starting

formulation.

the steps outlined

structure.

pro-

can be started

.I-, = 0. i =

Nest.

which

This

method.

we used a piping

any

eled by finite

program.

a stiffness

to illustrate,

both

by, using

elements.

the building using

though for

computer

the nonlinear

formulation.

an example

is performed

on the stiffness

treating

under

loads corresponding

and piping

are obtained Even

structure

loads and the unit

to the unknown

tion

point

.c, are decision

and (I,,

variables.

=

cl,,. The

Subject

to the

constraints

The correct ues given used

value

I,2

,...,

for

M.

for As, is the smallest

( 16) and (18). Equation

by

only

X=

state

variables

(18)

of the val-

(18) should

which

are also

be

slack

variables.

2

C,,.K,

a

fl,

+

i =

0.

, ,\I.

I, 2.

(17)

,=I

When

an increment

has caused where

( 12) represents

inequalities ties

a slack

L’,, and (1, are constants.

Thus

M linear

can instead

be expressed

with

nonnegativity

together

by introducing

M

M slack

variables

variable,

to become

These

ensure

that further

as ;\I equali-

iables,

do not force

inequalities.

conditions

i = I, 2.

I~_,.

. ,\I.

do that

X, and I~_,,.

C,,.\; + Cl,.

i =

, 2.

I

.

i = I. 2.

.t-,y_; 3 0.

differential

plicity

algorithm[j]

programming

and since

Gauss-Seidel

is used

problem

it is basically

method

ear equations,

with

( 13)

. .\I.

(14)

to solve

because

to become thus making

var-

negative. variable

To

by in-

X, a state var-

most

LLr,v_O

step is to write

I, 2.

engineers

. A’. By differentiating

cp,l.r,

+ CP,A.li

this

of the

simultaneous

the decision

= 2 /fi

we can write

or

linare fa-

Since .Y,~_,, is considered replacing formation

The first

form,

of its sim-

an extension

for solving which

, .I/.

miliar. i’;, i =

we need to

iable.

_\‘.VL, = ,F,

quadratic

then

X,

is also

in other decision

.c,\. _,, a decision

In an incremental

The

I~_;,

variable

.Y,~_,,, which

zero,

increments

vve make

terchanging

A.\-, in decision

the pth state variable

derivatives eqn

( I I ).

cision

the old x,. of variables

variables

variables

as a new decision

we can write between

variable

( 19) as a trans-

the set of old

de-

_L{_T},~~and the set of new decision

_I{_Y}“~~. as follows:

Treatment

of frictionless

contact

43

boundaries

I 1

I --CpI

i-

ICE2

...

1

-

-- C&V

Cpi

Cpi

Cpi

Cpi

..,

I

or let us write

it as

1 {.\.Lld = i-

Let change

us express of state

the

state

.A{.\.jnc,, = [Ply {.\.lncm.

‘..

!

PA PC ...

p,;

‘..;,

’

variables

set in an incremental

before form

-l{i’)

the

After

the change

straint

=

= [a] -1(X) (,,,,.

(24)

as After

l{.rI

(20)

I

p,,

the change

of the state

set. we have

(21)

la1 ~{.hl.

of state set. the equations

-l{?}

of con-

= [O’jl(.\.}“,,.

(25)

become where -I(s)

= [al

[PI A{.V)ns,r

= [a’]

-l(~r)nCu.

(??)

[a’] = [PIT [n] [PI = [PIT [Cl

where

ith column (~II

+ al, Pi11

(aI2 + al,P,2) ..' (ali + a:;P,I) (a:? + a,,p,z) ... =

:

pth -+

[a’]

P;,

p,:

ah P,, a3 P!,

“’

If we save the ith column

of [a] in a working

. a.viJthen the computations

al;,

(23) can be performed

column-wise.

When

array

..,

(Cl1

in eqn

(C2l

=

working

ith --j row

(23)

+ I: p,,

row

[a,;.

..' (a1.v+ al,P,.v) ..' (a?, + a2, P,.v)

+

Pi\

p,,c,,, cc,: + P,,CiJJ

“’

P,2C;,) (C2 + P,2C;2)

“’

a,,>;:

P.&

. (26)

.”

where

.

I

lth column 1 1 ((I,,

on the kth besides

column.

this working

+ n,,P,,J

the only array.

((111 +

information

Cll,P,2)

needed.

...

Thus.

Ul,Pi,

two

“’

+

tr,,p,,,

If we save the ith column In an incremental

form.

we write

eqn

( 15) as

be performed

to compute

[u’].

[C] and the second to compute

[Cl,,. (12;. . . . . c~.~;J. then

Deci.sion deri\~rtril~r.s

1

passes are required

first to compute

is P,r.

((I,\

columnwise.

of [(I] in a working the operations Also.

since

the [(I’]. array

in (27) can [pi,.

t3,2.

. p,\] is av,ailabie as a working array. the opin (26) can be performed columnwise.

erations

To handle problems of any size and keep generality. it was decided to store the matrix of coefficients of decision derivatives [n] and the matrix of coefficients of the constraints [a] on a disk file. It was found to be advantageous to store these matrices columnwise as shown in Fig. 3. The record size is selected as a multiple of (N + M) and can be as small or as large as necessary to tit the available core. When ~~~ is being incremented, all the coefficients needed to check decision derivatives, as well as state variables. are in the ith column. Thus, we proceed sequentially, column by column. A change of state set, when encountered, is handled as described earlier. Convergence

Convergence is achieved when all decision variables satisfy one of the following: 1. A decision variable is a slack variable which is zero and the corresponding decision derivative 2 0. 2. A decision variable is a slack variable > 0 and the corresponding decision derivative = 0. 3. A decision variable is not a slack variable and the corresponding derivative = 0. Thus, we see that for all decisions i, UiXi = 0. 5.

NUMERICAL EXAMPLES

Two examples are selected to illustrate the characteristics of this class of uroblems, namely: Conservatism. In the absence of friction, there is no energy loss. Path independence and reversibility. These follow from (a). The sequence of loading is immaterial, and what determines the state of deformation and stresses is the final loading configuration. Nonlinearity. Even for linear, elastic materials, the response can be highly nonlinear. Direct solution. Since the problem is nonlinear and elastic, we do not have to follow an incremental approach. For a given load, the solution is found in one shot. The results for both examples were obtained using

C. F. Braun’s proprietary computer program for the analysis of piping systems. It contains the solution technique outlined in this paper. Both examples involve piping systems. since in such systems. nonlinear boundaries such as lifting supports, gaps. and friction are quite common. Exumple 1

The piping system shown in Fig. 4 is selected with the following properties: E = 30 x lo6 psi.

Thermal

expansion

= 0.6 in/l00 ft.

The supports at 1 IO, 120, . , 180. 200 allow only upward motion and the support at 190 allows only downward motion. All supports have zero gaps. First, dead load and thermal were applied simultaneously. The computer program found the corresponding solution in one shot. The reactions are as shown in Fig. 5. We see that only the supports at 180 and 200 are active. Instead of obtaining the solution in one shot, we can proceed incrementally. For example, let us apply the full thermal load first, and then apply the dead load incrementally. The sequence of reactions corresponding to different levels of dead load will be as shown in Fig. 6. As seen, the final result is the same as the one obtained in one shot. Configuration 6(b) is the one when support 190 just becomes inactive; 6(c) is the one where support 180 is just becoming active. If we reverse the sequence of loading, such that the full dead load is applied first and then thermal ioatllt~ t;:-+:*::*+i lcm110 120 i3J

1011

14015016017u

480 190 x0 .

= 0

Fig. 4. Piping system.

Fig. 3. Columnwise

storage

of matrices.

Fig. 5. Simultaneous

application

of dead load and thermal

Treatment of frictionless

contact boundaries

45

seen, the final result of both load sequences is the same. l Two quantities are selected to illustrate the non2lb4 IQ1 Thrmol only linearity of the problem. Figure 8 shows the history t 1 of the horizontal displacement at point 100 and Fig. (b) Tturma~ + 892 D L 9 shows the history of the vertical reaction at the ? fixed base 50. 1 (c) Thermal + ,945 0 L 849 Looking at Fig. 9, we see one consequence of ? the nonlinearity. When thermal was applied by itt IdI Thwmal+D L 94 828 self, it produced a reaction of 316 Ibs. However, Fig. cj. The sequence of reactions corresponding to dif- when the same thermal was applied after the apferent levels of dead load. plication of D.L.. it produced a reaction of 1140 95 = 1045 Ibs. Also, when D.L. was applied by itself, it produced a reaction of 95 Ibs but the same 4 ?tt+T++T D.L. applied after the application of thermal pro283tR2022Ot20t199173362 175 D.L only duces 1140 - 316 = 824 Ibs. almost 9 times as ‘) 196173 362 175 DL + 0047T 28292l52ce much. Obviously. the principle of superposition can 4 + + not be applied to such problems. I 362 1 DL + 015T 366.2 220 When unloading takes place, it follows the same I +++a+ l 1 ti.6 694 I tj3.41 kl DL + 034T loading path. For example, if after applying both 2ot.3 3621 r l *** D.L. and thermal, we want to remove thermal, then 1 5ti.6 2%5~4;~~3k2 1 1is DL + 0688T path (B) will be followed. On the other hand, if D.L. 1 + + 4 were removed, then path (A) would be followed. 593 5 2t9 1753 DL + 1728T 73 A third loading path (C) is shown in Figs. 8 and +a2Q ? I 64a5 9 17t35 DL + 3048T 9. It is intended to demonstrate the case of nonzero 1 ? + gaps. All the properties already discussed are still 9211 2066 DL 4894T valid. The following nonzero gaps were used, which ? ? L could have resulted from settlement of supports. 94 628 DL +T 241

l

Fig. 7. The sequence of reactions with reverse sequence of loading.

Point

is applied incrementally, we would get the sequence shown in Fig. 7. Each step in Fig. 7 represents one more of the nonlinear supports becoming inactive. As can be

120 130 140 150 160

Hor~ronlol

dlsplocement

al

Gap (inches) -

0.0004 0.0008 0.00096 0.0008 0.0004

100

Fig. 8. The history of the horizontal displacement

at point 100

46

Vertical

reactton

ol 50, lb

Fig. 9. The history of the vertical reaction at the fixed base SO.

Example

2

Usually, when dealing with proportional loading, we do not expect unloading. For example, we do not expect the displacement at a point to initially increase and then, while the load is still increasing, to start decreasing. When nonzero gaps are present, this behavior can occur, even for proportional loading. Actually, nonzero gaps make the loading nonproportional. This example illustrates this phenomenon. The method of solution of this paper is not restricted to proportional loading, and in the absence of friction, the solution can still be obtained in one shot. The piping system shown in Fig. 10 is similar to

Fig. 10. The piping system

that in Example I; however, we now have only three nonlinear supports at points 110, 170, and 200, with nonzero gaps at 110 and 200. Supports at 110 and 200 can exert downward reaction only and that at 170 can exert only upward reaction. Full thermal loading corresponds to an expansion of 10 in/100 ft. and only thermal loading is applied.

E = 30

x

IO6psi.

Various percentages of the full thermal load were applied. and in each case, the computer program found the solution in one shot. The history of the vertical deflection at point 200 is shown in Fig. 11. We see that the deflection increases linearly until the gap at 200 closes (segment l-2). The gap remains closed while the load increases from 10% to about 42.5% (segment 2-3). Now, with further increase of the load. we see that the displacement decreases and even reverses direction (segment 3-1) until it reaches the extreme value at about 47.5%. With the load still increasing, upward deflection resumes (segment 4-5) until the gap is closed for the second time (point .5), at about 84%. The gap remains closed for the remainder of the load (segment 5-6). It is seen that the overall force-deflection curve is highly nonlinear. even though each segment is linear. Figure I? shows the vertical deflection at another point, 140. An example of force history is shown in Fig. 13 where the vertical reaction at the fixed base, 50, is plotted. We see the the first 40%

Treatment of frictionless

Vertical deflection 01 pt. 200,

contact boundaries

47

I”.

Fig. I I. The history of the vertical deflection at point 200.

of the load produced about 0.1% of the final reaction. The gradual deformation of the top (point lOO200) is shown in Fig. 14, where the various percentages of the full load are indicated on the curves. We see that the gap at 200 closes first, to be followed by the closing of the gap at 110. As soon as the gap at 110 closes, the gap at 200 opens and the segment 1 lo-200 starts to move downwards instead

N~rtrol

of upwards. This continues until the support at 170 becomes active, at which point the segment 170200 starts to move up again, until the gap at 200 closes for the second time.

6. SUMMARY

AND CONCLUSIONS

The problem of nonlinear boundaries such as resting and lifting supports and gap boundaries was

deflection

at

pt. WO, In

Fig. 12. The vertical deflection at point 140.

ANTONE

40

60

F.

and FR.ANK K. Tso

SAYEGH

a0

al

t20

MO

200

Vertical reaction 01 pt. 50, kips

Fig. 13. An example of force history. For an elastic material, and in the absence of friction, it was shown that this is a nonlinear, elastic problem which is path independent. The generalized principle of minimum complementary potential was the basis for a force formulation of a direct solution method. For a linear elastic material, this leads to a quadratic programming problem, which was solved using a differential algorithm which was described. Numerical examples were given to illustrate the treated.

capabilities of this direct solution method and the characteristics of this class of problems. The implementation of this solution technique in C. F. Braun’s computer program for the analysis of piping systems and its successful use for numerous problems establishes both the feasibility and the efficiency of the method. Acknowledgement-The authors wish to express their gratitude to K. H. Pang, R. C. Harbage, and R. B. Hill of

Fig. 14. The gradual deformation of the top (point 100-21X1), where the various percentages load are indicated on the curves.

of the full

Treatment of frictionless C. F. Braun and Company for their interest and encouragement in planning the development of this investigation; to S. G. Sharma and the members of the Stress Section, Piping Department of C. F. Braun and Company for pointing out the need for this development; and to the management of C. F. Braun and Company for providing the means for conducting the study.

REFERENCES 1.

B. Torstenfelt, Contact problems with friction in general purpose finite element computer programs. Comput. Srruct. 16, 487-493 (1983). 2. A. F. Sayegh, Elastic analysis with indeterminate

contact boundaries

49

boundary conditions. Proc. AXE, Engng Mech. Div. lOO(E!vll), 49-62 (1974). 3. P. G. Hodge. “Numerical Applications of Minimum Principles in Plasticity”, Engineering Plnsricity. Papers for a conference held in Cambridge, LMarch 1968. Cambridge U. P., Cambridge (1%8), pp. 237-256. Variational Principles in the Mathe4. D. C. Drucker, matical Theory of Plasticity, Proceedings of the Eighth Symposium in Applied Mathematics American Mothemnrical Society. Held April

of the

12-13, 1956 at the University of Chicago. McGraw-Hill, New York, (1958). pp. 1l-12. 5. D. J. Wilde and C. S. Beightler, Foundutions of Oprimiznrion. Prentice Hall. Englewood Cliffs, N.J. (1967). pp. 69-81.