Variational Equations
Definition 1 2
Let the space $X = \mathbb{R}^{n}$ and the function $f : X \to X$ be given such that the following vector field is represented by a differential equation: $$ \dot{x} = f(x) $$ For the Jacobian matrix $J$ of $f$, the following is called the variational equation: $$ \dot{Y} = J Y $$ Here, the initial condition of the matrix function $Y = Y(t) \in \mathbb{R}^{n \times n}$ is set to be the identity matrix $Y(0) = I$.
Explanation
Since the Jacobian matrix $J$ is a matrix function $J = J(t) = J \left( x(t) \right)$ that continuously changes according to the trajectory $x(t)$ of the original system, the variational equation is not a simple linear differential equation as it might appear.
Geometrically, $Y$ can be thought of as showing how that tangent vector itself acts while turning into $x(t)$ with a slight movement from $x(0)$ of the original system.
Lyapunov Spectrum
When calculating the Lyapunov spectrum in continuous systems, solving the variational equation involves repeatedly calculating $Y(t)$ and observing how the open ball is geometrically transformed, measuring the axis lengths. The corresponding RK4 is a method viewed as solving a linear system $\dot{Y} = J Y$ for a linear transformation $T_{J}$ as expressed by $$ \dot{U} = T_{J} (U) = J U $$ Thus, when a sufficiently small timestep $h$ is given for $U(t)$, the following computation is performed to find $U(t + h)$. $$ \begin{align*} U(t + h) \approx RK4(J, U, h) =& {\frac{ h }{ 6 }} \left( V_{1} + 2 V_{2} + 2 V_{3} + V_{4} \right) \\ V_{1} =& J U \\ V_{2} =& J \left( U + {\frac{ h }{ 2 }} V_{1} \right) \\ V_{3} =& J \left( U + {\frac{ h }{ 2 }} V_{2} \right) \\ V_{4} =& J \left( U + h V_{3} \right) \end{align*} $$ Below is Julia code implementing the discussed method for matrices.
function RK4(J::AbstractMatrix, U::AbstractMatrix, dt=1e-2)
V1 = J*U
V2 = J*(U + (dt/2)*V1)
V3 = J*(U + (dt/2)*V2)
V4 = J*(U + dt*V3)
return U + (dt/6)*(V1 + 2V2 + 2V3 + V4)
end
Yorke. (1996). CHAOS: An Introduction to Dynamical Systems: p382. ↩︎
Karlheinz Geist, Ulrich Parlitz, Werner Lauterborn, Comparison of Different Methods for Computing Lyapunov Exponents, Progress of Theoretical Physics, Volume 83, Issue 5, May 1990, Pages 875–893, https://doi.org/10.1143/PTP.83.875 ↩︎