Fantastic linear systems and how to solve them - Pleasingly parallel computation of resolution and uncertainty in seismic tomography

In classic global seismic tomography one ultimately needs to solve a large, typically sparse, linear system. Specific formulations of the tomographic problem then allow for computing the resulting model, as well as resolution and uncertainty. Although the latter two properties are crucial for the assessment of a given tomographic image, they are often not calculated because of high computational costs.

Based on an actual data set used for global tomography, we discuss the usual and an alternative tomographic system and investigate the possibility of embarrassingly parallel (or pleasingly, perfectly parallel) computations of the full model result including resolution and uncertainty. The necessary numerical strategies and computational concepts are briefly explained from a practical point of view. This includes possible iterative and direct numerical solvers like, for example, LSQR and QR decomposition, a discussion of multiprocessing vs. multithreading and first personal experience with GPU computing using the programming language Julia.