[Paper Review] Development of Massively Parallel Near Peak Performance Solvers for Three-Dimensional Geodynamic Modelling

Samuel Omlin, Räss, Ludovic+2 more|arXiv (Cornell University)|Jul 18, 2022

Drilling and Well EngineeringEngineering17 citations

TL;DR

This paper presents HPC.m, a MATLAB-based compiler that enables the development of massively parallel, near-peak performance solvers for 3D geodynamic simulations. It introduces parallel algorithms for modeling nonlinear porosity waves—both mechanical and reactive—demonstrating linear scaling on 80 GPUs and near-ideal scaling on 5,000 GPUs of Piz Daint, with soliton-like wave propagation observed across varying Damköhler numbers.

ABSTRACT

We address in this thesis the current need to design new parallel algorithms and tools that ease the development of geodynamic modelling applications that are suited for today's and tomorrow's hardware. We present (1) the MATLAB HPC compiler HPC.m, which greatly simplifies the building of parallel high performance applications and (2) parallel algorithms for the 3D simulation of strongly nonlinear processes as mechanical and reactive porosity waves. To simulate mechanical porosity waves we employ a massively parallel algorithm that permits to resolve the deformation of fluid-filled viscoelastic porous media in 3D. The utilized mathematical model is based on Biot's poroelastic theory, extended to account for viscous deformation and plastic yielding. The modelling results exhibit the impact of decompaction weakening on the formation of three-dimensional solitary-wave-like moving porosity channels. To simulate reactive porosity waves we use a solver for 3D deformation of fluid-filled reactive viscous porous media. The Damköhler number (Da) of the simulations is varied in order to estimate the respective roles of viscous deformation (low Da) and reaction (high Da) on wave propagation. 3D waves are found to propagate independently of their source at constant speed by going through each other for all the investigated Da. Soliton-like wave propagation as a result of metamorphic reaction provides an efficient mechanism for fluid flow in the Earth's crust. We illustrate the great performance and versatility of HPC.m by deploying it to generate solvers for a variety of physics across multiple Earth Science disciplines. All solvers run close to hardware's peak performance and were shown to scale linearly on a institute cluster with 80 GPUs. Moreover, our nonlinear poroviscoelastic two-phase flow solver scales close to ideally on Piz Daint's 5000 GPUs at the Swiss National Supercomputing Centre.

Motivation & Objective

• Address the growing need for high-performance, portable parallel geodynamic modeling tools tailored to modern heterogeneous hardware.
• Develop a domain-specific compiler (HPC.m) that simplifies the creation of high-performance parallel applications in MATLAB.
• Enable accurate, scalable simulation of complex 3D geodynamic processes such as mechanical and reactive porosity waves.
• Achieve near-peak performance and strong scaling on large GPU clusters for nonlinear poroviscoelastic flow problems.
• Investigate the role of viscous deformation versus chemical reaction in 3D porosity wave dynamics using variable Damköhler numbers.

Proposed method

• Implement HPC.m, a MATLAB HPC compiler that translates high-level MATLAB code into optimized, portable parallel code for GPU and distributed systems.
• Apply Biot's poroelastic theory extended with viscous deformation and plastic yielding to model fluid-filled viscoelastic porous media.
• Develop a nonlinear, two-phase flow solver for reactive porous media, incorporating reaction kinetics and variable Damköhler numbers (Da).
• Use a fully explicit, fully parallel finite difference scheme with domain decomposition for 3D simulations on distributed-memory GPU clusters.
• Employ a staggered grid discretization and iterative solvers optimized for memory bandwidth and load balancing.
• Validate and benchmark solvers on institutional GPU clusters and the Piz Daint supercomputer at CSCS, using up to 5,000 GPUs.

Experimental results

Research questions

• RQ1 Can a high-level compiler like HPC.m enable near-peak performance in 3D geodynamic simulations without low-level coding?
• RQ2 How do mechanical porosity waves evolve in 3D viscoelastic porous media under decompaction weakening?
• RQ3 What is the role of viscous deformation (low Da) versus chemical reaction (high Da) in the propagation of 3D reactive porosity waves?
• RQ4 Do 3D porosity waves exhibit soliton-like behavior, maintaining shape and speed when interacting?
• RQ5 To what extent can the solvers scale linearly on large GPU clusters, such as Piz Daint with 5,000 GPUs?

Key findings

• The HPC.m compiler enables the development of 3D geodynamic solvers that achieve near-peak performance on GPU hardware.
• The nonlinear poroviscoelastic solver scales linearly on an 80-GPU cluster, demonstrating strong portability and efficiency.
• On Piz Daint, the solver achieves near-ideal scaling across 5,000 GPUs, confirming its readiness for exascale geodynamic simulations.
• Mechanical porosity waves form 3D solitary-wave-like channels due to decompaction weakening, as predicted by the extended Biot model.
• Reactive porosity waves propagate at constant speed and maintain shape when crossing each other, indicating soliton-like behavior across all tested Damköhler numbers.
• The study confirms that metamorphic reactions can drive efficient, long-range fluid flow in the Earth's crust through stable, propagating porosity waves.