text:
# Particle Physics in Intense Electromagnetic Fields 11footnote 1Published in: IL Nuovo Cimento A, 1999, vol.112A, p.977-1000.

## 1 Introduction

Quantum field theory in the presence of an intense electromagnetic background has a long fruitful history. The interest promoting investigations in this area is based on several reasons. The first one being purely theoretical is associated with the self-consistency problem of the relativistic quantum theory at very high energies. It is well known that in quantum electrodynamics (QED) perturbative vacuum probably is not a true ground state of the theory because of the fictitious pole in photon propagator in the ultraviolet region : $`|k^2|m_e^2\mathrm{exp}(3\pi /\alpha )`$ . This challenge has inspired many studies dealing with QED of intense fields in which the convergence of perturbation series in fine-structure constant $`\alpha =e^2/4\pi `$ has been analyzed (see e.g. , ,  and references therein). The matter is that ultraviolet behaviour of the theory and the structure of divergencies determine in a uniform way the properties of Green functions in very intense electromagnetic fields and at large momentum transfer. In other words the physics of extremely strong background fields is intimately connected with the one at small distances ,. This can be illustrated by the effect of vacuum polarization in QED which modifies the Lagrange density of Maxwell electrodynamics. For weak fields one loop corrections have been calculated long ago by Heisenberg and Euler .

$$=+\frac{2\alpha ^2}{45m_e^4}\left(4^2+7𝒢^2\right),$$
(1)
where

$$=\frac{1}{4}F_{\mu \nu }F^{\mu \nu }=\frac{1}{2}\left(𝐁^2𝐄^2\right),𝒢=\frac{1}{4}F_{\mu \nu }\stackrel{~}{F}^{\mu \nu }=\mathrm{𝐄𝐁}.$$
(2)
However the above expression is valid only if the field strength is relatively small ($`,𝒢m_e^4/e^2`$). The domain of extremely strong fields ($`m_e^4/e^2;𝒢`$) is described by another formula  which in logarithmic approximation can be written as

$$=+\frac{\alpha }{6\pi }\mathrm{ln}\left(\frac{2e^2}{m_e^4}\right)$$
(3)
We see that asymptotically the effective Lagrangian (3) vanishes as far as the field strength $`B=\sqrt{2}`$ tends to the critical value $`m_e^2/e\mathrm{exp}(3\pi /\alpha )`$. So, the Moscow zero-charge problem  reveals itself as a nought of effective Lagrangian if we impose the following interrelation:

$$|k^2|eB=e\sqrt{2}$$
(4)
Two loop calculations  have confirmed that the analogy (4) is retained in all orders of perturbation. From this viewpoint the physics of intense electromagnetic fields is a powerful theoretical tool for searching the stable non-perturbative vacuum state in which the artifacts mentioned above can be removed.

Let us also note, that the problem of vacuum instability and phase transitions in a strong background field has received additional interest in the context of non-Abelian gauge theories . These investigations rely essentially on the philosophy and technique elaborated in the framework of QED with an intense background.

The second motivation for studying particles interactions with very strong electromagnetic fields originates from the quantum theory of synchrotron radiation . The process of bremsstrahlung by electrons in a homogeneous magnetic field is has long been known. Classical electrodynamics provides a satisfactory description of the radiation power and spectrum being observed in synchrotrons (for a recent review see ). At the same time particular phenomena require the apparatus of QED with an intense magnetic field to be applied for an adequate theoretical explanation. Systematic investigations of quantum processes in the presence of strong magnetic fields have been initiated in works of Sokolov and Ternov  and were continued by their collaborators . These investigations resulted in discovering many interesting theoretical predictions some of which have been confirmed experimentally. Here we can highlight the famous effect of self-polarization of electrons in storage rings due to the synchrotron radiation.

The third reason stimulating interest in physics with intense electromagnetic background is caused by recent CERN experiments involving SPS accelerator in which the particle interactions with single crystals were explored , . Notice that, because of regular and systematic structure, electromagnetic interactions in single crystals are substantially enhanced in comparison with the amorphous medium. This gives a unique possibility of creating a strong electromagnetic background for relativistic particles penetrating single crystal near axial or planar directions. Thus we have an opportunity for searching new physics in such unusual conditions because quantum phenomena in the presence of the intense electromagnetic field differ considerably from those occurring in a vacuum. Due to non-linear and non-perturbative influence of the external field one can observe absolutely new processes which as a rule are forbidden under normal conditions. The electromagnetic field open new channels of reactions taking away the bans for transitions between definite quantum states. The most wonderful feature is the possibility for very light particles to decay into heavy species capturing a lacking amount of energy from the intense electromagnetic environment. It would be pertinent to mention here the process of electron-positron pairs production by solitary photons incident along crystal axes. The data of experiments CERN WA-81, CERN NA-046  concerning with this reaction $`\gamma e^+e^{}`$ are in a satisfactory agreement with predictions of QED with the strong electromagnetic field  (see also ).

Thus, we see that particle physics in the presence of non-perturbative background fields is an impetuously developing branch of physical sciences having many interesting applications. Unfortunately it is impossible to itemize all exploitations of this theoretical technique because they are very numerous. Here I can only mention some other important employments such as the problem of vacuum in quantum chromodynamics , astrophysics of neutron stars , the cosmology of the early Universe  and etc.

Lastly, let us say some words about the method itself. The philosophy of calculations in strong electromagnetic fields originates from the Feynman non-perturbative approach in QED  which was elaborated in details by Furry . This technique has become widely known as the “Furry picture”. Classical examples of QED with external fields were given in works of Schwinger , . Further developments of the Feynman-Furry method were achieved in the research of Fradkin and Gitman with their collaborators who have generalized and modified the above approach for a wide class of electromagnetic fields with unstable vacuum (see  and references therein). It is worth noting that the procedure of second quantization in the intense background field is a separate, uneasy task having common features with the quantum field theory in a curved space-time . A particular difficulty is connected with appropriate treatment of the classical solutions of relativistic wave equations that do not admit one particle interpretation. For instance, in a case of a constant uniform electric field exact solutions of the Dirac equation have a form which cannot be splitted into positive and negative frequency modes corresponding to particle and antiparticle wave functions. The later circumstance is regarded as the instability of the physical vacuum with electric fields in relation with the processes of pairs production. However, as was shown by Schwinger  the probability of pairs creation from the vacuum becomes feasible only in very strong electric fields being comparable with the critical value $`E_{cr}=m_e^2/e=1.3210^{18}V/m`$.

Inasmuch as present values of the field strength being available in experiments are much less than $`E_{cr}`$ it is advantageous to employ an effective approximation which was proposed in works of Ritus and Nikishov . The main idea comes from a detailed analysis of probabilities of quantum phenomena in an arbitrary constant homogeneous electromagnetic field. Suppose that the initial quantum state is characterized by the only one particle with momentum $`p^\mu `$. Then, in general the probability of quantum transition to some final state $`P_{fi}`$ depends on the background field strength $`F_{\mu \nu }`$ through the following Lorentz-invariant dimensionless parameters:

$$a=\frac{e^2𝒢}{m_e^4}=\frac{e^2}{m_e^4}\mathrm{𝐄𝐁},b=\frac{e^2}{m_e^4}=\frac{e^2}{2m_e^4}\left(𝐁^2𝐄^2\right).$$
(5)
$$\chi =\frac{e}{m_e^3}\sqrt{(F_{\mu \nu }p^\nu )^2}=\frac{e}{m_e^3}\left[(p_0𝐄+𝐩\times 𝐁)^2(\mathrm{𝐩𝐄})^2\right]^{1/2}$$
(6)
The domain of electromagnetic fields that currently can be handled in a laboratory satisfies the conditions: $`a1,b1`$ (i.e. $`EE_{cr},BB_{cr}=m_e^2/e=4.4110^9T`$). Now, if the initial particle is relativistic ($`p_0m_e`$) then $`\chi ^2a;b`$ and one can expand the probability $`P_{fi}`$ into a series of these small quantities:

$$P_{fi}(\chi ,a,b)=P(\chi ,0,0)+a\frac{P}{a}(\chi ,0,0)+b\frac{P}{b}(\chi ,0,0)+\mathrm{}$$
(7)
The first term in eq.(7) represents the probability of quantum transition in a so-called crossed field, which is a combination of the orthogonal electric and magnetic fields having the same magnitude ($`=𝒢=0`$). So, the formula describing a reaction in the crossed field provides a good approximation for that in an arbitrary constant electromagnetic field if one can omit in eq.(7) the corrections proportional to small parameters $`a,b`$. The merits of the crossed field approximation lie in the most simple form of charged-particle wave functions and propagators which greatly facilitate the scheme of calculations. All what said above may be summarized as the following prescription: if you want to investigate a process in a constant homogeneous electromagnetic field, then the most straightforward way is to calculate the corresponding probability assuming that $`=𝒢=0`$ (2). The obtained result depending on the sole background field parameter $`\chi `$ (6) will simulate the exact one as the first term in the expansion (7). Note, that the crossed field is a particular case of the plane-wave electromagnetic field fitting the choice: $`A_\mu (\phi )=a_\mu \phi ;F_{\mu \nu }=k_\mu a_\nu k_\nu a_\mu `$ (see below).

The purpose of this paper is to give a pedagogical introduction to the particle physics in the presence of an intense electromagnetic background. I briefly reproduce the main milestones of the Feynman-Furry method allowing for non-perturbative interactions with very strong electromagnetic fields. Taking into account the spectrum of the known electrically charged particles, I consider the cases of scalar, spinor and vector quantum fields separately and obtain various forms for wave functions and propagators that are suitable for practical calculations (sections 3,4,5). Then I discuss some general features that are typical for quantum phenomena in the background electromagnetic field (section 6). And finally, I analyze the possibilities of searching new physics in the framework of Minimal Standard Model through the investigation of quantum effects induced by the strong electromagnetic background.

## 2 Scalar Particles in an Intense Electromagnetic Background

Let us consider the Lagrangian of scalar particles with charge $`e>0`$ interacting with a strong electromagnetic field.

$$=(_\mu \varphi ^+)(^\mu \varphi ^{})m^2\varphi ^+\varphi ^{}+ieA^\mu (\varphi ^+_\mu \varphi ^{}\varphi ^{}_\mu \varphi ^+)+e^2\varphi ^+\varphi ^{}(A^\mu A_\mu )$$
(8)
The electromagnetic potential $`A_\mu (x)`$ is regarded here as a given function of space-time coordinates which corresponds to the classical field strength tensor $`F_{\mu \nu }=_\mu A_\nu _\nu A_\mu `$. We wish to describe the impact of the external electromagnetic field in a non-perturbative manner employing exact solutions of the equation of motion which follows from the Lagrangian (8).

$$\left[(_\mu \pm ieA_\mu )^2+m^2\right]\varphi ^\pm (x)=0$$
(9)
The greatest embarrassment of this approach is that of solving analytically equation (9) and to obtain non-perturbative wave functions for an arbitrary predetermined configuration of the background field. At present exact solutions are known for the case of a constant uniform electromagnetic field, for the Coulomb field potential and for the so-called Redmond’s configuration (see e.g. ). Besides there is an explicit expression for the scalar particle wave function in the presence of a plane electromagnetic wave, which is characterized by some unspecified potential $`A_\mu (\phi )`$ depending on space-time coordinates $`x^\mu `$ only through the phase $`\phi =k_\mu x^\mu =\omega t\mathrm{𝐤𝐱}`$. This is the case we intend to discuss in more details.

The starting point of our analysis is to invoke the WKB method for the equation (9) which as it comes out will give us the exact result in the first approximation. Now, substituting $`\varphi ^+(x)=e^{iS(x)}`$ and imposing the Lorentz gauge for the vector potential, we get

$$(_\mu S+eA_\mu )^2m^2=i_\mu ^\mu S.$$
(10)
In the semiclassical limit being considered the right-hand side of eq. (10) can be neglected and it gives us the Hamilton-Yacobi equation for the classical action of a relativistic particle moving in a plane-wave electromagnetic field. The classical action as a function of coordinates $`x^\mu `$ and the vector of initial momentum $`p^\mu `$ is well known and it can be written in the following form:

$$S(x,p)=p_\mu x^\mu \underset{0}{\overset{kx}{}}\left(\frac{ep_\mu A^\mu (\phi )}{p_\nu k^\nu }\frac{e^2A_\mu (\phi )A^\mu (\phi )}{2p_\nu k^\nu }\right)𝑑\phi +S_0.$$
(11)
Here it is implied that the wave vector $`k^\mu `$ obeys the isotropic and transversal conditions

$$k^2=0,k^\mu A_\mu =0.$$
(12)
In order to obtain quantum corrections one must insert formula (11) into eq.(10) and take into account the right-hand side of this expression. However it is easy to see that classical action (11) makes both sides of eq. (10) vanish if the momentum vector $`p^\mu `$ is taken on the mass shell: $`p_0=\sqrt{𝐩^2+m^2}`$. So, we have found an exact classical solution of eq.(9) which should be treated as the scalar particle wave function in a background electromagnetic wave field. For practical use the wave function must be normalized in a proper way. Assuming that the configuration volume space is equal to $`V`$, we finally get

$$\varphi ^\pm (x,𝐩)=\frac{1}{\sqrt{2p_0V}}\mathrm{exp}\left[ipxi\underset{0}{\overset{kx}{}}\left(\pm \frac{epA}{pk}\frac{e^2A^2}{2pk}\right)𝑑\phi \right].$$
(13)
The next step is to realize the second quantization procedure. As usually, we decompose the scalar field into a series of the obtained-above solutions involving the Fock creation and annihilation operators.

