3.5 Approximations for stiffness matrices

Back

A stiffness matrix is exact if the strain energy product a(., .) is exact, if the

unit displacements ϕi are exact, and if the integration is done exactly:

a(ϕi, ϕj) =

_ l

0

MiMj

EI

dx . (3.118)

Green’s first identity explains how the strain energy product is defined. It

is more difficult to find the homogeneous solution of the differential equation.

The unit displacements are based on this solution. Approximate stiffness

matrices are usually based on approximate unit displacements, which are

substituted into the correct strain energy product.

3.5 Approximations for stiffness matrices 299

Fig. 3.22. Variable width

In a belly-shaped beam (Fig. 3.22) for example, the bending stiffness

E I(x) varies with x, so the beam equation becomes lengthy:

EI

__(x)w

__(x) + 2EI

_(x)w

___(x)

_ _ _

additional terms

+EI(x)wIV (x) = p(x) , (3.119)

while in the strain energy product little changes, only that I(x) now depends

on x:

a(w, ˆ w) =

_ l

0

M ˆM

EI(x) dx =

_ l

0

EI(x)w

__ ˆ w

__

dx . (3.120)

If the width h varies linearly,

EI(x) = E

bh3(x)

12 h(x) = a0 + a1 x (3.121)

it is possible to find the homogeneous solution of (3.119). But even Maple or

Mathematica will have difficulty finding the homogeneous solutions of more

complicated shapes.

Then the only way out is to substitute the unit displacements of the beam

with a constant bending stiffness EI into the exact strain energy product and

to work with this approximate matrix ˜K :

˜k ij =

_ l

0

EI(x) ϕ

__

i (x) ϕ

__

j (x) dx . (3.122)

Because the unit displacements are not exact, also the equivalent nodal forces

(or negative end-fixing forces p i) will not be consistent, i.e., they are only

approximations of the real fi:

˜ f i =

_ l

0

p ϕi(x) dx . (3.123)

In the case of a beam on an elastic foundation (see Fig. 3.23),

EIwIV (x) + cw(x) = p(x) , (3.124)

the exact unit displacements are well known, but program authors prefer

to use the unit displacements of the standard beam element, because this

300 3 Frames

Fig. 3.23. Beam on an elastic

foundation (Winkler)

facilitates program maintenance. That is, they calculate the elements of the

stiffness matrix

k ij =

_ l

0

[MiMj

EI

+ cϕi ϕj ] dx (3.125)

with the unit displacements ϕi(x) of the standard beam element. This yields

the approximation

˜K

= EI

l3

⎡

⎢⎢⎣

12 −6 l −12 −6 l

−6 l 4 l2 6 l 2 l2

−12 6 l 12 6 l

−6 l 2 l2 6 l 4 l2

⎤

⎥⎥⎦

+ c

420

⎡

⎢⎢⎣

156 l −22 l2 54 l 13 l2

−22 l2 4 l3 −13 l2 −3 l3

54 l 13 l2 156 l 22 l2

13 l2 −3 l3 22 l2 4 l3

⎤

⎥⎥⎦

.

(3.126)

Hence the FE program connects the end points of the deformed beam (end

displacements ui) with a curve w that is the sum of the four unit displacements

of the standard beam,

w(x) =

_4

i=1

ui ϕi(x) . (3.127)

This elastic curve deviates from the exact shape because it is not a homogeneous

solution of (3.124). Rather, a residual force appears,

EIwIV (x) + cw(x) = c[u1ϕ1(x) + u2ϕ2(x) + u3ϕ3(x) + u4ϕ4(x)] ,

(3.128)

which is just the distributed load p that must be applied to force the beam

into the shape w(x) (see Eq. (3.127)).

In second-order beam theory,

EIwIV (x) +P w

__(x) = p , (3.129)

the procedure is virtually the same. If the unit displacements ϕi(x) of the

standard beam are substituted into the exact strain energy product

a(w, ˆ w) =

_ l

0

[EIw

__ ˆ w

__ −P w

_ ˆ w

_] dx , (3.130)

3.5 Approximations for stiffness matrices 301

Fig. 3.24. Second-order

beam theory

the resulting approximate stiffness matrix

˜K

= EI

l3

⎡

⎢⎢⎣

12 −6l −12 −6l

−6l 4l2 6l 2l2

−12 6l 12 6l

−6l 2l2 6l 4l2

⎤

⎥⎥⎦

− P

30 l

⎡

⎢⎢⎣

36 −3 l −36 −3 l

−3 l 4l2 3 l −l2

−36 3l 36 3 l

−3 l −l2 3 l 4 l2

⎤

⎥⎥⎦

, (3.131)

is the sum of the first-order stiffness matrix and the so-called geometric matrix.

This procedure amounts to a Taylor expansion of the exact stiffness matrix

K = K(P) at the point P = 0.

If in a Timoshenko beam element [−1, 1] the deflection w(ξ) and the rotation

θ(ξ) are linear

_

w(ξ)

