3.5 Rigid Body Motion

Back

A first example of nonlinear multiple point constraints constitutes rigid body motion. Here,

nonlinearity arises because of large rotations. In what follows, rectangular coordinates are

assumed and the spatial frame coincides with the material frame.

3.5.1 Large rotations

Consider a vector θ = θn along an axis AB (Figure 3.9), and a vector r0. Now, the vector

r0 is rotated about the axis AB until the new vector r includes an angle θ = _θ_ with r0.

We would like to find an expression for r as a function of r0 and θ.

156 GEOMETRIC NONLINEAR EFFECTS

_

_

_

_

_

Start of increment

Update MPC

_

Update MPC

Uind

0

= U0

U

dep

0

= f (Uind

0 )

S0 = f (U0)

Fint

0

= f (S0,U0)

_

i = 1

_

i ++

Determine Fext

i1

_

K _Uind

i1

= Fext

i1

Fint

i1

_

Uind

i

= Uind

i1

+ _Uind

i1

U

dep

i

= f (Uind

i )

Si = f (Ui )

Fint

i

= f (Si ,Ui )

_

_Fext Fint_ < _?

_______

__ __ __ _

_______

__ __ __ _

_

_

_

_

_

_

_

_

End of increment

No

Yes

Figure 3.8 Stream chart of the nonlinear solution procedure

“dep” = dependent, “ind” = independent, “int” = internal, “ext” = external

GEOMETRIC NONLINEAR EFFECTS 157

For an infinitesimal angle dθ, the change dr of r is perpendicular to r and satisfies

dr = dθ(n × r) (3.59)

in component notation:

dri = dθeij knj rk. (3.60)

Defining the matrix S by

Sik := eij knj (3.61)

one finds

dr = dθS r (3.62)

or

dr

dθ

= S r. (3.63)

This is a linear homogeneous vector differential equation with the solution

r = eSθ r0 (3.64)

satisfying the initial condition r(0) = r0. Equation (3.64) can be expanded into

r = _I + θS + 1

2!

θ2S2 + 1

3!

θ3S3 + _ r0. (3.65)

A

B

r

o

r0

dr0

θ = _θ_

θn = θ

Figure 3.9 Large rotation about the axis AB

158 GEOMETRIC NONLINEAR EFFECTS

Since S r = n × r (Equations (3.59) and (3.62)) and a × (b × c) = (a c)b (a b)c,

one finds

S2 r = S (S r) = n × (n × r) = (n r)n r (3.66)

S3 r = S (S2 r) = n × [(n r)n r] = n × r = S r (3.67)

from which one finds

S3 = S. (3.68)

Accordingly, all powers of S exceeding 2 can be reduced to Ѓ}S or Ѓ}S2. Consequently,

eSθ = I + S _

θ 1

3!

θ3 + 1

5!

θ5 _

+ S2 _ 1

2!

θ2 1

4!

θ4 + 1

6!

θ6 _

= I + sin θS + (1 cos θ)S2. (3.69)

Hence,

r = _I + sin θS + (1 cos θ)S2_ r0. (3.70)

Since

S2 = n n I (3.71)

this also reduces to

r = [cos θI + sin θS + (1 cos θ)n n] r0. (3.72)

Defining

ˆ θ = θS (3.73)

finally yields

r = C r0 (3.74)

where

C = _cos θI + sin θ

θ

ˆ θ + (1 cos θ)

θ θ

θ2

_ (3.75)

or in component notation,

Cij = δij cos θ + sin θ

θ

eikj θk + _1 cos θ

θ2

_

θiθj . (3.76)

GEOMETRIC NONLINEAR EFFECTS 159

Notice that this is a nonlinear relation in θ. Therefore, only a truly nonlinear calculation can

take large rotations into account. In simple linear calculations, Equation (3.59) is sometimes

used for finite rotations, yielding

r = r0 + θ(n × r0). (3.77)

Using this relation amounts to the motion in Figure 3.10 and is only feasible for a small

θ. The true angle α satisfies

α = arctan θ θ θ3

3

+ (3.78)

and _r_ satisfies

