2020
Activity report
Project-Team
TONUS
RNSR: 201221129U
Research center
In partnership with:
CNRS, Université de Strasbourg
Team name:
TOkamaks and NUmerical Simulations
In collaboration with:
Institut de recherche mathématique avancée (IRMA)
Domain
Digital Health, Biology and Earth
Theme
Earth, Environmental and Energy Sciences
Creation of the Team: 2012 January 01, updated into Project-Team: 2019 May 01

# Keywords

• A6.1.1. Continuous Modeling (PDE, ODE)
• A6.1.4. Multiscale modeling
• A6.1.5. Multiphysics modeling
• A6.2.1. Numerical analysis of PDE and ODE
• A6.2.7. High performance computing
• A6.5.2. Fluid mechanics
• B4.2.2. Fusion

# 1 Team members, visitors, external collaborators

## Research Scientists

• Emmanuel Franck [Inria, Researcher]
• Victor Michel-Dansac [Inria, from Sep 2020, Starting Faculty Position]

## Faculty Members

• Philippe Helluy [Team leader, Univ de Strasbourg, Professor, HDR]
• Clementine Courtès [Univ de Strasbourg, Associate Professor]
• Michaël Gutnic [Univ de Strasbourg, Associate Professor]
• Laurent Navoret [Univ de Strasbourg, Associate Professor]
• Yannick Privat [Univ de Strasbourg, Professor, HDR]

## Post-Doctoral Fellow

• Pierre Gerhard [Univ de Strasbourg, from Sep 2020]

## PhD Students

• Mickael Bestard [Univ de Strasbourg, from Oct 2020]
• Leo Bois [Inria, from Sep 2020]
• Clement Flint [Univ de Strasbourg, from Nov 2020]
• Romane Helie [Univ de Strasbourg]
• Marie Houillon [Univ de Strasbourg, until Nov 2020]
• Guillaume Mestdagh [Univ de Strasbourg]
• Lucie Quibel [Univ de Strasbourg, until Sep 2020]

## Technical Staff

• Matthieu Boileau [CNRS, Engineer]
• Laura Mendoza [Inria, Engineer]

## Interns and Apprentices

• Mickael Bestard [Univ de Strasbourg, from Feb 2020 until Aug 2020]
• Leo Bois [Inria, from Feb 2020 until Aug 2020]

• Ouiza Herbi [Inria]

# 2 Overall objectives

TONUS started in January 2014. It is a team of the Inria Nancy-Grand Est center. It is located in the mathematics institute (IRMA) of the University of Strasbourg.

The International Thermonuclear Experimental Reactor (ITER) is a large-scale scientific experiment that aims to demonstrate that it is possible to produce energy from fusion, by confining a very hot hydrogen plasma inside a toroidal chamber, called tokamak. In addition to physics and technology research, tokamak design also requires mathematical modelling and numerical simulations on supercomputers.

The objective of the TONUS project is to deal with such mathematical and computing issues. We are mainly interested in kinetic, gyrokinetic and fluid simulations of tokamak plasmas. In the TONUS project-team we are working on the development of new numerical methods devoted to such simulations. We investigate several classical plasma models, study new reduced models and new numerical schemes adapted to these models. We implement our methods in two software projects: Selalib 1 and SCHNAPS 2 adapted to recent computer architectures.

We have strong relations with the CEA-IRFM team and participate in the development of their gyrokinetic simulation software GYSELA. We are involved in two Inria Project Labs, respectively devoted to tokamak mathematical modelling and high performance computing. The numerical tools developed from plasma physics can also be applied in other contexts. For instance, we collaborate with a small company in Strasbourg specialized in numerical software for applied electromagnetism. We also study kinetic acoustic models with the CEREMA and multiphase flows with EDF.

Finally, our topics of interest are at the interaction between mathematics, computer science, High Performance Computing, physics and practical applications.

# 3 Research program

## 3.1 Kinetic models for plasmas

The fundamental model for plasma physics is the coupled Vlasov-Maxwell kinetic model: the Vlasov equation describes the distribution function of particles (ions and electrons), while the Maxwell equations describe the electromagnetic field. In some applications, it may be necessary to take relativistic particles into account, which leads to consider the relativistic Vlasov equation, even if in general, tokamak plasmas are supposed to be non-relativistic. The distribution function of particles depends on seven variables (three for space, three for the velocity and one for time), which yields a huge amount of computation. To these equations we must add several types of source terms and boundary conditions for representing the walls of the tokamak, the applied electromagnetic field that confines the plasma, fuel injection, collision effects, etc.

Tokamak plasmas possess particular features, which require developing specialized theoretical and numerical tools.

Because the magnetic field is strong, the particle trajectories have a very fast rotation around the magnetic field lines. A full resolution would require a prohibitive amount of computation. It is necessary to develop reduced models for large magnetic fields in order to obtain tractable calculations. The resulting model is called a gyrokinetic model. It allows us to reduce the dimensionality of the problem. Such models are implemented in GYSELA and Selalib.

On the boundary of the plasma, the collisions can no more be neglected. Fluid models, such as MagnetoHydroDynamics (MHD) become again relevant. For the good operation of the tokamak, it is necessary to control MHD instabilities that arise at the plasma boundary. Computing these instabilities requires special implicit numerical discretizations with excellent long time behavior.

In addition to theoretical modelling tools, it is necessary to develop numerical schemes adapted to kinetic, gyrokinetic and fluid models. Three kinds of methods are studied in TONUS: Particle-In-Cell (PIC) methods, semi-Lagrangian and fully Eulerian approaches.

## 3.2 Gyrokinetic models: theory and approximation

In most phenomena where oscillations are present, we can establish a three-model hierarchy: $\left(i\right)$ the model parameterized by the oscillation period, $\left(ii\right)$ the limit model and $\left(iii\right)$ the two-scale model, possibly with its corrector. In a context where one wishes to simulate such a phenomenon where the oscillation period is small and the oscillation amplitude is not small, it is important to have numerical methods based on an approximation of the two-scale model. If the oscillation period varies significantly over the domain of simulation, it is important to have numerical methods that approximate properly and effectively the model parameterized by the oscillation period and the two-scale model. Implementing two-scale numerical methods (for instance by Frénod et al. 26) is based on a numerical approximation of the Two-Scale model. These are called of order 0. A Two-Scale Numerical Method is called of order 1 if it incorporates information from the corrector and from the equation of which this corrector is a solution. If the oscillation period varies between very small values and values of order 1, it is necessary to have new types of numerical schemes (Two-Scale Asymptotic Preserving Schemes of order 1 or TSAPS) that preserve the asymptotics between the model parameterized by the oscillation period and the Two-Scale model with its corrector. A first work in this direction has been initiated by Crouseilles et al. 23.

