Can you use Oseledet's theorem to numerically approximate the Lyapunov spectra? Let's say you have a set of first order differential equations with known Jacobian $J$. Let $x_0, x_1, ..., x_n$ be sampled points on the trajectory near the attractor.
Let $T_n = J(x_{n-1})J(x_{n-2})...J(x_0)$.
Oseledet says that
$$\frac{1}{n} \log(\mbox{the singular values of } T_n)$$
converge to the Lyapunov spectra as $n \to \infty$.
However, can this be used in practice to numerically approximate the Lyapunov spectra (say just to the 2nd decimal place)? Is it known how big $n$ must be to do so?
 A: Once you have a method to estimate the top Lyapunov exponent, you can use the action on wedge products to estimate the other exponents (This should be clear if you examine the proof of Oseledet's theorem.) Estimating the top exponent is difficult, in general. Doing it via the definition is very slow as you observed. More efficient approaches are usually based on Furstenberg's foundational paper [1] that preceded Oseledets' work. The key is to estimate the stationary measure on projective space constructed in [1]. See also [2].  Sometimes this measure is discrete which allows for efficient computation [3]. A sophisticated complex analysis approach is in [4].
[1] Furstenberg, Harry. "Noncommuting random products." Transactions of the American Mathematical Society 108, no. 3 (1963): 377-428.
[2] Gol'dsheid, I. Ya, and Grigorii Aleksandrovich Margulis. "Lyapunov indices of a product of random matrices." Russian mathematical surveys 44, no. 5 (1989): 11-71.
[3] Kenyon, Richard, and Yuval Peres. "Intersecting random translates of invariant Cantor sets." Inventiones mathematicae 104, no. 1 (1991): 601-629.
[4] Pollicott, Mark. "Maximal Lyapunov exponents for random matrix products." Inventiones mathematicae 181, no. 1 (2010): 209-226.
