Mixed precision iterative refinement with adaptive precision sparse approximate inverse preconditioning

This paper derives tradeoffs between three basic costs of a parallel algorithm: synchronization, data movement, and computational cost. These tradeoffs are lower bounds on the execution time of the algorithm which are independent of the number of processors, but dependent on the problem size. Therefore, they provide lower bounds on the parallel execution time of any algorithm computed by a system composed of any number of homogeneous components, each with associated computational, communication, and synchronization payloads. We employ a theoretical model counts the amount of work and data movement as a maximum of any execution path during the parallel computation. By considering this metric, rather than the total communication volume over the whole machine, we obtain new insights into the characteristics of parallel schedules for algorithms with non-trivial dependency structures. We also present reductions from BSP and LogP algorithms to our execution model, extending our lower bounds to these two models of parallel computation. We first develop our results for general dependency graphs and hypergraphs based on their expansion properties, then we apply the theorem to a number of specific algorithms in numerical linear algebra, namely triangular substitution, Gaussian elimination, and Krylov subspace methods. Our lower bound for LU factorization demonstrates the optimality of Tiskin's LU algorithm answering an open question posed in his paper, as well as of the 2.5D LU algorithm which has analogous costs. We treat the computations in a general manner by noting that the computations share a similar dependency hypergraph structure and analyzing the communication requirements of lattice hypergraph structures.

Matrix computations lie at the heart of many scientific computations. While sophisticated algorithms have been established for various numerical linear algebra problems such as the solution of linear simultaneous equations and eigenvalue problems, they require considerable modification with the advent of exaFLOPS- scale supercomputers, which are expected to have a huge number of computing cores, deep memory hierarchy, and increased probability of hardware errors. In this chapter, we discuss projected hardware characteristics of exaFLOPS machines and summarize the challenges to be faced by numerical linear algebra algorithms in the near future. Based on these preparations, we present a brief survey of recent research efforts in the field of numerical linear algebra targeted at meeting these challenges.

Eigenstructure of order-one-quasiseparable matrices. Three-term and two-term recurrence relations

We investigate randomized distributed algorithms for matrix computations over loosely coupled distributed systems, such as P2P networks or sensor networks. In this poster, we discuss orthogonalization methods and orthogonal iteration. These algorithms are very well understood in the sequential or in the classical parallel context, and they are important building blocks for many algorithms in numerical linear algebra.

Efficient approximate solution of sparse linear systems

The development of efficient, general-purpose software for the iterative solution of sparse linear systems on parallel MIMD computers depends on recent results from a wide variety of research areas. Parallel graph heuristics, convergence analysis, and basic linear algebra implementation issues must all be considered.

SummaryAs the computational power of high‐performance computing systems continues to increase by using a huge number of cores or specialized processing units, high‐performance computing applications are increasingly prone to faults. In this paper, we present a new class of numerical fault tolerance algorithms to cope with node crashes in parallel distributed environments. This new resilient scheme is designed at application level and does not require extra resources, that is, computational unit or computing time, when no fault occurs. In the framework of iterative methods for the solution of sparse linear systems, we present numerical algorithms to extract relevant information from available data after a fault, assuming a separate mechanism ensures the fault detection. After data extraction, a well‐chosen part of missing data is regenerated through interpolation strategies to constitute meaningful inputs to restart the iterative scheme. We have developed these methods, referred to as interpolation–restart techniques, for Krylov subspace linear solvers. After a fault, lost entries of the current iterate computed by the solver are interpolated to define a new initial guess to restart the Krylov method. A well‐suited initial guess is computed by using the entries of the faulty iterate available on surviving nodes. We present two interpolation policies that preserve key numerical properties of well‐known linear solvers, namely, the monotonic decrease of the A‐norm of the error of the conjugate gradient or the residual norm decrease of generalized minimal residual algorithm for solving. The qualitative numerical behavior of the resulting scheme has been validated with sequential simulations, when the number of faults and the amount of data losses are varied. Finally, the computational costs associated with the recovery mechanism have been evaluated through parallel experiments. Copyright © 2016 John Wiley & Sons, Ltd.

Singular value decomposition (SVD) is a key step in many algorithms in statistics, machine learning and numerical linear algebra. While classical singular value decomposition has been made efficient in terms of computational complexity, classical algorithms are not able to fully utilise modern computing environments. The goal of this work is to survey various implementations and applications of randomized algorithms for SVD. Algorithms are compared in terms of accuracy and time of execution. On example of robust principal component analysis (RPCA), it is shown that using randomized algorithms can yield a significant speedup for image processing and similar applications.

This survey highlights the recent advances in algorithms for numerical linear algebra that have come from the technique of linear sketching, whereby given a matrix, one first compresses it to a much smaller matrix by multiplying it by a (usually) random matrix with certain properties. Much of the expensive computation can then be performed on the smaller matrix, thereby accelerating the solution for the original problem. In this survey we consider least squares as well as robust regression problems, low rank approximation, and graph sparsification. We also discuss a number of variants of these problems. Finally, we discuss the limitations of sketching methods.

