Algorithms For Sparse Matrices Research Articles

Algorithm xxx: Faster Randomized SVD with Dynamic Shifts

Aiming to provide a faster and convenient truncated SVD algorithm for large sparse matrices from real applications (i.e. for computing a few of largest singular values and the corresponding singular vectors), a dynamically shifted power iteration technique is applied to improve the accuracy of the randomized SVD method. This results in a d yn a mic sh ifts based randomized SVD (dashSVD) algorithm, which also collaborates with the skills for handling sparse matrices. An accuracy-control mechanism is included in the dashSVD algorithm to approximately monitor the per vector error bound of computed singular vectors with negligible overhead. Experiments on real-world data validate that the dashSVD algorithm largely improves the accuracy of randomized SVD algorithm or attains same accuracy with fewer passes over the matrix, and provides an efficient accuracy-control mechanism to the randomized SVD computation, while demonstrating the advantages on runtime and parallel efficiency. A bound of the approximation error of the randomized SVD with the shifted power iteration is also proved.

Optimization of Invariant Ellipsoid Technique for Sparse Controllers Design

The paper deals with the method for the design of linear controllers with sparse state feedback matrices for control the plants under conditions of unknown and bounded disturbances. The importance of the sparsity property in feedback can be explained by two factors. First, by minimizing the columnar norm of the feedback matrix in the control it becomes pos- sible to use a minimum number of measuring devices. Secondly, by minimizing the row norm of the feedback matrix, the required number of executive (control) devices is minimized. Both properties, if they are achievable in the synthesis of the controller, reduce the cost of the system and improve the fault tolerance and quality of regulation by reducing the structural complexity. The search algorithm for sparse matrices is based on the method of invariant ellipsoids and is formulated as a solution to a system of linear matrix inequalities with additional constraints. A special set of optimization conditions is proposed which for a disturbed system minimizes overshoot and overshoots in transient processes of the disturbed closed- loop system simultaneously with minimizing errors in the steady state. The proposed method also assumes the possibility of minimizing both the row norm of the feedback matrix and the column one, while preserving the robustness properties, which makes it possible to solve the sparse control problem (a sparse control is understood as a linear controller with a sparse feedback matrix). The efficiency of the proposed control scheme is confirmed by the results of computer modeling and comparison with some existing ones.

MATRIX FACTORIZATION OF BIG DATA IN THE INDUSTRIAL SYSTEMS

The creation of new technologies and their implementation in various fields necessitated Big Data processing and storage. In industrial systems, modernization means the use of a large number of smart devices that perform specialized functions. Data from such devices are used to control the system and automate production processes. A change in the parameters of individual components of the manufacturing system may indicate the need to adjust the global management strategy. The intelligent industrial systems main characteristics were defined in the paper. The Industrial Internet of Things concept and the relevance of the modernization problem for manufacturing were analyzed. The problems of processing Big Data in Industrial Internet of Things systems were examined in the paper. The use of recommendation systems for quickly finding relationships between users and production services was considered. The advantages of Big Data analysis by recommendation systems, which have a favourable effect on industrial enterprise efficiency were given. The use of SVD and FunkSVD matrix factorization algorithms for processing sparse data matrices was analyzed. The possibility of optimizing arrays of information, choosing the most important, and rejecting redundancy with the help of the above algorithms was determined. The proposed algorithms were simulated. The advantages of FunkSVD for working with sparse data were assigned. It was found that the FunkSVD algorithm processes the data in a shorter time than SVD, but this does not affect the accuracy of the result. The SVD is also more difficult to implement and it requires more computing resources was established. It has been shown that FunkSVD uses a lot of data to determine the relationships between it and make recommendations about the products most likely to be of interest to users. To increase the efficiency of processing large sets of information the FunkSVD algorithm was improved in such a way that it uses fewer data to generate recommendations. Based on the results of the research, the modified method works faster than the non-modified one but retains high calculation accuracy, which is important for work in recommender systems. The possibility of providing recommendations to users of industrial systems in a shorter period, thus improving their relevance, was revealed. It was proposed to continue research for finding the optimal parameters of the FunkSVD algorithm for Big Data processing.

