Engineering Papers⌕ Search

SEARCH · Engineering Papers

Results for “sparse matrices”

Search indexed NASA NTRS and DOE OSTI research on propulsion, heat transfer, battery materials and energy systems. Follow report and document links to the original sources.

Quote a phrase for an exact phrase match. Source license links do not imply unrestricted reuse.

At least 163 records · Page 9

Sparse Approximate Multifrontal Factorization with Butterfly Compression for High-Frequency Wave Equations

In this work, we present a fast and approximate multifrontal solver for large-scale sparse linear systems arising from finite-difference, finite-volume or finite-element discretization of high-frequency wave equations. The proposed solver leverages the butterfly algorithm and its hierarchical matrix extension for compressing and factorizing large frontal matrices via graph-distance guided entry evaluation or randomized matrix-vector multiplication-based schemes. Complexity analysis and numerical experiments demonstrate $\mathcal{O}(N\log^2 N)$ computation and $\mathcal{O}(N)$ memory complexity when applied to an $N\times N$ sparse system arising from 3D high-frequency Helmholtz and Maxwell problems.

97 MATHEMATICS AND COMPUTING↗

GCoD: Graph Convolutional Network Acceleration via Dedicated Algorithm and Accelerator Co-Design

Graph Convolutional Networks (GCNs) have emerged as the state-of-the-art graph learning model. However, it remains notoriously challenging to inference GCNs over large graph datasets, limiting their application to large real-world graphs and hindering the exploration of deeper and more sophisticated GCN graphs. This is because real-world graphs can be extremely large and sparse. Furthermore, the node degree of GCNs tends to follow the power-law distribution and therefore have highly irregular adjacency matrices, resulting in prohibitive inefficiencies in both data processing and movement and thus substantially limiting the achievable GCN acceleration efficiency. To this end, this paper proposes the first GCN algorithm and accelerator Co-Design framework dubbed GCoD which can largely alleviate the aforementioned GCN irregularity and boost GCNs' inference efficiency. Specifically, on the algorithm level, GCoD integrates a divide and conquer GCN training strategy that polarizes the graphs to be either denser or sparser in local neighborhoods without compromising the model accuracy, resulting in graph adjacency matrices that (mostly) have merely two levels of workload and enjoys largely enhanced regularity and thus ease of acceleration. On the hardware level, we further develop a dedicated two-pronged accelerator with a separated engine to process each of the aforementioned workloads, further boosting the overall utilization and acceleration efficiency. Extensive experiments and ablation studies validate that our GCoD consistently outperforms state-of-the-art designs in terms of accelerator efficiency while maintaining or even improving the task accuracy. Additionally, we visualize GCoD trained graph adjacency matrices to better understand its advantages. All codes and pre-trained models will be released upon acceptance.

You, Haoran↗

A parallel-plate avalanche counter for the prompt fission neutron spectrum measurement

Neutrons, gamma’s, and fission fragments are among the prompt fission observables. Their precision measurements are fundamental to advance our understanding of fission and any application of fission. For our programmatic need, the χ matrix, which defines as the outgoing neutron spectrum as a function of incoming neutron energy, is sparsely populated and poorly measured for 233,235,238 U and 239,240 Pu. To improve the quality of those χ matrices, a joint LANLLLNL project was undertaken to measure their prompt fission neutron spectra at the LANSCE/WNR facility using neutron detector arrays and a charged-particle detector for the fission-fragment detection. A typical setup is shown in Fig. 1, where one of neutron detector array with a charged-particle detector is given.

73 NUCLEAR PHYSICS AND RADIATION PHYSICS↗

Highly parallel sparse Cholesky factorization

