Fast linear system solver for D? [closed]

2019-06-18 22:58发布

问题:

Where can I get a fast linear system solver written in D? It should be able to take a square matrix A and a vector b and solve the equation Ax = b for b and, ideally, also perform explicit inversion on A. I have one I wrote myself, but it's pretty slow, probably because it's completely cache-naive. However, for my use case, I need something with the following absolute, non-negotiable requirements, i.e. if it doesn't meet these, then I don't otherwise care how good it otherwise is:

  1. Must be licensed public domain, Boost license, or some similar permissive license. Ideally it should not require attribution in binaries (i.e. not BSD), though this point is somewhat negotiable.

  2. Must be written in pure D or easily translatable to pure D. Inscrutable Fortran code (i.e. LAPACK) is not a good answer no matter how fast it is.

  3. Must be optimized for large (i.e. n > 1000) systems. I don't want something that's designed for game programmers to solve 4x4 matrices really, really fast.

  4. Must not be inextricably linked to a huge library of stuff I don't need.

Edit: The reason for these seemingly insane requirements is that I need this code for a permissively licensed open source library that I don't want to have any third-party dependencies.

回答1:

If you don't like Fortran code, one reasonably fast C++ dense matrix library with modest multi-core support, well-written code and a good user-interface is Eigen. It should be straightforward to translate its code to D (or to take some algorithms from it).

And now my "think about your requirements": there is a reason why "everyone" (Mathematica, Matlab, Maple, SciPy, GSL, R, ...) uses ATLAS / LAPACK, UMFPACK, PARDISO, CHOLMOD etc. It is hard work to write fast, multi-threaded, memory-efficient, portable and numerically stable matrix solvers (trust me, I have tried). A lot of this hard work has gone into ATLAS and the rest.

So my approach would be to write bindings for the relevant library depending on your matrix type, and link from D against the C interfaces. Maybe the bindings in multiarray are enough (I haven't tried). Otherwise, I'd suggest looking at another C++ library, namely uBlas and the respective bindings for ideas.