Column Pivoting Research Articles

QR decomposition of dual matrices and its application

Dual number matrix decompositions have played an important role in fields such as kinematics and computer graphics in recent years. In this paper, we present a QR decomposition algorithm for dual number matrices. When dealing with large-scale problems, we present the thin QR decomposition of dual number matrices, along with its algorithm with column pivoting. In numerical experiments, we discuss the suitability of different QR algorithms when confronted with various large-scale dual matrices, providing their respective domains of applicability. Finally, we employ the QR decomposition of dual matrices to compute the DMPGI, attaining results of higher precision.

Machine Learning K-Means Clustering in Interpolative Separable Density Fitting Algorithm: Advancing Accurate and Efficient Cubic-Scaling Density Functional Perturbation Theory Calculations within Plane Waves.

Density functional perturbation theory (DFPT) is a crucial tool for accurately describing lattice dynamics. The adaptively compressed polarizability (ACP) method reduces the computational complexity of DFPT calculations from O(N4) to O(N3) by combining the interpolative separable density fitting (ISDF) algorithm. However, the conventional QR factorization with column pivoting (QRCP) algorithm, used for selecting the interpolation points in ISDF, not only incurs a high cubic-scaling computational cost, O(N3), but also leads to suboptimal convergence. This convergence issue is particularly pronounced when considering the complex interplay between the external potential and atomic displacement in ACP-based DFPT calculations. Here, we present a machine learning K-means clustering algorithm to select the interpolation points in ISDF, which offers a more efficient quadratic-scaling O(N2) alternative to the computationally intensive cubic-scaling O(N3) QRCP algorithm. We implement this efficient K-means-based ISDF algorithm to accelerate plane-wave DFPT calculations in KSSOLV, which is a MATLAB toolbox for performing Kohn-Sham density functional theory calculations within plane waves. We demonstrate that this K-means algorithm not only offers comparable accuracy to QRCP in ISDF but also achieves better convergence for ACP-based DFPT calculations. In particular, K-means can remarkably reduce the computational cost of selecting the interpolation points by nearly 2 orders of magnitude compared to QRCP in ISDF.

Real-time calculation method of transformer winding temperature field based on sparse sensor placement

The Real-time obtain of the transformer winding temperature field is of great significance for understanding the operation status of the transformer. Limited by the spatial sampling rate of the sensor measurement method and the computational efficiency of the Computational Fluid Dynamics (CFD) method, the existing techniques cannot complete the real-time obtain of the transformer winding temperature field well. In this study, an estimator is designed, and used as a bridge to establish the relationship between the sparse measurement data of the winding temperature and the winding temperature field, and the real-time obtain of the winding temperature field is realized. Firstly, the proper orthogonal decomposition (POD) modes are obtained based on the prior information on the winding temperature field. Then, a sparse placement strategy of sensors is designed based on the QR factorization with column pivoting (QR-pivot) method, the measurement matrix is formed, and the sparse measurement data of the winding temperature are obtained. Finally, an estimator is designed to complete the winding temperature field's real-time calculation using the winding temperature's sparse measurement data. The transformer temperature rise test proves the effectiveness of the method. This study can provide a new approach for the real-time obtain of transformer winding temperature field.

Numerical Considerations and a new implementation for invariant coordinate selection

Invariant coordinate selection (ICS) is a multivariate data transformation and a dimension reduction method that can be useful in many different contexts. It can be used for outlier detection or cluster identification, and can be seen as an independent component or a non-Gaussian component analysis method. The usual implementation of ICS is based on a joint diagonalization of two scatter matrices, and may be numerically unstable in some ill-conditioned situations. We focus on one-step M-scatter matrices and propose a new implementation of ICS based on a pivoted QR factorization of the centered data set. This factorization avoids the direct computation of the scatter matrices and their inverse and brings numerical stability to the algorithm. Furthermore, the row and column pivoting leads to a rank revealing procedure that allows computation of ICS when the scatter matrices are not full rank. Several artificial and real data sets illustrate the interest of using the new implementation compared to the original one.

Multiple color image cryptosystem based on coupled-logistic-map-biometric keys, QR decomposition with column pivoting and optical Fresnel transform