Several fine grained parallel algorithms were developed and compared to compute the Cholesky factorization of a sparse matrix. The experimental implementations are on the Connection Machine, a distributed memory SIMD machine whose programming model conceptually supplies one processor per data element. In contrast to special purpose algorithms in which the matrix structure conforms to the connection structure of the machine, the focus is on matrices with arbitrary sparsity structure. The most promising algorithm is one whose inner loop performs several dense factorizations simultaneously on a 2-D grid of processors. Virtually any massively parallel dense factorization algorithm can be used as the key subroutine. The sparse code attains execution rates comparable to those of the dense subroutine. Although at present architectural limitations prevent the dense factorization from realizing its potential efficiency, it is concluded that a regular data parallel architecture can be used efficiently to solve arbitrarily structured sparse problems. A performance model is also presented and it is used to analyze the algorithms.

Gilbert, John R.↗

Impact of Reordering on the LU Factorization Performance of Bordered Block-Diagonal Sparse Matrix

Power engineers rely on computer-based simulation tools to assess grid performance and ensure security. At the core of these tools are solvers for sparse linear equations. When transformed into a bordered block-diagonal (BBD) structure, part of the sparse linear equation solving can be parallelized. This work focuses on using the Schur-complement-based method for LU factorization on BBD matrices, specifically, Jacobian matrices from large-scale systems. Our findings show that the natural ordering method outperforms the default ordering method in computational performance for each block of the BBD matrix. This observation is validated using synthetic 25k-bus and 70k-bus cases, showing a speedup of up to 38% when using natural ordering without permutation. Additionally, the impact of the number of partitions is studied, and the result shows that computational performance improves with more, smaller partitions in the BBD matrices.

BBD matrix↗

Some iterative schemes for transonic potential flows

The minimal residual (MR) method for the numerical solution of transonic potential flows is closely related to the conjugate gradient method, which has found widespread use in the solution of large sparse, symmetric, and positive-definite linear equations. The primary advantage of the MR method is its applicability to both symmetric and nonsymmetric matrices.

Wong, Y. S.↗

An implementation of the look-ahead Lanczos algorithm for non-Hermitian matrices, part 2

It is shown how the look-ahead Lanczos process (combined with a quasi-minimal residual QMR) approach) can be used to develop a robust black box solver for large sparse non-Hermitian linear systems. Details of an implementation of the resulting QMR algorithm are presented. It is demonstrated that the QMR method is closely related to the biconjugate gradient (BCG) algorithm; however, unlike BCG, the QMR algorithm has smooth convergence curves and good numerical properties. We report numerical experiments with our implementation of the look-ahead Lanczos algorithm, both for eigenvalue problem and linear systems. Also, program listings of FORTRAN implementations of the look-ahead algorithm and the QMR method are included.

Freund, Roland W.↗

A hill-sliding strategy for initialization of Gaussian clusters in the multidimensional space

A hill sliding technique was devised to extract Gaussian clusters from the multivariate probability density estimate of sample data for the first step of iterative unsupervised classification. Each cluster was assumed to posses a unimodal normal distribution. A clustering function proposed distinguished elements of a cluster under formation from the rest in the feature space. Initial clusters were extracted one by one according to the hill sliding tactics. A dimensionless cluster compactness parameter was proposed as a universal measure of cluster goodness and used satisfactorily in test runs with LANDSAT multispectral scanner data. The normalized divergence, defined by the cluster divergence divided by the entropy of the entire sample data, was utilized as a general separability measure between clusters. An overall clustering objective function was set forth in terms of cluster covariance matrices, from which the cluster compactness measure could be deduced. Minimal improvement of initial data partitioning was evaluated by this objective function in eliminating scattered sparse data points. The hill sliding clustering technique developed herein has the potential applicability to decomposition any multivariate mixture distribution into a number of unimodal distributions when an appropriate distribution function to the data set is employed.

Park, J. K.↗

Sparse Linear Solvers for Large-scale Electromagnetic Transient Simulations

