| Russian Mathematical Surveys |
|
|
|
|
|
Boltzmann-type kinetic equations and discrete models A. V. Bobylevab a Keldysh Institute of Applied Mathematics of Russian Academy of Sciences b Peoples' Friendship University of Russia
Russian version:
Uspekhi Matematicheskikh Nauk, 2024, Volume 79, Issue 3(477), DOI: https://doi.org/10.4213/rm10161 
Bibliographic databases:
   1. IntroductionThe classical Boltzmann equation, which is the main mathematical tool for the description of rarefied gases, occupies a very specific place in mathematical physics. Speaking about equations of mathematical physics, we normally have in mind linear or non-linear equations in partial derivatives. Non-linear kinetic equations of Boltzmann type belong to a slightly different kind of equations, though they contain partial derivatives together with multiple integrals. These equations look too cumbersome at first glance. Perhaps, partly for this reason they are not included in the standard university courses for mathematicians. Nevertheless, the interest of mathematicians all over the world to this part of mathematical physics was growing fast in recent decades. It is very typical for the history of physics and mathematics that some parts of physics, which looked not very clear (even for physicists) after first fundamental discoveries, become rather a part of mathematics with time. This can be a slow process. Of course, this does not mean that their importance for physics has disappeared. The famous examples are Newton’s mechanics and Maxwell’s electrodynamics. Something like that happened and is still happening with the kinetic theory of gases founded by Maxwell [47], Boltzmann [16], and their predecessors approximately 150 years ago. The Boltzmann equation was first published in 1872 [16]. This equation has a very interesting history; we mention only one important point. Formally, this equation was supposed to describe more accurately — than the equations of gas dynamics known at that time — the behaviour of a rarefied gas of particles interacting in accordance with the laws of classical mechanics. The first conclusion made by Boltzmann in his paper cited above was his $H$-theorem, which formally proves the existence of a functional that is monotone decreasing it time on any solution of his equation. In fact, this was the discovery of the mechanical meaning of entropy, a result which is probably more important than the Boltzmann equation itself. The immediate consequence was that the solution of the Boltzmann equation cannot be invertible in time in contrast to the equations of Newtonian mechanics. This caused certain doubts in this equation, especially among mathematicians. The famous remark by Zermelo [61] based on Poincaré’s recurrence theorem was made at the end of the 19th century. This question and some others were clarified since then, and the validity of the Boltzmann equation has been justified rigorously, at least for short time intervals. We describe briefly below the progress in mathematical results for the Boltzmann equation. Only a few such results were obtained before the 1960s. These were, in particular, works by Hilbert [36], Carleman [19], Morgenstern [48], and Grad [34]. Of course, the great contribution of Bogolyubov [15], who proposed in 1946 his famous method of derivation of the kinetic equations from dynamics, should be mentioned. Some of the works cited above were written at the formal level of mathematical rigour, but this is inevitable at the early stage of the development of any mathematical theory related to physics. The Chapman–Enskog method [22], [27], invented by physicists for the transition to hydrodynamics from the Boltzmann equation, became a standard tool for the study of dynamical systems with small parameter. The number of rigorous mathematical results in kinetic theory began to grow faster in the 1960s–1970s. In particular, we mention: (i) the mathematically rigorous theory of the Cauchy problem for the linearized Boltzmann equation constructed in works by Arsen’ev [8] (for short-range intermolecular potentials) and Ellis and Pinsky [26] (for power-like potentials with Grad’s angular cut-off); (ii) the complete theory of the existence and uniqueness of solutions to the Cauchy problem for the spatially homogeneous Boltzmann equation for a wide class of potentials developed by Arkeryd [2]; (iii) the first global in time existence and uniqueness theorems for the Cauchy problem for the spatially inhomogeneous Boltzmann equation by Maslova and Firsov [46] and Ukai [55] under the assumption that the initial conditions are sufficiently close to an equilibrium; (iv) the proof by Lanford [42] of the first validation theorem for the Boltzmann equation in the Boltzmann–Grad limit for hard spheres on a time interval of an order of the particle mean free path. For the sake of brevity we do not discuss the related results obtained for other two famous classical non-linear kinetic equations introduced by Landau [40] and Vlasov [59], respectively. The kinetic theory and, in particular, the theory of the Boltzmann equation is gradually becoming a more-or-less regular part of mathematical physics. Mathematical conferences on this subject became more and more frequent in Europe, USA, and Japan during the last two decades of the 20th century. There was a sort of competition of pure mathematicians from different countries to prove a global existence theorem for the Boltzmann equation with initial data far from equilibrium. Finally, the result was obtained by DiPerna and Lions in 1989 [23]. The Fields Medal obtained by P.-L. Lions apparently attracted a great number of good young mathematicians to kinetic theory, especially in France. One of them, Villani, obtained another Fields Medal in 2010 for his works in this field of mathematical physics. He also wrote an excellent survey [58] of the mathematical results on the Boltzmann equation. This survey, together with its renewed on-line version and the book [21], contain most of the important mathematical results on the Boltzmann equation obtained before 2005. A review of some more recent results can be found, for example, in [10]. The development of kinetic theory in the last decades was connected with applications of its ideas to various unusual objects, which are far from traditional rarefied gases. It is also related to the rapid development of computers and numerical methods. Kinetic equations are now used for modelling traffic flows, the distribution of ‘active particles’ (viruses and so on) in biology, and socio-economic processes. Of course, many of these equations are perhaps ‘too young’ to become a subject of a rigorous mathematical study. On the other hand there are some classes of relatively ‘old’ kinetic equations, which have actively been used by physicists since the 1960s (see, for example, [30] and [60]). These are the so-called ‘wave kinetic equations’ used in the theory of weak turbulence. There are some mathematical results on these equations (see, for example, [28], [25], and the references there), but still many questions, like the long-time asymptotic behaviour of solutions, remain unclear. One goal of the present paper is to study these and similar equations (for example, the quantum Nordheim–Uehling–Uhlenbeck equation [49], [54]) from a unified point of view as particular cases of the general class of Boltzmann-type equations [12]. Our approach to this class of equations is, to a certain extent, close to the approach of Arkeryd in [5]. To this class of equations we apply some methods used for the classical Boltzmann equation. In particular, the properties of discrete models of Boltzmann-type equations are discussed. The main attention is paid to the wave kinetic equation (WKE): it is proved that solutions of any normal discrete model of the WKE tend to the equilibrium distribution as time tends to infinity. The consequences of this fact for solutions of the WKE are discussed in detail. It should be noted that we do not discuss in this paper important questions relating to the validation of various Boltzmann-type equations, that is, their derivation from some more general mathematical models. This is not an easy question. For example, it took more than a hundred years to find a mathematical proof (Lanford’s theorem [42]) of the connection between the classical Boltzmann equation and the system of $N$ hard spheres as $N$ tends to infinity. A recent mathematical result in [25] on the derivation of the WKE from the non-linear Schrödinger equation (in the presence of a random force field) looks very promising in this sense. The paper is organized as follows. Section 2 presents a concise introduction to key ideas and technique that lead from Hamiltonian mechanics to the classical Boltzmann equation. In principle, it allows one to read the paper without preliminary knowledge in kinetic theory. On the other hand important notation and some technical tools are introduced there. We begin with the $N$-particle system described by the Newtonian equations in the Hamiltonian form. Then we introduce the notion of an $N$-particle distribution function and Liouville’s equation. Following Grad’s version of the BBGKY-hierarchy, we finally present a formal derivation of the Boltzmann equation for hard spheres. A formal generalization of this equation to the case of any reasonable intermolecular potential is also performed. Then we study all general properties of the Boltzmann equation, including conservation laws and the $H$-theorem. The general Boltzmann-type kinetic equation is introduced at the beginning of § 3 as a natural generalization of the classical Boltzmann equation from § 2. This new equation depends on an arbitrary function $F(x_1,x_2;x_3,x_4)$ of four variables. For the Boltzmann equation we have $F=x_3 x_4-x_1 x_2$. Other kinetic equations mentioned above correspond to different polynomial forms of $F$. We show that all these equations can be studied from a unified point of view. Different forms of the general kinetic equation are constructed for the 3D-case (Proposition 3.1). The generalization to the $d$-dimensional case for $d \geqslant 2$ is also considered. The weak form of the general kinetic equation (the equation for average values) and conservation laws are also discussed. We define a class of functions $F$ that can lead to an analogue of Boltzmann’s $H$-theorem (Proposition 3.2). Then we introduce discrete kinetic models of the general Boltzmann-type kinetic equation by using the formal analogy with discrete velocity models of the Boltzmann equation. We define similarly the notion of a normal discrete model. Then we prove Theorem 3.1 on the main properties of normal discrete models which possess an analogue of the $H$-theorem. Section 4 is devoted to some properties of solutions to normal discrete kinetic models of WKE. The main result of this section is the proof of the convergence of any positive solution of this model to a unique equilibrium solution. This result is formulated at the beginning of § 4 (Theorem 4.1). The proof is given in the rest of that section. First we construct the solution for any positive initial data and prove its global in time existence and uniqueness (Lemma 4.1). Then we construct a positive stationary solution of the model and prove its uniqueness under fixed invariants (mass and energy) in Lemma 4.2. Then we improve some estimates for strictly positive initial data and complete the proof of convergence to an equilibrium by more or less standard methods of the theory of ODEs. At the end of § 4 we discuss a general problem of the approximation of the Boltzmann- type kinetic equation by a sequence of discrete models as the order of the model, that is, the number of its discrete points, tends to infinity. The results of the paper and related open problems are briefly discussed in § 5. 2. From particle dynamics to the Boltzmann equation2.1. $N$-particle dynamics and modelling of rarefied gasesWe consider $N \geqslant 1$ identical particles with mass $m=1$. This system is characterized by a $6N$-dimensional phase vector $Z_N=\{z_1,\dots,z_N\}$ with components $z_i=(x_i,v_i)$, where the $x_i \in \mathbb{R}^3$ and $v_i \in \mathbb{R}^3$ denote, respectively, the position and velocity of the $i\operatorname{th}$ particle, $i=1,\dots,N$. Usually, we assume below that particles interact via a given pair potential $\Phi(r)$, where $r>0$ denotes the distance between two interacting particles. We also assume that $\Phi(r) \to 0$ as $r \to\infty$. The equations of motion of the system have the following Hamiltonian form (see any textbook in classical mechanics, for example, [41]):
