I have a sparse $n\times n$ matrix $A$ with a pretty interesting structure. It has a block structure with a symmetric structure but asymmetric blocks. Expressed mathematically the block $A_{jk} = A_{kj}$ but $A_{jk} \neq A_{jk}^T$, so my matrix is not actually symmetric. Are there any algorithms that can take advantage of this kind of system?
$\begingroup$
$\endgroup$
1
-
$\begingroup$ There would a possibility of an advantage if various block submatrices commute. $\endgroup$– hardmath ♦Sep 11, 2022 at 20:54
Add a comment
|
2 Answers
$\begingroup$
$\endgroup$
The PARDISO solver from oneMKL supports structurally symmetric matrices. Quoting from the Intel oneMKL library documentation:
Structurally Symmetric Matrices
The solver first computes a symmetric fill-in reducing permutation P followed by the parallel numerical factorization of $PAP^T = QLU^T$. The solver uses partial pivoting in the supernodes and an approximation of $X$ is found by forward and backward substitution and optional iterative refinement.
The library also supports block-sparse storage formats.
I'd also have a look at the PANUA - PARDISO, which has improved threading and GPU capabilities. A license can be obtained free of charge for academic users.
$\begingroup$
$\endgroup$
4
First of all - your mathematical description seems wrong: $A_{jk} = A_{kj}, \forall\ k,j \in \{1,\dots,n\}$ implies that the matrix is symmetric. So it is not entirely clear what you mean. Also, $A^T_{jk}$ is exactly defined via $A^T_{jk} = A_{kj}$ (see e.g. Wikipedia). This would somehow be in conflict to your equations.
As to the actual question: While I am not sure what your structure actually is - I am also not aware of any method that might exploit a structure that resembles yours. For a nice overview on iterative methods I'd recommend the templates book. But the most "structure exploitation" there is with regard to symmetry. I hope this helps some..
-
2$\begingroup$ Is it possible that Flusslauf means that the locations of the non-zeros is symmetric but their values are different? There was some work done in the late '90s/early 2000s to show that ILU preconditioners which allowed fill in the locations that were structurally non-zero did better than the ILU preconditioners with levels of fill or thresholds, but I don't have that reference handy. These are very challenging to parallelize for distributed memory machines with MPI. $\endgroup$ Aug 17, 2021 at 15:08
-
3$\begingroup$ Ajk refers the the blocks not the individual entries which is why I talk about the transpose of the block Ajk. The transpose of an individual element would not make sense. $\endgroup$– EMPAug 17, 2021 at 16:18
-
$\begingroup$ @EMP, usually $A_{jk}$ means the $j,k^{\rm th}$ element of matrix $A$ not the block, so some people are confused by your question. There's really no standard notation that I'm aware of to denote a block vs. an element. $\endgroup$ Aug 17, 2021 at 16:32
-
1$\begingroup$ Agreed that there's not typically a good way to delineate block versus element, I have seen the use of indices to refer to blocks rather than elements in some papers. But the lack of clarity is why I had specified the problem beforehand discussing the block structure and also I think the use of the transpose notation makes it pretty clear I'm discussing the blocks. $\endgroup$– EMPAug 17, 2021 at 16:41