Math 331.1 •  Fall 2005

Damped harmonic oscillators

Copyright © 2005 Murray Eisenberg.  All rights reserved.

Overdamped harmonic oscillator 
(mass m = 1, spring constant k = 8, damping coefficient b = 6)

y''  +  6 y'  +  8 y  =  0

Equivalent system of 1st order ODEs:

\(TraditionalForm
      \`y'\  = \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ v, 

      v'\  = \ \ \ \(-8\)\ y\  - \ 6\ v\)

Needs["Graphics`PlotField`"]

Clear[A];  A  =  { {0,  1},  {-8,  -6} };  
    MatrixForm[A]

Clear [F];  F[x_,  v_]  :=  A  .  {x,  v};  
    F[x,  v]

field  =  
        PlotVectorField[
          F[x,  v],  {x,  -3,  3}, {v,  -3,  3},  
          Axes  ->  True,  
          AxesLabel  ->  {"x",  "y"},  
          ScaleFunction  ->  (0.25 &),  
          ScaleFactor  ->  0.275,  
          PlotPoints  ->  20];

Eigensystem[A]

yk[t_]  :=  k1 Exp[-4 t]  +  k2 Exp[-2 t]

Solve[{yk[0]  ==  1,  yk'[0]  ==  2},  {k1,  k2}]

y[t_]  :=  -2  Exp[-4 t]  +  
          3 Exp[-2 t];
v[t_]  :=  y'[t];  v[t]

phasecurve  =  
        ParametricPlot[Evaluate[{y[t],  v[t]}],  {t,  0,  5},  
          AspectRatio  ->  Automatic,  
          PlotRange  ->  {{-3,  3},  {-3,  3}},  
          PlotStyle  ->  {Thickness[0.005],  Red}];

Show[{field,  phasecurve}];

Plot[{y[t],  v[t]},  {t,  0,  5},  
        PlotStyle  ->  {RGBColor[0, 0, 0],  
            Dashing[{0.04,  0.03}]}];

Underdamped harmonic oscillator 
(mass m = 1, spring constant k = 5, damping coefficient b = 2)

y''  +  2 y'  +  5 y  =  0

Equivalent system of 1st order ODEs:

\(TraditionalForm\`y'\  = \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ v, 
v'\  = \ \ \ \(-5\)\ y\  - \ 2\ v\)

Needs["Graphics`PlotField`"]

Clear[A];  A  =  { {0,  1},  {-5,  -2} };  
    MatrixForm[A]

Clear [F];  F[x_,  v_]  :=  A  .  {x,  v};  
    F[x,  v]

field  =  
        PlotVectorField[
          F[x,  v],  {x,  -3,  3}, {v,  -3,  3},  
          Axes  ->  True,  
          AxesLabel  ->  {"x",  "y"},  
          ScaleFunction  ->  (0.25 &),  
          ScaleFactor  ->  0.275,  
          PlotPoints  ->  20];

Eigensystem[A]

yk[t_]  :=  
      k1 Exp[(-1  +  2  I) t]  +  
        k2 Exp[(-1  -  2  I) t]

Solve[{yk[0]  ==  1,  yk'[0]  ==  2},  {k1,  k2}]

y[
          t_]  :=  ((1/2) - ((3 I)/4)) Exp[(-1  +  2  I) t]  +  ((1/2) + ((3 I)/4)) 
            Exp[(-1  -  2  I) t];
v[t_]  :=  y'[t];

ComplexExpand[y[t]]
ComplexExpand[v[t]]

phasecurve  =  
        ParametricPlot[Evaluate[{y[t],  v[t]}],  {t,  0,  5},  
          AspectRatio  ->  Automatic,  
          PlotRange  ->  {{-3,  3},  {-3,  3}},  
          PlotStyle  ->  {Thickness[0.005],  Red}];

Show[{field,  phasecurve}];

Plot[{y[t],  v[t]},  {t,  0,  5},  
        PlotStyle  ->  {RGBColor[0, 0, 0],  
            Dashing[{0.04,  0.03}]}];

Example:  Critically damped harmonic oscillator 
(mass m = 1, spring constant k = 4, damping coefficient b = 4)

y''  +  + 4 y'  +  4 y  =  0

Equivalent system of 1st order ODEs:

\(TraditionalForm\`y'\  = \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ v, 
v'\  = \ \ \ \(-4\)\ y\  - \ \(\(4\)\(\ \)\(v\)\(\ \)\)\)

Needs["Graphics`PlotField`"]

Clear[A];  A  =  { {0,  1},  {-4,  -4} };  
    MatrixForm[A]

Clear [F];  F[x_,  v_]  :=  A  .  {x,  v};  
    F[x,  v]

field  =  
        PlotVectorField[
          F[x,  v],  {x,  -3,  3}, {v,  -3,  3},  
          Axes  ->  True,  
          AxesLabel  ->  {"x",  "y"},  
          ScaleFunction  ->  (0.25 &),  
          ScaleFactor  ->  0.275,  
          PlotPoints  ->  19];

Eigensystem[A]

yk[t_]  :=  k1 Exp[-2 t]  +  k2 t Exp[-2 t]

Solve[{yk[0]  ==  1,  yk'[0]  ==  2},  {k1,  k2}]

y[t_]  :=  Exp[-2 t]  + 4  t Exp[-2 t]
v[t_]  :=  y'[t];  v[t]

phasecurve  =  
        ParametricPlot[Evaluate[{y[t],  v[t]}],  {t,  0,  5},  
          AspectRatio  ->  Automatic,  
          PlotRange  ->  {{-3,  3},  {-3,  3}},  
          PlotStyle  ->  {Thickness[0.005],  Red}];

Show[{field,  phasecurve}];

Plot[{y[t],  v[t]},  {t,  0,  5},  
        PlotStyle  ->  {RGBColor[0, 0, 0],  
            Dashing[{0.04,  0.03}]}];

Critical damping vs. overdamping

Critically damped system

Clear[A];  A  =  { {0,  1},  {-4,  -4} };  
    Eigensystem[A]

Eigensystem[A]

yk[t_]  :=  k1 Exp[-2 t]  +  k2 t Exp[-2 t]

Solve[{yk[0]  ==  -1,  yk'[0]  ==  3},  {k1,  
        k2}]

y[t_]  :=  - Exp[-2 t]  + t Exp[-2 t];  
    v[t_]  :=  y'[t];  v[t]

criticallyDampedPhase  =  
        ParametricPlot[Evaluate[{y[t],  v[t]}],  {t,  0,  2.5},  
          AspectRatio  ->  Automatic,  
          PlotRange  ->  {{-0.3,  0.3},  {-0.3,  
                0.3}},   PlotLabel  ->  "Critically damped",  
          DisplayFunction  ->  Identity];

criticallyDampedFns  =  
        Plot[{y[t],  v[t]},  {t,  0,  2.5},  
          PlotStyle  ->  {RGBColor[0, 0, 0],  
              Dashing[{0.04,  0.03}]},  
          PlotLabel  ->  "Critically damped", 
          DisplayFunction  ->  Identity];

Overdamped system

Clear[A];  A  =  { {0,  1},  {-8,  -6} };  
    Eigensystem[A]

yk[t_]  :=  k1 Exp[-4 t]  +  k2 Exp[-2 t]

Solve[{yk[0]  ==  -1,  yk'[0]  ==  3},  {k1,  
        k2}]

y[t_]  :=  -(1/2)  Exp[-4 
              t]  +  -(1/2) Exp[-2 t]; 
    v[t_]  :=  y'[t];  v[t]

overdampedPhase  =  
        ParametricPlot[Evaluate[{y[t],  v[t]}],  {t,  0,  2.5},  
          AspectRatio  ->  Automatic,  
          PlotRange  ->  {{-0.3,  0.3},  {-0.3,  0.3}},  
          PlotLabel  ->  "Overdamped", 
          DisplayFunction  ->  Identity];

overdampedFns  =  
        Plot[{y[t],  v[t]},  {t,  0,  2.5},  
          PlotStyle  ->  {Black,  Dashing[{0.04,  0.03}]},  
          PlotLabel  ->  "Overdamped", 
          DisplayFunction  ->  Identity];

Comparison of solutions with same same initial conditions

Show[GraphicsArray[{criticallyDampedPhase,  overdampedPhase}],  
        DisplayFunction  ->  $DisplayFunction,  
        ImageSize  ->  72*7];

Show[GraphicsArray[{criticallyDampedFns,  overdampedFns}],  
        DisplayFunction  ->  $DisplayFunction,  
        ImageSize  ->  72*7];