Linear solvers form the basis for electromagnetic transient (EMT) simulations. There is a need to speed up EMT simulations as larger regions are analyzed using EMT simulations. For the same, the performance of linear solvers plays an important role. Exploiting the sparsity of the matrices generated in EMT simulations could assist with speed-up. Scalability is also crucial as power grids expand, demanding solutions capable of accommodating the increasing system size. Recent studies from the North American Electric Reliability Corporation (NERC) increasingly emphasize that EMT simulation models of the power grid will grow larger with the inclusion of power electronics components. Parallelisms in sparsity patterns exploit modern central processing units (CPUs), multi-core CPUs, and graphics processing units (GPUs) architectures in sparse solver designs. Therefore, this paper explores publicly available existing linear solvers and investigates their efficiency in large-scale power grid simulations. A large-scale power grid is developed by increasing the size of the IEEE 39 bus test system to up to 39000 bus systems.

Hsu, Kuan-Chieh↗

Block encoding bosons by signal processing

Block Encoding (BE) is a crucial subroutine in many modern quantum algorithms, including those with near-optimal scaling for simulating quantum many-body systems, which often rely on Quantum Signal Processing (QSP). Currently, the primary methods for constructing BEs are the Linear Combination of Unitaries (LCU) and the sparse oracle approach. In this work, we demonstrate that QSP-based techniques, such as Quantum Singular Value Transformation (QSVT) and Quantum Eigenvalue Transformation for Unitary Matrices (QETU), can themselves be efficiently utilized for BE implementation. Specifically, we present several examples of using QSVT and QETU algorithms, along with their combinations, to block encode Hamiltonians for lattice bosons, an essential ingredient in simulations of high-energy physics. We also introduce a straightforward approach to BE based on the exact implementation of Linear Operators Via Exponentiation and LCU (LOVE-LCU). We find that, while using QSVT for BE results in the best asymptotic gate count scaling with the number of qubits per site, LOVE-LCU outperforms all other methods for operators acting on up to qubits, highlighting the importance of concrete circuit constructions over mere comparisons of asymptotic scalings. Using LOVE-LCU to implement the BE, we simulate the time evolution of single-site and two-site systems in the lattice theory using the Generalized QSP algorithm and compare the gate counts to those required for Trotter simulation.

Kane, Christopher F↗

An Optimized Multicolor Point-Implicit Solver for Unstructured Grid Applications on Graphics Processing Units

In the field of computational fluid dynamics, the Navier-Stokes equations are often solved using an unstructuredgrid approach to accommodate geometric complexity. Implicit solution methodologies for such spatial discretizations generally require frequent solution of large tightly-coupled systems of block-sparse linear equations. The multicolor point-implicit solver used in the current work typically requires a significant fraction of the overall application run time. In this work, an efficient implementation of the solver for graphics processing units is proposed. Several factors present unique challenges to achieving an efficient implementation in this environment. These include the variable amount of parallelism available in different kernel calls, indirect memory access patterns, low arithmetic intensity, and the requirement to support variable block sizes. In this work, the solver is reformulated to use standard sparse and dense Basic Linear Algebra Subprograms (BLAS) functions. However, numerical experiments show that the performance of the BLAS functions available in existing CUDA libraries is suboptimal for matrices representative of those encountered in actual simulations. Instead, optimized versions of these functions are developed. Depending on block size, the new implementations show performance gains of up to 7x over the existing CUDA library functions.

Zubair, Mohammad↗

Fast Solution in Sparse LDA for Binary Classification

