The Physics of the Damped Harmonic Oscillator

Open Live Script

This example explores the physics of the damped harmonic oscillator by solving the equations of motion in the case of no driving forces. This example investigates the cases of under-, over-, and critical-damping.

Derive Equation of Motion
Solve the Equation of Motion (F = 0)
Underdamped Case ( $ζ < 1$ )
Overdamped Case ( $ζ > 1$ )
Critically Damped Case ( $ζ = 1$ )
Conclusion

1. Derive Equation of Motion

Consider a forced harmonic oscillator with damping shown below. Model the resistance force as proportional to the speed with which the oscillator moves.

Define the equation of motion where

$m$ is the mass
$c$ is the damping coefficient
$k$ is the spring constant
$F$ is a driving force

syms x(t) m c k F(t)
eq = m*diff(x,t,t) + c*diff(x,t) + k*x == F

eq(t) = 
 $m \frac{\partial^{2}}{\partial t^{2}} x (t) + c \frac{\partial}{\partial t} x (t) + k x (t) = F (t)$

Rewrite the equation using $c = m γ$ and $k = m ω_{0}^{2}$ .

syms gamma omega_0
eq = subs(eq, [c k], [m*gamma, m*omega_0^2])

eq(t) = 
 $m \frac{\partial^{2}}{\partial t^{2}} x (t) + γ m \frac{\partial}{\partial t} x (t) + m x (t) {ω_{0}}^{2} = F (t)$

Divide out the mass $m$ . Now we have the equation in a convenient form to analyze.

eq = collect(eq, m)/m

eq(t) = 
 $\frac{\partial^{2}}{\partial t^{2}} x (t) + γ \frac{\partial}{\partial t} x (t) + x (t) {ω_{0}}^{2} = \frac{F (t)}{m}$

2. Solve the Equation of Motion where F = 0

Solve the equation of motion using dsolve in the case of no external forces where $F = 0$ . Use the initial conditions of unit displacement and zero velocity.

vel = diff(x,t);
cond = [x(0) == 1, vel(0) == 0];
eq = subs(eq,F,0);
sol = dsolve(eq, cond)

sol = 
 $\begin{array}{l} \frac{e^{- t (\frac{γ}{2} - \frac{σ_{1}}{2})} (γ + σ_{1})}{2 σ_{1}} - \frac{e^{- t (\frac{γ}{2} + \frac{σ_{1}}{2})} (γ - σ_{1})}{2 σ_{1}} \\ where \\ σ_{1} = \sqrt{(γ - 2 ω_{0}) (γ + 2 ω_{0})} \end{array}$

Examine how to simplify the solution by expanding it.

sol = expand(sol)

sol = 
 $\begin{array}{l} \frac{σ_{1} σ_{2}}{2} + \frac{σ_{1} e^{\frac{t \sqrt{γ^{2} - 4 {ω_{0}}^{2}}}{2}}}{2} - \frac{γ σ_{1} σ_{2}}{2 \sqrt{γ^{2} - 4 {ω_{0}}^{2}}} + \frac{γ σ_{1} e^{\frac{t \sqrt{γ^{2} - 4 {ω_{0}}^{2}}}{2}}}{2 \sqrt{γ^{2} - 4 {ω_{0}}^{2}}} \\ where \\ σ_{1} = e^{- \frac{γ t}{2}} \\ σ_{2} = e^{- \frac{t \sqrt{γ^{2} - 4 {ω_{0}}^{2}}}{2}} \end{array}$

Notice that each term has a factor of $σ_{1}$ , or $e^{- γ t / 2}$ , use collect to gather these terms

sol = collect(sol, exp(-gamma*t/2))

sol = 
 $\begin{array}{l} (\frac{e^{- σ_{1}}}{2} + \frac{e^{σ_{1}}}{2} - \frac{γ e^{- σ_{1}}}{2 \sqrt{γ^{2} - 4 {ω_{0}}^{2}}} + \frac{γ e^{σ_{1}}}{2 \sqrt{γ^{2} - 4 {ω_{0}}^{2}}}) e^{- \frac{γ t}{2}} \\ where \\ σ_{1} = \frac{t \sqrt{γ^{2} - 4 {ω_{0}}^{2}}}{2} \end{array}$

The term $\sqrt{γ^{2} - 4 ω_{0}^{2}}$ appears in various parts of the solution. Rewrite it in a simpler form by introducing the damping ratio $ζ \equiv \frac{γ}{2 ω_{0}}$ .

Substituting ζ into the term above gives:

$\sqrt{γ^{2} - 4 ω_{0}^{2}} = 2 ω_{0} \sqrt{{(\frac{γ}{2 ω_{0}})}^{2} - 1} = 2 ω_{0} \sqrt{ζ^{2} - 1},$

syms zeta;
sol = subs(sol, ...
    sqrt(gamma^2 - 4*omega_0^2), ...
    2*omega_0*sqrt(zeta^2-1))

