[Camllist] Matrix libraries
[
Home
]
[ Index:
by date

by threads
]
[ Message by date: previous  next ] [ Message in thread: previous  next ] [ Thread: previous  next ]
[ Message by date: previous  next ] [ Message in thread: previous  next ] [ Thread: previous  next ]
Date:   (:) 
From:  Markus Mottl <markus@o...> 
Subject:  Re: [Camllist] Matrix libraries 
On Sat, 13 Dec 2003, Oleg Trott wrote: > 1.) I think a matrix library has to have >= 2 levels of userfrienliness. One > level should be easy to use and have a MATLAB feel to it, i.e. if someone > wants the eigenvalues of a real matrix, he should just be able to write > e_vals = eigenvalues(A) > or > (e_vals, e_vecs) = eigensystem(A) > without having to worry about details. This should be built on top of a > lowlevel interface similar to LAPACK functions themselves. E.g. DGEEV does > not even return the eigenvectors corresponding to complex eigenvalues, but > rather something from which they can be reconstructed (which is what > one might actually want sometimes) This is true for many functions in LAPACK. E.g. there is none for computing the inverse of a matrix directly, because this is inefficient and not very stable if you just want to solve linear equations. You'd have to take the detour over e.g. xgetrf+xgetri. > Right now, it appears to me that LaCaml is close to the low level, but not > quite there (because e.g. some LAPACK interfaces _return_ a matrix rather > then letting the user modify some piece of memory in place like in FORTRAN) Maybe this is just a misunderstanding: some functions do return matrices (and vectors), but it's always the argument that is overwritten anyway. There is probably no purpose to having these functions return anything if a nonoptional argument containing the result is overwritten. This will most likely just confuse people. I have changed this now. Thanks for pointing this out! In any case, LACAML is supposed to stay a lowlevel interface to BLAS/LAPACK. Others are working on higherlevel ones, e.g.: http://www.math.ucsb.edu/~lyons/camlFloat/index.html I am probably going to meet up with one of the authors, Bill Lyons, soon, who happens to come to a nearby conference, and we are intending to talk about coordinating our work. It might be interesting to factor out the lowlevel stuff using LACAML and leave higherlevel functions and operators to CamlFloat. > 2.) Currently, LaCaml interfaces only a fraction of LAPACK. LAPACK is big, and > interfacing all of it by hand is a huge job. On the other hand, MOST of the > information needed for the lowlevel interface can be parsed from function > declarations in *.f files (I kept thinking about this while using LAPACK in > my C++ programs). Maybe the way to interface all of LAPACK is > to parse those *.f files into some intermediate form, let humans edit it (e.g. > add minimum workspace size expressions, etc.  obvious from docs, but hard to > parse/NLP), and then autogenerate the whole interface from the result. This > will have the addtional benefit of being consistent (If you just let a bunch > of developers contribute interface code directly, they will bring their own > "style" to their parts of the interface) I'd also like to see something like that, but I am skeptical that this is realistically possible. Writing the Ccode for the interface is very little work. It's the OCamlcode around it that costs time to write, especially error handling. > 3.) Regarding "WORK" arguments. Why not have a shared workspace: [snip] > This way, one does not need to "hoist" workspace allocation out of loops > manually. OTOH performance does not suffer from allocating workspace > for each function call. I don't know how this would interoperate with > threading or SMP though. Yes, that's exactly the problem: what about SMPmachines and threading then? It's really not this difficult to wrap all LACAMLfunctions that use workspaces into ones that share a pool of them if threading is not an issue. > 4) Submatrices. For a function f that takes a matrix argument "a" of size m by > n, you normally use something like > > f ar ac a m n > > where ar and ac refer to the "beginning" of the submatrix. Even though they > can default to 1, this is kind of inconvenient compared to > > f a > > I would propose using the latter everywhere, and letting the "mat" type > include ar, ac, m and n. > > submatrix : int > int > int > int > mat > mat > (* submatrix x1 y1 x2 y2 a *) > > would provide "views" into matrix a. I had indeed thought about this, but that would have made it more inconvenient to people who want to keep accessing "a" directly using the Bigarraymodule and the .{}notation. CamlFloat, the higherlevel library, has such a matrix type (introduced by the constructor "Dense"). I suppose the developers are intending to add further ones for triangular, band, symmetric matrices, etc. > 5) I think a function that lets one view matrices as vectors vec_of_mat is > needed. They are all just ordered sets of numbers after all. This would collide with 4), wouldn't it? Vector arguments in LAPACK don't have anything like a "leading dimension" so this wouldn't help you unless you allow copying data. If, however, we keep using bigarrays as concrete representation for matrices, then you could as well just call "genarray_of_array2" followed by "reshape_1" in the Bigarraymodule. This is a bit inconvenient, I agree, but feel free to complain to INRIA to introduce some more direct conversions in the Bigarraymodule :) > 6) License. Have you thought of adding a static linking exception to > your LGPL license? Even INRIA did with Caml standard library. I might > even consider contributing then :) Ok, I have relaxed the conditions of the license. Please contribute now :) Regards, Markus  Markus Mottl http://www.oefai.at/~markus markus@oefai.at  To unsubscribe, mail camllistrequest@inria.fr Archives: http://caml.inria.fr Bug reports: http://caml.inria.fr/bin/camlbugs FAQ: http://caml.inria.fr/FAQ/ Beginner's list: http://groups.yahoo.com/group/ocaml_beginners