Stack Exchange Network

Stack Exchange network consists of 175 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.

1
vote
0answers
142 views

Numerical Precision in Matrix Inversion Routines

Let's say I have a matrix $\mathbf{A}\in\mathbb{R}^{n\times n}$, and its associated inverse, $\mathbf{A}^{-1}$. The elements of $\mathbf{A}$ are given in IEEE single precision, i.e. 23-bit mantissa, ...
2
votes
5answers
527 views

c++ Libraries for large linear system of equations

I am looking for libraries for solving large scale linear system (10e5) of equations using parallelization and shared memory. 1. Sparse, complex symmetric, SPD. 2. Suitable for higher order FEM ,...
4
votes
1answer
402 views

Poor SVD reconstruction of singular matrix

I am trying to calculate the singular value decomposition of this matrix using numpy.linalg.svd . However, reconstructing the matrix from the SVD gives a poor ...
1
vote
0answers
151 views

How to make LAPACK eigenvectors orthogonal like Matlab?

I'm using LAPACK zgeev to calculate eigenvectors of a symmetric complex matrix of high dimensions ($n \approx 2000$). I need these eigenvectors to satisfy $$\sum_{...
3
votes
1answer
1k views

Computing the Cholesky decomposition based of the QR decomposition

Let A be a n×n positive-definite Hermitian matrix. I already have the QR decomposition of A. Is there an efficient way to utilize this knowledge to speed up the Cholesky decomposition of A?
1
vote
2answers
106 views

Numerically stable computation of $F((ax+b)^k)$