sol = 
 $\begin{array}{l} e^{- \frac{γ t}{2}} (\frac{σ_{2}}{2} + \frac{σ_{1}}{2} + \frac{γ σ_{2}}{4 ω_{0} \sqrt{ζ^{2} - 1}} - \frac{γ σ_{1}}{4 ω_{0} \sqrt{ζ^{2} - 1}}) \\ where \\ σ_{1} = e^{- ω_{0} t \sqrt{ζ^{2} - 1}} \\ σ_{2} = e^{ω_{0} t \sqrt{ζ^{2} - 1}} \end{array}$

Further simplify the solution by substituting $γ$ in terms of $ω_{0}$ and $ζ$ ,

sol = subs(sol, gamma, 2*zeta*omega_0)

sol = 
 $\begin{array}{l} e^{- ω_{0} t ζ} (\frac{σ_{2}}{2} + \frac{σ_{1}}{2} + \frac{ζ σ_{2}}{2 \sqrt{ζ^{2} - 1}} - \frac{ζ σ_{1}}{2 \sqrt{ζ^{2} - 1}}) \\ where \\ σ_{1} = e^{- ω_{0} t \sqrt{ζ^{2} - 1}} \\ σ_{2} = e^{ω_{0} t \sqrt{ζ^{2} - 1}} \end{array}$

We have derived the general solution for the motion of the damped harmonic oscillator with no driving forces. Next, we'll explore three special cases of the damping ratio $ζ$ where the motion takes on simpler forms. These cases are called

underdamped $(ζ < 1)$ ,
overdamped $(ζ > 1)$ , and
critically damped $(ζ = 1)$ .

3. Underdamped Case ( $ζ < 1$ )

If $ζ < 1$ , then $\sqrt{ζ^{2} - 1} = i \sqrt{1 - ζ^{2}}$ is purely imaginary

solUnder = subs(sol, sqrt(zeta^2-1), 1i*sqrt(1-zeta^2))

solUnder = 
 $\begin{array}{l} e^{- ω_{0} t ζ} (\frac{σ_{1}}{2} + \frac{σ_{2}}{2} + \frac{ζ σ_{1} i}{2 \sqrt{1 - ζ^{2}}} - \frac{ζ σ_{2} i}{2 \sqrt{1 - ζ^{2}}}) \\ where \\ σ_{1} = e^{- ω_{0} t \sqrt{1 - ζ^{2}} i} \\ σ_{2} = e^{ω_{0} t \sqrt{1 - ζ^{2}} i} \end{array}$

Notice the terms $e^{i ω_{0} t \sqrt{ζ^{2} - 1}} \pm e^{- {i ω}_{0} t \sqrt{ζ^{2} - 1}}$ in the above equation and recall the identity $e^{i x} = \cos (x) + i \sin (x) .$

Rewrite the solution in terms of $\cos$ .

solUnder = coeffs(solUnder, zeta);
solUnder = solUnder(1);
c = exp(-omega_0 * zeta * t);
solUnder = c * rewrite(solUnder / c, 'cos')

solUnder =  $e^{- ω_{0} t ζ} \cos (ω_{0} t \sqrt{1 - ζ^{2}})$

solUnder(t, omega_0, zeta) = solUnder

solUnder(t, omega_0, zeta) =  $e^{- ω_{0} t ζ} \cos (ω_{0} t \sqrt{1 - ζ^{2}})$

The system oscillates at a natural frequency of $ω_{0} \sqrt{1 - ζ^{2}}$ and decays at an exponential rate of $1 / ω_{0} ζ$ .

Plot the solution with fplot as a function of $ω_{0} t$ and $ζ$ .

z = [0 1/4 1/2 3/4];
w = 1;
T = 4*pi;
lineStyle = {'-','--',':k','-.'};

fplot(@(t)solUnder(t, w, z(1)), [0 T], lineStyle{1});

hold on;
for k = 2:numel(z)
    fplot(@(t)solUnder(t, w, z(k)), [0 T], lineStyle{k});
end
hold off;
grid on;
xticks(T*linspace(0,1,5));
xticklabels({'0','\pi','2\pi','3\pi','4\pi'});
xlabel('t / \omega_0');
ylabel('amplitude');
lgd = legend('0','1/4','1/2','3/4');
title(lgd,'\zeta');
title('Underdamped');

4. Overdamped Case ( $ζ > 1$ )

If $ζ > 1$ , then $\sqrt{ζ^{2} - 1}$ is purely real and the solution can be rewritten as

solOver = sol

solOver = 
 $\begin{array}{l} e^{- ω_{0} t ζ} (\frac{σ_{2}}{2} + \frac{σ_{1}}{2} + \frac{ζ σ_{2}}{2 \sqrt{ζ^{2} - 1}} - \frac{ζ σ_{1}}{2 \sqrt{ζ^{2} - 1}}) \\ where \\ σ_{1} = e^{- ω_{0} t \sqrt{ζ^{2} - 1}} \\ σ_{2} = e^{ω_{0} t \sqrt{ζ^{2} - 1}} \end{array}$

