<<

. 49
( 54 .)



>>

& stifp)
implicit double precision (a-h,o-z)
dimension elbdb(ndime,nnodb+1,*),lbnod(nnodb,*),




& press(*),x(ndime,*),
& ldgof(ndime,*),tload(*),
& react(*),resid(*),kprof(*),
& stifd(*),stifp(*),dxis(3,2)
data ((dxis(i,j),i=1,3),j=1,2)
& /0.,0.,0.,0.,0.,0./
if(nbpel.eq.0) return
do ie=1,nbpel
epres=press(ie)
do ig=1,ngaub
34 COMPUTER IMPLEMENTATION



For cases where there are surface or line elements with pressure applied,
this routine ¬rst loops over these elements and then over their Gauss points.

bpress segment 2 “ Equivalent forces and sti¬ness component

if(ndime.eq.2) dxis(3,2)=-1.
do id=1,ndime
do jd=1,ndime-1
dxis(id,jd)=0.0
do in=1,nnodb
ip=lbnod(in,ie)
dxis(id,jd)=dxis(id,jd)+x(id,ip)*
& elbdb(jd+1,in,ig)
enddo
enddo
enddo
call pforce(ndime,nnodb,epres,xlamb,elbdb(1,1,ig),
& lbnod(1,ie),ldgof,tload,resid,react,
& dxis)
apres=epres*xlamb
call kpress(ndime,nnodb,apres,elbdb(1,1,ig),
& lbnod(1,ie),ldgof,dxis,kprof,stifd,
& stifp)
enddo
enddo
return
end



For each Gauss point, ¬rst the tangent vectors ‚x/‚ξ and ‚x/‚· are
evaluated. For 2-D cases the vector ‚x/‚· is set to [0, 0, ’1]T , so that when
the cross product with ‚x/‚ξ is taken the correct normal vector is obtained.
These tangent vectors are stored in the array dxis. Using these vectors
and the pressure value epress, the routine pface evaluates and assembles
into tload in a standard manner the equivalent nodal force components
corresponding to this Gauss point using Equation (7.15c). At the same time
these nodal forces are added to the residual force vector and to the reaction
vector. Finally, the routine kpress is called to evaluate the initial pressure
sti¬ness component.
35
8.12 EXAMPLES



y
1
1
1/4 1/4
1/4
3/4
1




1/2 1 x
2
1

FIGURE 8.5 Patch test.


8.12 EXAMPLES

A limited number of examples are described in this section in order to show
the capabilities of FLagSHyP and illustrate some of the di¬culties that may
arise in the analysis of highly nonlinear problems.


8.12.1 SIMPLE PATCH TEST
As a simple check on the code the nonlinear plane strain patch test example
shown in Figure 8.5 is studied. Two irregular six-noded elements making up
a square are employed and boundary displacements are prescribed as shown
in the ¬gure in order to obtain a uniform deformation gradient tensor given
by,
20
F= ; J = 3/2
0 3/4
Assuming Material 1 with values » = 100 and µ = 100, the program gives the
correct state of uniform stress that can be evaluated from Equation (5.29)
as,
µ » 227.03 0 40
b = FFT =
σ= (b ’ I) + (ln J)I = ;
0 ’2.136 0 9/16
J J



8.12.2 NONLINEAR TRUSS
The nonlinear truss example already described in Chapter 1, Section 1.3.2,
can now be re-examined. For this purpose, a single four-noded element is
36 COMPUTER IMPLEMENTATION



99 1
0.4
down
0.3
up
0.2
0.1




F/EA
y
99 0
’0.1
(not to scale) ’0.2
’0.4
1 x
’4 ’3 ’2 ’1 0 1 2 3 4
x/L
(a) (b)

FIGURE 8.6 Truss example: (a) geometry; (b) load-displacement curve.