$$\varphi ^+(x)=\underset{𝐩}{}\left[a_𝐩\varphi ^+(x,𝐩)+b_𝐩^{}\varphi ^{}(x,𝐩)\right]$$
(14)
$$\varphi ^{}(x)=\underset{𝐩}{}\left[a_𝐩^{}\varphi ^+(x,𝐩)+b_𝐩\varphi ^{}(x,𝐩)\right]$$
(15)
These operators must satisfy a standard set of commutation relations

$`[a_𝐩;a_𝐩^{}^{}]=\delta _{\mathrm{𝐩𝐩}^{}},[b_𝐩;b_𝐩^{}^{}]=\delta _{\mathrm{𝐩𝐩}^{}},`$
$`[a_𝐩;a_𝐩^{}]=[b_𝐩;b_𝐩^{}]=[a_𝐩;b_𝐩^{}]=[a_𝐩;b_𝐩^{}^{}]=0,`$ (16)
which enable us to reproduce the ordinary scheme of quantum field theory. All further reasoning almost exactly repeats the ones of the conventional approach. The negative frequency modes of eq.(9) are regarded as the scalar antiparticles with wave functions $`\varphi ^{}(x,𝐩)`$ and the opposite charge $`e<0`$. The physical vacuum is defined as a ground quantum state with no particles present: $`a_𝐩|0>=b_𝐩|0>=0`$, and the Feynman propagator is determined by the formula

$$D(x,x^{})=i<0|T\varphi ^+(x)\varphi ^{}(x^{})|0>$$
(17)
So, the essence of the Furry picture of interaction with a background electromagnetic field is to replace the ordinary de Broglie’s wave functions by some exact solutions which incorporate adequately all non-perturbative effects. As a calculation prescription, we present explicitly the modification of the scalar particle propagator (17) in a plane-wave field:

$$D(x,x^{})=\frac{d^4p}{(2\pi )^4}\frac{1}{m^2p^2i0}\mathrm{exp}\left[ip(xx^{})i\underset{kx^{}}{\overset{kx}{}}\left(\frac{epA}{pk}\frac{e^2A^2}{2pk}\right)𝑑\phi \right].$$
(18)
Notice that the propagator (18) is a Green function of the Klein-Gordon equation (9) id est it obeys inhomogeneous equation with the Dirac $`\delta `$-function in the right-hand side:

$$\left[(_\mu +ieA_\mu )^2+m^2\right]D(x,x^{})=\delta ^4(xx^{})$$
(19)
Another useful representation for the scalar propagator can be obtained from eq.(19) through the Fock-Schwinger proper time formalism , . A comprehensive expounding of this method can be found in many text books . For the sake of completeness let us introduce here the final expression for the propagator (17) in such an approach which is equivalent to eq.(18):

$`D(x,x^{})={\displaystyle \frac{\mathrm{\Phi }(x,x^{})}{(4\pi )^2}}{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{ds}{s^2}}\mathrm{exp}[ism^2i{\displaystyle \frac{(xx^{})^2}{4s}}+{\displaystyle \frac{ie^2s}{(kxkx^{})}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}A_\mu (\phi )A^\mu (\phi )d\phi `$
$`{\displaystyle \frac{ie^2s}{(kxkx^{})^2}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}d\phi _1{\displaystyle \underset{kx^{}}{\overset{kx}{}}}d\phi _2A_\mu (\phi _1)A^\mu (\phi _2)],`$ (20)
where the gauge factor $`\mathrm{\Phi }(x,x^{})`$ is determined by the integral along the straight line from the point $`x_\mu `$ to the point $`x_\mu ^{}`$:

$$\mathrm{\Phi }(x,x^{})=\mathrm{exp}\left[ie\underset{x^{}}{\overset{x}{}}A_\mu (z)𝑑z^\mu \right]=\mathrm{exp}\left[ie\frac{(xx^{})^\mu }{(kxkx^{})}\underset{kx^{}}{\overset{kx}{}}A_\mu (\phi )𝑑\phi \right].$$
(21)
## 3 Charged Fermions in an Intense Electromagnetic Field

The results of previous section can be easily generalized for the spin $`1/2`$ particles being governed by QED-like Lagrangian

$$=\frac{i}{2}\overline{\psi }\gamma ^\mu (_\mu \psi )\frac{i}{2}(_\mu \overline{\psi })\gamma ^\mu \psi m(\overline{\psi }\psi )eA_\mu (\overline{\psi }\gamma ^\mu \psi ),$$
(22)
which gives the Dirac equation for the four component bispinor $`\psi `$:

$$\left(i\gamma ^\mu _\mu e\gamma ^\mu A_\mu m\right)\psi (x)=0$$
(23)
Exact solutions in a plane electromagnetic wave have been obtained by Volkov  and they are deduced from the Klein-Gordon equation (9) by a standard substitution (see also ):

$$\psi (x)=\left(i\gamma ^\mu _\mu e\gamma ^\mu A_\mu +m\right)\phi (x).$$
(24)
The squared Dirac equation

$$\left[(_\mu +ieA_\mu )^2+m^2+\frac{e}{2}\sigma ^{\mu \nu }F_{\mu \nu }\right]\phi (x)=0$$
(25)
differs from that for scalar particles only in the last term arising because of $`\gamma `$-matrices commutator

$$\sigma ^{\mu \nu }=\frac{i}{2}(\gamma ^\mu \gamma ^\nu \gamma ^\nu \gamma ^\mu ).$$
(26)
Note that in a plane electromagnetic wave all powers of this term exceeding the first are equal to zero due to eq.(12). So the desired solution of equation (25) is a slightly modified form of expression (13).

$$\phi (x)=\mathrm{exp}\left[ipxi\underset{0}{\overset{kx}{}}\left(\frac{epA}{pk}\frac{e^2A^2}{2pk}\right)𝑑\phi \right]\left[1+\frac{e}{2pk}\widehat{k}\widehat{A}(kx)\right]v,$$
(27)
where

$$\widehat{k}=k_\mu \gamma ^\mu ,\widehat{A}=\gamma ^\mu A_\mu ,$$

and $`v`$ is an arbitrary constant bispinor. Now, returning to eq.(24), we derive the sought for non-perturbative wave function of a charged fermion in an intense electromagnetic background:

$$\mathrm{\Psi }^\pm (x,𝐩)=\mathrm{exp}\left[ipxi\underset{0}{\overset{kx}{}}\left(\pm \frac{epA}{pk}\frac{e^2A^2}{2pk}\right)𝑑\phi \right]\left[1\pm \frac{e}{2pk}\widehat{k}\widehat{A}(kx)\right]\frac{u_\sigma (p)}{\sqrt{2p_0V}}$$
(28)
Polarization effects are described here by an ordinary Dirac bispinor

$$u_\sigma (p)=\sqrt{2p_0V}(\gamma ^\mu p_\mu +m)v$$

which is subject to the following restrictions:

$$(\gamma ^\mu p_\mu m)u_\sigma (p)=0\overline{u}_\sigma (p)u_\sigma (p)=2m\overline{u}_\sigma (p)\gamma ^\mu u_\sigma (p)=2p^\mu $$
(29)
The second quantization is carried out by analogy with eq. (14),(15)

$$\psi (x)=\underset{𝐩\sigma }{}\left[a_{𝐩\sigma }\mathrm{\Psi }^+(x,𝐩)+b_{𝐩\sigma }^{}\mathrm{\Psi }^c(x,𝐩)\right]$$
(30)
$$\overline{\psi }(x)=\underset{𝐩\sigma }{}\left[a_{𝐩\sigma }^{}\overline{\mathrm{\Psi }^+}(x,𝐩)+b_{𝐩\sigma }\overline{\mathrm{\Psi }^c}(x,𝐩)\right]$$
(31)
and involves also the antifermion wave functions being obtained through the $`C`$-matrix of charge conjugation:

$`\mathrm{\Psi }^c(x,𝐩)=C\overline{\mathrm{\Psi }}^T(x,𝐩)=\gamma ^0C\mathrm{\Psi }^{}(x,𝐩)`$
$`\overline{\mathrm{\Psi }^c}(x,𝐩)=\mathrm{\Psi }^T(x,𝐩)C^{}`$ (32)
A specific form of this matrix is inessential provided that it satisfies the following properties:

$$C^{}=C^1C^T=CC^1\gamma ^\mu C=\gamma ^{\mu T}$$
(33)
For example we can choose $`C=i\gamma ^0\gamma ^2`$ in Weyl and Dirac $`\gamma `$-matrices representations.

Taking into account the spin-statistics relation, we must impose the anticommutation rules for the Fock operators $`a_{𝐩\sigma },b_{𝐩\sigma }`$ which have no difference with the conventional ones:

$`\{a_{𝐩\sigma };a_{𝐩^{}\sigma ^{}}^{}\}=\delta _{\mathrm{𝐩𝐩}^{}}\delta _{\sigma \sigma ^{}}\{b_{𝐩\sigma };b_{𝐩^{}\sigma ^{}}^{}\}=\delta _{\mathrm{𝐩𝐩}^{}}\delta _{\sigma \sigma ^{}}`$
$`\{a_{𝐩\sigma };a_{𝐩^{}\sigma ^{}}\}=\{b_{𝐩\sigma };b_{𝐩^{}\sigma ^{}}\}=\{a_{𝐩\sigma };b_{𝐩^{}\sigma ^{}}\}=\{a_{𝐩\sigma };b_{𝐩^{}\sigma ^{}}^{}\}=0`$ (34)
The corresponding Feynman propagator is defined as usual

$$S(x,x^{})=i<0|T\psi (x)\overline{\psi }(x^{})|0>$$
(35)
and has the most simple form in the basis of exact solutions (28):

$`S(x,x^{})={\displaystyle }{\displaystyle \frac{d^4p}{(2\pi )^4}}[1+{\displaystyle \frac{e\widehat{k}}{2pk}}\widehat{A}(kx)]{\displaystyle \frac{\gamma ^\mu p_\mu +m}{p^2m^2+i0}}[1{\displaystyle \frac{e\widehat{k}}{2pk}}\widehat{A}(kx^{})]`$
$`\mathrm{exp}\left[ip(xx^{})i{\displaystyle \underset{kx^{}}{\overset{kx}{}}}\left({\displaystyle \frac{epA}{pk}}{\displaystyle \frac{e^2A^2}{2pk}}\right)𝑑\phi \right].`$ (36)
It is often more convenient to use the Fock-Schwinger representation which enables us to rewrite the last formula as follows.

$`S(x,x^{})={\displaystyle \frac{\mathrm{\Phi }(x,x^{})}{(4\pi )^2}}{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{ds}{s^2}}\mathrm{exp}[ism^2i{\displaystyle \frac{(xx^{})^2}{4s}}+{\displaystyle \frac{ie^2s}{(kxkx^{})}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}A_\mu (\phi )A^\mu (\phi )d\phi `$
$`{\displaystyle \frac{ie^2s}{(kxkx^{})^2}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}d\phi _1{\displaystyle \underset{kx^{}}{\overset{kx}{}}}d\phi _2A_\mu (\phi _1)A^\mu (\phi _2)\left]\right\{m+{\displaystyle \frac{(\widehat{x}\widehat{x}^{})}{2s}}+mes\widehat{k}{\displaystyle \frac{\widehat{A}(kx)\widehat{A}(kx^{})}{(kxkx^{})}}+`$
$`+{\displaystyle \frac{e}{2(kxkx^{})}}\left[\widehat{k}\widehat{A}(kx)(\widehat{x}\widehat{x}^{})(\widehat{x}\widehat{x}^{})\widehat{k}\widehat{A}(kx^{})\right]+{\displaystyle \frac{e}{(kxkx^{})}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}\widehat{A}(\phi )𝑑\phi `$
$`e\widehat{k}{\displaystyle \frac{(xx)^\mu }{(kxkx^{})^2}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}A_\mu (\phi )𝑑\phi {\displaystyle \frac{e^2s\widehat{k}}{(kxkx^{})}}\widehat{A}(kx)\widehat{A}(kx^{})+{\displaystyle \frac{e^2s\widehat{k}}{(kxkx^{})^2}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}A_\mu (\phi )A^\mu (\phi )𝑑\phi +`$
$`+{\displaystyle \frac{e^2s\widehat{k}}{(kxkx^{})^2}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}[\widehat{A}(kx)\widehat{A}(\phi )+\widehat{A}(\phi )\widehat{A}(kx^{})]d\phi {\displaystyle \frac{2e^2s\widehat{k}}{(kxkx^{})^3}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}d\phi _1{\displaystyle \underset{kx^{}}{\overset{kx}{}}}d\phi _2A_\mu (\phi _1)A^\mu (\phi _2)\}`$ (37)
One can easily check that fermion propagator (3), (3) is a Green function of the Dirac equation (23) in the presence of the intense electromagnetic field.

$$[i\gamma ^\mu _\mu e\gamma ^\mu A_\mu m]S(x,x^{})=\delta ^4(xx^{}).$$
(38)
## 4 W-bosons in a Background Electromagnetic Field

Let us now concentrate on electromagnetic interactions of weak intermediate vector bosons $`W^\pm `$. The corresponding part of the Lagrangian of the Glashow-Weinberg-Salam model has the form (see e.g. )