θ(ξ)

_

=

1

2

_

1 − ξ 0 1+ξ 0

0 1− ξ 0 1+ξ

_

⎡

⎢⎢⎣

w1

θ2

w3

θ4

⎤

⎥⎥⎦

(3.132)

then the resulting stiffness matrix is only an approximation

K = EI

le

⎡

⎢⎢⎣

0 0 0 0

0 1 0 −1

0 0 0 0

0 −1 0 1

⎤

⎥⎥⎦

+ GAs

6 le

⎡

⎢⎢⎣

6 3le 6 3le

3 le 2 l2

e

−3 l2 l2

e

6 −3 le 6 −3 le

3 le l2

e

−3 l2 2 l2

e

⎤

⎥⎥⎦

(3.133)

because the linear shape functions do not satisfy the differential equations

(3.105).

Equilibrium ?

In first-order beam theory the load and the support reactions maintain equilibrium

with regard to the undeformed structure and in second-order beam

theory with regard to the deformed structure—correct? No. Second-order theory

actually is a mixture of both; see Fig. 3.24. The lateral deformation is

considered, but the longitudinal displacement of the axis is neglected.

Theoretically it is therefore not possible to check the equilibrium of a

frame, that was analyzed with second-order beam theory, because first-order

and second-order effects contribute to the movements of the joints. If we ignore

302 3 Frames

Fig. 3.25. Buckling

this and check the equilibrium conditions by employing the nodal displacements

of the deformed structure then the check will fail. The reason that it

does not fail in practice is that normally the displacements are so small, that

the deviations are considered to be rounding errors.

Buckling length

Many design codes allow to establish the stability of a structure by determining

the buckling loads of the individual members. The critical load Pcrit of a

frame member is expressed in the form

Pcrit = π2 EI

(K · l)2 (3.134)

where K is termed effective length factor . That is the load Pcrit is expressed in

terms of the critical load of an equivalent pin-ended frame member of length

sK = K · l.

The standard FE-approach instead is a global approach based on 2nd order

theory and usually also including possible imperfections of the structure. The

result of the analysis are the eigenvalues λ of the structure based on the elastic

stiffness matrix K and the geometric stiffness matrix KG

(K + λKG)u = 0 u       = 0 . (3.135)

3.5 Approximations for stiffness matrices 303

Fig. 3.26. Flag pole

However it is not always easy to determine the critical eigenvalues, sometimes

even negative eigenvalues are observed which pose a problem for most of the

algorithms, [215].

In the case of a multi-story building as in Fig. 3.25, with a slender antenna

on top of the roof, we clearly see that it is not the smallest eigenvalue which is

the critical eigenvalue but that we have to check a whole range of eigenvalues.

Another problem pose those structures where the length of the individual

members can change. In the case of the flag pole in Fig. 3.26 the maximum

bending moment occurs just above the horizontal support. The analysis of the

upper part of the structure is standard, but in the lower part, the short beam,

the buckling length is quite small, and it is not possible to find a correct stress

in that part. If the horizontal support also carries vertical loads the normal

force in the lower part is zero and then it is not possible to apply the buckling

length approach to the design of this part.

This holds also true for the structure in Fig. 3.27 where

N1 =F N2 = N3 = N4 = N5 = −S ε:= l

%

F

EI

. (3.136)

The work done by the two forces F and the moment Ma on acting through

the displacements in Fig. 3.27 b (= influence function for Ma) must be zero

δWe = (Ma + F · l) · 1 − S · L · (−0.5) · (−0.5) · 2 − S · L · 0.5 · 0.5 · 2

= Ma + F · l − S · L = 0. (3.137)

Given a horizontal force H, see Fig. 3.28, the bending moments are

MI

a = H ·l MII

a = Ma = H · l · tanh ε

ε

≤ H · l . (3.138)

304 3 Frames

Fig. 3.27. What is the critical load? a) structure; b) influence function for Ma,

displacements not to scale

Fig. 3.28. The line of action of R intersects

the frame element only once

The difference between MI

a and MII

a is due to the force F

MI

a

−MII

a = Ma · ( ε

tanh ε

− 1) = F · l · 1.0 (3.139)

where we assumed that—as in the influence function—the tip of the element

moves sideways by Δw = l · tan ψ = l · 1.0 units. So that with

Ma = F · l · 1.0

ε/ tanh ε − 1 S = − F

2 · sin α

(3.140)

it follows

tanh ε

ε

= cos(2 · α) . (3.141)

3.6 Cables 305

Fig. 3.29. Cable

Because of

tanh ε

ε

≤ 1 ⇒ cos(2 · α) ≤ 1 ⇒ α ≤ π

4

(3.142)

the system is stable if α > π/4 while for α ≤ π/4, the system will buckle if

F exceeds the critical load Fcrit = Fcrit(α). But because the line of action of

the resultant R intersects the buckled frame element only once (see Fig. 3.28)

the buckling length approach is not applicable to this problem, [210].

Навигация