Assessing the robustness and scalability of the accelerated pseudo-transient method towards exascale computing

The development of highly efficient, robust and scalable numerical algorithms lags behind the rapid increase in massive parallelism of modern hardware. We address this challenge with the accelerated pseudo-transient iterative method and present here a physically motivated derivation. We analytically determine optimal iteration parameters for a variety of basic physical processes and confirm the validity of theoretical predictions with numerical experiments. We provide an efficient numerical implementation of pseudo-transient solvers on graphical processing units (GPUs) using the Julia language. We achieve a parallel efficiency over 96 % on 2197 GPUs in distributed memory parallelisation weak scaling benchmarks. 2197 GPUs allow for unprecedented terascale solutions of 3D variable viscosity Stokes flow on 49953 grid cells involving over 1.2 trillion degrees of freedom. We verify the robustness of the method by handling contrasts up to 9 orders of magnitude in material parameters such as viscosity, and arbitrary distribution of viscous inclusions for different flow configurations. Moreover, we show that this method is well suited to tackle strongly nonlinear problems such as shear-banding in a visco-elasto-plastic medium. A GPU-based implementation can outperform CPU-based direct-iterative solvers in terms of wall-time even at relatively low resolution. We additionally motivate the accessibility of the method by its conciseness, flexibility, physically motivated derivation and ease of implementation. This solution strategy has thus a great potential for future high-performance computing applications, and for paving the road to exascale in the geosciences and beyond.

Achieving proportional fairness in WiFi networks via bandit convex optimization

In this paper, we revisit proportional fair channel allocation in IEEE 802.11 networks. Traditional approaches are either based on the explicit solution of the optimization problem or use iterative solvers to converge to the optimum. Instead, we propose an algorithm able to learn the optimal slot transmission probability only by monitoring the throughput of the network. We have evaluated this algorithm both (i) using the true value of the function to optimize and (ii) considering estimation errors. We provide a comprehensive performance evaluation that includes assessing the sensitivity of the algorithm to different learning and network parameters as well as its reaction to network dynamics. We also evaluate the effect of noisy estimates on the convergence rate and propose a method to alleviate them. We believe our approach is a practical solution to improve the performance of wireless networks as it does not require knowing the network parameters in advance. Yet, we conclude that the setting of the parameters of the algorithm is crucial to guarantee fast convergence.

A GPU-Accelerated Linear Solver for Massively Parallel Underground Simulations

Modern engineering applications require the solution of linear systems of millions or even billions of equations. The solution of the linear system takes most of the simulation for large scale simulations, and represent the bottleneck in developing scientific and technical software. Usually, preconditioned iterative solvers are preferred because of their low memory requirements and they can have a high level of parallelism. Approximate inverses have been proven to be robust and effective preconditioners in several contexts. In this communication, we present an adaptive Factorized Sparse Approximate Inverse (FSAI) preconditioner with a very high level of parallelism in both set-up and application. Its inherent parallelism makes FSAI an ideal candidate for a GPU-accelerated implementation, even if taking advantage of this hardware is not a trivial task, especially in the set-up stage. An extensive numerical experimentation has been performed on industrial underground applications. It is shown that the proposed approach outperforms more traditional preconditioners in challenging underground simulation, greatly reducing time-to-solution.

All-at-once multigrid approaches for one-dimensional space-fractional diffusion equations

We focus on a time-dependent one-dimensional space-fractional diffusion equation with constant diffusion coefficients. An all-at-once rephrasing of the discretized problem, obtained by considering the time as an additional dimension, yields a large block linear system and paves the way for parallelization. In particular, in case of uniform space–time meshes, the coefficient matrix shows a two-level Toeplitz structure, and such structure can be leveraged to build ad-hoc iterative solvers that aim at ensuring an overall computational cost independent of time. In this direction, we study the behavior of certain multigrid strategies with both semi- and full-coarsening that properly take into account the sources of anisotropy of the problem caused by the grid choice and the diffusion coefficients. The performances of the aforementioned multigrid methods reveal sensitive to the choice of the time discretization scheme. Many tests show that Crank–Nicolson prevents the multigrid to yield good convergence results, while second-order backward-difference scheme is shown to be unconditionally stable and that it allows good convergence under certain conditions on the grid and the diffusion coefficients. The effectiveness of our proposal is numerically confirmed in the case of variable coefficients too and a two-dimensional example is given.

Solution of a Linear Nondegenerate Matrix Equation Based on the Zero Divisor

New formulas were obtained to solve the linear non-degenerate matrix equations based on zero divisors of numerical matrices. Two theorems were formulated, and a proof to one of them is provided. It is noted that the proof of the second theorem is similar to the proof of the first one. The proved theorem substantiates new formula in solving the equation equivalent in the sense of the solution uniqueness to the known formulas. Its fundamental difference lies in the following: any explicit matrix inversion or determinant calculation is missing; solution is "based" not on the left, but on the right side of the matrix equation; zero divisor method is used (it was never used in classical formulas for solving a matrix equation); zero divisor calculation is reduced to simple operations of permutating the vector elements on the right-hand side of the matrix equation. Examples are provided of applying the proposed method for solving a nondegenerate matrix equation to the numerical matrix equations. High accuracy of the proposed formulas for solving the matrix equations is demonstrated in comparison with standard solvers used in the MATLAB environment. Similar problems arise in the synthesis of fast and ultrafast iterative solvers of linear matrix equations, as well as in nonparametric identification of abnormal (emergency) modes in complex technical systems, for example, in the power systems

Time-Domain Multidimensional Deconvolution: A Physically Reliable and Stable Preconditioned Implementation

Multidimensional deconvolution constitutes an essential operation in a variety of geophysical scenarios at different scales ranging from reservoir to crustal, as it appears in applications such as surface multiple elimination, target-oriented redatuming, and interferometric body-wave retrieval just to name a few. Depending on the use case, active, microseismic, or teleseismic signals are used to reconstruct the broadband response that would have been recorded between two observation points as if one were a virtual source. Reconstructing such a response relies on the the solution of an ill-conditioned linear inverse problem sensitive to noise and artifacts due to incomplete acquisition, limited sources, and band-limited data. Typically, this inversion is performed in the Fourier domain where the inverse problem is solved per frequency via direct or iterative solvers. While this inversion is in theory meant to remove spurious events from cross-correlation gathers and to correct amplitudes, difficulties arise in the estimation of optimal regularization parameters, which are worsened by the fact they must be estimated at each frequency independently. Here we show the benefits of formulating the problem in the time domain and introduce a number of physical constraints that naturally drive the inversion towards a reduced set of stable, meaningful solutions. By exploiting reciprocity, time causality, and frequency-wavenumber locality a set of preconditioners are included at minimal additional cost as a way to alleviate the dependency on an optimal damping parameter to stabilize the inversion. With an interferometric redatuming example, we demonstrate how our time domain implementation successfully reconstructs the overburden-free reflection response beneath a complex salt body from noise-contaminated up- and down-going transmission responses at the target level.