## 3.3 Semi-Lagrangian schemes

The Strasbourg team has a long and recognized experience in numerical methods for Vlasov-type equations. We are specialized in both particle and phase space solvers for the Vlasov equation: Particle-in-Cell (PIC) methods and semi-Lagrangian methods. We also have a long-standing collaboration with CEA Cadarache for the development of the GYSELA software for gyrokinetic tokamak plasmas.

The Vlasov and the gyrokinetic models are partial differential equations that express the transport of the distribution function in the phase space. In the original Vlasov case, the phase space is the six-dimension position-velocity space. For the gyrokinetic model, the phase space is five-dimensional because we consider only the parallel velocity in the direction of the magnetic field and the gyrokinetic angular velocity instead of three velocity components.

A few years ago, Eric Sonnendrücker and his collaborators introduced a new family of methods for solving transport equations in the phase space. This family of methods are the semi-Lagrangian methods. The principle of these methods is to solve the equation on a grid of the phase space. The grid points are transported with the flow of the transport equation for a time step and interpolated back periodically onto the initial grid. The method is then a mix of particle Lagrangian methods and Eulerian methods. The characteristics can be solved forward or backward in time leading to the Forward Semi-Lagrangian (FSL) or Backward Semi-Lagrangian (BSL) schemes. Conservative schemes based on this idea can be developed and are called Conservative Semi-Lagrangian (CSL).

GYSELA is a 5D full gyrokinetic code based on a classical backward semi-Lagrangian scheme (BSL) 32 for the simulation of core turbulence that has been developed at CEA Cadarache in collaboration with our team 27.

More recently, we have started to apply the semi-Lagrangian methods to more general kinetic equations. Indeed, most of the conservation laws of physics can be represented by a kinetic model with a small set of velocities Compressible fluids or MHD equations have such representations. Semi-Lagrangian methods then become a very appealing and efficient approach for solving these equations.

## 3.4 PIC methods

Historically PIC methods have been very popular for solving the Vlasov equations. They allow solving the equations in the phase space at a relatively low cost. The main disadvantage of this approach is that, due to its random aspect, it produces an important numerical noise that has to be controlled in some way, for instance by regularizations of the particles, or by divergence correction techniques in the Maxwell solver. We have a long-standing experience in PIC methods and we started implementing them in Selalib. An important aspect is to adapt the method to new multicore computers. See the work by Crestetto and Helluy 22.

## 3.5 Fluid and reduced kinetic models for plasmas

As already said, kinetic plasmas computer simulations are very intensive, because of the gyrokinetic turbulence. In some situations, it is possible to make assumptions on the shape of the distribution function that simplify the model. We obtain in this way a family of fluid or reduced models.

Assuming that the distribution function has a Maxwellian shape, for instance, we obtain the MagnetoHydroDynamic (MHD) model. It is physically valid only in some parts of the tokamak (at the edges for instance). The fluid model is generally obtained from the hypothesis that the collisions between particles are strong.

But the reduction is not necessarily a consequence of collisional effects. Indeed, even without collisions, the plasma may still relax to an equilibrium state over sufficiently long time scales (Landau damping effect).

In the fluid or reduced-kinetic regions, the approximation of the distribution function could require fewer data while still achieving a good representation, even in the collisionless regime.

Therefore, a fluid or a reduced model is a model where the explicit dependency on the velocity variable is removed. In a more mathematical way, we consider that in some regions of the plasma, it is possible to exhibit a (preferably small) set of parameters $\alpha$ that allows us to describe the main properties of the plasma with a generalized "Maxwellian” $M$. Then

In this case it is sufficient to solve for $\alpha \left(x,t\right)$. Generally, the vector $\alpha$ is the solution of a first order hyperbolic system.

Another way to reduce the model is to try to find an abstract kinetic representation with an as small as possible set of kinetic velocities. The kinetic approach has then only a mathematical meaning. It allows solving very efficiently many equations of physics.

## 3.6 Numerical schemes

As previously indicated, an efficient method for solving the reduced models is the Discontinuous Galerkin (DG) approach. It is possible to make it of arbitrary order. It requires limiters when it is applied to nonlinear PDEs occurring for instance in fluid mechanics. But the reduced models that we intend to write are essentially linear. The nonlinearity is concentrated in a few coupling source terms.

In addition, this method, when written in a special set of variables, called the entropy variables, has nice properties concerning the entropy dissipation of the model. It opens the door to constructing numerical schemes with good conservation properties and no entropy dissipation, as already used for other systems of PDEs 33, 20, 29, 28.

## 3.7 Matrix-free implicit schemes

In tokamaks, the reduced model generally involves many time scales. Among these time scales, many of them, associated to the fastest waves, are not relevant. In order to filter them out, it is necessary to adopt implicit solvers in time. When the reduced model is based on a kinetic interpretation, it is possible to construct implicit schemes that do not impose solving costly linear systems. In addition the resulting solver is stable even at a very high CFL (Courant Friedrichs Lax) number.

## 3.8 Electromagnetic solvers

Precise resolution of the electromagnetic fields is essential for proper plasma simulation. Thus it is important to use efficient solvers for the Maxwell systems and its asymptotics: Poisson equation and magnetostatics.

The proper coupling of the electromagnetic solver with the Vlasov solver is also crucial for ensuring conservation properties and stability of the simulation.

Finally, plasma physics implies very different time scales. It is thus very important to develop implicit Maxwell solvers and Asymptotic Preserving (AP) schemes in order to obtain good behavior on long time scales.

## 3.9 Coupling

The coupling of the Maxwell equations to the Vlasov solver requires some precautions. The most important one is to control the charge conservation errors, which are related to the divergence conditions on the electric and magnetic fields. We will generally use divergence correction tools for hyperbolic systems presented for instance in 19 (and the references therein).

