There's one more really important generalization we can add to the damped harmonic oscillator, which is a driving force. The damped oscillator is sort of a special case, because with no energy input the damping will quickly dispose of all the mechanical energy and cause the oscillator to stop. In many realistic applications (like a clock), we want an oscillator that will continue to oscillate, which means an external force has to be applied.
We'll assume the applied force depends on time, but not on the position of the oscillator, so it takes the general form \( F(t) \). Then the ODE becomes
\[ \begin{aligned} m \ddot{x} + b \dot{x} + kx = F(t). \end{aligned} \]
(Again, if you're familiar with circuits, this is like an LRC circuit adding in a source of electromotive force like a battery.)
With the addition of the driving force, our ODE is no longer homogeneous. This means that we no longer have the principle of superposition: we can't find our general solution by just taking any two independent solutions and adding them together with constants.
Instead, we can use an old technique that we've seen for first-order equations. Suppose that we find a function \( x_c(t) \), the complementary solution, that just gives zero when plugged in to the left-hand side of the ODE, i.e.
\[ \begin{aligned} m \ddot{x}_c + b \dot{x}_c + kx_c = 0. \end{aligned} \]
(This equation with the non-homogeneous part thrown out is called the complementary equation.) Next, we find a function \( x_p(t) \) called the particular solution, which satisfies the full equation:
\[ \begin{aligned} m\ddot{x}_p + b\dot{x}_p + kx_p = F(t). \end{aligned} \]
Although the particular solution works, since this is a linear second-order ODE, we know that the general solution requires two undetermined constants. We can't add any constants to \( x_p \) itself; for example, \( 2x_p(t) \) is not a solution since it gives \( 2F(t) \) on the right-hand side. But since the complementary solution gives zero, we can add that in with an unknown constant.
Thus, the full general solution is
\[ \begin{aligned} x(t) = C_1 x_{1,c}(t) + C_2 x_{2,c}(t) + x_p(t) \end{aligned} \]
The good news is that we've already found the complementary solutions, since the complementary equation is just the regular damped oscillator that we just finished solving.
How we find the particular solution depends strongly on what the non-homogeneous part of the equation looks like. Constant is the simplest possibility: let's do a quick example to demonstrate.
Example: non-homogeneous with constant driving
Let's warm up by finding the general solution of the following equation:
\[ \begin{aligned} \ddot{y} - \dot{y} - 6y = 4. \end{aligned} \]
We start with the complementary equation,
\[ \begin{aligned} \ddot{y}_c - \dot{y}_c - 6y_c = 0 \end{aligned} \]
We can solve this using the algebra of differential operators: rewriting the complementary equation gives
\[ \begin{aligned} (D^2 - D - 6) y_c = 0 \end{aligned} \]
which we can factorize into
\[ \begin{aligned} (D+2)(D-3)y_c = 0. \end{aligned} \]
As usual we can write a simple exponential solution for each monomial:
\[ \begin{aligned} y_{1,c}(t) = e^{-2t},\ \ y_{2,c}(t) = e^{+3t}. \end{aligned} \]
On to the particular solution. We always want to look for a particular solution that has a similar form to the driving. A constant particular solution is especially nice, since all of the derivatives vanish immediately. It should be easy to see right away that
\[ \begin{aligned} y_p(t) = -\frac{2}{3} \end{aligned} \]
will satisfy the equation. Thus, putting everything together, we have for the general solution
\[ \begin{aligned} y(t) = C_1 e^{-2t} + C_2 e^{+3t} - \frac{2}{3}. \end{aligned} \]
Coming back to the specific case of the oscillator, suppose that wehave a constant driving force \( F(t) = F_0 \). Rewriting the equation into our standard form,
\[ \begin{aligned} \ddot{x} + 2\beta \dot{x} + \omega_0^2 x = \frac{F_0}{m} \end{aligned} \]
so to our damped oscillator solutions from before, we simply have to add the particular solution
\[ \begin{aligned} x_p(t) = \frac{F_0}{m \omega_0^2}. \end{aligned} \]
This is just a constant offset to our \( x \) coordinate. So, for example, if we take an oscillator moving horizontally and reorient it so that it's hanging vertically subject to gravity, the motion will be exactly the same except that the equilibrium position will be displaced to
\[ \begin{aligned} x_p = \frac{g}{\omega_0^2} = \frac{mg}{k} \end{aligned} \]
(which is exactly the new equilibrium point where the spring force \( kx_p \) balances the gravitational force \( mg \).)
Let's do a more interesting example of a driving force.
Complex exponential driving
Now we'll assume a driving force of exponential form:
\[ \begin{aligned} F(t) = F e^{ct} \end{aligned} \]
where both \( F \) and \( c \) can be complex. This covers both regular exponential forces (if \( c \) is real), as well as the very important case of oscillating driving force (if \( c \) is imaginary.)
The obvious thing to try for our particular solution is exactly the same form: we'll guess that
\[ \begin{aligned} x_p(t) = C_0 e^{ct} \end{aligned} \]
will give us a solution, and then try it out. Note that although I put in an unknown constant \( C_0 \), as we've seen with other particular solutions this is not an arbitrary constant; we expect only a specific choice of \( C_0 \) to work given a certain driving force.
Let's plug our guessed solution back in and see what we get: from
\[ \begin{aligned} \ddot{x} + 2\beta \dot{x} + \omega_0^2 x = \frac{F}{m} e^{ct} \end{aligned} \]
we have
\[ \begin{aligned} \left[ c^2 + 2\beta c + \omega_0^2 \right] C_0 e^{ct} = \frac{F}{m} e^{ct} \end{aligned} \]
or cancelling off the exponentials,
\[ \begin{aligned} C_0 = \frac{F}{m(c^2 + 2\beta c + \omega_0^2)}. \end{aligned} \]
Once again, this is not arbitrary! Both \( c \) and \( F \) are given to us as part of the external driving force: they are literally numbers that are dialed in on whatever machine is providing the driving for our experiment.
Let's step back for a moment and look at the big picture. Remember that the complementary solutions take the form \( e^{r_{\pm} t} \), where \( r_{\pm} \) solve the equation
\[ \begin{aligned} r^2 + 2\beta r + \omega_0^2 = 0. \end{aligned} \]
The full solution is two complementary solutions with arbitrary coefficients, plus one particular:
\[ \begin{aligned} x(t) = C_+ e^{r_+ t} + C_- e^{r_- t} + C_0 e^{ct} \end{aligned} \]
with \( C_0 \) determined by \( F \) and \( c \) as we found above. All of these exponents can be complex, so we can have individual contributions that oscillate at one, two, or three different frequencies; we'll make more sense of this soon.
Notice that there's an important edge case here: the same quadratic equation for \( r \) appears again in the denominator of \( C0 \). So if we have \( c = r{\pm} \), then our formula for \( C_0 \) blows up. This is similar to what happened in the critical damping case: our would-be particular solution fails, because it's not independent from the complementary solution. The way to fix this, similarly, is to take a particular solution of the form \( C_0 te^{ct} \) or \( C_0 t^2 e^{ct} \) as needed (we only need the second one if the system is critically damped, so \( te^{ct} \) is already being used.)
Sinusoidal driving
So far this has been totally general, but real exponential driving forces are sort of an unusual special case; we won't really be interested in them. So let's specialize now to the very important case of an oscillating driving force, which we'll write as follows:
\[ \begin{aligned} F(t) = F e^{ct} = F_0 e^{-i\delta_0} e^{i \omega t} \end{aligned} \]
so that \( c = i \omega \), and we've split the complex force \( F \) into two real numbers \( F_0 \) and \( \delta_0 \). Plugging all of this back in to \( C_0 \) and splitting into real and imaginary parts in the denominator gives
\[ \begin{aligned} C_0 = \frac{F_0 e^{-i \delta_0}/m}{(\omega_0^2 - \omega^2) + 2i \beta \omega}. \end{aligned} \]
It will be useful to rewrite this in polar form, \( C_0 = Ae^{-i\Delta} \). To find the amplitude \( A \) first, we multiply by the complex conjugate:
\[ \begin{aligned} |C_0|^2 = C_0 C_0^\star \\ = \frac{F_0e^{-i\delta_0}/m}{(\omega_0^2 - \omega^2) - 2i \beta \omega} \frac{F_0e^{+i \delta_0}/m}{(\omega_0^2 - \omega^2) + 2i \beta \omega} \\ = \frac{F_0^2/m^2}{(\omega_0^2 - \omega^2)^2 + 4\beta^2 \omega^2}, \end{aligned} \]
so the amplitude is
\[ \begin{aligned} A = |C_0| = \frac{F_0/m}{\sqrt{(\omega_0^2 - \omega^2)^2 + 4\beta^2 \omega^2}}. \end{aligned} \]
Next, the phase \( \Delta \). We can simplify finding this by using the fact that the force is already in polar form, i.e. we can rewrite
\[ \begin{aligned} C_0 = Ae^{-i\Delta} = Ae^{-i(\delta_0 +\delta)} = Ae^{-i\delta_0} e^{-i\delta} \end{aligned} \]
defining the angle \( \delta \) to be the phase shift, i.e. the difference in phase between the input force phase \( \delta_0 \) and the particular solution phase \( \Delta \). This lets us pull out the \( \delta_0 \) from \( C_0 \) like so:
\[ \begin{aligned} C_0 e^{+i\delta_0} = Ae^{-i\delta} \\ = \frac{F_0 /m}{(\omega_0^2 - \omega^2) + 2i \beta \omega} \left[\frac{(\omega_0^2 - \omega^2) - 2i \beta \omega}{(\omega_0^2 - \omega^2) - 2i \beta \omega}\right] \\ = \frac{F_0/m}{(\omega_0^2 - \omega^2)^2 + 4\beta^2 \omega^2} \left[ (\omega_0^2 - \omega^2) - 2i \beta \omega \right]. \end{aligned} \]
where I've used the conjugate of the denominator to split this expression into real and imaginary parts.
Clicker Question
What is the phase angle \( \delta \) in this equation?
\[ \begin{aligned} Ae^{-i\delta} = \frac{F_0/m}{(\omega_0^2 - \omega^2)^2 + 4\beta^2 \omega^2} \left[ (\omega_0^2 - \omega^2) - 2i \beta \omega \right]. \end{aligned} \]
A. \( \tan \delta = \frac{2\beta \omega}{\omega_0^2 - \omega^2} \)
B. \( \tan \delta = \frac{-2\beta \omega}{\omega_0^2 - \omega^2} \)
C. \( \tan \delta = \frac{\omega_0^2 - \omega^2}{2\beta \omega} \)
D. \( \tan \delta = \frac{2\beta \omega}{\sqrt{\omega_0^2 - \omega^2 + 4\beta^2 \omega^2}} \)
E. \( \tan \delta = \frac{\omega_0^2 - \omega2}{\sqrt{\omega_02 - \omega^2 + 4\beta^2 \omega^2}} \)
_Answer:_ A
If we're interested in finding \( \delta \) the factor in front of the square brackets doesn't matter, since the angle only depends on the relative size of real and imaginary parts. The best way to see the answer is by making a quick sketch in the complex plane:
so we can read off that the complex phase is given by
\[ \begin{aligned} \tan \delta = \frac{2\beta \omega}{\omega_0^2 - \omega^2}. \end{aligned} \]
taking out a minus sign since we wrote this as \( Ae^{-i\delta} \) to begin with. (Or equivalently, the minus sign in \( Ae^{-i\delta} \) accounts for the fact that the imaginary part is negative, so the angle \( \delta \) in our sketch should be positive.)