Team Calvi

Overall Objectives
Scientific Foundations
Application Domains
New Results
Contracts and Grants with Industry

Section: New Results

Keywords : Vlasov, semi-Lagrangian, cubic splines, PIC, Maxwell, MPI.

Development of Vlasov solvers

Participants : Nicolas Besse, Jean-Philippe Braeunig, Nicolas Crouseilles, Alain Ghizzo, Michael Gutnic, Matthieu Haefele, Olivier Hoenen, Guillaume Latu, Michel Mehrenberger, Thomas Respaud, Stéphanie Salmon, Eric Sonnendrücker, Eric Violard.

Two-dimensional solvers

In [45] a new Forward semi-Lagrangian method has been developed and validated. The main differences with the classical Backward semi-Lagrangian method are twofold. First, since the characteristics curves along with the unknown is constant are followed forward in time, this method enables the use of classical high order time discretization (such as Runge-Kutta 4 algorithms). Second, the remapping (or the deposition) step which is based on cubic spline polynomials, enables to reconstruct the distribution function on a uniform mesh using the particles which have moved during one time step. The analogous step for the backward method involves an interpolation to update the unknown. Several steps of the algorithm are identical or very close to those using in PIC methods. Hence, some strategies used in PIC methods for example in order to enforce the exact conservation of charge at the discrete level should be straightforward to adapt.

Another strategy based on conservative semi-Lagrangian methods has been implemented and tested. Whereas semi-Lagrangian methods compute the distribution function on the grid points, conservative methods considers an average of the unknown on each cell. This approximation enables to solve multi-dimensional problems by a successive solution of one-dimensional equations with a time splitting method. This is not the case for non-conservative methods when the advection field is not constant. In such a case, the conservative methods exhibit a very good behavior on long time simulations compared to traditional semi-Lagrangian methods. The possibility to deal with one-dimensional problems has a huge advantage from a memory point of view or for parallel algorithms. Finally, we have proved that in the particular case of a constant advective field, the two methods are algebraically equivalent.

CALVI platform

The aim of the platform is to change the way numerical methods are implemented and tested. It has been initiated because most of the researchers of the CALVI project develop new numerical methods for almost the same equations. Until now, every researchers implemented their methods as stand-alone C or Fortran applications. So, each researcher, for each code, has to implement the validation process by himself without using previous implementation done by himself or another member of the project. The platform moves the implementation from stand-alone application to a module oriented one. Thanks to standardized application programing interfaces (API), the different numerical methods can be swapped between them and can be validated within a common skeleton. This common skeleton plus the standard API is actually the platform. A better reuse of existing modules is expected as well as an increased efficiency in numerical methods implementation. The basis of the CALVI platform has been built this year and the main work consisted in architecture and API design. A first release has been provided. It contains skeleton modules: initial functions (8 test-cases with 3 validating ones), silo exports (Python bindings to the silo library), time diagnostics, and 5 drivers for validation purpose. It contains also solvers which are codes that have been modified in order to respect the API: 3 Vlasov solvers, 2 Rho solvers, and 4 Poisson solvers.

Four-dimensional solvers

A four-dimensional cubic splines interpolation has been validated in the framework of the backward semi-Lagrangian method on the Vlasov-Poisson equations. This study was motivated to study the validity of the time splitting of the method when non-conservative advection terms are involved. The method is currently tested on more complex problems and appears to be competitive from a CPU time point of view. This approach benefits from the Local Splines strategy which enables a decomposition domain well suited for parallel implementation.

Otherwise, a basic uniform semi-Lagrangian Vlasov-Poisson solver in 2 D×2 D , named split4D which takes the benefit of total directional splitting has been written in C. In order to keep data locality, we have considered some local splines and quasi-interpolants. The code has been tested for the Landau damping. The benefit of using quasi-interpolants should be further developed. A straightforward parallelization of this uniform code has been implemented with calls to MPI and OpenMP directives. Some investigations are currently done about how to best use the memory hierarchy within this code.

Adaptive solvers

We contributed to the development of the 2D adaptive MPI parallel code (referred to as Yoda4D). This code targets distributed memory architectures and clusters. We added a fast geometric mesh partitioner and data redistribution is now overlapped with computations in order to reduce overheads during the dynamic load balancing phase. This balancing scheme and the efficiency of the code up to 32 processors has been validated in [36] .