$$
\begin{equation}
\begin{gathered} \, \partial_t x_i=\frac{\partial H_N}{\partial v_i}\,,\qquad \partial_t v_i=-\frac{\partial H_N}{\partial x_i}\,, \\ H_N=\frac{1}{2}\sum_{i=1}^N |v_i|^2+\sum_{1\leqslant i < j\leqslant N} \Phi(|x_i-x_j|), \\ x_i(0)=x_i^{(0)},\quad v_i(0)=v_i^{(0)},\qquad i=1,\dots,N. \end{gathered}
\end{equation}
\tag{2.1}
$$
Thus the temporal evolution of the system can be understood as the motion of the phase point $Z_N=\{z_1,\dots,z_N\}$ in the phase space $\mathbb{R}^{6N}$. This motion obeys the conservation laws of energy $E_N \in \mathbb{R}$ and momentum $P_N\in \mathbb{R}^3$,
$$
\begin{equation}
E_N=H_N[Z_N(t)]=\mathrm{const},\qquad P_N=\sum_{i=1}^N v_i(t)=\mathrm{const},
\end{equation}
\tag{2.2}
$$
respectively. The laws (2.2) follow directly from equations (2.1). In principle, the $N$-particle system described above can be used to model real gases or fluids. Then the main problem is that the number of particles $N$ is of order of $10^{23}$ (Avogadro’s number). It is intuitively clear that the case of a rarefied gas is easier for description than the case of a dense gas. Indeed, in the low density limit we obtain a free molecular flow, that is, each particle moves independently with constant velocity. Consequently our main assumption will be as follows: a typical distance between particles is much greater than the effective diameter $d$ of the potential. It is known that the typical ‘size’ $d$ of a molecule of the air is roughly equal to $3.7 \cdot 10^{-8}$ cm, whereas the number density $\bar{n}$ of the air is about $2.5 \cdot 10^{19}$ cm$^{-3}$ under normal conditions. The inequality is the well-known criterion for an ideal gas. Note that $\delta \approx 10^{-3}$ for air under normal conditions on the surface of the Earth. This parameter $\delta$ is decreasing with height. Therefore, kinetic equations (in particular, the classical Boltzmann equation for rarefied gases) are important for applications in space science and technology.The simplest model intermolecular potential $\Phi(r)$ corresponds to particles interacting like hard spheres of diameter $d$. Then we formally obtain
$$
\begin{equation}
\Phi_{\rm HS}(r)=\begin{cases} \infty & \text{if} \ 0< r \leqslant d, \\ 0 & \text{otherwise}. \end{cases}
\end{equation}
\tag{2.4}
$$
Another well-known model corresponds to power-like repulsive potentials
$$
\begin{equation}
\Phi(r)=\frac{\alpha}{r^n}\,,\qquad \alpha > 0,\quad n\geqslant 1,
\end{equation}
\tag{2.5}
$$
including the Coulomb case $n=1$ for any sign of $\alpha$. In the next subsection we discuss some probabilistic aspects of the kinetic theory of gases. 2.2. Distribution functions and Liouville’s equationWe introduce an important notion of one-particle distribution function $f(x,v,t)$ (the words ‘one- particle’ are usually omitted below for the sake of brevity). The physical meaning of this function is as follows: the average number of particles in any measurable set $\Delta \subseteq \mathbb{R}^3 \times \mathbb{R}^3$ is given by the equality In other words, $f(x,v,t)$ is the density of the number of particles in the phase space. Usually, we assume that the initial data are fixed. How can we find the distribution function $f(x,v,t)$ for $t>0$? This is, in a sense, the main problem of kinetic theory. It can be shown that for some special physical systems, like rarefied gases, the temporal evolution of distribution function $f(x,v,t)$ is described by a so-called ‘kinetic equation’ where $A(f)$ is a non-linear operator acting on $f$. We usually assume that the initial value problem (2.7), (2.8) has a unique solution $f(x,v,t)$ on some time interval $0 \leqslant t \leqslant T$.Let us consider the simplest kinetic equation connected with the system (2.1). Omitting the index $i=1$ we obtain the equations of free motion Note that the motion of one particle can also be described by the distribution function $f(x,v,t)$, which has in this case the meaning of a probability density, provided that
$$
\begin{equation}
\int_{\mathbb{R}^3\times\mathbb{R}^3} dx\,dv\,f^{(0)}(x,v)=1
\end{equation}
\tag{2.10}
$$
in the notation of (2.7). The solution of equations (2.9) is obvious: Therefore, $f(x,v,t)$ satisfying conditions (2.7) reads
$$
\begin{equation}
f(x,v,t)=f^{(0)}(x-vt, v)=\exp(-tv\cdot \partial_x)f^{(0)}(x,v).
\end{equation}
\tag{2.11}
$$
Here and below a dot denotes the scalar product in $\mathbb{R}^3$. We can check by differentiation that This is the simplest kinetic equation. Note that the kinetic equations (2.12) are exactly equivalent to the dynamical equations (2.9). The probabilistic description is caused only by an uncertainty in the initial conditions. Let us now extend these arguments to the case of $N$ non-interacting particles. We consider equations (2.1) for $\Phi(r)\equiv 0$ and obtain
$$
\begin{equation}
\begin{gathered} \, \partial_t x_i=v_i, \quad \partial_t v_i=0; \\ x_i(0)=x_i^{(0)},\quad v_i(0)=v_i^{(0)},\qquad i=1,\dots, N. \end{gathered}
\end{equation}
\tag{2.13}
$$
Thus we have $N$ independent vector equations for each particle. It is natural to introduce an $N$-particle distribution function
$$
\begin{equation*}
F_N(z_1,\dots,z_N;t),\qquad z_i=(x_i,v_i),\quad 1 \leqslant i \leqslant N,
\end{equation*}
\notag
$$
with the meaning of a probability density in the $N$-particle phase space $\mathbb{R}^{6N}$. The initial condition reads
$$
\begin{equation}
F_N\big|_{t=0}=F_N^{(0)}(z_1,\dots,z_N),\quad \int_{\mathbb{R}^{6}\times \dots \times \mathbb{R}^{6}}\,dz_1\cdots dz_n\, F_N^{(0)}(z_1,\dots,z_N)=1.
\end{equation}
\tag{2.14}
$$
Remark 2.1. Here and below we use notation like $F_N(z_1,\dots,z_N)$ (with upper- case $F$) for various multi-particle distribution functions which have the meaning of a probability density. These functions are always normalized by $1$ in the whole phase space. The notation like $f_N(z_1,\dots,z_N)$ will be used for a slightly different class of functions related to equality (2.6). The difference disappears in the trivial case $N=1$. Then it is easy to see that
$$
\begin{equation*}
\begin{gathered} \, F_N(z_1,\dots,z_N;t)=F_N^{(0)}[z_1(t),\dots,z_N(t)], \\ z_i(t)=(x_i-v_it,v_i), \quad i=1,\dots,N. \end{gathered}
\end{equation*}
\notag
$$
Note that
$$
\begin{equation}
\biggl(\partial_t+\sum_{i=1}^N v\cdot \partial_{x_i}\biggr) F_N(x_1,v_1,\dots,x_N,v_N;t)=0.
\end{equation}
\tag{2.15}
$$
Thus we obtain the simplest version of Liouville’s equation. We assume that all particles are identical and independently distributed at $t=0$, that is, Then a similar factorization holds for all $t>0$:
$$
\begin{equation}
F_N(z_1,\dots,z_N;t)=\prod_{i=1}^N f(z_i,t)= \prod_{i=1}^N f^{(0)}(x_i-v_i t,v_i).
\end{equation}
\tag{2.17}
$$
This property is known as the propagation of chaos [39]. It is self-evident for non-interacting particles, but it can also be proved as an asymptotic property of more complex multi-particle systems. What changes if we consider the Hamiltonian system (2.1) for a non-zero potential $\Phi(r)\ne 0$? Then we can still use the $N$-particle distribution function $F_N=\{z_1,\dots,z_N;t)$. We shall see below that the equation for $F_N$ reads
$$
\begin{equation}
\begin{gathered} \, \biggl[\partial_t+\sum_{i=1}^N\biggl(v_i\cdot\partial_{x_i}- \frac{\partial \Phi_N}{\partial x_i}\cdot \partial_{v_i}\biggr)\biggr] F_N(x_1,v_1;\dots;x_N,v_N;t)=0, \\ \Phi_N=\sum_{1 \leqslant i <j \leqslant N}\Phi(|x_i-x_j|),\qquad i=1,\dots,N. \end{gathered}
\end{equation}
\tag{2.18}
$$
This is the famous Liouville equation [41]. It can be derived very easily. For a moment we simplify our notation in the following way:
$$
\begin{equation*}
\begin{gathered} \, z=\{z_1,\dots,z_N\} \in \mathbb{R}^{6N},\qquad z_i=(x_i,v_i); \\ w(z)=\{w_1,\dots,w_N\}\in \mathbb{R}^{6N},\qquad w_i=\biggl(\frac{\partial H}{\partial v_i}\,, -\frac{\partial H}{\partial x_i}\biggr); \\ H=H_N=\frac{1}{2}\sum_{i=1}^N |v_i|^2+ \sum_{1 \leqslant i <j \leqslant N}\Phi(|x_i-x_j|); \\ F(z;t)=F_N(z_1,\dots,z_N;t). \end{gathered}
\end{equation*}
\notag
$$
Then we can treat $F(z,t)$ as a density of a fluid in $\mathbb{R}^{6N}$, which moves in accordance with a dynamical system where $w(z)$ is assumed to be a ‘nice’ function. Then the density $F(z;t)$ satisfies the continuity equation where $\operatorname{div}_z$ denotes divergence with respect to $z$. Simple calculations yield
$$
\begin{equation*}
\operatorname{div}_z w(z)=\sum_{i=1}^N\biggl(\frac{\partial}{\partial x_i}\cdot \frac{\partial H}{\partial v_i}-\frac{\partial}{\partial v_i}\cdot \frac{\partial H}{\partial x_i}\biggr)=0.
\end{equation*}
\notag
$$
Hence we obtain
$$
\begin{equation}
F_t+\sum_{i=1}^{N}\biggl(\frac{\partial H}{\partial v_i}\cdot \frac{\partial F}{\partial x_i}-\frac{\partial H}{\partial x_i}\cdot \frac{\partial F}{\partial v_i}\biggr)=0,
\end{equation}
\tag{2.21}
$$
which is Liouville’s equation (2.18) in slightly different notation. The sum in (2.21) is usually called the Poisson bracket $\{F,H\}$ (only for particles with unit mass $m= 1$); see [41]. Liouville’s equation is very important because it allows one (at least formally) to build a bridge between $N$-particle dynamics, kinetic theory, and hydrodynamics. 2.3. The BBGKY-hierarchyIn § 2.2 the function $f_N(z_1,\dots,z_N;t)$ was assumed to be integrable over the whole phase space ${\mathbb{R}}^{6N}$. This assumption is sometimes too strong. For example, it does not allow one to consider translationally invariant (in the physical space ${\mathbb{R}}^3$) systems with initial data like $f_2^{(0)}(x_1-x_2;v_1,v_2)$. Therefore, it is more convenient to consider an $N$-particle system confined to a bounded domain $\Omega \subset {\mathbb{R}^3}$, say, a box or a sphere with ‘large’ diameter $L$. The volume of $\Omega$ is denoted by $|\Omega|$. Then the phase state of the $i$th particle is described by the point
$$
\begin{equation}
z_i=(x_i,v_i) \in \Omega \times {\mathbb{R}^3},\qquad i=1,\dots,N.
\end{equation}
\tag{2.22}
$$
We can assume for simplicity that the walls of $\Omega$ are specularly reflecting. Another possibility is to assume that $\Omega$ is a periodic box. What is important is to keep unchanged the total number $N$ of particles inside $\Omega$. For the sake of brevity we use below the symbolic notation
$$
\begin{equation}
\begin{gathered} \, F_N(x_1,v_1;\dots;x_N,v_N;t)=F_N(1,\dots,N;t), \\ \int_\Omega dx_i\int_{\mathbb{R}^3}\,dv_i\ldots= \int_G di\dots\,, \qquad G=\Omega\times {\mathbb{R}^3},\quad i=1,\dots,N, \end{gathered}
\end{equation}
\tag{2.23}
$$
where $F_N$ is normalized by the equality
$$
\begin{equation}
\int_{G_N} d1\,d2 \cdots\,dN\,F_N(1,\dots,N;t)=1,\qquad G_N=G^N,
\end{equation}
\tag{2.24}
$$
and satisfies Liouville’s equation (2.8):
$$
\begin{equation}
\begin{gathered} \, \biggl( \partial_t+\sum_{i=1}^N A_i-\sum_{1 \leqslant i< j \leqslant N} B_{ij}\biggr)F_N=0, \\ A_i=v_i\cdot\partial_{x_i},\qquad i=1,\dots,N, \\ B_{ij}=\frac{\partial \Phi(|x_i-x_j|)}{\partial x_i}\cdot (\partial_{v_i}-\partial_{v_j}),\qquad i,j=1,\dots,N,\quad i \ne j. \end{gathered}
\end{equation}
\tag{2.25}
$$
The function $F_N(1,\dots,N;t)$ is assumed to be invariant under permutations of its arguments $(1,\dots,N)$ because all $N$ particles are identical. This property is preserved by Liouville’s equation, provided that it is valid at $t=0$. We introduce $k$-particle probability distributions by the equalities
$$
\begin{equation}
F_k^{\kern1pt(N)}(1,\dots,k)=\int_{G_{N-k}}\,d(k+1)\cdots dN\,F_N(1,\dots,N),\qquad 1\leqslant k \leqslant N-1.
\end{equation}
\tag{2.26}
$$
Here and below the argument $t$ is omitted in such cases when this does not cause any confusion. The next step is to obtain a set of evolution equations for $F_k^{\kern1pt(N)}(1,\dots,k;t)$. We take any $1\leqslant k \leqslant N-1$ and integrate equation (2.25) with respect to $d(k+1)\cdots\,dN$. The result reads
$$
\begin{equation}
(\partial_t+L_k)F_k^{\kern1pt(N)}=\Gamma_k^{\kern1pt(N)},\qquad L_k=\sum_{i=1}^k A_i-\sum_{1 \leqslant i < j \leqslant k}B_{ij},
\end{equation}
\tag{2.27}
$$
$$
\begin{equation}
\nonumber \Gamma_k^{\kern1pt(N)}= \int_{G_{N-k}}d(k+1)\cdots\,dN \biggl[-\sum_{i=k+1}^N A_i+\sum_{i=1}^k\,\sum_{j=k+1}^N B_{ij}
\end{equation}
\notag
$$
$$
\begin{equation}
\qquad\qquad\qquad\qquad\qquad\qquad\qquad\quad+\sum_{i=k+1}^{N-1}\,\sum_{j=i+1}^N B_{ij}\biggr]F_N(1,\dots,N).
\end{equation}
\tag{2.28}
$$
We assume that the boundary conditions guarantee that for each $1 \leqslant i \leqslant N$
$$
\begin{equation}
\int_{G}di\,(v_i \cdot \partial_{x_i})\,F_N(1,\dots,N;t)=0.
\end{equation}
\tag{2.29}
$$
Then the first sum in (2.28) disappears. The third sum in (2.28) also disappears under the natural assumption that
$$
\begin{equation}
F_N(1,\dots,N;t) \to 0 \quad \text{if}\ \ |v_i| \to \infty \ \ \text{for some}\ \ 1 \leqslant i \leqslant N.
\end{equation}
\tag{2.30}
$$
Then we obtain
$$
\begin{equation*}
\Gamma_k^{\kern1pt(N)}=\sum_{i=1}^k\Gamma_{ik}^{\kern1pt(N)}, \qquad \Gamma_{ik}^{\kern1pt(N)}=\sum_{j=k+1}^N \,\int_{G_{N-k}} d(k+1)\cdots\,dN \,B_{ij}\,F_N(1,\dots,N).
\end{equation*}
\notag
$$
Let us consider the first term of the sum for $\Gamma_{ik}^{\kern1pt(N)}$. Then $j=k+1$ and we obtain
$$
\begin{equation*}
\begin{aligned} \, &\int_{G} d(k+1) B_{i,k+1}\int_{G_{N-k-1}}d(k+2)\cdots\,dN\,F_N(1,\dots,N) \\ &\qquad=\int_{G}d(k+1)\,B_{i,k+1}\,F_{k+1}^{\kern1pt(N)}(1,\dots,k+1),\quad 1\leqslant k \leqslant N-1,\qquad F_N^{\kern1pt(N)}=F_N. \end{aligned}
\end{equation*}
\notag
$$
If we take another value of $j$ in the sum for $\Gamma_{ik}^{\kern1pt(N)}$, then the result is the same. This follows from the symmetry of $F_N(1,\dots,N)$ with respect to permutations. Moreover, the operator $B_{i,k+1}$ from (2.25) can be replaced by
$$
\begin{equation*}
\widetilde{B}_{i,k+1}=\frac{\partial\Phi(|x_i-x_{k+1}|)}{\partial x_i} \cdot\frac{\partial}{\partial v_i}
\end{equation*}
\notag
$$
because of conditions (2.30). Therefore, we obtain
$$
\begin{equation*}
\Gamma_k^{\kern1pt(N)}=(N-k)\sum_{i=1}^k \int_{G} d(k+1)\,\widetilde{B}_{i,k+1}\, F_{k+1}^{\kern1pt(N)}(1,\dots,k+1).
\end{equation*}
\notag
$$
We substitute this formula into equations (2.27) and obtain the following set of equations:
$$
\begin{equation}
(\partial_t+L_k) F_k^{\kern1pt(N)}=(N-k)C_{k+1}F_{k+1}^{\kern1pt(N)},\qquad 1 \leqslant k \leqslant N-1,
\end{equation}
\tag{2.31}
$$
where the $F_k^{\kern1pt(N)}(1,\dots,k)$ are given in (2.26),
$$
\begin{equation}
\begin{gathered} \, F_N^{\kern1pt(N)}(1,\dots,N)=F_N(1,\dots,N), \\ L_k=\sum_{i=1}^k v_i \cdot \partial_{ x_i}- \sum_{1 \leqslant i < j \leqslant k}\,\frac{\Phi(|x_i-x_{j}|)}{\partial x_i} \cdot(\partial_{v_i}-\partial_{ v_j}), \\ C_{k+1}=\sum_{i=1}^k C_{i,k+1}, \\ C_{i,k+1}\,F_{k+1}^{\kern1pt(N)}=\int_{G} d(k+1)\, \frac{\Phi(|x_i-x_{k+1}|)}{\partial x_i}\cdot \partial_{v_i} F_{k+1}^{\kern1pt(N)}(1,\dots,k+1). \end{gathered}
\end{equation}
\tag{2.32}
$$
The set of equations (2.31), (2.32) is called (provided that we ignore its trivial modifications) the BBGKY-hierarchy. The acronym BBGKY stands for Bogolyubov, Born, Green, Kirkwood, and Yvon, the names of the physicists who introduced independently this system of equations (see, for example, [9] and [15]). Equation (2.31) for $k=N$ also makes sense if we set $F_N^{\kern1pt(N)}=F_N$ and $F_N^{(N+1)}=0$. Then it becomes Liouville’s equation for $F_N$. Note that equations (2.31) (the BBGKY- hierarchy for $N$-particle dynamical system (2.1) in the box $\Omega$) can formally be considered as exact. The only relevant assumption is equality (2.29), related to interactions with boundaries of the box $\Omega$. This equality is fulfilled, in particular, for a periodic box or a box with specularly reflecting walls. The BBGKY-hierarchy is a starting point for all classical kinetic equations. For example, assume that the interaction between particles is weak and replace $\Phi(r)$ by $\varepsilon\bar{\Phi}(r)$ in equations (2.31) in the following form:
$$
\begin{equation*}
\biggl(\partial_t+\sum_{i=1}^N A_i- \varepsilon \sum_{1 \leqslant i <j \leqslant N} B_{i j}\biggr) F_{k+1}^{\kern1pt(N)}=\varepsilon(N-k)C_{k+1} F_{k+1}^{\kern1pt(N)},\qquad 1 \leqslant k \leqslant N-1,
\end{equation*}
\notag
$$
in the notation of (2.25), (2.32). Then it is natural to consider the limit
$$
\begin{equation*}
N \to \infty,\quad \varepsilon \to 0,\quad N\varepsilon=\mathrm{const}
\end{equation*}
\notag
$$
and to assume that $F_{k}^{\kern1pt(N)} \to F_k$ in that limit for all $k=1,2,\dots$ . Then we formally obtain the following infinite set of limiting equations for $F_{k}$:
$$
\begin{equation*}
\begin{aligned} \, &\biggl(\partial_t+\sum_{i=1}^k \,v_i\cdot \partial_{x_i}\biggr) F_{k}(z_1,\dots,z_k;t) \\ &\qquad=N\varepsilon \sum_{i=1}^k \int_G dz_{k+1}\, \frac{\Phi(|x_i-x_{k+1}|)}{\partial x_i}\cdot \partial_{v_i}F_{k+1}(z_1,\dots,z_{k+1};t); \\ &z_i=(x_i,v_i) \in G,\quad dz_i=dx_i\,dv_i,\quad i=1,\dots,k+1, \quad k=1,2,\dots\,. \end{aligned}
\end{equation*}
\notag
$$
It is easy to verify that these equations admit a class of solutions in the factorized form:
$$
\begin{equation*}
F_k(z_1,\dots,z_k;t)=\prod_{i=1}^k F(z_i,t),\qquad k=1,2,\dots,
\end{equation*}
\notag
$$
where $F(x,v,t)$ satisfies Vlasov’s equation [59]
$$
\begin{equation*}
F_t+v \cdot F_x-N\varepsilon F_v \cdot \partial_x\int_\Omega dy\,dw \,\Phi(|x-y|)F(y,w,t)=0.
\end{equation*}
\notag
$$
The most important applications of Vlasov’s equations are related to the Coulomb potential $\Phi(r)=\alpha/r$ (gravitational or electrostatic forces). On the other hand, in the important case of potentials $\Phi(r)$ with compact support (strong interaction at short distances) it is more convenient to use a different approach. This approach is considered in the next subsection. 2.4. Hard spheres and the Boltzmann–Grad limitWe begin with the case of $N$ hard spheres of diameter $d$. Then the Liouville equation (2.18) cannot be used directly because the potential $\Phi_{\rm HS}(r)$ (2.4) is too singular. The function $F_N(x_1,v_1;\dots;x_N,v_N;t)$ is defined in this case by the free flow equation
$$
\begin{equation}
\partial_t F_N+\sum_{i=1}^N v_i \cdot \partial_{x_i}F_N=0,
\end{equation}
\tag{2.33}
$$
valid in the domain
$$
\begin{equation}
B_N=\{x_i \in \mathbb{R}^3\colon |x_i-x_j| > d; \, i,j=1,\dots,N, \,i\ne j\}.
\end{equation}
\tag{2.34}
$$
We add to this equation the boundary conditions on each of $N(N-1)/2$ boundary surfaces $|x_i-x_j|=d$ ($i,j=1,\dots,N$, $i\ne j$) in $\mathbb{R}^{3N}$. Taking, for example, $i=1$ and $j=2$, and assuming the specular reflection law we obtain
$$
\begin{equation}
F_N(x_1,v_1;x_2,v_2;\dots;t)\big|_{|x_1-x_2|=d,\,x\cdot u>0}= F_N(x_1,v_1';x_2,v_2';\dots;t),
\end{equation}
\tag{2.35}
$$
where
$$
\begin{equation}
\begin{gathered} \, x=x_1-x_2= d\,n, \quad n \in S^2, \qquad u=v_1-v_2, \\ v_1'=v_1-n (u \cdot n),\quad v_2'=v_2+n (u \cdot n). \end{gathered}
\end{equation}
\tag{2.36}
$$
The surfaces that correspond to multiple collisions of $k \geqslant 3$ particles are described by at least $k-1 \geqslant 2$ equalities. For example, for $k=3$ we need to satisfy simultaneously two conditions like $|x_1-x_2|=d$ and $|x_2-x_3|=d$. These surfaces have a zero measure as compared with the case $k=2$ of pair collisions. Therefore, we ignore multiple collisions. For brevity we assume that
$$
\begin{equation}
F_N \in L(G_N), \qquad G_N=G^N,\quad G=\mathbb{R}^3\times \mathbb{R}^3,
\end{equation}
\tag{2.37}
$$
and if $|x_i-x_j| < d$ for at least one pair of indices $1 \leqslant i<j \leqslant N$. Here and below we use the symbolic notation (2.23) from § 2.3 when it does not cause any confusion. We also assume that $F_N \geqslant 0$ is the probability density in $G_N$ with usual normalization condition (2.24). Our aim is to construct the equation for the one-particle probability density $F_1^{\kern1pt(N)}(1)$ defined in (2.26) for $k=1$. Note that this function can equivalently be defined by the equality
$$
\begin{equation}
\begin{gathered} \, F_1^{\kern1pt(N)}(1)=\int_{G_{N-1}}d2 \cdots\,dN \,\Psi(1|2,\dots,N)F_N(1,\dots,N), \\ \Psi(1|2,\dots,N)=\prod_{k=2}^N \eta[d^2-|x_1-x_k|^2], \end{gathered}
\end{equation}
\tag{2.38}
$$
where $\eta(y)$ is the step function
$$
\begin{equation}
\eta(y)=\begin{cases} 1 & \text{for}\ y>0, \\ 0 & \text{otherwise}. \end{cases}
\end{equation}
\tag{2.39}
$$
We multiply equation (2.33) by $\Psi(1|2,\dots,N)$ and integrate over $G_{N-1}$. The result reads
$$
\begin{equation}
\partial_t F_1^{\kern1pt(N)}+I_1+\sum_{k=2}^N I_k =0,
\end{equation}
\tag{2.40}
$$
where
$$
\begin{equation}
\begin{gathered} \, I_k=\int_{G_{N-1}} d2\cdots\,dN\,\Psi(1|2,\dots,N)A_k F_N(1,\dots,N), \\ A_k=v_k \cdot \partial_{x_k}, \quad k=1,\dots,N. \end{gathered}
\end{equation}
\tag{2.41}
$$
We separated the term with $k=1$ in (2.40) because all other terms in the sum are equal to $I_2$. Indeed, we always assume that all particles are identical and therefore $F_N(1,\dots,N)$ is symmetric with respect to permutations of the arguments $(2,\dots,N)$. Hence $I_k=I_2$ for any $3 \leqslant k \leqslant N$, and therefore The integral $I_2$ can be written as
$$
\begin{equation}
I_2=\int_{|x_1-x_2|>d} dx_2\,\operatorname{div}_{x_2}\int_{\mathbb{R}^3} dv_2 \,v_2 \,\widetilde{F}_2^{\kern1pt(N)}(x_1,v_1;x_2,v_2),
\end{equation}
\tag{2.43}
$$
where
$$
\begin{equation}
\begin{aligned} \, \nonumber \widetilde{F}_2^{\kern1pt(N)}(x_1,v_1;x_2,v_2)&= \widetilde{F}_2^{\kern1pt(N)}(1,2) \\ &=\int_{G_{N-2}}d3 \cdots\,dN\,F_N(1,2,\dots,N)\prod_{k=3}^N \eta[d^2-|x_1-x_k|^2]. \end{aligned}
\end{equation}
\tag{2.44}
$$
We use the notation $\widetilde{F}_2^{\kern1pt(N)}$ because this function coincides formally with $F_2^{\kern1pt(N)}$ in (2.26) only for $d=0$. We apply Gauss’s theorem to the integral (2.43) and obtain after simple transformations
$$
\begin{equation*}
\begin{aligned} \, I_2&=-\int_{|y|>d}dy\,\operatorname{div}_y\int_{\mathbb{R}^3}dv_2\, v_2 \widetilde{F}_2^{\kern1pt(N)}(x_1,v_1;x_1-y,v_2) \\ &=d^2\int_{\mathbb{R}^3 \times S^2}dv_2\,dn\,(v_2 \cdot n) \widetilde{F}_2^{\kern1pt(N)}(x_1,v_1;x_1-dn,v_2), \end{aligned}
\end{equation*}
\notag
$$
where $n$ denotes the outward unit normal vector to the unit sphere $S^2 $. It remains to evaluate the integral $I_1$ in (2.41). We note that
$$
\begin{equation*}
\Psi(1|2,\dots,N)A_1 F_N(1,\dots,N)= A_1\Psi F_N -F_N A_1 \Psi,\qquad A_1=v_1 \cdot \partial_{x_i}.
\end{equation*}
\notag
$$
Since $\eta'(y)=\delta(y)$, where $\delta(y) $ denotes the Dirac delta function, we obtain
$$
\begin{equation*}
\begin{aligned} \, A\Psi&=v_1 \cdot \partial_{x_1}\prod_{i=2}^N\eta[d^2-|x_1-x_i|^2] \\ &=2\sum_{i=2}^N\,[v_1 \cdot(x_1-x_i)]\delta[|x_1-x_i|^2-d^2] \prod_{\substack{2 \leqslant j \leqslant N\\ j \ne i}}\eta[d^2-|x_1-x_j|^2]. \end{aligned}
\end{equation*}
\notag
$$
Then we perform integration in (2.41) and use the symmetry of $F_N$ and $\Psi$ again. The result reads
$$
\begin{equation*}
\begin{aligned} \, I_1&= v \cdot \partial_{x_1} F_1^{\kern1pt(N)}(x_1,v_1) -2(N-1)\int_{\mathbb{R}^3\times \mathbb{R}^3}dx_2\,dv_2\, \delta\bigl[|x_1-x_2|^2-d^2\bigr]v_1 \\ &\qquad \times (x_1-x_2) \widetilde{F}_2^{\kern1pt(N)} (x_1,v_1;x_2,v_2) \end{aligned}
\end{equation*}
\notag
$$
in the notation of (2.44). Note that
$$
\begin{equation*}
2\int_{\mathbb{R}^3}dy\,\delta[(x-y)^2-d^2]F(y)= d \int_{S^2} dn\,F(x-dn).
\end{equation*}
\notag
$$
Therefore, we obtain equation (2.40) in the following form:
$$
\begin{equation}
\begin{aligned} \, \nonumber &(\partial_t+v_1 \cdot \partial_{x_1}) F_1^{\kern1pt(N)}(x_1,v_1)=Q^{\kern1pt(N)} \\ &\qquad=(N-1)\, d^2\int_{ \mathbb{R}^3\times S^2}dv_2\,dn \, [(v_1-v_2) \cdot n]\widetilde{F}_2^{\kern1pt(N)}(x_1,v_1;x_1-dn,v_2), \end{aligned}
\end{equation}
\tag{2.45}
$$
where $\widetilde{F}_2^{\kern1pt(N)}(x_1,v_1;x_2,v_2)$ was defined in (2.44). We can split the integral over $S^2$ into two parts in the following way:
$$
\begin{equation*}
\begin{gathered} \, \int_{S^2} dn\,(u \cdot n)\Psi(n)=\int_{S_+^2} dn\,|u \cdot n|\Psi(n)- \int_{S_-^2} dn\,|u \cdot n| \Psi(n)\,, \\ S_+^2=\{n \in S^2\colon u \cdot n >0\},\qquad S_-^2=\{n \in S^2\colon u \cdot n <0\}, \end{gathered}
\end{equation*}
\notag
$$
where $u=v_1-v_2$ and $\Psi(n)$ is an arbitrary integrable function. It is clear from (2.35) and (2.36) that the function $\Psi(n)=\widetilde{F}_2^{\kern1pt(N)}(x_1,v_1; x_1-dn,v_2)$ in the integral over $S_+^2$ can be expressed in terms of $\Psi(n)$ in the integral over $S^2_-$. Then we obtain
$$
\begin{equation}
\begin{aligned} \, Q^{\kern1pt(N)}&=(N-1)\,d^2 \,\int_{ \mathbb{R}^3\times S_+^2} dv_2\,dn\,|u \cdot n| \\ &\qquad\times\bigl[\widetilde{F}_2^{\kern1pt(N)}( x_1,v_1';x_1-dn,v_2')- \widetilde{F}_2^{\kern1pt(N)}(x_1,v_1;x_1+dn,v_2)\bigr], \\ v_1'&=v_1-n(u \cdot n), \qquad u=v_1-v_2,\qquad v_2'=v_2+n(u \cdot n). \end{aligned}
\end{equation}
\tag{2.46}
$$
Note that equations (2.45) and (2.46) are formally exact for hard spheres. To our knowledge, they were published for the first time by Grad no later than 1957 [34]. These equations are very important as a starting point for the mathematically rigorous derivation of the Boltzmann equation. For our goals it is sufficient to introduce the ‘chaotic’ initial data
$$
\begin{equation}
\begin{gathered} \, F_N(1,2,\dots,N)\big|_{t=0}=c_N\prod_{k=1}^N F_0(k) \prod_{1 \leqslant i < j \leqslant N}\eta[d^2-|x_i-x_j|^2], \\ \int_{G}d1 \, F_0(1)=1, \end{gathered}
\end{equation}
\tag{2.47}
$$
where $c_N$ is the normalization constant, and to consider the formal limit of equations (2.45) and (2.46) under the conditions that
$$
\begin{equation}
N \to \infty, \quad d \to 0, \quad Nd^2=\mathrm{const}.
\end{equation}
\tag{2.48}
$$
It is the so-called Boltzmann–Grad limit [20], [34]. To be more precise, we assume that
$$
\begin{equation}
F_1^{\kern1pt(N)}(x_1,v_1,t) \to F(x,v,t),\qquad F_1^{\kern1pt(N)}(x_1,v_1;x_2,v_2;t) \to \prod_{i=1}^2 F(x_i,v_i,t)
\end{equation}
\tag{2.49}
$$
under conditions (2.48). Then from (2.45) and (2.46) we formally obtain the Boltzmann equation for hard spheres. In the simplified notation $x=x_1$, $v=v_1$, $w=v_2$ this equation reads
$$
\begin{equation}
(\partial_t+v \cdot \partial_{x})F(x,v,t)=(Nd^2) \widetilde{Q}(F,F),
\end{equation}
\tag{2.50}
$$
where
$$
\begin{equation}
\begin{gathered} \, \widetilde{Q}(F,F)=\frac{1}{2}\int_{\mathbb{R}^3\times S^2}dw\,dn\, |u \cdot n|\bigl[F(x,v',t)F(x,w',t)-F(x,v,t)F(x,w,t)\bigr], \\ u=v-w,\quad n \in S^2, \quad v'=v-(u \cdot n) n, \quad w'=w+(u \cdot n) \end{gathered}
\end{equation}
\tag{2.51}
$$
(the domain of integration $S_+^2$ in (2.46) is extended to the whole unit sphere $S^2$ in an obvious way). The limiting initial conditions follows formally from (2.47):
$$
\begin{equation}
F(x,v)\big|_{t=0}=F_0(x,v),\quad \int_G dx\,dv\,F_0(x,v)=1.
\end{equation}
\tag{2.52}
$$
We consider in this section the case of the whole space $G=\mathbb{R}^3 \times \mathbb{R}^3$ in order to simplify formal calculations. In fact, the same equations (2.50)–(2.52) can be derived rigorously from $N$-particle dynamics in the case when all particles are confined to a bounded domain $\Omega$ with reflecting walls (see [21] and the references there for details). In our formal derivation of the Boltzmann equation in this section we have followed Grad’s scheme from [34]. It should be pointed out that the first mathematically rigorous proof of the validity of the Boltzmann equation for hard spheres was given by Lanford [42] in 1975. The deep and rigorous presentation of the validation problem for the Boltzmann equation can be found in the book [21]. 2.5. The Boltzmann equation for hard spheres and other potentialsThe classical Boltzmann equation is usually considered not for the probability density $F(x,v,t)$, but for the distribution function The equation for $f(x,v,t)$ reads
$$
\begin{equation}
\begin{gathered} \, f_t+v \cdot f_x=Q(f,f)=\frac{d^2}{2}\int_{\mathbb{R}^3\times S^2} dw\,dn \,|u \cdot n|\,[f(v') f(w')-f(v) f(w)], \\ u=v-w,\quad n \in S^2, \quad v'=v-(u \cdot n) n, \quad w'=w+(u \cdot n), \end{gathered}
\end{equation}
\tag{2.54}
$$
where the irrelevant arguments $(x,t)$ of $f(x,v,t)$ have been omitted in the so-called Boltzmann collision integral $Q(f,f)$. Note that $Q(f,f)$ is an operator quadratic with respect to $f$ and acting only on the variable $v \in \mathbb{R}^3$. The advantage of equation (2.54) is that it does not contain the number of particles $N$. The Boltzmann collision integral can be presented in different forms. In particular, the following form of $Q(f,f)$ is very useful:
$$
\begin{equation}
\begin{gathered} \, Q(f,f)=\frac{d^2}{4}\int_{\mathbb{R}^3\times S^2} dw\,d\omega\,|u|\,[f(v')f(w')-f(v) f(w)], \\ \nonumber u=v-w,\quad \omega \in S^2, \quad v'=\frac{1}{2}(v+w+|u|\omega), \quad w'=\frac{1}{2}(v+w-|u|\omega). \end{gathered}
\end{equation}
\tag{2.55}
$$
In order to prove that this integral coincides with $Q(f,f)$ in (2.54) we consider the simpler integral
$$
\begin{equation}
I(F)=\int_{\mathbb{R}^3}\,dk\,\delta\biggl(k \cdot u+ \frac{|k|^2}{2}\biggr)F(k),
\end{equation}
\tag{2.56}
$$
where $u \in \mathbb{R}^3$ and $F(u)$ is a continuous function. Then the following lemma can easily be proved. Proof. We evaluate the integral $I(F)$ in the spherical coordinates with polar axis directed along $u \in \mathbb{R}^3$. Thus we set $k=rn$, $n \in S^2$, in (2.56) and obtain
$$
\begin{equation*}
I(F)=\int_{0}^\infty\,dr\, r^2\int_{S^2}\,dn \, \delta\biggl[rn\cdot u+\frac{r^2}{2}\biggr]F(rn).
\end{equation*}
\notag
$$
Since
$$
\begin{equation}
\delta(\alpha x)=\frac{\delta(x)}{\alpha}\,,\qquad \alpha > 0,\quad x \in \mathbb{R},
\end{equation}
\tag{2.58}
$$
we change the order of integration and obtain
$$
\begin{equation*}
\begin{aligned} \, I(F)&=\int_{S^2}\,dn\int_{0}^{\infty}\,dr\, r\delta \biggl[\frac{1}{2}(2 n \cdot u+r)\biggr]F(rn) \\ &=4\int_{S^2_-}\,dn \,|n \cdot u|F[-2(u \cdot n)n],\quad S^2_-=\{n \in S^2\colon u \cdot n < 0\}. \end{aligned}
\end{equation*}
\notag
$$
The integrand is an even function of $n \in S^2$, and thus the first equality (2.57) follows. The second equality is based on the change of variables $k=\widetilde{k}-u$ in the integral (2.56). Then we obtain and evaluate this integral in the same way as above. This completes the proof of the lemma. Now we can prove the transformation of $Q(f,f)$ from (2.54) into (2.55). We consider (2.54) and set
$$
\begin{equation}
F(k)=f\biggl(v+\frac{k}{2}\biggr)f\biggl(w-\frac{k}{2}\biggr)-f(v)f(w),
\end{equation}
\tag{2.59}
$$
considering $v$ and $w$ as fixed parameters. Then from (2.54) we obtain
$$
\begin{equation*}
Q(f,f)=\frac{d^2}{2}\int_{\mathbb{R}^3 \times S^2} dw \, dn \, |u \cdot n|F[-2(u \cdot n)n].
\end{equation*}
\notag
$$
It remains to use identity (2.57) and obtain the following result:
$$
\begin{equation*}
Q(f,f)= \frac{d^2}{4}\int_{\mathbb{R}^3 \times S^2} dw \, d \omega\, |u| F(|u|\omega-u)
\end{equation*}
\notag
$$
in the notation of (2.59). It is easy to check that this formula for $Q(f,f)$ coincides with (2.55). Hence the equivalence of (2.54) and (2.55) is proved. Note also that the same identity (2.57) leads to the third useful representation of the collision integral for hard spheres:
$$
\begin{equation}
\begin{aligned} \, \nonumber Q(f,f)&=\frac{d^2}{4} \int_{\mathbb{R}^3\times \mathbb{R}^3} dw\,dk \, \delta\biggl(k \cdot u+\frac{|k|^2}{2}\biggr) \\ &\qquad\times\biggl[f\biggl(v+\frac{k}{2}\biggr)f\biggl(w-\frac{k}{2}\biggr)- f(v)f(w)\biggr]. \end{aligned}
\end{equation}
\tag{2.60}
$$
The physical meaning of the Boltzmann equation can be understood better by considering $Q(f,f)$ in the form (2.55). We set
$$
\begin{equation}
\langle f,\psi\rangle=\int_{\mathbb{R}^3}\,dv\, f(v)\psi(v),
\end{equation}
\tag{2.61}
$$
where $\psi(v)$ is an arbitrary function of velocity $v \in \mathbb{R}^3$ for which the integral exists. Then from (2.54) we formally obtain
$$
\begin{equation}
\partial_t\langle f,\psi \rangle + \partial_x \cdot \langle f, v\psi \rangle=\langle \psi,Q(f,f) \rangle.
\end{equation}
\tag{2.62}
$$
The right-hand side with $Q(f,f)$ from (2.55) reads
$$
\begin{equation*}
\langle \psi,Q(f,f)\rangle=\frac{d^2}{4} \int_{\mathbb{R}^3 \times\mathbb{R}^3\times S^2}\,dv\,dw\,d\omega\, |u|\psi(v)[f(v')f(w')-f(v)f(w)]
\end{equation*}
\notag
$$
in the notation of (2.55). We denote the centre of mass variables by $U$ and $u$, where
$$
\begin{equation}
U=\frac{v+w}{2}\,, \ \ u=v-w\quad\Longleftrightarrow\quad v=U+\frac{u}{2}\,, \ w=U-\frac{u}{2}\,.
\end{equation}
\tag{2.63}
$$
Hence $dv\,dw=dU\,du$. Therefore, we obtain
$$
\begin{equation*}
\begin{gathered} \, \langle \psi,Q(f,f)\rangle=\frac{d^2}{4} \int_{\mathbb{R}^3 \times\mathbb{R}^3 \times S^2} dU\,du\,d\omega\, |u|\psi\biggl(U+\frac{u}{2}\biggr)[F(U,|u|\omega)-F(U,u)], \\ F(U,u)=f(v)f(w). \end{gathered}
\end{equation*}
\notag
$$
If we set $u=r\omega_0$, $\omega_0 \in S^2$, and write the integral with respect to $du$ as
$$
\begin{equation*}
\int_{\mathbb{R}^3}\,du\,\varphi(u)=\int_0^\infty\,dr \, r^2 \int_{S^2}\,d\omega_0\,\varphi(r\omega_0),
\end{equation*}
\notag
$$
then the inner integral with respect to $d\omega_0\,d\omega$ reads
$$
\begin{equation*}
I=\int_{{S}^2\times S^2} d\omega_0\,d\omega\, \psi\biggl(U+\frac{r}{2}\omega_0\biggr)[F(U,r\omega)-F(U,r\omega_0)],
\end{equation*}
\notag
$$
where $r=|u|$. Obviously, we can interchange the variables $\omega$ and $\omega_0$ in the first term and obtain
$$
\begin{equation*}
I=\int_{{S}^2\times S^2}d\omega_0\,d\omega\,F(U,r\omega_0) \biggl[\psi\biggl(U+\frac{r}{2}\,\omega\biggr)- \psi\biggl(U+\frac{r}{2}\,\omega_0\biggr)\biggr].
\end{equation*}
\notag
$$
Returning to the original variables we get that
$$
\begin{equation}
\langle \psi,Q(f,f)\rangle=\int_{\mathbb{R}^3 \times \mathbb{R}^3} dv\,dw\,f(v)f(w)|u| \sigma_{\rm tot}[\overline{\psi(v')-\psi(v)}],
\end{equation}
\tag{2.64}
$$
where
$$
\begin{equation}
\begin{gathered} \, \sigma_{\rm tot}=\int_{S^2}\,d\omega\,\sigma_{\rm diff}=\pi d^2, \\ \overline{\psi(v')-\psi(v)}=\frac{1}{\sigma_{\rm tot}} \int_{S^2}\, d\omega\,\sigma_{\rm diff}\biggl[\psi\biggl(U+\frac{|u|}{2}\omega\biggr)- \psi\biggl(U+\frac{u}{2}\biggr)\biggr], \end{gathered}
\end{equation}
\tag{2.65}
$$
and $\sigma_{\rm diff}=d^2/4$ denotes the so-called differential scattering cross-section for hard spheres with diameter $d$. A bar in (2.64) and (2.65) means actually averaging over the random impact parameters. The physical meaning of equation (2.62) becomes clear if we write it as
$$
\begin{equation}
\partial_t\langle f,\psi\rangle+\partial_x \cdot \langle f,v\psi\rangle= \int_{\mathbb{R}^3\times \mathbb{R}^3}dv\,dw\,f(v)f(w)|u|\sigma_{\rm tot} [\overline{\psi(v')-\psi(v)}]
\end{equation}
\tag{2.66}
$$
in the notation of (2.65). Indeed, the average total number of collisions per unit time is given by the following integral:
$$
\begin{equation}
\nu_{\rm tot}(f,f)=\int_{\mathbb{R}^3\times \mathbb{R}^3}dv\,dw\,f(v)f(w) |u|\sigma_{\rm tot}.
\end{equation}
\tag{2.67}
$$
On the other hand the average change of $\psi$ in the collision of particles with velocities $v$ and $w$ is equal to the average difference $[\overline{\psi(v')-\psi(v)}]$ given in (2.65). Hence the right-hand side of (2.66) defines correctly (at the intuitive level) the rate of change of $\langle f,\psi \rangle$ due to collisions. These considerations allow us to generalize formulas (2.64)–(2.67) to the case of general (repulsive) potential $\Phi(r)$ with finite radius of action $R_{\max}=d$. In that case we have the same total cross-section $\sigma_{\rm tot}=\pi d^2$ as for hard spheres with diameter $d$. However, the differential cross-section of scattering $(v,w) \to (v',w')$, such that is given for the general potential $\Phi(r)$ by the function
$$
\begin{equation}
\sigma_{\rm diff}=\sigma\biggl(|u|,\omega \cdot \frac{u}{|u|}\biggr),\qquad \sigma_{\rm tot}=2 \pi \int_{-1}^{1}\,d\mu\,\sigma(|u|,\mu)
\end{equation}
\tag{2.68}
$$
which is discussed in detail in textbooks on classical mechanics (see [41]). The connection of $\sigma_{\rm diff}(|u|,\mu)$ with the potential $\Phi(x)$ is briefly discussed below in § 2.8. If we fix the intermolecular potential $\Phi(r)$ and the corresponding differential scattering cross-section $\sigma(|u|,\mu)$, then we obtain the same equation (2.64), where
$$
\begin{equation}
\begin{gathered} \, \overline{\psi(v')-\psi(v)}=\frac{1}{\sigma_{\rm tot}}\int_{{S}^2} d\omega \, \sigma\biggl(|u|,\omega \cdot \frac{u}{|u|}\biggr) [\psi(v')-\psi(v)], \\ v'=\frac{1}{2}(v+w+|u|\omega),\quad u=v-w. \end{gathered}
\end{equation}
\tag{2.69}
$$
Finally, we note that the total cross-section $\sigma_{\rm tot}$ disappears after the substitution of (2.65) into (2.64). In addition, the differential cross-section $\sigma(|u|,\mu)$ is always multiplied by $|u|$ in the integrand in (2.64). Therefore, it is more convenient to introduce the new function
$$
\begin{equation}
g(|u|,\mu)=|u|\sigma(|u|,\mu),\qquad \mu \in [-1,1].
\end{equation}
\tag{2.70}
$$
Then the general equation (2.64) reads
$$
\begin{equation}
\begin{gathered} \, \langle \psi,Q(f,f)\rangle= \int_{\mathbb{R}^3 \times \mathbb{R}^3\times S^2} dv\,dw\,d\omega\, f(v) f(w) g\biggl(|u|,\omega \cdot \frac{u}{|u|}\biggr)[\psi(v')-\psi(v)], \\ u=v-w,\quad v'=\frac{1}{2} [v+w+|u| \omega], \quad \omega \in S^2. \end{gathered}
\end{equation}
\tag{2.71}
$$
Note that $\psi(v)$ is an arbitrary test function. Therefore, this expression is sometimes called ‘the weak form of the Boltzmann collision integral’. The corresponding strong form of $Q(f,f)$ is
$$
\begin{equation}
\begin{gathered} \, Q(f,f)=\int_{\mathbb{R}^3 \times S^2} dw \, d \omega\, g\biggl(|u|,\omega \cdot \frac{u}{|u|}\biggr)[f(v')f(w')-f(v)f(w)], \\ \omega \in S^2,\quad u=v-w,\quad v'=\frac{1}{2} (v+w+|u| \omega), \quad w'=\frac{1}{2} (v+w-|u| \omega), \end{gathered}
\end{equation}
\tag{2.72}
$$
where $g(|u|,\mu)$ is given in (2.70). We consider below the Boltzmann equation in the notation of (2.70) and (2.72). We shall usually consider $g(|u|,\mu)$ as a given function. 2.6. Basic properties of the Boltzmann equationIn applications to rarefied gas dynamics we are mainly interested not in the distribution function $f(x,v,t)$ itself, but in the (macroscopic) characteristics of the gas averaged over the velocity space. In accordance with the physical meaning of $f(x,v,t)$, the density of the gas or, equivalently, the number of particles per unit volume is defined by the equality
$$
\begin{equation}
\rho(x,t)=\langle f,1 \rangle=\int_{\mathbb{R}^3}dv\,f(x,v,t),\qquad x \in\mathbb{R}^3,\quad t \geqslant 0.
\end{equation}
\tag{2.74}
$$
Other important macroscopic characteristics of the gas are the bulk or mean velocity $u(x,t)$ (not to be confused with the notation $u$ for the relative velocity in the collision integral (2.72)) and the absolute temperature $T(x,t)$. These functions are defined by the equalities
$$
\begin{equation}
u(x,t)=\frac{1}{\rho}\langle f,v\rangle,\quad T(x,t)=\frac{m}{3 \rho}\,\langle f,|v-u|^2\rangle,
\end{equation}
\tag{2.75}
$$
in the notation of (2.55). Here $m$ is the molecular mass, while $T$ is expressed in energy units. Usually, we assume in this paper that $m=1$ unless a mixture of different gases is considered. For given values of $\rho(x,t)$, $u(x,t)$, and $T(x,t)$ the following distribution function is called the local Maxwell distribution (or Maxwellian):
$$
\begin{equation}
f_{\rm M}(x,v,t)=\rho(2\pi T)^{-3/2}\exp\biggl(-\frac{|v-u|^2}{2 T}\biggr),
\end{equation}
\tag{2.76}
$$
where it is assumed that $m=1$ in (2.75). The same function is called the ‘absolute Maxwellian’ if the parameters $\rho$, $u$, and $T$ are independent of $x$ and $t$. Returning to the Boltzmann equation (2.73) we can easily understand the importance of the Maxwellian distribution (2.76). Indeed, it follows from (2.72) that which are the conservation laws for the momentum and energy in each pair collision. Hence any function of the form
$$
\begin{equation}
f(v)=\exp(\alpha+\beta \cdot v-\gamma|v|^2),\qquad \gamma >0,
\end{equation}
\tag{2.78}
$$
with constant parameters $(\alpha,\beta,\gamma)$ satisfies the equations
$$
\begin{equation}
f(v')f(w')=f(v)f(w),\qquad v \in \mathbb{R}^3,\quad w \in \mathbb{R}^3,\quad \omega \in S^2,
\end{equation}
\tag{2.79}
$$
in the notation of (2.72). Hence for any local Maxwellian (2.76). Another important property of the Boltzmann equation is connected with the conservation laws for the mass, momentum, and energy. We consider identity (2.70) for a fixed test function $\psi(v)$ and transform the integral by interchanging the variables $v$ and $w$. Then we easily obtain
$$
\begin{equation}
\begin{aligned} \, \nonumber \langle \psi,Q(f,f)\rangle&=\frac{1}{2}\, \int_{\mathbb{R}^3 \times \mathbb{R}^3 \times S^2}\,dv\,dw\,d\omega\, f(v) f(w)g\biggl(|u|,\omega \cdot \frac{u}{|u|}\biggr)[\psi(v') \\ &\qquad+\psi(w')-\psi(v)-\psi(w)] \end{aligned}
\end{equation}
\tag{2.81}
$$
in the notation of (2.72). Hence
$$
\begin{equation}
\langle \psi,Q(f,f)\rangle=0 \quad\text{if} \ \ \psi(v)=a+b\cdot v+ c|v|^2
\end{equation}
\tag{2.82}
$$
for arbitrary constant parameters $a$, $b$, and $c$. This identity leads to conservation laws for the mass, momentum, and energy. Indeed, consider equation (2.62) for $\psi=1$, $\psi=v$, and $\psi=|v|^2$, respectively. Then we obtain
$$
\begin{equation}
\begin{gathered} \, \partial_t \rho+\operatorname{div} \rho u =0, \\ \partial_t \rho u_\alpha+\partial_{x_\beta} \langle f,v_\alpha v_\beta \rangle =0,\qquad \alpha,\beta=1,2,3; \\ \partial_t \langle f,|v|^2 \rangle+ \operatorname{div}\langle f,|v|^2 v \rangle=0. \end{gathered}
\end{equation}
\tag{2.83}
$$
These equations are very basic for the Boltzmann equation and the transition to hydrodynamics. Finally, we prove the famous Boltzmann $H$-theorem. This theorem is based on the following inequality: To prove it we need one more identity for $\langle \psi,Q(f,f)\rangle$, namely,
$$
\begin{equation}
\begin{aligned} \, \nonumber \langle \psi, Q(f,f)\rangle&=-\frac{1}{4} \int_{\mathbb{R}^3 \times \mathbb{R}^3\times S^2} dv\,dw\,d\omega\, g\biggl(|u|,\omega \cdot \frac{u}{|u|}\biggr) \\ &\qquad\times [f(v') f(w')-f(v) f(w)][\psi(v')+\psi(w')-\psi(v)-\psi(w)]. \end{aligned}
\end{equation}
\tag{2.85}
$$
It follows from (2.81) and the other general equality
$$
\begin{equation}
\int_{\mathbb{R}^3 \times \mathbb{R}^3\times S^2} dv\,dw\,d\omega\, [\Psi(v,w;v',w')-\Psi(v',w';v,w)]=0,
\end{equation}
\tag{2.86}
$$
which is valid for any function $\Psi(v_1,v_2;v_3,v_4)$ such that the integral is convergent. For the proof it is sufficient to pass to the variables $U$ and $u$ from (2.63) in the integrand and to repeat the considerations used for the proof of identity (2.81). For brevity we omit these straightforward calculations. To complete the proof of (2.85) we consider (2.81) and set
$$
\begin{equation*}
\begin{gathered} \, \Psi(v,w;v',w')=\frac{1}{2}f(v)f(w)g\biggl(|u|,u'\cdot\frac{u}{|u|^2}\biggr) [\psi(v')+\psi(w')-\psi(v)-\psi(w)], \\ u'=v'-w'=|u|\omega, \qquad u=v-w. \end{gathered}
\end{equation*}
\notag
$$
Then we apply (2.86) and obtain
$$
\begin{equation*}
\langle \psi, Q(f,f)\rangle=\frac{1}{2}\, \int_{\mathbb{R}^3 \times \mathbb{R}^3\times S^2} dv\,dw\,d\omega\, [\Psi(v,w;v',w')-\Psi(v',w';v,w)],
\end{equation*}
\notag
$$
which is identity (2.85). Now we can prove inequality (2.84) by substituting $\psi=\log f(v)$ into (2.85). Then we obtain
$$
\begin{equation}
\begin{aligned} \, \nonumber \langle \log f, Q(f,f)\rangle&=-\frac{1}{4} \int_{\mathbb{R}^3 \times \mathbb{R}^3\times S^2} dv\,dw\,d\omega\, g\biggl(|u|,\omega \cdot \frac{u}{|u|}\biggr) \\ &\qquad\times[f(v') f(w')-f(v) f(w)] \log\frac{f(v')f(w')} {f(v)f(w)}\leqslant 0. \end{aligned}
\end{equation}
\tag{2.87}
$$
This completes the proof of inequality (2.84). The main application of (2.84) is the proof of Boltzmann’s $H$-theorem. We introduce the Boltzmann $H$-functional
$$
\begin{equation}
H(f)(x,t)=\langle f,\log f\rangle= \int_{\mathbb{R}^3}\,dv\,f(x,v,t)\log f(x,v,t),
\end{equation}
\tag{2.88}
$$
where $f(x,v,t)$ is a solution of equation (2.73). Note that
$$
\begin{equation*}
(\partial_t+v \cdot \partial_x)f\log f=(1+\log f)(f_t+v \cdot f_x)= (1+\log f)Q(f,f).
\end{equation*}
\notag
$$
Hence integrating with respect to $v$ we obtain
$$
\begin{equation}
\partial_t \langle f,\log f \rangle+\operatorname{div} \langle f,v \log f\rangle=\langle \log f,Q(f,f)\rangle \leqslant 0.
\end{equation}
\tag{2.89}
$$
This inequality is known as Boltzmann’s $H$-theorem. Its importance can easily be understood in the spatially homogeneous case, considered in the next section. 2.7. The spatially homogeneous problemThe Boltzmann equation (2.73) admits a class of spatially homogeneous solutions $f(v,t)$. Usually, we consider the initial value problem in the notation of (2.72).The conservation laws (2.83) show that
$$
\begin{equation}
\begin{gathered} \, \rho =\langle f,1\rangle=\mathrm{const},\qquad u=\frac{1}{\rho} \langle f,v\rangle=\mathrm{const}, \\ T=\frac{1}{3 \rho} \langle f,|v-u|^2\rangle= \frac{1}{3 \rho}[\langle f,|v|^2 \rangle-\rho |u|^2]=\mathrm{const}. \end{gathered}
\end{equation}
\tag{2.91}
$$
We also note that the operator $Q(f,f)$ is invariant under the shifts $v \to v+v_0$, $v_0 \in \mathbb{R}^3$, in the $v$-space. Therefore, if $f(v,t)$ is a solution of the equation in (2.90), then $f(v+v_0,t)$ is also a solution for any $v_0 \in \mathbb{R}^3$. Hence we can always reduce problem (2.90) to the case
$$
\begin{equation}
u=\frac{1}{\rho} \langle f,v\rangle=0,\qquad t \geqslant 0.
\end{equation}
\tag{2.92}
$$
Moreover, if $f(v,t)$ is a solution of the spatially homogeneous Boltzmann equation, then so is the function $\widetilde{f}(v,t)=\alpha f(v,\alpha t)$ for any $\alpha >0$. This transformation allows one to reduce the general problem (2.90) to the case The corresponding Maxwell distribution (2.76) reads
$$
\begin{equation}
f_{\rm M}(v)=(2 \pi T)^{-3/2}\exp\biggl(-\frac{|v|^2}{2T}\biggr),
\end{equation}
\tag{2.94}
$$
where The $H$-theorem (2.89) shows that the functional
$$
\begin{equation}
H(f)=\langle f,\log f \rangle=\int_{\mathbb{R}^3}dv\,f(v,t)\log f(v,t)
\end{equation}
\tag{2.96}
$$
cannot increase with time on a solution of (2.90) because
$$
\begin{equation}
\partial_t H(f)(t)=\langle \log f,Q(f,f)\rangle \leqslant 0.
\end{equation}
\tag{2.97}
$$
If we consider the explicit formula (2.87), then it becomes clear that $\langle f, Q(f,f)\rangle=0$ if and only if for almost all values $(v,w, \omega) \in \mathbb{R}^3\times \mathbb{R}^3\times S^2$, provided the kernel $g(|u|,\mu)$ in (2.72) is positive almost everywhere. This functional equation was studied in various classes of functions by many authors, beginning with Boltzmann (see [21] and the references there). They proved the uniqueness of its well-known solution (2.78). On the other hand the only function (2.78) which satisfies the conservation laws discussed above is the Maxwellian $f_{\rm M}$ (2.94). Hence we can conclude at the formal level that $H(f)$ decreases monotonically in time unless $f=f_{\rm M}$. This conclusion can be confirmed by the general inequality
$$
\begin{equation}
\langle f_{\rm M},\log f_{\rm M} \rangle \leqslant \langle f,\log f\rangle
\end{equation}
\tag{2.98}
$$
in the notation of (2.74)–(2.76). Its proof is very simple. Note that Therefore, it is sufficient to prove that $\langle f,\log f-\log f_{\rm M}\rangle \geqslant 0$. This follows from the elementary inequality
$$
\begin{equation*}
G(z,y)=z(\log z-\log y)+y-z=zG_1\biggl(\frac{z}{y}\biggr) \geqslant 0,\qquad z>0,\quad y>0,
\end{equation*}
\notag
$$
where $G_1(t)=\log t+t^{-1}-1 \geqslant 0$. We set $z=f(v)$, $y=f_{\rm M}(v)$ and integrate the inequality $G(f,f_{\rm M}) \geqslant 0$ over $v \in \mathbb{R}^3$. This completes the proof of (2.98). Inequality (2.98) shows that the Maxwellian (2.74)–(2.76) is the minimizer of the $H$- functional $H(f)=\langle f,\log f \rangle$ in the class of distribution functions with fixed lower moments $(\rho,u,T)$. Returning to the initial value problem (2.90) and assuming conditions (2.92) and (2.93), we know that there is a unique positive stationary solution $f_{\rm M}$, given in (2.94) and (2.95), such that all conservation laws are satisfied. This stationary solution $f_{\rm M}$ minimizes the $H$-functional, and therefore we expect that the solution $f(v,t)$ converges (in some precise sense) to $f_{\rm M}$ for large values of time. This is the qualitative behaviour of solutions of problem (2.90) that we expect on the basis of the above formal considerations. The corresponding physical process is called ‘relaxation to equilibrium’. The rigorous mathematical theory of problem (2.90) is not simple. First steps in its development were made by Carleman in the 1930s (see [19]) for the model of hard spheres. Then a more general and detailed theory was presented by Arkeryd [2] in the early 1970s (see also [3]). To understand these and more recent results in this area we need to introduce a sort of classification of collisional kernels $g(|u|,\mu)$ in the Boltzmann integral (2.72). This is done in the next subsection. 2.8. Collisional kernelsRecall that the kernel $g(|u|,\mu)$ is equal to $|u|$ times the differential cross-section $\sigma(|u|,\mu)$ expressed as a function of $\mu=\cos \theta$. The scattering angle $\theta \in [0,\pi]$ is given in the form (see [41])
$$
\begin{equation}
\theta(b,|u|)=\pi-2b|u|\int_{r_{\min}}^\infty\frac{r}{r^2} \biggl[|u|^2\biggl(1-\frac{b^2}{|u|^2}\biggr)- \frac{2\Phi(r)}{m}\biggr]^{-1/2},
\end{equation}
\tag{2.99}
$$
where $b$ is the impact parameter, $\Phi(r)$ is the intermolecular potential, and $m$ is the reduced mass of colliding particles ($m=1/2$ for particles with unit mass). To find $\sigma(|u|,\mu)$ we need to construct the inverse function $b=b(|u|,\theta)$. Then we express this function as $b=\widetilde{ b}(|u|,\cos \theta)$ and finally obtain (see [41])
$$
\begin{equation}
\sigma(|u|,\mu)= \biggl|\frac{1}{2}\partial_\mu \widetilde{ b}^2(|u|,\mu)\biggr|.
\end{equation}
\tag{2.100}
$$
Generally speaking, this is a rather complicated calculation. Fortunately it leads to a simple explicit formula $\sigma=d^2/4$ in the important case, when particles are hard spheres with diameter $d$. It also can be shown that for power-like potentials $\Phi(r)=\alpha r^{-n}$, $\alpha > 0$, we obtain
$$
\begin{equation*}
\sigma(|u|,\mu)= \biggl(\frac{\alpha}{m|u|^2}\biggr)^{2/n}\widetilde{A}_n(\mu),\qquad \mu=\cos \theta,\quad n \geqslant 1,
\end{equation*}
\notag
$$
where a function of $\theta$ is expressed as the function $\widetilde{A}_n$ of $\mu$. Hence in the case of power-like potentials $\Phi(r)=\alpha r^{-n}$ the collisional kernel in (2.72) reads
$$
\begin{equation}
g(|u|,\mu)= |u|^{\gamma_n} g_n(\mu),\qquad \gamma_n=1-\frac{4}{n}\,.
\end{equation}
\tag{2.101}
$$
We can use the same formula for hard spheres by assuming that $n=\infty$ and $g_\infty=d^2/4$. There is, however, an important difference between hard spheres and power-like potentials. We consider the collision integral (2.72) again and split it formally into two parts: where
$$
\begin{equation}
\begin{gathered} \, Q^{\rm loss}(f,f)=f(v)\nu(v), \\ \nu(v)=\int_{\mathbb{R}^3}\,dw\,f(w)g_{\rm tot}(|v-w|), \\ g_{\rm tot}(|u|)=|u|\sigma_{\rm tot}(|u|)= 2\pi|u|\int_{-1}^{1}\,d\mu\,\sigma(|u|,\mu). \end{gathered}
\end{equation}
\tag{2.103}
$$
It was already discussed in § 2.5 that $\sigma_{\rm tot}=\pi R^2_{\max}$, where $R_{\max}$ denotes the radius of action of the potential. In the case of hard spheres of diameter $d$ or any potential with $R_{\max}=d$ we obtain the universal formula for the collision frequency:
$$
\begin{equation}
\nu(v)=\pi d^2 \int_{\mathbb{R}^3}\,dw\,f(w)|v-w|.
\end{equation}
\tag{2.104}
$$
However, if we consider the power-like potential for any $n > 0$, then $R_{\max}=\infty$, and therefore the integral $\nu(v)$ diverges. Hence the splitting (2.102) is impossible, though the ‘full’ collision integral (2.72) can be convergent. The matter is that the kernel $g(|u|,\mu)$ has a non-integrable singularity at $\mu=1$, that is, $\theta=0$, for long-range potentials with $R_{\max}=\infty$. At the same time $v'=v$ and $w'=w$ if $\mu=1$. Therefore, the second factor in the integrand is equal to zero at that point. It is easy to see that the integral (2.72) is convergent for a large class of functions $f(v)$, provided that in the case of power-like potentials. It can be shown that this condition is satisfied for all $n > 1$. The Coulomb case $n=1$ is always considered separately.There also were many publications in the last two decades that were related to the Boltzmann equation with long-range potentials, but we do not consider related problems below (see, for example, [1] and the references there). Our main goal is to discuss a general class of Boltzmann-type equations, which are used in applications by physicists. This class of equations is introduced in the next section. 3. Boltzmann-type kinetic equations and their discrete models3.1. A generalization of the Boltzmann equationIn this section we introduce a general class of kinetic equations that includes the spatially homogeneous Boltzmann equation (2.90) as a particular case. To this end we choose a function $F(x_1,x_2;x_3,x_4)$ of four real (or complex) variables and assume that
$$
\begin{equation}
F(x_1,x_2;x_3,x_4)=F(x_2,x_1;x_3,x_4)=F(x_1,x_2;x_4,x_3)=-F(x_3,x_4;x_1,x_2).
\end{equation}
\tag{3.1}
$$
We also introduce a non-negative function $R(u_1,u_2)$ of two vectors $u_1,u_2 \in \mathbb{R}^3$. It is usually assumed that $R(u_1,u_2)$ is invariant under rotations of $\mathbb{R}^3$. Then it is known that that is, such a function can be reduced to a function of three scalar variables. This property holds for $d$-dimensional vectors for any $d \geqslant 2$. The function $R(u_1,u_2)$ plays below the role of the kernel of a certain integral operator.Then we define the general Boltzmann-type kinetic equation for a function $f(v,t)$, where $v \in \mathbb{R}^3$ and $t \in \mathbb{R}_+$, by the equality where the general kinetic operator $K$ acts on the $v$-variable only. It is defined by the formula
$$
\begin{equation}
\begin{aligned} \, \nonumber K[f](v)&=\int_{\mathbb{R}^3 \times \mathbb{R}^3 \times \mathbb{R}^3} dv_2\,dv_3\,dv_4\,\delta[v+v_2-v_3-v_4]\delta[|v|^2+|v_2|^2-|v_3|^2-|v_4|^2] \\ &\qquad\times {R}(v-v_2,v_3-v_4)F[f(v),f(v_2);f(v_3),f(v_4)]. \end{aligned}
\end{equation}
\tag{3.4}
$$
Our closest goal is to simplify this integral for arbitrary $F(x_1,x_2;x_3,x_4)$ and, in particular, to show that $K[f]$ coincides with the Boltzmann collision integral if We note that for any $\alpha \ne 0$,
$$
\begin{equation*}
\delta(\alpha x)=| \alpha |^{-1}\delta(x)\quad\text{and}\quad |v|^2+|v_2|^2-|v_3|^2-|v+v_2-v_3|^2=-2 (v_3-v_2)\cdot(v_3-v).
\end{equation*}
\notag
$$
Hence
$$
\begin{equation}
\begin{aligned} \, K[f](v)&=\frac{1}{2}\,\int_{\mathbb{R}^3 \times \mathbb{R}^3} dv_2\, dv_3\, \delta[(v_3-v_2) \cdot (v_3-v)] \nonumber \\ &\qquad\times R(v-v_2,2v_3-v-v_2)F[f(v),f(v_2);f(v_3),f(v+v_2-v_3)]. \end{aligned}
\end{equation}
\tag{3.6}
$$
Changing the variables by the formulae we obtain
$$
\begin{equation}
\begin{aligned} \, K[f](v)&=\frac{1}{8}\, \int_{\mathbb{R}^3 \times \mathbb{R}^3} dw\, dk\, \delta\biggl(k \cdot u+\frac{|k|^2}{2}\biggr)R(u,u+k) \nonumber \\ &\qquad\times F\biggl[f(v),f(w);f\biggl(v+\frac{k}{2}\biggr), f\biggl(w-\frac{k}{2}\biggr)\biggr],\qquad u=v-w. \end{aligned}
\end{equation}
\tag{3.7}
$$
Then we use formulae (2.56) and (2.57) and obtain after simple calculations
$$
\begin{equation}
K[f](v)=\frac{1}{8}\int_{\mathbb{R}^3 \times S^2} dw\, d\omega\, |u| R(u, |u| \omega)F[f(v),f(w);f(v'),f(w')],
\end{equation}
\tag{3.8}
$$
where
$$
\begin{equation}
\omega \in S^2,\quad u=v-w,\quad v'=\frac{1}{2}(v+w +u'),\quad u'= |u| \omega, \quad w'=\frac{1}{2}(v+w -u').
\end{equation}
\tag{3.9}
$$
This formula for $K[f](v)$ is valid for any kernel $R(u_1,u_2)$. In the case of an isotropic kernel (invariant under rotations) we can use (3.2) and set
$$
\begin{equation}
R(u,|u|\omega)=8 |u|^{-1}g(|u|,\widehat{u}\cdot \omega), \qquad \widehat{u}=\frac{u}{|u|}\,.
\end{equation}
\tag{3.10}
$$
Then we obtain
$$
\begin{equation}
\begin{gathered} \, K[f](v)=\int_{\mathbb{R}^3 \times S^2} dw\, d\omega\, g(|u|, \widehat{u} \cdot \omega)F[f(v),f(w);f(v'),f(w')], \end{gathered}
\end{equation}
\tag{3.11}
$$
where the notation is as in (2.72). If, in addition, we assume that $F$ is given by (3.5), then $K[f](v)$ coincides with the Boltzmann collision integral (2.72). Note that the simplest case of a constant kernel $R(u_1,u_2)=\mathrm{const}$ in (3.4) corresponds, under condition (3.5), to the case of hard spheres for the Boltzmann equation. Thus the following statement is proved. Proposition 3.1. The equation (3.3), (3.4) with isotropic kernel (3.2) can be reduced by formal transformations to the Boltzmann-type equation (3.3), (3.11). The connection between the kernels of the corresponding integral operators is described by equality (3.10). If $F(x_1,x_2;x_3,x_4)$ in the operator $K$ in (3.4) is given by formula (3.5), then equation (3.3), (3.4) coincides with the spatially homogeneous Boltzmann equation (2.90), (2.72). Some authors (see, for example, [25]) consider a $d$-dimensional version of the integral (3.4), where $\mathbb{R}^3$ is replaced by $\mathbb{R}^d$, $d \geqslant 2$. Then all above transformations can be repeated for the integral $K[f](v)$, $v \in \mathbb{R}^d$, with minimal changes [12]. The $d$-dimensional analogue of equality (3.8) reads
$$
\begin{equation}
K[f](v)=2^{-d}\int_{\mathbb{R}^3 \times S^{d-1}}dw\, d\omega\, |u|^{d-2} R(u,|u|\omega)F[f(v),f(w);f(v'),f(w')],
\end{equation}
\tag{3.12}
$$
in the notation of (3.9), where $v \in \mathbb{R}^d$, $w \in \mathbb{R}^d$, and $\omega \in S^{d-1}$. If we assume that the kernel $R(u_1,u_2)$ is invariant under rotations of $\mathbb{R}^d$, then we can use the same identity (3.2) and set
$$
\begin{equation}
R(u,|u|\omega)=2^d|u|^{2-d}g(|u|,\widehat{u} \cdot \omega),\qquad \widehat{u}=\frac{u}{|u|}\,.
\end{equation}
\tag{3.13}
$$
The substitution of formula (3.13) into (3.12) leads to the $d$-dimensional ‘collision integral’ (3.11), where the domain of integration $\mathbb{R}^3 \times S^2$ is replaced by $\mathbb{R}^d \times S^{d-1}$. Hence Proposition 3.1 can formally be generalized to an arbitrary dimension $d \geqslant 2$. For the simplicity of presentation we mainly consider below the case $d=3$. We present without derivation two other useful forms of the integral (3.4) [12]. The first reads
$$
\begin{equation}
\begin{aligned} \, K[f](v)&=\frac{1}{8}\,\int_{\mathbb{R}^3 \times \mathbb{R}^3}dw\, dk\, \delta(k \cdot u) R\biggl(u-\frac{k}{2}\,,u+\frac{k}{2}\biggr) \nonumber \\ &\qquad\times F\biggl[f(v),f\biggl(w+\frac{k}{2}\biggr); f\biggl(v+\frac{k}{2}\biggr),f(w)\biggr],\qquad u=v-w. \end{aligned}
\end{equation}
\tag{3.14}
$$
It is clear that the integral over $k \in \mathbb{R}^3$ can be reduced to an integral over the plane orthogonal to $u \in \mathbb{R}^3$. For the Boltzmann equation a transformation of this kind was originally used by Carleman [19]. The second form of $K[f](v)$ reads
$$
\begin{equation}
K[f](v)=\frac{1}{4}\,\int_{\mathbb{R}^3 \times S^2}dw\, dn\, |u \cdot n| R(u,u')F[f(v),f(w);f(v'),f(w')],
\end{equation}
\tag{3.15}
$$
where
$$
\begin{equation}
u=v-w, \quad n \in S^{2}, \quad v'=v-(u \cdot n) n,\quad w'=w+(u \cdot n) n,\quad u'=v'-w'.
\end{equation}
\tag{3.16}
$$
Note that the notation for $v'$ and $w'$ in this formula differs from similar notation in (3.9). We do not use below representation (3.15), (3.16) for $K[f](v)$, but it should be mentioned because in the case (3.5) of the Boltzmann equation it is the most conventional form of the collision integral (see, for example, the books [20], [21], and [43]). If the kernel $R(u_1,u_2)$ is isotropic, then we can use equality (3.13) and obtain
$$
\begin{equation}
K[f](v)=\int_{\mathbb{R}^3 \times S^2} dw\, dn\,{B}(|u|, \widehat{u} \cdot n) F[f(v),f(w); f(v'), f(w')],
\end{equation}
\tag{3.17}
$$
in the notation of (3.16), where
$$
\begin{equation*}
{B}(|u|,\widehat{u} \cdot n)= 2|\widehat{u} \cdot n|g[|u|,1-2(\widehat{u} \cdot n)^2],\qquad \widehat{u}=\frac{u}{|u|}\,.
\end{equation*}
\notag
$$
In the next subsection we discuss some basic properties of the kinetic equation (3.3), (3.4). 3.2. Conservation laws and the generalized $H$-theoremFor any test function $h(v)$ we denote
$$
\begin{equation}
\langle f,h \rangle=\int_{\mathbb{R}^3}dv\, f(v)h(v)
\end{equation}
\tag{3.18}
$$
assuming that the integral exists. If $f(v,t)$ is a solution of (3.3), then we formally obtain
$$
\begin{equation}
\frac{d}{dt}\langle f,h \rangle=\langle K[f],h \rangle.
\end{equation}
\tag{3.19}
$$
After straightforward transformations by using $K[f](v)$, for example, in the form (3.11), we obtain
$$
\begin{equation}
\begin{aligned} \, \langle K[f],h \rangle&=-\frac{1}{4}\, \int_{\mathbb{R}^3 \times \mathbb{R}^3 \times S^2} dv\, dw\, d\omega\, g(|u|,\widehat{u} \cdot \omega)G(v,w;v',w') \nonumber \\ &\qquad \times [h(v')+h(w')-h(v)-h(w)], \end{aligned}
\end{equation}
\tag{3.20}
$$
where for any function $F(x_1,x_2;x_3,x_4)$ satisfying conditions (3.1). By considering the functional equation in the notation of (3.9) one can easily check that the two scalar functions $h_1=1$ and $h_3= |v|^2$, and also the vector-function $h_2=v \in {\mathbb{R}^3}$ are linearly independent solutions of this equation. The uniqueness of these solutions in different classes of functions has been proved by many authors (see the discussion of equation (3.5) in [21]). Hence we have the following conservation laws for equation (3.3), (3.4):
$$
\begin{equation}
\langle f,1 \rangle=\mathrm{const},\qquad \langle f,v \rangle=\mathrm{const},\qquad \langle f,|v|^2 \rangle=\mathrm{const},
\end{equation}
\tag{3.23}
$$
provided that conditions (3.1) are fulfilled for $F$ in (3.4). The corresponding integrals in (3.23) in the case (3.5) of the Boltzmann equation have the physical meaning of the total number of particles (gas molecules), total momentum, and total kinetic energy, respectively. These properties of the Boltzmann equation were discussed in § 2.6 above. Let us assume that there exists a function $p(x)$ such that
$$
\begin{equation}
F(x_1,x_2; x_3,x_4)[ p(x_3)+p(x_4)-p(x_1)-p(x_2)] \geqslant 0
\end{equation}
\tag{3.24}
$$
for almost all $x_i \geqslant 0$, $i=1,2,3,4$. Then we can formally introduce a generalized $H$-functional (see the end of § 2.6) on the set of non-negative solutions $f(v,t)$ of equation (3.3), (3.4) by the formula
$$
\begin{equation}
\widehat{H}[f(\,\cdot\,,t)]=\int_{\mathbb{R}^3}\,dv\,I[f(v,t)],\quad I(x)=\int_{0}^{x}\,dy\,p(y),
\end{equation}
\tag{3.25}
$$
assuming the convergence of integrals. Then formal differentiation yields
$$
\begin{equation*}
\frac{d}{dt}\widehat{H}[f(\,\cdot\,,t)]=\langle f_t,p(f) \rangle= \langle K(f),p(f)\rangle.
\end{equation*}
\notag
$$
We always assume that $F$ in (3.4) satisfies conditions (3.1). Therefore, we can apply identity (3.20) and conclude from inequality (3.24) that $\widehat{H}[f(\,\cdot\,,t)]$ cannot increase with time. In the case (3.5) of the Boltzmann equation inequality (3.24) holds for $p(x)=\log x$, and we obtain
$$
\begin{equation}
\widehat{H}(f)=\langle f(v,t),\log f(v,t)-1 \rangle.
\end{equation}
\tag{3.26}
$$
Note that $\langle f,1\rangle=\mathrm{const}$ because of the conservation laws (2.1). Therefore, the functional $\widehat{H}(f)$ is basically the same as the classical Boltzmann $H$-functional $H(f)=\langle f,\log f\rangle$ considered in § 2.6. The results of this subsection can be formulated as follows. Proposition 3.2. Equation (3.3), (3.11), where $F(x_1,x_2;x_3,x_4)$ satisfies conditions (3.1), has the same conservation laws (3.23) as the Boltzmann equation. If the function $F$ also satisfies inequality (3.24) for some function $p(x)$, $x \in \mathbb{R}_+$, then (at least formally) the functional $\widehat{H}[f(\,\cdot\,,t)]$ (3.25) cannot increase with time $t \geqslant 0$ on any solution $f(v,t)$ of equation (3.3), (3.11). Of course, all our considerations in §§ 3.1 and 3.2 were made at the formal level of mathematical rigour, since we did not specify the function $F(x_1,x_2;x_3,x_4)$ in (3.4). In the next subsection we consider some specific cases, which are different from the Boltzmann case (3.5), but are also interesting for applications. 3.3. The Nordheim–Uehling–Uhlenbeck equation and wave kinetic equationIt is clear that specific operators $K$ in (3.4) can have different functions $F(x_1,x_2;x_3,x_4)$. In all interesting applications the function $F$ can be represented as a difference of two functions:
$$
\begin{equation}
F(x_1,x_2;x_3,x_4)=P(x_3,x_4;x_1,x_2)-P(x_1,x_2;x_3,x_4).
\end{equation}
\tag{3.27}
$$
There are at least three cases of kinetic equations (3.3), (3.4) of interest to physics, for which $F$ has the structure (3.4) with different functions $P$. These are the following cases: (A) the classical Boltzmann kinetic equation (B) the quantum Nordheim–Uehling–Uhlenbeck equation for bosons and fermions [49], [54]
$$
\begin{equation}
P_{\rm NUU}(x_1,x_2;x_3,x_4)=x_1 x_2(1+\theta x_3)(1+\theta x_4),
\end{equation}
\tag{3.29}
$$
where $\theta=\pm 1$; (C) the wave kinetic equation (WKE) (see [28], [25], and the references there) Using similar notation for $F$ it is easy to verify that
$$
\begin{equation*}
F_{\rm NUU}(x_1,x_2;x_3,x_4)=F_{\rm B}(x_1,x_2;x_3,x_4)+ \theta F_{\rm W}(x_1,x_2;x_3,x_4),
\end{equation*}
\notag
$$
because $\theta^2=1$ and therefore terms of the fourth order in $F_{\rm NUU}$ vanish. A review of mathematical results for the NUU-equation can be found in [4] and [5] (see also [6] and [7]). An interesting formal generalization of this equation to the case of so-called anyons (quasi-particles with any fractional spin between 0 and 1) was also considered in these papers. This model corresponds to (3.27), with (D) anyons
$$
\begin{equation}
\begin{gathered} \, P(x_1,x_2 ; x_3,x_4)=x_1 x_2 \Phi(x_3) \Phi(x_4), \\ \Phi(x)=(1-\alpha x)^\alpha[1+(1-\alpha)x]^{1-\alpha},\qquad 0 < \alpha < 1. \end{gathered}
\end{equation}
\tag{3.31}
$$
The limiting values $\alpha=0,1$ correspond to the NUU-equation for bosons ($\alpha =0$) and fermions ($\alpha=1$). The existence of global solutions in $L^1\cap L^\infty$ of the Cauchy problem for equation (3.3), (3.4), where $F$ is given in (3.27) and (3.31), was proved in [4] under some restriction on the initial conditions and the kernel $R$ of the operator (3.4). We note that in all above cases (A)–(D) it is possible to find a function $p(x)$ that satisfies inequality (3.24). Indeed the function $F$ for cases (A), (B), and (D) can be written as
$$
\begin{equation}
\begin{aligned} \, \nonumber &F(x_1,x_2;x_3,x_4)=x_3 x_4\Phi(x_1)\Phi(x_2)-x_1 x_2\Phi(x_3)\Phi(x_4) \\ &\qquad=[\Psi(x_3)\Psi(x_4)-\Psi(x_1)\Psi(x_2)]\prod_{i=1}^{4}\Phi(x_i),\qquad \Psi(x)=\frac{x}{\Phi(x)}\,, \end{aligned}
\end{equation}
\tag{3.32}
$$
where $\Phi(x)=1$ in case (A), $\Phi(x)=(1+\theta x)^{-1}$ in case (B), and $\Phi(x)$ is given in (3.31) in case (D). Then it is easy to see that the function satisfies (3.24), provided that $\Phi(x) > 0$. The positivity condition for $\Phi(x)$ is fulfilled for all $x \geqslant 0$ in cases (A) and (B) with $\theta =1$. It is also fulfilled for $0 \leqslant x < 1$, for $\theta= -1$ in case (B) and for $0 \leqslant x < 1/\alpha$ in case (D). The known results on the existence of solutions to kinetic equations (3.3) with the corresponding operators (3.4) show that it is sufficient to satisfy these restrictions only at $t= 0$ [4]. Thus the kinetic equation (3.3) has in cases (A), (B), and (D) the monotone decreasing functionals
$$
\begin{equation*}
\widehat{H}[f(\,\cdot\,,t)]=\int_{\mathbb{R}^3}\,dv\,I[f(v,t)],\quad I(x)=\int_{0}^{x}\,dy\,[\log y-\log\Phi(y)],
\end{equation*}
\notag
$$
with the corresponding functions $\Phi(x)$. Inequality (3.24) in cases (A), (B), and (D) allows one to answer an important question about stationary solutions of equation (3.3). If $K[f^{\rm st}](v)=0$, then we can integrate this equality against any ‘nice’ function $h(v)$ and obtain the identity (with respect to $h(v)$) $\langle K(f^{\rm st}),h\rangle=0$ in the notation of (3.18). Then we take $h(v)=\log \Psi[f^{\rm st}]$ and use the transformation (3.20) and inequality (3.24). Since the resulting integral of a non-negative function over the set $\mathbb{R}^3 \times \mathbb{R}^3 \times S^2$ must be equal to zero, we conclude that this function is equal to zero almost everywhere in this set. This leads to the equation
$$
\begin{equation*}
h(v')+h(w')-h(v)-h(w)=0, \qquad h(v)=\log \Psi[f^{\rm st}],
\end{equation*}
\notag
$$
in the notation of (3.9) and (3.32). This equation holds almost everywhere in $\mathbb{R}^3 \times \mathbb{R}^3 \times S^2$. Then we obtain (see comments to equation (3.22) above)
$$
\begin{equation*}
\Psi[f^{\rm st}(v)]=\frac{f^{\rm st}(v)}{\Phi[f^{\rm st}(v)]}=M(v)= \exp(\alpha+\beta \cdot v+\gamma |v|^2),
\end{equation*}
\notag
$$
where $\alpha \in \mathbb{R}$, $\gamma \in \mathbb{R}$, and $\beta \in \mathbb{R}^3$ are arbitrary constant parameters. For brevity we do not discuss these known stationary solutions; see, for example, [5] for details. We stress that the above considerations just repeat the usual arguments in the proof of the uniqueness of the Maxwellian stationary solution to the Boltzmann equation; see, for example, [20]. Because of many similarities with the Boltzmann case (A), one can expect a similar behaviour of solutions to equation (3.3) in cases (B) and (D), in particular, convergence to stationary solutions discussed above for large values of time. The situation looks more complex in case (C) of the WKE. In this case we can also satisfy inequality (3.24) by choosing $p(x)=-1/x$. Then this inequality for $F$ from (3.27) and (3.30) reads
$$
\begin{equation}
x_1 x_2 x_3 x_4(x_1^{-1}+x_2^{-1}-x_3^{-1}-x_4^{-1})^2 \geqslant 0.
\end{equation}
\tag{3.33}
$$
However, an attempt to construct the $\widehat{H}$-functional (3.25) leads to divergent integral $I(x)$. It looks reasonable to replace this integral in (3.25) by $I(x)=-\log x$; then we formally obtain
$$
\begin{equation*}
\widehat{H}[f(\,\cdot\,,t)]=-\int_{\mathbb{R}^3}dv\,\log f(v,t).
\end{equation*}
\notag
$$
This integral is divergent for large $|v|$, because we always assume that $f(v,t) \to 0$ as $|v| \to \infty$. Below we will try to clarify the situation with the WKE by using discrete kinetic models introduced in the next subsection. 3.4. Discrete kinetic modelsThe idea of using discrete velocity models for qualitative description of solutions to the Boltzmann equation seems to be very natural. Implicitly, it was already used by Boltzmann in the first publication [16] of his famous equation. We also mention the first two toy models with a few velocities introduced by Carleman [19] and Broadwell [17], respectively. An important role in the development of this idea was played by Cabannes [18] and Gatignol [31] in the 1970 –1980s. Moreover, it was proved in the 1990s (see [13], [50], and the references there) that the Boltzmann equation can be approximated by its discrete velocity models as the number of velocities tends to infinity. These results show that discrete models can not only be used for a qualitative, but also for a quantitative description of solutions to the Boltzmann equation. A similar scheme of the construction of discrete models can be applied to the general kinetic equation (3.3), (3.4). We introduce the velocity space $V \subset \mathbb{R}^d$ that contains $n \geqslant 4$ points and replace the function $f(v,t)$ by a vector $f(t) \in \mathbb{R}^n$, where
$$
\begin{equation}
V=\{v_1,\dots,v_n\}\quad\text{and}\quad f(t)=\{f_1(t),\dots,f_n(t)\}.
\end{equation}
\tag{3.34}
$$
It is implicitly assumed here that $f_i(t)$ approximates for large $n$ the function $f(v,t)$ at the point $v=v_i \in \mathbb{R}^d$, $i=1,\dots,n$. Speaking about discrete models it is convenient to use an arbitrary dimension $d \geqslant 2$, as we will see below. The simplest and most transparent case is, of course, the planar case $d=2$. The kinetic equation (3.3), (3.4) in the $d$-dimensional case (3.12) is replaced by the following set of ordinary differential equations:
$$
\begin{equation}
\frac{df_i}{dt}=\sum_{j,k,l=1}^{n}\Gamma_{ij}^{kl} F_w(f_i,f_j ; f_k,f_l),\qquad \Gamma_{ij}^{kl}=\Gamma_{ji}^{kl}=\Gamma^{ij}_{kl},\quad 1 \leqslant i \leqslant n,
\end{equation}
\tag{3.35}
$$
where the constant (for the given set $V$) parameters $\Gamma_{kl}^{ij}$ depend only on $|v_i-v_j|=|v_k-v_l|$ and $(v_i-v_j) \cdot (v_k-v_l)$ for any integer values of the indices $1 \leqslant i,j,k,l \leqslant n$. The strict inequality $\Gamma^{kl}_{ij} > 0$ is possible only for
$$
\begin{equation}
v_i+v_j=v_k+v_l,\quad |v_i|^2+|v_j|^2=|v_k|^2+|v_l|^2.
\end{equation}
\tag{3.36}
$$
Note that equations (3.35) have a universal form for any dimension $d \geqslant 2$, though the coefficients $\Gamma^{kl}_{ij}$ can depend on $d$. The equalities (3.36) have a simple geometric meaning: the points $\{v_i,v_j,v_k,v_l\}$ form a rectangle, where the two pairs of points $\{v_i,v_j\}$ and $\{v_k,v_l\}$ belong to two different diagonals. Obviously, this geometrical meaning does not depend on the dimension. The simplest non-trivial example of the set $V \subset \mathbb{R}^2$ in (3.34) has just four ‘velocities’
$$
\begin{equation}
\begin{gathered} \, v_1=(1,0),\qquad v_2=(-1,0),\qquad v_3=(0,1),\qquad v_4=(0,-1); \\ f(v_i,t)=f_i(t),\quad i=1,\dots,4;\qquad \Gamma_{12}^{34}=\Gamma_{21}^{34}=\Gamma_{34}^{12}=1, \end{gathered}
\end{equation}
\tag{3.37}
$$
like in Broadwell’s planar model [17] of the Boltzmann equation. Equations (3.35) of this model read
$$
\begin{equation}
\frac{\partial f_1}{\partial t}=\frac{\partial f_2}{\partial t}= -\frac{\partial f_3}{\partial t}=-\frac{\partial f_4}{\partial t}= \Gamma^{34}_{12} F(f_1,f_2;f_3,f_4),
\end{equation}
\tag{3.38}
$$
where $\Gamma^{34}_{12} > 0$ is a constant. In the Boltzmann case (3.25) this equation can easily be reduced to a linear equation. The case (3.27), (3.30) of the WKE is a bit more complicated; it was discussed in [12] in more detail. In order to construct a discrete model (3.34), (3.35), which has all the relevant properties of the initial kinetic equation (3.3), (3.4), we need to impose some restrictions on the set $V$ and the coefficients of equations (3.35). Definition 3.1. The model (3.34), (3.35) is called normal if it satisfies the following conditions on the set $V$: (a) all of its $n$ elements are pairwise different and do not lie in a linear subspace of dimension $d' \leqslant d-1$ or on a sphere in $\mathbb{R}^d$; (b) the set $V$ does not have isolated points, that is, for any $1 \leqslant i \leqslant n$ in (3.35) there exist $1 \leqslant j,k,l \leqslant n$ such that $\Gamma_{ij}^{kl} > 0$; (c) if the functional equation is fulfilled for all indices $(i,j;k,l)$ for which $\Gamma_{ij}^{kl} > 0$, then there exist constants $\alpha,\gamma \in \mathbb{R}$ and $\beta \in \mathbb{R}^d$ such that $h(v)=\alpha+\beta \cdot v+\gamma |v|^2$.Methods of construction of normal models are well developed (see, in particular, [11], [14], [56], and [57]). It is easy to show that $n \geqslant 6$ for them. Most models that we consider below are assumed to be normal. The condition of normality is very important for the large-time asymptotic behaviour of solutions, as we will see below. 3.5. The properties of discrete modelsThe discrete kinetic models (3.34), (3.35) of the kinetic equation (3.3), (3.4) is uniquely defined by (a) the function $F(x_1,x_2;x_3,x_4)$ satisfying conditions (3.1); (b) the phase $V=\{v_1,\dots,v_n\} \subset \mathbb{R}^d$; (c) the set of coefficients $\Gamma=\{ \Gamma_{ij}^{kl} \geqslant 0, 1 \leqslant i, j, k, l \leqslant n\}$, where $\Gamma_{ij}^{kl}$ can depend on $|v_i-v_j |=|v_k-v_l |$ and $(v_i-v_j) \cdot (v_k-v_l)$. Recall that the inequality $\Gamma_{ij}^{kl} > 0$ is possible only under conditions (3.36). Moreover, the symmetry conditions from (3.35) are fulfilled for all elements of the set $\Gamma$. By using these symmetry conditions for $\Gamma_{ij}^{kl}$ and related conditions (3.1) for $F$ it is easy to derive from (3.35) the following identity:
$$
\begin{equation}
\frac{d}{dt} \sum_{i=1}^{n} f_i(t) h_i=-\frac{1}{4}\sum_{i,j,k,l=1}^{n} \Gamma_{ij}^{kl}F(f_i,f_j;f_k,f_l)(h_k+h_l-h_i-h_j ),
\end{equation}
\tag{3.40}
$$
where $h_1,\dots,h_n$ are constant numbers. Obviously, this is a discrete analogue of identity (3.20) for the kinetic equation (3.3), (3.4). Then conditions (3.36) lead to the following conservation laws:
$$
\begin{equation}
\sum_{i=1}^{n} f_i(t)=\mathrm{const},\qquad \sum_{i=1}^{n} f_i(t)v_i=\mathrm{const},\qquad \sum_{i=1}^{n} f_i(t) |v_i|^2=\mathrm{const},
\end{equation}
\tag{3.41}
$$
which are similar to integrals (3.23) for the kinetic equation. Note that both identity (3.40) and the conservation laws (3.41) are valid for any discrete kinetic model, not only for normal models which cannot have other linear conservation laws than the ones listed in (3.41) or linear combinations of these. However, this is not true for all normal models, but at least this holds for normal models with function $F$ satisfying inequality (3.24) for some function $p(x)$, $x > 0$. We can prove the following statement. Theorem 3.1. Assume that the model (3.34), (3.35) is normal and there exists a function $p(x)$ such that inequality (3.24) is satisfied for all $x_i > 0$, $i=1,2,3,4$. Also assume that Then
Proof. Let $I(x)$ be any function such that $I'(x)=p(x)$ for all $x > 0$. Then we set If $\{f_i(t)>0, i=1,\dots,n\}$ satisfy (3.35), then
$$
\begin{equation}
\begin{aligned} \, \nonumber &\frac{d}{dt} H[f_1(t),\dots,f_n(t)]=\sum_{i=1}^{n} p(f_i) \frac{d f_i}{dt} \\ &\qquad=-\frac{1}{4}\sum_{i,j,k,l=1}^{n} \Gamma_{ij}^{kl}F(f_i,f_j;f_k,f_l) [p(f_k)+p(f_l)-p(f_i)-p(f_j)] \leqslant 0 \end{aligned}
\end{equation}
\tag{3.46}
$$
as follows from (3.40) and (3.24). Hence (a) is proved.
To prove part (3.1-b) we use (3.46) and reduce (3.42) to the identity
$$
\begin{equation}
\sum_{i,j,k,l=1}^{n}\Gamma_{ij}^{kl}F(f_i,f_j;f_k,f_l)(h_k+h_l-h_i-h_j)=0,
\end{equation}
\tag{3.47}
$$
which is supposed to be valid for any $f_i > 0$, $i=1,\dots,n$. Assume that (b) is wrong and that in this identity $h=(h_1,\dots,h_n)$ is not a linear combination of the vectors $\varphi_1,\dots,\varphi_{d+2}$, which correspond to the conservation laws (3.41). Without loss of generality we can assume that
$$
\begin{equation}
\alpha \leqslant h_i \leqslant \beta, \qquad i=1,\dots,n,
\end{equation}
\tag{3.48}
$$
where $\alpha < \beta$ is any pair of given real numbers. Indeed, if $h=(h_1,\dots,h_n)$ satisfies (3.47), then so does
$$
\begin{equation*}
\widetilde{h}=\lambda h+\mu \varphi_1,\qquad \varphi_1=(1,1,\dots,1),
\end{equation*}
\notag
$$
where $\lambda$ and $\mu$ are arbitrary real numbers. We can always choose these numbers in such a way that conditions (3.48) for $\widetilde{h}=(\widetilde{h}_1,\dots,\widetilde{h}_n)$ are fulfilled. Tildes are omitted below. The numbers $\alpha$ and $\beta$ in (3.48) are chosen in the following way. It follows from assumption (1) of the theorem that the function $p(x)$ maps the interval $[a,b]$ to another interval, say, $[\alpha,\beta]$. Moreover, there is an inverse function $x(p)$, which maps any point $p \in [\alpha,\beta]$ to $x(p) \in [a,b]$. Then we can easily construct a counterexample to our assumption by substituting $f_i=x(h_i)$, $i=1,\dots,n$, into identity (3.47). We obtain a sum of non-negative terms and conclude that each term vanishes, that is, for any $1 \leqslant i,j,k,l \leqslant n$ such that $\Gamma_{ij}^{kl}>0$. Then we use assumption (2) of the theorem and Definition 3.1. This proves (3.1-b).
In order to prove (3.1-c) we consider identity (3.47) for a stationary solution $f^{\rm st}$. Then the left-hand side of (3.47) is equal to zero for an arbitrary vector $h=(h_1,\dots,h_n)$. We substitute $f_i=f^{\rm st}_i$ and $h_i=p(f^{\rm st}_i)$, $i=1,\dots,n$, into (3.47) and obtain again a sum of non-negative terms. The same considerations as above lead to the equalities
$$
\begin{equation*}
p(f^{\rm st}_k)+p(f^{\rm st}_l)-p(f^{\rm st}_i)-p(f^{\rm st}_j)=0
\end{equation*}
\notag
$$
for all $1 \leqslant i,j,k,l \leqslant n$ such that $\Gamma_{ij}^{kl} > 0$. Then we use Definition 3.1 again and prove (3.1-c). This completes the proof. 3.6. Some transformations of equations and initial dataWe return for a moment to the kinetic equation (3.3), (3.4) and note that this equation is invariant under rotations and translations of variable $v \in \mathbb{R}^3$ (or $v \in \mathbb{R}^d$, $d \geqslant 2$, in the general case; see (3.12)). Usually, we consider initial data $f(v,0) \geqslant 0$ for (3.3) such that The invariance of equation (3.3), (3.4) under translations means that if $f(v,t)$ is a solution of this equation, then so is $f_a(v,t)=f(v+a,t)$ for any $a \in \mathbb{R}^3$. Then
$$
\begin{equation}
\langle f_a(v,0),v \rangle=\langle f(v,0),v \rangle-a\langle f(v,0),1\rangle
\end{equation}
\tag{3.50}
$$
in the notation of (3.18). It is always assumed that $\langle f(v,0),1 \rangle \ne 0$. Therefore, we can always choose $a \in \mathbb{R}^3$ (or $a \in \mathbb{R}^d$ in the general case) in such a way that $\langle f_a(v,0),v \rangle=0$. This almost trivial observation allows us to consider only initial conditions for (3.3), (3.4) such that $\langle f(v,0),v \rangle=0$. This transition from the general kinetic equation (3.3), (3.4) to its discrete model (3.34), (3.35) preserves the translational symmetry of equations (3.35). Indeed, a shift $v \to v+a$ of the $v$-variables means the following transformation of the set $V$ in (3.34) for the discrete model:
$$
\begin{equation}
V=\{v_1,\dots,v_n\} \to V=\{v_1+a,\dots,v_n+a\}, \qquad a \in \mathbb{R}^d.
\end{equation}
\tag{3.51}
$$
Equations (3.35) of the model are connected with the set $V \subset \mathbb{R}^d$ only through the coefficients $\Gamma_{ij}^{kl}$. However, these coefficients depend (for fixed indices $1\leqslant i,j,k,l \ n$) only on differences $v_i-v_j$ and $v_k-v_l$. Therefore, each coefficient $\Gamma_{ij}^{kl}$ is invariant under translations (3.51) of the whole of $V$. It is also easy to check that Definition 3.1 of a normal discrete model is invariant under translations. In other words, if the pair $(V,\Gamma)$, where $\Gamma=\{\Gamma_{ij}^{kl}, \, 1 \leqslant i,j,k,l\}$, defines a normal model, then so does any pair $(V_a,\Gamma)$ in the notation of (3.51). 4. Convergence to equilibrium for discrete models of wave kinetic equations4.1. The statement of the problem and the formulation of resultsWe consider below the discrete models (3.34), (3.35), where
$$
\begin{equation}
F(x_1,x_2;x_3,x_4)=x_3 x_4(x_1+x_2)-x_1 x_2(x_3+x_4),
\end{equation}
\tag{4.1}
$$
that is, the models of WKE (3.3), (3.4) with function $F$ defined in (3.27) and (3.30). In simplified notation the set of ODEs (3.35) reads where
$$
\begin{equation}
\begin{gathered} \, f(t)=\{f_1(t),\dots,f_n(t)\}, \qquad Q^{\pm}(f)=\{Q^{\pm}_1(f),\dots,Q^{\pm}_n(f)\}, \\ Q^{+}(f)=\sum_{i,j,k,l=1}^{n} \Gamma_{ij}^{kl}f_k f_l(f_i+f_j),\qquad Q^{-}(f)=f_i B_i(f), \\ B_i(f)=\sum_{j,k,l=1 }^{n}\Gamma_{ij}^{kl}f_j (f_k+f_l), \qquad i=1,\dots, n, \\ \Gamma_{ij}^{kl}=\Gamma_{ji}^{kl}=\Gamma_{kl}^{ij}, \qquad 1 \leqslant i,j,k,l \leqslant n. \end{gathered}
\end{equation}
\tag{4.3}
$$
The constant coefficients $\Gamma_{ij}^{kl} \geqslant 0$ depend on the phase set
$$
\begin{equation}
V=\{v_i \in \mathbb{R}^d,i=1,\dots,n\},\qquad d \geqslant 2.
\end{equation}
\tag{4.4}
$$
We recall that the strict inequality $\Gamma_{ij}^{kl} > 0$ is possible only for indices such that
$$
\begin{equation*}
v_i+v_j=v_k+v_l,\quad |v_i|^2+|v_j|^2=|v_k|^2+|v_l|^2,\qquad 1 \leqslant i,j,k,l \leqslant n.
\end{equation*}
\notag
$$
Each coefficient $\Gamma_{ij}^{kl}$ depends on $|v_i-v_j|=|v_k-v_l|$ and $(v_i-v_j) \cdot (v_k-v_l)$. We consider the Cauchy problem for equations (4.2) and the initial conditions
$$
\begin{equation}
\begin{gathered} \, f\big|_{t=0}=f^{(0)}=\{f_1^{(0)},\dots,f_n^{(0)}\}, \qquad f_i^{(0)}>0, \quad i=1,\dots,n; \\ \sum_{i=1}^{n}f_i^{(0)}v_i=0. \end{gathered}
\end{equation}
\tag{4.5}
$$
The last restriction does not lead to a loss of generality, as explained in § 3.6. A strict positivity condition is needed in order to avoid some ‘non-typical’ solutions and simplify the proofs. For example, if $\Gamma_{12}^{34} > 0$ and $f_i^{(0)}=0$ for all $i \geqslant 5$, then $f_5(t)=\dots=f_n(t)=0$ for $t > 0$. In this case we obtain from (3.35) a simple set of four equations (3.38), which should be considered separately (see, for example, [12]). Neglecting such special cases does not look very important for the general qualitative behaviour of solutions of the model. The main result of Section 4 can be formulated in the following way. Theorem 4.1. Assume that the discrete model (4.2)–(4.4) is normal, that is, the set $V$ in (4.4) and the coefficients $\Gamma_{ij}^{kl}$, $1 \leqslant i,j,k,l \leqslant n$, satisfy Definition 3.1. Also assume that (1) $v_1=0$ in (4.4); (2) if $v_i \in V$, then $(-v_i) \in V$ for all $i=1,\dots,n$. Then the Cauchy problem for equations (4.2) with the initial conditions (4.5) has a unique solution $f(t)=\{f_1(t),\dots,f_n(t)\}$ for all $t > 0$. Moreover, (a) for all $1 \leqslant i \leqslant n$,
$$
\begin{equation}
0< f^{(0)}_{i}\exp(-c \rho_{0}^{2} t) \leqslant f_i(t) < \rho_0, \qquad \rho_0=\sum_{i=1}^{n} f^{(0)}_{i},
\end{equation}
\tag{4.6}
$$
where $c > 0$ is a constant independent of $f^{(0)}$; (b) for all $1 \leqslant i \leqslant n$, where the number $b >-M^{-1}$, $M=\max\{|v_i|^2, \, 1 \leqslant i \leqslant n\}$, is the maximum real root of the equationIt is easy to see that the function $T_0(b)$ defined by equality (4.8) decreases monotonically on the interval $-M^{-1} \leqslant b < \infty$ from its maximum value $T_0(-M^{-1})=M$ to zero as $b \to \infty$. Therefore, the root $b(T_0)$ defined in part (b) of the theorem is unique. The proof of Theorem 4.1 is given in §§ 4.2–4.5. It is based on the simple estimates (4.6), conservation laws, and on the fact that equations (4.2) have a Lyapunov function which decreases monotonically on positive solutions of (4.2). 4.2. The existence and uniqueness of global non-negative solutionsWe consider equations (4.2), (4.3) and note that
$$
\begin{equation}
\sum_{i=1 }^{n} Q_i(f)=\sum_{i=1}^{n}Q_i(f)|v_i|^2 =0,\qquad \sum_{i=1 }^{n} Q_i(f) v_i =0
\end{equation}
\tag{4.9}
$$
for each $ f(t)=\{f_1(t),\dots,f_n(t)\}$. Therefore, any solution $f(t)$ of problem (4.2), (4.5) satisfies the conservation laws
$$
\begin{equation}
\rho[f(t)]=\rho[f_0]=\rho_0,\qquad E[f(t)]=E(f_0)=\rho_0 T_0,\qquad \sum_{i=1}^{n} f_i(t) v_i=0,
\end{equation}
\tag{4.10}
$$
where
$$
\begin{equation}
\rho(f)=\sum_{i=1 }^{n} f_i,\qquad E(f)=\sum_{i=1 }^{n} f_i |v_i|^2.
\end{equation}
\tag{4.11}
$$
The third equality in (4.10) follows from (4.5); it is not used in this subsection. Note that the existence and uniqueness of local in time solutions to problem (4.2), (4.5) follow from the general theory of ODEs. However, we need to construct a global solution for the positive initial data (4.5). To do this we use a simple trick, which is more or less standard for the Boltzmann equation (see, for example, [21]). Namely, we modify equation (4.2) in the following way:
$$
\begin{equation}
\begin{gathered} \, \nonumber \frac{d\varphi}{dt}+\lambda \varphi=A(\varphi), \qquad \varphi=(\varphi_1,\dots,\varphi_n), \\ A(\varphi)=[A_1(\varphi),\dots,A_n(\varphi)], \qquad A_i(\varphi)=Q_i(\varphi)+g \varphi_i \rho^2 (\varphi), \quad 1 \leqslant i \leqslant n, \end{gathered}
\end{equation}
\tag{4.12}
$$
where
$$
\begin{equation}
\varphi\big|_{t=0}=f^{(0)},\quad \lambda=g \rho^2, \quad g=2\,\max_{1 \leqslant i,j,k \leqslant {n}}\,\sum_{l=1}^{n}\Gamma_{ij}^{kl},
\end{equation}
\tag{4.13}
$$
in the notation of (4.10) and (4.11). Note that $\lambda$ and $g$ are positive constants. It is easy to see that $A_i(\varphi) \geqslant Q^+_i(\varphi) \geqslant 0$ because
$$
\begin{equation*}
g\varphi_i\rho^2(\varphi)-Q_i^-(\varphi)= \varphi_i\sum_{j,k=1}^{n}a_{jk}^{i}\varphi_j\varphi_k,
\end{equation*}
\notag
$$
where
$$
\begin{equation*}
a_{jk}^{i}=g-2\sum_{l=1}^{n}\Gamma_{ij}^{kl}\geqslant 0, \qquad 1 \leqslant i,j,k \leqslant {n}.
\end{equation*}
\notag
$$
The terms $Q^{+}_{i}(\varphi)$ are also polynomials with non-negative coefficients for each $1 \leqslant i \leqslant n$. Therefore, for any two vectors
$$
\begin{equation}
\varphi=(\varphi_1,\dots,\varphi_n),\quad \psi=(\psi_1,\dots,\psi_n)
\end{equation}
\tag{4.14}
$$
with non-negative components $\varphi_i \geqslant 0$ and $\psi_i \geqslant 0$ such that $\varphi_i \geqslant \psi_i$ for all $i=1,\dots,n$ we obtain
$$
\begin{equation}
A_i(\varphi) \geqslant A_i(\psi)\geqslant 0, \qquad i=1,\dots,n.
\end{equation}
\tag{4.15}
$$
To construct the solution of problem (4.12), (4.13) we transform (4.12) into the integral equation
$$
\begin{equation}
\varphi(t)=f^{(0)} e^{-\lambda t}+\int_{0}^{t}\,d\tau\, e^{-\lambda(t-\tau)}A[\varphi(\tau)]
\end{equation}
\tag{4.16}
$$
and try to solve this equation using the iterations
$$
\begin{equation}
\varphi^{(k+1)}(t)=f^{(0)}e^{-\lambda t}+\int_{0}^{t}d\tau\, e^{-\lambda(t-\tau)}A[\varphi^{(k)}(\tau)],\qquad k=0,1,\dots\,,
\end{equation}
\tag{4.17}
$$
where $\varphi^{(0)}(t)=0 $. It follows from inequalities (4.15) that
$$
\begin{equation}
0 < f^{(0)}_{i} e^{-\lambda t} \leqslant \varphi^{(k)}_{i}(t) \leqslant \varphi^{(k+1)}_{i} (t), \qquad 1 \leqslant i \leqslant n, \quad k=1,2,\dots\,.
\end{equation}
\tag{4.18}
$$
Hence we obtain a monotone increasing sequence of positive functions. It remains to prove that it is bounded above. To this end we consider the sequence of sums (with some abuse of notation of (4.11))
$$
\begin{equation}
\rho^{(k) }(t)=\sum_{i= 1 }^{n} \varphi^{(k)}_{i}(t),\qquad n= 0,1,\dots\,.
\end{equation}
\tag{4.19}
$$
Using the first identity in (4.9) and the definition of $A(\varphi)$ in (4.12) we obtain
$$
\begin{equation}
\rho^{(k+1)}(t)=\rho^{(0)} e^{-\lambda t}+g\int_{0}^{t}\,d\tau\, e^{-\lambda(t-\tau)}[\rho^{(k)}(\tau)]^{3}, \qquad k \geqslant 0,
\end{equation}
\tag{4.20}
$$
where $\lambda=g\rho_{0}^{2}$. We can easily prove by induction that $\rho^{(k)} \leqslant \rho_{0}$ for all $k \geqslant 0$ because
$$
\begin{equation*}
\rho_{0}e^{-\lambda t}+g\rho_{0}^{3}\int_{0}^{t}\,d\tau\, e^{-\lambda \tau}=\rho_0[e^{-\lambda t}+(1-e^{-\lambda t})]=\rho_{0}.
\end{equation*}
\notag
$$
The sequence $\{\rho^{(k)}(t),\, k=0,1,\dots\}$ is obviously monotone increasing and bounded. Taking its limit and the limit of equations (4.20) as $k \to \infty$, we can easily show that The transition to the limit under the integral sign here and below is justified by Lebesgues’s dominated convergence theorem. For brevity we ignore sets of zero measure. On the other hand it follows from (2.15), (2.16) that $0 \leqslant \varphi^{(k)}_{i} \leqslant \rho_{0}$ for all $k \geqslant 0$ and $i=1,\dots,n$. Therefore,
$$
\begin{equation*}
\varphi(t)=\Bigl\{\varphi_{i}(t)=\lim_{k \to \infty}\varphi_{i}^{(k)}(t), \, i=1,\dots,n\Bigr\}.
\end{equation*}
\notag
$$
It follows from (4.18) and (4.19) that the function $\varphi(t)$ solves equation (4.16). Note that in the notation of (4.11). Therefore, the components of $ A(\varphi )=\{ A_i(\varphi ),\, i =1,\dots,n \} $ in equation (4.16) read (see (4.12))
$$
\begin{equation*}
A_i(\varphi)=Q_i(\varphi)+g\varphi_i(t)\rho_{0}^{2},\qquad i=1,\dots,n.
\end{equation*}
\notag
$$
Hence the equation (4.16) can be written as
$$
\begin{equation*}
\varphi(t)e^{\lambda t}=f^{(0)}+\int_{0}^{t}\,d\tau\,e^{\lambda\tau} \{Q[\varphi(\tau)]+\lambda \varphi(\tau)\},
\end{equation*}
\notag
$$
where $\lambda=g \rho_{0}^{2}$. Then we can prove by differentiation that $\varphi(t)$ solves the Cauchy problem
$$
\begin{equation}
\frac{d \varphi}{dt}=Q(\varphi), \quad \varphi\big|_{t=0}=f^{(0)}.
\end{equation}
\tag{4.22}
$$
The uniqueness of its solution follows from standard theorems for autonomous ODEs with polynomial right-hand side. Note that we have not used any connection of the coefficients $\Gamma_{ij}^{kl}$ in equation (4.3) with the set in (4.4). Thus the following lemma is almost proved. Lemma 4.1. Consider equations (4.2), (4.3). Then for any non-negative data there exists a unique global in time solution $f(t)$ of equations (4.2), (4.3). The functions $f_i(t)$ satisfy the inequalities Proof. To finish the proof it is sufficient to note that the problem (4.22) for $\varphi(t)$ coincides with the problem (4.2), (4.3), (4.23) for $f(t)$. Therefore, we just need to set $f(t)=\varphi(t)$, where $\varphi(t)$ was already constructed above. Inequalities (4.24) follow from (4.17) and (4.11). This completes the proof. It is clear that Lemma 4.1 proves the first part of Theorem 4.1 (without statement (b)) under much weaker conditions, independent of the set $V$ in (4.4) and the specific properties of normal discrete model of WKE. The role of these stronger conditions will be clear in the next subsection. 4.3. The existence of a unique stationary solutionOur goal in this section is to prove that the stationary equation (4.2), that is, the equation in the notation of (4.3), has a unique solution under the assumptions of Theorem 4.1 and the additional assumption that
$$
\begin{equation}
\sum_{i=1 }^{n} f_i=\rho, \qquad \sum_{i=1 }^{n} f_i |v_i|^2=E=\rho T, \qquad \sum_{i=1 }^{n} f_i v_i= 0,
\end{equation}
\tag{4.26}
$$
in the notation of (4.4). This fact follows from Theorem 3.1 (c) in § 3.5: we just need to check that all conditions of that theorem are satisfied. We note that the model (4.2)–(4.4) is normal by the assumptions of Theorem 4.1. Moreover, the function $p(x)=-1/x$ satisfies inequality (3.24) for $F(x_1,x_2;x_3,x_4)$ given in (4.1), as it was shown in (3.33). The function $p(x)=-1/x$ also satisfies assumptions (1) and (2) of Theorem 3.1 in § 3.5. Hence we can use part (c) of that theorem and conclude that any solution of (4.25) reads
$$
\begin{equation}
f=\{f_1,\dots,f_n\},\qquad f_i=(\alpha+\beta \cdot v_i+\gamma|v_i|^2)^{-1},\quad i=1,\dots,n,
\end{equation}
\tag{4.27}
$$
where $\alpha \in \mathbb{R}$, $\gamma \in \mathbb{R}$, and $\beta \in \mathbb{R}^d$ are free parameters. We obtain these parameters from conditions (4.26) and the assumptions of Theorem 4.1. Note that $v_1=0$, $v_i \ne v_j$ if $i \ne j$, and the set $V$ in (4.4) is invariant (perhaps with changed numbering) under the transformation $v_i \to (-v_i)$, $1 \leqslant i \leqslant n$, as it follows from these assumptions. Hence $\alpha \ne 0$ in (4.27), and we can change the notation of (4.27) to
$$
\begin{equation}
f_i=a(1+b|v_i|^2+\beta \cdot v_i)^{-1}, \qquad 1 \leqslant i \leqslant n.
\end{equation}
\tag{4.28}
$$
Then we apply the third condition from (4.26) and obtain
$$
\begin{equation}
\begin{aligned} \, \nonumber 0&=a\sum_{i=1 }^{n} v_i \psi(v_i)=\frac{a}{2} \sum_{i=1}^{n} v_i[\psi(v_i)-\psi(-v_i)] \\ &=a\sum_{i=1 }^{n}v_i(\beta \cdot v_i)\bigl[(1+b |v_i|^2)^2- (\beta \cdot v_i)^2\bigr]^{-1}, \quad \psi(v)=1+b|v|^2+\beta \cdot v. \end{aligned}
\end{equation}
\tag{4.29}
$$
We are interested in bounded non-negative functions $f_i$ in (4.28). Therefore,
$$
\begin{equation*}
(1+b|v_i|^2 )^2 > (\beta \cdot v_i)^2, \qquad 1 \leqslant i \leqslant n.
\end{equation*}
\notag
$$
Then we consider the scalar product of equality (4.29) and the constant vector $\beta \in \mathbb{R}^d$ and obtain a sum of non-negative numbers on the right-hand side. Hence $\beta=0$, since the equality $\beta \cdot v_i=0$ for all $1 \leqslant i \leqslant n$ is impossible for any normal model (see Definition 3.1). Then it follows from the algebraic equations (4.26) that
$$
\begin{equation}
a=\frac{\rho}{S_0(b)}\,, \qquad S_k(b)=\sum_{i=1}^{n}|v_i|^{2k}(1+b|v_i|^2)^{-1}, \quad k=0,1,\dots,
\end{equation}
\tag{4.30}
$$
where $b$ is the maximum real root of the algebraic equation The existence of such a root follows from simple considerations. If we consider (4.31) as the definition of the function $T(b)$, then we can compute its derivative $T'(b)$:
$$
\begin{equation*}
T'(b)=-\frac{1}{2}\sum_{i,j=1}^{n} \biggl[\frac{b(|v_i|^2-|v_j|^{2})}{(1+b|v_i|^2)(1+b|v_j|^2)}\biggr]^{2}>0.
\end{equation*}
\notag
$$
It is easy to see that $T(b)$ decreases monotonically from to zero for $b \to \infty$. Hence the inverse function $b(T)$ satisfying (4.31) is uniquely defined for all $0 < T < M$. The limiting value $b(M)=-M^{-1}$ means that the solution (4.28) becomes singular. This limit is irrelevant for Theorem 4.1 because we always have $T_0 < M$ for the initial data (4.5). The result of this subsection can be formulated as follows. Lemma 4.2. Assume that conditions of Theorem 4.1 are satisfied for equations (4.2), (4.3) and the set (4.4). Then the stationary equation (4.25) has a unique non-negative solution The proof of Lemma 4.2 has already been given above. We use the notation $f^{\rm st}$ in (4.32) in order to distinguish the stationary solution from the time-dependent solution $f(t)$ of equations (4.2), (4.3). In the next subsection we study some properties of the stationary solution (4.32). 4.4. The properties of the stationary solutionIt was shown in the beginning of the previous subsection that all conditions of Theorem 3.1 are satisfied for the function $p(x)=-1/x$ under the assumptions of Theorem 4.1. In particular, inequality (3.24) for $F(x_1,x_2;x_3,x_4)$ from (4.1) and $p(x)=-1/x$ is presented in (3.33). Hence we can apply part (a) of Theorem 3.1 from § 3.5 and equations (3.45), (3.46) from its proof. Obviously, we can choose
$$
\begin{equation}
I(x)=-\log x\quad\text{and} \quad H[f(t)]=-\sum_{i=1}^{n} \log f_i(t)
\end{equation}
\tag{4.33}
$$
and conclude that on any positive solution $f(t)=\{f_{1}(t)> 0,\dots,f_n(t)> 0\}$ of equations (4.2), (4.3). One can check directly that
$$
\begin{equation}
\frac{d}{dt} H[f(t)] =-\frac{1}{4} \sum_{i,j,k,l=1 }^{n} \Gamma_{ij}^{kl} f_i f_j f_k f_l ( f_i^{-1}+f_j^{-1}-f_k^{-1}- f_l^{-1} )^2 \leqslant 0
\end{equation}
\tag{4.34}
$$
in accordance with equation (3.33). Thus, $H(f)$ is the Lyapunov function for equations (4.2), (4.3). Let us try to minimize $H(f)$ in the domain $\Omega \subset \mathbb{R}^n$ such that
$$
\begin{equation}
\Omega(\rho,T)=\{f=(f_1,\dots,f_n)\colon\, (1)\ f_i > 0,\ i=1,\dots,n;\ (2)\ f\ \text{satisfies}\ (4.26)\},
\end{equation}
\tag{4.35}
$$
where the set $V=\{v_1,\dots,v_n\}$ in (4.4) satisfies conditions of Theorem 4.1. It is straightforward to see that the standard method of Lagrange multipliers shows that the point $f^{\rm st} \in \Omega$ in Lemma 4.1 is a unique point of extremum in $\Omega$. In fact, it is a point of minimum because
$$
\begin{equation*}
\frac{\partial^2 H}{\partial f_i\,\partial f_j}= \delta_{ij }f_i^{-2} \geqslant 0, \qquad i,j=1,\dots,n.
\end{equation*}
\notag
$$
Hence we can construct the modified Lyapunov function
$$
\begin{equation}
U(f)=H(f)-H(f^{\rm st}) \geqslant 0,\quad f \in \Omega(\rho,T),
\end{equation}
\tag{4.36}
$$
in the notation of Theorem 4.1. The point $f =f^{\rm st} \in \Omega$ is the only point in $\Omega$ where $U(f)=0$. The function $U(f)$ is used in the next subsection for the proof of Theorem 4.1. 4.5. The proof of convergence to equilibriumWe consider again the Cauchy problem in the notation of (4.2)–(4.5) and assume that all conditions of Theorem 4.1 are fulfilled. The unique solution of problem (4.37) was constructed in § 4.2. This solution obviously satisfies the conservation laws
$$
\begin{equation}
\begin{gathered} \, \sum_{i=1}^{n} f_i(t)=\sum_{i=1}^{n} f^{(0)}_{i}=\rho, \qquad \sum_{i=1}^{n} f_i(t)v_i=0, \\ \sum_{i=1}^{n}f_i(t)|v_i|^{2}=\sum_{i=1}^{n}f_i^{(0)}|v_i|^{2}=\rho T, \end{gathered}
\end{equation}
\tag{4.38}
$$
in the notation of (4.4). Moreover, $f_i(t) > 0$ for all $1 \leqslant i \leqslant n$ and all $t > 0$ because of the lower estimates in (4.24). Hence for each $t \geqslant 0$, in the notation of (4.35). Note that $\Omega \subset \mathbb{R}^n$ is a bounded domain for any values of the parameters $\rho >0$ and $0 < T < \max\{|v_i|^{2},\, i=1,\dots,n\}$. In particular, $\Omega \subset B(\rho \sqrt{d}\,)$, where $B(\rho \sqrt{d}\,)$ is the ball of radius $\rho \sqrt{d}$ centred at the origin. We have constructed the Lyapunov function $U(f)$ (see (4.36)) for equation (4.37) such that $U(f) \geqslant 0$ for any $f \in \Omega$. It was also shown that there is a unique point $f^{\rm st} \in \Omega$ such that $U(f^{\rm st})=0$. Note that the derivative
$$
\begin{equation}
\frac{d}{dt}\,U[f(t)]=\frac{d}{dt}\,H[f(t)] \leqslant 0
\end{equation}
\tag{4.40}
$$
was computed in (4.34). We can rewrite (4.40) as
$$
\begin{equation*}
\frac{d}{dt}\,U[f]=\operatorname{grad}_f U[f] \cdot Q[f]=-W(f),
\end{equation*}
\notag
$$
where the dot denotes the scalar product in $\mathbb{R}^n$ and
$$
\begin{equation}
W(f)=\frac{1}{4}\,\sum_{i,j,k,l=1}^{n}\Gamma_{ij}^{kl} f_i f_j f_k f_l(f_i^{-1}+f_j^{-1}-f_k^{-1}-f_l^{-1})^2.
\end{equation}
\tag{4.41}
$$
It is clear that under conditions of Theorem 4.1 the equation $ W(f)=0$ has a unique solution $f=f^{\rm st}$ (4.32) in $\Omega$. This fact was actually used in the proof of Lemma 4.2. Hence, similarly to $U(f)$ the function $W(f)$ has the following properties:
$$
\begin{equation}
{\rm(a)}\ \ W(f) > 0\ \ \text{if}\ \ f \in \Omega,\ \ f \ne f^{\rm st};\qquad {\rm(b)}\ \ W(f^{\rm st})=0,
\end{equation}
\tag{4.42}
$$
where $f^{\rm st} \in \Omega$ is given in (4.32). We are almost prepared to prove that on the basis of well-known facts from the theory of ODEs. To this end we begin with the following lemma.Lemma 4.3. For all $t \geqslant 0$ the solution (4.39) constructed above satisfies the inequality Proof. We consider the function $H(f)$ (4.33) and note that $H[f(t)] \leqslant H[f(0)]$ because of inequality (4.34) or, equivalently, Since $f_i(t) \leqslant \rho$ for any $i=1,\dots,n$, we obtain a lower estimate for each component of $f(t)$: This completes the proof. Hence for any trajectory $f(t)$ satisfying (4.39) we can introduce a closed bounded domain $\Omega_1 \subset \Omega$ such that
$$
\begin{equation}
\Omega_1=\biggl\{f=(f_1,\dots,f_n) \in \Omega\colon f_i \geqslant \rho^{-(n-1)}\prod_{j=1}^{n} f_j(0)\biggr\}.
\end{equation}
\tag{4.44}
$$
It follows from Lemma 4.3 that $f(t) \in \Omega_1$, $t \geqslant 0$, where $f(t)$ is the solution of problem (4.37). Note that the definition of $\Omega_1$ depends on the initial data not only through the sums (4.38), but also through the product of components of $f^{(0)}$. We recall some known applications of Lyapunov functions (see, for example, the textbook on ODEs [53]). With some abuse of notation, we consider the vector ODE (like (4.2))
$$
\begin{equation}
\frac{df}{dt}=Q(f),\qquad f \in \mathbb{R}^n, \quad Q(f) \in \mathbb{R}^n.
\end{equation}
\tag{4.45}
$$
It is assumed for simplicity that the components $Q_i(f)$, $1 \leqslant i \leqslant n$, are polynomials in the components of $f=(f_1,\dots,f_n)$. Of course, these polynomials can differ from those shown in (4.3). Let $D \subset \mathbb{R}^n$ be a closed bounded domain, and assume that there exists a solution $f(t)$ of equation (4.45) such that $f(t) \subset D$ for all $t \geqslant 0$. The following theorem is a simple modification of Theorem 5.3 in [53]. In fact, the first version of this theorem was proved by Lyapunov in 1892, in his thesis [45]. Theorem 4.2. Assume that equation (4.45) has a Lyapunov function $U(f)$ such that (a) $U(f) \geqslant 0$ for all $f \in D$, where $D \subset \mathbb{R}^n$ is closed bounded domain; (b) $W(f)=-\operatorname{grad}_f U \cdot Q(f) \leqslant 0$ for all $f \in D$; (c) both the functions $U(f)$ and $W(f)$ are continuous in $D$; (d) there exists a unique vector (a point in $\mathbb{R}^n$) $f^{\rm st} \in D$ such that Then any solution $f(t)$ of (4.45) such that $f(t) \in D$ for all $t \geqslant 0$ converges to $f^{\rm st}$, that is, Proof. The proof is simply a repetition of the proof of Theorem 5.3 in [53]. Therefore, we just outline the scheme of the proof. The first step is to prove that $U[f(t)] \to 0$ as $t \to \infty$. Assuming the opposite we obtain $U[f(t)] \geqslant \alpha > 0$, $t \geqslant 0$. Then we can prove by contradiction that $|f(t)-f^{\rm st}| \geqslant \beta$ and $W[f(t)] \geqslant \gamma$ for all $t \geqslant 0$. Hence we obtain that
$$
\begin{equation*}
\frac{d U[f(t)]}{dt} \leqslant-\gamma\quad\Longrightarrow\quad U[f(t)] \leqslant U[f(0)]-\gamma t, \qquad t \geqslant 0.
\end{equation*}
\notag
$$
This obviously contradicts the assumption that $f(t) \in D$, and therefore $U[f(t)] \geqslant 0$ for all $t \geqslant 0$. Hence $U[f(t)] \to 0$ as $t \to \infty$. The second step is to prove the limiting equality (4.47). Again, we assume the opposite. Then there exist $\varepsilon > 0$ and a sequence $\{t_k,k \geqslant 1\}$ such that $t_n \to \infty$ as $n \to \infty$ but for all $k \geqslant 1$. The sequence $\{f_k=f(t_k), k=1,2,\dots\} \subset D$ is bounded and therefore contains a convergent subsequence, which converges to a point $\overline{f} \in D$ (here we need the domain $D$ to be closed). The equality $\overline{f}=f^{\rm st}$ is impossible because of inequalities (4.48). If we substitute this convergent subsequence into the continuous Lyapunov function $U(f)$, then the corresponding sequence converges to $U(\overline{f}) \ne 0$. Hence we obtain a contradiction, and this completes the proof. The end of the proof of Theorem 4.1. It remains to apply Theorem 4.2 to the case of equations (4.2), (4.3) with initial data (4.5). We consider the solution $f \in \Omega_1$ constructed in Lemma 2.1. It is clear that all conditions of Theorem 4.1 are fulfilled. Indeed, we choose the domain $D=\Omega_1$ in the notation of (4.44) and the functions $U(f)$ and $W(f)$ in the notation of (4.27) and (4.41), respectively. The stationary solution $f^{\rm st} \in \Omega_1$ was constructed in Lemma 4.2 in the explicit form (4.32) under the assumptions of Theorem 4.1. Equalities (4.46) follow from formulae (4.36) for $U(f)$ and (4.41) for $W(f)$, under the same assumption. The uniqueness in $\Omega$ (and therefore in $\Omega_1 \subset \Omega$) of the root $f=f^{\rm st}$ of the equations $U(f)=W(f)=0$ was proved above: see (4.47) and the comments after (4.36). Hence the assumptions of Theorem 4.1 allow us to use Theorem 4.2 and prove the limiting equality (4.7). This completes the proof of Theorem 4.1. In the next subsection we discuss in greater detail the connection of Boltzmann-type kinetic equations with their discrete models. 4.6. On the approximation of Boltzmann-type equations by discrete kinetic modelsWe consider below the Boltzmann-type equation (3.1) for the distribution function $f(v,t)$, $v \in \mathbb{R}^d$, where $K[f]$ is written in the form (3.12). Omitting the irrelevant constant factor in (3.12) we obtain
$$
\begin{equation}
f_t(v,t)=K[f](v)=\int_{\mathbb{R}^d \times S^{d-1}}dw\, d\omega\, |u|^{d-2}R(u,u')G(v,w;v',w'),
\end{equation}
\tag{4.49}
$$
where $R(u,u')=R(u',u)$ and
$$
\begin{equation}
F(x_1,x_2;x_3,x_4)=F(x_2,x_1;x_3,x_4)=F(x_1,x_2;x_4,x_3)= -F(x_3,x_4;x_1,x_2),
\end{equation}
\tag{4.51}
$$
$$
\begin{equation}
\omega \in S^2,\quad u=v-w,\quad v'=\frac{v+w+u'}{2}\,,\quad u'=|u| \omega, \quad w'=\frac{v+w-u'}{2}\,.
\end{equation}
\tag{4.52}
$$
Suppose that we want to approximate the integral $K[f](v)$ by a quadrature formula on some discrete lattice in the $v$-space. To this end we introduce a regular grid in $\mathbb{R}^d$: Let $\widetilde{f}_{i}(t)$ be an approximation of $f(v_i,t)$. Then a discrete version of equation (4.49) reads
$$
\begin{equation*}
\frac{d \widetilde{f}_{i}}{dt}=\widetilde{K}_i(\widetilde{f}),\qquad \widetilde{f}=(\widetilde{f}_1,\dots,\widetilde{f}_n),\quad i=1,\dots,n,
\end{equation*}
\notag
$$
where $\widetilde{K}_i(\widetilde{f})$ denotes a quadrature formula for $K[f](v_i)$. In principle, we can use any such formula, but there are obvious advantages associated with approximations of the form
$$
\begin{equation*}
K_i(f)=\sum_{j,k,l=1}^{n}\Gamma_{ij}^{kl}F(f_i,f_j;f_k,f_l).
\end{equation*}
\notag
$$
Tildes are omitted here and below. It is assumed that the constant coefficients $\Gamma_{ij}^{kl}$ are non-negative and, moreover, $\Gamma_{ij}^{kl} > 0$ only if
$$
\begin{equation*}
v_i+v_j=v_k+v_l,\quad |v_i|^2+|v_j|^2=|v_k|^2+|v_l|^2.
\end{equation*}
\notag
$$
Then we obtain the general discrete kinetic model of equation (4.49):
$$
\begin{equation}
\frac{df_i}{dt}=\sum_{j,k,l=1}^{n}\Gamma_{ij}^{kl}F(f_i,f_j;f_k,f_l),\qquad i=1,\dots,n.
\end{equation}
\tag{4.53}
$$
Such models were discussed above, in §§ 3.4–3.6 and 4.1–4.5 – under the assumption that It was shown there that solutions of system (4.53) of ODEs inherit the properties of solutions of the kinetic equation (4.49), namely, the conservation laws and $H$-theorem (provided that the $H$-theorem is valid for equation (4.49)).Our goal is to explain some methods for constructing such approximations. Our presentation can be considered as a generalization of the papers [13] and [50], devoted to similar problems for the classical Boltzmann equation and its formal approximation by discrete velocity models proposed in [32]. Let $(k_1 h,\dots,k_d h)$ be our discrete set of points in $\mathbb{R}^d$, where $h$ is any positive number (mesh step) and $(k_1,\dots,k_d)$ are integers. We identify this grid with $\mathbb{Z}^d$ and call its points ‘integer points’. We also use below the term ‘even’ points if all the $n_i$ are even integers. A natural first step is to use the simplest rectangular formula
$$
\begin{equation}
K[f](v_i) \approx (2h)^d \sum_{v_j} \int_{S^{d-1}}\,d\omega\, R(u_{ij},|u_{ij}|\omega)G(v_i,v_j;v_i',v_j')
\end{equation}
\tag{4.54}
$$
in the notation of (4.49)–(4.52). Here $v_i$ is a given point of the lattice and the sum is taken over all points $v_j$ such that $u_{ij}=v_i-v_j$ (this choice of $v_j$ is explained below). Then the next step is to choose an approximation of the integrals over the unit sphere $S^{d-1}$
$$
\begin{equation}
I(v_i,v_j)=\int_{S^{d-1}}\,d\omega\,R(u_{ij},|u_{ij}|\omega) G(v_i,v_j;v_i',v_j').
\end{equation}
\tag{4.55}
$$
If two integer vectors $v_i$ and $v_j$ are chosen in the way explained above, then $u_{ij}=v_i-v_j$ is an even vector and $U_{ij}=(v_i+v_j)/2$ is an integer vector. We introduce the shorthand notation
$$
\begin{equation}
\varphi\biggl(\frac{|u_{ij}|^2}{4}\,,\omega\biggr)= |S^{d-1}|R(u_{ij},|u_{ij}|\omega)G(v_i,v_j;v_i',v_j'),
\end{equation}
\tag{4.56}
$$
where $v_i'$ and $v_j'$ are given in (4.52) and $|S^{d-1}|$ denotes the area of the unit sphere in $\mathbb{R}^d$. Then the integral (4.55) reads
$$
\begin{equation}
I(v_i,v_j)=\frac{1}{|S^{d-1}|}\int_{S^{d-1}}\,d\omega\, \varphi\biggl(\frac{|u_{ij}|^2}{4}\,,\omega\biggr),
\end{equation}
\tag{4.57}
$$
where $v_i \in h \mathbb{Z}^d$, $v_j \in \mathbb{Z}^d$, and $u_{ij}=v_i-v_j \in 2h\mathbb{Z}^d$ by construction. It is implicitly assumed that the function $\varphi(|u_{ij}|^2/4,\omega)$ is known only at points $\omega_k$ where $v_i'$ and $v_j'$ in (4.55) are integer vectors. Let us assume that
$$
\begin{equation}
u_{ij}=2h (n_1,\dots,n_d),\qquad |u_{ij}|^2=4 h^2 m,\qquad m=\sum_{l=1}^{d} n_l^2.
\end{equation}
\tag{4.58}
$$
Then it is clear from formula (4.52) for $v'$ and $w'$ that $v_i'$ and $v_j'$ in (4.56) are integer vectors if and only if
$$
\begin{equation*}
u_{ij}=v_i'-v_j'=2h(x_1,\dots,x_d) \in 2h\mathbb{Z}^d,
\end{equation*}
\notag
$$
where in the notation of (4.58). In our case a natural idea of the approximation of the average value (4.57) of the function $\varphi(h^2 m,\omega)$ over the unit sphere $S^{d-1}$ is to replace the integral (4.57) by the mean value over the ‘integer’ points (in the above sense) on the unit sphere. In other words, we assume that for large $m$,
$$
\begin{equation}
\begin{gathered} \, \frac{1}{|S^{d-1}|}\int_{S^{d-1}}\,d\omega\,\varphi(h^2 m,\omega) \approx \frac{1}{r_d(m)}\sum_{x \in V_d(m)} \varphi\biggl(h^2 m,\frac{x}{\sqrt{m}}\biggr), \\ V_d(m)=\{x=(x_1,\dots,x_d) \in \mathbb{Z}^d\colon x^2=m\}, \end{gathered}
\end{equation}
\tag{4.60}
$$
where $r_d(m)$ is the number of integer solutions of equation (4.59), that is, the cardinality of $V_d(m)$. The first argument $h^2m$ of $\varphi(h^2 m,\omega)$ is not important. Normally, in the general approximate formula (4.54) we consider the limit as $h \to 0$, $m \to \infty$, and $h^2 m=\mathrm{const}$. Thus, the approximate formula (4.60) is closely connected with the classical number-theoretical problem of finding integer solutions of equation (4.59). There is a large literature on this problem (see, for example, the books [44] and [35]). It is intuitively clear that the approximate equality (4.60) becomes exact for a continuous function of $\omega \in S^{d-1}$ in the limit as $m \to \infty$, provided that (1) the number $r(m) $ of solutions of (4.59) tends to infinity as $m \to \infty$, and (2) these solutions tend to be equidistributed on the sphere. This is true for $d \geqslant 4$; see [44] for details. However, the situation is less simple in the more practically interesting cases of $d=2$ and $d=3$. We discuss briefly these two cases below. All details can be found in [13] and [50] for $d=3$ and in [29] for $d=2$. We begin with the case $d=3$. The first difficulty is that equation (4.59) for $d=3$ does not have integer solutions if $m=4^{a}(8k+7)$, where $a$ and $k$ are integers. Fortunately this is irrelevant for our problem, since it is known from the equality for $m$ in (4.58) that there exists at least one integer solution of (4.59). Still, there are ‘bad’ sequences of the form $m_{a}=4^{a}m_0$, where $a\to \infty$, because $r(m_{a})=r(m_0)$ for any $a=0,1,\dots$ . It follows from elementary observations that $x_1$, $x_2$, and $x_3$ in (4.59) are even if and only if $m=0\pmod4$. If $m=4^{a}m_1$, where $m_1\ne 0\pmod4$ and $m_1 \ne 7\pmod8$, then $r_3(m) \to \infty$ as $m_{1} \to \infty$. Of course, this is not enough for a rigorous justification of the approximate formula (4.60). The more difficult problem of the equidistribution of solutions of (4.59) on the sphere of radius $\sqrt{m}$ in $\mathbb{R}^3$ was solved in [13] and [50] on the basis of rather complicated number-theoretical results due to Iwaniec [38] (see also [33], [24], and [52]). Without going into details, we cite here one of the results in [50] (see Corollary 4 on p. 1873 there). Proposition 4.1. Let $\varphi(\omega)$ be a continuous function on $S^2$. Then for every $m=1,2,3,5,6\pmod8$. Then we substitute the left-hand side of the approximate formula (4.60) in the notation of (4.56) into the approximate equality (4.54), and for $d \geqslant 3$ we obtain
$$
\begin{equation}
K[f](v_i) \approx K_h[f](v_i)=(2h)^d \sum_{v_j,v_k,v_l \in h \mathbb{Z}^d} \Gamma_{ij}^{kl}F(f_i,f_j;f_k,f_l),
\end{equation}
\tag{4.62}
$$
where $f_i=f(v_i)$, $v_i \in h \mathbb{Z}$, and
$$
\begin{equation}
\Gamma_{ij}^{kl}=\widetilde{R}(v_i-v_i,v_k-v_l) \frac{|S^{d-1}|}{r_d(|v_i-v_j|^2)/(4h^2)}\delta[(v_i-v_j)^2-(v_k-v_l)^2] \delta[v_i+v_i-v_k-v_l].
\end{equation}
\tag{4.63}
$$
Here the delta functions of integer vectors and squares of such vectors denote simply the corresponding Kronecker symbols. The kernel $\widetilde{R}(u,u')$ is defined by the equality
$$
\begin{equation*}
\widetilde{R}(u,u')=\begin{cases} \widetilde{R}(u,u') & \text{if}\ u,u' \in 2h \mathbb{Z}^d, \\ 0 & \text{otherwise}. \end{cases}
\end{equation*}
\notag
$$
The convergence was proved in [50] for the classical Boltzmann equation with $d=3$ and for non-negative continuous functions such that (see, for example, [19]) This is the result of Theorem 11 in [50]. It can easily be generalized to the case $d \geqslant 3$ for $\|f\|$ from (4.65) and to the whole class of Boltzmann-type operators $K[f](v)$ for a function $F(f_i(x_1),f_j(x_2);f_k(x_3),f_l(x_4))$ satisfying the inequality
$$
\begin{equation}
|F(f_i(x_1),f_j(x_2);f_k(x_3),f_l(x_4))| \leqslant C(R)(x_1 x_2+x_3 x_4)
\end{equation}
\tag{4.66}
$$
for all $0 \leqslant x_i \leqslant R$, $1 \leqslant i \leqslant 4$, and for any $R > 0$. We do not discuss here various estimates of the rate of convergence of the quadrature formula (4.62) as $h \to 0$, obtained for the Boltzmann equation in [50]. Instead we give some comments on the difficult case $d=2$. In fact, the quadrature formula (4.62) was also proved for the Boltzmann equation under some smoothness assumptions on the functions under the integral sign in (4.49). Roughly speaking, it was proved in [29] that the lattice points on circles are equidistributed on the average in the following sense: the exponential sums where $\theta_{u'} \in [0,2\pi]$ is the angular coordinate of $u' \in \mathbb{Z}^2$, converge to zero as $m$ tends to infinity. This is a key point in the proof in [29] of the quadrature formula (4.62) in the planar case $d=2$.There is also another way of approximation of the Boltzmann-type integrals (4.59) by infinite sums. The idea is to use Carleman’s form of the integral (see formula (3.14) in Section 3.1) with inner integral over a plane. Discrete models of the Boltzmann equation based on this idea were constructed and discussed in [51]. Independently, the convergence of certain infinite sums to the integral $K[f]$ for the WKE written in a similar form was discussed in [25]. For the sake of brevity we do not discuss here the problem of the approximation of infinite sums by finite sums. This is always possible to do for functions $f$ and $F$ satisfying conditions (4.65) and (4.66), respectively. 5. ConclusionsA large class of non-linear kinetic equations of Boltzmann type was considered in §§ 2–4 from a unified point of view. This class includes, in particular, such well-known equations as (a) the classical Boltzmann equation, (b) the quantum Nordheim–Uehling–Uhlenbeck equation, (c) the wave kinetic equation used in the theory of weak turbulence. It was shown that all these equations can naturally by considered as different forms of the general Boltzmann-type equation introduced in § 3. The general properties (conservation laws and monotone functionals) of this equation were also studied there. By analogy with discrete velocity models of the Boltzmann equation the class of discrete models of the general kinetic equation was introduced and the properties of models were studied. The long-time behaviour of solutions to discrete models of the WKE was investigated in detail in § 4. First we proved the existence of a unique global in time solution of the corresponding set of ODEs for any non-negative initial conditions. The Lyapunov function was constructed for any positive solution of the model and then used for the proof of convergence to equilibrium at the end of § 4. This result was established for so-called normal models, which do not have any spurious conservation law. Perhaps, similar results can also be proved for discrete models of the NUU-equation for fermions, but the case of WKE looks more interesting for some reasons. The matter is that it is natural to expect that the time-evolution of solutions to normal discrete kinetic models imitates, to a certain extent, the behaviour of the corresponding solutions to the kinetic equations. In principle, we can approximate with any prescribed accuracy the kinetic equation by a sequence of discrete models with sufficiently large number of discrete phase points, as it was shown in § 4.6. These arguments work very well in the case of the Boltzmann equation, for which the discrete models predict that a solution $f(v,t)$, $v \in \mathbb{R}^d$, with finite moments up to the second order tends, as $t \to \infty$, to a Maxwellian distribution of the form $M(v)=a\exp(-b|v|^2)$, under some irrelevant extra conditions. The positive parameters $a$ and $b$ are determined by conservation laws. This prediction is absolutely correct for the Boltzmann equation. On the contrary, the attractor for solutions to discrete models of the WKE has the form $f^{\rm st}(v)=a(1+b|v|^2)^{-1}$. This function is obviously not integrable in $\mathbb{R}^d$, $d \geqslant 2$, for any positive $a$ and $b$. Therefore, it cannot be an attractor for integrable solutions of the WKE. Hence a straightforward prediction of similar long-time behaviour is impossible in that case. At the same time the information about the long-time behaviour of solutions to discrete models of WKE can still be useful. We hope to return to this question in subsequent publications. AcknowledgementsI thank S. B. Kuksin for important discussions. I am also grateful to I. F. Potapenko for her help in the preparation of the manuscript. 
Citation in format AMSBIB
\Bibitem{Bob24} |
|
||||||||||||||||||||||||||||||||||||||||||||||||||||
|
Contact us: |
Terms of Use
|
Registration to the website |
Logotypes |
|