Friday, July 2, 2010

Research activity. Numerical solution of kinetic equations

The numerical solution of Boltzmann-type kinetic equations represents a major computational challenge in rarefied gas dynamics and related fields. Typically this is due to the high dimensionality of the problem and to the presence of different time and/or space scales in near-continuum regimes. The necessity of fast solvers for the kinetic integral operators is then an essential part of any numerical schemes for such problems.


My research interests can be summarized in the following topics:
  1. Spectral and fast methods for kinetic integral operators
  2. Exponential schemes for stiff nonlinear kinetic equations
  3. Monte Carlo methods
  4. Multiscale hybrid methods 
Read Less...» 


1. Spectral and fast methods for kinetic integral operators
Kinetic equations play a fundamental role in many field of applied sciences and have the form

where f = f(x,v,t) is a nonnegative function and Q(f,f) is a collision integral. The most well-known example is represented by the Boltzmann equation of rarefied gas dynamics. Besides other classical examples, like the Landau equation of plasma physics, kinetic equations appears in modelling granular gases, charged particles in semiconductors, neutron transport and quantum gases. More recently applications of kinetic equations have been considered for car traffic flows, chemotactical movements, tumor immune cells competition, coagulation-fragmentation processes, population dynamics, market economies, supply chains and many others (see also the research section on kinetic and mean field modelling).
The numerical solution of kinetic equations present several difficulties for numerical methods, mainly because of the high dimensionality of the problem. Even nowadays the deterministic numerical solution of the Boltzmann equation still represents a challenge for scientific computing.

  Figure 1: Growth of the support of f(v) at each time step

An important class of numerical methods is based on the use of spectral techniques in the velocity space. The methods were first derived in [1], inspired by previous works on the use of Fourier transform techniques (see A.V.Bobylev, Dokl. Akad. Nauk SSSR, 1975 [link] or L. Desvillettes, Riv. Mat. Univ. Parma, 2003 [pdf file]). The numerical method is based on approximating the distribution function by a periodic function in velocity space, and on its representation by Fourier series. The resulting scheme can be evaluated with a quadratic computational cost by avoiding the cost of the angular integration.
The method was further developed in [2,3] where evolution equations for the Fourier modes were explicitly derived and spectral accuracy of the method was proved. Strictly speaking these methods are not conservative, since they preserve mass, whereas momentum and energy are approximated with spectral accuracy. This trade off between accuracy and conservations seems to be an unavoidable compromise in the development of numerical schemes for the Boltzmann equation.
Figure 2: Computation of the "ghost effect" by spectral methods

We recall here that the spectral method has been applied also to the Landau equation [4,5] (follow this link for some movies), where fast algorithms can be readily derived, to the case of the grazing asymptotic limit [6] and to the case of granular gases [7] (follow this link for some movies). 
Recently in [9,10], using a suitable representation of the collision operator, the computational cost of spectral methods has been reduced from O(n2) to O(n log n) without loosing the spectral accuracy thus making the methods competitive with probabilistic Monte Carlo techniques.
Fast Boltzmann solvers based on direct discretizations of the collision integral in energy space have been  developed in [8] for quantum kinetic equations. Such schemes are able to satisfy most physical requirementes and to treat the challenging case of the Bose-Einstein condensation (follow this link for some movies)

