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

0
votes
1answer
31 views

Least squares fitting

I have the following equation I came across which was solved using least squares $x = \sum_{n=1}^{N} A_{n} y_{n}$ Where $x$ is a $m \times p$ matrix and $y$ would be of size $m \times p$ as well ...
2
votes
0answers
55 views

Inverse problem with a rank-1 update

I hope you can help me out with this. I have to find the solution x to an inverse system $$ x=A^{-1}b $$ This inverse problem is basically a least square problem with a rank-1 update. $$ ...
0
votes
1answer
51 views

LSA, SVD and the Frobenius norm

In Latent Semantic Analysis one uses the SVD to perform a dimensional reduction of the term-document matrix, via the Eckart-Young theorem. Now, the rank $k$ approximation obtained by E-Y is proven to ...
3
votes
3answers
90 views

Efficient solver for a symmetric tridiagonal system where the upper/lower diagonals are offset

I'm looking for an efficient way to solve a symmetric tridiagonal system $Mx = d$, where the upper and lower diagonals of $M$ are offset from the main diagonal by $k$ rows/columns: $$ \begin{bmatrix} ...
1
vote
0answers
41 views

Lanczos algorithm with thick restart on a dynamic matrix

According to a recommendation, this is a re-post of that. currently, I'm working on a way to compute the 2 biggest eigenvalues of a real, symmetric, huge and sparse matrix that changes a few entries ...
0
votes
0answers
66 views

How to solve singular non symmetric poisson equation with Neumann boundary condtions?

I am trying to solve 2D Poisson equations with Neumann boundary conditions. When the mesh is uniform, Poisson equation is singular and symmetric, so the method listed in Null Space Projection for ...
4
votes
2answers
121 views

Applying matrix square root inverse in matrix-free regime

Let $A$ be a large symmetric positive definite matrix, and suppose that we can efficiently apply $A$ and have a fast solver to apply $A^{-1}$, but we do not have access to the matrix entries for ...
3
votes
2answers
135 views

Methods for solving linear systems

This is such a basic topic but there are so many different methods proposed for solving a linear system of equations. I recently found a very good source but couldn't really make sense of all the ...
3
votes
3answers
154 views

How to solve a Linear Matrix Equation: AX-XA=B efficiently?

recently I have been working on solving some math problems using Fortran. There occurs to me that a linear matrix equation: $$ AX-XA=B $$ where $A$ and $B$ are known $n\times n$ matrices and $X$ is ...
2
votes
1answer
135 views

Good Finite Element Library for a small project

I'm currently working on this project and I have a basic structural analyzer that uses the finite element method. Essentially, I turn each block into a set of trusses, construct a stiffness matrix ...
5
votes
1answer
62 views

Computing eigendecomposition of a Hermitian matrix that is almost unitary

