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

learn more… | top users | synonyms

1
vote
0answers
41 views

Block preconditioners and condition number

I am working with block Jacobi like preconditioners which are very cheap for my problem. But I could not find much about the dynamics of basic preconditioners (block Jacobi, Gauss-Seidel, ILU etc). ...
1
vote
2answers
128 views

How to obtain a convergent solution iteratively for a linear system of equations?

I am working on a problem that requires an iterative procedure to solve a linear system of equations, the system of equations in matrix form is: $$\underbrace{\begin{bmatrix} r_{11} & r_{12} ...
2
votes
1answer
61 views

Effect of banded matrix on error

It is necessary for me to solve a Poisson problem with a numerical method on a square domain with two types of triangular mesh: uniform triangular mesh (using uniform distribution nodes on square) and ...
3
votes
2answers
110 views

Well-posedness of a linear elasticity problem and Navier-Cauchy equation

I read a master thesis on a topic I'm interested too. This work concern the solution of the displacement equation of motion for a homogeneous, elastic, isotropic material: $$\rho \ddot{\mathbf{u}} - ...
0
votes
0answers
93 views

Estimating the second largest eigenvalue

I am currently dealing with the following problem. I am given a matrix $A$ of order $n \times n$ where $n \leq 20.$ The principal $(n-1) \times (n-1)$ matrix of $A$ is symmetric and contains only ...
0
votes
1answer
38 views

compute change of phase along closed contour

The following image represents the phase of a wavefunction (in radians) on a square lattice, where $m$ and $n$ label the lattice sites. Computationally speaking, it is the density plot of a 41x41 real ...
1
vote
0answers
44 views

Library for calculating determinants with Kronecker products

I need to calculate a determinant consisting of vectors, using the Kronecker product as product. As an example I would need to be able to calculate: $\left| \begin{array}{cc} ...
2
votes
0answers
73 views

Using SVD to biorthogonalize left and right eigenvectors?

I have a set of left and right eigenvectors from an nonsymmetric eigenproblem, and I'd like to biorthogonalize them. I tried Gram-Schmidt, but this fails for most cases. I then read that the SVD is ...
1
vote
0answers
86 views

eigs routine in octave

I am using octave and observed a problem with the eigs-routine for non symmetric matrices. Using GNU octave version 3.8.1 the code below gives significant difference of eigenvalues although same ...
1
vote
0answers
62 views

Search Direction in Conjugate Gradient

Could you help me with a Conjugate Gradient question? In using CG to solve Ax=b, why is the search direction $p_{k+1}$ in CG chosen as a linear combination of the residual $r_k$ and previous direction ...
1
vote
1answer
134 views

Solve for a matrix given two vectors

I'm programming a beam finite element model by following a book (Nonlinear Finite Element Analysis of Solids and Structures Volume 2, in case you're wondering!). I've come across the following ...
2
votes
0answers
45 views

Algorithm code for Drazin and Bott-Duffin inverse (Matlab or C)

I could find the common Moore-Penrose algorithm, but I couldn't find the Drazin or the Bott-Duffin generalised inverse, except for some very specific cases, useless for my studying purposes. Is there ...
5
votes
0answers
58 views

Find the solution of linear equation using Wiedemann/ Krylov method

Let given $M =$ 1 0 1 0 1 1 1 1 1 and $b =$ 1 0 1 How to find the solution $x_3$ where ...
1
vote
1answer
64 views

Numerical eigenbasis for a unitary matrix

Do you know what numerical software computes an eigenvector basis for a unitary matrix? Say I have a unitary matrix $U$. If its eigenvalues are simple (no multiplicities), then for instance Matlab ...
2
votes
0answers
38 views

In-place QR update: deleting a column

Background I'm trying to do an update to a "thin" QR decomposition ($A = QR$, where $Q$ is $\mathbf{R}^{m,n}$, the first few columns (up to the matrix rank) of an orthogonal matrix and $R$ is upper ...
2
votes
2answers
214 views

Eigenvectors: MATLAB vs LAPACK DGGEV or DGGEVX

If we call LAPACK DGGEV or DGGEVX routines for two badly-conditioned matrices in a C++ code, will we get the same eigen-values ...
3
votes
1answer
194 views

Why is my MATLAB code for back-substitution slower than the backslash operator?

I wrote the code below to invert an upper triangular matrix, avoiding any possible multiplication/subtraction by zero. It just uses $\frac{1}{6}n^3+\ldots$ flops instead of $n^3+\ldots$ flops. ...
2
votes
3answers
125 views

