VIBRATION CONTROL OF CONTINUOUS SYSTEMS
USING BOUNDARY CONSTRAINT
By
Assaad Awad Abbass

A DISSERTATION
Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of
DOCTOR OF PHILOSOPHY
Mechanical Engineering
2012

ABSTRACT
VIBRATION CONTROL OF CONTINUOUS SYSTEMS
USING BOUNDARY CONSTRAINT
By
Assaad Awad Abbass
Vibration suppression in continuous systems, namely strings and membranes, using
boundary constraints was investigated. A problem of string vibration, where the string is
subjected to a fixed boundary constraint, is studied first. The boundary constraint is in the
form of a smooth obstacle and the string is assumed to wrap and unwrap around the obstacle
during vibration. Assuming that wrapping of the string around the obstacle results in loss of
kinetic energy due to inelastic collision and unwrapping conserves energy, the behavior of the
string was investigated for different modes of oscillation and different obstacle shapes. The
loss of energy is found to be greater for higher modes of oscillation and for obstacles that
induce greater length of wrapping. Next, the vibration of the string subjected to a moving
constraint near one boundary is investigated. The constraint is applied by a scabbard, which
moves a small distance along the mean position of the string. The constraint is removed by
moving the scabbard back to its original position and the change in energy of the string is
investigated for different values of scabbard travel distance and time of application of the
constraint. Unlike the fixed obstacle at the boundary, simulation results show that the energy
of the string can increase or decrease depending on the time of application of the constraint.
Based of this finding, a semi-active control strategy is proposed for vibration suppression.
The control strategy is verified experimentally by substituting the scabbard-like actuator by
a pair of solenoids. The experimental results show good match with results obtained through
simulations. The control strategy developed for the string is extended to a circular

membrane. It is assumed that the membrane is fixed at its outer boundary and a zero displacement circular areal constraint is sequentially applied and removed. To investigate the
effect of constraint application and removal on the energy of the membrane, the dynamics
of constrained membranes is investigated. For arbitrary size and location of the constraint,
the orthogonality of distinct modes is mathematically established and a procedure for accurate computation of the eigenfrequencies and mode shapes is presented. Assuming that the
constraint is applied and removed instantaneously at arbitrary time intervals, the change in
total energy is investigated for different sizes and locations of the constraint and for different
times of application of the constraint. Similar to the vibrating string with the scabbard-like
actuator, the results show that energy of the membrane can decrease or increase depending
on the time of application of the constraint. Three different semi-active control strategies
are presented to the vibration suppression of circular membrane. These control strategies
have the potential for use in large space structures where high dimensional tolerances are
required.

“ The wind gasps with the midday heat,
Like a nightmare in the late afternoon
And on the masts, it continues to fold, to spread for departure
The gulf is crowded with them–laborers roaming the seas
Barefoot, half-naked
And on the sand, by the gulf
A stranger sat–a baffled vision wanders the gulf
Destroying the pillars of light with the rising wail
Higher than the torrents roaring foam, than the clamor
A voice thunders in the abyss of my bereaved soul: Iraq
Like the crest rising, like a cloud, like tears to the eyes
The wind cries to me: Iraq
The wave howls at me: Iraq. Iraq. Nothing but Iraq
The sea is as wide as can be, and you are as distant
The sea is between you and me: Oh Iraq“
A Stranger by the Gulf: by Badr Shakir Al-Sayyab (1953)

To my family and my home country Iraq

iv

ACKNOWLEDGMENTS

In the name of Allah, the Most Gracious and the Most Merciful, all praises to Allah for the
strengths and His blessing in completing this work. I would like to express my most sincere
gratitude to my advisor, Prof. Ranjan Mukherjee for his excellent guidance, consistent
encouragement, I am greatly indebted to him for his invaluable inspiration. Special thanks
go to Prof. Alejandro R. Diaz , Prof. Hassan Khalil and Prof. Steven W. Shaw for their
willingness of being in my committee members despite their busy schedules. Also, I would
like to thank my friends in the control lab for their help and support.
My gratitude extends to my parents for their support and encouragement. Their prayers
for me during the course of this study have brought many blessings into my life. Special
gratitude and appreciation are expressed to my family for their patience and support during
my study. Last but not least, I also would like to thank my uncle Jawad Al-Sahlani and my
two uncles families in United States for their help and support during my study.

v

TABLE OF CONTENTS

List of Tables . . . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . . .

ix

. . . . . . . . . . . . . .

. . . . . . . . . . . . . . . . . .

xi

List of Figures

Chapter 1 Introduction . . . . . . . . . . . . . .
1.1 Motivation and Objectives . . . . . . . . . . .
1.2 Vibration Control of Strings . . . . . . . . . .
1.2.1 Passive Control of Vibrating Strings . .
1.2.2 Active Control of String vibration . . .
1.3 Dynamics and Vibration Control of Membrane
Chapter 2
2.1
2.2
2.3

2.4
2.5

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

Vibration of a String Wrapping and Unwrapping
Obstacle . . . . . . . . . . . . . . . . . . . . . . . . . .
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Problem Statement and Assumptions . . . . . . . . . . . . . .
Analytical Model . . . . . . . . . . . . . . . . . . . . . . . . .
2.3.1 Boundary conditions and general solution . . . . . . .
2.3.2 Wrapping of the string . . . . . . . . . . . . . . . . . .
2.3.3 Unwrapping of the string . . . . . . . . . . . . . . . . .
Numerical Simulations . . . . . . . . . . . . . . . . . . . . . .
Effect of Change in Slope of Obstacle . . . . . . . . . . . . . .

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

1
1
3
3
5
7

Around an
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .

12
12
13
16
16
20
22
23
29

Chapter 3 Vibration Control of a String Using A Scabbard-Like Actuator
3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2 Problem Statement and Assumptions . . . . . . . . . . . . . . . . . . . . . .
3.3 Effect of Applying Constraint on the Dynamics of the String . . . . . . . . .
3.3.1 Equation of motion after application of the constraint . . . . . . . . .
3.3.2 Energetics of constraint application . . . . . . . . . . . . . . . . . . .
3.4 Effect of Removing Constraint on the Dynamics of the String . . . . . . . . .
3.4.1 Equation of motion after removal of the constraint . . . . . . . . . . .
3.4.2 Energetics of constraint removal . . . . . . . . . . . . . . . . . . . . .
3.5 Simulation of One Cycle of Constraint Application and Removal . . . . . . .
3.6 Sequential Application and Removal of Constraints . . . . . . . . . . . . . .
vi

31
31
32
34
34
37
39
40
42
43
47

3.6.1
3.6.2
3.6.3

Dynamics of string with arbitrary initial conditions . . . . . . . . . .
Control strategy for vibration suppression . . . . . . . . . . . . . . .
Numerical simulations . . . . . . . . . . . . . . . . . . . . . . . . . .

Chapter 4
4.1
4.2

4.3

4.4
4.5

Vibration Control of a String Through Cyclic Application
Removal of Constraints: Experimental Verification . . . . .
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Cyclic Application and Removal of Constraint . . . . . . . . . . . . . .
4.2.1 Actuation and Sensing Mechanism . . . . . . . . . . . . . . . .
4.2.2 Mathematical Modeling . . . . . . . . . . . . . . . . . . . . . .
4.2.2.1 Linear model of vibrating string with damping . . . . .
4.2.2.2 Effect of constraint application . . . . . . . . . . . . .
4.2.2.3 Effect of constraint removal . . . . . . . . . . . . . . .
Experimental Setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.3.1 Coil Extension Spring - A Substitute for the String . . . . . . .
4.3.2 Internal Damping of the String . . . . . . . . . . . . . . . . . .
4.3.2.1 Estimation of Damping Ratio . . . . . . . . . . . . . .
4.3.2.2 Variation of Damping Coefficient with Length . . . . .
Simulation and Experimental Results on Vibration Suppression . . . .
Additional Investigations . . . . . . . . . . . . . . . . . . . . . . . . . .
4.5.1 Efficacy of Vibration Suppression - A Parametric Study . . . . .
4.5.2 Excitation of Higher Modes . . . . . . . . . . . . . . . . . . . .

Chapter 5
5.1
5.2

5.3
5.4

5.5

Dynamics of a Circular Membrane with an Eccentric
Areal Constraint: Analysis and Accurate Simulations
Problem Statement and Assumptions . . . . . . . . . . . . . . . .
Background - Modes of Vibration . . . . . . . . . . . . . . . . . .
5.2.1 General Solution . . . . . . . . . . . . . . . . . . . . . . .
5.2.2 Boundary Conditions . . . . . . . . . . . . . . . . . . . . .
5.2.3 Even or Symmetric Modes . . . . . . . . . . . . . . . . . .
5.2.4 Odd or Antisymmetric Modes . . . . . . . . . . . . . . . .
Orthogonality of Modes . . . . . . . . . . . . . . . . . . . . . . .
Accurate Computation of Eigenfrequencies and Mode Shapes . . .
5.4.1 Computational Algorithm . . . . . . . . . . . . . . . . . .
5.4.2 Numerical Results . . . . . . . . . . . . . . . . . . . . . .
Simulation of Constrained Membrane Dynamics . . . . . . . . . .
5.5.1 Computation of Modal Coefficients . . . . . . . . . . . . .
5.5.2 Numerical Results Verifying Orthogonality . . . . . . . . .
5.5.3 Simulation of Free Vibration . . . . . . . . . . . . . . . . .

Chapter 6
6.1
6.2

Vibration Control of a Circular Membrane
Constraint . . . . . . . . . . . . . . . . . . . . . .
Introduction . . . . . . . . . . . . . . . . . . . . . . . . .
Problem Statement and Assumptions . . . . . . . . . . .
vii

Using
. . . .
. . . .
. . . .

47
50
55

and
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .
. . .

59
59
60
60
62
62
64
67
68
68
71
71
72
74
78
78
83

Circular
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .
. . . . . .

85
86
88
88
90
91
94
96
101
101
103
113
113
117
118

An Areal
. . . . . . . 124
. . . . . . . 124
. . . . . . . 125

6.3

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

128
128
129
131
132
135
135
140
141
145
145
147
147
152
159
164

Chapter 7 Conclusion Remarks and Future Work . . . . . . . . . . .
7.1 Vibration of a String Wrapping and Unwrapping Around an Obstacle
7.2 Semi-Active Control for String Vibration . . . . . . . . . . . . . . .
7.3 Dynamics and Vibration Control of Membrane . . . . . . . . . . . . . . .

.
.
.
.

. . .
. . .
....
. . .

171
172
173
175

APPENDICES . . . . . . . . . . . . . . . . . . . . . . .
Appendix A Vibration of Circular Membrane: Preliminary . . . .
A.1 Bessel Functions . . . . . . . . . . . . . .. . . . . . . . .
A.2 Vibration of Circular Membranes . . . . .. . . . . . . . .
Appendix B Change of Summation Indices . . . . . . . . . . . .
Appendix C Energy of Circular Membrane . . . . . . . . . . .

.
.
.
.
.
.

.
.
.
.
.
.

179
179
179
183
188
190

6.4

6.5
6.6

Effect of Applying Constraint on the Dynamics of the Membrane
6.3.1 Equation of motion after application of the constraint . . .
6.3.2 Initial Conditions . . . . . . . . . . . . . . . . . . . . . . .
6.3.3 Membrane Energy Before Application of Constraint . . . .
6.3.4 Membrane Energy After Application of Constraint . . . .
Effect of Removing Constraint on the Dynamics of the Membrane
6.4.1 Equation of motion after removal of the constraint . . . . .
6.4.2 Membrane Energy After Removal of Constraint . . . . . .
Simulation of One Cycle of Constraint Application and Removal .
Sequential Application and Removal of Constraints . . . . . . . .
6.6.1 Dynamics of membrane with arbitrary initial conditions . .
6.6.2 Control strategies for vibration suppression . . . . . . . . .
6.6.3 Simple Application and Removal of Constraint . . . . . . .
6.6.4 Radial Constraint Application . . . . . . . . . . . . . . . .
6.6.5 Circumferential Constraint Application . . . . . . . . . .
6.6.6 Numerical simulations . . . . . . . . . . . . . . . . . . . .

.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

. .
..
. .
. .
. .
. .

.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.

.
.
.
.
.
.

.
.
.
.
.
.

Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 193

viii

LIST OF TABLES

Table 2.1

Percentage energy loss over one cycle of oscillation and xk for different
¯
values of E0 and three modes of oscillation, all with R = 1 m . . . .

23

Percentage energy loss over one cycle of oscillation and xk for obsta¯
cles of different shapes and sizes, all with E0 = 0.50 J . . . . . . . .

28

Percentage energy loss over one cycle of oscillation and xk for two
¯
modes of oscillation with different values of θ . . . . . . . . . . . . .

30

Table 3.1

Comparison of values of Ec and Er for x0 = 0.1l and tr = 0.8 sec .

44

Table 3.2

Comparison of percentage energy loss of the string for cases (A) and
(B) over 12 cycles of constraint application and removal . . . . . . .

56

Table 4.1

Specifications of coil extension spring . . . . . . . . . . . . . . . . .

69

Table 4.2

Experimentally determined natural frequencies of the spring/string
in rad/s . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

71

Effect of three different parameters: (x0 /ℓ), tr , and N on the settling
time . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

79

Table 2.2

Table 2.3

Table 4.3

Table 5.1

First eight non-dimensional frequencies of antisymmetric modes of
constrained membrane with a/R = 0.25 and d/a = 1.0. The results in
parenthesis are from [1], [2], and [3], and the star superscript indicates
frequencies with accurate mode shapes computed using our algorithm 105
ix

Table 5.2

First five non-dimensional frequencies of symmetric modes of constrained membranes and values of N used to compute the mode
shapes. For a given set of a/R and d/a values, the first row contains values that correspond to accurate mode shapes. The second
row contain values that are closest to those presented in the literature
- the eigenfrequencies are underlined to indicate that the corresponding mode shapes are inaccurate. The results from the literature [1],
[2], [3] are presented in the following rows with eigenfrequencies in
parenthesis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107

Table 5.3

Checking orthogonality between the accurate and inaccurate third
symmetric mode and the accurate fourth symmetric mode of constrained membranes with a/R = 0.5 and d/a = 0.6, 0.8, 1.0. The
inaccurate third symmetric mode is underlined . . . . . . . . . . . . 117

Table 5.4

Number of integration points n required to achieve orthogonality between symmetric-symmetric and symmetric-antisymmetric modes of
constrained membranes with a/R = 0.1 and d/a = 0, 4, 8 . . . . . . 118

Table 6.1

Comparison of values of Er with a fixed value of Ec = 2.4455 J for
constrained membrane with d/R = 0.88, a/R = 0.12 and tr = 0.4 sec 143

Table 6.2

Comparison of percentage energy loss of the membrane for the three
different constraint applications strategies (A) simple (B) radial (C)
circumferential, over 10 constraint applications . . . . . . . . . . . . 167

Table A.1

Some values for γmn which is zeros for Bessel function Jm (γnm ). . 185

x

LIST OF FIGURES

Figure 1.1

Figure 2.1

Figure 2.2

Figure 2.3

Figure 2.4

Figure 2.5

Figure 2.6
Figure 3.1

The spartan 207 inflatable antenna, NASA [4] “For interpretation
of the references to color in this and all other figures, the reader is
referred to the electronic version of this dissertation“ . . . . . . . . .

3

A string vibrating against an obstacle is shown in three configurations: (a) the string has both potential and kinetic energy and is not
in contact with the obstacle, (b) the string has kinetic energy but no
potential energy and is not in contact with the obstacle, (c) the string
has both potential and kinetic energy and has wrapped around the
obstacle. In a wrapped configuration, (¯, y ) denotes the coordinate
x ¯
where the string breaks contact with the obstacle . . . . . . . . . . .
Small string segment ∆x wraps around obstacle after perfectly inelastic collision. In the magnified image, AB denotes the small string
segment of length ∆x that wraps around the obstacle over the region
AC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

20

Plot of percentage energy content of the string over three consecutive
cycles of vibration in (a) Mode 1, and (b) Mode 2. For both cases,
the initial energy of the string was E0 = 0.5 J . . . . . . . . . . . .

25

Exponential decay in the energy of a string wrapping and unwrapping
around an obstacle. The plots show energy decay for single-mode
vibration in the first four modes with E0 = 0.5 J . . . . . . . . . . .

26

A string vibrating against a circular obstacle in (a) Mode 1 (b) Mode
2, and (c) Mode 3 . . . . . . . . . . . . . . . . . . . . . . . . . . . .

27

A string vibrating against an obstacle. The obstacle is identical to
the one in Fig. 2.1 but rotated clockwise by angle θ about point O .

29

(a) A vibrating string (b) A zero-displacement constraint is applied
to the string at time t = tc over the length segment x ∈ [0, x0 ) using
a scabbard . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

33

xi

15

Figure 3.2
Figure 3.3

Figure 3.4

Figure 3.5

Figure 3.6

Figure 4.1
Figure 4.2

Figure 4.3

Figure 4.4

Figure 4.5

Figure 4.6

Figure 4.7

Figure 4.8

At time t = tr , the scabbard is instantaneously moved back to its
original position . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
(a) Position of the string at different time instants (b) Percentage
change in energy of the string due to application of the constraint for
different values of x0 and tc . . . . . . . . . . . . . . . . . . . . . .

40

46

The sensor measures the displacement of the string at a distance x0
from the tip of the scabbard. When this displacement is zero, the
scabbard moves incrementally by distance x0 . . . . . . . . . . . . .

50

Displacement of the string and percentage energy loss at different
time instants during the first cycle of constraint application . . . . .

57

Displacement and percentage energy loss of the string at the end of
the first cycle and after completion of 12 cycles of constraint application and removal for cases (A) and (B) . . . . . . . . . . . . . . . .

58

Actuation and sensing mechanism for application and removal of constraint . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
(a) Experimental setup showing the string between two fixed supports
and the actuation mechanism, (b) A close-up view of the actuation
mechanism and sensor . . . . . . . . . . . . . . . . . . . . . . . . . .

62

70

Free vibration of a string: (a) Experimental results for the string
in Fig. 4.2 (b) Simulation results obtained using a damping ratio
computed from the experimental results in (a) . . . . . . . . . . . .

72

Effect of changing the length of the string on the product ζω for three
different strings . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

73

Initial displacement of the string used for both simulations and experiments. The guide for generating repeatable initial conditions in
experiments can be seen in Fig. 4.2 . . . . . . . . . . . . . . . . . .

75

Results of vibration suppression obtained from (a) simulations and
(b) experiments . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

76

Results of vibration suppression obtained from (a) simulations and
(b) experiments, with (x0 /ℓ) = 0.05, tr = 0.5 s, and N = 12 . . . . .

80

Results of vibration suppression obtained from (a) simulations and
(b) experiments, with (x0 /ℓ) = 0.10, tr = 1 s, and N = 12 . . . . .

81

xii

Figure 4.9

Results of vibration suppression obtained from (a) simulations and
(b) experiments, with (x0 /ℓ) = 0.15, tr = 1.5 s, and N = 12 . . . . .

82

Fast Fourier-Transform results showing multiple modes of vibration
of the (a) unconstrained string and (b) constrained string . . . . . .

83

Figure 5.1

A circular membrane with an inner circular areal constraint . . . . .

87

Figure 5.2

Shaded region showing the area of integration for the two cases where
the integration is done (a) first with respect to r and then with respect
to θ, and (b) first with respect to θ and then with respect to r . . .

99

Figure 4.10

Figure 5.3

The eigenfrequencies are the zero crossings of the determinant of the
matrix P c for symmetric modes, and P s for antisymmetric modes . 102

Figure 5.4

Algorithm for computing eigenfrequencies and mode shapes for symmetric modes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 104

Figure 5.5

c
Mode shapes k3 a for constrained membranes in Table 5.2 with a/R =
0.5 and (a) d/a = 0.6 (b) d/a = 0.8 (c) d/a = 1.0. For each of these
three cases, the inaccurate mode shapes are shown above the accurate
mode shapes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110

Figure 5.6

First symmetric mode and sixteenth symmetric mode of a constrained
membrane with a/R = 0.1 and d/a = 9.0. For both these cases,
the isometric view is placed above the top view and the value of N
required for accuracy is shown . . . . . . . . . . . . . . . . . . . . . 111

Figure 5.7

First five modes of a membrane with a/R = 0.1 and d/a = 2.0 (a)
Antisymmetric modes (b) Symmetric modes. For both, the isometric
view is placed next the top view . . . . . . . . . . . . . . . . . . . . 112
Algorithm for computing modal coefficients . . . . . . . . . . . . . . 116

Figure 5.8
Figure 5.9
Figure 5.10
Figure 5.11
Figure 5.12
Figure 5.13
Figure 6.1

Initial displacement of the membrane described in case 1 . . . . . .
Initial displacement of the membrane described in case 2 . . . . . .
Snapshots of the membrane described in case 1 during one second of
its vibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Snapshots of the membrane described in case 2 during one second of
its vibration . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Snapshots of the membrane described in section case 3 during one
cycle of its vibration . . . . . . . . . . . . . . . . . . . . . . . . . . .

119
119
121
122
123

Circular membrane with radius R . . . . . . . . . . . . . . . . . . . 126
xiii

Figure 6.2

Figure 6.3
Figure 6.4

(a) A vibrating membrane (b) A zero-displacement circular constraint
is applied to the membrane at time t = tc over a small region of the
membrane (c) Geometry of the membrane and the applied constraint
after constraint application . . . . . . . . . . . . . . . . . . . . . . 127
(a) A vibrating constrained membrane (b) The circular constraint is
removed from the membrane at time t ≃ tr . . . . . . . . . . . . . 136
Percentage change in total energy of the membrane due to application
of the constraint of different size of values of of d and to. (a) position
of the membrane at different time instants. (b) a/R = 0.1 (c) a/R =
0.12 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 144

Figure 6.5

Sensing mechanism for simple application and removal of constraint
(a) side view. (b) top view . . . . . . . . . . . . . . . . . . . . . . . 148

Figure 6.6

Flowchart for simple application and removal of constraint strategy . 151

Figure 6.7

Sensing mechanism and path of constraint motion for radial constraint motion (a) side view. (b) top view. . . . . . . . . . . . . . . 153

Figure 6.8

The coordinates of two different locations of a constraint moved from
location at d1 to location at d2 . . . . . . . . . . . . . . . . . . . . . 154

Figure 6.9

Flowchart for radial application of constraint strategy, for each p increment, new constraint location dp is defined and the eigenfrequencies and mode shapes are calculates at this location . . . . . . . . . 158

Figure 6.10

Sensing mechanism for circumferential constraint motion (a) back
side view. (b) top view . . . . . . . . . . . . . . . . . . . . . . . . . 160

Figure 6.11

The coordinates systems of two different locations of a constraint
moved in a circular path from location (d, 0) to (d, φ1 ) . . . . . . . 161

Figure 6.12

Flowchart for circumferential application of constraint strategy, for
each p increment, new constraint location (angle φp ) is defined and
the eigenfrequencies and mode shapes are calculates at this location

Figure 6.13
Figure 6.14

165

Snapshots of one cycle for simple constraint application and removal
strategy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 169
Displacement and percentage energy loss E% of the membrane at
the end of the first cycle, 6 cycles and after completion of 10 cycles
of constraint application and removal for cases; simple, radial and
circumferential applications of constraint . . . . . . . . . . . . . . . 170
xiv

Figure A.1

Bessel functions of the first kind Jn (x) . . . . . . . . . . . . . . . . 181

Figure A.2

Bessel functions of the second kind Yn (x) . . . . . . . . . . . . . . . 182

Figure A.3

Circular membrane with radius R . . . . . . . . . . . . . . . . . . . 183

Figure A.4

Fundamental modes of circular membrane with (a) First four redial
nodes and no angular nodes (b) First four redial nodes and one angular nodes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 186

Figure A.5

Fundamental modes of circular membrane with (a) First four redial
nodes and two angular nodes (b) First four redial nodes and three
angular nodes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 187

xv

Chapter 1
Introduction

1.1

Motivation and Objectives

The suppression of undesired vibration in industrial and engineering applications has been a
fertile field of research for decades. Vibration occurs in structures when they are subjected
to disturbances, regardless of the size of the structure. Some structures are huge such as
bridges and skyscrapers whereas others structures are tiny such as MEMS devices and nano
scale systems. To study vibration in mechanical systems, the systems are often modeled as
discrete systems with finite degrees of freedom. For some other applications, the mechanical
vibration occurs in systems that have infinite degrees of freedom, such systems are called
continuous systems. Examples of continuous systems are strings, beams, plates and membranes which are widely used in many engineering applications such as space structures. In
space applications, the weight and size play a crucial role in the cost of transportation and
deployment. For this reason, aerospace companies prefer to use light weight flexible continuous systems as a substitute to heavier and rigid elements. For instance, the regular space
1

telescope rigid mirrors are now substituted by thin deployable membrane mirrors that are
much lighter in weight. The same is true for satellite antennas where large and lightweight inflatable antennas, as shown in Fig. 1.1, are replacing regular antennas. The new lightweight
antennas are functionally efficient and can be folded into significantly smaller volume [4].
Although lightweight deployable structures have low cost transportation and deployment,
they are prone to disturbances.

In general, vibration in continuous systems can be suppressed passively by changing the
mass or stiffness in the systems, or actively, by sensing the states of the system and using
feedback control. Passive control methods are not very effective and are not well suited
for many applications. Active control methods are effective but they require many sensing
and actuators complicated sensing and actuating mechanisms and complex mathematical
models to compute the control effort. In this work, we present semi-control strategies to
suppress the vibration in continuous systems namely, strings and circular membranes using
boundary constraints. First, we investigate the energetics of string wrapping and unwrapping
around fixed obstacle located at one of the boundaries where the obstacle represents a passive
mechanism. Second, we present semi-active control strategy to suppress the vibration in
string using a scabbard-like actuator that interacts with string near one of its boundaries. We
present experimental verification to the control strategy based on the scabbard-like actuator
control strategy. Third, we extend our semi-active control strategy to vibrating circular
membranes. We use a small circular areal constraint that is applied sequentially on the
membrane to suppress the vibration. We present the complete mathematical solution of
the dynamic and energetics of circular membrane subjected to eccentric areal constraint and
different control strategies are presented to the vibration suppression of circular membrane.
2

S77E5022 1996 :05 :20 08 :04 :38

Figure 1.1: The spartan 207 inflatable antenna, NASA [4] “For interpretation of the references to color in this and all other figures, the reader is referred to the electronic version of
this dissertation“

1.2
1.2.1

Vibration Control of Strings
Passive Control of Vibrating Strings

The dynamics of vibrating strings has been a subject of study for a very long time but the
motion of strings vibrating against obstacles appeared in the technical literature relatively
recently. Early work on this problem can be credited to Citrini [5] who considered pointshaped obstacles. The element of string that comes in contact with point-shaped obstacles
can be assumed to be massless and hence the energy of the string, in the absence of damping,
was assumed to remain conserved. Amerio [6] investigated the motion of a string vibrating
against a rigid wall, parallel to the position of the string at rest. The motion of the string in
3

the presence of the unilateral constraint was posed as a problem in impact. The nature of
the impact was assumed to be elastic and the problem was formulated based on conservation
of energy of the string. A number of other researchers have also based their work on the
premise of energy conservation of the string. These include Schatzman [7], who investigated
the existence and uniqueness of solutions for concave obstacles and Haraux and Cabannes
[8], who established almost-periodic nature of solutions for straight and fixed obstacles.

In 1982, Burridge et al. [9] investigated the vibration of the sitar, an Indian stringed instrument. The sitar differs from the Western stringed instruments in that the bridge across
which the strings pass form a broad support, rather than a well-defined edge. During vibration, the sitar string wraps and unwraps around the gentle slope of the bridge and the length
of the vibrating part of the string varies during oscillation [9]. Burridge et al. [9] modeled
the impact of the string with the bridge as perfectly inelastic, discarding the assumption of
energy conservation of the string. Subsequently, Bamberger and Schatzman [10] proved the
existence of solutions which do not conserve energy with arbitrary obstacles and Ahn [11]
claimed energy loss of the string vibrating against flat obstacles. In conformity with earlier
work by Citrini [5], Ahn [11] also showed that energy remains conserved for highly-peaked
obstacles. Other work on string vibration against obstacles includes discretization [12] and
finite difference methods [13] for numerical simulation, and study of nonlinear effects of
varying amplitude and gravity [14] on extensible and non-extensible cables. In our work we
investigate the vibration of a string against an obstacle located at its boundary. Similar
to the work by Burridge, et al. [9], we assume the string to wrap and unwrap around the
obstacle during each oscillation. The impact of the string during wrapping is assumed to
be perfectly inelastic and the obstacle is implicitly assumed to be convex. The assumption
4

of convexity of the obstacle is both convenient and practical. Assuming that the string vibrates in a single mode at all times, it is shown that energy loss is higher for higher modes
of oscillation.

1.2.2

Active Control of String vibration

Several methods have been proposed in the literature for actively controlling string vibration. Some of these methods involve direct physical interaction with the string, using point
or distributed forces, for example. Other methods are less intrusive and are based on tension
variation or boundary control. In our work we propose to use a scabbard-like actuator that
imposes a zero-displacement constraint over some length of the string from one boundary.
Although this particular approach has not been proposed earlier, several papers in the literature have explored fundamentally similar ideas. Early work by Dutt and Ramakrishna
[15] investigated the application of a distributed control force to minimize the energy of a
part of the string while maintaining the energy of the uncontrolled part at a desired level.
The work was later extended [16] to investigate optimal positions and magnitudes of the
control forces. In our approach, the scabbard is intended to reduce the energy of the overall
system while dissipating the energy of a short length of the string. When the scabbard
is removed, the energy of the uncontrolled portion of the string is redistributed over the
entire length and thus cyclic application and removal of the scabbard can lead to gradual
dissipation of the total energy of the string. This idea is derived from earlier work on finite
degree-of-freedom systems [17, 18] and has similarities with the work by Hebrard and Henrot
[19], where exponential decay in the total energy is achieved by subjecting a portion of the
string to constant damping. The zero-displacement constraint enforced by the scabbard can
5

be viewed as impulsive actuation and in this regard our approach is similar to the approach
proposed by Myshkis [20].
The scabbard-like actuator proposed in this study essentially moves one boundary of
the string in the axial direction. Although vibration control using a scabbard has not been
proposed earlier, the effect of moving the boundary and varying the length of the string has
been investigated [21]. It has been shown that reducing the length of the string can increase
its energy through compression of the propagating waves. Many researchers have also investigated transverse boundary control of vibrating strings. Shahruz and Kurmaji [22] designed
a controller for a nonlinear model of axially moving strings using vertical control forces at
one support. Fung, et al. [23] proposed an exponentially stabilizing mass-damper-spring controller and Li, et al. [24] proposed boundary velocity feedback to stabilize a nonlinear model
obtained using Hamilton’s principle. Do and Pan [25] studied boundary control of flexible
marine risers; the controller was designed using Lyapunov’s direct method and backstepping
[26], and implemented using a hydraulic actuator. Zhang, et al. [27] designed boundary
controllers for a general class of nonlinear string-actuator systems; a nonlinear distributed
parameter model was used to account for large amplitude displacement and the associated
variation in tension. Other boundary control methods include, for example, wave control
[28, 29], wave cancellation [30, 31], and active sinks [32]. While wave control is a general
methodology, wave cancellation refers to control using a boundary force that prevents reflection of the wave and active sinks is a strategy to make vibration modes inactive by applying
a control force at some distance from the source of excitation. Our method is fundamentally
different in that a zero displacement constraint is applied at the boundary.
Although there are many theoretical studies on the active control of vibrating strings, exper6

imental results on control of string vibration have been scarce. This can be partly attributed
to the fact that traditional strings have high internal damping and do not lend themselves
well to vibration control experiments. Another challenge in performing vibration control experiments with strings relates to actuator placement. Unlike a beam, to which piezoelectric
actuators can be easily mounted, it is inconvenient to mount actuators on a string. Many papers in the literature present experimental results, but they primarily consider axially moving
strings; few papers consider strings that are not translating between supports. Baicu, et al.
[33] developed an active boundary controller for the nonlinear model of an elastic cable; the
angle of the cable at one boundary was sensed and its position was controlled to asymptotically decay modal vibration. A wave absorption strategy was developed by Saigo, et al. [29];
a voice-coil actuator was placed near a fixed boundary of the string to suppress vibration.
For the experimental verification of scabbard like actuator method, we present a control
methodology in which energy reduction is achieved through cyclic application and removal
of a constraint. In our method, a zero-displacement constraint is enforced at one point on
the string, close to one of its boundaries. This effectively results in two vibrating strings,
one of which is much shorter in length than the other. The vibration of the shorter length
string decays rapidly due to high internal damping and when the constraint is released, the
remaining energy of the string is redistributed over its entire length; this allows the cycle of
constraint application and removal to be repeated for vibration suppression.

1.3

Dynamics and Vibration Control of Membrane

The dynamics of a membrane is a classical problem in mechanical vibrations. Among the different types of membranes, circular membranes are the most commonly studied due to their
7

large number of applications. From the study of musical notes of percussion instruments,
circular membranes have been used to design diaphragms for condenser microphones, model
the dynamics of the human ear [34], understand the vibration characteristics of membrane
mirrors and gossamer structures [35], measure surface tension [36], and design ink-jet printers [37]. For a majority of these applications, unconstrained membrane models have been
used. A constrained membrane is one in which a portion of the membrane does not undergo
vibration, and an application where constrained membrane models have been used is waveguides. Early work on this subject dates back to the 1960’s and 70’s. Yee and Audeh [38], [1]
used the point-matching method to determine eigenfrequencies of waveguides with eccentric
cross-section. Davies and Muilwyk [39] and Steele [40] used the finite difference method to
compute eigenvalues of waveguides with arbitrary cross-section. The Helmholtz equation,
which describes the dynamics of waveguides, was solved by Arlett, et al. [41] and Gass [42]
using finite element methods, by Laura [43] and Hine [44] using the Galerkin method, and
by Bulley and Davies [45] using the Rayleigh-Ritz method.

