VIBRATION ISOLATION DESIGN FOR AEROSPACE

STRUCTURES USING MATLAB OPTIMIZATION

A colloquium report submitted by

DEVIKA SANTHOSH
M240360ME

SAREPALLI RAMHARAN M240549ME

In partial fulfilment for award of degree of

Master of Technology
in Machine Design

under the guidance of

Dr. JAGADEESHA

NATIONAL INSTITUTE OF TECHNOLOGY, CALICUT

NIT Campus P.O, Calicut,
Kerala, India, 673601

April 2025

i
 Department of Mechanical Engineering

NATIONAL INSTITUTE OF TECHNOLOGY, CALICUT

Certificate
This is to certify that the project report entitled “Vibration Isolation Design For
Aerospace Structures Using MATLAB Optimization.” submitted by Ms. Devika
Santhosh S (M240360ME)and Mr. Sarepalli Ramharan(M240549ME),National Institute
of Technology Calicut, towards the partial fulfillment of the requirements for the award of the
Degree of Master of Technology in Machine Design, is a bonafide record of the work carried
out by him under my supervision and guidance.

Dr. Jagadeesha T
Associate Professor
Department of Mechanical Engineering
National Institute of Technology, Calicut

ii
 DECLARATION
We hereby declare that this project report on “Vibration Isolation Design For Aerospace
Structures Using MATLAB Optimization.” is our own work and that, to the best of our
knowledge and belief, it contains no material previously published or written by another
person nor material which has been accepted for the award of any other degree or diploma of
the university or other institutes of higher learning, except where due acknowledgement has
been made in the text.

iii
 ACKNOWLEDGMENT
We sincerely express gratitude to all those who contributed to the successful completion of
this project. Guidance and support from mentors and peers were instrumental in navigating
challenges and achieving meaningful results.

We also acknowledge the availability of essential resources and infrastructure, which made
the simulation and analysis possible. Constructive feedback and encouragement from
colleagues played a significant role in refining the work and maintaining focus throughout the
process.

Finally, We deeply thankful for the support and understanding from everyone who provided
motivation and encouragement during this endeavor. This achievement reflects the collective
effort and shared commitment of all involved.

iv
 ABSTRACT
This study investigates the dynamic behavior and damping characteristics of a composite
cantilever beam consisting of a steel base layer and a variable-thickness rubber layer. The
analysis combines analytical modeling and numerical simulation to evaluate how damping
performance evolves with changes in the rubber layer's thickness. Equivalent material
properties—such as density and Young’s modulus—are calculated based on the composite
geometry. Using modal analysis, the first mode frequency is derived from the standard
cantilever beam characteristic equation, and the corresponding mode shape is obtained and
normalized.

To quantify energy dissipation, an effective damping coefficient is computed using a gradient-

based formulation involving the mode shape and rubber viscosity. The beam's timedomain
response is modeled as a single-degree-of-freedom (SDOF) damped oscillator and is solved
using the explicit Euler integration method. Key performance metrics—including peak time,
rise time, settling time, and steady-state error—are extracted from the simulated displacement
responses. A parametric sweep over rubber thickness (0–1 times the steel thickness) is
conducted, and the impact on vibrational behavior is analyzed. Results reveal how increasing
rubber thickness improves damping but affects response time and residual oscillation. These
insights are essential for optimizing layered beam structures for vibration control
applications.

TABLE OF CONTENTS

i
CHAPTER Title Page No

ABSTRACT I

LIST OF FIGURES III

1 INTRODUCTION 1

1.1 PROBLEM STATEMENT 2

1.2 OBJECTIVE 2

1.3 SCOPE OF WORK 2

2 METHODOLOGY 8

3 RESULTS 10

4 CONCLUSION 13

5 REFERENCE 14

6 MATLAB CODE 15

LIST OF FIGURES
Figure No. Title Page No

ii
3.1 Sample response of beam when 0.2mm thickness of 10
rubber is used
3.2 variation of peak time when rubber thickness ratio is 11
increasing
3.3 Peak time vs rubber thickness ratio 11

3.4 settling time vs rubber thickness ratio 12

3.5 steady-state error vs rubber thickness ratio 12

iii
 CHAPTER 1