References
  1. L.Pareschi, B.Perthame, A Fourier spectral method for homogeneous Boltzmann equations, Transport Theory and Statistical Physics, Vol.25, No.3-5, (1996), pp.369-383. 
  2. L.Pareschi, G.Russo, Numerical solution of the Boltzmann equation I: Spectrally accurate approximation of the collision operator, SIAM J. Numerical Analysis, Vol. 37, No. 4, pp. 1217-1245 (2000).  
  3. L.Pareschi, G.Russo, On the stability of spectral methods for the homogeneous Boltzmann equation, Transp. Theo. Stat. Phys. 29, 3-5, pp.431-447, (2000)
  4. L.Pareschi, G.Russo, G.Toscani, Fast spectral methods for the Fokker-Planck-Landau collision operator, J. Comp. Phys. 165, pp. 216-236, (2000).
  5. F.Filbet, L.Pareschi, A numerical method for the accurate solution of the Fokker-Planck-Landau equation in the non homogeneous case, Journal of Computational Physics, 179, 1-26 (2002). 
  6. L.Pareschi, G.Toscani, C.Villani, Spectral methods for the non cut-off Boltzmann equation and numerical grazing collision limit, Numerische Mathematik, 93, pp.527-548, (2003).  
  7. G.Naldi, L.Pareschi, G.Toscani, Spectral methods for one-dimensional kinetic models of granular flows and numerical quasi elastic limit, Mathematical Models and Numerical Analysis, 37, (2003), pp. 73-90.
  8. P.Markowich, L.Pareschi, Fast, conservative and entropic numerical methods for the Bosonic Boltzmann equation, Numerische Math. 99 (2005), pp.509--532.
  9. C.Mouhot, L.Pareschi, Fast algorithms for computing the Boltzmann collision operator, Math. Comp. 75 (2006) 1833-1852. 
  10. F.Filbet, C.Mouhot, L.Pareschi, Solving the Boltzmann equation in NlogN, SIAM J. Sci. Comput. 28, 1029 (2006)  
Read Less...»

2.  Exponential schemes for stiff nonlinear kinetic equations
The numerical solution of the Boltzmann equation close to continuum flow regimes represents a major computational challenge in rarefied gas dynamics. In such regimes, in fact, the Knudsen number (proportional to the mean free path) becomes very small and standard computational approaches loose their efficiency due to the necessity of using very small time steps in deterministic schemes or, equivalently, a large number of collisions in probabilistic approaches. The most standard techniques applied are based on domain-decomposition strategies and/or model reduction asymptotic methods. In fact, the direct time discretization of the Boltzmann equation it's a difficult task in such stiff regimes due to the high dimensionality and the nonlinearity of the collision operator which makes unpractical the use of implicit solvers.

Figure 1: Knudsen number regimes

Exponential methods based on the Wild sum expansion of the Boltzmann equation combined with splitting strategies have been introduced in [1,2] represent one possible way to overcome the problem. The main advantage of the approach proposed is that it works uniformly for a wide range of relaxation times and avoids the solution of nonlinear systems of equations even in stiff regimes. These discretizations, usually referred to as Time Relaxed (TR) methods, have been applied with success both in the context of spectral methods as well as Monte Carlo methods (see next paragraph on Monte Carlo methods). Extensions to non splitting situations have been also considered [3].

Figure 2: Convergence rates of third order exponential scheme

Recently the methods have been generalized using exponential Runge-Kutta techniques [4]. The schemes are based on a decomposition of the gain term of the collision operator into an equilibrium and a non equilibrium part which corresponds to a linearization of the collision operator around equilibrium (see S.Jin, F.Filbet, J. Comp. Phys. 2010 [pdf file]). They can be applied to any nonlinear stiff systems of ODEs which can be linearized around its unique equilibrium states that characterize its large time behavior.

References
  1. E.Gabetta, L.Pareschi, G.Toscani, Wild's sums and numerical approximation of nonlinear kinetic equations, Transport Theory and Statistical Physics, Vol.25, No.3-5, (1996), pp.515-531. 
  2. E.Gabetta, L.Pareschi, G.Toscani, Relaxation schemes for nonlinear kinetic equations, SIAM J. Numerical Analysis, Vol. 34, No. 6, pp. 2168-2194, (1997).
  3. L.Pareschi, Characteristic-based numerical schemes for hyperbolic systems with nonlinear relaxation, Rendiconti Circolo Matematico di Palermo, Serie II, Suppl. 57, pp. 375-380, (1998). 
  4. G.Dimarco, L.Pareschi, Exponential methods for kinetic equations, SIAM J. Num. Anal to appear (2011). 
