Gaussian Elimination Algorithm Research Articles

Gaussian Elimination: A case study in efficient genericity with MetaOCaml

The Gaussian Elimination algorithm is in fact an algorithm family — common implementations contain at least six (mostly independent) “design choices”. A generic implementation can easily be parametrized by all these design choices, but this usually leads to slow and bloated code. Using MetaOCaml’s staging facilities, we show how we can produce a natural and type-safe implementation of Gaussian Elimination which exposes its design choices at code-generation time, so that these choices can effectively be specialized away, and where the resulting code is quite efficient.

Anomalous dynamics in two- and three-dimensional Heisenberg-Mattis spin glasses

We investigate the spectral and localization properties of unmagnetized Heisenberg-Mattis spin glasses, in space dimensionalities $d=2$ and 3, at T=0. We use numerical transfer-matrix methods combined with finite-size scaling to calculate Lyapunov exponents, and eigenvalue-counting theorems, coupled with Gaussian elimination algorithms, to evaluate densities of states. In $d=2$ we find that all states are localized, with the localization length diverging as $\omega^{-1}$, as energy $\omega \to 0$. Logarithmic corrections to density of states behave in accordance with theoretical predictions. In $d=3$ the density-of-states dependence on energy is the same as for spin waves in pure antiferromagnets, again in agreement with theoretical predictions, though the corresponding amplitudes differ.

A class of novel parallel algorithms for the solution of tridiagonal systems

In this paper, a new class of parallel Gaussian elimination algorithms is presented for the solution of tridiagonal matrix systems. The new algorithms, called ACER (alternating cyclic elimination and reduction), combine the advantages of the well known cyclic elimination algorithm (which is fast) and the cyclic reduction algorithms (which requires fewer operations). The ACER algorithms are developed with the unifying graph model.

An optimal schedule for Gaussian elimination on an MIMD architecture

In this paper, we propose an algorithm for allocating the tasks of the well known Gaussian Elimination Algorithm on an MIMD architecture and prove that the schedule is optimal in order of magnitude, up to a polylog factor.

Multi-Variable Non-Singular BEM for 2-D Potential Problems

A multi-variable non-singular boundary element method (MNBEM) is presented for 2-D potential problems. This method is based on the coincident collocation of non-singular boundary integral equations (BIEs) of the potential and its derivatives, where the nodal potential derivatives are considered independent of the nodal potential and flux. The system equation is solved to determine the unknown boundary potentials and fluxes, with high accuracy boundary nodal potential derivatives obtained from the solution at the same time. A modified Gaussian elimination algorithm was developed to improve the solution efficiency of the final system equation. Numerical examples verify the validity of the proposed algorithm.

Weak theories of linear algebra

We investigate the theories Open image in new window of linear algebra, which were originally defined to study the question of whether commutativity of matrix inverses has polysize Frege proofs. We give sentences separating quantified versions of these theories, and define a fragment Open image in new window in which we can interpret a weak theory V1 of bounded arithmetic and carry out polynomial time reasoning about matrices - for example, we can formalize the Gaussian elimination algorithm. We show that, even if we restrict our language, Open image in new window proves the commutativity of inverses.

Graph Transformation and Designing Parallel Sparse Matrix Algorithms beyond Data Dependence Analysis

Algorithms are often parallelized based on data dependence analysis manually or by means of parallel compilers. Some vector/matrix computations such as the matrix‐vector products with simple data dependence structures (data parallelism) can be easily parallelized. For problems with more complicated data dependence structures, parallelization is less straightforward. The data dependence graph is a powerful means for designing and analyzing parallel algorithms. However, for sparse matrix computations, parallelization based on solely exploiting the existing parallelism in an algorithm does not always give satisfactory results. For example, the conventional Gaussian elimination algorithm for the solution of a tri‐diagonal system is inherently sequential, so algorithms specially for parallel computation has to be designed. After briefly reviewing different parallelization approaches, a powerful graph formalism for designing parallel algorithms is introduced. This formalism will be discussed using a tri‐diagonal system as an example. Its application to general matrix computations is also discussed. Its power in designing parallel algorithms beyond the ability of data dependence analysis is shown by means of a new algorithm called ACER (Alternating Cyclic Elimination and Reduction algorithm).

Cost-optimal Gaussian Elimination Algorithm on RCC-FULL

This paper presents a cost-optimal parallel algorithm to find the solution of n simultaneous linear equations on RCC-FULL λ(N, L) network. The algorithm is designed following the strategy of the sequential Gaussian Elimination algorithm. The algorithm utilizes strength of the RCC-FULL network which helps to reduce the time complexity and has the time complexity of O(n) with O(n 2) processors.

Impact of the implementation of MPI point-to-point communications on the performance of two general sparse solvers