solOver = coeffs(solOver, zeta);
solOver = solOver(1)

solOver = 
 $e^{- ω_{0} t ζ} (\frac{e^{ω_{0} t \sqrt{ζ^{2} - 1}}}{2} + \frac{e^{- ω_{0} t \sqrt{ζ^{2} - 1}}}{2})$

Notice the terms $\frac{(e^{ω_{0} t \sqrt{ζ^{2} - 1}} + e^{- ω_{0} t \sqrt{ζ^{2} - 1}})}{2}$ and recall the identity $\cosh (x) = \frac{e^{x} + e^{- x}}{2}$ .

Rewrite the expression in terms of $\cosh$ .

c = exp(-omega_0*t*zeta);
solOver = c*rewrite(solOver / c, 'cosh')

solOver =  $\cosh (ω_{0} t \sqrt{ζ^{2} - 1}) e^{- ω_{0} t ζ}$

solOver(t, omega_0, zeta) = solOver

solOver(t, omega_0, zeta) =  $\cosh (ω_{0} t \sqrt{ζ^{2} - 1}) e^{- ω_{0} t ζ}$

Plot the solution to see that it decays without oscillating.

z = 1 + [1/4 1/2 3/4 1];
w = 1;
T = 4*pi;
lineStyle = {'-','--',':k','-.'};

fplot(@(t)solOver(t, w, z(1)), [0 T], lineStyle{1});

hold on;
for k = 2:numel(z)
    fplot(@(t)solOver(t, w, z(k)), [0 T], lineStyle{k});
end
hold off;
grid on;
xticks(T*linspace(0,1,5));
xticklabels({'0','\pi','2\pi','3\pi','4\pi'});
xlabel('\omega_0 t');
ylabel('amplitude');
lgd = legend('1+1/4','1+1/2','1+3/4','2');
title(lgd,'\zeta');
title('Overdamped');

5. Critically Damped Case ( $ζ = 1$ )

If $ζ = 1$ , then the solution simplifies to

solCritical(t, omega_0) = limit(sol, zeta, 1)

solCritical(t, omega_0) =  $e^{- ω_{0} t} (ω_{0} t + 1)$

Plot the solution for the critically damped case.

w = 1;
T = 4*pi;

fplot(solCritical(t, w), [0 T])
xlabel('\omega_0 t');
ylabel('x');
title('Critically damped, \zeta = 1');
grid on;
xticks(T*linspace(0,1,5));
xticklabels({'0','\pi','2\pi','3\pi','4\pi'});

6. Conclusion

We have examined the different damping states for the harmonic oscillator by solving the ODEs which represents its motion using the damping ratio $ζ$ . Plot all three cases together to compare and contrast them.

zOver  = pi;
zUnder = 1/zOver;
w = 1;
T = 2*pi;
lineStyle = {'-','--',':k'};

fplot(@(t)solOver(t, w, zOver), [0 T], lineStyle{1},'LineWidth',2);
hold on;
fplot(solCritical(t, w), [0 T], lineStyle{2},'LineWidth',2)
fplot(@(t)solUnder(t, w, zUnder), [0 T], lineStyle{3},'LineWidth',2);
hold off;

textColor = lines(3);
text(3*pi/2, 0.3 , 'over-damped'      ,'Color',textColor(1,:));
text(pi*3/4, 0.05, 'critically-damped','Color',textColor(2,:));
text(pi/8  , -0.1, 'under-damped');

grid on;
xlabel('\omega_0 t');
ylabel('amplitude');
xticks(T*linspace(0,1,5));
xticklabels({'0','\pi/2','\pi','3\pi/2','2\pi'});
yticks((1/exp(1))*[-1 0 1 2 exp(1)]);
yticklabels({'-1/e','0','1/e','2/e','1'});
lgd = legend('\pi','1','1/\pi');
title(lgd,'\zeta');
title('Damped Harmonic Oscillator');

The Physics of the Damped Harmonic Oscillator

Contents

1. Derive Equation of Motion

2. Solve the Equation of Motion where F = 0

3. Underdamped Case ( $ζ < 1$ )

4. Overdamped Case ( $ζ > 1$ )

5. Critically Damped Case ( $ζ = 1$ )

6. Conclusion

See Also

Functions

The Physics of the Damped Harmonic Oscillator

Contents

1. Derive Equation of Motion

2. Solve the Equation of Motion where F = 0

3. Underdamped Case (ζ<1)

4. Overdamped Case (ζ>1)

5. Critically Damped Case (ζ=1)

6. Conclusion

See Also

Functions

3. Underdamped Case ( $ζ < 1$ )

4. Overdamped Case ( $ζ > 1$ )

5. Critically Damped Case ( $ζ = 1$ )