4.8: Example - the simple harmonic oscillator
- Page ID
- 219165
\( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)
\( \newcommand{\dsum}{\displaystyle\sum\limits} \)
\( \newcommand{\dint}{\displaystyle\int\limits} \)
\( \newcommand{\dlim}{\displaystyle\lim\limits} \)
\( \newcommand{\id}{\mathrm{id}}\) \( \newcommand{\Span}{\mathrm{span}}\)
( \newcommand{\kernel}{\mathrm{null}\,}\) \( \newcommand{\range}{\mathrm{range}\,}\)
\( \newcommand{\RealPart}{\mathrm{Re}}\) \( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)
\( \newcommand{\Argument}{\mathrm{Arg}}\) \( \newcommand{\norm}[1]{\| #1 \|}\)
\( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)
\( \newcommand{\Span}{\mathrm{span}}\)
\( \newcommand{\id}{\mathrm{id}}\)
\( \newcommand{\Span}{\mathrm{span}}\)
\( \newcommand{\kernel}{\mathrm{null}\,}\)
\( \newcommand{\range}{\mathrm{range}\,}\)
\( \newcommand{\RealPart}{\mathrm{Re}}\)
\( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)
\( \newcommand{\Argument}{\mathrm{Arg}}\)
\( \newcommand{\norm}[1]{\| #1 \|}\)
\( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)
\( \newcommand{\Span}{\mathrm{span}}\) \( \newcommand{\AA}{\unicode[.8,0]{x212B}}\)
\( \newcommand{\vectorA}[1]{\vec{#1}} % arrow\)
\( \newcommand{\vectorAt}[1]{\vec{\text{#1}}} % arrow\)
\( \newcommand{\vectorB}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\( \newcommand{\vectorC}[1]{\textbf{#1}} \)
\( \newcommand{\vectorD}[1]{\overrightarrow{#1}} \)
\( \newcommand{\vectorDt}[1]{\overrightarrow{\text{#1}}} \)
\( \newcommand{\vectE}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash{\mathbf {#1}}}} \)
\( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)
\(\newcommand{\longvect}{\overrightarrow}\)
\( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)
\(\newcommand{\ket}[1]{\left| #1 \right>}\)
\(\newcommand{\bra}[1]{\left< #1 \right|}\)
\(\newcommand{\braket}[2]{\left< #1 \vphantom{#2} \right| \left. #2 \vphantom{#1} \right>}\)
\(\newcommand{\braopket}[3]{\left< #1 \vphantom{#2}\vphantom{#3} \right| #2 \vphantom{#1}\vphantom{#3} \left| #3 \vphantom{#1}\vphantom{#2} \right>}\)
\(\newcommand{\qmvec}[1]{\mathbf{\vec{#1}}}\)
\(\newcommand{\op}[1]{\hat{\mathbf{#1}}}\)
\(\newcommand{\expect}[1]{\langle #1 \rangle}\)
\(\newcommand{\dfn}[1]{\emph{\textbf{#1}}}\)
Let’s explore using RK4 using our old friend the SHO. Once again broken down into two first order ODEs, the SHO is
\[\begin{aligned} \nonumber \frac{d u}{d t} &= -\frac{k}{m} v \\ \nonumber \frac{d v}{d t} &= u\end{aligned}\]
To use RK4 we write this in the general form
\[\begin{aligned} \frac{d u}{d t} &= f(t, u, v) \\ \frac{d v}{d t} &= g(t, u, v) \end{aligned} \label{eq:4.2}\]
with \(f(t, u, v) = -(k/m) v\) and \(g(t, u, v) = u\). To implement this method on the computer we treat all variables as vectors. That is, we gather the LHS expression as
\[\nonumber y = \begin{pmatrix} du/dt \\ dv/dt \end{pmatrix}\]
and the RHS as
\[\nonumber F = \begin{pmatrix} f(t, u, v) \\ g(t, u, v) \end{pmatrix}\]
so the vectorized system to integrate is
\[\nonumber \frac{dy}{dt} = F(t, y)\]
Then each of the intermediate steps \(k_i\) is a \(2 \times 1\) column vector. Writing the problem this way makes it easy to implement RK4 using vectorized expressions. Results obtained for integrating the SHO using RK4 are shown in 4.12 alongside the analytic result \(\cos(\omega t)\). The results look great – both \(h = 0.1\) and \(h=0.01\) produce cosine waves which look picture-perfect – at least over the scale and time period simulated. However, the results are not perfect – 4.13 shows the difference (error) between the computed and the analytic results. Although the difference is small, it is nonzero. Worse, the error grows exponentially with increasing time. This suggests the RK4 method – like forward and backward Euler – is unstable for the simple harmonic oscillator, albeit much better than the other methods.

