<<

. 41
( 54 .)



>>

K p,23 = ’ ° 0 0 4»
2
’1 ’4 0




7.4.5 TANGENT MATRIX
The linearized virtual work Equation (7.30) can now be discretized for el-
ement (e) linking nodes a and b (see Equation (7.32) and Figure 7.2(a)),
18 DISCRETIZATION AND SOLUTION



δva
δva δva
ub ub


a
a a
b b


(b) (c)
(a)

FIGURE 7.2 Assembly of linearized virtual work



in terms of the total substi¬ness matrix K ab obtained by combining Equa-
tions (7.36), (7.45), and (7.50) to give,
(e)
DδW (e) (φ, Na δv a )[Nb ub ] = δv a · K ab ub ;
(e) (e) (e) (e)
K ab = K c,ab + K σ,ab ’ K p,ab (7.51a,b)

The assembly of the total linearized virtual work can now be accomplished
by establishing (i) the contribution to node a from node b associated with
all elements (e) (1 to ma,b ) containing nodes a and b (see Figure 7.2(b)), (ii)
summing these contributions to node a from all nodes b = 1, na , where na is
the number of nodes connected to node a (see Figure 7.2(c)), (iii) summing
contributions from all nodes a = 1, N . This assembly process is summarized
as,
ma,b
DδW (e) (φ, Na δv a )[Nb ub ] (7.52a)
(i) DδW (φ, Na δv a )[Nb ub ] =
e=1
e a,b
na
(ii) DδW (φ, Na δv a )[u] = DδW (φ, Na δv a )[Nb ub ] (7.52b)
b=1
N
(iii) DδW (φ, δv)[u] = DδW (φ, Na δv a )[u] (7.52c)
a=1

This standard ¬nite element assembly procedure can alternatively be
expressed using the complete virtual velocity vector given in Equation (7.20)
together with the corresponding nodal displacements, uT = [uT , uT , . . . , uT ]
1 2 n
and the assembled tangent sti¬ness matrix K to yield,

DδW (φ, δv)[u] = δvT Ku (7.53a)

where the tangent sti¬ness matrix K is de¬ned by assembling the nodal
19
7.5 MEAN DILATATION METHOD FOR INCOMPRESSIBILITY



components as,

···
K 11 K 12 K 1n
® 

···
K K 22 K 2n 
 21
K= .

. .
..
°. . .»
.
. . .
···
K n1 K n2 K nn



7.5 MEAN DILATATION METHOD FOR
INCOMPRESSIBILITY

The standard discretization presented above is unfortunately not applica-
ble to simulations involving incompressible or nearly incompressible three-
dimensional or plain strain behavior. It is well known that without further
development such a formulation is kinematically overconstrained, resulting
in the oversti¬ phenomenon known as volumetric locking. These de¬cien-
cies in the standard formulation can be overcome using the three-¬eld Hu-
Washizu variational approach together with an appropriate distinction being
made between the discretization of distortional and volumetric components.
¯
The resulting independent volumetric variables p and J can now be interpo-
lated either continuously or discontinuously across element boundaries. In
the former case new nodal unknowns are introduced into the ¬nal solution
process, which leads to a cumbersome formulation that will not be pursued
¯
herein. In the latter case the volumetric variables p and J pertain only to an
element and can be eliminated at the element level. In such a situation the
¯
simplest discontinuous interpolation is to make p and J constant through-
out the element. This is the so-called mean dilatation technique discussed
in Section 6.6.5. Observe, however, that for simple constant stress elements
such as the linear triangle and tetrahedron the mean dilatation method co-
incides with the standard formulation and therefore su¬ers the detrimental
locking phenomenon.


7.5.1 IMPLEMENTATION OF THE MEAN DILATATION
METHOD
Recall from Chapter 6, Section 6.6.5, that the mean dilatation approach for
a given volume v leads to a constant pressure over the volume, as indicated
by Equations (6.50a“b). When this formulation is applied to each element
(e) in a ¬nite element mesh the pressure becomes constant over the element
volume. In particular, assuming for instance that the potential shown in
20 DISCRETIZATION AND SOLUTION