A lecture has been given by Martin Campos-Pinto [39] . A new class of adaptive semi-Lagrangian schemes-based on performing a semi-Lagrangian method on adaptive interpolation grids – in the context of solving Vlasov equations with underlying smooth flow, such as the one-dimensional Vlasov-Poisson system. After recalling the main features of the semi-Lagrangian method and its error analysis in a uniform setting, we describe two frameworks for implementing adaptive interpolations, namely multilevel finite elements and interpolating wavelets. For both discretizations we introduce a notion of good adaptivity to a given function, and show that it is preserved by a low-cost prediction algorithm which transport multilevel grids along any smooth flow. As a consequence, error estimates are established for the resulting “predict and re-adapt” schemes under the essential assumption that the flow underlying the transport equation, as well as its approximation, is a diffeomorphism. Some complexity results are stated in addition, together with a conjectured convergence rate for the overall adaptive scheme. As for the wavelet case, these results are new and also apply to high order interpolation.

Electromagnetic Particle In Cell (PIC) solvers

This project funded by ANR proposes to develop and compare Finite Element Time Domain (FETD) solvers based on the one hand on high order H(curl) conforming elements and on the other hand on high order Discontinuous Galerkin (DG) finite elements and investigate their coupling to the particles. These self consistent relativistic PIC solvers will be the first of this kind in this context and promise to have an impact for the simulation of realistic problems in accelerator and plasma physics.

We have developed a general mathematical framework for electromagnetic Particle-In-Cell (PIC) codes that can be used on structured as well as unstructured or hybrid grids. It can accommodate high order Maxwell finite element solvers conforming in H( curl ) that prevent to create spurious electromagnetic modes which could be amplified by the numerical noise of the PIC method (see [49] for high order Maxwell finite solvers and their properties). A 2D code, named HYBRID_FEM_PIC_2D, coupling a H( curl ) conforming finite element solver of any order 1, 2 or 3 on hybrid grid for the Maxwell equations and a PIC solver for the Vlasov equation was developed including an exact discrete continuity equation. We have also developed a 3D code on a regular grid (STRUCTURED_FEM_PIC_3D) including the same features and finally we are working on an unstructured one (UNSTRUCTURED_FEM_PIC_3D). One of the key points in the coupling between Maxwell and PIC solvers is the continuity equation. Indeed, the computation of the current density created by the particles and entering as the source for the electromagnetic solver must enforce the discrete charge conservation to ensure that the numerical solution is physical. This equation reads Im23 ${\mfrac {\#8706 \#961 }{\#8706 t}+divJ=0}$ , where $ \rho$ is the charge density and J the current density. We have shown that if the charge density is given by Im24 ${\#961 _h^n{(\#119857 )}:=\#8721 _{k=1}^N_pw_kS{(\#119857 -\#119857 _k^n)}\#916 t}$ , where S is a local smooth shape function such as a B-spline, and Im25 ${w_k,~\#119857 _k,~\#119855 _k}$ are respectively the weight, the position and velocity of the k -th particle, then the discrete continuity equation is guaranteed by the following formula for the current density:

Im26 ${\mover \#119817 ¯_h^{n+\mfrac 12}{(\#119857 )}:=\munderover \#8721 {k=1}N_pw_k\#119855 _k^{n+\mfrac 12}\mfrac 1{\#916 t}\#8747 _{t_n}^t_{n+1}S{(\#119857 -\#119857 _k{(t)})}dt.}$(3)

We have given a practical algorithm for implementing this new method, which requires the same computational tools as the classical PIC scheme. We also have shown in [47] that the lowest order version of our formulation (order 1 for the finite element) corresponds, on a rectangular grid, to the traditional Cloud-In-Cell (CIC) version of the PIC method coupled to a Yee Maxwell solver and the Villasenor-Buneman charge deposition scheme.

A parallel code (Spin) has been developed to solve bigger test cases using the same numerical method. Spin solves Maxwell equations on hybrid meshes (quadrangles and triangles) using edge finite elements. The work distribution on processors is expected to be balanced because the Maxwell solver on quadrangles is highly parallel.


Logo Inria