$`={\displaystyle \frac{1}{2}}(_\mu W_\nu ^+_\nu W_\mu ^+)(^\mu W^\nu ^\nu W^\mu )+m_W^2W_\mu ^+W^\mu ieF^{\mu \nu }W_\mu ^+W_\nu ^{}`$
$`ie(_\mu W_\nu ^+_\nu W_\mu ^+)W^\mu A^\nu +ie(_\mu W_\nu ^{}_\nu W_\mu ^{})W^{+\mu }A^\nu e^2(A_\mu A^\mu )(W_\nu ^+W^\nu )+`$
$`+e^2(A^\mu W_\mu ^+)(A^\nu W_\nu ^{}){\displaystyle \frac{1}{\xi }}(^\mu W_\mu ^++ieA^\mu W_\mu ^+)(^\nu W_\nu ^{}ieA^\nu W_\nu ^{})`$ (39)
The last summand in (4) arises from a non-linear gauge fixing term

$`_{gf}={\displaystyle \frac{1}{2\xi }}(^\mu A_\mu ^1gA_\mu ^2A^{3\mu }+\xi m_W\phi _1)^2{\displaystyle \frac{1}{2\xi }}(^\mu A_\mu ^2+gA_\mu ^1A^{3\mu }+\xi m_W\phi _2)^2`$
$`{\displaystyle \frac{1}{2\xi }}(^\mu A_\mu ^3+\xi m_Z\mathrm{cos}\theta _W\phi _3)^2{\displaystyle \frac{1}{2\xi }}(^\mu A_\mu ^0\xi m_Z\mathrm{sin}\theta _W\phi _3)^2`$ (40)
which generalizes the linear $`R_\xi `$-gauge and is more suitable for calculations in a background electromagnetic field . $`A_\mu ^a`$ and $`A_\mu ^0`$ denote here the $`SU(2)`$ and $`U(1)`$ gauge fields respectively whereas $`\phi _a`$ designates the Goldstone bosons. The Lagrangian (4) brings the equation of motion

$$\left[(_\nu \pm ieA_\nu )^2+m_W^2\right]W^{\pm \mu }+(1/\xi 1)(^\mu \pm ieA^\mu )(^\nu W_\nu ^\pm \pm ieA^\nu W_\nu ^\pm )\pm 2ieF^{\mu \nu }W_\nu ^\pm =0$$
(41)
describing charged massive particles with spin $`0`$ and $`1`$. The longitudinal mode of the vector field $`W^\pm `$

$$B^\pm (x)=\frac{1}{\xi m_W}(^\mu W_\mu ^\pm \pm ieA^\mu W_\mu ^\pm )$$
(42)
corresponds to the unphysical scalar with mass $`\xi m_W`$ and it obeys the Klein-Gordon equation

$$\left[(_\nu \pm ieA_\nu )^2+\xi m_W^2\right]B^\pm (x)=0.$$
(43)
The wave functions of $`W`$-bosons come from the transversal mode being associated with vector quantum states.

$$V_\mu ^\pm (x)=W_\mu ^\pm (x)+\frac{1}{\xi m_W^2}(_\mu \pm ieA_\mu )(^\nu W_\nu ^\pm \pm ieA^\nu W_\nu ^\pm )$$
(44)
The equation for the transversal mode $`V_\mu ^\pm `$ is derived from eq.(41) and it can be written in the following way.

$$\left[(_\nu \pm ieA_\nu )^2+m_W^2\right]V^{\pm \mu }(x)\pm 2ieF^{\mu \nu }V_\nu ^\pm (x)=0$$
(45)
Besides, due to eq.(43) the field $`V_\mu ^\pm `$ satisfies additional condition

$$(^\mu \pm ieA^\mu )V_\mu ^\pm (x)=0$$
(46)
which means that only $`3`$ components of the field $`V_\mu ^\pm `$ are independent (exactly as it should be for spin-1 particles). Thus, we have obtained a Lorentz-covariant decomposition of the field $`W_\mu ^\pm `$ into the scalar and vector parts:

$$W_\mu ^\pm (x)=V_\mu ^\pm (x)+\frac{1}{m_W}(_\mu \pm ieA_\mu )B^\pm (x)$$
(47)
Notice that all physically sensible information is contained in the vector field $`V_\mu ^\pm `$, whereas the scalar field $`B^\pm `$ plays an auxiliary role and its contribution to the S-matrix is canceled by that of the Goldstone bosons $`\phi ^\pm =\frac{1}{\sqrt{2}}(\phi _2\pm i\phi _1)`$.

To carry out the scheme of canonical quantization, one needs an explicit form of the $`W`$-boson wave functions in an intense electromagnetic background. This requires obtaining exact solutions of eq.(45) for a given configuration of the external field. If we restrict ourselves to the case of a plane-wave electromagnetic field with potential $`A_\mu (\phi )`$ then eq.(45) can be easily integrated and we get:

$`W_\mu ^\pm (x,𝐩)=\mathrm{exp}[ipxi{\displaystyle \underset{0}{\overset{kx}{}}}(\pm {\displaystyle \frac{epA}{pk}}{\displaystyle \frac{e^2A^2}{2pk}})d\phi ]`$
$`\left[g_{\mu \nu }\pm {\displaystyle \frac{e}{pk}}\left(k_\mu A_\nu (kx)k_\nu A_\mu (kx)\right){\displaystyle \frac{e^2}{2(pk)^2}}A^2(kx)k_\mu k_\nu \right]{\displaystyle \frac{v^{\pm \nu }(p,\sigma )}{\sqrt{2p_0V}}},`$ (48)
where $`v^{\pm \nu }(p,\sigma )`$ is a complex polarization vector which is subject to the following constraints:

$$v_\mu ^\pm p^\mu =0,v_\mu ^\pm v^{\pm \mu }=1,\underset{\sigma }{}v_\mu ^\pm (p,\sigma )v_\nu ^\pm (p,\sigma )=g_{\mu \nu }+p_\mu p_\nu /m_W^2.$$
(49)
Now the task is to represent the field $`W_\mu ^\pm `$ as the sum of one-particle quantum excitations. Using expression (47) we expand $`W_\mu ^\pm `$ into series of classical solutions (4), (13)

$$W_\mu ^+(x)=\underset{𝐩\sigma }{}\left[a_{𝐩\sigma }W_\mu ^+(x,𝐩)+b_{𝐩\sigma }^{}W_\mu ^{}(x,𝐩)\right]+\frac{(_\mu +ieA_\mu )}{m_W}\underset{𝐩}{}\left[a_{𝐩0}\varphi ^+(x,𝐩)+b_{𝐩0}^{}\varphi ^{}(x,𝐩)\right],$$
(50)
$$W_\mu ^{}(x)=\underset{𝐩\sigma }{}\left[a_{𝐩\sigma }^{}W_\mu ^+(x,𝐩)+b_{𝐩\sigma }W_\mu ^{}(x,𝐩)\right]+\frac{(_\mu ieA_\mu )}{m_W}\underset{𝐩}{}\left[a_{𝐩0}^{}\varphi ^+(x,𝐩)+b_{𝐩0}\varphi ^{}(x,𝐩)\right],$$
(51)
which incorporate the Fock operators of two kinds. The observables $`a_{𝐩\sigma },b_{𝐩\sigma }`$ correspond to the vector quantum states of $`W`$-bosons and satisfy the usual commutation relations (with a polarization subscript $`\sigma =1,2,3`$.)

$`[a_{𝐩\sigma };a_{𝐩^{}\sigma ^{}}^{}]=\delta _{\mathrm{𝐩𝐩}^{}}\delta _{\sigma \sigma ^{}}[b_{𝐩\sigma };b_{𝐩^{}\sigma ^{}}^{}]=\delta _{\mathrm{𝐩𝐩}^{}}\delta _{\sigma \sigma ^{}}`$
$`[a_{𝐩\sigma };a_{𝐩^{}\sigma ^{}}]=[b_{𝐩\sigma };b_{𝐩^{}\sigma ^{}}]=[a_{𝐩\sigma };b_{𝐩^{}\sigma ^{}}]=[a_{𝐩\sigma };b_{𝐩^{}\sigma ^{}}^{}]=0`$ (52)
At the same time the operators $`a_{𝐩0},b_{𝐩0}`$ connected with the scalar part (42) must obey different rules of quantization:

$`[a_{𝐩0};a_{𝐩^{}0}^{}]=\delta _{\mathrm{𝐩𝐩}^{}}[b_{𝐩0};b_{𝐩^{}0}^{}]=\delta _{\mathrm{𝐩𝐩}^{}}`$
$`[a_{𝐩0};a_{𝐩^{}0}]=[b_{𝐩0};b_{𝐩^{}0}]=[a_{𝐩0};b_{𝐩^{}0}]=[a_{𝐩0};b_{𝐩^{}0}^{}]=0`$ (53)
The odd minus in eq.(4) is a price that must be paid for the Lorentz-covariance. It immediately causes the appearance of indefinite metrics in Hilbert space $``$ which is typical for Gupta-Bleuler quantization. However, since longitudinal mode $`B^\pm `$ (42) makes no physical sense, one can eliminate it by considering the physical subspace $`_{phys}`$ without these scalar excitations. Similar arguments can be found in many textbooks on quantum field theory . Loop calculations in a background electromagnetic field involve Feynman propagator being defined by the formula:

$$D_{\mu \nu }(x,x^{})=i<0|TW_\mu ^+(x)W_\nu ^{}(x^{})|0>,$$
(54)
and which is a Green function of the equation (41)

$$\left[\left((_\alpha +ieA_\alpha )^2+m_W^2\right)g^{\mu \nu }+(1/\xi 1)(^\mu +ieA^\mu )(^\nu +ieA^\nu )+2ieF^{\mu \nu }\right]D_{\nu \lambda }(x,x^{})=\delta _\lambda ^\mu \delta ^4(xx^{})$$
(55)
Inserting expansions (50) and (51) into eq.(54) and employing the commutation rules (4), (4), one can obtain an explicit covariant form of the $`W`$-boson propagator:

$$D_{\mu \nu }(x,x^{})=\frac{d^4p}{(2\pi )^4}\frac{W_{\mu \alpha }(x,p)W_{\beta \nu }^{}(x^{},p)}{m_W^2p^2i0}\left[g^{\alpha \beta }+(\xi 1)\frac{p^\alpha p^\beta }{p^2\xi m_W^2+i0}\right].$$
(56)
This representation is based on the eigen functions of the differential operator corresponding to eq.(45)

$`W_{\mu \alpha }(x,p)=\mathrm{exp}[ipxi{\displaystyle \underset{0}{\overset{kx}{}}}({\displaystyle \frac{epA}{pk}}{\displaystyle \frac{e^2A^2}{2pk}})d\phi ]`$
$`\left\{g_{\mu \alpha }+{\displaystyle \frac{e}{pk}}\left(k_\mu A_\alpha (kx)k_\alpha A_\mu (kx)\right){\displaystyle \frac{e^2}{2(pk)^2}}A^2(kx)k_\mu k_\alpha \right\}`$ (57)
Let us note that classical solution (4) is simply a particular form of the eigen function (4) taken on the mass shell and contracted through the last index with the ordinary polarization vector. For an arbitrary value of gauge fixing parameter $`\xi `$ a proper time representation of the propagator (54) is rather awkward. However, practically one can use a diagonal gauge ($`\xi =1`$) in which expression (56) can be rewritten in the most elegant way:

$`D_{\mu \nu }(x,x^{})={\displaystyle \frac{\mathrm{\Phi }(x,x^{})}{(4\pi )^2}}{\displaystyle \underset{0}{\overset{\mathrm{}}{}}}{\displaystyle \frac{ds}{s^2}}\mathrm{exp}[ism^2i{\displaystyle \frac{(xx^{})^2}{4s}}+{\displaystyle \frac{ie^2s}{(kxkx^{})}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}A_\mu (\phi )A^\mu (\phi )d\phi `$
$`{\displaystyle \frac{ie^2s}{(kxkx^{})^2}}{\displaystyle \underset{kx^{}}{\overset{kx}{}}}d\phi _1{\displaystyle \underset{kx^{}}{\overset{kx}{}}}d\phi _2A_\mu (\phi _1)A^\mu (\phi _2)]\{g_{\mu \nu }+{\displaystyle \frac{2es}{(kxkx^{})}}[k_\mu (A_\nu (kx)A_\nu (kx^{}))`$
$`k_\nu (A_\mu (kx)A_\mu (kx^{}))]{\displaystyle \frac{2e^2s^2}{(kxkx^{})^2}}(A_\lambda (x)A_\lambda (x^{}))^2k_\mu k_\nu \}`$ (58)
## 5 General Features of Quantum Processes in Strong Electromagnetic Fields

Many distinctive properties of particle physics with intense electromagnetic background can be traced in the reaction of $`W`$-boson and neutrino production via the decay of a massive charged lepton $`\mathrm{}`$. It is well known that the process $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ should be forbidden in vacuum if the lepton mass is small in comparison with the ones of final particles. However, in the presence of electromagnetic fields the restriction $`m_{\mathrm{}}>m_W+m_\nu `$ can be removed and there appears a possibility to observe an effect analogous to quantum transmission of the potential barrier through tunneling. The matrix element of the reaction $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ can be calculated employing the wave functions of the charged lepton $`\mathrm{}`$ (28) and the one of W-boson (4):

$$S_{fi}=\frac{ig}{2\sqrt{2}}d^4xW_\mu ^{}(x,p^{})\overline{\nu _{\mathrm{}}}(x,p^{\prime \prime })\gamma ^\mu (1+\gamma ^5)\mathrm{\Psi }_{\mathrm{}}^{}(x,p)$$
(59)
In order to obtain the probability of the decay $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ per unit of time one must integrate $`|S_{fi}|^2`$ over the phase volume of final particles and sum up their polarizations. In the case of a constant crossed field ($`=𝒢=0`$) the final result can be written as :