A novel multiple color image cryptosystem based on coupled-logistic-map-biometric keys (CBK), QR decomposition with column pivoting and optical Fresnel transform (FT), is proposed. In this system, the four sub-bands of each primary color image are produced by the discrete wavelet transform (DWT). The LL sub-bands of four images are fused. The fused image is divided into R, G, and B channels. The two sets of channels (R and B) and (G and B) are independently encoded into complex images, scrambled with generalized Arnold transform (GAT) and then decomposed by QRD with column pivoting. A novel concept of single-channel image constructed by adding of Q1.R2 and Q2.R1 modulated with the R channel of CBK is introduced because the parameters of the coupled logistic map and uniqueness of iris image protect against unauthorized access. The single-channel image is modulated by G channel of CBK and then Fresnel transformed. The phase and amplitude truncation operations are applied to the resultant image to obtain the first cipher image and the first phase key, respectively. The first cipher image is modulated by B channel of CBK and then Fresnel transformed. Finally, the phase and amplitude truncation operations are performed to the transformed image to get the final cipher image and the second phase key, respectively. Two decryption keys produced by QR decomposition with column pivoting and two-phase keys are asymmetric keys to resist potential attacks. Moreover, the parameters of the coupled logistic map (in R, G, and B channels of CBK), GAT and FT supply a large key space of the cryptosystem. A hybrid opto-electronic system can be used to perform the proposed approach. Numerical simulations demonstrate the validity and viability of the proposed scheme.

Primal-Dual Path-Following Methods and the Trust-Region Updating Strategy for Linear Programming with Noisy Data

In this article, we consider the primal-dual path-following method and the trust-region updating strategy for the standard linear programming problem. For the rank-deficient problem with the small noisy data, we also give the preprocessing method based on the QR decomposition with column pivoting. Then, we prove the global convergence of the new method when the initial point is strictly primal-dual feasible. Finally, for some rank-deficient problems with or without the small noisy data from the NETLIB collection, we compare it with other two popular interior-point methods, i.e. the subroutine pathfollow.m and the built-in subroutine linprog.m of the MATLAB environment. Numerical results show that the new method is more robust than the other two methods for the rank-deficient problem with the small noise data.

Identification of Independent Patterns of COVID-19 Mortality in Brazil by a Functional QR Decomposition

The subset selection problem of linear algebra is applied to identify independent patterns of COVID-19 evolution within Brazil. The data consist of a set of mortality curves in states of Brazil. A subset of the most independent curves is selected by using a functional version of the QR matrix decomposition technique with column pivoting. The selected subset is used next as a basis to represent the remaining curves filtering out any data redundancy. For each independent curve, an associated epidemiological region of influence is defined. The results show two main independent curves with a similar two-peak pattern and a 50-day shift between the patterns. Two main epidemiological regions are next identified: one encompassing most of the country from the center and northeast states to the south, an another one containing the Amazonian region at the northwest.

Education

The Education section in this issue presents two contributions. The first paper, “$LU$ and $CR$ Elimination,” is written by Gilbert Strang and Cleve Moler. This is an expository paper which is helpful to instructors and students interested in linear algebra and matrix manipulations. The focus of the paper is the column-row factorization for any matrix $A$ of rank $r$. The matrix is represented as $A=CR$, where the matrix $C$ contains the first $r$ independent columns of $A$, and the matrix $R$ contains the nonzero rows of the reduced row echelon form of $A$. The authors explain how to obtain the factorization and contrast it with $LU$, $QR$ factorizations, the single-value decomposition, and another factorization, $A= CW^{-1}B$, which they call magic. In section 6, the authors point to an important timely application pertaining to randomized algorithms for a low-rank approximation of large sparse matrices. The paper is clearly and nicely written. It includes a discussion about relevant MATLAB functions, examples, and several references to randomized linear algebra literature. The second paper is “Newton's Method in Mixed Precision,” presented by C. T. Kelley. The paper discusses the impact of precision in calculations on the rate of convergence of Newton's method for solving a system of nonlinear equations $F(x)=0.$ The vector-function $F$ is differentiable and its Jacobian $F' (x)$ is nonsingular. It is assumed that the step $s$ in each iteration is determined by solving the linear equation $F' (x)s= - F(x)$ via Gaussian elimination with column pivoting. The paper starts with the classical result on the local convergence of the method under standard assumptions and then passes onto analyzing the effect of errors arising in calculation of the function and in approximation of the Jacobian. The author's derivations predict no significant difference in the convergence of the nonlinear iteration between a double-precision analytic Jacobian with double precision in the linear solver and a forward-difference approximate Jacobian with a single precision in the linear solver. An interesting part of the paper deals with the question of estimating the effect of the backward error in the solver. The author observes that the worst case estimates are too pessimistic and rarely seen in practice. Invoking recent techniques in probabilistic rounding analysis, he derives more realistic results. The theoretical statements are supplemented by a numerical illustration in section 3. An additional example demonstrates how the theory breaks down if the Jacobian is singular at the solution.

