caml-list - the Caml user's mailing list
 help / color / mirror / Atom feed
* [Caml-list] Matrix libraries
@ 2003-12-13 13:51 romildo
  2003-12-13 16:36 ` Markus Mottl
  0 siblings, 1 reply; 8+ messages in thread
From: romildo @ 2003-12-13 13:51 UTC (permalink / raw)
  To: caml-list

Hello.

I am looking for a good library for numerical
matrixes manipulation in OCaml. Can anybody help me?

Regards,

Romildo

-------------------
To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr
Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/
Beginner's list: http://groups.yahoo.com/group/ocaml_beginners


^ permalink raw reply	[flat|nested] 8+ messages in thread

* Re: [Caml-list] Matrix libraries
  2003-12-13 13:51 [Caml-list] Matrix libraries romildo
@ 2003-12-13 16:36 ` Markus Mottl
  2003-12-13 21:38   ` Oleg Trott
  0 siblings, 1 reply; 8+ messages in thread
From: Markus Mottl @ 2003-12-13 16:36 UTC (permalink / raw)
  To: romildo; +Cc: OCaml

On Sat, 13 Dec 2003, romildo@uber.com.br wrote:
> I am looking for a good library for numerical
> matrixes manipulation in OCaml. Can anybody help me?

Although it has already been on the web for a few weeks, I hadn't actually
announced it yet, waiting for comments of some co-developers:

Lacaml is available in a new major version. This library interfaces
the Fortran libraries BLAS and LAPACK for heavy-weight linear algebra
(i.e. matrix computations). New features include support for complex
transformations (complex numbers), and a convenient way of accessing
submatrices using labels. As usual, you can choose either single or
double precision computations. The computations can run in parallel
on SMP-machines.

You can download the library here:

  http://www.oefai.at/~markus/home/ocaml_sources.html#LACAML

I'd be grateful for feedback!

Best regards,
Markus Mottl

-- 
Markus Mottl          http://www.oefai.at/~markus          markus@oefai.at

-------------------
To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr
Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/
Beginner's list: http://groups.yahoo.com/group/ocaml_beginners


^ permalink raw reply	[flat|nested] 8+ messages in thread

* Re: [Caml-list] Matrix libraries
  2003-12-13 16:36 ` Markus Mottl
@ 2003-12-13 21:38   ` Oleg Trott
  2003-12-14 15:01     ` Markus Mottl
                       ` (2 more replies)
  0 siblings, 3 replies; 8+ messages in thread
From: Oleg Trott @ 2003-12-13 21:38 UTC (permalink / raw)
  To: Markus Mottl, romildo; +Cc: OCaml

On Saturday 13 December 2003 11:36 am, Markus Mottl wrote:
> On Sat, 13 Dec 2003, romildo@uber.com.br wrote:
> > I am looking for a good library for numerical
> > matrixes manipulation in OCaml. Can anybody help me?
>
> Although it has already been on the web for a few weeks, I hadn't actually
> announced it yet, waiting for comments of some co-developers:
>
> Lacaml is available in a new major version. This library interfaces
> the Fortran libraries BLAS and LAPACK for heavy-weight linear algebra
> (i.e. matrix computations). New features include support for complex
> transformations (complex numbers), and a convenient way of accessing
> submatrices using labels. As usual, you can choose either single or
> double precision computations. The computations can run in parallel
> on SMP-machines.
>
> You can download the library here:
>
>   http://www.oefai.at/~markus/home/ocaml_sources.html#LACAML
>
> I'd be grateful for feedback!

All right! I was actually thinking of giving it lately. First of all, I think 
the library is very interesting. Here are some of my assorted thoughts 
regarding improving it.

1.) I think a matrix library has to have >= 2 levels of user-frienliness. 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 
low-level 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)

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)

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 low-level 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 auto-generate 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)

3.) Regarding "WORK" arguments. Why not have a shared workspace:

(* module Workspace, untested *)

(* private *)
let workspace = ref (Array1.create float64 100)

(* public *)
let size () = Array1.dim !workspace

let resize n = workspace := Array1.create float64 n

let delete () = resize 0

let get_float64 n = if size () >= n then !workspace 
                                                else (resize n; !workspace)

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 inter-operate with 
threading or SMP though.

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.

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. 

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 :-)

-- 
Oleg Trott <oleg_trott@columbia.edu>

-------------------
To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr
Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/
Beginner's list: http://groups.yahoo.com/group/ocaml_beginners


^ permalink raw reply	[flat|nested] 8+ messages in thread

* Re: [Caml-list] Matrix libraries
  2003-12-13 21:38   ` Oleg Trott