_r_ = r0_θ2 + 1 r0

_1 + θ2

2

_

. (3.79)

3.5.2 Rigid body formulation

Defining a set of nodes to behave like a rigid body means that all degrees of freedom of

the set are reduced to six degrees of freedom: three translations w of a point A and three

rotations θ about point A. Point A can be the center of gravity of the node set, but this

does not have to be. Any point will do. Usually, we take an existing node belonging to

the rigid node set to be point A. However, we can also generate an additional fictitious

node to be point A. Hence, the motion u of a node at location p can be described as

(Figure 3.11)

u = w + [C(θ) I ] (p q) (3.80)

r

r0

θ(n × r0)

α

Figure 3.10 Linearized rotation

160 GEOMETRIC NONLINEAR EFFECTS

A

u

q

p

p q

p q

w

w

(C I ) (p q)

C (p q)

θ

Figure 3.11 Rigid body motion of p about A

where w represents the motion of point A and q its location. The first term on the right-hand

side represents the translation and the second represents the rotation. Equation (3.80) is a

nonlinear relationship since C(θ) is nonlinear (Equation (3.75)). Linearizing at (w0, θ0),

as described in Section 3.4, yields

u0 + I (u u0) = w0 + I (w w0) + [C(θ0) I ] (p q) (3.81)

+ _C

θ

(θ0) (θ θ0)

_ (p q)

or

(u u0) = I (w w0) + _C

θ

(θ0) (θ θ0)

_ (p q) (3.82)

+ w0 + [C(θ0) I ] (p q) u0.

In component notation, this reads

ui u0i = wi w0i + _C

θ

(θ0)

_

ij l

(θl θ0l)(p q)j

+ w0i +            Cij (θ0) δij  (p q)j u0i (3.83)

GEOMETRIC NONLINEAR EFFECTS 161

where, differentiating Equation (3.76),

_C

θ

(θ0)

_

ij l

= Cij

θl

(3.84)

= cos θ

θl

δij +

θl

_sin θ

θ

_

eikj θk

+ _sin θ

θ

_

eilj +

θl

_1 cos θ

θ2

_

θiθj

+ _1 cos θ

θ2

_

_δilθj + θiδjl_ . (3.85)

The first term in Equation (3.83) is linear in the translations wi , the second term is linear

in the rotations θl and the third term is constant. The derivatives in Equation (3.85) are

easily determined

cos θ

θl

= θl

sin θ

θ

(3.86)

θl

_sin θ

θ

_ = θl

θ3 (θ cos θ sin θ) (3.87)

θl

_1 cos θ

θ

_ = θl

θ4 (θ sin θ 2 + 2 cos θ). (3.88)

For small values of θ, these expressions are undetermined and the limit must be taken

lim

θ0

sin θ

θ

= 1 (3.89)

lim

θ0

1 cos θ

θ2

= 1

2

(3.90)

lim

θ0

θ cos θ sin θ

θ3

= 1

3

(3.91)

lim

θ0

θ sin θ 2 + 2 cos θ

θ4

= 1

12

. (3.92)

Equations (3.83) are the linearized rigid body multiple point constraints at (w0, θ0).

Whereas the translational degrees of freedom can be associated with an existing node,

this is not the case for the rotational degrees of freedom. The easiest solution is to assign

them to a new fictitious node, that is, the translational degrees of freedom of the new node

are interpreted as the rotational degrees of freedom of the rigid body.

The above procedure assumes that there is a one-to-one relationship between the motion

of the body and the translation and rotation expression given by (w, θ). If this is not the

case, additional measures must be taken. For instance, if the body consists of points lying

on a straight line, the rotation about this line is not uniquely determined. In that case, the

rotation about the line must be explicitly assigned. Assume that a is a unit vector along the

line, then, setting the rotation about the line to zero amounts to the linear multiple point

constraint a θ = 0.

162 GEOMETRIC NONLINEAR EFFECTS

3.5.3 Beam and shell elements

The present section looks into a three-dimensional expansion theory of beam and shell

elements. Beam and shell structures are characterized by small dimensions across their

