text:
# A Perturbative Solution for Gravitational Waves in Quadratic Gravity

## 1 Introduction

Gravitational Lagrangians which include quadratic curvature terms, have been proposed at classical level as an extension of Einstein’s theory . More recently, quadratic Lagrangians have been used to yield renormalizable theories of gravity coupled to matter . Higher derivative theories also arise as the low energy limit of string theories .

For the sake of definiteness we will deal with the quadratic (fourth-order) gravity generated by the action

$`S={\displaystyle \frac{1}{16\pi G}}{\displaystyle d^4x\sqrt{g}\left\{R+\alpha R^2+\beta R_{\mu \nu }R^{\mu \nu }+16\pi G_m\right\}}.`$ (1)
Here, $`R`$ is the Ricci scalar, $`R_{\mu \nu }`$ is the Ricci tensor, $`G`$ is the Newton’s gravitational constant, and $`\alpha `$, $`\beta `$ the coupling parameters of the quadratic interactions. We use throughout this work units in which $`\mathrm{}=c=1`$, metric signature $`(+)`$, Riemann tensor $`R_{}^{\alpha }{}_{\beta \mu \nu }{}^{}=\mathrm{\Gamma }_{\beta \mu ;\nu }^\alpha +\mathrm{}`$ and Ricci tensor $`R_{\mu \nu }=g^{\alpha \beta }R_{\alpha \mu \nu \beta }`$.

Let us consider the action $`S`$. For $`\alpha =\beta =0`$, $`S`$ is exactly reduced to the Einstein-Hilbert action. Then, we can expect that for $`\alpha `$ and $`\beta `$ sufficiently small the quadratic gravity can be approximately described by Einstein’s General Relativity second-order equations with small corrections due to quadratic curvature couplings. In this case approximate solutions to quadratic gravity can be constructed by adding successive perturbations to a solution of the Einstein’s gravity. This perturbative approach is develop from the concept of regular reduction of a system of partial differential equations  and is used by Lousto et al.  to obtain perturbative solutions to quadratic gravity for some non-radiative astrophysical scenarios. In the present work instead, we use a perturbative approach to build radiative solutions to quadratic gravity.

In section 2, we give a description of the scheme on how the perturbative field equations are obtained. In section 3, we use the perturbative method to build a solution to the linearized quadratic gravity. In section 4 we discuss the implications of the solution obtained to astrophysics and add some considerations concerning the detectability of the computed effect.

## 2 The general perturbative scheme

We give in this section, a description of the general scheme on how the perturbative equations are obtained. We then show how a solution to quadratic gravity, as a power series in $`\alpha `$ and $`\beta `$ parameters, can be built.

Let us now write the general field equations to the fourth-order gravity. Variations of the action $`S`$, with respect to the metric $`g_{\mu \nu }`$ leads to the following field equations:

$`G_{\mu \nu }+\alpha H_{\mu \nu }+\beta I_{\mu \nu }=8\pi GT_{\mu \nu },`$ (2)
where $`G_{\mu \nu }`$ is the Einstein’s tensor given by

$`G_{\mu \nu }=R_{\mu \nu }{\displaystyle \frac{1}{2}}g_{\mu \nu }R.`$ (3)
The terms in the fourth-order derivatives of the space-time metric are given by

$`H_{\mu \nu }=2R_{,\mu \nu }2g_{\mu \nu }\mathrm{}R+{\displaystyle \frac{1}{2}}g_{\mu \nu }R^22RR_{\mu \nu }`$ (4)
and

$`I_{\mu \nu }=R_{,\mu \nu }{\displaystyle \frac{1}{2}}g_{\mu \nu }\mathrm{}R\mathrm{}R_{\mu \nu }2R_{\mu }^{}{}_{}{}^{\alpha }R_{\alpha \nu }+{\displaystyle \frac{1}{2}}g_{\mu \nu }R_{\alpha \beta }R^{\alpha \beta }.`$ (5)
Let us now write the perturbative equations for the system (2), that is, the equations which give as solutions the contribution to the space-time metric of a given order in $`\alpha `$ and $`\beta `$ parameters.

In the zeroth-order ($`\alpha =\beta =0`$), the perturbative equation is by construction the ordinary Einstein’s gravity equation, which is written as

$`G_{\mu \nu }^{(0)}=8\pi GT_{\mu \nu }^{(0)},`$ (6)
where

$`G_{\mu \nu }^{(0)}=R_{\mu \nu }^{(0)}{\displaystyle \frac{1}{2}}g_{\mu \nu }R^{(0)}.`$ (7)
The solution of these equations is denoted by $`g_{\mu \nu }^{(0)}`$ and is independent of $`\alpha `$ and $`\beta `$ parameters. From the above relations we obtain

$`R^{(0)}=8\pi GT^{(0)}`$ (8)
and

$`R_{\mu \nu }^{(0)}=8\pi G\left(T_{\mu \nu }^{(0)}{\displaystyle \frac{1}{2}}\eta _{\mu \nu }T^{(0)}\right).`$ (9)
At the 1th-order we have the equations

$`G_{\mu \nu }^{(1)}=8\pi GT_{\mu \nu }^{(1)}\alpha H_{\mu \nu }^{(0)}\beta I_{\mu \nu }^{(0)},`$ (10)
where $`H_{\mu \nu }^{(0)}`$ and $`I_{\mu \nu }^{(0)}`$ can be written in terms of the energy-momentum tensor from (8) and (2). The solution of these equations at first order in the coupling parameters $`\alpha `$ and $`\beta `$ is denoted by $`g_{\mu \nu }^{(1)}`$.

By repeating this procedure to higher order we obtain at nth-order the following of perturbative equations:

$`G_{\mu \nu }^{(n)}=8\pi GT_{\mu \nu }^{(n)}\alpha H_{\mu \nu }^{(n1)}\beta I_{\mu \nu }^{(n1)}.`$ (11)
These equations show that at each order of perturbation scheme we have the Einstein’s General Relativity equations with a new effective energy-momentum source. Thus, in this perturbative scheme, the quadratic curvature couplings are left out of the dynamics, which is given by Einstein’s gravity equation. The higher order terms enter only as components of the source in the perturbative field equations.

The approximate solution to the quadratic gravity can be constructed by adding the contributions to the space-time metric order by order, that is,

$`g_{\mu \nu }=g_{\mu \nu }^{(0)}+g_{\mu \nu }^{(1)}+g_{\mu \nu }^{(2)}+\mathrm{}.`$ (12)
## 3 The perturbative scheme in linear approximation

The ordinary theory of gravitational emission and propagation is developed in the radiation zone, where the variations of the gravitational field can be treated linearly. Then, we apply the scheme delineated above to the linearized quadratic gravity. To deal with the linearized theory we expand the metric $`g_{\mu \nu }`$ as follows

$`g_{\mu \nu }=\eta _{\mu \nu }+\kappa h_{\mu \nu },`$ (13)
where $`\kappa 1`$, and $`\eta _{\mu \nu }`$ is the Minkowski space-time metric. Within this approximation, we can neglect the nonlinear curvature terms present in (4) and (5). Thus, the linearized version of the eq. (2) can be written as

$`G_{\mu \nu }^{(L)}=\kappa T_{\mu \nu }2{\displaystyle \frac{\alpha }{\kappa }}H_{\mu \nu }^{(L)}2{\displaystyle \frac{\beta }{\kappa }}I_{\mu \nu }^{(L)},`$ (14)
where $`\kappa =\sqrt{16\pi G}`$, $`G_{\mu \nu }^{(L)}`$ is the linearized Einstein tensor defined by

$`G_{\mu \nu }^{(L)}\mathrm{}h_{\mu \nu }h_{(\mu \alpha \nu )}^{,\alpha }+h_{,\mu \nu }\eta _{\mu \nu }(\mathrm{}hh_{}^{\alpha \beta }{}_{,\alpha \beta }{}^{})`$ (15)
and the terms which involve fourth-order derivatives of the metric reduce after linearization to

$`{\displaystyle \frac{H_{\mu \nu }^{(L)}}{\kappa }}=2\mathrm{}h_{,\mu \nu }2h_{}^{\alpha \beta }{}_{,\alpha \beta \mu \nu }{}^{}2\eta _{\mu \nu }(\mathrm{}^2h\mathrm{}h_{}^{\alpha \beta }{}_{,\alpha \beta }{}^{})`$ (16)
and

$`\begin{array}{cc}\hfill {\displaystyle \frac{I_{\mu \nu }^{(L)}}{\kappa }}& =\mathrm{}h_{,\mu \nu }h_{}^{\alpha \beta }{}_{,\alpha \beta \mu \nu }{}^{}{\displaystyle \frac{1}{2}}\eta _{\mu \nu }\mathrm{}(\mathrm{}hh_{}^{\alpha \beta }{}_{,\alpha \beta }{}^{})\hfill \\ & {\displaystyle \frac{1}{2}}\mathrm{}\left(\mathrm{}h_{\mu \nu }h_{(\mu \alpha \nu )}^{,\alpha }+h_{,\mu \nu }\right).\hfill \end{array}`$ (19)
Let us adopt the Hilbert gauge, namely

$`h_{\mu \nu }^{,\nu }={\displaystyle \frac{1}{2}}h_{,\mu }`$ (20)
and define the trace reverse of $`h_{\mu \nu }`$ by

$`\phi _{\mu \nu }h_{\mu \nu }{\displaystyle \frac{1}{2}}\eta _{\mu \nu }h.`$ (21)
In our approximation we deal with metric independent sources, therefore the energy-momentum tensor is independent of $`\alpha `$ and $`\beta `$ parameters. Thus,

$$\begin{array}{cc}\hfill T_{\mu \nu }^{(0)}=T_{\mu \nu }\mathrm{and}T_{\mu \nu }^{(n)}=0\mathrm{for}n1.& \end{array}$$
(22)
Now, we write the perturbative equations expanding the field variable $`\phi _{\mu \nu }`$, following the scheme delineated in the previous section.

The zeroth-order perturbative equation is the linearized Einstein’s gravity equation

$`\mathrm{}\phi _{\mu \nu }^{(0)}=\kappa T_{\mu \nu }^{(0)},`$ (23)
whose solutions in the radiation zone are given by the usual retarded gravitational potential

$`\phi _{\mu \nu }^{(0)}={\displaystyle \frac{\kappa }{4\pi r}}{\displaystyle d^3𝐱^{}T_{\mu \nu }(tr,𝐱^{})},`$ (24)
where $`r`$ denotes the distance from the source to the observer which is very large as compared to the dimensions of the source.

These are just the ordinary Einstein’s General Relativity solutions from which we know that they can be made transverse, traceless and have only spatial components.

By virtue of (22), and the traceless nature of $`\phi _{kl}^{(0)}`$ the 1-th order perturbative equations are independent of $`\alpha `$ and reduce to

$`\mathrm{}\phi _{\mu \nu }^{(1)}=\beta \kappa \mathrm{}T_{\mu \nu }.`$ (25)
These equations can be solved by applying the usual Green’s function method. Cancelling the surface integrals we obtain

$`\phi _{\mu \nu }^{(1)}={\displaystyle \frac{\beta \kappa }{4\pi r}}{\displaystyle \frac{^2}{t^2}}{\displaystyle d^3𝐱^{}T_{\mu \nu }(tr,𝐱^{})},`$ (26)
which can also be written in the form

$`\phi _{kl}^{(1)}=\beta {\displaystyle \frac{^2}{t^2}}\phi _{kl}^{(0)}.`$ (27)
Since $`\phi _{kl}^{(1)}`$ is proportional to $`\phi _{kl}^{(0)}`$ the former is also transverse and traceless.

The approximate solution up to 1th-order is given by

$`\phi _{kl}^{(1)}=\phi _{kl}^{(0)}+\beta {\displaystyle \frac{^2}{t^2}}\phi _{kl}^{(0)}.`$ (28)
Since $`\phi ^{(0)}=0`$, $`h_{kl}^{(1)}=\phi _{kl}^{(1)}`$. Thus we can compute the scalar amplitude, $`h_s=\sqrt{h_+^2+h_\times ^2}`$, for the first quadratic curvature order solution via equation (28). At 1th-order in $`\beta `$ we have

$$h_{s}^{}{}_{}{}^{(1)}=\mathrm{h}^{(0)}+\beta \left(\frac{h_{xx}^{(0)}\ddot{h}_{xx}^{(0)}+h_{xy}^{(0)}\ddot{h}_{xy}^{(0)}}{h^{(0)}}\right),$$
(29)
where

$$h_{s}^{}{}_{}{}^{(0)}=\sqrt{(h_+^{(0)})^2+(h_\times ^{(0)})^2}$$

is the scalar amplitude given by Einstein’s General Relativity and,

$$\ddot{f}=\frac{^2f}{t^2}.$$

For monochromatic signals of angular frequency $`\omega `$ we get the simple result

$$h_{s}^{}{}_{}{}^{(1)}=\left(1+\frac{\beta \omega ^2}{c^2}\right)h_{s}^{}{}_{}{}^{(0)},$$
(30)
whenever $`(\beta \omega ^2)/c^2<1`$, with $`c`$ being the velocity of light. The above result shows the amplitude of a single monochromatic gravitational wave in quadratic gravity to the first approximation in the higher curvature terms.

### 3.1 An exact solution to linearized quadratic gravity

By adding the contributions of higher-order perturbative terms, the solution at the n-th order can be written as

$$h_{kl}^{(n)}=\left(1+\frac{\beta \omega ^2}{c^2}+\frac{\beta ^2\omega ^4}{c^4}+\frac{\beta ^3\omega ^6}{c^6}+\mathrm{}+\frac{\beta ^n\omega ^{2n}}{c^{2n}}\right)h_{kl}^{(0)}.$$
(31)
The series in $`(\beta \omega ^2)/c^2`$ is a geometric one. This series can be summed to give the exact result $`(1\frac{\beta \omega ^2}{c^2})^1`$, whenever $`(\beta \omega ^2)/c^2<1`$. Therefore, an exact solution for a monochromatic wave in quadratic gravity in the linear weak field approximation is given by

$$h_{kl}=\left(1\frac{\beta \omega ^2}{c^2}\right)^1h_{kl}^{(0)},$$
(32)
where $`h_{kl}^{(0)}`$ is the Einstein’s General Relativity solution.

The scalar amplitude $`h_s`$ for this solution is related to Einstein’s gravity scalar amplitude $`\mathrm{h}^{(0)}`$ by the same factor, i. e.,

$$h_s=\left(1\frac{\beta \omega ^2}{c^2}\right)^1h_{s}^{}{}_{}{}^{(0)}.$$
(33)
We observe that this result is valid only for $`(\beta \omega ^2)/c^2<1`$ or, for waves with frequency $`0<\omega <c/\beta ^{1/2}`$.

The result derived above shows a correction to the amplitude of linear oscillations of the space-time metric due to the presence of a Ricci-squared term in gravitational Lagrangian. This result is remarkable simple and predicts a frequency dependent shift in the wave amplitude. We point out that, in the context of higher-order gravity, this is an interesting result since it enable us to know the effect of quadratic curvature terms in linear oscillations of space-time metric by avoiding the negative energy problems that can arise due to the effect of a massive ghost, associated to the Ricci squared term, on the linearized fourth-order field equations. Since the original fourth order field equations are reduced, by the perturbative scheme, to the linearized Einstein’s gravity equations in presence of an effective source, the energy density and the energy flux of gravitational waves can be computed from the same expressions for this quantities derived in Einstein’s linearized gravity.

## 4 Possible implications to astrophysics and cosmology and perspectives to detectability of the amplitude correction

Let us now discuss how the result derived in this paper can affect the predictions for the gravitational radiation generated by astrophysical and cosmological sources and discuss a possibility to detect the predicted quadratic curvature amplitude correction.

Firstly, we need an order of magnitude estimative for the amplitude shift. To first order in $`\beta `$ the shift is given by the product of Einstein’s linearized wave amplitude by the factor $`\beta \omega ^2/c^2`$. Observational constraints on $`\beta `$ parameter derived from the bending of light by a gravitational field in the fourth-order gravity indicates an upper bound of $`|\beta |10^4\mathrm{cm}^2`$ . The same bound is obtained form a translation of the results of the sub-millimeter tests of the inverse square law of gravitational attraction among massive bodies . Due to the small value of $`\beta `$ upper bound, the effect of the Ricci squared correction to the linearized wave amplitude will appear only for high frequency gravitational waves. Even for inspiraling compact binaries of neutron star or stellar black hole pairs, which experiments an intense growth of the frequency when the system approach the coalescence, the frequency of the gravitational oscillations cannot exceed $`10^4\mathrm{Hz}`$. In this case, the shift in the wave amplitude will be $`10^{16}h^{(0)}`$, where $`h^{(0)}`$ is the Einstein’s linearized wave amplitude. This shift is clearly non observable.

For the amplitude shift to be observable, taking the maximal value acceptable for $`\beta `$, the frequency of gravitational waves must lie in the range of $`10^{10}10^{12}\mathrm{Hz}`$. Compact binaries formed by pairs of mini black holes having masses $`M4\times 10^5M_{}`$ are able to emit gravitational waves with frequencies above $`10^{10}\mathrm{Hz}`$. As pointed out in , a significant fraction of the dark halo of our galaxy can consist of a large number ($`10^{19}`$) of sub-lunar mass primordial black holes with mass $`10^{12}M_{}M10^7M_{}`$. Therefore, it is natural to expect that most of these objects form binaries which can emmit very-high-frequency gravitational radiation before coalescing. Then, we can expect a background of gravitational radiation associated with these objects in the frequency band where the quadratic curvature effects might be relevant.

On the other hand, gravitational waves with frequencies in the GHz band are predicted by several models for stochastic backgrounds of gravitational waves of cosmological origin derived from inflation, string cosmology and cosmic strings scenarios . Then, we can expect that the result derived in the present work can affect the predictions for cosmological backgrounds of gravitational radiation, in the high (GHz) frequency spectral band. An indication of the effect of quadratic correction on the spectrum of a stochastic background, which is assumed to be isotropic, stationary and unpolarized can be easily computed. The energy density of the gravitational waves is given by

$$\rho _{gw}=\frac{c^2}{32\pi G}<\dot{h}_{kl}\dot{h}^{kl}>,$$
(34)
where $`\dot{h}_{kl}`$ is the time derivative of the transverse traceless gauge filed and $`<>`$ means spatial average over several wavelengths. If we insert the quadratic correction to gravitational waves in the equation (34) and compute the ensemble average for a stochastic background  we get the result