## 3.10 Implicit solvers

As already pointed out, in a tokamak, the plasma presents several different space and time scales. It is not possible in practice to solve the initial Vlasov-Maxwell model. It is first necessary to establish asymptotic models by letting some parameters (such as the Larmor frequency or the speed of light) tend to infinity. This is the case for the electromagnetic solver and this requires implementing implicit time solvers in order to efficiently capture the stationary state, the solution of the magnetic induction equation or the Poisson equation.

# 4 Application domains

## 4.1 Controlled fusion and ITER

The search for alternative energy sources is a major issue for the future. Among others, controlled thermonu-clear fusion in a hot hydrogen plasma is a promising possibility. The principle is to confine the plasma in a toroidal chamber, called a tokamak, and to attain the necessary temperatures to sustain nuclear fusion reactions.The International Thermonuclear Experimental Reactor (ITER) is a tokamak being constructed in Cadarache, France. This was the result of a joint decision by an international consortium including the European Union, Canada, USA, Japan, Russia, South Korea, India and China. ITER is a huge project. As of today, the budget is estimated at 20 billion euros. The first plasma shot is planned for 2025 and the first deuterium-tritium operation for 2027. Many technical and conceptual difficulties have to be overcome before the actual exploitation of fusion energy. Consequently, much research has been carried out around magnetically confined fusion. Among these studies, it is important to carry out computer simulations of the burning plasma. Thus, mathematicians and computer scientists are also needed in the design of ITER. The reliability and the precision of numerical simulations allow a better understanding of the physical phenomena and thus would lead to better designs. TONUS’s main involvement is in such research. The required temperatures to attain fusion are very high, of the order of a hundred million degrees. Thus it is imperative to prevent the plasma from touching the tokamak inner walls. This confinement is obtained thanks to intense magnetic fields. The magnetic field is created by poloidal coils, which generate the toroidal component of the field. The toroidal plasma current also induces a poloidal component of the magnetic field that twists the magnetic field lines. The twisting is very important for the stability of the plasma. The idea goes back to research by Tamm and Sakharov, two Russian physicists, in the 50’s. Other devices are essential for the proper operation of the tokamak: a divertor for collecting the escaping particles, microwave heating for reaching higher temperatures, a fuel injector for sustaining the fusion reactions, toroidal coils for controlling instabilities, etc.

## 4.2 Other applications

The software and numerical methods that we develop can also be applied to other fields of physics or of engineering.

• For instance, we have a collaboration with the company AxesSim in Strasbourg for the development of efficient Discontinuous Galerkin (DG) solvers on hybrid computers (CPU/GPU). The applications are electro-magnetic simulations for the conception of antennas, electronic devices or aircraft electromagnetic compatibility.
• The acoustic conception of large rooms requires huge numerical simulations. It is not always possible to solve the full wave equation and many reduced acoustic models have been developed. A popular model consists in considering "acoustic" particles moving at the speed of sound. The resulting Partial Differential Equation (PDE) is very similar to the Vlasov equation. The same modelling is used in radiation theory. We have started to work on the reduction of the acoustic particles model and realized that our reduction approach perfectly applies to this situation. A PhD with CEREMA (Centre d’études et d’expertise sur les risques, l’environnement, la mobilité et l’aménagement) has started in October 2015 (PhD of Pierre Gerhard). The objective was to investigate the model reduction and to implement the resulting acoustic model in our DG solver.
• In September 2017, we started a collaboration with EDF Chatou (PhD of Lucie Quibel) on the modelling of multiphase fluids with complex equations of state. The goal is to simulate the high temperature liquid-vapor flow occurring in a nuclear plant. Among others, we will apply our recent kinetic method for designing efficient implicit schemes for this kind of flows.

# 5 New software and platforms

## 5.1 New software

### 5.1.1 CLAC

• Name: Conservation Laws Approximation on many Cores
• Keywords: Discontinuous Galerkin, GPU, Opencl, Maxwell equations
• Scientific Description:

It is clear now that future computers will be made of a collection of thousands of interconnected multicore processors. Globally it appears as a classical distributed memory MIMD machine. But at a lower level, each of the multicore processors is itself made of a shared memory MIMD unit (a few classical CPU cores) and a SIMD unit (a GPU). When designing new algorithms, it is important to adapt them to this kind of architecture. Our philosophy will be to program our algorithms in such a way that they can be run efficiently on this kind of computers. Practically, we will use the MPI library for managing the coarse grain parallelism, while the OpenCL library will efficiently operate the fine grain parallelism.

We have invested for several years until now into scientific computing on GPUs, using the open standard OpenCL (Open Computing Language). We were recently awarded a prize in the international AMD OpenCL innovation challenge thanks to an OpenCL two-dimensional Vlasov-Maxwell solver that fully runs on a GPU. OpenCL is a very interesting tool because it is an open standard now available on almost all brands of multicore processors and GPUs. The same parallel program can run on a GPU or a multicore processor without modification.

Because of the envisaged applications of CLAC, which may be either academic or commercial, it is necessary to conceive a modular framework. The heart of the library is made of generic parallel algorithms for solving conservation laws. The parallelism can be both fine-grained (oriented towards GPUs and multicore processors) and coarse-grained (oriented towards GPU clusters). The separate modules allow managing the meshes and some specific applications. In this way, it is possible to isolate parts that should be protected for trade secret reasons.

• Functional Description: CLAC is a generic Discontinuous Galerkin solver, written in C/C++, based on the OpenCL and MPI frameworks.
• URL:
• Contact: Philippe Helluy
• Partner: AxesSim

### 5.1.2 SCHNAPS

• Name: Solver for Conservative Hyperbolic Nonlinear Applications for PlasmaS
• Keywords: Discontinuous Galerkin, StarPU, Kinetic scheme
• Functional Description: Generic systems of conservation laws. Specific models: fluids, Maxwell, Vlasov, acoustics (with kinetic representation). Multitasking with StarPU. Explicit solvers (RK2, RK3, RK4): accelerated with OpenCL Implicit solvers: through kinetic representations and palindromic time integration.
• URL:
• Contact: Philippe Helluy
• Participants: Philippe Helluy, Matthieu Boileau, Bérenger Bramas

### 5.1.3 Patapon

