LU decomposition

This paper addresses the problem of algorithm discov- ery, via evolutionary search, in the context of matrix multiplication. The traditional multiplication algorithm requires O(n 3) multiplications for square matrices of or- der n. Strassen (Strassen 1969) discovered a re, cursive matrix multiplication algorithm requiring only seven multiplications at each level, resulting in a runtime of O(n/x7), or O(n2"81). We have

Strassen's algorithm for fast matrix-matrix multiplication has been implemented for matrices of arbitrary shapes on the CRAY-2 and CRAY Y-MP supercomputers. Several techniques have been usCd to reduce the scratch space requirement for this algorithm while simultaneously preserving a high level of performance. When the resulting Strassen-based matrix multiply routine is combined with some routines from the new LAPACK library, LV decomposition can be performed with rates significantly higher than those achieved by conventional means. We succeeded in factoring a 2048 X 2048 matrix on the CRAY Y-MP at a rate equivalent to 325 MFLOPS.

LU decomposition is a fundamental in linear algebra. Numerous tools exists that provide this important factorization. The authors present the conditions for a matrix to have none, one, or infinitely many LU factorizations. In the case where no factorization exists, the authors illustrate how to approximate an LU decomposition by considering LU factorization of nearby matrices.

LU decomposition is a fundamental in linear algebra. Numerous tools exists that provide this important factorization. The authors present the conditions for a matrix to have none, one, or infinitely many LU factorizations. In the case where no factorization exists, the authors illustrate how to approximate an LU decomposition by considering LU factorization of nearby matrices.

An efficient implicit lower-upper symmetric Gauss-Seidel (LU-SGS) solution algorithm has been developed for a high order multi-domain spectral difference method on unstructured hexahedral grids. The LU-SGS solver is preconditioned by the block element matrix, and the system of equations is then solved with an exact LU decomposition approach. In addition, the effects of several parameters on the convergence rate in the implicit scheme have been investigated for both external and internal flows. The implicit scheme has shown a speed up factor of more than an order of magnitude relative to the multi-stage explicit Runge-Kutta scheme for several demonstration problems. he advantage of high-order methods (order of accuracy > 2) over first and second-order ones is well known in the CFD community. Generally speaking, with the same number of degrees-of-freedom (DOFs) or solution unknowns, high-order methods are capable of producing much more accurate results. For problems requiring very high accuracy, e.g., wave propagation problems in computational aeroacoustics, high-order methods have been the main choice. Many high-order methods were developed for structured grids, e.g., ENO/WENO methods 1, compact methods 2-3 , optimized methods 4 , to name just a few. In the last two decades, there have been intensive research efforts on high-order methods for unstructured grids since many real world applications have complex geometries. An incomplete list of notable examples includes the spectral element method 5 , multi-domain spectral method 6-7 , kexact finite volume method 8 , WENO methods 9 , discontinuous Galerkin (DG) method 10-12 , high-order residual distribution methods 13 , spectral volume (SV) 14-17 and spectral difference (SD) methods 18-21 . Among those methods, some are based on the weighted residual form of the governing equations, for instance the DG method 10-12 . Some are based on the integral form of the governing equations, e. g., the k-exact finite volume method 8 and SV methods 14-17 . Others, such as the staggered grid multi-domain spectral method 6-7 and the SD method 18-21 , are based on the differential form. When one chooses a particular method for three-dimensional applications, the cost and the complexity in implementing the method is often an important factor. It is obvious that methods based on the differential form are the easiest to implement since they do not involve surface or volume integrals. This is particularly true when highorder curved boundaries need to be dealt with. We recently developed a high order SD method 22 for the three dimensional Navier-Stokes equations on unstructured hexahedral grids. High order accuracy and spectral convergence are achieved for several benchmark problems. It was also shown that the wall boundaries must be approximated with high-order surfaces. An explicit Runge-Kutta time integration scheme was used in the implementation. Although the explicit scheme is easy to implement and has high-order accuracy in time, it suffered from slow convergence, especially for viscous grids which are clustered in the viscous boundary layer. It is wellknown that high-order methods are restricted to a smaller CFL number than low order ones. In addition, they also

This research introduces a row compression and nested product decomposition of an n × n hierarchical representation of a rank structured matrix A, which extends the compression and nested product decomposition of a quasiseparable matrix. The hierarchical parameter extraction algorithm of a quasiseparable matrix is efficient, requiring only O(nlog(n)) operations, and is proven backward stable. The row compression is comprised of a sequence of small Householder transformations that are formed from the low-rank, lower triangular, off-diagonal blocks of the hierarchical representation. The row compression forms a factorization of matrix A, where A = QC, Q is the product of the Householder transformations, and C preserves the low-rank structure in both the lower and upper triangular parts of matrix A. The nested product decomposition is accomplished by applying a sequence of orthogonal transformations to the low-rank, upper triangular, off-diagonal blocks of the compressed matrix C. Both...

Dedicated to my loving parents, John and Anne Hudachek, who are watching from above; and to my very patient husband, Dan Buswell, who is ecstatic to have his wife back. v ACKNOWLEDGEMENTS The monumental task of a doctoral degree can only be undertaken with the consistent support and encouragement of many individuals. First and foremost, I must thank Professor Michael Stewart for directing my research all these years. His instruction, patience, and support have made me the researcher I am today. Working with Professor Stewart has been a great honor, and I look forward to our continued collaboration. I wish to thank the Computer Science faculty at Georgia State University for their leadership in guiding me towards a doctoral degree in computing. Special thanks goes to my co-chair, Professor Belkasim for introducing me to image processing and encouraging me to research in this area. I express gratitude to my committee members, Professors Raj Sunderraman, Yi Pan and Jon Preston, for their valuable and constructive feedback on my research. Many thanks go to my mentor, Professor Catherine Aust, who is responsible for starting me down this doctoral path in computer science, and assisted me in navigating the process. I am grateful to my dear colleague, Professor Catherine Matos, for her steadfast support, and boundless efforts to discuss my research at all hours of the day. I also extend my thanks to Professor Rebecca Rizzo for graciously opening her office to me this last year. I would like to acknowledge Mr. Trey Olmsted for his technical skills in maintaining my research computer, and providing me with the freedom to access my research anywhere at anytime. My sincere thanks also go to Ms. Tina Breckenridge for reminding me to laugh during the last decade, and cheering me onward. Lastly, I thank my entire family. In the end, all you have is family; it is only through their love, care, and support that I have been able to finish this research marathon. Especially, my husband Dan for travailing this decade long journey with me, and to my brother Michael for his unwavering moral support and financial assistance. vi TABLE OF CONTENTS ACKNOWLEDGEMENTS .

In this paper, we postulate a new decomposition theorem of a matrix A into two matrices, namely, a lower triangular matrix M, in which all entries are determinants, and an upper triangular matrix U whose entries are also in determinant form. From a well-known theorem on the pivot elements of the Doolittle-Gauss elimination process, we deduce a corollary to obtain a diagonal matrix D. With it, we scale the elementary lower triangular matrix of the Doolittle- Gauss elimination process and deduce a new elementary lower triangular matrix. Applying this linear transformation to A by means of both minimum and complete pivoting strategies, we obtain the determinant of A as if it had been calculated by means of a Laplace expansion. If we apply this new linear transformation and the above pivot strategy to an augmented matrix (A|b), we obtain a Cramer's solution of the linear

Abstract. For saddle point problems in fluid dynamics, several popular preconditioners exploit the block struc-ture of the problem to construct block triangular preconditioners. The performance of such preconditioners depends on whether fast, approximate solvers for the linear ...

In this work, the solution of a large sparse linear system of equations with an arbitrary sparsity pattern is obtained by using LU-decomposition method as well as numerical structure approach. The LU-decomposition method is based on Doolittle's method while the numerical structure approach is based on Cramer's rule. The numerical structure approach produces direct solution without facing fill-in problems as encountered in LU-decomposition. In order to reduce the 'fill-ins' in the decomposition, the powers of a Boolean matrix, obtained from the coefficient matrix A are taken so that the 'fill-ins' in the structure of A can be known in advance. The position of fill-ins in A are thus determined in the best choice manner, that is, it is very effective and memory-wise cheap. We also outline a method by using numerical structure with reduced computation efforts. Finally, experiments are performed on eight examples to compare the efficiency of the proposed methods. The results obtained are reported in a table. It is found that the LU-decomposition method is much better than numerical structure. The usefulness of numerical structure approach is also discussed. (~) 2001 Elsevier Science Ltd. All rights reserved.

The robustness of the recently introduced circulant blockfactorization (CBF) preconditioners is studied in the case of finite element matrices arising from the discretization of the 2D Navier equations of elasticity. Conforming triangle finite elements are used for the numerical solution of the differential problem. The proposed preconditioner Me is constructed by CBF approximation of the block-diagonal part of the stiffness matrix. In other words, we implement in our algorithm the circulant block-factorization into the framework of the displacement decomposition technique. The estimate ~(M~IK)= 0 (~ is proved asymptotically on N, where N is the size of the discrete problem. Note, that the corresponding known estimate for the widely used incomplete factorization displacement decomposition preconditioner MILV is a(M~'~uK) = The theoretical estimate as well as the presented numerical tests show some significant advantages of this new approach for a PCG iterative solution of almost incompressible elastic problems, that is when the modified Poisson ratio v tends to the incompressible limit case v = 1.

Decomposition Approach for Inverse Matrix Calculation 113 of a corresponding inverse of a modified matrix A *-1 are derived in (Strassen, 1969). The components of the inverse matrix can be evaluated analytically. Finding the inverse matrix is related with a lot of calculations. Instead of direct finding an inverse matrix, it is worth to find analytical relations where lower dimensions inverse matrices components are available. Here analytical relations for inverse matrix calculation are derived and the corresponding MATLAB codes are illustrated.

The comprehensive LU decomposition of a parametric matrix consists of a case analysis of the LU factors for each specialization of the parameters. Special cases can be discontinuous with respect to the parameters, the discontinuities being triggered by zero pivots encountered during factorization. For polynomial matrices, we describe an implementation of comprehensive LU decomposition in Maple, using the RegularChains package.

The ScaLAPACK library for parallel dense matrix computations is built on top of the BLACS communications layer. In this work, we investigate the use of BSPlib as the basis for a communications layer. We examine the LU decomposition from ScaLAPACK and develop a BSP version which is signi cantly faster. The savings in communication time are typically 10-15%. The gain in overall execution time is less pronounced, but still signi cant. We present the main features of a new library, BSP2D, which w e propose to develop for porting the whole of ScaLAPACK.

The Cholesky decomposition plays an important role in finding the inverse of the correlation matrices. As it is a fast and numerically stable for linear system solving, inversion, and factorization compared to singular valued decomposition (SVD), QR factorization and LU decomposition. As different methods exist to find the Cholesky decomposition of a given matrix, this paper presents the comparative study of a proposed RChol algorithm with the conventional methods. The RChol algorithm is an explicit way to estimate the modified Cholesky factors of a dynamic correlation matrix. Cholesky decomposition is a fast and numerically stable for linear system solving, inversion, and factorization compared to singular valued decomposition (SVD), QR factorization and LU decomposition [1]. The wireless communication system is highly dependent on matrix inversion of the correlation matrix. Such system consists of a huge matrix inversion. An outdoor wireless communication has a time-varying channel which changes dynamically for mobile user. In case of narrowband channel, the channel is considered constant for a symbol duration, whereas for broadband, it is changing within a symbol period. Such time-varying channel forms the special structure of channel matrix and correlation matrix. To exploit such special structure, a novel modified recessive Cholesky (RChol) algorithm is introduced in [2]. Our proposed (RChol) algorithm is a computational efficient algorithm to compute the modified Cholesky factors of known as well as an unknown covariance matrix. In this paper, we present the comparative study of conventional Cholesky algorithm and the RChol algorithm to manifest the importance of the proposed algorithm in a highly dynamic wireless communication.

We investigate the performance and scalability of the randomized CX low-rank matrix factorization and demonstrate its applicability through the analysis of a 1TB mass spectrometry imaging (MSI) dataset, using Apache Spark on an Amazon EC2 cluster, a Cray XC40 system, and an experimental Cray cluster. We implemented this factorization both as a parallelized C implementation with hand-tuned optimizations and in Scala using the Apache Spark high-level cluster computing framework. We obtained consistent performance across the three platforms: using Spark we were able to process the 1TB size dataset in under 30 minutes with 960 cores on all systems, with the fastest times obtained on the experimental Cray cluster. In comparison, the C implementation processed the 1TB size dataset 21X faster on the Amazon EC2 system, due to careful cache optimizations, bandwidth-friendly access of matrices and vector computation using SIMD units. We report these results and their implications on the hardw...

Hereby, I, Ta Minh THANH, consciously assure that for the manuscript "Consideration of A Robust Watermarking Algorithm for Color Image Using Improved QR Decomposition" the following is fulfilled: 1) This material is the authors' own original work, which has not been previously published elsewhere. 2) The paper is not currently being considered for publication elsewhere. 3) The paper reflects the authors' own research and analysis in a truthful and complete manner. 4) The paper properly credits the meaningful contributions of co-authors and co-researchers. 5) The results are appropriately placed in the context of prior and existing research. 6) All sources used are properly disclosed (correct citation). Literally copying of text must be indicated as such by using quotation marks and giving proper reference. 7) All authors have been personally and actively involved in substantial work leading to the paper, and will take public responsibility for its content. I agree with the above statements and declare that this submission follows the policies of Solid State Ionics as outlined in the Guide for Authors and in the Ethical Statement.