We examine the send and receive mechanisms of MPI and show how to implement message passing robustly so that performance is not significantly affected by changes to the MPI system. We discuss this within the context of two different parallel algorithms for sparse Gaussian elimination: a multifrontal solver ( MUMPS), and a supernodal one ( SuperLU). The performance of our initial strategies based on simple MPI point-to-point communication primitives is very sensitive to the MPI system, particularly the way MPI buffers are used. Using nonblocking communication primitives improves the performance and robustness, but at the cost of increased code complexity.

A faster algorithm for solving linear algebraic equations on the star graph

The problem of solving a linear system of equations is widely encountered in many fields of science and engineering. In this paper, we present a parallel algorithm to solve the above problem on a star graph. The proposed solution (i) is based on a variant of the Gaussian elimination algorithm (GE) called the successive Gaussian elimination algorithm (SGE) (IEE Proc. Comput. Digit. Tech. 143 (4) (1996)) and (ii) supports partial pivoting to provide numerical stability. We present efficient matrix distribution techniques on the star graph. Our proposed parallel algorithm employs these techniques to reduce communication overhead during matrix operations on the star graph. We estimate the performance of our parallel algorithm and demonstrate its effectiveness by comparing it with a recent algorithm for the same problem on star graphs (IEEE Trans. Parallel Distrib. Systems 8 (8) (1997) 803).

Efficient and Stable Solution of M-Matrix Linear Systems of (Block) Hessenberg Form

The problem of solving large linear systems whose coefficient matrix is a sparse M-matrix in block Hessenberg form has recently received much attention, especially for applications in Markov chains and queueing theory. Stewart proposed a recursive algorithm which is shown to be backward stable. Although the theoretical derivation of such an algorithm is very simple, its efficient implementation is logically rather involved. An analysis of its computational cost in the case where the initial coefficient matrix satisfies quite general sparsity properties can be found in [P. Favati et al., Acta Tecnica Acad. Sci. Hungar., 108 (1997--1999), pp. 89--105]. In this paper we devise a different divide-and-conquer strategy for the solution of block Hessenberg linear systems. Our approach follows from a recursive application of the block Gaussian elimination algorithm. For dense matrices, the present method has a computational cost comparable to that of Stewart's algorithm; for sparse matrices it is more efficient and more robust.

Improved Symbolic and Numerical Factorization Algorithms for Unsymmetric Sparse Matrices

We present algorithms for the symbolic and numerical factorization phases in the direct solution of sparse unsymmetric systems of linear equations. We have modified a classical symbolic factorization algorithm for unsymmetric matrices to inexpensively compute minimal elimination structures. We give an efficient algorithm to compute a near-minimal data-dependency graph for unsymmetric multifrontal factorization that is valid irrespective of the amount of dynamic pivoting performed during factorization. Finally, we describe an unsymmetric-pattern multifrontal algorithm for Gaussian elimination with partial pivoting that uses the task- and data-dependency graphs computed during the symbolic phase. These algorithms have been implemented in WSMP---an industrial strength sparse solver package---and have enabled WSMP to significantly outperform other similar solvers. We present experimental results to demonstrate the merits of the new algorithms.

A Generalized Sylvester Identity and Fraction-free Random Gaussian Elimination

Sylvester’s identity is a well-known identity that can be used to prove that certain Gaussian elimination algorithms are fraction free. In this paper we will generalize Sylvester’s identity and use it to prove that certain random Gaussian elimination algorithms are fraction free. This can be used to yield fraction free algorithms for solving Ax=b(x≥ 0) and for the simplex method in linear programming.

A multi-region boundary element method for multigroup neutron diffusion calculations

For the analysis of a two-dimensional nuclear system consisting of a number of homogeneous regions (termed cells), first the cell matrices which depend solely on the material composition and geometrical dimension of the cell (hence on the cell type) are constructed using a boundary element formulation based on the multigroup boundary integral equation. For a particular nuclear system, the cell matrices are utilized in the assembly of the global system matrix in block-banded form using the newly introduced concept of virtual side. For criticality calculations, the classical fission source iteration is employed and linear system solutions are by the block Gaussian-elimination algorithm. The numerical applications show the validity of the proposed formulation both through comparison with analytical solutions and assessment of benchmark problem results against alternative methods.

PROCESSOR-TIME-OPTIMAL SYSTOLIC ARRAYS

Minimizing the amount of time and number of processors needed to perform an application reduces the application's fabrication cost and operation costs. A directed acyclic graph (dag) model of algorithms is used to define a time-minimal schedule and a processor-time-minimal schedule, We present a technique for finding a lower bound on the number of processors needed to achieve a given schedule of an algorithm. The application of this technique is illustrated with a tensor product computation. We then apply the technique to the free schedule of algorithms for matrix product, Gaussian elimination, and transitive closure. For each, we provide a time-minimal processor schedule that meets these processor lower bounds, including the one for tensor product.