$`{\displaystyle \frac{dP}{du}}(\mathrm{}^{}W^{}\nu _{\mathrm{}})={\displaystyle \frac{G_Fm_W^4}{4\pi ^2\sqrt{2}p_0}}\{[1{\displaystyle \frac{(m_{\mathrm{}}^2+m_\nu ^2)}{2m_W^2}}{\displaystyle \frac{(m_{\mathrm{}}^2m_\nu ^2)^2}{2m_W^4}}]\mathrm{\Phi }_1(z)+`$
$`+{\displaystyle \frac{2m_e^2(\chi ^2u)^{1/3}}{m_W^2(1u)^{4/3}}}[1+u^2+(1u)^2{\displaystyle \frac{(m_{\mathrm{}}^2+m_\nu ^2)}{2m_W^2}}]\mathrm{\Phi }^{}(z)\}`$ (60)
This formula represents the differential probability of the process with respect to the variable

$$u=\frac{(p_\mu ^{\prime \prime }k^\mu )}{(p_\lambda k^\lambda )}=\frac{(p_\mu ^{\prime \prime }F^{\mu \alpha }F_{\alpha \beta }p^\beta )}{(p_\mu F^{\mu \alpha }F_{\alpha \beta }p^\beta )},$$
(61)
where $`p_\mu ^{\prime \prime },p_\mu ^{}`$ are momenta of neutrino and $`W`$-boson, respectively, while $`p_\mu `$ stands for the one of initial lepton. Note that in a crossed electromagnetic field $`F_{\mu \lambda }`$ the wave vector $`k_\mu `$ must be substituted by the following expression:

$$k^\mu =\frac{e^2(pk)}{m_e^6\chi ^2}F^{\mu \alpha }F_{\alpha \beta }p^\beta (p_0𝐄^2𝐩(𝐄\times 𝐁);𝐄(\mathrm{𝐄𝐩})+𝐁(\mathrm{𝐁𝐩})𝐩(\mathrm{𝐁𝐁})+p_0(𝐄\times 𝐁)).$$
(62)
This interrelation meaning that relativistic particles perceive the electromagnetic background as an incident plane-wave have been noticed long ago by Williams and Weizsäcker  who proposed the famous equivalent-photon approximation (see also ). In the case of a constant homogeneous electromagnetic field being considered the kinematics of the reaction $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ resembles the one with real photons and is determined by the following constraint of 4-momentum conservation:

$$p_\mu +\frac{e^2M^2}{2m_e^6\chi ^2}F_{\mu \alpha }F^{\alpha \beta }p_\beta =p_\mu ^{}+p_\mu ^{\prime \prime },$$
(63)
where $`M^2`$ represents the energy deficit obstructing the decay $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ in a vacuum without fields.

$$M^2=(p^{}+p^{\prime \prime })^2p^2=m_W^2+m_\nu ^2m_{\mathrm{}}^2+2(p^{}p^{\prime \prime })(m_W+m_\nu )^2m_{\mathrm{}}^2$$
(64)
Now we see that for light lepton ($`m_{\mathrm{}}<m_W+m_\nu `$) $`M^2`$ is always positive and the reaction could proceed only due to the energy-momentum borrowed from the background field.

$$q^\mu =\frac{M^2}{2(pFFp)}F^{\mu \alpha }F_{\alpha \beta }p^\beta q^2=\frac{M^4}{4m_e^2}\left[\left(\frac{a}{\chi ^2}\right)^2\frac{2b}{\chi ^2}\right]0$$
(65)
If we suppose the lepton $`\mathrm{}`$ to be a hypothetical heavy particle of the fourth generation ($`m_{\mathrm{}}>m_W+m_\nu `$) then for some values of neutrino and $`W`$-boson momenta $`M^2`$ can become negative. This reflects another possibility when the background field takes some amount of energy from the initial lepton and absorbs it. Thus electromagnetic environment reveals itself as an additional virtual particle with momentum $`q^\mu `$ (65) that induces the decay $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ regardless whether it is forbidden or not when the external field is absent.

All said above is described by the formula (5) which involves special functions

$$\mathrm{\Phi }(z)=\underset{0}{\overset{\mathrm{}}{}}𝑑t\mathrm{cos}(zt+t^3/3),\mathrm{\Phi }_1(z)=\underset{z}{\overset{\mathrm{}}{}}𝑑t\mathrm{\Phi }(t),\mathrm{\Phi }^{}(z)=\frac{d\mathrm{\Phi }(z)}{dz}$$
(66)
which are known in mathematics as the functions of Airy (see the Appendix). The argument $`z`$ is proportional to the quantity $`M^2`$ (64) and in terms of spectral variable $`u`$ (61) it can be written as follows:

$$z=\frac{m_W^2u+m_\nu ^2(1u)m_{\mathrm{}}^2u(1u)}{m_e^2[\chi u^2(1u)]^{2/3}}$$
(67)
Note that due to the constraint (63) the momentum projection on the wave vector (62) is conserved

$$p_\mu F^{\mu \alpha }F_{\alpha \beta }p^\beta =(p_\mu ^{}+p_\mu ^{\prime \prime })F^{\mu \alpha }F_{\alpha \beta }p^\beta $$
(68)
that fixes the interval of variations for $`u`$: $`0u1`$. Now if the reaction $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ is forbidden in vacuum, then $`z>0`$ for all $`u[0;1]`$. In the opposite case: $`m_{\mathrm{}}>m_W+m_\nu `$ there is a domain of negative values

$`z<0foru_1<u<u_2,`$
$`u_{2,1}={\displaystyle \frac{1}{2m_{\mathrm{}}^2}}\left[m_{\mathrm{}}^2+m_\nu ^2m_W^2\pm \sqrt{[m_{\mathrm{}}^2(m_W+m_\nu )^2][m_{\mathrm{}}^2(m_Wm_\nu )^2]}\right],`$ (69)
which exactly coincide with the phase volume of the process at hand provided that external field is omitted.

The sign and magnitude of the argument $`z`$ exert primary control over the probability of any quantum process in the background electromagnetic field. This is caused by a particular behaviour of the Airy functions that is examined carefully in the Appendix. Let us only mention that for large positive values of $`z`$ these functions decrease exponentially as $`\mathrm{exp}(\frac{2}{3}z^{3/2})`$. So the probabilities of reactions being forbidden under normal conditions must contain a factor of exponential suppression. This exponent can be treated as the well-known semiclassical probability of quantum transmission below the potential energy barrier. The index of the power increasing proportionally with $`M^2`$ (64) characterizes how much energy is borrowed from the background field. The aforesaid arguments can be confirmed by the asymptotic estimate of the total probability of the decay $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ in the domain of weak electromagnetic fields. Integrating eq.(5) through the use of the saddle-point method one obtains.

$$P(\mathrm{}^{}W^{}\nu _{\mathrm{}})=\frac{2G_Fm_Wm_e^3\chi }{9\pi \sqrt{6}p_0}\mathrm{exp}\left(\frac{\sqrt{3}m_W^3}{\chi m_e^3}\right),for\chi \left(\frac{m_W}{m_e}\right)^3$$
(70)
Here it is implied that masses of particles involved in the reaction satisfy real phenomenological conditions ($`m_\nu m_{\mathrm{}}m_W`$). We see that all leptons being known up to now can produce $`W`$-bosons only if the field strength $`E`$ and lepton energy $`p_0`$ take the values that are comparable with the estimate:

$$\chi =\left(\frac{p_0}{m_e}\right)\left(\frac{E}{E_{cr}}\right)\sqrt{3}\left(\frac{m_W}{m_e}\right)^36.710^{15}$$
(71)
In other words this exotic process could occur perhaps only in the Early Universe when the energies $`p_010^{12}GeV`$ were accessible. This ultra-relativistic domain is described by the formula which is derived from eq.(5) by means of the expansions (99), (100) (see the Appendix).

$$P(\mathrm{}^{}W^{}\nu _{\mathrm{}})=\frac{3G_Fm_Wm_e^3\chi }{2\pi \sqrt{6}p_0},for\chi \left(\frac{m_W}{m_e}\right)^3$$
(72)
The above results can be summed up as the following statement. If the reaction proceeding by the background electromagnetic field is forbidden in vacuum, then for weak fields its probability is exponentially suppressed by the factor which has the lacking energy borrowed from the background field as the power index. This suppression becomes unimportant only in the domain of extremely strong fields and ultra-relativistic energies that practically are inaccessible in experiments, at least now.

The situation changes dramatically for the processes that could take place when the background field is absent. For example, if we consider the decay of a hypothetical heavy lepton of the fourth generation ($`m_{\mathrm{}}>m_W+m_\nu `$), then from eq.(5) it follows.

$`P(\mathrm{}^{}W^{}\nu _{\mathrm{}})={\displaystyle \frac{G_F\sqrt{2}}{16\pi p_0}}\sqrt{\left[1\left({\displaystyle \frac{m_W+m_\nu }{m_{\mathrm{}}}}\right)^2\right]\left[1\left({\displaystyle \frac{m_Wm_\nu }{m_{\mathrm{}}}}\right)^2\right]}\times `$
$`\left[(m_{\mathrm{}}^2m_\nu ^2)^2+m_W^2(m_{\mathrm{}}^2+m_\nu ^2)2m_W^4\right]+O(\chi ^2)`$ (73)
This expression exactly reproduces the probability of the reaction $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ which can be calculated through the standard quantum field theory technique. Note, that it arises from the formal substitution $`\mathrm{\Phi }_1(z)\pi \theta (z)`$ ($`\theta (x)`$ is the Heaviside step function) giving the first term in the asymptotic expansion of eq.(5) when $`\chi 0`$. Thus we see that the probability of the permitted in vacuum process acquires in the background field corrections proportional to $`\chi ^2`$ (the omitted term $`O(\chi ^2)`$ in eq.(5)). However these corrections which become essential only in the domain (71) can be neglected if $`\chi (m_W/m_e)^3`$.

## 6 New Processes in Background Electromagnetic Fields

Now I would like to discuss the possibilities for searching new physics in background electromagnetic fields. We have already noticed that the relativistic particle penetrating through the electromagnetic environment feels the background field as a beam of real photons which can supply it with substantial amount of energy. If the field strength is sufficiently high then new channels of reactions can be opened and the particle decays into more heavy species. Let us consider such a reaction when a light charged particle $`A^\pm `$ produces two others, one of which is charged too ($`B^\pm `$) while the other is neutral ($`C^0`$). We assume that the process $`A^\pm B^\pm C^0`$ is forbidden in a vacuum ($`m_A<m_B+m_C`$) and wish to estimate the field-energy domains where it could occur. As we have explained in the previous section, in general the probability of the decay $`A^\pm B^\pm C^0`$ having a form similar to eq.(5) can be represented through the Airy functions which depend on the following argument:

$$Z_\pm =\frac{m_B^2u+m_C^2(1u)m_A^2u(1u)}{m_e^2[\chi u^2(1u)]^{2/3}}$$
(74)
For weak electromagnetic fields that can be handled in experiments the above probability must contain the exponential factor which is connected with the minimal positive value of the argument $`Z_\pm `$ (74):

$$P(A^\pm B^\pm C^0)\mathrm{exp}(\frac{2}{3}Z_{\pm min}^{3/2})\mathrm{exp}(\gamma )$$
(75)
Using eq.(74) one can explicitly calculate the power index $`\gamma `$ for the arbitrary relationship among the masses of all particles involved in the reaction:

$$\gamma =\frac{\sqrt{6}}{16\chi }\left\{128\lambda _1^2\lambda _2+80\lambda _1(\lambda _1+\lambda _2\lambda _0)^2+\frac{(\lambda _1+\lambda _2\lambda _0)^4}{\lambda _2}\left[\left(1+\frac{32\lambda _1\lambda _2}{(\lambda _1+\lambda _2\lambda _0)^2}\right)^{3/2}1\right]\right\}^{1/2},$$
(76)
where

$$\lambda _0=(\frac{m_A}{m_e})^2,\lambda _1=(\frac{m_B}{m_e})^2,\lambda _2=(\frac{m_C}{m_e})^2$$
(77)
This quantity governing the range of exponential suppression makes a crucial impact on the probability and fixes a certain threshold of the reaction induced by the electromagnetic background. In the domain of the very small values of parameter $`\chi `$ the probability of the decay $`A^\pm B^\pm C^0`$ becomes negligible on the account of extremely large magnitude of the power index ($`\gamma 1`$). So the reaction could proceed with a feasible rate only if $`\gamma 1`$. I will refer to this condition as to the threshold of the reaction designating the corresponding value of $`\chi `$ by a zero subscript ($`\chi _0`$).

Let us analyze the consequences of the general expression (76). Suppose that neutral particle being produced in the decay $`A^\pm B^\pm C^0`$ is very light in comparison with the others. Then the power index $`\gamma `$ can be approximated as follows.

$$\gamma =\frac{\sqrt{3}}{\chi }\left(\frac{m_B}{m_e}\right)^3\left(1\frac{m_A^2}{m_B^2}\right),form_Cm_A,m_B.$$
(78)
A particular case of this formula have been employed in the previous section for the process of $`W`$-boson and neutrino production. The threshold of this reaction $`\mathrm{}^{}W^{}\nu _{\mathrm{}}`$ coinciding with the estimate (71) can be calculated without resort to the explicit form of the probability (70). The only thing that we need is the power index (78) being connected with the threshold value of $`\chi `$ by a simple equation $`\gamma =1`$. This reasoning can be extended to some other reactions being of practical interest. Substituting the values of electron, muon and pion masses in eq.(78) we obtain characteristic domains of the following processes induced by the external electromagnetic background.