An efficient technique for partitioning and programming linear algebra algorithms on concurrent architectures is described and applied to 2-D wavefront arrays. The mapping of the computational elements (processes) to processors is based on the concept of folding. The mapping pattern on the 2-D full-size mesh of processes is a composition of symmetric tiles of size 2J(N) x 2,/(N), N being the number of processors. The algorithm can be partitioned according to a globally sequential, locally parallel scheme. The code optimisation is performed by programming a few different types of tile, according to the algorithm. When the size of the problem is much larger than the size of the mesh of processors, a linear speed-up is achieved independently of the number of processors. Experimental results are presented for matrix multiplication, L U decomposition and the solution of triangular system equations on 2-D meshes of transputers programmed in Occam.

Recursive block decomposition algorithms (also known as quadtree algorithms when the blocks are all square) have been proposed to solve well-known problems such as matrix addition, multiplication, inversion, determinant computation, block LDU decomposition and Cholesky and QR factorization. Until now, such algorithms have been seen as impractical, since they require leading submatrices of the input matrix to be invertible (which is rarely guaranteed). We show how to randomize an input matrix to guarantee that submatrices meet these requirements, and to make recursive block decomposition methods practical on well-conditioned input matrices. The resulting algorithms are elegant, and we show the recursive programs can perform well for both dense and sparse matrices, although with randomization dense computations seem most practical. By ‘homogenizing’ the input, randomization provides a way to avoid degeneracy in numerical problems that permits simple recursive quadtree algorithms to so...

