Research

Numerical methods for fluids and plasma dynamics

CFL-less and parallel Discontinuous Galerkin solver

The Discontinuous Galerkin method is a general method for solving conservation laws. In 8, we described a parallel and quasi-explicit Discontinuous Galerkin (DG) kinetic scheme for solving systems of balance laws. The solver is unconditionally stable (i.e., the CFL number can be arbitrary) and has the complexity of an explicit scheme. It can be applied to any hyperbolic system of balance laws. In this work, we assessed the performance of the scheme in the particular case of the three-dimensional wave equation and of Maxwell’s equations. We measured the benefit of the unconditional stability by performing experiments with very large CFL numbers. In addition, we investigated how to parallelize this method.

A version of this solver (presented in section 5.1.1), called KOUGLOFV, was implemented in the RUST language. 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 shared-memory parallelization. This parallelization was tested in the code, and good efficiency results were obtained. In 2022, KOUGLOFV was enhanced with distributed-memory parallelization (via MPI) and with the capability to handle large-scale simulations of electromagnetic waves within the human body, see figure 1. This led to the preprints 18 as well as the short report 19; the research was performed in the context of the grant detailed in section 7.1.

Generic high-order well-balanced numerical schemes

In 5, we considered a general framework to build a linear high-order well-balanced scheme for systems of balance laws. Such a framework is suited to almost any equation describing the motion of a fluid or a plasma. This article deals with a well-known issue of high-order well-balanced schemes. Indeed, such high-order schemes are based on a polynomial reconstruction, which must preserve the steady states under consideration in order to get the required well-balancedness property. A priori, to capture such a steady state, one needs to solve some strongly nonlinear equations. Here, a very easy, linear correction is designed under the generic framework of a system of hyperbolic balance laws, which describe most of the fluid or plasma systems. This correction can be applied to any scheme of order greater than or equal to 2, such as a MUSCL-type scheme, and ensures that this scheme exactly preserves the steady solutions. The main discrepancy with usual techniques lies in never having to invert the nonlinear function governing the steady solutions.

Third-order asymptotic-preserving scheme for the isentropic Euler equations

In 11, we developed of first order total variation diminishing (TVD) implicit-explicit (IMEX) Runge-Kutta (RK) schemes as a basis of a Multidimensional Optimal Order detection (MOOD) approach to approximate the solution of hyperbolic multi-scale equations. A key feature of our newly proposed TVD schemes is that the resulting CFL condition does not depend on the fast waves of the considered model, as long as they are integrated implicitly. However, a result from Gottlieb et al. 26 gives a first order barrier for unconditionally stable implicit TVD-RK schemes and TVD-IMEX-RK schemes with scale-independent CFL conditions. Therefore, the goal of this work is to consistently improve the resolution of a first-order IMEX-RK scheme, while retaining its $L^1$ stability and TVD properties. In this work we present a new approach based on a convex combination between a first-order TVD IMEX Euler scheme and a potentially oscillatory high-order IMEX-RK scheme. We derive and analyse the TVD property for a scalar multi-scale equation. In addition, we numerically assess the performance of our TVD schemes compared to standard $L$ -stable and SSP IMEX RK schemes from the literature. Finally, the resulting TVD-MOOD schemes are applied to the isentropic Euler equations, with very convincing results, see figure 2.

Optimal design of stellarators

In 12, we have been interested, with M. Sigalotti and R. Robin, in the optimal design of stellarators, devices for the production of controlled nuclear fusion reactions, alternative to tokamaks. The confinement of the plasma is entirely achieved by a helical magnetic field created by the complex arrangement of coils powered by high currents around a toric domain. These coils describe a surface called “coil winding surface” (CWS). We modeled the design of the CWS as a shape optimization problem, so that the cost function reflects both the optimal plasma properties, through a least squares functional, and also manufacturability, through geometric terms involving the lateral surface and curvature of the CWS.

Other applications

While the main focus of the numerical tools we develop is plasma physics, they can also be used for other applications. We list below four such applications.

Reduced modeling and optimal control of epidemiological individual-based models with contact heterogeneity

Modeling epidemics using classical population-based models suffers from shortcomings that so-called individual-based models are able to overcome. They are able to take into account heterogeneity features, such as super-spreaders, and describe the dynamics involved in small clusters. In return, such models often involve large graphs which are expensive to simulate and difficult to optimize, both in theory and in practice. By combining the reinforcement learning philosophy with reduced models, we propose in 17 a numerical approach to determine optimal health policies for a stochastic individual-based model taking into account heterogeneity in the population. More precisely, we introduce a deterministic reduced population-based model involving a neural network, designed to faithfully mimic the local dynamics of the more complex individual-based model. Then the optimal control is determined by sequentially training the network until an optimal strategy for the population-based model succeeds in also containing the epidemic when simulated on the individual-based model. After describing the practical implementation of the method, several numerical tests are proposed to demonstrate its ability to determine controls for models with contact heterogeneity.

Optimal control for population dynamics

In collaboration with L. Almeida, J. Bellver Arnau, M. Duprez, G. Nadin, I. Mazari and N. Vauchelet, we pursued a series of works dedicated to the analysis and simulation of solutions of an optimal control problem motivated by population dynamics issues. In order to control the spread of mosquito-borne arboviruses, the sterile insect technique (SIT) consists in releasing mosquitoes infected with a bacterium called Wolbachia into the environment, which considerably reduces the transmission of the virus to humans. The goal is to effectively release the mosquitoes spatially so that the population of infected ones overwhelms the population of uninfected mosquitoes. Assuming very high fecundity rates, an asymptotic model on the proportion of infected mosquitoes is introduced, leading to an optimal control problem to determine the best spatial strategy to adopt.

We tackled the optimal strategy problem for SIT in 3 and 4, more specifically studying the robustness of optimal strategy for a general family of criteria and introduced an adapted optimization algorithm for computing optimal control strategies.

In addition, 10 is dedicated to the study of some qualitative properties of optimal control problems involving weighted eigenvalues and diffusion reaction systems used to describe the spatio-temporal dynamics of mosquito vectors of diseases such as dengue. This work, although fundamental and not dedicated to any particular application, have notably led to efficient optimization algorithms and analysis for the mentioned applied problems.

Gridless 3D Recovery of Image Sources from Room Impulse Responses

In 13, we were interested in the reconstruction of the shape of a parallelepipedic room from measurements made by microphones during a certain time. A difficulty of this work was the fact that only the low frequency part of the signal is measured. We introduced a model and an efficient reconstruction algorithm for such shapes. We now seek to generalize our approach to any shape of room.

Other applications in (shape) optimization and optimal control

In 6, algorithms specific to a particular application, epidemiology or quantum chemistry, were obtained using a precise analysis of the model and the optimality conditions.

Finally, 9 contains fine analysis results of the “geometric quantity” for the wave equation with an internal control, corresponding to the minimum time taken on average by a geodesic to intersect the control set.