An algorithm that performs sparse linear discriminant analysis (Sparse-LDA) finds near-optimal solutions in far less time than the prior art when specialized to binary classification (of 2 classes). Sparse-LDA is a type of feature- or variable- selection problem with numerous applications in statistics, machine learning, computer vision, computational finance, operations research, and bio-informatics. Because of its combinatorial nature, feature- or variable-selection problems are NP-hard or computationally intractable in cases involving more than 30 variables or features. Therefore, one typically seeks approximate solutions by means of greedy search algorithms. The prior Sparse-LDA algorithm was a greedy algorithm that considered the best variable or feature to add/ delete to/ from its subsets in order to maximally discriminate between multiple classes of data. The present algorithm is designed for the special but prevalent case of 2-class or binary classification (e.g. 1 vs. 0, functioning vs. malfunctioning, or change versus no change). The present algorithm provides near-optimal solutions on large real-world datasets having hundreds or even thousands of variables or features (e.g. selecting the fewest wavelength bands in a hyperspectral sensor to do terrain classification) and does so in typical computation times of minutes as compared to days or weeks as taken by the prior art. Sparse LDA requires solving generalized eigenvalue problems for a large number of variable subsets (represented by the submatrices of the input within-class and between-class covariance matrices). In the general (fullrank) case, the amount of computation scales at least cubically with the number of variables and thus the size of the problems that can be solved is limited accordingly. However, in binary classification, the principal eigenvalues can be found using a special analytic formula, without resorting to costly iterative techniques. The present algorithm exploits this analytic form along with the inherent sequential nature of greedy search itself. Together this enables the use of highly-efficient partitioned-matrix-inverse techniques that result in large speedups of computation in both the forward-selection and backward-elimination stages of greedy algorithms in general.

Moghaddam, Baback↗

A Performance and Energy Study of GPU-Resident Preconditioners for Conjugate Gradient Solvers: In the Context of Existing and Novel Approaches

Optimizing a particular subprogram out of the set of Basic (sparse) Linear Algebra Subprograms (BLAS) for a given architecture is a common topic of research. In applications, however, these BLAS functions rarely appear in isolation; usually, many of them are used together, in various combinations and with varying inputs. As the need to solve a large, sparse linear system is ubiquitous throughout HPC applications, linear solvers constitute a realistic, sufficiently complex and well-defined representative use case for composite BLAS routines. To this end, based on a representative set of matrices drawn from a diverse set of fields, we present a framework to study, from the performance and energy perspective, the efficacy of GPU- resident parallel Conjugate Gradient (CG) linear solver with different preconditioner options, including Gauss-Seidel, Jacobi, and incomplete Cholesky. We also propose a novel GPU-based preconditioner, in which the triangular solves are approximated by an iterative process. The development of this preconditioner was motivated by solving large graph Laplacian linear systems, for which the existing preconditioners either perform slow on GPU-based platforms or are not applicable. We compare the performance of these preconditioners on different hardware accelerator architectures, i.e., AMD MI250X, MI100, Nvidia A100, V100, and Jetson. Our experiments reveal performance trade-offs and provide information on how to select the best strategy for the given linear system, dictated by its properties, and the platform of interest. We demonstrate the application of our novel preconditioner for solving CG and graph Laplacian systems. Overall, the framework can be utilized as a benchmark to guide informed decisions in choosing a specific preconditioner, i.e., whether it is better to rely on the performance of a triangular solver or on the performance of sparse matrix-vector product. Finally, by considering power consumption to solve the linear systems, we report the energy footprint for the solvers.

Preconditioned Conjugate Gradient, GPUs, iterative↗

Multi-color incomplete Cholesky conjugate gradient methods for vector computers

In this research, we are concerned with the solution on vector computers of linear systems of equations, Ax = b, where A is a larger, sparse symmetric positive definite matrix. We solve the system using an iterative method, the incomplete Cholesky conjugate gradient method (ICCG). We apply a multi-color strategy to obtain p-color matrices for which a block-oriented ICCG method is implemented on the CYBER 205. (A p-colored matrix is a matrix which can be partitioned into a pXp block matrix where the diagonal blocks are diagonal matrices). This algorithm, which is based on a no-fill strategy, achieves O(N/p) length vector operations in both the decomposition of A and in the forward and back solves necessary at each iteration of the method. We discuss the natural ordering of the unknowns as an ordering that minimizes the number of diagonals in the matrix and define multi-color orderings in terms of disjoint sets of the unknowns. We give necessary and sufficient conditions to determine which multi-color orderings of the unknowns correpond to p-color matrices. A performance model is given which is used both to predict execution time for ICCG methods and also to compare an ICCG method to conjugate gradient without preconditioning or another ICCG method. Results are given from runs on the CYBER 205 at NASA's Langley Research Center for four model problems.