thickness. Therefore, simplified assumptions can be applied in the thickness direction,

leading to different formulations. In the simplest forms, straight fibers orthogonal to the

midplane in plates and shells and to the midline in beams are assumed to stay straight

and orthogonal during deformation. This leads to the Kirchhoff theory for plates and the

Bernoulli–Euler theory for beams. If the fibers do remain straight during deformation but

not necessarily orthogonal to the midplane/midline, the formulation is called the Mindlin

theory for plates and the Timoshenko theory for beams, see also (Zienkiewicz and Taylor

1989), (Graff 1975) and (Meirovitch 1967). The assumptions regarding the displacement

field across the thickness have the advantage that only the middle plane (midline) needs

to be modeled, while the changes across the thickness are covered by the introduction of

additional rotational degrees of freedom in the nodes. Accordingly, modeling needs are

basically reduced to the creation of a two-dimensional mesh of the (curved) surface (for

shells/plates) or a one-dimensional mesh of the beam axis. The price to be paid is the need

for the derivation of the material stiffness matrix specifically for shell and/or beam elements,

due to the special formulation in terms of rotational degrees of freedom. Therefore, the idea

of hybrid shell-solid and even pure-solid formulations has come up in different forms in

recent years, (Bischoff and Ramm 1999), (Flores and O˜nate 2001), (Wriggers et al. 1996),

(DЁuster et al. 2001) and (Sze et al. 2002).

In the present derivation, a new, pure-solid way is selected. The ease of modeling

is kept by reducing the shells and beams to their midplane and centerline respectively.

However, instead of introducing rotational degrees of freedom, 8-node quadratic shell or

plate elements and 3-node quadratic beam elements are expanded into 1 layer of 20-node

brick elements (with full or reduced integration). Quadratic elements are chosen because

of their intrinsically good properties: they are known to behave well for slender structures

and rarely exhibit locking or hourglassing. The way of expansion is shown in Figures 3.12

and 3.13.

As long as the plate, shell or beam is smooth, the expansion results in a threedimensional

connected continuum model. However, problems arise as soon as sharp kinks

need to be modeled, in areas where several shells and beams cross or the thickness of the

shells or beams changes discontinuously. At such locations, all nodes expanded from one

and the same node are considered to behave like a rigid body, and will be called a knot.

At a knot, the degrees of freedom are reduced to three translational and three rotational

degrees. All participating structures are expanded as stand-alone parts. Figure 3.14 shows

the expansion at a knot between shells and Figure 3.15 at a knot between beams. The

structures partially overlap.

A knot is also introduced between beams with a different offset and/or with different

cross section (Figure 3.16) and in composed shells and beams. The I-cross

section in Figure 3.17 consists of three beam elements with exactly the same nodes,

but with different cross section and different offset. Since the cross section is defined

as a rigid body, it will remain plane and no warping will occur. However, shear

deformation is possible since the cross section does not have to remain orthogonal

to the central axis. The expanded structure is a volume model and has no rotational

GEOMETRIC NONLINEAR EFFECTS 163

Thickness 1

Thickness 2

t

􀀀 􀀀 􀀀

_ _ _

􀀀 􀀀 􀀀

_ _ _

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀

􀀀_ _

􀀀

􀀀_ _

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀􀀀

􀀀􀀀

__

__

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀

􀀀

􀀀

_

_

_

􀀀 􀀀 􀀀

_ _ _

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀 􀀀 􀀀

_ _ _

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀 􀀀 􀀀

_ _ _

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀􀀀

􀀀􀀀

__

__

1

1

2

2

3

4

5 6

8 7

9

10

11

12

13

14

15

15

16

17 18

19

20

1

2

2 3

Nodes of the 1-D element

Nodes of the 3-D element

Figure 3.12 Expansion of the one-dimensional element

Thickness

􀀀􀀀

􀀀􀀀

__

__

􀀀

􀀀

_

_

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀􀀀

􀀀􀀀

__

__

􀀀

􀀀

_

_

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀

􀀀

􀀀

_

_

_

􀀀

􀀀

_

_

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀􀀀

􀀀􀀀

__

__

􀀀

