Sunday, 23 March 2025

Solving an Undamped System Under Harmonic Force Using Laplace Transform

Topic: Harmonically Excited Vibration
Subject: Mechanical Vibration
Tool: Euler Math Toolbox, QCAD, & Scilab

By: Gani Comia, Mar. 2025

  • Undamped System Under Harmonic Force - Equation of Motion

The undamped system under harmonic force is a dynamic system subjected to external force or excitation. Consider the harmonic excitations of the form

$$F(t) = F_0 \cos(\omega t) \tag{1}$$

Where:
\(F(t)\) – excitation of forcing function
\(F_0\) – constant static force
\(\omega\) – excitation angular frequency
\(t\) - time

If the force \(F(t)\) acts on an undamped spring-mass system, the equation of motion is

$$m\,y'' + k\,y = F(t) \tag{2}$$

or

$$m\,y'' + k\,y = F_0 \cos(\omega t) \tag{3}$$

Where:
\(m\) – mass
\(k\) – spring constant
\(y’’\) – acceleration of mass
\(y\) – displacement from rest position

Under a harmonic excitation, the response of the system will also be harmonic. If the frequency of excitation \(\omega\) coincides with the natural frequency \(\omega_n\) of the system, the response of the system will be very large. The condition, \(\omega = \omega_n\), is what we called resonance. Resonance is to be avoided to prevent system’s failure as in this case mechanical damage.

$$\omega = \omega_n \quad \left(\omega_n = \sqrt{k / m}\right) \tag{4} $$

Where:
\(\omega_n\) – natural angular frequency

  • Laplace Transform of Second-Order Linear Differential Equation

Laplace transform is one of the methods to solve for the response for the undamped system under harmonic force. Consider a general second-order linear differential equation

$$a\,\frac{d^2 y}{d t^2} + b\, \frac{dy}{dt} + c\,y = u(t) \quad \left(t \geq 0 \right) \tag{5}$$

Subject to initial conditions

$$t = 0, \quad y(0) = y_0, \quad \frac{dy}{dt}(0) = \nu_0 \tag{6}$$

Such differential equations may model the dynamics of some system for which the variable \(y(t)\) determines the response of the system to the forcing or excitation term \(u(t)\). The terms system input and system output are also frequently used for \(u(t)\) and \(y(t)\) respectively.

  • Laplace and Inverse Laplace Transform of an Undamped System Under Harmonic Force

Taking the Laplace transform of each term of Equation (3)

$$m\, \mathcal{L} \left\{ \frac{d^2 y}{d t^2} \right\} + k\, \mathcal{L} \{y\} = F_0 \, \mathcal{L} \{ \cos(\omega t)\} \tag{7}$$

With initial conditions

$$t = 0, \quad y(0) = \frac{F_0}{k - m\,\omega^2}, \quad \frac{dy}{dt}(0) = 0 \tag{8}$$

Leads to

$$m \, \left( s^2 \, Y(s) - \frac{s\,F_0}{k - m \, \omega^2} \right) + k\,Y(s) = \frac{s \, F_0}{\omega^2 + s^2} \tag{9}$$

By rearranging and solving for \(Y(s)\) gives

$$Y(s) = - \frac{s\,F_0}{m\,\omega^4 + \omega^2 \, (m\,s^2 - k) - k\,s^2} \tag{10}$$

The inverse Laplace transform solves the system output \(y(t)\) for \(\omega > 0\)

$$\mathcal{L^{-1}} \{Y(s)\} = \mathcal{L^{-1}} \left\{ - \frac{s\,F_0}{m\,\omega^4 + \omega^2 \, (m\,s^2 - k) - k\,s^2} \right\} \quad (\omega > 0) \tag{11}$$

The steady-state response \(y(t)\) is

$$y(t) = \frac{F_0 \, \cos(\omega \, t)}{m \, \omega^2 - k} \tag{12}$$

The response of the system can be identified to be of three types.

Case 1. When \(0 < {\omega}/{\omega_n} < 1\), the harmonic response of the system \(y(t)\) is said to be lagging the force.

Case 2. When \({\omega}/{\omega_n} > 1\), the response of the system to a harmonic force of very high frequency is close to zero. The response of the system is said to be leading the force.

Case 3. When \({\omega}/{\omega_n} = 1\), the response is called resonance and the amplitude of \(y(t)\) becomes infinite. This phenomenon can result on large displacements and stresses.

Consider a sample application illustrated in Figure 1, let us analyze the response \(y(t)\) in terms of excitation \(F(t)\) based on the Case 1 and 2 responses. Case 3 will not be possible to plot. For this let us assume an excitation frequency, \(\omega\), of less than and more than the natural frequency, \(\omega_n\).

