In linear algebra, an $m$x$n$ matrix maps vectors from space $R^m$ to $R^n$. Matrices are distinguished from multidimensional arrays, which is a general term for any rectangular array of information.
3
votes
1answer
56 views
$\mathbf{UDU}^\top$ decomposition routines in LAPACK/Eigen?
I would like to compute the decomposition of a real symmetric positive definite matrix $\mathbf{A} = \mathbf{UDU}^\top$.
LINPACK seems to have it as DSIFA, but I ...
1
vote
1answer
66 views
LU Decomposition with memory-mapped matrices
I have a ~4.12 Tb structured relatively-sparse matrix dataset (about 8% of the matrix entries are non-zero) that i want to apply an LU decomposition, however, given the size of it, loading it in ...
2
votes
0answers
59 views
SVD computation with “initial guess”
Suppose I have some matrix A whose SVD I know. Now, I am given B which is A plus some small ...
1
vote
0answers
28 views
On solution of a class of discrete-time Lyapunov equation for systems with multiplicaitve noise
Let's consider the following equation
$$X=F_{1}XF_{1}^{T}+...+F_{p}XF_{p}^{T}+C$$
where $p$ is an positive integer and $C$ is a known positive semidefinite matrix. If we augment $F=[F_{1}...F_{p}]$ ...
10
votes
2answers
213 views
How is the SVD of a matrix computed in practice
How does MATLAB, for instance, calculate the SVD of a given matrix? I assume the answer probably involves computing the eigenvectors and eigenvalues of A*A'. If ...
2
votes
2answers
110 views
Efficient computation of Markov chain transition probability matrix
Consider a continuous Markov chain $X=(X_t)$ on a finite state space and let $Q$ be the (given) transition rate matrix. This matrix is very sparse, with non-zero values on 3 diagonals only (so from ...
4
votes
1answer
120 views
Are there any specialized methods available for solving structurally symmetric sparse linear systems?
When solving $Ax=b$, prior knowledge about $A$'s structure can help in designing an efficient solver which exploits this information (e.g conjugate gradient method is to be used when $A$ is ...
6
votes
1answer
45 views
Bad scaling versus collinearity
I was trying to solve a linear system:
$$
\mathbf{A}\mathbf{x} = \mathbf{y}
$$
but the conditioning number was quite bad (around $10^{17}$).
I thought that the system was singular, but after scaling ...
2
votes
1answer
27 views
Trying to implement a simple/efficient combinations function in MATLAB
So, recently, I have found myself in the position of having to implement a combinations function in MATLAB. What I mean by this is the following: I simply need to list all possible combinations for an ...
5
votes
2answers
111 views
Rearrange an ordinary matrix to block diagonal form
Is there an algorithm to rearrange a matrix into block diagonal form, given that the matrix is block diagonal in nature but randomized with an unwise choice of basis?
In particular, are there any ...
4
votes
1answer
89 views
Parallel computation of big covariance matrices
We need to compute covariance matrices with sizes ranging from 10000x10000 to 100000x100000. We have access to GPUs and clusters, we wonder what is the best parallel approach for speeding up these ...
1
vote
0answers
53 views
Identifying unknown variables in a graph
I have graphs with just X and Y values. Is it possible to dissect a graph to find out number of parameters(Degrees of freedom) that fit the data-set?
For e.g I have a graph of Time Vs Distance and I ...
10
votes
4answers
166 views
Calculation of the sparsity structure for finite element matrices
Question: What methods are available to accurately and efficiently calculate the sparsity structure of a finite element matrix?
Info: I'm working on a Poisson Pressure Equation solver, using ...
4
votes
1answer
62 views
Sparse LU for block-sparse matrices
I frequently need to solve linear systems with sparse matrices of moderate dimension (say a few thousand). These matrices are composed entirely of small dense blocks (typically 5-10 in dimension), and ...
5
votes
2answers
223 views
Computing the pseudoinverse of a 3x3 matrix
I need to compute the (Moore-Penrose) pseudoinverse of fixed-size 3x3 matrices. I would prefer to have a simple method without bringing in the industrial strength machinery of Lapack. Are there any ...
5
votes
1answer
127 views
Eigenvalue Decomposition of Hermitian Matrix in Scala
I'm working on helping my friend create code to perform the Overlap Dirac Operator and have come across one part that I'm not sure what to do.
I need to compute the Eigenvalues and corresponding ...
3
votes
1answer
107 views
restriction and interpolation in multigrid method
I need detailed explanation of the formula below
A2=I1*A1*I2
I suppose this formula computes matrix A2 on a coarse grid and here A1 is original matrix on fine ...
1
vote
0answers
21 views
PCA performed on a configuration with scaled axes
Suppose a configuration $X\in\mathbb{R}^{n\times 2}$ is output of PCA on some high-dimensional data $Y\in\mathbb{R}^{n\times h}$. Note that this PCA is performed by $$X=Y\cdot U,$$ where columns of ...
2
votes
1answer
61 views
Optimality criterion of PCA via recovered distances
It is stated in
http://users.eecs.northwestern.edu/~yingwu/teaching/EECS510/Reading/Williams_NIPS01.pdf
that the PCA mapping from $h$-dimensional data to low $k$-dimensional space minimizes ...
3
votes
0answers
80 views
How to accurately decompose positive semidefinite matrix and use the lower triangular part in linear equations
I have $n$ arbitrary $p\times 1$ vectors $x_i$, and $p\times k$ matrices $A_i$, and $n$ $p \times p$ positive semidefinite matrices $S_i$, where some (often most) of the $S_i$'s are same (for example ...
9
votes
2answers
347 views
Row major versus Column major layout of matrices
In programming dense matrix computations, is there any reason
to choose a row-major layout of the over the column-major layout?
I know that depending on the layout of the matrix chosen, we need to ...
2
votes
1answer
151 views
How to apply a Galerkin finite element method to a linear, one-dimensional boundary value problem
I have the following boundary value problem:
$$-(\alpha u')' + \gamma u = f $$
in $\Omega = (a,b)$ with b.c. $u(a) = u(b) = 0$
and $\alpha > 0, \gamma ≥ 0$ and $f:(a,b) \to \Re $
The weak ...
4
votes
2answers
107 views
Imposing invertibility on a Matrix
I have a symmetric positive semidefinite covariance matrix $A$, which is approximately computed as the output of a quadratic regression. I then need to invert $A$, but often it is close to singular. ...
1
vote
1answer
42 views
Existence of a solution at a stationary equation for quadratics
Given a convex quadratic function $f(x)$, to obtain a solution for which $f(x)$ has minimal value one sets $\nabla f(x)=0$, and solves for $x$. Suppose that the result of differentiation of convex ...
4
votes
2answers
64 views
What is the algorithm for computing block reflectors in xDLARFB
The theory behind computing a single Householder reflector to zero out part of a column of a matrix is pretty well described in Matrix Computations by Golub and Van Loan. However, the blocked ...
5
votes
3answers
120 views
What is the criteria for switching between Strassen's and Regular matrix multiplication Algorithms
Strassen’s Matrix Multiplication algorithm has theoretical performance of $ O( n^{log_2 7}) $. Regular MM algorithm has performance of $ O( n^{3}) $.
At certain sizes of matrices (lets call it ...
2
votes
1answer
45 views
Configuration shift for determination of a true dimensionality
What would then be the way to determine a true dimensionality of a configuration of points $X\in\mathbb{R}^{n\times k}$ based on its Gram matrix $G=XX^T$? The "true" dimensionality refers to the ...
0
votes
2answers
92 views
Affecting the rank of a Gram matrix by configuration shift
Let certain configuration of $n$ points exist in $d-$dimensional space, $X\in\mathbb{R}^{n\times d}$, $d<<n$. Also, let the corresponding Gram matrix be defined as $G=XX^T$.
Since $X$ exists in ...
4
votes
1answer
149 views
How to get sparse complex matrices from my code to PETSc efficiently
What is the most efficient way to get a complex sparse matrix from my Fortran code to PETSc? I understand that this is problem dependent, so I tried to give as many relevant details as possible below.
...
5
votes
3answers
253 views
Fastest algorithm to compute the condition number of a large matrix in Matlab/Octave
From the definition of condition number it seems that a matrix inversion is needed to compute it, I'm wondering if for a generic square matrix (or better if symmetric positive definite) is possible to ...
8
votes
1answer
110 views
Algorithms for Large Sparse Integer Matrices
I'm looking for a library that performs matrix operations on large sparse matrices w/o sacrificing numerical stability. Matrices will be 1000+ by 1000+ and values of the matrix will be between 0 and ...
3
votes
1answer
259 views
Power Iteration on general matrices (with higher multiplicity of dominant eigenvalue)
To compute the eigenvector corresponding to a dominant eigenvalue of a matrix $A\in\mathbb{R}^{n\times n}$, one could apply the Power Iteration: $$v_1=\frac{Av_1}{\|Av_1\|}.$$
1) in case $A$ is ...
1
vote
1answer
82 views
Configuration shift to change the rank of a Gram matrix
Suppose a matrix $D\in\mathbb{R}^{n\times n}$ of Euclidean distances between $n$ points is given. To obtain a Gram matrix (matrix of inner-products of points that give rise to distances in $D$), one ...
3
votes
1answer
80 views
multiplications of graph adjacency matrix
Suppose $A$ is a directed graph adjacency matrix. Is there any good interpration of the $(i,j)-$entry of the matrix $(A^{32}\cdot (A^T)^{32})$ ?
3
votes
1answer
145 views
Manipulating a generalized eigenvector problem to plain eigenvector problem
Let $X\in\mathbb{R}^{n\times p}$ denote a matrix with $p$ linearly-independent columns, and let $L\in\mathbb{R}^{n\times n}$ denote a symmetric matrix. Furthermore, let $D\in\mathbb{R}^{n\times n}$ ...
5
votes
3answers
297 views
Eigenvectors with the Power Iteration
To compute the eigenvector corresponding to dominant eigenvalue of a symmetric matrix $A\in\mathbb{R}^{n\times n}$, one used Power Iteration, i.e., given some random initialization, ...
4
votes
2answers
234 views
Extracting a 2D matrix from a multidimensional array in C
In my c language program, I have to store multiple dense $m\times m$ matrices corresponding to gridpoints $x_i$ with $i=1,...,n$. I decided to create a three dimensional array $A\in R^{n\times ...
0
votes
3answers
131 views
Equivalence of linear systems, solving one instead of the other
This question is related to recently posted one, but I guess it deserves a separate attention.
Suppose a symmetric matrix $L\in\mathbb{R}^{n\times n}$ is given, and a rectangular matrix ...
2
votes
1answer
83 views
Recovering coordinates by eigendecomposition without double-centering
Suppose an Euclidean distance $D\in\mathbb{R}^{n\times n}$ matrix between a set of $n$ objects is given. To obtain inner-products (which will be further be used to recover coordinates), entries of $D$ ...
0
votes
2answers
167 views
Gram-Schmidt method to identify linearly dependent vectors
A method to orthogonalize a set of vectors (vectors of unit length that are mutually orthogonal) is the Gram-Schmidt process:
http://en.wikipedia.org/wiki/Gram%E2%80%93Schmidt_process
Note that the ...
1
vote
3answers
153 views
Convergence of the gradient descent and linear vs non-linear fixed point iteration
Suppose a system $$Ax=b$$ is given, with $A\in\mathbb{R}^{n\times n}$ being a symmetric positive-definite matrix, and some non-zero $b\in\mathbb{R}^n$. The gradient method with optimum step length can ...
11
votes
4answers
538 views
Dealing with the inverse of a positive definite symmetric (covariance) matrix?
In statistics and its various applications, we often calculate the covariance matrix, which is positive definite (in the cases considered) and symmetric, for various uses. Sometimes, we need the ...
0
votes
0answers
24 views
How do I construct cross-correlated event strings?
Suppose I have groups of event strings that I can easily self-correlate through having a specified pool size from which these event strings can be constructed from. For example,
I need 10 event ...
6
votes
1answer
89 views
Sudden drops in matrix multiplication performance
I've been reading about implementing dense matrix multiplication when the matrix doesn't fit in cache. One of the graphs I've seen (slide 9 from these slides) shows sudden drops in performance using ...
5
votes
1answer
102 views
Is there a way to inspect the graph of a sparse matrix with PETSc?
I am currently trying to code the CA-CG method within the PETSc framework. A mandatory step in this process is the implementation of the "matrix powers kernel" algorithm for a generic sparse matrix.
...
2
votes
1answer
103 views
Root Convergence rate of Iterative Scheme
I have an iterative sequence for optimizing an EM (Expectation Maximization) algorithm based loss function $L(X)$ with $t$ being the iteration number as:
$X_t=ABX_{t-1}+CX_{t-1}+X_{t-1}$ where $A$ is ...
2
votes
3answers
231 views
derivative of linsolve
Consider a vector $\mathbf{g} \in \mathbb{R}^{m}$ and a matrix $\mathbf{A} \equiv \mathbf{A(g)} \in \mathcal{M}_{p\times q} [\mathbb{R}]$, a function of $\mathbf{g}$.
Furthermore, let $\mathbf{S} ...
3
votes
1answer
48 views
Calculation of Multivariate Coherence
I trying to detecting whether a data set of time series has a global change in frequencies. Calculating the average (or median) pairwise coherence, I feel, misses the point because I am trying to get ...
3
votes
1answer
73 views
High-dimensional representation of arbitrary input
Given a symmetric matrix $A\in\mathbb{R}^{n\times n}$ with positive entries and zero diagonal, is it always possible to construct a high-dimensional configuration in Euclidean space, such that these ...
3
votes
2answers
156 views
Positive semi-definiteness of a (symmetric) matrix
Suppose a matrix $A\in\mathbb{R}^{n\times n}$ is given. Faced with a proof for $$x^TAx>0,$$
for a non-zero vector $x\in\mathbb{R}^{n}$, I was thinking to use the information of the spectrum of $A$ ...