$$\mathrm{\Omega }_{gw}(f)\frac{1}{\rho _c}\frac{d\rho _{gw}}{d(\mathrm{log}f)}=\left(1+8\pi ^2\frac{\beta f^2}{c^2}\right)\mathrm{\Omega }_{gw}^{(0)}(f),$$
(35)
where $`f`$ is the wave frequency ($`\omega =2\pi f`$), $`\mathrm{\Omega }_{gw}^{(0)}(f)`$ is the “spectral function” predicted by Einstein’s theory and $`\rho _c`$ is the critical energy density required to close the universe. We note that the $`\mathrm{\Omega }_{gw}`$ will differ from $`\mathrm{\Omega }_{gw}^{(0)}`$ only in the ($`10^{10}10^{12})\mathrm{Hz}`$ frequency range. Such very-high-frequency gravitational waves reflect the behaviour of the very early universe . We note that the expression (35) must be taken with care since the computation of $`\rho _{gw}`$ may be different for a given cosmological models. Thus, a deeper inspection of how the cosmological background is affected by the quadratic curvature corrections must be carried out.

Once very-high-frequency ($`10^{10}10^{12}\mathrm{Hz}`$) gravitational radiation could be detected, this detection could be used to set up new observational constraints on the $`\beta `$ parameter. The gravitational wave detectors such as LIGO, VIRGO and the mass resonant detectors are not adequate to observe such high frequency waves. The effect of quadratic curvature gravity on the linearized fluctuations of space-time metric will only be observable when a new generation of detectors with sufficient sensitivity to measure amplitudes of gravitational waves with frequencies in the range of $`(1010^3)\mathrm{GHz}`$ are available.

A gravitational wave detector which explores the effects on electromagnetic fields of the passage of gravitational waves has been proposed by Cruise . The proposed device can be designed to observe gravitational radiation with frequencies in the MHz and GHz range. Although the sensitivity of the Cruise detector does not reach the required sensitivity to detect gravitational wave signals, improvements in the technology might lead to more sensitive detectors in the near future. This detectors will be able to observe the gravitational radiation generated by cosmological sources at the very early moments in the universe. Therefore, the quadratic curvature corrections to gravitational waves might be detectable by a sufficiently sensitive generation of electromagnetic detectors. These detectors might also be able to observe the background radiation produced by a halo population of binary systems composed of sub-lunar mass primordial black holes, whereas these objects may actually exists.

E. C. de Rey Neto would like to thank the Brazilian agency FAPESP for financial support (grants 00/10374-5). O. D. Aguiar and J.C.N. de Araujo thanks the Brazilian agency CNPq for support (grants 300619/92-8, 380503/02-2, respectively).

## 5 References


format: text

text:
# Index Iteration Theory for Symplectic Paths with Applications to Nonlinear Hamiltonian Systems

## 1. Index function theory for symplectic paths

As usual we define the symplectic group by $`\mathrm{Sp}(2n)=\{M\mathrm{GL}(𝐑^{2n})|M^TJM`$ $`=J\}`$, where $`J=\left(\begin{array}{cc}0& I\\ I& 0\end{array}\right)`$, $`I`$ is the identity matrix on $`𝐑^n`$, and $`M^T`$ denotes the transpose of $`M`$. For $`\omega 𝐔`$, the unit circle in the complex plane $`𝐂`$, we define the $`\omega `$-singular subset in $`\mathrm{Sp}(2n)`$ by $`\mathrm{Sp}(2n)_\omega ^0=\{M\mathrm{Sp}(2n)|\omega ^ndet(\gamma (\tau )\omega I)=0\}`$. Here for any $`M\mathrm{Sp}(2n)_\omega ^0`$, we define the orientation of $`\mathrm{Sp}(2n)_\omega ^0`$ at $`M`$ by the positive direction $`\frac{d}{dt}M\mathrm{exp}(tJ)|_{t=0}`$. Since the fundamental solution of a general linear Hamiltonian system with continuous symmetric periodic coefficient $`2n\times 2n`$ matrix function $`B(t)`$,

$$\dot{x}(t)=JB(t)x(t),t𝐑,$$
(1.1)
is a path in $`\mathrm{Sp}(2n)`$ starting from the identity, for $`\tau >0`$ we define the set of symplectic matrix paths by $`𝒫_\tau (2n)=\{\gamma C([0,\tau ],\mathrm{Sp}(2n))|\gamma (0)=I\}`$. For any two path $`\xi `$ and $`\eta :[0,\tau ]\mathrm{Sp}(2n)`$ with $`\xi (\tau )=\eta (0)`$, as usual we define $`\eta \xi (t)`$ by $`\xi (2t)`$ if $`0t\tau /2`$, and $`\eta (2t\tau )`$ if $`\tau /2t\tau `$. We define a special path $`\zeta :[0,\tau ]\mathrm{Sp}(2n)`$ by

$$\zeta (t)=\mathrm{diag}(2\frac{t}{\tau },\mathrm{},2\frac{t}{\tau },(2\frac{t}{\tau })^1,\mathrm{},(2\frac{t}{\tau })^1),\mathrm{for}\mathrm{\hspace{0.33em}0}t\tau .$$

Definition 1. (cf. ) For any $`\tau >0`$, $`\omega 𝐔`$, and $`\gamma 𝒫_\tau (2n)`$, we define the $`\omega `$-nullity of $`\gamma `$ by

$$\nu _\omega (\gamma )=dim_𝐂\mathrm{ker}_𝐂(\gamma (\tau )\omega I).$$
(1.2)
If $`\gamma `$ is $`\omega `$ non-degenerate, i.e., $`\nu _\omega (\gamma )=0`$, we define the $`\omega `$-index of $`\gamma `$ by the intersection number

$$i_\omega (\gamma )=[\mathrm{Sp}(2n)_\omega ^0:\gamma \zeta ].$$
(1.3)
If $`\gamma `$ is $`\omega `$ degenerate, i.e., $`\nu _\omega (\gamma )>0`$, we let $`(\gamma )`$ be the set of all open neighborhoods of $`\gamma `$ in $`𝒫_\tau (2n)`$, and define

$$i_\omega (\gamma )=\underset{U(\gamma )}{sup}inf\{i_\omega (\beta )|\beta U,\nu _\omega (\beta )=0\}.$$
(1.4)
Then we call $`(i_\omega (\gamma ),\nu _\omega (\gamma ))𝐙\times \{0,1,\mathrm{},2n\}`$ the index function of $`\gamma `$ at $`\omega `$.

The relation of this index $`(i_1(\gamma ),\nu _1(\gamma ))`$ with the Morse index of $`\tau `$-periodic solutions of the problem (1) was proved by C. Conley, E. Zehnder, and the author in , , and  (cf. Theorem 6.1.1 of ).

## 2. Iteration theory of the index for symplectic paths

Given a path $`\gamma 𝒫_\tau (2n)`$, its iteration is defined inductively by $`\gamma (t+\tau )=\gamma (t)\gamma (\tau )`$ for $`t0`$, i.e.,

$$\gamma ^m(t)=\gamma (tj\tau )\gamma (\tau )^j,j\tau t(j+1)\tau ,j=0,1,\mathrm{},m1,$$
(2.1)
for any $`m`$ in the natural integer set $`𝐍`$. For our applications of this index theory to nonlinear Hamiltonian systems, we are facing two types of problems:

$`1`$ knowing the end point $`\gamma (\tau )`$ of a path $`\gamma 𝒫_\tau (2n)`$, the initial index $`(i_1(\gamma ),\nu _1(\gamma ))`$, and the iteration time $`m`$, want to find the index $`i_1(\gamma ^m)`$ of the $`m`$-th iterated path $`\gamma ^m`$;

$`2`$ knowing the end point $`\gamma (\tau )`$ of a path $`\gamma 𝒫_\tau (2n)`$, the initial index $`(i_1(\gamma ),\nu _1(\gamma ))`$, and the index $`(i_1(\gamma ^m),\nu _1(\gamma ^m))`$ of the $`m`$-th iterated path $`\gamma ^m`$, want to find the iteration time $`m`$.

To solve these problems, we first generalize Bott’s formula of the iterated Morse index for closed geodesics to the index theory for general symplectic paths:

Theorem 2 (cf. ). For any $`\tau >0`$, $`\gamma 𝒫_\tau (2n)`$, $`z𝐔`$, and $`m𝐍`$, there hold:

$$i_z(\gamma ^m)=\underset{\omega ^m=z}{}i_\omega (\gamma ),\nu _z(\gamma ^m)=\underset{\omega ^m=z}{}\nu _\omega (\gamma ).$$
(2.2)
By (2.2) it is easy to see that the mean index $`\widehat{i}(\gamma )=lim_{m+\mathrm{}}i_1(\gamma ^m)/m`$ for any $`\gamma 𝒫_\tau (2n)`$ is always a finite real number.

To further solve the problems $`1`$ and $`2`$, we need to go beyond the Bott-type formula (2.2). For a given path $`\gamma `$ we consider to deform it to a new path $`\eta `$ in $`𝒫_\tau (2n)`$ so that

$$i_1(\gamma ^m)=i_1(\eta ^m),\nu _1(\gamma ^m)=\nu _1(\eta ^m),m𝐍,$$
(2.3)
and that $`(i_1(\eta ^m),\nu _1(\eta ^m))`$ is easy enough to compute. This leads to finding homotopies $`\delta :[0,1]\times [0,\tau ]\mathrm{Sp}(2n)`$ starting from $`\gamma `$ in $`𝒫_\tau (2n)`$ and keeping the end points of the homotopy always stay in a certain suitably chosen maximal subset of $`\mathrm{Sp}(2n)`$ so that (2.3) always holds. By (2.2), this set is defined to be the path connected component $`\mathrm{\Omega }^0(M)`$ containing $`M=\gamma (\tau )`$ of the set

$`\mathrm{\Omega }(M)=\{N\mathrm{Sp}(2n)`$ $`|`$ $`\sigma (N)𝐔=\sigma (M)𝐔,\mathrm{and}`$ (2.4)
$`\nu _\lambda (N)=\nu _\lambda (M)\lambda \sigma (M)𝐔\}.`$
Here we call $`\mathrm{\Omega }^0(M)`$ the homotopy component of $`M`$ in $`\mathrm{Sp}(2n)`$.

Using normal forms of symplectic matrices (cf. , ), we then decompose $`\gamma (\tau )`$ within $`\mathrm{\Omega }^0(\gamma (\tau ))`$ into product of 10 special $`2\times 2`$ and $`4\times 4`$ symplectic normal form matrices, which we call basic normal forms. Correspondingly by the homotopy invariance and symplectic additivity of the index theory, the computations in (2.3) are reduced to iterations of those paths in $`\mathrm{Sp}(2)`$ or $`\mathrm{Sp}(4)`$ whose end points are one of the 10 basic normal form matrices. The study of the index for iterations of any symplectic paths is carried out for paths in $`\mathrm{Sp}(2)`$ via the $`𝐑^3`$-cylindrical coordinate representation of $`\mathrm{Sp}(2)`$, then for hyperbolic and elliptic paths in $`\mathrm{Sp}(2n)`$. This yields the precise iteration formula obtained in  of the index theory for any symplectic path $`\gamma 𝒫_\tau (2n)`$ in terms of the basic norm form decomposition of $`\gamma (\tau )`$, $`(i(\gamma ,1),\nu (\gamma ,1))`$, and the iteration time $`m`$.

For any $`M\mathrm{Sp}(2n)`$, its splitting numbers at an $`\omega 𝐔`$ is defined in  by

$$S_M^\pm (\omega )=\underset{ϵ0^+}{lim}i_{\omega \mathrm{exp}(\pm \sqrt{1}ϵ)}(\gamma )i_\omega (\gamma ),$$
(2.5)
via any $`\gamma 𝒫_\tau (2n)`$ satisfying $`\gamma (\tau )=M`$. Then it is proved that the splitting numbers of $`M`$ at $`\omega `$ can be characterized algebraically.

Motivated by the precise iteration formulae of , the following second index iteration formula of any symplectic path is established by C. Zhu and the author. Here we denote by $`(i(\gamma ,m),\nu (\gamma ,m))`$ $`=`$ $`(i_1(\gamma ^m),\nu _1(\gamma ^m))`$.

Theorem 3 (cf. ). For any $`\tau >0`$, $`\gamma 𝒫_\tau (2n)`$, and $`m𝐍`$, there holds:

$`i(\gamma ,m)`$ $`=`$ $`m(i(\gamma ,1)+S_M^+(1)C(M))`$ (2.6)
$`+2{\displaystyle \underset{\theta (0,2\pi )}{}}E({\displaystyle \frac{m\theta }{2\pi }})S_M^{}(e^{\sqrt{1}\theta })(S_M^+(1)+C(M)),`$
where $`M=\gamma (\tau )`$, $`C(M)=_{0<\theta <2\pi }S_M^{}(e^{\sqrt{1}\theta })`$, and $`E(a)=\mathrm{min}\{k𝐙|ka\}`$ for any $`a𝐑`$.

In order to solve problems on nonlinear Hamiltonian systems, various index iteration inequalities for any path $`\gamma 𝒫_\tau (2n)`$ and $`m𝐍`$ are proved by D. Dong, C. Liu, C. Zhu and the author in , , , and .

Theorem 4. For any $`\gamma 𝒫_\tau (2n)`$ and $`m𝐍`$, the following iteration inequalities always hold.

Estimate via mean index (cf. , ):

$$m\widehat{i}(\gamma )ni(\gamma ,m)m\widehat{i}(\gamma )+n\nu (\gamma ,m).$$
(2.7)
Estimate via initial index (cf. ):

$`m(i(\gamma ,1)`$ $`+`$ $`\nu (\gamma ,1)n)+n\nu (\gamma ,1)i(\gamma ,m)`$ (2.8)
$``$ $`m(i(\gamma ,1)+n)n(\nu (\gamma ,m)\nu (\gamma ,1)).`$
Successive index estimate (cf. ):

$`\nu (\gamma ,m){\displaystyle \frac{e(\gamma (\tau ))}{2}}`$ $``$ $`i(\gamma ,m+1)i(\gamma ,m)i(\gamma ,1)`$ (2.9)
$``$ $`\nu (\gamma ,1)\nu (\gamma ,m+1)+{\displaystyle \frac{e(\gamma (\tau ))}{2}}.`$
Here we define $`e(M)`$ to be the total multiplicity of eigenvalues of $`M`$ on $`𝐔`$ and call it the elliptic height of $`M`$.

A consequence of the iteration inequality (2.8) together with the necessary and sufficient conditions for any equality in (2.8) to hold for some $`m`$ yields a new proof of the following theorem of D. Dong and the author on controlling the iteration time $`m`$ via indices:

Theorem 5 (cf. ). For any $`\gamma 𝒫_\tau (2n)`$ and $`m𝐍`$, suppose $`i(\gamma ,m)n+1`$, $`i(\gamma ,1)n`$, and $`\nu (\gamma ,1)1`$. Then $`m=1`$.

Note also that the inequality (2.9) yields a way to estimate the ellipticity of solution orbits of Hamiltonian systems obtained by variational methods via their iterated indices.

In order to study the properties of solution orbits of the system (0.1) on a given energy hypersurface, when the number of orbits is finite, we need to study common properties of any given finite family of symplectic paths $`\gamma _j𝒫_{\tau _j}(2n)`$ with $`1jq`$. This leads to the following common index jump theorem of C. Zhu and the author proved in . For any $`\gamma 𝒫_\tau (2n)`$, its $`m`$-th index jump $`𝒢_m(\gamma )`$ is defined to be the open interval $`𝒢_m(\gamma )=(i(\gamma ,m)+\nu (\gamma ,m)1,i(\gamma ,m+2))`$.

Theorem 6 (cf. ). Let $`\gamma _j𝒫_{\tau _j}(2n)`$ with $`1jq`$ satisfying

$$\widehat{i}(\gamma _j)>0,i(\gamma _j,1)n,1jq.$$
(2.10)
Then there exist infinitely many positive integer tuples $`(N,m_1,\mathrm{},m_q)𝐍^{q+1}`$ such that

$$\mathrm{}[2N\kappa _1,2N+\kappa _2]\underset{j=1}{\overset{q}{}}𝒢_{2m_j1}(\gamma _j),$$
(2.11)
where $`\kappa _1=\mathrm{min}_{1jq}(i(\gamma _j,1)+2S_{\gamma _j(\tau _j)}^+(1)\nu (\gamma _j,1))`$ and $`\kappa _2=\mathrm{min}_{1jq}i(\gamma _j,1)1`$.

In order to prove this theorem, we need to make each index jump to be as big as possible, and to make their largest sizes happen simultaneously to guarantee the existence of a non-empty largest common intersection interval among them. By the term $`E(\frac{m\theta }{2\pi })`$ in the abstract iteration formula (2.6), such a problem is reduced to a dynamical system problem on a torus, and is solved by properties of closed additive subgroups of tori.

## 3. Applications to nonlinear Hamiltonian systems

So far, we have applied our index iteration theory to three important problems on periodic solutions of nonlinear Hamiltonian systems. Let $`T>0`$ and suppose $`x`$ is a non-constant $`T`$-periodic solution of the nonlinear Hamiltonian system (0.1). Suppose the minimal period of $`x`$ is $`\tau =T/k`$ for some $`k𝐍`$. We denote by $`\gamma _x𝒫_\tau (2n)`$ the fundamental solution of the linearized Hamiltonian system (1.1) at $`x`$ with $`B(t)=H^{\prime \prime }(x(t))`$, and the iterated index of $`x`$ by $`(i(x,m),\nu (x,m))=(i(\gamma _x,m),\nu (\gamma _x,m))`$ for all $`m𝐍`$.

### 3.1. Prescribed minimal period solution problem

In  of 1978, P. Rabinowitz posed a conjecture on whether the Hamiltonian system possesses periodic solutions with prescribed minimal period when the Hamiltonian function satisfies his superquadratic conditions. This conjecture is studied by D. Dong and the author as an application of our index iteration theory. Note that for a non-constant $`\tau `$-periodic solution $`x`$ of the autonomous system (0.1), the condition on the nullity in Theorem 5 always holds. Thus Theorem 5 yields:

Theorem 7 (cf. ). For any non-constant $`\tau `$-periodic solution $`x`$ of (0.1), denote its minimal period by $`\tau /m`$ for some $`m𝐍`$. Suppose $`i(x|_{[0,\tau ]},1)n+1`$ and $`ni(x|_{[0,\tau /m]},1)`$. Then $`m=1`$, i.e., $`\tau `$ is the minimal period of $`x`$.

Here the first estimate on the index holds if $`x`$ is obtained by minimax or minimization methods, and the second estimate on the index holds if the Hamiltonian function $`H`$ is convex in a certain weak sense along the orbit of $`x`$. This result reveals the intrinsic relationship between the minimal period of a periodic solution and its indices, and unifies all the results on Rabinowitz’s conjecture under various convexity conditions. Specially, it recovers the famous theorem of I. Ekeland and H. Hofer in 1985 (cf. ) who solved Rabinowitz’s conjecture for convex superquadratic Hamiltonian systems.

### 3.2. Periodic points of the Poincaré map of Lagrangian systems on tori

In 1984, C. Conley stated a conjecture on whether the Poincaré map of any $`1`$-periodic time dependent Hamiltonian system defined on the standard torus $`T^{2n}=𝐑^{2n}/𝐙^{2n}`$ always possesses infinitely many periodic points which are produced by contractible periodic solutions of the corresponding Hamiltonian system on $`T^{2n}`$. A celebrated partial answer to this conjecture was given by D. Salamon and E. Zehnder in 1992 (cf. ) for a large class of symplectic manifolds on which every contractible integer periodic solution of the Hamiltonian system has at least one Floquet multiplier not equal to $`1`$. So far Conley conjecture is still open and seems far from being completely understood.

In , we studied the Lagrangian system version of this conjecture. Consider

$$\frac{d}{dt}L_{\dot{x}}(t,x,\dot{x})L_x(t,x,\dot{x})=0,x𝐑^n,$$
(3.1)
where $`L_{\dot{x}}`$ and $`L_x`$ denote the gradients of $`L`$ with respect to $`\dot{x}`$ and $`x`$ respectively. The main result is the following:

Theorem 8 (cf. ). Suppose the Lagrangian function $`L`$ satisfies

(L1) $`L(t,x,p)=\frac{1}{2}A(t)pp+V(t,x)`$, where $`\frac{1}{2}A(t)pp\lambda |p|^2`$ for all $`(t,p)𝐑\times 𝐑^n`$ and some fixed constant $`\lambda >0`$.

(L2) $`AC^3(𝐑,_s(𝐑^n))`$, $`VC^3(𝐑\times 𝐑^n,𝐑)`$, both $`A`$ and $`V`$ are $`1`$-periodic in all of their variables, where $`_s(𝐑^n)`$ denotes the set of $`n\times n`$ real symmetric matrices.

Then the Poincaré map $`\mathrm{\Psi }`$ of the system (3.1) possesses infinitely many periodic points on $`TT^n`$ produced by contractible integer periodic solutions of the system (3.1).

In the proof of Theorem 8, the above inequality (2.7) plays a crucial role. By this inequality, at very high iteration level, a global homological injection map can be constructed which maps a generator of a certain non-trivial local critical group to a nontrivial homology class $`[\sigma ]`$ in a global homology group, if the number of contractible integer periodic solution towers of the system (3.1) is finite. But on the other hand, by a technique of V. Bangert and W. Klingenberg in , it is shown that this homology class $`[\sigma ]`$ must be trivial globally. This contradiction then yields the conclusion of Theorem 8.

### 3.3. Closed characteristics on convex compact hypersurfaces

Denote the set of all compact strictly convex $`C^2`$-hypersurfaces in $`𝐑^{2n}`$ by $`(2n)`$. For $`\mathrm{\Sigma }(2n)`$ and $`x\mathrm{\Sigma }`$, let $`N_\mathrm{\Sigma }(x)`$ be the outward normal unit vector at $`x`$ of $`\mathrm{\Sigma }`$. We consider the problem of finding $`\tau >0`$ and a curve $`xC^1([0,\tau ],𝐑^{2n})`$ such that

$$\{\begin{array}{ccc}\dot{x}(t)& =& JN_\mathrm{\Sigma }(x(t)),x(t)\mathrm{\Sigma },t𝐑,\hfill \\ x(\tau )& =& x(0).\hfill \end{array}$$
(3.1)
A solution $`(\tau ,x)`$ of the problem (3.1) is called a closed characteristic on $`\mathrm{\Sigma }`$. Two closed characteristics $`(\tau ,x)`$ and $`(\sigma ,y)`$ are geometrically distinct, if $`x(𝐑)y(𝐑)`$. We denote by $`𝒯(\mathrm{\Sigma })`$ the set of all geometrically distinct closed characteristics $`(\tau ,x)`$ on $`\mathrm{\Sigma }`$ with $`\tau `$ being the minimal period of $`x`$. Note that the problem (3.1) can be described in a Hamiltonian system version and solved by variational methods. A closed characteristic $`(\tau ,x)`$ is non-degenerate, if $`1`$ is a Floquet multiplier of $`x`$ of precisely algebraic multiplicity $`2`$, and is elliptic, if all the Floquet multipliers of $`x`$ are on $`𝐔`$. Let $`{}_{}{}^{\mathrm{\#}}A`$ denote the total number of elements in a set $`A`$.

This problem has been studied for more than 100 years since at least A. M. Liapunov in 1892. A long standing conjecture on the multiplicity of closed characteristics is whether

$${}_{}{}^{\mathrm{\#}}\stackrel{~}{𝒥}(\mathrm{\Sigma })n,\mathrm{\Sigma }(2n).$$
(3.2)
The first break through on this problem in the global sense was made by P. Rabinowitz  and A. Weinstein  in 1978. They proved $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })1`$ for all $`\mathrm{\Sigma }(2n)`$. Besides many results under pinching conditions, in 1987–1988, I. Ekeland-L. Lassoued, I. Ekeland-H. Hofer, and A, Szulkin proved $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })2`$ for all $`\mathrm{\Sigma }(2n)`$ and $`n2`$. In 1998, H. Hofer, K. Wysocki, and E. Zehnder proved in : $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })=2`$ or $`+\mathrm{}`$ for every $`\mathrm{\Sigma }(4)`$. In recent years C. Liu, C. Zhu, and the author gave the following answers to the conjecture (3.2):

Theorem 9 (cf. ). There holds

$${}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })[\frac{n}{2}]+1,\mathrm{\Sigma }(2n),$$
(3.3)
where $`[a]=\mathrm{max}\{k𝐙|ka\}`$ for any $`a𝐑`$. Moreover, if all the closed characteristics on $`\mathrm{\Sigma }`$ are non-degenerate, then $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })n`$.