Figure 1. Application of an Undamped System under Harmonic Force

Below is the Scilab script for plotting both the system input and output for the illustrated problem.

  • Scilab Script

// Copyright (C) 2025 - Gani Comia
// Date of creation: 22 Mar 2025
// Undamped System under Harmonic Force using Laplace Transform
clear;clc;
// (1) solution y(t) from m y''(t) + k y(t) = F0 cos(wt)
function y=position(t)
    y = -(F0.*cos(omega.*t))./((m.*omega.^2)-k)
endfunction

// (2) excitation of forcing function
function F=force(t)
    F = F0.*cos(omega.*t)
endfunction

// (3) primary parameters from m y''(t) + k y(t) = F0 cos(wt)
m = 100/9.8;                        // kg, mass (where W = 100 N)
k = 2000;                           // N/m, spring constant
F0 = 25;                            // N, constant static force

// (4) secondary parameters
omega_n = sqrt(k./m);               // Hz, natural angular frequency
mprintf("omega_n (Hz) = %3.1f\n", omega_n)

clf;
// (5) figure properties
fig = gcf()
fig.figure_size = [700,800];

// (6) for excitation or forcing function
omega = omega_n
omega_t = linspace(0,2*%pi,100);
t = omega_t./omega;

// (6-1) plot of excitation force
subplot(3,1,1)
F = force(t)
plot(t,F,"m-","linewidth",4)
title("$\Large\text{Excitation Function}\quad F(t)=F_0\cos(\omega t)$")
ylabel("$\Large \mathbf{F(t)}$")
xlabel("$\Large \text{time,}\; t$")
legend("$\LARGE\omega = \omega_n$",4,with_box=%F)
xgrid(3,0)
// (6-2) axes properties
ax = gca()
ax.x_location = 'origin' 
ax.data_bounds = [0 -30; 0.5 30]

// (7) for Case 1. 0 < omega < omega_n
omega = 0.95*omega_n
omega_t = linspace(0,2*%pi,100);
t = omega_t./omega;

// (7-1) plot of a sampe case 1
subplot(3,1,2)
y = position(t)
plot(t,y,"b-","linewidth",4)
title("$\Large\text{Displacement at}\quad \mathbf{\omega < \omega_n}$")
ylabel("$\Large \mathbf{y(t)}$")
xlabel("$\Large \text{time,}\; t$")
legend("$\LARGE \omega = \text{0.95} \; \omega_n$",4,with_box=%F)
xgrid(3,0)
ax = gca()
ax.x_location = 'origin'
ax.data_bounds = [0 -0.2; 0.5 0.2]

// (8) for Case 2. omega > omega_n
omega = 1.05*omega_n
omega_t = linspace(0,2*%pi,100);
t = omega_t./omega;

// (8-1) plot of a sampe case 2
subplot(3,1,3)
y = position(t)
plot(t,y,"r-","linewidth",4)
title("$\Large\text{Displacement at}\quad \mathbf{\omega > \omega_n}$")
ylabel("$\Large \mathbf{y(t)}$")
xlabel("$\Large \text{time,}\; t$")
legend("$\LARGE \omega = \text{1.05} \; \omega_n$",1,with_box=%F)
xstring(0.31,-0.2,"https://gani-mech-toolbox.blogspot.com")
xgrid(3,0)
ax = gca()
ax.x_location = 'origin'
ax.data_bounds = [0 -0.2; 0.5 0.2] 

  • Scilab Output (Figure 2)

Figure 2. Forcing Function and Responses to Undamped System with Harmonic Excitation

The harmonic steady-state response \(y(t)\) for case 1 is said to be lagging with the forcing function \(F(t)\). As for the case 2 the response \(y(t)\) is said to be leading with external harmonic force \(F(t)\).

Feel free to comment for inquiry, clarification, or suggestion for improvement. Drop your email to request the soft copy of the file.

References
  1. Singiresu S. Rao. Mechanical Vibrations. 2nd Ed. Addison-Wesley Publishing Company. 1990.
  2. G. James, D. Burley, P. Dyke, J. Searl, N. Steele, J. Wright. Advanced Modern Engineering Mathematics. Addison-Wesley Publishing Company, Inc. 1993.

No comments:

Post a Comment

C++ and Python for Numerical Solution to Spring-Mass-Damper Model

Topic: SDOF Spring-Mass-Damper System Subject: Numerical Methods & Mechanical Vibration Tools: C++, CodeBlocks, GitHub Copilot, LabP...