Principle of Virtual Displacements in Structural Dynamics

CEE 541. Structural Dynamics

Department of Civil and Environmental Engineering
Duke University
Henri P. Gavin
Fall 2020

1 An introductory example of potential energy in elastic structures and virtual work

Consider the two bar truss shown below. With the horizontal load F acting to the right, bar 1 is in
tension, T1 > 0 and bar 2 is in compression, T2 < 0. Suppose the loaded node moves to the right
Dx and also upward Dy . If the truss is made of a linear elastic material and if the displacements are
small compared to the overall size of the truss, then the displacements Dx and Dy increase linearly
with the force F . The work of the external force F increasing in proportion to the its collocated
displacement Dx is the area under the F − D line. This external work is simply.
1
W = F Dx
2
If the two-bar truss is made from elastic materials and there are no other energy dissipation mech-
anisms (e.g., friction) involved, then the external work of pushing the truss node to the right by a
displacement Dx is entirely stored within the bars of the truss. Knowing the bar tensions, the elastic
modulus of the truss bars (E1 and E2 ), their cross section areas (A1 and A2 ), and their lengths (L1
and L2 ), the stretch of each bar is
d1 = T1 L1 /(E1 A1 )
and
d2 = T2 L2 /(E2 A2 ).
The internal potential energy stored within the truss is
1 1
U = T1 d1 + T2 d2
2 2
Since the work of pushing the truss node to the right is entirely stored as potential energy within
the truss.
U =W
1 1 1
T1 d1 + T2 d2 = F Dx (1)
2 2 2
This is called the principle of real work.
Given numerical values for F , L1 , E1 , A1 , L2 , E2 , and A2 the principle of real work could be used to
compute the displacement Dx , collocated with the single applied force, F .
2 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

1 2

1 2

Figure 1. A two bar truss loaded first with a horizontal force F and subsequently with combined
horizontal force F and vertical force δF . The force F is held constant while the force δF is
added.

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 3

Now, with the force F held constant (and the related bar tensions T1 and T2 held constant), suppose
a vertical force δF is applied, resulting in additional displacements δDx and δDy , and additional bar
forces δT1 and δT2 . The displacements Dy and δDy are collocated with the applied force δF , and
so the work of δF increasing proportionally through the collocated displacement δDy is ( 21 δF δDy ).
This work is stored as potential energy within the bars of the truss with with tensions δT1 and δT2
proportional and bar stretches, δd1 and δd2 .
1 1 1
δT1 δd1 + δT2 δd2 = δF δDy (2)
2 2 2

Additionally, as δF increase, the constant force F moves through a displacement δDx , and the
constant bar tensions T1 and T2 move through displacements δd1 and δd2 . The work of the constant
force F moving through displacement δDx , (F δDx ), is called the virtual work of the external force.
And the work of the constant bar forces T1 and T2 moving through displacement δd1 and δd2 ,
(T1 δd1 + T2 δd2 ), is called the virtual work of internal forces.
With the total combined forces applied to the truss, the combined external work is
1 1
W = F Dx + δF δDy + F δDx
2 2
and the combined internal potential energy is
1 1 1 1
U = T1 d1 + T2 d2 + δT1 δd1 + δT2 δd2 + T1 δd1 + T2 δd2 .
2 2 2 2
Setting the external work of the combined forces equal to the internal potential energy of the combined
forces, and substituting (1) and (2) into this equality results in
T1 δd1 + T2 δd2 = F δDx . (3)
This is called the principle of virtual work, or in this specific example, the principle of virtual dis-
placements. In words, the principle of virtual displacements states that: For a “real” external force
(e.g., F ) in equilibrium with a set of “real” internal forces (e.g., T1 and T2 ), and a “virtual” external
displacement (e.g., δDx ) collocated with the external force (e.g., F ) and kinematically compatible
with “virtual” internal displacements (e.g., δd1 and δd2 ) collocated with internal forces (e.g, T1 and
T2 ), the external work of the “real” external force moving through a collocated “virtual” displacement
equals the internal work of the “real” internal forces moving through collocated “virtual” internal
displacements.
Given numerical values for δF , L1 , E1 , A1 , L2 , E2 , and A2 the principle of virtual displacements
could be used to compute the horizontal displacement δDx caused by the vertical force, δF . To do so,
one would solve for the bar forces T1 and T2 in equilibrium with F , then solve for the bar forces δT1
and δT2 in equilibrium with δF , compute δd1 and δd2 from δT1 and δT2 , and, finally solve equation
(3) for δDx , the horizontal displacement due to the vertical force δF .
• What would be the approach to find the vertical displacement Dy due to the horizontal force, F ?
• If F = δF , how are Dy and δDx related? Are you surprised by this result?
• If F = δF and Dy < 0, re-draw Figure 1.

cbnd H.P. Gavin October 10, 2022

4 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

2 Strain Energy in Elastic Solids

Consider an elastic solid object with external forces F and f in equilibrium with internal stresses σ
and τ .
Fi 1111111111111111
0000000000000000
0000000000000000
1111111111111111 F1
0000000000000000
1111111111111111
0000000000000000
1111111111111111 f R n Fn
Ri S1111111111111111
0000000000000000
0000000000000000
1111111111111111
00000000000000
11111111111111
0000000000000000
1111111111111111 R1
00000000000000
11111111111111
r
00000000000000
11111111111111

11111
00000
Ri Rj
000
111 00000
11111 00000
11111
111
000
000
111 00000
11111 11111
00000
00000
11111
000
111 00000
11111 00000
11111
000
111
000
111 Fi Fj
000
111
000
111
000
111
000
111
σ ε F1
Fi Fj Fn
R1 f
000000
111111
Ri Rj
111111
000000
000000
111111
000000
111111 R
r n
00000
11111
11111
00000 00000
11111
11111
00000 00000
11111
11111
00000
00000
11111 00000
11111 00000
11111

• F and f are real external forces in equilibrium, acting at points, or over a portion of a surface
S,

• R and r are real displacements, admissible with respect to the support conditions, collocated
with F and f ,

• σ are real internal stresses, distributed within the solid volume V , in equilibrium with F & f ,

• ϵ are real internal strains, distributed within the solid volume V , compatible with R and r,

2.1 External Work

The work of external forces increasing from 0 to F and f and pushing through displacements from
0 to R and r is Z R Z Z r
W = F (R̄) dR̄ + f (r̄) dr̄ dS (4)
0 S 0
where

• the forces F and f depend on displacements R and r

• R̄ and r̄ are dummy variables of integration

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 5

2.2 Internal Strain Energy