Theorem 10 (cf. ). For any $`\mathrm{\Sigma }(2n)`$, if $`\mathrm{\Sigma }`$ is symmetric with respect to the origin, i.e., $`x\mathrm{\Sigma }`$ implies $`x\mathrm{\Sigma }`$, then $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })n`$.

Very recently, Y. Dong and the author further proved the following result.

Theorem 11 (cf. ). Let $`\mathrm{\Sigma }(2n)`$ be $`P`$-symmetric with respect to the origin, i.e., $`x\mathrm{\Sigma }`$ implies $`Px\mathrm{\Sigma }`$, where $`P=\mathrm{diag}(I_{nk},I_k,I_{nk},I_k)`$ for some fixed integer $`k[0,n1]`$. Let $`\mathrm{\Sigma }(k)=\{(x,y)(𝐑^k)^2|(0,x,0,y)\mathrm{\Sigma }\}`$. Suppose $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma }(k))k`$ or $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma }(k))=+\mathrm{}`$ holds. Then $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })n2k`$.

Proof of Theorem 11 depends on a new index iteration theory for symplectic paths iterated by the formula $`\gamma (t+\tau )=P\gamma (t)P\gamma (\tau )`$ for $`t0`$.

The second long standing conjecture on closed characteristics is whether there always exists at least an elliptic closed characteristic on any $`\mathrm{\Sigma }(2n)`$. Up to the author’s knowledge, the existence of one elliptic closed characteristic on $`\mathrm{\Sigma }(2n)`$ was proved by I. Ekeland in 1990 when $`\mathrm{\Sigma }`$ is $`\sqrt{2}`$-pinched by two spheres, and by G.-F. Dell’Antonio, B. D’Onofrio, and I. Ekeland in 1992 when $`\mathrm{\Sigma }`$ is symmetric with respect to the origin. Recently using an enhanced version of the iteration estimate (2.9) on the elliptic height, based on results in  the following result was further proved by C. Zhu and the author.

Theorem 12 (cf. ). For $`\mathrm{\Sigma }(2n)`$, suppose $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })<+\mathrm{}`$. Then there exists at least an elliptic closed characteristic on $`\mathrm{\Sigma }`$. Moreover, suppose $`n2`$ and $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })2[n/2]`$. Then there exist at least two elliptic elements in $`𝒯(\mathrm{\Sigma })`$.

The main ingredient in the proofs of Theorems 9 to 12 is our index iteration theory mentioned above. To illustrate this method, we briefly describe below the main idea in the proof of (3.3) in Theorem 9. Because each closed characteristic on $`\mathrm{\Sigma }`$ corresponds to infinitely many critical values of the related dual action functional, our way to solve the problem is to study how the index intervals of iterated closed characteristics cover the set of integers $`2𝐍2+n`$ to count the number of closed characteristics on $`\mathrm{\Sigma }`$. Suppose $`q=^\mathrm{\#}\stackrel{~}{𝒥}(\mathrm{\Sigma })<+\mathrm{}`$. In the proof of the multiplicity claim (3.3) of Theorem 9, the most important ingredient is the following estimates:

$`q`$ $``$ $`{}_{}{}^{\mathrm{\#}}((2𝐍2+n)_{j=1}^q𝒢_{2m_j1}(\gamma _{x_j}))`$ (3.4)
$``$ $`{}_{}{}^{\mathrm{\#}}((2𝐍2+n)[2N\kappa _1,2N+\kappa _2])`$
$``$ $`[{\displaystyle \frac{n}{2}}]+1,`$
The first inequality in (3.4) is a new version of the Liusternik-Schnirelman theoretical argument at the iterated index level, which distinguishes solution orbits geometrically instead of critical points only as usual methods do. The second inequality in (3.4) uses the common index jump Theorem 6. The last inequality in (3.4) uses the Morse theoretical approach. Roughly speaking, the common index jump theorem picks up as many as possible points of $`2𝐍2+n`$ in the interval $`[2N\kappa _1,2N+\kappa _2]`$ $``$ $`_{j=1}^q𝒢_{2m_j1}(\gamma _{x_j})`$, which yields a lower bound for $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })`$.

As usual, a hypersurface $`\mathrm{\Sigma }𝐑^{2n}`$ is star-shaped if the tangent hyperplane at any $`x\mathrm{\Sigma }`$ does not intersect the origin. Closed characteristics on $`\mathrm{\Sigma }`$ can be defined by (3.1) too. In this case, the result $`{}_{}{}^{\mathrm{\#}}\stackrel{~}{𝒥}(\mathrm{\Sigma })1`$ was proved by P. Rabinowitz in  of 1978. Then multiplicity results were proved under certain pinching conditions on star-shaped $`\mathrm{\Sigma }`$. Recently, the following result for the free case was proved by X. Hu and the author:

Theorem 13 (cf. ). Let $`\mathrm{\Sigma }`$ be a star-shaped compact $`C^2`$-hypersurface in $`𝐑^{2n}`$. Suppose all the closed characteristics on $`\mathrm{\Sigma }`$ and all of their iterates are non-degenerate. Then $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })2`$. Moreover, if $`n=2`$ and $`{}_{}{}^{\mathrm{\#}}𝒯(\mathrm{\Sigma })<+\mathrm{}`$ further holds, then there exist at least two elliptic closed characteristics on $`\mathrm{\Sigma }`$.

Here the crucial point is to prove $`i(x,1)n`$ when $`(\tau ,x)`$ is the only geometrically distinct closed characteristic on $`\mathrm{\Sigma }`$. This conclusion is proved by using our index iteration theory and an identity of non-degenerate closed characteristics on $`\mathrm{\Sigma }`$ proved by C. Viterbo in 1989.

Because of Theorem 9 and other indications, we suspect that the following holds:

$$\{^\mathrm{\#}𝒯(\mathrm{\Sigma })|\mathrm{\Sigma }(2n)\}=\{k𝐙|[\frac{n}{2}]+1kn\}\{+\mathrm{}\}.$$
(3.5)
We also suspect that closed orbits of the Reeb field on a compact contact hypersurfaces in a symplectic manifold may have similar properties.

Many other problems related to iterations of periodic solution orbits are still open, for example, the Seifert conjecture on the existence of at least $`n`$ brake orbits for the given energy problem of classical Hamiltonian systems on $`𝐑^n`$ (cf. ,  and the references there in), and the conjecture on the existence of infinitely many geometrically distinct closed geodesics on every compact Riemannian manifold (cf.  and the solution for $`S^2`$ by J. Franks and V. Bangert). We believe that our index iteration theory for symplectic paths and the methods we developed to establish and apply it to nonlinear problems will have the potential to play more roles in the study on these problems and in other mathematical areas.

Acknowledgements. The author sincerely thanks the 973 Program of MOST, NNSF, MCME, RFDP, PMC Key Lab of MOE of China, S. S. Chern Foundation, CEC of Tianjin, and Qiu Shi Sci. Tech. Foundation of Hong Kong for their supports in recent years.


text:
# Massive Variability Surveys from Venezuela

## 1. Introduction

Recent surveys focusing on repeated observations of large expanses of the night sky (e.g. ROTSE, Akerlof et al. 2000) are transforming the study of variability phenomena one of the forefront research areas in modern day astrophysics. These searches have also demonstrated the unique value of small telescopes for this kind of research; in particular, the wide field Schmidt design, when equipped with mosaic CCD detectors, is a potentially ideal combination to collect multi-epoch data over extended areas of the sky in an efficient way. In this contribution I describe one of the variability surveys, spanning many hundreds of square degrees near the celestial equator, that are being carried at the Venezuela National Astronomical Observatory using an 8k x 8k CCD Mosaic Camera, optimized for drift-scanning, installed on a 1m Schmidt telescope.

## 2. A Variability Survey of the Orion OB1 Association

The study of how stars and planetary systems form and evolve is an area in which variability surveys are particularly useful as a means to search for pre-main sequence objects. Young (1-10 Myr) low mass ($`<2M_{}`$) stars are intrinsically variable at both optical (Herbst 1986) and near IR wavelengths (Carpenter, Hillenbrand & Strutskie 2001) by up to several magnitudes, thus, variability should help to single out these objects among the general field population. However, so far variability has been used mainly in follow-up studies rather than as a survey technique, to collect time-series data of young objects that were detected by some other method, like color-magnitude diagrams, X-ray and H$`\alpha `$ surveys.
In collaboration with Nuria Calvet (CfA), A.K. Vivas (CIDA) and Lee Hartmann (CfA), we have used the 8k $`\times `$ 8k CCD Mosaic Camera on the Venezuela 1m Schmidt telescope (see Baltay et al. 2002 for a description of the camera; also Vivas, A.K. in this volume) to conduct a large scale, systematic optical VRIH$`\alpha `$ survey spanning $`180`$ sqr.deg. in the Orion OB1 Association (Fig.1), one of the nearest ($`d400`$ pc, Brown et al. 1994) star forming regions. Though Orion has been extensively studied in the past, existing optical/IR studies in Orion have concentrated mostly on small regions such as the Orion Nebula Cluster (ONC; e.g., Hillenbrand 1997) and the surroundings of the star $`\sigma `$ Ori (Walter et al. 1998). A few large scale studies (Wiramihardja et al. 1991; Alcalá et al. 1996) have been done but did not find the more spread out, slightly older lower mass young stars. Only by identifying and studying somewhat older ($`10`$ Myr), low-mass stars in widely spread stellar populations, can we address fundamental questions like the duration of the planet formation phase in dusty disks around young stars, and the lifetimes of molecular clouds.

We have developed tools for identifying variable stars using differential photometry. With a $`\chi ^2`$ test and assuming a Gaussian distribution for the errors, we consider variable only those objects for which the probability that the observed distribution is a result of the random errors is $`<0.01`$%. The value of variability for picking out young stars is shown in Fig.2. Follow up spectra of a large number of candidate young stars indicate that our selection technique is very efficient: roughly 50-60% of our candidates turn out to be young, low mass stars. Our initial results have been described in Briceño et al. (2001).

## 3. The Equatorial Variability Survey

The exciting results in Orion, and those presented in this same conference by A.K. Vivas, have prompted us to extend our spatial coverage to span all the sky near the Celestial Equator, from 0-24h. We have now accumulated 4 years of data for a strip centered at $`\delta =1^o`$. Presently we are scanning a new strip centered at $`\delta =3^o`$, which we plan to observe repeatedly during the next 2 yrs. The actual temporal sampling is dictated by the number of nights available for our project, usually 8-10 nights/month, and the weather (Fig. 2c). Using the YALO telescope we have obtained data of selected fields along our survey strips for the photometric calibration of the data. The variability database will provide a unique research tool for our studies of nearby star forming regions, RR Lyrae in the Halo of the Milky Way, and Quasar variability, among others. In collaboration with E.A. Ponsot at Universidad de los Andes, Mérida, Venezuela, we are developing a variability database from the data collected in our surveys, together with web-based tools to mine the data. The variability database will provide a unique research tool for our studies of nearby star forming regions, RR Lyrae in the Halo of the Milky Way, and Quasar variability, among others. We expect to release an initial version of the database to the astronomical community in the near future.

### Acknowledgments.

This work is supported by grants S1-2001001144 of FONACYT, Venezuela, AST 9987367 of the National Science Foundation, and NASA Origins NGC5-10545, USA.

## References

Alcalá J.M. et al., 1996, A&AS, 199, 7

Akerlof, C. et al. 2000, AJ, 119, 1901

Baltay et al. 2002, PASP, 114, 780

Briceño, C. et al., 2001, Science, 291, No. 5501, 93

Carpenter, J.M., Hillenbrand, L.A. & Strutskie, M.F., 2001, AJ, 121, 3160

Brown, A.G.A., de Geus, E.J., & de Zeeuw, P.T. 1994, AA, 289, 101

Herbst, W. 1986, PASP, 98, 1088

Hillenbrand, L. 1997, AJ, 113, 1733

Walter, F.M., Wolk, S.J., Sherry, W., 1998, in Cool Stars 10, eds. R.A. Donahue and J.A. Bookbinder, ASP Conf. Ser., 154, 1793

Wiramihardja, S., Kogure, T., Yoshida, S., Ogura, K., & Nakano, M. 1993, PASJ, 45, 643


text:
# Fourier space intermittency of the small-scale turbulent dynamo

## I Introduction

In some astrophysical applications such as the interstellar medium and protogalactic plasmas the kinematic viscosity is greater than the magnetic field diffusivity by the factor $`10^{14}`$ to $`10^{22}`$ . In these situations, there is a vast scaling range where the magnetic field has a smaller characteristic length-scale than the velocity field. The dynamo process of the stochastic stretching and amplification of the magnetic field can be studied in this case in terms of the statistics of Lagrangian deformations . It can be pictured as a collection of magnetic field wavepackets, each moving along a fluid particle path and being distorted by the local strain.***The second and higher derivatives of the velocity field can be ignored in this case. Such a regime of smooth velocity fields is similar to the Batchelor’s regime found in the related problem of passive scalar advection. The Batchelor’s regime has a long history of study, some of the more recent advances can be found, for example, in  (and references therein). A common further simplification of the problem is to assume that the local strain matrix is a Gaussian white noise process, this is commonly known as the Kazantsev-Kraichnan model . However, this assumption is not always necessary, and some results have been shown to be universal for a broader class of stochastic flows .

The first analysis of the dynamo problem were based in Fourier space , with particular emphasis being placed on the second order moment, corresponding to the energy spectrum. It was realized, however, that the Fourier space and coordinate space equations have a similar structure. Recent studies of the turbulence intermittency (in both the passive scalar and dynamo problems) have focussed on coordinate space moments of second order or higher. This approach was motivated by a feeling that this description is more natural and can give more information than the $`k`$-space moments . As a result, the only $`k`$-space correlator that has been seriously studied to date is the energy spectrum and there is no theory describing the higher Fourier space correlators.

In the present paper, we turn our attention back to Fourier space and consider the one-point moments of the magnetic field Fourier transform. Although the problem under investigation here is the small-scale turbulent dynamo, we would like to give several reasons why the description of the Fourier space moments is important and why it should be developed in a broader turbulence context.

(i) The presence of singular structures in turbulence is known to affect the scaling of the structure functions. The most well-known examples of this being the $`\beta `$ and multi-fractal models . However, some coherent structures are singular in Fourier space rather than in coordinate space, and therefore can be detected by an investigation of the Fourier moments. A simple example of structures that are singular in Fourier space (and regular in coordinate space) is a sea of vortex filaments in 2D turbulence. The layered pattern of the vortices in coordinate space corresponds to a 1D curve in Fourier space.The Fourier transform should be taken over a local box, which is large enough to fit many layers of filaments, but is smaller than the large-scale vortices in this case. This will be discussed in the next section.

(ii) To date, the only Fourier moment that was studied in detail is the second order moment that describes the distribution of the turbulent energy over wavenumbers. However, there exist other Fourier moments that have a clear physical meaning and describe important properties of the turbulence. One of these objects is related to the fourth order moments and has a physical meaning of the mean turbulence polarization. It will be introduced and studied in this paper.

(iii) Phases of the Fourier modes can be dealt with directly, and therefore the validity of the random phase assumption can be examined.

(iv) Finally, in some cases a Fourier space analysis is the only way to have a treatable problem due to the introduction of non-localities, for example via a pressure or wave dispersion. Note that the dynamo and passive scalar systems do not fit into this class of problems. However, the methods developed in the present paper will be applied to the non-local Navier-Stokes equationSuch an equation is the basis of Rapid Distortion theory which will be developed for the case of a stochastic strain. (which involves a pressure term) in our next paper .

## II One-point and two-point correlators of Fourier amplitudes

To avoid any confusion, we should state that the Fourier transforms in this paper are performed over a finite box of a size much greater than typical length-scale of the magnetic turbulence, but much less than the scale of the advecting velocity field. It is only for finite box Fourier transforms that one-point correlators are well defined objects. Secondly, the center of each finite box is not fixed. The box center coordinate dependence of the Fourier transforms is a useful measure of the slow variability of the large-scale magnetic turbulence. In this paper we make each box move in unison with a fluid particle, and thus any coordinate dependence is replaced by a time dependence along a given fluid trajectory.

One should be aware that in this case the finite box size does not produce a discrete spectrum as the magnetic field is in general not periodic within the box. On the other hand, the lack of periodicity can produce unphysical power-law tails in the Fourier transforms at (large) wavenumbers that are greater than one over the box length. Slow decay of these tails can be undesirable in numerical simulations and other applications. However, this problem can be easily cured by using a smooth filter instead of the box, for example by using Gabor transforms instead . In this presentation the box (or filter) shape is actually not an issue. We will present the results in terms of finite box Fourier transforms. However, performing the analysis using Gabor transforms will in fact produce identical results.

Let us consider the general two-point correlator of the Fourier transformed magnetic field components

$$B_{i_1}(𝐤_1)B_{i_2}(𝐤_1)\mathrm{}B_{i_n}(𝐤_1)B_{j_1}(𝐤_2)B_{j_2}(𝐤_2)\mathrm{}B_{j_n}(𝐤_2),$$
(1)
where $`i_1,i_2,\mathrm{},i_n`$ and $`j_1,j_2,\mathrm{},j_n`$ take the values 1, 2 or 3 corresponding to the three components of the magnetic field, and $`n`$ is an arbitrary natural number. Let us assume the turbulence is quasi-homogeneous. That is, the correlator $`B_i(𝐱_1)B_j(𝐱_2)`$ (and similar) depend only on $`𝐱_1+𝐱_2`$ as slowly as the advecting strain field. On the other hand, such correlators decay in $`𝐱_1𝐱_2`$ at distances much less than the box size $`L`$. In this case

$`B_{i_1}(𝐤_1)`$ $`B_{i_2}(𝐤_1)\mathrm{}B_{i_n}(𝐤_1)B_{j_1}(𝐤_2)B_{j_2}(𝐤_2)\mathrm{}B_{j_n}(𝐤_2)=`$ (3)
$`B_{i_1}(𝐤_1)B_{i_2}(𝐤_1)\mathrm{}B_{i_n}(𝐤_1)B_{j_1}(𝐤_1)B_{j_2}(𝐤_1)\mathrm{}B_{j_n}(𝐤_1)f(𝐤_1+𝐤_2),`$
where $`f(𝐤)=\frac{1}{L^d}\mathrm{\Pi }_{l=1}^d\frac{\mathrm{sin}(Lk_l)}{k_l}`$ is the Fourier transform of the filter function and $`d`$ is the number of space dimensions. Here we see that any two-point correlator is fully determined by the one-point correlators

$$I_{i_1,i_2,\mathrm{},i_n;j_1,j_2,\mathrm{},j_n}(𝐤)=B_{i_1}(𝐤)B_{i_2}(𝐤)\mathrm{}B_{i_n}(𝐤)B_{j_1}(𝐤)B_{j_2}(𝐤)\mathrm{}B_{j_n}(𝐤).$$
(4)
This simple observation means we can avoid lengthy derivations dealing with two-point objects directly, and instead can obtain results from the one-point correlators that are themselves far easier to deal with. Turbulence isotropy and the divergence-free condition further narrow the class of possible one-point (and two-point) correlator tensors we need to consider (see for example )

$$I_{i_1,i_2,\mathrm{},i_n;j_1,j_2,\mathrm{},j_n}(𝐤)=\underset{\mathrm{pairs}}{}\mathrm{\Phi }_s^n(k)D_{l_1,l_2}D_{l_3,l_4}\mathrm{}D_{l_{n1},l_n},$$
(5)
where

$$D_{ij}\delta _{i,j}\frac{k_ik_j}{k^2},$$
(6)
and the summation is over the set of all possible pairings $`(l_1,l_2),(l_3,l_4),\mathrm{},(l_{n1},l_n)`$ of the indices $`i`$ and $`j`$ from the left hand-side of (5). For each pairing, the index $`s`$ is equal to the number of pairs that consist only of $`i`$’s, §§§ There are obviously going to be the same number of pairs that consist only of $`j`$’s. for example $`(i_1,i_3)`$ or $`(i_8,i_2)`$. One can see that any correlator of order $`2n`$ can be expressed in terms of $`N=\text{Max}(s)+1=\text{Int}(n/2)+1`$ independent functions $`\mathrm{\Phi }_s^n(k)`$. There is only one such function for the second order correlators, two functions for the fourth and sixth orders, and three functions at the eighth and tenth orders.

Instead of $`\mathrm{\Phi }_s^n(k)`$, we can choose another set of independent functions, in particular the following set of correlators

$$\mathrm{\Psi }_s^n=|𝐁(𝐤)|^{(2n4s)}|𝐁(𝐤)^2|^{2s}.$$
(7)
These correlators are linear combinations of $`\mathrm{\Phi }_s^n(k)`$’s, for example at the fourth order we have

$`\mathrm{\Psi }_0^{(2)}`$ $`=`$ $`d(d1)\mathrm{\Phi }_0^{(2)}+(d1)\mathrm{\Phi }_1^{(2)},`$ (8)
$`\mathrm{\Psi }_1^{(2)}`$ $`=`$ $`2(d1)\mathrm{\Phi }_0^{(2)}+(d1)^2\mathrm{\Phi }_1^{(2)}.`$ (9)
Below, we will derive a technique that allows us to deal with the fundamental set of correlators (7).

## III The model

Let us start with the equation for the magnetic field $`𝐁(𝐱,t)`$,

$$_t𝐁+(𝐮)𝐁=(𝐁)𝐮+\kappa \mathrm{\Delta }𝐁,$$
(10)
where $`𝐮`$ is the velocity field, which we assume has a much slower spatial variation than the magnetic field $`𝐁`$, and $`\kappa `$ is the magnetic diffusivity (determined by the fluid conductivity). We will make use of the previously discussed finite box Fourier transforms. Each box has sides of length $`L`$ which is chosen to lie between the length-scales associated with $`𝐁`$ and $`𝐮`$, namely $`L_B`$ and $`L_u`$. The box Fourier transform is defined as

$$\widehat{𝐁}(𝐤,𝐱,t)_{box}𝐁(𝐫,t)e^{i𝐤(𝐱𝐫)}𝑑𝐫,$$
(11)
where $`𝐱`$ is the coordinate of the box center. Applying this Fourier transform to (10) we have (with accuracy up to the first order of the scale separation parameter $`ϵ=L_B/L_u1`$)

$$_t\widehat{B}_m+u_i_i\widehat{B}_m=\sigma _{ij}k_i_j\widehat{B}_m+\sigma _{mi}\widehat{B}_i\kappa k^2\widehat{B}_m,$$
(12)
where $`\sigma _{ij}=_ju_i`$ is the strain matrix and the operators $`_i`$ and $`_i`$ correspond to derivatives with respect to $`x_i`$ and $`k_i`$ respectively. In considering a fluid path determined by $`\dot{𝐱}(t)=𝐮`$, and thus noting that $`\sigma _{ij}(𝐱,t)\sigma _{ij}(𝐱(t),t)`$ and $`\widehat{𝐁}(𝐤,𝐱,t)\widehat{𝐁}(𝐤,𝐱(t),t)`$, this equation becomesHereafter, we will the drop hat notation on $`\widehat{𝐁}`$ because only Fourier components will be considered. Also, we will not mention explicitly dependence on the fluid path and simply write $`𝐁𝐁(𝐤,t)`$.

$$_tB_m=\sigma _{ij}k_i_jB_m+\sigma _{mi}B_i\kappa k^2B_m.$$
(13)
It should also be noted that the strain $`\sigma _{ij}`$, taken along a fluid path, only enters this equation as a given function of time. To complete the model one has to specify this dependence or postulate the strain statistics. In what follows we choose the strain matrix to be Gaussian such that

$$\sigma _{ij}=\mathrm{\Omega }\left(A_{ij}\frac{A_{ll}}{d}\delta _{ij}\right),$$
(14)
where $`A_{ij}`$ is a matrix with statistically independent elements that are white in time

$$A_{ij}(t)A_{kl}(0)=\delta _{ij}\delta _{kl}\delta (t).$$
(15)
This choice of strain matrix ensures incompressibility and statistical isotropy. Indeed,

$$\sigma _{ij}(t)\sigma _{kl}(0)=\mathrm{\Omega }\left(\delta _{ik}\delta _{jl}\frac{1}{d}\delta _{ij}\delta _{kl}\right)\delta (t).$$
(16)
However, one must be aware that this is not the only way to satisfy the incompressibility and isotropy conditions; for example  have chosen

$$\sigma _{ij}(t)\sigma _{kl}(0)=\mathrm{\Omega }\left[(d+1)\delta _{ik}\delta _{jl}\delta _{il}\delta _{jk}\delta _{ij}\delta _{kl}\right]\delta (t).$$
(17)
More generally, there is an infinite one-parametric family of possible strains. The strain matrices that are Gaussian and white in time correspond to the Kazantsev-Kraichnan model, a natural starting point for an analytical analysis because of its simplicity. However, it should be noted that some results remain universal in the case of other smooth velocity fields found in a much wider class of statistical models . In our future work we will investigate the behaviour of the Fourier space correlators in the case of more general statistics.

## IV Generating Function

Let us consider the following generating function

$$Z(\lambda ,\alpha ,\beta ,k)=e^{\lambda |𝐁(𝐤)|^2+\alpha 𝐁^2(𝐤)+\beta \overline{𝐁}^2(𝐤)},$$
(18)
where the overline denotes the complex conjugation. This function allows one to obtain any of the fundamental one-point correlators (7) via differentiation with respect to $`\lambda ,\alpha `$ and $`\beta `$

$$\mathrm{\Psi }_s^n=|𝐁(𝐤)|^{(2n4s)}|𝐁^2(𝐤)|^{2s}=\left[_\lambda ^{(2n4s)}_\alpha ^s_\beta ^sZ\right]_{\lambda =\alpha =\beta =0}.$$
(19)
Differentiating (18) with respect to time and using the dynamical equation (13), we have

$`\dot{Z}=k_i_j\sigma _{ij}E+`$ $`\lambda \sigma _{ml}(\overline{B}_mB_l+\overline{B}_lB_m)E+2\alpha \sigma _{ml}B_mB_lE+2\beta \sigma _{ml}\overline{B}_m\overline{B}_lE`$ (21)
$`2\kappa k^2[\lambda |𝐁(𝐤)|^2+\alpha 𝐁^2(𝐤)+\beta \overline{𝐁}^2(𝐤)]E,`$
where

$$E=e^{\lambda |𝐁(𝐤)|^2+\alpha 𝐁^2(𝐤)+\beta \overline{𝐁}^2(𝐤)}.$$
(22)
To find the correlators on the right hand-side of (21), one needs to make use of the Gaussianity of the strain matrix $`\sigma _{ij}`$ and perform a Gaussian integration by parts. We then use the whiteness of the strain field to find the response function (functional derivative of $`B_l`$ with respect to $`\sigma _{ij}`$). Finally, one can use the statistical isotropy of the strain, so that the final equation involves only $`k=|𝐤|`$ and no angular dependence of the wavevector. This derivation is discussed in more detail in the Appendix. Here, we just write the final result

$`\dot{Z}`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Omega }}{2}}[(1{\displaystyle \frac{1}{d}})k^2Z_{kk}+{\displaystyle \frac{1}{d}}(4𝒟+d^21)kZ_k`$ (24)
$`+(2d6)𝒟Z+4(1{\displaystyle \frac{1}{d}})𝒟^2Z+2(\lambda ^24\alpha \beta )(Z_{\alpha \beta }Z_{\lambda \lambda })]2\kappa k^2𝒟Z,`$
where the $`k,\alpha ,\beta `$ and $`\lambda `$ subscripts in $`Z`$ denote differentiation with respect to $`k,\alpha ,\beta `$ and $`\lambda `$ respectively, and

$$𝒟=\lambda _\lambda +\alpha _\alpha +\beta _\beta .$$
(25)
Actually, one can reduce the number of independent variables in this equation by taking into account that $`Z`$ depends on $`\alpha `$ and $`\beta `$ only in the combination $`\eta =\alpha \beta `$. The easiest way to see this is to consider a Taylor series expansion of $`Z`$ from (18); any term that contains $`\alpha `$ and $`\beta `$ in a different combination will be zero because of the quasi-homogeneity of the turbulence. Thus, we can write

$`\dot{Z}`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Omega }}{2}}[(1{\displaystyle \frac{1}{d}})k^2Z_{kk}+{\displaystyle \frac{1}{d}}(4𝒟+d^21)kZ_k`$ (27)
$`+(2d6)𝒟Z+4(1{\displaystyle \frac{1}{d}})𝒟^2Z+2(\lambda ^24\eta )(Z_\eta +\eta Z_{\eta \eta }Z_{\lambda \lambda })]2\kappa k^2𝒟Z,`$
where

$$𝒟=\lambda _\lambda +2\eta _\eta .$$
(28)
In what follows, we will restrict our consideration to the 3D case ($`d=3`$) in which case the equation for $`Z`$ reduces to

$`\dot{Z}`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Omega }}{3}}[k^2Z_{kk}+(42𝒟)kZ_k`$ (30)
$`+4𝒟^2Z+3(\lambda ^24\eta )(Z_\eta +\eta Z_{\eta \eta }Z_{\lambda \lambda })]2\kappa k^2𝒟Z.`$
## V Energy spectrum

Let us now consider the energy spectrum of the magnetic turbulence, given by the second order correlator

$$E(k,t)=\mathrm{\Psi }_0^1=|𝐁(𝐤)|^2=\left[_\lambda Z\right]_{\lambda =\eta =0}.$$
(31)
Differentiating (30) with respect to $`\lambda `$ and taking the result at $`\lambda =\eta =0`$ we have

$$\dot{E}=\frac{\mathrm{\Omega }}{3}(k^2E_{kk}+2kE_k+4E)2\kappa k^2E.$$
(32)
This equation for the evolution of the energy spectrum was first obtained by Kazantsev  and Kraichnan and Nagarajan . Kazantsev analyzed an eigenvalue problem associated with this equation which allowed him to obtain the growth exponents of the total magnetic energy. Numerically, the energy spectrum was studied by Kulsrud and Anderson  who gave a detailed description of the $`k`$-space evolution of this spectrum. Recently, Schekochihin, Boldyrev and Kulsrud  presented the complete solution of (32) obtained by the use of Kontorovich-Lebedev transforms (KLT). Note that this integral transform approach has an advantage over the Kazantsev’s eigenvalue analysis in that it allows us to obtain not only the growth exponents, but also the power-law prefactors, of the large time asymptotic solutions.

As they will be of use later in this presentation, let us briefly review the previous results for the energy spectrum before we discuss the higher order moments. Using the following substitution

$$E=e^{5\mathrm{\Omega }t/4}k^{1/2}\varphi (k/k_{d_1},t)k_{d_n}=\sqrt{\frac{\mathrm{\Omega }}{6n\kappa }},$$
(33)
one can reduce (32) to

$$\frac{3}{\mathrm{\Omega }}\dot{\varphi }(p,t)=p^2\varphi _{pp}+p\varphi _pp^2\varphi ,$$
(34)
where $`p=k/k_{d_1}`$. At scales much greater than the dissipative one, $`p1`$, there is a perfect conductor regime for time $`t(\mathrm{ln}q)^2`$ (where $`q1`$ is the mean wavenumber of the initial condition). Thus, in this regime the last term in (32) can be neglected. By changing to logarithmic coordinates and a moving frame of reference one can transform this equation into a heat equation. For $`t1/\mathrm{\Omega }`$ the solution of which is just the Green’s function, which gives

$$\varphi =\text{const}t^{1/2}e^{\frac{3(\mathrm{ln}k/q)^2}{4\mathrm{\Omega }t}},$$
(35)
where the constant is fixed by the initial condition. This solution describes a spectrum with an expanding $`k^{1/2}`$ scaling range. At $`t(\mathrm{ln}q)^2`$ the front of this scaling range reaches the dissipative scales. To solve (32) in this case, we note that the right hand-side of this equation is just the modified Bessel operator and by using KLT one immediately obtains

$$\varphi (p,t)=\text{const}_0^{\mathrm{}}𝑑ss\text{sinh}(\pi s)K_{is}(p)K_{is}(q)e^{s^2t},$$
(36)
where $`K_{is}`$ is a MacDonald function of imaginary order. Again this solution is given by the Green’s function only because the condition that $`t1/\mathrm{\Omega }`$ is obviously satisfied if $`q1`$. For time $`t(\mathrm{ln}q)^2`$, the function $`e^{s^2t}`$ is strongly peaked at $`t=0`$ and integration of (36) gives

$$\varphi =\text{const}t^{3/2}K_0(p).$$
(37)
One should note that although the $`K_0(p)`$ shape is somehow predicted by the Kazantsev’s eigenmode analysis , the $`t^{3/2}`$ factor can only be obtained by solving the full initial value problem. For $`p1`$, $`K_0(p)\mathrm{ln}p`$, which means that at large times the scales far larger than the dissipative one are affected by diffusion via a logarithmic correction

$$E(k)k^{1/2}\mathrm{ln}(k_{d_1}/k).$$
(38)
The energy evolution gives an important, but incomplete picture, of the dynamo process. In particular, it does not capture the existence of small-scale intermittency and does not allow us to predict the type of coherent structures dominating the turbulence at large time. To deal with these issues one has to study higher order correlators. Higher one-point correlators in coordinate space were studied in  and used to predict the shape of the dominant structures. Below, we proceed to study the higher $`k`$-space correlators. In particular, this will lead to the discovery of a new quantity of interest, corresponding to the mean polarization of the magnetic turbulence.

## VI 4th-order correlators, turbulence polarization and flatness

There are two independent 4th order correlators

$$S(k,t)=\mathrm{\Psi }_0^{(2)}=|𝐁(𝐤)|^4=\left[Z_{\lambda \lambda }\right]_{\lambda =\eta =0},$$
(40)
$$T(k,t)=\mathrm{\Psi }_1^{(2)}=|𝐁^2(𝐤)|^2=\left[Z_\eta \right]_{\lambda =\eta =0}.$$
(41)
Differentiating (30) twice with respect to $`\lambda `$ and taking the result at $`\lambda =\eta =0`$ we have

$$\dot{S}=\frac{\mathrm{\Omega }}{3}(k^2S_{kk}+10S+6T)4\kappa k^2S.$$
(42)
Now, differentiating (30) with respect to $`\eta `$ and taking the result at $`\lambda =\eta =0`$ we get

$$\dot{T}=\frac{\mathrm{\Omega }}{3}(k^2T_{kk}+4T+12S)4\kappa k^2T.$$
(43)
Equations (42) and (43) make up a complete system for $`S`$ and $`T`$ and can be solved exactly in the general case. Observe that there is a closed equation for $`W=ST`$

$$\dot{W}=\frac{\mathrm{\Omega }}{3}(k^2W_{kk}2W)4\kappa k^2W.$$
(44)
Before solving this equation, let us examine the physical meaning of $`W`$ by writing it as

$$W=4\underset{jl}{\overset{3}{}}[\mathrm{}(B_j\overline{B}_l)]^2=4\underset{jl}{\overset{3}{}}|B_j|^2|B_l|^2\mathrm{sin}^2(\varphi _j\varphi _l)0,$$
(45)
where $`\mathrm{}`$ denotes the imaginary part and $`\varphi _j`$ and $`\varphi _l`$ are the phases of the components $`B_j`$ and $`B_l`$ respectively. We see therefore that $`W`$ contains information not only about the amplitudes but also about the phases of the Fourier modes. In particular, $`W0`$ corresponds to the case where all Fourier components of the magnetic field are plane polarization. If $`W0`$ then other polarizations (circular, elliptic) are present. This is the case for example for a Gaussian field where one finds $`W=E^2/2>0`$. On the other hand, the smallness of the phase differences in $`W`$ can be overpowered by large amplitudes. Therefore, a better measure of the mean polarization would be a normalized $`W`$, for example

$$P=W/S.$$
(46)
Defined in this way, the mean turbulence polarization is an example of an important physical quantity which can be obtained from the one-point Fourier correlators, and that is unavailable from the (one-point or two-point) coordinate space correlators.

In the perfect conductor regime, when the diffusivity term in the equations (42) – (44) can be ignored, the solution for $`W`$ can be obtained in a similar manner to the previous energy spectrum analysis. Namely, by reducing (44) into a heat equation via a logarithmic change of variables and passing into a moving frame of reference. The solution therefore is

$$W=W_0t^{1/2}e^{3\mathrm{\Omega }t/4}k^{1/2}e^{\frac{3(\mathrm{ln}k/q)^2}{4\mathrm{\Omega }t}},$$
(47)
where $`W_0`$ is a constant which can be found from the initial conditions. We see that $`W`$ develops a $`k^{1/2}`$ scaling range which is cut-off at low and high $`k`$ by exponentially propagating fronts. Within this scaling range, $`W`$ decays exponentially in time.

Given $`W`$, one can also find $`S`$ by representing it as $`S=V+cW`$ and choosing the constant $`c`$ such that the equation for $`V`$ is closed. This gives $`c=3/7`$ and

$`S`$ $`=`$ $`t^{1/2}k^{1/2}e^{\frac{3(\mathrm{ln}k/q)^2}{4\mathrm{\Omega }t}}\left(V_0e^{21\mathrm{\Omega }t/4}+{\displaystyle \frac{3}{7}}W_0e^{3\mathrm{\Omega }t/4}\right),`$ (48)
where $`V_0`$ is another constant which can be found from the initial conditions. For $`t1/\mathrm{\Omega }`$, the second term in the parenthesis should be neglected, and we have the following solution for the mean turbulence polarization

$$P=W/S=\frac{W_0}{V_0}e^{6\mathrm{\Omega }t}.$$
(49)
In the perfect conductor regime we see the mean polarization tends to a value that is independent of $`k`$ and decays exponentially in time. This means that all Fourier modes of the magnetic field eventually become plane polarized. Recall that such turbulence is very far from being Gaussian, remembering that the mean polarization of a Gaussian field remains finite.That is, elliptic and circular polarized modes are present.

It is also easy to obtain solutions to $`W`$ for the diffusive regime. Indeed, following the example of the energy spectrum we make use of the substitutions

$`W`$ $`=`$ $`e^{3\mathrm{\Omega }t/4}k^{1/2}\varphi (k/k_{d_2}),`$ (50)
$`V`$ $`=`$ $`e^{21\mathrm{\Omega }t/4}k^{1/2}\varphi (k/k_{d_2}),`$ (51)
to transform the governing equations for $`W`$ and $`V`$ into a form similar to (34). In each case we can solve the equation for $`\varphi `$ using KLT and find that for times $`t(\mathrm{ln}q)^2`$

$`W`$ $``$ $`W_0k^{1/2}t^{3/2}e^{3\mathrm{\Omega }t/4}K_0(k/k_{d_2}),`$ (52)
$`V`$ $``$ $`V_0k^{1/2}t^{3/2}e^{21\mathrm{\Omega }t/4}K_0(k/k_{d_2}).`$ (53)
Correspondingly, $`S`$ has the solution

$$Sk^{1/2}t^{3/2}\left[V_0e^{21\mathrm{\Omega }t/4}+\frac{3}{7}W_0e^{3\mathrm{\Omega }t/4}\right]K_0(k/k_{d_2}),$$
(54)
while importantly we find that the normalised polarisation $`P=W/S`$ behaves identically in both the diffusive and perfect conductor regimes.

Therefore, in the diffusive regime $`W`$ continues to decrease in time with the same exponential rate as perfect conductor case and has the same $`t^{3/2}`$ prefactor as the energy spectrum. Thus, by the time the diffusive regime is achieved $`W`$ can be essentially put equal to zero. The fact that the polarization becomes plane has quite a simple physical explanation. Indeed, a magnetic field wavepacket of arbitrary polarization will be strongly distorted by stretching. The stretching being strongest along the direction of the dominant eigenvector of the Lagrangian deformation matrix (corresponding to the greatest Lyapunov exponent). Such a stretching will make any initial “spiral” structure flat at large time, with the dominant field component lying in a plane passing through the eigenvector stretching and wavevector directions.\**\**\**Of course, the incompressibility of $`𝐁(𝐤)`$ ensures that it is perpendicular to $`𝐤`$.

Another measure of intermittency in turbulence is the flatness which can be defined in $`k`$-space as $`F=S/E^2`$. For large time, in the perfect conductor regime

$$Ft^{1/2}e^{11\mathrm{\Omega }t/4}k^{3/2}.$$
(55)
We see that the flatness grows both in time and in $`k`$ which indicates the presence of small-scale intermittency. Such an intermittency can be attributed to the presence of coherent structures in $`k`$-space. In the diffusive regime, in contrast to the polarisation the behaviour of the flatness is modified and at large time one finds

$$Fk^{3/2}t^{3/2}e^{11\mathrm{\Omega }t/4}\frac{K_0(k/k_{d_2})}{[K_0(k/k_{d_1})]^2}.$$
(56)
For small $`k`$, we again have a region of $`k^{3/2}`$ scaling but now with a log correction arising from the MacDonald functions. For large $`k`$, below the spectral cut-off, the additional MacDonald functions act to heighten the flatness. That is, the introduction of a finite diffusivity actually increases the small-scale intermittency. In the next section we will investigate the correlators of all orders.

## VII Large-time behavior of higher correlators

The observation in the last section, that there is a dominant field component, allows us to predict that for large time $`|𝐁|^4|𝐁^2|^2`$ in each realization, that is $`Z_{\lambda \lambda }Z_{\alpha \beta }`$. Therefore, the property $`Z_{\lambda \lambda }=Z_{\alpha \beta }`$, if valid initially, should be preserved by the equation for $`Z`$. Indeed, let us differentiate the equation (30) twice with respect to $`\lambda `$ and subtract it from the same equation differentiated with respect to $`\alpha `$ and $`\beta `$. This gives the following closed equation for the combination $`w=Z_{\lambda \lambda }Z_{\alpha \beta }`$

$`\dot{w}`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Omega }}{3}}[k^2w_{kk}2k𝒟w_k+4𝒟^2w`$ (58)
$`+16(𝒟w+w)+3(_{\lambda \lambda }_{\alpha \beta })\{(\lambda ^24\alpha \beta )w\}]2\kappa k^2(𝒟w+2w).`$
We see that if $`w0`$ at $`t=0`$, then it will remain identically zero for all time. Thus, we can consider a class of solutions of (30), corresponding to large-time asymptotics of the general solution, such that $`Z_{\lambda \lambda }=Z_{\alpha \beta }`$. Assuming this equality in (30) and putting $`\eta =0`$ we have

