Civil Engineering Reference
In-Depth Information
and subroutine
shape_der
delivers the derivatives with respect to
L
1
and
L
2
1
−
10
der
=
(3.57)
0
−
11
The nodal numbering is shown in Figure 3.9(b). Exactly the same sequence of oper-
ations (3.47) to (3.51) as was used for quadrilaterals places the required derivatives with
respect to
(x, y)
in
deriv
, finds the Jacobian determinant
det
, forms the
bee
matrix
and numerically integrates the terms of the stiffness matrix
km
. For this simple element,
only one integrating point at the element centroid is required (
nip=1
). For higher-order
elements more triangular integrating points would be required. For example, the 6-node
triangle would usually require
nip=3
for plane analysis. For integration over triangles,
the sampling points in local coordinates
(L
1
,L
2
)
are held in the array
points
and the
corresponding weighting coefficients in the array
weights
. As with quadrilaterals, both
of these items are provided by the subroutine
sample
. This subroutine allows the total
number of integrating points (
nip
) for triangles to take the values, 1, 3, 6, 7, 12, or 16. The
coding should be referred to in order to determine the sequence in which the integrating
points are sampled for
nip>1
.
3.7.5 Axisymmetric strain of elastic solids
The formation of the strain-displacement matrix follows a similar course to that described by
(3.47) to (3.49), however in this case
bee
must be augmented by a fourth row corresponding
to the “hoop” strain
θ
as shown in (2.76). The cylindrical coordinates
(r, z)
replace their
counterparts
(x, y)
. The stress-strain matrix is given by equation (2.77) and is returned
by subroutine
deemat
with
nst
, the number of stress and strain components now set
to 4.
In this case, the integrated element stiffness is given by (2.74), namely
[
B
]
T
[
D
][
B
]
r
d
r
d
z
[
k
m
]
=
(3.58)
where
r
is the radial coordinate given in the programs as
gc(1)
from the isoparametric
relationship,
(3.59)
gc = MATMUL(coord,fun)
where
gc
and hence
fun
are evaluated at the sampling points.
The numerical integration summation in axisymmetry is written as
nip
(3.60)
km =
det
i
*weights(i)*btdb
i
*gc(1)
i
i
=
1
By comparison with (3.51) it may be seen that when evaluated numerically, the algo-
rithms for axisymmetric and plane stiffness formation will be essentially the same, despite
the fact that they are algebraically quite different. This is very significant from the points
of view of programming effort and of program flexibility (e.g. Program 5.1).