Is it better to do normalization after all orthogonalization in Gram-Schmidt process?

In Gram-Schmidt process, is it better to do normalization after orthogonalization of all the vectors in a basis, or to normalize each new vector immediately after it is created, from computational ...
5
votes
1answer
92 views

Solving Generalization of Saddle point problem

I am interested in knowing if there is a generalization of the Uzawa iteration for the linear problems of the form $$\left[ \begin{array}{ccc}A& B^T&0\\ B&0&C^T\\ 0&C&0 ...
0
votes
1answer
85 views

Diagonalize a circulant-plus-rank-one matrix

Suppose $A=uv^T$ where $u$ and $v$ are non-zero column vectors in $\mathbb{R}^n, n\geq 3$. $\lambda=0$ is an eigenvalue of $A$ since $A$ is not of full rank. $\lambda=v^Tu$ is also an eigenvalue of ...
1
vote
1answer
100 views

writing linear system in matrix form

I have the following set of linear equations $$a_{m+1,n}+a_{m-1,n}+m(a_{m,n+1}+a_{m,n-1})+(m^{2}+n^{2})a_{m,n}=f_{m,n}$$ Here $m$ and $n$ run from 1 to $N$, so there are $N^2$ equations for the ...
3
votes
1answer
110 views

Get symmetric Finite Difference matrix in non Laplacian settings

I would like to solve a system of differential equations $u+\nabla(\nabla\cdot u)=f$ or in more detail $a+\partial_t^2a+\partial_t\partial_xb+\partial_t\partial_yc=f$ ...
5
votes
1answer
147 views

Least-squares for a diagonal matrix

This is a follow-up to a different question I asked with more detail. For $v\in\mathbb{R}^n$, denote $D_v\in\mathbb{R}^n$ as the diagonal matrix with elements in $v$. Given a "tall" matrix ...
3
votes
3answers
139 views

Large overdetermined system of linear equations

I'm looking for a method to solve a large overdetermined system of linear equations in a least squares sense. The matrix is dense. I'd like to use a method that works even with limited memory (we ...
4
votes
2answers
100 views

Solving “Hadamard systems”

Suppose we have two matrices $A$ and $B$ (we can assume they're symmetric; if absolutely necessary I think they may be positive definite). Then, is there any technique for solving $$(A\circ B)x=b,$$ ...
0
votes
2answers
89 views

How to determine the number of c points in algebraic multi grid

I am trying to write an algebraic multi-grid solver (in c++). At a given level I determine which nodes are c-points and which nodes are f-points (where the total number of c and f points equals the ...
2
votes
0answers
72 views

Does it matter if I use principal component analysis on the transpose instead of the original matrix?

My data set is a 60x10 matrix. I performed principal component analysis of this matrix with matlab using the princomp(AdjustedData) after I adjusting my original data set by subtracting the mean of ...
4
votes
2answers
114 views

Schur(QZ) Decomposition Differences

I am having issues understanding why different languages are producing different answers for the Schur(QZ) decomposition. I am working on writing some old stuff from Matlab into Julia and Python and ...
1
vote
0answers
28 views

linear solution of curve fitting on multiple linear functions differing by a multiplier

I am facing the following problem. I know nonlinear least squares can provide a solution but I am wondering if a linear way to solve this data fitting problem may exists. This is my input dataset: ...
1
vote
1answer
105 views

What is the algorithm that matlab used in its built-in function 'pca'?

Do anyone know what is the algorithm that MATLAB used in its built-in function "pca"? I have the following data set: 148.9820 55.8438 210.2150 149.3030 56.8891 208.4280 151.4400 ...
3
votes
1answer
173 views

Parallel Gram-Schmidt algorithms

I've heard that classical Gram-Schmidt is more amenable to parallelization than modified Gram-Schmidt; apparently the reason has something to do with "level 2 BLAS", which I'm not familiar with. ...
1
vote
1answer
142 views

Translate a 3D point along a heading

I need to translate a point (P1) in 3D a certain amount, call it stepSize, along a vector described by a heading composed of ...
0
votes
0answers
63 views

LU Factorization update when adding columns

I am looking for a way to update the LU factorization of a general $m \times n$ matrix after adding a column to the matrix. I have to iterate this procedure so I will begin with a matrix that is $m ...
0
votes
0answers
21 views

Updating the LU Factorization [duplicate]

I am looking for a way to update the LU factorization of a general m×n matrix after adding a column to the matrix. I have to iterate this procedure so I will begin with a matrix that is m×1 and go all ...
5
votes
2answers
195 views

Dimensionality reduction of the domain of f(x)

I'm wondering if there is something analogous to a PCA for data sets where there is a dependent variable. (Though I am interested in any method of dimensionality reduction, PCA is just an example.) ...
0
votes
0answers
51 views

Principal Components Analysis Not Behaving as I Expect it to

I have a bunch of points in $\mathbb{R}^3$ that I would like to translate and rotate so that their center is at the origin and the variance along the $x$ and $y$ axes are maximal (greedy, and in that ...
6
votes
1answer
151 views

Efficient RQ decomposition

I have an upper trapezoidal matrix stored in column major format. That is, my matrix looks like this: I'd like to RQ decompose it, and store Q in the rectangular part of my upper trapezoidal ...
2
votes
1answer
213 views

Comparing Eigenvectors, Mathematica vs. Matlab

I am trying to create the same out puts in Mathmatica and Matlab, however I am running into trouble aligning the eigenvectors with the eigenvalues, I think the Matlab is doing something slighly more ...
1
vote
3answers
89 views

rank-deficient NNLS

I want to find the minimum-norm solution to a rank-deficient least-squares problem, subject to positivity constraints, e.g. $$\min_x\ \|x\|^2 \quad s.t.\quad Ax = b,\ x \geq 0$$ where $A$ is large, ...
3
votes
1answer
99 views

How does matrix-matrix product scale with multiple CPUs?

These days, one can have 64 cores in a single node. I wonder how well the dense matrix-matrix product (SGEMM and DGEMM) scales ...
1
vote
1answer
109 views

Symbolic Computations with Block Matrices in Maple

Due to the derivation of an algorithm I tried to use Maple (9.5) to calculate some block matrix expression. Unfortunately Maple seems to ignore the assumption I set on the variable. Lets consider the ...
5
votes
2answers
233 views

Recommendations for symmetric preconditioner

Given symmetric, positive-definite matrix $A$ and its preconditioner $M^{-1}$. What's kind of preconditioner $M^{-1}$ which preserve the symmetry of $M^{-1}A$?
3
votes
1answer
73 views

State-of-the-art for active set optimization algorithms?

Given a problem like this: $$ \text{min } ||Ex-f|| \text{ s.t.}$$ $$ Gx \ge 0$$ $$ Cx = d $$ And assuming that the matrices are medium sized (dimensions in the low thousands) and dense, what's the ...
0
votes
0answers
58 views

Lapack : dgeev strange behaviour

I have a problem that is quite strange. I use c++ and dgeev Lapack to diagonalize and find the eigenvectors of a 36x36 real non symmetrical matrix. The documentation explains the following : if ...
5
votes
2answers
123 views

Sparse smallest eigenvalue problem on a linear subspace?

I am interested in methods for solving the optimization problem $$ \begin{array}{rl} \arg\min_x & x^T A x \\ \mathrm{s.t.} & x^T x = 1 \\ & Bx = 0 \end{array} $$ where $A$ is symmetric ...
7
votes
1answer
78 views

Do C++ matrix libraries translate compound vector operations to single loops?

I am trying to replicate Fortran90 array syntax using a C++ library. The libraries themselves are discussed at length in this question. They can all do something like this: ...
1
vote
1answer
60 views

Conditioning of matrix factorizations and square roots

For my application, I need factors $\tilde C$ so that $$ \tilde C{}^T \tilde C = C^TMC $$ where $C$ is a long matrix, i.e. $C$ has much more columns than rows, and $M$ is a small symmetric positive ...
5
votes
1answer
235 views

Testing if a matrix is positive semi-definite

I have a list ${\cal L}$ of symmetric matrices that I need to check for positive semi-definiteness (i.e their eigenvalues are non-negative.) The comment above implies that one could do it by ...
7
votes
1answer
171 views

Algorithm to calculate the exponential of an Hessenberg matrix

I am interested in computing the solution of a lage system of ODEs using a krylov method as in [1]. Such method involve functions related to the exponential (the so-called $\varphi$-functions). It ...
0
votes
1answer
115 views

solving underdetermined system of equations with a sparse matrix as input

I am using Matlab to solve Ax=b and my A is very large, sparse, binary and also rectangular. I saw the Matlab backlash \ operator help and it states that if A is rectangular then it will use the QR ...