• Name: Parallel Task in Python
• Keywords: Python, Parallel computing, High order time schemes
• Functional Description: Patapon is a code in PyOpenCL which allows to solve PDE like MHD using the vectorial Lattice Boltzmann method on Cartesian grids.
• Contact: Philippe Helluy
• Participant: Philippe Helluy

### 5.1.4 tofu

• Name: Tomography for Fusion
• Keywords: 3D, Data visualization, Visualization, Magnetic fusion, Tomography, Diagnostics, Plasma physics, Ray-tracing, Python
• Functional Description: tofu aims at providing the fusion and plasma community with an object-oriented, transparent and documented tool for designing tomography diagnostics, computing synthetic signal (direct problem) as well as tomographic inversions (inverse problem). It gives access to a full 3D description of the diagnostic geometry, thus reducing the impact of geometrical approximations on the direct and, most importantly, on the inverse problem.
• Release Contributions: Python 2.7 is not supported anymore. Python 3.6 and 3.7 are supported. Several changes to try and make installation easier (on clusters, windows, mac....) and less verbose for users. More explicit names for default saved configurations. Major bug fix in one of the methods for computing synthetic signal.Minor bug fixes in interactive figures. Minor bug fixes in Plasma2D interpolation. New configuration (ITER) available. First version of a class handling 2D XRay bragg spectrometers. First tools for magnetic field line tracing available on WEST. Better documentation, more resources. More informative error messages. Extra tools for computing LOS length, closest point to magnetic axis... Better PEP8 compliance
• URL:
• Contact: Laura Mendoza
• Partner: CEA

# 6 New results

## 6.1 Convolutional Neural Networks for PDEs

Participants: L. Bois, E. Franck, L. Navoret, V. Vigon (IRMA)

In this work, which was the object of a M2 internship, we have considered the Vlasov-Poisson equation for plasma physics. To reduce the CPU cost, it is usual to consider a model on the first moments in the velocity space (Euler equations). To obtain the final model it is important to add a closure (writing the last moment using the others). The classical closures are obtained by asymptotic analysis and are valid only on some regimes. In this work 13 we have proposed a closure based on a convolutional neural network call V-net and reference simulations obtained with a kinetic code. We have obtained a closure accuracy uniformly compared to the regime.

During two M1 internships, we have also tested the applications of Convolutional Neural Networks (CNNs) to the resolution of PDEs. The goal of the two M1 internships was to test the accuracy of convolutional networks like V-net to capture the solution of PDE problems for different applications. Regarding the M2 internship, its goal was to provide new closures for the Vlasov-Poisson equations.

In the first M1 internship, we proposed to construct a reduced model for tsunami propagation using data obtained solving the Shallow water equation in 1D. In the end, a first model, based on CNN, was able to predict the size and arrival time of the tsunami wave. In a second time, the CNN managed to predict the full wave which arrived on the beach. This work showed the interesting ability of CNNs to capture some physical effects and design cheaper qualitative models.

The second M1 internship deals with an inverse problem for the radiative transfer equations, which model the light propagation in a medium. In a first work, we consider a medium density constant and low with localized high density areas. Using simulations of the PDE, a first CNN model has been trained to localize the high density area using the results of light propagation only on the boundary. With a second neural network model, we have extended this method to the reconstruction of smooth densities. These two examples give a good idea of the CNN's ability to capture simplified physical solutions of the PDE.

## 6.2 Numerical methods for fluids and plasma dynamics

We have been working on building new numerical methods, as well as a simulation framework, for fluid and/or plasma dynamics. The development of numerical methods have been focused on the multi-scale aspect of the equations, and more specifically on the low-Mach and low-beta regimes. Namely, we constructed so-called asymptotic-preserving schemes, which behave well in every scale of the simulation.

### 6.2.1 Developement of a Python library for tomography diagnostics

Participants: L. Mendoza

The tofu (Tomography for Fusion) library (see section 5.1.4) aims at providing the fusion and plasma community with a tool for designing tomography diagnostics, computing synthetic signal (direct problem) as well as tomographic inversions (inverse problem). In order to compute the radiated power from the plasma on a certain volume received by a particle it is essential to accurately sample that volume. This procedure can be time- and memory-consuming. Most of the coding time was dedicated to this issue this year. The code was parallelized using OpenMP. Updates were made to use the latest and open source tools of CI/CD to ensure long-term maintenance and reliability. Finally, with the increasing number of users from different regions of the plasma community, we integrated more tokamak configurations to the code (DEMO, Asdex Upgrade, NSX, etc.)

### 6.2.2 Relaxation scheme for the Euler and MHD equations

Participants: F. Bouchut (CNRS & Gustave Eiffel University), C. Courtès, E. Franck, V. Michel-Dansac, L. Navoret

Previously, we have proposed implicit relaxation methods for fluid models that allow us to reverse a simple system. However, previous methods 21 were not very effective in the multi-scale regimes of interest. We have proposed in 2 a semi-implicit scheme based on a dynamic splitting and a relaxation allows to compute efficiently the low Mach limit of the Euler equations. During this year we have studied the MHD extension for plasmas simulations. First, we have written the relaxation model and proven the stability of the system.Second, we have proposed a spatial scheme for the relaxation which allows us to treat with a good accuracy the low Mach et low Beta flows. Currently we will try to write the splitting for the relaxation problem and the original equations.

### 6.2.3 Stable IMEX schemes for multi-scale equations

Participants: V. Michel-Dansac, A. Thomann (JGU Mainz, Germany)

When simulating the solutions to multi-scale problems, special care must be provided to the treatment of the different scales. This work aims at developing and applying IMEX (IMplicit-EXplicit) schemes to such problems, like the MHD system. The goal is twofold: developing generic schemes for a toy multi-scale problem, under physical constraints of stability and precision, and applying these schemes to complex systems. Following 25, 31, we have developed in 18 a new paradigm to build stable and accurate IMEX schemes for a toy scalar equation. The application to the Euler system, which represents flows in fluid dynamics, is underway and shows promising results.

### 6.2.4 Improvements and applications of the Discontinuous Galerkin method

Participants : Ph. Helluy, M. Houillon, P. Gerhard, L. Mendoza, V. Michel-Dansac, B. Bramas (CAMUS project-team)