This paper presents implementations of new methods solving LU factorization used in engineering applications. The implementations are done on the Alliant FX/80 minisupercomputer and use Level 3 Basic Linear Algebra Subprograms. Three ways of expressing the LU factorization in terms of blocked algorithms are considered. The performance of the blocked algorithms, using the parallel vector facilities, are compared to a noblock algorithm using only subprograms of leve l 1 a n d 2 B L A S .

The Gauss-Huard algorithm (the GHA) is a specialized version of Gauss-Jordan elimination for the solution of linear systems that, enhanced with column pivoting, exhibits numerical stability and computational cost close to those of the conventional solver based on the LU factorization with row pivoting. Furthermore, the GHA can be formulated as a procedure rich in matrix multiplications, so that high performance can be expected on current architectures with multi-layered memories. Unfortunately, in principle the GHA does not admit the introduction of look-ahead, a technique that has been demonstrated to be rather useful to improve the performance of the LU factorization on multi-threaded platforms with high levels of hardware concurrency. In this paper we analyze the effect of this drawback on the implementation of the GHA on systems accelerated with graphics processing units (GPUs), exposing the roles of the CPU-to-GPU and single precision-to-double precision performance ratios, as well as the contribution from the operations in the algorithm's critical path. Keywords: Linear systems of equations • Gauss-Huard algorithm • LU factorization • Multicore processors • Graphics processing units (GPUs) • Mixed precision • High performance All researchers acknowledge the support from the EHFARS project funded by the German Ministry of Education and Research BMBF. E.S. Quintana-Ortí was supported by the CICYT project TIN2014-53495-R of the Ministerio de Economía y Competitividad and FEDER.