Revisiting Co-Occurring Directions: Sharper Analysis and Efficient Algorithm for Sparse Matrices

We study the streaming model for approximate matrix multiplication (AMM). We are interested in the scenario that the algorithm can only take one pass over the data with limited memory. The state-of-the-art deterministic sketching algorithm for streaming AMM is the co-occurring directions (COD), which has much smaller approximation errors than randomized algorithms and outperforms other deterministic sketching methods empirically. In this paper, we provide a tighter error bound for COD whose leading term considers the potential approximate low-rank structure and the correlation of input matrices. We prove COD is space optimal with respect to our improved error bound. We also propose a variant of COD for sparse matrices with theoretical guarantees. The experiments on real-world sparse datasets show that the proposed algorithm is more efficient than baseline methods.

Fast Algorithms for Conducting Large-Scale GWAS of Age-at-Onset Traits Using Cox Mixed-Effects Models.

Age-at-onset is one of the critical traits in cohort studies of age-related diseases. Large-scale genome-wide association studies (GWAS) of age-at-onset traits can provide more insights into genetic effects on disease progression and transitions between stages. Moreover, proportional hazards (or Cox) regression models can achieve higher statistical power in a cohort study than a case-control trait using logistic regression. Although mixed-effects models are widely used in GWAS to correct for sample dependence, application of Cox mixed-effects models (CMEMs) to large-scale GWAS is so far hindered by intractable computational cost. In this work, we propose COXMEG, an efficient R package for conducting GWAS of age-at-onset traits using CMEMs. COXMEG introduces fast estimation algorithms for general sparse relatedness matrices including, but not limited to, block-diagonal pedigree-based matrices. COXMEG also introduces a fast and powerful score test for dense relatedness matrices, accounting for both population stratification and family structure. In addition, COXMEG generalizes existing algorithms to support positive semidefinite relatedness matrices, which are common in twin and family studies. Our simulation studies suggest that COXMEG, depending on the structure of the relatedness matrix, is orders of magnitude computationally more efficient than coxme and coxph with frailty for GWAS. We found that using sparse approximation of relatedness matrices yielded highly comparable results in controlling false-positive rate and retaining statistical power for an ethnically homogeneous family-based sample. By applying COXMEG to a study of Alzheimer's disease (AD) with a Late-Onset Alzheimer's Disease Family Study from the National Institute on Aging sample comprising 3456 non-Hispanic whites and 287 African Americans, we identified the APOE ε4 variant with strong statistical power (P = 1e-101), far more significant than that reported in a previous study using a transformed variable and a marginal Cox model. Furthermore, we identified novel SNP rs36051450 (P = 2e-9) near GRAMD1B, the minor allele of which significantly reduced the hazards of AD in both genders. These results demonstrated that COXMEG greatly facilitates the application of CMEMs in GWAS of age-at-onset traits.

GPU-Accelerated Parallel Sparse LU Factorization Method for Fast Circuit Analysis

