Nonlinear approximation theory for the homogeneous Boltzmann equation

Nonlinear approximation theory for the homogeneous Boltzmann equation

Minh-Binh Tran
Basque Center for Applied Mathematics
Mazarredo 14, 48009 Bilbao Spain

A challenging problem in solving the Boltzmann equation numerically is that the velocity space is approximated by a finite region. Therefore, most methods are based on a truncation technique and the computational cost is then very high if the velocity domain is large. Moreover, sometimes, non-physical conditions have to be imposed on the equation in order to keep the velocity domain bounded. In this paper, we introduce the first nonlinear approximation theory for the Boltzmann equation. Our nonlinear wavelet approximation is non-truncated and based on a nonlinear, adaptive spectral method associated with a new wavelet filtering technique and a new formulation of the equation. A complete and new theory to study the method is provided. The method is proved to converge and perfectly preserve most of the properties of the homogeneous Boltzmann equation. It could also be considered as a general frame work for approximating kinetic integral equations.

Keyword Boltzmann equation, wavelet, adaptive spectral method, Maxwell lower bound, propagation of polynomial moments, propagation of exponential moments, convergence to equilibrium, conservation laws, nonlinear approximation theory, numerical stability, convergence theory, wavelet filter.
MSC: 82C40, 65M70, 76P05, 41A46, 42C40.

1 Introduction

Numerical resolution methods for the Boltzmann equation plays a very important role in the practical an theoretical study of the theory of rarefied gas. The main difficulty in the approximation of the Boltzmann equation is due to the multidimensional structure of the Boltzmann collision operator. The best known numerical method for Boltzmann equation is the Direct Simulation Monte Carlo technique by Bird [4]. The method is very efficient and preserves the main physical property of the equation; however, it is quite expensive.
After the early work of Carleman ([12], [11]), Discrete Velocity Models - DVMs has been developed as one of the main classes of deterministic algorithms to resolve the Boltzmann equation numerically ([67], [8], [9], [6], [64], [7], [56], [10], [40]). They are based on a Cartesian grid in velocity and a discrete collision operator, which is a nonlinear system of conservation laws


where is the discrete collision operator. The velocity set is assumed to be a part of the regular grid

contained in a truncated set


However, in order to guarantee the convergence, the mesh size needs to be very small and the parameter needs to be large. DVMs are then very expensive, especially if we want to observe the behaviour of the solution for large velocities, which is a very important issue in the study of the Boltzmann equation. Even for small velocities, the methods are quite expensive. The models were proved to be consistent ([57], [31]), i.e. the discrete collision term could be seen as an approximation of the real collision operator. In [49], [58], [22] the approximate solutions are proved to converge weakly to the solution of the main equation when tends to and tends to infinity by DiPerna-Lions theory ([25]). However, because of this way of truncating the mesh from to , it is not easy to obtain an accuracy estimate of errors between the approximate solutions and the global solution on the entire non-truncated space.
The second deterministic approximation is the Fourier Spectral Methods - FSMs, which were first introduced in [60] inspired by spectral methods in fluid mechanics. The methods were later developed in several works, where a new way of accelerating the algorithms was also introduced ([61], [52], [63], [35], [53], [35], [52], [62], [59], [32], [36], [33], [44]). The analysis of the methods was provided in [34]. The idea of the methods is to truncate the Boltzmann equation on the velocity space and periodize the solution on this new bounded domain. To illustrate this idea, we consider the equation


where is the truncated collision operator and is periodic on . Since is periodic on , we can write an approximation of in terms of Fourier series

which leads to a system of ODEs