Strain energy is a kind of potential energy arising from stress and deformation of elastic solids. In
nonlinear elastic solids, the strain energy of stresses increasing from 0 to σ(x, y, z) and working
through strains from 0 to ϵ(x, y, z) is
Z Z ϵ
U= σ(ϵ̄) · dϵ̄ dV (5)
V 0

where

• V is the volume of the solid

n o
• σ(ϵ) = σxx σyy σzz τxy τyz τxz
n o
• ϵ= ϵxx ϵyy ϵzz γxy γyz γxz

• ϵ̄ is a dummy variable of integration ... so ...

Z Z ϵ
U = σ(ϵ̄) · dϵ̄ dV
ZV Z0ϵ
= (σxx (ϵ̄)dϵ̄xx + σyy (ϵ̄)dϵ̄yy + σzz (ϵ̄)dϵ̄zz + σxy (ϵ̄)dϵ̄xy + σyz (ϵ̄)dϵ̄yz + σxz (ϵ̄)dϵ̄xz ) dV
V 0

2.3 The Principle of Real Work

In an elastic solid, the work of external forces, W , is stored entirely as elastic strain energy, U , within
the solid.
U =W (6)
In linear elastic solids:

• Stresses σ increase linearly with strains ϵ,

σ = Eϵ ... and ... τ = Gγ

• Displacements D and rotations Θ increase linearly with forces F and moments M ,

F = kD ... and ... M = κΘ

• The work of an external force F acting through a displacement D on the solid is

1 1 1
W = F D = kD2 = F 2 /k
2 2 2

• The work of an external moment M acting through a rotation Θ on the solid is

1 1 1
W = M Θ = κΘ2 = M 2 /κ
2 2 2

cbnd H.P. Gavin October 10, 2022

6 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

111111111111111111111111111111111111
000000000000000000000000000000000000
000000000
111111111 1111111111111111111111111111111111111
0000000000000000000000000000000000000
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111
000000000
111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
000000000
111111111
σ(ε) 000000000000000000000000000000000000
111111111111111111111111111111111111
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
000000000
111111111
i 1111111111111111111111111111111111111
F(R’)0000000000000000000000000000000000000
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111
dU 000000000
111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
W
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111 000000000
111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111 000000000
111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111
000000000
111111111 000000000
111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111
000000000
111111111 d ε’ 0000000000000000000000000000000000000
1111111111111111111111111111111111111
000000000
111111111 dR’i
000000000000000000000000000000000000
111111111111111111111111111111111111
000000000
111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
000000000
111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
000000000
111111111
000000000000000000000000000000000000
111111111111111111111111111111111111 000000000
111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
0
000000000
111111111 0
000000000
111111111 Ri R’i
ε
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000 0000000000000000000000000000000000000
1111111111111111111111111111111111111
σ 1111111111111111111111111111111111111
0000000000000000000000000000000000000 F i 1111111111111111111111111111111111111
0000000000000000000000000000000000000
0000000000000000000000000000000000000
1111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
1
0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
1
(σ.ε) dV W= 2 Fi R i
0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111
dU=
0000000000000000000000000000000000000
1111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
2
0000000000000000000000000000000000000
1111111111111111111111111111111111111 0000000000000000000000000000000000000
1111111111111111111111111111111111111
0000000000000000000000000000000000000
1111111111111111111111111111111111111 0
0000000000000000000000000000000000000
1111111111111111111111111111111111111
R
0 ε i

2.4 Strain energy in slender structural elements

In slender structural elements (bars, beams, or shafts) the internal forces, moments, shears, and
torques can vary along the length of each element; so do the displacements and rotations.
The strain energy of spatially-varying internal forces F (x) acting through spatially-varying internal
displacements D(x) along a linear elastic prismatic solids is

1Z dD(x) 1Z
U= F (x) · dx = F (x)D′ (x) dx (7)
2 l dx 2 l

The strain energy of spatially-varying internal moments M (x) acting through spatially-varying in-
ternal rotations Θ(x) along linear elastic prismatic solids is

1Z dΘ(x) 1Z
U= M (x) · dx = M (x)Θ′ (x) dx (8)
2 l dx 2 l

In slender structural elements, the relation between internal forces and moments and internal dis-
placements and rotations depend on the kind of loading.

• Axial Nx (x) = E(x)A(x)u′ (x)

• Bending Mz (x) = E(x)I(x)v ′′ (x) ... assuming v ′′ (x) ≈ Θ′ (x)
• Shear Vy (x) = G(x)As (x)vs′ (x)
• Torsion Tx (x) = G(x)J( x)ϕ′ (x)

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 7

Inserting these expressions into the general expressions for internal strain energy above,

“force” deformation strain energy (U )

Z Z Z
Nx (x)2
Axial Nx (x) u (x)
′ 1
2
Nx (x)u (x)dx
′ 1
2 E(x)A(x)
dx 1
2
E(x)A(x)(u′ (x))2 dx
l l l
Z Z Z
Mz (x)2
Bending Mz (x) v (x)
′′ 1
2
Mz (x)v (x)dx
′′ 1
2 E(x)I(x)
dx 1
2
E(x)I(x)(v ′′ (x))2 dx
l l l
Z Z Z
Vy (x)2
Shear Vy (x) vs′ (x) 1
2
Vy (x)vs′ (x)dx 1
2 G(x)As (x)
dx 1
2
G(x)As (x) (vs′ (x))2 dx
l l l
Z Z Z
Tx (x)2
Torsion Tx (x) ϕ′ (x) 1
2
Tx (x)ϕ′ (x)dx 1
2 G(x)J(x)
dx 1
2
G(x)J(x)(ϕ′ (x))2 dx
l l l

E(x) is Young’s modulus

G(x) is the shear modulus

A(x) is the cross sectional area of a bar

I(x) is the bending moment of inertia of a beam
A(x)/α is the effective shear area of a beam
J(x) is the torsional moment of inertia of a shaft

Nx (x) is the axial force within a bar

Mz (x) is the bending moment within a beam
Vy (x) is the shear force within a beam
Tx (x) is the torque within a shaft

u(x) is the axial displacement along the bar

u′ (x) is the axial displacement per unit length, du(x)/dx, the axial strain

v(x) is the transverse bending displacement of the beam

v ′ (x) is the slope of the displacement of the beam
v ′′ (x) is the rotation per unit length, the curvature, approximately d2 v(x)/dx2

vs (x) is the transverse shear displacement of the beam

vs′ (x) is the transverse shear displacement per unit length, dvs (x)/dx

ϕ(x) is the torsional rotation (twist) of the shaft

ϕ′ (x) is the torsional rotation per unit length, dϕ(x)/dx

