2
$\begingroup$

I am interested in solving a large number of small linear systems of equations, $Ax=b$, with $A$ either $2\times2$ or $3\times3$. Assuming none of these systems are actually singular, is there anything to deter me from forming $A^{-1}$ explicitly using Cramer's rule and computing $x=A^{-1}b$? The only things that comes to mind is possible loss of precision in the calculation of $\mathrm{det}A$. Can I expect any improvement over this naive approach by using a canned solver, e.g., LAPACK's dgesv? Namely, would doing so reap numerical benefits that justify the overhead of calling out to LAPACK?

(I am open to methodologies in between these two extremes as well, e.g., hard-coding Gaussian elimination).

$\endgroup$
5
  • 2
    $\begingroup$ Did you check the answer to this question? $\endgroup$
    – nicoguaro
    Oct 16, 2018 at 17:18
  • $\begingroup$ Yes, I did. One answer (GertVdE's) does not actually address numerics. The other (Geoff's) does, but only in a cursory, qualitative way. I am hoping for more detail. $\endgroup$
    – Endulum
    Oct 16, 2018 at 19:16
  • $\begingroup$ mind sharing a bit more detail about your problem? if you're doing e.g. computer vision and use GPUs, here's a thought: group your many matrices into batches and use specialized code to process them all at once. this has been done. a couple of references to get the ball rolling: Lemaitre et al. - batch Cholesky, Dong - batch LU $\endgroup$
    – GoHokies
    Oct 16, 2018 at 19:42
  • $\begingroup$ Intel's MKL (v2018) supports batch operations through vectorized compact routines, see here. another optimization technique (again, MKL-specific) is the use of MKL_DIRECT_CALL's. AFAIK, what this does is to disable error checking. The Intel folks report decent DGEMM speedups over the standard (non-direct?) MKL lib functions. YMMV. $\endgroup$
    – GoHokies
    Oct 16, 2018 at 19:56
  • 1
    $\begingroup$ anyway, to quote from here: "libraries compiled for specific sizes, either via templating or just-in-time (JIT) compilation, provide a large performance boost for small (matrix) sizes; for instance, Intel’s libxsmm3 for extra-small matrix-multiply or batched BLAS for sets of small matrices. Thus, users with such small matrices are encouraged to use special purpose interfaces rather than try to optimize overheads in a general purpose interface." $\endgroup$
    – GoHokies
    Oct 16, 2018 at 20:03

1 Answer 1

2
$\begingroup$

Cramer's rule using a more stable determinant:

https://hal.archives-ouvertes.fr/hal-01500199/document

Implementation of the algorithm is consistently 2.5 times slower than LU factorization of matlab regardless of size of matrix.

$\endgroup$
1
  • 2
    $\begingroup$ To be precise, according to the data in the paper you linked, their Cramer's rule is faster than Matlab LU (whatever it's calling internally) for $N\lessapprox180$. However, I would be much more interested in comparison with an efficient LAPACK exactly for the matrix sizes of interest. $\endgroup$
    – Anton Menshov
    Oct 16, 2018 at 17:53

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.