Realizing Effective Cubic-Scaling Coulomb Hole Plus Screened Exchange Approximation in Periodic Systems via Interpolative Separable Density Fitting with a Plane-Wave Basis Set.

The GW approximation is an effective way to accurately describe the single-electron excitations of molecules and the quasiparticle energies of solids. However, a perceived drawback of the GW calculations is their high computational cost and large memory usage, which limit their applications to large systems. Herein, we demonstrate an accurate and effective low-rank approximation to accelerate non-self-consistent GW (G0W0) calculations under the static Coulomb hole plus screened exchange (COHSEX) approximation for periodic systems. Our approach is to adopt the interpolative separable density fitting (ISDF) decomposition and Cauchy's integral to construct low-rank representations of the dielectric matrix ϵ and self-energy matrix Σ. This approach reduces the number of floating point operations from to and requires a much smaller memory footprint. Two methods are used to select the interpolation points in ISDF, including the standard QR factorization with column pivoting (QRCP) procedure and the machine learning K-means clustering (K-means) algorithm. We demonstrate that these two methods can yield similar accuracy for both molecules and solids at much lower computational cost. In particular, K-means clustering can significantly reduce the computational cost of selecting the interpolation points by an order of magnitude compared to QRCP, resulting in an overall speedup factor of about ten times ISDF accelerated the static COHSEX calculations compared with conventional COHSEX approximation.

Machine Learning K-Means Clustering Algorithm for Interpolative Separable Density Fitting to Accelerate Hybrid Functional Calculations with Numerical Atomic Orbitals.

The interpolative separable density fitting (ISDF) is an efficient and accurate low-rank decomposition method to reduce the high computational cost and memory usage of the Hartree-Fock exchange (HFX) calculations with numerical atomic orbitals (NAOs). In this work, we present a machine learning K-means clustering algorithm to select the interpolation points in ISDF, which offers a much cheaper alternative to the expensive QR factorization with column pivoting (QRCP) procedure. We implement this K-means-based ISDF decomposition to accelerate hybrid functional calculations with NAOs in the HONPAS package. We demonstrate that this method can yield a similar accuracy for both molecules and solids at a much lower computational cost. In particular, K-means can remarkably reduce the computational cost of selecting the interpolation points by nearly two orders of magnitude compared to QRCP, resulting in a speedup of ∼10 times for ISDF-based HFX calculations.

Sparse Tucker Tensor Decomposition on a Hybrid FPGA–CPU Platform

Recommendation systems, social network analysis, medical imaging, and data mining often involve processing sparse high-dimensional data. Such high-dimensional data are naturally represented as tensors, and they cannot be efficiently processed by conventional matrix or vector computations. Sparse Tucker decomposition is an important algorithm for compressing and analyzing these sparse high-dimensional data sets. When energy efficiency and data privacy are major concerns, hardware accelerators on resource-constraint platforms become crucial for the deployment of tensor algorithms. In this work, we propose a hybrid computing framework containing CPU and FPGA to accelerate sparse Tucker factorization. This algorithm has three main modules: tensor-times-matrix (TTM), Kronecker products, and QR decomposition with column pivoting (QRP). In addition, we accelerate the former two modules on a Xilinx FPGA and the latter one on a CPU. Our hybrid platform achieves $23.6 \times \sim 1091\times$ speedup and over $93.519\% \sim 99.514 \%$ energy savings compared with CPU on the synthetic and real-world datasets.

A hybrid symmetry–PSO approach to finding the self-equilibrium configurations of prestressable pin-jointed assemblies

For pin-jointed assemblies with many members or self-stress states, the form-finding problem using conventional methods generally involves considerable computational complexities due to the large size of the solution spaces. Here, we propose an improved form-finding method for prestressable pin-jointed structures by combining symmetry-based qualitative analysis with particle swarm optimization. Expressed in the symmetry-adapted coordinate system, the nodal coordinate vectors of a structure with specific symmetry and topology are independently extracted from the key blocks of the small-sized force density matrices associated with rigid-body translations. Then, the first block of the equilibrium matrix is computed, in which the null space reveals integral self-stress states. Particle swarm optimization is introduced and adapted to find feasible prestress modes, where the uniformity and unilaterality conditions of the members are considered. Besides, the QR decomposition with column pivoting is adopted for efficient computations on the null space of these blocks. The QR decompositions of the small-sized blocks of the force density matrix and the equilibrium matrix are performed iteratively, to simultaneously find a stable self-equilibrium configuration and a feasible prestress mode. Representative examples show the presented method is computationally efficient and accurate for the form-finding of symmetric tensegrities and prestressed cable–strut structures.