cbnd H.P. Gavin October 10, 2022

8 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

3 Virtual Work in Elastic Solids — The Principle of Virtual Displacements

Now consider a second set of loads, δF , δf , in equilibrium and applied subsequently to the loads F
and f . The loads δF and δf give rise to displacements δR and δr collocated with forces F and f ,
and internal stresses δσ and strains δϵ. In other words, the displacements δR and δr are admissible
to the kinematic constraints.
Call δF and δf a set of any arbitrary “virtual” forces in equilibrium.
Call δR and δr a set of “virtual” displacements, collocated with forces F and f , and resulting from
forces δF and δf (and therefore kinematically admissible). The displacements δR and δr may also
be called variations of displacements, admissible to the constraints.
Forces F and f are held constant as loads δF and δf are applied. Stresses σ, in equilibrium with
forces F and f , are therefore also held constant as loads δF and δf are applied. Forces F and f do
not increase with displacements δR and δr. Strains δϵ increase as loads δF and δf are applied.
Fi
111111111111111
000000000000000
000000000000000
111111111111111
000000000000000
111111111111111 f
Ri S111111111111111
000000000000000
000000000000000
111111111111111
00000000000000
11111111111111
000000000000000
111111111111111
00000000000000
11111111111111
r
δf 00000000000000
11111111111111
000000000000000
111111111111111
0000000
1111111 000000000000000
111111111111111
1111111
0000000 δr
0000000
1111111 δRi δ Fj
0000000
1111111
0000000
1111111
0000000
1111111
0000000
1111111
000
111
000
111
111
000
000
111
000
111
000
111
000
111
000
111
000
111
000
111
δσ δε
00000
11111
11111
00000
00000
11111

111111
000000
000000
111111
000000
111111
000000
111111

000000000000000000000000
111111111111111111111111 1111111111111111111111111111111111111111
0000000000000000000000000000000000000000000000000000000000000000
111111111111111111111111
11111111111111111111111111111111111111111
00000000000000000000000000000000000000000
σ 000000000000000000000000
111111111111111111111111 0000000000000000000000000000000000000000
1111111111111111111111111111111111111111
Fi
000000000000000000000000
111111111111111111111111
00000000000000000000000000000000000000000
11111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111 0000000000000000000000000000000000000000
1111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111
00000000000000000000000000000000000000000
11111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111 0000000000000000000000000000000000000000
1111111111111111111111111111111111111111
0000000000000000000000000000000000000000111111111111111111111111
000000000000000000000000
00000000000000000000000000000000000000000111111111111111111111111
11111111111111111111111111111111111111111000000000000000000000000 1111111111111111111111111111111111111111
dV

00000000000000000000000000000000000000000
11111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111 0000000000000000000000000000000000000000
1111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111
δR

00000000000000000000000000000000000000000
11111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111 0000000000000000000000000000000000000000111111111111111111111111
1111111111111111111111111111111111111111000000000000000000000000
000000000000000000000000
111111111111111111111111
ε)

00000000000000000000000000000000000000000111111111111111111111111
11111111111111111111111111111111111111111000000000000000000000000 0000000000000000000000000000000000000000
1111111111111111111111111111111111111111
E=F

00000000000000000000000000000000000000000
11111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111 0000000000000000000000000000000000000000111111111111111111111111
1111111111111111111111111111111111111111000000000000000000000000
.δ

