Rewriting 2nd-order ODE as system of two 1st-order ODEs

The vast collection of ODE solvers (e.g., in Odespy) cannot be applied to $$u^{\prime\prime} + \omega^2 u = 0$$ unless we write this higher-order ODE as a system of 1st-order ODEs.

Introduce an auxiliary variable $v=u^{\prime}$: $$\begin{align} u^{\prime} &= v, \tag{1}\\ v^{\prime} &= -\omega^2 u \tag{2} \tp \end{align}$$

Initial conditions: $u(0)=I$ and $v(0)=0$.