SIGEST

The SIGEST article in this issue, “Randomized Projection for Rank-Revealing Matrix Factorizations and Low-Rank Approximations,” by Jed A. Duersch and Ming Gu, concerns the use of randomization to reduce the bottleneck in a central numerical linear algebra kernel. The authors consider QR factorization, which is a key component in algorithms for constructing stable bases, for solving least-squares problems, and generally for summarizing or extracting insights from data. Column pivoting is the standard technique for dealing with cases where the target matrix is not of full rank, or where the rank is in doubt. However, for large-scale problems column pivoting requires substantial communication between different levels of computer memory, which severely limits the overall efficiency. The authors therefore propose a strategy for choosing pivots based on randomization rather than exhaustive search. This approach reduces communication complexity while increasing uncertainty in the quality of the factorization in a controlled manner. The article makes use of clear pseudocode listings in order to summarize both the original deterministic strategies and the newly developed randomized approaches; see Algorithms 1 to 6. The authors also explain how a Bayesian framework can be set up to quantify the uncertainty introduced by the use of random projections. Extensive, large-scale tests are conducted using Fortran with OpenMP, an application programming interface that supports shared-memory multiprocessing. The original article appeared in the SIAM Journal on Scientific Computing in 2017, and the main idea has subsequently proved to be effective for a number of related tasks. In creating this highlighted article for the SIGEST section, the authors have revised the introductory material, providing additional context; reorganized the document structure and added new visualizations; updated and extended the computational experiments, while making this section more compact and easier to interpret; and discussed recent developments in the field and provided references to the growing literature.

Application of PCA Analysis and QR Decomposition to Address RFM's Ill-Posedness

Rational function model (<small>RFM</small>) is the most widely used sensor model in the remote sensing community. However, it suffers from ill-posedness, challenging its feasibility. This problem, i.e., ill-posedness, is mainly caused due to highly correlated coefficients of the <small>RFM</small>, which magnifies any small perturbations of observations, such as noise and instrumental error. This paper outlines a novel two-step method, called principal component analysis (<small>PCA</small>)-<small>RFM</small>, based on the integration of <small>PCA</small> and QR decomposition. In the first step, the <small>PCA</small>-<small>RFM</small> reduces the observational perturbations from the design matrix using the <small>PCA</small>. In the next step, the <small>RFM</small>'s coefficients are estimated using a <small>QR</small> decomposition with column pivoting and least square method. According to the results, the <small>PCA</small>-<small>RFM</small> is less sensitive than its rivals to the changes of the ground control point (<small>GCPs</small>) distribution. Geometrically speaking, in addition, <small>PCA</small>-<small>RFM</small> is more accurate than recently established methods even in the presence of the small number of <small>GCPs</small>.

Linear-time CUR approximation of BEM matrices

In this paper we propose linear-time CUR approximation algorithms for admissible matrices obtained from the hierarchical form of Boundary Element matrices. We propose a new approach called geometric sampling to obtain indices of most significant rows and columns using information from the domains where the problem is posed. Our strategy is tailored to Boundary Element Methods (BEM) since it uses directly and explicitly the cluster tree containing information from the problem geometry. Our CUR algorithm has precision comparable with low-rank approximations created with the truncated QR factorization with column pivoting (QRCP) and the Adaptive Cross Approximation (ACA) with full pivoting, which are quadratic-cost methods. When compared to the well-known linear-time algorithm ACA with partial pivoting, we show that our algorithm improves, in general, the convergence error and overcomes some cases where ACA fails. We provide a general relative error bound for CUR approximations created with geometrical sampling. Finally, we evaluate the performance of our algorithms on traditional BEM problems defined over different geometries.

An improved WLS-ENO method for solving hyperbolic conservation laws

Liu and Jiao proposed a WLS-ENO method (Weighted-Least-Squares based Essentially Non-Oscillatory). Comparing to the WENO method (Weighted ENO), the WLS-ENO method is more flexible for solving the hyperbolic conservation laws. Notice that, a WLS problem is generated for each WLS-ENO procedure, although the QR factorization with column pivoting is employed, the computational cost of solving WLS problem is still heavy. So in this paper, we will propose an improved WLS-ENO method for improving the efficiency. The basic idea of new method is that, let the reconstructed polynomial approximate the cell averages of some smooth cells in conservation fashion, and approximate other cell averages in WLS fashion, so that the resulted WLS problem only has one unknown. Because the resulted WLS problem can be solved by using the normal equation directly, the whole computational cost is small, thus the new method has better efficiency.