$`e^{}\mu ^{}\nu _e\stackrel{~}{\nu }_\mu ,e^+\mu ^+\stackrel{~}{\nu }_e\nu _\mu ,\chi _0=1.5310^7`$
$`e^{}\pi ^{}\nu _e,e^+\pi ^+\stackrel{~}{\nu }_e,\chi _0=3.5310^7`$
$`\mu ^{}\pi ^{}\nu _\mu ,\mu ^+\pi ^+\stackrel{~}{\nu }_\mu ,\chi _0=1.5110^7`$ (79)
On the other hand, if we assume that the mass of the charged particle emerged in the decay $`A^\pm B^\pm C^0`$ is the smallest parameter, then eq.(76) can be reduced to

$$\gamma =\frac{\sqrt{3}m_B}{\chi m_e}\left(\frac{m_C}{m_e}\right)^2\left(1\frac{m_A^2}{m_C^2}\right),form_Bm_A,m_C.$$
(80)
Now it is evident that the factor of exponential suppression grows up from the non-zero mass of the charged particle being produced, because the power index $`\gamma `$ (80) vanishes in the limit $`m_B0`$. By virtue of the fact that there are no charged particles with masses less than the one of electron, it is reasonable to look for the reactions with electrons and positrons in final states. For example, the processes of the inverse beta decay being forbidden in vacuum have the threshold that lies far below the ones of the lepton reactions (6). (This is also caused by a small mass difference between proton and neutron.)

$$p^+ne^+\nu _e,\stackrel{~}{p}^{}\stackrel{~}{n}e^{}\stackrel{~}{\nu }_e,\chi _0=1.6110^4$$
(81)
Another interesting possibility arises in connection with the process of neutral particle emission in the intense electromagnetic field $`A^\pm A^\pm C^0`$. The power index $`\gamma `$ that determines the rate of this reaction can be deduced from eq.(76) by a formal replacement $`m_B=m_A`$.

$$\gamma =\frac{\sqrt{6}}{16\chi }\left(\frac{m_C}{m_e}\right)^3[(1+32\frac{m_A^2}{m_C^2})^{3/2}1+80(\frac{m_A}{m_C})^2+128(\frac{m_A}{m_C})^4]^{1/2}.$$
(82)
Now provided that the neutral particle is much heavier than its fore-runner we can rewrite eq.(82) in the form.

$$\gamma =\frac{\sqrt{3}m_A}{\chi m_e}\left(\frac{m_C}{m_e}\right)^2,form_Am_C.$$
(83)
Here it is relevant to quote the reactions of pion radiation as an example governed by the last formula.

$$e^{}e^{}\pi ^0,e^+e^+\pi ^0,\chi _0=1.2110^5.$$
(84)
In the opposite case corresponding to the relatively small mass of the emitted particle the factor of exponential suppression appears with the following index.

$$\gamma =\frac{\sqrt{3}m_C}{\chi m_e}\left(\frac{m_A}{m_e}\right)^2,form_Am_C.$$
(85)
So the process of pion radiation from protons could occur only if the field strength and proton energy exceed significantly the ones of the reactions (84).

$$p^+p^+\pi ^0,\stackrel{~}{p}^{}\stackrel{~}{p}^{}\pi ^0,\chi _0=8.9010^8.$$
(86)
Next we intend to discuss the possibilities of charged particles production via the field induced decay $`A^0B^+C^{}`$. The argument of the Airy functions through which the probability of this reaction can be expressed has a form that slightly differs from eq.(74).

$$Z_0=\frac{m_B^2u+m_C^2(1u)m_A^2u(1u)}{m_e^2[\chi u(1u)]^{2/3}}$$
(87)
Following the reasoning given above it is straightforward to calculate the exponential factor which governs the probability rate in the domain of weak electromagnetic fields.

$$P(A^0B^+C^{})\mathrm{exp}(\frac{2}{3}Z_{0min}^{3/2})exp(\gamma )$$
(88)
At this stage, in order to obtain precise estimates of the threshold of the process $`A^0B^+C^{}`$ we need additional assumptions concerning with the masses of particles being implicated. If the initial neutral particle $`A^0`$ has a very small mass as compared to the charged ones ($`m_Am_B,m_C`$) then the power index of eq.(88) can be written as follows:

$$\gamma =\frac{\sqrt{2}}{6\chi }\left[(\lambda _1^2+14\lambda _1\lambda _2+\lambda _2^2)^{3/2}\lambda _1^3+33\lambda _1^2\lambda _2+33\lambda _1\lambda _2^2\lambda _2^3\right]^{1/2}$$
(89)
The last equation describes, for example, the reactions being caused by neutrinos moving throughout the intense electromagnetic field. This is especially the case when one considers the process which is cross symmetric to the muon decay.

$$\nu _\mu e^+\mu ^{}\nu _e,\stackrel{~}{\nu }_\mu e^{}\mu ^+\stackrel{~}{\nu }_e\chi _0=7.4110^4$$
(90)
Confronting the conditions appropriate to the above reaction with the threshold of another cross symmetric channel adduced in eq.(6) we see that reactions initiated by neutrinos require much less values of energy and background field strength than the ones with electrons. From this viewpoint it is more promising to look for new physics in external electromagnetic fields through the investigations of the processes induced by neutral particles. This conclusion acquires further confirmation from the asymptotic estimate of the power index (88) in the case when one of charged particles produced is extremely light.

$$\gamma =\frac{\sqrt{3}m_C(m_B^2m_A^2)}{\chi m_e^3},form_Cm_A,m_B.$$
(91)
Related conditions are characteristic for the following pion decays, which have the lowest thresholds among the all afore-mentioned reactions.

$$\pi ^0\pi ^+e^{}\stackrel{~}{\nu }_e,\pi ^0\pi ^{}e^+\nu _e\chi _0=8.3810^3,$$
(92)
It should also be noted that the power index (91) exactly coincides with the estimate (80) if the interchange $`m_Bm_C`$ is made. This fact explains once more the origin of the exponential suppression as a semiclassical probability of the field-induced quantum transmission of the electron from the Dirac sea to the upper continuum. The height of the potential barrier corresponding to the results (80), (91) can be easily calculated if one reproduces the reasoning which is relevant to the famous Klein paradox in quantum mechanics .

Finally, I would like to make some remarks about the processes of particle-antiparticle pairs production in the intense electromagnetic fields. The probability of the forbidden in a vacuum decay $`A^0B^+B^{}`$ is suppressed in the electromagnetic background by the factor which can be obtained from eq.(87), (88) under the assumption $`m_C=m_B`$ :

$$P(A^0B^+B^{})\mathrm{exp}\left[\frac{8m_B^3}{3\chi m_e^3}\left(1\frac{m_A^2}{4m_B^2}\right)^{3/2}\right]$$
(93)
This formula can be applied for the threshold evaluations related to the following reactions provoked by the external field.

$$\gamma e^+e^{},\nu _\mu \nu _\mu e^+e^{}\chi _0=2.67$$
(94)
$$\gamma \mu ^+\mu ^{},\nu _\mu \nu _\mu \mu ^+\mu ^{}\chi _0=2.3610^7.$$
(95)
We see that electron-positron pairs can be easily produced at the energies $`p_0>25GeV`$ if the background field strength is of the order $`E10^{10}V/m`$. Similar conditions have been created in CERN experiments  for the photons incident on the crystal of Germanium with a small angle to the axial direction. Neutrino fluxes of the like energies could give rise to the same effect caused by the coherent interactions with a regular arrangement of atoms in the single crystal. However due to the very small values of the neutrino couplings the process $`\nu _\mu \nu _\mu e^+e^{}`$ has an additional suppression by the factor $`(m_e/m_W)^410^{21}`$  which impedes its immediate experimental observation.

## 7 Conclusion

In this paper I have presented some methodological ideas providing the basis for practical calculations of quantum processes in very strong electromagnetic fields. The primary motivation of these studies is to look for new physics that could emerge due to the intense electromagnetic background. There is a hope that by combining modern achievements in the laser technique and in electromagnet construction with traditional methods of elementary-particle physics, one could obtain results that are inaccessible in other investigations.

Tracing some general features of quantum processes induced by strong electromagnetic fields I have estimated the energy-field domains that are needed in order to observe a number of reactions which are forbidden under normal conditions. At relatively small values of the field strength there remains the only one parameter $`\chi `$ (6) which has a crucial influence on the rate of the reaction and governs how the energy threshold moves with the increase of the field intensity. I drew a figure for illustrative purposes (see Fig.1) to display the characteristic domains relevant to the processes mentioned in the previous section. We see that present experimental situation is confined by the bound $`\chi 10`$. This conclusion comes from the estimates of electric fields extending over macroscopic distances along strings of atoms in single crystals . For example, the electric field along a $`<111>`$ axis in a crystal of Wolfram amounts to $`E510^{13}V/m`$ . For crystals of higher atomic number, the electric fields can be still larger than the value given above. However it is insufficient to reach the domain $`\chi 10^2÷10^5`$ as long as the energies of accelerated particles do not exceed the TeV scale. So we can expect that with the particle machines of the new generation there could appear a possibility for experimental observation of some reactions discussed in this paper.

## Acknowledgements

The author would like to thank the Theory Division of CERN for kind hospitality and financial support while part of this work was being done.

## Appendix

We have employed the Airy functions being connected with the ones described in mathematical handbooks in the following manner.

$$\mathrm{\Phi }(z)=\pi Ai(z),\mathrm{\Phi }_1(z)=\pi ^2\left[Ai(z)Gi^{}(z)Ai^{}(z)Gi(z)\right],\mathrm{\Phi }^{}(z)=\pi Ai^{}(z)$$
(96)
In accordance with standard mathematical notations these functions are defined through the integral representations

$$Ai(z)=\frac{1}{\pi }\underset{0}{\overset{\mathrm{}}{}}𝑑t\mathrm{cos}(zt+t^3/3),Gi(z)=\frac{1}{\pi }\underset{0}{\overset{\mathrm{}}{}}𝑑t\mathrm{sin}(zt+t^3/3)$$
(97)
and they obey the following differential equations

$$Ai^{\prime \prime }(z)zAi(z)=0,Gi^{\prime \prime }(z)zGi(z)=\frac{1}{\pi }.$$
(98)
The Airy functions that are the entire functions of the complex variable $`z`$ can be expanded in the Laurent series that converge everywhere ($`|z|<\mathrm{}`$)

$$Ai(z)=\frac{3^{2/3}}{\pi }\underset{n=0}{\overset{\mathrm{}}{}}\mathrm{\Gamma }\left(\frac{n+1}{3}\right)\mathrm{sin}\left(\frac{\pi }{3}\frac{2\pi n}{3}\right)\frac{(3^{1/3}z)^n}{n!},$$
(99)
$$Gi(z)=\frac{3^{2/3}}{\pi }\underset{n=0}{\overset{\mathrm{}}{}}\mathrm{\Gamma }\left(\frac{n+1}{3}\right)\mathrm{cos}\left(\frac{\pi }{3}\frac{2\pi n}{3}\right)\frac{(3^{1/3}z)^n}{n!}.$$
(100)
Many properties of the Airy functions can be deduced from their interrelation with the cylindrical (Bessel and McDonald) functions of the fractional index.

$$Ai(z)=\frac{1}{\pi }\sqrt{\frac{z}{3}}K_{1/3}\left(\frac{2}{3}z^{3/2}\right),Ai(z)=\frac{\sqrt{z}}{3}\left[J_{1/3}\left(\frac{2}{3}z^{3/2}\right)+J_{1/3}\left(\frac{2}{3}z^{3/2}\right)\right]$$
(101)
$$Ai^{}(z)=\frac{z}{\pi \sqrt{3}}K_{2/3}\left(\frac{2}{3}z^{3/2}\right),Ai^{}(z)=\frac{z}{3}\left[J_{2/3}\left(\frac{2}{3}z^{3/2}\right)J_{2/3}\left(\frac{2}{3}z^{3/2}\right)\right]$$
(102)
For large values of its arguments ($`|z|\mathrm{},|\mathrm{arg}z|<\pi `$) there is a possibility to approximate those by the following asymptotic estimates:

$$Ai(z)=\frac{z^{1/4}}{2\pi }\mathrm{exp}\left(\frac{2}{3}z^{3/2}\right)\underset{n=0}{\overset{\mathrm{}}{}}(1)^n\frac{\mathrm{\Gamma }(3n+1/2)}{(2n)!}\left(9z^{3/2}\right)^n$$
(103)
$$Gi(z)=\frac{1}{\pi z}\underset{n=0}{\overset{\mathrm{}}{}}\frac{(3n)!}{n!}\left(3z^3\right)^n$$
(104)


format: text

text:
# From Fermi-Pasta-Ulam Problem to Galaxies: The Quest for Relaxation

## 1 The Fermi-Pasta-Ulam system: Historical introduction

The numerical experiment by Fermi, Pasta and Ulam opened several new fields of research and appeared to be an invaluable source of inspiration for generations of physicists.

This work is outstanding at least by the following reasons:

1. It represents the first computer study of a nonlinear system;