The similarity between the differential equations of membranes and waveguides motivated
the study of circular membranes with constraints in the 1970’s and 80’s. The eigenvalues of
a circular membrane with an internal eccentric circular areal constraint were first computed
by Nagaya [2] using Fourier series. Another solution to this problem was presented by Lin
[46], who used Graf’s addition theorem to transform the Bessel functions from one set of
polar coordinates to another. Lin [46] missed a term in the final expression, which affected
the results of calculation. The missing term was identified by Singh and Kothari [47] and
Singh, et al. [3], who used Nagaya’s approach together with Graf’s addition theorem. Later,
Nagaya and Hai [48] studied composite membranes with arbitrarily shaped inner and outer
8

boundaries. S. Noga [49] presented numerical simulation for the free transversal vibrations
of a systems of two annular and circular membranes connected by a Winkler elastic layer.
Arokiasamy, A. [50] presented work for vibrating rectangular membrane by using computer
graphics. E. Alvarado-Anell and et al. [51] presented simulation and animation of circular,
rectangular and elliptical membranes by using computer MAPLE software. The circular
membrane is also used for the modelling of human tympanum, where computer simulations
are presented [52, 53].
The vibration of membranes remains an active subject of research with recent work
focusing on complex boundary conditions and shapes, inhomogeneities, accuracy of solutions, and applications. Investigations of membranes with complex boundary conditions and
shapes include the work by Laura and Gutierrez [54], Yu and Wang [55], Kang and Lee [56],
and Chen, et al. [57]. Investigations of non-homogeneous membranes include the work by
Cortinez and Laura [58] and Cap [59], and that of accuracy of solutions include the work by
Zhao and Stevens [60]. Applications of membranes are diverse and membrane vibrations have
been investigated in the context of micro-air vehicle wings [61], [62], [63], energy harvesting
[64, 65]. Wang, O. [66] presented a study about the wave propagation in a piezoelectric
coupled cylindrical membrane shell. Furthermore, membrane structures are used in many
different applications such as concentrators, planar configurations, solar sails and optical applications [67]. Although many theoretical works in the literature investigated the vibration
of membranes in general, only handful of studies are conducted to the experimental study
of vibration in membranes, a thorough literature review about membrane experiments can
be found in [67].
Some studies in the literature are conduct for the vibration control of membranes using
9

different active control schemes. Grosso, Ronald P. and Yellin, Martin [68] investigated the
vibration control of a membrane mirror assembly by using electrostatic actuators. Joo, I.
[69] and [70] presented a study on the control of circular and rectangular membranes by
reachability of movement states in fixed time. James D. Moore [71] presented a design and
testing of a one-meter membrane mirror with active boundary control. Mauro V. Aguanno
and et al. [72] presented an experimental work on investigation to demonstrate a full-field
laser vibrometer system that could replace electro-mechanical scanning with electronic scanning within a programmable stand-alone and relatively low-cost digital camera. Recently, a
study on the active vibration control of thin and flexible disc is presented by [73], where the
thin disc is equipped with two piezoelectric circular patches: one of them works as a sensor
and the other is used as an actuator for damping of the most vibrating modes in a specified
bandwidth. Another recent work on the vibration suppression of membrane is presented
by P. A. Tarazaga and et al. [74], who investigated the using of vibro-acoustic technique to
suppress the vibration in a pressurized optical membrane mirror.

In our work we use a small areal circular boundary constraint that interacts actively with the
membrane during its motion and suppresses the vibration through sequential applications
and removals of the constraint. Application and removal of constraint requires a complete
and accurate solution for dynamics of the membrane with and without constraint. As we
have seen, a thorough review of the literature indicates that investigations of constrained
membranes have been limited to computing the eigenfrequencies of vibration. Although
accurate computation of the eigenfrequencies may be sufficient for many applications, it
is not sufficient for accurate computation of the mode shapes and/or modal coefficients of
a constrained membrane. Simulation of the dynamics requires computation of the modal
10

coefficients using orthogonality property of the modes. The orthogonality property of modes
is well-known but it has not been established mathematically for constrained membranes. We
establish the orthogonality property of all modes of a circular membrane with an internal
circular areal constraint for the first time in this work. The computation of the modal
coefficients require accurate computation of the mode shapes and this requires proper choice
of the number of angular modes. We present an algorithm that determines the appropriate
number of angular modes for arbitrary size and location of the constraint. In addition to
this algorithm, which computes the mode shapes, we provide an algorithm for computation
of the modal coefficients for arbitrary initial conditions. Also, we investigate the effect of
application and removal of the constraint on the energetics of the vibrating membrane and
provide semi-active control strategies to suppress the vibrating in the membrane through
cyclic application and removal of the constraint.
Our entire work in this dissertation is organized as follows; in Chapter 2 we present
a passive control strategy for the string vibration by subjecting the string to an obstacle
located at one boundary. A semi-active control strategy for the string vibration is presented
in Chapter 3, where a scabbard-like actuator is used near one boundary to suppress the
string vibration. In Chapter 4 ,we present experimental verification to the control strategy
based on the scabbard-like actuator control strategy. Chapter 5 presents accurate dynamics
and simulation of constrained membrane under arbitrary initial conditions. In Chapter 6 we
use the analysis presented in Chapter 5 to investigate the energetics of the application and
removal of the areal boundary constraint and three different control strategies are presented.
Chapter 7 includes the concluding remarks and the future work.

11

Chapter 2
Vibration of a String Wrapping and
Unwrapping Around an Obstacle

2.1

Introduction

In this chapter we investigate the vibration of a string against an obstacle located at its
boundary. Similar to the work by Burridge, et al.[9], we assume the string to wrap and
unwrap around the obstacle during each oscillation. The impact of the string during wrapping is assumed to be perfectly inelastic and the obstacle is implicitly assumed to be convex.
The assumption of convexity of the obstacle is both convenient and practical. The obstacle
constrains the motion of the string and in this regard the mechanism for energy loss is a
continuous-system version of the energy dissipation methodology proposed for finite degreeof-freedom systems by Issa, et al.[18]. Since the energy of the string decreases even in the
absence of damping, the obstacle can be regarded as a passive mechanism for vibration suppression and control. This chapter is organized as follows. A formal problem statement
12

and a list of the assumptions made in our analysis is provided in Section 2.2. The analytical
model for computing the geometry of the string as it wraps and unwraps around the obstacle
during oscillation is presented in Section 2.3. In Section 2.4 we provide simulation results for
percentage energy loss and length of wrapping during each cycle of oscillation for different
modes with circular- and elliptic-shaped obstacles. Using numerical simulations, we show
in Section 2.5 that percentage energy loss can be increased significantly by changing the
orientation of the obstacle.

2.2

Problem Statement and Assumptions

Consider a string vibrating against an obstacle placed at one of its boundaries, as shown in
Fig. 2.1. We investigate energy dissipation in the string under the following assumptions:

A1. The obstacle is rigid and has the following geometry

y = f (x),

dy
=0
dx x=0

y(0) = 0,

(2.1)

A2. The string is homogeneous and has a constant mass per unit length denoted by ρ.
The tension in the string is equal to T and remains constant at all times. The string
undergoes transverse vibration in the xy plane and is not affected by gravity.

A3. The amplitude of oscillation of the string is small and therefore the equation of motion
of the string can be expressed by the standard relation [75]
1 ∂2y
∂2y
=
,
∂x2
c2 ∂t2
13

c

T /ρ

(2.2)

where y(x, t) is the displacement of the string at a distance x from the origin at time
t.

A4. The string wraps around the obstacle during vibration. Over each time step during
wrapping, a small element of the string comes to rest on the obstacle through perfectly
inelastic collisions. The wrapping process continues till the freely vibrating portion of
the string has no more kinetic energy.

A5. The surface of the obstacle is not sticky and the string unwraps from the obstacle
without any loss of energy.

A6. At the initial time t = 0, the string has no contact with the obstacle. It is in its mean
position with zero potential energy and kinetic energy equal to E0 .

A7. The string continues to vibrate in the mode in which it started its vibration at the
initial time. This implies that each point of the string, not in contact with the obstacle,
has the same frequency of vibration1 at any instant of time, and the number of nodes2
in the vibrating string remains constant.

A8. The string has no internal damping, i.e., the energy of the string will remain conserved
during free vibration.
1 The

frequency of the string is not constant; it varies with time and amplitude when the string
wraps and unwraps around the obstacle.
2 A node is a point with zero displacement. As the string wraps or unwraps around the obstacle,
the location of the node(s) change and therefore node(s) have nonzero horizontal velocity.

14

y

y = f (x)

x
y(x, t)
(a)
y

y = f (x)

y(x, t)

x
(b)

y

y = f (x)
y(x, t)
(¯, y )
x ¯
x
l
(c)

Figure 2.1: A string vibrating against an obstacle is shown in three configurations: (a) the
string has both potential and kinetic energy and is not in contact with the obstacle, (b) the
string has kinetic energy but no potential energy and is not in contact with the obstacle, (c)
the string has both potential and kinetic energy and has wrapped around the obstacle. In
a wrapped configuration, (¯, y ) denotes the coordinate where the string breaks contact with
x ¯
the obstacle
15

2.3
2.3.1

Analytical Model
Boundary conditions and general solution

A general solution to the partial differential equation in Eq. (2.2) can be written as [75]

y(x, t) = (α1 sin λx + α2 cos λx) (α3 sin ωt + α4 cos ωt)

(2.3)

where αi , i = 1, 2, 3, 4 are constants, ω is the circular frequency and λ is the wave length
and it is related to ω by the relation
ω

cλ

(2.4)

At time t = 0, the string is at the mean position, i.e. y(x, 0) ≡ 0, per assumption A6. This
implies α4 = 0. The solution in Eq. (2.3) can now be written as

y(x, t) = (A sin λx + B cos λx) sin ωt,

A = α1 α3 ,

B = α2 α3

(2.5)

At the right boundary, the string satisfies the relation y(l, t) = 0 for all t. Using Eq. (2.5)
we get
A sin λl + B cos λl = 0

=⇒

B = −A tan λl

(2.6)

Substitution of Eq. (2.6) into Eq. (2.5) gives the solution

y(x, t) = A (sin λx − tan λl cos λx) sin ωt

(2.7)

We now consider the boundary conditions at the contact break point. From Fig. 2.1 we
16

have
f (¯) = y(¯, t)
x
x

=⇒

f (¯) = A (sin λ¯ − tan λl cos λ¯) sin ωt
x
x
x

(2.8)

Also, the string is tangential to the obstacle at the contact break point x = x, i.e.,
¯

f ′ (¯) =
x

∂y(¯, t)
x
∂x

=⇒

f ′ (¯) = λA (cos λ¯ + tan λl sin λ¯) sin ωt
x
x
x

(2.9)

From Eqs. (2.8) and (2.9) we get

f (¯)
x
tan λ(l − x) = −λ ′
¯
f (¯)
x

(2.10)

which indicates that λ can be computed from the value of x. The solution of Eq. (2.10) is
¯
however not unique - each non-trivial value of λ corresponds to a mode of vibration of the
string. Since λ is an implicit function of x, we can rewrite Eq. (2.9) as follows
¯

A sin ωt = g(¯),
x

f ′ (x)
λ (cos λx + tan λl sin λx)

g(x)

(2.11)

Equation (2.11) can be used to compute t from the value of x. The existence of the solution,
¯
however, depends on the magnitude of A. We now discuss the procedure for computing A.
Let the total energy of the string at any time t be denoted by E. Then,

obs
vib
E = Epe + Eke = Epe + Epe + Eke

(2.12)

obs
vib
where Epe is the potential energy of the string wrapped around the obstacle, Epe is the
potential energy of the freely vibrating string, Epe is the total potential energy, and Eke is
17

the kinetic energy of the string. The total potential energy of the string is computed as the
product of the tension T (which is assumed constant) and elongation of the string [76]. The
elongation of the string is computed by integrating the strain of the string along the length
wrapped around the obstacle and along the length of string vibrating freely. Thus, the total
potential energy can be written as

Epe = T

dl = T

(

dx2 + dy 2 − dx) = T

(

l

1 + (dy/dx)2 − 1)dx

(

(ds − dx) = T

1 + (dy/dx)2 − 1)dx

0

(

1 + (dy/dx)2 − 1)dx

obs
vib
and Epe and Epe can be written as

obs
Epe = T

x
¯
0
l

vib
Epe = T

x
¯

Since the string conforms to the shape of the obstacle, (dy/dx) = f ′ (x) for x ∈ [0, x]. By
¯
expressing (dy/dx) = y ′ (x, t) for x ∈ [¯, l] and simplifying using Eqs. (2.7) and (2.11), we
x
get

obs
Epe = T
vib
Epe = T
≈

x
¯
0
l
x
¯

1 + f ′ (x) 2 − 1 dx
1 + y ′(x, t) 2 − 1 dx

2
T l ′
y (x, t) dx
2 x
¯

18

(2.13)

then
l
T 2 2 2
λ A sin ωt
[cos λx + tan λl sin λx]2 dx
2
x
¯
1
2 sin2 ωt sec2 λl 2λ(l − x) + sin[2λ(l − x)]
T λA
¯
¯
=
8
1
T λ sec2 λl 2λ(l − x) + sin[2λ(l − x)] g(¯) 2
¯
¯
x
=
8

vib
Epe =

(2.14)

An element of string of length dx has a mass of ρ dx and velocity is y(x, t). Thus, the kinetic
˙
energy of the freely vibrating string can be written and simplified as follows

Eke

1 l
ρ y(x, t) 2 dx
˙
2 x
¯
l
ρ 2 2
ω A cos2 ωt
=
[sin λx − tan λl cos λx]2 dx
2
x
¯
1
2 A2 cos2 ωt sec2 λl 2λ(l − x) − sin[2λ(l − x)]
ρω
¯
¯
=
8λ
1
ρ ω 2 sec2 λl 2λ(l − x) − sin[2λ(l − x)] A2 − g(¯) 2
¯
¯
x
=
8λ
=

(2.15)

From Eqs. (2.12), (2.13), (2.14) and (2.15) it is easy to verify that the energy expression has
the form
E = h(¯, A)
x

(2.16)

For a configuration in which the string is wrapped around the obstacle, the complete
solution can be determined from the values of x and E using the four-step algorithm below:
¯
1. Use Eq. (2.10) to determine the value of λ. Since Eq. (2.10) provides multiple nontrivial solutions that correspond to different modes of vibration, the solution corresponding to the initial mode of vibration should be chosen - see assumption A7.
19

y
(¯j+1 , yj+1 )
x
¯

y = f (x)

C
y(x, t)

(¯j , yj )
x ¯
xj
¯

Dx

x

B

A
Dx

xj+1
¯

Figure 2.2: Small string segment ∆x wraps around obstacle after perfectly inelastic collision.
In the magnified image, AB denotes the small string segment of length ∆x that wraps around
the obstacle over the region AC
2. Use Eq. (2.4) to compute ω.
3. Compute A from Eq. (2.16) using the values of x, E, λ and ω.
¯
4. Compute the time t from Eq. (2.11) by substituting in the values of x, A and λ.
¯

The complete solution can now be described using Eqs. (2.1) and (2.7) as follows:

y(x, t) =

2.3.2



 f (x)

: x ∈ [0, x]
¯

(2.17)


 A (sin λx − tan λl cos λx) sin ωt : x ∈ [¯, l ]
x

Wrapping of the string

From our discussion in the last section we know that the geometry of the string can be
determined from the values of x and E. In this section we discuss the method for computing
¯
these values at regular intervals of time. Let {¯i , Ei } denote the values of x and E at time
x
¯
t = ti , i = 0, 1, 2, · · · , k. We assume t0 = 0. Then, from assumption A6, x0 = 0 and the
¯
value of E0 is known. We will discuss the method for determining the value of tk which
denotes the time after which the string begins to unwrap.
Let us assume that for some i = j, {¯j , Ej } is known. We outline the method for
x
20

computing {¯j+1 , Ej+1 } from the values of {¯j , Ej }. Choose a small segment of the
x
x
vibrating string that is expected to wrap around the obstacle over a small interval of time.
Let the projection of this string segment AB on the x axis be ∆x as shown in Fig. 2.2. The
kinetic energy of this string segment, which will be lost due to inelastic collision, can be
computed from Eq. (2.7) as follows
¯
ρ xj +∆x
y(x, t) 2 dx
˙
Elost =
2 x
¯j
xj +∆x
¯
2
ρ 2 2
ωj Aj cos2 ωj tj
=
sin λj x − tan λj l cos λj x dx
2
xj
¯

(2.18)

where Aj , ωj , λj and tj denote values of A, ω, λ and t, respectively, derived for the pair
{¯j , Ej }. Using Eq. (2.18), Ej+1 can be computed as follows
x

Ej+1 = Ej − Elost ,

j = 1, 2, · · · , k − 1

(2.19)

To compute xj+1 , j = 1, 2, · · · , k − 1, we make the following general assumption:
¯

A9. With reference to Fig. 2.2, the potential energy of the vibrating string segment AB
at time tj is equal to the potential energy of the string segment AC wrapped on the
obstacle at time tj+1 .

Using Eqs. (2.13) and (2.14) assumption A9 can be mathematically expressed as follows
xj+1
¯
xj
¯
≈

1 + f ′ (x) 2 − 1 dx =

¯
2
1 xj +∆x ′
y (x, t) dx =⇒
2 x
¯j
21

xj +∆x
¯
xj
¯

1 + y ′ (x, t) 2 − 1 dx

¯
xj +∆x
¯
2
2
1
1 xj +∆x ′
y (x, t) dx = λ2 A2 sin2 ωj tj
cos λj x + tan λj l sin λj x dx
2 x
2 j j
¯j
xj
¯
(2.20)

Equation (2.20) can be used to determine xj+1 . The values of Aj+1 , ωj+1 , λj+1 and tj+1
¯
are computed from the values of xj+1 and Ej+1 . The iterative process is terminated when
¯
the kinetic energy of the vibrating string segment becomes approximately equal to zero. At
this time, which is denoted as tk , the string stops wrapping and begins to unwrap.

2.3.3

Unwrapping of the string

Similar to wrapping, the geometry of the string during unwrapping is computed from the
values of x and E. The string begins to unwrap at t = tk ; at this time the values of
¯
x = xk and E = Ek are known. Let {¯i , Ei } denote the values of x and E at time
¯
¯
x
¯
t = ti , i = k, k + 1, k + 2, · · · , l, where tl denotes the time when the string has unwrapped
completely. We outline the method for computing {¯j+1 , Ej+1 } from the values of {¯j , Ej }
x
x
for k ≤ j ≤ l − 1. One chooses a small segment of the string that is expected to unwrap over
a small interval of time. Then we let the projection of this string segment on the x axis be
∆x. Then,
xj+1 = xj − ∆x,
¯
¯

j = k, k + 1, · · · , l − 1

(2.21)

Since there is no loss of kinetic energy during unwrapping (see assumption A5), we have

Ej+1 = Ej ,

j = k, k + 1, · · · , l − 1

(2.22)

The values of Aj+1 , ωj+1 , λj+1 and tj+1 are computed iteratively from the values of
xj+1 and Ej+1 . The iterative process is terminated at t = tl when the potential energy of
¯
22

the string is equal to its value at the mean position.

2.4

Numerical Simulations

Consider a string with

T = 1 N,

ρ = 0.025 kg/m,

l=4m

(2.23)

The obstacle is assumed to be a circle of radius R and center coordinates (x, y) ≡ (0, R),
i.e.,
y = f (x) = R −

R2 − x2 ,

0≤x≤R

(2.24)

It can be verified that f (x) in Eq. (2.24) satisfies the boundary conditions in Eq. (2.1). For
R = 1 m and ∆x = 0.001 m, we compute the percentage loss of energy over one cycle of
string oscillation for three different values of initial energy E0 and for oscillation in the first,
second, and third modes, respectively. These values are shown in table 2.1 together with
the values of xk , which is a measure of the length of wrapping around the obstacle. For the
¯
special case of E0 = 0.5 J, we plot the percentage loss of energy for three consecutive cycles
of string vibration in the first two modes. These plots are shown in Fig. 2.3. Figure 2.6 plots

Table 2.1: Percentage energy loss over one cycle of oscillation and xk for different values of
¯
E0 and three modes of oscillation, all with R = 1 m
E0 = 1.00 J
E0 = 0.50 J
E0 = 0.25 J

Mode 1
0.491%, 0.729 m
0.256%, 0.596 m
0.113%, 0.461 m

Mode 2
1.598%, 0.696 m
0.861%, 0.569 m
0.402%, 0.445 m

23

Mode 3
2.752%, 0.654 m
1.547%, 0.537 m
0.766%, 0.424 m

the decay in energy as a function of time for vibration in the first four modes with E0 = 0.5
J. The following observations can be made from the plots in Figs. 2.3 and 2.6, and the data
in table 2.1
• For any mode of oscillation, it can be seen that the percentage energy loss is higher for
higher values of E0 . This is not surprising since higher values of E0 results in higher
kinetic energy and greater length of wrapping, as evident from the values of xk in table
¯
2.1, and consequently more energy loss through inelastic collision. The same argument
can explain the reduction in the percentage loss of energy over consecutive cycles of
vibration in Fig. 2.3.
• The percentage energy loss is higher for higher modes of oscillation for the same value
of E0 . This is true for the same number of cycles (see Fig. 2.3) as well as for the same
length of time (see Fig. 2.6) and is due to the fact that the velocities of the string
associated with higher frequencies are higher in higher modes, and as a consequence
the loss upon impact is higher. The value of xk is less for the higher modes but this
¯
does not have a significant effect on the percentage of energy loss.

24

100.0
0.256%

98.8

Percentage energy left

Percentage energy left

100.0

cycle-3
99.6

0.255%

99.4
cycle-1
0

0.254%

cycle-2
4

8

12

0.861%
cycle-3

99.0
0.853%
98.0
cycle-1
0

0.846%

cycle-2
4

2

Time (s)

Time (s)

(a)

6

(b)

Figure 2.3: Plot of percentage energy content of the string over three consecutive cycles of vibration in (a) Mode 1, and (b)
Mode 2. For both cases, the initial energy of the string was E0 = 0.5 J

25

100

Mode 1

Percentage energy left

80

60

Mode 2

40

Mode 3

20

Mode 4

0

4

8

12

16

20

Time (min)
Figure 2.4: Exponential decay in the energy of a string wrapping and unwrapping around
an obstacle. The plots show energy decay for single-mode vibration in the first four modes
with E0 = 0.5 J
The geometry of the string at different points in time during one cycle of oscillation is shown
in Fig. 2.5 for Mode 1 and Mode 2 with initial energy E0 = 0.5 J. It can be seen from these
plots that a fixed point on the string moves in the y direction only when the string is not
in contact with the obstacle but moves in both the x and y directions during wrapping and
unwrapping. From the plot for Mode 2, it is also clear that a node is not a fixed point on
the string. It is a point of zero displacement but has nonzero velocity during wrapping and
unwrapping.
26

2

2

2

1

1

1

0

0

0

−1

−1

−1

−2

0

2

4

−2

0

(a)

2
(b)

4

−2

0

2
(c)

Figure 2.5: A string vibrating against a circular obstacle in (a) Mode 1 (b) Mode 2, and (c) Mode 3

27

4

To study the effect of the shape of the obstacle on percentage energy loss, we fix the value
of the initial energy to E0 = 0.5 J and study the following four cases where the obstacle is:
(a) a circle with R = 0.5 m

(b) a circle with R = 1.0 m

(c) a circle with R = 1.5 m

(d) an ellipse with semi-major and semi-minor axes lengths of 1.2 m and 1.0 m, respectively, and with the major axis aligned with the x axis

and satisfy the boundary conditions in Eq. (2.1). The results are shown in Table 2.2. It is
clear from the results that for circular obstacles the percentage energy loss increases with
increase in radius and vice versa. This is in agreement with the results expected for the
limiting cases, namely, percentage energy loss is zero when the radius of the circle is zero
and is equal to 100% when the radius is infinity. The ellipse in case (d) circumscribes the
circle in case (b) and provides a lower slope for the wrapping curve. A comparison of the
data for cases (b) and (d) indicates that a slight decrease in slope of the obstacle results in
significantly higher percentage of energy loss.

Table 2.2: Percentage energy loss over one cycle of oscillation and xk for obstacles of different
¯
shapes and sizes, all with E0 = 0.50 J
Case
(a)
(b)
(c)
(d)

Mode 1
0.029%, 0.295
0.256%, 0.596
0.915%, 0.897
0.665%, 0.803

m
m
m
m

Mode 2
0.113%, 0.291
0.861%, 0.569
2.549%, 0.815
2.021%, 0.750

28

m
m
m
m

Mode 3
0.237%, 0.286
1.547%, 0.537
3.887%, 0.737
3.278%, 0.692

m
m
m
m

2.5

Effect of Change in Slope of Obstacle

We consider the obstacle in Fig. 2.6 where the curve y = g(x) is obtained by rotating the
curve y = f (x) in Fig. 2.1 clockwise by angle θ about point O. To deal with this problem,
we modify assumptions A1 and A6 as follows:

A1. The obstacle is rigid and has the following geometry

y = f (x),

y(0) = 0,

dy
= − tan θ
dx x=0

(2.25)

y
string wrapped around obstacle
y(x, t)
y = g(x)

(¯, y )
x ¯

O
θ

x

initial configuration of the string

Figure 2.6: A string vibrating against an obstacle. The obstacle is identical to the one in
Fig. 2.1 but rotated clockwise by angle θ about point O

A6. At the initial time t = 0, the string has no contact with the obstacle. It has zero kinetic
energy and potential energy equal to E0 . The displacement of the string at the initial
time corresponds to a single mode of free vibration as shown in Fig. 2.6.
29

The remaining assumptions, A2 through A5 and A6 through A9, are not changed. In table
2.3 we present simulation results for a string with

T = 1 N,

ρ = 0.025 kg/m,

l = 4 m,

E0 = 0.50 J

(2.26)

and a circular obstacle of radius R = 1 m. A comparison of the results indicates that
percentage energy loss is significantly higher for higher values of θ.
Table 2.3: Percentage energy loss over one cycle of oscillation and xk for two modes of
¯
oscillation with different values of θ
θ

Mode 1

Mode 2

0◦

0.256%, 0.596 m

0.861%, 0.569 m

15 ◦

1.486%, 0.866 m

3.478%, 0.841 m

30 ◦

9.233%, 1.095 m

13.72%, 1.080 m

In this chapter, the string was assumed to have no damping. In reality, the string will
have damping and this will enhance the rate of energy decay. The rate of energy decay will
however not be constant even if the damping ratio of the string is constant. As the string
wraps around the obstacle, its effective length decreases and frequency of vibration increases
- this will increase the rate of energy decay which depends on the product of damping ratio
and natural frequency.

30

Chapter 3
Vibration Control of a String Using A
Scabbard-Like Actuator

3.1

Introduction

In this chapter, we investigate the dynamics and energetics of vibrating string with fixed ends
subjected to moving constraint at one boundary. The constraint is applied by a scabbard
that moves a small distance along the mean position of the string. The scabbard is moved
instantaneously such that the position and velocity of the string outside the scabbard is
unaffected immediately after application of the constraint, whereas the length of the string
covered by the scabbard is brought to rest. The constraint is removed by moving the scabbard
back to its original position and the change in energy of the string is investigated for different
values of scabbard travel distance and time of application of the constraint. This chapter is
organized as follows. A formal problem statement and a list of assumptions are provided in
Section 3.2. Assuming that the string is vibrating in its fundamental mode, we use analytical
31

methods in Section 3.3 to study the dynamics of the string and change in its energy after
application of the constraint. The dynamics of the string and change in its energy after
removal of the constraint is studied in Section 3.4. In Section 3.5 we present numerical
simulation results for one cycle of constraint application and removal - these results indicate
that the energy of the string can increase or decrease depending on the time of application
of the constraint. In Section 3.6 we repeat the analyses of Sections 3.3 and 3.4 for arbitrary
initial conditions of the string and develop a control strategy to reduce the energy of the
string for every cycle of constraint application and removal. Numerical simulation results
are then presented to demonstrate the efficacy of the control strategy.

3.2

Problem Statement and Assumptions

Consider the vibrating string in Fig. 3.1 (a), that passes through a scabbard located at its
left boundary. At time t = tc , the scabbard is moved instantaneously to the right by distance
x0 along the mean position of the string. This is shown in Fig. 3.1 (b). At some future time
t = tr , tr > tc , the scabbard is moved back to its original position. To investigate the
effect of application and removal of the scabbard on the vibration of the string, we make the
following simplifying assumptions:

A1. The string is homogeneous and has a constant mass per unit length denoted by ρ.
The tension in the string is equal to T and remains constant at all times. The string
undergoes transverse vibration in the xy plane and is not affected by gravity.
A2. The string is initially vibrating in its fundamental mode. This assumption will be
removed and a general displacement profile of the string will be assumed in Section
32

3.6.

A3. The amplitude of oscillation of the string is small and therefore the equation of motion
of the string can be expressed by the standard relation [75]
∂2y
∂x2

=

1
c2

∂2y
∂t2

,

c

T /ρ

(3.1)

where y(x, t) is the displacement of the string at a distance x and time t.

A4. The string has no internal damping, i.e., the energy of the string will remain conserved
for free vibration.

t ≃ tc

y

x
l

scabbard

(a)
approximate deformation
t ≃ tc
y

expected real deformation
scabbard

x0

x

l
(b)

Figure 3.1: (a) A vibrating string (b) A zero-displacement constraint is applied to the string
at time t = tc over the length segment x ∈ [0, x0 ) using a scabbard
33

A5. At time t = tc , the scabbard is moved instantaneously to the right by a distance
x0 along the mean position of the string. The movement of the scabbard imposes a
zero-displacement constraint over the length interval x ∈ [0, x0 ). The displacement
and velocity of the string over the remaining interval x ∈ [x0 , l] remains unchanged
immediately after movement of the scabbard.
A6. The distance x0 is small compared to the length of the string. Consequently, the
discontinuity1 in the displacement of the string at x0 at time tc , shown in Fig. 3.1 (b),
will be small.
A7. At time t = tr , tr > tc , the scabbard is instantaneously moved back to its original
position, i.e., to the left by a distance x0 along the mean position of the string. The
displacement and velocity of the string over the interval x ∈ [x0 , l] remains unchanged
immediately after movement of the scabbard.

3.3

Effect of Applying Constraint on the Dynamics of
the String

3.3.1

Equation of motion after application of the constraint

The motion of the string prior to application of the constraint is shown in Fig. 3.1 (a) and
is described by the equation

y0 (x, t) = A0 sin λ0 x sin ω0 t,
1 In

λ0 = (π/l),

ω0 = cλ0

(3.2)

reality, there will be no discontinuity in the string displacement. The same is true for the
dynamic model since a finite number of Fourier coefficients will be used.

34

where A0 is the amplitude, ω0 is the circular frequency, and c is defined in Eq. (3.1).
To obtain the equation of motion of the string after the constraint has been applied, we
solve Eq. (3.1) subject to initial and boundary conditions consistent with the constraint. A
general solution to the partial differential equation in Eq. (3.1) is first written as [75]

y(x, t) = f (x) g(t)

(3.3)

where
f (x) = A sin λx + B cos λx

,

ω = cλ

(3.4)

g(t) = C sin ωt + D cos ωt
where A, B, C and D are constants, and ω is the circular frequency. The boundary conditions
for the constrained string, shown in Fig. 3.1 (b), give us the following relations

y(x0 , t) = 0
⇒ (A sin λx0 + B cos λx0 ) (C sin ωt + D cos ωt) = 0
⇒ B = −A tan λx0

(3.5)

y(l, t) = 0
⇒ (A sin λl + B cos λl) (C sin ωt + D cos ωt) = 0
⇒ B = −A tan λl

(3.6)

From Eqs. (3.5) and (3.6) we get

tan λx0 = tan λl

⇒

sin λ(l − x0 ) = 0
35

⇒

λn =

nπ
, n = 1, 2, · · ·
l − x0

(3.7)

By substituting the expressions for B and λn from Eqs. (3.6) and (3.7) in the expression for
f (x) in Eq. (3.4), we get

f (x) = A (sin λn x − tan λn l cos λn x) = −

A
sin λn (l − x)
cos λn l

(3.8)

The general solution to Eq. (3.1) is obtained by substituting Eqs. (3.4) and (3.8) into Eq. (3.3):

y(x, t) =











0

: x ∈ [0, x0 )

∞

(3.9)
sin λn (l − x) Cn sin cλn t + Dn cos cλn t

: x ∈ [x0 , l]

n=1

The constants Cn and Dn , n = 1, 2, · · · , are determined from initial conditions. Without
loss of generality, the time is first reset from t = tc to t = 0. Using assumption A5 in relation
to the displacement and velocity of the string over the interval x ∈ [x0 , l], we get

y(x, 0) = y0 (x, tc ) = A0 sin λ0 x sin ω0 tc , x ∈ [x0 , l]
∞
⇒
Dn sin λn (l − x) = A0 sin λ0 x sin ω0 tc , x ∈ [x0 , l]
n=1
l
2A0
sin ω0 tc
sin λ0 x sin λn (l − x) dx
(l − x0 )
x0
2A0 λn
sin λ0 x0 sin ω0 tc
= (−1)n
(l − x0 )(λ2 − λ2 )
n
0

