Solution (source code)

= Solution

Since $L^T=L$ and $B^T=-B$,
$$
[B,L]^T=(BL-LB)^T=L(-B)-(-B)L=[B,L],
$$
so the commutator is symmetric. If $Lv=\lambda v$ with $\lVert v\rVert=1$, symmetry permits a differentiable orthonormal eigenbasis locally, and
$$
\dot\lambda=v^T\dot Lv
=v^T(BL-LB)v
=\lambda v^TBv-\lambda v^TBv=0.
$$
Thus a <Lax pair> evolution is <Isospectral Lax equation>[isospectral].

For the harmonic oscillator,
$$
H(q,p)=\frac12(p^2+\omega^2q^2),\qquad
\dot q=p,\qquad \dot p=-\omega^2q.
$$
With the stated $L$, take
$$
B=\begin{pmatrix}0&-\omega/2\\ \omega/2&0\end{pmatrix}.
$$
Direct calculation gives
$$
[B,L]=
\begin{pmatrix}-\omega^2q&\omega p\\
\omega p&\omega^2q\end{pmatrix}
=\dot L.
$$
The eigenvalues of $L$ are
$$
\lambda_\pm=\pm\sqrt{p^2+\omega^2q^2}.
$$
Their constancy proves conservation of the oscillator energy $H=\lambda_+^2/2$.