Civil Engineering Reference
In-Depth Information
and since the Jacobian matrix is given by
∂y
∂ξ
∂x
∂ξ
[
J
]
=
(3.46)
∂x
∂η
∂y
∂η
it is clear that its terms can be obtained from (3.45), once the derivatives of the shape func-
tions with respect to the local coordinates have been provided by subroutine
shape_der
,
thus
CALL shape_der(der,points,i)
(3.47)
jac = MATMUL(der,coord)
det = determinant(jac)
The function
determinant
computes
det
, the determinant of the Jacobian matrix,
required later for the purposes of numerical integration.
In order to compute
deriv
we must invert
jac
using subroutine
invert
and finally
carry out the multiplication of this inverse by
der
to give
deriv
,
CALL invert(jac)
(3.48)
deriv = MATMUL(jac,der)
It should be noted that subroutine
invert
overwrites the original matrix by its inverse,
thus
jac
in fact holds
[
J
]
−
1
after the subroutine call.
The matrix
[
B
]
in (3.40) (called
bee
in program terminology) can now be assembled
as it consists of components of
deriv
, and this operation is performed by the call,
(3.49)
CALL beemat(bee,deriv)
The components of the integral of [
B
]
T
[
D
][
B
], at each of the
nip
integrating points, can
now be computed by transposing
bee
using the Fortran 95 intrinsic function
TRANSPOSE
,
and by forming the stress-strain matrix
dee
using the subroutine
deemat
.In2Danalysis
(ndim=2)
, subroutine
deemat
gives the plane strain stress-strain matrix (2.70). The
size of
dee
is given by
nst
, the number of components of stress and strain, which for 2D
elastic analysis is equal to 3. A plane stress analysis would be obtained by simply replacing
subroutine
deemat
by
fmdsig
.
The multiplication
(3.50)
btdb=MATMUL(MATMUL(TRANSPOSE(bee),dee),bee)
gives
btdb
, the quantity to be integrated numerically by
nip
(3.51)
km =
det
i
*weights(i)*btdb
i
i
=
1
where
weights(i)
are the numerical integration weighting coefficients from (3.8).