Let $F: \mathbb{R} \rightarrow \mathbb{R}$ be a linear map. I want to evaluate an expression of the type $$F((ax+b)^k)$$ in terms of $F(x)$ for some fixed value of $x$ (I already know $F(x)^r$ for $r=...
3
votes
2answers
316 views

Difference in performance of preconditioned GMRES and MINRES

I have two matrices $A, B$ coming from a finite element discretization of a system of partial differential equations. $A$ represents the system matrix and is symmetric and indefinite. $B$ is symmetric ...
13
votes
4answers
253 views

Rapidly determining whether or not a dense matrix is of low rank

In a software project that I'm working on, certain computations are vastly easier for dense low-rank matrices. Some problem instances involve dense low-rank matrices, but they're given to me in full, ...
0
votes
1answer
51 views

Maximize sum(AI) for matrix A and any permutation of identity matrix I

I have a random binary matrix $A$ $$ A=\left[\begin{array}{c c c c c}0&0&0&0&1\\0&1&0&1&0\\1&1&1&1&0\\0&1&1&0&0\\1&0&1&1&...
1
vote
1answer
153 views

How can I solve on a computer a large projection problem with redundant constraints?

This question is the essence of this one. After we remove all the cruft, we can recast it as follows: Problem: Given $b \in \mathbb{R}^n$, $C\in \mathbb{R}^{n\times m}$, and $g\in \mathrm{Range}(C^...
2
votes
0answers
55 views

How can I numerically solve a saddle point problem with repeated constraints?

I am interested in numerically solving the following constrained minimization problem; Find the value of $x\in \mathbb{R}^n$ that minimizes $f$ where $f\colon \mathbb{R}^n\to \mathbb{R}$ is defined ...
1
vote
1answer
365 views

Fast way to compute the diagonal elements of the inverse of covariance matrix

We have a variance-covariance matrix denoted with $X^TX$, where $X$ is the design matrix. In linear regression we can estimate beta coefficients with normal equations like $\hat{\beta} = (X^TX)^{-1}X^...
1
vote
0answers
154 views

Best algorithm for inversion of matrix spanning many orders of magnitude [duplicate]

I have a very similar problem to the one described in Calculate inverse of dense matrix with entries of very different magnitude. The reason why I am opening a new question is because as far as I ...
1
vote
1answer
68 views

What is the difference in the pivoting strategies between LAPACK's dpstrf and dpst2 and why?

dpstf and dpstrf sometimes give different pivot results. Of course I can read the source code, but I don't get the idea from it. Since pivoting is for stability of the Cholesky decomposition, one of ...
1
vote
1answer
490 views

Regularization vs constrained optimization of an ill posed tomography problem

I am trying to solve an ill-posed linear system of equations. The particular system has 160 equations and 400 variables. Moreover, the condition number of the left hand side matrix is of order $10^{16}...
0
votes
1answer
71 views

Invert a matrix only on a subset of variables / Compute the “equivalent circuit”

Suppose I have a complicatedly shaped conductor with inhomogeneous conductivity. The conductor is modeled using the Finite Element Method. The conductor has electrical contacts on both ends. Those ...
4
votes
0answers
102 views

How to stochastically estimate the trace of a matrix?

Specifically, the diagonal elements (can possibly both positive and negative) of the matrix can be computed efficiently but the total number is large ($\mathcal O(10^{18})$). My first thought about ...
3
votes
2answers
442 views

Quadratic programs with rank deficient positive semidefinite matrices

Let $A$ be a $n\times n$ square symmetric matrix. In addition, $A\succeq0$ and $\mathrm{rank}(A)<n$. This means that all eigenvalues are non-negative, but also that there are some zero eigenvalues. ...
0
votes
1answer
453 views

Set of linear ordinary differential equations with a mass matrix

What methods are known for efficiently solving a large set of linear homogeneous ordinary differential equations of the following form? \begin{equation} \mathbf{B} \frac{d\mathbf{y}}{dt} = \mathbf{A} \...
1
vote
1answer
150 views

Matrix Decomposition of Conics

I was reading about ellipse-ellipse intersection and I came across this article: https://math.stackexchange.com/questions/679622/intersection-between-conic-and-line-in-homogeneous-space/867428#...
1
vote
0answers
45 views

Any method to efficiently compute SVD of a perturbation of matrix $\bf A$ if the SVD of $\bf A$ is already known? [duplicate]

Suppose we know the SVD of matrix $\bf A$, and $\bf B$ is a slight perturbation of $A$ (e.g. $\|{\bf B}-{\bf A}\|_{\text F}$ is relatively small), then is there any method that can efficiently compute ...
0
votes
1answer
94 views

simple matrix multiplication?

can any one show the steps how I convert the summation to norm as below: $$ \sum_{j=1}^{M'} \big(\lambda_j (\mathbf{L}\mathbf{x})^2_j\big) = \|\mathbf{L}_\lambda \mathbf{x}\|^2, \qquad\text{where } \...
1
vote
0answers
45 views

How is the Gastner-Newman equation implemented to create value-by-area cartograms?

There is a paper called "Density-equalizing map projections: Diffusion-based algorithm and applications" by Michael T. Gastner and M. E. J. Newman, which explains their algorithm (which is based in ...
1
vote
1answer
73 views

Efficient computation of $BX=A$ when LU factorization of $A$ is given

First, $AX=B$ is solved, so I have the LU factorization of $A$ computed already. Now I need to solve $BX=A$. Is there any way to reuse this information (LAPACK ...
1
vote
1answer
1k views

How Jacobian matrix helps optimization faster?

I tried some python optimization functions and some of them needed Jacobian matrix prior for faster convergence. I understand Jacobians are basically transformation matrices that data from one space ...
2
votes
2answers
120 views

Efficient computation of AX=B where B has special structure (block-diagonal)

In case B(size ~ 2k, complex double) is block-diagonal, where block size is small(e.g. 2), is there any more efficient way to compute this other than Lapack gesv?
10
votes
2answers
392 views

Eigenvectors of a small norm adjustment

I have a dataset that is slowly changing, and I need to keep track of eigenvectors/eigenvalues of its covariance matrix. I've been using scipy.linalg.eigh, but it'...
1
vote
1answer
232 views

Why the product of symmetric-sparse matrices is not symmetric, or dense

For the equation $\frac{\partial u}{\partial t}=-(-\Delta)u$ with zero boundary condition using finite element method. Applying the corresponding weak formulation and taking $v=\phi_{j},j=1,2,...,N$, ...
1
vote
1answer
162 views

Fast counting of all submatrices of a binary matrix with a full column rank

I have a binary full-rank matrix of size, say, $25 \times 50$. I need to count how many subsets of its columns form matrices with a full column rank, i.e. all the columns in the subset are linearly ...
3
votes
2answers
471 views

Solving Ax = b with sparse A and sparse b

Let's suppose I'm numerically solving the Poisson equation for a delta function source: $$ \nabla^2 f(x) = \delta(x-x') $$ I can represent the Laplacian $\nabla^2$ using the finite difference method ...
2
votes
1answer
376 views

Distributed (MPI) matrix matrix multiplication

I perform matrix matrix multiplications (between rank-3 and rank-2 arrays) in fortran using following subroutine, ...
3
votes
1answer
2k views

Block-matrix SVD and rank bounds

Assume, we have an $m\times n$ block matrix $M=\left[\begin{array}{c c}A&C\\B&D \end{array}\right]$, where $A$ is an $m_1 \times n_1$ matrix of rank $k_A$. $B$ is an $m_2 \times n_1$ matrix ...
1
vote
0answers
123 views

Eigenvalues using QR iteration

I'm trying find the eigenvalues of a matrix A using QR iteration with Householder. I used this code which I found from Cornell University that decomposes QR with Householder. ...
2
votes
2answers
70 views

Using low rank property for maximal/minimal value search (or sorting)

I was thinking about the following problem: Suppose there is a positive semidefinite matrix $X$ of size $n$ (for example, a kernel). Suppose $X$ can be approximated as a low rank matrix, $X\approx ...
1
vote
0answers
56 views

Reduce large sparse linear operators to memory efficient loops?

I'm dealing a lot with large sparse linear operators these days and I'm quite new to them. A lot of the matrices I deal with originate with only a few unique integers, however, there are lots of them. ...
6
votes
1answer
200 views

Does mean removal increase accuracy of numerical differentiation?

I wish to compute the derivative of a vector through numerical differentiation. Let's say, we use a standard 2nd order central difference scheme, to arrive at a differentiation matrix, and apply it on ...
2
votes
2answers
349 views

Precision loss in Matrix-Vector product when applying Finite-Difference scheme

I am applying a 6th order Finite-Difference differentiation scheme as seen in http://www.scholarpedia.org/article/Method_of_lines/example_implementation/dss006 There seems to be severe numerical/...
3
votes
1answer
216 views

SLEPc eigensolvers take long time to converge for large sparse symmetric matrices

I am leveraging the SLEPc library for solving the first $k$ (where $k = 3$ or $4$) eigenvalues and their corresponding vectors for a matrix of size 200,000. The matrix is sparse and symmetric. I ...
1
vote
1answer
251 views

Preconditioning ARPACK eigenvalue solver

I am working on a generalized eigenvalue problem of the form $$ \boldsymbol{A}\cdot\boldsymbol{x}=\lambda\boldsymbol{B}\cdot\boldsymbol{x} $$ where $\boldsymbol{B}$ is not symmetric positive. ...
0
votes
1answer
96 views

Writing a non-square linear system in standard form $A\cdot{x}=b$

I have spend the last few days working my way through an interesting paper and I'm building a numerical model so I can apply the method. However, I am getting stuck at an "it can be shown" step. I am ...
6
votes
3answers
2k views

Method to check for positive definite matrices

I think it's already been asked, but I still can't figure out a way to do it computationally. I had to check for positive definiteness of an $n \times n$ matrix $A$. I know that for any nonzero ...
5
votes
2answers
145 views

Compute $x = B^{-1}(2A+I)(C^{-1}+A)b$ without calculating matrix inverses

If $A, B, C$ are $n \times n$ matrices, where both $B$ and $C$ are nonsingular, and $b$ is a vector of length $n$, how would you compute the following without computing any inverses? $$x = B^{-1}(2A+...
1
vote
1answer
97 views

Eigenvalue problem constrained with a penalty method

I am trying to constrain an eigenvalue problem. I am aware of the method utilizing the nullspace of the constraint vectors but I was wondering if it would be to use a penalty method for the same ...
2
votes
2answers
2k views

BLAS libraries for Octave or Matlab, preferrably with GPU support?

I just searched around a bit for BLAS implementations and was amazed by the sheer amount of libraries around. Does someone know of a benchmark or otherwise rating of the various libraries? How easy ...
2
votes
0answers
194 views

Connectivity and Clustering using Eigenvectors and the Fiedler Vector

Going off of the answer here: sorting adjacency matrix by the Fiedler vector So here, Jesse the answerer plotted the first 3 eigenvectors associated with non-zero eigenvalues of the Laplacian against ...
1
vote
2answers
223 views

Comparison between two matrices

I have a large dense matrix $A$. For simplicity in simulation purposes and application demand, I induced sparsity by replacing lower/insignificant values by zero and reordering it in block diagonal ...
1
vote
0answers
165 views

All eigenpairs of large sparse symmetric matrix

In advance I am sorry for my noobish question. I am a physics PHD student and basically I use python for my math/physics problems. But now I have a problem which requires more computing capacity and ...
7
votes
1answer
252 views

Numerically stable computation of the Characteristic Polynomial of a matrix for Cayley-Hamilton Theorem

I was wondering if there is any known way to compute the Charactaristic Polynomial P of a matrix A numerically stable in the sense of realizing the Cayley-Hamilton Theorem, i.e. that P(A)=0. I've ...
1
vote
2answers
467 views

Derivative of the inverse of the Right Cauchy-Green Deformation Tensor wrt itself

In continuum mechanics, we define the Right-Cauchy-Green Deformation Tensor as $\boldsymbol{C}=\boldsymbol{F}^T\boldsymbol{F}$ I want to compute $\frac{\partial \boldsymbol{C}^{-1}}{\partial \...
0
votes
3answers
1k views

What's the fastest implementation of elementwise vector multiplication in Fortran?

My fortran code contains lines like the following ...