New derivation of Lyapunov exponents in doc
This commit is contained in:
parent
acf1c27b1d
commit
328f3bd99c
@ -301,7 +301,7 @@ in which $\mathcal F^*$ is defined like $\mathcal F$ but with the opposite phase
|
||||
\sum_{p,q\in\mathcal K}
|
||||
\frac1{N_1N_2}
|
||||
\sum_{n\in\mathcal N}e^{-\frac{2i\pi}{N_1}n_1(p_1+q_1-k_1)-\frac{2i\pi}{N_2}n_2(p_2+q_2-k_2)}
|
||||
(q\cdot p^\perp)\frac{|q|}{|p|}\hat u_q\hat u_p
|
||||
(q\cdot p^\perp)\frac{|q|}{|p|}\hat u_p\hat u_q
|
||||
\end{equation}
|
||||
provided
|
||||
\begin{equation}
|
||||
@ -440,52 +440,98 @@ The enstrophy is defined as
|
||||
|
||||
\subsection{Lyapunov exponents}
|
||||
\indent
|
||||
We will consider several methods of computing the Lyapunov exponents.
|
||||
The Lyapunov are defined from the {\it tangent flow} of the dynamics.
|
||||
Consider an equation of the form
|
||||
\begin{equation}
|
||||
\dot u=f(t;u)
|
||||
.
|
||||
\end{equation}
|
||||
The tangent flow is given by
|
||||
\begin{equation}
|
||||
\dot\delta=Df(t;u(t))\delta
|
||||
\end{equation}
|
||||
where $Df$ is the Jacobian of $f$ with respect to $u$.
|
||||
The flow of this equation is denoted by
|
||||
\begin{equation}
|
||||
\varphi_{t_0,t_1}(\delta_0)
|
||||
\end{equation}
|
||||
and defined by
|
||||
\begin{equation}
|
||||
\frac d{dt}\varphi_{t_0,t}(\delta_0)=Df(t;\varphi_{t_0,t}(\delta_0))\varphi_{t_0,t}(\delta_0)
|
||||
,\quad
|
||||
\varphi_{t_0,t_0}(\delta_0)=\delta_0
|
||||
.
|
||||
\end{equation}
|
||||
The flow $\varphi_{t_0,t}(\delta_0)$ is linear in $\delta_0$, and so can be represented as a matrix.
|
||||
The Lyapnuov exponents are defined from the $QR$ decomposition of $\varphi_{0,t}$: if
|
||||
\begin{equation}
|
||||
\varphi_{0,t}=Q_tR_t
|
||||
\end{equation}
|
||||
where $Q_t$ is orthogonal, and $R_t$ is triangular with positive diagonal entries $(r_t^{(j)})_j$, then the Lyapunov exponents are
|
||||
\begin{equation}
|
||||
\lim_{t\to\infty}\frac 1t\log|r_t^{(j)}|
|
||||
.
|
||||
\end{equation}
|
||||
\bigskip
|
||||
|
||||
\subsubsection{Method 1: compute the Jacobian at every step}
|
||||
\indent
|
||||
To compute the Lyapunov exponents, we must first compute the Jacobian of $\hat u^{(n)}\mapsto\hat u^{(n+1)}$.
|
||||
This map is always of Runge-Kutta type, that is,
|
||||
One problem to compute the Lyapunov exponents numerically is that the spectrum depends exponentially on time, and so has a tendency to blow up or shrink very rapidly, which leads to large truncation errors.
|
||||
To avoid these, we compute the tangent flow {\it in parts}: we split time into intervals:
|
||||
\begin{equation}
|
||||
\hat u(t_{n+1})=\mathfrak F_{t_n}(\hat u(t_n))
|
||||
.
|
||||
\end{equation}
|
||||
Let $D\mathfrak F_{t_{n}}$ be the Jacobian of this map, in which we split the real and imaginary parts: if
|
||||
\begin{equation}
|
||||
\hat u_k(t_n)=:\rho_k+i\iota_k
|
||||
[0,t)=\bigcup_{i=0}^{N-1}[t_i,t_{i+1})
|
||||
,\quad
|
||||
\mathfrak F_{t_n}(\hat u(t_n))_k=:\phi_k+i\psi_k
|
||||
\end{equation}
|
||||
then
|
||||
\begin{equation}
|
||||
(D\mathfrak F_{t_n})_{k,p}:=\left(\begin{array}{cc}
|
||||
\partial_{\rho_p}\phi_k&\partial_{\iota_p}\phi_k\\
|
||||
\partial_{\rho_p}\psi_k&\partial_{\iota_p}\psi_k
|
||||
\end{array}\right)
|
||||
\end{equation}
|
||||
We compute this Jacobian numerically using a finite difference, by computing
|
||||
\begin{equation}
|
||||
(D\mathfrak F_{t_n})_{k,p}:=\frac1\epsilon
|
||||
\left(\begin{array}{cc}
|
||||
\phi_k(\hat u+\epsilon\delta_p)-\phi_k(\hat u)&\phi_k(\hat u+i\epsilon\delta_p)-\phi_k(\hat u)\\
|
||||
\psi_k(\hat u+\epsilon\delta_p)-\psi_k(\hat u)&\psi_k(\hat u+i\epsilon\delta_p)-\psi_k(\hat u)
|
||||
\end{array}\right)
|
||||
\varphi_{0,t}=\prod_{i=N-1}^0\varphi_{t_{i},t_{i+1}}
|
||||
.
|
||||
\end{equation}
|
||||
At the thresholds between these intervals, we perform a QR decomposition:
|
||||
\begin{equation}
|
||||
Q_0R_0:=\varphi_{0,t_0}
|
||||
,\quad
|
||||
Q_iR_i:=\varphi_{t_{i-1},t_i}Q_{i-1}
|
||||
\end{equation}
|
||||
where $Q_i$ is orthogonal and $R_i$ is upper triangular with $\geqslant 0$ diagonal entries (in addition, we assume these are not $0$).
|
||||
Doing so, we find
|
||||
\begin{equation}
|
||||
\varphi_{0,t}=Q_{N-1}\prod_{i=N-1}^0R_i
|
||||
\end{equation}
|
||||
and since the product of triangular matrices is triangular and the diagonal elements multiply, we find that the Lyapunov exponents are given by
|
||||
\begin{equation}
|
||||
\lim_{t_{N-1}\to\infty}\frac1{t_{N-1}}\sum_{i=N-1}^0\log|r_i^{(j)}|
|
||||
.
|
||||
\end{equation}
|
||||
To do the computation numerically, we drop the limit, and compute the logarithm of the product of the diagonal entries of $R_i$.
|
||||
\bigskip
|
||||
|
||||
\indent
|
||||
In practice, we approximate $\varphi_{t_{i-1},t_i}$ by running a Runge-Kutta algorithm for the tangent flow equation.
|
||||
To obtain the full matrix, we consider every element of the canonical basis as an initial condition $\delta_0$.
|
||||
We then iterate the Runge-Kutta algorithm until the time $t_0$ (chosen in one of two ways, see below), at which point we perform a QR decomposition, save the diagonal entries of $R$, replace the family of initial conditions with the columns of $Q$, and continue the flow from there.
|
||||
The choice of the times $t_i$ can be done either by fixed-length intervals, specified with the option {\tt lyapunov\_reset}, or the QR decomposition can be triggered whenever $\|\delta\|_1$ exceeds a threshold, specified in {\tt lyapunov\_maxu} (after all, the intervals are used to prevent $\delta$ from becoming too large).
|
||||
\bigskip
|
||||
|
||||
\indent
|
||||
To compute the Lyapunov exponents, we thus need the Jacobian of $f$.
|
||||
For the irreversible equation,
|
||||
\begin{equation}
|
||||
f(\hat u)=
|
||||
-\frac{4\pi^2}{L^2}\nu k^2\hat u_k+\hat g_k
|
||||
+\frac{4\pi^2}{L^2|k|}T(\hat u,k)
|
||||
\end{equation}
|
||||
and so
|
||||
\begin{equation}
|
||||
((D f(\hat u))\delta)_k
|
||||
=
|
||||
-\frac{4\pi^2}{L^2}\nu k^2\delta_k
|
||||
+\frac{4\pi^2}{L^2|k|}DT(\hat u,k)\delta
|
||||
\end{equation}
|
||||
where
|
||||
\begin{equation}
|
||||
DT(\hat u,k)\delta
|
||||
=
|
||||
\sum_{\displaystyle\mathop{\scriptstyle p,q\in\mathbb Z^2}_{p+q=k}}
|
||||
\left(\frac{(q\cdot p^\perp)|q|}{|p|}+\frac{(p\cdot q^\perp)|p|}{|q|}\right)\hat u_p\hat \delta_q
|
||||
.
|
||||
\end{equation}
|
||||
The parameter $\epsilon$ can be set using the parameter {\tt D\_epsilon}.
|
||||
%, so
|
||||
%\begin{equation}
|
||||
% D\hat u^{(n+1)}=\mathds 1+\delta\sum_{i=1}^q w_iD\mathfrak F(\hat u^{(n)})
|
||||
% .
|
||||
%\end{equation}
|
||||
%We then compute
|
||||
%\begin{equation}
|
||||
% (D\mathfrak F(\hat u))_{k,\ell}
|
||||
% =
|
||||
% -\frac{4\pi^2}{L^2}\nu k^2\delta_{k,\ell}
|
||||
% +\frac{4\pi^2}{L^2|k|}\partial_{\hat u_\ell}T(\hat u,k)
|
||||
%\end{equation}
|
||||
%and, by\-~(\ref{T}),
|
||||
%\begin{equation}
|
||||
% \partial_{\hat u_\ell}T(\hat u,k)
|
||||
@ -504,30 +550,115 @@ The parameter $\epsilon$ can be set using the parameter {\tt D\_epsilon}.
|
||||
% \right)\hat u_{k-\ell}
|
||||
% .
|
||||
%\end{equation}
|
||||
\bigskip
|
||||
For the reversible equation,
|
||||
\begin{equation}
|
||||
f(\hat u)=
|
||||
-\frac{4\pi^2}{L^2}\alpha(\hat u) k^2\hat u_k
|
||||
+\hat g_k
|
||||
+\frac{4\pi^2}{L^2|k|}T(\hat u,k)
|
||||
\end{equation}
|
||||
so
|
||||
\begin{equation}
|
||||
((D f(\hat u))\delta)_k
|
||||
=
|
||||
-\frac{4\pi^2}{L^2}\alpha(\hat u) k^2\delta_k
|
||||
-\frac{4\pi^2}{L^2}k^2\hat u_k D\alpha(\hat u)\delta
|
||||
+\frac{4\pi^2}{L^2|k|}DT(\hat u,k)\delta
|
||||
\end{equation}
|
||||
where
|
||||
\begin{equation}
|
||||
D\alpha(\hat u)\delta
|
||||
=
|
||||
\frac{\frac{L^2}{4\pi^2}\sum_k k^2\hat \delta_k^*\hat g_k+\sum_k|k|\hat (\delta_k^*T(\hat u,k)+\hat u_k^*DT(\hat u,k)\delta)}{\sum_kk^4|\hat u_k|^2}
|
||||
-\alpha(\hat u)\frac{2\sum_kk^4\delta_k^*\hat u_k}{\sum_kk^4|\hat u_k|^2}
|
||||
.
|
||||
\end{equation}
|
||||
|
||||
\indent
|
||||
The Lyapunov exponents are then defined as
|
||||
\begin{equation}
|
||||
\frac1{t_n}\log\mathrm{spec}(\pi_{t_n})
|
||||
,\quad
|
||||
\pi_{t_n}:=\prod_{i=1}^nD\hat u(t_i)
|
||||
.
|
||||
\end{equation}
|
||||
However, the product of $D\hat u$ may become quite large or quite small (if the exponents are not all 1).
|
||||
To avoid this, we periodically rescale the product.
|
||||
We set $\mathfrak L_r>0$ (set by adjusting the {\tt lyanpunov\_reset} parameter), and, when $t_n$ crosses a multiple of $\mathfrak L_r$, we rescale the eigenvalues of $\pi_i$ to 1.
|
||||
To do so, we perform a $QR$ decomposition:
|
||||
\begin{equation}
|
||||
\pi_{\alpha\mathfrak L_r}
|
||||
=R^{(\alpha)}Q^{(\alpha)}
|
||||
\end{equation}
|
||||
where $Q^{(\alpha)}$ is orthogonal and $R^{(\alpha)}$ is a diagonal matrix, and we divide by $R^{(\alpha)}$ (thus only keeping $Q^{(\alpha)}$).
|
||||
The Lyapunov exponents at time $\alpha\mathfrak L_r$ are then
|
||||
\begin{equation}
|
||||
\frac1{\alpha\mathfrak L_r}\sum_{\beta=1}^\alpha\log\mathrm{spec}(Q^{(\beta)})
|
||||
.
|
||||
\end{equation}
|
||||
|
||||
%\indent
|
||||
%To compute the Lyapunov exponents, we must first compute the Jacobian of $\hat u^{(n)}\mapsto\hat u^{(n+1)}$.
|
||||
%This map is always of Runge-Kutta type, that is,
|
||||
%\begin{equation}
|
||||
% \hat u(t_{n+1})=\mathfrak F_{t_n}(\hat u(t_n))
|
||||
% .
|
||||
%\end{equation}
|
||||
%Let $D\mathfrak F_{t_{n}}$ be the Jacobian of this map, in which we split the real and imaginary parts: if
|
||||
%\begin{equation}
|
||||
% \hat u_k(t_n)=:\rho_k+i\iota_k
|
||||
% ,\quad
|
||||
% \mathfrak F_{t_n}(\hat u(t_n))_k=:\phi_k+i\psi_k
|
||||
%\end{equation}
|
||||
%then
|
||||
%\begin{equation}
|
||||
% (D\mathfrak F_{t_n})_{k,p}:=\left(\begin{array}{cc}
|
||||
% \partial_{\rho_p}\phi_k&\partial_{\iota_p}\phi_k\\
|
||||
% \partial_{\rho_p}\psi_k&\partial_{\iota_p}\psi_k
|
||||
% \end{array}\right)
|
||||
%\end{equation}
|
||||
%We compute this Jacobian numerically using a finite difference, by computing
|
||||
%\begin{equation}
|
||||
% (D\mathfrak F_{t_n})_{k,p}:=\frac1\epsilon
|
||||
% \left(\begin{array}{cc}
|
||||
% \phi_k(\hat u+\epsilon\delta_p)-\phi_k(\hat u)&\phi_k(\hat u+i\epsilon\delta_p)-\phi_k(\hat u)\\
|
||||
% \psi_k(\hat u+\epsilon\delta_p)-\psi_k(\hat u)&\psi_k(\hat u+i\epsilon\delta_p)-\psi_k(\hat u)
|
||||
% \end{array}\right)
|
||||
% .
|
||||
%\end{equation}
|
||||
%The parameter $\epsilon$ can be set using the parameter {\tt D\_epsilon}.
|
||||
%%, so
|
||||
%%\begin{equation}
|
||||
%% D\hat u^{(n+1)}=\mathds 1+\delta\sum_{i=1}^q w_iD\mathfrak F(\hat u^{(n)})
|
||||
%% .
|
||||
%%\end{equation}
|
||||
%%We then compute
|
||||
%%\begin{equation}
|
||||
%% (D\mathfrak F(\hat u))_{k,\ell}
|
||||
%% =
|
||||
%% -\frac{4\pi^2}{L^2}\nu k^2\delta_{k,\ell}
|
||||
%% +\frac{4\pi^2}{L^2|k|}\partial_{\hat u_\ell}T(\hat u,k)
|
||||
%%\end{equation}
|
||||
%%and, by\-~(\ref{T}),
|
||||
%%\begin{equation}
|
||||
%% \partial_{\hat u_\ell}T(\hat u,k)
|
||||
%% =
|
||||
%% \sum_{\displaystyle\mathop{\scriptstyle q\in\mathbb Z^2}_{\ell+q=k}}
|
||||
%% \left(
|
||||
%% \frac{(q\cdot \ell^\perp)|q|}{|\ell|}
|
||||
%% +
|
||||
%% \frac{(\ell\cdot q^\perp)|\ell|}{|q|}
|
||||
%% \right)\hat u_q
|
||||
%% =
|
||||
%% (k\cdot \ell^\perp)\left(
|
||||
%% \frac{|k-\ell|}{|\ell|}
|
||||
%% -
|
||||
%% \frac{|\ell|}{|k-\ell|}
|
||||
%% \right)\hat u_{k-\ell}
|
||||
%% .
|
||||
%%\end{equation}
|
||||
%\bigskip
|
||||
%
|
||||
%\indent
|
||||
%The Lyapunov exponents are then defined as
|
||||
%\begin{equation}
|
||||
% \frac1{t_n}\log\mathrm{spec}(\pi_{t_n})
|
||||
% ,\quad
|
||||
% \pi_{t_n}:=\prod_{i=1}^nD\hat u(t_i)
|
||||
% .
|
||||
%\end{equation}
|
||||
%However, the product of $D\hat u$ may become quite large or quite small (if the exponents are not all 1).
|
||||
%To avoid this, we periodically rescale the product.
|
||||
%We set $\mathfrak L_r>0$ (set by adjusting the {\tt lyanpunov\_reset} parameter), and, when $t_n$ crosses a multiple of $\mathfrak L_r$, we rescale the eigenvalues of $\pi_i$ to 1.
|
||||
%To do so, we perform a $QR$ decomposition:
|
||||
%\begin{equation}
|
||||
% \pi_{\alpha\mathfrak L_r}
|
||||
% =R^{(\alpha)}Q^{(\alpha)}
|
||||
%\end{equation}
|
||||
%where $Q^{(\alpha)}$ is orthogonal and $R^{(\alpha)}$ is a diagonal matrix, and we divide by $R^{(\alpha)}$ (thus only keeping $Q^{(\alpha)}$).
|
||||
%The Lyapunov exponents at time $\alpha\mathfrak L_r$ are then
|
||||
%\begin{equation}
|
||||
% \frac1{\alpha\mathfrak L_r}\sum_{\beta=1}^\alpha\log\mathrm{spec}(Q^{(\beta)})
|
||||
% .
|
||||
%\end{equation}
|
||||
|
||||
|
||||
|
||||
|
||||
Loading…
x
Reference in New Issue
Block a user