⇒ Dn =

and,
36

(3.10)

y(x, 0) = y0 (x, tc ) = A0 ω0 sin λ0 x cos ω0 tc ,
˙
˙
x ∈ [x0 , l]
∞
⇒
cλn Cn sin λn (l − x) = A0 ω0 sin λ0 x cos ω0 tc , x ∈ [x0 , l]
n=1
l
2A0 ω0
cos ω0 tc
sin λ0 x sin λn (l − x) dx
ncπ
x0
2A0 ω0 λn
= (−1)n
sin λ0 x0 cos ω0 tc
ncπ(λ2 − λ2 )
n
0

⇒ Cn =

(3.11)

It is clear from Eqs. (3.10) and (3.11) that the values of the coefficients Cn and Dn , n =
1, 2, · · · , depend on the values of x0 and tc . The motion of the string after application of the
constraint is described by Eq. (3.9), and with the knowledge of these coefficients the motion
is completely defined.

3.3.2

Energetics of constraint application

In this section we compute the change in the energy of the string due to application of the
constraint. The energy of the string prior to application of the constraint is computed as
follows [75]
ρ l ∂y0 (x, t) 2
T l ∂y0 (x, t) 2
dx +
dx
E0 =
2 0
∂x
2 0
∂t

(3.12)

By differentiating the expression for y0 (x, t) in Eq. (3.2) with respect to x and t and substituting the relations in Eq. (3.12), we get

T l 2 2
ρl 2 2
E0 =
A0 λ0 sin2 ω0 t +
A ω cos2 ω0 t
22
22 0 0

(3.13)

Substituting the expressions for λ0 and ω0 from Eq. (3.2), λn from Eq. (3.7) and c from
37

Eq. (3.1) into Eq. (3.13), we get
E0 =

T π 2 A2
0

(3.14)

4l

After the constraint has been applied, the energy of the string is given by the relation

Ec =

ρ l ∂y(x, t) 2
T l ∂y(x, t) 2
dx +
dx
2 x
∂x
2 x
∂t
0
0

(3.15)

If we define the time functions pn , qn , n = 1, 2, · · · , as

pn (t) = Cn sin cλn t + Dn cos cλn t
qn (t) = Cn cos cλn t − Dn sin cλn t

(3.16)

Eq. (3.15) can be rewritten as follows

l



∞

2

T

pn λn cos λn (l − x) dx +
2 x
0 n=1


∞
l
T

p2 λ2 cos2 λn (l − x) dx +
=
n n
2 x
0 n=1


Ec =

+

+

l



∞

2

ρ

qn ωn sin λn (l − x) dx
2 x
0 n=1


∞
l
ρ
2 2

qn ωn sin2 λn (l − x) dx
2 x
0 n=1


∞ ∞

T l 


pm pn λm λn cos λm (l − x) cos λn (l − x) dx

2 x 

0 m=1 n=1
n=m


∞ ∞

ρ l 


qm qn ωm ωn sin λm (l − x) sin λn (l − x) dx

2 x 

0 m=1 n=1
n=m

(3.17)

2
2
2
Using the relation (p2 + qn ) = (Cn + Dn ), which can be easily shown from Eq. (3.16), the
n

38

expression for λn in Eq. (3.7), and the identities
l
2
x0 cos λn (l − x)dx = (l − x0 )/2,
l
2
x0 sin λn (l − x)dx = (l − x0 )/2,

l
x0 cos λn (l − x) cos λm (l − x)dx = 0
l
x0 sin λn (l − x) sin λm (l − x)dx = 0

Eq. (3.17) is simplified to the form

T
Ec = (l − x0 )
4

∞
2
2
2
2 (p2 + q 2 ) = T π
n2 (Cn + Dn )
λn n
n
4(l − x0 )
n=1
n=1
∞

(3.18)

The work done by the constraint is equal to the change in energy of the string, and is given
by the relation

Wc = (Ec − E0 ) =

T π2
4(l − x0 )

∞
n=1

T π 2 A2
2
2
0
n2 (Cn + Dn ) −
4l

(3.19)

From the expression for Wc in Eq. (3.19) and the expressions for Cn and Dn , n = 1, 2, · · · ,
in Eqs. 3.10) and (3.11), it is clear that the work done by the constraint depends on the
values of x0 and tc . Through numerical simulations, it will be shown in Section 3.5 that the
work done can be positive or negative, depending on the time of application of the constraint
tc .

3.4

Effect of Removing Constraint on the Dynamics of
the String

Assume the constraint is removed during the string vibration at an arbitrary time t = tr .
In this section we investigate the effect of the constraint removal on the dynamics of the
39

string and change in its energy. When the scabbard is removed instantaneously back, the
string configuration which represents the initial displacement and velocity, will have two
regions as shown in Fig. 3.2. The first region was covered by the scabbard and has zero
velocity and displacement, and the second region is the configuration of the string at the
instance of scabbard removal where the displacement and velocity are given by y(tr ) and
y(tr ) respectively.
˙

t = tr
y

x
scabbard

l

Figure 3.2: At time t = tr , the scabbard is instantaneously moved back to its original
position

3.4.1

Equation of motion after removal of the constraint

To obtain the equation of motion of the string after the constraint has been removed, we
start from the general solution to the partial differential equation in Eq. (3.1), namely

yr (x, t) = (α sin Λx + β cos Λx) (γ sin Ωt + δ cos Ωt),

Ω = cΛ

(3.20)

where x ∈ [0, l], α, β, γ and δ are constants, and Ω is the circular frequency. The boundary
40

conditions for the string, shown in Fig. 3.2, give us the following relations:

yr (0, t) = 0 ⇒ β = 0
yr (l, t) = 0

⇒ sin Λl = 0

(3.21)
⇒ Λk = (kπ/l), k = 1, 2, · · ·

Substituting the relations in Eq. (3.21) into Eq. (3.20), we get
∞
yr (x, t) =

sin Λk x γk sin cΛk t + δk cos cΛk t ,

x ∈ [0, l]

(3.22)

k=1
The constants γk and δk , k = 1, 2, · · · , are determined from initial conditions. Without loss
of generality, the time is first reset from t = tr to t = 0. Using assumption A7 in relation to
the displacement and velocity of the string over the interval x ∈ [0, l], we get

yr (x, 0) = y(x, tr )




0
∞

∞
⇒
δk sin Λk x =


sin λn (l − x) Cn sin cλn tr + Dn cos cλn tr

k=1

n=1

∞
2 l
⇒ δk =
sin Λk x sin λn (l − x) Cn sin cλn tr + Dn cos cλn tr dx
l x
0 n=1
∞
λn
2
(−1)n
= sin Λk x0
C sin cλn tr + Dn cos cλn tr
2 − λ2 n
l
Λ
n
n=1
k

and,

yr (x, 0) = y(x, tr ) ⇒
˙
˙

41

: x ∈ [0, x0 )
: x ∈ [x0 , l]

(3.23)

similarly,
∞

∞
cΛk γk sin Λk x =

k=1

cλn sin λn (l − x) Cn cos cλn tr − Dn sin cλn tr

: x ∈ [x0 , l]

n=1

l ∞
2
⇒ γk =
λn sin Λk x sin λn (l − x) Cn cos cλn tr − Dn sin cλn tr dx
lΛk x
0 n=1
∞
2
λ2
n
=
C cos cλn tr − Dn sin cλn tr
(3.24)
sin Λk x0
(−1)n
2 − λ2 n
kπ
Λ
n
n=1
k
The equation of motion of the string after removal of the constraint is given by Eq. (3.22),
where δk and γk are defined by Eqs. (3.23) and (3.24), respectively.

3.4.2

Energetics of constraint removal

The change in energy of the string due to removal of the constraint can be computed in the
same way as we computed the change in energy due to application of the constraint. The
energy of the string after removal of the constraint is given by the relation
ρ l ∂yr (x, t) 2
T l ∂yr (x, t) 2
dx +
dx
Er =
2 0
∂x
2 0
∂t

(3.25)

Using Eqs. 3.22), (3.23) and (3.24) to compute the derivatives of yr (x, t) with respect to x
and t, and substituting them in Eq. (3.25) and simplifying, we get the final expression
T π2
Er =
4l

∞
k=1
42

2
2
k 2 (γk + δk )

(3.26)

The change in energy due to removal of the constraint is computed as
T π2
Wr = (Er − Ec ) =
4

∞
k=1

2
2
2
2
2 (γk + δk ) − (Ck + Dk )
k
l
(l − x0 )

(3.27)

In the next section we will use numerical simulations to show that the change in energy
due to constraint removal is zero, i.e., Wr = 0. This result is intuitive because removal of
the constraint simply results in redistribution of the energy from its constrained length of
(l − x0 ) to its original length l.

3.5

Simulation of One Cycle of Constraint Application
and Removal

We consider a string of length l = 4 m, mass per unit length ρ = 0.25 kg/m, and tension
T = 1 N, vibrating in its first mode with unit amplitude A0 = 0.1 m. At t = 0 the string is
assumed to pass through its mean position, and the equation of motion of the string is

π
π
y0 (x, t) = 0.1 sin( x) sin( t)
4
2

(3.28)

The time period of oscillation of the string is 4 sec and Fig. 3.3 (a) shows the shape of
the string at different instants of time over the interval t ∈ [1.0, 3.0] sec. The string is in
its maximum potential energy configuration at t = 1.0 sec and t = 3.0 sec, and maximum
kinetic energy configuration at t = 2.0 sec.
We first present simulation results for percentage change in energy due to constraint
application for x0 ∈ {0.0l, 0.01l, 0.02l, · · · , 0.09l, 0.10l} and tc ∈ [1.7, 2.3] sec. The results,
43

shown in Fig. 3.3 (b), indicate that the energy of the string increases if the constraint is
applied when the string is far away from the mean position, and decreases if the constraint
is applied when the string is near its mean position, irrespective of the value of x0 . This
can be explained as follows. Upon application of the constraint, the potential energy of the
string increases. This increase is large when the string is far away from its mean position,
and although the constraint removes the kinetic energy of a portion of the string, there is
a net gain in energy. When the string is near its mean position, the increase in energy
due to change in potential energy is small compared to the loss of kinetic energy and as a
result the net change in energy is negative. When the string is at its mean position, there
is no change in potential energy upon application of the constraint and consequently the
energy loss is maximum in this configuration. For the same value of tc , a larger value of
x0 results in a larger increase in potential energy and a larger decrease in kinetic energy.
Consequently, a larger value of x0 results in a larger percentage change in the energy. This
can be verified from Fig. 3.3 (b). It should be noted that increase in energy due to movement
of the boundary that reduces the length of the string was also observed by Zhu and Zheng
[21].
Table 3.1: Comparison of values of Ec and Er for x0 = 0.1l and tr = 0.8 sec
No. of Fourier coefficients
50
100
1000

Ec (J)
0.6124
0.6127
0.6128

Er (J)
0.6096
0.6110
0.6123

error %
0.4723
0.2694
0.0853

We next investigate the change in energy of the string due to removal of the constraint.
Instead of computing the value of Wr in Eq. (3.27), we compute the values of Ec and Er ,
which equal the energy of the string before and after removal of the constraint. Using
44

Eqs. (3.18) and (3.26), Ec and Er are computed for varying number of Fourier coefficients.
From this data (see Table 3.1) it is clear that Ec and Er approach each other as the number of
Fourier coefficients increase. Although the data corresponds to the specific case of x0 = 0.1l
and tr = 0.8 sec, the same trends can be observed for all values of x0 and tr .
Although constraint removal does not change the energy of the system, it resets the
system for a new cycle of constraint application. In the next section we investigate the effect
of application of the constraint on the dynamics of the string with arbitrary initial conditions
such that sequential application and removal of constraints can be explored as a strategy for
vibration suppression.

45

2.0
x0 = 0.1l

1.0

t = 1.0
(Wc /E0 ) × 100

y(x, t)

0.1

0.0

0.0
x0 = 0.04l

−1.0

x0 = 0.03l

x0 = 0.02l

x0 = 0.01l
t = 3.0

−0.1
0

t = 2.0

1

2
x

3

4

−2.0
1.7

1.8

1.9

2.0

2.1

2.2

2.3

tc

(a)

(b)

Figure 3.3: (a) Position of the string at different time instants (b) Percentage change in energy of the string due to application
of the constraint for different values of x0 and tc

46

3.6

Sequential Application and Removal of Constraints

3.6.1

Dynamics of string with arbitrary initial conditions

In this section we relax assumption A2 and allow the string to have arbitrary initial conditions
and then repeat the analysis of sections 3.4.2 and 3.5. Prior to application of the constraint,
the equation of motion of the string has the form
∞
y0 (x, t) =

¯
sin Λk x γk sin cΛk t + δk cos cΛk t ,
¯

x ∈ [0, l]

(3.29)

k=1
which is similar to Eq. (3.22) and satisfies the boundary conditions y0 (0, t) = y0 (l, t) = 0.
In Eq. (3.29), the values of the coefficients γk , δk , k = 1, 2, · · · , are known. The energy of
¯ ¯
the string at the initial time is given by the relation
T π2
E0 =
4l

∞
k=1

k 2 (¯k + δk )
γ 2 ¯2

(3.30)

which is similar to Eq. (3.26). After application of the constraint, the equation of motion of
the string has the form

y(x, t) =











0
∞

: x ∈ [0, x0 )

¯
¯
sin λn (l − x) Cn sin cλn t + Dn cos cλn t

(3.31)
: x ∈ [x0 , l]

n=1

which is identical to Eq. (3.9) except that the coefficients Cn and Dn , n = 1, 2, · · · , have
¯
¯
¯
¯
been replaced by Cn and Dn . The constants Cn and Dn , n = 1, 2, · · · , are determined
from initial conditions by repeating the procedure used in section 3.3.1. Without loss of
generality, the time is reset from t = tc to t = 0, and the displacement and velocity of the
47

string given by Eqs. (3.29) and (3.31) are compared over the interval x ∈ [x0 , l]:
∞
y(x, 0) = y0 (x, tc ) =

¯
sin Λk x γk sin cΛk tc + δk cos cΛk tc ,
¯

k=1
∞
¯
¯ n sin λn (l − x) =
⇒
D
sin Λk x γk sin cΛk tc + δk cos cΛk tc ,
¯
n=1
k=1

x ∈ [x0 , l]

∞

x ∈ [x0 , l]

l ∞
2
¯
sin Λk x sin λn (l − x) γk sin cΛk tc + δk cos cΛk tc dx
¯
l − x0 x
0 k=1
∞
2
λn
¯
=
sin Λk x0 γk sin cΛk tc + δk cos cΛk tc
¯
(3.32)
(−1)n
2 − λ2
l − x0
Λ
n
k
k=1

¯
⇒ Dn =

and,
∞
y(x, 0) = y0 (x, tc ) =
˙
˙

¯
cΛk sin Λk x γk cos cΛk tc − δk sin cΛk tc ,
¯

x ∈ [x0 , l]

k=1
∞
¯
¯n sin λn (l − x) =
⇒
cλn C
cΛk sin Λk x γk cos cΛk tc − δk sin cΛk tc , x ∈ [x0 , l]
¯
n=1
k=1
∞
l
Λk
2
¯
¯
sin Λk x sin λn (l − x) γk cos cΛk tc − δk sin cΛk tc dx
¯
⇒ Cn =
l − x0 x
λn
0 k=1
∞
2
kλn
¯
=
(−1)n
sin Λk x0 γk cos cΛk tc − δk sin cΛk tc
¯
(3.33)
2 − λ2
nl
Λ
n
k
k=1
∞

The total energy of the system after application of the constraint is given by the expression
T π2
Ec =
4(l − x0 )

∞

¯2 ¯ 2
n2 (Cn + Dn )

(3.34)

n=1

which is similar to the expression in Eq. (3.18). Using Eqs. (3.30) and (3.34), the change in

48

energy due to constraint application can be expressed as
T π2
Wc = (Ec − E0 ) =
4

∞
n=1

¯2 ¯ 2
(Cn + Dn ) (¯n + δn )
γ 2 ¯2
n2
−
(l − x0 )
l

(3.35)

When the constraint is removed at t = tr > tc , the equation of motion of the string is
described by Eq. (3.22), namely
∞
yr (x, t) =

sin Λk x γk sin cΛk t + δk cos cΛk t ,

x ∈ [0, l]

(3.36)

k=1
where γk , δk , k = 1, 2, · · · , can be computed from Eqs. (3.23) and (3.24) after replacing Cn
¯
¯
and Dn , n = 1, 2, · · · , with Cn and Dn , respectively, i.e.,

2
sin Λk x0
δk =
l

∞

(−1)n

λn

¯
¯
C sin cλn tr + Dn cos cλn tr
2 − λ2 n
Λ
n
n=1
k
∞
2
λn
2
¯
¯
C cos cλn tr − Dn sin cλn tr
sin Λk x0
(−1)n
γk =
2 − λ2 n
kπ
Λ
n
n=1
k

(3.37)

The energy of the string after removal of the constraint is given by Eq. (3.26), namely
T π2
Er =
4l

∞
k=1

2
2
k 2 (γk + δk )

(3.38)

Similar to section 3.5, we rely on simulation results to claim Wr = (Er − Ec ) = 0. The net
change in energy due to one cycle of constraint application and removal is therefore equal to
T π2
∆W = (Wc + Wr ) = Wc =
4

∞

¯2 ¯ 2
(Cn + Dn ) (¯n + δn )
γ 2 ¯2
n2
−
(l − x0 )
l
n=1
49

(3.39)

The effect of sequential application and removal of constraints can be investigated by defining
¯
γk = γk , δk = δk , k = 1, 2, · · · , and y0 (x, t) = yr (x, t), and repeating the procedure outlined
¯
above.

3.6.2

Control strategy for vibration suppression

From our simulation results in section 3.5, we have seen that maximum reduction in energy
is achieved if the constraint is applied when the displacement of the string is zero over
the interval [0, x0 ). Since the string was assumed to be vibrating in its first mode, greater
reduction in energy was achieved by simply choosing larger values of x0 . For arbitrary initial
conditions, the displacement of the string can be uniformly zero only over a small length
segment at some given time, and therefore, energy reduction can be ensured by choosing a
small value of x0 . A small value of x0 will result in less energy reduction but greater reduction
in energy can be achieved by applying the constraint incrementally N times, for a total
distance of Nx0 . Figure 3.4 shows a sensing scheme that can be employed to determine the
time when the scabbard can move incrementally. The removal of the constraint will involve
moving the scabbard to its original configuration after N applications of the constraint.

y

sensor

scabbard
2x0
(p + 1)x0

x

x0
px0

Figure 3.4: The sensor measures the displacement of the string at a distance x0 from the
tip of the scabbard. When this displacement is zero, the scabbard moves incrementally by
distance x0
50

In the sequel we provide the equations needed to describe the vibration of the string
during the process of N applications and removal of the constraint. The initial conditions of
the string are assumed to be arbitrary and hence its equation of motion is given by Eq. (3.29),
namely
∞

¯
sin Λk x γk sin cΛk t + δk cos cΛk t ,
¯

y0 (x, t) =

x ∈ [0, l]

(3.40)

k=1
(1) (2)
(N)
We assume that the constraint is applied N times at t = tc , tc , · · · , tc , and removed
(1)
(2)
(N)
at t = tr , where tc < tc < · · · < tc
< tr . After the first application of the constraint,
the equation of motion of the string has the form

y1 (x, t) =











0
∞

: x ∈ [0, x0 )

(1)
(1)
(1)
¯ (1)
¯ (1)
sin λn (l − x) Cn sin cλn t + Dn cos cλn t

(3.41)
: x ∈ [x0 , l]

n=1

¯
¯
The above equation is identical to Eq. (3.31) with the only difference that Cn and Dn ,
n = 1, 2, · · · , and λn now have the superscript (1), which indicates that the constraint was
(1)
applied at t = tc . Using the same procedure as before, y0 (x, t) and y1 (x, t) in Eqs. (3.40)
¯ (1)
¯ (1)
and (3.41) are compared over the interval x ∈ [x0 , l] to determine Cn and Dn
¯ (1)
Dn =

2
l − x0

∞

(−1)n

k=1

(1)
λn
(1) ¯
(1)
sin Λk x0 γk sin cΛk tc + δk cos cΛk tc
¯
(1) 2
Λ2 − (λn )
k

(3.42)

and

2
¯ (1)
Cn =
nl

∞
k=1

(−1)n

(1)
kλn
(1) ¯
(1)
sin Λk x0 γk cos cΛk tc − δk sin cΛk tc
¯
(1)
Λ2 − (λn )2
k

51

(3.43)

where
(p)
λn =

nπ
,
l − px0

p = 1, 2, · · · , N

(3.44)

After the second application of the constraint, the equation of motion of the string has the
form

y2 (x, t) =











0
∞

: x ∈ [0, 2x0 )

(2)
(2)
(2)
¯ (2)
¯ (2)
sin λn (l − x) Cn sin cλn t + Dn cos cλn t

: x ∈ [2x0 , l]

n=1

(3.45)

By comparing y1 (x, t) and y2 (x, t) in Eqs. (3.41) and (3.45) over the interval x ∈ [2x0 , l],
¯ (2)
¯ (2)
Cn and Dn , n = 1, 2, · · · , are obtained as follows
¯ (2) =
Dn

1
l − 2x0

¯ (2) = 1
Cn
cnπ

∞

∞

¯
¯ (1) sin cλ(1) t(2) + D (1) cos cλ(1) t(2) U(k, 1)
C
k c
k c
k
k

(3.46)

(1) ¯ (1)
(1) (2) ¯ (1)
(1) (2)
C
cos cλ tc − D sin cλ tc U(k, 1)
k
k
k
k
k

(3.47)

k=1
cλ

k=1

where

U(k, i) =

(i+1)
(i+1)
(i)
(i)
sin[(λn
sin[(λn
+ λ )(l − (i + 1)x0 )]
− λ )(l − (i + 1)x0 )]
k
k
−
(i+1)
(i)
(i+1)
(i)
λn
−λ
λn
+λ
k
k
(3.48)

The behavior of the string hereafter can be iteratively described by the relation

yp (x, t) =











0
∞

: x ∈ [0, px0 )

(p)
(p)
(p)
¯ (p)
¯ (p)
sin λn (l − x) Cn sin cλn t + Dn cos cλn t

n=1

: x ∈ [px0 , l]
(3.49)

52

¯ (p)
¯ (p)
where p = 2, 3, · · · , N, and Cn and Dn , n = 1, 2, · · · , are given by the expressions
¯ (p) =
Dn

1
l − px0

¯ (p) = 1
Cn
cnπ

∞

¯ (p−1) sin cλ(p−1) t(p) + D (p−1) cos cλ(p−1) t(p) U(k, p − 1)
¯
C
c
c
k
k
k
k

k=1

∞
cλ
k=1

(p−1) ¯ (p−1)
(p−1) (p) ¯ (p−1)
(p−1) (p)
C
cos cλ
tc − D
sin cλ
tc U(k, p − 1)
k
k
k
k
k
(3.50)

The scabbard is moved back to its original configuration after it reaches its maximum travel
distance Nx0 , and this completes one cycle of constraint application and removal. The
equation of motion of the string after constraint removal is given by Eq. (3.36), namely
∞
sin Λk x γk sin cΛk t + δk cos cΛk t ,

yr (x, t) =

x ∈ [0, l]

(3.51)

k=1
¯ (N) and Dn as follows
¯ (N)
where γk , δk , k = 1, 2, · · · , are defined in terms of Cn
2
sin Λk Nx0
δk =
l

∞

(−1)n

n=1

(N)
λn
(N)
(N)
¯ (N)
¯ (N)
Cn sin cλn tr + Dn cos cλn tr
(N)
Λ2 − (λn )2
k
(3.52)

and

2
sin Λk Nx0
γk =
kπ

∞
n=1

(−1)n

(N)
(λn )2
(N)
(N)
¯ (N)
¯ (N)
Cn cos cλn tr − Dn sin cλn tr
(N)
Λ2 − (λn )2
k
(3.53)

Similar to Eq. (3.39), the work done in one cycle of constraint application and removal can
53

be shown to be

∆W = (Wc + Wr ) = Wc =

∞
T π2
4

n=1