􀀀

􀀀

_

_

_

􀀀 􀀀

􀀀

_

_

_

􀀀

􀀀

_

_

􀀀

􀀀

_

_

􀀀 􀀀

􀀀

_

_

_

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

􀀀

􀀀

_

_

􀀀

􀀀

􀀀

_

_

_

􀀀􀀀

􀀀􀀀

__

__

􀀀

􀀀

􀀀

_

_

_

􀀀􀀀

􀀀􀀀

􀀀􀀀

__

__

__

1 2

3

4

5 6

8 7

9

10

11

12

13

14

15

15

16

17 18

19

20

1 2

4 3

5

5

6

7

8

Nodes of the 2-D element

Nodes of the 3-D element

Figure 3.13 Expansion of the two-dimensional element

164 GEOMETRIC NONLINEAR EFFECTS

Figure 3.14 Knot between shells

Figure 3.15 Knot between beams

GEOMETRIC NONLINEAR EFFECTS 165

Figure 3.16 Knot between beams with different offset and different cross section

Figure 3.17 I-cross section composed of three simple beam elements each with a different

offset

166 GEOMETRIC NONLINEAR EFFECTS

x

A

B

R

F

F

Symmetry plane

Symmetry plane

18

0.004 R

One node fixed in z

Free

Free

Free

y

z

ν = 0.3

Figure 3.18 Hemispherical shell loaded by concentrated forces

degrees of freedom. Therefore, knots are also introduced at nodes where the user has

defined rotations.

The foregoing expansion can also be applied to plane stress, plane strain and axisymmetric

elements. Any mixing of these element types among each other or with beams

and shells is also taken care of by the knots. However, in the case of plane stress, plane

strain or axisymmetric elements, the rigid body definition is restricted to the nodes in the

midplane or along the centerline. Indeed, the off-center nodes in plane stress, plane strain

and axisymmetric elements are subject to additional conditions due to the z-symmetry or

axisymmetry, which would collide with the rigid body definition.

As an example, consider the thin hemispherical shell with a hole at the top and loaded

by concentrated forces (Figure 3.18). The shell is meshed in three different ways:

1. As a three-dimensional structure using genuine 20-node brick elements with full

integration. The 8 × 10 element mesh contains one element over the thickness (1872

degrees of freedom in total). The length to the thickness ratio of the elements is

about 40. All nodes in the x z plane are fixed in the y-direction, and the nodes in

the y z plane are fixed in the x-direction. This description contains translational

degrees of freedom only.

2. As a three-dimensional structure using genuine 20-node brick elements with reduced

integration. The same comments as under 1 also apply here.

3. As a shell structure meshed by 8 × 10 quadratic shell elements with reduced integration.

In the x z plane, the translational degrees of freedom in the y-direction

and the rotational degrees of freedom about the x-axis and z-axis are fixed, in the

y z plane the translational degrees of freedom in the x-direction and the rotational

GEOMETRIC NONLINEAR EFFECTS 167

Table 3.1 Displacements of nodes A and B.

Load ABAQUS20-node brick 20-node brick 8-node shell

109 F

ER2 4-node shell full integration red. integration red. integration

1

R

ux,A

1

R

uy,B

1

R

ux,A

1

R

uy,B

1

R

ux,A

1

R

uy,B

1

R

ux,A

1

R

uy,B

5.86 0.326 0.232 0.138 0.114 0.329 0.227 0.324 0.222

8.79 0.434 0.282 0.191 0.147 0.447 0.277 0.442 0.272

14.65 0.590 0.341 0.271 0.190 0.618 0.334 0.610 0.328

degrees of freedom about the y-axis and z-axis are fixed. The shell elements are

internally automatically expanded into 20-node brick elements with reduced integration.

Along x = 0 and y = 0, rigid knots are introduced to take care of the rotational

degrees of freedom.

The displacements of nodes A and B in x- and y-direction respectively, are listed in

Table 3.1 and compared with ABAQUSreference results. The 20-node brick elements

with full integration are clearly too stiff. However, the elements with reduced integration

show good agreement with the reference results even for highly nonlinear deformations.