@ 2003-12-14 15:01     ` Markus Mottl
  2003-12-14 23:10       ` Oleg Trott
  2003-12-15 17:05     ` [Caml-list] CamlFloat Shivkumar Chandrasekaran
  2003-12-15 17:14     ` [Caml-list] Matrix libraries Shivkumar Chandrasekaran
  2 siblings, 1 reply; 8+ messages in thread
From: Markus Mottl @ 2003-12-14 15:01 UTC (permalink / raw)
  To: Oleg Trott; +Cc: romildo, OCaml

On Sat, 13 Dec 2003, Oleg Trott wrote:
> 1.) I think a matrix library has to have >= 2 levels of user-frienliness. 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 
> low-level 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 non-optional 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 low-level interface to
BLAS/LAPACK. Others are working on higher-level 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 low-level stuff using LACAML and leave higher-level 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 low-level 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 auto-generate 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 C-code for the interface is very
little work. It's the OCaml-code 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 inter-operate with 
> threading or SMP though.

Yes, that's exactly the problem: what about SMP-machines and threading
then? It's really not this difficult to wrap all LACAML-functions 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 Bigarray-module and the .{}-notation. CamlFloat, the higher-level
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 Bigarray-module.

This is a bit inconvenient, I agree, but feel free to complain to INRIA
to introduce some more direct conversions in the Bigarray-module :-)

> 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 caml-list-request@inria.fr Archives: http://caml.inria.fr
Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/
Beginner's list: http://groups.yahoo.com/group/ocaml_beginners


^ permalink raw reply	[flat|nested] 8+ messages in thread

* Re: [Caml-list] Matrix libraries
  2003-12-14 15:01     ` Markus Mottl
@ 2003-12-14 23:10       ` Oleg Trott
  2003-12-15 10:01         ` Markus Mottl
  0 siblings, 1 reply; 8+ messages in thread
From: Oleg Trott @ 2003-12-14 23:10 UTC (permalink / raw)
  To: Markus Mottl; +Cc: OCaml

On Sunday 14 December 2003 10:01 am, Markus Mottl wrote:
> In any case, LACAML is supposed to stay a low-level interface to
> BLAS/LAPACK. 

The $1e6 quesion: do you want the library to be safe from user abuse, i.e.
no function input should result in corrupted memory ?

> Others are working on higher-level ones, e.g.:
>
>   http://www.math.ucsb.edu/~lyons/camlFloat/index.html

Thanks for the link! CamlFloat isn't even in Google yet. It must be new. Are 
there others?

> I'd also like to see something like that, but I am skeptical that this
> is realistically possible. Writing the C-code for the interface is very
> little work. It's the OCaml-code around it that costs time to write,
> especially error handling.

Some of this error-handling like checking that input matrices/vectors have 
compatible sizes seems tedious (and error-prone), and I think it can be 
auto-generated from parsing *.f files (including comments) But, yes, maybe 
it's too hard and not worth the effort.

> > 3.) Regarding "WORK" arguments. Why not have a shared workspace:
>
> Yes, that's exactly the problem: what about SMP-machines and threading
> then?

Each thread oviously needs its own workspace. I think this can be done using

(int * vec ref) list ref            (* int  = id (self ()) *)

association list (or hash table).  Now, "get", "resize", etc. could check if 
the thread has its workspace and return it, allocating if necessary.

I think there is a problem with this approach though: each thread's workspace
needs to be removed once the thread terminates. OCaml has at_exit but no 
at_thread_exit that I can find (Maybe it can be defined using 
Sys.set_signal's ?)

> 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 Bigarray-module and the .{}-notation. 

I think the inconvenience is minimal:

a.{...}  vs a.mat_data.{...}

(and it's just typing)

But it saves you from the very error-prone and boring task of having to 
remember which variables designate which dimensions ("Is it m x n or n x k, 
did I transpose that?"), etc.

OTOH submatrices/slices probably aren't the most frequently used features,
so I haven't made up my mind as to which is better.

> > 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. 

Just a runtime error if LD does not equal the number of rows would be good.
(If your matrix type includes info about both the number of rows and LD)

-- 
Oleg Trott <oleg_trott@columbia.edu>

-------------------
To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr
Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/
Beginner's list: http://groups.yahoo.com/group/ocaml_beginners


^ permalink raw reply	[flat|nested] 8+ messages in thread

* Re: [Caml-list] Matrix libraries
  2003-12-14 23:10       ` Oleg Trott
