Simple Harmonic Motion

We will start by reviewing simple harmonic motion (SHM) as it contains many of the important concepts that we will meet in waves. The most general form of the equation for a simple harmonic oscillator (SHO) including damping and driving forces can be written:

$$ m\frac{\partial^2 \psi}{\partial t^2} = -s\psi -b \frac{\partial \psi}{\partial t} + F_0 \cos \omega t, $$

where $\psi$ is the displacement, $m$ is the mass of the oscillator, $s$ is the constant giving the restoring force (e.g. a spring constant) sometimes known as the stiffness, $b$ is the damping coefficient or resistance and $F_0$ gives the magnitude of the driving force (which has angular frequency $\omega$).

SHM and Circular Motion

If we set the damping and driving coefficients to zero, we recover the original, SHM equation:

$$m\frac{\partial^2 \psi}{\partial t^2} = -s \psi,$$

which can be shown to be solved using sinusoidal motion: $$ \psi = A \sin \omega_0 t + B \cos \omega_0 t = C \cos \left(\omega t + \phi\right), $$ where $\omega_0 = \sqrt{s/m}$, $\phi$ is a constant phase and $A = -C\sin \phi$ and $B = C \cos \phi$. As you have seen in both PHAS1245 (Mathematical Methods I) and PHAS1247 (Classical Mechanics), we can use De Moivre's theorem to write the sinusoidal terms as a complex exponential ($e^{i\theta} = \cos \theta + i \sin \theta$ where $i = \sqrt{-1}$), giving: $$ \psi = Re\left[A e^{i(\omega_0 t + \phi)}\right] = Re\left[D e^{i\omega_0 t}\right] $$ where $D = Ae^{i\phi}$.

We often call the argument of a trigonometric function the phase, and $\phi$ is often called the phase difference when there is more than one oscillation; it represents the offset in the phase at $t=0$. Now that we have the representation of the oscillation in terms of the complex exponential (or the sum of a $\cos$ and $\sin$) we can see that there is an immediate link with circular motion: the amplitude of the oscillation is just the projection of circular motion onto the $x$-axis (or any other axis that is chosen). (Note that from now on I will not write the need to take the real part of $\psi$ explicitly, but will assume it.) All these ideas are illustrated in the figure below.

In [4]:
from IPython.display import Image
Image(filename='SHMCircularMotion.png',width=300) #Size is in pixels

We can also illustrate this more dynamically, with the python code below. The blue line shows the phasor at $t=0$ (with the angle relative to the x-axis being $\phi$. The red line shows the phasor at a later time $t$, with the real part shown as a dashed red line along the x-axis.

In [16]:
phi = 30.0*pi/180. # Phase in degrees
rad = 1.0
x = (0,rad*cos(phi))
y = (0,rad*sin(phi))
# Add arc to indicate phi
theta = linspace(0,phi,1000)
# Label
text(0.5, 0.18, r"$\phi$", fontsize=20, color="blue")
# Now draw full circle
theta = linspace(0,2.0*pi,1000)
# Plot the phasor
omega = 2.0
t = 2.5
psix = (0,rad*cos(omega*t+phi))
psiy = (0,rad*sin(omega*t+phi))
# Projection on the x-axis
# Label
text(0.4,-0.2, r"$\psi$", fontsize=20, color="red")
# Arc showing omega t
theta = linspace(phi,omega*t+phi,1000)
# Label
text(-0.5, 0.6, r"$\omega t$", fontsize=20, color="red")
<matplotlib.text.Text at 0x110c412d0>

The representation of a simple harmonic oscillator at a single point in time in the complex plane (e.g. figure above) is often called a phasor diagram, and the arrow representing the amplitude and phase of the oscillator relative to a fixed time or phase is called a phasor; note that for a phasor we just need the arrow, not the associated circle. We will see later that we can combine two oscillations using phasors (or using complex arithmetic - the two are completely equivalent).

A simple way to think about phase differences is to consider the velocity and acceleration of the oscillator:

\begin{eqnarray} \psi &=& A e^{i(\omega_0 t + \phi)}\\ \dot{\psi} &=& \frac{\partial \psi}{\partial t} = i\omega_{0} A e^{i(\omega_0 t + \phi)} = i\omega_{0} \psi\\ \ddot{\psi} &=& \frac{\partial^2 \psi}{\partial t^2} = -\omega_{0}^2 A e^{i(\omega_0 t + \phi)} = -\omega_{0}^2 \psi \end{eqnarray}

Notice that both of these quantities vary harmonically and with the same frequency as the oscillation, though with different amplitudes; more importantly, however, note that the velocity is a factor of $i$ different to the displacement, while the acceleration is a factor of $-1$ different. These are easier to understand when we write them in complex exponential form: $i = e^{i\pi/2}$ and $-1 = e^{i\pi}$. So we can write:

\begin{eqnarray} \psi &=& A e^{i(\omega_0 t + \phi)}\\ \dot{\psi} &=& \omega A e^{i(\omega_0 t + \phi + \pi/2)} = \omega_{0} \psi e^{i\pi/2}\\ \ddot{\psi} &=& \omega^2 A e^{i(\omega_0 t + \phi + \pi)} = \omega_{0}^2 \psi e^{i\pi} \end{eqnarray}

We see that the velocity leads the displacement by a phase factor of $\pi/2$ and the acceleration leads the velocity by a phase factor of $\pi/2$. This phase relationship is shown in the next figure.

In [3]:
from IPython.display import Image
Image(filename='SHMDispVelAcc.png',width=300) #Size is in pixels

The simple harmonic oscillator, as with all mechanical systems, has two forms of energy: kinetic and potential. If we have $W$ as the total energy and $T$ as kinetic and $V$ as potential, we can write:

\begin{eqnarray} T &=& \frac{1}{2} m\dot{\psi}^2\\ V &=& \frac{1}{2} s\psi^2\\ W &=& T + V \end{eqnarray}

The form of the potential energy follows directly from the form of the force (and is why the motion is called harmonic). As there is no form of dissipation in the motion, we can assert that the total energy does not change with time, just exchanging between potential (at a maximum when the displacement is at a maximum and the velocity is at a minimum) and kinetic (at a maximum when the displacement is at a minimum). They are out of phase - as we expect from the phase differences between displacement and velocity.

As the total energy is constant, we can write:

\begin{eqnarray} \frac{d W}{d t} = \frac{d T}{d t} + \frac{d V}{dt} &=& 0\\ \Rightarrow m\dot{\psi}\ddot{\psi} + s\psi\dot{\psi} &=&0\\ \Rightarrow m\ddot{\psi} &=& -s\psi \end{eqnarray}

which is of course just the original equation for simple harmonic motion.

Damping Oscillations

Restoring the damping term changes the equation and its solution. We will use $\gamma = b/2m$ in the equations below. These solutions were derived in PHAS1247 and are discussed in many textbooks, so I will not derive them here. We find:

\begin{eqnarray} m\frac{\partial^2 \psi}{\partial t^2} &=& -s \psi -b \frac{\partial \psi}{\partial t}\\ \psi &=& \begin{cases} A e^{-\gamma t}e^{i\omega t} &\gamma < \omega_0\\ \hskip 0.5cm\omega = \sqrt{\omega^2_0 - \gamma^2} &\, \\ A e^{-\mu_+ t} + B e^{-\mu_-t} & \gamma > \omega_0\\ \hskip 0.5cm\mu_\pm = \gamma\mp \sqrt{\gamma^2 - \omega_0^2} &\, \\ A(1+\omega_0 t)e^{-\omega_0 t} & \gamma = \omega_0 \end{cases} \end{eqnarray}

The cases for damped harmonic motion correspond to light (or underdamped, $\gamma < \omega_0$), heavy (or overdamped, $\gamma > \omega_0$) and critical ($\gamma=\omega_0$) damping respectively. The total energy is \emph{not} conserved in this case (contrast the undamped oscillator) as the damping force opposes the motion at all times. We write the energy as before, and find the change with time:

\begin{eqnarray} W &=& \frac{1}{2}m\dot{\psi}^2 + \frac{1}{2}s\psi^2\\ \frac{dW}{dt} &=& \frac{d W}{d \dot{\psi}}\frac{d \dot{\psi}}{d t} + \frac{d W}{d \psi}\frac{d \psi}{d t} = (m\ddot{\psi} + s\psi)\dot{\psi}\\ &=& -b\dot{\psi}^2 \end{eqnarray}

where the last line comes from using the SHM equation for damped motion. Notice that the change in energy is always less than zero (unless $\dot{\psi}=0$) and so the total energy of the system decreases with time (as we would expect).

We illustrate this behaviour in the interactive figure below: test the effect of changing b and k (and re-run the cell to render the figure).

In [27]:
t = linspace(0,10.0,1000)
m = 1.0
k = 2.0
b = 5.0
omega0 = sqrt(k/m)
gamma = b/(2.0*m)
print "omega0 and gamma are ",omega0,gamma
A = 1.0
if gamma<omega0 :
    print "Light damping"
    omega = sqrt(omega0*omega0 - gamma*gamma)
elif gamma>omega0 : 
    print "Heavy damping"
    mum = gamma + sqrt(gamma*gamma - omega0*omega0)
    mup = gamma - sqrt(gamma*gamma - omega0*omega0)
    B = 1.0
    plot(t,A*exp(-mum*t) + B*exp(-mup*t),label='sum')
omega0 and gamma are  1.41421356237 2.5
Heavy damping

Driving Oscillations

Now we will introduce a driving term to the system. As we have a harmonic oscillator, it will make sense to use a harmonic driving term (though we could, for instance, use impulses at regular or irregular intervals). Note that the driving term will have an angular frequency $\omega$ which is \emph{different} to the natural frequency of the system, $\omega_0 = \sqrt{s/m}$. Also remember that $\omega = 2\pi \nu$ where $\nu$ is the frequency. We will retain the damping term to give a damped, driven oscillator (again derived in PHAS1247 and textbooks).

\begin{eqnarray} m\frac{\partial^2 \psi}{\partial t^2} &=& -s \psi -b \frac{\partial \psi}{\partial t} + F_0 \cos \omega t\\ \psi &=& A\cos (\omega t + \phi)\\ A &=& \frac{F_0}{m} \left( \frac{1}{(\omega_0^2 - \omega^2)^2 + 4\gamma^2 \omega^2}\right)^{\frac{1}{2}}\\ \tan \phi &=& \frac{-2\gamma\omega}{\omega_0^2 - \omega^2} \end{eqnarray}

The response of the system (including both the phase and the amplitude) depends strongly on the frequency; we can consider three regimes (though there isn't really time to do this in depth).

  • When $\omega$ is small (i.e. $\omega \ll \omega_0$), the amplitude can be shown to be $A \simeq F_0/m\omega^2_0 = F_0/s$ and the motion is dominated by the spring constant (stiffness controlled).
  • When $\omega$ is large (i.e. $\omega \gg \omega_0$), then the amplitude is $A \simeq F_0/m\omega^2$ and the motion is dominated by the mass (mass controlled).
  • When $\omega \sim \omega_0$, we are at resonance, and the response is controlled by the drag term $\gamma$ (also known resistance limited).

The behaviour of $A$ and $\phi$ are plotted below; you should try changing the values of b (damping) and s (stiffness) to see the effects. In particular, you may need to increase or decrease F0 to see the graph on a sensible scale. (Remember to run the cell by pressing the play button in the toolbar above.)

In [30]:
# Set up parameters of oscillator
m = 1.0
k = 2.0
b = 0.2
gamma = b/(2.0*m)
# Find resonant frequency
omega0 = sqrt(k/m)
print "omega0 and gamma are ",omega0,gamma
F0 = 1.0
# Create array of values of omega from omega0-1 to omega0+1
omega = linspace(omega0-1.0,omega0+1.0,200)
# Calculate phi
phi = arctan(-2.0*gamma*omega/(omega0*omega0 - omega*omega))
# arctan in NumPy is defined to lie between -pi/2 and pi/2 so we correct
# by -pi to give physical reasonableness
phi[100:200] -=pi
# Calculate amplitude
fac = (omega0*omega0 - omega*omega)**2 + 4*gamma*gamma*omega*omega
A = F0/(sqrt(fac)*m)
omega0 and gamma are  1.41421356237 0.1
<matplotlib.legend.Legend at 0x110cba8d0>

Power absorption

Let's think a little about the power absorbed by the oscillator; to do that, we must consider the work done against the drag. If the displacement changes from $\psi$ to $\psi + \Delta \psi$ then the work done is $-F_d \Delta\psi$, where $F_d = -b\dot{\psi}$ is the work done against the drag. If that takes a time $\Delta t$ then the rate of work is $-F_d (\Delta \psi/\Delta t)$ which tends to $-F_d \dot{\psi}$ as $\Delta t \rightarrow 0$. So the instantaneous power adsorption becomes: \begin{eqnarray} P &=& -F_d \dot{\psi} = b\dot{\psi}^2 \end{eqnarray}

We are not actually interested in the instantaneous power, but the time averaged power (often written $\langle P\rangle$), which for a harmonic force can be shown to be: \begin{eqnarray} \langle P \rangle &=& b\langle \dot{\psi}^2 \rangle = \frac{1}{2}b\omega^2A^2\\ &=& \frac{F_0^2}{m^2}\frac{m\gamma\omega^2}{(\omega_0^2 - \omega^2)^2 + 4\gamma^2\omega^2} \end{eqnarray} where we have substituted $b = 2m\gamma$. This has a maximum value of $F^2_0/(4m\gamma)$ when $\omega = \omega_0$. Note that this is inversely proportional to the resistive force.

Finally we introduce the idea of impedance which is a measure of the resistance to the motion of the oscillator. It is defined as the amplitude of the driving force divided by the complex amplitude of the velocity, $Z(\omega) = F_0/\dot{\psi} = b + i(m\omega - s/\omega)$. At resonance, $Z(\omega) = b$; in the resonance region, it only departs very slightly from this value. Away from resonance, it includes a phase lag or lead, which reflects the relation of velocity to the driving force. We will encounter impedance again with waves.

We plot below the power absorbed, imaginary part of the impedance and the magnitude of the impedance, using the system given above. Note how the impedance has a minimum at the point that the maximum power is absorbed (and is purely real at that point). Away from that, the power is mainly dissipated (by the impedance).

In [35]:
# We use the values of b, m, k, omega and A defined above
power = 0.5*b*omega*omega*A*A
impi = m*omega - k/omega
Z = sqrt(b*b+impi*impi)
plot(omega,impi,label='Imaginary part of impedance')
plot(omega,Z,label=r'$\vert Z\vert$')
<matplotlib.legend.Legend at 0x110eac990>


Note that so far we have discussed a steady state solution: there is a transient behaviour at the beginning of the oscillation which takes the form of a damped oscillator added to a solution of the homogenous equation. So a real system will respond at its natural frequency, but damping will cause that solution to die off, while the driven oscillation persists into the steady state:

\begin{equation} \psi(t) = A_s \cos(\omega t +\phi) + A_t e^{-\gamma t} \cos(\omega_0 t + \phi_t) \end{equation} where the subscripts s and t stand for steady state and transient respectively. For lightly damped system near resonance we can approximate and find that $A_t \simeq A_s$ and $\phi_t \simeq \phi_s + \pi$, giving: \begin{equation} \psi(t) = A_s \left(\cos(\omega t +\phi) - e^{-\gamma t} \cos(\omega_0 t + \phi_s)\right) \end{equation}

We plot an example of this below, with the total response of the system over time in blue, and the transient in green: notice how the initial response of the system is at its natural frequency, $\omega_0$, and how this dies down and leads to the steady-state, driven response at $\omega$. You should try changing the value of the driving frequency omegaD (though you will also have to change the amplitude of the driving force, F0, and the initial response, F1).

In [44]:
t = linspace(0,50.0,1000)
# Driving frequency
omegaD = 7.0
# Amplitude of driving (F0) and response (F1)
F0 = 10.0
F1 = 0.10
fac = (omega0*omega0 - omegaD*omegaD)**2 + 4*gamma*gamma*omegaD*omegaD
A = F0/(sqrt(fac)*m)
# Define transient
transient = 2*(F1/m)*exp(-gamma*t)*cos(omega0*t+pi/2)
steady = A*cos(omegaD*t+pi/2)
[<matplotlib.lines.Line2D at 0x10da39e50>]

Combining Oscillations

Instead of combining a single driving force with a damped response, what happens if we have an oscillator with two driving forces ? The full equation looks like this:

$$ \begin{equation} m\frac{\partial^2 \psi}{\partial t^2} = -s \psi -b \frac{\partial \psi}{\partial t} + F_1 \cos (\omega_1 t + \phi_1) + F_2 \cos (\omega_2 t + \phi_2), \end{equation} $$

We know that the SHO equation is linear, so we can try solving this problem by summing the steady state solutions arising from each driving force independently. We will consider three specific cases (where it is easier to understand the behaviour) with the general case more complex:

  • Same amplitude and frequency ($F_1 = F_2 = F$, $\omega_1 = \omega_2 = \omega$), different phase ($\phi_1 \ne \phi_2$)
  • Same frequency ($\omega_1 = \omega_2$), different phase and amplitude ($F_1 \ne F_2, \phi_1 \ne \phi_2$)
  • Same amplitude and phase ($F_1 = F_2=F, \phi_1 = \phi_2$), different frequency ($\omega_1 \ne \omega_2$)

In all cases, we assume that we can write the displacement as a superposition of the responses to the individual oscillations: $$ \begin{equation} \psi = A_1 \cos(\omega_1 t + \phi_1) + A_2 \cos(\omega_2 t + \phi_2) \end{equation} $$

Same amplitude and frequency

For simplicity, we will set $\phi_1 = 0$ and set $\phi_2 = \phi$. Then we can perform the following manipulations, using the complex exponential form for simplicity: $$ \begin{eqnarray} \psi &=& Ae^{i\omega t} + Ae^{i(\omega t + \phi)}\\ &=& Ae^{i\omega t}\left(1 + e^{i\phi}\right)\\ &=& Ae^{i\omega t}\left(e^{i\phi/2}e^{-i\phi/2} + e^{i\phi/2}e^{i\phi/2}\right)\\ &=& Ae^{i\omega t}e^{i\phi/2}\left(e^{-i\phi/2} + e^{i\phi/2}\right)\\ &=& 2Ae^{i(\omega t+\phi/2)}\cos(\phi/2) = 2A\cos(\phi/2)e^{i(\omega t+\phi/2)} \end{eqnarray} $$ where the identification of $\cos (\phi/2)$ follows from De Moivre's theorem. So the resultant oscillation has magnitude $2A\cos(\phi/2)$ and phase $\phi/2$. The same result can be found using trigonometry and a phasor diagram

In [28]:
from IPython.display import Image

Same frequency

We can write a solution as before, by using superposition:

$$ \begin{equation} \psi = A_1 \cos(\omega t + \phi_1) + A_2 \cos(\omega t + \phi_2) = A \cos(\omega t + \theta) \end{equation} $$ This is illustrated in a phasor diagram below.

In [27]:
from IPython.display import Image

This is just an oscillation with frequency $\omega$ but a total amplitude and phase which depend on the amplitudes and phases of the two driving forces. If, for instance, the phase difference is $\phi_1

  • \phi_2 = \pi$ then the two driving forces are out-of-phase and we will get destructive interference. The use of phasors allows a simple visualisation of the resultant. Using trigonometry from the phasor diagram (or the cosine rule), we can write the amplitude of the resultant as: $$ \begin{equation} A^2 = A_1^2 + A^2_2 + 2A_1A_2\cos(\phi_2 - \phi_1) \end{equation} $$ Now we can see that the amplitude will vary between $A_1+A_2$ if $\phi_1=\phi_2$ and $|A_1 - A_2|$ if $|\phi_1 - \phi_2| = \pi$. If the driving forces are in phase then we have the maximum amplitude, while if they are in antiphase we have the minimum amplitude.

The phase is also found using trigonometry; by projecting onto thereal and imaginary axes, we can write: $$ \begin{equation} \tan \theta = \frac{A_1\sin \phi_1 + A_2\sin \phi_2}{A_1 \cos \phi_1 + A_2 \cos \phi_2} \end{equation} $$ We can get the same result using complex notation. Start by finding that, for two general complex numbers $z_1$ and $z_2$: $$ \begin{eqnarray} \left| z_1 + z_2\right|^2 &=& (z_1+z_2)(z_1 + z_2)^\star = (z_1+z_2)(z^\star_1+z^\star_2)\\ &=& |z_1|^2 + |z_2|^2 + (z_1z^\star_2 + z^\star_1z_2)\\ &=& |z_1|^2 + |z_2|^2 +2\mathrm{Re}(z_1z^\star_2) \end{eqnarray} $$ Now we have $z_1 = A_1 e^{i(\omega t +\phi_1)}$ and $z_2 = A_2 e^{i(\omega t +\phi_2)}$ with $A_1$ and $A_2$ real. This gives: $$ \begin{eqnarray} \left| z_1 + z_2\right|^2 &=& |z_1|^2 + |z_2|^2 +2\Re(z_1z^\star_2)\\ &=& A_1^2 + A_2^2 + 2\mathrm{Re}(A_1 e^{i(\omega t +\phi_1)}A_2 e^{-i(\omega t +\phi_2)})\\ &=& A_1^2 + A_2^2 + 2\mathrm{Re}(A_1A_2 e^{i(\phi_1 - \phi_2)})\\ &=& A_1^2 + A_2^2 + 2A_1A_2 \cos(\phi_2 - \phi_1) \end{eqnarray} $$ while we can find the phase from: $$ \begin{eqnarray} \mathrm{arg}(z_1+z_2) &=& \tan^{-1} \left[ \frac{\mathrm{Im}(z_1+z_2)}{\mathrm{Re}(z_1+z_2)} \right]\\ \Rightarrow \theta &=& \tan^{-1} \left[\frac{A_1\sin \phi_1 + A_2\sin \phi_2}{A_1 \cos \phi_1 + A_2 \cos \phi_2}\right] \end{eqnarray} $$

Same amplitude and phase: beats

For this case, we use the principle of superposition to write: $$ \begin{equation} \psi = A \cos \omega_1 t + A \cos \omega_2 t \end{equation} $$ Note that we have assumed that the common phase can be set to zero for simplicity. But we can rearrange this using a standard trigonometrical formula for the sum of two cosines: $$ \begin{eqnarray} \psi &=& A \left(\cos \omega_1 t + \cos \omega_2 t\right)\\ &=& 2A \cos\left(\frac{\omega_1 + \omega_2}{2} t\right) \cos\left(\frac{\omega_1 - \omega_2}{2} t\right) \end{eqnarray} $$

Now let's define two \emph{new} angular frequencies, and rewrite the solution: $$ \begin{eqnarray} \omega &=& \frac{1}{2}\left(\omega_1 + \omega_2\right)\\ \Delta \omega &=& \frac{1}{2}\left(\omega_1 - \omega_2\right)\\ \psi &=& 2A\cos (\omega t) \cos (\Delta\omega t) \end{eqnarray} $$ If the two angular frequencies are relatively close, then the form of the resulting oscillation is rather simple. There is an oscillation at the average frequency (the term $\cos \omega t$) whose amplitude is \emph{modulated} by a slow oscillation at the difference frequency (the term $\cos \Delta\omega t$ - also known as the envelope). This is a phenomenon known as \emph{beats}, and is illustrated in the figure below. It is important to understand that, while the angular frequency of the modulation is $\Delta \omega = \frac{1}{2}|\omega_1 - \omega_2|$, the angular frequency at which \emph{peaks} of activity occur is $|\omega_1 - \omega_2|$ (or equivalently zeroes of activity). The number of minima per second is $\Delta \omega/\pi$. The perceived effect (say for sound) will be of a sound at the average angular frequency with its amplitude varying according to the envelope.

In [48]:
fig = figure(figsize=[12,3])
# Add subplots: number in y, x, index number
ax = fig.add_subplot(121,autoscale_on=False,xlim=(0,65),ylim=(-2,2))
ax2 = fig.add_subplot(122,autoscale_on=False,xlim=(0,65),ylim=(-2,2))
t = linspace(0.0,65.0,1000)
o1 = 1.3
o2 = 1.0
omega = 0.5*(o1+o2)
domega = 0.5*(o1-o2)
# Plot the resulting pattern
# Show the envelope in a dashed red line - we need two lines
# On a separate plot, show the two individual oscillations
[<matplotlib.lines.Line2D at 0x10e206690>]

The result of two driving forces on one SHO with different frequencies but the same amplitude. Left: two forces shown together. Right: the resultant motion with the envelope superimposed in dashed lines. The original frequencies are $\omega = 1.0$ and $\omega=1.2 s^{-1}$; on ce how the peaks coincide at regular intervals, and this leads to the maxima in the envelope function.

The figure above shows an illustration of exactly this behaviour for the frequencies $\omega_1 = 1.2 s^{-1}$ and $\omega_2 = 1.0 s^{-1}$ which gives a resulting motion at frequency $\omega = 1.1 s^{-1}$ modulated by an envelope with frequency $\Delta \omega = 0.1 s^{-1}$.

The resulting motion will show a true periodicity if the ratio of the frequencies can be written as a ratio of integers (i.e. $\omega_1/\omega_2 = n_1/n_2$). The motion from two frequencies which are almost non-periodic is shown in the figure below (I used 25/3 and 7.1 here but we could have made it properly non-periodic if we'd tried harder).

Note that if there is a phase difference between the two oscillations then beats still result, but the phase difference will appear in both sum and difference phases. We can write $\cos(\omega_1 t + \phi_1) + \cos(\omega_2 t + \phi_2) = 2\cos(\omega t + \phi)\cos(\Delta \omega t + \Delta \phi)$ where $\omega$ and $\Delta \omega$ have the same meaning as before and $\phi = (\phi_1 + \phi_2)/2$ and $\Delta \phi = (\phi_1 - \phi_2)/2$.

In [24]:
fig = figure(figsize=[12,3])
# Add subplots: number in y, x, index number
xmax = 15.0
ax = fig.add_subplot(121,autoscale_on=False,xlim=(0,xmax),ylim=(-2,2))
ax2 = fig.add_subplot(122,autoscale_on=False,xlim=(0,xmax),ylim=(-2,2))
t = linspace(0.0,xmax,1000)
o1 = 25./3
o2 = 7.1
omega = 0.5*(o1+o2)
domega = 0.5*(o1-o2)
# Plot the resulting pattern
# Show the envelope in a dashed red line - we need two lines
# On a separate plot, show the two individual oscillations
[<matplotlib.lines.Line2D at 0x10d231810>]

General case

If the amplitudes and frequencies all differ, then the phasor diagram must be dynamic as illustrated below.

In [21]:
from IPython.display import Image
Image(filename='PhasorsDiffFreq2.png',width=200) #Size is in pixels

The tip of the resultant will trace out a shape in time called a cycloid. For instance, the path followed by the tip of the resultant vector in the complex plane for different amplitudes and frequencies differing by a factor of two is plotted below - you should try playing with the parameters, and seeing what effects they have.

In [15]:
# Define figure size to get square
fig = figure(figsize=[5,5])
# The range of t may need to change with angular frequencies
t = linspace(0,10*pi,200)
# Original frequencies: 1.2 and 0.6
o1 = 1.2
o2 = 0.6
# Original amplitudes: 1.5 and 1.3
A1 = 1.5
A2 = 1.3
# Calculate x and y projections of the two oscillations when added
x = A1*cos(o1*t) + A2*cos(o2*t)
y = A1*sin(o1*t) + A2*sin(o2*t)
# Plot path of the tip of the vector
[<matplotlib.lines.Line2D at 0x10cd57c10>]