High resolution propagation of seismic waves in space and time

January 10, 2017 - by Simone Ulmer

Geophysicists explore the Earth's interior by making use of use seismic waves, which may be caused by earthquakes or by artificial sources such as explosions. By using supercomputers to simulate the seismic waves’ behaviour in relation to their velocity of propagation, the scientists obtain a tomographic image of the Earth’s interior. They use this for exploring fault zones where earthquakes occur, as well as searching for mineral deposits.

Optimal time steps

The globe is diced up into a three-dimensional lattice for numerical calculation purposes. The chosen lattice spacing needs to be finer in places where there exist special spatial conditions such as major topographic differences, as well as to provide higher-resolution detail in geological fault zones. Computer scientist Maximilian Rietmann of ETH Zurich collaborated with geophysicist Daniel Peter, Assistant Professor at the King Abdullah University of Science and Technology (KAUST) in Saudi Arabia, Marcus Grote, Professor of Mathematics at the University of Basel and Olaf Schenk, Professor of Computational Science at Università della Svizzera italiana (USI) to develop a novel time stepping method. Their method runs on high-performance computers and is optimal for calculating wave propagation in Earth’s highly complex interior, say the researchers. The parallel-operating algorithm is the first to be capable of simulating seismic waves in differently sized lattice elements with varying, locally optimised time steps. This has sped up earthquake simulation on the CSCS supercomputer ‘Piz Daint’ by orders of magnitude: up to 100 times, depending on the model.

Earthquake simulation setup, with a Mw 4.6 earthquake in the Swiss-French Alps and two recording stations FUORN and SENIN. (Image: Rietmann et al.)

The researchers incorporated the new algorithm into what is known as SPECFEM code, which geophysicists have been using for around a decade to calculate tomographic images of Earth’s interior. The problem with most conventional codes for simulating wave propagation in the Earth's interior is that they allow selection of just one time step, regardless of the lattice spacing. For a stable calculation and correct simulation, the path traversed by the wave per time step must not exceed the lattice spacing. Therefore, if the lattice spacings differ, the maximum allowable time step must be chosen according to the smallest lattice spacing. So if the Earth is diced into millions of lattice elements including among them just a few thousand smaller ones, every one of the elements must nonetheless be calculated using the smallest time step. The simulation becomes correspondingly more time-consuming in result.

Flexible algorithm operation

Now the researchers have succeeded in constructing the algorithm such that simulation on a massively parallel supercomputer like ‘Piz Daint’ can use various time steps at once. So wave propagation in finely-spaced lattices can be calculated locally in small time steps (‘local time stepping’, LTS), while wave propagation in the surrounding coarse lattice elements is simulated in larger time steps. “The LTS method scales very well both on conventional CPUs and graphics processors”, says Rietmann. Rather than ten minutes, the calculation for a theoretical model takes just a few seconds using GPUs.

“We have redeveloped and optimised the wave propagation time scheme and made it very efficient in large simulations”, says geophysicist Daniel Peter. The researchers believe that many of the wave propagation codes currently operating with the traditional scheme stand to benefit. It is conceivable for the new algorithm to be used in any simulation of wave propagation in a medium that is aimed at creating a high-resolution spatial visualisation of its internal structure.

Switzerland example with localized refinement at source. (Image: Rietmann et al.)

The crucial breakthrough in this work is a first-ever solution to the load balancing problem by programming the algorithm to evenly distribute calculations for the differing lattice spacings among the processors, as Olaf Schenk and Marcus Grote point out.

Maximilian Rietmann earned his doctorate at the USI Institute for Computational Science, where he worked for around two years for the Platform for Advanced Scientific Computing (PASC) on the new algorithm and its implementation as part of the PASC GeoScale project ‘A framework for multiscale modelling and inversion’.

Further reading