I have a dense Hermitian matrix that is approximately unitary, so it has eigenvalues that are $\sim \pm1$. I would like to compute all the eigenvectors corresponding to the $+1$ eigenvalue (not ...
4
votes
0answers
48 views

Wanted: sequences of linear systems for recycling Krylov solver analysis

In the solution of sequences of linear systems $$A_ix_i=b_i\quad\text{for}\quad i=1,2,\dots$$ with Krylov subspace methods, data can be recycled from already solved linear systems in order to speed up ...
2
votes
1answer
70 views

Weighted Gauss-Seidel Algorithm

In Jacobi method's Wikipedia article there's a section that describes Weighted Jacobi method: http://en.wikipedia.org/wiki/Jacobi_method#Weighted_Jacobi_method. I need to implement the Weighted ...
1
vote
2answers
83 views

What sparse linear programming solver it is better to use?

I have the following LP problem: $$ \min \limits_{\varepsilon, x_{1}, \ldots, x_{n}}f(\varepsilon, x_{1}, \ldots, x_{n}) = \varepsilon \;\;\;\;\; s.t. \;\;C x \geq 0, \;\; x_{i}^{0} - \varepsilon ...
2
votes
3answers
91 views

Iterating through a 3D triangle

Given a triangular plane formed by three points in R3 space {p1, p2, p3}, I want to iterate through all points on the triangle plane by using two variables, x0,y0, something like in this example: ...
1
vote
0answers
39 views

Sparse matrix factorization of a rank deficient matrix by decomposition into linearly independent components

I've got a little conjecture I need to prove for a theoretical result related to causal Bayes net search with latent variables under sparsity constraints. If you're interested in the application ...
3
votes
0answers
32 views

BLACS context value and multiple MPI communicators

I am reposting here a question previously asked in stackoverflow: as suggested in a reply to my question, scicomp could be a better place to obtain some useful comments/suggestions for my problem. I ...
1
vote
2answers
89 views

How to find the smallest positive eigenvalue of a large general system if they are all in +/- pairs of real eigenvalues

I used MKL Lapack SBGVX because it can solve for "selected" eigenvalues/modes (both positive and negative), thinking it would be efficient, but it is extremely slow when compared with Bathe's subspace ...
3
votes
1answer
60 views

spectral decomposition in Numpy, sign difference

I am trying to follow along with an example from a book, but I get seemingly different answers depending on which spectral decomposition function I use in Numpy. I am trying to transform the Matrix G, ...
3
votes
1answer
53 views

Does the covariance matrix in Least Squares depend upon the input data?

I had always assumed that the covariance matrix depends upon the amount and quality of your input data, but I am finding out that this is not the case. Is this true? We want to fit $f(t) = ...
0
votes
1answer
43 views

Truncated SVD implementation in Java

I need the Truncated SVD implementation in java. I need to pass a matrix of doubles and an integer value representing the rank where to filter out noise. In output i need a filtered matrix of doubles. ...
1
vote
1answer
65 views

large symmetric positive band matrix

I use gpbsv command from Intel MKL to solve symmetric positive band system. But unfortunately when the system is large I get an error Access violating writing location in VisualStudio. Could someone ...
1
vote
3answers
107 views

Can numpy.linalg.solve use back substitution when possible?

The question is if Python Numpy library can use back subsitution to solve Ax=b if possible, that is, if A is lower triangular? Do numerical linear algebra packages do this? I would think Numpy would ...
0
votes
2answers
75 views

What to do with singular (non-invertible) rotation matrix

I have an orthotropic material with a (6x1) stress vector known in the global coordinate system and yield surfaces known in a local coordinate system. So far I have only needed to convert from local ...
4
votes
0answers
62 views

Tracking for two meshes

I'm dealing with physical simulation (position based dynamic). Now I'm trying to realize tracking for two meshes. In order to explain what does it mean, let's assume that we have two similar meshes: ...
7
votes
2answers
197 views

Which novel data structures are used in adaptive FEM?

A lot of adaptive FEM libraries use more advanced mesh data structures to handle adding/removing nodes, edges, triangles, tetrahedra, etc. For example, the p4est library uses octree data structures ...
0
votes
1answer
57 views

sorting adjacency matrix by the Fiedler vector

I have a fairly sparse adjacency matrix showing connections between approx 5,000 points in my dataset. I'm looking at various ways to analyze the relationships between the data points. This approach ...
1
vote
0answers
77 views

Algorithms to compute largest gap between smallest nonzero eigenvalues of sparse symmetric matrix

I am looking mainly for c/c++ implementations but also for theoretical algorithms to compute gaps between smallest positive eigenvalues of symmetric, singular matrix or real numbers. To be precise, I ...
0
votes
0answers
17 views

Minimum sparsity of two randomly populated matrices leading to performance increase when switchting from dense to sparse matrix multiplication

Given two 2D matrices A, B which are going to be multiplied by a matrix multiplication algorithm. Properties of matrix A The columns are grouped into blocks of size N. Within each block identified ...
10
votes
1answer
136 views

Is there any way to do “double preconditioning”

Question: Suppose that you have two different (factored) preconditioners for a symmetric positive definite matrix $A$: $$A \approx B^TB$$ and $$A \approx C^TC,$$ where the inverses of the factors $B, ...
0
votes
0answers
21 views

Markowitz Pivoting to reduce size of a dense integer system

I am dealing with a large sparse integer matrix that I need to find the nullspace of. I've seen Markowitz Pivoting come up in several places discussing similar problems such as here: ...
4
votes
0answers
64 views

Fast Eigenvalue and SVD Solver for Structured Matrices

I am looking for a fast Eigenvalue and SVD solver for small dense structured matrices (Hankel and Toeplitz). I have searched for efficient implementations in libraries like MKL but I am not able to ...
0
votes
0answers
101 views

discretizing $\frac{d^2}{dx^2}$

I am trying to discretize the differential operator $\frac{d^2}{dx^2}$ acting on $S^1 = [0,1]$ using finitely many points around a circle at $0, \frac{1}{N}, \frac{2}{N}, \dots, \frac{N-1}{N}$. Here ...
4
votes
1answer
93 views

Methods for calculating the symmetric part of a matrix

I am using a multigrid preconditioned GMRES method for a nonsymmetric matrix. The matrix is the discretisation of the derivative of a nonlinear operator. Since multigrid is not the best for ...
0
votes
0answers
42 views

Is there any backward error when applying SVD to solve linear systems

I checked on-line and only found some backward error bounds to apply LU or QR to solve linear systems. I hope any one can point me to some reference about SVD.
0
votes
0answers
19 views

MLLL algorithm for sparse, integer bases to find a nullspace

I am trying to find a suitable algorithm that can find a basis for the nullspace of a sparse, integer matrix. Reading A Course in Computational Algebraic Number Theory by Cohen, an algorithm based ...
1
vote
0answers
44 views

Efficient way to do congruent transformation using matrix inverse?

I know a square self-adjoint matrix $S_{vv}$ and I want to find: $S_{rr} = HS_{vv}H^{\dagger}$ where $\dagger$ denotes conjugate transpose. I do not know $H$ but I do know $H^{-1}$. What is the ...
2
votes
2answers
64 views

Algorithm for directly finding the leading eigenvector of an irreducible matrix

According to the Perron-Frobenius theorem, a real matrix with only positive entries (or one with non-negative entries with a property called irreducibility) will have a unique eigenvector that ...
2
votes
4answers
154 views

cartesian products in numPy

Have two arrays, in my case $X = \{1,2,\dots, n\}$ or X = np.arange(n). How do I get $Y = X \times X = \{ [i,j]: 1 \leq i,j \leq n \}$ as a 2D array in numPy? ...
4
votes
1answer
113 views

Nullspace algorithm for a sparse matrix

I am dealing with large, sparse, rational matrices that I need to determine the nullspace of. Currently, I have one that is about 12000x12000 (but not square), where one in every 2000ish elements is ...
2
votes
1answer
104 views

Question about extending Tikhonov regularization

I know that the Tikhonov regularization of a linear system has an analytical solution given by: \begin{equation} \hat{\mathbf{x}} = \mathrm{arg\;min}\left( \left| \mathbf{Ax} - \mathbf{b} \right|^{2} ...
2
votes
1answer
74 views

Algorithm for Eigenvalue Problem of a Real Symmetric nxn Matrix

I have a nxn covariance matrix (so, real, symmetric, dense, nxn... I mean positive semi-definite). 'n' may be very very very big! I'd like to solve partial (3 largest) eigenvalue (+eigenvectors) ...
4
votes
1answer
146 views

What is the best way to multiply a diagonal matrix (in fortran)

What is the best way to compute: $$ Y = D X $$ where $D \in \mathbb{R}^{m\times m}$ is diagonal and $X \in \mathbb{C}^{m \times n}$ is general. I am mostly interested in these two cases: $m ...
0
votes
0answers
34 views

polynomiality of a problem

I have a feasibility problem as follows: Does there exist a $p\geq 0$ such that $Ap+Gq=b$ for every non-negative $q\in D$. Here $A,G$ are adjacency matrices and $D$ is a convex set. Is this problem ...
2
votes
2answers
76 views

Advice on the regularisation of a linear problem

I'm numerically inverting an integral transform using a method suggested by a scicomp user from an earlier question. The problem is as follows: I wish to estimate $f(x)$ for a given $F(y)$, both of ...
1
vote
0answers
46 views

full rank update to cholesky decomposition for multivariate normal distribution

This question is a specialization of full rank update to cholesky decomposition, to which I hope to get a more positive answer. When calculating the minus log of the multivariate normal distribution, ...
2
votes
1answer
81 views

full rank update to cholesky decomposition

Let $A$ be a real, symmetric, positive definite matrix. It has at least 500 rows, possibly much more. I compute its Cholesky decomposition, which allows me to calculate $det(A)$ $A^{-1}X$ for some ...
0
votes
2answers
121 views

Algorithm to find non-negative integer solutions to x_1 + x_2 …=n

I know the number of solutions to the equation $$x_1+x_2+x_3+...+x_k=n$$ is given by $\binom{n+k-1}{k-1}$. Is there an algorithm to actually find all the solutions to this equation, without having to ...
9
votes
1answer
114 views

Purely rotational least squares match

Could anyone recommend a method for the following least squares problem: find $R \in \mathbb{R}^{3 \times 3}$ that minimizes : $\sum_{i=0}^N (Rx_i - b_i)^2 \rightarrow min$, where $R$ is unitary ...
3
votes
2answers
137 views

Algorithm for Complete Eigenvalue Problem of a Real Symmetric nxn Matrix

I have a nxn covariance matrix (so, real, symmetric, dense, nxn). 'n' may be very very very big! I'd like to solve complete eigenvalue (+eigenvectors) problem for this matrix. Could somebody tell me ...