How to best obtain the Lyapunov exponent for spatiotemporal chaos

Question

I want to investigate (numerically) whether the spatiotemporal dynamics of a PDE (reaction–diffusion) is chaotic. The Lyapunov exponent focuses on the dependence on initial conditions and is an indicator for chaos. I am not sure how to apply Lyapunov exponents to such a problem. I have three possibilities in mind:

1. Take a particular point in space and calculate the exponent for this point. I‘d stay that if the exponent is positive, a dependence exists. If not, however, this may only be true for that point.
2. Consider the spatial discretization as say 200 coupled ODEs and take the largest exponent.
3. Consider a spatial average.

What is the best way for spatiotemporal chaos? Why? Other ideas are also more than welcome.

Answer 1

Just so we are all on the same page: I assume that you want to evolve your dynamics for two slightly different initial conditions and see whether the difference of some observable diverges exponentially.

Unless some part of your system is uncoupled from the rest, all reasonable observables should yield the same Lyapunov exponent¹. However, the less sensitive an observable is to the phase local oscillations and similar, the less temporal averaging you need to determine the Lyapunov exponent robustly. Let’s reconsider your options in light of this:

1. take a particular point in space and calculate the exponent for this point. I‘d stay that if the exponent is positive, a dependence exists. If not, however, this may only be true for that point.

   Unless your point is somehow uncoupled from the rest of the system, the latter should not happen. For some weakly coupled point, it may however take some time till you see the effect of your perturbation.

2. consider the spatial discretization as say 200 coupled odes and take the largest exponent

   I presume that by this you mean that you consider the exponential growth of the difference $|x(t)-y(t)|$, where $x$ and $y$ are vectors containing all the 200 dynamical variables of your original and perturbed system, respectively.

   This would probably be the most feasible way since you take into account all information and average out the effect of the local phase of oscillations.

3. consider a spatial average

   I presume that by this you mean that you consider the exponential growth of the difference $\left| \bar{x}(t) - \bar{y}(t)\right|$ in the notation of the previous point, with $\bar{·}$ denoting the spacial average.

   While you should see an exponential divergence (with the same Lyapunov exponent) here, it should be rather small in absolute scales and probably short-lived, since the spacial average for most systems exhibiting spatiotemporal chaos should only exhibit small fluctuations – which vanish in the limit of an infinitely large space. It could even be that the spatial average is truly constant due to some inherent conservation laws of your system. I would therefore consider this highly infeasible.

---

¹ Notation:

• $x$ is the solution of your dynamics (comprising all dynamical variables),
• $y$ is the solution of your (infinitesimally) perturbed dynamics,
• $λ$ is the largest Lyapunov exponent,
• $v$ is the corresponding Lyapunov vector,
• $P$ is a linear projection representing the extraction of your observable,
• $≅$ is an equality that only holds if you consider an average over an infinite time or infinitely many initial conditions.

After discarding transients, you get:

$$ y(t) = x(t) + v(t) e^{λt} \qquad\text{with}\qquad |v(t)| ≙ 1.$$

And from this:

$$ \left| Py(t) - Px(t) \right| = |P v(t)| e^{λt} ≅ ce^{λt},$$

where the last equality holds since the Lyapunov vector $v$ is recurring and thus $|P v(t)| ≅ c$ for some constant $c$. The case where $P v(t) = 0 ~∀t$ would correspond to you observing some part of the dynamics that is uncoupled from the rest.

Answer 2

Lyapunov exponents are tied to a norm that will depend on what function space your dynamics are taking place in. In the full infinite dimensional case, your dynamics look something like $u:\mathbb{R}\to V$, where $V$ is some function space. Typically, this is an $L^p$ or a Sobolev space, each of which have associated norms. For numerical experiments, it is probably easier to discretize the PDE then for your norm, compute a discretization of the infinite dimensional norm, e.g. for the $L^2$ norm, compute $\sqrt{\sum_{i=1}^nw_iu_i^2}$, where $w_i$ are appropriate quadrature weights for your nodes.