INTRODUCTION
Vibration control is a critical aspect in the design of mechanical and structural systems,
especially in applications where oscillatory motion can lead to noise, fatigue, or failure. One
effective method for passive vibration suppression involves the integration of viscoelastic
damping layers, such as rubber, with high-stiffness structural materials like steel. In this study,
we investigate the damping characteristics of a composite cantilever beam composed of a steel
base and a rubber layer of varying thickness. The objective is to understand how the rubber
layer influences the system’s dynamic response characteristics, including peak time, rise time,
settling time, and steady-state error.

hr

Where h = thickness of steel bar hr

= thickness of rubber.

1
PROBLEM STATEMENT OBJECTIVE:

To optimize the vibration isolation performance of a composite cantilever beam by analyzing

the effect of rubber layer thickness on its dynamic response using MATLAB simulation.

Scope of work

• To model a composite cantilever beam with a steel core and variable-thickness rubber
layer.
• To compute equivalent mechanical properties for the composite section.
• To simulate the dynamic response using a single-degree-of-freedom (SDOF) model in
MATLAB.
• To evaluate key time-domain characteristics: peak time, rise time, settling time, and
steady-state error.
• To perform a parametric analysis over a range of rubber thicknesses.
• To identify the optimal rubber thickness range for effective vibration [Link]

Properties of beam and rubber used in the code

Beam

Length L = 0.25m

Breadth b = 0.03m

Height h_steel = 0.001m

Youngs modulus E_steel = 193 Gpa.

Density of steel = 8030 (Kg/m^3)

Rubber

Youngs modulus E_rubber = 5 Mpa

Density of rubber = 1200 (Kg/m^3)

Viscosity of rubber = 5000 pa.s

2
THEORETICAL FRAMEWORK

We have considered 1st mode of vibration of cantilever beam

For transverse vibration of a uniform beam according to Euler bernoulli beam theory

Taking only 1st mode w(x,t) = ϕ(x).q(t),

Normalisation of the mode shape,

Multiply the governing PDE by ϕ(x), integrate over the beam length [0,L] and use
orthogonality of modes:

This simplifies to the standard form;;

If damping is considered,

Modal damping is not derived from Euler Bernoulli equation;

3
 1. STRUCTURAL VIBRATION THEORY

The system is idealized as a single-degree-of-freedom (SDOF) damped oscillator to simplify

the time-domain simulation of beam response. The governing differential equation for the
damped SDOF system is:

mÿ(t) + c𝑦̇(t) + ky(t) = 0

Where:

• m: equivalent mass
• c: damping coefficient
• k: equivalent stiffness
• y(t): displacement
• 𝑦̇(t): velocity
• ÿ(t): acceleration

2. NUMERICAL SIMULATION

To facilitate numerical simulation, this second-order ODE is split into two first-order ODEs
and solved using the Euler integration method:

1. Compute acceleration:

2. Update velocity and displacement:

4
 3. COMPOSITE BEAM PROPERTIES

The beam is modelled with a total height ℎ𝑡𝑜𝑡𝑎𝑙 = ℎ𝑠𝑡𝑒𝑒𝑙 + ℎ𝑟𝑢𝑏𝑏𝑒𝑟 . The equivalent properties
for a composite cross-section are computed as:

• Area: A = b⋅ ℎ𝑡𝑜𝑡𝑎𝑙
𝑏ℎ3
• Moment of inertia: I =
12

• Equivalent density:

𝜌𝑒𝑞 = 𝜌𝑠𝑡𝑒𝑒𝑙 ℎ
𝑠𝑡𝑒𝑒𝑙 ℎ
+ 𝜌
𝑡𝑜𝑡𝑎𝑙 𝑟𝑢𝑏𝑏𝑒𝑟 ℎ𝑟𝑢𝑏𝑏𝑒𝑟

• Equivalent modulus:
𝐸𝑒𝑞 = 𝐸𝑠𝑡𝑒𝑒𝑙 𝑠𝑡𝑒𝑒𝑙
ℎ + 𝐸
ℎ 𝑡𝑜𝑡𝑎𝑙 𝑟𝑢𝑏𝑏𝑒𝑟 ℎ𝑟𝑢𝑏𝑏𝑒𝑟

Natural frequency and mode shape

For a cantilever beam, the fundamental frequency is determined by solving the characteristic
equation:

Cos(βL)Cosh(βL) + 1 = 0

Using the first non-trivial root (βL≈1.875), the natural frequency is:

Φ(x) = Cosh(βx)- Cos(βx) - Cosh(βL) +Cos(βL) (Sinh(βx) − Sin(βx)

Sinh(βL) +Sin(βL)

The first mode shape is used to estimate damping and displacement at any point along the
beam.

5
 4. VISCOELASTIC MODELING: KELVIN–VOIGT MODEL

The rubber layer in the beam is modeled as a Kelvin–Voigt material, which represents linear
elasticity in parallel with linear viscosity. The stress-strain relationship is given by:

σ(t)=E ε(t)+η ε̇ (t)

where

• σ(t) = stress
• ε(t) = strain
• ε˙(t) = strain rate
• E = elastic modulus
• η = viscous damping coefficient

This constitutive relation allows modeling of both instantaneous elastic response and
timedependent damping behavior. We also know that, shear strain in rubber, ε(x,t) =

𝛾(𝑥, 𝑡) = 𝜕 ∅

𝜕𝑥

Stress in the rubbere

The viscous component contributes to damping:

Effective damping co-efficient:

Using the Kelvin–Voigt model, the instantaneous power P(t)P(t)P(t) dissipated is:

6
 P(t) = Stress×Strain rate= ηr ⋅ (strain rate)2

In beam theory damping power per unit volume becomes:

The instantaneous power dissipated is:

Computing average power over one period T = 2𝜋

𝑤

Since the average 𝑠𝑖𝑛2𝑤𝑡 over one cycle is 0.5

Definition of effective damping co-efficient

Equating the two equations;

The above expression is the effective damping considering energy dissipation due to internal
friction.

Substituting the above in the governing equation:

M_eq = ρeq*A*L;

7
 𝐸𝑠𝑡𝑒𝑒𝑙 ℎ +𝐸
𝑠𝑡𝑒𝑒𝑙 𝑟𝑢𝑏𝑏𝑒𝑟 ℎ𝑟𝑢𝑏𝑏𝑒𝑟

𝐸𝑒𝑞 = ℎ𝑡𝑜𝑡𝑎𝑙

𝑐𝑒𝑓𝑓 = 𝜂𝑟 ∗ 𝑏 ∗ ℎ𝑟 ∗ ∫ (𝜓′(𝑥))2 𝑑𝑥
0

Keq = meq*𝑤𝑛2

FORMUALE FOR TIME CHARACTERISTICS: Used to analyze the

response of the system.

Settling time: Settling time is the time required for a system’s response to
remain within a specified percentage (commonly 2% or 5%) of its final
steadystate value after a disturbance or initial excitation. It indicates how
quickly the system stabilizes and stops oscillating significantly.

4
𝑇𝑠 =
𝜁 𝑤𝑛

Peak time: Peak time is the time it takes for the system response to reach its first
maximum (peak) value after an initial disturbance or input. It reflects how
quickly the system reacts and overshoots before settling down.

𝜋
𝑇𝑝 =
𝑤𝑑

Rise time: Rise time is the time taken by the system response to increase from a
low percentage (typically 10%) to a high percentage (typically 90%) of its first
peak value (or final value, depending on context). It indicates how fast the
system initially responds to an input or disturbance.

8
 Tr = 𝜋 −𝜃
where, 𝜃
𝑤𝑑 𝜁

Steady state error: It is the difference between the system’s final

(steadystate) response and the desired value or reference input as time
approaches infinity.

CHAPTER 2

9
 METHODOLOGY

1. Beam Model and Materials

A cantilever beam model was developed, consisting of a steel core with an

added layer of viscoelastic rubber. The geometric parameters such as beam
length, width, and steel thickness were defined. Material properties including
density, Young’s modulus, and viscosity were assigned for both steel and rubber.

2. Rubber Thickness Variation

To study the effect of damping, the rubber layer thickness was varied from 0 to
the same thickness as the steel. A total of 100 thickness values were analyzed in
a parametric sweep.

3. Calculation of Equivalent Properties

For each thickness value, the overall cross-sectional properties such as area,
moment of inertia, and equivalent material properties were computed based on
the contributions from both steel and rubber.

4. Natural Frequency and Mode Shape

The first bending natural frequency of the beam was calculated using standard
cantilever beam theory. The corresponding vibration mode shape was also
generated and normalized.

10
5. Damping Estimation

The damping effect of the rubber layer was estimated by integrating the mode
shape gradient, considering the rubber’s viscosity and its influence across the
beam.

6. Time-Domain Simulation

Settling time: It is calculated numerically by checking when the response stays

within 2% of the final value for a sustained duration.

Peak time: It is determined numerically by finding the time at which the second
peak occurs in the absolute value of the displacement response.

Rise time: Found as the time interval during which the absolute displacement
increases from 10% to 90% of the first peak value.

Steady state error: it is seen as the absolute value of the final displacement close
to zero.

7. Dynamic Response Analysis

For each rubber thickness: The damped displacement response was recorded.

Time-domain characteristics such as peak time, rise time, settling time, and
steady-state error were extracted using standard signal analysis techniques.

8. Visualization

Time responses for two selected rubber thickness values were plotted for
comparison. Finally, the variation of all time-domain response metrics with
respect to rubber thickness was displayed using subplots.

11
 CHAPTER 3

RESULTS

Figure 3.1 Sample response of beam when 0.2mm thickness of rubber is used

Time characteristics Time in (sec)

Peak time 0.076002

Rise time 0.004401

Settling time 9.8158

Steady state error 0.00015313

12
Figure 3.2 variation of peak time when rubber thickness ratoo is increasing

Figure 3.3 Peak time vs rubber thickness ratio

13
 Figure 3.4 settling time vs rubber thickness ratio

Figure 3.5 steady-state error vs rubber thickness ratio

Final result: The observation of settling time, peak time , rise time
and steady state error suggests that a rubber thickness of 0.3mm to
0.5 mm thickness is optimum value.

CHAPTER 4

14
 CONCLUSION
This study effectively demonstrates how integrating a viscoelastic
rubber layer with a steel cantilever beam significantly influences the
dynamic response of the structure. By modelling the system as a
damped single-degree-of-freedom oscillator and employing
MATLAB-based numerical simulations, critical time-domain
response metrics—such as peak time, rise time, settling time, and
steady-state error—were accurately evaluated for varying rubber
thicknesses.
The results clearly show that increasing rubber thickness enhances
damping behaviour, reducing steady-state error and improving
vibration isolation. However, this also leads to longer settling times,
indicating a trade-off between stability and responsiveness. From the
analysis, a rubber thickness range of 0.3 mm to 0.5 mm is identified
as optimal, offering a balance between effective damping and
acceptable dynamic response times. These findings are valuable for
the design and optimization of lightweight, vibration-resistant
aerospace structures using layered composite materials.

CHAPTER 5
REFERENCE
[1] Ehsan Omidi, S. Nima Mahmoodi, W. Steve Shepard Jr.,
‘Vibration reduction in aerospace structures via an optimized

15
modified positive velocity feedback control’ 1270-9638/© 2015
Elsevier Masson SAS.

CHAPTER 6
MATLAB CODE
clc; clear; close all;
%% Beam Geometry and Material Properties
L = 0.25; % Beam length (m) b=
0.03; % Width (m) h_steel = 0.001;
% Steel thickness (m)

16
 E_steel = 193e9; % Young's Modulus of steel (Pa)
rho_steel = 8030; % Density of steel (kg/m^3)

E_rubber = 5e6; % Young's modulus of rubber (Pa)

rho_rubber = 1200; % Density of rubber (kg/m^3) eta
= 5e3; % Viscosity of rubber (Pa.s)

%% Rubber Thickness Sweep

h_rubber_vec = linspace(0, 1, 100) * h_steel;

% Preallocate results t_peak_vec =

zeros(size(h_rubber_vec)); t_rise_vec =
zeros(size(h_rubber_vec)); t_settle_vec =
zeros(size(h_rubber_vec)); e_ss_vec =
zeros(size(h_rubber_vec));

for k = 1:length(h_rubber_vec)
h_rubber = h_rubber_vec(k); h_total
= h_steel + h_rubber;

% Composite Section Properties

A = b * h_total;
I = (b * h_total^3) / 12;
rho_eq = (rho_steel * h_steel + rho_rubber * h_rubber) / h_total;
E_eq = (E_steel * h_steel + E_rubber * h_rubber) / h_total;

17
 % Frequency Calculation f_eq = @(x)
cos(x).*cosh(x) + 1; bL = fzero(f_eq, 1.875); beta =
bL / L; omega_n_damped = beta^2 * sqrt(E_eq * I /
(rho_eq * A));

% Mode Shape x = linspace(0, L,

100); phi = cosh(beta * x) - cos(beta *
x) - ...
((cosh(bL) + cos(bL)) / (sinh(bL) + sin(bL))) * ...
(sinh(beta * x) - sin(beta * x));
phi = phi / max(abs(phi));

% Damping Coefficient c_damping = eta * b * h_rubber *

trapz(x, (gradient(phi, x)).^2); % Time and Observation obs_point =
0.1 * L;
[~, obs_index] = min(abs(x - obs_point));
t = linspace(0, 10, 50000); dt = t(2) - t(1);
amp = 1; phi_i = phi(obs_index);

m_eq = rho_eq * A * L / 2;
k_eq = m_eq * omega_n_damped^2;
c_eq = c_damping;

% Time Integration (Euler method)

y = zeros(size(t)); y_dot =
zeros(size(t)); y(1) = amp * phi_i;

18
 for n = 1:length(t)-1 y_ddot = (-c_eq *
y_dot(n) - k_eq * y(n)) / m_eq; y_dot(n+1) =
y_dot(n) + dt * y_ddot; y(n+1) = y(n) + dt *
y_dot(n+1); end

y_damped = y;
if k == 1 figure; hold
on;
plot(t, y_damped, 'r', 'DisplayName', sprintf('h_{rubber} = %.4f
m', h_rubber));
end

if k == 20 plot(t, y_damped, 'b', 'DisplayName',

sprintf('h_{rubber} = %.4f m', h_rubber));
end
% Time-Domain Characteristics
[peaks, locs] = findpeaks(abs(y_damped), t);

if length(peaks) >= 2
y_peak = peaks(2); t_peak
= locs(2);

% Rise Time: 10% to 90% of peak t_window =

[t_peak - 0.01, t_peak]; idx_window = find(t >=

19
t_window(1) & t <= t_window(2)); y_segment =
abs(y_damped(idx_window)); t_segment = t(idx_window);

y10 = 0.1 * y_peak;

y90 = 0.9 * y_peak;
idx10 = find(y_segment >=
y10, 1, 'first'); idx90 =
find(y_segment >= y90, 1,
'first');

if ~isempty(idx10) && ~isempty(idx90) && idx90 > idx10

t_rise = t_segment(idx90) - t_segment(idx10); else
t_rise = NaN; end else y_peak = NaN;
t_peak = NaN; t_rise = NaN; end xlabel('Time (s)');
ylabel('Displacement (m)'); title('response for thickness of
0.2mm thickness of rubber') grid on; legend('show');
% Settling Time (within 2% of final value) tolerance = 0.02
* abs(y_peak); window_length = round(0.01 * length(t)); %
window = 1% of total time samples

within_bounds = abs(y_damped - y_damped(end)) <= tolerance;

% Look for the first time index after which all remaining values are
within bounds for idx = 1:(length(t) - window_length) if
all(within_bounds(idx:idx + window_length)) t_settle = t(idx);
break; else t_settle = NaN; end end

20
 % Steady-state error
e_ss = abs(y_damped(end));

% Store results
t_peak_vec(k) = t_peak;
t_rise_vec(k) = t_rise;
t_settle_vec(k) = t_settle;
e_ss_vec(k) = e_ss; end

%% Plotting figure; subplot(2,2,1);

plot(h_rubber_vec*1e3, t_peak_vec, 'b');
xlabel('Rubber Thickness ratio');
ylabel('Peak Time (s)');

title('Peak Time vs Rubber Thickness'); grid on;

subplot(2,2,2); plot(h_rubber_vec*1e3, t_rise_vec,

'r'); xlabel('Rubber Thickness ratio'); ylabel('Rise Time
(s)'); title('Rise Time vs Rubber Thickness'); grid on;

subplot(2,2,3); plot(h_rubber_vec*1e3, t_settle_vec, 'g');

xlabel('Rubber Thickness ratio'); ylabel('Settling Time (s)');
title('Settling Time vs Rubber Thickness'); grid on;

subplot(2,2,4); plot(h_rubber_vec*1e3, e_ss_vec, 'k');

xlabel('Rubber Thickness ratio'); ylabel('Steady-State Error (m)');
title('Steady-State Error vs Rubber Thickness'); grid on;

21
 disp(['Peak time = ', num2str(t_peak_vec(5))]);
disp(['rise time = ', num2str(t_rise_vec(5))]);
disp(['settling time = ', num2str(t_settle_vec(5))]);
disp(['steady state error = ', num2str(e_ss_vec(5))]);

22