Recently, a lot of research has been done on the development of Neural Operators. The aim is to map infinite-dimensional functions to infinite-dimensional functions with neural networks. One approach is the Fourier Neural Operator [1]. We use this model to estimate the annual temperature based on the global land-sea-ice distribution. The input of the neural network is the earth’s geography, the output is the global temperature over one year. The training and testing data set is generated with the simulation framework “Klimakoffer” (for more information see Snapshot: Numerical simulations of earth’s climate and [2]).

The best results are achieved when using 4 Fourier integral operator layers, 16 Fourier modes, and a network width of 20. The model thus has 52 433 557 trainable parameters. The remaining specifications are as chosen by Li et al. [1]. It took 4 hours and 26 minutes to train the model on an NVIDIA GeForce RTX 3090 GPU and 0.77 seconds to evaluate it for one given example on a CPU. A relative L2 loss of approximately 0.0025 was achieved during training.

The first video shows one exemplary output. The second animation visualizes the difference between the temperature calculated by the FNO model and the “true” temperature determined by the framework “Klimakoffer”.

References:
[1] Li, Z., Kovachki, N., Azizzadenesheli, K., Liu, B., Bhattacharya, K., Stuart, A., Anandkumar, A. (2020). “Fourier Neural Operator for Parametric Partial Differential Equations”, arXiv, https://arxiv.org/abs/2010.08895
[2] Zhuang, K., North, G.R., Stevens, Mark J. (2017). “A NetCDF version of the two-dimensional energy balance model based on the full multigrid algorithm”, SoftwareX

To simulate elastic structures, we use a total Lagrangian SPH formulation [1].
We simulate a dam break against a flexible structure based on the setup given by [1].
For the fluid dynamics, a weakly-compressible Lagrangian SPH formulation with an artificial viscosity by Monaghan [2] and boundary force particles to model the tank [3] was used.
This animation shows the simulation in 0.1x real time with 8.7k fluid particles and 455 particles to model the flexible structure.

References:
[1] O’Connor, J., Rogers, B. D. (2021). “A fluid–structure interaction model for free-surface flows and flexible structures using smoothed particle hydrodynamics on a GPU.” In: Journal of Fluids and Structures, 104. https://doi.org/10.1016/J.JFLUIDSTRUCTS.2021.103312
[2] Monaghan, J. J. (1989). “On the Problem of Penetration in Particle Methods”. In: Journal of Computational Physics, 82(1), 1–15. https://doi.org/10.1016/0021-9991(89)90032-6
[3] Monaghan, J. J., & Kajtar, J. B. (2009). “SPH particle boundary forces for arbitrary boundaries”. In: Computer Physics Communications, 180(10), 1811–1820. https://doi.org/10.1016/j.cpc.2009.05.008

We solve the compressible Euler equations of gas dynamics with a fourth-order accurate entropy-stable discontinuous Galerkin (DG) method and combine it with a first-order accurate finite volume (FV) method at the node level to impose positivity of density and pressure [1,2].
The initial condition of this problem is given by:
\[
\begin{array}{rlrl}
\rho (t=0) &= \frac{1}{2}
+ \frac{3}{4} B,
&
p (t=0) &= 0.1, \\
v_1 (t=0) &= \frac{1}{2} \left( B-1 \right),
&
v_2 (t=0) &= \frac{1}{10} \sin(2 \pi x),
\end{array}
\]
with $B=\tanh \left( 15 y + 7.5 \right) – \tanh(15y-7.5)$, where $\rho$ is the density, $\vec{v}=(v_1,v_2)$ is the velocity, and $p$ is the pressure.

The video shows the evolution of the density in time using $1024^2$ and degrees of freedom. We use the entropy-conservative and kinetic energy preserving flux of Ranocha [3] for the volume fluxes and the Rusanov solver for the surface fluxes of the DG and FV methods. During this simulation, the FV method acts on average in 0.000086% of the computational domain, and never more than 0.002% of the computational domain at a specific time.

References
[1] A. M. Rueda-Ramírez, G. J. Gassner, A Subcell Finite Volume Positivity-Preserving Limiter for DGSEM Discretizations of the Euler Equations, WCCM-ECCOMAS2020, pp. 1–12.
[2] A. M. Rueda-Ramírez, W. Pazner, G. J. Gassner, Subcell limiting strategies for discontinuous galerkin spectral element methods, Computers & Fluids 247 (2022) 105627.
[3] H. Ranocha, Generalised summation-by-parts operators and entropy stability of numerical methods for hyperbolic balance laws, Cuvillier Verlag, 2018.

We solve the multi-component magneto-hydrodynamics (MHD) equations using a hybrid finite volume (FV)/discontinuous Galerkin (DG) method [1], which combines the two discretization approaches at the node level. To test the robustness of the hybrid FV/DG method, we compute a modification of the Orszag-Tang vortex problem [2], where we use two ion species and initialize the flow with a sharp interface at $y=0$:
\[
\begin{array}{rlrl}
\rho^1(x,y,t=0) &= \begin{cases}
{25}/{36 \pi} & y \ge 0 \\
10^{-8} & y < 0
\end{cases} ,
&\rho^2(x,y,t=0) &= \begin{cases}
10^{-8} & y < 0\\
{25}/{36 \pi} & y \ge 0
\end{cases} , \\
v_1 (x,y,t=0) &= – \sin (2 \pi y),
&v_2 (x,y,t=0) &= \sin (2 \pi x), \\
B_1 (x,y,t=0) &= -\frac{1}{\sqrt{4 \pi}} \sin (2 \pi y),
&B_2 (x,y,t=0) &= -\frac{1}{\sqrt{4 \pi}} \sin (4 \pi x). \\
p (x,y,t=0) &= \frac{5}{12 \pi},
& & \\
\end{array}
\]
where $\rho^1$ is the density of the first ion species, $\rho^2$ is the density of the second ion species,$\vec{v}=(v_1,v_2)$ is the plasma velocity, $\vec{B}=(B_1,B_2)$ is the magnetic field, and $p$ is the plasma pressure (without magnetic pressure).
 
