Questions tagged [matrix-factorization]

Decomposition of a matrix into a product of matrices with special properties. Common matrix factorizations include LU, QR, SVD, and Cholesky.

Filter by
Sorted by
Tagged with
0 votes
0 answers
36 views

Solving Symmetric/Hollow Matrix issue

I have a particular issue and need something creative or solution from calculus. I have Symmetric/Hollow Matrix, a numbers are % of mismatch between them. Ideally, all of them should be 0, but I have ...
user avatar
1 vote
1 answer
73 views

Fast algorithm to compute chi-square

I would like to evaluate the chi-square of the form $\chi^2=v^{T}C^{-1}v$ where $v$ is a column vector and $C$ is a covariance matrix. Both $v$ and $C$ are known and $C$ is a $740\times740$ matrix. ...
user avatar
  • 123
0 votes
1 answer
86 views

Incomplete Cholesky factorization algorithm

I want to implement incomplete Cholesky factorization to precondition, the algorithm I refer from incomplete Cholesky factorization, ...
user avatar
3 votes
1 answer
164 views

Incomplete Cholesky preconditioner for CG efficiency

I am currently solving the harmonic equation using a P1 FEM discretisation. The resulting matrix $A$ is SPD and fairly sparse so I use a preconditioned conjugate gradients (CG) solver to find a ...
user avatar
1 vote
1 answer
246 views

Givens rotation algorithm without matrix-matrix multiplication

I would like to implement a givenRotation algorithm without having matrix-matrix multiplication. Matrix-vector is fine or just for looping. I am to decompose a rectangular (m+1)xm Hessenberg matrix. I ...
user avatar
5 votes
2 answers
479 views

Computational method to compute both the (log) determinant and inverse of a matrix

Suppose I have a square matrix $\mathbf{A} \in \mathbb{R}^{n\times n}$ and a vector $\mathbf{b}\in\mathbb{R}^n$. In my application I need to accomplish two things. I need to find the solution of the ...
user avatar
6 votes
1 answer
198 views

Algorithm for solving systems which are nearly symmetric/adjoint?

I am familiar with Cholesky decomposition and LU factorization for solving systems of linear equations. I have a problem where I have large sparse matrices (say, 1000x1000 or larger) where only one or ...
user avatar
3 votes
0 answers
112 views

Invert a huge sparse operator;

please help me with this question, I want to invert a huge sparse (non-circulant) this below in a $Ax=y$ equation: $$(\lambda I+ \beta D+ \sigma C)x=y$$ where I is an Identity Matrix,D is a Diagonal ...
user avatar
1 vote
0 answers
51 views

Using the lower triangular portion of a matrix to return a Symmetric Positive Definite Matrix?

I've recently posted this question. To summarise, I'm dealing with supposedly Symmetric Positive Definite(SPD) matrices, but due to machine-precision they end up not being SPD. In a comment, a user ...
user avatar
3 votes
2 answers
127 views

Efficient schemes for solving the extended Saddle point problem

I am interested in knowing some efficient techniques for solving the following extended Saddle point problem. \begin{align} \begin{bmatrix} A & B^T & C^T \\ B & 0 & 0 \\ C & ...
user avatar
  • 565
4 votes
2 answers
173 views

Is there a simple way to add a sparse matrix to an LU decomposition of a dense matrix?

I am solving a parabolic equation in the form: $$ \left( {M \over\tau_j} + A \right) u^{j+1} = M f^j + {u^j \over \tau_j}, $$ where $A$ and $f$ are a dense stiffness matrix and the right hand side of ...
user avatar
8 votes
3 answers
734 views

Accurate Way to Calculate Matrix Powers and Matrix Exponential for Sparse Positive Semidefinite Matrices

I do need to numerically calculate the following forms for any $x\in\mathbb{R}^n$, possibly in python: $x^T M^k x$, where $M\in\mathbb{R^{n\times n}}$ is a PSD sparse matrix, $n$ can be quite large ...
user avatar
  • 183
1 vote
0 answers
52 views

Nonsymmetric permutations for LU factorisation of symmetric matrix

Let $A$ be a symmetric matrix. It is then well known that computing the LU factorisation of $PAP^T$ instead of $A$ for a suitably chosen permutation matrix $P$ can greatly reduce fill-in. My question ...
user avatar
  • 405
0 votes
1 answer
83 views

Solve for large array of PD matrices

I have N matrices that are positive definite, and I have to solve for a M vectors. As M is large in my case, doing all solves simultaneously using np.linalg.solve ...
user avatar
  • 109
3 votes
1 answer
199 views

Quick way to find a common basis of eigenvectors between 2 matrices : valid or not?

Following the advise of @Federico Polonion a previous post, one suggested, to find a basis of common eigen vectors between 2 matrices, to simply do : Generate 2 ...
user avatar
2 votes
1 answer
154 views

Efficiently compute a projection matrix from Householders reflectors

Let $A \in \mathbb{R}^{m \times n}$ where $m \geq n$. Let $B$ and $\tau$ be the result of applying LAPACK's dgeqrfp method (R on the upper right triangle, and the ...
user avatar
1 vote
2 answers
137 views

