Stack Exchange Network

Stack Exchange network consists of 174 Q&A communities including Stack Overflow, the largest, most trusted online community for developers to learn, share their knowledge, and build their careers.

Visit Stack Exchange

Questions tagged [linear-algebra]

Questions on the algorithmic/computational aspects of linear algebra, including the solution of linear systems, least squares problems, eigenproblems, and other such matters.

0
votes
1answer
70 views

Efficiently approximating sum of 2-norms

Suppose I have a real vector $\bf x$ of length $KN$, where $K<<N$. Let's say I break this vector $\bf x$ up into $N$ vectors each of length $K$, $${\bf x}_1,...{\bf x}_n.$$ Now I would like ...
1
vote
0answers
57 views

Efficiently computing the properties of a Chebyshev series

Suppose we have some function $f(x)$ defined as a Chebyshev series up to order $M$: $$ f(x) = \sum_{n=0}^{M} c_n T_n(x). $$ For a given coefficients vector $\mathbf{c}$, and $x \in [-1,1]$ I'm ...
1
vote
0answers
79 views

How to get the eigenvalues of Hamiltonian in an over complete basis

Let $|\psi_i\rangle$, $i=1...N+m$, be a set of overcomplete basis vector in a $N$-dim Hilbert space. The following are known: (Einstein's summation convention assumed) $$\hat{H}|\psi_i\rangle=H_{ji}|\...
1
vote
1answer
552 views

Eigen - Solving Complex Generalized Eigenvalue Problem

I've been using the Eigen C++ linear algebra library to solve various eigenvalue problems with complex matrices. I've recently had to use a generalized eigenvalue solution process, only to be ...
1
vote
0answers
107 views

Updating factorization of Laplacian (add/remove edges)

For a graph $G=(V,E)$, recall that the unweighted Laplacian is $L:=D^\top D$, where $D\in\{-1,0,1\}^{|E|\times|V|}$ is the graph "gradient" operator that subtracts adjacent vertex values onto edges. ...
2
votes
1answer
84 views

Solving system of related equations without completely recomputing LU decomposition for each equation

Let $\Sigma$ be positive definite and the $D_i$ positive diagonal. Let the $X_i$ be unknown square matrices. Consider the system of equations: $$(I+\Sigma D_i)X_i=\Sigma\hspace{5mm}\text{for}\...
8
votes
1answer
185 views

Solving two inverse problems with same solution

I've got two inverse problems, $$A_1 ~ x = b_1 \qquad A_2 ~ x = b_2$$ So far I've been solving them independently using Tikhonov Regularization and getting two estimates for $x$. However in my case $...
2
votes
1answer
282 views

Rectangular Cholesky Decomposition

Assume we have a positive semidefinite matrix $A$ with $A\in\mathbb{R^{n\times n}}$. Clearly a cholesky decomposition $A=B^TB$ exist with $B\in\mathbb{R^{n\times n}}$. For my research it would be ...
0
votes
1answer
161 views

Open Source Linear Algebra Library

I am making a code in C that requires the equivalent of Matlab's '\' command for a linear system of the form AX=B where A is an NxN matrix and X, B are Nx1 vectors- i.e a code that performs X=A\B that ...
3
votes
2answers
412 views

Python environments for AMG and Gauss Seidel as solvers instead of preconditioners

I am working on block preconditioning and seemingly it is common to write customised Krylov solvers for them. Within each solver, the individual block linear system with preconditioners are ...
3
votes
0answers
94 views

Eigenvalue decomposition of the sum: $AA^T$ + diag($u$)

Suppose $A\in\mathbb{R}^{n\times c}$,$u\in\mathbb{R}^n$,$n\gg c$. The time complexity of eigenvalue decomposing directly for matrix $AA^T+\text{diag}(u)$ is $O(n^3)$. And it is easy to avoid $O(n^3)$ ...
1
vote
0answers
26 views

Difference between Chebyshev first and second degree iterative methods

Consider linear equation $Au = f$. We want to solve it with iterative method (assuming $A$ is good). First order iterative method is: $$ u^{k+1} = u^k - \alpha_{k+1}(Au^k - f), $$ The second degree ...
7
votes
2answers
112 views

Lanczos algorithms for Hermitian system with Toeplitz kernel

Basically, I am trying to compute the SVD of a large Hermitian matrix $H$ using Lanczos iteration, while $H$ consists if a Toeplitz kernel $K$, which should be able to help speed up the matrix-vector ...
4
votes
3answers
1k views

Understanding Finite-Element Modal Analysis

I am teaching a basic course on computational physics and for the last part of the course I will introduce freshman physics undergraduates to finite-element modelling methods. I am preparing a COMSOL ...
1
vote
1answer
118 views

How to determine the truncation error with products and quotients

If I have an equation given by $$\displaystyle Y = \frac{a^2}{d^2}\frac{(1-c^2\frac{c}{a})}{(1-b^2)}$$ and I expand $a,b,c,d$ in a Taylor series, where $a$ is truncated at the $A^{th}$ order, $b$ is ...
5
votes
0answers
114 views

Rank-one Update to a Rank Revealing QR (RRQR) Factorization?

Suppose we are given an RRQR factorization for some matrix $A \in \mathbb{R}^{m \times n}$, $A\Pi = QR$ where $m > n$. Is there a cheap way to update $A' = A + uv^{\top}$ given this factorization?...
4
votes
3answers
163 views

Eigenvectors of Black-box matrix

$\DeclareMathOperator{\diag}{diag}$ Consider the generalized eigenproblem $A\mathbf{x}=\lambda B\mathbf{x}$. When solving PDEs numerically (specifically, I am interested on finding the Dirichlet ...
3
votes
2answers
467 views

Zero Eigenvalues in Lanczos Algorithm

I need to find the smallest few eigenvalues of a Hamiltonian (exact diagonalization) I use Python, and SciPy's built-in sparse eigenvalue solver. I notice, however, that for my small system (only a ...
0
votes
1answer
65 views

Show the symmetric Gauss-Seidel converges for any $x_0$

Let $A\in\mathbb{R}^{n\times n}$ is symmetric positive definite and consider solving linear system $Ax = b$. Show that the symmetric Gauss-Seidel iteration converges for any $x_0$. Solution - Since $...
0
votes
1answer
72 views

Linear stationary iteration method

Suppose you are attempting to solve $Ax = b$ using linear stationary iteration method defined by $$x_k = G x_{k-1} + f$$ that is consistent with $Ax = b$, i.e., for which $f = (I - G)A^{-1}b$. Suppose ...
3
votes
1answer
174 views

Efficient algorithm for a matrix product

Recall that a unit lower triangular matrix $L\in\mathbb{R}^{n\times n}$ is a lower triangular matrix with diagonal elements $e_i^{T}L e_i = \lambda_{ii} = 1$. An elementary unit lower triangular ...
3
votes
1answer
200 views

Determine a sufficient condition for a Hessenberg matrix to be nonsingular

Consider $A\in\mathbb{R}^{n\times n}$ whose nonzero elements are restricted to the main diagonal the strict upper triangular part, and the first subdiagonal. For $n = 8$ the locations that must be ...
0
votes
1answer
77 views

$LU$ Factorization of a nonsingular matrix with a particular pattern

Consider $S\in\mathbb{R}^{n\times n}$ whose nonzero elements have the following pattern for $n = 8$: $$\begin{pmatrix} 1 & 0 & 0 & 0 & \mu_1 & 0 & 0 & 0\\ 0 & 1 &...
2
votes
1answer
150 views

Schur complement of a matrix $A$

Let $A\in\mathbb{R}^{n\times n}$ and its inverse be partitioned $$A = \begin{pmatrix} A_{11} & A_{12}\\ A_{21} & A_{22}\\ \end{pmatrix},\:\: A^{-1} = \begin{pmatrix} \tilde{A_{11}} & \...
1
vote
1answer
340 views

Preconditioning of two step iteration for dense matrices

I would like to solve a dense linear system the form in python $$ L\left(\boldsymbol{x}\right):=\left[\gamma^+\left[\boldsymbol{A}+\frac{1}{2}\boldsymbol{B}^{-1}\right] +\gamma^-\left[\boldsymbol{A}-\...
0
votes
1answer
151 views

Parallel linear algebra without OpenMP

I have searched through the archives without success. Apparently, the question is simple: What linear algebra library can I use that is parallel (shared memory) but without OpenMP? As far as I've ...
1
vote
0answers
41 views

Optimal ordering in Jacobi SVD algorithm

In Jacobi SVD algorithm as given here every pair of columns of the matrix is orthogonalized until convergence. I want to know that how does the order of selection of the pair of columns affect the ...
6
votes
1answer
152 views

Is there an eigenvalue estimation method more accurate than Gershgorin's, which uses no multiplication?

Suppose I have a real symmetric matrix. I would like to tell wether it has at least $k$ strictly positive eigenvalues, but using only additions (no multiplications). Is there a method that I could use?...
1
vote
1answer
95 views

Numerical method for solving a system with positive definite blocks

I have a system with below coefficient matrix $$ C = \begin{pmatrix} A & B^T \\ B & D \end{pmatrix},$$ where, $A$ and $D$ are square and positive definite. Furthermore, if $B$ be square, ...
3
votes
2answers
1k views

A numerical GMRES example

I'm having trouble understanding how GMRES works. I've read the part in Saad's book and a few others but still I am confused. Can someone provide me a numerical example to understand it better? Or if ...
5
votes
1answer
129 views

Solve for $C$ such that $C^{T}AC$ is banded of given width

Given a symmetric matrix $A$, the Lanczos algorithm outputs $C$ such that $C^{T}AC$ is tridiagonal. Is there a generalization of this such that $C^{T}AC$ is banded of specific width $w$? Note that $C$...
5
votes
0answers
266 views

Compute sparsity pattern of $A^2$

Suppose we have a sparse matrix $A$. Is there any way to compute just the sparsity pattern of $A^2 = A \cdot A$ (I do not actually need to know what exactly the nonzero value are) faster than to ...
6
votes
1answer
231 views

Robust smoothers for geometric multigrid

I'm searching for robust smoothers for geometric multigrids. By robust I mean: Effective for high order approximations (say spectral element, spectral Discontinuous Galerkin), Parallel (suitable for ...
2
votes
0answers
257 views

Fixing a near singular covariance matrix

Given a near singular covariance matrix, the standard method of 'fixing' it seems to be to add a small damping coefficient $c>0$ to the diagonal, which serves to bump all the eigenvalues up by this ...
0
votes
2answers
42 views

Grouping Boolean vectors by similarity up to a rotation [duplicate]

I'm looking for an algorithm to take the entire vector space of length n Boolean vectors and partition it into vectors that are the same up to a rotation of the entries. For example if n=3 the ...
6
votes
2answers
250 views

Solving $A=B+AB$ without matrix inverse

I have a linear system of equations that can be expressed $$A=B+AB,$$ where $A$ and $B$ are real, symmetric matrices. I would like to solve for $A$ given $B$. At present, I solve for $A$ directly via $...
0
votes
1answer
74 views

Most efficient way to compute eigenvectors / values of this matrix?

I have a symmetric $ 3 \times 3 $ matrix $A$ and I need to compute the eigenvectors and eigenvalues of this. I know that I can use something like Lapack, but I also know that this can be computed ...
2
votes
1answer
190 views

Probability of reconstructing a word using c substrings from a random sample

Consider a voice recording split into it's phonemes as our sample $S=(s_1,...,s_k) \in \Omega = P^k$. The number of phonemes is $|P| = 40$. Then I have a word $w = (w_1,...,w_n) \in P^n$. I want to ...
1
vote
1answer
273 views

Solving a linear system whose matrix has imbalanced diagonal entries

I am trying to solve following set of equations: A(i,i-2)*u(i-2) + A(i,i-1)*u(i-1) + (A(i,i)+β(i) )*u(i) + A(i,i+1)*u(i+1) + A(i,i+2)*u(i+2)= B(i) + β(i) where i=1:1000000 If values of β varies(...
3
votes
0answers
223 views

Quasi Newton method for block diagonal Hessian

I am having an unconstrained optimization problem, where the Hessian is positive semidefinite and block diagonal. The function is strictly convex, hence, the curvature condition ($ s^{T}_{k}y_{k} > ...
1
vote
1answer
202 views

Efficiently rotate vector in 2D (and 3D)

I need to efficiently rotate a 2D (and 3D) vector in a CUDA kernel. I was thinking about generating random unitary rotation matrices. I don't need to know the angle, it just has to be randomly ...
5
votes
0answers
137 views

Galerkin FEM error when using even number of elements

Intro: I am developing a Galerkin FEM code in matlab - starting small with a simple 1D ODE. The equation I'm trying to solve is: $$a u_x = cos(x), x \in [0, 2\pi]$$ Which has a known exact solution $$...
2
votes
0answers
49 views

SPECT reconstrction using MLEM

In Single-Photon Emission Computerized Tomography (SPECT) parallel beam reconstruction using Maximum-Likelihood Expectation–Maximization(MLEM), is it sufficient to scan the object around 180 degree? ...
1
vote
0answers
27 views

Ordering of eigenvectors to maximise trace of diagonalising matrix

I asked a similar question on the Mathematics stack exchange here without much success, so I thought I'd ask it with a more practical bent here. Suppose we have a Hermitian matrix $H$ with (for the ...
1
vote
0answers
172 views

Preconditioned Steepest Descent

For my program assignment I need to write a preconditioned steepest descent algorithm. I have psedo-code from my professor which is here: From reading in my text on this method they say that $P$ is ...
2
votes
1answer
91 views

Condition Number of Rectangular Matrices

The 2-norm condition number can be easily extended to rectangular matrices. I'm wondering if the inequality for the product of matrices still holds in that case, i.e., $\operatorname{cond}(AB) \leq \...
3
votes
1answer
221 views

Methods for solving rectangular, full-rank systems of equations — which is best?

Suppose I have a large, sparse, $m \times n$ matrix $A$, with $m \gt n$ and $\text{rank}(A) = n$. I wish to solve $Ax=b$. Suppose I know that $A$ has the following characteristics: $A$ is somewhat ...
4
votes
1answer
157 views

Algorithm to decompose a sparse unitary matrix into a Kronecker product of smaller unitary matricies

Given some sparse unitary square matrix $A$ ($dim=2^n$ if it matters), is there an algorithm to decompose $A$ into a Kronecker/tensor product of smaller unitary matrices? In other words: decompose ...
4
votes
1answer
491 views

Power Iteration over Rayleigh Quotient Iteration?

It is a commonly known fact that the Rayleigh Quotient converges cubically (1), while the Power Iteration may converge slowly if the difference between the dominant and second-dominant eigenvalue is ...
8
votes
1answer
747 views

What is the state of the art algorithm for diagonalizing real symmetric matrices?

There are many methods for diagonalizing matrices, probably the most widely used is the combination of Householder transformations and the QR algorithm. Is there any superior method for diagonalizing ...