The Discontinuous Galerkin method is a general method for solving conservation laws. This year, the end of the thesis of Marie Houillon culminated with the realization of a challenging numerical simulation on the Jean Zay supercomputer, using several hundreds of GPUs. The objective is to simulate the electromagnetic waves generated by a BlueTooth antenna around a human body: https://www.inria.fr/fr/objets-connectes-un-pas-de-geant-dans-la-modelisation-de-leurs-interactions-electromagnetiques-avec. We have also described some optimizations obtained with the StarPU runtime on our Maxwell solver 3. We have also obtained interesting results with our CFL-less explicit kinetic DG scheme applied to electromagnetic problems. A version of the software has been written in the language RUST. It is a very reliable language that allows us to avoid most of the common memory bugs at compile time. It also provides nice tools for automatic and robust parallelization.

### 6.2.5 Applications to plasma physics and multiphase flows

Participants : Ph. Helluy, L. Navoret, E. Franck, M. Boileau, R. Hélie, L. Quibel, H. Baty (ObAS, Strasbourg), C. Klingenberg (JMU Würzburg, Germany)

This year, we continued to explore the wide possibilities of kinetic schemes, with applications to such areas as plasma physics, multiphase flows and multi-scale problems. We first applied kinetic schemes to plasma physics 6 and to multiphase flows 7. The kinetic scheme is extremely simple and allows us to implement very efficient software. Indeed, with a Lattice Boltzmann version of the scheme, we were able to reach 80% of the maximal bandwidth of high-end GPU for computing MHD flows. This allowed us to propose simulations on extremely fine meshes. In 7, in the context of multiphase flows, we propose the simulation of the explosion of a hot-water pressurized vessel. The flow is complex, with shock waves, phase transitions and high variations of physical quantities. The kinetic scheme was able to capture the phenomena with good robustness and accuracy.

## 6.3 Other applications

While developing numerical tools mainly for plasma physics, it appears that these tools can also be used for other applications. We list below two of these applications.

### 6.3.1 Optimal control for population dynamics

Participants: L. Almeida (Sorbonne University), M. Duprez (MIMESIS project-team), R. Hélie, Y. Privat and N. Vauchelet (Paris 13 University)

We are interested in the analysis and simulation of solutions to optimal control problems motivated by population dynamics issues. In order to control the spread of mosquito-borne arboviruses, the population replacement technique consists in releasing into the environment mosquitoes infected with the Wolbachia bacterium, which greatly reduces the transmission of the virus to the humans. Spatial releases are then sought in such a way that the infected mosquito population invades the uninfected mosquito population. In order to determine the best temporal strategy for releases, we first considered time-dependent dynamic models 12, 11. Having made good progress on this issue, we now seek to determine time-space strategies on concrete territories.

We have introduced an appropriate model on the proportion of infected mosquitoes in time and space, and then an optimal control problem to determine the best spatial strategy to achieve these releases. This problem has been theoretically analysed, including the optimality of natural candidates and we have carried out first numerical simulations in one dimension of space to illustrate the relevance of our approach. This first work having been done 14, we are now seeking to design powerful numerical algorithms, allowing us to address this problem in dimensions 2 and 3 of space.

In parallel with this work, in the same spirit, we are seeking to determine how to optimally allocate resources in a habitat, when these are in limited quantities, in order to optimize the survival of the population, see 17, 16, 15

### 6.3.2 Cell motion as a stochastic process controlled by focal contacts dynamics

Participants: S. Lo Vecchio (IGBMC, Strasbourg), R. Thiagarajan (IGBMC, Strasbourg), D. Caballero (IGBMC, Strasbourg), V. Vigon (IGBMC, Strasbourg), L. Navoret, R. Voituriez (Sorbonne University), D. Riveline (IGBMC, Strasbourg)

In collaboration with physicists, we are interested in the individual movements of cells and how migration can result from asymmetries in the cell environment. The experiment consists in cells moving over triangular adhesion zones arranged in a line. Focal contacts of cells on these adhesion zones were detected and analysed. These data allowed us to validate in 5 a geometric model that predicts the direction of migration from the geometry of the accessible adhesive zones to the cells.

# 7 Bilateral contracts and grants with industry

## 7.1 Bilateral grants with industry

• Postdoc of Pierre Gerhard:
• Title: CFL-less explicit solver for the Maxwell equations
• Support: ANR DGA RAPID
• Company: ONERA Toulouse, AxesSim
• Duration: 09/2020 – 08/2021
• Participants: Ph. Helluy, V. Michel-Dansac, P. Gerhard
• Abstract: The main objective of this project is to implement in an efficient way the CFL-less explicit scheme developed in Tonus for several years and applied it to electromagnetic simulations. Another objective is to design new methods for applying boundary conditions and coupling with other solvers from the partners.
• PhD of Lucie Quibel:
• Title: Simulation d'écoulements diphasiques eau-vapeur en présence d'incondensables
• Support: CIFRE
• Company: EDF Chatou
• Duration: 10/2017 – 09/2020
• Participants: Ph. Helluy, L. Quibel
• Abstract: The objective of this thesis was to create new mathematical models for the simulation of two-phase flows with phase transition under two constraints: physical accuracy and numerical robustness.

# 8 Partnerships and cooperations

## 8.1 European initiatives

### 8.1.1 FP7 & H2020 Projects

• Strengthening the non-linear MHD code JOREK for application to key questions of the fusion roadmap:
• Acronym: H2020-EUROFusion-8824 / WP19 ENR-MFE19.MPG-03-T002-D001
• Duration: 01/2019 – 12/2020
• Coordinator: M. Hoelzl, Max Planck Institute for Plasma Physics, Garching
• Participants: E. Franck
• Abstract: In this project, the non-linear MHD code JOREK was extended beyond the state of the art. Extensive interpretation and validation work was done with experiments in the fields of ELMs and disruptions as well predictive simulations were done. A comprehensive overview article has been compiled describing physics models, numerical methods, verification, and the large number of applications.
• Mathematics and Algorithms for GYrokinetic and Kinetic models:
• Acronym: MAGYK
• Duration: 01/2019 – 06/2021
• Coordinator: E. Sonnendrücker, Max Planck Institute for Plasma Physics, Garching
• Participants: L. Navoret
• Abstract: This proposal is aimed at developing new models and algorithms that will be instrumental in enabling the efficient and reliable simulation of the full tokamak including the edge and scrape-off layer up to the wall with gyrokinetic or full kinetic models. It is based on a collaboration between applied mathematicians and fusion physicists that has already been very successful in a previous enabling research project and brings new ideas and techniques into the magnetic fusion community.