The approximation was also used to derive DVMs for the Boltzmann equation ([54]) through Carleman’s representation of the equation. However, it is proved ([34]) that has a constant function equilibrium state, which is totally different from a normal equilibrium state of the Boltzmann equation. Therefore, solving does not give us the real solution of the Boltzmann equation.
The major problem with deterministic methods like DVMs and FSMs that use a fixed discretization in the velocity space is that the velocity space is approximated by a finite region. Physically, the velocity space is and even if the initial condition is compactly supported, the collision operator does not preserve this property. The collision operator indeed spreads out the supports by a factor (see [66]). Therefore in order to use both DVMs and FSMs, we have to impose nonphysical conditions to keep the supports of the solutions in the velocity space uniformly compact. For DVMs, we have to remove binary collisions which spread outside the bounded velocity space. This truncation breaks down the convolution structure of the collision operators. For FSMs, the convolution structure is perfectly preserved however we need to add nonphysical binary collisions by a periodized process. In [38], [39], Gamba and Tharkabhushanam proposed another class of FSMs, called Spectral-Lagrangian Methods (SLMs), to preserve the conservation of mass, momentum and energy on the numerical schemes. However, since these are truncated methods, we need to remove the values of the initial conditions with velocities lying outside of the truncated domain. Moreover, similar as the other Fourier-based algorithms, the positivity of the solution could not be automatically preserved. Indeed, the authors proved that if the computation domain is large enough, the negative parts of the approximate solutions are very small in the energy norms. When the computation domains are large, since we need to keep the mesh sizes small, we need to put more grid points and the methods become very expensive. Indeed, for the three deterministic approximations that we mention here DVMs and both classical FSMs and SLMs, the meshes are non-adaptive and therefore they are expensive to carry out computations on large domains. Another drawback of SLMs is that according to the theory, the initial conditions need to be regular enough to guarantee the convergence.
In order to be able to construct numerical schemes, it is natural that we require the computation domain to be bounded. To solve partial differential equations on unbounded domains, there is a famous method called Absorbing Boundary Conditions (ABCs) of Engquist and Majda ([27], [28]), where the PDE we need to solve is restricted onto a bounded domain and some artificial boundary conditions are introduced in order to guarantee that solving the equation on the bounded domain with the new boundary conditions and solving that equation on the whole space would give the same result. The construction of the ABCs is based on the partial differential structure of the PDEs and could not be used for integral equations like Boltzmann. Though it does not allow us to build ABCs, the integral structure of the Boltzmann equation gives us another advantage to build an equivalent strategy as ABCs for Boltzmann equation, which we will describe in the following: Consider the following change of variables from to

Apply this change of variables to the Boltzmann equation, we get a new formulation where the equation is considered on a bounded domain. Since there is no partial differential structure in the Boltzmann equation, the prize that we need to pay after using this change of variable is just the Jacobian of it, which is . This means that we need to introduce a weight or equivalently in all of the norms we consider for the solution of the new equation. Notice that is just the momentum with order , which appears quite a lot in the theory of Boltzmann equation. This new formulation of the Boltzmann equation is discussed in details in Section 2.
After having an equation on a bounded domain through a change of variables technique, we can construct a spectral algorithm similar as in [60]. However, different from [60], we do not use Fourier basis. Let us explain why. We recall some quantitative properties of the Boltzmann equation that we want to preserve on the numerical schemes. Notice that these properties could not be preserved with previous strategies.

  • Maxwellian lower bounds (Carleman [12], Pulvirenti and Wennberg[66]): if the initial condition satisfies



  • Production of polynomial moments (Povzner [65], Desvillettes [21], Wennberg [69], Mischler and Wennberg [50]): if the initial condition satisfies



  • Propagation of exponential moments (Bobylev, Gamba and Panferov [5], Gamba, Panferov and Villani [37], Alonso, Cañizo, Gamba and Mouhot[1]): Assume that the initial data satisfies for some

    then there are some constants such that



Suppose that we approximate by its truncated Fourier series


We can see that the approximate solution will never satisfy the properties that we mention above no matter how good is. The reason is that all components of the Fourier basis, i.e. the and functions are globally and smoothly defined on the whole interval and they encounter singular problems at the extremes and . This raises the need for a compactly supported wavelet basis and a new filtering technique. The idea of the technique is simple: we remove compactly supported wavelets which contain the singular points and . After having a good orthogonal basis based on this filtering technique, we can apply the normal spectral method to solve the equation. This filtering technique looks like a truncation technique, however it is more natural since we only need to remove some spectral components and different from classical approximations, the support of our approximate solutions spread to the whole space gradually after each approximate level . Moreover, it is designed to preserve properties , and , which are crucial in resolving the Boltzmann equation numerically. We preserve the good properties of both DVMs and FSMs: we are able to keep the convolution structure of the collision operators and do not have to impose a periodic boundary condition on the equation. In addition, our algorithm solves the entire, non-trucated problem with the complexity . The wavelet basis, the filtering technique and the spectral method will be presented in section 3. More precisely, our spectral equation is defined in . A comparison between Zuazua’s Fourier filtering technique ([70] and [71]) used to preserve the propagation, observation and control of waves and our wavelet filtering technique used to preserve the properties of propagation of polynomial and exponential moments will also be given in subsection 6.1.
We have provided our first point of view based on Absorbing Boundary Contions. In order to understand better the mechanism of our nonlinear, adaptive spectral method, we now provide a different point of view based on Nonlinear Approximation Theory ([23], [20], [24]). The fundamental problem of approximation theory is to resolve a complicated function, by simpler, easier to compute functions called ”the approximants”. The main idea of nonlinear approximation is that the approximants do not come from linear spaces but rather from nonlinear manifolds. An important application of nonlinear approximation is the adaptive finite element methods for elliptic equations originated in [3] and developed in [15], [14], [17]. These methods are based on the idea that fine meshes are put where the solutions are bad and coarse meshes are set where the solutions are good. Coming back to the Boltzmann equation, suppose that we use the Haar wavelet to solve the Botlzmann equation with the new variable on . As we see later from , and , solving the Boltzmann equation with on means that we need construct a mesh by to dividing into small cubes. To explain clearer our idea, suppose that we are in one dimension and we need to approximate the solution in a space spanned by the following orthogonal basis