Interpolative Separable Density Fitting through Centroidal Voronoi Tessellation with Applications to Hybrid Functional Electronic Structure Calculations.

The recently developed interpolative separable density fitting (ISDF) decomposition is a powerful way for compressing the redundant information in the set of orbital pairs and has been used to accelerate quantum chemistry calculations in a number of contexts. The key ingredient of the ISDF decomposition is to select a set of nonuniform grid points, so that the values of the orbital pairs evaluated at such grid points can be used to accurately interpolate those evaluated at all grid points. The set of nonuniform grid points, called the interpolation points, can be automatically selected by a QR factorization with column pivoting (QRCP) procedure. This is the computationally most expensive step in the construction of the ISDF decomposition. In this work, we propose a new approach to find the interpolation points based on the centroidal Voronoi tessellation (CVT) method, which offers a much less expensive alternative to the QRCP procedure when ISDF is used in the context of hybrid functional electronic structure calculations. The CVT method only uses information from the electron density and can be efficiently implemented using a K-Means algorithm. We find that this new method achieves comparable accuracy to the ISDF-QRCP method, at a cost that is negligible in the overall hybrid functional calculations. For instance, for a system containing 1000 silicon atoms simulated using the HSE06 hybrid functional on 2000 computational cores, the cost of the QRCP-based method for finding the interpolation points is 38.1 s, while the CVT procedure only takes 0.7 s. We also find that the ISDF-CVT method enhances the smoothness of the potential energy surface in the context of ab initio molecular dynamics (AIMD) simulations with hybrid functionals.

Analysis and comparison of numerical algorithms for finding the GCD of certain types of polynomials in the Chebyshev basis

This research investigates on the numerical methods for computing the greatest common divisors (GCD) of two polynomials in the orthogonal basis without having to convert to the power series form. Previous implementations were conducted using the Gauss Elimination with partial pivoting (GEPP) and QR Householder algorithms, respectively. This work proceeds to seek for a better approximate solution by comparing the results of the implementations with the QR with column pivoting (QRCP) algorithm. The results reveal that QRCP is as competent as the GEPP algorithm, up to a certain degree, giving a reasonably good approximate solution. It is also found that normalizing the columns of the associated coefficient matrix slightly reduces the condition number of the matrix but has no significant effect on the GCD solutions when applying the GEPP and QR Householder algorithms. However equilibration of the columns by computing its ∞-norm is capable to improve the solution when QRCP is applied. Comparing the three algorithms on some test problems, QR Householder outperforms the rest and is able to give a good approximate solution in the worst case condition when the smallest element of the matrix is 1, the entries ranging up to 15 digits integers.

Randomized QR with Column Pivoting

The dominant contribution to communication complexity in factorizing a matrix using QR with column pivoting is due to column-norm updates that are required to process pivot decisions. We use randomized sampling to approximate this process which dramatically reduces communication in column selection. We also introduce a sample update formula to reduce the cost of sampling trailing matrices. Using our column selection mechanism, we observe results that are comparable in quality to those obtained from the QRCP algorithm, but with performance near unpivoted QR. We also demonstrate strong parallel scalability on shared-memory multiple core systems using an implementation in Fortran with OpenMP. This work immediately extends to produce low-rank truncated approximations of large matrices. We propose a truncated QR factorization with column pivoting that avoids trailing matrix updates which are used in current implementations of level-3 BLAS QR and QRCP. Provided the truncation rank is small, avoiding trailing matrix updates reduces approximation time by nearly half. By using these techniques and employing a variation on Stewart's QLP algorithm, we develop an approximate truncated SVD that runs nearly as fast as truncated QR.

Householder QR Factorization With Randomization for Column Pivoting (HQRRP)

A fundamental problem when adding column pivoting to the Householder QR factorization is that only about half of the computation can be cast in terms of high performing matrix-matrix multiplications, which greatly limits the benefits that can be derived from so-called blocking of algorithms. This paper describes a technique for selecting groups of pivot vectors by means of randomized projections. It is demonstrated that the asymptotic flop count for the proposed method is $2mn^2 - (2/3)n^3$ for an $m\times n$ matrix, identical to that of the best classical unblocked Householder QR factorization algorithm (with or without pivoting). Experiments demonstrate acceleration in speed of close to an order of magnitude relative to the geqp3 function in LAPACK, when executed on a modern CPU with multiple cores. Further, experiments demonstrate that the quality of the randomized pivot selection strategy is roughly the same as that of classical column pivoting. The described algorithm is made available under open sou...