@ 2003-12-15 10:01         ` Markus Mottl
  0 siblings, 0 replies; 8+ messages in thread
From: Markus Mottl @ 2003-12-15 10:01 UTC (permalink / raw)
  To: Oleg Trott; +Cc: OCaml

On Sun, 14 Dec 2003, Oleg Trott wrote:
> On Sunday 14 December 2003 10:01 am, Markus Mottl wrote:
> > In any case, LACAML is supposed to stay a low-level interface to
> > BLAS/LAPACK. 
> 
> The $1e6 quesion: do you want the library to be safe from user abuse, i.e.
> no function input should result in corrupted memory ?

Yes, absolutely! That's why I think that checking all parameters before
passing them to LAPACK is a necessity even if it is a bit tedious
to implement considering the huge number of arguments that some
LAPACK-functions take.

> Thanks for the link! CamlFloat isn't even in Google yet. It must be new. Are 
> there others?

Yes, CamlFloat is fairly new. I don't know whether there are
others. Interfacing BLAS/LAPACK is not for the fainthearted: it has tons
of functions each with lots of parameters of various types so it's a
bit boring to do and very easy to make mistakes. Most people probably
don't have the patience for this.

CamlFloat definitely seems to aim at professional use so I expect that
it may become the library of choice for linear algebra, my library acting
as convenient low-level interface.

> Some of this error-handling like checking that input matrices/vectors have 
> compatible sizes seems tedious (and error-prone), and I think it can be 
> auto-generated from parsing *.f files (including comments) But, yes, maybe 
> it's too hard and not worth the effort.

By relying on as much abstraction as possible, the effort of parameter
checking can be reduced to a reasonable level. Using existing functions as
guideline, it most often shouldn't take longer than an hour to implement
and test new LAPACK-functions. Due to the macro system, this automatically
covers alls four kinds of functions (S/D/C/Z).

> Each thread oviously needs its own workspace. I think this can be done using
> 
> (int * vec ref) list ref            (* int  = id (self ()) *)
> 
> association list (or hash table).  Now, "get", "resize", etc. could check if 
> the thread has its workspace and return it, allocating if necessary.
> 
> I think there is a problem with this approach though: each thread's workspace
> needs to be removed once the thread terminates. OCaml has at_exit but no 
> at_thread_exit that I can find (Maybe it can be defined using 
> Sys.set_signal's ?)

Well, there is always the tension between convenience and simplicity. I
think that features like this should be added in a separate layer.

> > 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 Bigarray-module and the .{}-notation. 
> 
> I think the inconvenience is minimal:
> 
> a.{...}  vs a.mat_data.{...}
> 
> (and it's just typing)

And what about other matrix types like band, tridiagonal, etc.? You'd need
a discriminated union for this, which requires pattern matching.  Again,
I think this should be done on a higher level - e.g. as in CamlFloat.

> But it saves you from the very error-prone and boring task of having to 
> remember which variables designate which dimensions ("Is it m x n or n x k, 
> did I transpose that?"), etc.
> 
> OTOH submatrices/slices probably aren't the most frequently used features,
> so I haven't made up my mind as to which is better.

Exactly. I think that libraries should attempt to optimize for the likely
cases. And a low-level library, which is mainly used by other libraries
again, should IMHO not go overboard with convenience for end users.

Regards,
Markus

-- 
Markus Mottl          http://www.oefai.at/~markus          markus@oefai.at

-------------------
To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr
Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/
Beginner's list: http://groups.yahoo.com/group/ocaml_beginners


^ permalink raw reply	[flat|nested] 8+ messages in thread

* [Caml-list] CamlFloat
  2003-12-13 21:38   ` Oleg Trott
  2003-12-14 15:01     ` Markus Mottl
@ 2003-12-15 17:05     ` Shivkumar Chandrasekaran
  2003-12-15 17:14     ` [Caml-list] Matrix libraries Shivkumar Chandrasekaran
  2 siblings, 0 replies; 8+ messages in thread
From: Shivkumar Chandrasekaran @ 2003-12-15 17:05 UTC (permalink / raw)
  To: caml-list; +Cc: Bill Lyons

I was hoping to hold off on announcing this until the new year..... So 
please take this as just a "pre-announcement".

We plan to release CamlFloat, a "matlab-like" (hopefully better) 
interface to Lapack+Blas "soon". Unfortunately, it is *not* (yet) built 
on top of Lacaml.  It comes with a fully documented interface, a 
tutorial, and a block-sparse matrix module. The code has been 
extensively used in our own research.

For people who cannot wait, they can find a preliminary release here:

http://www.math.ucsb.edu/~lyons/camlFloat/

CleanFloat, a similar package for Clean (http://cs.kun.nl/~clean) will 
also be made available.

--shiv--


On Dec 13, 2003, at 1:38 PM, Oleg Trott wrote:

> On Saturday 13 December 2003 11:36 am, Markus Mottl wrote:
>> On Sat, 13 Dec 2003, romildo@uber.com.br wrote:
>>> I am looking for a good library for numerical
>>> matrixes manipulation in OCaml. Can anybody help me?
>>
>> Although it has already been on the web for a few weeks, I hadn't 
>> actually
>> announced it yet, waiting for comments of some co-developers:
>>
>> Lacaml is available in a new major version. This library interfaces
>> the Fortran libraries BLAS and LAPACK for heavy-weight linear algebra
>> (i.e. matrix computations). New features include support for complex
>> transformations (complex numbers), and a convenient way of accessing
>> submatrices using labels. As usual, you can choose either single or
>> double precision computations. The computations can run in parallel
>> on SMP-machines.
>>
>> You can download the library here:
>>
>>   http://www.oefai.at/~markus/home/ocaml_sources.html#LACAML
>>
>> I'd be grateful for feedback!
>
> All right! I was actually thinking of giving it lately. First of all, 
> I think
> the library is very interesting. Here are some of my assorted thoughts
> regarding improving it.
>
> 1.) I think a matrix library has to have >= 2 levels of 
> user-frienliness. 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
> low-level 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)
>
> 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)
>
> 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 low-level 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 auto-generate 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)
>
> 3.) Regarding "WORK" arguments. Why not have a shared workspace:
>
> (* module Workspace, untested *)
>
> (* private *)
> let workspace = ref (Array1.create float64 100)
>
> (* public *)
> let size () = Array1.dim !workspace
>
> let resize n = workspace := Array1.create float64 n
>
> let delete () = resize 0
>
> let get_float64 n = if size () >= n then !workspace
>                                                 else (resize n; 
> !workspace)
>
> 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 inter-operate with
> threading or SMP though.
>
> 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.
>
> 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.
>
> 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 :-)
>
> -- 
> Oleg Trott <oleg_trott@columbia.edu>
>
> -------------------
> To unsubscribe, mail caml-list-request@inria.fr Archives: 
> http://caml.inria.fr
> Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: 
> http://caml.inria.fr/FAQ/
> Beginner's list: http://groups.yahoo.com/group/ocaml_beginners

-------------------
To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr
Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/
Beginner's list: http://groups.yahoo.com/group/ocaml_beginners


^ permalink raw reply	[flat|nested] 8+ messages in thread

* Re: [Caml-list] Matrix libraries
  2003-12-13 21:38   ` Oleg Trott
  2003-12-14 15:01     ` Markus Mottl
  2003-12-15 17:05     ` [Caml-list] CamlFloat Shivkumar Chandrasekaran
@ 2003-12-15 17:14     ` Shivkumar Chandrasekaran
  2 siblings, 0 replies; 8+ messages in thread
From: Shivkumar Chandrasekaran @ 2003-12-15 17:14 UTC (permalink / raw)
  To: caml-list


On Dec 13, 2003, at 1:38 PM, Oleg Trott wrote:

>
> 3.) Regarding "WORK" arguments. Why not have a shared workspace:
>
> (* module Workspace, untested *)

I was thinking of using weak arrays for this, so that a one-time use of 
a large workspace does not consume memory for too long.

>
> (* private *)
> let workspace = ref (Array1.create float64 100)
>
> (* public *)
> let size () = Array1.dim !workspace
>
> let resize n = workspace := Array1.create float64 n
>
> let delete () = resize 0
>
> let get_float64 n = if size () >= n then !workspace
>                                                 else (resize n; 
> !workspace)
>
> 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 inter-operate with
> threading or SMP though.

Since the "naive" way of doing multiple OCaml-threads does not imply 
that multiple processors will be used (see old posts on this subject), 
I was hoping to write yet another interface to Lapack+Blas where each 
fortran call would run in its own C thread. That way the user can 
precisely control the threading. Anybody has any advice on this matter?

Apparently ATLAS BLAS already can exploit SMP, but it does not seem to 
matter for small matrix sizes, which is my case. I just happen to have 
lots of them!


>
> 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.

See my previous "pre-announcement" for CamlFloat. It provides such a 
mechanism, and tracks all the nasty details of unit diagonal 
lower-triangular matrices, etc., etc., for you. That is, all the usual 
things that make most non-trivial Lapack programs buggy.

>
> 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.
>
> 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 :-)
>
> -- 
> Oleg Trott <oleg_trott@columbia.edu>
>
> -------------------
> To unsubscribe, mail caml-list-request@inria.fr Archives: 
> http://caml.inria.fr
> Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: 
> http://caml.inria.fr/FAQ/
> Beginner's list: http://groups.yahoo.com/group/ocaml_beginners

-------------------
To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr
Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/
Beginner's list: http://groups.yahoo.com/group/ocaml_beginners


^ permalink raw reply	[flat|nested] 8+ messages in thread

end of thread, other threads:[~2003-12-15 21:15 UTC | newest]

Thread overview: 8+ messages (download: mbox.gz / follow: Atom feed)
-- links below jump to the message on this page --
2003-12-13 13:51 [Caml-list] Matrix libraries romildo
2003-12-13 16:36 ` Markus Mottl
2003-12-13 21:38   ` Oleg Trott
2003-12-14 15:01     ` Markus Mottl
2003-12-14 23:10       ` Oleg Trott
2003-12-15 10:01         ` Markus Mottl
2003-12-15 17:05     ` [Caml-list] CamlFloat Shivkumar Chandrasekaran
2003-12-15 17:14     ` [Caml-list] Matrix libraries Shivkumar Chandrasekaran

This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for NNTP newsgroup(s).