Equation (6.51) is used, the uniform element pressure is given as,
v (e) ’ V (e)
(e)
p =κ (7.54)
V (e)
where V (e) and v (e) are the initial and current element volumes.
The internal equivalent nodal forces for a typical element (e) are given
by Equation (7.15b) where now the Cauchy stress is evaluated from,
σ = σ + p(e) I (7.55)
and the deviatoric stress σ is evaluated using the appropriate constitutive
equation given by (5.51), or (5.103).
Continuing with the discretization, recall the modi¬ed linearized virtual
work Equation (6.60), which for an element (e) is,
(e)
σ : [( u)T
DδWint (φ, δv)[u] = δd : c : µ dv + δv] dv
v (e) v (e)

κv (e)
(e)
+ κv (·u)(·δv);
¯ κ = (e)
¯ (7.56)
V
ˆ ˆ
where the elasticity tensor is c = c +c p and c is the distortional component
that depends upon the material used, and c p is given by Equation (5.55b)
as,
c p = p(I — I ’ 2i ) (7.57)
The average divergences are now rede¬ned for an element (e) as,
n
1 1
·u = ·u dv = Na dv (7.7.58a,b)
ua ·
v (e) v (e)
v (e) v (e) a=1

n
1 1
·δv = ·δv dv = δv a · Na dv
(7.7.58c,d)
v (e) v (e)
(e) (e)
v v a=1

Discretization of the ¬rst two terms in Equation (7.56) is precisely as given
in the previous section, but the ¬nal dilatation term needs further atten-
tion. For element (e) the contribution to the linearized internal virtual work
related to the dilatation and associated, as before, with nodes a and b is,
(e)
DδWκ (φ, Na δv a )[Nb ub ]
κ
= (e) δv a · Na dv Nb dv
ub ·
V v (e) v (e)
κ
= δv a · Na dv — Nb dv ub
V (e) v (e) v (e)
(e)
= δv a · K κ,ab ub (7.59)
21
7.6 NEWTON“RAPHSON ITERATION AND SOLUTION PROCEDURE



where the dilatational tangent sti¬ness component is obtained in terms of
the average Cartesian derivatives of the shape functions* as,
1
(e)
K κ,ab = κv (e)
¯ Na — Nb ; Na = Na dv (7.60)
v (e) v (e)
The complete discretization of Equation (7.56) can now be written in terms
of the total element tangent substi¬ness matrix as,
(e) (e) (e) (e)
K ab = K c,ab + K σ,ab + K κ,ab (7.61)
(e)
into which the surface pressure component K p,ab (see Equation (7.51)), may,
if appropriate, be included. Assembly of the complete linearized virtual
work and hence the tangent matrix follows the procedure given in Equa-
tions (7.52“3).


7.6 NEWTON“RAPHSON ITERATION AND SOLUTION
PROCEDURE

7.6.1 NEWTON“RAPHSON SOLUTION ALGORITHM
In the previous sections it was shown that the equilibrium equation was dis-
cretized as δW (φ, δv) = δvT R, whereas the linearized virtual work term is
expressed in terms of the tangent matrix as DδW (φ, δv)[u] = δvT Ku. Con-
sequently, the Newton“Raphson equation δW (δv, φk )+DδW (φk , δv)[u] = 0
given in Equation (6.2) is expressed in a discretized form as,
δvT Ku = ’δvT R (7.62)
Because the nodal virtual velocities are arbitrary, a discretized Newton“
Raphson scheme is formulated as,
Ku = ’R(xk ); xk+1 = xk + u (7.63)
Although it is theoretically possible to achieve a direct solution for a
given load case, it is however more practical to consider the external load F
as being applied in a series of increments as,
l
F= ∆Fi (7.64)
i=1
where l is the total number of load increments. Clearly, the more increments
taken, the easier it becomes to ¬nd a converged solution for each individual

* The emergence of (v (e) )2 and the average Cartesian derivatives in Equation (7.60) is simply to
conform to expressions occurring in the literature.
22 DISCRETIZATION AND SOLUTION




BOX 7.1: Solution algorithm
r INPUT geometry, material properties and solution parameters
r INITIALIZE F = 0, x = X, R = 0
r LOOP over load increments
r FIND ∆F using (7.15c)
r SET F = F + ∆F
r SET R = R ’ ∆F
r DO WHILE ( R / F > tolerance )
FIND K using (7.51b)
r
SOLVE Ku = ’R
r
UPDATE x = x + u
r
FIND F (e) , b(e) and σ (e) using (7.5), (7.9d) and typically (5.29)
r
FIND T using (7.15b)
r
FIND R = T ’ F
r
r ENDDO
r ENDLOOP




load step. Observe that in the case of a hyperelastic material the ¬nal

<<

. 41
( 54 .)



>>