(N)
2 ¯2
¯n )2 + (Dn )2 } {¯n + δn }
¯ (N)
{(C
γ

n2 
−
(l − Nx0 )
l

(3.54)

Since the sensing scheme in Fig. 3.4 allows negative work to be done each time the scabbard
is moved incrementally, the net work done in each cycle of constraint application and removal
will be negative.
We define the work done after a number of constraint applications as the cost function, where
N
Ek+1 − Ek

J=

(3.55)

k=0
where, Ek+1 is total energy after k constraint application and Ek is the total energy before
k constraint application. And the cost is J = (E1 − E0 ) + (E2 − E1 ) + · · · + (EN−1 −
EN−2 ) + (EN − EN−1 ), then
J = EN − E0

(3.56)

To perform the control, we sense the states during the motion by the mean of a sensor. Where the sensor reads the displacement at x = (p + 1)x0 i.e y = y((p + 1)x0 , t) and
ˆ
˙
y = y((p + 1)x0 , t) as shown in Fig. 3.4, and the control strategy becomes
ˆ ˙


˙
 x : |ˆ | ≃ 0, |y | > 0, t = t
yk
ˆk
0
k
uk =

 0 : t=t +t
r
k

(3.57)

The condition of |ˆk | ≃ 0 guarantees the potential energy is minimum in the controlled
y
˙
segment, and |yk | > 0 guarantees the kinetic energy Ke > 0 where Ke will be the reduction
ˆ
54

in total system energy. tr is the time of constraint removal.

3.6.3

Numerical simulations

Similar to our example in section 3.5, we consider a string of length l = 4 m, mass per unit
length ρ = 0.25 kg/m, and tension T = 1 N, vibrating in its first mode with amplitude
A0 = 0.1 m. At t = 0 the string is assumed to pass through its mean position, and the
equation of motion of the string is therefore

π
π
y0 (x, t) = 0.1 sin( x) sin( t)
4
2

We simulate 12 cycles of constraint application and removal with constraint application
occurring incrementally in each cycle with x0 = 0.05l and N = 3, as discussed in section
3.6.2. In each cycle, the constraints are applied at the earliest instant of time, when the
string displacement is zero at the sensor2 location. For the first cycle, this is illustrated with
the help of Fig. 3.5, which shows the displacement and percentage energy loss of the string
at different time instances.
(3)
The constraint can be removed at any time tr , tr > tc , and there is flexibility in
choosing tr . In Table 3.2 and Fig. 3.6, we present simulation results for the two cases:
(A): the constraint is removed when the sensor measures zero displacement, and
(B): the constraint is removed when the sensor records the maximum displacement.
From the data in Table 3.2, it is clear that case (B) results in significantly higher reduction
in energy in comparison to case (A). Although constraint removal does not change the
2 The

sensor measures the displacement of the string at a distance x0 from the scabbard - see Fig.
3.4. Since x0 is small, a zero sensor reading implies that the string has almost zero displacement
over the entire interval x ∈ [0, x0 ).

55

overall energy of the string, the configuration of the string at the time of constraint removal
determines how energy of the constrained string is redistributed over the entire length;
and this has a significant effect on the rate of energy reduction in subsequent cycles of
constraint application. Beginning with the last plots in Fig. 3.5 and Fig. 3.6 shows the
configuration of the string at the instant of time when the constraint is removed in the first
cycle for cases (A) and (B). It also shows the maximum displacement of the string after the
completion of 12 cycles for the two cases. The maximum displacement plots corroborate
that the energy content of the string for case (B) is much less than that of case (A). The
plots in Fig. 3.6 can explain why case (B) results in higher reduction in energy. In case (B),
the constraint is removed when the string has the maximum deformation. This effectively
results in plucking of the string (without addition of energy), which excites the high-frequency
modes of the string and transfer of energy into these modes. When a significant fraction of
energy is funneled into the high-frequency modes, the scabbard is more effective in energy
reduction since these modes have higher energy density than the low-frequency modes for
the same amplitude of vibration. It can be verified from Fig. 3.6 that the string, after 12
cycles, predominantly vibrates in its first mode for case (A) but has many high-frequency
components for case (B).

Table 3.2: Comparison of percentage energy loss of the string for cases (A) and (B) over 12
cycles of constraint application and removal
No.
(A)
(B)

1
0.88
0.88

2
5.52
3.78

3
9.83
9.07

4
15.94
17.14

5
16.33
23.56

6
16.37
31.42

56

7
16.46
44.14

8
16.50
46.19

9
16.54
46.95

10
16.58
47.25

11
16.60
56.32

12
16.64
57.25

0.1

t = 0.0
y(m)

y(m)

0.1
0.0

energy loss = 0.00%

-0.1

0.0
scabbard

0.1

energy loss = 0.00%

0.1

(1)

t = tc = 2.0

-0.1

0.0
-0.1

y(m)

y(m)

0.1

t = 1.0

energy loss = 0.11%

t = 2.88

0.0

energy loss = 0.11%

-0.1
0.1

(2)

t = 4.76

y(m)

y(m)

t = tc = 3.9

0.0
scabbard

-0.1

0.1

(3)

t = tc = 5.75

0.0
-0.1

scabbard

0

1

energy loss = 0.22%

-0.1

y(m)

y(m)

0.1

energy loss = 0.22%

0.0

energy loss = 0.88%

2

3

4

x(m)

t = 6.61
energy loss = 0.88%

0.0
-0.1

0

1

2

3

x(m)

Figure 3.5: Displacement of the string and percentage energy loss at different time instants
during the first cycle of constraint application

57

4

y(m)

0.1

(A)

t = 6.61
energy loss = 0.88%

0.0
-0.1

scabbard

(B)
0.1

y(m)

y(m)

0.1
0.0

t = tr = 7.46 energy loss = 0.88%

-0.1

t = tr = 8.25 energy loss = 0.88%
-0.1
end of first cycle
(B)
end of 12 cycles
0.1

(A)

y(m)

y(m)

0.1
0.0
-0.1

0.0

t = 63.78
0

1

energy loss = 16.64%
2
x(m)

3

4

0.0
-0.1
0

t = 63.45
1

energy loss = 57.25%
2
x(m)

3

4

Figure 3.6: Displacement and percentage energy loss of the string at the end of the first
cycle and after completion of 12 cycles of constraint application and removal for cases (A)
and (B)

58

Chapter 4
Vibration Control of a String
Through Cyclic Application and
Removal of Constraints:
Experimental Verification

4.1

Introduction

In Chapter 3 we presented a control strategy for the vibration control of string by using
scabbard-like actuator. The scabbard-like actuator is difficult to implement in hardware
and an alternate actuation mechanism is therefore proposed in this chapter. In this chapter
we present an experimental control methodology that verifies the idea of the scabbardlike actuator method in which energy reduction is achieved through cyclic application and
removal of a constraint. In our method, a zero-displacement constraint is enforced at one
59

point on the string, close to one of its boundaries. This effectively results in two vibrating
strings, one of which is much shorter in length than the other. The vibration of the shorter
length string decays rapidly due to high internal damping and when the constraint is released,
the remaining energy of the string is redistributed over its entire length; this allows the cycle
of constraint application and removal to be repeated for vibration suppression. This chapter
is organized as follows. In Section 4.2 we present the mathematical model of the effect of
constraint application and removal, including change in energy of the string. This modeling
effort is similar to the one presented in Chapter 3 but includes damping; the addition of
damping allows us to compare simulation and experimental results. We use a stretched coil
extension spring as a string in our experiments. A coil extension spring can easily store
potential energy during stretching, and as a consequence has low damping - this makes it
well-suited for vibration control experiments. Although a stretched coil extension spring was
used earlier by Kashy, et al [77], we present experimental results to justify its use as a string.
These results are presented in Section 4.3. Section 4.3 also contains experimental results that
are used to identify the damping coefficient of the spring-string. Simulation and experimental
results of vibration suppression are presented in Section 4.4 - these results establish the
feasibility of vibration suppression through application and removal of constraints.

4.2
4.2.1

Cyclic Application and Removal of Constraint
Actuation and Sensing Mechanism

We propose the mechanism in Fig. 4.1 for constraint application and removal. It is comprised
of a pair of solenoids that can impose a zero velocity constraint at one point on the string,
60

located at a distance of x0 from the boundary. An optical sensor, located at a distance of η
from the solenoids, measures the distance of the string from its mean position. The solenoids
are activated when the sensor indicates that the string is passing through its mean position,
i.e. y(x0 +η) = 0. Since η is a small distance by design, the solenoids impose the zero velocity
constraint when y(x0 ) ≈ 0, i.e., the solenoids impose a zero displacement constraint. This
is necessary to ensure that the constraint does not alter the equilibrium configuration of the
string. The zero displacement constraint virtually results in two vibrating strings over the
intervals x ∈ [0, x0 ) and x ∈ (x0 , l]. The distance x0 is a small fraction of the total length l
by design and consequently the string vibrating over the interval x ∈ [0, x0 ) will have a much
higher natural frequency and damping. This implies that the energy of the shorter string
will dissipate quickly and naturally, and subsequent de-activation of the solenoids will result
in redistribution of the remaining energy of the string over its entire length. Each cycle of
solenoid activation and de-activation (constraint application and removal) will reduce the
energy of the string and eventually result in vibration suppression.

The actuation mechanism in Fig. 4.1 is fundamentally similar to the scabbard-like actuator in Chapter 3. By applying a zero displacement constraint over a small length of string
that is passing through its mean position, the scabbard proposes to remove kinetic energy
of the small string segment through impact. In contrast, the mechanism in Fig. 4.1 imposes
a zero displacement constraint at one point of the string. This point is proximal to one
fixed support and therefore the energy trapped in the string between the point of application of the constraint and the proximal support is dissipated naturally. The time required
for natural dissipation depends on the distance of the point of application of the constraint
and the proximal boundary - a smaller distance results in faster dissipation and vice versa.
61

The scabbard proposes to dissipate the energy of the small string segment instantaneously
through impact, but it is difficult to implement.

z

solenoid
top view
x
plunger
optical sensor
stator
vibrating string
x0

y

η

side view
x

distance measured

optical sensor
ℓ

Figure 4.1: Actuation and sensing mechanism for application and removal of constraint

4.2.2

Mathematical Modeling

4.2.2.1

Linear model of vibrating string with damping

In this section we consider a mathematical model of the string and compute the change
in energy of the string due to constraint application and constraint removal. To include
damping, we consider the mathematical model from [78]

T

∂2y
∂2 y
∂y
=ρ
+C
2
2
∂t
∂x
∂t
62

(4.1)

where T is the tension, ρ is the mass per unit length, and C is the damping coefficient of
the string. Through separation of variables y(x, t) = X(x) φ(t), we get
1 1 ∂2φ
∂φ
1 ∂2X
=
+K
= −λ2 ,
2
2 ∂t2
X ∂x
φc
∂t

K

C/ρ,

c

T /ρ

where λ is constant. The above equation gives

∂2φ
∂t2

∂2X
+ λ2 X = 0 ⇒ X(x) = {α cos λx + β sin λx}
2
∂x
∂φ
+K
+ c2 λ2 φ = 0 ⇒ φ(t) = e−ζωt A cos(ω 1 − ζ 2 )t + B sin(ω
∂t

1 − ζ 2 )t

and leads to the general solution

y(x, t) = e−ζωt {α cos λx + β sin λx} A cos(ω

1 − ζ 2 )t + B sin(ω

1 − ζ 2 )t

(4.2)

where
ω

cλ,

ζ

K/(2cλ)

Assuming the string to be of length ℓ and applying the boundary conditions y(0, t) = y(ℓ, t) =
0 or X(0) = X(ℓ) = 0, we get

X(x) = βn sin λn x,

λn =

nπ
,
ℓ

n = 1, 2, · · · , ∞

The general solution then takes the form
∞
y(x, t) =

φn (t) sin λn x
n=1
63

(4.3)

where
φn (t) = e−ζn ωn t An cos(ωn

2
1 − ζn )t

2
1 − ζn )t + Bn sin(ωn

and
ωn

cλn ,

ζn

K/(2cλn )

The total energy of the string is the sum of its potential and kinetic energy
ρ ℓ ∂y(x, t) 2
T ℓ ∂y(x, t) 2
E(t) =
dx +
dx
2 0
∂t
2 0
∂x

(4.4)

Substituting Eq. (4.3) into Eq. (4.4), we get

E(t) =

ℓ





∞



∞

ρ
˙
˙

φn (t) sin λn x 
φm (t) sin λm x dx
2 0
n=1
m=1



∞
ℓ ∞
T

+
λn φn (t) cos λn x 
λm φm (t) cos λm x dx
2 0
n=1
m=1

Using the orthogonality property of the modes and simplifying, we get

ρℓ
E(t) =
4

4.2.2.2

∞
n=1

Tℓ
˙
φ2 (t) +
n
4

∞

λ2 φ2 (t)
n n

(4.5)

n=1

Effect of constraint application

To obtain the equation of motion of the string after the constraint is applied, we solve
Eq. (4.1) subject to initial and boundary conditions consistent with the constraint. Assuming
¯
¯
the solution to be of the form yc (x, t) = X(x) φ(t), where

¯
¯
¯
¯
X(x) = α cos λx + β sin λx
¯
64

and
¯¯
¯
¯
ω
φ(t) = e−ζ ω t A cos(¯

¯
¯
1 − ζ 2 )t + B sin(¯
ω

¯
1 − ζ 2 )t

¯
¯
application of the boundary conditions yc (x0 , t) = yc (ℓ, t) = 0 or X(x0 ) = X(ℓ) = 0 yields

¯
¯
¯
X(x) = βm sin λm (ℓ − x),

¯
λm =

mπ
,
ℓ − x0

m = 1, 2, · · · , ∞

The constrained solution then takes the form
∞
yc (x, t) =

¯
¯
φm (t) sin λm (ℓ − x)

(4.6)

m=1
where

¯ ¯
¯
¯
ω
φm (t) = e−ζm ωm t Am cos(¯ m

¯2
¯
1 − ζm )t + Bm sin(¯ m
ω

¯2
1 − ζm )t

and
ωm
¯

¯
cλ m ,

¯
ζm

¯
K/(2cλm)

If the constraint is applied at time t = tc , we have the following conditions:

A1. yc (x, 0) = y(x, tc ). From Eqs. (4.3) and (4.6) we can write
∞

¯
¯
Am sin λm (ℓ − x) =

m=1

∞
φn (tc ) sin λn x

(4.7)

n=1

¯
where φn (t) was defined after Eq. (4.3). Multiplying both sides by sin λm (ℓ − x) and
65

integrating with respect to x from x0 to ℓ, we get

¯
Am =
=

2
ℓ − x0
2
ℓ − x0

∞
φn (tc )
n=1
∞

ℓ
x0

¯
sin λm (ℓ − x) sin λn x dx

(−1)m φn (tc )

n=1

¯
λm
sin λn x0
¯
λ2 − λ2
n
m

(4.8)

A2. yc (x, 0) = y(x, tc ). From Eqs. (4.3) and (4.6) we can again write
˙
˙
∞
ωm
¯

∞

¯2 ¯
¯ ¯
¯
1 − ζm Bm − ζm Am sin λm (ℓ − x) =

m=1

˙
φn (tc ) sin λn x

(4.9)

n=1

¯
Multiplying both sides by sin λm (ℓ − x), integrating with respect to x from x0 to ℓ,
we get

¯
Bm =
ωm
¯


¯m
λ
1
˙
¯ ¯ ¯
ζm ωm Am + 2
sin λn x0 
(−1)m φn (tc )
¯
ℓ − x0
λ2 − λ2
¯2
n
m
1 − ζm
n=1
(4.10)


∞

¯
where Am is defined by Eq. (4.8).

The energy of the string after application of the constraint is given by the relation

Ec (t) =

ρ ℓ ∂yc (x, t) 2
T ℓ ∂yc (x, t) 2
dx +
dx
2 x
∂t
2 x
∂x
0
0

(4.11)

Substituting Eq. (4.6) and simplifying using orthogonality property of the modes, we get

ρ(ℓ − x0 )
Ec (t) =
4

∞
m=1

T (ℓ − x0 )
˙
¯2
φm (t) +
4
66

∞
m=1

¯ ¯
λ2 φ2 (t)
m m

(4.12)

The work done by application of the constraint can be simply computed as

Wc = Ec (tc ) − E(tc )

4.2.2.3

(4.13)

Effect of constraint removal

After removal of the constraint, the equation of motion of the string will be described by the
relation
∞
yr (x, t) =

e−ζn ωn t µn cos(ωn

2
1 − ζn )t + νn sin(ωn

2
1 − ζn )t sin λn x

(4.14)

n=1
where
ωn

cλn ,

ζn

K/(2cλn )

If the constraint is applied at time t = tr , we have the following conditions:

A1. yr (x, 0) = yc (x, tr ). From Eqs. (4.6) and (4.14) we can write
∞

∞

¯
¯
φm (tr ) sin λm (ℓ − x)

µn sin λn x =

(4.15)

m=1

n=1

¯
where φm (t) was defined after Eq. (4.6). Multiplying both sides by sin λn x and integrating with respect to x from 0 to ℓ, we get

µn =
=

2
ℓ
2
ℓ

∞
m=1
∞

¯
φm (tr )

ℓ
x0

¯
sin λm (ℓ − x) sin λn x dx

¯
(−1)m φm (tr )

m=1
67

¯
λm
sin λn x0
¯
λ2 − λ2
m
n

(4.16)

A2. yr (x, 0) = yc (x, tr ). From Eqs. (4.6) and (4.14) we can again write
˙
˙
∞
ωn

2
1 − ζn νn − ζn µn sin λn x =

n=1

∞

˙
¯
¯
φm (tr ) sin λm (ℓ − x)

(4.17)

m=1

Multiplying both sides by sin λn x, integrating with respect to x from 0 to ℓ, we get

νn =
ωn


¯m
λ
1
˙
¯
ζn ωn µn + 2
sin λn x0 
(−1)m φm (tr )
2 − λ2
¯
ℓ
λn
2
m
1 − ζn
m=1


∞

(4.18)

where µn is defined by Eq. (4.16).
Once the values of µn and νn have been determined, we update the values of An and Bn
as follows
An = µn ,

Bn = νn ,

n = 1, 2, · · · , ∞

The change in the energy of the string due to constraint removal can now be computed as

Wr = E(tr ) − Ec (tr )

(4.19)

The update in the values of An and Bn also allow us to repeat the process of constraint
application and removal.

4.3
4.3.1

Experimental Setup
Coil Extension Spring - A Substitute for the String

A traditional string has high damping and is not well-suited for vibration control experiments.
As a substitute for the string, we used a stretched coil extension spring which has low
68

damping. The material and geometric specifications of the spring [79] used in experiments
are provided in Table 4.1:
Table 4.1: Specifications of coil extension spring
Material

Outside diameter

Inside diameter

Wire diameter

Steel

2.032 mm

1.016 mm

0.508 mm

A sample spring of free length 2.286 m (90 in) was found to have a mass of 0.034 kg. The
mass per unit length of the spring is therefore ρs = 0.0149 kg/m. Using a known load, the
stiffness of the spring in tension was computed to be ks = 150.42 N/m. The experimental
setup is shown in Fig. 4.2. In this setup, the distance between the fixed supports is ℓ = 1.64
m. The spring clamped between these supports had a free length of 1.525 m; the tension
in the spring is therefore T = ks ∆ = 150.42 × (1.64 − 1.525) = 17.29 N. The spring was
given an arbitrary initial displacement and the sensor was used to measure the displacement
of one point on the string - see Fig. 4.2. The sensor data was recorded for 100 sec at 1000
samples/sec. The natural frequencies of free vibration of the spring were computed from the
FFT of the sensor data and are shown in Table 4.2. The data indicates ω2 ≈ 2ω1 , ω3 ≈ 3ω1
and ω4 ≈ 4ω1 . This is characteristic of vibrating strings and therefore the coil extension
spring is a good substitute for the string.

69

fixed support actuation mechanism

guide for generating repeatable
initial condition
string
solenoids

(a)

optical sensor

solenoid

(b)

Figure 4.2: (a) Experimental setup showing the string between two fixed supports and the actuation mechanism, (b) A close-up
view of the actuation mechanism and sensor

70

Treating the spring as a string, the mass per unit length of the string was computed from
the first natural frequency in Table 4.2, using the relation

ρ=T

π 2
= 0.0164 kg/m
ℓ ω1

(4.20)

This is higher than the mass per unit length of the spring, ρs = 0.0149 kg/m. In our
simulations, we will use the value of mass per unit length of the equivalent string, ρ = 0.0164
kg/m.

Table 4.2: Experimentally determined natural frequencies of the spring/string in rad/s
Units

ω1

ω2

ω3

ω4

rad/s

62.07

124.15

186.29

248.37

Hz

9.88

19.76

29.65

39.53

4.3.2

Internal Damping of the String

4.3.2.1

Estimation of Damping Ratio

In order to compare simulation and experimental results on vibration suppression, it is necessary to determine the internal damping of the string. To this end, we collected experimental
data on the string vibrating freely in the first mode and computed the damping ratio ζ using
the method of logarithmic decrement
−1/2
ζ = 1 + (2π/δ)2
,
71

δ=

1
ln
k

xi
xk+i

y(x0 + η) (m)

0.004
1.194 × 10−3

0.0

-0.004

3.088 × 10−3

0

20

0

20

40

60

40

(a)

60

y(x0 + η) (m)

0.004

0.0

-0.004

time (sec)
(b)
Figure 4.3: Free vibration of a string: (a) Experimental results for the string in Fig. 4.2 (b)
Simulation results obtained using a damping ratio computed from the experimental results
in (a)
The experimental data is shown in Fig. 4.3 (a). For this data, we had xi = 3.088 × 10−3
m, xk+i = 1.194 × 10−3 m and k = 591; this gives ζ = 2.56 × 10−4 . Using this value of
ζ and values of T and ρ obtained experimentally, we simulated free vibration of the string
in its first mode. The results are shown in Fig. 4.3 (b); they match quite well with the
experimental results in Fig. 4.3 (a) and provide confidence in our mathematical model.

4.3.2.2

Variation of Damping Coefficient with Length

The damping ratio ζ = 2.56 × 10−4 was determined experimentally for the string vibrating
in its first mode. The product of the damping ratio and the natural frequency is therefore

ζω = 2.56 × 10−4 × 62.07 = 0.0159 rad/s
72

The product ζω was assumed to be constant in our mathematical model since it is proportional to the damping coefficient C (ζω = C/2ρ). To investigate the validity of this
assumption, we present experimental results on variation of the product ζω with change in
length of the string. The results, shown in Fig. 4.4, indicate that ζω initially remains constant but increases rapidly thereafter as the length is decreased. The trend is consistently
observed with three different strings, namely:

(1) string described in Section 3.1 with stiffness ks = 150.42 N/m and tension T = 17.29
N;

(2) string with spring coil diameter smaller than that of string (1) but of higher stiffness
ks = 284.90 N/m and with higher tension T = 34.50 N; and
(3) string with spring coil diameter similar to that of string (1) and same tension as that
of string (1) but of lower stiffness ks = 47.5 N/m.
0.20

ζω rad/s

0.15
0.10
(2)

(3)

0.05

(1)

0.00
20

60

40

80

100

Percentage of length ℓ (ℓ = 1.64 m)
Figure 4.4: Effect of changing the length of the string on the product ζω for three different
strings
73

4.4

Simulation and Experimental Results on Vibration
Suppression

For simulations, we used the mathematical model in Eq. (4.3) with the following parameter
values
ℓ = 1.64 m,

T = 17.29 N,

ρ = 0.0164 kg/m,

ζ = 2.56 × 10−4

which were obtained experimentally. To be consistent with our experimental setup, we used

x0 = 0.1 ℓ = 0.164 m,

η = 0.04 m

The string was given the displacement shown in Fig. 4.5 with d = ℓ/2 = 0.82 m and h =
0.01 m, and released from rest at the initial time. This initial condition, which can be
mathematically described as follows

y(x, 0) = f (x) =


 h
 − x

 d


0≤x≤d
,





 − h (ℓ − x)
ℓ−d

y(x, 0) = g(x) = 0
˙

d≤x≤ℓ

excites all modes of vibration. The modal coefficients can be determined from the expression
for y(x, t) in Eq. (4.3) and its derivative, as follows
∞
An sin λn x = f (x)
n=1
⇒

An = −

ℓ
nπ
h
2hℓ2
nπ
nπ
2 h d
sin
x dx +
x dx =
d
x sin
(ℓ − x) sin
ℓ d 0
ℓ
ℓ−d d
ℓ
ℓ
n2 π 2 d(ℓ − d)
74

y
actuator

side view
x

d

η

h

x0
optical sensor

guide for initial condition

Figure 4.5: Initial displacement of the string used for both simulations and experiments.
The guide for generating repeatable initial conditions in experiments can be seen in Fig. 4.2
∞
ωn
n=1
⇒

Bn =

2
1 − ζn Bn − ζn ωn An sin λn x = g(x) = 0
ζn
2
1 − ζn

An

(4.21)

The constraint was applied when the sensed point on the string passes through its mean
position, i.e., y(x0 + η) ≈ 0. The constraint was released exactly 0.5 sec after the constraint
was applied. The time interval of 0.5 sec is provided to allow vibration of the string over
the interval [0, x0 ) to dissipate naturally. In Fig. 4.6 (a), we present simulation results for
N = 12 cycles of constraint application and removal with a time interval of 0.5 sec between
cycles and with n = 3 modes. The experimental results are shown in Fig. 4.6 (b). The initial
condition used in our experiment was identical to that used in simulations and facilitated by
the guide for generating repeatable initial conditions shown in Fig. 4.2.
The following observations can be made from the simulation and experimental results in
Fig. 4.6:

• The time duration of active control is almost identical in simulations and experiments,
slightly greater than 12 sec. It includes 12 cycles of constraint application and removal
with a minimum duration of 1.0 sec per cycle.
75

y(x0 + η) (m)

0.004
free vibration
0.0
t ≈ 12.70 sec

-0.004
(a)
y(x0 + η) (m)

0.004
free vibration
0.0
t ≈ 12.25 sec

-0.004
0

10

(b)

20

30

time (sec)
Figure 4.6: Results of vibration suppression obtained from (a) simulations and (b) experiments
• For the period of time that the constraint is active, the amplitude of vibration is
significantly lower than when the constraint is not active. This is due to the fact that
constraint application results in a temporary boundary at x = x0 and Fig. 4.6 plots
the displacement of the string at x = x0 + η, i.e., at a small distance η from this
boundary.

• In simulations, the constrained string vibrates about the mean position of the unconstrained string - see plot (a), but in experiments it has a positive or a negative offset
- see plot (b). This is because it is possible to apply the constraint exactly when
the string is passing through its mean position in simulations; the same is difficult to
achieve in experiments due to reaction time of the solenoids. It may be noted that
the solenoids fail to keep the string completely constrained during the first cycle of
76

experiments. At this time the string has a large amount of energy and the friction
forces between the string and the solenoid plungers are insufficient to prevent some
slipping.

• At the initial time, we had d = 0.82 m, h = 0.01 m, x0 = 0.164 m and η = 0.04 m.
This corresponds to y(x0 + η) = 2.48 × 10−3 m - see Fig. 4.5. The initial amplitudes
of vibration in Fig. 4.6 matches well with this number.

• The ratio of the initial amplitude of vibration to the amplitude of vibration immediately
after active control is terminated is 4.29 for simulation results in plot (a) and 4.98 for
experimental results in plot (b). Assuming only the first mode of vibration to be
present, these numbers correspond to ≈ 94.5% and ≈ 95.9% of energy dissipation for
simulations and experiments, respectively. The large percentage of energy dissipation
is indicative of the effectiveness of the control strategy.

• The simulation results in plot (a) indicate the presence of higher modes after termination of active control. The experimental results in plot (b) indicate that only the first
mode is present and the higher modes have dissipated. This suggests that the damping ratio ζ may not be constant across modes; a higher damping ratio for the higher
modes may have contributed to higher energy dissipation in experiments compared to
simulations.

77

4.5
4.5.1

Additional Investigations
Efficacy of Vibration Suppression - A Parametric Study

The efficacy of vibration suppression through constraint application and removal depends on
several parameters. In this section, we present a parametric experimental study to illustrate
the effect of three specific parameters, namely:
A1. (x0 /ℓ): Location of constraint from proximal boundary
A2. tr : Duration of time for which the constraint is applied, and
A3. N: Number of times the constraint is applied and removed
on the settling time. We define the settling time as the time required for the amplitude of
vibration to decay to 20% of its initial value. For a string vibrating predominantly in its
first mode, this corresponds to the time required for dissipation of ≈ 96% of the energy.
The results of settling time are presented in Table 4.3. These results indicate that larger
number of cycles of constraint application and removal (higher value of N) results in smaller
settling time, i.e., faster vibration suppression, independent of the values of (x0 /ℓ) and tr .
Also, larger values of (x0 /ℓ) require proportionately longer time of application tr , for the
shortest settling time. It can be seen that independent of the value of N, tr = 0.5 results
in the the shortest settling time for (x0 /ℓ) = 0.05, tr = 1.0 results in the the shortest
settling time for (x0 /ℓ) = 0.10, and tr = 1.5 results in the the shortest settling time for
(x0 /ℓ) = 0.15. This because a larger value of x0 lowers the frequency of vibration of the
constrained string segment and increases the time required for natural dissipation of energy
in this segment. From the above discussion, it is clear that tr = 0.5 and N = 12 provides the
78

Table 4.3: Effect of three different parameters: (x0 /ℓ), tr , and N on the settling time
tr (sec)
0.5
1.0
1.5

(x0 /ℓ)=0.05
N = 6 N = 12
51.61
43.44
57.42
46.63
60.32
47.44

(x0 /ℓ)=0.10
N = 6 N = 12
64.52
38.12
54.70
30.43
63.66
43.22

(x0 /ℓ)=0.15
N = 6 N = 12
72.21
65.32
70.46
53.18
62.44
28.07

optimal settling time for (x0 /ℓ) = 0.05; tr = 1.0 and N = 12 provides the optimal settling
time for (x0 /ℓ) = 0.10; and tr = 1.5 and N = 12 provides the optimal settling time for
(x0 /ℓ) = 0.15. These three cases are highlighted in Table 4.3. Among them, (x0 /ℓ) = 0.15
has the shortest settling time. Simulation and experimental results for shaded results in
table 4.3 are presented in Figs. 4.7, 4.8 and 4.9 respectively. These figures show good match
in reduction of vibration amplitude and settling time ts .

79

y(x0 + η) (m)

0.004
free vibration
0.0
ts ≈ 46.7180 sec
−0.004
(a)

y(x0 + η) (m)

0.004
free vibration
0.0
ts ≈ 43.44 sec
−0.004
0

10

20

(b)

30

40

50

time (sec)
Figure 4.7: Results of vibration suppression obtained from (a) simulations and (b) experiments, with (x0 /ℓ) = 0.05, tr = 0.5
s, and N = 12

80

y(x0 + η) (m)

0.004
free vibration
0.0
ts ≈ 27.32 sec
−0.004
(a)

y(x0 + η) (m)

0.004
free vibration
0.0
ts ≈ 30.43 sec
−0.004
0

10

20

(b)

30

40

50

time (sec)
Figure 4.8: Results of vibration suppression obtained from (a) simulations and (b) experiments, with (x0 /ℓ) = 0.10, tr = 1 s,
and N = 12

81

y(x0 + η) (m)

0.004

free vibration

0.0
ts ≈ 26.66 sec
−0.004
(a)

y(x0 + η) (m)

0.004
free vibration
0.0
ts ≈ 28.07 sec
−0.004
0

10

20

(b)

30

40

50

time (sec)
Figure 4.9: Results of vibration suppression obtained from (a) simulations and (b) experiments, with (x0 /ℓ) = 0.15, tr = 1.5
s, and N = 12

82

4.5.2

Excitation of Higher Modes

The energy dissipation mechanism proposed in this chapter relies on trapping energy in a
short segment of the string, which is dissipated naturally due to high internal damping. The
energy of the string is also dissipated through excitation of high-frequency modes in the
longer segment of the string. The high-frequency modes are excited each time the constraint
is applied or removed; this can be explained mathematically by the mapping between Fourier
coefficients describing the motion of the string in the constrained and unconstrained states.
The high-frequency modes also have high rates of damping, and although they do not decay
out as quickly as the vibration in the short string segment, energy transfer to these modes
facilitate vibration suppression.
3.0
2.0

unconstrained string

1.0
0.0
(a)
0.40
0.30

constrained string

0.20
0.10
0.00
0

20

40

60

80

(b)
frequency (Hz)
Figure 4.10: Fast Fourier-Transform results showing multiple modes of vibration of the (a)
unconstrained string and (b) constrained string
83

We present results of fast Fourier-Transform (FFT) of the string vibration collected over
3 sec, both after application and after removal of the constraint. These plots shown in
Fig. 4.8 indicate that at least six higher modes are excited. The FFT plots indicate that
the frequencies of the constrained string are slightly higher than those of the unconstrained
string. This is because of the shorter length of the constrained string. The amplitudes of the
constrained string in the FFT plot are however lower than those of the unconstrained string.
This is because the sensor is located closer to the boundary when the string is constrained
than when it is unconstrained.

84

Chapter 5
Dynamics of a Circular Membrane
with an Eccentric Circular Areal
Constraint: Analysis and Accurate
Simulations
The studies on vibration of constrained circular membrane are restrict to the determination
of the eigenfrequencies where most of the results in the literature typically compute the first
few eigenfrequencies of vibration but they do not provide a method to simulate the dynamics
accurately. The simulation of the dynamics requires computation of modal coefficients using
orthogonality properties of the modes. Although orthogonality properties of modes are
well-known, they have not been established mathematically for constrained membranes. We
present orthogonality properties of all modes of a circular membrane with an internal circular
areal constraint for the first time in this chapter. The computation of the modal coefficients
85

require accurate computation of the mode shapes and this requires proper choice of the
number of angular modes. This chapter is organized as follows. A formal problem statement
and a list of assumptions are provided in Section 5.1. The expression for the symmetric and
antisymmetric modes are derived in Section 5.2 using Lin’s analysis [46]. In Section 5.3,
we establish orthogonality between the symmetric modes, the antisymmetric modes, and
between symmetric and antisymmetric modes. In Section 5.4 we provide the algorithm for
accurate computation of the mode shapes, and in Section 5.5 we provide the algorithm for
computation of modal coefficients. The dynamics of a constrained membrane is simulated
at the end of Section 5.5.

5.1

Problem Statement and Assumptions

Consider the circular membrane of radius R and outer boundary Γo , shown in Fig. 5.1. A
circular area of radius a and boundary Γi , located at a distance d from the center of the
membrane, has zero displacement at all times; this constrains the dynamics of the membrane. To study the vibration of the active region of the membrane, we make the following
assumptions:

A1. The membrane is homogeneous and has a constant mass per unit area equal to µ.
The tension in the membrane is equal to T and remains constant at all times. The
membrane undergoes transverse vibration and it is not affected by gravity.

A2. The amplitude of oscillation of the membrane is small and its motion can be expressed
86

by the standard relation in polar coordinates [75]:
∂ 2 Z 1 ∂Z
1 ∂2Z
1 ∂2Z
+
=
,
+
r ∂r
∂r 2
r 2 ∂θ2
c2 ∂t2

c=

T /µ

(5.1)

where Z = Z(r, θ, t) is the transverse displacement. The distance of any point P is
measured from the center of the inner circle O and is denoted by r; θ denotes the angle
that vector r makes with the x-axis.

A3. The transverse displacement of all points in the area enclosed by the inner circle is
zero, i.e., Z(r, θ, t) = 0 for r ∈ [0, a], θ ∈ [0, 2π]. The displacement of the membrane is
also zero along the outer boundary Γo .

y

P

R

r

ρ

Γi
ψ

θ
a

O

x

d
Γo

Figure 5.1: A circular membrane with an inner circular areal constraint

87

5.2
5.2.1

Background - Modes of Vibration
General Solution

A general solution to Eq. (5.1) is obtained by separation of variables. Assuming Z(r, θ, t) =
ϕ(r, θ) F (t) and rearranging terms, we get
1 ∂ 2 ϕ 1 ∂ϕ
1 1 ∂2F
1 ∂2ϕ
+
}=
{
+
ϕ ∂r 2
r ∂r
r 2 ∂θ2
c2 F ∂t2

(5.2)

Let −k 2 be the separation constant. Then, the right-hand side of Eq. (5.2) gives
1 1 ∂2F
= −k 2 ,
2 F ∂t2
c

If we define the natural frequency ω

or

∂2F
+ c2 k 2 F = 0
2
∂t

ck 1 , we get

F (t) = A cos ωt + B sin ωt

(5.3)

where A and B are constants whose values depend on initial conditions. Equating the
left-hand side of Eq. (5.2) to −k 2 and rearranging terms, we get

r2

∂ϕ
∂2ϕ
∂2ϕ
+r
+ k2 r2 ϕ = −
∂r
∂r 2
∂θ2

(5.4)

To separate the radial and angular terms in Eq. (5.4), we substitute ϕ(r, θ) = U(r) V (θ);
1 The

separation constant k is related to the natural frequency by the relation k = ω/c and has
the unit of rad/m. Since c is a constant, the value of k is proportional to the frequency. In the
sequel we will refer to k as the eigenfrequency.

88

this yields
1 2 ∂2U
∂U
1 ∂2V
+r
{r
+ k2 r2 } = −
U
∂r
V ∂θ2
∂r 2

(5.5)

By choosing m2 to be the separation constant, we get
1 ∂2V
= −m2 ,
2
V ∂θ

⇒

V (θ) = Cm cos mθ + Dm sin mθ

(5.6)

where Cm and Dm are constants whose values depend on initial conditions. Since we have
∂V
∂V
(π) =
(−π)
∂θ
∂θ

V (π) = V (−π),

it can be shown that m2 takes integer values 0, 1, 2, · · · . From the left-hand side of Eq. (5.5)
we get
r2

∂U
∂2U
+ (k 2 r 2 − m2 )U = 0
+r
2
∂r
∂r

(5.7)

For any integer value of m, the general solution to Eq. (5.7) is as follows:

U(r) = βm Jm (kr) + γm Ym (kr)

where βm and γm are constants, and Jm (kr) and Ym (kr) are Bessel functions [75] of the
first and second kind, respectively. The general solution to the equation of motion of the

2 The

integer m is the number of angular nodes. For an unconstrained membrane, it is equal to
the number of diametrical lines along which the membrane has zero displacement

89

circular membrane can now be written as
∞
Z(r, θ, t) =

{βm Jm (kr) + γm Ym (kr)} ×

m=0
{Cm cos mθ + Dm sin mθ} {A cos ωt + B sin ωt}

(5.8)

The solution in Eq. (5.8) can be rewritten as

Z(r, θ, t) = ϕs (r, θ) + ϕc (r, θ) {A cos ωt + B sin ωt}

(5.9)

where

ϕs (r, θ) =
ϕc (r, θ) =

∞
{β1m Jm (kr) + γ1m Ym (kr)} sin mθ
m=1
∞
{β2m Jm (kr) + γ2m Ym (kr)} cos mθ

(5.10)

m=0
and β1m = βm Dm , γ1m = γm Dm , β2m = βm Cm and γ1m = γm Cm . The constants
in Eq. (5.9), namely, k, A, B, m, β1m , γ1m , β2m , γ2m , can be determined by applying
boundary conditions and initial conditions.

5.2.2

Boundary Conditions

The circular membrane is assumed to be fixed from the outer boundary and with presence of
the internal circular constraint, the outer boundary of the constraint is assumed to be fixed.
The boundary conditions can be written as:
90

A1. The membrane has zero displacement for r = a. From Eq. (5.9) we get

Z(a, θ, t) = ϕs (a, θ) + ϕc (a, θ) {A cos ωt + B sin ωt} = 0

This implies {ϕs (a, θ) + ϕc (a, θ)} = 0, or

∞
{β1m Jm (ka) + γ1m Ym (ka)} sin mθ
m=1
∞
{β2m Jm (ka) + γ2m Ym (ka)} cos mθ = 0

+

(5.11)

m=0
A2. The membrane has zero displacement for all r values for which ρ = R. This implies
∞
β1m Jm (kr|ρ=R ) + γ1m Ym (kr|ρ=R ) sin mθ +
m=1
∞
β2m Jm (kr|ρ=R ) + γ2m Ym (kr|ρ=R ) cos mθ = 0

(5.12)

m=0
The solutions to Eqs. (5.11) and (5.12) were investigated by several researchers [2], [46],
[47], [3] and the eigenfrequencies were determined for the even (symmetric) and the odd
(antisymmetric) modes3 .

5.2.3

Even or Symmetric Modes

The symmetric modes are associated with the cosine term cos mθ in Eq. (5.10) and are
denoted by ϕc (r, θ). The superscript “c” is used to represent the cosine terms. Since the
3A

background about simple circular unconstrained membrane can be found in appendix A.2

91

symmetric modes satisfy the differential equation independently, the boundary conditions in
Eqs. (5.11) and (5.12) give
∞
{β2m Jm (ka) + γ2m Ym (ka)} cos mθ = 0
m=0

∞

β2m Jm (kr|ρ=R ) + γ2m Ym (kr|ρ=R ) cos mθ = 0

(5.13)

m=0
Equation (5.13) cannot be solved directly for k since the second equation contains r|ρ=R ,
which is a function of θ. From Fig. 5.1, we can show

r|ρ=R = r(θ) =

R2 + d2 − 2Rd cos ψ =

R2 − d2 sin2 θ − d cos θ

(5.14)

To solve Eq. (5.13) analytically, Lin [46] used Graf’s addition formula [80] to transform the
Bessel functions from (r, θ) coordinates to (ρ, ψ) coordinates. This transformation is given
below

Jm (kr) eimθ =
Ym (kr) eimθ =

∞

Jm+q (kρ) Jq (kd) ei(m+q)ψ

q=−∞
∞

Ym+q (kρ) Jq (kd) ei(m+q)ψ

(5.15)

q=−∞
Using the expressions in Eq. (5.15), Eq. (5.13) is rewritten as follows
∞
{β2m Jm (ka) + γ2m Ym (ka)} cos mθ = 0
∞

∞

m=0
β2m Jm+q (kR) + γ2m Ym+q (kR) Jq (kd) cos(m + q)ψ = 0

m=0 q=−∞

92

(5.16)

Expanding the inner summation and using the property of Bessel functions, we get4
∞
{β2n Jn (ka) + γ2n Yn (ka)} cos nθ = 0
∞

∞

n=0
ǫn {β2m Jn (kR) + γ2m Yn (kR)} Jn−m (kd) + (−1)m Jn+m (kd) cos nψ = 0

m=0 n=0
(5.17)

where ǫn is defined as



 1/2 if n = 0

ǫn =

 1

if n > 0

Expanding Eq. (5.17) from m = 0 to N horizontally and from n = 0 to N vertically, gives
the matrix equation
P cX c = 0

(5.18)

where P c ∈ R(2N+2)×(2N+2) and A ∈ R(2N+2) are given by the relations


J0 (ζ1 )
Y0 (ζ1 )
0
0


 1
 J0 (ζ2 )H c 1 Y0 (ζ2 )H c 1 J0 (ζ2 )H c 1 Y0 (ζ2 )H c
00 2
00 2
01 2
01
 2



0
0
J1 (ζ1 )
Y1 (ζ1 )



c
c
c
c
P c =  J1 (ζ2 )H10
Y1 (ζ2 )H11
J1 (ζ2 )H11
Y1 (ζ2 )H10


.
.
.
.

.
.
.
.

.
.
.
.




0
0
0
0


c
c
c
c
JN (ζ2 )HN 0 YN (ζ2 )HN 0 JN (ζ2 )HN 1 YN (ζ2 )HN 1

···

0

···

1
c
2 J0 (ζ2 )H0N

1
c
2 Y0 (ζ2 )H0N

···

0

0

···

c
JN (ζ2 )H1N

c
YN (ζ2 )H1N

..

.
.
.

.
.
.

JN (ζ1 )

YN (ζ1 )

.

···

c
· · · JN (ζ2 )HN N

T

Xc =
4 See

β20 γ20 β21 γ21 · · · β2N

appendix B

93

0

γ2N

c
YN (ζ2 )HN N

























c
and where ζ1 = ka, ζ2 = kR, and Hnm = Jn−m (kd) + (−1)m Jn+m (kd) . Equation
c
c
(5.18) is used to solve for the eigenfrequencies ki and eigenvectors Xi for i = 1, 2, · · · . Using
Eq. (5.10), the symmetric modes can be subsequently obtained as follows:

ϕc (r, θ) = β20,i
i

N

γ2m,i
β2m,i
c
c
Jm (ki r) +
Ym (ki r) cos mθ
β20,i
β20,i
m=0

(5.19)

where the terms (β2m,i /β20,i ), (γ2m,i /β20,i), m = 0, 1, · · · , N are entries of the i-th eigenc
vector Xi , and β20,i is a scaling factor for the i-th symmetric mode.

5.2.4

Odd or Antisymmetric Modes

The antisymmetric modes are associated with the sine term sin mθ in Eq. (5.10) and are
denoted by ϕs (r, θ); the superscript “s” is used to represent the sine terms. For the symmetric
modes, the boundary conditions in Eqs. (5.11) and (5.12) give
∞
{β1m Jm (ka) + γ1m Ym (ka)} sin mθ = 0
∞

m=1
β1m Jm (kr|ρ=R ) + γ1m Ym (kr|ρ=R ) sin mθ = 0

(5.20)

m=1
Using Graf’s addition formula [80], the boundary conditions in Eq. (5.20) can be rewritten
as
∞
{β1n Jn (ka) + γ1n Yn (ka)} sin nθ = 0
∞

∞

n=1
{β1m Jn (kR) + γ1m Yn (kR)} Jn−m (kd) − (−1)m Jn+m (kd) sin nψ = 0

m=1 n=1
(5.21)
94

Expanding Eq. (5.21) from m = 1 to M and n = 1 to M, M = (N + 1), gives the matrix
equation
P sX s = 0

(5.22)

where P s ∈ R2M×2M and B ∈ R2M are given by the relations



J1 (ζ1 )
Y1 (ζ1 )
0
0



s
s
s
 J1 (ζ2 )H s
Y1 (ζ2 )H12
J1 (ζ2 )H12
Y1 (ζ2 )H11
11




0
0
J2 (ζ1 )
Y2 (ζ1 )



s
s
s
s
P s =  J2 (ζ2 )H21
Y2 (ζ2 )H22
J2 (ζ2 )H22
Y2 (ζ2 )H21



.
.
.
.
.
.
.
.

.
.
.
.




0
0
0
0


s
s
s
s
JM (ζ2 )HM 1 YM (ζ2 )HM 1 JM (ζ2 )HM 2 YM (ζ2 )HM 2

···

0

0

···

s
J1 (ζ2 )H1M

s
Y1 (ζ2 )H1M

···

0

0

···

s
J2 (ζ2 )H2M

s
Y2 (ζ2 )H2M

..

.
.
.

.
.
.

JM (ζ1 )

YM (ζ1 )

.

···

s
· · · JM (ζ2 )HM M

s
YM (ζ2 )HM M

























T

Xs =

β11 γ11 β12 γ12 · · · β1M

γ1M

s
and where ζ1 = ka, ζ2 = kR, and Hnm = Jn−m (kd) − (−1)m Jn+m (kd) . Equation
s
s
(5.22) is used to solve for the eigenfrequencies ki and eigenvectors Xi for i = 1, 2, · · · . Using
Eq. (5.10), the symmetric modes can be subsequently obtained as follows:

ϕs (r, θ) = β11,i
i

M

γ1m,i
β1m,i
s
s
Jm (ki r) +
Ym (ki r) sin mθ
β11,i
β11,i
m=1

(5.23)

i
i
i
i
where the terms (β1m /β11 ), (γ1m /β11 ), m = 1, 2, · · · , M are entries of the i-th eignevector
s
i
Xi , and β11 is a scaling factor for the i-th antisymmetric mode.
From the expressions in Eqs. (5.19) and (5.23) it can be seen that each mode shape of a
95

constrained membrane ϕi is a summation of terms that are a function of m, i.e., each mode
shape depends on all the angular nodes. This is different from unconstrained membranes
where each mode shape is associated with a specific angular node only and is denoted by
ϕi,m , i.e., a pair of subscripts are used to identify the mode shape. Since each mode shape
is associated with an eigenfrequency, the eigenfrequencies of unconstrained membranes are
denoted by ki,m whereas the eigenfrequencies of constrained membranes are denoted by ki .
In the literature, the subscripts for frequencies and mode shapes of unconstrained membranes are numbered starting with zero, i.e., i, m = 0, 1, 2, · · · . In a slight deviation from this
tradition, the frequencies and mode shapes of constrained membranes have been numbered
starting with one in this chapter, i.e., i = 1, 2, 3, · · · .

5.3

Orthogonality of Modes

The differential equation in Eq. (5.2) is related to the eigenfrequency k as follows
1 1 ∂2F
1 ∂ 2 ϕ 1 ∂ϕ
1 ∂2ϕ
+
}=
= −k 2
{
+
2
2 ∂θ2
2 F ∂t2
ϕ ∂r
r ∂r
r
c

and can be rewritten in the form
∂ϕ
1 ∂2ϕ
1 ∂
+ k2 ϕ = 0
r
+
r ∂r ∂r
r 2 ∂θ2

(5.24)

Let ki and kj be the i-th and j-th eigenfrequency associated with the modes ϕi and ϕj ,
respectively. Each of the modes ϕi and ϕj may be symmetric or anti-symmetric. For these
96

two modes, Eq. (5.24) can be written as follows
∂ϕ
1 ∂
1 ∂ 2 ϕi
2
+ k i ϕi = 0
r i +
r ∂r
∂r
r 2 ∂θ2
2
∂ϕj
1 ∂ ϕj
1 ∂
2
r
+
+ k j ϕj = 0
2 ∂θ2
r ∂r
∂r
r

(5.25)
(5.26)

Multiplying Eq. (5.25) by ϕj and Eq. (5.26) by ϕi and subtracting one from the other, we
get
∂ϕj
∂ϕ
∂
∂
r i ϕj −
r
ϕi
∂r
∂r
∂r
∂r

1
r

∂ 2 ϕj
∂ 2 ϕi
ϕj −
ϕi
∂θ2
∂θ2

1
+
r2

2
2
+ k i − k j ϕi ϕj = 0
(5.27)

Integrating Eq. (5.27) over the active region of membrane now gives
2π
0
+
+

2π
0

r|ρ=R

∂ϕj
∂ϕ
∂
∂
r i ϕj −
r
ϕi
∂r
∂r
∂r
∂r

a
r|ρ=R 1
r
a

2
2
ki − kj

2π
0

∂ 2 ϕj
∂ 2 ϕi
ϕj −
ϕi
∂θ2
∂θ2
r|ρ=R
a

drdθ

(5.28)

drdθ

ϕi ϕj rdrdθ = 0

Integration of the first term by parts gives
2π

2π
0

r

∂ϕj
∂ϕi
ϕj − r
ϕ
∂r
∂r i

r|ρ=R

0
a
r|ρ=R
dθ −
r=a

∂ϕj
∂ϕ
∂
∂
r i ϕj −
r
ϕi
∂r
∂r
∂r
∂r
2π
0

drdθ =

r|ρ=R ∂ϕ ∂ϕj
∂ϕj ∂ϕ
i rdrdθ
i
−
∂r ∂r
∂r ∂r
a

Since both ϕi and ϕj are zero for r = a and r values where ρ = R, the first term on the
right-hand side of the above equation is zero. The second term on the right-hand side of the
97

above equation is zero trivially and this implies
2π
0

r|ρ=R 1
r
a

∂ 2 ϕj
∂ 2 ϕi
ϕj −
ϕi
∂θ2
∂θ2

2
2
drdθ + ki − kj

2π
0

r|ρ=R
a

ϕi ϕj rdrdθ = 0
(5.29)

The first term in Eq. (5.29) is given by the expression
2π
0

r|ρ=R 1
r
a

∂ 2 ϕj
∂ 2 ϕi
ϕj −
ϕi
∂θ2
∂θ2

drdθ

(5.30)

where the limits of integration correspond to the unconstrained area of the membrane in
Fig. 5.1. In Eq. (5.30), the integration is first carried out with respect to r and then with
respect to θ. For the first integration, the upper limit of r is expressed as a function of θ;
this functional dependence of r on θ is given by Eq. (5.14) and shown in Fig. 5.2 (a). We
evaluate the expression in Eq. (5.30) by changing the order of integration. This requires θ
to be expressed as a function of r, which is obtained from the geometry in Fig. 5.1:

θ = θ(r) = cos−1

R 2 − r 2 − d2
,
2rd

θ ∈ [0, π]

(5.31)

The area of integration is redrawn in Fig. 5.2 (b) and the integral in Eq. (5.30) is expressed
as the sum of the following three terms:
∂ 2 ϕj
∂ 2 ϕi
ϕj −
ϕi dθdr
∂θ2
∂θ2
a
R+d −θ(r) 1 ∂ 2 ϕ
∂ 2 ϕj
iϕ −
ϕi dθdr
r ∂θ2 j
∂θ2
R−d −π
R+d +π 1 ∂ 2 ϕ
∂ 2 ϕj
iϕ −
ϕi dθdr
j
∂θ2
R−d +θ(r) r ∂θ2
R−d

+
+

+π 1
−π r

98

(5.32)

where θ(r) is defined by Eq. (5.31). For each of the three terms in Eq. (5.32), the inner
r

θ
+π

r(θ)

+θ(r)
0
−θ(r)
a
0

2π

−π

θ

a

R+d

R−d

r

(b)

(a)

Figure 5.2: Shaded region showing the area of integration for the two cases where the
integration is done (a) first with respect to r and then with respect to θ, and (b) first with
respect to θ and then with respect to r
integration with respect to θ by parts, gives
∂ 2 ϕj
∂ 2 ϕi
ϕj −
ϕi dθ =
∂θ2
θ ∂θ2
=

θu
∂ϕj
∂ϕj ∂ϕ
∂ϕi
∂ϕi ∂ϕj
i dθ
ϕj −
ϕi
−
−
∂θ
∂θ
∂θ ∂θ
∂θ ∂θ
θ
θl
θu
∂ϕj
∂ϕi
ϕ −
ϕ
(5.33)
∂θ j
∂θ i
θl

where θl and θu denote the lower and upper limits of θ. Using Eq. (5.33), the three terms
in Eq. (5.32) can be simplified to the form
+π
R−d 1 ∂ϕ
R+d
∂ϕj
iϕ −
ϕi
dr +
r ∂θ j
∂θ
a
R−d
−π
R+d
+
R−d

−θ(r)
∂ϕj
1 ∂ϕi
ϕ −
ϕ
dr
r ∂θ j
∂θ i
−π
+π
∂ϕj
1 ∂ϕi
ϕ −
ϕ
dr
r ∂θ j
∂θ i
+θ(r)

Since θ = +θ(r) and θ = −θ(r) describe the outer boundary of the membrane where ϕi =
99

ϕj = 0 for all i and j, the above integrals simplify to the form

+

=

+π
R−d 1 ∂ϕ
∂ϕj
iϕ −
ϕ
dr
r ∂θ j
∂θ i
a
−π
+π
R+d 1 ∂ϕ
∂ϕj
iϕ −
ϕ
dr
j
∂θ i
R−d r ∂θ
−π
+π
R+d 1 ∂ϕ
∂ϕj
iϕ −
ϕ
dr
r ∂θ j
∂θ i
a
−π

(5.34)

Now consider the three separate cases:
A1. Modes ϕi and ϕj are both symmetric: Using Eq. (5.19), both (∂ϕi /∂θ) and (∂ϕj /∂θ)
can be expressed in the form
∂ϕi
=
∂θ

∞
(.) sin mθ,
m=0

∂ϕj
∂θ

∞
=

(.) sin mθ
m=0

Since their values are zero at θ = ±π, the integral in Eq. (5.34) is identically zero.
A2. Modes ϕi and ϕj are both anti-symmetric: It can be see from Eq. (5.23) that ϕi and
ϕj have the form
∞

∞
(.) sin mθ,

ϕi =

(.) sin mθ

ϕj =
m=0

m=0

Since their values are zero at θ = ±π, the integral in Eq. (5.34) is identically zero.
A3. Mode ϕi is symmetric and mode ϕj is anti-symmetric: Using Eqs. (5.19) and (5.23) it
can be shown that both the terms (∂ϕi /∂θ)ϕj and (∂ϕj /∂θ)ϕi are zero at θ = ±π,
and therefore the integral in Eq. (5.34) is identically zero.
100

This proves that the integral in Eq. (5.30) is zero for any two distinct modes ϕi and ϕj ,
irrespective of whether they are both symmetric, both antisymmetric, or one of them is
symmetric and the other is antisymmetric. From Eq. (5.29) we can now conclude
2π
0

r|ρ=R
a

ϕi ϕj rdrdθ = 0,

i=j

(5.35)

which is the orthogonality condition for a pair of distinct modes.

5.4

Accurate Computation of Eigenfrequencies and Mode
Shapes

5.4.1

Computational Algorithm

We present an algorithm for numerical computation of the eigenfrequencies and their corresponding mode shapes. The algorithm is identical for symmetric and antisymmetric modes
and therefore we discuss it here for the symmetric modes only. The algorithm first computes
c
all eigenfrequencies ki , i = 1, 2, · · · , that are less than a user-specified maximum eigenfrec
c
quency kmax . The interval [0, kmax ] is divided into small segments of length ∆k and the
eigenfrequencies within each segment are computed separately. The number of eigenfrequencies that show up within each segment depend on the value of N, which is related to the
size of the matrix P c in Eq. (5.18). Starting from zero, increasing the value of N typically
results in larger number of zero crossings of the determinant P c (see Fig. 5.3), and hence a
larger number of eigenfrequencies within a specific segment. Beyond a certain value of N,
the number of eigenfrequencies in a given segment will not change. At this point, increasing
101

|P |

large N
small N
ki

k

δ0

∆k

kmax
δf

Figure 5.3: The eigenfrequencies are the zero crossings of the determinant of the matrix P c
for symmetric modes, and P s for antisymmetric modes

the value of N produces very small changes in the values of the eigenfrequencies but such
changes can result in significant improvement or loss of accuracy of the associated mode
shapes. The flowchart shown in Fig. 5.8 provides the algorithm for computing all eigenc
frequencies in the interval [0, kmax ], and accurately computing the mode shape associated
with each eigenfrequency through proper choice of N. A brief summary of the algorithm is
provided next.
c
The interval [0, kmax ] is divided into small segments of length ∆k. For any given segment,
the boundaries of the segment are denoted by δ0 and δf . The value of k is changed in
increments of δ (small number) from δ0 to δf , and zero crossings of | P c | are stored as
potential eigenfrequencies. For a given segment, this procedure is repeated by increasing
N, N ∈ {0, 1, 2, · · · , Nmax }, till no additional zero crossings occur and the mode shapes
associated with the eigenfrequencies have been computed accurately. The mode shapes are
computed using Eq. (5.19) and are checked for accuracy by verifying ϕc = 0 along the inner
i
and outer boundary of the membrane. The procedure is repeated for each segment till the
102

maximum user-specified value of k = kmax is reached.
In the literature, several methods have been proposed for computing the eigenfrequency
and mode shapes of constrained membranes. Most of these methods are well suited to computing the first few eigenfrequencies and mode shapes only. For simulation of the dynamics
of a constrained membrane, it is necessary to accurately compute many frequencies and
mode shapes; this is difficult to compute using existing methods due to singularity. For
example, Nagaya [2], Singh and Kothari [47], and Singh, et al. [3] proposed to combine the
two equations for symmetric modes in Eq. (5.13) into the following single equation
∞
β2m

Jm (kr|ρ=R ) −

m=0

Jm (ka)
Ym (kr|ρ=R ) cos mθ = 0
Ym (ka)

(5.36)

This procedure reduces the dimension of the coefficient matrix by a factor of two and reduces
the computational burden but presence of Bessel function of the second kind, Ym (ka), in the
denominator makes it prone to singularity. This problem become critical during computation
of higher frequencies.
The eigenfrequencies and mode shapes of the antisymmetric modes can be computed
c
s
using the algorithm in Fig. 5.4 by replacing P c , ki , and ϕc with P s , ki , and ϕs , respectively,
i
i
and N with M. P s is defined in Eq. (5.22) and ϕs is defined in Eq. (5.23).
i

5.4.2

Numerical Results

We use the algorithm in Section 5.4.1 to accurately compute the eigenfrequencies and mode
shapes of a constrained membrane with R/a = 0.25 and d/a = 1.0. This problem was
considered in the literature [1], [2], [47] and the results were summarized by Singh, et al.
[3]. We compare these results with those obtained using our algorithm, implemented with
103

Input membrane data
R, a, d, T , ρ

i=1
δ0 = 0, δf = ∆k

Choose
Nmax , kmax , ∆k, δ

N =0

c
ki = δ0

Calculate
P c in Eq. (5.18)

δ0 = δf
δf = δf + ∆k

N =N +1

c
c
ki = ki + δ

No

δf ≥ kmax ?

No
Yes

N = Nmax ?

No

| P c |≈ 0 ?

No
Yes

c
ki = kf ?

Yes

No

Stop

Yes
Compute ϕc using
i
Eq. (5.19) along inner
and outer boundaries

ϕc ≈ 0 ?
i

i=i+1
Yes
Save
c
ki , ϕc , Ni = N
i

Yes

c
repeated ki ?

No
Figure 5.4: Algorithm for computing eigenfrequencies and mode shapes for symmetric modes
kmax = 10a and ∆k = kmax . Table 5.1 shows the non-dimensional fundamental antisymmetric frequencies for the constrained membrane obtained using our algorithm; the value of
104

Table 5.1: First eight non-dimensional frequencies of antisymmetric modes of constrained
membrane with a/R = 0.25 and d/a = 1.0. The results in parenthesis are from [1], [2],
and [3], and the star superscript indicates frequencies with accurate mode shapes computed
using our algorithm
M

c
k1 a

c
k2 a

c
k3 a

c
k4 a

c
k5 a

c
k6 a

c
k7 a

c
k8 a

1
2
3
4
5
6
7
8
9
10
11
12
13
40
-

1.1119
1.0636
1.0661
1.0660
1.0660*
1.0660*
1.0660*
1.0660*
1.0660*
1.0660*
1.0660*
1.0660*
1.0660*
1.0647
(1.0659)
(1.0500)

2.1342
1.4307
1.3811
1.3850
1.3848
1.3848
1.3848*
1.3848*
1.3848*
1.3848*
1.3848*
1.3848*
1.3848*
1.3822
(1.3848)
(1.3800)

1.9199
1.6875
1.6557
1.6591
1.6589
1.6589
1.6589*
1.6589*
1.6589*
1.6589*
1.6589*
1.6589*
1.6572
(1.6589)
-

1.9290
1.9288
1.9275
1.9284
1.9284
1.9284
1.9284*
1.9284*
1.9284*
1.9284*
1.9284*
1.9272
(1.9283)
-

2.3177
1.9473
1.9296
1.9312
1.9311
1.9311
1.9311*
1.9311*
1.9311*
1.9311*
1.9311*
1.9297
-

2.3280
2.2195
2.2083
2.2101
2.2099
2.2099
2.2099
2.2099*
2.2099*
2.2099*
2.2072
-

2.3292
2.3288
2.3288
2.3288
2.3288
2.3288*
2.3288*
2.3288*
2.3288*
2.3272
-

2.4982
2.4916
2.4928
2.4927
2.4927
2.4927
2.4927*
2.4927*
2.4922
-

M is gradually increased till all modes are accurately computed. As mentioned in Section
5.4.1, the accuracy of each mode shape is verified by computing the displacement of 100
points along the inner and outer membrane boundaries. For accuracy, the displacement of
the membrane at all these points must satisfy ǫ

10−4 . In Table 5.1, eigenfrequencies

associated with accurate mode shapes (this depends on correct choice of M) have a star superscript. The eigenfrequencies obtained from the literature are presented in separate rows
at the bottom of the table and are placed within parenthesis. The following observations
can be made from the results presented in Table 5.1:

A1. It is not possible to identify all the eigenfrequencies in an interval using a small value
of M. Increasing the value of M helps in identifying more eigenfrequencies and M ≥ 6
105

identifies all eight eigenfrequencies in this particular example.
A2. For M = 6, all eigenfrequencies are identified but only the first eigenfrequency has an
accurate mode shape. Increasing the value of M results in a larger number of accurate
mode shapes and M = 12 results in accurate computation of all mode shapes.
A3. A large value of M is not guaranteed to result in accurate mode shapes. For example,
M = 13 gives all accurate mode shapes but none of the mode shapes are accurate for
M = 40. This is because of the sensitivity of the Bessel functions. This justifies the
need for increasing M gradually rather than choosing a pre-specified large integer.
A4. Although all the mode shapes can be accurately computed using M = 12, the computational burden can be significantly reduced by using the lowest value of M for each
eigenfrequency that gives an accurate mode shape. For example, the mode shape for
c
c
k1 a should be computed using M = 5, the mode shape for k2 a should be computed
using M = 7, and so on.
A5. The eigenfrequencies in the literature have values that are very close to those obtained
by our algorithm. However, these small errors result in large errors in the shape of
the modes, which will violate boundary and orthogonality conditions. This will be
discussed further later.
We next present results for different sizes and locations of the constraint in the membrane.
The specific values of a/R and d/a used are taken from Singh, et al. [3]. Table 5.2 presents
the first five non-dimensional frequencies of symmetric modes of the constrained membranes
and values of N used to compute the corresponding mode shapes. For a given set of a/R
and d/a values, the first row contains eigenfrequencies and N values that correspond to
106

accurate mode shapes. The procedure for obtaining accurate eigenfrequencies and mode
shapes was discussed in Section 5.4.1. The second row contains eigenfrequencies that are

Table 5.2: First five non-dimensional frequencies of symmetric modes of constrained membranes and values of N used to compute the mode shapes. For a given set of a/R and d/a
values, the first row contains values that correspond to accurate mode shapes. The second
row contain values that are closest to those presented in the literature - the eigenfrequencies
are underlined to indicate that the corresponding mode shapes are inaccurate. The results
from the literature [1], [2], [3] are presented in the following rows with eigenfrequencies in
parenthesis
a/R

d/a
0.2

0.4
0.5
0.6

0.8

1.0

1

0.25

2

3

c
k1 a

N

c
k2 a

N

c
k3 a

N

c
k4 a

N

c
k5 a

N

2.7348
2.7348
(2.7348)
(2.7520)
2.4053
2.4053
(2.4053)
(2.4300)
2.1536
2.2081
(2.2800)
1.9570
1.9537
(1.9340)
1.8001
1.8001
(1.8000)
0.8681
0.8680
(0.8680)
(0.8720)
0.7456
0.7456
(0.7455)
(0.7430)
0.6674
0.6690
(0.6700)

5
3
5
4
6
1
6
2
7
4
5
2
8
3
9
3
-

3.2374
3.2373
(3.2373)
(3.2500)
3.0862
3.0862
(3.0862)
(3.0860)
2.9451
2.9429
(2.9310)
2.8256
2.8109
(2.8130)
2.7236
2.7239
(2.7240)
1.2305
1.2305
(1.2305)
(1.2200)
1.1967
1.1966
(1.1966)
(1.2040)
1.1024
1.1129
(1.1230)

6
4
7
5
8
4
8
3
9
5
6
4
10
5
11
3
-

3.6531
3.6531
(3.6530)
(3.6680)
3.6973
3.6973
(3.6972)
(3.6950)
3.6598
3.6505
(3.6620)
3.5975
3.5779
(3.5620)
3.4187
3.4183
(3.4180)
1.4817
1.4817
(1.4817)
(1.4400)
1.4521
1.4521
(1.4521)
(1.4620)
1.3636
1.3660
(1.3670)

7
5
9
6
9
5
11
5
9
5
10
5
8
6
12
4
-

3.9346
3.9346
(3.9346)
(3.9370)
4.2487
4.2487
(4.2487)
(4.2500)
4.1455

8
6
10
8
8

4.2482
4.2482
(4.2482)

8
6
-

4.6704
4.6704
(4.6704)

9
5
-

3.7439

11

11
3
12
6

3.5358

12

1.6789
1.6789
(1.6788)

8
6
-

4.3194
4.4738
(4.4320)
4.3140
4.4424
(4.4250)
4.2690
4.3044
(4.3750)
1.6882
1.6882
(1.6881)

1.5610
1.5610
(1.5609)

11
7
-

1.8897
1.8898
(1.8898)

12
8
-

1.5074

13

1.7571

12

107

13
7
8
6
-

closest to values presented in the literature; these eigenfrequencies are obtained by gradually
increasing N till the error between the computed value and the value in the literature is
minimum. In all cases, the error is minimum for a value of N for which the mode shape
is inaccurate - this is indicated by the underlined eigenfrequencies. The eigenfrequencies in
parenthesis are taken from the literature [3] and are presented in the third and fourth rows.
The following observations can be made from the results presented in Table 5.2 pertaining
to accurate eigenfrequencies and mode shapes:
A1. Lower values of N are required to compute lower eigenfrequencies, and vice versa. This
can also be observed from the data presented in Table 5.1.
A2. Lower values of N are required for lower values of d, i.e., constraints located farther
from the center of the membrane, and vice versa.
A3. Higher values of N are required for smaller values of a, i.e., constraints of smaller size,
and vice versa.
The effect of N on the accuracy of eigenfrequencies was first observed by Nagaya [2], who
used 8 ≤ N ≤ 12 to confine the error to ≈ 0.5%. Such choice of N may indeed result in
accurate computation of the eigenfrequency but there is no guarantee that the mode shape
will be accurate. The accuracy of the mode shape depends on the size and location of the
constraint, as discussed above, and incorrect choice of N can result in significant error in
the mode shape even when the eigenfrequency is computed accurately. This is illustrated
c
with the help of Fig. 5.5, which shows the accurate and inaccurate mode shapes for k3 a
for a/R = 0.5 and d/a = 0.6, 0.8, and 1.0. It should be noted that the eigenfrequencies
corresponding to these accurate and inaccurate mode shapes, shown in Table 5.2, indicate
negligible error. For example, for a/R = 0.5 and d/a = 0.8, the error in the eigenfrequency
108

is only
3.5975 − 3.5620
× 100 ≈ 1%
3.5975
but the mode shapes in Fig. 5.5 (b) are significantly different; and of course, the inaccurate
mode shape does not satisfy the boundary conditions. Similar observations can be made for
the cases a/R = 0.5, d/a = 0.6 and a/R = 0.5, d/a = 1.0. For these cases, the error in the
eigenfrequency is ≈ −0.06% and ≈ 0.02%, respectively but the corresponding mode shapes
are significantly different. We complete this section with numerical results for an example
where the size of the constraint is small and the location of the constraint is far away from
the center of the membrane, i.e., a/R is small and d/a is large. For such cases, it is difficult
to obtain the eigenfrequencies and mode shapes accurately because the Bessel functions are
sensitive to large values of d and small values of a. Despite this sensitivity, our algorithm
was able to compute the first twenty modes of a membrane with a/R = 0.1 and d/a = 9.0.
The first symmetric mode and the sixteenth symmetric mode of the membrane are shown
in Fig. 5.6 together with the non-dimensional eigenfrequencies and the required value of N.
Note that N = 18 was required for computing the sixteenth mode accurately. Figure 5.7
shows more mode shapes with d/a = 2.0.

109

c
k3 a = 3.6505, N = 5

c
k3 a = 3.6598, N = 9

c
k3 a = 3.5779, N = 5

c
k3 a = 3.4183, N = 5

c
k3 a = 3.5975, N = 11

c
k3 a = 3.4187, N = 9

(a)

(b)

(c)

c
Figure 5.5: Mode shapes k3 a for constrained membranes in Table 5.2 with a/R = 0.5 and (a) d/a = 0.6 (b) d/a = 0.8 (c)
d/a = 1.0. For each of these three cases, the inaccurate mode shapes are shown above the accurate mode shapes

110

c
k1 a = 0.24927, N = 6

c
k16 a = 1.17795, N = 18

(a)

(b)

Figure 5.6: First symmetric mode and sixteenth symmetric mode of a constrained membrane with a/R = 0.1 and d/a = 9.0.
For both these cases, the isometric view is placed above the top view and the value of N required for accuracy is shown

111

(1)

(1)

(2)

(2)

(3)

(3)

(4)

(4)

(5)

(5)

(a)

(b)

Figure 5.7: First five modes of a membrane with a/R = 0.1 and d/a = 2.0 (a) Antisymmetric
modes (b) Symmetric modes. For both, the isometric view is placed next the top view

112

5.5
5.5.1

Simulation of Constrained Membrane Dynamics
Computation of Modal Coefficients

From the general form of the dynamics of the membrane given in Eq. (5.9), the complete
solution can be written as follows
∞
Zc (r, θ, t) =

c
c
c
ϕc (ki r, θ)(A1i cos ωi t + B1i sin ωi t)
i

i=0
∞
s
s
s
ϕs (kj r, θ)(A2j cos ωi t + B2j sin ωj t)
+
j
j=1

(5.37)

c
c s
s
where ϕc , ϕs are defined in Eqs. (5.19) and (5.23), respectively, and ωi = cki , ωj = ckj ,
i i
i, j = 1, 2, · · · . The initial conditions are assumed to be

˙
Z(r, θ, 0) = H(r, θ)

Z(r, θ, 0) = G(r, θ),

Applying the displacement initial conditions to Eq. (5.37), we get
∞
i=1

Ac ϕc +
i i

∞
j=1

As ϕs = G(r, θ)
j j

(5.38)

Multiplying Eq. (5.38) by ϕc , integrating both sides over the unconstrained area of the
i
membrane, and using the orthogonality condition in Eq. (5.35), we get

Ac =
i

2π
0

r|ρ=R
a

G(r, θ) ϕc rdrdθ /
i

2π
0

r|ρ=R
a

(ϕc )2 r drdθ ,
i

i = 1, 2, · · ·
(5.39)

113

Multiplying Eq. (5.38) by ϕs , integrating both sides over the unconstrained area of the
i
membrane, and using the orthogonality condition in Eq. (5.35), we get

As =
j

2π

r|ρ=R

0

a

2π

G(r, θ) ϕs rdrdθ /
i

r|ρ=R

0

a

(ϕs )2 r drdθ ,
i

j = 1, 2, · · ·
(5.40)

Applying the velocity initial conditions to Eq. (5.37), we get
∞
i=1

c c
Bi ωi ϕc +
i

∞
j=1

s s
Bj ωj ϕs = H(r, θ)
j

(5.41)

By repeating the same procedure as above, the remaining modal coefficients can be obtained
as

1
c
Bi = c
ωi

2π
0

r|ρ=R
a

2π

H(r, θ) ϕc rdrdθ /
i

0

r|ρ=R
a

(ϕc )2 r drdθ ,
i

i = 1, 2, · · ·
(5.42)

1
s
Bj = s
ωi

2π
0

r|ρ=R
a

2π

H(r, θ) ϕs rdrdθ /
i

0

r|ρ=R
a

(ϕs )2 r drdθ ,
i

j = 1, 2, · · ·
(5.43)

The integrals in Eqs. (5.39), (5.40), (5.42) and (5.43) are of the form

2π
0
where g(θ) = r|ρ=R =

g(θ)
a

f (r, θ)drdθ

R2 − d2 sin2 θ − d cos θ. Through change of variables x = θ and
114

r = [(1 − y)a + yg(x)] the integrals are changed to the form
2π
0

1
0

f [(1 − y)a + yg(x), x] [g(x) − a]dydx =

2π
0

1
0

h(y, x) dy dx

The integrations are performed numerically using a n-point Gaussian quadrature rule as
follows

2π
0

n n
1
π
π 1
wi wj h(ζi , ηj )
h(ζ, η) dη dζ =
h(y, x) dy dx =
2 −1 −1
2
0
i=1 j=1
1

(5.44)

where x = π(η + 1), y = (ζ + 1)/2, wi and wj are the weights, and ζi and ηj are the integration points respectively. Since Bessel functions are very sensitive to numerical integration,
accurate computation of the modal coefficients requires correct choice of the number of integration points in Eq. (5.44). An algorithm for properly choosing n and accurately computing
the modal coefficients is described by the flow chart in Fig. 5.8, which is self-explanatory.

115

Stack all modes: symmetric and antisymmetric
ϕi , i = 1, 2, · · · , ℓ

Stop

i=1

Yes

No

Choose ϕi

i=i+1

i=ℓ?

Select n

Check for orthogonality
Iij =
ϕi ϕj rdrdθ

Compute modal coefficients
using Eqs. (5.40) and (5.43)
No

Increase n
Iij = 0 ?
No

∀j = i

Yes

Mode ϕi
symmetric ?

Yes

Compute modal coefficients
using Eqs. (5.39) and (5.42)

Figure 5.8: Algorithm for computing modal coefficients
116

5.5.2

Numerical Results Verifying Orthogonality

In Section 5.4.2 we showed that results in the literature provide reasonably accurate eigenfrequencies but the associated mode shapes are inaccurate and as such do not satisfy the
boundary conditions. In this section we show that these inaccurate modes do not satisfy the
orthogonality condition; the accurate modes satisfy the orthogonality conditions but requires
proper choice of the number of integration points n. The number of integration points can
be chosen properly by following the algorithm described by the flowchart in Fig. 5.8. We
consider constrained membranes with a/R = 0.5 and d/a = 0.6, 0.8, 1.0; these cases were
presented in Table 5.2. For this membrane, the third symmetric mode ϕc and fourth sym3
metric mode ϕc are used to check for orthogonality. In particular, we check the orthogonality
4
of ϕc (accurate mode shape) with both ϕc (accurate mode shape) and ϕc (inaccurate mode
4
3
3
shape). The results are presented in Table 5.3.
Table 5.3: Checking orthogonality between the accurate and inaccurate third symmetric
mode and the accurate fourth symmetric mode of constrained membranes with a/R = 0.5
and d/a = 0.6, 0.8, 1.0. The inaccurate third symmetric mode is underlined
d/a
n
4
8
12
16
20
24
28
32

0.6
c , ϕc > < ϕc , ϕc >
< ϕ3 4
3 4
0.0091
0.0149
-0.0366
-0.0821
0.0142
-0.0080
0.0056
0.0081
-0.0013
0.0011
-0.0015
0.0000
-0.0016
0.0000
-0.0016
0.0000

0.8
c , ϕc > < ϕc , ϕc >
< ϕ3 4
3 4
-0.0221
-0.0465
-0.0621
-0.1769
-0.0069
-0.0514
0.0187
0.0427
-0.0080
0.0184
-0.0061
0.0025
-0.0068
0.0002
-0.0067
0.0000

1.0
c , ϕc > < ϕc , ϕc >
< ϕ3 4
3 4
0.0560
-0.3005
-0.0171
0.1715
0.0107
0.0327
0.0026
0.0222
0.0030
0.0244
0.0019
0.0058
0.0021
0.0006
0.0022
0.0000

To show the orthogonality between symmetric-symmetric and symmetric-antisymmetric
modes, we consider a membrane with a/R = 0.1 and d/a = 0, 4, 8. The results are shown in
Table 5.4. It is clear that accurate computation of the mode shapes require proper choice of
117

the number of integration points n. The number of integration points depends on the size
and location of the constraint and can be chosen using the flowchart in Fig. 5.8. It may be
mentioned that symmetric-antisymmetric modes are orthogonal for small values of n and
higher values of n are required to satisfy orthogonality between symmetric-symmetric and
antisymmetric-antisymmetric modes.
Table 5.4: Number of integration points n required to achieve orthogonality between
symmetric-symmetric and symmetric-antisymmetric modes of constrained membranes with
a/R = 0.1 and d/a = 0, 4, 8
d/a
n
4
8
12
16
20
24
28
32
36

5.5.3

0
ϕc , ϕc
4 5
-0.0013
0.0013
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000

4
ϕc , ϕs
4 5
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000

ϕc , ϕc
4 5
-0.3009
0.2686
-0.1659
-0.0741
-0.0086
-0.0004
0.0000
0.0000
0.0000

8
ϕc , ϕs
4 5
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000

ϕc , ϕc
4 5
-0.0746
0.0215
0.0293
0.0182
0.0224
0.0100
0.0020
0.0002
0.0000

ϕc , ϕs
4 5
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000

Simulation of Free Vibration

In this section, we simulate the motion of a constrained membrane subjected to an initial
condition for different cases. The membrane is assumed to have the following specifications
R = 1 m,
µ = 0.25 kg.m−2 ,
T = 1 N.m−1

118

• Case 1 : a/R = 0.12 and d/R = 0.88 and the initial conditions are represented as
G(r, θ) = 0.1(r − a)8 (r − r(θ)),

H(r, θ) = 0,

where r(θ) =

r≥a

R2 − d2 sin2 θ − d cos θ
z

z

y
x

x

(a)

(b)

Figure 5.9: Initial displacement of the membrane described in case 1
z

z

y
x

x

(a)

(b)

Figure 5.10: Initial displacement of the membrane described in case 2

• Case 2 : a/R = 0.12 and d/R = 0.66 and the initial conditions are represented as

H(r, θ) = 0,

3
G(r, θ) = (r − a)(r − r(θ))8 ,
5

r≥a

• Case 3 : In this case, the membrane is assumed to be unconstrained and vibrating in
119

its first mode. As it passes through the mean position, the constraint is applied and
held fixed for all future time with a/R = 0.1 and d/R = 0.4. The position and velocity
of the unconstrained membrane are given by the expressions
c
c
Z(r, θ, t) = J0 (k1 ρ) cos ω1 t,

c
c
c
˙
Z(r, θ, t) = −ω1 J0 (k1 ρ) sin ω1 t

c
r 2 + d2 + 2rd cos θ is defined in Fig. 5.1, k1 R = 2.405 and is obtained
c
c
from the first zero of the Bessel function J0 , ω1
T /ρ k1 , and t = 0 denotes the
c
initial time when the membrane is released from rest. At t = π/2ω1 , the membrane
where ρ

passes through the mean position and the constraint is applied. The initial conditions
for the constrained membrane are therefore

G(r, θ) = 0,



 −ω c J (k c ρ) : a < r ≤ R2 − d2 sin2 θ − d cos θ
1 0 1
H(r, θ) =


0 : 0≤r≤a

The dynamics of the membrane was simulated using 16 symmetric and 16 antisymmetric
modes. Fig. 5.11 shows the evolution of the motion for the membrane in case 1, where the
isometric and side views are showed for each snapshot and the beginning of the motion is
showed in Fig. 5.11 (a) and continues in Fig. 5.11 (b). The simulation of membrane in case
2 is shown in Fig. 5.12, where the motion starts in Fig. 5.12 (a) and continues in Fig. 5.12
(b) and Fig. 5.12 (c) respectively. Fig. 5.13 shows the evolution in motion of the membrane
described in case 3 over one cycle of its motion . It is clear from the snapshots that both
internal and external boundary conditions are satisfied and this is indicative of the accuracy
of the simulations.

120

x
y

z

z
y

(a)

x
y

z

z
y

(b)

Figure 5.11: Snapshots of the membrane described in case 1 during one second of its vibration
121

z

z

z

x

x

x

y

y

y

(a)

(b)

(c)

Figure 5.12: Snapshots of the membrane described in case 2 during one second of its vibration
122

z
y
x

z

z
y

x

Figure 5.13: Snapshots of the membrane described in section case 3 during one cycle of its
vibration
123

Chapter 6
Vibration Control of a Circular
Membrane Using An Areal Constraint

6.1

Introduction

In this chapter we use the accurate dynamics analysis and simulation procedure presented in
Chapter 5 to investigate the effect of constraint application and removal on the energetics of
a vibrating membrane. The constraint is assumed to be circular and applied suddenly at an
arbitrary time during free vibration of the membrane, which is vibrating in its fundamental
mode. The energetics of the membrane due to constraint application is investigated for
different sizes and locations of the constraint. This chapter is organized as follows: A formal
problem statement and a list of assumptions are provided in Section 6.2. Assuming that the
membrane is vibrating in its fundamental mode, we use analytical methods in Section 6.3
to study the dynamics of the membrane and change in its energy after application of the
constraint. The dynamics of the membrane and change in its energy after removal of the
124

constraint is studied in Section 6.4. In Section 6.5 we present numerical simulation results for
one cycle of constraint application and removal - these results indicate that the energy of the
membrane can increase or decrease depending on the time of application of the constraint.
In Section 6.6 we repeat the analyses of Sections 6.3 and 6.4 for arbitrary initial conditions
of the membrane and develop different control strategies to suppress the vibration of the
membrane. Numerical simulation results are then presented to demonstrate the efficacy of
the control strategies.

6.2

Problem Statement and Assumptions

Consider a circular membrane of radius R initially vibrating freely in its fundamental mode
as shown in Fig. 6.2 (a). Assume a circular areal constraint of radius a to be applied instantaneously near the boundary of the membrane as shown Fig. 6.2 (b) at some arbitrary time
t = tc during the motion of the membrane. The constraint is assumed to be applied during
membrane motion at a distance d from the center of the membrane as shown in Fig. 6.2 (c).
We investigate the effect of applying the constraint on the total energy under the following
simplifying assumptions:
A1. The membrane is homogeneous and has a constant mass per unit area denoted by µ.
The tension in the membrane is equal to T and remains constant at all times.
A2. The membrane is initially vibrating in its fundamental mode. This assumption will
be removed and a general displacement profile of the membrane will be assumed in
Section 6.6.
A3. The amplitude of oscillation of the membrane is small and therefore the equation of
125

motion of the membrane can be expressed by the standard relation in polar coordinates
[75]
∂ 2 Z0 1 ∂Z0
1 ∂ 2 Z0
1 ∂ 2 Z0
+
=
+
ρ ∂ρ
∂ρ2
ρ2 ∂ψ 2
c2 ∂t2

(6.1)

where Z0 = Z0 (ρ, ψ, t) is the transversal displacement of the freely vibrating membrane
and c2 = T /µ. The distance ρ is measured from the center of the membrane, and ψ is
the angle that ρ makes with the x-axis as shown in Fig. 6.1.
y
P
ρ
R
ψ

x

Figure 6.1: Circular membrane with radius R .

A4. The membrane has no internal damping, i.e., the energy of the membrane will remain
conserved during free vibration.
A5. At time t = tc , the areal constrained is applied instantaneously as shown in Fig. 6.2
(b), at a distance d from the center of the membrane. The application of the constraint
imposes an areal zero-displacement over a small region of the vibrating membrane. The
displacement and velocity of the membrane over the remaining region of the membrane
remains unchanged immediately after application of the constraint.
126

A6. At time t = tr , tr > tc , the constraint is instantaneously moved back to its original
position, i.e., the constraint is removed to the outside of the membrane. The displacement and velocity of the membrane over the unconstrained region of the membrane
remains unchanged immediately after removal of the constraint.

z

y
x

(a)
t = tc
moving constraint

(b)
y
P
r

ρ

R

ψ

θ
a O

Γi
x

d
Γo
(c)
Figure 6.2: (a) A vibrating membrane (b) A zero-displacement circular constraint is applied
to the membrane at time t = tc over a small region of the membrane (c) Geometry of the
membrane and the applied constraint after constraint application
127

6.3

Effect of Applying Constraint on the Dynamics of
the Membrane

6.3.1

Equation of motion after application of the constraint

The geometry of the membrane and the constraint is shown in Fig. 6.2 (c). To investigate
the dynamics of the constrained membrane, we transform the coordinates of membrane from
(ρ, ψ) to (r, θ) where both of the coordinates systems are shown in Fig. 6.2 (c), and the
governing differential equation can be written as
∂ 2 Zc 1 ∂Zc
1 ∂ 2 Zc
1 ∂ 2 Zc
+
+
=
r ∂r
∂r 2
r 2 ∂θ2
c2 ∂t2

(6.2)

where Z = Z(r, θ, t) is the transversal displacement. The distance r is measured from the
center of the inner circle, O, and θ is the angle that r forms with the x-axis as shown in
Fig. 6.2. The accurate dynamics analysis of the constrained membrane in (r, θ) coordinates
is discussed in detail in Chapter 5 and the general solution is obtained in Eq. (5.37) as:
∞
Zc (r, θ, t) =

c
c
c
ϕc (ki r, θ)(A1i cos ωi t + B1i sin ωi t)
i

i=0
∞
s
s
s
+
ϕs (kj r, θ)(A2j cos ωi t + B2j sin ωj t)
j
j=1

(6.3)

c
s
where ϕc (ki r, θ) and ϕs (kj r, θ) are the symmetric and antisymmetric shape functions, which
i
j
c
s
are presented in Eqs. (5.19) and (5.23). Also, ki and kj are the symmetric and antisymmetric
eigenfrequencies of the constrained membrane. The model coefficients A1i , B1i , A2j and B2j
are to be determined by applying the initial conditions.
128

6.3.2

Initial Conditions

The membrane is initially vibrating in the its first mode as:

Z0 (ρ, ψ, t) = A0 J0 (k00 ρ) cos ω00 t

(6.4)

where A0 is the maximum amplitude of oscillation. The constraint is applied at an arbitrary
time tc and therefore the initial conditions are

A1. The displacement of the membrane can be written as follows


 A J (k ρ) cos ω tc : a < r ≤ R2 − d2 sin2 θ − d cos θ
0 0 00
00
Zc (r, θ, 0) =


0 : 0≤r≤a

(6.5)

Applying the initial conditions to Eqs. (6.3) and (6.4), we get
∞
i=0

c
A1i ϕc (ki r, θ) +
i

∞
j=1

s
A2j ϕs (kj r, θ) = A0 J0 (k00 ρ) cos ω00 tc
j

(6.6)

To obtain the coefficients A1j and A2j , we use the orthogonality conditions between
modes, where the orthogonality conditions between distinct modes are obtained in
Eq. (5.35). By multiplying Eq. (6.6) with ϕc and integrating both sides over the effeci
tive membrane surface and using the orthogonality condition presented in Eq. (5.35),
we get

A1i =

A0 cos ω00 tc 2π r(θ)
J0 (k00 ρ) ϕc r dr dθ
i
κc
a
0
129

(6.7)

where

κc =

2π
0

r(θ)
a

(ϕc )2 r dr dθ
i

(6.8)

Multiplying Eq. (6.6) with ϕs and integrating both sides over the effective membrane
j
surface and using the orthogonality conditions in Eq. (5.35), gives

A2j =

A0 cos ω00 tc 2π r(θ)
J0 (k00 ρ) ϕs r dr dθ
j
κs
0
a

(6.9)

where

κs =

2π
0

r(θ)
a

(ϕs )2 r dr dθ
j

(6.10)

A2. The velocity of the membrane can be written as follows

 ∂Z0 (ρ, ψ, tc )
: a < r ≤ R2 − d2 sin2 θ − d cos θ
∂Zc (r, θ, 0) 
∂t
=

∂t

0 : 0≤r≤a
∞
⇒
i=0

c
B1i ωi ϕc +
i

∞
j=1

s
B2j ωj ϕs = −ω00 A0 J0 (k00 ρ) sin ω00 tc
j

(6.11)

(6.12)

By multiplying Eq. (6.12) with ϕc and integrating both sides over the membrane suri
face and using the orthogonality conditions in Eq. (5.35), we get

B1i =

− ω00 A0 sin ω00 tc 2π r(θ)
J0 (k00 ρ) ϕc r dr dθ
c κc
i
ωi
0
a
130

(6.13)

By multiplying Eq. (6.12) with ϕs and integrating both sides over the effective memi
brane surface and using the orthogonality condition presented in Eq. (5.35), we get

B2j =

6.3.3

− ω00 A0 sin ω00 tc 2π r(θ)
J0 (k00 ρ) ϕs r dr dθ
s
j
ωj κs
0
a

(6.14)

Membrane Energy Before Application of Constraint

The total energy of the circular membrane in polar coordinates can be written as1
∂Z(r, θ, t) 2
r dr dθ
∂t

µ
E =
2
A
+

T
2
A

1
∂Z(r, θ, t) 2
+
∂r
r2

∂Z(r, θ, t) 2
∂θ

r dr dθ

(6.15)

Where µ is the mass per unit area of the membrane and T is the tension per unit length.
The initial energy of membrane can be given by rewriting Eq. (6.15) in (ρ, ψ) coordinates
and with integrating over the membrane area as follows:
µ 2π R ∂Z0 (ρ, ψ, t) 2
E0 =
ρ dρ dψ
2 0
∂t
0
∂Z0 (ρ, ψ, t) 2
1
T 2π R
+
+
2 0
∂ρ
ρ2
0

1 See

appendix C

131

∂Z0 (ρ, ψ, t) 2
∂ψ

ρ dρ dψ

(6.16)

By differentiating Eq. (6.4) with respect to t, ρ and ψ, we get

∂Z0 (ρ, ψ, t)
= −A0 ω0 J0 (ko ρ) sin ω0 t
∂t
∂Z0 (ρ, ψ, t)
= −A0 ko J1 (ko ρ) cos ω0 t
∂ρ
∂Z0 (ρ, ψ, t)
= 0
∂ψ

(6.17)

Inserting Eqs. (6.17) in Eq. (6.16) leads to
µ 2π R
E0 =
(−A0 ω0 J0 (ko ρ) sin ω0 t)2 ρ dρ dψ
2 0
0
2π R
T
(−A0 ko J1 (ko ρ) cos ω0 t)2 ρ dρ dψ
+
2 0
0
Note that
R
0

2
Jn (kρ) ρ dρ =

R2
2
Jn (kR) − Jn−1 (kR)Jn+1 (kR)
2

(6.18)

and J−1 (ko R) = −J1 (ko R), ω0 = ck0 , T = c2 µ, and J0 (ko R) = 0. Then the total energy
before the constraint application is

E0 =

6.3.4

Tπ
(A0 k0 J1 (ko R))2
2

(6.19)

Membrane Energy After Application of Constraint

In this section we investigate the energy of the membrane after the application of the constraint. The equation of membrane motion after the constraint application is obtained by
Zc (r, θ, t) in Eq. (6.3), hence, the total energy of the membrane after constraint application
132

can be given by using Eq. (6.15), and integrating over the constrained membrane area as
follows
1
T 2π r(θ)
∂Zc (r, θ, t) 2
+
2 0
∂r
r2
a
µ 2π r(θ) ∂Zc (r, θ, t) 2
+
r dr dθ
2 0
∂t
a

Ec =

∂Zc (r, θ, t) 2
∂θ

r dr dθ
(6.20)

By differentiating Eq. (6.4) with respect to t, r and θ, we get

∂Zc (r, θ, t)
=
∂t

∞
i=0
∞

+
j=1

c
c
c
c
ωi ϕc (ki r, θ)(−A1i sin ωi t + B1i cos ωi t)
i
s
s
s
s
ωj ϕs (kj r, θ)(−A2j sin ωi t + B2j cos ωj t)
j

∞ ∂ϕc (k c r, θ)
c
c
i i
(A1i cos ωi t + B1i sin ωi t)
∂r
i=0
∞ ϕs (k s r, θ)
j j
s
s
+
(A2j cos ωi t + B2j sin ωj t)
∂r
j=1

∂Zc (r, θ, t)
=
∂r

(6.21)

where,
c
∂ϕc (ki r, θ)
i
∂r
s
ϕs (kj r, θ)
j
∂r

N+1
=
m=0
N
=
m=1

i
β2m
j
β1m

c
∂Jm (ki r)
∂r
s
∂Jm (kj r)
∂r

133

c
i ∂Ym (ki r) cos mθ
+ γ2m
∂r
j
+ γ1m

s
∂Ym (kj r)
∂r

sin mθ

(6.22)

where, in general,

∂Jm (kr)
k
J
(kr) − Jm+1 (kr)
=
∂r
2 m−1
k
∂Ym (kr)
Y
(kr) − Ym+1 (kr)
=
∂r
2 m−1

(6.23)

Also,
∞ ∂ϕc (k c r, θ)
c
c
i i
(A1i cos ωi t + B1i sin ωi t)
∂θ
i=0
∞ ∂ϕs (k s r, θ)
j j
s
s
+
(A2j cos ωi t + B2j sin ωj t)
∂θ
j=1

∂Zc (r, θ, t)
=
∂θ

(6.24)

where
c
∂ϕc (ki r, θ)
i
∂θ

N
=

i
c
i
c
−m β2m Jm (ki r) + γ2m Ym (ki r) sin mθ

m=0
s (k s r, θ)
N+1
∂ϕj j
j
j
s
s
=
m β1m Jm (kj r) + γ1m Ym (kj r) cos mθ
∂θ
m=1

(6.25)

The total energy of the membrane after the application of the constraint Ec can be calculated
by inserting the equations (6.21-6.25) into Eq. (6.20) and integrating numerically. The work
done by the constraint Wd and the percentage change in energy are expressed as:

Wd = Ec − E0
%Change in Energy =

134

(6.26)

Ec − E0
× 100
E0

(6.27)

6.4

Effect of Removing Constraint on the Dynamics of
the Membrane

6.4.1

Equation of motion after removal of the constraint

We analyse the dynamics of the membrane after constraint removal. Assume the constraint is
removed at an arbitrary time t = tr after application of the constraint, such that tr > tc > 0.
When the constraint is removed, the membrane vibrates freely according to the equation of
motion of an unconstrained membrane. The equation of motion of simple unconstrained
¯
membrane can be obtained by solving the PDE in Eq. (6.1). Assume Zf (ρ, ψ, t) is the
transverse displacement of the membrane after the removal of the constraint, and the governing differential equation becomes
2¯
2¯
¯
¯
∂ 2 Zf
1 ∂ Zf
1 ∂ Zf
1 ∂ Zf
+
=
+
ρ ∂ρ
∂ρ2
ρ2 ∂ψ 2
c2 ∂t2

(6.28)

The general solution of Eq. (6.28) is given by using the separation of variables and it can be
written as2
∞ ∞
¯
Zf (ρ, ψ, t) =

Jq (kqi ρ)(Cqi cos ωqi t + Dqi sin ωqi t) cos qψ
q=0 i=0
∞ ∞

+

Jq (kqi ρ)(Eqi cos ωqi t + Fqi sin ωqi t) sin qψ

(6.29)

q=1 i=1
where kqi represents the eigenfrequencies of the unconstrained membrane and ωqi = c kqi
and the model coefficients Cqi , Dqi , Eqi and Fqi are constants to be determined by applying
2 See

Appendix A.2

135

the initial conditions. When the constraint is removed at t = tr , previous motion of the
constrained membrane represents the initial conditions for the upcoming membrane motion
after constraint removal. Thus, the initial conditions can be written as:
Zc (r, θ, t), t < tr

z
y
x

(a)
¯
Zf (ρ, ψ, t), t ≥ tr
t ≃ tr
constraint is removed

(b)
Figure 6.3: (a) A vibrating constrained membrane (b) The circular constraint is removed
from the membrane at time t ≃ tr

A1. The displacement of the membrane can be obtained by applying the initial conditions
to Eqs. (6.3) and (6.29) as follows


 Zc (r, θ, tr ) : a < r ≤ R2 − d2 sin2 θ − d cos θ
¯ (ρ, ψ, 0) =
Zf


0 : 0≤r≤a
136

(6.30)

then
∞ ∞
Zc (r, θ, tr ) = C00 J0 (k00 ρ) +

Jq (kqi ρ) Cqi cos qψ + Eqi sin qψ
q=1 i=1
(6.31)

Note that the orthogonality conditions for unconstrained membrane are:
R

Jq (kqi ρ)Jq (kqj ρ)ρdρ = 0, i = j
0
2π
2π
sin mψ sin qψ dψ =
cos mψ cos qψ dψ = 0,
0
0
2π
cos mψ sin qψ dψ = 0, ∀q, m
0

q=m
(6.32)

By multiplying Eqs. (6.31) by J0 (k00 ρ) and integrating both sides over the membrane
surface and using the orthogonality condition presented in Eq(6.32), we get

C00 =
R

where △00 =

△qi =

0
R
0

2π R
1
Zc (r, θ, tr )J0 (k00 ρ) ρ dρ dψ
2π △00 0
0

(6.33)

2
J0 (k00 ρ) ρ dρ, and note that from Eq. (6.18) we can write

2
Jn (kqi ρ) ρ dρ =

R2
2
Jn (kqi R) − Jn−1 (kqi R)Jn+1 (kqi R)
2

(6.34)

By repeating the procedure above after multiplying Eq. (6.31) by Jq (kqj ρ) cos mψ and

137

Jq (kqj ρ) sin mψ respectively, and integrating over the membrane area, we get
2π R
1
Zc (r, θ, tr )Jq (kqi ρ) cos qψ ρ dρ dψ
π △qi 0
0
2π R
1
Zc (r, θ, tr )Jq (kqi ρ) sin qψ ρ dρ dψ
Eqi =
π △qi 0
0

Cqi =

(6.35)

A2. The velocity of the membrane can be written as

 ∂Zc (r, θ, tr )
¯
: a < r ≤ R2 − d2 sin2 θ − d cos θ
∂ Zf (ρ, ψ, 0) 
∂t
=

∂t

0 : 0≤r≤a

(6.36)

this gives
∞ ∞
D00 ω00 J0 (k00 ρ) +

ωqi Jq (kqi ρ) Dqi cos qψ + Fqi sin qψ =
q=1 i=1

∂Zc (r, θ, tr )
∂t
(6.37)

Following the procedure in which the coefficients in Eqs. (6.33) and (6.35) were obtained, we get
2π

R ∂Zc (r, θ, tr )
J0 (k00 ρ) ρ dρ dψ
2πω00 △00 0
∂t
0
2π R ∂Zc (r, θ, tr )
1
Jq (kqi ρ) cos qψ ρ dρ dψ
Dqi =
πωqi △qi 0
∂t
0
2π R ∂Zc (r, θ, tr )
1
Fqi =
Jq (kqi ρ) sin qψ ρ dρ dψ
πωqi △qi 0
∂t
0
D00 =

1

(6.38)

One can note that the variables in the integrations in Eqs. (6.33), (6.35),and (6.38) are
different (ρ, ψ and r, θ ). In addition, the integrations should be performed by excluding the
138

area covered by the constraint because it has zero displacement and zero velocity. Excluding
the constrained area by using (ρ, ψ ) coordinates is quite challenging. Therefore, we change
the order of the integrations and the variables of the integration will be (r, θ) coordinates
system only. We can unify the integrations variables by using the relationships between (ρ, ψ
and r, θ ), where from Fig. 6.2 we can write

ρ=

d + r cos θ
ρ

ψ = cos−1

d2 + r 2 + 2 r d cos θ,

(6.39)

The area of the membrane is constant and can be obtained by

Area =

2π
0

r(θ)
0

2π

r dr dθ =

0

R
0

ρ dρ dψ

Using Eq. (6.39) and integrating over the area with (r dr θ) instead of (ρ dρ dψ) we can get
2π
0

R
0

h(r, θ, tr ) J0 (k00 ρ) ρ dρ dψ +

=

where h(r, θ, tr ) could be Zc (r, θ, tr ) or

0 J0 (k00 ρ) ρ dρ dψ
Ac
2π
0

r(θ)
a

∂Zc (r, θ, tr )
.
∂t

139

h(r, θ, tr ) J0 (k00 ρ(θ)) r dr dθ

6.4.2

Membrane Energy After Removal of Constraint

The total energy of the membrane after constraint removal can be obtained by using Eq. (6.15),
and integrating over the unconstrained membrane area

Er =


¯
R  ∂ Zf (ρ, ψ, t) 2
1
+
∂ρ
ρ2
0 

2π

T
2 0

µ 2π R
+
2 0
0

¯
∂ Zf (ρ, ψ, t) 2
∂t

ρ dρ dψ


¯ (ρ, ψ, t) 2 
∂ Zf
ρ dρ dψ

∂ψ

(6.40)

¯
The derivatives of Zf (ρ, ψ, t) with respect to ρ ,ψ and t are obtained from Eq. (6.29) as
follows
∞ ∞ ∂J (k ρ)
q qi
c
s
fqi (t) cos qψ + fqi (t) sin qψ
=
∂ρ
∂ρ
q=1 i=1
∞ ∞
¯
∂ Zf (ρ, ψ, t)
c
s
=
qJq (kqi ρ) −fqi (t) sin qψ + fqi (t) cos qψ
∂ψ
q=1 i=1
s
c
∞ ∞
¯
dfqi (t)
dfqi (t)
∂ Zf (ρ, ψ, t)
Jq (kqi ρ)
=
cos qψ +
sin qψ
∂t
dt
dt
q=1 i=1
¯
∂ Zf (ρ, ψ, t)

(6.41)

where
∂Jm (kqi ρ)
∂ρ

=

kqi
2

Jm−1 (kqi ρ) − Jm+1 (kqi ρ)

c
fqi (t) = Cqi cos ωqi t + Dqi sin ωqi t
s
fqi (t) = Eqi cos ωqi t + Fqi sin ωqi t
c
dfqi (t)
= ωqi −Cqi sin ωqi t + Dqi cos ωqi t
dt
s
dfqi (t)
= ωqi −Eqi sin ωqi t + Fqi cos ωqi t
dt
140

(6.42)

By inserting Equations (6.43) and (6.42) in Eq. (6.40), and using the numerical integration
described in Eq. (5.44), the total energy of the membrane after constraint removal can be
obtained.

6.5

Simulation of One Cycle of Constraint Application
and Removal

In this section we present numerical simulations for the circular vibrating membrane subjected to a sudden areal constraint for different time of constraint applications tc and constraint sizes and locations. Consider a circular membrane oscillating in the first mode as
follows:
Z0 (ρ, ψ, t) = A0 J0 (k00 ρ) cos(ω00 )t

(6.43)

where the specifications of the membrane are :

R = 1 m,

µ = 0.25 kg/m2 ,

T = 1 N/m

and k00 = 2.4048 m −1 , A0 = 0.1 m and ω00 = 4.8096 rad/sec
From ω00 we can obtain the time period which is 1.3064 sec and Fig. 6.3 (a) shows the
shape of the membrane at different instants of time over the interval t ∈ [0, 0.6532] sec. The
membrane is in its maximum potential energy configuration at t = 0 sec and t = 0.6532 sec,
and maximum kinetic energy configuration at t = 0.3266 sec when the membrane passes
through the mean position.
We first present simulation results for percentage change in energy due to constraint of
141

size a/R = 0.1 and locations d/R ∈ {0.0, 0.6, 0.9} and tc ∈ [0.24, 0.42] sec as shown in
Fig. 6.4 (b). The simulation is performed using 16 symmetric and 16 antisymmetric modes.
The results, shown in Fig. 6.4 (b), indicate that the energy of the membrane increases if the
constraint is applied when the membrane is far away from the mean position, and decreases
if the constraint is applied when the membrane is near its mean position, irrespective of
the value of d. Same result can be observed for different constraint size (a/R = 0.12)
and locations d/R ∈ {0.0, 0.44, 0.88} as shown in Fig 6.4 (c). The results also show that
larger constraint causes more change (increase/decrease) in the total energy, as expected.
Similar observation were made in our previous work on vibration suppression of string using
a scabbard-like actuator, presented in Chapter 3.

The results can be explained as follows. When the constraint is applied and an area
of the membrane is brought to rest, the kinetic energy of the seized membrane is removed.
However, the potential energy in the membrane is increased due to sudden change in the
slope adjacent to the constraint. The potential energy increases more as the constraint is
applied when the membrane is far away form the mean position, and hence, the total energy
is increased. When the constraint is applied as the membrane passes through the mean
position, the change in potential energy is almost zero and the kinetic of the area seized by
the constraint is lost, hence the total energy of the membrane is reduced. In addition, the
kinetic energy is proportional to the area of the membrane which is covered by the constraint,
and therefore, larger constraint size causes higher energy reduction when the constraint is
applied as the membrane passes through the mean position. That is because the velocity of
the constrained area is higher when the it is located close to the center of the membrane. For
the same value of tc , a smaller value of d results in a higher decrease in total energy when
142

the constraint is applied as the membrane passes through the mean position. The general
statement of energy reduction condition can be stated as: The maximum energy reduction
occurs when the constraint is applied as the membrane passes through its mean position .
The energy loss due to the loss of kinetic energy of the membrane area which is covered by
the constraint . We can use this condition to apply and remove the constraint sequentially,
or move the constraint over the membrane area to reduce the total energy of the membrane.
Next we investigate the change in energy of the membrane due to removal of the constraint. We compute the values of Ec and Er , which equal the energy of the membrane
before and after removal of the constraint respectively. Using Eqs. (6.20), Ec is computed
for 16 symmetric and 16 antisymmetric modes. We fix the value of the energy before the constraint removal Ec and we compute Er from Eq. (6.40) with varying the number of Fourier
coefficients. From this data (see Table 6.1) it is clear that Ec and Er approach each other
as the number of Fourier coefficients increase. Although the data corresponds to the specific
case of d/R = 0.88, a/R = 0.12 and tr = 0.4 sec, the same trends can be observed for all
values of d/R, a/R and tr . In the next section we investigate the effect of application of
the constraint on the dynamics of the membrane with arbitrary initial conditions such that
sequential application and removal of constraints can be explored as a strategy for vibration
suppression.
Table 6.1: Comparison of values of Er with a fixed value of Ec = 2.4455 J for constrained
membrane with d/R = 0.88, a/R = 0.12 and tr = 0.4 sec
No. of Fourier coefficients
10
30
64
80

143

Er (J)
2.4176
2.4335
2.4365
2.4391

error %
-1.1428
-0.4903
-0.3683
-0.2604

Z
0.1

t=0

0.0

t = 0.6532

t = 0.3266

−0.1
0.0

-1

x

1

(a)

E%

E%

0
d
R = 0.9

-8
0.24

0.28

0

-8
0.32

d
R =0
d
R = 0.44

8

d =0
R
d
R = 0.6

8

0.36

0.4

d
R = 0.88
0.24

0.28

0.32

Time (sec)

0.4

Time (sec)

(b)

0.36

(c)

Figure 6.4: Percentage change in total energy of the membrane due to application of the constraint of different size of values
of of d and to. (a) position of the membrane at different time instants. (b) a/R = 0.1 (c) a/R = 0.12
144

6.6
6.6.1

Sequential Application and Removal of Constraints
Dynamics of membrane with arbitrary initial conditions

The analysis in this chapter thus far is based on a membrane that initially vibrates in its
the fundamental mode. We extend the analysis to include arbitrary initial vibration of the
membrane. The equation of motion of unconstrained membrane that includes all possible
initial conditions can be written as
∞ ∞
Zo (ρ, ψ, t) =

Jq (kqi ρ)(αqi cos ωqi t + βqi sin ωqi t) cos qψ
q=0 i=0
∞ ∞

+

Jq (kqi ρ)(Γqi cos ωqi t + δqi sin ωqi t) sin qψ

(6.44)

q=1 i=1
which is similar to Eq(6.29) and satisfies the boundary conditions of unconstrained membrane, and the coefficients αqi , βqi , Γqi and δqi are known. By using Eq. (6.15), the initial
energy of the membrane can be written as
2π R
1
∂Zo (ρ, ψ, t) 2
¯o = T
+
E
2 0
∂ρ
ρ2
0
µ 2π R ∂Zo (ρ, ψ, t) 2
+
ρ dρ dψ
2 0
∂t
0

∂Zo (ρ, ψ, t) 2
∂ψ

ρ dρ dψ
(6.45)

where it is identical to Er in Eq. (6.40), and the same calculation procedure for Er can be
used to obtain Eo . The equation of motion of the membrane after constraint application
can be written as

¯
Zc (r, θ, t) =

∞
i=0

c
ϕc (ki r, θ) fi (t) +
¯i c
145

∞
j=1

s
ϕs (kj r, θ) fj (t)
¯j s

(6.46)

c
c
c
s
s
s
¯
¯
¯
¯
where fi (t) = A1i cos ωi t + B1i sin ωi t and fj (t) = A2j cos ωi t + B2j sin ωj t. Equation
(6.46) is identical to Eq. (6.3) except that ϕc , ϕs , A1i , A2i , B1j , and B2j are replaced by
i j
¯
¯
¯
¯
¯
¯
¯
ϕc , ϕs , A1i , A2i , B1j , and B2j . The coefficients A1i , A2i , B1j , and B2j are determined
¯i ¯j ¯
from the initial conditions by repeating the procedure used in Section 6.3.1. Without loss of
generality, the time is reset from t = tc to t = 0, and the displacement and velocity of the
membrane given by

¯
Zc (r, θ, 0) = Zo (ρ, ψ, tc )
¯
∂Zo (ρ, ψ, tc )
∂ Zc (r, θ, 0)
=
∂t
∂t

(6.47)

By using Eqs. (6.44), (6.46) and (6.47) and following the same procedure used in Section
6.3.1, we can write
2π r(θ)
1
¯
Zo (ρ, ψ, tc ) ϕc r dr dθ
¯i
A1i =
κc 0
a
2π r(θ)
1
¯
Zo (ρ, ψ, tc ) ϕs r dr dθ
¯j
A2j =
κs 0
a
2π r(θ) ∂Zo (ρ, ψ, tc )
1
¯
B1i =
ϕc r dr dθ
¯i
c κc
ωi
∂t
0
a
2π r(θ) ∂Zo (ρ, ψ, tc )
1
¯
B2j =
ϕs r dr dθ
¯j
s
ωj κs 0
∂t
a

(6.48)

The total energy of the membrane after constraint application is
¯
1
T 2π r(θ)
∂ Zc (r, θ, t) 2
¯
+
Ec =
2 0
∂r
r2
a
¯
µ 2π r(θ) ∂ Zc (r, θ, t) 2
+
r dr dθ
2 0
∂t
a

146

¯
∂ Zc (r, θ, t) 2
∂θ

r dr dθ
(6.49)

6.6.2

Control strategies for vibration suppression

Similar to the vibrating string subjected to a constraint applied by a scabbard-like (see
Chapter 3), the simulation results in Section 6.5 show that the maximum reduction in energy of the membrane is achieved if the constraint is applied when the displacement of the
membrane is zero over the area r ≤ a. For arbitrary initial conditions, the displacement of
the membrane can be uniformly zero over a small area adjacent to the constraint at some
given time, and therefore, energy reduction can be ensured by choosing a small distance
that the constraint moves. We investigate three different control strategies for the vibration
suppression of a circular membrane. These different control strategies are presented in next
three sections.

6.6.3

Simple Application and Removal of Constraint

In this section we investigate energy dissipation due to simple application and removal of the
areal constraint. The constraint is applied at a distance d from the center of the membrane
during its motion, and then removed to its original location outside the membrane as shown
in Fig 6.5 (b). We use the energy reduction condition from Section 6.5 to compute the
total energy of the vibrating membrane through sequential application and removal of the
constraint. Assume that the membrane is initially vibrating freely with arbitrary initial
conditions as in Eq. (6.44). A set of sensors are mounted on the mechanism applying the
¯
constraint to sense the displacements of the membrane of the circular area Ω to which the
¯
constraint will be applied, see Fig. 6.5 (b). When the displacement of the sensed area Ω is
approximately zero, the constraint is applied, and after some arbitrary time tr the constraint
is removed and the procedure will be repeated.

147

Z
X
vibrating membrane

(a)

sensors
application
¯
Ω

a

Y

d
R

removal

X

constraint

(b)
Figure 6.5: Sensing mechanism for simple application and removal of constraint (a) side view. (b) top view
148

The procedure of simple application and removal of the constraint can be summarized as
follows:

A1. The membrane initially vibrates freely and its equation of motion is given by Zp (ρ, ψ, t) =
Zo (ρ, ψ, t), where Zo (ρ, ψ, t) is presented in Eq. (6.44) and p denotes the number of
times the constraint has applied. The initial energy E0 of the membrane is computed
by using Eq. (6.16).

(p)
¯
A2. When the displacement of sensed area Ω is ≈ 0, then t = tc and the constraint is
applied.

A3. Reset the time t = 0 and find the equation of the motion of the membrane after con¯
straint application Zc (r, θ, t) which is presented in Eq. (6.46). The procedure presented
in Section 6.6.1 is used to obtain the model coefficients as follows :
2π r(θ)
1
(p)
¯(p)
A1i =
Zp (ρ, ψ, tc ) ϕc r dr dθ
¯i
κc 0
a
2π r(θ)
1
(p)
¯(p)
A2j =
Zp (ρ, ψ, tc ) ϕs r dr dθ
¯j
κs 0
a
2π r(θ) ∂Zp (ρ, ψ, t(p) )
1
c
¯ (p) =
ϕc r dr dθ
¯i
B1i
c κc
ωi
∂t
0
a
2π r(θ) ∂Zp (ρ, ψ, t(p) )
1
(p)
c
¯
B2j =
ϕs r dr dθ
¯j
sκ
ωj s 0
∂t
a
where κc and κs are presented in Eqs. (6.8) and (6.10) respectively.

¯ (p)
A4. Compute the total membrane energy after constraint application Ec by using Eq. (6.49)
149

and then compute the energy reduction by
¯ (p)
Eo − Ec
E% =
× 100
Eo

(6.50)

A5. Let the constrained membrane vibrate freely (starting from t = 0 ) and let the con(p)
straint be removed after some arbitrary time t = tr .
A6. Reset the time t = 0, and use the previous displacement and velocity of the con(p)
¯
strained membrane Zc (r, θ, tr ) as initial conditions to find the equation of motion
¯ (p+1) (ρ, ψ, t), which is presented in Eq. (6.29). The proceafter constraint removal Z
f
dure presented in Section 6.4.1 is used to calculate the modal coefficients, as follows
2π R (p)
1
(p)
¯
Zc (r, θ, tr ) J0 (k00 ρ) ρ dρ dψ
2π △00 0
0
2π R (p)
1
(p+1)
(p)
¯
=
Cqi
Zc (r, θ, tr ) Jq (kqi ρ) cos qψ ρ dρ dψ
π △qi 0
0
(p+1)
C00
=

2π R (p)
1
(p)
(p+1)
¯
Zc (r, θ, tr ) Jq (kqi ρ) sin qψ ρ dρ dψ
=
Eqi
π △qi 0
0
2π R ∂ Z (p) (r, θ, t(p) )
¯
1
(p+1)
c
r
D00
=
J0 (k00 ρ) ρ dρ dψ
2πω00 △00 0
∂t
0
2π R ∂ Z (p) (r, θ, t(p) )
¯c
1
(p+1)
r J (k ρ) cos qψ ρ dρ dψ
=
Dqi
q qi
πωqi △qi 0
∂t
0
2π R ∂ Z (p) (r, θ, t(p) )
¯
1
(p+1)
c
r J (k ρ) sin qψ ρ dρ dψ
=
Fqi
q qi
πωqi △qi 0
∂t
0

A7. Update the equation of the unconstrained membrane motion by setting Zp (ρ, ψ, t) =
¯ (p+1) (ρ, ψ, t) and increase the counter , where p = p + 1 and go to step 2.
Z
f

150

t = t + ∆t
Start

No

¯
Ω
Zo ≃ 0 ?

simulate Zo (ρ, ψ, t)

t=0
Zo (ρ, ψ, t) = Zf (ρ, ψ, t)

Yes
tc = t
remove constraint
set initial conditions
˙
Zo (ρ, ψ, tc ) and Zo (ρ, ψ, tc )

calculate model coefficient
of unconstrained membrane
=⇒ Zf (ρ, ψ, t)
set initial conditions
˙
Zc (r, θ, tr ) and Zc (r, θ, tr )

apply constraint
t = t + ∆t
No
calculate model coefficient
of constrained membrane
=⇒ Zc (r, θ, t)

tr = t

calculate E%
Yes
t=0

simulate Zc (r, θ, t)

t = tr ?

Figure 6.6: Flowchart for simple application and removal of constraint strategy
151

Note that simple constraint application and removal requires computing the eigenfrequencies
and mode shapes for only one time and storing it and using it offline during sequential
application and removal of constraint. A flowchart for the procedure described in this section
is shown in Fig 6.6.

6.6.4

Radial Constraint Application

In this section we study effect of moving the constraint radially on the energy of the membrane. The constraint is first applied near the membrane, boundary at a distance d1 from the
center of the membrane and then moved towards the center of the membrane at a distance
d2 , where d2 < d1 and so on till the constraint reaches a prescribed distance dp such that
dp < · · · < d2 < d1 as shown in Fig. 6.7 (b). Thereafter, the constraint is moved back to a
distance dp−1 and then dp−2 and so on till the constraint comes back to its original location
at a distance d1 . The sensing mechanism for this case is shown in Fig. 6.7, where the sensors
¯
record the displacement of the area Ω which lies on the path of the constraint movement (on
X-axis). The case of radial constraint motion is more complicated than the simple constraint
application and removal case. The change of the constraint location implies computation
of eigenfrequencies and mode shapes for each application distance d during the constraint
movement. When the constrained is moved from its location at d1 to its location at d2 , the
coordinate system changes, see Fig. 6.8. As a result, computation of modal coefficients of the
new location requires integration over the entire membrane area excluding both constraint
¯
¯
¯
areas Ω1 and Ω2 . This is because Ω1 represents an area of zero initial displacement and
velocity for the new constraint location at d2 and the upcoming motion is restricted to the
¯
area outside Ω2 .

152

Z
X

vibrating membrane

(a)

Y
X
application

Axis of constraint
motion

¯
Ω

sensors

a

constraint

d3
H

¯
Ω

d2
d1

(b)
Figure 6.7: Sensing mechanism and path of constraint motion for radial constraint motion (a) side view. (b) top view.
153

Y

P
R

ρ
ψ

r1

r2
θ2
¯
Ω2

θ1
¯
Ω1

X

d2
d1

Figure 6.8: The coordinates of two different locations of a constraint moved from location at d1 to location at d2 .

154

¯
¯
It is quite difficult to exclude both Ω1 and Ω2 mathematically, and therefore we use
an alternative method by storing the initial conditions of the previous motion (d1 location
with r1 , θ1 coordinates system) in the (ρ, ψ) coordinates system then transfer the initial
conditions into new coordinates r2 , θ2 which are associated with the new constraint location
d2 . The procedure of the constraint movement can be summarized as follows:
A1. The membrane is initially vibrates freely and the equation of motion Zo (ρ, ψ, t) is
presented in Eq. (6.44). The initial energy E0 of the membrane is computed by using
Eq. (6.16).
(1)
¯
A2. When the displacement of sensed area Ω1 is ≈ 0, then t = tc and the constraint is
applied.
A3. Reset the time t = 0 and find the equation of the motion of the membrane after
¯
constraint application Zc (r1 , θ1 , t) which is presented in Eq. (6.46) with replacing (r, θ)
by (r1 , θ1 ). The procedure presented in Section 6.6.1 is used to obtain the model
coefficients at constraint location d1 as follows
2π r(θ1 )
1
(1)
Zo (ρ, ψ, tc ) (ϕc )(1) r1 dr1 dθ1
¯i
(1) 0
a
κc
2π r(θ1 )
1
(1)
¯(1)
Zo (ρ, ψ, tc ) (ϕs )(1) r1 dr1 dθ1
¯j
A2j =
(1) 0
a
κs
2π r(θ1 ) ∂Zo (ρ, ψ, t(1) )
1
c
¯ (1) =
B1i
(ϕc )(1) r1 dr1 dθ1
¯i
(1) 0
∂t
c
a
(ωi )(1) κc
2π r(θ1 ) ∂Zo (ρ, ψ, t(1) )
1
c
¯ (1) =
(ϕs )(1) r1 dr1 dθ1
¯j
B2j
∂t
s )(1) κ(1) 0
a
(ωj
s
¯(1)
A1i =

Where the subscript (1) indicates that the calculations are associated with constraint
155

location d1 . Set the counter p = 1.
¯ (p)
A4. Compute the total membrane energy after constraint application Ec by using Eq. (6.49)
and then compute the energy reduction by using Eq. (6.50).
A5. Let the constrained membrane vibrate freely and when the displacement of sensed
(p)
¯
area Ωp+1 is ≈ 0, then t = tr and the constraint is removed. The initial conditions
are stored in the (ρ, ψ) coordinates system by setting the time t = 0, and finding
¯ (p)
the equation of motion after constraint removal Z (ρ, ψ, t) which is presented in
f
Eq. (6.29). By using similar procedure to that presented in Section 6.4.1, the modal
coefficients of the unconstrained motion can be obtained as
2π R (p)
1
(p)
¯
Zc (rp , θp , tr ) J0 (k00 ρ) ρ dρ dψ
2π △00 0
0
2π R (p)
1
(p)
(p)
¯
Cqi =
Zc (rp , θp , tr ) Jq (kqi ρ) cos qψ ρ dρ dψ
π △qi 0
0
(p)
C00 =

2π R (p)
1
(p)
(p)
¯
Zc (rp , θp , tr ) Jq (kqi ρ) sin qψ ρ dρ dψ
Eqi =
π △qi 0
0
¯ (p)
2π R ∂ Zc (rp , θp , t(p) )
1
(p)
r
J0 (k00 ρ) ρ dρ dψ
D00 =
2πω00 △00 0
∂t
0
¯ (p)
2π R ∂ Zc (rp , θp , t(p) )
1
(p)
r J (k ρ) cos qψ ρ dρ dψ
Dqi =
q qi
πωqi △qi 0
∂t
0
¯ (p)
2π R ∂ Zc (rp , θp , t(p) )
1
(p)
r J (k ρ) sin qψ ρ dρ dψ
Fqi =
q qi
πωqi △qi 0
∂t
0
where △qi is presented in Eq. (6.34).
A6. Transfer the initial conditions stored in (ρ, ψ) coordinates system (from step 5) to
(rp+1 , θp+1 ) which represents the coordinates system of the new constraint location
at dp+1 . Reset the time t = 0 and find the equation of the motion of the membrane
156

¯
after constraint application Zc (rp+1 , θp+1 , t) which is presented in Eq. (6.46) with
replacing (r, θ) by (rp+1 , θp+1 ). The procedure presented in Section 6.6.1 is used to
obtain the model coefficients at constraint location dp+1 , except the initial conditions
are considered at t = 0 instead of t = tc as follows
2π r(θp+1 ) (p)
1
¯
Z (ρ, ψ, 0) (ϕc )(p+1) rp+1 drp+1 dθp+1
¯i
f
(p+1) 0
a
κc
2π r(θp+1 ) (p)
1
¯
¯(p+1) =
Z (ρ, ψ, 0) (ϕs )(p+1) rp+1 drp+1 dθp+1
¯j
A2j
f
(p+1) 0
a
κs
¯ (p)
2π r(θp+1 ) ∂ Zf (ρ, ψ, 0)
1
(p+1)
¯
=
B1i
(ϕc )(p+1) rp+1 drp+1 dθp+1
¯i
Σc 0
∂t
a
¯ (p)
2π r(θp+1 ) ∂ Zf (ρ, ψ, 0)
1
(p+1)
¯
(ϕs )(p+1) rp+1 drp+1 dθp+1
¯j
=
B2j
Σs 0
∂t
a
¯(p+1) =
A1i

(p+1)
(p+1)
s
c
Where Σs = (ωj )(p+1) κs
and Σc = (ωi )(p+1) κc
. The subscript (p + 1)
indicates the calculations of mode shapes and eigenfrequencies are performed when the
constraint at dp+1 .
A7. Update last equation of the constrained membrane motion to be initial motion for the
¯
¯
next constraint location by setting Zc (rp , θp , t) = Zc (rp+1 , θp+1 , t) and define the
new constraint location dp+2 and go to step 4. Note that for each constraint location,
the integrations require the following relationship

ρ=

2
2
d2 + rp + 2 rp d2 cos θp ,
p
p

ψ = cos−1

dp + rp cos θp
ρ

A flowchart for the procedure described in this section is shown in Fig. 6.9.

157

t = t + ∆t
No
¯
Ω
Zo ≃ 0 ?

simulate Zo (ρ, ψ, t)

Start

t=0

Yes
tc = t

set initial conditions
˙
Zo (ρ, ψ, tc ) and Zo (ρ, ψ, tc )
(p)
(p+1)
Zc = Zc

calculate model coefficient
of constrained membrane
(p+1)
⇒ Zc
(rp+1 , θp+1 , t)
calculate model coefficient
of unconstrained membrane
=⇒ Zf (ρ, ψ, t)
remove constraint

calculate model coefficient
of constrained membrane
=⇒ Zc (r, θ, t)

apply constraint
p=p+1

apply constraint

t=0

(p)
simulate Zc

set initial conditions
˙
Zf (ρ, ψ, 0) and Zf (ρ, ψ, 0)
set initial conditions
(p)
˙ (p)
Zc (rp , θp , tr ) and Zc (rp , θp , tr )

calculate E%

tr = t

Yes

¯
Ω
Zc ≃ 0 ?

No
t = t + ∆t

Figure 6.9: Flowchart for radial application of constraint strategy, for each p increment, new constraint location dp is defined
and the eigenfrequencies and mode shapes are calculates at this location

158

6.6.5

Circumferential Constraint Application

In this section we study vibration suppression in the membrane due to circumferential motion
of the constraint. The constraint is first applied on the X-axis near the membrane boundary
at a distance d from the center of the membrane and then the constraint is moved along
a circular path around the membrane center at location of (d, φ1 ) from the center of the
membrane, where φ1 > 0. Thereafter, the constraint is moved further to new location on the
circular path at (d, φ2 ) such that φ2 > φ1 till the constraint reaches a prescribed location at
(d, φp ) such that φp > · · · > φ2 > φ1 as shown in Fig. 6.11. When the constraint location
is at (d, φ1 ), the eigenfrequencies and the mode shapes of symmetric and antisymmetric
modes should be computed according to the inclined coordinates system (x1 , y1 ) as shown
Fig. 6.11.
Similar to the case of radial constraint motion, the case of circumferential constraint motion
is more complicated than the simple constraint application and removal method. When the
constraint is moved from location (d, 0) to the new location (d, φ1 ),the coordinate system
changes as shown in Fig. 6.11. As a result and similar to radial constraint motion, computing
of modal coefficients of the new location requires integration over the entire membrane area
¯
¯
¯
while excluding the two constraints areas Ω0 and Ω1 . That is because Ω0 represents an
area of zero initial displacement and velocity for the new constraint location at (d, φ1 ) and
¯
the upcoming motion is restricted to the area outside Ω1 . Similar to the radial motion case,
we use an alternative method by storing the initial conditions of the previous motion (d, 0
location with r0 , θ0 coordinates system) in the (ρ, ψ) coordinates system then transfer it to
the new constraint location ( d, φ1 location with r1 , θ1 coordinates system).

159

Z
Y

(a)
vibrating membrane
Y
X
set of sensors

path of constraint motion

application
constraint

(b)
Figure 6.10: Sensing mechanism for circumferential constraint motion (a) back side view. (b) top view
160

Y

y1
P
r1

x1

θ1

ρ
¯
Ω1

R
ψ

d

r0
θ0

φ1

X
d

¯
Ω0

constraint motion

Figure 6.11: The coordinates systems of two different locations of a constraint moved in a circular path from location (d, 0)
to (d, φ1 )
161

The procedure of the constraint movement can be summarized as follows
A1. The membrane is initially vibrates freely and the equation of motion Zo (ρ, ψ, t) is
presented in Eq. (6.44). The initial energy E0 of the membrane is computed by using
Eq. (6.16).
(0)
¯
A2. When the displacement of sensed area Ω0 is ≈ 0, then t = tc and the constraint is
applied.
A3. Reset the time t = 0 and find the equation of the motion of the membrane after
¯
constraint application Zc (r0 , θ0 , t) which is presented in Eq. (6.46) with replacing (r, θ)
by (r0 , θ0 ). The procedure presented in Section 6.6.1 is used to obtain the model
coefficients at constraint location (d, φ0 )=(d, 0) as follows
¯(0)
A1i =

2π r(θ0 )
1
(0)
Zo (ρ, ψ, tc ) (ϕc )(0) r0 dr0 dθ0
¯i
(0) 0
a
κc

¯(0)
A2j =

2π r(θ0 )
1
(0)
¯j
Zo (ρ, ψ, tc ) (ϕs )(0) r0 dr0 dθ0
(0) 0
a
κs

¯ (0)
B1i =

¯ (0)
B2j =

1

2π

(0) 0
c
(ωi )(0) κc
1

2π

(0) 0
s
(ωj )(0) κs

r(θ0 ) ∂Zo (ρ, ψ, t(0) )
c
(ϕc )(0) r0 dr0 dθ0
¯i
∂t
a
r(θ0 ) ∂Zo (ρ, ψ, t(0) )
c
(ϕs )(0) r0 dr0 dθ0
¯j
∂t
a

Where the subscript (0) indicates that the calculations are performed at the location
(d, φ0 ) where φ0 = 0. Set the counter p = 1.
¯ (p)
A4. Compute the total membrane energy after constraint application Ec by using Eq. (6.49)
162

and then compute the energy reduction by using Eq. (6.50)
(p)
¯
A5. When the displacement of sensed area Ωp is ≈ 0, the constraint is removed and t = tr .
and the initial conditions are stored in (ρ, ψ) coordinates system by setting the time
¯ (p)
t = 0, and finding the equation of motion after constraint removal Z (ρ, ψ, t) which
f
is presented in Eq. (6.29). By using similar procedure to that presented in section
6.4.1, the modal coefficients of the unconstrained motion can be obtained as
A6. Transfer the initial conditions stored in (ρ, ψ) coordinates system (from step 5) to
(rp+1 , θp+1 ) which represents the coordinates system of the new constraint location
at (d, φp+1 ).
A7. Reset the time t = 0 and find the equation of the motion of the membrane after
¯
constraint application Zc (rp+1 , θp+1 , t) which is presented in Eq. (6.46) with replacing
(r, θ) by (rp+1 , θp+1 ). The procedure presented in Section 6.6.1 is used to obtain
the model coefficients at constraint location d, φp+1, except the initial conditions are
considered at t = 0 instead of t = tc as follows
2π r(θp+1 ) (p)
1
¯
Z (ρ, ψ, 0) (ϕc )(p+1) rp+1 drp+1 dθp+1
¯i
f
(p+1) 0
a
κc
2π r(θp+1 ) (p)
1
¯
¯(p+1) =
Z (ρ, ψ, 0) (ϕs )(p+1) rp+1 drp+1 dθp+1
¯j
A2j
f
(p+1) 0
a
κs
¯ (p)
2π r(θp+1 ) ∂ Zf (ρ, ψ, 0)
1
(p+1)
¯
=
B1i
(ϕc )(p+1) rp+1 drp+1 dθp+1
¯i
Σc 0
∂t
a
¯ (p)
2π r(θp+1 ) ∂ Zf (ρ, ψ, 0)
1
(p+1)
¯
(ϕs )(p+1) rp+1 drp+1 dθp+1
¯j
=
B2j
Σs 0
∂t
a
¯(p+1) =
A1i

(p+1)
(p+1)
s
c
where Σs = (ωj )(p+1) κs
and Σc = (ωi )(p+1) κc
.
163

A8. Update last equation of the constrained membrane motion to be initial motion for the
¯
¯
next constraint location by setting Zc (rp , θp , t) = Zc (rp+1 , θp+1 , t) and define the
new constraint location (d, φp+1 ) and go to step 4. A flowchart for the procedure
described in this section is shown in Fig 6.12. Note that for each constraint location,
the integrations require the following relationship

ρ=

6.6.6

ψ = cos−1

2
2
d2 + rp + 2 rp d2 cos θp ,

d + rp cos θp
ρ

+ φp

Numerical simulations

In this section we present numerical simulations for the three control strategies presented in
Sections (6.6.3), (6.6.4) and (6.6.5). For all of the three control strategies we assume that
the membrane is initially vibrating in its fundamental mode and the equation of motion is
Z0 (ρ, ψ, t) = A0 J0 (k00 ρ) cos(ω00 )t, where, k00 = 2.4048 m −1 , A0 = 0.1 m and ω00 =
4.8096 rad/sec and the specifications of the membrane and constraint are :

a
= 0.12,
R

R = 1 m,

µ = 0.25 kg/m2 ,

T = 1 N/m

We simulate 10 cycles of constraint application and removal for all of the three control
strategies with 16 symmetric and 16 antisymmetric modes. First, we simulate one cycle of
constraint application and removal for the simple application and removal of the constraint
strategy presented in Section 6.6.3 with d/R = 0.88, i.e the constraint is adjacent to the
¯
boundary. The constraint is applied when the sensed area Ωp passes through the mean
¯
position and removed when Ωp at its maximum amplitude.

164

t = t + ∆t
No
¯
Ω
Zo ≃ 0 ?

simulate Zo (ρ, ψ, t)

Start

t=0

Yes
tc = t

set initial conditions
˙
Zo (ρ, ψ, tc ) and Zo (ρ, ψ, tc )
(p)
(p+1)
Zc = Zc

calculate model coefficient
of constrained membrane
(p+1)
⇒ Zc
(rp+1 , θp+1 , t)
calculate model coefficient
of unconstrained membrane
=⇒ Zf (ρ, ψ, t)
remove constraint

calculate model coefficient
of constrained membrane
=⇒ Zc (r, θ, t)

apply constraint
p=p+1

apply constraint

t=0

(p)
simulate Zc

set initial conditions
˙
Zf (ρ, ψ, 0) and Zf (ρ, ψ, 0)
set initial conditions
(p)
˙ (p)
Zc (rp , θp , tr ) and Zc (rp , θp , tr )

calculate E%

tr = t

Yes

¯
Ω
Zc ≃ 0 ?

No
t = t + ∆t

Figure 6.12: Flowchart for circumferential application of constraint strategy, for each p increment, new constraint location
(angle φp ) is defined and the eigenfrequencies and mode shapes are calculates at this location

165

In the redial constraint application strategy discussed in Section 6.6.4, the constraint is
moved radially forward (toward the membrane center) and backward (toward the outer
boundary). At the beginning, the constraint is applied at the earliest instant of time adjacent
to the outer membrane boundary i.e. at a distance d/R = 0.88 when the forward sensors
record approximately zero displacements and the constrained membrane continues vibrating
freely. Next earliest instant of time when the forward sensors record zero displacements, the
constraint is moved forward to a distance d/R = 0.66 (the previous location is adjacent to
the next location) and process is repeated till the constraint reaches a distance d/R = 0.44
and then the constraint is moved backward gradually till it reaches the outer boundary
at a distance d/R = 0.88 in similar way by using the records of backward sensors. The
circumferential constraint application strategy discussed in Section 6.6.5 is simulated with
an angular increments of φ = 2 tan−1 (a/d) with d/R = 0.88 and a/R = 0.12, where the
constraint moves circumferentially adjacent to the outer boundary of the membrane. At
the beginning, the constraint is applied at the earliest instant of time adjacent to the outer
membrane boundary i.e. at a distance d/R = 0.88 then the constraint is moved counterclockwise incrementally, that is the next constraint location is adjacent to the previous
constraint location whereφ = 15.53◦ .

Figure 6.14 shows displacement and percentage energy loss E% of the membrane at the end of
the first cycle, 6 cycles and after completion of 10 cycles of constraint application and removal
for the three control strategies. Table 6.2 shows comparison between the three constraint
application strategies. It is clear form the data that the simple constraint application and
removal strategy results in significantly higher reduction in energy in comparison to radial
and circumferential constraint application. Although constraint removal does not change
166

the overall energy of the membrane, the configuration of the membrane at the time of
constraint removal determines how energy of the constrained membrane is redistributed over
the entire membrane area; and this has a significant effect on the rate of energy reduction in
subsequent cycles of constraint application. Figure 6.13 shows the snapshots of the first cycle
of constraint application and removal. It can be seen clearly that the removal of constraint as
the the membrane in its maximum displacement, generates a bump that implies excitation of
higher modes similar to our observation in the case of string vibration presented in Chapter
3. Exciting higher modes implies increasing of the velocities distribution over some small
areas that is implicitly included on the area of the membrane. The energy lost is proportional
to the kinetic energy of constrained area Ω, and the kinetic energy is proportional to the
velocity of Ω, therefore, exciting higher modes increases the velocity over the area Ω which
implies increasing in the kinetic energy and, hence, the energy lost is higher.
Table 6.2: Comparison of percentage energy loss of the membrane for the three different
constraint applications strategies (A) simple (B) radial (C) circumferential, over 10 constraint
applications
No.
(A)
(B)
(C)

1
0.305
0.305
0.305

2
2.247
1.203
0.498

3
4.425
2.236
0.606

4
6.000
2.945
0.844

5
8.110
3.289
0.867

6
10.278
5.479
1.251

7
12.102
6.692
1.158

8
14.035
7.182
1.295

9
16.834
9.292
1.3635

10
17.180
11.548
1.483

The energy reduction in radial constraint application is less than the simple application
strategy but higher than the circumferential motion strategy. This can be explained as
follows; in the radial and circumferential motion, the removal time tr can not be chosen
arbitrary like the simple application and removal case. Rather, in radial and circumferential
constraint motion the constraint is moved when the sensor record zero displacement on the
adjacent area, in fact each constraint motion step in these two cases represents application to
167

the new constraint location and removal for the previous constraint location simultaneously.
Hence, exciting higher modes is not guaranteed in these two cases, as a result the total
energy reduction will depend on the kinetic energy stored in constrained part. Furthermore,
the velocity over the area Ω increases inherently when the constraint approaches toward
the center of the membrane. Hence, in the radial constraint application case the energy is
reduced due to velocity increasing over the constrained area. In the circumferential motion
case the amount of kinetic energy removed by the constraint applications is small for two
reasons. First, the velocity over the area to be constrained Ω near the boundary is inherently
low, as a result, the amount of the kinetic energy of Ω is small. Second, the choosing of the
constraint removal time tr is not arbitrary like the simple case, i.e the constraint is moved
when the membrane is at the minimum deformation and therefore higher modes are not
excited.

168

Z

Z
t = 0 sec
X

X

Y

constraint application
t = 0.34 sec

t = 0.40 sec

t = 0.50 sec

constraint removal
t = 0.67 sec

t = 0.80 sec

t = 0.90 sec

t = 1.00 sec

constraint application
t = 1.03 sec

Figure 6.13: Snapshots of one cycle for simple constraint application and removal strategy

169

t = 0.0

E%= 0.0 s

Z
Y
X

E%= 0.305

t = 0.34 s

end of first cycle

simple application
E%= 10.278

t = 3.67 s

radial application
E%= 5.479

circumferential application

t = 4.03 s E%= 1.251

t = 4.11 s

end of 6 cycles
E%= 17.180

t = 6.25 s E%= 11.548

t = 5.42 s E%= 1.483

t = 6.10 s

end of 10 cycles
Figure 6.14: Displacement and percentage energy loss E% of the membrane at the end of the
first cycle, 6 cycles and after completion of 10 cycles of constraint application and removal
for cases; simple, radial and circumferential applications of constraint

170

Chapter 7
Conclusion Remarks and Future Work
In this dissertation, we investigated the vibration control of continuous systems, namely
string and circular membrane, using boundary constraint. The problem of vibrating string
subjected to fixed constraint represented by a rigid obstacle located at one boundary is investigated first. The energy dissipation in the string due to wrapping around the obstacle is
investigated and it is shown that the energy dissipation can be increased by changing the orientation of the obstacle, where the obstacle is considered as passive mechanism. Second, the
problem of vibrating string subjected to moving constraint is investigated. The constraint
is represented by a scabbard that is applied and removed at one of the string boundaries.
The effect of application and removal of the constraint on the total energy of the string is
investigated. The results show that the total energy can be reduced based on the time of
constraint application. Based on energy reduction condition, a semi-active control strategy
is presented to suppress the string vibration by applying and removing the constraint sequentially. The semi-active control strategy is verified experimentally and the results from
experiment and simulation matched well.
171

The vibration control of circular membrane using areal constraint is investigated third in
two parts. The first part includes investigating the accuracy of the dynamics of circular
membrane with the presence of eccentric areal constraint for arbitrary initial conditions.
The second part includes investigation of the effect of application and removal of the areal
constraint on the energy of the membrane. The results show that the total energy can be reduced by applying the constraint as the membrane passes through the mean position. Based
on energy reduction condition, three different semi-active control strategies are presented
and efficacy of each strategy is investigated. In the next sections we discuss the entire work
with details.

7.1

Vibration of a String Wrapping and Unwrapping
Around an Obstacle

In Chapter 2 we investigated the problem of a string vibrating against a smooth obstacle.
The obstacle is located at one of the boundaries and the string is assumed to wrap and
unwrap around the obstacle during vibration. Assuming linear behavior of the string, an
analytical model was developed for computing its geometry at each time step by bookkeeping
the energy. The wrapping of the obstacle is modeled by a series of perfectly inelastic collisions
between the obstacle and adjacent segments of the string and unwrapping is assumed to be
energy conserving. The geometry of the string is determined iteratively starting from an
initial configuration where the string is vibrating in a single mode and is not in contact with
the obstacle. The obstacle can be regarded as a passive mechanism for vibration suppression
in which the energy lost during each cycle of oscillation depends on the energy content of
172

the string at the beginning of the cycle. Numerical simulation results are provided for the
string vibrating in different modes for circular- and elliptic-shaped obstacles. The loss of
energy is found to be greater for higher modes of oscillation and for obstacles that induce
greater length of wrapping.
Since the energy dissipation is proportional to the length of string wrapping around the
obstacle, the energy dissipation can be increased by adding another obstacle at the boundary
of the string. Also, the obstacle geometry can have potential non-linear effect on the vibration
of the string. Another non-linearity effect can be presented due to the large amplitude of the
string vibration or due to inhomogeneity of the mass distribution along the string. Further
investigation on the non-linearity effect lies in the scope of our future work.

7.2

Semi-Active Control for String Vibration

The dynamics of a vibrating string subjected to a constraint at one boundary is investigated
in Chapter 3. The constraint is applied by a scabbard that moves a small distance along the
mean position of the string. The scabbard is moved instantaneously such that the position
and velocity of the string outside the scabbard is unaffected immediately after application
of the constraint, whereas the length of the string covered by the scabbard is brought to
rest. The constraint is removed by moving the scabbard back to its original position and the
change in energy of the string is investigated for different values of scabbard travel distance
and time of application of the constraint. Analytical and numerical simulation results are
first provided for the string vibrating in the first mode, and then for a more general case
where the string has arbitrary initial conditions. The results show that the energy of the
string can increase or decrease depending on the time of application of the constraint for
173

a given distance of travel of the scabbard. This provides the opportunity for semi-active
control of vibration of the string through direct physical interaction, using the scabbard
as an actuator. A simple feedback control strategy is proposed and numerical simulation
results are presented. These results indicate that although removal of the constraint does
not change the energy of the string, the effectiveness of the control strategy depends on the
time of removal of the constraint.

In Chapter 4 we presented an experimental verification to the idea of the scabbard-like
actuator method. The control strategy relies on sequential application and removal of a
zero displacement constraint at one point on the string near its boundary. Application of
the constraint results in rapid dissipation of vibration energy of the string segment that lies
between the point of application of the constraint and the proximal boundary. Removal of
the constraint redistributes the energy of the string over its entire length and allows the
cycle of constraint application and removal to be repeated for suppression of vibration. A
linear model of the string with damping was used to simulate the effect of constraint application and removal; the damping in the model was based on experimentally-determined
values. Simulation results showed that a few cycles of constraint application and removal
can suppress the vibration of the string rapidly. This result was corroborated by experiments
using an extension coil spring which behaves like a lightly damped string. The displacement
of the string was sensed using an optical sensor and the constraint was applied by a pair of
solenoids. Both simulation and experimental results establish the efficacy of our simple control strategy. This simple semi-active control strategies can be used for vibration suppression
of strings and tendons used in space applications such as tensegrity structures where rigid
pre-compressed members (usually bars) are connected by pre-tensioned cables or strings or
174

tendons [81].

7.3

Dynamics and Vibration Control of Membrane

Chapter 5 presented investigation on the dynamics of a circular membrane with an internal
eccentric circular areal constraint. The membrane is assumed to be fixed at its outer boundary and the areal constraint is assumed to impose zero displacement over the entire internal
circular area. This problem has been investigated by several researchers but prior efforts
have been limited to solving for the first few eigenfrequencies; the orthogonality property of
the modes have not been established and the procedure for computing the mode shapes or
simulating the dynamics has not been presented. In this chapter we establish the orthogonality property of the symmetric and antisymmetric modes of the constrained membrane
and provide a systematic method for computing the eigenfrequencies and mode shapes accurately. It is shown that the number of terms in the series expansion plays a critical role
in accurate computation of the mode shapes. While fewer terms result in truncation errors,
too many terms lead to numerical errors due to high sensitivity of the Bessel functions. An
algorithm is presented to choose the appropriate number of terms and the importance of the
algorithm is demonstrated through examples where fewer or more terms result in accurate
eigenfrequencies but inaccurate mode shapes. The orthogonality conditions established in
this chapter are used to determine modal coefficients for dynamics simulations. This however requires proper choice of the number of numerical integration points. An algorithm
is presented to determine the appropriate number of numerical integration points which is
shown to depend on the size and location of the areal constraint. Using the two algorithms
for determination of the mode shapes and computation of the modal coefficients, the dy175

namics of a constrained membrane is accurately simulated; the accuracy of the simulations
can be verified from the boundary conditions. The tools developed for simulation of constrained membrane dynamics, presented in this chapter, can be used to study the vibration
of membrane structures and the efficacy of control methodologies for vibration mitigation.
The energetics of circular membrane subjected to eccentric areal constraint is investigated in Chapter 6 . The constraint of radius a is assumed to be applied at a distance d
form the center of the membrane arbitrary time to during the membrane motion. Under the
assumptions of small displacements and linearity, it was shown that application of the constraint can increase or decrease the total energy of the membrane depending on the time of
application of the constraint. It was shown that the maximum energy reduction occurs when
the constraint is applied as the membrane passes through the mean position. Furthermore,
larger constraint size causes higher increasing/decreasing in total energy, and higher value of
eccentricity d leads to less energy reduction and vice versa. Cyclic application and removal
of constraint will result in vibration suppression when the constraint is always applied at
times when it removes energy from the membrane.
Since the problem of vibrating membrane is two dimensional problem, the constraint can
be moved in many different directions. Three different constraint movement directions are
investigated for the optimal vibration suppression for the membrane. The three directions
of constraint movement are ; (1) simple application and removal of constraint (2) radial
constraint motion (3) circumferential constraint motion. Among the three constraint motion strategies, the simple application and removal of constraint is the best strategy where
the energy is reduced significantly with comparison to radial and circumferential constraint
motion. Radial constraint motion reduces energy more than the circumferential constraint
176

motion. The methodology can be used for vibration control of membranes and other continuous structures through physical interaction with a relatively small section of the structure
near the boundary. The control mechanisms of our control strategies for the vibrating strings
and membranes are simple and do not require complicated feedback control design, thus, our
semi-active control strategies can be used in many different applications such as large membrane mirror and antennas used in space application. The semi-active control strategies
presented in this work can be modified to be used in some non-linear problems in vibration
of continuous systems and such studies are in the scope of out future work.

177

APPENDICES

178

Appendix A

Vibration of Circular Membrane:
Preliminary

A.1

Bessel Functions

One of the varieties of special functions which are encountered in the solution of physical
problems is the class of functions called Bessel functions. They are solutions to a very
important differential equation, the Bessel equation [80] :

d2 y
dy
x2
+ x + (x2 − n2 ) = 0
2
dx
dx

(A.1)

Where n could be any integer (the order of the Bessel function). The Bessel function
could have the following forms:

179

Bessel functions of the first kind Jn (x)
∞
Jn (x) =
m=0
Jn (x) =

(−1)m
1
( x)2m+n
m!Γ[m + n + 1] 2

1 π
cos(nτ − x sin(τ ))dτ.
π 0

(A.2)

(A.3)

Bessel functions of the second kind Yn (x)

Yv (x) =

Jv (x) cos(vπ) − J−v (x)
sin(vπ)

(A.4)

Where v is not an integer, and for v is an integer n we have

1
( x)−n n−1 (n − k − 1)! 1
2
1
( x2 )k + ln( x)Jn (x)
Yn (x) = − 2
π
k!
4
π
2
k=0
1
1
(− x2 )k
( x)n n−1
[ϕo (k + 1) + ϕo (n + k + 1)] 4
− 2
π
k!(n + k)!
k=0

(A.5)

Where ϕo is the digamma function.
The first and second kind Bessel function have the following properties

J−n (x) = (−1)n Jn (x)

(A.6)

Y−n (x) = (−1)n Yn (x)

(A.7)

180

1

J0

J1

J2

J3

J5

Jn (x)

J15

J10

0.5

0.0

−0.5
0

12

8

4

x

Figure A.1: Bessel functions of the first kind Jn (x)

181

16

20

1
Y0

Yn (x)

0.0

Y1
Y2

−1

Y15
Y3

Y5

−2
0

12

8

4

x
Figure A.2: Bessel functions of the second kind Yn (x)

182

16

20

A.2

Vibration of Circular Membranes

Consider a circular membrane of the radius R as shown in Fig. A.3

y

P
ρ
R
ψ
x

Figure A.3: Circular membrane with radius R .

The governing differential equation for the vibrating membrane is given as [75]
1 ∂2Z
1 ∂2Z
∂ 2 Z 1 ∂Z
+
=
+
ρ ∂ρ
∂ρ2
ρ2 ∂ψ 2
c2 ∂t2

(A.8)

Where Z = Z(ρ, ψ, t) is the transversal displacement (perpendicular to the page in Fig. A.3,
and c2 = T /µ. Where T is the tension in the membrane and µ is the mass density per unit
area.
General Solution
The solution for Eq. (A.8) can be obtained by applying separation of the variables.
183

∞
Z(ρ, ψ, t) =

{βm Jm (kρ) + γm Ym (kρ)} {Cm cos mψ + Dm sin mψ} {A cos ωt + B sin ωt}
m=0

(A.9)

Boundary conditions
The boundary conditions are associated with the spatial variables r and ψ, therefore, we
apply the boundary conditions to V (ρ, ψ) only.

A1. At r = 0, V (0, ψ) = {βJm (0) + γYm (0)}{Cm cos mψ + Dm sin mψ}, but Ym (0) = ∞
and it is not suitable because there is finite value for V (0, ψ), which requires that γ = 0.

A2. The membrane is fixed at r = R or V (R, ψ) = {βJm (kR)}{Cm cos mψ+Dm sin mψ} =
0
Now, the only term left to satisfy V (R, ψ) = 0 is the zeros of the Bessel function
Jm (kR). There will be many zeros (intersections with x-axis), let the value of the first
intersection of the Jm (γ) occurs at γm1 and the second intersection at γm2 and so
γ
on. Then knm R = γnm =⇒ knm = nm .
R
Then
∞

∞

V (ρ, ψ) =
m=0 n=1

ρ
Jm (γmn ){Cmn cos mψ + Dmn sin mψ}
R

(A.10)

and the form of Z(ρ, ψ, t) becomes
∞

∞

Z(ρ, ψ, t) =
m=0 n=1

ρ
c
s
Jm (γmn ){Tmn (t) cos mψ + Tmn (t) sin mψ}
R

(A.11)

where
c
Tmn (t) = Amn cos ωmn t + Bmn cos ωmn t
184

(A.12)

and
s
Tmn (t) = Cmn cos ωmn t + Dmn cos ωmn t

(A.13)

The coefficients Amn , Bmn ,Cmn and Dmn are to be determined by apply the initial conr
ditions and ωmn = cγmn . Equation (A.11) represents the eigenfunction of the membrane,
R
r
and the eigen frequencies ωnm in Eq. (A.12) can be obtained by ωnm = cknm = cγmn ,
R
where γmn can be obtained form the zeros of the Bessel function and table (A.1) shows
some values for γmn .
Figure A.4 shows some different possible shape modes. And it is clear that there are angular
Table A.1: Some values for γmn which is zeros for Bessel function Jm (γnm ).

\

n m
1
2
3
4

0
2.4048
5.5201
8.6537
11.7915

1
3.8317
7.0156
10.1735
13.3237

2
5.1356
8.4172
11.6198
14.7960

3
6.3802
9.7610
13.0152
16.2235

4
7.5883
11.0647
14.3725
17.6160

and redial nodes, where m refers to number of angular nodes and n refers to number of redial
nodes.

185

(n, m) = (1, 0)

(n, m) = (1, 1)

(n, m) = (2, 0)

(n, m) = (2, 1)

(n, m) = (3, 0)

(n, m) = (3, 1)

(n, m) = (4, 0)

(n, m) = (4, 1)

(a)

(b)

Figure A.4: Fundamental modes of circular membrane with (a) First four redial nodes and
no angular nodes (b) First four redial nodes and one angular nodes .

186

(n, m) = (1, 2)

(n, m) = (1, 3)

(n, m) = (2, 2)

(n, m) = (2, 3)

(n, m) = (3, 2)

(n, m) = (3, 3)

(n, m) = (4, 2)

(n, m) = (4, 3)

(a)

(b)

Figure A.5: Fundamental modes of circular membrane with (a) First four redial nodes and
two angular nodes (b) First four redial nodes and three angular nodes .

187

Appendix B

Change of Summation Indices

We have the following double summation
∞

∞

m=0 q=−∞

β2m Jm+q (kR) + γ2m Ym+q (kR) Jq (kd) cos(m + q)ψ = 0

(B.1)

we define n = m + q, the

∞

∞
{β2m Jn (kR) + γ2m Yn (kR)} Jn−m (kd) cos(n)ψ = 0

(B.2)

m=0 n=−∞

note that
∞

∞
Hn =

n=−∞

−1
Hn +

n=0

∞
Hn =

n=−∞
188

ǫn (Hn + H−n )
n=0

(B.3)

where ǫn is defined as
ǫn =

using Eq. (B.5) in Eq. (B.2) leads to
∞



 1/2 if n = 0

 1

if n > 0

∞
{β2m Jn (kR) + γ2m Yn (kR)} Jn−m (kd) cos(n)ψ

m=0 n=0
+ β2m J−n (kR) + γ2m Y−n (kR) J−n−m (kd) cos(−n)ψ

(B.4)

From the properties of Bessel functions, we have

J−ν (x) = (−1)ν Jν (x),

Y−ν (x) = (−1)ν Yν (x)

Then Eq. (B.1) finally becomes
∞

∞

ǫn {β2m Jn (kR) + γ2m Yn (kR)} Jn−m (kd) + (−1)m Jn+m (kd) cos nψ = 0

m=0 n=0
(B.5)

For the double summation with sine term, we can use the same procedure to prove that
∞

∞

m=1 q=−∞

β1m Jm+q (kR) + γ1m Ym+q (kR) Jq (kd) sin(m + q)ψ = 0

(B.6)

can be written as
∞

∞

{β1m Jn (kR) + γ1m Yn (kR)} Jn−m (kd) − (−1)m Jn+m (kd) sin nψ = 0

m=1 n=1
(B.7)
189

Appendix C

Energy of Circular Membrane

The total energy of the membrane in Cartesian coordinates is obtained as [82]

µ
2

E=

∂Z(x, y, t) 2 T
+
∂t
2

∂Z(x, y, t) 2 T
+
∂x
2

∂Z(x, y, t) 2
∂y

dx dy

(C.1)

A
We have

∂z
∂z ∂r ∂z
=
+
∂x
∂r ∂x ∂θ
∂z ∂r ∂z
∂z
=
+
∂y
∂r ∂y ∂θ

where r =

∂θ
∂z
∂z sin θ
=
cos θ −
∂x
∂r
∂θ r
∂z
∂z cos θ
∂θ
=
sin θ +
∂y
∂r
∂θ r

x2 + y 2 and

∂r
= cos x,
∂x

∂r
= cos y,
∂y

∂θ
sin θ
=−
,
∂x
r
190

∂θ
cos θ
=
∂y
r

(C.2)

By squaring the equations in Eq. (C.1), we get
∂z 2
=
∂x
∂z 2
=
∂y

∂z 2
1 ∂z ∂z
∂z 2 sin2 θ
cos2 θ −
cos θ sin θ +
∂r
r ∂r ∂θ
∂θ
r2
∂z 2
1 ∂z ∂z
∂z 2 cos2 θ
sin2 θ +
cos θ sin θ +
∂r
r ∂r ∂θ
∂θ
r2

(C.3)

By adding the equations in Eq. (C.3) we get
∂z 2
+
∂x

∂z 2
=
∂y

∂z 2
1
+
∂r
r2

∂z 2
∂θ

Note that dx dy = r dr dθ, then the energy expression in Eq.(C.1) becomes

E =

∂Z(r, θ, t) 2
r dr dθ
∂t

µ
2
A

+

T
2
A

1
∂Z(r, θ, t) 2
+
∂r
r2

191

∂Z(r, θ, t) 2
∂θ

r dr dθ

(C.4)

BIBLIOGRAPHY

192

BIBLIOGRAPHY
[1] H. Yee and N. Audeh, “Cutoff frequencies of eccentric waveguides,” IEEE Transactions
on Microwave Theory and Techniques, vol. 14, no. 10, pp. 487 – 493, 1966.
[2] N. K., “Vibration of a membrane having a circular outer boundary and an eccentric
circular inner boundary,” Journal of Sound and Vibration, vol. 50, no. 4, pp. 545 – 551,
1977.
[3] V. A. S. G. P. Singh, G. S. and L. S. Kothari, “Comments on free transverse vibrations
of uniform circular plates and membranes with eccentric holes,” Journal of Sound and
Vibration, vol. 100, no. 1, pp. 141–145, 1985.
[4] “The spartan 207 inflatable antenna http://www.nasaimages.org/.”
[5] C. Citrini, “The energy theorem in the impact of a string vibrating against a pointshaped
obstacle,” Atti Accad. Naz. Lincei Rend., vol. 59, pp. 143–148, 1977.
[6] L. Amerio, “Continuous solutions of the problem of a string vibrating against an obstacle,” Rendiconti del Seminario Mathematico della Universita di Padova, vol. 59, no. ,
pp. 67–96, 1978.
[7] Michelle and Schatzman, “A hyperbolic problem of second order with unilateral constraints: The vibrating string with a concave obstacle,” Journal of Mathematical Analysis and Applications, vol. 73, no. 1, pp. 138 – 191, 1980.
[8] A. Haraux and H. Cabannes, “Almost periodic motion of a string vibrating against a
straight fixed obstacle,” Nonlinear Analysis, Theory, Methods and Applications, vol. 7,
no. 2, pp. 129–141, 1983.
[9] R. Burridge, J. Kappraff, and C. Morshedi, “The sitar string, a vibrating string with a
one-sided inelastic constraint,” SIAM Journal on Applied Mathematics, vol. 42, no. 6,
pp. 1231–1251, 1982.
193

[10] A. Bamberger and M. Schatzman, “New results on the vibrating string with a continuous
obstacle,” SIAM Journal of Mathematical Analysis, vol. 14, no. 3, pp. 560–595, 1983.
[11] J. Ahn, “A vibrating string with frictionless dynamic impact,” Applied Numerical Mathematics, vol. 57, pp. 861–884, 2007.
[12] M. Frontini and L. Gotusso, “Numerical study of the motion of a string vibrating
against an obstacle by physical discretization,” Applied Mathematical Modeling, vol. 14,
pp. 489–494, 1990.
[13] A. Krishnaswamy and J. Smith, “Methods for simulating string collisions with rigid
spatial obstacles,” in IEEE Workshop on Applications of Signal Processing to Audio
and Acoustics, pp. 233 – 236, 2003.
[14] S. M. Han and M. A. Grosenbaugh, “Nonlinear free vibration of a cable against a
straight obstacle,” Journal of Sound and Vibration, vol. 273, pp. 337–361, 2004.
[15] D. N. Dutt and B. S. Ramakrishna, “Vibration control of stretched strings by an external
distributed force,” Journal of Sound and Vibration, vol. 57, no. 4, pp. 603–606, 1978.
[16] D. N. Dutt and B. S. Ramakrishna, “Application of optimization methods to vibration
control of stretched strings,” Journal of Sound and Vibration, vol. 57, no. 4, pp. 499–514,
1978.
[17] J. Issa, R. Mukherjee, A. R. Diaz, and S. W. Shaw, “Modal disparity and its experimental verification,” Journal of Sound and Vibration, vol. 311, no. 35, pp. 1465 – 1475,
2008.
[18] J. Issa, R. Mukherjee, and A. R. Diaz, “Energy dissipation in dynamical systems through
sequential application and removal of constraints,” Journal of Dynamic Systems, Measurement, and Control, vol. 131, no. 2, p. 021011, 2009.
[19] P. Hbrard and A. Henrott, “Optimal shape and position of the actuators for the stabilization of a string,” Systems and Control Letters, vol. 48, no. 34, pp. 199 – 209,
2003.
[20] A. D. Myshkis, “Vibrations of the string with energy dissipation and impulsive feedback
support,” Nonlinear Analysis, vol. 26, no. 7, pp. 1271 – 1278, 1996.
194

[21] W. D. Zhu and N. A. Zheng, “Exact response of a translating string with arbitrarily varying length under general excitation,” Applied Mathematics and Computation,
vol. 75, no. 3, 2008.
[22] S. Shahruz and D. Kurmaji, “Vibration suppression of a nonlinear axially moving string
by boundary control,” Journal of Sound and Vibration, vol. 201, no. 1, pp. 145 – 152,
1997.
[23] R.-F. Fung, J.-W. Wu, and S.-L. Wu, “Exponential stabilization of an axially moving
string by linear boundary feedback,” Automatica, vol. 35, no. 1, pp. 177 – 181, 1999.
[24] T. cheng Li, Z. chao Hou, and J. feng Li, “Stabilization analysis of a generalized nonlinear axially moving string by boundary velocity feedback,” Automatica, vol. 44, no. 2,
pp. 498 – 503, 2008.
[25] K. Do and J. Pan, “Boundary control of transverse motion of marine risers with actuator
dynamics,” Journal of Sound and Vibration, vol. 318, no. 45, pp. 768 – 791, 2008.
[26] H. K. Khalil, Nonlinear Systems, 3rd Ed. Prentice Hall, New York, NY., 2002.
[27] F. ZHANG, D. DAWSON, S. NAGARKATTI, and C. RAHN, “Boundary control for
a general class of non-linear string-actuator systems,” Journal of Sound and Vibration,
vol. 229, no. 1, pp. 113 – 132, 2000.
[28] J. Lu, M. Crocker, and P. Raju, “Active vibration control using wave control concepts,”
Journal of Sound and Vibration, vol. 134, no. 2, pp. 364 – 368, 1989.
[29] A. D. Myshkis, “Active wave control of string near fixed boundary,” ASME Journal of
Vibration and Acoustics, vol. 130, no. 1, 2008.
[30] C. TAN and S. YING, “Active wave control of the axially moving string: Theory and
experiment,” Journal of Sound and Vibration, vol. 236, no. 5, pp. 861 – 880, 2000.
[31] W. Zhang and L.-Q. Chen, “Vibration control of an axially moving string system: Wave
cancellation method,” Applied Mathematics and Computation, vol. 175, no. 1, pp. 851
– 863, 2006.
[32] N. Tanaka and Y. Kikushima, “Active wave control of a flexible beam fundamental
characteristics of an active-sink system and its verification,” Journal of Sound and
Vibration, vol. 35, no. 2, pp. 499–514, 1992.
195

[33] C. Baicu, C. Rahn, and B. Nibali, “Active boundary control of elastic cables: Theory
and experiment,” Journal of Sound and Vibration, vol. 198, no. 1, pp. 17 – 26, 1996.
[34] M. Sosa, A. Carneiro, J. Colafemina, and O. Baffa, “A new magnetic probe to study the
vibration of the tympanic membrane,” Journal of Magnetism and Magnetic Materials,
vol. 226230, no. 2, pp. 2067 – 2069, 2001.
[35] E. J. Ruggiero and D. J. Inman, “Modeling and vibration control of an active membrane
mirror,” Smart Materials and Structures, vol. 18, p. 095027, 2009.
[36] K. Miyano, “Surface tension measured by vibrating membranes: An application to
smectic-a and -b phases,” Phys. Rev. A, vol. 26, pp. 1820–1823, 1982.
[37] G. Percin, T. S. Lundgren, and B. T. Khuri-Yakub, “Controlled ink-jet printing and
deposition of organic polymers and solid particles,” Applied Physics Letters, vol. 73,
pp. 2375–2378, 1998.
[38] H. Y. Yee and N. F. Audeh, “Uniform waveguides with arbitrary cross sections considered by the point-matching method,” Institute of Electrical and Electronic Engineers
Transactions MTT13, pp. 847–851, 1965.
[39] J. Davies and C. Muilwyk, “Numerical solution of uniform hollow waveguides with
boundaries of arbitrary shape,” Proceedings of the Institution of Electrical Engineers,
vol. 113, no. 2, pp. 277 –284, 1966.
[40] C. W and Steele, “Numerical computation of electric and magnetic fields in a uniform
waveguide of arbitrary cross section,” Journal of Computational Physics, vol. 3, no. 1,
pp. 148 – 153, 1968.
[41] P. Arlett, A. Bahrani, and O. Zienkiewicz, “Application of finite elements to the solution
of helmholtz’s equation,” Proceedings of the Institution of Electrical Engineers, vol. 115,
no. 12, pp. 1762 –1766, 1968.
[42] N. Gass, “Solution of the three dimensional scalar helmholtz equation by a finite element
formulation,” ZAMM - Journal of Applied Mathematics and Mechanics / Zeitschrift fr
Angewandte Mathematik und Mechanik, vol. 57, no. 2, pp. 95–99, 1977.
[43] P. Laura, “A simple method for the determination of cutoff frequencies of waveguides
with arbitrary cross sections,” Proceedings of the IEEE, vol. 54, no. 10, pp. 1495 – 1497,
1966.
196

[44] M.J. and Hine, “Eigenvalues for a uniform fluid waveguide with an eccentric-annulus
cross-section,” Journal of Sound and Vibration, vol. 15, no. 3, pp. 295 – 305, 1971.
[45] R. Bulley and B. Davies, “Computation of approximate polynomial solutions to te
modes in an arbitrarily shaped waveguide,” IEEE Transactions on Microwave Theory
and Techniques, vol. 17, no. 8, pp. 440 – 446, 1969.
[46] L. W.H., “Free transverse vibrations of uniform circular plates and membranes with
eccentric holes,” Journal of Sound and Vibration, vol. 81, no. 3, pp. 425 – 435, 1982.
[47] G. S. Singh and L. S. Kothari, “On the solution of the two-dimensional helmholtz
equation,” Journal of Mathematical Physics, vol. 25, no. 3, pp. 810–811, 1984.
[48] K. Nagaya and Y. Hai, “Free vibration of composite membranes with arbitrary shape,”
Journal of Sound and Vibration, vol. 100, no. 1, pp. 123 – 134, 1985.
[49] S. Noga, “Free transverse vibration analysis of an elastically connected annular and
circular double-membrane compound system,” Journal of Sound and Vibration, vol. 329,
no. 9, pp. 1507 – 1522, 2010.
[50] A. Arokiasamy, “Simulation of vibrations of rectangular and square membranes using
computer graphics,” Simulation Modelling Practice and Theory, vol. 10, no. 12, pp. 3 –
12, 2002.
[51] E. Alvarado-Anell, M. Sosa, and M. A. Moreles, “Computational study of forced oscillations in a membrane,” REVISTA MEXICANA DE FISICA, vol. 51, no. 2, pp. 102–107,
2005.
[52] E. Alvarado-Anell, M. Sosa, and M. A. Moreles, “Numerical simulation of the dynamical
properties of the human tympanum,” REVISTA MEXICANA DE FISICA, vol. 54,
no. 2, pp. 135–140, 2008.
[53] A. BHIDE, “Computer simulation of human tympanic membrane,” The Journal of
Laryngology and Otology, vol. 106, pp. 878–881, 1992.
[54] P. A. A. Laura and R. H. Gutierrez, “Determination of the fundamental eigenvalue of
a rhombic membrane with a concentric circular hole,” Journal of Sound and Vibration,
vol. 153, no. 1, pp. 176–179, 1992.
[55] L. H. Wu and C. Y. Wang, “Fundamental frequencies of a circular membrane with a
centered strip,” Journal of Sound and Vibration, vol. 239, no. 2, pp. 363–368, 2001.
197

[56] S. W. Kang and J. M. Lee, “Free vibration analysis of composite rectangular membranes
with an oblique interface,” Journal of Sound and Vibration, vol. 251, no. 3, pp. 505–517,
2002.
[57] J. Chen, S. Kao, W. Lee, and Y. Lee, “Eigensolutions of the helmholtz equation
for a multiply connected domain with circular boundaries using the multipole trefftz
method,” Engineering Analysis with Boundary Elements, vol. 34, no. 5, pp. 463 – 470,
2010.
[58] V. H. Cortinez and P. A. A. Laura, “Vibrations of non-homogeneous rectangular membranes,” Journal of Sound and Vibration, vol. 156, no. 2, pp. 217–225, 1992.
[59] F. F. Cap, “Eigenfrequencies of membranes of arbitrary boundary and with varying
surface mass density,” Applied Mathematics and Computation, vol. 124, pp. 319–329,
2001.
[60] C. Zhao and G. P. Stevens, “A practical error estimator for finite element predicted
natural frequencies of membrane vibration problems,” Journal of Sound and Vibration,
vol. 195, no. 5, pp. 739–756, 1996.
[61] I. P. A. R. Stanford, B. and W. Shyy, “Fixed membrane wings for micro air vehicles: Experimental characterization, numerical modeling, and tailoring,” Progress in Aerospace
Sciences, vol. 44, pp. 258–294, 2008.
[62] T. X. I. E. G. R. B. K. S. S. Song, A. and B. K., “Aeromechanics of membrane wings
with implications for animal flight,” AIAA Journal, vol. 46, no. 8, 2008.
[63] G. M. S. W. Z. Rojratsirikul, P. and I. Gursul, “Flow-induced vibrations of low aspect
ratio rectangular membrane wings,” Journal of Fluids and Structures, vol. 27, no. 8,
pp. 1296–1309, 2011.
[64] D. A. Wang and H. H. Ko, “Piezoelectric energy harvesting from flow-induced vibration,” Journal of Micromechanics and Microengineering, vol. 20, no. 2, 2010.
[65] J. J. Allen and A. J. Smits, “Energy harvesting eel,” Journal of Fluids and Structures,
vol. 15, no. 3, pp. 629–640, 2011.
[66] Q. Wang, “Wave propagation in a piezoelectric coupled cylindrical membrane shell,”
International Journal of Solids and Structures, vol. 38, no. 4647, pp. 8207 – 8218, 2001.
198

[67] C. H. Jenkins and U. A. Korde, “Membrane vibration experiments: An historical review
and recent results,” Journal of Sound and Vibration, vol. 295, no. 35, pp. 602 – 613,
2006.
[68] R. P. Grosso and M. Yellin, “The membrane mirror as an adaptive optical element,”
The Optical Society of America, vol. 67, no. 3, pp. 399–406, 1977.
[69] I. Jo, “On the optimal control of circular membranes,” Acta Mathematica Hungarica,
vol. 58, pp. 365–376, 1991.
[70] I. Jo, “On the control of a rectangular membrane,” Acta Mathematica Hungarica, vol. 59,
pp. 25–32, 1992.
[71] D. M. James, B. G. Patrick, S. Chodimella, D. Marker, and B. deBlonk, “Design
and testing of a one-meter membrane mirror with active boundary control,” vol. 5899,
p. 58990Z, SPIE, 2005.
[72] M. V. Aguanno, F. Lakestani, M. P. Whelan, and M. J. Connelly, “Heterodyne speckle
interferometer for full-field velocity profile measurements of a vibrating membrane by
electronic scanning,” Optics and Lasers in Engineering, vol. 45, no. 6, pp. 677 – 683,
2007.
[73] S. Tliba, “Control of a vibrating axisymmetric membrane using piezoelectric tranducers,” vol. 18, 2011.
[74] J. M. E. Tarazaga, P. A. and D. J. Inman, “Vibro-acoustics of a pressurized optical
membrane,” Mechanical Systems and Signal Processing, vol. 30, pp. 373–392, 2012.
[75] E. Kreyzig, Advanced Engineering Mathematics 9th Ed. John Wiley and Sons, Hoboken,
NJ, 1993.
[76] L. Meirovitch, Elements of Vibration Analysis. McGraw Hill, New York, NY.
[77] E. Kashy, D. A. Johnson, J. McIntyre, and S. L. Wolfe, “Transverse standing waves in
a string with free ends,” American Journal of Physics, vol. 65, no. 4, pp. 310–313, 1997.
[78] A. W. Leissa and M. S. Qatu, Vibration of Continuous Systems, 1st Edition. Mc-Graw
Hill Professional, 2011.
[79] “http://www.drtempleman.com/spring-guides.htm.”
199

[80] G. N. Watson, A Treatise on the Theory of Bessel Functions. Cambridge University
Press, 1958.
[81] I. J. Oppenheim and W. O. Williams, “Vibration of an elastic tensegrity structure,”
European Journal of Mechanics - A/Solids, vol. 20, no. 6, pp. 1023 – 1031, 2001.
[82] L. N. Virgin, Vibration of Axially Loaded Structures. Cambridge University Press, 2007.

200