used with the support conditions shown in Figure 8.6. In order to mini-
mize the two-dimensional e¬ects a large length-to-thickness ratio of approx-
imately 100 is used. The initial angle is 45 √ degrees as in Figure 1.5 of
Chapter 1 and the initial thickness is set to 1/ 2 so that the initial area of
the truss is 1. Material 8 is used with µ = 1/3, which for ν = 0.5 implies
E = 1. The displacement of the top node is prescribed either upward or
downward in small increments of 0.5. The resulting vertical force is shown
in Figure 8.6(b) and is practically identical to that obtained in Section 1.3.2
using the logarithmic strain equation and shown in Figure 1.5.
It is worth noting that when the truss is pushed downward and reaches
the ’45 degree position the program appears to fail to converge. The reason
for this is simply that at this position the stresses in the truss, and hence
all the reactions, are zero to within machine accuracy. Consequently, when
the convergence norm is evaluated by dividing the norm of the residuals by
the norm of the reactions a practically random number is obtained. In order
to escape from this impasse, an arti¬cially large convergence norm must be
chosen for the increment that corresponds to this position (0.98 has been
used by the authors). The computation can then be continued using normal
convergence norms for the rest of the increments. Note that the restarting
facilities in FLagSHyP enable this arti¬ce to be easily performed.


8.12.3 STRIP WITH A HOLE
This is a well-known plane stress hyperelasticity example where an initial
6.5 — 6.5 — 0.079 mm3 strip with a hole 0.5 mm in diameter is stretched in
37
8.12 EXAMPLES



6




Dimensionless force F/Aµ
5
6.5 mm
4
3
F
6.5 mm




F
0.5 mm
2
1
0
1 2 3 4 5 6
Longitudinal stretch »
(b)
(a)




(c)

FIGURE 8.7 Strip with hole: (a) geometry; (b) load-stretch curve; (c) ¬nal mesh.



the horizontal direction and constrained in the vertical direction as shown in
Figure 8.7. Given the two planes of symmetry, 100 four-noded quadrilateral
elements are used to describe a quarter of the problem. A plane stress incom-
pressible neo-Hookean material model is used with µ = 0.4225 N/mm2 . The
strip is stretched to six times its horizontal length using only ¬ve increments.
The load-displacement curve is shown in Figure 8.7(b), and Figure 8.7(c)
displays the ¬nal mesh.


8.12.4 PLANE STRAIN NEARLY INCOMPRESSIBLE STRIP
This well-known example has been included in order to illustrate the di¬cul-
ties that can be encountered using the penalty type of nearly incompressible
plane strain or three-dimensional materials implemented in FLagSHyP in
conjunction with a displacement control process. In this case a 20 — 20 mm2
strip is clamped and stretched as shown in Figure 8.8. Because of the
symmetry only a quarter of the problem is actually modelled using 256
four-noded mean dilatation quadrilateral elements. Material 5 is used with
µ = 0.4225 N/mm2 and κ = 5 N/mm2 . The ¬nal mesh is shown in Fig-
ure 8.8(b), where for a total horizontal stretch of three a vertical stretch of
38 COMPUTER IMPLEMENTATION




20 mm




20 mm
(a)




(b)

FIGURE 8.8 Incompressible strip: (a) geometry; (b) ¬nal mesh.



0.3711 is observed (smaller values are obtained as κ is increased to reach the
incompressible limit.)
Although the material is slightly compressible (ν ≈ 0.46) the solution
requires 200 increments in order to ensure that the Newton“Raphson it-
erative process reaches convergence in a reasonable number of steps. In
fact if the value of κ is increased in order to approach the incompressible
limit (to say 50 or 500), an even larger number of increments is needed.
This is in clear contrast with the previous example, which could be run in
¬ve increments. Unfortunately, using large increments in the present case
leads to large changes in the volumes of those elements containing nodes
with prescribed displacements, which, because of the relatively large bulk
modulus, give extremely large internal pressures at these elements. These
unrealistic pressure values lead in turn to very large residual forces from
which the Newton“Raphson process fails to converge. This problem is typi-
cally solved using the augmented Lagrangian technique whereby equilibrium
is ¬rst reached using small initial values of κ, which are then progressively
increased, always maintaining equilibrium by further iterations if necessary,
until the desired κ/µ ratio is reached. This technique enables very large
¬nal values of κ/µ to be used without increasing the number of increments
needed. In order to keep the code as simple as possible, however, this tech-

<<

. 49
( 54 .)



>>