Lower upper (LU) factorization for sparse matrices is the most important computing step for circuit simulation problems. However, parallelizing LU factorization on the graphic processing units (GPUs) turns out to be a difficult problem due to intrinsic data dependence and irregular memory access, which diminish GPU computing power. In this paper, we propose a new sparse LU solver on GPUs for circuit simulation and more general scientific computing. The new method, which is called GPU accelerated LU factorization (GLU) solver (for GPU LU), is based on a hybrid right-looking LU factorization algorithm for sparse matrices. We show that more concurrency can be exploited in the right-looking method than the left-looking method, which is more popular for circuit analysis, on GPU platforms. At the same time, the GLU also preserves the benefit of column-based left-looking LU method, such as symbolic analysis and column-level concurrency. We show that the resulting new parallel GPU LU solver allows the parallelization of all three loops in the LU factorization on GPUs. While in contrast, the existing GPU-based left-looking LU factorization approach can only allow parallelization of two loops. Experimental results show that the proposed GLU solver can deliver $5.71\times $ and $1.46\times $ speedup over the single-threaded and the 16-threaded PARDISO solvers, respectively, $19.56\times $ speedup over the KLU solver, $47.13\times $ over the UMFPACK solver, and $1.47\times $ speedup over a recently proposed GPU-based left-looking LU solver on the set of typical circuit matrices from the University of Florida (UFL) sparse matrix collection. Furthermore, we also compare the proposed GLU solver on a set of general matrices from the UFL, GLU achieves $6.38\times $ and $1.12\times $ speedup over the single-threaded and the 16-threaded PARDISO solvers, respectively, $39.39\times $ speedup over the KLU solver, $24.04\times $ over the UMFPACK solver, and $2.35\times $ speedup over the same GPU-based left-looking LU solver. In addition, comparison on self-generated $RLC$ mesh networks shows a similar trend, which further validates the advantage of the proposed method over the existing sparse LU solvers.

An Explicit Link between Gaussian Fields and Gaussian Markov Random Fields: The Stochastic Partial Differential Equation Approach

SummaryContinuously indexed Gaussian fields (GFs) are the most important ingredient in spatial statistical modelling and geostatistics. The specification through the covariance function gives an intuitive interpretation of the field properties. On the computational side, GFs are hampered with the big n problem, since the cost of factorizing dense matrices is cubic in the dimension. Although computational power today is at an all time high, this fact seems still to be a computational bottleneck in many applications. Along with GFs, there is the class of Gaussian Markov random fields (GMRFs) which are discretely indexed. The Markov property makes the precision matrix involved sparse, which enables the use of numerical algorithms for sparse matrices, that for fields in ℝ2 only use the square root of the time required by general algorithms. The specification of a GMRF is through its full conditional distributions but its marginal properties are not transparent in such a parameterization. We show that, using an approximate stochastic weak solution to (linear) stochastic partial differential equations, we can, for some GFs in the Matérn class, provide an explicit link, for any triangulation of ℝd, between GFs and GMRFs, formulated as a basis function representation. The consequence is that we can take the best from the two worlds and do the modelling by using GFs but do the computations by using GMRFs. Perhaps more importantly, our approach generalizes to other covariance functions generated by SPDEs, including oscillating and non-stationary GFs, as well as GFs on manifolds. We illustrate our approach by analysing global temperature data with a non-stationary model defined on a sphere.

Cache friendly sparse matrix-vector multiplication

Sparse matrix-vector multiplication or SpMxV is an important kernel in scientific computing. For example, in the conjugate gradient method, where SpMxV is the main computational step. Though the total number of arithmetic operations in SpMxV is fixed, reducing the probability of cache misses per operation is still a challenging area of research. In this work, we present a new column ordering algorithm for sparse matrices. We analyze the cache complexity of SpMxV when A is ordered by our technique. The numerical experiments, with very large test matrices, clearly demonstrate the performance gains rendered by our proposed technique.

Fast algorithms for placing large entries along the diagonal of a sparse matrix

Solving a sparse system of linear equations A x = b is one of the most fundamental operations inside any circuit simulator. The equations/rows in the matrix A are often rearranged/permuted before factorization and applying direct or iterative methods to obtain the solution. Permuting the rows of the matrix A so that the entries with large absolute values lie on the diagonal has several advantages like better numerical stability for direct methods (e.g., Gaussian elimination) and faster convergence for indirect methods (such as the Jacobi method). Duff (2009) [3] has formulated this as a weighted bipartite matching problem (the MC64 algorithm). In this paper we improve the performance of the MC64 algorithm with a new labeling technique which improves the asymptotic complexity of updating dual variables from O ( | V | + | E | ) to O ( | V | ) , where | V | is the order of the matrix A and | E | is the number of non-zeros. Experimental results from using the new algorithm, when benchmarked with both industry benchmarks and UFL sparse matrix collection, are very promising. Our algorithm is more than 60 times faster (than Duff’s algorithm) for sparse matrices with at least a million non-zeros.