Read Less...» 

3. Monte Carlo methods
Since the early 1970’s the dominant method for computation for rarefied gas dynamics (RGD) has been the Direct Simulation Monte Carlo (DSMC) method pioneered by Graeme Bird (see DSMC resources from Graeme Bird), which moves particles according to their velocities and performs collisions between randomly chosen particles. This method has been tremendously successful in a wide range of applications (for example aerospace applications, materials processing and micro-electro-mechanical systems (MEMS)). However, there is an important flow regime in which the DSMC method loses its effectiveness: flow for which the Knudsen number is small enough that the collision rate is large, but not small enough that the flow is well described by fluid mechanics. In this near-continuum regime, the external length and time scales are nearly those for fluid mechanics, but the collisional length and time scales are quite small. Since accuracy of DSMC depends on resolution of the collisional length and time scales, it becomes slow and inaccurate in this regime, which is a significant limitation on its capability for many applications.

Figure 1: Sketch of the DSMC basics.

In [1-3] we introduced a new Monte Carlo method that will accelerate DSMC in the near-continuum regime, removing this limitation. The method, inspired by the Wild sum, is based on the Time Relaxed exponential schemes for the homogeneous Boltzmann equation (see previous paragraph on Exponential Schemes) and represents an improvement over existing acceleration methods for RGD in that it is a single uniform method, valid for the full range of Knudsen numbers, with the correct asymptotic behavior in the continuum and near-continuum regimes. Numerical simulations of Time Relaxed Monte Carlo (TRMC) for two-dimensional flows have been presented in [6] (follow this link for some simulations). The method has been improved to remove time discretization errors by using a suitable recursive approach combined with adaptive strategies [4,7]. Such recursive techniques are based on computing collisions in the McKean graphs originated from the Wild sum expansions.

Figure 2. McKean graphs

Using similar ideas Monte Carlo methods have been derived for conservation laws in [5] using a relaxation approximation. More recently we have considered the case of Coulomb intercations in plasma and we have extended classical DSMC and TRMC to treat this challenging case [8].

References
  1. R.E.Caflisch, L.Pareschi, An implicit Monte Carlo method for rarefied gas dynamics I: The space homogeneous case, J. Computational Physics, 154, pp. 90-116, (1999).  
  2. L.Pareschi, G.Russo, Asymptotic preserving Monte Carlo methods for the Boltzmann equation, Transp. Theo. Stat. Phys. 29, 3-5, pp.415-430, (2000). 
  3. L.Pareschi, G.Russo, Time Relaxed Monte Carlo methods for the Boltzmann equation, SIAM J. Sci. Comput. 23 (2001), no 4, 1253--1273 
  4. L.Pareschi, B.Wennberg, A recursive Monte Carlo algorithm for the Boltzmann equation in the Maxwellian case, Monte Carlo Methods and Applications, Vol. 7, no. 3-4, pp.~349-357, (2001).  
  5. L.Pareschi, M.Seaid, A New Monte Carlo Approach for Conservation Laws and Relaxation Systems, Lecture Notes in Computer Science, Vol. 3037, pp:281-288 (2004) 
  6. L.Pareschi, S.Trazzi, Numerical solution of the Boltzmann equation by Time Relaxed Monte Carlo (TRMC) methods, International Journal of Numerical Methods in Fluids, 48, (2005), pp. 947-983.  
  7. Trazzi, Stefano ; Pareschi, Lorenzo ; Wennberg, Bernt . Adaptive and recursive time relaxed Monte Carlo methods for rarefied gas dynamics. SIAM J. Sci. Comput. 31 (2008/09), no. 2, 1379--1398.
  8. Caflisch R., Dimarco G., Pareschi L., Direct simulation Monte Carlo schemes for Coulomb interactions in plasmas, Comm. App. Ind. Math., 1, (2010), pp. 72-91. 
Read Less...» 