$$\dot{Z}=\frac{\mathrm{\Omega }}{3}\left[k^2Z_{kk}+(2\lambda _\lambda +4)kZ_k+4(\lambda _\lambda )^2Z\right]2\kappa k^2\lambda _\lambda Z.$$
(59)
Let us consider a solution to this equation that is formally represented as a series in $`\lambda `$ (for example via a Taylor series):

$$Z=1+\underset{n=1}{\overset{\mathrm{}}{}}\frac{\lambda ^n}{n!}\mathrm{\Psi }^{(n)}(k,t).$$
(60)
Here, the function $`\mathrm{\Psi }^{(n)}`$ is the correlator of order $`2n`$. We omitted the lower index in $`\mathrm{\Psi }^{(n)}`$ because the correlators corresponding to different lower indices are identical in this case. Substituting (60) into (59) we have the following equation for these correlators

$$\dot{\mathrm{\Psi }}^{(n)}=\frac{\mathrm{\Omega }}{3}\left[k^2\mathrm{\Psi }_{kk}^{(n)}+(2n+4)k\mathrm{\Psi }_k^{(n)}+4n^2\mathrm{\Psi }^{(n)}\right]2\kappa k^2n\mathrm{\Psi }^{(n)}.$$
(61)
Note that for $`n=1`$ this equation agrees with the Kazantsev equation for the energy spectrum (32). Moreover, by the substitution

$$\mathrm{\Psi }^{(n)}=e^{(n1/2)(n+3/2)\mathrm{\Omega }t}k^{n3/2}\varphi (k/k_{d_n},t),$$
(62)
one can reduce (61) to an equation for $`\varphi `$, similar to (34), that is independent of $`n`$. We can therefore immediately write down the general solution for the correlators at any order $`n`$. In particular, in the perfect conductor regime $`1/\mathrm{\Omega }t(\mathrm{ln}q)^2`$ we have

$$\mathrm{\Psi }^{(n)}=\text{const}(n)t^{1/2}e^{(n1/2)(n+3/2)\mathrm{\Omega }t}e^{\frac{3(\mathrm{ln}k/q)^2}{4\mathrm{\Omega }t}}k^{n3/2},$$
(63)
and in the diffusive regime $`t(\mathrm{ln}q)^2`$

$$\mathrm{\Psi }^{(n)}=\text{const}(n)t^{3/2}e^{(n1/2)(n+3/2)\mathrm{\Omega }t}k^{n3/2}K_0(k/k_{d_n}).$$
(64)
We see that the main effect of the dissipation on all the moments (including the energy spectrum) is the prefactor changes from $`t^{1/2}t^{3/2}`$ (but without a change in the exponential growth rate), further the $`K_0(k/k_{d_n})`$ form-factor corresponds to a log-correction at $`k1`$ and an exponential cut-off at larger $`k`$.††††††Note that the higher the order of the correlator the earlier the spectral cut-off. It is this exponential cut-off that causes the exponential growth of the mean magnetic energy  and higher $`x`$-space moments of the magnetic field  to change. Indeed, simply integrating $`\mathrm{\Psi }^{(n)}`$ over $`k`$ with a cut-off at $`k=k_{d_n}`$ (and ignoring the prefactor and log corrections) one recovers the Kazantsev growth rate of the magnetic energy in the dissipative regime . Such an explanation was previously given in .