Solving a specific sparse linear system without dense materialization

I need to (computationally) solve a system of equations, for the purposes of an interior point method, of the form $$ \left[\begin{array}{cc}B & A^T \\ A & 0\end{array}\right] \left[\begin{...
user avatar
  • 113
2 votes
1 answer
176 views

Numerical Linear Algebra: When to use Direct methods versus iterative methods to solve a linear system - for PDEs in particular

I am reading the Chapra and Canale book on numerical methods, and was working through the chapters on solving linear systems. Now the book goes through direct methods including Gaussian Elimination, ...
user avatar
  • 255
2 votes
0 answers
71 views

Does incomplete LU preconditioning improve the asymptotic scaling of Krylov subspace methods?

It is well known that unpreconditioned Krylov subspace methods applied to the finite-difference-discretised Poisson equation with $n$ grid points per direction require $O(n \, |\log(\varepsilon)|)$ ...
user avatar
  • 405
1 vote
1 answer
192 views

Range of a matrix from its complete orthogonal decomposition

In this StackOverflow answer, @Gokul has shown how to get a basis of the kernel of a matrix with the help of the 'Eigen' function CompleteOrthogonalDecomposition. ...
user avatar
4 votes
0 answers
61 views

Efficient computation of marginalized multivariate normal likelihood

In general,if we know that the marginal Gaussian distribution for some variable $\textbf{x}$ and a conditional Gaussian distribution for some $\textbf{y}|\textbf{x}$ of the forms: $$p(\textbf{x}) = \...
user avatar
0 votes
2 answers
229 views

"WY" representation of QR factorization -- implementations?

I have a matrix $A \in \mathbb{R}^{m \times n}$ where $m \gg n$ and I want to compute the full QR decomposition $A = QR$. Where $Q$ is an orthogonal $m \times m$ matrix. Bishof & Van Loan (1987) ...
user avatar
  • 145
2 votes
1 answer
198 views

Methods to improve the efficiency and the memory requirement of LU factorization for complex symmetric system matrix

I want to solve a linear set of equations (Ax=b) using LU decomposition. My "A" matrix is a complex matrix which is ...
user avatar
  • 33
2 votes
0 answers
173 views

Computing Singular Value Decomposition of small ($4\times 4$) matrices

I need to compute the Singular Value Decomposition (SVD) of many $4 \times 4$ matrices. I'm looking for SVD algorithms specialized for small matrices. I've read that the ...
user avatar
  • 131
1 vote
0 answers
49 views

Binary data clustering by Matrix factorization [closed]

I have read an article talking about binary clustering using Matrix factorization(see attached), but i would like to understand some optimization concepts: Is it reasonable to use a Frobenius norm in ...
user avatar
1 vote
1 answer
303 views

Implementation of sparse matrix SVD for GPU

I have a sparse matrix $W$ which is almost-squared ($N+1 \times N$) and I would like to know the eigenvalues of $A = W^T W$. $A$ is Hermitian so the eigenvalues are real-positive valued. The usual ...
user avatar
  • 11
1 vote
1 answer
191 views

Efficient way to solve a set of linear equations $Ax=b$ when $A$ is sparse and some elements of $b$ are equal to zero

I have a set of linear equations, $Ax=b$. And about half of the elements in the right-hand side (vector $b$) are equal to zero. My system matrix $A$ is a sparse complex matrix. And $A$ is in the size ...
user avatar
  • 33
1 vote
1 answer
111 views

Functions from Scipy, Blas, or Lapack that compute only upper triangular matrix

My goal is to transform a matrix into upper triangular form in Python. I know the function scipy.linalg.lu will do LU decomposition and get both upper and lower ...
user avatar
  • 207
3 votes
1 answer
650 views

Why is 'scipy.sparse.linalg.spilu' less efficient than 'scipy.linalg.lu' for sparse matrix?

I have a matrix B which is sparse and try to utilize a function scipy.sparse.linalg.spilu specialized for sparse matrix to ...
user avatar
  • 207
5 votes
1 answer
177 views

Algorithm to factorize matrix whose many rows are already of upper triangular form?

I have a matrix whose many rows are already in the upper triangular form. $$\begin{bmatrix} x_{11} & x_{12} & x_{13} & x_{14} & x_{5} \\ 0 & x_{22} & x_{23} & x_{...
user avatar
  • 207
1 vote
0 answers
45 views

Triangle on top of diagonal least squares

I need to solve many least squares problems with the following matrices: $$ \pmatrix{ R \\ D_i } $$ where $R$ is upper triangular and $D_i$ is diagonal. $R$ is the same for all the problems, while $...
user avatar
  • 876
0 votes
2 answers
260 views

Reorder eigenvalues in Schur factorization in descending order

In this command: [US,TS] = ordschur(U,T,select) what should replace the select to rearrange the eigenvalues in descending ...
user avatar
  • 41
1 vote
1 answer
162 views

How can I get Cholesky decomposition from eigenvalue decomposition?

I have $$S = QLQ^T$$ I know $Q$, $L$, $Q^T$. How can I get the $R$ and $R^T$ for the Cholesky decomposition $S=R^TR$?
user avatar
3 votes
2 answers
261 views

Bareiss algorithm vs. LU-decomposition

I at the moment try to fully understand the Bareiss algorithm for calculating determinants. One question that came to my mind is the following: Why is LU-decomposition much more often used than the ...
user avatar
5 votes
0 answers
82 views

Why the two Gram-Schmidt algorithms produce different results for qr factorization?

For the qr factorization using classic Gram-Schmidt algorithm, I found the 2 different implementations below. The first one uses the for loop to compute the upper ...
user avatar
  • 951
1 vote
0 answers
64 views

Numerical methods. MDF (ILU) implementation

I am trying to implement Minimum Discarded Fill (MDF) Ordering algorithm for incomplete matrix factorization. The algorithm description is here on page 60 Preconditioning Techniques for a Newton–...
user avatar
  • 21
6 votes
2 answers
300 views

Block-matrix: optimal fill-in reduction for LU factorization

Consider a square $N \times N$ block-matrix $\mathbf{A}$, where each $n \times n$ block $\mathbf{A}_{ii}$ is either a dense block or a zero-block. So, $N$ denotes the number of blocks, $n$ denotes the ...
user avatar
  • 8,287
2 votes
1 answer
295 views

Reference for QR algorithm for complex matrix

I am trying to find out if the known QR algorithm to find the eigenvalues of a real matrix, which can be found in the book Fundamentals of Matrix Computations, can also be used for complex matrices ...
user avatar
  • 31
0 votes
0 answers
114 views

Why does the matlab command **chol(A)** slower than **chol(A,'lower')** for a large sparse SPD matrix?

For a SPD matrix A, there exists Cholesky factorization $A=LL^T$ or $A=R^TR$, where L, R are a lower and upper triangular matrix, respectively. Also in matlab, there has a command R = chol(A) which ...
user avatar
  • 951
7 votes
3 answers
430 views

Does a symmetric positive definite matrix also have $\mathbf{A} = \mathbf{L}^T\mathbf{L}$ (where $\mathbf{L}$ is a lower triangular matrix)?

As we know, for a symmetric positive definite (SPD) matrix $\mathbf{A}$, there is a theorem about the Cholesky factorization $\mathbf{A}= \mathbf{L}\mathbf{L}^T$, where $\mathbf{L}$ is a lower ...
user avatar
  • 951
4 votes
1 answer
204 views

Low rank update of QR of inverse

I am in a situation where as part of a sort of inverse power method scheme, I want to very often perform the following step: Apply a symmetric rank one update $uu^\top$ to my inverse matrix $A^{-1}$ ...
user avatar
5 votes
2 answers
310 views

Do most statistical packages and libraries in high-level programming languages rely on LAPACK for their matrix inversion operations?

Possible an open-ended question, but I am wondering if most statistical packages and libraries, for instance, Stata, R, Python's NumPy and MATLAB rely on LAPACK algorithms to perform matrix operations,...
user avatar
1 vote
1 answer
70 views

How to avoid unnecessary checks when inverting this LU decomposition

Background for the question I am currently working on a Matlab code in which the systems of linear equations $Ax_1 = b_1$, $Ax_2 = b_2$, ... have to be solved. As the matrix $A$ is constant during ...
user avatar
3 votes
1 answer
165 views

Weighted QR Implementation

Say I want a QR decomposition of matrix $A$, where orthogonality of $Q$ is with respect to a generic non-degenerate positive-definite bilinear form $\phi$ (in my case, $\phi$ is "defined" by a finite-...
user avatar
  • 379
3 votes
1 answer
1k views

Cholesky for ill-conditioned/singular covariance matrices

Can someone suggest a way to get Cholesky factorization of a singular covariance matrix? I need it to match Cholesky on full-rank matrices, ie coordinate order should be preserved. My attempt below ...
user avatar
1 vote
1 answer
271 views

Incomplete LU decomposition of sparse matrix

I have a sparse matrix stored in CSR format. For this matrix, I would like to get the incomplete LU decomposition. I tried to find algorithms which can utilize the CSR format but I could not find ...
user avatar
1 vote
1 answer
294 views

How to use QZ decomposition for single matrix in Matlab?

Can I use QZ decomposition on a single square matrix in Matlab? Like, [Aa,Q,Z]=qz(A);
user avatar
1 vote
1 answer
165 views

Classical vs. modified Gram-Schmidt

It is often said that modified Gram-Schmidt is more robust with respect to rounding errors than classical Gram-Schmidt, but it is very hard to find a good explanation / example of why this is so. Can ...
user avatar
  • 405
3 votes
1 answer
579 views

Givens rotation vs 2x2 Householder reflection

The usual story of Givens rotations vs Householder reflections is that Householder reflections are better if you want to map a long vector to $e_1$, while Givens is better if you want to map a 2-...
user avatar
  • 405
1 vote
1 answer
149 views

Software for parallel incomplete LU factorisation

I am looking for a software package to compute incomplete LU factorisations in parallel. Further considerations are: The package must allow for arbitrary level-of-fill or threshold-based truncation. ...
user avatar
  • 405