## 8.2 National initiatives

### 8.2.1 National projects

• Étude de stratégies multi-échelles pour la simulation en compatibilité électromagnétique:
• Acronym: ANR M2CEM
• Duration: 01/2021 – 12/2023
• Coordinator: X. Ferrières, ONERA Toulouse
• Partners: Tonus, IRMA, AxesSim, ONERA
• Participants: Ph. Helluy, V. Michel-Dansac, P. Gerhard
• Abstract: The objective of this project is to improve the numerical solvers of the partners for electromagnetic compatibility simulations. These tools are applied for the design of antennae or connected objects. Several aspects will be explored: code coupling, adaptative time stepping, CFL-less schemes, locally implicit schemes, innovative time integrators and mesh management.

# 9 Dissemination

## 9.1 Promoting scientific activities

### 9.1.1 Journal

#### Member of the editorial boards

• Editor of Computational and Applied Mathematics (Ph. Helluy)
• 05/2020 – now (Y. Privat): Member of the editorial board of the Journal of Optimization, Theory and Applications
• 11/2019 – now (Y. Privat): Member of the editorial board of the AIMS Applied Mathematics books series
• 04/2018 – now (Y. Privat): Member of the editorial board of the journal Evolution Equations and Control Theory

#### Reviewer - reviewing activities

• Proceedings of the Combustion Institute (M. Boileau)
• Engineering Computations (C. Courtès)
• European Journal of Mechanics / B Fluids (Ph. Helluy)
• Computers and Fluids (Ph. Helluy)
• Journal of Computational Physics (Ph. Helluy & L. Navoret)
• Shock Waves (SHOC) (Ph. Helluy)
• ESAIM proceedings 2019 (L. Navoret)
• Applied Mathematics and Computation (V. Michel-Dansac)
• Computers and Mathematics with Applications (V. Michel-Dansac)
• Journal of Fluid Mechanics (V. Michel-Dansac)
• zbMATH (V. Michel-Dansac)
• ESAIM Control Optimization and Calculus of Variations (Y. Privat)
• Journal de Mathématiques Pures et Appliquées (Y. Privat)
• SIAM Journal on Control and Optimization (Y. Privat)

### 9.1.2 Invited talks

• 12/2020: PDE and Applications Seminar, Poitiers, online (C. Courtès)
• 11/2020: PDE Seminar, Strasbourg (C. Courtès)
• 02/2020: 32nd CEA/GAMNI seminar on computational fluid dynamics, Paris (C. Courtès)
• 01/2020: Conférence itinérante du GDR EDP (C. Courtès)
• 06/2020: FVCA9, Bergen, Norway, online (E. Franck)
• “Journée scientifique du GdR MaNu 15 octobre 2020”, on Zoom...
• Oberwolfach Seminar: Structure-preserving Methods for Nonlinear Hyperbolic Problems, September 2020. Most of the talks were given on Zoom. (Ph. Helluy)
• 10/2020: PDE Seminar, Strasbourg (V. Michel-Dansac)
• 02/2020: 32nd CEA/GAMNI seminar on computational fluid dynamics, Paris (V. Michel-Dansac)

• Member of the managing committee of the “GDR Calcul” (M. Boileau)
• Member of the “Commission de la Recherche du Conseil Académique de l'Université de Strasbourg” (M. Boileau)
• Member of the CDT INRIA Nancy (M. Boileau)
• member of the jury (“comité de sélection”) for an associate professor position in Strasbourg (C. Courtès)
• Head of the mathematics department (IRMA) in Strasbourg University (M. Gutnic)
• Director of the CNRS Lab IRMA “Institut de Recherche Mathématique Avancée”, UMR 7501 (Ph. Helluy)
• 11/2019 – now: Member of CNU Section 26 (Y. Privat)
• 04/2019 – now: Member of the IRMA expert committee (Y. Privat)
• 10/2018 – now: Member of the IRMA scientific board (Y. Privat)

## 9.2 Teaching - Supervision - Juries

### 9.2.1 Teaching