Poole, E. L.↗

Krylov methods preconditioned with incompletely factored matrices on the CM-2

The performance is measured of the components of the key interative kernel of a preconditioned Krylov space interative linear system solver. In some sense, these numbers can be regarded as best case timings for these kernels. Sweeps were timed over meshes, sparse triangular solves, and inner products on a large 3-D model problem over a cube shaped domain discretized with a seven point template. The performance of the CM-2 is highly dependent on the use of very specialized programs. These programs mapped a regular problem domain onto the processor topology in a careful manner and used the optimized local NEWS communications network. The rather dramatic deterioration in performance was documented when these ideal conditions no longer apply. A synthetic workload generator was developed to produce and solve a parameterized family of increasingly irregular problems.

Berryman, Harry↗

A simultaneous coordinate relaxation algorithm for large, sparse matrix eigenvalue problems

An algorithm is proposed for a scheme of simultaneous coordinate relaxation. A variant of root-shifting coordinate relaxation, this procedure consists of iterating several vectors at the same time, instead of one at a time. Results of application of the algorithm to test matrices are discussed. For many matrix eigenvalue problems for which coordinate relaxation is a viable approach, the present algorithm is more effective than previous implementations of coordinate relaxation. Total central processor operations should be decreased due to significantly improved convergence.

Raffenetti, R. C.↗

Probing for the Trace Estimation of a Permuted Matrix Inverse Corresponding to a Lattice Displacement

We report thatpProbing is a general technique that is used to reduce the variance of the Hutchinson stochastic estimator for the trace of the inverse of a large, sparse matrix A. The variance of the estimator is the sum of the squares of the off-diagonal elements of A -1 . Therefore, this technique computes probing vectors that when used in the estimator annihilate the largest off-diagonal elements. For matrices that display decay of the magnitude of |A$^{-1}_{ij}$| with the graph distance between nodes i and j, this is achieved through graph coloring of increasing powers A k . Equivalently, when a matrix stems from a lattice discretization, it is computationally beneficial to find a distance-k coloring of the lattice. Previously, a hierarchical coloring was proposed so that k can be increased at runtime as needed without discarding previous work. In this work, we study probing for the more general problem of computing the trace of a permutation of A -1 , say PA -1 . The motivation comes from lattice quantum chromodynamics (QCD), where we need to construct “disconnected diagrams” to extract flavor-separated generalized parton functions. In lattice QCD, where the matrix has a four-dimensional toroidal lattice structure, these nonlocal operators correspond to a PA -1 , where P is the permutation relating to some displacement $\vec{p}$ in one or more dimensions. We focus on a single dimension displacement (p), but our methods are general. We show that probing on A k or (PA) k does not annihilate the largest magnitude elements. To resolve this issue, our displacement-based probing works on PA k using a new coloring scheme that works directly on appropriately displaced neighborhoods on the lattice. We prove lower bounds on the number of colors needed and study the effect of this scheme on variance reduction, both theoretically and experimentally on a real-world lattice QCD calculation. We achieve orders of magnitude speedup over the unprobed or the naively probed methods.

97 MATHEMATICS AND COMPUTING↗

Butterfly Factorization Via Randomized Matrix-Vector Multiplications

This paper presents an adaptive randomized algorithm for computing the butterfly factorization of an m × n matrix with m ≈ n provided that both the matrix and its transpose can be rapidly applied to arbitrary vectors. The resulting factorization is composed of O(log n) sparse factors, each containing O(n) nonzero entries. The factorization can be attained using O(n 3/2 log n) computation and O(n log n) memory resources. Furthermore, the proposed algorithm can be implemented in parallel and can apply to matrices with strong or weak admissibility conditions arising from surface integral equation solvers as well as multi-frontal-based finite-difference, finite-element, or finite-volume solvers. A distributed-memory parallel implementation of the algorithm demonstrates excellent scaling behavior.

97 MATHEMATICS AND COMPUTING↗