On Sat, 2007-08-25 at 20:45 +0200, "Sören Gebbert" wrote: > Hi Brad, Dear devs, > i am currently working on a ATLAS - BLAS/LAPACK implementation > within the gpde library. > The idea is to provide full access to ATLAS BLAS level 1-3 > functions and have multi-threaded BLAS 1-3 functions implemented in GRASS > as backup. I have almost finished the BLAS level 1 vector routines and > working on > level 2 the vector - matrix and level 3 matrix - matrix implementation.
It took me a little bit to get back on this. I needed to ensure that ATLAS supports all of the LAPACK functions we need (it only supports a subset). We also missed each other (probably by minutes) several times on IRC. :-) > The grass BLAS level 1-3 implementation is of course not as fast as > the ATLAS library. > ATLAS is cache optimized and uses recursive functions. > Additionally, ATLAS uses pthreads to provide multi threaded level 3 functions. Sounds good to me. I'll have to look into replacing autotools BLAS/LAPACK detection with ATLAS. It shouldn't be terribly difficult unless we want to support all derivatives of CBLAS/CLAPACK. I could use a little direction there as to which way to go. Glynn? Markus? > But for each specialised ATLAS function, a more general grass functions > will exists. Eg: there are level 3 functions for quadratic, symmetric and > hermitian > matrices in BLAS, grass will provide only one matrix-matrix function. > If the user compiles grass without ATLAS support, all the modules > which are using > BLAS functions will still work ... but not that fast. ;) That works for me. This prevents the current caveat of "install BLAS/LAPACK or these modules won't compile". Do all of the lower level functions in lib/gmath, then lib/gpde calls the gmath function, which automatically determins which version to run (actually, it'll be determined at compile time). > I try to keep the interface as simple as i can. The vectors are 1d > float or double arrays > and matrices are 1d float or double arrays with additional row > pointers (using the G_alloc_vector and > G_alloc_matrix functions). I use the same low level implementation > as the ATLAS interface. The user must take care of correct allocation > and row, column counting. > > Additionally i will try to create a wrapper to the LAPACK solver > provided in ATLAS. > IMHO the gpde lu decomposition with row pivoting > and the gpde bandwidth optimized cholesky solver are the backup routines > for most of the LAPACK sover. > > The krylov space solver (cg and bicgstab) will use the ATLAS BLASS 1-2 > functions. > > So the gpde interface is much more low level than the gmath LAPACK/BLAS > wrapper. And i think we can use the gpde routines instead of the native > lapack routines in the gmath warpper. But the matrix functions will work > in row major order. If you have lower level functions, then those should be moved to lib/gmath, IMO. > The blas functions are named like the (IMHO horrible) > netlib-blas naming convention, accept that i have put a N_ before the name. > Eg: > N_cblas_ddot, N_cblas_sdot ... . I'd rather use G_math_ddot (), G_math_sdot ()... I think it'll be simpler to determine the function operation at compile time rather than run time. -- 73, de Brad KB8UYR/6 <rez touchofmadness com> _______________________________________________ grass-dev mailing list [email protected] http://grass.itc.it/mailman/listinfo/grass-dev