2. The results contradicted the belief held since Poincare. Fermi had considered it ’a little discovery’ (as quoted by Ulam), thus immediately evaluating their extraordinary importance;

3. It was one of last Fermi’s works, completed after his ultimate death in 1954;

4. Remained unpublished during a decade;

5. Coincides in time with Kolmogorov’s theorem (1954), though FPU and Kolmogorov-Arnold-Moser (KAM) theory ’met’ only in 1966.

6. Inspired the discovery of solitons and numerous other studies;

7. Its results are not fully understood till now and the FPU model continues its inspiring mission also today, after half a century.

The circumstances around this work are now widely known, mainly due to the reminiscences of Fermi’s collaborator of Los Alamos period and coauthor of this and other works, renown mathematician Stanislav Ulam. We mention some aspects only.

The study was planned in the summer of 1952, during one of Fermi’s visits to Los Alamos, the main calculations were carried on in the summer of 1953, and final calculations in early 1955. It was documented as Los-Alamos internal report LA-1940 in May 1955 .

The calculations were performed on MANIAC, one of the first computers, together with the programming assistant John Pasta and with the help of Mary Tsingou.

Though Ulam had reported the results in several meetings and had published their short account , the report LA-1940 was fully published only in 1965, in the second volume of collected papers of Fermi . After that, Ulam had studied other systems as well, i.e. described by homogeneous quadratic transformations, a class of so-called Binary Reaction Systems in three and four variables and of cubic transformations in three variables, in collaboration with P.Stein and M.Menzel (Tsingou) .

What was the aim of FPU work?

According to Ulam :

”The ultimate aim was to discuss problems with more independent variables in the hope of obtaining material which would suggest some general features of behavior of systems with an infinite number of degrees of freedom, with nonlinear interaction terms as they occur between the oscillators in quantum theory or between the degrees of freedom in an electromagnetic field or a meson field. The mathematical possibilities concerning what might be called quasi-states were discussed with Fermi in this connection”.

FPU system represents a 1-dimensional system of N nonlinearly coupled oscillators given by a Hamiltonian

$$H(p,q)=\underset{j=1}{\overset{N+1}{}}[\frac{1}{2}p_j^2+\frac{1}{2}(q_jq_{j1})^2+\frac{\alpha }{3}(q_jq_{j1})^m],$$
(1)
at boundary conditions

$$q_j=q_{j+N}.$$

The calculations were performed for three type of perturbations, cubic, quadratic and of broken-linear function.

At $`m=3`$ ($`\alpha `$-FPU) and $`4`$ ($`\beta `$-FPU) and $`N=16,32,64`$ with zero velocities of particles at $`t=0`$, the system was run during time scales corresponding up to 80,000 cycles of the linear system.

Instead of expected thermalization, i.e. a trend towards equipartition over modes, a regular behavior has been observed, namely, ”mode 1 comes back to within one percent of its intimal value so that the system seems to be almost periodic” .

Ulam recounts :

The results were entirely different qualitatively from what even Fermi, with his great knowledge of wave motions, had expected.

Indeed, it was contradicting the widely held views based on the Poincare theorem, that such perturbed systems have to be chaotic.

By dramatic coincidence, at the same time that belief was shaken by a theorem published without a proof by Kolmogorov  in Moscow in 1954. Its proof was given almost a decade later, by Arnold and Moser in early 1960s, and Kolmogorov-Arnold-Moser (KAM) theory had to become an achievement of fundamental importance.

Fermi was well familiar with Poincare theorem. Already in 1923 he had applied Poincare’s method to prove that a system which has no other integrals than the energy one, has to be ergodic on the energy hypersurface .

It might be, though Ulam mentioned about this indirectly , that one of Fermi’s motivations for FPU study was to check his theorem via an experiment. Does the fact that Fermi immediately accepted the paradoxical result of the experiment, indicate his anticipation of new fundamental insights? What would be his reaction, if informed about Kolmogorov’s theorem?

Unfortunately these questions remained unanswered. Fermi died, even before FPU was completed, without delivering his scheduled Gibbs Lecture for 1955 in American Mathematical Society.

Kolmogorov’s theorem needed around 15 years to reach the wider physical audience in west.

Thus, FPU emerged both at the time of major breakthrough in the theory of perturbed nonlinear dynamical systems and at the moment of appearance of the first electronic computers.

Fermi’s intuition had clearly revealed both these aspects already at that time, as Ulam refers  to Fermi’s opinion on the importance of the ’understanding of nonlinear systems’ for the future fundamental theories, and the ’potentialities of the electronic computing machines’ and even mentions about Fermi’s learning of the actual coding (programming) during one summer.

The forthcoming decades indeed were marked by spectacular achievements both in the studies of nonlinear systems and in the development of the computer technique.

In 1958 Kolmogorov introduced the metric entropy, now known as Kolmogorov-Sinai (KS) entropy, as an invariant of the transformation, In 1959 he introduced K-systems, now known as Kolmogorov systems111Kolmogorov’s original notation for K-systems implied ’quasi-regular systems’..

Introducing the entropy as a metric invariant, Kolmogorov moved to probabilistic and geometric methods from the functional approach dominating since von Neumann and Hopf.

His works together with those of Sinai and Anosov in 1960s are considered to lead to the modern ergodic theory of smooth systems.

FPU work stimulated important studies, already before the publication in Fermi’s volumes. We mention only two. Tuck and Menzel (Tsingou) in 1961 observed a superperiod at more extensive studies of the same system. Their results remained unpublished until 1972, though Ulam mentioned them when commenting Fermi’s volume in 1965 . Still unpublished, FPU inspired the discovery of solitons by Zabusky and Kruskal  in 1965, ”as solitary wave pulses which pass through each other and do not sum in a usual linear-wave superposition.”

In forthcoming years FPU continued and still is continuing to serve as a reference point for numerous studies (for reviews and references see-). Already the fact, that the decade delay of its publication, which might be fatal for any other scientific result, had not affected its influence, indicates its extraordinary power.

In the present review, on the example of FPU, we will discuss the problem of thermalization of nonlinear systems associated with an important astrophysical problem, the relaxation of N-body gravitating systems. The structure of the paper is as follows.

First, we briefly discuss the association of FPU to the Poincare theorem and KAM theory. Then, we review the approach to the thermalization based on the Langevin equation, as of the simplest and most transparent version of the Brownian dynamics. Some results on the long-range interacting systems and iterated maps are presented thereafter and finally, we discuss the geometrical approach in stellar dynamics. Each topic is rather broad, and therefore we confine ourselves with brief accounts, without going into the details of applications.

## 2 Relaxation in Statistical Mechanics

Thus, the FPU puzzle concerned the thermalization, or relaxation towards equilibrium, one of the oldest problems in statistical physics.

On the one hand, the relaxation verifies the equilibrium statistical mechanics, since provides the basic equilibrium distribution functions — microcanonical, canonical and grand canonical — as the result of a dynamical relaxation process. On the other hand, the problem of thermalization constitutes the very subject of non-equilibrium statistical mechanics ,, since it is supposed to describe scenarios of relaxation and to give information on relaxation times scales.

Historically, the first step in this field was made by Boltzmann. His kinetic equation for weakly non-ideal dilute gases was the first physical example which derived thermalization starting from the microscopics dynamics,. In Boltzmann’s setup one considers a closed macroscopic system, and makes certain sensible kinetic assumptions on its dynamics which result in the kinetic equation. This is sufficient to derive the H-theorem, which states a monotonous increase of entropy, and as the main result one gets relaxation of the macroscopic system towards equilibrium (thermalization).

Unfortunately, this approach is not always easy to motivate and implement. In particular, the kinetic assumptions become mathematically rigorous only in rather artificial limits. For instance, it is known that a self-consistent description of thermalization for essentially non-ideal gases within Boltzmann’s approach meet serious problems .

The example of Fermi-Pasta-Ulam comes to verify that not all realistic systems need to satisfy conditions of Boltzmann’s kinetic equation.

Another approach to the thermalization problem was started by Einstein and Langevin . This is a drastic change in philosophy as compared with the Boltzmann’s approach, since one does not look for thermalization of the whole macroscopic system, but rather investigates its subsystem. This approach originated from the problem of Brownian motion, where one follows to the motion of one heavy particle immersed in a thermal bath.

In contrast to the Boltzmann’s one, in this approach the thermalization problem appears to be more transparent and more widely applicable, since after certain reasonable and rather universal assumptions on the structure of the bath, one derives thermalization of the Brownian particle via an exact dynamical analysis .

As a consequence, one gets relaxation towards Gibbs distribution. One of the essential advantages of this approach is that the extension to quantum-mechanical systems appears to be straightforward .

The third approach is provided by ergodic theory which classifies dynamical systems by so called statistical properties, including those of tending to the equilibrium state.

In particular, the property of mixing is defined via decay of time correlation functions, as necessary condition for relaxation, while the ergodicity alone is not sufficient for that.

Though this approach is known by difficulties of checking the sufficiency conditions for realistic physical problems, there are methods, among them are the geometrical ones, which enable to overpass those difficulties.

The geometric methods of ergodic theory have been advanced by Kolmogorov, Sinai, Anosov and others (see ), while in physical problems they were firstly applied by Krylov  based on the early results of Hedlund and Hopf.

In stellar dynamics, i.e. for the problem of N-body gravitating systems interacting with Newtonian potential, the Langevin equation approach has been used long time ago by Chandrasekhar and von Neumann . The geometrical methods of ergodic theory firstly have been applied in stellar dynamics in  (and their other papers quoted therein).

## 3 Poincare theorem, KAM and FPU

First, let us recall that N-dimensional system is considered as integrable if its first integrals $`I_1,\mathrm{},I_N`$ in involution are known, i.e. their Poisson brackets are zero. Then, as follows from the Liouville theorem if the set of levels

$$M_I=\{I_j(x)=I_j^0,j=1,\mathrm{},N\}$$

is compact and connected, then it is diffeomorphic to N-dimensional torus

$$T^N=\{(\theta _1,\mathrm{},\theta _N),modd2\pi \},$$

and the Hamiltonian system performs a conditional-periodic motion on $`M_I`$. Poincare theorem  states that for a system with perturbed Hamiltonian

$$H(I,\phi ,ϵ)=H_0(I)+ϵH_1(I,\phi ,ϵ),$$
(2)
where $`I,\phi `$ are action-angle coordinates, at small $`ϵ>0`$ no other integral exists besides the one of energy $`H=const`$, if $`H_0`$ fulfills the nondegeneracy condition,

$$det|\omega /I|0,$$
(3)
i.e. the functional independence of the frequencies $`\omega =H_0/I`$ of the torus over which the conditional-periodic winding is performed.

Though this theorem does not specify the behavior of the trajectories of the system on the energy hypersurface, up to 1950s it was widely believed that such perturbed systems have to be chaotic.

Kolmogorov’s theorem of 1954, now the main theorem of KAM theory, came to contradict to that belief, namely, formulating the conditions when the perturbed Hamiltonian systems can remain stable.

It states:

If the system (2) satisfies the nondegeneracy condition (3) and $`H_1`$ is an analytic function, then at enough small $`ϵ>0`$ most of non-resonant tori, i.e. tori with rationally independent frequencies satisfying the condition

$$n\omega _k0,$$
(4)
do not disappear and the measure of the complement of their union set $`\mu (M)0`$ at $`ϵ0`$.

The theorem says nothing about the limiting value of the perturbation $`ϵ`$. Later numerical studies showed that, for certain systems that theorem can remain valid even at not small perturbations. Also, the behavior in two-dimensional Hamiltonian systems when the torus separate the regions of three-dimensional system and those of higher dimensions have to be different, when the phase trajectories of the system can wander to other regions due to Arnold diffusion.

Izrailev and Chirikov in 1966 , were the first to match KAM and FPU, via the formulated resonance overlap (Chirikov) criterion. Though even earlier, in 1963, Ford and Waters quoted Kolmogorov’s address at the International Congress of Mathematicians of 1957:

”Kolmogorov claims to have shown that nonlinear systems, more general than the systems considered here, are not ergodic”, the importance of that theorem was not been recognized and the works of Arnold and Moser were unknown.

Now, after numerous studies it is clear that FPU does not satisfy the conditions of the KAM theorem, though its behavior is close to that of integrable systems and to KAM conditions. KAM-FPU association is still an open issue.

FPU does not satisfy Poincare theorem conditions either, and cannot be directly traced via Poincare recurrence time, as was attempted by some authors (see refs in ).

The rigorous understanding of the thermalization remains the challenge of FPU.

## 4 The Langevin equation

Let us first consider the derivation of Langevin’s equation within the Brownian motion (system-bath) approach. Instead of describing this derivation in its full generality, we will consider the simplest situation as to provide all essential mathematical details, and to make clear all necessary physical assumptions.

The Langevin equation is derived from the exact Hamiltonian description of a subsystem (Brownian particle) and a thermal bath, by tracing out the degrees of freedom of the bath. The influence of the particle on the bath is assumed to be sufficiently small. Thus, only the linear modes of the bath are excited, and the interaction of the particle with the bath is assumed to be linear. To be as pedagogic as possible, we first take a definite model for the bath, namely a collection of harmonic oscillators.

For the total Hamiltonian we thus assume