The multiplication of a sparse matrix by a dense vector (SpMV) is a centerpiece of scientific computing applications: it is the essential kernel for the solution of sparse linear systems and sparse eigenvalue problems by iterative methods. The efficient implementation of the sparse matrix-vector multiplication is therefore crucial and has been the subject of an immense amount of research, with interest renewed with every major new trend in high-performance computing architectures. The introduction of General-Purpose Graphics Processing Units (GPGPUs) is no exception, and many articles have been devoted to this problem.With this article, we provide a review of the techniques for implementing the SpMV kernel on GPGPUs that have appeared in the literature of the last few years. We discuss the issues and tradeoffs that have been encountered by the various researchers, and a list of solutions, organized in categories according to common features. We also provide a performance comparison across different GPGPU models and on a set of test matrices coming from various application domains.

Energy efficiency is a major concern in modern high-performance-computing. Still, few studies provide a deep insight into the power consumption of scientific applications. Especially for algorithms running on hybrid platforms equipped with hardware accelerators, like graphics processors, a detailed energy analysis is essential to identify the most costly parts, and to evaluate possible improvement strategies. In this paper we analyze the computational and power performance of iterative linear solvers applied to sparse systems arising in several scientific applications. We also study the gains yield by dynamic voltage/frequency scaling (DVFS), and illustrate that this technique alone cannot to reduce the energy cost to a considerable amount for iterative linear solvers. We then apply techniques that set the (multi-core processor in the) host system to a low-consuming state for the time that the GPU is executing. Our experiments conclusively reveal how the combination of these two techniques deliver a notable reduction of energy consumption without a noticeable impact on computational performance.

The numerical simulation of modern engineering problems via finite elements requires the solution of sparse linear systems of millions or even billions of unknowns. The algebraic multigrid (AMG) methods are the most common choice as linear solvers because of their fast convergence even for large-size problems. In this communication, we propose Chronos, a massively parallel implementation of a novel AMG framework, specifically designed to address complex problems by adapting its components, from the smoother, to the coarse grid correction and prolongation to the problem at hand. This work demonstrates not only the numerical performance of the proposed library, but also its robustness and adaptability to very challenging matrices, arising from different fields of application.

SummaryProper range and precision analysis play an important role in the development of fixed‐point algorithms for embedded system applications. Numerical linear algebra algorithms used to find singular value decomposition of symmetric matrices are suitable for signal and image‐processing applications. These algorithms have not been attempted much in fixed‐point arithmetic. The reason is wide dynamic range of data and vulnerability of the algorithms to round‐off errors. For any real‐time application, the range of the input matrix may change frequently. This poses difficulty for constant and variable fixed‐point formats to decide on integer wordlengths during float‐to‐fixed conversion process because these formats involve determination of integer wordlengths before the compilation of the program. Thus, these formats may not guarantee to avoid overflow for all ranges of input matrices. To circumvent this problem, a novel dynamic fixed‐point format has been proposed to compute integer wordlengths adaptively during runtime. Lanczos algorithm with partial orthogonalization, which is a tridiagonalization step in computation of singular value decomposition of symmetric matrices, has been taken up as a case study. The fixed‐point Lanczos algorithm is tested for matrices with different dimensions and condition numbers along with image covariance matrix. The accuracy of fixed‐point Lanczos algorithm in three different formats has been compared on the basis of signal‐to‐quantization‐noise‐ratio, number of accurate fractional bits, orthogonality and factorization errors. Results show that dynamic fixed‐point format either outperforms or performs on par with constant and variable formats. Determination of fractional wordlengths requires minimization of hardware cost subject to accuracy constraint. In this context, we propose an analytical framework for deriving mean‐square‐error or quantization noise power among Lanczos vectors, which can serve as an accuracy constraint for wordlength optimization. Error is found to propagate through different arithmetic operations and finally accumulate in the last Lanczos vector. It is observed that variable and dynamic fixed‐point formats produce vectors with lesser round‐off error than constant format. All the three fixed‐point formats of Lanczos algorithm have been synthesized on Virtex 7 field‐programmable gate array using Vivado high‐level synthesis design tool. A comparative study of resource usage and power consumption is carried out. Copyright © 2015 John Wiley & Sons, Ltd.

New features of our DSC system for distributing a symbolic computation task over a network of processors are described. A new scheduler sends parallel subtasks to those compute nodes that are best suited in handling the added load of CPU usage and memory. Furthermore, a subtask can communicate back to the process that spawned it by a co-routine style calling mechanism. Two large experiments are described in this improved setting. We have implemented an algorithm that can prove a number of more than 1,000 decimal digits prime in about 2 months elapsed time on some 20 computers. A parallel version of a sparse linear system solver is used to compute the solution of sparse linear systems over finite fields. We are able to find the solution of a 100,000 by 100,000 linear system with about 10.3 million non-zero entries over the Galois field with 2 elements using 3 computers in about 54 hours CPU time.

5 - Iterative schemes and convergence analysis