## VIII Physical interpretation of the scalings

Let us analyze the physical origins of the scalings obtained in the previous section because this can give us an indication to whether the same scalings should be expected for more general strain statistics.

Let us consider the expression (33) for the n-th order correlator and re-write it in the form $`C(n)exp((2n^2+2n)\mathrm{\Omega }t)k^nP(k,t)`$, where $`P(k,t)`$ is a universal function. The physical meaning of various terms in this expression can be easily analyzed. It follows from the Central Limit theorem that the large time statistics of the Lyapunov exponents is Gaussian with a dispersion $`D(t)\sqrt{t}`$. It follows from time reversal invariance that the average values of the Lyapunov exponents are $`\overline{\lambda },0,\overline{\lambda }`$. Hence, in the large time limit, $`B(t)exp(\lambda t)`$, where $`p.d.f.(t)exp(\frac{(\lambda \overline{\lambda })^2}{2\mathrm{\Delta }t^1})`$. Therefore, the statistics of the magnetic field is log-normal and $`B^{2n}exp((2\overline{\lambda }n^2+2\mathrm{\Delta }n)t)`$. For the statistics chosen in the present paper $`\overline{\lambda }=\mathrm{\Delta }=\mathrm{\Omega }`$. Substituting these values into the last expression for $`B^{2n}`$, we restore the correct $`n`$-dependence of the exponential growth rate of $`\mathrm{\Psi }^n`$. Note that terms of order $`n^2`$ in the growth rate are due to the Gaussian nature of the fluctuations of $`\lambda `$ around $`\overline{\lambda }`$, while the terms of order $`n`$ are due to the fact that $`\overline{\lambda }0`$. It is interesting that the magnetic field has log-normal statistics in both the perfect conductor and the dissipative regimes. For the coordinate space moments of $`B`$, persistence of the log-normality was emphasized in , although the $`e^{\lambda n^2t}`$ dependence was established earlier in . An equivalent result for the random matrices can be traced back to  (see also ).

The $`k`$-dependence of the magnetic field correlators is also very important because it gives us information about the dominant structures in wavenumber space. Suppose that initially the magnetic turbulence is isotropic and concentrated in a ball centered at the origin in wavenumber space. For each realization such a ball will stretch into an ellipsoid with one large, one short and one neutral dimension. One can visualize this ellipsoid as an elongated flat cactus leaf with thorns aligned to the direction of the magnetic field. Note that in this picture one component of the magnetic field (transverse to the cactus leaf) is dominant, this is captured by the fact that the polarization $`W`$ introduced in this paper tends to zero at large time. Another consequence of this picture is that the wavenumber space will be covered by the ellipsoids more and more sparsely at large $`k`$, this implies large intermittent fluctuations of the magnetic field Fourier transform. These fluctuations can be quantified by the flatness $`F`$ which was shown in (55) to grow as $`k^{3/2}`$, a clear indication of the small-scale intermittency.

To investigate the Kraichnan-Kazantsev model based dynamo problem further, a set of numerical experiments have been performed to, among other things, investigate the sensitivity of these analytical results to changes in the strain statistics. The details of this investigation can be found in . To help visualise the cactus leaf structures described above we will briefly include some snapshots of the configuration of a wavepacket ensemble in wavenumber space. As above each wavepacket is initially randomly distributed on a unit sphere in $`k`$-space with a randomly orientated magnetic field that lies in a plane perpendicular to the wavepacket’s wavevector, tangent to the sphere. Figure 1 shows the real magnetic fields of a set of $`500`$ wavepackets that have been subjected to two different realisations of the strain matrix.‡‡‡‡‡‡The imaginary magnetic fields are qualitatively similar and have therefore not been included. In the left-hand figure it is clear the magnetic field in this realisation is far from being plane polarised, the magnetic vectors still being predominantly random in their orientation at this given point in time. In contrast, in the right-hand figure, which has also been taken at the same point in time, we see the ellipsoid has become very elongated and the magnetic field appears plane polarised. In contrast, figure 2 shows the same two strain realisations but with a finite diffusivity. The right-hand figure in particular demonstrates why, in the diffusive regime, an ellipsoid will cover $`k`$-space more sparsely due to the decay of its magnetic field at its tips. This is the reason why spectral flatness increases in the diffusive regime.

## IX Conclusion

In this paper we introduced a description of the small-scale magnetic turbulence in terms of the one-point correlators of the Fourier amplitudes. From the classic works of  to more recent analysis (see for example  and references therein), the only correlator of this kind that was studied is the energy spectrum. The higher order correlators have only been considered in coordinate space . In this paper we have considered the Fourier correlators of all orders and showed that they contain important additional information about the turbulence, which is unavailable from a similar analysis of one-point or two-point coordinate space correlators. In particular, the fourth order Fourier correlators carry information about the mean polarization of the magnetic field modes. We showed that this polarization becomes plane for the Kraichnan-Kazantsev model. The $`k`$-scaling of the higher correlators allows us to determine the structures in Fourier space responsible for the intermittency, which for the Kraichnan-Kazantsev dynamo turns out to be elongated ellipsoids centered at the origin. The time scaling of the higher correlators allows one to conclude that the magnetic field has a log-normal statistics, although the same information is contained, and was established before, in analysis of the coordinate space correlators .

Finally, we would like to discuss an interesting connection between this work and the recent work of Schekochihin et al. . In particular the connection between the statistics of the magnetic field curvature studied by Schekochihin et al. and the magnetic polarization measure introduced in this present paper. Schekochihin et al. found that that the curvature of the magnetic field decreases, corresponding to folded and strongly stretched structures. This agrees with our results that the Fourier modes of the magnetic field tend to a state of plane polarization. However, the polarization gives more information than the curvature statistics. Indeed, zero curvature allows any structure that is constant along the magnetic field, in particular a set of magnetic filaments parallel to each other or a set of layered sheets, such that the magnetic field is constant on each sheet but its direction may change arbitrarily when passing from one layer to another. On the other hand, the wavenumber scalings obtained in this paper, (63) and (64), indicate that the magnetic field structures are layers in coordinate space and this rules out any filamentary structures. Further, our result about the plane polarization inhibit any “twists” of the magnetic field between layers, i.e. the magnetic field direction stays the same (or reverses) when passing from one layer to another. In fact, the presence of one neutral direction in the Lagrangian deformations tells us that these layers have a finite width in one direction and thus look like ribbons with the magnetic field directed along these ribbons.

##

Our aim here is to derive a closed equation for the generating function $`Z`$ (27), starting with equation (21). The last term in this equation is the easiest one

$$2\kappa k^2(\lambda |𝐁(𝐤)|^2+\alpha 𝐁^2(𝐤)+\beta \overline{𝐁}^2(𝐤))E=2\kappa k^2𝒟Z,$$
(65)
where $`𝒟`$ is a differential operator defined in (25). The correlators containing a factor of $`\sigma _{ij}`$ can be found using Gaussian integration by parts. In particular

$`\sigma _{ij}E=\mathrm{\Omega }{\displaystyle \frac{\delta E}{\delta \sigma _{ij}}}={\displaystyle \frac{\mathrm{\Omega }}{2}}`$ $`[\lambda (\mathrm{\Gamma }_{m,ij}\overline{B}_m+\overline{\mathrm{\Gamma }}_{m,ij}B_m)E`$ (67)
$`+2\alpha \mathrm{\Gamma }_{m,ij}B_mE+2\beta \overline{\mathrm{\Gamma }}_{m,ij}\overline{B}_mE],`$
where we have used the definition (22). Here, $`\mathrm{\Gamma }_{m,ij}`$ is a response function

$$\mathrm{\Gamma }_{m,ij}=\frac{\delta B_m}{\delta \sigma _{ij}}.$$
(68)
Differentiating (13) with respect to $`\sigma _{ij}`$ and using the statistical whiteness of the strain tensor we get

$$\mathrm{\Gamma }_{m,ij}=\left[k_i_j\frac{\delta _{ij}}{d}(1+k_l_l)\right]B_m+\delta _{mi}B_j.$$
(69)
In what follows we will make use of the isotropy of the turbulence, in particular, expressions of the type

$`(\overline{B}_iB_j+\overline{B}_jB_i)E`$ $`=`$ $`{\displaystyle \frac{2}{d1}}|B|^2E\left(\delta _{ij}{\displaystyle \frac{k_ik_j}{k^2}}\right),`$ (70)
$`B_iB_jE`$ $`=`$ $`{\displaystyle \frac{1}{d1}}B^2E\left(\delta _{ij}{\displaystyle \frac{k_ik_j}{k^2}}\right),`$ (71)
$`\overline{B}_i\overline{B}_jE`$ $`=`$ $`{\displaystyle \frac{1}{d1}}\overline{B}^2E\left(\delta _{ij}{\displaystyle \frac{k_ik_j}{k^2}}\right).`$ (72)
Substituting (69) into (67) and using the above isotropy relations we have

$$\sigma _{ij}E=\frac{\mathrm{\Omega }}{2}\left(k_i_j\frac{\delta _{ij}}{d}k_l_l\right)Z\frac{\mathrm{\Omega }}{d1}\left(\frac{k_ik_j}{k^2}\frac{\delta _{ij}}{d}\right)𝒟Z,$$
(73)
where $`𝒟`$ is a differential operator defined in (25). This allows us to find the first term on the right hand-side of (21)

$$k_i_j\sigma _{ij}E=\frac{\mathrm{\Omega }}{2}\left[\frac{d1}{d}k^2Z_{kk}+\frac{1}{d}\left(2𝒟+d^21\right)kZ_k2𝒟Z\right].$$
(74)
Similarly, the other three terms on the right hand-side of (21) can be obtained via Gaussian integration by parts, and the use of the response function (69) and isotropy condition. After some lengthy but straightforward algebra one obtains

$`\lambda \sigma _{ml}(\overline{B}_mB_l+\overline{B}_lB_m)E=\lambda \mathrm{\Omega }`$ $`[(d{\displaystyle \frac{2}{d}})Z_\lambda +2(1{\displaystyle \frac{1}{d}})𝒟Z_\lambda `$ (76)
$`+\lambda (Z_{\alpha \beta }Z_{\lambda \lambda }){\displaystyle \frac{1}{d}}k_i_iZ_\lambda ],`$
and

$`2\alpha \sigma _{ml}B_mB_lE=\alpha \mathrm{\Omega }`$ $`[(d{\displaystyle \frac{2}{d}})Z_\alpha +2(1{\displaystyle \frac{1}{d}})𝒟Z_\alpha `$ (78)
$`2\beta Z_{\alpha \beta }+2\beta Z_{\lambda \lambda }{\displaystyle \frac{1}{d}}k_i_iZ_\alpha ].`$
The fourth term can be obtained from (78) via interchanging $`\alpha `$ with $`\beta `$ and $`B`$ with $`\overline{B}`$

$`2\beta \sigma _{ml}\overline{B}_m\overline{B}_lE=\beta \mathrm{\Omega }`$ $`[(d{\displaystyle \frac{2}{d}})Z_\beta +2(1{\displaystyle \frac{1}{d}})𝒟Z_\beta `$ (80)
$`2\alpha Z_{\alpha \beta }+2\alpha Z_{\lambda \lambda }{\displaystyle \frac{1}{d}}k_i_iZ_\beta ].`$
Using the expressions (74), (76), (78), (80) and (65), we obtain the final equation

$`\dot{Z}`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Omega }}{2}}[(1{\displaystyle \frac{1}{d}})k^2Z_{kk}+{\displaystyle \frac{1}{d}}(4𝒟+d^21)kZ_k`$ (82)
$`+(2d6)𝒟Z+4(1{\displaystyle \frac{1}{d}})𝒟^2Z+2(\lambda ^24\alpha \beta )(Z_{\alpha \beta }Z_{\lambda \lambda })]2\kappa k^2𝒟Z.`$