On the choice of a multiplication algorithm for polynomials and polynomial matrices

We investigate multiplication algorithms for dense and sparse polynomials and polynomial matrices over different numerical domains and obtain expressions for the complexity of multiplication of polynomials and polynomial matrices understood as the expectation of the number of arithmetic operations. These expressions for a set of parameters of practical interest are tabulated. The results of experiments with the corresponding programs are presented. Bibliography: 8 titles.

Frequency-Domain Linearized Euler Model for Turbomachinery Noise Radiation Through Engine Exhaust

A numerical model for the exhaust noise radiation problem is presented. In the model, it is assumed that an incoming wave is propagating through the exhaust nozzle, or the fan duct, and radiating outside. The near-field propagation is based on the solution of the linearized Euler equations in the frequency domain: for each wave number, a linearized Euler problem is solved using a finite element method on unstructured grids for arbitrarily shaped axisymmetric geometries. The frequency-domain approach enables the suppression of the Kelvin-Helmholtz instability waves. Moreover, each single calculation, limited to a single frequency, is well suited to the exhaust noise radiation problem in which the incoming wave can be treated as a superposition of elementary duct modes. To reduce the memory requirements, a continuous Galerkin formulation with linear triangular and quadrangular elements is employed and the global matrix inversion is performed with a direct solver based on a parallel memory distributed multifrontal algorithm for sparse matrices. The acoustic near field is then radiated in the far field using the formulation of Ffowcs Williams and Hawkings. Numerical calculations for a validation test case, the Munt problem, and two turbomachinery configurations are compared with analytical solutions and experimental data.

Sheafhom: Software for Sparse Integer Matrices

Sheafhom is a free software package for large-scale computations in the category of finitely-generated modules over the integers and related rings. Its front end is a language for problems in algebraic topology and geometry. These problems come down to sparse systems of linear equations over the integers. Sheafhom's back end solves the sparse systems, with emphasis on avoiding fill-in and integer explosion. We survey and compare algorithms for integer sparse matrices, and we present implementation techniques in Common Lisp. The final section finds the quotient of the free abelian group on 26 letters by sets of words in the dictionary, in the spirit of [8].

Parallel Symbolic Factorization for Sparse LU with Static Pivoting

This paper presents the design and implementation of a memory scalable parallel symbolic factorization algorithm for general sparse unsymmetric matrices. Our parallel algorithm uses a graph partitioning approach, applied to the graph of $|A|+|A|^T$, to partition the matrix in such a way that is good for sparsity preservation as well as for parallel factorization. The partitioning yields a so-called separator tree which represents the dependencies among the computations. We use the separator tree to distribute the input matrix over the processors using a block cyclic approach and a subtree to subprocessor mapping. The parallel algorithm performs a bottom-up traversal of the separator tree. With a combination of right-looking and left-looking partial factorizations, the algorithm obtains one column structure of L and one row structure of U at each step. The algorithm is implemented in C and MPI. From a performance study on large matrices, we show that the parallel algorithm significantly reduces the memory requirement of the symbolic factorization step, as well as the overall memory requirement of the parallel solver. It also often reduces the runtime of the sequential algorithm, which is already relatively small. In general, the parallel algorithm prevents the symbolic factorization step from being a time or memory bottleneck of the parallel solver.

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.

Balancing sparse matrices for computing eigenvalues

Applying a permuted diagonal similarity transform DPAP T D −1 to a matrix A before calculating its eigenvalues can improve the speed and accuracy with which the eigenvalues are computed. This is often called balancing. This paper describes several balancing algorithms for sparse matrices and compares them against each other and the traditional dense algorithm. We first discuss our sparse implementation of the dense algorithm; our code is faster than the dense algorithm when the density of the matrix is no more than approximately .5, and is much faster for large, sparse matrices. We next describe a set of randomized balancing algorithms for matrices that are not given explicitly, i.e. given a vector x, we can compute only Ax and perhaps A T x. We motivate these Krylov-based algorithms using Perron–Frobenius theory. Results are given comparing the Krylov-based algorithms to each other and to the sparse and dense direct balancing algorithms, looking at norm reduction, running times, and the accuracy of eigenvalues computed after a matrix is balanced. We conclude that sparse balancing algorithms are efficient preconditioners for eigensolvers.