The concept of block{cyclic order elimination can be applied to out{of{ core LU and QR matrix factorizations on distributed memory architectures equipped with a parallel I/O system. This elimination scheme provides load balanced computation in both the factor and solve phases and further optimizes the use of the network bandwidth to perform I/O operations. Stability of LU factorization is enforced by full column pivoting. Performance results are presented for the Connection Machine system CM{5.

The Gauss-Huard algorithm (the GHA) is a specialized version of Gauss-Jordan elimination for the solution of linear systems that, enhanced with column pivoting, exhibits numerical stability and computational cost close to those of the conventional solver based on the LU factorization with row pivoting. Furthermore, the GHA can be formulated as a procedure rich in matrix multiplications, so that high performance can be expected on current architectures with multi-layered memories. Unfortunately, in principle the GHA does not admit the introduction of look-ahead, a technique that has been demonstrated to be rather useful to improve the performance of the LU factorization on multi-threaded platforms with high levels of hardware concurrency. In this paper we analyze the effect of this drawback on the implementation of the GHA on systems accelerated with graphics processing units (GPUs), exposing the roles of the CPU-to-GPU and single precision-to-double precision performance ratios, as well as the contribution from the operations in the algorithm's critical path. Keywords: Linear systems of equations • Gauss-Huard algorithm • LU factorization • Multicore processors • Graphics processing units (GPUs) • Mixed precision • High performance All researchers acknowledge the support from the EHFARS project funded by the German Ministry of Education and Research BMBF. E.S. Quintana-Ortí was supported by the CICYT project TIN2014-53495-R of the Ministerio de Economía y Competitividad and FEDER.

This paper presents a new approach for the solution of Linear Programming Problems with the help of LU Factorization Method of matrices. This method is based on the fact that a square matrix can be factorized into the product of unit lower triangular matrix and upper triangular matrix. In this method, we get direct solution without iteration. We also show that this method is better than simplex method.

Determining I/O lower bounds is a crucial step in obtaining communication-efficient parallel algorithms, both across the memory hierarchy and between processors. Current approaches either study specific algorithms individually, disallow programmatic motifs such as recomputation, or produce asymptotic bounds that exclude important constants. We propose a novel approach for obtaining precise I/O lower bounds on a general class of programs, which we call Simple Overlap Access Programs (SOAP). SOAP analysis covers a wide variety of algorithms, from ubiquitous computational kernels to full scientific computing applications. Using the red-blue pebble game and combinatorial methods, we are able to bound the I/O of the SOAP-induced Computational Directed Acyclic Graph (CDAG), taking into account multiple statements, input/output reuse, and optimal tiling. To deal with programs that are outside of our representation (e.g., non-injective access functions), we describe methods to approximate them with SOAP. To demonstrate our method, we analyze 38 different applications, including kernels from the Polybench benchmark suite, deep learning operators, and-for the first time-applications in unstructured physics simulations, numerical weather prediction stencil compositions, and full deep neural networks. We derive tight I/O bounds for several linear algebra kernels, such as Cholesky decomposition, improving the existing reported bounds by a factor of two. For stencil applications, we improve the existing bounds by a factor of up to 14. We implement our method as an open-source tool, which can derive lower bounds directly from provided C code.

Determining I/O lower bounds is a crucial step in obtaining communication-efficient parallel algorithms, both across the memory hierarchy and between processors. Current approaches either study specific algorithms individually, disallow programmatic motifs such as recomputation, or produce asymptotic bounds that exclude important constants. We propose a novel approach for obtaining precise I/O lower bounds on a general class of programs, which we call Simple Overlap Access Programs (SOAP). SOAP analysis covers a wide variety of algorithms, from ubiquitous computational kernels to full scientific computing applications. Using the red-blue pebble game and combinatorial methods, we are able to bound the I/O of the SOAP-induced Computational Directed Acyclic Graph (CDAG), taking into account multiple statements, input/output reuse, and optimal tiling. To deal with programs that are outside of our representation (e.g., non-injective access functions), we describe methods to approximate them with SOAP. To demonstrate our method, we analyze 38 different applications, including kernels from the Polybench benchmark suite, deep learning operators, and-for the first time-applications in unstructured physics simulations, numerical weather prediction stencil compositions, and full deep neural networks. We derive tight I/O bounds for several linear algebra kernels, such as Cholesky decomposition, improving the existing reported bounds by a factor of two. For stencil applications, we improve the existing bounds by a factor of up to 14. We implement our method as an open-source tool, which can derive lower bounds directly from provided C code.

Matrix factorizations are among the most important building blocks of scientific computing. However, state-of-the-art libraries are not communication-optimal, underutilizing current parallel architectures. We present novel algorithms for Cholesky and LU factorizations that utilize an asymptotically communication-optimal 2.5D decomposition. We first establish a theoretical framework for deriving parallel I/O lower bounds for linear algebra kernels, and then utilize its insights to derive Cholesky and LU schedules, both communicating 3 /(√) elements per processor, where M is the local memory size. The empirical results match our theoretical analysis: our implementations communicate significantly less than Intel MKL, SLATE, and the asymptotically communication-optimal CANDMC and CAPITAL libraries. Our code outperforms these state-of-the-art libraries in almost all tested scenarios, with matrix sizes ranging from 2,048 to 524,288 on up to 512 CPU nodes of the Piz Daint supercomputer, decreasing the time-to-solution by up to three times. Our code is ScaLAPACK-compatible and available as an open-source library.

A singular matrix A may have more than one LU factorizations. In this work the set of all LU factorizations of A is explicitly described when the lower triangular matrix L is nonsingular. To this purpose, a canonical form of A under left multiplication by unit lower triangular matrices is introduced. This canonical form allows us to characterize the matrices that have an LU factorization and to parametrize all possible LU factorizations. Formulae in terms of quotient of minors of A are presented for the entries of this canonical form.

SPEX Left LU is a software package for exactly solving unsymmetric sparse linear systems. As a component of the sparse exact (SPEX) software package, SPEX Left LU can be applied to any input matrix, A, whose entries are integral, rational, or decimal, and provides a solution to the system Ax = b which is either exact or accurate to user-speciied precision. SPEX Left LU preorders the matrix A with a user-speciied ill-reducing ordering and computes a left-looking LU factorization with the special property that each operation used to compute the L and U matrices is integral. Notable additional applications of this package include benchmarking the stability and accuracy of state-of-the-art linear solvers, and determining whether singular-to-double-precision matrices are indeed singular. Computationally, this paper evaluates the impact of several novel pivoting schemes in exact arithmetic, benchmarks the exact iterative solvers within Linbox, and benchmarks the accuracy of MATLAB sparse backslash. Most importantly, it is shown that SPEX Left LU outperforms the exact iterative solvers in run time on easy instances and in stability as the iterative solver fails on a sizeable subset of the tested (both easy and hard) instances. The SPEX Left LU package is written in ANSI C, comes with a MATLAB interface, and is distributed via GitHub, as a component of the SPEX software package, and as a component of SuiteSparse.

Gaussian Elimination is commonly used to solve dense linear systems in scientific models. In a large number of applications, a need arises to solve many small size problems, instead of few large linear systems. The size of each of these small linear systems depends, for example, on the number of the ordinary differential equations (ODEs) used in the model, and can be on the order of hundreds of unknowns. To efficiently exploit the computing power of modern accelerator hardware, these linear systems are processed in batches. To improve the numerical stability of the Gaussian Elimination, at least partial pivoting is required, most often accomplished with row pivoting. However, row pivoting can result in a severe performance penalty on GPUs because it brings in thread divergence and non-coalesced memory accesses. The state-of-the-art libraries for linear algebra that target GPUs, such as MAGMA, focus on large matrix sizes. They change the data layout by transposing the matrix to avoid these divergence and non-coalescing penalties. However, the data movement associated with transposition is very expensive for small matrices. In this paper, we propose a batched LU factorization for GPUs by using a multi-level blocked right looking algorithm that preserves the data layout but minimizes the penalty of partial pivoting. Our batched LU achieves up to 2.5-fold speedup when compared to the alternative CUBLAS solutions on a K40c GPU and 3.6-fold speedup over MKL on a node of the Titan supercomputer at ORNL in a nuclear reaction network simulation.

In this paper, we introduce novel fast matrix inversion algorithms that leverage triangular decomposition and recurrent formalism, incorporating Strassen's fast matrix multiplication. Our research places particular emphasis on triangular matrices, where we propose a novel computational approach based on combinatorial techniques for finding the inverse of a general non-singular triangular matrix. Unlike iterative methods, our combinatorial approach for (block) triangular-type matrices enables direct computation of the matrix inverse through a nonlinear combination of carefully selected combinatorial entries from the initial matrix. This unique characteristic makes our proposed method fully parallelizable, offering significant potential for efficient implementation on parallel computing architectures. While it is widely acknowledged that combinatorial algorithms typically suffer from exponential time complexity, thus limiting their practicality, our approach demonstrates intriguing features that allow the derivation of recurrent relations for constructing the matrix inverse. By combining the (block) combinatorial approach, with a recursive triangular split method for inverting triangular matrices, we develop potentially competitive algorithms that strike a balance between efficiency and accuracy. To establish the validity and effectiveness of our approach, we provide rigorous mathematical proofs of the newly presented method. Additionally, we conduct extensive numerical tests to showcase its applicability and efficiency. Furthermore, we propose several innovative numerical linear algebra algorithms that directly factorize the inverse of a given general matrix. These algorithms hold immense potential for offering preconditioners to accelerate Krylov subspace iterative methods and address large-scale systems of linear equations more efficiently. The comprehensive evaluation and experimental results presented in this paper confirm the practical utility of our proposed algorithms, demonstrating their superiority over classical approaches in terms of computational efficiency. Our research opens up new avenues for exploring advanced matrix inversion techniques, paving the way for improved numerical linear algebra algorithms and the development of effective preconditioners for various applications.

Sparse LU decomposition has been widely used to solve sparse linear systems of equations found in many scientific and engineering applications, such as circuit simulation, power system modeling and computer vision. However, it is considered a computationally expensive factorization tool. While parallel implementations have been explored to accelerate sparse LU decomposition, irregular sparsity patterns often limit their performance gains. Prior FPGA-based accelerators have been customized to domain-specific sparsity patterns of pre-ordered symmetric matrices. In this paper, we present an efficient architecture for sparse LU decomposition that supports both symmetric and asymmetric sparse matrices with arbitrary sparsity patterns. The control structure of our architecture parallelizes computation and pivoting operations. Also, on-chip resource utilization is configured based on properties of the matrices being processed. Our experimental results show a 1:6 to 14x speedup over an opti...

Results are given concerning the LU factorization of H-matrices, and Gaussian elimination with column-diagonaldominant pivoting is shown to be applicable to H-matrices. This algorithm, which uses a symmetric permutation to exchange the most diagonally dominant column of the unreduced submatrix into the pivotal position, is shown to be numerically stable by deriving an upper bound on the growth factor associated with the backward error analysis for Gaussian elimination. 1. INTRODUCTION An nXn realmatrix A=(a,.)isan M-matrix if aij<Oforall i# j and if ReXa for all x~o(A), & e spectrum of A. There are numerous *This author acknowledges support from the Graduate Research Council of Youngstown State University. 'This author acknowledges support from the Natural Sciences and Engineering Research Council under aant A-8214.

This paper presents a deep-pipelined FPGA implementation of real-time ellipse estimation for eye tracking. The system is constructed by the Starburst algorithm on a streamoriented architecture and the RANSAC algorithm without any external memories. In particular, the paper presents comparative results between three different hypothesis generators for the RANSAC algorithm based on Cramer's rule, Gauss-Jordan elimination and LU decomposition. The evaluation results showed that the Gauss-Jordan elimination achieved the highest throughput while the solver with Cramer's rule was the most compact and that our proposed architecture achieved a real-time throughput of 62.5 fps with a single FPGA chip without using any external memories.

In this paper we present new hybrid CPU-GPU routines to accelerate the solution of linear systems, with band coefficient matrix, by off-loading the major part of the computations to the GPU and leveraging highly tuned implementations of the BLAS for the graphics processor. Our experiments with an nVidia S2070 GPU report speed-ups up to 6× for the hybrid band solver based on the LU factorization over analogous CPU-only routines in Intel's MKL. As a practical demonstration of these benefits, we plug the new CPU-GPU codes into a sparse matrix Lyapunov equation solver, showing a 3× acceleration on the solution of a large-scale benchmark arising in model reduction.

An algorithm mainly consisting of a part of Divide and Conquer and the twisted factorization is proposed for bidiagonal SVD. The algorithm costs Oðn 2 Þ flops and is highly parallelizable when singular values are isolated. If strong clusters exist, the singular vector computation needs reorthgonalization. In such case, the cost of the algorithm increases to Oðn 2 þ nk 2 Þ flops and the parallelism may worsen depending on the distribution of singular values. Here k is the size of the largest cluster. The algorithm needs only OðnÞ working memory for every type of matrices.

At the heart of a frontal or multifrontal solver for the solution of sparse symmetric sets of linear equations, there is the need to partially factorize dense matrices (the frontal matrices) and to be able to use their factorizations in subsequent forward and backward substitutions. For a large problem, packing (holding only the lower or upper triangular part) is important to save memory. It has long been recognized that blocking is the key to efficiency and this has become particularly relevant on modern hardware. For stability in the indefinite case, the use of interchanges and 2 × 2 pivots as well as 1 × 1 pivots is equally well established. In this article, the challenge of using these three ideas (packing, blocking, and pivoting) together is addressed to achieve stable factorizations of large real-world symmetric indefinite problems with good execution speed. The ideas are not restricted to frontal and multifrontal solvers and are applicable whenever partial or complete factori...

We propose a hybrid sparse system solver for handling linear systems using algebraic domain decomposition-based techniques. The solver consists of several stages. The first stage uses a reordering scheme that brings as many of the largest matrix elements as possible closest to the main diagonal. This is followed by partitioning the coefficient matrix into a set of overlapped diagonal blocks that contain most of the largest elements of the coefficient matrix. The only constraint here is to minimize the size of each overlap. Separating these blocks into independent linear systems with the constraint of matching the solution parts of neighboring blocks that correspond to the overlaps, we obtain a balance system. This balance system is not formed explicitly and has a size that is much smaller than the original system. Our novel solver requires only a one-time factorization of each diagonal block, and in each outer iteration, obtaining only the upper and lower tips of a solution vector where the size of each tip is equal to that of the individual overlap. This scheme proves to be scalable on clusters of nodes in which each node has a multicore architecture. Numerical experiments comparing the scalability of our solver with direct and preconditioned iterative methods are also presented.

In this paper, we study the Singular Value Decomposition of an arbitrary matrix A  , especially its subspaces of activation, which leads in natural manner to the pseudo inverse of Moore-Bjenhammar-Penrose. Besides, we analyze the compatibility of linear systems and the uniqueness of the corresponding solution and our approach gives the Lanczos classification for these systems.

We present algorithms for the symbolic and numerical factorization phases in the direct solution of sparse unsymmetric systems of linear equations. We have modi ed a classical symbolic factorization algorithm for unsymmetric matrices to inexpensively compute minimal elimination structures. We give an e cient algorithm to compute a near-minimal data-dependency graph that is valid irrespective of the amount of dynamic pivoting performed during numerical 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 signi cantly outperform other similar solvers. We present experimental results to demonstrate the merits of the new algorithms.

This paper presents a new approach for the solution of Linear Programming Problems with the help of LU Factorization Method of matrices. This method is based on the fact that a square matrix can be factorized into the product of unit lower triangular matrix and upper triangular matrix. In this method, we get direct solution without iteration. We also show that this method is better than simplex method.

Realistic and accurate numerical simulations of electrostimulation of tissues and full-body biomodels have been developed and implemented. Typically, whole-body systems are very complex and consist of a multitude of tissues, organs, and subcomponents with diverse properties. From an electrical standpoint, these can be characterized in terms of separate conductivities and permittivities. Accuracy demands good spatial resolution; thus, the overall tissue/animal models need to be discretized into a finegrained mesh. This can lead to a large number of grid points (especially for a three-dimensional entity) and can place prohibitive requirements of memory storage and execution times on computing machines. Here, the authors include a simple yet fast and efficient numerical implementation. It is based on LU decomposition for execution on a cluster of computers running in parallel with distributed storage of the data in a sparse format. In this paper, the details of electrical tissue representation, the fast algorithm, the relevant biomodels, and specific applications to whole-animal studies of electrostimulation are discussed.

A method for solving systems of linear equations is presented based on direct decomposition of the coefficient matrix using the form LAX LB B′ = =. Elements of the reducing lower triangular matrix L can be determined using either row wise or column wise operations and are demonstrated to be sums of permutation products of the Gauss pivot row multipliers. These sums of permutation products can be constructed using a tree structure that can be easily memorized or alternatively computed using matrix products. The method requires only storage of the L matrix which is half in size compared to storage of the elements in the LU decomposition. Equivalence of the proposed method with both the Gauss elimination and LU decomposition is also shown in this paper.

In this paper, authors present their work on field-programmable gate array (FPGA) hardware implementation of proposed direction of arrival estimation algorithms employing LU factorization. Both L and U matrices were considered in computing the angle estimates. Hardware implementation was done on a Virtex-5 FPGA and its experimental verification was performed using National Instruments PXI platform which provides hardware modules for data acquisition, RF down-conversion, digitization, etc. A uniform linear array consisting of four antenna elements was deployed at the receiver. LabVIEW FPGA modules with high throughput math functions were used for implementing the proposed algorithms. MATLAB simulations of the proposed algorithms were also performed to validate the efficacy of the proposed algorithms prior to hardware implementation of the same. Both MATLAB simulation and experimental verification establish the superiority of the proposed methods over existing methods reported in the literature, such as QR decomposition-based implementations. FPGA compilation results report low resource usage and faster computation time compared with the QR-based hardware implementation. Performance comparison in terms of estimation accuracy, percentage resource utilization, and processing time is also presented for different data and matrix sizes. INDEX TERMS FPGAs, LU factorization, NI PXI platform, pipelined architecture.

In this paper, we postulate a new decomposition theorem of a matrix A into two matrices, namely, a lower triangular matrix M, in which all entries are determinants, and an upper triangular matrix U whose entries are also in determinant form. From a well-known theorem on the pivot elements of the Doolittle-Gauss elimination process, we deduce a corollary to obtain a diagonal matrix D. With it, we scale the elementary lower triangular matrix of the Doolittle- Gauss elimination process and deduce a new elementary lower triangular matrix. Applying this linear transformation to A by means of both minimum and complete pivoting strategies, we obtain the determinant of A as if it had been calculated by means of a Laplace expansion. If we apply this new linear transformation and the above pivot strategy to an augmented matrix (A|b), we obtain a Cramer&#39;s solution of the linear

Sparse parallel factorization is among the most complicated and irregular algorithms to analyze and optimize. Performance depends both on system characteristics such as the floating point rate, the memory hierarchy, and the interconnect performance, as well as input matrix characteristics such as such as the number and location of nonzeros. We present LUsim, a simulation framework for modeling the performance of sparse LU factorization. Our framework uses micro-benchmarks to calibrate the parameters of machine characteristics and additional tools to facilitate real-time performance modeling. We are using LUsim to analyze an existing parallel sparse LU factorization code, and to explore a latency tolerant variant. We developed and validated a model of the factorization in SuperLU DIST, then we modeled and implemented a new variant of SuperLU DIST, replacing a blocking collective communication phase with a non-blocking asynchronous point-to-point one. Our strategy realized a mean improvement of 11% over a suite of test matrices.

In this work, the solution of a large sparse linear system of equations with an arbitrary sparsity pattern is obtained by using LU-decomposition method as well as numerical structure approach. The LU-decomposition method is based on Doolittle's method while the numerical structure approach is based on Cramer's rule. The numerical structure approach produces direct solution without facing fill-in problems as encountered in LU-decomposition. In order to reduce the 'fill-ins' in the decomposition, the powers of a Boolean matrix, obtained from the coefficient matrix A are taken so that the 'fill-ins' in the structure of A can be known in advance. The position of fill-ins in A are thus determined in the best choice manner, that is, it is very effective and memory-wise cheap. We also outline a method by using numerical structure with reduced computation efforts. Finally, experiments are performed on eight examples to compare the efficiency of the proposed methods. The results obtained are reported in a table. It is found that the LU-decomposition method is much better than numerical structure. The usefulness of numerical structure approach is also discussed. (~

In this paper, we postulate a new decomposition theorem of a matrix A into two matrices, namely, a lower triangular matrix M, in which all entries are determinants, and an upper triangular matrix U whose entries are also in determinant form. From a well-known theorem on the pivot elements of the Doolittle-Gauss elimination process, we deduce a corollary to obtain a diagonal matrix D. With it, we scale the elementary lower triangular matrix of the Doolittle- Gauss elimination process and deduce a new elementary lower triangular matrix. Applying this linear transformation to A by means of both minimum and complete pivoting strategies, we obtain the determinant of A as if it had been calculated by means of a Laplace expansion. If we apply this new linear transformation and the above pivot strategy to an augmented matrix (A|b), we obtain a Cramer&#39;s solution of the linear

The comprehensive LU decomposition of a parametric matrix consists of a case analysis of the LU factors for each specialization of the parameters. Special cases can be discontinuous with respect to the parameters, the discontinuities being triggered by zero pivots encountered during factorization. For polynomial matrices, we describe an implementation of comprehensive LU decomposition in Maple, using the RegularChains package.