The fourth-order DG method delivers high accuracy, but fails to describe shocks and material interfaces. Therefore, we use a modal indicator [3] to detect the elements where the solution changes abruptly. In those elements, we combine the DG method with a robust lower order FV scheme locally (at the node level) to impose TVD-like conditions on the density of the ion species,
\[
\min_{k \in \mathcal{N} (j)} {\rho}^{i,\mathrm{FV}}_{k}
\le \rho^i_{j} \le
\max_{k \in \mathcal{N} (j)} {\rho}^{i,\mathrm{FV}}_{k},
\]
and a local minimum principle on the specific entropy,
\[
s({\mathbf{u}}_{j}) \ge \min_{k \in \mathcal{N} (j)} s(\mathbf{u}^{\mathrm{FV}}_{k}),.
\]
where $\mathcal{N} (j)$ denotes the collection of nodes in the low-order stencil of each node $j$.
 
The video shows the total density ($\rho=\rho^1+\rho^2$), the density of the first ion species ($\rho^1$), and the so-called blending coefficient ($\alpha$) for the combination of a fourth-order accurate entropy stable DG method with first- and second-order accurate FV methods (Rusanov solver) and $1024^2$ degrees of freedom. A blending coefficient $\alpha=0$ means that only the DG is being used, and a blending coefficient $\alpha=1$ means that only the FV method is being used. The example shows that the hybrid FV/DG scheme is able to capture the shocks of the simulation correctly, deal with near-vacuum conditions and sharp interfaces.

[1] Rueda-Ramírez, A. M.; Pazner, W., & Gassner, G. J. (2022). Subcell limiting strategies for discontinuous Galerkin spectral element methods. Computers & Fluids. arXiv:2202.00576.
[2] Orszag, S.; Tang, C (1979). Small-scale structure of two-dimensional magnetohydrodynamic turbulence. Journal of Fluid Mechanics.
[3] Persson, P. O.; Peraire, J. (2006). Sub-cell shock capturing for discontinuous Galerkin methods. In 44th AIAA Aerospace Sciences Meeting and Exhibit (p. 112).

To study the Smoothed Particle Hydrodynamics method, we simulate a dam break based on the setup given by [1].
We use a weakly-compressible Lagrangian SPH formulation with an artificial viscosity by Monaghan [2] and boundary force particles to model the tank [3].

This animation shows the simulation in 0.5x real time with 125k fluid particles and 7k boundary particles.

References:
[1] Marrone, S., Antuono, M., Colagrossi, A., Colicchio, G., le Touzé, D., & Graziani, G. (2011). “δ-SPH model for simulating violent impact flows”. In: Computer Methods in Applied Mechanics and Engineering, 200(13–16), 1526–1542. https://doi.org/10.1016/J.CMA.2010.12.016
[2] Monaghan, J. J. (1989). “On the Problem of Penetration in Particle Methods”. In: Journal of Computational Physics, 82(1), 1–15. https://doi.org/10.1016/0021-9991(89)90032-6
[3] Monaghan, J. J., & Kajtar, J. B. (2009). “SPH particle boundary forces for arbitrary boundaries”. In: Computer Physics Communications, 180(10), 1811–1820. https://doi.org/10.1016/j.cpc.2009.05.008

We investigate the Kelvin-Helmholtz instability setup of Fjordholm et al. [1], where the sharp interface in the initial condition is randomly perturbed. The simulation is done with an entropy stable DGSEM with polynomial degree N=3 and 32 x 32 grid cells. The final time is t = 2.0. Trixi.jl [2] is used with the setup found in examples/tree_2d_dgsem/elixir_euler_kelvin_helmholtz_instability_fjordholm_etal.jl.

The following gif is showing 100 simulation results of the density distribution for different random interface perturbations:

Here is the average density distribution of all 100 results:

[1] Ulrik S. Fjordholm, Roger Käppeli, Siddhartha Mishra, Eitan Tadmor. Construction of approximate entropy measure valued solutions for hyperbolic systems of conservation laws, 2014. (https://arxiv.org/abs/1402.0909)
[2] Trixi.jl, a numerical simulation framework for hyperbolic conservation laws written in Julia. (https://github.com/trixi-framework/Trixi.jl)

We show a detonation diffraction simulation with multiple obstacles in a reflective box using reactive multicomponent Euler equations with Trixi.jl. These kind of detonation waves may result in pressure and density drops close to zero which is numerically difficult not only but especially for high order schemes.

In this simulation we use a reacting model which consists of one reaction and two species. This example has been calculated using the high order DG-FV blending method in Trixi.jl with a HOHQMesh generated mesh whereas the chemical network is solved with KROME.

(Note: Not all of these features are currently available in the main code of Trixi.jl.)