00000000000000000000000000000000000000000111111111111111111111111
11111111111111111111111111111111111111111000000000000000000000000 0000000000000000000000000000000000000000
1111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111
(σ

0000000000000000000000000000000000000000111111111111111111111111
000000000000000000000000
δW

00000000000000000000000000000000000000000
11111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111 1111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111
00000000000000000000000000000000000000000111111111111111111111111
11111111111111111111111111111111111111111000000000000000000000000 0000000000000000000000000000000000000000
1111111111111111111111111111111111111111
=

00000000000000000000000000000000000000000
11111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111 0000000000000000000000000000000000000000111111111111111111111111
1111111111111111111111111111111111111111000000000000000000000000
δW

00000000000000000000000000000000000000000111111111111111111111111
000000000000000000000000 0000000000000000000000000000000000000000
1111111111111111111111111111111111111111000000000000000000000000
111111111111111111111111
I

11111111111111111111111111111111111111111
0
00000000000000000000000000000000000000000
11111111111111111111111111111111111111111
ε
000000000000000000000000
111111111111111111111111
ε+δε 0 R111111111111111111111111
0000000000000000000000000000000000000000
1111111111111111111111111111111111111111000000000000000000000000
i R + δR
i i

The principle of virtual displacements states that the virtual external work of real external forces (f
and F ) moving through collocated admissible virtual displacements (δr and δR) equals the internal
virtual work of real stresses (σ) in equilibrium with real forces (f and F ) with the virtual strains
(δϵ) compatible with the virtual displacements (δr and δR), integrated over the volume of the solid.

Z
δWI = δW
Z E
σ · δϵ dV = f · δr dS + (9)
X
Fi · δRi
V S i

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 9

3.1 Work of axial loads and transverse displacements in slender structural elements
In slender solid elements, nonuniform transverse displacements (dv(x) ̸= 0) induce longitudinal
shortening, de(x).

de+d δe
dx
δ v + d δv

v’+ δ v’

dx de de
δv
v’ v + dv v’ v + dv

v v

dx dx

Figure 2. Transverse deformation v ′ (x) and longitudinal shortening de(x).

A relation between dv and de can be derived from the Pythagorean theorem and is quadratic in dv
and de.
(dx − de)2 + (dv)2 = (dx)2
2(de)(dx) − (de)2 = (dv)2
!2
de 1 dv 1
≈ = (v ′ )2
dx 2 dx 2
With additional virtual displacements δv(x) a relation for the incremental virtual shortening dδe
may also be derived from the Pythagorean theorem.
(dx − de − dδe)2 + (dv + dδv)2 = (dx)2
2(de)(dx) − 2(de)(dδe) + 2(dδe)(dx) − (de)2 − (dδe)2 = (dv)2 + 2(dv)(dδv) + (dδv)2
Subtracting 2(de)(dx) − (de)2 = (dv)2 and dividing by (dx)2 leaves
!2
de dδe dδe dδe
−2 +2 − = 2(v ′ )(δv ′ ) + (δv ′ )2
dx dx dx dx
Neglecting higher order terms (assuming virtual displacements are infinitesimal), leaves
dδe
≈ (v ′ )(δv ′ ) (10)
dx
The virtual work of a distributed axial compression P (x) (applied externally, for example, by gravita-
tional acceleration) acting through virtual shortening displacements δe(x) integrated along a slender
element is, then, Z
dδe Z
δWG = P (x) dx = P (x) v ′ (x) δv ′ (x) dx (11)
l dx l
This result can also be obtained by integrating along the arc-length of the deformed element as is
done in Tedesco, McDougal, and Ross’s textbook, Structural Dynamics: Theory and Applications.

cbnd H.P. Gavin October 10, 2022

10 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

4 The Principle of Virtual Displacements for Dynamic Loading

The principle of virtual displacements applies to both static and dynamic forces. Elastic forces
k(x)r(x, t) are present in structural systems responding to static or dynamic loads. Forces arising
from dynamic effects only include viscous damping forces c(x)ṙ(x, t) and inertial forces m(x)r̈(x, t).
Elastic forces, viscous damping forces, and inertial forces can be developed within slender structural
elements in response to axial, bending, shear, and torsional deformations.

real virtual
“force” deformation internal virtual work (δWI )

Z Z
Axial Nx (x, t) δu (x, t)
′
Nx (x, t) δu (x, t) dx
′
EA(x) u′ (x, t) δu′ (x, t) dx
l Zl
ηa A(x) u̇′ (x, t) δu′ (x, t) dx
lZ
ρA(x) ü(x, t) δu(x, t) dx
l
Z Z
Bending Mz (x, t) δv ′′ (x, t) Mz (x, t) δv ′′ (x, t) dx EI(x) v ′′ (x, t) δv ′′ (x, t) dx
l Zl
ηa I(x) v̇ ′′ (x, t) δv ′′ (x, t) dx
l Z
ρA(x) v̈(x, t) δv(x, t) dx
l
Z Z
Shear Vy (x, t) δvs′ (x, t) Vy (x, t) δvs′ (x, t) dx GAs (x) vs′ (x, t) δvs′ (x, t) dx
l Zl
ηs As (x) v̇s′ (x, t) δvs′ (x, t) dx
l Z
ρA(x) v̈(x, t) δv(x, t) dx
l
Z Z
Torsion Tx (x, t) δϕ′ (x, t) Tx (x, t) δϕ′ (x, t) dx GJ(x) ϕ′ (x, t) δϕ′ (x, t) dx
l Zl
ηs J(x) ϕ̇′ (x, t) δϕ′ (x, t) dx
lZ
ρJ(x) ϕ̈(x, t) δϕ(x, t) dx
l
Z Z
Geometric P (x) δe(x, t) P (x) δe(x, t) dx P (x) v ′ (x, t) δv ′ (x, t) dx
l l

In this table:

• The internal virtual work of viscous effects is derived assuming linear viscous stress - strain-rate
relations: σ = ηa ϵ̇ and τ = ηs γ̇. As will be seen later in the course, the damping properties of
real structural materials are actually more complicated.
• Rotatory inertia effects are neglected in the virtual work of inertial forces in bending beams.

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 11

5 Generalized Coordinates
A dynamic response r(x, t) may be represented as an expansion of products of spatially dependent
quantities and time dependent quantities
r(x, t) = ψk (x) qk (t) (12)
X

The functions ψk (x) are called shape-functions, and the functions q(t) may be called generalized
coordinates. In order for the above expansion to yield realistic and accurate solutions, the shape
functions must at least satisfy the essential boundary conditions. (The shape functions must be
kinematically-admissible.) Shape functions which also satisfy the natural boundary conditions will
yield more accurate solutions. Also, if the shape functions are dimensionless, the generalized coordi-
nates have the same units as the response, which permits a useful interpretation of the generalized
coordinates. Further, if the shape functions are kinematically admissible, and the expansion (12) for
r is expressed in terms of q, but not q̇, then virtual displacements defined as variations in r(x, t) with
respect to the set of coordinates qk (t) are also kinematically admissible
∂r(x, t)
δr(x, t) = δqj (t) = ψj (x) δqj (t) ,
X X

j ∂qj (t) j

and the derivatives of r with respect to x and t are

r(x, t) = k qk (t)ψk (x) ṙ(x, t) = k q̇k (t)ψk (x) r̈(x, t) = k q̈k (t)ψk (x)
P P P

r′ (x, t) = k qk (t)ψk (x)

′
ṙ′ (x, t) = k q̇k (t)ψk (x)
′
r̈′ (x, t) = k q̈k (t)ψk (x)
′
P P P

r′′ (x, t) = k qk (t)ψk (x)

′′
ṙ′′ (x, t) = k q̇k (t)ψk (x)
′′
r̈′′ (x, t) = k q̈k (t)ψk (x)
′′
P P P

Internal virtual work can also be expressed in terms of generalized virtual displacements. For example
in the elastic bending of a beam, the work of moments (EIv ′′ (x, t)) moving through virtual rotations
(δv ′′ (x, t) dx) in terms of generalized coordinate displacements qk (t) and virtual displacements δqj (t)
is
Z
δWI = EI(x) v ′′ (x, t) δv ′′ (x, t) dx
Zl
= ψk′′ (x) qk (t) ψj′′ (x)δqj (t) dx
X X
EI(x)
l k j
X X Z 
= EI(x) ψj′′ (x) ψk′′ (x) dx qk (t) δqj (t) (13)
j k l

The work of transverse inertial forces (ρAv̈(x, t) dx) moving through virtual displacements (δv(x, t))
in terms of generalized coordinate accelerations q̈k (t) and virtual displacements δqj (t) is
Z
δWI = ρA(x) v̈(x, t) δv(x, t) dx
Zl
= ψk (x)q̈k (t) ψj (x)δqj (t) dx
X X
ρA(x)
l k j
X X Z 
= ρA(x) ψj (x) ψk (x) dx q̈k (t) δqj (t) (14)
j k l

cbnd H.P. Gavin October 10, 2022

12 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

The work of external forces f (x, t) and F (t) moving through collocated virtual displacements δv(x, t)
can be expressed in terms of virtual displacements of generalized coordinates, δqj (t).
Z
δWE = f (x, t) · δv(x, t) dx + Fi · δv(x, ti )
X
l i
Z
= f (x, t) · ψj (x) δqj (t) dx + ψj (xi ) δqj (t)
X X X
Fi ·
l j i j
" #
X Z 
= f (x, t) · ψj (x) dx δqj (t) + Fi · ψj (xi ) δqj (t) (15)
X X

j l j i

The external virtual work of axial compression P (x) moving through virtual end shortening (v ′ (x, t)δv ′ (x, t) dx)
in terms of generalized coordinate displacements qk (t) and virtual displacements δqj (t) is
Z
δWE = P (x) · v ′ (x, t) δv ′ (x, t) dx
Zl
= P (x) · ψk′ (x) qk (t) ψj′ (x) δqj (t) dx
X X
l k j
X X Z 
= P (x) ψj′ (x) ψk′ (x) dx qk (t) δqj (t) (16)
j k l

By setting the internal virtual work equal to the external virtual work, and factoring out the inde-
pendent and arbitrary variations δqj , equations (13), (14), (15), and (16), result in
   
[M ]q̈(t) + [KE ]q(t) − [KG ]q(t) − f (t) · δq(t) = 0
Noting that each variation ψj δqj is be arbitrary, and the set of variations j = 1, 2, ... must be
independent, not only must the dot product equal zero, but each term within the inner product
must be zero. Therefore, the term on the left of the inner product must evaluate to the zero-vector.
This is an important concept in the principle of virtual work and in the calculus of variations. It’s
application results in the matrix equations of motion,
[M ] q̈(t) + [KE ] q(t) − [KG ] q(t) = f (t)
where the j, k term of the mass matrix is,
Z
Mjk = ρA(x) ψj (x) ψk (x) dx
l
the j, k term of the elastic stiffness matrix is,
Z
KEjk = EI(x) ψj′′ (x) ψk′′ (x) dx
l
the j, k term of the geometric stiffness matrix is,
Z
KGjk = P (x) ψj′ (x) ψk′ (x) dx
l
and the j-th element of the forcing vector is the inner product of the forcing with the j-th shape
function, Z
fj = f (x) · ψj (x) dx + Fi · ψj (xi )
X
l i
From the above relations, it is clear that Mij = Mji (the mass matrix is symmetric), Kij = Kji (the
stiffness matrices is symmetric), and that [M ] and [K] are positive definite, provided that the set of
shape functions are linearly independent.

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 13

6 Choice of Shape Function

In the set of shape functions described by
2j − 1
 
ψj (x) = sin πx
2
ψj (0) = 0 and ψ ′ (1) = 0. The figure below shows a set of the first four shape functions (j ∈ (1, 2, 3, 4))
on the left and the set of the first seven shape functions on the right (j ∈ (1, 2, 3, 4, 5, 6, 7)).
The curves in black show examples of weighted sums of the basis functions.
v(x, t) = q1 (t)ψ1 (x) + q2 (t)ψ2 (x) + ... + qN (t)ψN (x)
in which the “weights” correspond to the time-dependent generalized coordinates qj (t). So one may
think of the black curves as snapshots of vibrational shapes taken at various instances in time.
Note that since all the shape functions satisfy ψj (0) = 0 and ψ ′ (1) = 0, then so must the weighted
sum of those shape functions.
Note also that the use of a larger set of shape functions permits more complicated vibrational shapes.

1 1

0.5 0.5
ψ (x)

ψ (x)

0 0
j

-0.5 -0.5

-1 -1
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1
x/L x/L

q(t) = [ -0.020 , -7.861 , 7.258 , -2.870 ] q(t) = [ 9.127 , -8.107 , -6.644 , 9.382 , 1.131 , 7.799 , 7.059 ]
20 30
v(x,t) = Σ ψ) (x) q (t)