• Licence 1, “Oral examinations”, 2h, Strasbourg University, France (C. Courtès)
• Licence 1, “Mathematics for biology”, 12h lectures + 23h exercise sessions, Strasbourg University, France (M. Gutnic)
• Licence 1, “Mathematics for health studies”, 36h, Strasbourg University, France (M. Gutnic)
• Licence 2, “Multivariable functions”, 12h, Strasbourg University, France (C. Courtès)
• Licence 2, “Applied numerical analysis”, 18h, Strasbourg University, France (C. Courtès)
• Licence 2, “Scientific computing”, 20h lectures + 17h computer labs, Strasbourg University, France (M. Gutnic)
• “Probability and Statistics”, 30h, CNAM, France (M. Gutnic)
• Licence 2, “Statistics for biologists”, 21.5h, Strasbourg University, France (L. Navoret)
• Licence 3, “Numerical analysis techniques 1”, 14 lectures (C. Courtès) + 18h computer labs (C. Courtès & E. Franck), Strasbourg University, France
• Licence 3, “Scientific computing”, 65h, Strasbourg University, France (C. Courtès)
• Licence 3, “Computer Science S6”, 20h lectures (E. Franck) + 34h computer labs (E. Franck & L. Navoret), Strasbourg University, France
• Licence 3, “Numerical analysis techniques 2”, 34h, Strasbourg University, France (Ph. Helluy)
• Licence 3, “Nonlinear Optimization”, 52h Strasbourg University, France (Y. Privat)
• Licence 3, “Introduction to C++”, 34h, Strasbourg University, France (V. Michel-Dansac)
• Master 1, “Parallel Computing”, 20h, Strasbourg University, France (M. Boileau)
• Master 1 Pure Mathematics, “Scientific computing”, 13h lectures (Ph. Helluy) + 14h computer labs (Ph. Helluy & V. Michel-Dansac), Strasbourg University, France
• Master 1 Applied Mathematics, “Scientific computing”, 35h, Strasbourg University, France (L. Navoret)
• Master 1, “Optimization”, 64h, Strasbourg University, France (Y. Privat)
• Master 2 Cell Physics, “Basics in mathematics”, 24h, Strasbourg University, France (L. Navoret)
• Master 2 Agrégation, “Scientific computing”, 20h, Strasbourg University, France (Ph. Helluy)
• Master 2 Agrégation, “Analysis”, 4.5h, Strasbourg University, France (L. Navoret)
• Master 2 Agrégation, “Scientific computing”, 66h, Strasbourg University, France (L. Navoret)
• Master 2 Agrégation, “Text study”, 30h, Strasbourg University, France (Y. Privat)
• Master 2, “Data Sciences”, 20h, Strasbourg University, France (M. Boileau)
• Master 2, “Numerical methods for hyperbolic PDEs”, 35h, Strasbourg University, France (Ph. Helluy)
• Master 2, “Shape optimization for fluids”, 15h, Sorbonne University, France (Y. Privat)
• Master 2, “Optimal control”, 64h, Strasbourg University, France (Y. Privat)
• Training Course, “Basics for scientific python”, 33h, Urfist (Unité Régionale de Formation à l'Information Scientifique et Technique), France (M. Boileau)

### 9.2.2 Supervision

• PhD in progress: M. Bestard, “Optimal control for numerical simulation in physics plasma.”, Strasbourg university. Beginning: September 2020. E. Franck, L. Navoret, Y. Privat.
• PhD in progress: L. Bois, “Machine learning for numericals methods. Application in plasma physics”, Strasbourg university. Beginning: September 2020. Ph. Helluy, E. Franck, L. Navoret, V. Vigon.
• PhD in progress: R. Hélie, “Relaxation methods for kinetic models in plasma physics”. Beginning: October 2019. Ph. Helluy, L. Navoret, E. Franck.
• PhD in progress: J. B. Arnau, “Stratégie de contrôle d’une population de moustiques pour la lutte contre les arbovirus.” Beginning: September 2019. Y. Privat, L. Almeida.
• PhD in progress: G. Mestdagh, “Suivi de tumeur en temps réel par des méthodes d’optimisation”, Strasbourg university. Beginning: September 2019. Y. Privat, S. Cotin.
• PhD in progress: A. Courtais, “Contrôle optimal de contacteurs à lit fixe par fabrication additive.” Beginning: September 2017. Y. Privat, F. Lesage and A. Latifi.
• PhD: A. Delyon, “Shape optimization problems around the geometry of Branchiopods eggs”, defended 29/10/2020, Y. Privat & A. Henrot, Lorraine University, 24
• PhD: P. Gerhard, “Réduction de modèles cinétiques et applications à l'acoustique du bâtiment ”, defended 28/01/2020, Ph. Helluy, 8
• PhD: M. Houillon, “Schémas Galerkin Discontinu optimisés pour les problèmes d'électromagnétisme avec des géométries complexes”, defended 19/11/2020, Ph. Helluy, 9
• PhD: I. Mazari, “Shape optimization and spatial heterogeneity in reaction-diffusion equations”, defended 06/07/2020, Y. Privat & G. Nadin, Sorbonne Université, 30
• PhD: L. Quibel, “Simulation d'écoulements diphasiques eau-vapeur en présence d'incondensables”, defended 18/09/2020, Ph. Helluy, 10

### 9.2.3 Juries

• president of the PhD jury of J.-B. Clément, University of Toulon (Ph. Helluy)
• member of the PhD jury of N. Dahmen, University of Toulouse (Ph. Helluy)
• reviewer in the PhD jury of T. Padioleau, University of Paris-Saclay (Ph. Helluy)
• reviewer of the PhD thesis of K. Mavreas, University of Nice (Y. Privat)
• member of the PhD jury of A. Léculier, Toulouse University (Y. Privat)
• member of the PhD jury of F. Vallet, Strasbourg University (Y. Privat)
• president of the PhD jury of J.N. Brunet, Lorraine University (Y. Privat)
• reviewer of the thesis of M. Boissier, École Polytechnique (Y. Privat)

## 9.3 Popularization

E. Franck and L. Navoret gave online talks for the “fête de la science”.

### 9.3.1 Articles and contents

• Article in the 2020 ECCOMAS newsletter: A Bernoulli's relation capturing scheme to simulate the 2011 Japan tsunami (V. Michel-Dansac)

# 10 Scientific production

## 10.1 Major publications

• 1 article ClémentineC. Courtès, DavidD. Coulette, EmmanuelE. Franck and LaurentL. Navoret. 'Vectorial kinetic relaxation model with central velocity. Application to implicit relaxations schemes'. Communications in Computational Physics 27 4 April 2020

## 10.2 Publications of the year

### International journals

• 2 articleFrançoisF. Bouchut, EmmanuelE. Franck and LaurentL. Navoret. 'A low cost semi-implicit low-Mach relaxation scheme for the full Euler equations'.Journal of Scientific Computing831April 2020, 24
• 3 articleBérengerB. Bramas, PhilippeP. Helluy, LauraL. Mendoza and BrunoB. Weber. 'Optimization of a discontinuous Galerkin solver with OpenCL and StarPU'.International Journal on Finite Volumes151January 2020, 1-19
• 4 article ClémentineC. Courtès, DavidD. Coulette, EmmanuelE. Franck and LaurentL. Navoret. 'Vectorial kinetic relaxation model with central velocity. Application to implicit relaxations schemes'. Communications in Computational Physics 27 4 April 2020
• 5 articleSimonS. Lo Vecchio, RaghavanR. Thiagarajan, DavidD. Caballero, VincentV. Vigon, LaurentL. Navoret, RaphaelR. Voituriez and DanielD. Riveline. 'Cell motion as a stochastic process controlled by focal contacts dynamics'.Cell Systems106June 2020, 535-542.e4

### Scientific book chapters

• 6 inbook RomaneR. Hélie, PhilippeP. Helluy, EmmanuelE. Franck and LaurentL. Navoret. 'Kinetic over-relaxation method for the convection equation with Fourier solver'. Finite Volumes for Complex Applications IX - Methods, Theoretical Aspects, Examples. FVCA 2020. June 2020
• 7 inbook PhilippeP. Helluy, OlivierO. Hurisse and LucieL. Quibel. 'Simulation of a liquid-vapour compressible flow by a Lattice Boltzmann Method'. Finite Volumes for Complex Applications IX - Methods, Theoretical Aspects, Examples. FVCA 2020. June 2020

### Doctoral dissertations and habilitation theses

• 8 thesis 'Reduction of kinetic models and their application to building acoustics'. Université de Strasbourg January 2020
• 9 thesis 'Optimized Discontinuous Galerkin Schemes for problems in electromagnetics with complex geometries'. Université de Strasbourg, IRMA UMR 7501 November 2020
• 10 thesis 'Simulation of water-vapor two-phase flows with non-condensable gas.'. Université de Strasbourg September 2020

### Reports & preprints

• 11 misc LuisL. Almeida, JesúsJ. Bellver Arnau, MichelM. Duprez and YannickY. Privat. 'Minimal cost-time strategies for mosquito population replacement'. November 2020
• 12 misc 'Optimal control strategies for the sterile mosquitoes technique'. 2020
• 13 misc LeoL. Bois, EmmanuelE. Franck, LaurentL. Navoret and VincentV. Vigon. 'A neural network closure for the Euler-Poisson system based on kinetic simulations'. October 2020
• 14 misc 'Optimization of spatial control strategies for population replacement, application to Wolbachia'. December 2020
• 15 misc AntoineA. Henrot, IdrissI. Mazari and YannickY. Privat. 'Shape optimization of a Dirichlet type energy for semilinear elliptic partial differential equations'. April 2020
• 16 misc IdrissI. Mazari, GrégoireG. Nadin and YannickY. Privat. 'Shape optimization of a weighted two-phase Dirichlet eigenvalue'. January 2020
• 17 misc IdrissI. Mazari, GrégoireG. Nadin and YannickY. Privat. 'Some challenging optimisation problems for logistic diffusive equations and numerical issues'. November 2020
• 18 misc VictorV. Michel-Dansac and AndreaA. Thomann. 'Large time step TVD IMEX Runge-Kutta schemes based on arbitrarily high order Butcher tableaux'. October 2020

## 10.3 Cited publications

• 19 inproceedings ChristophC. Altmann, ThomasT. Belat, MichaëlM. Gutnic, PhilippeP. Helluy, HélèneH. Mathis, EricE. Sonnendrücker, WilfredoW. Angulo and Jean-MarcJ.-M. Hérard. 'A local time-stepping Discontinuous Galerkin algorithm for the MHD system'. Modélisation et Simulation de Fluides Complexes - CEMRACS 2008 Marseille, France July 2009
• 20 articleT. Barth. 'On the role of involutions in the discontinous Galerkin discretization of Maxwell and magnetohydrodynamic systems'.IMA Vol. Math. Appl.1422006, 69–88
• 21 articleD. Coulette, E. Franck, Ph.P. Helluy, M. Mehrenberger and L. Navoret. 'High-order implicit palindromic discontinuous Galerkin method for kinetic-relaxation approximation'.Comput. & Fluids1902019, 485--502
• 22 inproceedingsAnaïsA. Crestetto and PhilippeP. Helluy. 'Resolution of the Vlasov-Maxwell system by PIC Discontinuous Galerkin method on GPU with OpenCL'.CEMRACS'1138FranceEDP Sciences2011, 257--274
• 23 articleNicolasN. Crouseilles, EmmanuelE. Frénod, Sever AdrianS. Hirstoaga and AlexandreA. Mouton. 'Two-Scale Macro-Micro decomposition of the Vlasov equation with a strong magnetic field'.Mathematical Models and Methods in Applied Sciences23082013, 1527--1559
• 24 phdthesis AlexandreA. Delyon. 'Shape optimisation problems around the geometry of Branchiopod eggs'. Université de Lorraine October 2020
• 25 articleG. Dimarco, R. Loubère, V. Michel-Dansac and M.-H. Vignal. 'Second-order implicit-explicit total variation diminishing schemes for the Euler system in the low Mach regime'.J. Comput. Phys.3722018, 178--201
• 26 articleEmmanuelE. Frénod, FrancescoF. Salvarani and EricE. Sonnendrücker. ' Long time simulation of a beam in a periodic focusing channel via a two-scale PIC-method'.Mathematical Models and Methods in Applied Sciences192ACM 82D10 35B27 76X052009, 175-197
• 27 articleVirginieV. Grandgirard, MauraM. Brunetti, PierreP. Bertrand, NicolasN. Besse, XavierX. Garbet, PhilippeP. Ghendrih, GiovanniG. Manfredi, YanickY. Sarazin, OlivierO. Sauter, EricE. Sonnendrücker, J. Vaclavik and LaurentL. Villard. 'A drift-kinetic Semi-Lagrangian 4D Vlasov code for ion turbulence simulation'.J. of Comput. Phys. 2172006, 395
• 28 articleC. Hauck and C.-D. Levermore. 'Convex Duality and Entropy-Based Moment Closures: Characterizing Degenerate Densities'.SIAM J. Control Optim.472008, 1977–2015
• 29 articleC.-D. Levermore. 'Entropy-based moment closures for kinetic equations'. Transport Theory Statist. Phys.264-51997, 591–606
• 30 phdthesis IdrissI. Mazari. 'Shape optimization and spatial heterogeneity in reaction-diffusion equations'. Sorbonne Université July 2020
• 31 incollection V. Michel-Dansac and A. Thomann. 'On high-precision ${L}^{}$-stable IMEX schemes for scalar hyperbolic multi-scale equations'. Proceedings of NumHyp 2019 SEMA SIMAI Springer Series Springer International Publishing 2019
• 32 articleEricE. Sonnendrücker, Jean-RodolpheJ.-R. Roche, PierreP. Bertrand and AlainA. Ghizzo. 'The semi-Lagrangian method for the numerical resolution of the Vlasov equation'.J. Comput. Phys.14921999, 201--220
• 33 incollectionE. Tadmor. 'Entropy conservative finite element schemes'.Numerical methods for Compressible Flows, Finite Difference Element and Volume TechniquesProc. Winter Annual Meeting, Amer. Soc. Mech. Eng, AMD- Vol. 781986, 149