Second Order Differential Equations

Return to the Math 2250 Overview Page.


Homogeneous Second Order Equations

Second-order differential equations relate a function to its first and second derivatives. Second-order differential equations are related to kinematic equations in physics and engineering. A general second-order differential equation is usually written as \[my''+\gamma y' +ky = f(t)\]
This form of the equation describes a damped harmonic oscillator with a forcing function. Think of it like a spring with a hanging mass on it that when excited will bounce up and down until the oscillation damps out. If used, the \(\textit{m}\) term is the mass on the spring, \(\gamma\) is the damping coefficient, and k is the spring constant. As with most kinematic expressions, \(\textit{y}\) is the position of the mass, \(\textit{y'}\) is the velocity, and \(\textit{y''}\) is the acceleration of the system. In a homogeneous system, \(f(t)=0\)

Charactersitic Polynomial

Given the second-order homogeneous differential equation

\[Ay''+By'+Cy=0\]

When solving for a first-order differential equation, we were able to assume a general solution of \(y_n=e^{st}\). Since our rate is not intrinsically visible in the equation, we will need to solve for it algebraically. We can do this by substituting our general solution and its derivatives into the differential equation.

\[y'=re^{rt}\]

\[y''=r^2e^{rt}\]

\[Ar^2e^{rt}+Bre^{rt}+Ce^{rt}=0\]

With this substitution, we are now able to divide both sides by \(e^{rt}\). This gives us the characteristic polynomial for the differential equation. We can then solve for \(\textit{r}\) using the quadratic formula.

\[Ar^2+Br+C=0\]

\[r=\frac{-B}{2A}\pm\frac{\sqrt{B^2-4AC}}{2A}\]

Second-order systems will always have two null solutions that satisfy the equation. To have a complete solution, both independent \(\textit{r}\) values must be included. The complete null solution will then come out to be

\[y_n=C_1e^{r_1t}+C_2e^{r_2t}\]

Solving for a second-order system also requires two initial conditions, \(y(0)\) and \(y'(0)\), to solve for the arbitrary constants on the differential equations. The derived rates may also turn out to be complex roots. If \(\textit{r}\) is complex then the solution will have an oscillating or sinusoidal behavior.

Undamped Harmonic Oscillator

An undamped harmonic oscillator \((\gamma=0)\) can be described using the differential equation

\[my''+ky=0\]

If we apply the characteristic polynomial to this solution, the rates come out to be purely imaginary.

\[mr^2+k\]

\[r=\pm\sqrt{\frac{-k}{m}}\]

In this case, our general solution can be adapted to a sinusoidal function using Euler's formula. Euler's formula relates complex exponentials to sine and cosine.

\[Re^{i\omega t}=R\cos(\omega t)+iR\sin(\omega t)\]

This proof will be discussed and applied in later sections. It can be used to simplify our general solution to the differential equation from a complex exponential function to real sine and cosine terms.

\[\omega=\sqrt{\frac{k}{m}}\]

\[y_n=C_1e^{i\omega t}+C_2e^{-i\omega t}=A\cos(\omega t)+B\sin(\omega t)\]

This term can be further simplified into a single cosine function. This is because sine functions can be converted to phase-shifted cosine functions. The substitution for these terms is given as the following

\[A\cos(\omega t) + B\sin(\omega t)= G\cos(\omega t + \phi)\]

\[G = \sqrt{A^2+B^2}\]

\[\phi=\tan^{-1}(\frac{B}{A})\]

These substitutions come from trigonometric substitution proofs. \(\textit{G}\) is the total amplitude of the function. \(\phi\) is the phase shift from the combination of the sine and cosine terms.

General Solution Behavior

The roots of the characteristic polynomial define the behavior of homogeneous solutions. The four potential behaviors of a harmonic oscillator are the following.

  • No Damping \(\gamma = 0\)
  • Under-damped \(\gamma^2 < 4mk\)
  • Over-damped \(\gamma^2 > 4mk\)
  • Critically Damped \(\gamma^2 = 4mk\)

Graphed Damped Harmonic Oscillators Functions
Graph of four general harmonic oscillating solutions

These behaviors apply only to positive \(\gamma\) values. Using the characteristic polynomial to solve for homogeneous exponential rates will always have decaying behavior over time. Using the general solution to second-order homogeneous differential equations will yield adequate results but can be translated into different forms using Euler's Formula. 

Under-damped System

Substituting the rates in under-damped systems into a general solution yields the following:

\[\sqrt{\gamma^2 - 4mk} = i\omega\]

\[y(t) = C_1e^{\frac{-\gamma}{2m} + i\omega t} + C_2e^{\frac{-\gamma}{2m} - i\omega t}\]

Factoring out a \(e^{\frac{-\gamma}{2m}t}\) allows the function to simplify even further. By inspection, we can see the terms in the parentheses match our undamped system. Using Euler's Formula, we can convert the imaginary exponential functions into sine and cosine.

\[y(t) = e^{\frac{-\gamma}{2m}t}(C_1e^{i\omega t} + C_2e^{-i\omega t})\]

\[y(t) = e^{\frac{-\gamma}{2m}t}(A\cos(\omega t) + B\sin(\omega t))\]

This form of the solution describes some oscillating function that has a decaying amplitude over time. This decay comes from the real term in the rate function. The oscillating behavior comes from the imaginary terms.

To visualize this behavior, imagine a car's suspension after it goes over a bump. The suspension of an under-damped system would cause the car to oscillate for a time until it reached an equilibrium point.

Over-damped System

In an over-damped harmonic oscillator, the system loses all oscillating behavior. The roots are purely real.

In this situation, the roots and solution to the homogeneous differential equation become

\[s = \frac{-\gamma}{2m} \pm \frac{\sqrt{\gamma^2 - 4mk}}{2m}\]

\[y_n = C_1 e^{s_1 t} + C_2 e^{s_2 t}\]

These roots in a harmonic oscillator will always decay. This only applies with all positive values of \(\textit{m}\), \(\textit{ \gamma }\), and \(\textit{k}\).

Non-Homogeneous Second-Order Equations

The general form of a non-homogeneous second-order differential equation is written as \[Ax''+Bx'+Cx=f(t)\]Where\(f(t)\) is some forcing functions. Common forcing functions include:

  •  \(f(t)=ae^{rt}\)
  •  \(f(t)=a_1\cos(\omega t)+a_2\sin(\omega t)\)
  •  \(f(t)=at^2+bt+c\)
  •  \(f(t)=at^ne^{rt}\)

As in first-order differential equations, a complete solution can be constructed from a combination of particular and homogeneous solutions. We solve for the particular solution by one of two means.

  1. Method of Undetermined Coefficients
  2. Convolution (Green's Function) with an impulse response

Method of Undetermined Coefficients

The method of undetermined coefficients is straightforward in solving particular solutions. To use this method, an assumed particular solution is made from the forcing function. An arbitrary coefficient, known as a transfer function is applied to the solution. This particular solution and its first and second derivatives are substituted into the original equation. Finally, the unknown coefficients are derived from the differential equation.

Method of Undetermined Coefficients: Exponentials

\[f(t) = a e^{rt}\]

Assuming \(r \ne s_1, s_2\), we are able to solve for the particular solution to this forcing function using the following undetermined form:

\[y_p = \mu e^{rt}\]

Where \(\mu\) is the transfer function, and \(r\) is the same rate in the forcing function. We then take derivatives and substitute into the differential equation:

\[y' = \mu r e^{rt}\]

\[y'' = \mu r^2 e^{rt}\]

\[A\mu r^2 e^{rt} + B\mu r e^{rt} + C\mu e^{rt} = a e^{rt}\]

Divide both sides by \(e^{rt}\) and solve for \(\mu\):

\[A\mu r^2 + B\mu r + C\mu = a\]

\[\mu(Ar^2 + Br + C) = a\]

\[\mu = \frac{a}{Ar^2 + Br + C}\]

When this \(\mu\) is applied to the assumed solution, we find the actual particular solution to the differential equation:

\[y_p = \frac{a}{Ar^2 + Br + C} e^{rt}\]

Method of Undetermined Coefficients: Sinusoids

\[f(t) = \cos(\omega t)\]

When the forcing function is any sinusoid, the general particular solution takes the following form:

\[y_p = M\cos(\omega t) + N\sin(\omega t)\]

This is similar to first-order differential equations with sinusoidal forcing functions. The first and second derivatives of the particular solution are as follows:

\[y' = -M\omega\sin(\omega t) + N\omega\cos(\omega t)\]

\[y'' = -M\omega^2\cos(\omega t) - N\omega^2\sin(\omega t)\]

Substituting the derivatives into the differential equation:

\[ A(-M\omega^2\cos(\omega t) - N\omega^2\sin(\omega t)) + B(-M\omega\sin(\omega t) + N\omega\cos(\omega t)) + C(M\cos(\omega t) + N\sin(\omega t)) = \cos(\omega t) \]

The constants \(M\) and \(N\) can be solved by grouping sine and cosine terms to create a linear system of equations:

\[1 = M(C - A\omega^2) + N B\omega\]

\[0 = N(C - A\omega^2) - M B\omega\]

Solving for \(M\) and \(N\) gives the coefficients of the particular solution:

\[M = \frac{C - A\omega^2}{B^2\omega^2 + (C - A\omega^2)^2}\]

\[N = \frac{B\omega}{B^2\omega^2 + (C - A\omega^2)^2}\]

\[ y_p = \frac{C - A\omega^2}{B^2\omega^2 + (C - A\omega^2)^2}\cos(\omega t) + \frac{B\omega}{B^2\omega^2 + (C - A\omega^2)^2}\sin(\omega t) \]

As this general solution is complex, the process of undetermined coefficients is used rather than always using the general form. This particular solution can also be found using complex exponentials by substituting an \(e^{i\omega t}\) in for cosine or sine and taking the real or imaginary part of the solution:

For cosine:

\[f(t) = \cos(\omega t) = \text{Re}[e^{i\omega t}]\]

\[y_p = \text{Re}[(M - iN)e^{i\omega t}]\]

For sine:

\[f(t) = \sin(\omega t) = \text{Im}[e^{i\omega t}]\]

\[y_p = \text{Im}[(M - iN)e^{i\omega t}]\]

Method of Undetermined Coefficients: Polynomials

\[f(t) = t^2\]

When the forcing term is a polynomial function, the particular solution takes the form of a polynomial to the same degree as the forcing function. From the general form, we can take the first and second derivatives and solve the created system of equations.

\[y_p = at^2 + bt + c\]

\[y' = 2at + b\]

\[y'' = 2a\]

As the function has similar coefficients to the general differential equation, we will use the general harmonic differential equation. Substituting the particular solution into the differential equation yields:

\[m(2a) + \gamma(2at + b) + k(at^2 + bt + c) = t^2\]

By grouping terms of the same power, we create a system of 3 equations:

\[ka = 1\]

\[2\gamma a + kb = 0\]

\[2ma + \gamma b + kc = 0\]

This system can be solved sequentially to give the value of the coefficients \(a\), \(b\), and \(c\):

\[a = \frac{1}{k}\]

\[b = \frac{-2\gamma}{k^2}\]

\[c = \frac{2\gamma^2}{k^3} - \frac{2m}{k^3}\]

These coefficients can be substituted back into the original particular solution.

Method of Undetermined Coefficients: Polynomials and Exponentials

\[f(t) = t^2e^{rt}\]

To solve for a forcing function with both exponential and polynomial terms, we combine the general solutions of each. The general particular solution can be written as:

\[y_p = (at^2 + bt + c)e^{rt}\]

When taking the first and second derivatives for this solution, we need to make sure not to neglect the product rule. Our derivatives become:

\[y' = r(at^2 + bt + c)e^{rt} + (2at + b)e^{rt}\]

\[y'' = r^2(at^2 + bt + c)e^{rt} + r(2at + b)e^{rt} + 2ae^{rt}\]

Substitute into the general solution:

\[m(r^2(at^2 + bt + c)e^{rt} + r(2at + b)e^{rt} + 2ae^{rt}) + \gamma(r(at^2 + bt + c)e^{rt} + (2at + b)e^{rt}) + k((at^2 + bt + c)e^{rt}) = t^2e^{rt}\]

Divide out \(e^{rt}\) and group terms based on the degree of the polynomial:

\[mr^2a + \gamma ra + ka = 1\]

\[m(r^2b + 2ra) + \gamma(rb + 2a) + kb = 0\]

\[m(r^2c + rb + 2a) + \gamma(rc + b) + kc = 0\]

Solving this system of equations to find \(a\), \(b\), and \(c\) will yield the correct coefficients to the particular solution that solves the system of equations.

Green's Function

A common equation to solve any particular solution to a differential equation comes from Green's function. This is also known as the convolution function. In the case of differential equations, the function that is convoluted with the forcing function \(f(t)\) is an impulse response function \(g(t)\).

\[y_p = \int \limits_0^t g(t - s)f(s)ds\]

In first-order differential equations, the impulse response was always some form of the null solution: \(g(t) = e^{at}\). In the case of second-order differential equations, the impulse response relies on both null solutions. This impulse response can be written generally as:

\[g(t) = \frac{e^{s_1t} - e^{s_2t}}{s_1 - s_2}\]

This impulse response can be rewritten in equivalent forms of sinusoidal functions when imaginary numbers are involved. A case to be aware of is when \(s_1 = s_2\), as the general equation falls apart. This double-root null solution uses an equivalent impulse response that can be found using L'Hôpital's Rule in the following limit expression:

\[g(t) = \lim_{s_1 \to s_2} \frac{e^{s_1t} - e^{s_2t}}{s_1 - s_2} = te^{s_1t}\]

The solution produced by Green's function may look different than the solution provided by the method of undetermined coefficients. The results of either solution produce the same results though as Green's function normally adds terms to the equation that can be accounted for in the null solution.