v(x,t) = Σ ψ) (x) q (t)

15 20
j

10
10
j

5
0
0
-10
-5
-10 -20
-15 -30
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1
x/L x/L

q(t) = [ 2.737 , -1.605 , 3.967 , 1.696 ] q(t) = [ -1.567 , 7.177 , -2.190 , -0.388 , -6.156 , 7.464 , -2.312 ]
8 20
v(x,t) = Σ ψ) (x) q (t)

v(x,t) = Σ ψ) (x) q (t)

6 10
j

4
0
j

2
-10
0
-2 -20

-4 -30
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1
x/L x/L

q(t) = [ 6.217 , 0.389 , -2.074 , -8.277 ] q(t) = [ 4.477 , -2.601 , -1.960 , -4.412 , -7.338 , 3.013 , 7.568 ]
15 20
v(x,t) = Σ ψ) (x) q (t)

v(x,t) = Σ ψ) (x) q (t)

10
j

10
5
j

0
0
-10
-5

-10 -20
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1
x/L x/L

cbnd H.P. Gavin October 10, 2022

14 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

It is essential that selected shape functions is kinematically admissible with respect to the essential
boundary conditions (the structural supports). The analytically “correct” shape function is both
kinematically admissible and satisfies equilibrium. Equations of motion resulting from the use of
kinematically admissible shape functions that do not satisfy equilibrium will provide approximate
solutions, which, in many cases are within the errors implied by other fundamental assumptions.
The true equations of motion for a particular system are unique. The principle of virtual displace-
ments provides a means to derive approximate equations of motion. The accuracy of the PVD
approximation depends on the set of shape functions used, [ψ1 (x), ..., ψN (x)].
Since v(x, t) has units of length, if ψ(x) is unitless, then the coordinate q(t) must have a unit of
length, and if ψ(x) has units of length, then q(t) is unitless (like a rotation).

7 Examples
7.1 Example 1: a single generalized coordinate, choice of two shape functions

Consider the vibration of a cantilever beam with a point end-mass (assuming that the rotatory inertia
of the end mass is negligible). And consider the choice between two similar shape functions,

3 x 1 x
2 3
πx
   
ψ(x) = − or ψ(x) = 1 − cos
2 L 2 L 2L
The cubic shape function is the static displacement of a cantilever beam with a concentrated tip load,
which would seem like a reasonable guess for the deformed shape in this problem. Note that the
static displacements satisfy equilibrium for a static load; they do not necessarily satisfy equilibrium
for a dynamic load.
The (1-cosine) shape function is a reasonable guess, since it is smooth and satisfies the essential
boundary conditions. But the (1-cosine) does not satisfy internal equilibrium for static or dynamic
loads.
Regardless of the choice of shape function, the internal virtual work is the work of inertial forces
moving through collocated virtual displacements plus the work of internal bending moments moving
through virtual rotations.
Z L Z L
δWINT = mv̈(x, t) δv(x, t) dx + M v̈(L, t) δv(L, t) + EIv ′′ (x, t) δv ′′ (x, t) dx
0 0
Z L Z L
= m (ψ(x)) dx q̈(t) δq(t) + M (ψ(L)) q̈(t) δq(t) + EI
2 2
(ψ ′′ (x)2 ) dx q(t) δq(t)
0 0

and the work of the external force moving through its collocated virtual rotation is

δWEXT = F (t)δv(L, t) = F (t) ψ(L) δq(t)

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 15