Let us make the change of variable .

We can see that solving the Boltzmann equation in on a uniform mesh in is equivalent with solving the Boltzmann equation in on a non-uniform mesh in . In other words, the role of the change of variables is to construct a new non-uniform mesh to approximate the Boltzmann equation. The non-uniform mesh has the following interesting property: the larger is the coarser the mesh is, and the smaller is the finer the mesh is. This is crucial, since properties , and play the role of a preconditioning analysis in our nonlinear approximation theory: the solution of the Boltzmann equation behaves like a Maxwellian as large, which means that if is large, we only need a coarse mesh to represent the value of . This is also the main difference between our approximation and classical ones. We can see from the spectral equations and that the mapping has a ”support-stretching” effect: it maps the wavelet basis supported in to a new ”nonlinear basis” supported in the whole space, which are ”the approximants” of our nonlinear approximation. Our method therefore gives a general frame work for solving kinetic integral equations (for example, the coagulation models [29], the quantum Boltzmann equations [30],…) numerically: Suppose that we need to solve the following problem

where is some bilinear form. We approximate as

and get the approximate equation on the unknown

Moreover, our approximation also provides a general view point for both DVMs and FSMs: FSMs and DVMs are special cases of our approximation using Fourier and Haar wavelet basis. If we take Haar wavelet basis as the spectral basis, our algorithm in this special case then gives an nonlinear, adaptive DVMs for Boltzmann equation, where no direct truncation as is imposed and the convolution structure of the collision operator is perfectly preserved. Our new adaptive DVMs is then cheap and it has a spectral accuracy. Therefore, both classical DVMs and FSMs could be seen as special linear and non-adaptive approximations in our theory. We will come back to this discussion at the end of subsection 3.2.
We also introduce a full new analysis to study theoretically our algorithm. Different with the periodized case ([34]) where the truncated Boltzmann collision operator is a bounded bilinear form and the projection of the collision operator onto the subspaces could be considered as a perturbation of with a small term , in our case, the analysis is much harder since the collision operator is unbounded. Since does not preserve the symmetry of , the first problem is how we could preserve the conservation laws with this approximation

Another problem is the preservation the ”coercivity” property of the gain part of the collision operator

Notice that this is one of the main advantages of our approximation: perfectly preserve the coercivity property of the gain part of the collision operator. Approximation strategies using Fourier basis could not preserve this coercivity structure because of the effect of the Gibbs phenomenon and the non-positivity of the projection using Fourier basis. We then introduce some new methods to resolve these problems. Since the methods are quite technical, we will leave this discussion for subsection 3.3. Based on these new techniques, we can construct a new method for studying our algorithm theoretically.

  • We approximate the projected operator by bounded operators and prove that the solutions produced by the bounded operators are uniformly bounded in and , moreover they are bounded from below by a Maxwellian. Notice that different from Fourier-based spectral algorithms, our approximate solutions are automatically positive because of the positivity of the wavelet projection.

  • We prove that converges to as tends to infinity. Moreover are uniformly bounded in , and they are bounded from below by a Maxwellian.

  • We preform a detailed analysis to prove that converges to which guarantees the convergence of the algorithm. Notice that different from [49], [58] our convergence proof is not based on averaging lemma techniques and gives a strong convergence result.

In order to preserve properties and , we need to overcome further difficulties. One difficulty is to obtain an estimate of the collision operator: how we could perform a Povzner’s inequality argument for the projected collision operator since the structure of the operator is totally different. Another difficulty should be: if we have

how could we have a uniform bound with respect to for

Due to the Gibbs phenomenon and the non-positivity of the projection , Fourier basis is not a good choice to preserve these properties. This is then another advantage of spectral approximations using wavelet basis. We introduce some new methods to overcome these difficulties and we will leave these discussions for sections 6 and 7. Our main results are:

  • We prove that the algorithm converges; the energy, mass and momentum of the approximate solution converge to that of the original equation; moreover the approximate solution is bounded from below by a Maxwellian. These are the results of theorems 4.1 and 4.2.

  • We prove that the polynomial moments of arbitrary orders of the approximate solutions are uniformly bounded. This is the result of theorem 6.1.

  • We also prove that the exponential moments of the approximate solutions are uniformly bounded. This is the result of theorem 7.1.

