# 22. Driven Damped Anharmonic Oscillators

Michael Fowler

## Introduction

Landau’ next sections (Chapter 6, sections 28,29) address nonlinear one-dimensional systems.  In particular, he focusses on driven damped oscillators with nonlinear, but small, added potential terms. Using ingenious semiquantitative techniques, he predicts some unexpected results: for example, a discontinuity in the oscillation amplitude on slowly varying the driving frequency at constant driving force (and constant damping). He also finds resonances when the driving frequency is a fraction, for example a third, of the oscillator’s natural frequency.

Fortunately, this system is easy to analyze numerically, and we have an applet to do just that. The parameters are set by sliders, and one can immediately find the large discontinuity in amplitude (factor of two or so) as the frequency is slightly changed.  At the end of this lecture, we show simple plots of amplitude response to a constant driving force as the frequency is varied.  These were found using the applet, the reader can easily check them, and venture into parts of the parameter space. The applet provides a measure of Landau’s (semiquantitative) accuracy, of course surprisingly good (of order 20% error or less) given the nature of the problem.

It should be added that this is one area where, thanks to computers, major advances have been made since Landau wrote the book, in particular the discovery for some systems of period doubling and chaos as the driving force is increased.  We’ve added a lecture (22a) on a particular system, the driven damped pendulum, a natural extension of Landau’s oscillator. This illustrates some of the novel features.  We will follow part of chapter 12 of Taylor’s excellent text, Classical Mechanics. Taylor provides many computer-generated graphs of the pendulum’s response as parameters are varied.  We provide applets that can generate these graphs.   The reader can easily use these applets to explore other parameter inputs.

In this lecture, to gain a bit of intuition about these nonlinear potentials, we’ll begin (following Landau) with no driving and no damping: just a particle oscillating in a potential that’s simple harmonic plus small ${x}^{3}$ and (positive) ${x}^{4}$ terms.  The basic questions are, how do these terms change the frequency of oscillation, and how does that frequency depend on the amplitude of oscillation?   The answers will guide us in understanding how a particle in such a potential will respond to a harmonic driving term, plus damping.

Next, we briefly review the driven damped linear oscillator (covered in detail in lecture 18, this is really just a reminder of the notation).  Then we add small cubic and quartic terms.  We present Landau’s argument  that--above a certain driving force--gradually increasing the driving frequency leads at a critical value to a discontinuous drop in the amplitude of the response, then use an applet to confirm and quantify his result.

## Frequency of Oscillation of a Particle is a Slightly Anharmonic Potential

See the applet illustrating this section.

Landau (para 28) considers a simple harmonic oscillator with added small potential energy terms $\frac{1}{3}m\alpha {x}^{3}+\frac{1}{4}m\beta {x}^{4}$.  In leading orders, these terms contribute separately, and differently, so it’s easier to treat them one at a time. We’ll first consider the quartic term, an equation of motion

$\stackrel{¨}{x}+{\omega }_{0}^{2}x=-\beta {x}^{3}.$