δ v(t,x)

v(t,x)
M

EI, m
F(t)

1
0.8
0.6
ψ(x)

0.4
cubic
0.2 1-cosine
0
0 0.2 0.4 0.6 0.8 1
x/L
3
2.5 cubic
2 1-cosine
ψ’’(x)

1.5
1
0.5
0
0 0.2 0.4 0.6 0.8 1
x/L
Figure 3. Two similar shape functions assumed for the displaced shape of a vibrating cantilever
beam with an end mass. Differences are clearer in the curvature ψ ′′ (x) (bending moments
M (x) = EIψ ′′ (x)) of the system. The cubic shape function corresponds to the triangular-
shaped bending moment of a cantilever beam with a static end-load. Note that the neither the
cubic shape function nor the 1-cosine shape function are exactly correct for this problem. The
true shape function would depend on the ratio mL/M .

cbnd H.P. Gavin October 10, 2022

16 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

Setting δWINT = δWEXT , factoring out the arbitrary virtual coordinate δq(t), and solving the integrals
gives for each of the candidate shape functions gives, for the cubic shape function,
33 3EI
 
mL + M q̈(t) + 3 q(t) = F (t)
140 L
and for the (1-cosine) shape function

3π − 8 π 4 EI
 
mL + M q̈(t) + q(t) = F (t)
2π 32L3
with natural frequency for the cubic shape function
v
3EI/L3
u
ωn =
u
t
33mL/140 + M

and for the (1-cosine) shape function

v
π 4 EI/(32L3 )
u
ωn =
u
t
(3π − 8)mL/(2π) + M

The shape function giving the lower natural frequency is more accurate.
As shown in the figure below, the cubic shape function gives a slightly more accurate dynamic model
as compared to the (1-cosine) function, by about one percent for mass ratios from 0 to 5, which could
be close to the difference of including or neglecting the rotational inertia of the end-mass.

1.8
cubic
1-cosine
1.7

1.6
ωn / (EI/(ML3))1/2

1.5

1.4

1.3

1.2

1.1
0 1 2 3 4 5
ratio of beam mass to end mass, mL / M

Figure 4. Natural frequencies of a beam with an end mass as a function of the ratio of the
beam mass to the end mass, using two choices for the shape function.

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 17

7.2 Example 2: a single generalized coordinate

In this example, the essential boundary conditions are v(t, 0) = 0 and v ′ (t, 0) = 0, so any shape
function used in this problem must also satisfy ψk (0) = 0 and ψk′ (0) = 0. In this first example, we
will consider a single (dimensionless) shape function, such as, ψ(x) = (x/L)2 , ψ(x) = (x/L)3 , or
ψ( x) = 1 − cos(πx/(2L)). Just to keep this simple for now, we choose ψ(x) = (x/L)3 . Forces and
associated virtual displacements are tabulated below.

δ v(t,x)

F(t) v(t,x)
EI, m M
P

c k f(t,x)

x=0 x=a x=b x=L

Element Real Internal Force Virtual Internal Displacement
M
M v̈(L, t) = M ψ(L)q̈(t) = M q̈(t) δv(L, t) = ψ(L)δq(t) = δq(t)

cv̇(a, t) = cψ(a)q̇(t) = c(a/L)3 q̇(t) δv(a, t) = ψ(a)δq(t) = (a/L)3 δq(t)

kv(b, t) = kψ(b)q(t) = k(b/L)3 q(t) δv(b, t) = ψ(b)δq(t) = (b/L)3 δq(t)

EI, m

EIv ′′ (x, t) = EIψ ′′ (x)q(t) = EI · 6x/L3 · q(t) δv ′′ (x, t) = ψ ′′ (x)δq(t) = 6x/L3 · δq(t)
EI, m

mv̈(x, t) = mψ(x)q̈(t) = m(x/L)3 · q̈(t) δv(x, t) = ψ(x)δq(t) = (x/L)3 δq(t)

Real External Force Virtual Displacement
F(t)
F (t) δv(a, t) = ψ(a)δq(t) = (a/L)3 δq(t)
f(t,x)

f (x, t) δv(x, t) = ψ(x)δq(t) = (x/L)3 δq(t)

P
P v ′ (x, t)δv ′ (x, t) = 9x4 /L6 q(t) δq(t)

Equating the work of real internal forces moving through internal virtual displacements, with real

cbnd H.P. Gavin October 10, 2022

18 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

external forces moving through collocated virtual displacements,

Z L Z L
M q̈ δq + c((a/L)3 )2 q̇ δq + k((b/L)3 )2 q δq + EI((6x/L3 ))2 dx q δq + m((x/L)3 )2 dx q̈ δq
0 0
Z L Z L
= F (t)(a/L) δq +
3
f (x, t)(x/L) dx δq +
3
P (9x4 /L6 ) dx q δq
b 0

Evaluating the definite integrals, factoring out the (arbitrary) virtual coordinate δq, specifying that
the distributed dynamic force is uniform with intensity fo , and grouping terms, the equation of
motion for this system is
 !6 
1 9  1 L4 − b4
6 3
a b EI a
   
M + mL q̈(t) + c q̇(t) + k + 12 3 − P q(t) = F (t) + fo (t)
6 L L L 5L L 4 L3

Note that this equation of motion is dimensionally homogeneous (as it should be).
The natural frequency of this system is
v
u  6
b
uk L
u
+ 12 EI
L3
− 9
5L
P
ωn = t
M + 61 mL

In this equation the term (9P q(t))/(5L) is moved to the left hand side of the equation, as it is a
function of position q(t). The coefficient (9P )/(5L) is called the geometric stiffness of this system.
The negative sign on this term shows that the axial compressive force P is destabilizing for this
system. Under the condition !6
b EI 9
k + 12 3 − P =0
L L 5L
the natural frequency would go to zero, and the system would buckle. So the critical axial buckling
load for the system is  
!6
b EI 5L
 
Pcr = k + 12 3
9
 
L L

Dynamical responses of complex systems require complex mathematical descriptions. The simple
approximation v(x, t) = (x/L)3 q(t) used here could be passable for a simple cantilever beam. But in
this example if the spring stiffness k were much higher than EI/L3 the dynamic response at x = b
would have a very small amplitude compared to responses the domains x < b and x > b. This kind
of response is not captured by the approximation ψ(x) = (x/L)3 . In fact, the nature of the free
dynamic response in systems such as the one in this example depends on the relative values of the
physical parameters, EI/L3 , M g/L, mg, P/L, k, etc. More complex mathematical models for v(x, t)
are required to describe the dynamic responses of complex systems such as this. The next example
shows an extension in this direction, in which v(x, t) is modeled by the superposition of two shape
functions, and two generalized coordinates.

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 19