Parallel direct method for linear equations

A simple direct method and its variants, based on matrix multiplications, to invert square non-singular matrices (need not be positive definite) and to solve systems of linear equations are developed. Since the method and its variants involve only matrix multiplications they are straightforward to parallelise and hence have an advantage over some existing well known direct methods. Theoretical background, analysis of the proposed method and the complexity of the algorithm for sparse and dense matrices are given. Two significantly different parallel algorithms for dense and sparse matrices are developed. In the case of the dense matrix standard parallelisation of matrix multiplications was undertaken. However, for sparse matrices the standard parallel matrix multiplication is not effective which led us to develop a parallel algorithm different from that for the dense matrix. The performances of our algorithms are tested via numerical examples.

Level 3 basic linear algebra subprograms for sparse matrices

This article proposes a set of Level 3 Basic Linear Algebra Subprograms and associated kernels for sparse matrices. A major goal is to design and develop a common framework to enable efficient, and portable, implementations of iterative algorithms for sparse matrices on high-performance computers. We have designed the routines to shield the developer of mathematical software from most of the complexities of the various data structures used for sparse matrices. We have kept the interface and suite of codes as simple as possible while at the same time including sufficient functionality to cover most of the requirements of iterative solvers and sufficient flexibility to cover most sparse matrix data structures. An important aspect of our framework is that it can be easily extended to incorporate new kernels if the need arises. We discuss the design, implementation, and use of subprograms for the multiplication of a fully matrix by a sparse one and for the solution of sparse triangular systems with one or more (full) right-hand sides. We include a routine for checking the input data, generating a new sparse data structure from the input, and scaling a sparse matrix. The new data structure for the transformation can be specified by the user or can be chosen automatically by vendors to be efficient on their machines. We also include a routine for permuting the columns of a sparse matrix and one for permuting the rows of a full matrix.

Overlapping Domain Decomposition Algorithms for General Sparse Matrices

Domain decomposition methods for finite element problems using a partition based on the underlying finite element mesh have been extensively studied. In this paper, we discuss algebraic extensions of the class of overlapping domain decomposition algorithms for general sparse matrices. The subproblems are created with an overlapping partition of the graph corresponding to the sparsity structure of the matrix. These algebraic domain decomposition methods are especially useful for unstructured mesh problems. We also discuss some difficulties encountered in the algebraic extension, particularly the issues related to the coarse solver.

Overlapping Domain Decomposition Algorithms for General Sparse Matrices

Domain decomposition methods for finite element problems using a partition based on the underlying finite element mesh have been extensively studied. In this paper, we discuss algebraic extensions of the class of overlapping domain decomposition algorithms for general sparse matrices. The subproblems are created with an overlapping partition of the graph corresponding to the sparsity structure of the matrix. These algebraic domain decomposition methods are especially useful for unstructured mesh problems. We also discuss some difficulties encountered in the algebraic extension, particularly the issues related to the coarse solver.

The Use of Pointwise Interpolation in Domain Decomposition Methods with Nonnested Meshes

In this paper, the author develops a new technique and a corresponding theory for Schwarz-type overlapping domain decomposition methods for solving large sparse linear systems which arise from finite element discretization of elliptic partial differential equations. The theory provides an optimal convergence of an additive Schwarz algorithm that is constructed with a nonnested coarse space and a not necessarily shape regular subdomain partitioning. The theory is also applicable to the graph partitioning algorithms recently developed by Cai and Saad [Overlapping domain decomposition algorithms for general sparse matrices, Preprint 93-027, Univ. of Minnesota] and Farhat and Lesoine [Internat. J. Numer. Methods Engrg., 36 (1993), pp. 745–764] for problems defined on unstructured meshes.