Modelling and analysis of communication overhead for parallel matrix algorithms

We develop and analyze novel algorithms that make efficient use of the communication system in distributed memory architectures with processing elements interconnected by a hypercube network. The algorithms studied here include the parallel Gauss-Jordan (GJ) matrix inversion algorithm and the Gaussian Elimination (GE) algorithm for LU factorization. We first propose a new broadcasting algorithm on the hypercube multiprocessor for the parallel GJ algorithm. This algorithm ensures that the data items are sent out from the source and arrive at the destinations at the earliest possible time. We then present a parallel GJ inversion algorithm using row partitioning. This algorithm exploits a compute-and-send-ahead strategy for achieving overlapping of communication and computation, and the resulting framework leads to rigorous analytical and model-based numerical performance analysis of our parallel algorithms. In particular, we prove a lower bound on the matrix size such that data transmission is fully overlapped by computation. We also prove that the message queue length in the input buffer of each processor is at most two. We next consider the GJ algorithm under submatrix partitioning, with or without pivoting. We show that when submatrix partitioning is used, even when communication is fully overlapped by computation, the communication overhead is larger than when using row partitioning. Thus, we show that by minimizing the communication overhead, the row partitioning scheme can indeed have better overall performance than the submatrix partitioning scheme. Finally, we extend the idea of overlapping communication and computation to the parallel LU factorization algorithm.

UNIFIED PARALLEL ALGORITHMS FOR GAUSSIAN ELIMINATION WITH BACKWARD SUBSTITUTION ON PRODUCT NETWORKS

The increasing interest in product networks (PNs) as a method of combining desirable properties of component networks, has prompted a need for the general study of the algorithmic issues related to this important class of interconnection networks. In this paper we present unified parallel algorithms for Gaussian elimination, with partial and complete pivoting, on product networks. A parallel algorithm for backward substitution is also presented. The proposed algorithms are network independent and are also independent of the matrix distribution methods employed. These algorithms can be used on a wide range of PNs including hypercube, mesh, and k-ary n-cube. Unified models for estimating computation time and interprocessor communication time are also presented. These models are then used to measure the performance of the proposed algorithms on several product networks

A systolic array for the parallel solution of block trid1agonal linear systems (Sats)

In this paper the parallel solution of block tridiagonal linear systems derived from the finite difference/element discretisation of 2D/3D elliptic partial differential equations which occur frequently in scientific and engineering applications is presented. The direct solution by the Gaussian elimination algorithm is used and a parallel architecture designed employing a systolic array (SATS) which comprises of triangular and rectangular systolic arrays to perform the block triangularisation of the system. The block triangular system is then solved by a pipeline using local memory systems. The results show that when the system order is > 10 then the speedup of the parallel system is significant

Scheduling algorithms for parallel Gaussian elimination with communication costs

We consider a graph theoretical model and study a parallel implementation of the well-known Gaussian elimination method on parallel distributed memory architectures, where the communication delay for the transmission of an elementary data is higher than the computation time of an elementary instruction. We propose and analyze two low-complexity algorithms for scheduling the tasks of the parallel Gaussian elimination on an unbounded number of completely connected processors. We compare these two algorithms with a higher-complexity general-purpose scheduling algorithm, the DSC heuristic, proposed by A. Gerasoulis and T. Yang (1993).

Symmetric Gaussian elimination for cauchy-type matrices with application to positive definite Toeplitz matrices

The problem of solving linear equations with a Toeplitz matrix \(T\) appears in many applications. Often \(T\) is positive definite but ill-conditioned with many small eigenvalues. In this case fast and superfast algorithms may show a very poor behavior or even break down. In recent papers the transformation of a Toeplitz matrix into a Cauchy-type matrix is proposed. The resulting new linear equations can be solved in \(O(n^2)\) operations using standard pivoting strategies which leads to very stable fast methods also for ill-conditioned systems. The basic tool is the formulation of Gaussian elimination for matrices with low displacement rank. In this paper, we will transform a Hermitian Toeplitz matrix \(T\) into a Cauchy-type matrix \(B=FTF^{\rm H}\) by applying the Fourier transform. We will prove some useful properties of \(B\) and formulate a symmetric Gaussian elimination algorithm for positive definite \(B\). Using the symmetry and persymmetry of \(T\) we can reduce the total costs of this algorithm compared with unsymmetric Gaussian elimination. For complex Hermitian \(T\), the complexity of the new algorithm is then nearly the same as for the Schur algorithm. Furthermore, it is possible to include some strategies for ill-conditioned positive definite matrices that are well-known in optimization. Numerical examples show that this new algorithm is fast and reliable.