7.3 Example 3: the same example with two shape functions and two generalized coordinates

In this example, the displaced shape is expressed as the sum of two (independent and kinematically
admissible) shape functions, ψ1 (x) and ψ2 (x)

3 x 1 x
" 2 3 # "  3 2 #
x x
  
v(x, t) = − q1 (t) + 8 −7 q2 (t)
2 L 2 L L L

1
ψ1(x/L) = (3/2)(x/L)2 - (1/2)(x/L)3
ψ2(x/L) = 8(x/L)3 - 7(x/L)2
shape functions: ψ1(x/L), ψ2(x/L)

0.5

-0.5

-1
0 0.2 0.4 0.6 0.8 1
x/L

Generalized coordinates associated with dimensionless shape functions have the same physical di-
mensions as the response variables, which is generally desirable. Shape functions that resemble the
actual dynamic responses correspond to more realistic dynamic models. Actual dynamic responses
must adhere to essential and natural boundary conditions. So as a first requirement, shape function
approximations must adhere to the essential boundary conditions. Shape functions that also adhere
to the natural boundary conditions correspond to more realistic models. Mass, and stiffness matrices
derived from sets of linearly independent shape functions are positive definite (assuming the system
has no rigid body modes). Mass and/or stiffness matrices derived from sets of mutually orthogo-
nal shape functions are numerically well conditioned. Because of this, models derived from sets of
mutually orthogonal shape functions are more precise over a broader frequency range.
In this example, ψ1 (x) corresponds to the static deflection of a cantilever beam with a point load at
x = L; ψ2 (x) has an inflection point and a zero-crossing.
The application of the principle of virtual displacements in which the responses are an expansion of
n (admissible and linearly independent) shape functions result in n dimensional matrix equations of
motion. Examples of mass and stiffness matrices for higher dimensional approximations are given in
equations (13), (14), (15), and (16). This problem is slightly more complex as it involves a spring, a
damper, and a concentrated mass.

cbnd H.P. Gavin October 10, 2022

20 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

Applying the principle of superposition, expressions for the internal and external virtual work corre-
sponding to each of these various components may be taken individually.
Work of the inertial force of the distributed mass of the beam, mv̈(x, t)dx, moving through virtual
displacements δv(x, t)
X X Z 
δWI = m(x) ψj (x) ψk (x) dx q̈k (t) δqj (t)
j k l

Work of the inertial force of the point mass of the beam, M v̈(L, t) moving through virtual displace-
ments δv(L, t) Z 
δWI = M δ(x − L) ψj (x) ψk (x) dx q̈k (t) δqj (t)
XX

j k l

Work of the bending moments distributed along the beam, EIv ′′ (x, t) moving through virtual rota-
tions distributed along the beam δv ′′ (x, t) dx
X X Z 
δWI = EI(x) ψj′′ (x) ψk′′ (x) dx qk (t) δqj (t)
j k l

Work of the spring force, kv(b, t), moving through virtual displacements δv(b, t)

X X Z 
δWI = kδ(x − b) ψj (x) ψk (x) dx qk (t) δqj (t)
j k l

Work of the damper force, cv̇(a, t), moving through virtual displacements δv(a, t)

X X Z 
δWI = cδ(x − a) ψj (x) ψk (x) dx q̇k (t) δqj (t)
j k l

Work of the dynamic point force, F (t), moving through virtual displacements δv(a, t)

X Z 
δWE = F (t)δ(x − a) ψj (x) dx δqj (t)
j l

Work of the dynamic distributed Force, f (t) dx, moving through virtual displacements δv(x, t)

" #
X Z L
δWE = f (x, t) ψj (x) dx δqj (t)
j b

Work of the uniform axial force, P , moving through distributed virtual end shortening v ′ (x, t)δv ′ (x, t)dx

X X Z 
δWE = P (x) ψj′ (x) ψk′ (x) dx qk (t) δqj (t)
j k l

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 21

Each j, k term within the square brackets corresponds to the j, k term of a mass, damping, or stiffness
matrix. In these derivations, δ(x − a) is the Dirac delta function, which has the defining property,
Z
g(x)δ(x − a) dx = g(a)
l

The evaluation of the associated derivatives and integrals can be easily carried out using symbolic
manipulation software like Mathematica, Maple, or Wolfram-α.

: login-teer-12 Sun Sep 02 14:11:35 ˜

: maple ## EVALUATE INTEGRALS ....................................
|\ˆ/| Maple 2017 (X86 64 LINUX)
._|\| |/|_. Copyright (c) Maplesoft, a division of Waterloo Maple Inc. 2017
\ MAPLE / All rights reserved. Maple is a trademark of
<____ ____> Waterloo Maple Inc.
| Type ? for help.
## MASS MATRIX TERMS .....................................
## INPUT THE SHAPE FUNCTION EQUATIONS ....................
> m11 := int(m*p1*p1,x=0..L);
> p1 := (3/2)*(x/L)ˆ2 - (1/2)*(x/L)ˆ3; 33 m L
2 3 m11 := ------
3 x x 140
p1 := ---- - ----
2 3 > m12 := int(m*p1*p2,x=0..L);
2 L 2 L 37 m L
m12 := - ------
> p2 := 8*(x/L)ˆ3 - 7*(x/L)ˆ2; 420
3 2
8 x 7 x > m22 := int(m*p2*p2,x=0..L);
p2 := ---- - ---- 29 m L
3 2 m22 := ------
L L 105

## EVALUATE DERIVITIVES .................................. > M11 := eval(M*p1*p1,x=L);

M11 := M

> dp1 := diff(p1,x); > M12 := eval(M*p1*p2,x=L);

2 M12 := M
3 x 3 x
dp1 := --- - ---- > M22 := eval(M*p2*p2,x=L);
2 3 M22 := M
L 2 L

> ddp1 := diff(dp1,x); ## STIFFNESS MATRIX TERMS ................................

3 3 x
ddp1 := ---- - --- > EI11 := int(EI*ddp1*ddp1,x=0..L);
2 3 3 EI
L L EI11 := ----
3
> dp2 := diff(p2,x); L
2
24 x 14 x > EI12 := int(EI*ddp1*ddp2,x=0..L);
dp2 := ----- - ---- 3 EI
3 2 EI12 := ----
L L 3
L
> ddp2 := diff(dp2,x);
48 x 14 > EI22 := int(EI*ddp2*ddp2,x=0..L);
ddp2 := ---- - ---- 292 EI
3 2 EI22 := ------
L L 3
L

