iraikov / chicken-lapack Goto Github PK

As per issue #3 test cases need to be added to ensure that the LAPACK functions used by ATLAS can be called correctly, and produce the expected result. Moreover, these test cases will serve as example documentation for users in the future to see how each of the functions can be used in practice.

Ideally I would like to use the test and maybe even the test-generative eggs to create these tests, as single static asserts are not always the most intuitive way to run tests. Moreover, the test runners in the test/test-generative eggs provide much more informative output.

Care should also be taken to ensure complex-number functions (prefixed with z or c) also work, as they can be subtly different than real-number functions (prefixed with s or d) as evidenced by issues #3 and #2.

This is a minor grievance, but I think it would make things more consistent if the lapack-wrap and lapack-wrapx macros were renamed to atlas-wrap and atlas-wrapx respectively.

I don't usually like submitting issues like this, but I think it would be valuable to distinguish between the macros that exist within atlas-lapack.scm and lapack-extras.scm, as the two macros are somewhat dissimilar at present (see https://github.com/ThatGeoGuy/chicken-lapack/tree/lapack-extras). Overall I think this will prevent anyone who might look at the code from confusing one as merely a copy of the other (though in many ways, they are very similar).

Since the lapack-extras module is not yet merged into master, I think this can at present wait until that time, but I thought it valuable to bring up beforehand.

This is somewhat of a feature request to complete some of the functionality made available via LAPACK. Currently, the egg only supports optimized LAPACK routines made available via ATLAS; however, most users will compile ATLAS with the complete LAPACK package using the following configure flag:

--with-netlib-lapack-tarfile=/path/to/netlib/lapack-x.x.x.tgz

where the xs in the file name correspond to the version of lapack compiled with ATLAS. Moreover, most Linux distributions (or package managers in general) likewise package the full lapack functionality with ATLAS. To see if you have all of the functions, look for dgeev_ or another driver routine that is in LAPACK but not provided via ATLAS in the symbols table:

$ nm /usr/lib/liblapack.a | grep dgeev
dgeev.o:
0000000000000000 T dgeev_
dgeevx.o:
0000000000000000 T dgeevx_

This means that these functions are available, but not fully provided via the ATLAS interface. You can compile a regular C program to take advantage of these functions as follows:

// Example C program
#include <stdio.h>
#include <stdlib.h>
#include <clapack.h>

int main() {
    dgeev_( /* args here */ );
    return 0;
}

Which is then compiled with:

$ gcc example.c -o example -llapack -lf77blas -lgfortran -lcblas -latlas -lm

So theoretically, this means that by perhaps adding libf77blas and libgfortran to the compiler flags in the .setup file, as well as making some changes to expose the interface to dgeev_ via CHICKEN scheme's bind egg, it should be plausible to expose additional LAPACK routines from this egg. The reason in particular I chose dgeev_ is because I have a vested interest in Eigenvalues / Eigenvectors, and I think beyond the general matrix solvers and the symmetric-positive-(semi)-definite solvers in LAPACK / ATLAS the Eigen solution functions are the most useful and versatile.

I don't yet know what changes to make to expose dgeev_ to scheme bindings, but I believe it will takes some changes to existing code, because currently the code works by prepending clapack_ to any of the called bindings (e.g. if you call dgesv in scheme it will expand and call to clapack_dgesv in ATLAS). This ought to be a good start to figuring out the problem in any case, any suggestions or advice are appreciated, as I am willing to help implement this feature.

As mentioned at #2 (comment), the atlas-lapack egg has a rather critical bug that should prevent any functions that operate over complex or double complex matrices from working.

The basic issue is in the lapack-wrapx macro, which defines the length operation over complex / double complex vectors as (lambda (v) (fx/ 2 (f32vector-length v))). In actuality, this makes no sense since we would be dividing the length of our vector v from the fixnum 2. Since this is fixnum division, this will almost always result in 2. The proposed fix is to change it such that:

(lambda (v) (fx/ (f32vector-length v) 2))

However, I think the fix should also include a relevant test-case. I think in general the library is lacking test-cases, and could use some basic identities to prove that the functions are being called appropriately. They would also serve as examples for people who may wish to use the egg in the future. I'll be blunt, I could really use the examples myself, even though I consider myself fairly familiar with the library as it stands.

I'm willing to make the fix and submit a pull request. Likewise, I will also try to construct some automated tests. Perhaps I should open a separate issue for adding tests?