$`H_{tot}=H+H_B+H_I={\displaystyle \frac{p^2}{2m}}+`$ $`V(x)+{\displaystyle \underset{i}{}}\left[{\displaystyle \frac{p_i^2}{2m_i}}+{\displaystyle \frac{m_i\omega _i^2}{2}}x_i^2\right]`$ (5)
$`+{\displaystyle \underset{i}{}}\left[c_ix_ix+{\displaystyle \frac{c_i^2}{2m_i\omega _i^2}}x^2\right],`$
where $`H`$, $`H_B`$ are the Hamiltonians of the particle and the bath, and $`H_I`$ is the interaction Hamiltonian. $`p`$, $`p_i`$, $`x`$, $`x_i`$ are the momenta and coordinate of the Brownian particle and the linear modes of the bath. $`c_i`$ are coupling constants with the bath; their concrete expressions will be given later on. $`V(x)`$ is the confining potential of the particle, and $`m`$ and $`m_i`$ are the corresponding masses.

The Brownian particle was taken to be one-dimensional. This is by no means a restriction of generality, because the generalization towards three-dimensional particle, as well as $`N`$ interacting Brownian particles, is straightforward.

Notice that the $`x_i`$-terms form a complete square, since it includes a self-interaction term proportional to $`x^2`$. This guarantees that the total Hamiltonian $`H_{tot}`$ will be positive definite.

At the initial time $`t=0`$ the bath and the particle are decoupled. The distribution function of the particle is not specified.

The bath, on the other hand, is assumed to be in equilibrium at temperature $`T=1/\beta `$ (we set Boltzmann’s constant $`k_B=1`$ in the body of this work), and described by its Gibbs distribution.

This means

$$\rho _{SB}(0)=\rho (0)\rho _B(0)=\rho _S(0)\frac{\mathrm{exp}[\beta H_B(0)]}{Z_B},$$
(6)
where $`\rho (0)`$, $`\rho _B(0)`$ are the initial distribution functions of the particle and the bath, and

$$Z_B=\underset{i}{}\mathrm{d}p_i\mathrm{d}x_i\mathrm{exp}[\beta H_B]$$

is the partition sum of bath.

The Hamilton equations of motion for the bath modes read

$`\dot{x}_i={\displaystyle \frac{1}{m_i}}p_i,`$ (7)
$`\dot{p}_i=x_im_i\omega _i^2+c_ix`$ (8)
These equations are solved readily:

$`x_i(t)`$ $`=`$ $`x_i(0)\mathrm{cos}\omega _it+{\displaystyle \frac{p_i(0)}{m_i\omega _i}}\mathrm{sin}\omega _it+{\displaystyle \frac{c_i}{m_i\omega _i}}{\displaystyle _0^t}ds\mathrm{sin}\omega _i(ts)x(s)`$ (9)
$`p_i(t)`$ $`=`$ $`m_i\omega _ix_i(0)\mathrm{sin}\omega _it+p_i(0)\mathrm{cos}\omega _it+c_i{\displaystyle _0^t}ds\mathrm{cos}\omega _i(ts)x(s)`$ (10)
Thus, the Heisenberg equations of motion for the Brownian particle read

$`\dot{x}={\displaystyle \frac{1}{m}}p,`$ (11)
$`\dot{p}=V^{}(x)+{\displaystyle \underset{i}{}}c_ix_ix{\displaystyle \underset{i}{}}{\displaystyle \frac{c_i^2}{m_i\omega _i^2}}`$ (12)
Combined with Eq. (9) the last equation becomes

$`\dot{p}=V^{}(x)x(0)\gamma (t){\displaystyle _0^t}ds\gamma (ts)\dot{x}(s)+\eta (t),`$ (13)
where

$`\eta (t)=`$ $`{\displaystyle \underset{i}{}}c_i[x_i(0)\mathrm{cos}\omega _it+{\displaystyle \frac{p_i(0)}{m_i\omega _i}}\mathrm{sin}\omega _it]`$
$`=`$ $`{\displaystyle \underset{i}{}}\sqrt{{\displaystyle \frac{\mathrm{}c_i^2}{2m_i\omega _i}}}[a_i^{}(0)e^{i\omega _it}+a_i(0)e^{i\omega _it}],`$ (14)
$$\gamma (t)=\underset{i}{}\frac{c_i^2}{m_i\omega _i^2}\mathrm{cos}(\omega _it),$$
(15)
are the noise related to the unperturbed bath, and the friction kernel, respectively. Notice that in this exact derivation the back-reaction of the bath on the particle has been taken into account.

It is described by the integrals in eqs. (9), (10), and brings the damping terms $`x(t)\gamma (0)x(0)\gamma (t)_0^tds\gamma (ts)\dot{x}(s)`$ in eq. (12) to yield (13).

### 4.1 Drude-Ullersma spectrum

For some, but not all, applications it is benefitable to consider a fully explicit case for the bath. The bath is assumed to have uniformly spaced modes

$$\omega _i=i\mathrm{\Delta }i=1,2,3,\mathrm{}$$
(16)
and for the couplings we choose the Drude-Ullersma spectrum

$$c_i=\sqrt{\frac{2\gamma m_i\omega _i^2\mathrm{\Delta }}{\pi }\frac{\mathrm{\Gamma }^2}{\omega _i^2+\mathrm{\Gamma }^2}}$$
(17)
Here $`\mathrm{\Gamma }`$ is the characteristic Debye cutoff frequency of the bath, and $`\gamma `$ stands for the coupling constant; it has dimension $`kg/s`$.

The thermodynamic limit for the bath is taken by sending $`\mathrm{\Delta }0`$. Notice that this creates an infinite “Poincare” timescale $`1/\mathrm{\Delta }`$, implying that in the remaining approach the limit of “large times” always means the quasi-stationary non-equilibrium state where time is still much less than $`1/\mathrm{\Delta }`$.

In the limit $`\mathrm{\Delta }0`$ each coupling $`c_i\sqrt{\mathrm{\Delta }}`$ is very weak. The fact that the bath has many modes nevertheless induces its non-trivial influence.

At finite but small $`\mathrm{\Delta }`$ the system has an initial relaxational behavior, which at times of order $`1/\mathrm{\Delta }`$ is changed in a chaotic behavior.

It is customary to define the spectral density

$$J(\omega )=\frac{\pi }{2}\underset{i}{}\frac{c_i^2}{m_i\omega _i}\delta (\omega \omega _i)=\frac{\gamma \omega \mathrm{\Gamma }^2}{\omega ^2+\mathrm{\Gamma }^2}$$
(18)
It has the Ohmic behavior $`J\gamma \omega `$ for $`\omega \mathrm{\Gamma }`$, and $`\gamma `$ is called the interaction strength or damping constant.

As $`J(\omega )`$ is cut off at the “Debye” frequency $`\mathrm{\Gamma }`$, it is called a quasi-Ohmic spectrum.

It can then be shown that the friction kernel (15) becomes

$$\gamma (t)=\frac{2\gamma }{\pi }_0^{\mathrm{}}d\omega \frac{\mathrm{\Gamma }^2}{\omega ^2+\mathrm{\Gamma }^2}\mathrm{cos}\omega t=\gamma \mathrm{\Gamma }e^{\mathrm{\Gamma }|t|}$$
(19)
It is non-local in time, but on timescales much larger than $`1/\mathrm{\Gamma }`$ it may be replaced by $`\gamma \delta _+(t)`$.

The resulting Langevin equation reads

$$\dot{p}+\frac{\gamma \mathrm{\Gamma }}{m}_0^tdt^{}e^{\mathrm{\Gamma }(tt^{})}p(t^{})+V^{}(x)=\gamma \mathrm{\Gamma }e^{\mathrm{\Gamma }t}x(0)+\eta (t),$$
(20)
It follows from (6) that the noise is stationary Gaussian, and has the following properties:

$$K(tt^{})=\eta (t)\eta (t^{})=\frac{2\gamma T}{\pi }_0^{\mathrm{}}d\omega \frac{\mathrm{cos}\omega (tt^{})}{1+(\omega /\mathrm{\Gamma })^2}=2\gamma T\frac{\mathrm{\Gamma }}{2}e^{\mathrm{\Gamma }|tt^{}|},$$
(21)
where, again, the average is taken over the initial state of the bath (over realizations of the noise).

The connection between properties of the noise and the friction kernel is the consequence of fluctuation-dissipation theorem .

The classical white noise situation $`K(t)2\gamma T\delta (t)`$ is recovered by taking the high-$`\mathrm{\Gamma }`$ limit (Ohmic spectrum of the bath frequencies).

In this Ohmic regime one gets finally for $`t>0`$:

$$\dot{p}+\frac{\gamma }{m}p(t)+V^{}(x)=\eta (t),\eta (t)\eta (t^{})=2\gamma T\delta (tt^{}).$$
(22)
Finally we wish to mention that there are alternative ways to derive the Langevin equation , since many of its properties are rigidly determined by general statements like the fluctuation-dissipation theorem .

Nevertheless, we choose to focus on concrete models, because they show in detail how the Langevin equation arises from first principles, and thus are better suited for pedagogical purposes.

### 4.2 Fokker-Planck equation

Once we are interested by one-time quantities, the state of the Brownian particle is described by the distribution function $`\rho (p,x,t)`$:

$$\rho (p,x,t)=dp_0dx_0\rho (p_0,x_0,0)\delta (p(t)p)\delta (x(t)x),$$
(23)
where the average is taken with respect to the noise, $`\rho (p,x,t)`$ and $`\rho (p_0,x_0,0)`$ are final and initial distribution functions, while $`p(t)`$, $`x(t)`$ are the solutions of (22) for the corresponding initial conditions, and for a particular realization of the Gaussian noise.

The Gaussian noise is distributed according to the functional

$$𝒫[\eta ]\mathrm{exp}\frac{1}{2}dtds\eta (t)K^1(ts)\eta (s).$$
(24)
It will be clear that the distribution function $`W(y_1,y_2,t)`$, and the propagator $`\delta (p(t)y_1)\delta (x(t)y_2)`$ satisfy the same differential equation although with different initial conditions.

Therefore, we will use them interchangeably, making additional specifications only when it will be necessary.

Differentiating $`\rho (y_1,y_2,t)`$ we get

$$\frac{\rho (y_1,y_2,t)}{t}=\underset{k=1}{\overset{2}{}}\frac{(v_k\rho )}{y_k}\frac{}{y_1}\delta (p(t)y_1)\delta (x(t)y_2)\eta (t),$$
(25)
where

$$v_1=\frac{\gamma }{m}p(t)V^{}(x),v_2=\frac{p}{m}.$$
(26)
are the damped Newtonian acceleration and the velocity, respectively.

The following two relations can be verified

$$\eta (t)𝒫[\eta ]=dt^{}K(t^{}t)\frac{\delta 𝒫[\eta ]}{\delta \eta (t^{})},$$
(27)
$$\left[\begin{array}{c}\hfill \delta p(t)/\delta \eta (t^{})\\ \hfill \delta x(t)/\delta \eta (t^{})\end{array}\right]=\theta (tt^{})\left\{\mathrm{exp}_t^{}^tduA(x(u))\right\}_+\left[\begin{array}{c}\hfill 1\\ \hfill 0\end{array}\right],$$
(28)
$$A(x)=\left(\begin{array}{cc}\hfill \gamma /m& \hfill V^{\prime \prime }(x)\\ \hfill 1/m& \hfill 0\end{array}\right),$$
(29)
where $`\{\mathrm{}\}_+`$ means the chronological ordering (time-ordering).

The first equality holds since $`𝒫[\eta ]`$ is Gaussian; the second follows directly from Eq. (22) because of the relations

$`{\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}{\displaystyle \frac{\delta p(t)}{\delta \eta (t^{})}}={\displaystyle \frac{\gamma }{m}}{\displaystyle \frac{\delta p(t)}{\delta \eta (t^{})}}V^{\prime \prime }(x){\displaystyle \frac{\delta x(t)}{\delta \eta (t^{})}}+\delta (tt^{}),`$ (30)
$`{\displaystyle \frac{\mathrm{d}}{\mathrm{d}t}}{\displaystyle \frac{\delta x(t)}{\delta \eta (t^{})}}={\displaystyle \frac{1}{m}}{\displaystyle \frac{\delta p(t)}{\delta \eta (t^{})}}`$ (31)
Using Eqs. (25, 27, 28, 29) we obtain

$$\frac{\rho (y_1,y_2,t)}{t}=\underset{k=1}{\overset{2}{}}\frac{}{y_k}\left\{v_k\rho +\frac{}{y_1}\delta (p(t)y_1)\delta (x(t)y_2)\mathrm{\Phi }_{k1}(\{x\},t)\right\},$$
(32)
where $`\mathrm{\Phi }`$ is the $`2\times 2`$ matrix

$$\mathrm{\Phi }(\{x\},t)=_0^tdt^{}K(t^{})\left\{\mathrm{exp}_{tt^{}}^tduA(x(u))\right\}_+$$
(33)
and $`\mathrm{\Phi }_{k1}`$ is the corresponding matrix element $`K(t^{})=2\gamma T\delta (t^{})`$, implies

$`\mathrm{\Phi }_{k1}=\gamma T\delta _{k1}.`$ (34)
The final result is that we obtain a diffusion-type equation (Kramers-Klein equation) for $`\rho `$ itself :

$$\frac{\rho (p,x,t)}{t}=\frac{p}{m}\frac{\rho }{x}+\frac{}{p}\left([\frac{\gamma }{m}p+V^{}(x)]\rho \right)+\gamma T\frac{^2\rho }{p^2}.$$
(35)
### 4.3 Thermalization and H-theorem

It is now possible to provide a general argument that Eq. (35) ensures thermalization. First one notices that its stationary solution, obtained from (35) by setting $`_t\rho =0`$, is given by the standard Gibbs distribution:

$`\rho _{st}={\displaystyle \frac{1}{Z}}\mathrm{exp}\left[\beta \left({\displaystyle \frac{p^2}{2m}}+V(x)\right)\right],`$ (36)
where $`Z`$ is the statistical sum. Second, one notices that the H-function

$`H(t)={\displaystyle dpdx\rho (p,x,t)\mathrm{ln}\frac{\rho (p,x,t)}{\rho _{st}(p,x,t)}},`$ (37)
where $`\rho (t)`$ is the actual time-dependent solution, monotonously decreases towards zero, which is its equilibrium value $`\rho =\rho _{st}`$.

Indeed, differentiating $`H(t)`$, using (35) and requiring that both $`\rho `$ and $`\rho _{st}`$ are zero at infinity, one gets:

$`\dot{H}(t)=\gamma T{\displaystyle dpdx\rho (p,x,t)\left[\frac{}{p}\mathrm{ln}\frac{\rho (p,x,t)}{\rho _{st}(p,x,t)}\right]^2}<0.`$ (38)
Thus, all solutions of (35) relax with time towards the stationary Gibbs distribution (36). Eq. (38) provides an information on the global relaxation time; it is proportional to the damping constant $`\gamma `$ and temperature $`T`$. Notice that the

H-function provides an information-theoretic distance between the distribution functions $`\rho _{st}`$ and $`\rho `$ .

Thus, the open system approach provides a rigorous justification of thermalization. Let us outline the main assumptions of the present setup.

1) The Brownian particle was assumed to be in contact with a macroscopic thermal bath at equilibrium.

2) The influence of the particle to the bath was assumed to be weak enough, so that the non-linear modes of the bath are not excited.

3) The choice of initial conditions (6) was special, and moreover, it has to be special. One can motivate that for any macroscopic system there are initial states which do display anti-relaxational behavior . For sufficiently simple cases one can show that those states are exceptional .

4) The thermalization was obtained via a rigorous approach; no a priori assumptions on ergodicity and mixing were used. As the main result Eq. (38) provides a general expression for the relaxation time.

## 5 Long-Range Interacting Systems

The problem of understanding of the observed thermalization in FPU system thus includes the explanation of the existence of large relaxation times $`\tau _R`$ and their dependence on the number of particles $`N`$ and the energy density $`ϵ=E/N`$ (see  and references therein).

We will mention some of the remarkable results obtained in the recent period. The following dependence was found

$$\tau _R\frac{\sqrt{N}}{ϵ},$$
(39)
for finite $`N`$ and fixed wave number interval for the $`\beta `$-FPU model, i.e. the divergence of the equipartition time scale with increase of $`N`$. On the other hand, in the thermodynamical limit the dependence on $`N`$ disappears

$$\tau _Rϵ^3.$$
(40)
The conclusion of these studies is that, ”the process of relaxation to equipartition in the FPU model is not regulated by the microscopic chaotic instability but by the typical time in which an orbit diffuses in phase space, which is determined by the interaction among the phonons” .

While moving from FPU to 1 and 2-dimensional models with long-range interaction, to mimic the gravitational systems, the problem becomes more complex since acquires new thermodynamics, negative specific heat, etc.

Indeed, different results have been obtained at the analysis of the chaotic properties on the following 1-dimensional

$$H=\underset{i=1}{\overset{N}{}}\frac{p^2}{2m}+\frac{m^2}{2N}\underset{i,j=1}{}N(1\mathrm{cos}(\theta _i\theta _j))$$
(41)
and 2-dimensional systems

$`H`$ $`={\displaystyle \underset{i=1}{\overset{N}{}}}{\displaystyle \frac{p_{i,x}^2+p_{i,y}^2}{2m}}`$
$`+{\displaystyle \frac{1}{2N}}{\displaystyle \underset{i,j=1}{}}N[3\mathrm{cos}(x_ix_j)\mathrm{cos}(y_iy_j)\mathrm{cos}(x_ix_j)\mathrm{cos}(y_iy_j)],`$
and, particularly of their equipartition time scales. As a result, ”few scaling laws have been found to be universal for all these models. In the gaseous phase the maximal Lyapunov exponent vanishes as

$$N^{1/3}$$
(43)
in analogy found in for some specific gravitational systems” .

In a similar model, the Konishi-Kaneko iterated map

$$p_i^{n+1}=p_i^n+k\underset{j=1}{\overset{N}{}}\mathrm{sin}2\pi (x_j^nx_i^n),$$
(44)
$$x_i^{n+1}=x_i^n+p_i^{n+1};(mod1).$$
(45)
which can be represented as a system with Hamiltonian

$$H=\underset{i=1}{}N[\frac{p_i^2}{2}+\frac{k}{4\pi }\underset{ij}{}Ncos(2\pi (x_ix_j))]$$
(46)
the appearance of cluster-type configurations, of power law correlations , and Feigenbaum period-doubling bifurcations have been observed . Particularly, the bifurcation scale $`\delta =8.72\mathrm{}`$ can be estimated at the following conditions for the bifurcation points

$$\underset{j}{\overset{N}{}}|x_j^{n+1}x_j^n|<ϵ,(2^1=2),$$

$$\underset{j}{\overset{N}{}}|x_j^{n+2}x_j^n|<ϵ,\underset{j}{\overset{N}{}}|x_j^{n+3}x_j^{n+1}|<ϵ,(2^2=4),$$

for each $`k_n,n=1,2,\mathrm{}`$, respectively, where $`ϵ`$ is the accuracy of the obtained values of $`k_n`$.

The use of iterated maps can enable one to avoid the principal difficulties associated with $`N`$-body systems – the non-compactness of the phase space and singularity of Newtonian interaction – which are not always crucial at the study of a particular property of a given stellar system. Also, as we saw for Konishi-Kaneko map, the iterated maps can be informative in revealing the mechanisms of developing of chaos.

## 6 Gravitating Systems

While moving to the gravitating systems first one must recall that, the Henon-Heiles  system being one the first demonstrations of chaos at numerical experiments, was motivated by the study of the motion in the axisymmetric Galactic potential.

Since then, developments in nonlinear dynamics had various impacts on the problem of N-body gravitating systems (for reviews see the volumes , ).

The chaotic dynamics had offered principally new insight on the properties and the evolution of the Solar system. In the works of Laskar , Tremaine, Wisdom and their collaborators it was shown that the inner planets can possess chaotic properties influenced by the perturbations of the large planets. Particularly, Mercury and Venus undergo chaotic variations of the eccentricity of their orbit, which can lead to overlapping of orbits with subsequent drastic change of the orbit of Mercury, up to its escape from the Solar system.

By means of the frequency map analysis technique developed by Laskar, the stabilizing role of the Moon in the chaotic variations of the obliquity of the Earth has been shown with direct consequences for the climate and hence for the evolution of life on the Earth. Chaotic variations of the obliquity were shown to be important for Mars.

The idea of the frequency map technique  is in the search of quasiperiodic approximations of the solutions of the perturbation problem by means of finite number of terms

$$x_j(t)=x_j^0+\underset{k=1}{\overset{N}{}}a(m_k;\nu )\mathrm{exp}[i<m_k,\nu >t]$$
(47)
where $`\nu _j`$ are the frequency vectors. Note, that the system is close to an integrable one, and possess KAM (invariant) tori in the phase space.

So, proper theoretical treatment on the perturbed systems together with advanced computations had led to such basic results.

While profound role of chaotic effects in the evolution of Solar system is established, their role in stellar dynamics is often remained underestimated222See e.g. , where the authors needed half a page to transfer from linear to nonlinear effects in describing the relaxation of stellar systems. even though the latter ones are not close to integrability and KAM conditions, and one might expect that chaos has to be of more importance for them than for planetary systems.

The problem of relaxation of stellar systems can illustrate the situation.

The two-body relaxation formula derived many decades ago via the linear sum of contributions of several two-body interactions, is still the basic reference for stellar system, even though it predicts time scales exceeding the age of such well-relaxed systems as the elliptical galaxies.

A priori it is not excluded that the two-body relaxation time scale can be the solution of the nonlinear problem even when the perturbations of other particles are included, however this does not follow from the ’linear’ derivation borrowed from plasma physics, involving the Debye screening and the Coulomb logarithm, i.e. ignoring the long range interaction of gravitating particles.

While even for 1-dimensional systems the relaxation is essentially different for short and long range interactions.

Geometrical methods appear as a useful tool for the study of nonlinear effects of gravitating systems. In , the K-mixing and exponential tending to an equilibrium state (quasi-equilibrium) was shown analytically for perturbed spherical N-body systems.

That was reached via the proof that while the two-dimensional curvature of the configuration space for N-body systems is sign-indefinite, for spherical systems in large N limit it is strongly negative.

The formula derived for the relaxation time scale for the parameters of real stellar systems differs from the two-body relaxation formula, and by now is supported by (a) alternative theoretical derivation as mentioned above; (b) numerical simulations and (c) observational data of star clusters (see   and references therein).

The importance of chaotic effects in relaxation of stellar systems was shown by Pfenniger using the Lyapunov exponent technique. That study is especially remarkable given the limited applicability of Lyapunov exponents for many-dimensional systems.

Back in 1964 at a remarkable computer study of his time Miller  had observed exponential growth of errors of the coordinates and velocities, $`\mathrm{\Delta }=(\delta x^2+\delta v^2)`$ for systems with up to 32 gravitating particles.

The result of Miller remained not well understood for quite a long time. Recently, in some papers (e.g. and refs therein) it is identified with the exponential instability (K-mixing) shown in , even though their exponential time-scales are quite different. Since such opinion is moving from one paper to another, we will consider this issue in more details.

As follows from the Maupertuis principle, only trajectories in the region of the configuration space $`M=[EU(r)>0]`$ defined by the Riemannian metric,

$$ds^2=(EU(r_{1,1}\mathrm{},r_{N,3})\underset{a=1}{\overset{N}{}}\underset{i=1}{\overset{3}{}}dr_{a,i}^2$$
(48)
are corresponding to the N-body Hamiltonian system

$$H(p,r)=\underset{a=1}{\overset{N}{}}\underset{i=1}{\overset{3}{}}\frac{p_{a,i}^2}{2m_a}+U(r),$$
(49)
where

$$U(r)=\underset{a<b}{}\frac{Gm_am_b)}{|r_ar_b|}.$$
(50)
The Maupertuis parameterization and hence Riemannian metric guarantee the conservation of total energy of the system

$$H=E=const,$$

and therefore enables to study the hyperbolicity and mixing of the system with exponentially deviating trajectories.

The works discussing the exponential instability can be divided into two groups:

1. Based on Maupertuis parameterization: by Krylov, Anosov and Arnold, and in stellar dynamics, ;

2. Without using Maupertuis parameterization: numerical experiments by Miller, repeated by Lecar and Standish in 1960s and again in 1990s in many papers by Kandrup, Goodman, Heggie, Hut and others (see  for references).

In the first group of works the exponent is linked with mixing and relaxation via theorems of ergodic theory. For the second group, no such link exists. The reason is as follows.

Without Maupertuis parameterization, two points on two phase trajectories can deviate while the trajectories themselves may not deviate at all. Such apparent exponential deviation will be observed for systems of any number of particles, even of very few, and of any configurations, say of disk type. In contrary, the exponential deviation of exists only for spherical systems and only at large N-limit, and not for disk or few-particle gravitating systems.

Therefore, a statement that, the second group papers also refer to chaos and mixing, is the same as to assign chaos to the two-body system, where the same apparent deviation can be observed.

This is the reason that, only the dynamical (crossing) time scale is observed at such numerical experiments.

Numerical simulations  taking into account the Maupertuis parameterization are in agreement with the relaxation time scale derived in .

At the same time, we are not aware of any numerical experiments which confirm the two-body relaxation time scale333Donald Lynden-Bell to whom we consulted (V.G., July,2001), was not aware of such experiments either..

The fact that Miller’s experiment does not imply mixing, by no means diminishes the pioneering character of his work. However, now the results of theory of dynamical systems provide its interpretation.

The geometrical approach enabled to arrive also to other general conclusions on stellar systems:

elliptical and spiral galaxies should have different origins.

The spiral galaxies particularly, are more regular than ellipticals, i.e. do not possess the property mixing . Due to the existence of fundamental frequencies they can be described by Laskar’s frequency map technique

At Geneva workshop in 1993 one of us (V.G.) presented a list of ”10 key problems of stellar dynamics”. Now, about 10 years after, one can note that though a considerable work has been done, especially a numerical one due to advances of computer technique , still there is a long way to go to fully evaluate the role of nonlinear effects for stellar systems.

## 7 Epilogue

The Fermi-Pasta-Ulam problem was able to address profound issues of nonlinear systems, though studied mainly from physical concepts, rather than from rigorous mathematical ones. It still lacks answers to basic questions. As mentioned Pesin : ”It is a very interesting question the relation between KAM and FPU and whether FPU can serve as an ”experiment” for KAM”. One can expect that the activity in FPU problem will shift towards investigation of its rigorous aspects.

The relaxation problem remains a central issue for many systems, including open and long-range interacting ones. In astrophysical context it concerns N-body gravitating systems, galaxies and star clusters. We briefly discussed two powerful approaches which have been applied in stellar dynamics, the Brownian particle approach and the geometrical methods of ergodic theory.

In stellar dynamics also one has to expect more rigorous studies of nonlinear effects, including in the performance and interpretation of numerical experiments, a lesson given to us by Fermi, Pasta and Ulam 50 years ago.