cbnd H.P. Gavin October 10, 2022

22 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

> k11 := eval(k*p1*p1,x=b); ## EXTERNAL FORCING TERMS ................................

> k11 := simplify(k11);

4 2 > F1 := eval(F*p1,x=a);
k b (3 L - b)
k11 := --------------- > F1 := simplify(F1);
6 2
4 L F a (3 L - a)
F1 := --------------
> k12 := eval(k*p1*p2,x=b); 3
2 L
> k12 := simplify(k12);
4 > F2 := eval(F*p2,x=a);
k b (3 L - b) (8 b - 7 L)
k12 := -------------------------- > F2 := simplify(F2);
6 2
2 L F a (8 a - 7 L)
F2 := ----------------
> k22 := eval(k*p2*p2,x=b); 3
L
> k22 := simplify(k22);
4 2 > f1 := int(fo*p1,x=b..L);
k b (7 L - 8 b)
k22 := ----------------- > f1 := simplify(f1);
6 4 3 4
L fo (3 L - 4 L b - b )
f1 := -----------------------
3
## DAMPING MATRIX TERMS ................................. 8 L

> c11 := eval(c*p1*p1,x=a);

> f2 := int(fo*p2,x=b..L);
> c11 := simplify(c11);
4 2 > f2 := simplify(f2);
c a (3 L - a) 4 3 4
c11 := --------------- fo (L - 7 L b + 6 b )
6 f2 := - ---------------------------
4 L 3
3 L
> c12 := eval(c*p1*p2,x=a);

> c12 := simplify(c12); ## GEOMETRIC STIFFNESS TERMS .............................

4
c a (3 L - a) (8 a - 7 L) > P11 := int(P*dp1*dp1,x=0..L);
c12 := -------------------------- 6 P
6 P11 := ---
2 L 5 L

> c22 := eval(c*p2*p2,x=a); > P12 := int(P*dp1*dp2,x=0..L);

41 P
> c22 := simplify(c22); P12 := ----
4 2 20 L
c a (8 a - 7 L)
c22 := ----------------- > P22 := int(P*dp2*dp2,x=0..L);
6 188 P
L P22 := -----
15 L

cbnd H.P. Gavin October 10, 2022

Virtual Displacements in Structural Dynamics 23

The resulting equations of motion in terms of generalized coordinates, q1 (t) and q2 (t) are
33
+ M − 420
37
mL + M q̈1 (t)
  
140
mL
  
  
37
− 420 mL + M 29
105
mL +M q̈2 (t)

 
(3L−a)2 (3L−a)(8a−7L)
q̇1 (t)

4
ca  4 2
+
 
 
L6   
(3L−a)(8a−7L)
2
(8a − 7L)2 q̇2 (t)

 
(3L−b)2 (3L−b)(8b−7L)
q1 (t)

4
kb  4 2
+ 6 
 

L   
(3L−b)(8b−7L)
2
(8b − 7L)2 q2 (t)

3 3 q1 (t)
  
EI 
+  
L3
  
3 292 q2 (t)

a2 (3L−a)
 
3L4 −4Lb3 −b4
q1 (t) F (t)
 6 41
  
P 2L3 8L3
 5 20
=  (17)
 
−    
L
   
q2 (t)
 
41
20
188
15
a2 (8a−7L) 4
− L −7Lb
3 +6b4
fo (t)
L3 3L3

With the scaling of the dimensionless shape functions, ψ1 (L) = ψ2 (L) = 1, q1 (t) and q2 (t) are the
values of v(L, t) corresponding to ψ1 (x) and ψ2 (x). With the dimensionless formulation of the shape
functions, every term in this equation has units of force.
This example is an introduction to methodologies that are invoked later in the course.
For more complex geometries, for example, beams with tapered sections, the derivation can become
very complex, and the analysis is more easily carried out numerically.
The same example, computed with matlab, is:

%% DEFINE NUMERICAL VALUES FOR CONSTANTS

EI = 1e7; % flexural rigidity N.mˆ2

k = 1e2; % concentrated stiffness N/m
m = 1e0; % distributed mass kg/m
M = 1e1; % lumped mass kg
c = 0.1; % spring damping rate N/m/s
L = 10; % overall length m
a = 3; % location of damper m
b = 5; % location of spring m

dx = 0.01; % increment of length along the beam

x = [ 0 : dx : L ]; % x-axis
xa = round(a/dx);
xb = round(b/dx);

cbnd H.P. Gavin October 10, 2022

24 CEE 541. Structural Dynamics – Duke University – Fall 2020 – H.P. Gavin

%% INPUT THE SHAPE FUNCTION EQUATIONS ....................

p1 = (3/2)*(x/L).ˆ2 - (1/2)*(x/L).ˆ3;

p2 = 8*(x/L).ˆ3 - 7*(x/L).ˆ2;

%% EVALUATE DERIVITIVES ..................................

dp1 = cdiff(p1)/dx;
ddp1 = cdiff(dp1)/dx;

dp2 = cdiff(p2)/dx;
ddp2 = cdiff(dp2)/dx;

%% EVALUATE INTEGRALS ....................................

%% MASS MATRIX TERMS .....................................

m11 = trapz(m*p1.*p1)*dx;
m12 = trapz(m*p1*p2)*dx;
m22 = trapz(m*p2*p2)*dx;

M11 = M*p1(end)*p1(end);
M12 = M*p1(end)*p2(end);
M22 = M*p2(end)*p2(end);

%% STIFFNESS MATRIX TERMS ................................

EI11 = trapz(EI*ddp1.*ddp1)*dx;
EI12 = trapz(EI*ddp1.*ddp2)*dx;
EI22 = trapz(EI*ddp2.*ddp2)*dx;

k11 = k*p1(xb)*p1(xb);
k12 = k*p1(xb)*p2(xb);
k22 = k*p2(xb)*p2(xb);

%% DAMPING MATRIX TERMS .................................

c11 = c*p1(xa)*p1(xa);
c12 = c*p1(xa)*p2(xa);
c22 = c*p2(xa)*p2(xa);

%% EXTERNAL FORCING TERMS ................................

F1 = F*p1(xa);
F2 = F*p2(xa);

f1 = trapz(fo*p1(xb:L))*dx;
f2 = trapz(fo*p2(xb:L))*dx;

%% GEOMETRIC STIFFNESS TERMS .............................

P11 = trapz(P*dp1.*dp1)*dx;
P12 = trapz(P*dp1.*dp2)*dx;
P22 = trapz(P*dp2.*dp2)*dx;

cbnd H.P. Gavin October 10, 2022