In other words, the algorithm is proved to converge and preserve all of the properties , and as well as the conservation laws. Moreover, we could prove that the approximate solutions belong to and converges to the main solution in for all (see remark 6.2). Since our nonlinear approximation preserves well the structure of the collision operator, we could expect that other properties of the solution could be reflected on the numerical scheme as well.
We also want to mention another important property of the Boltzmann equation: In the paper [51], Mouhot proved that the solution of the Boltzmann equation converges to its equilibrium with the rate . In theorem 5.1, we prove that this property is preserved by our approximation as well.

2 A new formulation of the Boltzmann equation

2.1 The Boltzmann equation

The Boltzmann equation describes the behaviour of a dilute gas of particles when the binary elastic collisions are the only interactions taken into account. It reads


where is the spacial domain and is the time-dependent particle distribution function for the phase space. The Boltzmann collision operator is a quadratic operator defined as


where , , , and

with .
In this work, we only assume that is locally integrable and


where and is a smooth function satisfying


and assumptions - in [55]


Under these assumptions, the collision operator could be split as



Formally, Boltzmann collision operator has the properties of conserving mass, momentum and energy

and it satisfies the Boltzmann’s H-theorem

in which is defined as the entropy of the solution. A consequence of the Boltzmann’s H-theorem is that any equilibrium distribution function has the form of a locally Maxwellian distribution

where , , are the density, macroscopic velocity and temperature of the gas

We suppose that the initial datum satisfies on and

We refer to [13] and [68] for further details and discussions on the Boltzmann equation. In this work, we only consider the equation in but the methodology would be exactly the same for other dimensions.

2.2 The new formulation

Different from [60], where a truncation technique is introduced in order to reduce the Boltzmann equation defined on the whole domain into an equation on a bounded domain, we introduce in this section a new formulation of the Boltzmann equation defined on based on a change of variables technique. Let us define the following change of variables mapping


where we restrict our attention to the norm with . The inverse mapping of reads

The idea of our technique is to replace the variable in by a new variable in through the mapping . Based on this idea, we define the new density function

where is the new variable in .
With the notice that the Jacobian of the change of variable is , we have

Therefore if belongs to with the weight , then belongs to with the weight . Notice that there are several one-to-one mappings that map to however the above property makes us choose to work on .
We now define


where , are real numbers. For further use, we also need

where , are some positive weights.
Moreover, we also need the notation

The Boltzmann equation for is now


which is our first new formulation of the Boltzmann equation.
Now define

which implies

This means if belongs to then belongs to . Notice that we define to make our proof simpler, however the theoretical results remain the same if with being any constant in , could be .
The Boltzmann equation for then reads






we get our second new formulation of the Boltzmann equation


The initial datum is now defined


Let us mention that though the two new formulations seem to be complicated, we only use them for theoretical purposes. Our spectral equation is based on the former formulation of the equation.

3 Approximating the homogeneous Boltzmann equation: an adaptive spectral method

In order to preserve the quantitative properties of the Boltzmann equation , and , as it is pointed out in the introduction, we cannot use the Fourier basis. We will construct a wavelet basis for in subsection 3.1. Our new spectral algorithm is defined in equation of subsection 3.2. In subsection 3.3 we discuss about the assumption that we need for the multiresolution analysis and the wavelet filtering technique.

3.1 Wavelets for

We first construct a wavelet multiresolution analysis for . Let be a positive scaling function which defines a multiresolution analysis, i.e., a ladder of embedded approximation subspaces of

such that constitutes an orthonormal basis for . The wavelet is built to characterize the missing details between two adjacent levels of approximation. More concretely, is an orthonormal basis of where

Multiresolution analysis is a frame work developed by Mallat [41] and Meyer [46], we refer to these two pioneering works or the books [19], [48] for more details, examples and proofs.
We now follow exactly the construction in [19, Section 9.3] to build the same ”periodized wavelets” for . Notice that there are other ways besides this way (see [47], [16]). Suppose that the scaling function and the wavelet have reasonable decays, for example . Define


Similar as [19, Note 6, Chapter 9] we have

which implies

These facts mean for are one dimensional spaces of constant functions. Moreover, similar as [19, Note 7, Chapter 9] we have

and for . As a consequence, we only need to consider the spaces and with . According to the property of the multiresolution analysis , , then , . We also have that and are orthogonal