Questions tagged [linear-solver]
Referring to methods for solving linear systems of equations.
313
questions
0
votes
1answer
65 views
Solving large sparse system
I am working on a problem with very large sparse matrices. I'd like to compute $A^{-1} B$, that is a crucial part of converting DAE to ODE (and there is no workaround). Here size of $A$ is 2E+5 x 2E+5 ...
1
vote
1answer
98 views
Efficient way to solve a set of linear equations $Ax=b$ when $A$ is sparse and some elements of $b$ are equal to zero
I have a set of linear equations, $Ax=b$. And about half of the elements in the right-hand side (vector $b$) are equal to zero. My system matrix $A$ is a sparse complex matrix. And $A$ is in the size ...
7
votes
2answers
243 views
Is there an iterative solver for dense matrices with possible zero diagonal entries?
Is there an iterative solver that can handle potentially zero entries on the central diagonal? I am implementing a polynomial fitting algorithm (up to $10^{th}$-order) and my matrix is a "...
0
votes
0answers
45 views
Solver for large dense BVP system in python
I have a large system of boundary value problems of the form
$$ \frac{d^2 y }{dt^2} = C(t) y + b(t), $$
where the variable $y$ is a vector that has anywhere from 50 to around 500 components, $C$ is a ...
0
votes
1answer
54 views
Solving a sparse linear system using transpose of lower triangular matrix without copying
I have a sparse lower-triangular matrix $L$, and a right-hand side $b$, and I'd like to solve the linear system
$$L^T x = b$$
but without explicitly creating $L^T$. Ideally, I could write something ...
2
votes
1answer
88 views
Solving an m x m symmetric linear system involving a matrix multiplication versus an (n+m) x (n+m) system
Suppose that $R$ and $D$ are an $n \times m$ and $m \times m$ matrices. Assume that $m \ll n$ and that $D$ is positive definite. We would like to solve the system $(R^T R + D) x = R^T b$. This ...
3
votes
1answer
112 views
Is there any reason to scale a matrix before (sparse) Cholesky decomposition?
I have a sparse symmetric positive-definite matrix $M$ and I expect the entries in some rows/columns to have very different orders of magnitude (up to a factor of $10^8$) than the entries in others.
...
1
vote
2answers
184 views
Why OpenFOAM uses its own data structures and linear solvers?
I wonder why OpenFOAM code has its own data structures Lists, HashTables, ... etc. when there is the STL in C++?
Another ...
3
votes
1answer
86 views
What happens when I use a conjugate gradient solver with a symmetric positive semi-definite matrix?
I have a symmetric positive semi-definite matrix, i.e., a laplacian and wonder what may happen when I use a CG solver, that is an algorithm for positive definite matrices.
What happens when the ...
1
vote
2answers
61 views
Solution of symmeric/non-symmetric linear system
I would like to understand what happens in the following:
I have a really simple Poisson problem, in 1D, with $u_0 = u_N = 0$. I assembled the stiffness matrix and the right-hand side, and I applied ...
4
votes
2answers
132 views
More stable method of back substitution?
I've been tinkering a little in Fortran (2008) and wrote the following to solve
$Rx=b$ for $R\in\mathbb{R}^{n\times n}$ upper-triangular, $x,b\in\mathbb{R}^n$.
My code looks like this:
...
1
vote
1answer
71 views
Solving triangular matrix equations on a GPU
Suppose I have these two $N\times N$ lower triangular banded matrices:
$A = \begin{bmatrix}
a_0 & & \\
a_1 & a_0 & \\
a_2 & a_1 & a_0 \\
a_3 & a_2 & a_1 & a_0 \\
&...
1
vote
1answer
69 views
Parallel solution of PDEs
Which is the best approach to solve a PDE in parallel:
1.To split the mesh the mesh in N parts and every processor works on its own part or
2.To take the global linear system Ax=b and solve it in ...
3
votes
3answers
225 views
Using matrix exponential to solve linear system
Consider the system of linear equations:
$$
Ax=b
\tag{1}
\label{eq1}
$$
where
$A\in\mathbb F^{n\times n}$, diagonalizable dense matrix, over the field $\mathbb F$ of real or complex numbers,
$x\...
2
votes
1answer
117 views
How to solve for f(A)x=b without GMRES?
How to solve for $f(A)x=b$?
For GMRES, an answer is given in this book chapter: http://link.springer.com/chapter/10.1007%2F978-3-642-58333-9_2.
Ungated version:
https://www.researchgate.net/profile/...
2
votes
0answers
61 views
Residual of Poisson equation with periodic boundaries
I am trying to write a multigrid solver for Poisson's equation, $-\Delta u=f$, on the unit square, $\Omega=(0,1)^2$ with periodic boundaries. My primary source has been Multigrid by Trottenberg, ...
2
votes
0answers
62 views
Solving a huge least squares system of equations when I can only evaluate Ax
I have a situation where I can generate a system of $M$ linear equations for $N$ variables ($N \ne M$). Implicitly this is of the form $Ax=b$ with $A \in \mathbb{R}^{M \times N}$, although I never ...
1
vote
2answers
91 views
How to understand the storage of the Hessenberg matrix of Krylov subspace matrix?
For the Krylov subspace method to solve the large sparse linear system, we first need to generate a subspace Km = span{v,Av,...A^{m-1}v}, which indeed a process ...
0
votes
0answers
49 views
How to implement the gmres method using Householder transformation instead of the Gram-Schmidt?
For Generalized Minimal Residual method GMRES, we usually use the Modified Gram-Schmidt MGS to generate an orthonormal basis of ...
3
votes
0answers
54 views
Levinson Recursion for Non Square Toeplitz Matrices
Given a rectangular Toeplitz Matrix $ H $, how could one solve:
$$ y = H x $$
For instance, $ H $ can be Linear Convolution Matrix of the filter $ h $:
$$ H = \begin{bmatrix}
{h}_{1} & 0 & ...
2
votes
1answer
107 views
Testing a block tridiagonal system of equations
In 1D problems, tridiagonal systems of equations are obtained when we use finite-difference or finite-volumes in a structured mesh. A wide solver is the TDMA algorithm here. In two-dimensional ...
4
votes
0answers
42 views
Nonlinear least squares optimized Jacobian calculation
I have a nonlinear least squares problem, in which I am trying to minimize residuals which can be divided into four classes:
$$
\min_x ||\epsilon(x)||^2 + ||\xi(x)||^2 + ||\delta(x)||^2 + ||s(x)||^2
$$...
4
votes
3answers
153 views
Solving for a vector in a linear system that is both left and right multiplied
I have a linear system where I am given 2 matrices, $A$ and $B$, and 2 vectors, $v$ and $c$, and I need to solve for the vector $x$. $A$ is $n\times n$, $B$ is $n \times n \times n$, and the vectors $...
2
votes
1answer
60 views
What is the standard, extrapolation, and modified version of Richardson iteration method?
I have been studying the iterative methods recently. For classical iterative methods solving $Ax=b$, I have seen that the most simplest iteration method is the so-called "Richardson iteration". But I ...
2
votes
1answer
74 views
What method to solve a sparse complex symmetric (non-Hermitian) system?
I have a sparse system (about 78% of zero entries) that is complex and symmetric (but not Hermitian). The following figure shows the structure of the problem. The off-diagonal blocks are incidence ...
2
votes
1answer
124 views
Why Krylov subspace iterative methods are faster than classical iteration?
This semester, I have been studying the most popular iterative methods, i.e., Krylov subspace iteration methods.
For a large sparse system linear
$$
Ax=b,
$$
where $A$ is nonsingular, I know that ...
1
vote
0answers
31 views
What kind of problem or matrices are suitable for multigrid method?
For Poisson or Convection-diffusion equation as follows:
$$
-\Delta u=f,\qquad u|_\Omega = g.
$$
or
$$
-\Delta u +\vec{w}.\nabla u=f,\qquad u|_\Omega = g.
$$
using FDM or FEM discretization, we can ...
0
votes
0answers
37 views
Why does the initial guess for linear system usually choose by zero vector?
For solving linear system
$$
Ax=b,
$$
using iterative mehods, we often use the terminate criterion as follows:
$$
\frac{\|r_k\|}{\|r_0\|}=\frac{\|b-Ax_k\|}{\|b-Ax_0\|}<eps.
$$where $x_0$ is the ...
2
votes
1answer
51 views
Partially Banded Matrix
I have a somewhat peculiar Jx=R system that I need to solve. The matrix J is 2N -by- 2N. The first N rows have all entries filled. The next N rows are banded in two places, i.e. for the (N+k)th row, ...
0
votes
1answer
60 views
Method to calculate solution of a linear equation system?
I am searching a solution method for the following equation system of equation systems:
Let $A, B \in \mathbb{R}^{n \times n}$ be s.p.d. Matrices and $O$ be the zero matrix of the same size. Further ...
5
votes
0answers
168 views
Symmetric sparse direct solvers in scipy
scipy.linalg.solve, in its newer versions, has a parameter assume_a that can be used to specify that the matrix $A$ is symmetric ...
1
vote
0answers
58 views
Hybrid Ellpack-Itpack (ELL) + COO Sparse Matrix Representation decomposition threshold
Hybrid ELL-COO sparse matrix representation can be done as in the picture, I was looking intensively, however I couldn't find out what is the threshold of decomposing the original matrix into ELL part ...
1
vote
0answers
61 views
implementation for coppersmith matrix multiplication
Is there any online implementation for the coppersmith matrix multiplication I have searched alot but can not find any? and if there is not any why is that Isn't this algotithm much faster than ...
0
votes
0answers
155 views
What algorithm do BLAS and ATLAS use for matrix multiplication?
I have searched and what I understood was that they use the naive one with several memory and cache optimizations.
However, I wanted to know whether they are using the Strassen or the Coppersmith-...
1
vote
1answer
100 views
Lapack symmetric update $B^{-1}AB^{-T}$
Does Lapack have a routine that, given symmetric $A=A^T$ and $B$, computes the symmetric matrix $B^{-1}AB^{-T}$ (while preserving symmetry exactly)?
It would be enough to have this routine for ...
0
votes
1answer
77 views
The proper way to assess the error of Jacobi iteration (for 2D Poisson equation)?
Motivation: I'm using 2D regular grid (it's actually a quadtree but I can still treat it as a finite difference thing if I weight-average the solution over smaller scale cells for the purpose of ...
2
votes
2answers
247 views
Why iterative method: AMG preconditioned PCG is slower than Matlab direct method 'A\b'?
Recently, I have met a question that
a saying goes that for large linear system: iterative methods are required because of memory problem of direct methods.
But when I implement some experiments ...
3
votes
1answer
73 views
Why are the round-off errors when solving the linear system $Ax = b$ of order $\varepsilon_\text{mach} x_j$?
I was reading a paper on arXiv where, in Section 2.4, the authors are discussing the error that arises in the solution of a linear system
$$Ax = b,$$
or, to match up better with the paper,
$$\Phi \...
1
vote
1answer
94 views
Sparse matrix inversion
I have the impedance matrix $Y$, formulated from an electrical network by augmented nodal analysis. The matrix $Y$ is shown as an image to illustrate its feature visually, where all the white blocks ...
0
votes
1answer
162 views
How to set an initial guess for the iterative solver in Comsol?
How to set the initial guess for the iterative solver GMRES or FGMRES for linear problems (Helmholtz equation of RF module) in Comsol?
5
votes
1answer
145 views
Complex differentiation of linear solvers
I have a linear system $$Ax=b$$ which I'm solving approximately, and I need to take the frechet derivative of x with respect to z. Were I solving the problem exactly (either analytically or to machine ...
7
votes
0answers
79 views
How to construct an effective preconditioner for this particular problem
A quick introduction to my problem
I am currently developing a method for simulation of water waves in three dimensions based on potential flow theory. The computational bottleneck of the method is ...
1
vote
0answers
45 views
Stability of SVD, Eigendecompositions, and pseudoinverse procedures in modern LAPACK routines
I have proposed an optimisation algorithm which I claim has improved upon the previous algorithm in a number of ways.
One of these claims is that my proposed solution requires no explicit SVD and ...
3
votes
1answer
478 views
A fast way to check if a Matrix is ill-conditioned, and turning it into well-conditioned
I'm running a simulation, and some linear solvers are returning a message of ill-conditioned matrix.
Hence, I'm looking for a fast, easy to implement, method to detect if a matrix is ill-conditioned, ...
6
votes
2answers
121 views
Efficient approach for solving matrix plus diagonal matrix system that varies in time
When solving a system of ODEs, as part of a preconditioner, I get the system
$(A + D(t))x = b(t)$
where $A$ is a sparse matrix and $D(t)$ is diagonal. I'm currently solving this by taking the LU-...
2
votes
1answer
88 views
Assessing numerical error in solving a least squares problem
I have a linear system of the type
$$Ax = b$$
I want to minimise $|b - Ax|^2$.
I know there are different approaches to directly solve the system (Normal equation + Cholesky, QR decomposition, SVD ...
2
votes
0answers
116 views
Richardson's Iteration, Gradient Method and Spectral Radius
Richardson's iteration introduce a scalar $\alpha$ to the update formula:
$$ \textbf{x}^{(k+1)} = \textbf{x}^{(k)} + \alpha \textbf{r}^{(k)} $$
And compute $\alpha$ by minimizing the spectral radius:...
2
votes
1answer
66 views
Right-preconditioning and fixed point linear iterations
Given a linear system $A\textbf{x}=\textbf{b}$, we can express it into the easier-to-solve right-preconditioned form:
$$ AM^{-1}\textbf{y}=\textbf{b}, \quad \textbf{y}= M\textbf{x} $$
On the other ...
9
votes
4answers
494 views
fastest linear system solve for small square matrices (10x10)
I am very interested in optimizing the hell out of linear system solving for small matrices (10x10), sometimes called tiny matrices. Is there a ready solution for this? The matrix can be assumed ...
3
votes
0answers
68 views
Numerical analysis, pivoting and incomplete LU decomposition
When doing LU decomposition, the algorithm will break down if any of the diagonal element $x_{ii}$ is zero. Therefore, we can use pivoting on the matrix such that $x_{ii}$ is no longer zero. That is ...