(We'll always take $\beta$ positive.)

Writing a perturbation theory expansion (following Landau):

$x={x}^{\left(1\right)}+{x}^{\left(2\right)}+\cdots \text{ }.$

(Standard practice in most books would be to write $x={x}^{\left(0\right)}+{x}^{\left(1\right)}+\dots$ with the superscript indicating the order of the perturbation--we're following Landau's notation, hopefully reducing confusion…)

We take as the leading term

${x}^{\left(1\right)}=a\mathrm{cos}\omega \text{\hspace{0.17em}}t$

with the exact value of $\omega$, $\omega ={\omega }_{\text{\hspace{0.17em}}0}+\Delta \omega$.  Of course, we don’t know the value of $\omega$ yet$—$this is what we’re trying to find!

And, as Landau points out, you can’t just write $\mathrm{cos}\left({\omega }_{\text{\hspace{0.17em}}0}+\Delta \omega \right)t=\mathrm{cos}{\omega }_{\text{\hspace{0.17em}}0}t-\left(\Delta \omega \right)t{\omega }_{0}\mathrm{sin}{\omega }_{0}t$ because that implies motion increasing in time, and our system is a particle oscillating in a fixed potential, with no energy supply. Furthermore, even if we did somehow have the value of $\omega$ exactly right, this expression would not be a full solution to the equation: the motion is certainly periodic with period $2\pi /\omega$, but the complete description of the motion is a Fourier series including frequencies $n\omega ,\text{ }n$ an integer, since the potential is no longer simple harmonic.

Anyway, putting this correct frequency into the equation of motion $\stackrel{¨}{x}+{\omega }_{0}^{2}x=-\beta {x}^{3}$ gives a nonzero left-hand side, so we rearrange.  We subtract $\left(1-\left({\omega }_{0}^{2}/{\omega }^{2}\right)\right)\stackrel{¨}{x}$ from both sides to get:

$\frac{{\omega }_{0}^{2}}{{\omega }^{2}}\stackrel{¨}{x}+{\omega }_{0}^{2}x=-\beta {x}^{3}-\left(1-\frac{{\omega }_{0}^{2}}{{\omega }^{2}}\right)\stackrel{¨}{x}.$

Now putting the leading term ${x}^{\left(1\right)}=a\mathrm{cos}\omega t$ into the left-hand side does give zero: if the equation had zero on the right hand side, this would just be a free (undamped) oscillator with natural frequency $\omega ,$ not ${\omega }_{\text{ }0}.$  This doesn’t look very promising, but keep reading.

The equation for the first-order correction ${x}^{\left(2\right)}$ is:

$\frac{{\omega }_{0}^{2}}{{\omega }^{2}}{\stackrel{¨}{x}}^{\left(2\right)}+{\omega }_{0}^{2}{x}^{\left(2\right)}=-\beta {\left({x}^{\left(1\right)}\right)}^{3}-\left(1-\frac{{\omega }_{0}^{2}}{{\omega }^{2}}\right){\stackrel{¨}{x}}^{\left(1\right)}.$

Notice that the second term on the right-hand side includes ${\stackrel{¨}{x}}^{\left(1\right)}=-{\omega }^{2}a\mathrm{cos}\omega \text{ }t$. This equation now represents a driving force on an undamped oscillator exactly at its resonant frequency, so would cause the amplitude to increase linearly, obviously an unphysical result, since we’re just modeling a particle sliding back and forth in a potential, with no energy being supplied from outside!

The key is that there is also a resonant driver in that first term $-\beta {\left({x}^{\left(1\right)}\right)}^{3}.$

Clearly these two driving terms have to cancel, and this requirement nails $\Delta \omega$: here’s how:

$-\beta {\left({x}^{\left(1\right)}\right)}^{3}=-\beta {a}^{3}{\mathrm{cos}}^{3}\omega t=-\beta {a}^{3}\left(\frac{3}{4}\mathrm{cos}\omega t+\frac{1}{4}\mathrm{cos}3\omega t\right)$

so the resonant driving terms cancel provided

$-\beta {a}^{3}\frac{3}{4}\mathrm{cos}\omega t-\left(1-\frac{{\omega }_{0}^{2}}{{\omega }^{2}}\right)\left(-{\omega }^{2}a\mathrm{cos}\omega t\right)=0.$

Remembering $\omega ={\omega }_{0}+\Delta \omega$, this gives (to this order)

$\Delta \omega =\frac{3\beta {a}^{2}}{8{\omega }_{0}}.$

(Strictly, ${\omega }_{0}+\frac{1}{2}\Delta \omega$ in the denominator, but that’s a higher order correction.)

Note that the frequency increases with amplitude:  the ${x}^{4}$ potential gives an increasingly stronger restoring force with amplitude than the harmonic well.  You can check this with the applet.

Now let’s consider the equation for a small cubic perturbation,

$\stackrel{¨}{x}+{\omega }_{0}^{2}x=-\alpha {x}^{2}.$

This represents an added potential $-\frac{1}{3}\alpha {x}^{3},$ which is an odd function, so to leading order it won’t change the period, speeding up one half of the oscillation and slowing the other half the same amount in leading order.  The first correction to the position as a function of time is the solution of

${\stackrel{¨}{x}}^{\left(2\right)}+{\omega }_{0}^{2}{x}^{\left(2\right)}=-\alpha {a}^{2}{\mathrm{cos}}^{2}\omega t=-\frac{1}{2}\alpha {a}^{2}-\frac{1}{2}\alpha {a}^{2}\mathrm{cos}2\omega t.$

The solution is

${x}^{\left(2\right)}=-\frac{\alpha {a}^{2}}{2{\omega }_{0}^{2}}+\frac{\alpha {a}^{2}}{6{\omega }_{0}^{2}}\mathrm{cos}2\omega t.$

Physically, adding this to the leading term, the particle is spending more of its time in the softer half of the potential, giving an amplitude-dependent correction to its average position.

To get the correction to the frequency, we need to go to the next order, $\omega ={\omega }_{0}+{\omega }^{\left(2\right)}.$ Dropping terms of higher order, the equation of motion for the next correction is

${\stackrel{¨}{x}}^{\left(3\right)}+{\omega }_{0}^{2}{x}^{\left(3\right)}=-2\alpha {x}^{\left(1\right)}{x}^{\left(2\right)}+2{\omega }_{0}{\omega }^{\left(2\right)}{x}^{\left(1\right)}$

and with $x={x}^{\left(1\right)}+{x}^{\left(2\right)}+{x}^{\left(3\right)},\omega ={\omega }_{0}+{\omega }^{\left(2\right)},$ following Landau,

${\stackrel{¨}{x}}^{\left(3\right)}+{\omega }_{0}^{2}{x}^{\left(3\right)}=-\frac{{\alpha }^{2}{a}^{3}}{6{\omega }_{0}^{2}}\mathrm{cos}3\omega t+a\left[2{\omega }_{0}{\omega }^{\left(2\right)}+\frac{5{a}^{2}{\alpha }^{2}}{6{\omega }_{0}^{2}}\right]\mathrm{cos}\omega t.$

Again, there cannot be a nonzero term driving the system at resonance, so the quantity in the square brackets must be zero, this gives us $\Delta \omega ={\omega }^{\left(2\right)}=-5{a}^{2}{\alpha }^{2}/12{\omega }_{0}^{3}.$

The total correction to frequency to leading order for the additional small potentials $\frac{1}{3}m\alpha {x}^{3}+\frac{1}{4}m\beta {x}^{4}$ is therefore (they add independently to this order)

$\Delta \omega =\left(\frac{3\beta }{8{\omega }_{0}}-\frac{5{\alpha }^{2}}{12{\omega }_{0}^{3}}\right){a}^{2}=\kappa {a}^{2},$

(the $\kappa$ here being a convenient notation Landau employs later).

### How Good Are These Approximations?

We have an applet that solves this equation numerically, so it is straightforward to check.

Beginning with the quartic perturbation potential $\frac{1}{4}m\beta {x}^{4},$ Landau finds a frequency correction $\Delta \omega =3\beta {a}^{2}/8{\omega }_{0}.$ Taking a rather large perturbation $a=\beta ={\omega }_{0}=1$ we find from the applet that $\Delta \omega =0.33$ whereas Landau’s perturbation theory predicts $\Delta \omega =\frac{3}{8}=0.375.$ However, if we correct Landau’s denominator (as mentioned above, he pointed out it should be $\omega ,$ but said that  was second-order) the error is very small.

Taking $\alpha =0.3,\text{ }\beta =0.1,\text{ }{\omega }_{0}=0,\text{ }a=1$ the formula gives $\Delta \omega =0.018,$ so less than two percent error, and for amplitude 0.2, the effort is less than 0.1%.

Explore with the applet here.

## Resonance in a Damped Driven Linear Oscillator: A Brief Review

This is just to remind you of what we covered in lecture 18, before we add anharmonic terms in the next section.

The linear damped driven oscillator has equation of motion:

$\stackrel{¨}{x}+2\lambda \stackrel{˙}{x}+{\omega }_{0}^{2}x=\left(f/m\right){e}^{i\gamma t}.$

(Following Landau’s notation here$—$note it means the actual frictional drag force is $2\lambda m\stackrel{˙}{x}$ )

Looking near resonance for steady state solutions at the driving frequency, with amplitude $b$, phase lag $\delta$, that is, $x\left(t\right)=b{e}^{i\left(\gamma t+\delta \right)}$, we find

$b{e}^{i\delta }\left(-{\gamma }^{2}+2i\lambda \gamma +{\omega }_{0}^{2}\right)=\left(f/m\right).$

For a near-resonant driving frequency

$\gamma ={\omega }_{0}+\epsilon$

and assuming the damping to be sufficiently small that we can drop the $\epsilon \lambda$ term along with ${\epsilon }^{2}$, the leading order terms give

$b{e}^{i\delta }=-f/2m\left(\epsilon -i\lambda \right){\omega }_{0},$

so the response, the dependence of amplitude $b$ on driving frequency $\Omega ={\omega }_{0}+\epsilon$ is to this accuracy

$b=\frac{f}{2m{\omega }_{0}\sqrt{{\left(\gamma -{\omega }_{0}\right)}^{2}+{\lambda }^{2}}}=\frac{f}{2m{\omega }_{0}\sqrt{{\epsilon }^{2}+{\lambda }^{2}}}.$

(Note also that the resonant frequency is itself lowered by the damping, another second-order effect we’ll ignore.) The rate of absorption of energy equals the frictional loss.  The friction force $2\lambda m\stackrel{˙}{x}$ on the mass moving at $\stackrel{˙}{x}$ is doing work at an average rate:

$2\lambda m\overline{{\stackrel{˙}{x}}^{2}}=\lambda m{b}^{2}{\gamma }^{2}.$

The half width of the resonance curve as a function of $\gamma$ is given by the damping.  The total area under the curve is independent of damping.

For future use, we’ll write the above equation for the amplitude $b$ in terms of deviation $\epsilon$ from the resonant frequency ${\omega }_{0},$

${b}^{2}\left({\epsilon }^{2}+{\lambda }^{2}\right)=\frac{{f}^{2}}{4{m}^{2}\omega {\text{ }}_{0}^{2}},\text{ }\epsilon =\gamma -{\omega }_{0}.$

## Damped Driven Nonlinear Oscillator: Qualitative Discussion

We now add to the damped driven linear oscillator a positive quartic potential term, giving equation of motion

$\stackrel{¨}{x}+2\lambda \stackrel{˙}{x}+{\omega }_{0}^{2}x=\left(f/m\right)\mathrm{cos}\gamma t-\beta {x}^{3}.$

As mentioned above, for a particle oscillating in this potential $\frac{1}{2}{\omega }_{0}^{2}{x}^{2}+\frac{1}{4}\beta {x}^{4}$ the frequency increases with amplitude: the bigger swings encounter a potential becoming stronger and stronger than the simple harmonic oscillator.

So if we drive the oscillator from rest at the frequency that resonates with its small amplitude oscillations (where the $\frac{1}{4}\beta {x}^{4}$ potential term has negligible effect), as the amplitude builds up, the oscillator frequency increases, and the driving force falls out of sync.

The way to keep the amplitude increasing is evidently to gradually increase the frequency of the driving force to match the natural frequency at the increased amplitude.  (Side note: this is the principle of the synchrocyclotron$—$except, in that case the frequency is lowered as the energy increases, because the particles go to bigger and bigger orbits as their mases increase relativistically.)  This way a small external driving force (enough to overcome frictional damping) can maintain a large amplitude oscillation at a frequency well above the frequency ${\omega }_{0}$ of small oscillations. But what if we apply this high frequency to a system initially at rest, rather than gradually ramping up in sync with the oscillations? Then for a small driving force, we can treat the system as a damped simple harmonic oscillator, and this off-resonance force will drive relatively small amplitude oscillations.

The bottom line is that for the same external driving force, with frequency in some range above ${\omega }_{0}$, there can be two possible steady state oscillation amplitudes, depending on the history of the system.

## Nonlinear Case: Landau’s Analysis

The equation of motion is:

$\stackrel{¨}{x}+2\lambda \stackrel{˙}{x}+{\omega }_{0}^{2}x=\left(f/m\right)\mathrm{cos}\gamma t-\beta {x}^{3}.$

We established earlier that the nonlinear quartic term brings in a correction to the oscillator’s frequency that depends on the amplitude $b$:

$\omega ={\omega }_{0}+\frac{3\beta {b}^{2}}{8{\omega }_{0}}={\omega }_{0}+\kappa {b}^{2}$

in Landau’s notation, $\kappa =3\beta /8{\omega }_{0}.$

The equation for the amplitude in the linear case (from the previous section) was, with $\epsilon =\gamma -{\omega }_{\text{ }0},$

${b}^{2}\left({\epsilon }^{2}+{\lambda }^{2}\right)=\frac{{f}^{2}}{4{m}^{2}{\omega }_{0}^{2}}.$

For the nonlinear case, the maximum amplitude will clearly be at the true (amplitude dependent!) resonance frequency $\omega \left(b\right)={\omega }_{0}+\kappa {b}^{2}$ so with $\epsilon =\gamma -{\omega }_{\text{ }0},$ as before, we now have a cubic equation for ${b}^{2}$:

${b}^{2}\left({\left(\epsilon -\kappa {b}^{2}\right)}^{2}+{\lambda }^{2}\right)=\frac{{f}^{2}}{4{m}^{2}{\omega }_{0}^{2}}$. Note that for small driving force $f\ll 2m{\omega }_{0}\lambda ,$ $b$ is small ( ${b}_{\mathrm{max}}^{2}\approx {f}^{2}/4{m}^{2}{\omega }_{0}^{2}{\lambda }^{2}$ ) but the center of the peak has shifted slightly upwards, to $\epsilon =\kappa {b}^{2},$ that is, at a driving frequency $\gamma ={\omega }_{0}+\kappa {b}^{2}.$  The cubic equation for ${b}^{2}$ has only this one real solution.

However, as the driving force is increased, the coefficients of the cubic equation change and at a critical force ${f}_{k}$ two more real roots appear.

The $b,\epsilon$ curve for driving force above ${f}_{k}$ looks like: So what’s going on here?  For a range of frequencies, including the vertical dashed red line in the figure, there appear to be three possible amplitudes of steady oscillation at one frequency.  However, it turns out that the middle one is unstable, so will exponentially deviate, going to one of the other two, both of which are stable.

If the oscillator is being driven at ${\omega }_{0}$, and the driving frequency is gradually increased, the amplitude will follow the upper curve to the point $C$, then drop discontinuously to the lower curve.  Further frequency increase (with the same strength driving force, of course) will give diminishing amplitude of oscillation$—$just as happens for the ordinary simple harmonic oscillator on going away from the resonant frequency.

If the frequency is now gradually lowered, the amplitude gradually will increase to point $D$, where it will jump discontinuously to the upper curve.  The overall response to driving frequency is sometimes called hysteresis, by analogy with the response of a magnetic material to a varying imposed external field.

To put in some numbers, the maximum amplitude for any of these curves is when  $db/d\epsilon =0,$ that is, at $\epsilon =\kappa {b}^{2},$ or

${b}_{\mathrm{max}}=f/2m{\omega }_{0}\lambda ,$

the same result as for small oscillations.

To find the critical value of the driving force for which the multiple solutions appear, in the graph above that’s when $C,D$ coincide. That is, $db/d\epsilon =\infty$ has coincident roots.

Differentiating the equation $b\left(\epsilon \right)$ for amplitude as a function of frequency (and of course this is at constant driving force $f$ )

$\frac{db}{d\epsilon }=-\frac{-\epsilon b+\kappa {b}^{2}}{{\epsilon }^{2}+{\lambda }^{2}-4\kappa \epsilon {b}^{2}+3{\kappa }^{2}{b}^{4}}.$

$C,D$ coincide when the discriminant in the denominator quadratic is zero, that is, at ${\kappa }^{2}{b}^{4}={\lambda }^{2},$ where $\epsilon =2\kappa {b}^{2}.$

Putting these values into the equation for $b\left(\epsilon \right)$ as a function of driving force $f,$ the critical driving force is

${f}_{\kappa }^{2}=\frac{8{m}^{2}{\omega }_{0}^{2}{\lambda }^{3}}{\left|\kappa \right|}.$

### Numerical Applet Results

The results above are all from Landau’s book, and are semiquantitative.  They can easily be checked using our online applet, which is accurate to one percent or better.  The curves below are plotted from applet results, and certainly exhibit the behavior predicted by Landau.

These plots are for ${\omega }_{0}^{2}=m=\beta =1,\text{ }\text{ }\alpha =0,\text{ }\text{ }2\lambda =0.34.$

For these values, Landau's ${f}_{k}$ is approximately 0.3. Ours looks a bit more. Note for $f=0.3,$ we show $\kappa {b}^{2}=0.38\cdot {\left(0.75\right)}^{2}=0.21,$ close to graph peak.     ## Frequency Multiples

The above analysis is for frequencies not very far from ${\omega }_{\text{​}\text{ }0}$. But nonlinear terms can cause resonance to occur at frequencies which are rational multiples of ${\omega }_{\text{ }0}$.  Landau shows that a small $\frac{1}{3}\alpha {x}^{3}$ in the potential  (so an additional force $\alpha {x}^{2}$ in the equation of motion) can generate a resonance near $\gamma =\frac{1}{2}{\omega }_{\text{ }0}$.  We’ve only considered a quartic addition to the potential, $\frac{1}{4}\beta {x}^{4}$, a force $\beta {x}^{3}$, we can show that gives a resonance near $\gamma =\frac{1}{3}{\omega }_{\text{ }0},$ and presumably this is the small bump near the beginning of the curves above for large driving strength.

We have $\stackrel{¨}{x}+2\lambda \stackrel{˙}{x}+{\omega }_{\text{ }0}^{2}x=\left(f/m\right)\mathrm{cos}\gamma t-\beta {x}^{3}.$

We’ll write $x={x}^{\left(0\right)}+{x}^{\left(1\right)}+\dots$

Let’s define ${x}^{\left(0\right)}$ by

${\stackrel{¨}{x}}^{\left(0\right)}+2\lambda {\stackrel{˙}{x}}^{\left(0\right)}+{\omega }_{\text{ }0}^{2}{x}^{\left(0\right)}=\left(f/m\right)\mathrm{cos}\gamma t$

So ${x}^{\left(0\right)}=b\mathrm{cos}\left(\gamma t+\delta \right)$.  Then

$\begin{array}{c}{\stackrel{¨}{x}}^{\left(1\right)}+2\lambda {\stackrel{˙}{x}}^{\left(1\right)}+{\omega }_{\text{ }0}^{2}{x}^{\left(1\right)}=-\beta {\left({x}^{\left(0\right)}\right)}^{3}\\ =-\beta {b}^{3}{\mathrm{cos}}^{3}\left(\gamma t+\delta \right)\\ =-\beta {b}^{3}\left[\frac{3}{4}\mathrm{cos}\left(\gamma t+\delta \right)+\frac{1}{4}\mathrm{cos}\left(3\gamma t+\delta \right)\right]\end{array}$

Then, for $\gamma =\frac{1}{3}{\omega }_{\text{ }0}$, the second term, $-\beta {b}^{3}\frac{1}{4}\mathrm{cos}\left(3\gamma t+\delta \right)=-\beta {b}^{3}\frac{1}{4}\mathrm{cos}\left({\omega }_{\text{ }0}t+\delta \right)$, will have a resonant response, although it is proportional to the (small) amplitude cubed.

Similar arguments work for other fractional frequencies.