4. Multiscale hybrid methods  
From a numerical viewpoint, the impact of multiscale analysis on scientific computing had an explosive growth in the recent years, with applications in many fields. A broad range of scientific problems involve multiple scales and multi-scale phenomena (material science, chemistry, fluid dynamics, biology...). These involve different physical laws which govern the processes at different scales. On the computational side, several important classes of numerical methods have been developed which address explicitly the multiscale nature of the solutions (wavelets, multigrid, domain decomposition, stiff solvers, adaptive mesh refinements...).
Multi-scale methods are closely related to multi-grid and multi-resolution (wavelet) methods. Still, there are some distinctive features since the latter aim at producing an accurate solution of the microscopic solution itself. In multi-scale methods, the focus is on the macro-scale, and the micro-scale is is interesting only insofar as it helps producing the right macro-scale solution.

Figure 1: The different space-time scales

For many problems, representation or solution on the fine-scale is often impossible because of the overwhelming costs. Couplings of atomistic or molecular, and more generally microscopic stochastic models, to macroscopic deterministic models based on ODEs and PDEs is highly desirable in many applications. Similar arguments apply also to numerical methods. In fact, in addition to multiscale modelling the schemes involve the development of heterogeneous numerical methods which hybridize different approaches of probabilistic and deterministic nature.
Clearly the details of the schemes are rather problem dependent. We quote the recent works by Weinan E and Bjorn Engquist for a general approach to heterogeneous multiscale methods in scientific computing (see E Weinan, Bjorn Engquist, Xiantao Li, Weiqing Ren, Eric Vanden-Eijnden, Communications in Computational Physics, 2007 [pdf-file]).
Figure 2: Hybrid representation of the solution

A classical field where mutiscale-methods play an important rule is that of kinetic equations. In such system the time scale is proportional to a relaxation time ε and a strong model (and dimension) reduction is obtained when ε → 0. The general method intoduced in [1, 2, 4, 5] is based on the introduction of an hybrid representation of the solution which evolves accordingly to different methods, a stochastic method for the microscopic part and a deterministic method for the macroscopic one. In order to use such hybrid approach, it is essential to identify a local equilibrium function, like the Maxwellian distribution for RGD, either analytically or numerically. This local equilibrium originates a model reduction from the microscale to the macroscale formulation and allows to ignore the details of microscopic interactions in terms of simplified equations (follow this link for some movies).
A related approach has been presented recently which aims at improve the accuracy on the microscopic scale by using corrective information collected from the macroscopic scale [6] (see also a description from Pierre Degond web page here).

References 
  1. R.E.Caflisch, L.Pareschi, Towards and hybrid Monte Carlo method for rarefied gasdynamics, Ben Abdallah, Naoufel (ed.) et al., Transport in transition regimes. New York, NY: Springer. IMA Vol. Math. Appl. 135, pp.57-73 (2004). 
  2. G.Dimarco, L.Pareschi, Hybrid methods for multiscale problems I: hyperbolic relaxation system, Commun. Math. Sci. 4 (2006), no. 1, 155--177.  
  3. E.Ferrari, L.Pareschi, Hybrid Monte Carlo schemes for the diffusion of impurities in a gas, Modelling and Numerics of Kinetic Dissipative Systems, Nova-Science, New York (2006). 
  4. Dimarco, Giacomo ; Pareschi, Lorenzo . Hybrid multiscale methods. II. Kinetic equations. Multiscale Model. Simul. 6 (2007/08), no. 4, 1169--1197. 
  5. Dimarco G., Pareschi L. Fluid Solver Independent Hybrid Methods For Multiscale Kinetic Equations. Siam Journal On Scientific Computing, Vol. 32; P. 603-634, (2010). 
  6. Degond P., Dimarco G., Pareschi L.. The Moment Guided Monte Carlo Method. International Journal For Numerical Methods In Fluids, (2010), DOI: 10.1002/fld.2345.  
Read Less...»
Etichette: